First version of the plymc vcg surface reconstructor

This commit is contained in:
Paolo Cignoni 2016-06-14 22:46:23 +02:00
parent cf308f9c2e
commit 78254e94af
8 changed files with 2318 additions and 0 deletions

10
apps/plymc/plymc.pro Executable file
View File

@ -0,0 +1,10 @@
TARGET = plymc
DEPENDPATH += .
INCLUDEPATH += ../..
CONFIG += console c++11
CONFIG -= app_bundle
TEMPLATE = app
SOURCES += ../../wrap/ply/plylib.cpp \
plymc_main.cpp

223
apps/plymc/plymc_main.cpp Normal file
View File

@ -0,0 +1,223 @@
/*
* plymc_main.cpp
* filter_plymc
*
* Created by Paolo Cignoni on 10/23/09.
* Copyright 2009 ISTI - CNR. All rights reserved.
*
*/
#include <vcg/complex/algorithms/create/plymc/plymc.h>
#include "simplemeshprovider.h"
#define _PLYMC_VER "4.0"
using namespace std;
using namespace vcg;
string MYbasename = "plymcout";
string VolumeBaseName;
string subtag=""; // la stringa da appendere al nome di file per generare i nomi di file relativi a vari sottopezzi
string alnfile;
/************ Command Line Parameters *******/
int saveMask=vcg::tri::io::Mask::IOM_ALL;
void usage()
{
printf(
"\nUsage:\n"
" plymc [options] filein.ply [filein.ply...]\n"
" plymc [options] filein.aln \n"
"Options: (no leading space before numeric values!) \n"
" -oname Set the base output name (default plymcout, without 'ply' extension)\n"
" -vname Enable the saving of the final volume with the specified filename\n"
" -C## Set numbers mesh that can be cached (default: 6)\n"
" -c## Set numbers of k cells (default: 10000)\n"
" -V# Set the required voxel size (override -c)\n"
" -s... Compute only a subvolume (specify 6 integers) \n"
" -S... Compute all the subvolumes of a partition (specify 3 int) \n"
" -X... Compute a range of the the subvolumes of a partition (specify 9 int)\n"
" -M Apply a 'safe' simplification step that removes only the unecessary triangles\n"
" -w# Set distance field Expansion factor in voxel (default 3)\n"
" -W# Set distance field Exp. as an absolute dist (override -w)\n"
" -a# Set angle threshold for distance field expansion (default 30)\n"
" -f# Set the fill threshold (default 12: all voxel having less \n"
" than 12 adjacent are not automaticall filled)\n"
" -L# Set Number of smoothing passes to be done after merging of all meshes\n"
" -R# Set Number of Refill passes to be done after merging of all meshes\n"
" -l# Make a single smoothing step after each expansion\n"
" -G# Disable Geodesic Quality\n"
" -F# Use per vertex quality defined in plyfile (geodesic is disabled)\n"
" -O# Set an Offset (<0!) threshold and build a double surface\n"
" -q# Set Quality threshold for smoothing. Only whose distance (in voxel)\n"
" is lower than the specified one are smoothed (default 3 voxel)\n"
" -Q# Same of above but expressed in absolute units.\n"
" -p use vertex splatting instead face rasterizing\n"
" -d# set <n> as verbose level (default 0)\n"
" -D# save <n> debug slices during processing\n"
"\nNotes:\n\n"
"The Quality threshold can be expressed in voxel unit or in absolute units.\n"
"It represents the geodetic distance from the mesh border.\n"
"I.e. -q3 means that all voxels that are within 3voxel from the mesh border\n"
"are smoothed.\n\n"
"A partition in the working volume is defined using 3 integers, that \n"
"specify the subdivision along the three axis.\n"
"To automatically compute ALL subvolumes of a given subdivision use '-S' \n"
"-S 1 1 1 default no subdivision at all\n"
"-S 2 2 2 compute all the octant of a octree-like subdivision\n"
"-S 1 1 4 compute four Z-slices\n\n"
"To work only on a SINGLE subvolume of the partition you have to specify \n"
"six integers: the first three ints specify the subdivision along the\n"
"three axis and the last three ones which subvolume you desire.\n"
"the parameter to be used is '-s' \n"
"-s 1 1 1 0 0 0 default no subdivision at all\n"
"-s 1 1 3 0 0 1 make three Z-slices and take the middle one \n"
"-s 2 2 2 1 1 1 the last octant in a octree-like subdivision\n\n"
"To START FROM a specific subvolume of the partition you have to specify \n"
"six integers: the first three ints specify the subdivision along the\n"
"three axis and the last three ones which subvolume you want to start\n"
"the process will go on using lexicographic order. Subvolumes are ordered\n"
"by Z, then by Y, then by X\n"
"The parameter to be used is '-K' \n"
"-K 4 4 4 0 0 0 a partition of 64 blocks, starting from the first one\n"
"-K 4 4 4 1 0 3 a partition of 64 blocks, starting from block 19 (index 1 0 3)\n\n"
"To work only on a specific subvolume range of the partition you have \n"
"to specify nine integers: the first three ints specify the subdivision\n"
"along the three axis and, the next three which is the starting subvolume\n"
"and the last three which is the last subvolume to be computed.\n"
"the process will compute all blocks with x,y,z index included in the interval\n"
"specified: Xstart<=X<=Xend Ystart<=Y<=Yend Zstart<=Z<=Zend\n"
"-X 3 3 3 0 0 0 2 2 2 three subdivision on each axis, all subvolumes\n"
"-X 2 2 2 1 0 0 1 1 1 three subdivision on each axis, only the 'right' part\n\n"
);
exit(-1);
}
int main(int argc, char *argv[])
{
Histogram<float> h;
tri::PlyMC<SMesh,SimpleMeshProvider<SMesh> > pmc;
tri::PlyMC<SMesh,SimpleMeshProvider<SMesh> >::Parameter &p = pmc.p;
// This line is required to be sure that the decimal separatore is ALWAYS the . and not the ,
// see the comment at the beginning of the file
setlocale(LC_ALL, "En_US");
printf( "\n PlyMC " _PLYMC_VER " (" __DATE__ ")\n"
" Copyright 2002-2016 Visual Computing Group I.S.T.I. C.N.R.\n"
" Paolo Cignoni (p.cignoni@isti.cnr.it)\n\n");
//// Parameters
int i=1;
if(argc<2) usage();
while(argv[i][0]=='-')
{
switch(argv[i][1])
{
case 'o' : p.basename=argv[i]+2;printf("Setting Basename to %s\n",MYbasename.c_str());break;
case 'C' : pmc.MP.setCacheSize(atoi(argv[i]+2));printf("Setting MaxSize of MeshCache to %i\n",atoi(argv[i]+2)); break;
case 'c' : p.NCell =atoi(argv[i]+2);printf("Setting NCell to %i\n",p.NCell); break;
case 'v' : p.SaveVolumeFlag=true; VolumeBaseName=argv[i]+2; printf("Saving Volume enabled: volume Basename to %s\n",VolumeBaseName.c_str());break;
case 'V' : p.VoxSize =atof(argv[i]+2);printf("Setting VoxSize to %f; overridden NCell\n",p.VoxSize);p.NCell=0;break;
case 'w' : p.WideNum =atoi(argv[i]+2);printf("Setting WideNum to %i\n",p.WideNum);break;
case 'W' : p.WideSize=atof(argv[i]+2);printf("Setting WideSize to %f;overridden WideNum\n",p.WideSize);break;
case 'L' : p.SmoothNum =atoi(argv[i]+2);printf("Setting Laplacian SmoothNum to %i\n",p.SmoothNum);break;
case 'R' : p.RefillNum =atoi(argv[i]+2);printf("Setting Refilling Num to %i\n",p.RefillNum);break;
case 'q' : p.QualitySmoothVox=atof(argv[i]+2);printf("Setting QualitySmoothThr to %f; \n",p.QualitySmoothVox);break;
case 'Q' : p.QualitySmoothAbs=atof(argv[i]+2);printf("Setting QualitySmoothAbsolute to %f; it will override the default %f voxel value\n",p.QualitySmoothAbs,p.QualitySmoothVox);break;
case 'l' : p.IntraSmoothFlag=true; printf("Setting Laplacian Smooth after expansion \n");break;
case 'G' : p.GeodesicQualityFlag=false; printf("Disabling Geodesic Quality\n");break;
case 'F' : p.PLYFileQualityFlag=true; p.GeodesicQualityFlag=false; printf("Enabling PlyFile (and disabling Geodesic) Quality\n");break;
case 'f' : p.FillThr=atoi(argv[i]+2);printf("Setting Fill Threshold to %i\n",p.FillThr);break;
case 'a' : p.ExpAngleDeg=atoi(argv[i]+2);printf("Setting Expanding Angle Threshold to %f Degree\n",p.ExpAngleDeg);break;
case 'O' : p.OffsetThr=atof(argv[i]+2);printf("Setting Offset Threshold to %f \n",p.OffsetThr);p.OffsetFlag=true;break;
case 's' :
p.IDiv[0]=atoi(argv[++i]); p.IDiv[1]=atoi(argv[++i]); p.IDiv[2]=atoi(argv[++i]);
p.IPosS[0]=atoi(argv[++i]);p.IPosS[1]=atoi(argv[++i]);p.IPosS[2]=atoi(argv[++i]);
p.IPosE[0]=p.IPosS[0]; p.IPosE[1]=p.IPosS[1]; p.IPosE[2]=p.IPosS[2];
if((p.IPosS[0]>=p.IDiv[0]) || (p.IPosS[1]>=p.IDiv[1]) || (p.IPosS[2]>=p.IDiv[2]))
{
printf("the subvolume you have requested is invalid (out of bounds)");
exit(-1);
}
printf("Computing ONLY subvolume [%i,%i,%i] on a %ix%ix%i\n",p.IPosS[0],p.IPosS[1],p.IPosS[2],p.IDiv[0],p.IDiv[1],p.IDiv[2]);
break;
case 'S' :
p.IDiv[0]=atoi(argv[++i]);p.IDiv[1]=atoi(argv[++i]);p.IDiv[2]=atoi(argv[++i]);
p.IPosS=Point3i(0,0,0);
p.IPosE[0]=p.IDiv[0]-1; p.IPosE[1]=p.IDiv[1]-1; p.IPosE[2]=p.IDiv[2]-1;
printf("Autocomputing ALL subvolumes on a %ix%ix%i\n",p.IDiv[0],p.IDiv[1],p.IDiv[2]);
break;
case 'K' :
p.IDiv[0]=atoi(argv[++i]); p.IDiv[1]=atoi(argv[++i]);p.IDiv[2]=atoi(argv[++i]);
p.IPosB[0]=atoi(argv[++i]);p.IPosB[1]=atoi(argv[++i]);p.IPosB[2]=atoi(argv[++i]);
p.IPosS=Point3i(0,0,0);
p.IPosE[0]=p.IDiv[0]-1; p.IPosE[1]=p.IDiv[1]-1; p.IPosE[2]=p.IDiv[2]-1;
if((p.IPosB[0]>=p.IDiv[0]) || (p.IPosB[1]>=p.IDiv[1]) || (p.IPosB[2]>=p.IDiv[2]))
{
printf("the start subvolume you have requested is invalid (out of bounds)");
exit(-1);
}
printf("Autocomputing ALL subvolumes FROM [%i,%i,%i] on a %ix%ix%i\n",p.IPosB[0],p.IPosB[1],p.IPosB[2],p.IDiv[0],p.IDiv[1],p.IDiv[2]);
break;
case 'X' :
p.IDiv[0]=atoi(argv[++i]); p.IDiv[1]=atoi(argv[++i]);p.IDiv[2]=atoi(argv[++i]);
p.IPosS[0]=atoi(argv[++i]);p.IPosS[1]=atoi(argv[++i]);p.IPosS[2]=atoi(argv[++i]);
p.IPosE[0]=atoi(argv[++i]);p.IPosE[1]=atoi(argv[++i]);p.IPosE[2]=atoi(argv[++i]);
// test if the interval is ok
int Istart,Iend;
Istart = p.IPosS[2] + (p.IPosS[1]*p.IDiv[2]) + (p.IPosS[0]*p.IDiv[2]*p.IDiv[1]);
Iend = p.IPosE[2] + (p.IPosE[1]*p.IDiv[2]) + (p.IPosE[0]*p.IDiv[2]*p.IDiv[1]);
if((Iend-Istart)<=0)
{
printf("the range you have requested is invalid (reversed or empty)");
exit(-1);
}
if((p.IPosS[0]>=p.IDiv[0]) || (p.IPosS[1]>=p.IDiv[1]) || (p.IPosS[2]>=p.IDiv[2]) ||
(p.IPosE[0]>=p.IDiv[0]) || (p.IPosE[1]>=p.IDiv[1]) || (p.IPosE[2]>=p.IDiv[2]))
{
printf("the subvolume you have requested is invalid (out of bounds)");
exit(-1);
}
printf("Autocomputing subvolumes FROM [%i,%i,%i] TO [%i,%i,%i] on a %ix%ix%i\n",p.IPosS[0],p.IPosS[1],p.IPosS[2],p.IPosE[0],p.IPosE[1],p.IPosE[2],p.IDiv[0],p.IDiv[1],p.IDiv[2]);
break;
// case 'B' : p.SafeBorder =atoi(argv[i]+2);printf("Setting SafeBorder among blocks to %i*%i (default 1)\n",p.SafeBorder,Volume<Voxelf>::BLOCKSIDE());break;
case 'p' : p.VertSplatFlag =true; printf("Enabling VertexSplatting instead of face rasterization\n");break;
case 'd' : p.VerboseLevel=atoi(argv[i]+2);printf("Enabling VerboseLevel= %i )\n",p.VerboseLevel);break;
case 'D' : p.VerboseLevel=1; p.SliceNum=atoi(argv[i]+2);printf("Enabling Debug Volume saving of %i slices (VerboseLevel=1)\n",p.SliceNum);break;
case 'M' : p.SimplificationFlag =true; printf("Enabling PostReconstruction simplification\n"); break;
default : {printf("Error unable to parse option '%s'\n",argv[i]); exit(0);}
}
++i;
}
Matrix44f Identity;
Identity.SetIdentity();
string alnfile;
while(i<argc)
{
if(strcmp(strrchr(argv[i],'.'),".aln")==0)
pmc.MP.openALN(argv[i]);
else
pmc.MP.AddSingleMesh(argv[i]);
++i;
}
if(pmc.MP.size()==0) usage();
printf("End Parsing\n\n");
pmc.Process();
return 0;
}

