Main Page | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members | Related Pages

EvtISGW2 Class Reference

#include <EvtISGW2.hh>

Inheritance diagram for EvtISGW2:

EvtDecayAmp EvtDecayBase List of all members.

Public Member Functions

void checkNArg (int a1, int a2=-1, int a3=-1, int a4=-1)
void checkNDaug (int d1, int d2=-1)
void checkQ ()
void checkSpinDaughter (int d1, EvtSpinType::spintype sp)
void checkSpinParent (EvtSpinType::spintype sp)
EvtDecayBase * clone ()
virtual void command (std::string cmd)
virtual std::string commandName ()
void decay (EvtParticle *p)
void disableCheckQ ()
 EvtISGW2 ()
double getArg (int j)
double * getArgs ()
std::string * getArgsStr ()
std::string getArgStr (int j)
double getBranchingFraction ()
EvtId getDaug (int i)
EvtId * getDaugs ()
int getDSum ()
std::string getModelName ()
void getName (std::string &name)
int getNArg ()
int getNDaug ()
EvtId getParentId ()
int getPHOTOS ()
double getProbMax (double prob)
void init ()
void initProbMax ()
void makeDecay (EvtParticle *p)
virtual bool matchingDecay (const EvtDecayBase &other) const
void noProbMax ()
virtual int nRealDaughters ()
void printSummary ()
double resetProbMax (double prob)
void saveDecayInfo (EvtId ipar, int ndaug, EvtId *daug, int narg, std::vector< std::string > &args, std::string name, double brfr)
void setPHOTOS ()
void setProbMax (double prbmx)
void setSummary ()
void setVerbose ()
void setWeight (double weight)
int summary ()
int verbose ()
void vertex (int *i1, const EvtComplex &amp)
void vertex (int i1, int i2, int i3, const EvtComplex &amp)
void vertex (int i1, int i2, const EvtComplex &amp)
void vertex (int i1, const EvtComplex &amp)
void vertex (const EvtComplex &amp)
virtual ~EvtISGW2 ()

Static Public Member Functions

void findMass (EvtParticle *p)
void findMasses (EvtParticle *p, int ndaugs, EvtId daugs[10], double masses[10])
double findMaxMass (EvtParticle *p)

Protected Member Functions

bool daugsDecayedByParentModel ()

Protected Attributes

EvtAmp _amp2
bool _daugsDecayedByParentModel

Private Attributes

EvtSemiLeptonicAmp * calcamp
EvtSemiLeptonicFF * isgw2ffmodel

Constructor & Destructor Documentation

EvtISGW2::EvtISGW2 (   ) 
  

00037                   :
00038    isgw2ffmodel(0)
00039   ,calcamp(0)
00040 {}

EvtISGW2::~EvtISGW2 (   )  [virtual]
  

00043                     {
00044   delete isgw2ffmodel;
00045   isgw2ffmodel=0;
00046   delete calcamp;
00047   calcamp=0;
00048 }

Member Function Documentation

void EvtDecayBase::checkNArg (  int  a1,
int  a2 = -1,
int  a3 = -1,
int  a4 = -1
)  [inherited]
  

00483                                                            {
00484 
00485   if ( _narg != a1 && _narg != a2 && _narg != a3 && _narg != a4 ) {
00486     report(ERROR,"EvtGen") << _modelname.c_str() << " generator expected "<<endl;
00487     report(ERROR,"EvtGen") << a1<<endl;; 
00488     if ( a2>-1) {
00489       report(ERROR,"EvtGen") << " or " << a2<<endl; 
00490     }
00491     if ( a3>-1) {
00492       report(ERROR,"EvtGen") << " or " << a3<<endl; 
00493     }
00494     if ( a4>-1) {
00495       report(ERROR,"EvtGen") << " or " << a4<<endl; 
00496     }
00497     report(ERROR,"EvtGen") << " arguments but found:"<< _narg << endl;
00498     printSummary();
00499     report(ERROR,"EvtGen") << "Will terminate execution!"<<endl;
00500     ::abort();
00501 
00502   } 
00503 
00504 }

void EvtDecayBase::checkNDaug (  int  d1,
int  d2 = -1
)  [inherited]
  

00505                                            {
00506 
00507   if ( _ndaug != d1 && _ndaug != d2 ) {
00508     report(ERROR,"EvtGen") << _modelname.c_str() << " generator expected ";
00509     report(ERROR,"EvtGen") << d1; 
00510     if ( d2>-1) {
00511       report(ERROR,"EvtGen") << " or " << d2; 
00512     }
00513     report(ERROR,"EvtGen") << " daughters but found:"<< _ndaug << endl;
00514     printSummary();
00515     report(ERROR,"EvtGen") << "Will terminate execution!"<<endl;
00516     ::abort();
00517   } 
00518 
00519 }

void EvtDecayBase::checkQ (   )  [inherited]
  

00035                           {
00036   int i;
00037   int q=0;
00038   int qpar;
00039 
00040   //If there are no daughters (jetset etc) then we do not
00041   //want to do this test.  Why?  Because sometimes the parent
00042   //will have a nonzero charge.
00043 
00044   if ( _ndaug != 0) {
00045     for(i=0; i<_ndaug; i++ ) {
00046       q += EvtPDL::chg3(_daug[i]);
00047     }
00048     qpar = EvtPDL::chg3(_parent);
00049 
00050     if ( q != qpar ) {
00051       report(ERROR,"EvtGen") <<_modelname.c_str()<< " generator expected "
00052                              << " charge to be conserved, found:"<<endl;
00053       report(ERROR,"EvtGen") << "Parent charge of "<<(qpar/3)<<endl;
00054       report(ERROR,"EvtGen") << "Sum of daughter charge of "<<(q/3)<<endl;
00055       report(ERROR,"EvtGen") << "The parent is "<< EvtPDL::name(_parent).c_str()<<endl;
00056       for(i=0; i<_ndaug; i++ ) {
00057       report(ERROR,"EvtGen") << "Daughter "<< EvtPDL::name(_daug[i]).c_str()<<endl;
00058       }
00059       report(ERROR,"EvtGen") << "Will terminate execution!"<<endl;
00060       
00061       ::abort();
00062     }
00063   }
00064 }

void EvtDecayBase::checkSpinDaughter (  int  d1,
EvtSpinType::spintype  sp
)  [inherited]
  

00534                                                                    {
00535 
00536   EvtSpinType::spintype parenttype = EvtPDL::getSpinType(getDaug(d1));
00537   if ( parenttype != sp ) {
00538     report(ERROR,"EvtGen") << _modelname.c_str() 
00539                            << " did not get the correct daughter spin d=" 
00540                            << d1 << endl;
00541     printSummary();
00542     report(ERROR,"EvtGen") << "Will terminate execution!"<<endl;
00543     ::abort();
00544   } 
00545 
00546 }

void EvtDecayBase::checkSpinParent (  EvtSpinType::spintype  sp  )  [inherited]
  

00521                                                          {
00522 
00523   EvtSpinType::spintype parenttype = EvtPDL::getSpinType(getParentId());
00524   if ( parenttype != sp ) {
00525     report(ERROR,"EvtGen") << _modelname.c_str() 
00526                            << " did not get the correct parent spin\n";
00527     printSummary();
00528     report(ERROR,"EvtGen") << "Will terminate execution!"<<endl;
00529     ::abort();
00530   } 
00531 
00532 }

EvtDecayBase * EvtISGW2::clone (   )  [virtual]
 

Implements EvtDecayBase.

00058                              {
00059 
00060   return new EvtISGW2;
00061 
00062 }

