Skip to content

Commit

Permalink
Merge pull request #41 from zxcvbnmike15/feature/area_correction
Browse files Browse the repository at this point in the history
Feature/area correction
  • Loading branch information
an-underwater-acoustician committed Feb 19, 2020
2 parents c6d5051 + ae649cc commit fca079c
Show file tree
Hide file tree
Showing 8 changed files with 283 additions and 31 deletions.
26 changes: 26 additions & 0 deletions examples/workground/reson_smith/ex_test_area.py
@@ -0,0 +1,26 @@
import logging

from pathlib import Path

from hyo2.openbst.lib.openbst import OpenBST
from hyo2.abc.lib.testing_paths import TestingPaths
from hyo2.openbst.lib.processing.process_methods.area_correction import AreaCorrectionEnum

logger = logging.getLogger(__name__)

testing = TestingPaths(root_folder=Path(__file__).parents[3].resolve())

bst = OpenBST(prj_name="test_area", force_new=True).prj

raw_path = testing.download_data_folder().joinpath('raw_reson', '20190730_144835.s7k')
bst.add_raw(raw_path)
bst.check_health()

# Test 1: Check we meet requirements = False
bst.area_correction()

# Test 2: Method: Flat Seafloor assumption
bst.raw_decode()
bst.interpolation()
bst.raytrace()
bst.area_correction()
Expand Up @@ -3,9 +3,8 @@
from pathlib import Path

from hyo2.openbst.lib.openbst import OpenBST
from hyo2.openbst.lib.processing.auxilaries.calibration import Calibration
from hyo2.abc.lib.testing_paths import TestingPaths
from hyo2.openbst.lib.processing.process_methods.radiation_pattern_compensation import RadiationPatternEnum
from hyo2.openbst.lib.processing.process_methods.calibration import CalibrationEnum


logger = logging.getLogger(__name__)
Expand All @@ -21,12 +20,11 @@
bst.add_raw(raw_path)
bst.check_health()
# Test 1: Check that we meet requirements = False
bst.parameters.calibration.method_type = RadiationPatternEnum.calibration_file
bst.parameters.calibration.method_type = CalibrationEnum.calibration_file
bst.parameters.calibration.fit_curve = True
bst.parameters.calibration.curve_order = 4
bst.calibration()

# Test 2: Method: Calibration File, Pass = No errors
bst.raw_decode()
bst.calibration()

11 changes: 8 additions & 3 deletions hyo2/openbst/lib/processing/parameters.py
@@ -1,8 +1,10 @@
import logging

from hyo2.openbst.lib.processing.process_methods.radiation_pattern_compensation import RadiationPatternParameters
from hyo2.openbst.lib.processing.process_methods.interpolation import InterpParameters
from hyo2.openbst.lib.processing.process_methods.dicts import ProcessMethods

from hyo2.openbst.lib.processing.process_methods.area_correction import AreaCorrectionParameters
from hyo2.openbst.lib.processing.process_methods.calibration import CalibrationParameters
from hyo2.openbst.lib.processing.process_methods.interpolation import InterpParameters
from hyo2.openbst.lib.processing.process_methods.raw_decoding import RawDecodeParameters
from hyo2.openbst.lib.processing.process_methods.raytracing import RayTraceParams
from hyo2.openbst.lib.processing.process_methods.source_level import SourceLevelParameters
Expand All @@ -18,7 +20,8 @@ class Parameters:
"""Class to store processing parameters"""

def __init__(self):
self.calibration = RadiationPatternParameters()
self.area_correction = AreaCorrectionParameters()
self.calibration = CalibrationParameters()
self.rawdecode = RawDecodeParameters()
self.static_gains = StaticGainParameters()
self.source_level = SourceLevelParameters()
Expand All @@ -31,6 +34,8 @@ def __init__(self):
def get_process_params(self, process_type: ProcessMethods):
if process_type == ProcessMethods.RAWDECODE:
return self.rawdecode
elif process_type == ProcessMethods.INSONIFIEDAREA:
return self.area_correction
elif process_type == ProcessMethods.CALIBRATION:
return self.calibration
elif process_type == ProcessMethods.INTERPOLATION:
Expand Down
31 changes: 21 additions & 10 deletions hyo2/openbst/lib/processing/process.py
Expand Up @@ -11,7 +11,9 @@
from hyo2.openbst.lib.processing.parameters import Parameters
from hyo2.openbst.lib.processing.process_management.process_manager import ProcessManager
from hyo2.openbst.lib.processing.process_methods.dicts import ProcessMethods
from hyo2.openbst.lib.processing.process_methods.radiation_pattern_compensation import RadiationPatternCorrection

