Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
mp104_processH1.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_multicore
3/// \notebook -nodraw
4/// Illustrate the usage of the multiproc to process the H1 analysis
5/// example.
6///
7/// \macro_code
8///
9/// \author Gerardo Ganis
10
11#include "TString.h"
12#include "TROOT.h"
13#include "TTree.h"
14#include "TH1F.h"
15#include "TH2F.h"
16#include "TEntryList.h"
17#include "TTreeReader.h"
18#include "TTreeReaderArray.h"
19#include "TTreeReaderValue.h"
20#include "TSystem.h"
21#include "TMath.h"
22#include "TCanvas.h"
23#include "TStyle.h"
24#include "TF1.h"
25#include "TLine.h"
26#include "TPaveStats.h"
27#include "TStopwatch.h"
29
30static std::string tutname = "mp104_processH1: ";
31static std::string logfile = "mp104_processH1.log";
32static RedirectHandle_t gRH;
33
34std::vector<std::string> files {"http://root.cern/files/h1/dstarmb.root",
35 "http://root.cern/files/h1/dstarp1a.root",
36 "http://root.cern/files/h1/dstarp1b.root",
37 "http://root.cern/files/h1/dstarp2.root"};
38
39int mp104_processH1()
40{
41
42 // MacOSX may generate connection to WindowServer errors
43 gROOT->SetBatch(kTRUE);
44
45 TStopwatch stp;
46
47 // Check and fit lambdas
48 #include "mp_H1_lambdas.C"
49
51
52 std::cout << tutname << "processing the H1 dataset with a lambda \n";
53
54 auto hListFun = pool.Process(files, doH1, "h42");
55
56 // Check the output
57 if (checkH1(hListFun) < 0)
58 return -1;
59
60 // Do the fit
61 if (doFit(hListFun, logfile.c_str()) < 0)
62 return -1;
63
64 stp.Print();
65 stp.Start();
66
67 // Run the analysis with a selector
68
69 TString selectorPath = gROOT->GetTutorialDir();
70 selectorPath += "/tree/h1analysisTreeReader.C+";
71 std::cout << tutname << "processing the H1 dataset with selector '" << selectorPath << "'\n";
72 auto sel = TSelector::GetSelector(selectorPath);
73
74 // In a second run we use sel
75 gSystem->RedirectOutput(logfile.c_str(), "w", &gRH);
76 auto hListSel = pool.Process(files, *sel, "h42");
77 gSystem->RedirectOutput(nullptr, nullptr, &gRH);
78
79 // Check the output
80 if (checkH1(hListSel) < 0)
81 return -1;
82
83 // Do the fit
84 if (doFit(hListSel, logfile.c_str()) < 0)
85 return -1;
86
87 stp.Print();
88 stp.Start();
89
90 return 0;
91}
constexpr Bool_t kTRUE
Definition RtypesCore.h:93
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h Atom_t Int_t ULong_t ULong_t unsigned char prop_list Atom_t sel
#define gROOT
Definition TROOT.h:406
R__EXTERN TSystem * gSystem
Definition TSystem.h:561
This class provides an interface to process a TTree dataset in parallel with multi-process technology...
static TSelector * GetSelector(const char *filename)
The code in filename is loaded (interpreted or compiled, see below), filename must contain a valid cl...
Stopwatch class.
Definition TStopwatch.h:28
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
void Print(Option_t *option="") const override
Print the real and cpu time passed between the start and stop events.
Basic string class.
Definition TString.h:139
virtual Int_t RedirectOutput(const char *name, const char *mode="a", RedirectHandle_t *h=nullptr)
Redirect standard output (stdout, stderr) to the specified file.
Definition TSystem.cxx:1715
<a href="https://nbviewer.jupyter.org/url/root.cern/doc/master/notebooks/mp_H1_lambdas....