void EvtDecayBase::command (  std::string  cmd  )  [virtual, inherited]
 

Reimplemented in EvtJetSet, EvtLunda, EvtLundCharm, EvtPythia, and EvtTauola.

00132                                        {
00133   report(ERROR,"EvtGen") << "Should never call EvtDecayBase::command"<<endl;
00134   ::abort();
00135 }

std::string EvtDecayBase::commandName (   )  [virtual, inherited]
 

Reimplemented in EvtJetSet, EvtLunda, EvtLundCharm, EvtPythia, and EvtTauola.

00129                                    {
00130   return std::string("");
00131 }

bool EvtDecayBase::daugsDecayedByParentModel (   )  [inline, protected, inherited]
  

00111 {return _daugsDecayedByParentModel;}

void EvtISGW2::decay (  EvtParticle *  p  )  [virtual]
 

Implements EvtDecayBase.

00064                                     {
00065 
00066   p->initializePhaseSpace(getNDaug(),getDaugs());
00067 
00068   calcamp->CalcAmp(p,_amp2,isgw2ffmodel);
00069 
00070 }

void EvtDecayBase::disableCheckQ (   )  [inline, inherited]
  

00062 {_chkCharge=0;};

void EvtDecayBase::findMass (  EvtParticle *  p  )  [static, inherited]
  

00347                                           {
00348 
00349   //Need to also check that this mass does not screw
00350   //up the parent
00351   //This code assumes that for the ith daughter, 0..i-1
00352   //already have a mass
00353   double maxOkMass=findMaxMass(p);
00354 
00355   int count=0;
00356   double mass;
00357   bool massOk=false;
00358   int i;
00359   while (!massOk) { 
00360     count++;
00361     if ( count > 10000 ) {
00362       report(INFO,"EvtGen") << "Can not find a valid mass for: " << EvtPDL::name(p->getId()).c_str() <<endl;
00363       report(INFO,"EvtGen") << "Now printing parent and/or grandparent tree\n";
00364       if ( p->getParent() ) {
00365         if ( p->getParent()->getParent() ) {
00366           p->getParent()->getParent()->printTree();
00367           report(INFO,"EvtGen") << p->getParent()->getParent()->mass() <<endl;
00368           report(INFO,"EvtGen") << p->getParent()->mass() <<endl;
00369         }
00370         else{
00371           p->getParent()->printTree();
00372           report(INFO,"EvtGen") << p->getParent()->mass() <<endl;
00373         }
00374       }
00375       else  p->printTree();
00376       report(INFO,"EvtGen") << "maxokmass=" << maxOkMass << " " << EvtPDL::getMinMass(p->getId()) << " " << EvtPDL::getMaxMass(p->getId())<<endl;
00377       if ( p->getNDaug() ) { 
00378         for (i=0; i<p->getNDaug(); i++) {
00379           report(INFO,"EvtGen") << p->getDaug(i)->mass()<<" ";
00380         }
00381         report(INFO,"EvtGen") << endl;
00382       }
00383       if ( maxOkMass >= EvtPDL::getMinMass(p->getId()) ) {
00384         report(INFO,"EvtGen") << "taking a default value\n";
00385         p->setMass(maxOkMass);
00386         return;
00387       } 
00388       assert(0);
00389     }
00390     mass = EvtPDL::getMass(p->getId());
00391     //Just need to check that this mass is > than
00392     //the mass of all daughters
00393     double massSum=0.;
00394     if ( p->getNDaug() ) { 
00395       for (i=0; i<p->getNDaug(); i++) {
00396         massSum+= p->getDaug(i)->mass();
00397       }
00398     }
00399     //some special cases are handled with 0 (stable) or 1 (k0->ks/kl) daughters
00400     if (p->getNDaug()<2)  massOk=true;
00401     if ( p->getParent() ) {
00402       if ( p->getParent()->getNDaug()==1 ) massOk=true;
00403     }
00404     if ( !massOk ) { 
00405       if (massSum < mass) massOk=true;
00406       if ( mass> maxOkMass) massOk=false;
00407     }
00408   }
00409 
00410   p->setMass(mass);
00411   
00412 }

void EvtDecayBase::findMasses (  EvtParticle *  p,
int  ndaugs,
EvtId  daugs[10],
double  masses[10]
)  [static, inherited]
  

00416                                                                      {
00417 
00418   int i;
00419   double mass_sum;
00420 
00421   int count=0;
00422 
00423   if (!( p->firstornot() )) {
00424     for (i = 0; i < ndaugs; i++ ) {
00425       masses[i] = p->getDaug(i)->mass();
00426     } //for
00427   } //if
00428   else {
00429     p->setFirstOrNot();
00430     // if only one daughter do it
00431 
00432     if (ndaugs==1) {
00433       masses[0]=p->mass();
00434       return;
00435     }
00436     
00437     //until we get a combo whose masses are less than _parent mass.
00438     do {
00439       mass_sum = 0.0;
00440 
00441       for (i = 0; i < ndaugs; i++ ) {
00442         masses[i] = EvtPDL::getMass(daugs[i]);
00443         mass_sum = mass_sum + masses[i];
00444       } 
00445 
00446       count++;
00447 
00448      
00449       if(count==10000) {
00450         report(ERROR,"EvtGen") <<"Decaying particle:"<<
00451           EvtPDL::name(p->getId()).c_str()<<" (m="<<p->mass()<<")"<<endl;
00452         report(ERROR,"EvtGen") <<"To the following daugthers"<<endl;
00453         for (i = 0; i < ndaugs; i++ ) {
00454           report(ERROR,"EvtGen") <<  
00455             EvtPDL::name(daugs[i]).c_str() << endl;
00456         } 
00457         report(ERROR,"EvtGen") << "Has been rejected "<<count
00458                                << " times, will now take minimal masses "
00459                                << " of daugthers"<<endl;
00460         
00461         mass_sum=0.;
00462         for (i = 0; i < ndaugs; i++ ) {
00463           masses[i] = EvtPDL::getMinMass(daugs[i]);
00464           mass_sum = mass_sum + masses[i];
00465         } 
00466         if (mass_sum > p->mass()){
00467           report(ERROR,"EvtGen") << "Parent mass="<<p->mass()
00468                                  << "to light for daugthers."<<endl
00469                                  << "Will throw the event away."<<endl;
00470           //dont terminate - start over on the event.
00471           EvtStatus::setRejectFlag();
00472           mass_sum=0.;
00473           //      ::abort();
00474         }
00475 
00476       }
00477     } while ( mass_sum > p->mass());
00478   } //else
00479   
00480   return;
00481 }

double EvtDecayBase::findMaxMass (  EvtParticle *  p  )  [static, inherited]
  

00312                                                {
00313 
00314   
00315   double maxOkMass=EvtPDL::getMaxMass(p->getId());
00316 
00317   //protect against vphotons
00318   if ( maxOkMass < 0.0000000001 ) return 10000000.;
00319   //and against already determined masses
00320   if ( p->hasValidP4() ) maxOkMass=p->mass();
00321 
00322   EvtParticle *par=p->getParent();
00323   if ( par ) {
00324     double maxParMass=findMaxMass(par);
00325     int i;
00326     double minDaugMass=0.;
00327     for(i=0;i<par->getNDaug();i++){
00328       EvtParticle *dau=par->getDaug(i);
00329       if ( dau!=p) {
00330         // it might already have a mass
00331         if ( dau->isInitialized() || dau->hasValidP4() )
00332           minDaugMass+=dau->mass();
00333         else
00334         //give it a bit of phase space 
00335           minDaugMass+=1.000001*EvtPDL::getMinMass(dau->getId());
00336       }
00337     }
00338     if ( maxOkMass>(maxParMass-minDaugMass)) maxOkMass=maxParMass-minDaugMass;
00339   }
00340   return maxOkMass;
00341 }

