Logo ROOT  
Reference Guide
mt102_readNtuplesFillHistosAndFit.C File Reference

Detailed Description

View in nbviewer Open in SWAN Read n-tuples in distinct workers, fill histograms, merge them and fit. Knowing that other facilities like TProcessExecutor might be more adequate for this operation, this tutorial complements mc101, reading and merging. We convey another message with this tutorial: the synergy of ROOT and STL algorithms is possible.

(int) 0
Int_t mt102_readNtuplesFillHistosAndFit()
{
// No nuisance for batch execution
gROOT->SetBatch();
// Perform the operation sequentially ---------------------------------------
TChain inputChain("multiCore");
inputChain.Add("mt101_multiCore_*.root");
TH1F outHisto("outHisto", "Random Numbers", 128, -4, 4);
inputChain.Draw("r >> outHisto");
outHisto.Fit("gaus");
// We now go MT! ------------------------------------------------------------
// The first, fundamental operation to be performed in order to make ROOT
// thread-aware.
// We adapt our parallelisation to the number of input files
const auto nFiles = inputChain.GetListOfFiles()->GetEntries();
// We define the histograms we'll fill
std::vector<TH1F> histograms;
auto workerIDs = ROOT::TSeqI(nFiles);
histograms.reserve(nFiles);
for (auto workerID : workerIDs) {
histograms.emplace_back(TH1F(Form("outHisto_%u", workerID), "Random Numbers", 128, -4, 4));
}
// We define our work item
auto workItem = [&histograms](UInt_t workerID) {
TFile f(Form("mt101_multiCore_%u.root", workerID));
auto ntuple = f.Get<TNtuple>("multiCore");
auto &histo = histograms.at(workerID);
for (auto index : ROOT::TSeqL(ntuple->GetEntriesFast())) {
ntuple->GetEntry(index);
histo.Fill(ntuple->GetArgs()[0]);
}
};
TH1F sumHistogram("SumHisto", "Random Numbers", 128, -4, 4);
// Create the collection which will hold the threads, our "pool"
std::vector<std::thread> workers;
// Spawn workers
// Fill the "pool" with workers
for (auto workerID : workerIDs) {
workers.emplace_back(workItem, workerID);
}
// Now join them
for (auto &&worker : workers)
worker.join();
// And reduce with a simple lambda
std::for_each(std::begin(histograms), std::end(histograms),
[&sumHistogram](const TH1F &h) { sumHistogram.Add(&h); });
sumHistogram.Fit("gaus", 0);
return 0;
}
Date
January 2016
Author
Danilo Piparo

Definition in file mt102_readNtuplesFillHistosAndFit.C.

f
#define f(i)
Definition: RSha256.hxx:122
TNtuple
Definition: TNtuple.h:28
ROOT::EnableThreadSafety
void EnableThreadSafety()
Enables the global mutex to make ROOT thread safe/aware.
Definition: TROOT.cxx:494
Form
char * Form(const char *fmt,...)
ROOT::TSeqI
TSeq< int > TSeqI
Definition: TSeq.hxx:194
Int_t
int Int_t
Definition: RtypesCore.h:45
TTree::GetEntriesFast
virtual Long64_t GetEntriesFast() const
Definition: TTree.h:460
TTree::GetEntry
virtual Int_t GetEntry(Long64_t entry=0, Int_t getall=0)
Read all branches of entry and return total number of bytes read.
Definition: TTree.cxx:5537
TNtuple::GetArgs
Float_t * GetArgs() const
Definition: TNtuple.h:56
h
#define h(i)
Definition: RSha256.hxx:124
UInt_t
unsigned int UInt_t
Definition: RtypesCore.h:46
TFile
Definition: TFile.h:54
TH1F
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:572
TChain
Definition: TChain.h:33
ROOT::TSeq
A pseudo container class which is a generator of indices.
Definition: TSeq.hxx:66
TNtuple::Fill
virtual Int_t Fill()
Fill a Ntuple with current values in fArgs.
Definition: TNtuple.cxx:169
gROOT
#define gROOT
Definition: TROOT.h:406