USML
proploss Class Reference
Collaboration diagram for proploss:

Detailed Description

Container for a list of targets and their associated propagation data.

Passing an object of this type to a wavefront object causes it to accumulates acoustic eigenrays at each location. After propagation is complete, the sum_eigenrays() method is used to collect the results into a phasor-summed propagation loss and phase at each target point.

Constructor & Destructor Documentation

proploss ( const wposition targets)

Initialize the acoustic propagation effects associated with each target.

Parameters
targetsGrid of targets to ensonify.
proploss ( const seq_vector frequencies,
const wposition1 source_pos,
const seq_vector source_de,
const seq_vector source_az,
double  time_step,
const wposition targets 
)

Initialize with references to wave front information.

Constructor Initialize the acoustic propagation effects associated with each target.

Parameters
frequenciesFrequencies over which to compute loss (Hz).
source_posLocation of the wavefront source.
source_deLaunch D/E angle at source (deg)
source_azLaunch AZ angle at source (deg)
time_stepPropagation step size (seconds).
targetsGrid of targets to ensonify.

Initialize with references to wave front information.

Member Function Documentation

bool addEigenray ( unsigned  targetRow,
unsigned  targetCol,
eigenray  pRay 
)
virtual

addEigenray - Adds an eigenray to the eigenray_list for the target specified.

Add eigenray via proplossListener.

implementation of the pure virtual method of proplossListener.

Parameters
targetRowRow number of the current target.
targetColColumn number of the current target.
pRayThe eigenray to add.
Returns
True on success, false on failure.

Implements proplossListener.

eigenray_list* eigenrays ( unsigned  t1,
unsigned  t2 
)
inline

Return eigenray list for a single target.

Parameters
t1Row number of the current target.
t2Column number of the current target.
Returns
eigenray_list Pointer to eigenray_list for single target.
const seq_vector* frequencies ( ) const
inline

Frequencies over which propagation is computed (Hz).

Returns
seq_vector Pointer to vector containing frequencies.
void initialize ( )
private

Initialize with references to wave front information.

_loss data structure.

wposition1 position ( unsigned  t1,
unsigned  t2 
)
inline

Position of a single target in the grid.

Parameters
t1Row number of the current target.
t2Column number of the current target.
Returns
wposition1 wposition1 of target.
unsigned size1 ( ) const
inline

Number of rows in target grid.

Returns
Number of rows.
unsigned size2 ( ) const
inline

Number of columns in target grid.

Returns
Number of columns.
void sum_eigenrays ( bool  coherent = true)

Compute propagation loss summed over all eigenrays.

Parameters
coherentCompute coherent propagation loss if true, and incoherent if false.
const eigenray* total ( unsigned  t1,
unsigned  t2 
)
inline

Propagation loss for a single target summed over eigenrays.

Includes eigenray element weighted averages.

Parameters
t1Row number of the current target.
t2Column number of the current target.
Returns
eigenray* Pointer to eigenray class.
void write_netcdf ( const char *  filename,
const char *  long_name = NULL 
)

Write proploss scenario data to a netCDF file using a ragged array structure.

Write proploss data to to netCDF file.

This ragged array concept (see reference) stores the proploss data in a one dimensional list and uses an externally defined index to lookup the appropriate elements for each target.

This ragged array concept is used to define the intensity, phase, source_de, source_az, target_de, target_az, surface, bottom, and caustic variables. The proploss_index variable defines the lookup index into these arrays for the summed proploss for each target. The eigenray_index variable defines a similar index for the beginning of the eigenray list. Subsequent eigenrays for this target immediately follow the 1st eigenray. The eigenray_number variable defines the number of eigenrays for each target.

The user is responsible for ensuring that sum_eigenrays() has been called prior to this routine.

