<div dir="ltr"><div>Hi,</div><div>Both your input volume and your stack of projections should be 3D. If you simulate only one projection, you can always make it 2D afterwards by slicing the stack. The third dimension of spacing and origin of the stack are never used in RTK but must be provided since it's a 3D image.<br></div><div>I might be wrong but I don't think SimpleITK and the RTK wrappings are compatible. I would suggest the use of the ITK wrappings described in their examples. See <a href="https://itkpythonpackage.readthedocs.io/en/master/Quick_start_guide.html">here</a> for an example.</div><div>Enclosed is a slightly modified version of your script which produces a projection stack. The geometry is probably not correct but I let you play with the script to set them as you wish.<br></div><div>Cheers,</div><div>Simon<br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Mon, Aug 9, 2021 at 9:55 AM Milou de Boer <<a href="mailto:m.m.deboer@student.utwente.nl">m.m.deboer@student.utwente.nl</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 style="overflow-wrap: break-word;" lang="NL"><div class="gmail-m_3052572566543848943WordSection1"><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal">Dear all,</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal">I’m quite new to Python (and the use of RTK and ITK), but I have a project for which I want to create digitally reconstructed radiographs (DRRs) of a CT. I therefore have successfully installed RTK for python and I have run the FirstReconstruction.py example. Now I want to edit this so that it is possible to create DRRs, thus forwardprojection. </p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal">I have read in the archives of RTK-users and people are mostly directed towards the C++ code <a href="https://github.com/SimonRit/RTK/blob/b32cffdc6f9d7a432c50023c370ed996a7985b69/applications/rtkforwardprojections/rtkforwardprojections.cxx" target="_blank">https://github.com/SimonRit/RTK/blob/b32cffdc6f9d7a432c50023c370ed996a7985b69/applications/rtkforwardprojections/rtkforwardprojections.cxx</a> for forward projections. However, I find this very hard to read since I don’t know C++ so I hope you can help. I have also read <a href="https://discourse.itk.org/t/itk-to-rtk-migration-forward-projection-issues-questions/2107" target="_blank">https://discourse.itk.org/t/itk-to-rtk-migration-forward-projection-issues-questions/2107</a> regarding this matter, but I cannot reproduce this.</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal">What I do in my code:</p><ul style="margin-top:0cm" type="disc"><li class="gmail-m_3052572566543848943MsoListParagraph" style="margin-left:0cm">Load 3D CT image and change it to float Image Type ( sitk.sitkFloat32)</li><li class="gmail-m_3052572566543848943MsoListParagraph" style="margin-left:0cm">Define the image type I want for DRR</li><li class="gmail-m_3052572566543848943MsoListParagraph" style="margin-left:0cm">Create a stack of empty projection images using ConstantImageSource, using ImageType of DRR</li><li class="gmail-m_3052572566543848943MsoListParagraph" style="margin-left:0cm">Define the RTK geometry object (ThreeDCircularProjectionGeometry)</li><li class="gmail-m_3052572566543848943MsoListParagraph" style="margin-left:0cm">Then create JosephForwardProjectionImageFilter</li></ul><p class="MsoNormal">Here, the problems arise.</p><p class="MsoNormal">The filter only seems to accept 3D imageTypes, so I have given it the [itk.F,3], [itk.F,3] input.</p><p class="MsoNormal">Then, if I try to setInput(0) with constantImageSource.GetOutput(), it throws the error: </p><div><pre style="word-break:break-all"><span><u></u> <u></u></span></pre><pre style="word-break:break-all"><span>Expecting argument of type itkImageF3 or itkImageSourceIF3.<u></u><u></u></span></pre><pre style="word-break:break-all"><span>Additional information:<u></u><u></u></span></pre><pre style="word-break:break-all"><span>Wrong number or type of arguments for overloaded function 'itkImageToImageFilterIF3IF3_SetInput'.<u></u><u></u></span></pre></div><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal">Probably since I have used the ImageType of the DRR (2D) for the constantImageSource. </p><p class="MsoNormal">Additionally, it seems that the setInput(1) = CT also doesn’t work, maybe because I have used sitk.sitkFloat32) to get it to a float image?</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal">Does anyone know what I have to do here to solve this? And/or is there some python script available which does forward projection so I can have a look how it is done in Python?</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal">I also don’t know what I would have to do after this step. Do I just do an FDK reconstruction etc. just as in the FirstReconstruction.py script?</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal">Thank you in advance for anyone who can answer (some) of my questions.</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal">M</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal">My code:</p><p class="MsoNormal">import sys</p><p class="MsoNormal">import itk</p><p class="MsoNormal">from itk import RTK as rtk</p><p class="MsoNormal">import SimpleITK as sitk</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"># Utility method that either downloads data from the Girder repository or</p><p class="MsoNormal"># if already downloaded returns the file name for reading from disk (cached data).</p><p class="MsoNormal">%run update_path_to_download_script</p><p class="MsoNormal">from downloaddata import fetch_data as fdata</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"># Always write output to a separate directory, we don't want to pollute the source directory. </p><p class="MsoNormal">import os</p><p class="MsoNormal">OUTPUT_DIR = 'Output_test'</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"># Loading 3D CT image</p><p class="MsoNormal">CT = sitk.ReadImage(fdata("training_001_ct.mha"), sitk.sitkFloat32)</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"># Defines the image type CT</p><p class="MsoNormal">Dimension_CT = 3</p><p class="MsoNormal">PixelType = itk.F</p><p class="MsoNormal">ImageType_CT = itk.Image[PixelType, Dimension_CT]</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"># Defines the image type DRR</p><p class="MsoNormal">Dimension_DRR = 2</p><p class="MsoNormal">ImageType_DRR = itk.Image[PixelType, Dimension_DRR]</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"># Create a stack of empty projection images</p><p class="MsoNormal">ConstantImageSourceType = rtk.ConstantImageSource[ImageType_DRR]</p><p class="MsoNormal">constantImageSource = ConstantImageSourceType.New()</p><p class="MsoNormal"># Define origin, sizeOutput and spacing (still need to change these)</p><p class="MsoNormal">origin = [ -127, -127]</p><p class="MsoNormal">sizeOutput = [ 512, 512]</p><p class="MsoNormal">spacing = [ 0.653595, 2.0]</p><p class="MsoNormal">constantImageSource.SetOrigin( origin )</p><p class="MsoNormal">constantImageSource.SetSpacing( spacing )</p><p class="MsoNormal">constantImageSource.SetSize( sizeOutput )</p><p class="MsoNormal">constantImageSource.SetConstant(0.)</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"># Defines the RTK geometry object</p><p class="MsoNormal">geometry = rtk.ThreeDCircularProjectionGeometry.New()</p><p class="MsoNormal">numberOfProjections = 360</p><p class="MsoNormal">firstAngle = 0.</p><p class="MsoNormal">angularArc = 360.</p><p class="MsoNormal">sid = 600 # source to isocenter distance</p><p class="MsoNormal">sdd = 1200 # source to detector distance</p><p class="MsoNormal">for x in range(0,numberOfProjections):</p><p class="MsoNormal">  angle = firstAngle + x * angularArc / numberOfProjections</p><p class="MsoNormal">  geometry.AddProjection(sid,sdd,angle)</p><p class="MsoNormal">    </p><p class="MsoNormal"># Writing the geometry to disk</p><p class="MsoNormal">xmlWriter = rtk.ThreeDCircularProjectionGeometryXMLFileWriter.New()</p><p class="MsoNormal">xmlWriter.SetFilename ( sys.argv[2] )</p><p class="MsoNormal">xmlWriter.SetObject ( geometry )</p><p class="MsoNormal">xmlWriter.WriteFile()</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal">REIType = rtk.JosephForwardProjectionImageFilter[ImageType_CT, ImageType_CT]</p><p class="MsoNormal">rei = REIType.New()</p><p class="MsoNormal">semiprincipalaxis = [ 50, 50, 50]</p><p class="MsoNormal">center = [ 0, 0, 10]</p><p class="MsoNormal"># Set GrayScale value, axes, center...</p><p class="MsoNormal">#rei.SetDensity(2)</p><p class="MsoNormal">#rei.SetAngle(0)</p><p class="MsoNormal">#rei.SetCenter(center)</p><p class="MsoNormal">#rei.SetAxis(semiprincipalaxis)</p><p class="MsoNormal">rei.SetGeometry( geometry )</p><p class="MsoNormal">rei.SetInput(0, constantImageSource.GetOutput())</p><p class="MsoNormal">rei.SetInput(1, CT)</p><p class="MsoNormal"><u></u> <u></u></p><p class="MsoNormal"><u></u> <u></u></p></div></div>
_______________________________________________<br>
Rtk-users mailing list<br>
<a href="mailto:Rtk-users@public.kitware.com" target="_blank">Rtk-users@public.kitware.com</a><br>
<a href="https://public.kitware.com/mailman/listinfo/rtk-users" rel="noreferrer" target="_blank">https://public.kitware.com/mailman/listinfo/rtk-users</a><br>
</blockquote></div>