laura is hosted by Hepforge, IPPP Durham
Laura++  v3r0p1
A maximum likelihood fitting package for performing Dalitz-plot analysis.
LauIsobarDynamics.hh
Go to the documentation of this file.
1 
2 // Copyright University of Warwick 2005 - 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 
19 #ifndef LAU_ISOBAR_DYNAMICS
20 #define LAU_ISOBAR_DYNAMICS
21 
22 #include <set>
23 #include <vector>
24 
25 #include "TString.h"
26 
27 #include "LauAbsResonance.hh"
28 #include "LauComplex.hh"
29 
30 class LauCacheData;
31 class LauDaughters;
32 class LauAbsEffModel;
33 class LauAbsIncohRes;
34 class LauFitDataTree;
37 class LauKinematics;
38 
40 
41  public:
43  typedef std::map<Int_t,LauAbsEffModel*> LauTagCatScfFractionModelMap;
44 
46  enum ToyMCStatus {
50  };
51 
53 
58  LauIsobarDynamics(LauDaughters* daughters, LauAbsEffModel* effModel, LauAbsEffModel* scfFractionModel = 0);
59 
61 
66  LauIsobarDynamics(LauDaughters* daughters, LauAbsEffModel* effModel, LauTagCatScfFractionModelMap scfFractionModel);
67 
69  virtual ~LauIsobarDynamics();
70 
72 
75  void initialise(const std::vector<LauComplex>& coeffs);
76 
79 
81 
84  inline void setIntFileName(const TString& fileName) {intFileName_ = fileName;}
85 
86  // Integration
88 
92  void setIntegralBinWidths(const Double_t m13BinWidth, const Double_t m23BinWidth);
93 
95 
101  void setNarrowResonanceThreshold(const Double_t narrowWidth) { narrowWidth_ = narrowWidth; }
102 
104 
110  void setIntegralBinningFactor(const Double_t binningFactor) { binningFactor_ = binningFactor; }
111 
113 
122 
124 
131  LauAbsResonance* addIncoherentResonance(const TString& resName, const Int_t resPairAmpInt, const LauAbsResonance::LauResonanceModel resType);
132 
134 
142  void defineKMatrixPropagator(const TString& propName, const TString& paramFileName,
143  Int_t resPairAmpInt, Int_t nChannels, Int_t nPoles, Int_t rowIndex = 1);
144 
146 
151  void addKMatrixProdPole(const TString& poleName, const TString& propName, Int_t poleIndex);
152 
154 
159  void addKMatrixProdSVP(const TString& SVPName, const TString& propName, Int_t channelIndex);
160 
162 
165  inline void setASqMaxValue(Double_t value) {aSqMaxSet_ = value;}
166 
168 
171  inline Double_t getASqMaxSetValue() const { return aSqMaxSet_; }
172 
174 
177  inline Double_t getASqMaxVarValue() const { return aSqMaxVar_; }
178 
180 
183  Bool_t generate();
184 
186 
191  ToyMCStatus checkToyMC(Bool_t printErrorMessages = kTRUE, Bool_t printInfoMessages = kFALSE);
192 
194 
197  inline Int_t maxGenIterations() const {return iterationsMax_;}
198 
200 
203  void calcLikelihoodInfo(const UInt_t iEvt);
204 
206 
210  void calcLikelihoodInfo(const Double_t m13Sq, const Double_t m23Sq);
211 
213 
219  void calcLikelihoodInfo(const Double_t m13Sq, const Double_t m23Sq, const Int_t tagCat);
220 
222 
225  void calcExtraInfo(const Bool_t init = kFALSE);
226 
228 
232  Bool_t gotReweightedEvent();
233 
235 
238  Double_t getEventWeight();
239 
241 
244  inline const LauComplex& getEvtDPAmp() const {return totAmp_;}
245 
247 
250  inline Double_t getEvtm13Sq() const {return m13Sq_;}
251 
253 
256  inline Double_t getEvtm23Sq() const {return m23Sq_;}
257 
259 
262  inline Double_t getEvtmPrime() const {return mPrime_;}
263 
265 
268  inline Double_t getEvtthPrime() const {return thPrime_;}
269 
271 
274  inline Double_t getEvtEff() const {return eff_;}
275 
277 
280  inline Double_t getEvtScfFraction() const {return scfFraction_;}
281 
283 
286  inline Double_t getEvtJacobian() const {return jacobian_;}
287 
289 
294  inline Double_t getEvtLikelihood() const {return evtLike_;}
295 
297 
301  inline LauComplex getDynamicAmp(Int_t resID) const {return ff_[resID].scale(fNorm_[resID]);}
302 
304 
307  inline const std::vector< std::vector<LauComplex> >& getFiFjSum() const {return fifjSum_;}
308 
310 
313  inline const std::vector< std::vector<LauComplex> >& getFiFjEffSum() const {return fifjEffSum_;}
314 
316 
319  inline const std::vector<Double_t>& getFNorm() const {return fNorm_;}
320 
322 
325  void fillDataTree(const LauFitDataTree& fitDataTree);
326 
328  void modifyDataTree();
329 
331 
335  Bool_t hasResonance(const TString& resName) const;
336 
338 
342  TString getConjResName(const TString& resName) const;
343 
345 
348  void updateCoeffs(const std::vector<LauComplex>& coeffs);
349 
351 
354  inline void flipHelicityForCPEigenstates(Bool_t boolean) {flipHelicity_ = boolean;}
355 
357 
360  inline LauParameter getMeanEff() const {return meanDPEff_;}
361 
363 
366  inline LauParameter getDPRate() const {return DPRate_;}
367 
369 
372  inline const LauParArray& getFitFractions() const {return fitFrac_;}
373 
375 
379 
381 
384  inline UInt_t getnTotAmp() const {return nAmp_+nIncohAmp_;}
385 
387 
390  inline UInt_t getnCohAmp() const {return nAmp_;}
391 
393 
396  inline UInt_t getnIncohAmp() const {return nIncohAmp_;}
397 
399 
402  inline Double_t getDPNorm() const {return DPNorm_;}
403 
405 
408  inline UInt_t nData() const {return data_.size();}
409 
411 
414  inline const std::vector<LauCacheData*>& getCacheData() const {return data_;}
415 
417 
421 
423 
427 
429 
433 
435 
439 
441 
444  inline std::map <Int_t,LauAbsEffModel*> getScfFractionModels() {return scfFractionModel_;}
445 
447 
450  inline Bool_t usingScfModel() { return ! scfFractionModel_.empty(); }
451 
453 
456  inline std::vector<LauParameter> getExtraParameters() {return extraParameters_;}
457 
459 
462  inline std::vector<LauParameter*>& getFloatingParameters() {return resonancePars_;}
463 
464  protected:
466  void initSummary();
467 
469  void initialiseVectors();
470 
472  void resetNormVectors();
473 
475  void calcDPNormalisation();
476 
479 
482 
484 
488 
490  void writeIntegralsFile();
491 
493 
498  void setFFTerm(const UInt_t index, const Double_t realPart, const Double_t imagPart);
499 
501 
505  void setIncohIntenTerm(const UInt_t index, const Double_t value);
506 
508  void calculateAmplitudes();
509 
511 
516  void calculateAmplitudes( LauDPPartialIntegralInfo* intInfo, const UInt_t m13Point, const UInt_t m23Point );
517 
519 
522  void addGridPointToIntegrals(const Double_t weight);
523 
525 
528  void calcTotalAmp(const Bool_t useEff);
529 
531 
534  Double_t retrieveEfficiency();
535 
537 
541  Double_t retrieveScfFraction(Int_t tagCat);
542 
544 
547  inline void setASqMaxVarValue(Double_t value) {aSqMaxVar_ = value;}
548 
550 
553  Double_t calcSigDPNorm();
554 
556 
559  void resAmp(const UInt_t index);
560 
562 
565  void incohResAmp(const UInt_t index);
566 
568 
571  void setDataEventNo(UInt_t iEvt);
572 
574 
578  LauAbsResonance* findResonance(const TString& name);
579 
581 
585  const LauAbsResonance* findResonance(const TString& name) const;
586 
588 
591  void removeCharge(TString& string) const;
592 
594 
599  Bool_t gotKMatrixMatch(UInt_t resAmpInt, const TString& propName) const;
600 
601  private:
604 
607 
609  typedef std::map<TString, LauKMatrixPropagator*> KMPropMap;
610 
612  typedef std::map<TString, TString> KMStringMap;
613 
616 
619 
622 
624 
629 
631  UInt_t nAmp_;
632 
634  UInt_t nIncohAmp_;
635 
637  std::vector<LauComplex> Amp_;
638 
640  Double_t DPNorm_;
641 
644 
647 
650 
653 
655  std::vector<LauCacheData*> data_;
656 
659 
661  std::vector<LauParameter> extraParameters_;
662 
664  std::vector<LauAbsResonance*> sigResonances_;
665 
667  std::vector<LauAbsIncohRes*> sigIncohResonances_;
668 
671 
674 
676  std::vector<TString> resTypAmp_;
677 
679  std::vector<Int_t> resPairAmp_;
680 
682  std::vector<TString> incohResTypAmp_;
683 
685  std::vector<Int_t> incohResPairAmp_;
686 
688  std::vector<Int_t> typDaug_;
689 
692 
695 
698 
701 
703  std::vector<LauDPPartialIntegralInfo*> dpPartialIntegralInfo_;
704 
706  TString intFileName_;
707 
709  Double_t m13BinWidth_;
710 
712  Double_t m23BinWidth_;
713 
715  Double_t narrowWidth_;
716 
718  Double_t binningFactor_;
719 
721  Double_t m13Sq_;
722 
724  Double_t m23Sq_;
725 
727  Double_t mPrime_;
728 
730  Double_t thPrime_;
731 
733  Int_t tagCat_;
734 
736  Double_t eff_;
737 
739  Double_t scfFraction_;
740 
742  Double_t jacobian_;
743 
745  Double_t ASq_;
746 
748  Double_t evtLike_;
749 
752 
754 
757  std::vector< std::vector<LauComplex> > fifjEffSum_;
758 
760 
763  std::vector< std::vector<LauComplex> > fifjSum_;
764 
766  std::vector<LauComplex> ff_;
767 
769  std::vector<Double_t> incohInten_;
770 
772  std::vector<Double_t> fSqSum_;
773 
775  std::vector<Double_t> fSqEffSum_;
776 
778  std::vector<Double_t> fNorm_;
779 
782 
785 
787  Double_t aSqMaxSet_;
788 
790  Double_t aSqMaxVar_;
791 
794 
797 
799  std::vector<LauParameter*> resonancePars_;
800 
802  std::vector<Double_t> resonanceParValues_;
803 
805  std::vector< std::vector<UInt_t> > resonanceParResIndex_;
806 
808  std::set<UInt_t> integralsToBeCalculated_;
809 
810  ClassDef(LauIsobarDynamics,0)
811 };
812 
813 #endif
KMStringMap kMatrixPropSet_
The names of the M-matrix components in the model mapped to their propagators.
std::map< Int_t, LauAbsEffModel * > getScfFractionModels()
Retrieve the model for the fraction of events that are poorly reconstructed (the self cross feed frac...
std::vector< LauComplex > Amp_
The complex coefficients for the amplitude components.
void modifyDataTree()
Recache the amplitude values for those that have changed.
Int_t nSigGenLoop_
The number of unsucessful attempts to generate an event so far.
Double_t getASqMaxVarValue() const
Retrieve the maximum of A squared that has been found while generating.
UInt_t nData() const
Retrieve the number of cached events.
Bool_t recalcNormalisation_
Flag to recalculate the normalisation.
virtual ~LauIsobarDynamics()
Destructor.
const std::vector< Double_t > & getFNorm() const
Retrieve the normalisation factors for the dynamic parts of the amplitudes for all of the amplitude c...
void updateCoeffs(const std::vector< LauComplex > &coeffs)
Update the complex coefficients for the resonances.
std::set< UInt_t > integralsToBeCalculated_
Resonance indices for which the amplitudes and integrals should be recalculated.
Double_t m13BinWidth_
The bin width to use when integrating over m13.
Bool_t flipHelicity_
The helicity flip flag for new amplitude components.
std::vector< Int_t > resPairAmp_
The index of the daughter not produced by the resonance for each amplitude component.
Bool_t hasResonance(const TString &resName) const
Check whether this model includes a named resonance.
std::vector< LauParameter > extraParameters_
any extra parameters/quantities (e.g. K-matrix total fit fractions)
Class that defines the particular 3-body decay under study.
Definition: LauDaughters.hh:33
Double_t getEvtJacobian() const
Retrieve the Jacobian, for the transformation into square DP coordinates, for the current event...
void defineKMatrixPropagator(const TString &propName, const TString &paramFileName, Int_t resPairAmpInt, Int_t nChannels, Int_t nPoles, Int_t rowIndex=1)
Define a new K-matrix Propagator.
Double_t getDPNorm() const
Retrieve the normalisation factor for the log-likelihood function.
std::vector< Double_t > fSqEffSum_
The event-by-event running total of the dynamical amplitude squared for each amplitude component...
void calcDPNormalisationScheme()
Calculate the Dalitz plot normalisation integrals across the whole Dalitz plot.
LauAbsEffModel * effModel_
The efficiency model across the Dalitz plot.
void calcExtraInfo(const Bool_t init=kFALSE)
Calculate the fit fractions, mean efficiency and total DP rate.
std::vector< Int_t > incohResPairAmp_
The index of the daughter not produced by the resonance for each incoherent amplitude component...
void setASqMaxValue(Double_t value)
Set the maximum value of A squared to be used in the accept/reject.
Double_t DPNorm_
The normalisation factor for the log-likelihood function.
LauTagCatScfFractionModelMap scfFractionModel_
The self cross feed fraction models across the Dalitz plot.
LauParArray fitFrac_
The fit fractions for the amplitude components.
Double_t aSqMaxVar_
The maximum value of A squared that has been seen so far while generating.
Double_t scfFraction_
The fraction of events that are poorly reconstructed (the self cross feed fraction) at the current po...
Abstract class for defining incoherent resonant amplitude models.
Double_t m13Sq_
The invariant mass squared of the first and third daughters.
Int_t iterationsMax_
The maximum allowed number of attempts when generating an event.
LauComplex totAmp_
The total amplitude for the current event.
Double_t getEvtthPrime() const
Retrieve the square Dalitz plot coordinate, theta&#39;, for the current event.
Double_t getASqMaxSetValue() const
Retrieve the maximum value of A squared to be used in the accept/reject.
TString getConjResName(const TString &resName) const
Retrieve the name of the charge conjugate of a named resonance.
std::vector< std::vector< LauParameter > > LauParArray
Type to define an array of parameters.
LauAbsEffModel * getEffModel()
Retrieve the model for the efficiency across the Dalitz plot.
Double_t aSqMaxSet_
The maximum allowed value of A squared.
std::vector< Double_t > resonanceParValues_
List of floating resonance parameter values from previous calculation.
Bool_t symmetricalDP_
Whether the Dalitz plot is symmetrical.
UInt_t getnIncohAmp() const
Retrieve the number of incoherent amplitude components.
Pure abstract base class for defining the efficiency description across the signal Dalitz plot...
void setIncohIntenTerm(const UInt_t index, const Double_t value)
Set the dynamic part of the intensity for a given incoherent amplitude component at the current point...
Class for defining (a section of) the Dalitz plot integration binning scheme.
Class to contain cached data relating to an event.
Definition: LauCacheData.hh:30
void calcDPNormalisation()
Calculate the Dalitz plot normalisation integrals across the whole Dalitz plot.
const LauParArray & getFitFractions() const
Retrieve the fit fractions for the amplitude components.
void recalculateNormalisation()
recalculate Normalization
std::vector< TString > incohResTypAmp_
The resonance types of all of the incoherent amplitude components.
std::vector< LauAbsIncohRes * > sigIncohResonances_
The incoherent resonances in the model.
void setASqMaxVarValue(Double_t value)
Set the maximum of A squared that has been found.
Double_t getEvtScfFraction() const
Retrieve the fraction of events that are poorly reconstructed (the self cross feed fraction) for the ...
void resAmp(const UInt_t index)
Calculate the dynamic part of the amplitude for a given component at the current point in the Dalitz ...
LauDaughters * getDaughters()
Retrieve the daughters.
LauParameter DPRate_
The overall Dalitz plot rate.
void calcLikelihoodInfo(const UInt_t iEvt)
Calculate the likelihood (and all associated information) for the given event number.
ToyMCStatus
The possible statuses for toy MC generation.
LauParameter getMeanEff() const
Retrieve the mean efficiency across the Dalitz plot.
KMPropMap kMatrixPropagators_
The K-matrix propagators.
std::vector< Double_t > fSqSum_
The event-by-event running total of the dynamical amplitude squared for each amplitude component...
Double_t getEvtmPrime() const
Retrieve the square Dalitz plot coordinate, m&#39;, for the current event.
UInt_t nAmp_
The number of amplitude components.
Double_t m23Sq_
The invariant mass squared of the second and third daughters.
Double_t thPrime_
The square Dalitz plot coordinate theta&#39;.
Double_t evtLike_
The normalised likelihood for the current event.
void setIntegralBinningFactor(const Double_t binningFactor)
Set the factor relating the width of a narrow resonance and the binning size in its integration regio...
std::vector< std::vector< LauComplex > > fifjEffSum_
The event-by-event running total of efficiency corrected amplitude cross terms for each pair of ampli...
void initSummary()
Print a summary of the model to be used.
Bool_t fullySymmetricDP_
Whether the Dalitz plot is fully symmetric.
void fillDataTree(const LauFitDataTree &fitDataTree)
Fill the internal data structure that caches the resonance dynamics.
LauAbsResonance * findResonance(const TString &name)
Retrieve the named resonance.
void resetNormVectors()
Zero the various values used to store integrals info.
Double_t ASq_
The value of A squared for the current event.
const std::vector< std::vector< LauComplex > > & getFiFjSum() const
Retrieve the event-by-event running totals of amplitude cross terms for all pairs of amplitude compon...
Double_t getEvtm23Sq() const
Retrieve the invariant mass squared of the second and third daughters in the current event...
Double_t m23BinWidth_
The bin width to use when integrating over m23.
TString intFileName_
The name of the file to save integrals to.
std::vector< LauParameter * > & getFloatingParameters()
Retrieve the floating parameters of the resonance models.
void setNarrowResonanceThreshold(const Double_t narrowWidth)
Set the value below which a resonance width is considered to be narrow.
void addKMatrixProdPole(const TString &poleName, const TString &propName, Int_t poleIndex)
Add a K-matrix production pole term to the model.
void calcTotalAmp(const Bool_t useEff)
Calculate the total Dalitz plot amplitude at the current point in the Dalitz plot.
const LauComplex & getEvtDPAmp() const
Retrieve the total amplitude for the current event.
Double_t retrieveScfFraction(Int_t tagCat)
Obtain the self cross feed fraction of the current event from the model.
std::vector< TString > resTypAmp_
The resonance types of all of the amplitude components.
LauKinematics * getKinematics()
Retrieve the Dalitz plot kinematics.
LauIsobarDynamics(LauDaughters *daughters, LauAbsEffModel *effModel, LauAbsEffModel *scfFractionModel=0)
Constructor.
Int_t tagCat_
The tagging category.
LauCacheData * currentEvent_
The cached data for the current event.
void initialiseVectors()
Initialise the internal storage for this model.
Double_t jacobian_
The Jacobian, for the transformation into square DP coordinates at the current point in the Dalitz pl...
LauDaughters * daughters_
The daughters of the decay.
std::map< Int_t, LauAbsEffModel * > LauTagCatScfFractionModelMap
The type used for containing multiple self cross feed fraction models for different categories (e...
LauParameter getDPRate() const
Retrieve the overall Dalitz plot rate.
Double_t narrowWidth_
The value below which a resonance width is considered to be narrow.
File containing declaration of LauComplex class.
void incohResAmp(const UInt_t index)
Calculate the dynamic part of the intensity for a given incoherent component at the current point in ...
const LauParArray & getFitFractionsEfficiencyUncorrected() const
Retrieve the fit fractions for the amplitude components.
std::vector< std::vector< LauComplex > > fifjSum_
The event-by-event running total of the amplitude cross terms for each pair of amplitude components...
Class for defining the fit parameter objects.
Definition: LauParameter.hh:34
UInt_t getnTotAmp() const
Retrieve the total number of amplitude components.
Double_t calcSigDPNorm()
Calculate the normalisation factor for the log-likelihood function.
std::vector< LauParameter > getExtraParameters()
Retrieve any extra parameters/quantities (e.g. K-matrix total fit fractions)
void addGridPointToIntegrals(const Double_t weight)
Add the amplitude values (with the appropriate weight) at the current grid point to the running integ...
Bool_t generate()
Generate a toy MC signal event.
LauAbsResonance * addResonance(const TString &resName, const Int_t resPairAmpInt, const LauAbsResonance::LauResonanceModel resType, const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory=LauBlattWeisskopfFactor::Default, const LauBlattWeisskopfFactor::BarrierType bwType=LauBlattWeisskopfFactor::BWPrimeBarrier)
Add a resonance to the Dalitz plot.
void removeCharge(TString &string) const
Remove the charge from the given particle name.
std::vector< Int_t > typDaug_
The PDG codes of the daughters.
LauResonanceModel
Define the allowed resonance types.
void flipHelicityForCPEigenstates(Bool_t boolean)
Set the helicity flip flag for new amplitude components.
Double_t eff_
The efficiency at the current point in the Dalitz plot.
void calculateAmplitudes()
Calculate the amplitudes for all resonances for the current kinematics.
BarrierType
Define the allowed types of barrier factors.
Int_t maxGenIterations() const
Retrieve the maximum number of iterations allowed when generating an event.
Bool_t gotReweightedEvent()
Calculates whether an event with the current kinematics should be accepted in order to produce a dist...
Abstract class for defining type for resonance amplitude models (Breit-Wigner, Flatte etc...
Double_t binningFactor_
The factor relating the width of the narrowest resonance and the binning size.
void addKMatrixProdSVP(const TString &SVPName, const TString &propName, Int_t channelIndex)
Add a K-matrix slowly-varying part (SVP) term to the model.
std::vector< LauCacheData * > data_
The cached data for all events.
void findIntegralsToBeRecalculated()
Determine which amplitudes and integrals need to be recalculated.
Class for defining signal dynamics using the isobar model.
Double_t getEventWeight()
Calculate the acceptance rate, for events with the current kinematics, when generating events accordi...
void setIntFileName(const TString &fileName)
Set the name of the file to which to save the results of the integrals.
LauAbsResonance * addIncoherentResonance(const TString &resName, const Int_t resPairAmpInt, const LauAbsResonance::LauResonanceModel resType)
Add an incoherent resonance to the Dalitz plot.
void writeIntegralsFile()
Write the results of the integrals (and related information) to a file.
File containing declaration of LauAbsResonance class.
void calcDPPartialIntegral(LauDPPartialIntegralInfo *intInfo)
Calculate the Dalitz plot normalisation integrals over a given range.
ToyMCStatus checkToyMC(Bool_t printErrorMessages=kTRUE, Bool_t printInfoMessages=kFALSE)
Check the status of the toy MC generation.
LauKinematics * kinematics_
The kinematics of the decay.
UInt_t getnCohAmp() const
Retrieve the number of coherent amplitude components.
LauIsobarDynamics & operator=(const LauIsobarDynamics &rhs)
Copy assignment operator (not implemented)
const std::vector< LauCacheData * > & getCacheData() const
Retrieve the cached data.
Class for defining a complex number.
Definition: LauComplex.hh:47
Double_t getEvtEff() const
Retrieve the efficiency for the current event.
UInt_t nIncohAmp_
The number of incoherent amplitude components.
Bool_t gotKMatrixMatch(UInt_t resAmpInt, const TString &propName) const
Check whether a resonance is a K-matrix component of a given propagator.
Double_t getEvtLikelihood() const
Retrieve the likelihood for the current event.
std::vector< LauDPPartialIntegralInfo * > dpPartialIntegralInfo_
The storage of the integration scheme.
std::map< TString, TString > KMStringMap
The type used for mapping K-matrix components to their propagators.
Bool_t usingScfModel()
Check whether a self cross feed fraction model is being used.
void setDataEventNo(UInt_t iEvt)
Load the data for a given event.
std::map< TString, LauKMatrixPropagator * > KMPropMap
The type used for containing the K-matrix propagators.
Class for calculating 3-body kinematic quantities.
void setIntegralBinWidths(const Double_t m13BinWidth, const Double_t m23BinWidth)
Set the widths of the bins to use when integrating across the Dalitz plot.
Bool_t normalizationSchemeDone_
Whether the scheme for the integration has been determined.
std::vector< LauParameter * > resonancePars_
List of floating resonance parameters.
std::vector< std::vector< UInt_t > > resonanceParResIndex_
Indices in sigResonances_ to point to the corresponding signal resonance(s) for each floating paramet...
BlattWeisskopfCategory
Define resonance categories that will share common barrier factor radii.
LauAbsEffModel * getScfFractionModel()
Retrieve the model for the fraction of events that are poorly reconstructed (the self cross feed frac...
LauParameter meanDPEff_
The mean efficiency across the Dalitz plot.
Bool_t integralsDone_
Whether the integrals have been performed.
std::vector< Double_t > incohInten_
The dynamic part of the intensity for each incoherent amplitude component at the current point in the...
LauComplex getDynamicAmp(Int_t resID) const
Retrieve the normalised dynamic part of the amplitude of the given amplitude component at the current...
std::vector< LauComplex > ff_
The dynamic part of the amplitude for each amplitude component at the current point in the Dalitz plo...
void initialise(const std::vector< LauComplex > &coeffs)
Initialise the Dalitz plot dynamics.
const std::vector< std::vector< LauComplex > > & getFiFjEffSum() const
Retrieve the event-by-event running totals of efficiency corrected amplitude cross terms for all pair...
Double_t getEvtm13Sq() const
Retrieve the invariant mass squared of the first and third daughters in the current event...
Double_t mPrime_
The square Dalitz plot coordinate, m&#39;.
std::vector< Double_t > fNorm_
The normalisation factors for the dynamic parts of the amplitude for each amplitude component...
Class to store the input fit variables.
std::vector< LauAbsResonance * > sigResonances_
The resonances in the model.
void setFFTerm(const UInt_t index, const Double_t realPart, const Double_t imagPart)
Set the dynamic part of the amplitude for a given amplitude component at the current point in the Dal...
Double_t retrieveEfficiency()
Obtain the efficiency of the current event from the model.
LauParArray fitFracEffUnCorr_
The efficiency-uncorrected fit fractions for the amplitude components.
Class for defining a K-matrix propagator.