[vtkusers] Re: Re: MarchingCubes

Sean McInerney seanm at nmr.mgh.harvard.edu
Thu Jul 29 13:56:10 EDT 2004


Goodwin,

   Thanks for pointing out differences in the two methods. I had 
originally overlooked this and it well explains the inaccuracies seen in 
my result.

   If you were inspired to implement a signed version of 
vtkImplicitModeller, more power to you! It would be immensely more 
useful with this addition ... I, for one, would be happy and impressed.

-Sean

Goodwin Lawlor wrote:
> Hi Sean
> 
> Thats true too... the advantage of vtkImplicitModeller is that it uses the
> dataset topology to calculate an accurate distance to the surface.
> vtkSurfaceReconstructionFilter estimates a surface tangent to calculate the
> distance since the dataset has no surface topology (if it did why would you
> be reconstructing it?).
> 
> vtkImplicitModeller could be modified to calculate a signed distance based
> on the normals of the input.
> 
> I had thought about this awhile ago as a method to create a hex FE mesh from
> a polygonal surface. I might give it a go when I get time (a couple of
> months).
> 
> Goodwin
> 
> ----- Original Message ----- 
> From: "Sean McInerney" <seanm at nmr.mgh.harvard.edu>
> To: "Goodwin Lawlor" <goodwin.lawlor at ucd.ie>
> Cc: <vtkusers at public.kitware.com>
> Sent: Thursday, July 29, 2004 6:32 PM
> Subject: Re: [vtkusers] Re: Re: MarchingCubes
> 
> 
> 
>>Right you are, Goodwin,
>>
>>   vtkVoxelModeller is the easiest choice of all and produces a true
>>binary volume. The resultant volume, however, will be "hollow". I was
>>experimenting with vtkSurfaceReconstructionFilter instead of
>>vtkImplicitModeller because it produces a signed distance from the
>>surface. vtkImplicitModeller has no concept of "inside" versus "outside".
>>
>>-Sean
>>
>>Goodwin Lawlor wrote:
>>
>>>You could have a go with any of the sub-classes of
> 
> vtkDataSetToImageFilter.
> 
>>>vtkImplicitModeller and and vtkVoxelModeller would be worth a look
>>>
>>>Goodwin
>>>
>>>"Sean McInerney" <seanm at nmr.mgh.harvard.edu> wrote in message
>>>news:41085C71.8090500 at nmr.mgh.harvard.edu...
>>>
>>>
>>>>Jonathan,
>>>>
>>>>  I'll check out your code.
>>>>  In the meantime, I had a quicker and dirtier idea that I decided to
>>>>try out using vtkSurfaceReconstructionFilter and vtkImageThreshold to
>>>>produce an binary ImageData volume from a PolyData surface. It's rough
>>>>but interesting. The output is a volume rendering.
>>>>
>>>>-Sean
>>>>
>>>>Jonathan Bailleul wrote:
>>>>
>>>>
>>>>>Sean McInerney wrote:
>>>>>
>>>>>
>>>>>
>>>>>>Jonathan,
>>>>>>
>>>>>> That's a cool idea!
>>>>>>
>>>>>> Foreach desired slice plane for the ImageData volume that you will
>>>>>>ultimately produce, do something like this:
>>>>>>
>>>>>>1.) Use vtkClipPolyData to get the 2D contour of the given slice.
>>>>>>
>>>>>>2.) Use vtkPolyDataToImageStencil to produce ImageStencilData
>>>>>>
>>>>>>3.) Use vtkImageStencil to produce the ImageData
>>>>>>   (need to do a little filtering to produce the desired binary
>>>
>>>image).
>>>
>>>
>>>>>>4.) Use vtkAppendFilter to fuse the 2D binary images into a 3D volume.
>>>>>>   Hopefully, the ImageStencilData retains a 3D extent.
>>>>>>
>>>>>>5.) Voila! Binary volume from PolyData ... maybe.
>>>>>>
>>>>>>-Sean
>>>>>
>>>>>
>>>>>
>>>>>Thanks for the help!
>>>>>
>>>>>If someone is inspirated to actually implement this (or has already
>>>>>completed this), please notify me!
>>>>>Anyway, here is a stupid but working solution. Code is sometimes ugly,
>>>>>but it works (only voxels at border of each slice are "colored").
> 
> Output
> 
>>>>>is AIR analyse format, recognized on ITK and usable on VTK if you know
>>>>>the parameters and consider the .img file as the raw volume.
>>>>>
>>>>>All personal includes in:
>>>>>http://www.greyc.ismra.fr/~bailleul/Ressources/summary.html (see gpl
>>>>>code section)
>>>>>for those who really want to end up testing the program.
>>>>>
>>>>>I really should have been more investigating on VTK abilities before
>>>>>doing all this by hand, but I was VTK newbie when I coded this and had
>>>>>most of the necessary code at hand on my personal libraries.
>>>>>
>>>>>
>>>>
>>>
>>--------------------------------------------------------------------------
> 
> --
> 
>>>----
>>>
>>>
>>>
>>>
>>>>// VTK Common
>>>>#include "vtkPolyData.h"
>>>>#include "vtkImageData.h"
>>>>#include "vtkPiecewiseFunction.h"
>>>>#include "vtkColorTransferFunction.h"
>>>>// VTK Imaging
>>>>#include "vtkSurfaceReconstructionFilter.h"
>>>>#include "vtkImageThreshold.h"
>>>>#include "vtkImageToStructuredPoints.h"
>>>>// VTK Graphics
>>>>#include "vtkOutlineFilter.h"
>>>>// VTK IO
>>>>#include "vtkBYUReader.h"
>>>>// VTK Rendering
>>>>#include "vtkPolyDataMapper.h"
>>>>#include "vtkProperty.h"
>>>>#include "vtkActor.h"
>>>>#include "vtkVolumeRayCastCompositeFunction.h"
>>>>#include "vtkVolumeRayCastMapper.h"
>>>>#include "vtkVolumeProperty.h"
>>>>#include "vtkVolume.h"
>>>>#include "vtkRenderer.h"
>>>>#include "vtkRenderWindow.h"
>>>>#include "vtkRenderWindowInteractor.h"
>>>>
>>>>int
>>>>main (int argc, char* argv[])
>>>>{
>>>> const char* dataRoot = getenv("VTK_DATA_ROOT");
>>>>
>>>> if (dataRoot == NULL)
>>>>   {
>>>>   cerr << "Must set VTK_DATA_ROOT!" << endl;
>>>>   return 0;
>>>>   }
>>>>
>>>> const char* subdir = "/Data/";
>>>> const char* file1 = "Viewpoint/cow.g";
>>>> char* path1 = new char
>>>
>>>[strlen(dataRoot)+strlen(subdir)+strlen(file1)+1];
>>>
>>>
>>>> strcpy(path1, dataRoot);
>>>> strcat(path1, subdir);
>>>> strcat(path1, file1);
>>>>
>>>> double length;
>>>> double range[2];
>>>> double bounds[6];
>>>> int extent[6];
>>>>
>>>> vtkRenderer* ren1 = vtkRenderer::New();
>>>> vtkRenderWindow* renWin = vtkRenderWindow::New();
>>>>   {
>>>>   renWin->AddRenderer(ren1);
>>>>   renWin->SetSize(300,300);
>>>>   ren1->Delete();
>>>>   }
>>>> vtkRenderWindowInteractor* iren = vtkRenderWindowInteractor::New();
>>>>   {
>>>>   iren->SetRenderWindow(renWin);
>>>>   renWin->Delete();
>>>>   }
>>>>
>>>> vtkImageData* imageData = vtkImageData::New();
>>>>   {
>>>>   vtkSurfaceReconstructionFilter* surfaceRecon =
>>>>     vtkSurfaceReconstructionFilter::New();
>>>>     {
>>>>     vtkBYUReader* reader = vtkBYUReader::New();
>>>>       {
>>>>       reader->SetGeometryFileName(path1);
>>>>       reader->Update();
>>>>       }
>>>>     surfaceRecon->SetInput(reader->GetOutput());
>>>>     surfaceRecon->SetNeighborhoodSize(10);
>>>>
>>>
>>>surfaceRecon->SetSampleSpacing(reader->GetOutput()->GetLength()/100.0);
>>>
>>>
>>>>     surfaceRecon->Update();
>>>>     reader->Delete();
>>>>     }
>>>>   imageData->DeepCopy(surfaceRecon->GetOutput());
>>>>   surfaceRecon->Delete();
>>>>   }
>>>>
>>>> vtkActor* outlineActor = vtkActor::New();
>>>>   {
>>>>   vtkPolyDataMapper* outlineMapper = vtkPolyDataMapper::New();
>>>>     {
>>>>     vtkOutlineFilter* outlineFilter = vtkOutlineFilter::New();
>>>>       {
>>>>       outlineFilter->SetInput(imageData);
>>>>       }
>>>>     outlineMapper->SetInput(outlineFilter->GetOutput());
>>>>     outlineFilter->Delete();
>>>>     }
>>>>   outlineActor->SetMapper(outlineMapper);
>>>>   outlineMapper->Delete();
>>>>   }
>>>> ren1->AddProp(outlineActor);
>>>> outlineActor->Delete();
>>>>
>>>> renWin->Render();
>>>>
>>>> vtkVolume* volume = vtkVolume::New();
>>>>   {
>>>>   vtkVolumeRayCastMapper* volumeMapper =
> 
> vtkVolumeRayCastMapper::New();
> 
>>>>     {
>>>>     vtkVolumeRayCastCompositeFunction* compositeFunction =
>>>>       vtkVolumeRayCastCompositeFunction::New();
>>>>     vtkImageThreshold* imageThreshold = vtkImageThreshold::New();
>>>>       {
>>>>       imageThreshold->SetInput(imageData);
>>>>       imageThreshold->ReplaceInOn();
>>>>       imageThreshold->ReplaceOutOn();
>>>>       imageThreshold->SetInValue(1.0);
>>>>       imageThreshold->SetOutValue(0.0);
>>>>       imageThreshold->ThresholdByLower(0.0);
>>>>       imageThreshold->SetOutputScalarTypeToUnsignedChar();
>>>>       }
>>>>     volumeMapper->SetVolumeRayCastFunction(compositeFunction);
>>>>     volumeMapper->SetInput(imageThreshold->GetOutput());
>>>>     imageThreshold->Delete();
>>>>     }
>>>>   vtkVolumeProperty* volumeProperty = vtkVolumeProperty::New();
>>>>     {
>>>>     // Create transfer mapping scalar value to opacity
>>>>     vtkPiecewiseFunction* opacityTF = vtkPiecewiseFunction::New();
>>>>       {
>>>>       opacityTF->AddPoint(0.0, 0.167);
>>>>       opacityTF->AddPoint(1.0, 0.667);
>>>>       opacityTF->ClampingOn();
>>>>       }
>>>>     // Create transfer mapping scalar value to color
>>>>     vtkColorTransferFunction* colorTF =
> 
> vtkColorTransferFunction::New();
> 
>>>>       {
>>>>       colorTF->AddHSVPoint(0.0, 0.00, 1.0, 1.0);
>>>>       colorTF->AddHSVPoint(1.0, 0.50, 1.0, 1.0);
>>>>       }
>>>>     volumeProperty->SetColor(colorTF);
>>>>     volumeProperty->SetScalarOpacity(opacityTF);
>>>>     volumeProperty->SetInterpolationTypeToLinear();
>>>>     opacityTF->Delete();
>>>>     colorTF->Delete();
>>>>     }
>>>>   volume->SetMapper(volumeMapper);
>>>>   volume->SetProperty(volumeProperty);
>>>>   volumeMapper->Delete();
>>>>   volumeProperty->Delete();
>>>>   }
>>>> ren1->AddProp(volume);
>>>> volume->Delete();
>>>>
>>>> renWin->Render();
>>>>
>>>> iren->Start();
>>>>
>>>> iren->Delete();
>>>> imageData->Delete();
>>>>
>>>> delete [] path1;
>>>>
>>>> return 0;
>>>>}
>>>>
>>>
>>>
>>>
>>--------------------------------------------------------------------------
> 
> --
> 
>>>----
>>>
>>>
>>>
>>>
>>>>_______________________________________________
>>>>This is the private VTK discussion list.
>>>>Please keep messages on-topic. Check the FAQ at:
>>>
>>><http://public.kitware.com/cgi-bin/vtkfaq>
>>>
>>>>Follow this link to subscribe/unsubscribe:
>>>>http://www.vtk.org/mailman/listinfo/vtkusers
>>>>
>>>
>>>
>>>
>>>
>>>_______________________________________________
>>>This is the private VTK discussion list.
>>>Please keep messages on-topic. Check the FAQ at:
> 
> <http://public.kitware.com/cgi-bin/vtkfaq>
> 
>>>Follow this link to subscribe/unsubscribe:
>>>http://www.vtk.org/mailman/listinfo/vtkusers
>>>
>>
> 
> _______________________________________________
> This is the private VTK discussion list. 
> Please keep messages on-topic. Check the FAQ at: <http://public.kitware.com/cgi-bin/vtkfaq>
> Follow this link to subscribe/unsubscribe:
> http://www.vtk.org/mailman/listinfo/vtkusers
> 



More information about the vtkusers mailing list