1/* 2Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho 3All rights reserved. 4 5Redistribution and use in source and binary forms, with or without modification, 6are permitted provided that the following conditions are met: 7 8Redistributions of source code must retain the above copyright notice, this list of 9conditions and the following disclaimer. Redistributions in binary form must reproduce 10the above copyright notice, this list of conditions and the following disclaimer 11in the documentation and/or other materials provided with the distribution. 12 13Neither the name of the Johns Hopkins University nor the names of its contributors 14may be used to endorse or promote products derived from this software without specific 15prior written permission. 16 17THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY 18EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES 19OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT 20SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, 21INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED 22TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR 23BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN 24CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN 25ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH 26DAMAGE. 27*/ 28 29#include <stdio.h> 30#include <stdlib.h> 31#include <math.h> 32#include <float.h> 33#ifdef _WIN32 34#include <Windows.h> 35#include <Psapi.h> 36#endif // _WIN32 37#include "MyTime.h" 38#include "MarchingCubes.h" 39#include "Octree.h" 40#include "SparseMatrix.h" 41#include "CmdLineParser.h" 42#include "PPolynomial.h" 43#include "Ply.h" 44#include "MemoryUsage.h" 45#ifdef _OPENMP 46#include "omp.h" 47#endif // _OPENMP 48void DumpOutput( const char* format , ... ); 49void DumpOutput2( std::vector< char* >& comments , const char* format , ... ); 50#include "MultiGridOctreeData.h" 51 52#define DEFAULT_FULL_DEPTH 5 53 54#define XSTR(x) STR(x) 55#define STR(x) #x 56#if DEFAULT_FULL_DEPTH 57#pragma message ( "[WARNING] Setting default full depth to " XSTR(DEFAULT_FULL_DEPTH) ) 58#endif // DEFAULT_FULL_DEPTH 59 60#include <stdarg.h> 61char* outputFile=NULL; 62int echoStdout=0; 63void DumpOutput( const char* format , ... ) 64{ 65 if( outputFile ) 66 { 67 FILE* fp = fopen( outputFile , "a" ); 68 va_list args; 69 va_start( args , format ); 70 vfprintf( fp , format , args ); 71 fclose( fp ); 72 va_end( args ); 73 } 74 if( echoStdout ) 75 { 76 va_list args; 77 va_start( args , format ); 78 vprintf( format , args ); 79 va_end( args ); 80 } 81} 82void DumpOutput2( std::vector< char* >& comments , const char* format , ... ) 83{ 84 if( outputFile ) 85 { 86 FILE* fp = fopen( outputFile , "a" ); 87 va_list args; 88 va_start( args , format ); 89 vfprintf( fp , format , args ); 90 fclose( fp ); 91 va_end( args ); 92 } 93 if( echoStdout ) 94 { 95 va_list args; 96 va_start( args , format ); 97 vprintf( format , args ); 98 va_end( args ); 99 } 100 comments.push_back( new char[1024] ); 101 char* str = comments.back(); 102 va_list args; 103 va_start( args , format ); 104 vsprintf( str , format , args ); 105 va_end( args ); 106 if( str[strlen(str)-1]=='\n' ) str[strlen(str)-1] = 0; 107} 108 109 110cmdLineString 111 In( "in" ) , 112 Out( "out" ) , 113 VoxelGrid( "voxel" ) , 114 XForm( "xForm" ); 115 116cmdLineReadable 117#ifdef _WIN32 118 Performance( "performance" ) , 119#endif // _WIN32 120 Complete( "complete" ) , 121 ShowResidual( "showResidual" ) , 122 NoComments( "noComments" ) , 123 PolygonMesh( "polygonMesh" ) , 124 Confidence( "confidence" ) , 125 NormalWeights( "nWeights" ) , 126 NonManifold( "nonManifold" ) , 127 Dirichlet( "dirichlet" ) , 128 ASCII( "ascii" ) , 129 Density( "density" ) , 130 LinearFit( "linearFit" ) , 131 PrimalVoxel( "primalVoxel" ) , 132 Verbose( "verbose" ) , 133 Double( "double" ); 134 135cmdLineInt 136 Degree( "degree" , 2 ) , 137 Depth( "depth" , 8 ) , 138 CGDepth( "cgDepth" , 0 ) , 139 KernelDepth( "kernelDepth" ) , 140 AdaptiveExponent( "adaptiveExp" , 1 ) , 141 Iters( "iters" , 8 ) , 142 VoxelDepth( "voxelDepth" , -1 ) , 143 FullDepth( "fullDepth" , DEFAULT_FULL_DEPTH ) , 144 MinDepth( "minDepth" , 0 ) , 145 MaxSolveDepth( "maxSolveDepth" ) , 146 Threads( "threads" , omp_get_num_procs() ); 147 148cmdLineFloat 149 Color( "color" , 16.f ) , 150 SamplesPerNode( "samplesPerNode" , 1.5f ) , 151 Scale( "scale" , 1.1f ) , 152 CSSolverAccuracy( "cgAccuracy" , float(1e-3) ) , 153 PointWeight( "pointWeight" , 4.f ); 154 155 156cmdLineReadable* params[] = 157{ 158 &In , &Degree , &Depth , &Out , &XForm , 159 &Scale , &Verbose , &CSSolverAccuracy , &NoComments , &Double , 160 &KernelDepth , &SamplesPerNode , &Confidence , &NormalWeights , &NonManifold , &PolygonMesh , &ASCII , &ShowResidual , &VoxelDepth , 161 &PointWeight , &VoxelGrid , &Threads , &MaxSolveDepth , 162 &AdaptiveExponent , &Dirichlet , 163 &Density , 164 &FullDepth , 165 &MinDepth , 166 &CGDepth , &Iters , 167 &Complete , 168 &Color , 169 &LinearFit , 170 &PrimalVoxel , 171#ifdef _WIN32 172 &Performance , 173#endif // _WIN32 174}; 175 176 177void ShowUsage(char* ex) 178{ 179 printf( "Usage: %s\n" , ex ); 180 printf( "\t --%s <input points>\n" , In.name ); 181 182 printf( "\t[--%s <ouput triangle mesh>]\n" , Out.name ); 183 printf( "\t[--%s <ouput voxel grid>]\n" , VoxelGrid.name ); 184 185 printf( "\t[--%s <b-spline degree>=%d]\n" , Degree.name , Degree.value ); 186 187 printf( "\t[--%s <maximum reconstruction depth>=%d]\n" , Depth.name , Depth.value ); 188 printf( "\t\t Running at depth d corresponds to solving on a 2^d x 2^d x 2^d\n" ); 189 printf( "\t\t voxel grid.\n" ); 190 191 printf( "\t[--%s <full depth>=%d]\n" , FullDepth.name , FullDepth.value ); 192 printf( "\t\t This flag specifies the depth up to which the octree should be complete.\n" ); 193 194 printf( "\t[--%s <depth at which to extract the voxel grid>=<%s>]\n" , VoxelDepth.name , Depth.name ); 195 196 printf( "\t[--%s <conjugate-gradients depth>=%d]\n" , CGDepth.name , CGDepth.value ); 197 printf( "\t\t The depth up to which a conjugate-gradients solver should be used.\n"); 198 199 printf( "\t[--%s <scale factor>=%f]\n" , Scale.name , Scale.value ); 200 printf( "\t\t Specifies the factor of the bounding cube that the input\n" ); 201 printf( "\t\t samples should fit into.\n" ); 202 203 printf( "\t[--%s <minimum number of samples per node>=%f]\n" , SamplesPerNode.name, SamplesPerNode.value ); 204 printf( "\t\t This parameter specifies the minimum number of points that\n" ); 205 printf( "\t\t should fall within an octree node.\n" ); 206 207 printf( "\t[--%s <interpolation weight>=%f]\n" , PointWeight.name , PointWeight.value ); 208 printf( "\t\t This value specifies the weight that point interpolation constraints are\n" ); 209 printf( "\t\t given when defining the (screened) Poisson system.\n" ); 210 211 printf( "\t[--%s <iterations>=%d]\n" , Iters.name , Iters.value ); 212 printf( "\t\t This flag specifies the (maximum if CG) number of solver iterations.\n" ); 213 214 printf( "\t[--%s <pull factor>]\n" , Color.name ); 215 printf( "\t\t This flag specifies the pull factor for color interpolation\n" ); 216 217#ifdef _OPENMP 218 printf( "\t[--%s <num threads>=%d]\n" , Threads.name , Threads.value ); 219 printf( "\t\t This parameter specifies the number of threads across which\n" ); 220 printf( "\t\t the solver should be parallelized.\n" ); 221#endif // _OPENMP 222 223 printf( "\t[--%s]\n" , Confidence.name ); 224 printf( "\t\t If this flag is enabled, the size of a sample's normals is\n" ); 225 printf( "\t\t used as a confidence value, affecting the sample's\n" ); 226 printf( "\t\t constribution to the reconstruction process.\n" ); 227 228 printf( "\t[--%s]\n" , NormalWeights.name ); 229 printf( "\t\t If this flag is enabled, the size of a sample's normals is\n" ); 230 printf( "\t\t used as to modulate the interpolation weight.\n" ); 231 232#if 0 233 printf( "\t[--%s]\n" , NonManifold.name ); 234 printf( "\t\t If this flag is enabled, the isosurface extraction does not add\n" ); 235 printf( "\t\t a planar polygon's barycenter in order to ensure that the output\n" ); 236 printf( "\t\t mesh is manifold.\n" ); 237#endif 238 239 printf( "\t[--%s]\n" , PolygonMesh.name); 240 printf( "\t\t If this flag is enabled, the isosurface extraction returns polygons\n" ); 241 printf( "\t\t rather than triangles.\n" ); 242 243#if 0 244 printf( "\t[--%s <minimum depth>=%d]\n" , MinDepth.name , MinDepth.value ); 245 printf( "\t\t This flag specifies the coarsest depth at which the system is to be solved.\n" ); 246 247 printf( "\t[--%s <cg solver accuracy>=%g]\n" , CSSolverAccuracy.name , CSSolverAccuracy.value ); 248 printf( "\t\t This flag specifies the accuracy cut-off to be used for CG.\n" ); 249 250 printf( "\t[--%s <adaptive weighting exponent>=%d]\n", AdaptiveExponent.name , AdaptiveExponent.value ); 251 printf( "\t\t This flag specifies the exponent scale for the adaptive weighting.\n" ); 252 253#ifdef _WIN32 254 printf( "\t[--%s]\n" , Performance.name ); 255 printf( "\t\t If this flag is enabled, the running time and peak memory usage\n" ); 256 printf( "\t\t is output after the reconstruction.\n" ); 257#endif // _WIN32 258#endif 259 printf( "\t[--%s]\n" , Dirichlet.name); 260 printf( "\t\t If this flag is enabled, Dirichlet boundary constraints are used for reconstruction.\n" ); 261 262 printf( "\t[--%s]\n" , Density.name ); 263 printf( "\t\t If this flag is enabled, the sampling density is written out with the vertices.\n" ); 264 265 printf( "\t[--%s]\n" , LinearFit.name ); 266 printf( "\t\t If this flag is enabled, the iso-surfacing will be performed using linear fitting.\n" ); 267 268 printf( "\t[--%s]\n" , PrimalVoxel.name ); 269 printf( "\t\t If this flag is enabled, voxel sampling is performed at corners rather than centers.\n" ); 270 271#if 0 272 printf( "\t[--%s]\n" , ASCII.name ); 273 printf( "\t\t If this flag is enabled, the output file is written out in ASCII format.\n" ); 274 275 printf( "\t[--%s]\n" , NoComments.name ); 276 printf( "\t\t If this flag is enabled, the output file will not include comments.\n" ); 277#endif 278 279 printf( "\t[--%s]\n" , Double.name ); 280 printf( "\t\t If this flag is enabled, the reconstruction will be performed with double-precision floats.\n" ); 281 282 printf( "\t[--%s]\n" , Verbose.name ); 283 printf( "\t\t If this flag is enabled, the progress of the reconstructor will be output to STDOUT.\n" ); 284} 285 286Point3D< unsigned char > ReadASCIIColor( FILE* fp ) 287{ 288 Point3D< unsigned char > c; 289 if( fscanf( fp , " %c %c %c " , &c[0] , &c[1] , &c[2] )!=3 ) fprintf( stderr , "[ERROR] Failed to read color\n" ) , exit( 0 ); 290 return c; 291} 292 293PlyProperty PlyColorProperties[]= 294{ 295 { "r" , PLY_UCHAR , PLY_UCHAR , int( offsetof( Point3D< unsigned char > , coords[0] ) ) , 0 , 0 , 0 , 0 } , 296 { "g" , PLY_UCHAR , PLY_UCHAR , int( offsetof( Point3D< unsigned char > , coords[1] ) ) , 0 , 0 , 0 , 0 } , 297 { "b" , PLY_UCHAR , PLY_UCHAR , int( offsetof( Point3D< unsigned char > , coords[2] ) ) , 0 , 0 , 0 , 0 } , 298 { "red" , PLY_UCHAR , PLY_UCHAR , int( offsetof( Point3D< unsigned char > , coords[0] ) ) , 0 , 0 , 0 , 0 } , 299 { "green" , PLY_UCHAR , PLY_UCHAR , int( offsetof( Point3D< unsigned char > , coords[1] ) ) , 0 , 0 , 0 , 0 } , 300 { "blue" , PLY_UCHAR , PLY_UCHAR , int( offsetof( Point3D< unsigned char > , coords[2] ) ) , 0 , 0 , 0 , 0 } 301}; 302 303bool ValidPlyColorProperties( const bool* props ){ return ( props[0] || props[3] ) && ( props[1] || props[4] ) && ( props[2] || props[5] ); } 304 305template< class Real , int Degree , class Vertex > 306int _Execute( int argc , char* argv[] ) 307{ 308 Reset< Real >(); 309 int paramNum = sizeof(params)/sizeof(cmdLineReadable*); 310 std::vector< char* > comments; 311 312 if( Verbose.set ) echoStdout=1; 313 314 XForm4x4< Real > xForm , iXForm; 315 if( XForm.set ) 316 { 317 FILE* fp = fopen( XForm.value , "r" ); 318 if( !fp ) 319 { 320 fprintf( stderr , "[WARNING] Could not read x-form from: %s\n" , XForm.value ); 321 xForm = XForm4x4< Real >::Identity(); 322 } 323 else 324 { 325 for( int i=0 ; i<4 ; i++ ) for( int j=0 ; j<4 ; j++ ) 326 { 327 float f; 328 if( fscanf( fp , " %f " , &f )!=1 ) fprintf( stderr , "[ERROR] Execute: Failed to read xform\n" ) , exit( 0 ); 329 xForm(i,j) = (Real)f; 330 } 331 fclose( fp ); 332 } 333 } 334 else xForm = XForm4x4< Real >::Identity(); 335 iXForm = xForm.inverse(); 336 337 DumpOutput2( comments , "Running Screened Poisson Reconstruction (Version 8.0)\n" ); 338 char str[1024]; 339 for( int i=0 ; i<paramNum ; i++ ) 340 if( params[i]->set ) 341 { 342 params[i]->writeValue( str ); 343 if( strlen( str ) ) DumpOutput2( comments , "\t--%s %s\n" , params[i]->name , str ); 344 else DumpOutput2( comments , "\t--%s\n" , params[i]->name ); 345 } 346 347 double t; 348 double tt=Time(); 349 Real isoValue = 0; 350 351 Octree< Real > tree; 352 tree.threads = Threads.value; 353 if( !In.set ) 354 { 355 ShowUsage( argv[0] ); 356 return 0; 357 } 358 if( !MaxSolveDepth.set ) MaxSolveDepth.value = Depth.value; 359 360 OctNode< TreeNodeData >::SetAllocator( MEMORY_ALLOCATOR_BLOCK_SIZE ); 361 362 t=Time(); 363 int kernelDepth = KernelDepth.set ? KernelDepth.value : Depth.value-2; 364 if( kernelDepth>Depth.value ) 365 { 366 fprintf( stderr,"[WARNING] %s can't be greater than %s: %d <= %d\n" , KernelDepth.name , Depth.name , KernelDepth.value , Depth.value ); 367 kernelDepth = Depth.value; 368 } 369 370 double maxMemoryUsage; 371 t=Time() , tree.maxMemoryUsage=0; 372 SparseNodeData< PointData< Real > , 0 >* pointInfo = new SparseNodeData< PointData< Real > , 0 >(); 373 SparseNodeData< Point3D< Real > , NORMAL_DEGREE >* normalInfo = new SparseNodeData< Point3D< Real > , NORMAL_DEGREE >(); 374 SparseNodeData< Real , WEIGHT_DEGREE >* densityWeights = new SparseNodeData< Real , WEIGHT_DEGREE >(); 375 SparseNodeData< Real , NORMAL_DEGREE >* nodeWeights = new SparseNodeData< Real , NORMAL_DEGREE >(); 376 int pointCount; 377 typedef typename Octree< Real >::template ProjectiveData< Point3D< Real > > ProjectiveColor; 378 SparseNodeData< ProjectiveColor , DATA_DEGREE >* colorData = NULL; 379 380 char* ext = GetFileExtension( In.value ); 381 if( Color.set && Color.value>0 ) 382 { 383 colorData = new SparseNodeData< ProjectiveColor , DATA_DEGREE >(); 384 OrientedPointStreamWithData< float , Point3D< unsigned char > >* pointStream; 385 if ( !strcasecmp( ext , "bnpts" ) ) pointStream = new BinaryOrientedPointStreamWithData< float , Point3D< unsigned char > >( In.value ); 386 else if( !strcasecmp( ext , "ply" ) ) pointStream = new PLYOrientedPointStreamWithData< float , Point3D< unsigned char > >( In.value , PlyColorProperties , 6 , ValidPlyColorProperties ); 387 else pointStream = new ASCIIOrientedPointStreamWithData< float , Point3D< unsigned char > >( In.value , ReadASCIIColor ); 388 pointCount = tree.template SetTree< float , NORMAL_DEGREE , WEIGHT_DEGREE , DATA_DEGREE , Point3D< unsigned char > >( pointStream , MinDepth.value , Depth.value , FullDepth.value , kernelDepth , Real(SamplesPerNode.value) , Scale.value , Confidence.set , NormalWeights.set , PointWeight.value , AdaptiveExponent.value , *densityWeights , *pointInfo , *normalInfo , *nodeWeights , colorData , xForm , Dirichlet.set , Complete.set ); 389 delete pointStream; 390 391 for( const OctNode< TreeNodeData >* n = tree.tree().nextNode() ; n!=NULL ; n=tree.tree().nextNode( n ) ) 392 { 393 int idx = colorData->index( n ); 394 if( idx>=0 ) colorData->data[idx] *= (Real)pow( Color.value , n->depth() ); 395 } 396 } 397 else 398 { 399 OrientedPointStream< float >* pointStream; 400 if ( !strcasecmp( ext , "bnpts" ) ) pointStream = new BinaryOrientedPointStream< float >( In.value ); 401 else if( !strcasecmp( ext , "ply" ) ) pointStream = new PLYOrientedPointStream< float >( In.value ); 402 else pointStream = new ASCIIOrientedPointStream< float >( In.value ); 403 pointCount = tree.template SetTree< float , NORMAL_DEGREE , WEIGHT_DEGREE , DATA_DEGREE , Point3D< unsigned char > >( pointStream , MinDepth.value , Depth.value , FullDepth.value , kernelDepth , Real(SamplesPerNode.value) , Scale.value , Confidence.set , NormalWeights.set , PointWeight.value , AdaptiveExponent.value , *densityWeights , *pointInfo , *normalInfo , *nodeWeights , colorData , xForm , Dirichlet.set , Complete.set ); 404 delete pointStream; 405 } 406 delete[] ext; 407 if( !Density.set ) delete densityWeights , densityWeights = NULL; 408 { 409 std::vector< int > indexMap; 410 if( NORMAL_DEGREE>Degree ) tree.template EnableMultigrid< NORMAL_DEGREE >( &indexMap ); 411 else tree.template EnableMultigrid< Degree >( &indexMap ); 412 if( pointInfo ) pointInfo->remapIndices( indexMap ); 413 if( normalInfo ) normalInfo->remapIndices( indexMap ); 414 if( densityWeights ) densityWeights->remapIndices( indexMap ); 415 if( nodeWeights ) nodeWeights->remapIndices( indexMap ); 416 if( colorData ) colorData->remapIndices( indexMap ); 417 } 418 419 DumpOutput2( comments , "# Tree set in: %9.1f (s), %9.1f (MB)\n" , Time()-t , tree.maxMemoryUsage ); 420 DumpOutput( "Input Points: %d\n" , pointCount ); 421 DumpOutput( "Leaves/Nodes: %d/%d\n" , tree.leaves() , tree.nodes() ); 422 DumpOutput( "Memory Usage: %.3f MB\n" , float( MemoryInfo::Usage() )/(1<<20) ); 423 424 maxMemoryUsage = tree.maxMemoryUsage; 425 t=Time() , tree.maxMemoryUsage=0; 426 DenseNodeData< Real , Degree > constraints = tree.template SetLaplacianConstraints< Degree >( *normalInfo ); 427 delete normalInfo; 428 DumpOutput2( comments , "# Constraints set in: %9.1f (s), %9.1f (MB)\n" , Time()-t , tree.maxMemoryUsage ); 429 DumpOutput( "Memory Usage: %.3f MB\n" , float( MemoryInfo::Usage())/(1<<20) ); 430 maxMemoryUsage = std::max< double >( maxMemoryUsage , tree.maxMemoryUsage ); 431 432 t=Time() , tree.maxMemoryUsage=0; 433 DenseNodeData< Real , Degree > solution = tree.SolveSystem( *pointInfo , constraints , ShowResidual.set , Iters.value , MaxSolveDepth.value , CGDepth.value , CSSolverAccuracy.value ); 434 delete pointInfo; 435 constraints.resize( 0 ); 436 437 DumpOutput2( comments , "# Linear system solved in: %9.1f (s), %9.1f (MB)\n" , Time()-t , tree.maxMemoryUsage ); 438 DumpOutput( "Memory Usage: %.3f MB\n" , float( MemoryInfo::Usage() )/(1<<20) ); 439 maxMemoryUsage = std::max< double >( maxMemoryUsage , tree.maxMemoryUsage ); 440 441 CoredFileMeshData< Vertex > mesh; 442 443 if( Verbose.set ) tree.maxMemoryUsage=0; 444 t=Time(); 445 isoValue = tree.GetIsoValue( solution , *nodeWeights ); 446 delete nodeWeights; 447 DumpOutput( "Got average in: %f\n" , Time()-t ); 448 DumpOutput( "Iso-Value: %e\n" , isoValue ); 449 450 if( VoxelGrid.set ) 451 { 452 double t = Time(); 453 FILE* fp = fopen( VoxelGrid.value , "wb" ); 454 if( !fp ) fprintf( stderr , "Failed to open voxel file for writing: %s\n" , VoxelGrid.value ); 455 else 456 { 457 int res = 0; 458 Pointer( Real ) values = tree.Evaluate( solution , res , isoValue , VoxelDepth.value , PrimalVoxel.set ); 459 fwrite( &res , sizeof(int) , 1 , fp ); 460 if( sizeof(Real)==sizeof(float) ) fwrite( values , sizeof(float) , res*res*res , fp ); 461 else 462 { 463 float *fValues = new float[res*res*res]; 464 for( int i=0 ; i<res*res*res ; i++ ) fValues[i] = float( values[i] ); 465 fwrite( fValues , sizeof(float) , res*res*res , fp ); 466 delete[] fValues; 467 } 468 fclose( fp ); 469 DeletePointer( values ); 470 } 471 DumpOutput( "Got voxel grid in: %f\n" , Time()-t ); 472 } 473 474 if( Out.set ) 475 { 476 t = Time() , tree.maxMemoryUsage = 0; 477 tree.template GetMCIsoSurface< Degree , WEIGHT_DEGREE , DATA_DEGREE >( densityWeights , colorData , solution , isoValue , mesh , !LinearFit.set , !NonManifold.set , PolygonMesh.set ); 478 if( PolygonMesh.set ) DumpOutput2( comments , "# Got polygons in: %9.1f (s), %9.1f (MB)\n" , Time()-t , tree.maxMemoryUsage ); 479 else DumpOutput2( comments , "# Got triangles in: %9.1f (s), %9.1f (MB)\n" , Time()-t , tree.maxMemoryUsage ); 480 maxMemoryUsage = std::max< double >( maxMemoryUsage , tree.maxMemoryUsage ); 481 DumpOutput2( comments , "# Total Solve: %9.1f (s), %9.1f (MB)\n" , Time()-tt , maxMemoryUsage ); 482 483 if( NoComments.set ) 484 { 485 if( ASCII.set ) PlyWritePolygons( Out.value , &mesh , PLY_ASCII , NULL , 0 , iXForm ); 486 else PlyWritePolygons( Out.value , &mesh , PLY_BINARY_NATIVE , NULL , 0 , iXForm ); 487 } 488 else 489 { 490 if( ASCII.set ) PlyWritePolygons( Out.value , &mesh , PLY_ASCII , &comments[0] , (int)comments.size() , iXForm ); 491 else PlyWritePolygons( Out.value , &mesh , PLY_BINARY_NATIVE , &comments[0] , (int)comments.size() , iXForm ); 492 } 493 DumpOutput( "Vertices / Polygons: %d / %d\n" , mesh.outOfCorePointCount()+mesh.inCorePoints.size() , mesh.polygonCount() ); 494 } 495 solution.resize( 0 ); 496 if( colorData ){ delete colorData ; colorData = NULL; } 497 return 1; 498} 499 500#ifdef _WIN32 501inline double to_seconds( const FILETIME& ft ) 502{ 503 const double low_to_sec=100e-9; // 100 nanoseconds 504 const double high_to_sec=low_to_sec*4294967296.0; 505 return ft.dwLowDateTime*low_to_sec+ft.dwHighDateTime*high_to_sec; 506} 507#endif // _WIN32 508 509template< class Real , class Vertex > 510int Execute( int argc , char* argv[] ) 511{ 512 switch( Degree.value ) 513 { 514 case 1: return _Execute< Real , 1 , Vertex >( argc , argv ); 515 case 2: return _Execute< Real , 2 , Vertex >( argc , argv ); 516 case 3: return _Execute< Real , 3 , Vertex >( argc , argv ); 517 case 4: return _Execute< Real , 4 , Vertex >( argc , argv ); 518 default: 519 fprintf( stderr , "[ERROR] Only B-Splines of degree 1 - 4 are supported" ); 520 return EXIT_FAILURE; 521 } 522} 523int main( int argc , char* argv[] ) 524{ 525#if defined(WIN32) && defined(MAX_MEMORY_GB) 526 if( MAX_MEMORY_GB>0 ) 527 { 528 SIZE_T peakMemory = 1; 529 peakMemory <<= 30; 530 peakMemory *= MAX_MEMORY_GB; 531 printf( "Limiting memory usage to %.2f GB\n" , float( peakMemory>>30 ) ); 532 HANDLE h = CreateJobObject( NULL , NULL ); 533 AssignProcessToJobObject( h , GetCurrentProcess() ); 534 535 JOBOBJECT_EXTENDED_LIMIT_INFORMATION jeli = { 0 }; 536 jeli.BasicLimitInformation.LimitFlags = JOB_OBJECT_LIMIT_JOB_MEMORY; 537 jeli.JobMemoryLimit = peakMemory; 538 if( !SetInformationJobObject( h , JobObjectExtendedLimitInformation , &jeli , sizeof( jeli ) ) ) 539 fprintf( stderr , "Failed to set memory limit\n" ); 540 } 541#endif // defined(WIN32) && defined(MAX_MEMORY_GB) 542 double t = Time(); 543 544 cmdLineParse( argc-1 , &argv[1] , sizeof(params)/sizeof(cmdLineReadable*) , params , 1 ); 545 if( Density.set ) 546 if( Color.set ) 547 if( Double.set ) Execute< double , PlyColorAndValueVertex< float > >( argc , argv ); 548 else Execute< float , PlyColorAndValueVertex< float > >( argc , argv ); 549 else 550 if( Double.set ) Execute< double , PlyValueVertex< float > >( argc , argv ); 551 else Execute< float , PlyValueVertex< float > >( argc , argv ); 552 else 553 if( Color.set ) 554 if( Double.set ) Execute< double , PlyColorVertex< float > >( argc , argv ); 555 else Execute< float , PlyColorVertex< float > >( argc , argv ); 556 else 557 if( Double.set ) Execute< double , PlyVertex< float > >( argc , argv ); 558 else Execute< float , PlyVertex< float > >( argc , argv ); 559#ifdef _WIN32 560 if( Performance.set ) 561 { 562 HANDLE cur_thread=GetCurrentThread(); 563 FILETIME tcreat, texit, tkernel, tuser; 564 if( GetThreadTimes( cur_thread , &tcreat , &texit , &tkernel , &tuser ) ) 565 printf( "Time (Wall/User/Kernel): %.2f / %.2f / %.2f\n" , Time()-t , to_seconds( tuser ) , to_seconds( tkernel ) ); 566 else printf( "Time: %.2f\n" , Time()-t ); 567 HANDLE h = GetCurrentProcess(); 568 PROCESS_MEMORY_COUNTERS pmc; 569 if( GetProcessMemoryInfo( h , &pmc , sizeof(pmc) ) ) printf( "Peak Memory (MB): %d\n" , pmc.PeakWorkingSetSize>>20 ); 570 } 571#endif // _WIN32 572 return EXIT_SUCCESS; 573} 574