Loading [MathJax]/extensions/tex2jax.js
1.7.3
Geostatistics & Machine Learning toolbox | https://gstlearn.org
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
ACov Class Referenceabstract

Class containing the Covariance part of the Model. More...

#include <ACov.hpp>

Inheritance diagram for ACov:
ASpaceObject ICloneable AStringable CorAniso CorGneiting CorMatern CovBase CovList CovProportional CovAnisoList CovAniso CovLMCAnamorphosis CovLMCConvolution CovLMCTapering CovLMGradient ACovGradient CovGradientFunctional CovGradientNumerical

Detailed Description

Class containing the Covariance part of the Model.

It is the uppermost class of the Covariance Tree and is conceived as simple as possible on purpose (in order to let the user defined its own version if necessary): it must simply be able to return its value between two end-point (see eval method).

It is mainly implemented in CovAniso.hpp or CovAnisoList.hpp

Public Member Functions

 ACov (const CovContext &ctxt=CovContext())
 
 ACov (const ACov &r)
 
ACovoperator= (const ACov &r)
 
virtual ~ACov ()
 
virtual int getNVar () const
 ACov Interface.
 
virtual bool isIndexable () const
 
bool isNoStat () const
 
virtual void loadInfoValues ()
 
const CovContextgetContext () const
 
void setContext (const CovContext &ctxt)
 
void updateFromContext ()
 
virtual void copyCovContext (const CovContext &ctxt)
 
void initFromContext ()
 
CovContext getContextCopy () const
 
