All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Classes | Public Types | Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
ana::Fitter Class Reference

Perform MINUIT fits in one or two dimensions. More...

#include <Fit.h>

Inheritance diagram for ana::Fitter:

Classes

struct  SeedPt
 

Public Types

enum  Verbosity { kQuiet, kVerbose }
 
enum  Precision {
  kFast = 0, kNormal = 1, kCareful = 2, kGradDesc = 3,
  kAlgoMask = 3, kIncludeSimplex = 4, kIncludeHesse = 8, kIncludeMinos = 16
}
 

Public Member Functions

void SetPrecision (Precision prec)
 
 Fitter (const IExperiment *expt, std::vector< const IFitVar * > vars, std::vector< const ISyst * > systs={}, Precision prec=kNormal)
 
double Fit (osc::IOscCalcAdjustable *seed, SystShifts &bestSysts=junkShifts, const std::map< const IFitVar *, std::vector< double >> &seedPts={}, const std::vector< SystShifts > &systSeedPts={}, Verbosity verb=kVerbose) const
 
double Fit (osc::IOscCalcAdjustable *seed, SystShifts &systSeed, Verbosity verb) const
 Variant with no seedPts. More...
 
double Fit (osc::IOscCalcAdjustable *seed, Verbosity verb) const
 Variant with no seedPts and no systematics result returned. More...
 
double Fit (SystShifts &systSeed, Verbosity verb=kVerbose) const
 Variant with no oscillations - useful for ND fits. More...
 
TMatrixDSym * GetCovariance ()
 Return the fit covariance. More...
 
int GetCovarianceStatus ()
 covariance matrix status (0 = not valid, 1 approximate, 2, full but made pos def, 3 accurate and not pos def) More...
 
std::vector< std::string > GetParamNames ()
 Return the fit names. More...
 
std::vector< double > GetPreFitValues ()
 Return the prefit values. More...
 
std::vector< double > GetPreFitErrors ()
 Return the prefit errors. More...
 
std::vector< double > GetPostFitValues ()
 Return the postfit values. More...
 
std::vector< double > GetPostFitErrors ()
 Return the postfit errors. More...
 
std::vector< std::pair< double,
double > > 
GetMinosErrors ()
 Return the minos errors. More...
 
int GetNFCN ()
 Return number of function calls. More...
 
double GetEDM ()
 Return edm form the fit. More...
 
bool GetIsValid ()
 Say whether the fit was good. More...
 
SystShifts GetSystShifts () const
 
virtual double DoEval (const double *pars) const override
 Evaluate the log-likelihood, as required by MINUT interface. More...
 
virtual unsigned int NDim () const override
 
FitterClone () const override
 

Protected Member Functions

std::vector< SeedPtExpandSeeds (const std::map< const IFitVar *, std::vector< double >> &seedPts, std::vector< SystShifts > systSeedPts) const
 
std::unique_ptr
< ROOT::Math::Minimizer > 
FitHelperSeeded (osc::IOscCalcAdjustable *seed, SystShifts &systSeed, Verbosity verb) const
 Helper for FitHelper. More...
 
double FitHelper (osc::IOscCalcAdjustable *seed, SystShifts &bestSysts, const std::map< const IFitVar *, std::vector< double >> &seedPts, std::vector< SystShifts > systSeedPts, Verbosity verb) const
 Helper for Fit. More...
 
void DecodePars (const double *pars) const
 Updates mutable fCalc and fShifts. More...
 

Protected Attributes

const IExperimentfExpt
 
std::vector< const IFitVar * > fVars
 
std::vector< const ISyst * > fSysts
 
Precision fPrec = kNormal
 
osc::IOscCalcAdjustablefCalc
 
SystShifts fShifts
 
int fNEval = 0
 
int fNEvalGrad = 0
 
int fNEvalFiniteDiff = 0
 
double fEdm = -1
 
bool fIsValid = false
 
TMatrixDSym * fCovar
 
bool fCovarStatus
 
std::vector< std::string > fParamNames
 
std::vector< double > fPreFitValues
 
std::vector< double > fPreFitErrors
 
std::vector< double > fPostFitValues
 
