00001 #ifndef C1DRECOBJECTIVE_H 00002 #define C1DRECOBJECTIVE_H 00003 #include "PlottableObjective.h" 00004 #include "RecCalc.h" 00005 #include "SeismicDataComp.h" 00006 #include <boost/shared_ptr.hpp> 00007 00008 //! Calculate the misfit between observed receiver function for a given 1D model by calculating a synthetic receiver function from that model 00009 /*! The constructor takes a few essential parameters that are expected not to change during the 00010 * program (in most cases inversion), e.g. the data to fit. See also #C1DRecObjective for the necessary 00011 * constructor parameters. 00012 * 00013 * There are two ways to calculate the Misfit for a given model, either calling GeneralObjective#CalcPerformance, 00014 * or, for parallel programs, calling the three function PreParallel, SafeParallel and PostParallel in that order. 00015 * In both cases you have to give the parameter #member, which is a ::ttranscribed of size 2n, where 00016 * n is the number of layers of the model. The first n entries are the layer thicknesses in km and the last n entries 00017 * are the S-wave velocities in km/s. Densities and P-Wave velocities are calculated from S-Velocities by hard-coded relationships. 00018 */ 00019 class C1DRecObjective : public PlottableObjective { 00020 private: 00021 //! We need to write some temporary files, this is the name 00022 std::string tempname; 00023 //! The object holding the calculated synthetic data 00024 SeismicDataComp RecSynthData; 00025 //! The slowness used in the calculation 00026 double slowness; 00027 //! The object holding the 1D model for forward calculation 00028 ResPkModel Model; 00029 //! The object used to calculate the synthetics 00030 RecCalc RecCalculator; 00031 //! Errorfloor 00032 double errorlevel; 00033 //! Poisson's ratio vp/vs, this might become part of the model 00034 double poisson; 00035 //! We might not want to use the whole receiver function for misfit calculation, so we can give a starttime 00036 double starttime; 00037 //! and an endtime 00038 double endtime; 00039 //! The measured data is stored in a shared pointer. 00040 /*! For parallel implementations we often work with several instances of this 00041 * class that all work with the same ObservedData. To save memory we use a shared_pointer 00042 * and to avoid collisions SeismicDataComp is const */ 00043 boost::shared_ptr<const SeismicDataComp> ObservedData; 00044 protected: 00045 const SeismicDataComp &GetObservedData(){return *ObservedData;} 00046 public: 00047 //! return a pointer to a copy of the current object 00048 virtual C1DRecObjective *clone() const { 00049 return new C1DRecObjective(*this); 00050 } 00051 //! Set the time window used for misfit calculations, start and end are seconds from start of receiver function 00052 void SetTimeWindow(const double start, const double end) { 00053 starttime = start; 00054 endtime=end; 00055 } 00056 //! Set the errorlevel for fit, this is relative to the maximum amplitude, not for each individual data point 00057 void SetErrorLevel(const double level) { 00058 errorlevel = level; 00059 } 00060 //! Set poisson's ratio, at the moment the same for all layers, used for calculating P-velocity 00061 void SetPoisson(const double ratio) { 00062 poisson = ratio; 00063 } 00064 //! Write the synthetic data to a sac file with name filename, makes only sense after calculating the misfit 00065 virtual void WriteData(const std::string &filename) { 00066 RecSynthData.WriteAsSac(filename); 00067 } 00068 //! Write the current model to ascii file for calculations 00069 void WriteModel(const std::string &filename) { 00070 Model.WriteData(filename); 00071 } 00072 //! Write the current model to ascii file for plotting 00073 void WritePlot(const std::string &filename) { 00074 Model.WritePlot(filename); 00075 } 00076 //! We have to write runfiles before parallel execution 00077 virtual void PreParallel(const ttranscribed &member); 00078 //! We also clean up files serially 00079 virtual double PostParallel(const ttranscribed &member); 00080 //! Calculate the misfit between the data calculated from model vector member and measured data given in the constructor. 00081 virtual void SafeParallel(const ttranscribed &member); 00082 C1DRecObjective(const C1DRecObjective &Old); 00083 C1DRecObjective& operator=(const C1DRecObjective& source); 00084 //! The constructor needs a few essential parameters 00085 /*! @param TheRecData Shared pointer to the measured receiver function 00086 * @param myshift the time shift used for calculating the measured receiver function 00087 * @param myomega omega used for calculating the measured receiver function 00088 * @param mysigma sigma used for calculating the measured receiver function 00089 * @param myc water level used for calculating the measured receiver function 00090 * @param myslowness slowness used for calculating the measured receiver function 00091 * @param method The method used to calculate the observed data. Can be specdiv or iterdecon 00092 */ 00093 C1DRecObjective(boost::shared_ptr<const SeismicDataComp> TheRecData, const int myshift, 00094 const double myomega, const double mysigma, const double myc, 00095 const double myslowness, 00096 const RecCalc::trfmethod method = RecCalc::specdiv); 00097 00098 virtual ~C1DRecObjective(); 00099 friend class AbsVelRecObjective; 00100 }; 00101 00102 #endif // C1DRECOBJECTIVE_H
1.5.5