diff --git a/examples/workground/reson_smith/ex_test_area.py b/examples/workground/reson_smith/ex_test_area.py new file mode 100644 index 0000000..80a4def --- /dev/null +++ b/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() diff --git a/examples/workground/reson_smith/ex_test_calibration_compensation.py b/examples/workground/reson_smith/ex_test_calibration_compensation.py index e9b977f..c532341 100644 --- a/examples/workground/reson_smith/ex_test_calibration_compensation.py +++ b/examples/workground/reson_smith/ex_test_calibration_compensation.py @@ -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__) @@ -21,7 +20,7 @@ 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() @@ -29,4 +28,3 @@ # Test 2: Method: Calibration File, Pass = No errors bst.raw_decode() bst.calibration() - diff --git a/hyo2/openbst/lib/processing/parameters.py b/hyo2/openbst/lib/processing/parameters.py index ad3cf6a..123981b 100644 --- a/hyo2/openbst/lib/processing/parameters.py +++ b/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 @@ -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() @@ -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: diff --git a/hyo2/openbst/lib/processing/process.py b/hyo2/openbst/lib/processing/process.py index 0e9b093..7e5c6b6 100644 --- a/hyo2/openbst/lib/processing/process.py +++ b/hyo2/openbst/lib/processing/process.py @@ -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 @@ -95,8 +97,14 @@ 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) @@ -104,12 +112,12 @@ def run_process(self, process_method: ProcessMethods, process_file_path: Path, r 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: @@ -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) diff --git a/hyo2/openbst/lib/processing/process_methods/area_correction.py b/hyo2/openbst/lib/processing/process_methods/area_correction.py new file mode 100644 index 0000000..93a0e45 --- /dev/null +++ b/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 + diff --git a/hyo2/openbst/lib/processing/process_methods/radiation_pattern_compensation.py b/hyo2/openbst/lib/processing/process_methods/calibration.py similarity index 89% rename from hyo2/openbst/lib/processing/process_methods/radiation_pattern_compensation.py rename to hyo2/openbst/lib/processing/process_methods/calibration.py index 023b872..9e00b84 100644 --- a/hyo2/openbst/lib/processing/process_methods/radiation_pattern_compensation.py +++ b/hyo2/openbst/lib/processing/process_methods/calibration.py @@ -11,30 +11,30 @@ # ## Calibration Enum and Dictionaries ## -class RadiationPatternEnum(Enum): +class CalibrationEnum(Enum): calibration_file = 0 custom_curve = 1 -radiation_pattern_title = { - RadiationPatternEnum.calibration_file: "Corrected Backscatter - Compensated Radiation Pattern using Calibration " +calibration_title = { + CalibrationEnum.calibration_file: "Corrected Backscatter - Compensated Radiation Pattern using Calibration " "File", - RadiationPatternEnum.custom_curve: "Corrected Backscatter - Compensated Radiation Pattern using a Fitted Curve" + CalibrationEnum.custom_curve: "Corrected Backscatter - Compensated Radiation Pattern using a Fitted Curve" } # ## Radiation Pattern Parameters Object ## -class RadiationPatternParameters: +class CalibrationParameters: process_name = "radiation_pattern_compensation" def __init__(self): - self.method_type = RadiationPatternEnum.calibration_file + self.method_type = CalibrationEnum.calibration_file self.fit_curve = True self.curve_order = 2 def nc_write_parameters(self, grp_process: Group): try: - grp_process.title = radiation_pattern_title[self.method_type] + grp_process.title = calibration_title[self.method_type] grp_process.method_type = self.method_type.name grp_process.fit_curve = str(self.fit_curve) grp_process.curve_order = self.curve_order @@ -43,7 +43,7 @@ def nc_write_parameters(self, grp_process: Group): return False def process_identifiers(self) -> list: - process_string = RadiationPatternParameters.process_name + process_string = CalibrationParameters.process_name parameter_string = str() for key, value in self.__dict__.items(): parameter_string += key + str(value) @@ -53,18 +53,18 @@ def process_identifiers(self) -> list: # Radiation Pattern Class and Methods ## -class RadiationPatternCorrection: +class Calibration: def __init__(self): pass @classmethod def radiation_pattern_correction(cls, ds_process: Dataset, ds_raw: Dataset, ds_aux: Dataset, - parent: str, parameters: RadiationPatternParameters, + parent: str, parameters: CalibrationParameters, calibration_list: Optional[list] = None): p_method_type = parameters.method_type - if p_method_type is RadiationPatternEnum.calibration_file: + if p_method_type is CalibrationEnum.calibration_file: # Check if there is a calibration file if len(calibration_list) < 1: @@ -97,7 +97,7 @@ def radiation_pattern_correction(cls, ds_process: Dataset, ds_raw: Dataset, ds_a cal_curve_angles=data_angles, cal_curve_values=data_cal, fit_curve=p_fit_curve, curve_order=p_curve_order) - elif p_method_type is RadiationPatternEnum.custom_curve: + elif p_method_type is CalibrationEnum.custom_curve: data_out = None pass else: @@ -175,4 +175,3 @@ def calibration_file(cls, backscatter_data: np.ma.MaskedArray, @classmethod def custom_curve(cls): pass - \ No newline at end of file diff --git a/hyo2/openbst/lib/processing/process_methods/dicts.py b/hyo2/openbst/lib/processing/process_methods/dicts.py index 28856e1..ca98727 100644 --- a/hyo2/openbst/lib/processing/process_methods/dicts.py +++ b/hyo2/openbst/lib/processing/process_methods/dicts.py @@ -35,7 +35,7 @@ class ProcessMethods(Enum): ProcessMethods.GEOLOCATION: ProcessMethods.RAYTRACE, ProcessMethods.GRIDDING: ProcessMethods.GEOLOCATION, ProcessMethods.TRANSMISSIONLOSS: ProcessMethods.RAYTRACE, - ProcessMethods.INSONIFIEDAREA: ProcessMethods.INTERPOLATION, + ProcessMethods.INSONIFIEDAREA: ProcessMethods.RAYTRACE, ProcessMethods.ANGULARDEPENDENCE: ProcessMethods.INSONIFIEDAREA, ProcessMethods.DESPECKLE: None, ProcessMethods.ANTIALIASING: None diff --git a/hyo2/openbst/lib/project.py b/hyo2/openbst/lib/project.py index 396ba19..93fb5f0 100644 --- a/hyo2/openbst/lib/project.py +++ b/hyo2/openbst/lib/project.py @@ -264,6 +264,19 @@ def remove_calibration(self, path: Path): # ### Process ### # TODO: Significant code duplication below, need to reduce this # TODO: Add error handling in the event we havent added a raw yet and are calling these functions + def area_correction(self): + for path_hash in self.raws.raws_list: + raw_file_path = self.raws_folder.joinpath(path_hash + self.raws.ext) + process_file_path = self.process_folder.joinpath(path_hash + self.process.ext) + + processed = self.process.run_process(process_method=self.process.process_method_types.INSONIFIEDAREA, + process_file_path=process_file_path, + raw_path=raw_file_path, + parameters=self.parameters) + if processed is True: + self.info.manage_parent(parent=self.process.proc_manager.parent_process) + print('File corrected for radiation pattern correction: %s' % process_file_path.resolve()) + def calibration(self): for path_hash in self.raws.raws_list: raw_file_path = self.raws_folder.joinpath(path_hash + self.raws.ext)