From 210ae4e2f4841f170aaa3219e86f8c506cd4f89f Mon Sep 17 00:00:00 2001 From: Martijn van den Ende Date: Wed, 14 Jan 2026 17:01:04 +0100 Subject: [PATCH] Implemented ASN ROI logic for distance interpolation --- xdas/io/asn.py | 38 ++++++++++++++++++++++++++++++++++---- 1 file changed, 34 insertions(+), 4 deletions(-) diff --git a/xdas/io/asn.py b/xdas/io/asn.py index 7dc4234..7b41ba4 100644 --- a/xdas/io/asn.py +++ b/xdas/io/asn.py @@ -1,4 +1,5 @@ import json +from bisect import bisect_left import h5py import numpy as np @@ -10,16 +11,45 @@ from ..virtual import VirtualSource -def read(fname): +def read(fname: str) -> DataArray: with h5py.File(fname, "r") as file: header = file["header"] + demod = file["demodSpec"] + t0 = np.datetime64(round(header["time"][()] * 1e9), "ns") dt = np.timedelta64(round(1e9 * header["dt"][()]), "ns") - dx = header["dx"][()] * np.median(np.diff(header["channels"])) + dx = float(header["dx"][()]) # Note: dx before (internal) downsampling! data = VirtualSource(file["data"]) - nt, nx = data.shape + + # Get the optical distance for all the recorded channels (after downsampling) + # Note that this vector is not continuous for more than one ROI + all_dists = file["cableSpec"]["sensorDistances"][...] + + # Buffer for the data index at which each ROI starts/stops + dist_tie_inds = [] + # Buffer for the optical distance at which each ROI starts/stops + dist_tie_vals = [] + + # Loop over ROIs, get the start/stop index before downsampling + for n_start, n_end in zip(demod["roiStart"], demod["roiEnd"]): + # Get the index where the ROI starts based on the position in the + # distance vector. This solves the issue of rounding during decimation + i = bisect_left(all_dists, n_start * dx) + # Append the data index and optical distance to the buffers + dist_tie_inds.append(i) + dist_tie_vals.append(float(all_dists[i])) + + # Repeat the procedure for the index/distance at which the ROI ends. + # A "discontinuity" in the interpolation scheme is created in the + # following way: n_roi = [start, stop-1, stop, start, stop-1, stop, ...] + i = bisect_left(all_dists, n_end * dx) + for j in reversed(range(2)): + dist_tie_inds.append(i-j) + dist_tie_vals.append(float(all_dists[i-j])) + + nt = data.shape[0] time = {"tie_indices": [0, nt - 1], "tie_values": [t0, t0 + (nt - 1) * dt]} - distance = {"tie_indices": [0, nx - 1], "tie_values": [0.0, (nx - 1) * dx]} + distance = {"tie_indices": dist_tie_inds, "tie_values": dist_tie_vals} return DataArray(data, {"time": time, "distance": distance})