double EvtDecayBase::getArg (  int  j  )  [inherited]
  

00565                                  {
00566 
00567   // Verify string
00568 
00569   const char* str = _args[j].c_str();
00570   int i = 0;
00571   while(str[i]!=0){
00572     if (isalpha(str[i]) && str[i]!='e') {
00573 
00574       report(INFO,"EvtGen") << "String " << str << " is not a number" << endl;
00575       assert(0);
00576     }
00577     i++;
00578   }
00579   
00580   char** tc=0; 
00581   return strtod(_args[j].c_str(),tc);
00582 }

double * EvtDecayBase::getArgs (   )  [inherited]
  

00548                               {
00549 
00550   if ( _argsD ) return _argsD;
00551   //The user has asked for a list of doubles - the arguments 
00552   //better all be doubles...
00553   if ( _narg==0 ) return _argsD;
00554 
00555   _argsD = new double[_narg];
00556 
00557   int i;
00558   char * tc;
00559   for(i=0;i<_narg;i++) { 
00560     _argsD[i] =  strtod(_args[i].c_str(),&tc);
00561   }
00562   return _argsD;
00563 }

std::string* EvtDecayBase::getArgsStr (   )  [inline, inherited]
  

00073 {return _args;}

std::string EvtDecayBase::getArgStr (  int  j  )  [inline, inherited]
  

00075 {return _args[j];}

double EvtDecayBase::getBranchingFraction (   )  [inline, inherited]
  

00061 {return _brfr;}

EvtId EvtDecayBase::getDaug (  int  i  )  [inline, inherited]
  

00066 {return _daug[i];}

EvtId* EvtDecayBase::getDaugs (   )  [inline, inherited]
  

00065 {return _daug;}

int EvtDecayBase::getDSum (   )  [inline, inherited]
  

00077 {return _dsum; }

std::string EvtDecayBase::getModelName (   )  [inline, inherited]
  

00076 {return _modelname; }

void EvtISGW2::getName (  std::string &  name  )  [virtual]
 

Implements EvtDecayBase.

00050                                            {
00051 
00052   model_name="ISGW2";     
00053 
00054 }

int EvtDecayBase::getNArg (   )  [inline, inherited]
  

00067 {return _narg;}

int EvtDecayBase::getNDaug (   )  [inline, inherited]
  

00064 {return _ndaug;}

EvtId EvtDecayBase::getParentId (   )  [inline, inherited]
  

00060 {return _parent;}

int EvtDecayBase::getPHOTOS (   )  [inline, inherited]
  

00068 {return _photos;}

double EvtDecayBase::getProbMax (  double  prob  )  [inherited]
  

00067                                              {
00068 
00069   int i;
00070 
00071   //diagnostics
00072   sum_prob+=prob;
00073   if (prob>max_prob) max_prob=prob;
00074 
00075 
00076   if ( defaultprobmax && ntimes_prob<=500 ) { 
00077     //We are building up probmax with this iteration
00078      ntimes_prob += 1;
00079      if ( prob > probmax ) { probmax = prob;}
00080      if (ntimes_prob==500) { 
00081        probmax*=1.2;
00082      }
00083      return 1000000.0*prob;
00084   }
00085 
00086   if ( prob> probmax*1.0001) {
00087 
00088     report(INFO,"EvtGen") << "prob > probmax:("<<prob<<">"<<probmax<<")";
00089     report(INFO,"") << "("<<_modelname.c_str()<<") ";
00090     report(INFO,"") << EvtPDL::name(_parent).c_str()<<" -> ";
00091     for(i=0;i<_ndaug;i++){
00092        report(INFO,"") << EvtPDL::name(_daug[i]).c_str() << " ";
00093     }
00094     report(INFO,"") << endl;
00095 
00096     if (defaultprobmax) probmax = prob;
00097 
00098   }
00099 
00100   ntimes_prob += 1;
00101 
00102 
00103   return probmax;
00104 
00105 } //getProbMax

void EvtISGW2::init (   )  [virtual]
 

Reimplemented from EvtDecayBase.

00779                    {
00780 
00781   checkNArg(0);
00782   checkNDaug(3);
00783 
00784   //We expect the parent to be a scalar 
00785   //and the daughters to be X lepton neutrino
00786 
00787   checkSpinParent(EvtSpinType::SCALAR);
00788   checkSpinDaughter(1,EvtSpinType::DIRAC);
00789   checkSpinDaughter(2,EvtSpinType::NEUTRINO);
00790 
00791   EvtSpinType::spintype mesontype=EvtPDL::getSpinType(getDaug(0));
00792 
00793   isgw2ffmodel = new EvtISGW2FF;
00794   
00795   if ( mesontype==EvtSpinType::SCALAR ) { 
00796     calcamp = new EvtSemiLeptonicScalarAmp; 
00797   }
00798   if ( mesontype==EvtSpinType::VECTOR ) { 
00799     calcamp = new EvtSemiLeptonicVectorAmp; 
00800   }
00801   if ( mesontype==EvtSpinType::TENSOR ) { 
00802     calcamp = new EvtSemiLeptonicTensorAmp; 
00803   }
00804   
00805 }

void EvtISGW2::initProbMax (   )  [virtual]
 

Reimplemented from EvtDecayBase.

