<div dir="ltr"><div>Dear Laura,</div><div>I don't remember what was NormalizedJoseph but it's been removed a while ago:</div><div><a href="https://github.com/RTKConsortium/RTK/pull/171/commits/892e5cc5c765a13a7205ad2cabeb6193a75631b5">https://github.com/RTKConsortium/RTK/pull/171/commits/892e5cc5c765a13a7205ad2cabeb6193a75631b5</a></div><div>To answer your questions:</div><div>1/ Possibly. The adjoint of CUDAForwardProjectionImageFilter is CudaRayCastBackProjectionImageFilter. The adjoint of JosephBackProjectionImageFilter is JosephFowardProjectionImageFilter.</div><div>2/ Good!<br></div><div>3/ This is not the expected behavior and I have opened an issue for this:</div><div><a href="https://github.com/RTKConsortium/RTK/issues/670">https://github.com/RTKConsortium/RTK/issues/670</a></div><div>It would be great if you have time to check what's going on. The forward projection code is here:</div><div><a href="https://github.com/RTKConsortium/RTK/blob/master/src/rtkCudaForwardProjectionImageFilter.cu#L61-L153">https://github.com/RTKConsortium/RTK/blob/master/src/rtkCudaForwardProjectionImageFilter.cu#L61-L153</a></div><div>The backprojection code is here:</div><div><a href="https://github.com/RTKConsortium/RTK/blob/master/src/rtkCudaRayCastBackProjectionImageFilter.cu#L63-L238">https://github.com/RTKConsortium/RTK/blob/master/src/rtkCudaRayCastBackProjectionImageFilter.cu#L63-L238</a></div><div>The trickiest part is that the wired interpolation in the forward projection has to be transposed to "splat" the values.</div><div>Let us know if you spot something wrong. We'll try to investigate if you don't manage to do it.</div><div>Thanks for the report in any case!</div><div>Simon<br></div><div><br></div></div><br><div class="gmail_quote gmail_quote_container"><div dir="ltr" class="gmail_attr">On Thu, Jan 9, 2025 at 1:19 PM Laura Hellwege <<a href="mailto:l.hellwege@uni-luebeck.de">l.hellwege@uni-luebeck.de</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div class="msg2340827923785718647">
<div lang="EN-GB" style="overflow-wrap: break-word;">
<div class="m_2340827923785718647WordSection1">
<p class="MsoNormal"><span style="font-size:11pt">Hello!<u></u><u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt"><u></u> <u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt">I have some questions related to this post from 2015:<u></u><u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt"><a href="https://public.kitware.com/pipermail/rtk-users/2015-February/009750.html" target="_blank">https://public.kitware.com/pipermail/rtk-users/2015-February/009750.html</a><u></u><u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt"><u></u> <u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt">First, the referenced website “<a href="http://www.openrtk.org/Doxygen/classrtk_1_1NormalizedJosephBackProjectionImageFilter.html" target="_blank">http://www.openrtk.org/Doxygen/classrtk_1_1NormalizedJosephBackProjectionImageFilter.html</a>”
seems to not exist anymore. <u></u><u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt"><u></u> <u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt">Second, I have some issues concerning backprojection and adjointness:<u></u><u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt">For the following I am using these parameters:
<u></u><u></u></span></p>
<p class="MsoNormal" style="margin-left:36pt"><span style="font-size:11pt">sid = 500<u></u><u></u></span></p>
<p class="MsoNormal" style="margin-left:36pt"><span style="font-size:11pt">sdd = 1000<u></u><u></u></span></p>
<p class="MsoNormal" style="margin-left:36pt"><span style="font-size:11pt">n_projections = 4<u></u><u></u></span></p>
<p class="MsoNormal" style="margin-left:36pt"><span style="font-size:11pt">angular_increment = 90.0<u></u><u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt">I am using the python version of ITK and tried forward projection (CUDAForwardProjectionImageFilter) followed by backprojection (JosephBackProjectionImageFilter/CudaBackProjectionImageFilter/CudaRayCastBackProjectionImageFilter)
of an image with a single pixel value set to one:<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">img = np.zeros(image_size)<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">img[image_size[0] // 2, image_size[1] // 2, image_size[2] // 2] = 1.0<u></u><u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt"><u></u> <u></u></span></p>
<ol style="margin-top:0cm" start="1" type="1">
<li class="m_2340827923785718647MsoListParagraph" style="margin-left:0cm"><span style="font-size:11pt">I do not get expected results for the 2D case:<u></u><u></u></span></li></ol>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">image_size = [3, 1, 3]<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">pixel_size_mm = 1.0<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">slice_thickness_mm = 1.0<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt"><u></u> <u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">projection_size = [5, 1]<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">projection_spacing = [1.0, 1.0]<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt"><u></u> <u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">After forward projection, the backprojected image by JosephBackProjectionImageFilter is all
<b>zeros</b>. Is this expected behavior?<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt"> <u></u><u></u></span></p>
<ol style="margin-top:0cm" start="2" type="1">
<li class="m_2340827923785718647MsoListParagraph" style="margin-left:0cm"><span style="font-size:11pt">According to the definition of adjoint operator, the following should hold:
<u></u><u></u></span></li></ol>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt"><Ax, y> = <x, adjoint(A)y> (with <.,.> indicating the scalar product)<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt"><u></u> <u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">This works for 3D <b>CPU</b> case (JosephBackProjectionImageFilter):<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">image_size = [3, 3, 3]<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">pixel_size_mm = 1.0<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">slice_thickness_mm = 1.0<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt"><u></u> <u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">projection_size = [5, 5]<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">projection_spacing = [1.0, 1.0]<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt"><u></u> <u></u></span></p>
<p class="MsoNormal" style="margin-left:36pt"><span style="font-size:11pt"><Ax, y> = <fp, fp> = 9.0<u></u><u></u></span></p>
<p class="MsoNormal" style="margin-left:36pt"><span style="font-size:11pt"><x, adjoint(A)y> = <img, bp> = 9.000001907348633<u></u><u></u></span></p>
<p class="MsoNormal" style="margin-left:36pt"><span style="font-size:11pt">This discrepancy is within tolerance.<u></u><u></u></span></p>
<p class="MsoNormal" style="margin-left:36pt"><span style="font-size:11pt"><u></u> <u></u></span></p>
<ol style="margin-top:0cm" start="3" type="1">
<li class="m_2340827923785718647MsoListParagraph" style="margin-left:0cm"><span style="font-size:11pt">However, in the 3D
<b>GPU</b> case:<u></u><u></u></span></li></ol>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt"><x, adjoint(A)y> = <img, bp> = 4.0 for CudaBackProjectionImageFilter<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt">and<u></u><u></u></span></p>
<p class="m_2340827923785718647MsoListParagraph"><span style="font-size:11pt"><x, adjoint(A)y> = <img, bp> = 0.33841705322265625 for CudaRayCastBackProjectionImageFilter<u></u><u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt"><u></u> <u></u></span></p>
<p class="MsoNormal" style="margin-left:36pt"><span style="font-size:11pt">What kind of filtering or normalization has been used on the projections/images? I could not verify adjointment with the aforementioned method of backprojecting a sinogram of ones
and applying division.<u></u><u></u></span></p>
<p class="MsoNormal" style="margin-left:36pt"><span style="font-size:11pt">Is there a GPU backprojection filter for which adjointness holds?<u></u><u></u></span></p>
<p class="MsoNormal" style="margin-left:36pt"><span style="font-size:11pt"><u></u> <u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt">Thank you in advance for your reply and suggestions.<u></u><u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt">Best regards, Laura<u></u><u></u></span></p>
<p class="MsoNormal"><span style="font-size:11pt"><u></u> <u></u></span></p>
</div>
</div>
_______________________________________________<br>
Rtk-users mailing list<br>
<a href="mailto:rtk-users@openrtk.org" target="_blank">rtk-users@openrtk.org</a><br>
<a href="https://www.creatis.insa-lyon.fr/mailman/listinfo/rtk-users" rel="noreferrer" target="_blank">https://www.creatis.insa-lyon.fr/mailman/listinfo/rtk-users</a><br>
</div></blockquote></div>