<br>Hi Somi,<br><br><br>The problem is that a Mesh Spatial object assumes<br>that "inside" means: <br><br> "in the 2D manifold of its surface"<br><br>So, the points of space that are "inside" your mesh<br>
are the ones "close enough to the mesh surface".<br><br><br>If you look at the source code implementation of this<br>class, the test is done by checking whether a given<br>point is inside any of the triangles of the Mesh.<br>
<br><br>--<br><br><br>You probably should use the itkEllipseSpatialObject<br>along with the SpatialObjectToImageFilter<br><br>as illustrated in the example:<br><br> Insight/Examples/Filtering/<br> SpatialObjectToImage1.cxx<br>
<br><br>This example shows you how to generate a rasterized<br>Sphere. (you probably would like to remove the cylinder<br>that is also used in that source code example).<br><br><br><br> Regards,<br><br><br> Luis<br>
<br><br>-------------------------------------------------------------------------<br><div class="gmail_quote">On Thu, Apr 8, 2010 at 6:32 PM, somi <span dir="ltr"><<a href="mailto:seesomi@gmail.com">seesomi@gmail.com</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;">Hi Luis,<br>Thanks for your reply. itkTubeSpatialObject worked for me, I feed stupid not to think of this before :)<br>
I took the following approach:<br>vtkConeSource->VtkPolyData -> ITKMesh -> ITKMeshSpatialObject->ITKSpatialObjectToImage<br>
<br>While working on this I got an unexpected behavior as described below:<br>I used the follow pipeline:<br><br>vtkShrereSource->VtkPolyData -> ITKMesh -> ITKMeshSpatialObject->ITKSpatialObjectToImage<br>
<br><br>I expected to get a solid spherical object in the final image, but I got a hollow sphere (see attached screenshot).<br>Is this expected behavior ?<br><br>I used the following code:<br>///////<br>//Code Starts<br>
////<br>#include <iostream><br>#include <iomanip><br>#include <itkBinaryMask3DMeshSource.h><br>#include <itkMesh.h><br>#include <itkMeshSpatialObject.h><br>#include <vtkPolyData.h><br>
#include <vtkCellArray.h><br>
#include <vtkSphereSource.h><br>#include <vtkCellLinks.h><br>#include <itkImage.h><br><br>#include <iostream><br>#include <iomanip><br>#include <vtkSmartPointer.h><br>#include <itkSpatialObjectToImageFilter.h><br>
#include <itkImageFileWriter.h><br><br>#if defined(_MSC_VER)<br>#pragma warning ( disable : 4786 )<br>#endif<br><br>#ifdef __BORLANDC__<br>#define ITK_LEAN_AND_MEAN<br>#endif<br><br><br><br><br>int main( int argc, char *argv[] )<br>
{<br><br> // Generate a Sphere in VTK<br> vtkSmartPointer<vtkSphereSource> polyData = vtkSmartPointer<vtkSphereSource>::New();<br> polyData->SetRadius(10);<br> polyData->SetThetaResolution(10);<br>
polyData->SetPhiResolution(10);<br> polyData->SetCenter(80,50,50);<br> polyData->Update();<br><br>#ifndef vtkFloatingPointType<br>#define vtkFloatingPointType float<br>#endif<br><br><br><br> //<br> // Convert VTK polydata to ITK Mesh<br>
//<br><br> const unsigned int PointDimension = 3;<br> const unsigned int MaxCellDimension = 2;<br><br> typedef itk::DefaultStaticMeshTraits<<br> vtkFloatingPointType,<br> PointDimension,<br>
MaxCellDimension,<br> vtkFloatingPointType,<br> vtkFloatingPointType > MeshTraits;<br><br><br> typedef itk::Mesh<<br> vtkFloatingPointType,<br> PointDimension,<br>
MeshTraits > MeshType;<br><br><br> MeshType::Pointer mesh = MeshType::New();<br><br> //<br> // Transfer the points from the vtkPolyData into the itk::Mesh<br> //<br> const unsigned int numberOfPoints = polyData->GetOutput()->GetNumberOfPoints();<br>
<br> vtkPoints * vtkpoints = polyData->GetOutput()->GetPoints();<br><br> mesh->GetPoints()->Reserve( numberOfPoints );<br><br> for(unsigned int p =0; p < numberOfPoints; p++)<br> {<br><br> vtkFloatingPointType * apoint = vtkpoints->GetPoint( p );<br>
<br> mesh->SetPoint( p, MeshType::PointType( apoint ));<br><br> }<br><br><br> //<br> // Transfer the cells from the vtkPolyData into the itk::Mesh<br> //<br> vtkCellArray * triangleStrips = polyData->GetOutput()->GetStrips();<br>
<br><br> vtkIdType * cellPoints;<br> vtkIdType numberOfCellPoints;<br><br> //<br> // First count the total number of triangles from all the triangle strips.<br> //<br> unsigned int numberOfTriangles = 0;<br>
<br> triangleStrips->InitTraversal();<br><br> while( triangleStrips->GetNextCell( numberOfCellPoints, cellPoints ) )<br> {<br> numberOfTriangles += numberOfCellPoints-2;<br> }<br><br><br> vtkCellArray * polygons = polyData->GetOutput()->GetPolys();<br>
<br> polygons->InitTraversal();<br><br> while( polygons->GetNextCell( numberOfCellPoints, cellPoints ) )<br> {<br> if( numberOfCellPoints == 3 )<br> {<br> numberOfTriangles ++;<br>
}<br>
}<br><br><br><br> //<br> // Reserve memory in the itk::Mesh for all those triangles<br> //<br> mesh->GetCells()->Reserve( numberOfTriangles );<br><br><br> //<br> // Copy the triangles from vtkPolyData into the itk::Mesh<br>
//<br> //<br><br> typedef MeshType::CellType CellType;<br><br> typedef itk::TriangleCell< CellType > TriangleCellType;<br><br> int cellId = 0;<br><br> // first copy the triangle strips<br> triangleStrips->InitTraversal();<br>
while( triangleStrips->GetNextCell( numberOfCellPoints, cellPoints ) )<br> {<br><br> unsigned int numberOfTrianglesInStrip = numberOfCellPoints - 2;<br><br> unsigned long pointIds[3];<br> pointIds[0] = cellPoints[0];<br>
pointIds[1] = cellPoints[1];<br> pointIds[2] = cellPoints[2];<br><br> for( unsigned int t=0; t < numberOfTrianglesInStrip; t++ )<br> {<br> MeshType::CellAutoPointer c;<br> TriangleCellType * tcell = new TriangleCellType;<br>
tcell->SetPointIds( pointIds );<br> c.TakeOwnership( tcell );<br> mesh->SetCell( cellId, c );<br> cellId++;<br> pointIds[0] = pointIds[1];<br> pointIds[1] = pointIds[2];<br>
pointIds[2] = cellPoints[t+3];<br> }<br><br><br> }<br><br><br> // then copy the normal triangles<br> polygons->InitTraversal();<br> while( polygons->GetNextCell( numberOfCellPoints, cellPoints ) )<br>
{<br> if( numberOfCellPoints !=3 ) // skip any non-triangle.<br> {<br> continue;<br> }<br> MeshType::CellAutoPointer c;<br> TriangleCellType * t = new TriangleCellType;<br>
t->SetPointIds( (unsigned long*)cellPoints );<br> c.TakeOwnership( t );<br> mesh->SetCell( cellId, c );<br> cellId++;<br> }<br><br><br><br> std::cout << "Mesh " << std::endl;<br>
std::cout << "Number of Points = " << mesh->GetNumberOfPoints() << std::endl;<br> std::cout << "Number of Cells = " << mesh->GetNumberOfCells() << std::endl;<br>
<br><br><br><br> //<br> // Convert ITK Mesh to ITK Mesh spatial object<br> //<br><br> typedef itk::MeshSpatialObject< MeshType > SpatialObjectType;<br> SpatialObjectType::Pointer meshSpa = SpatialObjectType::New();<br>
meshSpa->SetMesh( mesh);<br> meshSpa->Update();<br><br><br><br> typedef unsigned char PixelType;<br> const unsigned int Dimension = 3;<br> typedef itk::Image< PixelType, Dimension > ImageType;<br>
<br> typedef itk::SpatialObjectToImageFilter<<br> SpatialObjectType, ImageType > SpatialObjectToImageFilterType;<br><br> SpatialObjectToImageFilterType::Pointer imageFilter =<br> SpatialObjectToImageFilterType::New();<br>
ImageType::SizeType size;<br> size[ 0 ] = 100;<br> size[ 1 ] = 100;<br> size[ 2 ] = 100;<br> imageFilter->SetSize( size );<br> ImageType::SpacingType spacing;<br> spacing[0] = 100.0 / size[0];<br>
spacing[1] = 100.0 / size[1];<br> spacing[2] = 100.0 / size[2];<br><br> imageFilter->SetSpacing( spacing );<br> imageFilter->SetInput( meshSpa );<br><br> ImageType::PointType origin;<br> origin[0]=0;origin[1]=0;origin[2]=0;<br>
imageFilter->SetOrigin(origin);<br> const PixelType outValue = 0;<br> const PixelType inValue = 1;<br> meshSpa->SetDefaultInsideValue( inValue );<br> meshSpa->SetDefaultOutsideValue( outValue );<br>
<br> imageFilter->SetUseObjectValue( true );<br> imageFilter->SetOutsideValue( airHounsfieldUnits );<br><br> // Write the image<br> typedef itk::ImageFileWriter< ImageType > WriterType;<br> WriterType::Pointer writer = WriterType::New();<br>
<br> writer->SetFileName("Output.nii");<br> writer->SetInput( imageFilter->GetOutput() );<br> try<br> {<br> imageFilter->Update();<br> writer->Update();<br> }<br> catch( itk::ExceptionObject & excp )<br>
{<br> std::cerr << excp << std::endl;<br> return EXIT_FAILURE;<br> }<br><br><br><br> return EXIT_SUCCESS;<br>}<br>/////////<br>///code ends<br>/////<br><br>Thanks,<div><div></div><div class="h5">
<br><br><div class="gmail_quote">
On Thu, Apr 8, 2010 at 2:51 PM, Luis Ibanez <span dir="ltr"><<a href="mailto:luis.ibanez@kitware.com" target="_blank">luis.ibanez@kitware.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;">
<div>Hi Somi,<br>
<br>
<br>
Nope, there is no ConeSpatialObject.<br>
<br>
<br>
However, you could easily create this new class by reusing the<br>
code of the CylinderSpatialObject.<br>
<br>
or<br>
<br>
you could fake one by using the<br>
<br>
itkTubeSpatialObject<br>
<br>
with only two points, one of them having Radius == 0.<br>
<br>
<br>
Regards,<br>
<br>
<br>
Luis<br>
<br>
<br>
-------------------------------------------------------------<br>
On Mon, Apr 5, 2010 at 1:56 PM, somi <<a href="mailto:seesomi@gmail.com" target="_blank">seesomi@gmail.com</a>> wrote:<br>
</div><div><div></div><div>> Hi,<br>
> Does a cone spatial object exist in ITK ?<br>
> In my application a user interactively places a cone over an ITK image<br>
> volume.<br>
><br>
> I then have to voxelize the cone. However volexiling in VTK takes a lot of<br>
> time (I tried vtkImplicitModeller,vtkSelectEnclosedPoints they were too<br>
> slow).<br>
> So is it possible to pass the parameters of the cone to itk and generate it<br>
> using some spatial object ?<br>
><br>
> Regards,<br>
> Somesh<br>
><br>
><br>
</div></div><div><div></div><div>> _____________________________________<br>
> Powered by <a href="http://www.kitware.com" target="_blank">www.kitware.com</a><br>
><br>
> Visit other Kitware open-source projects at<br>
> <a href="http://www.kitware.com/opensource/opensource.html" target="_blank">http://www.kitware.com/opensource/opensource.html</a><br>
><br>
> Kitware offers ITK Training Courses, for more information visit:<br>
> <a href="http://www.kitware.com/products/protraining.html" target="_blank">http://www.kitware.com/products/protraining.html</a><br>
><br>
> Please keep messages on-topic and check the ITK FAQ at:<br>
> <a href="http://www.itk.org/Wiki/ITK_FAQ" target="_blank">http://www.itk.org/Wiki/ITK_FAQ</a><br>
><br>
> Follow this link to subscribe/unsubscribe:<br>
> <a href="http://www.itk.org/mailman/listinfo/insight-users" target="_blank">http://www.itk.org/mailman/listinfo/insight-users</a><br>
><br>
><br>
</div></div></blockquote></div><br>
</div></div></blockquote></div><br>