00072                            {
00073 
00074   //added by Lange Jan4,2000
00075   static EvtId EM=EvtPDL::getId("e-");
00076   static EvtId EP=EvtPDL::getId("e+");
00077   static EvtId MUM=EvtPDL::getId("mu-");
00078   static EvtId MUP=EvtPDL::getId("mu+");
00079   static EvtId TAUM=EvtPDL::getId("tau-");
00080   static EvtId TAUP=EvtPDL::getId("tau+");
00081 
00082   static EvtId BP=EvtPDL::getId("B+");
00083   static EvtId BM=EvtPDL::getId("B-");
00084   static EvtId B0=EvtPDL::getId("B0");
00085   static EvtId B0B=EvtPDL::getId("anti-B0");
00086   static EvtId BS0=EvtPDL::getId("B_s0");
00087   static EvtId BSB=EvtPDL::getId("anti-B_s0");
00088   static EvtId BCP=EvtPDL::getId("B_c+");
00089   static EvtId BCM=EvtPDL::getId("B_c-");
00090 
00091   static EvtId DST0=EvtPDL::getId("D*0");
00092   static EvtId DSTB=EvtPDL::getId("anti-D*0");
00093   static EvtId DSTP=EvtPDL::getId("D*+");
00094   static EvtId DSTM=EvtPDL::getId("D*-");
00095   static EvtId D0=EvtPDL::getId("D0");
00096   static EvtId D0B=EvtPDL::getId("anti-D0");
00097   static EvtId DP=EvtPDL::getId("D+");
00098   static EvtId DM=EvtPDL::getId("D-");
00099 
00100   static EvtId D1P1P=EvtPDL::getId("D_1+");
00101   static EvtId D1P1N=EvtPDL::getId("D_1-");
00102   static EvtId D1P10=EvtPDL::getId("D_10");
00103   static EvtId D1P1B=EvtPDL::getId("anti-D_10");
00104 
00105   static EvtId D3P2P=EvtPDL::getId("D_2*+");
00106   static EvtId D3P2N=EvtPDL::getId("D_2*-");
00107   static EvtId D3P20=EvtPDL::getId("D_2*0");
00108   static EvtId D3P2B=EvtPDL::getId("anti-D_2*0");
00109 
00110   static EvtId D3P1P=EvtPDL::getId("D'_1+");
00111   static EvtId D3P1N=EvtPDL::getId("D'_1-");
00112   static EvtId D3P10=EvtPDL::getId("D'_10");
00113   static EvtId D3P1B=EvtPDL::getId("anti-D'_10");
00114 
00115   static EvtId D3P0P=EvtPDL::getId("D_0*+");
00116   static EvtId D3P0N=EvtPDL::getId("D_0*-");
00117   static EvtId D3P00=EvtPDL::getId("D_0*0");
00118   static EvtId D3P0B=EvtPDL::getId("anti-D_0*0");
00119 
00120   static EvtId D21S0P=EvtPDL::getId("D(2S)+");
00121   static EvtId D21S0N=EvtPDL::getId("D(2S)-");
00122   static EvtId D21S00=EvtPDL::getId("D(2S)0");
00123   static EvtId D21S0B=EvtPDL::getId("anti-D(2S)0");
00124 
00125   static EvtId D23S1P=EvtPDL::getId("D*(2S)+");
00126   static EvtId D23S1N=EvtPDL::getId("D*(2S)-");
00127   static EvtId D23S10=EvtPDL::getId("D*(2S)0");
00128   static EvtId D23S1B=EvtPDL::getId("anti-D*(2S)0");
00129 
00130   static EvtId RHO2S0=EvtPDL::getId("rho(2S)0");
00131   static EvtId RHO2SP=EvtPDL::getId("rho(2S)+");
00132   static EvtId RHO2SM=EvtPDL::getId("rho(2S)-");
00133   static EvtId OMEG2S=EvtPDL::getId("omega(2S)");
00134   static EvtId ETA2S=EvtPDL::getId("eta(2S)");
00135 
00136   static EvtId PI2S0=EvtPDL::getId("pi(2S)0");
00137   static EvtId PI2SP=EvtPDL::getId("pi(2S)+");
00138   static EvtId PI2SM=EvtPDL::getId("pi(2S)-");
00139 
00140   static EvtId PIP=EvtPDL::getId("pi+");
00141   static EvtId PIM=EvtPDL::getId("pi-");
00142   static EvtId PI0=EvtPDL::getId("pi0");
00143 
00144   static EvtId RHOP=EvtPDL::getId("rho+");
00145   static EvtId RHOM=EvtPDL::getId("rho-");
00146   static EvtId RHO0=EvtPDL::getId("rho0");
00147 
00148   static EvtId A2P=EvtPDL::getId("a_2+");
00149   static EvtId A2M=EvtPDL::getId("a_2-");
00150   static EvtId A20=EvtPDL::getId("a_20");
00151 
00152   static EvtId A1P=EvtPDL::getId("a_1+");
00153   static EvtId A1M=EvtPDL::getId("a_1-");
00154   static EvtId A10=EvtPDL::getId("a_10");
00155 
00156   static EvtId A0P=EvtPDL::getId("a_0+");
00157   static EvtId A0M=EvtPDL::getId("a_0-");
00158   static EvtId A00=EvtPDL::getId("a_00");
00159 
00160   static EvtId B1P=EvtPDL::getId("b_1+");
00161   static EvtId B1M=EvtPDL::getId("b_1-");
00162   static EvtId B10=EvtPDL::getId("b_10");
00163 
00164   static EvtId H1=EvtPDL::getId("h_1");
00165   static EvtId H1PR=EvtPDL::getId("h'_1");
00166 
00167   static EvtId F1=EvtPDL::getId("f_1");
00168   static EvtId F1PR=EvtPDL::getId("f'_1");
00169   static EvtId F0=EvtPDL::getId("f_0");
00170   static EvtId F0PR=EvtPDL::getId("f'_0");
00171   static EvtId F2=EvtPDL::getId("f_2");
00172   static EvtId F2PR=EvtPDL::getId("f'_2");
00173 
00174   static EvtId ETA=EvtPDL::getId("eta");
00175   static EvtId ETAPR=EvtPDL::getId("eta'");
00176   static EvtId OMEG=EvtPDL::getId("omega");
00177 
00178   static EvtId KP=EvtPDL::getId("K+");
00179   static EvtId KM=EvtPDL::getId("K-");
00180   static EvtId K0=EvtPDL::getId("K0");
00181   static EvtId KB=EvtPDL::getId("anti-K0");
00182   static EvtId K0S=EvtPDL::getId("K_S0");
00183   static EvtId K0L=EvtPDL::getId("K_L0");
00184 
00185   static EvtId KSTP=EvtPDL::getId("K*+");
00186   static EvtId KSTM=EvtPDL::getId("K*-");
00187   static EvtId KST0=EvtPDL::getId("K*0");
00188   static EvtId KSTB=EvtPDL::getId("anti-K*0");
00189 
00190   static EvtId K1P=EvtPDL::getId("K_1+");
00191   static EvtId K1M=EvtPDL::getId("K_1-");
00192   static EvtId K10=EvtPDL::getId("K_10");
00193   static EvtId K1B=EvtPDL::getId("anti-K_10");
00194 
00195   static EvtId K1STP=EvtPDL::getId("K'_1+");
00196   static EvtId K1STM=EvtPDL::getId("K'_1-");
00197   static EvtId K1ST0=EvtPDL::getId("K'_10");
00198   static EvtId K1STB=EvtPDL::getId("anti-K'_10");
00199 
00200   static EvtId K2STP=EvtPDL::getId("K_2*+");
00201   static EvtId K2STM=EvtPDL::getId("K_2*-");
00202   static EvtId K2ST0=EvtPDL::getId("K_2*0");
00203   static EvtId K2STB=EvtPDL::getId("anti-K_2*0");
00204 
00205   static EvtId K0STP=EvtPDL::getId("K_0*+");
00206   static EvtId K0STM=EvtPDL::getId("K_0*-");
00207   static EvtId K0ST0=EvtPDL::getId("K_0*0");
00208   static EvtId K0STB=EvtPDL::getId("anti-K_0*0");
00209 
00210   static EvtId PHI=EvtPDL::getId("phi");
00211   static EvtId DSP=EvtPDL::getId("D_s+");
00212   static EvtId DSM=EvtPDL::getId("D_s-");
00213 
00214   static EvtId DSSTP=EvtPDL::getId("D_s*+");
00215   static EvtId DSSTM=EvtPDL::getId("D_s*-");
00216   static EvtId DS1P=EvtPDL::getId("D_s1+");
00217   static EvtId DS1M=EvtPDL::getId("D_s1-");
00218   static EvtId DS0STP=EvtPDL::getId("D_s0*+");
00219   static EvtId DS0STM=EvtPDL::getId("D_s0*-");
00220   static EvtId DPS1P=EvtPDL::getId("D'_s1+");
00221   static EvtId DPS1M=EvtPDL::getId("D'_s1-");
00222   static EvtId DS2STP=EvtPDL::getId("D_s2*+");
00223   static EvtId DS2STM=EvtPDL::getId("D_s2*-");
00224 
00225 EvtId parnum,mesnum,lnum;
00226 
00227 parnum = getParentId();
00228 mesnum = getDaug(0);
00229 lnum = getDaug(1);
00230 
00231 
00232 if ( parnum==BP||parnum==BM||parnum==B0||parnum==B0B||parnum==BS0||parnum==BSB ) {
00233 
00234  
00235   if ( mesnum==DST0||mesnum==DSTP||mesnum==DSTB||mesnum==DSTM||mesnum==DSSTP||mesnum==DSSTM) {
00236 
00237     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00238        setProbMax(10000.0);
00239        return;
00240     }
00241     if ( lnum==TAUP||lnum==TAUM ) {
00242        setProbMax(7000.0);
00243        return;
00244     }
00245   }
00246 
00247  
00248   if ( mesnum==D0||mesnum==DP||mesnum==D0B||mesnum==DM||mesnum==DSP||mesnum==DSM) {
00249 
00250     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00251        setProbMax(4000.0);
00252        return;
00253     }
00254     if ( lnum==TAUP||lnum==TAUM ) {
00255        setProbMax(3500.0);
00256        return;
00257     }
00258   }
00259 
00260 
00261   if ( mesnum==D1P1P||mesnum==D1P1N||mesnum==D1P10||mesnum==D1P1B||mesnum==DS1P||mesnum==DS1M) {
00262     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00263       setProbMax(1300.0);
00264        return;
00265     }
00266     if ( lnum==TAUP||lnum==TAUM ) {
00267       setProbMax(480.0);
00268        return;
00269     }
00270   }
00271 
00272   if ( mesnum==D3P1P||mesnum==D3P1N||mesnum==D3P10||mesnum==D3P1B||mesnum==DS0STP||mesnum==DS0STM) {
00273 
00274     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00275        setProbMax(450.0);
00276        return;
00277     }
00278     if ( lnum==TAUP||lnum==TAUM ) {
00279       setProbMax(73.0);//???
00280        return;
00281     }
00282   }
00283 
00284   if ( mesnum==D3P0P||mesnum==D3P0N||mesnum==D3P00||mesnum==D3P0B||mesnum==DPS1P||mesnum==DPS1M) {
00285 
00286     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00287       setProbMax(200.0);
00288        return;
00289     }
00290     if ( lnum==TAUP||lnum==TAUM ) {
00291       setProbMax(90.0);
00292        return;
00293     }
00294   }
00295   if ( mesnum==D3P2P||mesnum==D3P2N||mesnum==D3P20||mesnum==D3P2B||mesnum==DS2STP||mesnum==DS2STM) {
00296 
00297     if ( mesnum==DS2STP|| mesnum==DS2STM) {
00298       setProbMax(550.0);
00299       return;
00300     }
00301     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00302        setProbMax(400.0);
00303        return;
00304     }
00305     if ( lnum==TAUP||lnum==TAUM ) {
00306        setProbMax(220.0);
00307        return;
00308     }
00309   }
00310 
00311   if ( mesnum==D21S0P||mesnum==D21S0N||mesnum==D21S00||mesnum==D21S0B) {
00312 
00313     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00314        setProbMax(16.0);
00315        return;
00316     }
00317     if ( lnum==TAUP||lnum==TAUM ) {
00318        setProbMax(3.0);
00319        return;
00320     }
00321   }
00322 
00323   if ( mesnum==D23S1P||mesnum==D23S1N||mesnum==D23S10||mesnum==D23S1B) {
00324 
00325     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00326        setProbMax(500.0);
00327        return;
00328     }
00329     if ( lnum==TAUP||lnum==TAUM ) {
00330        setProbMax(250.0);
00331        return;
00332     }
00333   }
00334 
00335   if ( mesnum==RHOP||mesnum==RHOM||mesnum==RHO0) {
00336 
00337     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00338        setProbMax(6500.0);
00339        return;
00340     }
00341     if ( lnum==TAUP||lnum==TAUM ) {
00342        setProbMax(6000.0);
00343        return;
00344     }
00345   }
00346 
00347   if ( mesnum==OMEG) {
00348 
00349     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00350        setProbMax(6800.0);
00351        return;
00352     }
00353     if ( lnum==TAUP||lnum==TAUM ) {
00354        setProbMax(6000.0);
00355        return;
00356     }
00357   }
00358 
00359   if ( mesnum==PIP||mesnum==PIM||mesnum==PI0) {
00360 
00361     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00362        setProbMax(1200.0);
00363        return;
00364     }
00365     if ( lnum==TAUP||lnum==TAUM ) {
00366        setProbMax(1150.0);
00367        return;
00368     }
00369   }
00370 
00371   if ( mesnum==ETA) {
00372 
00373     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00374        setProbMax(1800.0);
00375        return;
00376     }
00377     if ( lnum==TAUP||lnum==TAUM ) {
00378        setProbMax(1900.0);
00379        return;
00380     }
00381   }
00382 
00383   if ( mesnum==ETAPR) {
00384 
00385     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00386        setProbMax(3000.0);
00387        return;
00388     } 
00389     if ( lnum==TAUP||lnum==TAUM ) {
00390        setProbMax(3000.0);
00391        return;
00392     }
00393   }
00394 
00395 
00396   if ( mesnum==B1P||mesnum==B1M||mesnum==B10) {
00397 
00398     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00399        setProbMax(2500.0);
00400        return;
00401     }
00402     if ( lnum==TAUP||lnum==TAUM ) {
00403        setProbMax(1700.0);
00404        return;
00405     }
00406   }
00407 
00408   if ( mesnum==A0P||mesnum==A0M||mesnum==A00) {
00409 
00410     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00411        setProbMax(80.0);
00412        return;
00413     }
00414     if ( lnum==TAUP||lnum==TAUM ) {
00415        setProbMax(62.0);
00416        return;
00417     }
00418   }
00419 
00420   if ( mesnum==A1P||mesnum==A1M||mesnum==A10) {
00421 
00422     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00423        setProbMax(4500.0);
00424        return;
00425     }
00426     if ( lnum==TAUP||lnum==TAUM ) {
00427        setProbMax(3500.0);
00428        return;
00429     }
00430   }
00431 
00432   if ( mesnum==A2P||mesnum==A2M||mesnum==A20) {
00433 
00434     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00435        setProbMax(1200.0);
00436        return;
00437     }
00438     if ( lnum==TAUP||lnum==TAUM ) {
00439        setProbMax(1000.0);
00440        return;
00441     }
00442   }
00443 
00444   if ( mesnum==H1) {
00445 
00446     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00447        setProbMax(2600.0);
00448        return;
00449     }
00450     if ( lnum==TAUP||lnum==TAUM ) {
00451        setProbMax(2900.0);
00452        return;
00453     }
00454   }
00455 
00456   if ( mesnum==H1PR) {
00457 
00458     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00459        setProbMax(1400.0);
00460        return;
00461     }
00462     if ( lnum==TAUP||lnum==TAUM ) {
00463        setProbMax(1500.0);
00464        return;
00465     }
00466   }
00467 
00468   if ( mesnum==F2) {
00469 
00470     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00471        setProbMax(1100.0);
00472        return;
00473     }
00474     if ( lnum==TAUP||lnum==TAUM ) {
00475        setProbMax(1100.0);
00476        return;
00477     }
00478   }
00479 
00480   if ( mesnum==F2PR) {
00481 
00482     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00483        setProbMax(804.0);
00484        return;
00485     }
00486     if ( lnum==TAUP||lnum==TAUM ) {
00487        setProbMax(600.0);
00488        return;
00489     }
00490   }
00491 
00492   if ( mesnum==F1) {
00493 
00494     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00495        setProbMax(2500.0);
00496        return;
00497     }
00498     if ( lnum==TAUP||lnum==TAUM ) {
00499        setProbMax(2000.0) ;
00500        return;
00501     }
00502   }
00503 
00504   if ( mesnum==F1PR) {
00505 
00506     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00507        setProbMax(2400.0);
00508        return;
00509     }
00510     if ( lnum==TAUP||lnum==TAUM ) {
00511        setProbMax(1700.0);
00512        return;
00513     }
00514   }
00515 
00516   if ( mesnum==F0) {
00517 
00518     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00519        setProbMax( 80.0);
00520        return;
00521     }
00522     if ( lnum==TAUP||lnum==TAUM ) {
00523        setProbMax(63.0);
00524        return;
00525     }
00526   }
00527 
00528   if ( mesnum==F0PR) {
00529 
00530     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00531        setProbMax(120.0);
00532        return;
00533     }
00534     if ( lnum==TAUP||lnum==TAUM ) {
00535        setProbMax(120.0);
00536        return;
00537     }
00538   }
00539 
00540 
00541   if ( mesnum==RHO2SP||mesnum==RHO2SM||mesnum==RHO2S0) {
00542 
00543     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00544        setProbMax( 2400.0);
00545        return;
00546     }
00547     if ( lnum==TAUP||lnum==TAUM ) {
00548        setProbMax(2000.0);
00549        return;
00550     }
00551   }
00552 
00553   if ( mesnum==OMEG2S) {
00554 
00555     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00556        setProbMax(1600.0);
00557        return;
00558     }
00559     if ( lnum==TAUP||lnum==TAUM ) {
00560        setProbMax(1400.0) ;
00561        return;
00562     }
00563   }
00564 
00565   if ( mesnum==PI2SP||mesnum==PI2SM||mesnum==PI2S0) {
00566 
00567     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00568        setProbMax( 500.0);
00569        return;
00570     }
00571     if ( lnum==TAUP||lnum==TAUM ) {
00572        setProbMax(300.0);
00573        return;
00574     }
00575   }
00576 
00577   if ( mesnum==ETA2S) {
00578 
00579     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00580        setProbMax(344.0);
00581        return;
00582     }
00583     if ( lnum==TAUP||lnum==TAUM ) {
00584        setProbMax(300.0);
00585        return;
00586     }
00587   }
00588 
00589   if ( mesnum==KP||mesnum==KM||
00590        mesnum==K1P||mesnum==K1M||mesnum==K1STP||mesnum==K1STM) {
00591 
00592     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00593       setProbMax(2000.0);
00594       return;
00595     }
00596     if ( lnum==TAUP||lnum==TAUM ) {
00597       setProbMax(1000.0);
00598       return;
00599     }
00600   }
00601 
00602   if ( mesnum==KSTP||mesnum==KSTM ) {
00603 
00604     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00605       setProbMax(10000.0);
00606       return;
00607     }
00608     if ( lnum==TAUP||lnum==TAUM ) {
00609       setProbMax(7000.0);
00610       return;
00611     }
00612   }
00613 
00614 
00615 }
00616 
00617 if ( parnum==D0||parnum==DP||parnum==DM||parnum==D0B ) {
00618 
00619 
00620   if ( mesnum==RHOP||mesnum==RHOM||mesnum==RHO0) {
00621 
00622     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00623        setProbMax(110.0);
00624        return;
00625     }
00626   }
00627 
00628   if ( mesnum==OMEG) {
00629 
00630     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00631        setProbMax(75.0);
00632        return;
00633     }
00634   }
00635 
00636   if ( mesnum==PIP||mesnum==PIM||mesnum==PI0) {
00637 
00638     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00639        setProbMax(40.0);
00640        return;
00641     }
00642   }
00643 
00644   if ( mesnum==ETA) {
00645 
00646     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00647        setProbMax( 65.0);
00648        return;
00649     }
00650   }
00651 
00652   if ( mesnum==ETAPR) {
00653 
00654     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00655        setProbMax( 60.0);
00656        return;
00657     }
00658   }
00659 
00660   if ( mesnum==KP||mesnum==KM||mesnum==K0||
00661        mesnum==K0S||mesnum==K0L||mesnum==KB) {
00662 
00663     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00664        setProbMax( 70.0);
00665        return;
00666     }
00667   }
00668 
00669   if ( mesnum==K1STP||mesnum==K1STM||mesnum==K1ST0||mesnum==K1STB) {
00670 
00671     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00672        setProbMax( 3.3);
00673        return;
00674     }
00675   }
00676 
00677   if ( mesnum==K1P||mesnum==K1M||mesnum==K10||mesnum==K1B) {
00678 
00679     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00680        setProbMax( 100.0);
00681        return;
00682     }
00683   }
00684 
00685   if ( mesnum==KSTP||mesnum==KSTM||mesnum==KST0||mesnum==KSTB) {
00686 
00687     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00688        setProbMax( 135.0);
00689        return;
00690     }
00691   }
00692 
00693   if ( mesnum==K2STP||mesnum==K2STM||mesnum==K2ST0||mesnum==K2STB) {
00694 
00695     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00696       //Lange - Oct 26,2001 - increasing from 0.75 to 
00697       //accomodate
00698       setProbMax( 9.0);
00699       // setProbMax( 0.75);
00700        return;
00701     }
00702   }
00703 
00704   if ( mesnum==F0) {
00705     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00706       setProbMax(1.0);
00707       return;
00708     }
00709   }
00710 
00711 }
00712 
00713 if ( parnum==DSP||parnum==DSM ) {
00714 
00715 
00716   if ( mesnum==PHI ) {
00717 
00718     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00719        setProbMax( 90.0 );
00720        return;
00721     }
00722   }
00723 
00724   if ( mesnum==ETA ) {
00725 
00726     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00727       setProbMax( 75.0 );
00728        return;
00729     }
00730   }
00731 
00732   if ( mesnum==ETAPR ) {
00733 
00734     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00735       setProbMax( 80.0) ;
00736        return;
00737     }
00738   }
00739 
00740   if ( mesnum==KST0||mesnum==KSTB ) {
00741 
00742     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00743       setProbMax( 100.0) ;
00744        return;  
00745     }
00746   }
00747 
00748 
00749   if ( mesnum==K0 || mesnum==KB || mesnum==K0S || mesnum==K0L ) {
00750 
00751     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00752        setProbMax( 45.0 );
00753        return;
00754     }
00755   }
00756 
00757   if ( mesnum==F0) {
00758     if ( lnum==EP||lnum==EM||lnum==MUP||lnum==MUM ) {
00759       setProbMax(1.0);
00760       return;
00761     }
00762   }
00763 
00764 
00765 }
00766 
00767  if ( parnum==BCP||parnum==BCM ) {
00768    setProbMax(1000.0 );
00769    return;
00770  }
00771 
00772 
00773 
00774 //This is a real cludge.. (ryd)
00775  setProbMax(0.0);
00776 
00777 }

