lib
/
VTK-5.0.0
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity
Cloned library of VTK-5.0.0 with extra build files for internal package management.
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
2 Commits
1 Branch
0 Tags
9.4 MiB
 
 
 
 
 
 
Tag: Branch: Tree: main
Branches Tags
${ item.name }
Create tag ${ searchTerm }
Create branch ${ searchTerm }
from 'main'
${ noResults }
VTK-5.0.0/Graphics/vtkSynchronizedTemplatesCut...

542 lines
17 KiB
Raw Permalink Blame History

/*=========================================================================
Program: Visualization Toolkit
Module: $RCSfile: vtkSynchronizedTemplatesCutter3D.cxx,v $
Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
All rights reserved.
See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
This software is distributed WITHOUT ANY WARRANTY; without even
the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
PURPOSE. See the above copyright notice for more information.
=========================================================================*/
#include "vtkSynchronizedTemplatesCutter3D.h"
#include "vtkCellArray.h"
#include "vtkCellData.h"
#include "vtkCharArray.h"
#include "vtkDoubleArray.h"
#include "vtkExtentTranslator.h"
#include "vtkFloatArray.h"
#include "vtkImplicitFunction.h"
#include "vtkInformation.h"
#include "vtkInformationVector.h"
#include "vtkIntArray.h"
#include "vtkLongArray.h"
#include "vtkMath.h"
#include "vtkObjectFactory.h"
#include "vtkPointData.h"
#include "vtkPolyData.h"
#include "vtkShortArray.h"
#include "vtkStreamingDemandDrivenPipeline.h"
#include "vtkStructuredPoints.h"
#include "vtkUnsignedCharArray.h"
#include "vtkUnsignedIntArray.h"
#include "vtkUnsignedLongArray.h"
#include "vtkUnsignedShortArray.h"
#include <math.h>
vtkCxxRevisionMacro(vtkSynchronizedTemplatesCutter3D, "$Revision: 1.5 $");
vtkStandardNewMacro(vtkSynchronizedTemplatesCutter3D);
vtkCxxSetObjectMacro(vtkSynchronizedTemplatesCutter3D,CutFunction,vtkImplicitFunction);
//----------------------------------------------------------------------------
// Description:
// Construct object with initial scalar range (0,1) and single contour value
// of 0.0. The ImageRange are set to extract the first k-plane.
vtkSynchronizedTemplatesCutter3D::vtkSynchronizedTemplatesCutter3D()
{
this->CutFunction = 0;
}
//----------------------------------------------------------------------------
vtkSynchronizedTemplatesCutter3D::~vtkSynchronizedTemplatesCutter3D()
{
this->SetCutFunction(NULL);
}
//----------------------------------------------------------------------------
void vtkSynchronizedTemplatesCutter3DInitializeOutput(
int *ext,vtkImageData *input, vtkPolyData *o)
{
vtkPoints *newPts;
vtkCellArray *newPolys;
long estimatedSize;
estimatedSize = (int) pow ((double)
((ext[1]-ext[0]+1)*(ext[3]-ext[2]+1)*(ext[5]-ext[4]+1)), .75);
if (estimatedSize < 1024)
{
estimatedSize = 1024;
}
newPts = vtkPoints::New();
newPts->Allocate(estimatedSize,estimatedSize);
newPolys = vtkCellArray::New();
newPolys->Allocate(newPolys->EstimateSize(estimatedSize,3));
o->GetPointData()->CopyAllOn();
o->GetPointData()->InterpolateAllocate(input->GetPointData(),
estimatedSize,estimatedSize/2);
o->GetCellData()->CopyAllocate(input->GetCellData(),
estimatedSize,estimatedSize/2);
o->SetPoints(newPts);
newPts->Delete();
o->SetPolys(newPolys);
newPolys->Delete();
}
//----------------------------------------------------------------------------
//
// Contouring filter specialized for images
//
template <class T>
void ContourImage(vtkSynchronizedTemplatesCutter3D *self, int *exExt,
vtkImageData *data, vtkPolyData *output, T *ptr)
{
int *inExt = data->GetExtent();
int xdim = exExt[1] - exExt[0] + 1;
int ydim = exExt[3] - exExt[2] + 1;
double *values = self->GetValues();
int numContours = self->GetNumberOfContours();
T *inPtrX, *inPtrY, *inPtrZ;
T *s0, *s1, *s2, *s3;
int xMin, xMax, yMin, yMax, zMin, zMax;
int xInc, yInc, zInc, scalarZInc;
double *origin = data->GetOrigin();
double *spacing = data->GetSpacing();
int *isect1Ptr, *isect2Ptr;
double y, z, t;
int i, j, k;
int zstep, yisectstep;
int offsets[12];
int *tablePtr;
int idx, vidx;
double x[3], xz[3];
int v0, v1, v2, v3;
vtkIdType ptIds[3];
double value;
// We need to know the edgePointId's for interpolating attributes.
int edgePtId, inCellId, outCellId;
vtkPointData *inPD = data->GetPointData();
vtkCellData *inCD = data->GetCellData();
vtkPointData *outPD = output->GetPointData();
vtkCellData *outCD = output->GetCellData();
// Use to be arguments
vtkPoints *newPts;
vtkCellArray *newPolys;
ptr += self->GetArrayComponent();
vtkSynchronizedTemplatesCutter3DInitializeOutput(exExt, data, output);
newPts = output->GetPoints();
newPolys = output->GetPolys();
// this is an exploded execute extent.
xMin = exExt[0];
xMax = exExt[1];
yMin = exExt[2];
yMax = exExt[3];
zMin = exExt[4];
zMax = exExt[5];
vtkImplicitFunction *func = self->GetCutFunction();
if (!func)
{
return;
}
// implicit functions are 1 component
xInc = 1;
yInc = xInc*(inExt[1]-inExt[0]+1);
zInc = yInc*(inExt[3]-inExt[2]+1);
scalarZInc = zInc;
// Kens increments, probably to do with edge array
zstep = xdim*ydim;
yisectstep = xdim*3;
// compute offsets probably how to get to the edges in the edge array.
offsets[0] = -xdim*3;
offsets[1] = -xdim*3 + 1;
offsets[2] = -xdim*3 + 2;
offsets[3] = -xdim*3 + 4;
offsets[4] = -xdim*3 + 5;
offsets[5] = 0;
offsets[6] = 2;
offsets[7] = 5;
offsets[8] = (zstep - xdim)*3;
offsets[9] = (zstep - xdim)*3 + 1;
offsets[10] = (zstep - xdim)*3 + 4;
offsets[11] = zstep*3;
// allocate storage array
int *isect1 = new int [xdim*ydim*3*2];
// set impossible edges to -1
for (i = 0; i < ydim; i++)
{
isect1[(i+1)*xdim*3-3] = -1;
isect1[(i+1)*xdim*3*2-3] = -1;
}
for (i = 0; i < xdim; i++)
{
isect1[((ydim-1)*xdim + i)*3 + 1] = -1;
isect1[((ydim-1)*xdim + i)*3*2 + 1] = -1;
}
// allocate scalar storage for two slices
T *scalars = 0;
T *scalars1 = 0;
T *scalars2 = 0;
scalars = new T [xdim*ydim*2];
scalars1 = scalars;
scalars2 = scalars + xdim*ydim;
// for each contour
for (vidx = 0; vidx < numContours; vidx++)
{
value = values[vidx];
inPtrZ = ptr;
// fill the first slice
z = origin[2] + spacing[2]*zMin;
x[2] = z;
T *scalarsTmp = scalars1;
scalars1 = scalars2;
scalars2 = scalarsTmp;
for (j = yMin; j <= yMax; j++)
{
x[1] = origin[1] + spacing[1]*j;
for (i = xMin; i <= xMax; i++)
{
x[0] = origin[0] + spacing[0]*i;
*scalarsTmp = func->FunctionValue(x);
scalarsTmp++;
}
}
scalarZInc = -scalarZInc;
//==================================================================
for (k = zMin; k <= zMax; k++)
{
self->UpdateProgress((double)vidx/numContours +
(k-zMin)/((zMax - zMin+1.0)*numContours));
inPtrY = inPtrZ;
// for each slice compute the scalars
z = origin[2] + spacing[2]*(k+1);
x[2] = z;
scalarsTmp = scalars1;
scalars1 = scalars2;
scalars2 = scalarsTmp;
// if not the last slice then get more scalars
if (k < zMax)
{
for (j = yMin; j <= yMax; j++)
{
x[1] = origin[1] + spacing[1]*j;
for (i = xMin; i <= xMax; i++)
{
x[0] = origin[0] + spacing[0]*i;
*scalarsTmp = func->FunctionValue(x);
scalarsTmp++;
}
}
}
inPtrY = scalars1;
scalarZInc = -scalarZInc;
z = origin[2] + spacing[2]*k;
x[2] = z;
// swap the buffers
if (k%2)
{
offsets[8] = (zstep - xdim)*3;
offsets[9] = (zstep - xdim)*3 + 1;
offsets[10] = (zstep - xdim)*3 + 4;
offsets[11] = zstep*3;
isect1Ptr = isect1;
isect2Ptr = isect1 + xdim*ydim*3;
}
else
{
offsets[8] = (-zstep - xdim)*3;
offsets[9] = (-zstep - xdim)*3 + 1;
offsets[10] = (-zstep - xdim)*3 + 4;
offsets[11] = -zstep*3;
isect1Ptr = isect1 + xdim*ydim*3;
isect2Ptr = isect1;
}
for (j = yMin; j <= yMax; j++)
{
// Should not impact perfomance here/
edgePtId = (j-inExt[2])*yInc + (k-inExt[4])*zInc;
// Increments are different for cells. Since the cells are not
// contoured until the second row of templates, subtract 1 from
// i,j,and k. Note: first cube is formed when i=0, j=1, and k=1.
inCellId =
(xMin-inExt[0]) + (inExt[1]-inExt[0])*
( (j-inExt[2]-1) + (k-inExt[4]-1)*(inExt[3]-inExt[2]) );
y = origin[1] + j*spacing[1];
xz[1] = y;
s1 = inPtrY;
v1 = (*s1 < value ? 0 : 1);
inPtrX = inPtrY;
for (i = xMin; i <= xMax; i++)
{
s0 = s1;
v0 = v1;
*isect2Ptr = -1;
*(isect2Ptr + 1) = -1;
*(isect2Ptr + 2) = -1;
if (i < xMax)
{
s1 = (inPtrX + xInc);
v1 = (*s1 < value ? 0 : 1);
if (v0 ^ v1)
{
// watch for degenerate points
if (*s0 == value)
{
if (i > xMin && *(isect2Ptr-3) > -1)
{
*isect2Ptr = *(isect2Ptr-3);
}
else if (j > yMin && *(isect2Ptr - yisectstep + 1) > -1)
{
*isect2Ptr = *(isect2Ptr - yisectstep + 1);
}
else if (k > zMin && *(isect1Ptr+2) > -1)
{
*isect2Ptr = *(isect1Ptr+2);
}
}
else if (*s1 == value)
{
if (j > yMin && *(isect2Ptr - yisectstep +4) > -1)
{
*isect2Ptr = *(isect2Ptr - yisectstep + 4);
}
else if (k > zMin && i < xMax && *(isect1Ptr + 5) > -1)
{
*isect2Ptr = *(isect1Ptr + 5);
}
}
// if the edge has not been set yet then it is a new point
if (*isect2Ptr == -1)
{
t = (value - (double)(*s0)) / ((double)(*s1) - (double)(*s0));
x[0] = origin[0] + spacing[0]*(i+t);
x[1] = y;
*isect2Ptr = newPts->InsertNextPoint(x);
outPD->InterpolateEdge(inPD, *isect2Ptr, edgePtId, edgePtId+1, t);
}
}
}
if (j < yMax)
{
s2 = (inPtrX + yInc);
v2 = (*s2 < value ? 0 : 1);
if (v0 ^ v2)
{
if (*s0 == value)
{
if (*isect2Ptr > -1)
{
*(isect2Ptr + 1) = *isect2Ptr;
}
else if (i > xMin && *(isect2Ptr-3) > -1)
{
*(isect2Ptr + 1) = *(isect2Ptr-3);
}
else if (j > yMin && *(isect2Ptr - yisectstep + 1) > -1)
{
*(isect2Ptr + 1) = *(isect2Ptr - yisectstep + 1);
}
else if (k > zMin && *(isect1Ptr+2) > -1)
{
*(isect2Ptr + 1) = *(isect1Ptr+2);
}
}
else if (*s2 == value && k > zMin && *(isect1Ptr + yisectstep + 2) > -1)
{
*(isect2Ptr+1) = *(isect1Ptr + yisectstep + 2);
}
// if the edge has not been set yet then it is a new point
if (*(isect2Ptr + 1) == -1)
{
t = (value - (double)(*s0)) / ((double)(*s2) - (double)(*s0));
x[0] = origin[0] + spacing[0]*i;
x[1] = y + spacing[1]*t;
*(isect2Ptr + 1) = newPts->InsertNextPoint(x);
outPD->InterpolateEdge(inPD, *(isect2Ptr+1), edgePtId, edgePtId+yInc, t);
}
}
}
if (k < zMax)
{
s3 = (inPtrX + scalarZInc);
v3 = (*s3 < value ? 0 : 1);
if (v0 ^ v3)
{
if (*s0 == value)
{
if (*isect2Ptr > -1)
{
*(isect2Ptr + 2) = *isect2Ptr;
}
else if (*(isect2Ptr+1) > -1)
{
*(isect2Ptr + 2) = *(isect2Ptr+1);
}
else if (i > xMin && *(isect2Ptr-3) > -1)
{
*(isect2Ptr + 2) = *(isect2Ptr-3);
}
else if (j > yMin && *(isect2Ptr - yisectstep + 1) > -1)
{
*(isect2Ptr + 2) = *(isect2Ptr - yisectstep + 1);
}
else if (k > zMin && *(isect1Ptr+2) > -1)
{
*(isect2Ptr + 2) = *(isect1Ptr+2);
}
}
if (*(isect2Ptr + 2) == -1)
{
t = (value - (double)(*s0)) / ((double)(*s3) - (double)(*s0));
xz[0] = origin[0] + spacing[0]*i;
xz[2] = z + spacing[2]*t;
*(isect2Ptr + 2) = newPts->InsertNextPoint(xz);
outPD->InterpolateEdge(inPD, *(isect2Ptr+2), edgePtId, edgePtId+zInc, t);
}
}
}
// To keep track of ids for interpolating attributes.
++edgePtId;
// now add any polys that need to be added
// basically look at the isect values,
// form an index and lookup the polys
if (j > yMin && i < xMax && k > zMin)
{
idx = (v0 ? 4096 : 0);
idx = idx + (*(isect1Ptr - yisectstep) > -1 ? 2048 : 0);
idx = idx + (*(isect1Ptr -yisectstep +1) > -1 ? 1024 : 0);
idx = idx + (*(isect1Ptr -yisectstep +2) > -1 ? 512 : 0);
idx = idx + (*(isect1Ptr -yisectstep +4) > -1 ? 256 : 0);
idx = idx + (*(isect1Ptr -yisectstep +5) > -1 ? 128 : 0);
idx = idx + (*(isect1Ptr) > -1 ? 64 : 0);
idx = idx + (*(isect1Ptr + 2) > -1 ? 32 : 0);
idx = idx + (*(isect1Ptr + 5) > -1 ? 16 : 0);
idx = idx + (*(isect2Ptr -yisectstep) > -1 ? 8 : 0);
idx = idx + (*(isect2Ptr -yisectstep +1) > -1 ? 4 : 0);
idx = idx + (*(isect2Ptr -yisectstep +4) > -1 ? 2 : 0);
idx = idx + (*(isect2Ptr) > -1 ? 1 : 0);
tablePtr = VTK_SYNCHRONIZED_TEMPLATES_3D_TABLE_2
+ VTK_SYNCHRONIZED_TEMPLATES_3D_TABLE_1[idx];
while (*tablePtr != -1)
{
ptIds[0] = *(isect1Ptr + offsets[*tablePtr]);
tablePtr++;
ptIds[1] = *(isect1Ptr + offsets[*tablePtr]);
tablePtr++;
ptIds[2] = *(isect1Ptr + offsets[*tablePtr]);
tablePtr++;
if (ptIds[0] != ptIds[1] &&
ptIds[0] != ptIds[2] &&
ptIds[1] != ptIds[2])
{
outCellId = newPolys->InsertNextCell(3,ptIds);
outCD->CopyData(inCD, inCellId, outCellId);
}
}
}
inPtrX += xInc;
isect2Ptr += 3;
isect1Ptr += 3;
// To keep track of ids for copying cell attributes..
++inCellId;
}
inPtrY += yInc;
}
inPtrZ += zInc;
}
}
delete [] isect1;
if (scalars)
{
delete [] scalars;
}
}
//----------------------------------------------------------------------------
//
// Contouring filter specialized for images (or slices from images)
//
void vtkSynchronizedTemplatesCutter3D::ThreadedExecute(vtkImageData *data,
vtkInformation *outInfo,
int *exExt, int)
{
vtkPolyData *output;
vtkDebugMacro(<< "Executing Cutter3D structured contour");
output = vtkPolyData::SafeDownCast(outInfo->Get(vtkDataObject::DATA_OBJECT()));
if ( exExt[0] >= exExt[1] || exExt[2] >= exExt[3] || exExt[4] >= exExt[5] )
{
vtkDebugMacro(<<"Cutter3D structured contours requires Cutter3D data");
return;
}
// Check data type and execute appropriate function
ContourImage(this, exExt, data, output, (double *)0);
}
//----------------------------------------------------------------------------
int vtkSynchronizedTemplatesCutter3D::RequestData(
vtkInformation *request,
vtkInformationVector **inputVector,
vtkInformationVector *outputVector)
{
// get the info objects
vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
vtkInformation *outInfo = outputVector->GetInformationObject(0);
// get the input and ouptut
vtkImageData *input = vtkImageData::SafeDownCast(
inInfo->Get(vtkDataObject::DATA_OBJECT()));
vtkPolyData *output = vtkPolyData::SafeDownCast(
outInfo->Get(vtkDataObject::DATA_OBJECT()));
// to be safe recompute the
this->RequestUpdateExtent(request,inputVector,outputVector);
// Just call the threaded execute directly.
this->ThreadedExecute(input, outInfo, this->ExecuteExtent, 0);
output->Squeeze();
return 1;
}
//----------------------------------------------------------------------------
void vtkSynchronizedTemplatesCutter3D::PrintSelf(ostream& os, vtkIndent indent)
{
this->Superclass::PrintSelf(os,indent);
os << indent << "Cut Function: " << this->CutFunction << "\n";
}