ufo-kit / tofu

Helper scripts for tomographic reconstruction using the ufo-core framework
GNU Lesser General Public License v3.0
18 stars 9 forks source link

Optimizing the TIE phase-retrieval #69

Closed sgasilov closed 5 years ago

sgasilov commented 6 years ago

Hello guys,

please advise on the best way to make TIE phase retrieval most efficient. IRL people apply phase retrieval to 80% of CT data thus it will be very useful but currently UFO user must write an intermediate data set on the HDD in order to make the TIE-PR correctly. The reason for extra write/read operation is that TIE has to be applied to flat-corrected frames before the "-log" is taken (that is the "absorptivity" parameter). I do not see any way to do it with existing tools so currently I first do flatcorrect without the absorptivity, and then in the second pipeline add a separate, home-made "absorptivity" filter after the phase-retrieval.

I suggest to add section "tofu PR" in the main tofu tree. As for the actual implementation, I have looked in the reco.py and preprocess.py and was wondering what is most correct way to do this: 1) make use of the existing pipelines: f.i. create a separate create_flat_correct_pipeline which bypasses the absorptivity or sets absorptivity.props.expression = '-log(v)' to one. Then establish pipeline consisting of three steps "flat correction" --> create_phase_retrieval_pipeline --> with the separate "absorptivity" filter. 2) initialize a separate pipeline purely dedicated to TIR PR. Once again IRL 80% of synchrotron data sets require phase retrieval (ANKA TOTO is a rare exception from this case because almost everything is done with the white beam) so it is really worth to investigate this question.

Cheers, Sergei.

tfarago commented 6 years ago

the phase retrieval should give you the phase itself which doesn't need any logarithm anymore or am I wrong? That's why I thought it would work just fine without specifying the --absorptivity option.

sgasilov commented 6 years ago

Hello Tomas,

could you refer to an equation? I apply the UFO TIE phase-retrieval in the sense of Paganin's filter ( equation (10) from the reference: J Microsc. 2002 Apr;206(Pt 1):33-40 Simultaneous phase and amplitude extraction from a single defocused image of a homogeneous object.) As you my see when logarithm is applied the result obtains a physical meaning. Paganin simply calls it "projected thickness" which in fact must mean the projected sum of "weighting factors" in front of different delta/beta ratios in your sample.

Certainly, the difference can be very large since logarithm is a nonlinear operation.

Kind regards, Sergei.

tfarago commented 6 years ago

Thanks for info. In this case I would link this to issue ufo-filters:145. So I will first check what the filter actually does with @rshkarin and then we can get back here.

sgasilov commented 6 years ago

Hello,

we've just reviewed the source of the preprocess.py and apparently there is already everything that I need. I will simply add the absorptivity in the end in the retrieve-phase pipeline: log_phase_retrieve = get_task('calculate', processing_node=processing_node) log_phase_retrieve.props.expression = '-log(v)' graph.connect_nodes(crop_phase_retrieve, log_phase_retrieve) and call tofu preprocessing with both flatcorrect and retrieve-phase parameters so that "if" on the line 294 returns True.

I know that hardcoding absorptivity is not the most flexible solution, but in practice in 99% of cases simplest TIE is used anyways. If someone will need other functionality of the phase-retrieval filter, it can be always provided using ufo-launch.

I would have asked you to close the question, but I'm going to use this opportunity to remind you one more time:

PLEASE confirm that spatial frequencies as defined in the phase-retrieval filter are identical to the spatial frequencies of the target image as your fft module returns them !!!

Have a good weekend, Serg.

tfarago commented 6 years ago

@rshkarin Can you check that please?

sgasilov commented 6 years ago

Hello again,

and many thanks! I've just realized that there is everything that I need because someone has created all that preprocess section in the past several weeks )) I have tested it - it work as you intended it to work, but I see several problems:

minor problem: --projection-filter ‘none’ must be explicitly specified when calling tofu preprocess otherwise “if” on line 288 is True by default and it makes a mess.

major: whatever I do with default parameters of the argument absorptivity, the ffc_task totally disregards it and treats it as true - I shall open a separate issue for this (also applies to the previous tofu version which i downloaded in November and to flatcorrect thing).

