laura is hosted by Hepforge, IPPP Durham
Laura++  v3r5
A maximum likelihood fitting package for performing Dalitz-plot analysis.
LauAbsFitModel.cc
Go to the documentation of this file.
1 
2 /*
3 Copyright 2004 University of Warwick
4 
5 Licensed under the Apache License, Version 2.0 (the "License");
6 you may not use this file except in compliance with the License.
7 You may obtain a copy of the License at
8 
9  http://www.apache.org/licenses/LICENSE-2.0
10 
11 Unless required by applicable law or agreed to in writing, software
12 distributed under the License is distributed on an "AS IS" BASIS,
13 WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 See the License for the specific language governing permissions and
15 limitations under the License.
16 */
17 
18 /*
19 Laura++ package authors:
20 John Back
21 Paul Harrison
22 Thomas Latham
23 */
24 
29 #include <iostream>
30 #include <limits>
31 #include <vector>
32 
33 #include "TMessage.h"
34 #include "TMonitor.h"
35 #include "TServerSocket.h"
36 #include "TSocket.h"
37 #include "TSystem.h"
38 #include "TVirtualFitter.h"
39 
40 #include "LauAbsFitModel.hh"
41 #include "LauAbsFitter.hh"
42 #include "LauAbsPdf.hh"
43 #include "LauComplex.hh"
44 #include "LauFitter.hh"
45 #include "LauFitDataTree.hh"
46 #include "LauGenNtuple.hh"
47 #include "LauParameter.hh"
48 #include "LauParamFixed.hh"
49 #include "LauPrint.hh"
50 #include "LauSPlot.hh"
51 
53 
54 
56  fixParams_(kFALSE),
57  compareFitData_(kFALSE),
58  savePDF_(kFALSE),
59  writeLatexTable_(kFALSE),
60  writeSPlotData_(kFALSE),
61  storeDPEff_(kFALSE),
62  randomFit_(kFALSE),
63  emlFit_(kFALSE),
64  poissonSmear_(kFALSE),
65  enableEmbedding_(kFALSE),
66  usingDP_(kTRUE),
67  pdfsDependOnDP_(kFALSE),
68  inputFitData_(0),
69  genNtuple_(0),
70  sPlotNtuple_(0),
71  nullString_(""),
72  doSFit_(kFALSE),
73  sWeightBranchName_(""),
74  sWeightScaleFactor_(1.0),
75  outputTableName_(""),
76  fitToyMCFileName_("fitToyMC.root"),
77  fitToyMCTableName_("fitToyMCTable"),
78  fitToyMCScale_(10),
79  fitToyMCPoissonSmear_(kFALSE),
80  sPlotFileName_(""),
81  sPlotTreeName_(""),
82  sPlotVerbosity_("")
83 {
84 }
85 
86 
88 {
89  delete inputFitData_; inputFitData_ = 0;
90  delete genNtuple_; genNtuple_ = 0;
91  delete sPlotNtuple_; sPlotNtuple_ = 0;
92 
93  // Remove the components created to apply constraints to fit parameters
94  for (std::vector<LauAbsRValue*>::iterator iter = conVars_.begin(); iter != conVars_.end(); ++iter){
95  if ( !(*iter)->isLValue() ){
96  delete (*iter);
97  (*iter) = 0;
98  }
99  }
100 }
101 
102 void LauAbsFitModel::run(const TString& applicationCode, const TString& dataFileName, const TString& dataTreeName,
103  const TString& histFileName, const TString& tableFileName)
104 {
105  // Chose whether you want to generate or fit events in the Dalitz plot.
106  // To generate events choose applicationCode = "gen", to fit events choose
107  // applicationCode = "fit".
108 
109  TString runCode(applicationCode);
110  runCode.ToLower();
111 
112  TString histFileNameCopy(histFileName);
113  TString tableFileNameCopy(tableFileName);
114  TString dataFileNameCopy(dataFileName);
115  TString dataTreeNameCopy(dataTreeName);
116 
117  // Initialise the fit par vectors. Each class that inherits from this one
118  // must implement this sensibly for all vectors specified in clearFitParVectors,
119  // i.e. specify parameter names, initial, min, max and fixed values
120  this->initialise();
121 
122  // Add variables to Gaussian constrain to a list
123  this->addConParameters();
124 
125  if (dataFileNameCopy == "") {dataFileNameCopy = "data.root";}
126  if (dataTreeNameCopy == "") {dataTreeNameCopy = "genResults";}
127 
128  if (runCode.Contains("gen")) {
129 
130  if (histFileNameCopy == "") {histFileNameCopy = "parInfo.root";}
131  if (tableFileNameCopy == "") {tableFileNameCopy = "genResults";}
132 
133  this->setGenValues();
134  this->generate(dataFileNameCopy, dataTreeNameCopy, histFileNameCopy, tableFileNameCopy);
135 
136  } else if (runCode.Contains("fit")) {
137 
138  if (histFileNameCopy == "") {histFileNameCopy = "parInfo.root";}
139  if (tableFileNameCopy == "") {tableFileNameCopy = "fitResults";}
140 
141  this->fit(dataFileNameCopy, dataTreeNameCopy, histFileNameCopy, tableFileNameCopy);
142 
143  } else if (runCode.Contains("plot")) {
144 
145  this->savePDFPlots("plot");
146 
147  } else if (runCode.Contains("weight")) {
148 
149  this->weightEvents(dataFileNameCopy, dataTreeNameCopy);
150  }
151 }
152 
153 void LauAbsFitModel::doSFit( const TString& sWeightBranchName, Double_t scaleFactor )
154 {
155  if ( sWeightBranchName == "" ) {
156  std::cerr << "WARNING in LauAbsFitModel::doSFit : sWeight branch name is empty string, not setting-up sFit." << std::endl;
157  return;
158  }
159 
160  doSFit_ = kTRUE;
161  sWeightBranchName_ = sWeightBranchName;
162  sWeightScaleFactor_ = scaleFactor;
163 }
164 
165 void LauAbsFitModel::setBkgndClassNames( const std::vector<TString>& names )
166 {
167  if ( !bkgndClassNames_.empty() ) {
168  std::cerr << "WARNING in LauAbsFitModel::setBkgndClassNames : Names already stored, not changing them." << std::endl;
169  return;
170  }
171 
172  UInt_t nBkgnds = names.size();
173  for ( UInt_t i(0); i < nBkgnds; ++i ) {
174  bkgndClassNames_.insert( std::make_pair( i, names[i] ) );
175  }
176 
177  this->setupBkgndVectors();
178 }
179 
180 Bool_t LauAbsFitModel::validBkgndClass( const TString& className ) const
181 {
182  if ( bkgndClassNames_.empty() ) {
183  return kFALSE;
184  }
185 
186  Bool_t found(kFALSE);
187  for ( LauBkgndClassMap::const_iterator iter = bkgndClassNames_.begin(); iter != bkgndClassNames_.end(); ++iter ) {
188  if ( iter->second == className ) {
189  found = kTRUE;
190  break;
191  }
192  }
193 
194  return found;
195 }
196 
197 UInt_t LauAbsFitModel::bkgndClassID( const TString& className ) const
198 {
199  if ( ! this->validBkgndClass( className ) ) {
200  std::cerr << "ERROR in LauAbsFitModel::bkgndClassID : Request for ID for invalid background class \"" << className << "\"." << std::endl;
201  return (bkgndClassNames_.size() + 1);
202  }
203 
204  UInt_t bgID(0);
205  for ( LauBkgndClassMap::const_iterator iter = bkgndClassNames_.begin(); iter != bkgndClassNames_.end(); ++iter ) {
206  if ( iter->second == className ) {
207  bgID = iter->first;
208  break;
209  }
210  }
211 
212  return bgID;
213 }
214 
215 const TString& LauAbsFitModel::bkgndClassName( UInt_t classID ) const
216 {
217  LauBkgndClassMap::const_iterator iter = bkgndClassNames_.find( classID );
218 
219  if ( iter == bkgndClassNames_.end() ) {
220  std::cerr << "ERROR in LauAbsFitModel::bkgndClassName : Request for name of invalid background class ID " << classID << "." << std::endl;
221  return nullString_;
222  }
223 
224  return iter->second;
225 }
226 
228 {
229  std::cout << "INFO in LauAbsFitModel::clearFitParVectors : Clearing fit variable vectors" << std::endl;
230 
231  // Remove the components created to apply constraints to fit parameters
232  for (std::vector<LauAbsRValue*>::iterator iter = conVars_.begin(); iter != conVars_.end(); ++iter){
233  if ( !(*iter)->isLValue() ){
234  delete (*iter);
235  (*iter) = 0;
236  }
237  }
238  conVars_.clear();
239  fitVars_.clear();
240 }
241 
243 {
244  std::cout << "INFO in LauAbsFitModel::clearExtraVarVectors : Clearing extra ntuple variable vectors" << std::endl;
245  extraVars_.clear();
246 }
247 
249 {
250  // makes sure each parameter holds its genValue as its current value
251  for (LauParameterPList::iterator iter = fitVars_.begin(); iter != fitVars_.end(); ++iter) {
252  (*iter)->value((*iter)->genValue());
253  }
254  this->propagateParUpdates();
255 }
256 
257 void LauAbsFitModel::writeSPlotData(const TString& fileName, const TString& treeName, Bool_t storeDPEfficiency, const TString& verbosity)
258 {
259  if (this->writeSPlotData()) {
260  std::cerr << "ERROR in LauAbsFitModel::writeSPlotData : Already have an sPlot ntuple setup, not doing it again." << std::endl;
261  return;
262  }
263  writeSPlotData_ = kTRUE;
264  sPlotFileName_ = fileName;
265  sPlotTreeName_ = treeName;
266  sPlotVerbosity_ = verbosity;
267  storeDPEff_ = storeDPEfficiency;
268 }
269 
270 // TODO : histFileName isn't used here at the moment but could be used for
271 // storing the values of the parameters used in the generation.
272 // These could then be read and used for setting the "true" values
273 // in a subsequent fit.
274 void LauAbsFitModel::generate(const TString& dataFileName, const TString& dataTreeName, const TString& /*histFileName*/, const TString& tableFileNameBase)
275 {
276  // Create the ntuple for storing the results
277  std::cout << "INFO in LauAbsFitModel::generate : Creating generation ntuple." << std::endl;
278  if (genNtuple_ != 0) {delete genNtuple_; genNtuple_ = 0;}
279  genNtuple_ = new LauGenNtuple(dataFileName,dataTreeName);
280 
281  // add branches for storing the experiment number and the number of
282  // the event within the current experiment
283  this->addGenNtupleIntegerBranch("iExpt");
284  this->addGenNtupleIntegerBranch("iEvtWithinExpt");
285  this->setupGenNtupleBranches();
286 
287  // Start the cumulative timer
288  cumulTimer_.Start();
289 
290  const UInt_t firstExp = this->firstExpt();
291  const UInt_t nExp = this->nExpt();
292 
293  Bool_t genOK(kTRUE);
294  do {
295  // Loop over the number of experiments
296  for (UInt_t iExp = firstExp; iExp < (firstExp+nExp); ++iExp) {
297 
298  // Start the timer to see how long each experiment takes to generate
299  timer_.Start();
300 
301  // Store the experiment number in the ntuple
302  this->setGenNtupleIntegerBranchValue("iExpt",iExp);
303 
304  // Do the generation for this experiment
305  std::cout << "INFO in LauAbsFitModel::generate : Generating experiment number " << iExp << std::endl;
306  genOK = this->genExpt();
307 
308  // Stop the timer and see how long the program took so far
309  timer_.Stop();
310  timer_.Print();
311 
312  if (!genOK) {
313  // delete and recreate an empty tree
315 
316  // then break out of the experiment loop
317  std::cerr << "WARNING in LauAbsFitModel::generate : Problem in toy MC generation. Starting again with updated parameters..." << std::endl;
318  break;
319  }
320 
321  if (this->writeLatexTable()) {
322  TString tableFileName(tableFileNameBase);
323  tableFileName += "_";
324  tableFileName += iExp;
325  tableFileName += ".tex";
326  this->writeOutTable(tableFileName);
327  }
328 
329  } // Loop over number of experiments
330  } while (!genOK);
331 
332  // Print out total timing info.
333  cumulTimer_.Stop();
334  std::cout << "INFO in LauAbsFitModel::generate : Finished generating all experiments." << std::endl;
335  std::cout << "INFO in LauAbsFitModel::generate : Cumulative timing:" << std::endl;
336  cumulTimer_.Print();
337 
338  // Build the event index
339  std::cout << "INFO in LauAbsFitModel::generate : Building experiment:event index." << std::endl;
340  // TODO - can test this return value?
341  //Int_t nIndexEntries =
342  genNtuple_->buildIndex("iExpt","iEvtWithinExpt");
343 
344  // Write out toy MC ntuple
345  std::cout << "INFO in LauAbsFitModel::generate : Writing data to file " << dataFileName << "." << std::endl;
347 }
348 
350 {
352 }
353 
355 {
357 }
358 
360 {
361  genNtuple_->setIntegerBranchValue(name,value);
362 }
363 
365 {
366  genNtuple_->setDoubleBranchValue(name,value);
367 }
368 
370 {
371  return genNtuple_->getIntegerBranchValue(name);
372 }
373 
375 {
376  return genNtuple_->getDoubleBranchValue(name);
377 }
378 
380 {
382 }
383 
385 {
387 }
388 
390 {
392 }
393 
395 {
396  sPlotNtuple_->setIntegerBranchValue(name,value);
397 }
398 
400 {
401  sPlotNtuple_->setDoubleBranchValue(name,value);
402 }
403 
405 {
407 }
408 
409 void LauAbsFitModel::fit(const TString& dataFileName, const TString& dataTreeName, const TString& histFileName, const TString& tableFileNameBase)
410 {
411  // Routine to perform the total fit.
412 
413  const UInt_t firstExp = this->firstExpt();
414  const UInt_t nExp = this->nExpt();
415 
416  std::cout << "INFO in LauAbsFitModel::fit : First experiment = " << firstExp << std::endl;
417  std::cout << "INFO in LauAbsFitModel::fit : Number of experiments = " << nExp << std::endl;
418 
419  // Start the cumulative timer
420  cumulTimer_.Start();
421 
422  this->resetFitCounters();
423 
424  // Create and setup the fit results ntuple
425  this->setupResultsOutputs( histFileName, tableFileNameBase );
426 
427  // Create and setup the sPlot ntuple
428  if (this->writeSPlotData()) {
429  std::cout << "INFO in LauAbsFitModel::fit : Creating sPlot ntuple." << std::endl;
430  if (sPlotNtuple_ != 0) {delete sPlotNtuple_; sPlotNtuple_ = 0;}
432  this->setupSPlotNtupleBranches();
433  }
434 
435  // This reads in the given dataFile and creates an input
436  // fit data tree that stores them for all events and experiments.
437  Bool_t dataOK = this->verifyFitData(dataFileName,dataTreeName);
438  if (!dataOK) {
439  std::cerr << "ERROR in LauAbsFitModel::fit : Problem caching the fit data." << std::endl;
440  gSystem->Exit(EXIT_FAILURE);
441  }
442 
443  // Loop over the number of experiments
444  for (UInt_t iExp = firstExp; iExp < (firstExp+nExp); ++iExp) {
445 
446  // Start the timer to see how long each fit takes
447  timer_.Start();
448 
449  this->setCurrentExperiment( iExp );
450 
451  UInt_t nEvents = this->readExperimentData();
452  if (nEvents < 1) {
453  std::cerr << "WARNING in LauAbsFitModel::fit : Zero events in experiment " << iExp << ", skipping..." << std::endl;
454  timer_.Stop();
455  continue;
456  }
457 
458  // Now the sub-classes must implement whatever they need to do
459  // to cache any more input fit data they need in order to
460  // calculate the likelihoods during the fit.
461  // They need to use the inputFitData_ tree as input. For example,
462  // inputFitData_ contains m13Sq and m23Sq. The appropriate fit model
463  // then caches the resonance dynamics for the signal model, as
464  // well as the background likelihood values in the Dalitz plot
465  this->cacheInputFitVars();
466  if ( this->doSFit() ) {
467  this->cacheInputSWeights();
468  }
469 
470  // Do the fit for this experiment
471  this->fitExpt();
472 
473  // Write the results into the ntuple
475 
476  // Stop the timer and see how long the program took so far
477  timer_.Stop();
478  timer_.Print();
479 
480  // Store the per-event likelihood values
481  if ( this->writeSPlotData() ) {
482  this->storePerEvtLlhds();
483  }
484 
485  // Create a toy MC sample using the fitted parameters so that
486  // the user can compare the fit to the data.
487  if (compareFitData_ == kTRUE && this->statusCode() == 3) {
489  }
490 
491  } // Loop over number of experiments
492 
493  // Print out total timing info.
494  cumulTimer_.Stop();
495  std::cout << "INFO in LauAbsFitModel::fit : Cumulative timing:" << std::endl;
496  cumulTimer_.Print();
497 
498  // Print out stats on OK fits.
499  const UInt_t nOKFits = this->numberOKFits();
500  const UInt_t nBadFits = this->numberBadFits();
501  std::cout << "INFO in LauAbsFitModel::fit : Number of OK Fits = " << nOKFits << std::endl;
502  std::cout << "INFO in LauAbsFitModel::fit : Number of Failed Fits = " << nBadFits << std::endl;
503  Double_t fitEff(0.0);
504  if (nExp != 0) {fitEff = nOKFits/(1.0*nExp);}
505  std::cout << "INFO in LauAbsFitModel::fit : Fit efficiency = " << fitEff*100.0 << "%." << std::endl;
506 
507  // Write out any fit results (ntuples etc...).
508  this->writeOutAllFitResults();
509  if ( this->writeSPlotData() ) {
510  this->calculateSPlotData();
511  }
512 }
513 
514 void LauAbsFitModel::setupResultsOutputs( const TString& histFileName, const TString& tableFileName )
515 {
516  this->LauSimFitSlave::setupResultsOutputs( histFileName, tableFileName );
517 
518  outputTableName_ = tableFileName;
519 }
520 
521 Bool_t LauAbsFitModel::verifyFitData(const TString& dataFileName, const TString& dataTreeName)
522 {
523  // From the input data stream, store the variables into the
524  // internal tree inputFitData_ that can be used by the sub-classes
525  // in calculating their likelihood functions for the fit
526  delete inputFitData_;
527  inputFitData_ = new LauFitDataTree(dataFileName,dataTreeName);
528  Bool_t dataOK = inputFitData_->findBranches();
529 
530  if (!dataOK) {
531  delete inputFitData_; inputFitData_ = 0;
532  }
533 
534  return dataOK;
535 }
536 
538 {
539  Bool_t hasBranch = inputFitData_->haveBranch( sWeightBranchName_ );
540  if ( ! hasBranch ) {
541  std::cerr << "ERROR in LauAbsFitModel::cacheInputSWeights : Input data does not contain variable \"" << sWeightBranchName_ << "\".\n";
542  std::cerr << " : Turning off sFit!" << std::endl;
543  doSFit_ = kFALSE;
544  return;
545  }
546 
547  UInt_t nEvents = this->eventsPerExpt();
548  sWeights_.clear();
549  sWeights_.reserve( nEvents );
550 
551  for (UInt_t iEvt = 0; iEvt < nEvents; ++iEvt) {
552 
553  const LauFitData& dataValues = inputFitData_->getData(iEvt);
554 
555  LauFitData::const_iterator iter = dataValues.find( sWeightBranchName_ );
556 
557  sWeights_.push_back( iter->second * sWeightScaleFactor_ );
558  }
559 }
560 
562 {
563  // Routine to perform the actual fit for the given experiment
564 
565  // Update initial fit parameters if required (e.g. if using random numbers).
566  this->checkInitFitParams();
567 
568  // Initialise the fitter
572 
573  this->startNewFit( LauFitter::fitter()->nParameters(), LauFitter::fitter()->nFreeParameters() );
574 
575  // Now ready for minimisation step
576  std::cout << "\nINFO in LauAbsFitModel::fitExpt : Start minimisation...\n";
578 
579  // If we're doing a two stage fit we can now release (i.e. float)
580  // the 2nd stage parameters and re-fit
581  if (this->twoStageFit()) {
582  if ( fitResult.status != 3 ) {
583  std::cerr << "WARNING in LauAbsFitModel:fitExpt : Not running second stage fit since first stage failed." << std::endl;
585  } else {
587  this->startNewFit( LauFitter::fitter()->nParameters(), LauFitter::fitter()->nFreeParameters() );
588  fitResult = LauFitter::fitter()->minimise();
589  }
590  }
591 
592  const TMatrixD& covMat = LauFitter::fitter()->covarianceMatrix();
593  this->storeFitStatus( fitResult, covMat );
594 
595  // Store the final fit results and errors into protected internal vectors that
596  // all sub-classes can use within their own finalFitResults implementation
597  // used below (e.g. putting them into an ntuple in a root file)
599 }
600 
602 {
603  if (sPlotNtuple_ != 0) {
606  LauSPlot splot(sPlotNtuple_->fileName(), sPlotNtuple_->treeName(), this->firstExpt(), this->nExpt(),
607  this->variableNames(), this->freeSpeciesNames(), this->fixdSpeciesNames(), this->twodimPDFs(),
608  this->splitSignal(), this->scfDPSmear());
610  splot.writeOutResults();
611  }
612 }
613 
614 void LauAbsFitModel::compareFitData(UInt_t toyMCScale, const TString& mcFileName, const TString& tableFileName, Bool_t poissonSmearing)
615 {
616  compareFitData_ = kTRUE;
617  fitToyMCScale_ = toyMCScale;
618  fitToyMCFileName_ = mcFileName;
619  fitToyMCTableName_ = tableFileName;
620  fitToyMCPoissonSmear_ = poissonSmearing;
621 }
622 
623 void LauAbsFitModel::createFitToyMC(const TString& mcFileName, const TString& tableFileName)
624 {
625  // Create a toy MC sample so that the user can compare the fitted
626  // result with the data.
627  // Generate more toy MC to reduce statistical fluctuations:
628  // - use the rescaling value fitToyMCScale_
629 
630  // Store the info on the number of experiments, first expt and current expt
631  const UInt_t oldNExpt(this->nExpt());
632  const UInt_t oldFirstExpt(this->firstExpt());
633  const UInt_t oldIExpt(this->iExpt());
634 
635  // Turn off Poisson smearing if required
636  const Bool_t poissonSmearing(this->doPoissonSmearing());
638 
639  // Turn off embedding, since we need toy MC, not reco'ed events
640  const Bool_t enableEmbeddingOrig(this->enableEmbedding());
641  this->enableEmbedding(kFALSE);
642 
643  // Need to make sure that the generation of the DP co-ordinates is
644  // switched on if any of our PDFs depend on it
645  const Bool_t origUseDP = this->useDP();
646  if ( !origUseDP && this->pdfsDependOnDP() ) {
647  this->useDP( kTRUE );
648  this->initialiseDPModels();
649  }
650 
651  // Construct a unique filename for this experiment
652  TString exptString("_expt");
653  exptString += oldIExpt;
654  TString fileName( mcFileName );
655  fileName.Insert( fileName.Last('.'), exptString );
656 
657  // Generate the toy MC
658  std::cout << "INFO in LauAbsFitModel::createFitToyMC : Generating toy MC in " << fileName << " to compare fit with data..." << std::endl;
659  std::cout << " : Number of experiments to generate = " << fitToyMCScale_ << "." << std::endl;
660  std::cout << " : This is to allow the toy MC to be made with reduced statistical fluctuations." << std::endl;
661 
662  // Set the genValue of each parameter to its current (fitted) value
663  // but first store the original genValues for restoring later
664  std::vector<Double_t> origGenValues; origGenValues.reserve(this->nTotParams());
665  Bool_t blind(kFALSE);
666  for (LauParameterPList::iterator iter = fitVars_.begin(); iter != fitVars_.end(); ++iter) {
667  origGenValues.push_back((*iter)->genValue());
668  (*iter)->genValue((*iter)->unblindValue());
669  if ( (*iter)->blind() ) {
670  blind = kTRUE;
671  }
672  }
673  if ( blind ) {
674  std::cerr << "WARNING in LauAbsFitModel::createFitToyMC : One or more parameters are blind but the toy will be created using the unblind values - use with caution!!" << std::endl;
675  }
676 
677  // If we're asked to generate more than 100 experiments then split it
678  // up into multiple files since otherwise can run into memory issues
679  // when building the index
680  // TODO - this obviously depends on the number of events per experiment as well, so should do this properly
681 
682  UInt_t totalExpts = fitToyMCScale_;
683  if ( totalExpts > 100 ) {
684  UInt_t nFiles = totalExpts/100;
685  if ( totalExpts%100 ) {
686  nFiles += 1;
687  }
688  TString fileNameBase {fileName};
689  for ( UInt_t iFile(0); iFile < nFiles; ++iFile ) {
690 
691  UInt_t firstExp( iFile*100 );
692 
693  // Set number of experiments and first experiment to generate
694  UInt_t nExp = ((firstExp + 100)>totalExpts) ? totalExpts-firstExp : 100;
695  this->setNExpts(nExp, firstExp);
696 
697  // Create a unique filename and generate the events
698  fileName = fileNameBase;
699  TString extraname = "_file";
700  extraname += iFile;
701  fileName.Insert( fileName.Last('.'), extraname );
702  this->generate(fileName, "genResults", "dummy.root", tableFileName);
703  }
704  } else {
705  // Set number of experiments to new value
706  this->setNExpts(fitToyMCScale_, 0);
707  // Generate the toy
708  this->generate(fileName, "genResults", "dummy.root", tableFileName);
709  }
710 
711  // Reset number of experiments to original value
712  this->setNExpts(oldNExpt, oldFirstExpt);
713  this->setCurrentExperiment(oldIExpt);
714 
715  // Restore the Poisson smearing to its former value
716  this->doPoissonSmearing(poissonSmearing);
717 
718  // Restore the embedding status
719  this->enableEmbedding(enableEmbeddingOrig);
720 
721  // Restore "useDP" to its former status
722  this->useDP( origUseDP );
723 
724  // Restore the original genValue to each parameter
725  for (UInt_t i(0); i<this->nTotParams(); ++i) {
726  fitVars_[i]->genValue(origGenValues[i]);
727  }
728 
729  std::cout << "INFO in LauAbsFitModel::createFitToyMC : Finished in createFitToyMC." << std::endl;
730 }
731 
733 {
734  // Calculate the total negative log-likelihood over all events.
735  // This function assumes that the fit parameters and data tree have
736  // already been set-up correctly.
737 
738  // Loop over the data points to calculate the log likelihood
739  Double_t logLike = this->getLogLikelihood( 0, this->eventsPerExpt() );
740 
741  // Include the Poisson term in the extended likelihood if required
742  if (this->doEMLFit()) {
743  logLike -= this->getEventSum();
744  }
745 
746  // Calculate any penalty terms from Gaussian constrained variables
747  if ( ! conVars_.empty() ){
748  logLike -= this->getLogLikelihoodPenalty();
749  }
750 
751  Double_t totNegLogLike = -logLike;
752  return totNegLogLike;
753 }
754 
756 {
757  Double_t penalty(0.0);
758 
759  for ( LauAbsRValuePList::const_iterator iter = conVars_.begin(); iter != conVars_.end(); ++iter ) {
760  Double_t val = (*iter)->unblindValue();
761  Double_t mean = (*iter)->constraintMean();
762  Double_t width = (*iter)->constraintWidth();
763 
764  Double_t term = ( val - mean )*( val - mean );
765  penalty += term/( 2*width*width );
766  }
767 
768  return penalty;
769 }
770 
771 Double_t LauAbsFitModel::getLogLikelihood( UInt_t iStart, UInt_t iEnd )
772 {
773  // Calculate the total negative log-likelihood over all events.
774  // This function assumes that the fit parameters and data tree have
775  // already been set-up correctly.
776 
777  // Loop over the data points to calculate the log likelihood
778  Double_t logLike(0.0);
779  const Double_t worstLL = this->worstLogLike();
780 
781  // Loop over the number of events in this experiment
782  Bool_t ok(kTRUE);
783  for (UInt_t iEvt = iStart; iEvt < iEnd; ++iEvt) {
784 
785  Double_t likelihood = this->getTotEvtLikelihood(iEvt);
786 
787  if (likelihood > std::numeric_limits<Double_t>::min()) { // Is the likelihood zero?
788  Double_t evtLogLike = TMath::Log(likelihood);
789  if ( doSFit_ ) {
790  evtLogLike *= sWeights_[iEvt];
791  }
792  logLike += evtLogLike;
793  } else {
794  ok = kFALSE;
795  std::cerr << "WARNING in LauAbsFitModel::getLogLikelihood : Strange likelihood value for event " << iEvt << ": " << likelihood << "\n";
796  this->printEventInfo(iEvt);
797  this->printVarsInfo(); //Write the values of the floated variables for which the likelihood is zero
798  break;
799  }
800  }
801 
802  if (!ok) {
803  std::cerr << " : Returning worst NLL found so far to force MINUIT out of this region." << std::endl;
804  logLike = worstLL;
805  } else if (logLike < worstLL) {
806  this->worstLogLike( logLike );
807  }
808 
809  return logLike;
810 }
811 
812 void LauAbsFitModel::setParsFromMinuit(Double_t* par, Int_t npar)
813 {
814  // This function sets the internal parameters based on the values
815  // that Minuit is using when trying to minimise the total likelihood function.
816 
817  // MINOS reports different numbers of free parameters depending on the
818  // situation, so disable this check
819  if ( ! this->withinAsymErrorCalc() ) {
820  const UInt_t nFreePars = this->nFreeParams();
821  if (static_cast<UInt_t>(npar) != nFreePars) {
822  std::cerr << "ERROR in LauAbsFitModel::setParsFromMinuit : Unexpected number of free parameters: " << npar << ".\n";
823  std::cerr << " Expected: " << nFreePars << ".\n" << std::endl;
824  gSystem->Exit(EXIT_FAILURE);
825  }
826  }
827 
828  // Despite npar being the number of free parameters
829  // the par array actually contains all the parameters,
830  // free and floating...
831 
832  // Update all the floating ones with their new values
833  // Also check if we have any parameters on which the DP integrals depend
834  // and whether they have changed since the last iteration
835  Bool_t recalcNorm(kFALSE);
836  const LauParameterPSet::const_iterator resVarsEnd = resVars_.end();
837  for (UInt_t i(0); i<this->nTotParams(); ++i) {
838  if (!fitVars_[i]->fixed()) {
839  if ( resVars_.find( fitVars_[i] ) != resVarsEnd ) {
840  if ( fitVars_[i]->value() != par[i] ) {
841  recalcNorm = kTRUE;
842  }
843  }
844  fitVars_[i]->value(par[i]);
845  }
846  }
847 
848  // If so, then recalculate the normalisation
849  if (recalcNorm) {
850  this->recalculateNormalisation();
851  }
852 
853  this->propagateParUpdates();
854 }
855 
857 {
858  UInt_t nParsAdded(0);
859  for (LauPdfList::iterator pdf_iter = pdfList.begin(); pdf_iter != pdfList.end(); ++pdf_iter) {
860  LauAbsPdf* pdf = (*pdf_iter);
861  if ( pdf->isDPDependent() ) {
862  this->pdfsDependOnDP( kTRUE );
863  }
864  LauAbsRValuePList& pars = pdf->getParameters();
865  for (LauAbsRValuePList::iterator pars_iter = pars.begin(); pars_iter != pars.end(); ++pars_iter) {
866  LauParameterPList params = (*pars_iter)->getPars();
867  for (LauParameterPList::iterator params_iter = params.begin(); params_iter != params.end(); ++params_iter) {
868  if ( !(*params_iter)->clone() && ( !(*params_iter)->fixed() || ( this->twoStageFit() && (*params_iter)->secondStage() ) ) ) {
869  fitVars_.push_back(*params_iter);
870  ++nParsAdded;
871  }
872  }
873  }
874  }
875  return nParsAdded;
876 }
877 
879 {
880  for ( LauParameterPList::const_iterator iter = fitVars_.begin(); iter != fitVars_.end(); ++iter ) {
881  if ( (*iter)->gaussConstraint() ) {
882  conVars_.push_back( *iter );
883  std::cout << "INFO in LauAbsFitModel::addConParameters : Added Gaussian constraint to parameter "<< (*iter)->name() << std::endl;
884  }
885  }
886 
887  // Add penalties from the constraints to fit parameters
888  const std::vector<StoreConstraints>& storeCon = this->constraintsStore();
889  for ( std::vector<StoreConstraints>::const_iterator iter = storeCon.begin(); iter != storeCon.end(); ++iter ) {
890  const std::vector<TString>& names = (*iter).conPars_;
891  std::vector<LauParameter*> params;
892  for ( std::vector<TString>::const_iterator iternames = names.begin(); iternames != names.end(); ++iternames ) {
893  for ( LauParameterPList::const_iterator iterfit = fitVars_.begin(); iterfit != fitVars_.end(); ++iterfit ) {
894  if ( (*iternames) == (*iterfit)->name() ){
895  params.push_back(*iterfit);
896  }
897  }
898  }
899 
900  // If the parameters are not found, skip it
901  if ( params.size() != (*iter).conPars_.size() ) {
902  std::cerr << "WARNING in LauAbsFitModel::addConParameters: Could not find parameters to constrain in the formula... skipping" << std::endl;
903  continue;
904  }
905 
906  LauFormulaPar* formPar = new LauFormulaPar( (*iter).formula_, (*iter).formula_, params );
907  formPar->addGaussianConstraint( (*iter).mean_, (*iter).width_ );
908  conVars_.push_back(formPar);
909 
910  std::cout << "INFO in LauAbsFitModel::addConParameters : Added Gaussian constraint to formula\n";
911  std::cout << " : Formula: " << (*iter).formula_ << std::endl;
912  for ( std::vector<LauParameter*>::iterator iterparam = params.begin(); iterparam != params.end(); ++iterparam ) {
913  std::cout << " : Parameter: " << (*iterparam)->name() << std::endl;
914  }
915  }
916 
917 }
918 
920 {
921  for (LauPdfList::iterator pdf_iter = pdfList.begin(); pdf_iter != pdfList.end(); ++pdf_iter) {
922  (*pdf_iter)->updatePulls();
923  }
924 }
925 
926 void LauAbsFitModel::printFitParameters(const LauPdfList& pdfList, std::ostream& fout) const
927 {
928  LauPrint print;
929  for (LauPdfList::const_iterator pdf_iter = pdfList.begin(); pdf_iter != pdfList.end(); ++pdf_iter) {
930  const LauAbsRValuePList& pars = (*pdf_iter)->getParameters();
931  for (LauAbsRValuePList::const_iterator pars_iter = pars.begin(); pars_iter != pars.end(); ++pars_iter) {
932  LauParameterPList params = (*pars_iter)->getPars();
933  for (LauParameterPList::iterator params_iter = params.begin(); params_iter != params.end(); ++params_iter) {
934  if (!(*params_iter)->clone()) {
935  fout << (*params_iter)->name() << " & $";
936  print.printFormat(fout, (*params_iter)->value());
937  if ((*params_iter)->fixed() == kTRUE) {
938  fout << "$ (fixed) \\\\";
939  } else {
940  fout << " \\pm ";
941  print.printFormat(fout, (*params_iter)->error());
942  fout << "$ \\\\" << std::endl;
943  }
944  }
945  }
946  }
947  }
948 }
949 
951 {
952  for (LauPdfList::iterator pdf_iter = pdfList.begin(); pdf_iter != pdfList.end(); ++pdf_iter) {
953  (*pdf_iter)->cacheInfo(theData);
954  }
955 }
956 
957 Double_t LauAbsFitModel::prodPdfValue(LauPdfList& pdfList, UInt_t iEvt)
958 {
959  Double_t pdfVal = 1.0;
960  for (LauPdfList::iterator pdf_iter = pdfList.begin(); pdf_iter != pdfList.end(); ++pdf_iter) {
961  (*pdf_iter)->calcLikelihoodInfo(iEvt);
962  pdfVal *= (*pdf_iter)->getLikelihood();
963  }
964  return pdfVal;
965 }
966 
967 void LauAbsFitModel::printEventInfo(UInt_t iEvt) const
968 {
969  const LauFitData& data = inputFitData_->getData(iEvt);
970  std::cerr << " : Input data values for this event:" << std::endl;
971  for (LauFitData::const_iterator iter = data.begin(); iter != data.end(); ++iter) {
972  std::cerr << " " << iter->first << " = " << iter->second << std::endl;
973  }
974 }
975 
977 {
978  std::cerr << " : Current values of fit parameters:" << std::endl;
979  for (UInt_t i(0); i<this->nTotParams(); ++i) {
980  std::cerr << " " << (fitVars_[i]->name()).Data() << " = " << fitVars_[i]->value() << std::endl;
981  }
982 }
983 
985 {
986  // Update initial fit parameters if required (e.g. if using random numbers).
987  this->checkInitFitParams();
988 
989  // Store the total number of parameters and the number of free parameters
990  UInt_t nPars = fitVars_.size();
991  UInt_t nFreePars = 0;
992 
993  // Send the fit parameters
994  for ( LauParameterPList::iterator iter = fitVars_.begin(); iter != fitVars_.end(); ++iter ) {
995  if ( ! (*iter)->fixed() ) {
996  ++nFreePars;
997  }
998  array.Add( *iter );
999  }
1000 
1001  this->startNewFit( nPars, nFreePars );
1002 }
1003 
1004 void LauAbsFitModel::finaliseExperiment( const LauAbsFitter::FitStatus& fitStat, const TObjArray* parsFromMaster, const TMatrixD* covMat, TObjArray& parsToMaster )
1005 {
1006  // Copy the fit status information
1007  this->storeFitStatus( fitStat, *covMat );
1008 
1009  // Now process the parameters
1010  const UInt_t nPars = this->nTotParams();
1011  UInt_t nParsFromMaster = parsFromMaster->GetEntries();
1012  if ( nParsFromMaster != nPars ) {
1013  std::cerr << "ERROR in LauAbsFitModel::finaliseExperiment : Unexpected number of parameters received from master" << std::endl;
1014  std::cerr << " : Received " << nParsFromMaster << " when expecting " << nPars << std::endl;
1015  gSystem->Exit( EXIT_FAILURE );
1016  }
1017 
1018  for ( UInt_t iPar(0); iPar < nParsFromMaster; ++iPar ) {
1019  LauParameter* parameter = dynamic_cast<LauParameter*>( (*parsFromMaster)[iPar] );
1020  if ( ! parameter ) {
1021  std::cerr << "ERROR in LauAbsFitModel::finaliseExperiment : Error reading parameter from master" << std::endl;
1022  gSystem->Exit( EXIT_FAILURE );
1023  }
1024 
1025  if ( parameter->name() != fitVars_[iPar]->name() ) {
1026  std::cerr << "ERROR in LauAbsFitModel::finaliseExperiment : Error reading parameter from master" << std::endl;
1027  gSystem->Exit( EXIT_FAILURE );
1028  }
1029 
1030  *(fitVars_[iPar]) = *parameter;
1031  }
1032 
1033  // Write the results into the ntuple
1035 
1036  // Store the per-event likelihood values
1037  if ( this->writeSPlotData() ) {
1038  this->storePerEvtLlhds();
1039  }
1040 
1041  // Create a toy MC sample using the fitted parameters so that
1042  // the user can compare the fit to the data.
1043  if (compareFitData_ == kTRUE && fitStat.status == 3) {
1045  }
1046 
1047  // Send the finalised fit parameters
1048  for ( LauParameterPList::iterator iter = fitVars_.begin(); iter != fitVars_.end(); ++iter ) {
1049  parsToMaster.Add( *iter );
1050  }
1051 }
1052 
1054 {
1055  // retrieve the data and find out how many events have been read
1056  const UInt_t exptIndex = this->iExpt();
1057  inputFitData_->readExperimentData( exptIndex );
1058  const UInt_t nEvent = inputFitData_->nEvents();
1059  this->eventsPerExpt( nEvent );
1060  return nEvent;
1061 }
1062 
1063 void LauAbsFitModel::setParametersFromFile(const TString& fileName, const TString& treeName, const Bool_t fix)
1064 {
1065  fixParamFileName_ = fileName;
1066  fixParamTreeName_ = treeName;
1067  fixParams_ = fix;
1068 }
1069 
1070 void LauAbsFitModel::setParametersFromMap(const std::map<TString, Double_t>& parameters, const Bool_t fix)
1071 {
1072  fixParamMap_ = parameters;
1073  fixParams_ = fix;
1074 }
1075 
1076 void LauAbsFitModel::setNamedParameters(const TString& fileName, const TString& treeName, const std::set<TString>& parameters, const Bool_t fix)
1077 {
1078  fixParamFileName_ = fileName;
1079  fixParamTreeName_ = treeName;
1080  fixParamNames_ = parameters;
1081  fixParams_ = fix;
1082 }
1083 
1084 void LauAbsFitModel::setParametersFileFallback(const TString& fileName, const TString& treeName, const std::map<TString, Double_t>& parameters, const Bool_t fix)
1085 {
1086  fixParamFileName_ = fileName;
1087  fixParamTreeName_ = treeName;
1088  fixParamMap_ = parameters;
1089  fixParams_ = fix;
1090 }
Class for defining the SPlot technique.
Definition: LauSPlot.hh:69
UInt_t nTotParams() const
Access the total number of fit parameters.
virtual const FitStatus & minimise()=0
Perform the minimisation of the fit function.
Int_t statusCode() const
Access the fit status code.
TString fitToyMCTableName_
The output table name for Toy MC.
virtual void addSPlotNtupleIntegerBranch(const TString &name)
Add a branch to the sPlot tree for storing an integer.
Bool_t findBranches()
Find all of the branches in the tree.
virtual Bool_t genExpt()=0
The method that actually generates the toy MC events for the given experiment.
void createFitToyMC(const TString &mcFileName, const TString &tableFileName)
Create a toy MC sample from the fitted parameters.
File containing declaration of LauFitDataTree class.
virtual Double_t getEventSum() const =0
Returns the sum of the expected events over all hypotheses; used in the EML fit scenario.
Bool_t fixed() const
Check whether the parameter is fixed or floated.
virtual void updateParameters()=0
Update the values and errors of the parameters based on the fit minimum.
void setGenValues()
Make sure all parameters hold their genValue as the current value.
TString sPlotFileName_
The name of the sPlot file.
virtual void setupResultsOutputs(const TString &histFileName, const TString &tableFileName)
Setup saving of fit results to ntuple/LaTeX table etc.
virtual void writeOutTable(const TString &outputFile)=0
Write the latex table.
virtual UInt_t readExperimentData()
Read in the data for the current experiment.
Bool_t writeLatexTable() const
Determine whether writing out of the latex table is enabled.
Bool_t doSFit_
Option to perfom the sFit.
Bool_t writeSPlotData() const
Determine whether the sPlot data is to be written out.
TStopwatch timer_
The fit timer.
File containing declaration of LauParamFixed class.
Bool_t withinAsymErrorCalc() const
Query whether the fit is calculating the asymmetric errors.
Definition: LauFitObject.hh:88
virtual Double_t getTotNegLogLikelihood()
Calculates the total negative log-likelihood.
LauAbsRValuePList conVars_
Internal vectors of Gaussian parameters.
Bool_t storeDPEff_
Option to store DP efficiencies in the sPlot ntuple.
void cacheInfo(LauPdfList &pdfList, const LauFitDataTree &theData)
Have all PDFs in the list cache the data.
void runCalculations(const TString &option="q")
Method to calculate the sWeights and cN coeffs.
Definition: LauSPlot.cc:640
virtual void setGenNtupleIntegerBranchValue(const TString &name, Int_t value)
Set the value of an integer branch in the gen tree.
ClassImp(LauAbsCoeffSet)
virtual void setupBkgndVectors()=0
Method to set up the storage for background-related quantities called by setBkgndClassNames.
LauGenNtuple * genNtuple_
The generated ntuple.
File containing declaration of LauAbsFitModel class.
virtual Bool_t verifyFitData(const TString &dataFileName, const TString &dataTreeName)
Open the input file and verify that all required variables are present.
virtual void writeOutAllFitResults()
Write out any fit results.
const TString & name() const
The parameter name.
Class for defining combinations of fit parameter objects.
virtual void setupGenNtupleBranches()=0
Setup the generation ntuple branches.
UInt_t addFitParameters(LauPdfList &pdfList)
Add parameters of the PDFs in the list to the list of all fit parameters.
virtual void generate(const TString &dataFileName, const TString &dataTreeName, const TString &histFileName, const TString &tableFileNameBase)
Generate toy MC.
virtual void setBkgndClassNames(const std::vector< TString > &names)
Setup the background class names.
void setNExpts(UInt_t nExperiments, UInt_t firstExperiment=0)
Set the number of experiments and the first experiment.
Definition: LauFitObject.hh:95
Bool_t useDP() const
Is the Dalitz plot term in the likelihood.
Int_t buildIndex(const TString &majorName, const TString &minorName="0")
Create an index table using leaves of the tree.
Double_t worstLogLike() const
Access the worst log likelihood found so far.
void compareFitData(UInt_t toyMCScale=10, const TString &mcFileName="fitToyMC.root", const TString &tableFileName="fitToyMCTable.tex", Bool_t poissonSmearing=kTRUE)
Specify that a toy MC sample should be created for a successful fit to an experiment.
const TString & treeName() const
Ntuple tree name.
Definition: LauGenNtuple.hh:70
virtual void finaliseExperiment(const LauAbsFitter::FitStatus &fitStat, const TObjArray *parsFromMaster, const TMatrixD *covMat, TObjArray &parsToMaster)
Perform all finalisation actions.
Int_t getIntegerBranchValue(const TString &name) const
Get value of an integer branch.
virtual void cacheInputSWeights()
Cache the value of the sWeights to be used in the sFit.
virtual Bool_t twoStageFit() const =0
Determine whether the two-stage fit is enabled.
virtual void setGenNtupleDoubleBranchValue(const TString &name, Double_t value)
Set the value of a double branch in the gen tree.
void deleteAndRecreateTree()
Delete and recreate tree.
File containing declaration of LauAbsFitter class.
void storeFitStatus(const LauAbsFitter::FitStatus &status, const TMatrixD &covMatrix)
Store fit status information.
Definition: LauFitObject.cc:70
void addFriendTree(const TString &rootFileName, const TString &rootTreeName)
Add a friend tree.
virtual void storePerEvtLlhds()=0
Store the per-event likelihood values.
virtual void printEventInfo(UInt_t iEvt) const
Prints the values of all the fit variables for the specified event - useful for diagnostics.
Bool_t pdfsDependOnDP() const
Do any of the PDFs have a dependence on the DP?
virtual ~LauAbsFitModel()
Destructor.
Bool_t compareFitData_
Option to make toy from 1st successful experiment.
File containing declaration of LauPrint class.
Double_t prodPdfValue(LauPdfList &pdfList, UInt_t iEvt)
Calculate the product of the per-event likelihoods of the PDFs in the list.
void setParametersFileFallback(const TString &fileName, const TString &treeName, const std::map< TString, Double_t > &parameters, const Bool_t fix)
Set named model parameters from a given std::map, with fallback to those from a file.
Bool_t twoStageFit() const
Report whether the two-stage fit is enabled.
Definition: LauFitObject.hh:73
Double_t getDoubleBranchValue(const TString &name) const
Get value of a double branch.
std::vector< LauAbsPdf * > LauPdfList
List of Pdfs.
std::map< TString, Double_t > LauFitData
Type for holding event data.
void fillBranches()
Fill branches in the ntuple.
virtual Int_t getGenNtupleIntegerBranchValue(const TString &name) const
Get the value of an integer branch in the gen tree.
virtual LauSPlot::NumbMap fixdSpeciesNames() const =0
Returns the names and yields of species that are fixed in the fit.
LauBkgndClassMap bkgndClassNames_
The background class names.
Class to define various output print commands.
Definition: LauPrint.hh:43
virtual void setupSPlotNtupleBranches()=0
Setup the branches of the sPlot tuple.
virtual void fillGenNtupleBranches()
Fill the gen tuple branches.
TString fixParamFileName_
Imported parameters file name.
virtual void setSPlotNtupleDoubleBranchValue(const TString &name, Double_t value)
Set the value of a double branch in the sPlot tree.
UInt_t nFreeParams() const
Access the total number of fit parameters.
void printFitParameters(const LauPdfList &pdfList, std::ostream &fout) const
Print the fit parameters for all PDFs in the list.
std::map< TString, Double_t > fixParamMap_
Map from imported parameter name to value.
void addGaussianConstraint(Double_t newGaussMean, Double_t newGaussWidth)
Add a Gaussian constraint (or modify an existing one)
virtual void initialise(LauFitObject *fitObj, const std::vector< LauParameter * > &parameters)=0
Initialise the fitter, setting the information on the parameters.
void readExperimentData(UInt_t iExpt)
Read events only for the given experiment.
std::vector< LauAbsRValue * > LauAbsRValuePList
List of parameter pointers.
void resetFitCounters()
Reset the good/bad fit counters.
Definition: LauFitObject.cc:53
virtual Double_t getTotEvtLikelihood(UInt_t iEvt)=0
Calculates the likelihood for a given event.
TStopwatch cumulTimer_
The total fit timer.
File containing declaration of LauFitter class.
virtual void propagateParUpdates()=0
This function (specific to each model) calculates anything that depends on the fit parameter values...
void addConParameters()
Add parameters to the list of Gaussian constrained parameters.
void setParametersFromMap(const std::map< TString, Double_t > &parameters, const Bool_t fix)
Set model parameters from a given std::map.
virtual void recalculateNormalisation()=0
Recalculate normalisation the signal DP model(s)
LauFitDataTree * inputFitData_
The input data.
Bool_t doSFit() const
Return the flag to store the status of using an sFit or not.
void run(const TString &applicationCode, const TString &dataFileName, const TString &dataTreeName, const TString &histFileName, const TString &tableFileName="")
Start the toy generation / fitting.
std::vector< LauParameter * > LauParameterPList
List of parameter pointers.
void setNamedParameters(const TString &fileName, const TString &treeName, const std::set< TString > &parameters, const Bool_t fix)
Set named model parameters from a file.
const TString & treeName() const
Retrieve the tree name.
LauParameterList extraVars_
Extra variables that aren&#39;t in the fit but are stored in the ntuple.
Bool_t blind() const
The blinding state.
void addIntegerBranch(const TString &name)
Add integer branch to tree.
Definition: LauGenNtuple.cc:87
const TString & fileName() const
Ntuple file name.
Definition: LauGenNtuple.hh:64
UInt_t eventsPerExpt() const
Obtain the total number of events in the current experiment.
void setCurrentExperiment(const UInt_t curExpt)
Set the ID of the current experiment.
Abstract interface to the fitting and toy MC model.
const std::vector< StoreConstraints > & constraintsStore() const
Const access to the constraints store.
Int_t status
The status code of the fit.
Definition: LauAbsFitter.hh:57
virtual void addGenNtupleIntegerBranch(const TString &name)
Add a branch to the gen tree for storing an integer.
Bool_t fixParams_
Whether to fix the loaded parameters (kTRUE) or leave them floating (kFALSE)
virtual void prepareInitialParArray(TObjArray &array)
Package the initial fit parameters for transmission to the master.
File containing declaration of LauParameter class.
Bool_t doPoissonSmearing() const
Determine whether Poisson smearing is enabled for the toy MC generation.
Struct to store fit status information.
Definition: LauAbsFitter.hh:55
virtual Double_t getGenNtupleDoubleBranchValue(const TString &name) const
Get the value of a double branch in the gen tree.
virtual void weightEvents(const TString &dataFileName, const TString &dataTreeName)=0
Weighting - allows e.g. MC events to be weighted by the DP model.
Bool_t doEMLFit() const
Determine whether an extended maximum likelihood fit it being performed.
File containing declaration of LauComplex class.
const TString & fileName() const
Retrieve the file name.
virtual Bool_t isDPDependent() const
Specifies whether or not the PDF is DP dependent.
Definition: LauAbsPdf.hh:125
virtual LauSPlot::NumbMap freeSpeciesNames() const =0
Returns the names and yields of species that are free in the fit.
Class to store the results from the toy MC generation into an ntuple.
Definition: LauGenNtuple.hh:46
File containing declaration of LauSPlot class.
virtual Bool_t splitSignal() const =0
Check if the signal is split into well-reconstructed and mis-reconstructed types. ...
virtual void releaseSecondStageParameters()=0
Release parameters marked as &quot;second stage&quot;.
TString sPlotTreeName_
The name of the sPlot tree.
virtual void initialiseDPModels()=0
Initialise the DP models.
void setParametersFromFile(const TString &fileName, const TString &treeName, const Bool_t fix)
Set model parameters from a file.
LauGenNtuple * sPlotNtuple_
The sPlot ntuple.
UInt_t nExpt() const
Obtain the number of experiments.
TString sPlotVerbosity_
Control the verbosity of the sFit.
const TString nullString_
An empty string.
void clearFitParVectors()
Clear the vectors containing fit parameters.
Class for defining the fit parameter objects.
Definition: LauParameter.hh:49
virtual void addSPlotNtupleDoubleBranch(const TString &name)
Add a branch to the sPlot tree for storing a double.
void clearExtraVarVectors()
Clear the vectors containing extra ntuple variables.
virtual void cacheInputFitVars()=0
Cache the input data values to calculate the likelihood during the fit.
void fitExpt()
Routine to perform the actual fit for a given experiment.
virtual const std::vector< LauAbsRValue * > & getParameters() const
Retrieve the parameters of the PDF, e.g. so that they can be loaded into a fit.
Definition: LauAbsPdf.hh:253
Double_t sWeightScaleFactor_
The sWeight scaling factor.
File containing declaration of LauAbsPdf class.
UInt_t bkgndClassID(const TString &className) const
The number assigned to a background class.
const LauFitData & getData(UInt_t iEvt) const
Retrieve the data for a given event.
virtual void finaliseFitResults(const TString &tableFileName)=0
Write the results of the fit into the ntuple.
static LauAbsFitter * fitter()
Method that provides access to the singleton fitter.
Definition: LauFitter.cc:51
void setIntegerBranchValue(const TString &name, Int_t value)
Set value of an integer branch.
virtual void checkInitFitParams()=0
Update initial fit parameters if required.
virtual void printVarsInfo() const
Same as printEventInfo, but printing out the values of the variables in the fit.
Bool_t writeSPlotData_
Option to write sPlot data.
Bool_t validBkgndClass(const TString &className) const
Check if the given background class is in the list.
Bool_t fitToyMCPoissonSmear_
Option to perform Poisson smearing.
TString sWeightBranchName_
The name of the sWeight branch.
LauParameterPList fitVars_
Internal vector of fit parameters.
Bool_t useAsymmFitErrors() const
Report whether or not calculation of asymmetric errors is enabled.
Definition: LauFitObject.hh:57
virtual Bool_t scfDPSmear() const =0
Check if the mis-reconstructed signal is to be smeared in the DP.
void fit(const TString &dataFileName, const TString &dataTreeName, const TString &histFileName, const TString &tableFileNameBase)
Perform the total fit.
UInt_t numberOKFits() const
Access the number of successful fits.
TString fixParamTreeName_
Imported parameters tree name.
UInt_t numberBadFits() const
Access the number of failed fits.
virtual void setSPlotNtupleIntegerBranchValue(const TString &name, Int_t value)
Set the value of an integer branch in the sPlot tree.
UInt_t fitToyMCScale_
The scaling factor (toy vs data statistics)
virtual void fillSPlotNtupleBranches()
Fill the sPlot tuple.
TString fitToyMCFileName_
The output file name for Toy MC.
void writeOutGenResults()
Write out the results from the generation.
void printFormat(std::ostream &stream, Double_t value) const
Method to choose the printing format to a specified level of precision.
Definition: LauPrint.cc:46
virtual LauSPlot::TwoDMap twodimPDFs() const =0
Returns the species and variables for all 2D PDFs in the fit.
Double_t getLogLikelihoodPenalty()
Calculate the penalty terms to the log likelihood from Gaussian constraints.
Bool_t haveBranch(const TString &name) const
Check if the named branch is stored.
virtual void savePDFPlots(const TString &label)=0
Save the pdf Plots for all the resonances of experiment number fitExp.
const TString & bkgndClassName(UInt_t classID) const
Get the name of a background class from the number.
void startNewFit(const UInt_t nPars, const UInt_t nFreePars)
Indicate the start of a new fit.
Definition: LauFitObject.cc:60
LauParameterPSet resVars_
Internal set of fit parameters upon which the DP normalisation depends.
virtual Bool_t useAsymmFitErrors() const =0
Determine whether calculation of asymmetric errors is enabled.
Bool_t enableEmbedding() const
Determine whether embedding of events is enabled in the generation.
Class for defining the abstract interface for PDF classes.
Definition: LauAbsPdf.hh:55
Double_t getLogLikelihood(UInt_t iStart, UInt_t iEnd)
Calculate the sum of the log-likelihood over the specified events.
virtual LauSPlot::NameSet variableNames() const =0
Returns the names of all variables in the fit.
Double_t value() const
The value of the parameter.
UInt_t nEvents() const
Retrieve the number of events.
virtual void addGenNtupleDoubleBranch(const TString &name)
Add a branch to the gen tree for storing a double.
void updateFitParameters(LauPdfList &pdfList)
Update the fit parameters for all PDFs in the list.
void addDoubleBranch(const TString &name)
Add double branch to tree.
Definition: LauGenNtuple.cc:96
std::set< TString > fixParamNames_
Imported parameter names.
UInt_t iExpt() const
Obtain the number of the current experiment.
TString outputTableName_
The output table name.
virtual void setParsFromMinuit(Double_t *par, Int_t npar)
This function sets the parameter values from Minuit.
virtual void initialise()=0
Initialise the fit par vectors.
UInt_t firstExpt() const
Obtain the number of the first experiment.
Class to store the input fit variables.
virtual const TMatrixD & covarianceMatrix() const =0
Retrieve the fit covariance matrix.
virtual void setupResultsOutputs(const TString &histFileName, const TString &tableFileName)
Setup saving of fit results to ntuple/LaTeX table etc.
virtual void calculateSPlotData()
Calculate the sPlot data.
File containing declaration of LauGenNtuple class.
void setDoubleBranchValue(const TString &name, Double_t value)
Set value of a double branch.
std::vector< Double_t > sWeights_
The vector of sWeights.