Add Trapezoidalization

codebook/8_Geometry/Trapezoidalization.cpp

@@ -0,0 +1,99 @@
+template<class T>
+struct SweepLine {
+  struct cmp {
+    cmp(const SweepLine &_swp): swp(_swp) {}
+    bool operator()(int a, int b) const {
+      if (abs(swp.get_y(a) - swp.get_y(b)) <= swp.eps)
+        return swp.slope_cmp(a, b);
+      return swp.get_y(a) + swp.eps < swp.get_y(b); 
+    }
+    const SweepLine &swp;
+  } _cmp;
+  T curTime, eps, curQ;
+  vector<Line> base;
+  multiset<int, cmp> sweep;
+  multiset<pair<T, int>> event;
+  vector<typename multiset<int, cmp>::iterator> its;
+  vector<typename multiset<pair<T, int>>::iterator> eits;
+  bool slope_cmp(int a, int b) const {
+    assert(a != -1);
+    if (b == -1) return 0;
+    return sign(cross(base[a].Y - base[a].X, base[b].Y - base[b].X)) < 0;
+  }
+  T get_y(int idx) const {
+    if (idx == -1) return curQ;
+    Line l = base[idx];
+    if (l.X.X == l.Y.X) return l.Y.Y;
+    return ((curTime - l.X.X) * l.Y.Y + (l.Y.X - curTime) * l.X.Y) / (l.Y.X - l.X.X);
+  }
+  void insert(int idx) {
+    its[idx] = sweep.insert(idx);
+    if (its[idx] != sweep.begin()) 
+      update_event(*prev(its[idx]));
+    update_event(idx);
+    event.emplace(base[idx].Y.X, idx + 2 * SZ(base));
+  }
+  void erase(int idx) {
+    assert(eits[idx] == event.end());
+    auto p = sweep.erase(its[idx]);
+    its[idx] = sweep.end();
+    if (p != sweep.begin())
+      update_event(*prev(p));
+  }
+  void update_event(int idx) {
+    if (eits[idx] != event.end())  
+      event.erase(eits[idx]);
+    eits[idx] = event.end();
+    auto nxt = next(its[idx]);
+    if (nxt == sweep.end() || !slope_cmp(idx, *nxt)) return;
+    auto t = intersect(base[idx].X, base[idx].Y, base[*nxt].X, base[*nxt].Y).X;
+    if (t + eps < curTime || t >= min(base[idx].Y.X, base[*nxt].Y.X)) return;
+    eits[idx] = event.emplace(t, idx + SZ(base));
+  }
+  void swp(int idx) {
+    assert(eits[idx] != event.end());
+    eits[idx] = event.end();
+    int nxt = *next(its[idx]);
+    swap((int&)*its[idx], (int&)*its[nxt]);
+    swap(its[idx], its[nxt]);
+    if (its[nxt] != sweep.begin())
+      update_event(*prev(its[nxt]));
+    update_event(idx);
+  }
+  // only expected to call the functions below
+  SweepLine(T t, T e, vector<Line> vec): _cmp(*this), curTime(t), eps(e), curQ(), base(vec), sweep(_cmp), event(), its(SZ(vec), sweep.end()), eits(SZ(vec), event.end()) {
+    for (int i = 0; i < SZ(base); ++i) {
+      auto &[p, q] = base[i];
+      if (p > q) swap(p, q);
+      if (p.X <= curTime && curTime <= q.X)
+        insert(i);
+      else if (curTime < p.X)
+        event.emplace(p.X, i);
+    }
+  }
+  void setTime(T t, bool ers = false) {
+    assert(t >= curTime);
+    while (!event.empty() && event.begin()->X <= t) {
+      auto [et, idx] = *event.begin();
+      int s = idx / SZ(base);
+      idx %= SZ(base);
+      if (abs(et - t) <= eps && s == 2 && !ers) break;
+      curTime = et;
+      event.erase(event.begin());
+      if (s == 2) erase(idx);
+      else if (s == 1) swp(idx);
+      else insert(idx); 
+    }
+    curTime = t;
+  }
+  T nextEvent() {
+    if (event.empty()) return INF;
+    return event.begin()->X;
+  }
+  int lower_bound(T y) {
+    curQ = y;
+    auto p = sweep.lower_bound(-1); 
+    if (p == sweep.end()) return -1;
+    return *p;
+  }
+};

codebook/content.tex

@@ -211,6 +211,8 @@ \subsection{PolyUnion*} % test by CF 101673 A
 \lstinputlisting{8_Geometry/PolyUnion.cpp}
 \subsection{PolyCut}
 \lstinputlisting{8_Geometry/PolyCut.cpp}
+\subsection{Trapezoidalization} % test by aizu CGL_3_C, TIOJ 1403
+\lstinputlisting{8_Geometry/Trapezoidalization.cpp}
 \subsection{Polar Angle Sort*} % test by NTUJ 2270
 \lstinputlisting{8_Geometry/Polar_Angle_Sort.cpp}
 \subsection{Half plane intersection*} % test by qoj 2162 (World Finals' problem)