[Rtk-users] unable to use CudaParkerShortScanImageFilter in python

Mon Jul 10 21:13:50 CEST 2023

Hi Simon,

Cuda PSSF is working fine with the wheel file you shared. But it looks like
something I missed from the reconstruction process and the result is not up
to mark as I expected.
The result looks like as below:
[image: image.png]
[image: image.png]
and the expected is more or less something like below:

[image: image.png]

The input projections and the parameter file (Which I used to fill the
input parameters) are available here:
https://drive.google.com/drive/folders/13EC2Lusyrx3HvkKJL80aSBHfled0FduT?usp=sharing

Can you please check and provide your valuable suggestions on how I can
improve the result.

Here is the updated code:

# Defines the image type
ImageType = itk.Image[itk.F,3]
GPUImageType = itk.CudaImage[itk.F,3]

# Defines the RTK geometry object
geometry = rtk.ThreeDCircularProjectionGeometry.New()
numberOfProjections = 512
firstAngle = -0.54853
angularArc = 270.66079712
sid = 537.09972 # source to isocenter distance
sdd = 908.00674 # source to detector distance
size_x = 1024
size_y = 1024
projOffsetX = -2.06667995
projOffsetY = -0.19235000
outOfPlaneAngle = 0
inPlaneAngle = 0
sourceOffsetX = 0
sourceOffsetY = 0

angle_path = os.path.join("/skull/","angle.bin")
angles2 = []
with open(angle_path, 'rb') as file:
    while True:
        data = file.read(4)  # Read 4 bytes (32 bits) at a time
        if not data:
            break  # Reached end of file
        value = struct.unpack('f', data)[0]  # Unpack binary data as a float
        angles2.append(value)
print("Read angles from bin completed")

for x in angles2:
    geometry.AddProjection(sid,sdd,x,projOffsetX, projOffsetY,
outOfPlaneAngle, inPlaneAngle, sourceOffsetX, sourceOffsetY)

tiffio = itk.TIFFImageIO.New()

extension = ".tif"
fileNameFormat = "/binary/raw" + "%04d" + extension
fnames = itk.NumericSeriesFileNames.New()
fnames.SetStartIndex(0)
fnames.SetEndIndex(numberOfProjections-1)
fnames.SetIncrementIndex(1)
fnames.SetSeriesFormat(fileNameFormat)

ProjectionsReaderType = rtk.ProjectionsReader[ImageType]
projectionsSource = ProjectionsReaderType.New()
projectionsSource.SetImageIO(tiffio)
projectionsSource.SetFileNames(fnames.GetFileNames())
projOrigin = [ -0.42*(size_x-1)/2, 0, -0.42*(size_y-1)/2 ] #input images
are size_x x size_y pixels with a 0.42mm pixel size
projSpacing = [ 0.42, 1, 0.42 ]
projectionsSource.SetOrigin( projOrigin )
projectionsSource.SetSpacing( projSpacing )

# Graft the projections to an itk::CudaImage
projections = GPUImageType.New()
projectionsSource.Update()
projections.SetPixelContainer(projectionsSource.GetOutput().GetPixelContainer())
projections.CopyInformation(projectionsSource.GetOutput())
projections.SetBufferedRegion(projectionsSource.GetOutput().GetBufferedRegion())
projections.SetRequestedRegion(projectionsSource.GetOutput().GetRequestedRegion())

ParkerShortScanImageFilterType = rtk.CudaParkerShortScanImageFilter
ParkerFilter = ParkerShortScanImageFilterType.New()
ParkerFilter.SetInput(projections)
ParkerFilter.SetGeometry(geometry)
ParkerFilter.InPlaceOff()

ConstantImageSourceType = rtk.ConstantImageSource[GPUImageType]
# Create reconstructed image
constantImageSource2 = ConstantImageSourceType.New()
render_x = 512
render_y = 512
sizeOutput = [ render_x, render_y, numberOfProjections ]
origin = [ -(render_x/2.0), 0 , -(render_y/2.0) ]
spacing = [ 1.0, 0.5, 1.0 ]
constantImageSource2.SetOrigin( origin )
constantImageSource2.SetSpacing( spacing )
constantImageSource2.SetSize( sizeOutput )
constantImageSource2.SetConstant(0.)

