00001 #include <iostream>
00002 #include "Multi1DRecObjective.h"
00003 #include "SurfaceWaveData.h"
00004 #include "ResPkModel.h"
00005 #include "SeismicDataComp.h"
00006 #include "CCalcRecConf.h"
00007 #include <string>
00008 #include <boost/bind.hpp>
00009 #include "SeisTools.h"
00010 using namespace std;
00011 int main()
00012 {
00013 string datafilename, modelfilename;
00014 cout << "Base for Datafiles: ";
00015 cin >> datafilename;
00016 cout << "Name of Modelfile: ";
00017 cin >> modelfilename;
00018 double slowness;
00019 cout << "Slowness: ";
00020 cin >> slowness;
00021 CCalcRecConf Config;
00022 Config.GetData("calcrec.conf");
00023
00024 SurfaceWaveData SurfData;
00025 boost::shared_ptr<const SeismicDataComp> RecData(
00026 new SeismicDataComp(datafilename + ".rec", SeismicDataComp::sac));
00027 Normalize(const_cast<SeismicDataComp*> (RecData.get())->GetData());
00028 SurfData.ReadAscii(datafilename + ".vel");
00029 ResPkModel Model;
00030 Model.GetData(modelfilename);
00031
00032 Multi1DRecObjective SeisObjective;
00033
00034 SeisObjective.AddRecFunction(RecData, Config.shift, Config.sigma,
00035 Config.cc, slowness, RecCalc::iterdecon, 0.02, true);
00036 SeisObjective.AddAbsVelFunction(RecData, SurfData, Config.shift,
00037 Config.sigma, Config.cc, slowness, RecCalc::iterdecon, 0.02, true, 0.0, 1.0);
00038
00039 SeisObjective.SetPoisson(1.7320508);
00040 SeisObjective.SetTimeWindow(0, 40);
00041 const unsigned int nlayers = Model.GetThickness().size();
00042 const unsigned int nparams = nlayers * 2;
00043 ttranscribed member(nparams);
00044 for (unsigned i = 0; i < nlayers; ++i)
00045 {
00046 member(i) = Model.SetThickness().at(i);
00047 member(i + nlayers) = Model.SetSVelocity().at(i);
00048 }
00049 cout << "Misfit: " << SeisObjective.CalcPerformance(member) << endl;
00050 SeisObjective.WriteData("obj");
00051 SeisObjective.WriteModel("obj.mod");
00052 }
