home *** CD-ROM | disk | FTP | other *** search
/ OS/2 Shareware BBS: 10 Tools / 10-Tools.zip / ledar34.zip / leda-r-3_4_tar / LEDA-3.4 / src / plane_alg / balaban.c < prev    next >
C/C++ Source or Header  |  1996-09-03  |  14KB  |  557 lines
  1. /*******************************************************************************
  2. +
  3. +  LEDA 3.4
  4. +
  5. +  balaban.c
  6. +
  7. +  This file is part of the LEDA research version (LEDA-R) that can be 
  8. +  used free of charge in academic research and teaching. Any commercial
  9. +  use of this software requires a license which is distributed by the
  10. +  LEDA Software GmbH, Postfach 151101, 66041 Saarbruecken, FRG
  11. +  (fax +49 681 31104).
  12. +
  13. +  Copyright (c) 1991-1996  by  Max-Planck-Institut fuer Informatik
  14. +  Im Stadtwald, 66123 Saarbruecken, Germany     
  15. +  All rights reserved.
  16. + 
  17. *******************************************************************************/
  18. #include <LEDA/list.h>
  19. #include <LEDA/slist.h>
  20. #include <LEDA/sortseq.h>
  21. #include <LEDA/map.h>
  22.  
  23. // Problem :
  24. // Find all $K$ intersecting pairs in a set $S$ consisting of $n$ straight-line 
  25. // segments in the plane.
  26.  
  27. // Algorithm :
  28. // Literature : Ivan J. Balanban ???? S.211 - 219
  29.  
  30. // The following implementation should make no assumptions about the input :
  31. // -  segments may have length zero                              
  32. // -  segments may be vertical                                   
  33. // -  end  and intersection points may have the same abscissae   
  34. // -  end  and intersection points may have the same coordinates 
  35.  
  36. typedef void (*rep_int_func) (const SEGMENT&, const SEGMENT&);
  37.  
  38.  
  39.  
  40. // |stairs| : pointer to an array that implements the "staircase recursion
  41. //            stack"
  42. //            $O(n)$ space
  43.  
  44. static SEGMENT*  stairs = 0;
  45.  
  46. // |StairCount| : integer that stores the number of stairs that currently 
  47. // stored in |stairs|
  48.  
  49. static int      StairCount = 0;
  50.  
  51.  
  52.  
  53. // Classes :
  54.  
  55. class STAIRCASE
  56. {  
  57.   int   bcount;
  58.   int   ecount;
  59.   
  60.   public :
  61.   
  62.   int   b;   
  63.   int   e;
  64.  
  65.   STAIRCASE(int p1, int p2) : bcount(StairCount), ecount(0), b(p1), e(p2) {}
  66.   
  67.   int    Count() const { return ecount; }
  68.   bool   empty() const { return ecount == 0; }
  69.  
  70.   SEGMENT& operator[] (int i) const { return stairs[bcount+i]; }
  71.   
  72.   void append(SEGMENT& s)
  73.   { stairs[++StairCount] = s;
  74.     ecount++;
  75.    }
  76.  
  77. };
  78.  
  79.  
  80.  
  81.  
  82.  
  83. // function prototyps :
  84.  
  85. void TreeSearch_nloghoch2n(slist<SEGMENT>& L, slist<SEGMENT>& I, int b, int e,
  86.                slist<SEGMENT>& R);
  87. void SearchInStrip(slist<SEGMENT>& L, slist<SEGMENT>& R, int b, int e);
  88.  
  89. void InnerIntPoints(const STAIRCASE& D, const slist<SEGMENT>& I);
  90. int  Loc(const STAIRCASE& D, const SEGMENT s);
  91. void ReportIntersection(const SEGMENT& s1, const SEGMENT& s2);
  92.  
  93. void Split(slist<SEGMENT>& L, STAIRCASE& D);
  94. void Lrs_from_Rls(slist<SEGMENT>& Rls, slist<SEGMENT>& Lrs, int c);
  95. void Merge(STAIRCASE& D, slist<SEGMENT>& R);
  96.  
  97. // GlobalVar
  98.  
  99. // |report_func| : function to report intersecting pairs
  100. static rep_int_func              report_func = 0;
  101.  
  102. // |original| : pointer to a map that associates each input segment with an 
  103. //              internally used segment, which is directed from left to rigth 
  104. //              or upwards
  105. //              $O(n)$ space 
  106. static map<SEGMENT,SEGMENT>*     original    = 0; 
  107.  
  108. // |init_points| : pointer to an array that stores the different endpoints
  109. //                 of the input segments ordered lexicographically
  110. //                 $O(n)$ space
  111. static POINT*    init_points = 0;
  112.  
  113. // |seg_queue| : list of internal segment sorted, using the ordering defined  
  114. //               by the function |compare2(const SEGMENT&, const SEGMENT&)| 
  115. //               $O(n)$ space
  116. static list<SEGMENT>    seg_queue;
  117.  
  118. static int compare2(const SEGMENT& s1, const SEGMENT& s2)
  119. {
  120.   int c = compare(s1.start(), s2.start());
  121.   if (c) return c;
  122.   return cmp_slopes(s1, s2);
  123. }
  124.  
  125. void balaend(const list<SEGMENT>&)
  126. { // free used memory
  127.   
  128.   if(init_points) delete[] init_points; init_points = 0;
  129.   if(stairs)      delete[] stairs;      stairs = 0;
  130.   if(original)    delete original;      original = 0;
  131.   
  132.   seg_queue.clear();
  133. }
  134.  
  135. void BALABAN_SEGMENTS(const list<SEGMENT>& S, rep_int_func rif)
  136. {
  137.   report_func  = rif;
  138.  
  139.   if ( S.length() )   // no segments => no intersecting pairs
  140.     {
  141.       // allocate space for |original|
  142.       original    = new map<SEGMENT,SEGMENT>;
  143.       if( !original )  return;
  144.       
  145.       // sort the endpoints of the segments in $S$ lexico
  146.       sortseq <POINT, int> ps_queue;
  147.       SEGMENT s;
  148.       seq_item sit1, sit2;
  149.       
  150.       forall(s, S)
  151.     {
  152.       sit1 = ps_queue.insert(s.start(), 0);
  153.       sit2 = ps_queue.insert(s.end(), 0); 
  154.  
  155.       SEGMENT slex( ps_queue.key(sit1), ps_queue.key(sit2) );
  156.       
  157.       // adjust segments from left to rigth or upwards
  158.       if ( compare(slex.start(), slex.end()) > 0 ) 
  159.         slex = SEGMENT( slex.end(), slex.start() );
  160.  
  161.       (*original)[slex] = s;
  162.       seg_queue.append(slex);
  163.     }
  164.     
  165.       // sort |seg_queue| using the ordering defined by the 
  166.       // function |compare2(const SEGMENT&, const SEGMENT&)|
  167.       seg_queue.sort(compare2);
  168.  
  169.       int N = ps_queue.size();
  170.       StairCount = 0;
  171.  
  172.       // allocate space for |init_points| and |stairs|
  173.       init_points = new POINT[N+1];
  174.       stairs      = new SEGMENT[S.length()+1];
  175.  
  176.       if(!init_points || !stairs )  { balaend(S); return; }
  177.  
  178.       int n = 1;
  179.  
  180.       // initialize |init_points|
  181.       while (!ps_queue.empty())
  182.     {
  183.       sit1 = ps_queue.min();
  184.       init_points[n++]  = ps_queue.key(sit1);
  185.       ps_queue.del_item(sit1);
  186.     }
  187.       
  188.       slist<SEGMENT> L; 
  189.       slist<SEGMENT> I;
  190.       slist<SEGMENT> R;
  191.  
  192.       // initialize L 
  193.       Lrs_from_Rls(R, L, 1);
  194.       
  195.       // initialize I
  196.       forall( s, seg_queue )
  197.     if ( !identical(s.end(), init_points[N]) )
  198.       I.append(s);
  199.       
  200.       // start $TreeSearch$
  201.       if( N > 1 )
  202.     {
  203.       TreeSearch_nloghoch2n(L, I, 1, N, R);
  204.       Lrs_from_Rls(R, L, N);
  205.     }
  206.     }
  207.   
  208.   balaend(S);
  209. }
  210.  
  211. void TreeSearch_nloghoch2n(slist<SEGMENT>& Lv, slist<SEGMENT>& Iv,
  212.                int b, int e, slist<SEGMENT>& Rv)
  213. {   
  214.   // 1. if e-b = 1 then SearchInStripbe(Lv, Rv); exit;
  215.   
  216.   if( e-b == 1 )  { SearchInStrip(Lv, Rv, b, e); return; }
  217.   
  218.   // 2. Splitbe(Lv, Qv, Lls); Dv := (Qv, <b,e>);
  219.   //
  220.   //        and
  221.   // 
  222.   // 3. Find Int(Dv, Lls);
  223.   
  224.   STAIRCASE     Dv(b, e);
  225.   
  226.   Split(Lv, Dv);
  227.  
  228.   // 10. For s ele Iv do find Loc(Dv, s) endfor
  229.   // 11. Find Int(Dv, Iv)
  230.   
  231.   InnerIntPoints(Dv, Iv);
  232.  
  233.   // 4. c := (b+e)/2;
  234.   
  235.   int   c = (b+e)/2;  
  236.   POINT pc = init_points[c];
  237.   
  238.   // 5. Place segments of Iv
  239.   //     inner for <b,c> into Ils
  240.   //     inner for <c,e> into Irs
  241.  
  242.   slist<SEGMENT> Ils;
  243.   slist<SEGMENT> Irs;
  244.   
  245.   while (!Iv.empty())
  246.     {
  247.       if( compare( Iv.head().end(), pc) < 0 ) // inner <b,c>
  248.     Ils.append(Iv.pop());
  249.       else
  250.     if( compare(Iv.head().start(), pc) > 0 ) // inner <c, pe>
  251.       Irs.conc(Iv);
  252.     else 
  253.       Iv.pop();
  254.     }
  255.   
  256.   // 6. TreeSearch(Lls, Ils, b, c, Rls);
  257.  
  258.   TreeSearch_nloghoch2n( Lv, Ils, b, c, Rv);
  259.  
  260.   // 7. If pc is left endpoint of s(c) Lrs := insert s(c) into Rls
  261.   //                           else    Lrs := remove s(c) from Rls
  262.   
  263.   Lrs_from_Rls(Rv, Lv, c); 
  264.  
  265.   // 8. TreeSearch(Lrs, Irs, c, e, Rrs);
  266.  
  267.   TreeSearch_nloghoch2n( Lv, Irs, c, e, Rv);
  268.     
  269.   //  9. Find Int(Dv, Rrs)
  270.   // 12. Rv := Mergee(Qv, Rrs)
  271.   
  272.   Merge(Dv, Rv);
  273.  
  274.   StairCount -= Dv.Count();  // forget these stairs
  275. }
  276.  
  277. void SearchInStrip(slist<SEGMENT>& L, slist<SEGMENT>& R, int b, int e)
  278. {
  279.   STAIRCASE     D(b, e);
  280.   
  281.   Split(L, D);
  282.   
  283.   if( !L.empty() ) 
  284.     if( D.empty() ) { L.clear(); return; }
  285.     else SearchInStrip(L, R, b, e);
  286.   
  287.   Merge(D, R);
  288.  
  289.   StairCount -= D.Count();
  290. }
  291.  
  292. void InnerIntPoints(const STAIRCASE& D, const slist<SEGMENT>& I)
  293. { // Precondition : all segments in I are inner from <D.b, D.e>
  294.   // Find intersections between segments in D.Q and I (I is not ordered)
  295.   // O( |I|log|D.Q| + |Int(D,I)|)-time
  296.   
  297.   if ( D.empty() || I.empty() ) return;
  298.   
  299.   SEGMENT s(I.head());
  300.   
  301.   forall(s, I)
  302.     {
  303.       int i = Loc(D, s);   int j = i + 1;
  304.       
  305.       // intersections between s and stairs below s.start()
  306.       while( i > 0 && ( orientation(D[i], s.end()) <= 0 || 
  307.                   !orientation(D[i], s.start()) ) )
  308.           ReportIntersection(D[i--], s);
  309.       
  310.       // intersections between s and stairs above s.start()
  311.       while ( j<=D.Count() && ( orientation(D[j], s.end()) >= 0 ||
  312.                   !orientation(D[j], s.start()) ) )
  313.     ReportIntersection(D[j++], s);
  314.     }
  315. }
  316.  
  317. int  Loc(const STAIRCASE& D, const SEGMENT s)
  318. { // Precondition : s is inner from <D.b, D.e>
  319.   //                D is not empty
  320.    
  321.   if( orientation (D[1], s.start()) < 0 )          return 0;
  322.   if( orientation (D[D.Count()], s.start()) > 0 )  return D.Count();
  323.   
  324.   int locmin = 1;
  325.   int locmax = D.Count(); 
  326.   
  327.   while(locmax-locmin > 1)
  328.     {
  329.       int i = (locmin+locmax)/2;
  330.       SEGMENT sd(D[i]);
  331.       int o = orientation( sd, s.start() ); 
  332.       
  333.       if( o > 0 ) locmin = i;
  334.       else
  335.     if( o < 0 ) locmax = i;
  336.     else { locmin = i; locmax = i; }
  337.     }
  338.   
  339.   return locmin;
  340. }
  341.  
  342. void ReportIntersection(const SEGMENT& s1, const SEGMENT& s2)
  343. {
  344.   if ( report_func ) report_func( (*original)[s1],(*original)[s2] );
  345. }
  346.  
  347. void Split(slist<SEGMENT>& L, STAIRCASE& D)
  348. { 
  349.   slist_item lit = L.first();                 if(!lit) return;
  350.   slist_item lits_prev = nil;
  351.   SEGMENT    s(L.head());
  352.   POINT      pb(init_points[D.b]);
  353.   POINT      pe(init_points[D.e]);
  354.   slist<int> locs;
  355.   int        loc=0, i;
  356.  
  357.   while ( lit )
  358.     {
  359.       s = L.inf( lit );
  360.       if ( compare(s.end(), pe) >= 0 ) // s is spanning $(|pb|,|pe|)$ 
  361.     if ( D.empty() || cmp_at_x(D[D.Count()], s, pe) <= 0 )
  362.       {  
  363.         D.append(s); 
  364.         lit = L.succ(lit);
  365.         if ( lits_prev )
  366.           L.del_succ_item(lits_prev);
  367.         else L.pop();
  368.         loc++; 
  369.       }
  370.     else { locs.append(-loc); lits_prev = lit; lit = L.succ(lit); }
  371.       else   { locs.append(loc); lits_prev = lit; lit = L.succ(lit); }
  372.     }  
  373.  
  374.   if ( D.empty() || L.empty() ) return;
  375.  
  376.   lit = locs.first();
  377.   forall(s, L)
  378.     {
  379.       loc = locs.inf(lit); lit = locs.succ(lit);
  380.       if( loc < 0 )
  381.     {
  382.       i = -loc;
  383.       // intersections between s and stairs below s
  384.       ReportIntersection(D[i--], s); // s would be a stair if not      
  385.       while( i > 0 &&  cmp_at_x(D[i], s, pe) > 0  )
  386.         ReportIntersection(D[i--], s);
  387.     }
  388.       else
  389.     {
  390.       i = loc;
  391.       // intersections between s and stairs below s
  392.       while( i > 0 && orientation(D[i], s.end()) <= 0 )
  393.         ReportIntersection(D[i--], s);
  394.       i = loc+1;
  395.       // intersections between s and stairs above s
  396.       while ( i<=D.Count() && orientation(D[i], s.end()) >= 0 )
  397.         ReportIntersection(D[i++], s);
  398.     }
  399.     }
  400.   locs.clear();
  401. }
  402.  
  403. void Lrs_from_Rls(slist<SEGMENT>& R, slist<SEGMENT>& L, int c)
  404. {    
  405.   SEGMENT         s; 
  406.   POINT           pc(init_points[c]);
  407.  
  408.   // put all segments from R to L
  409.   L.conc(R);
  410.   slist_item      lits_prev = nil;
  411.   slist_item      lit = L.first();
  412.   
  413.   while ( lit && orientation(R.inf(lit), pc) > 0  ) 
  414.     { lits_prev = lit; lit = L.succ(lit); }
  415.   
  416.   slist<SEGMENT>   endings;
  417.   slist<SEGMENT>   pass_start;
  418.   slist<SEGMENT>   nuller;
  419.   SEGMENT          si;
  420.   
  421.   while ( lit && !orientation(R.inf(lit), pc) ) 
  422.     { 
  423.       s = R.inf(lit);
  424.       if( identical( s.end(), pc ) ) // segments ending at pc
  425.     {
  426.       forall(si, endings)    ReportIntersection(s, si);
  427.       forall(si, pass_start) ReportIntersection(s, si);
  428.       endings.append(s);
  429.       lit = L.succ(lit);
  430.       if ( lits_prev )  L.del_succ_item(lits_prev);
  431.       else L.pop();
  432.     }
  433.       else                          // segments passing pc
  434.     {  
  435.       forall(si, endings)    ReportIntersection(s, si);
  436.       forall(si, pass_start)
  437.         if( orientation(s, si.start()) ) ReportIntersection(s, si);
  438.       pass_start.push(s); 
  439.       lit = L.succ(lit);
  440.       if ( lits_prev )  L.del_succ_item(lits_prev);
  441.       else L.pop(); 
  442.     }
  443.     }
  444.  
  445.   while(!seg_queue.empty() && identical(seg_queue.head().start(), pc))
  446.     { 
  447.       s = seg_queue.pop();
  448.       if ( identical( s.end(), pc ) ) // zero length segments
  449.     {
  450.       forall(si, endings)    ReportIntersection(s, si);
  451.       forall(si, pass_start) ReportIntersection(s, si);
  452.       forall(si, nuller)     ReportIntersection(s, si);
  453.       nuller.append(s);
  454.     }
  455.       else                            // segments starting at pc
  456.     { 
  457.       forall(si, endings) ReportIntersection(s, si);
  458.       forall(si, nuller)  ReportIntersection(s, si);
  459.       slist_item lit1 = pass_start.first();
  460.       slist_item lit2 = 0;
  461.       while(lit1)
  462.         {
  463.           si = pass_start.inf(lit1);
  464.           int o = orientation(si, s.end());
  465.           if(o) ReportIntersection(s, si);
  466.           if( o > 0) lit2 = lit1; 
  467.           lit1 = pass_start.succ(lit1);
  468.         }
  469.       if (!lit2) pass_start.push(s);
  470.       else       pass_start.insert(s, lit2);
  471.     }
  472.     }
  473.  
  474.   forall(s, pass_start) 
  475.     {
  476.     if(lits_prev) 
  477.       { 
  478.         L.insert(s, lits_prev); 
  479.         lits_prev = L.succ(lits_prev); 
  480.       }
  481.     else  lits_prev = L.push(s);
  482.     }
  483. }
  484.  
  485. void Merge(STAIRCASE& D, slist<SEGMENT>& R)
  486. { // Merge the stairs in D and the segments in R and return result in R.
  487.   // Precondition : D.Q, and R are orderd by <(D.e.xcoord)
  488.   
  489.   int        loc = 1, i, j;
  490.   POINT      pb(init_points[D.b]);
  491.   POINT      pe(init_points[D.e]);
  492.   slist_item lit = R.first();
  493.   slist_item lits_prev = nil;
  494.   
  495.   if( D.Count() == 0 )  return;
  496.  
  497.   while( lit || loc <= D.Count() )
  498.     {
  499.       if( !lit )      
  500.     { 
  501.       for(; loc <= D.Count();) R.append(D[loc++]);
  502.       return; 
  503.     }
  504.       
  505.       SEGMENT s( R.inf(lit) );
  506.       
  507.       int c_x = loc > D.Count() ? 1 : cmp_at_x(s, D[loc], pe);
  508.       
  509.       if ( c_x < 0 )
  510.     {
  511.       if( compare( s.start(), pb ) > 0 )
  512.         {
  513.           i = loc-1; j = loc;
  514.           
  515.           while( i > 0 && orientation(D[i], s.start()) <= 0 )
  516.         ReportIntersection(D[i--], s);
  517.           
  518.           while( j<=D.Count() && orientation(D[j], s.start()) >= 0)
  519.         ReportIntersection(D[j++], s);
  520.         }
  521.       lits_prev = lit; lit = R.succ(lit);
  522.     }
  523.       else
  524.     if ( !c_x ) 
  525.       {
  526.         do
  527.           {
  528.         lits_prev = lit; lit = R.succ(lit);
  529.           } while( lit && !cmp_at_x(R.inf(lit), D[loc], pe) );
  530.         
  531.         while( loc<D.Count() && !cmp_at_x(D[loc+1],D[loc],pe))
  532.           { R.insert(D[loc++], lits_prev); lits_prev = R.succ(lits_prev); }
  533.       }
  534.     else 
  535.       {
  536.         if(loc > D.Count())
  537.           {
  538.         if( compare( s.start(), pb ) > 0 )
  539.           {
  540.             i = loc-1;
  541.             
  542.             while( i > 0 && orientation(D[i], s.start()) <= 0 )
  543.               ReportIntersection(D[i--], s);
  544.           }
  545.         lits_prev = lit; lit = R.succ(lit);
  546.           }
  547.         else
  548.           if(lits_prev)
  549.         { 
  550.           R.insert(D[loc++], lits_prev); lits_prev = R.succ(lits_prev);
  551.         }
  552.           else  lits_prev = R.push(D[loc++]);
  553.       }
  554.     }
  555. }
  556.  
  557.