Ekhara Class Reference

#include <Ekhara.h>

Inheritance diagram for Ekhara:

Public Member Functions
	Ekhara (const std::string &name, ISvcLocator *pSvcLocator)
StatusCode	initialize ()
StatusCode	execute ()
StatusCode	finalize ()
StatusCode	checkAndPrintParameters ()

Detailed Description

Definition at line 15 of file Ekhara.h.

Constructor & Destructor Documentation

Ekhara()

Ekhara::Ekhara	(	const std::string &	name,
		ISvcLocator *	pSvcLocator )

• CMS-energy (GeV)

• # of events to determine the max xs

• handle UB underestim <0> = halt,
• <positive> = use new UB if overshoot,
• <negative> = ignore

• Enlarge factor for upper bound

• 0:pi0pi0, 1:pi+pi-, 2:pi0, 3:eta, 4:eta', 5:Chi_c0, 6:Chi_c1, 7:Chi_c2

• pi0,eta,eta': s(1); t(2); s+t(3); t:signal

• include vacuum polarization (true or false)

• NLO radiative corrections (true or false)

• return weighted events as output for NLO (true or false)

• phase space regulator eps_ph: mgamma=me*eps_ph

• pi0: WZWconst(1); rho pole(2);
  • LMD(3); LMD+V(4); LMD+Vnew(5);
  • 1-Octet(7); 2-Octet(8);
  • eta/eta': WZWconst(1); VMD(3);
  • 1-Octet(7); 2-Octet(8);

• e+ minimum angle

• e+ maximum angle

• e+ minimum Energy

• e+ maximum Energy

• e- minimum angle

• e- maximum angle

• e- minimum Energy

• e- maximum Energy

• pi0, eta, eta' minimum angle

• pi0, eta, eta' maximum angle

• pi0, eta, eta' minimum Energy

• pi0, eta, eta' maximum Energy

• the choice of generation algorithm: 0-v.1.0; 1-v.3.1

• pi+pi-: s(1); s+t(2); s+t+2g_Born(3); s+t+2g_Full(4)

• the choice of the pion form factor to be used: 0-KS; 1-GS; 2-GS new

• pion angle cut

• missing momentum angle cut

• Chi_cj minimum angle

• Chi_cj maximum angle

• Chi_cj minimum energy

• Chi_cj minimum energy

• Tagging angle (symmetric) in degrees

• minimum Q^2 of tagged lepton in GeV^2

• Tagging mode 0:no tagging; 1/11:untagged
• 2/12: single tagging; 3/13:double tagging

Definition at line 26 of file Ekhara.cxx.

    : Algorithm( name, pSvcLocator ) {
  declareProperty( "Ecm", m_cmsEnergy = 3.773 ); //!- CMS-energy (GeV)
  declareProperty( "InitialEvents",
                   m_numUpperBound = 50000 ); //!- # of events to determine the max xs
  declareProperty( "IgnoreUpperBound",
                   m_ignoreUpperBound = 0.0 ); //!- handle UB underestim <0> = halt,
                                               //!- <positive> = use new UB if overshoot,
                                               //!- <negative> = ignore
  declareProperty( "UpperBoundEnlarge",
                   m_upperBoundEnlarge = 1.2 );        //!- Enlarge factor for upper bound
  declareProperty( "FinalStateID", m_finalState = 2 ); //!- 0:pi0pi0, 1:pi+pi-, 2:pi0, 3:eta,
                                                       //! 4:eta', 5:Chi_c0, 6:Chi_c1, 7:Chi_c2
 
  declareProperty( "MesonProdAmplitudes",
                   m_mesonAmplitudes = 2 ); //!- pi0,eta,eta': s(1); t(2); s+t(3); t:signal
  declareProperty( "IncludeVacuumPolarization",
                   m_switchVP = false ); //!- include vacuum polarization (true or false)
  declareProperty( "ConsiderNLO",
                   m_switchNLO = false ); //!- NLO radiative corrections (true or false)
  declareProperty( "OutputWeightedNLOEvents",
                   m_nloWithWeights =
                       false ); //!- return weighted events as output for NLO (true or false)
  declareProperty( "PhaseSpaceRegulator",
                   m_eps_ph = 0.01 ); //!- phase space regulator eps_ph: mgamma=me*eps_ph
  declareProperty( "MesonFormfactor", m_TFF_ID = 8 ); //!- pi0: WZWconst(1); rho pole(2);
                                                      //!-      LMD(3); LMD+V(4); LMD+Vnew(5);
                                                      //!-      1-Octet(7); 2-Octet(8);
                                                      //!- eta/eta': WZWconst(1); VMD(3);
                                                      //!-           1-Octet(7); 2-Octet(8);
  declareProperty( "PositronThetaMin", m_posiThetaMin = 0.0 );     //!- e+ minimum angle
  declareProperty( "PositronThetaMax", m_posiThetaMax = 180.0 );   //!- e+ maximum angle
  declareProperty( "PositronEnergyMin", m_posiEnergyMin = 0.0 );   //!- e+ minimum Energy
  declareProperty( "PositronEnergyMax", m_posiEnergyMax = 110.0 ); //!- e+ maximum Energy
  declareProperty( "ElectronThetaMin", m_elecThetaMin = 0.0 );     //!- e- minimum angle
  declareProperty( "ElectronThetaMax", m_elecThetaMax = 180.0 );   //!- e- maximum angle
  declareProperty( "ElectronEnergyMin", m_elecEnergyMin = 0.0 );   //!- e- minimum Energy
  declareProperty( "ElectronEnergyMax", m_elecEnergyMax = 110.0 ); //!- e- maximum Energy
 
  declareProperty( "MesonThetaMin", m_mesonThetaMin = 0.0 ); //!- pi0, eta, eta' minimum angle
  declareProperty( "MesonThetaMax",
                   m_mesonThetaMax = 180.0 ); //!- pi0, eta, eta' maximum angle
  declareProperty( "MesonEnergyMin",
                   m_mesonEnergyMin = 0.0 ); //!- pi0, eta, eta' minimum Energy
  declareProperty( "MesonEnergyMax",
                   m_mesonEnergyMax = 100.0 ); //!- pi0, eta, eta' maximum Energy
 
  declareProperty( "TwoPionPhaseSpaceAlg",
                   m_twoPionPhspAlg =
                       1 ); //!- the choice of generation algorithm: 0-v.1.0; 1-v.3.1
  declareProperty( "TwoPionProdAmplitudes",
                   m_twoPionAmplitudes =
                       4 ); //!- pi+pi-:  s(1); s+t(2); s+t+2g_Born(3); s+t+2g_Full(4)
  declareProperty( "TwoPionFormFactor",
                   m_twoPionFormFactor = 2 ); //!- the choice of the pion form factor to be
                                              //! used: 0-KS; 1-GS; 2-GS new
  declareProperty( "TwoPionThetaMin", m_twoPionThetaMin = 0.0 ); //!- pion angle cut
  declareProperty( "TwoPionThetaMax", m_twoPionThetaMax = 180.0 );
  declareProperty( "TwoPionMissThetaMin",
                   m_twoPionMissThetaMin = 0.0 ); //!- missing momentum angle cut
  declareProperty( "TwoPionMissThetaMax", m_twoPionMissThetaMax = 180.0 );
 
  declareProperty( "Chi_cjThetaMin", m_chicjThetaMin = 0.0 );     //!- Chi_cj minimum angle
  declareProperty( "Chi_cjThetaMax", m_chicjThetaMax = 180.0 );   //!- Chi_cj maximum angle
  declareProperty( "Chi_cjEnergyMin", m_chicjEnergyMin = 0.0 );   //!- Chi_cj minimum energy
  declareProperty( "Chi_cjEnergyMax", m_chicjEnergyMax = 100.0 ); //!- Chi_cj minimum energy
 
  declareProperty( "TaggingAngle",
                   m_taggingAngle = 20. ); //!- Tagging angle (symmetric) in degrees
  declareProperty( "TaggingQsquare",
                   m_taggingQsquare = 0.3 ); //!- minimum Q^2 of tagged lepton in GeV^2
  declareProperty( "TaggingMode",
                   m_taggingMode =
                       0 ); //!- Tagging mode 0:no tagging; 1/11:untagged
                            //!-              2/12: single tagging; 3/13:double tagging
}

Member Function Documentation

checkAndPrintParameters()

StatusCode Ekhara::checkAndPrintParameters

(

)

• pi0pi0

• Check switches for NLO

• pi+pi-

• Check switches for NLO

• 1 - s, 2 - s+t, 3 - s+t+g_Born, 4 - s+t+g_Full

• 0 - KS, 1 - GS, 2 - new GS

-Single Pseudoscalar Meson

• Check switches for NLO and VacPol

• 1 - s, 2 - t, 3 - s+t

• 1:Untagged , 2:Single, 3: double

• Chi_c0, Chi_c1, Chi_c2

• Check switches for NLO and VacPol

Definition at line 553 of file Ekhara.cxx.

                                           {
  MsgStream log( msgSvc(), name() );
 
  log << MSG::INFO
      << "=============================================================" << endmsg;
  log << MSG::INFO << "This is EKHARA, Version 3.1" << endmsg;
  switch ( m_finalState )
  {
  case 0:
    log << MSG::INFO << "\tSimulating the process:\te+e- -> e+e- pi0pi0" << endmsg;
    break;
  case 1:
    log << MSG::INFO << "\tSimulating the process:\te+e- -> e+e- pi+pi-" << endmsg;
    break;
  case 2: log << MSG::INFO << "\tSimulating the process:\te+e- -> e+e- pi0" << endmsg; break;
  case 3: log << MSG::INFO << "\tSimulating the process:\te+e- -> e+e- eta" << endmsg; break;
  case 4:
    log << MSG::INFO << "\tSimulating the process:\te+e- -> e+e- etaPRIME" << endmsg;
    break;
  case 5:
    log << MSG::INFO << "\tSimulating the process:\te+e- -> e+e- Chi_c0" << endmsg;
    break;
  case 6:
    log << MSG::INFO << "\tSimulating the process:\te+e- -> e+e- Chi_c1" << endmsg;
    break;
  case 7:
    log << MSG::INFO << "\tSimulating the process:\te+e- -> e+e- Chi_c2" << endmsg;
    break;
  default:
    log << MSG::ERROR << "WRONG channel ID: Process not implemented!" << endmsg;
    return StatusCode::FAILURE;
    break;
  }
 
  if ( m_finalState == 0 )
  { //!- pi0pi0
    if ( m_switchNLO )
    { //!- Check switches for NLO
      log << MSG::ERROR << "NLO not implemented for this Process!" << endmsg;
      return StatusCode::FAILURE;
    }
    if ( m_twoPionAmplitudes != 4 )
    {
      log << MSG::ERROR
          << " Neutral Pion Pairs only with full TwoPionAmplitudes == 4 ! \t You used "
          << m_twoPionAmplitudes << endmsg;
      return StatusCode::FAILURE;
    }
    if ( m_twoPionPhspAlg != 1 )
    {
      log << MSG::ERROR
          << " Neutral Pion Pairs only with full TwoPionPhaseSpaceAlg == 1 ! \t You used "
          << m_twoPionPhspAlg << endmsg;
      return StatusCode::FAILURE;
    }
    log << MSG::INFO << "\tPhase space generation using algorithm: " << m_twoPionPhspAlg
        << endmsg;
    log << MSG::INFO << "\tMatrix element: \t|M_s + M_t + M_2g(full)|^2" << endmsg;
 
    log << MSG::INFO << "The following conditions are applied:" << endmsg;
    log << MSG::INFO << "\tPion Momentum:    " << m_twoPionThetaMin << " < Theta [deg] < "
        << m_twoPionThetaMax << endmsg;
    log << MSG::INFO << "\tMissing Momentum: " << m_twoPionMissThetaMin << " < Theta [deg] < "
        << m_twoPionMissThetaMax << endmsg;
  }
  else if ( m_finalState == 1 )
  { //!- pi+pi-
    if ( m_switchNLO )
    { //!- Check switches for NLO
      log << MSG::ERROR << "NLO not implemented for this Process!" << endmsg;
      return StatusCode::FAILURE;
    }
    if ( m_twoPionPhspAlg < 0 || m_twoPionPhspAlg > 1 )
    {
      log << MSG::ERROR
          << " Wrong choice of phase space generation algorithm: " << m_twoPionPhspAlg
          << endmsg;
      return StatusCode::FAILURE;
    }
 
    log << MSG::INFO << "\tPhase space generation using algorithm: " << m_twoPionPhspAlg
        << endmsg;
    log << MSG::INFO << "\tMatrix element: ";
    switch ( m_twoPionAmplitudes )
    { //!-  1 - s, 2 - s+t, 3 - s+t+g_Born, 4 - s+t+g_Full
    case 1: log << MSG::INFO << "\t|M_s|^2" << endmsg; break;
    case 2: log << MSG::INFO << "\t|M_s + M_t|^2" << endmsg; break;
    case 3: log << MSG::INFO << "\t|M_s + M_t + M_2g(Born)|^2" << endmsg; break;
    case 4: log << MSG::INFO << "\t|M_s + M_t + M_2g(full)|^2" << endmsg; break;
    default: break;
    }
 
    switch ( m_twoPionFormFactor )
    { //!-  0 - KS, 1 - GS, 2 - new GS
    case 0: log << MSG::INFO << "\tUsing KS pion form factor" << endmsg; break;
    case 1: log << MSG::INFO << "\tUsing GS pion form factor" << endmsg; break;
    case 2: log << MSG::INFO << "\tUsing new GS pion form factor" << endmsg; break;
    default:
      log << MSG::WARNING
          << "Undefined pion form factor switch! Using new GS pion form factor instead!"
          << endmsg;
      m_twoPionFormFactor = 2;
      break;
    }
 
    log << MSG::INFO << "The following conditions are applied:" << endmsg;
    log << MSG::INFO << "\tPion Momentum:    " << m_twoPionThetaMin << " < Theta [deg] < "
        << m_twoPionThetaMax << endmsg;
    log << MSG::INFO << "\tMissing Momentum: " << m_twoPionMissThetaMin << " < Theta [deg] < "
        << m_twoPionMissThetaMax << endmsg;
  }
  else if ( m_finalState >= 2 && m_finalState <= 4 )
  { //!-Single Pseudoscalar Meson
    if ( m_switchNLO )
    { //!- Check switches for NLO and VacPol
      log << MSG::INFO << "\tWith NLO corrections: mgamma = " << m_eps_ph * 0.51099906E-3
          << endmsg;
      if ( m_nloWithWeights )
      { log << MSG::WARNING << "Output will contain weighted events!" << endmsg; }
      else
      { log << MSG::WARNING << "Attention! Events will be unweighted for output!" << endmsg; }
      if ( m_switchVP )
      {
        log << MSG::INFO << "\tVacuum polarization included:" << endmsg;
        log << MSG::INFO << "\thttp://www-com.physik.hu-berlin.de/~fjeger/software.html"
            << endmsg;
      }
      else { log << MSG::INFO << "\tVacuum polarization NOT included:" << endmsg; }
    }
    else
    {
      log << MSG::INFO << "\tWithout NLO corrections!" << endmsg;
      if ( m_switchVP )
      {
        log << MSG::ERROR << "Vacuum polarization can ONLY be included for NLO!" << endmsg;
        return StatusCode::FAILURE;
      }
    }
 
    log << MSG::INFO << "The cross section will be calculated" << endmsg;
    log << MSG::INFO << "\tMatrix element: ";
    switch ( m_mesonAmplitudes )
    { //!-  1 - s, 2 - t, 3 - s+t
    case 1:
      log << MSG::INFO << "\t|M_s|^2" << endmsg;
      if ( m_switchNLO == 1 )
      {
        log << MSG::ERROR << " NLO not implemented for s-channel!" << endmsg;
        return StatusCode::FAILURE;
      }
      if ( m_TFF_ID == 7 || m_TFF_ID == 8 )
      {
        log << MSG::WARNING
            << "TFF Models 7 and 8 were not compared to data in the only-s-channel "
               "configuration!"
            << endmsg;
      }
      break;
    case 2: log << MSG::INFO << "\t|M_t|^2" << endmsg; break;
    case 3: log << MSG::INFO << "\t|M_s + M_t|^2" << endmsg; break;
    default: break;
    }
 
    log << MSG::INFO << "The following conditions are applied:" << endmsg;
    log << MSG::INFO << "\te+   : " << m_posiThetaMin << " < Theta [deg] < " << m_posiThetaMax
        << endmsg;
    log << MSG::INFO << "\t       " << m_posiEnergyMin << " < Energy [GeV] < "
        << m_posiEnergyMin << endmsg;
    log << MSG::INFO << "\te-   : " << m_elecThetaMin << " < Theta [deg] < " << m_elecThetaMax
        << endmsg;
    log << MSG::INFO << "\t       " << m_elecEnergyMin << " < Energy [GeV] < "
        << m_elecEnergyMin << endmsg;
 
    double cosTagAngleRad = fabs( cos( m_taggingAngle * TMath::DegToRad() ) );
    switch ( m_taggingMode )
    { //!-  1:Untagged , 2:Single, 3: double
    case 1:
      log << MSG::INFO << "\t Generating Untagged event configuration!" << endmsg;
      log << MSG::INFO << "\t Accepting only events with: |cos(Theta Lepton)| > "
          << cosTagAngleRad << "!" << endmsg;
      break;
    case 2:
      log << MSG::INFO << "\t Generating Single Tagged event configuration!" << endmsg;
      log << MSG::INFO << "\t Accepting only events with:" << endmsg;
      log << MSG::INFO << "\t\t |cos(Theta e+/-)| > " << cosTagAngleRad
          << " and |cos(Theta e-/+)| < " << cosTagAngleRad << endmsg;
      break;
    case 3:
      log << MSG::INFO << "\t Generating Double Tagged event configuration!" << endmsg;
      log << MSG::INFO << "\t Accepting only events with: |cos(Theta Lepton)| < "
          << cosTagAngleRad << "!" << endmsg;
      break;
    case 11:
      log << MSG::INFO << "\t Generating Untagged event configuration!" << endmsg;
      log << MSG::INFO << "\t Accepting only events with: Q^2 (Lepton)| < " << m_taggingQsquare
          << "!" << endmsg;
      break;
    case 12:
      log << MSG::INFO << "\t Generating Single Tagged event configuration!" << endmsg;
      log << MSG::INFO << "\t Accepting only events with:" << endmsg;
      log << MSG::INFO << "\t\t Q^2 (e+/-) > " << m_taggingQsquare << " and Q^2 (e-/+) < "
          << m_taggingQsquare << endmsg;
      break;
    case 13:
      log << MSG::INFO << "\t Generating Double Tagged event configuration!" << endmsg;
      log << MSG::INFO << "\t Accepting only events with: Q^2 (Lepton) > " << m_taggingQsquare
          << "!" << endmsg;
      break;
    default:
      if ( m_taggingMode > 0 )
      {
        log << MSG::INFO << "\t Unknown tagging mode selected!" << endmsg;
        log << MSG::INFO << "\t Generating events without tagging configuration!" << endmsg;
      }
      m_taggingMode = 0;
      break;
    }
 
    switch ( m_finalState )
    {
    case 2: log << MSG::INFO << "\tpi0  : "; break;
    case 3: log << MSG::INFO << "\teta  : "; break;
    case 4: log << MSG::INFO << "\teta' : "; break;
    default: break;
    }
    log << MSG::INFO << m_mesonThetaMin << " < Theta [deg] < " << m_mesonThetaMax << endmsg;
    log << MSG::INFO << "\t       " << m_mesonEnergyMin << " < Energy [GeV] < "
        << m_mesonEnergyMin << endmsg;
  }
  else if ( m_finalState >= 5 && m_finalState <= 7 )
  { //!- Chi_c0, Chi_c1, Chi_c2
    if ( m_switchNLO )
    { //!- Check switches for NLO and VacPol
      log << MSG::ERROR << "NLO not implemented for this Process!" << endmsg;
      return StatusCode::FAILURE;
    }
    log << MSG::INFO << "The following conditions are applied:" << endmsg;
    log << MSG::INFO << "\te+   : " << m_posiThetaMin << " < Theta [deg] < " << m_posiThetaMax
        << endmsg;
    log << MSG::INFO << "\t       " << m_posiEnergyMin << " < Energy [GeV] < "
        << m_posiEnergyMin << endmsg;
    log << MSG::INFO << "\te-   : " << m_elecThetaMin << " < Theta [deg] < " << m_elecThetaMax
        << endmsg;
    log << MSG::INFO << "\t       " << m_elecEnergyMin << " < Energy [GeV] < "
        << m_elecEnergyMin << endmsg;
    switch ( m_finalState )
    {
    case 5: log << MSG::INFO << "\tchi_c0  : "; break;
    case 6: log << MSG::INFO << "\tchi_c1  : "; break;
    case 7: log << MSG::INFO << "\tchi_c2 : "; break;
    default: break;
    }
    log << MSG::INFO << m_chicjThetaMin << " < Theta [deg] < " << m_chicjThetaMax << endmsg;
    log << MSG::INFO << "\t       " << m_chicjEnergyMin << " < Energy [GeV] < "
        << m_chicjEnergyMin << endmsg;
  }
  log << MSG::INFO
      << "=============================================================" << endmsg;
 
  return StatusCode::SUCCESS;
}

Referenced by initialize().

execute()

StatusCode Ekhara::execute

(

)

• Generate an event

-Store weigth as particle at rest:

Definition at line 179 of file Ekhara.cxx.

                           {
  MsgStream log( msgSvc(), name() );
  log << MSG::DEBUG << "Ekhara in execute()" << endmsg;
 
  //!- Generate an event
  if ( m_finalState >= 2 && m_finalState <= 4 && m_switchNLO && m_nloWithWeights )
  { EKHARA( 2 ); }
  else { EKHARA( 0 ); }
 
  double p1[4], p2[4], q1[4], q2[4];
  double pi1[4], pi2[4], qpion[4], qcj[4], kphp[4];
 
  for ( int i = 0; i < 4; i++ )
  {
    p1[i]    = 0.0; // four vector of incoming positron
    p2[i]    = 0.0; // four vector of incoming electron
    q1[i]    = 0.0; // four vector of outgoung positron
    q2[i]    = 0.0; // four vector of outgoing electron
    pi1[i]   = 0.0; // four vector of produced pi0,pi+ (m_finalState == 0,1)
    pi2[i]   = 0.0; // four vector of produced pi0,pi- (m_finalState == 0,1)
    qpion[i] = 0.0; // four vector of produced pi0,eta,eta' (m_finalState == 2,3,4)
    qcj[i]   = 0.0; // four vector of produced chi_cj (m_finalState == 5,6,7)
    kphp[i]  = 0.0; // four vector of radiative photon
  }
 
  // Fill event information
  GenEvent* evt = new GenEvent( 1, 1 );
 
  GenVertex* prod_vtx = new GenVertex();
  evt->add_vertex( prod_vtx );
 
  GET_FOURMOMENTA_LEPTONS( p1, p2, q1, q2 );
 
  // incoming beam e+
  GenParticle* posi =
      new GenParticle( CLHEP::HepLorentzVector( p1[1], p1[2], p1[3], p1[0] ), -11, 3 );
  posi->suggest_barcode( 1 );
  prod_vtx->add_particle_in( posi );
 
  // incoming beam e-
  GenParticle* elec =
      new GenParticle( CLHEP::HepLorentzVector( p2[1], p2[2], p2[3], p2[0] ), 11, 3 );
  elec->suggest_barcode( 2 );
  prod_vtx->add_particle_in( elec );
 
  evt->set_beam_particles( posi, elec );
 
  // outcoming beam e+
  GenParticle* p =
      new GenParticle( CLHEP::HepLorentzVector( q1[1], q1[2], q1[3], q1[0] ), -11, 1 );
  p->suggest_barcode( 3 );
  prod_vtx->add_particle_out( p );
 
  // outcoming beam e-
  p = new GenParticle( CLHEP::HepLorentzVector( q2[1], q2[2], q2[3], q2[0] ), 11, 1 );
  p->suggest_barcode( 4 );
  prod_vtx->add_particle_out( p );
 
  switch ( m_finalState )
  {
  case 0: // e+e-pi0 pi0
    GET_FOURMOMENTA_TWOPI( pi1, pi2 );
    p = new GenParticle( CLHEP::HepLorentzVector( pi1[1], pi1[2], pi1[3], pi1[0] ), 111, 1 );
    p->suggest_barcode( 5 );
    prod_vtx->add_particle_out( p );
    p = new GenParticle( CLHEP::HepLorentzVector( pi2[1], pi2[2], pi2[3], pi2[0] ), 111, 1 );
    p->suggest_barcode( 6 );
    prod_vtx->add_particle_out( p );
    break;
  case 1: // e+e-pi+ pi-
    GET_FOURMOMENTA_TWOPI( pi1, pi2 );
    p = new GenParticle( CLHEP::HepLorentzVector( pi1[1], pi1[2], pi1[3], pi1[0] ), 211, 1 );
    p->suggest_barcode( 5 );
    prod_vtx->add_particle_out( p );
    p = new GenParticle( CLHEP::HepLorentzVector( pi2[1], pi2[2], pi2[3], pi2[0] ), -211, 1 );
    p->suggest_barcode( 6 );
    prod_vtx->add_particle_out( p );
    double weights[6]; // event weight + five weights for the various contributions to e+e- -->
                       // e+e-pi+pi-
    for ( int i = 0; i < 6; i++ ) { weights[i] = 0.0; }
    GET_TWOPI_WEIGHTS( weights );
    for ( int ai = 0; ai < ( sizeof( weights ) / sizeof( double ) ); ai++ )
    { evt->weights().push_back( weights[ai] ); }
    p = new GenParticle( CLHEP::HepLorentzVector( 0, 0, 0, weights[0] ), 99, 1 );
    p->suggest_barcode( 7 );
    prod_vtx->add_particle_out( p );
 
    break;
  case 2: // e+e- pi0
    GET_FOURMOMENTA_PION( qpion );
    p = new GenParticle( CLHEP::HepLorentzVector( qpion[1], qpion[2], qpion[3], qpion[0] ),
                         111, 1 );
    p->suggest_barcode( 5 );
    prod_vtx->add_particle_out( p );
    break;
  case 3: // e+e- eta
    GET_FOURMOMENTA_PION( qpion );
    p = new GenParticle( CLHEP::HepLorentzVector( qpion[1], qpion[2], qpion[3], qpion[0] ),
                         221, 1 );
    p->suggest_barcode( 5 );
    prod_vtx->add_particle_out( p );
    break;
  case 4: // e+e- etaPrime
    GET_FOURMOMENTA_PION( qpion );
    p = new GenParticle( CLHEP::HepLorentzVector( qpion[1], qpion[2], qpion[3], qpion[0] ),
                         331, 1 );
    p->suggest_barcode( 5 );
    prod_vtx->add_particle_out( p );
    break;
  case 5: // e+e- chi_c0
    GET_FOURMOMENTA_CHICJ( qcj );
    p = new GenParticle( CLHEP::HepLorentzVector( qcj[1], qcj[2], qcj[3], qcj[0] ), 10441, 1 );
    p->suggest_barcode( 5 );
    prod_vtx->add_particle_out( p );
    break;
  case 6: // e+e- chi_c1
    GET_FOURMOMENTA_CHICJ( qcj );
    p = new GenParticle( CLHEP::HepLorentzVector( qcj[1], qcj[2], qcj[3], qcj[0] ), 20443, 1 );
    p->suggest_barcode( 5 );
    prod_vtx->add_particle_out( p );
    break;
  case 7: // e+e- chi_c2
    GET_FOURMOMENTA_CHICJ( qcj );
    p = new GenParticle( CLHEP::HepLorentzVector( qcj[1], qcj[2], qcj[3], qcj[0] ), 445, 1 );
    p->suggest_barcode( 5 );
    prod_vtx->add_particle_out( p );
    break;
  default: break;
  }
 
  if ( m_finalState >= 2 && m_finalState <= 4 && m_switchNLO )
  {
    if ( NLOTYPE.EvtPhTyp == 1 )
    {
      GET_FOURMOMENTA_PHOTON( kphp );
      p = new GenParticle( CLHEP::HepLorentzVector( kphp[1], kphp[2], kphp[3], kphp[0] ), 22,
                           1 );
      p->suggest_barcode( 6 );
      prod_vtx->add_particle_out( p );
    }
    if ( m_nloWithWeights )
    {
      double weight = GET_WEIGHT();
      evt->weights().push_back( weight );
      //!-Store weigth as particle at rest:
      if ( NLOTYPE.EvtPhTyp == 1 )
      {
        p = new GenParticle( CLHEP::HepLorentzVector( 0, 0, 0, weight ), 99, 1 );
        p->suggest_barcode( 7 );
        prod_vtx->add_particle_out( p );
      }
      else if ( NLOTYPE.EvtPhTyp == 0 )
      {
        p = new GenParticle( CLHEP::HepLorentzVector( 0, 0, 0, weight ), 98, 1 );
        p->suggest_barcode( 6 );
        prod_vtx->add_particle_out( p );
      }
    }
  }
 
  if ( log.level() < MSG::INFO ) { evt->print(); }
 
  // Check if the McCollection already exists
  SmartDataPtr<McGenEventCol> anMcCol( eventSvc(), "/Event/Gen" );
  if ( anMcCol != 0 )
  {
    log << MSG::DEBUG << "Adding McGenEvent to existing collection" << endmsg;
    McGenEvent* mcEvent = new McGenEvent( evt );
    anMcCol->push_back( mcEvent );
  }
  else
  {
    McGenEventCol* mcColl  = new McGenEventCol;
    McGenEvent*    mcEvent = new McGenEvent( evt );
    mcColl->push_back( mcEvent );
    StatusCode sc = eventSvc()->registerObject( "/Event/Gen", mcColl );
    if ( sc != StatusCode::SUCCESS )
    {
      log << MSG::ERROR << "Could not register McGenEvent" << endmsg;
      delete mcColl;
      delete evt;
      delete mcEvent;
      return StatusCode::FAILURE;
    }
    else { log << MSG::INFO << "McGenEventCol created " << endmsg; }
  }
 
  m_totalEvents++;
 
  log << MSG::DEBUG << "before execute() return" << endmsg;
  return StatusCode::SUCCESS;
}

finalize()

StatusCode Ekhara::finalize

(

)

e+e- --> e+e-pi0pi0

e+e- --> e+e-pi+pi-

Definition at line 372 of file Ekhara.cxx.

                            {
  MsgStream log( msgSvc(), name() );
  log << MSG::INFO << "Ekhara in finalize()" << endmsg;
 
  EKHARA( 1 );
 
  char xsect[100];
  if ( m_finalState == 0 )
  { //! e+e- --> e+e-pi0pi0
    double tnpfinpar[10];
    for ( int i = 0; i < 10; i++ ) { tnpfinpar[i] = 0.0; }
 
    GET_FINAL_MESON_INFO( 0, tnpfinpar );
 
    log << MSG::INFO << "Total number of iterations = " << tnpfinpar[0] << endmsg;
    log << MSG::INFO << "Number of MC trials = " << tnpfinpar[1] << endmsg;
    log << MSG::INFO << "ContribMax_MCloop  = " << tnpfinpar[2] << endmsg;
    log << MSG::INFO << "MC efficiency (evt/eff.iter) = " << m_totalEvents / tnpfinpar[1]
        << endmsg;
 
    log << MSG::INFO << "================================" << endmsg;
    log << MSG::INFO << "From all weights:" << endmsg;
    double CS = tnpfinpar[3] / tnpfinpar[0];
    double CSerr =
        sqrt( fabs( tnpfinpar[4] / tnpfinpar[0] - CS * CS ) / ( tnpfinpar[0] - 1. ) );
    log << MSG::INFO << " Integral Cross Section = " << CS << " +- " << CSerr << " [nb]"
        << endmsg;
    log << MSG::INFO << endmsg;
    log << MSG::INFO << "From UNWEIGHTED events:" << endmsg;
    if ( tnpfinpar[1] <= 0. )
    {
      CS    = 0.;
      CSerr = 0.;
    }
    else
    {
      CS    = tnpfinpar[5] * m_totalEvents / tnpfinpar[1];
      CSerr = sqrt( CS * ( tnpfinpar[5] - CS ) ) / sqrt( tnpfinpar[1] );
    }
    log << MSG::INFO << " Integral Cross Section = " << CS << " +- " << CSerr << " [nb]"
        << endmsg;
    log << MSG::INFO << "================================" << endmsg;
    if ( tnpfinpar[7] > 0. )
    {
      log << MSG::WARNING << "There were overshootings." << endmsg;
      log << MSG::WARNING << "The events, corresponding to them were ignored." << endmsg;
      log << MSG::WARNING << "The unweighted sample may be not reliable!" << endmsg;
      log << MSG::WARNING << "\nContribution from overshooted events:" << endmsg;
      sprintf( xsect, " Integral Cross Section = %1.5E +- %1.5E  [GeV^-2]", tnpfinpar[8],
               tnpfinpar[9] );
      log << MSG::WARNING << xsect << endmsg;
      sprintf( xsect, "                        = %1.5E +- %1.5E  [nb]",
               tnpfinpar[8] * tnpfinpar[6], tnpfinpar[9] * tnpfinpar[6] );
      log << MSG::WARNING << xsect << endmsg;
      log << MSG::INFO << "================================" << endmsg;
    }
  }
  else if ( m_finalState == 1 )
  { //! e+e- --> e+e-pi+pi-
    double pipifinpar[2];
    for ( int i = 0; i < 2; i++ ) { pipifinpar[i] = 0.0; }
    GET_FINAL_TWOPI_INFO( pipifinpar );
 
    double CS = pipifinpar[0] / m_totalEvents;
    double CSerr =
        sqrt( ( pipifinpar[1] / m_totalEvents - CS * CS ) / ( m_totalEvents - 1. ) );
    log << MSG::INFO << "================================" << endmsg;
    log << MSG::INFO << "Cross Section = " << CS << " +- " << CSerr << " [nb]" << endmsg;
    log << MSG::INFO << "================================" << endmsg;
  }
  else if ( m_finalState > 1 && m_finalState < 5 )
  {
    double mfp[10];
    for ( int i = 0; i < 10; i++ ) { mfp[i] = 0.0; }
    GET_FINAL_MESON_INFO( 0, mfp );
    if ( m_switchNLO )
    {
      log << MSG::INFO << "================================" << endmsg;
      log << MSG::INFO << "Events 0ph:" << endmsg;
      log << MSG::INFO << "-------------------------------- " << endmsg;
    }
    log << MSG::INFO << "Total number of iterations = " << mfp[0] << endmsg;
    if ( !m_switchNLO )
    {
      log << MSG::INFO << "Number of MC trials = " << mfp[1] << endmsg;
      log << MSG::INFO << "ContribMax_MCloop  = " << mfp[2] << endmsg;
      log << MSG::INFO << "MC efficiency (evt/eff.iter) = " << m_totalEvents / mfp[1]
          << endmsg;
    }
    log << MSG::INFO << "================================" << endmsg;
    log << MSG::INFO << "From all weights:" << endmsg;
    double CS       = mfp[3] / mfp[0];
    double CSerr    = sqrt( fabs( mfp[4] / mfp[0] - CS * CS ) / ( mfp[0] - 1. ) );
    double sumCS    = CS;
    double sumCSerr = CSerr * CSerr;
    log << MSG::INFO << " Integral Cross Section = " << CS << " +- " << CSerr << " [GeV^-2]"
        << endmsg;
    log << MSG::INFO << "                        = " << CS * mfp[6] << " +- " << CSerr * mfp[6]
        << " [nb]" << endmsg;
    log << MSG::INFO << endmsg;
    if ( m_switchNLO )
    {
      for ( int i = 0; i < 10; i++ ) { mfp[i] = 0.0; }
      GET_FINAL_MESON_INFO( 1, mfp );
      log << MSG::INFO << "================================" << endmsg;
      log << MSG::INFO << "Events 1ph:" << endmsg;
      log << MSG::INFO << "-------------------------------- " << endmsg;
      log << MSG::INFO << "Total number of iterations = " << mfp[0] << endmsg;
      log << MSG::INFO << "================================" << endmsg;
      log << MSG::INFO << "From all weights:" << endmsg;
      CS    = mfp[3] / mfp[0];
      CSerr = sqrt( fabs( mfp[4] / mfp[0] - CS * CS ) / ( mfp[0] - 1. ) );
      sumCS += CS;
      sumCSerr += CSerr * CSerr;
      log << MSG::INFO << " Integral Cross Section = " << CS << " +- " << CSerr << " [GeV^-2]"
          << endmsg;
      log << MSG::INFO << "                        = " << CS * mfp[6] << " +- "
          << CSerr * mfp[6] << " [nb]" << endmsg;
      log << MSG::INFO << endmsg;
      log << MSG::INFO << " -------------------------------- " << endmsg;
      log << MSG::INFO << " sigma (0ph+1ph events) = " << sumCS << " +- " << sqrt( sumCSerr )
          << " [GeV^-2]" << endmsg;
      log << MSG::INFO << "                        = " << sumCS * mfp[6] << " +- "
          << sqrt( sumCSerr ) * mfp[6] << " [nb]" << endmsg;
      log << MSG::INFO << "================================" << endmsg;
    }
    else
    {
      log << MSG::INFO << "From UNWEIGHTED events:" << endmsg;
      if ( mfp[1] <= 0. )
      {
        CS    = 0.;
        CSerr = 0.;
      }
      else
      {
        CS    = mfp[5] * m_totalEvents / mfp[1];
        CSerr = sqrt( CS * ( mfp[5] - CS ) ) / sqrt( mfp[1] );
      }
      sprintf( xsect, " Integral Cross Section = %1.5E +- %1.5E  [GeV^-2]", CS, CSerr );
      log << MSG::INFO << xsect << endmsg;
      sprintf( xsect, "                        = %1.5E +- %1.5E  [nb]", CS * mfp[6],
               CSerr * mfp[6] );
      log << MSG::INFO << xsect << endmsg;
      log << MSG::INFO << "================================" << endmsg;
 
      if ( mfp[7] > 0. )
      {
        log << MSG::WARNING << "There were overshootings." << endmsg;
        log << MSG::WARNING << "The events, corresponding to them were ignored." << endmsg;
        log << MSG::WARNING << "The unweighted sample may be not reliable!" << endmsg;
        log << MSG::WARNING << "\nContribution from overshooted events:" << endmsg;
        sprintf( xsect, " Integral Cross Section = %1.5E +- %1.5E  [GeV^-2]", mfp[8], mfp[9] );
        log << MSG::WARNING << xsect << endmsg;
        sprintf( xsect, "                        = %1.5E +- %1.5E  [nb]", mfp[8] * mfp[6],
                 mfp[9] * mfp[6] );
        log << MSG::WARNING << xsect << endmsg;
        log << MSG::INFO << "================================" << endmsg;
      }
    }
  }
  else if ( m_finalState > 4 && m_finalState < 8 )
  {
    double chicjfinpar[3];
    for ( int i = 0; i < 3; i++ ) { chicjfinpar[i] = 0.0; }
    GET_FINAL_CHICJ_INFO( chicjfinpar );
 
    double CS = chicjfinpar[1] / chicjfinpar[0];
    double CSerr =
        sqrt( ( chicjfinpar[2] / chicjfinpar[0] - CS * CS ) / ( chicjfinpar[0] - 1. ) );
    log << MSG::INFO << "================================" << endmsg;
    sprintf( xsect, " Cross Section = %1.5E +- %1.5E  [nb]", CS, CSerr );
    log << MSG::INFO << xsect << endmsg;
    log << MSG::INFO << "================================" << endmsg;
  }
 
  log << MSG::DEBUG << "Ekhara has terminated successfully" << endmsg;
 
  return StatusCode::SUCCESS;
}

initialize()

StatusCode Ekhara::initialize

(

)

• Set Bes unified random engine

• Perform consitency checks of jobOptions and print screen messages

• Transfer integer parameters to FORTRAN common blocks

• Prepare floating point parameters to be transfered to quadruple precision

• Transfer floating point parameters to FORTAN routines
• and initialize FORTAN routines

• Initialize event counter

Definition at line 103 of file Ekhara.cxx.

                              {
  MsgStream log( msgSvc(), name() );
  log << MSG::DEBUG << "Ekhara in initialize()" << endmsg;
 
  //!- Set Bes unified random engine
  static const bool CREATEIFNOTTHERE( true );
  StatusCode        RndmStatus = service( "BesRndmGenSvc", p_BesRndmGenSvc, CREATEIFNOTTHERE );
  if ( !RndmStatus.isSuccess() || 0 == p_BesRndmGenSvc )
  {
    log << MSG::ERROR << " Could not initialize Random Number Service" << endmsg;
    return RndmStatus;
  }
  HepRandomEngine* engine = p_BesRndmGenSvc->GetEngine( "Ekhara" );
  EkharaRandom::setRandomEngine( engine );
 
  //!- Perform consitency checks of jobOptions and print screen messages
  if ( !checkAndPrintParameters() ) return StatusCode::FAILURE;
 
  //!- Transfer integer parameters to FORTRAN common blocks
  CHANNELSEL.channel_id = m_finalState;
  CHANNELSEL.NLO1P      = 0;
  if ( m_switchNLO ) CHANNELSEL.NLO1P = 1;
  CHANNELSEL.VPSW = 0;
  if ( m_switchVP ) CHANNELSEL.VPSW = 1;
  CHANNELSEL.phsp_2pi = m_twoPionPhspAlg;
  SWDIAG.sw_2pi       = m_twoPionAmplitudes;
  SWDIAG.sw_1pi       = m_mesonAmplitudes;
  PIONFFSW.piggFFsw   = m_TFF_ID;
  TAGGINGMODE.tagmode = m_taggingMode;
  FFPARAMSET.FF_pion  = m_twoPionFormFactor;
  NLOTYPE.weighted    = 0;
  if ( m_nloWithWeights ) NLOTYPE.weighted = 1;
 
  //!- Prepare floating point parameters to be transfered to quadruple precision
  double xpar[27];
  xpar[0]  = m_cmsEnergy;
  xpar[1]  = m_numUpperBound;
  xpar[2]  = m_ignoreUpperBound;
  xpar[3]  = m_upperBoundEnlarge;
  xpar[4]  = m_eps_ph;
  xpar[5]  = m_posiThetaMin;
  xpar[6]  = m_posiThetaMax;
  xpar[7]  = m_posiEnergyMin;
  xpar[8]  = m_posiEnergyMax;
  xpar[9]  = m_elecThetaMin;
  xpar[10] = m_elecThetaMax;
  xpar[11] = m_elecEnergyMin;
  xpar[12] = m_elecEnergyMax;
  xpar[13] = m_mesonThetaMin;
  xpar[14] = m_mesonThetaMax;
  xpar[15] = m_mesonEnergyMin;
  xpar[16] = m_mesonEnergyMax;
  xpar[17] = m_twoPionThetaMin;
  xpar[18] = m_twoPionThetaMax;
  xpar[19] = m_twoPionMissThetaMin;
  xpar[20] = m_twoPionMissThetaMax;
  xpar[21] = m_chicjThetaMin;
  xpar[22] = m_chicjThetaMax;
  xpar[23] = m_chicjEnergyMin;
  xpar[24] = m_chicjEnergyMax;
  xpar[25] = m_taggingAngle;
  xpar[26] = m_taggingQsquare;
 
  //!- Transfer floating point parameters to FORTAN routines
  //!- and initialize FORTAN routines
  BOSS_INIT_EKHARA( xpar );
 
  if ( log.level() < MSG::INFO ) { DIAGNOSE(); }
 
  //!- Initialize event counter
  // should be replaced by request to ApplicationManager for EvtMax
  m_totalEvents = 0;
 
  return StatusCode::SUCCESS;
}

---

The documentation for this class was generated from the following files:

• Ekhara.h
• Ekhara.cxx