diff --git a/src/dtscalibration/datastore_utils.py b/src/dtscalibration/datastore_utils.py index 77b8747d..1e1d3a79 100644 --- a/src/dtscalibration/datastore_utils.py +++ b/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) diff --git a/tests/test_datastore.py b/tests/test_datastore.py index d27e01fa..5bafaa22 100644 --- a/tests/test_datastore.py +++ b/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}" +