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/dbDeepTexts.cc b/src/db/db/dbDeepTexts.cc index cabf79bfe..107a3f6ed 100644 --- a/src/db/db/dbDeepTexts.cc +++ b/src/db/db/dbDeepTexts.cc @@ -462,7 +462,7 @@ DeepTexts::apply_filter (const TextFilterBase &filter, bool with_true, bool with for (db::Shapes::shape_iterator si = s.begin (db::ShapeIterator::Texts); ! si.at_end (); ++si) { db::Text text; si->text (text); - if (filter.selected (text.transformed (*v), si->prop_id ())) { + if (filter.selected (text.transformed (tr), si->prop_id ())) { if (st_true) { st_true->insert (*si); } 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..45d1c0a53 --- /dev/null +++ b/src/db/db/dbQuadTree.h @@ -0,0 +1,825 @@ + +/* + + 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 "tlLog.h" +#include + +namespace db +{ + +/** + * @brief The quad tree node implementation class + */ +template +class quad_tree_node +{ +public: + typedef typename BC::box_type box_type; + typedef typename box_type::point_type point_type; + typedef typename box_type::vector_type vector_type; + typedef std::vector objects_vector; + typedef db::coord_traits coord_traits; + + quad_tree_node (const point_type ¢er) + : m_center (center) + { + init (true); + } + + ~quad_tree_node () + { + clear (); + } + + quad_tree_node (const quad_tree_node &other) + : m_center (center) + { + init (true); + operator= (other); + } + + quad_tree_node &operator= (const quad_tree_node &other) + { + if (this != &other) { + clear (); + m_center = other.m_center; + m_objects = other.m_objects; + if (! other.is_leaf ()) { + init (false); + for (unsigned int i = 0; i < 4; ++i) { + if (other.m_q[i]) { + m_q[i] = other.m_q[i]->clone (); + } + } + } + } + return *this; + } + + quad_tree_node (quad_tree_node &&other) + { + init (true); + swap (other); + } + + quad_tree_node &operator= (quad_tree_node &&other) + { + swap (other); + return *this; + } + + void swap (quad_tree_node &other) + { + if (this != &other) { + std::swap (m_center, other.m_center); + m_objects.swap (other.m_objects); + for (unsigned int i = 0; i < 4; ++i) { + std::swap (m_q[i], other.m_q[i]); + } + } + } + + quad_tree_node *clone () const + { + quad_tree_node *node = new quad_tree_node (m_center); + *node = *this; + return node; + } + + void clear () + { + m_objects.clear (); + if (! is_leaf ()) { + for (unsigned int i = 0; i < 4; ++i) { + if (m_q[i]) { + delete m_q[i]; + } + m_q[i] = 0; + } + init (true); + } + } + + const point_type ¢er () const + { + return m_center; + } + + void insert_top (const T &value, const box_type &total_box, const box_type &b) + { + insert (value, propose_ucenter (total_box), b); + } + + bool erase (const T &value, const box_type &b) + { + int n = quad_for (b); + + if (is_leaf () || n < 0) { + + for (auto i = m_objects.begin (); i != m_objects.end (); ++i) { + if (CMP () (*i, value)) { + m_objects.erase (i); + return true; + } + } + + } else if (m_q[n]) { + + if (m_q[n]->erase (value, b)) { + if (m_q[n]->empty ()) { + delete m_q[n]; + m_q[n] = 0; + } + return true; + } + + } + + return false; + } + + const objects_vector &objects () const + { + return m_objects; + } + + box_type q_box (unsigned int n) const + { + if (! is_leaf () && m_q[n]) { + return m_q[n]->box (m_center); + } else { + return box_type (); + } + } + + quad_tree_node *node (unsigned int n) const + { + tl_assert (! is_leaf ()); + return m_q [n]; + } + + bool empty () const + { + if (m_objects.empty ()) { + if (! is_leaf ()) { + 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 (); + if (! is_leaf ()) { + for (unsigned int n = 0; n < 4; ++n) { + if (m_q[n]) { + count += m_q[n]->size (); + } + } + } + return count; + } + + size_t levels () const + { + size_t l = 1; + if (! is_leaf ()) { + for (unsigned int n = 0; n < 4; ++n) { + if (m_q[n]) { + l = std::max (l, m_q[n]->levels () + 1); + } + } + } + return l; + } + + bool check_top (const box_type &total_box) const + { + return check (propose_ucenter (total_box)); + } + +private: + point_type m_center; + quad_tree_node *m_q [4]; + objects_vector m_objects; + + void init (bool is_leaf) + { + for (unsigned int i = 0; i < 4; ++i) { + m_q[i] = 0; + } + if (is_leaf) { + m_q[0] = reinterpret_cast (1); + } + } + + bool is_leaf () const + { + return m_q[0] == reinterpret_cast (1); + } + + int quad_for (const box_type &box) const + { + int sx = coord_traits::less (box.right (), m_center.x ()) ? 0 : (coord_traits::less (m_center.x (), box.left ()) ? 1 : -1); + int sy = coord_traits::less (box.top (), m_center.y ()) ? 0 : (coord_traits::less (m_center.y (), box.bottom ()) ? 2 : -1); + if (sx < 0 || sy < 0) { + return -1; + } else { + return sx + sy; + } + } + + 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) + { + init (false); + + objects_vector ov; + ov.swap (m_objects); + + for (auto o = ov.begin (); o != ov.end (); ++o) { + insert (*o, ucenter, BC () (*o)); + } + } + + void insert (const T &value, const point_type &ucenter, const box_type &b) + { + if (is_leaf () && m_objects.size () + 1 < thr) { + + m_objects.push_back (value); + + } else { + + if (is_leaf ()) { + split (ucenter); + } + + if (inside (b, box (ucenter))) { + + int n = quad_for (b); + if (n < 0) { + m_objects.push_back (value); + } else { + if (! m_q[n]) { + box_type bq = q (n, ucenter); + m_q[n] = new quad_tree_node (bq.center ()); + } + m_q[n]->insert (value, m_center, b); + } + + } else { + + tl_assert (m_q[0] || m_q[1] || m_q[2] || m_q[3]); + point_type new_ucenter = m_center - (m_center - ucenter) * 2.0; + grow (new_ucenter); + insert (value, new_ucenter, b); + + } + + } + } + + 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]->init (false); + m_q[i]->m_q[3 - i] = n; + } + } + } + + point_type propose_ucenter (const box_type &total_box) const + { + if (! is_leaf ()) { + for (unsigned int i = 0; i < 4; ++i) { + if (m_q[i]) { + return m_center - (m_center - m_q[i]->center ()) * 2.0; + } + } + } + + typename box_type::coord_type dx = std::max (std::abs (total_box.left () - m_center.x ()), std::abs (total_box.right () - m_center.y ())); + typename box_type::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 * 2, dy * 2); + } + + bool check (const point_type &ucenter) const + { + bool result = true; + + box_type bq = box (ucenter); + + for (auto i = m_objects.begin (); i != m_objects.end (); ++i) { + box_type b = BC () (*i); + if (! inside (b, bq)) { + tl::error << "Box " << b.to_string () << " not inside quad box " << bq.to_string (); + result = false; + } + point_type ucenter_opp = m_center + (m_center - ucenter); + if (coord_traits::equal (b.left (), ucenter.x ()) || coord_traits::equal (b.right (), ucenter.x ()) || + coord_traits::equal (b.bottom (), ucenter.y ()) || coord_traits::equal (b.top (), ucenter.y ()) || + coord_traits::equal (b.left (), ucenter_opp.x ()) || coord_traits::equal (b.right (), ucenter_opp.x ()) || + coord_traits::equal (b.bottom (), ucenter_opp.y ()) || coord_traits::equal (b.top (), ucenter_opp.y ())) { + tl::error << "Box " << b.to_string () << " touches quad boundary " << ucenter.to_string () << " .. " << ucenter_opp.to_string (); + result = false; + } + } + + if (! is_leaf ()) { + + for (auto i = m_objects.begin (); i != m_objects.end (); ++i) { + box_type b = BC () (*i); + int n = quad_for (b); + if (n >= 0) { + tl::error << "Box " << b.to_string () << " on quad level not overlapping multiple quads"; + result = false; + } + } + + for (unsigned int n = 0; n < 4; ++n) { + if (m_q[n]) { + if (! m_q[n]->check (m_center)) { + result = false; + } + box_type bbq = m_q[n]->box (m_center); + if (! bbq.equal (q (n, ucenter))) { + tl::error << "Quad not centered (quad box is " << bbq.to_string () << ", should be " << q (n, ucenter).to_string (); + result = false; + } + } + } + + } else { + + if (m_objects.size () > thr) { + tl::error << "Non-split object count exceeds threshold " << m_objects.size () << " > " << thr; + result = false; + } + + } + + return result; + } + + static bool inside (const box_type &box, const box_type &in) + { + if (box.empty () || in.empty ()) { + return false; + } else { + return coord_traits::less (in.left (), box.left ()) && coord_traits::less (box.right (), in.right ()) && + coord_traits::less (in.bottom (), box.bottom ()) && coord_traits::less (box.top (), in.top ()); + } + } +}; + +/** + * @brief The iterator implementation class + */ +template +class quad_tree_iterator +{ +public: + typedef quad_tree_node quad_tree_node_type; + typedef typename BC::box_type 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; + } + } + + m_stack.pop_back (); + + while (! m_stack.empty ()) { + + qn = m_stack.back ().first; + + 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; + } + } + + m_stack.pop_back (); + + } + } +}; + +/** + * @brief The selector for implementing the all-iterator + */ +template +class quad_tree_always_sel +{ +public: + typedef typename BC::box_type box_type; + + bool select (const T &) const + { + return true; + } + + bool select_quad (const box_type &) const + { + return true; + } +}; + +/** + * @brief The selector for implementing the touching iterator + */ +template +class quad_tree_touching_sel +{ +public: + typedef typename BC::box_type 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; +}; + +/** + * @brief The selector for implementing the overlapping iterator + */ +template +class quad_tree_overlapping_sel +{ +public: + typedef typename BC::box_type 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; +}; + +/** + * @brief The default compare function + */ +template +struct quad_tree_default_cmp +{ + bool operator() (const T &a, const T &b) const + { + return a == b; + } +}; + +/** + * @brief A generic quad tree implementation + * + * In contrast to the box_tree implementation, this is a self-sorting implementation + * which is more generic. + * + * @param T The value to be stored + * @param BC The box converter + * @param thr The number of items per leaf node before splitting + * @param CMP The compare function (equality) + * + * T needs to have a type member named "coord_type". + * BC is a function of T and delivers a db::box for T + * CMP is a function that delivers a bool value for a pair of T (equality) + */ +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 BC::box_type box_type; + typedef typename box_type::point_type point_type; + typedef typename box_type::vector_type vector_type; + typedef std::vector objects_vector; + + /** + * @brief Default constructor + * + * This creates an empty tree. + */ + quad_tree () + : m_root (point_type ()) + { + // .. nothing yet .. + } + + /** + * @brief Copy constructor + */ + quad_tree (const quad_tree &other) + : m_root (point_type ()) + { + operator= (other); + } + + /** + * @brief Assignment + */ + quad_tree &operator= (const quad_tree &other) + { + if (this != &other) { + m_root = other.m_root; + m_total_box = other.m_total_box; + } + return *this; + } + + /** + * @brief Move constructor + */ + quad_tree (quad_tree &&other) + : m_root (point_type ()) + { + swap (other); + } + + /** + * @brief Move assignment + */ + quad_tree &operator= (quad_tree &&other) + { + swap (other); + return *this; + } + + /** + * @brief Empties the tree + */ + void clear () + { + m_root.clear (); + m_total_box = box_type (); + } + + /** + * @brief Swaps the tree with another + */ + void swap (quad_tree &other) + { + if (this != &other) { + m_root.swap (other.m_root); + std::swap (m_total_box, other.m_total_box); + } + } + + /** + * @brief Returns a value indicating whether the tree is empty + */ + bool empty () const + { + return m_root.empty (); + } + + /** + * @brief Returns the number of items stored in the tree + */ + size_t size () const + { + return m_root.size (); + } + + /** + * @brief Returns the number of quad levels (for testing) + */ + size_t levels () const + { + return m_root.levels (); + } + + /** + * @brief Checks the tree for consistency (for testing) + */ + bool check () const + { + return m_root.check_top (m_total_box); + } + + /** + * @brief Inserts an object into the tree + */ + 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, b); + } + + /** + * @brief Erases the given element from the tree + * + * @return true, if the element was found and erased. + * + * If multiple elements of the same kind are stored, the + * first one is erased. + */ + bool erase (const T &value) + { + return m_root.erase (value, BC () (value)); + } + + /** + * @brief begin iterator for all elements + */ + quad_tree_flat_iterator begin () const + { + return quad_tree_flat_iterator (&m_root, quad_tree_always_sel ()); + } + + /** + * @brief begin iterator for all elements overlapping the given box + */ + 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 (); + } + } + + /** + * @brief begin iterator for all elements touching the given box + */ + 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..867187295 100644 --- a/src/db/db/dbTriangle.cc +++ b/src/db/db/dbTriangle.cc @@ -363,21 +363,33 @@ Triangle::Triangle (TriangleEdge *e1, TriangleEdge *e2, TriangleEdge *e3) } mp_v[2] = mp_e[1]->other (mp_v[1]); - // establish link to edges - for (int i = 0; i < 3; ++i) { - TriangleEdge *e = mp_e[i]; - int side_of = e->side_of (*mp_v[i == 0 ? 2 : i - 1]); - // NOTE: in the degenerated case, the triangle is not attached to an edge! - if (side_of < 0) { - e->set_left (this); - } else if (side_of > 0) { - e->set_right (this); - } + // enforce clockwise orientation + int s = db::vprod_sign (*mp_v[2] - *mp_v[0], *mp_v[1] - *mp_v[0]); + if (s < 0) { + std::swap (mp_v[2], mp_v[1]); + } else if (s == 0) { + // Triangle is not orientable + tl_assert (false); } - // enforce clockwise orientation - if (db::vprod_sign (*mp_v[2] - *mp_v[0], *mp_v[1] - *mp_v[0]) < 0) { - std::swap (mp_v[2], mp_v[1]); + // establish link to edges + for (int i = 0; i < 3; ++i) { + + TriangleEdge *e = mp_e[i]; + + unsigned int i1 = 0; + for ( ; e->v1 () != mp_v[i1] && i1 < 3; ++i1) + ; + unsigned int i2 = 0; + for ( ; e->v2 () != mp_v[i2] && i2 < 3; ++i2) + ; + + if ((i1 + 1) % 3 == i2) { + e->set_right (this); + } else { + e->set_left (this); + } + } } @@ -436,22 +448,40 @@ Triangle::bbox () const std::pair -Triangle::circumcircle () const +Triangle::circumcircle (bool *ok) const { - db::DVector v1 = *mp_v[0] - *mp_v[1]; - db::DVector v2 = *mp_v[0] - *mp_v[2]; - db::DVector n1 = db::DVector (v1.y (), -v1.x ()); - db::DVector n2 = db::DVector (v2.y (), -v2.x ()); + // see https://en.wikipedia.org/wiki/Circumcircle + // we set A=(0,0), so the formulas simplify - double p1s = v1.sq_length (); - double p2s = v2.sq_length (); + if (ok) { + *ok = true; + } - double s = db::vprod (v1, v2); - tl_assert (fabs (s) > db::epsilon); + db::DVector b = *mp_v[1] - *mp_v[0]; + db::DVector c = *mp_v[2] - *mp_v[0]; - db::DVector r = (n1 * p2s - n2 * p1s) * (0.5 / s); - db::DPoint center = *mp_v[0] + r; - double radius = r.length (); + double b2 = b.sq_length (); + double c2 = c.sq_length (); + + double sx = 0.5 * (b2 * c.y () - c2 * b.y ()); + double sy = 0.5 * (b.x () * c2 - c.x() * b2); + + double a1 = b.x() * c.y(); + double a2 = c.x() * b.y(); + double a = a1 - a2; + double a_abs = std::abs (a); + + if (a_abs < (std::abs (a1) + std::abs (a2)) * db::epsilon) { + if (ok) { + *ok = false; + return std::make_pair (db::DPoint (), 0.0); + } else { + tl_assert (false); + } + } + + double radius = sqrt (sx * sx + sy * sy) / a_abs; + db::DPoint center = *mp_v[0] + db::DVector (sx / a, sy / a); return std::make_pair (center, radius); } @@ -507,18 +537,28 @@ Triangle::common_edge (const Triangle *other) const int Triangle::contains (const db::DPoint &point) const { + auto c = *mp_v[2] - *mp_v[0]; + auto b = *mp_v[1] - *mp_v[0]; + + int vps = db::vprod_sign (c, b); + if (vps == 0) { + return db::vprod_sign (point - *mp_v[0], b) == 0 && db::vprod_sign (point - *mp_v[0], c) == 0 ? 0 : -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]; + int n = db::vprod_sign (point - *vl, *v - *vl) * vps; + if (n < 0) { return -1; + } else if (n == 0) { + res = 0; } vl = v; } + return res; } @@ -536,8 +576,9 @@ double Triangle::b () const { double lmin = min_edge_length (); - auto cr = circumcircle (); - return lmin / cr.second; + bool ok = false; + auto cr = circumcircle (&ok); + return ok ? lmin / cr.second : 0.0; } bool diff --git a/src/db/db/dbTriangle.h b/src/db/db/dbTriangle.h index dadcabc4a..29ff09abb 100644 --- a/src/db/db/dbTriangle.h +++ b/src/db/db/dbTriangle.h @@ -483,8 +483,10 @@ public: /** * @brief Gets the center point and radius of the circumcircle + * If ok is non-null, it will receive a boolean value indicating whether the circumcircle is valid. + * An invalid circumcircle is an indicator for a degenerated triangle with area 0 (or close to). */ - std::pair circumcircle () const; + std::pair circumcircle (bool *ok = 0) const; /** * @brief Gets the vertex opposite of the given edge diff --git a/src/db/db/dbTriangles.cc b/src/db/db/dbTriangles.cc index 9837f5104..23710b681 100644 --- a/src/db/db/dbTriangles.cc +++ b/src/db/db/dbTriangles.cc @@ -36,6 +36,12 @@ namespace db { +static inline bool is_equal (const db::DPoint &a, const db::DPoint &b) +{ + return std::abs (a.x () - b.x ()) < std::max (1.0, (std::abs (a.x ()) + std::abs (b.x ()))) * db::epsilon && + std::abs (a.y () - b.y ()) < std::max (1.0, (std::abs (a.y ()) + std::abs (b.y ()))) * db::epsilon; +} + Triangles::Triangles () : m_is_constrained (false), m_level (0), m_id (0), m_flips (0), m_hops (0) { @@ -44,9 +50,7 @@ Triangles::Triangles () Triangles::~Triangles () { - while (! mp_triangles.empty ()) { - remove_triangle (mp_triangles.begin ().operator-> ()); - } + clear (); } db::Vertex * @@ -88,6 +92,7 @@ Triangles::create_triangle (TriangleEdge *e1, TriangleEdge *e2, TriangleEdge *e3 db::Triangle *res = new db::Triangle (e1, e2, e3); res->set_id (++m_id); mp_triangles.push_back (res); + return res; } @@ -366,26 +371,35 @@ Triangles::insert (db::Vertex *vertex, std::list > *n // check, if the new vertex is on an edge (may be edge between triangles or edge on outside) std::vector on_edges; + std::vector on_vertex; for (int i = 0; i < 3; ++i) { db::TriangleEdge *e = tris.front ()->edge (i); if (e->side_of (*vertex) == 0) { - on_edges.push_back (e); + if (is_equal (*vertex, *e->v1 ()) || is_equal (*vertex, *e->v2 ())) { + on_vertex.push_back (e); + } else { + on_edges.push_back (e); + } } } - if (! on_edges.empty ()) { - if (on_edges.size () == size_t (1)) { - split_triangles_on_edge (tris, vertex, on_edges.front (), new_triangles); - return vertex; - } else { - // the vertex is already present - tl_assert (on_edges.size () == size_t (2)); - return on_edges.front ()->common_vertex (on_edges [1]); - } + if (! on_vertex.empty ()) { + + tl_assert (on_vertex.size () == size_t (2)); + return on_vertex.front ()->common_vertex (on_vertex [1]); + + } else if (! on_edges.empty ()) { + + tl_assert (on_edges.size () == size_t (1)); + split_triangles_on_edge (vertex, on_edges.front (), new_triangles); + return vertex; + } else if (tris.size () == size_t (1)) { + // the new vertex is inside one triangle split_triangle (tris.front (), vertex, new_triangles); return vertex; + } tl_assert (false); @@ -404,6 +418,7 @@ Triangles::find_triangle_for_point (const db::DPoint &point) } } } + return res; } @@ -642,7 +657,7 @@ Triangles::split_triangle (db::Triangle *t, db::Vertex *vertex, std::list &tris, db::Vertex *vertex, db::TriangleEdge *split_edge, std::list > *new_triangles_out) +Triangles::split_triangles_on_edge (db::Vertex *vertex, db::TriangleEdge *split_edge, std::list > *new_triangles_out) { TriangleEdge *s1 = create_edge (split_edge->v1 (), vertex); TriangleEdge *s2 = create_edge (split_edge->v2 (), vertex); @@ -651,6 +666,12 @@ Triangles::split_triangles_on_edge (const std::vector &tris, db: std::vector new_triangles; + std::vector tris; + tris.reserve (2); + for (auto t = split_edge->begin_triangles (); t != split_edge->end_triangles (); ++t) { + tris.push_back (t.operator-> ()); + } + for (auto t = tris.begin (); t != tris.end (); ++t) { (*t)->unlink (); @@ -931,13 +952,15 @@ Triangles::is_illegal_edge (db::TriangleEdge *edge) return false; } - auto lr = left->circumcircle (); - if (right->opposite (edge)->in_circle (lr.first, lr.second) > 0) { + bool ok = false; + + auto lr = left->circumcircle (&ok); + if (! ok || right->opposite (edge)->in_circle (lr.first, lr.second) > 0) { return true; } - auto rr = right->circumcircle(); - if (left->opposite (edge)->in_circle (rr.first, rr.second) > 0) { + auto rr = right->circumcircle(&ok); + if (! ok || left->opposite (edge)->in_circle (rr.first, rr.second) > 0) { return true; } @@ -1163,16 +1186,16 @@ Triangles::search_edges_crossing (Vertex *from, Vertex *to) } db::Vertex * -Triangles::find_vertex_for_point (const db::DPoint &pt) +Triangles::find_vertex_for_point (const db::DPoint &point) { - db::TriangleEdge *edge = find_closest_edge (pt); + db::TriangleEdge *edge = find_closest_edge (point); if (!edge) { return 0; } db::Vertex *v = 0; - if (edge->v1 ()->equal (pt)) { + if (is_equal (*edge->v1 (), point)) { v = edge->v1 (); - } else if (edge->v2 ()->equal (pt)) { + } else if (is_equal (*edge->v2 (), point)) { v = edge->v2 (); } return v; @@ -1186,7 +1209,7 @@ Triangles::find_edge_for_points (const db::DPoint &p1, const db::DPoint &p2) return 0; } for (auto e = v->begin_edges (); e != v->end_edges (); ++e) { - if ((*e)->other (v)->equal (p2)) { + if (is_equal (*(*e)->other (v), p2)) { return *e; } } @@ -1329,7 +1352,8 @@ 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); } } @@ -1632,7 +1656,28 @@ Triangles::refine (const TriangulateParameters ¶meters) auto cr = (*t)->circumcircle(); auto center = cr.first; - if ((*t)->contains (center) >= 0) { + int s = (*t)->contains (center); + if (s >= 0) { + + if (s > 0) { + + double snap = 1e-3; + + // Snap the center to a segment center if "close" to it. + // This avoids generating very skinny triangles that can't be fixed as the + // segment cannot be flipped. This a part of the issue #1996 problem. + for (unsigned int i = 0; i < 3; ++i) { + if ((*t)->edge (i)->is_segment ()) { + auto e = (*t)->edge (i)->edge (); + auto c = e.p1 () + e.d () * 0.5; + if (c.distance (center) < e.length () * 0.5 * snap - db::epsilon) { + center = c; + break; + } + } + } + + } if (tl::verbosity () >= parameters.base_verbosity + 20) { tl::info << "Inserting in-triangle center " << center.to_string () << " of " << (*t)->to_string (true); @@ -1642,7 +1687,7 @@ Triangles::refine (const TriangulateParameters ¶meters) } else { db::Vertex *vstart = 0; - for (int i = 0; i < 3; ++i) { + for (unsigned int i = 0; i < 3; ++i) { db::TriangleEdge *edge = (*t)->edge (i); vstart = (*t)->opposite (edge); if (edge->side_of (*vstart) * edge->side_of (center) < 0) { diff --git a/src/db/db/dbTriangles.h b/src/db/db/dbTriangles.h index f026ce820..31da3f48b 100644 --- a/src/db/db/dbTriangles.h +++ b/src/db/db/dbTriangles.h @@ -20,8 +20,6 @@ */ - - #ifndef HDR_dbTriangles #define HDR_dbTriangles @@ -309,6 +307,15 @@ protected: std::vector find_inside_circle (const db::DPoint ¢er, double radius) const; private: + struct TriangleBoxConvert + { + typedef db::DBox box_type; + box_type operator() (db::Triangle *t) const + { + return t ? t->bbox () : box_type (); + } + }; + tl::list mp_triangles; tl::stable_vector m_edges_heap; std::vector m_returned_edges; @@ -332,7 +339,7 @@ private: db::TriangleEdge *find_closest_edge (const db::DPoint &p, db::Vertex *vstart = 0, bool inside_only = false); db::Vertex *insert (db::Vertex *vertex, std::list > *new_triangles = 0); void split_triangle (db::Triangle *t, db::Vertex *vertex, std::list > *new_triangles_out); - void split_triangles_on_edge (const std::vector &tris, db::Vertex *vertex, db::TriangleEdge *split_edge, std::list > *new_triangles_out); + void split_triangles_on_edge (db::Vertex *vertex, db::TriangleEdge *split_edge, std::list > *new_triangles_out); void add_more_triangles (std::vector &new_triangles, db::TriangleEdge *incoming_edge, db::Vertex *from_vertex, db::Vertex *to_vertex, diff --git a/src/db/unit_tests/dbQuadTreeTests.cc b/src/db/unit_tests/dbQuadTreeTests.cc new file mode 100644 index 000000000..d92e07237 --- /dev/null +++ b/src/db/unit_tests/dbQuadTreeTests.cc @@ -0,0 +1,595 @@ + +/* + + 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" +#include "tlTimer.h" + +#include + +struct MyQuadTreeCMP +{ + bool operator() (const db::DBox &a, const db::DBox &b) const + { + return a.equal (b); + } +}; + +typedef db::quad_tree, size_t (1), MyQuadTreeCMP> MyQuadTree; + +std::string find_all (const MyQuadTree &qt) +{ + std::vector v; + auto i = qt.begin (); + while (! i.at_end ()) { + v.push_back (i->to_string ()); + ++i; + } + std::sort (v.begin (), v.end ()); + return tl::join (v, "/"); +} + +std::string find_touching (const MyQuadTree &qt, const db::DBox &box, bool report = false) +{ + std::vector v; + auto i = qt.begin_touching (box); + while (! i.at_end ()) { + v.push_back (i->to_string ()); + ++i; + } + if (report) { + tl::info << v.size () << " items found."; + } + std::sort (v.begin (), v.end ()); + return tl::join (v, "/"); +} + +std::string find_touching_from_all (const MyQuadTree &qt, const db::DBox &box) +{ + std::vector v; + auto i = qt.begin (); + while (! i.at_end ()) { + if (i->touches (box)) { + v.push_back (i->to_string ()); + } + ++i; + } + std::sort (v.begin (), v.end ()); + return tl::join (v, "/"); +} + +std::string find_overlapping (const MyQuadTree &qt, const db::DBox &box, bool report = false) +{ + std::vector v; + auto i = qt.begin_overlapping (box); + while (! i.at_end ()) { + v.push_back (i->to_string ()); + ++i; + } + if (report) { + tl::info << v.size () << " items found."; + } + std::sort (v.begin (), v.end ()); + return tl::join (v, "/"); +} + +std::string find_overlapping_from_all (const MyQuadTree &qt, const db::DBox &box) +{ + std::vector v; + auto i = qt.begin (); + while (! i.at_end ()) { + if (i->overlaps (box)) { + v.push_back (i->to_string ()); + } + ++i; + } + std::sort (v.begin (), v.end ()); + return tl::join (v, "/"); +} + +TEST(basic) +{ + MyQuadTree tree; + EXPECT_EQ (tree.empty (), true); + EXPECT_EQ (tree.size (), size_t (0)); + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (tree.levels (), size_t (1)); + + tree.insert (db::DBox ()); + EXPECT_EQ (tree.empty (), true); + EXPECT_EQ (tree.size (), size_t (0)); + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (tree.levels (), size_t (1)); + + tree.insert (db::DBox (-1, -2, 3, 4)); + EXPECT_EQ (tree.empty (), false); + EXPECT_EQ (tree.size (), size_t (1)); + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (tree.levels (), size_t (1)); + + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)"); + + db::DBox bx; + + bx = db::DBox (-2, 0, -1, 0); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, -3, -1, -2); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, -3, -1, -2.5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 4, -1, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 4.5, -1, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 3, -1, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 3, -1.5, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + + tree.insert (db::DBox (-1, -3, 3, 0)); + EXPECT_EQ (tree.empty (), false); + EXPECT_EQ (tree.size (), size_t (2)); + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (tree.levels (), size_t (1)); + + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;3,0)"); + + bx = db::DBox (-2, 0, -1, 0); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, -3, -1, -2); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, -3, -1, -2.5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 4, -1, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 4.5, -1, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 3, -1, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 3, -1.5, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + + tree.insert (db::DBox (-1, -3, -0.5, -2)); + EXPECT_EQ (tree.empty (), false); + EXPECT_EQ (tree.size (), size_t (3)); + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (tree.levels (), size_t (3)); + + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;3,0)"); + + bx = db::DBox (-2, 0, -1, 0); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, -3, -1, -2); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, -3, -1, -2.5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 4, -1, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 4.5, -1, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 3, -1, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 3, -1.5, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + + tree.insert (db::DBox (-1, -3, -0.5, 2)); + EXPECT_EQ (tree.empty (), false); + EXPECT_EQ (tree.size (), size_t (4)); + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (tree.levels (), size_t (3)); + + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + + bx = db::DBox (-2, 0, -1, 0); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, -3, -1, -2); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, -3, -1, -2.5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 4, -1, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 4.5, -1, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 3, -1, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); + bx = db::DBox (-2, 3, -1.5, 5); + EXPECT_EQ (find_touching (tree, bx), find_touching_from_all (tree, bx)); + EXPECT_EQ (find_overlapping (tree, bx), find_overlapping_from_all (tree, bx)); +} + +TEST(remove) +{ + MyQuadTree tree; + tree.insert (db::DBox (-1, -2, 3, 4)); + tree.insert (db::DBox (-1, -3, 3, 0)); + tree.insert (db::DBox (-1, -3, -0.5, -2)); + tree.insert (db::DBox (-1, -3, -0.5, 2)); + + EXPECT_EQ (tree.check (), true); + + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + + EXPECT_EQ (tree.erase (db::DBox (-1, -3, -0.5, -1)), false); + EXPECT_EQ (tree.erase (db::DBox (-1, -3, -0.5, -2)), true); + EXPECT_EQ (tree.check (), true); + + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + + while (! tree.empty ()) { + EXPECT_EQ (tree.erase (*tree.begin ()), true); + EXPECT_EQ (tree.check (), true); + } + + EXPECT_EQ (tree.size (), size_t (0)); + EXPECT_EQ (tree.levels (), size_t (1)); +} + +TEST(grow) +{ + MyQuadTree tree; + tree.insert (db::DBox (-1, -2, 3, 4)); + tree.insert (db::DBox (-1, -3, 3, 0)); + tree.insert (db::DBox (-1, -3, -0.5, -2)); + tree.insert (db::DBox (-1, -3, -0.5, 2)); + EXPECT_EQ (tree.levels (), size_t (3)); + tree.insert (db::DBox (-100, -3, -99, 2)); + EXPECT_EQ (tree.levels (), size_t (8)); + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)/(-100,-3;-99,2)"); + EXPECT_EQ (find_overlapping (tree, db::DBox (-100, -100, -90, 100)), "(-100,-3;-99,2)"); + + bool r = true; + while (r && ! tree.empty ()) { + r = tree.erase (*tree.begin ()); + EXPECT_EQ (r, true); + EXPECT_EQ (tree.check (), true); + } + + EXPECT_EQ (tree.size (), size_t (0)); + EXPECT_EQ (tree.levels (), size_t (1)); +} + +TEST(grow2) +{ + MyQuadTree tree; + tree.insert (db::DBox (-1, -2, 3, 4)); + tree.insert (db::DBox (-1, -3, 3, 0)); + tree.insert (db::DBox (-1, -3, -0.5, -2)); + tree.insert (db::DBox (-1, -3, -0.5, 2)); + EXPECT_EQ (tree.levels (), size_t (3)); + tree.insert (db::DBox (-100, -3, -99, -1)); + EXPECT_EQ (tree.levels (), size_t (8)); + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)/(-100,-3;-99,-1)"); + EXPECT_EQ (find_overlapping (tree, db::DBox (-100, -100, -90, 100)), "(-100,-3;-99,-1)"); + + bool r = true; + while (r && ! tree.empty ()) { + r = tree.erase (*tree.begin ()); + EXPECT_EQ (r, true); + EXPECT_EQ (tree.check (), true); + } + + EXPECT_EQ (tree.size (), size_t (0)); + EXPECT_EQ (tree.levels (), size_t (1)); +} + +TEST(clear) +{ + MyQuadTree tree; + tree.insert (db::DBox (-1, -2, 3, 4)); + tree.insert (db::DBox (-1, -3, 3, 0)); + tree.insert (db::DBox (-1, -3, -0.5, -2)); + tree.insert (db::DBox (-1, -3, -0.5, 2)); + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + + tree.clear (); + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (tree.empty (), true); + EXPECT_EQ (tree.size (), size_t (0)); + EXPECT_EQ (tree.levels (), size_t (1)); + EXPECT_EQ (find_all (tree), ""); +} + +TEST(copy) +{ + MyQuadTree tree; + tree.insert (db::DBox (-1, -2, 3, 4)); + tree.insert (db::DBox (-1, -3, 3, 0)); + tree.insert (db::DBox (-1, -3, -0.5, -2)); + tree.insert (db::DBox (-1, -3, -0.5, 2)); + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + EXPECT_EQ (tree.levels (), size_t (3)); + + MyQuadTree tree2 (tree); + tree.clear (); + + EXPECT_EQ (tree2.check (), true); + EXPECT_EQ (find_all (tree2), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + EXPECT_EQ (tree2.levels (), size_t (3)); +} + +TEST(assign) +{ + MyQuadTree tree; + tree.insert (db::DBox (-1, -2, 3, 4)); + tree.insert (db::DBox (-1, -3, 3, 0)); + tree.insert (db::DBox (-1, -3, -0.5, -2)); + tree.insert (db::DBox (-1, -3, -0.5, 2)); + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + EXPECT_EQ (tree.levels (), size_t (3)); + + MyQuadTree tree2; + tree2 = tree; + tree.clear (); + + EXPECT_EQ (tree2.check (), true); + EXPECT_EQ (find_all (tree2), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + EXPECT_EQ (tree2.levels (), size_t (3)); +} + +TEST(swap) +{ + MyQuadTree tree; + tree.insert (db::DBox (-1, -2, 3, 4)); + tree.insert (db::DBox (-1, -3, 3, 0)); + tree.insert (db::DBox (-1, -3, -0.5, -2)); + tree.insert (db::DBox (-1, -3, -0.5, 2)); + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + EXPECT_EQ (tree.levels (), size_t (3)); + + MyQuadTree tree2; + tree2.swap (tree); + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (tree.empty (), true); + EXPECT_EQ (find_all (tree), ""); + EXPECT_EQ (tree.levels (), size_t (1)); + + EXPECT_EQ (tree2.check (), true); + EXPECT_EQ (find_all (tree2), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + EXPECT_EQ (tree2.levels (), size_t (3)); +} + +TEST(move) +{ + MyQuadTree tree; + tree.insert (db::DBox (-1, -2, 3, 4)); + tree.insert (db::DBox (-1, -3, 3, 0)); + tree.insert (db::DBox (-1, -3, -0.5, -2)); + tree.insert (db::DBox (-1, -3, -0.5, 2)); + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + EXPECT_EQ (tree.levels (), size_t (3)); + + MyQuadTree tree2; + tree2 = std::move (tree); + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (tree.empty (), true); + EXPECT_EQ (find_all (tree), ""); + EXPECT_EQ (tree.levels (), size_t (1)); + + EXPECT_EQ (tree2.check (), true); + EXPECT_EQ (find_all (tree2), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + EXPECT_EQ (tree2.levels (), size_t (3)); +} + +TEST(move_ctor) +{ + MyQuadTree tree; + tree.insert (db::DBox (-1, -2, 3, 4)); + tree.insert (db::DBox (-1, -3, 3, 0)); + tree.insert (db::DBox (-1, -3, -0.5, -2)); + tree.insert (db::DBox (-1, -3, -0.5, 2)); + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (find_all (tree), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + EXPECT_EQ (tree.levels (), size_t (3)); + + MyQuadTree tree2 (std::move (tree)); + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (tree.empty (), true); + EXPECT_EQ (find_all (tree), ""); + EXPECT_EQ (tree.levels (), size_t (1)); + + EXPECT_EQ (tree2.check (), true); + EXPECT_EQ (find_all (tree2), "(-1,-2;3,4)/(-1,-3;-0.5,-2)/(-1,-3;-0.5,2)/(-1,-3;3,0)"); + EXPECT_EQ (tree2.levels (), size_t (3)); +} + +static double rvalue () +{ + return ((rand () % 1000000) - 5000) * 0.001; +} + +static db::DBox rbox () +{ + db::DBox box; + while ((box = db::DBox (db::DPoint (rvalue (), rvalue ()), db::DPoint (rvalue (), rvalue ()))).empty ()) { + ; + } + return box; +} + +static db::DBox rbox (double dim) +{ + db::DBox box; + db::DPoint c (rvalue (), rvalue ()); + return box = db::DBox (c, c).enlarged (db::DVector (dim * 0.5, dim * 0.5)); +} + +TEST(many) +{ + MyQuadTree tree; + + unsigned int n = 1000; + unsigned int ntests = 100; + + for (unsigned int i = 0; i < n; ++i) { + tree.insert (rbox ()); + } + + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (tree.size (), size_t (n)); + + bool report = false; + for (unsigned int i = 0; i < ntests; ++i) { + if (report) { + tl::info << "Test iteration " << i << " ..."; + } + auto bx = rbox (); + EXPECT_EQ (find_overlapping (tree, bx, report), find_overlapping_from_all (tree, bx)); + EXPECT_EQ (find_touching (tree, bx, report), find_touching_from_all (tree, bx)); + bx = db::DBox (bx.center (), bx.center ()); + EXPECT_EQ (find_touching (tree, bx, report), find_touching_from_all (tree, bx)); + } + + bool r = true; + while (r && ! tree.empty ()) { + r = tree.erase (*tree.begin ()); + EXPECT_EQ (r, true); + EXPECT_EQ (tree.check (), true); + } + + EXPECT_EQ (tree.empty (), true); + EXPECT_EQ (tree.check (), true); + EXPECT_EQ (tree.levels (), size_t (1)); + EXPECT_EQ (tree.size (), size_t (0)); +} + +TEST(timing_insert) +{ + MyQuadTree tree; + + { + unsigned int n = 1000000; + tl::SelfTimer timer (tl::sprintf ("%d inserts ..", int (n))); + for (unsigned int i = 0; i < n; ++i) { + tree.insert (rbox ()); + } + tl::info << "Quad levels: " << tree.levels (); + } + + tree.clear (); + + { + unsigned int n = 2000000; + tl::SelfTimer timer (tl::sprintf ("%d inserts ..", int (n))); + for (unsigned int i = 0; i < n; ++i) { + tree.insert (rbox ()); + } + tl::info << "Quad levels: " << tree.levels (); + } +} + +TEST(timing_lookup) +{ + test_is_long_runner (); + + MyQuadTree tree; + + unsigned int n = 1000000; + for (unsigned int i = 0; i < n; ++i) { + tree.insert (rbox (5.0)); + } + + unsigned int ntests = 1000; + + std::vector > > tests; + for (unsigned int i = 0; i < ntests; ++i) { + db::DBox bx = rbox (5.0); + tests.push_back (std::make_pair (bx, std::make_pair (size_t (0), size_t (0)))); + } + + { + tl::SelfTimer timer (tl::sprintf ("%d tests (lookup) ..", int (ntests))); + for (auto t = tests.begin (); t != tests.end (); ++t) { + size_t n = 0; + for (auto i = tree.begin_touching (t->first); ! i.at_end (); ++i) { + ++n; + } + t->second.first = n; + } + } + + { + tl::SelfTimer timer (tl::sprintf ("%d tests (brute force) ..", int (ntests))); + for (auto t = tests.begin (); t != tests.end (); ++t) { + size_t n = 0; + for (auto i = tree.begin (); ! i.at_end (); ++i) { + if (i->touches (t->first)) { + ++n; + } + } + t->second.second = n; + } + } + + for (auto t = tests.begin (); t != tests.end (); ++t) { + EXPECT_EQ (t->second.first, t->second.second); + } +} + 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..c3855ce29 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" @@ -699,6 +700,9 @@ TEST(triangulate_basic) EXPECT_GT (tri.num_triangles (), size_t (100)); EXPECT_LT (tri.num_triangles (), size_t (150)); + // for debugging: + // tri.dump ("debug.gds"); + param.min_b = 1.0; param.max_area = 0.1; @@ -913,7 +917,86 @@ TEST(triangulate_problematic) EXPECT_GE (t->b (), param.min_b); } - EXPECT_GT (tri.num_triangles (), size_t (470)); - EXPECT_LT (tri.num_triangles (), size_t (490)); + EXPECT_GT (tri.num_triangles (), size_t (540)); + EXPECT_LT (tri.num_triangles (), size_t (560)); } +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"); + + for (auto t = tri.begin (); t != tri.end (); ++t) { + EXPECT_GE (t->b (), param.min_b); + } + + EXPECT_GT (tri.num_triangles (), size_t (13000)); + EXPECT_LT (tri.num_triangles (), size_t (13200)); +} + +TEST(triangulate_issue1996) +{ + db::DPoint contour[] = { + db::DPoint (-8000, -8075), + db::DPoint (-8000, 8075), + db::DPoint (18000, 8075), + db::DPoint (18000, -8075) + }; + + db::DPolygon poly; + poly.assign_hull (contour + 0, contour + sizeof (contour) / sizeof (contour[0])); + + double dbu = 0.001; + + db::Triangles::TriangulateParameters param; + param.min_b = 0.5; + param.max_area = 5000.0 * dbu * 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"); + + for (auto t = tri.begin (); t != tri.end (); ++t) { + EXPECT_LE (t->area (), param.max_area); + EXPECT_GE (t->b (), param.min_b); + } + + EXPECT_GT (tri.num_triangles (), size_t (128000)); + EXPECT_LT (tri.num_triangles (), size_t (132000)); +} 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 \