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: EvtResonance.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/EvtResonance.hh"
27 : #include "EvtGenBase/EvtReport.hh"
28 : #include "EvtGenBase/EvtConst.hh"
29 : using std::endl;
30 :
31 0 : EvtResonance::~EvtResonance(){}
32 :
33 :
34 : EvtResonance& EvtResonance::operator = ( const EvtResonance &n)
35 : {
36 0 : if ( &n == this ) return *this;
37 0 : _p4_p = n._p4_p;
38 0 : _p4_d1 = n._p4_d1;
39 0 : _p4_d2 = n._p4_d2;
40 0 : _ampl = n._ampl;
41 0 : _theta = n._theta;
42 0 : _gamma = n._gamma;
43 0 : _spin = n._spin;
44 0 : _bwm = n._bwm;
45 0 : return *this;
46 0 : }
47 :
48 :
49 : EvtResonance::EvtResonance(const EvtVector4R& p4_p, const EvtVector4R& p4_d1,
50 : const EvtVector4R& p4_d2, double ampl,
51 : double theta, double gamma, double bwm, int spin):
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) {}
54 :
55 : EvtComplex EvtResonance::resAmpl() {
56 :
57 0 : double pi180inv = 1.0/EvtConst::radToDegrees;
58 :
59 0 : EvtComplex ampl;
60 : //EvtVector4R _p4_d3 = _p4_p-_p4_d1-_p4_d2;
61 :
62 : //get cos of the angle between the daughters from their 4-momenta
63 : //and the 4-momentum of the parent
64 :
65 : //in general, EvtDecayAngle(parent, part1+part2, part1) gives the angle
66 : //the missing particle (not listed in the arguments) makes
67 : //with part2 in the rest frame of both
68 : //listed particles (12)
69 :
70 : //angle 3 makes with 2 in rest frame of 12 (CS3)
71 0 : double cos_phi_0 = EvtDecayAngle(_p4_p, _p4_d1+_p4_d2, _p4_d1);
72 : //angle 3 makes with 1 in 12 is, of course, -cos_phi_0
73 :
74 0 : switch (_spin) {
75 :
76 : case 0 :
77 0 : ampl=(_ampl*EvtComplex(cos(_theta*pi180inv),sin(_theta*pi180inv))*
78 0 : sqrt(_gamma/EvtConst::twoPi)*
79 0 : (1.0/((_p4_d1+_p4_d2).mass()-_bwm-EvtComplex(0.0,0.5*_gamma))));
80 0 : break;
81 :
82 : case 1 :
83 0 : ampl=(_ampl*EvtComplex(cos(_theta*pi180inv),sin(_theta*pi180inv))*
84 0 : sqrt(_gamma/EvtConst::twoPi)*
85 0 : (cos_phi_0/((_p4_d1+_p4_d2).mass()-_bwm-EvtComplex(0.0,0.5*_gamma))));
86 0 : break;
87 :
88 : case 2:
89 0 : ampl=(_ampl*EvtComplex(cos(_theta*pi180inv),sin(_theta*pi180inv))*
90 0 : sqrt(_gamma/EvtConst::twoPi)*
91 0 : ((1.5*cos_phi_0*cos_phi_0-0.5)/((_p4_d1+_p4_d2).mass()-_bwm-EvtComplex(0.0, 0.5*_gamma))));
92 0 : break;
93 :
94 : case 3:
95 0 : ampl=(_ampl*EvtComplex(cos(_theta*pi180inv),sin(_theta*pi180inv))*
96 0 : sqrt(_gamma/EvtConst::twoPi)*
97 0 : ((2.5*cos_phi_0*cos_phi_0*cos_phi_0-1.5*cos_phi_0)/((_p4_d1+_p4_d2).mass()-_bwm-EvtComplex(0.0, 0.5*_gamma))));
98 0 : break;
99 :
100 : default:
101 0 : report(Severity::Debug,"EvtGen") << "EvtGen: wrong spin in EvtResonance" << endl;
102 0 : ampl = EvtComplex(0.0);
103 0 : break;
104 :
105 : }
106 :
107 : return ampl;
108 0 : }
109 :
110 : EvtComplex EvtResonance::relBrWig(int i) {
111 :
112 : //this function returns relativistic Breit-Wigner amplitude
113 : //for a given resonance (for P-wave decays of scalars only at the moment!)
114 :
115 0 : EvtComplex BW;
116 0 : EvtVector4R _p4_d3 = _p4_p-_p4_d1-_p4_d2;
117 0 : EvtVector4R _p4_12 = _p4_d1 + _p4_d2;
118 :
119 0 : double msq13 = (_p4_d1 + _p4_d3).mass2();
120 0 : double msq23 = (_p4_d2 + _p4_d3).mass2();
121 0 : double msqParent = _p4_p.mass2();
122 0 : double msq1 = _p4_d1.mass2();
123 0 : double msq2 = _p4_d2.mass2();
124 0 : double msq3 = _p4_d3.mass2();
125 :
126 : double M;
127 :
128 0 : double p2 = sqrt((_p4_12.mass2() - (_p4_d1.mass() + _p4_d2.mass())*(_p4_d1.mass() + _p4_d2.mass()))*(_p4_12.mass2() - (_p4_d1.mass() - _p4_d2.mass())*(_p4_d1.mass() - _p4_d2.mass())))/(2.0*_p4_12.mass());
129 :
130 0 : double p2R = sqrt((_bwm*_bwm - (_p4_d1.mass() + _p4_d2.mass())*(_p4_d1.mass() + _p4_d2.mass()))*(_bwm*_bwm - (_p4_d1.mass() - _p4_d2.mass())*(_p4_d1.mass() - _p4_d2.mass())))/(2.0*_bwm);
131 :
132 : double gam, R;
133 :
134 0 : if (i == 1) {
135 :
136 : R = 2.0/(0.197);
137 :
138 0 : }
139 : else R = 5.0/(0.197);
140 :
141 0 : gam = _gamma*(_bwm/_p4_12.mass())*(p2/p2R)*(p2/p2R)*(p2/p2R)*((1 + R*R*p2R*p2R)/(1 + R*R*p2*p2));
142 0 : M = (msq13 - msq23 - (msqParent - msq3)*(msq1 - msq2)/(_bwm*_bwm))*sqrt((1 + R*R*p2R*p2R)/(1 + R*R*p2*p2));
143 :
144 0 : BW = sqrt(_gamma)*M/((_bwm*_bwm - _p4_12.mass2()) - EvtComplex(0.0,1.0)*gam*_bwm);
145 :
146 : return BW;
147 :
148 0 : }
149 :
150 :
|