EvtBto2piCPiso.cc

Go to the documentation of this file.

//--------------------------------------------------------------------------
//
// Environment:
//      This software is part of the EvtGen package developed jointly
//      for the BaBar and CLEO collaborations.  If you use all or part
//      of it, please give an appropriate acknowledgement.
//
// Copyright Information: See EvtGen/COPYRIGHT
//      Copyright (C) 1998      Caltech, UCSB
//
// Module: EvtBto2piCPiso.cc
//
// Description: Routine to decay B -> pi pi with isospin amplitudes 
//
// Modification history:
//
//    RYD,NK     Febuary 7, 1998         Module created
//
//------------------------------------------------------------------------
//
#include "EvtGenBase/EvtPatches.hh"
#include <stdlib.h>
#include "EvtGenBase/EvtParticle.hh"
#include "EvtGenBase/EvtRandom.hh"
#include "EvtGenBase/EvtGenKine.hh"
#include "EvtGenBase/EvtCPUtil.hh"
#include "EvtGenBase/EvtPDL.hh"
#include "EvtGenModels/EvtBto2piCPiso.hh"
#include "EvtGenBase/EvtReport.hh"
#include "EvtGenBase/EvtId.hh"
#include <string>
#include "EvtGenBase/EvtConst.hh"
 
EvtBto2piCPiso::~EvtBto2piCPiso() {}
 
void EvtBto2piCPiso::getName(std::string& model_name){
 
  model_name="BTO2PI_CP_ISO";     
 
}
 
 
EvtDecayBase* EvtBto2piCPiso::clone(){
 
  return new EvtBto2piCPiso;
 
}
 
void EvtBto2piCPiso::init(){
 
  // check that there are 11 arguments
 
  checkNArg(11);
  checkNDaug(2);
 
  checkSpinParent(EvtSpinType::SCALAR);
 
  checkSpinDaughter(0,EvtSpinType::SCALAR);
  checkSpinDaughter(1,EvtSpinType::SCALAR);
 
}
 
 
void EvtBto2piCPiso::initProbMax() {
 
  //added by Lange Jan4,2000
  static EvtId PI0=EvtPDL::getId("pi0");
  static EvtId PIP=EvtPDL::getId("pi+");
  static EvtId PIM=EvtPDL::getId("pi-");
 
//this may need to be revised
 
if (((getDaugs()[0]==PIP) && (getDaugs()[1]==PIM)) || ((getDaugs()[0]==PIM) && (getDaugs()[1]==PIP))) {
   setProbMax(4.0*(getArg(2)*getArg(2)+getArg(4)*getArg(4)));
}
 
if ((getDaugs()[0]==PI0) && (getDaugs()[1]==PI0)) {
   setProbMax(2.0*(4.0*getArg(2)*getArg(2)+getArg(4)*getArg(4)));
}
 
if (((getDaugs()[0]==PIP) && (getDaugs()[1]==PI0)) || ((getDaugs()[0]==PI0) && (getDaugs()[1]==PIP))) {
   setProbMax(6.0*getArg(2)*getArg(2));
}
 
if (((getDaugs()[0]==PI0) && (getDaugs()[1]==PIM)) || ((getDaugs()[0]==PIM) && (getDaugs()[1]==PI0))) {
   setProbMax(6.0*getArg(4)*getArg(4));
}
 
}
 
void EvtBto2piCPiso::decay( EvtParticle *p ){
 
  p->initializePhaseSpace(getNDaug(),getDaugs());
 
  //added by Lange Jan4,2000
  static EvtId B0=EvtPDL::getId("B0");
  static EvtId B0B=EvtPDL::getId("anti-B0");
  static EvtId PI0=EvtPDL::getId("pi0");
  static EvtId PIP=EvtPDL::getId("pi+");
  static EvtId PIM=EvtPDL::getId("pi-");
 
  double t;
  EvtId other_b;
  int charged=0;
 
//randomly generate the tag (B0 or B0B) 
 
  double tag = EvtRandom::Flat(0.0,1.0);
  if (tag < 0.5) {
 
   EvtCPUtil::OtherB(p,t,other_b,1.0);
   other_b = B0;
  }
  else {
   
   EvtCPUtil::OtherB(p,t,other_b,0.0);
   other_b = B0B;
  }
 
  EvtComplex amp;
 
  EvtComplex A,Abar;
  EvtComplex A2, A2_bar, A0, A0_bar;
 
  A2 = EvtComplex(getArg(2)*cos(getArg(3)),getArg(2)*sin(getArg(3)));
  A2_bar = EvtComplex(getArg(4)*cos(getArg(5)),getArg(4)*sin(getArg(5)));
  
  A0 = EvtComplex(getArg(6)*cos(getArg(7)),getArg(6)*sin(getArg(7)));
  A0_bar = EvtComplex(getArg(8)*cos(getArg(9)),getArg(8)*sin(getArg(9)));
  
  //depending on what combination of pi pi we have, there will be different 
  //A and Abar
 
  if (((getDaugs()[0]==PIP) && (getDaugs()[1]==PI0)) || ((getDaugs()[0]==PI0) && (getDaugs()[1]==PIP))) {
 
    //pi+ pi0, so just A_2
  
    charged = 1;
    A = 3.0*A2;
    
  }
 
  if (((getDaugs()[0]==PI0) && (getDaugs()[1]==PIM)) || ((getDaugs()[0]==PIM) && (getDaugs()[1]==PI0))) {
 
    //pi- pi0, so just A2_bar
    
    charged = 1;  
    A = 3.0*A2_bar;
    
  }
 
  if (((getDaugs()[0]==PIP) && (getDaugs()[1]==PIM)) || ((getDaugs()[0]==PIM) && (getDaugs()[1]==PIP))) {
    
//pi+ pi-, so A_2 - A_0
 
    charged = 0;
    A=sqrt(2.0)*(A2-A0);
    Abar=sqrt(2.0)*(A2_bar-A0_bar);
  }
  
  if ((getDaugs()[0]==PI0) && (getDaugs()[1]==PI0)) {
    
    //pi0 pi0, so 2*A_2 + A_0
    
    charged = 0;
    A=2.0*A2 + A0;
    Abar=2.0*A2_bar + A0_bar;
  }
  
  if(charged == 0) {
    
    if (other_b==B0B){
      amp=A*cos(getArg(1)*t/(2*EvtConst::c))+
    EvtComplex(cos(-2.0*getArg(0)),sin(-2.0*getArg(0)))*
    EvtComplex(0.0,1.0)*Abar*sin(getArg(1)*t/(2*EvtConst::c));
    }
    if (other_b==B0){
      amp=A*EvtComplex(cos(2.0*getArg(0)),sin(2.0*getArg(0)))*
    EvtComplex(0.0,1.0)*sin(getArg(1)*t/(2*EvtConst::c))+       
    Abar*cos(getArg(1)*t/(2*EvtConst::c));
    }
  }
  else amp = A;
  
  vertex(amp);
  
  return ;
}