TEvtGen/EvtGenModels/EvtVVPIPI_WEIGHTED.cxx

   1 //--------------------------------------------------------------------------
   2 //
   3 // Environment:
   4 //      This software is part of the EvtGen package developed jointly
   5 //      for the BaBar and CLEO collaborations.  If you use all or part
   6 //      of it, please give an appropriate acknowledgement.
   7 //
   8 // Copyright Information: See EvtGen/COPYRIGHT
   9 //      Copyright (C) 1998      Caltech, UCSB
  10 //
  11 // Module: EvtVVSPwave.cc
  12 //
  13 // Description: Routine to decay vector-> vector pi pi where the
  14 //              decay is S-wave dominated.
  15 //
  16 // Modification history:
  17 //
  18 // Jim Hunt      June 4, 2008 Module Created
  19 //
  20 //------------------------------------------------------------------------
  21 //
  22 #include "EvtGenBase/EvtPatches.hh"
  23 #include <stdlib.h>
  24 #include "EvtGenBase/EvtParticle.hh"
  25 #include "EvtGenBase/EvtGenKine.hh"
  26 #include "EvtGenBase/EvtPDL.hh"
  27 #include "EvtGenBase/EvtVector4C.hh"
  28 #include "EvtGenBase/EvtTensor4C.hh"
  29 #include "EvtGenBase/EvtReport.hh"
  30 #include "EvtGenModels/EvtVVPIPI_WEIGHTED.hh"
  31 #include <string>
  32 using std::endl;
  33
  34 EvtVVPIPI_WEIGHTED::~EvtVVPIPI_WEIGHTED() {}
  35
  36 std::string EvtVVPIPI_WEIGHTED::getName(){
  37
  38   return "VVPIPI_WEIGHTED";
  39
  40 }
  41
  42
  43 EvtDecayBase* EvtVVPIPI_WEIGHTED::clone(){
  44
  45   return new EvtVVPIPI_WEIGHTED;
  46
  47 }
  48
  49 void EvtVVPIPI_WEIGHTED::init(){
  50
  51   static EvtId PIP=EvtPDL::getId("pi+");
  52   static EvtId PIM=EvtPDL::getId("pi-");
  53   static EvtId PI0=EvtPDL::getId("pi0");
  54
  55   // check that there are 0 arguments
  56   checkNArg(0);
  57   checkNDaug(3);
  58
  59   checkSpinParent(EvtSpinType::VECTOR);
  60   checkSpinDaughter(0,EvtSpinType::VECTOR);
  61
  62
  63
  64   if ((!(getDaug(1)==PIP&&getDaug(2)==PIM))&&
  65       (!(getDaug(1)==PI0&&getDaug(2)==PI0))) {
  66     report(ERROR,"EvtGen") << "EvtVVPIPI_WEIGHTED generator expected "
  67                            << " pi+ and pi- (or pi0 and pi0) "
  68                            << "as 2nd and 3rd daughter. "<<endl;
  69     report(ERROR,"EvtGen") << "Will terminate execution!"<<endl;
  70     ::abort();
  71   }
  72
  73 }
  74
  75 void EvtVVPIPI_WEIGHTED::initProbMax() {
  76   //Hard coded... should not be hard to calculate...
  77   setProbMax(0.08* 1.13 );
  78 }
  79
  80 double reweight_event(double pipi_mass)
  81 {
  82     pipi_mass *= 1000.0;
  83     return sqrt(
  84         -3.6911336508223251 +
  85         0.019119831948029617 * pipi_mass +
  86         -1.8962883732377376e-05 * pipi_mass * pipi_mass
  87         );
  88 }
  89
  90 void EvtVVPIPI_WEIGHTED::decay( EvtParticle *psi_prime){
  91
  92   psi_prime->initializePhaseSpace(getNDaug(),getDaugs());
  93
  94   EvtParticle *jpsi,*pi1,*pi2;
  95
  96   jpsi=psi_prime->getDaug(0);
  97   pi1=psi_prime->getDaug(1);
  98   pi2=psi_prime->getDaug(2);
  99
 100 //  Put phase space results into the daughters.
 101
 102   EvtVector4C ep0,ep1,ep2;
 103
 104   ep0=psi_prime->eps(0);
 105   ep1=psi_prime->eps(1);
 106   ep2=psi_prime->eps(2);
 107
 108   EvtVector4C e0,e1,e2;
 109
 110   e0 = jpsi->epsParent(0);
 111   e1 = jpsi->epsParent(1);
 112   e2 = jpsi->epsParent(2);
 113
 114   double mass2 = (pi1->getP4()+pi2->getP4()).mass2();
 115
 116   double fac = mass2-4*pi1->mass()*pi2->mass();
 117
 118   fac *= reweight_event(sqrt(mass2));
 119
 120   vertex(0,0,fac*(ep0*e0.conj()));
 121   vertex(0,1,fac*(ep0*e1.conj()));
 122   vertex(0,2,fac*(ep0*e2.conj()));
 123
 124   vertex(1,0,fac*(ep1*e0.conj()));
 125   vertex(1,1,fac*(ep1*e1.conj()));
 126   vertex(1,2,fac*(ep1*e2.conj()));
 127
 128   vertex(2,0,fac*(ep2*e0.conj()));
 129   vertex(2,1,fac*(ep2*e1.conj()));
 130   vertex(2,2,fac*(ep2*e2.conj()));
 131
 132   return ;
 133
 134 }
 135