EvtEtap2pipieta.cc

Go to the documentation of this file.

//--------------------------------------------------------------------------
//
// Environment:
//      This software is part of models developed at BES collaboration
//      based on the EvtGen framework.  If you use all or part
//      of it, please give an appropriate acknowledgement.
//
// Copyright Information: See EvtGen/BesCopyright
//      Copyright (A) 2006      Ping Rong-Gang @IHEP
//
// Module:  EvtEtap2pipieta.cc
//
// Description: etaprime -> pi+pi- eta
// see PHYSICAL REVIEW D 97, 012003 (2018)
//
// Modification history:
// Jun 03 00:42:10 2024  Liaoyuan Dong  Module created
//------------------------------------------------------------------------
 
#include "EvtGenBase/EvtPatches.hh"
#include <stdlib.h>
#include "EvtGenBase/EvtParticle.hh"
#include "EvtGenBase/EvtGenKine.hh"
#include "EvtGenBase/EvtPDL.hh"
#include "EvtGenBase/EvtVector4R.hh"
#include "EvtGenBase/EvtReport.hh"
#include "EvtGenModels/EvtEtap2pipieta.hh"
#include <string>
using namespace std;  //::endl;
 
EvtEtap2pipieta::~EvtEtap2pipieta() {}
 
void EvtEtap2pipieta::getName(std::string& model_name){
  model_name="Etap2pipieta";
}
 
EvtDecayBase* EvtEtap2pipieta::clone(){
  return new EvtEtap2pipieta;
}
 
void EvtEtap2pipieta::init(){
  checkNDaug(3);
  checkSpinParent(EvtSpinType::SCALAR);
  EvtSpinType::spintype parenttype = EvtPDL::getSpinType(getParentId());
  m_MetaOverMpip = 0.547862/0.13957018; // m_Meta/m_Mpi
  m_sqrt3=sqrt(3.0);
  par[0]=-0.127;
  par[1]=-0.106;
  par[2]=1.5e-03;
  par[3]=-0.082;
}
 
void EvtEtap2pipieta::initProbMax(){
  setProbMax(1.03805);
}
 
void EvtEtap2pipieta::decay( EvtParticle *p ){
/*
  double maxprob = 0.0;
  for(int ir=0;ir<=60000000;ir++){
     p->initializePhaseSpace(getNDaug(),getDaugs());
     _pd[0]=p->getDaug(0)->getP4();
     _pd[1]=p->getDaug(1)->getP4();
     _pd[2]=p->getDaug(2)->getP4();
     double Prob = AMPsq();
     if(Prob>maxprob) {
         maxprob=Prob;
         std::cout << "Max PDF = " << ir << " prob= " << Prob << std::endl;
     }
  }
  std::cout << "Max!!!!!!!!!!! " << maxprob<< std::endl;
*/
  p->initializePhaseSpace(getNDaug(),getDaugs());
  _pd[0]=p->getDaug(0)->getP4();
  _pd[1]=p->getDaug(1)->getP4();
  _pd[2]=p->getDaug(2)->getP4();
  double prob = AMPsq();
  setProb(prob);
  return;
}
 
double EvtEtap2pipieta::AMPsq(){
  EvtVector4R petap=_pd[0]+_pd[1]+_pd[2];
  EvtVector4R Betap(petap.get(0),-1.0*petap.get(1),-1.0*petap.get(2),-1.0*petap.get(3));
 
  EvtVector4R pip_etapcms = boostTo(_pd[0],Betap);
  EvtVector4R pim_etapcms = boostTo(_pd[1],Betap);
  EvtVector4R eta_etapcms = boostTo(_pd[2],Betap);
 
  double Q=Betap.mass()-pip_etapcms.mass()-pim_etapcms.mass()-eta_etapcms.mass();
 
  double Tp=pip_etapcms.get(0)-pip_etapcms.mass();
  double Tm=pim_etapcms.get(0)-pim_etapcms.mass();
  double Te=eta_etapcms.get(0)-eta_etapcms.mass();
 
  double X=m_sqrt3*(Tp-Tm)/Q;
  double Y=(2.0+m_MetaOverMpip)*Te/Q-1.0;
 
  double amp2=1.0+par[0]*Y+par[1]*Y*Y+par[2]*X+par[3]*X*X;
  return amp2;
}