USML
wave_queue.h
1 
5 #ifndef USML_WAVEQ3D_WAVE_QUEUE_H
6 #define USML_WAVEQ3D_WAVE_QUEUE_H
7 
8 #include <usml/ocean/ocean.h>
9 #include <usml/waveq3d/reverb_model.h>
10 #include <usml/waveq3d/wave_front.h>
11 #include <usml/waveq3d/proplossListener.h>
12 #include <netcdfcpp.h>
13 
14 namespace usml {
15 namespace waveq3d {
16 
17 using namespace usml::ocean ;
18 class reverb_model ; // forward references for friend declarations
19 class reflection_model ;
20 class spreading_model ;
21 class spreading_ray ;
22 class spreading_hybrid_gaussian ;
23 class proplossListener ;
24 
27 
63 class USML_DECLSPEC wave_queue {
64 
65  friend class reflection_model ;
66  friend class spreading_ray ;
68 
69  private:
70 
76 
82 
87 
94 
101 
103  double _time_step ;
104 
106  double _time ;
107 
112 
120  matrix<double> _targets_sin_theta ;
121 
124 
130 
137  double _intensity_threshold; //In dB
138 
148  wave_front *_past, *_prev, *_curr, *_next ;
149 
150 
156  std::vector<proplossListener *> _proplossListenerVec;
157 
165 
170  bool _de_branch ;
171 
172  public:
173 
177  typedef enum { CLASSIC_RAY, HYBRID_GAUSSIAN } spreading_type ;
178 
179  //**************************************************
180  // methods
181 
204  wave_queue(
205  ocean_model& ocean,
206  const seq_vector& freq,
207  const wposition1& pos,
208  const seq_vector& de, const seq_vector& az,
209  double time_step,
210  const wposition* targets=NULL,
211  spreading_type type=HYBRID_GAUSSIAN
212  ) ;
213 
215  virtual ~wave_queue() ;
216 
222  inline const wposition1& source_pos() {
223  return _source_pos ;
224  }
225 
233  inline double source_de( unsigned de ) {
234  return (*_source_de)(de) ;
235  }
236 
244  inline double source_az( unsigned az ) {
245  return (*_source_az)(az) ;
246  }
247 
251  inline double time() {
252  return _time ;
253  }
254 
258  inline const wave_front* next() {
259  return _next ;
260  }
261 
265  inline const wave_front* curr() {
266  return _curr ;
267  }
268 
272  inline const wave_front* prev() {
273  return _prev ;
274  }
275 
279  inline const wave_front* past() {
280  return _past ;
281  }
282 
286  inline unsigned num_de() const {
287  return _source_de->size() ;
288  }
289 
293  inline unsigned num_az() const {
294  return _source_az->size() ;
295  }
296 
302  inline void setIntensityThreshold(double dThreshold) {
303 
304  // Convert to absolute value for later comparison
305  // with the positive ray.intensity value.
306  _intensity_threshold = abs(dThreshold);
307  }
312  inline double getIntensityThreshold() {
313  return _intensity_threshold;
314  }
318  void set_bottom_reverb( reverb_model* model ) ;
319 
323  void set_surface_reverb( reverb_model* model ) ;
324 
328  bool addProplossListener(proplossListener* pListener);
329 
333  bool removeProplossListener(proplossListener* pListener);
334 
339  bool notifyProplossListeners(unsigned targetRow, unsigned targetCol, eigenray pEigenray);
340 
341  private:
342 
358  void init_wavefronts() ;
359 
360  public:
361 
383  void step() ;
384 
385  private:
386 
387  //**************************************************
388  // reflections and caustics
389 
407  void detect_reflections() ;
408 
420  bool detect_reflections_surface( unsigned de, unsigned az ) ;
421 
433  bool detect_reflections_bottom( unsigned de, unsigned az ) ;
434 
442  void detect_caustics() ;
443 
444  //**************************************************
445  // eigenray estimation routines
446 
453  void detect_eigenrays() ;
454 
485  bool is_closest_ray(
486  unsigned t1, unsigned t2,
487  unsigned de, unsigned az,
488  const double &center,
489  double distance2[3][3][3] ) ;
490 
542  void build_eigenray(
543  unsigned t1, unsigned t2,
544  unsigned de, unsigned az,
545  double distance2[3][3][3] ) ;
546 
592  static void compute_offsets(
593  const double distance2[3][3][3], const c_vector<double,3>& delta,
594  c_vector<double,3>& offset, c_vector<double,3>& distance,
595  bool& unstable ) ;
596 
615  static void make_taylor_coeff(
616  const double value[3][3][3], const c_vector<double,3>& delta,
617  double& center, c_vector<double,3>& gradient, c_matrix<double,3,3>& hessian,
618  bool diagonal_only = false ) ;
619 
620  //**************************************************
621  // wavefront_netcdf routines
622 
623  private:
624 
629  NcFile* _nc_file ;
630 
632  NcVar *_nc_time, *_nc_latitude, *_nc_longitude, *_nc_altitude,
633  *_nc_surface, *_nc_bottom, *_nc_caustic, *_nc_on_edge;
634 
636  int _nc_rec ;
637 
638  public:
639 
692  void init_netcdf( const char* filename, const char* long_name=NULL ) ;
693 
699  void save_netcdf() ;
700 
704  void close_netcdf() ;
705 
706 };
707 
709 } // end of namespace waveq3d
710 } // end of namespace usml
711 
712 #endif
reflection_model * _reflection_model
Reference to the reflection loss model component.
Definition: wave_queue.h:123
ocean_model & _ocean
Reference to the environmental parameters.
Definition: wave_queue.h:75
Definition: spreading_model.h:21
const wposition * _targets
List of acoustic targets.
Definition: wave_queue.h:111
spreading_model * _spreading_model
Reference to the spreading loss model component.
Definition: wave_queue.h:129
const seq_vector * _source_de
Initial depression/elevation angle (D/E) at the source location (degrees, positive is up)...
Definition: wave_queue.h:93
World location in geodetic earth coordinates (latitude, longitude, and altitude). ...
Definition: wposition.h:39
Wavefront characteristics at a specific point in time.
Definition: wave_front.h:68
unsigned num_de() const
Number of D/E angles in the ray fan.
Definition: wave_queue.h:286
std::vector< proplossListener * > _proplossListenerVec
Vector containing the references of objects that will be used to update classes that require eigenray...
Definition: wave_queue.h:156
const wposition1 _source_pos
Location of the wavefront source in spherical earth coordinates.
Definition: wave_queue.h:86
A single acoustic path between a source and target.
Definition: eigenray.h:24
Wavefront propagator for the WaveQ3D model.
Definition: wave_queue.h:63
void setIntensityThreshold(double dThreshold)
setIntensityThreshold
Definition: wave_queue.h:302
double time()
Elapsed time for the current element in the wavefront.
Definition: wave_queue.h:251
double _time
Time for current entry in the wave_front circular queue (seconds).
Definition: wave_queue.h:106
double getIntensityThreshold()
getIntensityThreshold
Definition: wave_queue.h:312
BOOST_UBLAS_INLINE matrix_unary1_traits< E, scalar_abs< typename E::value_type > >::result_type abs(const matrix_expression< E > &e)
Magnitude of a complex matrix.
Definition: matrix_math.h:257
const wave_front * next()
Return next element in the wavefront.
Definition: wave_queue.h:258
const wave_front * prev()
Return previous element in the wavefront.
Definition: wave_queue.h:272
wave_front * _prev
Definition: wave_queue.h:148
Combines the effects of surface, bottom, and profile into a single model.
Definition: ocean_model.h:20
Definition: reflection_model.h:49
double source_de(unsigned de)
Initial depression/elevation angle at the source location.
Definition: wave_queue.h:233
bool _az_boundary
Create an Azimuthal boundary loop condition upon initialization.
Definition: wave_queue.h:164
const wposition1 & source_pos()
Location of the wavefront source in spherical earth coordinates.
Definition: wave_queue.h:222
static const double time_step
Definition: eigenray_test.cc:20
matrix< double > _targets_sin_theta
Intermediate term: sin of colatitude for targets.
Definition: wave_queue.h:120
A reverberation model listens for interface collision callbacks from a wavefront. ...
Definition: reverb_model.h:24
double _time_step
Propagation step size (seconds).
Definition: wave_queue.h:103
bool _de_branch
Treat all targets that are slightly away from directly above the source as special cases...
Definition: wave_queue.h:170
const seq_vector * _source_az
Initial azimuthal angle (AZ) at the source location (degrees, clockwise from true north)...
Definition: wave_queue.h:100
Definition: spreading_ray.h:33
int _nc_rec
Current record number in netDCF file.
Definition: wave_queue.h:636
double _intensity_threshold
The value of the intensity threshold in dB Any eigenray intensity values that are weaker than this th...
Definition: wave_queue.h:137
double source_az(unsigned az)
Initial azimuthal angle at the source location.
Definition: wave_queue.h:244
NcVar * _nc_time
The netCDF variables used to record the wavefront log.
Definition: wave_queue.h:632
World location in geodetic earth coordinates (latitude, longitude, and altitude). ...
Definition: wposition1.h:23
This class is part of a Observer/Subject pattern for the wave_queue class and allows for multiple pro...
Definition: proplossListener.h:22
const wave_front * past()
Return past element in the wavefront.
Definition: wave_queue.h:279
const wave_front * curr()
Return current element in the wavefront.
Definition: wave_queue.h:265
const seq_vector * _frequencies
Frequencies over which to compute propagation loss (Hz).
Definition: wave_queue.h:81
NcFile * _nc_file
The netCDF file used to record the wavefront log.
Definition: wave_queue.h:629
static const seq_log freq(10e3, 10e3, 1)
A read-only, monotonic sequence of values.
Definition: seq_vector.h:36
unsigned num_az() const
Number of AZ angles in the ray fan.
Definition: wave_queue.h:293
Definition: spreading_hybrid_gaussian.h:44