std::vector< double > fPostFitErrors
 
std::vector< std::pair< double,
double > > 
fMinosErrors
 
std::vector< std::pair< double,
double > > 
fTempMinosErrors
 

Detailed Description

Perform MINUIT fits in one or two dimensions.

Definition at line 37 of file Fit.h.

Member Enumeration Documentation

Enumerator
kFast 
kNormal 
kCareful 
kGradDesc 
kAlgoMask 
kIncludeSimplex 
kIncludeHesse 
kIncludeMinos 

Definition at line 42 of file Fit.h.

42  {
43  // You must select one of these. The first three codes match the settings
44  // used by migrad. The fourth is a custom minimizer.
45  kFast = 0,
46  kNormal = 1,
47  kCareful = 2,
48  kGradDesc = 3,
49  // Allow bitmask operations to extract these first four options
50  kAlgoMask = 3,
51  // You may optionally specify this (eg kNormal | kIncludeSimplex) to
52  // improve the chances of escaping from invalid minima
53  kIncludeSimplex = 4,
54  // You may optionally specify these to improve the final error estimates
55  kIncludeHesse = 8,
56  kIncludeMinos = 16
57  };
Enumerator
kQuiet 
kVerbose 

Definition at line 40 of file Fit.h.

Constructor & Destructor Documentation

ana::Fitter::Fitter ( const IExperiment expt,
std::vector< const IFitVar * >  vars,
std::vector< const ISyst * >  systs = {},
Precision  prec = kNormal 
)

Definition at line 53 of file Fit.cxx.

57  : fExpt(expt), fVars(vars), fSysts(systs), fPrec(prec), fCalc(0),
58  fCovar(0), fCovarStatus(-1)
59  {
60  }
std::vector< const IFitVar * > fVars
Definition: Fit.h:173
Precision fPrec
Definition: Fit.h:175
TMatrixDSym * fCovar
Definition: Fit.h:186
bool fCovarStatus
Definition: Fit.h:187
osc::IOscCalcAdjustable * fCalc
Definition: Fit.h:176
std::vector< const ISyst * > fSysts
Definition: Fit.h:174
int prec
const IExperiment * fExpt
Definition: Fit.h:172

Member Function Documentation

Fitter* ana::Fitter::Clone ( ) const
inlineoverride

Definition at line 140 of file Fit.h.

141  {
142  std::cout << "Fitter::Clone() not implemented" << std::endl;
143  abort();
144  }
BEGIN_PROLOG could also be cout
void ana::Fitter::DecodePars ( const double *  pars) const
protected

Updates mutable fCalc and fShifts.

Definition at line 338 of file Fit.cxx.

339  {
340  for(unsigned int i = 0; i < fVars.size(); ++i){
341  const double val = pars[i];
342  if(isnan(val) || isinf(val)){
343  std::cout << "Fitter::DecodePars(). Minuit wants to set " << fVars[i]->ShortName() << " = " << val << ". Will ignore it and hope for the best." << std::endl;
344  // To debug this problem uncomment this abort and run in gdb
345  // abort();
346  }
347  else{
348  fVars[i]->SetValue(fCalc, val);
349  }
350  }
351  for(unsigned int j = 0; j < fSysts.size(); ++j){
352  const double val = pars[fVars.size()+j];
353  if(isnan(val) || isinf(val)){
354  std::cout << "Fitter::DecodePars(). Minuit wants to set " << fSysts[j]->ShortName() << " = " << val << ". Will ignore it and hope for the best." << std::endl;
355  // To debug this problem uncomment this abort and run in gdb
356  // abort();
357  }
358  else{
359  fShifts.SetShift(fSysts[j], val);
360  }
361  }
362  }
void SetShift(const ISyst *syst, double shift)
Shifts are 0=nominal, -1,+1 = 1 sigma shifts.
Definition: SystShifts.cxx:47
std::vector< const IFitVar * > fVars
Definition: Fit.h:173
osc::IOscCalcAdjustable * fCalc
Definition: Fit.h:176
SystShifts fShifts
Definition: Fit.h:177
std::vector< const ISyst * > fSysts
Definition: Fit.h:174
BEGIN_PROLOG could also be cout
double ana::Fitter::DoEval ( const double *  pars) const
overridevirtual

Evaluate the log-likelihood, as required by MINUT interface.

Definition at line 365 of file Fit.cxx.

366  {
367  ++fNEval;
368 
369  DecodePars(pars); // Updates fCalc and fShifts
370 
371  // Have to re-fetch the values because DecodePars() will have truncated
372  // values to the physical range where necessary.
373  double penalty = 0;
374  for(unsigned int i = 0; i < fVars.size(); ++i){
375  const double x = pars[i];
376  if(isinf(x) || isnan(x)) continue; // DecodePars() should have warned
377  penalty += fVars[i]->Penalty(x, fCalc);
378  }
379  for(unsigned int j = 0; j < fSysts.size(); ++j){
380  const double x = pars[fVars.size()+j];
381  if(isinf(x) || isnan(x)) continue; // DecodePars() should have warned
382  penalty += fSysts[j]->Penalty(x);
383  }
384 
385  if(fNEval%1000 == 0) std::cout << fNEval << ": EXPT chi2 = " << fExpt->ChiSq(fCalc, fShifts) << "; penalty = " << penalty << std::endl;
386  return fExpt->ChiSq(fCalc, fShifts) + penalty;
387  }
process_name opflash particleana ie x
virtual double ChiSq(osc::IOscCalcAdjustable *osc, const SystShifts &syst=SystShifts::Nominal()) const =0
int fNEval
Definition: Fit.h:179
std::vector< const IFitVar * > fVars
Definition: Fit.h:173
osc::IOscCalcAdjustable * fCalc
Definition: Fit.h:176
SystShifts fShifts
Definition: Fit.h:177
void DecodePars(const double *pars) const
Updates mutable fCalc and fShifts.
Definition: Fit.cxx:338
std::vector< const ISyst * > fSysts
Definition: Fit.h:174
const IExperiment * fExpt
Definition: Fit.h:172
BEGIN_PROLOG could also be cout
std::vector< Fitter::SeedPt > ana::Fitter::ExpandSeeds ( const std::map< const IFitVar *, std::vector< double >> &  seedPts,
std::vector< SystShifts systSeedPts 
) const
protected

Definition at line 64 of file Fit.cxx.

67  {
68  std::vector<SeedPt> ret;
69  ret.push_back(SeedPt());
70 
71  for(auto it: seedPts){
72  // For every variable, duplicate every entry in ret with the value set to
73  // each possibility.
74  const IFitVar* fv = it.first;
75  std::vector<SeedPt> newret;
76  for(double val: it.second){
77  for(SeedPt pt: ret){
78  pt.fitvars[fv] = val;
79  newret.push_back(pt);
80  }
81  } // end for val
82  ret = newret;
83  } // end for it
84 
85  // Now duplicate as many times as required for the syst seeds
86  if(!systSeedPts.empty()){
87  std::vector<SeedPt> newret;
88  for(const SystShifts& s: systSeedPts){
89  for(SeedPt pt: ret){
90  pt.shift = s;
91  newret.push_back(pt);
92  }
93  }
94  ret = newret;
95  }
96 
97  return ret;
98  }
then echo File list $list not found else cat $list while read file do echo $file sed s
Definition: file_to_url.sh:60
double ana::Fitter::Fit ( osc::IOscCalcAdjustable seed,
SystShifts bestSysts = junkShifts,
const std::map< const IFitVar *, std::vector< double >> &  seedPts = {},
const std::vector< SystShifts > &  systSeedPts = {},
Verbosity  verb = kVerbose 
) const
Parameters
[out]seedSeed parameter and output best-fit point
[out]bestSystsBest systematics result returned here
seedPtsSet each var to each of the values. Try all combinations. Beware of combinatorical explosion...
systSeedPtsIf non-empty, try fit starting at each of these
verbIf quiet, no printout
Returns
-2x the log-likelihood of the best-fit point

Definition at line 276 of file Fit.cxx.

