Interpolate additional points on lines. More...
#include <mitkInterpolateLinesFilter.h>
Public Types | |
typedef InterpolateLinesFilter | Self |
typedef SurfaceSource | Superclass |
typedef itk::SmartPointer< Self > | Pointer |
typedef itk::SmartPointer < const Self > | ConstPointer |
Public Member Functions | |
virtual const char * | GetClassName () const |
virtual void | GenerateOutputInformation () |
virtual void | GenerateData () |
const mitk::Mesh * | GetInput (void) |
virtual void | SetInput (const mitk::Mesh *image) |
virtual unsigned int | GetSplineResolution () |
Get spline resolution. | |
virtual void | SetSplineResolution (unsigned int _arg) |
Set spline resolution. | |
virtual const mitk::Geometry2D * | GetGeometryForInterpolation () |
Get/set geometry for interpolation. | |
virtual void | SetGeometryForInterpolation (mitk::Geometry2D *_arg) |
virtual ScalarType | GetLength () |
Get the overall length of the interpolated lines. | |
Static Public Member Functions | |
static Pointer | New () |
Protected Member Functions | |
InterpolateLinesFilter () | |
virtual | ~InterpolateLinesFilter () |
void | BuildPointAndVectorList (mitk::Mesh::CellType &cell, vtkPoints *points, vtkCellArray *cellarray) |
Protected Attributes | |
unsigned int | m_SplineResolution |
Spline resolution of created Surface. | |
mitk::Geometry2D::Pointer | m_GeometryForInterpolation |
vtkCardinalSpline * | m_SpX |
vtkCardinalSpline * | m_SpY |
vtkCardinalSpline * | m_SpZ |
ScalarType | m_Length |
Interpolate additional points on lines.
If a Geometry2D is set (by SetGeometryForInterpolation), we do an interpolation in the 2D-space of the Geometry2D: Map two neighboring original points on the Geometry2D, resulting in two 2D-points, interpolate in 2D between them, and map them back via the Geometry2D in the 3D-world.
Definition at line 46 of file mitkInterpolateLinesFilter.h.
typedef itk::SmartPointer<const Self> mitk::InterpolateLinesFilter::ConstPointer |
Definition at line 49 of file mitkInterpolateLinesFilter.h.
typedef itk::SmartPointer<Self> mitk::InterpolateLinesFilter::Pointer |
Definition at line 49 of file mitkInterpolateLinesFilter.h.
Definition at line 49 of file mitkInterpolateLinesFilter.h.
Definition at line 49 of file mitkInterpolateLinesFilter.h.
mitk::InterpolateLinesFilter::InterpolateLinesFilter | ( | ) | [protected] |
Definition at line 30 of file mitkInterpolateLinesFilter.cpp.
References m_SpX, m_SpY, m_SpZ, and New().
: m_SplineResolution(10), m_GeometryForInterpolation(NULL), m_Length(0.0) { m_SpX=vtkCardinalSpline::New(); m_SpY=vtkCardinalSpline::New(); m_SpZ=vtkCardinalSpline::New(); }
mitk::InterpolateLinesFilter::~InterpolateLinesFilter | ( | ) | [protected, virtual] |
Definition at line 39 of file mitkInterpolateLinesFilter.cpp.
void mitk::InterpolateLinesFilter::BuildPointAndVectorList | ( | mitk::Mesh::CellType & | cell, |
vtkPoints * | points, | ||
vtkCellArray * | cellarray | ||
) | [protected] |
Definition at line 94 of file mitkInterpolateLinesFilter.cpp.
References mitk::PointSet::GetPoint(), and QuadProgPP::t().
{ const mitk::Mesh* input = GetInput(); Mesh::PointIdIterator ptIt; Mesh::PointIdIterator ptEnd; ptEnd = cell.PointIdsEnd(); Point3D pt; int i, size=cell.GetNumberOfPoints(); int closed_loop_pre_load=0;//m_SplineResolution; if(m_GeometryForInterpolation.IsNull()) { //bei geschlossener Kontur: vor dem ersten Punkt die zwei letzten einfügen für glatten Übergang ptIt = ptEnd; ptIt-=closed_loop_pre_load+1; for(i=0;i<closed_loop_pre_load;++i, ++ptIt) { pt = input->GetPoint(*ptIt); m_SpX->AddPoint(i, pt[0]); m_SpY->AddPoint(i, pt[1]); m_SpZ->AddPoint(i, pt[2]); } //Punkte einfügen for(ptIt = cell.PointIdsBegin();i<size+closed_loop_pre_load;++i, ++ptIt) { pt = input->GetPoint(*ptIt); m_SpX->AddPoint(i, pt[0]); m_SpY->AddPoint(i, pt[1]); m_SpZ->AddPoint(i, pt[2]); } //bei geschlossener Kontur: nach dem letzten Punkt die zwei ersten einfügen für glatten Übergang int j; for(j=0,ptIt = cell.PointIdsBegin();j<closed_loop_pre_load;++j,++i, ++ptIt) { pt = input->GetPoint(*ptIt); m_SpX->AddPoint(i, pt[0]); m_SpY->AddPoint(i, pt[1]); m_SpZ->AddPoint(i, pt[2]); } //bool first = true; Point3D lastPt, firstPt; Vector3D vec; float t, step=1.0f/m_SplineResolution; size=(size-1)*m_SplineResolution; i=closed_loop_pre_load; cellarray->InsertNextCell(size); for(t=closed_loop_pre_load;i<size+closed_loop_pre_load;++i, t+=step) { float pt[3]; FillVector3D(pt, m_SpX->Evaluate(t), m_SpY->Evaluate(t), m_SpZ->Evaluate(t)); cellarray->InsertCellPoint(points->InsertNextPoint(pt)); } } else //m_GeometryForInterpolation!=NULL { Point2D pt2d; //bei geschlossener Kontur: vor dem ersten Punkt die zwei letzten einfügen für glatten Übergang ptIt = ptEnd; ptIt-=closed_loop_pre_load+1; for(i=0;i<closed_loop_pre_load;++i, ++ptIt) { m_GeometryForInterpolation->Map(input->GetPoint(*ptIt), pt2d); m_SpX->AddPoint(i, pt2d[0]); m_SpY->AddPoint(i, pt2d[1]); } //Punkte einfügen for(ptIt = cell.PointIdsBegin();i<size+closed_loop_pre_load;++i, ++ptIt) { m_GeometryForInterpolation->Map(input->GetPoint(*ptIt), pt2d); m_SpX->AddPoint(i, pt2d[0]); m_SpY->AddPoint(i, pt2d[1]); } //bei geschlossener Kontur: nach dem letzten Punkt die zwei ersten einfügen für glatten Übergang int j; for(j=0,ptIt = cell.PointIdsBegin();j<closed_loop_pre_load;++j,++i, ++ptIt) { m_GeometryForInterpolation->Map(input->GetPoint(*ptIt), pt2d); m_SpX->AddPoint(i, pt2d[0]); m_SpY->AddPoint(i, pt2d[1]); } bool first = true; Point3D lastPt; lastPt.Fill(0); Vector3D vec; float t, step=1.0f/m_SplineResolution; size=(size-1)*m_SplineResolution; i=closed_loop_pre_load; cellarray->InsertNextCell(size); for(t=closed_loop_pre_load;i<size+closed_loop_pre_load;++i, t+=step) { pt2d[0] = m_SpX->Evaluate(t); pt2d[1] = m_SpY->Evaluate(t); m_GeometryForInterpolation->Map(pt2d, pt); if(first==false) { vec=pt-lastPt; m_Length+=vec.GetNorm(); } first=false; float pvtk[3]; itk2vtk(pt, pvtk); cellarray->InsertCellPoint(points->InsertNextPoint(pvtk)); lastPt = pt; } } }
void mitk::InterpolateLinesFilter::GenerateData | ( | ) | [virtual] |
Definition at line 61 of file mitkInterpolateLinesFilter.cpp.
{ mitk::Mesh::ConstPointer input = this->GetInput(); mitk::Surface::Pointer output = this->GetOutput(); vtkPolyData *polyData = vtkPolyData::New(); vtkPoints *points = vtkPoints::New(); vtkCellArray *cellarray = vtkCellArray::New(); mitk::Mesh::PointType thisPoint; m_Length = 0.0; //iterate through all cells and build tubes Mesh::ConstCellIterator cellIt, cellEnd; cellEnd = input->GetMesh()->GetCells()->End(); for( cellIt = input->GetMesh()->GetCells()->Begin(); cellIt != cellEnd; ++cellIt ) { if(((*cellIt->Value()).GetType()==mitk::Mesh::CellType::POLYGON_CELL) && ((*cellIt->Value()).GetNumberOfPoints()>=2)) BuildPointAndVectorList(*cellIt->Value(), points, cellarray); } polyData->SetPoints( points ); points->Delete(); polyData->SetLines( cellarray ); cellarray->Delete(); output->SetVtkPolyData(polyData); polyData->Delete(); }
void mitk::InterpolateLinesFilter::GenerateOutputInformation | ( | void | ) | [virtual] |
Definition at line 46 of file mitkInterpolateLinesFilter.cpp.
{ mitk::Mesh::ConstPointer input = this->GetInput(); mitk::Surface::Pointer output = this->GetOutput(); itkDebugMacro(<<"GenerateOutputInformation()"); if(input.IsNull()) return; if(m_GeometryForInterpolation.IsNotNull()) output->SetGeometry(static_cast<Geometry3D*>(m_GeometryForInterpolation->Clone().GetPointer())); else output->SetGeometry(static_cast<Geometry3D*>(input->GetGeometry()->Clone().GetPointer())); }
virtual const char* mitk::InterpolateLinesFilter::GetClassName | ( | ) | const [virtual] |
virtual const mitk::Geometry2D* mitk::InterpolateLinesFilter::GetGeometryForInterpolation | ( | ) | [virtual] |
Get/set geometry for interpolation.
If this is set (not NULL), we do an interpolation in the 2D-space of the Geometry2D: Map two neighboring original points on the Geometry2D, resulting in two 2D-points, interpolate in 2D between them, and map them back via the Geometry2D in the 3D-world.
const mitk::Mesh * mitk::InterpolateLinesFilter::GetInput | ( | void | ) |
Definition at line 196 of file mitkInterpolateLinesFilter.cpp.
{ if (this->GetNumberOfInputs() < 1) { return 0; } return static_cast<const mitk::Mesh * > (this->ProcessObject::GetInput(0) ); }
virtual ScalarType mitk::InterpolateLinesFilter::GetLength | ( | ) | [virtual] |
Get the overall length of the interpolated lines.
virtual unsigned int mitk::InterpolateLinesFilter::GetSplineResolution | ( | ) | [virtual] |
Get spline resolution.
static Pointer mitk::InterpolateLinesFilter::New | ( | ) | [static] |
Reimplemented from mitk::SurfaceSource.
Referenced by InterpolateLinesFilter().
virtual void mitk::InterpolateLinesFilter::SetGeometryForInterpolation | ( | mitk::Geometry2D * | _arg ) | [virtual] |
void mitk::InterpolateLinesFilter::SetInput | ( | const mitk::Mesh * | image ) | [virtual] |
Definition at line 207 of file mitkInterpolateLinesFilter.cpp.
{
// Process object is not const-correct so the const_cast is required here
this->ProcessObject::SetNthInput(0,
const_cast< mitk::Mesh * >( input ) );
}
virtual void mitk::InterpolateLinesFilter::SetSplineResolution | ( | unsigned int | _arg ) | [virtual] |
Set spline resolution.
Definition at line 97 of file mitkInterpolateLinesFilter.h.
ScalarType mitk::InterpolateLinesFilter::m_Length [protected] |
Definition at line 100 of file mitkInterpolateLinesFilter.h.
unsigned int mitk::InterpolateLinesFilter::m_SplineResolution [protected] |
Spline resolution of created Surface.
Definition at line 95 of file mitkInterpolateLinesFilter.h.
vtkCardinalSpline* mitk::InterpolateLinesFilter::m_SpX [protected] |
Definition at line 98 of file mitkInterpolateLinesFilter.h.
Referenced by InterpolateLinesFilter().
vtkCardinalSpline * mitk::InterpolateLinesFilter::m_SpY [protected] |
Definition at line 98 of file mitkInterpolateLinesFilter.h.
Referenced by InterpolateLinesFilter().
vtkCardinalSpline * mitk::InterpolateLinesFilter::m_SpZ [protected] |
Definition at line 98 of file mitkInterpolateLinesFilter.h.
Referenced by InterpolateLinesFilter().