CGEM BOSS 6.6.5.h
BESIII Offline Software System
Loading...
Searching...
No Matches
EvtLundCharm Class Reference

#include <EvtLundCharm.hh>

+ Inheritance diagram for EvtLundCharm:

Public Member Functions

 EvtLundCharm ()
 
virtual ~EvtLundCharm ()
 
void getName (std::string &name)
 
EvtDecayBaseclone ()
 
void decay (EvtParticle *p)
 
std::string commandName ()
 
void command (std::string cmd)
 
void init ()
 
void initProbMax ()
 
int getTotalEvt ()
 
- Public Member Functions inherited from EvtDecayIncoherent
void makeDecay (EvtParticle *p)
 
virtual ~EvtDecayIncoherent ()
 
void setDaughterSpinDensity (int daughter)
 
int isDaughterSpinDensitySet (int daughter)
 
- Public Member Functions inherited from EvtDecayBase
double getProbMax (double prob)
 
double resetProbMax (double prob)
 
 EvtDecayBase ()
 
virtual ~EvtDecayBase ()
 
virtual bool matchingDecay (const EvtDecayBase &other) const
 
EvtId getParentId ()
 
double getBranchingFraction ()
 
void disableCheckQ ()
 
void checkQ ()
 
int getNDaug ()
 
EvtIdgetDaugs ()
 
EvtId getDaug (int i)
 
int getNArg ()
 
int getPHOTOS ()
 
void setPHOTOS ()
 
void setVerbose ()
 
void setSummary ()
 
double * getArgs ()
 
std::string * getArgsStr ()
 
double getArg (int j)
 
std::string getArgStr (int j)
 
std::string getModelName ()
 
int getDSum ()
 
int summary ()
 
int verbose ()
 
void saveDecayInfo (EvtId ipar, int ndaug, EvtId *daug, int narg, std::vector< std::string > &args, std::string name, double brfr)
 
void printSummary ()
 
void setProbMax (double prbmx)
 
void noProbMax ()
 
void checkNArg (int a1, int a2=-1, int a3=-1, int a4=-1)
 
void checkNDaug (int d1, int d2=-1)
 
void checkSpinParent (EvtSpinType::spintype sp)
 
void checkSpinDaughter (int d1, EvtSpinType::spintype sp)
 
virtual int nRealDaughters ()
 

Static Public Member Functions

static void LundcrmInit (int f)
 
- Static Public Member Functions inherited from EvtDecayBase
static void findMasses (EvtParticle *p, int ndaugs, EvtId daugs[10], double masses[10])
 
static void findMass (EvtParticle *p)
 
static double findMaxMass (EvtParticle *p)
 

Additional Inherited Members

- Protected Member Functions inherited from EvtDecayBase
bool daugsDecayedByParentModel ()
 
- Protected Attributes inherited from EvtDecayBase
bool _daugsDecayedByParentModel
 

Detailed Description

Definition at line 34 of file EvtLundCharm.hh.

Constructor & Destructor Documentation

◆ EvtLundCharm()

EvtLundCharm::EvtLundCharm ( )

Definition at line 72 of file EvtLundCharm.cc.

72{}

Referenced by clone().

◆ ~EvtLundCharm()

EvtLundCharm::~EvtLundCharm ( )
virtual

Definition at line 74 of file EvtLundCharm.cc.

74 {
75
76
77 int i;
78
79
80 //the deletion of commands is really uggly!
81
82 if (nlundcharmdecays==0) {
83 delete [] commands;
84 commands=0;
85 return;
86 }
87
88 for(i=0;i<nlundcharmdecays;i++){
89 if (lundcharmdecays[i]==this){
90 lundcharmdecays[i]=lundcharmdecays[nlundcharmdecays-1];
91 nlundcharmdecays--;
92 if (nlundcharmdecays==0) {
93 delete [] commands;
94 commands=0;
95 }
96 return;
97 }
98 }
99
100 report(ERROR,"EvtGen") << "Error in destroying LundCharm model!"<<endl;
101
102}
ostream & report(Severity severity, const char *facility)
Definition EvtReport.cc:36
@ ERROR
Definition EvtReport.hh:49

Member Function Documentation

◆ clone()

EvtDecayBase * EvtLundCharm::clone ( )
virtual

Implements EvtDecayBase.

Definition at line 111 of file EvtLundCharm.cc.

111 {
112
113 return new EvtLundCharm;
114
115}

◆ command()

void EvtLundCharm::command ( std::string cmd)
virtual

Reimplemented from EvtDecayBase.

Definition at line 154 of file EvtLundCharm.cc.

