Jpp 20.0.0-27-g39925593c-D
the software that should make you happy
Loading...
Searching...
No Matches
Functions
JMergeRunAnalyzer.cc File Reference
#include <string>
#include <iostream>
#include <iomanip>
#include <vector>
#include <map>
#include "TROOT.h"
#include "TFile.h"
#include "TH1D.h"
#include "TH2D.h"
#include "TKey.h"
#include "TString.h"
#include "TMath.h"
#include "TTimeStamp.h"
#include "km3net-dataformat/online/JDAQHeader.hh"
#include "JROOT/JRootFileReader.hh"
#include "JROOT/JRootFileWriter.hh"
#include "JSupport/JFileRecorder.hh"
#include "JSupport/JMultipleFileScanner.hh"
#include "JSupport/JSingleFileScanner.hh"
#include "JSupport/JTreeScanner.hh"
#include "JSupport/JMeta.hh"
#include "JTrigger/JTriggerInterface.hh"
#include "JTrigger/JTriggerBits.hh"
#include "JGizmo/JGizmoToolkit.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"
#include "JRunAnalyzer/JRunInfo.hh"

Go to the source code of this file.

Functions

void fill2DHistogram (int run, TH1D *h1, TH2D *h2, bool w_mc, double daq_livetime=1, double mc_livetime=1)
 
int main (int argc, char **argv)
 

Function Documentation

◆ fill2DHistogram()

void fill2DHistogram ( int run,
TH1D * h1,
TH2D * h2,
bool w_mc,
double daq_livetime = 1,
double mc_livetime = 1 )
inline

Definition at line 45 of file JMergeRunAnalyzer.cc.

45 {
46
47 if(h1 != NULL && h2 != NULL){
48
49 double w = 1;
50
51 for(int i = 1; i <= h1->GetNbinsX(); ++i){
52
53 if(w_mc){
54 w = daq_livetime / mc_livetime;
55 }
56
57 h2->Fill(run, h1->GetXaxis()->GetBinCenter(i), h1->GetBinContent(i) * w);
58 }
59 }
60}
JOSCPROB::JOscProbInterpolator
Template definition of a multi-dimensional oscillation probability interpolation table.
Definition JOscProbInterpolator.hh:71

◆ main()

int main ( int argc,
char ** argv )

Definition at line 63 of file JMergeRunAnalyzer.cc.

