Hi, I created another topic since the content is different from the previous one and people could land here with a different search in the archives.
While trying to reconstruct a set of images from 0 and 180° without having to add the fan angle, Simon recommended to use the conjugate gradient algorithm in RTK. So, I first tried to replicate in Python the conjugate gradient application written in C++ ( https://github.com/SimonRit/RTK/tree/master/applications/rtkconjugategradient). But whatever types of input I set, my python kernel dies when updating the filter without any error or debug message. I spent some time thinking about it but I can't find out what is wrong. I have attached my code at the end of this message if you have suggestions about what would be wrong. Thanks for your help, Clément ------------------------------------------------------------------------------------------------------------ # Defines the image type ImageType = itk.Image[itk.F,3] # Defines the RTK geometry object geometry = rtk.ThreeDCircularProjectionGeometry.New() numberOfProjections = 180 firstAngle = 90 sid = 600 # source to isocenter distance sdd = 1200 # source to detector distance angularArc = 180 for x in range(0,numberOfProjections): angle = firstAngle + x * angularArc / numberOfProjections geometry.AddProjection(sid,sdd,angle) # Create a stack of empty projection images ConstantImageSourceType = rtk.ConstantImageSource[ImageType] constantImageSource = ConstantImageSourceType.New() origin = [ -199, -199, 0. ] sizeOutput = [ 200, 200, 180 ] spacing = [ 2.0, 2.0, 2.0 ] constantImageSource.SetOrigin( origin ) constantImageSource.SetSpacing( spacing ) constantImageSource.SetSize( sizeOutput ) constantImageSource.SetConstant(0.) # Forward projection of input volume proj = rtk.JosephForwardProjectionImageFilter[ImageType, ImageType].New() proj.SetInput(constantImageSource) proj.SetInput(1, input) # input is ImageType from a reader proj.SetGeometry(geometry) # Create weights weightsSource = ConstantImageSourceType.New() weightsSource.SetInformationFromImage(proj.GetOutput()) weightsSource.SetConstant(1.0) # Create output like volume outputImageSource = ConstantImageSourceType.New() sizeOutput = [ 400, 400, 400] origin = [ -200, -200, -200 ] spacing = [ 1.0, 1.0, 1.0 ] outputImageSource.SetOrigin(origin) outputImageSource.SetSpacing(spacing) outputImageSource.SetSize(sizeOutput) outputImageSource.SetConstant(0.) # Create reconstructed image RCGType = rtk.ConjugateGradientConeBeamReconstructionFilter[ImageType] conjugateGradient = RCGType.New() conjugateGradient.SetInput(outputImageSource.GetOutput()) conjugateGradient.SetInput(1, proj.GetOutput()) conjugateGradient.SetInput(2, weightsSource.GetOutput()) conjugateGradient.SetGeometry(geometry) conjugateGradient.SetGamma(2.0) conjugateGradient.SetNumberOfIterations(30) conjugateGradient.SetDisableDisplacedDetectorFilter(True) conjugateGradient.Update()
_______________________________________________ Rtk-users mailing list Rtk-users@public.kitware.com https://public.kitware.com/mailman/listinfo/rtk-users
