Home  · Classes  · Annotated Classes  · Modules  · Members  · Namespaces  · Related Pages

ConvexHull2D.h (Maintainer: Marc Sturm)

Go to the documentation of this file.
00001 // -*- Mode: C++; tab-width: 2; -*-
00002 // vi: set ts=2:
00003 //
00004 // --------------------------------------------------------------------------
00005 //                   OpenMS Mass Spectrometry Framework 
00006 // --------------------------------------------------------------------------
00007 //  Copyright (C) 2003-2008 -- Oliver Kohlbacher, Knut Reinert
00008 //
00009 //  This library is free software; you can redistribute it and/or
00010 //  modify it under the terms of the GNU Lesser General Public
00011 //  License as published by the Free Software Foundation; either
00012 //  version 2.1 of the License, or (at your option) any later version.
00013 //
00014 //  This library is distributed in the hope that it will be useful,
00015 //  but WITHOUT ANY WARRANTY; without even the implied warranty of
00016 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00017 //  Lesser General Public License for more details.
00018 //
00019 //  You should have received a copy of the GNU Lesser General Public
00020 //  License along with this library; if not, write to the Free Software
00021 //  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00022 //
00023 // --------------------------------------------------------------------------
00024 // $Maintainer: Marc Sturm $
00025 // --------------------------------------------------------------------------
00026 
00027 #ifndef OPENMS_DATASTRUCTURES_CONVEXHULL2D_H
00028 #define OPENMS_DATASTRUCTURES_CONVEXHULL2D_H
00029 
00030 #include <OpenMS/config.h>
00031 #include <OpenMS/CONCEPT/Types.h>
00032 #include <OpenMS/CONCEPT/Exception.h>
00033 #include <OpenMS/DATASTRUCTURES/DBoundingBox.h>
00034 
00035 #include <vector>
00036 #include <map>
00037 #include <math.h>
00038 
00039 #include <CGAL/Cartesian.h>
00040 #include <CGAL/convex_hull_2.h>
00041 
00042 namespace OpenMS
00043 {
00049   class ConvexHull2D
00050   {
00051     public:
00052       typedef DPosition<2> PointType;
00053       typedef std::vector< PointType > PointArrayType;
00054       typedef PointArrayType::size_type SizeType;
00055       typedef PointArrayType::const_iterator PointArrayTypeConstIterator;
00057       ConvexHull2D()
00058         : points_()
00059       {
00060       }
00061 
00063       ConvexHull2D(const PointArrayType& points)
00064       {
00065         *this = points;
00066       }
00067       
00069       ConvexHull2D& operator=(const ConvexHull2D& rhs)
00070       {
00071         if (&rhs == this) return *this;
00072         
00073         points_ = rhs.points_;
00074         
00075         return *this;
00076       }
00077       
00078 
00080       ConvexHull2D& operator=(const PointArrayType& points)
00081       {
00082         //init
00083         points_.clear();
00084         if (points.size()==0) return *this;
00085         //convert input to cgal
00086         std::vector<Point_2> cgal_points;
00087         for (PointArrayTypeConstIterator it = points.begin(); it!=points.end(); ++it)
00088         {
00089           cgal_points.push_back( Point_2((*it)[0], (*it)[1]) );    
00090         }
00091         //calculate convex hull
00092         std::vector<Point_2> cgal_result;
00093         CGAL::convex_hull_2( cgal_points.begin(), cgal_points.end(), std::inserter(cgal_result, cgal_result.begin() ) );
00094         // add the points
00095         for (std::vector<Point_2>::const_iterator cit = cgal_result.begin(); cit != cgal_result.end(); ++cit)
00096         {
00097           points_.push_back( PointType(cit->x(),cit->y()) );      
00098         }   
00099         return *this;
00100       }
00101 
00103       bool operator==(const ConvexHull2D& rhs) const
00104       {
00105         // different size => return false
00106         if (points_.size()!= rhs.points_.size()) return false;
00107         
00108         //different points now => return false
00109         for (PointArrayTypeConstIterator it = rhs.points_.begin(); it !=  rhs.points_.end(); ++it)
00110         {
00111           if (find(points_.begin(),points_.end(),*it)==points_.end())
00112           {
00113             return false;
00114           }
00115         }
00116         return true;
00117       }
00118       
00120       void clear()
00121       {
00122         points_.clear();
00123       }
00124 
00126       const PointArrayType& getPoints() const
00127       {
00128         return points_;
00129       }
00130       
00132       DBoundingBox<2> getBoundingBox() const
00133       {
00134         DBoundingBox<2> bb;
00135         
00136         for (PointArrayTypeConstIterator it = points_.begin(); it!=points_.end(); ++it)
00137         {
00138           bb.enlarge(*it);
00139         }
00140         
00141         return bb;
00142       }
00143       
00145       bool addPoint(const PointType& point)
00146       {
00147         if (!encloses(point))
00148         {
00149           points_.push_back(point);
00150           return true;
00151         }
00152         return false;
00153       }
00154       
00163       bool encloses(const PointType& point) const
00164       {
00165         if (!getBoundingBox().encloses(point))
00166         {
00167           return false;
00168         }       
00169         // re-calculate convex hull
00170         std::vector<Point_2> cgal_points;
00171         std::vector<Point_2> cgal_result;
00172         
00173         for (PointArrayTypeConstIterator it = points_.begin(); it!=points_.end(); ++it)
00174         {
00175           cgal_points.push_back( Point_2((*it)[0], (*it)[1]) );    
00176         }
00177         points_.clear();
00178           
00179         CGAL::convex_hull_2( cgal_points.begin(), cgal_points.end(), std::inserter(cgal_result, cgal_result.begin() ) );
00180         // add the points
00181         for (std::vector<Point_2>::const_iterator cit = cgal_result.begin(); cit != cgal_result.end(); ++cit)
00182         {
00183           points_.push_back( PointType(cit->x(),cit->y()) );      
00184         }           
00185         
00186         //convert input to cgal
00187         cgal_points.clear();
00188         cgal_result.clear();
00189         for (PointArrayTypeConstIterator it = points_.begin(); it!=points_.end(); ++it)
00190         {
00191           cgal_points.push_back( Point_2((*it)[0], (*it)[1]) );    
00192         }
00193         // add point to be tested
00194         cgal_points.push_back( Point_2(point[0],point[1]) ); 
00195         
00196         //calculate convex hull
00197         CGAL::convex_hull_2( cgal_points.begin(), cgal_points.end(), std::inserter(cgal_result, cgal_result.begin() ) );
00198         
00199         //point added => return false
00200         if (cgal_result.size()!=points_.size()) return false;
00201         
00202         //different points now => return false
00203         for (std::vector<Point_2>::const_iterator cit = cgal_result.begin(); cit != cgal_result.end(); ++cit)
00204         {
00205           if (find(points_.begin(),points_.end(),PointType(cit->x(),cit->y()))==points_.end())
00206           {
00207             return false;
00208           }
00209         }
00210         return true;
00211       }
00212       
00213     protected:
00214       mutable PointArrayType points_;
00215       typedef CGAL::Cartesian<DoubleReal>::Point_2 Point_2;
00216   };
00217 } // namespace OPENMS
00218 
00219 #endif // OPENMS_DATASTRUCTURES_DCONVEXHULL_H

Generated Tue Apr 1 15:36:33 2008 -- using doxygen 1.5.4 OpenMS / TOPP 1.1