mirror of
https://github.com/lucaspalomodevelop/meshlab.git
synced 2026-03-20 03:16:10 +00:00
Original version of the Misha code (version 6.11)
This commit is contained in:
Paolo Cignoni cignoni
parent
e4d15a4966
commit
a7c1787fb9
114
src/plugins_experimental/filter_screened_poisson/Src/Geometry.cpp
Executable file
114
src/plugins_experimental/filter_screened_poisson/Src/Geometry.cpp
Executable file
@ -0,0 +1,114 @@
|
||||
/*
|
||||
Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
|
||||
All rights reserved.
|
||||
|
||||
Redistribution and use in source and binary forms, with or without modification,
|
||||
are permitted provided that the following conditions are met:
|
||||
|
||||
Redistributions of source code must retain the above copyright notice, this list of
|
||||
conditions and the following disclaimer. Redistributions in binary form must reproduce
|
||||
the above copyright notice, this list of conditions and the following disclaimer
|
||||
in the documentation and/or other materials provided with the distribution.
|
||||
|
||||
Neither the name of the Johns Hopkins University nor the names of its contributors
|
||||
may be used to endorse or promote products derived from this software without specific
|
||||
prior written permission.
|
||||
|
||||
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
|
||||
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
|
||||
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
|
||||
SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
|
||||
INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
|
||||
TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
|
||||
BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
|
||||
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
|
||||
ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
|
||||
DAMAGE.
|
||||
*/
|
||||
#include "Geometry.h"
|
||||
#include <stdio.h>
|
||||
#include <string.h>
|
||||
|
||||
///////////////////
|
||||
// CoredMeshData //
|
||||
///////////////////
|
||||
|
||||
TriangulationEdge::TriangulationEdge(void){pIndex[0]=pIndex[1]=tIndex[0]=tIndex[1]=-1;}
|
||||
TriangulationTriangle::TriangulationTriangle(void){eIndex[0]=eIndex[1]=eIndex[2]=-1;}
|
||||
|
||||
///////////////////////////
|
||||
// BufferedReadWriteFile //
|
||||
///////////////////////////
|
||||
BufferedReadWriteFile::BufferedReadWriteFile( char* fileName , int bufferSize )
|
||||
{
|
||||
_bufferIndex = 0;
|
||||
_bufferSize = bufferSize;
|
||||
if( fileName ) strcpy( _fileName , fileName ) , tempFile = false;
|
||||
#ifdef _WIN32
|
||||
else strcpy( _fileName , _tempnam( "." , "foo" ) ) , tempFile = true;
|
||||
#else // !_WIN32
|
||||
else strcpy( _fileName , tempnam( "." , "foo" ) ) , tempFile = true;
|
||||
#endif // _WIN32
|
||||
_fp = fopen( _fileName , "w+b" );
|
||||
if( !_fp ) fprintf( stderr , "[ERROR] Failed to open file: %s\n" , _fileName ) , exit( 0 );
|
||||
_buffer = (char*) malloc( _bufferSize );
|
||||
}
|
||||
BufferedReadWriteFile::~BufferedReadWriteFile( void )
|
||||
{
|
||||
free( _buffer );
|
||||
fclose( _fp );
|
||||
if( tempFile ) remove( _fileName );
|
||||
}
|
||||
void BufferedReadWriteFile::reset( void )
|
||||
{
|
||||
if( _bufferIndex ) fwrite( _buffer , 1 , _bufferIndex , _fp );
|
||||
_bufferIndex = 0;
|
||||
fseek( _fp , 0 , SEEK_SET );
|
||||
_bufferIndex = 0;
|
||||
_bufferSize = fread( _buffer , 1 , _bufferSize , _fp );
|
||||
}
|
||||
bool BufferedReadWriteFile::write( const void* data , size_t size )
|
||||
{
|
||||
if( !size ) return true;
|
||||
char* _data = (char*) data;
|
||||
size_t sz = _bufferSize - _bufferIndex;
|
||||
while( sz<=size )
|
||||
{
|
||||
memcpy( _buffer+_bufferIndex , _data , sz );
|
||||
fwrite( _buffer , 1 , _bufferSize , _fp );
|
||||
_data += sz;
|
||||
size -= sz;
|
||||
_bufferIndex = 0;
|
||||
sz = _bufferSize;
|
||||
}
|
||||
if( size )
|
||||
{
|
||||
memcpy( _buffer+_bufferIndex , _data , size );
|
||||
_bufferIndex += size;
|
||||
}
|
||||
return true;
|
||||
}
|
||||
bool BufferedReadWriteFile::read( void* data , size_t size )
|
||||
{
|
||||
if( !size ) return true;
|
||||
char *_data = (char*) data;
|
||||
size_t sz = _bufferSize - _bufferIndex;
|
||||
while( sz<=size )
|
||||
{
|
||||
if( size && !_bufferSize ) return false;
|
||||
memcpy( _data , _buffer+_bufferIndex , sz );
|
||||
_bufferSize = fread( _buffer , 1 , _bufferSize , _fp );
|
||||
_data += sz;
|
||||
size -= sz;
|
||||
_bufferIndex = 0;
|
||||
if( !size ) return true;
|
||||
sz = _bufferSize;
|
||||
}
|
||||
if( size )
|
||||
{
|
||||
if( !_bufferSize ) return false;
|
||||
memcpy( _data , _buffer+_bufferIndex , size );
|
||||
_bufferIndex += size;
|
||||
}
|
||||
return true;
|
||||
}
|
||||
213
src/plugins_experimental/filter_screened_poisson/Src/MAT.inl
Executable file
213
src/plugins_experimental/filter_screened_poisson/Src/MAT.inl
Executable file
@ -0,0 +1,213 @@
|
||||
/*
|
||||
Copyright (c) 2007, Michael Kazhdan
|
||||
All rights reserved.
|
||||
|
||||
Redistribution and use in source and binary forms, with or without modification,
|
||||
are permitted provided that the following conditions are met:
|
||||
|
||||
Redistributions of source code must retain the above copyright notice, this list of
|
||||
conditions and the following disclaimer. Redistributions in binary form must reproduce
|
||||
the above copyright notice, this list of conditions and the following disclaimer
|
||||
in the documentation and/or other materials provided with the distribution.
|
||||
|
||||
Neither the name of the Johns Hopkins University nor the names of its contributors
|
||||
may be used to endorse or promote products derived from this software without specific
|
||||
prior written permission.
|
||||
|
||||
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
|
||||
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
|
||||
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
|
||||
SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
|
||||
INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
|
||||
TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
|
||||
BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
|
||||
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
|
||||
ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
|
||||
DAMAGE.
|
||||
*/
|
||||
//////////////////////////////
|
||||
// MinimalAreaTriangulation //
|
||||
//////////////////////////////
|
||||
template <class Real>
|
||||
MinimalAreaTriangulation<Real>::MinimalAreaTriangulation(void)
|
||||
{
|
||||
bestTriangulation=NULL;
|
||||
midPoint=NULL;
|
||||
}
|
||||
template <class Real>
|
||||
MinimalAreaTriangulation<Real>::~MinimalAreaTriangulation(void)
|
||||
{
|
||||
if(bestTriangulation)
|
||||
delete[] bestTriangulation;
|
||||
bestTriangulation=NULL;
|
||||
if(midPoint)
|
||||
delete[] midPoint;
|
||||
midPoint=NULL;
|
||||
}
|
||||
template <class Real>
|
||||
void MinimalAreaTriangulation<Real>::GetTriangulation(const std::vector<Point3D<Real> >& vertices,std::vector<TriangleIndex>& triangles)
|
||||
{
|
||||
if(vertices.size()==3)
|
||||
{
|
||||
triangles.resize(1);
|
||||
triangles[0].idx[0]=0;
|
||||
triangles[0].idx[1]=1;
|
||||
triangles[0].idx[2]=2;
|
||||
return;
|
||||
}
|
||||
else if(vertices.size()==4)
|
||||
{
|
||||
TriangleIndex tIndex[2][2];
|
||||
Real area[2];
|
||||
|
||||
area[0]=area[1]=0;
|
||||
triangles.resize(2);
|
||||
|
||||
tIndex[0][0].idx[0]=0;
|
||||
tIndex[0][0].idx[1]=1;
|
||||
tIndex[0][0].idx[2]=2;
|
||||
tIndex[0][1].idx[0]=2;
|
||||
tIndex[0][1].idx[1]=3;
|
||||
tIndex[0][1].idx[2]=0;
|
||||
|
||||
tIndex[1][0].idx[0]=0;
|
||||
tIndex[1][0].idx[1]=1;
|
||||
tIndex[1][0].idx[2]=3;
|
||||
tIndex[1][1].idx[0]=3;
|
||||
tIndex[1][1].idx[1]=1;
|
||||
tIndex[1][1].idx[2]=2;
|
||||
|
||||
Point3D<Real> n,p1,p2;
|
||||
for(int i=0;i<2;i++)
|
||||
for(int j=0;j<2;j++)
|
||||
{
|
||||
p1=vertices[tIndex[i][j].idx[1]]-vertices[tIndex[i][j].idx[0]];
|
||||
p2=vertices[tIndex[i][j].idx[2]]-vertices[tIndex[i][j].idx[0]];
|
||||
CrossProduct(p1,p2,n);
|
||||
area[i] += Real( Length(n) );
|
||||
}
|
||||
if(area[0]>area[1])
|
||||
{
|
||||
triangles[0]=tIndex[1][0];
|
||||
triangles[1]=tIndex[1][1];
|
||||
}
|
||||
else
|
||||
{
|
||||
triangles[0]=tIndex[0][0];
|
||||
triangles[1]=tIndex[0][1];
|
||||
}
|
||||
return;
|
||||
}
|
||||
if(bestTriangulation)
|
||||
delete[] bestTriangulation;
|
||||
if(midPoint)
|
||||
delete[] midPoint;
|
||||
bestTriangulation=NULL;
|
||||
midPoint=NULL;
|
||||
size_t eCount=vertices.size();
|
||||
bestTriangulation=new Real[eCount*eCount];
|
||||
midPoint=new int[eCount*eCount];
|
||||
for(size_t i=0;i<eCount*eCount;i++)
|
||||
bestTriangulation[i]=-1;
|
||||
memset(midPoint,-1,sizeof(int)*eCount*eCount);
|
||||
GetArea(0,1,vertices);
|
||||
triangles.clear();
|
||||
GetTriangulation(0,1,vertices,triangles);
|
||||
}
|
||||
template <class Real>
|
||||
Real MinimalAreaTriangulation<Real>::GetArea(const std::vector<Point3D<Real> >& vertices)
|
||||
{
|
||||
if(bestTriangulation)
|
||||
delete[] bestTriangulation;
|
||||
if(midPoint)
|
||||
delete[] midPoint;
|
||||
bestTriangulation=NULL;
|
||||
midPoint=NULL;
|
||||
int eCount=vertices.size();
|
||||
bestTriangulation=new double[eCount*eCount];
|
||||
midPoint=new int[eCount*eCount];
|
||||
for(int i=0;i<eCount*eCount;i++)
|
||||
bestTriangulation[i]=-1;
|
||||
memset(midPoint,-1,sizeof(int)*eCount*eCount);
|
||||
return GetArea(0,1,vertices);
|
||||
}
|
||||
template<class Real>
|
||||
void MinimalAreaTriangulation<Real>::GetTriangulation(const size_t& i,const size_t& j,const std::vector<Point3D<Real> >& vertices,std::vector<TriangleIndex>& triangles)
|
||||
{
|
||||
TriangleIndex tIndex;
|
||||
size_t eCount=vertices.size();
|
||||
size_t ii=i;
|
||||
if(i<j)
|
||||
ii+=eCount;
|
||||
if(j+1>=ii)
|
||||
return;
|
||||
ii=midPoint[i*eCount+j];
|
||||
if(ii>=0)
|
||||
{
|
||||
tIndex.idx[0] = int( i );
|
||||
tIndex.idx[1] = int( j );
|
||||
tIndex.idx[2] = int( ii );
|
||||
triangles.push_back(tIndex);
|
||||
GetTriangulation(i,ii,vertices,triangles);
|
||||
GetTriangulation(ii,j,vertices,triangles);
|
||||
}
|
||||
}
|
||||
|
||||
template<class Real>
|
||||
Real MinimalAreaTriangulation<Real>::GetArea(const size_t& i,const size_t& j,const std::vector<Point3D<Real> >& vertices)
|
||||
{
|
||||
Real a=FLT_MAX,temp;
|
||||
size_t eCount=vertices.size();
|
||||
size_t idx=i*eCount+j;
|
||||
size_t ii=i;
|
||||
if(i<j)
|
||||
ii+=eCount;
|
||||
if(j+1>=ii)
|
||||
{
|
||||
bestTriangulation[idx]=0;
|
||||
return 0;
|
||||
}
|
||||
if(midPoint[idx]!=-1)
|
||||
return bestTriangulation[idx];
|
||||
int mid=-1;
|
||||
for(size_t r=j+1;r<ii;r++)
|
||||
{
|
||||
size_t rr=r%eCount;
|
||||
size_t idx1=i*eCount+rr,idx2=rr*eCount+j;
|
||||
Point3D<Real> p,p1,p2;
|
||||
p1=vertices[i]-vertices[rr];
|
||||
p2=vertices[j]-vertices[rr];
|
||||
CrossProduct(p1,p2,p);
|
||||
temp = Real( Length(p) );
|
||||
if(bestTriangulation[idx1]>=0)
|
||||
{
|
||||
temp+=bestTriangulation[idx1];
|
||||
if(temp>a)
|
||||
continue;
|
||||
if(bestTriangulation[idx2]>0)
|
||||
temp+=bestTriangulation[idx2];
|
||||
else
|
||||
temp+=GetArea(rr,j,vertices);
|
||||
}
|
||||
else
|
||||
{
|
||||
if(bestTriangulation[idx2]>=0)
|
||||
temp+=bestTriangulation[idx2];
|
||||
else
|
||||
temp+=GetArea(rr,j,vertices);
|
||||
if(temp>a)
|
||||
continue;
|
||||
temp+=GetArea(i,rr,vertices);
|
||||
}
|
||||
|
||||
if(temp<a)
|
||||
{
|
||||
a=temp;
|
||||
mid=int(rr);
|
||||
}
|
||||
}
|
||||
bestTriangulation[idx]=a;
|
||||
midPoint[idx]=mid;
|
||||
|
||||
return a;
|
||||
}
|
||||
200
src/plugins_experimental/filter_screened_poisson/Src/MemoryUsage.h
Executable file
200
src/plugins_experimental/filter_screened_poisson/Src/MemoryUsage.h
Executable file
@ -0,0 +1,200 @@
|
||||
/*
|
||||
Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
|
||||
All rights reserved.
|
||||
|
||||
Redistribution and use in source and binary forms, with or without modification,
|
||||
are permitted provided that the following conditions are met:
|
||||
|
||||
Redistributions of source code must retain the above copyright notice, this list of
|
||||
conditions and the following disclaimer. Redistributions in binary form must reproduce
|
||||
the above copyright notice, this list of conditions and the following disclaimer
|
||||
in the documentation and/or other materials provided with the distribution.
|
||||
|
||||
Neither the name of the Johns Hopkins University nor the names of its contributors
|
||||
may be used to endorse or promote products derived from this software without specific
|
||||
prior written permission.
|
||||
|
||||
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
|
||||
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
|
||||
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
|
||||
SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
|
||||
INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
|
||||
TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
|
||||
BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
|
||||
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
|
||||
ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
|
||||
DAMAGE.
|
||||
*/
|
||||
|
||||
#ifndef MEMORY_USAGE_INCLUDED
|
||||
#define MEMORY_USAGE_INCLUDED
|
||||
|
||||
#ifdef WIN32
|
||||
|
||||
#include <Windows.h>
|
||||
class MemoryInfo
|
||||
{
|
||||
public:
|
||||
size_t TotalPhysicalMemory;
|
||||
size_t FreePhysicalMemory;
|
||||
size_t TotalSwapSpace;
|
||||
size_t FreeSwapSpace;
|
||||
size_t TotalVirtualAddressSpace;
|
||||
size_t FreeVirtualAddressSpace;
|
||||
size_t PageSize;
|
||||
|
||||
void set(void){
|
||||
MEMORYSTATUSEX Mem;
|
||||
SYSTEM_INFO Info;
|
||||
ZeroMemory( &Mem, sizeof(Mem));
|
||||
ZeroMemory( &Info, sizeof(Info));
|
||||
Mem.dwLength = sizeof(Mem);
|
||||
::GlobalMemoryStatusEx( &Mem );
|
||||
::GetSystemInfo( &Info );
|
||||
|
||||
TotalPhysicalMemory = (size_t)Mem.ullTotalPhys;
|
||||
FreePhysicalMemory = (size_t)Mem.ullAvailPhys;
|
||||
TotalSwapSpace = (size_t)Mem.ullTotalPageFile;
|
||||
FreeSwapSpace = (size_t)Mem.ullAvailPageFile;
|
||||
TotalVirtualAddressSpace = (size_t)Mem.ullTotalVirtual;
|
||||
FreeVirtualAddressSpace = (size_t)Mem.ullAvailVirtual;
|
||||
PageSize = (size_t)Info.dwPageSize;
|
||||
}
|
||||
size_t usage(void) const {return TotalVirtualAddressSpace-FreeVirtualAddressSpace;}
|
||||
|
||||
static size_t Usage(void){
|
||||
MEMORY_BASIC_INFORMATION mbi;
|
||||
size_t dwMemUsed = 0;
|
||||
PVOID pvAddress = 0;
|
||||
|
||||
|
||||
memset(&mbi, 0, sizeof(MEMORY_BASIC_INFORMATION));
|
||||
while(VirtualQuery(pvAddress, &mbi, sizeof(MEMORY_BASIC_INFORMATION)) == sizeof(MEMORY_BASIC_INFORMATION)){
|
||||
if(mbi.State == MEM_COMMIT && mbi.Type == MEM_PRIVATE){dwMemUsed += mbi.RegionSize;}
|
||||
pvAddress = ((BYTE*)mbi.BaseAddress) + mbi.RegionSize;
|
||||
}
|
||||
return dwMemUsed;
|
||||
}
|
||||
};
|
||||
|
||||
#else // !WIN32
|
||||
|
||||
#ifndef __APPLE__ // Linux variants
|
||||
|
||||
#include <sys/time.h>
|
||||
#include <sys/resource.h>
|
||||
|
||||
class MemoryInfo
|
||||
{
|
||||
public:
|
||||
static size_t Usage(void)
|
||||
{
|
||||
FILE* f = fopen("/proc/self/stat","rb");
|
||||
|
||||
int d;
|
||||
long ld;
|
||||
unsigned long lu;
|
||||
unsigned long long llu;
|
||||
char* s;
|
||||
char c;
|
||||
|
||||
int pid;
|
||||
unsigned long vm;
|
||||
|
||||
int n = fscanf(f, "%d %s %c %d %d %d %d %d %lu %lu %lu %lu %lu %lu %lu %ld %ld %ld %ld %d %ld %llu %lu %ld %lu %lu %lu %lu %lu %lu %lu %lu %lu %lu %lu %lu %lu %d %d %lu %lu"
|
||||
,&pid ,&s ,&c ,&d ,&d ,&d ,&d ,&d ,&lu ,&lu ,&lu ,&lu ,&lu ,&lu ,&lu ,&ld ,&ld ,&ld ,&ld ,&d ,&ld ,&llu ,&vm ,&ld ,&lu ,&lu ,&lu ,&lu ,&lu ,&lu ,&lu ,&lu ,&lu ,&lu ,&lu ,&lu ,&lu ,&d ,&d ,&lu ,&lu );
|
||||
|
||||
fclose(f);
|
||||
/*
|
||||
pid %d
|
||||
comm %s
|
||||
state %c
|
||||
ppid %d
|
||||
pgrp %d
|
||||
session %d
|
||||
tty_nr %d
|
||||
tpgid %d
|
||||
flags %lu
|
||||
minflt %lu
|
||||
cminflt %lu
|
||||
majflt %lu
|
||||
cmajflt %lu
|
||||
utime %lu
|
||||
stime %lu
|
||||
cutime %ld
|
||||
cstime %ld
|
||||
priority %ld
|
||||
nice %ld
|
||||
0 %ld
|
||||
itrealvalue %ld
|
||||
starttime %lu
|
||||
vsize %lu
|
||||
rss %ld
|
||||
rlim %lu
|
||||
startcode %lu
|
||||
endcode %lu
|
||||
startstack %lu
|
||||
kstkesp %lu
|
||||
kstkeip %lu
|
||||
signal %lu
|
||||
blocked %lu
|
||||
sigignore %lu
|
||||
sigcatch %lu
|
||||
wchan %lu
|
||||
nswap %lu
|
||||
cnswap %lu
|
||||
exit_signal %d
|
||||
processor %d
|
||||
rt_priority %lu (since kernel 2.5.19)
|
||||
policy %lu (since kernel 2.5.19)
|
||||
*/
|
||||
return vm;
|
||||
}
|
||||
|
||||
};
|
||||
#else // __APPLE__: has no "/proc" pseudo-file system
|
||||
|
||||
// Thanks to David O'Gwynn for providing this fix.
|
||||
// This comes from a post by Michael Knight:
|
||||
//
|
||||
// http://miknight.blogspot.com/2005/11/resident-set-size-in-mac-os-x.html
|
||||
|
||||
#include <unistd.h>
|
||||
#include <stdio.h>
|
||||
#include <stdlib.h>
|
||||
#include <sys/types.h>
|
||||
#include <sys/sysctl.h>
|
||||
#include <mach/task.h>
|
||||
#include <mach/mach_init.h>
|
||||
|
||||
void getres(task_t task, unsigned long *rss, unsigned long *vs)
|
||||
{
|
||||
struct task_basic_info t_info;
|
||||
mach_msg_type_number_t t_info_count = TASK_BASIC_INFO_COUNT;
|
||||
|
||||
task_info(task, TASK_BASIC_INFO, (task_info_t)&t_info, &t_info_count);
|
||||
*rss = t_info.resident_size;
|
||||
*vs = t_info.virtual_size;
|
||||
}
|
||||
|
||||
class MemoryInfo
|
||||
{
|
||||
public:
|
||||
static size_t Usage(void)
|
||||
{
|
||||
unsigned long rss, vs, psize;
|
||||
task_t task = MACH_PORT_NULL;
|
||||
|
||||
if (task_for_pid(current_task(), getpid(), &task) != KERN_SUCCESS)
|
||||
abort();
|
||||
getres(task, &rss, &vs);
|
||||
return rss;
|
||||
}
|
||||
|
||||
};
|
||||
|
||||
#endif // !__APPLE__
|
||||
|
||||
#endif // WIN32
|
||||
|
||||
#endif // MEMORY_USAGE_INCLUDE
|
||||
File diff suppressed because it is too large
Load Diff
438
src/plugins_experimental/filter_screened_poisson/Src/MultiGridOctreeData.h
Executable file
438
src/plugins_experimental/filter_screened_poisson/Src/MultiGridOctreeData.h
Executable file
@ -0,0 +1,438 @@
|
||||
/*
|
||||
Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
|
||||
All rights reserved.
|
||||
|
||||
Redistribution and use in source and binary forms, with or without modification,
|
||||
are permitted provided that the following conditions are met:
|
||||
|
||||
Redistributions of source code must retain the above copyright notice, this list of
|
||||
conditions and the following disclaimer. Redistributions in binary form must reproduce
|
||||
the above copyright notice, this list of conditions and the following disclaimer
|
||||
in the documentation and/or other materials provided with the distribution.
|
||||
|
||||
Neither the name of the Johns Hopkins University nor the names of its contributors
|
||||
may be used to endorse or promote products derived from this software without specific
|
||||
prior written permission.
|
||||
|
||||
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
|
||||
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
|
||||
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
|
||||
SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
|
||||
INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
|
||||
TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
|
||||
BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
|
||||
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
|
||||
ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
|
||||
DAMAGE.
|
||||
*/
|
||||
|
||||
#ifndef MULTI_GRID_OCTREE_DATA_INCLUDED
|
||||
#define MULTI_GRID_OCTREE_DATA_INCLUDED
|
||||
|
||||
#define NEW_CODE 1
|
||||
|
||||
#define MAX_MEMORY_GB 15
|
||||
|
||||
#define GRADIENT_DOMAIN_SOLUTION 1 // Given the constraint vector-field V(p), there are two ways to solve for the coefficients, x, of the indicator function
|
||||
// with respect to the B-spline basis {B_i(p)}
|
||||
// 1] Find x minimizing:
|
||||
// || V(p) - \sum_i \nabla x_i B_i(p) ||^2
|
||||
// which is solved by the system A_1x = b_1 where:
|
||||
// A_1[i,j] = < \nabla B_i(p) , \nabla B_j(p) >
|
||||
// b_1[i] = < \nabla B_i(p) , V(p) >
|
||||
// 2] Formulate this as a Poisson equation:
|
||||
// \sum_i x_i \Delta B_i(p) = \nabla \cdot V(p)
|
||||
// which is solved by the system A_2x = b_2 where:
|
||||
// A_2[i,j] = - < \Delta B_i(p) , B_j(p) >
|
||||
// b_2[i] = - < B_i(p) , \nabla \cdot V(p) >
|
||||
// Although the two system matrices should be the same (assuming that the B_i satisfy dirichlet/neumann boundary conditions)
|
||||
// the constraint vectors can differ when V does not satisfy the Neumann boundary conditions:
|
||||
// A_1[i,j] = \int_R < \nabla B_i(p) , \nabla B_j(p) >
|
||||
// = \int_R [ \nabla \cdot ( B_i(p) \nabla B_j(p) ) - B_i(p) \Delta B_j(p) ]
|
||||
// = \int_dR < N(p) , B_i(p) \nabla B_j(p) > + A_2[i,j]
|
||||
// and the first integral is zero if either f_i is zero on the boundary dR or the derivative of B_i across the boundary is zero.
|
||||
// However, for the constraints we have:
|
||||
// b_1(i) = \int_R < \nabla B_i(p) , V(p) >
|
||||
// = \int_R [ \nabla \cdot ( B_i(p) V(p) ) - B_i(p) \nabla \cdot V(p) ]
|
||||
// = \int_dR < N(p) , B_i(p) V(p) > + b_2[i]
|
||||
// In particular, this implies that if the B_i satisfy the Neumann boundary conditions (rather than Dirichlet),
|
||||
// and V is not zero across the boundary, then the two constraints are different.
|
||||
// Forcing the < V(p) , N(p) > = 0 on the boundary, by killing off the component of the vector-field in the normal direction
|
||||
// (FORCE_NEUMANN_FIELD), makes the two systems equal, and the value of this flag should be immaterial.
|
||||
// Note that under interpretation 1, we have:
|
||||
// \sum_i b_1(i) = < \nabla \sum_ i B_i(p) , V(p) > = 0
|
||||
// because the B_i's sum to one. However, in general, we could have
|
||||
// \sum_i b_2(i) \neq 0.
|
||||
// This could cause trouble because the constant functions are in the kernel of the matrix A, so CG will misbehave if the constraint
|
||||
// has a non-zero DC term. (Again, forcing < V(p) , N(p) > = 0 along the boundary resolves this problem.)
|
||||
|
||||
#define FORCE_NEUMANN_FIELD 1 // This flag forces the normal component across the boundary of the integration domain to be zero.
|
||||
// This should be enabled if GRADIENT_DOMAIN_SOLUTION is not, so that CG doesn't run into trouble.
|
||||
|
||||
#define ROBERTO_TOLDO_FIX 1
|
||||
|
||||
#if !FORCE_NEUMANN_FIELD
|
||||
#pragma message( "[WARNING] Not zeroing out normal component on boundary" )
|
||||
#endif // !FORCE_NEUMANN_FIELD
|
||||
|
||||
#include "Hash.h"
|
||||
#include "BSplineData.h"
|
||||
|
||||
class TreeNodeData
|
||||
{
|
||||
public:
|
||||
static int NodeCount;
|
||||
int nodeIndex;
|
||||
|
||||
TreeNodeData( void );
|
||||
~TreeNodeData( void );
|
||||
};
|
||||
|
||||
class VertexData
|
||||
{
|
||||
typedef OctNode< TreeNodeData > TreeOctNode;
|
||||
public:
|
||||
static const int VERTEX_COORDINATE_SHIFT = ( sizeof( long long ) * 8 ) / 3;
|
||||
static long long EdgeIndex( const TreeOctNode* node , int eIndex , int maxDepth , int index[DIMENSION] );
|
||||
static long long EdgeIndex( const TreeOctNode* node , int eIndex , int maxDepth );
|
||||
static long long FaceIndex( const TreeOctNode* node , int fIndex , int maxDepth,int index[DIMENSION] );
|
||||
static long long FaceIndex( const TreeOctNode* node , int fIndex , int maxDepth );
|
||||
static long long CornerIndex( const TreeOctNode* node , int cIndex , int maxDepth , int index[DIMENSION] );
|
||||
static long long CornerIndex( const TreeOctNode* node , int cIndex , int maxDepth );
|
||||
static long long CenterIndex( const TreeOctNode* node , int maxDepth , int index[DIMENSION] );
|
||||
static long long CenterIndex( const TreeOctNode* node , int maxDepth );
|
||||
static long long CornerIndex( int depth , const int offSet[DIMENSION] , int cIndex , int maxDepth , int index[DIMENSION] );
|
||||
static long long CenterIndex( int depth , const int offSet[DIMENSION] , int maxDepth , int index[DIMENSION] );
|
||||
static long long CornerIndexKey( const int index[DIMENSION] );
|
||||
};
|
||||
|
||||
class SortedTreeNodes
|
||||
{
|
||||
typedef OctNode< TreeNodeData > TreeOctNode;
|
||||
protected:
|
||||
void _sortByZCoordinate( void );
|
||||
public:
|
||||
Pointer( TreeOctNode* ) treeNodes;
|
||||
int *nodeCount;
|
||||
int maxDepth;
|
||||
SortedTreeNodes( void );
|
||||
~SortedTreeNodes( void );
|
||||
void set( TreeOctNode& root , std::vector< int >* map );
|
||||
Pointer( Pointer( int ) ) sliceOffsets;
|
||||
static int Slices( int depth );
|
||||
std::pair< int , int > sliceSpan( int depth , int off , int d ) const;
|
||||
|
||||
template< int Indices >
|
||||
struct _Indices
|
||||
{
|
||||
int idx[Indices];
|
||||
_Indices( void ){ memset( idx , -1 , sizeof( int ) * Indices ); }
|
||||
int& operator[] ( int i ) { return idx[i]; }
|
||||
const int& operator[] ( int i ) const { return idx[i]; }
|
||||
};
|
||||
typedef _Indices< Square::CORNERS > SquareCornerIndices;
|
||||
typedef _Indices< Square::EDGES > SquareEdgeIndices;
|
||||
typedef _Indices< Square::FACES > SquareFaceIndices;
|
||||
|
||||
struct SliceTableData
|
||||
{
|
||||
std::vector< SquareCornerIndices > cTable;
|
||||
std::vector< SquareEdgeIndices > eTable;
|
||||
std::vector< SquareFaceIndices > fTable;
|
||||
int cCount , eCount , fCount , nodeOffset , nodeCount;
|
||||
SliceTableData( void ){ fCount , eCount = cCount = 0; }
|
||||
~SliceTableData( void ){ clear(); }
|
||||
void clear( void ) { cTable.clear() , eTable.clear() , fTable.clear() , fCount = eCount = cCount = 0; }
|
||||
SquareCornerIndices& cornerIndices( const TreeOctNode* node );
|
||||
SquareCornerIndices& cornerIndices( int idx );
|
||||
const SquareCornerIndices& cornerIndices( const TreeOctNode* node ) const;
|
||||
const SquareCornerIndices& cornerIndices( int idx ) const;
|
||||
SquareEdgeIndices& edgeIndices( const TreeOctNode* node );
|
||||
SquareEdgeIndices& edgeIndices( int idx );
|
||||
const SquareEdgeIndices& edgeIndices( const TreeOctNode* node ) const;
|
||||
const SquareEdgeIndices& edgeIndices( int idx ) const;
|
||||
SquareFaceIndices& faceIndices( const TreeOctNode* node );
|
||||
SquareFaceIndices& faceIndices( int idx );
|
||||
const SquareFaceIndices& faceIndices( const TreeOctNode* node ) const;
|
||||
const SquareFaceIndices& faceIndices( int idx ) const;
|
||||
protected:
|
||||
std::vector< int > _cMap , _eMap , _fMap;
|
||||
friend class SortedTreeNodes;
|
||||
};
|
||||
struct XSliceTableData
|
||||
{
|
||||
std::vector< SquareCornerIndices > eTable;
|
||||
std::vector< SquareEdgeIndices > fTable;
|
||||
int fCount , eCount , nodeOffset , nodeCount;
|
||||
XSliceTableData( void ){ fCount , eCount = 0; }
|
||||
~XSliceTableData( void ){ clear(); }
|
||||
void clear( void ) { fTable.clear() , eTable.clear() , fCount = eCount = 0; }
|
||||
SquareCornerIndices& edgeIndices( const TreeOctNode* node );
|
||||
SquareCornerIndices& edgeIndices( int idx );
|
||||
const SquareCornerIndices& edgeIndices( const TreeOctNode* node ) const;
|
||||
const SquareCornerIndices& edgeIndices( int idx ) const;
|
||||
SquareEdgeIndices& faceIndices( const TreeOctNode* node );
|
||||
SquareEdgeIndices& faceIndices( int idx );
|
||||
const SquareEdgeIndices& faceIndices( const TreeOctNode* node ) const;
|
||||
const SquareEdgeIndices& faceIndices( int idx ) const;
|
||||
protected:
|
||||
std::vector< int > _eMap , _fMap;
|
||||
friend class SortedTreeNodes;
|
||||
};
|
||||
void setSliceTableData ( SliceTableData& sData , int depth , int offset , int threads ) const;
|
||||
void setXSliceTableData( XSliceTableData& sData , int depth , int offset , int threads ) const;
|
||||
};
|
||||
|
||||
|
||||
template< class Real , int Degree >
|
||||
class Octree
|
||||
{
|
||||
typedef OctNode< TreeNodeData > TreeOctNode;
|
||||
struct _PointData
|
||||
{
|
||||
Point3D< Real > position;
|
||||
Real weightedCoarserValue;
|
||||
Real weight;
|
||||
_PointData( Point3D< Real > p=Point3D< Real >() , Real w=0 ) { position = p , weight = w , weightedCoarserValue = Real(0); }
|
||||
};
|
||||
public:
|
||||
struct NormalInfo
|
||||
{
|
||||
std::vector< int > normalIndices;
|
||||
std::vector< Point3D< Real > > normals;
|
||||
int normalIndex( const TreeOctNode* node ) const { return node->nodeData.nodeIndex>=normalIndices.size() ? -1 : normalIndices[ node->nodeData.nodeIndex ]; }
|
||||
};
|
||||
struct PointInfo
|
||||
{
|
||||
std::vector< int > pointIndices;
|
||||
std::vector< _PointData > points;
|
||||
int pointIndex( const TreeOctNode* node ) const { return node->nodeData.nodeIndex>=pointIndices.size() ? -1 : pointIndices[ node->nodeData.nodeIndex ]; }
|
||||
};
|
||||
protected:
|
||||
SortedTreeNodes _sNodes;
|
||||
Real _samplesPerNode;
|
||||
int _splatDepth;
|
||||
int _minDepth;
|
||||
int _fullDepth;
|
||||
bool _constrainValues;
|
||||
int _boundaryType;
|
||||
Real _scale;
|
||||
Point3D< Real > _center;
|
||||
std::vector< int > _pointCount;
|
||||
Real _normalSmooth;
|
||||
BSplineData< Degree > _fData;
|
||||
|
||||
bool _InBounds( Point3D< Real > ) const;
|
||||
|
||||
double GetLaplacian ( const typename BSplineData< Degree >::Integrator& integrator , int d , const int off1[3] , const int off2[3] , bool childParent ) const;
|
||||
double GetDivergence1( const typename BSplineData< Degree >::Integrator& integrator , int d , const int off1[3] , const int off2[3] , bool childParent , const Point3D< Real >& normal1 ) const;
|
||||
double GetDivergence2( const typename BSplineData< Degree >::Integrator& integrator , int d , const int off1[3] , const int off2[3] , bool childParent , const Point3D< Real >& normal2 ) const;
|
||||
Point3D< double > GetDivergence1( const typename BSplineData< Degree >::Integrator& integrator , int d , const int off1[3] , const int off2[3] , bool childParent ) const;
|
||||
Point3D< double > GetDivergence2( const typename BSplineData< Degree >::Integrator& integrator , int d , const int off1[3] , const int off2[3] , bool childParent ) const;
|
||||
|
||||
template< class C , int N > struct Stencil{ C values[N][N][N]; };
|
||||
struct CenterValueStencil
|
||||
{
|
||||
Stencil< double , 3 > stencil;
|
||||
Stencil< double , 3 > stencils[8];
|
||||
};
|
||||
struct CornerValueStencil
|
||||
{
|
||||
Stencil< double , 3 > stencil[8];
|
||||
Stencil< double , 3 > stencils[8][8];
|
||||
};
|
||||
struct CornerNormalStencil
|
||||
{
|
||||
Stencil< Point3D< double > , 3 > stencil[8];
|
||||
Stencil< Point3D< double > , 3 > stencils[8][8];
|
||||
};
|
||||
|
||||
void _setMultiColorIndices( int start , int end , std::vector< std::vector< int > >& indices ) const;
|
||||
int _SolveSystemGS( PointInfo& pointInfo , int depth , const typename BSplineData< Degree >::Integrator& integrator , const SortedTreeNodes& sNodes , Pointer( Real ) solution , Pointer( Real ) constraints , Pointer( Real ) metSolutionConstraints , int iters , bool coarseToFine , bool showResidual=false , double* bNorm2=NULL , double* inRNorm2=NULL , double* outRNorm2=NULL , bool forceSilent=false );
|
||||
int _SolveSystemCG( PointInfo& pointInfo , int depth , const typename BSplineData< Degree >::Integrator& integrator , const SortedTreeNodes& sNodes , Pointer( Real ) solution , Pointer( Real ) constraints , Pointer( Real ) metSolutionConstraints , int iters , bool coarseToFine , bool showResidual=false , double* bNorm2=NULL , double* inRNorm2=NULL , double* outRNorm2=NULL , double accuracy=0 );
|
||||
|
||||
int GetMatrixRowSize( const typename TreeOctNode::Neighbors5& neighbors5 , bool symmetric ) const;
|
||||
int SetMatrixRow( const PointInfo& pointInfo , const typename TreeOctNode::Neighbors5& neighbors5 , Pointer( MatrixEntry< Real > ) row , int offset , const typename BSplineData< Degree >::Integrator& integrator , const Stencil< double , 5 >& stencil , bool symmetric ) const;
|
||||
|
||||
void SetDivergenceStencil ( int depth , const typename BSplineData< Degree >::Integrator& integrator , Stencil< Point3D< double > , 5 >& stencil , bool scatter ) const;
|
||||
void SetDivergenceStencils( int depth , const typename BSplineData< Degree >::Integrator& integrator , Stencil< Point3D< double > , 5 > stencil[2][2][2] , bool scatter ) const;
|
||||
void SetLaplacianStencil ( int depth , const typename BSplineData< Degree >::Integrator& integrator , Stencil< double , 5 >& stencil ) const;
|
||||
void SetLaplacianStencils ( int depth , const typename BSplineData< Degree >::Integrator& integrator , Stencil< double , 5 > stencil[2][2][2] ) const;
|
||||
void SetCenterEvaluationStencil ( const typename BSplineData< Degree >::template CenterEvaluator< 1 >& evaluator , int depth , Stencil< double , 3 >& stencil ) const;
|
||||
void SetCenterEvaluationStencils( const typename BSplineData< Degree >::template CenterEvaluator< 1 >& evaluator , int depth , Stencil< double , 3 > stencil[8] ) const;
|
||||
void SetCornerEvaluationStencil ( const typename BSplineData< Degree >::template CornerEvaluator< 2 >& evaluator , int depth , Stencil< double , 3 > stencil [8] ) const;
|
||||
void SetCornerEvaluationStencils( const typename BSplineData< Degree >::template CornerEvaluator< 2 >& evaluator , int depth , Stencil< double , 3 > stencils[8][8] ) const;
|
||||
void SetCornerNormalEvaluationStencil ( const typename BSplineData< Degree >::template CornerEvaluator< 2 >& evaluator , int depth , Stencil< Point3D< double > , 3 > stencil [8] ) const;
|
||||
void SetCornerNormalEvaluationStencils( const typename BSplineData< Degree >::template CornerEvaluator< 2 >& evaluator , int depth , Stencil< Point3D< double > , 3 > stencils[8][8] ) const;
|
||||
void SetCornerNormalEvaluationStencil ( const typename BSplineData< Degree >::template CornerEvaluator< 2 >& evaluator , int depth , Stencil< Point3D< double > , 5 > stencil [8] ) const;
|
||||
void SetCornerNormalEvaluationStencils( const typename BSplineData< Degree >::template CornerEvaluator< 2 >& evaluator , int depth , Stencil< Point3D< double > , 5 > stencils[8][8] ) const;
|
||||
|
||||
static void UpdateCoarserSupportBounds( const TreeOctNode* node , int& startX , int& endX , int& startY , int& endY , int& startZ , int& endZ );
|
||||
|
||||
void UpdateConstraintsFromCoarser( const PointInfo& pointInfo , const typename TreeOctNode::Neighbors5& neighbors5 , const typename TreeOctNode::Neighbors5& pNeighbors5 , TreeOctNode* node , Pointer( Real ) constraints , ConstPointer( Real ) metSolution , const typename BSplineData< Degree >::Integrator& integrator , const Stencil< double , 5 >& stencil ) const;
|
||||
// Updates the constraints @(depth-1) based on the solution coefficients @(depth)
|
||||
void UpdateConstraintsFromFiner( const typename BSplineData< Degree >::Integrator& integrator , int depth , const SortedTreeNodes& sNodes , ConstPointer( Real ) fineSolution , Pointer( Real ) coarseConstraints ) const;
|
||||
// Evaluate the points @(depth) using coefficients @(depth-1)
|
||||
void SetPointValuesFromCoarser( PointInfo& pointInfo , int depth , const SortedTreeNodes& sNodes , ConstPointer( Real ) coarseCoefficients );
|
||||
// Evalutes the solution @(depth) at the points @(depth-1) and updates the met constraints @(depth-1)
|
||||
void SetPointConstraintsFromFiner( const PointInfo& pointInfo , int depth , const SortedTreeNodes& sNodes , ConstPointer( Real ) finerCoefficients , Pointer( Real ) metConstraints) const;
|
||||
Real _WeightedCoarserFunctionValue( const _PointData& pointData , const typename TreeOctNode::NeighborKey3& neighborKey3 , const TreeOctNode* node , ConstPointer( Real ) coarseCoefficients ) const;
|
||||
Real _WeightedFinerFunctionValue ( const _PointData& pointData , const typename TreeOctNode::NeighborKey3& neighborKey3 , const TreeOctNode* node , ConstPointer( Real ) finerCoefficients ) const;
|
||||
|
||||
// Down samples constraints @(depth) to constraints @(depth-1)
|
||||
template< class C > void DownSample( int depth , const SortedTreeNodes& sNodes , ConstPointer( C ) fineConstraints , Pointer( C ) coarseConstraints ) const;
|
||||
// Up samples solution @(depth-1) to solution @(depth)
|
||||
template< class C > void UpSample ( int depth , const SortedTreeNodes& sNodes , ConstPointer( C ) coarseCoefficients , Pointer( C ) fineCoefficients ) const;
|
||||
int GetSliceMatrixAndUpdateConstraints( const PointInfo& pointInfo , SparseMatrix< Real >& matrix , Pointer( Real ) constraints , const typename BSplineData< Degree >::Integrator& integrator , int depth , const SortedTreeNodes& sNodes , ConstPointer( Real ) metSolution , bool coarseToFine , int nStart , int nEnd );
|
||||
int GetMatrixAndUpdateConstraints( const PointInfo& pointInfo , SparseSymmetricMatrix< Real >& matrix , Pointer( Real ) constraints , const typename BSplineData< Degree >::Integrator& integrator , int depth , const SortedTreeNodes& sNodes , ConstPointer( Real ) metSolution , bool coarseToFine );
|
||||
|
||||
|
||||
int UpdateWeightContribution( std::vector< Real >& kernelDensityWeights , TreeOctNode* node , const Point3D<Real>& position , typename TreeOctNode::NeighborKey3& neighborKey , Real weight=Real(1.0) );
|
||||
Real GetSampleWeight( ConstPointer( Real ) kernelDensityWeight , const Point3D<Real>& position , typename TreeOctNode::NeighborKey3& neighborKey , int splatDepth );
|
||||
Real GetSampleWeight( ConstPointer( Real ) kernelDensityWeight , const TreeOctNode* node , const Point3D<Real>& position , typename TreeOctNode::ConstNeighborKey3& neighborKey );
|
||||
void GetSampleDepthAndWeight( ConstPointer( Real ) kernelDensityWeight , const TreeOctNode* node , const Point3D<Real>& position , typename TreeOctNode::ConstNeighborKey3& neighborKey , Real samplesPerNode , Real& depth , Real& weight );
|
||||
Real GetSampleWeight( ConstPointer( Real ) kernelDensityWeight , TreeOctNode* node , const Point3D<Real>& position , typename TreeOctNode::NeighborKey3& neighborKey );
|
||||
void GetSampleDepthAndWeight( ConstPointer( Real ) kernelDensityWeight , TreeOctNode* node , const Point3D<Real>& position , typename TreeOctNode::NeighborKey3& neighborKey , Real samplesPerNode , Real& depth , Real& weight );
|
||||
int SplatOrientedPoint( ConstPointer( Real ) kernelDensityWeights , TreeOctNode* node , const Point3D<Real>& point , const Point3D< Real >& normal , NormalInfo& normalInfo , typename TreeOctNode::NeighborKey3& neighborKey );
|
||||
Real SplatOrientedPoint( ConstPointer( Real ) kernelDensityWeights , const Point3D<Real>& point , const Point3D<Real>& normal , NormalInfo& normalInfo , typename TreeOctNode::NeighborKey3& neighborKey , int kernelDepth , Real samplesPerNode , int minDepth , int maxDepth );
|
||||
|
||||
int HasNormals( TreeOctNode* node , const NormalInfo& normalInfo );
|
||||
|
||||
///////////////////////////
|
||||
// Iso-Surfacing Methods //
|
||||
///////////////////////////
|
||||
struct IsoEdge
|
||||
{
|
||||
long long edges[2];
|
||||
IsoEdge( void ){ edges[0] = edges[1] = 0; }
|
||||
IsoEdge( long long v1 , long long v2 ){ edges[0] = v1 , edges[1] = v2; }
|
||||
long long& operator[]( int idx ){ return edges[idx]; }
|
||||
const long long& operator[]( int idx ) const { return edges[idx]; }
|
||||
};
|
||||
struct FaceEdges
|
||||
{
|
||||
IsoEdge edges[2];
|
||||
int count;
|
||||
};
|
||||
template< class Vertex >
|
||||
struct SliceValues
|
||||
{
|
||||
typename SortedTreeNodes::SliceTableData sliceData;
|
||||
Pointer( Real ) cornerValues ; Pointer( Point3D< Real > ) cornerNormals ; Pointer( char ) cornerSet;
|
||||
Pointer( long long ) edgeKeys ; Pointer( char ) edgeSet;
|
||||
Pointer( FaceEdges ) faceEdges ; Pointer( char ) faceSet;
|
||||
Pointer( char ) mcIndices;
|
||||
hash_map< long long , std::vector< IsoEdge > > faceEdgeMap;
|
||||
hash_map< long long , std::pair< int , Vertex > > edgeVertexMap;
|
||||
hash_map< long long , long long > vertexPairMap;
|
||||
|
||||
SliceValues( void );
|
||||
~SliceValues( void );
|
||||
void reset( bool nonLinearFit );
|
||||
protected:
|
||||
int _oldCCount , _oldECount , _oldFCount , _oldNCount;
|
||||
};
|
||||
template< class Vertex >
|
||||
struct XSliceValues
|
||||
{
|
||||
typename SortedTreeNodes::XSliceTableData xSliceData;
|
||||
Pointer( long long ) edgeKeys ; Pointer( char ) edgeSet;
|
||||
Pointer( FaceEdges ) faceEdges ; Pointer( char ) faceSet;
|
||||
hash_map< long long , std::vector< IsoEdge > > faceEdgeMap;
|
||||
hash_map< long long , std::pair< int , Vertex > > edgeVertexMap;
|
||||
hash_map< long long , long long > vertexPairMap;
|
||||
|
||||
XSliceValues( void );
|
||||
~XSliceValues( void );
|
||||
void reset( void );
|
||||
protected:
|
||||
int _oldECount , _oldFCount;
|
||||
};
|
||||
template< class Vertex >
|
||||
struct SlabValues
|
||||
{
|
||||
XSliceValues< Vertex > _xSliceValues[2];
|
||||
SliceValues< Vertex > _sliceValues[2];
|
||||
SliceValues< Vertex >& sliceValues( int idx ){ return _sliceValues[idx&1]; }
|
||||
const SliceValues< Vertex >& sliceValues( int idx ) const { return _sliceValues[idx&1]; }
|
||||
XSliceValues< Vertex >& xSliceValues( int idx ){ return _xSliceValues[idx&1]; }
|
||||
const XSliceValues< Vertex >& xSliceValues( int idx ) const { return _xSliceValues[idx&1]; }
|
||||
};
|
||||
template< class Vertex >
|
||||
void SetSliceIsoCorners( ConstPointer( Real ) solution , ConstPointer( Real ) coarseSolution , Real isoValue , int depth , int slice , std::vector< SlabValues< Vertex > >& sValues , const typename BSplineData< Degree >::template CornerEvaluator< 2 >& evaluator , const Stencil< double , 3 > stencil[8] , const Stencil< double , 3 > stencils[8][8] , const Stencil< Point3D< double > , 3 > nStencil[8] , const Stencil< Point3D< double > , 3 > nStencils[8][8] , int threads );
|
||||
template< class Vertex >
|
||||
void SetSliceIsoCorners( ConstPointer( Real ) solution , ConstPointer( Real ) coarseSolution , Real isoValue , int depth , int slice , int z , std::vector< SlabValues< Vertex > >& sValues , const typename BSplineData< Degree >::template CornerEvaluator< 2 >& evaluator , const Stencil< double , 3 > stencil[8] , const Stencil< double , 3 > stencils[8][8] , const Stencil< Point3D< double > , 3 > nStencil[8] , const Stencil< Point3D< double > , 3 > nStencils[8][8] , int threads );
|
||||
template< class Vertex >
|
||||
void SetSliceIsoVertices( ConstPointer( Real ) kernelDensityWeights , Real isoValue , int depth , int slice , int& vOffset , CoredMeshData< Vertex >& mesh , std::vector< SlabValues< Vertex > >& sValues , int threads );
|
||||
template< class Vertex >
|
||||
void SetSliceIsoVertices( ConstPointer( Real ) kernelDensityWeights , Real isoValue , int depth , int slice , int z , int& vOffset , CoredMeshData< Vertex >& mesh , std::vector< SlabValues< Vertex > >& sValues , int threads );
|
||||
template< class Vertex >
|
||||
void SetXSliceIsoVertices( ConstPointer( Real ) kernelDensityWeights , Real isoValue , int depth , int slab , int& vOffset , CoredMeshData< Vertex >& mesh , std::vector< SlabValues< Vertex > >& sValues , int threads );
|
||||
template< class Vertex >
|
||||
void CopyFinerSliceIsoEdgeKeys( int depth , int slice , std::vector< SlabValues< Vertex > >& sValues , int threads );
|
||||
template< class Vertex >
|
||||
void CopyFinerSliceIsoEdgeKeys( int depth , int slice , int z , std::vector< SlabValues< Vertex > >& sValues , int threads );
|
||||
template< class Vertex >
|
||||
void CopyFinerXSliceIsoEdgeKeys( int depth , int slab , std::vector< SlabValues< Vertex > >& sValues , int threads );
|
||||
template< class Vertex >
|
||||
void SetSliceIsoEdges( int depth , int slice , std::vector< SlabValues< Vertex > >& slabValues , int threads );
|
||||
template< class Vertex >
|
||||
void SetSliceIsoEdges( int depth , int slice , int z , std::vector< SlabValues< Vertex > >& slabValues , int threads );
|
||||
template< class Vertex >
|
||||
void SetXSliceIsoEdges( int depth , int slice , std::vector< SlabValues< Vertex > >& slabValues , int threads );
|
||||
|
||||
template< class Vertex >
|
||||
int SetIsoSurface( int depth , int offset , const SliceValues< Vertex >& bValues , const SliceValues< Vertex >& fValues , const XSliceValues< Vertex >& xValues , CoredMeshData< Vertex >& mesh , bool polygonMesh , bool addBarycenter , int& vOffset , int threads );
|
||||
|
||||
template< class Vertex >
|
||||
static int AddIsoPolygons( CoredMeshData< Vertex >& mesh , std::vector< std::pair< int , Vertex > >& polygon , bool polygonMesh , bool addBarycenter , int& vOffset );
|
||||
|
||||
template< class Vertex >
|
||||
bool GetIsoVertex( ConstPointer( Real ) kernelDensityWeights , Real isoValue , typename TreeOctNode::ConstNeighborKey3& neighborKey3 , const TreeOctNode* node , int edgeIndex , int z , const SliceValues< Vertex >& sValues , Vertex& vertex );
|
||||
template< class Vertex >
|
||||
bool GetIsoVertex( ConstPointer( Real ) kernelDensityWeights , Real isoValue , typename TreeOctNode::ConstNeighborKey3& neighborKey3 , const TreeOctNode* node , int cornerIndex , const SliceValues< Vertex >& bValues , const SliceValues< Vertex >& fValues , Vertex& vertex );
|
||||
|
||||
|
||||
////////////////////////
|
||||
// Evaluation Methods //
|
||||
////////////////////////
|
||||
Real getCornerValue( const typename TreeOctNode::ConstNeighborKey3& neighborKey3 , const TreeOctNode* node , int corner , ConstPointer( Real ) solution , ConstPointer( Real ) metSolution , const typename BSplineData< Degree >::template CornerEvaluator< 2 >& evaluator , const Stencil< double , 3 >& stencil , const Stencil< double , 3 > stencils[8] , bool isInterior ) const;
|
||||
Point3D< Real > getCornerNormal( const typename TreeOctNode::ConstNeighbors5& neighbors5 , const typename TreeOctNode::ConstNeighbors5& pNeighbors5 , const TreeOctNode* node , int corner , ConstPointer( Real ) solution , ConstPointer( Real ) metSolution , const typename BSplineData< Degree >::template CornerEvaluator< 2 >& evaluator , const Stencil< Point3D< double > , 5 >& nStencil , const Stencil< Point3D< double > , 5 > nStencils[8] , bool isInterior ) const;
|
||||
std::pair< Real , Point3D< Real > > getCornerValueAndNormal( const typename TreeOctNode::ConstNeighborKey3& neighborKey3 , const TreeOctNode* node , int corner , ConstPointer( Real ) solution , ConstPointer( Real ) metSolution , const typename BSplineData< Degree >::template CornerEvaluator< 2 >& evaluator , const Stencil< double , 3 >& vStencil , const Stencil< double , 3 > vStencils[8] , const Stencil< Point3D< double > , 3 >& nStencil , const Stencil< Point3D< double > , 3 > nStencils[8] , bool isInterior ) const;
|
||||
Real getCenterValue( const typename TreeOctNode::ConstNeighborKey3& neighborKey3 , const TreeOctNode* node , ConstPointer( Real ) solution , ConstPointer( Real ) metSolution , const typename BSplineData< Degree >::template CenterEvaluator< 1 >& evaluator , const Stencil< double , 3 >& stencil , const Stencil< double , 3 >& pStencil , bool isInterior ) const;
|
||||
|
||||
static bool _IsInset( const TreeOctNode* node );
|
||||
static bool _IsInsetSupported( const TreeOctNode* node );
|
||||
|
||||
void refineBoundary( std::vector< int >* map );
|
||||
public:
|
||||
int threads;
|
||||
static double maxMemoryUsage;
|
||||
TreeOctNode tree;
|
||||
|
||||
static double MemoryUsage( void );
|
||||
Octree( void );
|
||||
|
||||
void MakeComplete( std::vector< int >* map=NULL );
|
||||
void Finalize( std::vector< int >* map=NULL );
|
||||
void ClipTree( const NormalInfo& normalInfo );
|
||||
|
||||
Real Evaluate( ConstPointer( Real ) coefficients , Point3D< Real > p , const BSplineData< Degree >* fData=NULL ) const;
|
||||
Pointer( Real ) Evaluate( ConstPointer( Real ) coefficients , int& res , Real isoValue=0.f , int depth=-1 );
|
||||
template< class PointReal >
|
||||
int SetTree( char* fileName , int minDepth , int maxDepth , int fullDepth , int splatDepth , Real samplesPerNode ,
|
||||
Real scaleFactor , bool useConfidence , bool useNormalWeight , Real constraintWeight , int adaptiveExponent ,
|
||||
PointInfo& pointInfo , NormalInfo& normalInfo , std::vector< Real >& kernelDensityWeights , std::vector< Real >& centerWeights ,
|
||||
int boundaryType=BSplineElements< Degree >::NONE , XForm4x4< Real > xForm=XForm4x4< Real >::Identity , bool makeComplete=false );
|
||||
Pointer( Real ) SetLaplacianConstraints( const NormalInfo& normalInfo );
|
||||
Pointer( Real ) SolveSystem( PointInfo& pointInfo , Pointer( Real ) constraints , bool showResidual , int iters , int maxSolveDepth , int cgDepth=0 , double cgAccuracy=0 );
|
||||
|
||||
Real GetIsoValue( ConstPointer( Real ) solution , const std::vector< Real >& centerWeights );
|
||||
template< class Vertex >
|
||||
void GetMCIsoSurface( ConstPointer( Real ) kernelDensityWeights , ConstPointer( Real ) solution , Real isoValue , CoredMeshData< Vertex >& mesh , bool nonLinearFit=true , bool addBarycenter=false , bool polygonMesh=false );
|
||||
};
|
||||
|
||||
#include "MultiGridOctreeData.inl"
|
||||
#include "MultiGridOctreeData.SortedTreeNodes.inl"
|
||||
#include "MultiGridOctreeData.IsoSurface.inl"
|
||||
#endif // MULTI_GRID_OCTREE_DATA_INCLUDED
|
||||
113
src/plugins_experimental/filter_screened_poisson/Src/PPolynomial.h
Executable file
113
src/plugins_experimental/filter_screened_poisson/Src/PPolynomial.h
Executable file
@ -0,0 +1,113 @@
|
||||
/*
|
||||
Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
|
||||
All rights reserved.
|
||||
|
||||
Redistribution and use in source and binary forms, with or without modification,
|
||||
are permitted provided that the following conditions are met:
|
||||
|
||||
Redistributions of source code must retain the above copyright notice, this list of
|
||||
conditions and the following disclaimer. Redistributions in binary form must reproduce
|
||||
the above copyright notice, this list of conditions and the following disclaimer
|
||||
in the documentation and/or other materials provided with the distribution.
|
||||
|
||||
Neither the name of the Johns Hopkins University nor the names of its contributors
|
||||
may be used to endorse or promote products derived from this software without specific
|
||||
prior written permission.
|
||||
|
||||
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
|
||||
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
|
||||
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
|
||||
SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
|
||||
INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
|
||||
TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
|
||||
BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
|
||||
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
|
||||
ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
|
||||
DAMAGE.
|
||||
*/
|
||||
|
||||
#ifndef P_POLYNOMIAL_INCLUDED
|
||||
#define P_POLYNOMIAL_INCLUDED
|
||||
#include <vector>
|
||||
#include "Polynomial.h"
|
||||
|
||||
template<int Degree>
|
||||
class StartingPolynomial{
|
||||
public:
|
||||
Polynomial<Degree> p;
|
||||
double start;
|
||||
|
||||
template<int Degree2>
|
||||
StartingPolynomial<Degree+Degree2> operator * (const StartingPolynomial<Degree2>& p) const;
|
||||
StartingPolynomial scale(double s) const;
|
||||
StartingPolynomial shift(double t) const;
|
||||
int operator < (const StartingPolynomial& sp) const;
|
||||
static int Compare(const void* v1,const void* v2);
|
||||
};
|
||||
|
||||
template<int Degree>
|
||||
class PPolynomial
|
||||
{
|
||||
public:
|
||||
size_t polyCount;
|
||||
StartingPolynomial<Degree>* polys;
|
||||
|
||||
PPolynomial(void);
|
||||
PPolynomial(const PPolynomial<Degree>& p);
|
||||
~PPolynomial(void);
|
||||
|
||||
PPolynomial& operator = (const PPolynomial& p);
|
||||
|
||||
int size(void) const;
|
||||
|
||||
void set( size_t size );
|
||||
// Note: this method will sort the elements in sps
|
||||
void set( StartingPolynomial<Degree>* sps , int count );
|
||||
void reset( size_t newSize );
|
||||
|
||||
|
||||
double operator()( double t ) const;
|
||||
double integral( double tMin , double tMax ) const;
|
||||
double Integral( void ) const;
|
||||
|
||||
template<int Degree2>
|
||||
PPolynomial<Degree>& operator = (const PPolynomial<Degree2>& p);
|
||||
|
||||
PPolynomial operator + (const PPolynomial& p) const;
|
||||
PPolynomial operator - (const PPolynomial& p) const;
|
||||
|
||||
template<int Degree2>
|
||||
PPolynomial<Degree+Degree2> operator * (const Polynomial<Degree2>& p) const;
|
||||
|
||||
template<int Degree2>
|
||||
PPolynomial<Degree+Degree2> operator * (const PPolynomial<Degree2>& p) const;
|
||||
|
||||
|
||||
PPolynomial& operator += ( double s );
|
||||
PPolynomial& operator -= ( double s );
|
||||
PPolynomial& operator *= ( double s );
|
||||
PPolynomial& operator /= ( double s );
|
||||
PPolynomial operator + ( double s ) const;
|
||||
PPolynomial operator - ( double s ) const;
|
||||
PPolynomial operator * ( double s ) const;
|
||||
PPolynomial operator / ( double s ) const;
|
||||
|
||||
PPolynomial& addScaled(const PPolynomial& poly,double scale);
|
||||
|
||||
PPolynomial scale( double s ) const;
|
||||
PPolynomial shift( double t ) const;
|
||||
|
||||
PPolynomial< Degree-1 > derivative(void) const;
|
||||
PPolynomial< Degree+1 > integral(void) const;
|
||||
|
||||
void getSolutions(double c,std::vector<double>& roots,double EPS,double min=-DBL_MAX,double max=DBL_MAX) const;
|
||||
|
||||
void printnl( void ) const;
|
||||
|
||||
PPolynomial< Degree+1 > MovingAverage( double radius ) const;
|
||||
static PPolynomial BSpline( double radius=0.5 );
|
||||
|
||||
void write( FILE* fp , int samples , double min , double max ) const;
|
||||
};
|
||||
#include "PPolynomial.inl"
|
||||
#endif // P_POLYNOMIAL_INCLUDED
|
||||
431
src/plugins_experimental/filter_screened_poisson/Src/PPolynomial.inl
Executable file
431
src/plugins_experimental/filter_screened_poisson/Src/PPolynomial.inl
Executable file
@ -0,0 +1,431 @@
|
||||
/*
|
||||
Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
|
||||
All rights reserved.
|
||||
|
||||
Redistribution and use in source and binary forms, with or without modification,
|
||||
are permitted provided that the following conditions are met:
|
||||
|
||||
Redistributions of source code must retain the above copyright notice, this list of
|
||||
conditions and the following disclaimer. Redistributions in binary form must reproduce
|
||||
the above copyright notice, this list of conditions and the following disclaimer
|
||||
in the documentation and/or other materials provided with the distribution.
|
||||
|
||||
Neither the name of the Johns Hopkins University nor the names of its contributors
|
||||
may be used to endorse or promote products derived from this software without specific
|
||||
prior written permission.
|
||||
|
||||
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
|
||||
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
|
||||
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
|
||||
SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
|
||||
INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
|
||||
TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
|
||||
BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
|
||||
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
|
||||
ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
|
||||
DAMAGE.
|
||||
*/
|
||||
|
||||
#include "Factor.h"
|
||||
|
||||
////////////////////////
|
||||
// StartingPolynomial //
|
||||
////////////////////////
|
||||
template<int Degree>
|
||||
template<int Degree2>
|
||||
StartingPolynomial<Degree+Degree2> StartingPolynomial<Degree>::operator * (const StartingPolynomial<Degree2>& p) const{
|
||||
StartingPolynomial<Degree+Degree2> sp;
|
||||
if(start>p.start){sp.start=start;}
|
||||
else{sp.start=p.start;}
|
||||
sp.p=this->p*p.p;
|
||||
return sp;
|
||||
}
|
||||
template<int Degree>
|
||||
StartingPolynomial<Degree> StartingPolynomial<Degree>::scale(double s) const{
|
||||
StartingPolynomial q;
|
||||
q.start=start*s;
|
||||
q.p=p.scale(s);
|
||||
return q;
|
||||
}
|
||||
template<int Degree>
|
||||
StartingPolynomial<Degree> StartingPolynomial<Degree>::shift(double s) const{
|
||||
StartingPolynomial q;
|
||||
q.start=start+s;
|
||||
q.p=p.shift(s);
|
||||
return q;
|
||||
}
|
||||
|
||||
|
||||
template<int Degree>
|
||||
int StartingPolynomial<Degree>::operator < (const StartingPolynomial<Degree>& sp) const{
|
||||
if(start<sp.start){return 1;}
|
||||
else{return 0;}
|
||||
}
|
||||
template<int Degree>
|
||||
int StartingPolynomial<Degree>::Compare(const void* v1,const void* v2){
|
||||
double d=((StartingPolynomial*)(v1))->start-((StartingPolynomial*)(v2))->start;
|
||||
if (d<0) {return -1;}
|
||||
else if (d>0) {return 1;}
|
||||
else {return 0;}
|
||||
}
|
||||
|
||||
/////////////////
|
||||
// PPolynomial //
|
||||
/////////////////
|
||||
template<int Degree>
|
||||
PPolynomial<Degree>::PPolynomial(void){
|
||||
polyCount=0;
|
||||
polys=NULL;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree>::PPolynomial(const PPolynomial<Degree>& p){
|
||||
polyCount=0;
|
||||
polys=NULL;
|
||||
set(p.polyCount);
|
||||
memcpy(polys,p.polys,sizeof(StartingPolynomial<Degree>)*p.polyCount);
|
||||
}
|
||||
|
||||
template<int Degree>
|
||||
PPolynomial<Degree>::~PPolynomial(void){
|
||||
if(polyCount){free(polys);}
|
||||
polyCount=0;
|
||||
polys=NULL;
|
||||
}
|
||||
template<int Degree>
|
||||
void PPolynomial<Degree>::set(StartingPolynomial<Degree>* sps,int count){
|
||||
int i,c=0;
|
||||
set(count);
|
||||
qsort(sps,count,sizeof(StartingPolynomial<Degree>),StartingPolynomial<Degree>::Compare);
|
||||
for( i=0 ; i<count ; i++ )
|
||||
{
|
||||
if( !c || sps[i].start!=polys[c-1].start ) polys[c++] = sps[i];
|
||||
else{polys[c-1].p+=sps[i].p;}
|
||||
}
|
||||
reset( c );
|
||||
}
|
||||
template <int Degree>
|
||||
int PPolynomial<Degree>::size(void) const{return int(sizeof(StartingPolynomial<Degree>)*polyCount);}
|
||||
|
||||
template<int Degree>
|
||||
void PPolynomial<Degree>::set( size_t size )
|
||||
{
|
||||
if(polyCount){free(polys);}
|
||||
polyCount=0;
|
||||
polys=NULL;
|
||||
polyCount=size;
|
||||
if(size){
|
||||
polys=(StartingPolynomial<Degree>*)malloc(sizeof(StartingPolynomial<Degree>)*size);
|
||||
memset(polys,0,sizeof(StartingPolynomial<Degree>)*size);
|
||||
}
|
||||
}
|
||||
template<int Degree>
|
||||
void PPolynomial<Degree>::reset( size_t newSize )
|
||||
{
|
||||
polyCount=newSize;
|
||||
polys=(StartingPolynomial<Degree>*)realloc(polys,sizeof(StartingPolynomial<Degree>)*newSize);
|
||||
}
|
||||
|
||||
template<int Degree>
|
||||
PPolynomial<Degree>& PPolynomial<Degree>::operator = (const PPolynomial<Degree>& p){
|
||||
set(p.polyCount);
|
||||
memcpy(polys,p.polys,sizeof(StartingPolynomial<Degree>)*p.polyCount);
|
||||
return *this;
|
||||
}
|
||||
|
||||
template<int Degree>
|
||||
template<int Degree2>
|
||||
PPolynomial<Degree>& PPolynomial<Degree>::operator = (const PPolynomial<Degree2>& p){
|
||||
set(p.polyCount);
|
||||
for(int i=0;i<int(polyCount);i++){
|
||||
polys[i].start=p.polys[i].start;
|
||||
polys[i].p=p.polys[i].p;
|
||||
}
|
||||
return *this;
|
||||
}
|
||||
|
||||
template<int Degree>
|
||||
double PPolynomial<Degree>::operator ()( double t ) const
|
||||
{
|
||||
double v=0;
|
||||
for( int i=0 ; i<int(polyCount) && t>polys[i].start ; i++ ) v+=polys[i].p(t);
|
||||
return v;
|
||||
}
|
||||
|
||||
template<int Degree>
|
||||
double PPolynomial<Degree>::integral( double tMin , double tMax ) const
|
||||
{
|
||||
int m=1;
|
||||
double start,end,s,v=0;
|
||||
start=tMin;
|
||||
end=tMax;
|
||||
if(tMin>tMax){
|
||||
m=-1;
|
||||
start=tMax;
|
||||
end=tMin;
|
||||
}
|
||||
for(int i=0;i<int(polyCount) && polys[i].start<end;i++){
|
||||
if(start<polys[i].start){s=polys[i].start;}
|
||||
else{s=start;}
|
||||
v+=polys[i].p.integral(s,end);
|
||||
}
|
||||
return v*m;
|
||||
}
|
||||
template<int Degree>
|
||||
double PPolynomial<Degree>::Integral(void) const{return integral(polys[0].start,polys[polyCount-1].start);}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree> PPolynomial<Degree>::operator + (const PPolynomial<Degree>& p) const{
|
||||
PPolynomial q;
|
||||
int i,j;
|
||||
size_t idx=0;
|
||||
q.set(polyCount+p.polyCount);
|
||||
i=j=-1;
|
||||
|
||||
while(idx<q.polyCount){
|
||||
if (j>=int(p.polyCount)-1) {q.polys[idx]= polys[++i];}
|
||||
else if (i>=int( polyCount)-1) {q.polys[idx]=p.polys[++j];}
|
||||
else if(polys[i+1].start<p.polys[j+1].start){q.polys[idx]= polys[++i];}
|
||||
else {q.polys[idx]=p.polys[++j];}
|
||||
idx++;
|
||||
}
|
||||
return q;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree> PPolynomial<Degree>::operator - (const PPolynomial<Degree>& p) const{
|
||||
PPolynomial q;
|
||||
int i,j;
|
||||
size_t idx=0;
|
||||
q.set(polyCount+p.polyCount);
|
||||
i=j=-1;
|
||||
|
||||
while(idx<q.polyCount){
|
||||
if (j>=int(p.polyCount)-1) {q.polys[idx]= polys[++i];}
|
||||
else if (i>=int( polyCount)-1) {q.polys[idx].start=p.polys[++j].start;q.polys[idx].p=p.polys[j].p*(-1.0);}
|
||||
else if(polys[i+1].start<p.polys[j+1].start){q.polys[idx]= polys[++i];}
|
||||
else {q.polys[idx].start=p.polys[++j].start;q.polys[idx].p=p.polys[j].p*(-1.0);}
|
||||
idx++;
|
||||
}
|
||||
return q;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree>& PPolynomial<Degree>::addScaled(const PPolynomial<Degree>& p,double scale){
|
||||
int i,j;
|
||||
StartingPolynomial<Degree>* oldPolys=polys;
|
||||
size_t idx=0,cnt=0,oldPolyCount=polyCount;
|
||||
polyCount=0;
|
||||
polys=NULL;
|
||||
set(oldPolyCount+p.polyCount);
|
||||
i=j=-1;
|
||||
while(cnt<polyCount){
|
||||
if (j>=int( p.polyCount)-1) {polys[idx]=oldPolys[++i];}
|
||||
else if (i>=int(oldPolyCount)-1) {polys[idx].start= p.polys[++j].start;polys[idx].p=p.polys[j].p*scale;}
|
||||
else if (oldPolys[i+1].start<p.polys[j+1].start){polys[idx]=oldPolys[++i];}
|
||||
else {polys[idx].start= p.polys[++j].start;polys[idx].p=p.polys[j].p*scale;}
|
||||
if(idx && polys[idx].start==polys[idx-1].start) {polys[idx-1].p+=polys[idx].p;}
|
||||
else{idx++;}
|
||||
cnt++;
|
||||
}
|
||||
free(oldPolys);
|
||||
reset(idx);
|
||||
return *this;
|
||||
}
|
||||
template<int Degree>
|
||||
template<int Degree2>
|
||||
PPolynomial<Degree+Degree2> PPolynomial<Degree>::operator * (const PPolynomial<Degree2>& p) const{
|
||||
PPolynomial<Degree+Degree2> q;
|
||||
StartingPolynomial<Degree+Degree2> *sp;
|
||||
int i,j,spCount=int(polyCount*p.polyCount);
|
||||
|
||||
sp=(StartingPolynomial<Degree+Degree2>*)malloc(sizeof(StartingPolynomial<Degree+Degree2>)*spCount);
|
||||
for(i=0;i<int(polyCount);i++){
|
||||
for(j=0;j<int(p.polyCount);j++){
|
||||
sp[i*p.polyCount+j]=polys[i]*p.polys[j];
|
||||
}
|
||||
}
|
||||
q.set(sp,spCount);
|
||||
free(sp);
|
||||
return q;
|
||||
}
|
||||
template<int Degree>
|
||||
template<int Degree2>
|
||||
PPolynomial<Degree+Degree2> PPolynomial<Degree>::operator * (const Polynomial<Degree2>& p) const{
|
||||
PPolynomial<Degree+Degree2> q;
|
||||
q.set(polyCount);
|
||||
for(int i=0;i<int(polyCount);i++){
|
||||
q.polys[i].start=polys[i].start;
|
||||
q.polys[i].p=polys[i].p*p;
|
||||
}
|
||||
return q;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree> PPolynomial<Degree>::scale( double s ) const
|
||||
{
|
||||
PPolynomial q;
|
||||
q.set(polyCount);
|
||||
for(size_t i=0;i<polyCount;i++){q.polys[i]=polys[i].scale(s);}
|
||||
return q;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree> PPolynomial<Degree>::shift( double s ) const
|
||||
{
|
||||
PPolynomial q;
|
||||
q.set(polyCount);
|
||||
for(size_t i=0;i<polyCount;i++){q.polys[i]=polys[i].shift(s);}
|
||||
return q;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree-1> PPolynomial<Degree>::derivative(void) const{
|
||||
PPolynomial<Degree-1> q;
|
||||
q.set(polyCount);
|
||||
for(size_t i=0;i<polyCount;i++){
|
||||
q.polys[i].start=polys[i].start;
|
||||
q.polys[i].p=polys[i].p.derivative();
|
||||
}
|
||||
return q;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree+1> PPolynomial<Degree>::integral(void) const{
|
||||
int i;
|
||||
PPolynomial<Degree+1> q;
|
||||
q.set(polyCount);
|
||||
for(i=0;i<int(polyCount);i++){
|
||||
q.polys[i].start=polys[i].start;
|
||||
q.polys[i].p=polys[i].p.integral();
|
||||
q.polys[i].p-=q.polys[i].p(q.polys[i].start);
|
||||
}
|
||||
return q;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree>& PPolynomial<Degree>::operator += ( double s ) {polys[0].p+=s;}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree>& PPolynomial<Degree>::operator -= ( double s ) {polys[0].p-=s;}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree>& PPolynomial<Degree>::operator *= ( double s )
|
||||
{
|
||||
for(int i=0;i<int(polyCount);i++){polys[i].p*=s;}
|
||||
return *this;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree>& PPolynomial<Degree>::operator /= ( double s )
|
||||
{
|
||||
for(size_t i=0;i<polyCount;i++){polys[i].p/=s;}
|
||||
return *this;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree> PPolynomial<Degree>::operator + ( double s ) const
|
||||
{
|
||||
PPolynomial q=*this;
|
||||
q+=s;
|
||||
return q;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree> PPolynomial<Degree>::operator - ( double s ) const
|
||||
{
|
||||
PPolynomial q=*this;
|
||||
q-=s;
|
||||
return q;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree> PPolynomial<Degree>::operator * ( double s ) const
|
||||
{
|
||||
PPolynomial q=*this;
|
||||
q*=s;
|
||||
return q;
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree> PPolynomial<Degree>::operator / ( double s ) const
|
||||
{
|
||||
PPolynomial q=*this;
|
||||
q/=s;
|
||||
return q;
|
||||
}
|
||||
|
||||
template<int Degree>
|
||||
void PPolynomial<Degree>::printnl(void) const{
|
||||
Polynomial<Degree> p;
|
||||
|
||||
if(!polyCount){
|
||||
Polynomial<Degree> p;
|
||||
printf("[-Infinity,Infinity]\n");
|
||||
}
|
||||
else{
|
||||
for(size_t i=0;i<polyCount;i++){
|
||||
printf("[");
|
||||
if (polys[i ].start== DBL_MAX){printf("Infinity,");}
|
||||
else if (polys[i ].start==-DBL_MAX){printf("-Infinity,");}
|
||||
else {printf("%f,",polys[i].start);}
|
||||
if(i+1==polyCount) {printf("Infinity]\t");}
|
||||
else if (polys[i+1].start== DBL_MAX){printf("Infinity]\t");}
|
||||
else if (polys[i+1].start==-DBL_MAX){printf("-Infinity]\t");}
|
||||
else {printf("%f]\t",polys[i+1].start);}
|
||||
p=p+polys[i].p;
|
||||
p.printnl();
|
||||
}
|
||||
}
|
||||
printf("\n");
|
||||
}
|
||||
template< >
|
||||
PPolynomial< 0 > PPolynomial< 0 >::BSpline( double radius )
|
||||
{
|
||||
PPolynomial q;
|
||||
q.set(2);
|
||||
|
||||
q.polys[0].start=-radius;
|
||||
q.polys[1].start= radius;
|
||||
|
||||
q.polys[0].p.coefficients[0]= 1.0;
|
||||
q.polys[1].p.coefficients[0]=-1.0;
|
||||
return q;
|
||||
}
|
||||
template< int Degree >
|
||||
PPolynomial< Degree > PPolynomial<Degree>::BSpline( double radius )
|
||||
{
|
||||
return PPolynomial< Degree-1 >::BSpline().MovingAverage( radius );
|
||||
}
|
||||
template<int Degree>
|
||||
PPolynomial<Degree+1> PPolynomial<Degree>::MovingAverage( double radius ) const
|
||||
{
|
||||
PPolynomial<Degree+1> A;
|
||||
Polynomial<Degree+1> p;
|
||||
StartingPolynomial<Degree+1>* sps;
|
||||
|
||||
sps=(StartingPolynomial<Degree+1>*)malloc(sizeof(StartingPolynomial<Degree+1>)*polyCount*2);
|
||||
|
||||
for(int i=0;i<int(polyCount);i++){
|
||||
sps[2*i ].start=polys[i].start-radius;
|
||||
sps[2*i+1].start=polys[i].start+radius;
|
||||
p=polys[i].p.integral()-polys[i].p.integral()(polys[i].start);
|
||||
sps[2*i ].p=p.shift(-radius);
|
||||
sps[2*i+1].p=p.shift( radius)*-1;
|
||||
}
|
||||
A.set(sps,int(polyCount*2));
|
||||
free(sps);
|
||||
return A*1.0/(2*radius);
|
||||
}
|
||||
template<int Degree>
|
||||
void PPolynomial<Degree>::getSolutions(double c,std::vector<double>& roots,double EPS,double min,double max) const{
|
||||
Polynomial<Degree> p;
|
||||
std::vector<double> tempRoots;
|
||||
|
||||
p.setZero();
|
||||
for(size_t i=0;i<polyCount;i++){
|
||||
p+=polys[i].p;
|
||||
if(polys[i].start>max){break;}
|
||||
if(i<polyCount-1 && polys[i+1].start<min){continue;}
|
||||
p.getSolutions(c,tempRoots,EPS);
|
||||
for(size_t j=0;j<tempRoots.size();j++){
|
||||
if(tempRoots[j]>polys[i].start && (i+1==polyCount || tempRoots[j]<=polys[i+1].start)){
|
||||
if(tempRoots[j]>min && tempRoots[j]<max){roots.push_back(tempRoots[j]);}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
template<int Degree>
|
||||
void PPolynomial<Degree>::write(FILE* fp,int samples,double min,double max) const{
|
||||
fwrite(&samples,sizeof(int),1,fp);
|
||||
for(int i=0;i<samples;i++){
|
||||
double x=min+i*(max-min)/(samples-1);
|
||||
float v=(*this)(x);
|
||||
fwrite(&v,sizeof(float),1,fp);
|
||||
}
|
||||
}
|
||||
705
src/plugins_experimental/filter_screened_poisson/Src/Ply.h
Executable file
705
src/plugins_experimental/filter_screened_poisson/Src/Ply.h
Executable file
@ -0,0 +1,705 @@
|
||||
/*
|
||||
|
||||
Header for PLY polygon files.
|
||||
|
||||
- Greg Turk, March 1994
|
||||
|
||||
A PLY file contains a single polygonal _object_.
|
||||
|
||||
An object is composed of lists of _elements_. Typical elements are
|
||||
vertices, faces, edges and materials.
|
||||
|
||||
Each type of element for a given object has one or more _properties_
|
||||
associated with the element type. For instance, a vertex element may
|
||||
have as properties three floating-point values x,y,z and three unsigned
|
||||
chars for red, green and blue.
|
||||
|
||||
---------------------------------------------------------------
|
||||
|
||||
Copyright (c) 1994 The Board of Trustees of The Leland Stanford
|
||||
Junior University. All rights reserved.
|
||||
|
||||
Permission to use, copy, modify and distribute this software and its
|
||||
documentation for any purpose is hereby granted without fee, provided
|
||||
that the above copyright notice and this permission notice appear in
|
||||
all copies of this software and that you do not sell the software.
|
||||
|
||||
THE SOFTWARE IS PROVIDED "AS IS" AND WITHOUT WARRANTY OF ANY KIND,
|
||||
EXPRESS, IMPLIED OR OTHERWISE, INCLUDING WITHOUT LIMITATION, ANY
|
||||
WARRANTY OF MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE.
|
||||
|
||||
*/
|
||||
|
||||
#ifndef __PLY_H__
|
||||
#define __PLY_H__
|
||||
|
||||
#ifndef WIN32
|
||||
#define _strdup strdup
|
||||
#endif
|
||||
|
||||
#ifdef __cplusplus
|
||||
extern "C" {
|
||||
#endif
|
||||
|
||||
#include <stdlib.h>
|
||||
#include <stdio.h>
|
||||
#include <stddef.h>
|
||||
#include <string.h>
|
||||
|
||||
#define PLY_ASCII 1 /* ascii PLY file */
|
||||
#define PLY_BINARY_BE 2 /* binary PLY file, big endian */
|
||||
#define PLY_BINARY_LE 3 /* binary PLY file, little endian */
|
||||
#define PLY_BINARY_NATIVE 4 /* binary PLY file, same endianness as current architecture */
|
||||
|
||||
#define PLY_OKAY 0 /* ply routine worked okay */
|
||||
#define PLY_ERROR -1 /* error in ply routine */
|
||||
|
||||
/* scalar data types supported by PLY format */
|
||||
|
||||
#define PLY_START_TYPE 0
|
||||
#define PLY_CHAR 1
|
||||
#define PLY_SHORT 2
|
||||
#define PLY_INT 3
|
||||
#define PLY_UCHAR 4
|
||||
#define PLY_USHORT 5
|
||||
#define PLY_UINT 6
|
||||
#define PLY_FLOAT 7
|
||||
#define PLY_DOUBLE 8
|
||||
#define PLY_INT_8 9
|
||||
#define PLY_UINT_8 10
|
||||
#define PLY_INT_16 11
|
||||
#define PLY_UINT_16 12
|
||||
#define PLY_INT_32 13
|
||||
#define PLY_UINT_32 14
|
||||
#define PLY_FLOAT_32 15
|
||||
#define PLY_FLOAT_64 16
|
||||
|
||||
#define PLY_END_TYPE 17
|
||||
|
||||
#define PLY_SCALAR 0
|
||||
#define PLY_LIST 1
|
||||
|
||||
|
||||
typedef struct PlyProperty { /* description of a property */
|
||||
|
||||
char *name; /* property name */
|
||||
int external_type; /* file's data type */
|
||||
int internal_type; /* program's data type */
|
||||
int offset; /* offset bytes of prop in a struct */
|
||||
|
||||
int is_list; /* 1 = list, 0 = scalar */
|
||||
int count_external; /* file's count type */
|
||||
int count_internal; /* program's count type */
|
||||
int count_offset; /* offset byte for list count */
|
||||
|
||||
} PlyProperty;
|
||||
|
||||
typedef struct PlyElement { /* description of an element */
|
||||
char *name; /* element name */
|
||||
int num; /* number of elements in this object */
|
||||
int size; /* size of element (bytes) or -1 if variable */
|
||||
int nprops; /* number of properties for this element */
|
||||
PlyProperty **props; /* list of properties in the file */
|
||||
char *store_prop; /* flags: property wanted by user? */
|
||||
int other_offset; /* offset to un-asked-for props, or -1 if none*/
|
||||
int other_size; /* size of other_props structure */
|
||||
} PlyElement;
|
||||
|
||||
typedef struct PlyOtherProp { /* describes other properties in an element */
|
||||
char *name; /* element name */
|
||||
int size; /* size of other_props */
|
||||
int nprops; /* number of properties in other_props */
|
||||
PlyProperty **props; /* list of properties in other_props */
|
||||
} PlyOtherProp;
|
||||
|
||||
typedef struct OtherData { /* for storing other_props for an other element */
|
||||
void *other_props;
|
||||
} OtherData;
|
||||
|
||||
typedef struct OtherElem { /* data for one "other" element */
|
||||
char *elem_name; /* names of other elements */
|
||||
int elem_count; /* count of instances of each element */
|
||||
OtherData **other_data; /* actual property data for the elements */
|
||||
PlyOtherProp *other_props; /* description of the property data */
|
||||
} OtherElem;
|
||||
|
||||
typedef struct PlyOtherElems { /* "other" elements, not interpreted by user */
|
||||
int num_elems; /* number of other elements */
|
||||
OtherElem *other_list; /* list of data for other elements */
|
||||
} PlyOtherElems;
|
||||
|
||||
typedef struct PlyFile { /* description of PLY file */
|
||||
FILE *fp; /* file pointer */
|
||||
int file_type; /* ascii or binary */
|
||||
float version; /* version number of file */
|
||||
int nelems; /* number of elements of object */
|
||||
PlyElement **elems; /* list of elements */
|
||||
int num_comments; /* number of comments */
|
||||
char **comments; /* list of comments */
|
||||
int num_obj_info; /* number of items of object information */
|
||||
char **obj_info; /* list of object info items */
|
||||
PlyElement *which_elem; /* which element we're currently writing */
|
||||
PlyOtherElems *other_elems; /* "other" elements from a PLY file */
|
||||
} PlyFile;
|
||||
|
||||
/* memory allocation */
|
||||
extern char *my_alloc();
|
||||
#define myalloc(mem_size) my_alloc((mem_size), __LINE__, __FILE__)
|
||||
|
||||
#ifndef ALLOCN
|
||||
#define REALLOCN(PTR,TYPE,OLD_N,NEW_N) \
|
||||
{ \
|
||||
if ((OLD_N) == 0) \
|
||||
{ ALLOCN((PTR),TYPE,(NEW_N));} \
|
||||
else \
|
||||
{ \
|
||||
(PTR) = (TYPE *)realloc((PTR),(NEW_N)*sizeof(TYPE)); \
|
||||
if (((PTR) == NULL) && ((NEW_N) != 0)) \
|
||||
{ \
|
||||
fprintf(stderr, "Memory reallocation failed on line %d in %s\n", \
|
||||
__LINE__, __FILE__); \
|
||||
fprintf(stderr, " tried to reallocate %d->%d\n", \
|
||||
(OLD_N), (NEW_N)); \
|
||||
exit(-1); \
|
||||
} \
|
||||
if ((NEW_N)>(OLD_N)) \
|
||||
memset((char *)(PTR)+(OLD_N)*sizeof(TYPE), 0, \
|
||||
((NEW_N)-(OLD_N))*sizeof(TYPE)); \
|
||||
} \
|
||||
}
|
||||
|
||||
#define ALLOCN(PTR,TYPE,N) \
|
||||
{ (PTR) = (TYPE *) calloc(((unsigned)(N)),sizeof(TYPE));\
|
||||
if ((PTR) == NULL) { \
|
||||
fprintf(stderr, "Memory allocation failed on line %d in %s\n", \
|
||||
__LINE__, __FILE__); \
|
||||
exit(-1); \
|
||||
} \
|
||||
}
|
||||
|
||||
|
||||
#define FREE(PTR) { free((PTR)); (PTR) = NULL; }
|
||||
#endif
|
||||
|
||||
|
||||
/*** delcaration of routines ***/
|
||||
|
||||
extern PlyFile *ply_write(FILE *, int, const char **, int);
|
||||
extern PlyFile *ply_open_for_writing(char *, int, const char **, int, float *);
|
||||
extern void ply_describe_element(PlyFile *, char *, int, int, PlyProperty *);
|
||||
extern void ply_describe_property(PlyFile *, char *, PlyProperty *);
|
||||
extern void ply_element_count(PlyFile *, char *, int);
|
||||
extern void ply_header_complete(PlyFile *);
|
||||
extern void ply_put_element_setup(PlyFile *, char *);
|
||||
extern void ply_put_element(PlyFile *, void *);
|
||||
extern void ply_put_comment(PlyFile *, char *);
|
||||
extern void ply_put_obj_info(PlyFile *, char *);
|
||||
extern PlyFile *ply_read(FILE *, int *, char ***);
|
||||
extern PlyFile *ply_open_for_reading( char *, int *, char ***, int *, float *);
|
||||
extern PlyProperty **ply_get_element_description(PlyFile *, char *, int*, int*);
|
||||
extern void ply_get_element_setup( PlyFile *, char *, int, PlyProperty *);
|
||||
extern int ply_get_property(PlyFile *, char *, PlyProperty *);
|
||||
extern PlyOtherProp *ply_get_other_properties(PlyFile *, char *, int);
|
||||
extern void ply_get_element(PlyFile *, void *);
|
||||
extern char **ply_get_comments(PlyFile *, int *);
|
||||
extern char **ply_get_obj_info(PlyFile *, int *);
|
||||
extern void ply_close(PlyFile *);
|
||||
extern void ply_get_info(PlyFile *, float *, int *);
|
||||
extern PlyOtherElems *ply_get_other_element (PlyFile *, char *, int);
|
||||
extern void ply_describe_other_elements ( PlyFile *, PlyOtherElems *);
|
||||
extern void ply_put_other_elements (PlyFile *);
|
||||
extern void ply_free_other_elements (PlyOtherElems *);
|
||||
extern void ply_describe_other_properties(PlyFile *, PlyOtherProp *, int);
|
||||
|
||||
extern int equal_strings(const char *, const char *);
|
||||
|
||||
#ifdef __cplusplus
|
||||
}
|
||||
#endif
|
||||
#include "Geometry.h"
|
||||
#include <vector>
|
||||
|
||||
typedef struct PlyFace
|
||||
{
|
||||
unsigned char nr_vertices;
|
||||
int *vertices;
|
||||
int segment;
|
||||
} PlyFace;
|
||||
static PlyProperty face_props[] =
|
||||
{
|
||||
{ "vertex_indices" , PLY_INT , PLY_INT , offsetof(PlyFace,vertices) , 1 , PLY_UCHAR, PLY_UCHAR , offsetof(PlyFace,nr_vertices) },
|
||||
};
|
||||
|
||||
template< class Real >
|
||||
class PlyVertex
|
||||
{
|
||||
public:
|
||||
const static int Components=3;
|
||||
static PlyProperty Properties[];
|
||||
|
||||
Point3D< Real > point;
|
||||
|
||||
PlyVertex( void ) { ; }
|
||||
PlyVertex( Point3D< Real > p ) { point=p; }
|
||||
PlyVertex operator + ( PlyVertex p ) const { return PlyVertex( point+p.point ); }
|
||||
PlyVertex operator - ( PlyVertex p ) const { return PlyVertex( point-p.point ); }
|
||||
template< class _Real > PlyVertex operator * ( _Real s ) const { return PlyVertex( point*s ); }
|
||||
template< class _Real > PlyVertex operator / ( _Real s ) const { return PlyVertex( point/s ); }
|
||||
PlyVertex& operator += ( PlyVertex p ) { point += p.point ; return *this; }
|
||||
PlyVertex& operator -= ( PlyVertex p ) { point -= p.point ; return *this; }
|
||||
template< class _Real > PlyVertex& operator *= ( _Real s ) { point *= s ; return *this; }
|
||||
template< class _Real > PlyVertex& operator /= ( _Real s ) { point /= s ; return *this; }
|
||||
};
|
||||
template< class Real , class _Real > PlyVertex< Real > operator * ( XForm4x4< _Real > xForm , PlyVertex< Real > v ) { return PlyVertex< Real >( xForm * v.point ); }
|
||||
template<>
|
||||
PlyProperty PlyVertex< float >::Properties[]=
|
||||
{
|
||||
{"x", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyVertex,point.coords[0])), 0, 0, 0, 0},
|
||||
{"y", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyVertex,point.coords[1])), 0, 0, 0, 0},
|
||||
{"z", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyVertex,point.coords[2])), 0, 0, 0, 0}
|
||||
};
|
||||
template<>
|
||||
PlyProperty PlyVertex< double >::Properties[]=
|
||||
{
|
||||
{"x", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyVertex,point.coords[0])), 0, 0, 0, 0},
|
||||
{"y", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyVertex,point.coords[1])), 0, 0, 0, 0},
|
||||
{"z", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyVertex,point.coords[2])), 0, 0, 0, 0}
|
||||
};
|
||||
template< class Real >
|
||||
class PlyValueVertex
|
||||
{
|
||||
public:
|
||||
const static int Components=4;
|
||||
static PlyProperty Properties[];
|
||||
|
||||
Point3D<Real> point;
|
||||
Real value;
|
||||
|
||||
PlyValueVertex( void ) : value( Real(0) ) { ; }
|
||||
PlyValueVertex( Point3D< Real > p , Real v ) : point(p) , value(v) { ; }
|
||||
PlyValueVertex operator + ( PlyValueVertex p ) const { return PlyValueVertex( point+p.point , value+p.value ); }
|
||||
PlyValueVertex operator - ( PlyValueVertex p ) const { return PlyValueVertex( point-p.value , value-p.value ); }
|
||||
template< class _Real > PlyValueVertex operator * ( _Real s ) const { return PlyValueVertex( point*s , Real(value*s) ); }
|
||||
template< class _Real > PlyValueVertex operator / ( _Real s ) const { return PlyValueVertex( point/s , Real(value/s) ); }
|
||||
PlyValueVertex& operator += ( PlyValueVertex p ) { point += p.point , value += p.value ; return *this; }
|
||||
PlyValueVertex& operator -= ( PlyValueVertex p ) { point -= p.point , value -= p.value ; return *this; }
|
||||
template< class _Real > PlyValueVertex& operator *= ( _Real s ) { point *= s , value *= Real(s) ; return *this; }
|
||||
template< class _Real > PlyValueVertex& operator /= ( _Real s ) { point /= s , value /= Real(s) ; return *this; }
|
||||
};
|
||||
template< class Real , class _Real > PlyValueVertex< Real > operator * ( XForm4x4< _Real > xForm , PlyValueVertex< Real > v ) { return PlyValueVertex< Real >( xForm * v.point , v.value ); }
|
||||
template< >
|
||||
PlyProperty PlyValueVertex< float >::Properties[]=
|
||||
{
|
||||
{"x", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyValueVertex,point.coords[0])), 0, 0, 0, 0},
|
||||
{"y", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyValueVertex,point.coords[1])), 0, 0, 0, 0},
|
||||
{"z", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyValueVertex,point.coords[2])), 0, 0, 0, 0},
|
||||
{"value", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyValueVertex,value)), 0, 0, 0, 0}
|
||||
};
|
||||
template< >
|
||||
PlyProperty PlyValueVertex< double >::Properties[]=
|
||||
{
|
||||
{"x", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyValueVertex,point.coords[0])), 0, 0, 0, 0},
|
||||
{"y", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyValueVertex,point.coords[1])), 0, 0, 0, 0},
|
||||
{"z", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyValueVertex,point.coords[2])), 0, 0, 0, 0},
|
||||
{"value", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyValueVertex,value)), 0, 0, 0, 0}
|
||||
};
|
||||
template< class Real >
|
||||
class PlyOrientedVertex
|
||||
{
|
||||
public:
|
||||
const static int Components=6;
|
||||
static PlyProperty Properties[];
|
||||
|
||||
Point3D<Real> point , normal;
|
||||
|
||||
PlyOrientedVertex( void ) { ; }
|
||||
PlyOrientedVertex( Point3D< Real > p , Point3D< Real > n ) : point(p) , normal(n) { ; }
|
||||
PlyOrientedVertex operator + ( PlyOrientedVertex p ) const { return PlyOrientedVertex( point+p.point , normal+p.normal ); }
|
||||
PlyOrientedVertex operator - ( PlyOrientedVertex p ) const { return PlyOrientedVertex( point-p.value , normal-p.normal ); }
|
||||
template< class _Real > PlyOrientedVertex operator * ( _Real s ) const { return PlyOrientedVertex( point*s , normal*s ); }
|
||||
template< class _Real > PlyOrientedVertex operator / ( _Real s ) const { return PlyOrientedVertex( point/s , normal/s ); }
|
||||
PlyOrientedVertex& operator += ( PlyOrientedVertex p ) { point += p.point , normal += p.normal ; return *this; }
|
||||
PlyOrientedVertex& operator -= ( PlyOrientedVertex p ) { point -= p.point , normal -= p.normal ; return *this; }
|
||||
template< class _Real > PlyOrientedVertex& operator *= ( _Real s ) { point *= s , normal *= s ; return *this; }
|
||||
template< class _Real > PlyOrientedVertex& operator /= ( _Real s ) { point /= s , normal /= s ; return *this; }
|
||||
};
|
||||
template< class Real , class _Real > PlyOrientedVertex< Real > operator * ( XForm4x4< _Real > xForm , PlyOrientedVertex< Real > v ) { return PlyOrientedVertex< Real >( xForm * v.point , xForm.inverse().transpose() * v.normal ); }
|
||||
template<>
|
||||
PlyProperty PlyOrientedVertex< float >::Properties[]=
|
||||
{
|
||||
{"x", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyOrientedVertex,point.coords[0])), 0, 0, 0, 0},
|
||||
{"y", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyOrientedVertex,point.coords[1])), 0, 0, 0, 0},
|
||||
{"z", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyOrientedVertex,point.coords[2])), 0, 0, 0, 0},
|
||||
{"nx", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyOrientedVertex,normal.coords[0])), 0, 0, 0, 0},
|
||||
{"ny", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyOrientedVertex,normal.coords[1])), 0, 0, 0, 0},
|
||||
{"nz", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyOrientedVertex,normal.coords[2])), 0, 0, 0, 0}
|
||||
};
|
||||
template<>
|
||||
PlyProperty PlyOrientedVertex< double >::Properties[]=
|
||||
{
|
||||
{"x", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyOrientedVertex,point.coords[0])), 0, 0, 0, 0},
|
||||
{"y", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyOrientedVertex,point.coords[1])), 0, 0, 0, 0},
|
||||
{"z", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyOrientedVertex,point.coords[2])), 0, 0, 0, 0},
|
||||
{"nx", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyOrientedVertex,normal.coords[0])), 0, 0, 0, 0},
|
||||
{"ny", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyOrientedVertex,normal.coords[1])), 0, 0, 0, 0},
|
||||
{"nz", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyOrientedVertex,normal.coords[2])), 0, 0, 0, 0}
|
||||
};
|
||||
template< class Real >
|
||||
class PlyColorVertex
|
||||
{
|
||||
public:
|
||||
const static int Components=6;
|
||||
static PlyProperty Properties[];
|
||||
|
||||
Point3D<Real> point;
|
||||
unsigned char color[3];
|
||||
|
||||
operator Point3D<Real>& () {return point;}
|
||||
operator const Point3D<Real>& () const {return point;}
|
||||
PlyColorVertex(void) {point.coords[0]=point.coords[1]=point.coords[2]=0,color[0]=color[1]=color[2]=0;}
|
||||
PlyColorVertex(const Point3D<Real>& p) {point=p;}
|
||||
};
|
||||
template<>
|
||||
PlyProperty PlyColorVertex< float >::Properties[]=
|
||||
{
|
||||
{"x", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyColorVertex,point.coords[0])), 0, 0, 0, 0},
|
||||
{"y", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyColorVertex,point.coords[1])), 0, 0, 0, 0},
|
||||
{"z", PLY_FLOAT, PLY_FLOAT, int(offsetof(PlyColorVertex,point.coords[2])), 0, 0, 0, 0},
|
||||
{"red", PLY_UCHAR, PLY_UCHAR, int(offsetof(PlyColorVertex,color[0])), 0, 0, 0, 0},
|
||||
{"green", PLY_UCHAR, PLY_UCHAR, int(offsetof(PlyColorVertex,color[1])), 0, 0, 0, 0},
|
||||
{"blue", PLY_UCHAR, PLY_UCHAR, int(offsetof(PlyColorVertex,color[2])), 0, 0, 0, 0}
|
||||
};
|
||||
template<>
|
||||
PlyProperty PlyColorVertex< double >::Properties[]=
|
||||
{
|
||||
{"x", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyColorVertex,point.coords[0])), 0, 0, 0, 0},
|
||||
{"y", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyColorVertex,point.coords[1])), 0, 0, 0, 0},
|
||||
{"z", PLY_DOUBLE, PLY_DOUBLE, int(offsetof(PlyColorVertex,point.coords[2])), 0, 0, 0, 0},
|
||||
{"red", PLY_UCHAR, PLY_UCHAR, int(offsetof(PlyColorVertex,color[0])), 0, 0, 0, 0},
|
||||
{"green", PLY_UCHAR, PLY_UCHAR, int(offsetof(PlyColorVertex,color[1])), 0, 0, 0, 0},
|
||||
{"blue", PLY_UCHAR, PLY_UCHAR, int(offsetof(PlyColorVertex,color[2])), 0, 0, 0, 0}
|
||||
};
|
||||
|
||||
template< class Vertex , class Real >
|
||||
int PlyWritePolygons( char* fileName , CoredMeshData< Vertex >* mesh , int file_type , const Point3D< float >& translate , float scale , char** comments=NULL , int commentNum=0 , XForm4x4< Real > xForm=XForm4x4< Real >::Identity() );
|
||||
|
||||
template< class Vertex , class Real >
|
||||
int PlyWritePolygons( char* fileName , CoredMeshData< Vertex >* mesh , int file_type , char** comments=NULL , int commentNum=0 , XForm4x4< Real > xForm=XForm4x4< Real >::Identity() );
|
||||
|
||||
template<class Vertex>
|
||||
int PlyReadPolygons(char* fileName,
|
||||
std::vector<Vertex>& vertices,std::vector<std::vector<int> >& polygons,
|
||||
PlyProperty* properties,int propertyNum,
|
||||
int& file_type,
|
||||
char*** comments=NULL,int* commentNum=NULL , bool* readFlags=NULL );
|
||||
|
||||
template<class Vertex>
|
||||
int PlyWritePolygons(char* fileName,
|
||||
const std::vector<Vertex>& vertices,const std::vector<std::vector<int> >& polygons,
|
||||
PlyProperty* properties,int propertyNum,
|
||||
int file_type,
|
||||
char** comments=NULL,const int& commentNum=0);
|
||||
|
||||
template<class Vertex>
|
||||
int PlyWritePolygons(char* fileName,
|
||||
const std::vector<Vertex>& vertices , const std::vector< std::vector< int > >& polygons,
|
||||
PlyProperty* properties,int propertyNum,
|
||||
int file_type,
|
||||
char** comments,const int& commentNum)
|
||||
{
|
||||
int nr_vertices=int(vertices.size());
|
||||
int nr_faces=int(polygons.size());
|
||||
float version;
|
||||
const char *elem_names[] = { "vertex" , "face" };
|
||||
PlyFile *ply = ply_open_for_writing( fileName , 2 , elem_names , file_type , &version );
|
||||
if (!ply){return 0;}
|
||||
|
||||
//
|
||||
// describe vertex and face properties
|
||||
//
|
||||
ply_element_count(ply, "vertex", nr_vertices);
|
||||
for(int i=0;i<propertyNum;i++)
|
||||
ply_describe_property(ply, "vertex", &properties[i]);
|
||||
|
||||
ply_element_count(ply, "face", nr_faces);
|
||||
ply_describe_property(ply, "face", &face_props[0]);
|
||||
|
||||
// Write in the comments
|
||||
if(comments && commentNum)
|
||||
for(int i=0;i<commentNum;i++)
|
||||
ply_put_comment(ply,comments[i]);
|
||||
|
||||
ply_header_complete(ply);
|
||||
|
||||
// write vertices
|
||||
ply_put_element_setup(ply, "vertex");
|
||||
for (int i=0; i < int(vertices.size()); i++)
|
||||
ply_put_element(ply, (void *) &vertices[i]);
|
||||
|
||||
// write faces
|
||||
PlyFace ply_face;
|
||||
int maxFaceVerts=3;
|
||||
ply_face.nr_vertices = 3;
|
||||
ply_face.vertices = new int[3];
|
||||
|
||||
ply_put_element_setup(ply, "face");
|
||||
for (int i=0; i < nr_faces; i++)
|
||||
{
|
||||
if(int(polygons[i].size())>maxFaceVerts)
|
||||
{
|
||||
delete[] ply_face.vertices;
|
||||
maxFaceVerts=int(polygons[i].size());
|
||||
ply_face.vertices=new int[maxFaceVerts];
|
||||
}
|
||||
ply_face.nr_vertices=int(polygons[i].size());
|
||||
for(int j=0;j<ply_face.nr_vertices;j++)
|
||||
ply_face.vertices[j]=polygons[i][j];
|
||||
ply_put_element(ply, (void *) &ply_face);
|
||||
}
|
||||
|
||||
delete[] ply_face.vertices;
|
||||
ply_close(ply);
|
||||
return 1;
|
||||
}
|
||||
template<class Vertex>
|
||||
int PlyReadPolygons(char* fileName,
|
||||
std::vector<Vertex>& vertices , std::vector<std::vector<int> >& polygons ,
|
||||
PlyProperty* properties , int propertyNum ,
|
||||
int& file_type ,
|
||||
char*** comments , int* commentNum , bool* readFlags )
|
||||
{
|
||||
int nr_elems;
|
||||
char **elist;
|
||||
float version;
|
||||
int i,j,k;
|
||||
PlyFile* ply;
|
||||
char* elem_name;
|
||||
int num_elems;
|
||||
int nr_props;
|
||||
PlyProperty** plist;
|
||||
PlyFace ply_face;
|
||||
|
||||
ply = ply_open_for_reading(fileName, &nr_elems, &elist, &file_type, &version);
|
||||
if(!ply) return 0;
|
||||
|
||||
if(comments)
|
||||
{
|
||||
(*comments)=new char*[*commentNum+ply->num_comments];
|
||||
for(int i=0;i<ply->num_comments;i++)
|
||||
(*comments)[i]=_strdup(ply->comments[i]);
|
||||
*commentNum=ply->num_comments;
|
||||
}
|
||||
|
||||
for (i=0; i < nr_elems; i++) {
|
||||
elem_name = elist[i];
|
||||
plist = ply_get_element_description(ply, elem_name, &num_elems, &nr_props);
|
||||
if(!plist)
|
||||
{
|
||||
for(i=0;i<nr_elems;i++){
|
||||
free(ply->elems[i]->name);
|
||||
free(ply->elems[i]->store_prop);
|
||||
for(j=0;j<ply->elems[i]->nprops;j++){
|
||||
free(ply->elems[i]->props[j]->name);
|
||||
free(ply->elems[i]->props[j]);
|
||||
}
|
||||
free(ply->elems[i]->props);
|
||||
}
|
||||
for(i=0;i<nr_elems;i++){free(ply->elems[i]);}
|
||||
free(ply->elems);
|
||||
for(i=0;i<ply->num_comments;i++){free(ply->comments[i]);}
|
||||
free(ply->comments);
|
||||
for(i=0;i<ply->num_obj_info;i++){free(ply->obj_info[i]);}
|
||||
free(ply->obj_info);
|
||||
ply_free_other_elements (ply->other_elems);
|
||||
|
||||
for(i=0;i<nr_elems;i++){free(elist[i]);}
|
||||
free(elist);
|
||||
ply_close(ply);
|
||||
return 0;
|
||||
}
|
||||
if (equal_strings("vertex", elem_name))
|
||||
{
|
||||
for( int i=0 ; i<propertyNum ; i++)
|
||||
{
|
||||
int hasProperty = ply_get_property(ply,elem_name,&properties[i]);
|
||||
if( readFlags ) readFlags[i] = (hasProperty!=0);
|
||||
}
|
||||
vertices.resize(num_elems);
|
||||
for (j=0; j < num_elems; j++) ply_get_element (ply, (void *) &vertices[j]);
|
||||
}
|
||||
else if (equal_strings("face", elem_name))
|
||||
{
|
||||
ply_get_property (ply, elem_name, &face_props[0]);
|
||||
polygons.resize(num_elems);
|
||||
for (j=0; j < num_elems; j++)
|
||||
{
|
||||
ply_get_element (ply, (void *) &ply_face);
|
||||
polygons[j].resize(ply_face.nr_vertices);
|
||||
for(k=0;k<ply_face.nr_vertices;k++) polygons[j][k]=ply_face.vertices[k];
|
||||
delete[] ply_face.vertices;
|
||||
} // for, read faces
|
||||
} // if face
|
||||
else{ply_get_other_element (ply, elem_name, num_elems);}
|
||||
|
||||
for(j=0;j<nr_props;j++){
|
||||
free(plist[j]->name);
|
||||
free(plist[j]);
|
||||
}
|
||||
free(plist);
|
||||
} // for each type of element
|
||||
|
||||
for(i=0;i<nr_elems;i++){
|
||||
free(ply->elems[i]->name);
|
||||
free(ply->elems[i]->store_prop);
|
||||
for(j=0;j<ply->elems[i]->nprops;j++){
|
||||
free(ply->elems[i]->props[j]->name);
|
||||
free(ply->elems[i]->props[j]);
|
||||
}
|
||||
if(ply->elems[i]->props && ply->elems[i]->nprops){free(ply->elems[i]->props);}
|
||||
}
|
||||
for(i=0;i<nr_elems;i++){free(ply->elems[i]);}
|
||||
free(ply->elems);
|
||||
for(i=0;i<ply->num_comments;i++){free(ply->comments[i]);}
|
||||
free(ply->comments);
|
||||
for(i=0;i<ply->num_obj_info;i++){free(ply->obj_info[i]);}
|
||||
free(ply->obj_info);
|
||||
ply_free_other_elements (ply->other_elems);
|
||||
|
||||
|
||||
for(i=0;i<nr_elems;i++){free(elist[i]);}
|
||||
free(elist);
|
||||
ply_close(ply);
|
||||
return 1;
|
||||
}
|
||||
|
||||
template< class Vertex , class Real >
|
||||
int PlyWritePolygons( char* fileName , CoredMeshData< Vertex >* mesh , int file_type , const Point3D<float>& translate , float scale , char** comments , int commentNum , XForm4x4< Real > xForm )
|
||||
{
|
||||
int i;
|
||||
int nr_vertices=int(mesh->outOfCorePointCount()+mesh->inCorePoints.size());
|
||||
int nr_faces=mesh->polygonCount();
|
||||
float version;
|
||||
const char *elem_names[] = { "vertex" , "face" };
|
||||
PlyFile *ply = ply_open_for_writing( fileName , 2 , elem_names , file_type , &version );
|
||||
if( !ply ) return 0;
|
||||
|
||||
mesh->resetIterator();
|
||||
|
||||
//
|
||||
// describe vertex and face properties
|
||||
//
|
||||
ply_element_count( ply , "vertex" , nr_vertices );
|
||||
for( int i=0 ; i<Vertex::Components ; i++ ) ply_describe_property( ply , "vertex" , &Vertex::Properties[i] );
|
||||
|
||||
ply_element_count( ply , "face" , nr_faces );
|
||||
ply_describe_property( ply , "face" , &face_props[0] );
|
||||
|
||||
// Write in the comments
|
||||
for( i=0 ; i<commentNum ; i++ ) ply_put_comment( ply , comments[i] );
|
||||
|
||||
ply_header_complete( ply );
|
||||
|
||||
// write vertices
|
||||
ply_put_element_setup( ply , "vertex" );
|
||||
for( i=0 ; i<int( mesh->inCorePoints.size() ) ; i++ )
|
||||
{
|
||||
Vertex vertex = xForm * ( mesh->inCorePoints[i] * scale + translate );
|
||||
ply_put_element(ply, (void *) &vertex);
|
||||
}
|
||||
for( i=0; i<mesh->outOfCorePointCount() ; i++ )
|
||||
{
|
||||
Vertex vertex;
|
||||
mesh->nextOutOfCorePoint( vertex );
|
||||
vertex = xForm * ( vertex * scale +translate );
|
||||
ply_put_element(ply, (void *) &vertex);
|
||||
} // for, write vertices
|
||||
|
||||
// write faces
|
||||
std::vector< CoredVertexIndex > polygon;
|
||||
ply_put_element_setup( ply , "face" );
|
||||
for( i=0 ; i<nr_faces ; i++ )
|
||||
{
|
||||
//
|
||||
// create and fill a struct that the ply code can handle
|
||||
//
|
||||
PlyFace ply_face;
|
||||
mesh->nextPolygon( polygon );
|
||||
ply_face.nr_vertices = int( polygon.size() );
|
||||
ply_face.vertices = new int[ polygon.size() ];
|
||||
for( int i=0 ; i<int(polygon.size()) ; i++ )
|
||||
if( polygon[i].inCore ) ply_face.vertices[i] = polygon[i].idx;
|
||||
else ply_face.vertices[i] = polygon[i].idx + int( mesh->inCorePoints.size() );
|
||||
ply_put_element( ply, (void *) &ply_face );
|
||||
delete[] ply_face.vertices;
|
||||
} // for, write faces
|
||||
|
||||
ply_close( ply );
|
||||
return 1;
|
||||
}
|
||||
template< class Vertex , class Real >
|
||||
int PlyWritePolygons( char* fileName , CoredMeshData< Vertex >* mesh , int file_type , char** comments , int commentNum , XForm4x4< Real > xForm )
|
||||
{
|
||||
int i;
|
||||
int nr_vertices=int(mesh->outOfCorePointCount()+mesh->inCorePoints.size());
|
||||
int nr_faces=mesh->polygonCount();
|
||||
float version;
|
||||
const char *elem_names[] = { "vertex" , "face" };
|
||||
PlyFile *ply = ply_open_for_writing( fileName , 2 , elem_names , file_type , &version );
|
||||
if( !ply ) return 0;
|
||||
|
||||
mesh->resetIterator();
|
||||
|
||||
//
|
||||
// describe vertex and face properties
|
||||
//
|
||||
ply_element_count( ply , "vertex" , nr_vertices );
|
||||
for( int i=0 ; i<Vertex::Components ; i++ ) ply_describe_property( ply , "vertex" , &Vertex::Properties[i] );
|
||||
|
||||
ply_element_count( ply , "face" , nr_faces );
|
||||
ply_describe_property( ply , "face" , &face_props[0] );
|
||||
|
||||
// Write in the comments
|
||||
for( i=0 ; i<commentNum ; i++ ) ply_put_comment( ply , comments[i] );
|
||||
|
||||
ply_header_complete( ply );
|
||||
|
||||
// write vertices
|
||||
ply_put_element_setup( ply , "vertex" );
|
||||
for( i=0 ; i<int( mesh->inCorePoints.size() ) ; i++ )
|
||||
{
|
||||
Vertex vertex = xForm * mesh->inCorePoints[i];
|
||||
ply_put_element(ply, (void *) &vertex);
|
||||
}
|
||||
for( i=0; i<mesh->outOfCorePointCount() ; i++ )
|
||||
{
|
||||
Vertex vertex;
|
||||
mesh->nextOutOfCorePoint( vertex );
|
||||
vertex = xForm * ( vertex );
|
||||
ply_put_element(ply, (void *) &vertex);
|
||||
} // for, write vertices
|
||||
|
||||
// write faces
|
||||
std::vector< CoredVertexIndex > polygon;
|
||||
ply_put_element_setup( ply , "face" );
|
||||
for( i=0 ; i<nr_faces ; i++ )
|
||||
{
|
||||
//
|
||||
// create and fill a struct that the ply code can handle
|
||||
//
|
||||
PlyFace ply_face;
|
||||
mesh->nextPolygon( polygon );
|
||||
ply_face.nr_vertices = int( polygon.size() );
|
||||
ply_face.vertices = new int[ polygon.size() ];
|
||||
for( int i=0 ; i<int(polygon.size()) ; i++ )
|
||||
if( polygon[i].inCore ) ply_face.vertices[i] = polygon[i].idx;
|
||||
else ply_face.vertices[i] = polygon[i].idx + int( mesh->inCorePoints.size() );
|
||||
ply_put_element( ply, (void *) &ply_face );
|
||||
delete[] ply_face.vertices;
|
||||
} // for, write faces
|
||||
|
||||
ply_close( ply );
|
||||
return 1;
|
||||
}
|
||||
inline int PlyDefaultFileType(void){return PLY_ASCII;}
|
||||
|
||||
#endif /* !__PLY_H__ */
|
||||
2727
src/plugins_experimental/filter_screened_poisson/Src/PlyFile.cpp
Executable file
2727
src/plugins_experimental/filter_screened_poisson/Src/PlyFile.cpp
Executable file
File diff suppressed because it is too large
Load Diff
221
src/plugins_experimental/filter_screened_poisson/Src/PlyFile.h
Executable file
221
src/plugins_experimental/filter_screened_poisson/Src/PlyFile.h
Executable file
@ -0,0 +1,221 @@
|
||||
/*
|
||||
|
||||
Header for PLY polygon files.
|
||||
|
||||
- Greg Turk, March 1994
|
||||
|
||||
A PLY file contains a single polygonal _object_.
|
||||
|
||||
An object is composed of lists of _elements_. Typical elements are
|
||||
vertices, faces, edges and materials.
|
||||
|
||||
Each type of element for a given object has one or more _properties_
|
||||
associated with the element type. For instance, a vertex element may
|
||||
have as properties three floating-point values x,y,z and three unsigned
|
||||
chars for red, green and blue.
|
||||
|
||||
---------------------------------------------------------------
|
||||
|
||||
Copyright (c) 1994 The Board of Trustees of The Leland Stanford
|
||||
Junior University. All rights reserved.
|
||||
|
||||
Permission to use, copy, modify and distribute this software and its
|
||||
documentation for any purpose is hereby granted without fee, provided
|
||||
that the above copyright notice and this permission notice appear in
|
||||
all copies of this software and that you do not sell the software.
|
||||
|
||||
THE SOFTWARE IS PROVIDED "AS IS" AND WITHOUT WARRANTY OF ANY KIND,
|
||||
EXPRESS, IMPLIED OR OTHERWISE, INCLUDING WITHOUT LIMITATION, ANY
|
||||
WARRANTY OF MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE.
|
||||
|
||||
*/
|
||||
|
||||
#ifndef PLY_FILE_INCLUDED
|
||||
#define PLY_FILE_INCLUDED
|
||||
|
||||
|
||||
#ifndef WIN32
|
||||
#define _strdup strdup
|
||||
#endif
|
||||
|
||||
#ifdef __cplusplus
|
||||
extern "C" {
|
||||
#endif
|
||||
|
||||
#include <stdlib.h>
|
||||
#include <stdio.h>
|
||||
#include <stddef.h>
|
||||
#include <string.h>
|
||||
|
||||
#define PLY_ASCII 1 /* ascii PLY file */
|
||||
#define PLY_BINARY_BE 2 /* binary PLY file, big endian */
|
||||
#define PLY_BINARY_LE 3 /* binary PLY file, little endian */
|
||||
#define PLY_BINARY_NATIVE 4 /* binary PLY file, same endianness as current architecture */
|
||||
|
||||
#define PLY_OKAY 0 /* ply routine worked okay */
|
||||
#define PLY_ERROR -1 /* error in ply routine */
|
||||
|
||||
/* scalar data types supported by PLY format */
|
||||
|
||||
#define PLY_START_TYPE 0
|
||||
#define PLY_CHAR 1
|
||||
#define PLY_SHORT 2
|
||||
#define PLY_INT 3
|
||||
#define PLY_UCHAR 4
|
||||
#define PLY_USHORT 5
|
||||
#define PLY_UINT 6
|
||||
#define PLY_FLOAT 7
|
||||
#define PLY_DOUBLE 8
|
||||
#define PLY_INT_8 9
|
||||
#define PLY_UINT_8 10
|
||||
#define PLY_INT_16 11
|
||||
#define PLY_UINT_16 12
|
||||
#define PLY_INT_32 13
|
||||
#define PLY_UINT_32 14
|
||||
#define PLY_FLOAT_32 15
|
||||
#define PLY_FLOAT_64 16
|
||||
|
||||
#define PLY_END_TYPE 17
|
||||
|
||||
#define PLY_SCALAR 0
|
||||
#define PLY_LIST 1
|
||||
|
||||
#define PLY_STRIP_COMMENT_HEADER 0
|
||||
|
||||
typedef struct PlyProperty { /* description of a property */
|
||||
|
||||
char *name; /* property name */
|
||||
int external_type; /* file's data type */
|
||||
int internal_type; /* program's data type */
|
||||
int offset; /* offset bytes of prop in a struct */
|
||||
|
||||
int is_list; /* 1 = list, 0 = scalar */
|
||||
int count_external; /* file's count type */
|
||||
int count_internal; /* program's count type */
|
||||
int count_offset; /* offset byte for list count */
|
||||
|
||||
} PlyProperty;
|
||||
|
||||
typedef struct PlyElement { /* description of an element */
|
||||
char *name; /* element name */
|
||||
int num; /* number of elements in this object */
|
||||
int size; /* size of element (bytes) or -1 if variable */
|
||||
int nprops; /* number of properties for this element */
|
||||
PlyProperty **props; /* list of properties in the file */
|
||||
char *store_prop; /* flags: property wanted by user? */
|
||||
int other_offset; /* offset to un-asked-for props, or -1 if none*/
|
||||
int other_size; /* size of other_props structure */
|
||||
} PlyElement;
|
||||
|
||||
typedef struct PlyOtherProp { /* describes other properties in an element */
|
||||
char *name; /* element name */
|
||||
int size; /* size of other_props */
|
||||
int nprops; /* number of properties in other_props */
|
||||
PlyProperty **props; /* list of properties in other_props */
|
||||
} PlyOtherProp;
|
||||
|
||||
typedef struct OtherData { /* for storing other_props for an other element */
|
||||
void *other_props;
|
||||
} OtherData;
|
||||
|
||||
typedef struct OtherElem { /* data for one "other" element */
|
||||
char *elem_name; /* names of other elements */
|
||||
int elem_count; /* count of instances of each element */
|
||||
OtherData **other_data; /* actual property data for the elements */
|
||||
PlyOtherProp *other_props; /* description of the property data */
|
||||
} OtherElem;
|
||||
|
||||
typedef struct PlyOtherElems { /* "other" elements, not interpreted by user */
|
||||
int num_elems; /* number of other elements */
|
||||
OtherElem *other_list; /* list of data for other elements */
|
||||
} PlyOtherElems;
|
||||
|
||||
typedef struct PlyFile { /* description of PLY file */
|
||||
FILE *fp; /* file pointer */
|
||||
int file_type; /* ascii or binary */
|
||||
float version; /* version number of file */
|
||||
int nelems; /* number of elements of object */
|
||||
PlyElement **elems; /* list of elements */
|
||||
int num_comments; /* number of comments */
|
||||
char **comments; /* list of comments */
|
||||
int num_obj_info; /* number of items of object information */
|
||||
char **obj_info; /* list of object info items */
|
||||
PlyElement *which_elem; /* which element we're currently writing */
|
||||
PlyOtherElems *other_elems; /* "other" elements from a PLY file */
|
||||
} PlyFile;
|
||||
|
||||
/* memory allocation */
|
||||
extern char *my_alloc();
|
||||
#define myalloc(mem_size) my_alloc((mem_size), __LINE__, __FILE__)
|
||||
|
||||
#ifndef ALLOCN
|
||||
#define REALLOCN(PTR,TYPE,OLD_N,NEW_N) \
|
||||
{ \
|
||||
if ((OLD_N) == 0) \
|
||||
{ ALLOCN((PTR),TYPE,(NEW_N));} \
|
||||
else \
|
||||
{ \
|
||||
(PTR) = (TYPE *)realloc((PTR),(NEW_N)*sizeof(TYPE)); \
|
||||
if (((PTR) == NULL) && ((NEW_N) != 0)) \
|
||||
{ \
|
||||
fprintf(stderr, "Memory reallocation failed on line %d in %s\n", \
|
||||
__LINE__, __FILE__); \
|
||||
fprintf(stderr, " tried to reallocate %d->%d\n", \
|
||||
(OLD_N), (NEW_N)); \
|
||||
exit(-1); \
|
||||
} \
|
||||
if ((NEW_N)>(OLD_N)) \
|
||||
memset((char *)(PTR)+(OLD_N)*sizeof(TYPE), 0, \
|
||||
((NEW_N)-(OLD_N))*sizeof(TYPE)); \
|
||||
} \
|
||||
}
|
||||
|
||||
#define ALLOCN(PTR,TYPE,N) \
|
||||
{ (PTR) = (TYPE *) calloc(((unsigned)(N)),sizeof(TYPE));\
|
||||
if ((PTR) == NULL) { \
|
||||
fprintf(stderr, "Memory allocation failed on line %d in %s\n", \
|
||||
__LINE__, __FILE__); \
|
||||
exit(-1); \
|
||||
} \
|
||||
}
|
||||
|
||||
|
||||
#define FREE(PTR) { free((PTR)); (PTR) = NULL; }
|
||||
#endif
|
||||
|
||||
|
||||
/*** delcaration of routines ***/
|
||||
|
||||
extern PlyFile *ply_write(FILE *, int, const char **, int);
|
||||
extern PlyFile *ply_open_for_writing(char *, int, const char **, int, float *);
|
||||
extern void ply_describe_element(PlyFile *, char *, int, int, PlyProperty *);
|
||||
extern void ply_describe_property(PlyFile *, char *, PlyProperty *);
|
||||
extern void ply_element_count(PlyFile *, char *, int);
|
||||
extern void ply_header_complete(PlyFile *);
|
||||
extern void ply_put_element_setup(PlyFile *, char *);
|
||||
extern void ply_put_element(PlyFile *, void *);
|
||||
extern void ply_put_comment(PlyFile *, char *);
|
||||
extern void ply_put_obj_info(PlyFile *, char *);
|
||||
extern PlyFile *ply_read(FILE *, int *, char ***);
|
||||
extern PlyFile *ply_open_for_reading( char *, int *, char ***, int *, float *);
|
||||
extern PlyProperty **ply_get_element_description(PlyFile *, char *, int*, int*);
|
||||
extern void ply_get_element_setup( PlyFile *, char *, int, PlyProperty *);
|
||||
extern int ply_get_property(PlyFile *, char *, PlyProperty *);
|
||||
extern PlyOtherProp *ply_get_other_properties(PlyFile *, char *, int);
|
||||
extern void ply_get_element(PlyFile *, void *);
|
||||
extern char **ply_get_comments(PlyFile *, int *);
|
||||
extern char **ply_get_obj_info(PlyFile *, int *);
|
||||
extern void ply_close(PlyFile *);
|
||||
extern void ply_get_info(PlyFile *, float *, int *);
|
||||
extern PlyOtherElems *ply_get_other_element (PlyFile *, char *, int);
|
||||
extern void ply_describe_other_elements ( PlyFile *, PlyOtherElems *);
|
||||
extern void ply_put_other_elements (PlyFile *);
|
||||
extern void ply_free_other_elements (PlyOtherElems *);
|
||||
extern void ply_describe_other_properties(PlyFile *, PlyOtherProp *, int);
|
||||
|
||||
extern int equal_strings(const char *, const char *);
|
||||
|
||||
#ifdef __cplusplus
|
||||
}
|
||||
#endif
|
||||
#endif // PLY_FILE_INCLUDED
|
||||
183
src/plugins_experimental/filter_screened_poisson/Src/PointStream.inl
Executable file
183
src/plugins_experimental/filter_screened_poisson/Src/PointStream.inl
Executable file
@ -0,0 +1,183 @@
|
||||
/*
|
||||
Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
|
||||
All rights reserved.
|
||||
|
||||
Redistribution and use in source and binary forms, with or without modification,
|
||||
are permitted provided that the following conditions are met:
|
||||
|
||||
Redistributions of source code must retain the above copyright notice, this list of
|
||||
conditions and the following disclaimer. Redistributions in binary form must reproduce
|
||||
the above copyright notice, this list of conditions and the following disclaimer
|
||||
in the documentation and/or other materials provided with the distribution.
|
||||
|
||||
Neither the name of the Johns Hopkins University nor the names of its contributors
|
||||
may be used to endorse or promote products derived from this software without specific
|
||||
prior written permission.
|
||||
|
||||
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
|
||||
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
|
||||
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
|
||||
SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
|
||||
INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
|
||||
TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
|
||||
BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
|
||||
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
|
||||
ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
|
||||
DAMAGE.
|
||||
*/
|
||||
#include "Ply.h"
|
||||
|
||||
template< class Real >
|
||||
ASCIIPointStream< Real >::ASCIIPointStream( const char* fileName )
|
||||
{
|
||||
_fp = fopen( fileName , "r" );
|
||||
if( !_fp ) fprintf( stderr , "Failed to open file for reading: %s\n" , fileName ) , exit( 0 );
|
||||
}
|
||||
template< class Real >
|
||||
ASCIIPointStream< Real >::~ASCIIPointStream( void )
|
||||
{
|
||||
fclose( _fp );
|
||||
_fp = NULL;
|
||||
}
|
||||
template< class Real >
|
||||
void ASCIIPointStream< Real >::reset( void ) { fseek( _fp , SEEK_SET , 0 ); }
|
||||
template< class Real >
|
||||
bool ASCIIPointStream< Real >::nextPoint( Point3D< Real >& p , Point3D< Real >& n )
|
||||
{
|
||||
float c[2*DIMENSION];
|
||||
if( fscanf( _fp , " %f %f %f %f %f %f " , &c[0] , &c[1] , &c[2] , &c[3] , &c[4] , &c[5] )!=2*DIMENSION ) return false;
|
||||
p[0] = c[0] , p[1] = c[1] , p[2] = c[2];
|
||||
n[0] = c[3] , n[1] = c[4] , n[2] = c[5];
|
||||
return true;
|
||||
}
|
||||
template< class Real >
|
||||
BinaryPointStream< Real >::BinaryPointStream( const char* fileName )
|
||||
{
|
||||
_pointsInBuffer = _currentPointIndex = 0;
|
||||
_fp = fopen( fileName , "rb" );
|
||||
if( !_fp ) fprintf( stderr , "Failed to open file for reading: %s\n" , fileName ) , exit( 0 );
|
||||
}
|
||||
template< class Real >
|
||||
BinaryPointStream< Real >::~BinaryPointStream( void )
|
||||
{
|
||||
fclose( _fp );
|
||||
_fp = NULL;
|
||||
}
|
||||
template< class Real >
|
||||
void BinaryPointStream< Real >::reset( void )
|
||||
{
|
||||
fseek( _fp , SEEK_SET , 0 );
|
||||
_pointsInBuffer = _currentPointIndex = 0;
|
||||
}
|
||||
template< class Real >
|
||||
bool BinaryPointStream< Real >::nextPoint( Point3D< Real >& p , Point3D< Real >& n )
|
||||
{
|
||||
if( _currentPointIndex<_pointsInBuffer )
|
||||
{
|
||||
p[0] = _pointBuffer[ _currentPointIndex*6+0 ];
|
||||
p[1] = _pointBuffer[ _currentPointIndex*6+1 ];
|
||||
p[2] = _pointBuffer[ _currentPointIndex*6+2 ];
|
||||
n[0] = _pointBuffer[ _currentPointIndex*6+3 ];
|
||||
n[1] = _pointBuffer[ _currentPointIndex*6+4 ];
|
||||
n[2] = _pointBuffer[ _currentPointIndex*6+5 ];
|
||||
_currentPointIndex++;
|
||||
return true;
|
||||
}
|
||||
else
|
||||
{
|
||||
_currentPointIndex = 0;
|
||||
_pointsInBuffer = int( fread( _pointBuffer , sizeof( Real ) * 6 , POINT_BUFFER_SIZE , _fp ) );
|
||||
if( !_pointsInBuffer ) return false;
|
||||
else return nextPoint( p , n );
|
||||
}
|
||||
}
|
||||
|
||||
template< class Real >
|
||||
PLYPointStream< Real >::PLYPointStream( const char* fileName )
|
||||
{
|
||||
_fileName = new char[ strlen( fileName )+1 ];
|
||||
strcpy( _fileName , fileName );
|
||||
_ply = NULL;
|
||||
reset();
|
||||
}
|
||||
template< class Real >
|
||||
void PLYPointStream< Real >::reset( void )
|
||||
{
|
||||
int fileType;
|
||||
float version;
|
||||
PlyProperty** plist;
|
||||
if( _ply ) _free();
|
||||
_ply = ply_open_for_reading( _fileName, &_nr_elems, &_elist, &fileType, &version );
|
||||
if( !_ply )
|
||||
{
|
||||
fprintf( stderr, "[ERROR] Failed to open ply file for reading: %s\n" , _fileName );
|
||||
exit( 0 );
|
||||
}
|
||||
bool foundVertices = false;
|
||||
for( int i=0 ; i<_nr_elems ; i++ )
|
||||
{
|
||||
int num_elems;
|
||||
int nr_props;
|
||||
char* elem_name = _elist[i];
|
||||
plist = ply_get_element_description( _ply , elem_name , &num_elems , &nr_props );
|
||||
if( !plist )
|
||||
{
|
||||
fprintf( stderr , "[ERROR] Failed to get element description: %s\n" , elem_name );
|
||||
exit( 0 );
|
||||
}
|
||||
|
||||
if( equal_strings( "vertex" , elem_name ) )
|
||||
{
|
||||
foundVertices = true;
|
||||
_pCount = num_elems , _pIdx = 0;
|
||||
for( int i=0 ; i<PlyOrientedVertex< Real >::Components ; i++ )
|
||||
if( !ply_get_property( _ply , elem_name , &(PlyOrientedVertex< Real >::Properties[i]) ) )
|
||||
{
|
||||
fprintf( stderr , "[ERROR] Failed to find property in ply file: %s\n" , PlyOrientedVertex< Real >::Properties[i].name );
|
||||
exit( 0 );
|
||||
}
|
||||
}
|
||||
for( int j=0 ; j<nr_props ; j++ )
|
||||
{
|
||||
free( plist[j]->name );
|
||||
free( plist[j] );
|
||||
}
|
||||
free( plist );
|
||||
if( foundVertices ) break;
|
||||
}
|
||||
if( !foundVertices )
|
||||
{
|
||||
fprintf( stderr , "[ERROR] Could not find vertices in ply file\n" );
|
||||
exit( 0 );
|
||||
}
|
||||
}
|
||||
template< class Real >
|
||||
void PLYPointStream< Real >::_free( void )
|
||||
{
|
||||
if( _ply ) ply_close( _ply ) , _ply = NULL;
|
||||
if( _elist )
|
||||
{
|
||||
for( int i=0 ; i<_nr_elems ; i++ ) free( _elist[i] );
|
||||
free( _elist );
|
||||
}
|
||||
}
|
||||
template< class Real >
|
||||
PLYPointStream< Real >::~PLYPointStream( void )
|
||||
{
|
||||
_free();
|
||||
if( _fileName ) delete[] _fileName , _fileName = NULL;
|
||||
}
|
||||
template< class Real >
|
||||
bool PLYPointStream< Real >::nextPoint( Point3D< Real >& p , Point3D< Real >& n )
|
||||
{
|
||||
if( _pIdx<_pCount )
|
||||
{
|
||||
PlyOrientedVertex< Real > op;
|
||||
ply_get_element( _ply, (void *)&op );
|
||||
p = op.point;
|
||||
n = op.normal;
|
||||
_pIdx++;
|
||||
return true;
|
||||
}
|
||||
else return false;
|
||||
}
|
||||
218
src/plugins_experimental/filter_screened_poisson/Src/SparseMatrix.h
Executable file
218
src/plugins_experimental/filter_screened_poisson/Src/SparseMatrix.h
Executable file
@ -0,0 +1,218 @@
|
||||
/*
|
||||
Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
|
||||
All rights reserved.
|
||||
|
||||
Redistribution and use in source and binary forms, with or without modification,
|
||||
are permitted provided that the following conditions are met:
|
||||
|
||||
Redistributions of source code must retain the above copyright notice, this list of
|
||||
conditions and the following disclaimer. Redistributions in binary form must reproduce
|
||||
the above copyright notice, this list of conditions and the following disclaimer
|
||||
in the documentation and/or other materials provided with the distribution.
|
||||
|
||||
Neither the name of the Johns Hopkins University nor the names of its contributors
|
||||
may be used to endorse or promote products derived from this software without specific
|
||||
prior written permission.
|
||||
|
||||
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
|
||||
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
|
||||
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
|
||||
SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
|
||||
INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
|
||||
TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
|
||||
BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
|
||||
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
|
||||
ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
|
||||
DAMAGE.
|
||||
*/
|
||||
|
||||
#ifndef __SPARSEMATRIX_HPP
|
||||
#define __SPARSEMATRIX_HPP
|
||||
|
||||
#define ZERO_TESTING_JACOBI 1
|
||||
|
||||
|
||||
#include "Vector.h"
|
||||
#include "Array.h"
|
||||
|
||||
template <class T>
|
||||
struct MatrixEntry
|
||||
{
|
||||
MatrixEntry( void ) { N =-1; Value = 0; }
|
||||
MatrixEntry( int i ) { N = i; Value = 0; }
|
||||
MatrixEntry( int i , T v ) { N = i; Value = v; }
|
||||
int N;
|
||||
T Value;
|
||||
};
|
||||
|
||||
template<class T> class SparseMatrix
|
||||
{
|
||||
private:
|
||||
bool _contiguous;
|
||||
int _maxEntriesPerRow;
|
||||
void _init( void );
|
||||
public:
|
||||
int rows;
|
||||
Pointer( int ) rowSizes;
|
||||
Pointer( Pointer( MatrixEntry< T > ) ) m_ppElements;
|
||||
Pointer( MatrixEntry< T > ) operator[] ( int idx ) { return m_ppElements[idx]; }
|
||||
ConstPointer( MatrixEntry< T > ) operator[] ( int idx ) const { return m_ppElements[idx]; }
|
||||
|
||||
SparseMatrix( void );
|
||||
SparseMatrix( int rows );
|
||||
SparseMatrix( int rows , int maxEntriesPerRow );
|
||||
void Resize( int rows );
|
||||
void Resize( int rows , int maxEntriesPerRow );
|
||||
void SetRowSize( int row , int count );
|
||||
int Entries( void ) const;
|
||||
|
||||
SparseMatrix( const SparseMatrix& M );
|
||||
~SparseMatrix();
|
||||
|
||||
void SetZero();
|
||||
void SetIdentity();
|
||||
|
||||
SparseMatrix<T>& operator = (const SparseMatrix<T>& M);
|
||||
|
||||
SparseMatrix<T> operator * (const T& V) const;
|
||||
SparseMatrix<T>& operator *= (const T& V);
|
||||
|
||||
|
||||
SparseMatrix<T> operator * (const SparseMatrix<T>& M) const;
|
||||
SparseMatrix<T> Multiply( const SparseMatrix<T>& M ) const;
|
||||
SparseMatrix<T> MultiplyTranspose( const SparseMatrix<T>& Mt ) const;
|
||||
|
||||
template<class T2>
|
||||
Vector<T2> operator * (const Vector<T2>& V) const;
|
||||
template<class T2>
|
||||
Vector<T2> Multiply( const Vector<T2>& V ) const;
|
||||
template<class T2>
|
||||
void Multiply( const Vector<T2>& In , Vector<T2>& Out , int threads=1 ) const;
|
||||
|
||||
|
||||
SparseMatrix<T> Transpose() const;
|
||||
|
||||
static int Solve (const SparseMatrix<T>& M,const Vector<T>& b, int iters,Vector<T>& solution,const T eps=1e-8);
|
||||
|
||||
template<class T2>
|
||||
static int SolveSymmetric( const SparseMatrix<T>& M , const Vector<T2>& b , int iters , Vector<T2>& solution , const T2 eps=1e-8 , int reset=1 , int threads=1 );
|
||||
|
||||
bool write( FILE* fp ) const;
|
||||
bool write( const char* fileName ) const;
|
||||
bool read( FILE* fp );
|
||||
bool read( const char* fileName );
|
||||
|
||||
template< class T2 >
|
||||
static int SolveJacobi( const SparseMatrix<T>& M , const Vector<T2>& diagonal , const Vector<T2>& b , Vector<T2>& x , Vector<T2>& Mx , T2 sor , int threads=1 , int offset=0 );
|
||||
template< class T2 >
|
||||
static int SolveGS( const SparseMatrix<T>& M , const Vector<T2>& diagonal , const Vector<T2>& b , Vector<T2>& x , bool forward , int offset=0 );
|
||||
template< class T2 >
|
||||
static int SolveGS( const std::vector< std::vector< int > >& mcIndices , const SparseMatrix<T>& M , const Vector<T2>& diagonal , const Vector<T2>& b , Vector<T2>& x , bool forward , int threads=1 , int offset=0 );
|
||||
|
||||
template< class T2 >
|
||||
static int SolveJacobi( const SparseMatrix<T>& M , const Vector<T2>& b , Vector<T2>& x , Vector<T2>& Mx , T2 sor , int threads=1 , int offset=0 );
|
||||
template< class T2 >
|
||||
static int SolveGS( const SparseMatrix<T>& M , const Vector<T2>& b , Vector<T2>& x , bool forward , int offset=0 );
|
||||
template< class T2 >
|
||||
static int SolveGS( const std::vector< std::vector< int > >& mcIndices , const SparseMatrix<T>& M , const Vector<T2>& b , Vector<T2>& x , bool forward , int threads=1 , int offset=0 );
|
||||
|
||||
template< class T2 >
|
||||
void getDiagonal( Vector< T2 >& diagonal , int threads=1 , int offset=0 ) const;
|
||||
};
|
||||
|
||||
|
||||
template< class T2 >
|
||||
struct MapReduceVector
|
||||
{
|
||||
private:
|
||||
int _dim;
|
||||
public:
|
||||
std::vector< T2* > out;
|
||||
MapReduceVector( void ) { _dim = 0; }
|
||||
~MapReduceVector( void )
|
||||
{
|
||||
if( _dim ) for( int t=0 ; t<int(out.size()) ; t++ ) delete[] out[t];
|
||||
out.resize( 0 );
|
||||
}
|
||||
T2* operator[]( int t ) { return out[t]; }
|
||||
const T2* operator[]( int t ) const { return out[t]; }
|
||||
int threads( void ) const { return int( out.size() ); }
|
||||
void resize( int threads , int dim )
|
||||
{
|
||||
if( threads!=out.size() || _dim<dim )
|
||||
{
|
||||
for( int t=0 ; t<int(out.size()) ; t++ ) delete[] out[t];
|
||||
out.resize( threads );
|
||||
for( int t=0 ; t<int(out.size()) ; t++ ) out[t] = new T2[dim];
|
||||
_dim = dim;
|
||||
}
|
||||
}
|
||||
|
||||
};
|
||||
|
||||
template< class T >
|
||||
class SparseSymmetricMatrix : public SparseMatrix< T >
|
||||
{
|
||||
public:
|
||||
|
||||
template< class T2 >
|
||||
Vector< T2 > operator * ( const Vector<T2>& V ) const;
|
||||
|
||||
template< class T2 >
|
||||
Vector< T2 > Multiply( const Vector<T2>& V ) const;
|
||||
|
||||
template< class T2 >
|
||||
void Multiply( const Vector<T2>& In, Vector<T2>& Out , bool addDCTerm=false ) const;
|
||||
|
||||
template< class T2 >
|
||||
void Multiply( const Vector<T2>& In, Vector<T2>& Out , MapReduceVector< T2 >& OutScratch , bool addDCTerm=false ) const;
|
||||
|
||||
template< class T2 >
|
||||
void Multiply( const Vector<T2>& In, Vector<T2>& Out , std::vector< T2* >& OutScratch , const std::vector< int >& bounds ) const;
|
||||
|
||||
template< class T2 >
|
||||
static int SolveCG( const SparseSymmetricMatrix<T>& M , const Vector<T2>& b , int iters , Vector<T2>& x , T2 eps=1e-8 , int reset=1 , int threads=0 , bool addDCTerm=false , bool solveNormal=false );
|
||||
|
||||
template< class T2 >
|
||||
static int SolveCG( const SparseSymmetricMatrix<T>& M , const Vector<T2>& b , int iters , Vector<T2>& x , MapReduceVector<T2>& scratch , T2 eps=1e-8 , int reset=1 , bool addDCTerm=false , bool solveNormal=false );
|
||||
#ifdef WIN32
|
||||
template< class T2 >
|
||||
static int SolveCGAtomic( const SparseSymmetricMatrix<T>& M , const Vector<T2>& b , int iters , Vector<T2>& x , T2 eps=1e-8 , int reset=1 , int threads=0 , bool solveNormal=false );
|
||||
#endif // WIN32
|
||||
template< class T2 >
|
||||
static int SolveJacobi( const SparseSymmetricMatrix<T>& M , const Vector<T2>& diagonal , const Vector<T2>& b , Vector<T2>& x , MapReduceVector<T2>& scratch , Vector<T2>& Mx , T2 sor , int reset );
|
||||
template< class T2 >
|
||||
static int SolveJacobi( const SparseSymmetricMatrix<T>& M , const Vector<T2>& b , int iters , Vector<T2>& x , MapReduceVector<T2>& scratch , T2 sor=T2(1.) , int reset=1 );
|
||||
template< class T2 >
|
||||
static int SolveJacobi( const SparseSymmetricMatrix<T>& M , const Vector<T2>& diagonal , const Vector<T2>& b , Vector<T2>& x , Vector<T2>& Mx , T2 sor , int reset );
|
||||
template< class T2 >
|
||||
static int SolveJacobi( const SparseSymmetricMatrix<T>& M , const Vector<T2>& b , int iters , Vector<T2>& x , T2 sor=T2(1.) , int reset=1 );
|
||||
|
||||
enum
|
||||
{
|
||||
ORDERING_UPPER_TRIANGULAR ,
|
||||
ORDERING_LOWER_TRIANGULAR ,
|
||||
ORDERING_NONE
|
||||
};
|
||||
template< class T2 >
|
||||
static int SolveGS( const std::vector< std::vector< int > >& mcIndices , const SparseSymmetricMatrix<T>& M , const Vector<T2>& diagonal , const Vector<T2>& b , Vector<T2>& x , MapReduceVector<T2>& scratch , Vector<T2>& Mx , Vector<T2>& dx , bool forward , int reset );
|
||||
template< class T2 >
|
||||
static int SolveGS( const std::vector< std::vector< int > >& mcIndices , const SparseSymmetricMatrix<T>& M , const Vector<T2>& b , int iters , Vector<T2>& x , MapReduceVector<T2>& scratch , bool forward , int reset=1 );
|
||||
|
||||
template< class T2 >
|
||||
static int SolveGS( const SparseSymmetricMatrix<T>& M , const Vector<T2>& diagonal , const Vector<T2>& b , Vector<T2>& x , MapReduceVector<T2>& scratch , Vector<T2>& Mx , Vector<T2>& dx , bool forward , int reset , int ordering );
|
||||
template< class T2 >
|
||||
static int SolveGS( const SparseSymmetricMatrix<T>& M , const Vector<T2>& b , int iters , Vector<T2>& x , MapReduceVector<T2>& scratch , bool forward , int reset=1 , int ordering=ORDERING_NONE );
|
||||
template< class T2 >
|
||||
static int SolveGS( const SparseSymmetricMatrix<T>& M , const Vector<T2>& diagonal , const Vector<T2>& b , Vector<T2>& x , Vector<T2>& Mx , Vector<T2>& dx , bool forward , int reset , int ordering );
|
||||
template< class T2 >
|
||||
static int SolveGS( const SparseSymmetricMatrix<T>& M , const Vector<T2>& b , int iters , Vector<T2>& x , bool forward , int reset=1 , int ordering=ORDERING_NONE );
|
||||
|
||||
template< class T2 >
|
||||
void getDiagonal( Vector< T2 >& diagonal , int threads=1 ) const;
|
||||
};
|
||||
|
||||
#include "SparseMatrix.inl"
|
||||
|
||||
#endif
|
||||
|
||||
426
src/plugins_experimental/filter_screened_poisson/Src/SurfaceTrimmer.cpp
Executable file
426
src/plugins_experimental/filter_screened_poisson/Src/SurfaceTrimmer.cpp
Executable file
@ -0,0 +1,426 @@
|
||||
/*
|
||||
Copyright (c) 2013, Michael Kazhdan
|
||||
All rights reserved.
|
||||
|
||||
Redistribution and use in source and binary forms, with or without modification,
|
||||
are permitted provided that the following conditions are met:
|
||||
|
||||
Redistributions of source code must retain the above copyright notice, this list of
|
||||
conditions and the following disclaimer. Redistributions in binary form must reproduce
|
||||
the above copyright notice, this list of conditions and the following disclaimer
|
||||
in the documentation and/or other materials provided with the distribution.
|
||||
|
||||
Neither the name of the Johns Hopkins University nor the names of its contributors
|
||||
may be used to endorse or promote products derived from this software without specific
|
||||
prior written permission.
|
||||
|
||||
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
|
||||
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
|
||||
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
|
||||
SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
|
||||
INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
|
||||
TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
|
||||
BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
|
||||
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
|
||||
ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
|
||||
DAMAGE.
|
||||
*/
|
||||
|
||||
#include <stdio.h>
|
||||
#include <stdlib.h>
|
||||
#include <float.h>
|
||||
#include <omp.h>
|
||||
#include <algorithm>
|
||||
#include "CmdLineParser.h"
|
||||
#include "Geometry.h"
|
||||
#include "Ply.h"
|
||||
#include "MAT.h"
|
||||
#include "Time.h"
|
||||
|
||||
#define FOR_RELEASE 1
|
||||
|
||||
cmdLineString In( "in" ) , Out( "out" );
|
||||
cmdLineInt Smooth( "smooth" , 5 );
|
||||
cmdLineFloat Trim( "trim" ) , IslandAreaRatio( "aRatio" , 0.001f );
|
||||
cmdLineFloatArray< 2 > ColorRange( "color" );
|
||||
cmdLineReadable PolygonMesh( "polygonMesh" );
|
||||
|
||||
cmdLineReadable* params[] =
|
||||
{
|
||||
&In , &Out , &Trim , &PolygonMesh , &ColorRange , &Smooth , &IslandAreaRatio
|
||||
};
|
||||
|
||||
void ShowUsage( char* ex )
|
||||
{
|
||||
printf( "Usage: %s\n" , ex );
|
||||
printf( "\t --%s <input polygon mesh>\n" , In.name );
|
||||
printf( "\t[--%s <ouput polygon mesh>]\n" , Out.name );
|
||||
printf( "\t[--%s <smoothing iterations>=%d]\n" , Smooth.name , Smooth.value );
|
||||
printf( "\t[--%s <trimming value>]\n" , Trim.name );
|
||||
printf( "\t[--%s <relative area of islands>=%f]\n" , IslandAreaRatio.name , IslandAreaRatio.value );
|
||||
printf( "\t[--%s]\n" , PolygonMesh.name );
|
||||
#if !FOR_RELEASE
|
||||
printf( "\t[--%s <color range>]\n" , ColorRange.name );
|
||||
#endif // !FOR_RELEASE
|
||||
}
|
||||
|
||||
long long EdgeKey( int key1 , int key2 )
|
||||
{
|
||||
if( key1<key2 ) return ( ( (long long)key1 )<<32 ) | ( (long long)key2 );
|
||||
else return ( ( (long long)key2 )<<32 ) | ( (long long)key1 );
|
||||
}
|
||||
template< class Real >
|
||||
PlyValueVertex< Real > InterpolateVertices( const PlyValueVertex< Real >& v1 , const PlyValueVertex< Real >& v2 , const float& value )
|
||||
{
|
||||
if( v1.value==v2.value ) return (v1+v2)/Real(2.);
|
||||
PlyValueVertex< Real > v;
|
||||
|
||||
Real dx = (v1.value-value)/(v1.value-v2.value);
|
||||
for( int i=0 ; i<3 ; i++ ) v.point.coords[i]=v1.point.coords[i]*(1.f-dx)+v2.point.coords[i]*dx;
|
||||
v.value=v1.value*(1.f-dx)+v2.value*dx;
|
||||
return v;
|
||||
}
|
||||
|
||||
template< class Real >
|
||||
void ColorVertices( const std::vector< PlyValueVertex< Real > >& inVertices , std::vector< PlyColorVertex< Real > >& outVertices , float min , float max )
|
||||
{
|
||||
outVertices.resize( inVertices.size() );
|
||||
for( size_t i=0 ; i<inVertices.size() ; i++ )
|
||||
{
|
||||
outVertices[i].point = inVertices[i].point;
|
||||
float temp = ( inVertices[i].value-min ) / ( max - min );
|
||||
temp = std::max< float >( 0.f , std::min< float >( 1.f , temp ) );
|
||||
temp *= 255;
|
||||
outVertices[i].color[0] = outVertices[i].color[1] = outVertices[i].color[2] = (int)temp;
|
||||
}
|
||||
}
|
||||
template< class Real >
|
||||
void SmoothValues( std::vector< PlyValueVertex< Real > >& vertices , const std::vector< std::vector< int > >& polygons )
|
||||
{
|
||||
std::vector< int > count( vertices.size() );
|
||||
std::vector< Real > sums( vertices.size() , 0 );
|
||||
for( size_t i=0 ; i<polygons.size() ; i++ )
|
||||
{
|
||||
int sz = int(polygons[i].size());
|
||||
for( int j=0 ; j<sz ; j++ )
|
||||
{
|
||||
int j1 = j , j2 = (j+1)%sz;
|
||||
int v1 = polygons[i][j1] , v2 = polygons[i][j2];
|
||||
count[v1]++ , count[v2]++;
|
||||
sums[v1] += vertices[v2].value , sums[v2] += vertices[v1].value;
|
||||
}
|
||||
}
|
||||
for( size_t i=0 ; i<vertices.size() ; i++ ) vertices[i].value = ( sums[i] + vertices[i].value ) / ( count[i] + 1 );
|
||||
}
|
||||
template< class Real >
|
||||
void SmoothValues( std::vector< PlyValueVertex< Real > >& vertices , const std::vector< std::vector< int > >& polygons , Real min , Real max )
|
||||
{
|
||||
std::vector< int > count( vertices.size() );
|
||||
std::vector< Real > sums( vertices.size() , 0 );
|
||||
for( int i=0 ; i<polygons.size() ; i++ ) for( int j=0 ; j<polygons[i].size() ; j++ )
|
||||
{
|
||||
int j1 = j , j2 = (j+1)%polygons[i].size();
|
||||
int v1 = polygons[i][j1] , v2 = polygons[i][j2];
|
||||
if( vertices[v1].value>min && vertices[v1].value<max ) count[v1]++ , sums[v1] += vertices[v2].value;
|
||||
if( vertices[v2].value>min && vertices[v2].value<max ) count[v2]++ , sums[v2] += vertices[v1].value;
|
||||
}
|
||||
for( int i=0 ; i<vertices.size() ; i++ ) vertices[i].value = ( sums[i] + vertices[i].value ) / ( count[i] + 1 );
|
||||
}
|
||||
|
||||
template< class Real >
|
||||
void SplitPolygon
|
||||
(
|
||||
const std::vector< int >& polygon ,
|
||||
std::vector< PlyValueVertex< Real > >& vertices ,
|
||||
std::vector< std::vector< int > >* ltPolygons , std::vector< std::vector< int > >* gtPolygons ,
|
||||
std::vector< bool >* ltFlags , std::vector< bool >* gtFlags ,
|
||||
hash_map< long long , int >& vertexTable ,
|
||||
Real trimValue
|
||||
)
|
||||
{
|
||||
int sz = int( polygon.size() );
|
||||
std::vector< bool > gt( sz );
|
||||
int gtCount = 0;
|
||||
for( int j=0 ; j<sz ; j++ )
|
||||
{
|
||||
gt[j] = ( vertices[ polygon[j] ].value>trimValue );
|
||||
if( gt[j] ) gtCount++;
|
||||
}
|
||||
if ( gtCount==sz ){ if( gtPolygons ) gtPolygons->push_back( polygon ) ; if( gtFlags ) gtFlags->push_back( false ); }
|
||||
else if( gtCount==0 ){ if( ltPolygons ) ltPolygons->push_back( polygon ) ; if( ltFlags ) ltFlags->push_back( false ); }
|
||||
else
|
||||
{
|
||||
int start;
|
||||
for( start=0 ; start<sz ; start++ ) if( gt[start] && !gt[(start+sz-1)%sz] ) break;
|
||||
|
||||
bool gtFlag = true;
|
||||
std::vector< int > poly;
|
||||
|
||||
// Add the initial vertex
|
||||
{
|
||||
int j1 = (start+int(sz)-1)%sz , j2 = start;
|
||||
int v1 = polygon[j1] , v2 = polygon[j2];
|
||||
int vIdx;
|
||||
hash_map< long long , int >::iterator iter = vertexTable.find( EdgeKey( v1 , v2 ) );
|
||||
if( iter==vertexTable.end() )
|
||||
{
|
||||
vertexTable[ EdgeKey( v1 , v2 ) ] = vIdx = int( vertices.size() );
|
||||
vertices.push_back( InterpolateVertices( vertices[v1] , vertices[v2] , trimValue ) );
|
||||
}
|
||||
else vIdx = iter->second;
|
||||
poly.push_back( vIdx );
|
||||
}
|
||||
|
||||
for( int _j=0 ; _j<=sz ; _j++ )
|
||||
{
|
||||
int j1 = (_j+start+sz-1)%sz , j2 = (_j+start)%sz;
|
||||
int v1 = polygon[j1] , v2 = polygon[j2];
|
||||
if( gt[j2]==gtFlag ) poly.push_back( v2 );
|
||||
else
|
||||
{
|
||||
int vIdx;
|
||||
hash_map< long long , int >::iterator iter = vertexTable.find( EdgeKey( v1 , v2 ) );
|
||||
if( iter==vertexTable.end() )
|
||||
{
|
||||
vertexTable[ EdgeKey( v1 , v2 ) ] = vIdx = int( vertices.size() );
|
||||
vertices.push_back( InterpolateVertices( vertices[v1] , vertices[v2] , trimValue ) );
|
||||
}
|
||||
else vIdx = iter->second;
|
||||
poly.push_back( vIdx );
|
||||
if( gtFlag ){ if( gtPolygons ) gtPolygons->push_back( poly ) ; if( ltFlags ) ltFlags->push_back( true ); }
|
||||
else { if( ltPolygons ) ltPolygons->push_back( poly ) ; if( gtFlags ) gtFlags->push_back( true ); }
|
||||
poly.clear() , poly.push_back( vIdx ) , poly.push_back( v2 );
|
||||
gtFlag = !gtFlag;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
template< class Real >
|
||||
void Triangulate( const std::vector< PlyValueVertex< Real > >& vertices , const std::vector< std::vector< int > >& polygons , std::vector< std::vector< int > >& triangles )
|
||||
{
|
||||
triangles.clear();
|
||||
for( size_t i=0 ; i<polygons.size() ; i++ )
|
||||
if( polygons.size()>3 )
|
||||
{
|
||||
MinimalAreaTriangulation< Real > mat;
|
||||
std::vector< Point3D< Real > > _vertices( polygons[i].size() );
|
||||
std::vector< TriangleIndex > _triangles;
|
||||
for( int j=0 ; j<int( polygons[i].size() ) ; j++ ) _vertices[j] = vertices[ polygons[i][j] ].point;
|
||||
mat.GetTriangulation( _vertices , _triangles );
|
||||
|
||||
// Add the triangles to the mesh
|
||||
size_t idx = triangles.size();
|
||||
triangles.resize( idx+_triangles.size() );
|
||||
for( int j=0 ; j<int(_triangles.size()) ; j++ )
|
||||
{
|
||||
triangles[idx+j].resize(3);
|
||||
for( int k=0 ; k<3 ; k++ ) triangles[idx+j][k] = polygons[i][ _triangles[j].idx[k] ];
|
||||
}
|
||||
}
|
||||
else if( polygons[i].size()==3 ) triangles.push_back( polygons[i] );
|
||||
}
|
||||
template< class Vertex >
|
||||
void RemoveHangingVertices( std::vector< Vertex >& vertices , std::vector< std::vector< int > >& polygons )
|
||||
{
|
||||
hash_map< int , int > vMap;
|
||||
std::vector< bool > vertexFlags( vertices.size() , false );
|
||||
for( size_t i=0 ; i<polygons.size() ; i++ ) for( size_t j=0 ; j<polygons[i].size() ; j++ ) vertexFlags[ polygons[i][j] ] = true;
|
||||
int vCount = 0;
|
||||
for( int i=0 ; i<int(vertices.size()) ; i++ ) if( vertexFlags[i] ) vMap[i] = vCount++;
|
||||
for( size_t i=0 ; i<polygons.size() ; i++ ) for( size_t j=0 ; j<polygons[i].size() ; j++ ) polygons[i][j] = vMap[ polygons[i][j] ];
|
||||
|
||||
std::vector< Vertex > _vertices( vCount );
|
||||
for( int i=0 ; i<int(vertices.size()) ; i++ ) if( vertexFlags[i] ) _vertices[ vMap[i] ] = vertices[i];
|
||||
vertices = _vertices;
|
||||
}
|
||||
void SetConnectedComponents( const std::vector< std::vector< int > >& polygons , std::vector< std::vector< int > >& components )
|
||||
{
|
||||
std::vector< int > polygonRoots( polygons.size() );
|
||||
for( size_t i=0 ; i<polygons.size() ; i++ ) polygonRoots[i] = int(i);
|
||||
hash_map< long long , int > edgeTable;
|
||||
for( size_t i=0 ; i<polygons.size() ; i++ )
|
||||
{
|
||||
int sz = int( polygons[i].size() );
|
||||
for( int j=0 ; j<sz ; j++ )
|
||||
{
|
||||
int j1 = j , j2 = (j+1)%sz;
|
||||
int v1 = polygons[i][j1] , v2 = polygons[i][j2];
|
||||
long long eKey = EdgeKey( v1 , v2 );
|
||||
hash_map< long long , int >::iterator iter = edgeTable.find( eKey );
|
||||
if( iter==edgeTable.end() ) edgeTable[ eKey ] = int(i);
|
||||
else
|
||||
{
|
||||
int p = iter->second;
|
||||
while( polygonRoots[p]!=p )
|
||||
{
|
||||
int temp = polygonRoots[p];
|
||||
polygonRoots[p] = int(i);
|
||||
p = temp;
|
||||
}
|
||||
polygonRoots[p] = int(i);
|
||||
}
|
||||
}
|
||||
}
|
||||
for( size_t i=0 ; i<polygonRoots.size() ; i++ )
|
||||
{
|
||||
int p = int(i);
|
||||
while( polygonRoots[p]!=p ) p = polygonRoots[p];
|
||||
int root = p;
|
||||
p = int(i);
|
||||
while( polygonRoots[p]!=p )
|
||||
{
|
||||
int temp = polygonRoots[p];
|
||||
polygonRoots[p] = root;
|
||||
p = temp;
|
||||
}
|
||||
}
|
||||
int cCount = 0;
|
||||
hash_map< int , int > vMap;
|
||||
for( int i= 0 ; i<int(polygonRoots.size()) ; i++ ) if( polygonRoots[i]==i ) vMap[i] = cCount++;
|
||||
components.resize( cCount );
|
||||
for( int i=0 ; i<int(polygonRoots.size()) ; i++ ) components[ vMap[ polygonRoots[i] ] ].push_back(i);
|
||||
}
|
||||
template< class Real >
|
||||
inline Point3D< Real > CrossProduct( Point3D< Real > p1 , Point3D< Real > p2 ){ return Point3D< Real >( p1[1]*p2[2]-p1[2]*p2[1] , p1[2]*p2[0]-p1[0]*p2[2] , p1[0]*p1[1]-p1[1]*p2[0] ); }
|
||||
template< class Real >
|
||||
double TriangleArea( Point3D< Real > v1 , Point3D< Real > v2 , Point3D< Real > v3 )
|
||||
{
|
||||
Point3D< Real > n = CrossProduct( v2-v1 , v3-v1 );
|
||||
return sqrt( n[0]*n[0] + n[1]*n[1] + n[2]*n[2] ) / 2.;
|
||||
}
|
||||
template< class Real >
|
||||
double PolygonArea( const std::vector< PlyValueVertex< Real > >& vertices , const std::vector< int >& polygon )
|
||||
{
|
||||
if( polygon.size()<3 ) return 0.;
|
||||
else if( polygon.size()==3 ) return TriangleArea( vertices[polygon[0]].point , vertices[polygon[1]].point , vertices[polygon[2]].point );
|
||||
else
|
||||
{
|
||||
Point3D< Real > center;
|
||||
for( size_t i=0 ; i<polygon.size() ; i++ ) center += vertices[ polygon[i] ].point;
|
||||
center /= Real( polygon.size() );
|
||||
double area = 0;
|
||||
for( size_t i=0 ; i<polygon.size() ; i++ ) area += TriangleArea( center , vertices[ polygon[i] ].point , vertices[ polygon[ (i+1)%polygon.size() ] ].point );
|
||||
return area;
|
||||
}
|
||||
}
|
||||
int main( int argc , char* argv[] )
|
||||
{
|
||||
int paramNum = sizeof(params)/sizeof(cmdLineReadable*);
|
||||
cmdLineParse( argc-1 , &argv[1] , paramNum , params , 0 );
|
||||
|
||||
#if FOR_RELEASE
|
||||
if( !In.set || !Trim.set )
|
||||
{
|
||||
ShowUsage( argv[0] );
|
||||
return EXIT_FAILURE;
|
||||
}
|
||||
#else // !FOR_RELEASE
|
||||
if( !In.set )
|
||||
{
|
||||
ShowUsage( argv[0] );
|
||||
return EXIT_FAILURE;
|
||||
}
|
||||
#endif // FOR_RELEASE
|
||||
float min , max;
|
||||
std::vector< PlyValueVertex< float > > vertices;
|
||||
std::vector< std::vector< int > > polygons;
|
||||
|
||||
int ft , commentNum = paramNum+2;
|
||||
char** comments;
|
||||
bool readFlags[ PlyValueVertex< float >::Components ];
|
||||
PlyReadPolygons( In.value , vertices , polygons , PlyValueVertex< float >::Properties , PlyValueVertex< float >::Components , ft , &comments , &commentNum , readFlags );
|
||||
if( !readFlags[3] ){ fprintf( stderr , "[ERROR] vertices do not have value flag\n" ) ; return EXIT_FAILURE; }
|
||||
#if 0
|
||||
if( Trim.set ) for( int i=0 ; i<Smooth.value ; i++ ) SmoothValues( vertices , polygons , Trim.value-0.5f , Trim.value+0.5f );
|
||||
else for( int i=0 ; i<Smooth.value ; i++ ) SmoothValues( vertices , polygons );
|
||||
#else
|
||||
for( int i=0 ; i<Smooth.value ; i++ ) SmoothValues( vertices , polygons );
|
||||
#endif
|
||||
min = max = vertices[0].value;
|
||||
for( size_t i=0 ; i<vertices.size() ; i++ ) min = std::min< float >( min , vertices[i].value ) , max = std::max< float >( max , vertices[i].value );
|
||||
printf( "Value Range: [%f,%f]\n" , min , max );
|
||||
|
||||
|
||||
if( Trim.set )
|
||||
{
|
||||
hash_map< long long , int > vertexTable;
|
||||
std::vector< std::vector< int > > ltPolygons , gtPolygons;
|
||||
std::vector< bool > ltFlags , gtFlags;
|
||||
|
||||
for( int i=0 ; i<paramNum+2 ; i++ ) comments[i+commentNum]=new char[1024];
|
||||
sprintf( comments[commentNum++] , "Running Surface Trimmer (V5)" );
|
||||
if( In.set ) sprintf(comments[commentNum++],"\t--%s %s" , In.name , In.value );
|
||||
if( Out.set ) sprintf(comments[commentNum++],"\t--%s %s" , Out.name , Out.value );
|
||||
if( Trim.set ) sprintf(comments[commentNum++],"\t--%s %f" , Trim.name , Trim.value );
|
||||
if( Smooth.set ) sprintf(comments[commentNum++],"\t--%s %d" , Smooth.name , Smooth.value );
|
||||
if( IslandAreaRatio.set ) sprintf(comments[commentNum++],"\t--%s %f" , IslandAreaRatio.name , IslandAreaRatio.value );
|
||||
if( PolygonMesh.set ) sprintf(comments[commentNum++],"\t--%s" , PolygonMesh.name );
|
||||
|
||||
double t=Time();
|
||||
for( size_t i=0 ; i<polygons.size() ; i++ ) SplitPolygon( polygons[i] , vertices , <Polygons , >Polygons , <Flags , >Flags , vertexTable , Trim.value );
|
||||
if( IslandAreaRatio.value>0 )
|
||||
{
|
||||
std::vector< std::vector< int > > _ltPolygons , _gtPolygons;
|
||||
std::vector< std::vector< int > > ltComponents , gtComponents;
|
||||
SetConnectedComponents( ltPolygons , ltComponents );
|
||||
SetConnectedComponents( gtPolygons , gtComponents );
|
||||
std::vector< double > ltAreas( ltComponents.size() , 0. ) , gtAreas( gtComponents.size() , 0. );
|
||||
std::vector< bool > ltComponentFlags( ltComponents.size() , false ) , gtComponentFlags( gtComponents.size() , false );
|
||||
double area = 0.;
|
||||
for( size_t i=0 ; i<ltComponents.size() ; i++ )
|
||||
{
|
||||
for( size_t j=0 ; j<ltComponents[i].size() ; j++ )
|
||||
{
|
||||
ltAreas[i] += PolygonArea( vertices , ltPolygons[ ltComponents[i][j] ] );
|
||||
ltComponentFlags[i] = ( ltComponentFlags[i] || ltFlags[ ltComponents[i][j] ] );
|
||||
}
|
||||
area += ltAreas[i];
|
||||
}
|
||||
for( size_t i=0 ; i<gtComponents.size() ; i++ )
|
||||
{
|
||||
for( size_t j=0 ; j<gtComponents[i].size() ; j++ )
|
||||
{
|
||||
gtAreas[i] += PolygonArea( vertices , gtPolygons[ gtComponents[i][j] ] );
|
||||
gtComponentFlags[i] = ( gtComponentFlags[i] || gtFlags[ gtComponents[i][j] ] );
|
||||
}
|
||||
area += gtAreas[i];
|
||||
}
|
||||
for( size_t i=0 ; i<ltComponents.size() ; i++ )
|
||||
{
|
||||
if( ltAreas[i]<area*IslandAreaRatio.value && ltComponentFlags[i] ) for( size_t j=0 ; j<ltComponents[i].size() ; j++ ) _gtPolygons.push_back( ltPolygons[ ltComponents[i][j] ] );
|
||||
else for( size_t j=0 ; j<ltComponents[i].size() ; j++ ) _ltPolygons.push_back( ltPolygons[ ltComponents[i][j] ] );
|
||||
}
|
||||
for( size_t i=0 ; i<gtComponents.size() ; i++ )
|
||||
{
|
||||
if( gtAreas[i]<area*IslandAreaRatio.value && gtComponentFlags[i] ) for( size_t j=0 ; j<gtComponents[i].size() ; j++ ) _ltPolygons.push_back( gtPolygons[ gtComponents[i][j] ] );
|
||||
else for( size_t j=0 ; j<gtComponents[i].size() ; j++ ) _gtPolygons.push_back( gtPolygons[ gtComponents[i][j] ] );
|
||||
}
|
||||
ltPolygons = _ltPolygons , gtPolygons = _gtPolygons;
|
||||
}
|
||||
if( !PolygonMesh.set )
|
||||
{
|
||||
{
|
||||
std::vector< std::vector< int > > polys = ltPolygons;
|
||||
Triangulate( vertices , ltPolygons , polys ) , ltPolygons = polys;
|
||||
}
|
||||
{
|
||||
std::vector< std::vector< int > > polys = gtPolygons;
|
||||
Triangulate( vertices , gtPolygons , polys ) , gtPolygons = polys;
|
||||
}
|
||||
}
|
||||
|
||||
RemoveHangingVertices( vertices , gtPolygons );
|
||||
sprintf( comments[commentNum++] , "#Trimmed In: %9.1f (s)" , Time()-t );
|
||||
if( Out.set ) PlyWritePolygons( Out.value , vertices , gtPolygons , PlyValueVertex< float >::Properties , PlyValueVertex< float >::Components , ft , comments , commentNum );
|
||||
}
|
||||
else
|
||||
{
|
||||
if( ColorRange.set ) min = ColorRange.values[0] , max = ColorRange.values[1];
|
||||
std::vector< PlyColorVertex< float > > outVertices;
|
||||
ColorVertices( vertices , outVertices , min , max );
|
||||
if( Out.set ) PlyWritePolygons( Out.value , outVertices , polygons , PlyColorVertex< float >::Properties , PlyColorVertex< float >::Components , ft , comments , commentNum );
|
||||
}
|
||||
|
||||
return EXIT_SUCCESS;
|
||||
}
|
||||
|
||||
302
src/plugins_experimental/filter_screened_poisson/Src/Vector.inl
Executable file
302
src/plugins_experimental/filter_screened_poisson/Src/Vector.inl
Executable file
@ -0,0 +1,302 @@
|
||||
/*
|
||||
Copyright (c) 2006, Michael Kazhdan and Matthew Bolitho
|
||||
All rights reserved.
|
||||
|
||||
Redistribution and use in source and binary forms, with or without modification,
|
||||
are permitted provided that the following conditions are met:
|
||||
|
||||
Redistributions of source code must retain the above copyright notice, this list of
|
||||
conditions and the following disclaimer. Redistributions in binary form must reproduce
|
||||
the above copyright notice, this list of conditions and the following disclaimer
|
||||
in the documentation and/or other materials provided with the distribution.
|
||||
|
||||
Neither the name of the Johns Hopkins University nor the names of its contributors
|
||||
may be used to endorse or promote products derived from this software without specific
|
||||
prior written permission.
|
||||
|
||||
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
|
||||
EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO THE IMPLIED WARRANTIES
|
||||
OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT
|
||||
SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
|
||||
INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED
|
||||
TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
|
||||
BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
|
||||
CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
|
||||
ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
|
||||
DAMAGE.
|
||||
*/
|
||||
|
||||
#ifndef __VECTORIMPL_HPP
|
||||
#define __VECTORIMPL_HPP
|
||||
|
||||
////////////
|
||||
// Vector //
|
||||
////////////
|
||||
template<class T>
|
||||
Vector<T>::Vector( void )
|
||||
{
|
||||
m_N = 0;
|
||||
m_pV = NullPointer< T >();
|
||||
}
|
||||
template< class T >
|
||||
Vector< T >::Vector( const Vector<T>& V )
|
||||
{
|
||||
m_N = 0;
|
||||
m_pV = NullPointer< T >();
|
||||
Resize( V.m_N );
|
||||
memcpy( m_pV , V.m_pV , m_N*sizeof(T) );
|
||||
}
|
||||
template<class T>
|
||||
Vector<T>::Vector( size_t N )
|
||||
{
|
||||
m_N=0;
|
||||
m_pV = NullPointer< T >();
|
||||
Resize(N);
|
||||
}
|
||||
template<class T>
|
||||
void Vector<T>::Resize( size_t N )
|
||||
{
|
||||
if( m_N!=N )
|
||||
{
|
||||
if( m_N ) DeletePointer( m_pV );
|
||||
m_N = N;
|
||||
m_pV = NewPointer< T >( N );
|
||||
}
|
||||
if( N ) memset( m_pV , 0 , N*sizeof(T) );
|
||||
}
|
||||
|
||||
template<class T>
|
||||
Vector<T>::Vector( size_t N, ConstPointer( T ) pV )
|
||||
{
|
||||
Resize(N);
|
||||
memcpy( m_pV, pV, N*sizeof(T) );
|
||||
}
|
||||
template<class T>
|
||||
Vector<T>::~Vector(){Resize(0);}
|
||||
template<class T>
|
||||
Vector<T>& Vector<T>::operator = (const Vector& V)
|
||||
{
|
||||
Resize(V.m_N);
|
||||
memcpy( m_pV, V.m_pV, m_N*sizeof(T) );
|
||||
return *this;
|
||||
}
|
||||
template<class T>
|
||||
size_t Vector<T>::Dimensions() const{return m_N;}
|
||||
template<class T>
|
||||
void Vector<T>::SetZero(void){for (size_t i=0; i<m_N; i++){m_pV[i] = T(0);}}
|
||||
template<class T>
|
||||
const T& Vector<T>::operator () (size_t i) const
|
||||
{
|
||||
return m_pV[i];
|
||||
}
|
||||
template<class T>
|
||||
T& Vector<T>::operator () (size_t i)
|
||||
{
|
||||
return m_pV[i];
|
||||
}
|
||||
template<class T>
|
||||
const T& Vector<T>::operator [] (size_t i) const
|
||||
{
|
||||
return m_pV[i];
|
||||
}
|
||||
template<class T>
|
||||
T& Vector<T>::operator [] (size_t i)
|
||||
{
|
||||
return m_pV[i];
|
||||
}
|
||||
template<class T>
|
||||
Vector<T> Vector<T>::operator * (const T& A) const
|
||||
{
|
||||
Vector V(*this);
|
||||
for (size_t i=0; i<m_N; i++) V.m_pV[i] *= A;
|
||||
return V;
|
||||
}
|
||||
template<class T>
|
||||
Vector<T>& Vector<T>::operator *= (const T& A)
|
||||
{
|
||||
for (size_t i=0; i<m_N; i++) m_pV[i] *= A;
|
||||
return *this;
|
||||
}
|
||||
template<class T>
|
||||
Vector<T> Vector<T>::operator / (const T& A) const
|
||||
{
|
||||
Vector V(*this);
|
||||
for (size_t i=0; i<m_N; i++) V.m_pV[i] /= A;
|
||||
return V;
|
||||
}
|
||||
template<class T>
|
||||
Vector<T>& Vector<T>::operator /= (const T& A)
|
||||
{
|
||||
for (size_t i=0; i<m_N; i++) m_pV[i] /= A;
|
||||
return *this;
|
||||
}
|
||||
|
||||
template<class T>
|
||||
Vector<T> Vector<T>::operator + (const T& A) const
|
||||
{
|
||||
Vector V(*this);
|
||||
for (size_t i=0; i<m_N; i++) V.m_pV[i] += A;
|
||||
return V;
|
||||
}
|
||||
template<class T>
|
||||
Vector<T>& Vector<T>::operator += (const T& A)
|
||||
{
|
||||
for (size_t i=0; i<m_N; i++) m_pV[i] += A;
|
||||
return *this;
|
||||
}
|
||||
template<class T>
|
||||
Vector<T> Vector<T>::operator - (const T& A) const
|
||||
{
|
||||
Vector V(*this);
|
||||
for (size_t i=0; i<m_N; i++) V.m_pV[i] -= A;
|
||||
return V;
|
||||
}
|
||||
template<class T>
|
||||
Vector<T>& Vector<T>::operator -= (const T& A)
|
||||
{
|
||||
for (size_t i=0; i<m_N; i++) m_pV[i] -= A;
|
||||
return *this;
|
||||
}
|
||||
|
||||
template<class T>
|
||||
Vector<T> Vector<T>::operator + (const Vector<T>& V0) const
|
||||
{
|
||||
Vector<T> V(m_N);
|
||||
for (size_t i=0; i<m_N; i++) V.m_pV[i] = m_pV[i] + V0.m_pV[i];
|
||||
return V;
|
||||
}
|
||||
template<class T>
|
||||
Vector<T>& Vector<T>::operator += (const Vector<T>& V)
|
||||
{
|
||||
for ( size_t i=0 ; i<m_N ; i++ ) m_pV[i] += V.m_pV[i];
|
||||
return *this;
|
||||
}
|
||||
template<class T>
|
||||
Vector<T> Vector<T>::operator - (const Vector<T>& V0) const
|
||||
{
|
||||
Vector<T> V(m_N);
|
||||
for (size_t i=0; i<m_N; i++) V.m_pV[i] = m_pV[i] - V0.m_pV[i];
|
||||
return V;
|
||||
}
|
||||
template<class T>
|
||||
Vector<T>& Vector<T>::operator -= (const Vector<T>& V)
|
||||
{
|
||||
for (size_t i=0; i<m_N; i++) m_pV[i] -= V.m_pV[i];
|
||||
return *this;
|
||||
}
|
||||
|
||||
template<class T>
|
||||
Vector<T> Vector<T>::operator - (void) const
|
||||
{
|
||||
Vector<T> V(m_N);
|
||||
for (size_t i=0; i<m_N; i++) V.m_pV[i] = -m_pV[i];
|
||||
return V;
|
||||
}
|
||||
|
||||
template< class T >
|
||||
Vector< T >& Vector< T >::Add( const Vector< T >* V , int count )
|
||||
{
|
||||
for( int c=0 ; c<count ; c++ ) for( size_t i=0 ; i<m_N ; i++ ) m_pV[i] += V[c].m_pV[i];
|
||||
return *this;
|
||||
}
|
||||
|
||||
template< class T >
|
||||
Vector< T >& Vector< T >::AddScaled( const Vector<T>& V , const T& scale )
|
||||
{
|
||||
for (size_t i=0; i<m_N; i++) m_pV[i] += V.m_pV[i]*scale;
|
||||
return *this;
|
||||
}
|
||||
template<class T>
|
||||
Vector<T>& Vector<T>::SubtractScaled(const Vector<T>& V,const T& scale)
|
||||
{
|
||||
for (size_t i=0; i<m_N; i++) m_pV[i] -= V.m_pV[i]*scale;
|
||||
return *this;
|
||||
}
|
||||
|
||||
template<class T>
|
||||
void Vector<T>::Add( const Vector<T>& V1 , const T& scale1 , const Vector<T>& V2 , const T& scale2 , Vector<T>& Out )
|
||||
{
|
||||
for( size_t i=0 ; i<V1.m_N ; i++ ) Out.m_pV[i]=V1.m_pV[i]*scale1+V2.m_pV[i]*scale2;
|
||||
}
|
||||
template<class T>
|
||||
void Vector<T>::Add(const Vector<T>& V1,const T& scale1,const Vector<T>& V2,Vector<T>& Out)
|
||||
{
|
||||
for( size_t i=0 ; i<V1.m_N ; i++ ) Out.m_pV[i]=V1.m_pV[i]*scale1+V2.m_pV[i];
|
||||
}
|
||||
|
||||
template<class T>
|
||||
T Vector<T>::Norm( size_t Ln ) const
|
||||
{
|
||||
T N = T();
|
||||
for (size_t i = 0; i<m_N; i++)
|
||||
N += pow(m_pV[i], (T)Ln);
|
||||
return pow(N, (T)1.0/Ln);
|
||||
}
|
||||
template<class T>
|
||||
void Vector<T>::Normalize()
|
||||
{
|
||||
T N = 1.0f/Norm(2);
|
||||
for (size_t i = 0; i<m_N; i++)
|
||||
m_pV[i] *= N;
|
||||
}
|
||||
template< class T >
|
||||
T Vector< T >::Average( void ) const
|
||||
{
|
||||
T N = T();
|
||||
for( size_t i=0 ; i<m_N ; i++ ) N += m_pV[i];
|
||||
return N / m_N;
|
||||
}
|
||||
template<class T>
|
||||
T Vector<T>::Length() const
|
||||
{
|
||||
T N = T();
|
||||
for (size_t i = 0; i<m_N; i++)
|
||||
N += m_pV[i]*m_pV[i];
|
||||
return sqrt(N);
|
||||
}
|
||||
template<class T>
|
||||
T Vector<T>::Dot( const Vector<T>& V ) const
|
||||
{
|
||||
T V0 = T();
|
||||
for( size_t i=0 ; i<m_N ; i++ ) V0 += m_pV[i]*V.m_pV[i];
|
||||
return V0;
|
||||
}
|
||||
|
||||
template< class T >
|
||||
bool Vector< T >::read( const char* fileName )
|
||||
{
|
||||
FILE* fp = fopen( fileName , "rb" );
|
||||
if( !fp ) return false;
|
||||
bool ret = read( fp );
|
||||
fclose( fp );
|
||||
return ret;
|
||||
}
|
||||
template< class T >
|
||||
bool Vector< T >::write( const char* fileName ) const
|
||||
{
|
||||
FILE* fp = fopen( fileName , "wb" );
|
||||
if( !fp ) return false;
|
||||
bool ret = write( fp );
|
||||
fclose( fp );
|
||||
return ret;
|
||||
}
|
||||
template< class T >
|
||||
bool Vector< T >::read( FILE* fp )
|
||||
{
|
||||
int d;
|
||||
if( fread( &d , sizeof(int) , 1 , fp )!=1 ) return false;
|
||||
Resize( d );
|
||||
if( fread( &(*this)[0] , sizeof( T ) , d , fp )!=d ) return false;
|
||||
return true;
|
||||
}
|
||||
template< class T >
|
||||
bool Vector< T >::write( FILE* fp ) const
|
||||
{
|
||||
if( fwrite( &m_N , sizeof( int ) , 1 , fp )!=1 ) return false;
|
||||
if( fwrite( &(*this)[0] , sizeof( T ) , m_N , fp )!=m_N ) return false;
|
||||
return true;
|
||||
}
|
||||
|
||||
|
||||
|
||||
#endif
|
||||
Loading…
x
Reference in New Issue
Block a user