[Insight-users] fastest way to get image voxels(pixels) within an itkMesh ?
gregthom
gregthom99 at yahoo.com
Sat Jul 4 17:13:56 EDT 2009
Hi Luis,
Thank you for taking some time to look into this.
Yes you describe the problem accurately.
>> Can we assume that (B) and (C) share the same image grid ?
Answer: not necessarily.
>>It is not clear from your email,
>>if what you want at the end is an image in the grid of (B) with
>>the values of (C), or whether you want to simply have a list
>>(e.g. std::vector<> with the list of values from (C).
I need the in the end an image in the grid of (B) with the values of C as
you thought.
Can you suggest how to obtain this efficiently ? I will start to work on the
algorithm in the
pseudocode, in fact I had already been doing something similar.
Hope to hear from you
best wishes
Luis Ibanez wrote:
>
> Hi Gregthom,
>
> So you have:
>
>
> A) a Mesh
> B) a binary image that you get from rasterizing (A)
> C) a grayscale image
>
>
> You want to get the grayscales values of (C) that correspond
> to positions of (B) that are inside of (A).
>
> Can we assume that (B) and (C) share the same image grid ?
>
> That is,
> Do the images B and C have the same origin, spacing, orientation
> and number of pixels ?
>
> If yes,
> then you could do this faster by using the MaskImageFilter,
> and passing (B) and (C) as its inputs.
>
> If no,
> You will have to use code similar to what is used in the
> ImageMetrics, when collecting the values of the moving image.
>
> In pseudo code:
>
> Run Mesh rasterization filter
>
> Create interpolator and connect it to (C)
>
> For all pixels in the mask (B)
> if pixel is ON
> convert index to physical point
> Evaluate interpolator at point
>
>
>
>
> It is not clear from your email,
> if what you want at the end is an image in the grid of (B) with
> the values of (C), or whether you want to simply have a list
> (e.g. std::vector<> with the list of values from (C).
>
>
> Please clarify what kind of final output you want/need.
>
>
> Thanks
>
>
> Luis
>
>
>
>
> ----------------
> gregthom wrote:
>> Hello all,
>>
>> Since my last post, I have been able to concoct a piece of code that does
>> what I originally wanted. i'e
>>
>> take a mesh, binarize it, and then give me all voxels(physical) within my
>> binary volume.
>> I want to use these physical coordinates to interpolate in another 3D
>> volume
>> with scalar values, i.e for all voxels in my binary volume, what is the
>> scalar value interpolated from the 3D volume ? My main worry is about
>> speed.
>> For the first part, is this really the fastest way ? And what about the
>> interpolation part ? how to do that ?
>>
>> Thank you
>>
>>
>> Here is my working code
>> **********************************************************************************************************************************************
>> typedef unsigned char BinaryPixelType;
>> typedef unsigned short PixelType;
>> //typedef float PixelType;
>> typedef itk::Image<BinaryPixelType, 3> BinaryImageType;
>> typedef itk::Image<PixelType, 3> ImageType;
>> typedef itk::TriangleMeshToBinaryImageFilter<MeshType, BinaryImageType>
>> TriangleMeshToBinaryImageFilterType;
>> TriangleMeshToBinaryImageFilterType::Pointer triangleMeshToImage =
>> TriangleMeshToBinaryImageFilterType::New();
>> typedef itk::MeshSpatialObject<MeshType> MeshSpatialObjectType;
>>
>> //convert matlab meshS (i.e meshS.vertices, meshS.triangles) to itkmesh
>> MeshType::Pointer mymesh = meshSstructToItkMesh(prhs[0]);
>> triangleMeshToImage->SetInput(mymesh);
>> BinaryImageType::SizeType size;
>> size[0] = (unsigned int)((double *)(mxGetPr(prhs[3])))[0];
>> size[1] = (unsigned int)((double *)(mxGetPr(prhs[3])))[1];
>> size[2] = (unsigned int)((double *)(mxGetPr(prhs[3])))[2];
>>
>>
>> triangleMeshToImage->SetSpacing(spc);
>> triangleMeshToImage->SetOrigin(org);
>> triangleMeshToImage->SetSize (size);
>>
>> //typedef itk::ImageFileWriter<BinaryImageType> WriterType;
>> //WriterType::Pointer writer = WriterType::New();
>> //writer->SetFileName("testMeshBinaryImage.mha");
>> //writer->SetInput(triangleMeshToImage->GetOutput());
>>
>> try
>> {
>> triangleMeshToImage->Update();
>> //writer->Update();
>> }
>> catch( itk::ExceptionObject & excep )
>> {
>> std::cerr << "Exception caught !" << std::endl;
>> std::cerr << excep << std::endl;
>> }
>>
>> //create meshtoimage filter
>> BinaryImageType::ValueType m_MaskValue =
>> itk::NumericTraits<BinaryImageType::ValueType>::One;
>> typedef itk::ImageRegionConstIteratorWithIndex<BinaryImageType>
>> MaskIteratorType;
>> MaskIteratorType mit( triangleMeshToImage->GetOutput(),
>> triangleMeshToImage->GetOutput()->GetBufferedRegion() );
>> BinaryImageType::IndexType index;
>> typedef itk::PointSet< double, 3 > PointSetType;
>> typedef PointSetType::PointType InputPointType;
>>
>> // lets get number of voxels
>> int nvoxels;
>> for (mit.GoToBegin(); !mit.IsAtEnd(); ++mit)
>> {
>> index = mit.GetIndex();
>> InputPointType iPoint;
>> triangleMeshToImage->GetOutput()->TransformIndexToPhysicalPoint(
>> index,
>> iPoint );
>> if (mit.Get() == m_MaskValue)
>> {
>> nvoxels++;
>> }
>> else
>> {
>>
>> }
>> }//
>>
>> //create output variable for result
>> mwSize dims[] = {nvoxels,3};
>> plhs[0] = mxCreateNumericArray( 2, dims, mxDOUBLE_CLASS, mxREAL
>> );
>> double *pr = mxGetPr( plhs[0] );
>>
>> // second pass thru mask to get only values within the mask
>> unsigned int tcount;
>> for (mit.GoToBegin(); !mit.IsAtEnd(); ++mit)
>> {
>> index = mit.GetIndex();
>> InputPointType iPoint;
>> triangleMeshToImage->GetOutput()->TransformIndexToPhysicalPoint(
>> index,
>> iPoint );
>> if (mit.Get() == m_MaskValue)
>> {
>> //mexPrintf("Point in mesh: [%d %f %f
>> %f]\n",tcount,iPoint[0],iPoint[1],iPoint[2]);
>> pr[0*nvoxels + tcount] = (double)iPoint[0];
>> pr[1*nvoxels + tcount] = (double)iPoint[1];
>> pr[2*nvoxels + tcount] = (double)iPoint[2];
>>
>> tcount++;
>> }
>> else
>> {
>>
>> }
>> }//
>>
>>
>>
>>
>>
>>
>> Dan Mueller-2 wrote:
>>
>>>Hi Greg,
>>>
>>>You could try rasterizing your mesh, then iterating over the image
>>>considering only values under the raster mask. Further, you could
>>>constrict the iteration to the bounding box of the raster mask. If
>>>your mesh consists of triangles, you can use
>>>itkTriangleMeshToImageFilter to create a binary image mask.
>>>
>>>If rasterizing the mesh is already too slow, I'm not sure what else
>>>you could do...
>>>
>>>Some more information regarding your specific problem might also help
>>>discover solutions "outside the box".
>>>
>>>Hope this helps.
>>>
>>>Regards, Dan
>>>
>>>2009/6/25 gregthom <gregthom99 at yahoo.com>:
>>>
>>>>Greetings
>>>>
>>>>This is more of a design help question. I have an itkMesh and and
>>>>itkImage
>>>>in the same frame of reference. I need to get the image values of all
>>>>image
>>>>voxels(pixels) bounded by the itkMesh. Anybody have any thoughts on how
>>>>to
>>>>achieve this ? What is the fastest way to achieve this ?
>>>>
>>>>Thank you
>>>>
>>>>
>>>>greg
>>>>--
>>>>View this message in context:
>>>>http://www.nabble.com/fastest-way-to-get-image-voxels%28pixels%29-within-an-itkMesh---tp24199751p24199751.html
>>>>Sent from the ITK - Users mailing list archive at Nabble.com.
>>>>
>>>>_____________________________________
>>>>Powered by www.kitware.com
>>>>
>>>>Visit other Kitware open-source projects at
>>>>http://www.kitware.com/opensource/opensource.html
>>>>
>>>>Please keep messages on-topic and check the ITK FAQ at:
>>>>http://www.itk.org/Wiki/ITK_FAQ
>>>>
>>>>Follow this link to subscribe/unsubscribe:
>>>>http://www.itk.org/mailman/listinfo/insight-users
>>>>
>>>
>>>_____________________________________
>>>Powered by www.kitware.com
>>>
>>>Visit other Kitware open-source projects at
>>>http://www.kitware.com/opensource/opensource.html
>>>
>>>Please keep messages on-topic and check the ITK FAQ at:
>>>http://www.itk.org/Wiki/ITK_FAQ
>>>
>>>Follow this link to subscribe/unsubscribe:
>>>http://www.itk.org/mailman/listinfo/insight-users
>>>
>>>
>>
>>
> _____________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at
> http://www.kitware.com/opensource/opensource.html
>
> Please keep messages on-topic and check the ITK FAQ at:
> http://www.itk.org/Wiki/ITK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://www.itk.org/mailman/listinfo/insight-users
>
>
--
View this message in context: http://www.nabble.com/fastest-way-to-get-image-voxels%28pixels%29-within-an-itkMesh---tp24199751p24338074.html
Sent from the ITK - Users mailing list archive at Nabble.com.
More information about the Insight-users
mailing list