154 {
155
156 if (ncommand==lcommand){
157
158 lcommand=10+2*lcommand;
159
160 std::string* newcommands=new std::string[lcommand];
161
162 int i;
163
164 for(i=0;i<ncommand;i++){
165 newcommands[i]=commands[i];
166 }
167
168 delete [] commands;
169
170 commands=newcommands;
171
172 }
173
174 commands[ncommand]=cmd;
175
176 ncommand++;
177
178
179}

◆ commandName()

std::string EvtLundCharm::commandName ( )
virtual

Reimplemented from EvtDecayBase.

Definition at line 148 of file EvtLundCharm.cc.

148 {
149
150 return std::string("LundCharmPar");
151
152}

◆ decay()

void EvtLundCharm::decay ( EvtParticle * p)
virtual

Implements EvtDecayBase.

Definition at line 183 of file EvtLundCharm.cc.

183 {
184
185 static int iniflag=0;
186
187 static EvtId STRNG=EvtPDL::getId("string");
188
189 int istdheppar=EvtPDL::getStdHep(p->getId());
190
191/*
192 if (pycomp_(&istdheppar)==0){
193 report(ERROR,"EvtGen") << "LundCharm can not decay:"
194 <<EvtPDL::name(p->getId()).c_str()<<endl;
195 return;
196 }
197*/
198
199// std::cout<<"Lundcharm decaying "<<EvtPDL::name(p->getId())<<" mass: "<<p->getP4().mass()<<std::endl;
200
201// no eta_c(2S) in jetset74, so we don't include it in lundcharm
202 if(istdheppar != 443 && istdheppar != 100443 && istdheppar != 10441 &&istdheppar != 20443 &&istdheppar != 445 &&istdheppar != 10443 &&istdheppar != 441 &&istdheppar!= 30443){
203 std::cout<<"EvtGen: EvtLundCharm cann't not decay the particle pid= "<<istdheppar<<endl;
204 ::abort();
205 }
206
207 double mp=p->mass();
208 float xmp=mp;
209 double totEn=0;
210// std::cout<<"float xmp="<<xmp<<std::endl;
211
212 EvtVector4R p4[20];
213
214 int i,more, pflag;;
215 int ip=EvtPDL::getStdHep(p->getId());
216 int ndaugjs;
217 static int kf[100];
218 EvtId evtnumstable[100],evtnumparton[100];
219 int stableindex[100],partonindex[100];
220 int numstable;
221 int numparton;
222 static int km[100];
223 EvtId type[MAX_DAUG];
224
225 static double px[100],py[100],pz[100],e[100];
226 static int myflag;
227 if (iniflag==0) lundcrm_(&iniflag,&istdheppar,&xmp,&ndaugjs,kf,km,px,py,pz,e, &myflag);
228 LundcrmInit(0); // Allow user to set LundCharmPar in decay list
229
230 if ( p->getNDaug() != 0 ) { p->deleteDaughters(true);}
231
232 string name_parent = EvtPDL::name(p->getId());
233 double parityi=parityC::getC(name_parent);
234 int count=0;
235 do{
236 //report(INFO,"EvtGen") << "calling lundcharm " << ip<< " " << mp <<endl;
237 iniflag=iniflag+1; //to count the event number
238 lundcrm_(&iniflag,&istdheppar,&xmp,&ndaugjs,kf,km,px,py,pz,e, &myflag);
239 //-- change myflag to unsigned int
240
241 p->setGeneratorFlag(myflag);
242 // std::cout<<"EvtLundCharm::setGeneratorFalg= "<<myflag<<std::endl;
243 numstable=0;
244 numparton=0;
245 //report(INFO,"EvtGen") << "found some daughters " << ndaugjs << endl;
246 totEn=0;
247 double parityf=1;
248 for(i=0;i<ndaugjs;i++){
249 //std::cout<<"ndaugjs,kf,km,px,py,pz,e: "<<i<<", "<<km[i]<<", "<<kf[i]<<", "<<px[i]<<" ,"<<py[i]<<", "<<pz[i]<<", "<<e[i]<<std::endl; //for debugging
250 totEn +=e[i];
251 string name_daugi = EvtPDL::name( EvtPDL::evtIdFromStdHep(kf[i]) );
252 parityf = parityf*parityC::getC(name_daugi);
253
254 if (EvtPDL::evtIdFromStdHep(kf[i])==EvtId(-1,-1)) {
255 report(ERROR,"EvtGen") << "LundCharm returned particle:"<<kf[i]<<endl;
256 report(ERROR,"EvtGen") << "This can not be translated to evt number"<<endl;
257 report(ERROR,"EvtGen") << "and the decay will be rejected!"<<endl;
258 report(ERROR,"EvtGen") << "The decay was of particle:"<<ip<<endl;
259
260 }
261
262 //sort out the partons
263 if (abs(kf[i])<=6||kf[i]==21){
264 partonindex[numparton]=i;
265 evtnumparton[numparton]=EvtPDL::evtIdFromStdHep(kf[i]);
266 numparton++;
267 }
268 else{
269 stableindex[numstable]=i;
270 evtnumstable[numstable]=EvtPDL::evtIdFromStdHep(kf[i]);
271 numstable++;
272 }
273
274
275 // have to protect against negative mass^2 for massless particles
276 // i.e. neutrinos and photons.
277 // this is uggly but I need to fix it right now....
278
279 if (px[i]*px[i]+py[i]*py[i]+pz[i]*pz[i]>=e[i]*e[i]){
280
281 e[i]=sqrt(px[i]*px[i]+py[i]*py[i]+pz[i]*pz[i])+0.0000000001;
282
283 }
284
285 p4[i].set(e[i],px[i],py[i],pz[i]);
286
287
288 }
289
290 int channel=EvtDecayTable::inChannelList(p->getId(),numstable,evtnumstable);
291
292 // Test the branching fraction of lundcharm
293 // the specified decay mode is put as the 0-th channel with specifing mother particle
294 /*
295 if(ip==100443 && channel==0){
296 nevt++;
297 std::cout<<"nevt= "<<nevt<<std::endl;
298 channel=-1;
299 }
300 */
301 // std::cout<<"channel= "<<channel<<std::endl;
302 if(parityi!=0 && parityf!=0){
303 pflag=(parityi!=parityf);
304 }else{pflag=2;}
305
306 bool eck = fabs(xmp-totEn)>0.01;
307 // std::cout<<"eck= "<<eck<<", "<<fabs(xmp-totEn)<<std::endl;
308 more=(channel!=-1 || pflag ==1 || eck);
309 // more=(channel!=-1);
310
311 //---debugging
312 //std::cout<<"parentId= "<<istdheppar<<", pflag= "<<pflag<<std::endl;
313 //if(pflag==1) abort();
314
315
316 count++;
317
318 }while( more && (count<10000) );
319
320 /*
321 if(fabs(xmp-totEn)>0.01){
322 std::cout<<"Warning:LUNDCHARM generate incomplet final state, "<<mp<<" "<<totEn<<endl;
323 ::abort();
324 }
325 */
326
327 if (count>9999) {
328 report(INFO,"EvtGen") << "Too many loops in EvtLundCharm!!!"<<endl;
329 report(INFO,"EvtGen") << "Parent:"<<EvtPDL::name(getParentId()).c_str()<<endl;
330 for(i=0;i<numstable;i++){
331 report(INFO,"EvtGen") << "Daug("<<i<<")"<<EvtPDL::name(evtnumstable[i]).c_str()<<endl;
332 }
333
334 }
335
336
337
338 if (numparton==0){
339
340 p->makeDaughters(numstable,evtnumstable);
341 int ndaugFound=0;
342 for(i=0;i<numstable;i++){
343 p->getDaug(i)->init(evtnumstable[i],p4[stableindex[i]]);
344 ndaugFound++;
345 }
346 if ( ndaugFound == 0 ) {
347 report(ERROR,"EvtGen") << "Lundcharm has failed to do a decay ";
348 report(ERROR,"EvtGen") << EvtPDL::name(p->getId()).c_str() << " " << p->mass()<<endl;
349 assert(0);
350 }
351
352 fixPolarizations(p);
353
354 return ;
355
356 }
357 else{
358
359 //have partons in LUNDCHARM
360
361 EvtVector4R p4string(0.0,0.0,0.0,0.0);
362
363 for(i=0;i<numparton;i++){
364 p4string+=p4[partonindex[i]];
365 }
366
367 int nprimary=1;
368 type[0]=STRNG;
369 for(i=0;i<numstable;i++){
370 if (km[stableindex[i]]==0){
371 type[nprimary++]=evtnumstable[i];
372 }
373 }
374
375 p->makeDaughters(nprimary,type);
376
377 p->getDaug(0)->init(STRNG,p4string);
378
379 EvtVector4R p4partons[10];
380
381 for(i=0;i<numparton;i++){
382 p4partons[i]=p4[partonindex[i]];
383 }
384
385 ((EvtStringParticle*)p->getDaug(0))->initPartons(numparton,p4partons,evtnumparton);
386
387
388
389 nprimary=1;
390
391 for(i=0;i<numstable;i++){
392
393 if (km[stableindex[i]]==0){
394 p->getDaug(nprimary++)->init(evtnumstable[i],p4[stableindex[i]]);
395 }
396 }
397
398
399 int nsecond=0;
400 for(i=0;i<numstable;i++){
401 if (km[stableindex[i]]!=0){
402 type[nsecond++]=evtnumstable[i];
403 }
404 }
405
406
407 p->getDaug(0)->makeDaughters(nsecond,type);
408
409 EvtVector4R p4stringboost(p4string.get(0),-p4string.get(1),
410 -p4string.get(2),-p4string.get(3));
411
412 nsecond=0;
413 for(i=0;i<numstable;i++){
414 if (km[stableindex[i]]!=0){
415 p4[stableindex[i]]=boostTo(p4[stableindex[i]],p4stringboost);
416 p->getDaug(0)->getDaug(nsecond)->init(evtnumstable[i],p4[stableindex[i]]);
417 p->getDaug(0)->getDaug(nsecond)->setDiagonalSpinDensity();
418 p->getDaug(0)->getDaug(nsecond)->decay();
419 nsecond++;
420 }
421 }
422
423 if ( nsecond == 0 ) {
424 report(ERROR,"EvtGen") << "Jetset has failed to do a decay ";
425 report(ERROR,"EvtGen") << EvtPDL::name(p->getId()).c_str() << " " << p->mass() <<endl;
426 assert(0);
427 }
428
429 fixPolarizations(p);
430
431 return ;
432
433 }
434
435}
double abs(const EvtComplex &c)
EvtDiracSpinor boostTo(const EvtDiracSpinor &sp, const EvtVector4R p4)
void lundcrm_(int *, int *, float *, int *, int *, int *, double *, double *, double *, double *, int *)
const int MAX_DAUG
@ INFO
Definition EvtReport.hh:52
EvtId getParentId()
static int inChannelList(EvtId parent, int ndaug, EvtId *daugs)
Definition EvtId.hh:27
static void LundcrmInit(int f)
static int getStdHep(EvtId id)
Definition EvtPDL.hh:56
static EvtId evtIdFromStdHep(int stdhep)
Definition EvtPDL.cc:244
static std::string name(EvtId i)
Definition EvtPDL.hh:64
static EvtId getId(const std::string &name)
Definition EvtPDL.cc:287
void makeDaughters(int ndaug, EvtId *id)
virtual void init(EvtId part_n, const EvtVector4R &p4)=0
EvtId getId() const
void setGeneratorFlag(int flag)
void setDiagonalSpinDensity()
int getNDaug() const
EvtParticle * getDaug(int i)
void deleteDaughters(bool keepChannel=false)
double mass() const
void set(int i, double d)
static double getC(string parname)
Definition EvtParityC.cc:26
const double mp
uint32_t count(const node_t &list)
Definition node.cxx:42

