Jpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Functions
JMonitorMultiplicity.cc File Reference

Auxiliary program to analyse coincidence multiplicity;. More...

#include <iostream>
#include <string>
#include "TFile.h"
#include "TH1D.h"
#include "TH2D.h"
#include "TMath.h"
#include "TROOT.h"
#include "TVectorD.h"
#include "JAAnet/JHead.hh"
#include "JCalibrate/JCalibrateK40.hh"
#include "JDAQ/JDAQ.hh"
#include "JDAQ/JDAQEvaluator.hh"
#include "JDAQ/JDAQSummaryslice.hh"
#include "JDAQ/JDAQTimeslice.hh"
#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JDetector/JModuleRouter.hh"
#include "JDetector/JDetectorAddressMap.hh"
#include "JDetector/JDetectorAddressMapToolkit.hh"
#include "JDetector/JModuleAddressMap.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"
#include "JMath/JMathToolkit.hh"
#include "Jeep/JTimer.hh"
#include "JLang/JObjectMultiplexer.hh"
#include "JLang/JString.hh"
#include "JGizmo/JManager.hh"
#include "JGizmo/JGizmoToolkit.hh"
#include "JROOT/JROOTClassSelector.hh"
#include "JSupport/JMeta.hh"
#include "JSupport/JMonteCarloFileSupportkit.hh"
#include "JSupport/JMultipleFileScanner.hh"
#include "JSupport/JSupport.hh"
#include "JSupport/JTreeScanner.hh"
#include "JSupport/JTriggerParametersSupportkit.hh"
#include "JTrigger/JHit.hh"
#include "JTrigger/JHitR0.hh"
#include "JTrigger/JSummaryRouter.hh"
#include "JTrigger/JSuperFrame1D.hh"
#include "JTrigger/JSuperFrame2D.hh"
#include "JTools/JRange.hh"

Go to the source code of this file.

Functions

int main (int argc, char **argv)
 

Detailed Description

Auxiliary program to analyse coincidence multiplicity;.

Author
mlincetto, mkarel

Definition in file JMonitorMultiplicity.cc.

Function Documentation

int main ( int  argc,
char **  argv 
)

Definition at line 51 of file JMonitorMultiplicity.cc.

