diff --git a/src/db/db/db.pro b/src/db/db/db.pro index 499046cc7..bbae9671a 100644 --- a/src/db/db/db.pro +++ b/src/db/db/db.pro @@ -81,6 +81,7 @@ SOURCES = \ dbPolygonGenerators.cc \ dbPropertiesFilter.cc \ dbPropertiesRepository.cc \ + dbQuadTree.cc \ dbReader.cc \ dbRecursiveInstanceIterator.cc \ dbRecursiveShapeIterator.cc \ @@ -319,6 +320,7 @@ HEADERS = \ dbPropertiesFilter.h \ dbPropertiesRepository.h \ dbPropertyConstraint.h \ + dbQuadTree.h \ dbReader.h \ dbRecursiveInstanceIterator.h \ dbRecursiveShapeIterator.h \ diff --git a/src/db/db/dbQuadTree.cc b/src/db/db/dbQuadTree.cc new file mode 100644 index 000000000..7e9992fbc --- /dev/null +++ b/src/db/db/dbQuadTree.cc @@ -0,0 +1,24 @@ + +/* + + KLayout Layout Viewer + Copyright (C) 2006-2025 Matthias Koefferlein + + This program is free software; you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program; if not, write to the Free Software + Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +*/ + +#include "dbQuadTree.h" + diff --git a/src/db/db/dbQuadTree.h b/src/db/db/dbQuadTree.h new file mode 100644 index 000000000..e3ab61ef2 --- /dev/null +++ b/src/db/db/dbQuadTree.h @@ -0,0 +1,494 @@ + +/* + + KLayout Layout Viewer + Copyright (C) 2006-2025 Matthias Koefferlein + + This program is free software; you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program; if not, write to the Free Software + Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +*/ + +#ifndef HDR_dbQuadTree +#define HDR_dbQuadTree + +#include "dbBox.h" +#include + +namespace db +{ + +template +class quad_tree_node +{ +public: + typedef typename T::coord_type coord_type; + typedef db::box box_type; + typedef db::point point_type; + typedef db::vector vector_type; + typedef std::vector objects_vector; + + quad_tree_node (const point_type ¢er) + : m_split (false), m_center (center) + { + for (unsigned int i = 0; i < 4; ++i) { + m_q [i] = 0; + } + } + + ~quad_tree_node () + { + for (unsigned int i = 0; i < 4; ++i) { + delete m_q [i]; + m_q [i] = 0; + } + } + + const point_type ¢er () const + { + return m_center; + } + + void insert_top (const T &value, const box_type &total_box) + { + insert (value, propose_ucenter (total_box)); + } + + bool remove (const T &value) + { + box_type b = BC () (value); + + if (! m_split || b.contains (m_center)) { + for (auto i = m_objects.begin (); i != m_objects.end (); ++i) { + if (*i == value) { + m_objects.erase (i); + return true; + } + } + } + + for (unsigned int i = 0; i < 4; ++i) { + if (m_q[i] && b.inside (m_q[i]->box (m_center))) { + return b.remove (value); + } + } + + return false; + } + + const objects_vector &objects () const + { + return m_objects; + } + + box_type q_box (unsigned int n) const + { + if (m_q[n]) { + return m_q[n]->box (m_center); + } else { + return box_type (); + } + } + + quad_tree_node *node (unsigned int n) const + { + return m_q [n]; + } + + bool empty () const + { + if (m_objects.empty ()) { + for (unsigned int n = 0; n < 4; ++n) { + if (m_q[n] && ! m_q[n]->empty ()) { + return false; + } + } + return true; + } else { + return false; + } + } + + size_t size () const + { + size_t count = m_objects.size (); + for (unsigned int n = 0; n < 4; ++n) { + if (m_q[n]) { + count += m_q[n]->size (); + } + } + return count; + } + +private: + bool m_split; + point_type m_center; + quad_tree_node *m_q [4]; + objects_vector m_objects; + + box_type box (const point_type &ucenter) const + { + return box_type (ucenter, ucenter - (ucenter - m_center) * 2.0); + } + + box_type q (unsigned int n, const point_type &ucenter) const + { + // NOTE: with this definition the opposite quad index is 3 - n + + vector_type vx (std::abs (ucenter.x () - m_center.x ()), 0); + vector_type vy (0, std::abs (ucenter.y () - m_center.y ())); + switch (n) { + case 0: + return box_type (m_center - vx - vy, m_center); + case 1: + return box_type (m_center - vy, m_center + vx); + case 2: + return box_type (m_center - vx, m_center + vy); + default: + return box_type (m_center, m_center + vx + vy); + } + } + + void split (const point_type &ucenter) + { + m_split = true; + + objects_vector ov; + ov.swap (m_objects); + + for (auto o = ov.begin (); o != ov.end (); ++o) { + insert (*o, ucenter); + } + } + + void insert (const T &value, const point_type &ucenter) + { + if (! m_split && m_objects.size () + 1 < thr) { + + m_objects.push_back (value); + + } else { + + if (! m_split) { + split (ucenter); + } + + box_type b = BC () (value); + // @@@ should exclude m_center on box + if (b.contains (m_center)) { + m_objects.push_back (value); + return; + } + + for (unsigned int i = 0; i < 4; ++i) { + box_type bq = q (i, ucenter); + if (b.inside (bq)) { + if (! m_q[i]) { + m_q[i] = new quad_tree_node (bq.center ()); + } + m_q[i]->insert (value, m_center); + return; + } + } + + for (unsigned int i = 0; i < 4; ++i) { + if (m_q[i]) { + grow (m_center - (m_center - m_q[i]->center ()) * 2.0); + insert (value, ucenter); + return; + } + } + + tl_assert (false); + + } + } + + void grow (const point_type &ucenter) + { + for (unsigned int i = 0; i < 4; ++i) { + if (m_q[i]) { + quad_tree_node *n = m_q[i]; + m_q[i] = new quad_tree_node (q (i, ucenter).center ()); + m_q[i]->m_q[3 - i] = n; + } + } + } + + point_type propose_ucenter (const box_type &total_box) const + { + for (unsigned int i = 0; i < 4; ++i) { + if (m_q[i]) { + return m_center - (m_center - m_q[i]->center ()) * 2.0; + } + } + + coord_type dx = std::max (std::abs (total_box.left () - m_center.x ()), std::abs (total_box.right () - m_center.y ())); + coord_type dy = std::max (std::abs (total_box.bottom () - m_center.y ()), std::abs (total_box.top () - m_center.y ())); + return m_center - vector_type (dx, dy); + } +}; + +template +class quad_tree_iterator +{ +public: + typedef quad_tree_node quad_tree_node_type; + typedef typename T::coord_type coord_type; + typedef db::box box_type; + + quad_tree_iterator () + : m_s (), m_i (0) + { + // .. nothing yet .. + } + + quad_tree_iterator (const quad_tree_node_type *root, const S &s) + : m_s (s), m_i (0) + { + m_stack.push_back (std::make_pair (root, -1)); + validate (); + } + + bool at_end () const + { + return m_stack.empty (); + } + + quad_tree_iterator &operator++ () + { + ++m_i; + validate (); + return *this; + } + + const T &operator* () const + { + return m_stack.back ().first->objects () [m_i]; + } + + const T *operator-> () const + { + return (m_stack.back ().first->objects ().begin () + m_i).operator-> (); + } + +public: + S m_s; + std::vector > m_stack; + size_t m_i; + + void validate () + { + auto s = m_stack.end (); + tl_assert (s != m_stack.begin ()); + --s; + + const quad_tree_node_type *qn = s->first; + while (m_i < s->first->objects ().size () && ! m_s.select (s->first->objects () [m_i])) { + ++m_i; + } + + if (m_i < qn->objects ().size ()) { + return; + } + + m_i = 0; + + for (unsigned int n = 0; n < 4; ++n) { + box_type bq = qn->q_box (n); + if (! bq.empty () && m_s.select_quad (bq)) { + m_stack.back ().second = n; + m_stack.push_back (std::make_pair (qn->node (n), -1)); + validate (); + return; + } + } + + while (! m_stack.empty ()) { + + m_stack.pop_back (); + + int &n = m_stack.back ().second; + while (++n < 4) { + box_type bq = qn->q_box (n); + if (! bq.empty () && m_s.select_quad (bq)) { + m_stack.push_back (std::make_pair (qn->node (n), -1)); + validate (); + return; + } + } + + } + } +}; + +template +class quad_tree_always_sel +{ +public: + typedef typename T::coord_type coord_type; + typedef db::box box_type; + + bool select (const T &) const + { + return true; + } + + bool select_quad (const box_type &) const + { + return true; + } +}; + +template +class quad_tree_touching_sel +{ +public: + typedef typename T::coord_type coord_type; + typedef db::box box_type; + + quad_tree_touching_sel () + { + // .. nothing yet .. + } + + quad_tree_touching_sel (const box_type &box) + : m_box (box) + { + // .. nothing yet .. + } + + bool select (const T &value) const + { + return select_quad (BC () (value)); + } + + bool select_quad (const box_type &box) const + { + return m_box.touches (box); + } + +private: + box_type m_box; +}; + +template +class quad_tree_overlapping_sel +{ +public: + typedef typename T::coord_type coord_type; + typedef db::box box_type; + + quad_tree_overlapping_sel () + { + // .. nothing yet .. + } + + quad_tree_overlapping_sel (const box_type &box) + : m_box (box) + { + // .. nothing yet .. + } + + bool select (const T &value) const + { + return select_quad (BC () (value)); + } + + bool select_quad (const box_type &box) const + { + return m_box.overlaps (box); + } + +private: + box_type m_box; +}; + +template +class quad_tree +{ +public: + typedef quad_tree_node quad_tree_node_type; + typedef quad_tree_iterator > quad_tree_flat_iterator; + typedef quad_tree_iterator > quad_tree_touching_iterator; + typedef quad_tree_iterator > quad_tree_overlapping_iterator; + typedef typename T::coord_type coord_type; + typedef db::box box_type; + typedef db::point point_type; + typedef db::vector vector_type; + typedef std::vector objects_vector; + + quad_tree () + : m_root (point_type ()) + { + // .. nothing yet .. + } + + bool empty () const + { + return m_root.empty (); + } + + size_t size () const + { + return m_root.size (); + } + + void insert (const T &value) + { + box_type b = BC () (value); + if (b.empty ()) { + return; + } + + m_total_box += b; + m_root.insert_top (value, m_total_box); + } + + void erase (const T &value) + { + m_root.remove (value); + } + + quad_tree_flat_iterator begin () const + { + return quad_tree_flat_iterator (&m_root, quad_tree_always_sel ()); + } + + quad_tree_overlapping_iterator begin_overlapping (const box_type &box) const + { + if (m_total_box.overlaps (box)) { + return quad_tree_overlapping_iterator (&m_root, quad_tree_overlapping_sel (box)); + } else { + return quad_tree_overlapping_iterator (); + } + } + + quad_tree_touching_iterator begin_touching (const box_type &box) const + { + if (m_total_box.touches (box)) { + return quad_tree_touching_iterator (&m_root, quad_tree_touching_sel (box)); + } else { + return quad_tree_touching_iterator (); + } + } + +private: + box_type m_total_box; + quad_tree_node_type m_root; +}; + +} + +#endif diff --git a/src/db/db/dbTriangle.cc b/src/db/db/dbTriangle.cc index 1016ad3ec..0fa3018c9 100644 --- a/src/db/db/dbTriangle.cc +++ b/src/db/db/dbTriangle.cc @@ -507,15 +507,23 @@ Triangle::common_edge (const Triangle *other) const int Triangle::contains (const db::DPoint &point) const { + // the relative distance of the point from the edge + double alpha = db::epsilon; + + double d = db::vprod (*mp_v[2] - *mp_v[0], *mp_v[1] - *mp_v[0]); + if (fabs (d) < db::epsilon * db::epsilon) { + return -1; + } + int res = 1; - const Vertex *vl = mp_v[2];; + const Vertex *vl = mp_v[2]; for (int i = 0; i < 3; ++i) { - const Vertex *v = mp_v[i];; - int s = db::DEdge (*vl, *v).side_of (point); - if (s == 0) { - res = 0; - } else if (s > 0) { + const Vertex *v = mp_v[i]; + double n = db::vprod (point - *vl, *v - *vl) / d; + if (n < -alpha) { return -1; + } else if (n < alpha) { + res = 0; } vl = v; } diff --git a/src/db/db/dbTriangles.cc b/src/db/db/dbTriangles.cc index e6aff5ee5..5a0475a5d 100644 --- a/src/db/db/dbTriangles.cc +++ b/src/db/db/dbTriangles.cc @@ -355,6 +355,7 @@ Triangles::insert_point (db::DCoord x, db::DCoord y, std::list > *new_triangles) { + // @@@ printf("@@@ insert %d\n", (int)mp_triangles.size ()); fflush(stdout); std::vector tris = find_triangle_for_point (*vertex); // the new vertex is outside the domain @@ -412,7 +413,24 @@ Triangles::find_closest_edge (const db::DPoint &p, db::Vertex *vstart, bool insi { if (!vstart) { - if (! mp_triangles.empty ()) { + if (m_is_constrained && ! m_initial_segments.empty ()) { + + // Use the closest initial segment for the starting point + // TODO: m_initial_segments should be some search tree + + double d = -1; + for (auto i = m_initial_segments.begin (); i != m_initial_segments.end (); ++i) { + if (db::sprod_sign (p - i->first.p1 (), i->first.d ()) >= 0 && + db::sprod_sign (p - i->first.p2 (), i->first.d ()) < 0) { + double ds = std::abs (i->first.distance (p)); + if (d < 0.0 || ds < d) { + vstart = i->second; + d = ds; + } + } + } + + } else if (! mp_triangles.empty ()) { unsigned int ls = 0; size_t n = m_vertex_heap.size (); @@ -1319,6 +1337,7 @@ void Triangles::constrain (const std::vector > &contours) { tl_assert (! m_is_constrained); + m_initial_segments.clear (); std::vector > > resolved_edges; @@ -1329,8 +1348,10 @@ Triangles::constrain (const std::vector > &contours) if (vv == c->end ()) { vv = c->begin (); } - resolved_edges.push_back (std::make_pair (db::DEdge (**v, **vv), std::vector ())); + db::DEdge e (**v, **vv); + resolved_edges.push_back (std::make_pair (e, std::vector ())); resolved_edges.back ().second = ensure_edge (*v, *vv); + m_initial_segments.push_back (std::make_pair (e, *v)); } } @@ -1606,6 +1627,7 @@ Triangles::refine (const TriangulateParameters ¶meters) auto cr = (*t)->circumcircle(); auto center = cr.first; + // printf("@@@ %s -- %s\n", (*t)->to_string().c_str(), center.to_string ().c_str()); fflush(stdout); if ((*t)->contains (center) >= 0) { if (tl::verbosity () >= parameters.base_verbosity + 20) { diff --git a/src/db/db/dbTriangles.h b/src/db/db/dbTriangles.h index f86f77c1c..f15b2dc7a 100644 --- a/src/db/db/dbTriangles.h +++ b/src/db/db/dbTriangles.h @@ -311,6 +311,7 @@ private: std::vector m_returned_edges; tl::stable_vector m_vertex_heap; bool m_is_constrained; + std::vector > m_initial_segments; size_t m_level; size_t m_id; size_t m_flips, m_hops; diff --git a/src/db/unit_tests/dbQuadTreeTests.cc b/src/db/unit_tests/dbQuadTreeTests.cc new file mode 100644 index 000000000..fd74e164b --- /dev/null +++ b/src/db/unit_tests/dbQuadTreeTests.cc @@ -0,0 +1,77 @@ + +/* + + KLayout Layout Viewer + Copyright (C) 2006-2025 Matthias Koefferlein + + This program is free software; you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation; either version 2 of the License, or + (at your option) any later version. + + This program is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with this program; if not, write to the Free Software + Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +*/ + + +#include "dbQuadTree.h" +#include "dbBoxConvert.h" +#include "tlUnitTest.h" +#include "tlString.h" + +typedef db::quad_tree, size_t (1)> my_quad_tree; + +std::string find_touching (const my_quad_tree &qt, const db::DBox &box) +{ + std::vector v; + auto i = qt.begin_touching (box); + while (! i.at_end ()) { + v.push_back (i->to_string ()); + ++i; + } + std::sort (v.begin (), v.end ()); + return tl::join (v, "/"); +} + +std::string find_overlapping (const my_quad_tree &qt, const db::DBox &box) +{ + std::vector v; + auto i = qt.begin_overlapping (box); + while (! i.at_end ()) { + v.push_back (i->to_string ()); + ++i; + } + std::sort (v.begin (), v.end ()); + return tl::join (v, "/"); +} + +TEST(basic) +{ + my_quad_tree tree; + EXPECT_EQ (tree.empty (), true); + EXPECT_EQ (tree.size (), size_t (0)); + + tree.insert (db::DBox ()); + EXPECT_EQ (tree.empty (), true); + EXPECT_EQ (tree.size (), size_t (0)); + + tree.insert (db::DBox (-1, -2, 3, 4)); + EXPECT_EQ (tree.empty (), false); + EXPECT_EQ (tree.size (), size_t (1)); + + EXPECT_EQ (find_touching (tree, db::DBox (-2, 0, -1, 0)), "..."); + + tree.insert (db::DBox (-1, -3, 3, 0)); + EXPECT_EQ (tree.empty (), false); + EXPECT_EQ (tree.size (), size_t (1)); + + EXPECT_EQ (find_touching (tree, db::DBox (-2, 0, -1, 0)), "..."); +} + diff --git a/src/db/unit_tests/dbTriangleTests.cc b/src/db/unit_tests/dbTriangleTests.cc index c93c93b8d..439eb1c65 100644 --- a/src/db/unit_tests/dbTriangleTests.cc +++ b/src/db/unit_tests/dbTriangleTests.cc @@ -302,6 +302,42 @@ TEST(Triangle_contains) } } +TEST(Triangle_contains_small) +{ + db::Vertex v1; + db::Vertex v2 (0.001, 0.002); + db::Vertex v3 (0.002, 0.001); + + TestableTriangleEdge s1 (&v1, &v2); + TestableTriangleEdge s2 (&v2, &v3); + TestableTriangleEdge s3 (&v3, &v1); + + { + db::Triangle tri (&s1, &s2, &s3); + EXPECT_EQ (tri.contains (db::DPoint (0, 0)), 0); + EXPECT_EQ (tri.contains (db::DPoint (-0.001, -0.002)), -1); + EXPECT_EQ (tri.contains (db::DPoint (0.0005, 0.001)), 0); + EXPECT_EQ (tri.contains (db::DPoint (0.0005, 0.002)), -1); + EXPECT_EQ (tri.contains (db::DPoint (0.0025, 0.001)), -1); + EXPECT_EQ (tri.contains (db::DPoint (0.001, -0.001)), -1); + EXPECT_EQ (tri.contains (db::DPoint (0.001, 0.001)), 1); + } + + s1.reverse (); + s2.reverse (); + s3.reverse (); + + { + db::Triangle tri2 (&s3, &s2, &s1); + EXPECT_EQ (tri2.contains(db::DPoint(0, 0)), 0); + EXPECT_EQ (tri2.contains(db::DPoint(0.0005, 0.001)), 0); + EXPECT_EQ (tri2.contains(db::DPoint(0.0005, 0.002)), -1); + EXPECT_EQ (tri2.contains(db::DPoint(0.0025, 0.001)), -1); + EXPECT_EQ (tri2.contains(db::DPoint(0.001, -0.001)), -1); + EXPECT_EQ (tri2.contains(db::DPoint(0.001, 0.001)), 1); + } +} + TEST(Triangle_circumcircle) { db::Vertex v1; diff --git a/src/db/unit_tests/dbTrianglesTests.cc b/src/db/unit_tests/dbTrianglesTests.cc index 0bc00893e..9b7dd021c 100644 --- a/src/db/unit_tests/dbTrianglesTests.cc +++ b/src/db/unit_tests/dbTrianglesTests.cc @@ -23,6 +23,7 @@ #include "dbTriangles.h" #include "dbWriter.h" +#include "dbRegionProcessors.h" #include "tlUnitTest.h" #include "tlStream.h" #include "tlFileUtils.h" @@ -917,3 +918,39 @@ TEST(triangulate_problematic) EXPECT_LT (tri.num_triangles (), size_t (490)); } +TEST(triangulate_thin) +{ + db::DPoint contour[] = { + db::DPoint (18790, 58090), + db::DPoint (18790, 58940), + db::DPoint (29290, 58940), + db::DPoint (29290, 58090) + }; + + db::DPoint hole[] = { + db::DPoint (18791, 58091), + db::DPoint (29289, 58091), + db::DPoint (29289, 58939), + db::DPoint (18791, 58939) + }; + + db::DPolygon poly; + poly.assign_hull (contour + 0, contour + sizeof (contour) / sizeof (contour[0])); + poly.insert_hole (hole + 0, hole + sizeof (hole) / sizeof (hole[0])); + + double dbu = 0.001; + + db::Triangles::TriangulateParameters param; + param.min_b = 0.5; + param.max_area = 0.0; + param.min_length = 2 * dbu; + + TestableTriangles tri; + db::DCplxTrans trans = db::DCplxTrans (dbu) * db::DCplxTrans (db::DTrans (db::DPoint () - poly.box ().center ())); + tri.triangulate (trans * poly, param); + + EXPECT_EQ (tri.check (false), true); + + // for debugging: + tri.dump ("debug.gds"); +} diff --git a/src/db/unit_tests/unit_tests.pro b/src/db/unit_tests/unit_tests.pro index a24390237..ecb3c7ccb 100644 --- a/src/db/unit_tests/unit_tests.pro +++ b/src/db/unit_tests/unit_tests.pro @@ -14,6 +14,7 @@ SOURCES = \ dbObjectWithPropertiesTests.cc \ dbPolygonNeighborhoodTests.cc \ dbPropertiesFilterTests.cc \ + dbQuadTreeTests.cc \ dbRecursiveInstanceIteratorTests.cc \ dbRegionCheckUtilsTests.cc \ dbTriangleTests.cc \