Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes | List of all members
cmf::CovarianceFitHelper Class Reference

#include "/cvmfs/nova-development.opensciencegrid.org/novasoft/releases/N20-11-28/CovarianceMatrixFit/modules/CovarianceFitHelper.h"

Public Member Functions

void AddIteration (double const &chi2, std::vector< std::string > const &pars, std::vector< double > const &vals, cmf::ChiSqrCalc_t const &chiSqrCalc=cmf::kCovMat)
 
void Initialize (fhicl::ParameterSet const &helperPars, fhicl::ParameterSet const &chiSqrCalcPars, fhicl::ParameterSet const &manipulatorPars, fhicl::ParameterSet const &minimizerPars, fhicl::ParameterSet const &sawPars)
 
void Initialize (fhicl::ParameterSet const &chiSqrCalcPars, fhicl::ParameterSet const &sawPars)
 
void Initialize (fhicl::ParameterSet const &helperPars, fhicl::ParameterSet const &chiSqrCalcPars, fhicl::ParameterSet const &manipulatorPars, fhicl::ParameterSet const &minimizerPars, fhicl::ParameterSet const &sawPars, cmf::Spectrum const &dataSpectrum)
 
void FillDataSpectrum (cmf::EventListColl const &eventLists)
 
void FillDataSpectrum (std::vector< float > const &dataSpectrum)
 
void FillMCSpectrum (cmf::EventListColl const &eventLists)
 
void FindInitialGuess (std::vector< cmf::OscParamPoint > const &library)
 
void MakeResultPlots (std::string const &dirBaseName, std::string const &uniqueID=std::string(), bool fillFromLists=true)
 
void Minimize (cmf::ChiSqrCalc_t const &chiSqrCalc=cmf::kCovMat)
 
std::unique_ptr< cmf::PointResultResult (bool loadFromEventList=true)
 
cmf::Spectrum const & BestFitSpectrum (bool loadFromEventList=true)
 

Static Public Member Functions

static CovarianceFitHelperInstance ()
 

Private Member Functions

 CovarianceFitHelper ()
 
 ~CovarianceFitHelper ()
 
void InitializeEventLists (fhicl::ParameterSet const &manipulatorPars, cmf::Spectrum const &dataSpectrum)
 
void InitializeFitHelper (fhicl::ParameterSet const &helperPars)
 
void InitializeFit ()
 
void InitializeGlobalVars ()
 
void InitializeMinimizer (fhicl::ParameterSet const &minimizerPars)
 
void InitializeShifterAndWeighter (fhicl::ParameterSet const &sawPars)
 
void MinimizeCP (double const &dCPLow, double const &dCPHigh, std::vector< double > &fitVals)
 
void Make1DSpectra (art::TFileDirectory &tfd, cmf::EventListColl const &mcLists, std::string const &uniqueID)
 
void MakeCovarianceMatrixHistogram (art::TFileDirectory &tfd, std::string const &uniqueID)
 
void MakeDataMCCanv (std::map< long, TH1D * > &dataSpecMap, std::map< long, TH1D * > &mcSpecMap, cmf::SelectionType_t const &selVec, std::string const &canvName, art::TFileDirectory &tfd, std::string const &uniqueID)
 
void MakeIterationGraphs (art::TFileDirectory &tfd, std::string const &uniqueID)
 

Private Attributes

int fNumCPRanges
 number of ranges to divide up d_CP for seeding fits More...
 
fhicl::ParameterSet fMinimizerPars
 alternative algorithm to try if minimization fails More...
 
art::ServiceHandle< art::TFileServicefTFS
 TFileService. More...
 
std::map< int, double > fDataPOT
 Data POT values for each detector/beam/selection. More...
 
std::vector< float > fDataSpectrum
 data spectrum in bins More...
 
std::vector< float > fMCSpectrum
 MC spectrum in bins. More...
 
std::vector< float > fDataCount
 number of data events in each bin More...
 
std::vector< float > fMCCount
 number of MC events in each bin More...
 
std::map< cmf::ChiSqrCalc_t, std::vector< cmf::CovFitIteration > > fIterations
 keep track of the minimizer iterations More...
 

Detailed Description

Definition at line 42 of file CovarianceFitHelper.h.

Constructor & Destructor Documentation

cmf::CovarianceFitHelper::CovarianceFitHelper ( )
private

Definition at line 132 of file CovarianceFitHelper.cxx.

References fDataCount, fDataSpectrum, fMCCount, fMCSpectrum, and cmf::CovarianceBinUtility::Instance().

Referenced by Instance().

133  {
134  fDataSpectrum.resize(cmf::CovarianceBinUtility::Instance()->TotalBins(), 0.);
135  fMCSpectrum .resize(cmf::CovarianceBinUtility::Instance()->TotalBins(), 0.);
136  fDataCount .resize(cmf::CovarianceBinUtility::Instance()->TotalBins(), 0.);
137  fMCCount .resize(cmf::CovarianceBinUtility::Instance()->TotalBins(), 0.);
138  }
std::vector< float > fMCSpectrum
MC spectrum in bins.
std::vector< float > fDataCount
number of data events in each bin
std::vector< float > fMCCount
number of MC events in each bin
std::vector< float > fDataSpectrum
data spectrum in bins
static CovarianceBinUtility * Instance()
cmf::CovarianceFitHelper::~CovarianceFitHelper ( )
private

Definition at line 141 of file CovarianceFitHelper.cxx.

142  {
143  // fIterations .clear();
144  // fDataSpectrum.clear();
145  // fMCSpectrum .clear();
146  // fDataCount .clear();
147  // fMCCount .clear();
148  // fDataPOT .clear();
149 
150  // gLoc->FDLoc .clear();
151  // gLoc->NDLoc .clear();
152  // gLoc .reset();
153  // gExposureMap.reset();
154  // gFFVals .reset();
155  // gFFNames .reset();
156  }

Member Function Documentation

void cmf::CovarianceFitHelper::AddIteration ( double const &  chi2,
std::vector< std::string > const &  pars,
std::vector< double > const &  vals,
cmf::ChiSqrCalc_t const &  chiSqrCalc = cmf::kCovMat 
)
inline

Definition at line 46 of file CovarianceFitHelper.h.

References TMVAGlob::Initialize(), cmf::kCovMat, and string.

Referenced by FindInitialGuess(), cmf::FitFunction_ForROOT(), and Minimize().

49  { fIterations[chiSqrCalc].emplace_back(pars, vals, chi2); }
double chi2()
std::map< cmf::ChiSqrCalc_t, std::vector< cmf::CovFitIteration > > fIterations
keep track of the minimizer iterations
std::string pars("Th23Dmsq32")
cmf::Spectrum const & cmf::CovarianceFitHelper::BestFitSpectrum ( bool  loadFromEventList = true)

Definition at line 1115 of file CovarianceFitHelper.cxx.

References FillMCSpectrum(), fMCSpectrum, and Instance().

Referenced by Result().

1116  {
1117  // fill the MC spectrum
1119 
1120  return fMCSpectrum;
1121  }
std::vector< float > fMCSpectrum
MC spectrum in bins.
static CovarianceFitHelper * Instance()
void FillMCSpectrum(cmf::EventListColl const &eventLists)
std::unique_ptr< cmf::EventListColl > gMCEventLists
void cmf::CovarianceFitHelper::FillDataSpectrum ( cmf::EventListColl const &  eventLists)

Definition at line 411 of file CovarianceFitHelper.cxx.

References fDataCount, fDataSpectrum, cmf::FillSpectrum(), MECModelEnuComparisons::i, cmf::ChiSqrCalculator::Instance(), and cmf::ChiSqrCalculator::SetSpectrum().

Referenced by InitializeEventLists(), and cmf::FitFeldmanCousinsPoint::writeResults().

412  {
413  // now fill the data spectrum - set the values to 0
414  // for all the elements as cmf::FillSpectrum does not reset the vectors
415  for(size_t i = 0; i < fDataSpectrum.size(); ++i){
416  fDataSpectrum[i] = 0;
417  fDataCount[i] = 0.;
418  }
419 
420  cmf::FillSpectrum(eventLists,
421  *gExposureMap,
423  fDataCount);
424 
426  // for(auto const& itr : *fDataSpectrum)
427  // LOG_VERBATIM("CovarianceFitHelper")
428  // << itr;
429  }
std::vector< float > fDataCount
number of data events in each bin
std::unique_ptr< cmf::ExposureMap > gExposureMap
void FillSpectrum(cmf::EventListColl const &eventLists, cmf::ExposureMap const &exposureMap, std::vector< float > &spectrum, std::vector< float > &count, cmf::InteractionType_t intType, bool applyExposureNorm)
void SetSpectrum(std::vector< float > const &spectrum, bool isData)
static cmf::ChiSqrCalculator * Instance()
std::vector< float > fDataSpectrum
data spectrum in bins
void cmf::CovarianceFitHelper::FillDataSpectrum ( std::vector< float > const &  dataSpectrum)