View File

@ -0,0 +1,206 @@
#ifndef SIMPLEMESHPROVIDER_H
#define SIMPLEMESHPROVIDER_H
#include <list>
#include <vector>
#include <vcg/space/box3.h>
#include <wrap/ply/plystuff.h>
#include <wrap/io_trimesh/import.h>
using namespace std;
using namespace vcg;
template<class TriMeshType>
class MeshCache
{
class Pair
{
public:
Pair(){used=0;}
TriMeshType *M;
std::string Name;
int used; // 'data' dell'ultimo accesso. si butta fuori quello lru
};
std::list<Pair> MV;
public:
void clear();
MeshCache() {MaxSize=6;}
~MeshCache() {
typename std::list<Pair>::iterator mi;
for(mi=MV.begin();mi!=MV.end();++mi)
delete (*mi).M;
}
// Restituisce true se la mesh e' in cache;
// restituisce in ogni caso il puntatore dove sta (o dovrebbe stare) la mesh
// Gestione LRU
bool Find(std::string &name, TriMeshType * &sm)
{
typename std::list<Pair>::iterator mi;
typename std::list<Pair>::iterator oldest; // quello che e' piu' tempo che non viene acceduto.
int last;
last = std::numeric_limits<int>::max();
oldest = MV.begin();
for(mi=MV.begin();mi!=MV.end();++mi)
{
if((*mi).used<last)
{
last=(*mi).used;
oldest=mi;
}
if((*mi).Name==name) {
sm=(*mi).M;
(*mi).used++;
return true;
}
}
// we have not found the requested mesh
// either allocate a new mesh or give back a previous mesh.
if(MV.size()>MaxSize) {
sm=(*oldest).M;
(*oldest).used=0;
(*oldest).Name=name;
} else {
MV.push_back(Pair());
MV.back().Name=name;
MV.back().M=new TriMeshType();
sm=MV.back().M;
}
return false;
}
size_t MaxSize;
size_t size() const {return MV.size();}
};
template<class TriMeshType>
class SimpleMeshProvider
{
std::vector< std::string > meshnames;
std::vector<vcg::Matrix44f> TrV;
std::vector<float> WV; // vettore con i pesi da applicare alla mesh.
std::vector<vcg::Box3f> BBV; // vettore con i bbox trasformati delle mesh da scannare.
vcg::Box3f fullBBox;
MeshCache<TriMeshType> MC;
public:
int size() {return meshnames.size();}
int getCacheSize() {return MC.MaxSize;}
int setCacheSize(size_t newsize)
{
if(newsize == MC.MaxSize)
return MC.MaxSize;
if(newsize <= 0)
return MC.MaxSize;
MC.MaxSize = newsize;
return newsize;
}
bool openALN (const char* alnName)
{
// vector<RangeMap> rmaps;
// ALNParser::ParseALN(rmaps, alnName);
// for(size_t i=0; i<rmaps.size(); i++)
// AddSingleMesh(rmaps[i].filename.c_str(), rmaps[i].trasformation, rmaps[i].quality);
return true;
}
bool AddSingleMesh(const char* meshName, Matrix44f &tr, float meshWeight=1)
{
assert(WV.size()==meshnames.size() && TrV.size() == WV.size());
TrV.push_back(tr);
meshnames.push_back(meshName);
WV.push_back(meshWeight);
BBV.push_back(Box3f());
return true;
}
bool AddSingleMesh(const char* meshName)
{
Matrix44f identity; identity.SetIdentity();
return AddSingleMesh(meshName, identity);
}
vcg::Box3f bb(int i) {return BBV[i];}
vcg::Box3f fullBB(){ return fullBBox;}
vcg::Matrix44f Tr(int i) const {return TrV[i];}
std::string MeshName(int i) const {return meshnames[i];}
float W(int i) const {return WV[i];}
void Clear()
{
meshnames.clear();
TrV.clear();
WV.clear();
BBV.clear();
fullBBox.SetNull();
MC.clear();
}
bool Find(int i, TriMeshType * &sm)
{
return MC.Find(meshnames[i],sm);
}
bool InitBBox()
{
fullBBox.SetNull();
for(int i=0;i<int(meshnames.size());++i)
{
Box3d b;
bool ret;
Matrix44f mt;
Matrix44f Id; Id.SetIdentity();
mt.Import(TrV[i]);
printf("bbox scanning %4i/%i [%16s] \r",i+1,(int)meshnames.size(), meshnames[i].c_str());
if(tri::io::Importer<TriMeshType>::FileExtension(meshnames[i],"PLY") || tri::io::Importer<TriMeshType>::FileExtension(meshnames[i],"ply"))
{
if(!(TrV[i]==Id))
ret=ply::ScanBBox(meshnames[i].c_str(),BBV[i],mt,true,0);
else
ret=vcg::ply::ScanBBox(meshnames[i].c_str(),BBV[i]);
}
else
{ printf("Trying to import a non-ply file %s\n",meshnames[i].c_str());fflush(stdout);
TriMeshType m;
int retVal=tri::io::Importer<TriMeshType>::Open(m,meshnames[i].c_str());
ret = (retVal==0);
tri::UpdateBounding<TriMeshType>::Box(m);
BBV[i].Import(m.bbox);
}
if( ! ret)
{
printf("\n\nwarning:\n file '%s' not found\n",meshnames[i].c_str());fflush(stdout);
continue;
}
fullBBox.Add(BBV[i]);
}
return true;
}
};
class SVertex;
class SFace;
class SUsedTypes: public vcg::UsedTypes < vcg::Use<SVertex>::AsVertexType,
vcg::Use<SFace >::AsFaceType >{};
class SVertex : public Vertex< SUsedTypes, vertex::Coord3f, vertex::Normal3f,vertex::VFAdj, vertex::BitFlags, vertex::Color4b, vertex::Qualityf>{};
class SFace : public Face< SUsedTypes, face::VertexRef, face::Normal3f,face::Qualityf, face::VFAdj, face::BitFlags> {};
class SMesh : public vcg::tri::TriMesh< std::vector< SVertex>, std::vector< SFace > > {};
#endif // SIMPLEMESHPROVIDER_H

