Hi, Please find attached some example files generated by the new version of my MakeClass suite for the $ROOTSYS/tutorials/hsimple.root file. Some notes on how to use the generated skeleton analysis : 1. You can use the ntuple.[hc]xx directly in root ( as a CINT macro ) : root > .L ntuple.cxx root > ntuple *n = new ntuple() or, you can build the shared library ntuple.so ( an example Makefile for Linux is also attached ), and then : root > .L ntuple.so root > ntuple *n = new ntuple() The second method is preferred, as it gives much better performance. 2. You can "hardcode" your analysis into the ntuple::Loop() function, as in the example file ntuple.cxx, but this makes the problem, that whenever you need to make any small change in it, you need to "reload" the whole analysis ( including creating a new shared library, if you use this method ) : root > .L ntuple.so // or ".L ntuple.cxx" root > ntuple *n = new ntuple() root > n->Loop() 3. You can leave the ntuple::Loop() function EMPTY ( you won't need it later ) and create the ntuple.so ONCE. Then you can copy the ntuple.cxx into the loop.cxx, and modify the ntuple::Loop() skeleton ( in loop.cxx file ) to be an "Int_t loop(ntuple *)" function, and later you can : root > .L ntuple.so root > ntuple *n = new ntuple() root > .L loop.cxx root > loop(n) or simply : root > .L ntuple.so root > ntuple *n = new ntuple() root > .x loop.cxx(n) Note, whenever you change anything in your analysis "loop" you only need to reload the ".L loop.cxx", no changes to the shared library ntuple.so are required ( you get the full compiled code speed in accessing your ntuple mixed with CINT flexibility ). An example loop.cxx is also attached 4. Last, but not least, there is an automatically generated function : Int_t ntuple::Loop(Int_t (*analysis)(ntuple *)) which loops over all events and for every event calls the given "Int_t analysis(ntuple *)" function which is expected to analyze the current ( single ) event. And again, you can mix compiled ntuple.so ( created ONCE ) with interpreted "analysis" : root > .L ntuple.so root > ntuple *n = new ntuple() root > .L analysis.cxx root > n->Loop(analysis) Note again, whenever you change anything in your "analysis", you only need to reload the ".L analysis.cxx", no changes to the shared library ntuple.so are required. An example analysis.cxx is also attached. Warning : current CINT, as shiped with root 2.21/08 has some problems with this method, you probably need to wait till the new root 2.22 is out which has the new CINT in it. Due to current CINT limitations one cannot use code generated by this macro for multi-dimensional arrays directly in CINT ( as a CINT macro ), one needs to compile the analysis code and create the shared library. Have fun, Jacek. // # // # This class has been automatically generated // # (Thu May 20 20:23:25 1999 by ROOT version 2.21/08) // # from TTree ntuple/Demo ntuple // # found on file: hsimple.root // # // # *** The ntuple class interface *** #ifndef ntuple_hxx #define ntuple_hxx // # Required includes #ifndef __CINT__ #include <stdio.h> #include <stream.h> #include "TROOT.h" #include "TTree.h" #include "TFile.h" #include "Api.h" #else class TBranch; class TTree; class TFile; #endif // # class ntuple class ntuple { public: // # public variables TTree *fTree; // pointer to the analysed TTree Int_t fEvent; // current Event number // # declaration of branches Int_t branch_px; TBranch *branch_pointer_px; Int_t branch_py; TBranch *branch_pointer_py; Int_t branch_pz; TBranch *branch_pointer_pz; Int_t branch_random; TBranch *branch_pointer_random; Int_t branch_i; TBranch *branch_pointer_i; // # declaration of leaves Float_t leaf_px; Float_t &px() {if (branch_px) return leaf_px;branch_pointer_px->GetEvent(fEvent);branch_px=1;return leaf_px;} Float_t leaf_py; Float_t &py() {if (branch_py) return leaf_py;branch_pointer_py->GetEvent(fEvent);branch_py=1;return leaf_py;} Float_t leaf_pz; Float_t &pz() {if (branch_pz) return leaf_pz;branch_pointer_pz->GetEvent(fEvent);branch_pz=1;return leaf_pz;} Float_t leaf_random; Float_t &random() {if (branch_random) return leaf_random;branch_pointer_random->GetEvent(fEvent);branch_random=1;return leaf_random;} Float_t leaf_i; Float_t &i() {if (branch_i) return leaf_i;branch_pointer_i->GetEvent(fEvent);branch_i=1;return leaf_i;} // # public functions ntuple(Char_t *filename = 0); ntuple(TTree *tree) {if (!tree) ntuple();else Init(tree);} ~ntuple() {;} Stat_t GetEntries() {if (!fTree) return 0;return fTree->GetEntries();} Int_t GetEvent(Int_t event); void Init(TTree *tree); Int_t Loop(); Int_t Loop(Int_t (*analysis)(ntuple *)); void Show(Int_t event = -1); }; #endif // # *** The ntuple class implementation *** #ifdef ntuple_cxx // # ntuple::ntuple(Char_t *filename) ntuple::ntuple(Char_t *filename) { // # If parameter filename is not specified (or zero), connect the file // # used to generate this class ( hsimple.root ). if (!filename || !(*filename)) { filename = "hsimple.root"; } // # Find the required file in ROOT and, if not found, open it. TFile *f = (TFile*)gROOT->GetListOfFiles()->FindObject(filename); if (!f) { f = new TFile(filename); } // # Find the ntuple tree and initialize object. TTree *t = (TTree*)gDirectory->Get("ntuple"); Init(t); } // # Int_t ntuple::GetEvent(Int_t event) Int_t ntuple::GetEvent(Int_t event) { // # Prepare to read specified event from the Tree. if (!fTree) return 0; fEvent = event; // # Set branch statuses to 0 branch_px = 0; branch_py = 0; branch_pz = 0; branch_random = 0; branch_i = 0; return 1; } // # void ntuple::Init(TTree *tree) void ntuple::Init(TTree *tree) { // # Initialize public variables. fTree = tree; fEvent = -1; if (!tree) return; // # Set branch addresses. fTree->SetBranchAddress("px",&leaf_px); fTree->SetBranchAddress("py",&leaf_py); fTree->SetBranchAddress("pz",&leaf_pz); fTree->SetBranchAddress("random",&leaf_random); fTree->SetBranchAddress("i",&leaf_i); // # Set branch statuses to 0, initialize branch pointers fTree->SetBranchStatus("*",0); // disable all branches branch_px = 0; branch_pointer_px = fTree->GetBranch("px"); branch_py = 0; branch_pointer_py = fTree->GetBranch("py"); branch_pz = 0; branch_pointer_pz = fTree->GetBranch("pz"); branch_random = 0; branch_pointer_random = fTree->GetBranch("random"); branch_i = 0; branch_pointer_i = fTree->GetBranch("i"); } // # Int_t ntuple::Loop(Int_t (*analysis)(ntuple *)) Int_t ntuple::Loop(Int_t (*analysis)(ntuple *)) { // # Execute the analysis function for every event. In case the // # analysis function returns a non zero value, break the loop. if (!fTree || !analysis) return 0; // # Local variables ( common to compiled and interpreted code ). Stat_t nentries = GetEntries(); // Number of entries in the Tree. Int_t nevents = 0; // How many events were analysed. // # This is the event loop in compiled code. #ifndef __CINT__ // # First define some local variables and objects. char temp[64]; // INTERPRETEDFUNC long offset = 0; // INTERPRETEDFUNC G__ClassInfo globalscope; // INTERPRETEDFUNC G__CallFunc func; // INTERPRETEDFUNC, COMPILEDINTERFACEMETHOD, BYTECODEFUNC // # Then execute the loop. switch(G__isinterpretedp2f(((void*)analysis))) { // # using function call as string case G__INTERPRETEDFUNC: sprintf(temp,"(ntuple *)%p",(void*)this); func.SetFunc(&globalscope,(char*)analysis,temp,&offset); for (Int_t i=0; i<nentries; i++) { nevents += GetEvent(i); if (func.ExecInt((void*)NULL)) break; } break; // # using interface method case G__COMPILEDINTERFACEMETHOD: func.SetFunc((G__InterfaceMethod)analysis); func.SetArg(((long)this)); for (Int_t i=0; i<nentries; i++) { nevents += GetEvent(i); if (func.ExecInt((void*)NULL)) break; } break; // # bytecode version of interpreted func case G__BYTECODEFUNC: func.SetBytecode((struct G__bytecodefunc*)analysis); func.SetArg(((long)this)); for (Int_t i=0; i<nentries; i++) { nevents += GetEvent(i); if (func.ExecInt((void*)NULL)) break; } break; // # using true pointer to function case G__COMPILEDTRUEFUNC: // # pointer not in CINT global function table case G__UNKNOWNFUNC: for (Int_t i=0; i<nentries; i++) { nevents += GetEvent(i); if ((*analysis)(this)) break; } break; // # this should never happen ( unknown kind of pointer ) default: cerr << "Error : Unknown kind of pointer to function" << endl; break; } // # This is the event loop in interpreted code. #else // # current CINT cannot deal with this #if 0 for (Int_t i=0; i<nentries; i++) { nevents += GetEvent(i); if ((*analysis)(this)) break; // current CINT cannot deal with it } // # so we need to use this #else Int_t result; for (Int_t i=0; i<nentries; i++) { nevents += GetEvent(i); result = (*analysis)(this); if (result) break; } #endif #endif // # We are done, return. return nevents; } // # void ntuple::Show(Int_t event) void ntuple::Show(Int_t event) { // # Print contents of event ( all branches ). // # If event is not specified, print current event. if (!fTree) return; if (event>=0) fEvent = event; else event = fEvent; if (fEvent<0) return; fTree->SetBranchStatus("*",1); // enable all branches fTree->Show(fEvent); fTree->SetBranchStatus("*",0); // disable all branches // # Set branch statuses to 1 branch_px = 1; branch_py = 1; branch_pz = 1; branch_random = 1; branch_i = 1; } #endif // # End of file ntuple.hxx // # Include the ntuple class interface specification. #ifndef ntuple_cxx #define ntuple_cxx #endif #include "ntuple.hxx" #undef ntuple_cxx // # Place here all ROOT related includes that you need. #ifndef __CINT__ #include "TH1.h" #include "TH2.h" #include "TCanvas.h" #endif // # Int_t ntuple::Loop() Int_t ntuple::Loop() { // # Local variables. Stat_t nentries = GetEntries(); // Number of entries in the Tree. Int_t nevents = 0; // How many events were analysed. // # Create histograms. TH1F *h_px= new TH1F("px","px",100,-4,4); TH1F *h_py= new TH1F("py","py",100,-5,5); TH1F *h_pz= new TH1F("pz","pz",100,0,20); TH1F *h_random= new TH1F("random","random",100,0,1); TH1F *h_i= new TH1F("i","i",100,0,30000); TH2F *h_pxpz= new TH2F("pxpz","pxpz",100,0,20,100,-4,4); // # This is the loop skeleton. for (Int_t l=0; l<nentries; l++) { nevents += GetEvent(l); h_px->Fill(px()); h_py->Fill(py()); h_pz->Fill(pz()); h_random->Fill(random()); h_i->Fill(i()); h_pxpz->Fill(pz(),px()); } // # Show results. TCanvas *c = new TCanvas("c"); c->Divide(2,3); c->cd(1); h_px->Draw(); c->cd(2); h_py->Draw(); c->cd(3); h_pz->Draw(); c->cd(4); h_random->Draw(); c->cd(5); h_i->Draw(); c->cd(6); h_pxpz->Draw(); // # We are done, return. return nevents; } // # End of file ntuple.cxx SrcSuf = cxx IncSuf = hxx ObjSuf = o ExeSuf = DllSuf = so OutPutOpt = -o ROOTLIBS = -L$(ROOTSYS)/lib -lNew -lBase -lCint -lClib -lCont -lFunc \ -lGraf -lGraf3d -lHist -lHtml -lMatrix -lMeta -lMinuit -lNet \ -lPhysics -lPostscript -lProof -lTree -lUnix -lZip ROOTGLIBS = -lGpad -lGui -lGX11 -lX3d # Linux with egcs EGCS = /opt/egcs/pro CXX = g++ CXXFLAGS = -g -O2 -Wall -fno-rtti -fno-exceptions -fPIC -I$(ROOTSYS)/include LD = g++ LDFLAGS = -g -O2 -Wl,-rpath,$(EGCS)/lib:$(ROOTSYS)/lib SOFLAGS = -shared LIBS = $(ROOTLIBS) -lm -ldl -rdynamic GLIBS = $(ROOTLIBS) $(ROOTGLIBS) -L/usr/X11R6/lib \ -lXpm -lX11 -lm -ldl -rdynamic #------------------------------------------------------------------------------ NTUPLEO = ntuple.$(ObjSuf) ntupleDict.$(ObjSuf) NTUPLES = ntuple.$(SrcSuf) ntupleDict.$(SrcSuf) NTUPLESO = ntuple.$(DllSuf) OBJS = $(NTUPLEO) PROGRAMS = $(NTUPLESO) all: $(PROGRAMS) $(NTUPLESO): $(NTUPLEO) $(LD) $(SOFLAGS) $(LDFLAGS) $(NTUPLEO) $(OutPutOpt) $(NTUPLESO) clean: @rm -f $(OBJS) *Dict.* *~ .*~ core .SUFFIXES: .$(SrcSuf) ### ntuple.$(ObjSuf): ntuple.$(SrcSuf) ntuple.$(IncSuf) ntupleDict.$(SrcSuf): ntuple.$(IncSuf) @echo "Generating dictionary ntupleDict..." @$(ROOTSYS)/bin/rootcint -f ntupleDict.$(SrcSuf) -c ntuple.$(IncSuf) .$(SrcSuf).$(ObjSuf): $(CXX) $(CXXFLAGS) -c $< // # Place here all ROOT related includes that you need. #ifndef __CINT__ #include "ntuple.hxx" // Include the ntuple class interface specification. #include "TH1.h" #include "TH2.h" #include "TCanvas.h" #else class ntuple; #endif // # Int_t loop(ntuple *n) Int_t loop(ntuple *n) { // # Local variables. Stat_t nentries = n->GetEntries(); // Number of entries in the Tree. Int_t nevents = 0; // How many events were analysed. // # Create histograms. TH1F *h_px= new TH1F("px","px",100,-4,4); TH1F *h_py= new TH1F("py","py",100,-5,5); TH1F *h_pz= new TH1F("pz","pz",100,0,20); TH1F *h_random= new TH1F("random","random",100,0,1); TH1F *h_i= new TH1F("i","i",100,0,30000); TH2F *h_pxpz= new TH2F("pxpz","pxpz",100,0,20,100,-4,4); // # This is the loop skeleton. for (Int_t l=0; l<nentries; l++) { nevents += n->GetEvent(l); h_px->Fill(n->px()); h_py->Fill(n->py()); h_pz->Fill(n->pz()); h_random->Fill(n->random()); h_i->Fill(n->i()); h_pxpz->Fill(n->pz(),n->px()); } // # Show results. TCanvas *c = new TCanvas("c"); c->Divide(2,3); c->cd(1); h_px->Draw(); c->cd(2); h_py->Draw(); c->cd(3); h_pz->Draw(); c->cd(4); h_random->Draw(); c->cd(5); h_i->Draw(); c->cd(6); h_pxpz->Draw(); // # We are done, return. return nevents; } // # End of file loop.cxx // # Place here all ROOT related includes that you need. #ifndef __CINT__ #include "ntuple.hxx" // Include the ntuple class interface specification. #include "TH1.h" #include "TH2.h" #include "TCanvas.h" #else class ntuple; #endif // # Create histograms. TH1F *h_px= new TH1F("px","px",100,-4,4); TH1F *h_py= new TH1F("py","py",100,-5,5); TH1F *h_pz= new TH1F("pz","pz",100,0,20); TH1F *h_random= new TH1F("random","random",100,0,1); TH1F *h_i= new TH1F("i","i",100,0,30000); TH2F *h_pxpz= new TH2F("pxpz","pxpz",100,0,20,100,-4,4); // # Int_t analysis(ntuple *n) Int_t analysis(ntuple *n) { // # Analyse single event. h_px->Fill(n->px()); h_py->Fill(n->py()); h_pz->Fill(n->pz()); h_random->Fill(n->random()); h_i->Fill(n->i()); h_pxpz->Fill(n->pz(),n->px()); // # We are done, return. return 0; } // # End of file analysis.cxx
This archive was generated by hypermail 2b29 : Tue Jan 04 2000 - 00:43:33 MET