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.
369 lines
9.9 KiB
369 lines
9.9 KiB
/*=========================================================================
|
|
|
|
Program: Visualization Toolkit
|
|
Module: $RCSfile: vtkLoopSubdivisionFilter.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 "vtkLoopSubdivisionFilter.h"
|
|
|
|
#include "vtkCell.h"
|
|
#include "vtkCellArray.h"
|
|
#include "vtkEdgeTable.h"
|
|
#include "vtkIdList.h"
|
|
#include "vtkInformation.h"
|
|
#include "vtkInformationVector.h"
|
|
#include "vtkObjectFactory.h"
|
|
#include "vtkPointData.h"
|
|
#include "vtkPolyData.h"
|
|
#include "vtkStreamingDemandDrivenPipeline.h"
|
|
|
|
vtkCxxRevisionMacro(vtkLoopSubdivisionFilter, "$Revision: 1.22 $");
|
|
vtkStandardNewMacro(vtkLoopSubdivisionFilter);
|
|
|
|
static double LoopWeights[4] =
|
|
{.375, .375, .125, .125};
|
|
|
|
void vtkLoopSubdivisionFilter::GenerateSubdivisionPoints (vtkPolyData *inputDS,vtkIntArray *edgeData, vtkPoints *outputPts, vtkPointData *outputPD)
|
|
{
|
|
double *weights;
|
|
vtkIdType *pts = 0;
|
|
vtkIdType numPts, cellId, newId;
|
|
int edgeId;
|
|
vtkIdType npts;
|
|
vtkIdType p1, p2;
|
|
vtkCellArray *inputPolys=inputDS->GetPolys();
|
|
vtkEdgeTable *edgeTable;
|
|
vtkIdList *cellIds = vtkIdList::New();
|
|
vtkIdList *stencil = vtkIdList::New();
|
|
vtkPoints *inputPts=inputDS->GetPoints();
|
|
vtkPointData *inputPD=inputDS->GetPointData();
|
|
|
|
weights = new double[256];
|
|
|
|
// Create an edge table to keep track of which edges we've processed
|
|
edgeTable = vtkEdgeTable::New();
|
|
edgeTable->InitEdgeInsertion(inputDS->GetNumberOfPoints());
|
|
|
|
// Generate even points. these are derived from the old points
|
|
numPts = inputDS->GetNumberOfPoints();
|
|
for (vtkIdType ptId=0; ptId < numPts; ptId++)
|
|
{
|
|
this->GenerateEvenStencil (ptId, inputDS, stencil, weights);
|
|
this->InterpolatePosition (inputPts, outputPts, stencil, weights);
|
|
outputPD->InterpolatePoint (inputPD, ptId, stencil, weights);
|
|
}
|
|
|
|
// Generate odd points. These will be inserted into the new dataset
|
|
for (cellId=0, inputPolys->InitTraversal();
|
|
inputPolys->GetNextCell(npts, pts); cellId++)
|
|
{
|
|
if ( inputDS->GetCellType(cellId) != VTK_TRIANGLE )
|
|
{
|
|
continue;
|
|
}
|
|
|
|
// start with one edge
|
|
p1 = pts[2];
|
|
p2 = pts[0];
|
|
|
|
for (edgeId=0; edgeId < 3; edgeId++)
|
|
{
|
|
// Do we need to create a point on this edge?
|
|
if (edgeTable->IsEdge (p1, p2) == -1)
|
|
{
|
|
edgeTable->InsertEdge (p1, p2);
|
|
inputDS->GetCellEdgeNeighbors (-1, p1, p2, cellIds);
|
|
if (cellIds->GetNumberOfIds() == 1)
|
|
{
|
|
// Compute new Position and PointData using the same subdivision scheme
|
|
stencil->SetNumberOfIds(2);
|
|
stencil->SetId(0,p1);
|
|
stencil->SetId(1,p2);
|
|
weights[0] = .5; weights[1] = .5;
|
|
} // boundary edge
|
|
else
|
|
{
|
|
this->GenerateOddStencil (p1, p2,
|
|
inputDS, stencil, weights);
|
|
}
|
|
newId = this->InterpolatePosition (inputPts, outputPts,
|
|
stencil, weights);
|
|
outputPD->InterpolatePoint (inputPD, newId, stencil, weights);
|
|
}
|
|
else // we have already created a point on this edge. find it
|
|
{
|
|
newId = this->FindEdge (inputDS, cellId, p1, p2, edgeData, cellIds);
|
|
}
|
|
edgeData->InsertComponent(cellId,edgeId,newId);
|
|
p1 = p2;
|
|
if (edgeId < 2)
|
|
{
|
|
p2 = pts[edgeId + 1];
|
|
}
|
|
} // each interior edge
|
|
} // each cell
|
|
|
|
// cleanup
|
|
delete [] weights;
|
|
edgeTable->Delete();
|
|
stencil->Delete ();
|
|
cellIds->Delete();
|
|
}
|
|
|
|
void vtkLoopSubdivisionFilter::GenerateEvenStencil (vtkIdType p1,
|
|
vtkPolyData *polys,
|
|
vtkIdList *stencilIds,
|
|
double *weights)
|
|
{
|
|
vtkIdList *cellIds = vtkIdList::New();
|
|
vtkIdList *ptIds = vtkIdList::New();
|
|
vtkCell *cell;
|
|
|
|
int i, j;
|
|
int numCellsInLoop;
|
|
int startCell, nextCell;
|
|
vtkIdType p, p2;
|
|
vtkIdType bp1, bp2;
|
|
int K;
|
|
double beta, cosSQ;
|
|
|
|
// Get the cells that use this point
|
|
polys->GetPointCells (p1, cellIds);
|
|
numCellsInLoop = cellIds->GetNumberOfIds();
|
|
if (numCellsInLoop < 1)
|
|
{
|
|
vtkWarningMacro("numCellsInLoop < 1: " << numCellsInLoop);
|
|
stencilIds->Reset();
|
|
return;
|
|
}
|
|
// Find an edge to start with that contains p1
|
|
polys->GetCellPoints (cellIds->GetId(0), ptIds);
|
|
p2 = ptIds->GetId(0);
|
|
i = 1;
|
|
while (p1 == p2)
|
|
{
|
|
p2 = ptIds->GetId(i++);
|
|
}
|
|
polys->GetCellEdgeNeighbors (-1, p1, p2, cellIds);
|
|
|
|
nextCell = cellIds->GetId(0);
|
|
bp2 = -1;
|
|
bp1 = p2;
|
|
if (cellIds->GetNumberOfIds() == 1)
|
|
{
|
|
startCell = -1;
|
|
}
|
|
else
|
|
{
|
|
startCell = cellIds->GetId(1);
|
|
}
|
|
|
|
stencilIds->Reset();
|
|
stencilIds->InsertNextId(p2);
|
|
|
|
// walk around the loop counter-clockwise and get cells
|
|
for (j = 0; j < numCellsInLoop; j++)
|
|
{
|
|
cell = polys->GetCell(nextCell);
|
|
p = -1;
|
|
for (i = 0; i < 3; i++)
|
|
{
|
|
if ((p = cell->GetPointId(i)) != p1 && cell->GetPointId(i) != p2)
|
|
{
|
|
break;
|
|
}
|
|
}
|
|
p2 = p;
|
|
stencilIds->InsertNextId (p2);
|
|
polys->GetCellEdgeNeighbors (nextCell, p1, p2, cellIds);
|
|
if (cellIds->GetNumberOfIds() != 1)
|
|
{
|
|
bp2 = p2;
|
|
j++;
|
|
break;
|
|
}
|
|
nextCell = cellIds->GetId(0);
|
|
}
|
|
|
|
// now walk around the other way. this will only happen if there
|
|
// is a boundary cell left that we have not visited
|
|
nextCell = startCell;
|
|
p2 = bp1;
|
|
for (; j < numCellsInLoop && startCell != -1; j++)
|
|
{
|
|
cell = polys->GetCell(nextCell);
|
|
p = -1;
|
|
for (i = 0; i < 3; i++)
|
|
{
|
|
if ((p = cell->GetPointId(i)) != p1 && cell->GetPointId(i) != p2)
|
|
{
|
|
break;
|
|
}
|
|
}
|
|
p2 = p;
|
|
stencilIds->InsertNextId (p2);
|
|
polys->GetCellEdgeNeighbors (nextCell, p1, p2, cellIds);
|
|
if (cellIds->GetNumberOfIds() != 1)
|
|
{
|
|
bp1 = p2;
|
|
break;
|
|
}
|
|
nextCell = cellIds->GetId(0);
|
|
}
|
|
|
|
if (bp2 != -1) // boundary edge
|
|
{
|
|
stencilIds->SetNumberOfIds(3);
|
|
stencilIds->SetId(0,bp2);
|
|
stencilIds->SetId(1,bp1);
|
|
stencilIds->SetId(2,p1);
|
|
weights[0] = .125;
|
|
weights[1] = .125;
|
|
weights[2] = .75;
|
|
}
|
|
else
|
|
{
|
|
K = stencilIds->GetNumberOfIds();
|
|
// Remove last id. It's a duplicate of the first
|
|
K--;
|
|
if (K > 3)
|
|
{
|
|
// Generate weights
|
|
#define VTK_PI 3.14159265358979
|
|
cosSQ = .375 + .25 * cos (2.0 * VTK_PI / (double) K);
|
|
cosSQ = cosSQ * cosSQ;
|
|
beta = (.625 - cosSQ) / (double) K;
|
|
}
|
|
else
|
|
{
|
|
beta = 3.0 / 16.0;
|
|
}
|
|
for (j = 0; j < K; j++)
|
|
{
|
|
weights[j] = beta;
|
|
}
|
|
weights[K] = 1.0 - K * beta;
|
|
stencilIds->SetId (K,p1);
|
|
}
|
|
cellIds->Delete();
|
|
ptIds->Delete();
|
|
}
|
|
|
|
void vtkLoopSubdivisionFilter::GenerateOddStencil (vtkIdType p1, vtkIdType p2,
|
|
vtkPolyData *polys,
|
|
vtkIdList *stencilIds,
|
|
double *weights)
|
|
{
|
|
vtkIdList *cellIds = vtkIdList::New();
|
|
vtkCell *cell;
|
|
int i;
|
|
vtkIdType cell0, cell1;
|
|
vtkIdType p3=0, p4=0;
|
|
|
|
polys->GetCellEdgeNeighbors (-1, p1, p2, cellIds);
|
|
cell0 = cellIds->GetId(0);
|
|
cell1 = cellIds->GetId(1);
|
|
|
|
cell = polys->GetCell(cell0);
|
|
for (i = 0; i < 3; i++)
|
|
{
|
|
if ((p3 = cell->GetPointId(i)) != p1 && cell->GetPointId(i) != p2)
|
|
{
|
|
break;
|
|
}
|
|
}
|
|
cell = polys->GetCell(cell1);
|
|
for (i = 0; i < 3; i++)
|
|
{
|
|
if ((p4 = cell->GetPointId(i)) != p1 && cell->GetPointId(i) != p2)
|
|
{
|
|
break;
|
|
}
|
|
}
|
|
|
|
stencilIds->SetNumberOfIds (4);
|
|
stencilIds->SetId(0, p1);
|
|
stencilIds->SetId(1, p2);
|
|
stencilIds->SetId(2, p3);
|
|
stencilIds->SetId(3, p4);
|
|
|
|
for (i = 0; i < stencilIds->GetNumberOfIds (); i++)
|
|
{
|
|
weights[i] = LoopWeights[i];
|
|
}
|
|
cellIds->Delete();
|
|
}
|
|
|
|
int vtkLoopSubdivisionFilter::RequestUpdateExtent(
|
|
vtkInformation *request,
|
|
vtkInformationVector **inputVector,
|
|
vtkInformationVector *outputVector)
|
|
{
|
|
int numPieces, ghostLevel;
|
|
vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
|
|
vtkInformation *outInfo = outputVector->GetInformationObject(0);
|
|
|
|
if (!this->Superclass::RequestUpdateExtent(request, inputVector,
|
|
outputVector))
|
|
{
|
|
return 0;
|
|
}
|
|
|
|
numPieces =
|
|
outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_PIECES());
|
|
ghostLevel =
|
|
outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS());
|
|
|
|
if (numPieces > 1 && this->NumberOfSubdivisions > 0)
|
|
{
|
|
inInfo->Set(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS(),
|
|
ghostLevel + 1);
|
|
}
|
|
|
|
return 1;
|
|
}
|
|
|
|
int vtkLoopSubdivisionFilter::RequestData(
|
|
vtkInformation *request,
|
|
vtkInformationVector **inputVector,
|
|
vtkInformationVector *outputVector)
|
|
{
|
|
vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
|
|
|
|
vtkPolyData *input = vtkPolyData::SafeDownCast(
|
|
inInfo->Get(vtkDataObject::DATA_OBJECT()));
|
|
vtkCellArray *polys = input->GetPolys();
|
|
int hasTris = 0;
|
|
vtkIdType numPts = 0, *pts = 0;
|
|
|
|
input->BuildLinks();
|
|
|
|
polys->InitTraversal();
|
|
while (polys->GetNextCell(numPts, pts))
|
|
{
|
|
if (numPts == 3)
|
|
{
|
|
if (input->IsTriangle(pts[0], pts[1], pts[2]))
|
|
{
|
|
hasTris = 1;
|
|
break;
|
|
}
|
|
}
|
|
}
|
|
|
|
if (!hasTris)
|
|
{
|
|
vtkWarningMacro("vtkLoopSubdivisionFilter only operates on triangles, but this data set has no triangles to operate on.")
|
|
return 0;
|
|
}
|
|
|
|
return this->Superclass::RequestData(request, inputVector, outputVector);
|
|
}
|
|
|