laura is hosted by Hepforge, IPPP Durham
Laura++  v1r0
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 - 2013.
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 
47 #ifndef LAU_ABS_FIT_MODEL
48 #define LAU_ABS_FIT_MODEL
49 
50 #include "TString.h"
51 #include "TStopwatch.h"
52 #include "TVectorDfwd.h"
53 
54 #include <vector>
55 #include <iosfwd>
56 
57 // LauSPlot included to get LauSPlot::NameSet typedef
58 #include "LauSPlot.hh"
59 
60 class TMessage;
61 class TMonitor;
62 class TSocket;
63 class TTree;
64 class LauAbsPdf;
65 class LauComplex;
66 class LauFitDataTree;
67 class LauFitNtuple;
68 class LauGenNtuple;
69 class LauParameter;
70 
71 class LauAbsFitModel : public TObject {
72 
73  public:
76 
78  virtual ~LauAbsFitModel();
79 
81  Bool_t useDP() const { return usingDP_; }
82 
84 
87  void useDP(Bool_t usingDP) { usingDP_ = usingDP; }
88 
90  Bool_t doSFit() const { return doSFit_; }
91 
93 
97  void doSFit( const TString& sWeightBranchName, Double_t scaleFactor = 1.0 );
98 
100  Bool_t doEMLFit() const {return emlFit_;}
101 
103 
106  void doEMLFit(Bool_t emlFit) {emlFit_ = emlFit;}
107 
109  Bool_t doPoissonSmearing() const {return poissonSmear_;}
110 
112 
115  void doPoissonSmearing(Bool_t poissonSmear) {poissonSmear_ = poissonSmear;}
116 
118  Bool_t enableEmbedding() const {return enableEmbedding_;}
119 
121 
124  void enableEmbedding(Bool_t enable) {enableEmbedding_ = enable;}
125 
127  Bool_t twoStageFit() const {return twoStageFit_;}
128 
130 
136  void twoStageFit(Bool_t doTwoStageFit) {twoStageFit_ = doTwoStageFit;}
137 
139  Bool_t useAsymmFitErrors() const {return useAsymmFitErrors_;}
140 
142 
145  void useAsymmFitErrors(Bool_t useAsymmErrors) {useAsymmFitErrors_ = useAsymmErrors;}
146 
148  Bool_t writeLatexTable() const {return writeLatexTable_;}
149 
151 
154  void writeLatexTable(Bool_t writeTable) {writeLatexTable_ = writeTable;}
155 
157 
163  void writeSPlotData(const TString& fileName, const TString& treeName, Bool_t storeDPEfficiency, const TString& verbosity = "q");
164 
166  Bool_t writeSPlotData() const {return writeSPlotData_;}
167 
169  Bool_t storeDPEff() const {return storeDPEff_;}
170 
172  Bool_t useRandomInitFitPars() const {return randomFit_;}
173 
175  void useRandomInitFitPars(Bool_t boolean) {randomFit_ = boolean;}
176 
178 
182  void setNExpts(UInt_t nExperiments, UInt_t firstExperiment = 0) {
183  nExpt_ = nExperiments;
184  firstExpt_ = firstExperiment;
185  }
186 
188  UInt_t eventsPerExpt() const {return evtsPerExpt_;}
189 
191  UInt_t nExpt() const {return nExpt_;}
192 
194  UInt_t firstExpt() const {return firstExpt_;}
195 
197  UInt_t iExpt() const {return iExpt_;}
198 
200 
203  virtual void setBkgndClassNames( const std::vector<TString>& names );
204 
206  inline UInt_t nBkgndClasses() const {return bkgndClassNames_.size();}
207 
209 
213  virtual void setNSigEvents(Double_t nSigEvents, Bool_t fixNSigEvents) = 0;
214 
216 
222  virtual void setNBkgndEvents(const TString& bgClass, Double_t nBkgndEvents, Bool_t fixNBkgndEvents) = 0;
223 
225 
233  void compareFitData(UInt_t toyMCScale = 10, const TString& mcFileName = "fitToyMC.root",
234  const TString& tableFileName = "fitToyMCTable.tex", Bool_t poissonSmearing = kTRUE);
235 
237 
241  virtual void weightEvents( const TString& dataFileName, const TString& dataTreeName ) = 0;
242 
244 
251  void run(const TString& applicationCode, const TString& dataFileName = "", const TString& dataTreeName = "",
252  const TString& histFileName = "", const TString& tableFileName = "");
253 
255 
263  void runMaster(const TString& applicationCode, const TString& dataFileName = "", const TString& dataTreeName = "",
264  const TString& histFileName = "", const TString& tableFileName = "",
265  const UInt_t nSlaves = 0);
266 
268 
276  void runSlave(const TString& applicationCode, const TString& dataFileName = "", const TString& dataTreeName = "",
277  const TString& histFileName = "", const TString& tableFileName = "",
278  const TString& addressMaster = "");
279 
281 
288  virtual void setParsFromMinuit(Double_t* par, Int_t npar);
289 
291 
299  Double_t calculateLogLike( Int_t npar, Double_t *par, Int_t iflag );
300 
301  protected:
302 
303  // Some typedefs
304 
306  typedef std::vector<LauAbsPdf*> LauPdfList;
308  typedef std::vector<LauParameter*> LauParameterPList;
310  typedef std::vector<LauParameter> LauParameterList;
312  typedef std::map<UInt_t,TString> LauBkgndClassMap;
313 
315  void initSockets();
316 
318  void clearFitParVectors();
319 
321  void clearExtraVarVectors();
322 
324 
330  virtual void generate(const TString& dataFileName, const TString& dataTreeName, const TString& histFileName, const TString& tableFileNameBase);
331 
333 
336  virtual Bool_t genExpt() = 0;
337 
339 
345  void fit(const TString& dataFileName, const TString& dataTreeName, const TString& histFileName, const TString& tableFileNameBase);
346 
348 
354  void fitMaster(const TString& dataFileName, const TString& dataTreeName, const TString& histFileName, const TString& tableFileNameBase);
355 
357 
361  void fitSlave(const TString& dataFileName, const TString& dataTreeName);
362 
364  void fitExpt();
365 
367 
370  Bool_t runMinimisation();
371 
373 
377  void createFitToyMC(const TString& mcFileName, const TString& tableFileName);
378 
380 
384  Bool_t cacheFitData(const TString& dataFileName, const TString& dataTreeName);
385 
387  virtual void cacheInputFitVars() = 0;
388 
390  virtual void cacheInputSWeights();
391 
393 
399  virtual void initialise() = 0;
400 
402  virtual void initialiseDPModels() = 0;
403 
411  virtual void updateCoeffs() = 0;
412 
414  virtual void propagateParUpdates() = 0;
415 
417  Double_t getTotNegLogLikelihood();
418 
420 
424  Double_t getLogLikelihood( UInt_t iStart, UInt_t iEnd );
425 
427 
430  virtual Double_t getTotEvtLikelihood(UInt_t iEvt) = 0;
431 
433  virtual Double_t getEventSum() const = 0;
434 
436 
439  virtual void printEventInfo(UInt_t iEvt) const;
440 
442  virtual void printVarsInfo() const;
443 
445  virtual void checkInitFitParams() = 0;
446 
448 
451  virtual void finaliseFitResults(const TString& tableFileName) = 0;
452 
454 
457  virtual void writeOutTable(const TString& outputFile) = 0;
458 
460  virtual void storePerEvtLlhds() = 0;
461 
463  virtual void writeOutAllFitResults();
464 
466  virtual void calculateSPlotData();
467 
469  void setGenValues();
470 
472  virtual void setupBkgndVectors() = 0;
473 
475 
479  Bool_t validBkgndClass( const TString& className ) const;
480 
482 
486  UInt_t bkgndClassID( const TString& className ) const;
487 
489 
493  const TString& bkgndClassName( UInt_t classID ) const;
494 
496  void eventsPerExpt(UInt_t nEvents) {evtsPerExpt_ = nEvents;}
497 
499  virtual void setupGenNtupleBranches() = 0;
500 
502 
505  virtual void addGenNtupleIntegerBranch(const TString& name);
506 
508 
511  virtual void addGenNtupleDoubleBranch(const TString& name);
512 
514 
518  virtual void setGenNtupleIntegerBranchValue(const TString& name, Int_t value);
519 
521 
525  virtual void setGenNtupleDoubleBranchValue(const TString& name, Double_t value);
526 
528 
532  virtual Int_t getGenNtupleIntegerBranchValue(const TString& name) const;
533 
535 
539  virtual Double_t getGenNtupleDoubleBranchValue(const TString& name) const;
540 
542  virtual void fillGenNtupleBranches();
543 
545  virtual void setupSPlotNtupleBranches() = 0;
546 
548 
551  virtual void addSPlotNtupleIntegerBranch(const TString& name);
552 
554 
557  virtual void addSPlotNtupleDoubleBranch(const TString& name);
558 
560 
564  virtual void setSPlotNtupleIntegerBranchValue(const TString& name, Int_t value);
565 
567 
571  virtual void setSPlotNtupleDoubleBranchValue(const TString& name, Double_t value);
572 
574  virtual void fillSPlotNtupleBranches();
575 
577  virtual LauSPlot::NameSet variableNames() const = 0;
578 
580  virtual LauSPlot::NumbMap freeSpeciesNames() const = 0;
581 
583  virtual LauSPlot::NumbMap fixdSpeciesNames() const = 0;
584 
586  virtual LauSPlot::TwoDMap twodimPDFs() const = 0;
587 
589  virtual Bool_t splitSignal() const = 0;
590 
592  virtual Bool_t scfDPSmear() const = 0;
593 
595 
599  UInt_t addFitParameters(LauPdfList& pdfList);
600 
602 
606  void printFitParameters(const LauPdfList& pdfList, std::ostream& fout) const;
607 
609 
612  void updateFitParameters(LauPdfList& pdfList);
613 
615 
619  void cacheInfo(LauPdfList& pdfList, const LauFitDataTree& theData);
620 
622 
626  Double_t prodPdfValue(LauPdfList& pdfList, UInt_t iEvt);
627 
629 
632  Bool_t pdfsDependOnDP() const {return pdfsDependOnDP_;}
633 
635 
638  void pdfsDependOnDP(Bool_t dependOnDP) { pdfsDependOnDP_ = dependOnDP; }
639 
641  const LauParameterPList& fitPars() const {return fitVars_;}
643 
645  const LauParameterList& extraPars() const {return extraVars_;}
647 
649  const LauFitNtuple* fitNtuple() const {return fitNtuple_;}
651 
653  const LauGenNtuple* genNtuple() const {return genNtuple_;}
655 
657  const LauGenNtuple* sPlotNtuple() const {return sPlotNtuple_;}
659 
661  const LauFitDataTree* fitData() const {return inputFitData_;}
663 
665  Double_t nll() const {return NLL_;}
666 
668  Int_t fitStatus() const {return fitStatus_;}
669 
670  private:
671 
672  // Various control booleans
673 
683  Bool_t storeDPEff_;
685  Bool_t randomFit_;
687  Bool_t emlFit_;
689  Bool_t poissonSmear_;
693  Bool_t twoStageFit_;
695  Bool_t usingDP_;
698 
699  // Info on number of experiments and number of events
700 
702  UInt_t firstExpt_;
704  UInt_t nExpt_;
706  UInt_t evtsPerExpt_;
708  UInt_t iExpt_;
709 
712 
715 
716  // Input data and output ntuple
717 
726 
727  // Fit bookeeping variables
728 
730  Int_t fitStatus_;
732  Double_t NLL_;
734  UInt_t numberOKFits_;
736  UInt_t numberBadFits_;
738  UInt_t nParams_;
740  UInt_t nFreeParams_;
742  Double_t worstLogLike_;
744  Bool_t withinMINOS_;
745 
746  // Background class names
750  const TString nullString_;
751 
752  // sFit related variables
753 
755  Bool_t doSFit_;
759  std::vector<Double_t> sWeights_;
762 
763  // Fit timers
764 
766  TStopwatch timer_;
768  TStopwatch cumulTimer_;
769 
770  // Comparison toy MC related variables
771 
780 
781  // sPlot related variables
782 
784  TString sPlotFileName_;
786  TString sPlotTreeName_;
789 
790  // Parallel fitting variables
791 
793  TMonitor * socketMonitor_;
795  TSocket * sMaster_;
797  std::vector<TSocket*> sSlaves_;
799  TVectorD * vectorPar_;
801  TVectorD * vectorRes_;
803  TMessage * messageFromSlave_;
805  TMessage * messageFromMaster_;
807  UInt_t slaveId_;
809  UInt_t nSlaves_;
810 
811 
812  ClassDef(LauAbsFitModel,0) // Abstract interface to fit/toyMC model
813 };
814 
815 #endif
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.
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.
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.
void run(const TString &applicationCode, const TString &dataFileName="", const TString &dataTreeName="", const TString &histFileName="", const TString &tableFileName="")
Start the toy generation / fitting.
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.
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()
Retrieves the total negative log-likelihood.
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.
Bool_t twoStageFit() const
Determine whether the two-stage fit is enabled.
std::vector< TSocket * > sSlaves_
Sockets for each of the slaves to communicate with the master.
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.
TVectorD * vectorPar_
Parameters to send to the slaves.
LauGenNtuple * genNtuple_
The generated ntuple.
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.
TVectorD * vectorRes_
Parameters returned from the slaves.
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.
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.
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()
void runMaster(const TString &applicationCode, const TString &dataFileName="", const TString &dataTreeName="", const TString &histFileName="", const TString &tableFileName="", const UInt_t nSlaves=0)
Start the master process for the generation / fitting.
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 number of slaves.
Bool_t poissonSmear_
Option to perform Poisson smearing.
virtual void setNBkgndEvents(const TString &bgClass, Double_t nBkgndEvents, Bool_t fixNBkgndEvents)=0
Set the number of background events.
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.
TMonitor * socketMonitor_
Parallel setup monitor.
virtual void fillGenNtupleBranches()
Fill the gen tuple branches.
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.
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 initSockets()
Initialise socket connections for the slaves.
Bool_t withinMINOS_
Flag to indicate if MINOS is currently running.
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.
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.
Abstract interface to the fitting and toy MC model.
Double_t calculateLogLike(Int_t npar, Double_t *par, Int_t iflag)
Calculate the new value of the negative log likelihood.
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.
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.
void setNExpts(UInt_t nExperiments, UInt_t firstExperiment=0)
Set the number of experiments and the first experiment.
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.
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.
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:31
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.
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.
void fitSlave(const TString &dataFileName, const TString &dataTreeName)
Slaves required when performing a simultaneous fit.
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.
TMessage * messageFromSlave_
Message from slaves to the master.
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 twoStageFit(Bool_t doTwoStageFit)
Turn on or off the two stage fit.
void runSlave(const TString &applicationCode, const TString &dataFileName="", const TString &dataTreeName="", const TString &histFileName="", const TString &tableFileName="", const TString &addressMaster="")
Start the slave process for the generation / fitting.
void fit(const TString &dataFileName, const TString &dataTreeName, const TString &histFileName, const TString &tableFileNameBase)
Perform the total fit.
virtual void setNSigEvents(Double_t nSigEvents, Bool_t fixNSigEvents)=0
Set the number of signal events.
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.
void fitMaster(const TString &dataFileName, const TString &dataTreeName, const TString &histFileName, const TString &tableFileNameBase)
Master fit used when simultaneous fit is required.
UInt_t fitToyMCScale_
The scaling factor (toy vs data statistics)
Bool_t runMinimisation()
Routine to perform the minimisation.
virtual void fillSPlotNtupleBranches()
Fill the sPlot tuple.
virtual void updateCoeffs()=0
TString fitToyMCFileName_
The output file name for Toy MC.
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:41
virtual LauSPlot::TwoDMap twodimPDFs() const =0
Returns the species and variables for all 2D PDFs in the fit.
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:40
Double_t getLogLikelihood(UInt_t iStart, UInt_t iEnd)
Calculate the sum of the log-likelihood over the specified events.
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.
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
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.