This file structure is illustrated (for a single target with direct path, surface, and bottom eigenrays) in the netCDF sample below:

    netcdf eigenray_basic {
    dimensions:
        frequency = 1 ;
        rows = 1 ;
        cols = 1 ;
        eigenrays = 4 ;
        launch_de = 25 ;
        launch_az = 5 ;
    variables:
        double source_latitude ;
                source_latitude:units = "degrees_north" ;
        double source_longitude ;
                source_longitude:units = "degrees_east" ;
        double source_altitude ;
                source_altitude:units = "meters" ;
                source_altitude:positive = "up" ;
        double launch_de(launch_de) ;
                launch_de:units = "degrees" ;
                launch_de:positive = "up" ;
        double launch_az(launch_az) ;
                launch_az:units = "degrees_true" ;
                launch_az:positive = "clockwise" ;
        double time_step ;
                time_step:units = "seconds" ;
        double frequency(frequency) ;
                frequency:units = "hertz" ;
        double latitude(rows, cols) ;
                latitude:units = "degrees_north" ;
        double longitude(rows, cols) ;
                longitude:units = "degrees_east" ;
        double altitude(rows, cols) ;
                altitude:units = "meters" ;
                altitude:positive = "up" ;
        short proploss_index(rows, cols) ;
                proploss_index:units = "count" ;
        short eigenray_index(rows, cols) ;
                eigenray_index:units = "count" ;
        short eigenray_num(rows, cols) ;
                eigenray_num:units = "count" ;
        double intensity(eigenrays, frequency) ;
                intensity:units = "dB" ;
        double phase(eigenrays, frequency) ;
                phase:units = "radians" ;
        double travel_time(eigenrays) ;
                travel_time:units = "seconds" ;
        double source_de(eigenrays) ;
                source_de:units = "degrees" ;
                source_de:positive = "up" ;
        double source_az(eigenrays) ;
                source_az:units = "degrees_true" ;
                source_az:positive = "clockwise" ;
        double target_de(eigenrays) ;
                target_de:units = "degrees" ;
                target_de:positive = "up" ;
        double target_az(eigenrays) ;
                target_az:units = "degrees_true" ;
                target_az:positive = "clockwise" ;
        short surface(eigenrays) ;
                surface:units = "count" ;
        short bottom(eigenrays) ;
                bottom:units = "count" ;
        short caustic(eigenrays) ;
                caustic:units = "count" ;
    // global attributes:
                :long_name = "eigenray_basic test" ;
                :Conventions = "COARDS" ;
    data:
     source_latitude = 45 ;
     source_longitude = -45 ;
     source_altitude = -1000 ;
     launch_de = -60, -55, -50, -45, -40, -35, -30, -25, -20, -15, -10, -5, 0, 5, 
        10, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60 ;
     launch_az = -2, -1, 0, 1, 2 ;
     time_step = 0.1 ;
     frequency = 100000 ;
     latitude = 45.02 ;
     longitude = -45 ;
     altitude = -1000 ;
     proploss_index = 0 ;
     eigenray_index = 1 ;
     eigenray_num = 3 ;
     intensity = 64.1427331557955, 66.9694875322853, 69.5617366839523, 77.9679676720843 ;
     phase = 0.082693193293971, 0, -3.14159265358979, 0 ;
     travel_time = 1.88788322659623, 1.48401881064351, 1.99622641373789, 3.03542140949814 ;
     source_de = 6.86829906648868, -0.00988717967364023, 41.9270555318522, -61.0113452826929 ;
     source_az = 0, 0, 0, 0 ;
     target_de = -6.86819212315001, 0.0101128206074365, -41.9270291816411, 61.0112432784064 ;
     target_az = 0, 0, 0, 0 ;
     surface = 0, 0, 1, 0 ;
     bottom = 0, 0, 0, 1 ;
     caustic = 0, 0, 0, 0 ;
    }
Parameters
filenameName of the file to write to disk.
long_nameOptional global attribute for identifying data-set.
References:
"The NetCDF Users Guide - Data Model, Programming Interfaces, and Format for Self-Describing, Portable Data NetCDF", Version 3.6.3, Section 3.4, 7 June 2008.

Member Data Documentation

matrix< eigenray_list > _eigenrays
private

List of eigenrays associated with each target.

const seq_vector* _frequencies
private

Frequencies over which loss was computed (Hz).

Linked from wavefront object so we can use it in the sum_eigenrays() calculation.

matrix< eigenray > _loss
private

Propagation loss summed over all eigenrays.

Estimates of time and angle are averages weighted by the amplitude in linear (non-dB) space. The number of surface bounces, bottom bounces, and caustics are taken from the strongest path. If there is no path to a particular target, the number of surface bounces, bottom bounces, and caustics are all set to -1.

int _num_eigenrays
private

Total number of eigenrays.

Used by write_netcdf().

const seq_vector* _source_az
private

Initial azimuthal angle at the source location (degrees, clockwise from true north).

Ray fans that wrap around 360 deg must exclude 360 itself. Linked from wavefront object so we can write it to a netCDF file.

const seq_vector* _source_de
private

Initial depression/elevation angle at the source location (degrees, positive is up).

Linked from wavefront object so we can write it to a netCDF file.

const wposition1 _source_pos
private

Location of the wavefront source in spherical earth coordinates.

Linked from wavefront object so we can write it to a netCDF file.

const wposition* _targets
private

Matrix of target positions in world coordinates.

double _time_step
private

Propagation step size (seconds).

Linked from wavefront object so we can write it to a netCDF file.