Main Page | Namespace List | Class Hierarchy | Alphabetical List | Class List | Directories | File List | Namespace Members | Class Members | File Members | Related Pages

Bhlumi Class Reference

#include <Bhlumi.h>

List of all members.

Public Member Functions

 Bhlumi (const string &name, ISvcLocator *pSvcLocator)
 Bhlumi (const string &name, ISvcLocator *pSvcLocator)
StatusCode execute ()
StatusCode execute ()
StatusCode finalize ()
StatusCode finalize ()
StatusCode initialize ()
StatusCode initialize ()

Private Attributes

int m_angleMode
double m_cmEnergy
double m_infraredCut
std::vector< int > m_initSeed
std::vector< int > m_initSeed
double m_maxThetaAngle
double m_minThetaAngle
int npar [100]
double xpar [100]

Constructor & Destructor Documentation

Bhlumi::Bhlumi (  const string &  name,
ISvcLocator *  pSvcLocator
) 
  

00061                                                           :Algorithm( name, pSvcLocator )
00062 {
00063   declareProperty("CMEnergy", m_cmEnergy = 3.097); // 2*Ebeam [GeV]
00064   declareProperty("AngleMode", m_angleMode = 0);   //0: rad; 1: degree; 2: cos(theta);
00065                                                    //if 1 or 2, angle will be controlled absolutely by user
00066   declareProperty("MinThetaAngle", m_minThetaAngle = 0.24); // [rad]
00067   declareProperty("MaxThetaAngle", m_maxThetaAngle = 0.58); // [rad]
00068   declareProperty("SoftPhotonCut", m_infraredCut   = 1E-4); // dimensionless, Ephoton > m_infraredCut*sqrt(s)/2
00069   m_initSeed.clear();
00070   m_initSeed.push_back(54217137); m_initSeed.push_back(0); m_initSeed.push_back(0);
00071   declareProperty("InitializedSeed", m_initSeed);
00072 }

Bhlumi::Bhlumi (  const string &  name,
ISvcLocator *  pSvcLocator
) 
 

Member Function Documentation

StatusCode Bhlumi::execute (   ) 
 

StatusCode Bhlumi::execute (   ) 
  

00166 { 
00167   MsgStream log(messageService(), name());
00168   log << MSG::INFO << "Bhlumi executing" << endreq;
00169 
00170 
00171   BHLUMI( 0,xpar,npar);
00172 
00173   if( log.level() < MSG::INFO )
00174   {
00175     DUMPS(6);
00176    // dump output to file
00177    //  DUMPS(16);
00178   }
00179 
00180   int npart = 0;
00181   
00182   // Fill event information
00183   GenEvent* evt = new GenEvent(1,1);
00184 
00185   GenVertex* prod_vtx = new GenVertex();
00186 //  prod_vtx->add_particle_out( p );
00187   evt->add_vertex( prod_vtx );
00188                           
00189   // incoming beam e+
00190   GenParticle* p = new GenParticle( HepLorentzVector( MOMSET.p1[0],  MOMSET.p1[1], 
00191                                                       MOMSET.p1[2],  MOMSET.p1[3]),
00192                                                       -11, 3); 
00193   p->suggest_barcode( ++npart );
00194   prod_vtx->add_particle_in(p);
00195 
00196   // incoming beam e-
00197   p = new GenParticle( HepLorentzVector( MOMSET.q1[0],  MOMSET.q1[1], MOMSET.q1[2],  MOMSET.q1[3]),
00198                        11, 3); 
00199   p->suggest_barcode( ++npart );
00200   prod_vtx->add_particle_in(p);
00201   
00202   // scattered e+
00203   p = new GenParticle( HepLorentzVector( MOMSET.p2[0],  MOMSET.p2[1], 
00204                                          MOMSET.p2[2],  MOMSET.p2[3]),
00205                                          -11, 1); 
00206   p->suggest_barcode( ++npart );
00207   prod_vtx->add_particle_out(p);
00208 
00209   // scattered e-
00210   p = new GenParticle( HepLorentzVector( MOMSET.q2[0],  MOMSET.q2[1], MOMSET.q2[2],  MOMSET.q2[3]),
00211                        11, 1); 
00212   p->suggest_barcode( ++npart );
00213   prod_vtx->add_particle_out(p);
00214 
00215   int iphot=0;
00216   for (iphot=0; iphot<MOMSET.nphot; iphot++)
00217   {
00218     // gamma
00219     p = new GenParticle( HepLorentzVector( MOMSET.phot[0][iphot],  MOMSET.phot[1][iphot], 
00220                                            MOMSET.phot[2][iphot],  MOMSET.phot[3][iphot]),
00221                                            22, 1); 
00222     p->suggest_barcode( ++npart );
00223     prod_vtx->add_particle_out(p);
00224   }
00225   
00226   if( log.level() < MSG::INFO )
00227   {
00228     evt->print();  
00229   }
00230 
00231   // Check if the McCollection already exists
00232   SmartDataPtr<McGenEventCol> anMcCol(eventSvc(), "/Event/Gen");
00233   if (anMcCol!=0) 
00234   {
00235     // Add event to existing collection
00236     MsgStream log(messageService(), name());
00237     log << MSG::INFO << "Add McGenEvent to existing collection" << endreq;
00238     McGenEvent* mcEvent = new McGenEvent(evt);
00239     anMcCol->push_back(mcEvent);
00240   }
00241   else 
00242   {
00243     // Create Collection and add  to the transient store
00244     McGenEventCol *mcColl = new McGenEventCol;
00245     McGenEvent* mcEvent = new McGenEvent(evt);
00246     mcColl->push_back(mcEvent);
00247     StatusCode sc = eventSvc()->registerObject("/Event/Gen",mcColl);
00248     if (sc != StatusCode::SUCCESS) 
00249     {
00250       log << MSG::ERROR << "Could not register McGenEvent" << endreq;
00251       delete mcColl;
00252       delete evt;
00253       delete mcEvent;
00254       return StatusCode::FAILURE;
00255     }
00256     else 
00257     {
00258       log << MSG::INFO << "McGenEventCol created and " << npart <<" particles stored in McGenEvent" << endreq;
00259     }
00260   }
00261 
00262   // retrieve event from Transient Store (Storegate)
00263 /*  SmartDataPtr<McGenEventCol> McEvtColl(eventSvc(), "/Event/Gen");    
00264   if ( McEvtColl == 0 ) 
00265   {
00266     log << MSG::ERROR << "Could not retrieve McGenEventCollection" << endreq;
00267     return StatusCode::FAILURE;
00268   };
00269   
00270   McGenEventCol::iterator mcItr;
00271   for( mcItr = McEvtColl->begin(); mcItr != McEvtColl->end(); mcItr++ )   
00272   {
00273      HepMC::GenEvent* hEvt = (*mcItr)->getGenEvt();
00274     // MeVToGeV( hEvt );
00275 //     callEvtGen( hEvt );
00276     // GeVToMeV( hEvt );
00277   };
00278 */
00279      
00280   return StatusCode::SUCCESS; 
00281 }

StatusCode Bhlumi::finalize (   ) 
 

StatusCode Bhlumi::finalize (   ) 
  

00284 {
00285   MsgStream log(messageService(), name());
00286 
00287   BHLUMI( 2,xpar,npar);
00288 
00289   log << MSG::INFO << "Bhlumi finalized" << endreq;
00290 
00291   return StatusCode::SUCCESS;
00292 }

StatusCode Bhlumi::initialize (   ) 
 

StatusCode Bhlumi::initialize (   ) 
 

||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| ! User should cross-check the folowing two output cross sections ! which are calculated and printed at the very end of the output: ! Workshop95, Table14, BARE1 WW for zmin=0.5: KeyGen=3, KeyPia=0, KeyZet=0 ! Workshop95, Table18, CALO2 WW for zmin=0.5: KeyGen=3, KeyPia=2, KeyZet=1 !|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||

--------------------------------------------------------------------------- Input parameters for Bhlumi ---------------------------------------------------------------------- Technical parameters, nothing should depend on them:

Try both options for KeyWgt, result should be the same

-------------------------------------------------- Physics parameters:

--------------------------------------

2*Ebeam [GeV], as in Workshop95 

00074                              {
00075 
00076   MsgStream log(messageService(), name());
00077   
00078   log << MSG::INFO << "Bhlumi initialize" << endreq;
00079 
00080   GLIMIT(50000);
00081                                                                                             
00091 
00092 
00093 
00094 
00095   int KeyGen =   3;  // ! Multiphoton Bhlumi
00096   int KeyRem =   1;  // ! No remooval of photons below epsCM
00097   KeyRem =   0;  // ! Remooval of photons below epsCM, Necessary for Z!!!
00099   int KeyWgt =   2;  // ! weighted events, with t generation down to zero
00100   KeyWgt =   0;  // ! WT=1 events, for detector simulation
00101   int KeyRnd =   1;  // ! RANMAR random numbers
00102   int KeyOpt =1000*KeyGen +100*KeyRem +10*KeyWgt +KeyRnd;
00105   int KeyPia =   0;  // ! Vacuum polarization OFF
00106   int KeyMod =   2;  //   ! Matrix element default version 4.x
00107   KeyPia =   2;  //   ! Vacuum polarization ON
00108   int  KeyZet =   0;  //   ! Z contribution OFF
00109   KeyZet =   1;  //   ! Z contribution ON
00110   int KeyRad =1000*KeyZet +10*KeyMod +KeyPia;
00112   npar[0]=   KeyOpt;
00113   npar[1]=   KeyRad;
00114   double CmsEne =  m_cmEnergy;  
00115   xpar[0]=   CmsEne;
00116   double th1,th2,thmin,thmax;
00117   if(m_angleMode==0){
00118     th1    =  m_minThetaAngle;  //       ! Detector range ThetaMin [rad]
00119     th2    =  m_maxThetaAngle;  //       ! Detector range ThetaMax [rad]
00120     thmin    = 0.7*th1;  //   ! thmin has to be lower  than th1
00121     thmax    = 2.0*th2;  //   ! thmax has to be higher than th2
00122     if(thmin<0.||thmax>3.1415926) {
00123       log << MSG::WARNING << "input angles exceed range (0~pi), so effect will be unprospectable" << endreq;
00124       return StatusCode::FAILURE;
00125     }
00126   }
00127   else if(m_angleMode==1){
00128     th1 = m_minThetaAngle*3.1415926/180.;
00129     th2 = m_maxThetaAngle*3.1415926/180.;
00130     // not multiply  0.7(2.0) coefficient while large angle
00131     thmin    = th1;
00132     thmax    = th2;
00133   }
00134   else if(m_angleMode==2){
00135     th1 = acos(max(m_minThetaAngle,m_maxThetaAngle));
00136     th2 = acos(min(m_minThetaAngle,m_maxThetaAngle));
00137     // not multiply  0.7(2.0) coefficient while large angle
00138     thmin    = th1;
00139     thmax    = th2;
00140   }
00141   else{
00142     log << MSG::FATAL << "unknown angle mode!" << endreq;
00143     return StatusCode::FAILURE;
00144   }
00145   if(thmin<0.||thmax>3.1415926) {
00146     log << MSG::FATAL << "input angles exceed range (0~pi), unprospectable" << endreq;
00147     return StatusCode::FAILURE;
00148   }
00149   else if(thmin>thmax) {
00150     log << MSG::FATAL << "thmin>thmax, unprospectable" << endreq;
00151     return StatusCode::FAILURE;
00152   }
00153   if(KeyWgt == 2) thmin=th1;  //  ! Because generation below th1 is on!!!
00154   xpar[1]=   CmsEne*CmsEne*(1-cos(thmin))/2;  //  ! TransMin [GeV**2]
00155   xpar[2]=   CmsEne*CmsEne*(1-cos(thmax))/2;  //  ! TransMax [GeV**2]
00156   xpar[3]=   m_infraredCut;  //        ! Infrared cut on photon energy
00157 
00158   MARINI(m_initSeed[0], m_initSeed[1], m_initSeed[2]);
00159 
00160   BHLUMI(-1,xpar,npar);
00161 
00162   return StatusCode::SUCCESS;
00163 }

Member Data Documentation

int Bhlumi::m_angleMode [private]
 

double Bhlumi::m_cmEnergy [private]
 

double Bhlumi::m_infraredCut [private]
 

std::vector<int> Bhlumi::m_initSeed [private]
 

std::vector<int> Bhlumi::m_initSeed [private]
 

double Bhlumi::m_maxThetaAngle [private]
 

double Bhlumi::m_minThetaAngle [private]
 

int Bhlumi::npar [private]
 

double Bhlumi::xpar [private]
 

The documentation for this class was generated from the following files:
Generated on Wed Feb 2 15:55:30 2011 for BOSS6.5.5 by  doxygen 1.3.9.1