52 {
53  using namespace std;
54  using namespace KM3NETDAQ;
55  using namespace JPP;
56 
57  //-----------------------------------------------------------
58  // parameter interface
59  //-----------------------------------------------------------
60 
61  JMultipleFileScanner<JDAQTimesliceTypes_t> inputFile;
62  JLimit_t& numberOfEvents = inputFile.getLimit();
63  string outputFile;
64  string detectorFile;
65  Double_t TMax_ns;
66  JRange<double> rateVeto_Hz;
67  JRange<double> totVeto_ns;
68  JROOTClassSelector selector;
69  int debug;
70  int filterLevel;
71  int muteChannel;
72  bool twoDim;
73  bool monitorOccupancy;
74  bool notJoin;
75 
76  try {
77 
78  JParser<> zap("Auxiliary program to estimate PMT and hit multiplicities.");
79 
80  zap['f'] = make_field(inputFile);
81  zap['o'] = make_field(outputFile) = "monitormultiplicity.root";
82  zap['a'] = make_field(detectorFile);
83  zap['n'] = make_field(numberOfEvents) = JLimit::max();
84  zap['T'] = make_field(TMax_ns) = 10.0; // [ns]
85  zap['V'] = make_field(rateVeto_Hz) = JRange<double>(0, 1e5); // [Hz]
86  zap['t'] = make_field(totVeto_ns) = JRange<double>(0, 1e4); // min and max value of ToT of hits to be monitored [ns]
87  zap['C'] = make_field(selector) = getROOTClassSelection<JDAQTimesliceTypes_t>();
88  zap['d'] = make_field(debug) = 1;
89  zap['M'] = make_field(muteChannel) = -1;
90  zap['F'] = make_field(filterLevel, "0 = raw data; 1 = filtered data; 2+ = only clean frames") = 1;
91  zap['D'] = make_field(twoDim, "2D mode for background subtraction");
92  zap['O'] = make_field(monitorOccupancy, "produces 2D PMT vs multiplicity plots");
93  zap['j'] = make_field(notJoin, "do not join consecutive hits on same PMT (diagnostic purpose)");
94 
95  zap(argc, argv);
96  }
97  catch(const exception &error) {
98  FATAL(error.what() << endl);
99  }
100 
101  cout.tie(&cerr);
102 
103  //-----------------------------------------------------------
104  // check the input parameters
105  //-----------------------------------------------------------
106 
107  JTriggerParameters parameters;
108  bool hasTriggerParameters = true;
109 
110  try {
111  parameters = getTriggerParameters(inputFile);
112  NOTICE("Get trigger parameters from input." << endl);
113  }
114  catch(const JException& error) {
115  ERROR("No trigger parameters from input." << endl);
116  hasTriggerParameters = false;
117 
118  }
119 
120  // keep track of prescale factor for proper MC normalization.
121  int prescale = 1;
122 
123  if (hasTriggerParameters) {
124 
125  // for L1 (most common for sea analysis) check that time window is not larger than in trigger
126 
127  if (parameters.writeL1.prescale > 0 && (selector == "JDAQTimeslice" || selector == "JDAQTimesliceL1") &&
128  parameters.TMaxLocal_ns < TMax_ns) {
129  FATAL("TMax_ns (-T) " << TMax_ns << " ns is larger than in the trigger " << parameters.TMaxLocal_ns << " ns." << endl);
130  }
131 
132  if (selector == "JDAQTimeslice") {
133  // Jpp v8 legacy compatibility
134  prescale = (parameters.writeL1.prescale != 0) ? parameters.writeL1.prescale : parameters.writeL0.prescale;
135  }
136 
137  if (selector == "JDAQTimesliceL0") {
138  prescale = parameters.writeL0.prescale;
139  }
140 
141  if (selector == "JDAQTimesliceL1") {
142  prescale = parameters.writeL1.prescale;
143  }
144 
145  if (selector == "JDAQTimesliceSN") {
146  prescale = parameters.writeSN.prescale;
147  }
148 
149  }
150 
151  if (prescale == 0) {
152  FATAL("[R] According to trigger parameters, the " << selector << " stream should be empty.");
153  } else {
154  NOTICE("[R] Prescale factor is " << prescale << endl);
155  }
156 
157  if (!totVeto_ns.is_valid()) {
158  FATAL("Invalid time over threshold " << totVeto_ns << endl);
159  }
160 
161  //-----------------------------------------------------------
162  // load the detector file
163  //-----------------------------------------------------------
164 
165  JDetector detector;
166 
167  try {
168  load(detectorFile, detector);
169  }
170  catch(const JException& error) {
171  FATAL(error);
172  }
173 
174  if (detector.empty())
175  FATAL("Empty detector." << endl);
176 
177  const JModuleRouter router(detector);
178  JSummaryRouter summaryRouter;
179 
180  int detSize = getNumberOfModules(detector);
181 
182  int detID = detector.getID();
183 
184  JDetectorAddressMap& detectorAddressMap = getDetectorAddressMap(detID);
185 
186 
187  //-----------------------------------------------------------
188  // summary of input parameters
189  //-----------------------------------------------------------
190 
191  NOTICE("[R] Frame time [ms] " << getFrameTime() * 1e-6 << endl);
192  NOTICE("[R] Tmax [ns] " << TMax_ns << endl);
193  NOTICE("[R] Rate range [Hz] " << rateVeto_Hz << endl);
194  NOTICE("[R] ToT range [ns] " << totVeto_ns << endl);
195  NOTICE("[R] Timeslice stream " << selector << endl);
196  NOTICE("[R] Detector file has " << detSize << " DOMs." << endl);
197 
198  //-----------------------------------------------------------
199  // initialise histograms
200  //-----------------------------------------------------------
201 
202  // nominal livetime + effective livetime for each DOM
203  TH1D* h_livetime = new TH1D("LiveTime", "Livetime", 1 + detSize, 0.0, 1.0 + detSize);
204 
205  int ibin = 1;
206  h_livetime->GetXaxis()->SetBinLabel(ibin++, "Nominal");
207  for (JDetector::iterator module = detector.begin(); module != detector.end(); ++module) {
208  h_livetime->GetXaxis()->SetBinLabel(ibin++, TString(getModuleLabel(*module)));
209  }
210 
211  h_livetime->GetYaxis()->SetTitle("Livetime (s)");
212 
213 
214  // multiplicity histograms
215  const int nx = 1 + NUMBER_OF_PMTS;
216  const double xmin = -0.5;
217  const double xmax = nx - 0.5;
218 
219  typedef JManager<TString, TH1D> JManager_t;
220  typedef JManager<TString, TH2D> JManager2D_t;
221  JManager_t H(new TH1D("%_H", NULL, nx, xmin, xmax)); H->Sumw2();
222  JManager_t P(new TH1D("%_P", NULL, nx, xmin, xmax)); P->Sumw2();
223 
224  JManager2D_t HT(new TH2D("%_HT", NULL, nx, xmin, xmax, 100, -TMax_ns, +TMax_ns));
225 
226  TH2D* CO_proto = new TH2D("%_CO", NULL, NUMBER_OF_PMTS, 0.5, 0.5 + NUMBER_OF_PMTS, NUMBER_OF_PMTS, -0.5, -0.5 + NUMBER_OF_PMTS);
227 
228  if (monitorOccupancy) {
229  JModuleAddressMap memo = getDetectorAddressMap<JKM3NeT_t>().getDefaultModuleAddressMap();
230  setAxisLabels(CO_proto->GetYaxis(), memo);
231  }
232 
233  JManager2D_t CO(CO_proto);
234 
235  // benchmark histogram for time distribution of high-multiplicities
236  const int sn_mul_min = 6;
237  TH1D* time_distr = new TH1D("T_distr", "Coincidence (M > 6) distribution in timeslice", 100, 0, 1e8);
238 
239  //-----------------------------------------------------------
240  // output file and auxiliary data structures
241  //-----------------------------------------------------------
242 
243  // output file
244  TFile out(outputFile.c_str(), "recreate");
245 
246  // auxiliary counters
247 
248  int count_HRV = 0; // high-rate-veto counter
249  int count_FAF = 0; // fifo-almost-full counter
250  int count_URV = 0; // user-rate-veto counter
251 
252  int treeMismatchCount = 0;
253 
254  //-----------------------------------------------------------
255  // file pre-processing
256  //-----------------------------------------------------------
257 
258  JTreeScanner<JDAQSummaryslice, JDAQEvaluator> summaryFile(inputFile);
259 
260  JObjectMultiplexer<JDAQTimesliceTypes_t, JDAQTimeslice> in(inputFile, selector);
261 
262  vector<JHitR0> filteredData;
263 
264  typedef JSuperFrame2D<JHitR0> JSuperFrame2D_t;
265  typedef JSuperFrame1D<JHitR0> JSuperFrame1D_t;
266  const JMatch_t match(TMax_ns);
267 
268  counter_type counter = 0;
269 
270  // process Monte Carlo aanet header
271  bool hasMCHeader = true ;
272  bool isMupage = false;
273  double liveTimeMC = 0 ;
274  double liveTimeMCErr = 0 ;
275 
276  JHead mc_header;
277 
278  try {
279  mc_header = (JHead)getHeader(inputFile);
280  }
281  catch (const JLANG::JNullPointerException error) {
282  hasMCHeader = false;
283  }
284 
285  NOTICE("[R] File source is " << (hasMCHeader ? "MC" : "DAQ") << "." << endl);
286  if (hasMCHeader) {
287  liveTimeMC = mc_header.livetime.numberOfSeconds;
288  liveTimeMCErr = mc_header.livetime.errorOfSeconds;
289  if (liveTimeMC > 0) {
290  isMupage = true;
291  NOTICE("[R] mupage live time is (" << liveTimeMC << " +/- " << liveTimeMCErr << ") s." << endl);
292  } else {
293  NOTICE("[R] aanet header found but mupage live time is zero. Likely not a mupage file.");
294  }
295  }
296 
297  // retrieving and exporting run number to tree metadata;
298 
299  int runNumber = 0;
300 
301  if (summaryFile.hasNext()) {
302  runNumber = summaryFile.getEntry(0)->getRunNumber();
303  NOTICE("[R] Processing run #" << runNumber << endl);
304  summaryFile.rewind();
305  }
306 
307  //-----------------------------------------------------------
308  // file processing: loop over timeslices
309  //-----------------------------------------------------------
310 
311  for ( ; in.hasNext() && counter != inputFile.getLimit(); ++counter) {
312 
313  STATUS("Entry: " << setw(10) << counter << '\r');
314 
315  const JDAQTimeslice* timeslice = in.next();
316  const Long64_t index = summaryFile.find(*timeslice);
317  JDAQSummaryslice* summary = (index != -1 ? summaryFile.getEntry(index) : NULL);
318 
319  summaryRouter.update(summary);
320 
321  // check correspondence between summary data and timeslice data
322  if (summary != NULL) {
323  if (timeslice->getFrameIndex() != summary->getFrameIndex()) {
324  ++treeMismatchCount;
325  ERROR("[R>T] Frame indices do not match [counter=" << counter << ", timeslice=" << timeslice->getFrameIndex() << ", summaryslice=" << summary->getFrameIndex() << "]" << endl);
326  if (treeMismatchCount > 1) {
327  FATAL("ROOT TTrees not aligned. Run #" << runNumber << endl);
328  }
329  continue;
330  }
331  }
332 
333  //-----------------------------------------------------------
334  // timeslice processing: loop over superframes
335  //-----------------------------------------------------------
336  for (JDAQTimeslice::const_iterator super_frame = timeslice->begin(); super_frame != timeslice->end(); ++super_frame) {
337 
338  // empty superframes are not skipped because SN timeslices could have empty superframes even if the module was active;
339  // if DAQ error or FIFO full -> frame counts as dead time;
340  // if frame is missing but no errors -> frame counts as live time (SN stream has many empty frames);
341 
342  int moduleID = super_frame->getModuleID();
343 
344  if (!router.hasModule(moduleID)) {
345  // module is not mapped in detector file, nothing can be done
346  continue;
347  }
348 
349  JModuleAddressMap moduleAddressMap = detectorAddressMap.get(moduleID);
350 
351  //----------------------------------------------------------
352  // retrieval of summary data from superframe or summaryframe
353  //----------------------------------------------------------
354 
355  vector<double> rate_Hz(NUMBER_OF_PMTS, 0.0);
356  JDAQFrameStatus status;
357 
358  if (summary == NULL) {
359  WARNING("[R>T>F] Missing summary for timeslice." << endl);
360  // PMT rates estimated from timeslice data
361  for (JDAQSuperFrame::const_iterator i = super_frame->begin(); i != super_frame->end(); ++i) {
362  rate_Hz[i->getPMT()] += 1e9 / getFrameTime();
363  }
364 
365  // frame status retrieved from superframe
366  status = super_frame->getDAQFrameStatus();
367 
368  } else {
369 
370  // summary frame lookup
371 
372  JDAQSummaryFrame summary_frame;
373 
374  if (!summaryRouter.hasSummaryFrame(super_frame->getModuleIdentifier())) {
375  summaryRouter.print(cout);
376  FATAL("[R>D>F] Summary frame not found, module = " << super_frame->getModuleID() << endl);
377  } else {
378  summary_frame = summaryRouter.getSummaryFrame(super_frame->getModuleID());
379  }
380 
381  // PMT rates recovered from summaryslice data
382  for (int i = 0; i != NUMBER_OF_PMTS; ++i) {
383  rate_Hz[i] = summary_frame.getRate(i);
384  }
385 
386  // frame status retrieved from summaryframe
387  status = summary_frame.getDAQFrameStatus();
388 
389  }
390 
391  //---------------------------
392  // processing of summary data
393  //---------------------------
394 
395  // DAQ error (e.g. missing UDP trailer) -> discard frame
396  if (!status.testDAQStatus()) {
397  WARNING("[R>D>F] DAQ status not okay in timeslice #" << timeslice->getFrameIndex() << ", DOM=" << super_frame->getModuleID() << endl);
398  continue;
399  }
400 
401  vector<bool> veto(NUMBER_OF_PMTS, false);
402  bool moduleLiveTime = 0;
403 
404  // default status of module is alive if at least one PMT is active
405  for (vector<double>::const_iterator r = rate_Hz.begin(); !moduleLiveTime && r != rate_Hz.end(); ++r) {
406  if ((*r) > 0) {
407  moduleLiveTime = 1;
408  }
409  }
410 
411  int frame_HRV_count = status.countHighRateVeto();
412  int frame_FAF_count = status.countFIFOStatus();
413  int frame_URV_count = 0; // user rate veto
414 
415  for (int i = 0; i != NUMBER_OF_PMTS; ++i) {
416  if (veto[i] = (veto[i] || !rateVeto_Hz(rate_Hz[i]))) {
417  frame_URV_count++;
418  }
419  }
420 
421  count_HRV += frame_HRV_count;
422  count_FAF += frame_FAF_count;
423  count_URV += frame_URV_count;
424 
425  int frame_VETO_count = frame_HRV_count + frame_FAF_count + frame_URV_count;
426 
427 
428  if (filterLevel >= 2 && frame_VETO_count > 0) {
429  moduleLiveTime = 0;
430  fill(veto.begin(), veto.end(), true);
431  } else if (filterLevel >= 1 && frame_HRV_count + frame_FAF_count > 0) {
432  for (int i = 0; i != NUMBER_OF_PMTS; ++i) {
433  veto[i] = (veto[i] || (status.testFIFOStatus(i) || status.testHighRateVeto(i)));
434  }
435  }
436 
437  // recover DOM in case the dead channel corresponds to the muted channel
438 
439  if (muteChannel != -1) {
440 
441  // veto conditions are not exclusive
442  int mute_VETO_count =
443  status.testHighRateVeto(muteChannel) +
444  status.testFIFOStatus(muteChannel ) +
445  !rateVeto_Hz(rate_Hz[muteChannel] );
446 
447  if (mute_VETO_count == frame_VETO_count) {
448  moduleLiveTime = 1;
449  fill(veto.begin(), veto.end(), false);
450  }
451 
452  veto[muteChannel] = true;
453 
454  }
455 
456  // get address and module for current superframe DOM;
457  const JModuleAddress& address = router.getAddress(super_frame->getModuleID());
458  const JModule& module = detector.getModule(address);
459  const TString moduleLabel = getModuleLabel(module);
460 
461  // acknowledge DOM is live in this timeslice for analysis purpose;
462  if (moduleLiveTime) {
463  h_livetime->Fill(moduleLabel, getFrameTime() * 1e-9 * moduleLiveTime);
464  }
465 
466  //--------------------------------------------------------------------
467  // superframe processing: hit joining, calibration, filter, processing
468  //--------------------------------------------------------------------
469 
470  JSuperFrame2D_t& buffer = JSuperFrame2D_t::demultiplex(*super_frame, router.getModule(super_frame->getModuleID()));
471 
472  if (!notJoin) {
473  for (JSuperFrame2D_t::iterator i = buffer.begin(); i != buffer.end(); ++i) {
474  i->join(match);
475  }
476  }
477 
478  JSuperFrame1D_t& data = JSuperFrame1D_t::multiplex(buffer);
479 
480  data.pop_back(); // remove end marker
481 
482  filteredData.clear();
483 
484  for (vector<JHitR0>::const_iterator h = data.begin(); h != data.end(); ++h) {
485 
486  const int pmt = h->getPMT();
487 
488  // veto and calibration;
489  if (!veto[pmt] &&
490  rateVeto_Hz(rate_Hz[pmt]) &&
491  totVeto_ns(h->getToT()) ) {
492  filteredData.push_back(*h);
493  }
494  }
495 
496  // processing signal;
497 
498  if (filteredData.size() > 1) {
499 
500  TH1D* h_hit = H[moduleLabel];
501  TH1D* h_pmt = P[moduleLabel];
502 
503  TH2D* h2_hit = HT[moduleLabel];
504  TH2D* h2_co = CO[moduleLabel];
505 
506  sort(filteredData.begin(), filteredData.end());
507 
508  //-----------------------------------------------------------
509  // superframe processing: loop on calibrated hits
510  //-----------------------------------------------------------
511  for (vector<JHitR0>::const_iterator p = filteredData.begin(); p != filteredData.end(); ) {
512 
514 
515  // coincidence building using a fixed time window
516  while (++q != filteredData.end() && q->getT() - p->getT() <= TMax_ns ) {}
517 
518  // multiplicity estimation
519  int hit_multiplicity = distance(p,q);
520  set<int> pmthit;
521  for (vector<JHitR0>::const_iterator __p = p ; __p != q; ++__p ) { pmthit.insert(__p->getPMT()) ; }
522  int pmt_multiplicity = pmthit.size() ;
523 
524  // histogram filling
525  h_pmt->Fill(pmt_multiplicity);
526  h_hit->Fill(hit_multiplicity);
527 
528  if (twoDim && hit_multiplicity > 1) {
529  const double W = factorial(hit_multiplicity, 2);
530  for (vector<JHitR0>::const_iterator __p = p; __p != q; ++__p) {
531  for (vector<JHitR0>::const_iterator __q = __p; ++__q != q; ) {
532  double dt = JCombinatorics::getSign(__p->getPMT(), __q->getPMT()) * (__q->getT() - __p->getT());
533  h2_hit->Fill(hit_multiplicity, dt, 1.0/W);
534  }
535  }
536  }
537 
538  if (monitorOccupancy) {
539  for (set<int>::const_iterator r = pmthit.begin(); r != pmthit.end(); ++r) {
540  JPMTPhysicalAddress pmtAddress = moduleAddressMap.getAddressTranslator(*r);
541  TString pmtLabel(pmtAddress.toString());
542  h2_co->Fill(pmt_multiplicity, pmtLabel, 1.0/pmt_multiplicity);
543  }
544  }
545 
546  // benchmarking
547  if (hit_multiplicity >= sn_mul_min) {
548  time_distr->Fill(p->getT());
549  }
550 
551  p = q;
552 
553  }
554 
555  }
556 
557  }
558 
559  }
560 
561  NOTICE("[R] " << counter << " timeslices processed." << endl);
562 
563  if (isMupage) {
564  h_livetime->Scale(liveTimeMC / (getFrameTime() * 1e-9 * counter * prescale));
565  }
566 
567  // nominal live time;
568  double nominalLiveTime = (isMupage ? liveTimeMC : counter * getFrameTime() * 1e-9) / prescale;
569  h_livetime->Fill("Nominal", nominalLiveTime);
570 
571  out.cd();
572 
573  int nLiveDOMs = H.size();
574 
575  h_livetime->Write();
576 
577  // writing hit multiplicities is redundant given new options
578  // H.Write(out);
579 
580  P.Write(out);
581 
582  if (twoDim) {
583  HT.Write(out);
584  }
585 
586  if (monitorOccupancy) {
587  CO.Write(out);
588  }
589 
590  time_distr->Write();
591 
592  double rate_HRV = count_HRV / (1.0 * counter);
593  double rate_FAF = count_FAF / (1.0 * counter);
594  double rate_URV = count_URV / (1.0 * counter);
595  double biolumIndex = (rate_HRV + rate_FAF) / (nLiveDOMs * NUMBER_OF_PMTS);
596  NOTICE("[R] Average HRV count per timeslice: " << rate_HRV << endl);
597  NOTICE("[R] Average FIFO-almost-full count per timeslice: " << rate_FAF << endl);
598  NOTICE("[R] Average user rate veto count per timeslice: " << rate_URV << endl);
599  NOTICE("[R] " << H.size() << " DOMs were active in the run." << endl);
600  NOTICE("[R] Bioluminescence proxy index: " << biolumIndex << endl);
601 
602  putObject(&out, JMeta(argc,argv));
603 
604  out.Close();
605 
606  NOTICE(endl);
607 
608  return (counter ? 0 : 1);
609 
610 }
static const double H
Planck constant [eV s].
Definition: JConstants.hh:25
Auxiliary class for ROOT I/O of application specific meta data.
Definition: JMeta.hh:71
Utility class to parse command line options.
Definition: JParser.hh:1410
#define WARNING(A)
Definition: JMessage.hh:63
Data structure for all trigger parameters.
double getRate(const int tdc, const double factor=1.0) const
Get count rate.
int countFIFOStatus() const
Count FIFO status.
Auxiliary class to manage set of histograms.
#define STATUS(A)
Definition: JMessage.hh:61
Long64_t counter_type
Type definition for counter.
Definition: JCounter.hh:24
JSuperFrame2D< hit_type > JSuperFrame2D_t
Definition: JDataFilter.cc:93
string outputFile
long long int factorial(const long long int n)
Determine factorial.
Definition: JMathToolkit.hh:40
int getFrameIndex() const
Get frame index.
Head getHeader(const JMultipleFileScanner_t &file_list)
Get Monte Carlo header.
JLimit JLimit_t
Type definition of limit.
Definition: JLimit.hh:214
Exception for null pointer operation.
Definition: JException.hh:198
Hit data structure.
Definition: JDAQHit.hh:36
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1836
double getFrameTime()
Get frame time duration.
Definition: JDAQClock.hh:162
Data storage class for rate measurements of all PMTs in one module.
void setAxisLabels(TAxis *axis, const JModuleAddressMap &memo)
initialize axis with PMT address labels
#define NOTICE(A)
Definition: JMessage.hh:62
#define ERROR(A)
Definition: JMessage.hh:64
Data time slice.
void load(const JString &file_name, JDetector &detector)
Load detector from input file.
int debug
debug level
Definition: JSirene.cc:59
JDetectorAddressMap & getDetectorAddressMap()
Get detector address map.
Monte Carlo run header.
Definition: JHead.hh:727
#define FATAL(A)
Definition: JMessage.hh:65
bool testHighRateVeto() const
Test high-rate veto status.
std::string getModuleLabel(const JModuleLocation &location)
Get module label (DU-floor) for JMonitor applications.
int countHighRateVeto() const
Count high-rate veto status.
std::vector< frame_type >::iterator iterator
2-dimensional frame with time calibrated data from one optical module.
const JLimit & getLimit() const
Get limit.
Definition: JLimit.hh:72
static const int NUMBER_OF_PMTS
Total number of PMTs in module.
Definition: JDAQ.hh:26
bool testDAQStatus() const
Test DAQ status of packets.
const JDAQFrameStatus & getDAQFrameStatus() const
Get DAQ frame status.
int getNumberOfModules(const JDetector &detector)
Get number of modules.
JTriggerParameters getTriggerParameters(const JMultipleFileScanner_t &file_list)
Get trigger parameters.
JSuperFrame1D< hit_type > JSuperFrame1D_t
Definition: JDataFilter.cc:92
bool testFIFOStatus() const
Test FIFO status.
bool putObject(TDirectory *dir, const T &object)
Write object to ROOT directory.