addnoise.cpp

Go to the documentation of this file.

/*! \file 
 * This program adds gaussian random noise to the impedance estimates taken
 * from infilename and writes a new ".mtt" file with the noisy data and updated error
 * information. The output file has the same name as the input file with ".ns" appended.
 * The noise level is the relative noise for real and imaginary parts of the impedance at
 * each frequency.
 */

#include <iostream>
#include <complex>
#include <ctime>
#include "MTStation.h"
#include "Util.h"
#include <boost/random/lagged_fibonacci.hpp>
#include <boost/random/normal_distribution.hpp>
#include <boost/random/variate_generator.hpp>

using namespace std;


string version = "$Id: addnoise.cpp 1706 2008-06-06 07:14:04Z mmoorkamp $";

void AddNoise(std::complex<double> &impelement, double &noiseest,
                const double noiselevel, boost::lagged_fibonacci607 &generator) {
        complex<double> oldimp(impelement);
        impelement = boost::variate_generator<boost::lagged_fibonacci607&,
        boost::normal_distribution<> >(generator,
                        boost::normal_distribution<>(impelement.real(),
                                        fabs(impelement.real()* noiselevel)))();
        impelement += complex<double>(0.0,1.0)*boost::variate_generator<boost::lagged_fibonacci607&,
        boost::normal_distribution<> >(generator,
                        boost::normal_distribution<>(oldimp.imag(),
                                        fabs(oldimp.imag()* noiselevel)))();
        noiseest = noiselevel * abs(oldimp);
}

int main(int argc, char *argv[]) {
        string infilename, outfilename;
        double noiselevel;
        MTStation Data;
        boost::lagged_fibonacci607
                        generator(static_cast<unsigned int>(std::time(0)));
        complex<double> I(0.0, 1.0);
        cout << "This is addnoise: Add random noise to MT impedance estimates."
                        << endl<< endl;
        cout << " Usage: addnoise  infilename  noiselevel" << endl;
        cout
                        << " The output files will have the same name as the input files + .ns "
                        << endl << endl;
        cout << " This is Version: " << version << endl << endl;
        if (argc > 3) {
                infilename = argv[1];
                outfilename = argv[2];
                noiselevel = atof(argv[3]);
        }

        if (argc == 3) {
                infilename = argv[1];

                noiselevel = atof(argv[2]);
        } else {

                infilename = AskFilename("Infilename: ");
                cout << "Noiselevel:";
                cin >> noiselevel;
                outfilename = infilename + ".ns";
        }

        Data.GetData(infilename);
        double zxxerr, zxyerr, zyxerr, zyyerr;
        for (unsigned int i = 0; i < Data.GetMTData().size(); ++i) {
                AddNoise(Data.SetMTData().at(i).SetZxx(), zxxerr, noiselevel, generator);
                AddNoise(Data.SetMTData().at(i).SetZxy(), zxyerr, noiselevel, generator);
                AddNoise(Data.SetMTData().at(i).SetZyx(), zyxerr, noiselevel, generator);
                AddNoise(Data.SetMTData().at(i).SetZyy(), zyyerr, noiselevel, generator);
                Data.SetMTData().at(i).SetErrors(zxxerr, zxyerr, zyxerr, zyyerr);
        }
        Data.WriteAsMtt(outfilename);
}

---