Definition at line 400 of file CovarianceFitHelper.cxx.

References fDataSpectrum, cmf::ChiSqrCalculator::Instance(), and cmf::ChiSqrCalculator::SetSpectrum().

401  {
402  fDataSpectrum = dataSpectrum;
403 
405  // for(auto const& itr : *fDataSpectrum)
406  // LOG_VERBATIM("CovarianceFitHelper")
407  // << itr;
408  }
void SetSpectrum(std::vector< float > const &spectrum, bool isData)
static cmf::ChiSqrCalculator * Instance()
std::vector< float > fDataSpectrum
data spectrum in bins
void cmf::CovarianceFitHelper::FillMCSpectrum ( cmf::EventListColl const &  eventLists)

Definition at line 432 of file CovarianceFitHelper.cxx.

References cmf::FillSpectrum(), fMCCount, fMCSpectrum, MECModelEnuComparisons::i, cmf::ChiSqrCalculator::Instance(), and cmf::ChiSqrCalculator::SetSpectrum().

Referenced by BestFitSpectrum(), cmf::FitFunction_ForROOT(), InitializeFit(), and Minimize().

433  {
434  // now fill the MC spectrum with the current point - set the values to 0
435  // for all the elements as cmf::FillSpectrum does not reset the vectors
436  for(size_t i = 0; i < fMCSpectrum.size(); ++i){
437  fMCSpectrum[i] = 0;
438  fMCCount[i] = 0.;
439  }
440 
441  cmf::FillSpectrum(eventLists,
442  *gExposureMap,
443  fMCSpectrum,
444  fMCCount);
445 
447 
448  }
std::vector< float > fMCSpectrum
MC spectrum in bins.
std::vector< float > fMCCount
number of MC events in each bin
std::unique_ptr< cmf::ExposureMap > gExposureMap
void FillSpectrum(cmf::EventListColl const &eventLists, cmf::ExposureMap const &exposureMap, std::vector< float > &spectrum, std::vector< float > &count, cmf::InteractionType_t intType, bool applyExposureNorm)
void SetSpectrum(std::vector< float > const &spectrum, bool isData)
static cmf::ChiSqrCalculator * Instance()
void cmf::CovarianceFitHelper::FindInitialGuess ( std::vector< cmf::OscParamPoint > const &  library)

Definition at line 451 of file CovarianceFitHelper.cxx.

References AddIteration(), cmf::ChiSqrCalculator::ChiSqr(), fMCSpectrum, cmf::gMinChiSqr, cmf::ChiSqrCalculator::Instance(), Instance(), cmf::ShifterAndWeighter::Instance(), cmf::IsOscillationParameter(), cmf::kGarbageDouble, LOG_VERBATIM, cmf::ShifterAndWeighter::SetCurrentVals(), cmf::ChiSqrCalculator::SetSpectrum(), string, cmf::StringToOscParmType(), and art::to_string().

Referenced by cmf::CovarianceMatrixFitter::writeResults(), and cmf::FitFeldmanCousinsPoint::writeResults().

452  {
453  // loop over the points in the library to get the initial guess
455  double chiSqr;
456  cmf::OscillationParameterMap minOscPars;
457 
458 // TStopwatch ts;
459 // ts.Start();
460 
461  for(auto const& itr : library.front().OscPoint())
462  minOscPars.emplace(itr.first, itr.second);
463 
464  cmf::OscillationParameterMap oscPars = minOscPars;
465 
466  // The data spectrum was set when the CovarianceFitHelper was initialized
467  // so now for each point in space set the MC spectrum for the ChiSqrCalculator
468  for(auto const& prediction : library){
469 
470  for(auto const& oscItr : prediction.OscPoint())
471  oscPars[oscItr.first] = oscItr.second;
472 
473  cmf::ChiSqrCalculator::Instance()->SetSpectrum(prediction.PredictedSpectrum(),
474  false);
475  chiSqr = cmf::ChiSqrCalculator::Instance()->ChiSqr(gChiSqrCalc, oscPars);
476  if(chiSqr < gMinChiSqr){
477  gMinChiSqr = chiSqr;
478  minOscPars.swap(oscPars);
479  fMCSpectrum = cmf::Spectrum(prediction.PredictedSpectrum());
480  } // end if the minimum chiSqr
481  } // end loop over the predictions
482 
483  // now set the initial guess to the parameters at the minimum
484  // we only care about the oscillation parameters here
485  std::string initialGuess("initial guess from grid search - ");
486  cmf::OscParm_t oscPar;
487  for(size_t p = 0; p < gFFNames->size(); ++p){
488 
489  if(!cmf::IsOscillationParameter((*gFFNames)[p])) continue;
490 
491  oscPar = cmf::StringToOscParmType((*gFFNames)[p]);
492 
493  if(minOscPars.count(oscPar) > 0)
494  (*gFFVals)[p] = minOscPars.find(oscPar)->second;
495 
496  initialGuess += (*gFFNames)[p] + ": " + std::to_string((*gFFVals)[p]) + "\t";
497  }
498 
499  // Now update the ShifterAndWeighter to check the goodness of fit
501  gFFVals->data());
502 
503  // store this iteration of the fitter
505  *gFFNames,
506  *gFFVals,
507  gChiSqrCalc);
508 
509  LOG_VERBATIM("CovarianceFitHelper")
510  << initialGuess
511  << " with chi^2: "
512  << gMinChiSqr;
513 // << " took "
514 // << ts.RealTime()
515 // << " to go through "
516 // << library.size()
517 // << " library spectra";
518  }
enum cmf::osc_params OscParm_t
const char * p
Definition: xmltok.h:285
static const double kGarbageDouble
Definition: Constants.h:22
std::vector< float > fMCSpectrum
MC spectrum in bins.
std::map< cmf::OscParm_t, float > OscillationParameterMap
Definition: Constants.h:612
static ShifterAndWeighter * Instance()
void SetSpectrum(std::vector< float > const &spectrum, bool isData)
std::vector< float > Spectrum
Definition: Constants.h:610
static cmf::ChiSqrCalculator * Instance()
double ChiSqr(cmf::ChiSqrCalc_t const &chiSqrCalc, std::vector< std::string > const &pars, const double *vals)
static CovarianceFitHelper * Instance()
static bool IsOscillationParameter(std::string const &str)
void AddIteration(double const &chi2, std::vector< std::string > const &pars, std::vector< double > const &vals, cmf::ChiSqrCalc_t const &chiSqrCalc=cmf::kCovMat)
std::unique_ptr< std::vector< double > > gFFVals
#define LOG_VERBATIM(category)
std::string to_string(ModuleType mt)
Definition: ModuleType.h:32
double gMinChiSqr
cmf::ChiSqrCalc_t gChiSqrCalc
void SetCurrentVals(cmf::Location const &loc)
static cmf::OscParm_t StringToOscParmType(std::string const &str)
std::unique_ptr< std::vector< std::string > > gFFNames
enum BeamMode string
void cmf::CovarianceFitHelper::Initialize ( fhicl::ParameterSet const &  helperPars,
fhicl::ParameterSet const &  chiSqrCalcPars,
fhicl::ParameterSet const &  manipulatorPars,
fhicl::ParameterSet const &  minimizerPars,
fhicl::ParameterSet const &  sawPars 
)

Definition at line 159 of file CovarianceFitHelper.cxx.

References cmf::ChiSqrCalculator::InitializeCovarianceMatrix(), InitializeEventLists(), InitializeFit(), InitializeFitHelper(), InitializeGlobalVars(), InitializeMinimizer(), InitializeShifterAndWeighter(), and cmf::ChiSqrCalculator::Instance().

Referenced by cmf::FitFeldmanCousinsPoint::beginJob(), and cmf::CovarianceMatrixFitter::writeResults().

