vtkImagePolyDataSeedConnectivity Class Reference

#include <vtkImagePolyDataSeedConnectivity.h>

List of all members.

Public Member Functions
	vtkTypeRevisionMacro (vtkImagePolyDataSeedConnectivity, vtkStructuredPointsToPolyDataFilter)
void	PrintSelf (ostream &os, vtkIndent indent)
virtual void	SetAxis (vtkPolyData *)
	vtkGetObjectMacro (Axis, vtkPolyData)
	vtkSetMacro (ThresholdRatio, double)
	vtkGetMacro (ThresholdRatio, double)
	vtkGetObjectMacro (OuterMold, vtkPolyData)
Static Public Member Functions
static vtkImagePolyDataSeedConnectivity *	New ()
Protected Member Functions
	vtkImagePolyDataSeedConnectivity ()
	~vtkImagePolyDataSeedConnectivity ()
void	Execute ()
void	ClipImageWithAxis ()
Protected Attributes
double	ThresholdRatio
vtkPolyData *	Axis
vtkImageData *	ClipImageData
vtkPolyData *	OuterMold
Private Member Functions
	vtkImagePolyDataSeedConnectivity (const vtkImagePolyDataSeedConnectivity &)
void	operator= (const vtkImagePolyDataSeedConnectivity &)

---

Detailed Description

Definition at line 28 of file vtkImagePolyDataSeedConnectivity.h.

---

Constructor & Destructor Documentation

vtkImagePolyDataSeedConnectivity::vtkImagePolyDataSeedConnectivity

(

)

[protected]

Definition at line 43 of file vtkImagePolyDataSeedConnectivity.cxx.

References Axis, ClipImageData, New(), OuterMold, and ThresholdRatio.

{
  this->Axis = NULL;
  this->ClipImageData = vtkImageData::New();
  this->OuterMold = vtkPolyData::New();
  this->ThresholdRatio = 0.45;
}

Here is the call graph for this function:

vtkImagePolyDataSeedConnectivity::~vtkImagePolyDataSeedConnectivity

(

)

[protected]

Definition at line 51 of file vtkImagePolyDataSeedConnectivity.cxx.

References Axis, ClipImageData, and OuterMold.

{
  if(this->Axis)
  {
    this->Axis->Delete();
  }
  this->ClipImageData->Delete();
  this->ClipImageData = NULL;
  this->OuterMold->Delete();
  this->OuterMold = NULL;
}

vtkImagePolyDataSeedConnectivity::vtkImagePolyDataSeedConnectivity

(

const vtkImagePolyDataSeedConnectivity &

)

[private]

---

Member Function Documentation

void vtkImagePolyDataSeedConnectivity::ClipImageWithAxis

(

)

[protected]

Instead of working directly on leo we might remove the heart by a trick: If we could cut with a disk at place where the last(first) axis point is we could avoid the flood fill algo to spread in the heart...smart !

First we construct a virtual disk (=half sphere for now) See: Disk-Stencil.py

What we do here: Get first point (p2) and second point (p1) : assume it exists ! Then calculate vector p2->p1 normalize it use it to define intersection plane between sphere and plane

What we do here: Get last point (p2) and before last point (p1) : assume it exists ! Then calculate vector p2->p1 normalize it use it to define intersection plane between sphere and plane

Now use dgobbi class: Stencil stuff to remove the 'bad' heart part

Should I do my own vtkCappedCylinder ?? solution: no , cause ... I can't remember why... :p

Definition at line 285 of file vtkImagePolyDataSeedConnectivity.cxx.

References Axis, ClipImageData, and New().

Referenced by Execute().