#Displaced detector weighting
ddf = rtk.CudaDisplacedDetectorImageFilter.New()
ddf.SetInput(ParkerFilter.GetOutput())
ddf.SetGeometry(geometry)

# FDK reconstruction
print("Reconstructing...")
FDKGPUType = rtk.CudaFDKConeBeamReconstructionFilter
feldkamp = FDKGPUType.New()
feldkamp.SetInput(0, constantImageSource2.GetOutput())
feldkamp.SetInput(1, ddf.GetOutput())
feldkamp.SetGeometry(geometry)
feldkamp.GetRampFilter().SetTruncationCorrection(0.0)
feldkamp.GetRampFilter().SetHannCutFrequency(0.0)

# Field-of-view masking
FOVFilterType = rtk.FieldOfViewImageFilter[ImageType, ImageType]
fieldofview = FOVFilterType.New()
fieldofview.SetInput(0, feldkamp.GetOutput())
fieldofview.SetProjectionsStack(feldkamp.GetOutput())
fieldofview.SetGeometry(geometry)

# Writer
print("Writing output image...")
WriterType = rtk.ImageFileWriter[ImageType]
writer = WriterType.New()
writer.SetFileName("image_output_skull00_ver1.3.tiff")
writer.SetInput(fieldofview.GetOutput())
writer.Update()

Regards,
Akshara

On Thu, 6 Jul 2023 at 18:23, Akshara P K <akshara at advitech.in> wrote:

> Hi Simon,
>
> Thank you for the quick response.
> It works perfectly fine except some warnings as you mentioned.
> Thank you again for the support.
>
> Regards,
> Akshara
>
> On Wed, 5 Jul 2023 at 18:36, Simon Rit <simon.rit at creatis.insa-lyon.fr>
> wrote:
>
>> Hi,
>> Following a recent similar request (see here
>> <https://www.creatis.insa-lyon.fr/pipermail/rtk-users/2023-June/001979.html>),
>> I have added CudaParkerShortScanImageFilter with this PR
>> <https://www.creatis.insa-lyon.fr/pipermail/rtk-users/2023-June/001979.html>.
>> New Python packages are available as artifacts of the CI here
>> <https://github.com/RTKConsortium/RTK/actions/runs/5448263998>. I
>> haven't tested them, let me know whether they work... In any case, I have
>> noticed new warnings when using these wheels which must be fixed (see issue
>> #553 <https://github.com/RTKConsortium/RTK/issues/553>).
>> Simon
>>
>> On Wed, Jul 5, 2023 at 2:20 PM Akshara P K <akshara at advitech.in> wrote:
>>
>>> Hi all,
>>>
>>> I was trying to use FirstCudaReconstruction.py to reconstruct a set of
>>> tiff projections images using the python binary of itk-rtk. Since the angle
>>> coverage is around 268 degrees, I had to implement a
>>> ParkerShortScanImageFilter between the output of the ProjectionsReader and
>>> the FDK reconstruction and it is working fine with the CPU version of API.
>>> However, the use of CudaParkerShortScanImageFilter causes an error saying
>>> module 'RTK' has no attribute 'CudaParkerShortScanImageFilter' even though
>>> CudaParkerShortScanImageFilter is added and implemented in the cpp version.
>>> Can anyone provide an insight on how to solve this issue?
>>>
>>> Regards,
>>> Akshara
>>> _______________________________________________
>>> Rtk-users mailing list
>>> rtk-users at openrtk.org
>>> https://www.creatis.insa-lyon.fr/mailman/listinfo/rtk-users
>>>
>>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.creatis.insa-lyon.fr/pipermail/rtk-users/attachments/20230711/6aa75a7a/attachment-0001.htm>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: image.png
Type: image/png
Size: 185685 bytes
Desc: not available
URL: <http://www.creatis.insa-lyon.fr/pipermail/rtk-users/attachments/20230711/6aa75a7a/attachment-0003.png>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: image.png
Type: image/png
Size: 10665 bytes
Desc: not available
URL: <http://www.creatis.insa-lyon.fr/pipermail/rtk-users/attachments/20230711/6aa75a7a/attachment-0004.png>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: image.png
Type: image/png
Size: 187374 bytes
Desc: not available
URL: <http://www.creatis.insa-lyon.fr/pipermail/rtk-users/attachments/20230711/6aa75a7a/attachment-0005.png>