void EvtDecayAmp::makeDecay (  EvtParticle *  p  )  [virtual, inherited]
 

Implements EvtDecayBase.

00035                                          {
00036 
00037   int ntimes=10000;
00038 
00039   int i,more;
00040 
00041   EvtSpinDensity rho;
00042   double prob,prob_max;
00043 
00044   _amp2.init(p->getId(),getNDaug(),getDaugs());
00045 //  report(INFO,"EvtGen") << "Decaying " << EvtPDL::name(p->getId()) << endl;
00046   do{
00047 
00048     _daugsDecayedByParentModel=false;
00049     _weight = 1.0;
00050     decay(p);
00051 
00052     rho=_amp2.getSpinDensity();
00053     prob=p->getSpinDensityForward().NormalizedProb(rho);
00054     if (prob<0.0) {
00055 
00056       report(ERROR,"EvtGen")<<"Negative prob:"<<p->getId().getId()
00057                             <<" "<<p->getChannel()<<endl;
00058 
00059       report(ERROR,"EvtGen") << "rho_forward:"<<endl;
00060       report(ERROR,"EvtGen") << p->getSpinDensityForward();
00061       report(ERROR,"EvtGen") << "rho decay:"<<endl;
00062       report(ERROR,"EvtGen") << rho <<endl;
00063     }
00064 
00065     if (prob!=prob) {
00066 
00067       report(DEBUG,"EvtGen") << "Forward density matrix:"<<endl;
00068       report(ERROR,"EvtGen") << p->getSpinDensityForward();
00069 
00070       report(DEBUG,"EvtGen") << "Decay density matrix:"<<endl;
00071       report(ERROR,"EvtGen") << rho;
00072 
00073       report(DEBUG,"EvtGen") << "prob:"<<prob<<endl;
00074       
00075       report(DEBUG,"EvtGen") << "Particle:"
00076                              <<EvtPDL::name(p->getId()).c_str()<<endl;
00077       report(DEBUG,"EvtGen") << "channel        :"<<p->getChannel()<<endl;
00078       report(DEBUG,"EvtGen") << "Momentum:" << p->getP4() << " " << p->mass() << endl;
00079       if( p->getParent()!=0){
00080         report(DEBUG,"EvtGen") << "parent:"
00081                                <<EvtPDL::name(
00082                                 p->getParent()->getId()).c_str()<<endl;
00083         report(DEBUG,"EvtGen") << "parent channel        :"
00084                                <<p->getParent()->getChannel()<<endl;
00085 
00086         int i;
00087         report(DEBUG,"EvtGen") << "parent daughters  :";
00088         for (i=0;i<p->getParent()->getNDaug();i++){
00089           report(DEBUG,"") << EvtPDL::name(
00090                             p->getParent()->getDaug(i)->getId()).c_str()
00091                                  << " ";
00092         }
00093         report(DEBUG,"") << endl;
00094 
00095         report(DEBUG,"EvtGen") << "daughters  :";
00096         for (i=0;i<p->getNDaug();i++){
00097           report(DEBUG,"") << EvtPDL::name(
00098                             p->getDaug(i)->getId()).c_str()
00099                                  << " ";
00100         }
00101         report(DEBUG,"") << endl;
00102 
00103         report(DEBUG,"EvtGen") << "daughter momenta  :" << endl;;
00104         for (i=0;i<p->getNDaug();i++){
00105           report(DEBUG,"") << p->getDaug(i)->getP4() << " " << p->getDaug(i)->mass();
00106           report(DEBUG,"") << endl;
00107         }
00108 
00109       }
00110     }
00111 
00112 
00113     prob/=_weight;
00114 
00115 
00116     prob_max = getProbMax(prob);
00117     p->setDecayProb(prob/prob_max);
00118 
00119     //report(INFO,"EvtGen") << "Prob,prob_max,weight:"<<prob<<" "<<prob_max<<" "<<_weight<<endl;
00120 
00121     more=prob<EvtRandom::Flat(prob_max);
00122 
00123     ntimes--;
00124 
00125   }while(ntimes&&more);
00126   //report(INFO,"EvtGen") << "Done\n";
00127 
00128   if (ntimes==0){
00129     report(DEBUG,"EvtGen") << "Tried accept/reject:10000"
00130                            <<" times, and rejected all the times!"<<endl;
00131     report(DEBUG,"")<<p->getSpinDensityForward()<<endl;
00132     report(DEBUG,"EvtGen") << "Is therefore accepting the last event!"<<endl;
00133     report(DEBUG,"EvtGen") << "Decay of particle:"<<
00134       EvtPDL::name(p->getId()).c_str()<<"(channel:"<<
00135       p->getChannel()<<") with mass "<<p->mass()<<endl;
00136     
00137     int ii;
00138     for(ii=0;ii<p->getNDaug();ii++){
00139       report(DEBUG,"EvtGen") <<"Daughter "<<ii<<":"<<
00140         EvtPDL::name(p->getDaug(ii)->getId()).c_str()<<" with mass "<<
00141         p->getDaug(ii)->mass()<<endl;
00142     }                              
00143   }
00144 
00145 
00146   EvtSpinDensity rho_list[10];
00147 
00148   rho_list[0]=p->getSpinDensityForward();
00149    
00150 
00151   EvtAmp ampcont;
00152 
00153   if (_amp2._pstates!=1){
00154     ampcont=_amp2.contract(0,p->getSpinDensityForward());   // J2BB2 bugging here
00155 
00156   }
00157   else{
00158     ampcont=_amp2;
00159   }
00160 
00161           
00162 
00163   // it may be that the parent decay model has already
00164   // done the decay - this should be rare and the
00165   // model better know what it is doing..
00166   
00167   if ( !daugsDecayedByParentModel() ){
00168     
00169 
00170  //   report(INFO,"EvtGen") << "Found " << p->getNDaug() << " daughters\n";
00171     for(i=0;i<p->getNDaug();i++){
00172     
00173       rho.SetDim(_amp2.dstates[i]);
00174 
00175       if (_amp2.dstates[i]==1) {
00176         rho.Set(0,0,EvtComplex(1.0,0.0));
00177       }
00178       else{
00179         rho=ampcont.contract(_amp2._dnontrivial[i],_amp2);
00180       }
00181       
00182       if (!rho.Check()) {
00183         
00184         report(ERROR,"EvtGen") << "-------start error-------"<<endl;
00185         report(ERROR,"EvtGen")<<"forward rho failed Check:"<<
00186           EvtPDL::name(p->getId()).c_str()<<" "<<p->getChannel()<<" "<<i<<endl;
00187         
00188         report(ERROR,"EvtGen")<<"Parent:"<<EvtPDL::name(p->getParent()->getId()).c_str()<<endl;
00189         report(ERROR,"EvtGen")<<"GrandParent:"<<EvtPDL::name(p->getParent()->getParent()->getId()).c_str()<<endl;
00190         report(ERROR,"EvtGen")<<"GrandGrandParent:"<<EvtPDL::name(p->getParent()->getParent()->getParent()->getId()).c_str()<<endl;
00191         
00192         report(ERROR,"EvtGen") << rho;
00193         int ii; 
00194         _amp2.dump();
00195         
00196         for(ii=0;ii<i+1;ii++){
00197           report(ERROR,"EvtGen") << rho_list[ii];
00198         }
00199         report(ERROR,"EvtGen") << "-------Done with error-------"<<endl;  
00200       }
00201 
00202 
00203       
00204       p->getDaug(i)->setSpinDensityForward(rho);
00205 
00206 
00207      p->getDaug(i)->decay();
00208 
00209      rho_list[i+1]= p->getDaug(i)->getSpinDensityBackward();
00210 
00211       if (_amp2.dstates[i]!=1){
00212         ampcont=ampcont.contract(_amp2._dnontrivial[i],rho_list[i+1]);
00213       }
00214   
00215       
00216     }
00217 
00218     p->setSpinDensityBackward(_amp2.getBackwardSpinDensity(rho_list));
00219    
00220     
00221     if (!p->getSpinDensityBackward().Check()) {
00222       
00223       report(ERROR,"EvtGen")<<"rho_backward failed Check"<<
00224         p->getId().getId()<<" "<<p->getChannel()<<endl;
00225       
00226       report(ERROR,"EvtGen") << p->getSpinDensityBackward();
00227       
00228     }
00229   }
00230 
00231   if (getPHOTOS() || EvtRadCorr::alwaysRadCorr()) {
00232     EvtRadCorr::doRadCorr(p);
00233 
00234   }
00235 
00236 }

