Geophysical Inversion and Modelling Library  v1.3.1-2-g7599abf9
GIMLI::RInversion Class Reference
+ Inheritance diagram for GIMLI::RInversion:
+ Collaboration diagram for GIMLI::RInversion:

Public Types

typedef RVector Vec
typedef RVector ModelVector
- Public Types inherited from GIMLI::InversionBase< double >
typedef Vector< double > ModelVector

Public Member Functions

 RInversion (bool verbose=false, bool dosave=false)
 RInversion (ModellingBase &forward, bool verbose=false, bool dosave=false)
 RInversion (const Vec &data, ModellingBase &forward, bool verbose=false, bool dosave=false)
 RInversion (const Vec &data, ModellingBase &forward, Trans< Vec > &transData, bool verbose=true, bool dosave=false)
 RInversion (const Vec &data, ModellingBase &forward, Trans< Vec > &transData, Trans< Vec > &transModel, bool verbose=true, bool dosave=false)
virtual ~RInversion ()
void setData (const Vec &data)
const Vecdata () const
void setRelativeError (const Vec &e)
void setRelativeError (double relerr)
void setAbsoluteError (const Vec &abserr)
void setAbsoluteError (double abserr)
void setError (const Vec &err, bool isRelative=true)
void setError (double err, bool isRelative=true)
const Vecerror () const
void checkError ()
Vec errorDefault_ ()
void setForwardOperator (ModellingBase &forward)
ModellingBaseforwardOperator ()
ModellingBasefop ()
void setTransData (Trans< Vec > &tD)
Trans< Vec > & transData ()
void setTransModel (Trans< Vec > &tM)
Trans< Vec > & transModel ()
uint constraintsCount () const
uint dataCount () const
void setVerbose (bool verbose)
bool verbose () const
void setDoSave (bool d)
bool doSave () const
void saveModelHistory (bool doSaveModelHistory)
void setLineSearch (bool linesearch)
bool lineSearch () const
void setBlockyModel (bool isBlocky)
bool blockyModel () const
void setRobustData (bool isRobust)
bool robustData () const
void setLambda (double lambda)
double lambda () const
double getLambda () const
void setLambdaFactor (double lambdaFactor)
double lambdaFactor () const
void setLambdaMinimum (double l)
double lambdaMinimum () const
void setLocalRegularization (bool localReg)
bool localRegularization () const
void setOptimizeLambda (bool opt)
bool optimizeLambda () const
void setMaxIter (int maxiter)
int maxIter () const
void setMaxCGLSIter (int iter)
int maxCGLSIter () const
void setCGLSTolerance (double tol)
double maxCGLSTolerance () const
uint iter () const
bool isRunning () const
void abort ()
void stopAtChi1 (bool stopAtChi1)
void setDeltaPhiAbortPercent (double dPhi)
double deltaPhiAbortPercent () const
void setMarquardtScheme (double lambdaFactor=0.8)
void checkTransFunctions ()
void setRecalcJacobian (bool recalc)
bool recalcJacobian () const
void setBroydenUpdate (bool broydenUpdate)
bool broydenUpdate () const
void setModel (const Vec &model)
const ModelVectormodel () const
void setReferenceModel (const Vec &model)
void setResponse (const Vec &response)
void setConstraintsH (const Vec &constraintsH)
const Vecresponse () const
Vec getIRLS () const
void setCWeight (const Vec &cWeight)
const VeccWeight () const
void setMWeight (const Vec &mweight)
const VecmWeight () const
void robustWeighting ()
void constrainBlocky ()
void checkConstraints ()
void checkJacobian (bool force=false)
void echoStatus () const
void echoStatus (const Vec &response, const Vec &model, const std::string &xtra="") const
Vec modelCellResolution (int iModel)
Vec modelCellResolution (const RVector3 &pos)
RMatrix modelResolutionMatrix ()
RVector roughness (const RVector &model) const
RVector roughness () const
RVector pureRoughness (const RVector &model) const
double getPhiD (const Vec &response) const
double getPhiM (const Vec &model) const
double getPhi (const Vec &model, const Vec &response) const
double getPhiD () const
double getPhiM () const
double getPhi () const
double getChi2 () const
double getChi2 (const Vec &response) const
double chi2 () const
double absrms () const
double relrms () const
double linesearch (const Vec &modelNew, const Vec &responseNew) const
double linesearchQuad (const Vec &modelNew, const Vec &responseNew, const Vec &modelQuad, const Vec &responseQuad, double tauquad) const
const std::vector< RVector > & modelHistory () const
Vec invSubStep (const Vec &rhs)
Vec optLambda (const Vec &deltaData, const Vec &deltaModel0)
bool oneStep ()
 !! h-variante More...