◆ getName()

void EvtLundCharm::getName ( std::string & name)
virtual

Implements EvtDecayBase.

Definition at line 105 of file EvtLundCharm.cc.

105 {
106
107 model_name="LUNDCHARM";
108
109}

◆ getTotalEvt()

int EvtLundCharm::getTotalEvt ( )
inline

Definition at line 51 of file EvtLundCharm.hh.

51{return nevt;}

◆ init()

void EvtLundCharm::init ( )
virtual

Reimplemented from EvtDecayBase.

Definition at line 125 of file EvtLundCharm.cc.

125 {
126
127// checkNArg(1);
129
130 if (getParentId().isAlias()){
131
132 report(ERROR,"EvtGen") << "EvtLundCharm finds that you are decaying the"<<endl
133 << " aliased particle "
134 << EvtPDL::name(getParentId()).c_str()
135 << " with the LundCharm model"<<endl
136 << " this does not work, please modify decay table."
137 << endl;
138 report(ERROR,"EvtGen") << "Will terminate execution!"<<endl;
139 ::abort();
140
141 }
142
143 store(this);
144
145}
static void readParityC()
Definition EvtParityC.cc:6

◆ initProbMax()

void EvtLundCharm::initProbMax ( )
virtual

Reimplemented from EvtDecayBase.

Definition at line 118 of file EvtLundCharm.cc.

118 {
119
120 noProbMax();
121
122}

◆ LundcrmInit()

void EvtLundCharm::LundcrmInit ( int f)
static

Definition at line 499 of file EvtLundCharm.cc.

499 {
500
501 static int first=1;
502 if (first){
503
504 first=0;
505 for(int i=0;i<ncommand;i++)
506 lugive_(commands[i].c_str(),strlen(commands[i].c_str()));
507 }
508
509
510}
void lugive_(const char *cnfgstr, int length)
char * c_str(Index i)
Index first(Pair i)

Referenced by decay().


The documentation for this class was generated from the following files: