00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 #include "EvtGenBase/EvtPatches.hh"
00022
00023 #include <iostream>
00024 #include <iomanip>
00025 #include <fstream>
00026 #include <ctype.h>
00027 #include <stdlib.h>
00028 #include <string.h>
00029 #include "EvtGenBase/EvtParticle.hh"
00030 #include "EvtGenBase/EvtRandom.hh"
00031 #include "EvtGenBase/EvtDecayTable.hh"
00032 #include "EvtGenBase/EvtPDL.hh"
00033 #include "EvtGenBase/EvtDecayParm.hh"
00034 #include "EvtGenBase/EvtSymTable.hh"
00035 #include "EvtGenBase/EvtDecayBase.hh"
00036 #include "EvtGenBase/EvtModel.hh"
00037 #include "EvtGenBase/EvtParticleDecayList.hh"
00038 #include "EvtGenBase/EvtParser.hh"
00039 #include "EvtGenBase/EvtReport.hh"
00040 #include "EvtGenBase/EvtModelAlias.hh"
00041 #include "EvtGenBase/EvtRadCorr.hh"
00042 #include <vector>
00043 using std::endl;
00044 using std::fstream;
00045 using std::ifstream;
00046
00047 static std::vector<EvtParticleDecayList> decaytable;
00048
00049 int EvtDecayTable::getNMode(int ipar){
00050 return decaytable[ipar].getNMode();
00051 }
00052
00053 EvtDecayBase* getDecay(int ipar, int imode){
00054 return decaytable[ipar].getDecayModel(imode);
00055 }
00056
00057 EvtDecayBase* EvtDecayTable::gettheDecay(EvtId parent, int imode){
00058 int ipar=parent.getAlias();
00059 EvtDecayBase* thedecaymodel=decaytable[ipar].getDecay(imode).getDecayModel();
00060 return thedecaymodel;
00061 }
00062
00063 void EvtDecayTable::printSummary(){
00064
00065 int i;
00066
00067 for(i=0;i<EvtPDL::entries();i++){
00068 decaytable[i].printSummary();
00069 }
00070
00071 }
00072
00073 EvtDecayBase* EvtDecayTable::getDecayFunc(EvtParticle *p){
00074 int partnum;
00075
00076 partnum=p->getId().getAlias();
00077
00078 if ( decaytable[partnum].getNMode()==0 ) return 0;
00079 return decaytable[partnum].getDecayModel(p);
00080
00081 }
00082
00083 void EvtDecayTable::readDecayFile(const std::string dec_name){
00084
00085 if ( decaytable.size() < EvtPDL::entries() ) decaytable.resize(EvtPDL::entries());
00086 EvtModel &modelist=EvtModel::instance();
00087 int i;
00088
00089 report(INFO,"EvtGen") << "In readDecayFile, reading:"<<dec_name.c_str()<<endl;
00090
00091 ifstream fin;
00092
00093 fin.open(dec_name.c_str());
00094 if (!fin) {
00095 report(ERROR,"EvtGen") << "Could not open "<<dec_name.c_str()<<endl;
00096 }
00097 fin.close();
00098
00099 EvtParser parser;
00100 parser.Read(dec_name);
00101
00102 int itok;
00103
00104 int hasend=0;
00105
00106 std::string token;
00107
00108 for(itok=0;itok<parser.getNToken();itok++){
00109
00110 token=parser.getToken(itok);
00111
00112 if (token=="End") hasend=1;
00113
00114 }
00115
00116 if (!hasend){
00117 report(ERROR,"EvtGen") << "Could not find an 'End' in "<<dec_name.c_str()<<endl;
00118 report(ERROR,"EvtGen") << "Will terminate execution."<<endl;
00119 ::abort();
00120 }
00121
00122
00123
00124 std::string model,parent,sdaug;
00125
00126 EvtId ipar;
00127
00128 int n_daugh;
00129 EvtId daught[MAX_DAUG];
00130 double brfr;
00131
00132 int itoken=0;
00133
00134 std::vector<EvtModelAlias> modelAliasList;
00135
00136
00137 do{
00138
00139 token=parser.getToken(itoken++);
00140
00141
00142 if (token=="noPhotos"){
00143 EvtRadCorr::setNeverRadCorr();
00144 report(INFO,"EvtGen")
00145 << "As requested, PHOTOS will be turned off."<<endl;
00146 }
00147 else if (token=="yesPhotos"){
00148 EvtRadCorr::setAlwaysRadCorr();
00149 report(INFO,"EvtGen")
00150 << "As requested, PHOTOS will be turned on."<<endl;
00151 }
00152 else if (token=="normalPhotos"){
00153 EvtRadCorr::setNormalRadCorr();
00154 report(INFO,"EvtGen")
00155 << "As requested, PHOTOS will be turned on only when requested."<<endl;
00156 }
00157 else if (token=="Alias"){
00158
00159 std::string newname;
00160 std::string oldname;
00161
00162 newname=parser.getToken(itoken++);
00163 oldname=parser.getToken(itoken++);
00164
00165 EvtId id=EvtPDL::getId(oldname);
00166
00167 if (id==EvtId(-1,-1)) {
00168 report(ERROR,"EvtGen") <<"Unknown particle name:"<<oldname.c_str()
00169 <<" on line "<<parser.getLineofToken(itoken)<<endl;
00170 report(ERROR,"EvtGen") <<"Will terminate execution!"<<endl;
00171 ::abort();
00172 }
00173
00174 EvtPDL::alias(id,newname);
00175 if ( decaytable.size() < EvtPDL::entries() ) decaytable.resize(EvtPDL::entries());
00176
00177 } else if (token=="ModelAlias"){
00178 std::vector<std::string> modelArgList;
00179
00180 std::string aliasName=parser.getToken(itoken++);
00181 std::string modelName=parser.getToken(itoken++);
00182
00183 std::string nameTemp;
00184 do{
00185 nameTemp=parser.getToken(itoken++);
00186 if (nameTemp!=";") {
00187 modelArgList.push_back(nameTemp);
00188 }
00189 }while(nameTemp!=";");
00190 EvtModelAlias newAlias(aliasName,modelName,modelArgList);
00191 modelAliasList.push_back(newAlias);
00192 } else if (token=="ChargeConj"){
00193
00194 std::string aname;
00195 std::string abarname;
00196
00197 aname=parser.getToken(itoken++);
00198 abarname=parser.getToken(itoken++);
00199
00200 EvtId a=EvtPDL::getId(aname);
00201 EvtId abar=EvtPDL::getId(abarname);
00202
00203 if (a==EvtId(-1,-1)) {
00204 report(ERROR,"EvtGen") <<"Unknown particle name:"<<aname.c_str()
00205 <<" on line "<<parser.getLineofToken(itoken)<<endl;
00206 report(ERROR,"EvtGen") <<"Will terminate execution!"<<endl;
00207 ::abort();
00208 }
00209
00210 if (abar==EvtId(-1,-1)) {
00211 report(ERROR,"EvtGen") <<"Unknown particle name:"<<abarname.c_str()
00212 <<" on line "<<parser.getLineofToken(itoken)<<endl;
00213 report(ERROR,"EvtGen") <<"Will terminate execution!"<<endl;
00214 ::abort();
00215 }
00216
00217
00218 EvtPDL::aliasChgConj(a,abar);
00219
00220 } else if (modelist.isCommand(token)){
00221
00222 std::string cnfgstr;
00223
00224 cnfgstr=parser.getToken(itoken++);
00225
00226 modelist.storeCommand(token,cnfgstr);
00227
00228 } else if (token=="CDecay"){
00229
00230 std::string name;
00231
00232 name=parser.getToken(itoken++);
00233 ipar=EvtPDL::getId(name);
00234
00235 if (ipar==EvtId(-1,-1)) {
00236 report(ERROR,"EvtGen") <<"Unknown particle name:"<<name.c_str()
00237 <<" on line "
00238 <<parser.getLineofToken(itoken-1)<<endl;
00239 report(ERROR,"EvtGen") <<"Will terminate execution!"<<endl;
00240 ::abort();
00241 }
00242
00243 EvtId cipar=EvtPDL::chargeConj(ipar);
00244
00245 if (decaytable[ipar.getAlias()].getNMode()!=0) {
00246
00247 report(DEBUG,"EvtGen") <<
00248 "Redefined decay of "<<name.c_str()<<" in CDecay"<<endl;
00249
00250 decaytable[ipar.getAlias()].removeDecay();
00251 }
00252
00253
00254 decaytable[ipar.getAlias()].makeChargeConj(&decaytable[cipar.getAlias()]);
00255
00256 } else if (token=="Define"){
00257
00258 std::string name;
00259
00260 name=parser.getToken(itoken++);
00261
00262
00263 EvtSymTable::Define(name,parser.getToken(itoken++));
00264
00265
00266
00267
00268 } else if (token=="Particle"){
00269
00270 std::string pname;
00271 pname=parser.getToken(itoken++);
00272 report(INFO,"EvtGen") << pname.c_str() << endl;
00273
00274 double newMass=atof(parser.getToken(itoken++).c_str());
00275 EvtId thisPart = EvtPDL::getId(pname);
00276 double newWidth=EvtPDL::getMeanMass(thisPart);
00277 if ( parser.getNToken() > 3 ) newWidth=atof(parser.getToken(itoken++).c_str());
00278
00279
00280 EvtPDL::reSetMass(thisPart, newMass);
00281 EvtPDL::reSetWidth(thisPart, newWidth);
00282
00283 report(INFO,"EvtGen") << "Changing particle properties of " <<
00284 pname.c_str() << " Mass=" << newMass << " Width="<<newWidth<<endl;
00285
00286 } else if ( token=="ChangeMassMin") {
00287 std::string pname;
00288 pname=parser.getToken(itoken++);
00289 double tmass=atof(parser.getToken(itoken++).c_str());
00290
00291 EvtId thisPart = EvtPDL::getId(pname);
00292 EvtPDL::reSetMassMin(thisPart,tmass);
00293 report(DEBUG,"EvtGen") <<"Refined minimum mass for " << EvtPDL::name(thisPart).c_str() << " to be " << tmass << endl;
00294
00295 } else if ( token=="ChangeMassMax") {
00296 std::string pname;
00297 pname=parser.getToken(itoken++);
00298 double tmass=atof(parser.getToken(itoken++).c_str());
00299 EvtId thisPart = EvtPDL::getId(pname);
00300 EvtPDL::reSetMassMax(thisPart,tmass);
00301 report(DEBUG,"EvtGen") <<"Refined maximum mass for " << EvtPDL::name(thisPart).c_str() << " to be " << tmass << endl;
00302
00303 } else if ( token=="IncludeBirthFactor") {
00304 std::string pname;
00305 pname=parser.getToken(itoken++);
00306 bool yesno=false;
00307 if ( parser.getToken(itoken++).c_str()=="yes") yesno=true;
00308 EvtId thisPart = EvtPDL::getId(pname);
00309 EvtPDL::includeBirthFactor(thisPart,yesno);
00310 if ( yesno ) report(DEBUG,"EvtGen") <<"Include birth factor for " << EvtPDL::name(thisPart).c_str() <<endl;
00311 if ( !yesno ) report(DEBUG,"EvtGen") <<"No longer include birth factor for " << EvtPDL::name(thisPart).c_str() <<endl;
00312
00313
00314 } else if ( token=="IncludeDecayFactor") {
00315 std::string pname;
00316 pname=parser.getToken(itoken++);
00317 bool yesno=false;
00318 if ( parser.getToken(itoken++).c_str()=="yes") yesno=true;
00319 EvtId thisPart = EvtPDL::getId(pname);
00320 EvtPDL::includeDecayFactor(thisPart,yesno);
00321 if ( yesno ) report(DEBUG,"EvtGen") <<"Include decay factor for " << EvtPDL::name(thisPart).c_str() <<endl;
00322 if ( !yesno ) report(DEBUG,"EvtGen") <<"No longer include decay factor for " << EvtPDL::name(thisPart).c_str() <<endl;
00323
00324 }else if ( token=="AddFactorPn") {
00325 std::string pname;
00326 pname=parser.getToken(itoken++);
00327 double factor=atof(parser.getToken(itoken++).c_str());
00328 EvtId thisPart = EvtPDL::getId(pname);
00329 EvtPDL::addFactorPn(thisPart,factor);
00330 report(DEBUG,"EvtGen") <<"Include momentum factor Pn= "<<factor <<" for " << EvtPDL::name(thisPart).c_str() <<endl;
00331 }else if ( token=="LSNONRELBW") {
00332 std::string pname;
00333 pname=parser.getToken(itoken++);
00334 EvtId thisPart = EvtPDL::getId(pname);
00335 std::string tstr="NONRELBW";
00336 EvtPDL::changeLS(thisPart,tstr);
00337 report(DEBUG,"EvtGen") <<"Change lineshape to non-rel BW for " << EvtPDL::name(thisPart).c_str() <<endl;
00338 } else if ( token=="SP6LSFIX") {
00339 std::string pname;
00340 pname=parser.getToken(itoken++);
00341 EvtId thisPart = EvtPDL::getId(pname);
00342 EvtPDL::fixLSForSP8(thisPart);
00343 report(DEBUG,"EvtGen") <<"Fixed lineshape for SP6 --from M.Baak " << EvtPDL::name(thisPart).c_str() <<endl;
00344
00345 } else if ( token=="LSFLAT") {
00346 std::string pname;
00347 pname=parser.getToken(itoken++);
00348 EvtId thisPart = EvtPDL::getId(pname);
00349 std::string tstr="FLAT";
00350 EvtPDL::changeLS(thisPart,tstr);
00351 report(DEBUG,"EvtGen") <<"Change lineshape to flat for " << EvtPDL::name(thisPart).c_str() <<endl;
00352 } else if ( token=="LSMANYDELTAFUNC") {
00353 std::string pname;
00354 pname=parser.getToken(itoken++);
00355 EvtId thisPart = EvtPDL::getId(pname);
00356 std::string tstr="MANYDELTAFUNC";
00357 EvtPDL::changeLS(thisPart,tstr);
00358 report(DEBUG,"EvtGen") <<"Change lineshape to spikes for " << EvtPDL::name(thisPart).c_str() <<endl;
00359
00360 } else if ( token=="BlattWeisskopf") {
00361 std::string pname;
00362 pname=parser.getToken(itoken++);
00363 double tnum=atof(parser.getToken(itoken++).c_str());
00364 EvtId thisPart = EvtPDL::getId(pname);
00365 EvtPDL::reSetBlatt(thisPart,tnum);
00366 report(DEBUG,"EvtGen") <<"Redefined Blatt-Weisskopf factor " << EvtPDL::name(thisPart).c_str() << " to be " << tnum << endl;
00367 } else if ( token=="SetLineshapePW") {
00368 std::string pname;
00369 pname=parser.getToken(itoken++);
00370 EvtId thisPart = EvtPDL::getId(pname);
00371 std::string pnameD1=parser.getToken(itoken++);
00372 EvtId thisD1 = EvtPDL::getId(pnameD1);
00373 std::string pnameD2=parser.getToken(itoken++);
00374 EvtId thisD2 = EvtPDL::getId(pnameD2);
00375 int pw=atoi(parser.getToken(itoken++).c_str());
00376 report(DEBUG,"EvtGen") <<"Redefined Partial wave for " << pname.c_str() << " to " << pnameD1.c_str() << " " << pnameD2.c_str() << " ("<<pw<<")"<<endl;
00377 EvtPDL::setPWForDecay(thisPart,pw,thisD1,thisD2);
00378
00379 } else if (token=="Decay") {
00380
00381 std::string temp_fcn_new_model;
00382
00383 EvtDecayBase* temp_fcn_new;
00384
00385 double brfrsum=0.0;
00386
00387
00388
00389 parent=parser.getToken(itoken++);
00390 ipar=EvtPDL::getId(parent);
00391
00392 if (ipar==EvtId(-1,-1)) {
00393 report(ERROR,"EvtGen") <<"Unknown particle name:"<<parent.c_str()
00394 <<" on line "
00395 <<parser.getLineofToken(itoken-1)<<endl;
00396 report(ERROR,"EvtGen") <<"Will terminate execution!"<<endl;
00397 ::abort();
00398 }
00399
00400 if (decaytable[ipar.getAlias()].getNMode()!=0) {
00401 report(DEBUG,"EvtGen") <<"Redefined decay of "
00402 <<parent.c_str()<<endl;
00403 decaytable[ipar.getAlias()].removeDecay();
00404 }
00405
00406
00407 do{
00408
00409 token=parser.getToken(itoken++);
00410
00411 if (token!="Enddecay"){
00412
00413 i=0;
00414 while (token.c_str()[i++]!=0){
00415 if (isalpha(token.c_str()[i])){
00416 report(ERROR,"EvtGen") <<
00417 "Expected to find a branching fraction or Enddecay "<<
00418 "but found:"<<token.c_str()<<" on line "<<
00419 parser.getLineofToken(itoken-1)<<endl;
00420 report(ERROR,"EvtGen") << "Possibly to few arguments to model "<<
00421 "on previous line!"<<endl;
00422 report(ERROR,"EvtGen") << "Will terminate execution!"<<endl;
00423 ::abort();
00424 }
00425 }
00426
00427 brfr=atof(token.c_str());
00428
00429 int isname=EvtPDL::getId(parser.getToken(itoken)).getId()>=0;
00430 int ismodel=modelist.isModel(parser.getToken(itoken));
00431
00432 if (!(isname||ismodel)){
00433
00434 int iAlias;
00435 for(iAlias=0;iAlias<modelAliasList.size();iAlias++){
00436 if ( modelAliasList[iAlias].matchAlias(parser.getToken(itoken)) ) {
00437 ismodel=2;
00438 break;
00439 }
00440 }
00441 }
00442
00443 if (!(isname||ismodel)){
00444
00445 report(INFO,"EvtGen") << parser.getToken(itoken).c_str()
00446 << " is neither a particle name nor "
00447 << "the name of a model. "<<endl;
00448 report(INFO,"EvtGen") << "It was encountered on line "<<
00449 parser.getLineofToken(itoken)<<" of the decay file."<<endl;
00450 report(INFO,"EvtGen") << "Please fix it. Thank you."<<endl;
00451 report(INFO,"EvtGen") << "Be sure to check that the "
00452 << "correct case has been used. \n";
00453 report(INFO,"EvtGen") << "Terminating execution. \n";
00454 ::abort();
00455
00456 itoken++;
00457 }
00458
00459 n_daugh=0;
00460
00461 while(EvtPDL::getId(parser.getToken(itoken)).getId()>=0){
00462 sdaug=parser.getToken(itoken++);
00463 daught[n_daugh++]=EvtPDL::getId(sdaug);
00464 if (daught[n_daugh-1]==EvtId(-1,-1)) {
00465 report(ERROR,"EvtGen") <<"Unknown particle name:"<<sdaug.c_str()
00466 <<" on line "<<parser.getLineofToken(itoken)<<endl;
00467 report(ERROR,"EvtGen") <<"Will terminate execution!"<<endl;
00468 ::abort();
00469 }
00470 }
00471
00472
00473 model=parser.getToken(itoken++);
00474
00475
00476 int photos=0;
00477 int verbose=0;
00478 int summary=0;
00479
00480 do{
00481 if (model=="PHOTOS"){
00482 photos=1;
00483 model=parser.getToken(itoken++);
00484 }
00485 if (model=="VERBOSE"){
00486 verbose=1;
00487 model=parser.getToken(itoken++);
00488 }
00489 if (model=="SUMMARY"){
00490 summary=1;
00491 model=parser.getToken(itoken++);
00492 }
00493 }while(model=="PHOTOS"||
00494 model=="VERBOSE"||
00495 model=="SUMMARY");
00496
00497
00498 int iAlias;
00499 int foundAnAlias=-1;
00500 for(iAlias=0;iAlias<modelAliasList.size();iAlias++){
00501 if ( modelAliasList[iAlias].matchAlias(model) ) {
00502 foundAnAlias=iAlias;
00503 break;
00504 }
00505 }
00506
00507 if ( foundAnAlias==-1 ) {
00508 if(!modelist.isModel(model)){
00509 report(ERROR,"EvtGen") <<
00510 "Expected to find a model name,"<<
00511 "found:"<<model.c_str()<<" on line "<<
00512 parser.getLineofToken(itoken)<<endl;
00513 report(ERROR,"EvtGen") << "Will terminate execution!"<<endl;
00514 ::abort();
00515 }
00516 }
00517 else{
00518 model=modelAliasList[foundAnAlias].getName();
00519 }
00520
00521 temp_fcn_new_model=model;
00522 temp_fcn_new=modelist.getFcn(model);
00523
00524
00525 if (photos){
00526 temp_fcn_new->setPHOTOS();
00527 }
00528 if (verbose){
00529 temp_fcn_new->setVerbose();
00530 }
00531 if (summary){
00532 temp_fcn_new->setSummary();
00533 }
00534
00535
00536 std::vector<std::string> temp_fcn_new_args;
00537
00538 std::string name;
00539 int ierr;
00540
00541 if ( foundAnAlias==-1 ) {
00542 do{
00543 name=parser.getToken(itoken++);
00544 if (name!=";") {
00545 temp_fcn_new_args.push_back(EvtSymTable::Get(name,ierr));
00546 if (ierr) {
00547 report(ERROR,"EvtGen")
00548 <<"Reading arguments and found:"<<
00549 name.c_str()<<" on line:"<<
00550 parser.getLineofToken(itoken-1)<<endl;
00551 report(ERROR,"EvtGen")
00552 << "Will terminate execution!"<<endl;
00553 ::abort();
00554 }
00555 }
00556
00557 int ismodel=modelist.isModel(name);
00558 if (ismodel) {
00559 report(ERROR,"EvtGen")
00560 <<"Expected ';' but found:"<<
00561 name.c_str()<<" on line:"<<
00562 parser.getLineofToken(itoken-1)<<endl;
00563 report(ERROR,"EvtGen")
00564 << "Most probable error is omitted ';'."<<endl;
00565 report(ERROR,"EvtGen")
00566 << "Will terminate execution!"<<endl;
00567 ::abort();
00568 }
00569 }while(name!=";");
00570 }
00571 else{
00572 std::vector<std::string> copyMe=modelAliasList[foundAnAlias].getArgList();
00573 temp_fcn_new_args=copyMe;
00574 itoken++;
00575 }
00576
00577
00578 brfrsum+=brfr;
00579
00580 temp_fcn_new->saveDecayInfo(ipar,n_daugh,
00581 daught,
00582 temp_fcn_new_args.size(),
00583 temp_fcn_new_args,
00584 temp_fcn_new_model,
00585 brfr);
00586
00587 double massmin=0.0;
00588
00589
00590 for (i=0;i<temp_fcn_new->nRealDaughters();i++){
00591 if ( EvtPDL::getMinMass(daught[i])>0.0001 ){
00592 massmin+=EvtPDL::getMinMass(daught[i]);
00593 } else {
00594 massmin+=EvtPDL::getMeanMass(daught[i]);
00595 }
00596 }
00597
00598 decaytable[ipar.getAlias()].addMode(temp_fcn_new,brfrsum,massmin);
00599
00600
00601 }
00602 } while(token!="Enddecay");
00603
00604 decaytable[ipar.getAlias()].finalize();
00605
00606 }
00607
00608
00609 else if (token=="CopyDecay") {
00610 std::string newname;
00611 std::string oldname;
00612
00613 newname=parser.getToken(itoken++);
00614 oldname=parser.getToken(itoken++);
00615
00616 EvtId newipar=EvtPDL::getId(newname);
00617 EvtId oldipar=EvtPDL::getId(oldname);
00618
00619 if (oldipar==EvtId(-1,-1)) {
00620 report(ERROR,"EvtGen") <<"Unknown particle name:"<<oldname.c_str()
00621 <<" on line "<<parser.getLineofToken(itoken)<<endl;
00622 report(ERROR,"EvtGen") <<"Will terminate execution!"<<endl;
00623 ::abort();
00624 }
00625 if (newipar==EvtId(-1,-1)) {
00626 report(ERROR,"EvtGen") <<"Unknown particle name:"<<newname.c_str()
00627 <<" on line "<<parser.getLineofToken(itoken)<<endl;
00628 report(ERROR,"EvtGen") <<"Will terminate execution!"<<endl;
00629 ::abort();
00630 }
00631 if (decaytable[newipar.getAlias()].getNMode()!=0) {
00632 report(DEBUG,"EvtGen") <<"Redefining decay of "
00633 <<newname<<endl;
00634 decaytable[newipar.getAlias()].removeDecay();
00635 }
00636 decaytable[newipar.getAlias()] = decaytable[oldipar.getAlias()];
00637 }
00638
00639
00640 else if (token=="RemoveDecay") {
00641 parent = parser.getToken(itoken++);
00642 ipar = EvtPDL::getId(parent);
00643
00644 if (ipar==EvtId(-1,-1)) {
00645 report(ERROR,"EvtGen") <<"Unknown particle name:"<<parent.c_str()
00646 <<" on line "
00647 <<parser.getLineofToken(itoken-1)<<endl;
00648 report(ERROR,"EvtGen") <<"Will terminate execution!"<<endl;
00649 ::abort();
00650 }
00651
00652 if (decaytable[ipar.getAlias()].getNMode()==0) {
00653 report(DEBUG,"EvtGen") << "No decays to delete for "
00654 << parent.c_str() << endl;
00655 } else {
00656 report(DEBUG,"EvtGen") <<"Deleting selected decays of "
00657 <<parent.c_str()<<endl;
00658 }
00659
00660 do {
00661 token = parser.getToken(itoken);
00662
00663 if (token != "Enddecay") {
00664 n_daugh = 0;
00665 while (EvtPDL::getId(parser.getToken(itoken)).getId() >= 0) {
00666 sdaug = parser.getToken(itoken++);
00667 daught[n_daugh++] = EvtPDL::getId(sdaug);
00668 if (daught[n_daugh-1]==EvtId(-1,-1)) {
00669 report(ERROR,"EvtGen") <<"Unknown particle name:"<<sdaug.c_str()
00670 <<" on line "<<parser.getLineofToken(itoken)<<endl;
00671 report(ERROR,"EvtGen") <<"Will terminate execution!"<<endl;
00672 ::abort();
00673 }
00674 }
00675 token = parser.getToken(itoken);
00676 if (token != ";") {
00677 report(ERROR,"EvtGen")
00678 <<"Expected ';' but found:"<<
00679 token <<" on line:"<<
00680 parser.getLineofToken(itoken-1)<<endl;
00681 report(ERROR,"EvtGen")
00682 << "Most probable error is omitted ';'."<<endl;
00683 report(ERROR,"EvtGen")
00684 << "Will terminate execution!"<<endl;
00685 ::abort();
00686 }
00687 token = parser.getToken(itoken++);
00688 EvtDecayBase* temp_fcn_new = modelist.getFcn("PHSP");
00689 std::vector<std::string> temp_fcn_new_args;
00690 std::string temp_fcn_new_model("PHSP");
00691 temp_fcn_new->saveDecayInfo(ipar, n_daugh,
00692 daught,
00693 0,
00694 temp_fcn_new_args,
00695 temp_fcn_new_model,
00696 0.);
00697 decaytable[ipar.getAlias()].removeMode(temp_fcn_new);
00698 }
00699 } while (token != "Enddecay");
00700 itoken++;
00701 }
00702 else if (token!="End"){
00703
00704 report(ERROR,"EvtGen") << "Found unknown command:'"<<token.c_str()<<"' on line "
00705 <<parser.getLineofToken(itoken)<<endl;
00706 report(ERROR,"EvtGen") << "Will terminate execution!"<<endl;
00707 ::abort();
00708
00709 }
00710
00711 } while ((token!="End")&&itoken!=parser.getNToken());
00712
00713
00714
00715 int ii;
00716 for ( ii=0; ii<EvtPDL::entries(); ii++){
00717 EvtId temp(ii,ii);
00718 int nModTot=getNMode(ii);
00719
00720 if ( nModTot == 0 ) continue;
00721
00722 if ( EvtPDL::getWidth(temp) < 0.0000001 ) continue;
00723 int jj;
00724 double minMass=EvtPDL::getMaxMass(temp);
00725 for (jj=0; jj<nModTot; jj++) {
00726 double tmass=decaytable[ii].getDecay(jj).getMassMin();
00727 if ( tmass< minMass) minMass=tmass;
00728 }
00729 if ( minMass > EvtPDL::getMinMass(temp) ) {
00730
00731 report(INFO,"EvtGen") << "Given allowed decays, resetting minMass " << EvtPDL::name(temp).c_str() << " "
00732 << EvtPDL::getMinMass(temp) << " to " << minMass << endl;
00733 EvtPDL::reSetMassMin(temp,minMass);
00734 }
00735 }
00736
00737
00738 }
00739
00740 int EvtDecayTable::findChannel(EvtId parent, std::string model,
00741 int ndaug, EvtId *daugs,
00742 int narg, std::string *args){
00743
00744 int i,j,right;
00745 EvtId daugs_scratch[50];
00746 int nmatch,k;
00747
00748 for(i=0;i<decaytable[parent.getAlias()].getNMode();i++){
00749
00750 right=1;
00751
00752 right=right&&model==decaytable[parent.getAlias()].
00753 getDecay(i).getDecayModel()->getModelName();
00754 right=right&&(ndaug==decaytable[parent.getAlias()].
00755 getDecay(i).getDecayModel()->getNDaug());
00756 right=right&&(narg==decaytable[parent.getAlias()].
00757 getDecay(i).getDecayModel()->getNArg());
00758
00759 if ( right ){
00760
00761
00762
00763 for(j=0;j<ndaug;j++){
00764 daugs_scratch[j]=daugs[j];
00765 }
00766
00767 nmatch=0;
00768
00769 for(j=0;j<decaytable[parent.getAlias()].
00770 getDecay(i).getDecayModel()->getNDaug();j++){
00771
00772 for(k=0;k<ndaug;k++){
00773 if (daugs_scratch[k]==decaytable[parent.getAlias()].
00774 getDecay(i).getDecayModel()->getDaug(j)){
00775 daugs_scratch[k]=EvtId(-1,-1);
00776 nmatch++;
00777 break;
00778 }
00779 }
00780 }
00781
00782 right=right&&(nmatch==ndaug);
00783
00784 for(j=0;j<decaytable[parent.getAlias()].
00785 getDecay(i).getDecayModel()->getNArg();j++){
00786 right=right&&(args[j]==decaytable[parent.getAlias()].
00787 getDecay(i).getDecayModel()->getArgStr(j));
00788 }
00789 }
00790 if (right) return i;
00791 }
00792 return -1;
00793 }
00794
00795 int EvtDecayTable::inChannelList(EvtId parent, int ndaug, EvtId *daugs){
00796
00797 int i,j,k;
00798
00799 EvtId daugs_scratch[MAX_DAUG];
00800
00801 int dsum=0;
00802 for(i=0;i<ndaug;i++){
00803 dsum+=daugs[i].getAlias();
00804 }
00805
00806 int nmatch;
00807
00808 int ipar=parent.getAlias();
00809
00810 int nmode=decaytable[ipar].getNMode();
00811
00812 for(i=0;i<nmode;i++){
00813
00814 EvtDecayBase* thedecaymodel=decaytable[ipar].getDecay(i).getDecayModel();
00815
00816 if (thedecaymodel->getDSum()==dsum){
00817
00818 int nd=thedecaymodel->getNDaug();
00819
00820 if (ndaug==nd){
00821 for(j=0;j<ndaug;j++){
00822 daugs_scratch[j]=daugs[j];
00823 }
00824 nmatch=0;
00825 for(j=0;j<nd;j++){
00826 for(k=0;k<ndaug;k++){
00827 if (EvtId(daugs_scratch[k])==thedecaymodel->getDaug(j)){
00828 daugs_scratch[k]=EvtId(-1,-1);
00829 nmatch++;
00830 break;
00831 }
00832 }
00833 }
00834 if ((nmatch==ndaug)&&
00835 (!
00836 ((thedecaymodel->getModelName()=="JETSET")||
00837 (thedecaymodel->getModelName()=="LUNDCHARM")||
00838 (thedecaymodel->getModelName()=="PYTHIA")))){
00839 return i;
00840 }
00841 }
00842 }
00843 }
00844
00845 return -1;
00846 }
00847