diff --git a/Filters/Core/Testing/Data/Baseline/TestElevationFilter.png.md5 b/Filters/Core/Testing/Data/Baseline/TestElevationFilter.png.md5 new file mode 100644 index 0000000000000000000000000000000000000000..da67ef9b8f16129fb762e9e572d6a2b249080db7 --- /dev/null +++ b/Filters/Core/Testing/Data/Baseline/TestElevationFilter.png.md5 @@ -0,0 +1 @@ +3e99ef9b2f4b4438084874531b97b5f5 diff --git a/Filters/Core/Testing/Python/CMakeLists.txt b/Filters/Core/Testing/Python/CMakeLists.txt index 3075a8b62ad5fbe660e2c05b27f66c73b5c601d6..04f246dc706ae721416a9c83745b7e1a5895a55b 100644 --- a/Filters/Core/Testing/Python/CMakeLists.txt +++ b/Filters/Core/Testing/Python/CMakeLists.txt @@ -15,6 +15,7 @@ vtk_add_test_python( Delaunay3DAlphaTest.py QuadricDecimation.py StreamPolyData.py + TestElevationFilter.py TestGridSynchronizedTemplates3D.py TestMarchingSquares.py TestRectilinearSynchronizedTemplates.py diff --git a/Filters/Core/Testing/Python/TestElevationFilter.py b/Filters/Core/Testing/Python/TestElevationFilter.py new file mode 100755 index 0000000000000000000000000000000000000000..2b0676219e6f69e011bd37d6465433d8a8e0eaae --- /dev/null +++ b/Filters/Core/Testing/Python/TestElevationFilter.py @@ -0,0 +1,92 @@ +#!/usr/bin/env python + +import math +import vtk +import vtk.test.Testing +from vtk.util.misc import vtkGetDataRoot +VTK_DATA_ROOT = vtkGetDataRoot() + +res = 200 + +# Create the RenderWindow, Renderer and both Actors +# +ren1 = vtk.vtkRenderer() +ren2 = vtk.vtkRenderer() +renWin = vtk.vtkRenderWindow() +renWin.AddRenderer(ren1) +renWin.AddRenderer(ren2) + +# Create the data -- a plane +# +plane = vtk.vtkPlaneSource() +plane.SetXResolution(res) +plane.SetYResolution(res) +plane.SetOrigin(-10,-10,0) +plane.SetPoint1(10,-10,0) +plane.SetPoint2(-10,10,0) + +tf = vtk.vtkTriangleFilter() +tf.SetInputConnection(plane.GetOutputPort()) +tf.Update() + +pd = vtk.vtkPolyData() +pd.CopyStructure(tf.GetOutput()) +numPts = pd.GetNumberOfPoints() +oldPts = tf.GetOutput().GetPoints() +newPts = vtk.vtkPoints() +newPts.SetNumberOfPoints(numPts) +for i in range(0, numPts): + pt = oldPts.GetPoint(i) + r = math.sqrt(pt[0]*pt[0] + pt[1]*pt[1]) + z = 1.5*math.cos(2*r) + newPts.SetPoint(i, pt[0], pt[1], z) +pd.SetPoints(newPts) + +ele = vtk.vtkSimpleElevationFilter() +ele.SetInputData(pd) + +mapper = vtk.vtkPolyDataMapper() +mapper.SetInputConnection(ele.GetOutputPort()) + +actor = vtk.vtkActor() +actor.SetMapper(mapper) + +ele2 = vtk.vtkElevationFilter() +ele2.SetInputData(pd) +ele2.SetLowPoint(0,0,-1.5) +ele2.SetHighPoint(0,0,1.5) +ele2.SetScalarRange(-1.5,1.5) + +mapper2 = vtk.vtkPolyDataMapper() +mapper2.SetInputConnection(ele2.GetOutputPort()) + +actor2 = vtk.vtkActor() +actor2.SetMapper(mapper2) + +# Add the actors to the renderer, set the background and size +# +ren1.SetViewport(0, 0, .5, 1) +ren2.SetViewport(.5, 0, 1, 1) + +ren1.AddActor(actor) +ren2.AddActor(actor2) + +camera = vtk.vtkCamera() +camera.SetPosition(1,1,1) +ren1.SetActiveCamera(camera) +ren2.SetActiveCamera(camera) + +ren1.SetBackground(0, 0, 0) +ren2.SetBackground(0, 0, 0) + +renWin.SetSize(500, 250) + +# render and interact with data + +iRen = vtk.vtkRenderWindowInteractor() +iRen.SetRenderWindow(renWin); +ren1.ResetCamera() +renWin.Render() + +iRen.Initialize() +#iRen.Start() diff --git a/Filters/Core/vtkElevationFilter.cxx b/Filters/Core/vtkElevationFilter.cxx index 4846424c6fca45e3bbda2fdcad4c77b3f75bbf27..602be248d3d4308e20646a215e01945d6559107b 100644 --- a/Filters/Core/vtkElevationFilter.cxx +++ b/Filters/Core/vtkElevationFilter.cxx @@ -16,6 +16,7 @@ #include "vtkCellData.h" #include "vtkDataSet.h" +#include "vtkPointSet.h" #include "vtkFloatArray.h" #include "vtkInformation.h" #include "vtkInformationVector.h" @@ -23,10 +24,103 @@ #include "vtkObjectFactory.h" #include "vtkPointData.h" #include "vtkSmartPointer.h" +#include "vtkSMPTools.h" vtkStandardNewMacro(vtkElevationFilter); +// The heart of the algorithm plus interface to the SMP tools. Double templated +// over point and scalar types. +template <class TP> +class vtkElevationAlgorithm +{ +public: + vtkIdType NumPts; + double LowPoint[3]; + double HighPoint[3]; + double ScalarRange[2]; + const TP *Points; + float *Scalars; + const double *V; + double L2; + + // Contructor + vtkElevationAlgorithm(); + + // Interface between VTK and templated functions + static void Elevate(vtkElevationFilter *self, vtkIdType numPts, + double v[3], double l2, TP *points, float *scalars); + + // Interface implicit function computation to SMP tools. + template <class T> class ElevationOp + { + public: + ElevationOp(vtkElevationAlgorithm<T> *algo) + { this->Algo = algo;} + vtkElevationAlgorithm *Algo; + void operator() (vtkIdType k, vtkIdType end) + { + double ns, vec[3]; + const double *range = this->Algo->ScalarRange; + const double diffScalar = range[1] - range[0]; + const double *v = this->Algo->V; + const double l2 = this->Algo->L2; + const double *lp = this->Algo->LowPoint; + const TP *p = this->Algo->Points + 3*k; + float *s = this->Algo->Scalars + k; + for ( ; k < end; ++k) + { + vec[0] = p[0] - lp[0]; + vec[1] = p[1] - lp[1]; + vec[2] = p[2] - lp[2]; + ns = (vec[0]*v[0] + vec[1]*v[1] + vec[2]*v[2]) / l2; + ns = (ns < 0.0 ? 0.0 : ns > 1.0 ? 1.0 : ns); + + // Store the resulting scalar value. + *s = range[0] + ns*diffScalar; + + p+=3; + ++s; + } + } + }; +}; + +//---------------------------------------------------------------------------- +// Initialized mainly to eliminate compiler warnings. +template <class TP> vtkElevationAlgorithm<TP>:: +vtkElevationAlgorithm():Points(NULL),Scalars(NULL) +{ + this->LowPoint[0] = this->LowPoint[1] = this->LowPoint[2] = 0.0; + this->HighPoint[0] = this->HighPoint[1] = 0.0; + this->HighPoint[2] = 1.0; + this->ScalarRange[0] = 0.0; + this->ScalarRange[1] = 1.0; +} + //---------------------------------------------------------------------------- +// Templated class is glue between VTK and templated algorithms. +template <class TP> void vtkElevationAlgorithm<TP>:: +Elevate(vtkElevationFilter *self, vtkIdType numPts, + double *v, double l2, TP *points, float *scalars) +{ + // Populate data into local storage + vtkElevationAlgorithm<TP> algo; + algo.NumPts = numPts; + self->GetLowPoint(algo.LowPoint); + self->GetHighPoint(algo.HighPoint); + self->GetScalarRange(algo.ScalarRange); + algo.Points = points; + algo.Scalars = scalars; + algo.V = v; + algo.L2 = l2; + + // Okay now generate samples using SMP tools + ElevationOp<TP> values(&algo); + vtkSMPTools::For(0,algo.NumPts, values); +} + +//---------------------------------------------------------------------------- +// Begin the class proper vtkElevationFilter::vtkElevationFilter() { this->LowPoint[0] = 0.0; @@ -100,34 +194,56 @@ int vtkElevationFilter::RequestData(vtkInformation*, length2 = 1.0; } - // Support progress and abort. - vtkIdType tenth = (numPts >= 10? numPts/10 : 1); - double numPtsInv = 1.0/numPts; - int abort = 0; - - // Compute parametric coordinate and map into scalar range. - double diffScalar = this->ScalarRange[1] - this->ScalarRange[0]; vtkDebugMacro("Generating elevation scalars!"); - for(vtkIdType i=0; i < numPts && !abort; ++i) + + // Create a fast path for point set input + // + vtkPointSet *ps = vtkPointSet::SafeDownCast(input); + if ( ps ) { - // Periodically update progress and check for an abort request. - if(i % tenth == 0) + float *scalars = + static_cast<float*>(newScalars->GetVoidPointer(0)); + vtkPoints *points = ps->GetPoints(); + void *pts = points->GetData()->GetVoidPointer(0); + switch ( points->GetDataType() ) { - this->UpdateProgress((i+1)*numPtsInv); - abort = this->GetAbortExecute(); + vtkTemplateMacro( + vtkElevationAlgorithm<VTK_TT>::Elevate(this,numPts,diffVector,length2, + (VTK_TT *)pts,scalars)); } + }//fast path + + else + { + // Too bad, got to take the scenic route. + // Support progress and abort. + vtkIdType tenth = (numPts >= 10? numPts/10 : 1); + double numPtsInv = 1.0/numPts; + int abort = 0; - // Project this input point into the 1D system. - double x[3]; - input->GetPoint(i, x); - double v[3] = { x[0] - this->LowPoint[0], - x[1] - this->LowPoint[1], - x[2] - this->LowPoint[2] }; - double s = vtkMath::Dot(v, diffVector) / length2; - s = (s < 0.0 ? 0.0 : s > 1.0 ? 1.0 : s); - - // Store the resulting scalar value. - newScalars->SetValue(i, this->ScalarRange[0] + s*diffScalar); + // Compute parametric coordinate and map into scalar range. + double diffScalar = this->ScalarRange[1] - this->ScalarRange[0]; + for(vtkIdType i=0; i < numPts && !abort; ++i) + { + // Periodically update progress and check for an abort request. + if(i % tenth == 0) + { + this->UpdateProgress((i+1)*numPtsInv); + abort = this->GetAbortExecute(); + } + + // Project this input point into the 1D system. + double x[3]; + input->GetPoint(i, x); + double v[3] = { x[0] - this->LowPoint[0], + x[1] - this->LowPoint[1], + x[2] - this->LowPoint[2] }; + double s = vtkMath::Dot(v, diffVector) / length2; + s = (s < 0.0 ? 0.0 : s > 1.0 ? 1.0 : s); + + // Store the resulting scalar value. + newScalars->SetValue(i, this->ScalarRange[0] + s*diffScalar); + } } // Copy all the input geometry and data to the output. diff --git a/Filters/Core/vtkElevationFilter.h b/Filters/Core/vtkElevationFilter.h index 43740a59c6bd82e5503258b7dbfdbe97ddf48b1a..251991eaff4ffee91c98693f231d266f4b049348 100644 --- a/Filters/Core/vtkElevationFilter.h +++ b/Filters/Core/vtkElevationFilter.h @@ -20,6 +20,17 @@ // a line. The line can be oriented arbitrarily. A typical example is // to generate scalars based on elevation or height above a plane. +// .SECTION Caveats +// vtkSimpleElevationFilter may be easier to use in many cases; e.g., +// compute vertical elevation above zero z-point. +// +// This class has been threaded with vtkSMPTools. Using TBB or other +// non-sequential type (set in the CMake variable +// VTK_SMP_IMPLEMENTATION_TYPE) may improve performance significantly. + +// .SECTION See Also +// vtkSimpleElevationFilter + #ifndef vtkElevationFilter_h #define vtkElevationFilter_h @@ -61,6 +72,7 @@ protected: double LowPoint[3]; double HighPoint[3]; double ScalarRange[2]; + private: vtkElevationFilter(const vtkElevationFilter&); // Not implemented. void operator=(const vtkElevationFilter&); // Not implemented. diff --git a/Filters/Core/vtkSimpleElevationFilter.cxx b/Filters/Core/vtkSimpleElevationFilter.cxx index fda63faa0e35d7f4f230295caf552073221128ed..0e014f8f9d41f6d3d4a6a362a7109f4edce22f78 100644 --- a/Filters/Core/vtkSimpleElevationFilter.cxx +++ b/Filters/Core/vtkSimpleElevationFilter.cxx @@ -15,18 +15,87 @@ #include "vtkSimpleElevationFilter.h" #include "vtkCellData.h" -#include "vtkDataSet.h" +#include "vtkPointSet.h" #include "vtkFloatArray.h" #include "vtkInformation.h" #include "vtkInformationVector.h" #include "vtkMath.h" #include "vtkObjectFactory.h" #include "vtkPointData.h" +#include "vtkSMPTools.h" vtkStandardNewMacro(vtkSimpleElevationFilter); -// Construct object with LowPoint=(0,0,0) and HighPoint=(0,0,1). Scalar -// range is (0,1). +// The heart of the algorithm plus interface to the SMP tools. Double templated +// over point and scalar types. +template <class TP> +class vtkSimpleElevationAlgorithm +{ +public: + vtkIdType NumPts; + double Vector[3]; + const TP *Points; + float *Scalars; + + // Contructor + vtkSimpleElevationAlgorithm(); + + // Interface between VTK and templated functions + static void Elevate(vtkSimpleElevationFilter *self, vtkIdType numPts, + const TP *points, float *scalars); + + // Interface implicit function computation to SMP tools. + template <class T> class ElevationOp + { + public: + ElevationOp(vtkSimpleElevationAlgorithm<T> *algo) + { this->Algo = algo;} + vtkSimpleElevationAlgorithm *Algo; + void operator() (vtkIdType k, vtkIdType end) + { + const double *v = this->Algo->Vector; + const TP *p = this->Algo->Points + 3*k; + float *s = this->Algo->Scalars + k; + for ( ; k < end; ++k) + { + *s = v[0]*p[0] + v[1]*p[1] + v[2]*p[2]; + p+=3; + ++s; + } + } + }; +}; + +//---------------------------------------------------------------------------- +// Initialized mainly to eliminate compiler warnings. +template <class TP> vtkSimpleElevationAlgorithm<TP>:: +vtkSimpleElevationAlgorithm():Points(NULL),Scalars(NULL) +{ + this->Vector[0] = this->Vector[1] = this->Vector[2] = 0.0; +} + +//---------------------------------------------------------------------------- +// Templated class is glue between VTK and templated algorithms. +template <class TP> void vtkSimpleElevationAlgorithm<TP>:: +Elevate(vtkSimpleElevationFilter *self, vtkIdType numPts, + const TP *points, float *scalars) +{ + // Populate data into local storage + vtkSimpleElevationAlgorithm<TP> algo; + algo.NumPts = numPts; + self->GetVector(algo.Vector); + algo.Points = points; + algo.Scalars = scalars; + + // Okay now generate samples using SMP tools + ElevationOp<TP> values(&algo); + vtkSMPTools::For(0,algo.NumPts, values); +} + + +// Okay begin the class proper +//---------------------------------------------------------------------------- +// Construct object with Vector=(0,0,1). vtkSimpleElevationFilter::vtkSimpleElevationFilter() { this->Vector[0] = 0.0; @@ -34,7 +103,8 @@ vtkSimpleElevationFilter::vtkSimpleElevationFilter() this->Vector[2] = 1.0; } -// Convert position along ray into scalar value. Example use includes +//---------------------------------------------------------------------------- +// Convert position along the ray into scalar value. Example use includes // coloring terrain by elevation. // int vtkSimpleElevationFilter::RequestData( @@ -60,7 +130,7 @@ int vtkSimpleElevationFilter::RequestData( // vtkDebugMacro(<<"Generating elevation scalars!"); - // First, copy the input to the output as a starting point + // First, copy the input to the output as a starting point output->CopyStructure( input ); if ( ((numPts=input->GetNumberOfPoints()) < 1) ) @@ -74,7 +144,7 @@ int vtkSimpleElevationFilter::RequestData( newScalars = vtkFloatArray::New(); newScalars->SetNumberOfTuples(numPts); - // Set up 1D parametric system + // Ensure that there is a valid vector // if ( vtkMath::Dot(this->Vector,this->Vector) == 0.0) { @@ -82,21 +152,41 @@ int vtkSimpleElevationFilter::RequestData( this->Vector[0] = this->Vector[1] = 0.0; this->Vector[2] = 1.0; } - // Compute dot product + // Create a fast path for point set input // - int abort=0; - vtkIdType progressInterval=numPts/20 + 1; - for (i=0; i<numPts && !abort; i++) + vtkPointSet *ps = vtkPointSet::SafeDownCast(input); + if ( ps ) { - if ( ! (i % progressInterval) ) + float *scalars = + static_cast<float*>(newScalars->GetVoidPointer(0)); + vtkPoints *points = ps->GetPoints(); + void *pts = points->GetData()->GetVoidPointer(0); + switch ( points->GetDataType() ) { - this->UpdateProgress ((double)i/numPts); - abort = this->GetAbortExecute(); + vtkTemplateMacro(vtkSimpleElevationAlgorithm<VTK_TT>:: + Elevate(this,numPts,(VTK_TT *)pts,scalars)); } + }//fast path - input->GetPoint(i,x); - s = vtkMath::Dot(this->Vector,x); - newScalars->SetComponent(i,0,s); + else + { + // Too bad, got to take the scenic route. + // Compute dot product. + // + int abort=0; + vtkIdType progressInterval=numPts/20 + 1; + for (i=0; i<numPts && !abort; i++) + { + if ( ! (i % progressInterval) ) + { + this->UpdateProgress ((double)i/numPts); + abort = this->GetAbortExecute(); + } + + input->GetPoint(i,x); + s = vtkMath::Dot(this->Vector,x); + newScalars->SetComponent(i,0,s); + } } // Update self @@ -114,6 +204,7 @@ int vtkSimpleElevationFilter::RequestData( return 1; } +//---------------------------------------------------------------------------- void vtkSimpleElevationFilter::PrintSelf(ostream& os, vtkIndent indent) { this->Superclass::PrintSelf(os,indent); diff --git a/Filters/Core/vtkSimpleElevationFilter.h b/Filters/Core/vtkSimpleElevationFilter.h index a59c819b28b2f7ff69f0015feb265fac91e11235..d7d180244ae8a3f9e961c90a2a435b9cbc031420 100644 --- a/Filters/Core/vtkSimpleElevationFilter.h +++ b/Filters/Core/vtkSimpleElevationFilter.h @@ -18,6 +18,18 @@ // dataset. The scalar values are generated by dotting a user-specified // vector against a vector defined from the input dataset points to the // origin. + +// .SECTION Caveats +// This class has been threaded with vtkSMPTools. Using TBB or other +// non-sequential type (set in the CMake variable +// VTK_SMP_IMPLEMENTATION_TYPE) may improve performance significantly. Note +// however that the associated threading/templating implements a fast path, +// which will only kick in when the input has explicit point representation +// (e.g. a vtkPointSet). +// +// See also vtkElevationFilter provides more control over the operation, +// including clamping the output scalars within a range. + // .SECTION See Also // vtkElevationFilter @@ -34,11 +46,11 @@ public: void PrintSelf(ostream& os, vtkIndent indent); // Description: - // Construct object with Vector=(0,0,1); + // Construct object with Vector=(0,0,1). static vtkSimpleElevationFilter *New(); // Description: - // Define one end of the line (small scalar values). + // Define the vector with which to dot against. vtkSetVector3Macro(Vector,double); vtkGetVectorMacro(Vector,double,3); @@ -54,5 +66,3 @@ private: }; #endif - -