64{
65 using namespace std;
66 using namespace JPP;
67 using namespace KM3NETDAQ;
68
69 vector<string> inputFile;
70 JFileRecorder<JTYPELIST<run_info, JMeta, TH1D, TH2D>::typelist> outputFile;
71 bool weight_mc;
72 int debug;
73
74 try {
75
76 JParser<> zap("Auxiliary program to merge JRunAnalyzer histograms.");
77
78 zap['f'] = make_field(inputFile, "input file (output from JRunAnalyzer).");
79 zap['o'] = make_field(outputFile, "output file.") = "merge-jra.root";
80 zap['w'] = make_field(weight_mc, "weight mc events");
81 zap['d'] = make_field(debug, "debug.") = 1;
82
83 zap(argc, argv);
84 }
85 catch(const exception &error) {
86 FATAL(error.what() << endl);
87 }
88
89 int min_run = numeric_limits<int>::max(), max_run = 0;
90
91 size_t min_date = numeric_limits<int>::max(), max_date = 0;
92
93 run_info runinfo;
94
95 double total_mc_livetime = 0;
96 int run = 0;
97
98 for (vector<string>::const_iterator i = inputFile.begin(); i != inputFile.end(); ++i) {
99
100 {
101 JSingleFileScanner<run_info> in(*i);
102 if(in.hasNext()){
103 runinfo = *in.next();
104 if(in.hasNext()){
105 FATAL("MULTIPLE RUN INFO");
106 }
107
108 int run_id = runinfo.run_id;
109
110 if(run_id != run){
111
112 run = run_id;
113
114 total_mc_livetime = runinfo.mc_livetime;
115
116 } else {
117
118 total_mc_livetime += runinfo.mc_livetime;
119 }
120
121 } else {
122 FATAL("NO RUN INFO");
123 }
124 }
125
126 for (JRootFileReader<JDAQHeader> in(i->c_str()); in.hasNext(); ) {
127
128 const JDAQHeader* p = in.next();
129
130 if(p->getTimesliceStart().getUTCseconds() < min_date) min_date = p->getTimesliceStart().getUTCseconds();
131 if(p->getTimesliceStart().getUTCseconds() > max_date) max_date = p->getTimesliceStart().getUTCseconds();
132 if(p->getRunNumber() < min_run) min_run = p->getRunNumber();
133 if(p->getRunNumber() > max_run) max_run = p->getRunNumber();
134 }
135 }
136
137 time_t tmin = min_date, tmax = max_date;
138
139 TTimeStamp date_min(tmin, 0);
140 TTimeStamp date_max(tmax, 0);
141
142 cout << "Minimum RUN: " << min_run << ", "; date_min.Print("s");
143 cout << "Maximum RUN: " << max_run << ", "; date_max.Print("s");
144
145 int xbins = (max_run - min_run) + 1;
146 string title = string("Data taken from: ") + date_min.AsString("s") + " - " + date_max.AsString("s");
147
148 TH2D* h_pmt_rate_distribution = new TH2D("h_pmt_rate_distribution", string(title + "; RUN; PMT Rate Distribution [kHz]").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, JDAQRate::getN(), JDAQRate::getData(1.0e-4));
149 TH2D* h_Trigger_bit_hit = new TH2D("h_Trigger_bit_hit", string(title + "; RUN; Hits trigger bit").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, NUMBER_OF_TRIGGER_BITS, -0.5, NUMBER_OF_TRIGGER_BITS - 0.5);
150 TH2D* h_Triggered_hits = new TH2D("h_Triggered_hits", string(title + "; RUN; Number of triggered hits").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 50, 0, 4);
151 setLogarithmicY(h_Triggered_hits);
152 TH2D* h_Triggered_hits_3dmuon = new TH2D("h_Triggered_hits_3dmuon", string(title + "; RUN; Number of triggered hits - JTRIGGER3DMUON").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 50, 0, 3);
153 setLogarithmicY(h_Triggered_hits_3dmuon);
154 TH2D* h_Snapshot_hits = new TH2D("h_Snapshot_hits", string(title + "; RUN; Number of snapshot hits").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 50, 0, 4);
155 setLogarithmicY(h_Snapshot_hits);
156 TH2D* h_Triggered_over_Snapshot_hits = new TH2D("h_Triggered_over_Snapshot_hits", string(title + "; RUN; Triggered/Snapshot hits").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 100, 0, 0.5);
157 TH2D* h_Number_of_overlays = new TH2D("h_Number_of_overlays", string(title + ";RUN;Number of Overlays").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 1000, -0.5, 1000 - 0.5);
158 TH2D* h_pmt_distribution_triggered_hits = new TH2D("h_pmt_distribution_triggered_hits", string(title + ";RUN;Normalised PMT Distrib triggered hits (upper PMTs=[0-11])").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, NUMBER_OF_PMTS, -0.5 , NUMBER_OF_PMTS - 0.5);
159 TH2D* h_pmt_distribution_snapshot_hits = new TH2D("h_pmt_distribution_snapshot_hits", string(title + ";RUN; Normalised PMT Distrib snapshot hits (upper PMTs=[0-11])").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, NUMBER_OF_PMTS, -0.5 , NUMBER_OF_PMTS - 0.5);
160 TH2D* h_event_duration = new TH2D("h_event_duration", string(title + ";RUN;Event Duration [ns]").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 60, 1, 6);
161 setLogarithmicY(h_event_duration);
162 TH2D* h_tres = new TH2D("h_tres", string(title + ";RUN;Time residuals [ns]").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 100, -50, 150);
163 TH2D* h_likelihood = new TH2D ("h_likelihood", string(title + ";RUN;Likelihood").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 100, -1000, 1000);
164 TH2D* h_beta0 = new TH2D("h_beta0", string(title + ";RUN;beta0").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 20, 0, 1);
165 TH2D* h_energy = new TH2D ("h_energy", string(title + ";RUN;Energy [GeV];").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 65, 0, 9);
166 setLogarithmicY(h_energy);
167 TH2D* h_zenith = new TH2D("h_zenith", string(title + ";RUN;cosZenith").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 20, -1, 1);
168 TH2D* h_azimuth = new TH2D("h_azimuth", string(title + ";RUN;cosAzimuth").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 20, -1, 1);
169 TH2D* h_radial_position = new TH2D ("h_radial_position", string(title + ";RUN;Radial Position [m]").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 60, 0, 4.2);
170 setLogarithmicY(h_radial_position);
171 TH2D* h_z_position = new TH2D ("h_z_position", string(title + ";RUN;Z Position [m]").c_str(), xbins, (double) min_run-0.5, (double) max_run+0.5, 50, 0, 3.2);
172 setLogarithmicY(h_z_position);
173
174 for (vector<string>::const_iterator i = inputFile.begin(); i != inputFile.end(); ++i) {
175
176 DEBUG("Processing " << *i << endl) ;
177
178 {
179 JSingleFileScanner<run_info> in(*i);
180 if(in.hasNext()){
181 runinfo = *in.next();
182 if(in.hasNext()){
183 FATAL("MULTIPLE RUN INFO");
184 }
185 } else {
186 FATAL("NO RUN INFO");
187 }
188 }
189
190 TFile in(i->c_str(), "read");
191
192 cout << i->c_str() << endl;
193
194 int run = runinfo.run_id;
195
196 TH1D* h1d_prd = new TH1D();
197 TH1D* h1d_thb = new TH1D();
198 TH1D* h1d_th = new TH1D();
199 TH1D* h1d_th3D = new TH1D();
200 TH1D* h1d_sh = new TH1D();
201 TH1D* h1d_tosh = new TH1D();
202 TH1D* h1d_noo = new TH1D();
203 TH1D* h1d_pdth = new TH1D();
204 TH1D* h1d_pdsh = new TH1D();
205 TH1D* h1d_ed = new TH1D();
206 TH1D* h1d_tr = new TH1D();
207 TH1D* h1d_lk = new TH1D();
208 TH1D* h1d_b0 = new TH1D();
209 TH1D* h1d_en = new TH1D();
210 TH1D* h1d_zn = new TH1D();
211 TH1D* h1d_az = new TH1D();
212 TH1D* h1d_rp = new TH1D();
213 TH1D* h1d_zp = new TH1D();
214
215 in.GetDirectory("Detector")->GetObject("h_pmt_rate_distribution", h1d_prd);
216 fill2DHistogram(run, h1d_prd, h_pmt_rate_distribution, weight_mc, runinfo.daq_livetime, total_mc_livetime);
217
218 TDirectory* dir = in.GetDirectory("JDAQEvent");
219
220 dir->GetObject("h_Trigger_bit_hit", h1d_thb);
221 fill2DHistogram(run, h1d_thb, h_Trigger_bit_hit, weight_mc, runinfo.daq_livetime, total_mc_livetime);
222
223 dir->GetObject("h_Triggered_hits", h1d_th);
224 fill2DHistogram(run, h1d_th, h_Triggered_hits, weight_mc, runinfo.daq_livetime, total_mc_livetime);
225
226 dir->GetObject("h_Triggered_hits_3dmuon", h1d_th3D);
227 fill2DHistogram(run, h1d_th3D, h_Triggered_hits_3dmuon, weight_mc, runinfo.daq_livetime, total_mc_livetime);
228
229 dir->GetObject("h_Snapshot_hits", h1d_sh);
230 fill2DHistogram(run, h1d_sh, h_Snapshot_hits, weight_mc, runinfo.daq_livetime, total_mc_livetime);
231
232 dir->GetObject("h_Triggered_over_Snapshot_hits", h1d_tosh);
233 fill2DHistogram(run, h1d_tosh, h_Triggered_over_Snapshot_hits, weight_mc, runinfo.daq_livetime, total_mc_livetime);
234
235 dir->GetObject("h_Number_of_overlays", h1d_noo);
236 fill2DHistogram(run, h1d_noo, h_Number_of_overlays, weight_mc, runinfo.daq_livetime, total_mc_livetime);
237
238 dir->GetObject("h_pmt_distribution_triggered_hits", h1d_pdth);
239 fill2DHistogram(run, h1d_pdth, h_pmt_distribution_triggered_hits, weight_mc, runinfo.daq_livetime, total_mc_livetime);
240
241 dir->GetObject("h_pmt_distribution_snapshot_hits", h1d_pdsh);
242 fill2DHistogram(run, h1d_pdsh, h_pmt_distribution_snapshot_hits, weight_mc, runinfo.daq_livetime, total_mc_livetime);
243
244 dir->GetObject("h_event_duration", h1d_ed);
245 fill2DHistogram(run, h1d_ed, h_event_duration, weight_mc, runinfo.daq_livetime, total_mc_livetime);
246
247 if (in.GetKey("Reco")) {
248 TDirectory* dir_r = in.GetDirectory("Reco");
249
250 dir_r->GetObject("h_tres", h1d_tr);
251 fill2DHistogram(run, h1d_tr, h_tres, weight_mc, runinfo.daq_livetime, total_mc_livetime);
252
253 dir_r->GetObject("h_likelihood", h1d_lk);
254 fill2DHistogram(run, h1d_lk, h_likelihood, weight_mc, runinfo.daq_livetime, total_mc_livetime);
255
256 dir_r->GetObject("h_beta0", h1d_b0);
257 fill2DHistogram(run, h1d_b0, h_beta0, weight_mc, runinfo.daq_livetime, total_mc_livetime);
258
259 dir_r->GetObject("h_energy", h1d_en);
260 fill2DHistogram(run, h1d_en, h_energy, weight_mc, runinfo.daq_livetime, total_mc_livetime);
261
262 dir_r->GetObject("h_zenith", h1d_zn);
263 fill2DHistogram(run, h1d_zn, h_zenith, weight_mc, runinfo.daq_livetime, total_mc_livetime);
264
265 dir_r->GetObject("h_azimuth", h1d_az);
266 fill2DHistogram(run, h1d_az, h_azimuth, weight_mc, runinfo.daq_livetime, total_mc_livetime);
267
268 dir_r->GetObject("h_radial_position", h1d_rp);
269 fill2DHistogram(run, h1d_rp, h_radial_position, weight_mc, runinfo.daq_livetime, total_mc_livetime);
270
271 dir_r->GetObject("h_z_position", h1d_zp);
272 fill2DHistogram(run, h1d_zp, h_z_position, weight_mc, runinfo.daq_livetime, total_mc_livetime);
273 }
274 }
275
276 outputFile.open();
277
278 outputFile.put(*h_event_duration);
279 outputFile.put(*h_pmt_distribution_snapshot_hits);
280 outputFile.put(*h_pmt_distribution_triggered_hits);
281 outputFile.put(*h_Number_of_overlays);
282 outputFile.put(*h_Triggered_over_Snapshot_hits);
283 outputFile.put(*h_Snapshot_hits);
284 outputFile.put(*h_Triggered_hits_3dmuon);
285 outputFile.put(*h_Triggered_hits);
286 outputFile.put(*h_Trigger_bit_hit);
287 outputFile.put(*h_tres);
288 outputFile.put(*h_likelihood);
289 outputFile.put(*h_beta0);
290 outputFile.put(*h_energy);
291 outputFile.put(*h_zenith);
292 outputFile.put(*h_azimuth);
293 outputFile.put(*h_radial_position);
294 outputFile.put(*h_z_position);
295
296 JMultipleFileScanner<run_info> io(inputFile);
297
298 io >> outputFile;
299
300 outputFile.close();
301}
outputFile
string outputFile
Definition JDAQTimesliceSelector.cc:37
fill2DHistogram
void fill2DHistogram(int run, TH1D *h1, TH2D *h2, bool w_mc, double daq_livetime=1, double mc_livetime=1)
Definition JMergeRunAnalyzer.cc:45
DEBUG
#define DEBUG(A)
Message macros.
Definition JMessage.hh:62
FATAL
#define FATAL(A)
Definition JMessage.hh:67
debug
int debug
debug level
Definition JSirene.cc:72
make_field
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition JParser.hh:2142
KM3NETDAQ::JDAQChronometer::getRunNumber
int getRunNumber() const
Get run number.
Definition JDAQChronometer.hh:121
KM3NETDAQ::JDAQChronometer::getTimesliceStart
JDAQUTCExtended getTimesliceStart() const
Get start of timeslice.
Definition JDAQChronometer.hh:144
KM3NETDAQ::JDAQRate::getN
static int getN()
Get number of bins.
Definition JDAQSummaryFrame.hh:138
KM3NETDAQ::JDAQRate::getData
static const double * getData(const double factor=1.0)
Get abscissa values.
Definition JDAQSummaryFrame.hh:150
KM3NETDAQ::JDAQUTCExtended::getUTCseconds
uint32_t getUTCseconds() const
Get major time.
Definition JDAQUTCExtended.hh:87
JGIZMO::setLogarithmicY
void setLogarithmicY(TList *list)
Make y-axis of objects in list logarithmic (e.g. after using log10()).
Definition JGizmoToolkit.hh:1027
JPP
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Definition JAAnetToolkit.hh:45
KM3NETDAQ
KM3NeT DAQ data structures and auxiliaries.
Definition DataQueue.cc:39
run_info
Definition JRunInfo.hh:10
Generated by doxygen 1.10.0