From b95d43b0e48d0bf366e282027dcac0f9f4f533f5 Mon Sep 17 00:00:00 2001 From: Martyn Gigg <martyn.gigg@stfc.ac.uk> Date: Thu, 14 Jul 2011 14:37:02 +0000 Subject: [PATCH] Refs #3341. This should deal with small bins in a much better way in Rebin2D and general polygon intersections. There is probably still room for improvement on the speed front. --- .../Framework/Algorithms/src/Rebin2D.cpp | 14 +- Code/Mantid/Framework/Geometry/CMakeLists.txt | 12 +- .../inc/MantidGeometry/Math/ConvexPolygon.h | 34 ++- ...gonIntersection.h => LaszloIntersection.h} | 30 +-- .../inc/MantidGeometry/Math/PolygonEdge.h | 5 +- .../inc/MantidGeometry/Math/Vertex2D.h | 22 ++ .../Geometry/src/Math/ConvexPolygon.cpp | 177 +++++--------- .../Geometry/src/Math/LaszloIntersection.cpp | 221 ++++++++++++++++++ .../Geometry/src/Math/PolygonEdge.cpp | 30 +++ .../Geometry/src/Math/PolygonIntersection.cpp | 33 --- .../Geometry/test/ConvexPolygonTest.h | 117 ++++------ .../Geometry/test/LaszloIntersectionTest.h | 175 ++++++++++++++ .../Geometry/test/PolygonIntersectionTest.h | 106 --------- .../Framework/Geometry/test/Vertex2DTest.h | 79 ++++--- .../Kernel/inc/MantidKernel/Vertex2D.h | 38 --- 15 files changed, 649 insertions(+), 444 deletions(-) rename Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/{PolygonIntersection.h => LaszloIntersection.h} (64%) create mode 100644 Code/Mantid/Framework/Geometry/src/Math/LaszloIntersection.cpp delete mode 100644 Code/Mantid/Framework/Geometry/src/Math/PolygonIntersection.cpp create mode 100644 Code/Mantid/Framework/Geometry/test/LaszloIntersectionTest.h delete mode 100644 Code/Mantid/Framework/Geometry/test/PolygonIntersectionTest.h delete mode 100644 Code/Mantid/Framework/Kernel/inc/MantidKernel/Vertex2D.h diff --git a/Code/Mantid/Framework/Algorithms/src/Rebin2D.cpp b/Code/Mantid/Framework/Algorithms/src/Rebin2D.cpp index ed1f0759c51..7765ce311b5 100644 --- a/Code/Mantid/Framework/Algorithms/src/Rebin2D.cpp +++ b/Code/Mantid/Framework/Algorithms/src/Rebin2D.cpp @@ -8,7 +8,7 @@ #include "MantidAPI/WorkspaceProperty.h" #include "MantidAPI/NumericAxis.h" #include "MantidGeometry/Math/ConvexPolygon.h" -#include "MantidGeometry/Math/PolygonIntersection.h" +#include "MantidGeometry/Math/LaszloIntersection.h" namespace Mantid { @@ -20,7 +20,6 @@ namespace Mantid using namespace API; using Geometry::ConvexPolygon; - using Geometry::Vertex2DList; using Kernel::V2D; //-------------------------------------------------------------------------- @@ -135,7 +134,7 @@ namespace Mantid std::pair<double,double> binValues(0,0); if( inputWS->isDistribution() ) { - const double newWidth = newBin[1].X() - newBin[0].X(); // For distribution + const double newWidth = newBin[3].X() - newBin[0].X(); // For distribution binValues = calculateDistYE(inputWS, overlaps, newWidth); } else @@ -212,24 +211,23 @@ namespace Mantid overlaps.reserve(5); // Have a guess at a posible limit const size_t nxpoints(oldAxis1.size()-1), nypoints(oldAxis2.size()-1); - Vertex2DList vertices(4); for(size_t i = 0; i < nypoints; ++i) { const double yo_lo(oldAxis2[i]), yo_hi(oldAxis2[i+1]); // Check if there is a possibility of overlap - if( yo_hi < newPoly[0].Y() || yo_lo > newPoly[2].Y() ) continue; + if( yo_hi < newPoly[0].Y() || yo_lo > newPoly[1].Y() ) continue; for(size_t j = 0; j < nxpoints; ++j) { const double xo_lo(oldAxis1[j]), xo_hi(oldAxis1[j+1]); // Check if there is a possibility of overlap - if( xo_hi < newPoly[0].X() || xo_lo > newPoly[1].X() ) continue; + if( xo_hi < newPoly[0].X() || xo_lo > newPoly[3].X() ) continue; ConvexPolygon oldPoly(xo_lo, xo_hi, yo_lo, yo_hi); try { - ConvexPolygon overlap = intersection(newPoly, oldPoly, Geometry::PolygonIntersection::ORourke); + ConvexPolygon overlap = intersectionByLaszlo(newPoly, oldPoly); overlaps.push_back(BinWithWeight(i,j,overlap.area()/oldPoly.area())); } - catch(std::runtime_error &) + catch(Geometry::NoIntersectionException &) {} } } diff --git a/Code/Mantid/Framework/Geometry/CMakeLists.txt b/Code/Mantid/Framework/Geometry/CMakeLists.txt index bb3be0d8f6d..48cbdfb6fea 100644 --- a/Code/Mantid/Framework/Geometry/CMakeLists.txt +++ b/Code/Mantid/Framework/Geometry/CMakeLists.txt @@ -42,11 +42,9 @@ set ( SRC_FILES src/Math/Algebra.cpp src/Math/BnId.cpp src/Math/ConvexPolygon.cpp -# src/Math/LaszloIntersection.cpp - src/Math/ORourkeIntersection.cpp + src/Math/LaszloIntersection.cpp src/Math/PolyBase.cpp src/Math/PolygonEdge.cpp - src/Math/PolygonIntersection.cpp src/Math/RotCounter.cpp src/Math/Triple.cpp src/Math/Vertex2D.cpp @@ -141,12 +139,10 @@ set ( INC_FILES inc/MantidGeometry/Math/Algebra.h inc/MantidGeometry/Math/BnId.h inc/MantidGeometry/Math/ConvexPolygon.h -# inc/MantidGeometry/Math/LaszloIntersection.h + inc/MantidGeometry/Math/LaszloIntersection.h inc/MantidGeometry/Math/MapSupport.h - inc/MantidGeometry/Math/ORourkeIntersection.h inc/MantidGeometry/Math/PolyBase.h inc/MantidGeometry/Math/PolygonEdge.h - inc/MantidGeometry/Math/PolygonIntersection.h inc/MantidGeometry/Math/RotCounter.h inc/MantidGeometry/Math/Triple.h inc/MantidGeometry/Math/Vertex2D.h @@ -197,7 +193,7 @@ set ( TEST_FILES test/IndexingUtilsTest.h test/InstrumentRayTracerTest.h test/InstrumentTest.h -# test/LaszloIntersectionTest.h + test/LaszloIntersectionTest.h test/LineIntersectVisitTest.h test/LineTest.h test/MDCellTest.h @@ -213,7 +209,6 @@ set ( TEST_FILES test/MaterialTest.h test/MathSupportTest.h test/NearestNeighboursTest.h - test/ORourkeIntersectionTest.h test/ObjCompAssemblyTest.h test/ObjComponentTest.h test/ObjectTest.h @@ -228,7 +223,6 @@ set ( TEST_FILES test/ParametrizedComponentTest.h test/PlaneTest.h test/PolygonEdgeTest.h - test/PolygonIntersectionTest.h test/RectangularDetectorTest.h test/ReflectionConditionTest.h test/RotCounterTest.h diff --git a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/ConvexPolygon.h b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/ConvexPolygon.h index d6a3c324a6e..ce47cb86dc5 100644 --- a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/ConvexPolygon.h +++ b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/ConvexPolygon.h @@ -5,8 +5,8 @@ // Includes //----------------------------------------------------------------------------- #include "MantidGeometry/DllConfig.h" -#include "MantidGeometry/Math/Vertex2DList.h" #include "MantidGeometry/Math/PolygonEdge.h" +#include <vector> namespace Mantid { @@ -17,6 +17,10 @@ namespace Mantid //--------------------------------------------------------------------------- class Vertex2D; + //--------------------------------------------------------------------------- + // Typedefs + //--------------------------------------------------------------------------- + /** An implementation of a convex polygon. It contains a list of vertices that make up a convex polygon and the list is assumed to be ordered in an @@ -52,11 +56,11 @@ namespace Mantid */ class MANTID_GEOMETRY_DLL ConvexPolygon { + public: + public: /// Construct a polygon with a head vertex - ConvexPolygon(Vertex2D & head); - /// Construct a polygon with a collection of vertices - ConvexPolygon(const Vertex2DList & vertices); + ConvexPolygon(Vertex2D * head); /// Construct a rectangle as these will be quite common ConvexPolygon(const double x_lower, const double x_upper, const double y_lower, const double y_upper); @@ -64,42 +68,34 @@ namespace Mantid ConvexPolygon(const ConvexPolygon & rhs); /// Destructor ~ConvexPolygon(); - /// Access the current point - const Kernel::V2D & point() const; - /// Returns an edge on the polygon between the current vertex and the next - PolygonEdge edge() const; + /// Return a pointer to the head vertex + inline const Vertex2D * head() const { return m_head; } /// Index access. const Kernel::V2D& operator[](const size_t index) const; /// Return the number of vertices inline size_t numVertices() const { return m_numVertices; } - /// Advance the current vertex to the next in the chain - const Kernel::V2D& advance(); - + /// Is a point inside this polygon + bool contains(const Kernel::V2D & point) const; /// Compute the area of the polygon using triangulation double area() const; /// Compute the 'determinant' of the points double determinant() const; - /// Compute the orientation - int orientation() const; private: /// Default constructor ConvexPolygon(); - /// Test if the set of points is valid - void validate(const Vertex2DList & points) const; /// Test if a list of vertices is valid - void validate(const Vertex2D & head) const; + void validate(const Vertex2D * head) const; /// Compute the area of a triangle given by 3 points double triangleArea(const Kernel::V2D & a, const Kernel::V2D & b, const Kernel::V2D & c) const; - /// The current vertex - Vertex2D *m_currentVertex; /// The size of the polygon size_t m_numVertices; - /// Head vertex (enables indexing) + /// Head vertex Vertex2D *m_head; }; + /// Print a polygon to a stream MANTID_GEOMETRY_DLL std::ostream & operator<<(std::ostream & os, const ConvexPolygon & polygon); diff --git a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/PolygonIntersection.h b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/LaszloIntersection.h similarity index 64% rename from Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/PolygonIntersection.h rename to Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/LaszloIntersection.h index 3eb5262d56f..e498a552ac3 100644 --- a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/PolygonIntersection.h +++ b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/LaszloIntersection.h @@ -1,20 +1,23 @@ -#ifndef MANTID_GEOMETRY_POLYGONINTERSECTION_H_ -#define MANTID_GEOMETRY_POLYGONINTERSECTION_H_ +#ifndef MANTID_GEOMETRY_LASZLOINTERSECTION_H_ +#define MANTID_GEOMETRY_LASZLOINTERSECTION_H_ //----------------------------------------------------------------------------- // Includes //----------------------------------------------------------------------------- #include "MantidGeometry/DllConfig.h" #include "MantidGeometry/Math/ConvexPolygon.h" +#include <stdexcept> namespace Mantid { namespace Geometry { /** - This header defines various algorithms for defining polygon intersection + This header defines an implementation of the convex polygon intersection method + by Michael Laszlo - @author Martyn Gigg, Tessella plc + @author Martyn Gigg + @date 2011-07-12 Copyright © 2011 ISIS Rutherford Appleton Laboratory & NScD Oak Ridge National Laboratory @@ -37,18 +40,19 @@ namespace Mantid Code Documentation is available at: <http://doxygen.mantidproject.org> */ - /// Define an enum for the intersection algorithm to use - namespace PolygonIntersection + class MANTID_GEOMETRY_DLL NoIntersectionException : public std::runtime_error { - enum Method{Laszlo, ORourke}; + public: + NoIntersectionException() : std::runtime_error("No intersections found that form a complete convex polygon.") + { + } }; - /// Compute the intersection of two polygons using a chosen method + /// Compute the intersection of two ConvexPolygons MANTID_GEOMETRY_DLL - ConvexPolygon intersection(const ConvexPolygon &p, const ConvexPolygon &q, - PolygonIntersection::Method method = PolygonIntersection::ORourke); + ConvexPolygon intersectionByLaszlo(const ConvexPolygon &P, const ConvexPolygon &Q); - } //namespace Geometry -} //namespace Mantid + } // namespace Geometry +} // namespace Mantid -#endif //MANTID_GEOMETRY_POLYGONINTERSECTION_H_ +#endif /* MANTID_GEOMETRY_LASZLOINTERSECTION_H_ */ diff --git a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/PolygonEdge.h b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/PolygonEdge.h index 396e2c80b0a..f1b9203d571 100644 --- a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/PolygonEdge.h +++ b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/PolygonEdge.h @@ -56,6 +56,8 @@ namespace Geometry inline const Kernel::V2D & start() const { return m_start; } /// Access the end point inline const Kernel::V2D & end() const { return m_end; } + /// Return the direction + inline Kernel::V2D direction() const { return m_end - m_start; } /// Create a point a given fraction along this edge Kernel::V2D point(const double fraction) const; @@ -86,7 +88,8 @@ namespace Geometry MANTID_GEOMETRY_DLL PolygonEdge::Orientation orientation(const PolygonEdge & focusEdge, const PolygonEdge & refEdge, double & t); /// Calculate the crossing point of one edge with wrt another MANTID_GEOMETRY_DLL PolygonEdge::Orientation crossingPoint(const PolygonEdge & edgeOne, const PolygonEdge & edgeTwo, Kernel::V2D & crossPoint); - + /// Return if the edges aim at each other + MANTID_GEOMETRY_DLL bool edgeAimsAt(const PolygonEdge &a, const PolygonEdge &b, PointClassification pclass, PolygonEdge::Orientation crossType); } // namespace Geometry } // namespace Mantid diff --git a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/Vertex2D.h b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/Vertex2D.h index 59918d308db..922cf6aef52 100644 --- a/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/Vertex2D.h +++ b/Code/Mantid/Framework/Geometry/inc/MantidGeometry/Math/Vertex2D.h @@ -5,6 +5,7 @@ // Includes //----------------------------------------------------------------------------- #include "MantidGeometry/DllConfig.h" +#include "MantidGeometry/Math/PolygonEdge.h" #include "MantidKernel/V2D.h" namespace Mantid @@ -47,6 +48,8 @@ namespace Mantid Vertex2D(const double x, const double y); /// Constructor with a point Vertex2D(const Kernel::V2D & point); + /// Return the vertex as a point + inline const Kernel::V2D & point() const { return *this; } /// Insert a vertex so that it is next Vertex2D * insert(Vertex2D *vertex); /// Remove this node from the chain @@ -72,6 +75,25 @@ namespace Mantid Vertex2D * m_prev; }; + /** + * A small iterator type structure + */ + class Vertex2DIterator + { + public: + /// Constructor + Vertex2DIterator(const Vertex2D * start) + : m_vertex(start) {} + /// Advance the iterator + void advance() { m_vertex = m_vertex->next(); } + /// Get the point + const Kernel::V2D & point() { return m_vertex->point(); } + /// Get an edge between this and the next + PolygonEdge edge() { return PolygonEdge(*m_vertex, *(m_vertex->next()));} + private: + /// A pointer to the current vertex + Vertex2D const * m_vertex; + }; } // namespace Geometry } // namespace Mantid diff --git a/Code/Mantid/Framework/Geometry/src/Math/ConvexPolygon.cpp b/Code/Mantid/Framework/Geometry/src/Math/ConvexPolygon.cpp index ce77a917829..f70ef2e3cf4 100644 --- a/Code/Mantid/Framework/Geometry/src/Math/ConvexPolygon.cpp +++ b/Code/Mantid/Framework/Geometry/src/Math/ConvexPolygon.cpp @@ -18,43 +18,24 @@ namespace Mantid // Public functions //----------------------------------------------------------------------------- /** - * Constructor with a head vertex. It must be linked to at least 2 others to be considered valid - * @param head :: A reference to the head vertex + * Constructor with a head vertex. + * @param head :: A pointer the head vertex * @throws std::invalid_argument If the vertex list is invalid */ - ConvexPolygon::ConvexPolygon(Vertex2D & head) + ConvexPolygon::ConvexPolygon(Vertex2D * head) { validate(head); - m_currentVertex = m_head = &head; + m_head = head; m_numVertices = 1; - Vertex2D *nextVertex = m_currentVertex->next(); - while( nextVertex != m_head) + Vertex2D *nextVertex = head->next(); + // Count the vertices + while( nextVertex != head ) { ++m_numVertices; nextVertex = nextVertex->next(); } } - /** - * Constructor with a collection of vertices - * @param points :: The points forming the polygon, must be at least 3 - * @throws std::invalid_argument if the vertex list is does not meet - * the convex polygon requirement - */ - ConvexPolygon::ConvexPolygon(const Vertex2DList & points) - : m_currentVertex(NULL), m_numVertices(0) - { - validate(points); - m_numVertices = points.size(); - m_currentVertex = m_head = new Vertex2D(points[0]); - for( size_t i = 1; i < m_numVertices; ++i ) - { - m_currentVertex = m_currentVertex->insert(new Vertex2D(points[i])); - } - // Point it to the first one - m_currentVertex = m_currentVertex->next(); - } - /** * Construct a rectange * @param x_lower :: Lower x coordinate @@ -64,13 +45,13 @@ namespace Mantid */ ConvexPolygon::ConvexPolygon(const double x_lower, const double x_upper, const double y_lower, const double y_upper) - : m_currentVertex(new Vertex2D(x_lower, y_lower)), m_numVertices(4), m_head(m_currentVertex) + : m_numVertices(4), m_head(new Vertex2D(x_lower, y_lower)) { - m_currentVertex = m_currentVertex->insert(new Vertex2D(x_upper, y_lower)); // Bottom right - m_currentVertex = m_currentVertex->insert(new Vertex2D(x_upper, y_upper)); // Top right - m_currentVertex = m_currentVertex->insert(new Vertex2D(x_lower, y_upper)); // Top left - // Point it to the first one - m_currentVertex = m_currentVertex->next(); + // Iterating from the head now produces a clockwise pattern + m_head->insert(new Vertex2D(x_upper, y_lower)); // Bottom right + m_head->insert(new Vertex2D(x_upper, y_upper)); // Top right + m_head->insert(new Vertex2D(x_lower, y_upper)); // Top left + } /** @@ -82,18 +63,16 @@ namespace Mantid if( this != &rhs ) { m_numVertices = rhs.m_numVertices; - Vertex2D *rhsVertex = rhs.m_currentVertex; - m_currentVertex = new Vertex2D(rhsVertex->X(), rhsVertex->Y()); - for(size_t i = 1; i < m_numVertices; ++i) + m_head = new Vertex2D(rhs.m_head->X(), rhs.m_head->Y()); + /// Iterating next() around the other polygon produces the points in clockwise order + /// so we need to go backwards here to ensure that they remain in the same order + Vertex2D *rhsVertex = rhs.m_head->previous(); + // Count the vertices + while( rhsVertex != rhs.m_head ) { - rhsVertex = rhsVertex->next(); - m_currentVertex = m_currentVertex->insert(new Vertex2D(rhsVertex->X(), rhsVertex->Y())); - if( rhsVertex == rhs.m_head ) - { - m_head = m_currentVertex; - } + m_head->insert(new Vertex2D(rhsVertex->X(), rhsVertex->Y())); + rhsVertex = rhsVertex->previous(); } - m_currentVertex = m_currentVertex->next(); } } @@ -103,33 +82,13 @@ namespace Mantid ConvexPolygon::~ConvexPolygon() { // Ensure we delete the vertices - if( m_currentVertex ) + Vertex2D *nextVertex = m_head->next(); + while(nextVertex != m_head) { - Vertex2D *nextVertex = m_currentVertex->next(); - while(m_currentVertex != nextVertex) - { - delete nextVertex->remove(); - nextVertex = m_currentVertex->next(); - } - delete m_currentVertex; + delete nextVertex->remove(); + nextVertex = m_head->next(); } - } - - /** - * Returns the current point - */ - const Kernel::V2D & ConvexPolygon::point() const - { - return *m_currentVertex; - } - - /** - * Returns an edge on the polygon between the current vertex and the next - * @returns A PolygonEdge object - */ - PolygonEdge ConvexPolygon::edge() const - { - return PolygonEdge(*m_currentVertex, *(m_currentVertex->next())); + delete m_head; } /** @@ -155,13 +114,20 @@ namespace Mantid } /** - * Advance the current vertex to the next in the chain and return it - * @returns The next point in the chain + * Is a point inside this polygon + * @param point :: The point to test + * @returns True if the point is inside the polygon */ - const V2D & ConvexPolygon::advance() + bool ConvexPolygon::contains(const Kernel::V2D & point) const { - m_currentVertex = m_currentVertex->next(); - return *m_currentVertex; + Vertex2D *v = m_head; + while( v != m_head ) + { + PolygonEdge edge(v->point(), v->next()->point()); + if( classify(point, edge) == Left ) return false; + v = v->next(); + } + return true; } /** @@ -191,67 +157,40 @@ namespace Mantid // square matrices. We could fool it by putting extra zeroes but this // would increase the workload for no gain - // Loop over all and compute the individual elements - // The ends are handled by the vertex pointer linkage + // Loop over all and compute the individual elements. We assume + // that calling next() on the vertex takes us clockwise within + // the polygon. double lhs(0.0), rhs(0.0); - const Vertex2D *start = m_currentVertex; - Vertex2D *v_i = m_currentVertex; + Vertex2D *v_i = m_head; Vertex2D *v_ip1 = v_i->next(); do { - lhs += v_i->X()*v_ip1->Y(); - rhs += v_ip1->X()*v_i->Y(); + lhs += v_ip1->X()*v_i->Y(); + rhs += v_i->X()*v_ip1->Y(); v_i = v_i->next(); v_ip1 = v_i->next(); } - while( v_i != start ); + while( v_i != m_head ); return lhs - rhs; } - /** - * Compute the 'orientation' of the polygon. The possible values are: - * 1: 2*area > 0.5 - * -1: 2*area < -0.5 else - * 0. - * @returns An integer representing the orientation - */ - int ConvexPolygon::orientation() const - { - const double determinant = this->determinant(); - if( determinant > 0.5 ) return 1; - else if( determinant < -0.5 ) return -1; - else return 0; - } - - /** - * Check this is a valid convex polygon - * @param points :: Validate a list of points - * @throws std::invalid_argument if it is not - */ - void ConvexPolygon::validate(const Vertex2DList & points) const - { - if( points.size() < 3 ) - { - std::ostringstream os; - os << "Expected greater than 2 vertices when constructing a convex polygon, found " - << points.size(); - throw std::invalid_argument(os.str()); - } - } - /** * Check this is a valid polygon * @param head :: A pointer to the head vertex * @throws std::invalid_argument if it is not */ - void ConvexPolygon::validate(const Vertex2D & head) const + void ConvexPolygon::validate(const Vertex2D * head) const { + if( !head ) + { + throw std::invalid_argument("ConvexPolygon::validate - NULL pointer is an invalid head for a convex polygon"); + } // Must have at least two neighbours - if( head.next() == head.previous() ) + if( head->next() == head->previous() ) { std::ostringstream os; - os << "Expected 3 or more vertices when constructing a convex polygon, found "; - if( head.next() == &head ) os << "1"; + os << "ConvexPolygon::validate - Expected 3 or more vertices when constructing a convex polygon, found "; + if( head->next() == head ) os << "1"; else os << "2"; throw std::invalid_argument(os.str()); } @@ -284,11 +223,13 @@ namespace Mantid std::ostream & operator<<(std::ostream & os, const ConvexPolygon & polygon) { os << "ConvexPolygon("; - const size_t nverts(polygon.numVertices()); - for( size_t i = 0; i < nverts; ++i ) + const size_t numVertices(polygon.numVertices()); + Vertex2DIterator pIter(polygon.head()); + for( size_t i = 0 ; i < numVertices; ++i ) { - os << polygon[i]; - if( i < nverts - 1 ) os << ","; + os << pIter.point(); + if( i < numVertices - 1 ) os << ","; + pIter.advance(); } os << ")"; return os; diff --git a/Code/Mantid/Framework/Geometry/src/Math/LaszloIntersection.cpp b/Code/Mantid/Framework/Geometry/src/Math/LaszloIntersection.cpp new file mode 100644 index 00000000000..26826315f21 --- /dev/null +++ b/Code/Mantid/Framework/Geometry/src/Math/LaszloIntersection.cpp @@ -0,0 +1,221 @@ +//----------------------------------------------------------------------------- +// Includes +//----------------------------------------------------------------------------- +#include "MantidGeometry/Math/LaszloIntersection.h" +#include "MantidGeometry/Math/Vertex2D.h" +#include "MantidKernel/Exception.h" +#include <iostream> + +namespace Mantid +{ + namespace Geometry + { + using Kernel::V2D; + + namespace + { + + enum eEdgeIn + { + Unknown, /**< Which edge is inside the other is not known */ + PIsInside, /**< Edge P is inside of edge Q */ + QIsInside /**< Edge Q is inside of edge P */ + }; + + /** + * Advance the current vertex on the given polygon. If the inside flag is set + * record a point of intersection on the given vertex + * @param iter :: An iterator on the current polygon + * @param lastIntersect :: A pointer to the last vertex in the intersection list + * @param inside :: True if the current polygon point is inside another + */ + void advanceVertex(Vertex2DIterator & iter, Vertex2D *&lastIntersect, const bool inside) + { + iter.advance(); + const Kernel::V2D & curPolyPt = iter.point(); + if (inside && (lastIntersect->point() != curPolyPt)) + { + // Add an intersection as the point is inside the polygon + lastIntersect = lastIntersect->insert(new Vertex2D(curPolyPt)); +#ifdef VERBOSE + std::cout << "Advance adds cross pt: (" << curPolyPt.X() << "," << curPolyPt.Y() << ")" << std::endl; +#endif + } + } + + } // Anonymous namespace + + /** + * Compute the polygon that defines the intersection between two + * other concex polygons using the method of chasing edges implemented by Laszlo + * @param P :: A reference to the first polygon + * @param Q :: A reference to the second polygon + * @returns A new polygon defining the region of intersection + */ + ConvexPolygon intersectionByLaszlo(const ConvexPolygon &P, const ConvexPolygon &Q) + { + // The algorithm requires that the polygon with the greatest unsigned area + // be on the "Left" + if( P.determinant() < Q.determinant() ) return intersectionByLaszlo(Q,P); + + Vertex2D *curIntersection(NULL); + V2D iPnt, startPnt; + Vertex2DIterator pIter(P.head()), qIter(Q.head()); + eEdgeIn inflag = Unknown; + int phase = 1; + size_t maxItns = 2 * (P.numVertices() + Q.numVertices()); + for (size_t i = 1; i <= maxItns; ++i) + { +#ifdef VERBOSE + std::cout << "Iteration " << i << " Phase = " << phase << std::endl; +#endif + const PolygonEdge edgeP = pIter.edge(); + const PolygonEdge edgeQ = qIter.edge(); + PointClassification pclass = classify(edgeP.end(), edgeQ); +#ifdef VERBOSE + std::cout << "Class P Pt" << std::endl; + std::cout << "Class Pt: (" << edgeP.end().X() << "," << edgeP.end().Y() << ")"; + std::cout << std::endl; + std::cout << "Edge Orig Pt (" << edgeQ.start().X() << "," << edgeQ.start().Y() << ")"; + std::cout << std::endl; + std::cout << "Edge Dest Pt (" << edgeQ.end().X() << "," << edgeQ.end().Y() << ")"; + std::cout << std::endl; + std::cout << "P pt class: " << pclass << std::endl; +#endif + PointClassification qclass = classify(edgeQ.end(), edgeP); +#ifdef VERBOSE + std::cout << "Class Q Pt" << std::endl; + std::cout << "Class Pt: (" << edgeQ.end().X() << "," << edgeQ.end().Y() << ")"; + std::cout << std::endl; + std::cout << "Edge Orig Pt (" << edgeP.start().X() << "," << edgeP.start().Y() << ")"; + std::cout << std::endl; + std::cout << "Edge Dest Pt (" << edgeP.end().X() << "," << edgeP.end().Y() << ")"; + std::cout << std::endl; + std::cout << "Q pt class: " << qclass << std::endl; +#endif + PolygonEdge::Orientation crossType = crossingPoint(edgeP, edgeQ, iPnt); +#ifdef VERBOSE + std::cout << "PQ Orient: " << crossType << std::endl; +#endif + if (crossType == PolygonEdge::SkewCross) + { + if (phase == 1) + { + phase = 2; +#ifdef VERBOSE + std::cout << "Found a crossing pt: (" << iPnt.X() << ","; + std::cout << iPnt.Y() << ")" << std::endl; +#endif + + curIntersection = new Vertex2D(iPnt); + startPnt = iPnt; + } + else if (iPnt != *curIntersection) + { +#ifdef VERBOSE + std::cout << "Found a crossing pt: (" << iPnt.X() << ","; + std::cout << iPnt.Y() << ")" << std::endl; +#endif + if (iPnt != startPnt) + { + curIntersection = curIntersection->insert(new Vertex2D(iPnt)); + } + else // Back to the start, we're done + { + try + { + // Make the head vertex of the polygon the first one we found + return ConvexPolygon(curIntersection->next()); + } + catch(std::invalid_argument&) + { + throw NoIntersectionException(); + } + } + } + if (pclass == Right) + { + inflag = PIsInside; + } + else if (qclass == Right) + { + inflag = QIsInside; + } + else + { + inflag = Unknown; + } + } + else if ( (crossType == PolygonEdge::Collinear) && + (pclass != Behind) && + (qclass != Behind) ) + { + inflag = Unknown; + } +#ifdef VERBOSE + std::cout << "Current in flag: " << inflag << std::endl; +#endif + + bool pAIMSq = edgeAimsAt(edgeP, edgeQ, pclass, crossType); + bool qAIMSp = edgeAimsAt(edgeQ, edgeP, qclass, crossType); +#ifdef VERBOSE + std::cout << "P aims at Q:" << pAIMSq << std::endl; + std::cout << "Q aims at P:" << qAIMSp << std::endl; +#endif + if (pAIMSq && qAIMSp) + { + if( (inflag == QIsInside) || ((inflag == Unknown) && (pclass == Left)) ) + { +#ifdef VERBOSE + std::cout << "Move edge on P" << std::endl; +#endif + advanceVertex(pIter, curIntersection, false); + } + else + { +#ifdef VERBOSE + std::cout << "Move edge on Q" << std::endl; +#endif + advanceVertex(qIter, curIntersection, false); + } + } + else if (pAIMSq) + { +#ifdef VERBOSE + std::cout << "Move edge on P" << std::endl; +#endif + advanceVertex(pIter, curIntersection, inflag==PIsInside); + } + else if (qAIMSp) + { +#ifdef VERBOSE + std::cout << "Move edge on Q" << std::endl; +#endif + advanceVertex(qIter, curIntersection, inflag==QIsInside); + } + else + { + if ((inflag == QIsInside) || ((inflag == Unknown) && (pclass == Left))) + { +#ifdef VERBOSE + std::cout << "Move edge on P" << std::endl; +#endif + advanceVertex(pIter, curIntersection, false); + } + else + { +#ifdef VERBOSE + std::cout << "Move edge on Q" << std::endl; +#endif + advanceVertex(qIter, curIntersection, false); + } + } + } // end-for + + // Reaching here means no intersections have been found + throw NoIntersectionException(); + } + + } // namespace Geometry +} // namespace Mantid + diff --git a/Code/Mantid/Framework/Geometry/src/Math/PolygonEdge.cpp b/Code/Mantid/Framework/Geometry/src/Math/PolygonEdge.cpp index ae86295ffbb..13914d20f97 100644 --- a/Code/Mantid/Framework/Geometry/src/Math/PolygonEdge.cpp +++ b/Code/Mantid/Framework/Geometry/src/Math/PolygonEdge.cpp @@ -161,6 +161,36 @@ namespace Mantid } } + /** + * Return if the edges aim at each other + * @param a :: First edge + * @param b :: Second edge + * @param pclass :: The point classification of a point on edge a + */ + bool edgeAimsAt(const PolygonEdge &a, const PolygonEdge &b, PointClassification aclass, PolygonEdge::Orientation crossType) + { + V2D va = a.direction(); + V2D vb = b.direction(); + if (crossType != PolygonEdge::Collinear) + { + double ca = va.X() * vb.Y(); + double cb = vb.X() * va.Y(); + if( Kernel::gtEquals(ca, cb) ) + { + return (aclass != Right); + } + else + { + return (aclass != Left); + } + } + else + { + return (aclass != Beyond); + } + } + + } // namespace Geometry } // namespace Mantid diff --git a/Code/Mantid/Framework/Geometry/src/Math/PolygonIntersection.cpp b/Code/Mantid/Framework/Geometry/src/Math/PolygonIntersection.cpp deleted file mode 100644 index 6b05105d23e..00000000000 --- a/Code/Mantid/Framework/Geometry/src/Math/PolygonIntersection.cpp +++ /dev/null @@ -1,33 +0,0 @@ -//----------------------------------------------------------------------------- -// Includes -//----------------------------------------------------------------------------- -#include "MantidGeometry/Math/PolygonIntersection.h" -#include "MantidGeometry/Math/ORourkeIntersection.h" -//#include "MantidGeometry/Math/LaszloIntersection.h" -#include "MantidKernel/Exception.h" - -namespace Mantid -{ - namespace Geometry - { - /** - * Compute the intersection of two polygons via the chosen method - * @param p :: The first polygon - * @param q :: The second polygon - * @param method :: The chosen method (Default = Laszlo) - * @returns The overlap of the two polygons as a ConvexPolygon object - */ - ConvexPolygon intersection(const ConvexPolygon &p, const ConvexPolygon &q, - PolygonIntersection::Method method) - { - switch(method) - { - //case PolygonIntersection::Laszlo: return intersectionByLaszlo(p,q); - case PolygonIntersection::ORourke: return intersectionByORourke(p,q); - default: throw std::invalid_argument("Unknown type of polygon intersection algorithm requested."); - }; - } - - - } //namespace Geometry -} //namespace Mantid diff --git a/Code/Mantid/Framework/Geometry/test/ConvexPolygonTest.h b/Code/Mantid/Framework/Geometry/test/ConvexPolygonTest.h index e8c208eca41..98fbaa45530 100644 --- a/Code/Mantid/Framework/Geometry/test/ConvexPolygonTest.h +++ b/Code/Mantid/Framework/Geometry/test/ConvexPolygonTest.h @@ -11,7 +11,6 @@ using Mantid::Kernel::V2D; using Mantid::Geometry::ConvexPolygon; using Mantid::Geometry::Vertex2D; -using Mantid::Geometry::Vertex2DList; using Mantid::Geometry::PolygonEdge; class ConvexPolygonTest : public CxxTest::TestSuite @@ -23,30 +22,20 @@ public: TS_ASSERT_THROWS_NOTHING(makeEquilateralTriangle()); } - void test_Building_With_An_Too_Small_A_Set_Throws_Invalid_Arg() + void test_Building_With_An_Isolated_Vertex_Throws_Invalid_Arg() { - Vertex2DList vertices; - doExceptionCheck<std::invalid_argument>(vertices); - // Single vertex - vertices.insert(V2D()); - doExceptionCheck<std::invalid_argument>(vertices); - // Line - vertices.insert(V2D(1.,1.)); - doExceptionCheck<std::invalid_argument>(vertices); + Vertex2D *vertex = new Vertex2D(); + doExceptionCheck<std::invalid_argument>(vertex); + delete vertex; } - void test_Building_With_Non_Unique_Vertices_Only_Keeps_Unique() + void test_Building_With_A_Line_Throws_Invalid_Arg() { - Vertex2DList vertices; - vertices.insert(V2D()); - vertices.insert(V2D(2.0,0.0)); - vertices.insert(V2D(1.0,std::sqrt(3.0))); - vertices.insert(V2D(2.0,0.0)); - ConvexPolygon poly(vertices); - TS_ASSERT_EQUALS(poly.numVertices(), 3); - TS_ASSERT_EQUALS(poly[0], V2D()); - TS_ASSERT_EQUALS(poly[1], V2D(2.0,0.0)); - TS_ASSERT_EQUALS(poly[2], V2D(1.0,std::sqrt(3.0))); + Vertex2D *vertex = new Vertex2D(); + Vertex2D *second = vertex->insert(new Vertex2D(1.0,1.0)); + doExceptionCheck<std::invalid_argument>(vertex); + delete vertex; + delete second; } void test_Building_With_Head_Vertex_With_Two_Other_Points_Doesnt_Not_Throw() @@ -55,7 +44,7 @@ public: origin = origin->insert(new Vertex2D(1.0,1.0)); origin = origin->insert(new Vertex2D(0.0,1.0)); - TS_ASSERT_THROWS_NOTHING(new ConvexPolygon(*origin)); + TS_ASSERT_THROWS_NOTHING(new ConvexPolygon(origin)); } void test_Building_With_Head_Vertex_Gives_Correct_Number_Of_Vertices() @@ -63,17 +52,17 @@ public: Vertex2D *origin = new Vertex2D; origin = origin->insert(new Vertex2D(1.0,1.0)); origin = origin->insert(new Vertex2D(0.0,1.0)); - ConvexPolygon poly(*origin); + ConvexPolygon poly(origin); TS_ASSERT_EQUALS(poly.numVertices(), 3); } void test_Buildling_With_Head_Vertex_With_Less_Than_Two_Other_Points_Throws() { Vertex2D *origin = new Vertex2D; - TS_ASSERT_THROWS(new ConvexPolygon(*origin), std::invalid_argument); + TS_ASSERT_THROWS(new ConvexPolygon(origin), std::invalid_argument); Vertex2D *second = new Vertex2D(1.0,1.0); origin->insert(second); - TS_ASSERT_THROWS(new ConvexPolygon(*origin), std::invalid_argument); + TS_ASSERT_THROWS(new ConvexPolygon(origin), std::invalid_argument); delete second; delete origin; } @@ -81,46 +70,30 @@ public: void test_Copying_Preserves_Polygon() { // Returns by value but with some optimizations may not copy the object - ConvexPolygon triangle = makeEquilateralTriangle(); - ConvexPolygon copy(triangle); - - TS_ASSERT_EQUALS(copy.numVertices(), 3); - TS_ASSERT_EQUALS(copy.point(), V2D()); - copy.advance(); - TS_ASSERT_EQUALS(copy.point(), V2D(2.0, 0.0)); - copy.advance(); - TS_ASSERT_EQUALS(copy.point(), V2D(1.0,std::sqrt(3.0))); - } + ConvexPolygon rect = makeRectangle(); + TS_ASSERT_EQUALS(rect.numVertices(), 4); + TS_ASSERT_EQUALS(rect[0], V2D()); + TS_ASSERT_EQUALS(rect[1], V2D(0.0, 1.0)); + TS_ASSERT_EQUALS(rect[2], V2D(2.0, 1.0)); + TS_ASSERT_EQUALS(rect[3], V2D(2.0, 0.0)); - void test_Point_Access_Returns_Current_Point() - { - ConvexPolygon poly = makeRectangle(); - TS_ASSERT_EQUALS(poly.point(), V2D()); - poly.advance(); - TS_ASSERT_EQUALS(poly.point(), V2D(2.0, 0.0)); - poly.advance(); - TS_ASSERT_EQUALS(poly.point(), V2D(2.0, 1.0)); - poly.advance(); - TS_ASSERT_EQUALS(poly.point(), V2D(0.0, 1.0)); + ConvexPolygon copy(rect); + TS_ASSERT_EQUALS(copy[0], V2D()); + TS_ASSERT_EQUALS(copy[1], V2D(0.0, 1.0)); + TS_ASSERT_EQUALS(copy[2], V2D(2.0, 1.0)); + TS_ASSERT_EQUALS(copy[3], V2D(2.0, 0.0)); } - void test_Edge_Gives_Back_Current_Edge() + void test_Head_Returns_Correct_Vertex() { ConvexPolygon poly = makeRectangle(); - PolygonEdge edge1 = poly.edge(); - TS_ASSERT_EQUALS(edge1.start(), V2D()); - TS_ASSERT_EQUALS(edge1.end(), V2D(2.0,0.0)); - poly.advance(); - PolygonEdge edge2 = poly.edge(); - TS_ASSERT_EQUALS(edge2.start(), V2D(2.0,0.0)); - TS_ASSERT_EQUALS(edge2.end(), V2D(2.0,1.0)); - + TS_ASSERT_EQUALS(poly.head()->point(), V2D(0.0, 0.0)); } void test_Index_Access_Returns_Correct_Object_For_Valid_Index() { ConvexPolygon triangle = makeEquilateralTriangle(); - const V2D & apex = triangle[2]; + const V2D & apex = triangle[1]; TS_ASSERT_EQUALS(apex, V2D(1.0, std::sqrt(3.0))); } @@ -132,12 +105,12 @@ public: TS_ASSERT_THROWS(triangle[-1], IndexError); } - void test_Advance_Increments_The_Current_Vertex_And_Returns_The_Next() + void test_Point_Inside_Polygon_Returns_True() { - ConvexPolygon triangle = makeEquilateralTriangle(); - V2D nextPoint; - TS_ASSERT_THROWS_NOTHING(nextPoint = triangle.advance()); - TS_ASSERT_EQUALS(nextPoint, V2D(2.0,0.0)); + ConvexPolygon poly = makeRectangle(); + TS_ASSERT(poly.contains(V2D(1.0, 0.25))); + TS_ASSERT(poly.contains(V2D(1.0, 0.0))); + TS_ASSERT(poly.contains(poly.head()->point())); } void test_The_Determinant_For_A_Triangle() @@ -169,11 +142,10 @@ private: /// Side length 2 ConvexPolygon makeEquilateralTriangle() { - Vertex2DList vertices; - vertices.insert(V2D()); - vertices.insert(V2D(2.0,0.0)); - vertices.insert(V2D(1.0,std::sqrt(3.0))); - return ConvexPolygon(vertices); + Vertex2D *head = new Vertex2D(); + head->insert(new Vertex2D(2.0,0.0)); + head->insert(new Vertex2D(1.0,std::sqrt(3.0))); + return ConvexPolygon(head); } /// Short side 1, long side 2 @@ -185,22 +157,21 @@ private: /// Short side 2-1-2-1 ConvexPolygon makeParallelogram() { - Vertex2DList vertices; - vertices.insert(V2D()); - vertices.insert(V2D(2.0, 0.0)); - vertices.insert(V2D(2.0 + 0.5*std::sqrt(2.0), 0.5*std::sqrt(2.0))); - vertices.insert(V2D(0.5*std::sqrt(2.0), 0.5*std::sqrt(2.0))); - return ConvexPolygon(vertices); + Vertex2D *head = new Vertex2D(); + head->insert(new Vertex2D(2.0,0.0)); + head->insert(new Vertex2D(2.0 + 0.5*std::sqrt(2.0), 0.5*std::sqrt(2.0))); + head->insert(new Vertex2D(0.5*std::sqrt(2.0), 0.5*std::sqrt(2.0))); + return ConvexPolygon(head); } // If a class has no accessible default constructor we cannot use // TS_ASSERT_THROWS() template <typename ExceptionType> - void doExceptionCheck(const Vertex2DList & vertices) + void doExceptionCheck(Vertex2D * vertex) { try { - ConvexPolygon p(vertices); + ConvexPolygon p(vertex); } catch(ExceptionType &) { diff --git a/Code/Mantid/Framework/Geometry/test/LaszloIntersectionTest.h b/Code/Mantid/Framework/Geometry/test/LaszloIntersectionTest.h new file mode 100644 index 00000000000..2a1d9924c1f --- /dev/null +++ b/Code/Mantid/Framework/Geometry/test/LaszloIntersectionTest.h @@ -0,0 +1,175 @@ +#ifndef MANTID_GEOMETRY_LASZLOINTERSECTIONTEST_H_ +#define MANTID_GEOMETRY_LASZLOINTERSECTIONTEST_H_ + +//----------------------------------------------------------------------------- +// Includes +//----------------------------------------------------------------------------- +#include "MantidGeometry/Math/LaszloIntersection.h" +#include "MantidGeometry/Math/ConvexPolygon.h" +#include "MantidGeometry/Math/Vertex2D.h" +#include <cfloat> +#include <cxxtest/TestSuite.h> + +using Mantid::Kernel::V2D; +using Mantid::Geometry::Vertex2D; +using Mantid::Geometry::ConvexPolygon; + +class LaszloIntersectionTest : public CxxTest::TestSuite +{ +public: + + void test_Intersection_Of_Axis_Aligned_Squares() + { + // Define two squares that partially overlap + ConvexPolygon squareOne(0.0, 2.0, 0.0, 2.0); //2x2, bottom left-hand corner at origin + ConvexPolygon squareTwo(1.0, 3.0, 1.0, 3.0); //2x2, bottom left-hand corner at centre of first + + ConvexPolygon overlap = intersectionByLaszlo(squareOne,squareTwo); + TS_ASSERT_EQUALS(overlap.numVertices(), 4); + // Are they correct + TS_ASSERT_EQUALS(overlap[0], V2D(1,2)); + TS_ASSERT_EQUALS(overlap[1], V2D(2,2)); + TS_ASSERT_EQUALS(overlap[2], V2D(2,1)); + TS_ASSERT_EQUALS(overlap[3], V2D(1,1)); + + // Symmetry + ConvexPolygon overlapRev = intersectionByLaszlo(squareTwo,squareOne); + TS_ASSERT_EQUALS(overlapRev.numVertices(), 4); + // Are they correct + TS_ASSERT_EQUALS(overlapRev[0], V2D(1,2)); + TS_ASSERT_EQUALS(overlapRev[1], V2D(2,2)); + TS_ASSERT_EQUALS(overlapRev[2], V2D(2,1)); + TS_ASSERT_EQUALS(overlapRev[3], V2D(1,1)); + + } + + void test_House() + { + Vertex2D *head = new Vertex2D(); + head->insert(new Vertex2D(200,0)); + head->insert(new Vertex2D(200,100)); + head->insert(new Vertex2D(100,200)); + head->insert(new Vertex2D(0,100)); + ConvexPolygon house(head); + + head = new Vertex2D(100,100); + head->insert(new Vertex2D(300,100)); + head->insert(new Vertex2D(300,200)); + head->insert(new Vertex2D(100,200)); + ConvexPolygon rectangle(head); + + ConvexPolygon overlap = intersectionByLaszlo(house,rectangle); + TS_ASSERT_EQUALS(overlap.numVertices(), 3); + TS_ASSERT_EQUALS(overlap[0], V2D(100,200)); + TS_ASSERT_EQUALS(overlap[1], V2D(200,100)); + TS_ASSERT_EQUALS(overlap[2], V2D(100,100)); + + ConvexPolygon overlapRev = intersectionByLaszlo(rectangle, house); + TS_ASSERT_EQUALS(overlapRev.numVertices(), 3); + TS_ASSERT_EQUALS(overlapRev[0], V2D(100,200)); + TS_ASSERT_EQUALS(overlapRev[1], V2D(200,100)); + TS_ASSERT_EQUALS(overlapRev[2], V2D(100,100)); + } + + void test_Intersection_Of_Parallelogram_And_Square() + { + Vertex2D *head = new Vertex2D(100,50); + head->insert(new Vertex2D(175,50)); + head->insert(new Vertex2D(175,125)); + head->insert(new Vertex2D(100,125)); + ConvexPolygon square(head); + + head = new Vertex2D(); + head->insert(new Vertex2D(200,0)); + head->insert(new Vertex2D(300,100)); + head->insert(new Vertex2D(100,100)); + ConvexPolygon parallelogram(head); + + ConvexPolygon overlap = intersectionByLaszlo(square,parallelogram); + TS_ASSERT_EQUALS(overlap.numVertices(), 4); + // Are they correct + TS_ASSERT_EQUALS(overlap[0], V2D(100,100)); + TS_ASSERT_EQUALS(overlap[1], V2D(175,100)); + TS_ASSERT_EQUALS(overlap[2], V2D(175,50)); + TS_ASSERT_EQUALS(overlap[3], V2D(100,50)); + + //Symmetry + ConvexPolygon overlapRev = intersectionByLaszlo(parallelogram, square); + TS_ASSERT_EQUALS(overlapRev.numVertices(), 4); + // Are they correct + TS_ASSERT_EQUALS(overlapRev[0], V2D(100,100)); + TS_ASSERT_EQUALS(overlapRev[1], V2D(175,100)); + TS_ASSERT_EQUALS(overlapRev[2], V2D(175,50)); + TS_ASSERT_EQUALS(overlapRev[3], V2D(100,50)); + } + + void test_Intersection_With_Self() + { + ConvexPolygon squareOne(0.0, 2.0, 0.0, 2.0); //2x2, bottom left-hand corner at origin + ConvexPolygon overlap = intersectionByLaszlo(squareOne, squareOne); + TS_ASSERT_EQUALS(overlap.numVertices(), 4); + + // Are they correct + TS_ASSERT_EQUALS(overlap[0], V2D(0.0, 2.0)); + TS_ASSERT_EQUALS(overlap[1], V2D(2.0, 2.0)); + TS_ASSERT_EQUALS(overlap[2], V2D(2.0, 0.0)); + TS_ASSERT_EQUALS(overlap[3], V2D()); + } + + void test_Shapes_Sharing_A_Line_Throws() + { + Vertex2D *plist = new Vertex2D(-3.0, -3.0); + plist->insert(new Vertex2D(-0.5, -3.0)); + plist->insert(new Vertex2D(0.5, -1.0)); + plist->insert(new Vertex2D(-2.0, -1.0)); + ConvexPolygon parallelogram(plist); + + Vertex2D *s2list = new Vertex2D(1.0, -1.0); + s2list->insert(new Vertex2D(1.0, 3.0)); + s2list->insert(new Vertex2D(-4.0, 3.0)); + s2list->insert(new Vertex2D(-4.0, -1.0)); + ConvexPolygon rect2(s2list); + + /** + * The overlap here is a line-segment [-3,-1]->[0.5,-1] which is not + * valid polygon so this should throw + */ + TS_ASSERT_THROWS(intersectionByLaszlo(rect2, parallelogram), Mantid::Geometry::NoIntersectionException); + TS_ASSERT_THROWS(intersectionByLaszlo(parallelogram, rect2), Mantid::Geometry::NoIntersectionException); + + + } + +}; + +//------------------------------------------------------------------------ +// Performance Tests +//------------------------------------------------------------------------ + +class LaszloIntersectionTestPerformance : public CxxTest::TestSuite +{ +public: + void test_Intersection_Of_Large_Number() + { + const size_t niters(100000); + for( size_t i = 0; i < niters; ++i ) + { + // These are created each loop iteration to simulate a more real-life case + // of constructing polygons inside a loop and then testing their intersection + ConvexPolygon squareOne(0.0, 2.0, 0.0, 2.0); //2x2, bottom left-hand corner at origin + ConvexPolygon squareTwo(1.0, 3.0, 1.0, 3.0); //2x2, bottom left-hand corner at centre of first + try + { + intersectionByLaszlo(squareOne, squareTwo); + } + catch(Mantid::Geometry::NoIntersectionException&) + { + } + } + + } +}; + + +#endif /* MANTID_GEOMETRY_LASZLOINTERSECTIONTEST_H_ */ + diff --git a/Code/Mantid/Framework/Geometry/test/PolygonIntersectionTest.h b/Code/Mantid/Framework/Geometry/test/PolygonIntersectionTest.h deleted file mode 100644 index 2c379205e8a..00000000000 --- a/Code/Mantid/Framework/Geometry/test/PolygonIntersectionTest.h +++ /dev/null @@ -1,106 +0,0 @@ -#ifndef MANTID_GEOMETRY_POLYGONINTERSECTIONTEST_H_ -#define MANTID_GEOMETRY_POLYGONINTERSECTIONTEST_H_ - -//----------------------------------------------------------------------------- -// Includes -//----------------------------------------------------------------------------- -#include "MantidGeometry/Math/PolygonIntersection.h" -#include <cfloat> -#include <cxxtest/TestSuite.h> - -using Mantid::Kernel::V2D; -using Mantid::Geometry::Vertex2DList; -using Mantid::Geometry::ConvexPolygon; - -class PolygonIntersectionTest : public CxxTest::TestSuite -{ -public: - - void test_Stub() - { - } - - void xtest_Squares_With_Side_Length_Less_Than_One() - { - Vertex2DList vertices; - vertices.insert(V2D()); - vertices.insert(V2D(1,0)); - vertices.insert(V2D(1,1)); - vertices.insert(V2D(0,1)); - ConvexPolygon squareOne(vertices); - - vertices = Vertex2DList(); - vertices.insert(V2D(0,0.1)); - vertices.insert(V2D(0.2,0.1)); - vertices.insert(V2D(0.2,0.2)); - vertices.insert(V2D(0,0.2)); - ConvexPolygon squareTwo(vertices); - - ConvexPolygon overlap = intersection(squareOne,squareTwo); - TS_ASSERT_EQUALS(overlap.numVertices(), 4); - // Are they correct - TS_ASSERT_EQUALS(overlap[0], V2D(0,0.2)); - TS_ASSERT_EQUALS(overlap[1], V2D(0,0.1)); - TS_ASSERT_EQUALS(overlap[2], V2D(0.2,0.1)); - TS_ASSERT_EQUALS(overlap[3], V2D(0.2,0.2)); - TS_ASSERT_DELTA(overlap.area(), squareTwo.area(), DBL_EPSILON); - } - - void xtest_Squares_With_One_Larger() - { - Vertex2DList vertices; - vertices.insert(V2D(5.0, -0.5)); - vertices.insert(V2D(5.2, -0.5)); - vertices.insert(V2D(5.2, 0.5)); - vertices.insert(V2D(5.0, 0.5)); - ConvexPolygon squareOne(vertices); - - vertices = Vertex2DList(); - vertices.insert(V2D(5.0, -0.5)); - vertices.insert(V2D(5.1, -0.5)); - vertices.insert(V2D(5.1, 0.5)); - vertices.insert(V2D(5.0, 0.5)); - ConvexPolygon squareTwo(vertices); - - intersection(squareOne, squareTwo); - - ConvexPolygon overlap = intersection(squareOne, squareTwo); - TS_ASSERT_EQUALS(overlap.numVertices(), 4); - // Are they correct - TS_ASSERT_DELTA(overlap.area(), squareTwo.area(), DBL_EPSILON); - } -}; - - -//------------------------------------------------------------------------ -// Performance Tests -//------------------------------------------------------------------------ - -class PolygonIntersectionTestPerformance : public CxxTest::TestSuite -{ -public: - void test_Intersection_Of_Large_Number() - { - const size_t niters(100000); - for( size_t i = 0; i < niters; ++i ) - { - // These are created each loop iteration to simulate a more real-life case - // of constructing polygons inside a loop and then testing their intersection - ConvexPolygon squareOne(0.0, 2.0, 0.0, 2.0); //2x2, bottom left-hand corner at origin - ConvexPolygon squareTwo(1.0, 3.0, 1.0, 3.0); //2x2, bottom left-hand corner at centre of first - try - { - intersection(squareOne, squareTwo); - } - catch(std::runtime_error&) - { - } - } - - } - - -}; - -#endif //MANTID_GEOMETRY_POLYGONINTERSECTIONTEST_H_ - diff --git a/Code/Mantid/Framework/Geometry/test/Vertex2DTest.h b/Code/Mantid/Framework/Geometry/test/Vertex2DTest.h index 5d0c3996b93..27b64fb9896 100644 --- a/Code/Mantid/Framework/Geometry/test/Vertex2DTest.h +++ b/Code/Mantid/Framework/Geometry/test/Vertex2DTest.h @@ -5,6 +5,7 @@ #include "MantidGeometry/Math/Vertex2D.h" using Mantid::Geometry::Vertex2D; +using Mantid::Geometry::Vertex2DIterator; using Mantid::Kernel::V2D; class Vertex2DTest : public CxxTest::TestSuite @@ -47,6 +48,15 @@ public: TS_ASSERT_EQUALS(vertexPt.previous(), &vertexPt); } + void test_Vertex_As_Pt_Returns_Correct_Value() + { + Vertex2D vertex(5.1, 10.9); + TS_ASSERT_EQUALS(vertex.point(), V2D(5.1,10.9)); + Vertex2D *vertex2 = new Vertex2D(5.1, 10.9); + TS_ASSERT_EQUALS(vertex2->point(), V2D(5.1,10.9)); + delete vertex2; + } + void test_Insert_Yields_Next_As_Inserted_Vertex() { makeThreeVertexChain(true, false); @@ -57,49 +67,66 @@ public: makeThreeVertexChain(false, true); } + void test_Iteration_Advances_Correctly() + { + Vertex2D * start = makeThreeVertexChain(false,false); + Vertex2DIterator pIter(start); + TS_ASSERT_EQUALS(pIter.point(), V2D()); + pIter.advance(); + TS_ASSERT_EQUALS(pIter.point(), V2D(0.0,1.0)); + pIter.advance(); + TS_ASSERT_EQUALS(pIter.point(), V2D(1.0,1.0)); + pIter.advance(); //Back to the start + TS_ASSERT_EQUALS(pIter.point(), V2D()); + } + private: - void makeThreeVertexChain(const bool doInsertTests, const bool doRemoveTests) + Vertex2D * makeThreeVertexChain(const bool doInsertTests, const bool doRemoveTests) { - Vertex2D origin; - Vertex2D two(0.0,1.0); - Vertex2D *vertexTwo = origin.insert(&two); + Vertex2D *origin = new Vertex2D; + Vertex2D *two = new Vertex2D(0.0,1.0); + Vertex2D *vertexTwo = origin->insert(two); if( doInsertTests ) { - TS_ASSERT_EQUALS(vertexTwo, &two); + TS_ASSERT_EQUALS(vertexTwo, two); - TS_ASSERT_EQUALS(origin.next(), &two); - TS_ASSERT_EQUALS(origin.previous(), &two); - TS_ASSERT_EQUALS(vertexTwo->previous(), &origin); - TS_ASSERT_EQUALS(vertexTwo->next(), &origin); + TS_ASSERT_EQUALS(origin->next(), two); + TS_ASSERT_EQUALS(origin->previous(), two); + TS_ASSERT_EQUALS(vertexTwo->previous(), origin); + TS_ASSERT_EQUALS(vertexTwo->next(), origin); } //and a third - Vertex2D third(1.0, 1.0); - Vertex2D *vertexThree = two.insert(&third); + Vertex2D *third = new Vertex2D(1.0, 1.0); + Vertex2D *vertexThree = two->insert(third); if( doInsertTests ) { - TS_ASSERT_EQUALS(vertexThree, &third); - - TS_ASSERT_EQUALS(origin.next(), &two); - TS_ASSERT_EQUALS(origin.previous(), &third); - TS_ASSERT_EQUALS(vertexTwo->previous(), &origin); - TS_ASSERT_EQUALS(vertexTwo->next(), &third); - TS_ASSERT_EQUALS(vertexThree->previous(), &two); - TS_ASSERT_EQUALS(vertexThree->next(), &origin); + TS_ASSERT_EQUALS(vertexThree, third); + + TS_ASSERT_EQUALS(origin->next(), two); + TS_ASSERT_EQUALS(origin->previous(), third); + TS_ASSERT_EQUALS(vertexTwo->previous(), origin); + TS_ASSERT_EQUALS(vertexTwo->next(), third); + TS_ASSERT_EQUALS(vertexThree->previous(), two); + TS_ASSERT_EQUALS(vertexThree->next(), origin); } if( doRemoveTests ) { Vertex2D *removedOne = vertexThree->remove(); TS_ASSERT_EQUALS(removedOne, vertexThree); - TS_ASSERT_EQUALS(origin.next(), &two); - TS_ASSERT_EQUALS(origin.previous(), &two); - TS_ASSERT_EQUALS(vertexTwo->previous(), &origin); - TS_ASSERT_EQUALS(vertexTwo->next(), &origin); + TS_ASSERT_EQUALS(origin->next(), two); + TS_ASSERT_EQUALS(origin->previous(), two); + TS_ASSERT_EQUALS(vertexTwo->previous(), origin); + TS_ASSERT_EQUALS(vertexTwo->next(), origin); Vertex2D *removedTwo = vertexTwo->remove(); - TS_ASSERT( removedTwo); - TS_ASSERT_EQUALS(origin.next(), &origin); - TS_ASSERT_EQUALS(origin.previous(), &origin); + TS_ASSERT_EQUALS(origin->next(), origin); + TS_ASSERT_EQUALS(origin->previous(), origin); + return NULL; + } + else + { + return origin; } } diff --git a/Code/Mantid/Framework/Kernel/inc/MantidKernel/Vertex2D.h b/Code/Mantid/Framework/Kernel/inc/MantidKernel/Vertex2D.h deleted file mode 100644 index 7edcc10bdf6..00000000000 --- a/Code/Mantid/Framework/Kernel/inc/MantidKernel/Vertex2D.h +++ /dev/null @@ -1,38 +0,0 @@ -#ifndef MANTID_KERNEL_VERTEX2D_H_ -#define MANTID_KERNEL_VERTEX2D_H_ - -//----------------------------------------------------------------------------- -// Includes -//----------------------------------------------------------------------------- -#include "MantidKernel/V2D.h" - -namespace Mantid -{ - namespace Kernel - { - /** - Implements a vertex in a 2D plane embedded in a 3D space - - @author Martyn Gigg, Tessella plc - - Copyright © 2011 ISIS Rutherford Appleton Laboratory & NScD Oak Ridge National Laboratory - - This file is part of Mantid. - - Mantid 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 3 of the License, or - (at your option) any later version. - - Mantid 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, see <http://www.gnu.org/licenses/>. - - File change history is stored at: <https://svn.mantidproject.org/mantid/trunk/Code/Mantid>. - Code Documentation is available at: <http://doxygen.mantidproject.org> - */ - typedef V2D Vertex2D; -- GitLab