Thanks for getting back on this. I really need to update the pypi packages but I'm fighting a bit with the CI, it's no longer working (see #395 <https://github.com/SimonRit/RTK/pull/395>). I'll keep everyone posted on the mailing list of course.
On Fri, Dec 18, 2020 at 12:17 PM Clément Mirabel <clement.mira...@gmail.com> wrote: > Hi, > > Thank you for your answer. I still had the issue with your updated code > but I noticed I was still using itk 5.0.1 because of an issue > <https://tinyurl.com/y3dm3h86> with itk 5.1.2 due to numpy 1.9.4 with > Windows 10 recent update. I downgraded numpy to 1.9.3 so itk could install > 5.1.2 and then installed the artifact you pointed for RTK. > Everything works fine now, time to play with the parameters! > > Clément > > > > Le ven. 18 déc. 2020 à 08:00, Simon Rit <simon....@creatis.insa-lyon.fr> > a écrit : > >> Hi, >> I slightly modified your script and it worked for me. What version are >> you using? I used the latest ITK v5.1.2 with, for RTK, an artifact >> published here <https://github.com/SimonRit/RTK/actions/runs/344449575> >> (visible only if you're logged in github). >> Simon >> >> PS: I don't think you need the displaced detector, it is centered on the >> central ray in your simulation. >> >> import itk >> from itk import RTK as rtk >> >> # 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 = [ 100, 100, 180 ] >> spacing = [ 4.0, 4.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, itk.imread('test.mha')) # input is ImageType from a >> reader >> proj.SetGeometry(geometry) >> proj.Update() >> >> # Create weights >> weightsSource = ConstantImageSourceType.New() >> weightsSource.SetInformationFromImage(proj.GetOutput()) >> weightsSource.SetConstant(1.0) >> >> # Create output like volume >> outputImageSource = ConstantImageSourceType.New() >> sizeOutput = [ 100, 100, 100] >> origin = [ -200, -200, -200 ] >> spacing = [ 4.0, 4.0, 4.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(2) >> conjugateGradient.SetDisableDisplacedDetectorFilter(False) >> conjugateGradient.Update() >> >> itk.imwrite(conjugateGradient.GetOutput(), 'output.mha') >> >> On Wed, Dec 16, 2020 at 10:36 AM Clément Mirabel < >> clement.mira...@gmail.com> wrote: >> >>> 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 >>> >>
_______________________________________________ Rtk-users mailing list Rtk-users@public.kitware.com https://public.kitware.com/mailman/listinfo/rtk-users