164  {
166 
167  this->InitializeGlobalVars();
168  this->InitializeFitHelper(helperPars);
169  this->InitializeShifterAndWeighter(sawPars);
170  this->InitializeEventLists(manipulatorPars,
171  cmf::Spectrum());
172  this->InitializeMinimizer(minimizerPars);
173  this->InitializeFit();
174  }
void InitializeCovarianceMatrix(fhicl::ParameterSet const &pset)
void InitializeMinimizer(fhicl::ParameterSet const &minimizerPars)
void InitializeEventLists(fhicl::ParameterSet const &manipulatorPars, cmf::Spectrum const &dataSpectrum)
void InitializeFitHelper(fhicl::ParameterSet const &helperPars)
void InitializeShifterAndWeighter(fhicl::ParameterSet const &sawPars)
std::vector< float > Spectrum
Definition: Constants.h:610
static cmf::ChiSqrCalculator * Instance()
void cmf::CovarianceFitHelper::Initialize ( fhicl::ParameterSet const &  chiSqrCalcPars,
fhicl::ParameterSet const &  sawPars 
)

Definition at line 177 of file CovarianceFitHelper.cxx.

References cmf::ChiSqrCalculator::InitializeCovarianceMatrix(), InitializeGlobalVars(), InitializeShifterAndWeighter(), and cmf::ChiSqrCalculator::Instance().

179  {
181  this->InitializeGlobalVars();
182  this->InitializeShifterAndWeighter(sawPars);
183  }
void InitializeCovarianceMatrix(fhicl::ParameterSet const &pset)
void InitializeShifterAndWeighter(fhicl::ParameterSet const &sawPars)
static cmf::ChiSqrCalculator * Instance()
void cmf::CovarianceFitHelper::Initialize ( fhicl::ParameterSet const &  helperPars,
fhicl::ParameterSet const &  chiSqrCalcPars,
fhicl::ParameterSet const &  manipulatorPars,
fhicl::ParameterSet const &  minimizerPars,
fhicl::ParameterSet const &  sawPars,
cmf::Spectrum const &  dataSpectrum 
)

Definition at line 186 of file CovarianceFitHelper.cxx.

References cmf::ChiSqrCalculator::InitializeCovarianceMatrix(), InitializeEventLists(), InitializeFit(), InitializeFitHelper(), InitializeGlobalVars(), InitializeMinimizer(), InitializeShifterAndWeighter(), and cmf::ChiSqrCalculator::Instance().

192  {
194 
195  this->InitializeGlobalVars();
196  this->InitializeFitHelper(helperPars);
197  this->InitializeShifterAndWeighter(sawPars);
198  this->InitializeEventLists(manipulatorPars,
199  dataSpectrum);
200  this->InitializeMinimizer(minimizerPars);
201  this->InitializeFit();
202  }
void InitializeCovarianceMatrix(fhicl::ParameterSet const &pset)
void InitializeMinimizer(fhicl::ParameterSet const &minimizerPars)
void InitializeEventLists(fhicl::ParameterSet const &manipulatorPars, cmf::Spectrum const &dataSpectrum)
void InitializeFitHelper(fhicl::ParameterSet const &helperPars)
void InitializeShifterAndWeighter(fhicl::ParameterSet const &sawPars)
static cmf::ChiSqrCalculator * Instance()
void cmf::CovarianceFitHelper::InitializeEventLists ( fhicl::ParameterSet const &  manipulatorPars,
cmf::Spectrum const &  dataSpectrum 
)
private

Definition at line 217 of file CovarianceFitHelper.cxx.

References cmf::EventListManipulator::Deserialize(), cmf::EventListManipulator::ExposureMap(), fDataPOT, FillDataSpectrum(), cmf::kMC, and LOG_DEBUG.

Referenced by Initialize().

219  {
220 
221  if(!dataSpectrum.empty()){
222  this->FillDataSpectrum(dataSpectrum);
223  }
224 
225  // get the manipulator to deserialize the data and MC
226  // we'll keep the MC lists around, but will drop the data lists as
227  // soon as we leave this method
228  cmf::EventListManipulator elm(manipulatorPars);
229  gExposureMap = std::make_unique<cmf::ExposureMap>(elm.ExposureMap());
230  gMCEventLists = std::make_unique<cmf::EventListColl>();
231 
232  LOG_DEBUG("CovarianceFitHelper")
233  << "Exposure map has "
234  << gExposureMap->size()
235  << " entries";
236 
237  for(auto const& itr : *gExposureMap){
238  LOG_DEBUG("CovarianceFitHelper")
239  << itr.first.ToString()
240  << " has "
241  << itr.second
242  << " POT";
243 
244  fDataPOT[itr.first.DetectorBeamSelectionKey() + itr.first.PeriodKey()] = itr.second.goodPOT;
245  }
246 
247  // deserialize the ND and FD
248  elm.Deserialize(*gMCEventLists, cmf::kMC);
249  }
#define LOG_DEBUG(stream)
Definition: Messenger.h:149
std::unique_ptr< cmf::ExposureMap > gExposureMap
std::map< int, double > fDataPOT
Data POT values for each detector/beam/selection.
std::unique_ptr< cmf::EventListColl > gMCEventLists
void FillDataSpectrum(cmf::EventListColl const &eventLists)
void cmf::CovarianceFitHelper::InitializeFit ( )
private

Definition at line 252 of file CovarianceFitHelper.cxx.

References FillMCSpectrum().

Referenced by Initialize().

253  {
254  // initial fill of the MC spectrum
256 
257  // for(size_t b = 0; b< fMCSpectrum.size(); ++b)
258  // LOG_VERBATIM("CovarianceFitHelper")
259  // << "Initial MC Spectrum "
260  // << b
261  // << " "
262  // << fMCSpectrum[b];
263 
264  }
void FillMCSpectrum(cmf::EventListColl const &eventLists)
std::unique_ptr< cmf::EventListColl > gMCEventLists
void cmf::CovarianceFitHelper::InitializeFitHelper ( fhicl::ParameterSet const &  helperPars)
private

Definition at line 267 of file CovarianceFitHelper.cxx.

References fNumCPRanges, fhicl::ParameterSet::get(), LOG_DEBUG, and fhicl::ParameterSet::to_indented_string().

Referenced by Initialize().

268  {
269  LOG_DEBUG("CovarianceFitHelper")
270  << helperPars.to_indented_string();
271 
272  fNumCPRanges = helperPars.get<int >("NumCPRanges", 3);
273  }
#define LOG_DEBUG(stream)
Definition: Messenger.h:149
int fNumCPRanges
number of ranges to divide up d_CP for seeding fits
void cmf::CovarianceFitHelper::InitializeGlobalVars ( )
private

Definition at line 276 of file CovarianceFitHelper.cxx.

References confusionMatrixTree::count, cmf::ParameterUtility::Instance(), cmf::kGarbageDouble, LOG_DEBUG, cet::sqlite::max(), and cmf::ParameterUtility::ParametersAsLocation().

Referenced by Initialize().

277  {
278  gFFVals = std::make_unique<std::vector<double> >();
279  gFFNames = std::make_unique<std::vector<std::string>>();
282  gLoc = std::make_unique<cmf::Location>(cmf::ParameterUtility::Instance()->ParametersAsLocation());
283 
284  // Add all FD entry keys to the psl_key_list
285  // Add ND entry keys only if they don't show up on FD list
286  // E.g., FD entry will override ND entry if common
287  // Don't add fixed parameters or systematic uncertainties because
288  // those are never fit
289  std::list<std::string> psl_key_list_fd;
290  for(auto const& psl : gLoc->FDLocation()){
291  if(psl.second.IsFixed() ||
292  !psl.second.IsOscPar() ) continue;
293  LOG_DEBUG("ROOTFactoryFit")
294  << "FD Adding "
295  << psl.second;
296  gFFVals ->emplace_back(psl.second.Value());
297  gFFNames->emplace_back(psl.first );
298 
299  psl_key_list_fd.emplace_back(psl.first);
300  }
301 
302  // The ND list will contain an entry only if it occurs solely for the ND
303  // and not for the FD
304  for(auto const& psl : gLoc->NDLocation() ){
305  if(psl.second.IsFixed() ||
306  !psl.second.IsOscPar() ||
307  std::count(psl_key_list_fd.begin(), psl_key_list_fd.end(), psl.first) > 0 ) continue;
308 
309  LOG_DEBUG("ROOTFactoryFit")
310  << "ND Adding "
311  << psl.first;
312  gFFVals ->emplace_back(psl.second.Value());
313  gFFNames->emplace_back(psl.first );
314  }
315 
316  }
#define LOG_DEBUG(stream)
Definition: Messenger.h:149
static const double kGarbageDouble
Definition: Constants.h:22
static ParameterUtility * Instance()
std::unique_ptr< cmf::Location > gLoc
std::unique_ptr< std::vector< double > > gFFVals
cmf::Location ParametersAsLocation()
double gMinChiSqr
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
Definition: statistics.h:68
std::unique_ptr< std::vector< std::string > > gFFNames
void cmf::CovarianceFitHelper::InitializeMinimizer ( fhicl::ParameterSet const &  minimizerPars)
private