from hyo2.openbst.lib.processing.process_methods.calibration import Calibration
from hyo2.openbst.lib.processing.process_methods.area_correction import AreaCompensation
from hyo2.openbst.lib.processing.process_methods.interpolation import Interpolation
from hyo2.openbst.lib.processing.process_methods.raw_decoding import RawDecoding
from hyo2.openbst.lib.processing.process_methods.raytracing import RayTrace
Expand Down Expand Up @@ -95,21 +97,27 @@ def run_process(self, process_method: ProcessMethods, process_file_path: Path, r
# Run the process
method_parameters = parameters.get_process_params(process_type=process_method)
if process_method is ProcessMethods.RAWDECODE:
data_out = RawDecoding.decode(ds_raw=ds_raw,
parameters=method_parameters)
data_out = RawDecoding.decode(ds_raw=ds_raw, parameters=method_parameters)

elif process_method is ProcessMethods.INSONIFIEDAREA:
data_out = AreaCompensation.area_correction(ds_process=ds_process,
ds_raw=ds_raw,
parent=self.proc_manager.parent_process,
parameters=method_parameters)

elif process_method is ProcessMethods.INTERPOLATION:
data_out = Interpolation.interpolate(ds_raw=ds_raw, parameters=method_parameters)

elif process_method is ProcessMethods.CALIBRATION:
cal_list = self.auxiliary_files.calibration_list
ds_aux = Dataset(filename=self.auxiliary_files.path, mode='r')

data_out = RadiationPatternCorrection.radiation_pattern_correction(ds_process=ds_process,
ds_raw=ds_raw,
ds_aux=ds_aux,
parent=self.proc_manager.parent_process,
parameters=method_parameters,
calibration_list=cal_list)
data_out = Calibration.radiation_pattern_correction(ds_process=ds_process,
ds_raw=ds_raw,
ds_aux=ds_aux,
parent=self.proc_manager.parent_process,
parameters=method_parameters,
calibration_list=cal_list)
ds_aux.close()

elif process_method is ProcessMethods.STATICGAIN:
Expand Down Expand Up @@ -180,8 +188,11 @@ def store_process(self, process_method: ProcessMethods, nc_process_file: Path, p
if process_method is ProcessMethods.RAWDECODE:
process_written = RawDecoding.write_data_to_nc(data_dict=data, grp_process=grp_process)

elif process_method is ProcessMethods.INSONIFIEDAREA:
process_written = AreaCompensation.write_data_to_nc(data_dict=data, grp_process=grp_process)

elif process_method is ProcessMethods.CALIBRATION:
process_written = RadiationPatternCorrection.write_data_to_nc(data_dict=data, grp_process=grp_process)
process_written = Calibration.write_data_to_nc(data_dict=data, grp_process=grp_process)

elif process_method is ProcessMethods.INTERPOLATION:
process_written = Interpolation.write_data_to_nc(data_dict=data, grp_process=grp_process)
Expand Down
200 changes: 200 additions & 0 deletions hyo2/openbst/lib/processing/process_methods/area_correction.py
@@ -0,0 +1,200 @@
import logging
import numpy as np

from enum import Enum
from netCDF4 import Dataset, Group
from typing import Optional

from hyo2.openbst.lib.nc_helper import NetCDFHelper
from hyo2.openbst.lib.processing.process_methods.raytracing import RayTraceParams
logger = logging.getLogger(__name__)


# ## Area correction enum and dictionaries ##
class AreaCorrectionEnum(Enum):
flat_seafloor = 0
local_slope = 1


area_correction_title = {
AreaCorrectionEnum.flat_seafloor: "Corrected Backscatter - Area Correction Using Flat Sealoor Assumption",
AreaCorrectionEnum.local_slope: "Corrected Backscatter - Area Correction Using Local Slope"
}


# ## Area Correction Parameters Object ##
class AreaCorrectionParameters:
process_name = "area_correction"

def __init__(self):
self.method_type = AreaCorrectionEnum.flat_seafloor

def nc_write_parameters(self, grp_process: Group):
try:
grp_process.title = area_correction_title[self.method_type]
grp_process.method_type = self.method_type.name
return True
except TypeError:
return False

def process_identifiers(self) -> list:
process_string = self.process_name
parameter_string = str()
for key, value in self.__dict__.items():
parameter_string += key + str(value)
parameter_hash = NetCDFHelper.hash_string(input_str=parameter_string)
process_ids = [process_string, parameter_hash]
return process_ids


# ## Area Compensation Class and Methods ##
class AreaCompensation:

def __init__(self):
pass

@classmethod
def area_correction(cls, ds_process: Dataset, ds_raw: Dataset,
parent: str, parameters: AreaCorrectionParameters):

p_method_type = parameters.method_type

grp_beam_geo = ds_raw.groups['beam_geometry']
var_rx_bw_across = grp_beam_geo.variables['across_beamwidth']
data_rx_bw_across = var_rx_bw_across[:]

grp_runtime = ds_raw.groups['runtime_settings']
var_tx_bw_along = grp_runtime.variables['tx_along_beam_width']
data_tx_be_along = var_tx_bw_along[:]
var_sound_speed = grp_runtime.variables['sound_velocity']
data_sound_speed = var_sound_speed[:]
var_pulse_length = grp_runtime.variables['tx_pulse_width']
data_pulse_length = var_pulse_length[:]

grp_raw = ds_raw.groups['raw_bathymetry_data']
var_rx_angle = grp_raw.variables['rx_angle']
data_rx_angle = var_rx_angle[:]

data_backscatter = cls.find_bs_data(ds_process=ds_process, parent=parent)
data_raypath = cls.find_ray_data(ds_process=ds_process, parent=parent)

if data_backscatter is None:
logger.error("Step not computed: Backscatter data not found")
return False

if data_raypath is None:
logger.error("Step not computed: Raypath data not found")
return False

if p_method_type is AreaCorrectionEnum.flat_seafloor:
data_out = cls.flat_seafloor(backscatter=data_backscatter,
rx_angle=data_rx_angle, raypath=data_raypath,
rx_beamwidth_across=data_rx_bw_across,
tx_beamwidth_along=data_tx_be_along[:, np.newaxis],
pulse_len=data_pulse_length[:, np.newaxis],
sound_speed=data_sound_speed[:, np.newaxis])

elif p_method_type is AreaCorrectionEnum.local_slope:
data_out = None
pass
else:
raise TypeError("Urecognized Method Type: %s" % p_method_type)
return data_out

@classmethod
def write_data_to_nc(cls, data_dict: dict, grp_process: Group):
try:
grp_process.createDimension(dimname='ping', size=None)
grp_process.createDimension(dimname='beam', size=None)
grp_process.createDimension(dimname='angle', size=None)

for data_name, data in data_dict.items():
if data_name == 'backscatter_data':
var_bs_data = grp_process.createVariable(varname='backscatter_data',
datatype='f8',
dimensions=('ping', 'beam'))
var_bs_data[:] = data

elif data_name == 'area_correction':
var_area_corr = grp_process.createVariable(varname='area_correction',
datatype='f8',
dimensions=('ping', 'beam'))
var_area_corr[:] = data

elif data_name == 'area_corr_beam_limited':
var_area_bl = grp_process.createVariable(varname='beam_limited_area_correction',
datatype='f8',
dimensions=('ping', 'beam'))
var_area_bl[:] = data
elif data_name == 'area_corr_pulse_limited':
var_area_pl = grp_process.createVariable(varname='pulse_limited_area_correction',
datatype='f8',
dimensions=('ping', 'beam'))
var_area_pl[:] = data

return True
except RuntimeError:
return False

@classmethod
def flat_seafloor(cls,
backscatter: np.ma.MaskedArray, rx_angle: np.ma.MaskedArray, raypath: np.ma.MaskedArray,
rx_beamwidth_across: np.ndarray, tx_beamwidth_along: np.ndarray,
pulse_len: np.ndarray, sound_speed: np.ndarray):

rx_beamwidth_across = np.deg2rad(rx_beamwidth_across)
tx_beamwidth_along = np.deg2rad(tx_beamwidth_along)
rx_angle = np.deg2rad(rx_angle)

area_beam_lim = rx_beamwidth_across * tx_beamwidth_along * raypath ** 2
area_pulse_lim = ((sound_speed * pulse_len) / (2 * np.sin(np.abs(rx_angle)))) * (tx_beamwidth_along * raypath)

area_correction = 10 * np.ma.log10(np.minimum(area_beam_lim, area_pulse_lim))

bs_corrected = backscatter - area_correction

data_out = {
'backscatter_data': bs_corrected,
'area_correction': area_correction,
'area_corr_beam_limited': 10 * np.ma.log10(area_beam_lim),
'area_corr_pulse_limited': 10 * np.ma.log10(area_pulse_lim)
}

return data_out

@classmethod
def local_slope(cls):
pass

@classmethod
def find_bs_data(cls, ds_process: Dataset, parent: str):
grp_parent = ds_process.groups[parent]
try:
var_backscatter = grp_parent.variables['backscatter_data']
data_backscatter = var_backscatter[:]
return data_backscatter
except KeyError:
parent = grp_parent.parent_process
if parent == 'ROOT':
return None

data_backscatter = cls.find_bs_data(ds_process=ds_process, parent=parent)
return data_backscatter

@classmethod
def find_ray_data(cls, ds_process: Dataset, parent: str):
grp_parent = ds_process.groups[parent]

if RayTraceParams.process_name in parent:
var_ray_path = grp_parent.variables['path_length']
data_ray_path = var_ray_path[:]
return data_ray_path

elif grp_parent.parent_process == 'ROOT':
return None

else:
parent = grp_parent.parent_process
data_ray_path = cls.find_ray_data(ds_process=ds_process, parent=parent)
return data_ray_path

0 comments on commit fca079c

Please sign in to comment.