View File

@ -54,6 +54,8 @@
#include <wrap/ply/plystuff.h>
#include <vcg/complex/algorithms/create/marching_cubes.h>
#include <vcg/complex/algorithms/create/mc_trivial_walker.h>
// local optimization
#include <vcg/complex/algorithms/local_optimization.h>
#include <vcg/complex/algorithms/local_optimization/tri_edge_collapse.h>

View File

@ -0,0 +1,174 @@
/****************************************************************************
* VCGLib o o *
* Visual and Computer Graphics Library o o *
* _ O _ *
* Copyright(C) 2004-2016 \/)\/ *
* Visual Computing Lab /\/| *
* ISTI - Italian National Research Council | *
* \ *
* All rights reserved. *
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) *
* for more details. *
* *
****************************************************************************/
#ifndef __SVOXEL_H__
#define __SVOXEL_H__
template<class SCALAR_TYPE=float>
class SVoxel
{
private:
short v; // distanza dalla superficie espressa in centesimi di voxel;
short q; // come distanza dal bordo espressa in centesimi di voxel; // nota che questo implica che non si potrebbe rasterizzare una singola mesh in un volume di lato > 3000^3
char n[3];
unsigned char cnt;
// se != 0 b deve essere true;
// b e'tenuto implicitamente usando cnt;
// se cnt == 0 || cnt>0 -> b=false
// se cnt == 255 -> b=true
// b==false cnt==0 totalmente non inzializzato (Zero)
// b==false cnt >0 da normalizzare
// b==true cnt==0 gia' normalizzato
// b==true cnt >0 Errore!!!
public:
typedef SCALAR_TYPE scalar;
SVoxel(SCALAR_TYPE vv, bool /*bb*/, Point3<scalar> &nn, scalar qq) {SetV(vv);SetN(nn);SetQ(qq);}
SVoxel(SCALAR_TYPE vv, const Point3<scalar> &nn, scalar qq) {SetV(vv);SetN(nn);SetQ(qq);cnt=255;}
bool B() const {return cnt==255;} // puo' essere a true solo se cnt == 0; (il che significa che e' stato gia' normalizzato
void SetB(bool val) {
assert( val == (cnt==255 || cnt==0) );
if(val) cnt=255;
else if(cnt==255) cnt=0;
}
int Cnt() {
if(cnt==255) return 0;
else return int(cnt);
}
void SetCnt(int val) { cnt=(unsigned char)val;}
Point3<scalar> N() const {
return Point3<scalar>(scalar(n[0])/127.0f,scalar(n[1])/127.0f,scalar(n[2])/127.0f);
}
scalar N(const int i) const
{
return scalar(n[i])/127.0f;
}
void SetN(const Point3<scalar> &nn)
{
n[0]=char( nn[0]*127.0f );
n[1]=char( nn[1]*127.0f );
n[2]=char( nn[2]*127.0f );
}
scalar V() const
{
return scalar(v)/100.0f;
}
inline void Blend( SVoxel const & vx, scalar w)
{
float w1=1.0-w;
SetV(V()*w1+vx.V()*w);
SetQ(Q()*w1+vx.Q()*w);
SetN(N()*w1+vx.N()*w);
//return *this;
}
void SetV(const float &vv)
{
v= short(vv*100.0f);
if(v==0) v=1;
}
scalar Q() const
{
return scalar(q)/20.0f;
}
void SetQ(const float &qq)
{
int qi = qq * 20.0f;
if (qi>32767) qi =32767;
q = qi;
}
void Merge(const SVoxel &VOX)
{
v=(v*Q()+VOX.Q()*VOX.v)/(Q()+VOX.Q());
SetQ(Q()+VOX.Q());
}
void Set(const SVoxel &VOX)
{
v=VOX.v;
n[0]=VOX.n[0];
n[1]=VOX.n[1];
n[2]=VOX.n[2];
q=VOX.q;
}
inline SVoxel & operator += ( SVoxel const & vx)
{
if(cnt==0)
{
v=vx.v;
q=vx.q;
n[0]=vx.n[0];
n[1]=vx.n[1];
n[2]=vx.n[2];
cnt=1;
}
else
{
const int cnt1=cnt+1;
v+=vx.v;
q = (q*cnt+vx.q)/(cnt1);
n[0]=(vx.n[0]*int(cnt) +vx.n[0])/(cnt1) ;
n[1]=(vx.n[1]*int(cnt) +vx.n[1])/(cnt1) ;
n[2]=(vx.n[2]*int(cnt) +vx.n[2])/(cnt1) ;
if(cnt==255) cnt=1;
else ++cnt;
}
return *this;
}
inline bool Normalize(int thr)
{
assert(cnt>0);
if(cnt<thr)
{
(*this) = Zero();
return false;
}
v= v/cnt; // gli altri membri sono gia' normalizzati
cnt=255;
// b=true; inutile!
return true;
}
bool IsZero() const { return v==0; }
static const SVoxel &Zero() {
static SVoxel tt(0,false,Point3f(0,0,0),0);
return tt;
}
};
#endif