281  {
282  // Fits with no osc calc shouldn't be trying to optimize any
283  // oscilation parameters...
284  assert(seed || fVars.empty());
285 
286  for(const auto& it: seedPts){
287  if(std::find(fVars.begin(), fVars.end(), it.first) == fVars.end()){
288  std::cout << "ERROR Fitter::Fit() trying to seed '"
289  << it.first->ShortName()
290  << "' which is not part of the fit." << std::endl;
291  abort();
292  }
293  }
294  for(const auto& it: systSeedPts){
295  for(const ISyst* s: it.ActiveSysts()){
296  if(std::find(fSysts.begin(), fSysts.end(), s) == fSysts.end()){
297  std::cout << "ERROR Fitter::Fit() trying to seed '"
298  << s->ShortName()
299  << "' which is not part of the fit." << std::endl;
300  abort();
301  }
302  }
303  }
304 
305  if(verb == kVerbose){
306  std::cout << "Finding best fit for";
307  for(const IFitVar* v: fVars) std::cout << " " << v->ShortName();
308  for(const ISyst* s: fSysts) std::cout << " " << s->ShortName();
309  std::cout << "..." << std::endl;
310  }
311 
312  // Do all the actual work. This wrapper function is just so we can have
313  // better control of printing.
314  const double chi = FitHelper(seed, bestSysts, seedPts, systSeedPts, verb);
315 
316  if(verb == kVerbose){
317  std::cout << "Best fit";
318  for(const IFitVar* v: fVars){
319  std::cout << ", " << v->ShortName() << " = " << v->GetValue(seed);
320  }
321  for(const ISyst* s: fSysts){
322  std::cout << ", " << s->ShortName() << " = " << bestSysts.GetShift(s);
323  }
324  std::cout << ", LL = " << chi << std::endl;
325 
326  std::cout << " found with " << fNEval << " evaluations of the likelihood";
327  if(fNEvalFiniteDiff > 0)
328  std::cout << " (" << fNEvalFiniteDiff << " to evaluate gradients numerically)";
329  if(fNEvalGrad > 0)
330  std::cout << " and " << fNEvalGrad << " evaluations of the gradient";
331  std::cout << std::endl;
332  }
333 
334  return chi;
335  }
double FitHelper(osc::IOscCalcAdjustable *seed, SystShifts &bestSysts, const std::map< const IFitVar *, std::vector< double >> &seedPts, std::vector< SystShifts > systSeedPts, Verbosity verb) const
Helper for Fit.
Definition: Fit.cxx:176
int fNEvalGrad
Definition: Fit.h:180
int fNEval
Definition: Fit.h:179
int fNEvalFiniteDiff
Definition: Fit.h:181
std::vector< const IFitVar * > fVars
Definition: Fit.h:173
std::vector< const ISyst * > fSysts
Definition: Fit.h:174
then echo File list $list not found else cat $list while read file do echo $file sed s
Definition: file_to_url.sh:60
BEGIN_PROLOG could also be cout
double ana::Fitter::Fit ( osc::IOscCalcAdjustable seed,
SystShifts systSeed,
Verbosity  verb 
) const
inline

Variant with no seedPts.

Definition at line 79 of file Fit.h.

82  {
83  return Fit(seed, systSeed, {}, std::vector<SystShifts>(1, systSeed), verb);
84  }
double Fit(osc::IOscCalcAdjustable *seed, SystShifts &bestSysts=junkShifts, const std::map< const IFitVar *, std::vector< double >> &seedPts={}, const std::vector< SystShifts > &systSeedPts={}, Verbosity verb=kVerbose) const
Definition: Fit.cxx:276
double ana::Fitter::Fit ( osc::IOscCalcAdjustable seed,
Verbosity  verb 
) const
inline

Variant with no seedPts and no systematics result returned.

Definition at line 87 of file Fit.h.

89  {
90  return Fit(seed, junkShifts, {}, {}, verb);
91  }
double Fit(osc::IOscCalcAdjustable *seed, SystShifts &bestSysts=junkShifts, const std::map< const IFitVar *, std::vector< double >> &seedPts={}, const std::vector< SystShifts > &systSeedPts={}, Verbosity verb=kVerbose) const
Definition: Fit.cxx:276
double ana::Fitter::Fit ( SystShifts systSeed,
Verbosity  verb = kVerbose 
) const
inline

