Hi Cyril. I performed a SIRT reconstruction with 100 iterations as you suggested. After a whopping 80 hours of calculation time on 64 cores, the result was a bit disappointing (in the dropbox: https://www.dropbox.com/sh/7a4bkyjg43zjmy7/AAAp4tJrk9HedMEEuIKsGYDSa?dl=0). If anything, it looks worse than 3 iterations. Do you think we need even more iterations? Another question, you said the conjugate gradient method minimizes the same cost function as the SIRT method, but in the doxygen documentation I find that the conjugate gradient is similar to SART (not SIRT). Could you elaborate a bit on this please? Thanks again!
On 18-07-17 14:28, Cyril Mory wrote: > The weights are a way to adapt the cost function if you know that some > pixels of your projection data are more reliable than others (you then > give them a larger weight than the others). This is called Weighted > Least Squares (WLS), you can easily find a lot of theory on it. IYou > should only do this if you have a reliable way to estimate the > variance of the noise on your projection data, and use the inverse > variances as weights. There is one weight per pixel, so the weights > image should be exactly the same size as the full projection dataset. > Note that this input is mandatory for the > rtk::ConjugateGradientConeBeamReconstructionFilter, so if you do not > want to use it, you still have to generate an image full of ones and > pass it as the weights. Clearly that behavior is not optimal, but that > is the way it works at the moment. > > The support mask enforces a support constraint. It should be a binary > volume (ones in the ROI you are interested in reconstructing, zeros > outside). If it is provided, over the course of the iterations, the > iterates are multiplied by the mask, leaving the voxels inside the ROI > unchanged, and setting the ones outside the ROI to zero. > Mathematically, it is a special case of preconditioner (a diagonal, > binary preconditioner). Not passing a mask yields the same result (but > faster) as passing a mask full of ones. > > Hope this helps, > Cyril > > PS: If at some point you decide to use the "Weights" input of the > conjugate gradient filter, in order to perform WLS reconstructions, do > write another email on this list: I have a lot of advice to share on > this topic. > > On 18/07/2017 14:02, Lotte Schyns wrote: >> Thanks! I will try both the SIRT with 100 iterations and the conjugate >> gradient method. What exactly are the weights (input 2) and input >> support mask in rtkConjugateGradientConeBeamReconstructionFilter? >> >> On 18-07-17 11:31, Cyril Mory wrote: >>> Hi Lotte, >>> >>> I had a very quick look at what you are doing. It looks fine to me, >>> you just need to perform way more iterations in the SIRT case. 100 >>> would be a good start. SIRT is more stable than SART when there are >>> inconsistencies in the projection data, but converges slowly. >>> >>> An alternative to SIRT, which minimizes the same cost function with a >>> faster algorithm, is the conjugate gradient algorithm. You should >>> obtain nice results with something like 30-40 iterations (look for >>> rtk::ConjugateGradientConeBeamReconstructionFilter if you want to give >>> it a try). >>> >>> Best regards, >>> Cyril >>> >>> >>> On 18/07/2017 11:17, Lotte Schyns wrote: >>>> Hello, >>>> >>>> We are having problems with the SIRT reconstructions. They seem very >>>> strange and blurry. However, other reconstructions (SART, FDK, >>>> iterativeFDK) look perfect. I uploaded an example of a SART and a SIRT >>>> reconstruction (same parameters) to >>>> https://www.dropbox.com/sh/7a4bkyjg43zjmy7/AAAp4tJrk9HedMEEuIKsGYDSa?dl=0. >>>> >>>> >>>> I also uploaded a minimalistic version of the code that I used (the >>>> paths in CMakeLists.txt probably need to be adapted to your >>>> system). You >>>> can alternate between SART and SIRT by changing line 143. The raw data >>>> is also available on dropbox. I didn't use >>>> rtkXRadRawToAttenuationImageFilter, because our projections are >>>> already >>>> corrected for the dark field and flood field, so >>>> (signal-dark)/(flood-dark). I just take the natural logarithm and >>>> multiply by -1. Do you know what could be the problem with the SIRT >>>> reconstructions? Are we using wrong parameters? Thanks for your time. >>>> >>>> Lotte >>>> _______________________________________________ >>>> Rtk-users mailing list >>>> [email protected] >>>> http://public.kitware.com/mailman/listinfo/rtk-users >> _______________________________________________ >> Rtk-users mailing list >> [email protected] >> http://public.kitware.com/mailman/listinfo/rtk-users > _______________________________________________ Rtk-users mailing list [email protected] http://public.kitware.com/mailman/listinfo/rtk-users
