laura is hosted by Hepforge, IPPP Durham
Laura++  3.6.0
A maximum likelihood fitting package for performing Dalitz-plot analysis.

The coordinator process for simultaneous/combined fits. More...

#include <LauSimFitCoordinator.hh>

Inheritance diagram for LauSimFitCoordinator:
LauFitObject

Public Member Functions

 LauSimFitCoordinator (UInt_t numTasks, UInt_t port=9090)
 Constructor. More...
 
virtual ~LauSimFitCoordinator ()
 Destructor.
 
void runSimFit (const TString &fitNtupleFileName, const Bool_t useAsymmErrors=kFALSE, const Bool_t doTwoStageFit=kFALSE)
 Run the fit. More...
 
virtual void withinAsymErrorCalc (const Bool_t inAsymErrCalc)
 Mark that the fit is calculating asymmetric errors. More...
 
virtual void setParsFromMinuit (Double_t *par, Int_t npar)
 This function sets the parameter values from Minuit. More...
 
virtual Double_t getTotNegLogLikelihood ()
 Calculate the new value of the negative log likelihood. More...
 
virtual void withinAsymErrorCalc (const Bool_t inAsymErrCalc)
 Mark that the fit is calculating asymmetric errors. More...
 
Bool_t withinAsymErrorCalc () const
 Query whether the fit is calculating the asymmetric errors. More...
 
- Public Member Functions inherited from LauFitObject
virtual ~LauFitObject ()=default
 Destructor.
 
void useAsymmFitErrors (Bool_t useAsymmErrors)
 Turn on or off the computation of asymmetric errors (e.g. MINOS routine in Minuit) More...
 
Bool_t useAsymmFitErrors () const
 Report whether or not calculation of asymmetric errors is enabled.
 
void twoStageFit (Bool_t doTwoStageFit)
 Turn on or off the two stage fit. More...
 
Bool_t twoStageFit () const
 Report whether the two-stage fit is enabled.
 
Bool_t withinAsymErrorCalc () const
 Query whether the fit is calculating the asymmetric errors. More...
 
void setNExpts (UInt_t nExperiments, UInt_t firstExperiment, Bool_t toyExpts)
 Set the number of experiments, the first experiment, and whether this is toy. More...
 
UInt_t eventsPerExpt () const
 Obtain the total number of events in the current experiment.
 
UInt_t nExpt () const
 Obtain the number of experiments.
 
UInt_t firstExpt () const
 Obtain the number of the first experiment.
 
UInt_t iExpt () const
 Obtain the number of the current experiment.
 
Bool_t toyExpts () const
 Obtain whether this is toy.
 
void addConstraint (const TString &formula, const std::vector< TString > &pars, const Double_t mean, const Double_t width)
 Store constraint information for fit parameters. More...
 
void addFormulaConstraint (const TString &formula, const std::vector< TString > &pars, const Double_t mean, const Double_t width)
 Store constraint information for fit parameters. More...
 
void addMultiDimConstraint (const std::vector< TString > &pars, const TVectorD &means, const TMatrixD &covMat)
 Store n-dimensional constraint information for fit parameters. More...
 

Protected Member Functions

void printParInfo () const
 Print information on the parameters.
 
void initialise ()
 Initialise.
 
void initSockets ()
 Initialise socket connections for the tasks.
 
void getParametersFromTasks ()
 Determine/update the parameter initial values from all tasks.
 
void getParametersFromTasksFirstTime ()
 Determine the parameter names and initial values from all tasks.
 
void updateParametersFromTasks ()
 Update and verify the parameter initial values from all tasks.
 
void checkParameter (const LauParameter *param, UInt_t index) const
 Check for compatibility between two same-named parameters, which should therefore be identical.
 
void addConParameters ()
 Add parameters to the list of Gaussian constrained parameters.
 
Double_t getLogLikelihoodPenalty ()
 Calculate the penalty terms to the log likelihood from Gaussian constraints.
 
Bool_t readData ()
 Instruct the tasks to read the input data for the given experiment. More...
 
Bool_t cacheInputData ()
 Instruct the tasks to perform the caching. More...
 
void fitExpt ()
 Perform the fit for the current experiment.
 
void checkInitFitParams ()
 Instruct the tasks to update the initial fit parameter values, if required.
 
