RSTT  3.1.0
Regional Seismic Travel Time
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 "Uncertainty.h"
53 #include "UncertaintyPathDep.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 
111  string getVersion() { return SlbmVersion; };
112 
124  void loadVelocityModel(const string& modelPath);
125 
170  void saveVelocityModel(const string& modelFileName, const int& format=4);
171 
173 
190  void loadVelocityModelBinary(const string& modelPath)
191  { loadVelocityModel(modelPath); };
192 
209  void specifyOutputDirectory(const string& directoryName);
210 
220  void saveVelocityModelBinary();
221 
230  void loadVelocityModelBinary(util::DataBuffer& buffer);
231 
237  void saveVelocityModelBinary(util::DataBuffer& buffer);
238 
240 
246  int getBufferSize() const;
247 
255  void setInterpolatorType(const string& interpolatorType);
256 
262  string getInterpolatorType();
263 
279  void createGreatCircle(const string& phase,
280  const double& sourceLat,
281  const double& sourceLon,
282  const double& sourceDepth,
283  const double& receiverLat,
284  const double& receiverLon,
285  const double& receiverDepth);
286 
302  void createGreatCircle(const int& phase,
303  const double& sourceLat,
304  const double& sourceLon,
305  const double& sourceDepth,
306  const double& receiverLat,
307  const double& receiverLon,
308  const double& receiverDepth);
309 
329  void clear();
330 
336  bool isValid() { return valid; };
337 
341  string getPhase() { return greatCircle->getPhaseString(); };
342 
348  void getDistance(double& distance) { distance = getDistance(); };
349 
354  double getDistance();
355 
363  void getSourceDistance(double& dist);
364 
372  void getReceiverDistance(double& dist);
373 
385  void getHeadwaveDistance(double& dist);
386 
399  void getHeadwaveDistanceKm(double& dist);
400 
409  void getTravelTime(double& travelTime);
410 
425  void getTravelTimeComponents(double& tTotal, double& tSource, double& tReceiver,
426  double& tHeadwave, double& tGradient);
427 
433  void getSlowness(double& slowness);
434 
440  void get_dtt_ddist(double& slowness) { getSlowness(slowness); };
441  double get_dtt_ddist() { double x; getSlowness(x); return x; };
442 
449  void get_dtt_dlat(double& dtt_dlat);
450  double get_dtt_dlat() { double x; get_dtt_dlat(x); return x; };
451 
458  void get_dtt_dlon(double& dtt_dlon);
459  double get_dtt_dlon() { double x; get_dtt_dlon(x); return x; };
460 
467  void get_dtt_ddepth(double& dtt_ddepth);
468  double get_dtt_ddepth() { double x; get_dtt_ddepth(x); return x; };
469 
476  //void get_dtt_dlat_fast(double& dtt_dlat);
477 
484  //void get_dtt_dlon_fast(double& dtt_dlon);
485 
486 // //! \brief Retrieve the derivative of horizontal slowness wrt to source-receiver distance,
487 // //! in seconds/radian^2.
488 // //!
489 // //! Retrieve the derivative of horizontal slowness wrt to source-receiver distance,
490 // //! in seconds/radian^2.
491 // void get_dsh_ddist(double& dsh_ddist);
492 // double get_dsh_ddist() { double x; get_dsh_ddist(x); return x; };
493 //
494 // //! \brief Retrieve the derivative of horizontal slowness wrt to source latitude,
495 // //! in seconds/radian^2.
496 // //!
497 // //! Retrieve the derivative of horizontal slowness wrt to source latitude,
498 // //! in seconds/radian^2.
499 // //! @param dsh_dlat the derivative of horizontal slowness wrt to source latitude.
500 // void get_dsh_dlat(double& dsh_dlat);
501 //
502 // //! \brief Retrieve the derivative of horizontal slowness wrt to source longitude,
503 // //! in seconds/radian^2.
504 // //!
505 // //! Retrieve the derivative of horizontal slowness wrt to source longitude,
506 // //! in seconds/radian^2.
507 // //! @param dsh_dlon the derivative of horizontal slowness wrt to source longitude.
508 // void get_dsh_dlon(double& dsh_dlon);
509 //
510 // //! \brief Retrieve the derivative of horizontal slowness wrt to source depth,
511 // //! in seconds/radian-km.
512 // //!
513 // //! Retrieve the derivative of horizontal slowness wrt to source depth,
514 // //! in seconds/radian-km.
515 // //! @param dsh_ddepth the derivative of horizontal slowness wrt to source depth.
516 // void get_dsh_ddepth(double& dsh_ddepth);
517 // double get_dsh_ddepth() { double x; get_dsh_ddepth(x); return x; }
518 
556  void getWeights(int nodeId[], double weight[], int& nWeights);
557 
596  void getActiveNodeWeights(int nodeId[], double weight[], int& nWeights);
597 
631  void getWeights(vector<int>& nodeId, vector<double>& weight);
632 
665  void getActiveNodeWeights(vector<int>& nodeId, vector<double>& weight);
666 
680  void getWeightsSource(int nodeids[], double weights[], int& nWeights);
681 
695  void getActiveNodeWeightsSource(int nodeids[], double weights[], int& nWeights);
696 
707  void getWeightsReceiver(int nodeids[], double weights[], int& nWeights);
708 
719  void getActiveNodeWeightsReceiver(int nodeids[], double weights[], int& nWeights);
720 
742  string toString(const int& verbosity);
743 
747  void getNGridNodes(int& n);
748 
752  int getNGridNodes();
753 
763  void getNHeadWavePoints(int& nHeadWavePoints);
764 
785  const int& nodeId,
786  double& latitude,
787  double& longitude,
788  double depth[NLAYERS],
789  double pvelocity[NLAYERS],
790  double svelocity[NLAYERS],
791  double gradient[2]);
792 
813  const int& nodeId,
814  double& latitude,
815  double& longitude,
816  double depth[NLAYERS],
817  double pvelocity[NLAYERS],
818  double svelocity[NLAYERS],
819  double gradient[2]);
820 
837  const int& nodeId,
838  double depths[NLAYERS],
839  double pvelocity[NLAYERS],
840  double svelocity[NLAYERS],
841  double gradient[2]);
842 
859  const int& nodeId,
860  double depths[NLAYERS],
861  double pvelocity[NLAYERS],
862  double svelocity[NLAYERS],
863  double gradient[2]);
864 
901  void getGreatCircleData(
902  string& phase,
903  double& actual_path_increment,
904  double sourceDepth[NLAYERS],
905  double sourceVelocity[NLAYERS],
906  double receiverDepth[NLAYERS],
907  double receiverVelocity[NLAYERS],
908  int& npoints,
909  double headWaveVelocity[],
910  double gradient[]
911  );
912 
926  void getGreatCircleLocations(double lat[], double lon[], double depth[], int& npoints);
927 
955  void getGreatCircleNodeInfo(
956  int** neighbors,
957  double** coefficients,
958  const int& maxpoints,
959  const int& maxnodes,
960  int& npoints,
961  int* nnodes
962  );
963 
987  void getInterpolatedPoint(
988  const double& lat,
989  const double& lon,
990  int* nodeId,
991  double* coefficients,
992  int& nnodes,
993  double depth[NLAYERS],
994  double pvelocity[NLAYERS],
995  double svelocity[NLAYERS],
996  double& pgradient,
997  double& sgradient
998  );
999 
1025  void getInterpolatedTransect(
1026  double lat[],
1027  double lon[],
1028  const int& nLatLon,
1029  int** neighbors,
1030  double** coefficients,
1031  int* nNeighbors,
1032  double depth[][NLAYERS],
1033  double pvelocity[][NLAYERS],
1034  double svelocity[][NLAYERS],
1035  double pgradient[NLAYERS],
1036  double sgradient[NLAYERS],
1037  int& nInvalid
1038  );
1039 
1052  void initializeActiveNodes(const double& latmin, const double& lonmin,
1053  const double& latmax, const double& lonmax);
1054 
1071  void initializeActiveNodes(const string& polygonFileName);
1072 
1085  void initializeActiveNodes(GeoTessPolygon* polygon);
1086 
1108  void initializeActiveNodes(double* lat, double* lon, const int& npoints, const bool& inDegrees=true);
1109 
1125  void initializeActiveNodes(vector<double*>& unitVectors);
1126 
1130  int getNActiveNodes();
1131 
1134  void clearActiveNodes();
1135 
1141  int getGridNodeId(int activeNodeId);
1142 
1148  int getActiveNodeId(int gridNodeId);
1149 
1158  void getNodeHitCount(const int& nodeId, int& hitCount);
1159 
1165  void clearNodeHitCount();
1166 
1174  void getNodeNeighbors(const int& nid, int neighbors[], int& nNeighbors);
1175 
1183  void getActiveNodeNeighbors(const int& nid, int neighbors[], int& nNeighbors);
1184 
1189  void getNodeNeighbors(const int& nid, vector<int>& neighbors);
1190 
1195  void getActiveNodeNeighbors(const int& nid, vector<int>& neighbors);
1196 
1204  void getNodeNeighborInfo(const int& nid, int neighbors[], double distance[],
1205  double azimuth[], int& nNeighbors);
1206 
1214  void getActiveNodeNeighborInfo(const int& nid, int neighbors[], double distance[],
1215  double azimuth[], int& nNeighbors);
1216 
1224  void getActiveNodeNeighborInfo(const int& nid,
1225  vector<int>& neighbors, vector<double>& distance, vector<double>& azimuth);
1226 
1230  void getNodeSeparation(const int& node1, const int& node2, double& distance);
1231 
1235  void getNodeAzimuth(const int& node1, const int& node2, double& azimuth);
1236 
1245  GreatCircle* getGreatCircleObject();
1246 
1251  Grid* getGridObject() { return grid; };
1252 
1261  void getTravelTimeUncertainty( const int& phase, const double& distance, double& uncert );
1262 
1271  void getTravelTimeUncertainty(double& travelTimeUncertainty, bool calcRandomError = false);
1272 
1279  void getTravelTimeUncertainty1D(double& travelTimeUncertainty);
1280 
1289  void getSlownessUncertainty( const int& phase, const double& distance, double& uncert );
1290 
1297  void getSlownessUncertainty(double& slownessUncertainty);
1298 
1299  string getUncertaintyTable(const int& attribute, const int& phase);
1300 
1301  string getUncertaintyFileFormat(const int& attribute, const int& phase);
1302 
1324  void getZhaoParameters(double& Vm, double& Gm, double& H, double& C, double& Cm, int& udSign);
1325 
1342  void getPgLgComponents(double& tTotal,
1343  double& tTaup, double& tHeadwave,
1344  double& pTaup, double& pHeadwave,
1345  double& trTaup, double& trHeadwave);
1346 
1355  void static setCHMax(const double& chMax);
1356 
1365  void static getCHMax(double& chMax);
1366 
1376  void getAverageMantleVelocity(const int& type, double& velocity);
1377 
1389  void setAverageMantleVelocity(const int& type, const double& velocity);
1390 
1394  void getTessId(string& tessId);
1395 
1401  void getFractionActive(double& fractionActive);
1402 
1410  void static setMaxDistance(const double& maxDistance);
1411 
1417  void static getMaxDistance(double& maxDistance);
1418 
1426  void static setMaxDepth(const double& maxDepth);
1427 
1431  void static getMaxDepth(double& maxDepth);
1432 
1439  string getClassCount();
1440 
1445  const string& getModelPath() const;
1446 
1462  void getDistAz(const double& aLat, const double& aLon,
1463  const double& bLat, const double& bLon,
1464  double& distance, double& azimuth, const double& naValue);
1465 
1478  void movePoint(const double& aLat, const double& aLon,
1479  const double& distance, const double& azimuth,
1480  double& bLat, double& bLon);
1481 
1490  void getPiercePointSource(double& lat, double& lon, double& depth);
1491 
1500  void getPiercePointReceiver(double& lat, double& lon, double& depth);
1501 
1518  void getGreatCirclePoints(
1519  const double& aLat,
1520  const double& aLon,
1521  const double& bLat,
1522  const double& bLon,
1523  const int& npoints,
1524  double latitude[],
1525  double longitude[]);
1526 
1541  void getGreatCirclePointsOnCenters(
1542  const double& aLat,
1543  const double& aLon,
1544  const double& bLat,
1545  const double& bLon,
1546  const int& npoints,
1547  double latitude[],
1548  double longitude[]);
1549 
1555  void setDelDistance(const double& del_distance);
1556 
1562  void getDelDistance(double& del_distance);
1563 
1569  void setDelDepth(const double& del_depth);
1570 
1575  void getDelDepth(double& del_depth);
1576 
1580  void getRayParameter(double& ray_parameter);
1581 
1585  void getTurningRadius(double& turning_radius);
1586 
1600  void setPathIncrement(const double& pathIncrement);
1601 
1614  void getPathIncrement(double& pathIncrement);
1615 
1629 
1630  string getModelString() { return grid->toString(); }
1631 
1632 protected:
1633 
1638 
1643 
1647  static double CH_MAX;
1648 
1654  bool valid;
1655 
1659 
1660 private:
1661 
1662  // copies of the source and receiver locations that were
1663  // passed to createGreatCircle the last time it was called.
1664  // Used to compute distance in instances when createGreateCircle
1665  // fails. Units are radians.
1666  double srcLat, srcLon, srcDep, rcvLat, rcvLon, rcvDep;
1667 };
1668 
1669 //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1670 //
1671 // INLINE FUNCTIONS
1672 //
1673 //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1674 
1675 inline void SlbmInterface::createGreatCircle(
1676  const string& p,
1677  const double& sourceLat,
1678  const double& sourceLon,
1679  const double& sourceDepth,
1680  const double& receiverLat,
1681  const double& receiverLon,
1682  const double& receiverDepth)
1683 {
1684  int iphase = (p=="Pn" ? Pn : (p=="Sn" ? Sn : (p=="Pg" ? Pg : (p=="Lg" ? Lg : -1))));
1685  if (iphase == -1)
1686  {
1687  ostringstream os;
1688  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1689  os << endl << "ERROR in SlbmInterface::createGreatCircle" << endl
1690  << p << " is not a recognized phase. Must be one of Pn, Sn, Pg, Lg" << endl
1691  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1692  throw SLBMException(os.str(),112);
1693  }
1694 
1695  createGreatCircle(
1696  iphase,
1697  sourceLat,
1698  sourceLon,
1699  sourceDepth,
1700  receiverLat,
1701  receiverLon,
1702  receiverDepth
1703  );
1704 }
1705 
1706 inline GreatCircle* SlbmInterface::getGreatCircleObject()
1707 {
1708  return greatCircle;
1709 }
1710 
1711 inline void SlbmInterface::get_dtt_dlat(double& value)
1712 {
1713  if (!isValid())
1714  {
1715  value = NA_VALUE;
1716  ostringstream os;
1717  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1718  os << endl << "ERROR in SlbmInterface::get_dtt_dlat" << endl
1719  << "GreatCircle is invalid." << endl
1720  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1721  throw SLBMException(os.str(),113);
1722  }
1723  greatCircle->get_dtt_dlat(value);
1724 }
1725 
1726 inline void SlbmInterface::get_dtt_dlon(double& value)
1727 {
1728  if (!isValid())
1729  {
1730  value = NA_VALUE;
1731  ostringstream os;
1732  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1733  os << endl << "ERROR in SlbmInterface::get_dtt_dlon" << endl
1734  << "GreatCircle is invalid." << endl
1735  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1736  throw SLBMException(os.str(),113);
1737  }
1738  greatCircle->get_dtt_dlon(value);
1739 }
1740 
1741 inline void SlbmInterface::get_dtt_ddepth(double& value)
1742 {
1743  if (!isValid())
1744  {
1745  value = NA_VALUE;
1746  ostringstream os;
1747  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1748  os << endl << "ERROR in SlbmInterface::get_dtt_ddepth" << endl
1749  << "GreatCircle is invalid." << endl
1750  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1751  throw SLBMException(os.str(),113);
1752  }
1753  greatCircle->get_dtt_ddepth(value);
1754 }
1755 
1756 //inline void SlbmInterface::get_dtt_dlat_fast(double& value)
1757 //{
1758 // if (!isValid())
1759 // {
1760 // value = NA_VALUE;
1761 // ostringstream os;
1762 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1763 // os << endl << "ERROR in SlbmInterface::get_dtt_dlat" << endl
1764 // << "GreatCircle is invalid." << endl
1765 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1766 // throw SLBMException(os.str(),113);
1767 // }
1768 // greatCircle->get_dtt_dlat_fast(value);
1769 //}
1770 //
1771 //inline void SlbmInterface::get_dtt_dlon_fast(double& value)
1772 //{
1773 // if (!isValid())
1774 // {
1775 // value = NA_VALUE;
1776 // ostringstream os;
1777 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1778 // os << endl << "ERROR in SlbmInterface::get_dtt_dlon" << endl
1779 // << "GreatCircle is invalid." << endl
1780 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1781 // throw SLBMException(os.str(),113);
1782 // }
1783 // greatCircle->get_dtt_dlon_fast(value);
1784 //}
1785 //
1786 //inline void SlbmInterface::get_dsh_ddist(double& value)
1787 //{
1788 // if (!isValid())
1789 // {
1790 // value = NA_VALUE;
1791 // ostringstream os;
1792 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1793 // os << endl << "ERROR in SlbmInterface::get_dsh_ddist" << endl
1794 // << "GreatCircle is invalid." << endl
1795 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1796 // throw SLBMException(os.str(),113);
1797 // }
1798 // greatCircle->get_dsh_ddist(value);
1799 //}
1800 //
1801 //inline void SlbmInterface::get_dsh_dlat(double& value)
1802 //{
1803 // if (!isValid())
1804 // {
1805 // value = NA_VALUE;
1806 // ostringstream os;
1807 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1808 // os << endl << "ERROR in SlbmInterface::get_dsh_dlat" << endl
1809 // << "GreatCircle is invalid." << endl
1810 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1811 // throw SLBMException(os.str(),113);
1812 // }
1813 // greatCircle->get_dsh_dlat(value);
1814 //}
1815 //
1816 //inline void SlbmInterface::get_dsh_dlon(double& value)
1817 //{
1818 // if (!isValid())
1819 // {
1820 // value = NA_VALUE;
1821 // ostringstream os;
1822 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1823 // os << endl << "ERROR in SlbmInterface::get_dsh_dlon" << endl
1824 // << "GreatCircle is invalid." << endl
1825 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1826 // throw SLBMException(os.str(),113);
1827 // }
1828 // greatCircle->get_dsh_dlon(value);
1829 //}
1830 //
1831 //inline void SlbmInterface::get_dsh_ddepth(double& value)
1832 //{
1833 // if (!isValid())
1834 // {
1835 // value = NA_VALUE;
1836 // ostringstream os;
1837 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1838 // os << endl << "ERROR in SlbmInterface::get_dsh_ddepth" << endl
1839 // << "GreatCircle is invalid." << endl
1840 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1841 // throw SLBMException(os.str(),113);
1842 // }
1843 // greatCircle->get_dsh_ddepth(value);
1844 //}
1845 
1846 inline void SlbmInterface::getSlowness(double& value)
1847 {
1848  if (!isValid())
1849  {
1850  value = NA_VALUE;
1851  ostringstream os;
1852  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1853  os << endl << "ERROR in SlbmInterface::getSlowness" << endl
1854  << "GreatCircle is invalid." << endl
1855  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1856  throw SLBMException(os.str(),113);
1857  }
1858  greatCircle->get_dtt_ddist(value);
1859 }
1860 
1861 inline double SlbmInterface::getDistance()
1862 {
1863  if (greatCircle == NULL)
1864  {
1865  Location s(srcLat, srcLon, 0.);
1866  Location r(rcvLat, rcvLon, 0.);
1867  return s.distance(r);
1868  }
1869  else
1870  return greatCircle->getDistance();
1871 }
1872 
1873 inline void SlbmInterface::getSourceDistance(double& distance)
1874 {
1875  if (!isValid())
1876  {
1877  distance = NA_VALUE;
1878  ostringstream os;
1879  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1880  os << endl << "ERROR in SlbmInterface::getSourceDistance" << endl
1881  << "GreatCircle is invalid." << endl
1882  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1883  throw SLBMException(os.str(),113);
1884  }
1885  distance = greatCircle->getSourceDistance();
1886 }
1887 
1888 inline void SlbmInterface::getReceiverDistance(double& distance)
1889 {
1890  if (!isValid())
1891  {
1892  distance = NA_VALUE;
1893  ostringstream os;
1894  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1895  os << endl << "ERROR in SlbmInterface::getReceiverDistance" << endl
1896  << "GreatCircle is invalid." << endl
1897  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1898  throw SLBMException(os.str(),113);
1899  }
1900  distance = greatCircle->getReceiverDistance();
1901 }
1902 
1903 inline void SlbmInterface::getHeadwaveDistance(double& distance)
1904 {
1905  if (!isValid())
1906  {
1907  distance = NA_VALUE;
1908  ostringstream os;
1909  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1910  os << endl << "ERROR in SlbmInterface::getHeadwaveDistance" << endl
1911  << "GreatCircle is invalid." << endl
1912  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1913  throw SLBMException(os.str(),113);
1914  }
1915  distance = greatCircle->getHeadwaveDistance();
1916 }
1917 
1918 inline void SlbmInterface::getHeadwaveDistanceKm(double& distance)
1919 {
1920  if (!isValid())
1921  {
1922  distance = NA_VALUE;
1923  ostringstream os;
1924  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1925  os << endl << "ERROR in SlbmInterface::getHeadwaveDistanceKm" << endl
1926  << "GreatCircle is invalid." << endl
1927  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1928  throw SLBMException(os.str(),113);
1929  }
1930  distance = greatCircle->getHeadwaveDistanceKm();
1931 }
1932 
1933 inline void SlbmInterface::getTravelTime(double& tTotal)
1934 {
1935  if (!isValid())
1936  {
1937  tTotal = NA_VALUE;
1938  ostringstream os;
1939  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1940  os << endl << "ERROR in SlbmInterface::getTravelTime" << endl
1941  << "GreatCircle is invalid." << endl
1942  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1943  throw SLBMException(os.str(),113);
1944  }
1945  tTotal=greatCircle->getTravelTime();
1946 }
1947 
1948 inline void SlbmInterface::getTravelTimeComponents(
1949  double& tTotal, double& tSource, double& tReceiver,
1950  double& tHeadwave, double& tGradient)
1951 {
1952  if (!isValid())
1953  {
1954  tTotal = NA_VALUE;
1955  tSource = NA_VALUE;
1956  tReceiver = NA_VALUE;
1957  tHeadwave = NA_VALUE;
1958  tGradient = NA_VALUE;
1959  ostringstream os;
1960  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1961  os << endl << "ERROR in SlbmInterface::getTravelTimeComponents" << endl
1962  << "GreatCircle is invalid." << endl
1963  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1964  throw SLBMException(os.str(),113);
1965  }
1966  greatCircle->getTravelTime(tTotal, tSource, tReceiver,
1967  tHeadwave, tGradient);
1968 }
1969 
1970 inline void SlbmInterface::getWeights(int nodeId[], double weight[], int& nWeights)
1971 {
1972  if (!isValid())
1973  {
1974  nWeights=-1;
1975  ostringstream os;
1976  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
1977  os << endl << "ERROR in SlbmInterface::getWeights" << endl
1978  << "GreatCircle is invalid." << endl
1979  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
1980  throw SLBMException(os.str(),113);
1981  }
1982  greatCircle->getWeights(nodeId, weight, nWeights);
1983 }
1984 
1985 inline void SlbmInterface::getActiveNodeWeights(int nodeId[], double weight[], int& nWeights)
1986 {
1987  getWeights(nodeId, weight, nWeights);
1988  for (int i=0; i<nWeights; ++i)
1989  nodeId[i] = grid->getActiveNodeId(nodeId[i]);
1990 }
1991 
1992 inline void SlbmInterface::getWeights(
1993  vector<int>& nodeId, vector<double>& weight)
1994 {
1995  if (!isValid())
1996  {
1997  nodeId.clear();
1998  weight.clear();
1999  ostringstream os;
2000  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2001  os << endl << "ERROR in SlbmInterface::getWeights" << endl
2002  << "GreatCircle is invalid." << endl
2003  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2004  throw SLBMException(os.str(),113);
2005  }
2006  greatCircle->getWeights(nodeId, weight);
2007 }
2008 
2009 inline void SlbmInterface::getActiveNodeWeights(
2010  vector<int>& nodeId, vector<double>& weight)
2011 {
2012  getWeights(nodeId, weight);
2013  for (int i=0; i<(int)nodeId.size(); ++i)
2014  nodeId[i] = grid->getActiveNodeId(nodeId[i]);
2015 }
2016 
2017 inline void SlbmInterface::getWeightsSource(int nodeids[], double weights[], int& nWeights)
2018 {
2019  if (!isValid())
2020  {
2021  ostringstream os;
2022  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2023  os << endl << "ERROR in SlbmInterface::getWeightsSource" << endl
2024  << "GreatCircle is invalid." << endl
2025  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2026  throw SLBMException(os.str(),113);
2027  }
2028  greatCircle->getSourceProfile()->getWeights(nodeids, weights, nWeights);
2029 }
2030 
2031 inline void SlbmInterface::getActiveNodeWeightsSource(int nodeids[], double weights[], int& nWeights)
2032 {
2033  getWeightsSource(nodeids, weights, nWeights);
2034  for (int i=0; i<nWeights; ++i)
2035  nodeids[i] = grid->getActiveNodeId(nodeids[i]);
2036 }
2037 
2038 inline void SlbmInterface::getWeightsReceiver(int nodeids[], double weights[], int& nWeights)
2039 {
2040  if (!isValid())
2041  {
2042  ostringstream os;
2043  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2044  os << endl << "ERROR in SlbmInterface::getWeightsReceiver" << endl
2045  << "GreatCircle is invalid." << endl
2046  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2047  throw SLBMException(os.str(),113);
2048  }
2049  greatCircle->getReceiverProfile()->getWeights(nodeids, weights, nWeights);
2050 }
2051 
2052 inline void SlbmInterface::getActiveNodeWeightsReceiver(int nodeids[], double weights[], int& nWeights)
2053 {
2054  getWeightsReceiver(nodeids, weights, nWeights);
2055  for (int i=0; i<nWeights; ++i)
2056  nodeids[i] = grid->getActiveNodeId(nodeids[i]);
2057 }
2058 
2059 inline string SlbmInterface::toString(const int& verbosity)
2060 {
2061  if (!isValid())
2062  {
2063  ostringstream os;
2064  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2065  os << endl << "ERROR in SlbmInterface::toString" << endl
2066  << "GreatCircle is invalid." << endl
2067  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2068  throw SLBMException(os.str(),113);
2069  }
2070 
2071  ostringstream os;
2072  if (verbosity > 0)
2073  {
2074  os << endl
2075  << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << endl
2076  << "Great Circle " << endl << endl
2077  << greatCircle->toString(verbosity);
2078 
2079  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint)
2080  << setprecision(4);
2081  double slowness;
2082  getSlowness(slowness);
2083  os << endl << "Horizontal slowness = "
2084  << slowness << " sec/radian"
2085  << endl << endl;
2086  }
2087 
2088  return os.str();
2089 }
2090 
2091 inline void SlbmInterface::getGreatCircleData(
2092  string& phase,
2093  double& actual_path_increment,
2094  double sourceDepth[NLAYERS],
2095  double sourceVelocity[NLAYERS],
2096  double receiverDepth[NLAYERS],
2097  double receiverVelocity[NLAYERS],
2098  int& npoints,
2099  double headWaveVelocity[],
2100  double gradient[]
2101  )
2102 {
2103  if (!isValid())
2104  {
2105  phase = "";
2106  actual_path_increment = NA_VALUE;
2107  ostringstream os;
2108  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2109  os << endl << "ERROR in SlbmInterface::getGreatCircleData" << endl
2110  << "GreatCircle is invalid." << endl
2111  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2112  throw SLBMException(os.str(),113);
2113  }
2114 
2115  int p;
2116  greatCircle->getData(p, actual_path_increment,
2117  sourceDepth, sourceVelocity, receiverDepth, receiverVelocity,
2118  npoints, headWaveVelocity, gradient);
2119  phase = (p==Pn ? "Pn" : (p==Sn ? "Sn" : (p==Pg ? "Pg" : (p==Lg ? "Lg"
2120  : "unknown phase"))));
2121 }
2122 
2123 inline void SlbmInterface::getGreatCircleNodeInfo(
2124  int** neighbors,
2125  double** coefficients,
2126  const int& maxpoints,
2127  const int& maxnodes,
2128  int& npoints,
2129  int* nnodes
2130  )
2131 {
2132  if (!isValid())
2133  {
2134  ostringstream os;
2135  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2136  os << endl << "ERROR in SlbmInterface::getGreatCircleNodeInfo" << endl
2137  << "GreatCircle is invalid." << endl
2138  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2139  throw SLBMException(os.str(),113);
2140  }
2141 
2142  greatCircle->getNodeInfo(neighbors, coefficients, maxpoints, maxnodes, npoints, nnodes);
2143 }
2144 
2145 inline void SlbmInterface::getGreatCircleLocations(
2146  double lat[],
2147  double lon[],
2148  double depth[],
2149  int& npoints
2150  )
2151 {
2152  if (!isValid())
2153  {
2154  ostringstream os;
2155  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2156  os << endl << "ERROR in SlbmInterface::getGreatCircleData" << endl
2157  << "GreatCircle is invalid." << endl
2158  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2159  throw SLBMException(os.str(),113);
2160  }
2161 
2162  npoints = greatCircle->getNProfiles();
2163  Location loc;
2164  for (int i=0; i<greatCircle->getNProfiles(); i++)
2165  {
2166  // get the location of this layer profile
2167  greatCircle->getLayerProfileLocation(i, loc);
2168  lat[i] = loc.getLat();
2169  lon[i] = loc.getLon();
2170  depth[i] = loc.getDepth();
2171  }
2172 }
2173 
2174 
2175 inline void SlbmInterface::getInterpolatedPoint(
2176  const double& lat,
2177  const double& lon,
2178  int* nodeId,
2179  double* coefficients,
2180  int& nNeighbors,
2181  double depth[NLAYERS],
2182  double pvelocity[NLAYERS],
2183  double svelocity[NLAYERS],
2184  double& pgradient,
2185  double& sgradient
2186  )
2187 {
2188  Location location(lat, lon, 0.);
2189  QueryProfile* profile = grid->getQueryProfile(location);
2190  profile->getData(nodeId, coefficients, nNeighbors, depth, pvelocity, svelocity, pgradient, sgradient);
2191  delete profile;
2192 }
2193 
2194 inline void SlbmInterface::getInterpolatedTransect(
2195  double lat[],
2196  double lon[],
2197  const int& nLatLon,
2198  int** nodeId,
2199  double** coefficients,
2200  int* nNeighbors,
2201  double depth[][NLAYERS],
2202  double pvelocity[][NLAYERS],
2203  double svelocity[][NLAYERS],
2204  double pgradient[NLAYERS],
2205  double sgradient[NLAYERS],
2206  int& nInvalid
2207  )
2208 {
2209  nInvalid = 0;
2210  for (int i=0; i<nLatLon; i++)
2211  {
2212  try
2213  {
2214  getInterpolatedPoint(lat[i], lon[i], nodeId[i],
2215  coefficients[i], nNeighbors[i], depth[i], pvelocity[i], svelocity[i],
2216  pgradient[i], sgradient[i]);
2217  }
2218  catch (SLBMException ex)
2219  {
2220  for (int j=0; j<nNeighbors[i]; j++)
2221  {
2222  nodeId[i][j] = -1;
2223  coefficients[i][j] = NA_VALUE;
2224  }
2225 
2226  for (int j=0; j<NLAYERS; j++)
2227  {
2228  depth[i][j] = NA_VALUE;
2229  pvelocity[i][j] = NA_VALUE;
2230  svelocity[i][j] = NA_VALUE;
2231  }
2232 
2233  pgradient[i] = NA_VALUE;
2234  sgradient[i] = NA_VALUE;
2235 
2236  ++nInvalid;
2237  }
2238  }
2239 }
2240 
2241 inline void SlbmInterface::getNGridNodes(int& n)
2242 {
2243  if (!grid)
2244  {
2245  n = -1;
2246  ostringstream os;
2247  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2248  os << endl << "ERROR in SlbmInterface::getNGridNodes" << endl
2249  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2250  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2251  throw SLBMException(os.str(),114);
2252  }
2253  n = grid->getNNodes();
2254 }
2255 
2256 inline int SlbmInterface::getNGridNodes()
2257 {
2258  if (!grid)
2259  {
2260  ostringstream os;
2261  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2262  os << endl << "ERROR in SlbmInterface::getNGridNodes" << endl
2263  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2264  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2265  throw SLBMException(os.str(),114);
2266  }
2267  return grid->getNNodes();
2268 }
2269 
2270 inline void SlbmInterface::getNHeadWavePoints(int& nHeadWavePoints)
2271 {
2272  if (!greatCircle)
2273  {
2274  nHeadWavePoints = -1;
2275  ostringstream os;
2276  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2277  os << endl << "ERROR in SlbmInterface::getNHeadWavePoints" << endl
2278  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2279  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2280  throw SLBMException(os.str(),113);
2281  }
2282  nHeadWavePoints = greatCircle->getNProfiles();
2283 }
2284 
2285 //inline void SlbmInterface::getMaxInterpotionNodes(int& maxNodes)
2286 //{
2287 // if (!greatCircle)
2288 // {
2289 // nHeadWavePoints = -1;
2290 // ostringstream os;
2291 // os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2292 // os << endl << "ERROR in SlbmInterface::getMaxInterpotionNodes" << endl
2293 // << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2294 // << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2295 // throw SLBMException(os.str(),113);
2296 // }
2297 //
2298 // maxNodes = 1;
2299 // for (int i=0; i<greatCircle->getNProfiles(); i++)
2300 // if (greatCircle->profiles[i]->getNodes().size() > maxNodes)
2301 // maxNodes = greatCircle->profiles[i]->getNodes().size();
2302 //}
2303 
2304 inline void SlbmInterface::getSlownessUncertainty( const int& phase, const double& distance, double& uncert )
2305 {
2306  if (!grid)
2307  {
2308  ostringstream os;
2309  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2310  os << endl << "ERROR in SlbmInterface::getSlownessUncertainty" << endl
2311  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2312  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2313  throw SLBMException(os.str(),114);
2314  }
2315 
2316  if (grid->getUncertainty()[phase][SH] == NULL)
2317  {
2318  ostringstream os;
2319  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2320  os << endl << "ERROR in SlbmInterface::getSlownessUncertainty" << endl
2321  << "Uncertainty object is invalid.." << endl
2322  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2323  throw SLBMException(os.str(),602);
2324  }
2325  else
2326  {
2327  // retrieve slowness uncertainty in sec/radian
2328  uncert = grid->getUncertainty()[phase][SH]->getUncertainty( distance );
2329  }
2330 }
2331 
2332 inline void SlbmInterface::getSlownessUncertainty(double& slownessUncertainty)
2333 {
2334  if (!grid)
2335  {
2336  ostringstream os;
2337  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2338  os << endl << "ERROR in SlbmInterface::getSlownessUncertainty" << endl
2339  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2340  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2341  throw SLBMException(os.str(),114);
2342  }
2343 
2344  getSlownessUncertainty(greatCircle->getPhase(),
2345  greatCircle->getDistance(), slownessUncertainty);
2346 }
2347 
2348 inline string SlbmInterface::getUncertaintyTable(const int& phase, const int& attribute)
2349 {
2350  if (!grid)
2351  {
2352  ostringstream os;
2353  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2354  os << endl << "ERROR in SlbmInterface::getUncertaintyTable" << endl
2355  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2356  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2357  throw SLBMException(os.str(),114);
2358  }
2359 
2360  if (grid->getUncertainty()[phase][attribute] == NULL)
2361  {
2362  ostringstream os;
2363  os << "No uncertainty information is available for phase " <<
2364  Uncertainty::getPhase(phase) << " attribute " <<
2365  Uncertainty::getAttribute(attribute) << endl;
2366  return os.str();
2367  }
2368  return grid->getUncertainty()[phase][attribute]->toStringTable();
2369 }
2370 
2371 inline string SlbmInterface::getUncertaintyFileFormat(const int& phase, const int& attribute)
2372 {
2373  if (!grid)
2374  {
2375  ostringstream os;
2376  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2377  os << endl << "ERROR in SlbmInterface::getUncertaintyFileFormat" << endl
2378  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2379  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2380  throw SLBMException(os.str(),114);
2381  }
2382 
2383  if (grid->getUncertainty()[phase][attribute] == NULL)
2384  {
2385  ostringstream os;
2386  os << "No uncertainty information is available for phase " << phase << " attribute " << attribute << endl;
2387  return os.str();
2388  }
2389  return grid->getUncertainty()[phase][attribute]->toStringFile();
2390 }
2391 
2392 inline void SlbmInterface::getZhaoParameters(double& Vm, double& Gm, double& H, double& C, double& Cm, int& udSign)
2393 {
2394  if (!greatCircle)
2395  {
2396  ostringstream os;
2397  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2398  os << endl << "ERROR in SlbmInterface::getZhaoParameters" << endl
2399  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2400  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2401  throw SLBMException(os.str(),113);
2402  }
2403  greatCircle->getZhaoParameters(Vm, Gm, H, C, Cm, udSign);
2404 }
2405 
2406 inline void SlbmInterface::getPgLgComponents(double& tTotal,
2407  double& tTaup, double& tHeadwave,
2408  double& pTaup, double& pHeadwave,
2409  double& trTaup, double& trHeadwave)
2410 {
2411  if (!greatCircle)
2412  {
2413  ostringstream os;
2414  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2415  os << endl << "ERROR in SlbmInterface::getPgLgComponents" << endl
2416  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2417  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2418  throw SLBMException(os.str(),113);
2419  }
2420  greatCircle->getPgLgComponents(tTotal, tTaup, tHeadwave,
2421  pTaup, pHeadwave, trTaup, trHeadwave);
2422 }
2423 
2424 
2425 inline void SlbmInterface::getNodeNeighbors(const int& nid, int neighbors[], int& nNeighbors)
2426 {
2427  if (!grid)
2428  {
2429  ostringstream os;
2430  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2431  os << endl << "ERROR in SlbmInterface::getNodeNeighbors" << endl
2432  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2433  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2434  throw SLBMException(os.str(),114);
2435  }
2436  grid->getNodeNeighbors(nid, neighbors, nNeighbors);
2437 }
2438 
2439 inline void SlbmInterface::getActiveNodeNeighbors(const int& nid, int neighbors[], int& nNeighbors)
2440 {
2441  if (!grid)
2442  {
2443  ostringstream os;
2444  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2445  os << endl << "ERROR in SlbmInterface::getActiveNodeNeighbors" << endl
2446  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2447  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2448  throw SLBMException(os.str(),114);
2449  }
2450  grid->getActiveNodeNeighbors(nid, neighbors, nNeighbors);
2451 }
2452 
2453 inline void SlbmInterface::getNodeNeighbors(const int& nid, vector<int>& neighbors)
2454 {
2455  if (!grid)
2456  {
2457  ostringstream os;
2458  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2459  os << endl << "ERROR in SlbmInterface::getNodeNeighbors" << endl
2460  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2461  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2462  throw SLBMException(os.str(),114);
2463  }
2464  grid->getNodeNeighbors(nid, neighbors);
2465 }
2466 
2467 inline void SlbmInterface::getActiveNodeNeighbors(const int& nid, vector<int>& neighbors)
2468 {
2469  if (!grid)
2470  {
2471  ostringstream os;
2472  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2473  os << endl << "ERROR in SlbmInterface::getActiveNodeNeighbors" << endl
2474  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2475  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2476  throw SLBMException(os.str(),114);
2477  }
2478  grid->getActiveNodeNeighbors(nid, neighbors);
2479 }
2480 
2481 inline void SlbmInterface::getNodeNeighborInfo(const int& nid, int neighbors[],
2482  double distance[], double azimuth[], int& nNeighbors)
2483 {
2484  if (!grid)
2485  {
2486  ostringstream os;
2487  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2488  os << endl << "ERROR in SlbmInterface::getNodeNeighborInfo" << endl
2489  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2490  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2491  throw SLBMException(os.str(),114);
2492  }
2493  grid->getNodeNeighborInfo(nid, neighbors, distance, azimuth, nNeighbors);
2494 }
2495 
2496 inline void SlbmInterface::getActiveNodeNeighborInfo(const int& nid, int neighbors[],
2497  double distance[], double azimuth[], int& nNeighbors)
2498 {
2499  if (!grid)
2500  {
2501  ostringstream os;
2502  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2503  os << endl << "ERROR in SlbmInterface::getActiveNodeNeighborInfo" << endl
2504  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2505  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2506  throw SLBMException(os.str(),114);
2507  }
2508  grid->getActiveNodeNeighborInfo(nid, neighbors, distance, azimuth, nNeighbors);
2509 }
2510 
2511 inline void SlbmInterface::getActiveNodeNeighborInfo(const int& nid,
2512  vector<int>& neighbors, vector<double>& distance, vector<double>& azimuth)
2513 {
2514  if (!grid)
2515  {
2516  ostringstream os;
2517  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2518  os << endl << "ERROR in SlbmInterface::getActiveNodeNeighborInfo" << endl
2519  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2520  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2521  throw SLBMException(os.str(),114);
2522  }
2523  grid->getActiveNodeNeighborInfo(nid, neighbors, distance, azimuth);
2524 }
2525 
2526 inline void SlbmInterface::getNodeSeparation(const int& node1, const int& node2, double& distance)
2527 {
2528  if (!grid)
2529  {
2530  ostringstream os;
2531  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2532  os << endl << "ERROR in SlbmInterface::getNodeSeparation" << endl
2533  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2534  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2535  throw SLBMException(os.str(),114);
2536  }
2537  grid->getNodeSeparation(node1, node2, distance);
2538 }
2539 
2540 inline void SlbmInterface::getNodeAzimuth(const int& node1, const int& node2, double& azimuth)
2541 {
2542  if (!grid)
2543  {
2544  ostringstream os;
2545  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2546  os << endl << "ERROR in SlbmInterface::getNodeAzimuth" << endl
2547  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2548  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2549  throw SLBMException(os.str(),114);
2550  }
2551  grid->getNodeAzimuth(node1, node2, azimuth);
2552 }
2553 
2554 inline void SlbmInterface::initializeActiveNodes(const double& latmin, const double& lonmin,
2555  const double& latmax, const double& lonmax)
2556 {
2557  if (!grid)
2558  {
2559  ostringstream os;
2560  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2561  os << endl << "ERROR in SlbmInterface::initializeActiveNodes" << endl
2562  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2563  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2564  throw SLBMException(os.str(),114);
2565  }
2566  grid->initializeActiveNodes(latmin, lonmin, latmax, lonmax);
2567 }
2568 
2569 inline void SlbmInterface::initializeActiveNodes(double* lat, double* lon, const int& npoints, const bool& degrees)
2570 {
2571  vector<double*> v;
2572  v.reserve(npoints);
2573  if (degrees)
2574  for (int i=0; i<npoints; ++i)
2575  v.push_back(GeoTessUtils::getVectorDegrees(lat[i], lon[i]));
2576  else
2577  for (int i=0; i<npoints; ++i)
2578  v.push_back(GeoTessUtils::getVector(lat[i], lon[i]));
2579  initializeActiveNodes(v);
2580 }
2581 
2582 inline void SlbmInterface::initializeActiveNodes(vector<double*>& unitVectors)
2583 {
2584  initializeActiveNodes(new GeoTessPolygon(unitVectors));
2585 }
2586 
2587 inline void SlbmInterface::initializeActiveNodes(GeoTessPolygon* polygon)
2588 {
2589  grid->initializeActiveNodes(polygon);
2590 }
2591 
2592 inline void SlbmInterface::initializeActiveNodes(const string& fileName)
2593 {
2594  grid->initializeActiveNodes(new GeoTessPolygon(fileName));
2595 }
2596 
2597 inline void SlbmInterface::clearActiveNodes()
2598 {
2599  grid->clearActiveNodes();
2600 }
2601 
2602 
2603 inline int SlbmInterface::getNActiveNodes()
2604 {
2605  if (!grid)
2606  {
2607  ostringstream os;
2608  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2609  os << endl << "ERROR in SlbmInterface::nextActiveNode" << endl
2610  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2611  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2612  throw SLBMException(os.str(),114);
2613  }
2614  return grid->getNActiveNodes();
2615 }
2616 
2617 inline int SlbmInterface::getGridNodeId(int activeNodeId)
2618 {
2619  if (!grid)
2620  {
2621  ostringstream os;
2622  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2623  os << endl << "ERROR in SlbmInterface::getGridNodeId" << endl
2624  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2625  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2626  throw SLBMException(os.str(),114);
2627  }
2628  return grid->getGridNodeId(activeNodeId);
2629 }
2630 
2631 inline int SlbmInterface::getActiveNodeId(int gridNodeId)
2632 {
2633  if (!grid)
2634  {
2635  ostringstream os;
2636  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2637  os << endl << "ERROR in SlbmInterface::getActiveNodeId" << endl
2638  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2639  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2640  throw SLBMException(os.str(),114);
2641  }
2642  return grid->getActiveNodeId(gridNodeId);
2643 }
2644 
2645 inline void SlbmInterface::getNodeHitCount(const int& nodeId, int& hitCount)
2646 {
2647  if (!grid)
2648  {
2649  ostringstream os;
2650  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2651  os << endl << "ERROR in SlbmInterface::getNodeHitCount" << endl
2652  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2653  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2654  throw SLBMException(os.str(),114);
2655  }
2656  grid->getNodeHitCount(nodeId, hitCount);
2657 }
2658 
2659 inline void SlbmInterface::clearNodeHitCount()
2660 {
2661  if (!grid)
2662  {
2663  ostringstream os;
2664  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2665  os << endl << "ERROR in SlbmInterface::clearNodeHitCount" << endl
2666  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2667  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2668  throw SLBMException(os.str(),114);
2669  }
2670  grid->clearNodeHitCount();
2671 }
2672 
2673 inline void SlbmInterface::getAverageMantleVelocity(const int& type, double& velocity)
2674 {
2675  if (!grid)
2676  {
2677  ostringstream os;
2678  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2679  os << endl << "ERROR in SlbmInterface::setAverageMantleVelocity" << endl
2680  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2681  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2682  throw SLBMException(os.str(),114);
2683  }
2684  velocity = grid->getAverageMantleVelocity(type);
2685 }
2686 
2687 inline void SlbmInterface::setAverageMantleVelocity(const int& type,
2688  const double& velocity)
2689 {
2690  if (!grid)
2691  {
2692  ostringstream os;
2693  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2694  os << endl << "ERROR in SlbmInterface::setAverageMantleVelocity" << endl
2695  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2696  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2697  throw SLBMException(os.str(),114);
2698  }
2699  grid->setAverageMantleVelocity(type, velocity);
2700 }
2701 
2702 inline void SlbmInterface::getTessId(string& tessId)
2703 {
2704  if (!grid)
2705  {
2706  ostringstream os;
2707  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2708  os << endl << "ERROR in SlbmInterface::getTessId" << endl
2709  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2710  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2711  throw SLBMException(os.str(),114);
2712  }
2713  tessId = grid->getTessId();
2714 }
2715 
2716 inline void SlbmInterface::getFractionActive(double& fractionActive)
2717 {
2718  if (!isValid())
2719  {
2720  fractionActive = NA_VALUE;
2721  ostringstream os;
2722  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2723  os << endl << "ERROR in SlbmInterface::getFractionActive" << endl
2724  << "GreatCircle is invalid." << endl
2725  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2726  throw SLBMException(os.str(),113);
2727  }
2728  fractionActive = greatCircle->getFractionActive();
2729 }
2730 
2731 inline string SlbmInterface::getClassCount()
2732 {
2733  ostringstream os;
2734  os << "Class counts:" << endl;
2735  os << "GreatCircle = " << GreatCircle::getClassCount() << endl;
2736  os << "GridProfile = " << GridProfile::getClassCount() << endl;
2737  os << "GeoStack = " << GeoStack::getClassCount() << endl;
2738  os << "InterpolatedProfile = " << InterpolatedProfile::getClassCount() << endl;
2739  os << "CrustalProfile = " << CrustalProfile::getClassCount() << endl;
2740  os << "LayerProfile = " << LayerProfile::getClassCount() << endl;
2741  os << "QueryProfile = " << QueryProfile::getClassCount() << endl;
2742  os << "Location = " << Location::getClassCount() << endl;
2743  return os.str();
2744 }
2745 
2746 inline const string& SlbmInterface::getModelPath() const { return grid->getModelPath(); }
2747 
2748 inline void SlbmInterface::getPiercePointSource(double& lat, double& lon, double& depth)
2749 {
2750  if (!isValid())
2751  {
2752  lat = lon = NA_VALUE;
2753  ostringstream os;
2754  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2755  os << endl << "ERROR in SlbmInterface::getPiercePointSource" << endl
2756  << "GreatCircle is invalid." << endl
2757  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2758  throw SLBMException(os.str(),113);
2759  }
2760 
2761 
2762  if (greatCircle->getPhase() == Pg || greatCircle->getPhase() == Lg)
2763  {
2764  lat = lon = NA_VALUE;
2765  ostringstream os;
2766  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2767  os << endl << "ERROR in SlbmInterface::getPiercePointSource" << endl
2768  << "Cannot compute moho pierce points for " << getPhase() << endl
2769  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2770  throw SLBMException(os.str(),113);
2771  }
2772 
2773  if (!greatCircle->getSourceProfile()->isInCrust())
2774  {
2775  lat = lon = NA_VALUE;
2776  ostringstream os;
2777  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2778  os << endl << "ERROR in SlbmInterface::getPiercePointSource" << endl
2779  << "Cannot compute moho pierce point for source in the mantle." << endl
2780  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2781  throw SLBMException(os.str(),113);
2782  }
2783 
2784  Location loc;
2785  greatCircle->getGreatCircleLocation(greatCircle->getSourceDistance(), loc);
2786 
2787  lat = loc.getLat();
2788  lon = loc.getLon();
2789 
2790  QueryProfile* profile = NULL;
2791  profile = grid->getQueryProfile(loc);
2792  depth = profile->getDepth()[MANTLE];
2793  delete profile;
2794 }
2795 
2796 inline void SlbmInterface::getPiercePointReceiver(double& lat, double& lon, double& depth)
2797 {
2798  if (!isValid())
2799  {
2800  lat = lon = NA_VALUE;
2801  ostringstream os;
2802  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2803  os << endl << "ERROR in SlbmInterface::getPiercePointReceiver" << endl
2804  << "GreatCircle is invalid." << endl
2805  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2806  throw SLBMException(os.str(),113);
2807  }
2808 
2809 
2810  if (greatCircle->getPhase() == Pg || greatCircle->getPhase() == Lg)
2811  {
2812  lat = lon = NA_VALUE;
2813  ostringstream os;
2814  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2815  os << endl << "ERROR in SlbmInterface::getPiercePointReceiver" << endl
2816  << "Cannot compute moho pierce points for " << getPhase() << endl
2817  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2818  throw SLBMException(os.str(),113);
2819  }
2820 
2821  if (!greatCircle->getReceiverProfile()->isInCrust())
2822  {
2823  lat = lon = NA_VALUE;
2824  ostringstream os;
2825  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2826  os << endl << "ERROR in SlbmInterface::getPiercePointReceiver" << endl
2827  << "Cannot compute moho pierce point for receiver in the mantle." << endl
2828  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2829  throw SLBMException(os.str(),113);
2830  }
2831 
2832  Location loc;
2833  greatCircle->getGreatCircleLocation(
2834  greatCircle->getDistance()-greatCircle->getReceiverDistance(), loc);
2835 
2836  lat = loc.getLat();
2837  lon = loc.getLon();
2838 
2839  QueryProfile* profile = NULL;
2840  profile = grid->getQueryProfile(loc);
2841  depth = profile->getDepth()[MANTLE];
2842  delete profile;
2843 }
2844 
2845 inline void SlbmInterface::getGreatCirclePoints(
2846  const double& aLat,
2847  const double& aLon,
2848  const double& bLat,
2849  const double& bLon,
2850  const int& npoints,
2851  double latitude[],
2852  double longitude[])
2853 {
2854  Location a(aLat, aLon, 0.);
2855  Location b(bLat, bLon, 0.);
2856  double dx = a.distance(b)/(npoints-1);
2857 
2858  double moveDirection[3];
2859  a.vectorTripleProduct(b, moveDirection);
2860  for (int i=0; i<npoints; i++)
2861  {
2862  a.move(moveDirection, i*dx, b);
2863  latitude[i] = b.getLat();
2864  longitude[i] = b.getLon();
2865  }
2866 }
2867 
2868 inline void SlbmInterface::getGreatCirclePointsOnCenters(
2869  const double& aLat,
2870  const double& aLon,
2871  const double& bLat,
2872  const double& bLon,
2873  const int& npoints,
2874  double latitude[],
2875  double longitude[])
2876 {
2877  Location a(aLat, aLon, 0.);
2878  Location b(bLat, bLon, 0.);
2879  double dx = a.distance(b)/npoints;
2880 
2881  double moveDirection[3];
2882  a.vectorTripleProduct(b, moveDirection);
2883  for (int i=0; i<npoints; i++)
2884  {
2885  a.move(moveDirection, (i+0.5)*dx, b);
2886  latitude[i] = b.getLat();
2887  longitude[i] = b.getLon();
2888  }
2889 }
2890 
2891 inline void SlbmInterface::getDistAz(const double& aLat, const double& aLon,
2892  const double& bLat, const double& bLon,
2893  double& distance, double& azimuth, const double& naValue)
2894 {
2895  Location ptA(aLat, aLon);
2896  Location ptB(bLat, bLon);
2897  distance = ptA.distance(ptB);
2898  azimuth = ptA.azimuth(ptB, naValue);
2899 }
2900 
2901 inline void SlbmInterface::movePoint(const double& aLat, const double& aLon,
2902  const double& distance, const double& azimuth,
2903  double& bLat, double& bLon)
2904 {
2905  Location ptA(aLat, aLon);
2906  Location ptB;
2907  ptA.move(azimuth, distance, ptB);
2908  bLat = ptB.getLat();
2909  bLon = ptB.getLon();
2910 }
2911 
2912 inline void SlbmInterface::setInterpolatorType(const string& interpolatorType)
2913 {
2914  if (!grid)
2915  {
2916  ostringstream os;
2917  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2918  os << endl << "ERROR in SlbmInterface::setInterpolatorType" << endl
2919  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2920  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2921  throw SLBMException(os.str(),114);
2922  }
2923  grid->setInterpolatorType(interpolatorType);
2924 }
2925 
2926 inline string SlbmInterface::getInterpolatorType()
2927 {
2928  if (!grid)
2929  {
2930  ostringstream os;
2931  os << setiosflags(ios::fixed) << setiosflags(ios::showpoint) << setprecision(9);
2932  os << endl << "ERROR in SlbmInterface::getInterpolatorType" << endl
2933  << "Grid is invalid. Has the earth model been loaded with call to loadVelocityModel()?" << endl
2934  << "Version " << SlbmVersion << " File " << __FILE__ << " line " << __LINE__ << endl << endl;
2935  throw SLBMException(os.str(),114);
2936  }
2937  return grid->getInterpolatorType();
2938 }
2939 
2940 
2941 
2942 } // end slbm namespace
2943 
2944 #endif // SlbmInterface_H
slbm::SlbmInterface::CH_MAX
static double CH_MAX
Definition: SlbmInterface.h:1647
slbm::SlbmInterface::getDelDistance
void getDelDistance(double &del_distance)
Retrieve the value of step change in distance used to compute horizontal derivatives (radians)
slbm::QueryProfile
An InterpolatedProfile object that also has information about the the P and S wave velocity as a func...
Definition: QueryProfile.h:68
slbm::SlbmInterface::getTravelTimeUncertainty1D
void getTravelTimeUncertainty1D(double &travelTimeUncertainty)
Retrieve travel time uncertainty in sec using the phase and distance specified in last call to getGre...
GreatCircle.h
slbm::QueryProfile::getDepth
double * getDepth()
Retrieve the depth of the top of the k'th interval, in km.
Definition: QueryProfile.h:129
slbm::SlbmInterface::setGridData
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.
slbm::SlbmInterface::getTravelTimeUncertainty
void getTravelTimeUncertainty(const int &phase, const double &distance, double &uncert)
Retrieve the travel time uncertainty in sec for specified phase, distance (in radians).
slbm::SlbmInterface::setPathIncrement
void setPathIncrement(const double &pathIncrement)
Set the desired spacing of great circle nodes along the head wave interface, in radians.
slbm::Location
The Location Class manages a single point in/on the Earth, which is described by the GRS80 ellipsoid.
Definition: Location.h:78
slbm::SlbmInterface::setCHMax
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...
slbm::SlbmInterface::getGridData
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...
SLBMException.h
slbm::SlbmInterface::valid
bool valid
true if the current GreatCirlce object has been instantiated and is ready to be interrogated.
Definition: SlbmInterface.h:1654
slbm::SlbmInterface::createGreatCircle
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.
slbm::SlbmInterface::clearGreatCircles
void clearGreatCircles()
GridGeoTess.h
UncertaintyPathDep.h
slbm::Location::vectorTripleProduct
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
slbm::SlbmInterface::getBufferSize
int getBufferSize() const
Returns the size of a DataBuffer object required to store this SLBMInterface objects model data.
GreatCircleFactory.h
slbm::SlbmInterface::get_dtt_dlat
double get_dtt_dlat()
Definition: SlbmInterface.h:450
slbm::SlbmInterface::getActiveNodeData
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...
slbm::SlbmInterface::setMaxDistance
static void setMaxDistance(const double &maxDistance)
Set the maximum source-receiver separation for Pn/Sn phase, in radians.
slbm::QueryProfile::getData
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
slbm::Location::getLat
double getLat() const
Retrieve the geographic latitude of this Location, radians.
Definition: Location.h:542
slbm::SlbmInterface::isValid
bool isValid()
Returns true if the current GreatCirlce object has been instantiated and is ready to be interrogated.
Definition: SlbmInterface.h:336
slbm::Location::getLon
double getLon() const
Retrieve the longitude of this Location. Value will be between -PI and PI radians.
Definition: Location.h:548
Uncertainty.h
slbm::SlbmInterface::setActiveNodeData
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...
slbm::SlbmInterface::getModelString
string getModelString()
Definition: SlbmInterface.h:1630
slbm::SlbmInterface::getRayParameter
void getRayParameter(double &ray_parameter)
Retrieve the ray parameter.
slbm::SlbmInterface::get_dtt_dlon
double get_dtt_dlon()
Definition: SlbmInterface.h:459
slbm::SlbmInterface::getPathIncrement
void getPathIncrement(double &pathIncrement)
Retrieve the current value of the spacing of great circle nodes along the head wave interface,...
SLBM_EXP
#define SLBM_EXP
Definition: SLBMGlobals.h:181
slbm::SlbmInterface::getDistance
void getDistance(double &distance)
Retrieve the source-receiver separation, in radians.
Definition: SlbmInterface.h:348
slbm::SlbmInterface::getCHMax
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 ...
slbm::SlbmInterface::getMaxDistance
static void getMaxDistance(double &maxDistance)
Retrieve the current value for the maximum source-receiver separation, in radians.
slbm::SlbmInterface::setDelDepth
void setDelDepth(const double &del_depth)
Change the value of step change in depth used to compute depth derivatives (km)
slbm::SlbmInterface::setMaxDepth
static void setMaxDepth(const double &maxDepth)
Set the maximum source depth for Pn/Sn phase, in km.
slbm::Grid
A 2 dimensional, horizontal grid of GirdProfile objects.
Definition: Grid.h:90
slbm::SlbmInterface::~SlbmInterface
virtual ~SlbmInterface()
Destructor.
slbm::SlbmInterface::getGridObject
Grid * getGridObject()
Retrieve a pointer to the Grid object.
Definition: SlbmInterface.h:1251
slbm::SlbmInterface::getTurningRadius
void getTurningRadius(double &turning_radius)
Retrieve the turning radius of the ray.
slbm::SlbmInterface::getVersion
string getVersion()
Retrieve the SLBM Version number.
Definition: SlbmInterface.h:111
slbm::SlbmInterface::clear
void clear()
Delete the current GreatCircle object from memory and clear the pool of stored CrustalProfile objects...
slbm::SlbmInterface::get_dtt_ddist
void get_dtt_ddist(double &slowness)
Retrieve the horizontal slowness, i.e., the derivative of travel time wrt to receiver-source distance...
Definition: SlbmInterface.h:440
slbm::SlbmInterface::greatCircle
GreatCircle * greatCircle
The most recently requested GreatCircle object.
Definition: SlbmInterface.h:1642
slbm::SlbmInterface::getTravelTimeUncertainty
void getTravelTimeUncertainty(double &travelTimeUncertainty, bool calcRandomError=false)
Retrieve travel time uncertainty in sec using the phase and distance specified in last call to getGre...
slbm::SlbmInterface::grid
Grid * grid
The Grid object that stores the velocity model.
Definition: SlbmInterface.h:1637
slbm::Location::azimuth
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
slbm::SlbmInterface
The primary interface to the SLBM library, providing access to all supported functionality.
Definition: SlbmInterface.h:77
slbm::SlbmInterface::getMaxDepth
static void getMaxDepth(double &maxDepth)
Retrieve the current value for the maximum source depth, in km.
slbm::SlbmInterface::setDelDistance
void setDelDistance(const double &del_distance)
Change the value of step change in distance used to compute horizontal derivatives(in radians).
slbm::SLBMException
An Exception class for Grid and related objects.
Definition: SLBMException.h:55
slbm::Location::move
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
slbm::SlbmInterface::get_dtt_ddist
double get_dtt_ddist()
Definition: SlbmInterface.h:441
slbm::SlbmInterface::SlbmInterface
SlbmInterface()
Default constructor. Instantiates an SlbmInterface object based on an ellipsoidal earth.
slbm::SlbmInterface::SlbmInterface
SlbmInterface(const double &earthRadius)
Parameterized constructor. Instantiates an SlbmInterface object that is only partly based on an ellip...
slbm::SlbmInterface::getPathIncrement
double getPathIncrement()
Retrieve the current value of the spacing of great circle nodes along the head wave interface,...
slbm::Location::distance
double distance(const Location &other) const
Definition: Location.h:583
slbm::GreatCircle
The GreatCircle class manages information related to a great circle path between two Locations on the...
Definition: GreatCircle.h:115
slbm::SlbmInterface::loadVelocityModel
void loadVelocityModel(const string &modelPath)
Load the velocity model into memory from the specified file or directory. This method automatically d...
slbm::Location::getDepth
double getDepth() const
Definition: Location.h:573
slbm
Definition: CrustalProfile.h:59
slbm::SlbmInterface::get_dtt_ddepth
double get_dtt_ddepth()
Definition: SlbmInterface.h:468
util::DataBuffer
A byte array container used to hold binary data in the same manner as disk based file system.
Definition: DataBuffer.h:81
slbm::SlbmInterface::getDelDepth
void getDelDepth(double &del_depth)
Retrieve the value of step change in depth used to compute depth derivatives (km)
slbm::SlbmInterface::saveVelocityModel
void saveVelocityModel(const string &modelFileName, const int &format=4)
Save the velocity model currently in memory to the specified file.
slbm::SlbmInterface::getPhase
string getPhase()
Retrieve the phase specified in last call to createGreatCircle().
Definition: SlbmInterface.h:341