{
  //no input need we'll work with member data:
  //Now we only need the region growing part:
  double normal[3];

  this->Axis->Update();
  vtkPoints *points = this->Axis->GetPoints();
  int numPts = points->GetNumberOfPoints();
  double p1[3];
  double p2[3];
  
  points->GetPoint( 0, p2 );
  points->GetPoint( 1, p1 );
  
  normal[0] = p1[0] - p2[0];
  normal[1] = p1[1] - p2[1];
  normal[2] = p1[2] - p2[2];
  vtkMath::Normalize( normal );

  vtkPlane *plane1 = vtkPlane::New();
  plane1->SetOrigin( p2[ 0 ], p2[ 1 ], p2[ 2 ] );
  plane1->SetNormal( normal[0], normal[1], normal[2]);

  vtkDataArray* scalars = this->Axis->GetPointData()->GetScalars();
  double tuple[2];
  scalars->GetTuple(0, tuple);
  double radius = tuple[0];

   vtkSphere *sphere1 = vtkSphere::New();
   sphere1->SetCenter( p2[ 0 ], p2[ 1 ], p2[ 2 ] );
   sphere1->SetRadius( radius + 1);

   vtkImplicitBoolean *disk1 = vtkImplicitBoolean::New();
   disk1->SetOperationTypeToIntersection();
   disk1->AddFunction(sphere1);
   disk1->AddFunction(plane1);

  //ax->getControlPoint(numPts - 1, p2, p2+3);
  points->GetPoint( numPts - 1, p2);
  //ax->getControlPoint(numPts - 2, p1, p1+3); //assume it exists
  points->GetPoint( numPts - 2, p1);
  normal[0] = p1[0] - p2[0];
  normal[1] = p1[1] - p2[1];
  normal[2] = p1[2] - p2[2];
  vtkMath::Normalize( normal );

  vtkPlane *plane2 = vtkPlane::New();
  plane2->SetOrigin( p2[ 0 ], p2[ 1 ], p2[ 2 ] );
  plane2->SetNormal( normal[0], normal[1], normal[2]);

  vtkSphere *sphere2 = vtkSphere::New();
  sphere2->SetCenter( p2[ 0 ], p2[ 1 ], p2[ 2 ] );
  sphere2->SetRadius( radius + 1 );

  vtkImplicitBoolean *disk2 = vtkImplicitBoolean::New();
  disk2->SetOperationTypeToIntersection();
  disk2->AddFunction(sphere2);
  disk2->AddFunction(plane2);

  vtkImplicitBoolean *disk = vtkImplicitBoolean::New();
  disk->SetOperationTypeToUnion();
  disk->AddFunction(disk1);
  disk->AddFunction(disk2);

  vtkImplicitFunctionToImageStencil *functionToStencil = vtkImplicitFunctionToImageStencil::New ();
  functionToStencil->SetInput ( disk );

  vtkImageStencil *stencil = vtkImageStencil::New();
  stencil->SetInput( this->GetInput() );
  stencil->SetStencil( functionToStencil->GetOutput() );
  stencil->ReverseStencilOn();
  stencil->SetBackgroundValue (0);
  //don't forget to call ->Update before using an iterator !!!
  stencil->Update();
  
  this->ClipImageData->ShallowCopy( stencil->GetOutput() );
  
  //Delete stuff:
  plane1                        -> Delete();
  sphere1                       -> Delete();
  disk1                         -> Delete();
  plane2                        -> Delete();
  sphere2                       -> Delete();
  disk2                         -> Delete();
  disk                          -> Delete();
  functionToStencil     -> Delete();
  stencil                       -> Delete();
}

Here is the call graph for this function:

Here is the caller graph for this function:

void vtkImagePolyDataSeedConnectivity::Execute

(

)

[protected]

Instead of working on the real image, we might take advantage of class Stencil made by dgobbi

Steps:

• Create a fake disk (black) of radius ~ Optimale sphere and normal orientation = last point first eigen vector See: Disk-Stencil.py for more info

