]>
Commit | Line | Data |
---|---|---|
da0e9ce3 | 1 | #include "EvtGenBase/EvtPatches.hh" |
2 | /******************************************************************************* | |
3 | * Project: BaBar detector at the SLAC PEP-II B-factory | |
4 | * Package: EvtGenBase | |
5 | * File: $Id: EvtPto3PAmp.cc,v 1.21 2009/02/18 03:31:38 ryd Exp $ | |
6 | * Author: Alexei Dvoretskii, dvoretsk@slac.stanford.edu, 2001-2002 | |
7 | * | |
8 | * Copyright (C) 2002 Caltech | |
9 | * | |
10 | * Modified: May 30, 2005, Denis Dujmic (ddujmic@slac.stanford.edu): flip sign | |
11 | * for vector resonances (-i -> i at resonance mass). | |
12 | * Introduce Flatte lineshape. | |
13 | *******************************************************************************/ | |
14 | ||
15 | #include <assert.h> | |
16 | #include <math.h> | |
17 | #include <iostream> | |
18 | #include "EvtGenBase/EvtComplex.hh" | |
19 | #include "EvtGenBase/EvtPto3PAmp.hh" | |
20 | #include "EvtGenBase/EvtDalitzCoord.hh" | |
21 | #include "EvtGenBase/EvtdFunction.hh" | |
22 | #include "EvtGenBase/EvtCyclic3.hh" | |
23 | #include "EvtGenBase/EvtReport.hh" | |
24 | ||
25 | using std::endl; | |
26 | using EvtCyclic3::Index; | |
27 | using EvtCyclic3::Pair; | |
28 | ||
29 | ||
30 | EvtPto3PAmp::EvtPto3PAmp(EvtDalitzPlot dp, Pair pairAng, Pair pairRes, | |
31 | EvtSpinType::spintype spin, | |
32 | const EvtPropagator& prop, NumType typeN) | |
33 | : EvtAmplitude<EvtDalitzPoint>(), | |
34 | _pairAng(pairAng), _pairRes(pairRes), | |
35 | _spin(spin), | |
36 | _typeN(typeN), | |
37 | _prop((EvtPropagator*) prop.clone()), | |
38 | _g0(prop.g0()), | |
39 | _min(0), | |
40 | _max(0), | |
41 | _vb(prop.m0(),dp.m(EvtCyclic3::other(pairRes)),dp.bigM(),spin), | |
42 | _vd(dp.m(EvtCyclic3::first(pairRes)),dp.m(EvtCyclic3::second(pairRes)),prop.m0(),spin) | |
43 | {} | |
44 | ||
45 | ||
46 | ||
47 | EvtPto3PAmp::EvtPto3PAmp(const EvtPto3PAmp& other) | |
48 | : EvtAmplitude<EvtDalitzPoint>(other), | |
49 | _pairAng(other._pairAng), | |
50 | _pairRes(other._pairRes), | |
51 | _spin(other._spin), | |
52 | _typeN(other._typeN), | |
53 | _prop( (other._prop) ? (EvtPropagator*) other._prop->clone() : 0), | |
54 | _g0(other._g0), | |
55 | _min(other._min), | |
56 | _max(other._max), | |
57 | _vb(other._vb), _vd(other._vd) | |
58 | {} | |
59 | ||
60 | ||
61 | EvtPto3PAmp::~EvtPto3PAmp() | |
62 | { | |
63 | if(_prop) delete _prop; | |
64 | } | |
65 | ||
66 | ||
67 | void EvtPto3PAmp::set_fd(double R) | |
68 | { | |
69 | _vd.set_f(R); | |
70 | } | |
71 | ||
72 | void EvtPto3PAmp::set_fb(double R) | |
73 | { | |
74 | _vb.set_f(R); | |
75 | } | |
76 | ||
77 | ||
78 | EvtComplex EvtPto3PAmp::amplitude(const EvtDalitzPoint& x) const | |
79 | { | |
80 | EvtComplex amp(1.0,0.0); | |
81 | ||
82 | double m = sqrt(x.q(_pairRes)); | |
83 | ||
84 | if ((_max>0 && m > _max) || (_min>0 && m < _min)) return EvtComplex(0.0,0.0); | |
85 | ||
86 | EvtTwoBodyKine vd(x.m(EvtCyclic3::first(_pairRes)), | |
87 | x.m(EvtCyclic3::second(_pairRes)),m); | |
88 | EvtTwoBodyKine vb(m,x.m(EvtCyclic3::other(_pairRes)),x.bigM()); | |
89 | ||
90 | ||
91 | // Compute mass-dependent width for relativistic propagators | |
92 | ||
93 | if(_typeN!=NBW && _typeN!=FLATTE) { | |
94 | ||
95 | _prop->set_g0(_g0*_vd.widthFactor(vd)); | |
96 | } | |
97 | ||
98 | ||
99 | // Compute propagator | |
100 | ||
101 | amp *= evalPropagator(m); | |
102 | ||
103 | ||
104 | ||
105 | // Compute form-factors | |
106 | ||
107 | amp *= _vd.formFactor(vd); | |
108 | amp *= _vb.formFactor(vb); | |
109 | ||
110 | amp *= numerator(x); | |
111 | ||
112 | return amp; | |
113 | } | |
114 | ||
115 | ||
116 | EvtComplex EvtPto3PAmp::numerator(const EvtDalitzPoint& x) const | |
117 | { | |
118 | EvtComplex ret(0.,0.); | |
119 | double m = sqrt(x.q(_pairRes)); | |
120 | EvtTwoBodyKine vd(x.m(EvtCyclic3::first(_pairRes)), | |
121 | x.m(EvtCyclic3::second(_pairRes)),m); | |
122 | EvtTwoBodyKine vb(m,x.m(EvtCyclic3::other(_pairRes)),x.bigM()); | |
123 | ||
124 | // Non-relativistic Breit-Wigner | |
125 | ||
126 | if(NBW == _typeN) { | |
127 | ||
128 | ret = angDep(x); | |
129 | } | |
130 | ||
131 | // Standard relativistic Zemach propagator | |
132 | ||
133 | else if(RBW_ZEMACH == _typeN) { | |
134 | ||
135 | ret = _vd.phaseSpaceFactor(vd,EvtTwoBodyKine::AB)*angDep(x); | |
136 | } | |
137 | ||
138 | // Kuehn-Santamaria normalization: | |
139 | ||
140 | else if(RBW_KUEHN == _typeN) { | |
141 | ||
142 | ret = _prop->m0()*_prop->m0() * angDep(x); | |
143 | } | |
144 | ||
145 | ||
146 | // CLEO amplitude is not factorizable | |
147 | // | |
148 | // The CLEO amplitude numerator is proportional to: | |
149 | // | |
150 | // m2_AC - m2_BC + (m2_D - m2_C)(m2_B - m2_A)/m2_0 | |
151 | // | |
152 | // m2_AC = (eA + eC)^2 + (P - P_C cosTh(BC))^2 | |
153 | // m2_BC = (eB + eC)^2 + (P + P_C cosTh(BC))^2 | |
154 | // | |
155 | // The first term m2_AB-m2_BC is therefore a p-wave term | |
156 | // - 4PP_C cosTh(BC) | |
157 | // The second term is an s-wave, the amplitude | |
158 | // does not factorize! | |
159 | // | |
160 | // The first term is just Zemach. However, the sign is flipped! | |
161 | // Let's consistently use the convention in which the amplitude | |
162 | // is proportional to +cosTh(BC). In the CLEO expressions, I will | |
163 | // therefore exchange AB to get rid of the sign flip. | |
164 | ||
165 | ||
166 | if(RBW_CLEO == _typeN || FLATTE == _typeN || GS == _typeN) { | |
167 | ||
168 | Index iA = EvtCyclic3::other(_pairAng); // A = other(BC) | |
169 | Index iB = EvtCyclic3::common(_pairRes,_pairAng); // B = common(AB,BC) | |
170 | Index iC = EvtCyclic3::other(_pairRes); // C = other(AB) | |
171 | ||
172 | double M = x.bigM(); | |
173 | double mA = x.m(iA); | |
174 | double mB = x.m(iB); | |
175 | double mC = x.m(iC); | |
176 | double qAB = x.q(EvtCyclic3::combine(iA,iB)); | |
177 | double qBC = x.q(EvtCyclic3::combine(iB,iC)); | |
178 | double qCA = x.q(EvtCyclic3::combine(iC,iA)); | |
179 | ||
180 | //double m0 = _prop->m0(); | |
181 | ||
182 | if(_spin == EvtSpinType::SCALAR) ret = EvtComplex(1.,0.); | |
183 | else | |
184 | if(_spin == EvtSpinType::VECTOR) { | |
185 | ||
186 | //ret = qCA - qBC - (M*M - mC*mC)*(mA*mA - mB*mB)/m0/m0; | |
187 | ret = qCA - qBC - (M*M - mC*mC)*(mA*mA - mB*mB)/qAB; | |
188 | } | |
189 | else | |
190 | if(_spin == EvtSpinType::TENSOR) { | |
191 | ||
192 | //double x1 = qBC - qCA + (M*M - mC*mC)*(mA*mA - mB*mB)/m0/m0; | |
193 | double x1 = qBC - qCA + (M*M - mC*mC)*(mA*mA - mB*mB)/qAB; | |
194 | double x2 = M*M - mC*mC; | |
195 | //double x3 = qAB - 2*M*M - 2*mC*mC + x2*x2/m0/m0; | |
196 | double x3 = qAB - 2*M*M - 2*mC*mC + x2*x2/qAB; | |
197 | double x4 = mB*mB - mA*mA; | |
198 | //double x5 = qAB - 2*mB*mB - 2*mA*mA + x4*x4/m0/m0; | |
199 | double x5 = qAB - 2*mB*mB - 2*mA*mA + x4*x4/qAB; | |
200 | ret = (x1*x1 - 1./3.*x3*x5); | |
201 | } | |
202 | else assert(0); | |
203 | } | |
204 | ||
205 | return ret; | |
206 | } | |
207 | ||
208 | ||
209 | double EvtPto3PAmp::angDep(const EvtDalitzPoint& x) const | |
210 | { | |
211 | // Angular dependece for factorizable amplitudes | |
212 | // unphysical cosines indicate we are in big trouble | |
213 | ||
214 | double cosTh = x.cosTh(_pairAng,_pairRes); | |
215 | if(fabs(cosTh) > 1.) { | |
216 | ||
217 | report(INFO,"EvtGen") << "cosTh " << cosTh << endl; | |
218 | assert(0); | |
219 | } | |
220 | ||
221 | // in units of half-spin | |
222 | ||
223 | return EvtdFunction::d(EvtSpinType::getSpin2(_spin),2*0,2*0,acos(cosTh)); | |
224 | } | |
225 |