[ITK] [ITK-users] Extract slice from ResampleFilter

[Dženan Zukić]

Hi Sidharta,

the older versions of ITK-SNAP had index-based overlay, ignoring different
spatial position which stem from different origin or direction vectors. I
don't know whether the latest (3.6) breaks with that.

Slicer always uses spatial position and orientation, so for your
experiments you should use Slicer to verify whether your code does what you
want.

This difference in size setting should not matter, I think your calculation
of the direction vectors is more likely cause of error. You could check
that by setting different sizes in between the two extremes you trying
right now. In case there is some bug in ITK (unlikely), you will flush it
out with this. More likely this will make it more apparent to you what you
are doing wrong.

Regards,
Dženan Zukić, PhD, Senior R&D Engineer, Kitware (Carrboro, N.C.)

On Thu, May 11, 2017 at 12:00 PM, sidharta <sidharta.gupta93 at gmail.com>
wrote:

> Dear Dženan,
>
>
> Thank you for your swift reply.
>
> This is the entire function which I used:
>
> int ResampleImage(std::vector<Vector3d> trajectory_points,
> std::vector<Vector3d> coordinate_system, InputImageType::Pointer image){
>
>         Vector3d traj_start, traj_end, third_point;
>         traj_start = trajectory_points.at(0);
>         traj_end = trajectory_points.back();
>
>         third_point = traj_end + coordinate_system.at(3) * 2 /
> image->GetSpacing()[0];
>
>         Vector3d ST, SP, normal_to_plane;
>         ST = traj_end - traj_start;
>         SP = third_point - traj_start;
>         normal_to_plane = ST.cross(SP);
>
>         InputImageType::DirectionType imageDirection, trajectoryDirection;
>         imageDirection = image->GetDirection();
>         std::cout << "Image direction cosines : \n" << imageDirection <<
> std::endl;
>
>         trajectoryDirection[0][0] = ST.normalized()[0];
>         trajectoryDirection[0][1] = ST.normalized()[1];
>         trajectoryDirection[0][2] = ST.normalized()[2];
>
>         Vector3d second = normal_to_plane.cross(ST);
>         Vector3d normalized_second = second.normalized();
>         trajectoryDirection[1][0] = normalized_second[0];
>         trajectoryDirection[1][1] = normalized_second[1];
>         trajectoryDirection[1][2] = normalized_second[2];
>
>         Vector3d normalized_normal = normal_to_plane.normalized();
>         trajectoryDirection[2][0] = normalized_normal[0];
>         trajectoryDirection[2][1] = normalized_normal[1];
>         trajectoryDirection[2][2] = normalized_normal[2];
>
>         std::cout << "Trajectory direction cosines : \n" <<
> trajectoryDirection <<
> std::endl;
>
>         InputImageType::IndexType voxelIndex;
>         InputImageType::PointType point;
>
>         point[0] = traj_start[0], point[1] = traj_start[1], point[2] =
> traj_start[2];
>         image->TransformPhysicalPointToIndex(point, voxelIndex);
>         std::cout << "Index for point " << point << " is " << voxelIndex <<
> std::endl;
>
>         point[0] = traj_end[0], point[1] = traj_end[1], point[2] =
> traj_end[2];
>         image->TransformPhysicalPointToIndex(point, voxelIndex);
>         std::cout << "Index for point " << point << " is " << voxelIndex <<
> std::endl;
>
>         double mpr_start[3] = { traj_start[0], traj_start[1],
> traj_start[2] };
>
>         typedef itk::ResampleImageFilter<InputImageType, InputImageType>
> resampleFilterType;
>         resampleFilterType::Pointer resampleFilter =
> resampleFilterType::New();
>
>         typedef itk::LinearInterpolateImageFunction<InputImageType,
> InputPixelType>
> interpolatorType;
>         interpolatorType::Pointer linearInterpolator =
> interpolatorType::New();
>
>         resampleFilter->SetInterpolator(linearInterpolator);
>         resampleFilter->SetDefaultPixelValue(-2000);
>         resampleFilter->SetOutputSpacing(image->GetSpacing());
>         resampleFilter->SetOutputOrigin(mpr_start);
>         resampleFilter->SetOutputDirection(trajectoryDirection);
>         InputImageType::SizeType size;
>         size[0] = trajectory_points.size();
>         size[1] = 1.0 / image->GetSpacing()[1];
>         size[2] = 1.0 / image->GetSpacing()[2];
>         resampleFilter->SetSize(size); // This is for the trajectory, this
> corresponds to a size of [137, 8, 8] which is essentially the volume I want
> to extract.
>         //resampleFilter->SetSize(image->GetLargestPossibleRegion().GetSize());
> //
> This is for the entire resampled image
>         resampleFilter->SetInput(image);
>         //resampleFilter->Update();
>
>         //InputImageType::Pointer resampled_image = InputImageType::New();
>         //resampled_image = resampleFilter->GetOutput();
>
>         FileWriter(resampleFilter->GetOutput(),
> "C:/Users/api/Desktop/resampled_volume_1.mhd");
>         return EXIT_SUCCESS;
>
> }
>
> I should also add that when I open the original image with the resampled
> image in ITK-SNAP (by adding another image), the overlay is exactly the
> same. Additionally, when I opened this in Slicer, only the top part of the
> volume seems to be rotated and the rest isn't. Kindly let me know if you
> need more information.
> I apologize for asking too many questions. We don't have an experienced ITK
> user in our lab, although vector understanding is sound :)
>
> Thank you for your time.
>
>
>
> --
> View this message in context: http://itk-users.7.n7.nabble.
> com/Extract-slice-from-ResampleFilter-tp38074p38230.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
>
> Kitware offers ITK Training Courses, for more information visit:
> http://www.kitware.com/products/protraining.php
>
> 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://public.kitware.com/mailman/listinfo/insight-users
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://public.kitware.com/pipermail/community/attachments/20170512/385a00d8/attachment.html>
-------------- next part --------------
_____________________________________
Powered by www.kitware.com

Visit other Kitware open-source projects at
http://www.kitware.com/opensource/opensource.html

Kitware offers ITK Training Courses, for more information visit:
http://www.kitware.com/products/protraining.php

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://public.kitware.com/mailman/listinfo/insight-users