View File

@ -0,0 +1,129 @@
/****************************************************************************
* VCGLib o o *
* Visual and Computer Graphics Library o o *
* _ O _ *
* Copyright(C) 2004-2016 \/)\/ *
* Visual Computing Lab /\/| *
* ISTI - Italian National Research Council | *
* \ *
* All rights reserved. *
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) *
* for more details. *
* *
****************************************************************************/
#ifndef __TRI_EDGE_COLLAPSE_MC__
#define __TRI_EDGE_COLLAPSE_MC__
#include<vcg/simplex/face/topology.h>
class TriEdgeCollapseMCParameter : public BaseParameterClass
{
public:
Box3f bb;
bool preserveBBox;
float areaThr;
void SetDefaultParams()
{
bb.SetNull();
preserveBBox=true;
areaThr=0;
}
TriEdgeCollapseMCParameter() {SetDefaultParams();}
};
template<class MCTriMesh, class VertexPair, class MYTYPE >
// Specialized Simplification classes for removing all small pieces on meshes.
class MCTriEdgeCollapse: public tri::TriEdgeCollapse< MCTriMesh, VertexPair, MYTYPE> {
public:
typedef typename MCTriMesh::VertexPointer VertexPointer;
typedef typename MCTriMesh::FaceType FaceType;
typedef typename MCTriMesh::VertexType::CoordType CoordType;
typedef typename MCTriMesh::VertexType::ScalarType ScalarType;
inline MCTriEdgeCollapse( const VertexPair &p, int mark, BaseParameterClass *pp)
{
this->localMark = mark;
this->pos=p;
this->_priority = ComputePriority(pp);
}
// The priority is simply the edge lenght,
// but we consider collapsing edges that lies on one plane of the MC grid.
virtual inline ScalarType ComputePriority(BaseParameterClass *_pp)
{
TriEdgeCollapseMCParameter *pp=(TriEdgeCollapseMCParameter *)_pp;
const CoordType & p0 = this->pos.V(0)->cP();
const CoordType & p1 = this->pos.V(1)->cP();
int diffCnt=0;
if( (p0[0] != p1[0]) ) diffCnt++;
if( (p0[1] != p1[1]) ) diffCnt++;
if( (p0[2] != p1[2]) ) diffCnt++;
// non MC plane collapse return highest cost
// if(diffCnt>2) return this->_priority=std::numeric_limits<float>::max() ;
if(pp->preserveBBox)
{
const Box3f &bb=pp->bb;
// collapse on the bbox border. Avoid it.
if(p0[0]==bb.min[0] || p0[0]==bb.max[0]) return this->_priority=std::numeric_limits<float>::max() ;
if(p0[1]==bb.min[1] || p0[1]==bb.max[1]) return this->_priority=std::numeric_limits<float>::max() ;
if(p0[2]==bb.min[2] || p0[2]==bb.max[2]) return this->_priority=std::numeric_limits<float>::max() ;
if(p1[0]==bb.min[0] || p1[0]==bb.max[0]) return this->_priority=std::numeric_limits<float>::max() ;
if(p1[1]==bb.min[1] || p1[1]==bb.max[1]) return this->_priority=std::numeric_limits<float>::max() ;
if(p1[2]==bb.min[2] || p1[2]==bb.max[2]) return this->_priority=std::numeric_limits<float>::max() ;
}
return this->_priority=Distance(p0,p1);
}
inline void Execute(MCTriMesh &m,BaseParameterClass *)
{
const CoordType & p0 = this->pos.V(0)->cP();
const CoordType & p1 = this->pos.V(1)->cP();
std::vector<VertexPointer> starVec0;
std::vector<VertexPointer> starVec1;
// int count0=0,count1=0;
vcg::face::VVStarVF<FaceType>(this->pos.V(0),starVec0);
// for(size_t i=0;i<starVec.size();++i)
// {
// if( (p0[0]==starVec[i]->cP()[0]) ) count0++;
// if( (p0[1]==starVec[i]->cP()[1]) ) count0++;
// if( (p0[2]==starVec[i]->cP()[2]) ) count0++;
// }
vcg::face::VVStarVF<FaceType>(this->pos.V(1),starVec1);
// for(size_t i=0;i<starVec.size();++i)
// {
// if( (p1[0]==starVec[i]->cP()[0]) ) count1++;
// if( (p1[1]==starVec[i]->cP()[1]) ) count1++;
// if( (p1[2]==starVec[i]->cP()[2]) ) count1++;
// }
CoordType MidPoint=(p0+p1)/2.0;
if(starVec0.size()>starVec1.size()) MidPoint=p0;
if(starVec0.size()<starVec1.size()) MidPoint=p1;
// assert( (p0[0]==p1[0]) ||
// (p0[1]==p1[1]) ||
// (p0[2]==p1[2]) );
// DoCollapse(m, this->pos, MidPoint);
vcg::tri::EdgeCollapser<MCTriMesh,VertexPair>::Do(m, this->pos, MidPoint);
}
};
#endif

