-
Notifications
You must be signed in to change notification settings - Fork 1
Normalization step in the YMIR image workflow. #42
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Conversation
| """ | ||
| sample_images = _select_images_by_key(image_stacks, ImageKey.SAMPLE) | ||
| dark_current = _select_image_by_key(bg_images, ImageKey.DARK_CURRENT) | ||
| return CleansedSampleImages(sample_images - dark_current) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This part is very slow... something must be wrong here.
src/ess/imaging/normalize.py
Outdated
| def calculate_d0(background: BackgroundImage) -> D0: | ||
| """Calculate the D0 value from background image stack. | ||
| :math:`D0 = mean(background counts of all pixels)` | ||
| """ | ||
| return D0(sc.mean(background)) | ||
|
|
||
|
|
||
| def calculate_d(samples: CleansedSampleImages) -> D: | ||
| """Calculate the D value from the sample image stack. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This part is so slow too.
e228af4 to
6f4bc9f
Compare
84922f8 to
dab573f
Compare
ca2c9ac to
b20ae4f
Compare
35c3f41 to
9fb1e48
Compare
a11c081 to
ffabecb
Compare
| OpenBeam = mean(open_beam, 'time') | ||
| """ | ||
| return OpenBeamImage(sc.mean(open_beam, dim=TIME_COORD_NAME)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Question: Is the time coordinate here time during the measurement, or is it time-of-flight?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
it's timestamp of each image, when it was taken.
It's an optical 2-d detector so it doesn't have time-of-flight.
src/ess/imaging/normalize.py
Outdated
| Threshold for the background pixel values. | ||
| Any pixel values less than ``background_threshold`` | ||
| are replaced with ``background_threshold``. | ||
| Default is 1.0[counts]. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Is this done to avoid nans when doing the division in the normalization?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes. There is another python project called pymureh and I mostly followed their method,
but I thought it was a bit weird to hard-code the number so I exposed it as an argument.
src/ess/imaging/normalize.py
Outdated
| .. math:: | ||
| Background = mean(OpenBeam, 'time') - mean(DarkCurrent, 'time') | ||
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Can you also add (in words or in the math) the fact that you are applying a threshold after the operation?
src/ess/imaging/normalize.py
Outdated
| .. math:: | ||
| CleansedSample_{i} = Sample_{i} - mean(DarkCurrent, dim='time') |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Same here: also add (in words or in the math) the fact that you are applying a threshold after the operation?
src/ess/imaging/normalize.py
Outdated
| and returned negative values. | ||
| """ | ||
| return AverageSamplePixelCounts( | ||
| _mean_all_dims(sample_images.data) - dark_current.data.mean() |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I am wondering if sums should be used instead of means for both the AverageBackgroundPixelCounts and the AverageSamplePixelCounts?
Usually when we normalize by e.g. a monitor, we sum the counts?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
But we don't take same number of dark-current/open-beam images as the sample images.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Probably best to discuss in person / over a call, as I am wondering if the reduction over the time dimension should also be sums instead of means?
If I understand correctly, it's to make sure everything has been scaled to the same integration time, before performing the division, right?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
If I understand correctly, it's to make sure everything has been scaled to the same integration time, before performing the division, right?
Yes, that was the idea.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Or maybe it should be more like:
(sample_image.sum('time') / sample_image_integration_time - dark_frame.sum('time') / dark_frame_integration_time) / (openbeam_image.sum('time') / openbeam_image_integration_time - dark_frame.sum('time') / dark_frame_integration_time)
?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Number of images is same as the integration_time since the shutter(?) speed is fixed.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ah so you're saying that doing a mean is essentially doing image.sum(time) / integration_time since you divide by the number of images (there is just a factor of time_shutter_is_open which cancels out everywhere)?
I guess if the shutter speed is adjustable this should be verified?
Are the dark frames and open beam always taken just before the sample is, or is a dark frame sometimes taken weeks in advance and used for several experiments in a row?
If it's the former, then I guess we are fine with means.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ah so you're saying that doing a mean is essentially doing image.sum(time) / integration_time since you divide by the number of images (there is just a factor of time_shutter_is_open which cancels out everywhere)?
Yes...! That's what I meant.
But maybe I should check with Robin once more and mention this in the docstring.
I guess if the shutter speed is adjustable this should be verified?
Yes. Totally.
Are the dark frames and open beam always taken just before the sample is, or is a dark frame sometimes taken weeks in advanced and used for several experiments in a row?
They will be taken again every time just before or after the sample images.
Background images(dark-current and open-beam) takes much less time than samples and
Odin team wanted every types of images in one file.
If it's the former, then I guess we are fine with means.
Yes. In either case, I'll update the docstring.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I double checked with Søren and opened an issue as a reminder.
So the exposure time is constant within a single file, but may vary in different files.
And it gives warning in the average/normalization functions about exposure time.
I chose warning over docstring because I think we should use exposure-time as you said,
and it should be avilable in the file.
src/ess/imaging/normalize.py
Outdated
| def calculate_white_beam_background( | ||
| open_beam: OpenBeamImage, | ||
| dark_current: DarkCurrentImage, | ||
| background_threshold: BackgroundPixelThreshold = DEFAULT_BACKGROUND_THRESHOLD, |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'm thinking that BackgroundPixelThreshold should not have a default value here, but instead a default param should be set on the workflow/pipeline ?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
It does have a default param in the workflow/pipeline.
The default value in the function signature is not used by the sciline.Pipeline.
I thought it'll make it easier to use the provider itself, i.e. for testing...
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I would say if you want to use it for testing, then it's good that the parameter is required, so that you understand exactly what is being done?
An alternative would be to split it up into 2 providers: one that does the subtraction, and one that applies the threshold afterwards, with a very obvious name like ReplaceZerosAndNegativesWithOne?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I would say if you want to use it for testing, then it's good that the parameter is required, so that you understand exactly what is being done?
You can still test it with the default argument value though.
An alternative would be to split it up into 2 providers: one that does the subtraction, and one that applies the threshold afterwards, with a very obvious name like ReplaceZerosAndNegativesWithOne?
I didn't want to hard-code the numbers so the name should be different but
I can split this up into two steps if that makes it more interpretable.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Second thoughts... I think rounding up should be done within the normalization step since it's actually where it's needed.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@nvaytet I ended up inserting threshold applying steps instead of merging them into the normalization step
since they are used for calculating scale factor as well...
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
But I still don't see why we can't use the default argument value though
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'm saying that because I think in esssans, we don't have default arguments in providers, only default parameters on the workflow (unless I missed something somewhere?).
So I thought we should try to be consistent across techniques.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@neil Okay I fixed it accordingly...!
src/ess/imaging/normalize.py
Outdated
| def cleanse_sample_images( | ||
| sample_images: SampleImageStacks, | ||
| dark_current: DarkCurrentImage, | ||
| sample_threshold: SamplePixelThreshold = DEFAULT_SAMPLE_THRESHOLD, |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Same as above: set the parameter on the workflow?
| SamplePixelThreshold: DEFAULT_SAMPLE_THRESHOLD, | ||
| FileLock: DEFAULT_FILE_LOCK, | ||
| }, | ||
| ) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Tests will come from different PR
Unless this change is urgent (e.g. for the STAP demo), I would rather see the unit tests as part of the same PR?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yeah I wanted to release it before Christian leaves for vacation.
But I think I can just test them by myself.... I'll add tests here.
e2f173d to
1c1edb7
Compare
Co-authored-by: Neil Vaytet <39047984+nvaytet@users.noreply.github.com>
1c1edb7 to
3f5b477
Compare
| assert isinstance(normalized, sc.DataArray) | ||
| assert normalized.sizes['time'] == 2 | ||
| assert normalized.unit == "dimensionless" | ||
| assert_allclose(normalized, expected_normalized_sample_images) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Can you also add some tests on the functions that apply the thresholds on the sample and dark images?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@nvaytet Done...!
Co-authored-by: Neil Vaytet <39047984+nvaytet@users.noreply.github.com>
Normalization method is implemented.
Tests will come from different PRTest included in this PR.