Obsidian
Probabilistic Geophysical Joint Inversion
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
obsidian::fwd Namespace Reference

Namespaces

 detail
 

Classes

struct  GravmagInterpolatorParams
 Structure containing gravity and magnetic interpolation parameters. More...
 

Functions

template<>
ContactPointCache generateCache< ForwardModel::CONTACTPOINT > (const std::vector< world::InterpolatorSpec > &boundaryInterpolation, const WorldSpec &worldSpec, const ContactPointSpec &spec)
 Generate a cache object for a contact point forward model. More...
 
template<>
ContactPointResults forwardModel< ForwardModel::CONTACTPOINT > (const ContactPointSpec &spec, const ContactPointCache &cache, const WorldParams &world)
 Run a contact point forward model. More...
 
template<ForwardModel f>
Types< f >::Cache generateCache (const std::vector< world::InterpolatorSpec > &boundaryInterpolation, const WorldSpec &worldSpec, const typename Types< f >::Spec &spec)
 Generate a cache object for a specific forward model. Cache objects contain repeatly used information that only needs to be computed once by the forward model. More...
 
template<ForwardModel f>
Types< f >::Results forwardModel (const typename Types< f >::Spec &spec, const typename Types< f >::Cache &cache, const WorldParams &world)
 Run a particular forward model. More...
 
GlobalCache generateGlobalCache (const std::vector< world::InterpolatorSpec > &boundaryInterpolation, const GlobalSpec &spec, const std::set< ForwardModel > &enabled)
 Generate a global cache object for all forward models. More...
 
GlobalResults forwardModelAll (const GlobalSpec &spec, const GlobalCache &cache, const GlobalParams &params, const std::set< ForwardModel > &enabled)
 Run all the forward models. More...
 
template<>
GravCache generateCache< ForwardModel::GRAVITY > (const std::vector< world::InterpolatorSpec > &boundaryInterpolation, const WorldSpec &worldSpec, const GravSpec &gravSpec)
 Generate a cache object for a gravity forward model. More...
 
template<>
GravResults forwardModel< ForwardModel::GRAVITY > (const GravSpec &spec, const GravCache &cache, const WorldParams &world)
 Run a gravity forward model. More...
 
Eigen::MatrixXd gravSens (const Eigen::VectorXd &xEdges, const Eigen::VectorXd &yEdges, const Eigen::VectorXd &zEdges, const Eigen::MatrixXd &locations)
 Compute the gravity sensitivity matrix. More...
 
GravmagInterpolatorParams makeInterpParams (const VoxelSpec &voxelisation, const Eigen::MatrixXd &locations, const WorldSpec &worldSpec)
 Create a GravmagInterpolatorParams object. More...
 
Eigen::VectorXd computeField (const Eigen::MatrixXd &sens, const Eigen::MatrixXi sensorIndices, const Eigen::MatrixXd sensorWeights, const Eigen::VectorXd &properties)
 Computes the field values for either gravity or magnetic. More...
 
template<>
MagCache generateCache< ForwardModel::MAGNETICS > (const std::vector< world::InterpolatorSpec > &boundaryInterpolation, const WorldSpec &worldSpec, const MagSpec &magSpec)
 Generate a cache object for a magnetic forward model. More...
 
template<>
MagResults forwardModel< ForwardModel::MAGNETICS > (const MagSpec &spec, const MagCache &cache, const WorldParams &world)
 Run a magnetic forward model. More...
 
Eigen::MatrixXd magSens (const Eigen::VectorXd &xEdges, const Eigen::VectorXd &yEdges, const Eigen::VectorXd &zEdges, const Eigen::MatrixXd &locations, const double &bX, const double &bY, const double &bZ)
 Compute the magnetic sensitivity matrix. More...
 
Eigen::MatrixXd magSens (const Eigen::VectorXd &xEdges, const Eigen::VectorXd &yEdges, const Eigen::VectorXd &zEdges, const Eigen::MatrixXd &locations, const Eigen::VectorXd &bX, const Eigen::VectorXd &bY, const Eigen::VectorXd &bZ)
 Compute the magnetic sensitivity matrix. More...
 
Eigen::MatrixX4d phaseTensor1d (const Eigen::MatrixX4cd &Z)
 Calculate the phase tensor given the impedance matrix. More...
 
Eigen::VectorXd alpha1d (const Eigen::MatrixX4d &tensor)
 Computes the alpha value of the phase tensor, which is a measure of its rotations from the major axis. More...
 
