Page Menu
Home
HEPForge
Search
Configure Global Search
Log In
Files
F19300901
D34.1759167517.diff
No One
Temporary
Actions
View File
Edit File
Delete File
View Transforms
Subscribe
Award Token
Flag For Later
Size
12 KB
Referenced Files
None
Subscribers
None
D34.1759167517.diff
View Options
diff --git a/EvtGenModels/EvtEtaLLPiPi.hh b/EvtGenModels/EvtEtaLLPiPi.hh
new file mode 100644
--- /dev/null
+++ b/EvtGenModels/EvtEtaLLPiPi.hh
@@ -0,0 +1,81 @@
+
+/***********************************************************************
+* Copyright 1998-2020 CERN for the benefit of the EvtGen authors *
+* *
+* This file is part of EvtGen. *
+* *
+* EvtGen is free software: you can redistribute it and/or modify *
+* it under the terms of the GNU General Public License as published by *
+* the Free Software Foundation, either version 3 of the License, or *
+* (at your option) any later version. *
+* *
+* EvtGen is distributed in the hope that it will be useful, *
+* but WITHOUT ANY WARRANTY; without even the implied warranty of *
+* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+* GNU General Public License for more details. *
+* *
+* You should have received a copy of the GNU General Public License *
+* along with EvtGen. If not, see <https://www.gnu.org/licenses/>. *
+***********************************************************************/
+
+#ifndef EVT_ETALLPIPI_HH
+#define EVT_ETALLPIPI_HH
+
+#include "EvtGenBase/EvtDecayProb.hh"
+#include <string>
+
+class EvtParticle;
+
+class EvtEtaLLPiPi : public EvtDecayProb {
+public:
+
+ EvtEtaLLPiPi();
+
+ void init() override;
+ void initProbMax() override;
+
+ std::string getName() override;
+ EvtDecayBase *clone() override;
+
+ void decay( EvtParticle *p ) override;
+
+private:
+
+ void updateMassPars( double mLep, double mPi );
+
+ double rhoWidth( double s, double m ) const;
+
+ double F0( double sLL, double sPiPi ) const;
+
+ double lambda( double a, double b, double c ) const;
+
+ double ampSquared( EvtParticle *p ) const;
+
+ double m_alpha;
+ double m_eSq;
+ double m_fPi;
+ double m_f8;
+ double m_f0;
+ double m_thetaMix;
+ double m_mixSq;
+ double m_c1;
+ double m_c2;
+ double m_c3;
+ double m_par1;
+ double m_parLL;
+ double m_parPiPi;
+
+ double m_rhoMass;
+ double m_rhoMassSq;
+ double m_rhoGamma;
+
+ double m_lepMass;
+ double m_lepMassSq;
+ double m_piMass;
+ double m_piMassSq;
+ double m_4LepMassSq;
+ double m_4PiMassSq;
+
+};
+
+#endif
diff --git a/History.txt b/History.txt
--- a/History.txt
+++ b/History.txt
@@ -4,6 +4,15 @@
//
//===========================================================================
+21st Aug 2020 John Back
+ Add EvtEtaLLPiPi model for eta' -> l l pi pi decays (l = e or mu),
+ courtesy of Aleksei Luchinsky (LHCb).
+
+15th May 2020 Michal Kreps
+ Add EvtThreeBodyPhsp (rectangular Dalitz plot selection) and
+ EvtPhspDecaytimeCut (minimum decay time requirement) models.
+ Fix initialisation of constants in EvtBTo3hCP model.
+
R02-00-00-------------------------------------------------------------------
24th Apr 2020 Michal Kreps
Update particle properties file evt.pdl to 2019 version of RPP by PDG.
diff --git a/src/EvtGenModels/EvtEtaLLPiPi.cpp b/src/EvtGenModels/EvtEtaLLPiPi.cpp
new file mode 100644
--- /dev/null
+++ b/src/EvtGenModels/EvtEtaLLPiPi.cpp
@@ -0,0 +1,250 @@
+
+/***********************************************************************
+* Copyright 1998-2020 CERN for the benefit of the EvtGen authors *
+* *
+* This file is part of EvtGen. *
+* *
+* EvtGen is free software: you can redistribute it and/or modify *
+* it under the terms of the GNU General Public License as published by *
+* the Free Software Foundation, either version 3 of the License, or *
+* (at your option) any later version. *
+* *
+* EvtGen is distributed in the hope that it will be useful, *
+* but WITHOUT ANY WARRANTY; without even the implied warranty of *
+* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+* GNU General Public License for more details. *
+* *
+* You should have received a copy of the GNU General Public License *
+* along with EvtGen. If not, see <https://www.gnu.org/licenses/>. *
+***********************************************************************/
+
+#include "EvtGenModels/EvtEtaLLPiPi.hh"
+
+#include "EvtGenBase/EvtConst.hh"
+#include "EvtGenBase/EvtId.hh"
+#include "EvtGenBase/EvtKine.hh"
+#include "EvtGenBase/EvtParticle.hh"
+#include "EvtGenBase/EvtPDL.hh"
+#include "EvtGenBase/EvtSpinType.hh"
+#include "EvtGenBase/EvtVector4R.hh"
+
+#include <cmath>
+
+EvtEtaLLPiPi::EvtEtaLLPiPi() :
+ m_alpha( 1.0 / 137.0 ),
+ m_eSq( 4.0 * EvtConst::pi * m_alpha ),
+ m_fPi( 0.0924 ),
+ m_f8( 1.3 * m_fPi ),
+ m_f0( 1.04 * m_fPi ),
+ m_thetaMix( 20.0 * EvtConst::pi / 180.0 ),
+ m_mixSq( 0.0 ),
+ m_c1( 1.0 ),
+ m_c2( 0.0 ),
+ m_c3( m_c1 - m_c2 ), // Eq 9
+ m_par1( 1.0 - ( 3.0 * ( m_c1 - m_c2 + m_c3 ) / 4.0 ) ),
+ m_parLL( 3.0 * ( m_c1 - m_c2 - m_c3 ) / 4.0 ),
+ m_parPiPi( 3.0 * m_c3 / 2.0 ),
+ m_rhoMass( 0.775 ), // updated in init()
+ m_rhoMassSq( m_rhoMass * m_rhoMass ),
+ m_rhoGamma( 0.149 ), // updated in init()
+ m_lepMass( 0.106 ), // modified in updateMassPars()
+ m_lepMassSq ( m_lepMass * m_lepMass ),
+ m_piMass( 0.140 ), // modified in updateMassPars()
+ m_piMassSq( m_piMass * m_piMass ),
+ m_4LepMassSq( 4.0 * m_lepMassSq ),
+ m_4PiMassSq( 4.0 * m_piMassSq )
+{
+ // eta' -> mu+ mu- pi+ pi- or e+ e- pi+ pi-
+ // From Zhang Zhen-Yu et al, Chinese Phys. C 36, p926, 2012
+}
+
+void EvtEtaLLPiPi::init() {
+
+ // Check for 0 or 1 (maxProb) arguments
+ checkNArg( 0, 1 );
+
+ // Check particle types
+ checkSpinParent( EvtSpinType::SCALAR );
+ checkSpinDaughter( 0, EvtSpinType::DIRAC );
+ checkSpinDaughter( 1, EvtSpinType::DIRAC );
+ checkSpinDaughter( 2, EvtSpinType::SCALAR );
+ checkSpinDaughter( 3, EvtSpinType::SCALAR );
+
+ // Mass and width of rho0 from particle properties file
+ m_rhoMass = EvtPDL::getMeanMass( EvtPDL::getId("rho0") );
+ m_rhoMassSq = m_rhoMass * m_rhoMass;
+ m_rhoGamma = EvtPDL::getWidth( EvtPDL::getId( "rho0" ) );
+
+ // Mixing parameter squared, using Eq 6
+ const double denom = 8.0 * pow( EvtConst::pi * m_fPi, 2 );
+ const double factor = m_eSq / ( denom * denom * 3.0 );
+ const double fTerm8 = sin( m_thetaMix ) / m_f8;
+ const double fTerm0 = 2.0 * sqrt( 2.0 ) * cos( m_thetaMix ) / m_f0;
+ m_mixSq = factor * pow( fTerm8 + fTerm0, 2 );
+
+}
+
+void EvtEtaLLPiPi::initProbMax() {
+
+ if ( getNArg() == 1 ) {
+ setProbMax( getArg( 0 ) );
+
+ } else {
+
+ int lepId = getDaug(0).getId();
+ if ( lepId == EvtPDL::getId( "e-" ).getId() ||
+ lepId == EvtPDL::getId( "e+" ).getId() ) {
+ setProbMax( 27500.0 );
+
+ } else if ( lepId == EvtPDL::getId( "mu-" ).getId() ||
+ lepId == EvtPDL::getId( "mu+" ).getId() ) {
+ setProbMax( 20.0 );
+ }
+
+ }
+
+}
+
+std::string EvtEtaLLPiPi::getName() {
+ return "ETA_LLPIPI";
+}
+
+EvtDecayBase *EvtEtaLLPiPi::clone() {
+ return new EvtEtaLLPiPi();
+}
+
+void EvtEtaLLPiPi::decay( EvtParticle *p ) {
+
+ p->initializePhaseSpace( getNDaug(), getDaugs() );
+
+ const double mLep = p->getDaug( 0 )->mass();
+ const double mPi = p->getDaug( 2 )->mass();
+
+ updateMassPars( mLep, mPi );
+
+ const double prob = ampSquared( p );
+ setProb( prob );
+
+}
+
+void EvtEtaLLPiPi::updateMassPars( double mLep, double mPi ) {
+
+ // Update mass parameters used in various functions
+ m_lepMass = mLep;
+ m_lepMassSq = mLep*mLep;
+ m_4LepMassSq = 4.0*m_lepMassSq;
+
+ m_piMass = mPi;
+ m_piMassSq = mPi*mPi;
+ m_4PiMassSq = 4.0*m_piMassSq;
+
+}
+
+double EvtEtaLLPiPi::rhoWidth( double s, double m ) const {
+
+ // Define width of rho using modified vector meson dynamics, Eq 8
+ double gamma( 0.0 );
+
+ if ( s >= m_4PiMassSq ) {
+
+ const double mSq = m*m;
+ const double num = 1.0 - ( 4.0 * mSq / s );
+ const double denom = 1.0 - ( 4.0 * mSq / m_rhoMassSq );
+ const double ratio = denom > 0.0 ? num/denom : 0.0;
+ gamma = m_rhoGamma * ( s / m_rhoMassSq ) * pow ( ratio, 1.5 );
+
+ }
+
+ return gamma;
+}
+
+double EvtEtaLLPiPi::F0( double sLL, double sPiPi ) const {
+
+ // Equation 7
+ double ampSq( 0.0 );
+
+ const double rhoWidthL = rhoWidth( sLL, m_lepMass );
+ const double rhoWidthPi = rhoWidth( sPiPi, m_piMass );
+
+ const double mSqDiffL = m_rhoMassSq - sLL;
+ const double mRhoWidthL = m_rhoMass * rhoWidthL;
+
+ const double mSqDiffPi = m_rhoMassSq - sPiPi;
+ const double mRhoWidthPi = m_rhoMass * rhoWidthPi;
+
+ const double denomLL = mSqDiffL * mSqDiffL + mRhoWidthL * mRhoWidthL;
+ const double denomPiPi = mSqDiffPi * mSqDiffPi + mRhoWidthPi * mRhoWidthPi;
+
+ if ( denomLL > 0.0 && denomPiPi > 0.0 ) {
+
+ const double mRho4 = m_rhoMassSq * m_rhoMassSq;
+ const double denomProd = denomLL * denomPiPi;
+
+ double realAmp = m_par1 + m_parLL * ( m_rhoMassSq * mSqDiffL / denomLL );
+ realAmp += m_parPiPi * mRho4 *( ( mSqDiffPi * mSqDiffL ) - mRhoWidthL * mRhoWidthPi ) / denomProd;
+
+ double imagAmp = m_parLL * ( m_rhoMassSq * mRhoWidthL / denomLL );
+ imagAmp += m_parPiPi * mRho4 * ( mRhoWidthPi * mSqDiffL + mRhoWidthL * mSqDiffPi ) / denomProd;
+
+ ampSq = realAmp * realAmp + imagAmp * imagAmp;
+
+ }
+
+ return ampSq;
+
+}
+
+double EvtEtaLLPiPi::lambda( double a, double b, double c ) const {
+
+ const double sumSq = a * a + b * b + c * c;
+ const double prod = a * b + b * c + c * a;
+ const double L = sumSq - 2.0 * prod;
+ return L;
+
+}
+
+double EvtEtaLLPiPi::ampSquared( EvtParticle *p ) const {
+
+ // Equation 3
+ const double zeroProb( 0.0 );
+
+ // Mass of eta' meson
+ const double mEta = p->mass();
+
+ const EvtVector4R pL1 = p->getDaug(0)->getP4();
+ const EvtVector4R pL2 = p->getDaug(1)->getP4();
+ const EvtVector4R pPi1 = p->getDaug(2)->getP4();
+ const EvtVector4R pPi2 = p->getDaug(3)->getP4();
+
+ const EvtVector4R pLL = pL1 + pL2;
+ const double sLL = pLL.mass2();
+ const EvtVector4R pPiPi = pPi1 + pPi2;
+ const double sPiPi = pPiPi.mass2();
+
+ if ( sLL < 1e-4 || sPiPi < m_4PiMassSq || sLL < m_4LepMassSq ) {
+ // To avoid negative square roots etc
+ return zeroProb;
+ }
+
+ // Angles theta_p, theta_k and phi defined in Fig 1
+ const EvtVector4R pSum = pLL + pPiPi;
+ // Helicity angle of first lepton
+ const double cosThp = -EvtDecayAngle( pSum, pLL, pL1 );
+ const double sinThp = sqrt( 1.0 - cosThp * cosThp );
+ // Helicity angle of first pion
+ const double cosThk = -EvtDecayAngle( pSum, pPiPi, pPi2 );
+ const double sinThk = sqrt( 1.0 - cosThk * cosThk );
+ // Angle between the dilepton and dipion decay planes
+ const double phi = EvtDecayAngleChi( pSum, pL1, pL2, pPi1, pPi2 );
+ const double sinPhi = sin( phi );
+
+ const double betaLL = sqrt( 1.0 - ( m_4LepMassSq / sLL ) );
+ const double betaPiPi = sqrt( 1.0 - ( m_4PiMassSq / sPiPi ) );
+
+ const double betaProd = ( 1.0 - pow( betaLL * sinThp * sinPhi, 2) ) * sPiPi * pow( betaPiPi * sinThk, 2 );
+ const double L = lambda( mEta*mEta, sLL, sPiPi );
+ const double ampSq = m_eSq * F0( sLL, sPiPi ) * m_mixSq * L * betaProd / ( 8.0 * sLL );
+
+ return ampSq;
+
+}
diff --git a/src/EvtGenModels/EvtModelReg.cpp b/src/EvtGenModels/EvtModelReg.cpp
--- a/src/EvtGenModels/EvtModelReg.cpp
+++ b/src/EvtGenModels/EvtModelReg.cpp
@@ -160,6 +160,7 @@
#include "EvtGenModels/EvtbsToLLLL.hh"
#include "EvtGenModels/EvtbsToLLLLHyperCP.hh"
#include "EvtGenModels/EvtThreeBodyPhsp.hh"
+#include "EvtGenModels/EvtEtaLLPiPi.hh"
#include <assert.h>
#include <ctype.h>
@@ -341,4 +342,5 @@
modelist.registerModel( new EvtPsi2JpsiPiPi );
modelist.registerModel( new EvtThreeBodyPhsp );
+ modelist.registerModel( new EvtEtaLLPiPi );
}
File Metadata
Details
Attached
Mime Type
text/plain
Expires
Mon, Sep 29, 6:38 PM (19 h, 22 s)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
6571010
Default Alt Text
D34.1759167517.diff (12 KB)
Attached To
Mode
D34: Add EvtEtaLLPiPi model.
Attached
Detach File
Event Timeline
Log In to Comment