TEvtGen/EvtGenModels/EvtMultibody.cpp

   1 #include "EvtGenBase/EvtPatches.hh"
   2 #include "EvtGenBase/EvtGenKine.hh"
   3 #include "EvtGenBase/EvtPDL.hh"
   4 #include "EvtGenBase/EvtReport.hh"
   5 #include "EvtGenBase/EvtResonance.hh"
   6 #include "EvtGenBase/EvtResonance2.hh"
   7 #include "EvtGenModels/EvtMultibody.hh"
   8 #include "EvtGenBase/EvtConst.hh"
   9 #include "EvtGenBase/EvtdFunction.hh"
  10 #include "EvtGenBase/EvtKine.hh"
  11 #include "EvtGenBase/EvtParticle.hh"
  12
  13 EvtMultibody::~EvtMultibody()
  14 {
  15     if( _decayTree != NULL ) delete _decayTree;
  16     _decayTree=NULL;
  17     if( _ilist != NULL ) delete [] _ilist;
  18     _ilist=NULL;
  19 }
  20
  21 std::string EvtMultibody::getName()
  22 {
  23     return "D_MULTIBODY";
  24 }
  25
  26 EvtDecayBase* EvtMultibody::clone()
  27 {
  28     return new EvtMultibody;
  29 }
  30
  31 void EvtMultibody::init()
  32 {
  33     int N = getNArg();
  34
  35     _decayTree = new EvtMTree( getDaugs(), getNDaug() );
  36     _ilist = new int[getNDaug()+1];
  37
  38     for(int i=0; i<N-1; ++i) {
  39         if(getArgStr( i )=="RESONANCE") {
  40             _decayTree->addtree( getArgStr( ++i ) );
  41         } else {
  42             report(ERROR,"EvtGen")
  43                 << "Syntax error at " << getArgStr( i ) << std::endl;
  44             ::abort();
  45         }
  46     }
  47 }
  48
  49 // Set the maximum probability amplitude - if function is left blank then the
  50 // program will search for it.  This however is not deterministic and therefore
  51 // in the release cannot be in place.
  52 void EvtMultibody::initProbMax()
  53 {
  54     // setProbMax(1.0);
  55 }
  56
  57 void EvtMultibody::decay( EvtParticle *p )
  58 {
  59     // Initialize the phase space before doing anything else!
  60     p->initializePhaseSpace(getNDaug(),getDaugs());
  61
  62     EvtSpinAmp amp = _decayTree->amplitude( p );
  63
  64     vector<int> index = amp.iterallowedinit();
  65     vector<unsigned int> spins = amp.dims();
  66
  67     do {
  68         for( size_t i=0; i<index.size(); ++i ) {
  69             _ilist[i]=index[i]+spins[i];
  70         }
  71
  72         vertex( _ilist, amp( index ) );
  73     } while( amp.iterateallowed( index ) );
  74
  75 }