fix: match qlsi 2d with new 3d implementation

qpretrieve/interfere/if_qlsi.py

@@ -175,8 +175,8 @@ def run_pipeline(self, **pipeline_kws):
         # need to do this along the z axis, as skimage `unwrap_3d` does not
         # work for our use-case
         # todo: maybe use np.unwrap for the xy axes instead
-        px = np.zeros_like(hx)
-        py = np.zeros_like(hy)
+        px = np.zeros_like(hx, dtype=float)
+        py = np.zeros_like(hy, dtype=float)
         for i, (_hx, _hy) in enumerate(zip(hx, hy)):
             px[i] = unwrap_phase(np.angle(_hx))
             py[i] = unwrap_phase(np.angle(_hy))
@@ -210,10 +210,10 @@ def run_pipeline(self, **pipeline_kws):
                         copy=False)
         # Compute the frequencies that correspond to the frequencies of the
         # Fourier-transformed image.
-        fx = np.fft.fftfreq(rfft.shape[-1]).reshape(rfft.shape[0], -1, 1)
-        fy = np.fft.fftfreq(rfft.shape[-2]).reshape(rfft.shape[0], 1, -1)
+        fx = np.fft.fftfreq(rfft.shape[-2]).reshape(rfft.shape[0], -1, 1)
+        fy = np.fft.fftfreq(rfft.shape[-1]).reshape(rfft.shape[0], 1, -1)
         fxy = -2 * np.pi * 1j * (fx + 1j * fy)
-        fxy[0, 0] = 1
+        fxy[:, 0, 0] = 1
 
         # The wavefront is the real part of the inverse Fourier transform
         # of the filtered (divided by frequencies) data.

tests/test_fourier_base.py

@@ -137,7 +137,7 @@ def test_scale_to_filter_qlsi():
 
     ifh = interfere.QLSInterferogram(image, **pipeline_kws)
     raw_wavefront = ifh.run_pipeline()
-    assert raw_wavefront.shape == (720, 720)
+    assert raw_wavefront.shape == (1, 720, 720)
     assert ifh.phase.shape == (1, 720, 720)
     assert ifh.amplitude.shape == (1, 720, 720)
     assert ifh.field.shape == (1, 720, 720)
@@ -163,6 +163,6 @@ def test_scale_to_filter_qlsi():
     ifr.run_pipeline(**pipeline_kws_scale)
     phase_scaled = unwrap_phase(ifh.phase - ifr.phase)
 
-    assert phase_scaled.shape == (126, 126)
+    assert phase_scaled.shape == (1, 126, 126)
 
     assert np.allclose(phase_scaled.mean(), 0.1257080793074251, atol=1e-6)

tests/test_qlsi.py

@@ -118,8 +118,16 @@ def test_qlsi_rotate_2d_3d(hologram):
         axes=(-1, -2),  # the y and x axes
         reshape=False,
     )
+    rot_3d_2 = qpretrieve.interfere.if_qlsi.rotate_noreshape(
+        data_3d,
+        angle=2,
+        axes=(-2, -1),  # the y and x axes
+        reshape=False,
+    )
 
+    assert rot_2d.dtype == rot_3d.dtype
     assert np.array_equal(rot_2d, rot_3d[0])
+    assert np.array_equal(rot_2d, rot_3d_2[0])
 
 
 def test_qlsi_pad_2d_3d(hologram):
@@ -134,4 +142,27 @@ def test_qlsi_pad_2d_3d(hologram):
         data_3d, ((0, 0), (sx // 2, sx // 2), (sy // 2, sy // 2)),
         mode="constant", constant_values=0)
 
+    assert gradpad_2d.dtype == gradpad_3d.dtype
     assert np.array_equal(gradpad_2d, gradpad_3d[0])
+
+
+def test_fxy_complex_mul(hologram):
+    data_2d = hologram
+    data_3d, _ = qpretrieve.data_input._convert_2d_to_3d(data_2d)
+
+    assert np.array_equal(data_2d, data_3d[0])
+
+    # 2d
+    fx_2d = data_2d.reshape(-1, 1)
+    fy_2d = data_2d.reshape(1, -1)
+    fxy_2d = -2 * np.pi * 1j * (fx_2d + 1j * fy_2d)
+    fxy_2d[0, 0] = 1
+
+    # 3d
+    fx_3d = data_3d.reshape(data_3d.shape[0], -1, 1)
+    fy_3d = data_3d.reshape(data_3d.shape[0], 1, -1)
+    fxy_3d = -2 * np.pi * 1j * (fx_3d + 1j * fy_3d)
+    fxy_3d[:, 0, 0] = 1
+
+    assert np.array_equal(fx_2d, fx_3d[0])
+    assert np.array_equal(fxy_2d, fxy_3d[0])