[prev in list] [next in list] [prev in thread] [next in thread]
List: vtkusers
Subject: [vtkusers] DataSetSurfaceFilter supporting blanking with
From: "Dr. Helmut Schumacher" <hel.schumacher () fz-juelich ! de>
Date: 2010-09-29 12:52:53
Message-ID: 002401cb5fd5$3b19cb10$b14d6130$ () schumacher () fz-juelich ! de
[Download RAW message or body]
[Attachment #2 (multipart/mixed)]
Hello,
I found a very strange solution for my problem.
DataSetMapper uses DataSetSurfaceFilter.
In VTK 5.4.2 and 5.6.1 DataSetSurfaceFilter does NOT support blanking =
for a
StructuredGrid.
But the VTK-Version for ParaView has a method to handle this (for both =
Cell-
and PointBlanking).
Unfortunately the controlling "if"-statement selects this branch for
CellBlanking only.
So I copied the ParaView-Version of DataSetSurfaceFilter to my =
VTK/Graphics.
I changed line 161 of vtkDataSetSurfaceFilter.cxx from
if( grid->GetCellBlanking() )
to
if( grid->GetCellBlanking() || grid->GetPointBlanking() )
and compiled the changed VTK.
As the result I have a VTK with DataSetSurfaceFilter and therefore also
DataSetMapper supporting Cell- and PointBlanking.
I add my version to this mail and hope you will include it in the next
version of VTK.
Thanks in advance
Dr. Helmut Schumacher
Forschungszentrum J=FClich
J=FClich Supercompting Centre
Leo-Brandt-Stra=DFe
52425 J=FClich
Tel.: +49-2461/61-2482
Fax: +49-2461/61-6656
["vtkDataSetSurfaceFilter.cxx" (text/plain)]
/*=========================================================================
Program: Visualization Toolkit
Module: $RCSfile: vtkDataSetSurfaceFilter.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.
=========================================================================
Portions of this software is Copyright © 2008, Pittsburgh
Supercomputing Center (PSC).
Permission to use, copy, and modify this software and its
documentation without fee for personal use or non-commercial use
within your organization is hereby granted, provided that the above
copyright notice is preserved in all copies and that the copyright and
this permission notice appear in supporting documentation. PSC makes
no representations about the suitability of this software for any
purpose. It is provided "as is" without express or implied warranty.
=========================================================================*/
#include "vtkDataSetSurfaceFilter.h"
#include "vtkCellArray.h"
#include "vtkCellData.h"
#include "vtkGenericCell.h"
#include "vtkHexahedron.h"
#include "vtkInformation.h"
#include "vtkInformationVector.h"
#include "vtkMergePoints.h"
#include "vtkObjectFactory.h"
#include "vtkPointData.h"
#include "vtkPolyData.h"
#include "vtkPyramid.h"
#include "vtkRectilinearGrid.h"
#include "vtkStructuredGrid.h"
#include "vtkUniformGrid.h"
#include "vtkStreamingDemandDrivenPipeline.h"
#include "vtkStructuredGridGeometryFilter.h"
#include "vtkStructuredPoints.h"
#include "vtkTetra.h"
#include "vtkUnsignedCharArray.h"
#include "vtkUnstructuredGrid.h"
#include "vtkVoxel.h"
#include "vtkWedge.h"
#include "vtkIdTypeArray.h"
// Helper structure for hashing faces.
struct vtkFastGeomQuadStruct
{
vtkIdType p0;
vtkIdType p1;
vtkIdType p2;
vtkIdType p3;
vtkIdType p4;
vtkIdType p5;
vtkIdType SourceId;
struct vtkFastGeomQuadStruct *Next;
};
vtkCxxRevisionMacro(vtkDataSetSurfaceFilter, "$Revision: 1.10 $");
vtkStandardNewMacro(vtkDataSetSurfaceFilter);
//----------------------------------------------------------------------------
vtkDataSetSurfaceFilter::vtkDataSetSurfaceFilter()
{
this->QuadHash = NULL;
this->PointMap = NULL;
this->QuadHashLength = 0;
this->UseStrips = 0;
this->NumberOfNewCells = 0;
// Quad allocation stuff.
this->FastGeomQuadArrayLength = 0;
this->NumberOfFastGeomQuadArrays = 0;
this->FastGeomQuadArrays = NULL;
this->NextArrayIndex = 0;
this->NextQuadIndex = 0;
this->PieceInvariant = 0;
this->PassThroughCellIds = 0;
this->PassThroughPointIds = 0;
this->OriginalCellIds = NULL;
this->OriginalPointIds = NULL;
}
//----------------------------------------------------------------------------
vtkDataSetSurfaceFilter::~vtkDataSetSurfaceFilter()
{
if (this->QuadHash)
{
this->DeleteQuadHash();
}
if (this->OriginalCellIds != NULL)
{
this->OriginalCellIds->Delete();
this->OriginalCellIds = NULL;
}
}
//----------------------------------------------------------------------------
int vtkDataSetSurfaceFilter::RequestData(
vtkInformation *vtkNotUsed(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
vtkDataSet *input = vtkDataSet::SafeDownCast(
inInfo->Get(vtkDataObject::DATA_OBJECT()));
vtkPolyData *output = vtkPolyData::SafeDownCast(
outInfo->Get(vtkDataObject::DATA_OBJECT()));
vtkIdType numCells = input->GetNumberOfCells();
vtkIdType ext[6];
if (input->CheckAttributes())
{
return 1;
}
if (numCells == 0)
{
return 1;
}
switch (input->GetDataObjectType())
{
case VTK_UNSTRUCTURED_GRID:
{
if (!this->UnstructuredGridExecute(input, output))
{
return 1;
}
output->CheckAttributes();
return 1;
}
case VTK_RECTILINEAR_GRID:
{
vtkRectilinearGrid *grid = vtkRectilinearGrid::SafeDownCast(input);
int* tmpext = grid->GetExtent();
ext[0] = tmpext[0]; ext[1] = tmpext[1];
ext[2] = tmpext[2]; ext[3] = tmpext[3];
ext[4] = tmpext[4]; ext[5] = tmpext[5];
return this->StructuredExecute(grid, output, ext, inInfo);
}
case VTK_STRUCTURED_GRID:
{
vtkStructuredGrid *grid = vtkStructuredGrid::SafeDownCast(input);
if( grid->GetCellBlanking() || grid->GetPointBlanking() )
return this->BlankedStructuredExecute( grid, output );
else
return this->StructuredExecute( grid, output, ext, inInfo );
int* tmpext = grid->GetExtent();
ext[0] = tmpext[0]; ext[1] = tmpext[1];
ext[2] = tmpext[2]; ext[3] = tmpext[3];
ext[4] = tmpext[4]; ext[5] = tmpext[5];
}
case VTK_UNIFORM_GRID:
{
vtkUniformGrid *grid = vtkUniformGrid::SafeDownCast(input);
int* tmpext = grid->GetExtent();
ext[0] = tmpext[0]; ext[1] = tmpext[1];
ext[2] = tmpext[2]; ext[3] = tmpext[3];
ext[4] = tmpext[4]; ext[5] = tmpext[5];
return this->StructuredExecute(grid, output, ext, inInfo);
}
case VTK_STRUCTURED_POINTS:
{
vtkStructuredPoints *image = vtkStructuredPoints::SafeDownCast(input);
int* tmpext = image->GetExtent();
ext[0] = tmpext[0]; ext[1] = tmpext[1];
ext[2] = tmpext[2]; ext[3] = tmpext[3];
ext[4] = tmpext[4]; ext[5] = tmpext[5];
return this->StructuredExecute(image, output, ext, inInfo);
}
case VTK_IMAGE_DATA:
{
vtkImageData *image = vtkImageData::SafeDownCast(input);
int* tmpext = image->GetExtent();
ext[0] = tmpext[0]; ext[1] = tmpext[1];
ext[2] = tmpext[2]; ext[3] = tmpext[3];
ext[4] = tmpext[4]; ext[5] = tmpext[5];
return this->StructuredExecute(image, output, ext, inInfo);
}
case VTK_POLY_DATA:
{
vtkPolyData *inPd = vtkPolyData::SafeDownCast(input);
output->ShallowCopy(inPd);
if (this->PassThroughCellIds)
{
//make a 1:1 mapping
this->OriginalCellIds = vtkIdTypeArray::New();
this->OriginalCellIds->SetName("vtkOriginalCellIds");
this->OriginalCellIds->SetNumberOfComponents(1);
vtkCellData *outputCD = output->GetCellData();
outputCD->AddArray(this->OriginalCellIds);
vtkIdType numTup = output->GetNumberOfCells();
this->OriginalCellIds->SetNumberOfValues(numTup);
for (vtkIdType cId = 0; cId < numTup; cId++)
{
this->OriginalCellIds->SetValue(cId, cId);
}
this->OriginalCellIds->Delete();
this->OriginalCellIds = NULL;
}
if (this->PassThroughPointIds)
{
//make a 1:1 mapping
this->OriginalPointIds = vtkIdTypeArray::New();
this->OriginalPointIds->SetName("vtkOriginalPointIds");
this->OriginalPointIds->SetNumberOfComponents(1);
vtkPointData *outputPD = output->GetPointData();
outputPD->AddArray(this->OriginalPointIds);
vtkIdType numTup = output->GetNumberOfPoints();
this->OriginalPointIds->SetNumberOfValues(numTup);
for (vtkIdType cId = 0; cId < numTup; cId++)
{
this->OriginalPointIds->SetValue(cId, cId);
}
this->OriginalPointIds->Delete();
this->OriginalPointIds = NULL;
}
return 1;
}
default:
return this->DataSetExecute(input, output);
}
}
/*------------------------------------------------------------------------------
vtkDataSetSurfaceFilter::BlankedStructuredExecute
Create outer surface for a structured grid data set that uses blanking.
Note: points are called nodes.
Use node and cell visibility arrays directly to determine which quads
are needed.
This version ignores node and cell IDs; the request to keep the
original cell ids; the request to generate triangle strips, and the
update extent. Maybe later.
------------------------------------------------------------------------------*/
int vtkDataSetSurfaceFilter::BlankedStructuredExecute(
vtkStructuredGrid *SG,
vtkPolyData *PD
)
{
int i, iCV, *dim, nQuad, estNumQuads;
int iXN, iYN, iZN, nXN, nYN, nZN, nXYZN;
int iXC, iYC, iZC, nXC, nYC, nZC, nXYC, nXYZC;
vtkIdType nodes[4];
vtkCellArray *Q;
vtkCellData *inCD, *outCD;
unsigned char *nv, *cv, *NV, *CV;
vtkUnsignedCharArray *NVA, *CVA;
static int neighborCellIndexX[8] = {0,-1,-1,0,0,-1,-1,0};
static int neighborCellIndexY[8] = {0,0,-1,-1,0,0,-1,-1};
static int neighborCellIndexZ[8] = {0,0,0,0,-1,-1,-1,-1};
NV = CV = 0;
////
//// Prepare dimensions and number of things
////
dim = SG->GetDimensions();
nXN = dim[0]; // Nodes
nYN = dim[1];
nZN = dim[2];
nXYZN = nXN * nYN * nZN;
nXC = nXN - 1; // Cells
nYC = nYN - 1;
nZC = nZN - 1;
nXYC = nXC * nYC;
nXYZC = nXYC * nZC;
estNumQuads = 2*nXN*nYN + 2*nXN*nZN + 2*nYN*nZN;
////
//// Input nodes and node data are reused by the output
////
PD->SetPoints( SG->GetPoints() );
PD->GetPointData()->ShallowCopy( SG->GetPointData() );
////
//// Prepare to connect cell data, if any, to new quads
////
inCD = SG->GetCellData();
if( inCD->GetNumberOfArrays() ) {
outCD = PD->GetCellData();
inCD->CopyAllOn();
outCD->CopyAllOn();
outCD->CopyAllocate( inCD, (vtkIdType)estNumQuads );
}
else {
inCD = outCD = 0;
}
////
//// Prepare storage for output's list of quads
////
Q = vtkCellArray::New();
Q->Allocate( Q->EstimateSize( estNumQuads, 4 ) );
PD->SetPolys(Q);
Q->Delete();
////
//// Get copy of input's visibility flags
////
NVA = SG->GetPointVisibilityArray();
if( NVA ) {
NV = new unsigned char [nXYZN];
memcpy( NV, NVA->GetPointer(0), nXYZN );
}
CVA = SG->GetCellVisibilityArray();
CV = new unsigned char [nXYZC];
if( CVA ) {
memcpy( CV, CVA->GetPointer(0), nXYZC );
}
else {
memset( CV, (unsigned char)1, nXYZC );
}
////
//// Ensure self-consistent visibility arrays
//// If node is invisible ensure its cells are also invisible
////
if( NV ) {
nv = NV;
for( iZN=0; iZN<nZN; ++iZN ) {
for( iYN=0; iYN<nYN; ++iYN ) {
for( iXN=0; iXN<nXN; ++iXN, ++nv ) {
if( *nv ) continue;
for( i=0; i<8; ++i ) {
iXC = iXN + neighborCellIndexX[i];
iYC = iYN + neighborCellIndexY[i];
iZC = iZN + neighborCellIndexZ[i];
if( (iXC<0) || (iXC>=nXC) || (iYC<0) ) continue;
if( (iYC>=nYC) || (iZC<0) || (iZC>=nZC) ) continue;
CV[ iXC + nXC*( iYC + nYC*iZC ) ] = 0;
}
}
}
}
}
////
//// Look at each cell to see if it might contribute boundary quads
////
iCV = 0;
nQuad = 0;
for( iZC=0; iZC<nZC; ++iZC ) {
for( iYC=0; iYC<nYC; ++iYC ) {
for( iXC=0; iXC<nXC; ++iXC, ++iCV ) {
//// Only visible cells can contribute quads
if( !CV[ iCV ] ) continue;
//// Check -x side
if( !iXC || !CV[ iCV - 1 ] ) {
nodes[0] = iXC + nXN*( iYC + nYN*iZC );
nodes[1] = iXC + nXN*( iYC + nYN*(iZC+1) );
nodes[2] = iXC + nXN*( iYC+1 + nYN*(iZC+1) );
nodes[3] = iXC + nXN*( iYC+1 + nYN*iZC );
Q->InsertNextCell( 4, nodes );
if( inCD ) outCD->CopyData( inCD, iCV, nQuad );
++nQuad;
}
//// Check +x side
else if( ( iXC == nXC-1 ) || !CV[ iCV + 1 ] ) {
nodes[0] = iXC+1 + nXN*( iYC + nYN*iZC );
nodes[1] = iXC+1 + nXN*( iYC+1 + nYN*iZC );
nodes[2] = iXC+1 + nXN*( iYC+1 + nYN*(iZC+1) );
nodes[3] = iXC+1 + nXN*( iYC + nYN*(iZC+1) );
Q->InsertNextCell( 4, nodes );
if( inCD ) outCD->CopyData( inCD, iCV, nQuad );
++nQuad;
}
//// Check -y side
if( !iYC || !CV[ iCV - nXC ] ) {
nodes[0] = iXC+1 + nXN*( iYC + nYN*iZC );
nodes[1] = iXC+1 + nXN*( iYC + nYN*(iZC+1) );
nodes[2] = iXC + nXN*( iYC + nYN*(iZC+1) );
nodes[3] = iXC + nXN*( iYC + nYN*iZC );
Q->InsertNextCell( 4, nodes );
if( inCD ) outCD->CopyData( inCD, iCV, nQuad );
++nQuad;
}
//// Check +y side
else if( ( iYC == nYC-1 ) || !CV[ iCV + nXC ] ) {
nodes[0] = iXC+1 + nXN*( iYC+1 + nYN*iZC );
nodes[1] = iXC + nXN*( iYC+1 + nYN*iZC );
nodes[2] = iXC + nXN*( iYC+1 + nYN*(iZC+1) );
nodes[3] = iXC+1 + nXN*( iYC+1 + nYN*(iZC+1) );
Q->InsertNextCell( 4, nodes );
if( inCD ) outCD->CopyData( inCD, iCV, nQuad );
++nQuad;
}
//// Check -z side
if( !iZC || !CV[ iCV - nXYC ] ) {
nodes[0] = iXC+1 + nXN*( iYC+1 + nYN*iZC );
nodes[1] = iXC+1 + nXN*( iYC + nYN*iZC );
nodes[2] = iXC + nXN*( iYC + nYN*iZC );
nodes[3] = iXC + nXN*( iYC+1 + nYN*iZC );
Q->InsertNextCell( 4, nodes );
if( inCD ) outCD->CopyData( inCD, iCV, nQuad );
++nQuad;
}
//// Check +z side
else if( ( iZC == nZC-1 ) || !CV[ iCV + nXYC ] ) {
nodes[0] = iXC+1 + nXN*( iYC+1 + nYN*(iZC+1) );
nodes[1] = iXC + nXN*( iYC+1 + nYN*(iZC+1) );
nodes[2] = iXC + nXN*( iYC + nYN*(iZC+1) );
nodes[3] = iXC+1 + nXN*( iYC + nYN*(iZC+1) );
Q->InsertNextCell( 4, nodes );
if( inCD ) outCD->CopyData( inCD, iCV, nQuad );
++nQuad;
}
}
}
}
////
//// Wrapup
////
wrapup:
Q->Squeeze();
if( outCD ) outCD->Squeeze();
if( NV ) delete [] NV;
if( CV ) delete [] CV;
return(1);
}
//----------------------------------------------------------------------------
// It is a pain that structured data sets do not share a common super class
// other than data set, and data set does not allow access to extent!
int vtkDataSetSurfaceFilter::StructuredExecute(vtkDataSet *input,
vtkPolyData *output,
vtkIdType *ext,
vtkInformation *inInfo)
{
const int *wholeExt32;
vtkIdType wholeExt[6];
vtkIdType numPoints, cellArraySize;
vtkCellArray *outStrips;
vtkCellArray *outPolys;
vtkPoints *outPoints;
// Cell Array Size is a pretty good estimate.
// Does not consider direction of strip.
wholeExt32 = inInfo->Get(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT());
wholeExt[0] = wholeExt32[0];
wholeExt[1] = wholeExt32[1];
wholeExt[2] = wholeExt32[2];
wholeExt[3] = wholeExt32[3];
wholeExt[4] = wholeExt32[4];
wholeExt[5] = wholeExt32[5];
// Lets figure out how many cells and points we are going to have.
// It may be overkill comptuing the exact amount, but we can do it, so ...
cellArraySize = numPoints = 0;
// xMin face
if (ext[0] == wholeExt[0] && ext[2] != ext[3] && ext[4] != ext[5] && ext[0] != ext[1])
{
cellArraySize += 2*(ext[3]-ext[2]+1)*(ext[5]-ext[4]+1);
numPoints += (ext[3]-ext[2]+1)*(ext[5]-ext[4]+1);
}
// xMax face
if (ext[1] == wholeExt[1] && ext[2] != ext[3] && ext[4] != ext[5])
{
cellArraySize += 2*(ext[3]-ext[2]+1)*(ext[5]-ext[4]+1);
numPoints += (ext[3]-ext[2]+1)*(ext[5]-ext[4]+1);
}
// yMin face
if (ext[2] == wholeExt[2] && ext[0] != ext[1] && ext[4] != ext[5] && ext[2] != ext[3])
{
cellArraySize += 2*(ext[1]-ext[0]+1)*(ext[5]-ext[4]+1);
numPoints += (ext[1]-ext[0]+1)*(ext[5]-ext[4]+1);
}
// yMax face
if (ext[3] == wholeExt[3] && ext[0] != ext[1] && ext[4] != ext[5])
{
cellArraySize += 2*(ext[1]-ext[0]+1)*(ext[5]-ext[4]+1);
numPoints += (ext[1]-ext[0]+1)*(ext[5]-ext[4]+1);
}
// zMin face
if (ext[4] == wholeExt[4] && ext[0] != ext[1] && ext[2] != ext[3] && ext[4] != ext[5])
{
cellArraySize += 2*(ext[1]-ext[0]+1)*(ext[3]-ext[2]+1);
numPoints += (ext[1]-ext[0]+1)*(ext[3]-ext[2]+1);
}
// zMax face
if (ext[5] == wholeExt[5] && ext[0] != ext[1] && ext[2] != ext[3])
{
cellArraySize += 2*(ext[1]-ext[0]+1)*(ext[3]-ext[2]+1);
numPoints += (ext[1]-ext[0]+1)*(ext[3]-ext[2]+1);
}
int originalPassThroughCellIds = this->PassThroughCellIds;
if (this->UseStrips)
{
outStrips = vtkCellArray::New();
outStrips->Allocate(cellArraySize);
output->SetStrips(outStrips);
outStrips->Delete();
// disable cell ids passing since we are using tstrips.
this->PassThroughCellIds = 0;
}
else
{
outPolys = vtkCellArray::New();
outPolys->Allocate(cellArraySize);
output->SetPolys(outPolys);
outPolys->Delete();
}
outPoints = vtkPoints::New();
outPoints->Allocate(numPoints);
output->SetPoints(outPoints);
outPoints->Delete();
// Allocate attributes for copying.
output->GetPointData()->CopyGlobalIdsOn();
output->GetPointData()->CopyAllocate(input->GetPointData());
output->GetCellData()->CopyGlobalIdsOn();
output->GetCellData()->CopyAllocate(input->GetCellData());
if (this->PassThroughCellIds)
{
this->OriginalCellIds = vtkIdTypeArray::New();
this->OriginalCellIds->SetName("vtkOriginalCellIds");
this->OriginalCellIds->SetNumberOfComponents(1);
output->GetCellData()->AddArray(this->OriginalCellIds);
}
if (this->PassThroughPointIds)
{
this->OriginalPointIds = vtkIdTypeArray::New();
this->OriginalPointIds->SetName("vtkOriginalPointIds");
this->OriginalPointIds->SetNumberOfComponents(1);
output->GetPointData()->AddArray(this->OriginalPointIds);
}
if (this->UseStrips)
{
// xMin face
this->ExecuteFaceStrips(input, output, 0, ext, 0,1,2, inInfo);
// xMax face
this->ExecuteFaceStrips(input, output, 1, ext, 0,2,1, inInfo);
// yMin face
this->ExecuteFaceStrips(input, output, 0, ext, 1,2,0, inInfo);
// yMax face
this->ExecuteFaceStrips(input, output, 1, ext, 1,0,2, inInfo);
// zMin face
this->ExecuteFaceStrips(input, output, 0, ext, 2,0,1, inInfo);
// zMax face
this->ExecuteFaceStrips(input, output, 1, ext, 2,1,0, inInfo);
}
else
{
// xMin face
this->ExecuteFaceQuads(input, output, 0, ext, 0,1,2, inInfo);
// xMax face
this->ExecuteFaceQuads(input, output, 1, ext, 0,2,1, inInfo);
// yMin face
this->ExecuteFaceQuads(input, output, 0, ext, 1,2,0, inInfo);
// yMax face
this->ExecuteFaceQuads(input, output, 1, ext, 1,0,2, inInfo);
// zMin face
this->ExecuteFaceQuads(input, output, 0, ext, 2,0,1, inInfo);
// zMax face
this->ExecuteFaceQuads(input, output, 1, ext, 2,1,0, inInfo);
}
output->Squeeze();
if (this->OriginalCellIds != NULL)
{
this->OriginalCellIds->Delete();
this->OriginalCellIds = NULL;
}
if (this->OriginalPointIds != NULL)
{
this->OriginalPointIds->Delete();
this->OriginalPointIds = NULL;
}
this->PassThroughCellIds = originalPassThroughCellIds;
return 1;
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::ExecuteFaceStrips(vtkDataSet *input,
vtkPolyData *output,
int maxFlag, vtkIdType *ext,
int aAxis, int bAxis, int cAxis,
vtkInformation *inInfo)
{
vtkPoints *outPts;
vtkCellArray *outStrips;
vtkPointData *inPD, *outPD;
int *wholeExt;
vtkIdType pInc[3];
vtkIdType qInc[3];
vtkIdType ptCInc[3];
vtkIdType cOutInc;
double pt[3];
vtkIdType inStartPtId;
vtkIdType outStartPtId;
vtkIdType outPtId;
vtkIdType inId, outId;
vtkIdType ib, ic;
int aA2, bA2, cA2;
int rotatedFlag;
vtkIdType *stripArray;
vtkIdType stripArrayIdx;
outPts = output->GetPoints();
outPD = output->GetPointData();
inPD = input->GetPointData();
wholeExt = inInfo->Get(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT());
pInc[0] = 1;
pInc[1] = (ext[1]-ext[0]+1);
pInc[2] = (ext[3]-ext[2]+1) * pInc[1];
// quad increments (cell incraments, but cInc could be confused with c axis).
qInc[0] = 1;
qInc[1] = ext[1]-ext[0];
qInc[2] = (ext[3]-ext[2]) * qInc[1];
ptCInc[0] = 1;
ptCInc[1] = ext[1]-ext[0];
if (ptCInc[1] == 0)
{
ptCInc[1] = 1;
}
ptCInc[2] = (ext[3]-ext[2]);
if (ptCInc[2] == 0)
{
ptCInc[2] = 1;
}
ptCInc[2] = ptCInc[2] * ptCInc[1];
// Tempoprary variables to avoid many multiplications.
aA2 = aAxis * 2;
bA2 = bAxis * 2;
cA2 = cAxis * 2;
// We might as well put the test for this face here.
if (ext[bA2] == ext[bA2+1] || ext[cA2] == ext[cA2+1])
{
return;
}
if (maxFlag)
{ // max faces have a slightly different condition to avoid coincident faces.
if (ext[aA2] == ext[aA2+1] || ext[aA2+1] < wholeExt[aA2+1])
{
return;
}
}
else
{
if (ext[aA2] > wholeExt[aA2])
{
return;
}
}
// Lets rotate the image to make b the longest axis.
// This will make the tri strips longer.
rotatedFlag = 0;
if (ext[bA2+1]-ext[bA2] < ext[cA2+1]-ext[cA2])
{
int tmp;
rotatedFlag = 1;
tmp = cAxis;
cAxis = bAxis;
bAxis = tmp;
bA2 = bAxis * 2;
cA2 = cAxis * 2;
}
// Assuming no ghost cells ...
inStartPtId = 0;
if( maxFlag)
{
inStartPtId = pInc[aAxis]*(ext[aA2+1]-ext[aA2]);
}
vtkIdType outCellId = 0;
vtkIdType inStartCellId = 0;
vtkIdType inCellId = 0;
if (this->PassThroughCellIds)
{
outCellId = this->OriginalCellIds->GetNumberOfTuples();
if (maxFlag && ext[aA2] < ext[1+aA2])
{
inStartCellId = qInc[aAxis]*(ext[aA2+1]-ext[aA2]-1);
}
}
outStartPtId = outPts->GetNumberOfPoints();
// Make the points for this face.
for (ic = ext[cA2]; ic <= ext[cA2+1]; ++ic)
{
for (ib = ext[bA2]; ib <= ext[bA2+1]; ++ib)
{
inId = inStartPtId + (ib-ext[bA2])*pInc[bAxis]
+ (ic-ext[cA2])*pInc[cAxis];
input->GetPoint(inId, pt);
outId = outPts->InsertNextPoint(pt);
// Copy point data.
outPD->CopyData(inPD,inId,outId);
this->RecordOrigPointId(outId, inId);
}
}
// Do the cells.
cOutInc = ext[bA2+1] - ext[bA2] + 1;
// Tri Strips (no cell data ...).
// Allocate the temporary array used to create the tri strips.
stripArray = new vtkIdType[2*(ext[bA2+1]-ext[bA2]+1)];
// Make the cells for this face.
outStrips = output->GetStrips();
for (ic = ext[cA2]; ic < ext[cA2+1]; ++ic)
{
// Fill in the array describing the strips.
stripArrayIdx = 0;
outPtId = outStartPtId + (ic-ext[cA2])*cOutInc;
if (rotatedFlag)
{
for (ib = ext[bA2]; ib <= ext[bA2+1]; ++ib)
{
stripArray[stripArrayIdx++] = outPtId+cOutInc;
stripArray[stripArrayIdx++] = outPtId;
++outPtId;
if (this->PassThroughCellIds && ib != ext[bA2])
{
//Record the two triangular output cells just defined
//both belong to the same input quad cell
inCellId =
inStartCellId +
(ib-ext[bA2]-1)*ptCInc[bAxis] +
(ic-ext[cA2])*ptCInc[cAxis];
this->RecordOrigCellId(outCellId++, inCellId);
this->RecordOrigCellId(outCellId++, inCellId);
}
}
}
else
{ // Faster to justto dupicate the inner most loop.
for (ib = ext[bA2]; ib <= ext[bA2+1]; ++ib)
{
stripArray[stripArrayIdx++] = outPtId;
stripArray[stripArrayIdx++] = outPtId+cOutInc;
++outPtId;
if (this->PassThroughCellIds && ib != ext[bA2])
{
//Record the two triangular output cells just defined
//both belong to the same input quad cell
inCellId =
inStartCellId +
(ib-ext[bA2]-1)*ptCInc[bAxis] +
(ic-ext[cA2])*ptCInc[cAxis];
this->RecordOrigCellId(outCellId++, inCellId);
this->RecordOrigCellId(outCellId++, inCellId);
}
}
}
outStrips->InsertNextCell(stripArrayIdx, stripArray);
}
delete [] stripArray;
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::ExecuteFaceQuads(vtkDataSet *input,
vtkPolyData *output,
int maxFlag, vtkIdType *ext,
int aAxis, int bAxis, int cAxis,
vtkInformation *inInfo)
{
vtkPoints *outPts;
vtkCellArray *outPolys;
vtkPointData *inPD, *outPD;
vtkCellData *inCD, *outCD;
const vtkTypeInt32 *wholeExt32;
vtkIdType wholeExt[6];
vtkIdType pInc[3];
vtkIdType qInc[3];
vtkIdType cOutInc;
double pt[3];
vtkIdType inStartPtId;
vtkIdType inStartCellId;
vtkIdType outStartPtId;
vtkIdType outPtId;
vtkIdType inId, outId;
vtkIdType ib, ic;
int aA2, bA2, cA2;
outPts = output->GetPoints();
outPD = output->GetPointData();
inPD = input->GetPointData();
outCD = output->GetCellData();
inCD = input->GetCellData();
wholeExt32 = inInfo->Get(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT());
wholeExt[0] = wholeExt32[0]; wholeExt[1] = wholeExt32[1]; wholeExt[2] = wholeExt32[2];
wholeExt[3] = wholeExt32[3]; wholeExt[4] = wholeExt32[4]; wholeExt[5] = wholeExt32[5];
pInc[0] = 1;
pInc[1] = (ext[1]-ext[0]+1);
pInc[2] = (ext[3]-ext[2]+1) * pInc[1];
// quad increments (cell incraments, but cInc could be confused with c axis).
qInc[0] = 1;
qInc[1] = ext[1]-ext[0];
// The conditions are for when we have one or more degenerate axes (2d or 1d cells).
if (qInc[1] == 0)
{
qInc[1] = 1;
}
qInc[2] = (ext[3]-ext[2]) * qInc[1];
if (qInc[2] == 0)
{
qInc[2] = qInc[1];
}
// Temporary variables to avoid many multiplications.
aA2 = aAxis * 2;
bA2 = bAxis * 2;
cA2 = cAxis * 2;
// We might as well put the test for this face here.
if (ext[bA2] == ext[bA2+1] || ext[cA2] == ext[cA2+1])
{
return;
}
if (maxFlag)
{
if (ext[aA2+1] < wholeExt[aA2+1])
{
return;
}
}
else
{ // min faces have a slightly different condition to avoid coincident faces.
if (ext[aA2] == ext[aA2+1] || ext[aA2] > wholeExt[aA2])
{
return;
}
}
// Assuming no ghost cells ...
inStartPtId = inStartCellId = 0;
// I put this confusing conditional to fix a regression test.
// If we are creating a maximum face, then we indeed have to offset
// the input cell Ids. However, vtkGeometryFilter created a 2d
// image as a max face, but the cells are copied as a min face (no
// offset). Hence maxFlag = 1 and there should be no offset.
if (maxFlag && ext[aA2] < ext[1+aA2])
{
inStartPtId = pInc[aAxis]*(ext[aA2+1]-ext[aA2]);
inStartCellId = qInc[aAxis]*(ext[aA2+1]-ext[aA2]-1);
}
outStartPtId = outPts->GetNumberOfPoints();
// Make the points for this face.
for (ic = ext[cA2]; ic <= ext[cA2+1]; ++ic)
{
for (ib = ext[bA2]; ib <= ext[bA2+1]; ++ib)
{
inId = inStartPtId + (ib-ext[bA2])*pInc[bAxis]
+ (ic-ext[cA2])*pInc[cAxis];
input->GetPoint(inId, pt);
outId = outPts->InsertNextPoint(pt);
// Copy point data.
outPD->CopyData(inPD,inId,outId);
this->RecordOrigPointId(outId, inId);
}
}
// Do the cells.
cOutInc = ext[bA2+1] - ext[bA2] + 1;
outPolys = output->GetPolys();
// Old method for creating quads (needed for cell data.).
for (ic = ext[cA2]; ic < ext[cA2+1]; ++ic)
{
for (ib = ext[bA2]; ib < ext[bA2+1]; ++ib)
{
outPtId = outStartPtId + (ib-ext[bA2]) + (ic-ext[cA2])*cOutInc;
inId = inStartCellId + (ib-ext[bA2])*qInc[bAxis] + (ic-ext[cA2])*qInc[cAxis];
outId = outPolys->InsertNextCell(4);
outPolys->InsertCellPoint(outPtId);
outPolys->InsertCellPoint(outPtId+cOutInc);
outPolys->InsertCellPoint(outPtId+cOutInc+1);
outPolys->InsertCellPoint(outPtId+1);
// Copy cell data.
outCD->CopyData(inCD,inId,outId);
this->RecordOrigCellId(outId, inId);
}
}
}
//----------------------------------------------------------------------------
int vtkDataSetSurfaceFilter::DataSetExecute(vtkDataSet *input,
vtkPolyData *output)
{
vtkIdType cellId, newCellId;
int i, j;
vtkIdType numPts=input->GetNumberOfPoints();
vtkIdType numCells=input->GetNumberOfCells();
vtkGenericCell *cell;
vtkCell *face;
double x[3];
vtkIdList *cellIds;
vtkIdList *pts;
vtkPoints *newPts;
vtkIdType ptId, pt;
int npts;
vtkPointData *pd = input->GetPointData();
vtkCellData *cd = input->GetCellData();
vtkPointData *outputPD = output->GetPointData();
vtkCellData *outputCD = output->GetCellData();
if (numCells == 0)
{
return 1;
}
cellIds = vtkIdList::New();
pts = vtkIdList::New();
vtkDebugMacro(<<"Executing geometry filter");
cell = vtkGenericCell::New();
// Allocate
//
newPts = vtkPoints::New();
newPts->Allocate(numPts,numPts/2);
output->Allocate(4*numCells,numCells/2);
outputPD->CopyGlobalIdsOn();
outputPD->CopyAllocate(pd,numPts,numPts/2);
outputCD->CopyGlobalIdsOn();
outputCD->CopyAllocate(cd,numCells,numCells/2);
// Traverse cells to extract geometry
//
int abort=0;
vtkIdType progressInterval = numCells/20 + 1;
for(cellId=0; cellId < numCells && !abort; cellId++)
{
//Progress and abort method support
if ( !(cellId % progressInterval) )
{
vtkDebugMacro(<<"Process cell #" << cellId);
this->UpdateProgress (static_cast<double>(cellId)/numCells);
abort = this->GetAbortExecute();
}
input->GetCell(cellId,cell);
switch (cell->GetCellDimension())
{
// create new points and then cell
case 0: case 1: case 2:
npts = cell->GetNumberOfPoints();
pts->Reset();
for ( i=0; i < npts; i++)
{
ptId = cell->GetPointId(i);
input->GetPoint(ptId, x);
pt = newPts->InsertNextPoint(x);
outputPD->CopyData(pd,ptId,pt);
this->RecordOrigPointId(pt, ptId);
pts->InsertId(i,pt);
}
newCellId = output->InsertNextCell(cell->GetCellType(), pts);
outputCD->CopyData(cd,cellId,newCellId);
this->RecordOrigCellId(newCellId, cellId);
break;
case 3:
for (j=0; j < cell->GetNumberOfFaces(); j++)
{
face = cell->GetFace(j);
input->GetCellNeighbors(cellId, face->PointIds, cellIds);
if ( cellIds->GetNumberOfIds() <= 0)
{
npts = face->GetNumberOfPoints();
pts->Reset();
for ( i=0; i < npts; i++)
{
ptId = face->GetPointId(i);
input->GetPoint(ptId, x);
pt = newPts->InsertNextPoint(x);
outputPD->CopyData(pd,ptId,pt);
this->RecordOrigPointId(pt, ptId);
pts->InsertId(i,pt);
}
newCellId = output->InsertNextCell(face->GetCellType(), pts);
outputCD->CopyData(cd,cellId,newCellId);
this->RecordOrigCellId(newCellId, cellId);
}
}
break;
} //switch
} //for all cells
vtkDebugMacro(<<"Extracted " << newPts->GetNumberOfPoints() << " points,"
<< output->GetNumberOfCells() << " cells.");
// Update ourselves and release memory
//
cell->Delete();
output->SetPoints(newPts);
newPts->Delete();
//free storage
output->Squeeze();
cellIds->Delete();
pts->Delete();
return 1;
}
//----------------------------------------------------------------------------
int vtkDataSetSurfaceFilter::RequestUpdateExtent(
vtkInformation *vtkNotUsed(request),
vtkInformationVector **inputVector,
vtkInformationVector *outputVector)
{
// get the info objects
vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
vtkInformation *outInfo = outputVector->GetInformationObject(0);
int piece, numPieces, ghostLevels;
piece = outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_PIECE_NUMBER());
numPieces =
outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_PIECES());
ghostLevels =
outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS());
if (numPieces > 1 && this->PieceInvariant)
{
// The special execute for structured data handle boundaries internally.
// PolyData does not need any ghost levels.
vtkDataObject* dobj = inInfo->Get(vtkDataObject::DATA_OBJECT());
if (dobj && !strcmp(dobj->GetClassName(), "vtkUnstructuredGrid"))
{ // Processing does nothing fo ghost levels yet so ...
// Be careful to set output ghost level value one less than default
// when they are implemented. I had trouble with multiple executes.
++ghostLevels;
}
}
inInfo->Set(vtkStreamingDemandDrivenPipeline::UPDATE_PIECE_NUMBER(), piece);
inInfo->Set(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_PIECES(),
numPieces);
inInfo->Set(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS(),
ghostLevels);
inInfo->Set(vtkStreamingDemandDrivenPipeline::EXACT_EXTENT(), 1);
return 1;
}
//----------------------------------------------------------------------------
int vtkDataSetSurfaceFilter::FillInputPortInformation(int, vtkInformation *info)
{
info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkDataSet");
return 1;
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::PrintSelf(ostream& os, vtkIndent indent)
{
this->Superclass::PrintSelf(os,indent);
if (this->UseStrips)
{
os << indent << "UseStripsOn\n";
}
else
{
os << indent << "UseStripsOff\n";
}
os << indent << "PieceInvariant: " << this->PieceInvariant << endl;
os << indent << "PassThroughCellIds: " << (this->PassThroughCellIds ? "On\n" : "Off\n");
os << indent << "PassThroughPointIds: " << (this->PassThroughPointIds ? "On\n" : "Off\n");
}
//========================================================================
// Tris are now degenerate quads so we only need one hash table.
// We might want to change the method names from QuadHash to just Hash.
//----------------------------------------------------------------------------
int vtkDataSetSurfaceFilter::UnstructuredGridExecute(vtkDataSet *dataSetInput,
vtkPolyData *output)
{
vtkCellArray *newVerts;
vtkCellArray *newLines;
vtkCellArray *newPolys;
vtkPoints *newPts;
vtkIdType *ids;
int progressCount;
vtkIdType cellId;
int i, j;
vtkIdType *cellPointer;
int cellType;
vtkUnstructuredGrid *input = vtkUnstructuredGrid::SafeDownCast(dataSetInput);
vtkIdType numPts=input->GetNumberOfPoints();
vtkIdType numCells=input->GetNumberOfCells();
vtkGenericCell *cell;
int numFacePts, numCellPts;
vtkIdType inPtId, outPtId;
vtkPointData *inputPD = input->GetPointData();
vtkCellData *inputCD = input->GetCellData();
vtkCellData *cd = input->GetCellData();
vtkPointData *outputPD = output->GetPointData();
vtkCellData *outputCD = output->GetCellData();
vtkIdType outPts[6];
vtkFastGeomQuad *q;
unsigned char* cellTypes = input->GetCellTypesArray()->GetPointer(0);
// These are for the default case/
vtkIdList *pts;
vtkPoints *coords;
vtkCell *face;
int flag2D = 0;
pts = vtkIdList::New();
coords = vtkPoints::New();
cell = vtkGenericCell::New();
this->NumberOfNewCells = 0;
this->InitializeQuadHash(numPts);
// Allocate
//
newPts = vtkPoints::New();
newPts->Allocate(numPts);
newPolys = vtkCellArray::New();
newPolys->Allocate(4*numCells,numCells/2);
newVerts = vtkCellArray::New();
newLines = vtkCellArray::New();
outputPD->CopyGlobalIdsOn();
outputPD->CopyAllocate(inputPD, numPts, numPts/2);
outputCD->CopyGlobalIdsOn();
outputCD->CopyAllocate(inputCD, numCells, numCells/2);
if (this->PassThroughCellIds)
{
this->OriginalCellIds = vtkIdTypeArray::New();
this->OriginalCellIds->SetName("vtkOriginalCellIds");
this->OriginalCellIds->SetNumberOfComponents(1);
}
if (this->PassThroughPointIds)
{
this->OriginalPointIds = vtkIdTypeArray::New();
this->OriginalPointIds->SetName("vtkOriginalPointIds");
this->OriginalPointIds->SetNumberOfComponents(1);
}
// First insert all points. Points have to come first in poly data.
cellPointer = input->GetCells()->GetPointer();
for(cellId=0; cellId < numCells; cellId++)
{
// Direct access to cells.
cellType = cellTypes[cellId];
numCellPts = cellPointer[0];
ids = cellPointer+1;
// Move to the next cell.
cellPointer += (1 + *cellPointer);
// A couple of common cases to see if things go faster.
if (cellType == VTK_VERTEX || cellType == VTK_POLY_VERTEX)
{
newVerts->InsertNextCell(numCellPts);
for (i = 0; i < numCellPts; ++i)
{
inPtId = ids[i];
outPtId = this->GetOutputPointId(inPtId, input, newPts, outputPD);
newVerts->InsertCellPoint(outPtId);
}
this->RecordOrigCellId(this->NumberOfNewCells, cellId);
outputCD->CopyData(cd, cellId, this->NumberOfNewCells++);
}
}
// Traverse cells to extract geometry
//
progressCount = 0;
int abort=0;
vtkIdType progressInterval = numCells/20 + 1;
// First insert all points lines in output and 3D geometry in hash.
// Save 2D geometry for second pass.
// initialize the pointer to the cells for fast traversal.
cellPointer = input->GetCells()->GetPointer();
for(cellId=0; cellId < numCells && !abort; cellId++)
{
//Progress and abort method support
if ( progressCount >= progressInterval )
{
vtkDebugMacro(<<"Process cell #" << cellId);
this->UpdateProgress (static_cast<double>(cellId)/numCells);
abort = this->GetAbortExecute();
progressCount = 0;
}
progressCount++;
// Direct access to cells.
cellType = cellTypes[cellId];
numCellPts = cellPointer[0];
ids = cellPointer+1;
// Move to the next cell.
cellPointer += (1 + *cellPointer);
// A couple of common cases to see if things go faster.
if (cellType == VTK_VERTEX || cellType == VTK_POLY_VERTEX)
{
// Do nothing. This case was handled in the previous loop.
}
else if (cellType == VTK_LINE || cellType == VTK_POLY_LINE)
{
newLines->InsertNextCell(numCellPts);
for (i = 0; i < numCellPts; ++i)
{
inPtId = ids[i];
outPtId = this->GetOutputPointId(inPtId, input, newPts, outputPD);
newLines->InsertCellPoint(outPtId);
}
this->RecordOrigCellId(this->NumberOfNewCells, cellId);
outputCD->CopyData(cd, cellId, this->NumberOfNewCells++);
}
else if (cellType == VTK_HEXAHEDRON)
{
this->InsertQuadInHash(ids[0], ids[1], ids[5], ids[4], cellId);
this->InsertQuadInHash(ids[0], ids[3], ids[2], ids[1], cellId);
this->InsertQuadInHash(ids[0], ids[4], ids[7], ids[3], cellId);
this->InsertQuadInHash(ids[1], ids[2], ids[6], ids[5], cellId);
this->InsertQuadInHash(ids[2], ids[3], ids[7], ids[6], cellId);
this->InsertQuadInHash(ids[4], ids[5], ids[6], ids[7], cellId);
}
else if (cellType == VTK_VOXEL)
{
this->InsertQuadInHash(ids[0], ids[1], ids[5], ids[4], cellId);
this->InsertQuadInHash(ids[0], ids[2], ids[3], ids[1], cellId);
this->InsertQuadInHash(ids[0], ids[4], ids[6], ids[2], cellId);
this->InsertQuadInHash(ids[1], ids[3], ids[7], ids[5], cellId);
this->InsertQuadInHash(ids[2], ids[6], ids[7], ids[3], cellId);
this->InsertQuadInHash(ids[4], ids[5], ids[7], ids[6], cellId);
}
else if (cellType == VTK_TETRA)
{
this->InsertTriInHash(ids[0], ids[1], ids[3], cellId);
this->InsertTriInHash(ids[0], ids[2], ids[1], cellId);
this->InsertTriInHash(ids[0], ids[3], ids[2], cellId);
this->InsertTriInHash(ids[1], ids[2], ids[3], cellId);
}
else if (cellType == VTK_PENTAGONAL_PRISM)
{
// The quads :
this->InsertQuadInHash (ids[0], ids[1], ids[6], ids[5], cellId);
this->InsertQuadInHash (ids[1], ids[2], ids[7], ids[6], cellId);
this->InsertQuadInHash (ids[2], ids[3], ids[8], ids[7], cellId);
this->InsertQuadInHash (ids[3], ids[4], ids[9], ids[8], cellId);
this->InsertQuadInHash (ids[4], ids[0], ids[5], ids[9], cellId);
this->InsertPentaInHash(ids[0], ids[1], ids[2], ids[3], ids[4], cellId);
this->InsertPentaInHash(ids[5], ids[6], ids[7], ids[8], ids[9], cellId);
}
else if (cellType == VTK_HEXAGONAL_PRISM)
{
// The quads :
this->InsertQuadInHash(ids[0], ids[1], ids[7], ids[6], cellId);
this->InsertQuadInHash(ids[1], ids[2], ids[8], ids[7], cellId);
this->InsertQuadInHash(ids[2], ids[3], ids[9], ids[8], cellId);
this->InsertQuadInHash(ids[3], ids[4], ids[10], ids[9], cellId);
this->InsertQuadInHash(ids[4], ids[5], ids[11], ids[10], cellId);
this->InsertQuadInHash(ids[5], ids[0], ids[6], ids[11], cellId);
this->InsertHexInHash (ids[0], ids[1], ids[2], ids[3], ids[4], ids[5], cellId);
this->InsertHexInHash (ids[6], ids[7], ids[8], ids[9], ids[10], ids[11], cellId);
}
else if (cellType == VTK_PIXEL || cellType == VTK_QUAD ||
cellType == VTK_TRIANGLE || cellType == VTK_POLYGON ||
cellType == VTK_TRIANGLE_STRIP ||
cellType == VTK_QUADRATIC_TRIANGLE ||
cellType == VTK_QUADRATIC_QUAD ||
cellType == VTK_QUADRATIC_LINEAR_QUAD ||
cellType == VTK_BIQUADRATIC_QUAD )
{ // save 2D cells for second pass
flag2D = 1;
}
else
// Default way of getting faces. Differentiates between linear
// and higher order cells.
{
input->GetCell(cellId,cell);
if ( cell->IsLinear() )
{
if (cell->GetCellDimension() == 3)
{
int numFaces = cell->GetNumberOfFaces();
for (j=0; j < numFaces; j++)
{
face = cell->GetFace(j);
numFacePts = face->GetNumberOfPoints();
if (numFacePts == 4)
{
this->InsertQuadInHash(face->PointIds->GetId(0),
face->PointIds->GetId(1),
face->PointIds->GetId(2),
face->PointIds->GetId(3), cellId);
}
else if (numFacePts == 3)
{
this->InsertTriInHash(face->PointIds->GetId(0),
face->PointIds->GetId(1),
face->PointIds->GetId(2), cellId);
}
else
{
vtkWarningMacro(<< "I cannot deal with faces with " << numFacePts
<< " points.");
}
} // for all cell faces
} // if 3D
else
{
vtkDebugMacro("Missing cell type.");
}
} // a linear cell type
else //process nonlinear cells via triangulation
{
if ( cell->GetCellDimension() == 1 )
{
cell->Triangulate(0,pts,coords);
for (i=0; i < pts->GetNumberOfIds(); i+=2)
{
newLines->InsertNextCell(2);
inPtId = pts->GetId(i);
this->RecordOrigCellId(this->NumberOfNewCells, cellId);
outputCD->CopyData( cd, cellId, this->NumberOfNewCells++ );
outPtId = this->GetOutputPointId(inPtId, input, newPts, outputPD);
newLines->InsertCellPoint(outPtId);
inPtId = pts->GetId(i+1);
outPtId = this->GetOutputPointId(inPtId, input, newPts, outputPD);
newLines->InsertCellPoint(outPtId);
}
}
else if ( cell->GetCellDimension() == 2 )
{
vtkWarningMacro(<< "2-D nonlinear cells must be processed with all other 2-D cells.");
}
else //3D nonlinear cell
{
vtkIdList *cellIds = vtkIdList::New();
int numFaces = cell->GetNumberOfFaces();
for (j=0; j < numFaces; j++)
{
face = cell->GetFace(j);
input->GetCellNeighbors(cellId, face->PointIds, cellIds);
if ( cellIds->GetNumberOfIds() <= 0)
{
// FIXME: Face could not be consistent. vtkOrderedTriangulator is a better option
face->Triangulate(0,pts,coords);
for (i=0; i < pts->GetNumberOfIds(); i+=3)
{
this->InsertTriInHash(pts->GetId(i), pts->GetId(i+1),
pts->GetId(i+2), cellId);
}
}
}
cellIds->Delete();
} //3d cell
} //nonlinear cell
} // Cell type else.
} // for all cells.
// It would be possible to add these (except for polygons with 5+ sides)
// to the hashes. Alternatively, the higher order 2d cells could be handled
// in the following loop.
// Now insert 2DCells. Because of poly datas (cell data) ordering,
// the 2D cells have to come after points and lines.
// initialize the pointer to the cells for fast traversal.
cellPointer = input->GetCells()->GetPointer();
for(cellId=0; cellId < numCells && !abort && flag2D; cellId++)
{
// Direct acces to cells.
cellType = input->GetCellType(cellId);
numCellPts = cellPointer[0];
ids = cellPointer+1;
// Move to the next cell.
cellPointer += (1 + *cellPointer);
// A couple of common cases to see if things go faster.
if (cellType == VTK_PIXEL)
{ // Do we really want to insert the 2D cells into a hash?
pts->Reset();
pts->InsertId(0, this->GetOutputPointId(ids[0], input, newPts, outputPD));
pts->InsertId(1, this->GetOutputPointId(ids[1], input, newPts, outputPD));
pts->InsertId(2, this->GetOutputPointId(ids[3], input, newPts, outputPD));
pts->InsertId(3, this->GetOutputPointId(ids[2], input, newPts, outputPD));
newPolys->InsertNextCell(pts);
this->RecordOrigCellId(this->NumberOfNewCells, cellId);
outputCD->CopyData(cd, cellId, this->NumberOfNewCells++);
}
else if (cellType == VTK_POLYGON || cellType == VTK_TRIANGLE || cellType == VTK_QUAD)
{
pts->Reset();
for ( i=0; i < numCellPts; i++)
{
inPtId = ids[i];
outPtId = this->GetOutputPointId(inPtId, input, newPts, outputPD);
pts->InsertId(i, outPtId);
}
newPolys->InsertNextCell(pts);
this->RecordOrigCellId(this->NumberOfNewCells, cellId);
outputCD->CopyData(cd, cellId, this->NumberOfNewCells++);
}
else if (cellType == VTK_TRIANGLE_STRIP)
{
// Change strips to triangles so we do not have to worry about order.
int toggle = 0;
vtkIdType ptIds[3];
// This check is not really necessary. It was put here because of another (now fixed) bug.
if (numCellPts > 1)
{
ptIds[0] = this->GetOutputPointId(ids[0], input, newPts, outputPD);
ptIds[1] = this->GetOutputPointId(ids[1], input, newPts, outputPD);
for (i = 2; i < numCellPts; ++i)
{
ptIds[2] = this->GetOutputPointId(ids[i], input, newPts, outputPD);
newPolys->InsertNextCell(3, ptIds);
this->RecordOrigCellId(this->NumberOfNewCells, cellId);
outputCD->CopyData(cd, cellId, this->NumberOfNewCells++);
ptIds[toggle] = ptIds[2];
toggle = !toggle;
}
}
}
else if ( cellType == VTK_QUADRATIC_TRIANGLE
|| cellType == VTK_QUADRATIC_QUAD
|| cellType == VTK_BIQUADRATIC_QUAD
|| cellType == VTK_QUADRATIC_LINEAR_QUAD)
{
input->GetCell( cellId, cell );
cell->Triangulate( 0, pts, coords );
for ( i=0; i < pts->GetNumberOfIds(); i+=3 )
{
outPts[0] = this->GetOutputPointId( pts->GetId(i), input, newPts, outputPD );
outPts[1] = this->GetOutputPointId( pts->GetId(i+1), input, newPts, outputPD );
outPts[2] = this->GetOutputPointId( pts->GetId(i+2), input, newPts, outputPD );
newPolys->InsertNextCell( 3, outPts );
this->RecordOrigCellId(this->NumberOfNewCells, cellId);
outputCD->CopyData(cd, cellId, this->NumberOfNewCells++);
}
}
} // for all cells.
// Now transfer geometry from hash to output (only triangles and quads).
this->InitQuadHashTraversal();
while ( (q = this->GetNextVisibleQuadFromHash()) )
{
// If Tri or Quad:
if( q->p4 == -1 && q->p5 == -1)
{
if ( q->p0 == q->p3)
{
outPts[0] = this->GetOutputPointId(q->p0, input, newPts, outputPD);
outPts[1] = this->GetOutputPointId(q->p1, input, newPts, outputPD);
outPts[2] = this->GetOutputPointId(q->p2, input, newPts, outputPD);
newPolys->InsertNextCell(3, outPts);
this->RecordOrigCellId(this->NumberOfNewCells, q->SourceId);
outputCD->CopyData(inputCD, q->SourceId, this->NumberOfNewCells++);
}
else
{
outPts[0] = this->GetOutputPointId(q->p0, input, newPts, outputPD);
outPts[1] = this->GetOutputPointId(q->p1, input, newPts, outputPD);
outPts[2] = this->GetOutputPointId(q->p2, input, newPts, outputPD);
outPts[3] = this->GetOutputPointId(q->p3, input, newPts, outputPD);
newPolys->InsertNextCell(4, outPts);
this->RecordOrigCellId(this->NumberOfNewCells, q->SourceId);
outputCD->CopyData(inputCD, q->SourceId, this->NumberOfNewCells++);
}
}
else // Penta or Hex
{
if( q->p4 != -1 && q->p5 == -1)
{
outPts[0] = this->GetOutputPointId(q->p0, input, newPts, outputPD);
outPts[1] = this->GetOutputPointId(q->p1, input, newPts, outputPD);
outPts[2] = this->GetOutputPointId(q->p2, input, newPts, outputPD);
outPts[3] = this->GetOutputPointId(q->p3, input, newPts, outputPD);
outPts[4] = this->GetOutputPointId(q->p4, input, newPts, outputPD);
newPolys->InsertNextCell(5, outPts);
this->RecordOrigCellId(this->NumberOfNewCells, q->SourceId);
outputCD->CopyData(inputCD, q->SourceId, this->NumberOfNewCells++);
}
else
{
outPts[0] = this->GetOutputPointId(q->p0, input, newPts, outputPD);
outPts[1] = this->GetOutputPointId(q->p1, input, newPts, outputPD);
outPts[2] = this->GetOutputPointId(q->p2, input, newPts, outputPD);
outPts[3] = this->GetOutputPointId(q->p3, input, newPts, outputPD);
outPts[4] = this->GetOutputPointId(q->p4, input, newPts, outputPD);
outPts[5] = this->GetOutputPointId(q->p5, input, newPts, outputPD);
newPolys->InsertNextCell(6, outPts);
this->RecordOrigCellId(this->NumberOfNewCells, q->SourceId);
outputCD->CopyData(inputCD, q->SourceId, this->NumberOfNewCells++);
}
}
}
if (this->PassThroughCellIds)
{
outputCD->AddArray(this->OriginalCellIds);
}
if (this->PassThroughPointIds)
{
outputPD->AddArray(this->OriginalPointIds);
}
// Update ourselves and release memory
//
cell->Delete();
coords->Delete();
pts->Delete();
output->SetPoints(newPts);
newPts->Delete();
output->SetPolys(newPolys);
newPolys->Delete();
if (newVerts->GetNumberOfCells() > 0)
{
output->SetVerts(newVerts);
}
newVerts->Delete();
newVerts = NULL;
if (newLines->GetNumberOfCells() > 0)
{
output->SetLines(newLines);
}
newLines->Delete();
//free storage
output->Squeeze();
if (this->OriginalCellIds != NULL)
{
this->OriginalCellIds->Delete();
this->OriginalCellIds = NULL;
}
if (this->OriginalPointIds != NULL)
{
this->OriginalPointIds->Delete();
this->OriginalPointIds = NULL;
}
int ghostLevels = output->GetUpdateGhostLevel();
if (this->PieceInvariant)
{
output->RemoveGhostCells(ghostLevels+1);
}
this->DeleteQuadHash();
return 1;
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::InitializeQuadHash(vtkIdType numPoints)
{
vtkIdType i;
if (this->QuadHash)
{
this->DeleteQuadHash();
}
// Prepare our special quad allocator (for efficiency).
this->InitFastGeomQuadAllocation(numPoints);
this->QuadHash = new vtkFastGeomQuad*[numPoints];
this->QuadHashLength = numPoints;
this->PointMap = new vtkIdType[numPoints];
for (i = 0; i < numPoints; ++i)
{
this->QuadHash[i] = NULL;
this->PointMap[i] = -1;
}
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::DeleteQuadHash()
{
vtkIdType i;
this->DeleteAllFastGeomQuads();
for (i = 0; i < this->QuadHashLength; ++i)
{
this->QuadHash[i] = NULL;
}
delete [] this->QuadHash;
this->QuadHash = NULL;
this->QuadHashLength = 0;
delete [] this->PointMap;
this->PointMap = NULL;
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::InsertQuadInHash(vtkIdType a, vtkIdType b,
vtkIdType c, vtkIdType d,
vtkIdType sourceId)
{
int tmp;
vtkFastGeomQuad *quad, **end;
// Reorder to get smallest id in a.
if (b < a && b < c && b < d)
{
tmp = a;
a = b;
b = c;
c = d;
d = tmp;
}
else if (c < a && c < b && c < d)
{
tmp = a;
a = c;
c = tmp;
tmp = b;
b = d;
d = tmp;
}
else if (d < a && d < b && d < c)
{
tmp = a;
a = d;
d = c;
c = b;
b = tmp;
}
// Look for existing quad in the hash;
end = this->QuadHash + a;
quad = *end;
while (quad)
{
end = &(quad->Next);
// a has to match in this bin.
// c should be independant of point order.
if (c == quad->p2)
{
// Check boh orders for b and d.
if ((b == quad->p1 && d == quad->p3) || (b == quad->p3 && d == quad->p1))
{
// We have a match.
quad->SourceId = -1;
// That is all we need to do. Hide any quad shared by two or more cells.
return;
}
}
quad = *end;
}
// Create a new quad and add it to the hash.
quad = this->NewFastGeomQuad();
quad->Next = NULL;
quad->SourceId = sourceId;
quad->p0 = a;
quad->p1 = b;
quad->p2 = c;
quad->p3 = d;
quad->p4 = -1;
quad->p5 = -1;
*end = quad;
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::InsertTriInHash(vtkIdType a, vtkIdType b,
vtkIdType c, vtkIdType sourceId)
{
int tmp;
vtkFastGeomQuad *quad, **end;
// Reorder to get smallest id in a.
if (b < a && b < c)
{
tmp = a;
a = b;
b = c;
c = tmp;
}
else if (c < a && c < b)
{
tmp = a;
a = c;
c = b;
b = tmp;
}
// We can't put the second smnallest in b because it might change the order
// of the verticies in the final triangle.
// Look for existing tri in the hash;
end = this->QuadHash + a;
quad = *end;
while (quad)
{
end = &(quad->Next);
// a has to match in this bin.
if (quad->p0 == quad->p3)
{
if ((b == quad->p1 && c == quad->p2) || (b == quad->p2 && c == quad->p1))
{
// We have a match.
quad->SourceId = -1;
// That is all we need to do. Hide any tri shared by two or more cells.
return;
}
}
quad = *end;
}
// Create a new quad and add it to the hash.
quad = this->NewFastGeomQuad();
quad->Next = NULL;
quad->SourceId = sourceId;
quad->p0 = a;
quad->p1 = b;
quad->p2 = c;
quad->p3 = a;
quad->p4 = -1;
quad->p5 = -1;
*end = quad;
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::InsertPentaInHash(vtkIdType a, vtkIdType b,
vtkIdType c, vtkIdType d,
vtkIdType e, vtkIdType sourceId)
{
int idxmin = 0;
vtkIdType min = a;
vtkFastGeomQuad *quad, **end;
// Copy in an array
vtkIdType tab[5] = {a, b, c, d, e};
// Reorder to get smallest id in a.
for(int i=0;i<5;i++)
{
if( tab[i] < min )
{
min = tab[i];
idxmin = i;
}
}
a = tab[idxmin];
b = tab[(idxmin+1)%5];
c = tab[(idxmin+2)%5];
d = tab[(idxmin+3)%5];
e = tab[(idxmin+4)%5];
// Look for existing hex in the hash;
end = this->QuadHash + a;
quad = *end;
while (quad)
{
end = &(quad->Next);
// a has to match in this bin.
if (quad->p5 == -1)
{
if ((b == quad->p1 && c == quad->p2 && d == quad->p3 && e == quad->p4)
|| (b == quad->p4 && c == quad->p3 && d == quad->p2 && e == quad->p1))
{
// We have a match.
quad->SourceId = -1;
// That is all we need to do. Hide any tri shared by two or more cells.
return;
}
}
quad = *end;
}
// Create a new quad and add it to the hash.
quad = this->NewFastGeomQuad();
quad->Next = NULL;
quad->SourceId = sourceId;
quad->p0 = a;
quad->p1 = b;
quad->p2 = c;
quad->p3 = d;
quad->p4 = e;
quad->p5 = -1;
*end = quad;
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::InsertHexInHash(vtkIdType a, vtkIdType b,
vtkIdType c, vtkIdType d,
vtkIdType e, vtkIdType f,
vtkIdType sourceId)
{
int idxmin = 0;
vtkIdType min = a;
vtkFastGeomQuad *quad, **end;
// Copy in an array
vtkIdType tab[6] = {a, b, c, d, e, f};
// Reorder to get smallest id in a.
for(int i=0;i<6;i++)
{
if( tab[i] < min )
{
min = tab[i];
idxmin = i;
}
}
a = tab[idxmin];
b = tab[(idxmin+1)%6];
c = tab[(idxmin+2)%6];
d = tab[(idxmin+3)%6];
e = tab[(idxmin+4)%6];
f = tab[(idxmin+5)%6];
// Look for existing hex in the hash;
end = this->QuadHash + a;
quad = *end;
while (quad)
{
end = &(quad->Next);
// a has to match in this bin.
// d should be independant of point order.
if (d == quad->p3)
{
if ((b == quad->p1 && c == quad->p2 && e == quad->p4 && f == quad->p5)
|| (b == quad->p5 && c == quad->p4 && e == quad->p2 && f == quad->p1))
{
// We have a match.
quad->SourceId = -1;
// That is all we need to do. Hide any tri shared by two or more cells.
return;
}
}
quad = *end;
}
// Create a new quad and add it to the hash.
quad = this->NewFastGeomQuad();
quad->Next = NULL;
quad->SourceId = sourceId;
quad->p0 = a;
quad->p1 = b;
quad->p2 = c;
quad->p3 = d;
quad->p4 = e;
quad->p5 = f;
*end = quad;
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::InitFastGeomQuadAllocation(int numberOfCells)
{
int idx;
this->DeleteAllFastGeomQuads();
// Allocate 100 pointers to arrays.
// This should be plenty (unless we have riangle strips) ...
this->NumberOfFastGeomQuadArrays = 100;
this->FastGeomQuadArrays = new vtkFastGeomQuad*[100];
// Initalize all to NULL;
for (idx = 0; idx < 100; ++idx)
{
this->FastGeomQuadArrays[idx] = NULL;
}
// Set pointer to the begining.
this->NextArrayIndex = 0;
this->NextQuadIndex = 0;
// Lets keep the chunk size relatively small.
if (numberOfCells < 100)
{
this->FastGeomQuadArrayLength = 50;
}
else
{
this->FastGeomQuadArrayLength = numberOfCells / 2;
}
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::DeleteAllFastGeomQuads()
{
int idx;
for (idx = 0; idx < this->NumberOfFastGeomQuadArrays; ++idx)
{
if (this->FastGeomQuadArrays[idx])
{
delete [] this->FastGeomQuadArrays[idx];
this->FastGeomQuadArrays[idx] = NULL;
}
}
if (this->FastGeomQuadArrays)
{
delete [] this->FastGeomQuadArrays;
this->FastGeomQuadArrays = NULL;
}
this->FastGeomQuadArrayLength = 0;
this->NumberOfFastGeomQuadArrays = 0;
this->NextArrayIndex = 0;
this->NextQuadIndex = 0;
}
//----------------------------------------------------------------------------
vtkFastGeomQuad* vtkDataSetSurfaceFilter::NewFastGeomQuad()
{
if (this->FastGeomQuadArrayLength == 0)
{
vtkErrorMacro("Face hash allocation has not been initialized.");
return NULL;
}
// Although this should not happen often, check first.
if (this->NextArrayIndex >= this->NumberOfFastGeomQuadArrays)
{
int idx, num;
vtkFastGeomQuad** newArrays;
num = this->FastGeomQuadArrayLength * 2;
newArrays = new vtkFastGeomQuad*[num];
for (idx = 0; idx < num; ++idx)
{
newArrays[idx] = NULL;
if (idx < this->NumberOfFastGeomQuadArrays)
{
newArrays[idx] = this->FastGeomQuadArrays[idx];
}
}
delete [] this->FastGeomQuadArrays;
this->FastGeomQuadArrays = newArrays;
this->FastGeomQuadArrayLength = num;
}
// Next: allocate a new array if necessary.
if (this->FastGeomQuadArrays[this->NextArrayIndex] == NULL)
{
// We are allocating a whole bunch at a time to avoid
// allocation of many small pieces of memory.
this->FastGeomQuadArrays[this->NextArrayIndex]
= new vtkFastGeomQuad[this->FastGeomQuadArrayLength];
}
vtkFastGeomQuad* q = (this->FastGeomQuadArrays[this->NextArrayIndex]) + this->NextQuadIndex;
// Move to the next quad.
++(this->NextQuadIndex);
if (this->NextQuadIndex >= this->FastGeomQuadArrayLength)
{
++(this->NextArrayIndex);
this->NextQuadIndex = 0;
}
return q;
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::InitQuadHashTraversal()
{
this->QuadHashTraversalIndex = 0;
this->QuadHashTraversal = this->QuadHash[0];
}
//----------------------------------------------------------------------------
vtkFastGeomQuad *vtkDataSetSurfaceFilter::GetNextVisibleQuadFromHash()
{
vtkFastGeomQuad *quad;
quad = this->QuadHashTraversal;
// Move till traversal until we have a quad to return.
// Note: the current traversal has not been returned yet.
while (quad == NULL || quad->SourceId == -1)
{
if (quad)
{ // The quad must be hidden. Move to the next.
quad = quad->Next;
}
else
{ // must be the end of the linked list. Move to the next bin.
this->QuadHashTraversalIndex += 1;
if ( this->QuadHashTraversalIndex >= this->QuadHashLength)
{ // There are no more bins.
this->QuadHashTraversal = NULL;
return NULL;
}
quad = this->QuadHash[this->QuadHashTraversalIndex];
}
}
// Now we have a quad to return. Set the traversal to the next entry.
this->QuadHashTraversal = quad->Next;
return quad;
}
//----------------------------------------------------------------------------
vtkIdType vtkDataSetSurfaceFilter::GetOutputPointId(vtkIdType inPtId,
vtkDataSet *input,
vtkPoints *outPts,
vtkPointData *outPD)
{
vtkIdType outPtId;
outPtId = this->PointMap[inPtId];
if (outPtId == -1)
{
outPtId = outPts->InsertNextPoint(input->GetPoint(inPtId));
outPD->CopyData(input->GetPointData(), inPtId, outPtId);
this->PointMap[inPtId] = outPtId;
this->RecordOrigPointId(outPtId, inPtId);
}
return outPtId;
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::RecordOrigCellId(vtkIdType destIndex,
vtkIdType originalId)
{
if (this->OriginalCellIds != NULL)
{
this->OriginalCellIds->InsertValue(destIndex, originalId);
}
}
//----------------------------------------------------------------------------
void vtkDataSetSurfaceFilter::RecordOrigPointId(vtkIdType destIndex,
vtkIdType originalId)
{
if (this->OriginalPointIds != NULL)
{
this->OriginalPointIds->InsertValue(destIndex, originalId);
}
}
["smime.p7s" (application/x-pkcs7-signature)]
_______________________________________________
Powered by www.kitware.com
Visit other Kitware open-source projects at http://www.kitware.com/opensource/opensource.html
Please keep messages on-topic and check the VTK FAQ at: http://www.vtk.org/Wiki/VTK_FAQ
Follow this link to subscribe/unsubscribe:
http://www.vtk.org/mailman/listinfo/vtkusers
--===============0287992244==--
[prev in list] [next in list] [prev in thread] [next in thread]
Configure |
About |
News |
Add a list |
Sponsored by KoreLogic