Variant with no oscillations - useful for ND fits.

Definition at line 94 of file Fit.h.

95  {
96  return Fit(0, systSeed, {}, std::vector<SystShifts>(1, systSeed), verb);
97  }
double Fit(osc::IOscCalcAdjustable *seed, SystShifts &bestSysts=junkShifts, const std::map< const IFitVar *, std::vector< double >> &seedPts={}, const std::vector< SystShifts > &systSeedPts={}, Verbosity verb=kVerbose) const
Definition: Fit.cxx:276
double ana::Fitter::FitHelper ( osc::IOscCalcAdjustable seed,
SystShifts bestSysts,
const std::map< const IFitVar *, std::vector< double >> &  seedPts,
std::vector< SystShifts systSeedPts,
Verbosity  verb 
) const
protected

Helper for Fit.

Definition at line 176 of file Fit.cxx.

181  {
182  const std::vector<SeedPt> pts = ExpandSeeds(seedPts, systSeedPts);
183  double minchi = 1e10;
184  std::vector<double> bestFitPars, bestSystPars;
185 
186  for(const SeedPt& pt: pts){
187  osc::IOscCalcAdjustable *seed = initseed->Copy();
188 
189  for(auto it: pt.fitvars) it.first->SetValue(seed, it.second);
190 
191  // Need to deal with parameters that are not fit values!
192  SystShifts shift = pt.shift;
193 
194  // Need to copy over syst values into the seed for this fit
195  // that were in the input syst shifts, but are not being fit for,
196  // and therefore not part of the seedPts list
197  for(auto s: bestSysts.ActiveSysts()) {
198  auto fit_systs = shift.ActiveSysts();
199  if(std::find(fit_systs.begin(), fit_systs.end(), s) == fit_systs.end()) {
200  shift.SetShift(s, bestSysts.GetShift(s));
201  }
202  }
203 
204  std::unique_ptr<ROOT::Math::Minimizer> thisMin =
205  FitHelperSeeded(seed, shift, verb);
206 
207  // Check whether this is the best minimum we've found so far
208  if (thisMin->MinValue() < minchi) {
209  minchi = thisMin->MinValue();
210 
211  // Need to do Prefit here too...
212  fParamNames .clear();
213  fPreFitValues .clear();
214  fPreFitErrors .clear();
215  fPostFitValues .clear();
216  fPostFitErrors .clear();
217  fMinosErrors .clear();
219 
220  // Save pre-fit info
221  // In principle, these can change with the seed, so have to save them here...
222  for(const IFitVar* v: fVars){
223  const double val = v->GetValue(seed);
224  fParamNames .push_back(v->ShortName());
225  fPreFitValues.push_back(val);
226  fPreFitErrors.push_back(val ? val/2 : .1);
227  }
228  for(const ISyst* s: fSysts){
229  const double val = shift.GetShift(s);
230  fParamNames .push_back(s->ShortName());
231  fPreFitValues.push_back(val);
232  fPreFitErrors.push_back(1);
233  }
234 
235  // Now save postfit
236  fPostFitValues = std::vector<double>(thisMin->X(), thisMin->X()+thisMin->NDim());
237  fPostFitErrors = std::vector<double>(thisMin->Errors(), thisMin->Errors()+thisMin->NDim());
238 
239  fEdm = thisMin->Edm();
240  fIsValid = !thisMin->Status();
241 
242  delete fCovar;
243  fCovar = new TMatrixDSym(thisMin->NDim());
244 
245  for(unsigned int row = 0; row < thisMin->NDim(); ++row){
246  for(unsigned int col = 0; col < thisMin->NDim(); ++col){
247  (*fCovar)(row, col) = thisMin->CovMatrix(row, col);
248  }
249  }
250 
251  // Store the best fit values of all the parameters we know are being varied.
252  for(unsigned int i = 0; i < fVars.size(); ++i)
253  bestFitPars.push_back(fPostFitValues[i]);
254  for(unsigned int j = 0; j < fSysts.size(); ++j)
255  bestSystPars.push_back(fPostFitValues[fVars.size()+j]);
256  }
257  delete seed;
258  }
259 
260  // Stuff the results of the actual best fit back into the seeds
261  for(unsigned int i = 0; i < fVars.size(); ++i)
262  fVars[i]->SetValue(initseed, bestFitPars[i]);
263  for(unsigned int i = 0; i < fSysts.size(); ++i)
264  bestSysts.SetShift(fSysts[i], bestSystPars[i]);
265 
266  return minchi;
267  }
std::vector< std::pair< double, double > > fMinosErrors
Definition: Fit.h:193
std::vector< std::pair< double, double > > fTempMinosErrors
Definition: Fit.h:194
std::vector< SeedPt > ExpandSeeds(const std::map< const IFitVar *, std::vector< double >> &seedPts, std::vector< SystShifts > systSeedPts) const
Definition: Fit.cxx:64
double fEdm
Definition: Fit.h:184
shift
Definition: fcl_checks.sh:26
std::vector< double > fPreFitValues
Definition: Fit.h:189
std::vector< const IFitVar * > fVars
Definition: Fit.h:173
bool fIsValid
Definition: Fit.h:185
TMatrixDSym * fCovar
Definition: Fit.h:186
unsigned int seed
std::unique_ptr< ROOT::Math::Minimizer > FitHelperSeeded(osc::IOscCalcAdjustable *seed, SystShifts &systSeed, Verbosity verb) const
Helper for FitHelper.
Definition: Fit.cxx:102
std::vector< double > fPostFitValues
Definition: Fit.h:191
std::vector< double > fPreFitErrors
Definition: Fit.h:190
std::vector< const ISyst * > fSysts
Definition: Fit.h:174
std::vector< std::string > fParamNames
Definition: Fit.h:188
then echo File list $list not found else cat $list while read file do echo $file sed s
Definition: file_to_url.sh:60
then echo ***************************************echo Variable FHICL_FILE_PATH not found echo You porbably haven t set up larsoft echo Try setup uboonecode vXX_XX_XX q e10
Definition: find_fhicl.sh:6
std::vector< double > fPostFitErrors
Definition: Fit.h:192
std::unique_ptr< ROOT::Math::Minimizer > ana::Fitter::FitHelperSeeded ( osc::IOscCalcAdjustable seed,
SystShifts systSeed,
Verbosity  verb 
) const
protected

