#include <iostream>#include <iomanip>#include <fstream>#include <algorithm>#include <numeric>#include <sstream>#include <string>#include <boost/bind.hpp>#include <boost/shared_ptr.hpp>#include <boost/filesystem.hpp>#include <boost/date_time/posix_time/posix_time.hpp>#include "GA.h"#include "Iso1DMTObjective.h"#include "MTStation.h"#include "AbsVelRecObjective.h"#include "Multi1DRecObjective.h"#include "SeismicDataComp.h"#include "CombinedRoughness.h"#include "SeismicModelDiff.h"#include "ResPkModel.h"#include "SurfaceWaveObjective.h"#include "SurfaceWaveData.h"#include "Adaptors.h"#include "MTFitSetup.h"#include "Util.h"#include "MTRoughness.h"#include "FileAbort.h"#include "SeisTools.h"#include "IsoJointConf.h"#include "MTInvConf.h"#include "SurfInvConf.h"#include "RecInvConf.h"#include "GAConf.h"Go to the source code of this file.
Typedefs | |
| typedef boost::shared_ptr < GeneralObjective > | pCGeneralObjective |
Enumerations | |
| enum | tgatype { pareto, anneal, pareto, anneal, pareto, anneal } |
Functions | |
| void | SetupAnnealingGA (boost::shared_ptr< GeneralGA > &GA, const double inittemp, const int annealinggeneration, const double coolingratio) |
| Set the annealing GA temperature parameter, we should only call this function if we are SURE that it is an AnnealinGA. More... | |
| void | SetupRecObjective (const string recinfofile, Multi1DRecObjective &RecObjective, const bool normrec, RecCalc::trfmethod rfmethod) |
| Use the information stored in recinforfile to setup the multi receiver function objective function. More... | |
| int | main (int argc, char *argv[]) |
Detailed Description
The program 1dinvga is the core joint inversion program for MT, receiver functions and surface wave data. It reads in its settings from the configuration file 1dinvga.conf , described below, and produces several log files and the final data, model(s) and plot file(s) for each inverted dataset.
General Configuration in 1dinvga.conf
The options in 1dinvga.conf are as follows:
General program options
| Name | Default | Description |
| verbose | false | Print detailed information about all models used in the invbersion, only useful for debugging |
| outputbase | The start of the name of all logfiles. This can also be specified on the command lin.e | |
| threads | 1 | If the program has been compiled with openmp support, this specifies the number of simultaneous threads. |
Genetic algorithm options
| Name | Default | Description |
| gatype | Can be either anneal or pareto. If set to anneal we use an annealing type GA that minimizes a weighted sum of the objective functions with the weights specified below. For this type of GA we also have to specify inittemp and coolingratio. If set to pareto, we use NSGA-II as described by Deb et al. (2002). In this case inittemp and coolingratio are ignored. | |
| inittemp | (For GA type anneal)The initial temperature in the cooling schedule. High values mean that there will be no big differences in probability between well fitting models and bad fitting models. Low values make those differences much more severe. If set to negative values the temperature will be determined from the average misfit values in the first iteration. | |
| annealinggeneration | (For GA type anneal) The generation number when the initial tempature is divided by coolingratio to enforce stronger discrimination between good and bad fitting models. | |
| coolingratio | (For GA type anneal) The value by which inittemp is divided at generation annealinggeneration Higher values mean stronger discrimination./TR> | |
| popsize | The population size for the genetic algorithm. Typical values 50...1000 | |
| generations | The number of generations for the GA to run.Typical values 50...500. /TR> | |
| mutationprob | The probability for a mutation within one member, typically 0.1...0.2. The probability for mutation of individual bits is calculated from this number depending on the length of the encoding specified below. | |
| crossoverprob | The probability for two model to perform crossover. Typically 0.4...0.7. | |
| elitist | Do we want to keep the best model(s) from the last generation (true), or fully depend on crossover and mutation with the risk of loosing good models (false). |
Data information
| Name | Default | Description |
| recinfofile | The file containing information about the receiver function data. This file is decribed below. | |
| mtinputdata | The name of the file containing the MT impedances | |
| dispdata | The name of the file containing the Rayleigh wave dispersion data |
Error floor information
| Name | Default | Description |
| tensorerror | 0.02 | Error floor for MT impedances, if selected below |
| reserror | 0.02 | Error floor for MT apparent resistivities, if selected below |
| phaseerror | 0.02 | Error floor for MT phase, if selected below |
| surferror | 0.02 | Error floor for Rayleigh wave dispersion data |
Global receiver function parameters
| Name | Default | Description |
| poisson | Poissons ratio, we calculate Vp from Vs using this ratio | |
| normrec | true | Are we inverting receiver functions that have been normalized to a maximum amplitude of 1 ? |
| starttime | Start of the fit window in seconds. Any data in the receiver function before this time is ignored. | |
| endtime | End of the fit window in seconds. Any data in the receiver function after this time is ignored. | |
| recmethod | Method used to calculate the measured data. Can be specdiv for spectral division, or iterdecon for iterative deconvolution. |
MT data parameters
| Name | Default | Description |
| mode | Which component of the tensor to fit. Can be xy or yx. | |
| mtfit | Which aspect of the MT data to fit. Possible values are shown in the table below. |
Possible mtfit values
| Name | Description |
| appres | Fit only apparent resistivity of selected mode |
| phase | Fit only phase of selected mode |
| resphase | Fit apparent resistivity and phase of selected mode |
| berd | Fit Berdichevskiy invariant. This is independent of selected mode |
| det | Fit determinant of impedance tensor. This is independent of selected mode. |
Regularization and fit information
| Name | Default | Description |
| usevrefmodel | false | Do we want to regularize the seismic parameters by the difference to a reference model (true/false) ? |
| vrefmodel | If true above, the name of the reference model file. | |
| mtfitexponent | 2 | The exponent for the difference between measured and calculated data. Higher exponents penalize any departure from the measured data. This ensures a tighter fit, but also increases the influence of noise. |
| recfitexponent | 2 | Same for receiver function data. |
| surffitexponent | 2 | Same for dispersion data. |
Parameter encoding
The number of layers is determined by the number of entries, it has to be equal for all encoding parameters.
| Name | Default | Description |
| thickbase | Minimum thickness for each layer in km | |
| thickstep | Stepsize for each layer in km | |
| thicksizes | Number of bits used to encode the thickness of each layer | |
| resbase | Base 10 logarithm of the minimum resistivity, i.e. 0 corresponds to 1 Ohm.m | |
| resstep | Stepsize for base 10 logarithm of resistivity | |
| ressizes | Number of bits used to encode resistivity | |
| svelbase | Minimum S-Wave velocity for each layer in km/s | |
| svelstep | Stepsize for S_Wave velocity for each layer | |
| svelsizes | Number of bits used to encode Vs |
Weighting of datasets
| Name | Default | Description |
| weights | These are the weights for the different objective functions. The current version needs exactly 5 weights. They correspond to MT misfit, Receiver function misfit, Surface Wave Misfit, Seismic Regularization and MT regularization, respectively. If a weight is set to zero the calculation of the misfit is completely disabled and the inversion thus runs faster. For pareto type inversion the value of weight is irrelevant as long as it is different from zero. For annealing it determines the influence of the corresponding dataset. |
Recinfo file
The recinfo file contains information about the different receiver functions we are trying to fit. Each line in the file corresponds to one receiver function and requires the following entries
| Ray parameter in s/km |  for gaussian filter | Waterlevel | Relative error | Shift in s | Type of data | Filename |
Type of data is a single character. If it is 'a' we are also fitting absolute velocity information as calculated by the program rfvel and described in Svenningsen, GJI 170, 2007. In this case the following additional fields are needed
| Weight for RF | Weight for absolute velocities | Filename for absolute velocities |
If data type is any other character we only fit receiver function information and the absolute velocity information cannot be there.
Other features
- You can gracefully stop the inversion any time by creating a file called "abort" in the working directory. The inversion will stop after the current iteration and will write out the usual information about best models etc. This is useful when the chosen number of generations was too high and the algorithm has already reached an acceptable solution. The file is automatically deleted.
Definition in file 1dinvga.cpp.
Typedef Documentation
| typedef boost::shared_ptr<GeneralObjective> pCGeneralObjective |
Definition at line 187 of file 1dinvga.cpp.
Enumeration Type Documentation
| enum tgatype |
| Enumerator | |
|---|---|
| pareto | |
| anneal | |
| pareto | |
| anneal | |
| pareto | |
| anneal | |
Definition at line 182 of file 1dinvga.cpp.
Function Documentation
| int main | ( | int | argc, |
| char * | argv[] | ||
| ) |
Definition at line 267 of file 1dinvga.cpp.
References anneal, gplib::GAConf::annealinggeneration, gplib::GAConf::coolingratio, gplib::GAConf::crossoverprob, gplib::SurfInvConf::dispdata, gplib::GAConf::elitist, gplib::RecInvConf::endtime, f(), gplib::UniquePop::Find(), gplib::GAConf::gatype, gplib::GAConf::generations, gplib::SurfInvConf::GetData(), gplib::RecInvConf::GetData(), gplib::MTInvConf::GetData(), gplib::GAConf::GetData(), gplib::IsoJointConf::GetData(), gplib::MTStation::GetData(), gplib::GrayTranscribe::GetValues(), h(), gplib::GAConf::inittemp, gplib::UniquePop::Insert(), MTData, gplib::MTInvConf::mtfitexponent, gplib::MTInvConf::mtinputdata, MTObjective, gplib::GAConf::mutationprob, gplib::RecInvConf::normrec, gplib::IsoJointConf::outputbase, pareto, gplib::IsoJointConf::poisson, gplib::GAConf::popsize, gplib::UniquePop::PrintAll(), gplib::SurfaceWaveData::ReadFile(), gplib::ResPkModel::ReadModel(), gplib::RecInvConf::recinfofile, gplib::RecInvConf::recmethod, RecObjective(), gplib::IsoJointConf::resbase, gplib::IsoJointConf::ressizes, gplib::IsoJointConf::resstep, RFAbsVel, gplib::GeneralPropagation::SetParams(), SetupAnnealingGA(), gplib::SetupMTFitParameters(), SetupRecObjective(), gplib::RecInvConf::starttime, gplib::SurfInvConf::surferror, gplib::SurfInvConf::surffitexponent, gplib::IsoJointConf::svelbase, gplib::IsoJointConf::svelsizes, gplib::IsoJointConf::svelstep, gplib::IsoJointConf::thickbase, gplib::IsoJointConf::thicksizes, gplib::IsoJointConf::thickstep, gplib::IsoJointConf::usevrefmodel, gplib::IsoJointConf::verbose, version, gplib::IsoJointConf::vrefmodel, and gplib::IsoJointConf::weights.
| void SetupAnnealingGA | ( | boost::shared_ptr< GeneralGA > & | GA, |
| const double | inittemp, | ||
| const int | annealinggeneration, | ||
| const double | coolingratio | ||
| ) |
Set the annealing GA temperature parameter, we should only call this function if we are SURE that it is an AnnealinGA.
Definition at line 190 of file 1dinvga.cpp.
References gplib::AnnealingGA::SetParams().
Referenced by main().
| void SetupRecObjective | ( | const string | recinfofile, |
| Multi1DRecObjective & | RecObjective, | ||
| const bool | normrec, | ||
| RecCalc::trfmethod | rfmethod | ||
| ) |
Use the information stored in recinforfile to setup the multi receiver function objective function.
Definition at line 212 of file 1dinvga.cpp.
References gplib::Multi1DRecObjective::AddAbsVelFunction(), gplib::Multi1DRecObjective::AddRecFunction(), gplib::SurfaceWaveData::ReadAscii(), and RecData.
Referenced by main().
