laura is hosted by Hepforge, IPPP Durham
Laura++  v3r3
A maximum likelihood fitting package for performing Dalitz-plot analysis.
LauKinematics.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 
23 #ifndef LAU_KINEMATICS
24 #define LAU_KINEMATICS
25 
26 #include <vector>
27 
28 #include "TMath.h"
29 
30 
32 
33  public:
35 
44  LauKinematics(const Double_t m1, const Double_t m2, const Double_t m3, const Double_t mParent,
45  const Bool_t calcSquareDPCoords = kFALSE, const Bool_t symmetricalDP = kFALSE, const Bool_t fullySymmetricDP = kFALSE);
46 
48  virtual ~LauKinematics();
49 
51 
54  inline void squareDP( const Bool_t calcSquareDPCoords ) { squareDP_ = calcSquareDPCoords; }
55 
57 
60  inline Bool_t squareDP() const { return squareDP_; }
61 
63 
66  inline Bool_t gotSymmetricalDP() const { return (symmetricalDP_ && ! fullySymmetricDP_); }
67 
69 
72  inline Bool_t gotFullySymmetricDP() const { return fullySymmetricDP_; }
73 
75  inline void warningMessages(const Bool_t boolean) { warnings_ = boolean; }
76 
78 
84  void updateKinematics(const Double_t m13Sq, const Double_t m23Sq);
85 
87 
93  void updateSqDPKinematics(const Double_t mPrime, const Double_t thetaPrime);
94 
96 
100  void updateKinematicsFrom23(const Double_t m23, const Double_t c23);
101 
103 
107  void updateKinematicsFrom13(const Double_t m13, const Double_t c13);
108 
110 
114  void updateKinematicsFrom12(const Double_t m12, const Double_t c12);
115 
117 
122  Double_t calcSqDPJacobian(const Double_t mPrime, const Double_t thPrime) const;
123 
125 
128  Double_t calcSqDPJacobian() const;
129 
131 
135  void genFlatPhaseSpace(Double_t& m13Sq, Double_t& m23Sq) const;
136 
138 
142  void genFlatSqDP(Double_t& mPrime, Double_t& thetaPrime) const;
143 
145 
155  Bool_t withinDPLimits(const Double_t m13Sq, const Double_t m23Sq) const;
156 
158 
168  Bool_t withinDPLimits2(const Double_t m13Sq, const Double_t m23Sq) const;
169 
171 
176  Bool_t withinSqDPLimits(const Double_t mPrime, const Double_t thetaPrime) const;
177 
179 
184  Double_t calcThirdMassSq(const Double_t firstMassSq, const Double_t secondMassSq) const;
185 
187 
190  Double_t distanceFromDPCentre() const;
191 
193 
196  Double_t distanceFromDPCentre(const Double_t m13Sq, const Double_t m23Sq) const;
197 
199 
202  inline Double_t getm12() const {return m12_;}
204 
207  inline Double_t getm23() const {return m23_;}
209 
212  inline Double_t getm13() const {return m13_;}
213 
215 
218  inline Double_t getm12Sq() const {return m12Sq_;}
220 
223  inline Double_t getm23Sq() const {return m23Sq_;}
225 
228  inline Double_t getm13Sq() const {return m13Sq_;}
229 
231 
235  inline Double_t getc12() const {return c12_;}
237 
241  inline Double_t getc23() const {return c23_;}
243 
247  inline Double_t getc13() const {return c13_;}
248 
250 
253  inline Double_t getmPrime() const {return mPrime_;}
255 
258  inline Double_t getThetaPrime() const {return thetaPrime_;}
259 
261 
264  inline Double_t getmParent() const {return mParent_;}
266 
269  inline Double_t getmParentSq() const {return mParentSq_;}
270 
272 
278  inline Double_t getDPBoxArea() const {return (mSqMax_[1] - mSqMin_[1])*(mSqMax_[0] - mSqMin_[0]);}
279 
281 
285 
287 
291 
293 
296  inline Double_t getm1() const {return m1_;}
298 
301  inline Double_t getm2() const {return m2_;}
303 
306  inline Double_t getm3() const {return m3_;}
307 
309 
312  inline Double_t getm23Min() const {return TMath::Sqrt(mSqMin_[0]);}
314 
317  inline Double_t getm13Min() const {return TMath::Sqrt(mSqMin_[1]);}
319 
322  inline Double_t getm12Min() const {return TMath::Sqrt(mSqMin_[2]);}
323 
325 
328  inline Double_t getm23Max() const {return TMath::Sqrt(mSqMax_[0]);}
330 
333  inline Double_t getm13Max() const {return TMath::Sqrt(mSqMax_[1]);}
335 
338  inline Double_t getm12Max() const {return TMath::Sqrt(mSqMax_[2]);}
339 
341 
344  inline Double_t getm23SqMin() const {return mSqMin_[0];}
346 
349  inline Double_t getm13SqMin() const {return mSqMin_[1];}
351 
354  inline Double_t getm12SqMin() const {return mSqMin_[2];}
355 
357 
360  inline Double_t getm23SqMax() const {return mSqMax_[0];}
362 
365  inline Double_t getm13SqMax() const {return mSqMax_[1];}
367 
370  inline Double_t getm12SqMax() const {return mSqMax_[2];}
371 
373 
376  inline Double_t getp1_12() const {return p1_12_;}
378 
381  inline Double_t getp3_12() const {return p3_12_;}
383 
386  inline Double_t getp2_23() const {return p2_23_;}
388 
391  inline Double_t getp1_23() const {return p1_23_;}
393 
396  inline Double_t getp1_13() const {return p1_13_;}
398 
401  inline Double_t getp2_13() const {return p2_13_;}
402 
404 
407  inline Double_t getp1_Parent() const {return p1_Parent_;}
409 
412  inline Double_t getp2_Parent() const {return p2_Parent_;}
414 
417  inline Double_t getp3_Parent() const {return p3_Parent_;}
418 
420 
424  void drawDPContour(const Int_t orientation = 1323, const Int_t nbins = 100);
425 
427 
430  inline Double_t getcov12() const {return (mParentSq_ + m12Sq_ - m3Sq_)/(2.*mParent_*m12_);}
431 
433 
436  inline Double_t getcov13() const {return (mParentSq_ + m13Sq_ - m2Sq_)/(2.*mParent_*m13_);}
437 
439 
442  inline Double_t getcov23() const {return (mParentSq_ + m23Sq_ - m1Sq_)/(2.*mParent_*m23_);}
443 
444  protected:
446 
450  void updateMassSq_m12(const Double_t m12, const Double_t c12);
451 
453 
457  void updateMassSq_m23(const Double_t m23, const Double_t c23);
458 
460 
464  void updateMassSq_m13(const Double_t m13, const Double_t c13);
465 
467 
473  void updateMassSquares(const Double_t m13Sq, const Double_t m23Sq);
474 
476 
482  void updateSqDPMassSquares(const Double_t mPrime, const Double_t thetaPrime);
483 
485 
494  Double_t cFromM(const Double_t mijSq, const Double_t mikSq, const Double_t mij, const Int_t i, const Int_t j, const Int_t k) const;
495 
497 
506  Double_t mFromC(const Double_t mijSq, const Double_t cij, const Double_t mij, const Int_t i, const Int_t j, const Int_t k) const;
507 
509 
514  Double_t pCalc(const Double_t energy, const Double_t massSq) const;
515 
517 
520  Double_t genm13Sq() const;
522 
525  Double_t genm23Sq() const;
527 
530  Double_t genm12Sq() const;
531 
533  void calcm12Sq();
534 
536  void calcHelicities();
537 
539  void calcSqDPVars();
540 
542  void calcParentFrameMomenta();
543 
544  private:
546  LauKinematics(const LauKinematics& rhs);
547 
550 
552  const Bool_t symmetricalDP_;
554  const Bool_t fullySymmetricDP_;
555 
557  const Double_t m1_;
559  const Double_t m2_;
561  const Double_t m3_;
563  const Double_t mParent_;
564 
566  const Double_t m1Sq_;
568  const Double_t m2Sq_;
570  const Double_t m3Sq_;
572  const Double_t mParentSq_;
573 
575  std::vector<Double_t> mass_;
577  std::vector<Double_t> mMin_;
579  std::vector<Double_t> mMax_;
581  std::vector<Double_t> mDiff_;
582 
584  std::vector<Double_t> mSq_;
586  std::vector<Double_t> mSqMin_;
588  std::vector<Double_t> mSqMax_;
590  std::vector<Double_t> mSqDiff_;
591 
593  const Double_t mDTot_;
595  const Double_t massInt_;
597  const Double_t mSqDTot_;
598 
600  Double_t m12_;
602  Double_t m23_;
604  Double_t m13_;
605 
607  Double_t m12Sq_;
609  Double_t m23Sq_;
611  Double_t m13Sq_;
612 
614  Double_t c12_;
616  Double_t c23_;
618  Double_t c13_;
619 
621  Double_t mPrime_;
623  Double_t thetaPrime_;
624 
626  mutable Double_t qi_;
628  mutable Double_t qk_;
629 
631  Double_t p1_12_;
633  Double_t p3_12_;
635  Double_t p2_23_;
637  Double_t p1_23_;
639  Double_t p1_13_;
641  Double_t p2_13_;
642 
644  Double_t p1_Parent_;
646  Double_t p2_Parent_;
648  Double_t p3_Parent_;
649 
651  Bool_t squareDP_;
653  Bool_t warnings_;
654 
655  ClassDef(LauKinematics,0)
656 };
657 
658 Double_t dal(Double_t* x, Double_t* par);
659 
660 #endif
Bool_t withinDPLimits(const Double_t m13Sq, const Double_t m23Sq) const
Check whether a given (m13Sq,m23Sq) point is within the kinematic limits of the Dalitz plot...
void updateMassSq_m12(const Double_t m12, const Double_t c12)
Update the variables m23Sq_ and m13Sq_ given the invariant mass m12 and the cosine of the helicity an...
Double_t p1_23_
Momentum of track 1 in 2-3 rest frame.
Double_t getm23Max() const
Get the m23 maximum defined as (mParent - m1)
Double_t getm13SqMin() const
Get the m13Sq minimum, (m1 + m3)^2.
const Double_t mParentSq_
Mass of parent particle squared.
Bool_t squareDP() const
Are the square Dalitz plot co-ordinates being calculated?
Bool_t gotSymmetricalDP() const
Is the DP symmetric?
Double_t getc23() const
Get the cosine of the helicity angle theta23.
Double_t getm23SqMin() const
Get the m23Sq minimum, (m2 + m3)^2.
void updateMassSq_m13(const Double_t m13, const Double_t c13)
Update the variables m12Sq_ and m23Sq_ given the invariant mass m13 and the cosine of the helicity an...
Double_t qi_
Momentum q of particle i.
void warningMessages(const Bool_t boolean)
Enable/disable warning messages.
Double_t getp2_Parent() const
Get the momentum of the track 2 in parent rest frame.
std::vector< Double_t > mSqMin_
Vector of the minimum mijSq values.
LauKinematics & operator=(const LauKinematics &rhs)
Copy assignment operator (not implemented)
Double_t getcov12() const
Get covariant factor in 12 axis.
void squareDP(const Bool_t calcSquareDPCoords)
Enable/disable the calculation of square Dalitz plot co-ordinates.
Double_t getm23() const
Get the m23 invariant mass.
Double_t m23_
Invariant mass m23.
Double_t getcov23() const
Get covariant factor in 23 axis.
Double_t p2_Parent_
Momentum of track 2 in parent rest frame.
Double_t getmPrime() const
Get m&#39; value.
Bool_t gotFullySymmetricDP() const
Is the DP fully symmetric?
std::vector< Double_t > mass_
Vector of daughter particles masses.
Double_t getc13() const
Get the cosine of the helicity angle theta13.
const Double_t m1Sq_
Mass of particle 1 squared.
Double_t getm3() const
Get the m3 mass.
Double_t cFromM(const Double_t mijSq, const Double_t mikSq, const Double_t mij, const Int_t i, const Int_t j, const Int_t k) const
General method to calculate the cos(helicity) variables from the masses of the particles.
Double_t getp2_13() const
Get the momentum of the track 2 in 13 rest frame.
Double_t getm1() const
Get the m1 mass.
Double_t getm23SqMax() const
Get the m23Sq maximum, (mParent - m1)^2.
Double_t m13_
Invariant mass m13.
Double_t p1_13_
Momentum of track 1 in 1-3 rest frame.
Bool_t warnings_
Enable/disable warning messages.
Double_t getp1_12() const
Get the momentum of the track 1 in 12 rest frame.
void updateSqDPKinematics(const Double_t mPrime, const Double_t thetaPrime)
Update all kinematic quantities based on the square DP co-ordinates m&#39; and theta&#39;.
Double_t getm13Min() const
Get the m13 minimum defined as (m1 + m3)
Double_t getm13SqMax() const
Get the m13Sq maximum, (mParent - m2)^2.
Double_t getm23Sq() const
Get the m23 invariant mass square.
Double_t distanceFromDPCentre() const
Calculate the distance from the currently set (m13Sq, m23Sq) point to the centre of the Dalitz plot (...
void updateMassSq_m23(const Double_t m23, const Double_t c23)
Update the variables m12Sq_ and m13Sq_ given the invariant mass m23 and the cosine of the helicity an...
const Double_t mDTot_
Sum of the daughter masses.
const Double_t mSqDTot_
Sum of the squares of the daughter masses.
Double_t getm12SqMin() const
Get the m12Sq minimum, (m1 + m2)^2.
Double_t mFromC(const Double_t mijSq, const Double_t cij, const Double_t mij, const Int_t i, const Int_t j, const Int_t k) const
General method to calculate mikSq given mijSq and cosHel_ij.
std::vector< Double_t > mMin_
Vector of the minimum mij values.
Double_t getm13Max() const
Get the m13 maximum defined as (mParent - m2)
Double_t genm23Sq() const
Randomly generate the invariant mass squared m23Sq.
Double_t getp1_Parent() const
Get the momentum of the track 1 in parent rest frame.
Double_t getm12Min() const
Get the m12 minimum defined as (m1 + m2)
void calcm12Sq()
Calculate m12Sq from m13Sq and m23Sq.
Double_t p2_23_
Momentum of track 2 in 2-3 rest frame.
Double_t thetaPrime_
theta&#39; co-ordinate
std::vector< Double_t > mSqDiff_
Vector of the difference between the mSqMax and mSqMin.
Double_t calcSqDPJacobian() const
Calculate the Jacobian for the transformation m23^2, m13^2 -&gt; m&#39;, theta&#39; (square DP) at the currently ...
void updateKinematicsFrom23(const Double_t m23, const Double_t c23)
Update all kinematic quantities based on m23 and cos(theta23)
Double_t getm2() const
Get the m2 mass.
Double_t getcov13() const
Get covariant factor in 13 axis.
Bool_t withinDPLimits2(const Double_t m13Sq, const Double_t m23Sq) const
Check whether a given (m13Sq,m23Sq) point is within the kinematic limits of the Dalitz plot (alternat...
Double_t p2_13_
Momentum of track 2 in 1-3 rest frame.
void genFlatSqDP(Double_t &mPrime, Double_t &thetaPrime) const
Routine to generate events flat in the square Dalitz plot.
Double_t p1_Parent_
Momentum of track 1 in parent rest frame.
Double_t calcThirdMassSq(const Double_t firstMassSq, const Double_t secondMassSq) const
Calculate the third invariant mass square from the two provided (e.g. mjkSq from mijSq and mikSq) ...
Double_t getm12Max() const
Get the m12 maximum defined as (mParent - m3)
Double_t getp1_23() const
Get the momentum of the track 1 in 23 rest frame.
const Double_t m3Sq_
Mass of particle 3 squared.
Double_t getc12() const
Get the cosine of the helicity angle theta12.
Double_t getm13() const
Get the m13 invariant mass.
void flipAndUpdateKinematics()
Flips the DP variables m13^2 &lt;-&gt; m23^2 and recalculates all kinematic quantities. ...
Double_t getp1_13() const
Get the momentum of the track 1 in 13 rest frame.
Double_t getm23Min() const
Get the m23 minimum defined as (m2 + m3)
std::vector< Double_t > mMax_
Vector of the maximum mij values.
Double_t getmParent() const
Get parent mass.
Double_t m13Sq_
Invariant mass m13 squared.
void updateKinematics(const Double_t m13Sq, const Double_t m23Sq)
Update all kinematic quantities based on the DP co-ordinates m13Sq and m23Sq.
Double_t m12Sq_
Invariant mass m12 squared.
void updateSqDPMassSquares(const Double_t mPrime, const Double_t thetaPrime)
Update some kinematic quantities based on the square DP co-ordinates m&#39; and theta&#39;.
Double_t qk_
Momentum q of particle k.
Double_t p3_Parent_
Momentum of track 3 in parent rest frame.
void rotateAndUpdateKinematics()
Cyclically rotates the DP variables (m12 -&gt; m23, m23 -&gt; m13, m13 -&gt; m12) and recalculates all kinemat...
const Double_t m3_
Mass of particle 3.
Double_t genm13Sq() const
Randomly generate the invariant mass squared m13Sq.
Double_t getm12Sq() const
Get the m12 invariant mass square.
Double_t getp3_12() const
Get the momentum of the track 3 in 12 rest frame.
std::vector< Double_t > mSqMax_
Vector of the maximum mijSq values.
Double_t c23_
Cosine of the helicity angle theta23, which is defined as the angle between 1&amp;2 in the rest frame of ...
const Double_t massInt_
Mass difference between the parent particle and the sum of the daughter particles.
Double_t dal(Double_t *x, Double_t *par)
Double_t getm13Sq() const
Get the m13 invariant mass square.
void calcParentFrameMomenta()
Calculate the momenta of each daughter in the parent rest frame.
void updateMassSquares(const Double_t m13Sq, const Double_t m23Sq)
Update some kinematic quantities based on the DP co-ordinates m13Sq and m23Sq.
const Bool_t fullySymmetricDP_
Fully-symmetrical DP.
std::vector< Double_t > mDiff_
Vector of the difference between the mMax and mMin.
Double_t getm12() const
Get the m12 invariant mass.
Double_t getDPBoxArea() const
Get the box area defined from the kinematic bounds.
Bool_t squareDP_
Should we calculate the square DP co-ordinates or not?
Double_t m12_
Invariant mass m12.
Double_t getp2_23() const
Get the momentum of the track 2 in 23 rest frame.
Double_t c12_
Cosine of the helicity angle theta12, which is defined as the angle between 1&amp;3 in the rest frame of ...
const Double_t mParent_
Mass of parent particle.
void drawDPContour(const Int_t orientation=1323, const Int_t nbins=100)
Method to draw the Dalitz plot contours on the top of the histo previously drawn. ...
Double_t p3_12_
Momentum of track 3 in 1-2 rest frame.
Double_t getThetaPrime() const
Get theta&#39; value.
Double_t m23Sq_
Invariant mass m23 squared.
Double_t getm12SqMax() const
Get the m12Sq maximum, (mParent - m3)^2.
void genFlatPhaseSpace(Double_t &m13Sq, Double_t &m23Sq) const
Routine to generate events flat in phase-space.
virtual ~LauKinematics()
Destructor.
const Bool_t symmetricalDP_
Symmetrical DP.
Class for calculating 3-body kinematic quantities.
LauKinematics(const Double_t m1, const Double_t m2, const Double_t m3, const Double_t mParent, const Bool_t calcSquareDPCoords=kFALSE, const Bool_t symmetricalDP=kFALSE, const Bool_t fullySymmetricDP=kFALSE)
Constructor.
Double_t genm12Sq() const
Randomly generate the invariant mass squared m12Sq.
const Double_t m2_
Mass of particle 2.
Double_t c13_
Cosine of the helicity angle theta13, which is defined as the angle between 1&amp;2 in the rest frame of ...
void updateKinematicsFrom13(const Double_t m13, const Double_t c13)
Update all kinematic quantities based on m13 and cos(theta13)
void calcHelicities()
Calculate cosines of the helicity angles, momenta of daughters and bachelor in various ij rest frames...
const Double_t m2Sq_
Mass of particle 2 squared.
void updateKinematicsFrom12(const Double_t m12, const Double_t c12)
Update all kinematic quantities based on m12 and cos(theta12)
Double_t p1_12_
Momentum of track 1 in 1-2 rest frame.
Double_t pCalc(const Double_t energy, const Double_t massSq) const
General method to calculate the momentum of a particle, given its energy and invariant mass squared...
Double_t mPrime_
m&#39; co-ordinate
Bool_t withinSqDPLimits(const Double_t mPrime, const Double_t thetaPrime) const
Check whether a given (m&#39;,theta&#39;) point is within the kinematic limits of the Dalitz plot...
Double_t getmParentSq() const
Get parent mass squared.
Double_t getp3_Parent() const
Get the momentum of the track 3 in parent rest frame.
std::vector< Double_t > mSq_
Vector of daughter particles masses squared.
void calcSqDPVars()
Calculate the m&#39; and theta&#39; variables for the square Dalitz plot.
const Double_t m1_
Mass of particle 1.