[Rtk-users] Problem with reconstruction on Python

[Simon Rit]

Hi Ruslan,
You're calling GetOrigin and GetSpacing before calling
projectionsSource.Update() which is meaningless. If the values are still
incorrect after calling Update, this is a problem because the geometry
heavily relies on this. You need to set the spacing to the correct value
(possible using ChangeInformationImageFilter
<https://itk.org/Doxygen/html/classitk_1_1ChangeInformationImageFilter.html>).
For the origin, it strongly depends on the offsets, whether they are wrt
the projection center or corner. In the latter case, you can set it to 0.
Is the former case, you need to calculate the offset as
-0.5*(size-1)*spacing (in each coordinate).
You will also want to check if the output of the projectionsSource are line
integrals as expected by RTK reconstruction filter. If not, you need to
fine a way to convert the raw data to line integrals (e.g. taking the log
of the values divided by the value in air).
Best regards,
Simon

On Wed, Dec 6, 2023 at 2:33 PM Джаркинов, Руслан Бекзатович <
r.dzharkinov at g.nsu.ru> wrote:

> Hello for everyone,
>
> I have some problems with reconstruction of real data. If I use the code
> below, then the resulting reconstructions are of very poor quality. It
> seems to me that I don't understand something fundamental.
>
> Some info about my data:
> I have 7201 projections in .tiff format, angle between projections — 0.05.
> Detector pixel size — 85.0 micrometer.
> I also already have a reconstruction of this data (1490x1490x1408) from a
> different company, but I also want to try to reconstruct and compare the
> quality of resulting images.
>
> Here is the basic python code I am trying to use (modified from the
> example), if anyone can see what I might be doing wrong I would appreciate
> it:
>
> (some code of import modules)
> tiffio = itk.TIFFImageIO.New()
> ImageType = itk.Image[itk.F,3]
>
> ProjectionsReaderType = rtk.ProjectionsReader[ImageType]
> projectionsSource = ProjectionsReaderType.New()
> projectionsSource.SetImageIO(tiffio)
> projectionsSource.SetFileNames(data_names)
> projectionsSource.GetOrigin() # output — itkVectorD3 ([1.79769e+308,
> 1.79769e+308, 1.79769e+308])
> projectionsSource.GetSpacing() # output — itkPointD3 ([1.79769e+308,
> 1.79769e+308, 1.79769e+308])
> geometry = rtk.ThreeDCircularProjectionGeometry.New()
>
> projOffsetX=0.0375
> projOffsetY=-27.965
> numberOfProjections = 7201
> firstAngle = 0.
> angularArc = 360.
> sid = 17.2 # source to isocenter distance
> sdd = 795 # source to detector distance
> for x in range(0,numberOfProjections):
> angle = firstAngle + x * angularArc / numberOfProjections
> geometry.AddProjection(sid,sdd,angle, projOffsetX, projOffsetY)
>
> ConstantImageSourceType = rtk.ConstantImageSource[ImageType]
> constantImageSource2 = ConstantImageSourceType.New()
> origin = [-0.003677988*(1490-1)/2, -0.003677988*(1490-1)/2, 0]
> sizeOutput = [1490, 1490, 1410]
> spacing = [ 0.003677988, 0.003677988, 0.003677988] # Image pixel size —
> 3.677988 micrometers, so I convert to mm
> constantImageSource2.SetOrigin(origin)
> constantImageSource2.SetSpacing( spacing )
> constantImageSource2.SetSize( sizeOutput )
> constantImageSource2.SetConstant(0.)
>
> print("Reconstructing...")
> FDKCPUType = rtk.FDKConeBeamReconstructionFilter[ImageType]
> feldkamp = FDKCPUType.New()
> feldkamp.SetInput(0, constantImageSource2.GetOutput())
> feldkamp.SetInput(1, projectionsSource.GetOutput())
> feldkamp.SetGeometry(geometry)
> feldkamp.GetRampFilter().SetTruncationCorrection(0.0)
> feldkamp.GetRampFilter().SetHannCutFrequency(0.0)
>
> feldkamp.Update()
> output = feldkamp.GetOutput()
>
> print("Writing output image...")
> WriterType = rtk.ImageFileWriter[ImageType]
> writer = WriterType.New()
> writer.SetFileName("rec.tif")
> writer.SetInput(feldkamp.GetOutput())
> writer.Update()
>
> Regards, Ruslan
> _______________________________________________
> 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/20231206/dc0faa08/attachment-0001.htm>