Line data    Source code

       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             : //                    2011      University of Warwick, UK
      11             : // Module: EvtPythia.cc
      12             : //
      13             : // Description: Routine to decay a particle according th phase space
      14             : //
      15             : // Modification history:
      16             : //
      17             : //    RYD       January 8, 1997       Module created
      18             : //    JJB       April 2011            Modified to use new Pythia8 interface
      19             : //
      20             : //------------------------------------------------------------------------
      21             : 
      22             : #include "EvtGenBase/EvtPatches.hh"
      23             : 
      24             : #include "EvtGenBase/EvtParticle.hh"
      25             : #include "EvtGenBase/EvtId.hh"
      26             : #include "EvtGenBase/EvtPDL.hh"
      27             : #include "EvtGenBase/EvtSpinDensity.hh"
      28             : 
      29             : #include "EvtGenExternal/EvtPythia.hh"
      30             : 
      31             : #include "EvtGenExternal/EvtExternalGenFactory.hh"
      32             : #include "EvtGenModels/EvtAbsExternalGen.hh"
      33             : #include "EvtGenBase/EvtDecayBase.hh"
      34             : 
      35             : #include <iostream>
      36             : #include <cmath>
      37             : 
      38           0 : EvtPythia::EvtPythia() {
      39             : 
      40             :   // Set the Pythia engine to a null pointer at first.
      41             :   // When we do the decay, we retrieve the pointer to the Pythia engine
      42             :   // and use that for all decays. All clones will use the same Pythia engine.
      43           0 :   _pythiaEngine = 0;
      44             : 
      45           0 : }
      46             : 
      47           0 : EvtPythia::~EvtPythia() {
      48           0 :   _commandList.clear();
      49           0 : }
      50             : 
      51             : std::string EvtPythia::getName(){
      52             : 
      53           0 :   return "PYTHIA";
      54             : 
      55             : }
      56             : 
      57             : EvtDecayBase* EvtPythia::clone(){
      58             : 
      59           0 :   return new EvtPythia();
      60             : 
      61           0 : }
      62             : 
      63             : 
      64             : void EvtPythia::init(){
      65             : 
      66             :   // Do not check for any arguments. The PythiaEngine will check
      67             :   // to see if there is an integer specifying the decay physics,
      68             :   // otherwise it just uses phase-space.
      69             : 
      70           0 : }
      71             : 
      72             : void EvtPythia::initProbMax(){
      73             : 
      74           0 :   noProbMax();
      75             : 
      76           0 : }
      77             : 
      78             : void EvtPythia::decay( EvtParticle *p ){
      79             : 
      80             :   // We have to initialise the Pythia engine after the decay.dec files have been read in,
      81             :   // since we will be modifying Pythia data tables, and that is only possible once we have
      82             :   // defined all Pythia-type decays we want to use. 
      83             :   // We check to see if the engine has been created before doing the decay.
      84             :   // This should only create the full Pythia engine once, and all clones will point to the same engine.
      85             : 
      86           0 :   if (_pythiaEngine == 0) {
      87           0 :     _pythiaEngine = EvtExternalGenFactory::getInstance()->getGenerator(EvtExternalGenFactory::PythiaGenId);
      88           0 :   }
      89             :     
      90           0 :   if (_pythiaEngine != 0) {
      91           0 :     _pythiaEngine->doDecay(p);
      92           0 :   }
      93             : 
      94           0 :   this->fixPolarisations(p);
      95             : 
      96           0 : }
      97             : 
      98             : void EvtPythia::fixPolarisations(EvtParticle *p) {
      99             : 
     100             :   // Special case to handle the J/psi polarisation
     101             : 
     102           0 :   if (p == 0) {return;}
     103             : 
     104           0 :   int nDaug = p->getNDaug();  
     105             :   int i(0);
     106             : 
     107           0 :   static EvtId Jpsi = EvtPDL::getId("J/psi");
     108             : 
     109           0 :   for (i = 0; i < nDaug; i++){
     110             : 
     111           0 :     EvtParticle* theDaug = p->getDaug(i);
     112             : 
     113           0 :     if (theDaug != 0) {
     114             : 
     115           0 :       if (theDaug->getId() == Jpsi) {
     116             :   
     117           0 :         EvtSpinDensity rho;
     118             :       
     119           0 :         rho.setDim(3);
     120           0 :         rho.set(0,0,0.5);
     121           0 :         rho.set(0,1,0.0);
     122           0 :         rho.set(0,2,0.0);
     123             : 
     124           0 :         rho.set(1,0,0.0);
     125           0 :         rho.set(1,1,1.0);
     126           0 :         rho.set(1,2,0.0);
     127             : 
     128           0 :         rho.set(2,0,0.0);
     129           0 :         rho.set(2,1,0.0);
     130           0 :         rho.set(2,2,0.5);
     131             : 
     132           0 :         EvtVector4R p4Psi = theDaug->getP4();
     133             : 
     134           0 :         double alpha = atan2(p4Psi.get(2),p4Psi.get(1));
     135           0 :         double beta = acos(p4Psi.get(3)/p4Psi.d3mag());
     136             : 
     137           0 :         theDaug->setSpinDensityForwardHelicityBasis(rho,alpha,beta,0.0);
     138           0 :         setDaughterSpinDensity(i);
     139             : 
     140           0 :       }
     141             :     }
     142             :   }
     143           0 : }
     144             : 
     145             : std::string EvtPythia::commandName() {
     146             : 
     147             :   // Allow backward compatibility for decay.dec files
     148             :   // having JetSetPar parameters. They are obsolete for Pythia 8, 
     149             :   // since the JetSet-type array variables do not exist.
     150             :   // Need to think about including user defined parameters in 
     151             :   // EvtPythiaEngine::updatePhysicsParameters().
     152           0 :   return std::string("JetSetPar");
     153             :   
     154             : }
     155             : 
     156             : void EvtPythia::command(std::string cmd) {
     157             : 
     158             :   // Locally store commands in a vector
     159           0 :   _commandList.push_back(cmd);
     160             : 
     161           0 : }