Definition at line 319 of file CovarianceFitHelper.cxx.

References MakeMiniprodValidationCuts::f, cmf::FitFunction_ForROOT(), fMinimizerPars, fhicl::ParameterSet::get(), MECModelEnuComparisons::i, cmf::ParameterUtility::Instance(), cmf::Parameter::IsConstrained(), LOG_VERBATIM, cmf::ParameterUtility::ParameterInfo(), and string.

Referenced by Initialize().

320  {
321  fMinimizerPars = minimizerPars;
322 
323  // Look for minName and algoName (following ROOT factory nomenclature,
324  // see NumericalMinimization.C in ROOT fit documentation,
325  // root.cern.ch/root/html/tutorials/fit/NumericalMinimization.C.html
326  // Default here is set to Minuit2 && Migrad
327  auto minName = fMinimizerPars.get<std::string >("MinName", "Minuit2");
328  auto algoName = fMinimizerPars.get<std::string >("AlgoName", "Migrad");
329  auto maxFunctionCalls = fMinimizerPars.get<unsigned int >("MaxFunctionCalls", 1000);
330  auto maxIterations = fMinimizerPars.get<unsigned int >("MaxIterations", 1000);
331  auto tolerance = fMinimizerPars.get<double >("Tolerance", 0);
332  auto precision = fMinimizerPars.get<double >("Precision", 0);
333  auto strategyCode = fMinimizerPars.get<int >("StrategyCode", 0);
334  auto printLevel = fMinimizerPars.get<int >("MINUITPrintLevel", 3);
335 
336  // make the ROOT minimizer
337  gFFMin = ROOT::Math::Factory::CreateMinimizer(minName,
338  algoName);
339 
340  // Set up the minimizer particulars
341  gFFMin->SetMaxFunctionCalls(maxFunctionCalls);
342  gFFMin->SetMaxIterations (maxIterations );
343  gFFMin->SetStrategy (strategyCode );
344  gFFMin->SetPrintLevel (printLevel );
345  if(tolerance > 0) gFFMin->SetTolerance(tolerance);
346  if(precision > 0) gFFMin->SetPrecision(precision);
347 
348  for(size_t ipar = 0; ipar < gFFNames->size(); ++ipar){
350 
351  gFFMin->SetLimitedVariable(ipar,
352  (*gFFNames)[ipar],
353  (*gFFVals)[ipar],
354  0.01,
355  cmf::ParameterUtility::Instance()->ParameterInfo((*gFFNames)[ipar]).LowerBound(),
356  cmf::ParameterUtility::Instance()->ParameterInfo((*gFFNames)[ipar]).UpperBound());
357  }
358  else{
359  gFFMin->SetVariable(ipar,
360  (*gFFNames)[ipar],
361  (*gFFVals)[ipar],
362  0.01);
363  } // end if unconstrained parameter
364  } // end loop to set the initial parameters
365 
366  LOG_VERBATIM("CovarianceMatrixFit")
367  << "\nNum floatables: "
368  << gFFNames->size()
369  << "\nMax function calls: "
370  << maxFunctionCalls
371  << "\nMax iterations: "
372  << maxIterations
373  << "\nTolerance: "
374  << tolerance
375  << "\nInitial guess: ";
376 
377  ROOT::Fit::ParameterSettings parSet;
378  for(size_t i = 0; i < gFFNames->size(); ++i){
379  gFFMin->GetVariableSettings(i, parSet);
380  LOG_VERBATIM("CovarianceMatrixFit")
381  << parSet.Name()
382  << " : "
383  << parSet.Value()
384  << " lower limit ("
385  << parSet.HasLowerLimit()
386  << ") "
387  << parSet.LowerLimit()
388  << " upper limit ("
389  << parSet.HasUpperLimit()
390  << ") "
391  << parSet.UpperLimit();
392  }
393 
394  // Link in the function to be minimized
395  auto * f = new ROOT::Math::Functor(&FitFunction_ForROOT, gFFNames->size());
396  gFFMin->SetFunction(*f);
397  }
cmf::Parameter const & ParameterInfo(std::string const &parName)
bool const & IsConstrained() const
Definition: Parameter.h:109
double FitFunction_ForROOT(const double *params)
fhicl::ParameterSet fMinimizerPars
alternative algorithm to try if minimization fails
static ParameterUtility * Instance()
T get(std::string const &key) const
Definition: ParameterSet.h:231
ROOT::Math::Minimizer * gFFMin
std::unique_ptr< std::vector< double > > gFFVals
#define LOG_VERBATIM(category)
std::unique_ptr< std::vector< std::string > > gFFNames
enum BeamMode string
void cmf::CovarianceFitHelper::InitializeShifterAndWeighter ( fhicl::ParameterSet const &  sawPars)
private

Definition at line 205 of file CovarianceFitHelper.cxx.

References cmf::ShifterAndWeighter::InitShiftsAndWeightsToUse(), cmf::ParameterUtility::Instance(), cmf::ShifterAndWeighter::Instance(), PandAna.Demos.demo0::loc, cmf::ParameterUtility::ParametersAsLocation(), and cmf::ShifterAndWeighter::SetCurrentVals().

Referenced by Initialize().

206  {
208 
210  sawPars);
211 
213  //cmf::ShifterAndWeighter::Instance()->ReportCurrentVals();
214  }
static ParameterUtility * Instance()
static ShifterAndWeighter * Instance()
void InitShiftsAndWeightsToUse(cmf::Location const &loc, fhicl::ParameterSet const &parset)
cmf::Location ParametersAsLocation()
void SetCurrentVals(cmf::Location const &loc)
CovarianceFitHelper * cmf::CovarianceFitHelper::Instance ( )
static
void cmf::CovarianceFitHelper::Make1DSpectra ( art::TFileDirectory tfd,
cmf::EventListColl const &  mcLists,
std::string const &  uniqueID 
)
private

Definition at line 616 of file CovarianceFitHelper.cxx.

References b, bins, cmf::cBeamType_LatexStrings, cmf::cBeamType_Strings, cmf::cFileTypeStrings, cmf::cSelectionType_LatexStrings, cmf::cSelectionType_Strings, evt, fDataPOT, fDataSpectrum, fMCSpectrum, compareCafs::histName, cmf::CovarianceBinUtility::Instance(), cmf::PlotUtilities::Instance(), cmf::ShifterAndWeighter::Instance(), cmf::kNu_RecoE, art::TFileDirectory::make(), fetch_tb_beamline_files::md, norm, cmf::PlotUtilities::NormalizeBinContents(), cmf::CovarianceBinUtility::SelectionHistBinning(), string, and cmf::ShifterAndWeighter::Weight().

Referenced by MakeResultPlots().