const Vecinvert (const Vec &data)
const Vecstart ()
const Vecrun ()
Vec runChi1 (double acc=0.01, int maxiter=50)
const RVectormodelWeight () const
const RVectormodelRef () const
const RVectordataWeight () const
const RVectordeltaDataIter () const
const RVectordeltaModelIter () const
void reset ()
- Public Member Functions inherited from GIMLI::InversionBase< double >
 InversionBase ()
virtual ~InversionBase ()
virtual const ModelVectormodel () const=0
virtual const ModelVectorcWeight () const=0
virtual uint iter () const=0
virtual double chi2 () const=0
virtual bool isRunning () const=0

Protected Member Functions

void init_ ()

Protected Attributes

Vec data_
Trans< Vec > * tD_
Trans< Vec > * tM_
Trans< Vec > * transModelDefault_
Trans< Vec > * transDataDefault_
bool verbose_
bool dosave_
bool saveModelHistory_
Vec error_
Vec response_
Vec model_
Vec modelRef_
Vec constraintsH_
Vec constraintWeights_
Vec modelWeight_
Vec dataWeight_
Vec deltaDataIter_
Vec deltaModelIter_
int maxiter_
int iter_
int maxCGLSIter_
double lambda_
double lambdaFactor_
double lambdaMin_
double dPhiAbortPercent_
double CGLStol_
bool isBlocky_
bool isRobust_
bool isRunning_
bool useLinesearch_
bool optimizeLambda_
bool abort_
bool stopAtChi1_
bool doBroydenUpdate_
bool localRegularization_
bool haveReferenceModel_
bool recalcJacobian_
bool jacobiNeedRecalc_
bool activateFillConstraintWeights_
bool fixError_
std::vector< RVectormodelHist_

Detailed Description

