RSTT  3.2.0
Regional Seismic Travel Time
All Classes Namespaces Files Functions Variables Typedefs Friends Macros
SlbmInterface.h
Go to the documentation of this file.
1 //- ****************************************************************************
2 //-
3 //- Copyright 2009 National Technology & Engineering Solutions of Sandia, LLC
4 //- (NTESS). Under the terms of Contract DE-NA0003525 with NTESS, the U.S.
5 //- Government retains certain rights in this software.
6 //-
7 //- BSD Open Source License
8 //- All rights reserved.
9 //-
10 //- Redistribution and use in source and binary forms, with or without
11 //- modification, are permitted provided that the following conditions are met:
12 //-
13 //- 1. Redistributions of source code must retain the above copyright notice,
14 //- this list of conditions and the following disclaimer.
15 //-
16 //- 2. Redistributions in binary form must reproduce the above copyright
17 //- notice, this list of conditions and the following disclaimer in the
18 //- documentation and/or other materials provided with the distribution.
19 //-
20 //- 3. Neither the name of the copyright holder nor the names of its
21 //- contributors may be used to endorse or promote products derived from
22 //- this software without specific prior written permission.
23 //-
24 //- THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
25 //- AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
26 //- IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
27 //- ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
28 //- LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
29 //- CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
30 //- SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
31 //- INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
32 //- CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
33 //- ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
34 //- POSSIBILITY OF SUCH DAMAGE.
35 //-
36 //- ****************************************************************************
37 
38 #ifndef SlbmInterface_H
39 #define SlbmInterface_H
40 
41 // **** _SYSTEM INCLUDES_ ******************************************************
42 #include <map>
43 #include <cmath>
44 
45 using namespace std;
46 
47 // **** _LOCAL INCLUDES_ *******************************************************
48 
49 #include "GreatCircle.h"
50 #include "GreatCircleFactory.h"
51 #include "GridGeoTess.h"
52 #include "UncertaintyPIU.h"
53 #include "UncertaintyPDU.h"
54 #include "SLBMException.h"
55 
56 // **** _BEGIN SLBM NAMESPACE_ **************************************************
57 
58 namespace slbm {
59 
77 {
78 
79 public:
80 
89 
99  SlbmInterface(const double& earthRadius);
100 
105  virtual ~SlbmInterface();
106 
112  bool isEqual(SlbmInterface *other);
113 
114  bool operator == (SlbmInterface &other) { return isEqual(&other); };
115  bool operator != (SlbmInterface &other) { return !(*this == other); };
116 
117 
123  static bool modelsEqual(const string modelPath1, const string modelPath2);
124 
129  string getVersion() { return SlbmVersion; };
130 
142  void loadVelocityModel(const string& modelPath);
143 
188  void saveVelocityModel(const string& modelFileName, const int& format=4);
189 
191 
208  void loadVelocityModelBinary(const string& modelPath)
209  { loadVelocityModel(modelPath); };
210 
227  void specifyOutputDirectory(const string& directoryName);
228 
238  void saveVelocityModelBinary();
239 
248  void loadVelocityModelBinary(util::DataBuffer& buffer);
249 
255  void saveVelocityModelBinary(util::DataBuffer& buffer);
256 
258 
264  int getBufferSize() const;
265 
273  void setInterpolatorType(const string& interpolatorType);
274 
280  string getInterpolatorType();
281 
297  void createGreatCircle(const string& phase,
298  const double& sourceLat,
299  const double& sourceLon,
300  const double& sourceDepth,
301  const double& receiverLat,
302  const double& receiverLon,
303  const double& receiverDepth);
304 
320  void createGreatCircle(const int& phase,
321  const double& sourceLat,
322  const double& sourceLon,
323  const double& sourceDepth,
324  const double& receiverLat,
325  const double& receiverLon,
326  const double& receiverDepth);
327 
347  void clear();
348 
354  bool isValid() { return valid; };
355 
359  string getPhase() { return greatCircle->getPhaseString(); };
360 
366  void getDistance(double& distance) { distance = getDistance(); };
367 
372  double getDistance();
373 
381  void getSourceDistance(double& dist);
382  double getSourceDistance() { double x; getSourceDistance(x); return x; };
383 
391  void getReceiverDistance(double& dist);
392  double getReceiverDistance() { double x; getReceiverDistance(x); return x; };
393 
405  void getHeadwaveDistance(double& dist);
406  double getHeadwaveDistance() { double x; getHeadwaveDistance(x); return x; };
407 
420  void getHeadwaveDistanceKm(double& dist);
421  double getHeadwaveDistanceKm() { double x; getHeadwaveDistanceKm(x); return x; };
422 
431  void getTravelTime(double& travelTime);
432  double getTravelTime() { double x; getTravelTime(x); return x; };
433 
448  void getTravelTimeComponents(double& tTotal, double& tSource, double& tReceiver,
449  double& tHeadwave, double& tGradient);
450 
456  void getSlowness(double& slowness);
457  double getSlowness() { double x; getSlowness(x); return x; };
458 
464  void get_dtt_ddist(double& slowness) { getSlowness(slowness); };
465  double get_dtt_ddist() { double x; getSlowness(x); return x; };
466 
473  void get_dtt_dlat(double& dtt_dlat);
474  double get_dtt_dlat() { double x; get_dtt_dlat(x); return x; };
475 
482  void get_dtt_dlon(double& dtt_dlon);
483  double get_dtt_dlon() { double x; get_dtt_dlon(x); return x; };
484 
491  void get_dtt_ddepth(double& dtt_ddepth);
492  double get_dtt_ddepth() { double x; get_dtt_ddepth(x); return x; };
493 
500  //void get_dtt_dlat_fast(double& dtt_dlat);
501 
508  //void get_dtt_dlon_fast(double& dtt_dlon);
509 
510 // //! \brief Retrieve the derivative of horizontal slowness wrt to source-receiver distance,
511 // //! in seconds/radian^2.
512 // //!
513 // //! Retrieve the derivative of horizontal slowness wrt to source-receiver distance,
514 // //! in seconds/radian^2.
515 // void get_dsh_ddist(double& dsh_ddist);
516 // double get_dsh_ddist() { double x; get_dsh_ddist(x); return x; };
517 //
518 // //! \brief Retrieve the derivative of horizontal slowness wrt to source latitude,
519 // //! in seconds/radian^2.
520 // //!
521 // //! Retrieve the derivative of horizontal slowness wrt to source latitude,
522 // //! in seconds/radian^2.
523 // //! @param dsh_dlat the derivative of horizontal slowness wrt to source latitude.
524 // void get_dsh_dlat(double& dsh_dlat);
525 //
526 // //! \brief Retrieve the derivative of horizontal slowness wrt to source longitude,
527 // //! in seconds/radian^2.
528 // //!
529 // //! Retrieve the derivative of horizontal slowness wrt to source longitude,
530 // //! in seconds/radian^2.
531 // //! @param dsh_dlon the derivative of horizontal slowness wrt to source longitude.
532 // void get_dsh_dlon(double& dsh_dlon);
533 //
534 // //! \brief Retrieve the derivative of horizontal slowness wrt to source depth,
535 // //! in seconds/radian-km.
536 // //!
537 // //! Retrieve the derivative of horizontal slowness wrt to source depth,
538 // //! in seconds/radian-km.
539 // //! @param dsh_ddepth the derivative of horizontal slowness wrt to source depth.
540 // void get_dsh_ddepth(double& dsh_ddepth);
541 // double get_dsh_ddepth() { double x; get_dsh_ddepth(x); return x; }
542 
580  void getWeights(int nodeId[], double weight[], int& nWeights);
581 
620  void getActiveNodeWeights(int nodeId[], double weight[], int& nWeights);
621 
655  void getWeights(vector<int>& nodeId, vector<double>& weight);
656 
689  void getActiveNodeWeights(vector<int>& nodeId, vector<double>& weight);
690 
704  void getWeightsSource(int nodeids[], double weights[], int& nWeights);
705  void getWeightsSource(vector<int> &nodeids, vector<double> &weights);
706 
720  void getActiveNodeWeightsSource(int nodeids[], double weights[], int& nWeights);
721  void getActiveNodeWeightsSource(vector<int> &nodeids, vector<double> &weights);
722 
733  void getWeightsReceiver(int nodeids[], double weights[], int& nWeights);
734  void getWeightsReceiver(vector<int> &nodeids, vector<double> &weights);
735 
746  void getActiveNodeWeightsReceiver(int nodeids[], double weights[], int& nWeights);
747  void getActiveNodeWeightsReceiver(vector<int> &nodeids, vector<double> &weights);
748 
770  string toString(const int& verbosity);
771 
775  void getNGridNodes(int& n);
776 
780  int getNGridNodes();
781 
791  void getNHeadWavePoints(int& nHeadWavePoints);
792  int getNHeadWavePoints() { int x; getNHeadWavePoints(x); return x; };
793 
814  const int& nodeId,
815  double& latitude,
816  double& longitude,
817  double depth[NLAYERS],
818  double pvelocity[NLAYERS],
819  double svelocity[NLAYERS],
820  double gradient[2]);
821 
842  const int& nodeId,
843  double& latitude,
844  double& longitude,
845  double depth[NLAYERS],
846  double pvelocity[NLAYERS],
847  double svelocity[NLAYERS],
848  double gradient[2]);
849 
866  const int& nodeId,
867  double depths[NLAYERS],
868  double pvelocity[NLAYERS],
869  double svelocity[NLAYERS],
870  double gradient[2]);
871 
888  const int& nodeId,
889  double depths[NLAYERS],
890  double pvelocity[NLAYERS],
891  double svelocity[NLAYERS],
892  double gradient[2]);
893 
930  void getGreatCircleData(
931  string& phase,
932  double& actual_path_increment,
933  double sourceDepth[NLAYERS],
934  double sourceVelocity[NLAYERS],
935  double receiverDepth[NLAYERS],
936  double receiverVelocity[NLAYERS],
937  int& npoints,
938  double headWaveVelocity[],
939  double gradient[]
940  );
941  void getGreatCircleData(
942  string &phase,
943  double &actual_path_increment,
944  vector<double> &sourceDepth,
945  vector<double> &sourceVelocity,
946  vector<double> &receiverDepth,
947  vector<double> &receiverVelocity,
948  int &npoints,
949  vector<double> &headWaveVelocity,
950  vector<double> &gradient
951  );
952 
966  void getGreatCircleLocations(double lat[], double lon[], double depth[], int& npoints);
967  void getGreatCircleLocations(vector<double> &lat, vector<double> &lon, vector<double> &depth);
968 
996  void getGreatCircleNodeInfo(
997  int** neighbors,
998  double** coefficients,
999  const int& maxpoints,
1000  const int& maxnodes,
1001  int& npoints,
1002  int* nnodes
1003  );
1004  void getGreatCircleNodeInfo(vector<vector<int> >& neighbors, vector<vector<double> >& coefficients);
1005 
1029  void getInterpolatedPoint(
1030  const double& lat,
1031  const double& lon,
1032  int* nodeIds,
1033  double* coefficients,
1034  int& nWeights,
1035  double depth[NLAYERS],
1036  double pvelocity[NLAYERS],
1037  double svelocity[NLAYERS],
1038  double& pgradient,
1039  double& sgradient
1040  );
1041 
1064  void getInterpolatedPoint(
1065  const double &lat,
1066  const double &lon,
1067  vector<int> &nodeId,
1068  vector<double> &coefficients,
1069  vector<double> &depth,
1070  vector<double> &pvelocity,
1071  vector<double> &svelocity,
1072  double &pgradient,
1073  double &sgradient
1074  );
1075 
1101  void getInterpolatedTransect(
1102  double lat[],
1103  double lon[],
1104  const int& nLatLon,
1105  int** neighbors,
1106  double** coefficients,
1107  int* nNeighbors,
1108  double depth[][NLAYERS],
1109  double pvelocity[][NLAYERS],
1110  double svelocity[][NLAYERS],
1111  double pgradient[NLAYERS],
1112  double sgradient[NLAYERS],
1113  int& nInvalid
1114  );
1115 
1135  void getInterpolatedTransect(
1136  vector<double> lat,
1137  vector<double> lon,
1138  vector<vector<int> > &nodeId,
1139  vector<vector<double> > &coefficients,
1140  vector<vector<double> > &depth,
1141  vector<vector<double> > &pvelocity,
1142  vector<vector<double> > &svelocity,
1143  vector<double> &pgradient,
1144  vector<double> &sgradient,
1145  int& nInvalid
1146  );
1147 
1160  void initializeActiveNodes(const double& latmin, const double& lonmin,
1161  const double& latmax, const double& lonmax);
1162 
1179  void initializeActiveNodes(const string& polygonFileName);
1180 
1193  void initializeActiveNodes(GeoTessPolygon* polygon);
1194 
1216  void initializeActiveNodes(double* lat, double* lon, const int& npoints, const bool& inDegrees=true);
1217 
1237  void initializeActiveNodes(const vector<double> lat, const vector<double> lon, const bool& inDegrees=true);
1238 
1254  void initializeActiveNodes(vector<double*>& unitVectors);
1255 
1259  int getNActiveNodes();
1260 
1263  void clearActiveNodes();
1264 
1270  int getGridNodeId(int activeNodeId);
1271 
1277  int getActiveNodeId(int gridNodeId);
1278 
1287  void getNodeHitCount(const int& nodeId, int& hitCount);
1288  int getNodeHitCount(const int &nodeId) { int x; getNodeHitCount(nodeId, x); return x; };
1289 
1290 
1296  void clearNodeHitCount();
1297 
1305  void getNodeNeighbors(const int& nid, int neighbors[], int& nNeighbors);
1306 
1314  void getActiveNodeNeighbors(const int& nid, int neighbors[], int& nNeighbors);
1315 
1320  void getNodeNeighbors(const int& nid, vector<int>& neighbors);
1321 
1326  void getActiveNodeNeighbors(const int& nid, vector<int>& neighbors);
1327 
1335  void getNodeNeighborInfo(const int& nid, int neighbors[], double distance[],
1336  double azimuth[], int& nNeighbors);
1337  void getNodeNeighborInfo(const int nid,
1338  vector<int> &neighbors, vector<double> &distance, vector<double> &azimuth);
1339 
1347  void getActiveNodeNeighborInfo(const int& nid, int neighbors[], double distance[],
1348  double azimuth[], int& nNeighbors);
1349 
1357  void getActiveNodeNeighborInfo(const int& nid,
1358  vector<int>& neighbors, vector<double>& distance, vector<double>& azimuth);
1359 
1363  void getNodeSeparation(const int& node1, const int& node2, double& distance);
1364 
1368  void getNodeAzimuth(const int& node1, const int& node2, double& azimuth);
1369 
1378  GreatCircle* getGreatCircleObject();
1379 
1384  Grid* getGridObject() { return grid; };
1385 
1390  GeoTessModelSLBM* getModelObject() { return grid->getModel(); }
1391 
1400  void getTravelTimeUncertainty( const int& phase, const double& distance, double& uncert );
1401 
1410  void getTravelTimeUncertainty(double& travelTimeUncertainty, bool calcRandomError = false);
1411 
1418  void getTravelTimeUncertainty1D(double& travelTimeUncertainty);
1419 
1428  void getSlownessUncertainty( const int& phase, const double& distance, double& uncert );
1429 
1436  void getSlownessUncertainty(double& slownessUncertainty);
1437 
1438  string getUncertaintyTable(const int& phase, const int& attribute);
1439 
1440  string getUncertaintyFileFormat(const int& phase, const int& attribute);
1441 
1463  void getZhaoParameters(double& Vm, double& Gm, double& H, double& C, double& Cm, int& udSign);
1464 
1481  void getPgLgComponents(double& tTotal, double& tTaup, double& tHeadwave,
1482  double& pTaup, double& pHeadwave, double& trTaup, double& trHeadwave);
1483 
1492  void static setCHMax(const double& chMax);
1493 
1502  void static getCHMax(double& chMax);
1503  double static getCHMax() { double x; getCHMax(x); return x; }
1504 
1514  void getAverageMantleVelocity(const int& type, double& velocity);
1515  double getAverageMantleVelocity(const int &type);
1516 
1528  void setAverageMantleVelocity(const int& type, const double& velocity);
1529 
1533  void getTessId(string& tessId);
1534  string getTessId() { string s; getTessId(s); return s; };
1535 
1541  void getFractionActive(double& fractionActive);
1542  double getFractionActive() { double x; getFractionActive(x); return x; }
1543 
1551  void static setMaxDistance(const double& maxDistance);
1552 
1558  void static getMaxDistance(double& maxDistance);
1559  double static getMaxDistance() { double x; getMaxDistance(x); return x; }
1560 
1568  void static setMaxDepth(const double& maxDepth);
1569 
1573  void static getMaxDepth(double& maxDepth);
1574  double static getMaxDepth() { double x; getMaxDepth(x); return x; }
1575 
1582  string getClassCount();
1583 
1588  const string& getModelPath() const;
1589 
1605  void getDistAz(const double& aLat, const double& aLon,
1606  const double& bLat, const double& bLon,
1607  double& distance, double& azimuth, const double& naValue);
1608 
1621  void movePoint(const double& aLat, const double& aLon,
1622  const double& distance, const double& azimuth,
1623  double& bLat, double& bLon);
1624 
1633  void getPiercePointSource(double& lat, double& lon, double& depth);
1634 
1643  void getPiercePointReceiver(double& lat, double& lon, double& depth);
1644 
1661  void getGreatCirclePoints(
1662  const double& aLat,
1663  const double& aLon,
1664  const double& bLat,
1665  const double& bLon,
1666  const int& npoints,
1667  double latitude[],
1668  double longitude[]);
1669  void getGreatCirclePoints(
1670  const double& aLat,
1671  const double& aLon,
1672  const double& bLat,
1673  const double& bLon,
1674  const int& npoints,
1675  vector<double> &latitude,
1676  vector<double> &longitude);
1677 
1692  void getGreatCirclePointsOnCenters(
1693  const double& aLat,
1694  const double& aLon,
1695  const double& bLat,
1696  const double& bLon,
1697  const int& npoints,
1698  double latitude[],
1699  double longitude[]);
1700  void getGreatCirclePointsOnCenters(
1701  const double& aLat,
1702  const double& aLon,
1703  const double& bLat,
1704  const double& bLon,
1705  const int& npoints,
1706  vector<double> &latitude,
1707  vector<double> &longitude);
1708 
1714  void setDelDistance(const double& del_distance);
1715 
1721  void getDelDistance(double& del_distance);
1722  double getDelDistance() { double x; getDelDistance(x); return x; };
1723 
1729  void setDelDepth(const double& del_depth);
1730 
1735  void getDelDepth(double& del_depth);
1736  double getDelDepth() { double x; getDelDepth(x); return x; };
1737 
1741  void getRayParameter(double& ray_parameter);
1742  double getRayParameter() { double x; getRayParameter(x); return x; };
1743 
1747  void getTurningRadius(double& turning_radius);
1748  double getTurningRadius() { double x; getTurningRadius(x); return x; };
1749 
1763  void setPathIncrement(const double& pathIncrement);
1764 
1777  void getPathIncrement(double& pathIncrement);
1778 
1792 
1793  string getModelString() { return grid->toString(); }
1794 
1797  double getSrcLat() { return srcLat; }
1798 
1801  double getSrcLon() { return srcLon; }
1802 
1805  double getSrcDep() { return srcDep; }
1806 
1809  double getRcvLat() { return rcvLat; }
1810 
1813  double getRcvLon() { return rcvLon; }
1814 
1817  double getRcvDep() { return rcvDep; }
1818 
1819 
1820 protected:
1821 
1826 
1831 
1835  static double CH_MAX;
1836 
1842  bool valid;
1843 
1847 
1848 private:
1849 
1850  // copies of the source and receiver locations that were
1851  // passed to createGreatCircle the last time it was called.
1852  // Used to compute distance in instances when createGreateCircle
1853  // fails. Units are radians.
1854  double srcLat, srcLon, srcDep, rcvLat, rcvLon, rcvDep;
1855 };
1856 
1857 //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1858 //
1859 // INLINE FUNCTIONS
1860 //
1861 //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1862 
1863 inline void SlbmInterface::createGreatCircle(
1864  const string& p,
1865  const double& sourceLat,
1866  const double& sourceLon,
1867  const double& sourceDepth,
1868  const double& receiverLat,
1869  const double& receiverLon,
1870  const double& receiverDepth)
1871 {
1872  int iphase = (p=="Pn" ? Pn : (p=="Sn" ? Sn : (p=="Pg" ? Pg : (p=="Lg" ? Lg : -1))));
1873  if (iphase == -1)
1874  {
1875  ostringstream os;
1876  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1877  os << endl << "ERROR in SlbmInterface::createGreatCircle" << endl
1878  << p << " is not a recognized phase. Must be one of Pn, Sn, Pg, Lg" << endl
1879  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1880  throw SLBMException(os.str(),112);
1881  }
1882 
1883  createGreatCircle(
1884  iphase,
1885  sourceLat,
1886  sourceLon,
1887  sourceDepth,
1888  receiverLat,
1889  receiverLon,
1890  receiverDepth
1891  );
1892 }
1893 
1894 inline GreatCircle* SlbmInterface::getGreatCircleObject()
1895 {
1896  return greatCircle;
1897 }
1898 
1899 inline void SlbmInterface::get_dtt_dlat(double& value)
1900 {
1901  if (!isValid())
1902  {
1903  value = NA_VALUE;
1904  ostringstream os;
1905  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1906  os << endl << "ERROR in SlbmInterface::get_dtt_dlat" << endl
1907  << "GreatCircle is invalid." << endl
1908  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1909  throw SLBMException(os.str(),113);
1910  }
1911  greatCircle->get_dtt_dlat(value);
1912 }
1913 
1914 inline void SlbmInterface::get_dtt_dlon(double& value)
1915 {
1916  if (!isValid())
1917  {
1918  value = NA_VALUE;
1919  ostringstream os;
1920  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1921  os << endl << "ERROR in SlbmInterface::get_dtt_dlon" << endl
1922  << "GreatCircle is invalid." << endl
1923  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1924  throw SLBMException(os.str(),113);
1925  }
1926  greatCircle->get_dtt_dlon(value);
1927 }
1928 
1929 inline void SlbmInterface::get_dtt_ddepth(double& value)
1930 {
1931  if (!isValid())
1932  {
1933  value = NA_VALUE;
1934  ostringstream os;
1935  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1936  os << endl << "ERROR in SlbmInterface::get_dtt_ddepth" << endl
1937  << "GreatCircle is invalid." << endl
1938  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1939  throw SLBMException(os.str(),113);
1940  }
1941  greatCircle->get_dtt_ddepth(value);
1942 }
1943 
1944 //inline void SlbmInterface::get_dtt_dlat_fast(double& value)
1945 //{
1946 // if (!isValid())
1947 // {
1948 // value = NA_VALUE;
1949 // ostringstream os;
1950 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1951 // os << endl << "ERROR in SlbmInterface::get_dtt_dlat" << endl
1952 // << "GreatCircle is invalid." << endl
1953 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1954 // throw SLBMException(os.str(),113);
1955 // }
1956 // greatCircle->get_dtt_dlat_fast(value);
1957 //}
1958 //
1959 //inline void SlbmInterface::get_dtt_dlon_fast(double& value)
1960 //{
1961 // if (!isValid())
1962 // {
1963 // value = NA_VALUE;
1964 // ostringstream os;
1965 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1966 // os << endl << "ERROR in SlbmInterface::get_dtt_dlon" << endl
1967 // << "GreatCircle is invalid." << endl
1968 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1969 // throw SLBMException(os.str(),113);
1970 // }
1971 // greatCircle->get_dtt_dlon_fast(value);
1972 //}
1973 //
1974 //inline void SlbmInterface::get_dsh_ddist(double& value)
1975 //{
1976 // if (!isValid())
1977 // {
1978 // value = NA_VALUE;
1979 // ostringstream os;
1980 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1981 // os << endl << "ERROR in SlbmInterface::get_dsh_ddist" << endl
1982 // << "GreatCircle is invalid." << endl
1983 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1984 // throw SLBMException(os.str(),113);
1985 // }
1986 // greatCircle->get_dsh_ddist(value);
1987 //}
1988 //
1989 //inline void SlbmInterface::get_dsh_dlat(double& value)
1990 //{
1991 // if (!isValid())
1992 // {
1993 // value = NA_VALUE;
1994 // ostringstream os;
1995 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1996 // os << endl << "ERROR in SlbmInterface::get_dsh_dlat" << endl
1997 // << "GreatCircle is invalid." << endl
1998 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1999 // throw SLBMException(os.str(),113);
2000 // }
2001 // greatCircle->get_dsh_dlat(value);
2002 //}
2003 //
2004 //inline void SlbmInterface::get_dsh_dlon(double& value)
2005 //{
2006 // if (!isValid())
2007 // {
2008 // value = NA_VALUE;
2009 // ostringstream os;
2010 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2011 // os << endl << "ERROR in SlbmInterface::get_dsh_dlon" << endl
2012 // << "GreatCircle is invalid." << endl
2013 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2014 // throw SLBMException(os.str(),113);
2015 // }
2016 // greatCircle->get_dsh_dlon(value);
2017 //}
2018 //
2019 //inline void SlbmInterface::get_dsh_ddepth(double& value)
2020 //{
2021 // if (!isValid())
2022 // {
2023 // value = NA_VALUE;
2024 // ostringstream os;
2025 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2026 // os << endl << "ERROR in SlbmInterface::get_dsh_ddepth" << endl
2027 // << "GreatCircle is invalid." << endl
2028 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2029 // throw SLBMException(os.str(),113);
2030 // }
2031 // greatCircle->get_dsh_ddepth(value);
2032 //}
2033 
2034 inline void SlbmInterface::getSlowness(double& value)
2035 {
2036  if (!isValid())
2037  {
2038  value = NA_VALUE;
2039  ostringstream os;
2040  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2041  os << endl << "ERROR in SlbmInterface::getSlowness" << endl
2042  << "GreatCircle is invalid." << endl
2043  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2044  throw SLBMException(os.str(),113);
2045  }
2046  greatCircle->get_dtt_ddist(value);
2047 }
2048 
2049 inline double SlbmInterface::getDistance()
2050 {
2051  if (greatCircle == NULL)
2052  {
2053  Location s(srcLat, srcLon, 0.);
2054  Location r(rcvLat, rcvLon, 0.);
2055  return s.distance(r);
2056  }
2057  else
2058  return greatCircle->getDistance();
2059 }
2060 
2061 inline void SlbmInterface::getSourceDistance(double& distance)
2062 {
2063  if (!isValid())
2064  {
2065  distance = NA_VALUE;
2066  ostringstream os;
2067  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2068  os << endl << "ERROR in SlbmInterface::getSourceDistance" << endl
2069  << "GreatCircle is invalid." << endl
2070  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2071  throw SLBMException(os.str(),113);
2072  }
2073  distance = greatCircle->getSourceDistance();
2074 }
2075 
2076 inline void SlbmInterface::getReceiverDistance(double& distance)
2077 {
2078  if (!isValid())
2079  {
2080  distance = NA_VALUE;
2081  ostringstream os;
2082  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2083  os << endl << "ERROR in SlbmInterface::getReceiverDistance" << endl
2084  << "GreatCircle is invalid." << endl
2085  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2086  throw SLBMException(os.str(),113);
2087  }
2088  distance = greatCircle->getReceiverDistance();
2089 }
2090 
2091 inline void SlbmInterface::getHeadwaveDistance(double& distance)
2092 {
2093  if (!isValid())
2094  {
2095  distance = NA_VALUE;
2096  ostringstream os;
2097  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2098  os << endl << "ERROR in SlbmInterface::getHeadwaveDistance" << endl
2099  << "GreatCircle is invalid." << endl
2100  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2101  throw SLBMException(os.str(),113);
2102  }
2103  distance = greatCircle->getHeadwaveDistance();
2104 }
2105 
2106 inline void SlbmInterface::getHeadwaveDistanceKm(double& distance)
2107 {
2108  if (!isValid())
2109  {
2110  distance = NA_VALUE;
2111  ostringstream os;
2112  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2113  os << endl << "ERROR in SlbmInterface::getHeadwaveDistanceKm" << endl
2114  << "GreatCircle is invalid." << endl
2115  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2116  throw SLBMException(os.str(),113);
2117  }
2118  distance = greatCircle->getHeadwaveDistanceKm();
2119 }
2120 
2121 inline void SlbmInterface::getTravelTime(double& tTotal)
2122 {
2123  if (!isValid())
2124  {
2125  tTotal = NA_VALUE;
2126  ostringstream os;
2127  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2128  os << endl << "ERROR in SlbmInterface::getTravelTime" << endl
2129  << "GreatCircle is invalid." << endl
2130  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2131  throw SLBMException(os.str(),113);
2132  }
2133  tTotal=greatCircle->getTravelTime();
2134 }
2135 
2136 inline void SlbmInterface::getTravelTimeComponents(
2137  double& tTotal, double& tSource, double& tReceiver,
2138  double& tHeadwave, double& tGradient)
2139 {
2140  if (!isValid())
2141  {
2142  tTotal = NA_VALUE;
2143  tSource = NA_VALUE;
2144  tReceiver = NA_VALUE;
2145  tHeadwave = NA_VALUE;
2146  tGradient = NA_VALUE;
2147  ostringstream os;
2148  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2149  os << endl << "ERROR in SlbmInterface::getTravelTimeComponents" << endl
2150  << "GreatCircle is invalid." << endl
2151  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2152  throw SLBMException(os.str(),113);
2153  }
2154  greatCircle->getTravelTime(tTotal, tSource, tReceiver,
2155  tHeadwave, tGradient);
2156 }
2157 
2158 inline void SlbmInterface::getWeights(int nodeId[], double weight[], int& nWeights)
2159 {
2160  if (!isValid())
2161  {
2162  nWeights=-1;
2163  ostringstream os;
2164  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2165  os << endl << "ERROR in SlbmInterface::getWeights" << endl
2166  << "GreatCircle is invalid." << endl
2167  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2168  throw SLBMException(os.str(),113);
2169  }
2170  greatCircle->getWeights(nodeId, weight, nWeights);
2171 }
2172 
2173 inline void SlbmInterface::getActiveNodeWeights(int nodeId[], double weight[], int& nWeights)
2174 {
2175  getWeights(nodeId, weight, nWeights);
2176  for (int i=0; i<nWeights; ++i)
2177  nodeId[i] = grid->getActiveNodeId(nodeId[i]);
2178 }
2179 
2180 inline void SlbmInterface::getWeights(vector<int>& nodeId, vector<double>& weight)
2181 {
2182  if (!isValid())
2183  {
2184  nodeId.clear();
2185  weight.clear();
2186  ostringstream os;
2187  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2188  os << endl << "ERROR in SlbmInterface::getWeights" << endl
2189  << "GreatCircle is invalid." << endl
2190  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2191  throw SLBMException(os.str(),113);
2192  }
2193  greatCircle->getWeights(nodeId, weight);
2194 }
2195 
2196 inline void SlbmInterface::getActiveNodeWeights(vector<int>& nodeId, vector<double>& weight)
2197 {
2198  getWeights(nodeId, weight);
2199  for (int i=0; i<(int)nodeId.size(); ++i)
2200  nodeId[i] = grid->getActiveNodeId(nodeId[i]);
2201 }
2202 
2203 inline void SlbmInterface::getWeightsSource(int nodeids[], double weights[], int& nWeights)
2204 {
2205  if (!isValid())
2206  {
2207  ostringstream os;
2208  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2209  os << endl << "ERROR in SlbmInterface::getWeightsSource" << endl
2210  << "GreatCircle is invalid." << endl
2211  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2212  throw SLBMException(os.str(),113);
2213  }
2214  greatCircle->getSourceProfile()->getWeights(nodeids, weights, nWeights);
2215 }
2216 
2217 inline void SlbmInterface::getWeightsSource(vector<int> &nodeids, vector<double> &weights)
2218 {
2219  // get the number of weights
2220  int nWeights = greatCircle->getSourceProfile()->getNCoefficients();
2221 
2222  // resize arrays to proper size
2223  nodeids.resize(nWeights);
2224  weights.resize(nWeights);
2225 
2226  // call getWeightsSource (v.data() sends a pointer similar to a c++ array)
2227  getWeightsSource(nodeids.data(), weights.data(), nWeights);
2228 }
2229 
2230 inline void SlbmInterface::getActiveNodeWeightsSource(int nodeids[], double weights[], int& nWeights)
2231 {
2232  getWeightsSource(nodeids, weights, nWeights);
2233  for (int i=0; i<nWeights; ++i)
2234  nodeids[i] = grid->getActiveNodeId(nodeids[i]);
2235 }
2236 
2237 inline void SlbmInterface::getActiveNodeWeightsSource(vector<int> &nodeids, vector<double> &weights)
2238 {
2239  // get the number of weights
2240  int nWeights = greatCircle->getSourceProfile()->getNCoefficients();
2241 
2242  // resize arrays to proper size
2243  nodeids.resize(nWeights);
2244  weights.resize(nWeights);
2245 
2246  // call getActiveNodeWeightsSource (v.data() sends a pointer similar to a c++ array)
2247  getActiveNodeWeightsSource(nodeids.data(), weights.data(), nWeights);
2248 }
2249 
2250 inline void SlbmInterface::getWeightsReceiver(int nodeids[], double weights[], int& nWeights)
2251 {
2252  if (!isValid())
2253  {
2254  ostringstream os;
2255  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2256  os << endl << "ERROR in SlbmInterface::getWeightsReceiver" << endl
2257  << "GreatCircle is invalid." << endl
2258  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2259  throw SLBMException(os.str(),113);
2260  }
2261  greatCircle->getReceiverProfile()->getWeights(nodeids, weights, nWeights);
2262 }
2263 
2264 inline void SlbmInterface::getWeightsReceiver(vector<int> &nodeids, vector<double> &weights)
2265 {
2266  // get the number of weights
2267  int nWeights = greatCircle->getReceiverProfile()->getNCoefficients();
2268 
2269  // resize arrays to proper size
2270  nodeids.resize(nWeights);
2271  weights.resize(nWeights);
2272 
2273  // call getActiveNodeWeightsSource (v.data() sends a pointer similar to a c++ array)
2274  getWeightsReceiver(nodeids.data(), weights.data(), nWeights);
2275 }
2276 
2277 inline void SlbmInterface::getActiveNodeWeightsReceiver(int nodeids[], double weights[], int& nWeights)
2278 {
2279  getWeightsReceiver(nodeids, weights, nWeights);
2280  for (int i=0; i<nWeights; ++i)
2281  nodeids[i] = grid->getActiveNodeId(nodeids[i]);
2282 }
2283 
2284 inline void SlbmInterface::getActiveNodeWeightsReceiver(vector<int> &nodeids, vector<double> &weights)
2285 {
2286  // get the number of weights
2287  int nWeights = greatCircle->getReceiverProfile()->getNCoefficients();
2288 
2289  // resize arrays to proper size
2290  nodeids.resize(nWeights);
2291  weights.resize(nWeights);
2292 
2293  // call getActiveNodeWeightsSource (v.data() sends a pointer similar to a c++ array)
2294  getActiveNodeWeightsReceiver(nodeids.data(), weights.data(), nWeights);
2295 }
2296 
2297 inline string SlbmInterface::toString(const int& verbosity)
2298 {
2299  if (!isValid())
2300  {
2301  ostringstream os;
2302  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2303  os << endl << "ERROR in SlbmInterface::toString" << endl
2304  << "GreatCircle is invalid." << endl
2305  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2306  throw SLBMException(os.str(),113);
2307  }
2308 
2309  ostringstream os;
2310  if (verbosity > 0)
2311  {
2312  os << endl
2313  << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << endl
2314  << "Great Circle " << endl << endl
2315  << greatCircle->toString(verbosity);
2316 
2317  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint)
2318  << setprecision(4);
2319  double slowness;
2320  getSlowness(slowness);
2321  os << endl << "Horizontal slowness = "
2322  << slowness << " sec/radian"
2323  << endl << endl;
2324  }
2325 
2326  return os.str();
2327 }
2328 
2329 inline void SlbmInterface::getGreatCircleData(
2330  string& phase,
2331  double& actual_path_increment,
2332  double sourceDepth[NLAYERS],
2333  double sourceVelocity[NLAYERS],
2334  double receiverDepth[NLAYERS],
2335  double receiverVelocity[NLAYERS],
2336  int& npoints,
2337  double headWaveVelocity[],
2338  double gradient[]
2339  )
2340 {
2341  if (!isValid())
2342  {
2343  phase = "";
2344  actual_path_increment = NA_VALUE;
2345  ostringstream os;
2346  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2347  os << endl << "ERROR in SlbmInterface::getGreatCircleData" << endl
2348  << "GreatCircle is invalid." << endl
2349  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2350  throw SLBMException(os.str(),113);
2351  }
2352 
2353  int p;
2354  greatCircle->getData(p, actual_path_increment,
2355  sourceDepth, sourceVelocity, receiverDepth, receiverVelocity,
2356  npoints, headWaveVelocity, gradient);
2357  phase = (p==Pn ? "Pn" : (p==Sn ? "Sn" : (p==Pg ? "Pg" : (p==Lg ? "Lg"
2358  : "unknown phase"))));
2359 }
2360 
2361 inline void SlbmInterface::getGreatCircleData(
2362  string &phase,
2363  double &actual_path_increment,
2364  vector<double> &sourceDepth,
2365  vector<double> &sourceVelocity,
2366  vector<double> &receiverDepth,
2367  vector<double> &receiverVelocity,
2368  int &npoints,
2369  vector<double> &headWaveVelocity,
2370  vector<double> &gradient
2371  )
2372 {
2373  if (!isValid())
2374  {
2375  phase = "";
2376  actual_path_increment = NA_VALUE;
2377  ostringstream os;
2378  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2379  os << endl << "ERROR in SlbmInterface::getGreatCircleData" << endl
2380  << "GreatCircle is invalid." << endl
2381  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2382  throw SLBMException(os.str(),113);
2383  }
2384 
2385  // figure out how many points there are
2386  getNHeadWavePoints(npoints);
2387 
2388  // make some room in memory
2389  sourceDepth.resize(NLAYERS);
2390  sourceVelocity.resize(NLAYERS);
2391  receiverDepth.resize(NLAYERS);
2392  receiverVelocity.resize(NLAYERS);
2393  headWaveVelocity.resize(npoints);
2394  gradient.resize(npoints);
2395 
2396  int p;
2397  greatCircle->getData(p, actual_path_increment,
2398  sourceDepth.data(), sourceVelocity.data(), receiverDepth.data(),
2399  receiverVelocity.data(), npoints, headWaveVelocity.data(), gradient.data());
2400  phase = (p==Pn ? "Pn" : (p==Sn ? "Sn" : (p==Pg ? "Pg" : (p==Lg ? "Lg"
2401  : "unknown phase"))));
2402 }
2403 
2404 inline void SlbmInterface::getGreatCircleNodeInfo(
2405  int** neighbors,
2406  double** coefficients,
2407  const int& maxpoints,
2408  const int& maxnodes,
2409  int& npoints,
2410  int* nnodes
2411  )
2412 {
2413  if (!isValid())
2414  {
2415  ostringstream os;
2416  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2417  os << endl << "ERROR in SlbmInterface::getGreatCircleNodeInfo" << endl
2418  << "GreatCircle is invalid." << endl
2419  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2420  throw SLBMException(os.str(),113);
2421  }
2422 
2423  greatCircle->getNodeInfo(neighbors, coefficients, maxpoints, maxnodes, npoints, nnodes);
2424 }
2425 inline void SlbmInterface::getGreatCircleNodeInfo(
2426  vector<vector<int> >& neighbors, vector<vector<double> >& coefficients)
2427  {
2428  if (!isValid())
2429  {
2430  ostringstream os;
2431  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2432  os << endl << "ERROR in SlbmInterface::getGreatCircleNodeInfo" << endl
2433  << "GreatCircle is invalid." << endl
2434  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2435  throw SLBMException(os.str(),113);
2436  }
2437 
2438  greatCircle->getNodeInfo(neighbors, coefficients);
2439  }
2440 
2441 inline void SlbmInterface::getGreatCircleLocations(
2442  double lat[],
2443  double lon[],
2444  double depth[],
2445  int& npoints
2446  )
2447 {
2448  if (!isValid())
2449  {
2450  ostringstream os;
2451  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2452  os << endl << "ERROR in SlbmInterface::getGreatCircleData" << endl
2453  << "GreatCircle is invalid." << endl
2454  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2455  throw SLBMException(os.str(),113);
2456  }
2457 
2458  npoints = greatCircle->getNProfiles();
2459  Location loc;
2460  for (int i=0; i<greatCircle->getNProfiles(); i++)
2461  {
2462  // get the location of this layer profile
2463  greatCircle->getLayerProfileLocation(i, loc);
2464  lat[i] = loc.getLat();
2465  lon[i] = loc.getLon();
2466  depth[i] = loc.getDepth();
2467  }
2468 }
2469 
2470 inline void SlbmInterface::getGreatCircleLocations(
2471  vector<double> &lat,
2472  vector<double> &lon,
2473  vector<double> &depth
2474  )
2475 {
2476  if (!isValid())
2477  {
2478  ostringstream os;
2479  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2480  os << endl << "ERROR in SlbmInterface::getGreatCircleData" << endl
2481  << "GreatCircle is invalid." << endl
2482  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2483  throw SLBMException(os.str(),113);
2484  }
2485  int npoints = getNHeadWavePoints();
2486  lat.resize(npoints);
2487  lon.resize(npoints);
2488  depth.resize(npoints);
2489  getGreatCircleLocations(lat.data(), lon.data(), depth.data(), npoints);
2490 }
2491 
2492 
2493 inline void SlbmInterface::getInterpolatedPoint(
2494  const double& lat,
2495  const double& lon,
2496  int* nodeIds,
2497  double* coefficients,
2498  int& nWeights,
2499  double depth[NLAYERS],
2500  double pvelocity[NLAYERS],
2501  double svelocity[NLAYERS],
2502  double& pgradient,
2503  double& sgradient
2504  )
2505 {
2506  Location location(lat, lon, 0.);
2507  QueryProfile* profile = grid->getQueryProfile(location);
2508  profile->getData(nodeIds, coefficients, nWeights, depth, pvelocity, svelocity, pgradient, sgradient);
2509  delete profile;
2510 }
2511 
2512 inline void SlbmInterface::getInterpolatedPoint(
2513  const double &lat,
2514  const double &lon,
2515  vector<int> &nodeIds,
2516  vector<double> &coefficients,
2517  vector<double> &depth,
2518  vector<double> &pvelocity,
2519  vector<double> &svelocity,
2520  double &pgradient,
2521  double &sgradient
2522  )
2523 {
2524  // create a new location at the lat/lon
2525  Location location(lat, lon, 0.);
2526 
2527  // interpolate a profile at that location
2528  QueryProfile* profile = grid->getQueryProfile(location);
2529 
2530  // get the number of points used for that interpolation
2531  int nWeights = profile->getNCoefficients();
2532 
2533  // allocate some memory
2534  nodeIds.resize(nWeights);
2535  coefficients.resize(nWeights);
2536  depth.resize(NLAYERS);
2537  pvelocity.resize(NLAYERS);
2538  svelocity.resize(NLAYERS);
2539 
2540  // call the method
2541  profile->getData(nodeIds.data(), coefficients.data(), nWeights,
2542  depth.data(), pvelocity.data(), svelocity.data(), pgradient, sgradient);
2543  delete profile;
2544 }
2545 
2546 inline void SlbmInterface::getInterpolatedTransect(
2547  double lat[],
2548  double lon[],
2549  const int& nLatLon,
2550  int** nodeId,
2551  double** coefficients,
2552  int* nNeighbors,
2553  double depth[][NLAYERS],
2554  double pvelocity[][NLAYERS],
2555  double svelocity[][NLAYERS],
2556  double pgradient[NLAYERS],
2557  double sgradient[NLAYERS],
2558  int& nInvalid
2559  )
2560 {
2561  nInvalid = 0;
2562  for (int i=0; i<nLatLon; i++)
2563  {
2564  try
2565  {
2566  getInterpolatedPoint(lat[i], lon[i], nodeId[i],
2567  coefficients[i], nNeighbors[i], depth[i], pvelocity[i], svelocity[i],
2568  pgradient[i], sgradient[i]);
2569  }
2570  catch (SLBMException ex)
2571  {
2572  for (int j=0; j<nNeighbors[i]; j++)
2573  {
2574  nodeId[i][j] = -1;
2575  coefficients[i][j] = NA_VALUE;
2576  }
2577 
2578  for (int j=0; j<NLAYERS; j++)
2579  {
2580  depth[i][j] = NA_VALUE;
2581  pvelocity[i][j] = NA_VALUE;
2582  svelocity[i][j] = NA_VALUE;
2583  }
2584 
2585  pgradient[i] = NA_VALUE;
2586  sgradient[i] = NA_VALUE;
2587 
2588  ++nInvalid;
2589  }
2590  }
2591 }
2592 inline void SlbmInterface::getInterpolatedTransect(
2593  vector<double> lat,
2594  vector<double> lon,
2595  vector<vector<int> > &nodeId,
2596  vector<vector<double> > &coefficients,
2597  vector<vector<double> > &depth,
2598  vector<vector<double> > &pvelocity,
2599  vector<vector<double> > &svelocity,
2600  vector<double> &pgradient,
2601  vector<double> &sgradient,
2602  int &nInvalid
2603  )
2604 {
2605  nInvalid = 0;
2606 
2607  // check lat/lon sizes match
2608  int nLatLon;
2609  if (lat.size() == lon.size())
2610  {
2611  nLatLon = lat.size();
2612  }
2613  else
2614  {
2615  ostringstream os;
2616  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2617  os << endl << "ERROR in SlbmInterface::getInterpolatedTransect" << endl
2618  << "Size of lat and lon vectors do not match." << endl
2619  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2620  throw SLBMException(os.str(),114);
2621  }
2622 
2623  // resize vectors
2624  nodeId.resize(nLatLon);
2625  coefficients.resize(nLatLon);
2626  depth.resize(nLatLon);
2627  pvelocity.resize(nLatLon);
2628  svelocity.resize(nLatLon);
2629  pgradient.resize(nLatLon);
2630  sgradient.resize(nLatLon);
2631  for (int i=0; i<nLatLon; i++)
2632  {
2633  depth[i].resize(NLAYERS);
2634  pvelocity[i].resize(NLAYERS);
2635  svelocity[i].resize(NLAYERS);
2636  }
2637 
2638 
2639  for (int i=0; i<nLatLon; i++)
2640  {
2641  try
2642  {
2643  // create a new Location at the input lat/lon
2644  Location location(lat[i], lon[i], 0.);
2645 
2646  // interpolate a profile at that location
2647  QueryProfile profile(*grid, location); // create a query profile at that location
2648 
2649  // get the number of points used for that interpolation
2650  int nWeights = profile.getNCoefficients();
2651 
2652  // resize vectors
2653  nodeId[i].resize(nWeights);
2654  coefficients[i].resize(nWeights);
2655 
2656  getInterpolatedPoint(lat[i], lon[i], nodeId[i].data(), coefficients[i].data(), nWeights,
2657  depth[i].data(), pvelocity[i].data(), svelocity[i].data(), pgradient[i], sgradient[i]);
2658 
2659  }
2660  catch (SLBMException ex)
2661  {
2662  nodeId[i].push_back(-1);
2663  coefficients[i].push_back(NA_VALUE);
2664  pgradient[i] = NA_VALUE;
2665  sgradient[i] = NA_VALUE;
2666 
2667  for (int j=0; j<NLAYERS; j++)
2668  {
2669  depth[i][j] = NA_VALUE;
2670  pvelocity[i][j] = NA_VALUE;
2671  svelocity[i][j] = NA_VALUE;
2672  }
2673 
2674  ++nInvalid;
2675  }
2676  }
2677 }
2678 
2679 inline void SlbmInterface::getNGridNodes(int& n)
2680 {
2681  if (!grid)
2682  {
2683  n = -1;
2684  ostringstream os;
2685  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2686  os << endl << "ERROR in SlbmInterface::getNGridNodes" << endl
2687  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2688  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2689  throw SLBMException(os.str(),114);
2690  }
2691  n = grid->getNNodes();
2692 }
2693 
2694 inline int SlbmInterface::getNGridNodes()
2695 {
2696  if (!grid)
2697  {
2698  ostringstream os;
2699  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2700  os << endl << "ERROR in SlbmInterface::getNGridNodes" << endl
2701  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2702  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2703  throw SLBMException(os.str(),114);
2704  }
2705  return grid->getNNodes();
2706 }
2707 
2708 inline void SlbmInterface::getNHeadWavePoints(int& nHeadWavePoints)
2709 {
2710  if (!greatCircle)
2711  {
2712  nHeadWavePoints = -1;
2713  ostringstream os;
2714  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2715  os << endl << "ERROR in SlbmInterface::getNHeadWavePoints" << endl
2716  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2717  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2718  throw SLBMException(os.str(),113);
2719  }
2720  nHeadWavePoints = greatCircle->getNProfiles();
2721 }
2722 
2723 //inline void SlbmInterface::getMaxInterpotionNodes(int& maxNodes)
2724 //{
2725 // if (!greatCircle)
2726 // {
2727 // nHeadWavePoints = -1;
2728 // ostringstream os;
2729 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2730 // os << endl << "ERROR in SlbmInterface::getMaxInterpotionNodes" << endl
2731 // << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2732 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2733 // throw SLBMException(os.str(),113);
2734 // }
2735 //
2736 // maxNodes = 1;
2737 // for (int i=0; i<greatCircle->getNProfiles(); i++)
2738 // if (greatCircle->profiles[i]->getNodes().size() > maxNodes)
2739 // maxNodes = greatCircle->profiles[i]->getNodes().size();
2740 //}
2741 
2742 inline void SlbmInterface::getSlownessUncertainty( const int& phase, const double& distance, double& uncert )
2743 {
2744  if (!grid)
2745  {
2746  ostringstream os;
2747  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2748  os << endl << "ERROR in SlbmInterface::getSlownessUncertainty" << endl
2749  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2750  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2751  throw SLBMException(os.str(),114);
2752  }
2753 
2754  if (grid->getUncertainty()[phase][SH] == NULL)
2755  {
2756  ostringstream os;
2757  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2758  os << endl << "ERROR in SlbmInterface::getSlownessUncertainty" << endl
2759  << "Uncertainty object is invalid.." << endl
2760  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2761  throw SLBMException(os.str(),602);
2762  }
2763  else
2764  {
2765  // retrieve slowness uncertainty in sec/radian
2766  uncert = grid->getUncertainty()[phase][SH]->getUncertainty( distance );
2767  }
2768 }
2769 
2770 inline void SlbmInterface::getSlownessUncertainty(double& slownessUncertainty)
2771 {
2772  if (!grid)
2773  {
2774  ostringstream os;
2775  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2776  os << endl << "ERROR in SlbmInterface::getSlownessUncertainty" << endl
2777  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2778  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2779  throw SLBMException(os.str(),114);
2780  }
2781 
2782  getSlownessUncertainty(greatCircle->getPhase(),
2783  greatCircle->getDistance(), slownessUncertainty);
2784 }
2785 
2786 inline string SlbmInterface::getUncertaintyTable(const int& phase, const int& attribute)
2787 {
2788  if (!grid)
2789  {
2790  ostringstream os;
2791  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2792  os << endl << "ERROR in SlbmInterface::getUncertaintyTable" << endl
2793  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2794  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2795  throw SLBMException(os.str(),114);
2796  }
2797 
2798  if (grid->getUncertainty()[phase][attribute] == NULL)
2799  {
2800  ostringstream os;
2801  os << "No uncertainty information is available for phase " <<
2802  UncertaintyPIU::getPhase(phase) << " attribute " <<
2803  UncertaintyPIU::getAttribute(attribute) << endl;
2804  return os.str();
2805  }
2806  return grid->getUncertainty()[phase][attribute]->toStringTable();
2807 }
2808 
2809 inline string SlbmInterface::getUncertaintyFileFormat(const int& phase, const int& attribute)
2810 {
2811  if (!grid)
2812  {
2813  ostringstream os;
2814  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2815  os << endl << "ERROR in SlbmInterface::getUncertaintyFileFormat" << endl
2816  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2817  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2818  throw SLBMException(os.str(),114);
2819  }
2820 
2821  if (grid->getUncertainty()[phase][attribute] == NULL)
2822  {
2823  ostringstream os;
2824  os << "No uncertainty information is available for phase " << phase << " attribute " << attribute << endl;
2825  return os.str();
2826  }
2827  return grid->getUncertainty()[phase][attribute]->toStringFile();
2828 }
2829 
2830 inline void SlbmInterface::getZhaoParameters(double& Vm, double& Gm, double& H, double& C, double& Cm, int& udSign)
2831 {
2832  if (!greatCircle)
2833  {
2834  ostringstream os;
2835  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2836  os << endl << "ERROR in SlbmInterface::getZhaoParameters" << endl
2837  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2838  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2839  throw SLBMException(os.str(),113);
2840  }
2841  greatCircle->getZhaoParameters(Vm, Gm, H, C, Cm, udSign);
2842 }
2843 
2844 inline void SlbmInterface::getPgLgComponents(double& tTotal,
2845  double& tTaup, double& tHeadwave,
2846  double& pTaup, double& pHeadwave,
2847  double& trTaup, double& trHeadwave)
2848 {
2849  if (!greatCircle)
2850  {
2851  ostringstream os;
2852  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2853  os << endl << "ERROR in SlbmInterface::getPgLgComponents" << endl
2854  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2855  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2856  throw SLBMException(os.str(),113);
2857  }
2858  greatCircle->getPgLgComponents(tTotal, tTaup, tHeadwave,
2859  pTaup, pHeadwave, trTaup, trHeadwave);
2860 }
2861 
2862 
2863 inline void SlbmInterface::getNodeNeighbors(const int& nid, int neighbors[], int& nNeighbors)
2864 {
2865  if (!grid)
2866  {
2867  ostringstream os;
2868  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2869  os << endl << "ERROR in SlbmInterface::getNodeNeighbors" << endl
2870  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2871  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2872  throw SLBMException(os.str(),114);
2873  }
2874  grid->getNodeNeighbors(nid, neighbors, nNeighbors);
2875 }
2876 
2877 inline void SlbmInterface::getActiveNodeNeighbors(const int& nid, int neighbors[], int& nNeighbors)
2878 {
2879  if (!grid)
2880  {
2881  ostringstream os;
2882  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2883  os << endl << "ERROR in SlbmInterface::getActiveNodeNeighbors" << endl
2884  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2885  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2886  throw SLBMException(os.str(),114);
2887  }
2888  grid->getActiveNodeNeighbors(nid, neighbors, nNeighbors);
2889 }
2890 
2891 inline void SlbmInterface::getNodeNeighbors(const int& nid, vector<int>& neighbors)
2892 {
2893  if (!grid)
2894  {
2895  ostringstream os;
2896  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2897  os << endl << "ERROR in SlbmInterface::getNodeNeighbors" << endl
2898  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2899  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2900  throw SLBMException(os.str(),114);
2901  }
2902  grid->getNodeNeighbors(nid, neighbors);
2903 }
2904 
2905 inline void SlbmInterface::getActiveNodeNeighbors(const int& nid, vector<int>& neighbors)
2906 {
2907  if (!grid)
2908  {
2909  ostringstream os;
2910  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2911  os << endl << "ERROR in SlbmInterface::getActiveNodeNeighbors" << endl
2912  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2913  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2914  throw SLBMException(os.str(),114);
2915  }
2916  grid->getActiveNodeNeighbors(nid, neighbors);
2917 }
2918 
2919 inline void SlbmInterface::getNodeNeighborInfo(const int& nid, int neighbors[],
2920  double distance[], double azimuth[], int& nNeighbors)
2921 {
2922  if (!grid)
2923  {
2924  ostringstream os;
2925  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2926  os << endl << "ERROR in SlbmInterface::getNodeNeighborInfo" << endl
2927  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2928  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2929  throw SLBMException(os.str(),114);
2930  }
2931  grid->getNodeNeighborInfo(nid, neighbors, distance, azimuth, nNeighbors);
2932 }
2933 
2934 inline void SlbmInterface::getNodeNeighborInfo(const int nid,
2935  vector<int> &neighbors, vector<double> &distance, vector<double> &azimuth)
2936 {
2937  if (!grid)
2938  {
2939  ostringstream os;
2940  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2941  os << endl << "ERROR in SlbmInterface::getNodeNeighborInfo" << endl
2942  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2943  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2944  throw SLBMException(os.str(),114);
2945  }
2946  grid->getNodeNeighborInfo(nid, neighbors, distance, azimuth);
2947 }
2948 
2949 inline void SlbmInterface::getActiveNodeNeighborInfo(const int& nid, int neighbors[],
2950  double distance[], double azimuth[], int& nNeighbors)
2951 {
2952  if (!grid)
2953  {
2954  ostringstream os;
2955  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2956  os << endl << "ERROR in SlbmInterface::getActiveNodeNeighborInfo" << endl
2957  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2958  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2959  throw SLBMException(os.str(),114);
2960  }
2961  grid->getActiveNodeNeighborInfo(nid, neighbors, distance, azimuth, nNeighbors);
2962 }
2963 
2964 inline void SlbmInterface::getActiveNodeNeighborInfo(const int& nid,
2965  vector<int>& neighbors, vector<double>& distance, vector<double>& azimuth)
2966 {
2967  if (!grid)
2968  {
2969  ostringstream os;
2970  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2971  os << endl << "ERROR in SlbmInterface::getActiveNodeNeighborInfo" << endl
2972  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2973  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2974  throw SLBMException(os.str(),114);
2975  }
2976  grid->getActiveNodeNeighborInfo(nid, neighbors, distance, azimuth);
2977 }
2978 
2979 inline void SlbmInterface::getNodeSeparation(const int& node1, const int& node2, double& distance)
2980 {
2981  if (!grid)
2982  {
2983  ostringstream os;
2984  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2985  os << endl << "ERROR in SlbmInterface::getNodeSeparation" << endl
2986  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2987  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2988  throw SLBMException(os.str(),114);
2989  }
2990  grid->getNodeSeparation(node1, node2, distance);
2991 }
2992 
2993 inline void SlbmInterface::getNodeAzimuth(const int& node1, const int& node2, double& azimuth)
2994 {
2995  if (!grid)
2996  {
2997  ostringstream os;
2998  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2999  os << endl << "ERROR in SlbmInterface::getNodeAzimuth" << endl
3000  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3001  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3002  throw SLBMException(os.str(),114);
3003  }
3004  grid->getNodeAzimuth(node1, node2, azimuth);
3005 }
3006 
3007 inline void SlbmInterface::initializeActiveNodes(const double& latmin, const double& lonmin,
3008  const double& latmax, const double& lonmax)
3009 {
3010  if (!grid)
3011  {
3012  ostringstream os;
3013  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3014  os << endl << "ERROR in SlbmInterface::initializeActiveNodes" << endl
3015  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3016  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3017  throw SLBMException(os.str(),114);
3018  }
3019  grid->initializeActiveNodes(latmin, lonmin, latmax, lonmax);
3020 }
3021 
3022 inline void SlbmInterface::initializeActiveNodes(double* lat, double* lon, const int& npoints, const bool& inDegrees)
3023 {
3024  vector<double*> v;
3025  v.reserve(npoints);
3026  if (inDegrees)
3027  for (int i=0; i<npoints; ++i)
3028  v.push_back(GeoTessUtils::getVectorDegrees(lat[i], lon[i]));
3029  else
3030  for (int i=0; i<npoints; ++i)
3031  v.push_back(GeoTessUtils::getVector(lat[i], lon[i]));
3032  initializeActiveNodes(v);
3033 }
3034 inline void SlbmInterface::initializeActiveNodes(const vector<double> lat, const vector<double> lon, const bool& inDegrees)
3035 {
3036  int npoints;
3037  if (lat.size() == lon.size())
3038  {
3039  npoints = lat.size();
3040  }
3041  else
3042  {
3043  ostringstream os;
3044  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3045  os << endl << "ERROR in SlbmInterface::getInterpolatedTransect" << endl
3046  << "Size of lat and lon vectors do not match." << endl
3047  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3048  throw SLBMException(os.str(),114);
3049  }
3050  vector<double*> v;
3051  v.reserve(npoints);
3052  if (inDegrees)
3053  for (int i=0; i<npoints; ++i)
3054  v.push_back(GeoTessUtils::getVectorDegrees(lat[i], lon[i]));
3055  else
3056  for (int i=0; i<npoints; ++i)
3057  v.push_back(GeoTessUtils::getVector(lat[i], lon[i]));
3058  initializeActiveNodes(v);
3059 }
3060 
3061 inline void SlbmInterface::initializeActiveNodes(vector<double*>& unitVectors)
3062 {
3063  initializeActiveNodes(new GeoTessPolygon(unitVectors));
3064 }
3065 
3066 inline void SlbmInterface::initializeActiveNodes(GeoTessPolygon* polygon)
3067 {
3068  grid->initializeActiveNodes(polygon);
3069 }
3070 
3071 inline void SlbmInterface::initializeActiveNodes(const string& fileName)
3072 {
3073  grid->initializeActiveNodes(new GeoTessPolygon(fileName));
3074 }
3075 
3076 inline void SlbmInterface::clearActiveNodes()
3077 {
3078  grid->clearActiveNodes();
3079 }
3080 
3081 
3082 inline int SlbmInterface::getNActiveNodes()
3083 {
3084  if (!grid)
3085  {
3086  ostringstream os;
3087  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3088  os << endl << "ERROR in SlbmInterface::nextActiveNode" << endl
3089  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3090  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3091  throw SLBMException(os.str(),114);
3092  }
3093  return grid->getNActiveNodes();
3094 }
3095 
3096 inline int SlbmInterface::getGridNodeId(int activeNodeId)
3097 {
3098  if (!grid)
3099  {
3100  ostringstream os;
3101  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3102  os << endl << "ERROR in SlbmInterface::getGridNodeId" << endl
3103  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3104  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3105  throw SLBMException(os.str(),114);
3106  }
3107  return grid->getGridNodeId(activeNodeId);
3108 }
3109 
3110 inline int SlbmInterface::getActiveNodeId(int gridNodeId)
3111 {
3112  if (!grid)
3113  {
3114  ostringstream os;
3115  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3116  os << endl << "ERROR in SlbmInterface::getActiveNodeId" << endl
3117  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3118  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3119  throw SLBMException(os.str(),114);
3120  }
3121  return grid->getActiveNodeId(gridNodeId);
3122 }
3123 
3124 inline void SlbmInterface::getNodeHitCount(const int& nodeId, int& hitCount)
3125 {
3126  if (!grid)
3127  {
3128  ostringstream os;
3129  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3130  os << endl << "ERROR in SlbmInterface::getNodeHitCount" << endl
3131  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3132  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3133  throw SLBMException(os.str(),114);
3134  }
3135  grid->getNodeHitCount(nodeId, hitCount);
3136 }
3137 
3138 inline void SlbmInterface::clearNodeHitCount()
3139 {
3140  if (!grid)
3141  {
3142  ostringstream os;
3143  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3144  os << endl << "ERROR in SlbmInterface::clearNodeHitCount" << endl
3145  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3146  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3147  throw SLBMException(os.str(),114);
3148  }
3149  grid->clearNodeHitCount();
3150 }
3151 
3152 inline void SlbmInterface::getAverageMantleVelocity(const int& type, double& velocity)
3153 {
3154  if (!grid)
3155  {
3156  ostringstream os;
3157  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3158  os << endl << "ERROR in SlbmInterface::setAverageMantleVelocity" << endl
3159  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3160  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3161  throw SLBMException(os.str(),114);
3162  }
3163  velocity = grid->getAverageMantleVelocity(type);
3164 }
3165 
3166 inline double SlbmInterface::getAverageMantleVelocity(const int &type)
3167 {
3168  if (!grid)
3169  {
3170  ostringstream os;
3171  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3172  os << endl << "ERROR in SlbmInterface::setAverageMantleVelocity" << endl
3173  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3174  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3175  throw SLBMException(os.str(),114);
3176  }
3177  return grid->getAverageMantleVelocity(type);
3178 }
3179 
3180 inline void SlbmInterface::setAverageMantleVelocity(const int& type,
3181  const double& velocity)
3182 {
3183  if (!grid)
3184  {
3185  ostringstream os;
3186  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3187  os << endl << "ERROR in SlbmInterface::setAverageMantleVelocity" << endl
3188  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3189  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3190  throw SLBMException(os.str(),114);
3191  }
3192  grid->setAverageMantleVelocity(type, velocity);
3193 }
3194 
3195 inline void SlbmInterface::getTessId(string& tessId)
3196 {
3197  if (!grid)
3198  {
3199  ostringstream os;
3200  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3201  os << endl << "ERROR in SlbmInterface::getTessId" << endl
3202  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3203  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3204  throw SLBMException(os.str(),114);
3205  }
3206  tessId = grid->getTessId();
3207 }
3208 
3209 inline void SlbmInterface::getFractionActive(double& fractionActive)
3210 {
3211  if (!isValid())
3212  {
3213  fractionActive = NA_VALUE;
3214  ostringstream os;
3215  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3216  os << endl << "ERROR in SlbmInterface::getFractionActive" << endl
3217  << "GreatCircle is invalid." << endl
3218  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3219  throw SLBMException(os.str(),113);
3220  }
3221  fractionActive = greatCircle->getFractionActive();
3222 }
3223 
3224 inline string SlbmInterface::getClassCount()
3225 {
3226  ostringstream os;
3227  os << "Class counts:" << endl;
3228  os << "GreatCircle = " << GreatCircle::getClassCount() << endl;
3229  os << "GridProfile = " << GridProfile::getClassCount() << endl;
3230  os << "GeoStack = " << GeoStack::getClassCount() << endl;
3231  os << "InterpolatedProfile = " << InterpolatedProfile::getClassCount() << endl;
3232  os << "CrustalProfile = " << CrustalProfile::getClassCount() << endl;
3233  os << "LayerProfile = " << LayerProfile::getClassCount() << endl;
3234  os << "QueryProfile = " << QueryProfile::getClassCount() << endl;
3235  os << "Location = " << Location::getClassCount() << endl;
3236  return os.str();
3237 }
3238 
3239 inline const string& SlbmInterface::getModelPath() const { return grid->getModelPath(); }
3240 
3241 inline void SlbmInterface::getPiercePointSource(double& lat, double& lon, double& depth)
3242 {
3243  if (!isValid())
3244  {
3245  lat = lon = NA_VALUE;
3246  ostringstream os;
3247  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3248  os << endl << "ERROR in SlbmInterface::getPiercePointSource" << endl
3249  << "GreatCircle is invalid." << endl
3250  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3251  throw SLBMException(os.str(),113);
3252  }
3253 
3254 
3255  if ((greatCircle->getPhase() == Pg || greatCircle->getPhase() == Lg) && !greatCircle->getSourceProfile()->isInCrust())
3256  {
3257  lat = lon = NA_VALUE;
3258  ostringstream os;
3259  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3260  os << endl << "ERROR in SlbmInterface::getPiercePointSource" << endl
3261  << "Cannot compute moho pierce points for crustal phase (" << getPhase() << ") when source is in the mantle." << endl
3262  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3263  throw SLBMException(os.str(),113);
3264  }
3265 
3266  Location loc;
3267  greatCircle->getGreatCircleLocation(greatCircle->getSourceDistance(), loc);
3268 
3269  lat = loc.getLat();
3270  lon = loc.getLon();
3271 
3272  QueryProfile* profile = NULL;
3273  profile = grid->getQueryProfile(loc);
3274  if (greatCircle->getPhase() == Pg || greatCircle->getPhase() == Lg)
3275  depth = profile->getDepth()[MIDDLE_CRUST_G];
3276  else
3277  depth = profile->getDepth()[MANTLE];
3278  delete profile;
3279 }
3280 
3281 inline void SlbmInterface::getPiercePointReceiver(double& lat, double& lon, double& depth)
3282 {
3283  if (!isValid())
3284  {
3285  lat = lon = NA_VALUE;
3286  ostringstream os;
3287  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3288  os << endl << "ERROR in SlbmInterface::getPiercePointReceiver" << endl
3289  << "GreatCircle is invalid." << endl
3290  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3291  throw SLBMException(os.str(),113);
3292  }
3293 
3294 
3295  if ((greatCircle->getPhase() == Pg || greatCircle->getPhase() == Lg) && !greatCircle->getReceiverProfile()->isInCrust())
3296  {
3297  lat = lon = NA_VALUE;
3298  ostringstream os;
3299  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3300  os << endl << "ERROR in SlbmInterface::getPiercePointReceiver" << endl
3301  << "Cannot compute moho pierce points for crustal phase (" << getPhase() << ") when receiver is in the mantle." << endl
3302  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3303  throw SLBMException(os.str(),113);
3304  }
3305 
3306  Location loc;
3307  greatCircle->getGreatCircleLocation(
3308  greatCircle->getDistance()-greatCircle->getReceiverDistance(), loc);
3309 
3310  lat = loc.getLat();
3311  lon = loc.getLon();
3312 
3313  QueryProfile* profile = NULL;
3314  profile = grid->getQueryProfile(loc);
3315  if (greatCircle->getPhase() == Pg || greatCircle->getPhase() == Lg)
3316  depth = profile->getDepth()[MIDDLE_CRUST_G];
3317  else
3318  depth = profile->getDepth()[MANTLE];
3319  delete profile;
3320 }
3321 
3322 inline void SlbmInterface::getGreatCirclePoints(
3323  const double& aLat,
3324  const double& aLon,
3325  const double& bLat,
3326  const double& bLon,
3327  const int& npoints,
3328  double latitude[],
3329  double longitude[])
3330 {
3331  Location a(aLat, aLon, 0.);
3332  Location b(bLat, bLon, 0.);
3333  double dx = a.distance(b)/(npoints-1);
3334 
3335  double moveDirection[3];
3336  a.vectorTripleProduct(b, moveDirection);
3337  for (int i=0; i<npoints; i++)
3338  {
3339  a.move(moveDirection, i*dx, b);
3340  latitude[i] = b.getLat();
3341  longitude[i] = b.getLon();
3342  }
3343 }
3344 
3345 inline void SlbmInterface::getGreatCirclePoints(
3346  const double& aLat,
3347  const double& aLon,
3348  const double& bLat,
3349  const double& bLon,
3350  const int& npoints,
3351  vector<double> &latitude,
3352  vector<double> &longitude)
3353 {
3354  latitude.resize(npoints);
3355  longitude.resize(npoints);
3356 
3357  Location a(aLat, aLon, 0.);
3358  Location b(bLat, bLon, 0.);
3359  double dx = a.distance(b)/(npoints-1);
3360 
3361  double moveDirection[3];
3362  a.vectorTripleProduct(b, moveDirection);
3363  for (int i=0; i<npoints; i++)
3364  {
3365  a.move(moveDirection, i*dx, b);
3366  latitude[i] = b.getLat();
3367  longitude[i] = b.getLon();
3368  }
3369 }
3370 
3371 inline void SlbmInterface::getGreatCirclePointsOnCenters(
3372  const double& aLat,
3373  const double& aLon,
3374  const double& bLat,
3375  const double& bLon,
3376  const int& npoints,
3377  double latitude[],
3378  double longitude[])
3379 {
3380  Location a(aLat, aLon, 0.);
3381  Location b(bLat, bLon, 0.);
3382  double dx = a.distance(b)/npoints;
3383 
3384  double moveDirection[3];
3385  a.vectorTripleProduct(b, moveDirection);
3386  for (int i=0; i<npoints; i++)
3387  {
3388  a.move(moveDirection, (i+0.5)*dx, b);
3389  latitude[i] = b.getLat();
3390  longitude[i] = b.getLon();
3391  }
3392 }
3393 
3394 inline void SlbmInterface::getGreatCirclePointsOnCenters(
3395  const double& aLat,
3396  const double& aLon,
3397  const double& bLat,
3398  const double& bLon,
3399  const int& npoints,
3400  vector<double> &latitude,
3401  vector<double> &longitude)
3402 {
3403  latitude.resize(npoints);
3404  longitude.resize(npoints);
3405 
3406  Location a(aLat, aLon, 0.);
3407  Location b(bLat, bLon, 0.);
3408  double dx = a.distance(b)/npoints;
3409 
3410  double moveDirection[3];
3411  a.vectorTripleProduct(b, moveDirection);
3412  for (int i=0; i<npoints; i++)
3413  {
3414  a.move(moveDirection, (i+0.5)*dx, b);
3415  latitude[i] = b.getLat();
3416  longitude[i] = b.getLon();
3417  }
3418 }
3419 
3420 inline void SlbmInterface::getDistAz(const double& aLat, const double& aLon,
3421  const double& bLat, const double& bLon,
3422  double& distance, double& azimuth, const double& naValue)
3423 {
3424  Location ptA(aLat, aLon);
3425  Location ptB(bLat, bLon);
3426  distance = ptA.distance(ptB);
3427  azimuth = ptA.azimuth(ptB, naValue);
3428 }
3429 
3430 inline void SlbmInterface::movePoint(const double& aLat, const double& aLon,
3431  const double& distance, const double& azimuth,
3432  double& bLat, double& bLon)
3433 {
3434  Location ptA(aLat, aLon);
3435  Location ptB;
3436  ptA.move(azimuth, distance, ptB);
3437  bLat = ptB.getLat();
3438  bLon = ptB.getLon();
3439 }
3440 
3441 inline void SlbmInterface::setInterpolatorType(const string& interpolatorType)
3442 {
3443  if (!grid)
3444  {
3445  ostringstream os;
3446  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3447  os << endl << "ERROR in SlbmInterface::setInterpolatorType" << endl
3448  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3449  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3450  throw SLBMException(os.str(),114);
3451  }
3452  grid->setInterpolatorType(interpolatorType);
3453 }
3454 
3455 inline string SlbmInterface::getInterpolatorType()
3456 {
3457  if (!grid)
3458  {
3459  ostringstream os;
3460  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
3461  os << endl << "ERROR in SlbmInterface::getInterpolatorType" << endl
3462  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
3463  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
3464  throw SLBMException(os.str(),114);
3465  }
3466  return grid->getInterpolatorType();
3467 }
3468 
3469 
3470 
3471 } // end slbm namespace
3472 
3473 #endif // SlbmInterface_H
#define SLBM_EXP
Definition: SLBMGlobals.h:181
The GreatCircle class manages information related to a great circle path between two Locations on the...
Definition: GreatCircle.h:114
A 2 dimensional, horizontal grid of GirdProfile objects.
Definition: Grid.h:91
The Location Class manages a single point in/on the Earth, which is described by the GRS80 ellipsoid.
Definition: Location.h:78
double getLon() const
Retrieve the longitude of this Location. Value will be between -PI and PI radians.
Definition: Location.h:548
double distance(const Location &other) const
Definition: Location.h:583
double azimuth(const Location &other) const
Find the azimuth from this Location to some other Location. Result will be between 0 and 2*PI radians...
Definition: Location.h:593
bool vectorTripleProduct(const Location &other, double vtp[]) const
Compute the vector triple product (this x other) x this, normalized to unit length....
Definition: Location.h:691
void move(const double &azimuth, const double &distance, Location &loc) const
Retrieve a Location that is a specified distance away from this Location, in a specified direction.
Definition: Location.h:638
double getDepth() const
Definition: Location.h:573
double getLat() const
Retrieve the geographic latitude of this Location, radians.
Definition: Location.h:542
An InterpolatedProfile object that also has information about the the P and S wave velocity as a func...
Definition: QueryProfile.h:68
void getData(int *nodeIds, double *coefficients, int &nNeighbors, double *depths, double *pvelocities, double *svelocities, double &pgradient, double &sgradient)
Retrieve all the interval depth and velocity information contained in this QueryProfile object.
Definition: QueryProfile.h:185
double * getDepth()
Retrieve the depth of the top of the k'th interval, in km.
Definition: QueryProfile.h:129
An Exception class for Grid and related objects.
Definition: SLBMException.h:55
The primary interface to the SLBM library, providing access to all supported functionality.
Definition: SlbmInterface.h:77
SlbmInterface(const double &earthRadius)
Parameterized constructor. Instantiates an SlbmInterface object that is only partly based on an ellip...
double getRcvLon()
getter for read-only attributes property getter for read-only attributes property
double getSrcLon()
getter for read-only attributes property getter for read-only attributes property
void getActiveNodeData(const int &nodeId, double &latitude, double &longitude, double depth[NLAYERS], double pvelocity[NLAYERS], double svelocity[NLAYERS], double gradient[2])
Retrieve the lat (radians), lon (radians), interface depths (km), P and S wave interval velocities (k...
double getHeadwaveDistance()
void setDelDepth(const double &del_depth)
Change the value of step change in depth used to compute depth derivatives (km)
void setDelDistance(const double &del_distance)
Change the value of step change in distance used to compute horizontal derivatives(in radians).
bool isEqual(SlbmInterface *other)
Check if the model and/or greatCircle attached to the current SlbmInterface is equal that in another ...
double getSourceDistance()
int getNodeHitCount(const int &nodeId)
void getRayParameter(double &ray_parameter)
Retrieve the ray parameter.
void setGridData(const int &nodeId, double depths[NLAYERS], double pvelocity[NLAYERS], double svelocity[NLAYERS], double gradient[2])
Modify the velocity and gradient information associated with a specified node in the Grid.
GreatCircle * greatCircle
The most recently requested GreatCircle object.
void saveVelocityModel(const string &modelFileName, const int &format=4)
Save the velocity model currently in memory to the specified file.
SlbmInterface()
Default constructor. Instantiates an SlbmInterface object based on an ellipsoidal earth.
bool isValid()
Returns true if the current GreatCirlce object has been instantiated and is ready to be interrogated.
void setPathIncrement(const double &pathIncrement)
Set the desired spacing of great circle nodes along the head wave interface, in radians.
void getDelDistance(double &del_distance)
Retrieve the value of step change in distance used to compute horizontal derivatives (radians)
void setActiveNodeData(const int &nodeId, double depths[NLAYERS], double pvelocity[NLAYERS], double svelocity[NLAYERS], double gradient[2])
Modify the depth, velocity and gradient information associated with a specified active node in the Gr...
void createGreatCircle(const int &phase, const double &sourceLat, const double &sourceLon, const double &sourceDepth, const double &receiverLat, const double &receiverLon, const double &receiverDepth)
Instantiate a new GreatCircle object between two locations.
Grid * grid
The Grid object that stores the velocity model.
void getTravelTimeUncertainty1D(double &travelTimeUncertainty)
Retrieve travel time uncertainty in sec using the phase and distance specified in last call to getGre...
Grid * getGridObject()
Retrieve a pointer to the Grid object.
static bool modelsEqual(const string modelPath1, const string modelPath2)
Check if two models are equal.
void loadVelocityModel(const string &modelPath)
Load the velocity model into memory from the specified file or directory. This method automatically d...
void getTurningRadius(double &turning_radius)
Retrieve the turning radius of the ray.
static void getCHMax(double &chMax)
Retrieve the current value of chMax. c is the zhao c parameter and h is the turning depth of the ray ...
static void getMaxDepth(double &maxDepth)
Retrieve the current value for the maximum source depth, in km.
static void getMaxDistance(double &maxDistance)
Retrieve the current value for the maximum source-receiver separation, in radians.
GeoTessModelSLBM * getModelObject()
Retrieve a pointer to the GeoTessModelSLBM object.
void clear()
Delete the current GreatCircle object from memory and clear the pool of stored CrustalProfile objects...
static double getMaxDistance()
void getDistance(double &distance)
Retrieve the source-receiver separation, in radians.
int getBufferSize() const
Returns the size of a DataBuffer object required to store this SLBMInterface objects model data.
double getSrcLat()
getter for read-only attributes property getter for read-only attributes property
static double CH_MAX
double getRcvLat()
getter for read-only attributes property getter for read-only attributes property
static double getMaxDepth()
void getTravelTimeUncertainty(const int &phase, const double &distance, double &uncert)
Retrieve the travel time uncertainty in sec for specified phase, distance (in radians).
double getSrcDep()
getter for read-only attributes property getter for read-only attributes property
double getRcvDep()
getter for read-only attributes property getter for read-only attributes property
virtual ~SlbmInterface()
Destructor.
string getVersion()
Retrieve the SLBM Version number.
void getTravelTimeUncertainty(double &travelTimeUncertainty, bool calcRandomError=false)
Retrieve travel time uncertainty in sec using the phase and distance specified in last call to getGre...
void getPathIncrement(double &pathIncrement)
Retrieve the current value of the spacing of great circle nodes along the head wave interface,...
double getPathIncrement()
Retrieve the current value of the spacing of great circle nodes along the head wave interface,...
static double getCHMax()
void getDelDepth(double &del_depth)
Retrieve the value of step change in depth used to compute depth derivatives (km)
void getGridData(const int &nodeId, double &latitude, double &longitude, double depth[NLAYERS], double pvelocity[NLAYERS], double svelocity[NLAYERS], double gradient[2])
Retrieve the lat (radians), lon (radians), interface depths (km), P and S wave interval velocities (k...
void get_dtt_ddist(double &slowness)
Retrieve the horizontal slowness, i.e., the derivative of travel time wrt to receiver-source distance...
static void setMaxDistance(const double &maxDistance)
Set the maximum source-receiver separation for Pn/Sn phase, in radians.
static void setCHMax(const double &chMax)
Set the value of chMax. c is the zhao c parameter and h is the turning depth of the ray below the moh...
static void setMaxDepth(const double &maxDepth)
Set the maximum source depth for Pn/Sn phase, in km.
bool valid
true if the current GreatCirlce object has been instantiated and is ready to be interrogated.
string getPhase()
Retrieve the phase specified in last call to createGreatCircle().
double getReceiverDistance()
double getHeadwaveDistanceKm()
A byte array container used to hold binary data in the same manner as disk based file system.
Definition: DataBuffer.h:81