Bool_t finalise ()
 Return the final parameters to the tasks and instruct them to perform their finalisation.
 
Bool_t writeOutResults ()
 Instruct the tasks to write out the fit results.
 
- Protected Member Functions inherited from LauFitObject
 LauFitObject ()
 Constructor.
 
Bool_t checkRepetition (const std::vector< TString > &names, const ConstraintType conType)
 Check if parameters names for constraints have already been used elsewhere. More...
 
void generateConstraintMeans (std::vector< LauAbsRValue * > &conVars)
 Generate per-experiment mean for each Gaussian constraint. More...
 
const std::vector< FormulaConstraint > & formulaConstraints () const
 Const access to the formula constraints store.
 
std::vector< FormulaConstraint > & formulaConstraints ()
 Access to the formula constraints store.
 
const std::vector< MultiDimConstraint > & multiDimConstraints () const
 Const access to the ND constraints store.
 
std::vector< MultiDimConstraint > & multiDimConstraints ()
 Access to the ND constraints store.
 
const std::set< TString > & formulaConstrainedPars () const
 Const access to the parameter names used in formula constraints.
 
std::set< TString > & formulaConstrainedPars ()
 Access to the parameter names used in formula constraints.
 
const std::set< TString > & multiDimConstrainedPars () const
 Const access to the parameter names used in ND constraints.
 
std::set< TString > & multiDimConstrainedPars ()
 Access to the parameter names used in ND constraints.
 
void resetFitCounters ()
 Reset the good/bad fit counters.
 
void setCurrentExperiment (const UInt_t curExpt)
 Set the ID of the current experiment. More...
 
void startNewFit (const UInt_t nPars, const UInt_t nFreePars)
 Indicate the start of a new fit. More...
 
void eventsPerExpt (UInt_t nEvents)
 Set the number of events in the current experiment.
 
Double_t worstLogLike () const
 Access the worst log likelihood found so far.
 
void worstLogLike (const Double_t newWorstLogLike)
 Set a new value for the worst log likelihood. More...
 
void storeFitStatus (const LauAbsFitter::FitStatus &status, const TMatrixD &covMatrix)
 Store fit status information. More...
 
UInt_t nTotParams () const
 Access the total number of fit parameters.
 
UInt_t nFreeParams () const
 Access the total number of fit parameters.
 
const LauAbsFitter::FitStatusfitStatus () const
 Access the fit status information.
 
Double_t nll () const
 Access the current NLL value.
 
Double_t edm () const
 Access the current EDM value.
 
Int_t statusCode () const
 Access the fit status code.
 
const TMatrixD & covarianceMatrix () const
 Access the fit covariance matrix.
 
UInt_t numberOKFits () const
 Access the number of successful fits.
 
UInt_t numberBadFits () const
 Access the number of failed fits.
 

Private Member Functions

 LauSimFitCoordinator (const LauSimFitCoordinator &rhs)
 Copy constructor (not implemented)
 
LauSimFitCoordinatoroperator= (const LauSimFitCoordinator &rhs)
 Copy assignment operator (not implemented)
 

Private Attributes

const UInt_t nTasks_
 The number of tasks.
 
const UInt_t reqPort_
 The requested port.
 
std::vector< TMatrixD > covMatrices_
 The covariance sub-matrices for each task.
 
TMonitor * socketMonitor_
 Parallel setup monitor.
 
std::vector< TSocket * > socketTasks_
 Sockets for each of the tasks.
 
std::vector< TMessage * > messagesToTasks_
 Messages to tasks.
 
TMessage * messageFromTask_
 Message from tasks to the coordinator.
 
std::map< TString, UInt_t > parIndices_
 Map of parameter names to index in the values vector.
 
std::map< UInt_t, TString > parNames_
 Reverse map of index in the values vector to parameter names.
 
std::vector< LauParameter * > params_
 Parameters.
 
std::vector< LauAbsRValue * > conVars_
 Gaussian constraints.
 
std::vector< Double_t > parValues_
 Parameter values.
 
std::vector< std::vector< UInt_t > > taskIndices_
 Lists of indices for each task.
 
std::vector< std::vector< UInt_t > > taskFreeIndices_
 Lists of indices of free parameters for each task.
 
std::vector< Double_t * > vectorPar_
 Parameter values to send to the tasks.
 