Eigen::VectorXd beta1d (const Eigen::MatrixX4d &tensor)
 Computes the beta value of the phase tensor, which is a measure of its rotations from the z-axis. More...
 
Eigen::MatrixX4cd impedenceAniso1d (const Eigen::VectorXd &thicknesses, const Eigen::VectorXd &freqs, const Eigen::VectorXd &resx, const Eigen::VectorXd &resy, const Eigen::VectorXd &phases)
 Calculate the impedence matrix for anisotropic 1D MT. More...
 
template<>
MtAnisoCache generateCache< ForwardModel::MTANISO > (const std::vector< world::InterpolatorSpec > &boundaryInterpolation, const WorldSpec &worldSpec, const MtAnisoSpec &mtSpec)
 Generate a cache object for a MT forward model. More...
 
template<>
MtAnisoResults forwardModel< ForwardModel::MTANISO > (const MtAnisoSpec &spec, const MtAnisoCache &cache, const WorldParams &world)
 Run a MT forward model. More...
 
template<>
Seismic1dCache generateCache< ForwardModel::SEISMIC1D > (const std::vector< world::InterpolatorSpec > &boundaryInterpolation, const WorldSpec &worldSpec, const Seismic1dSpec &spec)
 Generate a cache object for a seismic forward model. More...
 
template<>
Seismic1dResults forwardModel< ForwardModel::SEISMIC1D > (const Seismic1dSpec &spec, const Seismic1dCache &cache, const WorldParams &world)
 Run a seismic forward model. More...
 
boost::multi_array< double, 3 > fillAndPad (const Eigen::VectorXd &values, const ThermalSpec &spec)
 
boost::multi_array< double, 3 > eastWest (const boost::multi_array< double, 3 > &pad)
 
boost::multi_array< double, 3 > cells (const boost::multi_array< double, 3 > &pad, uint axis)
 
boost::multi_array< double, 3 > isocells (const boost::multi_array< double, 3 > &pad)
 
boost::multi_array< double, 3 > temp (const boost::multi_array< double, 3 > &eastwest, const boost::multi_array< double, 3 > &northsouth, const boost::multi_array< double, 3 > &updown, const boost::multi_array< double, 3 > &sCells, const Eigen::MatrixXd &tempZ0, const Eigen::MatrixXd &zLowBound, bool zLowBoundIsHeatFlow, double xSize, double ySize, double zSize)
 
uint clip (int in, uint maxX)
 
Eigen::VectorXd evalAtLocations (const boost::multi_array< double, 3 > &tempVox, const Eigen::MatrixXd &locations, const std::pair< double, double > &xSize, const std::pair< double, double > &ySize, const std::pair< double, double > &zSize)
 
template<>
ThermalCache generateCache< ForwardModel::THERMAL > (const std::vector< world::InterpolatorSpec > &boundaryInterpolation, const WorldSpec &worldSpec, const ThermalSpec &thermSpec)
 
template<>
ThermalResults forwardModel< ForwardModel::THERMAL > (const ThermalSpec &spec, const ThermalCache &cache, const WorldParams &world)
 

Detailed Description

Namespace for all the forward models.

Function Documentation

Eigen::VectorXd obsidian::fwd::alpha1d ( const Eigen::MatrixX4d &  tensor)

Computes the alpha value of the phase tensor, which is a measure of its rotations from the major axis.

Parameters
tensorThe phase tensor computed by phaseTensor1d()
Eigen::VectorXd obsidian::fwd::beta1d ( const Eigen::MatrixX4d &  tensor)

Computes the beta value of the phase tensor, which is a measure of its rotations from the z-axis.

Parameters
tensorThe phase tensor computed by phaseTensor1d()
Eigen::VectorXd obsidian::fwd::computeField ( const Eigen::MatrixXd &  sens,
const Eigen::MatrixXi  sensorIndices,
const Eigen::MatrixXd  sensorWeights,
const Eigen::VectorXd &  properties 
)

Computes the field values for either gravity or magnetic.

Parameters
sensThe gravity or magnetic sensitivity matrix.
sensorIndicesThe indices of each of the sensors.
sensorWeightsThe weights of each of the sensors.
propertiesThe rock property for the sensor type.
template<ForwardModel f>
Types<f>::Results obsidian::fwd::forwardModel ( const typename Types< f >::Spec &  spec,
const typename Types< f >::Cache &  cache,
const WorldParams &  world 
)

Run a particular forward model.

Parameters
specThe forward model specification.
cacheThe forward model cache generated by generateCache().
worldThe world model parameters.
Returns
Forward model results.
template<>
ContactPointResults obsidian::fwd::forwardModel< ForwardModel::CONTACTPOINT > ( const ContactPointSpec &  spec,
const ContactPointCache &  cache,
const WorldParams &  world 
)

Run a contact point forward model.

Parameters
specThe forward model specification.
cacheThe forward model cache generated by generateCache().
worldThe world model parameters.
Returns
Forward model results.
template<>
GravResults obsidian::fwd::forwardModel< ForwardModel::GRAVITY > ( const GravSpec &  spec,
const GravCache &  cache,
const WorldParams &  world 
)

Run a gravity forward model.

Parameters
specThe forward model specification.
cacheThe forward model cache generated by generateCache().
worldThe world model parameters.
Returns
Forward model results.
template<>
MagResults obsidian::fwd::forwardModel< ForwardModel::MAGNETICS > ( const MagSpec &  spec,
const MagCache &  cache,
const WorldParams &  world 
)

Run a magnetic forward model.

Parameters
specThe forward model specification.
cacheThe forward model cache generated by generateCache().
worldThe world model parameters.
Returns
Forward model results.
template<>
MtAnisoResults obsidian::fwd::forwardModel< ForwardModel::MTANISO > ( const MtAnisoSpec &  spec,
const MtAnisoCache &  cache,
const WorldParams &  world 
)

Run a MT forward model.

Parameters
specThe forward model specification.
cacheThe forward model cache generated by generateCache().
worldThe world model parameters.
Returns
Forward model results.
template<>
Seismic1dResults obsidian::fwd::forwardModel< ForwardModel::SEISMIC1D > ( const Seismic1dSpec &  spec,
const Seismic1dCache &  cache,
const WorldParams &  world 
)

Run a seismic forward model.

Parameters
specThe forward model specification.
cacheThe forward model cache generated by generateCache().
worldThe world model parameters.
Returns
Forward model results.
GlobalResults obsidian::fwd::forwardModelAll ( const GlobalSpec &  spec,
const GlobalCache &  cache,
const GlobalParams &  params,
const std::set< ForwardModel > &  enabled 
)

Run all the forward models.

Parameters
specThe global world specifications.
cacheThe global forward model cache generated by generateCacheGlobal().
paramsThe global world parameters.
Returns
Results of all the forward models.
template<ForwardModel f>
Types<f>::Cache obsidian::fwd::generateCache ( const std::vector< world::InterpolatorSpec > &  boundaryInterpolation,
const WorldSpec &  worldSpec,
const typename Types< f >::Spec &  spec 
)

Generate a cache object for a specific forward model. Cache objects contain repeatly used information that only needs to be computed once by the forward model.

Parameters
boundaryInterpolationThe world model interpolation parameters.
worldSpecThe world model specification.
specThe forward model specification.
Returns
Forward model cache object.
template<>
ContactPointCache obsidian::fwd::generateCache< ForwardModel::CONTACTPOINT > ( const std::vector< world::InterpolatorSpec > &  boundaryInterpolation,
const WorldSpec &  worldSpec,
const ContactPointSpec &  spec 
)

Generate a cache object for a contact point forward model.

Parameters
boundaryInterpolationThe world model interpolation parameters.
worldSpecThe world model specification.
specThe forward model specification.
Returns
Forward model cache object.
template<>
GravCache obsidian::fwd::generateCache< ForwardModel::GRAVITY > ( const std::vector< world::InterpolatorSpec > &  boundaryInterpolation,
const WorldSpec &  worldSpec,
const GravSpec &  gravSpec 
)

Generate a cache object for a gravity forward model.

Parameters
boundaryInterpolationThe world model interpolation parameters.
worldSpecThe world model specification.
gravSpecThe forward model specification.
Returns
Forward model cache object.
template<>
MagCache obsidian::fwd::generateCache< ForwardModel::MAGNETICS > ( const std::vector< world::InterpolatorSpec > &  boundaryInterpolation,
const WorldSpec &  worldSpec,
const MagSpec &  magSpec 
)

Generate a cache object for a magnetic forward model.

Parameters
boundaryInterpolationThe world model interpolation parameters.
worldSpecThe world model specification.
magSpecThe forward model specification.
Returns
Forward model cache object.
template<>
MtAnisoCache obsidian::fwd::generateCache< ForwardModel::MTANISO > ( const std::vector< world::InterpolatorSpec > &  boundaryInterpolation,
const WorldSpec &  worldSpec,
const MtAnisoSpec &  mtSpec 
)

