laura is hosted by Hepforge, IPPP Durham
Laura++  v2r2p1
A maximum likelihood fitting package for performing Dalitz-plot analysis.
LauAbsFitModel.hh
Go to the documentation of this file.
1 
2 // Copyright University of Warwick 2004 - 2014.
3 // Distributed under the Boost Software License, Version 1.0.
4 // (See accompanying file LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
5 
6 // Authors:
7 // Thomas Latham
8 // John Back
9 // Paul Harrison
10 
48 #ifndef LAU_ABS_FIT_MODEL
49 #define LAU_ABS_FIT_MODEL
50 
51 #include "TMatrixD.h"
52 #include "TString.h"
53 #include "TStopwatch.h"
54 #include "TVectorDfwd.h"
55 
56 #include <vector>
57 #include <iosfwd>
58 
59 #include "LauFitObject.hh"
60 #include "LauFormulaPar.hh"
61 // LauSPlot included to get LauSPlot::NameSet typedef
62 #include "LauSPlot.hh"
63 
64 class TMessage;
65 class TMonitor;
66 class TSocket;
67 class TTree;
68 class LauAbsCoeffSet;
69 class LauAbsPdf;
70 class LauComplex;
71 class LauFitDataTree;
72 class LauFitNtuple;
73 class LauGenNtuple;
74 class LauAbsRValue;
75 class LauParameter;
76 
77 class LauAbsFitModel : public LauFitObject {
78 
79  public:
82 
84  virtual ~LauAbsFitModel();
85 
87  Bool_t useDP() const { return usingDP_; }
88 
90 
93  void useDP(Bool_t usingDP) { usingDP_ = usingDP; }
94 
96  Bool_t doSFit() const { return doSFit_; }
97 
99 
103  void doSFit( const TString& sWeightBranchName, Double_t scaleFactor = 1.0 );
104 
106  Bool_t doEMLFit() const {return emlFit_;}
107 
109 
112  void doEMLFit(Bool_t emlFit) {emlFit_ = emlFit;}
113 
115  virtual Bool_t twoStageFit() const {return twoStageFit_;}
116 
118 
128  virtual void twoStageFit(Bool_t doTwoStageFit) {twoStageFit_ = doTwoStageFit;}
129 
131  Bool_t useAsymmFitErrors() const {return useAsymmFitErrors_;}
132 
134 
137  void useAsymmFitErrors(Bool_t useAsymmErrors) {useAsymmFitErrors_ = useAsymmErrors;}
138 
140  Bool_t doPoissonSmearing() const {return poissonSmear_;}
141 
143 
146  void doPoissonSmearing(Bool_t poissonSmear) {poissonSmear_ = poissonSmear;}
147 
149  Bool_t enableEmbedding() const {return enableEmbedding_;}
150 
152 
155  void enableEmbedding(Bool_t enable) {enableEmbedding_ = enable;}
156 
158 
164  virtual void withinAsymErrorCalc(Bool_t inAsymErrCalc) {withinAsymErrorCalc_ = inAsymErrCalc;}
165 
167  Bool_t writeLatexTable() const {return writeLatexTable_;}
168 
170 
173  void writeLatexTable(Bool_t writeTable) {writeLatexTable_ = writeTable;}
174 
176 
182  void writeSPlotData(const TString& fileName, const TString& treeName, Bool_t storeDPEfficiency, const TString& verbosity = "q");
183 
185  Bool_t writeSPlotData() const {return writeSPlotData_;}
186 
188  Bool_t storeDPEff() const {return storeDPEff_;}
189 
191  Bool_t useRandomInitFitPars() const {return randomFit_;}
192 
194  void useRandomInitFitPars(Bool_t boolean) {randomFit_ = boolean;}
195 
197 
201  void setNExpts(UInt_t nExperiments, UInt_t firstExperiment = 0) {
202  nExpt_ = nExperiments;
203  firstExpt_ = firstExperiment;
204  }
205 
207  UInt_t eventsPerExpt() const {return evtsPerExpt_;}
208 
210  UInt_t nExpt() const {return nExpt_;}
211 
213  UInt_t firstExpt() const {return firstExpt_;}
214 
216  UInt_t iExpt() const {return iExpt_;}
217 
219 
222  virtual void setBkgndClassNames( const std::vector<TString>& names );
223 
225  inline UInt_t nBkgndClasses() const {return bkgndClassNames_.size();}
226 
228 
231  virtual void setNSigEvents(LauParameter* nSigEvents) = 0;
232 
234 
240  virtual void setNBkgndEvents(LauParameter* nBkgndEvents) = 0;
241 
243 
246  virtual void setAmpCoeffSet(LauAbsCoeffSet* coeffSet) = 0;
247 
249 
257  void compareFitData(UInt_t toyMCScale = 10, const TString& mcFileName = "fitToyMC.root",
258  const TString& tableFileName = "fitToyMCTable.tex", Bool_t poissonSmearing = kTRUE);
259 
261 
265  virtual void weightEvents( const TString& dataFileName, const TString& dataTreeName ) = 0;
266 
268 
275  void run(const TString& applicationCode, const TString& dataFileName, const TString& dataTreeName,
276  const TString& histFileName, const TString& tableFileName = "");
277 
279 
287  void runSlave(const TString& dataFileName, const TString& dataTreeName,
288  const TString& histFileName, const TString& tableFileName = "",
289  const TString& addressMaster = "localhost", const UInt_t portMaster = 9090);
290 
292 
299  virtual void setParsFromMinuit(Double_t* par, Int_t npar);
300 
302 
306  Double_t getTotNegLogLikelihood();
307 
309 
315  virtual void addConstraint(const TString& formula, const std::vector<TString>& pars, const Double_t mean, const Double_t width);
316 
317  protected:
318 
319  // Some typedefs
320 
322  typedef std::vector<LauAbsPdf*> LauPdfList;
324  typedef std::vector<LauParameter*> LauParameterPList;
326  typedef std::vector<LauAbsRValue*> LauAbsRValuePList;
328  typedef std::vector<LauParameter> LauParameterList;
330  typedef std::map<UInt_t,TString> LauBkgndClassMap;
331 
333  void initSockets();
334 
336  void clearFitParVectors();
337 
339  void clearExtraVarVectors();
340 
342 
348  virtual void generate(const TString& dataFileName, const TString& dataTreeName, const TString& histFileName, const TString& tableFileNameBase);
349 
351 
354  virtual Bool_t genExpt() = 0;
355 
357 
363  void fit(const TString& dataFileName, const TString& dataTreeName, const TString& histFileName, const TString& tableFileNameBase);
364 
366 
372  void fitSlave(const TString& dataFileName, const TString& dataTreeName, const TString& histFileName, const TString& tableFileNameBase);
373 
375  void fitExpt();
376 
378 
381  Bool_t runMinimisation();
382 
384 
388  void createFitToyMC(const TString& mcFileName, const TString& tableFileName);
389 
391 
395  Bool_t cacheFitData(const TString& dataFileName, const TString& dataTreeName);
396 
398  virtual void cacheInputFitVars() = 0;
399 
401  virtual void cacheInputSWeights();
402 
404 
410  virtual void initialise() = 0;
411 
413  virtual void initialiseDPModels() = 0;
414 
422  virtual void updateCoeffs() = 0;
423 
425  virtual void propagateParUpdates() = 0;
426 
428 
432  Double_t getLogLikelihood( UInt_t iStart, UInt_t iEnd );
433 
435  Double_t getLogLikelihoodPenalty();
436 
438 
441  virtual Double_t getTotEvtLikelihood(UInt_t iEvt) = 0;
442 
444  virtual Double_t getEventSum() const = 0;
445 
447 
450  virtual void printEventInfo(UInt_t iEvt) const;
451 
453  virtual void printVarsInfo() const;
454 
456  virtual void checkInitFitParams() = 0;
457 
459 
462  virtual void finaliseFitResults(const TString& tableFileName) = 0;
463 
465 
468  virtual void writeOutTable(const TString& outputFile) = 0;
469 
471  virtual void storePerEvtLlhds() = 0;
472 
474  virtual void writeOutAllFitResults();
475 
477  virtual void calculateSPlotData();
478 
480  void setGenValues();
481 
483  virtual void setupBkgndVectors() = 0;
484 
486 
490  Bool_t validBkgndClass( const TString& className ) const;
491 
493 
497  UInt_t bkgndClassID( const TString& className ) const;
498 
500 
504  const TString& bkgndClassName( UInt_t classID ) const;
505 
507  void eventsPerExpt(UInt_t nEvents) {evtsPerExpt_ = nEvents;}
508 
510  virtual void setupGenNtupleBranches() = 0;
511 
513 
516  virtual void addGenNtupleIntegerBranch(const TString& name);
517 
519 
522  virtual void addGenNtupleDoubleBranch(const TString& name);
523 
525 
529  virtual void setGenNtupleIntegerBranchValue(const TString& name, Int_t value);
530 
532 
536  virtual void setGenNtupleDoubleBranchValue(const TString& name, Double_t value);
537 
539 
543  virtual Int_t getGenNtupleIntegerBranchValue(const TString& name) const;
544 
546 
550  virtual Double_t getGenNtupleDoubleBranchValue(const TString& name) const;
551 
553  virtual void fillGenNtupleBranches();
554 
556  virtual void setupSPlotNtupleBranches() = 0;
557 
559 
562  virtual void addSPlotNtupleIntegerBranch(const TString& name);
563 
565 
568  virtual void addSPlotNtupleDoubleBranch(const TString& name);
569 
571 
575  virtual void setSPlotNtupleIntegerBranchValue(const TString& name, Int_t value);
576 
578 
582  virtual void setSPlotNtupleDoubleBranchValue(const TString& name, Double_t value);
583 
585  virtual void fillSPlotNtupleBranches();
586 
588  virtual LauSPlot::NameSet variableNames() const = 0;
589 
591  virtual LauSPlot::NumbMap freeSpeciesNames() const = 0;
592 
594  virtual LauSPlot::NumbMap fixdSpeciesNames() const = 0;
595 
597  virtual LauSPlot::TwoDMap twodimPDFs() const = 0;
598 
600  virtual Bool_t splitSignal() const = 0;
601 
603  virtual Bool_t scfDPSmear() const = 0;
604 
606 
610  UInt_t addFitParameters(LauPdfList& pdfList);
611 
613  void addConParameters();
614 
616 
620  void printFitParameters(const LauPdfList& pdfList, std::ostream& fout) const;
621 
623 
626  void updateFitParameters(LauPdfList& pdfList);
627 
629 
633  void cacheInfo(LauPdfList& pdfList, const LauFitDataTree& theData);
634 
636 
640  Double_t prodPdfValue(LauPdfList& pdfList, UInt_t iEvt);
641 
643 
646  Bool_t pdfsDependOnDP() const {return pdfsDependOnDP_;}
647 
649 
652  void pdfsDependOnDP(Bool_t dependOnDP) { pdfsDependOnDP_ = dependOnDP; }
653 
655  const LauParameterPList& fitPars() const {return fitVars_;}
657 
659  const LauParameterList& extraPars() const {return extraVars_;}
661 
663  const LauAbsRValuePList& conPars() const {return conVars_;}
665 
667  const LauFitNtuple* fitNtuple() const {return fitNtuple_;}
669 
671  const LauGenNtuple* genNtuple() const {return genNtuple_;}
673 
675  const LauGenNtuple* sPlotNtuple() const {return sPlotNtuple_;}
677 
679  const LauFitDataTree* fitData() const {return inputFitData_;}
681 
683  Double_t nll() const {return NLL_;}
684 
686  Int_t fitStatus() const {return fitStatus_;}
687 
689  const TMatrixD& covarianceMatrix() const {return covMatrix_;}
690 
691  private:
692  // Setup a struct to store information on constrained fit parameters
699  TString formula_;
701  std::vector<TString> conPars_;
703  Double_t mean_;
705  Double_t width_;
706  };
707 
709  std::vector<StoreConstraints> storeCon_;
710 
711  // Various control booleans
712 
714  Bool_t twoStageFit_;
724  Bool_t storeDPEff_;
726  Bool_t randomFit_;
728  Bool_t emlFit_;
730  Bool_t poissonSmear_;
734  Bool_t usingDP_;
737 
738  // Info on number of experiments and number of events
739 
741  UInt_t firstExpt_;
743  UInt_t nExpt_;
745  UInt_t evtsPerExpt_;
747  UInt_t iExpt_;
748 
751 
754 
757 
758  // Input data and output ntuple
759 
768 
769  // Fit bookeeping variables
770 
772  Int_t fitStatus_;
774  Double_t NLL_;
776  TMatrixD covMatrix_;
778  UInt_t numberOKFits_;
780  UInt_t numberBadFits_;
782  UInt_t nParams_;
784  UInt_t nFreeParams_;
786  Double_t worstLogLike_;
789 
790  // Background class names
794  const TString nullString_;
795 
796  // sFit related variables
797 
799  Bool_t doSFit_;
803  std::vector<Double_t> sWeights_;
806 
807  // Fit timers
808 
810  TStopwatch timer_;
812  TStopwatch cumulTimer_;
813 
814  // Comparison toy MC related variables
815 
824 
825  // sPlot related variables
826 
828  TString sPlotFileName_;
830  TString sPlotTreeName_;
833 
834  // Parallel fitting variables
835 
837  TSocket* sMaster_;
841  UInt_t slaveId_;
843  UInt_t nSlaves_;
845  Double_t* parValues_;
846 
847  ClassDef(LauAbsFitModel,0) // Abstract interface to fit/toyMC model
848 };
849 
850 #endif
TString fitToyMCTableName_
The output table name for Toy MC.
Double_t mean_
The mean value of the Gaussian constraint to be applied.
virtual void addSPlotNtupleIntegerBranch(const TString &name)
Add a branch to the sPlot tree for storing an integer.
virtual Bool_t twoStageFit() const
Determine whether the two-stage fit is enabled.
virtual Bool_t genExpt()=0
The method that actually generates the toy MC events for the given experiment.
File containing declaration of LauFormulaPar class.
void createFitToyMC(const TString &mcFileName, const TString &tableFileName)
Create a toy MC sample from the fitted parameters.
UInt_t firstExpt() const
Obtain the number of the first experiment.
std::vector< LauParameter > LauParameterList
List of parameters.
virtual Double_t getEventSum() const =0
Returns the sum of the expected events over all hypotheses; used in the EML fit scenario.
UInt_t eventsPerExpt() const
Obtain the total number of events in the current experiment.
void setGenValues()
Make sure all parameters hold their genValue as the current value.
TString sPlotFileName_
The name of the sPlot file.
Double_t * parValues_
Parameter values array (for reading from the master)
virtual void writeOutTable(const TString &outputFile)=0
Write the latex table.
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.
UInt_t nParams_
The number of fit parameters.
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.
virtual void setGenNtupleIntegerBranchValue(const TString &name, Int_t value)
Set the value of an integer branch in the gen tree.
void writeLatexTable(Bool_t writeTable)
Turn on or off the writing out of the latex table.
virtual void setupBkgndVectors()=0
Method to set up the storage for background-related quantities called by setBkgndClassNames.
LauGenNtuple * genNtuple_
The generated ntuple.
Bool_t withinAsymErrorCalc_
Flag to indicate if the asymmetric error calculation (e.g. MINOS) is currently running.
TString formula_
The formula to be used in the LauFormulaPar.
Bool_t useAsymmFitErrors() const
Determine whether calculation of asymmetric errors is enabled.
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.
Bool_t useDP() const
Is the Dalitz plot term in the likelihood.
void doPoissonSmearing(Bool_t poissonSmear)
Turn Poisson smearing (for the toy MC generation) on or off.
std::multimap< TString, std::pair< TString, TString > > TwoDMap
Type to associate the name of the species that have 2D PDFs with the names of the two variables invol...
Definition: LauSPlot.hh:68
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.
Bool_t randomFit_
Option to randomise the initial values of the fit parameters.
std::vector< TString > conPars_
The list of LauParameter names to be used in the LauFormulaPar.
virtual void cacheInputSWeights()
Cache the value of the sWeights to be used in the sFit.
virtual void setGenNtupleDoubleBranchValue(const TString &name, Double_t value)
Set the value of a double branch in the gen tree.
virtual void setNSigEvents(LauParameter *nSigEvents)=0
Set the number of signal events.
void eventsPerExpt(UInt_t nEvents)
Set the number of events in the current experiment.
TSocket * sMaster_
A socket to enable parallel setup.
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.
Int_t fitStatus() const
Access the fit status information.
Bool_t compareFitData_
Option to make toy from 1st successful experiment.
Bool_t writeLatexTable_
Option to output a Latex format table.
Double_t prodPdfValue(LauPdfList &pdfList, UInt_t iEvt)
Calculate the product of the per-event likelihoods of the PDFs in the list.
Bool_t emlFit_
Option to perform an extended ML fit.
LauFitDataTree * fitData()
virtual void addConstraint(const TString &formula, const std::vector< TString > &pars, const Double_t mean, const Double_t width)
Store constraint information for fit parameters.
std::vector< LauAbsPdf * > LauPdfList
List of Pdfs.
virtual Int_t getGenNtupleIntegerBranchValue(const TString &name) const
Get the value of an integer branch in the gen tree.
UInt_t nSlaves_
The total number of slaves.
Bool_t poissonSmear_
Option to perform Poisson smearing.
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.
TMessage * messageFromMaster_
Message from master to the slaves.
Bool_t twoStageFit_
Option to perform a two stage fit.
virtual void setupSPlotNtupleBranches()=0
Setup the branches of the sPlot tuple.
virtual void fillGenNtupleBranches()
Fill the gen tuple branches.
void fitSlave(const TString &dataFileName, const TString &dataTreeName, const TString &histFileName, const TString &tableFileNameBase)
Slaves required when performing a simultaneous fit.
Bool_t useAsymmFitErrors_
Option to use asymmetric errors.
virtual void setSPlotNtupleDoubleBranchValue(const TString &name, Double_t value)
Set the value of a double branch in the sPlot tree.
void printFitParameters(const LauPdfList &pdfList, std::ostream &fout) const
Print the fit parameters for all PDFs in the list.
Bool_t enableEmbedding_
Option to enable embedding.
std::vector< LauAbsRValue * > LauAbsRValuePList
List of parameter pointers.
UInt_t slaveId_
Slave id number.
LauParameterList & extraPars()
LauFitNtuple * fitNtuple()
virtual Double_t getTotEvtLikelihood(UInt_t iEvt)=0
Calculates the likelihood for a given event.
TStopwatch cumulTimer_
The total fit timer.
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 initSockets()
Initialise socket connections for the slaves.
LauFitDataTree * inputFitData_
The input data.
const LauParameterList & extraPars() const
Access the extra variables.
void doEMLFit(Bool_t emlFit)
Choice to perform an extended maximum likelihood fit.
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.
LauParameterList extraVars_
Extra variables that aren&#39;t in the fit but are stored in the ntuple.
std::vector< StoreConstraints > storeCon_
Store the constraints for fit parameters until initialisation is complete.
TMatrixD covMatrix_
The fit covariance matrix.
Abstract interface to the fitting and toy MC model.
UInt_t iExpt() const
Obtain the number of the current experiment.
UInt_t numberBadFits_
The number of bad fits.
virtual void addGenNtupleIntegerBranch(const TString &name)
Add a branch to the gen tree for storing an integer.
const LauFitNtuple * fitNtuple() const
Access the fit ntuple.
UInt_t nExpt() const
Obtain the number of experiments.
virtual void withinAsymErrorCalc(Bool_t inAsymErrCalc)
Mark that the fit is calculating asymmetric errors.
Bool_t doPoissonSmearing() const
Determine whether Poisson smearing is enabled for the toy MC generation.
virtual Double_t getGenNtupleDoubleBranchValue(const TString &name) const
Get the value of a double branch in the gen tree.
const LauAbsRValuePList & conPars() const
Access the Gaussian constrained variables.
void setNExpts(UInt_t nExperiments, UInt_t firstExperiment=0)
Set the number of experiments and the first experiment.
virtual void setAmpCoeffSet(LauAbsCoeffSet *coeffSet)=0
Set the DP amplitude coefficients.
virtual void weightEvents(const TString &dataFileName, const TString &dataTreeName)=0
Reweighting - allows e.g. MC events to be weighted by the DP model.
std::map< TString, Double_t > NumbMap
Type to associate a category name with a double precision number, e.g. a yield or PDF value for a giv...
Definition: LauSPlot.hh:62
Bool_t doEMLFit() const
Determine whether an extended maximum likelihood fit it being performed.
std::map< UInt_t, TString > LauBkgndClassMap
A type to store background classes.
virtual void twoStageFit(Bool_t doTwoStageFit)
Turn on or off the two stage fit.
UInt_t numberOKFits_
The number of good fits.
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:32
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. ...
LauParameterPList & fitPars()
Double_t worstLogLike_
The worst LL value found so far.
TString sPlotTreeName_
The name of the sPlot tree.
virtual void initialiseDPModels()=0
Initialise the DP models.
Class for defining the abstract interface for complex coefficient classes.
Double_t nll() const
Access the current NLL value.
LauAbsFitModel()
Constructor.
LauGenNtuple * sPlotNtuple_
The sPlot ntuple.
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:33
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.
Struct to store constraint information until the fit is run.
virtual void setNBkgndEvents(LauParameter *nBkgndEvents)=0
Set the number of background events.
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.
const LauGenNtuple * genNtuple() const
Access the gen ntuple.
LauGenNtuple * sPlotNtuple()
LauGenNtuple * genNtuple()
Double_t sWeightScaleFactor_
The sWeight scaling factor.
Bool_t storeDPEff() const
Determine whether the efficiency information should be stored in the sPlot ntuple.
const LauGenNtuple * sPlotNtuple() const
Access the sPlot ntuple.
UInt_t bkgndClassID(const TString &className) const
The number assigned to a background class.
virtual void writeOutAllFitResults()
Write out any fit results.
Bool_t useRandomInitFitPars() const
Determine whether the initial values of the fit parameters, in particular the isobar coefficient para...
virtual void finaliseFitResults(const TString &tableFileName)=0
Write the results of the fit into the ntuple.
const LauFitDataTree * fitData() const
Access the data store.
UInt_t firstExpt_
The number of the first experiment to consider.
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.
UInt_t iExpt_
The number of the current experiment.
LauParameterPList fitVars_
Internal vectors of fit parameters.
virtual Bool_t scfDPSmear() const =0
Check if the mis-reconstructed signal is to be smeared in the DP.
Bool_t cacheFitData(const TString &dataFileName, const TString &dataTreeName)
Store variables from the input file into the internal data storage.
void fit(const TString &dataFileName, const TString &dataTreeName, const TString &histFileName, const TString &tableFileNameBase)
Perform the total fit.
void enableEmbedding(Bool_t enable)
Turn on or off embedding of events in the generation.
void pdfsDependOnDP(Bool_t dependOnDP)
Do any of the PDFs have a dependence on the DP?
virtual void setSPlotNtupleIntegerBranchValue(const TString &name, Int_t value)
Set the value of an integer branch in the sPlot tree.
LauAbsRValuePList & conPars()
UInt_t fitToyMCScale_
The scaling factor (toy vs data statistics)
Bool_t runMinimisation()
Routine to perform the minimisation.
File containing declaration of LauFitObject class.
virtual void fillSPlotNtupleBranches()
Fill the sPlot tuple.
virtual void updateCoeffs()=0
TString fitToyMCFileName_
The output file name for Toy MC.
Double_t width_
The width of the Gaussian constraint to be applied.
Bool_t usingDP_
Option to include the DP as part of the fit.
Class to store the results from the fit into an ntuple.
Definition: LauFitNtuple.hh:42
void runSlave(const TString &dataFileName, const TString &dataTreeName, const TString &histFileName, const TString &tableFileName="", const TString &addressMaster="localhost", const UInt_t portMaster=9090)
Start the slave process for simultaneous fitting.
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.
Class for defining a complex number.
Definition: LauComplex.hh:47
UInt_t nBkgndClasses() const
Returns the number of background classes.
const TString & bkgndClassName(UInt_t classID) const
Get the name of a background class from the number.
void useRandomInitFitPars(Bool_t boolean)
Randomise the initial values of the fit parameters, in particular the isobar coefficient parameters...
void useDP(Bool_t usingDP)
Switch on/off the Dalitz plot term in the Likelihood (allows fits to other quantities, e.g. B mass)
Bool_t enableEmbedding() const
Determine whether embedding of events is enabled in the generation.
Double_t NLL_
The negative log-likelihood.
Class for defining the abstract interface for PDF classes.
Definition: LauAbsPdf.hh:41
Double_t getLogLikelihood(UInt_t iStart, UInt_t iEnd)
Calculate the sum of the log-likelihood over the specified events.
const TMatrixD & covarianceMatrix() const
Access the fit covariance matrix.
const LauParameterPList & fitPars() const
Access the fit variables.
UInt_t nExpt_
The number of experiments to consider.
LauFitNtuple * fitNtuple_
The fit ntuple.
virtual LauSPlot::NameSet variableNames() const =0
Returns the names of all variables in the fit.
UInt_t evtsPerExpt_
The number of events per experiment.
virtual void addGenNtupleDoubleBranch(const TString &name)
Add a branch to the gen tree for storing a double.
void useAsymmFitErrors(Bool_t useAsymmErrors)
Turn on or off the computation of asymmetric errors (e.g. MINOS routine in Minuit) ...
void updateFitParameters(LauPdfList &pdfList)
Update the fit parameters for all PDFs in the list.
Bool_t pdfsDependOnDP_
Option to state if pdfs depend on DP position.
virtual void setParsFromMinuit(Double_t *par, Int_t npar)
This function sets the parameter values from Minuit.
Pure abstract base class for defining a parameter containing an R value.
Definition: LauAbsRValue.hh:29
Int_t fitStatus_
The status of the fit.
virtual void initialise()=0
Initialise the fit par vectors.
Class to store the input fit variables.
std::set< TString > NameSet
Type to store names, e.g. of the discriminating/control variables.
Definition: LauSPlot.hh:59
The abstract interface for the objects that control the calculation of the likelihood.
Definition: LauFitObject.hh:26
virtual void calculateSPlotData()
Calculate the sPlot data.
std::vector< Double_t > sWeights_
The vector of sWeights.
UInt_t nFreeParams_
The number of free fit parameters.