Implemented the computation of the VBI on a common grid of 'r' values for the two selected frequencies. Exporting to a .csv file and updating the sand SSC plots still need work.

dev-bjvincent-fix-UBSediFlow
Bjarne Vincent 2025-10-21 18:27:15 +02:00
parent 5e6f2e0d69
commit fc438ee3a6
3 changed files with 105 additions and 85 deletions

View File

@ -1718,6 +1718,9 @@ class AcousticDataTab(QWidget):
stg.kt2D.append(np.array([]))
stg.kt3D.append(np.array([]))
stg.J_cross_section.append([np.array([]), np.array([])])
stg.J_cross_section_interp.append( [ np.array([]), np.array([]) ] )
stg.r2D_interp.append(np.array([]))
stg.t2D_interp.append(np.array([]))
stg.VBI_cross_section.append(np.array([]))
stg.SSC_fine.append(np.array([]))
stg.SSC_sand.append(np.array([]))

View File

@ -25,6 +25,7 @@ import gc
import time
import logging
import numpy as np
import scipy as sp
import pandas as pd
from math import isinf
@ -529,19 +530,85 @@ class AcousticInversionTab(QWidget):
self.plot_measured_vs_inverted_SSC_sand()
def compute_VBI(self):
'''
Compute the VBI (Volume Backscatter Index) for the selected dataset in the combobox.
If the data has been recorded by an Ubertone UB SediFlow, the data for the two
calibration frequencies is first interpolated on the same spatial and temporal
grids.
'''
# Get the dataset ID (integer relating the dataset selected in the combobox
# to the corresponding dataset loaded into AcouSed):
data_id = self.combobox_acoustic_data_choice.currentIndex()
# Get the spatial and temporal grids for the two selected frequencies, and create common
# grids to be used for the VBI and SSC computations (reminder: the grids have
# the same number of points, but the intervals are different). Here, we assume that the
# time delay between the channels (i.e., frequencies) is negligible, so that the i-th
# time stamp is taken as the average between the i-th time stamp of the first frequency
# and the i-th time stamp of the second frequency:
rGrid_freq1 = stg.depth_2D[data_id][ stg.frequencies_for_calibration[0][1] ]
rGrid_freq2 = stg.depth_2D[data_id][ stg.frequencies_for_calibration[1][1] ]
rVect_freq1 = rGrid_freq1[:,0] #Vector of vertical coordinates (m) for the first frequency
rVect_freq2 = rGrid_freq2[:,0] #Vector of vertical coordinates (m) for the second frequency
tVect_freq1 = stg.time[data_id][ stg.frequencies_for_calibration[0][1] ]
tVect_freq2 = stg.time[data_id][ stg.frequencies_for_calibration[1][1] ]
rMin = np.min( np.array( [ np.min(rVect_freq1), np.min(rVect_freq2) ] ) ) #Minimum vertical position (m)
rMax = np.max( np.array( [ np.max(rVect_freq1), np.max(rVect_freq2) ] ) ) #Maximum vertical position (m)
rVect_interp = np.linspace( rMin, rMax, rVect_freq1.shape[0] ) #Vector of vertical coordinates (m) used for interpolation
tVect_interp = ( tVect_freq1 + tVect_freq2 ) / 2 #Vector of time stamps (s) used for interpolation
rGrid_interp = np.zeros( ( rVect_interp.shape[0], rGrid_freq1.shape[1] ) )
tGrid_interp = np.zeros( rGrid_interp.shape )
for j in range( rGrid_interp.shape[1] ):
rGrid_interp[:, j] = rVect_interp
tGrid_interp[:, j] = tVect_interp[j] * np.ones( rGrid_interp.shape[0] )
stg.r2D_interp[data_id] = rGrid_interp
stg.t2D_interp[data_id] = tGrid_interp
# Interpolate on the same grid the quantity 'J' that has been previously
# computed for both frequencies:
j_interpolator_freq1 = sp.interpolate.RegularGridInterpolator( ( rVect_freq1, tVect_interp ),
stg.J_cross_section[data_id][0],
method="linear", fill_value=np.nan )
j_interpolator_freq2 = sp.interpolate.RegularGridInterpolator( ( rVect_freq2, tVect_interp ),
stg.J_cross_section[data_id][1],
method="linear", fill_value=np.nan )
j_interp_freq1 = j_interpolator_freq1( ( stg.r2D_interp[data_id], stg.t2D_interp[data_id] ) )
j_interp_freq2 = j_interpolator_freq2( ( stg.r2D_interp[data_id], stg.t2D_interp[data_id] ) )
stg.J_cross_section_interp[data_id] = [ j_interp_freq1, j_interp_freq2 ]
# Compute the Volume Backscatter Index (VBI) on the common grid:
stg.VBI_cross_section[data_id] = np.array([])
stg.VBI_cross_section[data_id] = self.inv_hc.VBI_cross_section(
freq1=stg.frequencies_for_calibration[0][0],
freq2=stg.frequencies_for_calibration[1][0],
zeta_freq1=stg.zeta[0], zeta_freq2=stg.zeta[1],
j_cross_section_freq1=stg.J_cross_section[data_id][0],
j_cross_section_freq2=stg.J_cross_section[data_id][1],
r2D=stg.depth_2D[data_id][
stg.frequency_for_inversion[1]
],
j_cross_section_freq1=stg.J_cross_section_interp[data_id][0],
j_cross_section_freq2=stg.J_cross_section_interp[data_id][1],
r2D=stg.r2D_interp[data_id],
water_attenuation_freq1=stg.alpha_s[0],
water_attenuation_freq2=stg.alpha_s[1],
X=stg.X_exponent[0]
@ -555,39 +622,24 @@ class AcousticInversionTab(QWidget):
'''
stg.SSC_fine[data_id] = self.inv_hc.SSC_fine(
zeta=stg.zeta[0],
r2D=stg.depth_2D[data_id][stg.frequencies_for_calibration[0][1]],
r2D=stg.r2D_interp[data_id],
VBI=stg.VBI_cross_section[data_id],
freq=stg.frequencies_for_calibration[0][0],
X=stg.X_exponent[0],
j_cross_section=stg.J_cross_section[data_id][0],
alpha_w=np.full(
shape=stg.depth_2D[data_id][
stg.frequencies_for_calibration[0][1]
].shape,
fill_value=stg.water_attenuation[data_id][
stg.frequencies_for_calibration[0][1]
]
)
alpha_w=stg.water_attenuation[data_id][ stg.frequencies_for_calibration[0][1] ]
) #Inversion using the first frequency
'''
stg.SSC_fine[data_id] = self.inv_hc.SSC_fine(
zeta=stg.zeta[1],
r2D=stg.depth_2D[data_id][stg.frequency_for_inversion[1]],
r2D=stg.r2D_interp[data_id],
VBI=stg.VBI_cross_section[data_id],
freq=stg.frequencies_for_calibration[1][0],
X=stg.X_exponent[0],
j_cross_section=stg.J_cross_section[data_id][1],
alpha_w=np.full(
shape=stg.depth_2D[data_id][
stg.frequency_for_inversion[1]
].shape,
fill_value=stg.water_attenuation[data_id][
stg.frequency_for_inversion[1]
]
)
j_cross_section=stg.J_cross_section_interp[data_id][1],
alpha_w=stg.water_attenuation[data_id][ stg.frequency_for_inversion[1] ]
) #Inversion using the second frequency
@ -658,9 +710,9 @@ class AcousticInversionTab(QWidget):
time_data = stg.time
depth_data = stg.depth
pcm_SSC_fine = self.axis_SSC_fine.pcolormesh(
time_data[data_id][stg.frequency_for_inversion[1]],
-depth_data[data_id][stg.frequency_for_inversion[1]],
stg.t2D_interp[data_id], -stg.r2D_interp[data_id],
stg.SSC_fine[data_id],
cmap='rainbow', norm=LogNorm(vmin=1e0, vmax=15),
shading='gouraud'
@ -668,21 +720,15 @@ class AcousticInversionTab(QWidget):
if stg.depth_bottom[data_id].shape != (0,):
self.axis_SSC_fine.plot(
time_data[data_id][stg.frequency_for_inversion[1]],
stg.t2D_interp[data_id][0,:],
-stg.depth_bottom[data_id],
color='black', linewidth=1, linestyle="solid"
)
self.pcm_SSC_fine_vertical_line, = self.axis_SSC_fine.plot(
time_data[data_id][stg.frequency_for_inversion[1],
self.slider_fine.value() - 1]
* np.ones(
depth_data[data_id][
stg.frequency_for_inversion[1]
].shape
),
-depth_data[data_id][stg.frequency_for_inversion[1]],
linestyle="solid", color='r', linewidth=2
stg.t2D_interp[data_id][0, self.slider_fine.value() - 1]
* np.ones( ( stg.r2D_interp[data_id].shape[0], 1 ) ),
- stg.r2D_interp[data_id][:,0], linestyle="solid", color='r', linewidth=2
)
# --- Plot samples of fine sediments ---
@ -789,36 +835,19 @@ class AcousticInversionTab(QWidget):
self.plot_fine , = self.axis_vertical_profile_SSC_fine.plot(
stg.SSC_fine[data_id][:, self.slider_fine.value() -1],
-depth_data[data_id][stg.frequency_for_inversion[1]],
-stg.r2D_interp[data_id][:,0],
linestyle="solid", linewidth=1, color="k"
)
self.pcm_SSC_fine_vertical_line.set_data(
time_data[data_id][
stg.frequency_for_inversion[1],
self.slider_fine.value() - 1
]
* np.ones(
depth_data[data_id][
stg.frequency_for_inversion[1]
].shape
),
-depth_data[data_id][
stg.frequency_for_inversion[1]
]
)
stg.t2D_interp[data_id][0, self.slider_fine.value() - 1]
* np.ones( ( stg.r2D_interp[data_id].shape[0], 1 ) ),
-stg.r2D_interp[data_id][:,0] )
self.figure_SSC_fine.canvas.draw_idle()
self.axis_vertical_profile_SSC_fine.set_ylim(
[
-np.nanmax(
depth_data[data_id][stg.frequency_for_inversion[1]]
),
-np.nanmin(
depth_data[data_id][stg.frequency_for_inversion[1]]
)
]
)
[ - np.nanmax( stg.r2D_interp[data_id][:,0] ), -np.nanmin( stg.r2D_interp[data_id][:,0] ) ] )
self.axis_vertical_profile_SSC_fine.set_xlim(0, 10)
self.axis_vertical_profile_SSC_fine.set_xlabel("Inverted Fine SSC (g/L)")
@ -845,30 +874,19 @@ class AcousticInversionTab(QWidget):
self.axis_vertical_profile_SSC_fine.plot(
stg.SSC_fine[data_id][:, self.slider_fine.value() - 1],
-depth_data[data_id][stg.frequency_for_inversion[1]],
-stg.r2D_interp[data_id][:,0],
linestyle="solid", linewidth=1, color="k"
)
self.pcm_SSC_fine_vertical_line.set_data(
time_data[data_id][stg.frequency_for_inversion[1],
self.slider_fine.value() - 1]
* np.ones(depth_data[data_id][
stg.frequency_for_inversion[1]
].shape),
-depth_data[data_id][stg.frequency_for_inversion[1]]
stg.t2D_interp[data_id][0, self.slider_fine.value() - 1]
* np.ones( ( stg.r2D_interp[data_id].shape[0], 1 ) ),
-stg.r2D_interp[data_id][:,0]
)
self.figure_SSC_fine.canvas.draw_idle()
self.axis_vertical_profile_SSC_fine.set_ylim(
[
-np.nanmax(
depth_data[data_id][stg.frequency_for_inversion[1]]
),
-np.nanmin(
depth_data[data_id][stg.frequency_for_inversion[1]]
)
]
)
[ - np.nanmax( stg.r2D_interp[data_id][:,0] ), -np.nanmin( stg.r2D_interp[data_id][:,0] ) ] )
for f in stg.fine_sample_profile:
time_fine_calibration = (
@ -1218,8 +1236,7 @@ class AcousticInversionTab(QWidget):
depth_data = stg.depth
pcm_SSC_sand = self.axis_SSC_sand.pcolormesh(
time_data[data_id][stg.frequency_for_inversion[1]],
-depth_data[data_id][stg.frequency_for_inversion[1]],
stg.t2D_interp[data_id], -stg.r2D_interp[data_id],
stg.SSC_sand[data_id],
cmap='rainbow', norm=LogNorm(vmin=1e0, vmax=10),
shading='gouraud'
@ -1227,18 +1244,15 @@ class AcousticInversionTab(QWidget):
if stg.depth_bottom[data_id].shape != (0,):
self.axis_SSC_sand.plot(
time_data[data_id][stg.frequency_for_inversion[1]],
stg.t2D_interp[data_id][0,:],
-stg.depth_bottom[data_id],
color='black', linewidth=1, linestyle="solid"
)
self.pcm_SSC_sand_vertical_line, = self.axis_SSC_sand.plot(
time_data[data_id][stg.frequency_for_inversion[1],
self.slider_sand.value() - 1]
* np.ones(
depth_data[data_id][stg.frequency_for_inversion[1]].shape
),
-depth_data[data_id][stg.frequency_for_inversion[1]],
stg.t2D_interp[data_id][0, self.slider_sand.value() - 1]
* np.ones( ( stg.r2D_interp[data_id].shape[0], 1 ) ),
-stg.r2D_interp[data_id][:,0],
linestyle="solid", color='r', linewidth=2
)

View File

@ -260,6 +260,9 @@ lin_reg = [] # (lin_reg_compute.slope, lin_reg_compute.interc
# Variables names # Description # Type
r2D_interp = [] # Spatial grid on which the BS (through 'J') is interpolated # List of 2D arrays
t2D_interp = [] # Temporal grid on which the BS (through 'J') is inteprolated # List of 2D arrays
J_cross_section_interp = [] # List of 'J_cross_section' interpolated on 'r2D_interp' # List of list of arrays
VBI_cross_section = [] # Volume Backscattering Index # List of 2D array
fine_sample_position = [] # Fine samples indexes position for time and depth [(time_index, depth_index)] # List of tuples
sand_sample_position = [] # Sand samples indexes position for time and depth [(time_index, depth_index)] # List of tuples