Generate a cache object for a MT forward model.

Parameters
boundaryInterpolationThe world model interpolation parameters.
worldSpecThe world model specification.
mtSpecThe forward model specification.
Returns
Forward model cache object.
template<>
Seismic1dCache obsidian::fwd::generateCache< ForwardModel::SEISMIC1D > ( const std::vector< world::InterpolatorSpec > &  boundaryInterpolation,
const WorldSpec &  worldSpec,
const Seismic1dSpec &  spec 
)

Generate a cache object for a seismic forward model.

Parameters
boundaryInterpolationThe world model interpolation parameters.
worldSpecThe world model specification.
specThe forward model specification.
Returns
Forward model cache object.
GlobalCache obsidian::fwd::generateGlobalCache ( const std::vector< world::InterpolatorSpec > &  boundaryInterpolation,
const GlobalSpec &  spec,
const std::set< ForwardModel > &  enabled 
)

Generate a global cache object for all forward models.

Parameters
boundaryInterpolationThe world model interpolation parameters.
specThe global world specifications.
Returns
Cache object containing the caches of all the forward models.
Eigen::MatrixXd obsidian::fwd::gravSens ( const Eigen::VectorXd &  xEdges,
const Eigen::VectorXd &  yEdges,
const Eigen::VectorXd &  zEdges,
const Eigen::MatrixXd &  locations 
)

Compute the gravity sensitivity matrix.

Parameters
xEdgesThe x coordinates of the mesh grid.
yEdgesThe y coordinates of the mesh grid.
zEdgesThe z coordinates of the mesh grid.
locationsA Nx3 matrix containing the coordinates of the sensor locations.
Eigen::MatrixX4cd obsidian::fwd::impedenceAniso1d ( const Eigen::VectorXd &  thicknesses,
const Eigen::VectorXd &  freqs,
const Eigen::VectorXd &  resx,
const Eigen::VectorXd &  resy,
const Eigen::VectorXd &  phases 
)

Calculate the impedence matrix for anisotropic 1D MT.

Parameters
thicknessesA vector containing the thicknesses for each layer.
freqsA vector containing the MT frequencies.
resxA vector containing the resistivities in x direction (not logarithmic).
resyA vector containing the resistivities in y direction (not logarithmic).
phasesA vector containing the phase of each layer.
Returns
The impedance matrix for anisotropic MT.
Eigen::MatrixXd obsidian::fwd::magSens ( const Eigen::VectorXd &  xEdges,
const Eigen::VectorXd &  yEdges,
const Eigen::VectorXd &  zEdges,
const Eigen::MatrixXd &  locations,
const double &  bX,
const double &  bY,
const double &  bZ 
)

Compute the magnetic sensitivity matrix.

Parameters
xEdgesThe x coordinates of the mesh grid.
yEdgesThe y coordinates of the mesh grid.
zEdgesThe z coordinates of the mesh grid.
locationsA Nx3 matrix containing the coordinates of the sensor locations.
bX,bY,bZThe field values along each axis. Assuming the field values are constant in each axis.
Eigen::MatrixXd obsidian::fwd::magSens ( const Eigen::VectorXd &  xEdges,
const Eigen::VectorXd &  yEdges,
const Eigen::VectorXd &  zEdges,
const Eigen::MatrixXd &  locations,
const Eigen::VectorXd &  bX,
const Eigen::VectorXd &  bY,
const Eigen::VectorXd &  bZ 
)

Compute the magnetic sensitivity matrix.

Parameters
xEdgesThe x coordinates of the mesh grid.
yEdgesThe y coordinates of the mesh grid.
zEdgesThe z coordinates of the mesh grid.
locationsA Nx3 matrix containing the coordinates of the sensor locations.
bX,bY,bZThe field values at each of the mesh grid points.
GravmagInterpolatorParams obsidian::fwd::makeInterpParams ( const VoxelSpec &  voxelisation,
const Eigen::MatrixXd &  locations,
const WorldSpec &  worldSpec 
)

Create a GravmagInterpolatorParams object.

Parameters
voxelisationThe voxelised world.
locationsThe sensor locations.
worldSpecThe world specifications.
Eigen::MatrixX4d obsidian::fwd::phaseTensor1d ( const Eigen::MatrixX4cd &  Z)

Calculate the phase tensor given the impedance matrix.

Parameters
ZThe impedance matrix.
Returns
The phase tensor.