RSTT  3.1.0
Regional Seismic Travel Time
slbm_C_shell.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 
118 
119 
120 #ifndef SLBM_C_SHELL_H
121 #define SLBM_C_SHELL_H
122 
123 #if defined(_WIN32) || defined(WIN32)
124 
125  #ifndef SLBM_C_EXPORT
126  #define SLBM_LIB __declspec(dllimport)
127  #else
128  #define SLBM_LIB __declspec(dllexport)
129  #endif
130 
131 #else
132 
133  #define SLBM_LIB
134 
135 #endif
136 
137 #ifdef __cplusplus
138 extern "C"
139 {
140 #endif
141 
142 //=================================================================================================
143 // ****NOTICE****
144 //
145 // The "C" shell implementation of the SLBM package has changed due to a package reorganization.
146 // The file "SLBMGlobals.h" has now been wrapped into a namespace, it is no longer possible to #include
147 // "SLBMGlobals.h" and compile in "C". The solution, so as not to require any changes to code already
148 // developed, is to hard-code the constants from "SLBMGlobals.h" into this header file. While it does
149 // satisfy that no changes are required to existing code, it is a liability if any of those constants
150 // do change. In the event (unlikely) the constants do change, the change will need to be made in two
151 // places.
152 //
153 
154 //static const char* SlbmVersion = "3.1.0";
155 
159 static const int PWAVE = 0;
160 
164 static const int SWAVE = 1;
165 
169 static const int Pn = 0;
170 
174 static const int Sn = 1;
175 
179 static const int Pg = 2;
180 
184 static const int Lg = 3;
185 
186 #define WATER 0
187 #define SEDIMENT1 1
188 #define SEDIMENT2 2
189 #define SEDIMENT3 3
190 #define UPPER_CRUST 4
191 #define MIDDLE_CRUST_N 5
192 #define MIDDLE_CRUST_G 6
193 #define LOWER_CRUST 7
194 #define MANTLE 8
195 #define NLAYERS 9
196 // MANTLE has to be last layer and must == NLAYERS-1
197 
198 // Parameter TOP_LAYER controls behavior when sources/receivers
199 // are in layer WATER.
200 // If TOP_LAYER = 0, and phase = Pn/Pg then rays that travel
201 // through the WATER layer travel at the velocity of water.
202 // If TOP_LAYER = 0, and phase = Sn/Lg then rays that travel
203 // through the WATER layer cause an exception to be thrown.
204 // If TOP_LAYER = 1 and the depth of a source or receiver is
205 // less than the depth of the top of layer SEDIMENT1, then
206 // the code will search downward through the layers, starting
207 // at layer SEDIMENT1, until if finds a layer with non-zero
208 // thickness. It will then extend the top of that layer up
209 // to the depth of the source/receiver. This essentially
210 // assumes that all sources/receivers in the ocean are actually
211 // located on islands that are too small to have been explicitly
212 // defined in the model.
213 //
214 // Bottom line is: to include water set TOP_LAYER = 0.
215 // To ignore water, set TOP_LAYER = 1;
216 static const int TOP_LAYER = 1;
217 
224 SLBM_LIB int slbm_shell_getVersion ( char* str );
225 
236 
242 
249 
257 
275 SLBM_LIB int slbm_shell_loadVelocityModelBinary ( const char* modelDirectory );
276 
294 SLBM_LIB int slbm_shell_specifyOutputDirectory ( const char* directoryName );
295 
307 
320 SLBM_LIB int slbm_shell_loadVelocityModel ( const char* modelPath );
321 
365 SLBM_LIB int slbm_shell_saveVelocityModelFormat ( const char* modelFileName, int format );
366 
379 SLBM_LIB int slbm_shell_saveVelocityModel ( const char* modelFileName );
380 
397 SLBM_LIB int slbm_shell_createGreatCircle ( char* phase, double* sourceLat, double* sourceLon, double* sourceDepth,
398  double* receiverLat, double* receiverLon, double* receiverDepth );
406 
428 
436 SLBM_LIB int slbm_shell_getDistance ( double* dist );
437 
447 
457 
471 
486 
496 SLBM_LIB int slbm_shell_getTravelTime ( double* travelTime );
497 
513 SLBM_LIB int slbm_shell_getTravelTimeComponents ( double* tTotal, double* tSource, double* tReceiver, double* tHeadwave, double* tGradient );
514 
552 SLBM_LIB int slbm_shell_getWeights ( int nodeId[], double weight[], int* nWeights );
553 
565 // TODO: test this
566 SLBM_LIB int slbm_shell_getWeightsSource ( int nodeids[], double weights[], int* nWeights );
567 
580 // TODO: test this
581 SLBM_LIB int slbm_shell_getWeightsReceiver ( int nodeids[], double weights[], int* nWeights );
582 
604 SLBM_LIB int slbm_shell_toString ( char* str, int verbosity );
605 
611 SLBM_LIB int slbm_shell_getNGridNodes ( int* numGridNodes );
612 
625 
646 SLBM_LIB int slbm_shell_getGridData ( int* nodeId, double* latitude, double* longitude, double* depth,
647  double* pvelocity, double* svelocity, double* gradient );
648 
664 SLBM_LIB int slbm_shell_setGridData ( int* nodeId, double* depth, double* pvelocity, double* svelocity, double* gradient );
665 
696 // TODO: test this
697 SLBM_LIB int slbm_shell_getGreatCircleData ( char* phase, double* path_increment, double sourceDepth[], double sourceVelocity[],
698  double receiverDepth[], double receiverVelocity[], int* npoints, double headWaveVelocity[], double gradient[]);
699 
728 // TODO: test this
729 SLBM_LIB int slbm_shell_getGreatCircleNodeInfo( int** neighbors, double** coefficients, const int* maxpoints,
730  const int* maxnodes, int* npoints, int* nnodes );
731 
756 SLBM_LIB int slbm_shell_getInterpolatedPoint ( double* lat, double* lon,
757  int* nodeIds, double* coefficients, int* nnodes, double* depth, double* pvelocity,
758  double* svelocity, double* pgradient, double* sgradient );
759 
786 // TODO: test this
787 SLBM_LIB int slbm_shell_getInterpolatedTransect ( double lat[], double lon[], int* nLatLon, int** nodeId,
788  double** coefficients, int* nNodes, double depth[][NLAYERS], double pvelocity[][NLAYERS],
789  double svelocity[][NLAYERS], double pgradient[NLAYERS], double sgradient[NLAYERS], int* nInvalid );
790 
803 SLBM_LIB int slbm_shell_initializeActiveNodes ( double* latmin, double* lonmin, double* latmax, double* lonmax );
804 
821 SLBM_LIB int slbm_shell_initActiveNodesFile(char* polygonFileName);
822 
844 SLBM_LIB int slbm_shell_initActiveNodesPoints(double* lat, double* lon, int* npoints, int* inDegrees);
845 
850 
854 
860 SLBM_LIB int slbm_shell_getGridNodeId ( int activeNodeId, int* gridNodeId );
861 
867 SLBM_LIB int slbm_shell_getActiveNodeId ( int gridNodeId );
868 
874 SLBM_LIB int slbm_shell_getNodeHitCount ( int* nodeId, int* hitCount );
875 
883 SLBM_LIB int slbm_shell_getNodeNeighbors ( int* nid, int neighbors[], int* nNeighbors );
884 
892 SLBM_LIB int slbm_shell_getNodeNeighborInfo ( int* nid, int neighbors[], double distance[], double azimuth[], int* nNeighbors );
893 
897 SLBM_LIB int slbm_shell_getNodeSeparation ( int* node1, int* node2, double* distance );
898 
902 SLBM_LIB int slbm_shell_getNodeAzimuth ( int* node1, int* node2, double* azimuth );
903 
913 SLBM_LIB int slbm_shell_getTravelTimeUncertainty ( int* phase, double* distance, double* uncertainty );
914 
922 SLBM_LIB int slbm_shell_getTTUncertainty ( double* uncertainty );
923 
934 
945 SLBM_LIB int slbm_shell_getTTUncertainty1D ( double* uncertainty );
946 
956 SLBM_LIB int slbm_shell_getSlownessUncertainty( int* phase, double* distance, double* uncert );
957 
965 SLBM_LIB int slbm_shell_getSHUncertainty(double* slownessUncertainty);
966 
989 SLBM_LIB int slbm_shell_getZhaoParameters ( double* Vm, double* Gm, double* H, double* C, double* Cm, int* udSign );
990 
1030 SLBM_LIB int slbm_shell_getActiveNodeWeights ( int nodeId[], double weight[], int* nWeights );
1031 
1039 // TODO: test this
1040 SLBM_LIB int slbm_shell_getActiveNodeWeightsSource ( int nodeids[], double weights[], int* nWeights );
1041 
1049 // TODO: test this
1050 SLBM_LIB int slbm_shell_getActiveNodeWeightsReceiver ( int nodeids[], double weights[], int* nWeights );
1051 
1052 
1058 SLBM_LIB int slbm_shell_getActiveNodeNeighbors ( int* nid, int neighbors[], int* nNeighbors );
1059 
1068 SLBM_LIB int slbm_shell_getActiveNodeNeighborInfo ( int* nid, int neighbors[], double distance[],
1069  double azimuth[], int* nNeighbors );
1070 
1092  double* latitude,
1093  double* longitude,
1094  double depth[NLAYERS],
1095  double pvelocity[NLAYERS],
1096  double svelocity[NLAYERS],
1097  double gradient[2] );
1098 
1116  double depth[NLAYERS],
1117  double pvelocity[NLAYERS],
1118  double svelocity[NLAYERS],
1119  double gradient[2] );
1120 
1130 SLBM_LIB int slbm_shell_setCHMax ( double* chMax );
1131 
1141 SLBM_LIB int slbm_shell_getCHMax ( double* chMax );
1142 
1153 SLBM_LIB int slbm_shell_getAverageMantleVelocity ( int* type, double* velocity );
1154 
1167 SLBM_LIB int slbm_shell_setAverageMantleVelocity ( int* type, double* velocity );
1168 
1173 SLBM_LIB int slbm_shell_getTessId ( char* tessId );
1174 
1181 SLBM_LIB int slbm_shell_getFractionActive ( double* fractionActive );
1182 
1191 SLBM_LIB int slbm_shell_setMaxDistance ( const double* maxDistance );
1192 
1199 SLBM_LIB int slbm_shell_getMaxDistance ( double* maxDistance );
1200 
1209 SLBM_LIB int slbm_shell_setMaxDepth ( const double* maxDepth );
1210 
1215 SLBM_LIB int slbm_shell_getMaxDepth ( double* maxDepth );
1216 
1234 SLBM_LIB int slbm_shell_getPgLgComponents(double* tTotal, double* tTaup,
1235  double* tHeadwave, double* pTaup, double* pHeadwave, double* trTaup, double* trHeadwave );
1236 
1243 SLBM_LIB int slbm_shell_getSlowness(double* slowness);
1244 
1245 // new gtb 16Jan2009
1246 //-------------------------------------------------------------
1254 SLBM_LIB int slbm_shell_get_dtt_dlat(double* dtt_dlat);
1255 
1263 SLBM_LIB int slbm_shell_get_dtt_dlon(double* dtt_dlon);
1264 
1272 SLBM_LIB int slbm_shell_get_dtt_ddepth(double* dtt_ddepth);
1273 
1282 //int slbm_shell_get_dsh_ddist(double* dsh_ddist);
1283 //
1291 //int slbm_shell_get_dsh_dlat(double* dsh_dlat);
1292 //
1300 //int slbm_shell_get_dsh_dlon(double* dsh_dlon);
1301 //
1309 //int slbm_shell_get_dsh_ddepth(double* dsh_ddepth);
1310 
1311 // new sb 2011/07/27 version 2.7.0
1312 //-------------------------------------------------------------------------
1313 
1330 SLBM_LIB int slbm_shell_getDistAz(double aLat, double aLon, double bLat, double bLon,
1331  double* distance, double* azimuth, double naValue);
1332 
1345 SLBM_LIB int slbm_shell_movePoint(double aLat, double aLon, double distance, double azimuth,
1346  double* bLat, double* bLon);
1347 
1361 SLBM_LIB int slbm_shell_getGreatCircleLocations ( double latitude[], double longitude[] ,
1362  double depth[], int* npoints);
1363 
1381 SLBM_LIB int slbm_shell_getGreatCirclePoints(double aLat, double aLon, double bLat, double bLon,
1382  int npoints, double latitude[], double longitude[]);
1383 
1399 SLBM_LIB int slbm_shell_getGreatCirclePointsOnCenters(double aLat, double aLon, double bLat, double bLon,
1400  int npoints, double latitude[], double longitude[]);
1401 
1411 SLBM_LIB int slbm_shell_getPiercePointSource(double* lat, double* lon, double* depth);
1412 
1422 SLBM_LIB int slbm_shell_getPiercePointReceiver(double* lat, double* lon, double* depth);
1423 
1431 SLBM_LIB int slbm_shell_getDelDistance(double* delDistance);
1432 
1440 SLBM_LIB int slbm_shell_setDelDistance( double delDistance);
1441 
1448 SLBM_LIB int slbm_shell_getDelDepth(double* delDepth);
1449 
1456 SLBM_LIB int slbm_shell_setDelDepth( double delDepth);
1457 
1463 SLBM_LIB int slbm_shell_getRayParameter(double* rayParameter);
1464 
1470 SLBM_LIB int slbm_shell_getTurningRadius(double* turningRadius);
1471 
1485 SLBM_LIB int slbm_shell_setPathIncrement(double pathIncrement);
1486 
1499 SLBM_LIB int slbm_shell_getPathIncrement(double* pathIncrement);
1500 
1508 SLBM_LIB int slbm_shell_setInterpolatorType(char* interpolatorType);
1509 
1515 SLBM_LIB int slbm_shell_getInterpolatorType(char* interpolatorType);
1516 
1522 SLBM_LIB int slbm_shell_getModelString(char* modelString, int* allocatedSize);
1523 
1537 SLBM_LIB int slbm_shell_getUncertaintyTable(int* phaseIndex, int* attributeIndex, char* uncertaintyTable,
1538  int* allocatedSize);
1539 
1553 SLBM_LIB int slbm_shell_getUncertaintyTableFileFormat(int* phaseIndex, int* attributeIndex, char* uncertaintyTable,
1554  int* allocatedSize);
1555 
1556 //-------------------------------------------------------------
1557 
1558 #ifdef __cplusplus
1559 }
1560 #endif
1561 
1562 #endif /* _SLBM_C_SHELL_H */
slbm_shell_getTTUncertainty1D
SLBM_LIB int slbm_shell_getTTUncertainty1D(double *uncertainty)
Retrieve travel time uncertainty in sec. This function will return the non-path-dependent 1D uncertai...
slbm_shell_getWeightsReceiver
SLBM_LIB int slbm_shell_getWeightsReceiver(int nodeids[], double weights[], int *nWeights)
Retrieve the node IDs and the interpolation coefficients for the receiver CrustalProfile.
slbm_shell_getSourceDistance
SLBM_LIB int slbm_shell_getSourceDistance(double *dist)
Retrieve horizontal offset below the source, in radians.
slbm_shell_getTravelTime
SLBM_LIB int slbm_shell_getTravelTime(double *travelTime)
Retrieve the total travel time for the GreatCircle, in seconds.
slbm_shell_getDistance
SLBM_LIB int slbm_shell_getDistance(double *dist)
Retrieve the source-receiver separation, in radians.
slbm_shell_getUncertaintyTableFileFormat
SLBM_LIB int slbm_shell_getUncertaintyTableFileFormat(int *phaseIndex, int *attributeIndex, char *uncertaintyTable, int *allocatedSize)
Retrieve an inconveniently formated table of the uncertainty values for the specified phaseIndex and ...
slbm_shell_getZhaoParameters
SLBM_LIB int slbm_shell_getZhaoParameters(double *Vm, double *Gm, double *H, double *C, double *Cm, int *udSign)
Retrieve some of the parameters that contribute to the calculation of of total travel time using the ...
slbm_shell_getRayParameter
SLBM_LIB int slbm_shell_getRayParameter(double *rayParameter)
Retrieve the ray parameter in sec/km.
slbm_shell_getNActiveNodes
SLBM_LIB int slbm_shell_getNActiveNodes(int *nNodes)
Retrieve the number of active nodes in the Grid.
slbm_shell_getTravelTimeUncertainty
SLBM_LIB int slbm_shell_getTravelTimeUncertainty(int *phase, double *distance, double *uncertainty)
Retrieve the travel time uncertainty in sec for specified phase, distance.
slbm_shell_getActiveNodeWeightsSource
SLBM_LIB int slbm_shell_getActiveNodeWeightsSource(int nodeids[], double weights[], int *nWeights)
Retrieve the active node IDs and the interpolation coefficients for the source CrustalProfile.
slbm_shell_getDelDepth
SLBM_LIB int slbm_shell_getDelDepth(double *delDepth)
Retrieve del_depth in km.
slbm_shell_getActiveNodeNeighbors
SLBM_LIB int slbm_shell_getActiveNodeNeighbors(int *nid, int neighbors[], int *nNeighbors)
Retrieve the node IDs of the nodes that surround the specified node.
slbm_shell_setActiveNodeData
SLBM_LIB int slbm_shell_setActiveNodeData(int *nodeId, double depth[NLAYERS], double pvelocity[NLAYERS], double svelocity[NLAYERS], double gradient[2])
Modify the velocity and gradient information associated with a specified active node in the Grid.
slbm_shell_getGreatCirclePointsOnCenters
SLBM_LIB int slbm_shell_getGreatCirclePointsOnCenters(double aLat, double aLon, double bLat, double bLon, int npoints, double latitude[], double longitude[])
Retrieve an array of lat, lon points along a great circle path between two specified points,...
slbm_shell_getSHUncertainty
SLBM_LIB int slbm_shell_getSHUncertainty(double *slownessUncertainty)
Retrieve uncertainty of the horizontal slowness, in seconds/radian using the phase and distance speci...
SLBM_LIB
#define SLBM_LIB
Definition: slbm_C_shell.h:133
slbm_shell_setMaxDistance
SLBM_LIB int slbm_shell_setMaxDistance(const double *maxDistance)
Set the maximum source-receiver separation for Pn/Sn phase, in radians.
slbm_shell_getInterpolatedTransect
SLBM_LIB int slbm_shell_getInterpolatedTransect(double lat[], double lon[], int *nLatLon, int **nodeId, double **coefficients, int *nNodes, double depth[][NLAYERS], double pvelocity[][NLAYERS], double svelocity[][NLAYERS], double pgradient[NLAYERS], double sgradient[NLAYERS], int *nInvalid)
Retrieve interpolated data from the earth model along a transect defined by equal sized,...
slbm_shell_setMaxDepth
SLBM_LIB int slbm_shell_setMaxDepth(const double *maxDepth)
Set the maximum source depth for Pn/Sn phase, in km.
slbm_shell_initActiveNodesFile
SLBM_LIB int slbm_shell_initActiveNodesFile(char *polygonFileName)
Specify the name of a file that contains a list of points that define a polygon that enclose the set ...
slbm_shell_getDelDistance
SLBM_LIB int slbm_shell_getDelDistance(double *delDistance)
Retrieve del_distance in radians.
slbm_shell_getPgLgComponents
SLBM_LIB int slbm_shell_getPgLgComponents(double *tTotal, double *tTaup, double *tHeadwave, double *pTaup, double *pHeadwave, double *trTaup, double *trHeadwave)
Retrieve information about Pg/Lg travel time calculations.
slbm_shell_getModelString
SLBM_LIB int slbm_shell_getModelString(char *modelString, int *allocatedSize)
Retrieve a string describing the contents of the model. Retrieve a string describing the contents of ...
slbm_shell_getFractionActive
SLBM_LIB int slbm_shell_getFractionActive(double *fractionActive)
Retrieve the fraction of the path length of the current GreatCircle object that is within the current...
slbm_shell_getUncertaintyTable
SLBM_LIB int slbm_shell_getUncertaintyTable(int *phaseIndex, int *attributeIndex, char *uncertaintyTable, int *allocatedSize)
Retrieve a conveniently formated table of the uncertainty values for the specified phaseIndex and att...
slbm_shell_getMaxDistance
SLBM_LIB int slbm_shell_getMaxDistance(double *maxDistance)
Retrieve the current value for the maximum source-receiver separation, in radians.
slbm_shell_setPathIncrement
SLBM_LIB int slbm_shell_setPathIncrement(double pathIncrement)
Set the desired spacing of great circle nodes along the head wave interface, in radians.
slbm_shell_getSlownessUncertainty
SLBM_LIB int slbm_shell_getSlownessUncertainty(int *phase, double *distance, double *uncert)
Retrieve the slowness uncertainty in sec/radian for specified phase, distance.
slbm_shell_getAverageMantleVelocity
SLBM_LIB int slbm_shell_getAverageMantleVelocity(int *type, double *velocity)
Retrieve the average P or S wave mantle velocity that is specified in the model input file,...
slbm_shell_getNodeNeighborInfo
SLBM_LIB int slbm_shell_getNodeNeighborInfo(int *nid, int neighbors[], double distance[], double azimuth[], int *nNeighbors)
Retrieve the node IDs of the nodes that surround the specified node.
slbm_shell_getNHeadWavePoints
SLBM_LIB int slbm_shell_getNHeadWavePoints(int *npoints)
Retrieve the number of LayerProfile objects positioned along the head wave interface.
slbm_shell_getMaxDepth
SLBM_LIB int slbm_shell_getMaxDepth(double *maxDepth)
Retrieve the current value for the maximum source depth, in km.
slbm_shell_getTessId
SLBM_LIB int slbm_shell_getTessId(char *tessId)
Retrieve the tessellation ID of the model currently in memory.
slbm_shell_loadVelocityModel
SLBM_LIB int slbm_shell_loadVelocityModel(const char *modelPath)
Load the velocity model into memory from the specified file or directory. This method automatically d...
slbm_shell_getWeights
SLBM_LIB int slbm_shell_getWeights(int nodeId[], double weight[], int *nWeights)
Retrieve the weight assigned to each grid node that was touched by the GreatCircle.
slbm_shell_setGridData
SLBM_LIB int slbm_shell_setGridData(int *nodeId, double *depth, double *pvelocity, double *svelocity, double *gradient)
Modify the velocity and gradient information associated with a specified node in the Grid.
slbm_shell_setInterpolatorType
SLBM_LIB int slbm_shell_setInterpolatorType(char *interpolatorType)
Specify the interpolation type to use, either 'linear' or 'natural_neighbor'.
slbm_shell_setDelDepth
SLBM_LIB int slbm_shell_setDelDepth(double delDepth)
Set the value of del_depth, in km.
slbm_shell_clearActiveNodes
SLBM_LIB int slbm_shell_clearActiveNodes()
Clear all active nodes. Clear all active nodes.
slbm_shell_movePoint
SLBM_LIB int slbm_shell_movePoint(double aLat, double aLon, double distance, double azimuth, double *bLat, double *bLon)
Find point B that is the specified distance and azimuth from point A. All quantities are in radians.
slbm_shell_getCHMax
SLBM_LIB int slbm_shell_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_shell_getGridData
SLBM_LIB int slbm_shell_getGridData(int *nodeId, double *latitude, double *longitude, double *depth, double *pvelocity, double *svelocity, double *gradient)
Retrieve the lat (radians), lon (radians), interface depths (km), P and S wave interval velocities (k...
slbm_shell_getTravelTimeComponents
SLBM_LIB int slbm_shell_getTravelTimeComponents(double *tTotal, double *tSource, double *tReceiver, double *tHeadwave, double *tGradient)
Retrieve the total travel time and the 4 components that contribute to it for the current GreatCircle...
slbm_shell_getErrorMessage
SLBM_LIB int slbm_shell_getErrorMessage(char *str)
Retrieve last error message.
slbm_shell_getActiveNodeWeights
SLBM_LIB int slbm_shell_getActiveNodeWeights(int nodeId[], double weight[], int *nWeights)
Retrieve the weight assigned to each active node that was touched by the GreatCircle.
slbm_shell_getGreatCircleLocations
SLBM_LIB int slbm_shell_getGreatCircleLocations(double latitude[], double longitude[], double depth[], int *npoints)
Retrieve the latitudes, longitudes and depths of all the profile positions along the moho.
slbm_shell_isValid
SLBM_LIB int slbm_shell_isValid()
Returns true if the current GreatCirlce object has been instantiated and is ready to be interrogated.
slbm_shell_getActiveNodeData
SLBM_LIB int slbm_shell_getActiveNodeData(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_shell_toString
SLBM_LIB int slbm_shell_toString(char *str, int verbosity)
Returns a human-readable string representation of the GreatCircle object.
slbm_shell_saveVelocityModelFormat
SLBM_LIB int slbm_shell_saveVelocityModelFormat(const char *modelFileName, int format)
Save the velocity model currently in memory to the specified file or directory.
slbm_shell_getDistAz
SLBM_LIB int slbm_shell_getDistAz(double aLat, double aLon, double bLat, double bLon, double *distance, double *azimuth, double naValue)
compute distance and azimuth between two points, A and B (all quantities are in radians).
slbm_shell_initializeActiveNodes
SLBM_LIB int slbm_shell_initializeActiveNodes(double *latmin, double *lonmin, double *latmax, double *lonmax)
Specify the latitude and longitude range in radians for active nodes.
slbm_shell_setCHMax
SLBM_LIB int slbm_shell_setCHMax(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_shell_getTTUncertainty
SLBM_LIB int slbm_shell_getTTUncertainty(double *uncertainty)
Retrieve travel time uncertainty in sec. This function will call either the path-dependent or 1D unce...
slbm_shell_loadVelocityModelBinary
SLBM_LIB int slbm_shell_loadVelocityModelBinary(const char *modelDirectory)
Load the velocity model into memory from the specified file or directory. This method automatically d...
slbm_shell_getHeadwaveDistanceKm
SLBM_LIB int slbm_shell_getHeadwaveDistanceKm(double *dist)
Retrieve horizontal distance traveled by the ray below the headwave interface, in radians.
slbm_shell_saveVelocityModel
SLBM_LIB int slbm_shell_saveVelocityModel(const char *modelFileName)
Save the velocity model currently in memory to the specified file.
slbm_shell_getReceiverDistance
SLBM_LIB int slbm_shell_getReceiverDistance(double *dist)
Retrieve horizontal offset below the receiver, in radians.
slbm_shell_get_dtt_dlon
SLBM_LIB int slbm_shell_get_dtt_dlon(double *dtt_dlon)
Retrieve the derivative of travel time wrt to source longitude, in seconds/radian.
slbm_shell_create_fixedEarthRadius
SLBM_LIB int slbm_shell_create_fixedEarthRadius(double *radius)
Instantiate a SLBM Interface object with fixed earth radius.
slbm_shell_getWeightsSource
SLBM_LIB int slbm_shell_getWeightsSource(int nodeids[], double weights[], int *nWeights)
Retrieve the node IDs and the interpolation coefficients for the source CrustalProfile.
slbm_shell_getInterpolatorType
SLBM_LIB int slbm_shell_getInterpolatorType(char *interpolatorType)
Retrieve the type of interpolator currently in use; either "LINEAR" or "NATUTAL_NEIGHBOR".
slbm_shell_getNodeSeparation
SLBM_LIB int slbm_shell_getNodeSeparation(int *node1, int *node2, double *distance)
Retrieve the angular separation of two grid nodes, in radians.
slbm_shell_getNodeNeighbors
SLBM_LIB int slbm_shell_getNodeNeighbors(int *nid, int neighbors[], int *nNeighbors)
Retrieve the node IDs of the nodes that surround the specified node.
slbm_shell_getGreatCircleData
SLBM_LIB int slbm_shell_getGreatCircleData(char *phase, double *path_increment, double sourceDepth[], double sourceVelocity[], double receiverDepth[], double receiverVelocity[], int *npoints, double headWaveVelocity[], double gradient[])
Retrieve the data required for input to the travel time calculation, for the current GreatCircle obje...
slbm_shell_getActiveNodeId
SLBM_LIB int slbm_shell_getActiveNodeId(int gridNodeId)
Retrieve the active node ID that corresponds to a specified grid node ID.
slbm_shell_getGreatCircleNodeInfo
SLBM_LIB int slbm_shell_getGreatCircleNodeInfo(int **neighbors, double **coefficients, const int *maxpoints, const int *maxnodes, int *npoints, int *nnodes)
Retrieve information about the interpolated points along the headwave path, including the number of p...
slbm_shell_getGridNodeId
SLBM_LIB int slbm_shell_getGridNodeId(int activeNodeId, int *gridNodeId)
Retrieve the grid node ID that corresponds to a specified active node ID.
slbm_shell_setDelDistance
SLBM_LIB int slbm_shell_setDelDistance(double delDistance)
Set the value of del_distance, radians.
slbm_shell_get_dtt_dlat
SLBM_LIB int slbm_shell_get_dtt_dlat(double *dtt_dlat)
Retrieve the derivative of travel time wrt to source latitude, in seconds/radian.
slbm_shell_getVersion
SLBM_LIB int slbm_shell_getVersion(char *str)
Retrieve the SLBM Version number.
slbm_shell_getPiercePointSource
SLBM_LIB int slbm_shell_getPiercePointSource(double *lat, double *lon, double *depth)
Retrieve the latitude and longitude of the moho pierce point below the source, in radians.
slbm_shell_get_dtt_ddepth
SLBM_LIB int slbm_shell_get_dtt_ddepth(double *dtt_ddepth)
Retrieve the derivative of travel time wrt to source depth, in seconds/km.
slbm_shell_getNodeAzimuth
SLBM_LIB int slbm_shell_getNodeAzimuth(int *node1, int *node2, double *azimuth)
Retrieve the azimuth from grid node1 to grid node2, radians.
slbm_shell_specifyOutputDirectory
SLBM_LIB int slbm_shell_specifyOutputDirectory(const char *directoryName)
Deprecated. Use saveVelocityModel() instead. Specify the directory into which the model that is curre...
slbm_shell_getTurningRadius
SLBM_LIB int slbm_shell_getTurningRadius(double *turningRadius)
Retrieve turning radius in km.
slbm_shell_getGreatCirclePoints
SLBM_LIB int slbm_shell_getGreatCirclePoints(double aLat, double aLon, double bLat, double bLon, int npoints, double latitude[], double longitude[])
Retrieve an array of lat, lon points along a great circle path between two specified points,...
slbm_shell_getActiveNodeNeighborInfo
SLBM_LIB int slbm_shell_getActiveNodeNeighborInfo(int *nid, int neighbors[], double distance[], double azimuth[], int *nNeighbors)
Retrieve the node IDs of the nodes that surround the specified node.
slbm_shell_getPathIncrement
SLBM_LIB int slbm_shell_getPathIncrement(double *pathIncrement)
Retrieve the current value of the spacing of great circle nodes along the head wave interface,...
NLAYERS
#define NLAYERS
Definition: slbm_C_shell.h:195
slbm_shell_getNGridNodes
SLBM_LIB int slbm_shell_getNGridNodes(int *numGridNodes)
Retrieve the number of Grid nodes in the Earth model.
slbm_shell_getInterpolatedPoint
SLBM_LIB int slbm_shell_getInterpolatedPoint(double *lat, double *lon, int *nodeIds, double *coefficients, int *nnodes, double *depth, double *pvelocity, double *svelocity, double *pgradient, double *sgradient)
Retrieve interpolated data from the earth model at a single specified latitude, longitude.
slbm_shell_getTTUncertainty_useRandErr
SLBM_LIB int slbm_shell_getTTUncertainty_useRandErr(double *uncertainty)
Retrieve travel time uncertainty in sec. This function will call either the path-dependent or 1D unce...
slbm_shell_getHeadwaveDistance
SLBM_LIB int slbm_shell_getHeadwaveDistance(double *dist)
Retrieve angular distance traveled by the ray below the headwave interface, in radians.
slbm_shell_getPiercePointReceiver
SLBM_LIB int slbm_shell_getPiercePointReceiver(double *lat, double *lon, double *depth)
Retrieve the latitude and longitude of the moho pierce point below the receiver, in radians.
slbm_shell_getActiveNodeWeightsReceiver
SLBM_LIB int slbm_shell_getActiveNodeWeightsReceiver(int nodeids[], double weights[], int *nWeights)
Retrieve the active node IDs and the interpolation coefficients for the receiver CrustalProfile.
slbm_shell_saveVelocityModelBinary
SLBM_LIB int slbm_shell_saveVelocityModelBinary()
Deprecated. Use saveVelocityModel() instead. Write the model in format 3 to directory previously spec...
slbm_shell_create
SLBM_LIB int slbm_shell_create()
Instantiate a SLBM Interface object.
slbm_shell_clear
SLBM_LIB int slbm_shell_clear()
Delete the current GreatCircle object from memory and clear the pool of stored CrustalProfile objects...
slbm_shell_getSlowness
SLBM_LIB int slbm_shell_getSlowness(double *slowness)
Retrieve the horizontal slowness, i.e., the derivative of travel time wrt to receiver-source distance...
slbm_shell_getNodeHitCount
SLBM_LIB int slbm_shell_getNodeHitCount(int *nodeId, int *hitCount)
Retrieve the number of times that node has been 'touched' by a GreatCircle object.
slbm_shell_initActiveNodesPoints
SLBM_LIB int slbm_shell_initActiveNodesPoints(double *lat, double *lon, int *npoints, int *inDegrees)
Specify a list of points that define a polygon that enclose the set of grid nodes that are to be cons...
slbm_shell_delete
SLBM_LIB int slbm_shell_delete()
Deletes the SlbmInterface object instantiated with the call slbm_shell_create().
slbm_shell_createGreatCircle
SLBM_LIB int slbm_shell_createGreatCircle(char *phase, double *sourceLat, double *sourceLon, double *sourceDepth, double *receiverLat, double *receiverLon, double *receiverDepth)
Instantiate a new GreatCircle object between two locations.
slbm_shell_setAverageMantleVelocity
SLBM_LIB int slbm_shell_setAverageMantleVelocity(int *type, double *velocity)
Set the average P or S wave mantle velocity that is recorded in the model input file,...