TEvtGen/EvtGenModels/EvtBToKpipiCP.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: EvtBToKpipiCP.cc
  12 //
  13 // Description: Routine to decay B->K pi pi
  14 //              and has CP violation.
  15 //       --- This is the routine to be called by the Main generator
  16 //          to get the decay of B0    -->-- K+ pi- pi0
  17 //          The decay proceeeds through three channels:
  18 //          a) B0 -->-- K*+ pi-  ; K*+    -->-- K+ pi0
  19 //          b)          K*0 pi0  ; K*0bar -->-- K+ pi-
  20 //          c)          K-  rho+ ; rho+   -->-- pi+ pi0
  21 //         It provides at the same time the CP conjugate decay
  22 //                              B0bar -->-- K- pi+ pi0
  23 //
  24 // Modification history:
  25 //
  26 //    Versille     September, 1997         Module created
  27 //
  28 //------------------------------------------------------------------------
  29 //
  30 #include "EvtGenBase/EvtPatches.hh"
  31 #include <stdlib.h>
  32 #include "EvtGenBase/EvtParticle.hh"
  33 #include "EvtGenBase/EvtGenKine.hh"
  34 #include "EvtGenBase/EvtCPUtil.hh"
  35 #include "EvtGenBase/EvtPDL.hh"
  36 #include "EvtGenBase/EvtReport.hh"
  37 #include "EvtGenModels/EvtBToKpipiCP.hh"
  38 #include "EvtGenBase/EvtId.hh"
  39 #include <string>
  40
  41 #ifdef WIN32
  42 extern "C" {
  43   extern void __stdcall EVTKPIPI(double *, double *, int *,double *,
  44                                  double *,double *,double *,double *,
  45                                  double *,double *,double *);
  46 }
  47 #else
  48 extern "C" {
  49   extern void evtkpipi_(double *, double *, int *,double *,
  50                         double *,double *,double *,double *,
  51                         double *,double *,double *);
  52 }
  53 #endif
  54
  55 EvtBToKpipiCP::~EvtBToKpipiCP() {}
  56
  57
  58 std::string EvtBToKpipiCP::getName(){
  59
  60   return "BTOKPIPI_CP";
  61
  62 }
  63
  64
  65 EvtDecayBase* EvtBToKpipiCP::clone(){
  66
  67   return new EvtBToKpipiCP;
  68
  69 }
  70
  71 void EvtBToKpipiCP::init(){
  72
  73   // check that there are 3 arguments
  74   checkNArg(3);
  75   checkNDaug(3);
  76
  77   checkSpinParent(EvtSpinType::SCALAR);
  78
  79   checkSpinDaughter(0,EvtSpinType::SCALAR);
  80   checkSpinDaughter(1,EvtSpinType::SCALAR);
  81   checkSpinDaughter(2,EvtSpinType::SCALAR);
  82
  83 }
  84
  85
  86 void EvtBToKpipiCP::decay( EvtParticle *p){
  87
  88   //added by Lange Jan4,2000
  89   static EvtId B0=EvtPDL::getId("B0");
  90   static EvtId B0B=EvtPDL::getId("anti-B0");
  91
  92   double t;
  93   EvtId other_b;
  94
  95   EvtCPUtil::OtherB(p,t,other_b);
  96
  97   EvtParticle *Kp,*pim,*pi0;
  98
  99   p->makeDaughters(getNDaug(),getDaugs());
 100   Kp=p->getDaug(0);
 101   pim=p->getDaug(1);
 102   pi0=p->getDaug(2);
 103
 104   EvtVector4R p4[3];
 105
 106   //double dm=getArg(0);
 107   double alpha=getArg(1);
 108   double beta=getArg(2);
 109   int iset;
 110
 111   static int first=1;
 112
 113   if (first==1) {
 114     iset=10000;
 115     first=0;
 116   }
 117   else{
 118     iset=0;
 119   }
 120
 121   double p4Kplus[4],p4piminus[4],p4gamm1[4],p4gamm2[4];
 122
 123   double realA,imgA,realbarA,imgbarA;
 124
 125 #ifdef WIN32
 126     EVTKPIPI(&alpha,&beta,&iset,p4Kplus,p4piminus,p4gamm1,p4gamm2,
 127              &realA,&imgA,&realbarA,&imgbarA);
 128 #else
 129   evtkpipi_(&alpha,&beta,&iset,p4Kplus,p4piminus,p4gamm1,p4gamm2,
 130             &realA,&imgA,&realbarA,&imgbarA);
 131 #endif
 132
 133   p4[0].set(p4Kplus[3],p4Kplus[0],p4Kplus[1],p4Kplus[2]);
 134   p4[1].set(p4piminus[3],p4piminus[0],p4piminus[1],p4piminus[2]);
 135   p4[2].set(p4gamm1[3]+p4gamm2[3],p4gamm1[0]+p4gamm2[0],
 136             p4gamm1[1]+p4gamm2[1],p4gamm1[2]+p4gamm2[2]);
 137
 138    Kp->init( getDaug(0), p4[0] );
 139    pim->init( getDaug(1), p4[1] );
 140    pi0->init( getDaug(2), p4[2] );
 141
 142    EvtComplex amp;
 143
 144    EvtComplex A(realA,imgA);
 145    EvtComplex Abar(realbarA,imgbarA);
 146
 147    if (other_b==B0B){
 148      amp=Abar;
 149    }
 150    if (other_b==B0){
 151      amp=A;
 152    }
 153
 154    vertex(amp);
 155
 156   return ;
 157 }
 158