Skip to content
GitLab
Commit 882c2091 authored by myron's avatar myron Committed by Pierre Paleo
Browse files

refined the selection of the best ROI

parent b2794cf5
Expand all Hide whitespace changes
Inline Side-by-side
nabu/preproc/alignment.py
View file @ 882c2091
This diff is collapsed.
nabu/preproc/tests/test_alignment.py
View file @ 882c2091
......@@ -48,8 +48,18 @@ def bootstrap_fcs(request):
cls.abs_tol_dist = 1e-2
cls.abs_tol_tilt = 2.5e-4
cls.data, cls.img_pos, cls.pixel_size, (calib_data_std, calib_data_angle) = get_focus_data("test_alignment_focus.h5")
cls.angle_best_ind, cls.angle_best_pos, cls.angle_tilt_v, cls.angle_tilt_h = calib_data_angle
(
cls.data,
cls.img_pos,
cls.pixel_size,
(calib_data_std, calib_data_angle),
) = get_focus_data("test_alignment_focus.h5")
(
cls.angle_best_ind,
cls.angle_best_pos,
cls.angle_tilt_v,
cls.angle_tilt_h,
) = calib_data_angle
cls.std_best_ind, cls.std_best_pos = calib_data_std
......@@ -64,10 +74,14 @@ def get_cor_data_h5(*dataset_path):
with h5py.File(dataset_downloaded_path, "r") as hf:
data = hf["/entry/instrument/detector/data"][()]
cor_global_pix = hf['/calibration/alignment/global/x_rotation_axis_pixel_position'][()]
cor_global_pix = hf[
"/calibration/alignment/global/x_rotation_axis_pixel_position"
][()]
cor_highlow_pix = hf['/calibration/alignment/highlow/x_rotation_axis_pixel_position'][()]
tilt_deg = hf['/calibration/alignment/highlow/z_camera_tilt'][()]
cor_highlow_pix = hf[
"/calibration/alignment/highlow/x_rotation_axis_pixel_position"
][()]
tilt_deg = hf["/calibration/alignment/highlow/z_camera_tilt"][()]
return data, (cor_global_pix, cor_highlow_pix, tilt_deg)
......@@ -109,10 +123,10 @@ def get_alignxc_data(*dataset_path):
img_pos = hf["/entry/instrument/detector/distance"][()]
unit_length_shifts_vh = [
hf['/calibration/alignxc/y_pixel_shift_unit'][()],
hf['/calibration/alignxc/x_pixel_shift_unit'][()]
hf["/calibration/alignxc/y_pixel_shift_unit"][()],
hf["/calibration/alignxc/x_pixel_shift_unit"][()],
]
all_shifts_vh = hf['/calibration/alignxc/yx_pixel_offsets'][()]
all_shifts_vh = hf["/calibration/alignxc/yx_pixel_offsets"][()]
return data, img_pos, (unit_length_shifts_vh, all_shifts_vh)
......@@ -129,18 +143,20 @@ def get_focus_data(*dataset_path):
data = hf["/entry/instrument/detector/data"][()]
img_pos = hf["/entry/instrument/detector/distance"][()]
pixel_size = np.mean([
hf['/entry/instrument/detector/x_pixel_size'][()],
hf['/entry/instrument/detector/y_pixel_size'][()]
])
pixel_size = np.mean(
[
hf["/entry/instrument/detector/x_pixel_size"][()],
hf["/entry/instrument/detector/y_pixel_size"][()],
]
)
angle_best_ind = hf['/calibration/focus/angle/best_img'][()]
angle_best_pos = hf['/calibration/focus/angle/best_pos'][()]
angle_tilt_v = hf['/calibration/focus/angle/tilt_v_rad'][()]
angle_tilt_h = hf['/calibration/focus/angle/tilt_h_rad'][()]
angle_best_ind = hf["/calibration/focus/angle/best_img"][()]
angle_best_pos = hf["/calibration/focus/angle/best_pos"][()]
angle_tilt_v = hf["/calibration/focus/angle/tilt_v_rad"][()]
angle_tilt_h = hf["/calibration/focus/angle/tilt_h_rad"][()]
std_best_ind = hf['/calibration/focus/std/best_img'][()]
std_best_pos = hf['/calibration/focus/std/best_pos'][()]
std_best_ind = hf["/calibration/focus/std/best_img"][()]
std_best_pos = hf["/calibration/focus/std/best_pos"][()]
calib_data_angle = (angle_best_ind, angle_best_pos, angle_tilt_v, angle_tilt_h)
calib_data_std = (std_best_ind, std_best_pos)
......@@ -221,9 +237,13 @@ class TestAlignmentBase(object):
yy, xx = np.meshgrid(fy, fx, indexing="ij")
peak_pos_yx = np.random.rand(2) * 1.6 - 0.8
f_vals = np.exp(-((yy - peak_pos_yx[0]) ** 2 + (xx - peak_pos_yx[1]) ** 2) / 100)
f_vals = np.exp(
-((yy - peak_pos_yx[0]) ** 2 + (xx - peak_pos_yx[1]) ** 2) / 100
)
fitted_peak_pos_yx = alignment.AlignmentBase.refine_max_position_2d(f_vals, fy, fx)
fitted_peak_pos_yx = alignment.AlignmentBase.refine_max_position_2d(
f_vals, fy, fx
)
message = (
"Computed peak position: (%f, %f) " % (*fitted_peak_pos_yx,)
......@@ -231,7 +251,9 @@ class TestAlignmentBase(object):
+ " Difference: (%f, %f)," % (*(fitted_peak_pos_yx - peak_pos_yx),)
+ " tolerance: %f" % self.abs_tol
)
assert np.all(np.isclose(peak_pos_yx, fitted_peak_pos_yx, atol=self.abs_tol)), message
assert np.all(
np.isclose(peak_pos_yx, fitted_peak_pos_yx, atol=self.abs_tol)
), message
def test_peak_fitting_2d_error_checking(self):
# Fit a 3 x 3 grid
......@@ -240,7 +262,9 @@ class TestAlignmentBase(object):
yy, xx = np.meshgrid(fy, fx, indexing="ij")
peak_pos_yx = np.random.rand(2) + 1.5
f_vals = np.exp(-((yy - peak_pos_yx[0]) ** 2 + (xx - peak_pos_yx[1]) ** 2) / 100)
f_vals = np.exp(
-((yy - peak_pos_yx[0]) ** 2 + (xx - peak_pos_yx[1]) ** 2) / 100
)
with pytest.raises(ValueError) as ex:
alignment.AlignmentBase.refine_max_position_2d(f_vals, fy, fx)
......@@ -259,10 +283,15 @@ class TestAlignmentBase(object):
cc_coords = np.arange(0, 8)
(found_peaks_val, found_peaks_pos) = alignment.AlignmentBase.extract_peak_regions_1d(img, axis=-1, cc_coords=cc_coords)
message = "The found peak positions do not correspond to the expected peak positions:\n Expected: %s\n Found: %s" % (
peaks_pos,
found_peaks_pos[1, :],
(
found_peaks_val,
found_peaks_pos,
) = alignment.AlignmentBase.extract_peak_regions_1d(
img, axis=-1, cc_coords=cc_coords
)
message = (
"The found peak positions do not correspond to the expected peak positions:\n Expected: %s\n Found: %s"
% (peaks_pos, found_peaks_pos[1, :],)
)
assert np.all(peaks_val == found_peaks_val[1, :]), message
......@@ -276,7 +305,10 @@ class TestCor(object):
CoR_calc = alignment.CenterOfRotation()
cor_position = CoR_calc.find_shift(radio1, radio2)
message = "Computed CoR %f " % cor_position + " and real CoR %f do not coincide" % self.cor_gl_pix
message = (
"Computed CoR %f " % cor_position
+ " and real CoR %f do not coincide" % self.cor_gl_pix
)
assert np.isclose(self.cor_gl_pix, cor_position, atol=self.abs_tol), message
def test_noisy_cor_posx(self):
......@@ -289,7 +321,10 @@ class TestCor(object):
CoR_calc = alignment.CenterOfRotation()
cor_position = CoR_calc.find_shift(radio1, radio2, median_filt_shape=(3, 3))
message = "Computed CoR %f " % cor_position + " and real CoR %f do not coincide" % self.cor_gl_pix
message = (
"Computed CoR %f " % cor_position
+ " and real CoR %f do not coincide" % self.cor_gl_pix
)
assert np.isclose(self.cor_gl_pix, cor_position, atol=self.abs_tol), message
@pytest.mark.skipif(not (__has_scipy__), reason="need scipy for this test")
......@@ -314,7 +349,10 @@ class TestCor(object):
# cor_position = CoR_calc.find_shift(radio1, radio2)
cor_position = CoR_calc.find_shift(radio1, radio2, low_pass=(6.0, 0.3))
message = "Computed CoR %f " % cor_position + " and real CoR %f do not coincide" % self.cor_gl_pix
message = (
"Computed CoR %f " % cor_position
+ " and real CoR %f do not coincide" % self.cor_gl_pix
)
assert np.isclose(self.cor_gl_pix, cor_position, atol=self.abs_tol), message
def test_half_tomo_cor_Err(self):
......@@ -330,7 +368,8 @@ class TestCor(object):
message = (
"Computed CoR %f " % cor_position
+ " and real CoR %f should not coincide when using the standard algorithm with half tomo data" % cor_pos
+ " and real CoR %f should not coincide when using the standard algorithm with half tomo data"
% cor_pos
)
assert not np.isclose(cor_pos, cor_position, atol=self.abs_tol), message
......@@ -342,11 +381,14 @@ class TestCor(object):
cor_pos = self.ht_cor
CoR_calc = alignment.CenterOfRotation()
cor_position = CoR_calc.find_shift(radio1, radio2, low_pass=1, high_pass=20, half_tomo_cor_guess=cor_pos - 10.0)
cor_position = CoR_calc.find_shift(
radio1, radio2, low_pass=1, high_pass=20, half_tomo_cor_guess=cor_pos - 10.0
)
message = (
"Computed CoR %f " % cor_position
+ " and real CoR %f should coincide when using the halftomo algorithm with hald tomo data" % cor_pos
+ " and real CoR %f should coincide when using the halftomo algorithm with hald tomo data"
% cor_pos
)
assert np.isclose(cor_pos, cor_position, atol=self.abs_tol), message
......@@ -360,19 +402,22 @@ class TestCor(object):
radio2 = np.fliplr(radio2)
cor_pos = 983.038
cor_pos = 983.107
CoR_calc = alignment.CenterOfRotation()
cor_position = CoR_calc.find_shift(radio1, radio2, low_pass=1, high_pass=20, global_search=True)
cor_position = CoR_calc.find_shift(
radio1, radio2, low_pass=1, high_pass=20, global_search=True
)
print("Found cor_position", cor_position)
message = (
"Computed CoR %f " % cor_position
+ " and real CoR %f should coincide when using the halftomo algorithm with hald tomo data" % cor_pos
+ " and real CoR %f should coincide when using the halftomo algorithm with hald tomo data"
% cor_pos
)
assert np.isclose(cor_pos, cor_position, atol=0.1), message
assert np.isclose(cor_pos, cor_position, atol=self.abs_tol), message
def test_cor_posx_linear(self):
radio1 = self.data[0, :, :]
......@@ -381,7 +426,10 @@ class TestCor(object):
CoR_calc = alignment.CenterOfRotation()
cor_position = CoR_calc.find_shift(radio1, radio2, padding_mode="edge")
message = "Computed CoR %f " % cor_position + " and real CoR %f do not coincide" % self.cor_gl_pix
message = (
"Computed CoR %f " % cor_position
+ " and real CoR %f do not coincide" % self.cor_gl_pix
)
assert np.isclose(self.cor_gl_pix, cor_position, atol=self.abs_tol), message
def test_error_checking_001(self):
......@@ -393,8 +441,13 @@ class TestCor(object):
with pytest.raises(ValueError) as ex:
CoR_calc.find_shift(radio1, radio2)
message = "Error should have been raised about img #1 shape, other error raised instead:\n%s" % str(ex.value)
assert "Images need to be 2-dimensional. Shape of image #1" in str(ex.value), message
message = (
"Error should have been raised about img #1 shape, other error raised instead:\n%s"
% str(ex.value)
)
assert "Images need to be 2-dimensional. Shape of image #1" in str(
ex.value
), message
def test_error_checking_002(self):
CoR_calc = alignment.CenterOfRotation()
......@@ -405,8 +458,13 @@ class TestCor(object):
with pytest.raises(ValueError) as ex:
CoR_calc.find_shift(radio1, radio2)
message = "Error should have been raised about img #2 shape, other error raised instead:\n%s" % str(ex.value)
assert "Images need to be 2-dimensional. Shape of image #2" in str(ex.value), message
message = (
"Error should have been raised about img #2 shape, other error raised instead:\n%s"
% str(ex.value)
)
assert "Images need to be 2-dimensional. Shape of image #2" in str(
ex.value
), message
def test_error_checking_003(self):
CoR_calc = alignment.CenterOfRotation()
......@@ -417,8 +475,9 @@ class TestCor(object):
with pytest.raises(ValueError) as ex:
CoR_calc.find_shift(radio1, radio2)
message = "Error should have been raised about different image shapes, " + "other error raised instead:\n%s" % str(
ex.value
message = (
"Error should have been raised about different image shapes, "
+ "other error raised instead:\n%s" % str(ex.value)
)
assert "Images need to be of the same shape" in str(ex.value), message
......@@ -429,10 +488,16 @@ class TestCor(object):
tilt_calc = alignment.CameraTilt()
cor_position, camera_tilt = tilt_calc.compute_angle(radio1, radio2)
message = "Computed tilt %f " % camera_tilt + " and real tilt %f do not coincide" % self.tilt_deg
message = (
"Computed tilt %f " % camera_tilt
+ " and real tilt %f do not coincide" % self.tilt_deg
)
assert np.isclose(self.tilt_deg, camera_tilt, atol=self.abs_tol), message
message = "Computed CoR %f " % cor_position + " and real CoR %f do not coincide" % self.cor_hl_pix
message = (
"Computed CoR %f " % cor_position
+ " and real CoR %f do not coincide" % self.cor_hl_pix
)
assert np.isclose(self.cor_gl_pix, cor_position, atol=self.abs_tol), message
......@@ -441,16 +506,25 @@ class TestDetectorTranslation(object):
def test_alignxc(self):
T_calc = alignment.DetectorTranslationAlongBeam()
shifts_v, shifts_h, found_shifts_list = T_calc.find_shift(self.data, self.img_pos, return_shifts=True)
shifts_v, shifts_h, found_shifts_list = T_calc.find_shift(
self.data, self.img_pos, return_shifts=True
)
message = "Computed shifts coefficients %s and expected %s do not coincide" % (
(shifts_v, shifts_h),
self.expected_shifts_vh,
)
assert np.all(np.isclose(self.expected_shifts_vh, [shifts_v, shifts_h], atol=self.abs_tol)), message
assert np.all(
np.isclose(self.expected_shifts_vh, [shifts_v, shifts_h], atol=self.abs_tol)
), message
message = "Computed shifts %s and expected %s do not coincide" % (found_shifts_list, self.all_shifts_vh)
assert np.all(np.isclose(found_shifts_list, self.all_shifts_vh, atol=self.abs_tol)), message
message = "Computed shifts %s and expected %s do not coincide" % (
found_shifts_list,
self.all_shifts_vh,
)
assert np.all(
np.isclose(found_shifts_list, self.all_shifts_vh, atol=self.abs_tol)
), message
@pytest.mark.skipif(not (__has_scipy__), reason="need scipy for this test")
def test_alignxc_synth(self):
......@@ -469,42 +543,77 @@ class TestDetectorTranslation(object):
stack, np.array([0.0, 1, 2, 3]), high_pass=1.0, return_shifts=True
)
message = "Found shifts per units %s and reference %s do not coincide" % ((shifts_v, shifts_h), (-1.234 * 2, -1.234))
assert np.all(np.isclose((shifts_v, shifts_h), (-1.234 * 2, -1.234), atol=self.abs_tol)), message
message = "Found shifts per units %s and reference %s do not coincide" % (
(shifts_v, shifts_h),
(-1.234 * 2, -1.234),
)
assert np.all(
np.isclose((shifts_v, shifts_h), (-1.234 * 2, -1.234), atol=self.abs_tol)
), message
@pytest.mark.skipif(not (__do_long_tests__), reason="need environment variable NABU_LONG_TESTS=1")
@pytest.mark.skipif(
not (__do_long_tests__), reason="need environment variable NABU_LONG_TESTS=1"
)
@pytest.mark.usefixtures("bootstrap_fcs")
class TestFocus(object):
def test_find_distance(self):
focus_calc = alignment.CameraFocus()
focus_pos, focus_ind = focus_calc.find_distance(self.data, self.img_pos)
message = "Computed focus motor position %f " % focus_pos + " and expected %f do not coincide" % self.std_best_pos
message = (
"Computed focus motor position %f " % focus_pos
+ " and expected %f do not coincide" % self.std_best_pos
)
assert np.isclose(self.std_best_pos, focus_pos, atol=self.abs_tol_dist), message
message = "Computed focus image index %f " % focus_ind + " and expected %f do not coincide" % self.std_best_ind
message = (
"Computed focus image index %f " % focus_ind
+ " and expected %f do not coincide" % self.std_best_ind
)
assert np.isclose(self.std_best_ind, focus_ind, atol=self.abs_tol_dist), message
def test_find_scintillator_tilt(self):
focus_calc = alignment.CameraFocus()
focus_pos, focus_ind, tilts_vh = focus_calc.find_scintillator_tilt(self.data, self.img_pos)
focus_pos, focus_ind, tilts_vh = focus_calc.find_scintillator_tilt(
self.data, self.img_pos
)
message = "Computed focus motor position %f " % focus_pos + " and expected %f do not coincide" % self.angle_best_pos
assert np.isclose(self.angle_best_pos, focus_pos, atol=self.abs_tol_dist), message
message = (
"Computed focus motor position %f " % focus_pos
+ " and expected %f do not coincide" % self.angle_best_pos
)
assert np.isclose(
self.angle_best_pos, focus_pos, atol=self.abs_tol_dist
), message
message = "Computed focus image index %f " % focus_ind + " and expected %f do not coincide" % self.angle_best_ind
assert np.isclose(self.angle_best_ind, focus_ind, atol=self.abs_tol_dist), message
message = (
"Computed focus image index %f " % focus_ind
+ " and expected %f do not coincide" % self.angle_best_ind
)
assert np.isclose(
self.angle_best_ind, focus_ind, atol=self.abs_tol_dist
), message
expected_tilts_vh = np.squeeze(np.array([self.angle_tilt_v, self.angle_tilt_h]))
computed_tilts_vh = -tilts_vh / (self.pixel_size / 1000)
message = "Computed tilts %s and expected %s do not coincide" % (computed_tilts_vh, expected_tilts_vh)
assert np.all(np.isclose(computed_tilts_vh, expected_tilts_vh, atol=self.abs_tol_tilt)), message
message = "Computed tilts %s and expected %s do not coincide" % (
computed_tilts_vh,
expected_tilts_vh,
)
assert np.all(
np.isclose(computed_tilts_vh, expected_tilts_vh, atol=self.abs_tol_tilt)
), message
def test_size_determination(self):
inp_shape = [2162, 2560]
exp_shape = np.array([2160, 2160])
new_shape = alignment.CameraFocus._check_img_block_size(inp_shape, 4, suggest_new_shape=True)
new_shape = alignment.CameraFocus._check_img_block_size(
inp_shape, 4, suggest_new_shape=True
)
message = "New suggested shape: %s and expected: %s do not coincide" % (new_shape, exp_shape)
message = "New suggested shape: %s and expected: %s do not coincide" % (
new_shape,
exp_shape,
)
assert np.all(new_shape == exp_shape), message
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment