vcglib: refine.h Source File

00001 /***********F*****************************************************************
00002 * VCGLib                                                            o o     *
00003 * Visual and Computer Graphics Library                            o     o   *
00004 *                                                                _   O  _   *
00005 * Copyright(C) 2004                                                \/)\/    *
00006 * Visual Computing Lab                                            /\/|      *
00007 * ISTI - Italian National Research Council                           |      *
00008 *                                                                    \      *
00009 * All rights reserved.                                                      *
00010 *                                                                           *
00011 * This program is free software; you can redistribute it and/or modify      *   
00012 * it under the terms of the GNU General Public License as published by      *
00013 * the Free Software Foundation; either version 2 of the License, or         *
00014 * (at your option) any later version.                                       *
00015 *                                                                           *
00016 * This program is distributed in the hope that it will be useful,           *
00017 * but WITHOUT ANY WARRANTY; without even the implied warranty of            *
00018 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the             *
00019 * GNU General Public License (http://www.gnu.org/licenses/gpl.txt)          *
00020 * for more details.                                                         *
00021 *                                                                           *
00022 ****************************************************************************/
00023 
00024 #ifndef __VCGLIB_REFINE
00025 #define __VCGLIB_REFINE
00026 
00027 #include <functional>
00028 #include <map>
00029 #include <vector>
00030 #include <vcg/space/sphere3.h>
00031 #include <vcg/space/plane3.h>
00032 #include <vcg/space/texcoord2.h>
00033 #include <vcg/space/color4.h>
00034 #include <vcg/simplex/face/pos.h>
00035 #include<vcg/complex/trimesh/allocate.h>
00036 #include<vcg/complex/trimesh/update/topology.h>
00037 #include<wrap/callback.h>
00038 #include <vcg/complex/trimesh/base.h>
00039 #include <vcg/space/triangle3.h>
00040 
00041 namespace vcg{
00042         
00043 /* A very short intro about the generic refinement framework,
00044         the main fuction is the 
00045         
00046  template<class MESH_TYPE,class MIDPOINT, class EDGEPRED>
00047  bool RefineE(MESH_TYPE &m, MIDPOINT mid, EDGEPRED ep,bool RefineSelected=false, CallBackPos *cb = 0)
00048  
00049  You have to provide two functor objects to this, one for deciding what edge has to be spltted and one to decide position and new values for the attributes of the new point.
00050         
00051  for example the minimal EDGEPRED is
00052  
00053  template <class MESH_TYPE, class FLT> class EdgeLen
00054  {
00055    public: 
00056          FLT thr2;
00057          bool operator()(face::Pos<typename MESH_TYPE::FaceType> ep) const
00058          {
00059                         return SquaredDistance(ep.f->V(ep.z)->P(), ep.f->V1(ep.z)->P())>thr2;
00060          }
00061  };
00062   
00063  With a bit of patience you can customize to make also slicing operation.
00064  
00065 */
00066         
00067 
00068 /* The table which encodes how to subdivide a triangle depending 
00069    on the splitted edges is organized as such:
00070 
00071     TriNum (the first number):    encodes the number of triangles
00072     TV (the following 4 triples): encodes the resulting triangles where
00073           0, 1, 2 are the original vertices of the triangles and 3, 4, 5 
00074           (mp01, mp12, mp20) are the midpoints of the three edges.
00075 
00076    In the case two edges are splitted the triangle has 2 possible splittings:
00077 we need to choose a diagonal of the resulting trapezoid.
00078 'swap' encodes the two diagonals to test: if diag1 < diag2 we swap the diagonal
00079 like this (140, 504 -> 150, 514) (the second vertex of each triangles is replaced 
00080  by the first vertex of the other one).
00081             2
00082            / \
00083           5---4
00084          /     \
00085         0-------1
00086 
00087 */
00088 
00089 class Split {
00090 public:
00091         int TriNum;                     // number of triangles
00092         int TV[4][3];   // The triangles coded as the following convention 
00093                                                                         //     0..2 vertici originali del triangolo 
00094                                                                         //     3..5 mp01, mp12, mp20 midpoints of the three edges
00095         int swap[2][2]; // the two diagonals to test for swapping
00096         int TE[4][3];   // the edge-edge correspondence between refined triangles and the old one
00097                                                                         //      (3) means the edge of the new triangle is internal;
00098 };
00099 
00100 const Split SplitTab[8]={
00101 /* m20 m12 m01 */
00102 /*  0   0   0 */        {1, {{0,1,2},{0,0,0},{0,0,0},{0,0,0}}, {{0,0},{0,0}},  {{0,1,2},{0,0,0},{0,0,0},{0,0,0}} },
00103 /*  0   0   1 */        {2, {{0,3,2},{3,1,2},{0,0,0},{0,0,0}}, {{0,0},{0,0}},  {{0,3,2},{0,1,3},{0,0,0},{0,0,0}} },
00104 /*  0   1   0 */        {2, {{0,1,4},{0,4,2},{0,0,0},{0,0,0}}, {{0,0},{0,0}},  {{0,1,3},{3,1,2},{0,0,0},{0,0,0}} },
00105 /*  0   1   1 */        {3, {{3,1,4},{0,3,2},{4,2,3},{0,0,0}}, {{0,4},{3,2}},  {{0,1,3},{0,3,2},{1,3,3},{0,0,0}} },
00106 /*  1   0   0 */        {2, {{0,1,5},{5,1,2},{0,0,0},{0,0,0}}, {{0,0},{0,0}},  {{0,3,2},{3,1,2},{0,0,0},{0,0,0}} },
00107 /*  1   0   1 */        {3, {{0,3,5},{3,1,5},{2,5,1},{0,0,0}}, {{3,2},{5,1}},  {{0,3,2},{0,3,3},{2,3,1},{0,0,0}} },
00108 /*  1   1   0 */        {3, {{2,5,4},{0,1,5},{4,5,1},{0,0,0}}, {{0,4},{5,1}},  {{2,3,1},{0,3,2},{3,3,1},{0,0,0}} },
00109 /*  1   1   1 */        //{4, {{0,3,5},{3,1,4},{5,4,2},{3,4,5}}, {{0,0},{0,0}},  {{0,3,2},{0,1,3},{3,1,2},{3,3,3}} },
00110 /*  1   1   1 */        {4, {{3,4,5},{0,3,5},{3,1,4},{5,4,2}}, {{0,0},{0,0}},  {{3,3,3},{0,3,2},{0,1,3},{3,1,2}} },
00111 };
00112 
00113 // Basic subdivision class
00114 // This class must provide methods for finding the position of the newly created vertices
00115 // In this implemenation we simply put the new vertex in the MidPoint position.
00116 // Color and TexCoords are interpolated accordingly.
00117 template<class MESH_TYPE>
00118 struct MidPoint : public   std::unary_function<face::Pos<typename MESH_TYPE::FaceType> ,  typename MESH_TYPE::CoordType >
00119 {
00120          MidPoint(MESH_TYPE *_mp) { mp=_mp; }
00121          
00122          MESH_TYPE *mp;
00123          
00124         void operator()(typename MESH_TYPE::VertexType &nv, face::Pos<typename MESH_TYPE::FaceType>  ep){
00125                 assert(mp);
00126                 nv.P()=   (ep.f->V(ep.z)->P()+ep.f->V1(ep.z)->P())/2.0;
00127 
00128                 if( MESH_TYPE::HasPerVertexNormal())
00129                         nv.N()= (ep.f->V(ep.z)->N()+ep.f->V1(ep.z)->N()).normalized();
00130 
00131                 if( MESH_TYPE::HasPerVertexColor())
00132                         nv.C().lerp(ep.f->V(ep.z)->C(),ep.f->V1(ep.z)->C(),.5f);
00133                 
00134                 if( MESH_TYPE::HasPerVertexQuality())
00135                         nv.Q() = ((ep.f->V(ep.z)->Q()+ep.f->V1(ep.z)->Q())) / 2.0;
00136 
00137                 if( tri::HasPerVertexTexCoord(*mp))
00138                         nv.T().P() = ((ep.f->V(ep.z)->T().P()+ep.f->V1(ep.z)->T().P())) / 2.0;
00139         }
00140 
00141         Color4<typename MESH_TYPE::ScalarType> WedgeInterp(Color4<typename MESH_TYPE::ScalarType> &c0, Color4<typename MESH_TYPE::ScalarType> &c1)
00142         {
00143                 Color4<typename MESH_TYPE::ScalarType> cc;
00144                 return cc.lerp(c0,c1,0.5f);
00145         }
00146 
00147         template<class FL_TYPE>
00148         TexCoord2<FL_TYPE,1> WedgeInterp(TexCoord2<FL_TYPE,1> &t0, TexCoord2<FL_TYPE,1> &t1)
00149         {
00150                 TexCoord2<FL_TYPE,1> tmp;
00151                 assert(t0.n()== t1.n());
00152                 tmp.n()=t0.n(); 
00153                 tmp.t()=(t0.t()+t1.t())/2.0;
00154                 return tmp;
00155         }
00156 };
00157 
00158 
00159 
00160 template<class MESH_TYPE>
00161 struct MidPointArc : public std::unary_function<face::Pos<typename MESH_TYPE::FaceType> ,  typename MESH_TYPE::CoordType>
00162 {
00163         void operator()(typename MESH_TYPE::VertexType &nv, face::Pos<typename MESH_TYPE::FaceType> ep)
00164         {
00165                 const typename MESH_TYPE::ScalarType EPS =1e-10;
00166                 typename MESH_TYPE::CoordType vp = (ep.f->V(ep.z)->P()+ep.f->V1(ep.z)->P())/2.0;
00167                 typename MESH_TYPE::CoordType  n = (ep.f->V(ep.z)->N()+ep.f->V1(ep.z)->N())/2.0;
00168                 typename MESH_TYPE::ScalarType w =n.Norm();
00169                 if(w<EPS) { nv.P()=(ep.f->V(ep.z)->P()+ep.f->V1(ep.z)->P())/2.0; return;}
00170                 n/=w;
00171                 typename MESH_TYPE::CoordType d0 = ep.f->V(ep.z)->P() - vp;
00172                 typename MESH_TYPE::CoordType d1 = ep.f->V1(ep.z)->P()- vp;
00173                 typename MESH_TYPE::ScalarType d=Distance(ep.f->V(ep.z)->P(),ep.f->V1(ep.z)->P())/2.0;
00174 
00175                 typename MESH_TYPE::CoordType  nn = ep.f->V1(ep.z)->N() ^ ep.f->V(ep.z)->N();
00176                 typename MESH_TYPE::CoordType  np = n ^ d0; //vector perpendicular to the edge plane, normal is interpolated
00177                 np.Normalize();
00178                 double sign=1;
00179                 if(np*nn<0) sign=-1; // se le normali non divergono sposta il punto nella direzione opposta
00180     
00181                 typename MESH_TYPE::CoordType n0=ep.f->V(ep.z)->N() -np*(ep.f->V(ep.z)->N()*np);
00182                 n0.Normalize();
00183                 typename MESH_TYPE::CoordType n1=ep.f->V1(ep.z)->N()-np*(ep.f->V1(ep.z)->N()*np);
00184                 assert(n1.Norm()>EPS);
00185                 n1.Normalize();
00186                 typename MESH_TYPE::ScalarType cosa0=n0*n;
00187                 typename MESH_TYPE::ScalarType cosa1=n1*n;
00188                 if(2-cosa0-cosa1<EPS) {nv.P()=(ep.f->V(ep.z)->P()+ep.f->V1(ep.z)->P())/2.0;return;}
00189                 typename MESH_TYPE::ScalarType cosb0=(d0*n)/d;
00190                 typename MESH_TYPE::ScalarType cosb1=(d1*n)/d;
00191                 assert(1+cosa0>EPS);
00192                 assert(1+cosa1>EPS);
00193                 typename MESH_TYPE::ScalarType delta0=d*(cosb0 +sqrt( ((1-cosb0*cosb0)*(1-cosa0))/(1+cosa0)) );
00194                 typename MESH_TYPE::ScalarType delta1=d*(cosb1 +sqrt( ((1-cosb1*cosb1)*(1-cosa1))/(1+cosa1)) );
00195                 assert(delta0+delta1<2*d);
00196                 nv.P()=vp+n*sign*(delta0+delta1)/2.0;
00197                 return ;
00198         }
00199 
00200         // Aggiunte in modo grezzo le due wedgeinterp
00201         Color4<typename MESH_TYPE::ScalarType> WedgeInterp(Color4<typename MESH_TYPE::ScalarType> &c0, Color4<typename MESH_TYPE::ScalarType> &c1)
00202         {
00203                 Color4<typename MESH_TYPE::ScalarType> cc;
00204                 return cc.lerp(c0,c1,0.5f);
00205         }
00206 
00207         template<class FL_TYPE>
00208         TexCoord2<FL_TYPE,1> WedgeInterp(TexCoord2<FL_TYPE,1> &t0, TexCoord2<FL_TYPE,1> &t1)
00209         {
00210                 TexCoord2<FL_TYPE,1> tmp;
00211                 assert(t0.n()== t1.n());
00212                 tmp.n()=t0.n(); 
00213                 tmp.t()=(t0.t()+t1.t())/2.0;
00214                 return tmp;
00215         }
00216 
00217 };
00218 
00219 /*
00220 Versione Della Midpoint basata sul paper:
00221 S. Karbacher, S. Seeger, G. Hausler
00222 A non linear subdivision scheme for triangle meshes
00223 
00224         Non funziona!
00225         Almeno due problemi:
00226         1) il verso delle normali influenza il risultato (e.g. si funziona solo se le normali divergono)
00227                  Risolvibile controllando se le normali divergono
00228   2) gli archi vanno calcolati sul piano definito dalla normale interpolata e l'edge.
00229                  funziona molto meglio nelle zone di sella e non semplici.
00230 
00231 */
00232 template<class MESH_TYPE>
00233 struct MidPointArcNaive : public std::unary_function< face::Pos<typename MESH_TYPE::FaceType> , typename MESH_TYPE::CoordType>
00234 {
00235         typename MESH_TYPE::CoordType operator()(face::Pos<typename MESH_TYPE::FaceType>  ep)
00236         {
00237 
00238                 typename MESH_TYPE::CoordType vp = (ep.f->V(ep.z)->P()+ep.f->V1(ep.z)->P())/2.0;
00239                 typename MESH_TYPE::CoordType  n = (ep.f->V(ep.z)->N()+ep.f->V1(ep.z)->N())/2.0;
00240                 n.Normalize();
00241                 typename MESH_TYPE::CoordType d0 = ep.f->V(ep.z)->P() - vp;
00242                 typename MESH_TYPE::CoordType d1 = ep.f->V1(ep.z)->P()- vp;
00243                 typename MESH_TYPE::ScalarType d=Distance(ep.f->V(ep.z)->P(),ep.f->V1(ep.z)->P())/2.0;
00244 
00245                 typename MESH_TYPE::ScalarType cosa0=ep.f->V(ep.z)->N()*n;
00246                 typename MESH_TYPE::ScalarType cosa1=ep.f->V1(ep.z)->N()*n;
00247                 typename MESH_TYPE::ScalarType cosb0=(d0*n)/d;
00248                 typename MESH_TYPE::ScalarType cosb1=(d1*n)/d;
00249 
00250                 typename MESH_TYPE::ScalarType delta0=d*(cosb0 +sqrt( ((1-cosb0*cosb0)*(1-cosa0))/(1+cosa0)) );
00251                 typename MESH_TYPE::ScalarType delta1=d*(cosb1 +sqrt( ((1-cosb1*cosb1)*(1-cosa1))/(1+cosa1)) );
00252                 
00253                 return vp+n*(delta0+delta1)/2.0;
00254         }
00255 };
00256 
00257 
00258 // Basic Predicate that tells if a given edge must be splitted.
00259 // the constructure requires the threshold. 
00260 // VERY IMPORTANT REQUIREMENT: this function must be symmetric
00261 // e.g. it must return the same value if the Pos is VFlipped.
00262 // If this function is not symmetric the Refine can crash.
00263 
00264 template <class MESH_TYPE, class FLT>
00265 class EdgeLen
00266 {
00267     FLT squaredThr;
00268 public:
00269         EdgeLen(){}; 
00270         EdgeLen(FLT threshold) {setThr(threshold);}
00271         void setThr(FLT threshold) {squaredThr = threshold*threshold; }
00272         bool operator()(face::Pos<typename MESH_TYPE::FaceType> ep) const
00273         {
00274                 return SquaredDistance(ep.V()->P(), ep.VFlip()->P())>squaredThr;
00275         }
00276 };
00277 
00278 /*********************************************************/
00279 /*********************************************************
00280 
00281 Given a mesh the following function refines it according to two functor objects:
00282 
00283 - a predicate that tells if a given edge must be splitted
00284 
00285 - a functor that gives you the new poistion of the created vertices (starting from an edge)
00286 
00287 If RefineSelected is true only selected faces are taken into account for being splitted.
00288 
00289 Requirement: FF Adjacency and Manifoldness
00290 
00291 **********************************************************/
00292 /*********************************************************/
00293 template <class VertexPointer> 
00294 class RefinedFaceData   
00295         {
00296                 public:
00297                 RefinedFaceData(){
00298                         ep[0]=0;ep[1]=0;ep[2]=0;
00299                         vp[0]=0;vp[1]=0;vp[2]=0;                
00300                 }
00301                 bool ep[3];
00302                 VertexPointer vp[3];
00303         };
00304 
00305 template<class MESH_TYPE,class MIDPOINT, class EDGEPRED>
00306 bool RefineE(MESH_TYPE &m, MIDPOINT mid, EDGEPRED ep,bool RefineSelected=false, CallBackPos *cb = 0)
00307 {
00308         // common typenames
00309         typedef typename MESH_TYPE::VertexIterator VertexIterator;
00310         typedef typename MESH_TYPE::FaceIterator FaceIterator;
00311         typedef typename MESH_TYPE::VertexPointer VertexPointer;
00312         typedef typename MESH_TYPE::FacePointer FacePointer;
00313         typedef typename MESH_TYPE::FaceType FaceType;  
00314         typedef typename MESH_TYPE::FaceType::TexCoordType TexCoordType;
00315         
00316         typedef face::Pos<FaceType>  PosType;
00317 
00318         int j,NewVertNum=0,NewFaceNum=0;
00319 
00320         typedef RefinedFaceData<VertexPointer> RFD;
00321         typedef typename MESH_TYPE :: template PerFaceAttributeHandle<RFD> HandleType;
00322         HandleType RD  = tri::Allocator<MESH_TYPE>:: template AddPerFaceAttribute<RFD> (m,std::string("RefineData"));
00323 
00324         // Callback stuff
00325         int step=0;
00326         int PercStep=std::max(1,m.fn/33);
00327         
00328         // First Loop: We analyze the mesh to compute the number of the new faces and new vertices 
00329         FaceIterator fi;
00330   for(fi=m.face.begin(),j=0;fi!=m.face.end();++fi) if(!(*fi).IsD())
00331         {
00332                 if(cb && (++step%PercStep)==0) (*cb)(step/PercStep,"Refining...");
00333                 // skip unselected faces if necessary
00334                 if(RefineSelected && !(*fi).IsS()) continue;
00335                 
00336                 for(j=0;j<3;j++)
00337                         {
00338                                 if(RD[fi].ep[j]) continue;
00339                                 
00340                                 PosType edgeCur(&*fi,j);
00341                                 if(RefineSelected && ! edgeCur.FFlip()->IsS()) continue;
00342                                 if(!ep(edgeCur)) continue;
00343                                                                 
00344                                 RD[edgeCur.F()].ep[edgeCur.E()]=true;
00345                                 ++NewFaceNum;
00346                                 ++NewVertNum;
00347                                 assert(edgeCur.IsManifold());
00348                                 if(!edgeCur.IsBorder()) 
00349                                 {
00350                                         edgeCur.FlipF();
00351                                         edgeCur.F()->SetV();
00352                                         RD[edgeCur.F()].ep[edgeCur.E()]=true;
00353                                         ++NewFaceNum;
00354                                 }
00355                         }
00356                 
00357   } // end face loop
00358         
00359         if(NewVertNum ==0 ) 
00360                 {       
00361                         tri::Allocator<MESH_TYPE> :: template DeletePerFaceAttribute<RefinedFaceData<VertexPointer> >  (m,RD);
00362                         return false;
00363                 }
00364         VertexIterator lastv = tri::Allocator<MESH_TYPE>::AddVertices(m,NewVertNum);
00365         
00366         // Secondo loop: We initialize a edge->vertex map 
00367         
00368         for(fi=m.face.begin();fi!=m.face.end();++fi) if(!(*fi).IsD())
00369   {
00370            if(cb && (++step%PercStep)==0)(*cb)(step/PercStep,"Refining...");
00371      for(j=0;j<3;j++)
00372                  {
00373                                 // skip unselected faces if necessary
00374                                 if(RefineSelected && !(*fi).IsS()) continue;
00375                                 for(j=0;j<3;j++)
00376                                 {
00377                                         PosType edgeCur(&*fi,j);
00378                                         if(RefineSelected && ! edgeCur.FFlip()->IsS()) continue;
00379                                         
00380                                         if( RD[edgeCur.F()].ep[edgeCur.E()] &&  RD[edgeCur.F()].vp[edgeCur.E()] ==0 )
00381                                         {
00382                                                 RD[edgeCur.F()].vp[edgeCur.E()] = &*lastv;
00383                                                 mid(*lastv,edgeCur);
00384                                                 if(!edgeCur.IsBorder()) 
00385                                                 {
00386                                                         edgeCur.FlipF();
00387                                                         assert(RD[edgeCur.F()].ep[edgeCur.E()]);
00388                                                         RD[edgeCur.F()].vp[edgeCur.E()] = &*lastv;
00389                                                 }
00390                                                 ++lastv;
00391                                         }                               
00392                                 }
00393                  }
00394   }
00395         
00396         assert(lastv==m.vert.end()); // critical assert: we MUST have used all the vertex that we forecasted we need
00397         
00398         FaceIterator lastf = tri::Allocator<MESH_TYPE>::AddFaces(m,NewFaceNum);
00399         FaceIterator oldendf = lastf; 
00400         
00401 /*
00402                 v0
00403 
00404    
00405                 f0
00406 
00407                                 mp01    f3      mp02
00408              
00409                                         
00410          f1             f2
00411 
00412  v1            mp12                v2
00413 
00414 */
00415 
00416         VertexPointer vv[6];    // The six vertices that arise in the single triangle splitting 
00417                                                                                                 //     0..2 Original triangle vertices  
00418                                                                                                 //     3..5 mp01, mp12, mp20 midpoints of the three edges
00419         FacePointer nf[4];   // The (up to) four faces that are created.
00420 
00421   TexCoordType wtt[6];  // per ogni faccia sono al piu' tre i nuovi valori 
00422                                                                                                                                                                                          // di texture per wedge (uno per ogni edge) 
00423   
00424         int fca=0,fcn =0;
00425         for(fi=m.face.begin();fi!=oldendf;++fi) if(!(*fi).IsD())
00426                 {
00427                                 if(cb && (++step%PercStep)==0)(*cb)(step/PercStep,"Refining...");
00428         fcn++;
00429                                 vv[0]=(*fi).V(0);
00430                                 vv[1]=(*fi).V(1);
00431                                 vv[2]=(*fi).V(2);
00432         vv[3] = RD[fi].vp[0];
00433         vv[4] = RD[fi].vp[1];
00434         vv[5] = RD[fi].vp[2];
00435 
00436                                 int ind=((&*vv[3])?1:0)+((&*vv[4])?2:0)+((&*vv[5])?4:0);
00437                                 
00438                                 nf[0]=&*fi;
00439                                 int i;
00440                                 for(i=1;i<SplitTab[ind].TriNum;++i){
00441                                                 nf[i]=&*lastf; ++lastf; fca++;
00442                                                 if(RefineSelected || (*fi).IsS()) (*nf[i]).SetS();
00443                                                 if(tri::HasPerFaceColor(m))
00444                                                                                         nf[i]->C()=(*fi).cC();
00445                                 }
00446         
00447                                                                         
00448         if(tri::HasPerWedgeTexCoord(m))
00449                                         for(i=0;i<3;++i)        {
00450                                                 wtt[i]=(*fi).WT(i);
00451                                                 wtt[3+i]=mid.WedgeInterp((*fi).WT(i),(*fi).WT((i+1)%3));
00452                                         }
00453 
00454                                 int orgflag=    (*fi).UberFlags();
00455                                 for(i=0;i<SplitTab[ind].TriNum;++i)
00456                                         for(j=0;j<3;++j){
00457                                                 (*nf[i]).V(j)=&*vv[SplitTab[ind].TV[i][j]];
00458                                                 
00459                                                 if(tri::HasPerWedgeTexCoord(m)) //analogo ai vertici...
00460                                                                         (*nf[i]).WT(j)=wtt[SplitTab[ind].TV[i][j]];
00461 
00462                                                 assert((*nf[i]).V(j)!=0);
00463                                                 if(SplitTab[ind].TE[i][j]!=3){
00464                                                         if(orgflag & (MESH_TYPE::FaceType::BORDER0<<(SplitTab[ind].TE[i][j]))) 
00465                                                                 (*nf[i]).SetB(j); 
00466                                                         else
00467                                                                 (*nf[i]).ClearB(j);
00468                                                 } 
00469                                                 else (*nf[i]).ClearB(j);                                                
00470                                         }
00471 
00472                 if(SplitTab[ind].TriNum==3 && 
00473                                                         SquaredDistance(vv[SplitTab[ind].swap[0][0]]->P(),vv[SplitTab[ind].swap[0][1]]->P()) < 
00474                                                         SquaredDistance(vv[SplitTab[ind].swap[1][0]]->P(),vv[SplitTab[ind].swap[1][1]]->P()) )
00475                                                         { // swap the last two triangles
00476                                                                 (*nf[2]).V(1)=(*nf[1]).V(0);
00477                                                                 (*nf[1]).V(1)=(*nf[2]).V(0);
00478                                                                 if(tri::HasPerWedgeTexCoord(m)){ //swap also textures coordinates
00479                                                                         (*nf[2]).WT(1)=(*nf[1]).WT(0);
00480                                                                         (*nf[1]).WT(1)=(*nf[2]).WT(0);
00481                                                                 }
00482                                                                 
00483                                                                 if((*nf[1]).IsB(0)) (*nf[2]).SetB(1); else (*nf[2]).ClearB(1);
00484                                                                 if((*nf[2]).IsB(0)) (*nf[1]).SetB(1); else (*nf[1]).ClearB(1);
00485                                                                 (*nf[1]).ClearB(0);
00486                                                                 (*nf[2]).ClearB(0);
00487                                                         }
00488                 }
00489         
00490         assert(lastf==m.face.end());     // critical assert: we MUST have used all the faces that we forecasted we need and that we previously allocated.
00491         assert(!m.vert.empty());
00492          for(fi=m.face.begin();fi!=m.face.end();++fi) if(!(*fi).IsD()){
00493                         assert((*fi).V(0)>=&*m.vert.begin() && (*fi).V(0)<=&m.vert.back() );
00494                         assert((*fi).V(1)>=&*m.vert.begin() && (*fi).V(1)<=&m.vert.back() );
00495                         assert((*fi).V(2)>=&*m.vert.begin() && (*fi).V(2)<=&m.vert.back() );
00496          }
00497         tri::UpdateTopology<MESH_TYPE>::FaceFace(m);    
00498         
00499         tri::Allocator<MESH_TYPE> :: template DeletePerFaceAttribute<RefinedFaceData<VertexPointer> >  (m,RD);
00500 
00501         return true;
00502 }
00503 
00504 /*************************************************************************/
00505 // simple wrapper of the base refine for lazy coder that do not need a edge predicate
00506 
00507 template<class MESH_TYPE,class MIDPOINT>
00508 bool Refine(MESH_TYPE &m, MIDPOINT mid, typename MESH_TYPE::ScalarType thr=0,bool RefineSelected=false, CallBackPos *cb = 0)
00509 {
00510         EdgeLen <MESH_TYPE, typename MESH_TYPE::ScalarType> ep(thr);
00511   return RefineE(m,mid,ep,RefineSelected,cb);
00512 }
00513 /*************************************************************************/
00514 
00515 /*
00516 Modified Butterfly interpolation scheme, 
00517 as presented in 
00518 Zorin, Schroeder
00519 Subdivision for modeling and animation
00520 Siggraph 2000 Course Notes
00521 */
00522 
00523 /*
00524 
00525     vul-------vu--------vur
00526                   \      /  \      /
00527                          \    /    \    /
00528         \  /  fu  \  /
00529          vl--------vr
00530         /  \  fd  /  \
00531        /    \    /    \
00532       /      \  /      \
00533     vdl-------vd--------vdr
00534 
00535 */
00536 
00537 template<class MESH_TYPE>
00538 struct MidPointButterfly : public std::unary_function<face::Pos<typename MESH_TYPE::FaceType> , typename MESH_TYPE::CoordType>
00539 {
00540         void operator()(typename MESH_TYPE::VertexType &nv, face::Pos<typename MESH_TYPE::FaceType>  ep)
00541         {       
00542                 face::Pos<typename MESH_TYPE::FaceType> he(ep.f,ep.z,ep.f->V(ep.z));
00543                 typename MESH_TYPE::CoordType *vl,*vr;
00544                 typename MESH_TYPE::CoordType *vl0,*vr0; 
00545                 typename MESH_TYPE::CoordType *vu,*vd,*vul,*vur,*vdl,*vdr;
00546                 vl=&he.v->P();
00547                 he.FlipV();
00548                 vr=&he.v->P();
00549                 
00550                 if( MESH_TYPE::HasPerVertexColor())
00551                         nv.C().lerp(ep.f->V(ep.z)->C(),ep.f->V1(ep.z)->C(),.5f);
00552 
00553                 if(he.IsBorder())
00554                 {
00555                         he.NextB();
00556                         vr0=&he.v->P();
00557                         he.FlipV();
00558                         he.NextB();
00559                         assert(&he.v->P()==vl);
00560                         he.NextB();
00561                         vl0=&he.v->P();
00562                         nv.P()=((*vl)+(*vr))*(9.0/16.0)-((*vl0)+(*vr0))/16.0 ;
00563                 }
00564                 else 
00565                 {
00566                         he.FlipE();he.FlipV();
00567                         vu=&he.v->P();
00568                         he.FlipF();he.FlipE();he.FlipV(); 
00569                         vur=&he.v->P();
00570                         he.FlipV();he.FlipE();he.FlipF(); assert(&he.v->P()==vu); // back to vu (on the right)
00571                         he.FlipE();
00572                         he.FlipF();he.FlipE();he.FlipV();
00573                         vul=&he.v->P();
00574                         he.FlipV();he.FlipE();he.FlipF(); assert(&he.v->P()==vu); // back to vu (on the left)
00575                         he.FlipV();he.FlipE();he.FlipF(); assert(&he.v->P()==vl);// again on vl (but under the edge)
00576                         he.FlipE();he.FlipV();
00577                         vd=&he.v->P();
00578                         he.FlipF();he.FlipE();he.FlipV();
00579                         vdl=&he.v->P();
00580                         he.FlipV();he.FlipE();he.FlipF(); assert(&he.v->P()==vd);// back to vd (on the right)
00581                         he.FlipE();
00582                         he.FlipF();he.FlipE();he.FlipV();
00583                         vdr=&he.v->P();
00584 
00585                         nv.P()=((*vl)+(*vr))/2.0+((*vu)+(*vd))/8.0 - ((*vul)+(*vur)+(*vdl)+(*vdr))/16.0;
00586                         }
00587         }
00588 
00590         Color4<typename MESH_TYPE::ScalarType> WedgeInterp(Color4<typename MESH_TYPE::ScalarType> &c0, Color4<typename MESH_TYPE::ScalarType> &c1)
00591         {
00592                 Color4<typename MESH_TYPE::ScalarType> cc;
00593                 return cc.lerp(c0,c1,0.5f);
00594         }
00595 
00596         template<class FL_TYPE>
00597         TexCoord2<FL_TYPE,1> WedgeInterp(TexCoord2<FL_TYPE,1> &t0, TexCoord2<FL_TYPE,1> &t1)
00598         {
00599                 TexCoord2<FL_TYPE,1> tmp;
00600                 assert(t0.n()== t1.n());
00601                 tmp.n()=t0.n(); 
00602                 tmp.t()=(t0.t()+t1.t())/2.0;
00603                 return tmp;
00604         }
00605 };
00606 
00607 
00608 #if 0
00609                         int rule=0;
00610                         if(vr==vul) rule+=1;
00611                         if(vl==vur) rule+=2;
00612                         if(vl==vdr) rule+=4;
00613                         if(vr==vdl) rule+=8;
00614                         switch(rule){
00615 /*      */
00616 /*      */                      case  0 :       return ((*vl)+(*vr))/2.0+((*vu)+(*vd))/8.0 - ((*vul)+(*vur)+(*vdl)+(*vdr))/16.0;
00617 /* ul   */              case  1 : return (*vl*6 + *vr*10 + *vu + *vd*3 - *vur - *vdl -*vdr*2 )/16.0; 
00618 /* ur   */              case  2 : return (*vr*6 + *vl*10 + *vu + *vd*3 - *vul - *vdr -*vdl*2 )/16.0; 
00619 /* dr   */              case  4 : return (*vr*6 + *vl*10 + *vd + *vu*3 - *vdl - *vur -*vul*2 )/16.0; 
00620 /* dl   */              case  8 : return (*vl*6 + *vr*10 + *vd + *vu*3 - *vdr - *vul -*vur*2 )/16.0; 
00621 /* ul,ur */             case  3 : return (*vl*4 + *vr*4 + *vd*2 + - *vdr - *vdl )/8.0; 
00622 /* dl,dr */             case 12 : return (*vl*4 + *vr*4 + *vu*2 + - *vur - *vul )/8.0; 
00623 
00624 /* ul,dr */             case  5 :
00625 /* ur,dl */             case 10 :                                               
00626                                                                 default:                                
00627                                                                         return (*vl+ *vr)/2.0; 
00628                         }
00629 
00630 
00631 
00632 #endif
00633 /*
00634     vul-------vu--------vur
00635                   \      /  \      /
00636                          \    /    \    /
00637         \  /  fu  \  /
00638          vl--------vr
00639         /  \  fd  /  \
00640        /    \    /    \
00641       /      \  /      \
00642     vdl-------vd--------vdr
00643 
00644 */
00645 
00646 // Versione modificata per tenere di conto in meniara corretta dei vertici con valenza alta
00647 
00648 template<class MESH_TYPE>
00649 struct MidPointButterfly2 : public std::unary_function<face::Pos<typename MESH_TYPE::FaceType> , typename MESH_TYPE::CoordType>
00650 {
00651         typename MESH_TYPE::CoordType operator()(face::Pos<typename MESH_TYPE::FaceType>  ep)
00652         {       
00653 double Rules[11][10] = 
00654 {
00655         {.0}, // valenza 0 
00656         {.0}, // valenza 1 
00657         {.0}, // valenza 2 
00658         {  .4166666667, -.08333333333 , -.08333333333  }, // valenza 3 
00659         {  .375       ,  .0           ,  -0.125        ,  .0          }, // valenza 4                                                                                   
00660         {  .35        ,  .03090169945 ,  -.08090169945 , -.08090169945,  .03090169945   }, // valenza 5 
00661         {  .5         ,  .125         ,  -0.0625       ,  .0          ,  -0.0625      , 0.125       }, // valenza 6     
00662         {  .25        ,  .1088899050  , -.06042933822  , -.04846056675, -.04846056675, -.06042933822,  .1088899050  }, // valenza 7       
00663         {  .21875     ,  .1196383476  , -.03125        , -.05713834763, -.03125      , -.05713834763, -.03125      ,.1196383476  }, // valenza 8          
00664         {  .1944444444,  .1225409480  , -.00513312590  , -.05555555556, -.03407448880, -.03407448880, -.05555555556, -.00513312590, .1225409480  }, // valenza 9          
00665         {  .175       ,  .1213525492  ,  .01545084973  , -.04635254918, -.04045084973, -.025        , -.04045084973, -.04635254918,  .01545084973,  .1213525492  } // valenza 10          
00666 };
00667 
00668 face::Pos<typename MESH_TYPE::FaceType> he(ep.f,ep.z,ep.f->V(ep.z));
00669         typename MESH_TYPE::CoordType *vl,*vr;
00670         vl=&he.v->P();
00671         vr=&he.VFlip()->P();
00672         if(he.IsBorder())
00673                 {he.FlipV();
00674                 typename MESH_TYPE::CoordType *vl0,*vr0; 
00675                         he.NextB();
00676                         vr0=&he.v->P();
00677                         he.FlipV();
00678                         he.NextB();
00679                         assert(&he.v->P()==vl);
00680                         he.NextB();
00681                         vl0=&he.v->P();
00682                         return ((*vl)+(*vr))*(9.0/16.0)-((*vl0)+(*vr0))/16.0 ;
00683                 }
00684 
00685         int kl=0,kr=0; // valence of left and right vertices
00686         bool bl=false,br=false; // if left and right vertices are of border
00687   face::Pos<typename MESH_TYPE::FaceType> heStart=he;assert(he.v->P()==*vl);
00688         do { // compute valence of left vertex
00689                 he.FlipE();he.FlipF();
00690                 if(he.IsBorder()) bl=true; 
00691                 ++kl;
00692         }       while(he!=heStart);
00693         
00694         he.FlipV();heStart=he;assert(he.v->P()==*vr);
00695         do { // compute valence of right vertex
00696                 he.FlipE();he.FlipF();
00697                 if(he.IsBorder()) br=true; 
00698                 ++kr;
00699         }       while(he!=heStart);
00700   if(br||bl) return MidPointButterfly<MESH_TYPE>()( ep );
00701         if(kr==6 && kl==6) return MidPointButterfly<MESH_TYPE>()( ep );
00702         // TRACE("odd vertex among valences of %i %i\n",kl,kr);
00703         typename MESH_TYPE::CoordType newposl=*vl*.75, newposr=*vr*.75;
00704         he.FlipV();heStart=he; assert(he.v->P()==*vl);
00705         int i=0;
00706         if(kl!=6) 
00707         do { // compute position  of left vertex
00708                 newposl+= he.VFlip()->P() * Rules[kl][i];
00709                 he.FlipE();he.FlipF();
00710                 ++i;
00711         }       while(he!=heStart);
00712         i=0;he.FlipV();heStart=he;assert(he.v->P()==*vr);
00713         if(kr!=6)
00714         do { // compute position of right vertex
00715                 newposr+=he.VFlip()->P()* Rules[kr][i];
00716                 he.FlipE();he.FlipF();
00717                 ++i;
00718         }       while(he!=heStart);
00719         if(kr==6) return newposl;
00720         if(kl==6) return newposr;
00721         return newposl+newposr;
00722         }
00723 };
00724 
00725 /*
00726 // Nuovi punti (e.g. midpoint)
00727 template<class MESH_TYPE>
00728 struct OddPointLoop : public std::unary_function<face::Pos<typename MESH_TYPE::FaceType> , typename MESH_TYPE::CoordType>
00729 {
00730 
00731 
00732 }
00733 
00734 // vecchi punti
00735 template<class MESH_TYPE>
00736 struct EvenPointLoop : public std::unary_function<face::Pos<typename MESH_TYPE::FaceType> , typename MESH_TYPE::CoordType>
00737 {
00738 }
00739 */
00740 
00741 template<class MESH_TYPE>
00742 struct MidPointPlane : public std::unary_function<face::Pos<typename MESH_TYPE::FaceType> , typename MESH_TYPE::CoordType>
00743 {
00744         Plane3<typename MESH_TYPE::ScalarType> pl;
00745         typedef Point3<typename MESH_TYPE::ScalarType> Point3x;
00746         
00747         void operator()(typename MESH_TYPE::VertexType &nv, face::Pos<typename MESH_TYPE::FaceType> ep){
00748                 Point3x &p0=ep.f->V0(ep.z)->P();
00749                 Point3x &p1=ep.f->V1(ep.z)->P();
00750                 double pp= Distance(p0,pl)/(Distance(p0,pl) - Distance(p1,pl));
00751                 
00752     nv.P()=p1*pp + p0*(1.0-pp);
00753         }
00754 
00755         Color4<typename MESH_TYPE::ScalarType> WedgeInterp(Color4<typename MESH_TYPE::ScalarType> &c0, Color4<typename MESH_TYPE::ScalarType> &c1)
00756         {
00757                 Color4<typename MESH_TYPE::ScalarType> cc;
00758                 return cc.lerp(c0,c1,0.5f);
00759         }
00760 
00761         template<class FL_TYPE>
00762         TexCoord2<FL_TYPE,1> WedgeInterp(TexCoord2<FL_TYPE,1> &t0, TexCoord2<FL_TYPE,1> &t1)
00763         {
00764                 TexCoord2<FL_TYPE,1> tmp;
00765                 assert(t0.n()== t1.n());
00766                 tmp.n()=t0.n(); 
00767                 tmp.t()=(t0.t()+t1.t())/2.0;
00768                 return tmp;
00769         }
00770 };
00771 
00772 
00773 template <class FLT>
00774 class EdgeSplPlane
00775 {
00776         public:
00777   Plane3<FLT> pl;
00778         bool operator()(const Point3<FLT> &p0, const  Point3<FLT> &p1) const
00779         {
00780                 if(Distance(pl,p0)>0) {
00781                         if(Distance(pl,p1)>0) return false;
00782                         else return true;
00783                 }
00784                 else if(Distance(pl,p1)<=0) return false;
00785                 return true;
00786         }
00787 };
00788 
00789 
00790 template<class MESH_TYPE>
00791 struct MidPointSphere : public std::unary_function<face::Pos<typename MESH_TYPE::FaceType> , typename MESH_TYPE::CoordType>
00792 {
00793         Sphere3<typename MESH_TYPE::ScalarType> sph;
00794         typedef Point3<typename MESH_TYPE::ScalarType> Point3x;
00795         
00796         void operator()(typename MESH_TYPE::VertexType &nv, face::Pos<typename MESH_TYPE::FaceType>  ep){
00797                 Point3x &p0=ep.f->V0(ep.z)->P();
00798                 Point3x &p1=ep.f->V1(ep.z)->P();
00799     nv.P()= sph.c+((p0+p1)/2.0 - sph.c ).Normalize();
00800         }
00801 
00802         Color4<typename MESH_TYPE::ScalarType> WedgeInterp(Color4<typename MESH_TYPE::ScalarType> &c0, Color4<typename MESH_TYPE::ScalarType> &c1)
00803         {
00804                 Color4<typename MESH_TYPE::ScalarType> cc;
00805                 return cc.lerp(c0,c1,0.5f);
00806         }
00807 
00808         template<class FL_TYPE>
00809         TexCoord2<FL_TYPE,1> WedgeInterp(TexCoord2<FL_TYPE,1> &t0, TexCoord2<FL_TYPE,1> &t1)
00810         {
00811                 TexCoord2<FL_TYPE,1> tmp;
00812                 assert(t0.n()== t1.n());
00813                 tmp.n()=t0.n(); 
00814                 tmp.t()=(t0.t()+t1.t())/2.0;
00815                 return tmp;
00816         }
00817 };
00818 
00819 
00820 template <class FLT>
00821 class EdgeSplSphere
00822 {
00823         public:
00824   Sphere3<FLT> sph;
00825         bool operator()(const Point3<FLT> &p0, const  Point3<FLT> &p1) const
00826         {
00827                 if(Distance(sph,p0)>0) {
00828                         if(Distance(sph,p1)>0) return false;
00829                         else return true;
00830                 }
00831                 else if(Distance(sph,p1)<=0) return false;
00832                 return true;
00833         }
00834 };
00835 
00840 template<class TRIMESH_TYPE>
00841 struct CenterPoint : public std::unary_function<typename TRIMESH_TYPE::FacePointer, typename TRIMESH_TYPE::CoordType>
00842 {
00843         typename TRIMESH_TYPE::CoordType operator()(typename TRIMESH_TYPE::FacePointer f){
00844                 return vcg::Barycenter<typename TRIMESH_TYPE::FaceType>(*f);
00845         }
00846 };
00847 
00848 template<class TRIMESH_TYPE, class CenterPoint>
00849 void TriSplit(typename TRIMESH_TYPE::FacePointer f,
00850                                                         typename TRIMESH_TYPE::FacePointer f1,typename TRIMESH_TYPE::FacePointer f2,
00851                                                         typename TRIMESH_TYPE::VertexPointer vB, CenterPoint    Center)
00852 {
00853         vB->P() = Center(f);
00854 
00855         //i tre vertici della faccia da dividere
00856         typename TRIMESH_TYPE::VertexType* V0,*V1,*V2;
00857         V0 = f->V(0);
00858         V1 = f->V(1);
00859         V2 = f->V(2);
00860 
00861         //risistemo la faccia di partenza
00862         (*f).V(2) = &(*vB);
00863         //Faccia nuova #1
00864         (*f1).V(0) = &(*vB);
00865         (*f1).V(1) = V1;
00866         (*f1).V(2) = V2;
00867         //Faccia nuova #2
00868         (*f2).V(0) = V0;
00869         (*f2).V(1) = &(*vB);
00870         (*f2).V(2) = V2;
00871 
00872         if(f->HasFFAdjacency())
00873         {
00874                 //adiacenza delle facce adiacenti a quelle aggiunte
00875                 f->FFp(1)->FFp(f->FFi(1)) = f1;
00876                 f->FFp(2)->FFp(f->FFi(2)) = f2;
00877 
00878                 //adiacenza ff
00879                 typename TRIMESH_TYPE::FacePointer FF0,FF1,FF2;
00880                 FF0 = f->FFp(0);
00881                 FF1 = f->FFp(1);
00882                 FF2 = f->FFp(2);
00883 
00884                 //Indici di adiacenza ff
00885                 char FFi0,FFi1,FFi2;
00886                 FFi0 = f->FFi(0);
00887                 FFi1 = f->FFi(1);
00888                 FFi2 = f->FFi(2);
00889 
00890                 //adiacenza della faccia di partenza
00891                 (*f).FFp(1) = &(*f1);
00892                 (*f).FFi(1) = 0;
00893                 (*f).FFp(2) = &(*f2);
00894                 (*f).FFi(2) = 0;
00895 
00896                 //adiacenza della faccia #1
00897                 (*f1).FFp(0) = f;
00898                 (*f1).FFi(0) = 1;
00899 
00900                 (*f1).FFp(1) = FF1;
00901                 (*f1).FFi(1) = FFi1;
00902 
00903                 (*f1).FFp(2) = &(*f2);
00904                 (*f1).FFi(2) = 1;
00905 
00906                 //adiacenza della faccia #2
00907                 (*f2).FFp(0) = f;
00908                 (*f2).FFi(0) = 2;
00909 
00910                 (*f2).FFp(1) = &(*f1);
00911                 (*f2).FFi(1) = 2;
00912 
00913                 (*f2).FFp(2) = FF2;
00914                 (*f2).FFi(2) = FFi2;
00915         }
00916 }
00917 
00918 
00919 } // namespace vcg
00920 
00921 
00922 
00923 
00924 #endif