619  {
620  size_t numBins = fDataSpectrum.size();
621 
622  if(numBins != fMCSpectrum.size())
623  throw cet::exception("CovarianceFitHelper")
624  << "data and MC spectra sizes are different: data - "
625  << numBins
626  << " mc - "
627  << fMCSpectrum.size();
628 
629  auto data1DHist = tfd.make<TH1D>(("dataSpectrumBins" + uniqueID).c_str(),
630  ";Logical Bin;Events",
631  numBins,
632  0,
633  numBins);
634 
635  auto mc1DHist = tfd.make<TH1D>(("mcSpectrumBins" + uniqueID).c_str(),
636  ";Logical Bin;Events",
637  numBins,
638  0,
639  numBins);
640 
641  for(size_t b = 0; b < numBins; ++b){
642  data1DHist->SetBinContent(b+1, fDataSpectrum[b]);
643  mc1DHist ->SetBinContent(b+1, fMCSpectrum[b]);
644  }
645 
646  // fill the energy spectra for each event selection, file type, beam and detector
647  // data and MC. The MC can be made from the event lists. We dropped the data
648  // from the event lists to save memory, so use the bin vector to make the data
649  // plots
650 
652  std::string histTitle;
653  std::map<std::string, TH1D*> nameToSpectrum;
654  std::vector<double> bins;
655  double norm;
656  TH1D *spectrum;
657 
658  for(auto const& itr : mcLists){
659  auto const& md = itr.ListMetaData();
660 
661  bins.clear();
663 
664  histName = ("SpectraEnergy" +
665  uniqueID +
666  cmf::cBeamType_Strings[md.BeamType()] +
667  md.DetectorString() +
668  cmf::cFileTypeStrings[md.fileType] +
669  cmf::cSelectionType_Strings[md.selectionType]);
670 
671  histTitle = (cmf::cBeamType_LatexStrings[md.BeamType()] + " " +
672  md.DetectorString() + " " +
673  cmf::cSelectionType_LatexStrings[md.selectionType] +
674  ";Energy (GeV);Events / 0.1 GeV");
675 
676 // LOG_DEBUG("CovarianceFitHelper")
677 // << "1D spectrum hist is "
678 // << histName;
679 
680 // for(auto const& binitr : bins)
681 // LOG_VERBATIM("CovarianceFitHelper")
682 // << "bin edge: "
683 // << binitr;
684 
685  if(nameToSpectrum.count(histName) < 1)
686  spectrum = nameToSpectrum.emplace(histName, tfd.make<TH1D>(histName.c_str(),
687  histTitle.c_str(),
688  bins.size() - 1,
689  &bins[0])).first->second;
690  else
691  spectrum = nameToSpectrum.find(histName)->second;
692  // get the normalization to use
693 
694  norm = 1.;
695  if(md.isMC){
696  norm = fDataPOT[md.DetectorBeamSelectionKey() +
697  md.PeriodKey()]/itr.ListSpillSummary().goodPOT;
698  }
699 
700  // fill the histogram with the proper weight and normalization
701  for(auto const& evt : itr){
702  spectrum->Fill(evt->DataVals().val_at(cmf::kNu_RecoE, md),
704  } // end loop over events
705 
706  // now normalize the bin contents if desired to the desired width
707  if(md.IsNuESelected())
709  else if(md.IsNuMuSelected())
711 
712  } // end loop over event lists
713 
714  }
void SelectionHistBinning(cmf::MetaData const &md, std::vector< double > &bins)
const std::string cSelectionType_Strings[12]
Definition: Constants.h:79
::xsd::cxx::tree::exception< char > exception
Definition: Database.h:225
std::vector< float > fMCSpectrum
MC spectrum in bins.
double Weight(cmf::EventPtr const &curEvent, cmf::MetaData const &md)
const std::string cFileTypeStrings[10]
Definition: Constants.h:163
int evt
const Binning bins
static ShifterAndWeighter * Instance()
std::map< int, double > fDataPOT
Data POT values for each detector/beam/selection.
T * make(ARGS...args) const
Float_t norm
std::vector< float > fDataSpectrum
data spectrum in bins
const hit & b
Definition: hits.cxx:21
static PlotUtilities * Instance()
void NormalizeBinContents(TH1D *hist, double normToVal)
const std::string cBeamType_LatexStrings[4]
Definition: Constants.h:40
const std::string cBeamType_Strings[4]
Definition: Constants.h:35
const std::string cSelectionType_LatexStrings[12]
Definition: Constants.h:92
static CovarianceBinUtility * Instance()
enum BeamMode string
void cmf::CovarianceFitHelper::MakeCovarianceMatrixHistogram ( art::TFileDirectory tfd,
std::string const &  uniqueID 
)
private

Definition at line 717 of file CovarianceFitHelper.cxx.

References plot_validation_datamc::c, cmf::ChiSqrCalculator::ChiSqr(), fDataSpectrum, fMCSpectrum, cmf::ChiSqrCalculator::Instance(), cmf::kCovMat, cmf::kStatMatrix, cmf::kSystMatrix, cmf::kTotalMatrix, art::TFileDirectory::make(), r(), and cmf::ChiSqrCalculator::SetSpectrum().

Referenced by MakeResultPlots().

719  {
720  // we need to fill the covariance matrix so that we have the correct
721  // output histogram for it. This method is called after the minimizations
722  // finish and we have the final fit values.
723  // Set the MC spectra in the ChiSqrCalculator to be sure it is at the best fit,
724  // ie final iteration
725  // then find the chi^2 to ensure the matrix is properly loaded
726  // just pass an empty nuisance parametermap, we aren't really trying to
727  // find the chi^2 anyway here
730  *gFFNames,
731  gFFVals->data());
732 
733 
734  auto matrixHistStat = tfd.make<TH2D>(("FitHelperStatCovarianceMatrix" + uniqueID).c_str(),
735  ";Logical Bin;Logical Bin",
736  fDataSpectrum.size(),
737  0,
738  fDataSpectrum.size(),
739  fDataSpectrum.size(),
740  0,
741  fDataSpectrum.size());
742 
743  auto matrixHistSyst = tfd.make<TH2D>(("FitHelperSystCovarianceMatrix" + uniqueID).c_str(),
744  ";Logical Bin;Logical Bin",
745  fDataSpectrum.size(),
746  0,
747  fDataSpectrum.size(),
748  fDataSpectrum.size(),
749  0,
750  fDataSpectrum.size());
751 
752  auto matrixHistTot = tfd.make<TH2D>(("FitHelperTotCovarianceMatrix" + uniqueID).c_str(),
753  ";Logical Bin;Logical Bin",
754  fDataSpectrum.size(),
755  0,
756  fDataSpectrum.size(),
757  fDataSpectrum.size(),
758  0,
759  fDataSpectrum.size());
760 
761  for(size_t r = 0; r < fDataSpectrum.size(); ++r){
762  for(size_t c = 0; c < fDataSpectrum.size(); ++c){
763  //if(this->ValidVectPosition(fCovarianceMatrixTot.size(), r * fDataSpectrum.size() + c))
764  matrixHistTot ->SetBinContent(r+1, c+1, cmf::ChiSqrCalculator::Instance()->CovarianceMatrixEntry(r * fDataSpectrum.size() + c, cmf::kStatMatrix));
765  matrixHistSyst->SetBinContent(r+1, c+1, cmf::ChiSqrCalculator::Instance()->CovarianceMatrixEntry(r * fDataSpectrum.size() + c, cmf::kSystMatrix));
766  matrixHistStat->SetBinContent(r+1, c+1, cmf::ChiSqrCalculator::Instance()->CovarianceMatrixEntry(r * fDataSpectrum.size() + c, cmf::kTotalMatrix));
767  }
768  }
769  }
std::vector< float > fMCSpectrum
MC spectrum in bins.
void SetSpectrum(std::vector< float > const &spectrum, bool isData)
static cmf::ChiSqrCalculator * Instance()
double ChiSqr(cmf::ChiSqrCalc_t const &chiSqrCalc, std::vector< std::string > const &pars, const double *vals)
T * make(ARGS...args) const
std::vector< float > fDataSpectrum
data spectrum in bins
std::unique_ptr< std::vector< double > > gFFVals
TRandom3 r(0)
std::unique_ptr< std::vector< std::string > > gFFNames
void cmf::CovarianceFitHelper::MakeDataMCCanv ( std::map< long, TH1D * > &  dataSpecMap,
std::map< long, TH1D * > &  mcSpecMap,
cmf::SelectionType_t const &  selVec,
std::string const &  canvName,
art::TFileDirectory tfd,
std::string const &  uniqueID 
)
private

Definition at line 774 of file CovarianceFitHelper.cxx.

References bins, canv, cmf::MetaData::detector, cmf::DetectorBeamSelectionTypesToKey(), cmf::CovarianceBinUtility::Instance(), cmf::SelectionUtility::Instance(), findDuplicateFiles::key, kRed, LOG_DEBUG, art::TFileDirectory::make(), cet::sqlite::max(), std::max(), fetch_tb_beamline_files::md, cmf::SelectionUtility::SelectionsUsedOfType(), cmf::MetaData::selectionType, std::sqrt(), and cmf::MetaData::ToString().

Referenced by MakeResultPlots().

