Skip to content
Snippets Groups Projects
Commit e445b1dd authored by Carl Philipp Klemm's avatar Carl Philipp Klemm
Browse files

add the ability to add random noise to output spectra

parent b6e7112b
Branches
No related tags found
No related merge requests found
#include "basicmath.h" #include "basicmath.h"
#include <algorithm> #include <algorithm>
#include <random>
#include "eistype.h" #include "eistype.h"
...@@ -115,3 +116,18 @@ std::vector<eis::DataPoint> eis::rescale(const std::vector<eis::DataPoint>& data ...@@ -115,3 +116,18 @@ std::vector<eis::DataPoint> eis::rescale(const std::vector<eis::DataPoint>& data
} }
return output; return output;
} }
void eis::noise(std::vector<eis::DataPoint>& data, double amplitude, bool relative)
{
std::random_device rd;
std::default_random_engine rande(rd());
for(eis::DataPoint& point : data)
{
double realNoise = (static_cast<double>(rande() >> 1) / (rande.max() >> 1))*amplitude;
point.im.real(relative ? point.im.real()+realNoise/point.im.real() : point.im.real()+realNoise);
double imgNoise = (static_cast<double>(rande() >> 1) / (rande.max() >> 1))*amplitude;
point.im.imag(relative ? point.im.imag()+imgNoise/point.im.imag() : point.im.imag()+imgNoise);
}
}
...@@ -12,4 +12,5 @@ namespace eis ...@@ -12,4 +12,5 @@ namespace eis
fvalue median(std::vector<fvalue> data); fvalue median(std::vector<fvalue> data);
std::complex<fvalue> median(const std::vector<eis::DataPoint>& data); std::complex<fvalue> median(const std::vector<eis::DataPoint>& data);
std::vector<eis::DataPoint> rescale(const std::vector<eis::DataPoint>& data, size_t outputSize); std::vector<eis::DataPoint> rescale(const std::vector<eis::DataPoint>& data, size_t outputSize);
void noise(std::vector<eis::DataPoint>& data, double amplitude, bool relative);
} }
...@@ -3,6 +3,7 @@ ...@@ -3,6 +3,7 @@
#include <chrono> #include <chrono>
#include <cmath> #include <cmath>
#include "basicmath.h"
#include "model.h" #include "model.h"
#include "log.h" #include "log.h"
#include "options.h" #include "options.h"
...@@ -36,7 +37,8 @@ static void paramSweepCb(std::vector<eis::DataPoint>& data, const std::vector<fv ...@@ -36,7 +37,8 @@ static void paramSweepCb(std::vector<eis::DataPoint>& data, const std::vector<fv
std::cout<<'.'<<std::flush; std::cout<<'.'<<std::flush;
} }
static void runSweep(const std::string& modelString, eis::Range omega, bool normalize = false, bool reduce = false, bool hertz = false, bool invert = false) static void runSweep(const std::string& modelString, eis::Range omega, bool normalize = false,
bool reduce = false, bool hertz = false, bool invert = false, double noise = 0)
{ {
std::vector<eis::DataPoint> results; std::vector<eis::DataPoint> results;
...@@ -64,6 +66,9 @@ static void runSweep(const std::string& modelString, eis::Range omega, bool norm ...@@ -64,6 +66,9 @@ static void runSweep(const std::string& modelString, eis::Range omega, bool norm
eis::Log(eis::Log::INFO)<<"results:"; eis::Log(eis::Log::INFO)<<"results:";
} }
if(noise > 0)
eis::noise(results, noise, false);
std::cout<<(hertz ? "freqency" : "omega")<<",real,im\n"; std::cout<<(hertz ? "freqency" : "omega")<<",real,im\n";
for(const eis::DataPoint& res : results) for(const eis::DataPoint& res : results)
...@@ -107,7 +112,7 @@ int main(int argc, char** argv) ...@@ -107,7 +112,7 @@ int main(int argc, char** argv)
switch(config.mode) switch(config.mode)
{ {
case MODE_SWEEP: case MODE_SWEEP:
runSweep(config.modelStr, config.omegaRange, config.normalize, config.reduce, config.hertz, config.invert); runSweep(config.modelStr, config.omegaRange, config.normalize, config.reduce, config.hertz, config.invert, config.noise);
break; break;
case MODE_PARAM_SWEEP: case MODE_PARAM_SWEEP:
runParamSweep(); runParamSweep();
......
...@@ -24,6 +24,7 @@ static struct argp_option options[] = ...@@ -24,6 +24,7 @@ static struct argp_option options[] =
{"reduce", 'r', 0, 0, "reduce values to \"interesting\" range" }, {"reduce", 'r', 0, 0, "reduce values to \"interesting\" range" },
{"hz", 'h', 0, 0, "freqency values as temporal frequency instead of angular frequency"}, {"hz", 'h', 0, 0, "freqency values as temporal frequency instead of angular frequency"},
{"invert", 'i', 0, 0, "inverts the imaginary axis"}, {"invert", 'i', 0, 0, "inverts the imaginary axis"},
{"noise", 'x', "[AMPLITUDE]", 0, "add noise to output"},
{ 0 } { 0 }
}; };
...@@ -42,6 +43,7 @@ struct Config ...@@ -42,6 +43,7 @@ struct Config
bool reduce = false; bool reduce = false;
bool hertz = false; bool hertz = false;
bool invert = false; bool invert = false;
double noise = 0;
Config(): omegaRange(1, 1001, 1, true) Config(): omegaRange(1, 1001, 1, true)
{} {}
...@@ -123,6 +125,9 @@ parse_opt (int key, char *arg, struct argp_state *state) ...@@ -123,6 +125,9 @@ parse_opt (int key, char *arg, struct argp_state *state)
} }
break; break;
} }
case 'x':
config->noise = std::stod(std::string(arg));
break;
default: default:
return ARGP_ERR_UNKNOWN; return ARGP_ERR_UNKNOWN;
} }
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment