CMS 3D CMS Logo

ZdcTestAnalysis Class Reference

#include <SimG4CMS/Forward/interface/ZdcTestAnalysis.h>

Inheritance diagram for ZdcTestAnalysis:

SimWatcher Observer< const BeginOfJob * > Observer< const BeginOfRun * > Observer< const EndOfRun * > Observer< const BeginOfEvent * > Observer< const EndOfEvent * > Observer< const G4Step * >

List of all members.

Public Member Functions

 ZdcTestAnalysis (const edm::ParameterSet &p)
virtual ~ZdcTestAnalysis ()

Private Member Functions

void finish ()
void update (const G4Step *step)
 This routine will be called when the appropriate signal arrives.
void update (const EndOfEvent *evt)
 This routine will be called when the appropriate signal arrives.
void update (const BeginOfEvent *evt)
 This routine will be called when the appropriate signal arrives.
void update (const EndOfRun *run)
 This routine will be called when the appropriate signal arrives.
void update (const BeginOfRun *run)
 This routine will be called when the appropriate signal arrives.
void update (const BeginOfJob *run)
 This routine will be called when the appropriate signal arrives.

Private Attributes

int doNTzdcevent
int doNTzdcstep
int eventIndex
std::string eventNtFileName
int stepIndex
std::string stepNtFileName
int verbosity
Float_t zdceventarray [15]
TNtuple * zdceventntuple
TFile * zdcOutputEventFile
TFile * zdcOutputStepFile
Float_t zdcsteparray [18]
TNtuple * zdcstepntuple


Detailed Description

Definition at line 65 of file ZdcTestAnalysis.h.


Constructor & Destructor Documentation

ZdcTestAnalysis::ZdcTestAnalysis ( const edm::ParameterSet p  ) 

Definition at line 25 of file ZdcTestAnalysis.cc.

References GenMuonPlsPt100GeV_cfg::cout, doNTzdcevent, doNTzdcstep, lat::endl(), eventNtFileName, edm::ParameterSet::getParameter(), stepNtFileName, verbosity, zdceventntuple, and zdcstepntuple.

00025                                                         {
00026   //constructor
00027   edm::ParameterSet m_Anal = p.getParameter<edm::ParameterSet>("ZdcTestAnalysis");
00028   verbosity    = m_Anal.getParameter<int>("Verbosity");
00029   doNTzdcstep  = m_Anal.getParameter<int>("StepNtupleFlag");
00030   doNTzdcevent  = m_Anal.getParameter<int>("EventNtupleFlag");
00031   stepNtFileName = m_Anal.getParameter<std::string>("StepNtupleFileName");
00032   eventNtFileName = m_Anal.getParameter<std::string>("EventNtupleFileName");
00033    
00034    if (verbosity > 0)
00035    std::cout<<std::endl;
00036    std::cout<<"============================================================================"<<std::endl;
00037    std::cout << "ZdcTestAnalysis:: Initialized as observer"<< std::endl;
00038    if (doNTzdcstep  > 0){
00039      std::cout <<" Step Ntuple will be created"<< std::endl;
00040      std::cout <<" Step Ntuple file: "<<stepNtFileName<<std::endl;
00041    }
00042    if (doNTzdcevent > 0){
00043      std::cout <<" Event Ntuple will be created"<< std::endl;
00044      std::cout <<" Step Ntuple file: "<<stepNtFileName<<std::endl;
00045    }
00046    std::cout<<"============================================================================"<<std::endl;
00047    std::cout<<std::endl;
00048 
00049    if (doNTzdcstep  > 0)
00050      zdcstepntuple = 
00051        new TNtuple("NTzdcstep","NTzdcstep",
00052                    "evt:trackid:charge:pdgcode:x:y:z:stepl:stepe:eta:phi:vpx:vpy:vpz:idx:idl:pvtype:ncherphot");
00053 
00054    if (doNTzdcevent  > 0)
00055      zdceventntuple = 
00056        new TNtuple("NTzdcevent","NTzdcevent",
00057                    "evt:ihit:fiberid:zside:subdet:layer:fiber:channel:enem:enhad:hitenergy:x:y:z:time:etot");
00058 
00059    //theZdcSD = new ZdcSD("ZDCHITSB", new ZdcNumberingScheme());
00060 }

ZdcTestAnalysis::~ZdcTestAnalysis (  )  [virtual]

Definition at line 62 of file ZdcTestAnalysis.cc.

References GenMuonPlsPt100GeV_cfg::cout, lat::endl(), finish(), and verbosity.

00062                                   {
00063   // destructor
00064   finish();
00065   if (verbosity > 0) {
00066     std::cout << std::endl << "ZdcTestAnalysis Dextructor  -------->  End of ZdcTestAnalysis : "
00067       << std::endl << std::endl; 
00068   }
00069 
00070   //if (doNTzdcstep  > 0)delete zdcstepntuple;
00071   //if (doNTzdcevent > 0)delete zdceventntuple;
00072 
00073   std::cout<<"ZdcTestAnalysis: End of process"<<std::endl;
00074 }


Member Function Documentation

void ZdcTestAnalysis::finish ( void   )  [private]

Definition at line 397 of file ZdcTestAnalysis.cc.

References GenMuonPlsPt100GeV_cfg::cout, doNTzdcevent, doNTzdcstep, lat::endl(), eventIndex, zdceventntuple, zdcOutputEventFile, zdcOutputStepFile, and zdcstepntuple.

Referenced by ~ZdcTestAnalysis().

00397                             {
00398   if (doNTzdcstep) {
00399     zdcOutputStepFile->cd();
00400     zdcstepntuple->Write();
00401     std::cout << "ZdcTestAnalysis: Ntuple step  written for event: "<<eventIndex<<std::endl;
00402     zdcOutputStepFile->Close();
00403     std::cout << "ZdcTestAnalysis: Step file closed" << std::endl;
00404   }
00405   
00406 
00407  if (doNTzdcevent) {
00408    zdcOutputEventFile->cd();
00409    zdceventntuple->Write("",TObject::kOverwrite);
00410    std::cout << "ZdcTestAnalysis: Ntuple event written for event: "<<eventIndex<<std::endl;   
00411    zdcOutputEventFile->Close();
00412    std::cout << "ZdcTestAnalysis: Event file closed" << std::endl;
00413  }
00414  
00415 }

void ZdcTestAnalysis::update ( const G4Step *   )  [private, virtual]

This routine will be called when the appropriate signal arrives.

Implements Observer< const G4Step * >.

Definition at line 115 of file ZdcTestAnalysis.cc.

References GenMuonPlsPt100GeV_cfg::cout, doNTzdcstep, lat::endl(), eta, eventIndex, funct::log(), ntzdcs_charge, ntzdcs_eta, ntzdcs_evt, ntzdcs_idl, ntzdcs_idx, ntzdcs_ncherphot, ntzdcs_pdgcode, ntzdcs_phi, ntzdcs_pvtype, ntzdcs_stepe, ntzdcs_stepl, ntzdcs_trackid, ntzdcs_vpx, ntzdcs_vpy, ntzdcs_vpz, ntzdcs_x, ntzdcs_y, ntzdcs_z, phi, stepIndex, verbosity, zdcsteparray, and zdcstepntuple.

00115                                                  {
00116   //step;
00117   stepIndex++;
00118 
00119   if (doNTzdcstep) {
00120     G4StepPoint * preStepPoint = aStep->GetPreStepPoint();
00121     // G4StepPoint * postStepPoint= aStep->GetPostStepPoint();
00122     G4double stepL = aStep->GetStepLength();
00123     G4double stepE = aStep->GetTotalEnergyDeposit();
00124     
00125     if (verbosity >= 2)
00126       std::cout << "Step " << stepL << "," << stepE <<std::endl;
00127     
00128     G4Track * theTrack    = aStep->GetTrack();
00129     G4int theTrackID      = theTrack->GetTrackID();
00130     G4double theCharge    = theTrack->GetDynamicParticle()->GetCharge();
00131     G4String particleType = theTrack->GetDefinition()->GetParticleName();
00132     G4int pdgcode         = theTrack->GetDefinition()->GetPDGEncoding();
00133     
00134     G4ThreeVector vert_mom = theTrack->GetVertexMomentumDirection();
00135     G4double vpx = vert_mom.x();
00136     G4double vpy = vert_mom.y();
00137     G4double vpz = vert_mom.z();
00138     double eta = 0.5 * log( (1.+vpz) / (1.-vpz) );
00139     double phi = atan2(vpy,vpx);
00140     
00141     G4ThreeVector hitPoint = preStepPoint->GetPosition();
00142     G4ThreeVector localPoint = theTrack->GetTouchable()->GetHistory()->
00143     GetTopTransform().TransformPoint(hitPoint);
00144     
00145     const G4VTouchable* touch = aStep->GetPreStepPoint()->GetTouchable();
00146     int idx = touch->GetReplicaNumber(0);
00147     int idLayer = -1;
00148     int thePVtype = -1;
00149     
00150     int historyDepth = touch->GetHistoryDepth();
00151 
00152     if (historyDepth > 0) {
00153       std::vector<int>                theReplicaNumbers(historyDepth);
00154       std::vector<G4VPhysicalVolume*> thePhysicalVolumes(historyDepth);
00155       std::vector<G4String>           thePVnames(historyDepth);
00156       std::vector<G4LogicalVolume*>   theLogicalVolumes(historyDepth);
00157       std::vector<G4String>           theLVnames(historyDepth);
00158       std::vector<G4Material*>        theMaterials(historyDepth);
00159       std::vector<G4String>           theMaterialNames(historyDepth);
00160       
00161       for (int jj = 0; jj < historyDepth; jj++) {
00162         theReplicaNumbers[jj] = touch->GetReplicaNumber(jj);
00163         thePhysicalVolumes[jj] = touch->GetVolume(jj);
00164         thePVnames[jj] = thePhysicalVolumes[jj]->GetName();
00165         theLogicalVolumes[jj] = thePhysicalVolumes[jj]->GetLogicalVolume();
00166         theLVnames[jj] = theLogicalVolumes[jj]->GetName();
00167         theMaterials[jj] = theLogicalVolumes[jj]->GetMaterial();
00168         theMaterialNames[jj] = theMaterials[jj]->GetName();
00169         if (verbosity >= 2)
00170           std::cout << " GHD " << jj << ": " << theReplicaNumbers[jj] << ","
00171                        << thePVnames[jj] << "," << theLVnames[jj] << ","
00172                     << theMaterialNames[jj]  << std::endl;
00173       }
00174 
00175       idLayer = theReplicaNumbers[1];
00176       if (thePVnames[0] == "ZDC_EMLayer")
00177         thePVtype = 1;
00178       else if (thePVnames[0] == "ZDC_EMAbsorber")
00179         thePVtype = 2;
00180       else if (thePVnames[0] == "ZDC_EMFiber")
00181         thePVtype = 3;
00182       else if (thePVnames[0] == "ZDC_HadLayer")
00183         thePVtype = 7;
00184       else if (thePVnames[0] == "ZDC_HadAbsorber")
00185         thePVtype = 8;
00186       else if (thePVnames[0] == "ZDC_HadFiber")
00187         thePVtype = 9;
00188       else if (thePVnames[0] == "ZDC_PhobosLayer")
00189         thePVtype = 11;
00190       else if (thePVnames[0] == "ZDC_PhobosAbsorber")
00191         thePVtype = 12;
00192       else if (thePVnames[0] == "ZDC_PhobosFiber")
00193       thePVtype = 13;
00194       else {
00195         thePVtype = 0;
00196         if (verbosity >= 2)
00197           std::cout << " pvtype=0 hd=" << historyDepth << " " << theReplicaNumbers[0] << ","
00198                     << thePVnames[0] << "," << theLVnames[0] << "," << theMaterialNames[0] << std::endl;
00199       }
00200     }    
00201     else if (historyDepth == 0) { 
00202       int theReplicaNumber = touch->GetReplicaNumber(0);
00203       G4VPhysicalVolume* thePhysicalVolume = touch->GetVolume(0);
00204       G4String thePVname = thePhysicalVolume->GetName();
00205       G4LogicalVolume * theLogicalVolume = thePhysicalVolume->GetLogicalVolume();
00206       G4String theLVname = theLogicalVolume->GetName();
00207       G4Material * theMaterial = theLogicalVolume->GetMaterial();
00208       G4String theMaterialName = theMaterial->GetName();
00209       if (verbosity >= 2)
00210         std::cout << " hd=0 " << theReplicaNumber << "," 
00211                   << thePVname << "," << theLVname << "," 
00212                   << theMaterialName << std::endl;
00213     }
00214     else {
00215       std::cout << " hd<0:  hd=" << historyDepth << std::endl;
00216     }
00217 
00218     
00219     double NCherPhot = -1;
00220     zdcsteparray[ntzdcs_evt] = (float)eventIndex;
00221     zdcsteparray[ntzdcs_trackid] = (float)theTrackID;
00222     zdcsteparray[ntzdcs_charge] = theCharge;
00223     zdcsteparray[ntzdcs_pdgcode] = (float)pdgcode;
00224     zdcsteparray[ntzdcs_x] = hitPoint.x();
00225     zdcsteparray[ntzdcs_y] = hitPoint.y();
00226     zdcsteparray[ntzdcs_z] = hitPoint.z();
00227     zdcsteparray[ntzdcs_stepl] = stepL;
00228     zdcsteparray[ntzdcs_stepe] = stepE;
00229     zdcsteparray[ntzdcs_eta] = eta;
00230     zdcsteparray[ntzdcs_phi] = phi;
00231     zdcsteparray[ntzdcs_vpx] = vpx;
00232     zdcsteparray[ntzdcs_vpy] = vpy;
00233     zdcsteparray[ntzdcs_vpz] = vpz;
00234     zdcsteparray[ntzdcs_idx] = (float)idx;
00235     zdcsteparray[ntzdcs_idl] = (float)idLayer;
00236     zdcsteparray[ntzdcs_pvtype] = thePVtype;
00237     zdcsteparray[ntzdcs_ncherphot] = NCherPhot;
00238     zdcstepntuple->Fill(zdcsteparray);
00239 
00240   }
00241 }

void ZdcTestAnalysis::update ( const EndOfEvent  )  [private, virtual]

This routine will be called when the appropriate signal arrives.

Implements Observer< const EndOfEvent * >.

Definition at line 243 of file ZdcTestAnalysis.cc.

References GenMuonPlsPt100GeV_cfg::cout, doNTzdcevent, lat::endl(), eta, eventIndex, CaloG4Hit::getEM(), CaloG4Hit::getEnergyDeposit(), CaloG4Hit::getHadr(), CaloG4Hit::getPosition(), CaloG4Hit::getTimeSliceID(), CaloG4Hit::getTrackID(), CaloG4Hit::getUnitID(), i, int, funct::log(), python::multivaluedict::map(), max, min, npart, ntzdce_channel, ntzdce_enem, ntzdce_enhad, ntzdce_etot, ntzdce_evt, ntzdce_fiber, ntzdce_fiberid, ntzdce_hitenergy, ntzdce_ihit, ntzdce_layer, ntzdce_subdet, ntzdce_time, ntzdce_x, ntzdce_y, ntzdce_z, ntzdce_zside, phi, funct::pow(), funct::sqrt(), stepIndex, funct::tan(), theta, ZdcNumberingScheme::unpackZdcIndex(), verbosity, zdceventarray, and zdceventntuple.

00243                                                    {
00244   //end of event 
00245    
00246   // Look for the Hit Collection
00247   std::cout  << std::endl << "ZdcTest::upDate(const EndOfEvent * evt) - event #" << (*evt)()->GetEventID()
00248              << std::endl << "  # of aSteps followed in event = " << stepIndex << std::endl;
00249   
00250   // access to the G4 hit collections
00251   G4HCofThisEvent* allHC = (*evt)()->GetHCofThisEvent();
00252   std::cout << "  accessed all HC";
00253   
00254   int theZDCHCid = G4SDManager::GetSDMpointer()->GetCollectionID("ZDCHITS");
00255   std::cout << " - theZDCHCid = " << theZDCHCid;
00256   
00257   CaloG4HitCollection* theZDCHC = (CaloG4HitCollection*) allHC->GetHC(theZDCHCid);
00258   std::cout << " - theZDCHC = " << theZDCHC << std::endl;
00259   
00260   ZdcNumberingScheme * theZdcNumScheme = new ZdcNumberingScheme(1);
00261   
00262   float ETot=0., SEnergy=0.;
00263   int maxTime=0;
00264   int fiberID=0;
00265   unsigned int unsignedfiberID=0;
00266   std::map<int,float,std::less<int> > energyInFibers;
00267   std::map<int,float,std::less<int> > primaries;
00268   float totalEnergy = 0;
00269   int nentries = theZDCHC->entries();
00270   std::cout << "  theZDCHC has " << nentries << " entries" << std::endl;
00271 
00272   if (doNTzdcevent) {
00273     if (nentries > 0) {
00274       for (int ihit = 0; ihit < nentries; ihit++) {
00275         CaloG4Hit* caloHit = (*theZDCHC)[ihit];
00276         totalEnergy += caloHit->getEnergyDeposit();
00277       }
00278 
00279       for (int ihit = 0; ihit < nentries; ihit++) {
00280         CaloG4Hit* aHit = (*theZDCHC)[ihit];
00281             fiberID = aHit->getUnitID();
00282             unsignedfiberID = aHit->getUnitID();
00283             double enEm = aHit->getEM();
00284             double enHad = aHit->getHadr();
00285             math::XYZPoint hitPoint = aHit->getPosition();
00286             double hitEnergy = aHit->getEnergyDeposit();
00287             if (verbosity >= 1)
00288               std::cout << " entry #" << ihit << ": fiberID=0x" << std::hex 
00289                         << fiberID << std::dec << "; enEm=" << enEm 
00290                         << "; enHad=" << enHad << "; hitEnergy=" << hitEnergy
00291                         << "z=" << hitPoint.z() << std::endl;
00292             energyInFibers[fiberID]+= enEm + enHad;
00293             primaries[aHit->getTrackID()]+= enEm + enHad;
00294             float time = aHit->getTimeSliceID();
00295             if (time > maxTime) maxTime=(int)time;
00296             
00297             int thesubdet, thelayer, thefiber, thechannel, thez;
00298             theZdcNumScheme->unpackZdcIndex(fiberID, thesubdet, thelayer, thefiber, thechannel, thez);
00299             int unsignedsubdet, unsignedlayer, unsignedfiber, unsignedchannel, unsignedz;
00300             theZdcNumScheme->unpackZdcIndex(unsignedfiberID, unsignedsubdet, 
00301                                             unsignedlayer, unsignedfiber, unsignedchannel, unsignedz);
00302             
00303             // unsigned int packidx1 = packZdcIndex(thesubdet, thelayer, thefiber, thechannel, thez);
00304             // unsigned int packidx1 = packZdcIndex(thesubdet, thelayer, thefiber, thechannel, thez);
00305             // unsigned int packidx1 = packZdcIndex(thesubdet, thelayer, thefiber, thechannel, thez);
00306             // unsigned int packidx1 = packZdcIndex(thesubdet, thelayer, thefiber, thechannel, thez);
00307 
00308             zdceventarray[ntzdce_evt] = (float)eventIndex;
00309             zdceventarray[ntzdce_ihit] = (float)ihit;
00310             zdceventarray[ntzdce_fiberid] = (float)fiberID;
00311             zdceventarray[ntzdce_zside] = (float)thez;
00312             zdceventarray[ntzdce_subdet] = (float)thesubdet;
00313             zdceventarray[ntzdce_layer] = (float)thelayer;
00314             zdceventarray[ntzdce_fiber] = (float)thefiber;
00315             zdceventarray[ntzdce_channel] = (float)thechannel;
00316             zdceventarray[ntzdce_enem] = enEm;
00317             zdceventarray[ntzdce_enhad] = enHad;
00318             zdceventarray[ntzdce_hitenergy] = hitEnergy;
00319             zdceventarray[ntzdce_x] = hitPoint.x();
00320             zdceventarray[ntzdce_y] = hitPoint.y();
00321             zdceventarray[ntzdce_z] = hitPoint.z();
00322             zdceventarray[ntzdce_time] = time;
00323             zdceventarray[ntzdce_etot] = totalEnergy;
00324             zdceventntuple->Fill(zdceventarray);
00325       }
00326       
00327       for (std::map<int,float,std::less<int> >::iterator is = energyInFibers.begin();
00328            is!= energyInFibers.end(); is++)
00329         {
00330           ETot = (*is).second;
00331           SEnergy += ETot;
00332         }
00333 
00334       // Find Primary info:
00335       int trackID = 0;
00336       int particleType = 0;
00337       G4PrimaryParticle* thePrim=0;
00338       G4int nvertex = (*evt)()->GetNumberOfPrimaryVertex();
00339       std::cout << "Event has " << nvertex << " vertex" << std::endl;
00340       if (nvertex==0)
00341         std::cout << "ZdcTest End Of Event  ERROR: no vertex" << std::endl;
00342       
00343       for (int i = 0 ; i<nvertex; i++) {
00344         
00345         G4PrimaryVertex* avertex = (*evt)()->GetPrimaryVertex(i);
00346         if (avertex == 0)
00347           std::cout << "ZdcTest End Of Event ERR: pointer to vertex = 0"
00348                        << std::endl;
00349         std::cout << "Vertex number :" <<i << std::endl;
00350         int npart = avertex->GetNumberOfParticle();
00351         if (npart ==0)
00352           std::cout << "ZdcTest End Of Event ERR: no primary!" << std::endl;
00353         if (thePrim==0) thePrim=avertex->GetPrimary(trackID);
00354       }
00355       
00356       double px=0.,py=0.,pz=0.;
00357       double eta = 0., phi = 0., pInit = 0.;
00358       
00359       if (thePrim != 0) {
00360         px = thePrim->GetPx();
00361         py = thePrim->GetPy();
00362         pz = thePrim->GetPz();
00363         pInit = sqrt(pow(px,2.)+pow(py,2.)+pow(pz,2.));
00364         if (pInit==0) {
00365           std::cout << "ZdcTest End Of Event  ERR: primary has p=0 " << std::endl;
00366         } else {
00367           float costheta = pz/pInit;
00368           float theta = acos(std::min(std::max(costheta,float(-1.)),float(1.)));
00369           eta = -log(tan(theta/2));
00370           
00371           if (px != 0) phi = atan(py/px);
00372         }
00373         particleType      = thePrim->GetPDGcode();
00374       } else {
00375         std::cout << "ZdcTest End Of Event ERR: could not find primary "
00376                      << std::endl;
00377       }
00378       
00379     } // nentries > 0
00380     
00381  
00382  } // createNTzdcevent
00383   
00384   int iEvt = (*evt)()->GetEventID();
00385   if (iEvt < 10)
00386     std::cout << " ZdcTest Event " << iEvt << std::endl;
00387   else if ((iEvt < 100) && (iEvt%10 == 0))
00388     std::cout << " ZdcTest Event " << iEvt << std::endl;
00389   else if ((iEvt < 1000) && (iEvt%100 == 0))
00390     std::cout << " ZdcTest Event " << iEvt << std::endl;
00391   else if ((iEvt < 10000) && (iEvt%1000 == 0))
00392     std::cout << " ZdcTest Event " << iEvt << std::endl;
00393 }

void ZdcTestAnalysis::update ( const BeginOfEvent  )  [private, virtual]

This routine will be called when the appropriate signal arrives.

Implements Observer< const BeginOfEvent * >.

Definition at line 107 of file ZdcTestAnalysis.cc.

References GenMuonPlsPt100GeV_cfg::cout, lat::endl(), eventIndex, and stepIndex.

00107                                                      {
00108   //event
00109   std::cout << "ZdcTest: Processing Event Number: "<<eventIndex<< std::endl;
00110   eventIndex++;
00111   stepIndex = 0;
00112 }

void ZdcTestAnalysis::update ( const EndOfRun  )  [private, virtual]

This routine will be called when the appropriate signal arrives.

Implements Observer< const EndOfRun * >.

Definition at line 395 of file ZdcTestAnalysis.cc.

00395 {;}

void ZdcTestAnalysis::update ( const BeginOfRun  )  [private, virtual]

This routine will be called when the appropriate signal arrives.

Implements Observer< const BeginOfRun * >.

Definition at line 84 of file ZdcTestAnalysis.cc.

References GenMuonPlsPt100GeV_cfg::cout, doNTzdcevent, doNTzdcstep, lat::endl(), eventIndex, eventNtFileName, stepNtFileName, zdcOutputEventFile, and zdcOutputStepFile.

00084                                                    {
00085   //run
00086 
00087  std::cout << std::endl << "ZdcTestAnalysis: Begining of Run"<< std::endl; 
00088   if (doNTzdcstep) { 
00089     std::cout << "ZDCTestAnalysis: output step file created"<< std::endl;
00090     TString stepfilename = stepNtFileName;
00091     zdcOutputStepFile = new TFile(stepfilename,"RECREATE");
00092 
00093   }
00094   
00095   if (doNTzdcevent) {
00096     std::cout << "ZDCTestAnalysis: output event file created"<< std::endl;
00097     TString stepfilename = eventNtFileName;
00098     zdcOutputEventFile = new TFile(stepfilename,"RECREATE");
00099   }
00100 
00101   eventIndex = 0;
00102 }

void ZdcTestAnalysis::update ( const BeginOfJob  )  [private, virtual]

This routine will be called when the appropriate signal arrives.

Implements Observer< const BeginOfJob * >.

Definition at line 77 of file ZdcTestAnalysis.cc.

References GenMuonPlsPt100GeV_cfg::cout, and lat::endl().

00077                                                    {
00078   //job
00079   std::cout<<"beggining of job"<<std::endl;;
00080 }


Member Data Documentation

int ZdcTestAnalysis::doNTzdcevent [private]

Definition at line 92 of file ZdcTestAnalysis.h.

Referenced by finish(), update(), and ZdcTestAnalysis().

int ZdcTestAnalysis::doNTzdcstep [private]

Definition at line 91 of file ZdcTestAnalysis.h.

Referenced by finish(), update(), and ZdcTestAnalysis().

int ZdcTestAnalysis::eventIndex [private]

Definition at line 102 of file ZdcTestAnalysis.h.

Referenced by finish(), and update().

std::string ZdcTestAnalysis::eventNtFileName [private]

Definition at line 94 of file ZdcTestAnalysis.h.

Referenced by update(), and ZdcTestAnalysis().

int ZdcTestAnalysis::stepIndex [private]

Definition at line 103 of file ZdcTestAnalysis.h.

Referenced by update().

std::string ZdcTestAnalysis::stepNtFileName [private]

Definition at line 93 of file ZdcTestAnalysis.h.

Referenced by update(), and ZdcTestAnalysis().

int ZdcTestAnalysis::verbosity [private]

Definition at line 90 of file ZdcTestAnalysis.h.

Referenced by update(), ZdcTestAnalysis(), and ~ZdcTestAnalysis().

Float_t ZdcTestAnalysis::zdceventarray[15] [private]

Definition at line 106 of file ZdcTestAnalysis.h.

Referenced by update().

TNtuple* ZdcTestAnalysis::zdceventntuple [private]

Definition at line 100 of file ZdcTestAnalysis.h.

Referenced by finish(), update(), and ZdcTestAnalysis().

TFile* ZdcTestAnalysis::zdcOutputEventFile [private]

Definition at line 96 of file ZdcTestAnalysis.h.

Referenced by finish(), and update().

TFile* ZdcTestAnalysis::zdcOutputStepFile [private]

Definition at line 97 of file ZdcTestAnalysis.h.

Referenced by finish(), and update().

Float_t ZdcTestAnalysis::zdcsteparray[18] [private]

Definition at line 105 of file ZdcTestAnalysis.h.

Referenced by update().

TNtuple* ZdcTestAnalysis::zdcstepntuple [private]

Definition at line 99 of file ZdcTestAnalysis.h.

Referenced by finish(), update(), and ZdcTestAnalysis().


The documentation for this class was generated from the following files:
Generated on Tue Jun 9 18:36:03 2009 for CMSSW by  doxygen 1.5.4