202 {
203
204
205
209
211 double theta=angles.getHelAng(1);
212 double phi =angles.getHelAng(2);
213
214
215
217
218 int ia,ib,ic;
219
220 double prob1=0.0;
221
222 for(ia=0;ia<_nA;ia++){
223 for(ib=0;ib<_nB;ib++){
224 for(ic=0;ic<_nC;ic++){
225 _amp[ia][ib][ic]=0.0;
226 if (
abs(_lambdaB2[ib]-_lambdaC2[ic])<=_JA2) {
228 _lambdaB2[ib]-_lambdaC2[ic],theta);
229
230 _amp[ia][ib][ic]=c*_HBC[ib][ic]*
232 _lambdaC2[ic])))*dfun;
233 }
234 prob1+=
real(_amp[ia][ib][ic]*
conj(_amp[ia][ib][ic]));
235 }
236 }
237 }
238
239 setUpRotationMatrices(p,theta,phi);
240
241 applyRotationMatrices();
242
243 double prob2=0.0;
244
245 for(ia=0;ia<_nA;ia++){
246 for(ib=0;ib<_nB;ib++){
247 for(ic=0;ic<_nC;ic++){
248 prob2+=
real(_amp[ia][ib][ic]*
conj(_amp[ia][ib][ic]));
249 if (_nA==1){
250 if (_nB==1){
251 if (_nC==1){
252 amp.
vertex(_amp[ia][ib][ic]);
253 }
254 else{
255 amp.
vertex(ic,_amp[ia][ib][ic]);
256 }
257 }
258 else{
259 if (_nC==1){
260 amp.
vertex(ib,_amp[ia][ib][ic]);
261 }
262 else{
263 amp.
vertex(ib,ic,_amp[ia][ib][ic]);
264
265 }
266 }
267 }else{
268 if (_nB==1){
269 if (_nC==1){
270 amp.
vertex(ia,_amp[ia][ib][ic]);
271 }
272 else{
273 amp.
vertex(ia,ic,_amp[ia][ib][ic]);
274 }
275 }
276 else{
277 if (_nC==1){
278 amp.
vertex(ia,ib,_amp[ia][ib][ic]);
279 }
280 else{
281 amp.
vertex(ia,ib,ic,_amp[ia][ib][ic]);
282 }
283 }
284 }
285 }
286 }
287 }
288
289 if (fabs(prob1-prob2)>0.000001*prob1){
290 report(
INFO,
"EvtGen") <<
"prob1,prob2:"<<prob1<<
" "<<prob2<<endl;
291 ::abort();
292 }
293
294 return ;
295
296}
Evt3Rank3C conj(const Evt3Rank3C &t2)
EvtComplex exp(const EvtComplex &c)
ostream & report(Severity severity, const char *facility)
void vertex(const EvtComplex &)
const EvtVector4R & getP4() const
EvtParticle * getDaug(int i)
static double d(int j, int m1, int m2, double theta)