1 /*
2 Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
3 All rights reserved.
4 
5 Redistribution and use in source and binary forms, with or without modification,
6 are permitted provided that the following conditions are met:
7 
8 Redistributions of source code must retain the above copyright notice, this list of
9 conditions and the following disclaimer. Redistributions in binary form must reproduce
10 the above copyright notice, this list of conditions and the following disclaimer
11 in the documentation and/or other materials provided with the distribution.
12 
13 Neither the name of the Johns Hopkins University nor the names of its contributors
14 may be used to endorse or promote products derived from this software without specific
15 prior written permission.
16 
17 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
18 EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
19 OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
20 SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
21 INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
22 TO, PROCUREMENT OF SUBSTITUTE  GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
23 BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
24 CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
25 ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
26 DAMAGE.
27 */
28 
29 #include <stdio.h>
30 #include <stdlib.h>
31 #include <math.h>
32 #include <float.h>
33 #include "MarchingCubes.h"
34 #include "Octree.h"
35 #include "SparseMatrix.h"
36 #include "FunctionData.h"
37 #include "PPolynomial.h"
38 #include "PoissonParam.h"
39 
40 #include <wrap/callback.h>
41 
42 #define SCALE 1.25
43 
44 #include <stdarg.h>
45 
46 #include "MultiGridOctreeData.h"
47 
ShowUsage(char * ex)48 void ShowUsage(char* ex)
49 {
50 	printf("Usage: %s\n",ex);
51 	printf("\t--in  <input points>\n");
52 
53 	printf("\t--out <output triangle mesh>\n");
54 
55 	printf("\t[--depth <maximum reconstruction depth>]\n");
56 	printf("\t\t Running at depth d corresponds to solving on a 2^d x 2^d x 2^d\n");
57 	printf("\t\t voxel grid.\n");
58 
59 	printf("\t[--scale <scale factor>]\n");
60 	printf("\t\t Specifies the factor of the bounding cube that the input\n");
61 	printf("\t\t samples should fit into.\n");
62 
63 	printf("\t[--binary]\n");
64 	printf("\t\t If this flag is enabled, the point set is read in in\n");
65 	printf("\t\t binary format.\n");
66 
67 	printf("\t[--solverDivide <subdivision depth>]\n");
68 	printf("\t\t The depth at which a block Gauss-Seidel solver is used\n");
69 	printf("\t\t to solve the Laplacian.\n");
70 
71 	printf("\t[--samplesPerNode <minimum number of samples per node>\n");
72 	printf("\t\t This parameter specifies the minimum number of points that\n");
73 	printf("\t\t should fall within an octree node.\n");
74 
75 	printf("\t[--confidence]\n");
76 	printf("\t\t If this flag is enabled, the size of a sample's normals is\n");
77 	printf("\t\t used as a confidence value, affecting the sample's\n");
78 	printf("\t\t constribution to the reconstruction process.\n");
79 
80 	printf("\t[--verbose]\n");
81 }
82 template<int Degree>
Execute(PoissonParam & Par,std::vector<Point3D<Real>> Pts,std::vector<Point3D<Real>> Nor,CoredVectorMeshData & mesh,Point3D<Real> & newCenter,Real & newScale,vcg::CallBackPos * cb)83 int Execute(PoissonParam &Par, std::vector<Point3D<Real> > Pts, std::vector<Point3D<Real> > Nor, 	CoredVectorMeshData &mesh, Point3D<Real> &newCenter, Real &newScale, vcg::CallBackPos *cb)
84 {
85 		int i;
86 //	int paramNum=sizeof(paramNames)/sizeof(char*);
87 	//int commentNum=0;
88 	//char **comments;
89 
90 	//comments=new char*[paramNum+7];
91 	//for(i=0;i<=paramNum+7;i++){comments[i]=new char[1024];}
92 
93 //	const char* Rev = "Rev: V2 ";
94 //	const char* Date = "Date: 2006-11-09 (Thur, 09 Nov 2006) ";
95 
96 	// cmdLineParse(argc-1,&argv[1],paramNames,paramNum,params,0);
97 
98 	double t;
99 	Point3D<float> center;
100 	Real scale=1.0;
101 	Real isoValue=0;
102 	Octree<Degree> tree;
103 	PPolynomial<Degree> ReconstructionFunction=PPolynomial<Degree>::GaussianApproximation();
104 
105 	center.coords[0]=center.coords[1]=center.coords[2]=0;
106 
107 	TreeOctNode::SetAllocator(MEMORY_ALLOCATOR_BLOCK_SIZE);
108 
109 	int kernelDepth=Par.Depth-2;
110 	if(Par.KernelDepth>=0){kernelDepth=Par.KernelDepth;}
111 
112 	tree.setFunctionData(ReconstructionFunction,Par.Depth,0,Real(1.0)/(1<<Par.Depth));
113 //	DumpOutput("Memory Usage: %.3f MB\n",float(MemoryInfo::Usage())/(1<<20));
114 	if(kernelDepth>Par.Depth){
115 		fprintf(stderr,"KernelDepth can't be greater than Depth: %d <= %d\n",kernelDepth,Par.Depth);
116 		return EXIT_FAILURE;
117 	}
118 
119 
120 #if 1
121 	tree.setTree(Pts,Nor,Par.Depth,kernelDepth,Real(Par.SamplesPerNode),Par.Scale,center,scale,!Par.NoResetSamples,Par.Confidence);
122 #else
123 if(Confidence.set){
124 	tree.setTree(Pts,Nor,Depth.value,kernelDepth,Real(SamplesPerNode.value),Scale.value,center,scale,!NoResetSamples.set,0,1);
125 }
126 else{
127 	tree.setTree(Pts,Nor,Depth.value,kernelDepth,Real(SamplesPerNode.value),Scale.value,center,scale,!NoResetSamples.set,0,0);
128 }
129 #endif
130 	printf("Leaves/Nodes: %d/%d\n",tree.tree.leaves(),tree.tree.nodes());
131 	printf("   Tree Size: %.3f MB\n",float(sizeof(TreeOctNode)*tree.tree.nodes())/(1<<20));
132 
133 	if(!Par.NoClipTree){
134 		tree.ClipTree();
135 		printf("Leaves/Nodes: %d/%d\n",tree.tree.leaves(),tree.tree.nodes());
136 	}
137 
138 	tree.finalize1(Par.Refine);
139 	printf("Leaves/Nodes: %d/%d\n",tree.tree.leaves(),tree.tree.nodes());
140 
141 	tree.maxMemoryUsage=0;
142 	tree.SetLaplacianWeights();
143 
144 	tree.finalize2(Par.Refine);
145 
146 	tree.maxMemoryUsage=0;
147 	tree.LaplacianMatrixIteration(Par.SolverDivide);
148 
149 	tree.maxMemoryUsage=0;
150 	isoValue=tree.GetIsoValue();
151 	printf("IsoValue is %f \n",isoValue);
152   isoValue = isoValue * Par.Offset;
153 	printf("IsoValue is %f \n",isoValue);
154 
155 	if(Par.IsoDivide){tree.GetMCIsoTriangles(isoValue,Par.IsoDivide,&mesh);}
156 	else{tree.GetMCIsoTriangles(isoValue,&mesh);}
157 //	PlyWriteTriangles(Out.value,&mesh,PLY_BINARY_NATIVE,center,scale,comments,commentNum);
158 newCenter=center;
159 newScale=scale;
160 	return 1;
161 }
162 
Execute2(PoissonParam & Par,std::vector<Point3D<Real>> Pts,std::vector<Point3D<Real>> Nor,CoredVectorMeshData & mesh,Point3D<Real> & newCenter,Real & newScale,vcg::CallBackPos * cb)163 int Execute2(PoissonParam &Par, std::vector<Point3D<Real> > Pts, std::vector<Point3D<Real> > Nor, 	CoredVectorMeshData &mesh, Point3D<Real> &newCenter, Real &newScale, vcg::CallBackPos *cb )
164 {
165 	return Execute<2>(Par,Pts,Nor,mesh,newCenter,newScale,cb);
166 }
167