portelli/LatAnalyze
1
0
Fork 0
mirror of https://github.com/aportelli/LatAnalyze.git synced 2024-11-10 00:45:36 +00:00
Code Releases Activity

Added cotanh fit function

Browse Source
This commit is contained in:
Andrew Zhen Ning Yong 2019-05-28 07:38:43 +01:00
parent c49da2ec2d
commit 4cab49cd92
1 changed files with 40 additions and 4 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

44
physics/2pt-fit.cpp
View File

@ -16,7 +16,7 @@ int main(int argc, char *argv[])
// parse arguments /////////////////////////////////////////////////////////
OptParser opt;
bool parsed, doPlot, doHeatmap, doCorr, fold;
string corrFileName, model, outFileName, outFmt;
string corrFileName, model, outFileName, outFmt, save;
Index ti, tf, shift, nPar, thinning;
double svdTol;
Minimizer::Verbosity verbosity;
@ -30,7 +30,7 @@ int main(int argc, char *argv[])
opt.addOption("s", "shift" , OptParser::OptType::value , true,
"time variable shift", "0");
opt.addOption("m", "model" , OptParser::OptType::value , true,
"fit model (exp|exp2|exp3|cosh|cosh2|cosh3|explin|const|<interpreter code>)", "cosh");
"fit model (exp|exp2|exp3|sinh|cosh|cosh2|cosh3|cotanh|explin|const|<interpreter code>)", "cosh");
opt.addOption("" , "nPar" , OptParser::OptType::value , true,
"number of model parameters for custom models "
"(-1 if irrelevant)", "-1");
@ -48,6 +48,8 @@ int main(int argc, char *argv[])
"show the fit plot");
opt.addOption("h", "heatmap" , OptParser::OptType::trigger, true,
"show the fit correlation heatmap");
opt.addOption("s", "save", OptParser::OptType::value, true,
"saves the source and .pdf", "");
opt.addOption("", "help" , OptParser::OptType::trigger, true,
"show this help message and exit");
parsed = opt.parse(argc, argv);
@ -71,6 +73,7 @@ int main(int argc, char *argv[])
fold = opt.gotOption("fold");
doPlot = opt.gotOption("p");
doHeatmap = opt.gotOption("h");
save = opt.optionValue("s");
switch (opt.optionValue<unsigned int>("v"))
{
case 0:
@ -117,7 +120,7 @@ int main(int argc, char *argv[])
// make model //////////////////////////////////////////////////////////////
DoubleModel mod;
bool coshModel = false, linearModel = false, constModel = false;
bool sinhModel = false, coshModel = false, cotanhModel = false, linearModel = false, constModel = false;
if ((model == "exp") or (model == "exp1"))
{
@ -144,6 +147,15 @@ int main(int argc, char *argv[])
+ p[5]*exp(-p[4]*x[0]);
}, 1, nPar);
}
else if (model == "sinh")
{
sinhModel = true;
nPar = 2;
mod.setFunction([nt](const double *x, const double *p)
{
return p[1]*(exp(-p[0]*x[0])-exp(-p[0]*(nt-x[0])));
}, 1, nPar);
}
else if ((model == "cosh") or (model == "cosh1"))
{
coshModel = true;
@ -174,6 +186,15 @@ int main(int argc, char *argv[])
+ p[5]*(exp(-p[2]*x[0])+exp(-p[4]*(nt-x[0])));
}, 1, nPar);
}
else if (model == "cotanh")
{
cotanhModel = true;
nPar = 2;
mod.setFunction([nt](const double *x, const double *p)
{
return p[0] + p[1]/tanh(0.1892*(nt/2-x[0])); // 0.1892 is the meson eff mass in lattice unit from 2pt fit
}, 1, nPar);
}
else if (model == "explin")
{
linearModel = true;
@ -247,6 +268,11 @@ int main(int argc, char *argv[])
init(0) = data.y(nt/4, 0)[central];
}
else if(cotanhModel)
{
init(0) = 0.5*(data.y(0, 0)[central]+ data.y(nt-1, 0)[central]);
init(1) = 0.5*(data.y(0, 0)[central]- data.y(nt-1, 0)[central]);
}
else
{
init(0) = log(data.y(nt/4, 0)[central]/data.y(nt/4 + 1, 0)[central]);
@ -269,9 +295,14 @@ int main(int argc, char *argv[])
{
globMin.setLowLimit(p, -10*fabs(init(0)));
locMin.setLowLimit(p, -10*fabs(init(0)));
// cout << "Suppressing low limits" << endl;
globMin.setHighLimit(p, 10*fabs(init(0)));
}
else if(cotanhModel)
{
globMin.setLowLimit(p,-0.5);
locMin.setLowLimit(p, -0.3);
globMin.setHighLimit(p,1);
}
else
{
globMin.setLowLimit(p, 0.);
@ -328,6 +359,11 @@ int main(int argc, char *argv[])
p << Color("rgb 'blue'") << PlotFunction(fit.getModel(), 0, nt - 1);
p << Color("rgb 'red'") << PlotData(data.getData());
p.display();
if(save != "")
{
cout << "Saving plot and source code to " << save << endl;
p.save(save);
}
// effective mass plot //////////////////////////////////////////////////////
if (!constModel)
{