780  {
781 
783  cmf::SelectionUtility::Instance()->SelectionsUsedOfType(selType, subsetOfSels);
784 
785  // make sure we have selections of the type we want
786  if(subsetOfSels.empty()) return;
787 
788  // figure out how many selection types we have and how many beams
789  size_t numSels = 0;
790  int cnt;
791  std::map<cmf::SelectionType_t, int> selToNum;
792  std::set<cmf::BeamType_t> beams;
793  for(auto const& itr : subsetOfSels){
794  beams.insert(itr.BeamType());
795  if(itr.Selections().size() > numSels){
796  numSels = itr.Selections().size();
797  cnt = 0;
798  for(auto const& selItr : itr.Selections()){
799  selToNum.emplace(selItr, cnt);
800  ++cnt;
801  }
802  }
803  }
804 
805  size_t numBeams = beams.size();
806 
807  // make 1 canvas each for the general selections, numu, nue and NC data/MC comparison
808  // The spectra maps have what we want and we can use the keys
809  // to determine how to make the canvases
810 
811  // for now we need pads for the ND and FD, FHC and RHC
812  // arrange them ND above FD, FHC grouped above RHC
813 
814  // We need a pad for each selection in each detector for each beam type.
815  // Arrange them ND above FD, FHC grouped above RHC.
816  // The total number of canvases is then 2 * 2 * selVec.size()
817  // We want selVec.size() columns of divisions and 2 * 2 rows
818  // Since we want an even number of divisions, the TCanvas::Divide function
819  // works for us here
820  auto canv = tfd.make<TCanvas>((canvName + uniqueID).c_str(),
821  canvName.c_str(),
822  1200,
823  1200);
824  canv->Divide(numSels, 2 * numBeams);
825 
827  TH1D *dataHist = nullptr;
828  TH1D *mcHist = nullptr;
829  long canvDiv;
830  long key;
831  std::vector<double> bins;
832  TVirtualPad *curPad = nullptr;
833  TLegend *curLeg = nullptr;
834 
835  for(auto const& dbsItr : subsetOfSels){
836  for(auto const& selItr : dbsItr.Selections()){
837 
838  md.selectionType = selItr;
839  md.detector = dbsItr.Detector();
840 
841  // make vectors of the bin boundaries for the new histogram
842  bins.clear();
843  bins.push_back(0.);
844  for(auto const& itr : cmf::CovarianceBinUtility::Instance()->SelectionHighEdges(md)){
845  bins.push_back(itr.first);
846  LOG_DEBUG("CovarianceFitHelper")
847  << md.ToString()
848  << " bins for total MC and ratios "
849  << bins.back()
850  << " "
851  << bins.size();
852  }
853 
854  canvDiv = (1 + selToNum.find(selItr)->second) + (numBeams * dbsItr.BeamType() * numSels) + (dbsItr.Detector() - 1) * numSels;
855 
856  curPad = canv->cd(canvDiv);
857 
858  key = cmf::DetectorBeamSelectionTypesToKey(dbsItr.Detector(), dbsItr.BeamType(), selItr);
859 
860  // grab the data and MC spectra for this key
861  dataHist = dataSpecMap.find(key)->second;
862  mcHist = mcSpecMap.find(key)->second;
863 
864  double max = std::max(dataHist->GetMaximum(), mcHist->GetMaximum());
865  mcHist->SetMaximum(1.2 * (max + std::sqrt(max)));
866 
867  mcHist->SetLineColor(kRed);
868 
869  mcHist->Draw("hist");
870  dataHist->Draw("psame");
871 
872  if(selItr == *(dbsItr.Selections().rbegin())){
873  curLeg = curPad->BuildLegend();
874  curLeg->SetBorderSize(0);
875  curLeg->SetFillStyle(0);
876  }
877 
878  canv->Update();
879 
880  } // end loop over selections
881  } // end loop over DetBeamSelSet
882 
883  // write out the canvas
884  canv->Write();
885 
886  delete curLeg;
887  }
T max(const caf::Proxy< T > &a, T b)
#define LOG_DEBUG(stream)
Definition: Messenger.h:149
enum BeamMode kRed
cmf::DetType_t detector
Definition: Structs.h:114
T sqrt(T number)
Definition: d0nt_math.hpp:156
static long DetectorBeamSelectionTypesToKey(cmf::DetType_t const &det, cmf::BeamType_t const &bt, cmf::SelectionType_t const &sel)
Definition: StaticFuncs.h:62
static SelectionUtility * Instance()
TCanvas * canv
cmf::SelectionType_t selectionType
Definition: Structs.h:116
const Binning bins
std::string ToString() const
Definition: Structs.cxx:133
T * make(ARGS...args) const
void SelectionsUsedOfType(cmf::SelectionType_t const &sel, DetBeamSelSet &sels)
std::set< cmf::SelectionUtility::DetBeamSels > DetBeamSelSet
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
Definition: statistics.h:68
static CovarianceBinUtility * Instance()
void cmf::CovarianceFitHelper::MakeIterationGraphs ( art::TFileDirectory tfd,
std::string const &  uniqueID 
)
private

Definition at line 890 of file CovarianceFitHelper.cxx.

References allTimeWatchdog::can, cmf::cChiSqrCalcType_Strings, fIterations, cmf::IsOscillationParameter(), art::TFileDirectory::make(), allInOneTrainingPlots::pad, and string.

Referenced by MakeResultPlots().

892  {
893  // loop over the different possible chi^2 calculators
894  std::string calcName;
895  for(auto const& itr : fIterations){
896 
897  calcName = cmf::cChiSqrCalcType_Strings[itr.first];
898  std::vector<cmf::CovFitIteration> const& cfiVec = itr.second;
899 
900  // create histograms for the fit parameters and the chi^2 to see
901  // how they changed during the minimization
902 
903  // make a collection of TGraphs for the fit parameters
904  std::map<std::string, TGraph *> grMap;
905  for(auto const& itr : cfiVec[0].parNames){
906 
907  if(!cmf::IsOscillationParameter(itr)) continue;
908 
909  grMap[itr] = tfd.make<TGraph>(cfiVec.size());
910  grMap[itr]->SetName((itr + calcName + "IterGraph" + uniqueID).c_str());
911  grMap[itr]->SetTitle(itr.c_str());
912  grMap[itr]->SetMarkerStyle(20);
913  }
914 
915  // Make the TGraph for the chi^2 values
916  auto chisq_graph = tfd.make<TGraph>(cfiVec.size());
917  chisq_graph->SetName((calcName + "IterGraphChiSqr" + uniqueID).c_str());
918  chisq_graph->SetTitle("#chi^{2}");
919  chisq_graph->SetMarkerStyle(20);
920 
921  //loop over the iterations to fill the graphs
922  // the last value in fIterations is not the minimum, it is
923  // the value that Minuit used to confirm it had found a minimum,
924  // so don't use it in the graphs
925  double parVal;
926  int ctr = 0;
927  for(auto const& itr : cfiVec){
928 
929  // get the value for each parameter we care about and put it into its
930  // graph
931  for(size_t p = 0; p < itr.parNames.size(); ++p){
932  parVal = itr.parVals[p];
933 
934  grMap.find(itr.parNames[p])->second->SetPoint(ctr, ctr * 1., parVal);
935  }
936 
937  // fill the chi^2 graph
938  chisq_graph->SetPoint(ctr, ctr * 1., itr.chiSqr);
939  ++ctr;
940  } // end loop over iterations
941 
942  // now loop over the graphs and create a canvas for each
943  std::string canvasname = calcName + uniqueID + "_ParameterIterationCanv";
944  TCanvas *can = tfd.make<TCanvas>(canvasname.data(), canvasname.data(), 800, 800);
945  can->Divide((grMap.size() / 2) + 1, (grMap.size() / 2) + 1);
946 
947  size_t pad = 1;
948  for(auto mapItr : grMap){
949 
950  can->cd(pad);
951  mapItr.second->Draw("ACP");
952  ++pad;
953  }//end of loop over parameter graphs
954 
955  // now do the chi^2 values for the iterations
956  can->cd(pad);
957  chisq_graph->Draw("ACP");
958 
959  // Save the canvas to file
960  can->Update();
961  can->Write();
962 
963  } // end loop over the different chi^2 calculators
964  }
const std::string cChiSqrCalcType_Strings[5]
Definition: Constants.h:55
const char * p
Definition: xmltok.h:285
std::map< cmf::ChiSqrCalc_t, std::vector< cmf::CovFitIteration > > fIterations
keep track of the minimizer iterations
static bool IsOscillationParameter(std::string const &str)
T * make(ARGS...args) const
enum BeamMode string
void cmf::CovarianceFitHelper::MakeResultPlots ( std::string const &  dirBaseName,
std::string const &  uniqueID = std::string(),
bool  fillFromLists = true 
)

Definition at line 521 of file CovarianceFitHelper.cxx.

References fDataCount, fDataSpectrum, cmf::ChiSqrCalculator::FillChiSqrHistograms(), fMCCount, fMCSpectrum, fTFS, cmf::CovarianceBinUtility::Instance(), cmf::ChiSqrCalculator::Instance(), cmf::PlotUtilities::Instance(), cmf::kNCSelection, cmf::kNuESelection, cmf::kNumOscParams, cmf::kNuMuSelection, LOG_DEBUG, Make1DSpectra(), MakeCovarianceMatrixHistogram(), MakeDataMCCanv(), cmf::PlotUtilities::MakeEnergySpectraFromBins(), MakeIterationGraphs(), art::TFileDirectory::mkdir(), and cmf::CovarianceBinUtility::TotalBins().

Referenced by cmf::CovarianceMatrixFitter::writeResults(), and cmf::FitFeldmanCousinsPoint::writeResults().

524  {
525 
526  auto tfd = fTFS->mkdir(dirBaseName + uniqueID);
527 
528  this->MakeIterationGraphs (tfd, uniqueID);
529  this->MakeCovarianceMatrixHistogram(tfd, uniqueID);
530  if(fillFromLists)
531  this->Make1DSpectra(tfd, *gMCEventLists, uniqueID);
532 
533  std::map<long, TH1D*> totMCSpectra;
534  std::map<long, TH1D*> dataSpectra;
535 
536  cmf::PlotUtilities::Instance()->MakeEnergySpectraFromBins(fDataSpectrum, fDataCount, dataSpectra, "DataEnergySpectrum", tfd, uniqueID);
537  cmf::PlotUtilities::Instance()->MakeEnergySpectraFromBins(fMCSpectrum, fMCCount, totMCSpectra, "TotMCEnergySpectrum", tfd, uniqueID);
538 
539  // make the chi^2 histograms
540  auto binByBin = tfd.make<TH1D>("binByBinChiSqr",
541  ";Logical Bin;#chi^{2} per Bin;",
543  0,
545 
546  auto cumulativeChiSqr = tfd.make<TH1D>("cumulativeChiSqr",
547  ";Logical Bin;Cumulative #chi^{2}",
549  0,
551 
554  fMCSpectrum,
555  binByBin,
556  cumulativeChiSqr);
557 
558 
559  // make 1 canvas each for the general selections, numu, nue and NC data/MC comparison
560  // The spectra maps have what we want and we can use the keys
561  // to determine how to make the canvases
562 
563  // if a desired selection type is not available from the SelectionUtility,
564  // the call to MakeDataMCCanv is a no-op
565  std::set<std::pair<cmf::SelectionType_t, std::string>> selToString({{cmf::kNuMuSelection, "NuMuResultsCanv"},
566  {cmf::kNuESelection, "NuEResultsCanv"},
567  {cmf::kNCSelection, "NCResultsCanv"},
568  });
569 
570  for(auto const& itr : selToString){
572  totMCSpectra,
573  itr.first,
574  itr.second,
575  tfd,
576  uniqueID);
577  }
578 
579  // Finally create a TTree with the output of the oscillation parameters at
580  // the best fit
581 
582  auto *fitResTree = tfd.make<TTree>("FitResultTree", "Best Fit Values");
583 
584  std::vector<double> parVals (cmf::kNumOscParams, 0.);
585  std::vector<int> parFixed(cmf::kNumOscParams, 1);
586 
587  // now fill the tree
588  fitResTree->Branch("ChiSqr", &gMinChiSqr);
589 
590  for(auto const& itr : gLoc->FDLocation()){
591  if(!itr.second.IsOscPar()) continue;
592  parVals [itr.second.Key()] = itr.second.Value();
593  parFixed[itr.second.Key()] = itr.second.IsFixed();
594 
595  fitResTree->Branch((itr.first + "Val").c_str(), &parVals [itr.second.Key()]);
596  fitResTree->Branch((itr.first + "Fix").c_str(), &parFixed[itr.second.Key()]);
597 
598  LOG_DEBUG("CovarianceFitHelper")
599  << "added branches "
600  << itr.first + "Val"
601  << " "
602  << parVals[itr.second.Key()]
603  << " "
604  << itr.first + "Fix"
605  << " "
606  << parFixed[itr.second.Key()]
607  << " to tree";
608  }
609 
610  fitResTree->Fill();
611  }
#define LOG_DEBUG(stream)
Definition: Messenger.h:149
void MakeEnergySpectraFromBins(std::vector< float > const &binSpec, std::vector< float > const &binCount, std::map< long, TH1D * > &energySpecMap, std::string const &name, art::TFileDirectory &tfd, std::string const &uniqueID)
TFileDirectory mkdir(std::string const &dir, std::string const &descr="")
std::vector< float > fMCSpectrum
MC spectrum in bins.
std::vector< float > fDataCount
number of data events in each bin
art::ServiceHandle< art::TFileService > fTFS
TFileService.
void MakeCovarianceMatrixHistogram(art::TFileDirectory &tfd, std::string const &uniqueID)
std::vector< float > fMCCount
number of MC events in each bin
void Make1DSpectra(art::TFileDirectory &tfd, cmf::EventListColl const &mcLists, std::string const &uniqueID)
void FillChiSqrHistograms(cmf::ChiSqrCalc_t const &chiSqrCalc, std::vector< float > const &dataSpectrum, std::vector< float > const &mcSpectrum, TH1D *binByBin, TH1D *cumulative)
static cmf::ChiSqrCalculator * Instance()
size_t TotalBins(bool allSels=false)
void MakeDataMCCanv(std::map< long, TH1D * > &dataSpecMap, std::map< long, TH1D * > &mcSpecMap, cmf::SelectionType_t const &selVec, std::string const &canvName, art::TFileDirectory &tfd, std::string const &uniqueID)
std::unique_ptr< cmf::EventListColl > gMCEventLists
std::unique_ptr< cmf::Location > gLoc
std::vector< float > fDataSpectrum
data spectrum in bins
static PlotUtilities * Instance()
double gMinChiSqr
cmf::ChiSqrCalc_t gChiSqrCalc
void MakeIterationGraphs(art::TFileDirectory &tfd, std::string const &uniqueID)
static CovarianceBinUtility * Instance()
void cmf::CovarianceFitHelper::Minimize ( cmf::ChiSqrCalc_t const &  chiSqrCalc = cmf::kCovMat)

Definition at line 967 of file CovarianceFitHelper.cxx.

References AddIteration(), cmf::ChiSqrCalculator::ChiSqr(), FillMCSpectrum(), fNumCPRanges, cmf::gMinChiSqr, MECModelEnuComparisons::i, cmf::ChiSqrCalculator::Instance(), cmf::ShifterAndWeighter::Instance(), cmf::kGarbageDouble, LOG_VERBATIM, LOG_WARNING, MinimizeCP(), cmf::ShifterAndWeighter::SetCurrentVals(), string, and registry_explorer::v.

Referenced by cmf::CovarianceMatrixFitter::writeResults(), and cmf::FitFeldmanCousinsPoint::writeResults().

968  {
969  gChiSqrCalc = chiSqrCalc;
970 
971  // get the initial chi^2
974  *gFFNames,
975  gFFVals->data());
976 
977  LOG_VERBATIM("CovarianceMatrixFitter")
978  << "initial chi^2: "
979  << gMinChiSqr;
980 
981  // store this iteration of the fitter
982  this->AddIteration(gMinChiSqr,
983  *gFFNames,
984  *gFFVals,
985  gChiSqrCalc);
986 
987  std::vector<double> fitVals(gFFMin->NDim(), 0.);
988 
989  // if we are minimizing d_CP we have just minimized allowing the
990  // value of d_CP to vary between 0 < d_CP < pi, now try it
991  // with pi < d_CP < 2pi
992  if(!(gLoc->FDLocation().find(std::string("dCP"))->second.IsFixed())){
993 
994  double dCPRange = 2. * gPI / (1. * fNumCPRanges);
995 
996  for(int i = 0; i < fNumCPRanges; ++i){
997  this->MinimizeCP(i * dCPRange,
998  (i + 1) * dCPRange,
999  fitVals);
1000  }
1001  } // end if we have to minimize over ranges of d_CP
1002  else{
1003  // the minimizer returns true if we converge
1004  if( !gFFMin->Minimize() ){
1005  LOG_WARNING("CovarianceMatrixFit")
1006  << "Fit failed to converge with the "
1007  << gFFMin->Options().MinimizerAlgorithm()
1008  << " algorithm";
1009  }
1010 
1011  // grab the results of the minimization
1012  gMinChiSqr = gFFMin->MinValue();
1013  gMinStatus = gFFMin->Status();
1014  for(size_t v = 0; v < gFFVals->size(); ++v) (*gFFVals)[v] = gFFMin->X()[v];
1015  } // end if dCP is fixed
1016 
1017  // Now update the ShifterAndWeighter to check the goodness of fit
1019 
1020  LOG_VERBATIM("CovarianceMatrixFitter")
1021  << "fill the MC spectrum with the final fit values\n";
1022  for(size_t i = 0; i < gFFNames->size(); ++i){
1023  LOG_VERBATIM("CovarianceMatrixFitter")
1024  << (*gFFNames)[i]
1025  << " "
1026  << (*gFFVals)[i];
1027  }
1028 
1029  // fill the MC spectrum
1030  this->FillMCSpectrum(*gMCEventLists);
1031  }
static const double kGarbageDouble
Definition: Constants.h:22
void MinimizeCP(double const &dCPLow, double const &dCPHigh, std::vector< double > &fitVals)
ROOT::Math::Minimizer * gFFMin
static ShifterAndWeighter * Instance()
static cmf::ChiSqrCalculator * Instance()
double ChiSqr(cmf::ChiSqrCalc_t const &chiSqrCalc, std::vector< std::string > const &pars, const double *vals)
void FillMCSpectrum(cmf::EventListColl const &eventLists)
#define LOG_WARNING(category)
std::unique_ptr< cmf::EventListColl > gMCEventLists
void AddIteration(double const &chi2, std::vector< std::string > const &pars, std::vector< double > const &vals, cmf::ChiSqrCalc_t const &chiSqrCalc=cmf::kCovMat)
std::unique_ptr< cmf::Location > gLoc
static const double gPI
std::unique_ptr< std::vector< double > > gFFVals
int fNumCPRanges
number of ranges to divide up d_CP for seeding fits
#define LOG_VERBATIM(category)
double gMinChiSqr
cmf::ChiSqrCalc_t gChiSqrCalc
void SetCurrentVals(cmf::Location const &loc)
std::unique_ptr< std::vector< std::string > > gFFNames
enum BeamMode string
void cmf::CovarianceFitHelper::MinimizeCP ( double const &  dCPLow,
double const &  dCPHigh,
std::vector< double > &  fitVals 
)
private