Helper for FitHelper.

Definition at line 102 of file Fit.cxx.

104  {
105  // Why, when this is called for each seed?
106  fCalc = seed;
107  fShifts = systSeed;
108 
109  if(fPrec & kIncludeSimplex) std::cout << "Simplex option specified but not implemented. Ignored." << std::endl;
110  if((fPrec & kAlgoMask) == kGradDesc){
111  std::cout << "GradDesc option specified but not implemented. Abort." << std::endl;
112  abort();
113  }
114 
115  std::unique_ptr<ROOT::Math::Minimizer> mnMin(
116  ROOT::Math::Factory::CreateMinimizer("Minuit2", "Combined"));
117  mnMin->SetStrategy(int(fPrec & kAlgoMask));
118  mnMin->SetMaxFunctionCalls(1E8);
119  mnMin->SetMaxIterations(1E8);
120 
121  for(const IFitVar* v: fVars){
122  const double val = v->GetValue(seed);
123  // name, value, error
124  mnMin->SetVariable(mnMin->NFree(), v->ShortName(), val, val ? fabs(val/2) : .1);
125  fParamNames.push_back(v->ShortName());
126  fPreFitValues.push_back(val);
127  fPreFitErrors.push_back(val ? val/2 : .1);
128  }
129  // One way this can go wrong is if two variables have the same ShortName
130  assert(mnMin->NFree() == fVars.size());
131  for(const ISyst*s: fSysts){
132  const double val = systSeed.GetShift(s);
133  // name, value, error
134  mnMin->SetVariable(mnMin->NFree(), s->ShortName(), val, 1);
135  fParamNames.push_back(s->ShortName());
136  fPreFitValues.push_back(val);
137  fPreFitErrors.push_back(1);
138  }
139  // One way this can go wrong is if two variables have the same ShortName
140  assert(mnMin->NFree() == fVars.size()+fSysts.size());
141 
142  mnMin->SetFunction(*this);
143 
144  if(verb == Verbosity::kQuiet) mnMin->SetPrintLevel(0);
145 
146  if(!mnMin->Minimize()){
147  std::cout << "*** ERROR: minimum is not valid ***" << std::endl;
148  std::cout << "*** Precision: " << mnMin->Precision() << std::endl;
149 
150  std::cout << "-- Stopped at: \n";
151  for(uint i = 0; i < mnMin->NFree(); ++i){
152  std::cout << "\t" << mnMin->VariableName(i) << " = " << mnMin->X()[i] << "\n";
153  }
154  }
155 
156  if (fPrec & kIncludeHesse){
157  std::cout << "It's Hesse o'clock" << std::endl;
158  mnMin->Hesse();
159  }
160 
161  if (fPrec & kIncludeMinos){
162  std::cout << "It's minos time" << std::endl;
163  fTempMinosErrors.clear();
164  for (uint i = 0; i < mnMin->NDim(); ++i){
165  double errLow = 0, errHigh = 0;
166  mnMin->GetMinosError(i, errLow, errHigh);
167  std::cout << i << "/" << mnMin->NDim() << " " << fParamNames[i] << ": " << errLow << ", +" << errHigh << " (" << mnMin->Errors()[i] << ")" << std::endl;
168  fTempMinosErrors.push_back(std::make_pair(errLow,errHigh));
169  }
170  }
171 
172  return mnMin;
173  }
std::vector< std::pair< double, double > > fTempMinosErrors
Definition: Fit.h:194
std::vector< double > fPreFitValues
Definition: Fit.h:189
std::vector< const IFitVar * > fVars
Definition: Fit.h:173
Precision fPrec
Definition: Fit.h:175
unsigned int seed
std::vector< double > fPreFitErrors
Definition: Fit.h:190
osc::IOscCalcAdjustable * fCalc
Definition: Fit.h:176
SystShifts fShifts
Definition: Fit.h:177
std::vector< const ISyst * > fSysts
Definition: Fit.h:174
std::vector< std::string > fParamNames
Definition: Fit.h:188
then echo File list $list not found else cat $list while read file do echo $file sed s
Definition: file_to_url.sh:60
BEGIN_PROLOG could also be cout
TMatrixDSym* ana::Fitter::GetCovariance ( )
inline

