mirror of https://github.com/CGAL/cgal
Added PCA for 3d. Also added tags for dimensions 0,1,2 and 3.
A Principal_component_analysis/include/CGAL/linear_least_squares_fitting_tetrahedrons_3.h A Principal_component_analysis/include/CGAL/linear_least_squares_fitting_spheres_3.h M Principal_component_analysis/include/CGAL/linear_least_squares_fitting_points_2.h M Principal_component_analysis/include/CGAL/centroid.h A Principal_component_analysis/include/CGAL/linear_least_squares_fitting_points_3.h M Principal_component_analysis/include/CGAL/linear_least_squares_fitting_2.h M Principal_component_analysis/include/CGAL/linear_least_squares_fitting_3.h M Principal_component_analysis/include/CGAL/util.h M Principal_component_analysis/include/CGAL/linear_least_squares_fitting_circles_2.h M Principal_component_analysis/include/CGAL/linear_least_squares_fitting_segments_2.h A Principal_component_analysis/include/CGAL/linear_least_squares_fitting_segments_3.h M Principal_component_analysis/include/CGAL/linear_least_squares_fitting_rectangles_2.h A Principal_component_analysis/include/CGAL/PCA_tags.h M Principal_component_analysis/include/CGAL/linear_least_squares_fitting_triangles_2.h A Principal_component_analysis/include/CGAL/linear_least_squares_fitting_triangles_3.h A Principal_component_analysis/include/CGAL/linear_least_squares_fitting_cuboids_3.h D Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_segments.cpp D Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_rectangles.cpp D Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_triangles.cpp A Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_spheres_3.cpp A Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_points_2.cpp M Principal_component_analysis/examples/Principal_component_analysis/centroid.cpp A Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_points_3.cpp A Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_circles_2.cpp A Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_ankit.cpp A Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_segments_2.cpp A Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_segments_3.cpp A Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_rectangles_2.cpp A Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_triangles_2.cpp A Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_cuboids_3.cpp A Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_triangles_3.cpp D Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting.cpp D Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_circles.cpp A Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_tetrahedrons_3.cpp
This commit is contained in:
Ankit Gupta
parent
9f430ac2e9
commit
9584e4e0af
|
|
@ -1896,15 +1896,30 @@ Principal_component_analysis/demo/Principal_component_analysis/windows/3d/res/id
|
|||
Principal_component_analysis/doc_tex/Principal_component_analysis/fit.eps -text svneol=unset#application/postscript
|
||||
Principal_component_analysis/doc_tex/Principal_component_analysis/fit.png -text svneol=unset#image/png
|
||||
Principal_component_analysis/doc_tex/Principal_component_analysis/teaserLeastSquaresFitting.png -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_circles.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_rectangles.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_segments.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_triangles.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_ankit.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_circles_2.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_cuboids_3.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_points_2.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_points_3.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_rectangles_2.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_segments_2.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_segments_3.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_spheres_3.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_tetrahedrons_3.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_triangles_2.cpp -text
|
||||
Principal_component_analysis/examples/Principal_component_analysis/linear_least_squares_fitting_triangles_3.cpp -text
|
||||
Principal_component_analysis/include/CGAL/PCA_tags.h -text
|
||||
Principal_component_analysis/include/CGAL/linear_least_squares_fitting_circles_2.h -text
|
||||
Principal_component_analysis/include/CGAL/linear_least_squares_fitting_cuboids_3.h -text
|
||||
Principal_component_analysis/include/CGAL/linear_least_squares_fitting_points_2.h -text
|
||||
Principal_component_analysis/include/CGAL/linear_least_squares_fitting_points_3.h -text
|
||||
Principal_component_analysis/include/CGAL/linear_least_squares_fitting_rectangles_2.h -text
|
||||
Principal_component_analysis/include/CGAL/linear_least_squares_fitting_segments_2.h -text
|
||||
Principal_component_analysis/include/CGAL/linear_least_squares_fitting_segments_3.h -text
|
||||
Principal_component_analysis/include/CGAL/linear_least_squares_fitting_spheres_3.h -text
|
||||
Principal_component_analysis/include/CGAL/linear_least_squares_fitting_tetrahedrons_3.h -text
|
||||
Principal_component_analysis/include/CGAL/linear_least_squares_fitting_triangles_2.h -text
|
||||
Principal_component_analysis/include/CGAL/linear_least_squares_fitting_triangles_3.h -text
|
||||
Principal_component_analysis/include/CGAL/util.h -text
|
||||
Principal_component_analysis/test/Principal_component_analysis/fitting_segments_2.cpp -text
|
||||
QP_solver/doc_tex/QP_solver/PkgDescription.tex -text
|
||||
|
|
|
|||
|
|
@ -85,41 +85,41 @@ BOOL CMeshDoc::OnOpenDocument(LPCTSTR lpszPathName)
|
|||
// off extension
|
||||
if(extension == ".off")
|
||||
{
|
||||
// read from stream
|
||||
std::ifstream stream(lpszPathName);
|
||||
if(!stream)
|
||||
{
|
||||
AfxMessageBox("Unable to open file");
|
||||
return false;
|
||||
}
|
||||
stream >> m_mesh;
|
||||
|
||||
// add mesh points to point set
|
||||
for(Mesh::Point_iterator it = m_mesh.points_begin();
|
||||
it != m_mesh.points_end();
|
||||
it++)
|
||||
m_points.push_back(*it);
|
||||
|
||||
// add mesh triangles to triangle set
|
||||
for(Mesh::Facet_iterator f = m_mesh.facets_begin();
|
||||
f != m_mesh.facets_end();
|
||||
f++)
|
||||
{
|
||||
const Point& a = f->halfedge()->vertex()->point();
|
||||
const Point& b = f->halfedge()->next()->vertex()->point();
|
||||
const Point& c = f->halfedge()->next()->next()->vertex()->point();
|
||||
m_triangles.push_back(Triangle(a,b,c));
|
||||
}
|
||||
// read from stream
|
||||
std::ifstream stream(lpszPathName);
|
||||
if(!stream)
|
||||
{
|
||||
AfxMessageBox("Unable to open file");
|
||||
return false;
|
||||
}
|
||||
stream >> m_mesh;
|
||||
|
||||
// add mesh points to point set
|
||||
for(Mesh::Point_iterator it = m_mesh.points_begin();
|
||||
it != m_mesh.points_end();
|
||||
it++)
|
||||
m_points.push_back(*it);
|
||||
|
||||
// add mesh triangles to triangle set
|
||||
for(Mesh::Facet_iterator f = m_mesh.facets_begin();
|
||||
f != m_mesh.facets_end();
|
||||
f++)
|
||||
{
|
||||
const Point& a = f->halfedge()->vertex()->point();
|
||||
const Point& b = f->halfedge()->next()->vertex()->point();
|
||||
const Point& c = f->halfedge()->next()->next()->vertex()->point();
|
||||
m_triangles.push_back(Triangle(a,b,c));
|
||||
}
|
||||
}
|
||||
else
|
||||
{
|
||||
AfxMessageBox("Unknown extension");
|
||||
return false;
|
||||
}
|
||||
|
||||
OnFitFitpointset();
|
||||
|
||||
{
|
||||
AfxMessageBox("Unknown extension");
|
||||
return false;
|
||||
}
|
||||
|
||||
OnFitFitpointset();
|
||||
|
||||
|
||||
return TRUE;
|
||||
}
|
||||
|
||||
|
|
|
|||
|
|
@ -19,7 +19,7 @@ int main()
|
|||
points_2.push_back(Point_2(2.0, 2.0));
|
||||
points_2.push_back(Point_2(3.0, 5.0));
|
||||
|
||||
Point_2 c2 = CGAL::centroid(points_2.begin(), points_2.end());
|
||||
Point_2 c2 = CGAL::centroid(points_2.begin(), points_2.end(),CGAL::PCA_dimension_0_tag());
|
||||
std::cout << c2 << std::endl;
|
||||
|
||||
// centroid of 3D points
|
||||
|
|
@ -28,7 +28,7 @@ int main()
|
|||
points_3.push_back(Point_3(2.0, 2.0, 1.2));
|
||||
points_3.push_back(Point_3(3.0, 5.0, 4.5));
|
||||
|
||||
Point_3 c3 = CGAL::centroid(points_3.begin(), points_3.end());
|
||||
Point_3 c3 = CGAL::centroid(points_3.begin(), points_3.end(),CGAL::PCA_dimension_0_tag());
|
||||
std::cout << c3 << std::endl;
|
||||
|
||||
return 0;
|
||||
|
|
|
|||
|
|
@ -0,0 +1,27 @@
|
|||
// Example program for the linear_least_square_fitting function
|
||||
|
||||
#include <CGAL/Cartesian.h>
|
||||
#include <CGAL/linear_least_squares_fitting_2.h>
|
||||
|
||||
#include <list>
|
||||
typedef double FT;
|
||||
typedef CGAL::Cartesian<FT> K;
|
||||
typedef K::Line_2 Line_2;
|
||||
typedef K::Point_2 Point_2;
|
||||
|
||||
int main(int argc, char** argv)
|
||||
{
|
||||
std::list<Point_2> points;
|
||||
points.push_back(Point_2(1.0,0.0));
|
||||
points.push_back(Point_2(0.0,1.0));
|
||||
points.push_back(Point_2(1.0,1.0));
|
||||
points.push_back(Point_2(-1.0,1.0));
|
||||
points.push_back(Point_2(-1.0,0.0));
|
||||
// points.push_back(Point_2(0.0,1.0));
|
||||
|
||||
Line_2 line;
|
||||
linear_least_squares_fitting_2(points.begin(),points.end(),line,CGAL::PCA_dimension_0_tag());
|
||||
|
||||
std::cout<<line<<std::endl;
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -15,10 +15,11 @@ int main()
|
|||
{
|
||||
std::list<Circle_2> circles;
|
||||
circles.push_back(Circle_2(Point_2(0.0,0.0),9));
|
||||
circles.push_back(Circle_2(Point_2(1.0,1.0),25));
|
||||
circles.push_back(Circle_2(Point_2(0.0,10.0),49));
|
||||
circles.push_back(Circle_2(Point_2(10.0,0.0),49));
|
||||
|
||||
Line_2 line;
|
||||
FT i = linear_least_squares_fitting_2(circles.begin(),circles.end(),line);
|
||||
FT i = linear_least_squares_fitting_2(circles.begin(),circles.end(),line,CGAL::PCA_dimension_1_tag());
|
||||
|
||||
std::cout<<"accuracy: "<<i<<std::endl;
|
||||
std::cout<<line<<std::endl;
|
||||
|
|
@ -0,0 +1,35 @@
|
|||
// Example program for the linear_least_square_fitting function
|
||||
|
||||
#include <CGAL/Cartesian.h>
|
||||
#include <CGAL/linear_least_squares_fitting_3.h>
|
||||
|
||||
#include <list>
|
||||
|
||||
typedef double FT;
|
||||
typedef CGAL::Cartesian<FT> K;
|
||||
typedef K::Line_3 Line_3;
|
||||
typedef K::Plane_3 Plane_3;
|
||||
typedef K::Point_3 Point_3;
|
||||
typedef K::Iso_cuboid_3 Iso_cuboid_3;
|
||||
|
||||
int main()
|
||||
{
|
||||
std::list<Iso_cuboid_3> Iso_cuboids;
|
||||
Iso_cuboids.push_back(Iso_cuboid_3(Point_3(0.0,0.0,0.0),Point_3(1.0,1.0,1.0)));
|
||||
// Iso_cuboids.push_back(Iso_cuboid_3(Point_3(4.0,0.0),Point_3(5.0,1.0)));
|
||||
|
||||
Line_3 line;
|
||||
FT i = linear_least_squares_fitting_3(Iso_cuboids.begin(),Iso_cuboids.end(),line,CGAL::PCA_dimension_3_tag());
|
||||
|
||||
std::cout<<"Line's accuracy: "<<i<<std::endl;
|
||||
std::cout<<"Line direction: "<<line.to_vector()<<std::endl;
|
||||
std::cout<<"Line's point: "<<line.point(0)<<std::endl;
|
||||
|
||||
Plane_3 plane;
|
||||
FT j = linear_least_squares_fitting_3(Iso_cuboids.begin(),Iso_cuboids.end(),plane,CGAL::PCA_dimension_3_tag());
|
||||
|
||||
std::cout<<"Plane's accuracy: "<<j<<std::endl;
|
||||
std::cout<<"Plane: "<<plane<<std::endl;
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -18,7 +18,7 @@ int main()
|
|||
points.push_back(Point_2(3.0,0.0));
|
||||
|
||||
Line_2 line;
|
||||
linear_least_squares_fitting_2(points.begin(),points.end(),line);
|
||||
linear_least_squares_fitting_2(points.begin(),points.end(),line,CGAL::PCA_dimension_0_tag());
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -0,0 +1,24 @@
|
|||
// Example program for the linear_least_square_fitting function
|
||||
|
||||
#include <CGAL/Cartesian.h>
|
||||
#include <CGAL/linear_least_squares_fitting_3.h>
|
||||
|
||||
#include <list>
|
||||
|
||||
typedef double FT;
|
||||
typedef CGAL::Cartesian<FT> K;
|
||||
typedef K::Line_3 Line_3;
|
||||
typedef K::Point_3 Point_3;
|
||||
|
||||
int main()
|
||||
{
|
||||
std::list<Point_3> points;
|
||||
points.push_back(Point_3(1.0,0.0,0.0));
|
||||
points.push_back(Point_3(2.0,0.0,0.0));
|
||||
points.push_back(Point_3(3.0,0.0,0.0));
|
||||
|
||||
Line_3 line;
|
||||
linear_least_squares_fitting_3(points.begin(),points.end(),line,CGAL::PCA_dimension_0_tag());
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -14,11 +14,11 @@ typedef K::Iso_rectangle_2 Iso_rectangle_2;
|
|||
int main()
|
||||
{
|
||||
std::list<Iso_rectangle_2> Iso_rectangles;
|
||||
Iso_rectangles.push_back(Iso_rectangle_2(Point_2(0.0,0.0),Point_2(1.0,1.0)));
|
||||
Iso_rectangles.push_back(Iso_rectangle_2(Point_2(0.0,0.0),Point_2(4.0,4.0)));
|
||||
Iso_rectangles.push_back(Iso_rectangle_2(Point_2(4.0,0.0),Point_2(5.0,1.0)));
|
||||
|
||||
Line_2 line;
|
||||
FT i = linear_least_squares_fitting_2(Iso_rectangles.begin(),Iso_rectangles.end(),line);
|
||||
FT i = linear_least_squares_fitting_2(Iso_rectangles.begin(),Iso_rectangles.end(),line,CGAL::PCA_dimension_0_tag());
|
||||
|
||||
std::cout<<"accuracy: "<<i<<std::endl;
|
||||
std::cout<<line<<std::endl;
|
||||
|
|
@ -14,13 +14,16 @@ typedef K::Segment_2 Segment_2;
|
|||
int main()
|
||||
{
|
||||
std::list<Segment_2> segments;
|
||||
segments.push_back(Segment_2(Point_2(0.0,0.0),Point_2(0.0,1.0)));
|
||||
segments.push_back(Segment_2(Point_2(0.0,2.0),Point_2(0.0,10.0)));
|
||||
// segments.push_back(Segment_2(Point_2(0.0,0.0),Point_2(0.0,1.0)));
|
||||
segments.push_back(Segment_2(Point_2(0.0,1.0),Point_2(1.0,0.0)));
|
||||
// segments.push_back(Segment_2(Point_2(0.0,0.0),Point_2(1.0,0.0)));
|
||||
|
||||
Line_2 line;
|
||||
FT i = linear_least_squares_fitting_2(segments.begin(),segments.end(),line);
|
||||
Point_2 c;
|
||||
FT i = linear_least_squares_fitting_2(segments.begin(),segments.end(),line,c,CGAL::PCA_dimension_1_tag());
|
||||
|
||||
std::cout<<"accuracy: "<<i<<std::endl;
|
||||
std::cout<<line<<std::endl;
|
||||
std::cout<<"centroid: "<<c<<std::endl;
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -0,0 +1,37 @@
|
|||
// Example program for the linear_least_square_fitting function
|
||||
|
||||
#include <CGAL/Cartesian.h>
|
||||
#include <CGAL/linear_least_squares_fitting_3.h>
|
||||
|
||||
#include <list>
|
||||
|
||||
typedef double FT;
|
||||
typedef CGAL::Cartesian<FT> K;
|
||||
typedef K::Line_3 Line_3;
|
||||
typedef K::Plane_3 Plane_3;
|
||||
typedef K::Point_3 Point_3;
|
||||
typedef K::Segment_3 Segment_3;
|
||||
|
||||
int main()
|
||||
{
|
||||
std::list<Segment_3> segments;
|
||||
segments.push_back(Segment_3(Point_3(1.0,1.0,1.0),Point_3(2.0,2.0,2.0)));
|
||||
// segments.push_back(Segment_3(Point_3(0.0,1.0),Point_3(1.0,0.0)));
|
||||
// segments.push_back(Segment_3(Point_3(0.0,0.0),Point_3(1.0,0.0)));
|
||||
|
||||
Line_3 line;
|
||||
Point_3 c;
|
||||
FT i = linear_least_squares_fitting_3(segments.begin(),segments.end(),line,c,CGAL::PCA_dimension_1_tag());
|
||||
|
||||
std::cout<<"accuracy: "<<i<<std::endl;
|
||||
std::cout<<line<<std::endl;
|
||||
std::cout<<"centroid: "<<c<<std::endl;
|
||||
|
||||
Plane_3 plane;
|
||||
i = linear_least_squares_fitting_3(segments.begin(),segments.end(),plane,c,CGAL::PCA_dimension_1_tag());
|
||||
|
||||
std::cout<<"accuracy: "<<i<<std::endl;
|
||||
std::cout<<plane<<std::endl;
|
||||
std::cout<<"centroid: "<<c<<std::endl;
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -0,0 +1,33 @@
|
|||
// Example program for the linear_least_square_fitting function
|
||||
|
||||
#include <CGAL/Cartesian.h>
|
||||
#include <CGAL/linear_least_squares_fitting_3.h>
|
||||
|
||||
#include <list>
|
||||
|
||||
typedef double FT;
|
||||
typedef CGAL::Cartesian<FT> K;
|
||||
typedef K::Line_3 Line_3;
|
||||
typedef K::Plane_3 Plane_3;
|
||||
typedef K::Point_3 Point_3;
|
||||
typedef K::Sphere_3 Sphere_3;
|
||||
|
||||
int main()
|
||||
{
|
||||
std::list<Sphere_3> spheres;
|
||||
spheres.push_back(Sphere_3(Point_3(0.0,0.0,0.0),9));
|
||||
spheres.push_back(Sphere_3(Point_3(0.0,10.0,0.0),25));
|
||||
|
||||
Line_3 line;
|
||||
FT i = linear_least_squares_fitting_3(spheres.begin(),spheres.end(),line,CGAL::PCA_dimension_3_tag());
|
||||
|
||||
std::cout<<"Line's accuracy: "<<i<<std::endl;
|
||||
std::cout<<"Line: "<<line<<std::endl;
|
||||
|
||||
Plane_3 plane;
|
||||
FT j = linear_least_squares_fitting_3(spheres.begin(),spheres.end(),plane,CGAL::PCA_dimension_3_tag());
|
||||
|
||||
std::cout<<"Plane's accuracy: "<<j<<std::endl;
|
||||
std::cout<<"Plane: "<<plane<<std::endl;
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -0,0 +1,37 @@
|
|||
// Example program for the linear_least_square_fitting function
|
||||
|
||||
#include <CGAL/Cartesian.h>
|
||||
#include <CGAL/linear_least_squares_fitting_3.h>
|
||||
|
||||
#include <list>
|
||||
|
||||
typedef double FT;
|
||||
typedef CGAL::Cartesian<FT> K;
|
||||
typedef K::Line_3 Line_3;
|
||||
typedef K::Plane_3 Plane_3;
|
||||
typedef K::Point_3 Point_3;
|
||||
typedef K::Tetrahedron_3 Tetrahedron_3;
|
||||
|
||||
int main()
|
||||
{
|
||||
std::list<Tetrahedron_3> tetrahedrons;
|
||||
tetrahedrons.push_back(Tetrahedron_3(Point_3(0.0,0.0,0.0),Point_3(1.0,0.0,0.0),Point_3(0.0,1.0,0.0),Point_3(0.0,0.0,1.0)));
|
||||
|
||||
Line_3 line;
|
||||
FT i = linear_least_squares_fitting_3(tetrahedrons.begin(),tetrahedrons.end(),line,CGAL::PCA_dimension_3_tag());
|
||||
|
||||
std::cout<<"Line's accuracy: "<<i<<std::endl;
|
||||
std::cout<<"Line direction: "<<line.to_vector()<<std::endl;
|
||||
std::cout<<"Line's point: "<<line.point(0)<<std::endl;
|
||||
|
||||
if(line.has_on(Point_3(2/3.0,0,0)) && line.has_on(Point_3(0,2/3.0,0)))
|
||||
std::cout<<"YES!!"<<std::endl;
|
||||
|
||||
Plane_3 plane;
|
||||
FT j = linear_least_squares_fitting_3(tetrahedrons.begin(),tetrahedrons.end(),plane,CGAL::PCA_dimension_3_tag());
|
||||
|
||||
std::cout<<"Plane's accuracy: "<<j<<std::endl;
|
||||
std::cout<<"Plane: "<<plane<<std::endl;
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -14,13 +14,15 @@ typedef K::Triangle_2 Triangle_2;
|
|||
int main()
|
||||
{
|
||||
std::list<Triangle_2> triangles;
|
||||
triangles.push_back(Triangle_2(Point_2(1.0,0.0),Point_2(0.0,1.0),Point_2(1.0,1.0)));
|
||||
triangles.push_back(Triangle_2(Point_2(-1.0,0.0),Point_2(0.0,1.0),Point_2(-1.0,1.0)));
|
||||
Point_2 c;
|
||||
triangles.push_back(Triangle_2(Point_2(0.0,0.0),Point_2(0.0,1.0),Point_2(1.0,0.0)));
|
||||
// triangles.push_back(Triangle_2(Point_2(-1.0,0.0),Point_2(0.0,1.0),Point_2(-1.0,1.0)));
|
||||
|
||||
Line_2 line;
|
||||
FT i = linear_least_squares_fitting_2(triangles.begin(),triangles.end(),line);
|
||||
FT i = linear_least_squares_fitting_2(triangles.begin(),triangles.end(),line,c,CGAL::PCA_dimension_0_tag());
|
||||
|
||||
std::cout<<"accuracy: "<<i<<std::endl;
|
||||
std::cout<<line<<std::endl;
|
||||
std::cout<<"centroid: "<<c<<std::endl;
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -0,0 +1,37 @@
|
|||
// Example program for the linear_least_square_fitting function
|
||||
|
||||
#include <CGAL/Cartesian.h>
|
||||
#include <CGAL/linear_least_squares_fitting_3.h>
|
||||
|
||||
#include <list>
|
||||
|
||||
typedef double FT;
|
||||
typedef CGAL::Cartesian<FT> K;
|
||||
typedef K::Line_3 Line_3;
|
||||
typedef K::Plane_3 Plane_3;
|
||||
typedef K::Point_3 Point_3;
|
||||
typedef K::Triangle_3 Triangle_3;
|
||||
|
||||
int main()
|
||||
{
|
||||
std::list<Triangle_3> triangles;
|
||||
triangles.push_back(Triangle_3(Point_3(1.0,0.0,0.0),Point_3(0.0,1.0,0.0),Point_3(0.0,0.0,0.0)));
|
||||
|
||||
Line_3 line;
|
||||
FT i = linear_least_squares_fitting_3(triangles.begin(),triangles.end(),line,CGAL::PCA_dimension_2_tag());
|
||||
|
||||
std::cout<<"Line's accuracy: "<<i<<std::endl;
|
||||
std::cout<<"Line direction: "<<line.to_vector()<<std::endl;
|
||||
std::cout<<"Line's point: "<<line.point(0)<<std::endl;
|
||||
|
||||
if(line.has_on(Point_3(2/3.0,0,0)) && line.has_on(Point_3(0,2/3.0,0)))
|
||||
std::cout<<"YES!!"<<std::endl;
|
||||
|
||||
Plane_3 plane;
|
||||
FT j = linear_least_squares_fitting_3(triangles.begin(),triangles.end(),plane,CGAL::PCA_dimension_2_tag());
|
||||
|
||||
std::cout<<"Plane's accuracy: "<<j<<std::endl;
|
||||
std::cout<<"Plane: "<<plane<<std::endl;
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
|
@ -0,0 +1,47 @@
|
|||
// Copyright (c) 2005 INRIA Sophia-Antipolis (France).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org); you may redistribute it under
|
||||
// the terms of the Q Public License version 1.0.
|
||||
// See the file LICENSE.QPL distributed with CGAL.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL: svn+ssh://gankit@scm.gforge.inria.fr/svn/cgal/trunk/Principal_component_analysis/include/CGAL/tags.h $
|
||||
// $Id: tags.h 37882 2007-04-03 15:15:30Z spion $
|
||||
//
|
||||
// Author(s) : Pierre Alliez and Sylvain Pion and Ankit Gupta
|
||||
|
||||
#ifndef CGAL_LINEAR_LEAST_SQUARES_FITTING_TAGS_H
|
||||
#define CGAL_LINEAR_LEAST_SQUARES_FITTING_TAGS_H
|
||||
|
||||
CGAL_BEGIN_NAMESPACE
|
||||
|
||||
//::::::::::::::::::::::::::::::::::::::::::::::::::::::
|
||||
//::These Tags are meant to represent the type of
|
||||
//::geometry we are dealing with. The CGAL kernel
|
||||
//::provides only for circles, triangles, spheres,
|
||||
//: rectangles and so on. But these can be solid,
|
||||
//::hollow, just outlined etc. These tags resolve
|
||||
//::between these different types of the same geometry
|
||||
//::::::::::::::::::::::::::::::::::::::::::::::::::::::
|
||||
|
||||
// Fully filled geometry object, such as disks in 2d and balls in 3d
|
||||
struct PCA_dimension_3_tag {};
|
||||
|
||||
// Empty geometry objects, such as triangles in 2d and hollow cubes with just the planar facets in 3d
|
||||
struct PCA_dimension_2_tag {};
|
||||
|
||||
// Only relevant for 3d objects, outlines of 3d geometry objects formed by just lines joing appropriate vetrices, such as outline of a cuboid.
|
||||
struct PCA_dimension_1_tag {};
|
||||
|
||||
// For the vertices of objects.
|
||||
struct PCA_dimension_0_tag {};
|
||||
|
||||
CGAL_END_NAMESPACE
|
||||
|
||||
#endif // CGAL_LINEAR_LEAST_SQUARES_FITTING_TAGS_H
|
||||
|
|
@ -25,7 +25,9 @@
|
|||
#include <iterator>
|
||||
#include <CGAL/Kernel_traits.h>
|
||||
#include <CGAL/Kernel/Dimension_utils.h>
|
||||
#include <CGAL/PCA_tags.h>
|
||||
|
||||
#include <list>
|
||||
// Functions to compute the centroid of N points.
|
||||
// Works in 2D and 3D.
|
||||
|
||||
|
|
@ -38,6 +40,8 @@ CGAL_BEGIN_NAMESPACE
|
|||
|
||||
namespace CGALi {
|
||||
|
||||
//:::::::::: 2D Objects :::::::::::::::::::
|
||||
|
||||
// computes the centroid of a 2D point set
|
||||
// takes an iterator range over K::Point_2
|
||||
template < typename InputIterator,
|
||||
|
|
@ -46,7 +50,8 @@ typename K::Point_2
|
|||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K&,
|
||||
const typename K::Point_2*)
|
||||
const typename K::Point_2*,
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::Vector_2 Vector;
|
||||
typedef typename K::FT FT;
|
||||
|
|
@ -63,40 +68,45 @@ centroid(InputIterator begin,
|
|||
return ORIGIN + v / (FT)nb_pts;
|
||||
}// end centroid of a 2D point set
|
||||
|
||||
// computes the centroid of a 3D point set
|
||||
// takes an iterator range over K::Point_3
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Point_3*)
|
||||
{
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::FT FT;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
Vector v = NULL_VECTOR;
|
||||
unsigned int nb_pts = 0;
|
||||
while (begin != end)
|
||||
{
|
||||
v = v + (*begin++ - ORIGIN);
|
||||
nb_pts++;
|
||||
}
|
||||
return ORIGIN + v / (FT)nb_pts;
|
||||
}// end centroid of a 3D point set
|
||||
|
||||
// computes the centroid of a 2D segment set
|
||||
// takes an iterator range over K::Segment_2
|
||||
|
||||
// centriod for 2D segment set with 0D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_2
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& k,
|
||||
const typename K::Segment_2*,
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::Point_2 Point;
|
||||
typedef typename K::Segment_2 Segment;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Segment& s = *it;
|
||||
points.push_back(s[0]);
|
||||
points.push_back(s[1]);
|
||||
}
|
||||
return centroid(points.begin(),points.end(),k,(Point*)NULL,tag);
|
||||
}// end centriod for 2D segment set with 0D tag
|
||||
|
||||
// centriod for 2D segment set with 1D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_2
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Segment_2*)
|
||||
const typename K::Segment_2*,
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Vector_2 Vector;
|
||||
|
|
@ -119,17 +129,78 @@ centroid(InputIterator begin,
|
|||
}
|
||||
CGAL_assertion(sum_lengths != 0.0);
|
||||
return ORIGIN + v / sum_lengths;
|
||||
} // end centroid of a 2D segment set
|
||||
} // end centroid of a 2D segment set with 1D tag
|
||||
|
||||
// computes the centroid of a 2D triangle set
|
||||
// takes an iterator range over K::Triangle_2
|
||||
|
||||
// centriod for 2D triangle set with 0D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_2
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& k,
|
||||
const typename K::Triangle_2*,
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::Triangle_2 Triangle;
|
||||
typedef typename K::Point_2 Point;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Triangle& triangle = *it;
|
||||
points.push_back(triangle[0]);
|
||||
points.push_back(triangle[1]);
|
||||
points.push_back(triangle[2]);
|
||||
}
|
||||
return centroid(points.begin(),points.end(),k,(Point*)NULL,tag);
|
||||
|
||||
} // end centroid of a 2D triangle set with 0D tag
|
||||
|
||||
// centriod for 2D triangle set with 1D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_2
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& k,
|
||||
const typename K::Triangle_2*,
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::Triangle_2 Triangle;
|
||||
typedef typename K::Segment_2 Segment;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
std::list<Segment> segments;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Triangle& triangle = *it;
|
||||
segments.push_back(triangle[0],triangle[1]);
|
||||
segments.push_back(triangle[1],triangle[2]);
|
||||
segments.push_back(triangle[2],triangle[0]);
|
||||
}
|
||||
return centroid(segments.begin(),segments.end(),k,(Segment*)NULL,tag);
|
||||
|
||||
} // end centroid of a 2D triangle set with 1D tag
|
||||
|
||||
// centriod for 2D triangle set with 2D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_2
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Triangle_2*)
|
||||
const typename K::Triangle_2*,
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Vector_2 Vector;
|
||||
|
|
@ -152,17 +223,53 @@ centroid(InputIterator begin,
|
|||
}
|
||||
CGAL_assertion(sum_areas != 0.0);
|
||||
return ORIGIN + v / sum_areas;
|
||||
} // end centroid of a 2D triangle set
|
||||
} // end centroid of a 2D triangle set with 2D tag
|
||||
|
||||
// computes the centroid of a 2D circle set
|
||||
// takes an iterator range over K::Circle_2
|
||||
|
||||
// centriod for 2D circle set with 1D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_2
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Circle_2*)
|
||||
const typename K::Circle_2*,
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Vector_2 Vector;
|
||||
typedef typename K::Point_2 Point;
|
||||
typedef typename K::Circle_2 Circle;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
Vector v = NULL_VECTOR;
|
||||
FT sum_lengths = 0;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Circle& s = *it;
|
||||
FT radius = std::sqrt(s.squared_radius());
|
||||
Point c = s.center();
|
||||
v = v + radius * (c - ORIGIN);
|
||||
sum_lengths += radius;
|
||||
}
|
||||
CGAL_assertion(sum_lengths != 0.0);
|
||||
return ORIGIN + v / sum_lengths;
|
||||
} // end centroid of a 2D circle set with 1D tag
|
||||
|
||||
// centriod for 2D circle set with 2D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_2
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Circle_2*,
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Vector_2 Vector;
|
||||
|
|
@ -179,24 +286,86 @@ centroid(InputIterator begin,
|
|||
{
|
||||
const Circle& s = *it;
|
||||
FT sq_radius = s.squared_radius();
|
||||
// Point c = K().construct_midpoint_2_object()(s[0],s[1]);
|
||||
Point c = s.center();
|
||||
v = v + sq_radius * (c - ORIGIN);
|
||||
sum_areas += sq_radius;
|
||||
}
|
||||
CGAL_assertion(sum_areas != 0.0);
|
||||
return ORIGIN + v / sum_areas;
|
||||
} // end centroid of a 2D circle set
|
||||
} // end centroid of a 2D circle set with 2D tag
|
||||
|
||||
// computes the centroid of a 2D rectangle set
|
||||
// takes an iterator range over K::Iso_Rectangle_2
|
||||
|
||||
// centriod for 2D rectangle set with 0D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_2
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& k,
|
||||
const typename K::Iso_rectangle_2*,
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::Iso_rectangle_2 Iso_rectangle;
|
||||
typedef typename K::Point_2 Point;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Iso_rectangle& r = *it;
|
||||
points.push_back(r[0]);
|
||||
points.push_back(r[1]);
|
||||
points.push_back(r[2]);
|
||||
points.push_back(r[3]);
|
||||
}
|
||||
return centroid(points.begin(),points.end(),k,(Point*)NULL,tag);
|
||||
|
||||
} // end centroid of a 2D rectangle set with 0D tag
|
||||
|
||||
// centriod for 2D rectangle set with 1D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_2
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& k,
|
||||
const typename K::Iso_rectangle_2*,
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::Iso_rectangle_2 Iso_rectangle;
|
||||
typedef typename K::Segment_2 Segment;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
std::list<Segment> segments;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Iso_rectangle& r = *it;
|
||||
segments.push_back(r[0],r[1]);
|
||||
segments.push_back(r[1],r[2]);
|
||||
segments.push_back(r[2],r[3]);
|
||||
segments.push_back(r[3],r[0]);
|
||||
}
|
||||
return centroid(segments.begin(),segments.end(),k,(Segment*)NULL,tag);
|
||||
|
||||
} // end centroid of a 2D rectangle set with 1D tag
|
||||
|
||||
// centriod for 2D rectangle set with 2D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_2
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Iso_rectangle_2*)
|
||||
const typename K::Iso_rectangle_2*,
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Vector_2 Vector;
|
||||
|
|
@ -219,17 +388,141 @@ centroid(InputIterator begin,
|
|||
}
|
||||
CGAL_assertion(sum_areas != 0.0);
|
||||
return ORIGIN + v / sum_areas;
|
||||
} // end centroid of a 2D rectangle set
|
||||
} // end centroid of a 2D rectangle set with 2D tag
|
||||
|
||||
// computes the centroid of a 3D triangle set
|
||||
// takes an iterator range over K::Triangle_3
|
||||
//::::::::::: 3D Objects ::::::::::::::::::::::::::
|
||||
|
||||
// computes the centroid of a 3D point set
|
||||
// takes an iterator range over K::Point_3
|
||||
|
||||
// centriod for 3D point set with 0D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K&,
|
||||
const typename K::Point_3*,
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::FT FT;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
Vector v = NULL_VECTOR;
|
||||
unsigned int nb_pts = 0;
|
||||
while (begin != end)
|
||||
{
|
||||
v = v + (*begin++ - ORIGIN);
|
||||
nb_pts++;
|
||||
}
|
||||
return ORIGIN + v / (FT)nb_pts;
|
||||
}// end centroid of a 3D point set with 0D tag
|
||||
|
||||
// centriod for 3D segment set with 1D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Triangle_3*)
|
||||
const typename K::Segment_3*,
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
Vector v = NULL_VECTOR;
|
||||
FT sum_lengths = 0;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Segment& s = *it;
|
||||
FT length = std::sqrt(std::abs(s.squared_length()));
|
||||
Point c = K().construct_midpoint_3_object()(s[0],s[1]);
|
||||
v = v + length * (c - ORIGIN);
|
||||
sum_lengths += length;
|
||||
}
|
||||
CGAL_assertion(sum_lengths != 0.0);
|
||||
return ORIGIN + v / sum_lengths;
|
||||
} // end centroid of a 3D segment set with 1D tag
|
||||
|
||||
// computes the centroid of a 3D triangle set
|
||||
// takes an iterator range over K::Triangle_3
|
||||
|
||||
// centriod for 3D triangle set with 0D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& k,
|
||||
const typename K::Triangle_3*,
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
typedef typename K::Point_3 Point;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Triangle& triangle = *it;
|
||||
points.push_back(triangle[0]);
|
||||
points.push_back(triangle[1]);
|
||||
points.push_back(triangle[2]);
|
||||
}
|
||||
return centroid(points.begin(),points.end(),k,(Point*)NULL,tag);
|
||||
|
||||
} // end centroid of a 3D triangle set with 0D tag
|
||||
|
||||
// centriod for 3D triangle set with 1D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& k,
|
||||
const typename K::Triangle_3*,
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
std::list<Segment> segments;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Triangle& triangle = *it;
|
||||
segments.push_back(triangle[0],triangle[1]);
|
||||
segments.push_back(triangle[1],triangle[2]);
|
||||
segments.push_back(triangle[2],triangle[0]);
|
||||
}
|
||||
return centroid(segments.begin(),segments.end(),k,(Segment*)NULL,tag);
|
||||
|
||||
} // end centroid of a 3D triangle set with 1D tag
|
||||
|
||||
// centriod for 3D triangle set with 2D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Triangle_3*,
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
|
|
@ -252,39 +545,284 @@ centroid(InputIterator begin,
|
|||
}
|
||||
CGAL_assertion(sum_areas != 0.0);
|
||||
return ORIGIN + v / sum_areas;
|
||||
} // end centroid of a 3D triangle set
|
||||
} // end centroid of a 3D triangle set with 2D tag
|
||||
|
||||
// computes the centroid of a 3D sphere set
|
||||
// takes an iterator range over K::Sphere_3
|
||||
|
||||
// centriod for 3D sphere set with 2D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Sphere_3*,
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Sphere_3 Sphere;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
Vector v = NULL_VECTOR;
|
||||
FT sum_areas = 0;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Sphere& sphere = *it;
|
||||
FT unsigned_area = sphere.squared_radius();
|
||||
Point c = sphere.center();
|
||||
v = v + unsigned_area * (c - ORIGIN);
|
||||
sum_areas += unsigned_area;
|
||||
}
|
||||
CGAL_assertion(sum_areas != 0.0);
|
||||
return ORIGIN + v / sum_areas;
|
||||
} // end centroid of a 3D sphere set with 2D tag
|
||||
|
||||
// centriod for 3D sphere set with 3D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Sphere_3*,
|
||||
const CGAL::PCA_dimension_3_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Sphere_3 Sphere;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
Vector v = NULL_VECTOR;
|
||||
FT sum_volumes = 0;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Sphere& sphere = *it;
|
||||
FT unsigned_volume = sphere.squared_radius() * std::sqrt(sphere.squared_radius());
|
||||
Point c = sphere.center();
|
||||
v = v + unsigned_volume * (c - ORIGIN);
|
||||
sum_volumes += unsigned_volume;
|
||||
}
|
||||
CGAL_assertion(sum_volumes != 0.0);
|
||||
return ORIGIN + v / sum_volumes;
|
||||
} // end centroid of a 3D sphere set with 3 tag
|
||||
|
||||
// computes the centroid of a 3D cuboid set
|
||||
// takes an iterator range over K::Iso_cuboid_3
|
||||
|
||||
// centriod for 3D cuboid set with 0D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& k,
|
||||
const typename K::Iso_cuboid_3*,
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
typedef typename K::Point_3 Point;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Iso_cuboid& cuboid = *it;
|
||||
points.push_back(cuboid[0]);
|
||||
points.push_back(cuboid[1]);
|
||||
points.push_back(cuboid[2]);
|
||||
points.push_back(cuboid[3]);
|
||||
points.push_back(cuboid[4]);
|
||||
points.push_back(cuboid[5]);
|
||||
points.push_back(cuboid[6]);
|
||||
points.push_back(cuboid[7]);
|
||||
}
|
||||
return centroid(points.begin(),points.end(),k,(Point*)NULL,tag);
|
||||
|
||||
} // end centroid of a 3D cuboid set with 0D tag
|
||||
|
||||
// centriod for 3D cuboid set with 1D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& k,
|
||||
const typename K::Iso_cuboid_3*,
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
std::list<Segment> segments;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Iso_cuboid& cuboid = *it;
|
||||
segments.push_back(cuboid[0],cuboid[1]);
|
||||
segments.push_back(cuboid[1],cuboid[2]);
|
||||
segments.push_back(cuboid[2],cuboid[3]);
|
||||
segments.push_back(cuboid[3],cuboid[0]);
|
||||
segments.push_back(cuboid[0],cuboid[5]);
|
||||
segments.push_back(cuboid[5],cuboid[4]);
|
||||
segments.push_back(cuboid[4],cuboid[3]);
|
||||
segments.push_back(cuboid[1],cuboid[6]);
|
||||
segments.push_back(cuboid[6],cuboid[7]);
|
||||
segments.push_back(cuboid[7],cuboid[2]);
|
||||
segments.push_back(cuboid[4],cuboid[7]);
|
||||
segments.push_back(cuboid[5],cuboid[6]);
|
||||
}
|
||||
return centroid(segments.begin(),segments.end(),k,(Segment*)NULL,tag);
|
||||
|
||||
} // end centroid of a 3D cuboid set with 1D tag
|
||||
|
||||
// centriod for 3D cuboid set with 2D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Iso_cuboid_3*,
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
Vector v = NULL_VECTOR;
|
||||
FT sum_areas = 0;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Iso_cuboid& cuboid = *it;
|
||||
FT unsigned_area = 2 * ((cuboid.xmax()-cuboid.xmin())*(cuboid.ymax()-cuboid.ymin()) + (cuboid.xmax()-cuboid.xmin())*(cuboid.zmax()-cuboid.zmin()) + (cuboid.ymax()-cuboid.ymin())*(cuboid.zmax()-cuboid.zmin()));
|
||||
Point c = K().construct_centroid_3_object()(cuboid[0],cuboid[1],cuboid[3],cuboid[5]);
|
||||
v = v + unsigned_area * (c - ORIGIN);
|
||||
sum_areas += unsigned_area;
|
||||
}
|
||||
CGAL_assertion(sum_areas != 0.0);
|
||||
return ORIGIN + v / sum_areas;
|
||||
} // end centroid of a 3D cuboid set with 2D tag
|
||||
|
||||
// centriod for 3D cuboid set with 3D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Iso_cuboid_3*,
|
||||
const CGAL::PCA_dimension_3_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
Vector v = NULL_VECTOR;
|
||||
FT sum_volumes = 0;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Iso_cuboid& cuboid = *it;
|
||||
FT unsigned_volume = cuboid.volume();
|
||||
Point c = K().construct_centroid_3_object()(cuboid[0],cuboid[1],cuboid[3],cuboid[5]);
|
||||
v = v + unsigned_volume * (c - ORIGIN);
|
||||
sum_volumes += unsigned_volume;
|
||||
}
|
||||
CGAL_assertion(sum_volumes != 0.0);
|
||||
return ORIGIN + v / sum_volumes;
|
||||
} // end centroid of a 3D cuboid set with 3D tag
|
||||
|
||||
// centriod for 3D Tetrahedron set with 3D tag
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::Point_3
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& ,
|
||||
const typename K::Tetrahedron_3*,
|
||||
const CGAL::PCA_dimension_3_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Tetrahedron_3 Tetrahedron;
|
||||
|
||||
CGAL_precondition(begin != end);
|
||||
|
||||
Vector v = NULL_VECTOR;
|
||||
FT sum_volumes = 0;
|
||||
for(InputIterator it = begin;
|
||||
it != end;
|
||||
it++)
|
||||
{
|
||||
const Tetrahedron& Tetrahedron = *it;
|
||||
FT unsigned_volume = Tetrahedron.volume();
|
||||
Point c = K().construct_centroid_3_object()(Tetrahedron[0],Tetrahedron[1],Tetrahedron[2],Tetrahedron[3]);
|
||||
v = v + unsigned_volume * (c - ORIGIN);
|
||||
sum_volumes += unsigned_volume;
|
||||
}
|
||||
CGAL_assertion(sum_volumes != 0.0);
|
||||
return ORIGIN + v / sum_volumes;
|
||||
} // end centroid of a 3D Tetrahedron set with 3D tag
|
||||
|
||||
} // namespace CGALi
|
||||
|
||||
// computes the centroid of a set of kernel objects
|
||||
// takes an iterator range over kernel objects
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K,
|
||||
typename tag >
|
||||
inline
|
||||
typename Point<Dimension<typename std::iterator_traits<InputIterator>::value_type, K>::value,
|
||||
K
|
||||
>::type
|
||||
centroid(InputIterator begin,
|
||||
InputIterator end,
|
||||
const K& k)
|
||||
const K& k,
|
||||
const tag& t)
|
||||
{
|
||||
typedef typename std::iterator_traits<InputIterator>::value_type Value_type;
|
||||
return CGALi::centroid(begin, end, k,(Value_type*) NULL);
|
||||
return CGALi::centroid(begin, end, k,(Value_type*) NULL, t);
|
||||
}
|
||||
|
||||
// this one takes an iterator range over kernel objects
|
||||
// and uses Kernel_traits<> to find out its kernel.
|
||||
template < typename InputIterator >
|
||||
template < typename InputIterator, typename tag >
|
||||
inline
|
||||
typename Point<Dimension<
|
||||
typename std::iterator_traits<InputIterator>::value_type,
|
||||
typename Kernel_traits<typename std::iterator_traits<InputIterator>::value_type>::Kernel >::value,
|
||||
typename Kernel_traits<typename std::iterator_traits<InputIterator>::value_type>::Kernel >::type
|
||||
centroid(InputIterator begin, InputIterator end)
|
||||
centroid(InputIterator begin, InputIterator end, const tag& t)
|
||||
{
|
||||
typedef typename std::iterator_traits<InputIterator>::value_type Point;
|
||||
typedef typename Kernel_traits<Point>::Kernel K;
|
||||
return CGAL::centroid(begin, end, K());
|
||||
return CGAL::centroid(begin, end, K(), t);
|
||||
}
|
||||
|
||||
CGAL_END_NAMESPACE
|
||||
|
|
|
|||
|
|
@ -21,19 +21,21 @@
|
|||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Object.h>
|
||||
#include <CGAL/centroid.h>
|
||||
#include <CGAL/Algebraic_structure_traits.h>
|
||||
#include <CGAL/IO/io.h>
|
||||
#include <CGAL/linear_least_squares_fitting_points_2.h>
|
||||
#include <CGAL/linear_least_squares_fitting_segments_2.h>
|
||||
#include <CGAL/linear_least_squares_fitting_triangles_2.h>
|
||||
#include <CGAL/linear_least_squares_fitting_circles_2.h>
|
||||
#include <CGAL/linear_least_squares_fitting_rectangles_2.h>
|
||||
#include <CGAL/PCA_tags.h>
|
||||
|
||||
#include <iterator>
|
||||
|
||||
CGAL_BEGIN_NAMESPACE
|
||||
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K , typename tag>
|
||||
inline
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_2(InputIterator first,
|
||||
|
|
@ -41,59 +43,63 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
typename K::Line_2& line,
|
||||
typename K::Point_2& centroid,
|
||||
const K& k,
|
||||
const bool non_standard_geometry = false)
|
||||
const tag& t)
|
||||
{
|
||||
typedef typename std::iterator_traits<InputIterator>::value_type Value_type;
|
||||
// BOOST_STATIC_ASSERT((boost::is_same<typename CGAL::Algebraic_structure_traits<Value_type>::Algebraic_category,CGAL::Field_with_sqrt_tag>::value));
|
||||
return CGALi::linear_least_squares_fitting_2(first, beyond, line,
|
||||
centroid, k, (Value_type*) NULL, non_standard_geometry);
|
||||
centroid, k, (Value_type*) NULL, t);
|
||||
}
|
||||
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K, typename tag >
|
||||
inline
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_2(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_2& line,
|
||||
const K& k,
|
||||
bool non_standard_geometry = false)
|
||||
const tag& t)
|
||||
{
|
||||
typedef typename std::iterator_traits<InputIterator>::value_type Value_type;
|
||||
// BOOST_STATIC_ASSERT((boost::is_same<typename CGAL::Algebraic_structure_traits<Value_type>::Algebraic_category,CGAL::Field_with_sqrt_tag>::value));
|
||||
typename K::Point_2 centroid;
|
||||
return CGALi::linear_least_squares_fitting_2(first, beyond, line,
|
||||
centroid, k,(Value_type*) NULL, non_standard_geometry);
|
||||
centroid, k,(Value_type*) NULL, t);
|
||||
}
|
||||
|
||||
// deduces the kernel from the points in container.
|
||||
template < typename InputIterator,
|
||||
typename Line,
|
||||
typename Point>
|
||||
typename Point, typename tag>
|
||||
inline
|
||||
typename Kernel_traits<Line>::Kernel::FT
|
||||
linear_least_squares_fitting_2(InputIterator first,
|
||||
InputIterator beyond,
|
||||
Line& line,
|
||||
Point& centroid,
|
||||
bool non_standard_geometry = false)
|
||||
const tag& t)
|
||||
{
|
||||
typedef typename std::iterator_traits<InputIterator>::value_type Value_type;
|
||||
// BOOST_STATIC_ASSERT((boost::is_same<typename CGAL::Algebraic_structure_traits<Value_type>::Algebraic_category,CGAL::Field_with_sqrt_tag>::value));
|
||||
typedef typename Kernel_traits<Value_type>::Kernel K;
|
||||
return CGAL::linear_least_squares_fitting_2(first,beyond,line,centroid,K(), non_standard_geometry);
|
||||
return CGAL::linear_least_squares_fitting_2(first,beyond,line,centroid,K(), t);
|
||||
}
|
||||
|
||||
// does not return the centroid and deduces the kernel as well.
|
||||
template < typename InputIterator,
|
||||
typename Line >
|
||||
typename Line, typename tag >
|
||||
inline
|
||||
typename Kernel_traits<Line>::Kernel::FT
|
||||
linear_least_squares_fitting_2(InputIterator first,
|
||||
InputIterator beyond,
|
||||
Line& line,
|
||||
bool non_standard_geometry = false)
|
||||
const tag& t)
|
||||
{
|
||||
typedef typename std::iterator_traits<InputIterator>::value_type Value_type;
|
||||
typedef typename Kernel_traits<Value_type>::Kernel K;
|
||||
return CGAL::linear_least_squares_fitting_2(first,beyond,line,K(), non_standard_geometry);
|
||||
// BOOST_STATIC_ASSERT((boost::is_same<typename CGAL::Algebraic_structure_traits<Value_type>::Algebraic_category,CGAL::Field_with_sqrt_tag>::value));
|
||||
return CGAL::linear_least_squares_fitting_2(first,beyond,line,K(), t);
|
||||
}
|
||||
|
||||
CGAL_END_NAMESPACE
|
||||
|
|
|
|||
|
|
@ -14,15 +14,22 @@
|
|||
// $URL$
|
||||
// $Id$
|
||||
//
|
||||
// Author(s) : Pierre Alliez and Sylvain Pion
|
||||
// Author(s) : Pierre Alliez and Sylvain Pion and Ankit Gupta
|
||||
|
||||
#ifndef CGAL_LINEAR_LEAST_SQUARES_FITTING_3_H
|
||||
#define CGAL_LINEAR_LEAST_SQUARES_FITTING_3_H
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Object.h>
|
||||
#include <CGAL/centroid.h>
|
||||
#include <CGAL/eigen.h>
|
||||
#include <CGAL/Algebraic_structure_traits.h>
|
||||
#include <CGAL/IO/io.h>
|
||||
#include <CGAL/linear_least_squares_fitting_points_3.h>
|
||||
#include <CGAL/linear_least_squares_fitting_segments_3.h>
|
||||
#include <CGAL/linear_least_squares_fitting_triangles_3.h>
|
||||
#include <CGAL/linear_least_squares_fitting_spheres_3.h>
|
||||
#include <CGAL/linear_least_squares_fitting_cuboids_3.h>
|
||||
#include <CGAL/linear_least_squares_fitting_tetrahedrons_3.h>
|
||||
#include <CGAL/PCA_tags.h>
|
||||
|
||||
#include <iterator>
|
||||
#include <list>
|
||||
|
|
@ -30,413 +37,106 @@
|
|||
|
||||
CGAL_BEGIN_NAMESPACE
|
||||
|
||||
namespace CGALi {
|
||||
|
||||
|
||||
// assemble covariance matrix from a point set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
void
|
||||
assemble_covariance_matrix_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
const K& , // kernel
|
||||
const typename K::Point_3*) // used for indirection
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
// 1 2
|
||||
// 3 4 5
|
||||
covariance[0] = covariance[1] = covariance[2] =
|
||||
covariance[3] = covariance[4] = covariance[5] = (FT)0.0;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Point& p = *it;
|
||||
Vector d = p - c;
|
||||
covariance[0] += d.x() * d.x();
|
||||
covariance[1] += d.x() * d.y();
|
||||
covariance[2] += d.y() * d.y();
|
||||
covariance[3] += d.x() * d.z();
|
||||
covariance[4] += d.y() * d.z();
|
||||
covariance[5] += d.z() * d.z();
|
||||
}
|
||||
}
|
||||
|
||||
// assemble covariance matrix from a triangle set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
void
|
||||
assemble_covariance_matrix_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
const K& , // kernel
|
||||
const typename K::Triangle_3*)// used for indirection
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
// 1 2
|
||||
// 3 4 5
|
||||
covariance[0] = covariance[1] = covariance[2] =
|
||||
covariance[3] = covariance[4] = covariance[5] = (FT)0.0;
|
||||
FT sum_areas = 0.0;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Triangle& triangle = *it;
|
||||
FT area = std::sqrt(triangle.squared_area());
|
||||
Point c_t = K().construct_centroid_3_object()(triangle[0],triangle[1],triangle[2]);
|
||||
sum_areas += area;
|
||||
|
||||
// e1 = ab, e2 = ac
|
||||
Vector e1 = Vector(triangle[0],triangle[1]);
|
||||
Vector e2 = Vector(triangle[0],triangle[2]);
|
||||
|
||||
FT c1 = (FT)(2.0 * area * 10.0 / 72.0);
|
||||
FT c2 = (FT)(2.0 * area * 7.0 / 72.0);
|
||||
|
||||
covariance[0] += c1*(e1[0]*e1[0] + e2[0]*e2[0]) + (FT)2.0*c2*e1[0]*e2[0];
|
||||
covariance[1] += c1*(e1[1]*e1[0] + e2[1]*e2[0]) + c2*(e1[1]*e2[0] + e1[0]*e2[1]);
|
||||
covariance[2] += c1*(e1[1]*e1[1] + e2[1]*e2[1]) + (FT)2.0*c2*e1[1]*e2[1];
|
||||
covariance[3] += c1*(e1[2]*e1[0] + e2[2]*e2[0]) + c2*(e1[2]*e2[0] + e1[0]*e2[2]);
|
||||
covariance[4] += c1*(e1[2]*e1[1] + e2[2]*e2[1]) + c2*(e1[2]*e2[1] + e1[1]*e2[2]);
|
||||
covariance[5] += c1*(e1[2]*e1[2] + e2[2]*e2[2]) + (FT)2.0*c2*e1[2]*e2[2];
|
||||
|
||||
// add area(t) c(t) * transpose(c(t))
|
||||
covariance[0] += area * c_t.x() * c_t.x();
|
||||
covariance[1] += area * c_t.y() * c_t.x();
|
||||
covariance[2] += area * c_t.y() * c_t.y();
|
||||
covariance[3] += area * c_t.z() * c_t.x();
|
||||
covariance[4] += area * c_t.z() * c_t.y();
|
||||
covariance[5] += area * c_t.z() * c_t.z();
|
||||
}
|
||||
|
||||
// remove sum(area) * (c * transpose(c))
|
||||
covariance[0] -= sum_areas * c.x() * c.x();
|
||||
covariance[1] -= sum_areas * c.y() * c.x();
|
||||
covariance[2] -= sum_areas * c.y() * c.y();
|
||||
covariance[3] -= sum_areas * c.z() * c.x();
|
||||
covariance[4] -= sum_areas * c.z() * c.y();
|
||||
covariance[5] -= sum_areas * c.z() * c.z();
|
||||
}
|
||||
|
||||
|
||||
// compute the eigen values and vectors of the covariance
|
||||
// matrix and deduces the best linear fitting plane.
|
||||
// returns fitting quality
|
||||
template < typename K >
|
||||
typename K::FT
|
||||
fitting_plane_3(const typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
const K& ) // kernel
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Plane_3 Plane;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
|
||||
// solve for eigenvalues and eigenvectors.
|
||||
// eigen values are sorted in descending order,
|
||||
// eigen vectors are sorted in accordance.
|
||||
FT eigen_values[3];
|
||||
FT eigen_vectors[9];
|
||||
eigen_symmetric<FT>(covariance,3,eigen_vectors,eigen_values);
|
||||
|
||||
// check unicity and build fitting line accordingly
|
||||
if(eigen_values[0] != eigen_values[1] &&
|
||||
eigen_values[0] != eigen_values[2])
|
||||
{
|
||||
// regular case
|
||||
Vector normal(eigen_vectors[6],
|
||||
eigen_vectors[7],
|
||||
eigen_vectors[8]);
|
||||
plane = Plane(c,normal);
|
||||
return (FT)1.0 - eigen_values[2] / eigen_values[0];
|
||||
} // end regular case
|
||||
else
|
||||
{
|
||||
// isotropic case (infinite number of directions)
|
||||
// by default: assemble a horizontal plane that goes
|
||||
// through the centroid.
|
||||
plane = Plane(c,Vector(0.0,0.0,1.0));
|
||||
return (FT)0.0;
|
||||
}
|
||||
}
|
||||
|
||||
// compute the eigen values and vectors of the covariance
|
||||
// matrix and deduces the best linear fitting line
|
||||
// (this is an internal function)
|
||||
// returns fitting quality
|
||||
template < typename K >
|
||||
typename K::FT
|
||||
fitting_line_3(const typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
typename K::Line_3& line, // best fit line
|
||||
const K& ) // kernel
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Line_3 Line;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
|
||||
// solve for eigenvalues and eigenvectors.
|
||||
// eigen values are sorted in descending order,
|
||||
// eigen vectors are sorted in accordance.
|
||||
FT eigen_values[3];
|
||||
FT eigen_vectors[9];
|
||||
eigen_symmetric<FT>(covariance,3,eigen_vectors,eigen_values);
|
||||
|
||||
// check unicity and build fitting line accordingly
|
||||
if(eigen_values[0] != eigen_values[1])
|
||||
{
|
||||
// regular case
|
||||
Vector direction(eigen_vectors[0],eigen_vectors[1],eigen_vectors[2]);
|
||||
line = Line(c,direction);
|
||||
return (FT)1.0 - eigen_values[1] / eigen_values[0];
|
||||
} // end regular case
|
||||
else
|
||||
{
|
||||
// isotropic case (infinite number of directions)
|
||||
// by default: assemble a horizontal plane that goes
|
||||
// through the centroid.
|
||||
line = Line(c,Vector(1.0,0.0,0.0));
|
||||
return (FT)0.0;
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
|
||||
// fits a plane to a 3D point set
|
||||
// returns a fitting quality (1 - lambda_min/lambda_max):
|
||||
// 1 is best (zero variance orthogonally to the fitting line)
|
||||
// 0 is worst (isotropic case, returns a plane with default direction)
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Point_3*) // used for indirection
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Plane_3 Plane;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K());
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6];
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Point*) NULL);
|
||||
|
||||
// compute fitting plane
|
||||
return fitting_plane_3(covariance,c,plane,k);
|
||||
} // end fit plane to point set
|
||||
|
||||
// fits a line to a 3D point set
|
||||
// returns a fitting quality (1 - lambda_min/lambda_max):
|
||||
// 1 is best (zero variance orthogonally to the fitting line)
|
||||
// 0 is worst (isotropic case, returns a line along x axis)
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Point_3*)// used for indirection
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Line_3 Line;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K());
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6];
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Point*) NULL);
|
||||
|
||||
// compute fitting line
|
||||
return fitting_line_3(covariance,c,line,k);
|
||||
} // end fit line to point set
|
||||
|
||||
// fits a plane to a 3D triangle set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Triangle_3*)// used for indirection
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K());
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6];
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Triangle*) NULL);
|
||||
|
||||
// compute fitting plane
|
||||
return fitting_plane_3(covariance,c,plane,k);
|
||||
} // end linear_least_squares_fitting_3
|
||||
|
||||
// fits a line to a 3D triangle set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Triangle_3*)// used for indirection
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
typedef typename K::Line_3 Line;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K());
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6];
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Triangle*) NULL);
|
||||
|
||||
// compute fitting plane
|
||||
return fitting_line_3(covariance,c,line,k);
|
||||
} // end linear_least_squares_fitting_3
|
||||
|
||||
} // end namespace CGALi
|
||||
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K , typename tag>
|
||||
inline
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane,
|
||||
typename K::Point_3& centroid,
|
||||
const K& k)
|
||||
const K& k,
|
||||
const tag& t)
|
||||
{
|
||||
typedef typename std::iterator_traits<InputIterator>::value_type Value_type;
|
||||
// BOOST_STATIC_ASSERT((boost::is_same<typename CGAL::Algebraic_structure_traits<Value_type>::Algebraic_category,CGAL::Field_with_sqrt_tag>::value));
|
||||
return CGALi::linear_least_squares_fitting_3(first, beyond, plane,
|
||||
centroid, k, (Value_type*) NULL);
|
||||
centroid, k, (Value_type*) NULL, t);
|
||||
}
|
||||
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K, typename tag >
|
||||
inline
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line,
|
||||
typename K::Point_3& centroid,
|
||||
const K& k)
|
||||
const K& k,
|
||||
const tag& t)
|
||||
{
|
||||
typedef typename std::iterator_traits<InputIterator>::value_type Value_type;
|
||||
// BOOST_STATIC_ASSERT((boost::is_same<typename CGAL::Algebraic_structure_traits<Value_type>::Algebraic_category,CGAL::Field_with_sqrt_tag>::value));
|
||||
return CGALi::linear_least_squares_fitting_3(first, beyond, line,
|
||||
centroid, k, (Value_type*) NULL);
|
||||
centroid, k, (Value_type*) NULL, t);
|
||||
}
|
||||
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K, typename tag >
|
||||
inline
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane,
|
||||
const K& k)
|
||||
const K& k,
|
||||
const tag& t)
|
||||
{
|
||||
typedef typename std::iterator_traits<InputIterator>::value_type Value_type;
|
||||
// BOOST_STATIC_ASSERT((boost::is_same<typename CGAL::Algebraic_structure_traits<Value_type>::Algebraic_category,CGAL::Field_with_sqrt_tag>::value));
|
||||
typename K::Point_3 centroid;
|
||||
return CGALi::linear_least_squares_fitting_3(first, beyond, plane,
|
||||
centroid, k,(Value_type*) NULL);
|
||||
centroid, k,(Value_type*) NULL, t);
|
||||
}
|
||||
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K, typename tag >
|
||||
inline
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line,
|
||||
const K& k)
|
||||
const K& k,
|
||||
const tag& t)
|
||||
{
|
||||
typedef typename std::iterator_traits<InputIterator>::value_type Value_type;
|
||||
// BOOST_STATIC_ASSERT((boost::is_same<typename CGAL::Algebraic_structure_traits<Value_type>::Algebraic_category,CGAL::Field_with_sqrt_tag>::value));
|
||||
typename K::Point_3 centroid;
|
||||
return CGALi::linear_least_squares_fitting_3(first, beyond, line,
|
||||
centroid, k,(Value_type*) NULL);
|
||||
centroid, k,(Value_type*) NULL, t);
|
||||
}
|
||||
|
||||
// deduces the kernel from the points in container.
|
||||
template < typename InputIterator,
|
||||
typename Object,
|
||||
typename Point>
|
||||
typename Point, typename tag>
|
||||
inline
|
||||
typename Kernel_traits<Object>::Kernel::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
Object& object,
|
||||
Point& centroid)
|
||||
Point& centroid,
|
||||
const tag& t)
|
||||
{
|
||||
typedef typename std::iterator_traits<InputIterator>::value_type Value_type;
|
||||
// BOOST_STATIC_ASSERT((boost::is_same<typename CGAL::Algebraic_structure_traits<Value_type>::Algebraic_category,CGAL::Field_with_sqrt_tag>::value));
|
||||
typedef typename Kernel_traits<Value_type>::Kernel K;
|
||||
return CGAL::linear_least_squares_fitting_3(first,beyond,object,centroid,K());
|
||||
return CGAL::linear_least_squares_fitting_3(first,beyond,object,centroid,K(), t);
|
||||
}
|
||||
|
||||
// does not return the centroid and deduces the kernel as well.
|
||||
template < typename InputIterator,
|
||||
typename Object >
|
||||
typename Object, typename tag >
|
||||
inline
|
||||
typename Kernel_traits<Object>::Kernel::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
Object& object)
|
||||
Object& object,
|
||||
const tag& t)
|
||||
{
|
||||
typedef typename std::iterator_traits<InputIterator>::value_type Value_type;
|
||||
// BOOST_STATIC_ASSERT((boost::is_same<typename CGAL::Algebraic_structure_traits<Value_type>::Algebraic_category,CGAL::Field_with_sqrt_tag>::value));
|
||||
typedef typename Kernel_traits<Value_type>::Kernel K;
|
||||
return CGAL::linear_least_squares_fitting_3(first,beyond,object,K());
|
||||
return CGAL::linear_least_squares_fitting_3(first,beyond,object,K(), t);
|
||||
}
|
||||
|
||||
CGAL_END_NAMESPACE
|
||||
|
|
|
|||
|
|
@ -19,8 +19,6 @@
|
|||
#ifndef CGAL_LINEAR_LEAST_SQUARES_FITTING_CIRCLES_2_H
|
||||
#define CGAL_LINEAR_LEAST_SQUARES_FITTING_CIRCLES_2_H
|
||||
|
||||
#define PI 3.1415926535
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Object.h>
|
||||
#include <CGAL/centroid.h>
|
||||
|
|
@ -50,7 +48,7 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
typename K::Point_2& c, // centroid
|
||||
const K&, // kernel
|
||||
const typename K::Circle_2*,// used for indirection
|
||||
const bool non_standard_geometry) // true means it is a disk
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
// types
|
||||
typedef typename K::FT FT;
|
||||
|
|
@ -63,98 +61,10 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
if(!non_standard_geometry) {
|
||||
|
||||
// ::::::::::DISK:::::::::::::::
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K());
|
||||
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
// 1 2
|
||||
//Final combined covariance matrix for all circles and their combined mass
|
||||
FT mass = 0.0;
|
||||
FT covariance[3] = {0.0,0.0,0.0};
|
||||
|
||||
// assemble 2nd order moment about the origin.
|
||||
FT temp[4] = {1.0, 0.0,
|
||||
0.0, 1.0};
|
||||
Matrix moment = init_Matrix<K>(2,temp);
|
||||
// Matrix moment = Matrix(2,true,PI);
|
||||
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
// Now for each circle, construct the 2nd order moment about the origin.
|
||||
// assemble the transformation matrix.
|
||||
const Circle& t = *it;
|
||||
|
||||
// defined for convenience.
|
||||
// FT example = CGAL::to_double(t[0].x());
|
||||
FT radius = std::sqrt(t.squared_radius());
|
||||
FT delta[4] = {radius, 0.0,
|
||||
0.0, radius};
|
||||
Matrix transformation = init_Matrix<K>(2,delta);
|
||||
FT area = t.squared_radius();
|
||||
CGAL_assertion(area != 0.0);
|
||||
|
||||
// Find the 2nd order moment for the circle wrt to the origin by an affine transformation.
|
||||
|
||||
// Transform the standard 2nd order moment using the transformation matrix
|
||||
transformation = area * transformation * moment * LA::transpose(transformation);
|
||||
|
||||
// Translate the 2nd order moment to the center of the circle.
|
||||
FT x0 = t.center().x();
|
||||
FT y0 = t.center().y();
|
||||
|
||||
// and add to covariance matrix
|
||||
covariance[0] += transformation[0][0] + area * x0*x0;
|
||||
covariance[1] += transformation[0][1] + area * x0*y0;
|
||||
covariance[2] += transformation[1][1] + area * y0*y0;
|
||||
|
||||
mass += area;
|
||||
}
|
||||
|
||||
// Translate the 2nd order moment calculated about the origin to
|
||||
// the center of mass to get the covariance.
|
||||
covariance[0] += mass * (-1.0 * c.x() * c.x());
|
||||
covariance[1] += mass * (-1.0 * c.x() * c.y());
|
||||
covariance[2] += mass * (-1.0 * c.y() * c.y());
|
||||
|
||||
// to remove later
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<std::endl;
|
||||
|
||||
// solve for eigenvalues and eigenvectors.
|
||||
// eigen values are sorted in descending order,
|
||||
// eigen vectors are sorted in accordance.
|
||||
std::pair<FT,FT> eigen_values;
|
||||
std::pair<Vector,Vector> eigen_vectors;
|
||||
// CGALi::eigen_symmetric_2<K>(final_cov, eigen_vectors, eigen_values);
|
||||
FT eigen_vectors1[4];
|
||||
FT eigen_values1[2];
|
||||
eigen_symmetric<FT>(covariance,2, eigen_vectors1, eigen_values1);
|
||||
eigen_values = std::make_pair(eigen_values1[0],eigen_values1[1]);
|
||||
eigen_vectors = std::make_pair(Vector(eigen_vectors1[0],eigen_vectors1[1]),Vector(eigen_vectors1[2],eigen_vectors1[3]));
|
||||
// check unicity and build fitting line accordingly
|
||||
if(eigen_values.first != eigen_values.second)
|
||||
{
|
||||
// regular case
|
||||
line = Line(c, eigen_vectors.first);
|
||||
return (FT)1.0 - eigen_values.second / eigen_values.first;
|
||||
}
|
||||
else
|
||||
{
|
||||
// isotropic case (infinite number of directions)
|
||||
// by default: assemble a line that goes through
|
||||
// the centroid and with a default horizontal vector.
|
||||
line = Line(c, Vector(1.0, 0.0));
|
||||
return (FT)0.0;
|
||||
}
|
||||
}
|
||||
else {
|
||||
// ::::::::::DISK:::::::::::::::
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K());
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
|
|
@ -239,8 +149,117 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
line = Line(c, Vector(1.0, 0.0));
|
||||
return (FT)0.0;
|
||||
}
|
||||
|
||||
} // end linear_least_squares_fitting_2 for circle set with 2D tag
|
||||
|
||||
template < typename InputIterator, typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_2(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_2& line, // best fit line
|
||||
typename K::Point_2& c, // centroid
|
||||
const K&, // kernel
|
||||
const typename K::Circle_2*,// used for indirection
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
// types
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Line_2 Line;
|
||||
typedef typename K::Point_2 Point;
|
||||
typedef typename K::Vector_2 Vector;
|
||||
typedef typename K::Circle_2 Circle;
|
||||
typedef typename CGAL::Linear_algebraCd<FT> LA;
|
||||
typedef typename LA::Matrix Matrix;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
// 1 2
|
||||
//Final combined covariance matrix for all circles and their combined mass
|
||||
FT mass = 0.0;
|
||||
FT covariance[3] = {0.0,0.0,0.0};
|
||||
|
||||
// assemble 2nd order moment about the origin.
|
||||
FT temp[4] = {1.0, 0.0,
|
||||
0.0, 1.0};
|
||||
Matrix moment = init_Matrix<K>(2,temp);
|
||||
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
// Now for each circle, construct the 2nd order moment about the origin.
|
||||
// assemble the transformation matrix.
|
||||
const Circle& t = *it;
|
||||
|
||||
// defined for convenience.
|
||||
// FT example = CGAL::to_double(t[0].x());
|
||||
FT radius = std::sqrt(t.squared_radius());
|
||||
FT delta[4] = {radius, 0.0,
|
||||
0.0, radius};
|
||||
Matrix transformation = init_Matrix<K>(2,delta);
|
||||
FT length = 2 * radius;
|
||||
CGAL_assertion(length != 0.0);
|
||||
|
||||
// Find the 2nd order moment for the circle wrt to the origin by an affine transformation.
|
||||
|
||||
// Transform the standard 2nd order moment using the transformation matrix
|
||||
transformation = 0.5 * length * transformation * moment * LA::transpose(transformation);
|
||||
|
||||
// Translate the 2nd order moment to the center of the circle.
|
||||
FT x0 = t.center().x();
|
||||
FT y0 = t.center().y();
|
||||
|
||||
// and add to covariance matrix
|
||||
covariance[0] += transformation[0][0] + length * x0*x0;
|
||||
covariance[1] += transformation[0][1] + length * x0*y0;
|
||||
covariance[2] += transformation[1][1] + length * y0*y0;
|
||||
|
||||
mass += length;
|
||||
}
|
||||
} // end linear_least_squares_fitting_2 for circle set
|
||||
|
||||
// Translate the 2nd order moment calculated about the origin to
|
||||
// the center of mass to get the covariance.
|
||||
covariance[0] += mass * (-1.0 * c.x() * c.x());
|
||||
covariance[1] += mass * (-1.0 * c.x() * c.y());
|
||||
covariance[2] += mass * (-1.0 * c.y() * c.y());
|
||||
|
||||
// to remove later
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<std::endl;
|
||||
|
||||
// solve for eigenvalues and eigenvectors.
|
||||
// eigen values are sorted in descending order,
|
||||
// eigen vectors are sorted in accordance.
|
||||
std::pair<FT,FT> eigen_values;
|
||||
std::pair<Vector,Vector> eigen_vectors;
|
||||
// CGALi::eigen_symmetric_2<K>(final_cov, eigen_vectors, eigen_values);
|
||||
FT eigen_vectors1[4];
|
||||
FT eigen_values1[2];
|
||||
eigen_symmetric<FT>(covariance,2, eigen_vectors1, eigen_values1);
|
||||
eigen_values = std::make_pair(eigen_values1[0],eigen_values1[1]);
|
||||
eigen_vectors = std::make_pair(Vector(eigen_vectors1[0],eigen_vectors1[1]),Vector(eigen_vectors1[2],eigen_vectors1[3]));
|
||||
// check unicity and build fitting line accordingly
|
||||
if(eigen_values.first != eigen_values.second)
|
||||
{
|
||||
// regular case
|
||||
line = Line(c, eigen_vectors.first);
|
||||
return (FT)1.0 - eigen_values.second / eigen_values.first;
|
||||
}
|
||||
else
|
||||
{
|
||||
// isotropic case (infinite number of directions)
|
||||
// by default: assemble a line that goes through
|
||||
// the centroid and with a default horizontal vector.
|
||||
line = Line(c, Vector(1.0, 0.0));
|
||||
return (FT)0.0;
|
||||
}
|
||||
} // end linear_least_squares_fitting_2 for circle set with 1D tag
|
||||
|
||||
} // end namespace CGALi
|
||||
|
||||
|
|
|
|||
|
|
@ -0,0 +1,330 @@
|
|||
// Copyright (c) 2005 INRIA Sophia-Antipolis (France).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org); you may redistribute it under
|
||||
// the terms of the Q Public License version 1.0.
|
||||
// See the file LICENSE.QPL distributed with CGAL.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL: svn+ssh://gankit@scm.gforge.inria.fr/svn/cgal/trunk/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_cuboids_3.h $
|
||||
// $Id: linear_least_squares_fitting_2.h 37882 2007-04-03 15:15:30Z spion $
|
||||
//
|
||||
// Author(s) : Pierre Alliez and Sylvain Pion and Ankit Gupta
|
||||
|
||||
#ifndef CGAL_LINEAR_LEAST_SQUARES_FITTING_CUBOIDS_3_H
|
||||
#define CGAL_LINEAR_LEAST_SQUARES_FITTING_CUBOIDS_3_H
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Object.h>
|
||||
#include <CGAL/centroid.h>
|
||||
#include <CGAL/eigen.h>
|
||||
#include <CGAL/util.h>
|
||||
|
||||
#include <iterator>
|
||||
|
||||
CGAL_BEGIN_NAMESPACE
|
||||
|
||||
namespace CGALi {
|
||||
|
||||
// fits a plane to a 3D cuboid set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Iso_cuboid_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_3_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Iso_cuboid*) NULL,tag);
|
||||
|
||||
// to remove later
|
||||
// std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
|
||||
// compute fitting plane
|
||||
return fitting_plane_3(covariance,c,plane,k);
|
||||
|
||||
|
||||
} // end linear_least_squares_fitting_cuboids_3
|
||||
|
||||
// fits a plane to a 3D cuboid set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Iso_cuboid_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Iso_cuboid*) NULL,tag);
|
||||
|
||||
// to remove later
|
||||
// std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
|
||||
// compute fitting plane
|
||||
return fitting_plane_3(covariance,c,plane,k);
|
||||
|
||||
|
||||
} // end linear_least_squares_fitting_cuboids_3
|
||||
|
||||
// fits a plane to a 3D cuboid set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Segment_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Iso_cuboid_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Segment> segments;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Iso_cuboid& t = *it;
|
||||
segments.push_back(t[0],t[1]);
|
||||
segments.push_back(t[1],t[2]);
|
||||
segments.push_back(t[2],t[3]);
|
||||
segments.push_back(t[3],t[0]);
|
||||
segments.push_back(t[4],t[5]);
|
||||
segments.push_back(t[5],t[6]);
|
||||
segments.push_back(t[6],t[7]);
|
||||
segments.push_back(t[7],t[4]);
|
||||
segments.push_back(t[5],t[0]);
|
||||
segments.push_back(t[1],t[6]);
|
||||
segments.push_back(t[2],t[7]);
|
||||
segments.push_back(t[3],t[4]);
|
||||
}
|
||||
|
||||
// compute fitting plane
|
||||
return linear_least_squares_fitting_3(segments.begin(),segments.end(),plane,c,k,(Segment*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_cuboids_3
|
||||
|
||||
// fits a plane to a 3D cuboid set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Iso_cuboid_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
typedef typename K::Point_3 Point;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Iso_cuboid& t = *it;
|
||||
points.push_back(t[0]);
|
||||
points.push_back(t[1]);
|
||||
points.push_back(t[2]);
|
||||
points.push_back(t[3]);
|
||||
}
|
||||
|
||||
// compute fitting plane
|
||||
return linear_least_squares_fitting_3(points.begin(),points.end(),plane,c,k,(Point*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_cuboids_3
|
||||
|
||||
// fits a line to a 3D cuboid set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Iso_cuboid_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_3_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Iso_cuboid*) NULL,tag);
|
||||
|
||||
// compute fitting line
|
||||
return fitting_line_3(covariance,c,line,k);
|
||||
|
||||
} // end linear_least_squares_fitting_cuboids_3
|
||||
|
||||
// fits a line to a 3D cuboid set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Iso_cuboid_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Iso_cuboid*) NULL,tag);
|
||||
|
||||
// to remove later
|
||||
// std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
|
||||
// compute fitting line
|
||||
return fitting_line_3(covariance,c,line,k);
|
||||
|
||||
|
||||
} // end linear_least_squares_fitting_cuboids_3
|
||||
|
||||
// fits a line to a 3D cuboid set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Segment_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Iso_cuboid_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Segment> segments;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Iso_cuboid& t = *it;
|
||||
segments.push_back(t[0],t[1]);
|
||||
segments.push_back(t[1],t[2]);
|
||||
segments.push_back(t[2],t[3]);
|
||||
segments.push_back(t[3],t[0]);
|
||||
segments.push_back(t[4],t[5]);
|
||||
segments.push_back(t[5],t[6]);
|
||||
segments.push_back(t[6],t[7]);
|
||||
segments.push_back(t[7],t[4]);
|
||||
segments.push_back(t[5],t[0]);
|
||||
segments.push_back(t[1],t[6]);
|
||||
segments.push_back(t[2],t[7]);
|
||||
segments.push_back(t[3],t[4]);
|
||||
}
|
||||
|
||||
// compute fitting line
|
||||
return linear_least_squares_fitting_3(segments.begin(),segments.end(),line,c,k,(Segment*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_cuboids_3
|
||||
|
||||
// fits a line to a 3D cuboid set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Iso_cuboid_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
typedef typename K::Point_3 Point;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Iso_cuboid& t = *it;
|
||||
points.push_back(t[0]);
|
||||
points.push_back(t[1]);
|
||||
points.push_back(t[2]);
|
||||
points.push_back(t[3]);
|
||||
}
|
||||
|
||||
// compute fitting line
|
||||
return linear_least_squares_fitting_3(points.begin(),points.end(),line,c,k,(Point*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_cuboids_3
|
||||
|
||||
} // end namespace CGALi
|
||||
|
||||
CGAL_END_NAMESPACE
|
||||
|
||||
#endif // CGAL_LINEAR_LEAST_SQUARES_FITTING_CUBOIDS_3_H
|
||||
|
|
@ -45,7 +45,7 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
typename K::Point_2& c, // centroid
|
||||
const K&, // kernel
|
||||
const typename K::Point_2*,// used for indirection
|
||||
const bool non_standard_geometry) // not useful
|
||||
const CGAL::PCA_dimension_0_tag& t)
|
||||
{
|
||||
// types
|
||||
typedef typename K::FT FT;
|
||||
|
|
@ -59,7 +59,7 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K());
|
||||
c = centroid(first,beyond,K(),t);
|
||||
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
|
|
|
|||
|
|
@ -0,0 +1,105 @@
|
|||
// Copyright (c) 2005 INRIA Sophia-Antipolis (France).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org); you may redistribute it under
|
||||
// the terms of the Q Public License version 1.0.
|
||||
// See the file LICENSE.QPL distributed with CGAL.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL: svn+ssh://gankit@scm.gforge.inria.fr/svn/cgal/trunk/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_points_3.h $
|
||||
// $Id: linear_least_squares_fitting_2.h 37882 2007-04-03 15:15:30Z spion $
|
||||
//
|
||||
// Author(s) : Pierre Alliez and Sylvain Pion and Ankit Gupta
|
||||
|
||||
#ifndef CGAL_LINEAR_LEAST_SQUARES_FITTING_POINTS_3_H
|
||||
#define CGAL_LINEAR_LEAST_SQUARES_FITTING_POINTS_3_H
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Object.h>
|
||||
#include <CGAL/centroid.h>
|
||||
#include <CGAL/eigen.h>
|
||||
#include <CGAL/util.h>
|
||||
|
||||
#include <iterator>
|
||||
|
||||
CGAL_BEGIN_NAMESPACE
|
||||
|
||||
namespace CGALi {
|
||||
// fits a plane to a 3D point set
|
||||
// returns a fitting quality (1 - lambda_min/lambda_max):
|
||||
// 1 is best (zero variance orthogonally to the fitting line)
|
||||
// 0 is worst (isotropic case, returns a plane with default direction)
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Point_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Plane_3 Plane;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6];
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Point*) NULL,tag);
|
||||
|
||||
// compute fitting plane
|
||||
return fitting_plane_3(covariance,c,plane,k);
|
||||
} // end fit plane to point set
|
||||
|
||||
// fits a line to a 3D point set
|
||||
// returns a fitting quality (1 - lambda_min/lambda_max):
|
||||
// 1 is best (zero variance orthogonally to the fitting line)
|
||||
// 0 is worst (isotropic case, returns a line along x axis)
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Point_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Line_3 Line;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6];
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Point*) NULL,tag);
|
||||
|
||||
// compute fitting line
|
||||
return fitting_line_3(covariance,c,line,k);
|
||||
} // end fit line to point set
|
||||
|
||||
} // end namespace CGALi
|
||||
|
||||
CGAL_END_NAMESPACE
|
||||
|
||||
#endif // CGAL_LINEAR_LEAST_SQUARES_FITTING_POINTS_3_H
|
||||
|
|
@ -48,7 +48,7 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
typename K::Point_2& c, // centroid
|
||||
const K&, // kernel
|
||||
const typename K::Iso_rectangle_2*,// used for indirection
|
||||
const bool non_standard_geometry) // true means it is an empty rectangle
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
// types
|
||||
typedef typename K::FT FT;
|
||||
|
|
@ -62,9 +62,9 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
if(!non_standard_geometry) {
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K());
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
|
|
@ -125,7 +125,7 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
covariance[2] += mass * (-1.0 * c.y() * c.y());
|
||||
|
||||
// to remove later
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<std::endl;
|
||||
// std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<std::endl;
|
||||
|
||||
// solve for eigenvalues and eigenvectors.
|
||||
// eigen values are sorted in descending order,
|
||||
|
|
@ -153,24 +153,75 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
line = Line(c, Vector(1.0, 0.0));
|
||||
return (FT)0.0;
|
||||
}
|
||||
}
|
||||
else {
|
||||
std::list<Segment_2> segments;
|
||||
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Iso_rectangle& t = *it;
|
||||
segments.push_back(Segment_2(t[0],t[1]));
|
||||
segments.push_back(Segment_2(t[1],t[2]));
|
||||
segments.push_back(Segment_2(t[2],t[3]));
|
||||
segments.push_back(Segment_2(t[3],t[0]));
|
||||
}
|
||||
} // end linear_least_squares_fitting_2 for rectangle set with 2D tag
|
||||
|
||||
return linear_least_squares_fitting_2(segments.begin(),segments.end(),line,c,K());
|
||||
}
|
||||
} // end linear_least_squares_fitting_2 for rectangle set
|
||||
template < typename InputIterator, typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_2(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_2& line, // best fit line
|
||||
typename K::Point_2& c, // centroid
|
||||
const K&, // kernel
|
||||
const typename K::Iso_rectangle_2*,// used for indirection
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
// types
|
||||
typedef typename K::Iso_rectangle_2 Iso_rectangle;
|
||||
typedef typename K::Segment_2 Segment_2;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Segment_2> segments;
|
||||
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Iso_rectangle& t = *it;
|
||||
segments.push_back(Segment_2(t[0],t[1]));
|
||||
segments.push_back(Segment_2(t[1],t[2]));
|
||||
segments.push_back(Segment_2(t[2],t[3]));
|
||||
segments.push_back(Segment_2(t[3],t[0]));
|
||||
}
|
||||
|
||||
return linear_least_squares_fitting_2(segments.begin(),segments.end(),line,c,K(),tag);
|
||||
|
||||
} // end linear_least_squares_fitting_2 for rectangle set with 1D tag
|
||||
|
||||
template < typename InputIterator, typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_2(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_2& line, // best fit line
|
||||
typename K::Point_2& c, // centroid
|
||||
const K&, // kernel
|
||||
const typename K::Iso_rectangle_2*,// used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
// types
|
||||
typedef typename K::Iso_rectangle_2 Iso_rectangle;
|
||||
typedef typename K::Point_2 Point_2;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Point_2> points;
|
||||
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Iso_rectangle& t = *it;
|
||||
points.push_back(Point_2(t[0]));
|
||||
points.push_back(Point_2(t[1]));
|
||||
points.push_back(Point_2(t[2]));
|
||||
points.push_back(Point_2(t[3]));
|
||||
}
|
||||
|
||||
return linear_least_squares_fitting_2(points.begin(),points.end(),line,c,K(),tag);
|
||||
|
||||
} // end linear_least_squares_fitting_2 for rectangle set with 0D tag
|
||||
|
||||
} // end namespace CGALi
|
||||
|
||||
|
|
|
|||
|
|
@ -48,7 +48,7 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
typename K::Point_2& c, // centroid
|
||||
const K&, // kernel
|
||||
const typename K::Segment_2*,// used for indirection
|
||||
const bool non_standard_geometry) // not useful
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
// types
|
||||
typedef typename K::FT FT;
|
||||
|
|
@ -63,7 +63,7 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K());
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
|
|
@ -111,6 +111,7 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
covariance[0] += mass * (-1.0 * c.x() * c.x());
|
||||
covariance[1] += mass * (-1.0 * c.x() * c.y());
|
||||
covariance[2] += mass * (-1.0 * c.y() * c.y());
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<(std::sqrt(2)/3.0)<<std::endl;
|
||||
|
||||
// to remove later
|
||||
// std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<(std::sqrt(2)/3.0)<<std::endl;
|
||||
|
|
@ -141,7 +142,37 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
line = Line(c, Vector(1.0, 0.0));
|
||||
return (FT)0.0;
|
||||
}
|
||||
} // end linear_least_squares_fitting_2 for segment set
|
||||
} // end linear_least_squares_fitting_2 for segment set with 1D tag
|
||||
|
||||
template < typename InputIterator, typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_2(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_2& line, // best fit line
|
||||
typename K::Point_2& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Segment_2*,// used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
// types
|
||||
typedef typename K::Point_2 Point;
|
||||
typedef typename K::Segment_2 Segment;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Segment& s = *it;
|
||||
points.push_back(s[0]);
|
||||
points.push_back(s[1]);
|
||||
}
|
||||
return linear_least_squares_fitting_2(points.begin(),points.end(),line,c,k,(Point*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_2 for segment set with 1D tag
|
||||
|
||||
} // end namespace CGALi
|
||||
|
||||
|
|
|
|||
|
|
@ -0,0 +1,172 @@
|
|||
// Copyright (c) 2005 INRIA Sophia-Antipolis (France).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org); you may redistribute it under
|
||||
// the terms of the Q Public License version 1.0.
|
||||
// See the file LICENSE.QPL distributed with CGAL.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL: svn+ssh://gankit@scm.gforge.inria.fr/svn/cgal/trunk/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_segments_3.h $
|
||||
// $Id: linear_least_squares_fitting_2.h 37882 2007-04-03 15:15:30Z spion $
|
||||
//
|
||||
// Author(s) : Pierre Alliez and Sylvain Pion and Ankit Gupta
|
||||
|
||||
#ifndef CGAL_LINEAR_LEAST_SQUARES_FITTING_SEGMENTS_3_H
|
||||
#define CGAL_LINEAR_LEAST_SQUARES_FITTING_SEGMENTS_3_H
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Object.h>
|
||||
#include <CGAL/centroid.h>
|
||||
#include <CGAL/eigen.h>
|
||||
#include <CGAL/util.h>
|
||||
|
||||
#include <iterator>
|
||||
|
||||
CGAL_BEGIN_NAMESPACE
|
||||
|
||||
namespace CGALi {
|
||||
|
||||
// fits a plane to a 3D segment set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Segment_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Segment*) NULL,tag);
|
||||
|
||||
// to remove later
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
|
||||
// compute fitting plane
|
||||
return fitting_plane_3(covariance,c,plane,k);
|
||||
|
||||
} // end linear_least_squares_fitting_segments_3
|
||||
|
||||
// fits a plane to a 3D segment set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Segment_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
typedef typename K::Point_3 Point;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Segment& t = *it;
|
||||
points.push_back(t[0]);
|
||||
points.push_back(t[1]);
|
||||
}
|
||||
|
||||
// compute fitting plane
|
||||
return linear_least_squares_fitting_3(points.begin(),points.end(),plane,c,k,(Point*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_segments_3
|
||||
|
||||
// fits a line to a 3D segment set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Segment_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Segment*) NULL,tag);
|
||||
|
||||
// to remove later
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
|
||||
// compute fitting line
|
||||
return fitting_line_3(covariance,c,line,k);
|
||||
|
||||
} // end linear_least_squares_fitting_segments_3
|
||||
|
||||
// fits a plane to a 3D segment set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Segment_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
typedef typename K::Point_3 Point;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Segment& t = *it;
|
||||
points.push_back(t[0]);
|
||||
points.push_back(t[1]);
|
||||
}
|
||||
|
||||
// compute fitting plane
|
||||
return linear_least_squares_fitting_3(points.begin(),points.end(),line,c,k,(Point*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_segments_3
|
||||
|
||||
} // end namespace CGALi
|
||||
|
||||
CGAL_END_NAMESPACE
|
||||
|
||||
#endif // CGAL_LINEAR_LEAST_SQUARES_FITTING_SEGMENTS_3_H
|
||||
|
|
@ -0,0 +1,164 @@
|
|||
// Copyright (c) 2005 INRIA Sophia-Antipolis (France).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org); you may redistribute it under
|
||||
// the terms of the Q Public License version 1.0.
|
||||
// See the file LICENSE.QPL distributed with CGAL.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL: svn+ssh://gankit@scm.gforge.inria.fr/svn/cgal/trunk/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_spheres_3.h $
|
||||
// $Id: linear_least_squares_fitting_2.h 37882 2007-04-03 15:15:30Z spion $
|
||||
//
|
||||
// Author(s) : Pierre Alliez and Sylvain Pion and Ankit Gupta
|
||||
|
||||
#ifndef CGAL_LINEAR_LEAST_SQUARES_FITTING_SPHERES_3_H
|
||||
#define CGAL_LINEAR_LEAST_SQUARES_FITTING_SPHERES_3_H
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Object.h>
|
||||
#include <CGAL/centroid.h>
|
||||
#include <CGAL/eigen.h>
|
||||
#include <CGAL/util.h>
|
||||
|
||||
#include <iterator>
|
||||
|
||||
CGAL_BEGIN_NAMESPACE
|
||||
|
||||
namespace CGALi {
|
||||
|
||||
// fits a plane to a 3D sphere set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Sphere_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_3_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Sphere_3 Sphere;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Sphere*) NULL,tag);
|
||||
|
||||
// compute fitting plane
|
||||
return fitting_plane_3(covariance,c,plane,k);
|
||||
|
||||
} // end linear_least_squares_fitting_spheres_3
|
||||
|
||||
// fits a plane to a 3D sphere set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Sphere_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Sphere_3 Sphere;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Sphere*) NULL,tag);
|
||||
|
||||
// compute fitting plane
|
||||
return fitting_plane_3(covariance,c,plane,k);
|
||||
|
||||
} // end linear_least_squares_fitting_spheres_3
|
||||
|
||||
// fits a line to a 3D sphere set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Sphere_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_3_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Sphere_3 Sphere;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Sphere*) NULL,tag);
|
||||
|
||||
// to remove later
|
||||
// std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
|
||||
// compute fitting line
|
||||
return fitting_line_3(covariance,c,line,k);
|
||||
|
||||
} // end linear_least_squares_fitting_spheres_3
|
||||
|
||||
// fits a line to a 3D sphere set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Sphere_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Sphere_3 Sphere;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Sphere*) NULL,tag);
|
||||
|
||||
// to remove later
|
||||
// std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
|
||||
// compute fitting line
|
||||
return fitting_line_3(covariance,c,line,k);
|
||||
|
||||
} // end linear_least_squares_fitting_spheres_3
|
||||
|
||||
} // end namespace CGALi
|
||||
|
||||
CGAL_END_NAMESPACE
|
||||
|
||||
#endif // CGAL_LINEAR_LEAST_SQUARES_FITTING_SPHERES_3_H
|
||||
|
|
@ -0,0 +1,325 @@
|
|||
// Copyright (c) 2005 INRIA Sophia-Antipolis (France).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (woidww.cgal.org); you may redistribute it under
|
||||
// the terms of the Q Public License version 1.0.
|
||||
// See the file LICENSE.QPL distributed with CGAL.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL: svn+ssh://gankit@scm.gforge.inria.fr/svn/cgal/trunk/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_tetrahedrons_3.h $
|
||||
// $Id: linear_least_squares_fitting_2.h 37882 2007-04-03 15:15:30Z spion $
|
||||
//
|
||||
// Author(s) : Pierre Alliez and Sylvain Pion and Ankit Gupta
|
||||
|
||||
#ifndef CGAL_LINEAR_LEAST_SQUARES_FITTING_TETRAHEDRONS_3_H
|
||||
#define CGAL_LINEAR_LEAST_SQUARES_FITTING_TETRAHEDRONS_3_H
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Object.h>
|
||||
#include <CGAL/centroid.h>
|
||||
#include <CGAL/eigen.h>
|
||||
#include <CGAL/util.h>
|
||||
|
||||
#include <iterator>
|
||||
|
||||
CGAL_BEGIN_NAMESPACE
|
||||
|
||||
namespace CGALi {
|
||||
|
||||
// fits a plane to a 3D tetrahedron set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Tetrahedron_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_3_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Tetrahedron_3 Tetrahedron;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Tetrahedron*) NULL,tag);
|
||||
|
||||
// to remove later
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
|
||||
// compute fitting plane
|
||||
return fitting_plane_3(covariance,c,plane,k);
|
||||
|
||||
|
||||
} // end linear_least_squares_fitting_tetrahedrons_3
|
||||
|
||||
// fits a plane to a 3D tetrahedron set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Tetrahedron_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Tetrahedron_3 Tetrahedron;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Triangle> triangles;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Tetrahedron& t = *it;
|
||||
triangles.push_back(t[0],t[1],t[2]);
|
||||
triangles.push_back(t[0],t[2],t[3]);
|
||||
triangles.push_back(t[0],t[3],t[1]);
|
||||
triangles.push_back(t[3],t[1],t[2]);
|
||||
}
|
||||
|
||||
// compute fitting plane
|
||||
return linear_least_squares_fitting_3(triangles.begin(),triangles.end(),plane,c,k,(Triangle*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_tetrahedrons_3
|
||||
|
||||
// fits a plane to a 3D tetrahedron set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Tetrahedron_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Tetrahedron_3 Tetrahedron;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Segment> segments;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Tetrahedron& t = *it;
|
||||
segments.push_back(t[0],t[1]);
|
||||
segments.push_back(t[1],t[2]);
|
||||
segments.push_back(t[1],t[3]);
|
||||
segments.push_back(t[2],t[3]);
|
||||
segments.push_back(t[0],t[2]);
|
||||
segments.push_back(t[0],t[3]);
|
||||
}
|
||||
|
||||
// compute fitting plane
|
||||
return linear_least_squares_fitting_3(segments.begin(),segments.end(),plane,c,k,(Segment*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_tetrahedrons_3
|
||||
|
||||
// fits a plane to a 3D tetrahedron set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Tetrahedron_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Tetrahedron_3 Tetrahedron;
|
||||
typedef typename K::Point_3 Point;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Tetrahedron& t = *it;
|
||||
points.push_back(t[0]);
|
||||
points.push_back(t[1]);
|
||||
points.push_back(t[2]);
|
||||
points.push_back(t[3]);
|
||||
}
|
||||
|
||||
// compute fitting plane
|
||||
return linear_least_squares_fitting_3(points.begin(),points.end(),plane,c,k,(Point*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_tetrahedrons_3
|
||||
|
||||
// fits a line to a 3D tetrahedron set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Tetrahedron_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_3_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Tetrahedron_3 Tetrahedron;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Tetrahedron*) NULL,tag);
|
||||
|
||||
// to remove later
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
|
||||
// compute fitting line
|
||||
return fitting_line_3(covariance,c,line,k);
|
||||
|
||||
} // end linear_least_squares_fitting_tetrahedrons_3
|
||||
|
||||
// fits a line to a 3D tetrahedron set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Tetrahedron_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Tetrahedron_3 Tetrahedron;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Triangle> triangles;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Tetrahedron& t = *it;
|
||||
triangles.push_back(t[0],t[1],t[2]);
|
||||
triangles.push_back(t[0],t[2],t[3]);
|
||||
triangles.push_back(t[0],t[3],t[1]);
|
||||
triangles.push_back(t[3],t[1],t[2]);
|
||||
}
|
||||
|
||||
// compute fitting line
|
||||
return linear_least_squares_fitting_3(triangles.begin(),triangles.end(),line,c,k,(Triangle*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_tetrahedrons_3
|
||||
|
||||
// fits a line to a 3D tetrahedron set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Tetrahedron_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Tetrahedron_3 Tetrahedron;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Segment> segments;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Tetrahedron& t = *it;
|
||||
segments.push_back(t[0],t[1]);
|
||||
segments.push_back(t[1],t[2]);
|
||||
segments.push_back(t[1],t[3]);
|
||||
segments.push_back(t[2],t[3]);
|
||||
segments.push_back(t[0],t[2]);
|
||||
segments.push_back(t[0],t[3]);
|
||||
}
|
||||
|
||||
// compute fitting line
|
||||
return linear_least_squares_fitting_3(segments.begin(),segments.end(),line,c,k,(Segment*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_tetrahedrons_3
|
||||
|
||||
// fits a line to a 3D tetrahedron set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Tetrahedron_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Tetrahedron_3 Tetrahedron;
|
||||
typedef typename K::Point_3 Point;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Tetrahedron& t = *it;
|
||||
points.push_back(t[0]);
|
||||
points.push_back(t[1]);
|
||||
points.push_back(t[2]);
|
||||
points.push_back(t[3]);
|
||||
}
|
||||
|
||||
// compute fitting line
|
||||
return linear_least_squares_fitting_3(points.begin(),points.end(),line,c,k,(Point*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_tetrahedrons_3
|
||||
|
||||
} // end namespace CGALi
|
||||
|
||||
CGAL_END_NAMESPACE
|
||||
|
||||
#endif // CGAL_LINEAR_LEAST_SQUARES_FITTING_TETRAHEDRONS_3_H
|
||||
|
|
@ -49,7 +49,7 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
typename K::Point_2& c, // centroid
|
||||
const K&, // kernel
|
||||
const typename K::Triangle_2*,// used for indirection
|
||||
const bool non_standard_geometry) // true means it is a hollow triangle
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
// types
|
||||
typedef typename K::FT FT;
|
||||
|
|
@ -64,9 +64,8 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
if(!non_standard_geometry) {
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K());
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
|
|
@ -97,24 +96,24 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
t[1].y() - y0, t[2].y() - y0};
|
||||
|
||||
Matrix transformation = init_Matrix<K>(2,delta);
|
||||
FT detA = fabs(LA::determinant(transformation));
|
||||
CGAL_assertion(detA!=0);
|
||||
FT area = 0.5 * std::abs(LA::determinant(transformation));
|
||||
CGAL_assertion(area!=0);
|
||||
|
||||
// Find the 2nd order moment for the triangle wrt to the origin by an affine transformation.
|
||||
|
||||
// Transform the standard 2nd order moment using the transformation matrix
|
||||
transformation = detA * transformation * moment * LA::transpose(transformation);
|
||||
transformation = 2 * area * transformation * moment * LA::transpose(transformation);
|
||||
|
||||
// Translate the 2nd order moment to (x0,y0).
|
||||
FT xav0 = (delta[0]+delta[1])/3.0;
|
||||
FT yav0 = (delta[2]+delta[3])/3.0;
|
||||
|
||||
// and add to the covariance matrix
|
||||
covariance[0] += transformation[0][0] + 0.5 * detA * (x0*xav0*2 + x0*x0);
|
||||
covariance[1] += transformation[0][1] + 0.5 * detA * (x0*yav0 + xav0*y0 + x0*y0);
|
||||
covariance[2] += transformation[1][1] + 0.5 * detA * (y0*yav0*2 + y0*y0);
|
||||
covariance[0] += transformation[0][0] + area * (x0*xav0*2 + x0*x0);
|
||||
covariance[1] += transformation[0][1] + area * (x0*yav0 + xav0*y0 + x0*y0);
|
||||
covariance[2] += transformation[1][1] + area * (y0*yav0*2 + y0*y0);
|
||||
|
||||
mass += 0.5 *detA;
|
||||
mass += area;
|
||||
}
|
||||
|
||||
// Translate the 2nd order moment calculated about the origin to
|
||||
|
|
@ -122,6 +121,8 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
covariance[0] += mass * (-1.0 * c.x() * c.x());
|
||||
covariance[1] += mass * (-1.0 * c.x() * c.y());
|
||||
covariance[2] += mass * (-1.0 * c.y() * c.y());
|
||||
// to remove later
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<std::endl;
|
||||
|
||||
// solve for eigenvalues and eigenvectors.
|
||||
// eigen values are sorted in descending order,
|
||||
|
|
@ -149,23 +150,73 @@ linear_least_squares_fitting_2(InputIterator first,
|
|||
line = Line(c, Vector(1.0, 0.0));
|
||||
return (FT)0.0;
|
||||
}
|
||||
}
|
||||
else {
|
||||
std::list<Segment_2> segments;
|
||||
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Triangle& t = *it;
|
||||
segments.push_back(Segment_2(t[0],t[1]));
|
||||
segments.push_back(Segment_2(t[1],t[2]));
|
||||
segments.push_back(Segment_2(t[2],t[0]));
|
||||
}
|
||||
} // end linear_least_squares_fitting_2 for triangle set with 2D tag
|
||||
|
||||
return linear_least_squares_fitting_2(segments.begin(),segments.end(),line,c,K());
|
||||
}
|
||||
} // end linear_least_squares_fitting_2 for triangle set
|
||||
template < typename InputIterator, typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_2(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_2& line, // best fit line
|
||||
typename K::Point_2& c, // centroid
|
||||
const K&, // kernel
|
||||
const typename K::Triangle_2*,// used for indirection
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
// types
|
||||
|
||||
typedef typename K::Triangle_2 Triangle;
|
||||
typedef typename K::Segment_2 Segment_2;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Segment_2> segments;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Triangle& t = *it;
|
||||
segments.push_back(Segment_2(t[0],t[1]));
|
||||
segments.push_back(Segment_2(t[1],t[2]));
|
||||
segments.push_back(Segment_2(t[2],t[0]));
|
||||
}
|
||||
|
||||
return linear_least_squares_fitting_2(segments.begin(),segments.end(),line,c,K(),tag);
|
||||
|
||||
} // end linear_least_squares_fitting_2 for triangle set with 1D tag
|
||||
|
||||
template < typename InputIterator, typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_2(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_2& line, // best fit line
|
||||
typename K::Point_2& c, // centroid
|
||||
const K&, // kernel
|
||||
const typename K::Triangle_2*,// used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
// types
|
||||
|
||||
typedef typename K::Triangle_2 Triangle;
|
||||
typedef typename K::Point_2 Point_2;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Point_2> points;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Triangle& t = *it;
|
||||
points.push_back(Point_2(t[0]));
|
||||
points.push_back(Point_2(t[1]));
|
||||
points.push_back(Point_2(t[2]));
|
||||
}
|
||||
|
||||
return linear_least_squares_fitting_2(points.begin(),points.end(),line,c,K(),tag);
|
||||
|
||||
} // end linear_least_squares_fitting_2 for triangle set with 0D tag
|
||||
|
||||
} // end namespace CGALi
|
||||
|
||||
|
|
|
|||
|
|
@ -0,0 +1,244 @@
|
|||
// Copyright (c) 2005 INRIA Sophia-Antipolis (France).
|
||||
// All rights reserved.
|
||||
//
|
||||
// This file is part of CGAL (www.cgal.org); you may redistribute it under
|
||||
// the terms of the Q Public License version 1.0.
|
||||
// See the file LICENSE.QPL distributed with CGAL.
|
||||
//
|
||||
// Licensees holding a valid commercial license may use this file in
|
||||
// accordance with the commercial license agreement provided with the software.
|
||||
//
|
||||
// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
|
||||
// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
|
||||
//
|
||||
// $URL: svn+ssh://gankit@scm.gforge.inria.fr/svn/cgal/trunk/Principal_component_analysis/include/CGAL/linear_least_squares_fitting_triangles_3.h $
|
||||
// $Id: linear_least_squares_fitting_2.h 37882 2007-04-03 15:15:30Z spion $
|
||||
//
|
||||
// Author(s) : Pierre Alliez and Sylvain Pion and Ankit Gupta
|
||||
|
||||
#ifndef CGAL_LINEAR_LEAST_SQUARES_FITTING_TRIANGLES_3_H
|
||||
#define CGAL_LINEAR_LEAST_SQUARES_FITTING_TRIANGLES_3_H
|
||||
|
||||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Object.h>
|
||||
#include <CGAL/centroid.h>
|
||||
#include <CGAL/eigen.h>
|
||||
#include <CGAL/util.h>
|
||||
|
||||
#include <iterator>
|
||||
|
||||
CGAL_BEGIN_NAMESPACE
|
||||
|
||||
namespace CGALi {
|
||||
|
||||
// fits a plane to a 3D triangle set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Triangle_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Triangle*) NULL,tag);
|
||||
|
||||
// to remove later
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
|
||||
// compute fitting plane
|
||||
return fitting_plane_3(covariance,c,plane,k);
|
||||
|
||||
} // end linear_least_squares_fitting_triangles_3
|
||||
|
||||
// fits a plane to a 3D triangle set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Segment_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Triangle_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Segment> segments;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Triangle& t = *it;
|
||||
segments.push_back(t[0],t[1]);
|
||||
segments.push_back(t[1],t[2]);
|
||||
segments.push_back(t[2],t[0]);
|
||||
}
|
||||
|
||||
// compute fitting plane
|
||||
return linear_least_squares_fitting_3(segments.begin(),segments.end(),plane,c,k,(Segment*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_triangles_3
|
||||
|
||||
// fits a plane to a 3D triangle set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Triangle_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
typedef typename K::Point_3 Point;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Triangle& t = *it;
|
||||
points.push_back(t[0]);
|
||||
points.push_back(t[1]);
|
||||
points.push_back(t[2]);
|
||||
}
|
||||
|
||||
// compute fitting plane
|
||||
return linear_least_squares_fitting_3(points.begin(),points.end(),plane,c,k,(Point*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_triangles_3
|
||||
|
||||
// fits a line to a 3D triangle set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Triangle_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
// compute centroid
|
||||
c = centroid(first,beyond,K(),tag);
|
||||
|
||||
// assemble covariance matrix
|
||||
FT covariance[6] = {0.0,0.0,0.0,0.0,0.0,0.0};
|
||||
assemble_covariance_matrix_3(first,beyond,covariance,c,k,(Triangle*) NULL,tag);
|
||||
|
||||
// to remove later
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
|
||||
// compute fitting line
|
||||
return fitting_line_3(covariance,c,line,k);
|
||||
|
||||
} // end linear_least_squares_fitting_triangles_3
|
||||
|
||||
// fits a line to a 3D triangle set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Segment_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Triangle_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Segment> segments;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Triangle& t = *it;
|
||||
segments.push_back(t[0],t[1]);
|
||||
segments.push_back(t[1],t[2]);
|
||||
segments.push_back(t[2],t[0]);
|
||||
}
|
||||
|
||||
// compute fitting line
|
||||
return linear_least_squares_fitting_3(segments.begin(),segments.end(),line,c,k,(Segment*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_triangles_3
|
||||
|
||||
// fits a line to a 3D triangle set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
typename K::FT
|
||||
linear_least_squares_fitting_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::Line_3& line, // best fit line
|
||||
typename K::Point_3& c, // centroid
|
||||
const K& k, // kernel
|
||||
const typename K::Triangle_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
typedef typename K::Point_3 Point;
|
||||
|
||||
// precondition: at least one element in the container.
|
||||
CGAL_precondition(first != beyond);
|
||||
|
||||
std::list<Point> points;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Triangle& t = *it;
|
||||
points.push_back(t[0]);
|
||||
points.push_back(t[1]);
|
||||
points.push_back(t[2]);
|
||||
}
|
||||
|
||||
// compute fitting line
|
||||
return linear_least_squares_fitting_3(points.begin(),points.end(),line,c,k,(Point*)NULL,tag);
|
||||
|
||||
} // end linear_least_squares_fitting_triangles_3
|
||||
|
||||
} // end namespace CGALi
|
||||
|
||||
CGAL_END_NAMESPACE
|
||||
|
||||
#endif // CGAL_LINEAR_LEAST_SQUARES_FITTING_TRIANGLES_3_H
|
||||
|
|
@ -22,6 +22,7 @@
|
|||
#include <CGAL/basic.h>
|
||||
#include <CGAL/Object.h>
|
||||
#include <CGAL/Linear_algebraCd.h>
|
||||
#include <CGAL/PCA_tags.h>
|
||||
|
||||
CGAL_BEGIN_NAMESPACE
|
||||
|
||||
|
|
@ -42,6 +43,685 @@ init_Matrix(int n, typename K::FT entries[]) {
|
|||
return ret;
|
||||
} // end initialization of matrix
|
||||
|
||||
// assemble covariance matrix from a point set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
void
|
||||
assemble_covariance_matrix_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
const K& , // kernel
|
||||
const typename K::Point_3*, // used for indirection
|
||||
const CGAL::PCA_dimension_0_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
// 1 2
|
||||
// 3 4 5
|
||||
covariance[0] = covariance[1] = covariance[2] =
|
||||
covariance[3] = covariance[4] = covariance[5] = (FT)0.0;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Point& p = *it;
|
||||
Vector d = p - c;
|
||||
covariance[0] += d.x() * d.x();
|
||||
covariance[1] += d.x() * d.y();
|
||||
covariance[2] += d.y() * d.y();
|
||||
covariance[3] += d.x() * d.z();
|
||||
covariance[4] += d.y() * d.z();
|
||||
covariance[5] += d.z() * d.z();
|
||||
}
|
||||
}
|
||||
|
||||
// assemble covariance matrix from a triangle set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
void
|
||||
assemble_covariance_matrix_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
const K& , // kernel
|
||||
const typename K::Triangle_3*,// used for indirection
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
typedef typename CGAL::Linear_algebraCd<FT> LA;
|
||||
typedef typename LA::Matrix Matrix;
|
||||
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
// 1 2
|
||||
// 3 4 5
|
||||
//Final combined covariance matrix for all triangles and their combined mass
|
||||
FT mass = 0.0;
|
||||
|
||||
// assemble 2nd order moment about the origin.
|
||||
FT temp[9] = {1.0/12.0, 1.0/24.0, 1.0/24.0,
|
||||
1.0/24.0, 1.0/12.0, 1.0/24.0,
|
||||
1.0/24.0, 1.0/24.0, 1.0/12.0};
|
||||
Matrix moment = init_Matrix<K>(3,temp);
|
||||
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
// Now for each triangle, construct the 2nd order moment about the origin.
|
||||
// assemble the transformation matrix.
|
||||
const Triangle& t = *it;
|
||||
|
||||
// defined for convenience.
|
||||
// FT example = CGAL::to_double(t[0].x());
|
||||
FT delta[9] = {t[0].x(), t[1].x(), t[2].x(),
|
||||
t[0].y(), t[1].y(), t[2].y(),
|
||||
t[0].z(), t[1].z(), t[2].z()};
|
||||
Matrix transformation = init_Matrix<K>(3,delta);
|
||||
FT area = std::sqrt(t.squared_area());
|
||||
CGAL_assertion(area != 0.0);
|
||||
|
||||
// Find the 2nd order moment for the triangle wrt to the origin by an affine transformation.
|
||||
|
||||
// Transform the standard 2nd order moment using the transformation matrix
|
||||
transformation = 2 * area * transformation * moment * LA::transpose(transformation);
|
||||
|
||||
// and add to covariance matrix
|
||||
covariance[0] += transformation[0][0];
|
||||
covariance[1] += transformation[1][0];
|
||||
covariance[2] += transformation[1][1];
|
||||
covariance[3] += transformation[2][0];
|
||||
covariance[4] += transformation[2][1];
|
||||
covariance[5] += transformation[2][2];
|
||||
|
||||
mass += area;
|
||||
}
|
||||
|
||||
// Translate the 2nd order moment calculated about the origin to
|
||||
// the center of mass to get the covariance.
|
||||
covariance[0] += mass * (-1.0 * c.x() * c.x());
|
||||
covariance[1] += mass * (-1.0 * c.x() * c.y());
|
||||
covariance[2] += mass * (-1.0 * c.y() * c.y());
|
||||
covariance[3] += mass * (-1.0 * c.z() * c.x());
|
||||
covariance[4] += mass * (-1.0 * c.z() * c.y());
|
||||
covariance[5] += mass * (-1.0 * c.z() * c.z());
|
||||
|
||||
}
|
||||
|
||||
// assemble covariance matrix from a cuboid set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
void
|
||||
assemble_covariance_matrix_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
const K& , // kernel
|
||||
const typename K::Iso_cuboid_3*,// used for indirection
|
||||
const CGAL::PCA_dimension_3_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
typedef typename CGAL::Linear_algebraCd<FT> LA;
|
||||
typedef typename LA::Matrix Matrix;
|
||||
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
// 1 2
|
||||
// 3 4 5
|
||||
//Final combined covariance matrix for all cuboids and their combined mass
|
||||
FT mass = 0.0;
|
||||
|
||||
// assemble 2nd order moment about the origin.
|
||||
FT temp[9] = {1.0/3.0, 1.0/4.0, 1.0/4.0,
|
||||
1.0/4.0, 1.0/3.0, 1.0/4.0,
|
||||
1.0/4.0, 1.0/4.0, 1.0/3.0};
|
||||
Matrix moment = init_Matrix<K>(3,temp);
|
||||
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
// Now for each cuboid, construct the 2nd order moment about the origin.
|
||||
// assemble the transformation matrix.
|
||||
const Iso_cuboid& t = *it;
|
||||
|
||||
// defined for convenience.
|
||||
// FT example = CGAL::to_double(t[0].x());
|
||||
FT x0 = t[0].x();
|
||||
FT y0 = t[0].y();
|
||||
FT z0 = t[0].z();
|
||||
FT delta[9] = {t[1].x()-x0, t[3].x()-x0, t[5].x()-x0,
|
||||
t[1].y()-y0, t[3].y()-y0, t[5].y()-y0,
|
||||
t[1].z()-z0, t[3].z()-z0, t[5].z()-z0};
|
||||
Matrix transformation = init_Matrix<K>(3,delta);
|
||||
FT volume = t.volume();
|
||||
CGAL_assertion(volume != 0.0);
|
||||
|
||||
// Find the 2nd order moment for the cuboid wrt to the origin by an affine transformation.
|
||||
|
||||
// Transform the standard 2nd order moment using the transformation matrix
|
||||
transformation = volume * transformation * moment * LA::transpose(transformation);
|
||||
|
||||
// Translate the 2nd order moment to the minimum corner (x0,y0,z0) of the cuboid.
|
||||
FT xav0 = (delta[0] + delta[1] + delta[2])/4.0;
|
||||
FT yav0 = (delta[3] + delta[4] + delta[5])/4.0;
|
||||
FT zav0 = (delta[6] + delta[7] + delta[8])/4.0;
|
||||
|
||||
// and add to covariance matrix
|
||||
covariance[0] += transformation[0][0] + volume * (2*x0*xav0 + x0*x0);
|
||||
covariance[1] += transformation[1][0] + volume * (xav0*y0 + yav0*x0 + x0*y0);
|
||||
covariance[2] += transformation[1][1] + volume * (2*y0*yav0 + y0*y0);
|
||||
covariance[3] += transformation[2][0] + volume * (x0*zav0 + xav0*z0 + x0*z0);
|
||||
covariance[4] += transformation[2][1] + volume * (yav0*z0 + y0*zav0 + z0*y0);
|
||||
covariance[5] += transformation[2][2] + volume * (2*zav0*z0 + z0*z0);
|
||||
|
||||
mass += volume;
|
||||
}
|
||||
|
||||
// Translate the 2nd order moment calculated about the origin to
|
||||
// the center of mass to get the covariance.
|
||||
covariance[0] += mass * (-1.0 * c.x() * c.x());
|
||||
covariance[1] += mass * (-1.0 * c.x() * c.y());
|
||||
covariance[2] += mass * (-1.0 * c.y() * c.y());
|
||||
covariance[3] += mass * (-1.0 * c.z() * c.x());
|
||||
covariance[4] += mass * (-1.0 * c.z() * c.y());
|
||||
covariance[5] += mass * (-1.0 * c.z() * c.z());
|
||||
|
||||
}
|
||||
|
||||
// assemble covariance matrix from a cuboid set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
void
|
||||
assemble_covariance_matrix_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
const K& , // kernel
|
||||
const typename K::Iso_cuboid_3*,// used for indirection
|
||||
const CGAL::PCA_dimension_2_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Iso_cuboid_3 Iso_cuboid;
|
||||
typedef typename CGAL::Linear_algebraCd<FT> LA;
|
||||
typedef typename LA::Matrix Matrix;
|
||||
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
// 1 2
|
||||
// 3 4 5
|
||||
//Final combined covariance matrix for all cuboids and their combined mass
|
||||
FT mass = 0.0;
|
||||
|
||||
// assemble 2nd order moment about the origin.
|
||||
FT temp[9] = {7.0/3.0, 1.5, 1.5,
|
||||
1.5, 7.0/3.0, 1.5,
|
||||
1.5, 1.5, 7.0/3.0};
|
||||
Matrix moment = init_Matrix<K>(3,temp);
|
||||
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
// Now for each cuboid, construct the 2nd order moment about the origin.
|
||||
// assemble the transformation matrix.
|
||||
const Iso_cuboid& t = *it;
|
||||
|
||||
// defined for convenience.
|
||||
// FT example = CGAL::to_double(t[0].x());
|
||||
FT x0 = t[0].x();
|
||||
FT y0 = t[0].y();
|
||||
FT z0 = t[0].z();
|
||||
FT delta[9] = {t[1].x()-x0, t[3].x()-x0, t[5].x()-x0,
|
||||
t[1].y()-y0, t[3].y()-y0, t[5].y()-y0,
|
||||
t[1].z()-z0, t[3].z()-z0, t[5].z()-z0};
|
||||
Matrix transformation = init_Matrix<K>(3,delta);
|
||||
FT volume = t.volume();
|
||||
CGAL_assertion(volume != 0.0);
|
||||
|
||||
// Find the 2nd order moment for the cuboid wrt to the origin by an affine transformation.
|
||||
|
||||
// Transform the standard 2nd order moment using the transformation matrix
|
||||
transformation = volume * transformation * moment * LA::transpose(transformation);
|
||||
|
||||
// Translate the 2nd order moment to the minimum corner (x0,y0,z0) of the cuboid.
|
||||
FT xav0 = (delta[0] + delta[1] + delta[2])/4.0;
|
||||
FT yav0 = (delta[3] + delta[4] + delta[5])/4.0;
|
||||
FT zav0 = (delta[6] + delta[7] + delta[8])/4.0;
|
||||
|
||||
// and add to covariance matrix
|
||||
covariance[0] += transformation[0][0] + volume * (2*x0*xav0 + x0*x0);
|
||||
covariance[1] += transformation[1][0] + volume * (xav0*y0 + yav0*x0 + x0*y0);
|
||||
covariance[2] += transformation[1][1] + volume * (2*y0*yav0 + y0*y0);
|
||||
covariance[3] += transformation[2][0] + volume * (x0*zav0 + xav0*z0 + x0*z0);
|
||||
covariance[4] += transformation[2][1] + volume * (yav0*z0 + y0*zav0 + z0*y0);
|
||||
covariance[5] += transformation[2][2] + volume * (2*zav0*z0 + z0*z0);
|
||||
|
||||
mass += volume;
|
||||
}
|
||||
|
||||
// Translate the 2nd order moment calculated about the origin to
|
||||
// the center of mass to get the covariance.
|
||||
covariance[0] += mass * (-1.0 * c.x() * c.x());
|
||||
covariance[1] += mass * (-1.0 * c.x() * c.y());
|
||||
covariance[2] += mass * (-1.0 * c.y() * c.y());
|
||||
covariance[3] += mass * (-1.0 * c.z() * c.x());
|
||||
covariance[4] += mass * (-1.0 * c.z() * c.y());
|
||||
covariance[5] += mass * (-1.0 * c.z() * c.z());
|
||||
|
||||
}
|
||||
|
||||
// assemble covariance matrix from a sphere set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
void
|
||||
assemble_covariance_matrix_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
const K& , // kernel
|
||||
const typename K::Sphere_3*,// used for indirection
|
||||
const CGAL::PCA_dimension_3_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Sphere_3 Sphere;
|
||||
typedef typename CGAL::Linear_algebraCd<FT> LA;
|
||||
typedef typename LA::Matrix Matrix;
|
||||
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
// 1 2
|
||||
// 3 4 5
|
||||
//Final combined covariance matrix for all spheres and their combined mass
|
||||
FT mass = 0.0;
|
||||
|
||||
// assemble 2nd order moment about the origin.
|
||||
FT temp[9] = {4.0/15.0, 0.0, 0.0,
|
||||
0.0, 4.0/15.0, 0.0,
|
||||
0.0, 0.0, 4.0/15.0};
|
||||
Matrix moment = init_Matrix<K>(3,temp);
|
||||
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
// Now for each sphere, construct the 2nd order moment about the origin.
|
||||
// assemble the transformation matrix.
|
||||
const Sphere& t = *it;
|
||||
|
||||
// defined for convenience.
|
||||
// FT example = CGAL::to_double(t[0].x());
|
||||
FT radius = std::sqrt(t.squared_radius());
|
||||
FT delta[9] = {radius, 0.0, 0.0,
|
||||
0.0, radius, 0.0,
|
||||
0.0, 0.0, radius};
|
||||
Matrix transformation = init_Matrix<K>(3,delta);
|
||||
FT area = 4/3.0 * radius*t.squared_radius();
|
||||
CGAL_assertion(area != 0.0);
|
||||
|
||||
// Find the 2nd order moment for the sphere wrt to the origin by an affine transformation.
|
||||
|
||||
// Transform the standard 2nd order moment using the transformation matrix
|
||||
transformation = (3.0/4.0) * area * transformation * moment * LA::transpose(transformation);
|
||||
|
||||
// Translate the 2nd order moment to the center of the sphere.
|
||||
FT x0 = t.center().x();
|
||||
FT y0 = t.center().y();
|
||||
FT z0 = t.center().z();
|
||||
|
||||
// and add to covariance matrix
|
||||
covariance[0] += transformation[0][0] + area * x0*x0;
|
||||
covariance[1] += transformation[1][0] + area * x0*y0;
|
||||
covariance[2] += transformation[1][1] + area * y0*y0;
|
||||
covariance[3] += transformation[2][0] + area * x0*z0;
|
||||
covariance[4] += transformation[2][1] + area * z0*y0;
|
||||
covariance[5] += transformation[2][2] + area * z0*z0;
|
||||
|
||||
mass += area;
|
||||
}
|
||||
|
||||
// Translate the 2nd order moment calculated about the origin to
|
||||
// the center of mass to get the covariance.
|
||||
covariance[0] += mass * (-1.0 * c.x() * c.x());
|
||||
covariance[1] += mass * (-1.0 * c.x() * c.y());
|
||||
covariance[2] += mass * (-1.0 * c.y() * c.y());
|
||||
covariance[3] += mass * (-1.0 * c.z() * c.x());
|
||||
covariance[4] += mass * (-1.0 * c.z() * c.y());
|
||||
covariance[5] += mass * (-1.0 * c.z() * c.z());
|
||||
|
||||
}
|
||||
|
||||
// assemble covariance matrix from a tetrahedron set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
void
|
||||
assemble_covariance_matrix_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
const K& , // kernel
|
||||
const typename K::Tetrahedron_3*,// used for indirection
|
||||
const CGAL::PCA_dimension_3_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Tetrahedron_3 Tetrahedron;
|
||||
typedef typename CGAL::Linear_algebraCd<FT> LA;
|
||||
typedef typename LA::Matrix Matrix;
|
||||
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
// 1 2
|
||||
// 3 4 5
|
||||
//Final combined covariance matrix for all tetrahedrons and their combined mass
|
||||
FT mass = 0.0;
|
||||
|
||||
// assemble 2nd order moment about the origin.
|
||||
FT temp[9] = {1.0/60.0, 1.0/120.0, 1.0/120.0,
|
||||
1.0/120.0, 1.0/60.0, 1.0/120.0,
|
||||
1.0/120.0, 1.0/120.0, 1.0/60.0};
|
||||
Matrix moment = init_Matrix<K>(3,temp);
|
||||
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
// Now for each tetrahedron, construct the 2nd order moment about the origin.
|
||||
// assemble the transformation matrix.
|
||||
const Tetrahedron& t = *it;
|
||||
|
||||
// defined for convenience.
|
||||
// FT example = CGAL::to_double(t[0].x());
|
||||
FT x0 = t[0].x();
|
||||
FT y0 = t[0].y();
|
||||
FT z0 = t[0].z();
|
||||
|
||||
FT delta[9] = {t[1].x()-x0, t[2].x()-x0, t[3].x()-x0,
|
||||
t[1].y()-y0, t[2].y()-y0, t[3].y()-y0,
|
||||
t[1].z()-z0, t[2].z()-z0, t[3].z()-z0};
|
||||
Matrix transformation = init_Matrix<K>(3,delta);
|
||||
FT volume = t.volume();
|
||||
CGAL_assertion(volume != 0.0);
|
||||
|
||||
// Find the 2nd order moment for the tetrahedron wrt to the origin by an affine transformation.
|
||||
|
||||
// Transform the standard 2nd order moment using the transformation matrix
|
||||
transformation = 6 * volume * transformation * moment * LA::transpose(transformation);
|
||||
|
||||
// Translate the 2nd order moment to the center of the tetrahedron.
|
||||
FT xav0 = (delta[0]+delta[1]+delta[2])/4.0;
|
||||
FT yav0 = (delta[3]+delta[4]+delta[5])/4.0;
|
||||
FT zav0 = (delta[6]+delta[7]+delta[8])/4.0;
|
||||
|
||||
// and add to covariance matrix
|
||||
covariance[0] += transformation[0][0] + volume * (2*x0*xav0 + x0*x0);
|
||||
covariance[1] += transformation[1][0] + volume * (xav0*y0 + yav0*x0 + x0*y0);
|
||||
covariance[2] += transformation[1][1] + volume * (2*y0*yav0 + y0*y0);
|
||||
covariance[3] += transformation[2][0] + volume * (x0*zav0 + xav0*z0 + x0*z0);
|
||||
covariance[4] += transformation[2][1] + volume * (yav0*z0 + y0*zav0 + z0*y0);
|
||||
covariance[5] += transformation[2][2] + volume * (2*zav0*z0 + z0*z0);
|
||||
|
||||
mass += volume;
|
||||
}
|
||||
|
||||
// Translate the 2nd order moment calculated about the origin to
|
||||
// the center of mass to get the covariance.
|
||||
covariance[0] += mass * (-1.0 * c.x() * c.x());
|
||||
covariance[1] += mass * (-1.0 * c.x() * c.y());
|
||||
covariance[2] += mass * (-1.0 * c.y() * c.y());
|
||||
covariance[3] += mass * (-1.0 * c.z() * c.x());
|
||||
covariance[4] += mass * (-1.0 * c.z() * c.y());
|
||||
covariance[5] += mass * (-1.0 * c.z() * c.z());
|
||||
|
||||
}
|
||||
|
||||
// assemble covariance matrix from a segment set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
void
|
||||
assemble_covariance_matrix_3(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
const K& , // kernel
|
||||
const typename K::Segment_3*,// used for indirection
|
||||
const CGAL::PCA_dimension_1_tag& tag)
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Segment_3 Segment;
|
||||
typedef typename CGAL::Linear_algebraCd<FT> LA;
|
||||
typedef typename LA::Matrix Matrix;
|
||||
|
||||
// assemble covariance matrix as a semi-definite matrix.
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
// 1 2
|
||||
// 3 4 5
|
||||
//Final combined covariance matrix for all segments and their combined mass
|
||||
FT mass = 0.0;
|
||||
|
||||
// assemble 2nd order moment about the origin.
|
||||
FT temp[9] = {1.0, 0.5, 0.0,
|
||||
0.5, 1.0, 0.0,
|
||||
0.0, 0.0, 0.0};
|
||||
Matrix moment = 1.0/3.0 * init_Matrix<K>(3,temp);
|
||||
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
// Now for each segment, construct the 2nd order moment about the origin.
|
||||
// assemble the transformation matrix.
|
||||
const Segment& t = *it;
|
||||
|
||||
// defined for convenience.
|
||||
// FT example = CGAL::to_double(t[0].x());
|
||||
FT delta[9] = {t[0].x(), t[1].x(), 0.0,
|
||||
t[0].y(), t[1].y(), 0.0,
|
||||
t[0].z(), t[1].z(), 1.0};
|
||||
Matrix transformation = init_Matrix<K>(3,delta);
|
||||
FT length = std::sqrt(t.squared_length());
|
||||
CGAL_assertion(length != 0.0);
|
||||
|
||||
// Find the 2nd order moment for the segment wrt to the origin by an affine transformation.
|
||||
|
||||
// Transform the standard 2nd order moment using the transformation matrix
|
||||
transformation = length * transformation * moment * LA::transpose(transformation);
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
// and add to covariance matrix
|
||||
covariance[0] += transformation[0][0];
|
||||
covariance[1] += transformation[1][0];
|
||||
covariance[2] += transformation[1][1];
|
||||
covariance[3] += transformation[2][0];
|
||||
covariance[4] += transformation[2][1];
|
||||
covariance[5] += transformation[2][2];
|
||||
std::cout<<covariance[0]<<" "<<covariance[1]<<" "<<covariance[2]<<" "<<covariance[3]<<" "<<covariance[4]<<" "<<covariance[5]<<std::endl;
|
||||
mass += length;
|
||||
}
|
||||
|
||||
// Translate the 2nd order moment calculated about the origin to
|
||||
// the center of mass to get the covariance.
|
||||
covariance[0] += mass * (-1.0 * c.x() * c.x());
|
||||
covariance[1] += mass * (-1.0 * c.x() * c.y());
|
||||
covariance[2] += mass * (-1.0 * c.y() * c.y());
|
||||
covariance[3] += mass * (-1.0 * c.z() * c.x());
|
||||
covariance[4] += mass * (-1.0 * c.z() * c.y());
|
||||
covariance[5] += mass * (-1.0 * c.z() * c.z());
|
||||
|
||||
}
|
||||
|
||||
/*
|
||||
// assemble covariance matrix from a triangle set
|
||||
template < typename InputIterator,
|
||||
typename K >
|
||||
void
|
||||
assemble_covariance_matrix_3_pierre(InputIterator first,
|
||||
InputIterator beyond,
|
||||
typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
const K& , // kernel
|
||||
const typename K::Triangle_3*)// used for indirection
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
typedef typename K::Triangle_3 Triangle;
|
||||
|
||||
// Matrix numbering:
|
||||
// 0
|
||||
// 1 2
|
||||
// 3 4 5
|
||||
covariance[0] = covariance[1] = covariance[2] =
|
||||
covariance[3] = covariance[4] = covariance[5] = (FT)0.0;
|
||||
FT sum_areas = 0.0;
|
||||
for(InputIterator it = first;
|
||||
it != beyond;
|
||||
it++)
|
||||
{
|
||||
const Triangle& triangle = *it;
|
||||
FT area = std::sqrt(triangle.squared_area());
|
||||
Point c_t = K().construct_centroid_3_object()(triangle[0],triangle[1],triangle[2]);
|
||||
sum_areas += area;
|
||||
|
||||
// e1 = ab, e2 = ac
|
||||
Vector e1 = Vector(triangle[0],triangle[1]);
|
||||
Vector e2 = Vector(triangle[0],triangle[2]);
|
||||
|
||||
FT c1 = (FT)(2.0 * area * 10.0 / 72.0);
|
||||
FT c2 = (FT)(2.0 * area * 7.0 / 72.0);
|
||||
|
||||
covariance[0] += c1*(e1[0]*e1[0] + e2[0]*e2[0]) + (FT)2.0*c2*e1[0]*e2[0];
|
||||
covariance[1] += c1*(e1[1]*e1[0] + e2[1]*e2[0]) + c2*(e1[1]*e2[0] + e1[0]*e2[1]);
|
||||
covariance[2] += c1*(e1[1]*e1[1] + e2[1]*e2[1]) + (FT)2.0*c2*e1[1]*e2[1];
|
||||
covariance[3] += c1*(e1[2]*e1[0] + e2[2]*e2[0]) + c2*(e1[2]*e2[0] + e1[0]*e2[2]);
|
||||
covariance[4] += c1*(e1[2]*e1[1] + e2[2]*e2[1]) + c2*(e1[2]*e2[1] + e1[1]*e2[2]);
|
||||
covariance[5] += c1*(e1[2]*e1[2] + e2[2]*e2[2]) + (FT)2.0*c2*e1[2]*e2[2];
|
||||
|
||||
// add area(t) c(t) * transpose(c(t))
|
||||
covariance[0] += area * c_t.x() * c_t.x();
|
||||
covariance[1] += area * c_t.y() * c_t.x();
|
||||
covariance[2] += area * c_t.y() * c_t.y();
|
||||
covariance[3] += area * c_t.z() * c_t.x();
|
||||
covariance[4] += area * c_t.z() * c_t.y();
|
||||
covariance[5] += area * c_t.z() * c_t.z();
|
||||
}
|
||||
|
||||
// remove sum(area) * (c * transpose(c))
|
||||
covariance[0] -= sum_areas * c.x() * c.x();
|
||||
covariance[1] -= sum_areas * c.y() * c.x();
|
||||
covariance[2] -= sum_areas * c.y() * c.y();
|
||||
covariance[3] -= sum_areas * c.z() * c.x();
|
||||
covariance[4] -= sum_areas * c.z() * c.y();
|
||||
covariance[5] -= sum_areas * c.z() * c.z();
|
||||
}
|
||||
*/
|
||||
// compute the eigen values and vectors of the covariance
|
||||
// matrix and deduces the best linear fitting plane.
|
||||
// returns fitting quality
|
||||
template < typename K >
|
||||
typename K::FT
|
||||
fitting_plane_3(const typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
typename K::Plane_3& plane, // best fit plane
|
||||
const K& ) // kernel
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Plane_3 Plane;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
|
||||
// solve for eigenvalues and eigenvectors.
|
||||
// eigen values are sorted in descending order,
|
||||
// eigen vectors are sorted in accordance.
|
||||
FT eigen_values[3];
|
||||
FT eigen_vectors[9];
|
||||
eigen_symmetric<FT>(covariance,3,eigen_vectors,eigen_values);
|
||||
|
||||
// check unicity and build fitting line accordingly
|
||||
if(eigen_values[0] != eigen_values[1] &&
|
||||
eigen_values[0] != eigen_values[2])
|
||||
{
|
||||
// regular case
|
||||
Vector normal(eigen_vectors[6],
|
||||
eigen_vectors[7],
|
||||
eigen_vectors[8]);
|
||||
plane = Plane(c,normal);
|
||||
return (FT)1.0 - eigen_values[2] / eigen_values[0];
|
||||
} // end regular case
|
||||
else
|
||||
{
|
||||
// isotropic case (infinite number of directions)
|
||||
// by default: assemble a horizontal plane that goes
|
||||
// through the centroid.
|
||||
plane = Plane(c,Vector(0.0,0.0,1.0));
|
||||
return (FT)0.0;
|
||||
}
|
||||
}
|
||||
|
||||
// compute the eigen values and vectors of the covariance
|
||||
// matrix and deduces the best linear fitting line
|
||||
// (this is an internal function)
|
||||
// returns fitting quality
|
||||
template < typename K >
|
||||
typename K::FT
|
||||
fitting_line_3(const typename K::FT covariance[6], // covariance matrix
|
||||
const typename K::Point_3& c, // centroid
|
||||
typename K::Line_3& line, // best fit line
|
||||
const K& ) // kernel
|
||||
{
|
||||
typedef typename K::FT FT;
|
||||
typedef typename K::Point_3 Point;
|
||||
typedef typename K::Line_3 Line;
|
||||
typedef typename K::Vector_3 Vector;
|
||||
|
||||
// solve for eigenvalues and eigenvectors.
|
||||
// eigen values are sorted in descending order,
|
||||
// eigen vectors are sorted in accordance.
|
||||
FT eigen_values[3];
|
||||
FT eigen_vectors[9];
|
||||
eigen_symmetric<FT>(covariance,3,eigen_vectors,eigen_values);
|
||||
|
||||
// check unicity and build fitting line accordingly
|
||||
if(eigen_values[0] != eigen_values[1])
|
||||
{
|
||||
// regular case
|
||||
Vector direction(eigen_vectors[0],eigen_vectors[1],eigen_vectors[2]);
|
||||
line = Line(c,direction);
|
||||
return (FT)1.0 - eigen_values[1] / eigen_values[0];
|
||||
} // end regular case
|
||||
else
|
||||
{
|
||||
// isotropic case (infinite number of directions)
|
||||
// by default: assemble a horizontal plane that goes
|
||||
// through the centroid.
|
||||
line = Line(c,Vector(1.0,0.0,0.0));
|
||||
return (FT)0.0;
|
||||
}
|
||||
}
|
||||
|
||||
} // end namespace CGALi
|
||||
|
||||
CGAL_END_NAMESPACE
|
||||
|
|
|
|||
Loading…
Reference in New Issue