Line data Source code
1 : //--------------------------------------------------------------------------
2 : //
3 : // Environment:
4 : // This software is part of the EvtGen package developed jointly
5 : // for the BaBar and CLEO collaborations. If you use all or part
6 : // of it, please give an appropriate acknowledgement.
7 : //
8 : // Copyright Information: See EvtGen/COPYRIGHT
9 : // Copyright (C) 1998 Caltech, UCSB
10 : //
11 : // Module: EvtResonance2.cc
12 : //
13 : // Description: resonance-defining class
14 : //
15 : // Modification history:
16 : //
17 : // NK September 4, 1997 Module created
18 : //
19 : //------------------------------------------------------------------------
20 : //
21 : #include "EvtGenBase/EvtPatches.hh"
22 : #include <math.h>
23 : #include "EvtGenBase/EvtVector4R.hh"
24 : #include "EvtGenBase/EvtKine.hh"
25 : #include "EvtGenBase/EvtComplex.hh"
26 : #include "EvtGenBase/EvtResonance2.hh"
27 : #include "EvtGenBase/EvtReport.hh"
28 : #include "EvtGenBase/EvtConst.hh"
29 :
30 0 : EvtResonance2::~EvtResonance2(){}
31 :
32 : EvtResonance2& EvtResonance2::operator = ( const EvtResonance2 &n)
33 : {
34 0 : if ( &n == this ) return *this;
35 0 : _p4_p = n._p4_p;
36 0 : _p4_d1 = n._p4_d1;
37 0 : _p4_d2 = n._p4_d2;
38 0 : _ampl = n._ampl;
39 0 : _theta = n._theta;
40 0 : _gamma = n._gamma;
41 0 : _spin = n._spin;
42 0 : _bwm = n._bwm;
43 0 : _invmass_angdenom = n._invmass_angdenom;
44 0 : return *this;
45 0 : }
46 :
47 :
48 : EvtResonance2::EvtResonance2(const EvtVector4R& p4_p, const EvtVector4R& p4_d1,
49 : const EvtVector4R& p4_d2, double ampl,
50 : double theta, double gamma, double bwm, int spin,
51 : bool invmass_angdenom):
52 0 : _p4_p(p4_p),_p4_d1(p4_d1), _p4_d2(p4_d2),_ampl(ampl), _theta(theta),
53 0 : _gamma(gamma), _bwm(bwm), _spin(spin), _invmass_angdenom(invmass_angdenom) {}
54 :
55 :
56 : EvtComplex EvtResonance2::resAmpl() {
57 :
58 0 : double pi180inv = 1.0/EvtConst::radToDegrees;
59 :
60 0 : EvtComplex ampl;
61 0 : EvtVector4R p4_d3 = _p4_p-_p4_d1-_p4_d2;
62 :
63 : //get cos of the angle between the daughters from their 4-momenta
64 : //and the 4-momentum of the parent
65 :
66 : //in general, EvtDecayAngle(parent, part1+part2, part1) gives the angle
67 : //the missing particle (not listed in the arguments) makes
68 : //with part2 in the rest frame of both
69 : //listed particles (12)
70 :
71 : //angle 3 makes with 2 in rest frame of 12 (CS3)
72 : //double cos_phi_0 = EvtDecayAngle(_p4_p, _p4_d1+_p4_d2, _p4_d1);
73 : //angle 3 makes with 1 in 12 is, of course, -cos_phi_0
74 :
75 : //first compute several quantities...follow CLEO preprint 00-23
76 :
77 0 : double mAB=(_p4_d1+_p4_d2).mass();
78 0 : double mBC=(_p4_d2+p4_d3).mass();
79 0 : double mAC=(_p4_d1+p4_d3).mass();
80 0 : double mA=_p4_d1.mass();
81 0 : double mB=_p4_d2.mass();
82 0 : double mD=_p4_p.mass();
83 0 : double mC=p4_d3.mass();
84 :
85 0 : double mR=_bwm;
86 0 : double gammaR=_gamma;
87 0 : double mdenom = _invmass_angdenom ? mAB : mR;
88 0 : double pAB=sqrt( (((mAB*mAB-mA*mA-mB*mB)*(mAB*mAB-mA*mA-mB*mB)/4.0) -
89 0 : mA*mA*mB*mB)/(mAB*mAB));
90 0 : double pR=sqrt( (((mR*mR-mA*mA-mB*mB)*(mR*mR-mA*mA-mB*mB)/4.0) -
91 0 : mA*mA*mB*mB)/(mR*mR));
92 :
93 0 : double pD= (((mD*mD-mR*mR-mC*mC)*(mD*mD-mR*mR-mC*mC)/4.0) -
94 0 : mR*mR*mC*mC)/(mD*mD);
95 0 : if ( pD>0 ) { pD=sqrt(pD); } else {pD=0;}
96 0 : double pDAB=sqrt( (((mD*mD-mAB*mAB-mC*mC)*(mD*mD-mAB*mAB-mC*mC)/4.0) -
97 0 : mAB*mAB*mC*mC)/(mD*mD));
98 :
99 :
100 :
101 : double fR=1;
102 : double fD=1;
103 : int power=0;
104 0 : switch (_spin) {
105 : case 0:
106 : fR=1.0;
107 : fD=1.0;
108 : power=1;
109 0 : break;
110 : case 1:
111 0 : fR=sqrt(1.0+1.5*1.5*pR*pR)/sqrt(1.0+1.5*1.5*pAB*pAB);
112 0 : fD=sqrt(1.0+5.0*5.0*pD*pD)/sqrt(1.0+5.0*5.0*pDAB*pDAB);
113 : power=3;
114 0 : break;
115 : case 2:
116 0 : fR = sqrt( (9+3*pow((1.5*pR),2)+pow((1.5*pR),4))/(9+3*pow((1.5*pAB),2)+pow((1.5*pAB),4)) );
117 0 : fD = sqrt( (9+3*pow((5.0*pD),2)+pow((5.0*pD),4))/(9+3*pow((5.0*pDAB),2)+pow((5.0*pDAB),4)) );
118 : power=5;
119 0 : break;
120 : default:
121 0 : report(Severity::Info,"EvtGen") << "Incorrect spin in EvtResonance22.cc\n";
122 0 : }
123 :
124 0 : double gammaAB= gammaR*pow(pAB/pR,power)*(mR/mAB)*fR*fR;
125 0 : switch (_spin) {
126 : case 0:
127 0 : ampl=_ampl*EvtComplex(cos(_theta*pi180inv),sin(_theta*pi180inv))*
128 0 : fR*fD/(mR*mR-mAB*mAB-EvtComplex(0.0,mR*gammaAB));
129 0 : break;
130 : case 1:
131 0 : ampl=_ampl*EvtComplex(cos(_theta*pi180inv),sin(_theta*pi180inv))*
132 0 : (fR*fD*(mAC*mAC-mBC*mBC+((mD*mD-mC*mC)*(mB*mB-mA*mA)/(mdenom*mdenom)))/
133 0 : (mR*mR-mAB*mAB-EvtComplex(0.0,mR*gammaAB)));
134 0 : break;
135 : case 2:
136 0 : ampl=_ampl*EvtComplex(cos(_theta*pi180inv),sin(_theta*pi180inv))*
137 0 : fR*fD/(mR*mR-mAB*mAB-EvtComplex(0.0,mR*gammaAB))*
138 0 : (pow((mBC*mBC-mAC*mAC+(mD*mD-mC*mC)*(mA*mA-mB*mB)/(mdenom*mdenom)),2)-
139 0 : (1.0/3.0)*(mAB*mAB-2*mD*mD-2*mC*mC+pow((mD*mD- mC*mC)/mdenom, 2))*
140 0 : (mAB*mAB-2*mA*mA-2*mB*mB+pow((mA*mA-mB*mB)/mdenom,2)));
141 0 : break;
142 :
143 : default:
144 0 : report(Severity::Info,"EvtGen") << "Incorrect spin in EvtResonance22.cc\n";
145 0 : }
146 :
147 : return ampl;
148 0 : }
149 :
150 :
151 :
|