Inversion template using a Generalized Minimization Approach, atm fixed to Gauss-Newton solver Inversion(bool verbose, bool dosave Inversion(RVector data, FOP f, bool verbose, bool dosave Inversion(RVector data, FOP f, transData, transModel, bool verbose, bool dosave

Constructor & Destructor Documentation

◆ RInversion() [1/5]

GIMLI::RInversion::RInversion ( bool  verbose = false,
bool  dosave = false 

Minimal constructor. verbose – gives some output and status information about the inversion progress; dosave – advanced debuging, saves alot of stuff and temporary data

◆ RInversion() [2/5]

GIMLI::RInversion::RInversion ( ModellingBase forward,
bool  verbose = false,
bool  dosave = false 

Contructor for forward operator only

◆ RInversion() [3/5]

GIMLI::RInversion::RInversion ( const Vec data,
ModellingBase forward,
bool  verbose = false,
bool  dosave = false 

Usual constructor. data – vector of the given data; forward – the associated forward operator

◆ RInversion() [4/5]

GIMLI::RInversion::RInversion ( const Vec data,
ModellingBase forward,
Trans< Vec > &  transData,
bool  verbose = true,
bool  dosave = false 

Intermediate constructor including transforms and data transformation function

◆ RInversion() [5/5]

GIMLI::RInversion::RInversion ( const Vec data,
ModellingBase forward,
Trans< Vec > &  transData,
Trans< Vec > &  transModel,
bool  verbose = true,
bool  dosave = false 

Full constructor including transforms. transData – data transformation function; transModel – model transformation function

◆ ~RInversion()

virtual GIMLI::RInversion::~RInversion ( )

Destructor. Frees allocated memory

Member Function Documentation

◆ abort()

void GIMLI::RInversion::abort ( )

Force abort at next iteration

Implements GIMLI::InversionBase< double >.

◆ absrms()

double GIMLI::RInversion::absrms ( ) const

Return last absolute RMS misfit

◆ checkConstraints()

void GIMLI::RInversion::checkConstraints ( )

Create constraints, check and compare size of constraint matrix with model/boundary control

References GIMLI::Vector< ValueType >::resize().

◆ checkError()

void GIMLI::RInversion::checkError ( )

Validate the data error. This check will be called at every new run or change of the error; Check size and content unequal to zero. If problems found, warn and set the data error to errorDefault_()

◆ checkJacobian()

void GIMLI::RInversion::checkJacobian ( bool  force = false)

Check size of Jacobian matrix against data and model number

◆ checkTransFunctions()

void GIMLI::RInversion::checkTransFunctions ( )

Check if transfunctions are valid, set default if no transfunctions are given or override given functions by regionManager.transfunctions if available

◆ chi2()

double GIMLI::RInversion::chi2 ( ) const

Return last chi-squared data misfit

◆ constrainBlocky()

void GIMLI::RInversion::constrainBlocky ( )

Apply blocky model constraints (see also setBlockyModel)

◆ constraintsCount()

uint GIMLI::RInversion::constraintsCount ( ) const

Return number of constraints

◆ cWeight()

const Vec& GIMLI::RInversion::cWeight ( ) const

Return reference to the current constraints weight vector

◆ data()

const Vec& GIMLI::RInversion::data ( ) const

Get data vector

◆ dataCount()

uint GIMLI::RInversion::dataCount ( ) const

Return number of data

◆ echoStatus() [1/2]

void GIMLI::RInversion::echoStatus ( ) const

Shortcut for echoStatus(response_, model_).

◆ echoStatus() [2/2]

void GIMLI::RInversion::echoStatus ( const Vec response,
const Vec model,
const std::string &  xtra = "" 
) const

Echo PhiD/PhiM/Phi and Chi^2 for given model and response

◆ error()

const Vec& GIMLI::RInversion::error ( ) const

Return the used data error

◆ errorDefault_()

Vec GIMLI::RInversion::errorDefault_ ( )

Return a copy of the default data error array, [ Relative data error of 1% ]

◆ fop()

ModellingBase* GIMLI::RInversion::fop ( )

Return forward operator. Shortcut for forwardOperator()

◆ forwardOperator()

ModellingBase* GIMLI::RInversion::forwardOperator ( )

Return forward operator.

Implements GIMLI::InversionBase< double >.

◆ getChi2() [1/2]

double GIMLI::RInversion::getChi2 ( ) const

Return the chi-squared data misfit DEPRECATED wrong nameing scheme

◆ getChi2() [2/2]

double GIMLI::RInversion::getChi2 ( const Vec response) const

Return the chi-squared data misfit for given forward response. DEPRECATED wrong nameing scheme

◆ getIRLS()

Vec GIMLI::RInversion::getIRLS ( ) const

Return IRLS function of roughness vector

◆ getPhi() [1/2]

double GIMLI::RInversion::getPhi ( ) const

Shortcut for getPhi(model_, response_) necessary ? DEPRECATED wrong nameing scheme

◆ getPhi() [2/2]

double GIMLI::RInversion::getPhi ( const Vec model,
const Vec response 
) const

Return total objective function (data OF plus lambda times model OF) DEPRECATED wrong nameing scheme

◆ getPhiD() [1/2]

double GIMLI::RInversion::getPhiD ( ) const

Shortcut for getPhiD(response_) necessary ? DEPRECATED wrong nameing scheme

◆ getPhiD() [2/2]

double GIMLI::RInversion::getPhiD ( const Vec response) const

Return data objective function (sum of squared data-weighted misfit)

◆ getPhiM() [1/2]

double GIMLI::RInversion::getPhiM ( ) const

Shortcut for getPhiM(model_) necessary ? DEPRECATED wrong nameing scheme

◆ getPhiM() [2/2]

double GIMLI::RInversion::getPhiM ( const Vec model) const

Return model objective function (squared model roughness)

◆ init_()

void GIMLI::RInversion::init_ ( )

Internal initialization function, which is called from constructor. Set default paramater and allocate required memory

◆ invert()

const RVector & GIMLI::RInversion::invert ( const Vec data)

Start the inversion with specific data.

◆ invSubStep()

Vec GIMLI::RInversion::invSubStep ( const Vec rhs)

Compute model update by solving one inverse sub-step

rhsThe right-hand-side vector of the system of equation

!! h-variante

◆ isRunning()

bool GIMLI::RInversion::isRunning ( ) const

Return true if inversion is running

◆ iter()

uint GIMLI::RInversion::iter ( ) const

Return curent iteration number

◆ lambdaMinimum()

double GIMLI::RInversion::lambdaMinimum ( ) const

Return the minimum lambda value that can be reached if lambda factor is set.

◆ linesearch()

double GIMLI::RInversion::linesearch ( const Vec modelNew,
const Vec responseNew 
) const

Start with linear interpolation, followed by quadratic fit if linesearch parameter tau is lower than 0.03. Tries to return values between 0.03 and 1

rigorous minimization of the total objective function

parabolic line search using step 0, 1 and tauquad to fit parabola

too large

negative or nearly zero (stucked) -> use small step instead

only if not nearly one (which saves a forward run

◆ linesearchQuad()

double GIMLI::RInversion::linesearchQuad ( const Vec modelNew,
const Vec responseNew,
const Vec modelQuad,
const Vec responseQuad,
double  tauquad 
) const

Compute objective function for old (tau=0), new (tau=1) and another model

◆ localRegularization()

bool GIMLI::RInversion::localRegularization ( ) const

Return whether regularization is global or local (e.g. Marquardt method)

◆ model()

const ModelVector& GIMLI::RInversion::model ( ) const

Return a const reference to the current model vector

◆ modelCellResolution() [1/2]

Vec GIMLI::RInversion::modelCellResolution ( const RVector3 pos)

Compute cell resolution for closest cell to given position

References GIMLI::Matrix< ValueType >::push_back().

◆ modelCellResolution() [2/2]

Vec GIMLI::RInversion::modelCellResolution ( int  iModel)

Compute model cell resolution (specific column of resolution matrix) by an LSCG solver (Guenther, 2004)

◆ modelHistory()

const std::vector< RVector >& GIMLI::RInversion::modelHistory ( ) const

Return the single models for each iteration. For debugging.

◆ modelResolutionMatrix()

RMatrix GIMLI::RInversion::modelResolutionMatrix ( )

Compute whole resolution matrix by single cell resolutions

◆ mWeight()

const Vec& GIMLI::RInversion::mWeight ( ) const

Return reference to the current model weight vector

References GIMLI::getIRLSWeights().

◆ oneStep()

bool GIMLI::RInversion::oneStep ( )

!! h-variante

One iteration step. Return true if the step can be calculated successfully else false is returned.

full step possible;

normal line search parameter between 0.03 and 0.94

** temporary stuff

References roughness().

◆ optLambda()

RVector GIMLI::RInversion::optLambda ( const Vec deltaData,
const Vec deltaModel0 

Optimization of regularization parameter by L-curve

◆ pureRoughness()

RVector GIMLI::RInversion::pureRoughness ( const RVector model) const

Return (C * m) , i.e. the pure (unweighted) roughness

◆ relrms()

double GIMLI::RInversion::relrms ( ) const

Return last relative RMS misfit

◆ reset()

void GIMLI::RInversion::reset ( )

Resets this inversion to the given startmodel.

◆ response()

const Vec& GIMLI::RInversion::response ( ) const

Return a reference to the current response vector

◆ robustWeighting()

void GIMLI::RInversion::robustWeighting ( )

Call robust date reweighting (see also setRobustData)

◆ roughness() [1/2]

RVector GIMLI::RInversion::roughness ( ) const

Shortcut for roughness for the current model vector

Referenced by oneStep().

◆ roughness() [2/2]

RVector GIMLI::RInversion::roughness ( const RVector model) const

Return (C * m * m_w) * c_w

◆ run()

const RVector & GIMLI::RInversion::run ( )

Start the inversion procedure from the last model and return the final model vector.

Run inversion with current model.

calculation of initial modelresponse

() clear the model history

validate and rebuild the data error if necessary

validate and rebuild the constraint matrix if necessary

compute roughness constraint and correct it for inter-region constraints

validate and rebuild the jacobian if necessary

Implements GIMLI::InversionBase< double >.

◆ runChi1()

Vec GIMLI::RInversion::runChi1 ( double  acc = 0.01,
int  maxiter = 50 

Specialized run function that tries to reach a datafit chi^2=1 by varying the regularization paramater lambda

◆ saveModelHistory()

void GIMLI::RInversion::saveModelHistory ( bool  doSaveModelHistory)

Set and get verbose behaviour

◆ setAbsoluteError() [1/2]

void GIMLI::RInversion::setAbsoluteError ( const Vec abserr)

Set the absolute data error to the vector abserr

◆ setAbsoluteError() [2/2]

void GIMLI::RInversion::setAbsoluteError ( double  abserr)

Set absolute data error to the scalar value abserr

◆ setBlockyModel()

void GIMLI::RInversion::setBlockyModel ( bool  isBlocky)

Set and get blocky model behaviour (by L1 reweighting of constraints)

◆ setBroydenUpdate()

void GIMLI::RInversion::setBroydenUpdate ( bool  broydenUpdate)

Enable/disable broyden update (which enforces scaling by transform function derivatives

◆ setCGLSTolerance()

void GIMLI::RInversion::setCGLSTolerance ( double  tol)

Set and get abort tolerance for cgls solver. -1 means default scaling [default]

◆ setConstraintsH()

void GIMLI::RInversion::setConstraintsH ( const Vec constraintsH)

Set constraint right-hand side by hand (very special cases, so be careful)

◆ setCWeight()

void GIMLI::RInversion::setCWeight ( const Vec cWeight)

Set the constraint weight (boundary control) vector

◆ setData()

void GIMLI::RInversion::setData ( const Vec data)

Set data vector

Implements GIMLI::InversionBase< double >.

◆ setDeltaPhiAbortPercent()

void GIMLI::RInversion::setDeltaPhiAbortPercent ( double  dPhi)

Set and get relative objective function decrease

◆ setDoSave()

void GIMLI::RInversion::setDoSave ( bool  d)

Set debug output

◆ setForwardOperator()

void GIMLI::RInversion::setForwardOperator ( ModellingBase forward)

Set the forward operator and the model-transform-function derived from fop-regionManager if not set before.

why is this so strictly necessary???

Always use a region manager

Implements GIMLI::InversionBase< double >.

◆ setLambda()

void GIMLI::RInversion::setLambda ( double  lambda)

Set and get regularization parameter and its change over iteration

Implements GIMLI::InversionBase< double >.

◆ setLambdaMinimum()

void GIMLI::RInversion::setLambdaMinimum ( double  l)

Set the minimum lambda value that can be reached if lambda factor is set. Default is 1.

◆ setLineSearch()

void GIMLI::RInversion::setLineSearch ( bool  linesearch)

Set and get line search

◆ setLocalRegularization()

void GIMLI::RInversion::setLocalRegularization ( bool  localReg)

Set whether regularization is global or local (e.g. Marquardt method)

◆ setMarquardtScheme()

void GIMLI::RInversion::setMarquardtScheme ( double  lambdaFactor = 0.8)

Marquardt scheme (local damping with decreasing regularization strength

no contribution of regularization to objective function

lambda is decreased towards zero

let the solution fully converge (important for statistics)

◆ setMaxCGLSIter()

void GIMLI::RInversion::setMaxCGLSIter ( int  iter)

Set and get maximum iteration number for inverse sub step

◆ setMaxIter()

void GIMLI::RInversion::setMaxIter ( int  maxiter)

Set and get maximum iteration number

Implements GIMLI::InversionBase< double >.

◆ setModel()

void GIMLI::RInversion::setModel ( const Vec model)

Set model vector . If you call run() the inversion starts with this model, otherwise it will start with fop.startModel().

Implements GIMLI::InversionBase< double >.

◆ setMWeight()

void GIMLI::RInversion::setMWeight ( const Vec mweight)

Set the model weight vector

◆ setOptimizeLambda()

void GIMLI::RInversion::setOptimizeLambda ( bool  opt)

Set and get whether lambda is being optimized by Lcurve

◆ setRecalcJacobian()

void GIMLI::RInversion::setRecalcJacobian ( bool  recalc)

Define and find out whether Jacobian is recalculated in each iteration

◆ setReferenceModel()

void GIMLI::RInversion::setReferenceModel ( const Vec model)

Set reference model vector

Implements GIMLI::InversionBase< double >.

References GIMLI::getIRLSWeights().

◆ setRelativeError() [1/2]

void GIMLI::RInversion::setRelativeError ( const Vec e)

Set the relative data error to the vector err.

◆ setRelativeError() [2/2]

void GIMLI::RInversion::setRelativeError ( double  relerr)

Set relative data error to scalar error value relerr. If you force relerr == 0 here. Data will not corrected or weighted.

◆ setResponse()

void GIMLI::RInversion::setResponse ( const Vec response)

Set current model response (e.g. to avoid time-consuming forward calculation, but be careful)

◆ setRobustData()

void GIMLI::RInversion::setRobustData ( bool  isRobust)

Set and get robust data weighting by L1 reweighting scheme

◆ setTransData()

void GIMLI::RInversion::setTransData ( Trans< Vec > &  tD)

Set and get data transform. WARNING! we just store a reference to the trans .. U have to ensure that the trans is and stayes valid. TODO change this!!!

Implements GIMLI::InversionBase< double >.

◆ setTransModel()

void GIMLI::RInversion::setTransModel ( Trans< Vec > &  tM)

Set and get Model transform WARNING! we just store a reference to the trans .. U have to ensure that the trans is and stayes valid. TODO change this!!!

Implements GIMLI::InversionBase< double >.

◆ setVerbose()

void GIMLI::RInversion::setVerbose ( bool  verbose)

Set verbose output

◆ start()

const RVector & GIMLI::RInversion::start ( )

Start the inversion procedure from starting model.

Start inversion from starting model.

◆ stopAtChi1()

void GIMLI::RInversion::stopAtChi1 ( bool  stopAtChi1)

Define whether minimization stops at chi^2=1

Member Data Documentation

◆ fixError_

bool GIMLI::RInversion::fixError_

Set this to zero if u want to use absolute errors == zero

◆ modelHist_

std::vector< RVector > GIMLI::RInversion::modelHist_

Hold old models, for debuging