From f2ce7bdf696e1731a28f7934f45ebf82beca8554 Mon Sep 17 00:00:00 2001
From: nico <nico.pietroni@gmail.com>
Date: Fri, 15 May 2020 10:02:30 +1000
Subject: [PATCH] added one more parameter in SmoothPCA to ignore some faces in
 the computation

---
 vcg/complex/algorithms/polygonal_algorithms.h | 15 +++++++++++++--
 1 file changed, 13 insertions(+), 2 deletions(-)

diff --git a/vcg/complex/algorithms/polygonal_algorithms.h b/vcg/complex/algorithms/polygonal_algorithms.h
index 496f21b1..bda99570 100644
--- a/vcg/complex/algorithms/polygonal_algorithms.h
+++ b/vcg/complex/algorithms/polygonal_algorithms.h
@@ -507,7 +507,8 @@ public:
                           bool isotropic=true,
                           ScalarType smoothTerm=0.1,
                           bool fixB=true,
-                          bool WeightByQuality=false)
+                          bool WeightByQuality=false,
+                          const std::vector<bool> *IgnoreF=NULL)
     {
         (void)isotropic;
         typedef typename PolyMeshType::FaceType PolygonType;
@@ -524,6 +525,7 @@ public:
         if (WeightByQuality)
             UpdateQuality(poly_m,QTemplate);
 
+        if (IgnoreF!=NULL){assert((*IgnoreF).size()==poly_m.face.size());}
         for (size_t s=0;s<(size_t)relax_step;s++)
         {
             //initialize the accumulation vector
@@ -533,6 +535,7 @@ public:
 
             for (size_t i=0;i<poly_m.face.size();i++)
             {
+                if ((IgnoreF!=NULL)&&((*IgnoreF)[i]))continue;
                 std::vector<typename PolygonType::CoordType> TemplatePos;
                 GetRotatedTemplatePos(poly_m.face[i],TemplatePos);
                 //then cumulate the position per vertex
@@ -566,7 +569,15 @@ public:
                 //               if (alpha<0)alpha=0;
                 //               if (alpha>1)alpha=1;
                 //               if (isnan(alpha))alpha=1;
-                CoordType newP=avgPos[i]/weightSum[i];
+
+                CoordType newP=poly_m.vert[i].P();
+                //safety checks
+                if (weightSum[i]>0)
+                    newP=avgPos[i]/weightSum[i];
+                if (isnan(newP.X())||isnan(newP.Y())||isnan(newP.Z()))
+                     newP=poly_m.vert[i].P();
+                if ((newP-poly_m.vert[i].P()).Norm()>poly_m.bbox.Diag())
+                    newP=poly_m.vert[i].P();
                 //std::cout<<"W "<<weightSum[i]<<std::endl;
                 newP=newP*(1-alpha)+AvVert[i]*alpha;
                 //newP=AvVert[i];