Skip to content
Snippets Groups Projects
Commit 1406914e authored by Bertram Kopf's avatar Bertram Kopf
Browse files

introduced ratio for parameters to be modified in evo minimizer

parent 0abe589d
No related branches found
No related tags found
No related merge requests found
......@@ -52,6 +52,7 @@ ParserBase::ParserBase(int argc,char **argv)
, _ratioMcToData(100000)
, _evoIterations(100)
, _evoPopulation(20)
, _evoRatioOfModParams(1.)
, _cacheAmps(true)
, _calcContributionError(false)
, _saveContributionHistos(false)
......@@ -108,6 +109,7 @@ ParserBase::ParserBase(int argc,char **argv)
("ratioMcToData",po::value<int>(&_ratioMcToData), "number of MC events defined by ratio #MCs/#Data")
("evoPopulation",po::value<int>(&_evoPopulation), "iteration population for evo minimizer")
("evoIterations",po::value<int>(&_evoIterations), "number of iterations for evo minimizer")
("evoRatioOfModParams",po::value<double>(&_evoRatioOfModParams), "chosen (avereaged) ratio of fit parameters to be changed for each population (value between 0. and 1.")
("cacheAmps",po::value<bool>(&_cacheAmps), "cache amplitudes")
("contributionError",po::value<bool>(&_calcContributionError), "calculate the wave contribution error")
("saveContributionHistos",po::value<bool>(&_saveContributionHistos), "creates a histogram root-file for each contribution of Option: calcContribution")
......@@ -215,6 +217,11 @@ bool ParserBase::parseCommandLine(int argc, char **argv)
Warning << "ErrorLogger not (properly) set -> Use mode 'DEBUG' " ; // << endmsg;
}
if (_evoRatioOfModParams<=0. || _evoRatioOfModParams>1.){
Alert << "_evoRatioOfModParams = " << _evoRatioOfModParams << " not possible\n"
<< "value must be between 0. and 1. !!!!" << endmsg;
exit(1);
}
if(_verbose){
std::cout << "\nRunning with the following options using " << _configFile << ":\n\n"
......
......@@ -73,6 +73,7 @@ public:
const int ratioMcToData() const {return _ratioMcToData;}
const int evoPopulation() const {return _evoPopulation;}
const int evoIterations() const {return _evoIterations;}
const double evoRatioOfModParams() const {return _evoRatioOfModParams;}
const bool cacheAmps() const {return _cacheAmps;}
const bool calcContributionError() const {return _calcContributionError;}
const bool saveContributionHistos() const {return _saveContributionHistos;}
......@@ -130,6 +131,7 @@ protected:
int _ratioMcToData;
int _evoIterations;
int _evoPopulation;
double _evoRatioOfModParams;
bool _cacheAmps;
bool _calcContributionError;
bool _saveContributionHistos;
......
......@@ -35,6 +35,7 @@ noOfClients = 2
noOfThreads = 2
evoIterations = 50
evoPopulation = 160
evoRatioOfModParams = 0.05
# Productions
......
......@@ -32,6 +32,7 @@
#include "FitParams/FitParColBase.hh"
#include "ErrLogger/ErrLogger.hh"
#include "PwaUtils/GlobalEnv.hh"
#include "ConfigParser/ParserBase.hh"
const double EvoMinimizer::DECREASESIGMAFACTOR = 0.9;
const double EvoMinimizer::INCREASESIGMAFACTOR = 1.1;
......@@ -46,10 +47,19 @@ EvoMinimizer::EvoMinimizer(std::shared_ptr<AbsFcn> theAbsFcnPtr, std::shared_ptr
AbsPawianMinimizer(theAbsFcnPtr, upar)
,_population(population)
, _iterations(iterations)
, _evoRatioOfModParams(GlobalEnv::instance()->parser()->evoRatioOfModParams())
, _currentBestParams(theAbsFcnPtr->defaultFitValParms())
, _defaultFitErrParms(theAbsFcnPtr->defaultFitErrParms())
, _currentResultFileName("currentEvoResult"+GlobalEnv::instance()->outputFileNameSuffix()+".dat")
{
if (_evoRatioOfModParams <= 0. || _evoRatioOfModParams>1.){
Alert << "_evoRatioOfModParams = " << _evoRatioOfModParams << " not possible\n"
<< "value must be set between 0. and 1. !!!!" << endmsg;
exit(1);
}
Info << "population: " << _population
<<"\niterations: " << _iterations
<< "\nratio of parameters to be modified: " << _evoRatioOfModParams << endmsg;
}
......@@ -156,46 +166,53 @@ void EvoMinimizer::ShuffleParams(){
typedef boost::variate_generator<RandomGenerator&, NormalDistribution> GaussianGenerator;
static RandomGenerator rng(static_cast<unsigned> (time(0)));
for(unsigned int i=0; i<_tmpParams->Params().size(); i++){
// Don't touch fixed parameters
if(_tmpParams->IsFixed(i)){
bool acceptNewParams=false;
while(!acceptNewParams){
for(unsigned int i=0; i<_tmpParams->Params().size(); i++){
// Don't touch fixed parameters
if(_tmpParams->IsFixed(i)){
continue;
}
// Decide whether parameter is increased or decreased
boost::random::uniform_int_distribution<> coin(0,1);
bool c = coin(rng);
// Initialize gaussian width as parameter error
double sigma = _tmpParams->Error(i);
// If gaussian collides with parameter limit, reduce width
if(_tmpParams->HasLimits(i)){
}
boost::random::uniform_real_distribution<> coinReal(0.,1.);
double c01 = coinReal(rng);
if (c01 > _evoRatioOfModParams) continue; //only 10% of the parameters changed
acceptNewParams=true;
// Decide whether parameter is increased or decreased
boost::random::uniform_int_distribution<> coin(0,1);
bool c = coin(rng);
// Initialize gaussian width as parameter error
double sigma = _tmpParams->Error(i);
// If gaussian collides with parameter limit, reduce width
if(_tmpParams->HasLimits(i)){
if(c && (2*sigma > (_tmpParams->UpperLimit(i) - _tmpParams->Value(i)))){
sigma = (_tmpParams->UpperLimit(i) - _tmpParams->Value(i)) / 2.0;
sigma = (_tmpParams->UpperLimit(i) - _tmpParams->Value(i)) / 2.0;
}
else if(!c && (2*sigma > (_tmpParams->Value(i) - _tmpParams->LowerLimit(i)))){
sigma = (_tmpParams->Value(i) - _tmpParams->LowerLimit(i)) / 2.0;
sigma = (_tmpParams->Value(i) - _tmpParams->LowerLimit(i)) / 2.0;
}
}
// Get random number and set new parameter
NormalDistribution gaussian_dist(0, sigma);
GaussianGenerator generator(rng, gaussian_dist);
double val = fabs(generator());
if(c) _tmpParams->SetValue(i, _tmpParams->Value(i) + val);
else _tmpParams->SetValue(i, _tmpParams->Value(i) - val);
// Check for limits
if(_tmpParams->HasLimits(i)){
if(_tmpParams->Value(i) < _tmpParams->LowerLimit(i))
_tmpParams->SetValue(i, _tmpParams->LowerLimit(i));
if(_tmpParams->Value(i) > _tmpParams->UpperLimit(i))
_tmpParams->SetValue(i, _tmpParams->UpperLimit(i));
}
}
// Get random number and set new parameter
NormalDistribution gaussian_dist(0, sigma);
GaussianGenerator generator(rng, gaussian_dist);
double val = fabs(generator());
if(c) _tmpParams->SetValue(i, _tmpParams->Value(i) + val);
else _tmpParams->SetValue(i, _tmpParams->Value(i) - val);
// Check for limits
if(_tmpParams->HasLimits(i)){
if(_tmpParams->Value(i) < _tmpParams->LowerLimit(i))
_tmpParams->SetValue(i, _tmpParams->LowerLimit(i));
if(_tmpParams->Value(i) > _tmpParams->UpperLimit(i))
_tmpParams->SetValue(i, _tmpParams->UpperLimit(i));
}
}
}
}
......
......@@ -47,6 +47,8 @@ public:
private:
int _population;
int _iterations;
double _evoRatioOfModParams;
fitParCol _currentBestParams;
fitParCol _defaultFitErrParms;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment