BOSS 8.0.0
BESIII Offline Software System
Loading...
Searching...
No Matches
EvtPhokhara_Lambda.cc
Go to the documentation of this file.
1//--------------------------------------------------------------------------
2//
3// Environment:
4// This software is part of models developed at BES collaboration
5// based on the EvtGen framework. If you use all or part
6// of it, please give an appropriate acknowledgement.
7//
8// Copyright Information: See EvtGen/BesCopyright
9// Copyright (A) 2006 Ping Rong-Gang @IHEP
10//
11// Module: EvtPhokhara.cc
12// the necessary file: phokharar.F
13// fist.inc,gen.inc mix.inc stdhep.inc
14// Description: Modified Lund model at tau-charm energy level, see
15// PHYSICAL REVIEW D, VOLUME 62, 034003
16// Modification history:
17//
18// Ping R.-G. Jan.25, 2010 Module created
19// The random engine RLU0 is unified with RLU for BesEvtGen
20//------------------------------------------------------------------------
21//
22#include "EvtPhokhara_Lambda.hh"
23#include "../EvtGenBase/EvtDecayTable.hh"
24#include "../EvtGenBase/EvtId.hh"
25#include "../EvtGenBase/EvtPDL.hh"
26#include "../EvtGenBase/EvtParticle.hh"
27#include "../EvtGenBase/EvtPatches.hh"
28#include "../EvtGenBase/EvtRandom.hh" //new added for using the same random seed with BesRndGenSvc 2022.3.17
29#include "../EvtGenBase/EvtReport.hh"
30#include "../EvtGenBase/EvtStringParticle.hh"
31#include "EvtPhokharaDef.hh"
32#include <fstream>
33#include <iomanip>
34#include <iostream>
35#include <stdio.h>
36#include <stdlib.h>
37#include <string.h>
38#include <string>
39#include <unistd.h>
40
41#include "BesRndmGenSvc/IBesRndmGenSvc.h"
42#include "CLHEP/Random/RandomEngine.h"
43#include "GeneratorObject/McGenEvent.h"
44#include "cfortran/cfortran.h"
45
46using namespace std;
47using namespace CLHEP;
48
49using std::endl;
50using std::fstream;
51using std::ios;
52using std::ofstream;
53using std::resetiosflags;
54using std::setiosflags;
55using std::setw;
56
57int EvtPhokhara_Lambda::nevtgen = 0;
58int EvtPhokhara_Lambda::nphokharadecays = 0;
59EvtDecayBasePtr* EvtPhokhara_Lambda::phokharadecays = 0;
60int EvtPhokhara_Lambda::ntable = 0;
61
62int EvtPhokhara_Lambda::ncommand = 0;
63int EvtPhokhara_Lambda::lcommand = 0;
64std::string* EvtPhokhara_Lambda::commands = 0;
65int EvtPhokhara_Lambda::nevt = 0;
66
67EvtPhokhara_Lambda::EvtPhokhara_Lambda() {}
68EvtPhokhara_Lambda::~EvtPhokhara_Lambda() {
69 int i;
70 // the deletion of commands is really uggly!
71
72 if ( nphokharadecays == 0 )
73 {
74 delete[] commands;
75 commands = 0;
76 return;
77 }
78
79 for ( i = 0; i < nphokharadecays; i++ )
80 {
81 if ( phokharadecays[i] == this )
82 {
83 phokharadecays[i] = phokharadecays[nphokharadecays - 1];
84 nphokharadecays--;
85 if ( nphokharadecays == 0 )
86 {
87 delete[] commands;
88 commands = 0;
89 }
90 return;
91 }
92 }
93
94 report( ERROR, "EvtGen" ) << "Error in destroying Phokhara model!" << endl;
95}
96
97void EvtPhokhara_Lambda::getName( std::string& model_name ) { model_name = "PHOKHARA_LAMBDA"; }
98
99EvtDecayBase* EvtPhokhara_Lambda::clone() { return new EvtPhokhara_Lambda; }
100
101void EvtPhokhara_Lambda::initProbMax() { noProbMax(); }
102
103void EvtPhokhara_Lambda::init_mode( EvtParticle* p ) {
104 m_pion = 9;
105 // mu+mu-(0),pi+pi-(1),2pi0pi+pi-(2),
106 // 2pi+2pi-(3),ppbar(4),nnbar(5),
107 // K+K-(6),K0K0bar(7),pi+pi-pi0(8),
108 // Lamb Lambbar->pi-pi+ppbar(9)
109#include "Phokhara_init_mode.txt"
110}
111
112void EvtPhokhara_Lambda::init() {
113 checkNArg( 0 );
114
115 std::string locvp = getenv( "BESEVTGENROOT" );
116 system( "cat $BESEVTGENROOT/share/phokhara_10.0.param>phokhara_10.0.param" );
117 system( "cat $BESEVTGENROOT/share/phokhara_10.0.fferr>phokhara_10.0.fferr" );
118 system( "cat $BESEVTGENROOT/share/phokhara_10.0.ffwarn>phokhara_10.0.ffwarn" );
119
120 if ( getParentId().isAlias() )
121 {
122
123 report( ERROR, "EvtGen" ) << "EvtPhokhara finds that you are decaying the" << endl
124 << " aliased particle " << EvtPDL::name( getParentId() ).c_str()
125 << " with the Phokhara model" << endl
126 << " this does not work, please modify decay table." << endl;
127 report( ERROR, "EvtGen" ) << "Will terminate execution!" << endl;
128 ::abort();
129 }
130
131 store( this );
132}
133
134std::string EvtPhokhara_Lambda::commandName() { return std::string( "PhokharaPar" ); }
135
136void EvtPhokhara_Lambda::command( std::string cmd ) {
137
138 if ( ncommand == lcommand )
139 {
140
141 lcommand = 10 + 2 * lcommand;
142
143 std::string* newcommands = new std::string[lcommand];
144
145 int i;
146
147 for ( i = 0; i < ncommand; i++ ) { newcommands[i] = commands[i]; }
148
149 delete[] commands;
150
151 commands = newcommands;
152 }
153
154 commands[ncommand] = cmd;
155
156 ncommand++;
157}
158
159void EvtPhokhara_Lambda::decay( EvtParticle* p ) {
160 bool debug = 0;
161 EvtId myvpho = EvtPDL::getId( "vpho" );
162 if ( p->getId() != myvpho )
163 {
164 std::cout << "Parent particle is required to be vpho for Phokhara model" << std::endl;
165 abort();
166 }
167 if ( nevtgen == 0 )
168 {
169 init_mode( p );
170 std::cout << "PHOKHARA : Lambda anti-Lambda mode " << std::endl;
171 }
172 else
173 {
174 init_evt( p );
175 if ( debug ) cout << "init_evt works!" << endl;
176 }
177 // else{init_mode(p);}
178 if ( debug )
179 {
180 std::cout << "flags "
181 ":ph0,nlo,pion,fsr,fsrnlo,ivac,FF_pion,f0_model,FF_kaon,narr_res,FF_pp,chi_"
182 "sw,chi_pion,FF_Pgg,nlo2 "
183 << std::endl;
184 std::cout << "= " << flags_.ph0 << "," << flags_.nlo << "," << flags_.pion << ","
185 << flags_.fsr << "," << flags_.fsrnlo << "," << flags_.ivac << ","
186 << flags_.FF_pion << "," << flags_.f0_model << "," << flags_.FF_kaon << ","
187 << flags_.narr_res << "," << flags_.FF_pp << "," << flags_.chi_sw << ","
188 << flags_.be_r << "," << flags_.FF_Pgg << "," << flags_.nlo2 << std::endl;
189 std::cout << "ctes: Sp = " << ctes_.Sp << std::endl;
190 std::cout << "cuts : w,q2min,q2_min_c,gmin,phot1cut,phot2cut,pi1cut,pi2cut" << std::endl
191 << "= " << cuts_.w << "," << cuts_.q2min << "," << cuts_.q2_min_c << ","
192 << cuts_.gmin << "," << cuts_.phot1cut << "," << cuts_.phot2cut << ","
193 << cuts_.pi1cut << "," << cuts_.pi2cut << std::endl;
194 }
195
196 int istdheppar = EvtPDL::getStdHep( p->getId() );
197 int ntrials = 0;
198 int tr_old[3];
199 tr_old[0] = (int)maxima_.tr[0];
200 tr_old[1] = (int)maxima_.tr[1];
201 tr_old[2] = (int)maxima_.tr[2];
202
203 while ( ntrials < 1000000 )
204 {
205 ievent++;
206 RANLXDF( Ar_r, 1 );
207 Ar[1] = Ar_r[0];
208
209 if ( Ar[1] <=
210 ( maxima_.Mmax[0] / ( maxima_.Mmax[0] + maxima_.Mmax[1] + maxima_.Mmax[2] ) ) )
211 {
212 maxima_.count[0] = maxima_.count[0] + 1.0;
213 GEN_0PH( 2, qqmin, ctes_.Sp, cos3min, cos3max );
214 }
215 else if ( Ar[1] <= ( ( maxima_.Mmax[0] + maxima_.Mmax[1] ) /
216 ( maxima_.Mmax[0] + maxima_.Mmax[1] + maxima_.Mmax[2] ) ) )
217 {
218 maxima_.count[1] = maxima_.count[1] + 1.0;
219 GEN_1PH( 2, qqmin, qqmax, cos1min, cos1max, cos3min, cos3max );
220 }
221 else
222 {
223 maxima_.count[2] = maxima_.count[2] + 1.0;
224 GEN_2PH( 2, qqmin, cos1min, cos1max, cos2min, cos2max, cos3min, cos3max );
225 }
226
227 if ( ( (int)maxima_.tr[0] + (int)maxima_.tr[1] + (int)maxima_.tr[2] ) >
228 ( tr_old[0] + tr_old[1] + tr_old[2] ) ) // event accepted after cuts
229 { goto storedEvents; }
230 ntrials++;
231 }
232 std::cout << "FATAL: Could not satisfy cuts after " << ntrials << "trials. Terminate."
233 << std::endl;
234 //----
235storedEvents:
236 int more = 0;
237 int numstable = 0;
238 int numparton = 0;
239 EvtId evtnumstable[100]; //
240 EvtVector4R p4[20];
241
242 // except ISR photos, products depending on channel
243 if ( flags_.pion == 0 )
244 { // mu+ mu-
245 // mu+
246 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( -13 );
247 p4[numstable].set( ctes_.momenta[0][5], ctes_.momenta[1][5], ctes_.momenta[2][5],
248 ctes_.momenta[3][5] );
249 numstable++;
250 // mu -
251 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 13 );
252 p4[numstable].set( ctes_.momenta[0][6], ctes_.momenta[1][6], ctes_.momenta[2][6],
253 ctes_.momenta[3][6] );
254 numstable++;
255 }
256 if ( flags_.pion == 1 )
257 { // pi+ pi-
258 // pi+
259 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 211 );
260 p4[numstable].set( ctes_.momenta[0][5], ctes_.momenta[1][5], ctes_.momenta[2][5],
261 ctes_.momenta[3][5] );
262 numstable++;
263 // pi -
264 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( -211 );
265 p4[numstable].set( ctes_.momenta[0][6], ctes_.momenta[1][6], ctes_.momenta[2][6],
266 ctes_.momenta[3][6] );
267 numstable++;
268 }
269 if ( flags_.pion == 2 )
270 { // pi+ pi-2pi0
271 // pi0
272 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 111 );
273 p4[numstable].set( ctes_.momenta[0][5], ctes_.momenta[1][5], ctes_.momenta[2][5],
274 ctes_.momenta[3][5] );
275 numstable++;
276 // pi0
277 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 111 );
278 p4[numstable].set( ctes_.momenta[0][6], ctes_.momenta[1][6], ctes_.momenta[2][6],
279 ctes_.momenta[3][6] );
280 numstable++;
281 // pi-
282 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( -211 );
283 p4[numstable].set( ctes_.momenta[0][7], ctes_.momenta[1][7], ctes_.momenta[2][7],
284 ctes_.momenta[3][7] );
285 numstable++;
286 // pi +
287 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 211 );
288 p4[numstable].set( ctes_.momenta[0][8], ctes_.momenta[1][8], ctes_.momenta[2][8],
289 ctes_.momenta[3][8] );
290 numstable++;
291 }
292 if ( flags_.pion == 3 )
293 { // 2(pi+ pi-)
294 // pi+
295 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 211 );
296 p4[numstable].set( ctes_.momenta[0][5], ctes_.momenta[1][5], ctes_.momenta[2][5],
297 ctes_.momenta[3][5] );
298 numstable++;
299 // pi-
300 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( -211 );
301 p4[numstable].set( ctes_.momenta[0][6], ctes_.momenta[1][6], ctes_.momenta[2][6],
302 ctes_.momenta[3][6] );
303 numstable++;
304 // pi+
305 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( -211 );
306 p4[numstable].set( ctes_.momenta[0][7], ctes_.momenta[1][7], ctes_.momenta[2][7],
307 ctes_.momenta[3][7] );
308 numstable++;
309 // pi -
310 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 211 );
311 p4[numstable].set( ctes_.momenta[0][8], ctes_.momenta[1][8], ctes_.momenta[2][8],
312 ctes_.momenta[3][8] );
313 numstable++;
314 }
315 if ( flags_.pion == 4 )
316 { // ppbar
317 // pbar
318 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( -2212 );
319 p4[numstable].set( ctes_.momenta[0][5], ctes_.momenta[1][5], ctes_.momenta[2][5],
320 ctes_.momenta[3][5] );
321 numstable++;
322 // p
323 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 2212 );
324 p4[numstable].set( ctes_.momenta[0][6], ctes_.momenta[1][6], ctes_.momenta[2][6],
325 ctes_.momenta[3][6] );
326 numstable++;
327 }
328 if ( flags_.pion == 5 )
329 { // nnbar
330 // pbar
331 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( -2112 );
332 p4[numstable].set( ctes_.momenta[0][5], ctes_.momenta[1][5], ctes_.momenta[2][5],
333 ctes_.momenta[3][5] );
334 numstable++;
335 // p
336 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 2112 );
337 p4[numstable].set( ctes_.momenta[0][6], ctes_.momenta[1][6], ctes_.momenta[2][6],
338 ctes_.momenta[3][6] );
339 numstable++;
340 }
341 if ( flags_.pion == 6 )
342 { // K+ K-
343 // K+
344 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 321 );
345 p4[numstable].set( ctes_.momenta[0][5], ctes_.momenta[1][5], ctes_.momenta[2][5],
346 ctes_.momenta[3][5] );
347 numstable++;
348 // K -
349 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( -321 );
350 p4[numstable].set( ctes_.momenta[0][6], ctes_.momenta[1][6], ctes_.momenta[2][6],
351 ctes_.momenta[3][6] );
352 numstable++;
353 }
354 if ( flags_.pion == 7 )
355 { // K0K0bar
356 // Kbar
357 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 311 );
358 p4[numstable].set( ctes_.momenta[0][5], ctes_.momenta[1][5], ctes_.momenta[2][5],
359 ctes_.momenta[3][5] );
360 numstable++;
361 // K0
362 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( -311 );
363 p4[numstable].set( ctes_.momenta[0][6], ctes_.momenta[1][6], ctes_.momenta[2][6],
364 ctes_.momenta[3][6] );
365 numstable++;
366 }
367 if ( flags_.pion == 8 )
368 { // pi+ pi-pi0
369 // pi+
370 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 211 );
371 p4[numstable].set( ctes_.momenta[0][5], ctes_.momenta[1][5], ctes_.momenta[2][5],
372 ctes_.momenta[3][5] );
373 numstable++;
374 // pi-
375 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( -211 );
376 p4[numstable].set( ctes_.momenta[0][6], ctes_.momenta[1][6], ctes_.momenta[2][6],
377 ctes_.momenta[3][6] );
378 numstable++;
379 // pi0
380 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 111 );
381 p4[numstable].set( ctes_.momenta[0][7], ctes_.momenta[1][7], ctes_.momenta[2][7],
382 ctes_.momenta[3][7] );
383 numstable++;
384 }
385 /*if (flags_.pion == 9) { //Lambda Lambdabar-> pi+ pi- ppbar
386 // pi+
387 evtnumstable[numstable]=EvtPDL::evtIdFromStdHep(211);
388 p4[numstable].set(ctes_.momenta[0][7],ctes_.momenta[1][7], ctes_.momenta[2][7],
389 ctes_.momenta[3][7]); numstable++;
390 // pbar
391 evtnumstable[numstable]=EvtPDL::evtIdFromStdHep(-2212);
392 p4[numstable].set(ctes_.momenta[0][8],ctes_.momenta[1][8], ctes_.momenta[2][8],
393 ctes_.momenta[3][8]); numstable++;
394 // pi-
395 evtnumstable[numstable]=EvtPDL::evtIdFromStdHep(-211);
396 p4[numstable].set(ctes_.momenta[0][9],ctes_.momenta[1][9], ctes_.momenta[2][9],
397 ctes_.momenta[3][9]); numstable++;
398 // p
399 evtnumstable[numstable]=EvtPDL::evtIdFromStdHep(2212);
400 p4[numstable].set(ctes_.momenta[0][10],ctes_.momenta[1][10], ctes_.momenta[2][10],
401 ctes_.momenta[3][10]); numstable++;
402 }*/
403 int LambdaDaus = 0;
404 EvtId LambdaMode[100]; //
405 EvtVector4R Lambdap4[20];
406 EvtVector4R ZERO( 0, 0, 0, 0 );
407 EvtVector4R tmp( 0, 0, 0, 0 );
408
409 if ( flags_.pion == 9 )
410 { // Lambda Lambdabar-> pi+ pi- ppbar
411 // anti-Lambda0
412 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( -3122 );
413 p4[numstable].set( ctes_.momenta[0][5], ctes_.momenta[1][5], ctes_.momenta[2][5],
414 ctes_.momenta[3][5] );
415 if ( debug )
416 std::cout << "Phokhara_Lambda: anti-Lambda0 p4[numstable] = " << p4[numstable]
417 << std::endl;
418 numstable++;
419 // Lambda0
420 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 3122 );
421 p4[numstable].set( ctes_.momenta[0][6], ctes_.momenta[1][6], ctes_.momenta[2][6],
422 ctes_.momenta[3][6] );
423 if ( debug )
424 std::cout << "Phokhara_Lambda: Lambda0 p4[numstable] = " << p4[numstable] << std::endl;
425 numstable++;
426 /////Lambda Decay
427 // pi+
428 LambdaMode[LambdaDaus] = EvtPDL::evtIdFromStdHep( 211 );
429 Lambdap4[LambdaDaus].set( ctes_.momenta[0][7], ctes_.momenta[1][7], ctes_.momenta[2][7],
430 ctes_.momenta[3][7] );
431 if ( debug )
432 std::cout << "Phokhara_Lambda: Pi+ p4[numstable] = " << Lambdap4[LambdaDaus]
433 << std::endl;
434 ////////boost pi+ to antiLambda0 frame
435 // tmp = ZERO-p4[0];wrong calculation
436 tmp.set( p4[0].get( 0 ), -p4[0].get( 1 ), -p4[0].get( 2 ), -p4[0].get( 3 ) );
437 Lambdap4[LambdaDaus] = boostTo( Lambdap4[LambdaDaus], tmp );
438 if ( debug )
439 std::cout << "Phokhara_Lambda:Boosted Pi+ p4[numstable] = " << Lambdap4[LambdaDaus]
440 << std::endl;
441 LambdaDaus++;
442 // pbar
443 LambdaMode[LambdaDaus] = EvtPDL::evtIdFromStdHep( -2212 );
444 Lambdap4[LambdaDaus].set( ctes_.momenta[0][8], ctes_.momenta[1][8], ctes_.momenta[2][8],
445 ctes_.momenta[3][8] );
446 if ( debug )
447 std::cout << "Phokhara_Lambda: pbar p4[numstable] = " << Lambdap4[LambdaDaus]
448 << std::endl;
449 ////////boost pbar to antiLambda0 frame
450 // tmp = ZERO-p4[0];
451 // tmp = p4[0];
452 tmp.set( p4[0].get( 0 ), -p4[0].get( 1 ), -p4[0].get( 2 ), -p4[0].get( 3 ) );
453 Lambdap4[LambdaDaus] = boostTo( Lambdap4[LambdaDaus], tmp );
454 if ( debug )
455 std::cout << "Phokhara_Lambda:Boosted pbar p4[numstable] = " << Lambdap4[LambdaDaus]
456 << std::endl;
457 LambdaDaus++;
458 // pi-
459 LambdaMode[LambdaDaus] = EvtPDL::evtIdFromStdHep( -211 );
460 Lambdap4[LambdaDaus].set( ctes_.momenta[0][9], ctes_.momenta[1][9], ctes_.momenta[2][9],
461 ctes_.momenta[3][9] );
462 if ( debug )
463 std::cout << "Phokhara_Lambda: pi- p4[numstable] = " << Lambdap4[LambdaDaus]
464 << std::endl;
465 ////////boost pi- to Lambda0 frame
466 // tmp = ZERO-p4[1];
467 // tmp = p4[1];
468 tmp.set( p4[1].get( 0 ), -p4[1].get( 1 ), -p4[1].get( 2 ), -p4[1].get( 3 ) );
469 Lambdap4[LambdaDaus] = boostTo( Lambdap4[LambdaDaus], tmp );
470 if ( debug )
471 std::cout << "Phokhara_Lambda:Boosted pi- p4[numstable] = " << Lambdap4[LambdaDaus]
472 << std::endl;
473 LambdaDaus++;
474 // p
475 LambdaMode[LambdaDaus] = EvtPDL::evtIdFromStdHep( 2212 );
476 Lambdap4[LambdaDaus].set( ctes_.momenta[0][10], ctes_.momenta[1][10], ctes_.momenta[2][10],
477 ctes_.momenta[3][10] );
478 if ( debug )
479 std::cout << "Phokhara_Lambda: p p4[numstable] = " << Lambdap4[LambdaDaus] << std::endl;
480 ////////boost p+ to Lambda0 frame
481 // tmp = ZERO-p4[1];
482 // tmp = p4[1];
483 tmp.set( p4[1].get( 0 ), -p4[1].get( 1 ), -p4[1].get( 2 ), -p4[1].get( 3 ) );
484 Lambdap4[LambdaDaus] = boostTo( Lambdap4[LambdaDaus], tmp );
485 if ( debug )
486 std::cout << "Phokhara_Lambda:Boosted p p4[numstable] = " << Lambdap4[LambdaDaus]
487 << std::endl;
488 LambdaDaus++;
489 }
490
491 // ISR gamma
492 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 22 );
493 p4[numstable].set( ctes_.momenta[0][2], ctes_.momenta[1][2], ctes_.momenta[2][2],
494 ctes_.momenta[3][2] );
495 if ( debug )
496 std::cout << "Phokhara_Lambda: first gamma p4[numstable] = " << p4[numstable] << std::endl;
497 numstable++;
498 if ( ctes_.momenta[0][3] != 0 ) // second photon exists
499 {
500 evtnumstable[numstable] = EvtPDL::evtIdFromStdHep( 22 );
501 p4[numstable].set( ctes_.momenta[0][3], ctes_.momenta[1][3], ctes_.momenta[2][3],
502 ctes_.momenta[3][3] );
503 if ( debug )
504 std::cout << "Phokhara_Lambda: second gamma p4[numstable] = " << p4[numstable]
505 << std::endl;
506 numstable++;
507 }
508
509 int channel = EvtDecayTable::inChannelList( p->getId(), numstable, evtnumstable );
510 more = ( channel != -1 );
511 if ( more )
512 {
513 std::cout << "EvtPhokhara_Lambda:Existence of mode " << channel
514 << " in exclusive decay list has the same final state as this one" << std::endl;
515 abort();
516 }
517
518 p->makeDaughters( numstable, evtnumstable );
519 // double weight = p->initializePhaseSpace(getNDaug(),getDaugs());
520 if ( EvtParticle::_NextLevelDauNum == 0 )
521 {
522 // std::cout<<"NextLevelDauNum==0,LambdaDaus = "<<LambdaDaus<<std::endl;
523 EvtParticle::_NextLevelDauNum = LambdaDaus;
524 for ( int i = 0; i < LambdaDaus; i++ )
525 {
526 EvtParticle::_NextLevelId[i] = LambdaMode[i];
527 EvtParticle::_NextLevelP4[i] = Lambdap4[i];
528 }
529 }
530 int ndaugFound = 0;
531 EvtVector4R SUMP4( 0, 0, 0, 0 );
532 for ( int i = 0; i < numstable; i++ )
533 {
534 p->getDaug( i )->init( evtnumstable[i], p4[i] );
535 ndaugFound++;
536 }
537 if ( ndaugFound == 0 )
538 {
539 report( ERROR, "EvtGen" ) << "Phokhara has failed to do a decay ";
540 report( ERROR, "EvtGen" ) << EvtPDL::name( p->getId() ).c_str() << " " << p->mass()
541 << endl;
542 assert( 0 );
543 }
544 if ( debug )
545 std::cout << "EvtPhokhara_Lambda SUMMARY: part p4" << p->getP4Lab() << std::endl;
546 if ( debug )
547 std::cout << "EvtPhokhara_Lambda SUMMARY: Daug0 p4" << p->getDaug( 0 )->getP4Lab()
548 << std::endl;
549 if ( debug )
550 std::cout << "EvtPhokhara_Lambda SUMMARY: Daug1 p4" << p->getDaug( 1 )->getP4Lab()
551 << std::endl;
552 if ( debug )
553 std::cout << "EvtPhokhara_Lambda SUMMARY: Daug2 p4" << p->getDaug( 2 )->getP4Lab()
554 << std::endl;
555
556 nevtgen++;
557 return;
558}
559
560void EvtPhokhara_Lambda::store( EvtDecayBase* jsdecay ) {
561
562 if ( nphokharadecays == ntable )
563 {
564
565 EvtDecayBasePtr* newphokharadecays = new EvtDecayBasePtr[2 * ntable + 10];
566 int i;
567 for ( i = 0; i < ntable; i++ ) { newphokharadecays[i] = phokharadecays[i]; }
568 ntable = 2 * ntable + 10;
569 delete[] phokharadecays;
570 phokharadecays = newphokharadecays;
571 }
572
573 phokharadecays[nphokharadecays++] = jsdecay;
574}
575
576void EvtPhokhara_Lambda::PhokharaInit( int dummy ) {
577 static int first = 1;
578 if ( first )
579 {
580
581 first = 0;
582 // for(int i=0;i<ncommand;i++)
583 // lugive0_(commands[i].c_str(),strlen(commands[i].c_str()));
584 }
585}
586
587void EvtPhokhara_Lambda::init_evt( EvtParticle* p ) {
588 m_pion = 9;
589 // mu+mu-(0),pi+pi-(1),2pi0pi+pi-(2),
590 // 2pi+2pi-(3),ppbar(4),nnbar(5),
591 // K+K-(6),K0K0bar(7),pi+pi-pi0(8),
592 // Lamb Lambbar->pi-pi+ppbar(9)
593#include "Phokhara_init_evt.txt"
594}
EvtDecayBasePtr
EvtDecayBase * EvtDecayBasePtr
Definition EvtDecayTable.hh:32
boostTo
EvtDiracSpinor boostTo(const EvtDiracSpinor &sp, const EvtVector4R p4)
Definition EvtDiracSpinor.cc:67
ctes_
struct @053254170326070136226344307237142165176240334330 ctes_
maxima_
struct @027003056066344010031101102046265032310161340072 maxima_
GEN_1PH
#define GEN_1PH(I, QQMIN, QQMAX, COS1MIN, COS1MAX, COS3MIN, COS3MAX)
Definition EvtPhokharaDef.hh:112
RANLXDF
#define RANLXDF(AR, VAL)
Definition EvtPhokharaDef.hh:109
GEN_0PH
#define GEN_0PH(I, QQMIN, SP, COS3MIN, COS3MAX)
GEN_2PH
#define GEN_2PH(I, QQMIN, COS1MIN, COS1MAX, COS2MIN, COS2MAX, COS3MIN, COS3MAX)
flags_
struct @366025114004024134344043225357106161032107165361 flags_
cuts_
struct @276363222274272223263152166363125067340201005217 cuts_
EvtPhokhara_Lambda.hh
report
ostream & report(Severity severity, const char *facility)
Definition EvtReport.cc:34
ERROR
@ ERROR
Definition EvtReport.hh:49
EvtDecayBase::getParentId
EvtId getParentId()
Definition EvtDecayBase.hh:63
EvtDecayBase::checkNArg
void checkNArg(int a1, int a2=-1, int a3=-1, int a4=-1)
Definition EvtDecayBase.cc:449
EvtDecayTable::inChannelList
static int inChannelList(EvtId parent, int ndaug, EvtId *daugs)
Definition EvtDecayTable.cc:831
EvtId
Definition EvtId.hh:27
EvtPDL::getStdHep
static int getStdHep(EvtId id)
Definition EvtPDL.hh:61
EvtPDL::evtIdFromStdHep
static EvtId evtIdFromStdHep(int stdhep)
Definition EvtPDL.cc:232
EvtPDL::name
static std::string name(EvtId i)
Definition EvtPDL.hh:70
EvtPDL::getId
static EvtId getId(const std::string &name)
Definition EvtPDL.cc:272
EvtParticle::makeDaughters
void makeDaughters(int ndaug, EvtId *id)
Definition EvtParticle.cc:1150
EvtParticle::init
virtual void init(EvtId part_n, const EvtVector4R &p4)=0
EvtParticle::getP4Lab
EvtVector4R getP4Lab()
Definition EvtParticle.cc:694
EvtParticle::getId
EvtId getId() const
Definition EvtParticle.cc:100
EvtParticle::_NextLevelId
static EvtId _NextLevelId[20]
Definition EvtParticle.hh:144
EvtParticle::_NextLevelDauNum
static int _NextLevelDauNum
Definition EvtParticle.hh:143
EvtParticle::getDaug
EvtParticle * getDaug(int i)
Definition EvtParticle.cc:80
EvtParticle::_NextLevelP4
static EvtVector4R _NextLevelP4[20]
Definition EvtParticle.hh:145
EvtParticle::mass
double mass() const
Definition EvtParticle.cc:114
EvtPhokhara_Lambda::clone
EvtDecayBase * clone()
Definition EvtPhokhara_Lambda.cc:99
EvtPhokhara_Lambda::init
void init()
Definition EvtPhokhara_Lambda.cc:112
EvtPhokhara_Lambda::PhokharaInit
void PhokharaInit(int dummy)
Definition EvtPhokhara_Lambda.cc:576
EvtPhokhara_Lambda::commandName
std::string commandName()
Definition EvtPhokhara_Lambda.cc:134
EvtPhokhara_Lambda::getName
void getName(std::string &name)
Definition EvtPhokhara_Lambda.cc:97
EvtPhokhara_Lambda::init_evt
void init_evt(EvtParticle *p)
Definition EvtPhokhara_Lambda.cc:587
EvtPhokhara_Lambda::~EvtPhokhara_Lambda
virtual ~EvtPhokhara_Lambda()
Definition EvtPhokhara_Lambda.cc:68
EvtPhokhara_Lambda::initProbMax
void initProbMax()
Definition EvtPhokhara_Lambda.cc:101
EvtPhokhara_Lambda::EvtPhokhara_Lambda
EvtPhokhara_Lambda()
Definition EvtPhokhara_Lambda.cc:67
EvtPhokhara_Lambda::decay
void decay(EvtParticle *p)
Definition EvtPhokhara_Lambda.cc:159
EvtPhokhara_Lambda::command
void command(std::string cmd)
Definition EvtPhokhara_Lambda.cc:136
EvtPhokhara_Lambda::init_mode
void init_mode(EvtParticle *p)
Definition EvtPhokhara_Lambda.cc:103
EvtVector4R::set
void set(int i, double d)
Definition EvtVector4R.hh:173