Internal function that take an input of a b&w or grayscale image (should be unsigned short/ unsigned char and apply a dilatation (see math morpho for info). NB: dilation o gaussian = gaussian o dilatation

Definition at line 71 of file vtkImagePolyDataSeedConnectivity.cxx.

References Axis, ClipImageData, ClipImageWithAxis(), New(), OuterMold, PSC_MAX, PSC_MIN, and ThresholdRatio.

{
  vtkImageData *input = this->GetInput();
  vtkPolyData *output = this->GetOutput();
  int intervalle = 2;
  int dim[3];

  vtkDebugMacro(<< "Executing polydata seed connectivity...");

  if (input == NULL)
    {
    vtkErrorMacro(<<"Input is NULL");
    return;
    }

  // just deal with volumes
  if ( input->GetDataDimension() != 3 )
    {
    vtkErrorMacro("Bad input: only treats 3D structured point datasets");
    return;
    }

  input->GetDimensions(dim);

  //Allocate a temporary image that will be filled in b&w
  //this is the base of the algothm, from a graysccale input + an axis
  //we produce a b&w image the the marching cubes can contour:

  vtkImageData *imagedata = vtkImageData::New();
  imagedata->SetOrigin( input->GetOrigin() );
  imagedata->SetSpacing( input->GetSpacing() );
  imagedata->SetExtent( input->GetExtent() );
  imagedata->SetScalarTypeToUnsignedChar (); //for connectityseed + reduce mem usage
  imagedata->AllocateScalars();

  /*
  Result: piece of crap !!
  There was some leakage from another vessel that was near enough from the heart.
  */

  vtkImageSeedConnectivity *seed = vtkImageSeedConnectivity::New();
  seed->SetInput( imagedata );
  seed->SetInputConnectValue( 255 ); //If 0 -> inverse video, 255 -> rien !

  // Before starting the loop we should make sure the image has been clipped
  // with the polydata:
  this->ClipImageWithAxis();

  vtkPoints *points = this->Axis->GetPoints();
  vtkDataArray* scalars = this->Axis->GetPointData()->GetScalars();
  
  int numPts = points->GetNumberOfPoints();
  double p1[3], p2[3], p3[3], tuple[2];
  double radius, signalvalue;
  int region[6];
  int i_p2; //allow saving of p2 index
  unsigned short* inSI, *inSIEnd;
  unsigned char*  outSI,*outSIEnd;
    
  for(int i=0; i<numPts; i+=intervalle)
  {
    i_p2 = i+intervalle/2;
    points->GetPoint(i, p1);
    if(numPts <= i+2*intervalle)
    {
      i_p2 = i+intervalle;
      points->GetPoint(numPts - 1, p3);
      i = numPts;       //no more iteration
    }
    else
    {
      points->GetPoint(i+intervalle, p3);
    }
    points->GetPoint(i_p2, p2);
    
    //Now we can add p2 point to the see list:
    //seed->AddSeed( 3, p2); //too bad :(
    seed->AddSeed( (int)p2[0], (int)p2[1], (int)p2[2] );

    //this parameter seems to be ok for most image
    //furthermore it doesn't make any change to make it bigger.
    scalars->GetTuple(i_p2, tuple);
    radius              = tuple[0]; 
        signalvalue = tuple[1];
    radius *= 4;                //multiply for a factor 4 for aneurism case

    region[0] = (int)PSC_MAX(p2[0] - radius, 0);
    region[2] = (int)PSC_MAX(p2[1] - radius, 0);
    region[4] = (int)PSC_MAX(p2[2] - radius, 0);

    region[1] = (int)PSC_MIN(p2[0] + radius, dim[0] - 1);
    region[3] = (int)PSC_MIN(p2[1] + radius, dim[1] - 1);
    region[5] = (int)PSC_MIN(p2[2] + radius, dim[2] - 1);

    // multiply signal by a factor this allow user to accuratley find the best
    // range
    signalvalue *= ThresholdRatio;  

      //Solution: this is too much a trouble: rather do it once on the entry (=first  axis point)

      vtkImageIterator<unsigned short> inIt(this->ClipImageData, region);//image_woheart
      vtkImageIterator<unsigned char> outIt(imagedata, region);
      // vtkImageProgressIterator should be nicer in conjonction with a wxGauge

      // Loop through ouput pixels
      while (!inIt.IsAtEnd())
      {
        inSI = inIt.BeginSpan(); inSIEnd = inIt.EndSpan();
        outSI = outIt.BeginSpan(); outSIEnd = outIt.EndSpan();

        // Pixel operation
        while (inSI != inSIEnd)
          *outSI++ = (signalvalue <= *inSI++) ? 255 : 0;

        inIt.NextSpan(); outIt.NextSpan();
      }
   }
   
  //NEW: introduced by Maciek/Marcela: use a vtkImageGaussian to get rid 
  //of the crenelage/steps introduced by MarchingCubes
  
  vtkImageGaussianSmooth *gauss = vtkImageGaussianSmooth::New();
  gauss->SetInput( seed->GetOutput() );
  gauss->SetDimensionality(3);
  //gauss->SetStandardDeviation(1.0);
  //gauss->SetRadiusFactor(3.);

  vtkMarchingCubes *contour = vtkMarchingCubes::New();
  contour->SetInput ( gauss->GetOutput() );
  contour->SetNumberOfContours( 1 );
  contour->SetValue(0, 128);
  contour->ComputeGradientsOn ();
  contour->ComputeScalarsOn ();

  vtkCleanPolyData *clean = vtkCleanPolyData::New();
  clean->SetInput ( contour->GetOutput() );

  vtkTriangleFilter *tf = vtkTriangleFilter::New();
  tf->SetInput( clean->GetOutput() );
  tf->Update();

  //Need to do a deep copy as we'll do a ShallowCopy afterwards
  output->DeepCopy( tf->GetOutput() );

  //output->Squeeze();  //usefull ?

  // Now we do a dilatation of seed->GetOutput() in order to have inner (=output)
  // and outer mold (=OuterMold):
  //First apply a morpho math: dilatation of 2mm:
  vtkImageContinuousDilate3D *dilate = vtkImageContinuousDilate3D ::New();
  dilate->SetInput ( seed->GetOutput() );
  //TODO: determine real kernel based on scale = cm/pixel
  //See: http://public.kitware.com/pipermail/vtkusers/2000-September/004210.html
  //[Spacing is x and y is the distance between the samples. The units of measure
  //are your choice. Often medical images use millimeters. Let's say your medical
  //data resolution is 512 x 512 and your study has a 250 mm field of view. Then
  //the spacing in x and y would be 250/512 or .488 mm.]
  //therefore we have here:

  int kernelsize[3];    
  double spacing[3];
  seed->GetOutput()->GetSpacing( spacing );  //beware input should be Update()
  for(int j=0; j<3; j++)
  {
    //2mm thickness, we assume spacing units is in mm !!
    kernelsize[j] = 3; //(int) (2. / spacing[i]);
    //std::cout << "kernel=" << kernelsize[i] << ",et:" << (int)(2./0.3334) << ",";
    //beware if spacing=0.3334 -> 2/0.3334 = 5 !!!! bad bad +/-0.5 to see
  }
  dilate->SetKernelSize (kernelsize[0], kernelsize[1], kernelsize[2]);

  //VTK filters rules !!!
  gauss->SetInput( dilate->GetOutput());
  tf->Update();
  // shallow copy result into output.
  // Could be re-written, as we know there is no vert neither strips...
  OuterMold->ShallowCopy( tf->GetOutput() );
  
  //OuterMold->Squeeze();  //usefull ?

  //Delete temp stuff:
  imagedata             -> Delete();
  seed                  -> Delete();
  gauss                 -> Delete();
  contour               -> Delete();
  clean                 -> Delete();
  tf                    -> Delete();
  dilate                -> Delete();
  
  //TODO: if this filter is used many times when working on different image size
  //I guess we need to find something brighter for this->ClipImageData
}

Here is the call graph for this function:

static vtkImagePolyDataSeedConnectivity* vtkImagePolyDataSeedConnectivity::New

(

)

[static]

Referenced by ClipImageWithAxis(), Execute(), and vtkImagePolyDataSeedConnectivity().

Here is the caller graph for this function:

void vtkImagePolyDataSeedConnectivity::operator=

(

const vtkImagePolyDataSeedConnectivity &

)

[private]

void vtkImagePolyDataSeedConnectivity::PrintSelf	(	ostream &	os,
		vtkIndent	indent
	)

Definition at line 409 of file vtkImagePolyDataSeedConnectivity.cxx.

References Axis.

{
  this->Superclass::PrintSelf(os,indent);
  if ( this->Axis )
    {
    os << indent << "Axis of " << this->Axis->GetNumberOfPoints()
       << "points defined\n";
    }
  else
    {
    os << indent << "Axis not defined\n";
    }
}

virtual void vtkImagePolyDataSeedConnectivity::SetAxis

(

vtkPolyData *

)

[virtual]

Referenced by vtk3DSurfaceSTLWidget::ConstructVessel().

Here is the caller graph for this function:

vtkImagePolyDataSeedConnectivity::vtkGetMacro	(	ThresholdRatio	,
		double
	)

vtkImagePolyDataSeedConnectivity::vtkGetObjectMacro	(	OuterMold	,
		vtkPolyData
	)

vtkImagePolyDataSeedConnectivity::vtkGetObjectMacro	(	Axis	,
		vtkPolyData
	)

vtkImagePolyDataSeedConnectivity::vtkSetMacro	(	ThresholdRatio	,
		double
	)

vtkImagePolyDataSeedConnectivity::vtkTypeRevisionMacro	(	vtkImagePolyDataSeedConnectivity	,
		vtkStructuredPointsToPolyDataFilter
	)

---

Member Data Documentation

vtkPolyData* vtkImagePolyDataSeedConnectivity::Axis [protected]

Definition at line 58 of file vtkImagePolyDataSeedConnectivity.h.

Referenced by ClipImageWithAxis(), Execute(), PrintSelf(), vtkImagePolyDataSeedConnectivity(), and ~vtkImagePolyDataSeedConnectivity().

vtkImageData* vtkImagePolyDataSeedConnectivity::ClipImageData [protected]

Definition at line 60 of file vtkImagePolyDataSeedConnectivity.h.

Referenced by ClipImageWithAxis(), Execute(), vtkImagePolyDataSeedConnectivity(), and ~vtkImagePolyDataSeedConnectivity().

vtkPolyData* vtkImagePolyDataSeedConnectivity::OuterMold [protected]

Definition at line 61 of file vtkImagePolyDataSeedConnectivity.h.

Referenced by Execute(), vtkImagePolyDataSeedConnectivity(), and ~vtkImagePolyDataSeedConnectivity().

double vtkImagePolyDataSeedConnectivity::ThresholdRatio [protected]

Definition at line 56 of file vtkImagePolyDataSeedConnectivity.h.

Referenced by Execute(), and vtkImagePolyDataSeedConnectivity().

---

The documentation for this class was generated from the following files:

• /home/davila/Creatis/creaMaracasVisu/lib/maracasVisuLib/src/interface/wxWindows/widgets/include/vtkImagePolyDataSeedConnectivity.h
• /home/davila/Creatis/creaMaracasVisu/lib/maracasVisuLib/src/interface/wxWindows/widgets/include/vtkImagePolyDataSeedConnectivity.cxx