Hello,
I am trying to write a macro.C to replicate in ROOT a code I wrote in PyRoot using TChain and MakeClass as follows
TFileCollection *fc=new TFileCollection("fc","fc","myfiles.txt")
TChain *c=new TChain("T")
c->AddFileInfoList(fc->GetList())
c->MakeClass("mymacro")
I am a beginner, so I am struggling with something probably quite simple. I attach below my macro.h and mymacro.C:
//////////////////////////////////////////////////////////
// This class has been automatically generated on
// Tue Feb 28 16:31:07 2017 by ROOT version 6.06/06
// from TChain Tree/
//////////////////////////////////////////////////////////
#ifndef mymacro_h
#define mymacro_h
#include <TROOT.h>
#include <TChain.h>
#include <TFile.h>
// Header file for the classes stored in the TTree if any.
class mymacro {
public :
TTree *fChain; //!pointer to the analyzed TTree or TChain
Int_t fCurrent; //!current Tree number in a TChain
// Fixed size dimensions of array or collections stored in the TTree if any.
// Declaration of leaf types
Int_t run;
Int_t evt;
Int_t bunchid;
Int_t bunchtype;
Int_t trighlt1;
Int_t trigmu;
Int_t l0dim;
Int_t l0m;
Int_t stripphi;
Int_t nlong;
Int_t nvelo;
Int_t down;
Int_t ntk;
Int_t nspd;
Int_t ntotphot;
Int_t nout;
Int_t nphotons;
Float_t etprev1;
Float_t etprev2;
Int_t tt[10]; //[nout]
Float_t tx[10]; //[nout]
Float_t ty[10]; //[nout]
Float_t tz[10]; //[nout]
Float_t px[10]; //[nout]
Float_t py[10]; //[nout]
Float_t pz[10]; //[nout]
Float_t pt[10]; //[nout]
Int_t charge[10]; //[nout]
Int_t muon[10]; //[nout]
Int_t muacc[10]; //[nout]
Float_t pide[10]; //[nout]
Float_t pidk[10]; //[nout]
Float_t pidpi[10]; //[nout]
Float_t pidp[10]; //[nout]
Float_t pidm[10]; //[nout]
Float_t ecal[10]; //[nout]
Float_t hcal[10]; //[nout]
Float_t eoverp[10]; //[nout]
Float_t phx[6]; //[nphotons]
Float_t phy[6]; //[nphotons]
Float_t phz[6]; //[nphotons]
Int_t spd[6]; //[nphotons]
Int_t prs[6]; //[nphotons]
Int_t nl0;
Float_t l0et[3]; //[nl0]
Float_t l0phi[3]; //[nl0]
Float_t l0theta[3]; //[nl0]
Int_t nher;
Int_t herchan[30]; //[nher]
Int_t heradc[30]; //[nher]
// List of branches
TBranch *b_run; //!
TBranch *b_evt; //!
TBranch *b_bunchid; //!
TBranch *b_bunchtype; //!
TBranch *b_trighlt1; //!
TBranch *b_trigmu; //!
TBranch *b_l0dim; //!
TBranch *b_l0m; //!
TBranch *b_stripphi; //!
TBranch *b_nlong; //!
TBranch *b_nvelo; //!
TBranch *b_ndown; //!
TBranch *b_ntk; //!
TBranch *b_nspd; //!
TBranch *b_ntotphot; //!
TBranch *b_nout; //!
TBranch *b_nphotons; //!
TBranch *b_etprev1; //!
TBranch *b_etprev2; //!
TBranch *b_tt; //!
TBranch *b_tx; //!
TBranch *b_ty; //!
TBranch *b_tz; //!
TBranch *b_px; //!
TBranch *b_py; //!
TBranch *b_pz; //!
TBranch *b_pt; //!
TBranch *b_charge; //!
TBranch *b_muon; //!
TBranch *b_muacc; //!
TBranch *b_pide; //!
TBranch *b_pidk; //!
TBranch *b_pidpi; //!
TBranch *b_pidp; //!
TBranch *b_pidm; //!
TBranch *b_ecal; //!
TBranch *b_hcal; //!
TBranch *b_eoverp; //!
TBranch *b_phx; //!
TBranch *b_phy; //!
TBranch *b_phz; //!
TBranch *b_spd; //!
TBranch *b_prs; //!
TBranch *b_nl0; //!
TBranch *b_l0et; //!
TBranch *b_l0phi; //!
TBranch *b_l0theta; //!
TBranch *b_nher; //!
TBranch *b_herchan; //!
TBranch *b_heradc; //!
mymacro(TTree *tree=0);
virtual ~mymacro();
virtual Int_t Cut(Long64_t entry);
virtual Int_t GetEntry(Long64_t entry);
virtual Long64_t LoadTree(Long64_t entry);
virtual void Init(TTree *tree);
virtual void Loop();
virtual Bool_t Notify();
virtual void Show(Long64_t entry = -1);
};
#endif
#ifdef mymacro_cxx
mymacro::mymacro(TTree *tree) : fChain(0)
{
// if parameter tree is not specified (or zero), connect the file
// used to generate this class and read the Tree.
if (tree == 0) {
#ifdef SINGLE_TREE
// The following code should be used if you want this class to access
// a single tree instead of a chain
TFile *f = (TFile*)gROOT->GetListOfFiles()->FindObject("Memory Directory");
if (!f || !f->IsOpen()) {
f = new TFile("Memory Directory");
}
f->GetObject("Tree",tree);
#else // SINGLE_TREE
// The following code should be used if you want this class to access a chain
// of trees.
TChain * chain = new TChain("Tree","");
chain->Add("file:///afs/cern.ch/work/b/bldelane/private/DST/Analysis/temp_root/365/EW26_selection.root/Tree");
tree = chain;
#endif // SINGLE_TREE
}
Init(tree);
}
mymacro::~mymacro()
{
if (!fChain) return;
delete fChain->GetCurrentFile();
}
Int_t mymacro::GetEntry(Long64_t entry)
{
// Read contents of entry.
if (!fChain) return 0;
return fChain->GetEntry(entry);
}
Long64_t mymacro::LoadTree(Long64_t entry)
{
// Set the environment to read one entry
if (!fChain) return -5;
Long64_t centry = fChain->LoadTree(entry);
if (centry < 0) return centry;
if (fChain->GetTreeNumber() != fCurrent) {
fCurrent = fChain->GetTreeNumber();
Notify();
}
return centry;
}
void mymacro::Init(TTree *tree)
{
// The Init() function is called when the selector needs to initialize
// a new tree or chain. Typically here the branch addresses and branch
// pointers of the tree will be set.
// It is normally not necessary to make changes to the generated
// code, but the routine can be extended by the user if needed.
// Init() will be called many times when running on PROOF
// (once per file to be processed).
// Set branch addresses and branch pointers
if (!tree) return;
fChain = tree;
fCurrent = -1;
fChain->SetMakeClass(1);
fChain->SetBranchAddress("run", &run, &b_run);
fChain->SetBranchAddress("evt", &evt, &b_evt);
fChain->SetBranchAddress("bunchid", &bunchid, &b_bunchid);
fChain->SetBranchAddress("bunchtype", &bunchtype, &b_bunchtype);
fChain->SetBranchAddress("trighlt1", &trighlt1, &b_trighlt1);
fChain->SetBranchAddress("trigmu", &trigmu, &b_trigmu);
fChain->SetBranchAddress("l0dim", &l0dim, &b_l0dim);
fChain->SetBranchAddress("l0m", &l0m, &b_l0m);
fChain->SetBranchAddress("stripphi", &stripphi, &b_stripphi);
fChain->SetBranchAddress("nlong", &nlong, &b_nlong);
fChain->SetBranchAddress("nvelo", &nvelo, &b_nvelo);
fChain->SetBranchAddress("down", &down, &b_ndown);
fChain->SetBranchAddress("ntk", &ntk, &b_ntk);
fChain->SetBranchAddress("nspd", &nspd, &b_nspd);
fChain->SetBranchAddress("ntotphot", &ntotphot, &b_ntotphot);
fChain->SetBranchAddress("nout", &nout, &b_nout);
fChain->SetBranchAddress("nphotons", &nphotons, &b_nphotons);
fChain->SetBranchAddress("etprev1", &etprev1, &b_etprev1);
fChain->SetBranchAddress("etprev2", &etprev2, &b_etprev2);
fChain->SetBranchAddress("tt", tt, &b_tt);
fChain->SetBranchAddress("tx", tx, &b_tx);
fChain->SetBranchAddress("ty", ty, &b_ty);
fChain->SetBranchAddress("tz", tz, &b_tz);
fChain->SetBranchAddress("px", px, &b_px);
fChain->SetBranchAddress("py", py, &b_py);
fChain->SetBranchAddress("pz", pz, &b_pz);
fChain->SetBranchAddress("pt", pt, &b_pt);
fChain->SetBranchAddress("charge", charge, &b_charge);
fChain->SetBranchAddress("muon", muon, &b_muon);
fChain->SetBranchAddress("muacc", muacc, &b_muacc);
fChain->SetBranchAddress("pide", pide, &b_pide);
fChain->SetBranchAddress("pidk", pidk, &b_pidk);
fChain->SetBranchAddress("pidpi", pidpi, &b_pidpi);
fChain->SetBranchAddress("pidp", pidp, &b_pidp);
fChain->SetBranchAddress("pidm", pidm, &b_pidm);
fChain->SetBranchAddress("ecal", ecal, &b_ecal);
fChain->SetBranchAddress("hcal", hcal, &b_hcal);
fChain->SetBranchAddress("eoverp", eoverp, &b_eoverp);
fChain->SetBranchAddress("phx", phx, &b_phx);
fChain->SetBranchAddress("phy", phy, &b_phy);
fChain->SetBranchAddress("phz", phz, &b_phz);
fChain->SetBranchAddress("spd", spd, &b_spd);
fChain->SetBranchAddress("prs", prs, &b_prs);
fChain->SetBranchAddress("nl0", &nl0, &b_nl0);
fChain->SetBranchAddress("l0et", l0et, &b_l0et);
fChain->SetBranchAddress("l0phi", l0phi, &b_l0phi);
fChain->SetBranchAddress("l0theta", l0theta, &b_l0theta);
fChain->SetBranchAddress("nher", &nher, &b_nher);
fChain->SetBranchAddress("herchan", herchan, &b_herchan);
fChain->SetBranchAddress("heradc", heradc, &b_heradc);
Notify();
}
Bool_t mymacro::Notify()
{
// The Notify() function is called when a new file is opened. This
// can be either for a new TTree in a TChain or when when a new TTree
// is started when using PROOF. It is normally not necessary to make changes
// to the generated code, but the routine can be extended by the
// user if needed. The return value is currently not used.
return kTRUE;
}
void mymacro::Show(Long64_t entry)
{
// Print contents of entry.
// If entry is not specified, print current entry
if (!fChain) return;
fChain->Show(entry);
}
Int_t mymacro::Cut(Long64_t entry)
{
// This function may be called from Loop.
// returns 1 if entry is accepted.
// returns -1 otherwise.
return 1;
}
#endif // #ifdef mymacro_cxx
and
#define mymacro_cxx
#include "mymacro.h"
#include <TH2.h>
#include <TStyle.h>
#include <TCanvas.h>
void mymacro::Loop()
{
// In a ROOT session, you can do:
// root> .L mymacro.C
// root> mymacro t
// root> t.GetEntry(12); // Fill t data members with entry number 12
// root> t.Show(); // Show values of entry 12
// root> t.Show(16); // Read and show values of entry 16
// root> t.Loop(); // Loop on all entries
//
// This is the loop skeleton where:
// jentry is the global entry number in the chain
// ientry is the entry number in the current Tree
// Note that the argument to GetEntry must be:
// jentry for TChain::GetEntry
// ientry for TTree::GetEntry and TBranch::GetEntry
//
// To read only selected branches, Insert statements like:
// METHOD1:
// fChain->SetBranchStatus("*",0); // disable all branches
// fChain->SetBranchStatus("branchname",1); // activate branchname
// METHOD2: replace line
// fChain->GetEntry(jentry); //read all branches
//by b_branchname->GetEntry(ientry); //read only this branch
TCanvas *c1= new TCanvas ("c1","c1",800,800);
if (fChain == 0) return;
fChain->SetBranchStatus("*",0); // disable all branches
fChain->SetBranchStatus("px",1); // activate branchname
fChain->SetBranchStatus("py",1); // activate branchname
fChain->SetBranchStatus("pz",1); // activate branchname
fChain->SetBranchStatus("muon",1); // activate branchname
fChain->SetBranchStatus("charge",1); // activate branchname
Long64_t nentries = fChain->GetEntriesFast();
TH1F *histo=new TH1F("histo","histo",100,3000,4000);//book histogram
Long64_t nbytes = 0, nb = 0;
for (Long64_t jentry=0; jentry<nentries;jentry++) {
Long64_t ientry = LoadTree(jentry);
if (ientry < 0) break;
nb = fChain->GetEntry(jentry); nbytes += nb;
// if (Cut(ientry) < 0) continue;
}
}
What I would like to understand is how to get the number of tracks per event, so to be able to insert something like:
double pxtot,pytot,pztot,etot;
for(int i=0;i<ntracks;i++){
pxtot+=px[i];
pytot+=py[i];
pztot+=pz[i];
etot+=sqrt(pow(px[i],2)+pow(py[i].2)+pow(pz[i],2));
}
double mass=sqrt(pow(etot,2)-pow(pxtot,2)-pow(pytot,2)-pow(pztot,2))l;
histo->Fill(mass);
in the jentry loop. Previously, in PyRoot I circumvented the issue by working with something like
for entry in tree:
#access the momenta and other vars
buf_px = np.array(entry.px)
and then looping with
#for every track in the event
for i in range(len(bufpx)):
I would appreciate any help on the matter.
Thanks,
bldelane