← Back to branch summary

~ubuntu-branches/ubuntu/lucid/meshlab/lucid

« back to all changes in this revision

Viewing changes to meshlab/src/meshlabplugins/filter_isoparametrization/local_parametrization.h

  • Committer: Bazaar Package Importer
  • Author(s): Teemu Ikonen
  • Date: 2009-10-08 16:40:41 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20091008164041-0c2ealqv8b8uc20c
Tags: 1.2.2-1
http://bugs.debian.org/539544
* New upstream version
* Do not build filter_isoparametrization because liblevmar dependency
  is not (yet) in Debian
* Fix compilation with gcc-4.4, thanks to Jonathan Liu for the patch
  (closes: #539544)
* rules: Add compiler variables to the qmake call (for testing with new
  GCC versions)
* io_3ds.pro: Make LIBS and INCLUDEPATH point to Debian version of lib3ds
* io_epoch.pro: Make LIBS point to Debian version of libbz2
* control:
  - Move Homepage URL to the source package section
  - Update to standards-version 3.8.3

Show diffs side-by-side

added added

removed removed

Lines of Context:
meshlab/src/meshlabplugins/filter_isoparametrization/local_parametrization.h
 
1
 
 
2
 
 
3
#include "defines.h"
 
4
///fitting
 
5
#include <vcg/space/fitting3.h>
 
6
#include <vcg/math/matrix33.h>
 
7
#include <vcg/space/triangle2.h>
 
8
#include "texcoord_optimization.h"
 
9
#include "mesh_operators.h"
 
10
 
 
11
#ifndef LOCAL_PARAMETRIZATION
 
12
#define LOCAL_PARAMETRIZATION
 
13
 
 
14
//#include <vcg/complex/trimesh/point_sampling.h>
 
15
 
 
16
//#define samples_area 80
 
17
 
 
18
template <class MeshType>
 
19
void ParametrizeExternal(MeshType &to_parametrize)
 
20
{
 
21
        typedef typename MeshType::FaceType FaceType;
 
22
        typedef typename MeshType::CoordType CoordType;
 
23
        typedef typename MeshType::ScalarType ScalarType;
 
24
        typedef typename MeshType::VertexType VertexType;
 
25
 
 
26
        std::vector<VertexType*> vertices;
 
27
        
 
28
        ///find first border vertex
 
29
        VertexType* Start=NULL;
 
30
        typename MeshType::VertexIterator Vi=to_parametrize.vert.begin();
 
31
        
 
32
        while ((Start==NULL)&&(Vi<to_parametrize.vert.end()))
 
33
        {
 
34
                if (((*Vi).IsB())&&(!(*Vi).IsD()))
 
35
                        Start=&(*Vi);
 
36
                Vi++;
 
37
        }
 
38
        if (Vi==to_parametrize.vert.end())
 
39
        {
 
40
                assert(0);
 
41
        }
 
42
        ///get sorted border vertices
 
43
        
 
44
        FindSortedBorderVertices<MeshType>(to_parametrize,Start,vertices);
 
45
 
 
46
        //assert(vertices.size()>=4);
 
47
        
 
48
        ///find perimeter
 
49
        ScalarType perimeter=0;
 
50
        int size=vertices.size();
 
51
        for (int i=0;i<size;i++)
 
52
                perimeter+=(vertices[i]->P()-vertices[(i+1)%size]->P()).Norm();
 
53
        
 
54
        ///find scaling factor
 
55
        /*ScalarType Sperimeter=(2.0*M_PI)/perimeter;*/
 
56
        
 
57
        ///set default texCoords
 
58
        for (Vi=to_parametrize.vert.begin();Vi!=to_parametrize.vert.end();Vi++)
 
59
        {
 
60
                (*Vi).T().U()=-2;
 
61
                (*Vi).T().V()=-2;
 
62
        }
 
63
 
 
64
        ///set border vertices
 
65
        typename std::vector<VertexType*>::iterator iteV;
 
66
        /*ScalarType curr_perim=0;*/
 
67
        ScalarType curr_angle=0;
 
68
        vertices[0]->T().U()=cos(curr_angle);
 
69
        vertices[0]->T().V()=sin(curr_angle);
 
70
        //for (int i=1;i<vertices.size();i++)
 
71
        //{
 
72
        //      curr_perim+=(vertices[i]->P()-vertices[(i-1)]->P()).Norm();
 
73
        //      //curr_perim+=perimeter/(ScalarType)size;
 
74
        //      curr_angle=curr_perim*Sperimeter;
 
75
        //      vertices[i]->T().U()=cos(curr_angle);
 
76
        //      vertices[i]->T().V()=sin(curr_angle);
 
77
        //      assert((vertices[i]->T().U()>=-1)&&(vertices[i]->T().U()<=1));
 
78
        //      assert((vertices[i]->T().V()>=-1)&&(vertices[i]->T().V()<=1));
 
79
        //}
 
80
        ScalarType anglediv=(2.0*M_PI)/(ScalarType)(vertices.size());
 
81
        curr_angle=0;
 
82
        for (unsigned int i=1;i<vertices.size();i++)
 
83
        {
 
84
                curr_angle+=anglediv;
 
85
                vertices[i]->T().U()=cos(curr_angle);
 
86
                vertices[i]->T().V()=sin(curr_angle);
 
87
                assert((vertices[i]->T().U()>=-1)&&(vertices[i]->T().U()<=1));
 
88
                assert((vertices[i]->T().V()>=-1)&&(vertices[i]->T().V()<=1));
 
89
        }
 
90
}
 
91
 
 
92
template <class MeshType>
 
93
void ParametrizeInternal(MeshType &to_parametrize)
 
94
{
 
95
        typedef typename MeshType::ScalarType ScalarType;
 
96
        const ScalarType Eps=(ScalarType)0.0001;
 
97
        ///set internal vertices
 
98
        for (typename MeshType::VertexIterator Vi=to_parametrize.vert.begin();Vi!=to_parametrize.vert.end();Vi++)
 
99
        {
 
100
                //assert(!Vi->IsD());
 
101
                if ((!Vi->IsB())&&(!Vi->IsD()))
 
102
                {
 
103
                        ///find kernel
 
104
                        std::vector<typename MeshType::VertexType*> star;
 
105
                        getVertexStar<MeshType>(&(*Vi),star);
 
106
                        ScalarType kernel=0;
 
107
                        for (unsigned int k=0;k<star.size();k++)
 
108
                                if (star[k]->IsB())
 
109
                                {
 
110
                                        ScalarType dist=((*Vi).P()-star[k]->P()).Norm();
 
111
                                        if (dist<Eps)
 
112
                                                dist=Eps;
 
113
                                        //ScalarType peso=1.0/(dist);
 
114
                                        ScalarType peso=(dist)/star.size();
 
115
                                        kernel+=(peso);//*dist);
 
116
                                }
 
117
                        assert(kernel>0);
 
118
                        ///then find factor
 
119
                        kernel=1.0/kernel;       
 
120
 
 
121
                        (*Vi).T().U()=0;
 
122
                        (*Vi).T().V()=0;
 
123
                        
 
124
                        int num=0;
 
125
                        ///find weighted media
 
126
                        for (unsigned int k=0;k<star.size();k++)
 
127
                                if (star[k]->IsB())
 
128
                                {
 
129
                                        ScalarType dist=((*Vi).P()-star[k]->P()).Norm();
 
130
                                        if (dist<Eps)
 
131
                                                dist=Eps;
 
132
                                        //ScalarType peso=1.0/(dist);
 
133
                                        ScalarType peso=(dist)/star.size();
 
134
                                        ScalarType kval=(peso)*kernel;
 
135
                                        assert(kval>0);
 
136
                                        (*Vi).T().U()+=kval*star[k]->T().U();
 
137
                                        (*Vi).T().V()+=kval*star[k]->T().V();
 
138
                                        num++;
 
139
                                }
 
140
                        ////on border case 2 neighbors
 
141
                        ///go next to the center
 
142
                        /*if (num==2)
 
143
                        {
 
144
                                (*Vi).T().U()/=2.0;
 
145
                                (*Vi).T().V()/=2.0;
 
146
                        }*/
 
147
                        /*ScalarType u=(*Vi).T().U();
 
148
                        ScalarType v=(*Vi).T().V();*/
 
149
                        assert(((*Vi).T().U()>=-1)&&((*Vi).T().U()<=1));
 
150
                        assert(((*Vi).T().V()>=-1)&&((*Vi).T().V()<=1));
 
151
                }
 
152
        }
 
153
        ///smoothing of txcoords
 
154
        InitDampRestUV(to_parametrize);
 
155
        for (typename MeshType::VertexIterator Vi=to_parametrize.vert.begin();Vi!=to_parametrize.vert.end();Vi++)
 
156
        {
 
157
                if ((!Vi->IsB())&&(!Vi->IsD()))
 
158
                {
 
159
                        std::vector<typename MeshType::VertexType*> star;
 
160
                        getVertexStar<MeshType>(&(*Vi),star);
 
161
                        vcg::Point2<ScalarType> UV=vcg::Point2<ScalarType>(0,0);
 
162
                        for (unsigned int k=0;k<star.size();k++)
 
163
                                UV+=star[k]->RestUV;
 
164
                        UV/=(ScalarType)star.size();
 
165
                        (*Vi).T().P()=UV;
 
166
                }
 
167
        }
 
168
}
 
169
 
 
170
 
 
171
template <class FaceType>
 
172
typename FaceType::CoordType InterpolatePos
 
173
(FaceType* f, const typename FaceType::CoordType &bary)
 
174
{return (f->V(0)->P()*bary.X()+f->V(1)->P()*bary.Y()+f->V(2)->P()*bary.Z());}
 
175
 
 
176
template <class FaceType>
 
177
typename FaceType::CoordType InterpolateRPos
 
178
(FaceType* f,const typename FaceType::CoordType &bary)
 
179
{       
 
180
        return (f->V(0)->RPos*bary.X()+f->V(1)->RPos*bary.Y()+f->V(2)->RPos*bary.Z());
 
181
}
 
182
 
 
183
template <class FaceType>
 
184
typename FaceType::CoordType InterpolateNorm
 
185
(FaceType* f, const typename FaceType::CoordType &bary)
 
186
{       
 
187
        typedef typename FaceType::CoordType CoordType;
 
188
        CoordType n0=f->V(0)->N();
 
189
        CoordType n1=f->V(1)->N();
 
190
        CoordType n2=f->V(2)->N();
 
191
        return (n0*bary.X()+n1*bary.Y()+n2*bary.Z());
 
192
}
 
193
 
 
194
template <class ScalarType>
 
195
int Approx(const ScalarType &value)
 
196
{
 
197
        ScalarType val0=floor(value);
 
198
        ScalarType val1=ceil(value);
 
199
        if (fabs(val0-value)<fabs(val1-value))
 
200
                return ((int)val0);
 
201
        else
 
202
                return ((int)val1);
 
203
}
 
204
 
 
205
template <class FaceType>
 
206
vcg::Point3i InterpolateColor
 
207
(FaceType* f,const typename FaceType::CoordType &bary)
 
208
{       
 
209
        typedef typename FaceType::ScalarType ScalarType;
 
210
        vcg::Color4b c0=f->V(0)->C();
 
211
        vcg::Color4b c1=f->V(1)->C();
 
212
        vcg::Color4b c2=f->V(2)->C();
 
213
        double R=(ScalarType)c0.X()*bary.X()+(ScalarType)c1.X()*bary.Y()+(ScalarType)c2.X()*bary.Z();
 
214
        double G=(ScalarType)c0.Y()*bary.X()+(ScalarType)c1.Y()*bary.Y()+(ScalarType)c2.Y()*bary.Z();
 
215
        double B=(ScalarType)c0.Z()*bary.X()+(ScalarType)c1.Z()*bary.Y()+(ScalarType)c2.Z()*bary.Z();
 
216
 
 
217
        vcg::Point3i p=vcg::Point3i(Approx(R),Approx(G),Approx(B));
 
218
 
 
219
        assert(p.X()<=255);
 
220
        assert(p.Y()<=255);
 
221
        assert(p.Z()<=255);
 
222
 
 
223
        return (p);
 
224
}
 
225
 
 
226
template <class VertexType>
 
227
typename VertexType::CoordType ProjectPos(const VertexType &v)
 
228
{
 
229
        typedef typename VertexType::FaceType FaceType;
 
230
        typedef typename VertexType::CoordType CoordType;
 
231
 
 
232
        FaceType *f=v.father;
 
233
        CoordType b=v.Bary;
 
234
        return (InterpolatePos<FaceType>(f,b));
 
235
}
 
236
 
 
237
template <class VertexType>
 
238
typename VertexType::CoordType Warp(const VertexType* v)
 
239
{
 
240
        typename VertexType::FaceType * father=v->father;
 
241
        typename VertexType::CoordType proj=father->P(0)*v->Bary.X()+father->P(1)*v->Bary.Y()+father->P(2)*v->Bary.Z();
 
242
        return proj;
 
243
}
 
244
 
 
245
template <class VertexType>
 
246
typename VertexType::CoordType WarpRpos(const VertexType* v)
 
247
{
 
248
        typename VertexType::FaceType * father=v->father;
 
249
        typename VertexType::CoordType proj=father->V(0)->RPos*v->Bary.X()+father->V(1)->RPos*v->Bary.Y()+father->V(2)->RPos*v->Bary.Z();
 
250
        return proj;
 
251
}
 
252
 
 
253
template <class MeshType>
 
254
typename MeshType::ScalarType EstimateLenghtByParam
 
255
                                                                (const typename MeshType::VertexType* v0,
 
256
                                                                 const typename MeshType::VertexType* v1,
 
257
                                                                 typename MeshType::FaceType* on_edge[2])
 
258
{
 
259
        typedef typename MeshType::FaceType FaceType;
 
260
        typedef typename MeshType::CoordType CoordType;
 
261
        typedef typename MeshType::VertexType VertexType;
 
262
        typedef typename MeshType::ScalarType ScalarType;
 
263
 
 
264
//      assert((on_edge.size()==0)||(on_edge.size()==1));
 
265
        ScalarType estimated[2]={0,0};
 
266
        int num[2]={0,0};
 
267
        
 
268
        for (int i=0;i<2;i++)
 
269
        {
 
270
                FaceType *test_face=on_edge[i];
 
271
 
 
272
                int edge_index=EdgeIndex(test_face,v0,v1);
 
273
                FaceType *Fopp=test_face->FFp(edge_index);
 
274
 
 
275
                if (test_face->vertices_bary.size()<2)
 
276
                {
 
277
                        ScalarType dist=Distance(v0->RPos,v1->RPos);
 
278
                        //#pragma omp atomic
 
279
                        estimated[i]+=dist;
 
280
                        num[i]=0;
 
281
                        continue;
 
282
                }
 
283
 
 
284
                ///collect vertices
 
285
                std::vector<VertexType*> vertices;
 
286
                vertices.reserve(test_face->vertices_bary.size());
 
287
                for (unsigned int k=0;k<test_face->vertices_bary.size();k++)
 
288
                        vertices.push_back(test_face->vertices_bary[k].first);
 
289
 
 
290
                
 
291
                ///collect faces
 
292
                std::vector<FaceType*> faces;
 
293
                getSharedFace<MeshType>(vertices,faces);
 
294
 
 
295
                ///get border edges
 
296
                std::vector<std::pair<VertexType*,VertexType*> > edges;
 
297
                for (unsigned int j=0;j<faces.size();j++)
 
298
                {
 
299
                        FaceType*f=faces[j];
 
300
                        ///find if there is an on-edge edge
 
301
                        bool found=false;
 
302
                        int k=0;
 
303
                        while  ((k<3)&&(!found))
 
304
                        {
 
305
                                if ((f->V0(k)->father==test_face)&&
 
306
                                        (f->V1(k)->father==test_face)&&
 
307
                                        (f->V2(k)->father==Fopp))
 
308
                                {
 
309
                                        edges.push_back(std::pair<VertexType*,VertexType*>(f->V0(k),f->V1(k)));
 
310
                                        found=true;
 
311
                                }
 
312
                                k++;
 
313
                        }
 
314
                }
 
315
                
 
316
                ///find if there's ant edge return inital lenght
 
317
                if (edges.size()==0)
 
318
                {
 
319
                        estimated[i]+=(Distance(v0->RPos,v1->RPos));
 
320
                        num[i]=0;
 
321
                        continue;
 
322
                }
 
323
                else
 
324
                {
 
325
                        //get edge direction
 
326
                        ///store the two nearest for each vertex
 
327
                        /*VertexType *n0=edges[0].first;
 
328
                        VertexType *n1=edges[0].second;
 
329
                        ScalarType d0=(Warp(n0)-v0->P()).Norm();
 
330
                        ScalarType d1=(Warp(n1)-v1->P()).Norm();*/
 
331
 
 
332
                        //CoordType edgedir=v0->cP()-v1->cP();
 
333
                        CoordType edgedir=v0->RPos-v1->RPos;
 
334
                        edgedir.Normalize();
 
335
                        num[i]=edges.size();
 
336
                        for (unsigned int e=0;e<edges.size();e++)
 
337
                        {
 
338
                                VertexType *vH0=edges[e].first;
 
339
                                VertexType *vH1=edges[e].second;
 
340
 
 
341
                                ///project points over the plane 
 
342
                                /*CoordType proj0=Warp(vH0);
 
343
                                CoordType proj1=Warp(vH1);*/
 
344
                                CoordType proj0=WarpRpos(vH0);
 
345
                                CoordType proj1=WarpRpos(vH1);
 
346
                                
 
347
                                CoordType dirproj=proj0-proj1;
 
348
                                dirproj.Normalize();
 
349
                                //estimated[i]+=fabs(dirproj*edgedir)*((vH0->P()-vH1->P()).Norm());
 
350
                                //#pragma omp atomic
 
351
                                estimated[i]+=fabs(dirproj*edgedir)*((vH0->RPos-vH1->RPos).Norm());
 
352
 
 
353
                        }
 
354
                }
 
355
        }
 
356
 
 
357
        ///media of estimated values
 
358
        ScalarType alpha0,alpha1;
 
359
        ScalarType max_num=abstraction_num;
 
360
 
 
361
        if (num[0]>=max_num)
 
362
                alpha0=1;
 
363
        else
 
364
                alpha0=num[0]/max_num;
 
365
 
 
366
        if (num[1]>=max_num)
 
367
                alpha1=1;
 
368
        else
 
369
                alpha1=num[1]/max_num;
 
370
 
 
371
        estimated[0]=alpha0*estimated[0]+(1.0-alpha0)*(Distance(v0->RPos,v1->RPos));
 
372
        estimated[1]=alpha1*estimated[1]+(1.0-alpha1)*(Distance(v0->RPos,v1->RPos));
 
373
        return((estimated[0]+estimated[1])/2.0);
 
374
}
 
375
 
 
376
template <class MeshType>
 
377
void MeanVal(const std::vector<vcg::Point2<typename MeshType::ScalarType> > &Points,
 
378
                         std::vector<typename MeshType::ScalarType> &Lamda,
 
379
                         typename MeshType::CoordType &p)
 
380
{
 
381
        typedef typename MeshType::FaceType FaceType;
 
382
        typedef typename MeshType::VertexType VertexType;
 
383
        typedef typename MeshType::ScalarType ScalarType;
 
384
        typedef typename MeshType::CoordType CoordType;
 
385
 
 
386
        int size=Points.size();
 
387
        Lamda.resize(size);
 
388
 
 
389
        ScalarType sum=0;       
 
390
        for (int i=0;i<size;i++)
 
391
        {
 
392
                int size=Points.size()-1;
 
393
                vcg::Point2<ScalarType> Pcurr=Points[i];
 
394
                vcg::Point2<ScalarType> Pprev=Points[(i+(size-1))%size];
 
395
                vcg::Point2<ScalarType> Pnext=Points[(i+1)%size];
 
396
 
 
397
                CoordType v0=Pprev-p;
 
398
                CoordType v1=Pcurr-p;
 
399
                CoordType v2=Pnext-p;
 
400
                ScalarType l=v1.Norm();
 
401
                v0.Normalize();
 
402
                v1.Normalize();
 
403
                v2.Normalize();
 
404
 
 
405
                ScalarType Alpha0=acos(v0*v1);
 
406
                ScalarType Alpha1=acos(v1*v2);
 
407
 
 
408
                Lamda[i]=(tan(Alpha0/2.0)+tan(Alpha1/2.0))/l;
 
409
                sum+=Lamda[i];
 
410
        }
 
411
 
 
412
        ///normalization
 
413
        for (int i=0;i<size;i++)
 
414
                Lamda[i]/=sum;
 
415
}
 
416
 
 
417
template <class FaceType>
 
418
typename FaceType::ScalarType EstimateAreaByParam(const FaceType* f)
 
419
{
 
420
        typedef typename FaceType::VertexType VertexType;
 
421
        typedef typename FaceType::ScalarType ScalarType;
 
422
        
 
423
        int num=0;
 
424
        ScalarType estimated=0;
 
425
        for (unsigned int k=0;k<f->vertices_bary.size();k++)
 
426
        {
 
427
                VertexType *HresVert=f->vertices_bary[k].first;
 
428
                estimated+=HresVert->area;
 
429
                num++;
 
430
        }
 
431
 
 
432
        ///media of estimated values
 
433
        ScalarType alpha;
 
434
        ScalarType max_num=abstraction_num;
 
435
 
 
436
        if (num>=max_num)
 
437
                alpha=1;
 
438
        else
 
439
                alpha=num/max_num;
 
440
 
 
441
        ScalarType Rarea=((f->V(1)->RPos-f->V(0)->RPos)^(f->V(2)->RPos-f->V(0)->RPos)).Norm()/2.0;
 
442
        estimated=alpha*estimated+(1.0-alpha)*Rarea;
 
443
 
 
444
        return(estimated);
 
445
}
 
446
 
 
447
template <class MeshType>
 
448
typename MeshType::ScalarType EstimateAreaByParam
 
449
                                                          (const typename MeshType::VertexType* v0,
 
450
                                                           const typename MeshType::VertexType* v1,
 
451
                                                           typename MeshType::FaceType* on_edge[2])
 
452
{
 
453
        typedef typename MeshType::FaceType FaceType;
 
454
        typedef typename MeshType::VertexType VertexType;
 
455
        typedef typename MeshType::ScalarType ScalarType;
 
456
 
 
457
        //MeshType::PerVertexAttributeHandle<AuxiliaryVertData> handle = vcg::tri::Allocator<MeshType>::GetPerVertexAttribute<AuxiliaryVertData>(mesh,"AuxiliaryVertData");
 
458
 
 
459
        ScalarType estimated[2]={0,0};
 
460
        int num[2]={0,0};
 
461
        VertexType *v2[2];
 
462
        for (int i=0;i<2;i++)
 
463
        {
 
464
                FaceType *test_face=on_edge[i];
 
465
 
 
466
                for (int k=0;k<3;k++)
 
467
                        if ((test_face->V(k)!=v0)&&(test_face->V(k)!=v1))
 
468
                                v2[i]=test_face->V(k);
 
469
                
 
470
                for (unsigned int k=0;k<test_face->vertices_bary.size();k++)
 
471
                {
 
472
                        VertexType *brother=test_face->vertices_bary[k].first;
 
473
                        estimated[i]+=brother->area;
 
474
                        //estimated[i]+=handle[brother].area;
 
475
                        num[i]++;
 
476
                }
 
477
        }
 
478
 
 
479
        ///media of estimated values
 
480
        ScalarType alpha0,alpha1;
 
481
        ScalarType max_num=abstraction_num;//20
 
482
 
 
483
        if (num[0]>=max_num)
 
484
                alpha0=1;
 
485
        else
 
486
                alpha0=num[0]/max_num;
 
487
 
 
488
        if (num[1]>=max_num)
 
489
                alpha1=1;
 
490
        else
 
491
                alpha1=num[1]/max_num;
 
492
        ScalarType Rarea0=((on_edge[0]->V(1)->RPos-on_edge[0]->V(0)->RPos)^(on_edge[0]->V(2)->RPos-on_edge[0]->V(0)->RPos)).Norm()/2.0;
 
493
        ScalarType Rarea1=((on_edge[1]->V(1)->RPos-on_edge[1]->V(0)->RPos)^(on_edge[1]->V(2)->RPos-on_edge[1]->V(0)->RPos)).Norm()/2.0;
 
494
        estimated[0]=alpha0*estimated[0]+(1.0-alpha0)*Rarea0;
 
495
        estimated[1]=alpha1*estimated[1]+(1.0-alpha1)*Rarea1;
 
496
        return((estimated[0]+estimated[1])/2.0);
 
497
}
 
498
 
 
499
///template class used to sample surface
 
500
template <class FaceType>
 
501
class VertexSampler{
 
502
        typedef typename FaceType::CoordType CoordType;
 
503
public: 
 
504
        std::vector<CoordType> points;
 
505
 
 
506
        void AddFace(const FaceType &f,const CoordType & bary)
 
507
        {points.push_back(f.P(0)*bary.X()+f.P(1)*bary.Y()+f.P(2)*bary.Z());}
 
508
};
 
509
 
 
510
 
 
511
///sample 3d vertex possible's position
 
512
///using area criterion
 
513
//template <class MeshType>
 
514
//void SamplingPoints(MeshType &mesh,
 
515
//                                      std::vector<typename MeshType::CoordType> &pos)
 
516
//{
 
517
//      typedef typename MeshType::CoordType CoordType;
 
518
//      typedef VertexSampler<MeshType::FaceType> Sampler;
 
519
//      pos.reserve(samples_area);
 
520
//      Sampler ps;
 
521
//      ps.points.reserve(samples_area);
 
522
//
 
523
//      vcg::tri::SurfaceSampling<MeshType,Sampler>::Montecarlo(mesh,ps,samples_area);
 
524
//      pos=std::vector<CoordType>(ps.points.begin(),ps.points.end());
 
525
//}
 
526
 
 
527
template <class MeshType> 
 
528
void InitDampRestUV(MeshType &m)
 
529
{
 
530
        for (unsigned int i=0;i<m.vert.size();i++)
 
531
                m.vert[i].RestUV=m.vert[i].T().P();
 
532
}
 
533
 
 
534
 
 
535
template <class MeshType> 
 
536
void RestoreRestUV(MeshType &m)
 
537
{
 
538
        for (unsigned int i=0;i<m.vert.size();i++)
 
539
                m.vert[i].T().P()=m.vert[i].RestUV;
 
540
}
 
541
 
 
542
 
 
543
///parametrize a submesh from trinagles that are incident on vertices 
 
544
template <class MeshType>
 
545
void ParametrizeLocally(MeshType &parametrized,
 
546
                                                bool fix_boundary=true,
 
547
                                                bool init_border=true)
 
548
{
 
549
        typedef typename MeshType::FaceType FaceType;
 
550
        typedef typename MeshType::ScalarType ScalarType;
 
551
        typedef typename MeshType::CoordType CoordType;
 
552
 
 
553
        //const ScalarType epsilon=(ScalarType)0.0001;
 
554
        
 
555
 
 
556
        ///save old positions
 
557
 
 
558
        std::vector<CoordType> positions;
 
559
        positions.resize(parametrized.vert.size());
 
560
        ///set rest position
 
561
        for (unsigned int i=0;i<parametrized.vert.size();i++)
 
562
        {
 
563
                positions[i]=parametrized.vert[i].P();
 
564
                parametrized.vert[i].P()=parametrized.vert[i].RPos;
 
565
        }
 
566
 
 
567
        UpdateTopologies(&parametrized);
 
568
 
 
569
        if (init_border)
 
570
                ParametrizeExternal(parametrized);
 
571
 
 
572
        ParametrizeInternal(parametrized);
 
573
        
 
574
        //for (int i=0;i<parametrized.vert.size();i++)
 
575
        //      parametrized.vert[i].RPos=parametrized.vert[i].P();
 
576
 
 
577
        vcg::tri::MeanValueTexCoordOptimization<MeshType> opt(parametrized);
 
578
        //vcg::tri::WachspressTexCoordOptimization<MeshType> opt(parametrized);
 
579
        vcg::tri::AreaPreservingTexCoordOptimization<MeshType> opt1(parametrized);
 
580
        opt.SetSpeed((ScalarType)0.0005);
 
581
        InitDampRestUV(parametrized);
 
582
        if (fix_boundary)
 
583
        {
 
584
                opt.TargetEquilateralGeometry();
 
585
                //opt.TargetCurrentGeometry();
 
586
                opt.SetBorderAsFixed();
 
587
                opt.IterateUntilConvergence((ScalarType)0.000001,100);
 
588
                /*opt.Iterate();
 
589
                opt.Iterate();*/
 
590
        }
 
591
        else
 
592
        {
 
593
                opt1.TargetCurrentGeometry();
 
594
                //opt1.SetBorderAsFixed();
 
595
                opt1.IterateUntilConvergence((ScalarType)0.000001,200);
 
596
        }
 
597
 
 
598
        ///assert parametrization
 
599
        for (unsigned int i=0;i<parametrized.face.size();i++)
 
600
        {
 
601
                FaceType *f=&parametrized.face[i];
 
602
                vcg::Point2<ScalarType> tex0=vcg::Point2<ScalarType>(f->V(0)->T().U(),f->V(0)->T().V());
 
603
                vcg::Point2<ScalarType> tex1=vcg::Point2<ScalarType>(f->V(1)->T().U(),f->V(1)->T().V());
 
604
                vcg::Point2<ScalarType> tex2=vcg::Point2<ScalarType>(f->V(2)->T().U(),f->V(2)->T().V());
 
605
                vcg::Triangle2<typename MeshType::ScalarType> t2d=vcg::Triangle2<typename MeshType::ScalarType>(tex0,tex1,tex2);
 
606
#ifndef NDEBUG
 
607
                ScalarType area=(tex1-tex0)^(tex2-tex0);
 
608
                assert(area>0);
 
609
#endif
 
610
//#ifndef _MESHLAB
 
611
//              if (area<0){
 
612
//                      vcg::tri::io::ExporterPLY<BaseMesh>::Save(parametrized,"case0.ply");
 
613
//
 
614
//                      for (int j=0;j<parametrized.vert.size();j++)
 
615
//                      {
 
616
//                              parametrized.vert[j].P().V(0)=parametrized.vert[j].T().U();
 
617
//                              parametrized.vert[j].P().V(1)=parametrized.vert[j].T().V();
 
618
//                              parametrized.vert[j].P().V(2)=0;
 
619
//                      }
 
620
//                      vcg::tri::io::ExporterPLY<BaseMesh>::Save(parametrized,"case1.ply");
 
621
//                      assert(0);
 
622
//              }
 
623
//#endif
 
624
        }
 
625
        ///restore position
 
626
        for (unsigned int i=0;i<parametrized.vert.size();i++)   
 
627
                parametrized.vert[i].P()=positions[i];
 
628
}
 
629
 
 
630
template <class MeshType>
 
631
void ForceInParam(vcg::Point2<typename MeshType::ScalarType> &UV,MeshType &domain)
 
632
{
 
633
        typedef typename MeshType::FaceType FaceType;
 
634
        typedef typename MeshType::CoordType CoordType;
 
635
        typedef typename MeshType::ScalarType ScalarType;
 
636
        typedef typename MeshType::VertexType VertexType;
 
637
 
 
638
        ScalarType minDist=(ScalarType)1000.0;
 
639
        vcg::Point2<ScalarType> closest;
 
640
        vcg::Point2<ScalarType> center=vcg::Point2<ScalarType>(0,0);
 
641
        for (unsigned int i=0;i<domain.face.size();i++) 
 
642
        {
 
643
                FaceType *f=&domain.face[i];
 
644
                vcg::Point2<ScalarType> tex0=vcg::Point2<ScalarType>(f->V(0)->T().U(),f->V(0)->T().V());
 
645
                vcg::Point2<ScalarType> tex1=vcg::Point2<ScalarType>(f->V(1)->T().U(),f->V(1)->T().V());
 
646
                vcg::Point2<ScalarType> tex2=vcg::Point2<ScalarType>(f->V(2)->T().U(),f->V(2)->T().V());
 
647
                center+=tex0;
 
648
                center+=tex1;
 
649
                center+=tex2;
 
650
                vcg::Triangle2<ScalarType> t2d=vcg::Triangle2<ScalarType>(tex0,tex1,tex2);
 
651
                ScalarType dist;
 
652
                vcg::Point2<ScalarType> temp;
 
653
                t2d.PointDistance(UV,dist,temp);
 
654
                if (dist<minDist)
 
655
                {
 
656
                        minDist=dist;
 
657
                        closest=temp;
 
658
                }
 
659
        }
 
660
        center/=(ScalarType)(domain.face.size()*3);
 
661
        UV=closest*(ScalarType)0.95+center*(ScalarType)0.05;
 
662
}
 
663
 
 
664
 
 
665
template <class VertexType>
 
666
bool testParamCoords(VertexType *v)
 
667
{
 
668
        typedef typename VertexType::ScalarType ScalarType;
 
669
        ScalarType eps=(ScalarType)0.00001;
 
670
        if (!(((v->T().P().X()>=-1-eps)&&(v->T().P().X()<=1+eps)&&
 
671
                   (v->T().P().Y()>=-1-eps)&&(v->T().P().Y()<=1+eps))))
 
672
                        return (false);
 
673
 
 
674
        return true;
 
675
}
 
676
 
 
677
template <class MeshType>
 
678
bool testParamCoords(MeshType &domain)
 
679
{
 
680
        for (unsigned int i=0;i<domain.vert.size();i++)
 
681
        {
 
682
                typename MeshType::VertexType *v=&domain.vert[i];
 
683
                bool b=testParamCoords<typename MeshType::VertexType>(v);
 
684
                if (!b)
 
685
                {
 
686
                        #ifndef _MESHLAB
 
687
                        printf("\n position %lf,%lf \n",v->T().U(),v->T().V());
 
688
                        #endif
 
689
                        return false;
 
690
                }
 
691
        }
 
692
        return true;
 
693
}
 
694
 
 
695
template <class CoordType>
 
696
bool testBaryCoords(CoordType &bary)
 
697
{
 
698
        typedef typename CoordType::ScalarType ScalarType;
 
699
        ///test
 
700
        float eps=(ScalarType)0.0001;
 
701
        if(!(fabs(bary.X()+bary.Y()+bary.Z()-1.0)<eps))
 
702
                return false;
 
703
        if(!((bary.X()<=1.0)&&(bary.X()>=0)&&(bary.Y()<=1.0)&&(bary.Y()>=0)&&(bary.Z()<=1.0)&&(bary.Z()>=0)))
 
704
                return false;
 
705
        return true;
 
706
}
 
707
 
 
708
 
 
709
template <class MeshType>
 
710
bool testParametrization(MeshType &domain,
 
711
                                                 MeshType &Hlev)
 
712
{
 
713
        typedef typename MeshType::FaceType FaceType;
 
714
        typedef typename MeshType::CoordType CoordType;
 
715
        typedef typename MeshType::ScalarType ScalarType;
 
716
        typedef typename MeshType::VertexType VertexType;
 
717
        bool is_good=true;
 
718
        for (unsigned int i=0;i<Hlev.vert.size();i++)
 
719
        {
 
720
                VertexType *v=&Hlev.vert[i];
 
721
                if (v->father==NULL)
 
722
                {
 
723
                        printf("\n PAR ERROR : father NULL\n");
 
724
                        is_good=false;
 
725
                }
 
726
                if (v->father->IsD())
 
727
                {
 
728
                        printf("\n PAR ERROR : father DELETED \n");
 
729
                        is_good=false;
 
730
                }
 
731
                if (!(((v->Bary.X()>=0)&&(v->Bary.X()<=1))&&
 
732
                ((v->Bary.Y()>=0)&&(v->Bary.Y()<=1))&&
 
733
                ((v->Bary.Z()>=0)&&(v->Bary.Z()<=1))))
 
734
                {
 
735
                        printf("\n PAR ERROR : bary coords exceeds: %f,%f,%f \n",v->Bary.X(),v->Bary.Y(),v->Bary.Z());
 
736
                        is_good=false;
 
737
                }
 
738
        }
 
739
        for (unsigned int i=0;i<domain.face.size();i++)
 
740
        {
 
741
                FaceType *face=&domain.face[i];
 
742
                if (!face->IsD())
 
743
                {
 
744
                        for (unsigned int j=0;j<face->vertices_bary.size();j++)
 
745
                        {
 
746
                                VertexType *v=face->vertices_bary[j].first;
 
747
                                if (v->father!=face)
 
748
                                {
 
749
                                        printf("\n PAR ERROR : Father<->son \n");
 
750
                                        is_good=false;
 
751
                                }
 
752
                        }
 
753
                }
 
754
        }
 
755
        return (is_good);
 
756
}
 
757
 
 
758
template <class MeshType>
 
759
bool NonFolded(MeshType &parametrized)
 
760
{
 
761
        //const ScalarType epsilon=0.00001;
 
762
        typedef typename MeshType::FaceType FaceType;
 
763
        typedef typename MeshType::CoordType CoordType;
 
764
        typedef typename MeshType::ScalarType ScalarType;
 
765
        typedef typename MeshType::VertexType VertexType;
 
766
        ///assert parametrization
 
767
        for (unsigned int i=0;i<parametrized.face.size();i++)
 
768
        {
 
769
                FaceType *f=&parametrized.face[i];
 
770
                if (!((f->V(0)->IsB())&&(f->V(1)->IsB())&&(f->V(2)->IsB())))
 
771
                {
 
772
 
 
773
                        vcg::Point2<ScalarType> tex0=vcg::Point2<ScalarType>(f->V(0)->T().U(),f->V(0)->T().V());
 
774
                        vcg::Point2<ScalarType> tex1=vcg::Point2<ScalarType>(f->V(1)->T().U(),f->V(1)->T().V());
 
775
                        vcg::Point2<ScalarType> tex2=vcg::Point2<ScalarType>(f->V(2)->T().U(),f->V(2)->T().V());
 
776
                        vcg::Triangle2<ScalarType> t2d=vcg::Triangle2<ScalarType>(tex0,tex1,tex2);
 
777
                        ScalarType area=(tex1-tex0)^(tex2-tex0);
 
778
                        if (area<=0)
 
779
                                return false;
 
780
                }
 
781
        }
 
782
        return true;
 
783
}
 
784
 
 
785
template <class MeshType>
 
786
bool NonFolded(MeshType &parametrized,std::vector<typename MeshType::FaceType*> &folded)
 
787
{
 
788
        
 
789
        typedef typename MeshType::FaceType FaceType;
 
790
        typedef typename MeshType::CoordType CoordType;
 
791
        typedef typename MeshType::ScalarType ScalarType;
 
792
        typedef typename MeshType::VertexType VertexType;
 
793
 
 
794
        const ScalarType epsilon=(ScalarType)0.00001;
 
795
        folded.resize(0);
 
796
        ///assert parametrization
 
797
        for (unsigned int i=0;i<parametrized.face.size();i++)
 
798
        {
 
799
                FaceType *f=&parametrized.face[i];
 
800
                if (!((f->V(0)->IsB())&&(f->V(1)->IsB())&&(f->V(2)->IsB())))
 
801
                {
 
802
 
 
803
                        vcg::Point2<ScalarType> tex0=vcg::Point2<ScalarType>(f->V(0)->T().U(),f->V(0)->T().V());
 
804
                        vcg::Point2<ScalarType> tex1=vcg::Point2<ScalarType>(f->V(1)->T().U(),f->V(1)->T().V());
 
805
                        vcg::Point2<ScalarType> tex2=vcg::Point2<ScalarType>(f->V(2)->T().U(),f->V(2)->T().V());
 
806
                        vcg::Triangle2<ScalarType> t2d=vcg::Triangle2<ScalarType>(tex0,tex1,tex2);
 
807
                        ScalarType area=(tex1-tex0)^(tex2-tex0);
 
808
                        if (area<=epsilon)
 
809
                                folded.push_back(f);
 
810
                }
 
811
        }
 
812
        return (folded.size()==0);
 
813
}
 
814
//getFoldedFaces(std::vector)
 
815
///parametrize a submesh from trinagles that are incident on vertices with equi-area subdivision
 
816
template <class MeshType>
 
817
void ParametrizeStarEquilateral(MeshType &parametrized,
 
818
                                                                const typename MeshType::ScalarType &radius=1)
 
819
{
 
820
        typedef typename MeshType::FaceType FaceType;
 
821
        typedef typename MeshType::CoordType CoordType;
 
822
        typedef typename MeshType::ScalarType ScalarType;
 
823
        typedef typename MeshType::VertexType VertexType;
 
824
        
 
825
        UpdateTopologies(&parametrized);
 
826
 
 
827
        //set borders
 
828
 
 
829
        ///find first border & non border vertex
 
830
        std::vector<VertexType*> non_border;
 
831
        VertexType* Start=NULL;
 
832
        for (unsigned int i=0;i<parametrized.vert.size();i++)
 
833
        {
 
834
                VertexType* vert=&parametrized.vert[i];
 
835
                if ((Start==NULL)&&(vert->IsB()))
 
836
                        Start=vert;
 
837
 
 
838
                if (!vert->IsB())
 
839
                        non_border.push_back(vert);
 
840
        }
 
841
        assert(non_border.size()>0);
 
842
 
 
843
        ///get sorted border vertices
 
844
        std::vector<VertexType*> vertices;
 
845
        FindSortedBorderVertices<MeshType>(parametrized,Start,vertices);
 
846
        
 
847
        ///set border vertices
 
848
        int num=vertices.size();
 
849
        typename std::vector<VertexType*>::iterator iteV;
 
850
        ScalarType curr_angle=0;
 
851
        vertices[0]->T().U()=cos(curr_angle)*radius;
 
852
        vertices[0]->T().V()=sin(curr_angle)*radius;
 
853
        ScalarType division=(2*M_PI)/(ScalarType)num;
 
854
        ///set border
 
855
        for (unsigned int i=1;i<vertices.size();i++)
 
856
        {
 
857
                curr_angle+=division;
 
858
                vertices[i]->T().U()=radius*cos(curr_angle);
 
859
                vertices[i]->T().V()=radius*sin(curr_angle);
 
860
        }
 
861
 
 
862
        if (non_border.size()==1)
 
863
        {
 
864
                ///if non-border vertex is one then set it to zero otherwise 
 
865
                ///set it to the average of neighbors
 
866
                non_border[0]->T().P()=vcg::Point2<ScalarType>(0,0);
 
867
        }
 
868
        else
 
869
        {
 
870
                ///set media of star vertices
 
871
                assert(non_border.size()==2);
 
872
                for (unsigned int i=0;i<non_border.size();i++)
 
873
                {
 
874
                        VertexType *v=non_border[i];
 
875
                        v->T().P()=vcg::Point2<ScalarType>(0,0);
 
876
                        int ariety_vert=0;
 
877
                        std::vector<VertexType*> star;
 
878
                        getVertexStar<MeshType>(v,star);
 
879
                        for (unsigned int k=0;k<star.size();k++)
 
880
                        {
 
881
                                if ((!star[k]->IsD())&&(star[k]->IsB()))
 
882
                                {
 
883
                                        v->T().P()+=star[k]->T().P();
 
884
                                        ariety_vert++;
 
885
                                }
 
886
                        }
 
887
                        v->T().P()/=(ScalarType)ariety_vert;
 
888
                }
 
889
                ///test particular cases
 
890
                if (!NonFolded(parametrized))
 
891
                {
 
892
                        std::vector<VertexType*> shared;
 
893
                        getSharedVertexStar<MeshType>(non_border[0],non_border[1],shared);
 
894
                
 
895
                        assert(shared.size()==2);
 
896
                        assert(shared[0]->IsB());
 
897
                        assert(shared[1]->IsB());
 
898
                        assert(shared[0]!=shared[1]);
 
899
                
 
900
                        //ScalarType epsilon=(ScalarType)0.001;
 
901
                        ///then get the media of two shared vertices
 
902
                        vcg::Point2<ScalarType> uvAve=shared[0]->T().P()+shared[1]->T().P();
 
903
                        assert(uvAve.Norm()>(ScalarType)0.001);
 
904
                        uvAve.Normalize();
 
905
                        vcg::Point2<ScalarType> p0=uvAve*(ScalarType)0.3;
 
906
                        vcg::Point2<ScalarType> p1=uvAve*(ScalarType)(-0.3);
 
907
                        ///then test and set right assignement
 
908
                        non_border[0]->T().P()=p0;
 
909
                        non_border[1]->T().P()=p1;
 
910
                        if (!NonFolded(parametrized)){
 
911
                                non_border[0]->T().P()=p1;
 
912
                                non_border[1]->T().P()=p0;
 
913
                        }
 
914
                }
 
915
        }
 
916
#ifndef _MESHLAB
 
917
        ///final assert parametrization
 
918
        if (!NonFolded(parametrized))
 
919
        {
 
920
                printf("ERROR FOLDED PARAMETRIZATION\n");
 
921
                vcg::tri::io::ExporterPLY<MeshType>::Save(parametrized,"case0.ply");
 
922
 
 
923
                for (int j=0;j<parametrized.vert.size();j++)
 
924
                {
 
925
                        parametrized.vert[j].P().V(0)=parametrized.vert[j].T().U();
 
926
                        parametrized.vert[j].P().V(1)=parametrized.vert[j].T().V();
 
927
                        parametrized.vert[j].P().V(2)=0;
 
928
                }
 
929
                vcg::tri::io::ExporterPLY<MeshType>::Save(parametrized,"case1.ply");
 
930
                assert(0);
 
931
        }
 
932
#endif
 
933
}
 
934
 
 
935
///given the mesh and the two edges (same) seen from face[0] and face[1] of the mesh construct
 
936
///a diamond parametrization using equilateral triangles of edge edge_len
 
937
template <class MeshType>
 
938
void ParametrizeDiamondEquilateral(MeshType &parametrized,
 
939
                                                                   const int &edge0,const int &edge1,
 
940
                                                                   const typename MeshType::ScalarType &edge_len=1)
 
941
{
 
942
 
 
943
        typedef typename MeshType::FaceType FaceType;
 
944
        typedef typename FaceType::CoordType CoordType;
 
945
        typedef typename FaceType::ScalarType ScalarType;
 
946
        typedef typename FaceType::VertexType VertexType;
 
947
        
 
948
        ScalarType h=(sqrt(3.0)/2.0)*edge_len;
 
949
 
 
950
        FaceType *fd0=&parametrized.face[0];
 
951
#ifndef NDEBUG
 
952
        FaceType *fd1=&parametrized.face[1];
 
953
#endif
 
954
        assert(fd0->FFp(edge0)==fd1);
 
955
        assert(fd1->FFp(edge1)==fd0);
 
956
 
 
957
        ///get 2 vertex on the edge
 
958
        VertexType *v0=fd0->V(edge0);
 
959
        VertexType *v1=fd0->V((edge0+1)%3);
 
960
 
 
961
#ifndef NDEBUG
 
962
        VertexType *vtest0=fd1->V(edge1);
 
963
        VertexType *vtest1=fd1->V((edge1+1)%3);
 
964
 
 
965
        assert(v0!=v1);
 
966
        assert(vtest0!=vtest1);
 
967
        assert(((v0==vtest0)&&(v1==vtest1))||((v1==vtest0)&&(v0==vtest1)));
 
968
#endif
 
969
 
 
970
        ///other 2 vertex
 
971
        VertexType *v2=parametrized.face[0].V((edge0+2)%3);
 
972
        VertexType *v3=parametrized.face[1].V((edge1+2)%3);
 
973
        assert((v2!=v3)&&(v0!=v2)&&(v0!=v3)&&(v1!=v2)&&(v1!=v3));
 
974
 
 
975
        ///assing texcoords
 
976
        v0->T().P()=vcg::Point2<ScalarType>(0,-edge_len/2.0);
 
977
        v1->T().P()=vcg::Point2<ScalarType>(0,edge_len/2.0);
 
978
        v2->T().P()=vcg::Point2<ScalarType>(-h,0);
 
979
        v3->T().P()=vcg::Point2<ScalarType>(h,0);
 
980
 
 
981
        ///test
 
982
        assert(NonFolded(parametrized));
 
983
}
 
984
 
 
985
///given the mesh and the two edges (same) seen from face[0] and face[1] of the mesh construct
 
986
///a diamond parametrization using equilateral triangles of edge edge_len
 
987
template <class MeshType>
 
988
void ParametrizeFaceEquilateral(MeshType &parametrized,
 
989
                                                                   const typename MeshType::ScalarType &edge_len=1)
 
990
{
 
991
        typedef typename MeshType::FaceType FaceType;
 
992
        typedef typename FaceType::CoordType CoordType;
 
993
        typedef typename FaceType::ScalarType ScalarType;
 
994
        typedef typename FaceType::VertexType VertexType;
 
995
        
 
996
        ScalarType h=(sqrt(3.0)/2.0)*edge_len;
 
997
 
 
998
        FaceType *f_param=&(parametrized.face[0]);
 
999
                                
 
1000
        f_param->V(0)->T().P()=vcg::Point2<ScalarType>(edge_len/2.0,0);
 
1001
        f_param->V(1)->T().P()=vcg::Point2<ScalarType>(0,h);
 
1002
        f_param->V(2)->T().P()=vcg::Point2<ScalarType>(-edge_len/2.0,0);
 
1003
}
 
1004
 
 
1005
 
 
1006
///parametrize and create a submesh from trinagles that are incident on
 
1007
/// vertices .... seturn a vetor of original faces
 
1008
template <class MeshType>
 
1009
void ParametrizeLocally(MeshType &parametrized,
 
1010
                                                const std::vector<typename MeshType::VertexType*> &subset,
 
1011
                                                std::vector<typename MeshType::FaceType*> &orderedFaces,
 
1012
                                                std::vector<typename MeshType::VertexType*> &orderedVertex)
 
1013
{
 
1014
        typedef typename MeshType::FaceType FaceType;
 
1015
        typedef typename FaceType::CoordType CoordType;
 
1016
        typedef typename FaceType::ScalarType ScalarType;
 
1017
        typedef typename FaceType::VertexType VertexType;
 
1018
 
 
1019
        orderedFaces.clear();
 
1020
        std::vector<VertexType*> vertices;
 
1021
        
 
1022
        
 
1023
        ///get faces referenced by vertices
 
1024
        getSharedFace<MeshType>(subset,orderedFaces);
 
1025
        
 
1026
        ///do a first copy of the mesh
 
1027
        ///and parametrize it
 
1028
        ///NB: order of faces is mantained
 
1029
        CopyMeshFromFaces<MeshType>(orderedFaces,orderedVertex,parametrized);
 
1030
        
 
1031
        //CreateMeshVertexStar(subset,orderedFaces,parametrized);
 
1032
        ParametrizeLocally(parametrized);       
 
1033
        
 
1034
}
 
1035
 
 
1036
 
 
1037
 
 
1038
template <class MeshType>
 
1039
void GetUV(const typename MeshType::FaceType* f,
 
1040
                   const typename MeshType::CoordType &bary,
 
1041
                   typename MeshType::ScalarType &U,
 
1042
                   typename MeshType::ScalarType &V)
 
1043
{
 
1044
        U=bary.X()*f->V(0)->T().U()+bary.Y()*f->V(1)->T().U()+bary.Z()*f->V(2)->T().U();
 
1045
        V=bary.X()*f->V(0)->T().V()+bary.Y()*f->V(1)->T().V()+bary.Z()*f->V(2)->T().V();
 
1046
        /*if ((!((U>=-1)&&(U<=1)))||(!((V>=-1)&&(V<=1))))
 
1047
        {
 
1048
                printf("Bary:%f,%f,%f \n",bary.X(),bary.Y(),bary.Z());
 
1049
                printf("texCoord:%f,%f \n",U,V);
 
1050
                assert(0);
 
1051
        }*/
 
1052
        //assert ((U>=-1)&&(U<=1));
 
1053
        //assert ((V>=-1)&&(V<=1));
 
1054
}
 
1055
 
 
1056
template <class MeshType>
 
1057
bool GetBaryFaceFromUV(const MeshType &m,
 
1058
                                       const typename MeshType::ScalarType &U,
 
1059
                                           const typename MeshType::ScalarType &V,
 
1060
                                       typename MeshType::CoordType &bary,
 
1061
                                       int &index)
 
1062
{
 
1063
        typedef typename MeshType::ScalarType ScalarType;
 
1064
        const ScalarType _EPSILON = ScalarType(0.0000001);
 
1065
        /*assert ((U>=-1)&&(U<=1));
 
1066
        assert ((V>=-1)&&(V<=1));*/
 
1067
        for (unsigned int i=0;i<m.face.size();i++)
 
1068
        {
 
1069
                const typename  MeshType::FaceType *f=&m.face[i];
 
1070
                vcg::Point2<ScalarType> tex0=vcg::Point2<ScalarType>(f->V(0)->T().U(),f->V(0)->T().V());
 
1071
                vcg::Point2<ScalarType> tex1=vcg::Point2<ScalarType>(f->V(1)->T().U(),f->V(1)->T().V());
 
1072
                vcg::Point2<ScalarType> tex2=vcg::Point2<ScalarType>(f->V(2)->T().U(),f->V(2)->T().V());
 
1073
 
 
1074
                vcg::Point2<ScalarType> test=vcg::Point2<ScalarType>(U,V);
 
1075
                vcg::Triangle2<ScalarType> t2d=vcg::Triangle2<ScalarType>(tex0,tex1,tex2);
 
1076
                ScalarType area=(tex1-tex0)^(tex2-tex0);
 
1077
                //assert(area>-_EPSILON);
 
1078
                ///then find if the point 2d falls inside
 
1079
                if ((area>_EPSILON)&&(t2d.InterpolationParameters(test,bary.X(),bary.Y(),bary.Z())))
 
1080
                {
 
1081
                        index=i;
 
1082
                        ///approximation errors
 
1083
                        ScalarType sum=0;
 
1084
                        for (int x=0;x<3;x++)
 
1085
                        {
 
1086
                                if (((bary[x])<=0)&&(bary[x]>=-_EPSILON))
 
1087
                                        bary[x]=0;
 
1088
                                else
 
1089
                                if (((bary[x])>=1)&&(bary[x]<=1+_EPSILON))
 
1090
                                        bary[x]=1;
 
1091
                                sum+=bary[x];
 
1092
                        }
 
1093
                        if (sum==0)
 
1094
                                printf("error SUM %f \n",sum);
 
1095
 
 
1096
                        bary/=sum;
 
1097
                        return true;
 
1098
                }
 
1099
        }
 
1100
 
 
1101
        return (false);
 
1102
}
 
1103
 
 
1104
template <class FaceType>
 
1105
bool GetBaryFaceFromUV(std::vector<FaceType*> faces,
 
1106
                                       const typename FaceType::ScalarType &U,
 
1107
                                           const typename FaceType::ScalarType &V,
 
1108
                                       typename FaceType::CoordType &bary,
 
1109
                                       int &index)
 
1110
{
 
1111
        typedef typename FaceType::CoordType CoordType;
 
1112
        typedef typename FaceType::ScalarType ScalarType;
 
1113
        typedef typename FaceType::VertexType VertexType;
 
1114
        typedef typename FaceType::ScalarType ScalarType;
 
1115
        const ScalarType _EPSILON = ScalarType(0.0000001);
 
1116
        /*assert ((U>=-1)&&(U<=1));
 
1117
        assert ((V>=-1)&&(V<=1));*/
 
1118
        for (unsigned int i=0;i<faces.size();i++)
 
1119
        {
 
1120
                FaceType *f=faces[i];
 
1121
                vcg::Point2<ScalarType> tex0=vcg::Point2<ScalarType>(f->V(0)->T().U(),f->V(0)->T().V());
 
1122
                vcg::Point2<ScalarType> tex1=vcg::Point2<ScalarType>(f->V(1)->T().U(),f->V(1)->T().V());
 
1123
                vcg::Point2<ScalarType> tex2=vcg::Point2<ScalarType>(f->V(2)->T().U(),f->V(2)->T().V());
 
1124
 
 
1125
                vcg::Point2<ScalarType> test=vcg::Point2<ScalarType>(U,V);
 
1126
                vcg::Triangle2<ScalarType> t2d=vcg::Triangle2<ScalarType>(tex0,tex1,tex2);
 
1127
                ScalarType area=fabs((tex1-tex0)^(tex2-tex0));
 
1128
                //assert(area>-_EPSILON);
 
1129
                ///then find if the point 2d falls inside
 
1130
                if ((area>_EPSILON)&&(t2d.InterpolationParameters(test,bary.X(),bary.Y(),bary.Z())))
 
1131
                {
 
1132
                        index=i;
 
1133
                        
 
1134
                        ///approximation errors
 
1135
                        ScalarType sum=0;
 
1136
                        for (int x=0;x<3;x++)
 
1137
                        {
 
1138
                                if (((bary[x])<=0)&&(bary[x]>=-_EPSILON))
 
1139
                                        bary[x]=0;
 
1140
                                else
 
1141
                                if (((bary[x])>=1)&&(bary[x]<=1+_EPSILON))
 
1142
                                        bary[x]=1;
 
1143
                                sum+=fabs(bary[x]);
 
1144
                        }
 
1145
                        if (sum==0)
 
1146
                                printf("error SUM %f \n",sum);
 
1147
                        
 
1148
                        bary/=sum;
 
1149
                        /*if (!((bary.X()>=0)&& (bary.X()<=1)))
 
1150
                                printf("error %f \n",bary.X());*/
 
1151
                        /*ScalarType diff=(1.0-bary.X()-bary.Y()-bary.Z());
 
1152
                        bary.X()+=diff;*/
 
1153
                        return true;
 
1154
                }
 
1155
        }
 
1156
 
 
1157
        return (false);
 
1158
}
 
1159
 
 
1160
template <class MeshType>
 
1161
bool GetBaryFaceFromUV(const MeshType &m,
 
1162
                                       const typename MeshType::ScalarType &U,
 
1163
                                           const typename MeshType::ScalarType &V,
 
1164
                                       const std::vector<typename MeshType::FaceType*> &orderedFaces,
 
1165
                                       typename MeshType::CoordType &bary,
 
1166
                                       typename MeshType::FaceType* &chosen)
 
1167
{
 
1168
        int index;
 
1169
        bool found=GetBaryFaceFromUV(m,U,V,bary,index);
 
1170
        if(!found) {
 
1171
                        chosen=0;
 
1172
                        return false;
 
1173
                }
 
1174
        chosen=orderedFaces[index];
 
1175
        return true;
 
1176
}
 
1177
 
 
1178
template <class MeshType>
 
1179
bool GetCoordFromUV(const MeshType &m,
 
1180
                                    const typename MeshType::ScalarType &U,
 
1181
                                        const typename MeshType::ScalarType &V,
 
1182
                                    typename MeshType::CoordType &val,
 
1183
                                        bool rpos=false)
 
1184
{
 
1185
        typedef typename MeshType::ScalarType ScalarType;
 
1186
        const ScalarType _EPSILON = (ScalarType)0.00001;
 
1187
        for (unsigned int i=0;i<m.face.size();i++)
 
1188
        {
 
1189
                const typename MeshType::FaceType *f=&m.face[i];
 
1190
                vcg::Point2<ScalarType> tex0=vcg::Point2<ScalarType>(f->V(0)->T().U(),f->V(0)->T().V());
 
1191
                vcg::Point2<ScalarType> tex1=vcg::Point2<ScalarType>(f->V(1)->T().U(),f->V(1)->T().V());
 
1192
                vcg::Point2<ScalarType> tex2=vcg::Point2<ScalarType>(f->V(2)->T().U(),f->V(2)->T().V());
 
1193
 
 
1194
                vcg::Point2<ScalarType> test=vcg::Point2<ScalarType>(U,V);
 
1195
                vcg::Triangle2<ScalarType> t2d=vcg::Triangle2<ScalarType>(tex0,tex1,tex2);
 
1196
                ScalarType area=(tex1-tex0)^(tex2-tex0);
 
1197
                ///then find if the point 2d falls inside
 
1198
                typename MeshType::CoordType bary;
 
1199
                if ((area>_EPSILON)&&(t2d.InterpolationParameters(test,bary.X(),bary.Y(),bary.Z())))
 
1200
                {
 
1201
                        ///approximation errors
 
1202
                        for (int x=0;x<3;x++)
 
1203
                        {
 
1204
                                if (((bary[x])<=0)&&(bary[x]>=-_EPSILON))
 
1205
                                        bary[x]=0;
 
1206
                                else
 
1207
                                if (((bary[x])>=1)&&(bary[x]<=1+_EPSILON))
 
1208
                                        bary[x]=1;
 
1209
                        }       
 
1210
                        ScalarType diff=(1.0-bary.X()-bary.Y()-bary.Z());
 
1211
                        bary.X()+=diff;
 
1212
                        if (!rpos)
 
1213
                                val=f->P(0)*bary.X()+f->P(1)*bary.Y()+f->P(0)*bary.Z();
 
1214
                        else
 
1215
                                val=f->V(0)->RPos*bary.X()+f->V(1)->RPos*bary.Y()+f->V(2)->RPos*bary.Z();
 
1216
                        
 
1217
                        return true;
 
1218
                }
 
1219
        }
 
1220
        return false;
 
1221
}
 
1222
 
 
1223
template <class MeshType>
 
1224
typename MeshType::ScalarType GetSmallestUVEdgeSize(const MeshType &m)
 
1225
{
 
1226
        typedef typename MeshType::ScalarType ScalarType;
 
1227
 
 
1228
        ScalarType smallest=100.f;
 
1229
        assert(m.fn>0);
 
1230
        for (int i=0;i<m.face.size();i++)
 
1231
        {
 
1232
 
 
1233
                ///approximation errors
 
1234
                        for (int j=0;j<3;j++)
 
1235
                        {
 
1236
 
 
1237
                                vcg::Point2<ScalarType> uv0=m.face[i].V(j)->T().P();
 
1238
                                vcg::Point2<ScalarType> uv1=m.face[i].V((j+1)%3)->T().P();
 
1239
                                ScalarType test=(uv0-uv1).Norm();
 
1240
                                if (test<smallest)
 
1241
                                        smallest=test;
 
1242
                        }
 
1243
        }
 
1244
        return smallest;
 
1245
}
 
1246
 
 
1247
template <class MeshType>
 
1248
typename MeshType::ScalarType GetSmallestUVHeight(const MeshType &m)
 
1249
{
 
1250
        typedef typename MeshType::ScalarType ScalarType;
 
1251
        ScalarType smallest=(ScalarType)100.0;
 
1252
        ScalarType eps=(ScalarType)0.0001;
 
1253
        assert(m.fn>0);
 
1254
        for (unsigned int i=0;i<m.face.size();i++)
 
1255
        {
 
1256
                const typename MeshType::FaceType *f=&m.face[i];
 
1257
                ///approximation errors
 
1258
                for (int j=0;j<3;j++)
 
1259
                {
 
1260
                        vcg::Point2<ScalarType> uv0=f->V(j)->cT().P();
 
1261
                        vcg::Point2<ScalarType> uv1=f->V1(j)->cT().P();
 
1262
                        vcg::Point2<ScalarType> uv2=f->V2(j)->cT().P();
 
1263
                        ScalarType area=fabs((uv1-uv0)^(uv2-uv0));
 
1264
                        ScalarType base=(uv1-uv2).Norm();
 
1265
                        ScalarType h_test=area/base;
 
1266
                        if (h_test<smallest)
 
1267
                                        smallest=h_test;
 
1268
                        }
 
1269
        }
 
1270
        if (smallest<eps)
 
1271
                smallest=(ScalarType)eps;
 
1272
        if (smallest>(ScalarType)0.2)
 
1273
                smallest=(ScalarType)0.2;
 
1274
        return smallest;
 
1275
}
 
1276
 
 
1277
template <class MeshType>
 
1278
void ParametrizeStarEquilateral(typename MeshType::VertexType *center,
 
1279
                                                                bool /*subvertices=true*/)
 
1280
{
 
1281
        ///initialize domain
 
1282
        typedef typename MeshType::VertexType VertexType;
 
1283
        typedef typename MeshType::FaceType FaceType;
 
1284
        typedef typename MeshType::CoordType CoordType;
 
1285
 
 
1286
        MeshType parametrized;
 
1287
        std::vector<VertexType*> vertices,ordVert;
 
1288
        std::vector<VertexType*> HresVert;
 
1289
        std::vector<FaceType*> faces;
 
1290
        vertices.push_back(center);
 
1291
        getSharedFace<MeshType>(vertices,faces);
 
1292
        CopyMeshFromFaces<MeshType>(faces,ordVert,parametrized);
 
1293
 
 
1294
        ///parametrize and then copy back
 
1295
        ParametrizeStarEquilateral<MeshType>(parametrized);
 
1296
        for (unsigned int i=0;i<ordVert.size();i++)
 
1297
                ordVert[i]->T().P()=parametrized.vert[i].T().P();
 
1298
 
 
1299
        ///initialize sub-vertices
 
1300
        getHresVertex<FaceType>(faces,HresVert);
 
1301
        for (unsigned int i=0;i<HresVert.size();i++)
 
1302
        {
 
1303
                FaceType *father=HresVert[i]->father;
 
1304
                CoordType Bary=HresVert[i]->Bary;
 
1305
                GetUV<MeshType>(father,Bary,HresVert[i]->T().U(),HresVert[i]->T().V());
 
1306
        }
 
1307
}
 
1308
 
 
1309
#endif