From 3a2303e66d0ac33b3fed7d4b2f61ccf760189e2b Mon Sep 17 00:00:00 2001 From: Paolo Cignoni cignoni Date: Sun, 3 Apr 2011 23:44:39 +0000 Subject: [PATCH] Corrections of include paths to comply the new folder arrangement of the VCG library (a first set of the many filters; removing useless includes) --- .../filter_plymc/filter_plymc.cpp | 2 +- .../filter_plymc/filter_plymc.pro | 6 +- src/meshlabplugins/filter_plymc/plymc.h | 28 +- .../filter_poisson/filter_poisson.cpp | 2 +- .../filter_sampling/filter_sampling.cpp | 12 +- .../filter_sampling/filter_sampling.pro | 2 +- src/meshlabplugins/filter_sdf/filter_sdf+.cpp | 16 +- src/meshlabplugins/filter_sdf/filter_sdf.cpp | 17 +- .../filter_select/filter_select.pro | 2 +- .../filter_select/meshselect.cpp | 4 +- .../filter_slice/filter_slice.cpp | 3 +- .../filter_slice/filter_slice.h | 8 +- src/meshlabplugins/filter_slice/kdtree.h | 396 +- .../filter_texture/filter_texture.h | 4 +- .../filter_texture/filter_texture_old.cpp | 2 +- src/meshlabplugins/filter_texture/rastering.h | 2 +- .../filter_unsharp/filter_unsharp.cpp | 9 +- .../filter_unsharp/filter_unsharp.pro | 2 +- .../filter_vsa/region_growing.h | 970 ++--- src/meshlabplugins/filter_vsa/vsa.cpp | 7 +- .../filter_zippering/filter_zippering.cpp | 3404 ++++++++--------- .../filter_zippering/filter_zippering.h | 940 ++--- 22 files changed, 2908 insertions(+), 2930 deletions(-) diff --git a/src/meshlabplugins/filter_plymc/filter_plymc.cpp b/src/meshlabplugins/filter_plymc/filter_plymc.cpp index 6ef337fa1..9c6ca40a8 100644 --- a/src/meshlabplugins/filter_plymc/filter_plymc.cpp +++ b/src/meshlabplugins/filter_plymc/filter_plymc.cpp @@ -22,7 +22,7 @@ ****************************************************************************/ #include "filter_plymc.h" -#include +#include #include #include "plymc.h" diff --git a/src/meshlabplugins/filter_plymc/filter_plymc.pro b/src/meshlabplugins/filter_plymc/filter_plymc.pro index 68476ee01..9adb6d354 100644 --- a/src/meshlabplugins/filter_plymc/filter_plymc.pro +++ b/src/meshlabplugins/filter_plymc/filter_plymc.pro @@ -5,12 +5,12 @@ HEADERS += filter_plymc.h \ voxel.h \ plymc.h \ simplemeshprovider.h \ -tri_edge_collapse_mc.h \ - $$VCGDIR/vcg/complex/local_optimization/tri_edge_collapse.h + tri_edge_collapse_mc.h \ + $$VCGDIR/vcg/algorithms/local_optimization/tri_edge_collapse.h SOURCES += filter_plymc.cpp \ -../../../../vcglib/wrap/ply/plylib.cpp + $$VCGDIR/wrap/ply/plylib.cpp TARGET = filter_plymc QT += opengl diff --git a/src/meshlabplugins/filter_plymc/plymc.h b/src/meshlabplugins/filter_plymc/plymc.h index 8c35cf019..571afb17e 100644 --- a/src/meshlabplugins/filter_plymc/plymc.h +++ b/src/meshlabplugins/filter_plymc/plymc.h @@ -30,28 +30,28 @@ #include #include #include -#include +#include -#include -#include -#include -#include -#include -#include -#include +#include +#include +#include +#include +#include +#include +#include #include -#include +#include #include #include #include -#include -#include +#include +#include #include "trivial_walker.h" // local optimization -#include -#include -#include +#include +#include +#include #include #include diff --git a/src/meshlabplugins/filter_poisson/filter_poisson.cpp b/src/meshlabplugins/filter_poisson/filter_poisson.cpp index 7debd7ef6..f58496658 100644 --- a/src/meshlabplugins/filter_poisson/filter_poisson.cpp +++ b/src/meshlabplugins/filter_poisson/filter_poisson.cpp @@ -41,7 +41,7 @@ add samplefilter #include #include -#include +#include #include "filter_poisson.h" #include "src/Geometry.h" diff --git a/src/meshlabplugins/filter_sampling/filter_sampling.cpp b/src/meshlabplugins/filter_sampling/filter_sampling.cpp index 2e02ecc65..7aba89395 100644 --- a/src/meshlabplugins/filter_sampling/filter_sampling.cpp +++ b/src/meshlabplugins/filter_sampling/filter_sampling.cpp @@ -34,12 +34,12 @@ $Log: samplefilter.cpp,v $ #include "filter_sampling.h" -#include -#include -#include -#include +#include +#include +#include +#include #include -#include +#include #include #include "voronoi_clustering.h" @@ -54,7 +54,7 @@ using namespace std; For example, you can compute Hausdorff distance (that is a sampler) using various sampling strategies (montecarlo, stratified etc). - For further details see vcg/complex/trimesh/point_sampling.h header file. + For further details see vcg/complex/algorithms/point_sampling.h header file. */ diff --git a/src/meshlabplugins/filter_sampling/filter_sampling.pro b/src/meshlabplugins/filter_sampling/filter_sampling.pro index 6f802048e..f749ae07c 100644 --- a/src/meshlabplugins/filter_sampling/filter_sampling.pro +++ b/src/meshlabplugins/filter_sampling/filter_sampling.pro @@ -1,7 +1,7 @@ include (../../shared.pri) HEADERS += filter_sampling.h \ - $$VCGDIR/vcg/complex/trimesh/point_sampling.h + $$VCGDIR/vcg/complex/algorithms/point_sampling.h SOURCES += filter_sampling.cpp TARGET = filter_sampling diff --git a/src/meshlabplugins/filter_sdf/filter_sdf+.cpp b/src/meshlabplugins/filter_sdf/filter_sdf+.cpp index c88323f7d..e191f373b 100644 --- a/src/meshlabplugins/filter_sdf/filter_sdf+.cpp +++ b/src/meshlabplugins/filter_sdf/filter_sdf+.cpp @@ -1,11 +1,11 @@ -#include -#include -#include -#include -#include -#include -#include -#include +#include +#include +#include +#include +#include +#include +#include +#include #include #include #include diff --git a/src/meshlabplugins/filter_sdf/filter_sdf.cpp b/src/meshlabplugins/filter_sdf/filter_sdf.cpp index 22cc9d489..c943b4685 100644 --- a/src/meshlabplugins/filter_sdf/filter_sdf.cpp +++ b/src/meshlabplugins/filter_sdf/filter_sdf.cpp @@ -1,15 +1,6 @@ -#include -#include -#include -#include -#include -#include -#include -#include -#include +#include #include #include -#include #include // vcg::GridStaticPtr #include // vcg::SpatialHashTable @@ -107,7 +98,7 @@ bool SdfPlugin::applyFilter(MeshDocument& md, RichParameterSet& pars, vcg::CallB float maxDist=m.bbox.Diag(); // This is a small number to avoid self-intersection during ray // casting. It's a very common trick - float epsilon = maxDist / 10000.0; + float epsilon = maxDist / 1000.0; //--- Ray casting vector cone; @@ -118,7 +109,7 @@ bool SdfPlugin::applyFilter(MeshDocument& md, RichParameterSet& pars, vcg::CallB for(int i=0; iloperc && coneSdf[i]=loperc && coneSdf[i]<=hiperc ){ totVal += coneSdf[i]; totCnt += 1; } diff --git a/src/meshlabplugins/filter_select/filter_select.pro b/src/meshlabplugins/filter_select/filter_select.pro index 3ccbc4398..72436036c 100644 --- a/src/meshlabplugins/filter_select/filter_select.pro +++ b/src/meshlabplugins/filter_select/filter_select.pro @@ -1,6 +1,6 @@ include (../../shared.pri) -HEADERS += $$VCGDIR/vcg/complex/trimesh/clean.h\ +HEADERS += $$VCGDIR/vcg/complex/algorithms/clean.h\ meshselect.h SOURCES += meshselect.cpp diff --git a/src/meshlabplugins/filter_select/meshselect.cpp b/src/meshlabplugins/filter_select/meshselect.cpp index a717683f7..638ca3f30 100644 --- a/src/meshlabplugins/filter_select/meshselect.cpp +++ b/src/meshlabplugins/filter_select/meshselect.cpp @@ -26,8 +26,8 @@ #include #include #include "meshselect.h" -#include -#include +#include +#include using namespace vcg; diff --git a/src/meshlabplugins/filter_slice/filter_slice.cpp b/src/meshlabplugins/filter_slice/filter_slice.cpp index 8bbf752e9..561ce188f 100644 --- a/src/meshlabplugins/filter_slice/filter_slice.cpp +++ b/src/meshlabplugins/filter_slice/filter_slice.cpp @@ -23,13 +23,12 @@ #include "filter_slice.h" //#include -#include +#include #include #include "filter_slice_functors.h" #include -#include #include "kdtree.h" diff --git a/src/meshlabplugins/filter_slice/filter_slice.h b/src/meshlabplugins/filter_slice/filter_slice.h index c394d65f3..b14795e68 100644 --- a/src/meshlabplugins/filter_slice/filter_slice.h +++ b/src/meshlabplugins/filter_slice/filter_slice.h @@ -30,12 +30,8 @@ #include #include -#include -#include -#include -#include -#include -#include +#include +#include #include #include diff --git a/src/meshlabplugins/filter_slice/kdtree.h b/src/meshlabplugins/filter_slice/kdtree.h index 236f7c202..5c4252fa1 100644 --- a/src/meshlabplugins/filter_slice/kdtree.h +++ b/src/meshlabplugins/filter_slice/kdtree.h @@ -1,198 +1,198 @@ -/**************************************************************************** -* MeshLab o o * -* An extendible mesh processor o o * -* _ O _ * -* Copyright(C) 2005, 2009 \/)\/ * -* Visual Computing Lab /\/| * -* ISTI - Italian National Research Council | * -* \ * -* All rights reserved. * -* * -* This program is free software; you can redistribute it and/or modify * -* it under the terms of the GNU General Public License as published by * -* the Free Software Foundation; either version 2 of the License, or * -* (at your option) any later version. * -* * -* This program is distributed in the hope that it will be useful, * -* but WITHOUT ANY WARRANTY; without even the implied warranty of * -* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * -* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) * -* for more details. * -* * -****************************************************************************/ - -#ifndef SLICE_KD_TREE -#define SLICE_KD_TREE -#include -#include -#include "filter_slice_functors.h" -#include "filter_slice.h" - -template -class KDTree -{ - protected: - KDTree *leftChild, *rightChild; - MeshType *slice; - MeshModel* mm; - Point3f planeAxis; - MeshDocument *m; - int axisIndex; - float eps; - QString name; - public: - bool delfather; - enum{X,Y,Z}; - KDTree(MeshDocument *_m, MeshModel* _mm, float _eps, int _axisIndex=X) - { - name=""; - m=_m; - mm=_mm; - eps=_eps; - planeAxis=Point3f(0,0,0); - leftChild=0; - rightChild=0; - slice=0; - axisIndex=_axisIndex; - assert(axisIndex>=0 && axisIndex<3); - planeAxis[axisIndex]=1; - delfather=false; - } - ~KDTree(){} - KDTree* L(){return leftChild;} - KDTree* R(){return rightChild;} - MeshType* S(){return slice;} - void Slice(vcg::CallBackPos *cb) - { - Plane3f slicingPlane; - vcg::tri::UpdateBounding::Box(mm->cm); - slicingPlane.Init(mm->cm.bbox.Center(),planeAxis); - actual_slice(slicingPlane,cb); - } - - void actual_slice(Plane3f slicingPlane,vcg::CallBackPos *cb=0) - { - if (mm->cm.vn<=0) - return; - mm->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); - if ( tri::Clean::CountNonManifoldEdgeFF(mm->cm)>0 || tri::Clean::CountNonManifoldVertexFF(mm->cm,false) != 0) - return; - - if (rightChild!=0 && leftChild!=0) - { - //find a plane - Plane3f slicingPlane; - vcg::Box3f leftbox=leftChild->mm->cm.bbox; - vcg::Box3f rightbox=rightChild->mm->cm.bbox; - leftbox.Translate(planeAxis*((leftbox.Dim()*.2)*planeAxis)); - rightbox.Translate(planeAxis*(-1)*((rightbox.Dim()*.2)*planeAxis)); - leftbox.Intersect(rightbox); - slicingPlane.Init(leftbox.Center(),leftChild->planeAxis); - - leftChild->actual_slice(slicingPlane,cb); - rightChild->actual_slice(slicingPlane,cb); - return; - } - //bool oriented; - //bool orientable; - //tolto per controllare le normali - //tri::Clean::IsOrientedMesh(mm->cm, oriented,orientable); - //actual slicing - - - //vcg::tri::UpdateFlags::FaceBorderFromNone(mm->cm); - - SlicedEdge slicededge(slicingPlane); - SlicingFunction slicingfunc(slicingPlane); - //after the RefineE call, the mesh will be half vertices selected - vcg::RefineE, SlicedEdge > - (mm->cm, slicingfunc, slicededge, false, cb); - //vcg::tri::UpdateTopology::FaceFace(mm->cm); - vcg::tri::UpdateNormals::PerVertexPerFace(mm->cm); - //vcg::tri::Clean::RemoveDuplicateVertex(mm->cm); - vcg::tri::UpdateTopology::FaceFace(mm->cm); - - - - - QString layername; - layername=name+"L.ply"; - MeshModel *slice1= m->addNewMesh("",qPrintable(layername)); - m->meshList.push_back(slice1); // mesh name - slice1->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); - vcg::tri::UpdateSelection::VertexFromQualityRange(mm->cm,VERTEX_LEFT,VERTEX_LEFT); - vcg::tri::UpdateSelection::FaceFromVertexLoose(mm->cm); - //if (hideSlices) - //slice1->visible=false; - //createSlice(mm->cm,slice1); - tri::Append::Mesh(slice1->cm, mm->cm, true); - tri::UpdateTopology::FaceFace(slice1->cm); - tri::UpdateBounding::Box(slice1->cm); // updates bounding box - slice1->cm.Tr = (mm->cm).Tr; // copy transformation - vcg::tri::UpdateFlags::FaceBorderFromNone(slice1->cm); - layername=name+"_slice.ply"; - MeshModel* cap= m->addNewMesh("",qPrintable(layername)); - m->meshList.push_back(cap); - cap->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); - ExtraFilter_SlicePlugin::capHole(slice1,cap); - - if (eps!=0) - { - layername=name+"_extr.ply"; - MeshModel* dup = m->addNewMesh("",qPrintable(layername)); - m->meshList.push_back(dup); - dup->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); - ExtraFilter_SlicePlugin::extrude(m,cap, dup, eps, planeAxis); - } - - tri::Append::Mesh(slice1->cm, cap->cm); - tri::Clean::RemoveDuplicateVertex(slice1->cm); - vcg::tri::UpdateTopology::FaceFace(slice1->cm); - vcg::tri::UpdateNormals::PerVertexPerFace(slice1->cm); - layername=name+"R.ply"; - MeshModel* slice2= m->addNewMesh("",qPrintable(layername)); - m->meshList.push_back(slice2); - - slice2->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); - vcg::tri::UpdateSelection::VertexFromQualityRange(mm->cm,VERTEX_RIGHT,VERTEX_RIGHT); - vcg::tri::UpdateSelection::FaceFromVertexLoose(mm->cm); - //createSlice(mesh,slice2); - tri::Append::Mesh(slice2->cm, mm->cm, true); - tri::UpdateTopology::FaceFace(slice2->cm); - tri::UpdateBounding::Box(slice2->cm); // updates bounding box - slice2->cm.Tr = (mm->cm).Tr; // copy transformation - tri::Clean::FlipMesh(cap->cm); - vcg::tri::UpdateFlags::FaceBorderFromNone(slice2->cm); - tri::Append::Mesh(slice2->cm, cap->cm); - tri::Clean::RemoveDuplicateVertex(slice2->cm); - vcg::tri::UpdateTopology::FaceFace(slice2->cm); - vcg::tri::UpdateNormals::PerVertexPerFace(slice2->cm); - - leftChild=new KDTree(m,slice1,eps,(axisIndex+1)%3); - leftChild->delfather=delfather; - layername=name+"L"; - leftChild->name=layername.toStdString().c_str(); - rightChild=new KDTree(m,slice2,eps,(axisIndex+1)%3); - rightChild->delfather=delfather; - layername=name+"R"; - rightChild->name=layername.toStdString().c_str(); - - slice1->visible=false; - slice1->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); - vcg::tri::UpdateTopology::FaceFace(slice1->cm); - vcg::tri::UpdateBounding::Box(slice1->cm); - - slice2->visible=false; - slice2->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); - vcg::tri::UpdateTopology::FaceFace(slice2->cm); - vcg::tri::UpdateBounding::Box(slice2->cm); - - cap->visible=false; - //if (delfather) - //m->delMesh(mm); - - } - -}; - -#endif +/**************************************************************************** +* MeshLab o o * +* An extendible mesh processor o o * +* _ O _ * +* Copyright(C) 2005, 2009 \/)\/ * +* Visual Computing Lab /\/| * +* ISTI - Italian National Research Council | * +* \ * +* All rights reserved. * +* * +* This program is free software; you can redistribute it and/or modify * +* it under the terms of the GNU General Public License as published by * +* the Free Software Foundation; either version 2 of the License, or * +* (at your option) any later version. * +* * +* This program is distributed in the hope that it will be useful, * +* but WITHOUT ANY WARRANTY; without even the implied warranty of * +* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * +* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) * +* for more details. * +* * +****************************************************************************/ + +#ifndef SLICE_KD_TREE +#define SLICE_KD_TREE +#include +#include +#include "filter_slice_functors.h" +#include "filter_slice.h" + +template +class KDTree +{ + protected: + KDTree *leftChild, *rightChild; + MeshType *slice; + MeshModel* mm; + Point3f planeAxis; + MeshDocument *m; + int axisIndex; + float eps; + QString name; + public: + bool delfather; + enum{X,Y,Z}; + KDTree(MeshDocument *_m, MeshModel* _mm, float _eps, int _axisIndex=X) + { + name=""; + m=_m; + mm=_mm; + eps=_eps; + planeAxis=Point3f(0,0,0); + leftChild=0; + rightChild=0; + slice=0; + axisIndex=_axisIndex; + assert(axisIndex>=0 && axisIndex<3); + planeAxis[axisIndex]=1; + delfather=false; + } + ~KDTree(){} + KDTree* L(){return leftChild;} + KDTree* R(){return rightChild;} + MeshType* S(){return slice;} + void Slice(vcg::CallBackPos *cb) + { + Plane3f slicingPlane; + vcg::tri::UpdateBounding::Box(mm->cm); + slicingPlane.Init(mm->cm.bbox.Center(),planeAxis); + actual_slice(slicingPlane,cb); + } + + void actual_slice(Plane3f slicingPlane,vcg::CallBackPos *cb=0) + { + if (mm->cm.vn<=0) + return; + mm->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); + if ( tri::Clean::CountNonManifoldEdgeFF(mm->cm)>0 || tri::Clean::CountNonManifoldVertexFF(mm->cm,false) != 0) + return; + + if (rightChild!=0 && leftChild!=0) + { + //find a plane + Plane3f slicingPlane; + vcg::Box3f leftbox=leftChild->mm->cm.bbox; + vcg::Box3f rightbox=rightChild->mm->cm.bbox; + leftbox.Translate(planeAxis*((leftbox.Dim()*.2)*planeAxis)); + rightbox.Translate(planeAxis*(-1)*((rightbox.Dim()*.2)*planeAxis)); + leftbox.Intersect(rightbox); + slicingPlane.Init(leftbox.Center(),leftChild->planeAxis); + + leftChild->actual_slice(slicingPlane,cb); + rightChild->actual_slice(slicingPlane,cb); + return; + } + //bool oriented; + //bool orientable; + //tolto per controllare le normali + //tri::Clean::IsOrientedMesh(mm->cm, oriented,orientable); + //actual slicing + + + //vcg::tri::UpdateFlags::FaceBorderFromNone(mm->cm); + + SlicedEdge slicededge(slicingPlane); + SlicingFunction slicingfunc(slicingPlane); + //after the RefineE call, the mesh will be half vertices selected + vcg::RefineE, SlicedEdge > + (mm->cm, slicingfunc, slicededge, false, cb); + //vcg::tri::UpdateTopology::FaceFace(mm->cm); + vcg::tri::UpdateNormals::PerVertexPerFace(mm->cm); + //vcg::tri::Clean::RemoveDuplicateVertex(mm->cm); + vcg::tri::UpdateTopology::FaceFace(mm->cm); + + + + + QString layername; + layername=name+"L.ply"; + MeshModel *slice1= m->addNewMesh("",qPrintable(layername)); + m->meshList.push_back(slice1); // mesh name + slice1->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); + vcg::tri::UpdateSelection::VertexFromQualityRange(mm->cm,VERTEX_LEFT,VERTEX_LEFT); + vcg::tri::UpdateSelection::FaceFromVertexLoose(mm->cm); + //if (hideSlices) + //slice1->visible=false; + //createSlice(mm->cm,slice1); + tri::Append::Mesh(slice1->cm, mm->cm, true); + tri::UpdateTopology::FaceFace(slice1->cm); + tri::UpdateBounding::Box(slice1->cm); // updates bounding box + slice1->cm.Tr = (mm->cm).Tr; // copy transformation + vcg::tri::UpdateFlags::FaceBorderFromNone(slice1->cm); + layername=name+"_slice.ply"; + MeshModel* cap= m->addNewMesh("",qPrintable(layername)); + m->meshList.push_back(cap); + cap->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); + ExtraFilter_SlicePlugin::capHole(slice1,cap); + + if (eps!=0) + { + layername=name+"_extr.ply"; + MeshModel* dup = m->addNewMesh("",qPrintable(layername)); + m->meshList.push_back(dup); + dup->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); + ExtraFilter_SlicePlugin::extrude(m,cap, dup, eps, planeAxis); + } + + tri::Append::Mesh(slice1->cm, cap->cm); + tri::Clean::RemoveDuplicateVertex(slice1->cm); + vcg::tri::UpdateTopology::FaceFace(slice1->cm); + vcg::tri::UpdateNormals::PerVertexPerFace(slice1->cm); + layername=name+"R.ply"; + MeshModel* slice2= m->addNewMesh("",qPrintable(layername)); + m->meshList.push_back(slice2); + + slice2->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); + vcg::tri::UpdateSelection::VertexFromQualityRange(mm->cm,VERTEX_RIGHT,VERTEX_RIGHT); + vcg::tri::UpdateSelection::FaceFromVertexLoose(mm->cm); + //createSlice(mesh,slice2); + tri::Append::Mesh(slice2->cm, mm->cm, true); + tri::UpdateTopology::FaceFace(slice2->cm); + tri::UpdateBounding::Box(slice2->cm); // updates bounding box + slice2->cm.Tr = (mm->cm).Tr; // copy transformation + tri::Clean::FlipMesh(cap->cm); + vcg::tri::UpdateFlags::FaceBorderFromNone(slice2->cm); + tri::Append::Mesh(slice2->cm, cap->cm); + tri::Clean::RemoveDuplicateVertex(slice2->cm); + vcg::tri::UpdateTopology::FaceFace(slice2->cm); + vcg::tri::UpdateNormals::PerVertexPerFace(slice2->cm); + + leftChild=new KDTree(m,slice1,eps,(axisIndex+1)%3); + leftChild->delfather=delfather; + layername=name+"L"; + leftChild->name=layername.toStdString().c_str(); + rightChild=new KDTree(m,slice2,eps,(axisIndex+1)%3); + rightChild->delfather=delfather; + layername=name+"R"; + rightChild->name=layername.toStdString().c_str(); + + slice1->visible=false; + slice1->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); + vcg::tri::UpdateTopology::FaceFace(slice1->cm); + vcg::tri::UpdateBounding::Box(slice1->cm); + + slice2->visible=false; + slice2->updateDataMask(MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEFLAGBORDER); + vcg::tri::UpdateTopology::FaceFace(slice2->cm); + vcg::tri::UpdateBounding::Box(slice2->cm); + + cap->visible=false; + //if (delfather) + //m->delMesh(mm); + + } + +}; + +#endif diff --git a/src/meshlabplugins/filter_texture/filter_texture.h b/src/meshlabplugins/filter_texture/filter_texture.h index fb7e8d866..b928d32d1 100644 --- a/src/meshlabplugins/filter_texture/filter_texture.h +++ b/src/meshlabplugins/filter_texture/filter_texture.h @@ -27,8 +27,8 @@ #include #include -#include -#include +#include +#include #include #include diff --git a/src/meshlabplugins/filter_texture/filter_texture_old.cpp b/src/meshlabplugins/filter_texture/filter_texture_old.cpp index 6edb302bf..206b5b53b 100644 --- a/src/meshlabplugins/filter_texture/filter_texture_old.cpp +++ b/src/meshlabplugins/filter_texture/filter_texture_old.cpp @@ -33,7 +33,7 @@ $Log: samplefilter.cpp,v $ #include #include -#include +#include #include "filter_texture.h" diff --git a/src/meshlabplugins/filter_texture/rastering.h b/src/meshlabplugins/filter_texture/rastering.h index a278d96df..6640a9da4 100644 --- a/src/meshlabplugins/filter_texture/rastering.h +++ b/src/meshlabplugins/filter_texture/rastering.h @@ -26,7 +26,7 @@ #include #include -#include +#include #include class VertexSampler diff --git a/src/meshlabplugins/filter_unsharp/filter_unsharp.cpp b/src/meshlabplugins/filter_unsharp/filter_unsharp.cpp index a8c42a893..fa57413dc 100644 --- a/src/meshlabplugins/filter_unsharp/filter_unsharp.cpp +++ b/src/meshlabplugins/filter_unsharp/filter_unsharp.cpp @@ -24,12 +24,9 @@ #include #include "filter_unsharp.h" -#include -#include -#include -#include -#include -#include +#include +#include +#include using namespace vcg; diff --git a/src/meshlabplugins/filter_unsharp/filter_unsharp.pro b/src/meshlabplugins/filter_unsharp/filter_unsharp.pro index cc06b8d9b..19e5491ed 100644 --- a/src/meshlabplugins/filter_unsharp/filter_unsharp.pro +++ b/src/meshlabplugins/filter_unsharp/filter_unsharp.pro @@ -1,7 +1,7 @@ include (../../shared.pri) HEADERS += filter_unsharp.h \ - $$VCGDIR/vcg/complex/trimesh/crease_cut.h + $$VCGDIR/vcg/complex/algorithms/crease_cut.h SOURCES += filter_unsharp.cpp diff --git a/src/meshlabplugins/filter_vsa/region_growing.h b/src/meshlabplugins/filter_vsa/region_growing.h index 5bc59b69a..48a6c9ca5 100755 --- a/src/meshlabplugins/filter_vsa/region_growing.h +++ b/src/meshlabplugins/filter_vsa/region_growing.h @@ -1,484 +1,484 @@ -#ifndef VCGLIB_REGION_GROWING -#define VCGLIB_REGION_GROWING - -#include -#include - -#include - -template -struct FaceError{ - FaceType * f; - float val; - FaceError(){}; - FaceError(FaceType * _f,double _val):f(_f),val(_val){} - const bool operator <(const FaceError & o) const {return val < o.val;} -}; - - -template -struct Region{ - Region():isd(false){} - typedef Region RegionType; - typedef typename MeshType MeshType; - typedef typename MeshType::FaceType FaceType; - typedef typename MeshType::FaceType::CoordType CoordType; - typedef typename MeshType::FaceType::ScalarType ScalarType; - typedef typename std::vector::iterator FaceIterator; - typedef typename std::list::iterator AdjIterator; - typedef FaceError FaceError; - - // adjacent regions - std::list adj,nx_adj; - - - // some flag - bool isd; - - // evaluate the gain if the triangle is added to the region - ScalarType Evaluate( FaceType & f); - - // update the approximation error of the region - void UpdateError( ); - - // check if two regions are mergeable in one - // i.e. if they have the same planes and coincide i non extreme - ScalarType Mergeable( RegionType & tr); - - // set the region as adjacent to this one - void Connect( RegionType * tr); - - // compute the faces which are on the border - void ComputeOnBorder( std::vector & onBorder); - - // clean the list of adjacencies - void Clean(); - - // refit the plane - void Refit(); - - // keep the best % percentile - void CutOff( FaceType*f); - - // restart from the best fitting face - void Restart(); - - void Decorate(); - -}; - - -template < class RegionType> -class RegionGrower{ -public: - RegionGrower():lastAdded(NULL),ave_error(-1),ave_var(-1),n_steps(0){} - typedef typename RegionType RegionType; - typedef typename RegionType::MeshType MeshType; - typedef typename RegionType::FaceType FaceType; - typedef typename FaceType::ScalarType ScalarType; - typedef typename std::list ::iterator TriRegIterator; - typedef typename std::list::iterator AdjIterator; - typedef typename RegionType::FaceError FaceError; - - std::list regions; - std::vector workingset; - int n_faces,target_max_regions; - FaceType * lastAdded; - ScalarType ave_error // average error (over single regions' error) - ,ave_var// average variance (over single regions' error) - ,changed // faces that have changed from previous step ([0..1)) - ,err - ,target_error; // taget error - - FaceType * worst; - ScalarType worst_err ; - MeshType * m; - - typename MeshType:: template PerFaceAttributeHandle attr_r; - typename MeshType:: template PerFaceAttributeHandle attr_r_old; - - unsigned int n_steps;// number of steps done - - struct CandiFace{ - FaceType * f; - float val; - RegionType * r; - char inPlane; // a quale piano di r - CandiFace(){}; - CandiFace(FaceType * _f,double _val,RegionType * _r):f(_f),val(_val),r(_r){} - const bool operator <(const CandiFace & o) const {return val < o.val;} - }; - - std::vector facesheap; - std::vector faceserr; - - struct ErrorEval { - - void Add(const float & v){ - if(!ns){ - if(v < samples[i_min]) i_min = 0; else {++i_min; - if(v > samples[i_max]) i_max = 0;else ++i_max; - } - if(i_min == samples.size()) {i_min = 0; for(int i= 0; i < samples.size(); ++i) if(samples[i]samples[i_max]) i_max = i; } - } - samples.pop_back();samples.push_front(v); - boxes.pop_back(); boxes.push_front(vcg::Point2f(samples[i_min],samples[i_max])); - ++ns ; - } - - float BoxOverlap(){ - float maxsize = std::max( boxes.back()[1]-boxes.back()[0], (*boxes.begin())[1]-(*boxes.begin())[0]); - float overlap = std::max(0.f, std::min(boxes.back()[1],(*boxes.begin())[1])-std::max(boxes.back()[0],(*boxes.begin())[0])); - assert(overlap <= maxsize); - return (maxsize > 0.f)?overlap / maxsize:0.0; - } - float RelativeDecrease(){ - if (ns<2) return std::numeric_limits::max(); - return (samples[0]<10e-22)?0.0:(samples[1]-samples[0])/samples[0]; - } - - void Init(int n ){ - - samples.clear(); - boxes.clear(); - for(int i = 0 ; i < n; ++i) samples.push_front(std::numeric_limits::max()); - for(int i = 0 ; i < n; ++i) boxes.push_front(vcg::Point2f(n,n-i)); - i_max = i_min = 0; - ns = 0; - } - - private: - int i_min,i_max; // index of min and max element in the queue - std::deque samples; - std::deque boxes; - int ns; - - }; - - ErrorEval erroreval; - // init - void Init(MeshType & mesh, int n_seeds,int max_regions, float max_err){ - erroreval.Init(10); - m = &mesh; - vcg::tri::Allocator::CompactFaceVector(*m); - vcg::tri::UpdateTopology::FaceFace(*m); - - float area = 0.f; - for(typename MeshType::FaceIterator fi = m->face.begin(); fi != m->face.end(); ++fi) - area += vcg::DoubleArea(*fi); - area/=2.f; - target_error = area*max_err*max_err; - target_max_regions = max_regions; - - // tte an attibute that will store the address in ocme for the vertex - attr_r = vcg::tri::Allocator::template GetPerFaceAttribute (*m,"r"); - if(!vcg::tri::Allocator::IsValidHandle(*m,attr_r)) - attr_r = vcg::tri::Allocator::template AddPerFaceAttribute (*m,"r"); - - attr_r_old = vcg::tri::Allocator::template GetPerFaceAttribute (*m,"r_old"); - if(!vcg::tri::Allocator::IsValidHandle(*m,attr_r_old)) - attr_r_old = vcg::tri::Allocator::template AddPerFaceAttribute (*m,"r_old"); - - regions.clear(); - - for(int i = 0; i < m->face.size(); ++i){ - attr_r[i] = NULL; - attr_r_old[i] = NULL; - if( (i%(m->fn/n_seeds))==0) - CreateRegion(&m->face[i]); - - } - } - - // add a region - void AddRegion(const RegionType & r){regions.push_back(r);} - - // remove a region - void DeleteRegion(const typename std::list::iterator & ri){std::remove(ri);} - - // initialize a region - void CreateRegion(FaceType * fi){ - AddRegion(RegionType()); - RegionType & tr =regions.back(); - AddFaceToRegion(tr,fi); - tr.Refit(); - tr.color = vcg::Color4b::Scatter(2000,(int)regions.size()); - } - void AddFaceToRegion( RegionType & r, FaceType * f){ - r.face.push_back(f); - attr_r[*f] = (int*) &r; - if(attr_r[*f]!=attr_r_old[*f]) ++r.changed; - attr_r_old[*f] = attr_r[*f]; - ++r.size; - } - - void Merge(RegionType & r0,RegionType & r1){ - assert(!r1.isd); - typename RegionType::FaceIterator fi; - AdjIterator ai; - for(fi = r1.face.begin();fi != r1.face.end(); ++fi) - AddFaceToRegion(r0,(*fi)); - for(ai= r1.adj.begin(); ai != r1.adj.end();++ai) - if( !(*ai)->isd && (*ai)!=&r0) - r0.nx_adj.push_back(*ai); - - r1.face.clear(); - r1.isd = true; - r0.Refit(); - } - - -void PushHeap(std::vector & candi, RegionType & r){ - typename std::vector::iterator ci; - for(ci = candi.begin(); ci != candi.end(); ++ci) - { - facesheap.push_back(CandiFace( *ci,-r.Evaluate(*(*ci)), &r)); - push_heap(facesheap.begin(),facesheap.end()); - } -} - -// the two regions are adjacent - void Connect(RegionType * r1,RegionType * r2){ - assert(r1!=r2); - r1->Connect(r2); - r2->Connect(r1); - } - -// for each region take the candidates and fill in facesheap - void Refill( ){ - facesheap.clear(); - - typename std::list::iterator ri; - - std::vector candi; - for(ri = regions.begin(); ri != regions.end(); ++ri) if(!(*ri).isd) - { - candi.clear(); - Candidates((*ri),candi); - PushHeap(candi,*ri); - } - std::make_heap(facesheap.begin(),facesheap.end()); - } - - - void Candidates(RegionType & r, std::vector< typename RegionType::FaceType*> & c){ - typename RegionType::FaceIterator fi; - - for(fi = r.face.begin(); fi!= r.face.end(); ++fi) - for(int i = 0; i < 3; ++i) - if( ((*fi)->FFp(i) != (*fi) ) && - (attr_r[(*fi)->FFp(i)] != (int*) &r) ) - c.push_back((*fi)->FFp(i)); - } - - bool IsRelaxed(){ - ScalarType _ave_error=0; - ScalarType _ave_var= 0; - ScalarType _changed = 0.0; - int nr=0; - - typename std::list::iterator ri; - worst=NULL;; - worst_err = -1; - qDebug("working set size: %d\n",workingset.size()); - for(ri = regions.begin(); ri != regions.end(); ++ri) if(!(*ri).isd){ - ++nr; - (*ri).UpdateError(); - _ave_error+=(*ri).approx_err; - _ave_var+=(*ri).approx_var; - _changed+=(*ri).changed; - (*ri).changed=0; - if((*ri).worst.val*(*ri).size > worst_err){ - worst = (*ri).worst.f; - worst_err = (*ri).worst.val*(*ri).size; - } - } - - - _ave_error/=nr; - _ave_var/=nr; - _changed/=n_faces; - - erroreval.Add(_ave_error); - qDebug("Err:%f ov: %f-------",_ave_error,erroreval.BoxOverlap()); - - - return (erroreval.BoxOverlap() > 0.5) || (erroreval.RelativeDecrease() < 0.1); - } - -void GrowStepOnce(){ - CandiFace cf; - typename std::list::iterator ri; - for(ri = regions.begin(); ri != regions.end(); ++ri) if(!(*ri).isd) - (*ri).Clean(); - - if(!facesheap.empty()){ - - std::vector toAdd; - std::pop_heap(facesheap.begin(),facesheap.end()); - cf = facesheap.back(); -// qDebug("err:%f\n",cf.val); - facesheap.pop_back(); - if (attr_r[*cf.f]==NULL){ - lastAdded = cf.f; - AddFaceToRegion(*cf.r,cf.f); // adds to region - for(int i=0; i < 3;++i) // put the adjacent in the set of faces to possibly add - if( (attr_r[cf.f->FFp(i)] == NULL) ) - toAdd.push_back(cf.f->FFp(i)); - PushHeap(toAdd,*cf.r); - } - else - { - int * aa = attr_r[*cf.f] ; - if ( attr_r[*cf.f] != (int*)(cf.r) ) - Connect((RegionType*)attr_r[*cf.f],(RegionType*)cf.r); - } - } - int h = (int)facesheap.size(); - // qDebug("----> %d\n",h); - } - - void GrowStep(){ - CandiFace cf; - typename std::list::iterator ri; - - n_steps++; - for(ri = regions.begin(); ri != regions.end(); ++ri) if(!(*ri).isd) - (*ri).Clean(); - - while(!facesheap.empty() ){ - std::vector toAdd; - std::pop_heap(facesheap.begin(),facesheap.end()); - cf = facesheap.back(); - //printf("err:%f\n",cf.val); - facesheap.pop_back(); - if (attr_r[*cf.f]==NULL){ - AddFaceToRegion( *cf.r,cf.f); // ads to region - lastAdded = &*cf.f; - for(int i=0; i < 3;++i) // put the adjacent in the set of faces to possibly add - if( attr_r[ *cf.f->FFp(i)] == NULL ) - toAdd.push_back(cf.f->FFp(i)); - else - if(attr_r[*cf.f->FFp(i)] != attr_r[*cf.f]) - Connect((RegionType*)attr_r[*cf.f->FFp(i)],(RegionType*)cf.r); - - PushHeap(toAdd,*cf.r); - } - else - { - if ( attr_r[*cf.f] != (int*)(cf.r) ) - Connect((RegionType*)attr_r[*cf.f],(RegionType*)cf.r); - } - } - int h = (int) facesheap.size(); - printf("----> %d\n",h); - } - - unsigned int MergeStep(){ - TriRegIterator tri; - unsigned int merged = 0; - typename RegionType::AdjIterator ai; - - for(tri = regions.begin(); tri != regions.end(); ++tri)if(!(*tri).isd) - for(ai = (*tri).adj.begin(); ai != (*tri).adj.end(); ++ai) if(!(*ai)->isd) - if((*tri).Mergeable(*(*ai))){ - Merge((*tri),*(*ai)); - merged++; - } - - - for(tri = regions.begin(); tri != regions.end(); ++tri){ - for(ai = (*tri).nx_adj.begin(); ai != (*tri).nx_adj.end();++ai) - if(!(*ai)->isd) - (*tri).adj.push_back(*ai); - (*tri).adj.sort(); - (*tri).adj.unique(); - } - return merged; - } - - //void TeleportStep(){ - // TriRegIterator tri; - // typename RegionType::AdjIterator ai; - - // for(tri = regions.begin(); tri != regions.end(); ++tri)if(!(*tri).isd){ - // std::vector & onBorder - // ComputeOnBorder(onBorder); - // if(onBorder.empty - // } - // - //} - - bool Restart(){ - std::vector candi; - TriRegIterator ri; - facesheap.clear(); - - if(IsRelaxed()){ - if( (worst_err <= target_error) || (regions.size() >= target_max_regions)) - return false; - else - { - erroreval.Init(10); - FaceError wrs; - wrs.f = worst; - wrs.val = worst_err; - printf("worst triangle error %f\n",worst_err); - - CreateRegion(wrs.f);// CreateRegion changes wr->r - - // reset state variables - ave_error=-1; - ave_var=-1; - err=0.0; - changed=0; - } - } - - for(ri = regions.begin(); ri != regions.end(); ) - if((*ri).isd) - ri = regions.erase(ri); - else - ++ri; - - for(ri = regions.begin(); ri != regions.end(); ++ri) - { - candi.clear(); - (*ri).Refit(); // fit a plane to the region - Restart(*ri); // clear stuff in the region, move the seed to the best fitting to the plabne - Candidates(*ri,candi); // take the (three) faces candidatees - PushHeap(candi,(*ri)); // put the faces on to the heap - } - return true; - } - - - void Restart(RegionType &r){ - if(!r.face.empty()){ - r.size=0; - float b_err = r.Evaluate(*(*r.face.begin())),err; - FaceType* b_face =(*r.face.begin()); - typename RegionType::FaceIterator fi; - - for( fi = r.face.begin(); fi != r.face.end(); ++fi) - { - err = r.Evaluate(**fi); - if(err < b_err) - { - b_err = err; - b_face = *fi; - } - } - for( fi = r.face.begin(); fi != r.face.end(); ++fi) attr_r[(*fi)] = NULL; - r.face.clear(); - r.adj.clear(); - AddFaceToRegion(r,b_face); - r.Refit(); - } - } - - void MakeCharts(){ +#ifndef VCGLIB_REGION_GROWING +#define VCGLIB_REGION_GROWING + +#include +#include + +#include + +template +struct FaceError{ + FaceType * f; + float val; + FaceError(){}; + FaceError(FaceType * _f,double _val):f(_f),val(_val){} + const bool operator <(const FaceError & o) const {return val < o.val;} +}; + + +template +struct Region{ + Region():isd(false){} + typedef Region RegionType; + typedef typename MeshType MeshType; + typedef typename MeshType::FaceType FaceType; + typedef typename MeshType::FaceType::CoordType CoordType; + typedef typename MeshType::FaceType::ScalarType ScalarType; + typedef typename std::vector::iterator FaceIterator; + typedef typename std::list::iterator AdjIterator; + typedef FaceError FaceError; + + // adjacent regions + std::list adj,nx_adj; + + + // some flag + bool isd; + + // evaluate the gain if the triangle is added to the region + ScalarType Evaluate( FaceType & f); + + // update the approximation error of the region + void UpdateError( ); + + // check if two regions are mergeable in one + // i.e. if they have the same planes and coincide i non extreme + ScalarType Mergeable( RegionType & tr); + + // set the region as adjacent to this one + void Connect( RegionType * tr); + + // compute the faces which are on the border + void ComputeOnBorder( std::vector & onBorder); + + // clean the list of adjacencies + void Clean(); + + // refit the plane + void Refit(); + + // keep the best % percentile + void CutOff( FaceType*f); + + // restart from the best fitting face + void Restart(); + + void Decorate(); + +}; + + +template < class RegionType> +class RegionGrower{ +public: + RegionGrower():lastAdded(NULL),ave_error(-1),ave_var(-1),n_steps(0){} + typedef typename RegionType RegionType; + typedef typename RegionType::MeshType MeshType; + typedef typename RegionType::FaceType FaceType; + typedef typename FaceType::ScalarType ScalarType; + typedef typename std::list ::iterator TriRegIterator; + typedef typename std::list::iterator AdjIterator; + typedef typename RegionType::FaceError FaceError; + + std::list regions; + std::vector workingset; + int n_faces,target_max_regions; + FaceType * lastAdded; + ScalarType ave_error // average error (over single regions' error) + ,ave_var// average variance (over single regions' error) + ,changed // faces that have changed from previous step ([0..1)) + ,err + ,target_error; // taget error + + FaceType * worst; + ScalarType worst_err ; + MeshType * m; + + typename MeshType:: template PerFaceAttributeHandle attr_r; + typename MeshType:: template PerFaceAttributeHandle attr_r_old; + + unsigned int n_steps;// number of steps done + + struct CandiFace{ + FaceType * f; + float val; + RegionType * r; + char inPlane; // a quale piano di r + CandiFace(){}; + CandiFace(FaceType * _f,double _val,RegionType * _r):f(_f),val(_val),r(_r){} + const bool operator <(const CandiFace & o) const {return val < o.val;} + }; + + std::vector facesheap; + std::vector faceserr; + + struct ErrorEval { + + void Add(const float & v){ + if(!ns){ + if(v < samples[i_min]) i_min = 0; else {++i_min; + if(v > samples[i_max]) i_max = 0;else ++i_max; + } + if(i_min == samples.size()) {i_min = 0; for(int i= 0; i < samples.size(); ++i) if(samples[i]samples[i_max]) i_max = i; } + } + samples.pop_back();samples.push_front(v); + boxes.pop_back(); boxes.push_front(vcg::Point2f(samples[i_min],samples[i_max])); + ++ns ; + } + + float BoxOverlap(){ + float maxsize = std::max( boxes.back()[1]-boxes.back()[0], (*boxes.begin())[1]-(*boxes.begin())[0]); + float overlap = std::max(0.f, std::min(boxes.back()[1],(*boxes.begin())[1])-std::max(boxes.back()[0],(*boxes.begin())[0])); + assert(overlap <= maxsize); + return (maxsize > 0.f)?overlap / maxsize:0.0; + } + float RelativeDecrease(){ + if (ns<2) return std::numeric_limits::max(); + return (samples[0]<10e-22)?0.0:(samples[1]-samples[0])/samples[0]; + } + + void Init(int n ){ + + samples.clear(); + boxes.clear(); + for(int i = 0 ; i < n; ++i) samples.push_front(std::numeric_limits::max()); + for(int i = 0 ; i < n; ++i) boxes.push_front(vcg::Point2f(n,n-i)); + i_max = i_min = 0; + ns = 0; + } + + private: + int i_min,i_max; // index of min and max element in the queue + std::deque samples; + std::deque boxes; + int ns; + + }; + + ErrorEval erroreval; + // init + void Init(MeshType & mesh, int n_seeds,int max_regions, float max_err){ + erroreval.Init(10); + m = &mesh; + vcg::tri::Allocator::CompactFaceVector(*m); + vcg::tri::UpdateTopology::FaceFace(*m); + + float area = 0.f; + for(typename MeshType::FaceIterator fi = m->face.begin(); fi != m->face.end(); ++fi) + area += vcg::DoubleArea(*fi); + area/=2.f; + target_error = area*max_err*max_err; + target_max_regions = max_regions; + + // tte an attibute that will store the address in ocme for the vertex + attr_r = vcg::tri::Allocator::template GetPerFaceAttribute (*m,"r"); + if(!vcg::tri::Allocator::IsValidHandle(*m,attr_r)) + attr_r = vcg::tri::Allocator::template AddPerFaceAttribute (*m,"r"); + + attr_r_old = vcg::tri::Allocator::template GetPerFaceAttribute (*m,"r_old"); + if(!vcg::tri::Allocator::IsValidHandle(*m,attr_r_old)) + attr_r_old = vcg::tri::Allocator::template AddPerFaceAttribute (*m,"r_old"); + + regions.clear(); + + for(int i = 0; i < m->face.size(); ++i){ + attr_r[i] = NULL; + attr_r_old[i] = NULL; + if( (i%(m->fn/n_seeds))==0) + CreateRegion(&m->face[i]); + + } + } + + // add a region + void AddRegion(const RegionType & r){regions.push_back(r);} + + // remove a region + void DeleteRegion(const typename std::list::iterator & ri){std::remove(ri);} + + // initialize a region + void CreateRegion(FaceType * fi){ + AddRegion(RegionType()); + RegionType & tr =regions.back(); + AddFaceToRegion(tr,fi); + tr.Refit(); + tr.color = vcg::Color4b::Scatter(2000,(int)regions.size()); + } + void AddFaceToRegion( RegionType & r, FaceType * f){ + r.face.push_back(f); + attr_r[*f] = (int*) &r; + if(attr_r[*f]!=attr_r_old[*f]) ++r.changed; + attr_r_old[*f] = attr_r[*f]; + ++r.size; + } + + void Merge(RegionType & r0,RegionType & r1){ + assert(!r1.isd); + typename RegionType::FaceIterator fi; + AdjIterator ai; + for(fi = r1.face.begin();fi != r1.face.end(); ++fi) + AddFaceToRegion(r0,(*fi)); + for(ai= r1.adj.begin(); ai != r1.adj.end();++ai) + if( !(*ai)->isd && (*ai)!=&r0) + r0.nx_adj.push_back(*ai); + + r1.face.clear(); + r1.isd = true; + r0.Refit(); + } + + +void PushHeap(std::vector & candi, RegionType & r){ + typename std::vector::iterator ci; + for(ci = candi.begin(); ci != candi.end(); ++ci) + { + facesheap.push_back(CandiFace( *ci,-r.Evaluate(*(*ci)), &r)); + push_heap(facesheap.begin(),facesheap.end()); + } +} + +// the two regions are adjacent + void Connect(RegionType * r1,RegionType * r2){ + assert(r1!=r2); + r1->Connect(r2); + r2->Connect(r1); + } + +// for each region take the candidates and fill in facesheap + void Refill( ){ + facesheap.clear(); + + typename std::list::iterator ri; + + std::vector candi; + for(ri = regions.begin(); ri != regions.end(); ++ri) if(!(*ri).isd) + { + candi.clear(); + Candidates((*ri),candi); + PushHeap(candi,*ri); + } + std::make_heap(facesheap.begin(),facesheap.end()); + } + + + void Candidates(RegionType & r, std::vector< typename RegionType::FaceType*> & c){ + typename RegionType::FaceIterator fi; + + for(fi = r.face.begin(); fi!= r.face.end(); ++fi) + for(int i = 0; i < 3; ++i) + if( ((*fi)->FFp(i) != (*fi) ) && + (attr_r[(*fi)->FFp(i)] != (int*) &r) ) + c.push_back((*fi)->FFp(i)); + } + + bool IsRelaxed(){ + ScalarType _ave_error=0; + ScalarType _ave_var= 0; + ScalarType _changed = 0.0; + int nr=0; + + typename std::list::iterator ri; + worst=NULL;; + worst_err = -1; + qDebug("working set size: %d\n",workingset.size()); + for(ri = regions.begin(); ri != regions.end(); ++ri) if(!(*ri).isd){ + ++nr; + (*ri).UpdateError(); + _ave_error+=(*ri).approx_err; + _ave_var+=(*ri).approx_var; + _changed+=(*ri).changed; + (*ri).changed=0; + if((*ri).worst.val*(*ri).size > worst_err){ + worst = (*ri).worst.f; + worst_err = (*ri).worst.val*(*ri).size; + } + } + + + _ave_error/=nr; + _ave_var/=nr; + _changed/=n_faces; + + erroreval.Add(_ave_error); + qDebug("Err:%f ov: %f-------",_ave_error,erroreval.BoxOverlap()); + + + return (erroreval.BoxOverlap() > 0.5) || (erroreval.RelativeDecrease() < 0.1); + } + +void GrowStepOnce(){ + CandiFace cf; + typename std::list::iterator ri; + for(ri = regions.begin(); ri != regions.end(); ++ri) if(!(*ri).isd) + (*ri).Clean(); + + if(!facesheap.empty()){ + + std::vector toAdd; + std::pop_heap(facesheap.begin(),facesheap.end()); + cf = facesheap.back(); +// qDebug("err:%f\n",cf.val); + facesheap.pop_back(); + if (attr_r[*cf.f]==NULL){ + lastAdded = cf.f; + AddFaceToRegion(*cf.r,cf.f); // adds to region + for(int i=0; i < 3;++i) // put the adjacent in the set of faces to possibly add + if( (attr_r[cf.f->FFp(i)] == NULL) ) + toAdd.push_back(cf.f->FFp(i)); + PushHeap(toAdd,*cf.r); + } + else + { + int * aa = attr_r[*cf.f] ; + if ( attr_r[*cf.f] != (int*)(cf.r) ) + Connect((RegionType*)attr_r[*cf.f],(RegionType*)cf.r); + } + } + int h = (int)facesheap.size(); + // qDebug("----> %d\n",h); + } + + void GrowStep(){ + CandiFace cf; + typename std::list::iterator ri; + + n_steps++; + for(ri = regions.begin(); ri != regions.end(); ++ri) if(!(*ri).isd) + (*ri).Clean(); + + while(!facesheap.empty() ){ + std::vector toAdd; + std::pop_heap(facesheap.begin(),facesheap.end()); + cf = facesheap.back(); + //printf("err:%f\n",cf.val); + facesheap.pop_back(); + if (attr_r[*cf.f]==NULL){ + AddFaceToRegion( *cf.r,cf.f); // ads to region + lastAdded = &*cf.f; + for(int i=0; i < 3;++i) // put the adjacent in the set of faces to possibly add + if( attr_r[ *cf.f->FFp(i)] == NULL ) + toAdd.push_back(cf.f->FFp(i)); + else + if(attr_r[*cf.f->FFp(i)] != attr_r[*cf.f]) + Connect((RegionType*)attr_r[*cf.f->FFp(i)],(RegionType*)cf.r); + + PushHeap(toAdd,*cf.r); + } + else + { + if ( attr_r[*cf.f] != (int*)(cf.r) ) + Connect((RegionType*)attr_r[*cf.f],(RegionType*)cf.r); + } + } + int h = (int) facesheap.size(); + printf("----> %d\n",h); + } + + unsigned int MergeStep(){ + TriRegIterator tri; + unsigned int merged = 0; + typename RegionType::AdjIterator ai; + + for(tri = regions.begin(); tri != regions.end(); ++tri)if(!(*tri).isd) + for(ai = (*tri).adj.begin(); ai != (*tri).adj.end(); ++ai) if(!(*ai)->isd) + if((*tri).Mergeable(*(*ai))){ + Merge((*tri),*(*ai)); + merged++; + } + + + for(tri = regions.begin(); tri != regions.end(); ++tri){ + for(ai = (*tri).nx_adj.begin(); ai != (*tri).nx_adj.end();++ai) + if(!(*ai)->isd) + (*tri).adj.push_back(*ai); + (*tri).adj.sort(); + (*tri).adj.unique(); + } + return merged; + } + + //void TeleportStep(){ + // TriRegIterator tri; + // typename RegionType::AdjIterator ai; + + // for(tri = regions.begin(); tri != regions.end(); ++tri)if(!(*tri).isd){ + // std::vector & onBorder + // ComputeOnBorder(onBorder); + // if(onBorder.empty + // } + // + //} + + bool Restart(){ + std::vector candi; + TriRegIterator ri; + facesheap.clear(); + + if(IsRelaxed()){ + if( (worst_err <= target_error) || (regions.size() >= target_max_regions)) + return false; + else + { + erroreval.Init(10); + FaceError wrs; + wrs.f = worst; + wrs.val = worst_err; + printf("worst triangle error %f\n",worst_err); + + CreateRegion(wrs.f);// CreateRegion changes wr->r + + // reset state variables + ave_error=-1; + ave_var=-1; + err=0.0; + changed=0; + } + } + + for(ri = regions.begin(); ri != regions.end(); ) + if((*ri).isd) + ri = regions.erase(ri); + else + ++ri; + + for(ri = regions.begin(); ri != regions.end(); ++ri) + { + candi.clear(); + (*ri).Refit(); // fit a plane to the region + Restart(*ri); // clear stuff in the region, move the seed to the best fitting to the plabne + Candidates(*ri,candi); // take the (three) faces candidatees + PushHeap(candi,(*ri)); // put the faces on to the heap + } + return true; + } + + + void Restart(RegionType &r){ + if(!r.face.empty()){ + r.size=0; + float b_err = r.Evaluate(*(*r.face.begin())),err; + FaceType* b_face =(*r.face.begin()); + typename RegionType::FaceIterator fi; + + for( fi = r.face.begin(); fi != r.face.end(); ++fi) + { + err = r.Evaluate(**fi); + if(err < b_err) + { + b_err = err; + b_face = *fi; + } + } + for( fi = r.face.begin(); fi != r.face.end(); ++fi) attr_r[(*fi)] = NULL; + r.face.clear(); + r.adj.clear(); + AddFaceToRegion(r,b_face); + r.Refit(); + } + } + + void MakeCharts(){ this->Refill(); while(this->Restart()){ //do{ @@ -487,7 +487,7 @@ void GrowStepOnce(){ } //while(this->MergeStep()); } - -}; - -#endif + +}; + +#endif diff --git a/src/meshlabplugins/filter_vsa/vsa.cpp b/src/meshlabplugins/filter_vsa/vsa.cpp index b5c65b0b6..55c433dfd 100755 --- a/src/meshlabplugins/filter_vsa/vsa.cpp +++ b/src/meshlabplugins/filter_vsa/vsa.cpp @@ -35,12 +35,7 @@ $Log: samplefilter.cpp,v $ #include "region_growing.h" #include "planar_region.h" -#include -#include -#include -#include -#include -#include +#include using namespace vcg; diff --git a/src/meshlabplugins/filter_zippering/filter_zippering.cpp b/src/meshlabplugins/filter_zippering/filter_zippering.cpp index 072c3fa82..8db071ccd 100644 --- a/src/meshlabplugins/filter_zippering/filter_zippering.cpp +++ b/src/meshlabplugins/filter_zippering/filter_zippering.cpp @@ -1,1702 +1,1702 @@ -/**************************************************************************** -* MeshLab o o * -* A versatile mesh processing toolbox o o * -* _ O _ * -* Copyright(C) 2005 \/)\/ * -* Visual Computing Lab /\/| * -* ISTI - Italian National Research Council | * -* \ * -* All rights reserved. * -* * -* This program is free software; you can redistribute it and/or modify * -* it under the terms of the GNU General Public License as published by * -* the Free Software Foundation; either version 2 of the License, or * -* (at your option) any later version. * -* * -* This program is distributed in the hope that it will be useful, * -* but WITHOUT ANY WARRANTY; without even the implied warranty of * -* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * -* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) * -* for more details. * -* * -****************************************************************************/ -/**************************************************************************** - History - - -****************************************************************************/ - -#include -#include "filter_zippering.h" - -#include -#include -#include -#include -#include -#include -#include -#include - - -#include - -using namespace vcg; -using namespace std; - -#define MAX_LOOP 150 -//function added to map call sto distance -template -ScalarType FilterZippering::SquaredDistance( vcg::Segment3 &s, vcg::Point3 &p) -{ - Point3 clos; - ScalarType dist; - vcg::SegmentPointDistance(s, p, clos,dist ); - return dist; -} - -template -vcg::Point3 FilterZippering::ClosestPoint( vcg::Segment3 &s, vcg::Point3 &p) -{ - Point3 clos; - ScalarType dist; - vcg::SegmentPointDistance(s, p, clos,dist ); - return clos; -} - -// Constructor usually performs only two simple tasks of filling the two lists -// - typeList: with all the possible id of the filtering actions -// - actionList with the corresponding actions. If you want to add icons to your filtering actions you can do here by construction the QActions accordingly - -FilterZippering::FilterZippering() -{ - typeList << FP_REDUNDANCY << FP_ZIPPERING; - - foreach(FilterIDType tt , types()) actionList << new QAction(filterName(tt), this); -} - -// ST() must return the very short string describing each filtering action -// (this string is used also to define the menu entry) - QString FilterZippering::filterName(FilterIDType filterId) const -{ - switch(filterId) { - case FP_REDUNDANCY : return QString("Select Redundant Faces"); - case FP_ZIPPERING : return QString("Zippering"); - default : assert(0); - } -} - -// Info() must return the longer string describing each filtering action -// (this string is used in the About plugin dialog) - QString FilterZippering::filterInfo(FilterIDType filterId) const -{ - switch(filterId) { - case FP_REDUNDANCY: return QString("Remove redundant faces from one mesh or from both of them, starting from borders."); - case FP_ZIPPERING : return QString("Merge two triangle meshes into a single one. This method doesn't provide check on redundancy. Based on Controlledand Adaptive Mesh Zippering, by S.Marras, F.Ganovelli, P.Cignoni."); - default : assert(0); - } - return QString("Unknown Filter"); -} - - int FilterZippering::getRequirements(QAction *action) -{ - switch(ID(action)) - { - case FP_REDUNDANCY: return MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEMARK | MeshModel::MM_VERTCOLOR | MeshModel::MM_VERTQUALITY; - case FP_ZIPPERING : return MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEMARK | MeshModel::MM_VERTCOLOR | MeshModel::MM_VERTQUALITY; - default: assert(0); - } - return 0; -} - -// The FilterClass describes in which generic class of filters it fits. -// This choice affect the submenu in which each filter will be placed -// More than a single class can be choosen. - FilterZippering::FilterClass FilterZippering::getClass(QAction *a) -{ - switch(ID(a)) - { - case FP_REDUNDANCY: return MeshFilterInterface::Selection; - case FP_ZIPPERING : return MeshFilterInterface::Remeshing; - default : assert(0); - } - return MeshFilterInterface::Generic; -} - -// This function define the needed parameters for each filter. Return true if the filter has some parameters -// it is called every time, so you can set the default value of parameters according to the mesh -// For each parmeter you need to define, -// - the name of the parameter, -// - the string shown in the dialog -// - the default value -// - a possibly long string describing the meaning of that parameter (shown as a popup help in the dialog) -void FilterZippering::initParameterSet(QAction *action, MeshDocument &md, RichParameterSet & parlst) -{ - MeshModel *target; - float maxVal = 0.0; - switch(ID(action)) { - case FP_REDUNDANCY: - //get diagonals - foreach (target, md.meshList) { - if ( target->cm.bbox.Diag() > maxVal ) maxVal = target->cm.bbox.Diag(); - if (target != md.mm()) break; - } - parlst.addParam( new RichMesh("FirstMesh", md.mm(), &md, "Source Mesh", "The mesh with holes.") ); - parlst.addParam( new RichMesh("SecondMesh", md.mm(), &md, "Target Mesh", "The mesh that will be used as patch.") ); - parlst.addParam( new RichAbsPerc("distance", maxVal*0.01, 0, maxVal, "Max distance", "Max distance between mesh and path") ); - parlst.addParam( new RichBool("UseQuality", false, "Use quality to select redundant face", "If selected, previously computed face quality will be used in order to select redundant faces.") ); - //parlst.addParam( new RichBool("FastErosion", false, "Use fast erosion algorithm", "If selected, improves the speed-up of algorithm (results may be not accurate). Useful for large meshes.") ); - parlst.addParam( new RichBool("FullProcessing", false, "Process the whole Target Mesh", "If selected, redundancy test is performed over the whole surface of the mesh") ); - break; - - case FP_ZIPPERING : - //get diagonals - foreach (target, md.meshList) { - if ( target->cm.bbox.Diag() > maxVal ) maxVal = target->cm.bbox.Diag(); - if (target != md.mm()) break; - } - parlst.addParam( new RichMesh("FirstMesh", md.mm(), &md, "Mesh (with holes)", "The mesh with holes.") ); - parlst.addParam( new RichMesh("SecondMesh", target, &md, "Patch", "The mesh that will be used as patch.") ); - parlst.addParam( new RichAbsPerc("distance", maxVal*0.01, 0, maxVal, "Max distance", "Max distance between mesh and path") ); - break; - default: break; // do not add any parameter for the other filters - } -} - -/* Given the mesh m (let's assume it's mesh with holes) and a face on another mesh (let's assume it's - * the patch), check if the face can be completely projected on to mesh m. "Completely projected" means - * that edges' projection completely lies on one or more faces of m, and does not intersect edge - * border of mesh m. - * @param face The query face (from patch) - * @param m The mesh with holes - * @param grid A face-grid created using the faces of mesh m - * @param max_dist Max Distance allowed between m and patch; if distance between face and m is higher than max_dist, face will be discarded - * @return true if face can be completely projected on m (redundant), false otherwise. - */ -bool FilterZippering::checkRedundancy( CMeshO::FacePointer face, - MeshModel *m, - MeshFaceGrid &grid, - CMeshO::ScalarType max_dist ) -{ - // Step1: check if border edge can be projected on m - - //search for border edge - //edge adjacent to a selected face can be considered as border edge - int i; - for ( i = 0; i < 3; i ++ ) { - if ( face::IsBorder(*face, i) ) break; - if ( face->FFp(i)->IsS() ) break; - } - //no border edge find; check edge 0 - if ( i == 3 ) i = 0; - - size_t samplePerEdge = SAMPLES_PER_EDGE; - - //samples edge in uniform way - vector< Point3< CMeshO::ScalarType > > edge_samples; - Point3< CMeshO::ScalarType > edge_dir = face->P1(i) - face->P(i); edge_dir.Normalize(); - float step = 1.0/(samplePerEdge+1); //step length - for( size_t j = 0; j <= samplePerEdge; j++ ) { - edge_samples.push_back( face->P(i) + edge_dir * (j * step) ); - } - - for ( unsigned int j = 0; j < edge_samples.size(); j ++ ) { - CMeshO::FacePointer nearestF = 0; - tri::FaceTmark markerFunctor; markerFunctor.SetMesh(&m->cm); tri::UnMarkAll(m->cm); - face::PointDistanceBaseFunctor PDistFunct; - MeshFaceGrid::ScalarType dist = max_dist; MeshFaceGrid::CoordType closest; - //Search closest point on A - nearestF = grid.GetClosest(PDistFunct, markerFunctor, edge_samples[j], max_dist, dist, closest); - - if ( nearestF == 0 ) return false; //no face within given range - if ( isOnBorder(closest, nearestF ) ) return false; //closest point on border - if ( nearestF->IsD() ) return false; //face is deleted - if ( nearestF->IsS() ) return false; //face is selected (will be deleted) - } - - //check if V2(i) has a closest point on border of m - CMeshO::FacePointer nearestF = 0; - tri::FaceTmark markerFunctor; markerFunctor.SetMesh(&m->cm); - face::PointDistanceBaseFunctor PDistFunct; - MeshFaceGrid::ScalarType dist = max_dist; MeshFaceGrid::CoordType closest; - nearestF = grid.GetClosest(PDistFunct, markerFunctor, face->P2(i), max_dist, dist, closest); - if ( nearestF == 0 ) return false; //no face within given range - if ( isOnBorder( closest, nearestF ) ) return false; //closest point on border - if ( nearestF->IsD() ) return false; //face has been removed - if ( nearestF->IsS() ) return false; //face has been selected - - //check if edges are completely projectable on m - for ( int j = (i+1)%3; j != i; j = (j+1)%3 ) { - samplePerEdge = SAMPLES_PER_EDGE; - edge_samples.clear(); edge_dir = face->P1(j) - face->P(j); edge_dir.Normalize(); - for( size_t k = 0; k <= samplePerEdge; k++ ) edge_samples.push_back( face->P(j) + (face->P1(j) - face->P(j)) * (k * step) ); //campionamento - // samples on A - for ( size_t k = 0; k < edge_samples.size(); k ++ ) { - CMeshO::FacePointer nearestF = 0; - tri::FaceTmark markerFunctor; markerFunctor.SetMesh(&m->cm); tri::UnMarkAll(m->cm); - face::PointDistanceBaseFunctor PDistFunct; - MeshFaceGrid::ScalarType dist = max_dist; MeshFaceGrid::CoordType closest; - //Search closest point on A - nearestF = grid.GetClosest(PDistFunct, markerFunctor, edge_samples[k], max_dist, dist, closest); - if ( nearestF == 0 ) return false; //no face within given range - if ( isOnBorder(closest, nearestF ) ) return false; //closest point on border - if ( nearestF->IsD() ) return false; //face has been removed - if ( nearestF->IsS() ) return false; //face has been selected - } - } - // redundant - return true; -} - -/* - * Perform a simple test on a face in order to check redundancy (F.Ganovelli improvement of redundancy test) - * CheckRedundancy procedure uses a number of sampled points; instead, this function performs redundancy test - * using only the barycenter of face f. Then, after that nearest face is found, function verifies if - * vertices' distance from border is greater than max edge length of f; in this case, f is considered redundant. - * Alternatively, it can be used in order to check if face is too distance from border of M. - * @param face The query face (from patch) - * @param m The mesh with holes (Note that distance from border must be previously calculated) - * @param grid A face-grid created using the faces of mesh m - * @param max_dist Max Distance allowed between m and patch; if distance between face and m is higher than max_dist, face will be discarded - * @param test parameted used to determine the type of test - * @return true if face can be considered redundant, false otherwise -*/ -bool FilterZippering::simpleCheckRedundancy( CMeshO::FacePointer f, //face - MeshModel *m, //mesh A - MeshFaceGrid &grid, //grid A - CMeshO::ScalarType max_dist, - bool test) { //Max search distance - - Point3f qp = Barycenter(*f); //f barycenter - //search for max_edge - float max_edge = max( Distance(f->P(0),f->P(1)), max( Distance(f->P(1),f->P(2)), Distance(f->P(2),f->P(0)) ) ); - float dist = max_dist; CMeshO::FacePointer nearestF = 0; Point3f closest; - tri::FaceTmark markerFunctor; markerFunctor.SetMesh(&m->cm); UnMarkAll(m->cm); - face::PointDistanceBaseFunctor PDistFunct; - nearestF = grid.GetClosest(PDistFunct, markerFunctor, qp, max_dist, dist, closest); - if (nearestF == 0) return false; //too far away - float min_q = min( nearestF->V(0)->Q(), min( nearestF->V(1)->Q(), nearestF->V(2)->Q() ) ); //min distance of nearestF's vertices from M-border - float max_q = max( f->V(0)->Q(), max( f->V(1)->Q(), f->V(2)->Q() ) ); //max distance of F's vertices from A-border - if ( min_q <= max_edge ) return false; - if (test) if ( min_q <= max_q ) return false; - return true; -} - -/* Check if point is on border of face f. - * @param point The query-point - * @param f Face containing point - * @return true if point lies on a border edge or vertex of f, false otherwise - */ -bool FilterZippering::isOnBorder( Point3f point, CMeshO::FacePointer f ) { - // for each edge, calculates distance point-edge - if ( f == 0 ) return false; //null face - - //compute barycentric coords - float bc[3]; - InterpolationParameters( *f, f->N(), point, bc[0], bc[1], bc[2] ); - //search for max and min - int min_el = min_element(bc, bc+3) - bc; - int max_el = max_element(bc, bc+3) - bc; - //coords of max el = 1.0f -> check vertex - if ( bc[max_el] >= 1.0f - eps ) return isBorderVert( f, max_el ); - //coords of min_el = 0.0f -> check edge - if ( bc[min_el] <= 0.0f + eps ) return ( face::IsBorder( *f, (min_el+1)%3 ) || f->FFp((min_el+1)%3)->IsS() ); - - return false; -} - -/* Check if vertex i on face f belong to a border edge. - * @param f face - * @param i index of vertex - * @return true if f->V(i) lies on border-edge, false otherwise - */ -bool FilterZippering::isBorderVert( CMeshO::FacePointer f, int i ) { - - face::Pos p( f, i, f->V(i) ); - //loop - do { - //border - if ( face::IsBorder( *p.F(), p.E() ) ) return true; - //adj to selected face - if ( p.F()->FFp(p.E())->IsS() ) return true; - p.FlipE(); p.FlipF(); - } while(p.F() != f); - return false; -} - -/* Check adjiacency between two faces - * @param f1 first face - * @param f2 second face - * return true if f1 is adjacent to f2, false otherwise - */ -bool FilterZippering::isAdjacent( CMeshO::FacePointer f1, CMeshO::FacePointer f2 ) { - if ( f1 == f2 ) return false; - return ((f1 == f2->FFp(0)) || (f1 == f2->FFp(1)) || (f1 == f2->FFp(2))); -} - -/* Using auxiliar informations, replace a face with a new triangulation, discarding part of the originale face - * and tesselating the other part. - * @param info Auxiliar information from face - * @param N Face normal - * @param coords Output coords of vertices - * @param pointers Output vertex indices, will be used for creation of new faces - */ -void FilterZippering::handleBorder( aux_info &info, //Auxiliar information for triangulatio - Point3f N, //face normal (useful for proiection) - vector &/*coords*/, //output coords - vector &pointers ) { //output triangles - // rotation matrix (will be used for projection on plane) - Matrix44f rot_matrix; - rot_matrix.SetRotateRad( Angle( N, Point3f(0.0, 0.0, 1.0) ), N ^ Point3f(0.0, 0.0, 1.0) ); - //border refinement - brand new method? - //Remove intersecating border - for (size_t i = 0; i < info.border.size(); i ++) { - for ( size_t j = 0; j < info.border[i].edges.size(); j++) { - //project one border segment on face - Segment2f s( Point2f((rot_matrix * info.border[i].edges[j].P0()).X(), (rot_matrix * info.border[i].edges[j].P0()).Y()), - Point2f((rot_matrix * info.border[i].edges[j].P1()).X(), (rot_matrix * info.border[i].edges[j].P1()).Y()) ); //projects edge on plane - Point2f x; - for ( size_t k = i+1; k < info.border.size(); k++) { - for ( size_t h = 0; h < info.border[k].edges.size(); h++) { - Segment2f t( Point2f((rot_matrix * info.border[k].edges[h].P0()).X(), (rot_matrix * info.border[k].edges[h].P0()).Y()), - Point2f((rot_matrix * info.border[k].edges[h].P1()).X(), (rot_matrix * info.border[k].edges[h].P1()).Y()) ); //projects edge on plane - if ( SegmentSegmentIntersection( s, t, x ) ) { - h = info.border[k].edges.size(); - info.border.erase(info.border.begin() + k); - k--; - } - } - } - } - } - - //Split border if necessary - for (size_t i = 0; i < info.border.size(); i ++) { - for ( size_t j = 1; j < info.border[i].edges.size(); j++ ) { - //check if P0 lies on trash edge, then split border (skip first edge) - for ( size_t k = 0; k < info.trash[0].edges.size(); k++ ) { - float EPSILON=info.trash[0].edges[k].Length()/100.f; - vcg::Point3 clos; - float dist; - vcg::SegmentPointDistance( info.trash[0].edges[k],info.border[i].edges[j].P0(),clos,dist); - if (dist( info.trash[0].edges[k], info.border[i].edges[j].P0() ) == 0.0f ) { //approxim - //Split is needed - polyline newborder; newborder.edges.insert( newborder.edges.begin(), info.border[i].edges.begin() + j, info.border[i].edges.end() ); - newborder.verts.insert( newborder.verts.begin(), info.border[i].verts.begin() + j, info.border[i].verts.end() ); - info.border[i].edges.erase( info.border[i].edges.begin() + j, info.border[i].edges.end() ); - info.border[i].verts.erase( info.border[i].verts.begin() + j, info.border[i].verts.end() ); - info.border.push_back( newborder ); - break; - } - } - } - } - - // For each border... - for (size_t i = 0; i < info.border.size(); i ++) { - //search for component intersecated by border and split it into two or more components - bool conn = true; int c = searchComponent( info, info.border[i].edges.front().P0(), info.border[i].edges.back().P1(), conn ); polyline current; - if ( c == -1 ) //no intersection - continue; - if ( conn ) current = info.conn[c]; else current = info.trash[c]; - info.AddCComponent( cutComponent( current, info.border[i], rot_matrix ) ); - polyline rev_border = info.border[i]; reverse( rev_border.edges.begin(), rev_border.edges.end() ); - for ( size_t k = 0; k < rev_border.edges.size(); k ++) rev_border.edges[k].Flip(); - info.AddTComponent( cutComponent( current, rev_border, rot_matrix ) ); - if ( conn ) info.RemoveCComponent( c ); else info.RemoveTComponent( c ); - } - //No border; triangulation of the whole face - if ( info.border.size() == 0 ) { - info.conn = info.trash; - for ( size_t i = 0; i < info.conn.size(); i ++ ) { - reverse( info.conn[i].edges.begin(), info.conn[i].edges.end() ); - for ( size_t j = 0; j < info.conn[i].edges.size(); j ++ ) info.conn[i].edges[j].Flip(); - reverse( info.conn[i].verts.begin(), info.conn[i].verts.end() ); - for ( size_t j = 0; j < info.conn[i].verts.size(); j ++ ) info.conn[i].verts[j] = make_pair(info.conn[i].verts[j].second, info.conn[i].verts[j].first); - } - } - //triangulation of Ccomponent - for ( int i = 0; i < info.nCComponent(); i ++ ) { - vector< Point3f > points; //coords vector - vector< int > vertices; //vertices vector - for ( size_t j = 0; j < info.conn[i].edges.size(); j ++ ) { - points.push_back( info.conn[i].edges[j].P0() ); - vertices.push_back( info.conn[i].verts[j].first ); - } - if ( points.size() < 3 ) continue; - vector< int > indices; int iters = 0; - vector< vector< Point3f > > outlines; outlines.push_back( points ); - while ( indices.size() == 0 && ++iters < MAX_LOOP ) { - glu_tesselator::tesselate( outlines, indices ); //glu tessellator - if ( indices.size() == 0 ) - for ( size_t k = 0; k < outlines[0].size(); k ++ ) outlines[0][k] = outlines[0][k] * 10.0; //glu tessellator doesn't work properly for close points, so we scale coords in order to obtain a triangulation, if needed - } - for ( size_t k = 0; k < indices.size(); k += 3 ) { - if ( (vertices[indices[k]] != vertices[indices[k+1]]) && - (vertices[indices[k]] != vertices[indices[k+2]]) && - (vertices[indices[k+1]] != vertices[indices[k+2]] ) ) { - pointers.push_back( vertices[indices[k]] ); //save indices, in order to create new faces - pointers.push_back( vertices[indices[k+1]] ); - pointers.push_back( vertices[indices[k+2]] ); - } - } - } -} - -/* Split a component into two component, using border as separating line. Discard component external to border polyline. - * @param comp To-be-split component - * @param border Border edges - * @param rot_mat Rotation matrix, needed to project comp points on z=k plane - */ - -polyline FilterZippering::cutComponent( polyline comp, //Component to be cut - polyline border, //border - Matrix44f rot_mat ) { //Rotation matrix - - Point3f startpoint = border.edges.front().P0(); - Point3f endpoint = border.edges.back().P1(); - Point2 startpoint2D ( (rot_mat * startpoint).X(), (rot_mat * startpoint).Y() ); - Point2 endpoint2D ( (rot_mat * endpoint).X(), (rot_mat * endpoint).Y() ); - //int startedge = 0, endedge = 0; float min_dist_s = SquaredDistance( comp.edges[0], startpoint ), min_dist_e = SquaredDistance( comp.edges[0], endpoint ); - int startedge = 0, endedge = 0; - Point3f clos; - float dist; - float min_dist_s; vcg::SegmentPointDistance( comp.edges[0], startpoint, clos,min_dist_s ); - float min_dist_e; vcg::SegmentPointDistance( comp.edges[0], endpoint, clos,min_dist_e ); - bool v_start = false, v_end = false; - // search where startpoint and endpoint lie - for ( size_t i = 0; i < comp.edges.size(); i ++ ) { - if ( !v_start && SquaredDistance( comp.edges[i], startpoint ) <= min_dist_s ) { startedge = i; min_dist_s = SquaredDistance( comp.edges[i], startpoint ); } - if ( !v_end && SquaredDistance( comp.edges[i], endpoint ) <= min_dist_e ) { endedge = i; min_dist_e = SquaredDistance( comp.edges[i], endpoint ); } - if ( comp.edges[i].P1() == startpoint ) { startedge = i; v_start = true; } //lies on a vertex - if ( comp.edges[i].P0() == endpoint ) { endedge = i; v_end = true; } //lies on a vertex - } - polyline p; - // border edges will be edges of new comp - p.edges.insert( p.edges.begin(), border.edges.begin(), border.edges.end() ); - p.verts.insert( p.verts.begin(), border.verts.begin(), border.verts.end() ); - // startedge == endedge - if ( startedge == endedge && !Convex( startpoint2D, Point2 ( (rot_mat * border.edges.front().P1()).X(), (rot_mat * border.edges.front().P1()).Y() ) , endpoint2D ) ) { - Segment3 join( endpoint, startpoint ); - p.edges.push_back( join ); p.verts.push_back( make_pair( border.verts.back().second, border.verts.front().first ) ); //Vertex pointers - return p; - } - - // startedge!=endedge - // search point on the right, create oriented segment and go on - int step = -1; - Point3f c0 = border.edges.back().P0(); - vector< Segment3 >::iterator edge_it = border.edges.end(); edge_it--; - //too short segment; not reliable - while ( Distance( rot_mat * c0, rot_mat * endpoint ) <= 5.0 * eps ) { - //previous - c0 = (*edge_it).P0(); - if (edge_it != border.edges.begin()) edge_it--; else break; - } - if ( v_end ) { - if ( !Convex( Point2 ( (rot_mat * c0).X(), (rot_mat * c0).Y() ), endpoint2D, - Point2 ( (rot_mat * comp.edges[(endedge+comp.edges.size()-1)%comp.edges.size()].P0()).X(), (rot_mat * comp.edges[(endedge+comp.edges.size()-1)%comp.edges.size()].P0()).Y() ) ) ) { - step = comp.edges.size() - 1; - } - else { - step = comp.edges.size() + 1; - Segment3 s( comp.edges[endedge].P0(), comp.edges[endedge].P1() ); - p.edges.push_back( s ); step = comp.edges.size() + 1; - p.verts.push_back( make_pair( comp.verts[endedge].first, comp.verts[endedge].second ) ); - } - } - else { - if ( !Convex( Point2 ( (rot_mat * c0).X(), (rot_mat * c0).Y() ), endpoint2D, - Point2 ( (rot_mat * comp.edges[endedge].P0()).X(), (rot_mat * comp.edges[endedge].P0()).Y() ) ) ) { - Segment3 s( endpoint, comp.edges[endedge].P0() ); - p.edges.push_back( s ); step = comp.edges.size() - 1; - p.verts.push_back( make_pair(border.verts.back().second, comp.verts[endedge].first ) ); - } - else { - Segment3 s( endpoint, comp.edges[endedge].P1() ); - p.edges.push_back( s ); step = comp.edges.size() + 1; - p.verts.push_back( make_pair(border.verts.back().second, comp.verts[endedge].second ) ); - } - } - for ( int i = (endedge + step)%(comp.edges.size()); i != startedge; i = (i + step)%(comp.edges.size()) ) { - p.edges.push_back( comp.edges[i] ); - pair vs( comp.verts[i] ); - if ( (p.edges[p.edges.size()-2].P0() == p.edges.back().P0()) || (p.edges[p.edges.size()-2].P1() == p.edges.back().P1()) ) { - p.edges.back().Flip(); //change direction - vs = make_pair( comp.verts[i].second, comp.verts[i].first ); //change direction - } - p.verts.push_back( vs ); - } - - //last segment - if ( v_start ) { - if ( p.edges.back().P1() == comp.edges[startedge].P0() ) { - Segment3 s( comp.edges[startedge].P0() , comp.edges[startedge].P1() ); - p.edges.push_back( s ); - p.verts.push_back( make_pair ( comp.verts[startedge].first , comp.verts[startedge].second ) ); - } - } - else { - Segment3 s( p.edges.back().P1() , startpoint ); - p.edges.push_back( s ); - p.verts.push_back( make_pair ( p.verts.back().second, border.verts.front().first ) ); - } - //comp - return p; -} - -/* Search for component intersected by points P0 and P1, which are respectively first and last point of a border polyline. - * @param info auxiliar information containing data about components - * @param P0 first point of border - * @param conn output parameter, true if resulting component is CC, false if it is trashC - * @return index of component - */ -int FilterZippering::searchComponent( aux_info &info, //Auxiliar info - Point3f P0, //Start border point - Point3f P1, //End border point - bool &conn ) { - int nearestC = -1; int nearestT = -1; - float distanceC = 100000*eps; float distanceT = 100000*eps; - - for ( int i = 0; i < info.nCComponent(); i ++ ) { - //for each ccon search for edges nearest to P0 and P1 - float distP0 = 200000*eps; float distP1 = 200000*eps; - for ( size_t j = 0; j < info.conn[i].edges.size(); j ++ ) { - //if ( SquaredDistance( info.conn[i].edges[j], P0 ) < distP0 ) distP0 = SquaredDistance( info.conn[i].edges[j], P0 ); - float dist_test; - vcg::Point3f clos; - vcg::SegmentPointSquaredDistance(info.conn[i].edges[j],P0,clos,dist_test); - if (dist_test( info.conn[i].edges[j], P1 ) < distP1 ) distP1 = SquaredDistance( info.conn[i].edges[j], P1 ); - vcg::SegmentPointSquaredDistance(info.conn[i].edges[j],P1,clos,dist_test); - if (dist_test( info.trash[i].edges[j], P0 ) < distP0 ) distP0 = SquaredDistance( info.trash[i].edges[j], P0 ); - if ( SquaredDistance( info.trash[i].edges[j], P1 ) < distP1 ) distP1 = SquaredDistance( info.trash[i].edges[j], P1 ); - } - if ( distP0 + distP1 < distanceT ) { distanceT = distP0 + distP1; nearestT = i; } - } - - if ( distanceC <= distanceT ) { conn = true; return nearestC; } - conn = false; return nearestT; -} - -/* Check if face f1 and f2 have a common vertex. - * @param f1 face#1 - * @param f2 face#2 - * @return true if faces share a vertex, false otherwise - */ -int FilterZippering::sharesVertex( CMeshO::FacePointer f1, CMeshO::FacePointer f2 ) { - for (int i = 0; i < 3; i ++) { - for (int k = 0; k < 3; k ++) - if ( f1->V(i) == f2->V(k) ) return i; - } - return -1; -} - -/* Given a face and an edge, check if edge projection intersect one of face edge. - * last_split edge is not checked. - * @param currentF Query face - * @param edge Query edge - * @param last_split Discard this edge during check - * @param splitted_edge Contains index of splitted edge - * @param hit Approximated intersection point - * @return true if there's intersection, false otherwise - */ -bool FilterZippering::findIntersection( CMeshO::FacePointer currentF, //face - Segment3 edge, //edge - int last_split, //previously splitted edge - int &splitted_edge, //currently splitted edge - Point3f &hit ) { //approximate intersection point - if ( currentF == NULL ) return false; - splitted_edge = -1; - Plane3 plane; plane.Init( currentF->P(0), currentF->N() ); //projection plane - Matrix44f rot_m; Point2f pt; //matrix - rot_m.SetRotateRad( Angle( currentF->N(), Point3f(0.0, 0.0, 1.0) ), currentF->N() ^ Point3f(0.0, 0.0, 1.0) ); - Segment2f s( Point2f((rot_m * plane.Projection(edge.P0())).X(), (rot_m * plane.Projection(edge.P0())).Y()), - Point2f((rot_m * plane.Projection(edge.P1())).X(), (rot_m * plane.Projection(edge.P1())).Y()) ); //projects edge on plane - for ( int e = 0; e < 3; e ++ ) { - if ( e != last_split && SegmentSegmentIntersection( s, Segment2f( Point2f( (rot_m * currentF->P(e)).X(), (rot_m * currentF->P(e)).Y() ), - Point2f( (rot_m * currentF->P1(e)).X(), (rot_m * currentF->P1(e)).Y() ) ), pt ) ) { - splitted_edge = e; break; - } - } - if (splitted_edge == -1) return false; //No intersection! - // search intersection point (approximation) - Segment3 b_edge( currentF->P(splitted_edge), currentF->P1(splitted_edge) ); - int sampleNum = SAMPLES_PER_EDGE; float step = 1.0 / (sampleNum + 1); - Point3f closest; float min_dist = b_edge.Length(); - for ( int k = 0; k <= sampleNum; k ++ ) { - Point3f currentP = edge.P0() + (edge.P1() - edge.P0())*(k*step); - if ( SquaredDistance( b_edge, currentP ) < min_dist ) { - closest = currentP; min_dist = SquaredDistance( b_edge, closest ); - } - } - if ( min_dist >= b_edge.Length() ) return false; //point not found - hit = ClosestPoint(b_edge, closest); //projection on edge - return true; -} - -/** - * Initialize queue for the selection of redundant faces. - * @param queue The queue (unsorted) - * @param m00 First mesh - * @param m01 Second mesh - * @param fullProcess If true, insert all the faces in the queue - */ -bool FilterZippering::Init_q( vector< pair >& queue, - MeshModel* a, - MeshModel* b, - bool fullProcess ) { - //full process mode: store all faces of target mesh in the in the queue - if ( fullProcess ) { - //store all the faces from B - for ( CMeshO::FaceIterator fi = b->cm.face.begin(); fi != b->cm.face.end(); ++fi ) - queue.push_back( std::make_pair(&*fi, 'B') ); - return true; - } - - //normal mode: store only border faces - vector< tri::Hole::Info > a_border; - vector< tri::Hole::Info > b_border; - //get information about border of the mesh - tri::Hole::GetInfo( a->cm, false, a_border ); - tri::Hole::GetInfo( b->cm, false, b_border ); - - if ( a_border.empty() && b_border.empty() ) { - Log( "No border face, exiting" ); - return false; - } - //face from A-border - for ( size_t i = 0; i < a_border.size(); i ++ ) { - face::Pos p = a_border[i].p; - if ( p.F()->IsD() ) continue; - do { - if ( !p.F()->IsD()) queue.push_back( make_pair(p.F(),'A') ); //label 'A' - p.NextB(); - } while ( p.F() != a_border[i].p.F() ); //visit the whole border - } - //face from B-border - for ( size_t i = 0; i < b_border.size(); i ++ ) { - face::Pos p = b_border[i].p; - if ( p.F()->IsD() ) continue; - do { - if ( !p.F()->IsD() ) queue.push_back( make_pair(p.F(),'B') ); //label 'B' - p.NextB(); - } while ( p.F() != b_border[i].p.F() ); - } - //check if queue is empty - if (queue.empty()) return false; - //queue is not empty, return true - return true; -} - -/** - * Initialize queue for the selection of redundant faces (overload for priority queue) - * @param queue The queue (unsorted) - * @param m00 First mesh - * @param m01 Second mesh - * @param fullProcess If true, insert all the faces in the queue - */ -bool FilterZippering::Init_pq( std::priority_queue< std::pair, std::vector< std::pair >, compareFaceQuality >& queue, - MeshModel* a, - MeshModel* b, - bool fullProcess ) { - //full process mode: store all faces in the queue - if ( fullProcess ) { - // all the faces from B - for ( CMeshO::FaceIterator fi = b->cm.face.begin(); fi != b->cm.face.end(); ++fi ) - queue.push( std::make_pair(&*fi, 'B') ); - return true; - } - - //normal mode: store only border faces - vector< tri::Hole::Info > a_border; - vector< tri::Hole::Info > b_border; - //get information about border of the mesh - tri::Hole::GetInfo( a->cm, false, a_border ); - tri::Hole::GetInfo( b->cm, false, b_border ); - - if ( a_border.empty() && b_border.empty() ) { - Log( "No border face, exiting" ); - return false; - } - //face from A-border - for ( size_t i = 0; i < a_border.size(); i ++ ) { - face::Pos p = a_border[i].p; - if ( p.F()->IsD() ) continue; - do { - if ( !p.F()->IsD()) queue.push( make_pair(p.F(),'A') ); //label 'A' - p.NextB(); - } while ( p.F() != a_border[i].p.F() ); //visit the whole border - } - //face from B-border - for ( size_t i = 0; i < b_border.size(); i ++ ) { - face::Pos p = b_border[i].p; - if ( p.F()->IsD() ) continue; - do { - if ( !p.F()->IsD() ) queue.push( make_pair(p.F(),'B') ); //label 'B' - p.NextB(); - } while ( p.F() != b_border[i].p.F() ); - } - //check if queue is empty - if (queue.empty()) return false; - //queue is not empty, return true - return true; -} - -int FilterZippering::preProcess (vector< std::pair >& queue, //queue - MeshModel* a, - MeshModel* b, - MeshFaceGrid grid_a, //grid on A - MeshFaceGrid grid_b, //grid on A - float max_dist ) { //max dist search - - //face count - int fc = 0; - - // update distance from border (A) - a->updateDataMask(MeshModel::MM_VERTFACETOPO); - vcg::tri::UpdateTopology::VertexFace(a->cm); - vcg::tri::UpdateFlags::FaceBorderFromVF(a->cm); - vcg::tri::UpdateQuality::VertexGeodesicFromBorder(a->cm); - a->updateDataMask(-MeshModel::MM_VERTFACETOPO); //is that correct? - vcg::tri::UpdateTopology::FaceFace(a->cm); - // update distance from border (B) - b->updateDataMask(MeshModel::MM_VERTFACETOPO); - vcg::tri::UpdateTopology::VertexFace(b->cm); - vcg::tri::UpdateFlags::FaceBorderFromVF(b->cm); - vcg::tri::UpdateQuality::VertexGeodesicFromBorder(b->cm); - b->updateDataMask(-MeshModel::MM_VERTFACETOPO); //is that correct? - vcg::tri::UpdateTopology::FaceFace(b->cm); - //perform simpleCheckRedundancy on faces of the queue - for ( int i = 0; i < queue.size(); i ++ ) { - if ( queue[i].second == 'B' ) { - if ( simpleCheckRedundancy( queue[i].first, a, grid_a, max_dist, true ) ) { - queue[i].first->SetS(); fc++; - } - } - if ( queue[i].second == 'A' ) { - if ( simpleCheckRedundancy( queue[i].first, b, grid_b, max_dist, true ) ) { - queue[i].first->SetS(); fc++; - } - } - } - return fc; -} - -int FilterZippering::preProcess_pq ( std::priority_queue< std::pair, std::vector< std::pair >, compareFaceQuality >& queue, //the queue - MeshModel* a, - MeshModel* b, - MeshFaceGrid grid_a, //grid on A - MeshFaceGrid grid_b, //grid on A - float max_dist ) { //max dist search - - //face count - int fc = 0; - - // update distance from border (A) - a->updateDataMask(MeshModel::MM_VERTFACETOPO); - vcg::tri::UpdateTopology::VertexFace(a->cm); - vcg::tri::UpdateFlags::FaceBorderFromVF(a->cm); - vcg::tri::UpdateQuality::VertexGeodesicFromBorder(a->cm); - a->updateDataMask(-MeshModel::MM_VERTFACETOPO); //is that correct? - vcg::tri::UpdateTopology::FaceFace(a->cm); - // update distance from border (B) - b->updateDataMask(MeshModel::MM_VERTFACETOPO); - vcg::tri::UpdateTopology::VertexFace(b->cm); - vcg::tri::UpdateFlags::FaceBorderFromVF(b->cm); - vcg::tri::UpdateQuality::VertexGeodesicFromBorder(b->cm); - b->updateDataMask(-MeshModel::MM_VERTFACETOPO); //is that correct? - vcg::tri::UpdateTopology::FaceFace(b->cm); - - //tmp queue - vector< std::pair > tmp_queue; - //copy data in a support vector - while (!queue.empty()) { - tmp_queue.push_back( queue.top() ); - queue.pop(); - } - //perform simpleCheckRedundancy on faces of the vector - for ( int i = 0; i < tmp_queue.size(); i ++ ) { - if ( tmp_queue[i].second == 'B' ) { - if ( simpleCheckRedundancy( tmp_queue[i].first, a, grid_a, max_dist, true ) ) { - tmp_queue[i].first->SetS(); fc++; - } - //store non-redundant faces for future check - else queue.push( tmp_queue[i] ); - } - if ( tmp_queue[i].second == 'A' ) { - if ( simpleCheckRedundancy( tmp_queue[i].first, b, grid_b, max_dist, true ) ) { - tmp_queue[i].first->SetS(); fc++; - } - else queue.push( tmp_queue[i] ); - } - } - - return fc; -} - -/** - * Select redundant faces from meshes A and B. A face is said to be redundant if - * a number of samples of the face project on the surface of the other mesh. - * @param queue Unsorted queue containing face-pointers from both meshes - * @param a, b the meshes involved in the process - * @param epsilon Maximum search distance - */ -int FilterZippering::selectRedundant( std::vector< std::pair >& queue, //queue - MeshModel* a, //mesh A - MeshModel* b, //mesh B - float epsilon ) { //max search distance - //create grid on the meshes (needed for nearest point search) - MeshFaceGrid grid_a; grid_a.Set( a->cm.face.begin(), a->cm.face.end() ); //Grid A - MeshFaceGrid grid_b; grid_b.Set( b->cm.face.begin(), b->cm.face.end() ); //Grid B - //clear selection on both meshes - tri::UpdateSelection::Clear( a->cm ); - tri::UpdateSelection::Clear( b->cm ); - //count selected faces - int sf = 0; - //fast pre processing - sf = preProcess( queue, a, b, grid_a, grid_b, epsilon ); - - //process face once at the time until queue is not empty - while ( !queue.empty() ) { - //extract face from the queue - CMeshO::FacePointer currentF = queue.back().first; char choose = queue.back().second; queue.pop_back(); - if ( currentF->IsD() || currentF->IsS() ) continue; //no op if face is deleted or selected (already tested) - //face from mesh A, test redundancy with respect to B - if (choose == 'A') { - if ( checkRedundancy( currentF, b, grid_b, epsilon ) ) { - //if face is redundant, remove it from A and put new border faces at the top of the queue - //in order to guarantee that A and B will be tested alternatively - currentF->SetS(); sf++; - //insert adjacent faces at the beginning of the queue - queue.push_back( make_pair(currentF->FFp(0),'A') ); - queue.push_back( make_pair(currentF->FFp(1),'A') ); - queue.push_back( make_pair(currentF->FFp(2),'A') ); - } - } - //face is from mesh B, test redundancy with respect to A - else { - if ( checkRedundancy( currentF, a, grid_a, epsilon ) ) { - //if face is redundant, remove it from B and put new border faces at the top of the queue - //in order to guarantee that A and B will be tested alternatively - currentF->SetS(); sf++; - //insert adjacent faces at the beginning of the queue - queue.push_back( make_pair(currentF->FFp(0),'B') ); - queue.push_back( make_pair(currentF->FFp(1),'B') ); - queue.push_back( make_pair(currentF->FFp(2),'B') ); - } - } - } - //return number of selected faces - return sf; -} - - -/** - * Select redundant faces from meshes A and B. A face is said to be redundant if - * a number of samples of the face project on the surface of the other mesh. - * @param queue priority queue containing face-pointers from both meshes ordered by quality - * @param a, b the meshes involved in the process - * @param epsilon Maximum search distance - */ -int FilterZippering::selectRedundant_pq( std::priority_queue< std::pair, std::vector< std::pair >, compareFaceQuality >& queue, //the queue - MeshModel* a, //mesh A - MeshModel* b, //mesh B - float epsilon ) { //max search distance - //create grid on the meshes (needed for nearest point search) - MeshFaceGrid grid_a; grid_a.Set( a->cm.face.begin(), a->cm.face.end() ); //Grid A - MeshFaceGrid grid_b; grid_b.Set( b->cm.face.begin(), b->cm.face.end() ); //Grid B - //clear selection on both meshes - tri::UpdateSelection::Clear( a->cm ); - tri::UpdateSelection::Clear( b->cm ); - //count selected faces - int sf = 0; - //fast pre processing - sf = preProcess_pq( queue, a, b, grid_a, grid_b, epsilon ); - - //process face once at the time until queue is not empty - while ( !queue.empty() ) { - //extract face from the queue - CMeshO::FacePointer currentF = queue.top().first; char choose = queue.top().second; queue.pop(); - if ( currentF->IsD() || currentF->IsS() ) continue; //no op if face is deleted or selected (already tested) - //face from mesh A, test redundancy with respect to B - if (choose == 'A') { - if ( checkRedundancy( currentF, b, grid_b, epsilon ) ) { - //if face is redundant, set is as Selectedand put new border faces in of the queue - currentF->SetS(); sf++; - //insert adjacent faces at the beginning of the queue - queue.push( make_pair(currentF->FFp(0),'A') ); - queue.push( make_pair(currentF->FFp(1),'A') ); - queue.push( make_pair(currentF->FFp(2),'A') ); - } - } - //face is from mesh B, test redundancy with respect to A - else { - if ( checkRedundancy( currentF, a, grid_a, epsilon ) ) { - //if face is redundant, remove it from B and put new border faces at the top of the queue - //in order to guarantee that A and B will be tested alternatively - currentF->SetS(); sf++; - //insert adjacent faces at the beginning of the queue - queue.push( make_pair(currentF->FFp(0),'B') ); - queue.push( make_pair(currentF->FFp(1),'B') ); - queue.push( make_pair(currentF->FFp(2),'B') ); - } - } - } - //return number of selected faces - return sf; -} - - -/** - * Search border-faces of the mesh having two or more border edges, then split them - * if they have only two border edges or delete them if they have three border edges. - */ -int FilterZippering::refineBorder( MeshModel* m ) { - //clear flags - tri::UpdateFlags::FaceClear(m->cm); - //recover information about the border of the mesh - vector< tri::Hole::Info > border; - tri::Hole::GetInfo( m->cm, false, border ); - //pos container - vector > border_pos; - //counter - int sf = 0; - //store pos information (index of face and of edge) for each border of the mesh - for ( size_t i = 0; i < border.size(); i ++ ) - border_pos.push_back( make_pair( tri::Index(m->cm, border[i].p.F()), border[i].p.E() ) ); - //explore the border - for ( size_t i = 0; i < border_pos.size(); i ++ ) { - //create and set a pos for the current border - face::Pos pos; - pos.Set( &(m->cm.face[border_pos[i].first]), //face pointer - border_pos[i].second, //edge index - m->cm.face[border_pos[i].first].V(border_pos[i].second) ); //vertex pointer - CMeshO::FacePointer start = pos.F(); //keep pointer to the first - //first case: start has three border face (spurious triangle). Remove it from m and step over - if ( face::BorderCount(*start) == 3 ) { - tri::Allocator::DeleteFace( m->cm, *start ); - continue; - } - //second case: start has ore or two border edges, so we visit the whole border - do { - //current face has not been deleted and has two border edges - if ( !pos.F()->IsD() && face::BorderCount(*pos.F()) >= 2 ) { - sf++; - //pointer updater - tri::Allocator::PointerUpdater vpu; - tri::Allocator::PointerUpdater fpu; - //add 4 new faces that will replace p.F() and an adjacent face - CMeshO::FaceIterator f = tri::Allocator::AddFaces( m->cm, 4, fpu ); - //update face pointers if necessary - if ( fpu.NeedUpdate() ) { fpu.Update( pos.F() ); fpu.Update( start ); } - //add 1 new vertex - CMeshO::VertexIterator v = tri::Allocator::AddVertices( m->cm, 1, vpu ); - //update vertex pointer if necessary - if ( vpu.NeedUpdate() ) - vpu.Update( pos.V() ); - //search non-border edge; - int j; for (j=0; j<3 && face::IsBorder(*pos.F(), j); j++); assert( j < 3 ); - //opposite face (vill be removed) and half-edge - CMeshO::FacePointer opp_f = pos.F()->FFp(j); int opp_j = pos.F()->FFi(j); - //place new vertex in the midpoint of the non-border edge - (*v).P() = (pos.F()->P(j) + pos.F()->P1(j))/2.0; - //create new faces - CMeshO::FacePointer f1 = &*f; f1->V(0) = pos.F()->V(j); f1->V(1) = &(*v); f1->V(2) = pos.F()->V2(j); ++f; - CMeshO::FacePointer f2 = &*f; f2->V(0) = pos.F()->V2(j); f2->V(1) = &(*v); f2->V(2) = pos.F()->V1(j); ++f; - CMeshO::FacePointer f3 = &*f; f3->V(0) = pos.F()->V(j); f3->V(1) = opp_f->V2(opp_j); f3->V(2) = &(*v); ++f; - CMeshO::FacePointer f4 = &*f; f4->V(0) = opp_f->V2(opp_j); f4->V(1) = pos.F()->V1(j); f4->V(2) = &(*v); - //attach face f1 - Not working, why? - //face::Attach( f1, 0, f3, 2 ); - //face::Attach( f1, 1, f2, 0 ); - //face::Attach( f1, 2, f1, 2 ); - //create topology for the new faces - f1->FFp(0) = f3; f1->FFp(1) = f2; f1->FFp(2) = f1; f1->FFi(0) = 2; f1->FFi(1) = 0; f1->FFi(2) = 2; - f2->FFp(0) = f1; f2->FFp(1) = f4; f2->FFp(2) = f2; f2->FFi(0) = 1; f2->FFi(1) = 1; f2->FFi(2) = 2; - f3->FFp(0) = face::IsBorder(*opp_f, (opp_j+1)%3)? f3 : opp_f->FFp((opp_j+1)%3); f3->FFp(1) = f4; f3->FFp(2) = f1; f3->FFi(0) = face::IsBorder(*opp_f, (opp_j+1)%3)? 0 : opp_f->FFi((opp_j+1)%3); f3->FFi(1) = 2; f3->FFi(2) = 0; - f4->FFp(0) = face::IsBorder(*opp_f, (opp_j+2)%3)? f4 : opp_f->FFp((opp_j+2)%3); f4->FFp(1) = f2; f4->FFp(2) = f3; f4->FFi(0) = face::IsBorder(*opp_f, (opp_j+2)%3)? 0 : opp_f->FFi((opp_j+2)%3); f4->FFi(1) = 1; f4->FFi(2) = 1; - //f3 is border-face on the 0 edge? - if( !face::IsBorder(*f3, 0) ) { - f3->FFp(0)->FFp(f3->FFi(0)) = f3; f3->FFp(0)->FFi(f3->FFi(0)) = 0; - } - //f4 is border-face on the 0 edge? - if( !face::IsBorder(*f4, 0) ) { - f4->FFp(0)->FFp(f4->FFi(0)) = f4; f4->FFp(0)->FFi(f4->FFi(0)) = 0; - } - //store pointers to face that must be deleted - CMeshO::FacePointer del_face_00 = pos.F(); - CMeshO::FacePointer del_face_01 = pos.F()->FFp(j); - //troubleshooting - //if we split start face, we have to set a new start in order to avoid infinite loop - if ( pos.F() == start ) - start = f1; - if ( (pos.F()->FFp(j)) == start && face::IsBorder(*f3,0) ) - start = f3; - if ( (pos.F()->FFp(j)) == start && face::IsBorder(*f4,0) ) - start = f4; - //it's a square, break the loop and go on - bool square = face::BorderCount(*(pos.F())) == 2 && face::BorderCount(*(pos.F()->FFp(j))) == 2; - //we deleted two faces, so we set a new pos using new faces - pos.Set( f1, 2, f1->V(0) ); - //finally remove faces - tri::Allocator::DeleteFace( m->cm, *del_face_00 ); - tri::Allocator::DeleteFace( m->cm, *del_face_01 ); - //it's a square, break the loop and go on - if (square) break; - } - pos.NextB(); - } while( pos.F() != start ); - - } - //return number of splitted faces - return sf; -} - - -/** - * Project a face from the mesh B on the surface of the mesh A. In order to project, we need - * to find nearest points of border-vertices of the face (we use the grid). - * Faces that must be retriangulated or deleted are stored in two separate vectors. - */ -void FilterZippering::projectFace( CMeshO::FacePointer f, //pointer to the face that will be projected - MeshModel* a, //mesh A - MeshFaceGrid grid_a, //grid on A - float max_dist, //max dist search - map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - vector< CMeshO::FacePointer >& tbt_faces, //vector to-be-triangulated faces - vector< CMeshO::FacePointer >& tbr_faces, //vector to-be-removed faces - vector< int >& verts ) { //indices of vertices that will define new faces - - assert( face::BorderCount( *f ) > 0 ); - int e; for ( e = 0; e < 3; e ++ ) if ( face::IsBorder(*f, e) ) break; - assert( face::IsBorder( *f, e ) ); //Check border correctness - //iteration counter - int cnt = 0; - //stack of border edges - vector < pair< int, int > > stack; - aux_info dummy; //dummy info - //store current border edge - stack.push_back( make_pair( tri::Index( a->cm, f->V(e) ), - tri::Index( a->cm, f->V1(e) ) ) ); //indices of border vertices - //while there are border edges... - while ( !stack.empty() ) { - - //avoid infinite loop - if ( cnt++ > 2*MAX_LOOP ) { - Log(GLLogStream::DEBUG, "Loop"); - stack.clear(); - continue; - } - - pair< int, int > current_edge = stack.back(); stack.pop_back(); //vertex indices - assert( current_edge.first != current_edge.second ); - tri::FaceTmark markerFunctor; markerFunctor.SetMesh(&a->cm); - face::PointDistanceBaseFunctor PDistFunct; - MeshFaceGrid::ScalarType dist = 2*max_dist; MeshFaceGrid::CoordType closestStart, closestEnd; - //search for nearest face of vertex e - CMeshO::FacePointer startF = grid_a.GetClosest(PDistFunct, markerFunctor, a->cm.vert[current_edge.first].P(), max_dist, dist, closestStart); - if ( fabs(dist) >= fabs(max_dist) ) startF = 0; - dist = 2*max_dist; - //search for nearest face of vertex e+1 - CMeshO::FacePointer endF = grid_a.GetClosest(PDistFunct, markerFunctor, a->cm.vert[current_edge.second].P(), max_dist, dist, closestEnd); - if ( fabs(dist) >= fabs(max_dist) ) endF = 0; - - //check if current edge projection fits together with an edge of the face - if ( startF != 0 && endF != 0 ) { - bool fit = false; - int vert = -1; - for ( int i = 0; i < 3; i ++ ) { - if ( vcg::Distance(closestStart, startF->P(i)) < eps ) - vert = i; - } - if ( vert != -1 ) { - //search around vert - face::JumpingPos p; - p.Set( startF, vert, startF->V(vert) ); - - do { - if ( vcg::Distance(closestStart, p.V()->P()) < eps && - vcg::Distance(closestEnd, p.F()->P1(p.E())) < eps && - face::IsBorder( *(p.F()), p.E() ) ) - fit = true; - p.FlipF(); - p.FlipE(); - }while (p.F() != startF); - //we found the face! - if (fit) - continue; - } - } - - //case 00: startF and endF are null faces: no op - if ( startF == 0 && endF == 0 ) continue; - - //case 04 (special) : vertices project on border edges. - if ( ( isOnBorder( closestStart, startF ) && isOnBorder( closestEnd, endF ) ) || - ( isOnBorder( closestStart, startF ) && endF == 0 ) || - ( isOnBorder( closestEnd, endF ) && startF == 0 ) ) { - //if closest point is really close to the border, then split border face - if ( startF != 0 && Distance(a->cm.vert[current_edge.first].P(), closestStart) <= eps ) { - map_info.insert( make_pair(startF, dummy) ); - map_info[startF].addVertex( &(a->cm.vert[current_edge.first]), current_edge.first ); - tbt_faces.push_back( startF ); - } - //check endF and eventually split it as seen before - if ( endF != 0 && Distance(a->cm.vert[current_edge.second].P(), closestEnd) <= eps ) { - map_info.insert( make_pair(endF, dummy) ); - map_info[endF].addVertex( &(a->cm.vert[current_edge.second]), current_edge.second ); - tbt_faces.push_back( endF ); - } - //check if the current_edge projects completely on border - //if true, no operation is needed - if ( handleBorderEdgeBB( current_edge, a, grid_a, max_dist, startF, endF, f, map_info, stack, verts ) ) - continue; - //false: current edge will be tested later - tbr_faces.push_back( f ); - } - - //case 05 (special): one of the border vertices doesn't project on the surface - if ( ( startF == 0 || isOnBorder( closestStart, startF ) ) || - ( endF == 0 || isOnBorder( closestEnd, endF ) ) ) { - //if closest point is really closest to the border, split border face - if ( startF != 0 && Distance(a->cm.vert[current_edge.first].P(), closestStart) <= eps ) { - map_info[startF].SetEps( eps ); - map_info[startF].Init( *startF, tri::Index( a->cm, startF->V(0) ), tri::Index( a->cm, startF->V(1) ), tri::Index( a->cm, startF->V(2) ) ); - map_info.insert( make_pair(startF, dummy) ); - map_info[startF].addVertex( &(a->cm.vert[current_edge.first]), current_edge.first ); - tbt_faces.push_back( startF ); - } - //if closest point is really closest to the border, split border face - if ( endF != 0 && Distance(a->cm.vert[current_edge.second].P(), closestEnd) <= eps ) { - map_info[endF].SetEps( eps ); - map_info[endF].Init( *endF, tri::Index( a->cm, endF->V(0) ), tri::Index( a->cm, endF->V(1) ), tri::Index( a->cm, endF->V(2) ) ); - map_info.insert( make_pair(endF, dummy) ); - map_info[endF].addVertex( &(a->cm.vert[current_edge.second]), current_edge.second ); - tbt_faces.push_back( endF ); - } - - //project one of the two vertices on the mesh - if ( startF != 0 ) a->cm.vert[current_edge.first].P() = closestStart; - if ( endF != 0 )a->cm.vert[current_edge.second].P() = closestEnd; - //we're going to split current face - tbr_faces.push_back( f ); - - //choose direction of splitting (from start to end or viceversa) according to which - //vertex doesn't project - int dir = 0; - if ( startF == 0 || isOnBorder( closestStart, startF ) ) dir = 1; - //split face - handleBorderEdgeOB( current_edge, dir, a, grid_a, max_dist, startF, endF, f, map_info, stack, tbt_faces, verts ); - - continue; - } - //end of special cases - - //initialize info on faces (if info is already initialized, no operation is done) - map_info.insert( make_pair(startF, dummy) ); map_info[startF].SetEps( eps ); map_info[startF].Init( *startF, tri::Index( a->cm, startF->V(0) ), tri::Index( a->cm, startF->V(1) ), tri::Index( a->cm, startF->V(2) ) ); - map_info.insert( make_pair(endF, dummy) ); map_info[endF].SetEps( eps ); map_info[endF].Init( *endF, tri::Index( a->cm, endF->V(0) ), tri::Index( a->cm, endF->V(1) ), tri::Index( a->cm, endF->V(2) ) ); - - //project points on the surface of the mesh - a->cm.vert[current_edge.first].P() = closestStart; - a->cm.vert[current_edge.second].P() = closestEnd; - - //case 01 : startF and endF are the same face - if ( startF == endF ) { - //start face will be retriangulated - tbt_faces.push_back( startF ); - //add information to startF - handleBorderEdgeSF( current_edge, a, startF, endF, f, map_info, stack, verts ); - continue; - } - - //case 02: startF is adjacent to endF - if ( isAdjacent( startF, endF ) ) { - //start face and end face will be retriangulated - tbt_faces.push_back( startF ); - tbt_faces.push_back( endF ); - //also, current face will be splitted and replaced by other two faces - tbr_faces.push_back( f ); - //split face and update information - handleBorderEdgeAF( current_edge, a, startF, endF, f, map_info, stack, verts ); - continue; - - } - - //case 03: startF is not adjacent to endF - if ( !isAdjacent( startF, endF ) ) { - //start face and end face will be retriangulated - tbt_faces.push_back( startF ); - tbt_faces.push_back( endF ); - //also, current face will be splitted and replaced by other two faces - tbr_faces.push_back( f ); - //split face and add edges to the stack - handleBorderEdgeNF( current_edge, a, startF, endF, f, map_info, stack, verts ); - continue; - } - - } - - -} - -//case 01: vertices of border edge project on the same face -void FilterZippering::handleBorderEdgeSF ( pair< int, int >& current_edge, //current border edge - MeshModel* a, //mesh A - CMeshO::FacePointer startF, //face where first vertex lies - CMeshO::FacePointer endF, //face where second vertex lies (=startF) - CMeshO::FacePointer splittingF, //splitting face - map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - vector < pair< int, int > >& stack, //stack containing border edges - vector< int >& verts ) { //vector of indices - - //border edge - int e; for ( e = 0; e < 3; e ++ ) if ( face::IsBorder(*splittingF, e) ) break; - - //add information to startF - if( map_info[startF].AddToBorder( Segment3 ( a->cm.vert[current_edge.first].P(), a->cm.vert[current_edge.second].P() ), - make_pair( current_edge.first, current_edge.second ) ) ) { - //avoid duplicate faces - if ( current_edge.first == tri::Index( a->cm, splittingF->V(e) ) && - current_edge.second == tri::Index( a->cm, splittingF->V1(e) ) ) return; - //store indices of vertices of the new face - verts.push_back( current_edge.first ); - verts.push_back( current_edge.second ); - verts.push_back( tri::Index( a->cm, splittingF->V2(e) ) ); - } -} - - - -//case 02: vertices of border edge project on adjacent faces -void FilterZippering::handleBorderEdgeAF ( pair< int, int >& current_edge, //current border edge - MeshModel* a, //mesh A - CMeshO::FacePointer startF, //face where first vertex lies - CMeshO::FacePointer endF, //face where second vertex lies - CMeshO::FacePointer splittingF, //splitting face - map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - vector < pair< int, int > >& stack, //stack containing border edges - vector< int >& verts ) { //vector of indices - - //calc. intersection point (approximate) and split edge - int shared; for ( int k = 0; k < 3; k ++ ) if ( startF->FFp(k) == endF ) shared = k; - Segment3 shared_edge( startF->P(shared), startF->P1(shared) ); - int sampleNum = SAMPLES_PER_EDGE; float step = 1.0/(sampleNum+1); - Point3f closest; float min_dist = shared_edge.Length(); - //subsample the current border edge and search for point which is closest to the edge shared by startF and endF - for ( int k = 0; k <= sampleNum; k ++ ) { - Point3f currentP = a->cm.vert[current_edge.first].P() + ( a->cm.vert[current_edge.second].P() - a->cm.vert[current_edge.first].P() ) * (k*step); - if ( SquaredDistance( shared_edge, currentP ) < min_dist ) { - closest = currentP; - min_dist = SquaredDistance( shared_edge, closest ); - } - } - //have we found something? - assert( SquaredDistance( shared_edge, closest ) < shared_edge.Length() ); - - //closest (pseudo-intersection) point - closest = ClosestPoint(shared_edge, closest); - - //border edge - int e; for ( e = 0; e < 3; e ++ ) if ( face::IsBorder(*splittingF, e) ) break; - - //no close vertices, add information to faces - tri::Allocator::PointerUpdater vpu; - CMeshO::VertexIterator v = tri::Allocator::AddVertices( a->cm, 1, vpu ); (*v).P() = closest; - if ( map_info[startF].AddToBorder( Segment3 ( a->cm.vert[current_edge.first].P(), (*v).P() ), - make_pair( current_edge.first, v - a->cm.vert.begin() ) ) ) { - verts.push_back( v - a->cm.vert.begin() ); - verts.push_back( tri::Index( a->cm, splittingF->V2(e) ) ); - verts.push_back( current_edge.first ); //new triangle - } - if ( map_info[endF].AddToBorder( Segment3 ( (*v).P(), a->cm.vert[current_edge.second].P() ), - make_pair( v - a->cm.vert.begin(), current_edge.second ) ) ) { - verts.push_back(v - a->cm.vert.begin()); - verts.push_back( current_edge.second ); - verts.push_back( tri::Index( a->cm, splittingF->V2(e) ) ); //new triangle - } -} - -//case 03: vertices of border edge project on non-adjacent faces -void FilterZippering::handleBorderEdgeNF ( pair< int, int >& current_edge, //current border edge - MeshModel* a, //mesh A - CMeshO::FacePointer startF, //face where first vertex lies - CMeshO::FacePointer endF, //face where second vertex lies - CMeshO::FacePointer splittingF, //splitting face - map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - vector < pair< int, int > >& stack, //stack containing border edges - vector< int >& verts ) { - //border edge - int e; for ( e = 0; e < 3; e ++ ) if ( face::IsBorder(*splittingF, e) ) break; - - //startF and endF are not adjacent, but they could share a vertex... - int w = sharesVertex( startF, endF ); - int cnt = 0; - Segment3 s( a->cm.vert[current_edge.first].P(), a->cm.vert[current_edge.second].P() ); - //if they share a vertex and border edge pass trough the vertex, we split current face using the vertex - if ( w != -1 && (cnt++ == MAX_LOOP || SquaredDistance( s, startF->P(w) ) <= eps) ) { - //too short and it's a vertex, do nothing - if ( s.Length() < eps && vcg::Distance( s.P0(), startF->P(w) ) < eps ) - return; - tri::Allocator::PointerUpdater vpu; - CMeshO::VertexIterator v = tri::Allocator::AddVertices( a->cm, 1, vpu ); (*v).P() = startF->P(w); - //add information to startF - if ( map_info[startF].AddToBorder( Segment3 ( a->cm.vert[current_edge.first].P(), (*v).P() ), - make_pair( current_edge.first, v - a->cm.vert.begin() ) ) ) { - verts.push_back( current_edge.first ); - verts.push_back( v - a->cm.vert.begin() ); - verts.push_back( tri::Index( a->cm, splittingF->V2(e) ) ); - } - //add information to endF - if ( map_info[endF].AddToBorder( Segment3 ( (*v).P(), a->cm.vert[current_edge.second].P() ), - make_pair( v - a->cm.vert.begin(), current_edge.second ) ) ) { - verts.push_back( v - a->cm.vert.begin() ); - verts.push_back( current_edge.second ); - verts.push_back( tri::Index( a->cm, splittingF->V2(e) ) ); - } - //stop - return; - } - - //add new vertex on the mid-point - tri::Allocator::PointerUpdater vpu; - CMeshO::VertexIterator v = tri::Allocator::AddVertices( a->cm, 1, vpu ); - (*v).P() = (a->cm.vert[current_edge.first].P() + a->cm.vert[current_edge.second].P())/2.00; - stack.push_back( make_pair( current_edge.first, v - a->cm.vert.begin() ) ); - stack.push_back( make_pair( v - a->cm.vert.begin(), current_edge.second ) ); -} - -//case 04: both vertices of border edge project on border edge -//return true if the whole current_edge project on border edge (no split needed) -bool FilterZippering::handleBorderEdgeBB ( std::pair< int, int >& current_edge, //current border edge - MeshModel* a, //mesh A - MeshFaceGrid grid_a, //grid on A (needed for sampling) - float max_dist, //max search dist (needed for sampling) - CMeshO::FacePointer startF, //face where first vertex lies - CMeshO::FacePointer endF, //face where second vertex lies - CMeshO::FacePointer splittingF, //splitting face - std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - std::vector < std::pair< int, int > >& stack, //stack containing border edges - std::vector< int >& verts ) { //vector of indices - - //Verify if the whole segment is on border - tri::FaceTmark markerFunctor; markerFunctor.SetMesh(&a->cm); - face::PointDistanceBaseFunctor PDistFunct; - int sampleNum = SAMPLES_PER_EDGE; float step = 1.0/(sampleNum+1); bool border = true; - Point3f closestP; float dist = 2*max_dist; - for ( int k = 0; k <= sampleNum; k ++ ) { - Point3f currentP = a->cm.vert[current_edge.first].P() + ( a->cm.vert[current_edge.second].P() - a->cm.vert[current_edge.first].P() ) * (k*step); - CMeshO::FacePointer closestFace = grid_a.GetClosest(PDistFunct, markerFunctor, currentP, 2*max_dist, dist, closestP); //closest point on mesh - if ( !isOnBorder( closestP, closestFace ) ) return false; //not completely on border; will be splitted later - } - - //border edge of splittingF - int e; for ( e = 0; e < 3; e ++ ) if ( face::IsBorder(*splittingF, e) ) break; - - //actualFace: no operation needed - if ( a->cm.vert[current_edge.first].P() == splittingF->P(e) && - a->cm.vert[current_edge.second].P() == splittingF->P1(e) ) {} - else {//create new triangle to file small hole - verts.push_back( current_edge.second ); - verts.push_back( tri::Index( a->cm, splittingF->V2(e) ) ); - verts.push_back( current_edge.first ); //new triangle - } - - return true; -} - - - -//case 05: one of the vertices doesn't project on the surface of the mesh -void FilterZippering::handleBorderEdgeOB ( std::pair< int, int >& current_edge, //current border edge - int direction, //split direction (1 from start to end, 0 from end to start) - MeshModel* a, //mesh A - MeshFaceGrid grid_a, //grid on A (needed for sampling) - float max_dist, //max search dist (needed for sampling) - CMeshO::FacePointer startF, //face where first vertex lies - CMeshO::FacePointer endF, //face where second vertex lies - CMeshO::FacePointer splittingF, //splitting face - map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - vector < pair< int, int > >& stack, //stack containing border edges - vector< CMeshO::FacePointer >& tbt_faces, //stack containing pointers to face that wille be retriangulated - vector< int >& verts ) { //vector of indices - - assert( direction == 1 || direction == 0 ); - //index of last splitted edge - int last_split = -1; - //face currently examinated - CMeshO::FacePointer currentF = direction? startF : endF; - //stop flag - bool stop = false; - - //border edge of splittingF - int e; for ( e = 0; e < 3; e ++ ) if ( face::IsBorder(*splittingF, e) ) break; - - //indices of vertices involved in the process - int startV = direction ? current_edge.first : current_edge.second; - int endV = direction ? current_edge.second : current_edge.first; - int thirdV = tri::Index( a->cm, splittingF->V2(e) ); - - int cnt = 0; - aux_info dummy; - - do { - int tosplit; - Point3f closest; - cnt++; - if (!findIntersection( currentF, Segment3(a->cm.vert[current_edge.first].P(),a->cm.vert[current_edge.second].P()), last_split, tosplit, closest )) { - stop = true; //no op - } - else { - //we found an intersection, so split the face using intersection point - last_split = currentF->FFi( tosplit ); - tri::Allocator::PointerUpdater vpu; - //new vertex on the splitting point - CMeshO::VertexIterator v = tri::Allocator::AddVertices( a->cm, 1, vpu ); (*v).P() = closest; - //init currentF aux_information - map_info.insert( make_pair(currentF, dummy) ); - map_info[currentF].Init( *currentF, tri::Index(a->cm, currentF->V(0)), tri::Index(a->cm, currentF->V(1)), tri::Index(a->cm, currentF->V(2)) ); - //add information to currentF - if ( map_info[currentF].AddToBorder( Segment3 ( (*v).P(), a->cm.vert[startV].P() ), make_pair( v - a->cm.vert.begin(), startV ) ) ) { - tbt_faces.push_back( currentF ); - verts.push_back( startV ); - verts.push_back( thirdV ); - verts.push_back( v - a->cm.vert.begin() ); - - } - //startV = v - a->cm.vert.begin(); //why? - if ( face::IsBorder( *currentF, tosplit ) ) { //border reached - verts.push_back( thirdV ); - verts.push_back( endV ); - verts.push_back( v - a->cm.vert.begin() ); - - stack.push_back( make_pair( endV, v - a->cm.vert.begin() ) ); - - stop = true; - } - else currentF = currentF->FFp(tosplit); - } - } while (!stop && cnt < MAX_LOOP); -} - - -/* Zippering of two meshes - * Given two mesh, a mesh with one or more holes (A) and a second mesh, a patch (B), fill a hole onto m surface - * using faces of patch. Algorithm const of three steps: - * - CheckRedundancy: delete face of patch that can be projected completely on A's surface - * - PatchRefinement: patch vertices are projected on A, and border edges are splitted; - * - MeshRefinement: faces where patch vertices lie are re-triangulated - * Based on Controlled and Adaptive mesh zippering - */ -bool FilterZippering::applyFilter(QAction *filter, MeshDocument &md, RichParameterSet & par, CallBackPos /**cb*/) -{ - //get the meshes - MeshModel* a = par.getMesh("FirstMesh"); - MeshModel* b = par.getMesh("SecondMesh"); - - if ( a == b ) { - errorMessage="Please add a second mesh"; - return false; - } - - //enable FF adjacency, mark, compute normals for face (both on A and B) - a->updateDataMask(MeshModel::MM_FACEFACETOPO + MeshModel::MM_FACEMARK + MeshModel::MM_FACEQUALITY + MeshModel::MM_FACECOLOR + MeshModel::MM_VERTQUALITY); - b->updateDataMask(MeshModel::MM_FACEFACETOPO + MeshModel::MM_FACEMARK + MeshModel::MM_FACEQUALITY + MeshModel::MM_FACECOLOR + MeshModel::MM_VERTQUALITY); - tri::UnMarkAll(a->cm); - tri::UnMarkAll(b->cm); - - tri::UpdateNormals::PerFaceNormalized(a->cm); - tri::UpdateFlags::FaceProjection(a->cm); - tri::UpdateNormals::PerVertexNormalized(a->cm); - tri::UpdateNormals::PerFaceNormalized(b->cm); - tri::UpdateFlags::FaceProjection(b->cm); - tri::UpdateNormals::PerVertexNormalized(b->cm); - //fixed eps - eps = 0.00001f; - - clock_t t1, t2; - t1 = clock(); - - /** - * Filter Redundancy: Select redundant faces from the surface of the meshes. - */ - if (ID(filter) == FP_REDUNDANCY) { - //queues - vector< pair > generic_queue; //unsorted queue - priority_queue< pair, vector< pair >, compareFaceQuality > priority_queue; //priority queue - - //if user chooses to usequality, initialize priority queue... - if ( par.getBool("UseQuality") ) { - if ( !Init_pq( priority_queue, a, b, par.getBool("FullProcessing") ) ) { - Log("Target mesh has no border faces - Please select Full Processing option"); - return false; - } - } - //...else initialize unsorted queue - else { - if ( !Init_q( generic_queue, a, b, par.getBool("FullProcessing") ) ) { - Log("Target mesh has no border faces - Please select Full Processing option"); - return false; - } - } - Log( "elapsed time: %f", (clock()-t1)/CLOCKS_PER_SEC ); - CMeshO::ScalarType epsilon = par.getAbsPerc("distance"); - int sf; - //if ( par.getBool("FastErosion") && par.getBool("UseQuality") ) selectRedundantFast_pq(); - //if ( par.getBool("FastErosion") && !par.getBool("UseQuality") ) selectRedundantFast(); - if ( /*!par.getBool("FastErosion") &&*/ par.getBool("UseQuality") ) { - Log( "Quality-based Redundancy" ); - sf = selectRedundant_pq( priority_queue, a, b, epsilon ); - } - if ( /*!par.getBool("FastErosion") &&*/ !par.getBool("UseQuality") ) { - Log( "Standard Redunancy" ); - sf = selectRedundant( generic_queue, a, b, epsilon ); - Log( "elapsed time: %f", (t2-t1)/CLOCKS_PER_SEC ); - } - - Log( "Selected %i redundant faces", sf ); - return true; - } //end FP_REDUNDANCY - - /** - * Filter zippering - * Assuming that redundant faces has already been removed, - * we project the border of mesh B on the surface of mesh A. - */ - if (ID(filter) == FP_ZIPPERING) { - //pre-processing step: visit the border of B and split faces with two border edges - refineBorder( b ); - //clean flags - tri::UpdateFlags::FaceClear(a->cm); - tri::UpdateFlags::FaceClear(b->cm); - //before append, clean mesh A - tri::Clean::RemoveUnreferencedVertex( a->cm ); - tri::Clean::RemoveDuplicateVertex( a->cm ); - //before append, clean mesh B - tri::Clean::RemoveUnreferencedVertex( b->cm ); - tri::Clean::RemoveDuplicateVertex( b->cm ); - //store face number - size_t fn_limit = a->cm.fn; - //append B to A and update flags - tri::UpdateTopology::FaceFace(a->cm); - tri::UpdateTopology::FaceFace(b->cm); - tri::Append::Mesh( a->cm, b->cm ); - tri::UpdateTopology::FaceFace(a->cm); - tri::UpdateFlags::FaceClear(a->cm); - tri::UpdateNormals::PerFaceNormalized(a->cm); - tri::UpdateFlags::FaceProjection(a->cm); - tri::UpdateNormals::PerVertexNormalized(a->cm); - - //create grid on mesh A - MeshFaceGrid grid; - grid.Set( a->cm.face.begin(), a->cm.face.begin() + fn_limit ); //compute grid on the original part of A - //recover information about border of B - vector< tri::Hole::Info > border; - tri::Hole::GetInfo( a->cm, false, border ); - - //Add optional attribute (true if face has been visited) - CMeshO::PerFaceAttributeHandle visited = tri::Allocator::AddPerFaceAttribute (a->cm); //check for already visited face - for ( size_t i = 0; i < a->cm.face.size(); i ++) visited[i] = false; //no face initially visited - - //create map between faces and auxiliar information - map map_info; - vector< CMeshO::FacePointer > tbt_faces; //To Be Triangulated - vector< CMeshO::FacePointer > tbr_faces; //To Be Removed - vector< int > verts; //vector containing indices of vertices of the new faces - //for each border of the mesh B, select a face and projected it on the surface of A - //add information to the faces of A when needed - for ( size_t c = 0; c < border.size(); c ++ ) { - //create pos on the surface - face::Pos p = border[c].p; p.FlipV();//CCW order - //face is from A border, go over - if ( tri::Index( a->cm, p.F() ) < fn_limit ) continue; - //face has already been deleted or visiter,go over - if ( visited[p.F()] || p.F()->IsD() ) continue; - do { - //Already deleted or visited, step over - if ( visited[p.F()] || p.F()->IsD()) { p.NextB(); continue; } - visited[p.F()] = true; - //store vertex index (avoid crash) - int v_ind = tri::Index( a->cm, p.V() ); - //project current face on the surface of A - projectFace( p.F(), a, grid, par.getFloat("distance"), map_info, tbt_faces, tbr_faces, verts ); - //restore vertex pointer - p.V() = &a->cm.vert[v_ind]; - p.NextB(); - } while ( p.F() != border[c].p.F() ); - } - - //remove splitted faces - for ( size_t i = 0; i < tbr_faces.size(); i++) { - if (!tbr_faces[i]->IsD()) tri::Allocator::DeleteFace(a->cm, *tbr_faces[i]); - } - - //remove duplicates from vector to-be-triangulated - sort( tbt_faces.begin(), tbt_faces.end() ); int k; - vector::iterator itr = unique( tbt_faces.begin(), tbt_faces.end() ); - tbt_faces.resize(itr - tbt_faces.begin() ); - //retriangulation of the faces and removal of the remaining part - vector< Point3f > coords; size_t patch_verts = verts.size(); - for ( size_t i = 0; i < tbt_faces.size(); i ++ ) { - if ( !tbt_faces[i]->IsD() ) { - handleBorder( map_info[tbt_faces[i]], tbt_faces[i]->N(), coords, verts ); - tri::Allocator::DeleteFace(a->cm, *tbt_faces[i]); - } - } - - //remove user-defined attribute - int debcnt = 0; - tri::Allocator::DeletePerFaceAttribute(a->cm, visited); - //add new faces and remove old ones - CMeshO::FaceIterator fn = tri::Allocator::AddFaces( a->cm, verts.size()/3 ); - for ( size_t k = 0; k < verts.size(); k += 3 ) { - CMeshO::VertexPointer v0 = &(a->cm.vert[verts[k]]); - CMeshO::VertexPointer v1 = &(a->cm.vert[verts[k+1]]); - CMeshO::VertexPointer v2 = &(a->cm.vert[verts[k+2]]); - if ( v0 == v1 || v1 == v2 || v0 == v2 ) - continue; - //correct orientation of face normals - if ( k < patch_verts ) { - (*fn).V(0) = v0; (*fn).V(1) = v1; (*fn).V(2) = v2; (*fn).N() = ( (*fn).P(0) - (*fn).P(2) )^( (*fn).P(1)-(*fn).P(2) ); - } - else { - (*fn).V(0) = v0; (*fn).V(2) = v1; (*fn).V(1) = v2; (*fn).N() = ( (*fn).P(0) - (*fn).P(2) )^( (*fn).P(1)-(*fn).P(2) ); - } - (*fn).SetS(); - ++fn; - } - - //remove degenerate faces - tri::Clean::RemoveDegenerateFace(a->cm); - tri::Clean::RemoveDuplicateVertex(a->cm); - //compact vectors and update topology - tri::Allocator::CompactFaceVector( a->cm ); - tri::Allocator::CompactVertexVector( a->cm ); - tri::UpdateTopology::FaceFace(a->cm); - tri::UpdateSelection::VertexFromFaceLoose(a->cm); - tri::UpdateSelection::FaceFromVertexLoose(a->cm); - tri::UpdateBounding::Box( a->cm ); - - return true; - } //end FP_ZIPPERING - - assert(0); - Log(GLLogStream::DEBUG, "**********************" ); - - return true; -} - -Q_EXPORT_PLUGIN(FilterZippering) +/**************************************************************************** +* MeshLab o o * +* A versatile mesh processing toolbox o o * +* _ O _ * +* Copyright(C) 2005 \/)\/ * +* Visual Computing Lab /\/| * +* ISTI - Italian National Research Council | * +* \ * +* All rights reserved. * +* * +* This program is free software; you can redistribute it and/or modify * +* it under the terms of the GNU General Public License as published by * +* the Free Software Foundation; either version 2 of the License, or * +* (at your option) any later version. * +* * +* This program is distributed in the hope that it will be useful, * +* but WITHOUT ANY WARRANTY; without even the implied warranty of * +* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * +* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) * +* for more details. * +* * +****************************************************************************/ +/**************************************************************************** + History + + +****************************************************************************/ + +#include +#include "filter_zippering.h" + +#include +#include +#include +#include +#include +#include +#include +#include + + +#include + +using namespace vcg; +using namespace std; + +#define MAX_LOOP 150 +//function added to map call sto distance +template +ScalarType FilterZippering::SquaredDistance( vcg::Segment3 &s, vcg::Point3 &p) +{ + Point3 clos; + ScalarType dist; + vcg::SegmentPointDistance(s, p, clos,dist ); + return dist; +} + +template +vcg::Point3 FilterZippering::ClosestPoint( vcg::Segment3 &s, vcg::Point3 &p) +{ + Point3 clos; + ScalarType dist; + vcg::SegmentPointDistance(s, p, clos,dist ); + return clos; +} + +// Constructor usually performs only two simple tasks of filling the two lists +// - typeList: with all the possible id of the filtering actions +// - actionList with the corresponding actions. If you want to add icons to your filtering actions you can do here by construction the QActions accordingly + +FilterZippering::FilterZippering() +{ + typeList << FP_REDUNDANCY << FP_ZIPPERING; + + foreach(FilterIDType tt , types()) actionList << new QAction(filterName(tt), this); +} + +// ST() must return the very short string describing each filtering action +// (this string is used also to define the menu entry) + QString FilterZippering::filterName(FilterIDType filterId) const +{ + switch(filterId) { + case FP_REDUNDANCY : return QString("Select Redundant Faces"); + case FP_ZIPPERING : return QString("Zippering"); + default : assert(0); + } +} + +// Info() must return the longer string describing each filtering action +// (this string is used in the About plugin dialog) + QString FilterZippering::filterInfo(FilterIDType filterId) const +{ + switch(filterId) { + case FP_REDUNDANCY: return QString("Remove redundant faces from one mesh or from both of them, starting from borders."); + case FP_ZIPPERING : return QString("Merge two triangle meshes into a single one. This method doesn't provide check on redundancy. Based on Controlledand Adaptive Mesh Zippering, by S.Marras, F.Ganovelli, P.Cignoni."); + default : assert(0); + } + return QString("Unknown Filter"); +} + + int FilterZippering::getRequirements(QAction *action) +{ + switch(ID(action)) + { + case FP_REDUNDANCY: return MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEMARK | MeshModel::MM_VERTCOLOR | MeshModel::MM_VERTQUALITY; + case FP_ZIPPERING : return MeshModel::MM_FACEFACETOPO | MeshModel::MM_FACEMARK | MeshModel::MM_VERTCOLOR | MeshModel::MM_VERTQUALITY; + default: assert(0); + } + return 0; +} + +// The FilterClass describes in which generic class of filters it fits. +// This choice affect the submenu in which each filter will be placed +// More than a single class can be choosen. + FilterZippering::FilterClass FilterZippering::getClass(QAction *a) +{ + switch(ID(a)) + { + case FP_REDUNDANCY: return MeshFilterInterface::Selection; + case FP_ZIPPERING : return MeshFilterInterface::Remeshing; + default : assert(0); + } + return MeshFilterInterface::Generic; +} + +// This function define the needed parameters for each filter. Return true if the filter has some parameters +// it is called every time, so you can set the default value of parameters according to the mesh +// For each parmeter you need to define, +// - the name of the parameter, +// - the string shown in the dialog +// - the default value +// - a possibly long string describing the meaning of that parameter (shown as a popup help in the dialog) +void FilterZippering::initParameterSet(QAction *action, MeshDocument &md, RichParameterSet & parlst) +{ + MeshModel *target; + float maxVal = 0.0; + switch(ID(action)) { + case FP_REDUNDANCY: + //get diagonals + foreach (target, md.meshList) { + if ( target->cm.bbox.Diag() > maxVal ) maxVal = target->cm.bbox.Diag(); + if (target != md.mm()) break; + } + parlst.addParam( new RichMesh("FirstMesh", md.mm(), &md, "Source Mesh", "The mesh with holes.") ); + parlst.addParam( new RichMesh("SecondMesh", md.mm(), &md, "Target Mesh", "The mesh that will be used as patch.") ); + parlst.addParam( new RichAbsPerc("distance", maxVal*0.01, 0, maxVal, "Max distance", "Max distance between mesh and path") ); + parlst.addParam( new RichBool("UseQuality", false, "Use quality to select redundant face", "If selected, previously computed face quality will be used in order to select redundant faces.") ); + //parlst.addParam( new RichBool("FastErosion", false, "Use fast erosion algorithm", "If selected, improves the speed-up of algorithm (results may be not accurate). Useful for large meshes.") ); + parlst.addParam( new RichBool("FullProcessing", false, "Process the whole Target Mesh", "If selected, redundancy test is performed over the whole surface of the mesh") ); + break; + + case FP_ZIPPERING : + //get diagonals + foreach (target, md.meshList) { + if ( target->cm.bbox.Diag() > maxVal ) maxVal = target->cm.bbox.Diag(); + if (target != md.mm()) break; + } + parlst.addParam( new RichMesh("FirstMesh", md.mm(), &md, "Mesh (with holes)", "The mesh with holes.") ); + parlst.addParam( new RichMesh("SecondMesh", target, &md, "Patch", "The mesh that will be used as patch.") ); + parlst.addParam( new RichAbsPerc("distance", maxVal*0.01, 0, maxVal, "Max distance", "Max distance between mesh and path") ); + break; + default: break; // do not add any parameter for the other filters + } +} + +/* Given the mesh m (let's assume it's mesh with holes) and a face on another mesh (let's assume it's + * the patch), check if the face can be completely projected on to mesh m. "Completely projected" means + * that edges' projection completely lies on one or more faces of m, and does not intersect edge + * border of mesh m. + * @param face The query face (from patch) + * @param m The mesh with holes + * @param grid A face-grid created using the faces of mesh m + * @param max_dist Max Distance allowed between m and patch; if distance between face and m is higher than max_dist, face will be discarded + * @return true if face can be completely projected on m (redundant), false otherwise. + */ +bool FilterZippering::checkRedundancy( CMeshO::FacePointer face, + MeshModel *m, + MeshFaceGrid &grid, + CMeshO::ScalarType max_dist ) +{ + // Step1: check if border edge can be projected on m + + //search for border edge + //edge adjacent to a selected face can be considered as border edge + int i; + for ( i = 0; i < 3; i ++ ) { + if ( face::IsBorder(*face, i) ) break; + if ( face->FFp(i)->IsS() ) break; + } + //no border edge find; check edge 0 + if ( i == 3 ) i = 0; + + size_t samplePerEdge = SAMPLES_PER_EDGE; + + //samples edge in uniform way + vector< Point3< CMeshO::ScalarType > > edge_samples; + Point3< CMeshO::ScalarType > edge_dir = face->P1(i) - face->P(i); edge_dir.Normalize(); + float step = 1.0/(samplePerEdge+1); //step length + for( size_t j = 0; j <= samplePerEdge; j++ ) { + edge_samples.push_back( face->P(i) + edge_dir * (j * step) ); + } + + for ( unsigned int j = 0; j < edge_samples.size(); j ++ ) { + CMeshO::FacePointer nearestF = 0; + tri::FaceTmark markerFunctor; markerFunctor.SetMesh(&m->cm); tri::UnMarkAll(m->cm); + face::PointDistanceBaseFunctor PDistFunct; + MeshFaceGrid::ScalarType dist = max_dist; MeshFaceGrid::CoordType closest; + //Search closest point on A + nearestF = grid.GetClosest(PDistFunct, markerFunctor, edge_samples[j], max_dist, dist, closest); + + if ( nearestF == 0 ) return false; //no face within given range + if ( isOnBorder(closest, nearestF ) ) return false; //closest point on border + if ( nearestF->IsD() ) return false; //face is deleted + if ( nearestF->IsS() ) return false; //face is selected (will be deleted) + } + + //check if V2(i) has a closest point on border of m + CMeshO::FacePointer nearestF = 0; + tri::FaceTmark markerFunctor; markerFunctor.SetMesh(&m->cm); + face::PointDistanceBaseFunctor PDistFunct; + MeshFaceGrid::ScalarType dist = max_dist; MeshFaceGrid::CoordType closest; + nearestF = grid.GetClosest(PDistFunct, markerFunctor, face->P2(i), max_dist, dist, closest); + if ( nearestF == 0 ) return false; //no face within given range + if ( isOnBorder( closest, nearestF ) ) return false; //closest point on border + if ( nearestF->IsD() ) return false; //face has been removed + if ( nearestF->IsS() ) return false; //face has been selected + + //check if edges are completely projectable on m + for ( int j = (i+1)%3; j != i; j = (j+1)%3 ) { + samplePerEdge = SAMPLES_PER_EDGE; + edge_samples.clear(); edge_dir = face->P1(j) - face->P(j); edge_dir.Normalize(); + for( size_t k = 0; k <= samplePerEdge; k++ ) edge_samples.push_back( face->P(j) + (face->P1(j) - face->P(j)) * (k * step) ); //campionamento + // samples on A + for ( size_t k = 0; k < edge_samples.size(); k ++ ) { + CMeshO::FacePointer nearestF = 0; + tri::FaceTmark markerFunctor; markerFunctor.SetMesh(&m->cm); tri::UnMarkAll(m->cm); + face::PointDistanceBaseFunctor PDistFunct; + MeshFaceGrid::ScalarType dist = max_dist; MeshFaceGrid::CoordType closest; + //Search closest point on A + nearestF = grid.GetClosest(PDistFunct, markerFunctor, edge_samples[k], max_dist, dist, closest); + if ( nearestF == 0 ) return false; //no face within given range + if ( isOnBorder(closest, nearestF ) ) return false; //closest point on border + if ( nearestF->IsD() ) return false; //face has been removed + if ( nearestF->IsS() ) return false; //face has been selected + } + } + // redundant + return true; +} + +/* + * Perform a simple test on a face in order to check redundancy (F.Ganovelli improvement of redundancy test) + * CheckRedundancy procedure uses a number of sampled points; instead, this function performs redundancy test + * using only the barycenter of face f. Then, after that nearest face is found, function verifies if + * vertices' distance from border is greater than max edge length of f; in this case, f is considered redundant. + * Alternatively, it can be used in order to check if face is too distance from border of M. + * @param face The query face (from patch) + * @param m The mesh with holes (Note that distance from border must be previously calculated) + * @param grid A face-grid created using the faces of mesh m + * @param max_dist Max Distance allowed between m and patch; if distance between face and m is higher than max_dist, face will be discarded + * @param test parameted used to determine the type of test + * @return true if face can be considered redundant, false otherwise +*/ +bool FilterZippering::simpleCheckRedundancy( CMeshO::FacePointer f, //face + MeshModel *m, //mesh A + MeshFaceGrid &grid, //grid A + CMeshO::ScalarType max_dist, + bool test) { //Max search distance + + Point3f qp = Barycenter(*f); //f barycenter + //search for max_edge + float max_edge = max( Distance(f->P(0),f->P(1)), max( Distance(f->P(1),f->P(2)), Distance(f->P(2),f->P(0)) ) ); + float dist = max_dist; CMeshO::FacePointer nearestF = 0; Point3f closest; + tri::FaceTmark markerFunctor; markerFunctor.SetMesh(&m->cm); UnMarkAll(m->cm); + face::PointDistanceBaseFunctor PDistFunct; + nearestF = grid.GetClosest(PDistFunct, markerFunctor, qp, max_dist, dist, closest); + if (nearestF == 0) return false; //too far away + float min_q = min( nearestF->V(0)->Q(), min( nearestF->V(1)->Q(), nearestF->V(2)->Q() ) ); //min distance of nearestF's vertices from M-border + float max_q = max( f->V(0)->Q(), max( f->V(1)->Q(), f->V(2)->Q() ) ); //max distance of F's vertices from A-border + if ( min_q <= max_edge ) return false; + if (test) if ( min_q <= max_q ) return false; + return true; +} + +/* Check if point is on border of face f. + * @param point The query-point + * @param f Face containing point + * @return true if point lies on a border edge or vertex of f, false otherwise + */ +bool FilterZippering::isOnBorder( Point3f point, CMeshO::FacePointer f ) { + // for each edge, calculates distance point-edge + if ( f == 0 ) return false; //null face + + //compute barycentric coords + float bc[3]; + InterpolationParameters( *f, f->N(), point, bc[0], bc[1], bc[2] ); + //search for max and min + int min_el = min_element(bc, bc+3) - bc; + int max_el = max_element(bc, bc+3) - bc; + //coords of max el = 1.0f -> check vertex + if ( bc[max_el] >= 1.0f - eps ) return isBorderVert( f, max_el ); + //coords of min_el = 0.0f -> check edge + if ( bc[min_el] <= 0.0f + eps ) return ( face::IsBorder( *f, (min_el+1)%3 ) || f->FFp((min_el+1)%3)->IsS() ); + + return false; +} + +/* Check if vertex i on face f belong to a border edge. + * @param f face + * @param i index of vertex + * @return true if f->V(i) lies on border-edge, false otherwise + */ +bool FilterZippering::isBorderVert( CMeshO::FacePointer f, int i ) { + + face::Pos p( f, i, f->V(i) ); + //loop + do { + //border + if ( face::IsBorder( *p.F(), p.E() ) ) return true; + //adj to selected face + if ( p.F()->FFp(p.E())->IsS() ) return true; + p.FlipE(); p.FlipF(); + } while(p.F() != f); + return false; +} + +/* Check adjiacency between two faces + * @param f1 first face + * @param f2 second face + * return true if f1 is adjacent to f2, false otherwise + */ +bool FilterZippering::isAdjacent( CMeshO::FacePointer f1, CMeshO::FacePointer f2 ) { + if ( f1 == f2 ) return false; + return ((f1 == f2->FFp(0)) || (f1 == f2->FFp(1)) || (f1 == f2->FFp(2))); +} + +/* Using auxiliar informations, replace a face with a new triangulation, discarding part of the originale face + * and tesselating the other part. + * @param info Auxiliar information from face + * @param N Face normal + * @param coords Output coords of vertices + * @param pointers Output vertex indices, will be used for creation of new faces + */ +void FilterZippering::handleBorder( aux_info &info, //Auxiliar information for triangulatio + Point3f N, //face normal (useful for proiection) + vector &/*coords*/, //output coords + vector &pointers ) { //output triangles + // rotation matrix (will be used for projection on plane) + Matrix44f rot_matrix; + rot_matrix.SetRotateRad( Angle( N, Point3f(0.0, 0.0, 1.0) ), N ^ Point3f(0.0, 0.0, 1.0) ); + //border refinement - brand new method? + //Remove intersecating border + for (size_t i = 0; i < info.border.size(); i ++) { + for ( size_t j = 0; j < info.border[i].edges.size(); j++) { + //project one border segment on face + Segment2f s( Point2f((rot_matrix * info.border[i].edges[j].P0()).X(), (rot_matrix * info.border[i].edges[j].P0()).Y()), + Point2f((rot_matrix * info.border[i].edges[j].P1()).X(), (rot_matrix * info.border[i].edges[j].P1()).Y()) ); //projects edge on plane + Point2f x; + for ( size_t k = i+1; k < info.border.size(); k++) { + for ( size_t h = 0; h < info.border[k].edges.size(); h++) { + Segment2f t( Point2f((rot_matrix * info.border[k].edges[h].P0()).X(), (rot_matrix * info.border[k].edges[h].P0()).Y()), + Point2f((rot_matrix * info.border[k].edges[h].P1()).X(), (rot_matrix * info.border[k].edges[h].P1()).Y()) ); //projects edge on plane + if ( SegmentSegmentIntersection( s, t, x ) ) { + h = info.border[k].edges.size(); + info.border.erase(info.border.begin() + k); + k--; + } + } + } + } + } + + //Split border if necessary + for (size_t i = 0; i < info.border.size(); i ++) { + for ( size_t j = 1; j < info.border[i].edges.size(); j++ ) { + //check if P0 lies on trash edge, then split border (skip first edge) + for ( size_t k = 0; k < info.trash[0].edges.size(); k++ ) { + float EPSILON=info.trash[0].edges[k].Length()/100.f; + vcg::Point3 clos; + float dist; + vcg::SegmentPointDistance( info.trash[0].edges[k],info.border[i].edges[j].P0(),clos,dist); + if (dist( info.trash[0].edges[k], info.border[i].edges[j].P0() ) == 0.0f ) { //approxim + //Split is needed + polyline newborder; newborder.edges.insert( newborder.edges.begin(), info.border[i].edges.begin() + j, info.border[i].edges.end() ); + newborder.verts.insert( newborder.verts.begin(), info.border[i].verts.begin() + j, info.border[i].verts.end() ); + info.border[i].edges.erase( info.border[i].edges.begin() + j, info.border[i].edges.end() ); + info.border[i].verts.erase( info.border[i].verts.begin() + j, info.border[i].verts.end() ); + info.border.push_back( newborder ); + break; + } + } + } + } + + // For each border... + for (size_t i = 0; i < info.border.size(); i ++) { + //search for component intersecated by border and split it into two or more components + bool conn = true; int c = searchComponent( info, info.border[i].edges.front().P0(), info.border[i].edges.back().P1(), conn ); polyline current; + if ( c == -1 ) //no intersection + continue; + if ( conn ) current = info.conn[c]; else current = info.trash[c]; + info.AddCComponent( cutComponent( current, info.border[i], rot_matrix ) ); + polyline rev_border = info.border[i]; reverse( rev_border.edges.begin(), rev_border.edges.end() ); + for ( size_t k = 0; k < rev_border.edges.size(); k ++) rev_border.edges[k].Flip(); + info.AddTComponent( cutComponent( current, rev_border, rot_matrix ) ); + if ( conn ) info.RemoveCComponent( c ); else info.RemoveTComponent( c ); + } + //No border; triangulation of the whole face + if ( info.border.size() == 0 ) { + info.conn = info.trash; + for ( size_t i = 0; i < info.conn.size(); i ++ ) { + reverse( info.conn[i].edges.begin(), info.conn[i].edges.end() ); + for ( size_t j = 0; j < info.conn[i].edges.size(); j ++ ) info.conn[i].edges[j].Flip(); + reverse( info.conn[i].verts.begin(), info.conn[i].verts.end() ); + for ( size_t j = 0; j < info.conn[i].verts.size(); j ++ ) info.conn[i].verts[j] = make_pair(info.conn[i].verts[j].second, info.conn[i].verts[j].first); + } + } + //triangulation of Ccomponent + for ( int i = 0; i < info.nCComponent(); i ++ ) { + vector< Point3f > points; //coords vector + vector< int > vertices; //vertices vector + for ( size_t j = 0; j < info.conn[i].edges.size(); j ++ ) { + points.push_back( info.conn[i].edges[j].P0() ); + vertices.push_back( info.conn[i].verts[j].first ); + } + if ( points.size() < 3 ) continue; + vector< int > indices; int iters = 0; + vector< vector< Point3f > > outlines; outlines.push_back( points ); + while ( indices.size() == 0 && ++iters < MAX_LOOP ) { + glu_tesselator::tesselate( outlines, indices ); //glu tessellator + if ( indices.size() == 0 ) + for ( size_t k = 0; k < outlines[0].size(); k ++ ) outlines[0][k] = outlines[0][k] * 10.0; //glu tessellator doesn't work properly for close points, so we scale coords in order to obtain a triangulation, if needed + } + for ( size_t k = 0; k < indices.size(); k += 3 ) { + if ( (vertices[indices[k]] != vertices[indices[k+1]]) && + (vertices[indices[k]] != vertices[indices[k+2]]) && + (vertices[indices[k+1]] != vertices[indices[k+2]] ) ) { + pointers.push_back( vertices[indices[k]] ); //save indices, in order to create new faces + pointers.push_back( vertices[indices[k+1]] ); + pointers.push_back( vertices[indices[k+2]] ); + } + } + } +} + +/* Split a component into two component, using border as separating line. Discard component external to border polyline. + * @param comp To-be-split component + * @param border Border edges + * @param rot_mat Rotation matrix, needed to project comp points on z=k plane + */ + +polyline FilterZippering::cutComponent( polyline comp, //Component to be cut + polyline border, //border + Matrix44f rot_mat ) { //Rotation matrix + + Point3f startpoint = border.edges.front().P0(); + Point3f endpoint = border.edges.back().P1(); + Point2 startpoint2D ( (rot_mat * startpoint).X(), (rot_mat * startpoint).Y() ); + Point2 endpoint2D ( (rot_mat * endpoint).X(), (rot_mat * endpoint).Y() ); + //int startedge = 0, endedge = 0; float min_dist_s = SquaredDistance( comp.edges[0], startpoint ), min_dist_e = SquaredDistance( comp.edges[0], endpoint ); + int startedge = 0, endedge = 0; + Point3f clos; + float dist; + float min_dist_s; vcg::SegmentPointDistance( comp.edges[0], startpoint, clos,min_dist_s ); + float min_dist_e; vcg::SegmentPointDistance( comp.edges[0], endpoint, clos,min_dist_e ); + bool v_start = false, v_end = false; + // search where startpoint and endpoint lie + for ( size_t i = 0; i < comp.edges.size(); i ++ ) { + if ( !v_start && SquaredDistance( comp.edges[i], startpoint ) <= min_dist_s ) { startedge = i; min_dist_s = SquaredDistance( comp.edges[i], startpoint ); } + if ( !v_end && SquaredDistance( comp.edges[i], endpoint ) <= min_dist_e ) { endedge = i; min_dist_e = SquaredDistance( comp.edges[i], endpoint ); } + if ( comp.edges[i].P1() == startpoint ) { startedge = i; v_start = true; } //lies on a vertex + if ( comp.edges[i].P0() == endpoint ) { endedge = i; v_end = true; } //lies on a vertex + } + polyline p; + // border edges will be edges of new comp + p.edges.insert( p.edges.begin(), border.edges.begin(), border.edges.end() ); + p.verts.insert( p.verts.begin(), border.verts.begin(), border.verts.end() ); + // startedge == endedge + if ( startedge == endedge && !Convex( startpoint2D, Point2 ( (rot_mat * border.edges.front().P1()).X(), (rot_mat * border.edges.front().P1()).Y() ) , endpoint2D ) ) { + Segment3 join( endpoint, startpoint ); + p.edges.push_back( join ); p.verts.push_back( make_pair( border.verts.back().second, border.verts.front().first ) ); //Vertex pointers + return p; + } + + // startedge!=endedge + // search point on the right, create oriented segment and go on + int step = -1; + Point3f c0 = border.edges.back().P0(); + vector< Segment3 >::iterator edge_it = border.edges.end(); edge_it--; + //too short segment; not reliable + while ( Distance( rot_mat * c0, rot_mat * endpoint ) <= 5.0 * eps ) { + //previous + c0 = (*edge_it).P0(); + if (edge_it != border.edges.begin()) edge_it--; else break; + } + if ( v_end ) { + if ( !Convex( Point2 ( (rot_mat * c0).X(), (rot_mat * c0).Y() ), endpoint2D, + Point2 ( (rot_mat * comp.edges[(endedge+comp.edges.size()-1)%comp.edges.size()].P0()).X(), (rot_mat * comp.edges[(endedge+comp.edges.size()-1)%comp.edges.size()].P0()).Y() ) ) ) { + step = comp.edges.size() - 1; + } + else { + step = comp.edges.size() + 1; + Segment3 s( comp.edges[endedge].P0(), comp.edges[endedge].P1() ); + p.edges.push_back( s ); step = comp.edges.size() + 1; + p.verts.push_back( make_pair( comp.verts[endedge].first, comp.verts[endedge].second ) ); + } + } + else { + if ( !Convex( Point2 ( (rot_mat * c0).X(), (rot_mat * c0).Y() ), endpoint2D, + Point2 ( (rot_mat * comp.edges[endedge].P0()).X(), (rot_mat * comp.edges[endedge].P0()).Y() ) ) ) { + Segment3 s( endpoint, comp.edges[endedge].P0() ); + p.edges.push_back( s ); step = comp.edges.size() - 1; + p.verts.push_back( make_pair(border.verts.back().second, comp.verts[endedge].first ) ); + } + else { + Segment3 s( endpoint, comp.edges[endedge].P1() ); + p.edges.push_back( s ); step = comp.edges.size() + 1; + p.verts.push_back( make_pair(border.verts.back().second, comp.verts[endedge].second ) ); + } + } + for ( int i = (endedge + step)%(comp.edges.size()); i != startedge; i = (i + step)%(comp.edges.size()) ) { + p.edges.push_back( comp.edges[i] ); + pair vs( comp.verts[i] ); + if ( (p.edges[p.edges.size()-2].P0() == p.edges.back().P0()) || (p.edges[p.edges.size()-2].P1() == p.edges.back().P1()) ) { + p.edges.back().Flip(); //change direction + vs = make_pair( comp.verts[i].second, comp.verts[i].first ); //change direction + } + p.verts.push_back( vs ); + } + + //last segment + if ( v_start ) { + if ( p.edges.back().P1() == comp.edges[startedge].P0() ) { + Segment3 s( comp.edges[startedge].P0() , comp.edges[startedge].P1() ); + p.edges.push_back( s ); + p.verts.push_back( make_pair ( comp.verts[startedge].first , comp.verts[startedge].second ) ); + } + } + else { + Segment3 s( p.edges.back().P1() , startpoint ); + p.edges.push_back( s ); + p.verts.push_back( make_pair ( p.verts.back().second, border.verts.front().first ) ); + } + //comp + return p; +} + +/* Search for component intersected by points P0 and P1, which are respectively first and last point of a border polyline. + * @param info auxiliar information containing data about components + * @param P0 first point of border + * @param conn output parameter, true if resulting component is CC, false if it is trashC + * @return index of component + */ +int FilterZippering::searchComponent( aux_info &info, //Auxiliar info + Point3f P0, //Start border point + Point3f P1, //End border point + bool &conn ) { + int nearestC = -1; int nearestT = -1; + float distanceC = 100000*eps; float distanceT = 100000*eps; + + for ( int i = 0; i < info.nCComponent(); i ++ ) { + //for each ccon search for edges nearest to P0 and P1 + float distP0 = 200000*eps; float distP1 = 200000*eps; + for ( size_t j = 0; j < info.conn[i].edges.size(); j ++ ) { + //if ( SquaredDistance( info.conn[i].edges[j], P0 ) < distP0 ) distP0 = SquaredDistance( info.conn[i].edges[j], P0 ); + float dist_test; + vcg::Point3f clos; + vcg::SegmentPointSquaredDistance(info.conn[i].edges[j],P0,clos,dist_test); + if (dist_test( info.conn[i].edges[j], P1 ) < distP1 ) distP1 = SquaredDistance( info.conn[i].edges[j], P1 ); + vcg::SegmentPointSquaredDistance(info.conn[i].edges[j],P1,clos,dist_test); + if (dist_test( info.trash[i].edges[j], P0 ) < distP0 ) distP0 = SquaredDistance( info.trash[i].edges[j], P0 ); + if ( SquaredDistance( info.trash[i].edges[j], P1 ) < distP1 ) distP1 = SquaredDistance( info.trash[i].edges[j], P1 ); + } + if ( distP0 + distP1 < distanceT ) { distanceT = distP0 + distP1; nearestT = i; } + } + + if ( distanceC <= distanceT ) { conn = true; return nearestC; } + conn = false; return nearestT; +} + +/* Check if face f1 and f2 have a common vertex. + * @param f1 face#1 + * @param f2 face#2 + * @return true if faces share a vertex, false otherwise + */ +int FilterZippering::sharesVertex( CMeshO::FacePointer f1, CMeshO::FacePointer f2 ) { + for (int i = 0; i < 3; i ++) { + for (int k = 0; k < 3; k ++) + if ( f1->V(i) == f2->V(k) ) return i; + } + return -1; +} + +/* Given a face and an edge, check if edge projection intersect one of face edge. + * last_split edge is not checked. + * @param currentF Query face + * @param edge Query edge + * @param last_split Discard this edge during check + * @param splitted_edge Contains index of splitted edge + * @param hit Approximated intersection point + * @return true if there's intersection, false otherwise + */ +bool FilterZippering::findIntersection( CMeshO::FacePointer currentF, //face + Segment3 edge, //edge + int last_split, //previously splitted edge + int &splitted_edge, //currently splitted edge + Point3f &hit ) { //approximate intersection point + if ( currentF == NULL ) return false; + splitted_edge = -1; + Plane3 plane; plane.Init( currentF->P(0), currentF->N() ); //projection plane + Matrix44f rot_m; Point2f pt; //matrix + rot_m.SetRotateRad( Angle( currentF->N(), Point3f(0.0, 0.0, 1.0) ), currentF->N() ^ Point3f(0.0, 0.0, 1.0) ); + Segment2f s( Point2f((rot_m * plane.Projection(edge.P0())).X(), (rot_m * plane.Projection(edge.P0())).Y()), + Point2f((rot_m * plane.Projection(edge.P1())).X(), (rot_m * plane.Projection(edge.P1())).Y()) ); //projects edge on plane + for ( int e = 0; e < 3; e ++ ) { + if ( e != last_split && SegmentSegmentIntersection( s, Segment2f( Point2f( (rot_m * currentF->P(e)).X(), (rot_m * currentF->P(e)).Y() ), + Point2f( (rot_m * currentF->P1(e)).X(), (rot_m * currentF->P1(e)).Y() ) ), pt ) ) { + splitted_edge = e; break; + } + } + if (splitted_edge == -1) return false; //No intersection! + // search intersection point (approximation) + Segment3 b_edge( currentF->P(splitted_edge), currentF->P1(splitted_edge) ); + int sampleNum = SAMPLES_PER_EDGE; float step = 1.0 / (sampleNum + 1); + Point3f closest; float min_dist = b_edge.Length(); + for ( int k = 0; k <= sampleNum; k ++ ) { + Point3f currentP = edge.P0() + (edge.P1() - edge.P0())*(k*step); + if ( SquaredDistance( b_edge, currentP ) < min_dist ) { + closest = currentP; min_dist = SquaredDistance( b_edge, closest ); + } + } + if ( min_dist >= b_edge.Length() ) return false; //point not found + hit = ClosestPoint(b_edge, closest); //projection on edge + return true; +} + +/** + * Initialize queue for the selection of redundant faces. + * @param queue The queue (unsorted) + * @param m00 First mesh + * @param m01 Second mesh + * @param fullProcess If true, insert all the faces in the queue + */ +bool FilterZippering::Init_q( vector< pair >& queue, + MeshModel* a, + MeshModel* b, + bool fullProcess ) { + //full process mode: store all faces of target mesh in the in the queue + if ( fullProcess ) { + //store all the faces from B + for ( CMeshO::FaceIterator fi = b->cm.face.begin(); fi != b->cm.face.end(); ++fi ) + queue.push_back( std::make_pair(&*fi, 'B') ); + return true; + } + + //normal mode: store only border faces + vector< tri::Hole::Info > a_border; + vector< tri::Hole::Info > b_border; + //get information about border of the mesh + tri::Hole::GetInfo( a->cm, false, a_border ); + tri::Hole::GetInfo( b->cm, false, b_border ); + + if ( a_border.empty() && b_border.empty() ) { + Log( "No border face, exiting" ); + return false; + } + //face from A-border + for ( size_t i = 0; i < a_border.size(); i ++ ) { + face::Pos p = a_border[i].p; + if ( p.F()->IsD() ) continue; + do { + if ( !p.F()->IsD()) queue.push_back( make_pair(p.F(),'A') ); //label 'A' + p.NextB(); + } while ( p.F() != a_border[i].p.F() ); //visit the whole border + } + //face from B-border + for ( size_t i = 0; i < b_border.size(); i ++ ) { + face::Pos p = b_border[i].p; + if ( p.F()->IsD() ) continue; + do { + if ( !p.F()->IsD() ) queue.push_back( make_pair(p.F(),'B') ); //label 'B' + p.NextB(); + } while ( p.F() != b_border[i].p.F() ); + } + //check if queue is empty + if (queue.empty()) return false; + //queue is not empty, return true + return true; +} + +/** + * Initialize queue for the selection of redundant faces (overload for priority queue) + * @param queue The queue (unsorted) + * @param m00 First mesh + * @param m01 Second mesh + * @param fullProcess If true, insert all the faces in the queue + */ +bool FilterZippering::Init_pq( std::priority_queue< std::pair, std::vector< std::pair >, compareFaceQuality >& queue, + MeshModel* a, + MeshModel* b, + bool fullProcess ) { + //full process mode: store all faces in the queue + if ( fullProcess ) { + // all the faces from B + for ( CMeshO::FaceIterator fi = b->cm.face.begin(); fi != b->cm.face.end(); ++fi ) + queue.push( std::make_pair(&*fi, 'B') ); + return true; + } + + //normal mode: store only border faces + vector< tri::Hole::Info > a_border; + vector< tri::Hole::Info > b_border; + //get information about border of the mesh + tri::Hole::GetInfo( a->cm, false, a_border ); + tri::Hole::GetInfo( b->cm, false, b_border ); + + if ( a_border.empty() && b_border.empty() ) { + Log( "No border face, exiting" ); + return false; + } + //face from A-border + for ( size_t i = 0; i < a_border.size(); i ++ ) { + face::Pos p = a_border[i].p; + if ( p.F()->IsD() ) continue; + do { + if ( !p.F()->IsD()) queue.push( make_pair(p.F(),'A') ); //label 'A' + p.NextB(); + } while ( p.F() != a_border[i].p.F() ); //visit the whole border + } + //face from B-border + for ( size_t i = 0; i < b_border.size(); i ++ ) { + face::Pos p = b_border[i].p; + if ( p.F()->IsD() ) continue; + do { + if ( !p.F()->IsD() ) queue.push( make_pair(p.F(),'B') ); //label 'B' + p.NextB(); + } while ( p.F() != b_border[i].p.F() ); + } + //check if queue is empty + if (queue.empty()) return false; + //queue is not empty, return true + return true; +} + +int FilterZippering::preProcess (vector< std::pair >& queue, //queue + MeshModel* a, + MeshModel* b, + MeshFaceGrid grid_a, //grid on A + MeshFaceGrid grid_b, //grid on A + float max_dist ) { //max dist search + + //face count + int fc = 0; + + // update distance from border (A) + a->updateDataMask(MeshModel::MM_VERTFACETOPO); + vcg::tri::UpdateTopology::VertexFace(a->cm); + vcg::tri::UpdateFlags::FaceBorderFromVF(a->cm); + vcg::tri::UpdateQuality::VertexGeodesicFromBorder(a->cm); + a->updateDataMask(-MeshModel::MM_VERTFACETOPO); //is that correct? + vcg::tri::UpdateTopology::FaceFace(a->cm); + // update distance from border (B) + b->updateDataMask(MeshModel::MM_VERTFACETOPO); + vcg::tri::UpdateTopology::VertexFace(b->cm); + vcg::tri::UpdateFlags::FaceBorderFromVF(b->cm); + vcg::tri::UpdateQuality::VertexGeodesicFromBorder(b->cm); + b->updateDataMask(-MeshModel::MM_VERTFACETOPO); //is that correct? + vcg::tri::UpdateTopology::FaceFace(b->cm); + //perform simpleCheckRedundancy on faces of the queue + for ( int i = 0; i < queue.size(); i ++ ) { + if ( queue[i].second == 'B' ) { + if ( simpleCheckRedundancy( queue[i].first, a, grid_a, max_dist, true ) ) { + queue[i].first->SetS(); fc++; + } + } + if ( queue[i].second == 'A' ) { + if ( simpleCheckRedundancy( queue[i].first, b, grid_b, max_dist, true ) ) { + queue[i].first->SetS(); fc++; + } + } + } + return fc; +} + +int FilterZippering::preProcess_pq ( std::priority_queue< std::pair, std::vector< std::pair >, compareFaceQuality >& queue, //the queue + MeshModel* a, + MeshModel* b, + MeshFaceGrid grid_a, //grid on A + MeshFaceGrid grid_b, //grid on A + float max_dist ) { //max dist search + + //face count + int fc = 0; + + // update distance from border (A) + a->updateDataMask(MeshModel::MM_VERTFACETOPO); + vcg::tri::UpdateTopology::VertexFace(a->cm); + vcg::tri::UpdateFlags::FaceBorderFromVF(a->cm); + vcg::tri::UpdateQuality::VertexGeodesicFromBorder(a->cm); + a->updateDataMask(-MeshModel::MM_VERTFACETOPO); //is that correct? + vcg::tri::UpdateTopology::FaceFace(a->cm); + // update distance from border (B) + b->updateDataMask(MeshModel::MM_VERTFACETOPO); + vcg::tri::UpdateTopology::VertexFace(b->cm); + vcg::tri::UpdateFlags::FaceBorderFromVF(b->cm); + vcg::tri::UpdateQuality::VertexGeodesicFromBorder(b->cm); + b->updateDataMask(-MeshModel::MM_VERTFACETOPO); //is that correct? + vcg::tri::UpdateTopology::FaceFace(b->cm); + + //tmp queue + vector< std::pair > tmp_queue; + //copy data in a support vector + while (!queue.empty()) { + tmp_queue.push_back( queue.top() ); + queue.pop(); + } + //perform simpleCheckRedundancy on faces of the vector + for ( int i = 0; i < tmp_queue.size(); i ++ ) { + if ( tmp_queue[i].second == 'B' ) { + if ( simpleCheckRedundancy( tmp_queue[i].first, a, grid_a, max_dist, true ) ) { + tmp_queue[i].first->SetS(); fc++; + } + //store non-redundant faces for future check + else queue.push( tmp_queue[i] ); + } + if ( tmp_queue[i].second == 'A' ) { + if ( simpleCheckRedundancy( tmp_queue[i].first, b, grid_b, max_dist, true ) ) { + tmp_queue[i].first->SetS(); fc++; + } + else queue.push( tmp_queue[i] ); + } + } + + return fc; +} + +/** + * Select redundant faces from meshes A and B. A face is said to be redundant if + * a number of samples of the face project on the surface of the other mesh. + * @param queue Unsorted queue containing face-pointers from both meshes + * @param a, b the meshes involved in the process + * @param epsilon Maximum search distance + */ +int FilterZippering::selectRedundant( std::vector< std::pair >& queue, //queue + MeshModel* a, //mesh A + MeshModel* b, //mesh B + float epsilon ) { //max search distance + //create grid on the meshes (needed for nearest point search) + MeshFaceGrid grid_a; grid_a.Set( a->cm.face.begin(), a->cm.face.end() ); //Grid A + MeshFaceGrid grid_b; grid_b.Set( b->cm.face.begin(), b->cm.face.end() ); //Grid B + //clear selection on both meshes + tri::UpdateSelection::Clear( a->cm ); + tri::UpdateSelection::Clear( b->cm ); + //count selected faces + int sf = 0; + //fast pre processing + sf = preProcess( queue, a, b, grid_a, grid_b, epsilon ); + + //process face once at the time until queue is not empty + while ( !queue.empty() ) { + //extract face from the queue + CMeshO::FacePointer currentF = queue.back().first; char choose = queue.back().second; queue.pop_back(); + if ( currentF->IsD() || currentF->IsS() ) continue; //no op if face is deleted or selected (already tested) + //face from mesh A, test redundancy with respect to B + if (choose == 'A') { + if ( checkRedundancy( currentF, b, grid_b, epsilon ) ) { + //if face is redundant, remove it from A and put new border faces at the top of the queue + //in order to guarantee that A and B will be tested alternatively + currentF->SetS(); sf++; + //insert adjacent faces at the beginning of the queue + queue.push_back( make_pair(currentF->FFp(0),'A') ); + queue.push_back( make_pair(currentF->FFp(1),'A') ); + queue.push_back( make_pair(currentF->FFp(2),'A') ); + } + } + //face is from mesh B, test redundancy with respect to A + else { + if ( checkRedundancy( currentF, a, grid_a, epsilon ) ) { + //if face is redundant, remove it from B and put new border faces at the top of the queue + //in order to guarantee that A and B will be tested alternatively + currentF->SetS(); sf++; + //insert adjacent faces at the beginning of the queue + queue.push_back( make_pair(currentF->FFp(0),'B') ); + queue.push_back( make_pair(currentF->FFp(1),'B') ); + queue.push_back( make_pair(currentF->FFp(2),'B') ); + } + } + } + //return number of selected faces + return sf; +} + + +/** + * Select redundant faces from meshes A and B. A face is said to be redundant if + * a number of samples of the face project on the surface of the other mesh. + * @param queue priority queue containing face-pointers from both meshes ordered by quality + * @param a, b the meshes involved in the process + * @param epsilon Maximum search distance + */ +int FilterZippering::selectRedundant_pq( std::priority_queue< std::pair, std::vector< std::pair >, compareFaceQuality >& queue, //the queue + MeshModel* a, //mesh A + MeshModel* b, //mesh B + float epsilon ) { //max search distance + //create grid on the meshes (needed for nearest point search) + MeshFaceGrid grid_a; grid_a.Set( a->cm.face.begin(), a->cm.face.end() ); //Grid A + MeshFaceGrid grid_b; grid_b.Set( b->cm.face.begin(), b->cm.face.end() ); //Grid B + //clear selection on both meshes + tri::UpdateSelection::Clear( a->cm ); + tri::UpdateSelection::Clear( b->cm ); + //count selected faces + int sf = 0; + //fast pre processing + sf = preProcess_pq( queue, a, b, grid_a, grid_b, epsilon ); + + //process face once at the time until queue is not empty + while ( !queue.empty() ) { + //extract face from the queue + CMeshO::FacePointer currentF = queue.top().first; char choose = queue.top().second; queue.pop(); + if ( currentF->IsD() || currentF->IsS() ) continue; //no op if face is deleted or selected (already tested) + //face from mesh A, test redundancy with respect to B + if (choose == 'A') { + if ( checkRedundancy( currentF, b, grid_b, epsilon ) ) { + //if face is redundant, set is as Selectedand put new border faces in of the queue + currentF->SetS(); sf++; + //insert adjacent faces at the beginning of the queue + queue.push( make_pair(currentF->FFp(0),'A') ); + queue.push( make_pair(currentF->FFp(1),'A') ); + queue.push( make_pair(currentF->FFp(2),'A') ); + } + } + //face is from mesh B, test redundancy with respect to A + else { + if ( checkRedundancy( currentF, a, grid_a, epsilon ) ) { + //if face is redundant, remove it from B and put new border faces at the top of the queue + //in order to guarantee that A and B will be tested alternatively + currentF->SetS(); sf++; + //insert adjacent faces at the beginning of the queue + queue.push( make_pair(currentF->FFp(0),'B') ); + queue.push( make_pair(currentF->FFp(1),'B') ); + queue.push( make_pair(currentF->FFp(2),'B') ); + } + } + } + //return number of selected faces + return sf; +} + + +/** + * Search border-faces of the mesh having two or more border edges, then split them + * if they have only two border edges or delete them if they have three border edges. + */ +int FilterZippering::refineBorder( MeshModel* m ) { + //clear flags + tri::UpdateFlags::FaceClear(m->cm); + //recover information about the border of the mesh + vector< tri::Hole::Info > border; + tri::Hole::GetInfo( m->cm, false, border ); + //pos container + vector > border_pos; + //counter + int sf = 0; + //store pos information (index of face and of edge) for each border of the mesh + for ( size_t i = 0; i < border.size(); i ++ ) + border_pos.push_back( make_pair( tri::Index(m->cm, border[i].p.F()), border[i].p.E() ) ); + //explore the border + for ( size_t i = 0; i < border_pos.size(); i ++ ) { + //create and set a pos for the current border + face::Pos pos; + pos.Set( &(m->cm.face[border_pos[i].first]), //face pointer + border_pos[i].second, //edge index + m->cm.face[border_pos[i].first].V(border_pos[i].second) ); //vertex pointer + CMeshO::FacePointer start = pos.F(); //keep pointer to the first + //first case: start has three border face (spurious triangle). Remove it from m and step over + if ( face::BorderCount(*start) == 3 ) { + tri::Allocator::DeleteFace( m->cm, *start ); + continue; + } + //second case: start has ore or two border edges, so we visit the whole border + do { + //current face has not been deleted and has two border edges + if ( !pos.F()->IsD() && face::BorderCount(*pos.F()) >= 2 ) { + sf++; + //pointer updater + tri::Allocator::PointerUpdater vpu; + tri::Allocator::PointerUpdater fpu; + //add 4 new faces that will replace p.F() and an adjacent face + CMeshO::FaceIterator f = tri::Allocator::AddFaces( m->cm, 4, fpu ); + //update face pointers if necessary + if ( fpu.NeedUpdate() ) { fpu.Update( pos.F() ); fpu.Update( start ); } + //add 1 new vertex + CMeshO::VertexIterator v = tri::Allocator::AddVertices( m->cm, 1, vpu ); + //update vertex pointer if necessary + if ( vpu.NeedUpdate() ) + vpu.Update( pos.V() ); + //search non-border edge; + int j; for (j=0; j<3 && face::IsBorder(*pos.F(), j); j++); assert( j < 3 ); + //opposite face (vill be removed) and half-edge + CMeshO::FacePointer opp_f = pos.F()->FFp(j); int opp_j = pos.F()->FFi(j); + //place new vertex in the midpoint of the non-border edge + (*v).P() = (pos.F()->P(j) + pos.F()->P1(j))/2.0; + //create new faces + CMeshO::FacePointer f1 = &*f; f1->V(0) = pos.F()->V(j); f1->V(1) = &(*v); f1->V(2) = pos.F()->V2(j); ++f; + CMeshO::FacePointer f2 = &*f; f2->V(0) = pos.F()->V2(j); f2->V(1) = &(*v); f2->V(2) = pos.F()->V1(j); ++f; + CMeshO::FacePointer f3 = &*f; f3->V(0) = pos.F()->V(j); f3->V(1) = opp_f->V2(opp_j); f3->V(2) = &(*v); ++f; + CMeshO::FacePointer f4 = &*f; f4->V(0) = opp_f->V2(opp_j); f4->V(1) = pos.F()->V1(j); f4->V(2) = &(*v); + //attach face f1 - Not working, why? + //face::Attach( f1, 0, f3, 2 ); + //face::Attach( f1, 1, f2, 0 ); + //face::Attach( f1, 2, f1, 2 ); + //create topology for the new faces + f1->FFp(0) = f3; f1->FFp(1) = f2; f1->FFp(2) = f1; f1->FFi(0) = 2; f1->FFi(1) = 0; f1->FFi(2) = 2; + f2->FFp(0) = f1; f2->FFp(1) = f4; f2->FFp(2) = f2; f2->FFi(0) = 1; f2->FFi(1) = 1; f2->FFi(2) = 2; + f3->FFp(0) = face::IsBorder(*opp_f, (opp_j+1)%3)? f3 : opp_f->FFp((opp_j+1)%3); f3->FFp(1) = f4; f3->FFp(2) = f1; f3->FFi(0) = face::IsBorder(*opp_f, (opp_j+1)%3)? 0 : opp_f->FFi((opp_j+1)%3); f3->FFi(1) = 2; f3->FFi(2) = 0; + f4->FFp(0) = face::IsBorder(*opp_f, (opp_j+2)%3)? f4 : opp_f->FFp((opp_j+2)%3); f4->FFp(1) = f2; f4->FFp(2) = f3; f4->FFi(0) = face::IsBorder(*opp_f, (opp_j+2)%3)? 0 : opp_f->FFi((opp_j+2)%3); f4->FFi(1) = 1; f4->FFi(2) = 1; + //f3 is border-face on the 0 edge? + if( !face::IsBorder(*f3, 0) ) { + f3->FFp(0)->FFp(f3->FFi(0)) = f3; f3->FFp(0)->FFi(f3->FFi(0)) = 0; + } + //f4 is border-face on the 0 edge? + if( !face::IsBorder(*f4, 0) ) { + f4->FFp(0)->FFp(f4->FFi(0)) = f4; f4->FFp(0)->FFi(f4->FFi(0)) = 0; + } + //store pointers to face that must be deleted + CMeshO::FacePointer del_face_00 = pos.F(); + CMeshO::FacePointer del_face_01 = pos.F()->FFp(j); + //troubleshooting + //if we split start face, we have to set a new start in order to avoid infinite loop + if ( pos.F() == start ) + start = f1; + if ( (pos.F()->FFp(j)) == start && face::IsBorder(*f3,0) ) + start = f3; + if ( (pos.F()->FFp(j)) == start && face::IsBorder(*f4,0) ) + start = f4; + //it's a square, break the loop and go on + bool square = face::BorderCount(*(pos.F())) == 2 && face::BorderCount(*(pos.F()->FFp(j))) == 2; + //we deleted two faces, so we set a new pos using new faces + pos.Set( f1, 2, f1->V(0) ); + //finally remove faces + tri::Allocator::DeleteFace( m->cm, *del_face_00 ); + tri::Allocator::DeleteFace( m->cm, *del_face_01 ); + //it's a square, break the loop and go on + if (square) break; + } + pos.NextB(); + } while( pos.F() != start ); + + } + //return number of splitted faces + return sf; +} + + +/** + * Project a face from the mesh B on the surface of the mesh A. In order to project, we need + * to find nearest points of border-vertices of the face (we use the grid). + * Faces that must be retriangulated or deleted are stored in two separate vectors. + */ +void FilterZippering::projectFace( CMeshO::FacePointer f, //pointer to the face that will be projected + MeshModel* a, //mesh A + MeshFaceGrid grid_a, //grid on A + float max_dist, //max dist search + map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + vector< CMeshO::FacePointer >& tbt_faces, //vector to-be-triangulated faces + vector< CMeshO::FacePointer >& tbr_faces, //vector to-be-removed faces + vector< int >& verts ) { //indices of vertices that will define new faces + + assert( face::BorderCount( *f ) > 0 ); + int e; for ( e = 0; e < 3; e ++ ) if ( face::IsBorder(*f, e) ) break; + assert( face::IsBorder( *f, e ) ); //Check border correctness + //iteration counter + int cnt = 0; + //stack of border edges + vector < pair< int, int > > stack; + aux_info dummy; //dummy info + //store current border edge + stack.push_back( make_pair( tri::Index( a->cm, f->V(e) ), + tri::Index( a->cm, f->V1(e) ) ) ); //indices of border vertices + //while there are border edges... + while ( !stack.empty() ) { + + //avoid infinite loop + if ( cnt++ > 2*MAX_LOOP ) { + Log(GLLogStream::DEBUG, "Loop"); + stack.clear(); + continue; + } + + pair< int, int > current_edge = stack.back(); stack.pop_back(); //vertex indices + assert( current_edge.first != current_edge.second ); + tri::FaceTmark markerFunctor; markerFunctor.SetMesh(&a->cm); + face::PointDistanceBaseFunctor PDistFunct; + MeshFaceGrid::ScalarType dist = 2*max_dist; MeshFaceGrid::CoordType closestStart, closestEnd; + //search for nearest face of vertex e + CMeshO::FacePointer startF = grid_a.GetClosest(PDistFunct, markerFunctor, a->cm.vert[current_edge.first].P(), max_dist, dist, closestStart); + if ( fabs(dist) >= fabs(max_dist) ) startF = 0; + dist = 2*max_dist; + //search for nearest face of vertex e+1 + CMeshO::FacePointer endF = grid_a.GetClosest(PDistFunct, markerFunctor, a->cm.vert[current_edge.second].P(), max_dist, dist, closestEnd); + if ( fabs(dist) >= fabs(max_dist) ) endF = 0; + + //check if current edge projection fits together with an edge of the face + if ( startF != 0 && endF != 0 ) { + bool fit = false; + int vert = -1; + for ( int i = 0; i < 3; i ++ ) { + if ( vcg::Distance(closestStart, startF->P(i)) < eps ) + vert = i; + } + if ( vert != -1 ) { + //search around vert + face::JumpingPos p; + p.Set( startF, vert, startF->V(vert) ); + + do { + if ( vcg::Distance(closestStart, p.V()->P()) < eps && + vcg::Distance(closestEnd, p.F()->P1(p.E())) < eps && + face::IsBorder( *(p.F()), p.E() ) ) + fit = true; + p.FlipF(); + p.FlipE(); + }while (p.F() != startF); + //we found the face! + if (fit) + continue; + } + } + + //case 00: startF and endF are null faces: no op + if ( startF == 0 && endF == 0 ) continue; + + //case 04 (special) : vertices project on border edges. + if ( ( isOnBorder( closestStart, startF ) && isOnBorder( closestEnd, endF ) ) || + ( isOnBorder( closestStart, startF ) && endF == 0 ) || + ( isOnBorder( closestEnd, endF ) && startF == 0 ) ) { + //if closest point is really close to the border, then split border face + if ( startF != 0 && Distance(a->cm.vert[current_edge.first].P(), closestStart) <= eps ) { + map_info.insert( make_pair(startF, dummy) ); + map_info[startF].addVertex( &(a->cm.vert[current_edge.first]), current_edge.first ); + tbt_faces.push_back( startF ); + } + //check endF and eventually split it as seen before + if ( endF != 0 && Distance(a->cm.vert[current_edge.second].P(), closestEnd) <= eps ) { + map_info.insert( make_pair(endF, dummy) ); + map_info[endF].addVertex( &(a->cm.vert[current_edge.second]), current_edge.second ); + tbt_faces.push_back( endF ); + } + //check if the current_edge projects completely on border + //if true, no operation is needed + if ( handleBorderEdgeBB( current_edge, a, grid_a, max_dist, startF, endF, f, map_info, stack, verts ) ) + continue; + //false: current edge will be tested later + tbr_faces.push_back( f ); + } + + //case 05 (special): one of the border vertices doesn't project on the surface + if ( ( startF == 0 || isOnBorder( closestStart, startF ) ) || + ( endF == 0 || isOnBorder( closestEnd, endF ) ) ) { + //if closest point is really closest to the border, split border face + if ( startF != 0 && Distance(a->cm.vert[current_edge.first].P(), closestStart) <= eps ) { + map_info[startF].SetEps( eps ); + map_info[startF].Init( *startF, tri::Index( a->cm, startF->V(0) ), tri::Index( a->cm, startF->V(1) ), tri::Index( a->cm, startF->V(2) ) ); + map_info.insert( make_pair(startF, dummy) ); + map_info[startF].addVertex( &(a->cm.vert[current_edge.first]), current_edge.first ); + tbt_faces.push_back( startF ); + } + //if closest point is really closest to the border, split border face + if ( endF != 0 && Distance(a->cm.vert[current_edge.second].P(), closestEnd) <= eps ) { + map_info[endF].SetEps( eps ); + map_info[endF].Init( *endF, tri::Index( a->cm, endF->V(0) ), tri::Index( a->cm, endF->V(1) ), tri::Index( a->cm, endF->V(2) ) ); + map_info.insert( make_pair(endF, dummy) ); + map_info[endF].addVertex( &(a->cm.vert[current_edge.second]), current_edge.second ); + tbt_faces.push_back( endF ); + } + + //project one of the two vertices on the mesh + if ( startF != 0 ) a->cm.vert[current_edge.first].P() = closestStart; + if ( endF != 0 )a->cm.vert[current_edge.second].P() = closestEnd; + //we're going to split current face + tbr_faces.push_back( f ); + + //choose direction of splitting (from start to end or viceversa) according to which + //vertex doesn't project + int dir = 0; + if ( startF == 0 || isOnBorder( closestStart, startF ) ) dir = 1; + //split face + handleBorderEdgeOB( current_edge, dir, a, grid_a, max_dist, startF, endF, f, map_info, stack, tbt_faces, verts ); + + continue; + } + //end of special cases + + //initialize info on faces (if info is already initialized, no operation is done) + map_info.insert( make_pair(startF, dummy) ); map_info[startF].SetEps( eps ); map_info[startF].Init( *startF, tri::Index( a->cm, startF->V(0) ), tri::Index( a->cm, startF->V(1) ), tri::Index( a->cm, startF->V(2) ) ); + map_info.insert( make_pair(endF, dummy) ); map_info[endF].SetEps( eps ); map_info[endF].Init( *endF, tri::Index( a->cm, endF->V(0) ), tri::Index( a->cm, endF->V(1) ), tri::Index( a->cm, endF->V(2) ) ); + + //project points on the surface of the mesh + a->cm.vert[current_edge.first].P() = closestStart; + a->cm.vert[current_edge.second].P() = closestEnd; + + //case 01 : startF and endF are the same face + if ( startF == endF ) { + //start face will be retriangulated + tbt_faces.push_back( startF ); + //add information to startF + handleBorderEdgeSF( current_edge, a, startF, endF, f, map_info, stack, verts ); + continue; + } + + //case 02: startF is adjacent to endF + if ( isAdjacent( startF, endF ) ) { + //start face and end face will be retriangulated + tbt_faces.push_back( startF ); + tbt_faces.push_back( endF ); + //also, current face will be splitted and replaced by other two faces + tbr_faces.push_back( f ); + //split face and update information + handleBorderEdgeAF( current_edge, a, startF, endF, f, map_info, stack, verts ); + continue; + + } + + //case 03: startF is not adjacent to endF + if ( !isAdjacent( startF, endF ) ) { + //start face and end face will be retriangulated + tbt_faces.push_back( startF ); + tbt_faces.push_back( endF ); + //also, current face will be splitted and replaced by other two faces + tbr_faces.push_back( f ); + //split face and add edges to the stack + handleBorderEdgeNF( current_edge, a, startF, endF, f, map_info, stack, verts ); + continue; + } + + } + + +} + +//case 01: vertices of border edge project on the same face +void FilterZippering::handleBorderEdgeSF ( pair< int, int >& current_edge, //current border edge + MeshModel* a, //mesh A + CMeshO::FacePointer startF, //face where first vertex lies + CMeshO::FacePointer endF, //face where second vertex lies (=startF) + CMeshO::FacePointer splittingF, //splitting face + map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + vector < pair< int, int > >& stack, //stack containing border edges + vector< int >& verts ) { //vector of indices + + //border edge + int e; for ( e = 0; e < 3; e ++ ) if ( face::IsBorder(*splittingF, e) ) break; + + //add information to startF + if( map_info[startF].AddToBorder( Segment3 ( a->cm.vert[current_edge.first].P(), a->cm.vert[current_edge.second].P() ), + make_pair( current_edge.first, current_edge.second ) ) ) { + //avoid duplicate faces + if ( current_edge.first == tri::Index( a->cm, splittingF->V(e) ) && + current_edge.second == tri::Index( a->cm, splittingF->V1(e) ) ) return; + //store indices of vertices of the new face + verts.push_back( current_edge.first ); + verts.push_back( current_edge.second ); + verts.push_back( tri::Index( a->cm, splittingF->V2(e) ) ); + } +} + + + +//case 02: vertices of border edge project on adjacent faces +void FilterZippering::handleBorderEdgeAF ( pair< int, int >& current_edge, //current border edge + MeshModel* a, //mesh A + CMeshO::FacePointer startF, //face where first vertex lies + CMeshO::FacePointer endF, //face where second vertex lies + CMeshO::FacePointer splittingF, //splitting face + map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + vector < pair< int, int > >& stack, //stack containing border edges + vector< int >& verts ) { //vector of indices + + //calc. intersection point (approximate) and split edge + int shared; for ( int k = 0; k < 3; k ++ ) if ( startF->FFp(k) == endF ) shared = k; + Segment3 shared_edge( startF->P(shared), startF->P1(shared) ); + int sampleNum = SAMPLES_PER_EDGE; float step = 1.0/(sampleNum+1); + Point3f closest; float min_dist = shared_edge.Length(); + //subsample the current border edge and search for point which is closest to the edge shared by startF and endF + for ( int k = 0; k <= sampleNum; k ++ ) { + Point3f currentP = a->cm.vert[current_edge.first].P() + ( a->cm.vert[current_edge.second].P() - a->cm.vert[current_edge.first].P() ) * (k*step); + if ( SquaredDistance( shared_edge, currentP ) < min_dist ) { + closest = currentP; + min_dist = SquaredDistance( shared_edge, closest ); + } + } + //have we found something? + assert( SquaredDistance( shared_edge, closest ) < shared_edge.Length() ); + + //closest (pseudo-intersection) point + closest = ClosestPoint(shared_edge, closest); + + //border edge + int e; for ( e = 0; e < 3; e ++ ) if ( face::IsBorder(*splittingF, e) ) break; + + //no close vertices, add information to faces + tri::Allocator::PointerUpdater vpu; + CMeshO::VertexIterator v = tri::Allocator::AddVertices( a->cm, 1, vpu ); (*v).P() = closest; + if ( map_info[startF].AddToBorder( Segment3 ( a->cm.vert[current_edge.first].P(), (*v).P() ), + make_pair( current_edge.first, v - a->cm.vert.begin() ) ) ) { + verts.push_back( v - a->cm.vert.begin() ); + verts.push_back( tri::Index( a->cm, splittingF->V2(e) ) ); + verts.push_back( current_edge.first ); //new triangle + } + if ( map_info[endF].AddToBorder( Segment3 ( (*v).P(), a->cm.vert[current_edge.second].P() ), + make_pair( v - a->cm.vert.begin(), current_edge.second ) ) ) { + verts.push_back(v - a->cm.vert.begin()); + verts.push_back( current_edge.second ); + verts.push_back( tri::Index( a->cm, splittingF->V2(e) ) ); //new triangle + } +} + +//case 03: vertices of border edge project on non-adjacent faces +void FilterZippering::handleBorderEdgeNF ( pair< int, int >& current_edge, //current border edge + MeshModel* a, //mesh A + CMeshO::FacePointer startF, //face where first vertex lies + CMeshO::FacePointer endF, //face where second vertex lies + CMeshO::FacePointer splittingF, //splitting face + map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + vector < pair< int, int > >& stack, //stack containing border edges + vector< int >& verts ) { + //border edge + int e; for ( e = 0; e < 3; e ++ ) if ( face::IsBorder(*splittingF, e) ) break; + + //startF and endF are not adjacent, but they could share a vertex... + int w = sharesVertex( startF, endF ); + int cnt = 0; + Segment3 s( a->cm.vert[current_edge.first].P(), a->cm.vert[current_edge.second].P() ); + //if they share a vertex and border edge pass trough the vertex, we split current face using the vertex + if ( w != -1 && (cnt++ == MAX_LOOP || SquaredDistance( s, startF->P(w) ) <= eps) ) { + //too short and it's a vertex, do nothing + if ( s.Length() < eps && vcg::Distance( s.P0(), startF->P(w) ) < eps ) + return; + tri::Allocator::PointerUpdater vpu; + CMeshO::VertexIterator v = tri::Allocator::AddVertices( a->cm, 1, vpu ); (*v).P() = startF->P(w); + //add information to startF + if ( map_info[startF].AddToBorder( Segment3 ( a->cm.vert[current_edge.first].P(), (*v).P() ), + make_pair( current_edge.first, v - a->cm.vert.begin() ) ) ) { + verts.push_back( current_edge.first ); + verts.push_back( v - a->cm.vert.begin() ); + verts.push_back( tri::Index( a->cm, splittingF->V2(e) ) ); + } + //add information to endF + if ( map_info[endF].AddToBorder( Segment3 ( (*v).P(), a->cm.vert[current_edge.second].P() ), + make_pair( v - a->cm.vert.begin(), current_edge.second ) ) ) { + verts.push_back( v - a->cm.vert.begin() ); + verts.push_back( current_edge.second ); + verts.push_back( tri::Index( a->cm, splittingF->V2(e) ) ); + } + //stop + return; + } + + //add new vertex on the mid-point + tri::Allocator::PointerUpdater vpu; + CMeshO::VertexIterator v = tri::Allocator::AddVertices( a->cm, 1, vpu ); + (*v).P() = (a->cm.vert[current_edge.first].P() + a->cm.vert[current_edge.second].P())/2.00; + stack.push_back( make_pair( current_edge.first, v - a->cm.vert.begin() ) ); + stack.push_back( make_pair( v - a->cm.vert.begin(), current_edge.second ) ); +} + +//case 04: both vertices of border edge project on border edge +//return true if the whole current_edge project on border edge (no split needed) +bool FilterZippering::handleBorderEdgeBB ( std::pair< int, int >& current_edge, //current border edge + MeshModel* a, //mesh A + MeshFaceGrid grid_a, //grid on A (needed for sampling) + float max_dist, //max search dist (needed for sampling) + CMeshO::FacePointer startF, //face where first vertex lies + CMeshO::FacePointer endF, //face where second vertex lies + CMeshO::FacePointer splittingF, //splitting face + std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + std::vector < std::pair< int, int > >& stack, //stack containing border edges + std::vector< int >& verts ) { //vector of indices + + //Verify if the whole segment is on border + tri::FaceTmark markerFunctor; markerFunctor.SetMesh(&a->cm); + face::PointDistanceBaseFunctor PDistFunct; + int sampleNum = SAMPLES_PER_EDGE; float step = 1.0/(sampleNum+1); bool border = true; + Point3f closestP; float dist = 2*max_dist; + for ( int k = 0; k <= sampleNum; k ++ ) { + Point3f currentP = a->cm.vert[current_edge.first].P() + ( a->cm.vert[current_edge.second].P() - a->cm.vert[current_edge.first].P() ) * (k*step); + CMeshO::FacePointer closestFace = grid_a.GetClosest(PDistFunct, markerFunctor, currentP, 2*max_dist, dist, closestP); //closest point on mesh + if ( !isOnBorder( closestP, closestFace ) ) return false; //not completely on border; will be splitted later + } + + //border edge of splittingF + int e; for ( e = 0; e < 3; e ++ ) if ( face::IsBorder(*splittingF, e) ) break; + + //actualFace: no operation needed + if ( a->cm.vert[current_edge.first].P() == splittingF->P(e) && + a->cm.vert[current_edge.second].P() == splittingF->P1(e) ) {} + else {//create new triangle to file small hole + verts.push_back( current_edge.second ); + verts.push_back( tri::Index( a->cm, splittingF->V2(e) ) ); + verts.push_back( current_edge.first ); //new triangle + } + + return true; +} + + + +//case 05: one of the vertices doesn't project on the surface of the mesh +void FilterZippering::handleBorderEdgeOB ( std::pair< int, int >& current_edge, //current border edge + int direction, //split direction (1 from start to end, 0 from end to start) + MeshModel* a, //mesh A + MeshFaceGrid grid_a, //grid on A (needed for sampling) + float max_dist, //max search dist (needed for sampling) + CMeshO::FacePointer startF, //face where first vertex lies + CMeshO::FacePointer endF, //face where second vertex lies + CMeshO::FacePointer splittingF, //splitting face + map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + vector < pair< int, int > >& stack, //stack containing border edges + vector< CMeshO::FacePointer >& tbt_faces, //stack containing pointers to face that wille be retriangulated + vector< int >& verts ) { //vector of indices + + assert( direction == 1 || direction == 0 ); + //index of last splitted edge + int last_split = -1; + //face currently examinated + CMeshO::FacePointer currentF = direction? startF : endF; + //stop flag + bool stop = false; + + //border edge of splittingF + int e; for ( e = 0; e < 3; e ++ ) if ( face::IsBorder(*splittingF, e) ) break; + + //indices of vertices involved in the process + int startV = direction ? current_edge.first : current_edge.second; + int endV = direction ? current_edge.second : current_edge.first; + int thirdV = tri::Index( a->cm, splittingF->V2(e) ); + + int cnt = 0; + aux_info dummy; + + do { + int tosplit; + Point3f closest; + cnt++; + if (!findIntersection( currentF, Segment3(a->cm.vert[current_edge.first].P(),a->cm.vert[current_edge.second].P()), last_split, tosplit, closest )) { + stop = true; //no op + } + else { + //we found an intersection, so split the face using intersection point + last_split = currentF->FFi( tosplit ); + tri::Allocator::PointerUpdater vpu; + //new vertex on the splitting point + CMeshO::VertexIterator v = tri::Allocator::AddVertices( a->cm, 1, vpu ); (*v).P() = closest; + //init currentF aux_information + map_info.insert( make_pair(currentF, dummy) ); + map_info[currentF].Init( *currentF, tri::Index(a->cm, currentF->V(0)), tri::Index(a->cm, currentF->V(1)), tri::Index(a->cm, currentF->V(2)) ); + //add information to currentF + if ( map_info[currentF].AddToBorder( Segment3 ( (*v).P(), a->cm.vert[startV].P() ), make_pair( v - a->cm.vert.begin(), startV ) ) ) { + tbt_faces.push_back( currentF ); + verts.push_back( startV ); + verts.push_back( thirdV ); + verts.push_back( v - a->cm.vert.begin() ); + + } + //startV = v - a->cm.vert.begin(); //why? + if ( face::IsBorder( *currentF, tosplit ) ) { //border reached + verts.push_back( thirdV ); + verts.push_back( endV ); + verts.push_back( v - a->cm.vert.begin() ); + + stack.push_back( make_pair( endV, v - a->cm.vert.begin() ) ); + + stop = true; + } + else currentF = currentF->FFp(tosplit); + } + } while (!stop && cnt < MAX_LOOP); +} + + +/* Zippering of two meshes + * Given two mesh, a mesh with one or more holes (A) and a second mesh, a patch (B), fill a hole onto m surface + * using faces of patch. Algorithm const of three steps: + * - CheckRedundancy: delete face of patch that can be projected completely on A's surface + * - PatchRefinement: patch vertices are projected on A, and border edges are splitted; + * - MeshRefinement: faces where patch vertices lie are re-triangulated + * Based on Controlled and Adaptive mesh zippering + */ +bool FilterZippering::applyFilter(QAction *filter, MeshDocument &md, RichParameterSet & par, CallBackPos /**cb*/) +{ + //get the meshes + MeshModel* a = par.getMesh("FirstMesh"); + MeshModel* b = par.getMesh("SecondMesh"); + + if ( a == b ) { + errorMessage="Please add a second mesh"; + return false; + } + + //enable FF adjacency, mark, compute normals for face (both on A and B) + a->updateDataMask(MeshModel::MM_FACEFACETOPO + MeshModel::MM_FACEMARK + MeshModel::MM_FACEQUALITY + MeshModel::MM_FACECOLOR + MeshModel::MM_VERTQUALITY); + b->updateDataMask(MeshModel::MM_FACEFACETOPO + MeshModel::MM_FACEMARK + MeshModel::MM_FACEQUALITY + MeshModel::MM_FACECOLOR + MeshModel::MM_VERTQUALITY); + tri::UnMarkAll(a->cm); + tri::UnMarkAll(b->cm); + + tri::UpdateNormals::PerFaceNormalized(a->cm); + tri::UpdateFlags::FaceProjection(a->cm); + tri::UpdateNormals::PerVertexNormalized(a->cm); + tri::UpdateNormals::PerFaceNormalized(b->cm); + tri::UpdateFlags::FaceProjection(b->cm); + tri::UpdateNormals::PerVertexNormalized(b->cm); + //fixed eps + eps = 0.00001f; + + clock_t t1, t2; + t1 = clock(); + + /** + * Filter Redundancy: Select redundant faces from the surface of the meshes. + */ + if (ID(filter) == FP_REDUNDANCY) { + //queues + vector< pair > generic_queue; //unsorted queue + priority_queue< pair, vector< pair >, compareFaceQuality > priority_queue; //priority queue + + //if user chooses to usequality, initialize priority queue... + if ( par.getBool("UseQuality") ) { + if ( !Init_pq( priority_queue, a, b, par.getBool("FullProcessing") ) ) { + Log("Target mesh has no border faces - Please select Full Processing option"); + return false; + } + } + //...else initialize unsorted queue + else { + if ( !Init_q( generic_queue, a, b, par.getBool("FullProcessing") ) ) { + Log("Target mesh has no border faces - Please select Full Processing option"); + return false; + } + } + Log( "elapsed time: %f", (clock()-t1)/CLOCKS_PER_SEC ); + CMeshO::ScalarType epsilon = par.getAbsPerc("distance"); + int sf; + //if ( par.getBool("FastErosion") && par.getBool("UseQuality") ) selectRedundantFast_pq(); + //if ( par.getBool("FastErosion") && !par.getBool("UseQuality") ) selectRedundantFast(); + if ( /*!par.getBool("FastErosion") &&*/ par.getBool("UseQuality") ) { + Log( "Quality-based Redundancy" ); + sf = selectRedundant_pq( priority_queue, a, b, epsilon ); + } + if ( /*!par.getBool("FastErosion") &&*/ !par.getBool("UseQuality") ) { + Log( "Standard Redunancy" ); + sf = selectRedundant( generic_queue, a, b, epsilon ); + Log( "elapsed time: %f", (t2-t1)/CLOCKS_PER_SEC ); + } + + Log( "Selected %i redundant faces", sf ); + return true; + } //end FP_REDUNDANCY + + /** + * Filter zippering + * Assuming that redundant faces has already been removed, + * we project the border of mesh B on the surface of mesh A. + */ + if (ID(filter) == FP_ZIPPERING) { + //pre-processing step: visit the border of B and split faces with two border edges + refineBorder( b ); + //clean flags + tri::UpdateFlags::FaceClear(a->cm); + tri::UpdateFlags::FaceClear(b->cm); + //before append, clean mesh A + tri::Clean::RemoveUnreferencedVertex( a->cm ); + tri::Clean::RemoveDuplicateVertex( a->cm ); + //before append, clean mesh B + tri::Clean::RemoveUnreferencedVertex( b->cm ); + tri::Clean::RemoveDuplicateVertex( b->cm ); + //store face number + size_t fn_limit = a->cm.fn; + //append B to A and update flags + tri::UpdateTopology::FaceFace(a->cm); + tri::UpdateTopology::FaceFace(b->cm); + tri::Append::Mesh( a->cm, b->cm ); + tri::UpdateTopology::FaceFace(a->cm); + tri::UpdateFlags::FaceClear(a->cm); + tri::UpdateNormals::PerFaceNormalized(a->cm); + tri::UpdateFlags::FaceProjection(a->cm); + tri::UpdateNormals::PerVertexNormalized(a->cm); + + //create grid on mesh A + MeshFaceGrid grid; + grid.Set( a->cm.face.begin(), a->cm.face.begin() + fn_limit ); //compute grid on the original part of A + //recover information about border of B + vector< tri::Hole::Info > border; + tri::Hole::GetInfo( a->cm, false, border ); + + //Add optional attribute (true if face has been visited) + CMeshO::PerFaceAttributeHandle visited = tri::Allocator::AddPerFaceAttribute (a->cm); //check for already visited face + for ( size_t i = 0; i < a->cm.face.size(); i ++) visited[i] = false; //no face initially visited + + //create map between faces and auxiliar information + map map_info; + vector< CMeshO::FacePointer > tbt_faces; //To Be Triangulated + vector< CMeshO::FacePointer > tbr_faces; //To Be Removed + vector< int > verts; //vector containing indices of vertices of the new faces + //for each border of the mesh B, select a face and projected it on the surface of A + //add information to the faces of A when needed + for ( size_t c = 0; c < border.size(); c ++ ) { + //create pos on the surface + face::Pos p = border[c].p; p.FlipV();//CCW order + //face is from A border, go over + if ( tri::Index( a->cm, p.F() ) < fn_limit ) continue; + //face has already been deleted or visiter,go over + if ( visited[p.F()] || p.F()->IsD() ) continue; + do { + //Already deleted or visited, step over + if ( visited[p.F()] || p.F()->IsD()) { p.NextB(); continue; } + visited[p.F()] = true; + //store vertex index (avoid crash) + int v_ind = tri::Index( a->cm, p.V() ); + //project current face on the surface of A + projectFace( p.F(), a, grid, par.getFloat("distance"), map_info, tbt_faces, tbr_faces, verts ); + //restore vertex pointer + p.V() = &a->cm.vert[v_ind]; + p.NextB(); + } while ( p.F() != border[c].p.F() ); + } + + //remove splitted faces + for ( size_t i = 0; i < tbr_faces.size(); i++) { + if (!tbr_faces[i]->IsD()) tri::Allocator::DeleteFace(a->cm, *tbr_faces[i]); + } + + //remove duplicates from vector to-be-triangulated + sort( tbt_faces.begin(), tbt_faces.end() ); int k; + vector::iterator itr = unique( tbt_faces.begin(), tbt_faces.end() ); + tbt_faces.resize(itr - tbt_faces.begin() ); + //retriangulation of the faces and removal of the remaining part + vector< Point3f > coords; size_t patch_verts = verts.size(); + for ( size_t i = 0; i < tbt_faces.size(); i ++ ) { + if ( !tbt_faces[i]->IsD() ) { + handleBorder( map_info[tbt_faces[i]], tbt_faces[i]->N(), coords, verts ); + tri::Allocator::DeleteFace(a->cm, *tbt_faces[i]); + } + } + + //remove user-defined attribute + int debcnt = 0; + tri::Allocator::DeletePerFaceAttribute(a->cm, visited); + //add new faces and remove old ones + CMeshO::FaceIterator fn = tri::Allocator::AddFaces( a->cm, verts.size()/3 ); + for ( size_t k = 0; k < verts.size(); k += 3 ) { + CMeshO::VertexPointer v0 = &(a->cm.vert[verts[k]]); + CMeshO::VertexPointer v1 = &(a->cm.vert[verts[k+1]]); + CMeshO::VertexPointer v2 = &(a->cm.vert[verts[k+2]]); + if ( v0 == v1 || v1 == v2 || v0 == v2 ) + continue; + //correct orientation of face normals + if ( k < patch_verts ) { + (*fn).V(0) = v0; (*fn).V(1) = v1; (*fn).V(2) = v2; (*fn).N() = ( (*fn).P(0) - (*fn).P(2) )^( (*fn).P(1)-(*fn).P(2) ); + } + else { + (*fn).V(0) = v0; (*fn).V(2) = v1; (*fn).V(1) = v2; (*fn).N() = ( (*fn).P(0) - (*fn).P(2) )^( (*fn).P(1)-(*fn).P(2) ); + } + (*fn).SetS(); + ++fn; + } + + //remove degenerate faces + tri::Clean::RemoveDegenerateFace(a->cm); + tri::Clean::RemoveDuplicateVertex(a->cm); + //compact vectors and update topology + tri::Allocator::CompactFaceVector( a->cm ); + tri::Allocator::CompactVertexVector( a->cm ); + tri::UpdateTopology::FaceFace(a->cm); + tri::UpdateSelection::VertexFromFaceLoose(a->cm); + tri::UpdateSelection::FaceFromVertexLoose(a->cm); + tri::UpdateBounding::Box( a->cm ); + + return true; + } //end FP_ZIPPERING + + assert(0); + Log(GLLogStream::DEBUG, "**********************" ); + + return true; +} + +Q_EXPORT_PLUGIN(FilterZippering) diff --git a/src/meshlabplugins/filter_zippering/filter_zippering.h b/src/meshlabplugins/filter_zippering/filter_zippering.h index 74e50f215..23bb6354c 100644 --- a/src/meshlabplugins/filter_zippering/filter_zippering.h +++ b/src/meshlabplugins/filter_zippering/filter_zippering.h @@ -1,470 +1,470 @@ -/**************************************************************************** -* MeshLab o o * -* A versatile mesh processing toolbox o o * -* _ O _ * -* Copyright(C) 2005 \/)\/ * -* Visual Computing Lab /\/| * -* ISTI - Italian National Research Council | * -* \ * -* All rights reserved. * -* * -* This program is free software; you can redistribute it and/or modify * -* it under the terms of the GNU General Public License as published by * -* the Free Software Foundation; either version 2 of the License, or * -* (at your option) any later version. * -* * -* This program is distributed in the hope that it will be useful, * -* but WITHOUT ANY WARRANTY; without even the implied warranty of * -* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * -* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) * -* for more details. * -* * -****************************************************************************/ -/**************************************************************************** - -****************************************************************************/ - -#ifndef FILTERZIPPERING_H -#define FILTERZIPPERING_H - -#include - -#include -#include -#include -#include - -#define SAMPLES_PER_EDGE 5 //modificare, length/epsilon - -// Polyline (set of consecutive segments) -struct polyline { - std::vector< vcg::Segment3 > edges; //polyline edges - std::vector< std::pair > verts; -}; - -//Auxiliar info for triangulation -struct aux_info { - std::vector< polyline > conn; //Close components (to be triangulated) - std::vector< polyline > trash; //Close components (to be deleted) - std::vector< polyline > border; //Segment intersecting components - float eps; //epsilon - //Add segment c to border - virtual bool AddToBorder( vcg::Segment3 c, std::pair v ) { - /****Insert new segment****/ - //Search for segment S having S.P0() == c.P1 or S.P1 == c.P0() - if ( v.first == 4013 && v.second == 4015 ) - int stop = 3; - if ( v.first == v.second ) return false; - if ( c.Length() < eps ) { - c.P0() = c.P1(); - v.first = v.second; - //if segment is too short (it's basically a point) check if it's one of the vertices - for ( size_t i = 0; i < trash.size(); i ++ ) { - for ( size_t j = 0; j < trash[i].edges.size(); j ++ ) { //Only one trash component - if ( vcg::Distance( trash[i].edges[j].P0(), c.P0() ) < eps ) return false; - } - } - } - - if ( v.first == 4013 && v.second == 4013 ) - int stop = 3; - //check if c doesn't lie on existing edges - //if it does, split edge - for ( size_t i = 0; i < trash.size(); i ++ ) { - for ( size_t j = 0; j < trash[i].edges.size(); j ++ ) { //Only one trash component - vcg::Segment3 a, b; - a = trash[i].edges[j]; b = c; - float tx0 = (float) (b.P0().X() - a.P0().X())/(a.P1().X() - a.P0().X()); - float ty0 = (float) (b.P0().Y() - a.P0().Y())/(a.P1().Y() - a.P0().Y()); - float tz0 = (float) (b.P0().Z() - a.P0().Z())/(a.P1().Z() - a.P0().Z()); - if ( (fabs(tx0 - ty0) < eps ) && (fabs(ty0 - tz0) < eps) ) { - float tx1 = (float) (b.P1().X() - a.P0().X())/(a.P1().X() - a.P0().X()); - float ty1 = (float) (b.P1().Y() - a.P0().Y())/(a.P1().Y() - a.P0().Y()); - float tz1 = (float) (b.P1().Z() - a.P0().Z())/(a.P1().Z() - a.P0().Z()); - if ( (fabs(tx1 - ty1) < eps) && (fabs(ty1 - tz1) < eps) ) { - if ( (tx0 < 0.0f + eps) && (tx1 > 1.0f - eps) ) return false; //external - if ( (tx1 < 0.0f + eps) && (tx0 > 1.0f - eps) ) return false; //external - if ( (tx0 >= 0.0f + eps) && (tx0 <= 1.0f - eps) && (tx1 >= 0.0f + eps) && (tx1 <= 1.0f - eps) ) { - //double split - if ( tx0 > tx1 ) { c.Flip(); v = std::make_pair( v.second, v.first ); } - //insert new edges - trash[i].edges.insert( trash[i].edges.begin() + j + 1, vcg::Segment3( trash[i].edges[j].P0(), c.P0() ) ); - trash[i].verts.insert( trash[i].verts.begin() + j + 1, std::make_pair( trash[i].verts[j].first, v.first ) ); - trash[i].edges.insert( trash[i].edges.begin() + j + 2, c ); - trash[i].verts.insert( trash[i].verts.begin() + j + 2, v ); - trash[i].edges.insert( trash[i].edges.begin() + j + 3, vcg::Segment3( c.P1(), trash[i].edges[j].P1() ) ); - trash[i].verts.insert( trash[i].verts.begin() + j + 3, std::make_pair( v.second, trash[i].verts[j].second ) ); - //erase old one - trash[i].edges.erase( trash[i].edges.begin() + j ); - trash[i].verts.erase( trash[i].verts.begin() + j ); - return true; - } - if ( (tx0 >= 0.0f + eps) && (tx0 <= 1.0f - eps) ) { - //single split in P0 - //insert new edges - trash[i].edges.insert( trash[i].edges.begin() + j + 1, vcg::Segment3( trash[i].edges[j].P0(), c.P0() ) ); - trash[i].verts.insert( trash[i].verts.begin() + j + 1, std::make_pair( trash[i].verts[j].first, v.first ) ); - trash[i].edges.insert( trash[i].edges.begin() + j + 2, vcg::Segment3( c.P0(), trash[i].edges[j].P1() ) ); - trash[i].verts.insert( trash[i].verts.begin() + j + 2, std::make_pair( v.first, trash[i].verts[j].second ) ); - //erase old one - trash[i].edges.erase( trash[i].edges.begin() + j ); - trash[i].verts.erase( trash[i].verts.begin() + j ); - return true; - } - if ( (tx1 >= 0.0f + eps) && (tx1 <= 1.0f - eps) ) { - //single split in P1 - //insert new edges - trash[i].edges.insert( trash[i].edges.begin() + j + 1, vcg::Segment3( trash[i].edges[j].P0(), c.P1() ) ); - trash[i].verts.insert( trash[i].verts.begin() + j + 1, std::make_pair( trash[i].verts[j].first, v.second ) ); - trash[i].edges.insert( trash[i].edges.begin() + j + 2, vcg::Segment3( c.P1(), trash[i].edges[j].P1() ) ); - trash[i].verts.insert( trash[i].verts.begin() + j + 2, std::make_pair( v.second, trash[i].verts[j].second ) ); - //erase old one - trash[i].edges.erase( trash[i].edges.begin() + j ); - trash[i].verts.erase( trash[i].verts.begin() + j ); - return true; - } - } - } - } - } - - bool found = false; - for ( unsigned int j = 0; j < border.size(); j ++ ) { - for ( size_t i = 0; i < border[j].verts.size() && !found; i ++ ) { - if ( border[j].verts[i].first == v.second ) { found = true; border[j].edges.insert( border[j].edges.begin() + i, c ); border[j].verts.insert( border[j].verts.begin() + i, v ); } //insert before i-th element - else - if ( border[j].verts[i].second == v.first ) { found = true; border[j].edges.insert( border[j].edges.begin() + i + 1, c ); border[j].verts.insert( border[j].verts.begin() + i + 1, v ); } //insert after i-th element - } - } - if (!found) { //Create a new polyline ad add it to the border list - polyline nwpoly; nwpoly.edges.push_back( c ); nwpoly.verts.push_back( v ); border.push_back( nwpoly ); - } else { - //Merge consecutive polylines into a single one - for ( unsigned int j = 0; j < border.size(); j ++ ) - for ( unsigned int i = j+1; i < border.size(); i ++ ) { - if ( border[j].verts.front().first == border[i].verts.back().second ) { - border[j].edges.insert( border[j].edges.begin(), border[i].edges.begin(), border[i].edges.end() ); - border[j].verts.insert( border[j].verts.begin(), border[i].verts.begin(), border[i].verts.end() ); - border.erase(border.begin() + i); - } - else if ( border[j].verts.back().second == border[i].verts.front().first ) { - border[j].edges.insert( border[j].edges.end(), border[i].edges.begin(), border[i].edges.end() ); - border[j].verts.insert( border[j].verts.end(), border[i].verts.begin(), border[i].verts.end() ); - border.erase(border.begin() + i); - } - } - }//end if (!found) - - for ( size_t k = 0; k < border.size(); k ++) { - for ( size_t i = 0; i < trash.size(); i ++ ) { - for ( size_t j = 0; j < trash[i].verts.size(); j ++ ) { //Only one trash component - if ( trash[i].verts[j].first == border[k].verts.back().second ) { - trash[i].edges[j].P0() = border[k].edges.back().P1(); - //trash[i].edges.erase( trash[i].edges.begin() + j ); trash[i].verts.erase( trash[i].verts.begin() + j ) - } else - if ( trash[i].verts[j].first == border[k].verts.front().first ) { - trash[i].edges[j].P0() = border[k].edges.front().P0(); - //trash[i].edges.erase( trash[i].edges.begin() + j ); trash[i].verts.erase( trash[i].verts.begin() + j ) - } - } - } - } - return true; - }//end AddToBorder - - // Add c.component - virtual void AddCComponent( polyline c ) { - conn.push_back(c); - } - // Add t.component - virtual void AddTComponent( polyline t ) { - trash.push_back(t); - } - - //Set eps - virtual void SetEps( float e ) { - eps = e; - } - - // Set initial t.component - virtual void Init( CMeshO::FaceType f, int a, int b, int c ) { - if (!trash.empty()) return; - polyline tri; tri.edges.push_back( vcg::Segment3(f.P(0), f.P(1)) ); - tri.edges.push_back( vcg::Segment3(f.P(1), f.P(2)) ); - tri.edges.push_back( vcg::Segment3(f.P(2), f.P(0)) ); - tri.verts.push_back( std::make_pair(a, b) ); - tri.verts.push_back( std::make_pair(b, c) ); - tri.verts.push_back( std::make_pair(c, a) ); - AddTComponent( tri ); - } - - // Remove c.component - virtual void RemoveCComponent( int i ) { - conn.erase( conn.begin() + i ); - } - - // Remove t.component - virtual void RemoveTComponent( int i ) { - trash.erase( trash.begin() + i ); - } - - // Number of c.component - virtual int nCComponent( ) { - return conn.size(); - } - - // Number of t.component - virtual size_t nTComponent( ) { - return trash.size(); - } - - // Add vertex in original triangle - virtual bool addVertex( CMeshO::VertexPointer v, int v_index ) { - int cnt = 0; int split = -1; - for ( size_t i = 0; i < trash.size(); i ++ ) { //one component only - - for ( size_t j = 0; j < trash[i].verts.size(); j ++ ) { //search for closest edge - if ( trash[i].verts[j].first == v_index ) return false; - } - - for ( size_t j = 0; j < trash[i].edges.size(); j ++ ) { //search for closest edge - //if ( vcg::SquaredDistance( trash[i].edges[j], v->P() ) <= eps ) { - float dist; - vcg::Point3f clos; - vcg::SegmentPointSquaredDistance(trash[i].edges[j], v->P(),clos,dist); - if (dist <= eps ) { - cnt++; split = j; - } - } - } - if (!cnt) return false; - - if ( cnt == 1 ) { //one edge only -> split the edge - vcg::Segment3 splitting_edge = trash[0].edges[split]; - std::pair splitting_edge_v = trash[0].verts[split]; - - trash[0].edges.erase(trash[0].edges.begin()+split); //remove edge - trash[0].verts.erase(trash[0].verts.begin()+split); - //replace edge using two new edges - trash[0].edges.insert(trash[0].edges.begin()+split, vcg::Segment3( splitting_edge.P0(), v->P() ) ); - trash[0].edges.insert(trash[0].edges.begin()+split+1, vcg::Segment3( v->P(), splitting_edge.P1() ) ); - trash[0].verts.insert(trash[0].verts.begin()+split, std::make_pair( splitting_edge_v.first, v_index ) ); - trash[0].verts.insert(trash[0].verts.begin()+split+1, std::make_pair( v_index, splitting_edge_v.second ) ); - } - - if ( cnt == 2 ) { // search for closest vertex and copy vertex coords - for ( size_t i = 0; i < trash.size(); i ++ ) { //one component only - for ( size_t j = 0; j < trash[i].edges.size(); j ++ ) { //search for closest edge - if ( vcg::Distance( trash[i].edges[j].P0(), v->P() ) <= eps ) { - v->P() = trash[i].edges[j].P0(); - } - } - } - } - - return true; - } - -};//end struct - -class compareFaceQuality { - public: - compareFaceQuality() { }; - - bool operator () (const std::pair f1, const std::pair f2) { - //quality f1 < quality f2 return true - return ( f1.first->Q() > f2.first->Q() ); - } -}; - -class FilterZippering : public QObject, public MeshFilterInterface -{ - Q_OBJECT - Q_INTERFACES(MeshFilterInterface) - - typedef vcg::GridStaticPtr MeshFaceGrid; - typedef vcg::GridStaticPtr MeshVertGrid; - -public: - //Different operations in different plugins - enum { FP_REDUNDANCY, - FP_ZIPPERING }; - - FilterZippering(); - - virtual QString filterName(FilterIDType filter) const; - virtual QString filterInfo(FilterIDType filter) const; - virtual void initParameterSet(QAction *,MeshDocument &/*m*/, RichParameterSet & /*parent*/); - int getRequirements(QAction *action); - virtual bool applyFilter(QAction *filter, MeshDocument &md, RichParameterSet & /*parent*/, vcg::CallBackPos * cb) ; - FilterClass getClass(QAction *a); - virtual int postCondition( QAction */*a*/ ) const { return MeshModel::MM_FACEFACETOPO|MeshModel::MM_VERTNORMAL; } - - - -private: - template - ScalarType SquaredDistance( vcg::Segment3 &s, vcg::Point3 &p); - - template - vcg::Point3 ClosestPoint( vcg::Segment3 &s, vcg::Point3 &p); - - bool checkRedundancy( CMeshO::FacePointer f, //face - MeshModel *a, //mesh A - MeshFaceGrid &grid, //grid A - CMeshO::ScalarType max_dist ); //Max search distance - bool simpleCheckRedundancy( CMeshO::FacePointer f, //face - MeshModel *a, //mesh A - MeshFaceGrid &grid, //grid A - CMeshO::ScalarType max_dist,//Max search distance - bool test ); - bool isBorderVert( CMeshO::FacePointer f, int i); - bool isOnBorder( CMeshO::CoordType point, CMeshO::FacePointer f ); - bool isOnEdge( CMeshO::CoordType point, CMeshO::FacePointer f ); - bool isAdjacent( CMeshO::FacePointer f1, CMeshO::FacePointer f2 ); - int sharesVertex( CMeshO::FacePointer f1, CMeshO::FacePointer f2 ); - void handleBorder( aux_info &info, //Auxiliar information for triangulatio - vcg::Point3 N, //face normal (useful for proiection) - std::vector &coords, //output coords - std::vector &output ); //output v. pointers - polyline cutComponent( polyline comp, //Component to be cut - polyline border, //border - vcg::Matrix44 rot_mat ); //Rotation matrix - int searchComponent( aux_info &info, //Auxiliar info - vcg::Point3 P0, //Start border point - vcg::Point3 P1, //End border point - bool &conn ); - bool findIntersection( CMeshO::FacePointer currentF, //face - vcg::Segment3 edge, //edge - int last_split, //last splitted edge - int &splitted_edge, //currently splitted edge - vcg::Point3 &hit ); //approximate intersection point - - - - //init unsorted queue - bool Init_q( std::vector< std::pair >& queue, //the queue - MeshModel* a, //mesh A - MeshModel* b, //mesh B - bool fullProcess ); //fullProcess flag - //init priority queue (overload) - bool Init_pq( std::priority_queue< std::pair, std::vector< std::pair >, compareFaceQuality >& queue, //the queue - MeshModel* a, //mesh A - MeshModel* b, //mesh B - bool fullProcess ); - //select redundant face (normal method, unsorted queue) - int selectRedundant( std::vector< std::pair >& queue, //queue - MeshModel* a, //mesh A - MeshModel* b, //mesh B - float epsilon ); //max search distance - //select redundant face (normal method, priority queue) - int selectRedundant_pq( std::priority_queue< std::pair, std::vector< std::pair >, compareFaceQuality >& queue, //the queue - MeshModel* a, //mesh A - MeshModel* b, //mesh B - float epsilon ); //max search distance - //refine border of a mesh, splitting faces having two border edges - int refineBorder( MeshModel* m ); - //project face of B on the surface of A - void projectFace( CMeshO::FacePointer f, //pointer to the face that will be projected - MeshModel* a, //mesh A - MeshFaceGrid grid_a, //grid on A - float max_dist, //max dist search - std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - std::vector< CMeshO::FacePointer >& tbt_faces, //vector to-be-triangulated faces - std::vector< CMeshO::FacePointer >& tbr_faces, //vector to-be-removed faces - std::vector< int >& verts ); //vector of indices - - //case 01: vertices of border edge project on the same face - void handleBorderEdgeSF ( std::pair< int, int >& current_edge, //current border edge - MeshModel* a, //mesh A - CMeshO::FacePointer startF, //face where first vertex lies - CMeshO::FacePointer endF, //face where second vertex lies - CMeshO::FacePointer splittingF, //splitting face - std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - std::vector < std::pair< int, int > >& stack, //stack containing border edges - std::vector< int >& verts ); //vector of indices - - //case 02: vertices of border edge project on adjacent faces - void handleBorderEdgeAF ( std::pair< int, int >& current_edge, //current border edge - MeshModel* a, //mesh A - CMeshO::FacePointer startF, //face where first vertex lies - CMeshO::FacePointer endF, //face where second vertex lies - CMeshO::FacePointer splittingF, //splitting face - std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - std::vector < std::pair< int, int > >& stack, //stack containing border edges - std::vector< int >& verts ); //vector of indices - - //case 03: vertices of border edge project on non-adjacent faces - void handleBorderEdgeNF ( std::pair< int, int >& current_edge, //current border edge - MeshModel* a, //mesh A - CMeshO::FacePointer startF, //face where first vertex lies - CMeshO::FacePointer endF, //face where second vertex lies - CMeshO::FacePointer splittingF, //splitting face - std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - std::vector < std::pair< int, int > >& stack, //stack containing border edges - std::vector< int >& verts ); //vector of indices - - //case 04: both vertices of current_edge project on another border edge - //return true if the whole current_edge project on border edge - bool handleBorderEdgeBB ( std::pair< int, int >& current_edge, //current border edge - MeshModel* a, //mesh A - MeshFaceGrid grid_a, //grid on A (needed for sampling) - float max_dist, //max search dist (needed for sampling) - vcg::Point3 closestStart, //closest point on startF - vcg::Point3 closestEnd, //closest point on endF - CMeshO::FacePointer startF, //face where first vertex lies - CMeshO::FacePointer endF, //face where second vertex lies - CMeshO::FacePointer splittingF, //splitting face - std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - std::vector < std::pair< int, int > >& stack, //stack containing border edges - std::vector< int >& verts ); //vector of indices - - - //case 04: both vertices of current_edge project on another border edge - //return true if the whole current_edge project on border edge - bool handleBorderEdgeBB ( std::pair< int, int >& current_edge, //current border edge - MeshModel* a, //mesh A - MeshFaceGrid grid_a, //grid on A (needed for sampling) - float max_dist, //max search dist (needed for sampling) - CMeshO::FacePointer startF, //face where first vertex lies - CMeshO::FacePointer endF, //face where second vertex lies - CMeshO::FacePointer splittingF, //splitting face - std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - std::vector < std::pair< int, int > >& stack, //stack containing border edges - std::vector< int >& verts ); //vector of indices - - //case 05: one of the vertices doesn't project on the surface of the mesh - void handleBorderEdgeOB ( std::pair< int, int >& current_edge, //current border edge - int direction, //splitting direction (1 from start to end, 0 from end to start) - MeshModel* a, //mesh A - MeshFaceGrid grid_a, //grid on A (needed for sampling) - float max_dist, //max search dist (needed for sampling) - CMeshO::FacePointer startF, //face where first vertex lies - CMeshO::FacePointer endF, //face where second vertex lies - CMeshO::FacePointer splittingF, //splitting face - std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information - std::vector < std::pair< int, int > >& stack, //stack containing border edges - std::vector< CMeshO::FacePointer >& tbt_faces, //stack containing pointers to face that wille be retriangulated - std::vector< int >& verts ); //vector of indices - - int preProcess ( std::vector< std::pair >& queue, //queue - MeshModel* a, - MeshModel* b, - MeshFaceGrid grid_a, //grid on A - MeshFaceGrid grid_b, //grid on A - float max_dist ); //max dist search - - int preProcess_pq ( std::priority_queue< std::pair, std::vector< std::pair >, compareFaceQuality >& queue, //the queue - MeshModel* a, - MeshModel* b, - MeshFaceGrid grid_a, //grid on A - MeshFaceGrid grid_b, //grid on A - float max_dist ); //max dist search - - - float eps; - - int dbg_cnt; -}; - -#endif +/**************************************************************************** +* MeshLab o o * +* A versatile mesh processing toolbox o o * +* _ O _ * +* Copyright(C) 2005 \/)\/ * +* Visual Computing Lab /\/| * +* ISTI - Italian National Research Council | * +* \ * +* All rights reserved. * +* * +* This program is free software; you can redistribute it and/or modify * +* it under the terms of the GNU General Public License as published by * +* the Free Software Foundation; either version 2 of the License, or * +* (at your option) any later version. * +* * +* This program is distributed in the hope that it will be useful, * +* but WITHOUT ANY WARRANTY; without even the implied warranty of * +* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * +* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) * +* for more details. * +* * +****************************************************************************/ +/**************************************************************************** + +****************************************************************************/ + +#ifndef FILTERZIPPERING_H +#define FILTERZIPPERING_H + +#include + +#include +#include +#include +#include + +#define SAMPLES_PER_EDGE 5 //modificare, length/epsilon + +// Polyline (set of consecutive segments) +struct polyline { + std::vector< vcg::Segment3 > edges; //polyline edges + std::vector< std::pair > verts; +}; + +//Auxiliar info for triangulation +struct aux_info { + std::vector< polyline > conn; //Close components (to be triangulated) + std::vector< polyline > trash; //Close components (to be deleted) + std::vector< polyline > border; //Segment intersecting components + float eps; //epsilon + //Add segment c to border + virtual bool AddToBorder( vcg::Segment3 c, std::pair v ) { + /****Insert new segment****/ + //Search for segment S having S.P0() == c.P1 or S.P1 == c.P0() + if ( v.first == 4013 && v.second == 4015 ) + int stop = 3; + if ( v.first == v.second ) return false; + if ( c.Length() < eps ) { + c.P0() = c.P1(); + v.first = v.second; + //if segment is too short (it's basically a point) check if it's one of the vertices + for ( size_t i = 0; i < trash.size(); i ++ ) { + for ( size_t j = 0; j < trash[i].edges.size(); j ++ ) { //Only one trash component + if ( vcg::Distance( trash[i].edges[j].P0(), c.P0() ) < eps ) return false; + } + } + } + + if ( v.first == 4013 && v.second == 4013 ) + int stop = 3; + //check if c doesn't lie on existing edges + //if it does, split edge + for ( size_t i = 0; i < trash.size(); i ++ ) { + for ( size_t j = 0; j < trash[i].edges.size(); j ++ ) { //Only one trash component + vcg::Segment3 a, b; + a = trash[i].edges[j]; b = c; + float tx0 = (float) (b.P0().X() - a.P0().X())/(a.P1().X() - a.P0().X()); + float ty0 = (float) (b.P0().Y() - a.P0().Y())/(a.P1().Y() - a.P0().Y()); + float tz0 = (float) (b.P0().Z() - a.P0().Z())/(a.P1().Z() - a.P0().Z()); + if ( (fabs(tx0 - ty0) < eps ) && (fabs(ty0 - tz0) < eps) ) { + float tx1 = (float) (b.P1().X() - a.P0().X())/(a.P1().X() - a.P0().X()); + float ty1 = (float) (b.P1().Y() - a.P0().Y())/(a.P1().Y() - a.P0().Y()); + float tz1 = (float) (b.P1().Z() - a.P0().Z())/(a.P1().Z() - a.P0().Z()); + if ( (fabs(tx1 - ty1) < eps) && (fabs(ty1 - tz1) < eps) ) { + if ( (tx0 < 0.0f + eps) && (tx1 > 1.0f - eps) ) return false; //external + if ( (tx1 < 0.0f + eps) && (tx0 > 1.0f - eps) ) return false; //external + if ( (tx0 >= 0.0f + eps) && (tx0 <= 1.0f - eps) && (tx1 >= 0.0f + eps) && (tx1 <= 1.0f - eps) ) { + //double split + if ( tx0 > tx1 ) { c.Flip(); v = std::make_pair( v.second, v.first ); } + //insert new edges + trash[i].edges.insert( trash[i].edges.begin() + j + 1, vcg::Segment3( trash[i].edges[j].P0(), c.P0() ) ); + trash[i].verts.insert( trash[i].verts.begin() + j + 1, std::make_pair( trash[i].verts[j].first, v.first ) ); + trash[i].edges.insert( trash[i].edges.begin() + j + 2, c ); + trash[i].verts.insert( trash[i].verts.begin() + j + 2, v ); + trash[i].edges.insert( trash[i].edges.begin() + j + 3, vcg::Segment3( c.P1(), trash[i].edges[j].P1() ) ); + trash[i].verts.insert( trash[i].verts.begin() + j + 3, std::make_pair( v.second, trash[i].verts[j].second ) ); + //erase old one + trash[i].edges.erase( trash[i].edges.begin() + j ); + trash[i].verts.erase( trash[i].verts.begin() + j ); + return true; + } + if ( (tx0 >= 0.0f + eps) && (tx0 <= 1.0f - eps) ) { + //single split in P0 + //insert new edges + trash[i].edges.insert( trash[i].edges.begin() + j + 1, vcg::Segment3( trash[i].edges[j].P0(), c.P0() ) ); + trash[i].verts.insert( trash[i].verts.begin() + j + 1, std::make_pair( trash[i].verts[j].first, v.first ) ); + trash[i].edges.insert( trash[i].edges.begin() + j + 2, vcg::Segment3( c.P0(), trash[i].edges[j].P1() ) ); + trash[i].verts.insert( trash[i].verts.begin() + j + 2, std::make_pair( v.first, trash[i].verts[j].second ) ); + //erase old one + trash[i].edges.erase( trash[i].edges.begin() + j ); + trash[i].verts.erase( trash[i].verts.begin() + j ); + return true; + } + if ( (tx1 >= 0.0f + eps) && (tx1 <= 1.0f - eps) ) { + //single split in P1 + //insert new edges + trash[i].edges.insert( trash[i].edges.begin() + j + 1, vcg::Segment3( trash[i].edges[j].P0(), c.P1() ) ); + trash[i].verts.insert( trash[i].verts.begin() + j + 1, std::make_pair( trash[i].verts[j].first, v.second ) ); + trash[i].edges.insert( trash[i].edges.begin() + j + 2, vcg::Segment3( c.P1(), trash[i].edges[j].P1() ) ); + trash[i].verts.insert( trash[i].verts.begin() + j + 2, std::make_pair( v.second, trash[i].verts[j].second ) ); + //erase old one + trash[i].edges.erase( trash[i].edges.begin() + j ); + trash[i].verts.erase( trash[i].verts.begin() + j ); + return true; + } + } + } + } + } + + bool found = false; + for ( unsigned int j = 0; j < border.size(); j ++ ) { + for ( size_t i = 0; i < border[j].verts.size() && !found; i ++ ) { + if ( border[j].verts[i].first == v.second ) { found = true; border[j].edges.insert( border[j].edges.begin() + i, c ); border[j].verts.insert( border[j].verts.begin() + i, v ); } //insert before i-th element + else + if ( border[j].verts[i].second == v.first ) { found = true; border[j].edges.insert( border[j].edges.begin() + i + 1, c ); border[j].verts.insert( border[j].verts.begin() + i + 1, v ); } //insert after i-th element + } + } + if (!found) { //Create a new polyline ad add it to the border list + polyline nwpoly; nwpoly.edges.push_back( c ); nwpoly.verts.push_back( v ); border.push_back( nwpoly ); + } else { + //Merge consecutive polylines into a single one + for ( unsigned int j = 0; j < border.size(); j ++ ) + for ( unsigned int i = j+1; i < border.size(); i ++ ) { + if ( border[j].verts.front().first == border[i].verts.back().second ) { + border[j].edges.insert( border[j].edges.begin(), border[i].edges.begin(), border[i].edges.end() ); + border[j].verts.insert( border[j].verts.begin(), border[i].verts.begin(), border[i].verts.end() ); + border.erase(border.begin() + i); + } + else if ( border[j].verts.back().second == border[i].verts.front().first ) { + border[j].edges.insert( border[j].edges.end(), border[i].edges.begin(), border[i].edges.end() ); + border[j].verts.insert( border[j].verts.end(), border[i].verts.begin(), border[i].verts.end() ); + border.erase(border.begin() + i); + } + } + }//end if (!found) + + for ( size_t k = 0; k < border.size(); k ++) { + for ( size_t i = 0; i < trash.size(); i ++ ) { + for ( size_t j = 0; j < trash[i].verts.size(); j ++ ) { //Only one trash component + if ( trash[i].verts[j].first == border[k].verts.back().second ) { + trash[i].edges[j].P0() = border[k].edges.back().P1(); + //trash[i].edges.erase( trash[i].edges.begin() + j ); trash[i].verts.erase( trash[i].verts.begin() + j ) + } else + if ( trash[i].verts[j].first == border[k].verts.front().first ) { + trash[i].edges[j].P0() = border[k].edges.front().P0(); + //trash[i].edges.erase( trash[i].edges.begin() + j ); trash[i].verts.erase( trash[i].verts.begin() + j ) + } + } + } + } + return true; + }//end AddToBorder + + // Add c.component + virtual void AddCComponent( polyline c ) { + conn.push_back(c); + } + // Add t.component + virtual void AddTComponent( polyline t ) { + trash.push_back(t); + } + + //Set eps + virtual void SetEps( float e ) { + eps = e; + } + + // Set initial t.component + virtual void Init( CMeshO::FaceType f, int a, int b, int c ) { + if (!trash.empty()) return; + polyline tri; tri.edges.push_back( vcg::Segment3(f.P(0), f.P(1)) ); + tri.edges.push_back( vcg::Segment3(f.P(1), f.P(2)) ); + tri.edges.push_back( vcg::Segment3(f.P(2), f.P(0)) ); + tri.verts.push_back( std::make_pair(a, b) ); + tri.verts.push_back( std::make_pair(b, c) ); + tri.verts.push_back( std::make_pair(c, a) ); + AddTComponent( tri ); + } + + // Remove c.component + virtual void RemoveCComponent( int i ) { + conn.erase( conn.begin() + i ); + } + + // Remove t.component + virtual void RemoveTComponent( int i ) { + trash.erase( trash.begin() + i ); + } + + // Number of c.component + virtual int nCComponent( ) { + return conn.size(); + } + + // Number of t.component + virtual size_t nTComponent( ) { + return trash.size(); + } + + // Add vertex in original triangle + virtual bool addVertex( CMeshO::VertexPointer v, int v_index ) { + int cnt = 0; int split = -1; + for ( size_t i = 0; i < trash.size(); i ++ ) { //one component only + + for ( size_t j = 0; j < trash[i].verts.size(); j ++ ) { //search for closest edge + if ( trash[i].verts[j].first == v_index ) return false; + } + + for ( size_t j = 0; j < trash[i].edges.size(); j ++ ) { //search for closest edge + //if ( vcg::SquaredDistance( trash[i].edges[j], v->P() ) <= eps ) { + float dist; + vcg::Point3f clos; + vcg::SegmentPointSquaredDistance(trash[i].edges[j], v->P(),clos,dist); + if (dist <= eps ) { + cnt++; split = j; + } + } + } + if (!cnt) return false; + + if ( cnt == 1 ) { //one edge only -> split the edge + vcg::Segment3 splitting_edge = trash[0].edges[split]; + std::pair splitting_edge_v = trash[0].verts[split]; + + trash[0].edges.erase(trash[0].edges.begin()+split); //remove edge + trash[0].verts.erase(trash[0].verts.begin()+split); + //replace edge using two new edges + trash[0].edges.insert(trash[0].edges.begin()+split, vcg::Segment3( splitting_edge.P0(), v->P() ) ); + trash[0].edges.insert(trash[0].edges.begin()+split+1, vcg::Segment3( v->P(), splitting_edge.P1() ) ); + trash[0].verts.insert(trash[0].verts.begin()+split, std::make_pair( splitting_edge_v.first, v_index ) ); + trash[0].verts.insert(trash[0].verts.begin()+split+1, std::make_pair( v_index, splitting_edge_v.second ) ); + } + + if ( cnt == 2 ) { // search for closest vertex and copy vertex coords + for ( size_t i = 0; i < trash.size(); i ++ ) { //one component only + for ( size_t j = 0; j < trash[i].edges.size(); j ++ ) { //search for closest edge + if ( vcg::Distance( trash[i].edges[j].P0(), v->P() ) <= eps ) { + v->P() = trash[i].edges[j].P0(); + } + } + } + } + + return true; + } + +};//end struct + +class compareFaceQuality { + public: + compareFaceQuality() { }; + + bool operator () (const std::pair f1, const std::pair f2) { + //quality f1 < quality f2 return true + return ( f1.first->Q() > f2.first->Q() ); + } +}; + +class FilterZippering : public QObject, public MeshFilterInterface +{ + Q_OBJECT + Q_INTERFACES(MeshFilterInterface) + + typedef vcg::GridStaticPtr MeshFaceGrid; + typedef vcg::GridStaticPtr MeshVertGrid; + +public: + //Different operations in different plugins + enum { FP_REDUNDANCY, + FP_ZIPPERING }; + + FilterZippering(); + + virtual QString filterName(FilterIDType filter) const; + virtual QString filterInfo(FilterIDType filter) const; + virtual void initParameterSet(QAction *,MeshDocument &/*m*/, RichParameterSet & /*parent*/); + int getRequirements(QAction *action); + virtual bool applyFilter(QAction *filter, MeshDocument &md, RichParameterSet & /*parent*/, vcg::CallBackPos * cb) ; + FilterClass getClass(QAction *a); + virtual int postCondition( QAction */*a*/ ) const { return MeshModel::MM_FACEFACETOPO|MeshModel::MM_VERTNORMAL; } + + + +private: + template + ScalarType SquaredDistance( vcg::Segment3 &s, vcg::Point3 &p); + + template + vcg::Point3 ClosestPoint( vcg::Segment3 &s, vcg::Point3 &p); + + bool checkRedundancy( CMeshO::FacePointer f, //face + MeshModel *a, //mesh A + MeshFaceGrid &grid, //grid A + CMeshO::ScalarType max_dist ); //Max search distance + bool simpleCheckRedundancy( CMeshO::FacePointer f, //face + MeshModel *a, //mesh A + MeshFaceGrid &grid, //grid A + CMeshO::ScalarType max_dist,//Max search distance + bool test ); + bool isBorderVert( CMeshO::FacePointer f, int i); + bool isOnBorder( CMeshO::CoordType point, CMeshO::FacePointer f ); + bool isOnEdge( CMeshO::CoordType point, CMeshO::FacePointer f ); + bool isAdjacent( CMeshO::FacePointer f1, CMeshO::FacePointer f2 ); + int sharesVertex( CMeshO::FacePointer f1, CMeshO::FacePointer f2 ); + void handleBorder( aux_info &info, //Auxiliar information for triangulatio + vcg::Point3 N, //face normal (useful for proiection) + std::vector &coords, //output coords + std::vector &output ); //output v. pointers + polyline cutComponent( polyline comp, //Component to be cut + polyline border, //border + vcg::Matrix44 rot_mat ); //Rotation matrix + int searchComponent( aux_info &info, //Auxiliar info + vcg::Point3 P0, //Start border point + vcg::Point3 P1, //End border point + bool &conn ); + bool findIntersection( CMeshO::FacePointer currentF, //face + vcg::Segment3 edge, //edge + int last_split, //last splitted edge + int &splitted_edge, //currently splitted edge + vcg::Point3 &hit ); //approximate intersection point + + + + //init unsorted queue + bool Init_q( std::vector< std::pair >& queue, //the queue + MeshModel* a, //mesh A + MeshModel* b, //mesh B + bool fullProcess ); //fullProcess flag + //init priority queue (overload) + bool Init_pq( std::priority_queue< std::pair, std::vector< std::pair >, compareFaceQuality >& queue, //the queue + MeshModel* a, //mesh A + MeshModel* b, //mesh B + bool fullProcess ); + //select redundant face (normal method, unsorted queue) + int selectRedundant( std::vector< std::pair >& queue, //queue + MeshModel* a, //mesh A + MeshModel* b, //mesh B + float epsilon ); //max search distance + //select redundant face (normal method, priority queue) + int selectRedundant_pq( std::priority_queue< std::pair, std::vector< std::pair >, compareFaceQuality >& queue, //the queue + MeshModel* a, //mesh A + MeshModel* b, //mesh B + float epsilon ); //max search distance + //refine border of a mesh, splitting faces having two border edges + int refineBorder( MeshModel* m ); + //project face of B on the surface of A + void projectFace( CMeshO::FacePointer f, //pointer to the face that will be projected + MeshModel* a, //mesh A + MeshFaceGrid grid_a, //grid on A + float max_dist, //max dist search + std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + std::vector< CMeshO::FacePointer >& tbt_faces, //vector to-be-triangulated faces + std::vector< CMeshO::FacePointer >& tbr_faces, //vector to-be-removed faces + std::vector< int >& verts ); //vector of indices + + //case 01: vertices of border edge project on the same face + void handleBorderEdgeSF ( std::pair< int, int >& current_edge, //current border edge + MeshModel* a, //mesh A + CMeshO::FacePointer startF, //face where first vertex lies + CMeshO::FacePointer endF, //face where second vertex lies + CMeshO::FacePointer splittingF, //splitting face + std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + std::vector < std::pair< int, int > >& stack, //stack containing border edges + std::vector< int >& verts ); //vector of indices + + //case 02: vertices of border edge project on adjacent faces + void handleBorderEdgeAF ( std::pair< int, int >& current_edge, //current border edge + MeshModel* a, //mesh A + CMeshO::FacePointer startF, //face where first vertex lies + CMeshO::FacePointer endF, //face where second vertex lies + CMeshO::FacePointer splittingF, //splitting face + std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + std::vector < std::pair< int, int > >& stack, //stack containing border edges + std::vector< int >& verts ); //vector of indices + + //case 03: vertices of border edge project on non-adjacent faces + void handleBorderEdgeNF ( std::pair< int, int >& current_edge, //current border edge + MeshModel* a, //mesh A + CMeshO::FacePointer startF, //face where first vertex lies + CMeshO::FacePointer endF, //face where second vertex lies + CMeshO::FacePointer splittingF, //splitting face + std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + std::vector < std::pair< int, int > >& stack, //stack containing border edges + std::vector< int >& verts ); //vector of indices + + //case 04: both vertices of current_edge project on another border edge + //return true if the whole current_edge project on border edge + bool handleBorderEdgeBB ( std::pair< int, int >& current_edge, //current border edge + MeshModel* a, //mesh A + MeshFaceGrid grid_a, //grid on A (needed for sampling) + float max_dist, //max search dist (needed for sampling) + vcg::Point3 closestStart, //closest point on startF + vcg::Point3 closestEnd, //closest point on endF + CMeshO::FacePointer startF, //face where first vertex lies + CMeshO::FacePointer endF, //face where second vertex lies + CMeshO::FacePointer splittingF, //splitting face + std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + std::vector < std::pair< int, int > >& stack, //stack containing border edges + std::vector< int >& verts ); //vector of indices + + + //case 04: both vertices of current_edge project on another border edge + //return true if the whole current_edge project on border edge + bool handleBorderEdgeBB ( std::pair< int, int >& current_edge, //current border edge + MeshModel* a, //mesh A + MeshFaceGrid grid_a, //grid on A (needed for sampling) + float max_dist, //max search dist (needed for sampling) + CMeshO::FacePointer startF, //face where first vertex lies + CMeshO::FacePointer endF, //face where second vertex lies + CMeshO::FacePointer splittingF, //splitting face + std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + std::vector < std::pair< int, int > >& stack, //stack containing border edges + std::vector< int >& verts ); //vector of indices + + //case 05: one of the vertices doesn't project on the surface of the mesh + void handleBorderEdgeOB ( std::pair< int, int >& current_edge, //current border edge + int direction, //splitting direction (1 from start to end, 0 from end to start) + MeshModel* a, //mesh A + MeshFaceGrid grid_a, //grid on A (needed for sampling) + float max_dist, //max search dist (needed for sampling) + CMeshO::FacePointer startF, //face where first vertex lies + CMeshO::FacePointer endF, //face where second vertex lies + CMeshO::FacePointer splittingF, //splitting face + std::map< CMeshO::FacePointer, aux_info >& map_info, //map with auxiliar information + std::vector < std::pair< int, int > >& stack, //stack containing border edges + std::vector< CMeshO::FacePointer >& tbt_faces, //stack containing pointers to face that wille be retriangulated + std::vector< int >& verts ); //vector of indices + + int preProcess ( std::vector< std::pair >& queue, //queue + MeshModel* a, + MeshModel* b, + MeshFaceGrid grid_a, //grid on A + MeshFaceGrid grid_b, //grid on A + float max_dist ); //max dist search + + int preProcess_pq ( std::priority_queue< std::pair, std::vector< std::pair >, compareFaceQuality >& queue, //the queue + MeshModel* a, + MeshModel* b, + MeshFaceGrid grid_a, //grid on A + MeshFaceGrid grid_b, //grid on A + float max_dist ); //max dist search + + + float eps; + + int dbg_cnt; +}; + +#endif