Return the fit covariance.

Definition at line 100 of file Fit.h.

100 { return this->fCovar;}
TMatrixDSym * fCovar
Definition: Fit.h:186
int ana::Fitter::GetCovarianceStatus ( )
inline

covariance matrix status (0 = not valid, 1 approximate, 2, full but made pos def, 3 accurate and not pos def)

Definition at line 103 of file Fit.h.

103 {return this->fCovarStatus;}
bool fCovarStatus
Definition: Fit.h:187
double ana::Fitter::GetEDM ( )
inline

Return edm form the fit.

Definition at line 127 of file Fit.h.

127 {return this->fEdm;}
double fEdm
Definition: Fit.h:184
bool ana::Fitter::GetIsValid ( )
inline

Say whether the fit was good.

Definition at line 130 of file Fit.h.

130 {return this->fIsValid;}
bool fIsValid
Definition: Fit.h:185
std::vector<std::pair<double,double> > ana::Fitter::GetMinosErrors ( )
inline

Return the minos errors.

Definition at line 121 of file Fit.h.

121 { return this->fMinosErrors;}
std::vector< std::pair< double, double > > fMinosErrors
Definition: Fit.h:193
int ana::Fitter::GetNFCN ( )
inline

Return number of function calls.

Definition at line 124 of file Fit.h.

124 {return this->fNEval;}
int fNEval
Definition: Fit.h:179
std::vector<std::string> ana::Fitter::GetParamNames ( )
inline

Return the fit names.

Definition at line 106 of file Fit.h.

106 { return this->fParamNames;}
std::vector< std::string > fParamNames
Definition: Fit.h:188
std::vector<double> ana::Fitter::GetPostFitErrors ( )
inline

