Added test to merge_double_ended_times

and added verbose kwarg to merge_double_ended

src/dtscalibration/datastore_utils.py

@@ -93,7 +93,7 @@ def check_timestep_allclose(ds, eps=0.01):
                                               'for all time steps'
 
 
-def merge_double_ended(ds_fw, ds_bw, cable_length, plot_result=True):
+def merge_double_ended(ds_fw, ds_bw, cable_length, plot_result=True, verbose=True):
     """
     Some measurements are not set up on the DTS-device as double-ended
     meausurements. This means that the two channels have to be merged manually.
@@ -123,7 +123,7 @@ def merge_double_ended(ds_fw, ds_bw, cable_length, plot_result=True):
             and ds_bw.attrs['isDoubleEnded'] == '0'), \
         "(one of the) input DataStores is already double ended"
 
-    ds_fw, ds_bw = merge_double_ended_times(ds_fw, ds_bw)
+    ds_fw, ds_bw = merge_double_ended_times(ds_fw, ds_bw, verbose=verbose)
 
     ds = ds_fw.copy()
     ds_bw = ds_bw.copy()
@@ -156,7 +156,7 @@ def merge_double_ended(ds_fw, ds_bw, cable_length, plot_result=True):
     return ds
 
 
-def merge_double_ended_times(ds_fw, ds_bw):
+def merge_double_ended_times(ds_fw, ds_bw, verbose=True):
     """Helper for `merge_double_ended()` to deal with missing measurements. The
     number of measurements of the forward and backward channels might get out
     of sync if the device shuts down before the measurement of the last channel
@@ -225,10 +225,12 @@ def merge_double_ended_times(ds_fw, ds_bw):
             pass
 
         elif direction == "fw" and direction_next == "fw":
-            print(f"Missing backward measurement beween {ds_fw.time.values[ind]} and {ds_fw.time.values[ind_next]}")
+            if verbose:
+                print(f"Missing backward measurement beween {ds_fw.time.values[ind]} and {ds_fw.time.values[ind_next]}")
 
         elif direction == "bw" and direction_next == "bw":
-            print(f"Missing forward measurement beween {ds_bw.time.values[ind]} and {ds_bw.time.values[ind_next]}")
+            if verbose:
+                print(f"Missing forward measurement beween {ds_bw.time.values[ind]} and {ds_bw.time.values[ind_next]}")
 
     return ds_fw.isel(time=iuse_chfw), ds_bw.isel(time=iuse_chbw)
 

tests/test_datastore.py

@@ -637,3 +637,60 @@ def test_merge_double_ended():
     result = (ds.isel(time=0).st - ds.isel(time=0).rst).sum().values
 
     np.testing.assert_approx_equal(result, -3712866.0382, significant=10)
+
+
+def test_merge_double_ended_times():
+    """
+    If all measurements are recorded: fw_t0, bw_t0, fw_t1, bw_t1, fw_t2, bw_t2, ..
+        > all are passed
+
+    If some are missing the accompanying measurement is skipped:
+     - fw_t0, bw_t0, bw_t1, fw_t2, bw_t2, .. > fw_t0, bw_t0, fw_t2, bw_t2, ..
+     - fw_t0, bw_t0, fw_t1, fw_t2, bw_t2, .. > fw_t0, bw_t0, fw_t2, bw_t2, ..
+     - fw_t0, bw_t0, bw_t1, fw_t2, fw_t3, bw_t3, .. > fw_t0, bw_t0, fw_t3, bw_t3,
+
+    Mixing forward and backward channels can be problematic when there is a pause
+    after measuring all channels. This function is not perfect as the following
+    situation is not caught:
+     - fw_t0, bw_t0, fw_t1, bw_t2, fw_t3, bw_t3, ..
+        > fw_t0, bw_t0, fw_t1, bw_t2, fw_t3, bw_t3, ..
+
+    This routine checks that the lowest channel
+    number is measured first (aka the forward channel), but it doesn't catch the
+    last case as it doesn't know that fw_t1 and bw_t2 belong to different cycles.
+    """
+    filepath_fw = data_dir_double_single_ch1
+    filepath_bw = data_dir_double_single_ch2
+    cable_length = 2017.7
+
+    ds_fw = read_silixa_files(directory=filepath_fw)
+    ds_bw = read_silixa_files(directory=filepath_bw)
+
+    # set stokes to varify proper time alignment
+    ds_fw.st.values = np.tile(np.arange(ds_fw.time.size)[None], (ds_fw.x.size, 1))
+    ds_bw.st.values = np.tile(np.arange(ds_bw.time.size)[None], (ds_bw.x.size, 1))
+
+    # time index that is not included in: fw, bw, merged
+    cases_iisnotin = [
+        ([], [], []),
+        ([1], [], [1]),
+        ([], [1], [1]),
+        ([1], [1], [1]),
+        ([1, 2], [], [1, 2]),
+        ([], [1, 2], [1, 2]),
+        ([1], [2], [1, 2]),
+    ]
+    # The following is not caught:
+    # ([2], [1], [1, 2]),
+    cases_iisin = [tuple(list(i for i in range(6) if i not in a) for a in b) for b in cases_iisnotin]
+
+    for ifw, ibw, iout in cases_iisin:
+        ds = merge_double_ended(
+            ds_fw.isel(time=ifw),
+            ds_bw.isel(time=ibw),
+            cable_length=cable_length,
+            plot_result=False,
+            verbose=False)
+        assert ds.time.size == len(iout) and np.all(ds.st.isel(x=0) == iout), \
+            f"FW:{ifw} & BW:{ibw} should lead to {iout} but instead leads to {ds.st.isel(x=0).values}"
+