Skip to content

GitLab

Commit 87012e51 authored by Will Schroeder's avatar Will Schroeder
Browse files
Options

First incarnation of Flying Edges

parent 1a965eb9
Hide whitespace changes
Inline Side-by-side
Showing with 55 additions and 177 deletions
Filters/Core/Testing/Python/TestFlyingEdges2D.py
View file @ 87012e51
......@@ -13,88 +13,28 @@ renWin.AddRenderer(ren1)
iren = vtk.vtkRenderWindowInteractor()
iren.SetRenderWindow(renWin)
# Create image manually for testing purposes
xSze = 5
ySze = 12
xSze = 4000
ySze = 6000
isoValue = 750
maxValue = 100
imageData = vtk.vtkImageData()
imageData.SetDimensions(xSze,ySze,1)
scalars = vtk.vtkUnsignedCharArray()
scalars.SetNumberOfTuples(xSze*ySze)
for i in range (0,xSze*ySze):
scalars.SetValue(i,0)
scalars.SetValue(9,maxValue)
scalars.SetValue(17,maxValue)
scalars.SetValue(18,maxValue)
scalars.SetValue(36,maxValue)
scalars.SetValue(42,maxValue)
scalars.SetValue(1000000,maxValue)
imageData.GetPointData().SetScalars(scalars)
# pipeline
reader = vtk.vtkStructuredPointsReader()
reader.SetFileName("C:\D\gitVTK\Bugs\FlyingEdges\StructuredPoints-Test.vtk")
reader2 = vtk.vtkPNGReader()
reader2.SetFileName("" + str(VTK_DATA_ROOT) + "/Data/fullhead15.png")
# Pipeline
reader = vtk.vtkPNGReader()
reader.SetFileName("" + str(VTK_DATA_ROOT) + "/Data/fullhead15.png")
iso = vtk.vtkFlyingEdges2D()
iso.SetInputData(imageData)
#iso.SetInputConnection(reader2.GetOutputPort())
iso.SetValue(0, isoValue)
#iso.GenerateValues(12,500,1150)
iso2 = vtk.vtkSynchronizedTemplates2D()
iso2.SetInputData(imageData)
#iso2.SetInputConnection(reader2.GetOutputPort())
iso2.SetValue(0, isoValue)
#iso2.GenerateValues(12,500,1150)
iso.SetInputConnection(reader.GetOutputPort())
iso.GenerateValues(12,500,1150)
isoMapper = vtk.vtkPolyDataMapper()
isoMapper.SetInputConnection(iso.GetOutputPort())
#isoMapper.SetInputConnection(iso2.GetOutputPort())
isoMapper.ScalarVisibilityOff()
isoActor = vtk.vtkActor()
isoActor.SetMapper(isoMapper)
isoActor.GetProperty().SetColor(1,1,1)
outline = vtk.vtkOutlineFilter()
outline.SetInputData(imageData)
#outline.SetInputConnection(reader2.GetOutputPort())
outline.SetInputConnection(reader.GetOutputPort())
outlineMapper = vtk.vtkPolyDataMapper()
outlineMapper.SetInputConnection(outline.GetOutputPort())
outlineActor = vtk.vtkActor()
outlineActor.SetMapper(outlineMapper)
outlineProp = outlineActor.GetProperty()
#eval $outlineProp SetColor 0 0 0
# Time the execution of the filter
timer = vtk.vtkExecutionTimer()
timer.SetFilter(iso)
iso.Update()
FEwallClock = timer.GetElapsedWallClockTime()
FEcpuClock = timer.GetElapsedCPUTime()
print ("FE:", FEwallClock, FEcpuClock)
timer2 = vtk.vtkExecutionTimer()
timer2.SetFilter(iso2)
iso2.Update()
wallClock = timer2.GetElapsedWallClockTime()
cpuClock = timer2.GetElapsedCPUTime()
print ("ST:", wallClock, cpuClock)
if FEwallClock > 0:
print ("Speedup:", (wallClock/FEwallClock))
#write output
writer = vtk.vtkPolyDataWriter()
writer.SetInputConnection(iso.GetOutputPort())
writer.SetFileName("C:\D\gitVTK\Bugs\FlyingEdges\FE-output.vtk")
#if xSze*ySze <= 100:
writer.Write()
# Add the actors to the renderer, set the background and size
#
......@@ -105,6 +45,5 @@ renWin.SetSize(400,400)
ren1.ResetCamera()
iren.Initialize()
iren.Start()
# prevent the tk window from showing up then start the event loop
renWin.Render()
# --- end of script --
Filters/Core/Testing/Python/TestFlyingEdges3D.py
View file @ 87012e51
......@@ -13,34 +13,7 @@ renWin.AddRenderer(ren1)
iren = vtk.vtkRenderWindowInteractor()
iren.SetRenderWindow(renWin)
# Create image manually for testing purposes
xSze = 6
ySze = 5
zSze = 8
#xSze = 150
#ySze = 175
#zSze = 750
isoValue = 200
maxValue = 250
imageData = vtk.vtkImageData()
imageData.SetDimensions(xSze,ySze,zSze)
scalars = vtk.vtkIntArray()
scalars.SetNumberOfTuples(xSze*ySze*zSze)
for i in range (0,xSze*ySze*zSze):
# scalars.SetValue(i,i)
scalars.SetValue(i,0)
scalars.SetValue(222,maxValue)
imageData.GetPointData().SetScalars(scalars)
# pipeline
reader = vtk.vtkImageReader()
reader.SetDataByteOrderToLittleEndian()
reader.SetDataExtent(0,63,0,63,1,93)
reader.SetDataSpacing(3.2,3.2,1.5)
reader.SetFilePrefix("" + str(VTK_DATA_ROOT) + "/Data/headsq/quarter")
reader.SetDataMask(0x7fff)
# another source
# Create a synthetic source
sphere = vtk.vtkSphere()
sphere.SetCenter( 0.0,0.0,0.0)
sphere.SetRadius(0.25)
......@@ -49,33 +22,16 @@ sphere.SetRadius(0.25)
sample = vtk.vtkSampleFunction()
sample.SetImplicitFunction(sphere)
sample.SetModelBounds(-0.5,0.5, -0.5,0.5, -0.5,0.5)
sample.SetSampleDimensions(200,200,200)
sample.SetSampleDimensions(250,350,400)
iso = vtk.vtkFlyingEdges3D()
#iso.SetInputData(imageData)
#iso.SetInputConnection(reader.GetOutputPort())
iso.SetInputConnection(sample.GetOutputPort())
#iso.SetValue(0, isoValue)
#iso.SetValue(0, 750)
#iso.SetValue(0, 0.0)
#iso.GenerateValues(12,500,1150)
iso.GenerateValues(3,-.11,.11)
iso.ComputeNormalsOn()
iso.ComputeGradientsOff()
iso2 = vtk.vtkSynchronizedTemplates3D()
#iso2.SetInputData(imageData)
#iso2.SetInputConnection(reader.GetOutputPort())
iso2.SetInputConnection(sample.GetOutputPort())
#iso2.SetValue(0, isoValue)
#iso2.SetValue(0, 750)
#iso2.SetValue(0, 0.0)
#iso2.GenerateValues(12,500,1150)
iso2.GenerateValues(3,-.11,.11)
iso.ComputeGradientsOn()
isoMapper = vtk.vtkPolyDataMapper()
isoMapper.SetInputConnection(iso.GetOutputPort())
#isoMapper.SetInputConnection(iso2.GetOutputPort())
isoMapper.ScalarVisibilityOff()
isoActor = vtk.vtkActor()
isoActor.SetMapper(isoMapper)
......@@ -83,40 +39,12 @@ isoActor.GetProperty().SetColor(1,1,1)
isoActor.GetProperty().SetOpacity(0.5)
outline = vtk.vtkOutlineFilter()
#outline.SetInputData(imageData)
#outline.SetInputConnection(reader.GetOutputPort())
outline.SetInputConnection(sample.GetOutputPort())
outlineMapper = vtk.vtkPolyDataMapper()
outlineMapper.SetInputConnection(outline.GetOutputPort())
outlineActor = vtk.vtkActor()
outlineActor.SetMapper(outlineMapper)
outlineProp = outlineActor.GetProperty()
#eval $outlineProp SetColor 0 0 0
# Time the execution of the filter
timer = vtk.vtkExecutionTimer()
timer.SetFilter(iso)
iso.Update()
FEwallClock = timer.GetElapsedWallClockTime()
FEcpuClock = timer.GetElapsedCPUTime()
print ("FE:", FEwallClock, FEcpuClock)
timer2 = vtk.vtkExecutionTimer()
timer2.SetFilter(iso2)
iso2.Update()
wallClock = timer2.GetElapsedWallClockTime()
cpuClock = timer2.GetElapsedCPUTime()
print ("ST:", wallClock, cpuClock)
if FEwallClock > 0:
print ("Speedup:", (wallClock/FEwallClock))
#write output
writer = vtk.vtkPolyDataWriter()
writer.SetInputConnection(iso.GetOutputPort())
writer.SetFileName("C:\D\gitVTK\Bugs\FlyingEdges\FE-output.vtk")
#if xSze*ySze <= 100:
#writer.Write()
# Add the actors to the renderer, set the background and size
#
......@@ -127,6 +55,5 @@ renWin.SetSize(400,400)
ren1.ResetCamera()
iren.Initialize()
iren.Start()
# prevent the tk window from showing up then start the event loop
renWin.Render()
# --- end of script --
Filters/Core/vtkFlyingEdges2D.h
View file @ 87012e51
......@@ -14,8 +14,8 @@
=========================================================================*/
// .NAME vtkFlyingEdges2D - generate isoline(s) from a structured points set
// .SECTION Description
// vtkFlyingEdges2D is a serial, reference implementation of the 2D version
// of the flying edges algorithm. It is designed to be highly scalable (i.e.,
// vtkFlyingEdges2D is a reference implementation of the 2D version of the
// flying edges algorithm. It is designed to be highly scalable (i.e.,
// parallelizable) for large data. It implements certain performance
// optimizations including computational trimming to rapidly eliminate
// processing of data regions, packed bit representation of case table
......@@ -127,7 +127,8 @@ protected:
vtkFlyingEdges2D();
~vtkFlyingEdges2D();
virtual int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *);
virtual int RequestData(vtkInformation *, vtkInformationVector **,
vtkInformationVector *);
virtual int FillInputPortInformation(int port, vtkInformation *info);
vtkContourValues *ContourValues;
......
Filters/Core/vtkFlyingEdges3D.cxx
View file @ 87012e51
......@@ -171,9 +171,9 @@ public:
{
if ( loc == Interior )
{
g[0] = ( *(s+this->Inc0) - *(s-this->Inc0) / this->Spacing[0] );
g[1] = ( *(s+this->Inc1) - *(s-this->Inc1) / this->Spacing[1] );
g[2] = ( *(s+this->Inc2) - *(s-this->Inc2) / this->Spacing[2] );
g[0] = 0.5*( (*(s+this->Inc0) - *(s-this->Inc0)) / this->Spacing[0] );
g[1] = 0.5*( (*(s+this->Inc1) - *(s-this->Inc1)) / this->Spacing[1] );
g[2] = 0.5*( (*(s+this->Inc2) - *(s-this->Inc2)) / this->Spacing[2] );
}
else
{
......@@ -346,11 +346,16 @@ vtkFlyingEdges3DAlgorithm():XCases(NULL),EdgeMetaData(NULL),NewScalars(NULL),
{
edgeCase = this->EdgeCases[eCase];
*edgeCase++ = numTris;
for ( edge = triCase->edges; edge[0] > -1; edge += 3 )
{// build new case table
*edgeCase++ = this->EdgeMap[edge[0]];
*edgeCase++ = this->EdgeMap[edge[1]];
*edgeCase++ = this->EdgeMap[edge[2]];
for ( edge = triCase->edges; edge[0] > -1; edge += 3, edgeCase+=3 )
{
// Build new case table. You're probably wondering why the
// crazy (0,2,1) edge order below. Simple: as originally
// presented the MC algorithm used a left-handed coordinate
// system, so we have to reverse the ordering of the triangle
// to make it consistent with any generated normals.
edgeCase[0] = this->EdgeMap[edge[0]];
edgeCase[2] = this->EdgeMap[edge[1]];
edgeCase[1] = this->EdgeMap[edge[2]];
}
}
}//x-edges
......@@ -477,7 +482,7 @@ ComputeBoundaryGradient(vtkIdType ijk[3], T *s, float g[3])
//----------------------------------------------------------------------------
// Interpolate a new point along a boundary edge. Make sure to consider
// proximity to boundary when computing gradients, etc.
// proximity to the boundary when computing gradients, etc.
template <class T> void vtkFlyingEdges3DAlgorithm<T>::
InterpolateEdge(double value, vtkIdType ijk[3], T *s, float x[3],
unsigned char edgeNum, unsigned char edgeUses[12],
......@@ -549,34 +554,34 @@ GeneratePoints(double value, unsigned char loc, vtkIdType ijk[3], T *sPtr,
{
// Create a slightly faster path for voxel axes interior to the volume.
float g0[3], x1[3];
vtkIdType offset[3];
vtkIdType ijk1[3];
if ( this->NeedGradients )
{
this->ComputeGradient(loc,ijk,sPtr,g0);
}
if ( edgeUses[0] ) //x axes edge
{
x1[0] = x[0] + this->Spacing[0]; offset[0] = ijk[0] + 1;
x1[1] = x[1]; offset[1] = 0;
x1[2] = x[2]; offset[2] = 0;
x1[0] = x[0] + this->Spacing[0]; ijk1[0] = ijk[0] + 1;
x1[1] = x[1]; ijk1[1] = ijk[1];
x1[2] = x[2]; ijk1[2] = ijk[2];
this->InterpolateAxesEdge(value, loc, sPtr, x, sPtr+this->Inc0,
x1, eIds[0], offset, g0);
x1, eIds[0], ijk1, g0);
}
if ( edgeUses[4] ) //y axes edge
{
x1[0] = x[0]; offset[0] = 0;
x1[1] = x[1] + this->Spacing[1]; offset[1] = ijk[1] + 1;
x1[2] = x[2]; offset[2] = 0;
x1[0] = x[0]; ijk1[0] = ijk[0];
x1[1] = x[1] + this->Spacing[1]; ijk1[1] = ijk[1] + 1;
x1[2] = x[2]; ijk1[2] = ijk[2];
this->InterpolateAxesEdge(value, loc, sPtr, x, sPtr+this->Inc1,
x1, eIds[4], offset, g0);
x1, eIds[4], ijk1, g0);
}
if ( edgeUses[8] ) //z axes edge
{
x1[0] = x[0]; offset[0] = 0;
x1[1] = x[1]; offset[1] = 0;
x1[2] = x[2] + this->Spacing[2]; offset[2] = ijk[2] + 1;
x1[0] = x[0]; ijk1[0] = ijk[0];
x1[1] = x[1]; ijk1[1] = ijk[1];
x1[2] = x[2] + this->Spacing[2]; ijk1[2] = ijk[2] + 1;
this->InterpolateAxesEdge(value, loc, sPtr, x, sPtr+this->Inc2,
x1, eIds[8], offset, g0);
x1, eIds[8], ijk1, g0);
}
// Otherwise do more general gyrations. These are boundary situations where
......@@ -1191,11 +1196,13 @@ int vtkFlyingEdges3D::RequestData(
{
newNormals = vtkFloatArray::New();
newNormals->SetNumberOfComponents(3);
newNormals->SetName("Gradients");
}
if (this->ComputeGradients)
{
newGradients = vtkFloatArray::New();
newGradients->SetNumberOfComponents(3);
newGradients->SetName("Gradients");
}
void *ptr = input->GetArrayPointerForExtent(inScalars, exExt);
......
Filters/Core/vtkFlyingEdges3D.h
View file @ 87012e51
......@@ -14,8 +14,8 @@
=========================================================================*/
// .NAME vtkFlyingEdges3D - generate isosurface from 3D image data (volume)
// .SECTION Description
// vtkFlyingEdges3D is a serial, reference implementation of the 3D version
// of the flying edges algorithm. It is designed to be highly scalable (i.e.,
// vtkFlyingEdges3D is a reference implementation of the 3D version of the
// flying edges algorithm. It is designed to be highly scalable (i.e.,
// parallelizable) for large data. It implements certain performance
// optimizations including computational trimming to rapidly eliminate
// processing of data regions, packed bit representation of case table
......@@ -28,9 +28,9 @@
// This is a three-pass algorithm. The first pass processes all x-edges and
// builds x-edge case values (which, when the four x-edges defining a voxel
// are combined, are equivalent to vertex-based case table except edge-based
// approaches are separable to parallel computing). Next x-voxel rows are
// processed to gather information from yz-edges (basically to count the
// number of edge intersections and triangles generated). Finally in the
// approaches are separable in support of parallel computing). Next x-voxel
// rows are processed to gather information from yz-edges (basically to count
// the number of edge intersections and triangles generated). Finally in the
// third pass output primitives are generated into pre-allocated arrays. This
// implementation uses voxel cell axes (a x-y-z triad located at the voxel
// origin) to ensure that each edge is intersected at most one time. Note
......@@ -152,8 +152,10 @@ protected:
int ArrayComponent;
vtkContourValues *ContourValues;
virtual int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *);
virtual int RequestUpdateExtent(vtkInformation *, vtkInformationVector **, vtkInformationVector *);
virtual int RequestData(vtkInformation *, vtkInformationVector **,
vtkInformationVector *);
virtual int RequestUpdateExtent(vtkInformation *, vtkInformationVector **,
vtkInformationVector *);
virtual int FillInputPortInformation(int port, vtkInformation *info);
private:
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment