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 / _polygon.c < prev    next >
C/C++ Source or Header  |  1996-09-03  |  10KB  |  453 lines
  1. /*******************************************************************************
  2. +
  3. +  LEDA 3.4
  4. +
  5. +  _polygon.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/polygon.h>
  19. #include <LEDA/graph.h>
  20. #include <math.h>
  21. #include <assert.h>
  22.  
  23. extern void SWEEP_SEGMENTS(const list<segment>&,GRAPH<point,segment>&);
  24.  
  25. //------------------------------------------------------------------------------
  26. // polygons
  27. //
  28. // by S. Naeher (1995)
  29. //------------------------------------------------------------------------------
  30.  
  31.  
  32. polygon_rep::polygon_rep(const list<segment>& sl) :  seg_list(sl)
  33. { segment s;
  34.   forall(s,sl) pt_list.append(s.start());
  35.  }
  36.  
  37.  
  38.  
  39. double polygon::compute_area(const list<segment>& seg_list) const
  40. {
  41.  
  42.   if (seg_list.length() < 3) return 0;
  43.  
  44.   list_item it = seg_list.item(1);
  45.   point     p  = seg_list[it].source();
  46.  
  47.   it = seg_list.succ(it);
  48.  
  49.   double    A  = 0;
  50.  
  51.   while (it)
  52.   { segment s = seg_list[it];
  53.     A += ::area(p,s.source(),s.target());
  54.     it = seg_list.succ(it);
  55.    }
  56.  
  57.   return A;
  58. }
  59.  
  60.  
  61. static void check_simplicity(const list<segment>& seg_list)
  62. { GRAPH<point,segment> G;
  63.   SWEEP_SEGMENTS(seg_list,G);
  64.   node v;
  65.   forall_nodes(v,G)
  66.     if (G.degree(v) != 2)
  67.       error_handler(1,"polygon: polygon must be simple");
  68. }
  69.  
  70.  
  71. polygon::polygon(const list<segment>& sl)
  72. { PTR = new polygon_rep(sl); }
  73.  
  74.  
  75. polygon::polygon(const list<point>& pl)
  76. { 
  77.   list<segment> seglist;
  78.  
  79.   for(list_item it = pl.first(); it; it = pl.succ(it))
  80.      seglist.append(segment(pl[it],pl[pl.cyclic_succ(it)]));
  81.  
  82.   if (compute_area(seglist) < 0)
  83.   { // reverse edge list
  84.     seglist.clear();
  85.     for(list_item it = pl.last(); it; it = pl.pred(it))
  86.       seglist.append(segment(pl[it],pl[pl.cyclic_pred(it)]));
  87.    }
  88.   check_simplicity(seglist);
  89.  
  90.   PTR = new polygon_rep(seglist);
  91. }
  92.  
  93.  
  94. polygon polygon::translate(double alpha, double d) const
  95. { list<segment> L;
  96.   segment s;
  97.   forall(s,ptr()->seg_list) L.append(s.translate(alpha,d));
  98.   return polygon(L);
  99. }
  100.  
  101.  
  102. polygon polygon::translate(const vector& v) const
  103. { list<segment> L;
  104.   segment s;
  105.   forall(s,ptr()->seg_list) L.append(s.translate(v));
  106.   return polygon(L);
  107. }
  108.  
  109.  
  110. polygon polygon::rotate(const point& p, double alpha) const
  111. { list<segment> L;
  112.   segment s;
  113.   forall(s,ptr()->seg_list) L.append(s.rotate(p,alpha));
  114.   return polygon(L);
  115. }
  116.  
  117. polygon polygon::rotate(double alpha) const
  118. { return rotate(point(0,0),alpha); }
  119.  
  120.  
  121. polygon polygon::rotate90(const point& p) const
  122. { list<segment> L;
  123.   segment s;
  124.   forall(s,ptr()->seg_list) L.append(s.rotate90(p));
  125.   return polygon(L);
  126. }
  127.  
  128.  
  129. polygon polygon::reflect(const point& p, const point& q) const
  130. { list<segment> L;
  131.   segment s;
  132.   forall(s,ptr()->seg_list) L.append(s.reflect(p,q));
  133.   return polygon(L);
  134. }
  135.  
  136.  
  137. polygon polygon::reflect(const point& p) const
  138. { list<segment> L;
  139.   segment s;
  140.   forall(s,ptr()->seg_list) L.append(s.reflect(p));
  141.   return polygon(L);
  142. }
  143.  
  144.  
  145. bool polygon::inside(const point& p) const
  146. {
  147.   list<segment>& seglist = ptr()->seg_list;
  148.  
  149.   int count = 0;
  150.  
  151.   double px = p.xcoord();
  152.  
  153.   list_item it0 = seglist.first();
  154.   list_item it1 = seglist.first();
  155.  
  156.   double x0 = seglist[it0].xcoord2();
  157.   double x1 = seglist[it1].xcoord2();
  158.  
  159.   list_item it;
  160.  
  161.   forall_items(it,seglist)
  162.   { segment s = seglist[it];
  163.     if (s.xcoord2() < x0) 
  164.     { it0 = it;
  165.       x0 = s.xcoord2();
  166.      }
  167.     if (s.xcoord2() > x1) 
  168.     { it1 = it1;
  169.       x1 = s.xcoord2();
  170.      }
  171.    }
  172.  
  173.   if (px <= x0 || x1 <= px) return false;
  174.  
  175.   while (seglist[it0].xcoord2() <= px) it0 = seglist.cyclic_succ(it0);
  176.  
  177.   it = it0;
  178.   do
  179.   { while (seglist[it].xcoord2() >= px) it = seglist.cyclic_succ(it);
  180.     if (orientation(seglist[it],p) > 0) count++;
  181.     while (seglist[it].xcoord2() <= px) it = seglist.cyclic_succ(it);
  182.     if (orientation(seglist[it],p) < 0) count++;
  183.   } while (it != it0);
  184.  
  185.   return count % 2;
  186.  
  187. }
  188.  
  189.  
  190.  
  191. bool polygon::outside(const point& p) const { return !inside(p); }
  192.  
  193.  
  194. list<point> polygon::intersection(const segment& s) const
  195. { list<point> result;
  196.   segment t;
  197.   point p;
  198.  
  199.   forall(t,ptr()->seg_list) 
  200.     if (s.intersection(t,p))
  201.      if (result.empty()) result.append(p);
  202.      else if (p != result.tail() ) result.append(p);
  203.  
  204.   return result;
  205. }
  206.  
  207.  
  208. list<point> polygon::intersection(const line& l) const
  209. { list<point> result;
  210.   segment t;
  211.   point p;
  212.  
  213.   forall(t,ptr()->seg_list) 
  214.     if (l.intersection(t,p))
  215.      if (result.empty()) result.append(p);
  216.      else if (p != result.tail() ) result.append(p);
  217.  
  218.   return result;
  219. }
  220.  
  221.  
  222. // intersection or union with polygon
  223.  
  224. static bool test_edge(GRAPH<point,segment>& G, edge i1, int mode)
  225. { 
  226.   node v  = G.target(i1);
  227.   edge i2 = G.cyclic_in_succ(i1);
  228.  
  229.   node u1 = G.source(i1);
  230.   node u2 = G.source(i2);
  231.  
  232.   if (i1 == i2) return false;
  233.  
  234.   // parallel edges
  235.   if (u1 == u2) return i1 == G.first_in_edge(v);
  236.  
  237.   edge o1 = G.first_adj_edge(v);
  238.   edge o2 = G.last_adj_edge(v);
  239.  
  240.   // anti-parallel edges
  241.   if (target(o1) == u1 || target(o2) == u1) return false;
  242.  
  243.  
  244.   point p1 = G[o1].target();
  245.   point p2 = G[o2].target();
  246.  
  247.   segment si1 = G[i1];
  248.   segment si2 = G[i2];
  249.  
  250.   if (mode == 0) // intersection
  251.   { if (orientation(si1,si2.source()) > 0)
  252.       return orientation(si1,p1) > 0 && orientation(si2,p1) < 0 && 
  253.              orientation(si1,p2) > 0 && orientation(si2,p2) < 0;
  254.    else
  255.       return (orientation(si1,p1) > 0 || orientation(si2,p1) < 0) && 
  256.              (orientation(si1,p2) > 0 || orientation(si2,p2) < 0);
  257.   }
  258.   else // union
  259.   { if (orientation(si1,si2.source()) < 0)
  260.       return orientation(si1,p1) <= 0 && orientation(si2,p1) >= 0 ||
  261.              orientation(si1,p2) <= 0 && orientation(si2,p2) >= 0;
  262.    else
  263.       return orientation(si1,p1) <= 0 || orientation(si2,p1) >= 0 || 
  264.              orientation(si1,p2) <= 0 || orientation(si2,p2) >= 0;
  265.    }
  266. }
  267.  
  268.  
  269.  
  270. static edge next_edge(GRAPH<point,segment>& G, edge i1, int dir)
  271. { 
  272.   // if dir = +1 turn left
  273.   // if dir = -1 turn right
  274.  
  275.   node v = G.target(i1);
  276.  
  277.   edge o1 = G.first_adj_edge(v);
  278.   edge o2 = G.last_adj_edge(v);
  279.  
  280.   if (o1 == o2) return o1;
  281.  
  282.   node w1 = G.target(o1);
  283.   node w2 = G.target(o2);
  284.  
  285.   if (w1 == w2) return (o1 == G.first_in_edge(w1)) ? o1 : o2;
  286.  
  287.   segment si1 = G[i1];
  288.   segment so1 = G[o1];
  289.   segment so2 = G[o2];
  290.  
  291.   int orient1 = orientation(si1,so1.target());
  292.   int orient2 = orientation(si1,so2.target());
  293.  
  294.   if (orient1 == orient2)
  295.      return (orientation(so1,so2.target()) == dir) ? o2 : o1;
  296.   else
  297.      if (orient1 < orient2)
  298.         return (dir == -1) ? o1 : o2;
  299.      else
  300.         return (dir == -1) ? o2 : o1;
  301. }
  302.  
  303.  
  304.  
  305. list_polygon_ polygon::intersection(const polygon& P) const
  306. {
  307.   list<polygon> result;
  308.   list<segment> seg_list;
  309.   GRAPH<point,segment> G;
  310.  
  311.   segment s;
  312.   forall(s,ptr()->seg_list) seg_list.append(s);
  313.   forall(s,P.ptr()->seg_list) seg_list.append(s);
  314.  
  315.   SWEEP_SEGMENTS(seg_list,G);
  316.  
  317.   bool borders_intersect = false;
  318.  
  319.   node v;
  320.   forall_nodes(v,G) 
  321.   { int deg = G.degree(v);
  322.     assert(deg == 2 || deg == 4);
  323.     if (deg == 4) borders_intersect = true;
  324.    }
  325.     
  326.  
  327.   if ( ! borders_intersect )
  328.   { // no intersections between edges of (*this) and P
  329.  
  330.     // check for inclusion
  331.  
  332.     segment s1 = ptr()->seg_list.head();
  333.     segment s2 = P.ptr()->seg_list.head();
  334.  
  335.     if ( P.inside(s1.start()) ) result.append(*this);
  336.     if ( inside(s2.start()) ) result.append(P);                   
  337.  
  338.     return result;
  339.  
  340.    }
  341.  
  342.   edge_array<bool> marked(G,false);
  343.  
  344.   edge e;
  345.   forall_edges(e,G) 
  346.   { node start = source(e);
  347.     if ( !marked[e] && test_edge(G,e,0) )
  348.     { list<segment> pol;
  349.       do { node v = source(e);
  350.            node w = target(e);
  351.            marked[e] = true;
  352.            pol.append(segment(G[v],G[w]));
  353.            e = next_edge(G,e,+1);
  354.           } while (source(e) != start);
  355.       result.append(polygon(pol));
  356.      }
  357.    }
  358.  
  359.  return result;
  360.  
  361. }
  362.  
  363.  
  364. list_polygon_ polygon::unite(const polygon& P) const
  365. {
  366.   list<polygon> result;
  367.   list<segment> seg_list;
  368.   GRAPH<point,segment> G;
  369.  
  370.   segment s;
  371.   forall(s,ptr()->seg_list) seg_list.append(s);
  372.   forall(s,P.ptr()->seg_list) seg_list.append(s);
  373.  
  374.   SWEEP_SEGMENTS(seg_list,G);
  375.  
  376.   bool borders_intersect = false;
  377.  
  378.   node v;
  379.   forall_nodes(v,G) 
  380.   { int deg = G.degree(v);
  381.     assert(deg == 2 || deg == 4);
  382.     if (deg == 4) borders_intersect = true;
  383.    }
  384.  
  385.   if (!borders_intersect)
  386.   { 
  387.     // check for inclusion
  388.  
  389.     segment s1 = ptr()->seg_list.head();
  390.     segment s2 = P.ptr()->seg_list.head();
  391.  
  392.     if ( ! P.inside(s1.start())) result.append(*this);
  393.     if ( ! inside(s2.start())) result.append(P);                   
  394.  
  395.     return result;
  396.    }
  397.  
  398.   edge_array<bool> marked(G,false);
  399.  
  400.   edge e;
  401.   forall_edges(e,G) 
  402.   { node start = source(e);
  403.     if ( !marked[e] && test_edge(G,e,1) )
  404.     { list<segment> pol;
  405.       do { node v = source(e);
  406.            node w = target(e);
  407.            marked[e] = true;
  408.            pol.append(segment(G[v],G[w]));
  409.            e = next_edge(G,e,-1);
  410.           } while (source(e) != start);
  411.       result.append(polygon(pol));
  412.      }
  413.    }
  414.  
  415.   if (result.length() > 1)
  416.   { // find outer polygon and move it to front of result
  417.  
  418.     list_item it;
  419.     segment s;
  420.     list_item min_it = nil;
  421.     double    min_x = MAXDOUBLE;
  422.  
  423.     forall_items(it,result)
  424.       forall(s,result[it].ptr()->seg_list) 
  425.         if (s.xcoord1() < min_x)
  426.         { min_x = s.xcoord1(); 
  427.           min_it = it;
  428.          }
  429.   
  430.     result.move_to_front(min_it);
  431.   }
  432.  
  433.  return result;
  434. }
  435.  
  436.  
  437.  
  438.  
  439.  
  440.  
  441. ostream& operator<<(ostream& out, const polygon& p) 
  442. { p.vertices().print(out);
  443.   return out << endl;
  444.  } 
  445.  
  446. istream& operator>>(istream& in,  polygon& p) 
  447. { list<point> L; 
  448.   L.read(in,'\n'); 
  449.   p = polygon(L); 
  450.   return in;
  451. }
  452.  
  453.