bool EvtDecayBase::matchingDecay (  const EvtDecayBase &  other  )  const [virtual, inherited]
  

00588                                                                 {
00589 
00590   if ( _ndaug != other._ndaug) return false;
00591   if ( _parent != other._parent) return false;
00592   
00593   std::vector<int> useDs;
00594   for ( unsigned int i=0; i<_ndaug; i++) useDs.push_back(0);
00595 
00596   for ( unsigned int i=0; i<_ndaug; i++) {
00597     bool foundIt=false;
00598     for ( unsigned int j=0; j<_ndaug; j++) {
00599       if ( useDs[j] == 1 ) continue;
00600       if ( _daug[i] == other._daug[j] && _daug[i].getAlias() == other._daug[j].getAlias()) {
00601         foundIt=true;
00602         useDs[j]=1;
00603         break;
00604       }
00605     }
00606     if ( foundIt==false) return false;
00607   }
00608   for ( unsigned int i=0; i<_ndaug; i++) if ( useDs[i]==0) return false;
00609 
00610   return true;
00611 
00612 }

void EvtDecayBase::noProbMax (   )  [inherited]
  

00305                             {
00306 
00307   defaultprobmax=0;
00308 
00309 }

virtual int EvtDecayBase::nRealDaughters (   )  [inline, virtual, inherited]
 