File diff suppressed because it is too large Load Diff

View File

@ -0,0 +1,195 @@
/****************************************************************************
* VCGLib o o *
* Visual and Computer Graphics Library o o *
* _ O _ *
* Copyright(C) 2004-2016 \/)\/ *
* Visual Computing Lab /\/| *
* ISTI - Italian National Research Council | *
* \ *
* All rights reserved. *
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License (http://www.gnu.org/licenses/gpl.txt) *
* for more details. *
* *
****************************************************************************/
#ifndef __VOXEL_H__
#define __VOXEL_H__
template<class SCALAR_TYPE=float>
class Voxel
{
public:
typedef SCALAR_TYPE scalar;
Voxel(SCALAR_TYPE vv, bool bb, Point3<scalar> nn, short _cnt) {v=vv;b=bb;n=nn;cnt=_cnt;}
Voxel(SCALAR_TYPE vv, Point3<scalar> nn, scalar qq) {v=vv;b=true;n=nn;cnt=0;q=qq;}
const scalar &N(const int i) const { return n[i]; }
const Point3<scalar> &N() const { return n; }
void SetN(const Point3<scalar> &nn) { n=nn; }
const scalar &V() const { return v; }
void SetV(const scalar &vv) { v=vv; }
const scalar &Q() const { return q; }
void SetQ(const scalar &qq) { q=qq; }
bool B() const {return b;};
void SetB(bool val) {b=val;}
int Cnt() const {return cnt;}
void SetCnt(int val) {cnt=val;}
inline void Blend( Voxel const & vx, scalar w)
{
float w1=1.0-w;
v=v*w1+vx.v*w;
q=q*w1+vx.q*w;
n=n*w1+vx.n*w;
//return *this;
}
inline Voxel & operator += ( Voxel const & vx)
{
if(cnt==0)
{
assert(!b);
v=vx.v;
q=vx.q;
n=vx.n;
cnt=1;
b=false;
}
else
{
assert(!b);
v+=vx.v;
q+=vx.q;
n+=vx.n;
++cnt;
}
return *this;
}
inline bool Normalize(int thr)
{
assert(cnt>0);
assert(!B());
if(cnt<thr)
{
(*this) = Zero();
return false;
}
v/=cnt;
q/=cnt;
n/=cnt;
cnt=0;
b=true;
return true;
}
static const Voxel &Zero() {
static Voxel tt(0,false,Point3f(0,0,0),0);
return tt;
}
void Merge(const Voxel &VOX)
{
v=(v*q+VOX.Q()*VOX.v)/(q+VOX.Q());
n=(n*q+VOX.n*VOX.Q())/(q+VOX.Q());
q=q+VOX.Q();
}
void Set(const Voxel &VOX)
{
v=VOX.v;
n=VOX.n;
q=VOX.q;
}
protected:
bool b;
short cnt;
scalar v;
scalar q;
Point3<SCALAR_TYPE> n;
};
class Voxelfc :public Voxel<float>
{
public:
Voxelfc(float vv, bool bb, Point3f nn, short _cnt) :Voxel<float>(vv,bb,nn,_cnt){}
Voxelfc(float vv, Point3f nn, scalar qq) :Voxel<float>(vv,nn,qq) {}
Voxelfc(float vv, Point3f nn, scalar qq,Color4b cc) :Voxel<float>(vv,nn,qq)
{
c[0]=cc[0];
c[1]=cc[1];
c[2]=cc[2];
}
inline bool Normalize(int thr)
{
if(cnt>=thr) c/=cnt;
return Voxel<float>::Normalize(thr);
}
static const Voxelfc &Zero() {
static Voxelfc tt(0,false,Point3f(0,0,0),0);
return tt;
}
void Merge(const Voxelfc &VOX)
{
c=( c*q + VOX.C()*VOX.Q() )/(q+VOX.Q());
Voxel<float>::Merge(VOX);
}
void Set(const Voxelfc &VOX)
{
Voxel<float>::Set(VOX);
c=VOX.c;
}
const float &C(const int i) const { return c[i]; }
const Point3f &C() const { return c; }
void SetC(const Point3f &cc) { c=cc; }
Color4b C4b() const
{
static Color4b cc;
cc=Color4b(c[0],c[1],c[2],255);
return cc;
}
inline void Blend( Voxelfc const & vx, scalar w)
{
float w1=1.0-w;
v=v*w1+vx.v*w;
q=q*w1+vx.q*w;
n=n*w1+vx.n*w;
c=c*w1+vx.c*w;
//return *this;
}
inline Voxelfc & operator += ( Voxelfc const & vx)
{
Voxel<float>::operator +=(vx);
if(cnt==1) c =vx.c;
else c+=vx.c;
return *this;
}
private:
Point3f c;
};
#endif