[Rtk-users] Understanding projection matrices management in rtkCudaBackProjectionImageFilter

[Simon Rit]

Hi Thibault,
You started digging in the code so that gets complicated. First thing
first, yes, the link you provide is still valid.
The offset in X is not random but computed in
BackProjectionImageFilter::GenerateInputRequestedRegion
<https://github.com/SimonRit/RTK/blob/master/code/rtkBackProjectionImageFilter.txx>.
The idea is that we request only the part of the projections we need to do
the backprojection. This optimizes the streaming (we only read from disk
what we need). If you use a huge volume (the test I'd suggest to validate
that I'm correct and there is nothing wrong with your code), then you'll
need the full projection and you will always get 0 in X as in Y. If you use
a small volume and do a computation projection per projection, you will get
a varying index depending on your input volume and the geometry for that
projection.
I think you understood the rest pretty well. Does that make sense or am I
missing something else?
Simon

On Tue, Nov 10, 2015 at 1:39 PM, Notargiacomo Thibault <gnthibault at gmail.com
> wrote:

> Hi all,
>
> I am currently trying to use the "AddProjectionInRadians" method from the
> ThreeDCircularProjectionGeometry class in RTK to import my geometry.
>
> First things first, I would like to know if the information on this pdf
> are still valid: http://www.openrtk.org/Doxygen/geometry.pdf
>
> Because I used them to deduce the source offset and Projection offset from
> my own matrix decomposition.
>
> Notice that I don't want to use geometric informations coming from the
> projection images themselves, so their geometric origin is assumed to be
> (0,0,0) in my case.
>
> In order to check if my geometry import method is correct, I print out one
> of the intermediate form of the projection matrix used in the
> backprojection process, in order to check if it is the same as in my own
> implementation.
>
> To do so, I print out one of the intermediate matrix that is generated in
> rtkCudaBackProjectionImageFilter while launching a SART reconstruction for
> instance.
>
> I was surpised by the fact that, the matrix allowing the transformation
> from 2D geometric world to 2D image index world, namely matrixIdxProj in
> rtkCudaBackProjectionImageFilte.cxx
> <http://www.openrtk.org/Doxygen/rtkCudaBackProjectionImageFilter_8cxx_source.html>
> is not always the same along the reconstruction process ?
>
> Here is the line that generates the matrix I don't understand:
>
> for(unsigned int i=0; i<2; i++)
>  //SR: 0.5 for 2D texture
>   matrixIdxProj[i][2] =
> -1*(this->GetInput(1)->GetBufferedRegion().GetIndex()[i])+0.5;
>
> What I understand here, is that, we are generating a 2D translation matrix
> embedded in a 3*3 matrix that allows to take into account two things
> related to the 2D projection:
>
> -If the internal behaviour of the imageContainer, allows buffering of only
> a part of the image, beginning at index (i,j), this information will be
> embeded in the projection matrix
>
> -the ITK floating point addressing definition seems to differ from the one
> from cuda (and openCV):
> When I take a look at this page
> <http://docs.nvidia.com/cuda/cuda-c-programming-guide/index.html#texture-fetching>
> I understand that the floating point indexing (0.5,0.5) in Cuda points to
> the center of the pixel (0,0) hence a linear contribution of 100%, while it
> seems that we must add a translation of (0.5,0.5) to the ITK geometric
> framework in order to get the same behaviour.
>
> In practice matrixIdxProj is based on a 3*3 identity matrix, where the
> inner 2D translation (ie 2 first term of the last column) in my case are:
> -in X, a random number between 10 and 150 plus 0.5
> -in Y: always 0.5
>
> I understand the 0.5 part, but why is there always a 0 offset in Y but a
> random offset in X remains a mystery to me.
>
> At first, I believed that each projection was splitted into multiple part,
> each having a different X and Y offset.
> But reconstructing using only one projection generates only one
> backprojection operation, with a non zero offset in X, 123 in my case, for
> an image of size 780 in X.
>
> Do you think there is something wrong whith the values of matrixIdxProjI
> am experiencing ?
>
> Thank you in advance.
>
> Kind regards
>
> Thibault Notargiacomo
>
> _______________________________________________
> Rtk-users mailing list
> Rtk-users at public.kitware.com
> http://public.kitware.com/mailman/listinfo/rtk-users
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://public.kitware.com/pipermail/rtk-users/attachments/20151110/8b9b86e6/attachment-0010.html>