Reimplemented in EvtBtoKD3P, and EvtVSSBMixCPT.

00105 { return _ndaug;}

void EvtDecayBase::printSummary (   )  [inherited]
  

00259                                 {
00260 
00261   int i;
00262 
00263   if (ntimes_prob>0) {
00264 
00265     report(INFO,"EvtGen") << "Calls="<<ntimes_prob<<" eff:"<<
00266       sum_prob/(probmax*ntimes_prob)<<" frac. max:"<<max_prob/probmax;
00267     report(INFO,"") <<" probmax:"<<probmax<<" max:"<<max_prob<<" : ";
00268   }
00269 
00270   report(INFO,"") << EvtPDL::name(_parent).c_str()<<" -> ";
00271   for(i=0;i<_ndaug;i++){
00272     report(INFO,"") << EvtPDL::name(_daug[i]).c_str() << " ";
00273   }
00274   report(INFO,"") << " ("<<_modelname.c_str()<<"):"<< endl;
00275   
00276   
00277   
00278 }

double EvtDecayBase::resetProbMax (  double  prob  )  [inherited]
  

00108                                              {
00109   
00110   report(INFO,"EvtGen") << "Reseting prob max\n"; 
00111   report(INFO,"EvtGen") << "prob > probmax:("<<prob<<">"<<probmax<<")";
00112   report(INFO,"") << "("<<_modelname.c_str()<<")";
00113   report(INFO,"") << EvtPDL::getStdHep(_parent)<<"->";
00114   
00115   for( int i=0;i<_ndaug;i++){
00116     report(INFO,"") << EvtPDL::getStdHep(_daug[i]) << " ";
00117   }
00118   report(INFO,"") << endl;
00119   
00120   probmax = 0.0;
00121   defaultprobmax = 0;
00122   ntimes_prob = 0;
00123   
00124   return prob;
00125 
00126 }

