home *** CD-ROM | disk | FTP | other *** search
/ Power-Programmierung / CD2.mdf / c / library / dos / diverses / leda / src / plane / _sweep_s.c < prev    next >
Encoding:
C/C++ Source or Header  |  1991-11-15  |  13.2 KB  |  586 lines
  1. /*******************************************************************************
  2. +
  3. +  LEDA  2.1.1                                                 11-15-1991
  4. +
  5. +
  6. +  _sweep_segments.c
  7. +
  8. +
  9. +  Copyright (c) 1991  by  Max-Planck-Institut fuer Informatik
  10. +  Im Stadtwald, 6600 Saarbruecken, FRG     
  11. +  All rights reserved.
  12. + 
  13. *******************************************************************************/
  14.  
  15.  
  16. //------------------------------------------------------------------------------
  17. //
  18. // SWEEP_SEGMENTS
  19. //
  20. // a plane sweep algorithm for computing line segment intersections
  21. //
  22. // SWEEP_SEGMENTS(list(segment)& L1, list(segment)& L2, GRAPH(point,int) G);
  23. //
  24. // computes the planar subdivision created by the (directed) straight
  25. // line segments of L1 and L2
  26. // 
  27. // G = (V,E) with
  28. //
  29. // 1. each node v in V contains a point of intersectoin between two segments
  30. //    (from L1 or L2)
  31. //
  32. // 2. each edge e = (v,w) corresponds to a
  33. //
  34. //
  35. //------------------------------------------------------------------------------
  36.  
  37.  
  38. #include <LEDA/plane.h>
  39. #include <LEDA/plane_alg.h>
  40.  
  41. #include <LEDA/sortseq.h>
  42. #include <LEDA/prio.h>
  43.  
  44. #include <math.h>
  45.  
  46. #define EPS  0.00001
  47. #define EPS2 0.0000000001
  48.  
  49. class S_POINT;
  50. class S_SEGMENT;
  51. typedef S_POINT* S_point;
  52. typedef S_SEGMENT* S_segment;
  53.  
  54. enum S_point_type {Cross=0,Rightend=1,Leftend=2}; 
  55.  
  56. static class S_POINT
  57. {
  58.   friend class S_SEGMENT;
  59.  
  60.   S_segment seg;
  61.   int     kind;
  62.   real  x;
  63.   real  y;
  64.  
  65. //size = 4 words
  66.  
  67.   public:
  68.  
  69.   S_POINT(real a,real b)  
  70.   { 
  71.     x=a; y=b; seg=0; kind=Cross;
  72.    }
  73.  
  74.  
  75.  
  76.   S_POINT(point p)         
  77.   { 
  78.     x=p.xcoord();y=p.ycoord();seg=0;kind=Cross;
  79.    }
  80.  
  81.  
  82. OPERATOR_NEW(4)
  83. OPERATOR_DEL(4)
  84.  
  85.   friend real    get_x(S_point p)      { return p->x; }
  86.   friend real    get_y(S_point p)      { return p->y; }
  87.   friend int     get_kind(S_point p)   { return p->kind; }
  88.   friend S_segment get_seg(S_point p)    { return p->seg; }   
  89.  
  90.   friend bool intersection(S_segment, S_segment, S_point&);
  91. };
  92.  
  93. static int compare(S_point& p1,S_point& p2)
  94. { if (p1==p2) return 0;
  95.  
  96.   real diffx = get_x(p1) - get_x(p2);
  97.   if (fabs(diffx) > EPS2 ) return (diffx > 0.0) ? 1 : -1;
  98.  
  99.   int    diffk = get_kind(p1)-get_kind(p2);
  100.   if (diffk != 0) return diffk;
  101.  
  102.   real diffy = get_y(p1) - get_y(p2);
  103.   if (fabs(diffy) > EPS2 ) return (diffy > 0.0) ? 1 : -1;
  104.  
  105.   return 0;
  106.  
  107. }
  108.  
  109.  
  110. static class S_SEGMENT
  111. {
  112.   S_point startpoint;
  113.   S_point endpoint;
  114.   real  slope;
  115.   real  yshift;
  116.   node  left_node;
  117.   int   orient;
  118.   int   color;
  119.   int   name;
  120.  
  121.  
  122. // size = 8 words
  123.  
  124.   public:
  125.  
  126.   S_SEGMENT(S_point, S_point,int,int);     
  127.  
  128.  ~S_SEGMENT() { delete startpoint; delete endpoint; }     
  129.  
  130. OPERATOR_NEW(8)
  131. OPERATOR_DEL(8)
  132.  
  133.  
  134.   friend S_point get_startpoint(S_segment seg)       { return seg->startpoint; }
  135.   friend S_point get_endpoint(S_segment seg)         { return seg->endpoint; }
  136.   friend real get_slope(S_segment seg)             { return seg->slope; }
  137.   friend real get_yshift(S_segment seg)            { return seg->yshift; }
  138.   friend int  get_orient(S_segment seg)            { return seg->orient; }
  139.   friend int  get_color(S_segment seg)             { return seg->color; }
  140.   friend int  get_name(S_segment seg)              { return seg->name; }
  141.   friend node get_left_node(S_segment seg)         { return seg->left_node; }
  142.   friend void set_left_node(S_segment seg, node v) { seg->left_node = v; }
  143.  
  144.   friend bool intersection(S_segment, S_segment, S_point&);
  145. };
  146.  
  147.  
  148. S_SEGMENT::S_SEGMENT(S_point p1,S_point p2,int c, int n)    
  149.   {
  150.     left_node  = nil;
  151.     color      = c;
  152.     name       = n;
  153.  
  154.     if (compare(p1,p2) < 0)
  155.      { startpoint = p1; 
  156.        endpoint = p2; 
  157.        orient = 0;
  158.       }
  159.     else
  160.      { startpoint = p2; 
  161.        endpoint = p1; 
  162.        orient = 1;
  163.       }
  164.  
  165.     startpoint->kind = Leftend; 
  166.     endpoint->kind = Rightend; 
  167.     startpoint->seg = this; 
  168.     endpoint->seg = this;
  169.  
  170.     if (endpoint->x != startpoint->x)
  171.     {
  172.       slope = (endpoint->y - startpoint->y)/(endpoint->x - startpoint->x);
  173.       yshift = startpoint->y - slope * startpoint->x;
  174.  
  175.       startpoint->x -= EPS;
  176.       startpoint->y -= EPS * slope;
  177.       endpoint->x += EPS;
  178.       endpoint->y += EPS * slope;
  179.     }
  180.     else //vertical segment
  181.     { startpoint->y -= EPS;
  182.       endpoint->y   += EPS;
  183.      }
  184.   }
  185.  
  186.  
  187. real x_sweep;
  188. real y_sweep;
  189.  
  190.  
  191. static int compare(S_segment& s1,S_segment& s2)
  192. {
  193.   real y1 = get_slope(s1)*x_sweep+get_yshift(s1);
  194.   real y2 = get_slope(s2)*x_sweep+get_yshift(s2);
  195.  
  196.   real diff = y1-y2;
  197.  
  198.   if (fabs(diff) > EPS2) return (diff > 0.0) ? 1 : -1;
  199.  
  200.   if (get_slope(s1) == get_slope(s2)) 
  201.         return compare(get_x(get_startpoint(s1)), get_x(get_startpoint(s2)));
  202.  
  203.   if (y1 <= y_sweep+EPS2)
  204.         return compare(get_slope(s1),get_slope(s2));
  205.   else
  206.         return compare(get_slope(s2),get_slope(s1));
  207.  
  208. }
  209.  
  210.  
  211. void Print(S_segment& x) 
  212. { S_point s = get_startpoint(x);
  213.   S_point e = get_endpoint(x);
  214.   cout << 
  215.     string("[(%f,%f)----(%f,%f)]",~get_x(s),~get_y(s), ~get_x(e),~get_y(e));
  216. }
  217.  
  218.  
  219. declare2(priority_queue,seq_item,S_point);
  220. priority_queue(seq_item,S_point) X_structure;
  221.  
  222. declare2(sortseq,S_segment,pq_item);
  223. sortseq(S_segment,pq_item) Y_structure;
  224.  
  225.  
  226. bool intersection(S_segment seg1,S_segment seg2, S_point& inter)
  227. {
  228.   if (seg1->slope == seg2->slope)
  229.     return false;
  230.   else
  231.   {
  232.     //x-coordinate  of the intersection
  233.     real Cross_x = (seg2->yshift - seg1->yshift) / (seg1->slope - seg2->slope);
  234.  
  235.     if (Cross_x <= x_sweep) return false;
  236.  
  237.     real s1 = seg1->startpoint->x;
  238.     real s2 = seg2->startpoint->x;
  239.     real e1 = seg1->endpoint->x;
  240.     real e2 = seg2->endpoint->x;
  241.  
  242.     if (s2>Cross_x || s1>Cross_x || Cross_x>e1 || Cross_x>e2) return false;
  243.  
  244.     //y-coordinate of the intersection
  245.     real Cross_y = (seg1->slope * Cross_x + seg1->yshift);
  246.     inter =  new S_POINT(Cross_x,Cross_y);
  247.  
  248.     return true;
  249.   }
  250. }
  251.  
  252.  
  253. static pq_item Xinsert(seq_item i, S_point p) 
  254. { 
  255.   return X_structure.insert(i,p);
  256. }
  257.  
  258. static S_point Xdelete(pq_item i) 
  259. {
  260.   S_point p = X_structure.inf(i);
  261.   X_structure.del_item(i);
  262.   return p;
  263. }
  264.  
  265.  
  266. static node New_Node(GRAPH(point,int)& G,real x, real y )
  267. { return G.new_node(point(x,y)); }
  268.  
  269. static void New_Edge(GRAPH(point,int)& G,node v, node w, S_segment l )
  270. { if (get_orient(l)==0)
  271.        G.new_edge(v,w,get_color(l));
  272.   else G.new_edge(w,v,get_color(l));
  273. }
  274.  
  275.  
  276. void handle_vertical_segment(GRAPH(point,int)& SUB, S_segment l)
  277. { 
  278.   S_point p = new S_POINT(get_x(get_startpoint(l)),get_y(get_startpoint(l)));
  279.   S_point q = new S_POINT(get_x(get_endpoint(l)),get_y(get_endpoint(l)));
  280.  
  281.   S_point r = new S_POINT(get_x(p)+1,get_y(p));
  282.   S_point s = new S_POINT(get_x(q)+1,get_y(q));
  283.  
  284.   S_segment bot = new S_SEGMENT(p,r,0,0);
  285.   S_segment top = new S_SEGMENT(q,s,0,0);
  286.  
  287.   seq_item bot_it = Y_structure.insert(bot,nil);
  288.   seq_item top_it = Y_structure.insert(top,nil);
  289.   seq_item sit;
  290.  
  291.   node u,v,w;
  292.   S_segment seg;
  293.   
  294.  
  295.   for(sit=Y_structure.succ(bot_it); sit != top_it; sit=Y_structure.succ(sit))
  296.   { seg = Y_structure.key(sit);
  297.  
  298.     real cross_y = (get_slope(seg) * get_x(p) + get_yshift(seg));
  299.  
  300.     v = get_left_node(seg);
  301.     if (v==nil)
  302.     { w = New_Node(SUB,get_x(p),cross_y);
  303.       set_left_node(seg,w);
  304.      }
  305.     else
  306.     { real vx = SUB[v].xcoord();
  307.       if ( vx < get_x(p)-EPS2) 
  308.       { w = New_Node(SUB,get_x(p),cross_y);
  309.         New_Edge(SUB,v,w,seg);
  310.         set_left_node(seg,w);
  311.        }
  312.       else w = v;
  313.      }
  314.  
  315.     u = get_left_node(l);
  316.     if (u!=nil && u!=w) New_Edge(SUB,u,w,l);
  317.     set_left_node(l,w);
  318.  
  319.    }
  320.   
  321.   delete l;
  322.   delete top;
  323.   delete bot;
  324.     
  325.   Y_structure.del_item(bot_it);
  326.   Y_structure.del_item(top_it);
  327.  
  328.  }
  329.  
  330.  
  331. void SWEEP_SEGMENTS(list(segment)& S, list(point)& P)
  332. { GRAPH(point,int) G;
  333.   list(segment) S0;
  334.  
  335.   SWEEP_SEGMENTS(S,S0,G);
  336.  
  337.   node v;
  338.   forall_nodes(v,G) P.append(G[v]);
  339.  
  340. }
  341.  
  342. void SWEEP_SEGMENTS(list(segment)& L1, list(segment)& L2, GRAPH(point,int)& SUB)
  343. {
  344.   S_point    p,inter;
  345.   S_segment  seg, l,lsit,lpred,lsucc,lpredpred;
  346.   pq_item  pqit,pxmin;
  347.   seq_item sitmin,sit,sitpred,sitsucc,sitpredpred;
  348.  
  349.  
  350.   int count=1;
  351.  
  352.   //initialization of the X-structure
  353.  
  354.   segment s;
  355.  
  356.   forall(s,L1) 
  357.    { S_point p = new S_POINT(s.start());
  358.      S_point q = new S_POINT(s.end());
  359.      seg = new S_SEGMENT(p,q,0,count++);
  360.      Xinsert(nil,get_startpoint(seg));
  361.    }
  362.  
  363.   count = -1;
  364.  
  365.   forall(s,L2) 
  366.    { S_point p = new S_POINT(s.start());
  367.      S_point q = new S_POINT(s.end());
  368.      seg = new S_SEGMENT(p,q,1,count--);
  369.      Xinsert(nil,get_startpoint(seg));
  370.    }
  371.  
  372.  
  373.   count = 0;
  374.  
  375.   x_sweep = -Infinity;
  376.   y_sweep = -Infinity;
  377.  
  378.  
  379.   while(!X_structure.empty())
  380.   {
  381.     pxmin = X_structure.find_min();
  382.     p = X_structure.inf(pxmin);
  383.  
  384.     sitmin = X_structure.key(pxmin);
  385.  
  386.     Xdelete(pxmin);
  387.  
  388.  
  389.  
  390.     if (get_kind(p) == Leftend)
  391.  
  392.     //left endpoint
  393.     { 
  394.  
  395.       l = get_seg(p); 
  396.  
  397.       x_sweep = get_x(p);
  398.       y_sweep = get_y(p);
  399.  
  400.  
  401.       if (get_x(p) == get_x(get_endpoint(l)))
  402.         handle_vertical_segment(SUB,l);
  403.       else
  404.       {
  405.  
  406. /*
  407.       sit = Y_structure.lookup(l);
  408.       if (sit!=nil)  
  409.            error_handler(1,"plane sweep: sorry, overlapping segments");
  410. */
  411.  
  412.  
  413.       sit = Y_structure.insert(l,nil);
  414.  
  415.       Xinsert(sit,get_endpoint(l));
  416.  
  417.       sitpred = Y_structure.pred(sit);
  418.       sitsucc = Y_structure.succ(sit);
  419.  
  420.       if (sitpred != nil) 
  421.       { if ((pqit = Y_structure.inf(sitpred)) != nil)
  422.           delete Xdelete(pqit);
  423.  
  424.         lpred = Y_structure.key(sitpred);
  425.  
  426.         Y_structure.change_inf(sitpred,nil);
  427.  
  428.         if (intersection(lpred,l,inter))
  429.             Y_structure.change_inf(sitpred,Xinsert(sitpred,inter));
  430.       }
  431.  
  432.  
  433.       if (sitsucc != nil)
  434.       { lsucc = Y_structure.key(sitsucc);
  435.         if (intersection(lsucc,l,inter))
  436.            Y_structure.change_inf(sit,Xinsert(sit,inter));
  437.       }
  438.      } /* else if vertical */
  439.  
  440.     }
  441.     else if (get_kind(p) == Rightend)
  442.          //right endpoint
  443.          { 
  444.            x_sweep = get_x(p);
  445.            y_sweep = get_y(p);
  446.  
  447.            sit = sitmin;
  448.  
  449.            sitpred = Y_structure.pred(sit);
  450.            sitsucc = Y_structure.succ(sit);
  451.  
  452.            S_segment seg = Y_structure.key(sit);
  453.  
  454.            Y_structure.del_item(sit);
  455.  
  456.            delete seg;
  457.  
  458.            if((sitpred != nil)&&(sitsucc != nil))
  459.            {
  460.              lpred = Y_structure.key(sitpred);
  461.              lsucc = Y_structure.key(sitsucc);
  462.              if (intersection(lsucc,lpred,inter))
  463.                 Y_structure.change_inf(sitpred,Xinsert(sitpred,inter));
  464.            }
  465.          }
  466.          else
  467.          /*point of intersection*/
  468.          { 
  469.            node w = New_Node(SUB,get_x(p),get_y(p));
  470.  
  471.            count++;
  472.  
  473.            /* Let L = list of all lines intersecting in p 
  474.  
  475.               we compute sit     = L.head();
  476.               and        sitpred = L.tail();
  477.  
  478.               by scanning the Y_structure in both directions 
  479.               starting at sitmin;
  480.  
  481.            */
  482.  
  483.            /* search for sitpred upwards from sitmin: */
  484.  
  485.            Y_structure.change_inf(sitmin,nil);
  486.  
  487.            sitpred = Y_structure.succ(sitmin);
  488.  
  489.  
  490.            while ((pqit=Y_structure.inf(sitpred)) != nil)
  491.            { S_point q = X_structure.inf(pqit);
  492.              if (compare(p,q) != 0) break; 
  493.              X_structure.del_item(pqit);
  494.              Y_structure.change_inf(sitpred,nil);
  495.              sitpred = Y_structure.succ(sitpred);
  496.             }
  497.  
  498.  
  499.            /* search for sit downwards from sitmin: */
  500.  
  501.            sit = sitmin;
  502.  
  503.            seq_item sit1;
  504.            
  505.            while ((sit1=Y_structure.pred(sit)) != nil)
  506.            { pqit = Y_structure.inf(sit1);
  507.              if (pqit == nil) break;
  508.              S_point q = X_structure.inf(pqit);
  509.              if (compare(p,q) != 0) break; 
  510.              X_structure.del_item(pqit);
  511.              Y_structure.change_inf(sit1,nil);
  512.              sit = sit1;
  513.             }
  514.  
  515.  
  516.  
  517.            // insert edges to p for all S_segments in sit, ..., sitpred into SUB
  518.            // and set left node to w 
  519.  
  520.            lsit = Y_structure.key(sitpred);
  521.  
  522.            node v = get_left_node(lsit);
  523.            if (v!=nil) New_Edge(SUB,v,w,lsit);
  524.            set_left_node(lsit,w);
  525.  
  526.            for(sit1=sit; sit1!=sitpred; sit1 = Y_structure.succ(sit1))
  527.            { lsit = Y_structure.key(sit1);
  528.  
  529.              v = get_left_node(lsit);
  530.              if (v!=nil) New_Edge(SUB,v,w,lsit);
  531.              set_left_node(lsit,w);
  532.             }
  533.  
  534.            lsit = Y_structure.key(sit);
  535.            lpred=Y_structure.key(sitpred);
  536.            sitpredpred = Y_structure.pred(sit);
  537.            sitsucc=Y_structure.succ(sitpred);
  538.  
  539.  
  540.            if (sitpredpred != nil)
  541.             { 
  542.               lpredpred=Y_structure.key(sitpredpred);
  543.  
  544.               if ((pqit = Y_structure.inf(sitpredpred)) != nil)
  545.                 delete Xdelete(pqit);
  546.  
  547.               Y_structure.change_inf(sitpredpred,nil);
  548.  
  549.  
  550.               if (intersection(lpred,lpredpred,inter))
  551.                 Y_structure.change_inf(sitpredpred,
  552.                                        Xinsert(sitpredpred,inter));
  553.              }
  554.  
  555.  
  556.            if (sitsucc != nil)
  557.             {
  558.               lsucc=Y_structure.key(sitsucc);
  559.  
  560.               if ((pqit = Y_structure.inf(sitpred)) != nil)
  561.                 delete Xdelete(pqit);
  562.                  
  563.               Y_structure.change_inf(sitpred,nil);
  564.  
  565.               if (intersection(lsucc,lsit,inter))
  566.                   Y_structure.change_inf(sit,Xinsert(sit,inter));
  567.              }
  568.  
  569.  
  570. // reverse the subsequence sit, ... ,sitpred  in the Y_structure
  571.  
  572.            x_sweep = get_x(p);
  573.            y_sweep = get_y(p);
  574.  
  575.            Y_structure.reverse_items(sit,sitpred);
  576.  
  577.           delete p;
  578.  
  579.          } // intersection
  580.  
  581.   }
  582.  
  583.   X_structure.clear();
  584.  
  585. }
  586.