From 831e3f9036755b97c8ec00413a12042d37148827 Mon Sep 17 00:00:00 2001 From: cignoni Date: Thu, 4 Feb 2016 18:18:53 +0000 Subject: [PATCH] fine tuning of the voronoi scaffolding engine --- .../algorithms/voronoi_volume_sampling.h | 77 ++++++++++++------- 1 file changed, 49 insertions(+), 28 deletions(-) diff --git a/vcg/complex/algorithms/voronoi_volume_sampling.h b/vcg/complex/algorithms/voronoi_volume_sampling.h index 2eadd3b7..6add323a 100644 --- a/vcg/complex/algorithms/voronoi_volume_sampling.h +++ b/vcg/complex/algorithms/voronoi_volume_sampling.h @@ -50,8 +50,8 @@ public: typedef typename vcg::tri::TrivialWalker MyWalker; typedef typename vcg::tri::MarchingCubes MyMarchingCubes; - VoronoiVolumeSampling(MeshType &_baseMesh, MeshType &_seedMesh) - :surfTree(0),seedTree(0),baseMesh(_baseMesh),seedMesh(_seedMesh),cb(0),restrictedRelaxationFlag(false) + VoronoiVolumeSampling(MeshType &_baseMesh) + :surfTree(0),seedTree(0),baseMesh(_baseMesh),cb(0),restrictedRelaxationFlag(false) { tri::RequirePerFaceMark(baseMesh); tri::UpdateBounding::Box(baseMesh); @@ -69,12 +69,13 @@ public: vcg::face::PointDistanceBaseFunctor PDistFunct; MeshType &baseMesh; - MeshType &seedMesh; + MeshType seedMesh; MeshType poissonSurfaceMesh; ScalarType poissonRadiusSurface; MeshType montecarloVolumeMesh; // we use this mesh as volume evaluator MeshType seedDomainMesh; // where we choose the seeds (by default is the montecarlo volume mesh) vcg::CallBackPos *cb; + math::MarsenneTwisterRNG rng; bool restrictedRelaxationFlag; @@ -91,12 +92,14 @@ public: ScalarType meshArea = Stat::ComputeMeshArea(baseMesh); int MontecarloSurfSampleNum = 10 * meshArea / (poissonRadiusSurface*poissonRadiusSurface); tri::MeshSampler sampler(montecarloSurfaceMesh); + tri::SurfaceSampling >::SamplingRandomGenerator()=rng; tri::SurfaceSampling >::Montecarlo(baseMesh, sampler, MontecarloSurfSampleNum); montecarloSurfaceMesh.bbox = baseMesh.bbox; // we want the same bounding box poissonSurfaceMesh.Clear(); tri::MeshSampler mps(poissonSurfaceMesh); typename tri::SurfaceSampling >::PoissonDiskParam pp; pp.geodesicDistanceFlag=false; + tri::SurfaceSampling >::PoissonDiskPruning(mps, montecarloSurfaceMesh, poissonRadiusSurface,pp); vcg::tri::UpdateBounding::Box(poissonSurfaceMesh); @@ -144,7 +147,7 @@ ScalarType DistanceFromVoronoiSeed(CoordType p_point) { ScalarType squaredDist; unsigned int ind; - surfTree->doQueryClosest(p_point,ind,squaredDist); + seedTree->doQueryClosest(p_point,ind,squaredDist); return math::Sqrt(squaredDist); } @@ -173,7 +176,7 @@ ScalarType DistanceFromVoronoiFace(CoordType p_point) * returns: distance between the point P and the line R */ - ScalarType DistanceFromVoronoiEdge(CoordType p_point) +ScalarType DistanceFromVoronoiEdge(CoordType p_point) { seedTree->doQueryK(p_point,3,pq); @@ -196,6 +199,30 @@ ScalarType DistanceFromVoronoiFace(CoordType p_point) return closestDist; } +ScalarType DistanceFromVoronoiCorner(CoordType p_point) +{ + + seedTree->doQueryK(p_point,4,pq); + std::vector > closeSeedVec; + CoordType p0= this->seedMesh.vert[pq.getIndex(0)].P(); + CoordType p1= this->seedMesh.vert[pq.getIndex(1)].P(); + CoordType p2= this->seedMesh.vert[pq.getIndex(2)].P(); + CoordType p3= this->seedMesh.vert[pq.getIndex(3)].P(); + + Plane3 pl01; pl01.Init((p0+p1)/2.0f,p0-p1); + Plane3 pl02; pl02.Init((p0+p2)/2.0f,p0-p2); + Plane3 pl03; pl03.Init((p0+p3)/2.0f,p0-p3); + Line3 voroLine; + + // Calculating the line R that intersect the planes pl01 and pl02 + vcg::IntersectionPlanePlane(pl01,pl02,voroLine); + CoordType intersectionPt; + bool ret = vcg::IntersectionLinePlane(voroLine,pl03,intersectionPt); + + return vcg::Distance(p_point,intersectionPt); +} + + void BarycentricRelaxVoronoiSamples(int relaxStep) { bool changed=false; @@ -319,26 +346,20 @@ void QuadricRelaxVoronoiSamples(int relaxStep) * PruningPoisson: mesh of inside and surface points, it's the voronoi3d diagram * n_voxel: number of voxels for the greater side */ - void BuildScaffoldingMesh(MeshType &scaffoldingMesh, int volumeSide, ScalarType isoThr,int elemEnum, bool surfFlag) + void BuildScaffoldingMesh(MeshType &scaffoldingMesh, float voxelSide, ScalarType isoThr,int elemEnum, bool surfFlag) { - printf("Scaffolding of the mesh \n"); MyVolume volume; - ScalarType max = math::Max(baseMesh.bbox.DimX(),baseMesh.bbox.DimY(),baseMesh.bbox.DimZ()); - ScalarType voxel = max / volumeSide; - int sizeX = (baseMesh.bbox.DimX() / voxel)+1; - int sizeY = (baseMesh.bbox.DimY() / voxel)+1; - int sizeZ = (baseMesh.bbox.DimZ() / voxel)+1; - - // Kdtree -// seedTree->setMaxNofNeighbors(4); + int sizeX = (baseMesh.bbox.DimX() / voxelSide)+1; + int sizeY = (baseMesh.bbox.DimY() / voxelSide)+1; + int sizeZ = (baseMesh.bbox.DimZ() / voxelSide)+1; + + printf("Scaffolding of the mesh %i %i %i\n",sizeX,sizeY,sizeZ); BoxType bb = BoxType::Construct(baseMesh.bbox); - bb.Offset(baseMesh.bbox.Diag()*0.04f); + bb.Offset(voxelSide+isoThr*2.0f); volume.Init(Point3i(sizeX,sizeY,sizeZ),bb); -// qDebug("Init Volume of %i %i %i",sizeX,sizeY,sizeZ); int cnt=0; - ScalarType offset= volume.voxel.Norm()*isoThr; for(ScalarType i=0;i0) this->BuildMontecarloSampling(montecarloSampleNum); @@ -486,10 +507,10 @@ void QuadricRelaxVoronoiSamples(int relaxStep) tri::Append::MeshCopy(seedDomainMesh,montecarloVolumeMesh); vector pruningVec; - if(poissonRadius ==0 && seedNum!=0) - tri::PoissonPruningExact(seedDomainMesh,pruningVec,poissonRadius,seedNum); + if(poissonRadius ==0 && poissonSampleNum!=0) + tri::PoissonPruningExact(seedDomainMesh,pruningVec,poissonRadius,poissonSampleNum,0.04,10,randSeed); else - tri::PoissonPruning(seedDomainMesh,pruningVec,poissonRadius,seedNum); + tri::PoissonPruning(seedDomainMesh,pruningVec,poissonRadius,randSeed); std::vector seedPts(pruningVec.size()); for(size_t i=0;i