So returning to FFC+TIE I have done two corrections and I suggest you to consider making them default:

(1) I've changed “return” line 161 to take negative log of TIE filtered images: if args.retrieval_method == 'tie': log_phase_retrieve = get_task('calculate', processing_node=processing_node) log_phase_retrieve.props.expression = \ '( isnan(v) || isinf(v) || (v<=0) ) ? 0.0 : -log(v)' graph.connect_nodes(crop_phase_retrieve, log_phase_retrieve) crop_phase_retrieve = log_phase_retrieve return (pad_phase_retrieve, crop_phase_retrieve)

(2) I’ve simply commented lines 288-293 because I do not see use for them now. Selection of filters is currently limited to tomographic filters, which I would never apply to projections anyways.

Many thanks again, And happy holidays to you, Serg.

tfarago commented 6 years ago

minor problem: --projection-filter ‘none’ must be explicitly specified when calling tofu preprocess otherwise “if” on line 288 is True by default and it makes a mess.

I left the filtering in place because 99.9 % of use cases for us is the reconstruction, so it is nice to have it as default. And it was actually a user request to enable filtering by preprocessing, so that one can run the reconstruction multiple times with different parameters afterward without the need to filter the data all the time. What does @matze think?

(1) I've changed “return” line 161 to take negative log of TIE filtered images:

@rshkarin we could really use your help here. I checked the formula and it is not the same as the equation (10) in the paper mentioned above. Can you shed some light on this? I think the logarithm should be inherent to the PR kernel.

sgasilov commented 6 years ago

I'm confused, should I then simly consider filtering functionality of preprocess as alternative to ufo-launch? I understand that you may apply any filter from the ufo collection, not only one of the CT filters. This is OK if you only need one filter. But if user need more then one he inevitably has to use ufo-launch to make a pipeline, is it right? Really, now when absorptivity and nan-and-inf-filters are included in ufo-filters I can create a pipeline which applies 4 different filters and does TIE phase-retrieval - and everything also in a single read/write cycle.

Anyways, returning to the original question, it still does not make sense that the default value of preprocess_filter is ramp_from_real -> you would never really apply this filter to projection before CT reconstruction... And even if you would like to apply tone of CT kernels to a certain type of input images, then you need backprojection anyways.

As for the TIE, if you rewrite the equation you may see that (1) log is missing, but the filter is the same as in Paganin paper, but, again but, its true that naming is very confusing like delta/beta ratio is called regularization rate, it is defined as power of 10, etc.
(2) the equation is never used correctly because literally almost no one applies Paganing filter to a phase contrast image normalized by the contact image. People simply take flatcorrected phase contrast image.

Cheers, Sergei.

sgasilov commented 6 years ago

Just another important remark on filtering and flatcorrection in tofu preprocess. Providing that someone decides to apply one of the legit ufo filters to images using preprocess, he must be aware of the fact that the filter will be always applied after the flat field correction is done. While strictly speaking in many cases filter must be applied to raw frames before the flat-field correction in order to eliminate possible artifacts as early as possible.

tfarago commented 6 years ago

You mean the tomographic filter, like ramp applied separately to darks, flats and radios?

sgasilov commented 6 years ago

Hi Tomas, no, certainly I do not mean this. I do not understand how preprocess module is connected to tomo module - neither logically nor in the python code of preprocess.py python code - except that for all modules share the same input parameters defined in config.py. I only thought that the preprocess module is designed to make flatcorrection+TIE PR in one read/write cycles - thus it is already fast and must be used for batch processing after all parameters are optimized. Certainly for optimization of parameters one needs to process one data set step-by-step and should probably save separately results of filtering, flat-correction and the phase-retrieval.

We are probably already beyond the original subject but just to sum up: For me the working sequence and corresponding tools are: ufo-launch for preliminary filtering of raw projections, preprocess for flatcorrects+TIE --> tomo for reconstruction. Is it right or not?

Also, maybe you let me know when I can call you, seems like typing does not bring us anywhere in this case ))

Cheers, Serg.