Jpp 19.3.0-rc.3
the software that should make you happy
Loading...
Searching...
No Matches
examples/JDAQ/JDAQTimeslice.cc
Go to the documentation of this file.
1
2#include <string>
3#include <iostream>
4#include <iomanip>
5#include <limits>
6
7#include "TROOT.h"
8#include "TFile.h"
9#include "TH1D.h"
10
13
15#include "JSupport/JSupport.hh"
16#include "JROOT/JManager.hh"
17#include "Jeep/JeepToolkit.hh"
18
20
21#include "Jeep/JParser.hh"
22#include "Jeep/JMessage.hh"
23
24
25/**
26 * \file
27 *
28 * Example program to histogram KM3NETDAQ::JDAQTimeslice.
29 * \author mdejong
30 */
31int main(int argc, char **argv)
32{
33 using namespace std;
34 using namespace JPP;
35 using namespace KM3NETDAQ;
36
38 JLimit_t& numberOfEvents = inputFile.getLimit();
39 string outputFile;
41 int debug;
42
43 try {
44
45 JParser<> zap("Example program to histogram timeslice data.");
46
47 zap['f'] = make_field(inputFile);
48 zap['o'] = make_field(outputFile) = "timeslice.root";
49 zap['n'] = make_field(numberOfEvents) = JLimit::max();
50 zap['P'] = make_field(PMT) = JDAQPMTIdentifier(-1, -1);
51 zap['d'] = make_field(debug) = 2;
52
53 zap(argc, argv);
54 }
55 catch(const exception& error) {
56 FATAL(error.what() << endl);
57 }
58
59
60 const bool zoom = (PMT.getModuleID() != -1 &&
61 PMT.getPMTAddress() != -1);
62
63 JManager<string, TH1D> H0(new TH1D("h0[%]", NULL, numeric_limits<JDAQHit::JPMT_t>::max(), -0.5, numeric_limits<JDAQHit::JPMT_t>::max() - 0.5));
64 JManager<string, TH1D> H1(new TH1D("h1[%]", NULL, numeric_limits<JDAQHit::JTOT_t>::max(), -0.5, numeric_limits<JDAQHit::JTOT_t>::max() - 0.5));
65 JManager<string, TH1D> H2(new TH1D("h2[%]", NULL, (zoom ? 20000000 : 1000), -0.5, getFrameTime() + 0.5));
66 JManager<string, TH1D> H3(new TH1D("h3[%]", NULL, 100, 0.0, 7.0));
67
68 H0->Sumw2();
69 H1->Sumw2();
70 H2->Sumw2();
71
73
74 for (counter_type counter = 0; in.hasNext(); ++counter) {
75
76 STATUS("event: " << setw(10) << counter << '\r'); DEBUG(endl);
77
78 JDAQTimeslice* timeslice = in.next();
79
80 TH1D* h0 = H0[getClassname(timeslice->GetName())];
81 TH1D* h1 = H1[getClassname(timeslice->GetName())];
82 TH1D* h2 = H2[getClassname(timeslice->GetName())];
83 TH1D* h3 = H3[getClassname(timeslice->GetName())];
84
85 for (JDAQTimeslice::const_iterator frame = timeslice->begin(); frame != timeslice->end(); ++frame) {
86
87 h3->Fill(log10((double) frame->size()));
88
89 for (JDAQSuperFrame::const_iterator hit = frame->begin(); hit != frame->end(); ++hit) {
90
91 if (JDAQPMTIdentifier::compare(PMT, JDAQPMTIdentifier(frame->getModuleID(), hit->getPMT()))) {
92 h0->Fill(hit->getPMT());
93 h1->Fill(hit->getToT());
94 h2->Fill(hit->getT());
95 }
96 }
97 }
98 }
99 STATUS(endl);
100
101 TFile out(outputFile.c_str(), "recreate");
102
103 out << H0 << H1 << H2 << H3;
104
105 out.Write();
106 out.Close();
107}
JDAQPMTIdentifier PMT
Command line options.
string outputFile
Dynamic ROOT object management.
General purpose messaging.
#define DEBUG(A)
Message macros.
Definition JMessage.hh:62
#define STATUS(A)
Definition JMessage.hh:63
#define FATAL(A)
Definition JMessage.hh:67
int debug
debug level
Definition JSirene.cc:72
Scanning of objects from multiple files according a format that follows from the extension of each fi...
Utility class to parse command line options.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition JParser.hh:2142
ROOT TTree parameter settings of various packages.
Auxiliary methods for handling file names, type names and environment.
Auxiliary class for multiplexing object iterators.
virtual bool hasNext() override
Check availability of next element.
virtual const pointer_type & next() override
Get next element.
Utility class to parse command line options.
Definition JParser.hh:1698
Auxiliary class to manage set of compatible ROOT objects (e.g. histograms) using unique keys.
Definition JManager.hh:47
void Write(TDirectory &out, const bool wm=false)
Write objects to file.
Definition JManager.hh:304
General purpose class for object reading from a list of file names.
Hit data structure.
Definition JDAQHit.hh:35
int getModuleID() const
Get module identifier.
int getPMTAddress() const
Get PMT identifier.
static bool compare(const JDAQPMTIdentifier &first, const JDAQPMTIdentifier &second)
Compare PMT identifiers.
int main(int argc, char **argv)
std::string getClassname(const std::string &type_name)
Get type name, i.e. part after JEEP::TYPENAME_SEPARATOR.
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Long64_t counter_type
Type definition for counter.
KM3NeT DAQ data structures and auxiliaries.
Definition DataQueue.cc:39
double getFrameTime()
Get frame time duration.
Definition JDAQClock.hh:162
Auxiliary class for defining the range of iterations of objects.
Definition JLimit.hh:45
static counter_type max()
Get maximum counter value.
Definition JLimit.hh:128