Jpp  15.0.1-rc.1-highQE
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JDAQProfile.cc
Go to the documentation of this file.
1 
2 #include <string>
3 #include <iostream>
4 #include <iomanip>
5 
6 #include "TROOT.h"
7 #include "TFile.h"
8 #include "TH1D.h"
9 #include "TProfile.h"
10 
12 #include "JDAQ/JDAQEventIO.hh"
14 #include "JDAQ/JDAQEvaluator.hh"
15 
17 #include "JSupport/JTreeScanner.hh"
19 #include "JSupport/JSupport.hh"
20 
21 #include "JROOT/JManager.hh"
22 
23 #include "Jeep/JParser.hh"
24 #include "Jeep/JMessage.hh"
25 
26 
27 /**
28  * \file
29  *
30  * Example program to histogram various data profiles.
31  * \author mdejong
32  */
33 int main(int argc, char **argv)
34 {
35  using namespace std;
36  using namespace JPP;
37  using namespace KM3NETDAQ;
38 
39  string inputFile;
40  JLimit_t numberOfEvents;
41  string outputFile;
42  Long64_t prescale;
43  int debug;
44 
45  try {
46 
47  JParser<> zap("Example program to histogram various data profiles.");
48 
49  zap['f'] = make_field(inputFile);
50  zap['o'] = make_field(outputFile) = "profile.root";
51  zap['n'] = make_field(numberOfEvents) = JLimit::max();
52  zap['P'] = make_field(prescale) = 1;
53  zap['d'] = make_field(debug) = 1;
54 
55  zap(argc, argv);
56  }
57  catch(const exception& error) {
58  FATAL(error.what() << endl);
59  }
60 
61 
62  cout.tie(&cerr);
63 
64 
65  NOTICE("Determine frame index range... " << flush);
66 
67  JFrameIndexRange frame_index = getFrameIndexRange<JDAQSummaryslice>(inputFile);
68 
69  NOTICE("= (" << frame_index.first << "," << frame_index.second << ")" << endl);
70 
71  const Long64_t NX = frame_index.second - frame_index.first + 1;
72  const Long64_t nx = (NX + prescale - 1) / prescale;
73  const Double_t xmin = (Double_t) frame_index.first - 0.5;
74  const Double_t xmax = (Double_t) frame_index.second + 0.5;
75 
76 
77  TFile out(outputFile.c_str(), "recreate");
78 
79  typedef JManager<int, TH1D> JManager_t;
80 
81  JManager_t map_summary(new TH1D("Summary[%]", NULL, NX, xmin, xmax));
82  JManager_t map_trigger(new TH1D("Trigger[%]", NULL, nx, xmin, xmax));
83  JManager_t map_status (new TH1D("Status[%]", NULL, NUMBER_OF_PMTS, -0.5, NUMBER_OF_PMTS - 0.5));
84 
85  TH1D hn("N", NULL, NX, xmin, xmax);
86  TH1D h0("L0", NULL, NX, xmin, xmax);
87  TH1D h1("WR", NULL, NX, xmin, xmax);
88  TH1D h2("trigger", NULL, nx, xmin, xmax);
89 
90  h2.Sumw2();
91 
92  for (JMultipleFileScanner<JDAQSummaryslice> in(inputFile, numberOfEvents); in.hasNext(); ) {
93 
94  STATUS("event: " << setw(10) << in.getCounter() << '\r'); DEBUG(endl);
95 
96  JDAQSummaryslice* summaryslice = in.next();
97 
98  const Double_t x = summaryslice->getFrameIndex();
99 
100 
101  for (JDAQSummaryslice::const_iterator frame = summaryslice->begin(); frame != summaryslice->end(); ++frame) {
102 
103  if (frame->testHighRateVeto()) {
104 
105  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
106 
107  if (frame->testHighRateVeto(pmt)) {
108  map_status[frame->getModuleID()]->Fill((Double_t) pmt);
109  }
110  }
111  }
112 
113  Double_t y = 0.0;
114 
115  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
116  if (!frame->testHighRateVeto(pmt)) {
117  y += frame->getRate(pmt) * 1e-3; // kHz
118  }
119  }
120 
121  map_summary[frame->getModuleID()]->Fill(x, y);
122  }
123 
124 
125  Double_t y = 0.0;
126 
127  for (JDAQSummaryslice::const_iterator frame = summaryslice->begin(); frame != summaryslice->end(); ++frame) {
128 
129  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
130  y += frame->getRate(pmt) * 1e-3; // kHz
131  }
132 
133  h1.Fill(x, (frame->testWhiteRabbitStatus() ? 1.0 : 0.0));
134  }
135 
136  hn.Fill(x, summaryslice->size());
137  h0.Fill(x, y);
138  }
139  STATUS(endl);
140 
141 
142  for (JMultipleFileScanner<JDAQEvent> in(inputFile, numberOfEvents); in.hasNext(); ) {
143 
144  STATUS("event: " << setw(10) << in.getCounter() << '\r'); DEBUG(endl);
145 
146  JDAQEvent* event = in.next();
147 
148  const Double_t x = event->getFrameIndex();
149  const Double_t y = 1.0e9 / getFrameTime() / prescale;
150 
151  h2.Fill(x, y);
152 
153  typedef JDAQTriggeredHit JHit_t;
154  //typedef JDAQSnapshotHit JHit_t;
155 
156  for (JDAQEvent::const_iterator<JHit_t> hit = event->begin<JHit_t>(); hit != event->end<JHit_t>(); ++hit) {
157  map_trigger[hit->getModuleID()]->Fill(x);
158  }
159  }
160  STATUS(endl);
161 
162  map_summary.Write(out);
163  map_trigger.Write(out);
164  map_status .Write(out);
165 
166  out.Write();
167  out.Close();
168 }
Auxiliary class to set-up Hit.
Definition: JHit_t.hh:25
Utility class to parse command line options.
Definition: JParser.hh:1500
int main(int argc, char *argv[])
Definition: Main.cc:15
ROOT TTree parameter settings of various packages.
#define STATUS(A)
Definition: JMessage.hh:63
Template const_iterator.
Definition: JDAQEvent.hh:68
then for HISTOGRAM in h0 h1
Definition: JMatrixNZ.sh:71
Dynamic ROOT object management.
string outputFile
int getFrameIndex() const
Get frame index.
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:41
Auxiliary class to manage set of compatible ROOT objects (e.g. histograms) using unique keys...
Definition: JManager.hh:43
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1961
double getFrameTime()
Get frame time duration.
Definition: JDAQClock.hh:162
#define NOTICE(A)
Definition: JMessage.hh:64
Support methods.
int debug
debug level
Definition: JSirene.cc:63
General purpose messaging.
#define FATAL(A)
Definition: JMessage.hh:67
Scanning of objects from multiple files according a format that follows from the extension of each fi...
Utility class to parse command line options.
then usage $script< input_file >< detector_file > fi set_variable OUTPUT_DIR set_variable SELECTOR JDAQTimesliceL1 set_variable DEBUG case set_variable DEBUG
static const int NUMBER_OF_PMTS
Total number of PMTs in module.
Definition: JDAQ.hh:26
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable INPUT_FILE $argv[2] eval JPrintDetector a $DETECTOR O IDENTIFIER eval JPrintDetector a $DETECTOR O SUMMARY source JAcoustics sh $DETECTOR_ID CHECK_EXIT_CODE typeset A TRIPODS get_tripods $WORKDIR tripod txt TRIPODS for EMITTER in
Definition: JCanberra.sh:41