virtual double eval0 (int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 Calculate the covariance between two variables for 0-distance (stationary case)
 
double evalCov (const SpacePoint &p1, const SpacePoint &p2, int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 Calculate the covariance between two variables and two points (general case)
 
virtual double evalCovOnSphere (double alpha, int degree=50, bool flagScaleDistance=false, const CovCalcMode *mode=nullptr) const
 
virtual VectorDouble evalSpectrumOnSphere (int n, bool flagNormDistance=false, bool flagCumul=false) const
 
virtual double evalSpectrum (const VectorDouble &freq, int ivar, int jvar) const
 
virtual void updateCovByPoints (int icas1, int iech1, int icas2, int iech2) const
 
void attachNoStatDb (const Db *db)
 
virtual bool isOptimEnabled () const
 Functions linked to Optimization during Covariance calculations.
 
void optimizationPreProcess (int mode, const std::vector< SpacePoint > &ps) const
 
SpacePointoptimizationLoadInPlace (int iech, int mode, int rank) const
 
void optimizationPostProcess () const
 
void optimizationSetTarget (SpacePoint &pt) const
 
VectorDouble eval (const std::vector< SpacePoint > &vec_p1, const std::vector< SpacePoint > &vec_p2, int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 Functions for evaluating Covariances.
 
MatrixSymmetric eval0Mat (const CovCalcMode *mode=nullptr) const
 
MatrixSymmetric evalCovMat0 (const Db *db, int iech, const KrigOpt &krigopt=KrigOpt()) const
 Functions for evaluating Covariance Matrices either in place or not.
 
MatrixDense evalCovMat (const Db *db1, const Db *db2=nullptr, int ivar0=-1, int jvar0=-1, const VectorInt &nbgh1=VectorInt(), const VectorInt &nbgh2=VectorInt(), const CovCalcMode *mode=nullptr, bool cleanOptim=true) const
 
MatrixSymmetric evalCovMatSym (const Db *db1, const VectorInt &nbgh1=VectorInt(), int ivar0=-1, const CovCalcMode *mode=nullptr, bool cleanOptim=true) const
 
MatrixSparseevalCovMatSparse (const Db *db1_arg, const Db *db2_arg=nullptr, int ivar0=-1, int jvar0=-1, const VectorInt &nbgh1=VectorInt(), const VectorInt &nbgh2=VectorInt(), const CovCalcMode *mode=nullptr, bool cleanOptim=true, double eps=EPSILON3) const
 
int evalCovMat0InPlace (MatrixSymmetric &mat, const Db *db, int iech, const KrigOpt &krigopt=KrigOpt()) const
 
int evalCovMatInPlace (MatrixDense &mat, const Db *db1, const Db *db2=nullptr, int ivar0=-1, int jvar0=-1, const VectorInt &nbgh1=VectorInt(), const VectorInt &nbgh2=VectorInt(), const CovCalcMode *mode=nullptr, bool cleanOptim=true) const
 
int evalCovMatSymInPlace (MatrixSymmetric &mat, const Db *db1, const VectorInt &nbgh1=VectorInt(), int ivar0=-1, const CovCalcMode *mode=nullptr, bool cleanOptim=true) const
 
int evalCovMatInPlaceFromIdx (MatrixDense &mat, const Db *db1, const Db *db2, const VectorVectorInt &index1, const VectorVectorInt &index2, const VectorInt &nbgh2=VectorInt(), const CovCalcMode *mode=nullptr, bool cleanOptim=true) const
 
int evalCovMatSymInPlaceFromIdx (MatrixSymmetric &mat, const Db *db1, const VectorVectorInt &index1, const CovCalcMode *mode=nullptr, bool cleanOptim=true) const
 
int evalCovMatRHSInPlaceFromIdx (MatrixDense &mat, const Db *db1, const Db *db2, const VectorVectorInt &index1, const int iech2=-1, const KrigOpt &krigopt=KrigOpt(), bool cleanOptim=true) const
 
int evalCovVecRHSInPlace (vect vect, const Db *db2, const VectorInt &index1, int iech2, const KrigOpt &krigopt, SpacePoint &pin, SpacePoint &pout, VectorDouble &tabwork, double lambda=1.) const
 
virtual int addEvalCovVecRHSInPlace (vect vect, const VectorInt &index1, const int iech2, const KrigOpt &krigopt, SpacePoint &pin, SpacePoint &pout, VectorDouble &tabwork, double lambda=1.) const
 
void eval0CovMatBiPointInPlace (MatrixSymmetric &mat, const CovCalcMode *mode) const
 
double evalIvarIpas (double step, const VectorDouble &dir=VectorDouble(), int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
double evalIvarIpasIncr (const VectorDouble &dincr, int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
VectorDouble evalIvarNlag (const VectorDouble &vec_step, const VectorDouble &dir=VectorDouble(), int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
MatrixSquare evalNvarIpas (double step, const VectorDouble &dir=VectorDouble(), const CovCalcMode *mode=nullptr) const
 
MatrixSquare evalNvarIpasIncr (const VectorDouble &dincr, const CovCalcMode *mode=nullptr) const
 
double evalIsoIvarIpas (double step, int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
VectorDouble evalIsoIvarNlag (const VectorDouble &vec_step, int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
MatrixSquare evalIsoNvarIpas (double step, const CovCalcMode *mode=nullptr) const
 
double evalCvv (const VectorDouble &ext, const VectorInt &ndisc, const VectorDouble &angles=VectorDouble(), int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
double evalCvvShift (const VectorDouble &ext, const VectorInt &ndisc, const VectorDouble &shift, const VectorDouble &angles=VectorDouble(), int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
MatrixSquare evalCvvM (const VectorDouble &ext, const VectorInt &ndisc, const VectorDouble &angles=VectorDouble(), const CovCalcMode *mode=nullptr) const
 
double evalCxv (const SpacePoint &p1, const VectorDouble &ext, const VectorInt &ndisc, const VectorDouble &angles=VectorDouble(), const VectorDouble &x0=VectorDouble(), int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
double evalCxv (const Db *db, const VectorDouble &ext, const VectorInt &ndisc, const VectorDouble &angles=VectorDouble(), const VectorDouble &x0=VectorDouble(), int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
MatrixSquare evalCxvM (const SpacePoint &p1, const VectorDouble &ext, const VectorInt &ndisc, const VectorDouble &angles=VectorDouble(), const VectorDouble &x0=VectorDouble(), const CovCalcMode *mode=nullptr) const
 
void evalPointToDb (VectorDouble &values, const SpacePoint &p1, const Db *db2, int ivar=0, int jvar=0, bool useSel=true, const VectorInt &nbgh2=VectorInt(), const CovCalcMode *mode=nullptr) const
 
void evalPointToDbAsSP (VectorDouble &values, const std::vector< SpacePoint > &p1s, const SpacePoint &p2, int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
double evalAverageDbToDb (const Db *db1, const Db *db2, int ivar=0, int jvar=0, double eps=0., int seed=434132, const CovCalcMode *mode=nullptr) const
 
double evalAverageIncrToIncr (const VectorVectorDouble &d1, const VectorVectorDouble &d2, int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
double evalAveragePointToDb (const SpacePoint &p1, const Db *db2, int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
double extensionVariance (const Db *db, const VectorDouble &ext, const VectorInt &ndisc, const VectorDouble &angles=VectorDouble(), const VectorDouble &x0=VectorDouble(), int ivar=0, int jvar=0) const
 
double samplingDensityVariance (const Db *db, const VectorDouble &ext, const VectorInt &ndisc, const VectorDouble &angles=VectorDouble(), const VectorDouble &x0=VectorDouble(), int ivar=0, int jvar=0) const
 
double specificVolume (const Db *db, double mean, const VectorDouble &ext, const VectorInt &ndisc, const VectorDouble &angles=VectorDouble(), const VectorDouble &x0=VectorDouble(), int ivar=0, int jvar=0) const
 
double coefficientOfVariation (const Db *db, double volume, double mean, const VectorDouble &ext, const VectorInt &ndisc, const VectorDouble &angles=VectorDouble(), const VectorDouble &x0=VectorDouble(), int ivar=0, int jvar=0) const
 
double specificVolumeFromCoV (Db *db, double cov, double mean, const VectorDouble &ext, const VectorInt &ndisc, const VectorDouble &angles=VectorDouble(), const VectorDouble &x0=VectorDouble(), int ivar=0, int jvar=0) const
 
double evaluateOneGeneric (const CovInternal *covint, const VectorDouble &d1=VectorDouble(), double weight=1., const CovCalcMode *mode=nullptr) const
 
double calculateStDev (Db *db1, int iech1, Db *db2, int iech2, bool verbose=false, double factor=1., const CovCalcMode *mode=nullptr) const
 
void evaluateMatInPlace (const CovInternal *covint, const VectorDouble &d1, MatrixSquare &covtab, bool flag_init=false, double weight=1., const CovCalcMode *mode=nullptr) const
 
VectorDouble evaluateFromDb (Db *db, int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
double evaluateOneIncr (double hh, const VectorDouble &codir=VectorDouble(), int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr) const
 
VectorDouble sample (const VectorDouble &h, const VectorDouble &codir=VectorDouble(), int ivar=0, int jvar=0, const CovCalcMode *mode=nullptr, const CovInternal *covint=nullptr) const
 
VectorDouble sampleUnitary (const VectorDouble &hh, int ivar=0, int jvar=0, VectorDouble codir=VectorDouble(), const CovCalcMode *mode=nullptr) const
 
VectorDouble envelop (const VectorDouble &hh, int ivar=0, int jvar=0, int isign=1, VectorDouble codir=VectorDouble(), const CovCalcMode *mode=nullptr) const
 
int buildVmapOnDbGrid (DbGrid *dbgrid, const NamingConvention &namconv=NamingConvention("VMAP")) const
 
double gofToVario (const Vario *vario, bool verbose=true) const
 
void manage (const Db *db1, const Db *db2) const
 
void load (const SpacePoint &p, bool case1) const
 
virtual void updateCovByMesh (int imesh, bool aniso=true) const
 
virtual double getValue (const EConsElem &econs, int iv1, int iv2) const
 
void makeStationary ()
 
virtual int makeElemNoStat (const EConsElem &econs, int iv1, int iv2, const AFunctional *func=nullptr, const Db *db=nullptr, const String &namecol=String())
 
void createNoStatTab ()
 
void informMeshByMesh (const AMesh *amesh) const
 
void informMeshByApex (const AMesh *amesh) const
 
VectorDouble informCoords (const VectorVectorDouble &coords, const EConsElem &econs, int iv1=0, int iv2=0) const
 
void informDbIn (const Db *dbin) const
 
void informDbOut (const Db *dbout) const
 
virtual void updateCovByPoints (int icas1, int iech1, int icas2, int iech2)
 
int getNDim (int ispace=-1) const
 
void optimizationPreProcessForData (const Db *db1=nullptr) const
 
virtual void setOptimEnabled (bool enabled) const
 
bool checkAndManageNoStatDb (const Db *db, const String &namecol)
 
std::shared_ptr< const DbgetDbNoStat () const
 
const DbgetDbNoStatRaw () const
 
void setNoStatDbIfNecessary (const Db *db)
 
void setNoStatDbIfNecessary (std::shared_ptr< const Db > &db)
 
- Public Member Functions inherited from ASpaceObject
 ASpaceObject (const ASpaceSharedPtr &space=ASpaceSharedPtr())
 
 ASpaceObject (const ASpaceObject &r)
 
ASpaceObjectoperator= (const ASpaceObject &r)
 
virtual ~ASpaceObject ()
 
virtual String toString (const AStringFormat *strfmt=nullptr) const override
 AStringable interface.
 
ASpaceSharedPtr getSpace () const
 Accessor to the current object space context.
 
bool isConsistent () const
 Indicate if I am consistent with my current space context.
 
VectorDouble getUnitaryVector () const
 Return unitary vector for the current space context.
 
bool isConsistent (const ASpaceSharedPtr &space) const
 Indicate if I am consistent with the provided space.
 
virtual bool isConsistent (const ASpace *space) const =0
 
unsigned int getNDim (int ispace=-1) const
 Shortcuts to ASpace methods.
 
const VectorDoublegetOrigin (int ispace=-1) const
 Return the current space context origin coordinates.
 
double getDistance (const SpacePoint &p1, const SpacePoint &p2, int ispace=0) const
 Return the distance between two space points for the current space context.
 
VectorDouble getDistances (const SpacePoint &p1, const SpacePoint &p2) const
 Return all the distances (space composits) between two space points for the current space context.
 
VectorDouble getIncrement (const SpacePoint &p1, const SpacePoint &p2, int ispace=0) const
 Return the increment vector between two space points for the current space context.
 
- Public Member Functions inherited from AStringable
 AStringable ()
 
 AStringable (const AStringable &r)
 
AStringableoperator= (const AStringable &r)
 
virtual ~AStringable ()
 
virtual void display (const AStringFormat *strfmt=nullptr) const final
 
virtual void display (int level) const final
 
- Public Member Functions inherited from ICloneable
 ICloneable ()
 
virtual ~ICloneable ()
 
virtual ICloneableclone () const =0
 

Static Public Member Functions

static void gofDisplay (double gof, bool byValue=true, const VectorDouble &thresholds={2., 5., 10., 100})
 

Constructor & Destructor Documentation

◆ ACov() [1/2]

ACov::ACov ( const CovContext ctxt = CovContext())

◆ ACov() [2/2]

ACov::ACov ( const ACov r)

◆ ~ACov()

ACov::~ACov ( )
virtual

Member Function Documentation

◆ addEvalCovVecRHSInPlace()

int ACov::addEvalCovVecRHSInPlace ( vect  vect,
const VectorInt index1,
const int  iech2,
const KrigOpt krigopt,
SpacePoint pin,
SpacePoint pout,
VectorDouble tabwork,
double  lambda = 1. 
) const
virtual

Reimplemented in CovList, CorAniso, and CovBase.

◆ attachNoStatDb()

void ACov::attachNoStatDb ( const Db db)

◆ buildVmapOnDbGrid()

int ACov::buildVmapOnDbGrid ( DbGrid dbgrid,
const NamingConvention namconv = NamingConvention("VMAP") 
) const

Calculate the variogram map from a Model (presented as Variogram, not Covariance)

Returns
Error return code
Parameters
[in]dbgridGrid structure
[in]namconvNaming convention

◆ calculateStDev()

double ACov::calculateStDev ( Db db1,
int  iech1,
Db db2,
int  iech2,
bool  verbose = false,
double  factor = 1.,
const CovCalcMode mode = nullptr 
) const

Returns the standard deviation at a given increment for a given model between two samples of two Dbs

Parameters
[in]db1First Db
[in]iech1Rank in the first Db
[in]db2Second Db
[in]iech2Rank in the second Db
[in]verboseVerbose flag
[in]factorMultiplicative factor for standard deviation
[in]modeCovCalcMode structure

◆ checkAndManageNoStatDb()

bool ACov::checkAndManageNoStatDb ( const Db db,
const String namecol 
)

◆ coefficientOfVariation()

double ACov::coefficientOfVariation ( const Db db,
double  volume,
double  mean,
const VectorDouble ext,
const VectorInt ndisc,
const VectorDouble angles = VectorDouble(),
const VectorDouble x0 = VectorDouble(),
int  ivar = 0,
int  jvar = 0 
) const

Calculate the Coefficient of Variation

Parameters
dbSet of data points
volumeSpecific production volume
meanValue of the Mean
extTarget Block extension
ndiscVector of discretization
anglesOptional rotation angle for block
x0Optional origin of the Block
ivarRank of the first variable
jvarRank of the second variable
Returns

◆ copyCovContext()

virtual void ACov::copyCovContext ( const CovContext ctxt)
inlinevirtual

Reimplemented in CovList.

◆ createNoStatTab()

void ACov::createNoStatTab ( )

◆ envelop()

VectorDouble ACov::envelop ( const VectorDouble hh,
int  ivar = 0,
int  jvar = 0,
int  isign = 1,
VectorDouble  codir = VectorDouble(),
const CovCalcMode mode = nullptr 
) const

◆ eval()

VectorDouble ACov::eval ( const std::vector< SpacePoint > &  vec_p1,
const std::vector< SpacePoint > &  vec_p2,
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

Functions for evaluating Covariances.

◆ eval0()

double ACov::eval0 ( int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const
virtual

Calculate the covariance between two variables for 0-distance (stationary case)

Reimplemented in CovAniso, CovAnisoList, CovGradientNumerical, CovList, CovLMCAnamorphosis, CovLMCConvolution, and CovLMCTapering.

◆ eval0CovMatBiPointInPlace()

void ACov::eval0CovMatBiPointInPlace ( MatrixSymmetric mat,
const CovCalcMode mode 
) const

Calculate the Matrix of covariance for zero distance

Parameters
matCovariance matrix (Dimension: nvar * nvar)
modeCalculation Options
Remarks
: Matrix 'mat' should be dimensioned and initialized beforehand

◆ eval0Mat()

MatrixSymmetric ACov::eval0Mat ( const CovCalcMode mode = nullptr) const

◆ evalAverageDbToDb()

double ACov::evalAverageDbToDb ( const Db db1,
const Db db2,
int  ivar = 0,
int  jvar = 0,
double  eps = 0.,
int  seed = 434132,
const CovCalcMode mode = nullptr 
) const

Calculate the (weighted) average Covariance between samples of two Dbs, for a pair of variables

Parameters
db1Pointer to the first Db
db2Pointer to the second Db
ivarRank of the first variables
jvarRank of the second variable
epsEpsilon used for randomization in calculation of CVV (optional)
seedSeed for the randomization
modeCovCalcMode structure
Returns

◆ evalAverageIncrToIncr()

double ACov::evalAverageIncrToIncr ( const VectorVectorDouble d1,
const VectorVectorDouble d2,
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

◆ evalAveragePointToDb()

double ACov::evalAveragePointToDb ( const SpacePoint p1,
const Db db2,
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

Calculate the (weighted) average Covariance between a point and a Db for a pair of variables

Parameters
p1Point location
db2Pointer to the second Db
ivarRank of the first variables
jvarRank of the second variable
modeCovCalcMode structure
Returns

◆ evalCov()

double ACov::evalCov ( const SpacePoint p1,
const SpacePoint p2,
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

Calculate the covariance between two variables and two points (general case)

◆ evalCovMat()

MatrixDense ACov::evalCovMat ( const Db db1,
const Db db2 = nullptr,
int  ivar0 = -1,
int  jvar0 = -1,
const VectorInt nbgh1 = VectorInt(),
const VectorInt nbgh2 = VectorInt(),
const CovCalcMode mode = nullptr,
bool  cleanOptim = true 
) const

Establish the covariance matrix between two Dbs Takes into account selection and heterotopy

Returns
Dense matrix containing the covariance matrix
Parameters
[in]db1First Db
[in]db2Second Db (= db1 if absent)
[in]ivar0Rank of the first variable (-1 for all variables)
[in]jvar0Rank of the second variable (-1 for all variables)
[in]nbgh1Vector of indices of active samples in db1 (optional)
[in]nbgh2Vector of indices of active samples in db2 (optional)
[in]modeCovCalcMode structure
[in]cleanOptimWhen True, clean optimization internal when ended
Remarks
If a Db does not contain any Z-variable defined, the covariance
cannot treat possible heterotopy and therefore uses all samples
The returned matrix if dimension to nrows * ncols where
each term is the product of the number of active samples
by the number of samples where the variable is defined
Note
'dbin' and 'dbout' cannot be made 'const' as they can be updated
due to the presence of 'nostat'

◆ evalCovMat0()

MatrixSymmetric ACov::evalCovMat0 ( const Db db,
int  iech,
const KrigOpt krigopt = KrigOpt() 
) const

Functions for evaluating Covariance Matrices either in place or not.

◆ evalCovMat0InPlace()

int ACov::evalCovMat0InPlace ( MatrixSymmetric mat,
const Db db,
int  iech,
const KrigOpt krigopt = KrigOpt() 
) const

◆ evalCovMatInPlace()

int ACov::evalCovMatInPlace ( MatrixDense mat,
const Db db1,
const Db db2 = nullptr,
int  ivar0 = -1,
int  jvar0 = -1,
const VectorInt nbgh1 = VectorInt(),
const VectorInt nbgh2 = VectorInt(),
const CovCalcMode mode = nullptr,
bool  cleanOptim = true 
) const

◆ evalCovMatInPlaceFromIdx()

int ACov::evalCovMatInPlaceFromIdx ( MatrixDense mat,
const Db db1,
const Db db2,
const VectorVectorInt index1,
const VectorVectorInt index2,
const VectorInt nbgh2 = VectorInt(),
const CovCalcMode mode = nullptr,
bool  cleanOptim = true 
) const

◆ evalCovMatRHSInPlaceFromIdx()

int ACov::evalCovMatRHSInPlaceFromIdx ( MatrixDense mat,
const Db db1,
const Db db2,
const VectorVectorInt index1,
const int  iech2 = -1,
const KrigOpt krigopt = KrigOpt(),
bool  cleanOptim = true 
) const

Establish covariance matrix between one Db and one sample of a Target Db

Returns
Dense matrix containing the covariance matrix
Parameters
[in]matMatrix (possibly resized)
[in]db1First Db
[in]db2Second Db
[in]index1Vector of vector indices of active samples in db1
[in]iech2Sample rank within db2
[in]krigoptKrigOpt structure
[in]cleanOptimWhen True, clean optimization internal when ended
Remarks
If a Db does not contain any Z-variable defined, the covariance
cannot treat possible heterotopy and therefore uses all samples
The returned matrix if dimension to nrows * 1 where
each 'nrows' is the number of active samples
by the number of samples where the variable is defined
Note
'dbin' and 'dbout' cannot be made 'const' as they can be updated
due to the presence of 'nostat'

◆ evalCovMatSparse()

MatrixSparse * ACov::evalCovMatSparse ( const Db db1,
const Db db2 = nullptr,
int  ivar0 = -1,
int  jvar0 = -1,
const VectorInt nbgh1 = VectorInt(),
const VectorInt nbgh2 = VectorInt(),
const CovCalcMode mode = nullptr,
bool  cleanOptim = true,
double  eps = EPSILON3 
) const

Establish the covariance matrix between two Dbs where samples are selected by ranks The output is stored in a Sparse Matrix

Returns
Sparse matrix containing the covariance matrix
Parameters
[in]db1First Db
[in]db2Second Db
[in]ivar0Rank of the first variable (-1: all variables)
[in]jvar0Rank of the second variable (-1: all variables)
[in]nbgh1Array giving ranks of selected samples (optional)
[in]nbgh2Array giving ranks of selected samples (optional)
[in]modeCovCalcMode structure
[in]cleanOptimWhen True, clean optimization internal when ended
[in]epsTolerance for discarding a covariance value
Remarks
The covariance matrix (returned) must be freed by calling routine
The covariance matrix is established for the first variable
and returned as a covariance
As the ranks are used, no test is performed on any selection
but only ranks positive or null are considered

◆ evalCovMatSym()

MatrixSymmetric ACov::evalCovMatSym ( const Db db1,
const VectorInt nbgh1 = VectorInt(),
int  ivar0 = -1,
const CovCalcMode mode = nullptr,
bool  cleanOptim = true 
) const

Establish the covariance matrix within a Db Takes into account selection and heterotopy This method takes advantage of calculating covariance between a Db and itself

Returns
Dense matrix containing the covariance matrix
Parameters
[in]db1First Db
[in]ivar0Rank of the first variable (-1 for all variables)
[in]nbgh1Vector of indices of active samples in db1 (optional)
[in]modeCovCalcMode structure
[in]cleanOptimWhen True, clean optimization internal arrays at end
Remarks
If a Db does not contain any Z-variable defined, the covariance
cannot treat possible heterotopy and therefore uses all samples
The returned matrix if dimension to nrows * ncols where
each term is the product of the number of active samples
by the number of samples where the variable is defined

◆ evalCovMatSymInPlace()

int ACov::evalCovMatSymInPlace ( MatrixSymmetric mat,
const Db db1,
const VectorInt nbgh1 = VectorInt(),
int  ivar0 = -1,
const CovCalcMode mode = nullptr,
bool  cleanOptim = true 
) const

◆ evalCovMatSymInPlaceFromIdx()

int ACov::evalCovMatSymInPlaceFromIdx ( MatrixSymmetric mat,
const Db db1,
const VectorVectorInt index1,
const CovCalcMode mode = nullptr,
bool  cleanOptim = true 
) const

◆ evalCovOnSphere()

virtual double ACov::evalCovOnSphere ( double  alpha,
int  degree = 50,
bool  flagScaleDistance = false,
const CovCalcMode mode = nullptr 
) const
inlinevirtual

Reimplemented in CorAniso, and CovAniso.

◆ evalCovVecRHSInPlace()

int ACov::evalCovVecRHSInPlace ( vect  vect,
const Db db2,
const VectorInt index1,
int  iech2,
const KrigOpt krigopt,
SpacePoint pin,
SpacePoint pout,
VectorDouble tabwork,
double  lambda = 1. 
) const

◆ evalCvv()

double ACov::evalCvv ( const VectorDouble ext,
const VectorInt ndisc,
const VectorDouble angles = VectorDouble(),
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

Average covariance over a block

Parameters
extVector of Block extensions
ndiscVector of Block discretization
anglesVector of rotation angles
ivarRank of the first variable
jvarRank of the second variable
modeCovCalcMode structure
Returns

◆ evalCvvM()

MatrixSquare ACov::evalCvvM ( const VectorDouble ext,
const VectorInt ndisc,
const VectorDouble angles = VectorDouble(),
const CovCalcMode mode = nullptr 
) const

◆ evalCvvShift()

double ACov::evalCvvShift ( const VectorDouble ext,
const VectorInt ndisc,
const VectorDouble shift,
const VectorDouble angles = VectorDouble(),
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

Average covariance between a block and the same block shifted

Parameters
extVector of Block extensions
ndiscVector of Block discretization
anglesVector of rotation angles
shiftShift between the two blocks
ivarRank of the first variable
jvarRank of the second variable
modeCovCalcMode structure
Returns

◆ evalCxv() [1/2]

double ACov::evalCxv ( const Db db,
const VectorDouble ext,
const VectorInt ndisc,
const VectorDouble angles = VectorDouble(),
const VectorDouble x0 = VectorDouble(),
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

◆ evalCxv() [2/2]

double ACov::evalCxv ( const SpacePoint p1,
const VectorDouble ext,
const VectorInt ndisc,
const VectorDouble angles = VectorDouble(),
const VectorDouble x0 = VectorDouble(),
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

Average covariance over a block

Parameters
p1Point location
extVector of Block extensions
ndiscVector of Block discretization
anglesVector of rotation angles
x0Vector for origin of block
ivarRank of the first variable
jvarRank of the second variable
modeCovCalcMode structure
Returns

◆ evalCxvM()

MatrixSquare ACov::evalCxvM ( const SpacePoint p1,
const VectorDouble ext,
const VectorInt ndisc,
const VectorDouble angles = VectorDouble(),
const VectorDouble x0 = VectorDouble(),
const CovCalcMode mode = nullptr 
) const

◆ evalIsoIvarIpas()

double ACov::evalIsoIvarIpas ( double  step,
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

Covariance for a given unit global distance (without anisotropy) for a pair of variables and a single step

Parameters
ivarRank of the first variable
jvarRank of the second variable
stepStep value
modeCovCalcMode structure
Returns

TODO : Not true whatever the space

◆ evalIsoIvarNlag()

VectorDouble ACov::evalIsoIvarNlag ( const VectorDouble vec_step,
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

Covariance for a given unit global distance (without anisotropy) for a pair of variables and a set of steps

Parameters
ivar
jvar
vec_step
mode
Returns

◆ evalIsoNvarIpas()

MatrixSquare ACov::evalIsoNvarIpas ( double  step,
const CovCalcMode mode = nullptr 
) const

Covariance for a given unit global distance (without anisotropy) for a set of variables and a single step

Parameters
stepStep value
modeCovCalcMode structure
Returns

◆ evalIvarIpas()

double ACov::evalIvarIpas ( double  step,
const VectorDouble dir = VectorDouble(),
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

Covariance from a given point (center) in a given direction (dir *step)

Parameters
ivarRank of the first variable
jvarRank of the second variable
stepStep value
dirDirection definition
modeCovCalcMode structure
Returns

◆ evalIvarIpasIncr()

double ACov::evalIvarIpasIncr ( const VectorDouble dincr,
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

◆ evalIvarNlag()

VectorDouble ACov::evalIvarNlag ( const VectorDouble vec_step,
const VectorDouble dir = VectorDouble(),
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

Covariance vector from a given point (center) in a given direction (dir * steps) for a pair of variables and a set of steps

Parameters
ivarRank of the first variable
jvarRank of the second variable
vec_stepVector of step values
dirDirection definition
modeCovCalcMode structure
Returns

◆ evalNvarIpas()

MatrixSquare ACov::evalNvarIpas ( double  step,
const VectorDouble dir = VectorDouble(),
const CovCalcMode mode = nullptr 
) const

Covariance Matrix from a given point (center) in a given direction (dir * step) for a set of variables and a given step

Parameters
stepStep value
dirDirection definition
modeCovCalcMode structure
Returns

◆ evalNvarIpasIncr()

MatrixSquare ACov::evalNvarIpasIncr ( const VectorDouble dincr,
const CovCalcMode mode = nullptr 
) const

◆ evalPointToDb()

void ACov::evalPointToDb ( VectorDouble values,
const SpacePoint p1,
const Db db2,
int  ivar = 0,
int  jvar = 0,
bool  useSel = true,
const VectorInt nbgh2 = VectorInt(),
const CovCalcMode mode = nullptr 
) const

Calculate the Covariance vector between a Point and all the samples of a Db, for a pair of variables

Parameters
valuesArray of returned values (possible resized)
p1Point location
db2Pointer to the second Db
ivarRank of the first variables
jvarRank of the second variable
useSelWhen TRUE, the returned vector is reduced to active samples Otherwise, returns TEST for masked samples
nbgh2Vector of indices of active samples in db2 (optional)
modeCovCalcMode structure

◆ evalPointToDbAsSP()

void ACov::evalPointToDbAsSP ( VectorDouble values,
const std::vector< SpacePoint > &  p1s,
const SpacePoint p2,
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

◆ evalSpectrum()

virtual double ACov::evalSpectrum ( const VectorDouble freq,
int  ivar,
int  jvar 
) const
inlinevirtual

Reimplemented in CorAniso, and CovAniso.

◆ evalSpectrumOnSphere()

virtual VectorDouble ACov::evalSpectrumOnSphere ( int  n,
bool  flagNormDistance = false,
bool  flagCumul = false 
) const
inlinevirtual

Reimplemented in CorAniso, and CovAniso.

◆ evaluateFromDb()

VectorDouble ACov::evaluateFromDb ( Db db,
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

Evaluate the model on a Db

Parameters
[in]dbDb structure
[in]ivarRank of the first variable
[in]jvarRank of the second variable
[in]modeCovCalcMode structure

◆ evaluateMatInPlace()

void ACov::evaluateMatInPlace ( const CovInternal covint,
const VectorDouble d1,
MatrixSquare covtab,
bool  flag_init = false,
double  weight = 1.,
const CovCalcMode mode = nullptr 
) const

Returns the covariances for an increment This is the generic internal function It can be called for stationary or non-stationary case

Parameters
[in]covintInternal structure for non-stationarityAddress for the next term after the drift or NULL (for stationary case)
[in]modeCovCalcMode structure
[in]flag_initInitialize the array beforehand
[in]weightMultiplicative weight
[in]d1Distance vector
[out]covtabCovariance array

◆ evaluateOneGeneric()

double ACov::evaluateOneGeneric ( const CovInternal covint,
const VectorDouble d1 = VectorDouble(),
double  weight = 1.,
const CovCalcMode mode = nullptr 
) const

Returns the covariance for an increment This is the generic internal function It can be called for stationary or non-stationary case

Parameters
[in]covintInternal structure for non-stationarityAddress for the next term after the drift or NULL (for stationary case)
[in]modeCovCalcMode structure
[in]weightMultiplicative weight
[in]d1Distance vector

◆ evaluateOneIncr()

double ACov::evaluateOneIncr ( double  hh,
const VectorDouble codir = VectorDouble(),
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr 
) const

Calculate the value of the model for a set of distances

Returns
The model value
Parameters
[in]ivarRank of the first variable
[in]jvarRank of the second variable
[in]modeCovCalcMode structure
[in]codirArray giving the direction coefficients (optional)
[in]hhVector of increments

◆ extensionVariance()

double ACov::extensionVariance ( const Db db,
const VectorDouble ext,
const VectorInt ndisc,
const VectorDouble angles = VectorDouble(),
const VectorDouble x0 = VectorDouble(),
int  ivar = 0,
int  jvar = 0 
) const

Variance of Extension of a set of points and the block

Parameters
dbReference Data Base
extVector giving the extensions of the target block
ndiscVector giving the discretization
anglesVector for the rotation angles of the block (optional)
x0Optional origin of the Block
ivarRank of the first variable
jvarRank of the second variable
Returns

◆ getContext()

const CovContext & ACov::getContext ( ) const
inline

◆ getContextCopy()

CovContext ACov::getContextCopy ( ) const
inline

◆ getDbNoStat()

std::shared_ptr< const Db > ACov::getDbNoStat ( ) const

◆ getDbNoStatRaw()

const Db * ACov::getDbNoStatRaw ( ) const

◆ getNDim()

int ACov::getNDim ( int  ispace = -1) const
inline

◆ getNVar()

virtual int ACov::getNVar ( ) const
inlinevirtual

ACov Interface.

Reimplemented in CorAniso, CorGneiting, CorMatern, CovAnisoList, CovBase, and CovList.

◆ getValue()

virtual double ACov::getValue ( const EConsElem &  econs,
int  iv1,
int  iv2 
) const
inlinevirtual

Reimplemented in CorAniso, and CovBase.

◆ gofDisplay()

void ACov::gofDisplay ( double  gof,
bool  byValue = true,
const VectorDouble thresholds = {2., 5., 10., 100} 
)
static

Printout of statement concerning the Quality of the GOF

Parameters
gofValue of the Gof
byValuetrue: display GOF value; false: print its quality level
thresholdsVector giving the Quality thresholds

◆ gofToVario()

double ACov::gofToVario ( const Vario vario,
bool  verbose = true 
) const

Evaluate the Goodness-of_fit of the Model on the Experimental Variogram It is expressed as the average departure between Model and Variogram scaled to the variance. As this variance may be poorly calculated (< gmax / 5), it may be replaced by the largest value (gmax) divided by 2 (highly non_stationary cases).

Parameters
varioExperimental variogram
verboseVerbose flag
Returns
Value for the Goodness-of_fit (as percentage of the total sill)

◆ informCoords()

VectorDouble ACov::informCoords ( const VectorVectorDouble coords,
const EConsElem &  econs,
int  iv1 = 0,
int  iv2 = 0 
) const

◆ informDbIn()

void ACov::informDbIn ( const Db dbin) const

◆ informDbOut()

void ACov::informDbOut ( const Db dbout) const

◆ informMeshByApex()

void ACov::informMeshByApex ( const AMesh amesh) const

◆ informMeshByMesh()

void ACov::informMeshByMesh ( const AMesh amesh) const

◆ initFromContext()

void ACov::initFromContext ( )
inline

◆ isIndexable()

virtual bool ACov::isIndexable ( ) const
inlinevirtual

Reimplemented in CovAnisoList, and CovList.

◆ isNoStat()

bool ACov::isNoStat ( ) const
inline

◆ isOptimEnabled()

virtual bool ACov::isOptimEnabled ( ) const
inlinevirtual

Functions linked to Optimization during Covariance calculations.

◆ load()

void ACov::load ( const SpacePoint p,
bool  case1 
) const

◆ loadInfoValues()

virtual void ACov::loadInfoValues ( )
inlinevirtual

Reimplemented in CovBase.

◆ makeElemNoStat()

int ACov::makeElemNoStat ( const EConsElem &  econs,
int  iv1,
int  iv2,
const AFunctional func = nullptr,
const Db db = nullptr,
const String namecol = String() 
)
virtual

Reimplemented in CovBase, and CovList.

◆ makeStationary()

void ACov::makeStationary ( )

◆ manage()

void ACov::manage ( const Db db1,
const Db db2 
) const
inline

◆ operator=()

ACov & ACov::operator= ( const ACov r)

◆ optimizationLoadInPlace()

SpacePoint & ACov::optimizationLoadInPlace ( int  iech,
int  mode,
int  rank 
) const

◆ optimizationPostProcess()

void ACov::optimizationPostProcess ( ) const

◆ optimizationPreProcess()

void ACov::optimizationPreProcess ( int  mode,
const std::vector< SpacePoint > &  ps 
) const

◆ optimizationPreProcessForData()

void ACov::optimizationPreProcessForData ( const Db db1 = nullptr) const

◆ optimizationSetTarget()

void ACov::optimizationSetTarget ( SpacePoint pt) const

◆ sample()

VectorDouble ACov::sample ( const VectorDouble h,
const VectorDouble codir = VectorDouble(),
int  ivar = 0,
int  jvar = 0,
const CovCalcMode mode = nullptr,
const CovInternal covint = nullptr 
) const

Calculate the value of the model for a set of distances

Returns
Array containing the model values
Parameters
[in]ivarRank of the first variable
[in]jvarRank of the second variable
[in]codirArray giving the direction coefficients (optional)
[in]hVector of increments
[in]modeCovCalcMode structure
[in]covintNon-stationary parameters

◆ sampleUnitary()

VectorDouble ACov::sampleUnitary ( const VectorDouble hh,
int  ivar = 0,
int  jvar = 0,
VectorDouble  codir = VectorDouble(),
const CovCalcMode mode = nullptr 
) const

Returns the value of the normalized covariance (by the variance/covariance value) for a given pair of variables

Parameters
hhVector of distances
ivarRank of the first variable
jvarRank of the second variable
codirDirection coefficients
modeCovCalcMode structure
Returns

◆ samplingDensityVariance()

double ACov::samplingDensityVariance ( const Db db,
const VectorDouble ext,
const VectorInt ndisc,
const VectorDouble angles = VectorDouble(),
const VectorDouble x0 = VectorDouble(),
int  ivar = 0,
int  jvar = 0 
) const

Calculate the Sampling Density Variance

Parameters
dbSet of data points
extBlock extension
ndiscDiscretization
anglesOptional rotation angles for the Block
x0Optional origin of the block
ivarRank of the first variable
jvarRank of the second variable
Returns

◆ setContext()

void ACov::setContext ( const CovContext ctxt)

◆ setNoStatDbIfNecessary() [1/2]

void ACov::setNoStatDbIfNecessary ( const Db db)

◆ setNoStatDbIfNecessary() [2/2]

void ACov::setNoStatDbIfNecessary ( std::shared_ptr< const Db > &  db)

◆ setOptimEnabled()

virtual void ACov::setOptimEnabled ( bool  enabled) const
inlinevirtual

Reimplemented in CovBase, and CovList.

◆ specificVolume()

double ACov::specificVolume ( const Db db,
double  mean,
const VectorDouble ext,
const VectorInt ndisc,
const VectorDouble angles = VectorDouble(),
const VectorDouble x0 = VectorDouble(),
int  ivar = 0,
int  jvar = 0 
) const

Calculate the Specific Volume

Parameters
dbSet of data points
meanValue of the Mean
extTarget Block extension
ndiscVector of discretization
anglesOptional rotation angle for block
x0Optional origin of the Block
ivarRank of the first variable
jvarRank of the second variable
Returns

◆ specificVolumeFromCoV()

double ACov::specificVolumeFromCoV ( Db db,
double  cov,
double  mean,
const VectorDouble ext,
const VectorInt ndisc,
const VectorDouble angles = VectorDouble(),
const VectorDouble x0 = VectorDouble(),
int  ivar = 0,
int  jvar = 0 
) const

Derive the Specific volume for a given CoV

Parameters
dbSet of data points
covTarget Coefficient of Variation
meanValue of the Mean
extTarget Block extension
ndiscVector of discretization
anglesOptional rotation angle for block
x0Optional origin of the Block
ivarRank of the first variable
jvarRank of the second variable
Returns

◆ updateCovByMesh()

virtual void ACov::updateCovByMesh ( int  imesh,
bool  aniso = true 
) const
inlinevirtual

Reimplemented in CorAniso, and CovBase.

◆ updateCovByPoints() [1/2]

virtual void ACov::updateCovByPoints ( int  icas1,
int  iech1,
int  icas2,
int  iech2 
)
inlinevirtual

Reimplemented in CorAniso.

◆ updateCovByPoints() [2/2]

virtual void ACov::updateCovByPoints ( int  icas1,
int  iech1,
int  icas2,
int  iech2 
) const
inlinevirtual

Reimplemented in CovBase, and CovList.

◆ updateFromContext()

void ACov::updateFromContext ( )
inline

The documentation for this class was generated from the following files: