#include "IO/LCReader.h"
#include "IOIMPL/LCFactory.h"
#include "EVENT/MCParticle.h"
#include "EVENT/LCCollection.h"
#include "IMPL/LCEventImpl.h"
#include "UTIL/LCTOOLS.h"
Include dependency graph for readEventTree.C:
Go to the source code of this file.
Functions | |
void | readEventTree (const char *FILEN) |
void readEventTree | ( | const char * | FILEN | ) |
Example script for testing the ROOT LCIO dictionary.
readEventTree: read LCEvents from tree with branch holding complete event e.g. created with writeEventTree.C and fill a histogram from MCParticle::getEnergy()
Definition at line 22 of file readEventTree.C.
00022 { 00023 00024 //just in case this script is executed multiple times 00025 delete gROOT->GetListOfFiles()->FindObject( FILEN ); 00026 delete gROOT->GetListOfCanvases()->FindObject("c1"); 00027 00028 if (!TClassTable::GetDict("IMPL::ReconstructedParticleImpl")) { 00029 unsigned res ; 00030 00031 res = gSystem->Load("$LCIO/lib/liblcio.so"); 00032 res = gSystem->Load("$LCIO/lib/liblcioDict.so"); 00033 } 00034 00035 00036 std::cout << " loaded LCIO library and dictionary ... " << std::endl ; 00037 00038 00039 //======== open file and get tree and branch ===================== 00040 00041 TFile* f = new TFile( FILEN, "READ") ; 00042 if (!f) { return; } 00043 00044 TTree *t; f->GetObject("LCIO",t); 00045 00046 IMPL::LCEventImpl* evt=0 ; 00047 TBranch* bevt = t->GetBranch("LCEvent") ; 00048 00049 if( bevt ){ 00050 00051 bevt->SetAddress( &evt ) ; 00052 00053 } else { 00054 00055 std::cout << " --- branch 'LCEvent' not found ;-( " << std::endl ; 00056 } 00057 00058 00059 // define a simple histogram 00060 00061 TH1F* h_mcpEner = new TH1F("h_mcpEner","MCParticles E(GeV) ", 100, 0. , 100. ) ; 00062 00063 std::string mcpName("MCParticlesSkimmed") ; 00064 00065 00066 int nEvents = 0 ; 00067 00068 int nevt = t->GetEntries(); 00069 00070 for (Int_t i = 0; i < nevt ; i++) { 00071 00072 Long64_t tentry = t->LoadTree(i); 00073 00074 int nbyte = bevt->GetEntry(tentry); 00075 00076 00077 //======== LCIO event loop below this line ============================ 00078 00079 00080 // UTIL::LCTOOLS::dumpEventDetailed( evt ) ; 00081 UTIL::LCTOOLS::dumpEvent( evt ) ; 00082 00083 00084 EVENT::LCCollection* col = evt->getCollection( mcpName ) ; 00085 00086 int nMcp = col->getNumberOfElements() ; 00087 for( int j = 0 ; j < nMcp ; ++j ) { 00088 00089 EVENT::MCParticle* mcp = (EVENT::MCParticle*) col->getElementAt(j) ; 00090 00091 h_mcpEner->Fill( mcp->getEnergy() ) ; 00092 00093 } 00094 00095 nEvents ++ ; 00096 00097 } 00098 00099 // -------- end of event loop ----------- 00100 00101 std::cout << std::endl 00102 << " " << nEvents 00103 << " events read from file: " 00104 << FILEN << std::endl ; 00105 00106 00107 00108 00109 //======== draw histogram(s) ============================ 00110 00111 // ------------------------- 00112 TCanvas *c1 = new TCanvas("c1","LCIO root example",200,10,700,500); 00113 // c1->Divide(2,2); 00114 c1->cd(1) ; 00115 h_mcpEner->Draw(); 00116 00117 }