Jpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JMonitorL1dt.cc
Go to the documentation of this file.
6 #include "JDAQ/JDAQEvaluator.hh"
7 #include "JDetector/JDetector.hh"
10 #include "JTrigger/JHitL1.hh"
11 #include "JTrigger/JBuildL2.hh"
12 #include "JTrigger/JPMTSelector.hh"
14 #include "JSupport/JTreeScanner.hh"
15 #include "JSupport/JSupport.hh"
18 #include "Jeep/JParser.hh"
19 
20 #include "TFile.h"
21 #include "TH1D.h"
22 #include "TH2D.h"
23 #include "TMath.h"
24 
25 #include <string>
26 #include <iostream>
27 
28 
29 namespace MONITORL1DT {
30 
31  /**
32  * Data structure for hit time and DOM identifier.
33  */
34  class JElement {
35  public:
37  id(0),
38  t (0.0)
39  {}
40 
41  JElement(const int __id,
42  const double __t) :
43  id(__id),
44  t (__t)
45  {}
46 
47  int id;
48  double t;
49  };
50 
51 
52  /**
53  * Sort JElement in time.
54  *
55  * \param first first element
56  * \param second second element
57  * \return true if second later than first; else false
58  */
59  inline bool operator<(const JElement& first, const JElement& second)
60  {
61  return first.t < second.t;
62  }
63 
64 
65  /**
66  * Auxiliary data structure for histogram management.
67  */
68  struct JHistogram {
69  /**
70  * Default constructor.
71  */
73  h2s(NULL),
74  h1l(NULL)
75  {}
76 
77 
78  /**
79  * Constructor.
80  *
81  * \param __h2s 2D histogram for signal
82  * \param __h1l 1D histogram for background
83  */
84  JHistogram(TH2D* __h2s,
85  TH1D* __h1l) :
86  h2s(__h2s),
87  h1l(__h1l)
88  {
89  h2s->Sumw2();
90  h1l->Sumw2();
91  }
92 
93  TH2D* h2s;
94  TH1D* h1l;
95  };
96 }
97 
98 /**
99  * \file
100  *
101  * Auxiliary program to plot L1 time difference correlations between DOMs.
102  *
103  * \author kmelis, lnauta
104  */
105 int main(int argc, char **argv)
106 {
107  using namespace std;
108  using namespace MONITORL1DT;
109  using namespace JPP;
110 
112  string outputFile;
113  string detectorFile;
114  Double_t Timewindow_ns;
115  double TmaxL1_ns;
116  JROOTClassSelector selector;
117  unsigned int multiplicity;
118  double livetime_s;
119  int debug;
120 
121  try {
122 
123  JParser<> zap("Program to create L1 hit time difference histograms from raw data.");
124 
125  zap['f'] = make_field(inputFile, "input file");
126  zap['o'] = make_field(outputFile, "output file") = "monitor.root";
127  zap['a'] = make_field(detectorFile, "detector file");
128  zap['t'] = make_field(TmaxL1_ns, "max time between L1 hits [ns]") = 1000.0;
129  zap['T'] = make_field(Timewindow_ns, "time window around t=0 [ns]") = 2400.0;
130  zap['C'] = make_field(selector, "datastream selector") = getROOTClassSelection<JDAQTimesliceTypes_t>();
131  zap['m'] = make_field(multiplicity, "minimal multiplicity of the L1 hits") = 2;
132  zap['L'] = make_field(livetime_s, "livetime of the data, set to positive value") = -1.0; // set to positive value if you want to set the livetime manually (i.e. MC)
133  zap['d'] = make_field(debug) = 1;
134 
135  if (zap.read(argc, argv) != 0)
136  return 1;
137  }
138  catch(const exception &error) {
139  FATAL(error.what() << endl);
140  }
141 
142  cout.tie(&cerr);
143 
145 
146  try {
147  load(detectorFile, detector);
148  }
149  catch(const JException& error) {
150  FATAL(error);
151  }
152 
153  if (detector.empty()) {
154  FATAL("Empty detector." << endl);
155  }
156 
157  const JModuleRouter router(detector);
158 
159  typedef vector<JHitL1> JFrameL1_t;
160 
161  const double ctMin = -1; // full angle over the DOM
162  const JBuildL2<JHitL1> buildL2(JL2Parameters(multiplicity, TmaxL1_ns, ctMin));
163 
164  vector<JHistogram> zmap(detector.size());
165 
166  const int nx = detector.size();
167  const double xmin = -0.5;
168  const double xmax = nx - 0.5;
169 
170  const double ymin = -floor(Timewindow_ns) + 0.5;
171  const double ymax = +floor(Timewindow_ns) + 0.5; // changed - to + to get an even number for Nbins, rebin needs numbers that have divisors.
172  const int ny = (int) ((ymax - ymin) / 1.0);
173 
174  NOTICE("Running JMonitorL1dt: Monitoring L1 time differences MonitorL1dt creating histograms.");
175  for (JDetector::iterator module = detector.begin(); module != detector.end(); ++module) {
176 
177  const JModuleAddress& address = router.getAddress(module->getID());
178 
179  STATUS("Booking histograms for module " << module->getID() << endl);
180 
181  const JString title(module->getID());
182  JString titleString1D;
183  JString titleString2D;
184  titleString1D = title + ".1L";
185  titleString2D = title + ".2S";
186 
187  zmap[address.first] = JHistogram(new TH2D((titleString2D).c_str(), NULL, nx, xmin, xmax, ny, ymin, ymax),
188  new TH1D((titleString1D).c_str(), NULL, nx, xmin, xmax));
189 
190  for (JDetector::iterator mod = detector.begin(); mod != detector.end(); ++mod) {
191  zmap[address.first].h2s->GetXaxis()->SetBinLabel(distance(detector.begin(), mod)+1, Form("%i", mod->getID()));
192  zmap[address.first].h1l->GetXaxis()->SetBinLabel(distance(detector.begin(), mod)+1, Form("%i", mod->getID()));
193  }
194  }
195 
196 
197  JTreeScanner<JDAQSummaryslice, JDAQEvaluator> summaryFile(inputFile);
198 
200 
201  int counter = 0;
202 
203  for ( ; in.hasNext() && counter != inputFile.getLimit(); ++counter) {
204 
205  STATUS("event: " << setw(10) << counter << '\r'); DEBUG(endl);
206 
207  const JDAQTimeslice* timeslice = in.next();
208 
209  JFrameL1_t frameL1;
211  vector<bool> DOM_OK(detector.size(), true);
212 
213  for (JDAQTimeslice::const_iterator super_frame = timeslice->begin(); super_frame != timeslice->end(); ++super_frame) {
214  if (router.hasModule(super_frame->getModuleID()) && !super_frame->empty()) {
215 
216  const JModuleAddress& address = router.getAddress(super_frame->getModuleID());
217  const JModule& module = detector.getModule(address);
218 
219  frameL1.clear();
220 
221  buildL2(*super_frame, module, back_inserter(frameL1));
222 
223  for (JFrameL1_t::iterator L1hit = frameL1.begin(); L1hit != frameL1.end(); ++L1hit) {
224  buffer.push_back(JElement(address.first, L1hit->begin()->getT()));
225  }
226  }
227  }
228 
229  for (vector<JHistogram>::iterator h1 = zmap.begin(); h1 != zmap.end(); ++h1) {
230  if (!DOM_OK[distance(zmap.begin(), h1)]) {
231  continue;
232  }
233  for (unsigned int i = 0; i < detector.size(); ++i) {
234  if (DOM_OK[i]) { h1->h1l->Fill(i, getFrameTime() * 1e-9); } // fill the th1 with the frametime (how long did it measure) in seconds, all weight will give total measured time
235  }
236  }
237 
238  // fill histograms with correlations
239  sort(buffer.begin(), buffer.end());
240 
241  for (vector<JElement>::const_iterator p = buffer.begin(); p != buffer.end(); ) {
243 
244  while (++q != buffer.end() && q->t - p->t <= Timewindow_ns ) {
245  // if (p->id == q->id) { continue; }
246 
247  zmap[p->id].h2s->Fill(q->id, q->t - p->t);
248  zmap[q->id].h2s->Fill(p->id, p->t - q->t);
249  }
250  p++; // move over to the next L1
251  }
252  buffer.clear();
253  }
254  STATUS(endl);
255 
256  // livetime in case set manually
257  if (livetime_s > 0.0) {
258  for (vector<JHistogram>::iterator i = zmap.begin(); i != zmap.end(); ++i) {
259  TH1D* hl = i->h1l;
260  for (int ibin = 1; ibin <= hl->GetNbinsX(); ++ibin) {
261  hl->SetBinContent(ibin, livetime_s);
262  hl->SetBinError(ibin, 0.0000001);
263  }
264  }
265  }
266 
267  TFile out(outputFile.c_str(), "recreate");
268 
269  for (vector<JHistogram>::iterator i = zmap.begin(); i != zmap.end(); ++i) {
270  i->h2s->Write();
271  i->h1l->Write();
272  }
273 
274  out.Write();
275  out.Close();
276 }
Utility class to parse command line options.
Definition: JParser.hh:1493
General exception.
Definition: JException.hh:23
ROOT TTree parameter settings.
Wrapper class around STL string class.
Definition: JString.hh:27
Data structure for a composite optical module.
Definition: JModule.hh:50
std::vector< T >::difference_type distance(typename std::vector< T >::const_iterator first, typename PhysicsEvent::const_iterator< T > second)
Specialisation of STL distance.
bool operator<(const JElement &first, const JElement &second)
Sort JElement in time.
Definition: JMonitorL1dt.cc:59
#define STATUS(A)
Definition: JMessage.hh:63
Detector data structure.
Definition: JDetector.hh:80
Auxiliary class to select ROOT class based on class name.
Router for direct addressing of module data in detector data structure.
JHistogram(TH2D *__h2s, TH1D *__h1l)
Constructor.
Definition: JMonitorL1dt.cc:84
esac print_variable DETECTOR INPUT_FILE OUTPUT_FILE CDF for TYPE in
Definition: JSirene.sh:45
then for HISTOGRAM in h0 h1
Definition: JMatrixNZ.sh:69
Auxiliary class for multiplexing object iterators.
string outputFile
Data structure for detector geometry and calibration.
then echo The file $DIR KM3NeT_00000001_00000000 root already please rename or remove it first
int first
index of module in detector data structure
Template L2 builder.
Definition: JBuildL2.hh:45
Detector file.
Definition: JHead.hh:130
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1954
Data structure for hit time and DOM identifier.
Definition: JMonitorL1dt.cc:34
double getFrameTime()
Get frame time duration.
Definition: JDAQClock.hh:162
#define NOTICE(A)
Definition: JMessage.hh:64
Data time slice.
virtual const pointer_type & next()
Get next element.
JHistogram()
Default constructor.
Definition: JMonitorL1dt.cc:72
Address of module in detector data structure.
void load(const JString &file_name, JDetector &detector)
Load detector from input file.
int debug
debug level
Definition: JSirene.cc:61
#define FATAL(A)
Definition: JMessage.hh:67
Scanning of objects from multiple files according a format that follows from the extension of each fi...
Direct access to module in detector data structure.
Data structure for L2 parameters.
JElement(const int __id, const double __t)
Definition: JMonitorL1dt.cc:41
General purpose class for object reading from a list of file names.
Utility class to parse command line options.
const JModuleAddress & getAddress(const JObjectID &id) const
Get address of module.
bool hasModule(const JObjectID &id) const
Has module.
virtual bool hasNext()
Check availability of next element.
KM3NeT DAQ constants, bit handling, etc.
Auxiliary data structure for histogram management.
Definition: JMonitorL1dt.cc:68
Basic data structure for L1 hit.
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62
Template definition of histogram object interface.
Definition: JHistogram.hh:27
int main(int argc, char *argv[])
Definition: Main.cpp:15