Definition at line 1037 of file CovarianceFitHelper.cxx.

References dCP, cmf::gMinChiSqr, cmf::ParameterUtility::Instance(), cmf::kFARDET, cmf::kNEARDET, LOG_VERBATIM, LOG_WARNING, par, cmf::ParameterUtility::ParameterInfo(), string, and registry_explorer::v.

Referenced by Minimize().

1040  {
1041  // set the starting guess for dCP to be midway in the range
1042  std::string dCP("dCP");
1043  gLoc->SetParameterValue(dCP, 0.25 * dCPHigh + 0.75 * dCPLow, cmf::kFARDET);
1044  gLoc->SetParameterValue(dCP, 0.25 * dCPHigh + 0.75 * dCPLow, cmf::kNEARDET);
1045 
1046  // Set initial values for all parameterSpaceLocations and
1047  // tell the minimizer what it is fitting
1048  double lowCon;
1049  double highCon;
1050 
1051  ROOT::Fit::ParameterSettings parSet;
1052 
1053  for(size_t ipar = 0; ipar < gFFNames->size(); ++ipar){
1054  auto const& par = cmf::ParameterUtility::Instance()->ParameterInfo((*gFFNames)[ipar]);
1055 
1056  if(par.IsConstrained()){
1057  lowCon = par.LowerBound();
1058  highCon = par.UpperBound();
1059 
1060  if((*gFFNames)[ipar] == "dCP"){
1061  lowCon = dCPLow;
1062  highCon = dCPHigh;
1063  }
1064 
1065  gFFMin->SetLimitedVariable(ipar,
1066  (*gFFNames)[ipar],
1067  (*gFFVals)[ipar],
1068  0.01,
1069  lowCon,
1070  highCon);
1071  }
1072  else{
1073  gFFMin->SetVariable(ipar,
1074  (*gFFNames)[ipar],
1075  (*gFFVals)[ipar],
1076  0.01);
1077  } // end if unconstrained parameter
1078 
1079  gFFMin->GetVariableSettings(ipar, parSet);
1080  LOG_VERBATIM("CovarianceMatrixFit")
1081  << parSet.Name()
1082  << " : "
1083  << parSet.Value()
1084  << " lower limit ("
1085  << parSet.HasLowerLimit()
1086  << ") "
1087  << parSet.LowerLimit()
1088  << " upper limit ("
1089  << parSet.HasUpperLimit()
1090  << ") "
1091  << parSet.UpperLimit();
1092 
1093  } // end loop over parameters
1094 
1095  // minimize
1096  if( !gFFMin->Minimize() ){
1097  LOG_WARNING("CovarianceMatrixFit")
1098  << "Fit failed to converge d_CP with the "
1099  << gFFMin->Options().MinimizerAlgorithm()
1100  << " algorithm";
1101  }
1102  else{
1103  // the fit converged so compare to the previous minimization and
1104  // keep the lowest chi^2 result
1105  if(gFFMin->MinValue() < gMinChiSqr){
1106  gMinChiSqr = gFFMin->MinValue();
1107  gMinStatus = gFFMin->Status();
1108  for(size_t v = 0; v < gFFVals->size(); ++v) (*gFFVals)[v] = gFFMin->X()[v];
1109  }
1110  } // end if the fit converged
1111 
1112  }
cmf::Parameter const & ParameterInfo(std::string const &parName)
Int_t par
Definition: SimpleIterate.C:24
static ParameterUtility * Instance()
double dCP
ROOT::Math::Minimizer * gFFMin
#define LOG_WARNING(category)
std::unique_ptr< cmf::Location > gLoc
std::unique_ptr< std::vector< double > > gFFVals
#define LOG_VERBATIM(category)
double gMinChiSqr
std::unique_ptr< std::vector< std::string > > gFFNames
enum BeamMode string
std::unique_ptr< cmf::PointResult > cmf::CovarianceFitHelper::Result ( bool  loadFromEventList = true)

Definition at line 1124 of file CovarianceFitHelper.cxx.

References BestFitSpectrum(), fIterations, cmf::gMinChiSqr, cmf::gMinStatus, and cmf::StringToOscParmType().

Referenced by cmf::CovarianceMatrixFitter::writeResults(), and cmf::FitFeldmanCousinsPoint::writeResults().

1125  {
1126  // get the fit parameters, we want the last iteration in the vector
1127  // because that is where the fit ended
1129  for(size_t p = 0; p < fIterations.begin()->second.back().parNames.size(); ++p){
1130 
1131  parMap.emplace(cmf::StringToOscParmType(fIterations.begin()->second.back().parNames[p]),
1132  fIterations.begin()->second.back().parVals[p]);
1133  }
1134 
1135  return std::make_unique<cmf::PointResult>(gMinChiSqr,
1136  gMinStatus,
1137  parMap,
1138  this->BestFitSpectrum(loadFromEventLists));
1139  }
cmf::Spectrum const & BestFitSpectrum(bool loadFromEventList=true)
const char * p
Definition: xmltok.h:285
std::map< cmf::ChiSqrCalc_t, std::vector< cmf::CovFitIteration > > fIterations
keep track of the minimizer iterations
std::map< cmf::OscParm_t, float > OscillationParameterMap
Definition: Constants.h:612
double gMinChiSqr
static cmf::OscParm_t StringToOscParmType(std::string const &str)

Member Data Documentation

std::vector<float> cmf::CovarianceFitHelper::fDataCount
private

number of data events in each bin

Definition at line 121 of file CovarianceFitHelper.h.

Referenced by CovarianceFitHelper(), FillDataSpectrum(), and MakeResultPlots().

std::map<int, double> cmf::CovarianceFitHelper::fDataPOT
private

Data POT values for each detector/beam/selection.

Definition at line 118 of file CovarianceFitHelper.h.

Referenced by InitializeEventLists(), and Make1DSpectra().

std::vector<float> cmf::CovarianceFitHelper::fDataSpectrum
private
std::map<cmf::ChiSqrCalc_t, std::vector<cmf::CovFitIteration> > cmf::CovarianceFitHelper::fIterations
private

keep track of the minimizer iterations

Definition at line 123 of file CovarianceFitHelper.h.

Referenced by MakeIterationGraphs(), and Result().

std::vector<float> cmf::CovarianceFitHelper::fMCCount
private

number of MC events in each bin

Definition at line 122 of file CovarianceFitHelper.h.

Referenced by CovarianceFitHelper(), FillMCSpectrum(), and MakeResultPlots().

std::vector<float> cmf::CovarianceFitHelper::fMCSpectrum
private
fhicl::ParameterSet cmf::CovarianceFitHelper::fMinimizerPars
private

alternative algorithm to try if minimization fails

Definition at line 116 of file CovarianceFitHelper.h.

Referenced by InitializeMinimizer().

int cmf::CovarianceFitHelper::fNumCPRanges
private

number of ranges to divide up d_CP for seeding fits

Definition at line 115 of file CovarianceFitHelper.h.

Referenced by InitializeFitHelper(), and Minimize().

art::ServiceHandle<art::TFileService> cmf::CovarianceFitHelper::fTFS
private

TFileService.

Definition at line 117 of file CovarianceFitHelper.h.

Referenced by MakeResultPlots().


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