std::vector< Double_t > vectorRes_
 Likelihood values returned from the tasks.
 
TStopwatch timer_
 The fit timer.
 
TStopwatch cumulTimer_
 The total fit timer.
 
LauFitNtuplefitNtuple_
 The fit results ntuple.
 

Additional Inherited Members

- Protected Types inherited from LauFitObject
enum  ConstraintType { ConstraintType::Formula, ConstraintType::MultDim }
 Enumeration of the different types of constraint. More...
 

Detailed Description

The coordinator process for simultaneous/combined fits.

Implementation of the JFit method described in arXiv:1409.5080 [physics.data-an].

This class acts as the interface between the task processes and the minimiser.

Definition at line 56 of file LauSimFitCoordinator.hh.

Constructor & Destructor Documentation

◆ LauSimFitCoordinator()

LauSimFitCoordinator::LauSimFitCoordinator ( UInt_t  numTasks,
UInt_t  port = 9090 
)

Constructor.

Parameters
[in]numTasksthe number of tasks processes to expect connections from
[in]portthe port on which to listen for connections from the tasks

Definition at line 52 of file LauSimFitCoordinator.cc.

Member Function Documentation

◆ cacheInputData()

Bool_t LauSimFitCoordinator::cacheInputData ( )
protected

Instruct the tasks to perform the caching.

Returns
success/failure of the caching operations

Definition at line 646 of file LauSimFitCoordinator.cc.

◆ getTotNegLogLikelihood()

Double_t LauSimFitCoordinator::getTotNegLogLikelihood ( )
virtual

Calculate the new value of the negative log likelihood.

This function has to be public since it is called from the global FCN. It should not be called otherwise!

Implements LauFitObject.

Definition at line 775 of file LauSimFitCoordinator.cc.

◆ readData()

Bool_t LauSimFitCoordinator::readData ( )
protected

Instruct the tasks to read the input data for the given experiment.

Returns
success/failure of the reading operations

Definition at line 596 of file LauSimFitCoordinator.cc.

◆ runSimFit()

void LauSimFitCoordinator::runSimFit ( const TString &  fitNtupleFileName,
const Bool_t  useAsymmErrors = kFALSE,
const Bool_t  doTwoStageFit = kFALSE 
)

Run the fit.

Parameters
[in]fitNtupleFileNamethe file to which the ntuple containing the fit results should be written
[in]useAsymmErrorsshould asymmetric errors be calculated or not
[in]doTwoStageFitshould the fit be performed in two stages or not

Definition at line 462 of file LauSimFitCoordinator.cc.

◆ setParsFromMinuit()

void LauSimFitCoordinator::setParsFromMinuit ( Double_t *  par,
Int_t  npar 
)
virtual

This function sets the parameter values from Minuit.

This function has to be public since it is called from the global FCN. It should not be called otherwise!

Parameters
[in]paran array storing the various parameter values
[in]nparthe number of free parameters

Implements LauFitObject.

Definition at line 742 of file LauSimFitCoordinator.cc.

◆ withinAsymErrorCalc() [1/3]

Bool_t LauFitObject::withinAsymErrorCalc
inline

Query whether the fit is calculating the asymmetric errors.

Returns
kTRUE if the fit is calculating the asymmetric errors, kFALSE otherwise

Definition at line 94 of file LauFitObject.hh.

◆ withinAsymErrorCalc() [2/3]

void LauSimFitCoordinator::withinAsymErrorCalc ( const Bool_t  inAsymErrCalc)
virtual

Mark that the fit is calculating asymmetric errors.

This function has to be public since it is called by the fitter interface to mark when entering and exiting the asymmetric error calculation. It should not be called otherwise!

Parameters
[in]inAsymErrCalcboolean marking that the fit is calculating the asymmetric errors

Reimplemented from LauFitObject.

Definition at line 551 of file LauSimFitCoordinator.cc.

◆ withinAsymErrorCalc() [3/3]

virtual void LauFitObject::withinAsymErrorCalc
inline

Mark that the fit is calculating asymmetric errors.

This is called by the fitter interface to mark when entering and exiting the asymmetric error calculation.

Parameters
[in]inAsymErrCalcboolean marking that the fit is calculating the asymmetric errors

Definition at line 85 of file LauFitObject.hh.


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