void EvtDecayBase::saveDecayInfo (  EvtId  ipar,
int  ndaug,
EvtId *  daug,
int  narg,
std::vector< std::string > &  args,
std::string  name,
double  brfr
)  [inherited]
  

00170                                               {
00171 
00172   int i;
00173 
00174   _brfr=brfr;
00175   _ndaug=ndaug;
00176   _narg=narg;
00177   _parent=ipar; 
00178 
00179   _dsum=0;
00180 
00181   if (_ndaug>0) {
00182     _daug=new EvtId [_ndaug];
00183     for(i=0;i<_ndaug;i++){
00184       _daug[i]=daug[i];
00185       _dsum+=daug[i].getAlias();
00186     }
00187   }
00188   else{
00189     _daug=0;
00190   }
00191 
00192   if (_narg>0) {
00193     _args=new std::string[_narg+1];
00194     for(i=0;i<_narg;i++){
00195       _args[i]=args[i];
00196     }
00197   }
00198   else{
00199      _args = 0;
00200   }
00201 
00202   _modelname=name;
00203 
00204   this->init();
00205   this->initProbMax();
00206 
00207   if (_chkCharge){
00208     this->checkQ();
00209   }
00210 
00211 
00212   if (defaultprobmax){
00213     report(INFO,"EvtGen") << "No default probmax for ";
00214     report(INFO,"") << "("<<_modelname.c_str()<<") ";
00215     report(INFO,"") << EvtPDL::name(_parent).c_str()<<" -> ";
00216     for(i=0;i<_ndaug;i++){
00217       report(INFO,"") << EvtPDL::name(_daug[i]).c_str() << " ";
00218     }
00219     report(INFO,"") << endl;
00220     report(INFO,"") << "This is fine for development, but must be provided for production."<<endl;
00221     report(INFO,"EvtGen") << "Never fear though - the decay will use the \n";
00222     report(INFO,"EvtGen") << "500 iterations to build up a good probmax \n";
00223     report(INFO,"EvtGen") << "before accepting a decay. "<<endl;
00224   }
00225 
00226 }

void EvtDecayBase::setPHOTOS (   )  [inline, inherited]
  

00069 {_photos=1;}

void EvtDecayBase::setProbMax (  double  prbmx  )  [inherited]
  

00298                                          {
00299 
00300   defaultprobmax=0;
00301   probmax=prbmx;
00302 
00303 }

void EvtDecayBase::setSummary (   )  [inline, inherited]
  

00071 {_summary=1;}

void EvtDecayBase::setVerbose (   )  [inline, inherited]
  

00070 {_verbose=1;}

void EvtDecayAmp::setWeight (  double  weight  )  [inline, inherited]
  

00032 {_weight=weight;}

int EvtDecayBase::summary (   )  [inline, inherited]
  

00078 {return _summary; }

int EvtDecayBase::verbose (   )  [inline, inherited]
  

00079 {return _verbose; }

void EvtDecayAmp::vertex (  int *  i1,
const EvtComplex &  amp
)  [inline, inherited]
 

sets the amplitudes calculated in the decay objects 

00061   { _amp2.vertex(i1,amp);}

void EvtDecayAmp::vertex (  int  i1,
int  i2,
int  i3,
const EvtComplex &  amp
)  [inline, inherited]
 

sets the amplitudes calculated in the decay objects 

00055   {_amp2.vertex(i1,i2,i3,amp);}

void EvtDecayAmp::vertex (  int  i1,
int  i2,
const EvtComplex &  amp
)  [inline, inherited]
 

sets the amplitudes calculated in the decay objects 

00048   {_amp2.vertex(i1,i2,amp);}

void EvtDecayAmp::vertex (  int  i1,
const EvtComplex &  amp
)  [inline, inherited]
 

sets the amplitudes calculated in the decay objects 

00042 {_amp2.vertex(i1,amp);}

void EvtDecayAmp::vertex (  const EvtComplex &  amp  )  [inline, inherited]
 

sets the amplitudes calculated in the decay objects 

00037 {_amp2.vertex(amp);}

Member Data Documentation

EvtAmp EvtDecayAmp::_amp2 [protected, inherited]
 

bool EvtDecayBase::_daugsDecayedByParentModel [protected, inherited]
 

EvtSemiLeptonicAmp* EvtISGW2::calcamp [private]
 

EvtSemiLeptonicFF* EvtISGW2::isgw2ffmodel [private]
 

The documentation for this class was generated from the following files:
Generated on Wed Feb 2 16:06:12 2011 for BOSS6.5.5 by  doxygen 1.3.9.1