manolas
/
vcglib
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

added VertexBorder sampling algorithm that simply collect all the vertexes on the boundary.

This commit is contained in:
Paolo Cignoni 2013-12-20 02:27:09 +00:00
parent 2acd02f102
commit d1a5d53a89
1 changed files with 314 additions and 304 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

618
vcg/complex/algorithms/point_sampling.h
View File

@ -185,14 +185,14 @@ public:
static math::MarsenneTwisterRNG &SamplingRandomGenerator() static math::MarsenneTwisterRNG &SamplingRandomGenerator()
{ {
static math::MarsenneTwisterRNG rnd; static math::MarsenneTwisterRNG rnd;
return rnd; return rnd;
} }
// Returns an integer random number in the [0,i-1] interval using the improve Marsenne-Twister method. // Returns an integer random number in the [0,i-1] interval using the improve Marsenne-Twister method.
static unsigned int RandomInt(unsigned int i) static unsigned int RandomInt(unsigned int i)
{ {
return (SamplingRandomGenerator().generate(0) % i); return (SamplingRandomGenerator().generate(0) % i);
} }
// Returns a random number in the [0,1) real interval using the improved Marsenne-Twister method. // Returns a random number in the [0,1) real interval using the improved Marsenne-Twister method.
@ -319,14 +319,14 @@ static int Poisson(double lambda)
static void AllVertex(MetroMesh & m, VertexSampler &ps) static void AllVertex(MetroMesh & m, VertexSampler &ps)
{ {
VertexIterator vi; VertexIterator vi;
for(vi=m.vert.begin();vi!=m.vert.end();++vi) for(vi=m.vert.begin();vi!=m.vert.end();++vi)
{ {
if(!(*vi).IsD()) if(!(*vi).IsD())
{ {
ps.AddVert(*vi); ps.AddVert(*vi);
} }
} }
} }
/// Sample the vertices in a weighted way. Each vertex has a probability of being chosen /// Sample the vertices in a weighted way. Each vertex has a probability of being chosen
@ -339,103 +339,103 @@ static void AllVertex(MetroMesh & m, VertexSampler &ps)
static void VertexWeighted(MetroMesh & m, VertexSampler &ps, int sampleNum) static void VertexWeighted(MetroMesh & m, VertexSampler &ps, int sampleNum)
{ {
ScalarType qSum = 0; ScalarType qSum = 0;
VertexIterator vi; VertexIterator vi;
for(vi = m.vert.begin(); vi != m.vert.end(); ++vi) for(vi = m.vert.begin(); vi != m.vert.end(); ++vi)
if(!(*vi).IsD()) if(!(*vi).IsD())
qSum += (*vi).Q(); qSum += (*vi).Q();
ScalarType samplePerUnit = sampleNum/qSum; ScalarType samplePerUnit = sampleNum/qSum;
ScalarType floatSampleNum =0; ScalarType floatSampleNum =0;
std::vector<VertexPointer> vertVec; std::vector<VertexPointer> vertVec;
FillAndShuffleVertexPointerVector(m,vertVec); FillAndShuffleVertexPointerVector(m,vertVec);
std::vector<bool> vertUsed(m.vn,false); std::vector<bool> vertUsed(m.vn,false);
int i=0; int cnt=0; int i=0; int cnt=0;
while(cnt < sampleNum) while(cnt < sampleNum)
{ {
if(vertUsed[i]) if(vertUsed[i])
{ {
floatSampleNum += vertVec[i]->Q() * samplePerUnit; floatSampleNum += vertVec[i]->Q() * samplePerUnit;
int vertSampleNum = (int) floatSampleNum; int vertSampleNum = (int) floatSampleNum;
floatSampleNum -= (float) vertSampleNum; floatSampleNum -= (float) vertSampleNum;
// for every sample p_i in T... // for every sample p_i in T...
if(vertSampleNum > 1) if(vertSampleNum > 1)
{ {
ps.AddVert(*vertVec[i]); ps.AddVert(*vertVec[i]);
cnt++; cnt++;
vertUsed[i]=true; vertUsed[i]=true;
} }
} }
i = (i+1)%m.vn; i = (i+1)%m.vn;
} }
} }
/// Sample the vertices in a uniform way. Each vertex has a probability of being chosen /// Sample the vertices in a uniform way. Each vertex has a probability of being chosen
/// that is proportional to the area it represent. /// that is proportional to the area it represent.
static void VertexAreaUniform(MetroMesh & m, VertexSampler &ps, int sampleNum) static void VertexAreaUniform(MetroMesh & m, VertexSampler &ps, int sampleNum)
{ {
VertexIterator vi; VertexIterator vi;
for(vi = m.vert.begin(); vi != m.vert.end(); ++vi) for(vi = m.vert.begin(); vi != m.vert.end(); ++vi)
if(!(*vi).IsD()) if(!(*vi).IsD())
(*vi).Q() = 0; (*vi).Q() = 0;
FaceIterator fi; FaceIterator fi;
for(fi = m.face.begin(); fi != m.face.end(); ++fi) for(fi = m.face.begin(); fi != m.face.end(); ++fi)
if(!(*fi).IsD()) if(!(*fi).IsD())
{ {
ScalarType areaThird = DoubleArea(*fi)/6.0; ScalarType areaThird = DoubleArea(*fi)/6.0;
(*fi).V(0)->Q()+=areaThird; (*fi).V(0)->Q()+=areaThird;
(*fi).V(1)->Q()+=areaThird; (*fi).V(1)->Q()+=areaThird;
(*fi).V(2)->Q()+=areaThird; (*fi).V(2)->Q()+=areaThird;
} }
VertexWeighted(m,ps,sampleNum); VertexWeighted(m,ps,sampleNum);
} }
static void FillAndShuffleFacePointerVector(MetroMesh & m, std::vector<FacePointer> &faceVec) static void FillAndShuffleFacePointerVector(MetroMesh & m, std::vector<FacePointer> &faceVec)
{ {
FaceIterator fi; FaceIterator fi;
for(fi=m.face.begin();fi!=m.face.end();++fi) for(fi=m.face.begin();fi!=m.face.end();++fi)
if(!(*fi).IsD()) faceVec.push_back(&*fi); if(!(*fi).IsD()) faceVec.push_back(&*fi);
assert((int)faceVec.size()==m.fn); assert((int)faceVec.size()==m.fn);
unsigned int (*p_myrandom)(unsigned int) = RandomInt; unsigned int (*p_myrandom)(unsigned int) = RandomInt;
std::random_shuffle(faceVec.begin(),faceVec.end(), p_myrandom); std::random_shuffle(faceVec.begin(),faceVec.end(), p_myrandom);
} }
static void FillAndShuffleVertexPointerVector(MetroMesh & m, std::vector<VertexPointer> &vertVec) static void FillAndShuffleVertexPointerVector(MetroMesh & m, std::vector<VertexPointer> &vertVec)
{ {
VertexIterator vi; VertexIterator vi;
for(vi=m.vert.begin();vi!=m.vert.end();++vi) for(vi=m.vert.begin();vi!=m.vert.end();++vi)
if(!(*vi).IsD()) vertVec.push_back(&*vi); if(!(*vi).IsD()) vertVec.push_back(&*vi);
assert((int)vertVec.size()==m.vn); assert((int)vertVec.size()==m.vn);
unsigned int (*p_myrandom)(unsigned int) = RandomInt; unsigned int (*p_myrandom)(unsigned int) = RandomInt;
std::random_shuffle(vertVec.begin(),vertVec.end(), p_myrandom); std::random_shuffle(vertVec.begin(),vertVec.end(), p_myrandom);
} }
/// Sample the vertices in a uniform way. Each vertex has the same probabiltiy of being chosen. /// Sample the vertices in a uniform way. Each vertex has the same probabiltiy of being chosen.
static void VertexUniform(MetroMesh & m, VertexSampler &ps, int sampleNum) static void VertexUniform(MetroMesh & m, VertexSampler &ps, int sampleNum)
{ {
if(sampleNum>=m.vn) { if(sampleNum>=m.vn) {
AllVertex(m,ps); AllVertex(m,ps);
return; return;
} }
std::vector<VertexPointer> vertVec; std::vector<VertexPointer> vertVec;
FillAndShuffleVertexPointerVector(m,vertVec); FillAndShuffleVertexPointerVector(m,vertVec);
for(int i =0; i< sampleNum; ++i) for(int i =0; i< sampleNum; ++i)
ps.AddVert(*vertVec[i]); ps.AddVert(*vertVec[i]);
} }
/// \brief Sample all the border corner vertices /// \brief Sample all the border corner vertices
/// ///
/// It assumes that the border flag have been set over the mesh. /// It assumes that the border flag have been set over the mesh both for vertex and for faces.
/// All the vertices on the border where the surface forms an angle smaller than the given threshold are sampled. /// All the vertices on the border where the surface forms an angle smaller than the given threshold are sampled.
/// ///
static void VertexBorderCorner(MetroMesh & m, VertexSampler &ps, float angleRad) static void VertexBorderCorner(MetroMesh & m, VertexSampler &ps, float angleRad)
@ -455,13 +455,23 @@ static void VertexBorderCorner(MetroMesh & m, VertexSampler &ps, float angleRad)
for(VertexIterator vi=m.vert.begin();vi!=m.vert.end();++vi) for(VertexIterator vi=m.vert.begin();vi!=m.vert.end();++vi)
{ {
if(angleSumH[vi]<angleRad) if(angleSumH[vi]<angleRad && vi->IsB())
ps.AddVert(*vi); ps.AddVert(*vi);
} }
tri::Allocator<MetroMesh>:: template DeletePerVertexAttribute<float> (m,angleSumH); tri::Allocator<MetroMesh>:: template DeletePerVertexAttribute<float> (m,angleSumH);
} }
/// \brief Sample all the border vertices
///
/// It assumes that the border flag have been set over the mesh.
/// All the vertices on the border are sampled.
///
static void VertexBorder(MetroMesh & m, VertexSampler &ps)
{
VertexBorderCorner(m,ps,std::numeric_limits<ScalarType>::max());
}
/// Sample all the crease vertices. /// Sample all the crease vertices.
/// It assumes that the crease edges had been marked as non-faux edges /// It assumes that the crease edges had been marked as non-faux edges
/// for example by using /// for example by using
@ -493,26 +503,26 @@ static void VertexCrease(MetroMesh & m, VertexSampler &ps)
static void FaceUniform(MetroMesh & m, VertexSampler &ps, int sampleNum) static void FaceUniform(MetroMesh & m, VertexSampler &ps, int sampleNum)
{ {
if(sampleNum>=m.fn) { if(sampleNum>=m.fn) {
AllFace(m,ps); AllFace(m,ps);
return; return;
} }
std::vector<FacePointer> faceVec; std::vector<FacePointer> faceVec;
FillAndShuffleFacePointerVector(m,faceVec); FillAndShuffleFacePointerVector(m,faceVec);
for(int i =0; i< sampleNum; ++i) for(int i =0; i< sampleNum; ++i)
ps.AddFace(*faceVec[i],Barycenter(*faceVec[i])); ps.AddFace(*faceVec[i],Barycenter(*faceVec[i]));
} }
static void AllFace(MetroMesh & m, VertexSampler &ps) static void AllFace(MetroMesh & m, VertexSampler &ps)
{ {
FaceIterator fi; FaceIterator fi;
for(fi=m.face.begin();fi!=m.face.end();++fi) for(fi=m.face.begin();fi!=m.face.end();++fi)
if(!(*fi).IsD()) if(!(*fi).IsD())
{ {
ps.AddFace(*fi,Barycenter(*fi)); ps.AddFace(*fi,Barycenter(*fi));
} }
} }
@ -534,31 +544,31 @@ static void AllEdge(MetroMesh & m, VertexSampler &ps)
static void EdgeUniform(MetroMesh & m, VertexSampler &ps,int sampleNum, bool sampleFauxEdge=true) static void EdgeUniform(MetroMesh & m, VertexSampler &ps,int sampleNum, bool sampleFauxEdge=true)
{ {
typedef typename UpdateTopology<MetroMesh>::PEdge SimpleEdge; typedef typename UpdateTopology<MetroMesh>::PEdge SimpleEdge;
std::vector< SimpleEdge > Edges; std::vector< SimpleEdge > Edges;
UpdateTopology<MetroMesh>::FillUniqueEdgeVector(m,Edges,sampleFauxEdge); UpdateTopology<MetroMesh>::FillUniqueEdgeVector(m,Edges,sampleFauxEdge);
// First loop compute total edge length; // First loop compute total edge length;
float edgeSum=0; float edgeSum=0;
typename std::vector< SimpleEdge >::iterator ei; typename std::vector< SimpleEdge >::iterator ei;
for(ei=Edges.begin(); ei!=Edges.end(); ++ei) for(ei=Edges.begin(); ei!=Edges.end(); ++ei)
edgeSum+=Distance((*ei).v[0]->P(),(*ei).v[1]->P()); edgeSum+=Distance((*ei).v[0]->P(),(*ei).v[1]->P());
float sampleLen = edgeSum/sampleNum; float sampleLen = edgeSum/sampleNum;
float rest=0; float rest=0;
for(ei=Edges.begin(); ei!=Edges.end(); ++ei) for(ei=Edges.begin(); ei!=Edges.end(); ++ei)
{ {
float len = Distance((*ei).v[0]->P(),(*ei).v[1]->P()); float len = Distance((*ei).v[0]->P(),(*ei).v[1]->P());
float samplePerEdge = floor((len+rest)/sampleLen); float samplePerEdge = floor((len+rest)/sampleLen);
rest = (len+rest) - samplePerEdge * sampleLen; rest = (len+rest) - samplePerEdge * sampleLen;
float step = 1.0/(samplePerEdge+1); float step = 1.0/(samplePerEdge+1);
for(int i=0;i<samplePerEdge;++i) for(int i=0;i<samplePerEdge;++i)
{ {
Point3f interp(0,0,0); Point3f interp(0,0,0);
interp[ (*ei).z ]=step*(i+1); interp[ (*ei).z ]=step*(i+1);
interp[((*ei).z+1)%3]=1.0-step*(i+1); interp[((*ei).z+1)%3]=1.0-step*(i+1);
ps.AddFace(*(*ei).f,interp); ps.AddFace(*(*ei).f,interp);
} }
} }
} }
// Generate the barycentric coords of a random point over a single face, // Generate the barycentric coords of a random point over a single face,
@ -572,30 +582,30 @@ static CoordType RandomBarycentric()
// Given a triangle return a random point over it // Given a triangle return a random point over it
static CoordType RandomPointInTriangle(const FaceType &f) static CoordType RandomPointInTriangle(const FaceType &f)
{ {
CoordType u = RandomBarycentric(); CoordType u = RandomBarycentric();
return f.cP(0)*u[0] + f.cP(1)*u[1] + f.cP(2)*u[2]; return f.cP(0)*u[0] + f.cP(1)*u[1] + f.cP(2)*u[2];
} }
static void StratifiedMontecarlo(MetroMesh & m, VertexSampler &ps,int sampleNum) static void StratifiedMontecarlo(MetroMesh & m, VertexSampler &ps,int sampleNum)
{ {
ScalarType area = Stat<MetroMesh>::ComputeMeshArea(m); ScalarType area = Stat<MetroMesh>::ComputeMeshArea(m);
ScalarType samplePerAreaUnit = sampleNum/area; ScalarType samplePerAreaUnit = sampleNum/area;
// Montecarlo sampling. // Montecarlo sampling.
double floatSampleNum = 0.0; double floatSampleNum = 0.0;
FaceIterator fi; FaceIterator fi;
for(fi=m.face.begin(); fi != m.face.end(); fi++) for(fi=m.face.begin(); fi != m.face.end(); fi++)
if(!(*fi).IsD()) if(!(*fi).IsD())
{ {
// compute # samples in the current face (taking into account of the remainders) // compute # samples in the current face (taking into account of the remainders)
floatSampleNum += 0.5*DoubleArea(*fi) * samplePerAreaUnit; floatSampleNum += 0.5*DoubleArea(*fi) * samplePerAreaUnit;
int faceSampleNum = (int) floatSampleNum; int faceSampleNum = (int) floatSampleNum;
// for every sample p_i in T... // for every sample p_i in T...
for(int i=0; i < faceSampleNum; i++) for(int i=0; i < faceSampleNum; i++)
ps.AddFace(*fi,RandomBarycentric()); ps.AddFace(*fi,RandomBarycentric());
floatSampleNum -= (double) faceSampleNum; floatSampleNum -= (double) faceSampleNum;
} }
} }
/** /**
@ -681,37 +691,37 @@ static void EdgeMontecarlo(MetroMesh & m, VertexSampler &ps, int sampleNum, bool
static void Montecarlo(MetroMesh & m, VertexSampler &ps,int sampleNum) static void Montecarlo(MetroMesh & m, VertexSampler &ps,int sampleNum)
{ {
typedef std::pair<ScalarType, FacePointer> IntervalType; typedef std::pair<ScalarType, FacePointer> IntervalType;
std::vector< IntervalType > intervals (m.fn+1); std::vector< IntervalType > intervals (m.fn+1);
FaceIterator fi; FaceIterator fi;
int i=0; int i=0;
intervals[i]=std::make_pair(0,FacePointer(0)); intervals[i]=std::make_pair(0,FacePointer(0));
// First loop: build a sequence of consecutive segments proportional to the triangle areas. // First loop: build a sequence of consecutive segments proportional to the triangle areas.
for(fi=m.face.begin(); fi != m.face.end(); fi++) for(fi=m.face.begin(); fi != m.face.end(); fi++)
if(!(*fi).IsD()) if(!(*fi).IsD())
{ {
intervals[i+1]=std::make_pair(intervals[i].first+0.5*DoubleArea(*fi), &*fi); intervals[i+1]=std::make_pair(intervals[i].first+0.5*DoubleArea(*fi), &*fi);
++i; ++i;
} }
ScalarType meshArea = intervals.back().first; ScalarType meshArea = intervals.back().first;
for(i=0;i<sampleNum;++i) for(i=0;i<sampleNum;++i)
{ {
ScalarType val = meshArea * RandomDouble01(); ScalarType val = meshArea * RandomDouble01();
// lower_bound returns the furthermost iterator i in [first, last) such that, for every iterator j in [first, i), *j < value. // lower_bound returns the furthermost iterator i in [first, last) such that, for every iterator j in [first, i), *j < value.
// E.g. An iterator pointing to the first element "not less than" val, or end() if every element is less than val. // E.g. An iterator pointing to the first element "not less than" val, or end() if every element is less than val.
typename std::vector<IntervalType>::iterator it = lower_bound(intervals.begin(),intervals.end(),std::make_pair(val,FacePointer(0)) ); typename std::vector<IntervalType>::iterator it = lower_bound(intervals.begin(),intervals.end(),std::make_pair(val,FacePointer(0)) );
assert(it != intervals.end()); assert(it != intervals.end());
assert(it != intervals.begin()); assert(it != intervals.begin());
assert( (*(it-1)).first <val ); assert( (*(it-1)).first <val );
assert( (*(it)).first >= val); assert( (*(it)).first >= val);
ps.AddFace( *(*it).second, RandomBarycentric() ); ps.AddFace( *(*it).second, RandomBarycentric() );
} }
} }
static ScalarType WeightedArea(FaceType f) static ScalarType WeightedArea(FaceType f)
{ {
ScalarType averageQ = ( f.V(0)->Q() + f.V(1)->Q() + f.V(2)->Q() ) /3.0; ScalarType averageQ = ( f.V(0)->Q() + f.V(1)->Q() + f.V(2)->Q() ) /3.0;
return DoubleArea(f)*averageQ/2.0; return DoubleArea(f)*averageQ/2.0;
} }
/// Compute a sampling of the surface that is weighted by the quality /// Compute a sampling of the surface that is weighted by the quality
@ -719,27 +729,27 @@ static ScalarType WeightedArea(FaceType f)
/// So the a face with a zero quality on all its vertices is never sampled and a face with average quality 2 get twice the samples of a face with the same area but with an average quality of 1; /// So the a face with a zero quality on all its vertices is never sampled and a face with average quality 2 get twice the samples of a face with the same area but with an average quality of 1;
static void WeightedMontecarlo(MetroMesh & m, VertexSampler &ps, int sampleNum) static void WeightedMontecarlo(MetroMesh & m, VertexSampler &ps, int sampleNum)
{ {
assert(tri::HasPerVertexQuality(m)); assert(tri::HasPerVertexQuality(m));
ScalarType weightedArea = 0; ScalarType weightedArea = 0;
FaceIterator fi; FaceIterator fi;
for(fi = m.face.begin(); fi != m.face.end(); ++fi) for(fi = m.face.begin(); fi != m.face.end(); ++fi)
if(!(*fi).IsD()) if(!(*fi).IsD())
weightedArea += WeightedArea(*fi); weightedArea += WeightedArea(*fi);
ScalarType samplePerAreaUnit = sampleNum/weightedArea; ScalarType samplePerAreaUnit = sampleNum/weightedArea;
// Montecarlo sampling. // Montecarlo sampling.
double floatSampleNum = 0.0; double floatSampleNum = 0.0;
for(fi=m.face.begin(); fi != m.face.end(); fi++) for(fi=m.face.begin(); fi != m.face.end(); fi++)
if(!(*fi).IsD()) if(!(*fi).IsD())
{ {
// compute # samples in the current face (taking into account of the remainders) // compute # samples in the current face (taking into account of the remainders)
floatSampleNum += WeightedArea(*fi) * samplePerAreaUnit; floatSampleNum += WeightedArea(*fi) * samplePerAreaUnit;
int faceSampleNum = (int) floatSampleNum; int faceSampleNum = (int) floatSampleNum;
// for every sample p_i in T... // for every sample p_i in T...
for(int i=0; i < faceSampleNum; i++) for(int i=0; i < faceSampleNum; i++)
ps.AddFace(*fi,RandomBarycentric()); ps.AddFace(*fi,RandomBarycentric());
floatSampleNum -= (double) faceSampleNum; floatSampleNum -= (double) faceSampleNum;
} }
@ -812,27 +822,27 @@ static int SingleFaceSubdivision(int sampleNum, const CoordType & v0, const Coor
static void FaceSubdivision(MetroMesh & m, VertexSampler &ps,int sampleNum, bool randSample) static void FaceSubdivision(MetroMesh & m, VertexSampler &ps,int sampleNum, bool randSample)
{ {
ScalarType area = Stat<MetroMesh>::ComputeMeshArea(m); ScalarType area = Stat<MetroMesh>::ComputeMeshArea(m);
ScalarType samplePerAreaUnit = sampleNum/area; ScalarType samplePerAreaUnit = sampleNum/area;
std::vector<FacePointer> faceVec; std::vector<FacePointer> faceVec;
FillAndShuffleFacePointerVector(m,faceVec); FillAndShuffleFacePointerVector(m,faceVec);
vcg::tri::UpdateNormal<MetroMesh>::PerFaceNormalized(m); vcg::tri::UpdateNormal<MetroMesh>::PerFaceNormalized(m);
double floatSampleNum = 0.0; double floatSampleNum = 0.0;
int faceSampleNum; int faceSampleNum;
// Subdivision sampling. // Subdivision sampling.
typename std::vector<FacePointer>::iterator fi; typename std::vector<FacePointer>::iterator fi;
for(fi=faceVec.begin(); fi!=faceVec.end(); fi++) for(fi=faceVec.begin(); fi!=faceVec.end(); fi++)
{ {
const CoordType b0(1.0, 0.0, 0.0); const CoordType b0(1.0, 0.0, 0.0);
const CoordType b1(0.0, 1.0, 0.0); const CoordType b1(0.0, 1.0, 0.0);
const CoordType b2(0.0, 0.0, 1.0); const CoordType b2(0.0, 0.0, 1.0);
// compute # samples in the current face. // compute # samples in the current face.
floatSampleNum += 0.5*DoubleArea(**fi) * samplePerAreaUnit; floatSampleNum += 0.5*DoubleArea(**fi) * samplePerAreaUnit;
faceSampleNum = (int) floatSampleNum; faceSampleNum = (int) floatSampleNum;
if(faceSampleNum>0) if(faceSampleNum>0)
faceSampleNum = SingleFaceSubdivision(faceSampleNum,b0,b1,b2,ps,*fi,randSample); faceSampleNum = SingleFaceSubdivision(faceSampleNum,b0,b1,b2,ps,*fi,randSample);
floatSampleNum -= (double) faceSampleNum; floatSampleNum -= (double) faceSampleNum;
} }
} }
//--------- //---------
// Subdivision sampling of a single face. // Subdivision sampling of a single face.
@ -967,17 +977,17 @@ static int SingleFaceSimilarDual(FacePointer fp, VertexSampler &ps, int n_sample
ps.AddFace(*fp, V0*rb[0]+V1*rb[1]+V2*rb[2]); ps.AddFace(*fp, V0*rb[0]+V1*rb[1]+V2*rb[2]);
} else ps.AddFace(*fp,(V0+V1+V2)/3.0); } else ps.AddFace(*fp,(V0+V1+V2)/3.0);
if( j < n_samples_per_edge-i-2 ) if( j < n_samples_per_edge-i-2 )
{ {
CoordType V3((i+1)*segmentLen,(j+1)*segmentLen, 1.0 - ((i+1)*segmentLen+(j+1)*segmentLen) ) ; CoordType V3((i+1)*segmentLen,(j+1)*segmentLen, 1.0 - ((i+1)*segmentLen+(j+1)*segmentLen) ) ;
n_samples++; n_samples++;
if(randomFlag) { if(randomFlag) {
CoordType rb=RandomBarycentric(); CoordType rb=RandomBarycentric();
ps.AddFace(*fp, V3*rb[0]+V1*rb[1]+V2*rb[2]); ps.AddFace(*fp, V3*rb[0]+V1*rb[1]+V2*rb[2]);
} else ps.AddFace(*fp,(V3+V1+V2)/3.0); } else ps.AddFace(*fp,(V3+V1+V2)/3.0);
} }
} }
return n_samples; return n_samples;
} }
// Similar sampling // Similar sampling
@ -1011,8 +1021,8 @@ static int SingleFaceSimilarDual(FacePointer fp, VertexSampler &ps, int n_sample
//void Sampling<MetroMesh>::SimilarFaceSampling() //void Sampling<MetroMesh>::SimilarFaceSampling()
static void FaceSimilar(MetroMesh & m, VertexSampler &ps,int sampleNum, bool dualFlag, bool randomFlag) static void FaceSimilar(MetroMesh & m, VertexSampler &ps,int sampleNum, bool dualFlag, bool randomFlag)
{ {
ScalarType area = Stat<MetroMesh>::ComputeMeshArea(m); ScalarType area = Stat<MetroMesh>::ComputeMeshArea(m);
ScalarType samplePerAreaUnit = sampleNum/area; ScalarType samplePerAreaUnit = sampleNum/area;
// Similar Triangles sampling. // Similar Triangles sampling.
int n_samples_per_edge; int n_samples_per_edge;
@ -1041,11 +1051,11 @@ static void FaceSimilar(MetroMesh & m, VertexSampler &ps,int sampleNum, bool dua
} }
// Rasterization fuction // Rasterization fuction
// Take a triangle // Take a triangle
// T deve essere una classe funzionale che ha l'operatore () // T deve essere una classe funzionale che ha l'operatore ()
// con due parametri x,y di tipo S esempio: // con due parametri x,y di tipo S esempio:
// class Foo { public void operator()(int x, int y ) { ??? } }; // class Foo { public void operator()(int x, int y ) { ??? } };
// This function does rasterization with a safety buffer area, thus accounting some points actually outside triangle area // This function does rasterization with a safety buffer area, thus accounting some points actually outside triangle area
// The safety area samples are generated according to face flag BORDER which should be true for texture space border edges // The safety area samples are generated according to face flag BORDER which should be true for texture space border edges
@ -1064,10 +1074,10 @@ static void FaceSimilar(MetroMesh & m, VertexSampler &ps,int sampleNum, bool dua
bboxf.Add(v1); bboxf.Add(v1);
bboxf.Add(v2); bboxf.Add(v2);
bbox.min[0] = floor(bboxf.min[0]); bbox.min[0] = floor(bboxf.min[0]);
bbox.min[1] = floor(bboxf.min[1]); bbox.min[1] = floor(bboxf.min[1]);
bbox.max[0] = ceil(bboxf.max[0]); bbox.max[0] = ceil(bboxf.max[0]);
bbox.max[1] = ceil(bboxf.max[1]); bbox.max[1] = ceil(bboxf.max[1]);
// Calcolo versori degli spigoli // Calcolo versori degli spigoli
Point2<S> d10 = v1 - v0; Point2<S> d10 = v1 - v0;
@ -1191,13 +1201,13 @@ static void FaceSimilar(MetroMesh & m, VertexSampler &ps,int sampleNum, bool dua
// check the radius constrain // check the radius constrain
static bool checkPoissonDisk(SampleSHT & sht, const Point3<ScalarType> & p, ScalarType radius) static bool checkPoissonDisk(SampleSHT & sht, const Point3<ScalarType> & p, ScalarType radius)
{ {
// get the samples closest to the given one // get the samples closest to the given one
std::vector<VertexType*> closests; std::vector<VertexType*> closests;
typedef VertTmark<MetroMesh> MarkerVert; typedef VertTmark<MetroMesh> MarkerVert;
static MarkerVert mv; static MarkerVert mv;
Box3f bb(p-Point3f(radius,radius,radius),p+Point3f(radius,radius,radius)); Box3f bb(p-Point3f(radius,radius,radius),p+Point3f(radius,radius,radius));
GridGetInBox(sht, mv, bb, closests); GridGetInBox(sht, mv, bb, closests);
ScalarType r2 = radius*radius; ScalarType r2 = radius*radius;
for(int i=0; i<closests.size(); ++i) for(int i=0; i<closests.size(); ++i)
@ -1253,9 +1263,9 @@ struct PoissonDiskParam
// It always return a point. // It always return a point.
static VertexPointer getSampleFromCell(Point3i &cell, MontecarloSHT & samplepool) static VertexPointer getSampleFromCell(Point3i &cell, MontecarloSHT & samplepool)
{ {
MontecarloSHTIterator cellBegin, cellEnd; MontecarloSHTIterator cellBegin, cellEnd;
samplepool.Grid(cell, cellBegin, cellEnd); samplepool.Grid(cell, cellBegin, cellEnd);
return *cellBegin; return *cellBegin;
} }
// Given a cell of the grid it search the point that remove the minimum number of other samples // Given a cell of the grid it search the point that remove the minimum number of other samples
@ -1286,24 +1296,24 @@ static VertexPointer getBestPrecomputedMontecarloSample(Point3i &cell, Montecarl
static ScalarType ComputePoissonDiskRadius(MetroMesh &origMesh, int sampleNum) static ScalarType ComputePoissonDiskRadius(MetroMesh &origMesh, int sampleNum)
{ {
ScalarType meshArea = Stat<MetroMesh>::ComputeMeshArea(origMesh); ScalarType meshArea = Stat<MetroMesh>::ComputeMeshArea(origMesh);
// Manage approximately the PointCloud Case, use the half a area of the bbox. // Manage approximately the PointCloud Case, use the half a area of the bbox.
// TODO: If you had the radius a much better approximation could be done. // TODO: If you had the radius a much better approximation could be done.
if(meshArea ==0) if(meshArea ==0)
{ {
meshArea = (origMesh.bbox.DimX()*origMesh.bbox.DimY() + meshArea = (origMesh.bbox.DimX()*origMesh.bbox.DimY() +
origMesh.bbox.DimX()*origMesh.bbox.DimZ() + origMesh.bbox.DimX()*origMesh.bbox.DimZ() +
origMesh.bbox.DimY()*origMesh.bbox.DimZ()); origMesh.bbox.DimY()*origMesh.bbox.DimZ());
} }
ScalarType diskRadius = sqrt(meshArea / (0.7 * M_PI * sampleNum)); // 0.7 is a density factor ScalarType diskRadius = sqrt(meshArea / (0.7 * M_PI * sampleNum)); // 0.7 is a density factor
return diskRadius; return diskRadius;
} }
static int ComputePoissonSampleNum(MetroMesh &origMesh, ScalarType diskRadius) static int ComputePoissonSampleNum(MetroMesh &origMesh, ScalarType diskRadius)
{ {
ScalarType meshArea = Stat<MetroMesh>::ComputeMeshArea(origMesh); ScalarType meshArea = Stat<MetroMesh>::ComputeMeshArea(origMesh);
int sampleNum = meshArea / (diskRadius*diskRadius *M_PI *0.7) ; // 0.7 is a density factor int sampleNum = meshArea / (diskRadius*diskRadius *M_PI *0.7) ; // 0.7 is a density factor
return sampleNum; return sampleNum;
} }
// Note that this function actually CHANGE the quality of the montecarlo samples so that it represents the expected radius of each sample // Note that this function actually CHANGE the quality of the montecarlo samples so that it represents the expected radius of each sample
@ -1311,16 +1321,16 @@ static int ComputePoissonSampleNum(MetroMesh &origMesh, ScalarType diskRadius)
static void ComputePoissonSampleRadii(MetroMesh &sampleMesh, ScalarType diskRadius, ScalarType radiusVariance, bool invert) static void ComputePoissonSampleRadii(MetroMesh &sampleMesh, ScalarType diskRadius, ScalarType radiusVariance, bool invert)
{ {
VertexIterator vi; VertexIterator vi;
std::pair<float,float> minmax = tri::Stat<MetroMesh>::ComputePerVertexQualityMinMax( sampleMesh); std::pair<float,float> minmax = tri::Stat<MetroMesh>::ComputePerVertexQualityMinMax( sampleMesh);
float minRad = diskRadius / radiusVariance; float minRad = diskRadius / radiusVariance;
float maxRad = diskRadius * radiusVariance; float maxRad = diskRadius * radiusVariance;
float deltaQ = minmax.second-minmax.first; float deltaQ = minmax.second-minmax.first;
float deltaRad = maxRad-minRad; float deltaRad = maxRad-minRad;
for (vi = sampleMesh.vert.begin(); vi != sampleMesh.vert.end(); vi++) for (vi = sampleMesh.vert.begin(); vi != sampleMesh.vert.end(); vi++)
{ {
(*vi).Q() = minRad + deltaRad*((invert ? minmax.second - (*vi).Q() : (*vi).Q() - minmax.first )/deltaQ); (*vi).Q() = minRad + deltaRad*((invert ? minmax.second - (*vi).Q() : (*vi).Q() - minmax.first )/deltaQ);
} }
} }
// initialize spatial hash table for searching // initialize spatial hash table for searching
@ -1447,9 +1457,9 @@ static void HierarchicalPoissonDisk(MetroMesh &origMesh, VertexSampler &ps, Metr
// This radius implies that when we pick a sample in a cell all that cell will not be touched again. // This radius implies that when we pick a sample in a cell all that cell will not be touched again.
ScalarType cellsize = 2.0f* diskRadius / sqrt(3.0); ScalarType cellsize = 2.0f* diskRadius / sqrt(3.0);
// inflating // inflating
BoxType bb=origMesh.bbox; BoxType bb=origMesh.bbox;
bb.Offset(cellsize); bb.Offset(cellsize);
int sizeX = std::max(1.0f,bb.DimX() / cellsize); int sizeX = std::max(1.0f,bb.DimX() / cellsize);
int sizeY = std::max(1.0f,bb.DimY() / cellsize); int sizeY = std::max(1.0f,bb.DimY() / cellsize);
@ -1461,17 +1471,17 @@ static void HierarchicalPoissonDisk(MetroMesh &origMesh, VertexSampler &ps, Metr
checkSHT.InitEmpty(bb, gridsize); checkSHT.InitEmpty(bb, gridsize);
// sampling algorithm // sampling algorithm
// ------------------ // ------------------
// //
// - generate millions of samples using montecarlo algorithm // - generate millions of samples using montecarlo algorithm
// - extract a cell (C) from the active cell list (with probability proportional to cell's volume) // - extract a cell (C) from the active cell list (with probability proportional to cell's volume)
// - generate a sample inside C by choosing one of the contained pre-generated samples // - generate a sample inside C by choosing one of the contained pre-generated samples
// - if the sample violates the radius constrain discard it, and add the cell to the cells-to-subdivide list // - if the sample violates the radius constrain discard it, and add the cell to the cells-to-subdivide list
// - iterate until the active cell list is empty or a pre-defined number of subdivisions is reached // - iterate until the active cell list is empty or a pre-defined number of subdivisions is reached
// //
int level = 0; int level = 0;
// initialize spatial hash to index pre-generated samples // initialize spatial hash to index pre-generated samples
montecarloSHTVec[0].InitEmpty(bb, gridsize); montecarloSHTVec[0].InitEmpty(bb, gridsize);
@ -1484,9 +1494,9 @@ static void HierarchicalPoissonDisk(MetroMesh &origMesh, VertexSampler &ps, Metr
if(pp.adaptiveRadiusFlag) if(pp.adaptiveRadiusFlag)
ComputePoissonSampleRadii(montecarloMesh, diskRadius, pp.radiusVariance, pp.invertQuality); ComputePoissonSampleRadii(montecarloMesh, diskRadius, pp.radiusVariance, pp.invertQuality);
do do
{ {
MontecarloSHT &montecarloSHT = montecarloSHTVec[level]; MontecarloSHT &montecarloSHT = montecarloSHTVec[level];
if(level>0) if(level>0)
{// initialize spatial hash with the remaining points {// initialize spatial hash with the remaining points
@ -1500,15 +1510,15 @@ static void HierarchicalPoissonDisk(MetroMesh &origMesh, VertexSampler &ps, Metr
unsigned int (*p_myrandom)(unsigned int) = RandomInt; unsigned int (*p_myrandom)(unsigned int) = RandomInt;
std::random_shuffle(montecarloSHT.AllocatedCells.begin(),montecarloSHT.AllocatedCells.end(), p_myrandom); std::random_shuffle(montecarloSHT.AllocatedCells.begin(),montecarloSHT.AllocatedCells.end(), p_myrandom);
// generate a sample inside C by choosing one of the contained pre-generated samples // generate a sample inside C by choosing one of the contained pre-generated samples
////////////////////////////////////////////////////////////////////////////////////////// //////////////////////////////////////////////////////////////////////////////////////////
int removedCnt=montecarloSHT.hash_table.size(); int removedCnt=montecarloSHT.hash_table.size();
int addedCnt=checkSHT.hash_table.size(); int addedCnt=checkSHT.hash_table.size();
for (int i = 0; i < montecarloSHT.AllocatedCells.size(); i++) for (int i = 0; i < montecarloSHT.AllocatedCells.size(); i++)
{ {
for(int j=0;j<4;j++) for(int j=0;j<4;j++)
{ {
if( montecarloSHT.EmptyCell(montecarloSHT.AllocatedCells[i]) ) continue; if( montecarloSHT.EmptyCell(montecarloSHT.AllocatedCells[i]) ) continue;
// generate a sample chosen from the pre-generated one // generate a sample chosen from the pre-generated one
typename MontecarloSHT::HashIterator hi = montecarloSHT.hash_table.find(montecarloSHT.AllocatedCells[i]); typename MontecarloSHT::HashIterator hi = montecarloSHT.hash_table.find(montecarloSHT.AllocatedCells[i]);
@ -1536,9 +1546,9 @@ static void HierarchicalPoissonDisk(MetroMesh &origMesh, VertexSampler &ps, Metr
// increase grid resolution // increase grid resolution
gridsize *= 2; gridsize *= 2;
// //
level++; level++;
} while(level < 5); } while(level < 5);
} }
//template <class MetroMesh> //template <class MetroMesh>
@ -1556,16 +1566,16 @@ static void HierarchicalPoissonDisk(MetroMesh &origMesh, VertexSampler &ps, Metr
// vcg::tri::UpdateFlags<Mesh>::FaceBorderFromFF(m); // vcg::tri::UpdateFlags<Mesh>::FaceBorderFromFF(m);
static void Texture(MetroMesh & m, VertexSampler &ps, int textureWidth, int textureHeight, bool correctSafePointsBaryCoords=true) static void Texture(MetroMesh & m, VertexSampler &ps, int textureWidth, int textureHeight, bool correctSafePointsBaryCoords=true)
{ {
FaceIterator fi; FaceIterator fi;
printf("Similar Triangles face sampling\n"); printf("Similar Triangles face sampling\n");
for(fi=m.face.begin(); fi != m.face.end(); fi++) for(fi=m.face.begin(); fi != m.face.end(); fi++)
if (!fi->IsD()) if (!fi->IsD())
{ {
Point2f ti[3]; Point2f ti[3];
for(int i=0;i<3;++i) for(int i=0;i<3;++i)
ti[i]=Point2f((*fi).WT(i).U() * textureWidth - 0.5, (*fi).WT(i).V() * textureHeight - 0.5); ti[i]=Point2f((*fi).WT(i).U() * textureWidth - 0.5, (*fi).WT(i).V() * textureHeight - 0.5);
// - 0.5 constants are used to obtain correct texture mapping // - 0.5 constants are used to obtain correct texture mapping
SingleFaceRaster(*fi, ps, ti[0],ti[1],ti[2], correctSafePointsBaryCoords); SingleFaceRaster(*fi, ps, ti[0],ti[1],ti[2], correctSafePointsBaryCoords);
} }
@ -1584,33 +1594,33 @@ TriMeshGrid gM;
static void RegularRecursiveOffset(MetroMesh & m, std::vector<Point3f> &pvec, ScalarType offset, float minDiag) static void RegularRecursiveOffset(MetroMesh & m, std::vector<Point3f> &pvec, ScalarType offset, float minDiag)
{ {
Box3<ScalarType> bb=m.bbox; Box3<ScalarType> bb=m.bbox;
bb.Offset(offset*2.0); bb.Offset(offset*2.0);
RRParam rrp; RRParam rrp;
rrp.markerFunctor.SetMesh(&m); rrp.markerFunctor.SetMesh(&m);
rrp.gM.Set(m.face.begin(),m.face.end(),bb); rrp.gM.Set(m.face.begin(),m.face.end(),bb);
rrp.offset=offset; rrp.offset=offset;
rrp.minDiag=minDiag; rrp.minDiag=minDiag;
SubdivideAndSample(m, pvec, bb, rrp, bb.Diag()); SubdivideAndSample(m, pvec, bb, rrp, bb.Diag());
} }
static void SubdivideAndSample(MetroMesh & m, std::vector<Point3f> &pvec, const Box3<ScalarType> bb, RRParam &rrp, float curDiag) static void SubdivideAndSample(MetroMesh & m, std::vector<Point3f> &pvec, const Box3<ScalarType> bb, RRParam &rrp, float curDiag)
{ {
Point3f startPt = bb.Center(); Point3f startPt = bb.Center();
ScalarType dist; ScalarType dist;
// Compute mesh point nearest to bb center // Compute mesh point nearest to bb center
FaceType *nearestF=0; FaceType *nearestF=0;
float dist_upper_bound = curDiag+rrp.offset; float dist_upper_bound = curDiag+rrp.offset;
Point3f closestPt; Point3f closestPt;
vcg::face::PointDistanceBaseFunctor<ScalarType> PDistFunct; vcg::face::PointDistanceBaseFunctor<ScalarType> PDistFunct;
dist=dist_upper_bound; dist=dist_upper_bound;
nearestF = rrp.gM.GetClosest(PDistFunct,rrp.markerFunctor,startPt,dist_upper_bound,dist,closestPt); nearestF = rrp.gM.GetClosest(PDistFunct,rrp.markerFunctor,startPt,dist_upper_bound,dist,closestPt);
curDiag /=2; curDiag /=2;
if(dist < dist_upper_bound) if(dist < dist_upper_bound)
{ {