Return the postfit errors.

Definition at line 118 of file Fit.h.

118 { return this->fPostFitErrors;}
std::vector< double > fPostFitErrors
Definition: Fit.h:192
std::vector<double> ana::Fitter::GetPostFitValues ( )
inline

Return the postfit values.

Definition at line 115 of file Fit.h.

115 { return this->fPostFitValues;}
std::vector< double > fPostFitValues
Definition: Fit.h:191
std::vector<double> ana::Fitter::GetPreFitErrors ( )
inline

Return the prefit errors.

Definition at line 112 of file Fit.h.

112 { return this->fPreFitErrors;}
std::vector< double > fPreFitErrors
Definition: Fit.h:190
std::vector<double> ana::Fitter::GetPreFitValues ( )
inline

Return the prefit values.

Definition at line 109 of file Fit.h.

109 { return this->fPreFitValues;}
std::vector< double > fPreFitValues
Definition: Fit.h:189
SystShifts ana::Fitter::GetSystShifts ( ) const
inline

Definition at line 132 of file Fit.h.

132 {return fShifts;}
SystShifts fShifts
Definition: Fit.h:177
virtual unsigned int ana::Fitter::NDim ( ) const
inlineoverridevirtual

Definition at line 138 of file Fit.h.

138 {return fVars.size()+fSysts.size();}
std::vector< const IFitVar * > fVars
Definition: Fit.h:173
std::vector< const ISyst * > fSysts
Definition: Fit.h:174
void ana::Fitter::SetPrecision ( Precision  prec)

Definition at line 270 of file Fit.cxx.

271  {
272  fPrec = prec;
273  }
Precision fPrec
Definition: Fit.h:175
int prec

Member Data Documentation

osc::IOscCalcAdjustable* ana::Fitter::fCalc
mutableprotected

Definition at line 176 of file Fit.h.

TMatrixDSym* ana::Fitter::fCovar
mutableprotected

Definition at line 186 of file Fit.h.

bool ana::Fitter::fCovarStatus
mutableprotected

Definition at line 187 of file Fit.h.

double ana::Fitter::fEdm = -1
mutableprotected

Definition at line 184 of file Fit.h.

const IExperiment* ana::Fitter::fExpt
protected

Definition at line 172 of file Fit.h.

bool ana::Fitter::fIsValid = false
mutableprotected

Definition at line 185 of file Fit.h.

std::vector<std::pair<double, double> > ana::Fitter::fMinosErrors
mutableprotected

Definition at line 193 of file Fit.h.

int ana::Fitter::fNEval = 0
mutableprotected

Definition at line 179 of file Fit.h.

int ana::Fitter::fNEvalFiniteDiff = 0
mutableprotected

Definition at line 181 of file Fit.h.

int ana::Fitter::fNEvalGrad = 0
mutableprotected

Definition at line 180 of file Fit.h.

std::vector<std::string> ana::Fitter::fParamNames
mutableprotected

Definition at line 188 of file Fit.h.

std::vector<double> ana::Fitter::fPostFitErrors
mutableprotected

Definition at line 192 of file Fit.h.

std::vector<double> ana::Fitter::fPostFitValues
mutableprotected

Definition at line 191 of file Fit.h.

Precision ana::Fitter::fPrec = kNormal
protected

Definition at line 175 of file Fit.h.

std::vector<double> ana::Fitter::fPreFitErrors
mutableprotected

Definition at line 190 of file Fit.h.

std::vector<double> ana::Fitter::fPreFitValues
mutableprotected

Definition at line 189 of file Fit.h.

SystShifts ana::Fitter::fShifts
mutableprotected

Definition at line 177 of file Fit.h.

std::vector<const ISyst*> ana::Fitter::fSysts
protected

Definition at line 174 of file Fit.h.

std::vector<std::pair<double, double> > ana::Fitter::fTempMinosErrors
mutableprotected

Definition at line 194 of file Fit.h.

std::vector<const IFitVar*> ana::Fitter::fVars
protected

Definition at line 173 of file Fit.h.


The documentation for this class was generated from the following files: