home *** CD-ROM | disk | FTP | other *** search
/ Resource Library: Graphics / graphics-16000.iso / general / raytrace / renaisnc / delaunay.lha / Delaunay / k-Dtree.c < prev    next >
C/C++ Source or Header  |  1992-01-04  |  9KB  |  354 lines
  1. /**
  2.  **  Maintain k-D tree
  3.  **/
  4.  
  5. #include <stdio.h>
  6. #include <math.h>
  7. #include <assert.h>
  8. #include "k-Dtree.h"
  9.  
  10. #ifndef TRUE
  11. #  define TRUE 1
  12. #  define FALSE 0
  13. #endif
  14.  
  15. /*
  16. ** Test to see if two rectangles overlap.  Return TRUE if the do, FALSE
  17. ** otherwise.
  18. */
  19. #define RectOverlap(Rect1,Rect2)  \
  20.  ( ! ( (Rect1)->xH < (Rect2)->xL || (Rect2)->xH < (Rect1)->xL ||     \
  21.        (Rect1)->yH < (Rect2)->yL || (Rect2)->yH < (Rect1)->yL    ))  
  22.  
  23.  
  24. /*------------------------------------
  25.   Check to see if the point X,Y is in 
  26.   the rectangle Rect.
  27.  ------------------------------------*/
  28. #define PointInRectangle(Xs,Ys,Rect)                     \
  29.     ( ((Xs) >= (Rect)->xL)  &&  ((Xs) < (Rect)->xH)  &&  \
  30.       ((Ys) >= (Rect)->yL)  &&  ((Ys) < (Rect)->yH) )
  31.  
  32.  
  33.  
  34. /*  Number of nodes/samples to allocate at a time */
  35. #define GLOB_SIZE 1024
  36.  
  37.  
  38.  
  39. /*----------------------------------------
  40.   Subdivision level
  41.   --------------------------------------*/
  42. static
  43. int level;
  44.  
  45. int AdaptiveSplits = FALSE;
  46. int Animating = FALSE;
  47.  
  48. HierarchicalRegion Root;
  49. RectType RootBoundary;
  50.  
  51. extern double atof();
  52. extern  HierarchicalRegion *PropagateSample();
  53.  
  54.  
  55. /* --------------------------------------------------------------------------
  56.    Split a node boundary rectangle into the rectangle for the Left 
  57.    (if LeftOrRight == TRUE) or Right (if LeftOrRight == FALSE) subnode.
  58.    ---------------------------------------------------------------------------*/
  59. SplitBoundary(LeftOrRight, splitDirection, splitValue,
  60.           OldBoundary, NewBoundary)
  61.      int LeftOrRight;
  62.      double splitValue;
  63.      RectType *OldBoundary, *NewBoundary;
  64. {
  65.   double Xs, Ys;
  66.   
  67.   *NewBoundary = *OldBoundary;
  68.   if ( splitDirection )  /* Odd levels split in X */
  69.     {  /* Split in X */
  70.       Xs = splitValue;
  71.       if (LeftOrRight) {    /* Left Child */
  72.     NewBoundary->xH = Xs;
  73.       } else {            /* Right child */
  74.     NewBoundary->xL = Xs; 
  75.       }
  76.     }
  77.   else
  78.     {  /* Split in Y */
  79.       Ys = splitValue;
  80.       if (LeftOrRight) {    /* Left Child  */
  81.     NewBoundary->yH = Ys;
  82.       } else {            /* Right Child */
  83.     NewBoundary->yL = Ys; 
  84.       }
  85.     }
  86. }
  87.  
  88.  
  89.  
  90. /*-----------------------
  91.   allocate a node
  92.   ----------------------*/
  93.  
  94. static
  95.   HierarchicalRegion *
  96.   CreateNode()
  97. {
  98.   char *calloc();
  99.   static HierarchicalRegion *FreeNodes;
  100.   static nFreeNodes = 0;
  101.  
  102.  
  103.   if (nFreeNodes <= 0) {
  104.     FreeNodes = (HierarchicalRegion *) 
  105.                   calloc( GLOB_SIZE, sizeof(HierarchicalRegion));
  106.     if(FreeNodes == NULL)
  107.       {
  108.     fprintf(stderr,"Fatal error - Can't allocate nodes\n");
  109.     abort();
  110.       }
  111.     nFreeNodes = GLOB_SIZE;
  112.   }
  113.   return &(FreeNodes[--nFreeNodes]);
  114. }
  115.  
  116.  
  117. static SampleData *SampleFreeList = NULL;
  118.  
  119. FreeSample( Sample )
  120. SampleData *Sample;
  121. {
  122.   /* Link the sample onto the free list */
  123.   *( (SampleData **) Sample) = SampleFreeList;
  124.   SampleFreeList = Sample;
  125. }
  126.  
  127. SampleData *CreateSample()
  128. {
  129.   char *calloc();
  130.   static SampleData *FreeSamples, *Sample;
  131.   static int nFreeSamples = 0;
  132.  
  133.   /* Get one off the free list if possible */
  134.   if (SampleFreeList != NULL) {
  135.     Sample = SampleFreeList;
  136.     SampleFreeList = *( (SampleData **) SampleFreeList );
  137.     bzero ( Sample, sizeof(SampleData) );
  138.     return Sample;
  139.   }
  140.  
  141.   if (nFreeSamples <= 0) {
  142.     FreeSamples = (SampleData *) calloc( GLOB_SIZE, sizeof(SampleData));
  143.     if(FreeSamples == NULL)
  144.       {
  145.     fprintf(stderr,"Fatal error - Can't allocate sample data\n");
  146.     abort();
  147.       }
  148.     nFreeSamples = GLOB_SIZE;
  149.   }
  150.   return & (FreeSamples[--nFreeSamples]);
  151. }
  152.  
  153. /* Split a node into two subnodes */
  154. SplitNode(Node)
  155.      HierarchicalRegion *Node;
  156. {
  157.   Node->Left = CreateNode();
  158.   Node->Right = CreateNode();
  159. }
  160.  
  161.  
  162. /*-----------------------------
  163.   propagate the current node's
  164.   sample to one of the subnodes
  165.   return the EMPTY SUBNODE
  166.   ---------------------------*/
  167. HierarchicalRegion *PropagateSample(Node, Boundary, NewSample)
  168. HierarchicalRegion *Node;
  169. RectType *Boundary;
  170. SampleData *NewSample;
  171. {
  172.   RectType LeftBoundary;
  173.   double dx, dy;
  174.  
  175.   if (AdaptiveSplits) {
  176.  
  177.     /** Choose a split direction which separates the new sample
  178.      ** from the old and maintains an aspect ratio as close to square
  179.      ** as possible for the new sub cells.
  180.      */
  181.     double xsplit, ysplit, x1, x2, y1, y2, xratio, xratio2, yratio, yratio2;
  182.     double width, height;
  183.     
  184.  
  185.     /* Split at the mid point between samples */
  186.     xsplit = 0.5*(Node->Sample->Xs + NewSample->Xs);
  187.     ysplit = 0.5*(Node->Sample->Ys + NewSample->Ys);
  188.  
  189.     /* These are the side lengths of the sub cells */
  190.     x1 = xsplit - Boundary->xL;
  191.     x2 = Boundary->xH - xsplit;
  192.     y1 = ysplit - Boundary->yL;
  193.     y2 = Boundary->yH - ysplit;
  194.  
  195.     /* These are the side lengths of the whole cell */
  196.     width  = (Boundary->xH - Boundary->xL);
  197.     height = (Boundary->yH - Boundary->yL);
  198.  
  199.     /* These are the aspect ratios if we split in x */
  200.     if (x1 > height)
  201.       xratio =  (height < 1.0E-20) ? 1.0E20 : x1 / height;
  202.     else
  203.       xratio =  (x1     < 1.0E-20) ? 1.0E20 : height / x1;
  204.     if (x2 > height)
  205.       xratio2 =  (height < 1.0E-20) ? 1.0E20 : x2 / height;
  206.     else
  207.       xratio2 =  (x2     < 1.0E-20) ? 1.0E20 : height / x2;
  208.  
  209.     
  210.     /* Use the worst of the two for comparison */
  211.     if (xratio2 > xratio)  xratio = xratio2;
  212.  
  213.     /* These are the aspect ratios if we split in y */
  214.     if (y1 > width)
  215.       yratio =  (width < 1.0E-20) ? 1.0E20 : y1 / width;
  216.     else
  217.       yratio =  (y1     < 1.0E-20) ? 1.0E20 : width / y1;
  218.     if (y2 > width)
  219.       yratio2 =  (width < 1.0E-20) ? 1.0E20 : y2 / width;
  220.     else
  221.       yratio2 =  (y2     < 1.0E-20) ? 1.0E20 : width / y2;
  222.  
  223.  
  224.     /* Use the worst of the two for comparison */
  225.     if (yratio2 > yratio)  yratio = yratio2;
  226.  
  227.     if (xratio < yratio) {  /* X split is better */
  228.       Node->splitDirection = 1;    /* Split in X */
  229.       Node->splitValue = 0.5*(Node->Sample->Xs + NewSample->Xs);
  230.     } else  {
  231.       Node->splitDirection = 0;    /* Split in Y */
  232.       Node->splitValue = 0.5*(Node->Sample->Ys + NewSample->Ys);
  233.     }
  234.   } else {
  235.     Node->splitDirection = (level & 1);  /* Odd levels split in X */
  236.     if (Node->splitDirection)
  237.       Node->splitValue = 0.5*(Boundary->xL + Boundary->xH);
  238.     else
  239.       Node->splitValue = 0.5*(Boundary->yL + Boundary->yH);
  240.   }
  241.     
  242.   SplitBoundary( TRUE, Node->splitDirection, Node->splitValue,
  243.           Boundary, &LeftBoundary );
  244.   
  245.   if( (Node->Sample->Xs >=   LeftBoundary.xL)  &&
  246.       (Node->Sample->Xs <    LeftBoundary.xH)  &&
  247.       (Node->Sample->Ys >=   LeftBoundary.yL)  &&
  248.       (Node->Sample->Ys <    LeftBoundary.yH) )
  249.     {
  250.       Node->Left->Sample    = Node->Sample;
  251.       return Node->Right; 
  252.     }
  253.   else
  254.     {
  255.       Node->Right->Sample    = Node->Sample;
  256.       return Node->Left;
  257.     } 
  258. }
  259.  
  260.  
  261. /*---------------------------------
  262.   Insert the new node in the tree.
  263.  ----------------------------------*/
  264. InsertNode(NewSample, Root, RootBoundary, InterestingBoundary)
  265.      SampleData *NewSample;
  266.      HierarchicalRegion *Root;
  267.      RectType *RootBoundary, *InterestingBoundary;
  268. {
  269.   HierarchicalRegion *Node = Root, *EmptyNode;
  270.   RectType Boundarys[2], *Boundary, *NewBoundary;
  271.   int old, new;
  272.   
  273.  
  274.   Boundarys[0] = *RootBoundary;
  275.   old = 0; new = 1;
  276.   Boundary = Boundarys+0;
  277.   NewBoundary = Boundarys+1;
  278.  
  279.   /* Search down to the leaf level for the cell this sample belongs in */
  280.   level = 0;
  281.   while (!IsLeaf(Node))
  282.     {
  283.       if (!RectOverlap( InterestingBoundary, Boundary )) {
  284.     /* Not in an interesting part of the image, ignore it */
  285.     FreeSample( NewSample );
  286.     return;
  287.       }
  288.       SplitBoundary( TRUE, Node->splitDirection, Node->splitValue,
  289.             Boundary, NewBoundary );
  290.       if (PointInRectangle( NewSample->Xs, NewSample->Ys, NewBoundary )) {
  291.     Node = Node->Left;
  292.       } else {
  293.     SplitBoundary( FALSE, Node->splitDirection, Node->splitValue,
  294.               Boundary, NewBoundary );
  295.     Node = Node->Right;
  296.       }
  297.       old = (old +1) %2;
  298.       new = (new +1) %2;
  299.       Boundary    = Boundarys+old;
  300.       NewBoundary = Boundarys+new;
  301.       level++;
  302.     }
  303.   
  304.   /* Split the Node */
  305.   SplitNode( Node );
  306.  
  307.   /* Fill in the two new cells */
  308.   EmptyNode = PropagateSample(Node, Boundary, NewSample);
  309.   EmptyNode->Sample = NewSample;
  310.  
  311.   AnimateShadeCell( Boundary, NewSample );
  312. }
  313.  
  314.  
  315. /*-----------------------
  316.   Initialize samples
  317.   ----------------------*/
  318. InitializeTree( xL, yL, xH, yH )
  319. float xL, yL, xH, yH;
  320. {
  321.   RootBoundary.xL = xL;
  322.   RootBoundary.yL = yL;
  323.   RootBoundary.xH = xH;
  324.   RootBoundary.yH = yH;
  325.   Root.Left = Root.Right = (HierarchicalRegion *) NULL;
  326.   Root.Sample = NULL;
  327. }
  328.  
  329. AddSampleToTree( x, y, index )
  330. float x,y;
  331. long index;
  332. {
  333.   if (Root.Sample == NULL)
  334.     {
  335.       Root.Sample = CreateSample();
  336.       Root.Sample->Xs = x;
  337.       Root.Sample->Ys = y;
  338.       Root.Sample->index = index;
  339.     }
  340.   else
  341.     {
  342.       SampleData *Sample = CreateSample();
  343.       Sample->Xs = x;
  344.       Sample->Ys = y;
  345.       Sample->index = index;
  346.       /* Insert this node in the tree */
  347.       InsertNode( Sample, &Root, &RootBoundary, &RootBoundary );
  348.     }
  349. }
  350.  
  351.  
  352.  
  353.  
  354.