Jpp 19.3.0-rc.3
the software that should make you happy
Loading...
Searching...
No Matches
JRandomTimesliceWriter.cc File Reference

Example program to test generation of KM3NETDAQ::JDAQTimeslice with random data. More...

#include <string>
#include <iostream>
#include <iomanip>
#include <vector>
#include "TROOT.h"
#include "TFile.h"
#include "TH1D.h"
#include "TH2D.h"
#include "TRandom3.h"
#include "km3net-dataformat/online/JDAQ.hh"
#include "km3net-dataformat/online/JDAQClock.hh"
#include "km3net-dataformat/online/JDAQTimeslice.hh"
#include "JTimeslice/JRandomTimeslice.hh"
#include "JROOT/JManager.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"

Go to the source code of this file.

Functions

int main (int argc, char **argv)
 

Detailed Description

Example program to test generation of KM3NETDAQ::JDAQTimeslice with random data.

Author
mdejong

Definition in file examples/JTimeslice/JRandomTimesliceWriter.cc.

Function Documentation

◆ main()

int main ( int argc,
char ** argv )

Definition at line 29 of file examples/JTimeslice/JRandomTimesliceWriter.cc.

30{
31 using namespace std;
32 using namespace JPP;
33 using namespace KM3NETDAQ;
34
35 string outputFile;
36 Long64_t numberOfSlices;
37 double rate_Hz;
38 pair<size_t, double> recycling;
39 ULong_t seed;
40 int debug;
41
42 try {
43
44 JParser<> zap("Example program to test generation of time slices with random data.");
45
46 zap['o'] = make_field(outputFile, "output file") = "timeslice.root";
47 zap['n'] = make_field(numberOfSlices);
48 zap['B'] = make_field(rate_Hz, "background rate [Hz]");
49 zap['N'] = make_field(recycling, "number of recycles / time interval for sampling data") = make_pair(0, 0.0);
50 zap['S'] = make_field(seed, "seed") = 0;
51 zap['d'] = make_field(debug, "debug") = 0;
52
53 zap(argc, argv);
54 }
55 catch(const exception &error) {
56 FATAL(error.what() << endl);
57 }
58
59
60 gRandom->SetSeed(seed);
61
62 const JDAQHit::JPMT_t PMT_L0[] = { 0, 1 };
63
64 JManager<size_t, TH1D> H0(new TH1D("h0[%]", NULL, 100, 0.0, getFrameTime()));
65
66 for (Long64_t counter = 0; counter != numberOfSlices; ++counter) {
67
68 STATUS("slice: " << setw(10) << counter << '\r'); DEBUG(endl);
69
70 JRandomTimeslice timeslice;
71
72 timeslice.resize(1);
73
74
75 vector<JDAQHit> buffer;
76
77 if (rate_Hz > 0.0) {
78
79 double period_ns = 1.0e9 / (NUMBER_OF_PMTS * rate_Hz);
80
81 for (double t1 = 0.0 * getFrameTime() + gRandom->Exp(period_ns); t1 < 0.5 * getFrameTime(); t1 += gRandom->Exp(period_ns)) {
82 buffer.push_back(JDAQHit(PMT_L0[0], t1, 0));
83 }
84
85 period_ns *= 2.0;
86
87 for (double t1 = 0.5 * getFrameTime() + gRandom->Exp(period_ns); t1 < 1.0 * getFrameTime(); t1 += gRandom->Exp(period_ns)) {
88 buffer.push_back(JDAQHit(PMT_L0[1], t1, 0));
89 }
90 }
91
92 timeslice[0].add(buffer.size(), buffer.data());
93
94 TH1D* h0 = H0[0];
95
96 for (JDAQSuperFrame::const_iterator hit = timeslice[0].begin(); hit != timeslice[0].end(); ++hit) {
97 h0->Fill((Double_t) hit->getT());
98 }
99
100 for (size_t i = 1; i <= recycling.first; ++i) {
101
102 timeslice.recycle(recycling.second);
103
104 TH1D* h0 = H0[i];
105
106 for (JDAQSuperFrame::const_iterator hit = timeslice[0].begin(); hit != timeslice[0].end(); ++hit) {
107 h0->Fill((Double_t) hit->getT());
108 }
109 }
110 }
111 STATUS(endl);
112
113
114 TFile out(outputFile.c_str(), "recreate");
115
116 out << H0;
117
118 out.Write();
119 out.Close();
120}
string outputFile
#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
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition JParser.hh:2142
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
Hit data structure.
Definition JDAQHit.hh:35
unsigned char JPMT_t
PMT channel in FPGA.
Definition JDAQHit.hh:38
JDAQTimeslice & add(const JDAQTimeslice &timeslice)
Add another timeslice.
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
KM3NeT DAQ data structures and auxiliaries.
Definition DataQueue.cc:39
double getFrameTime()
Get frame time duration.
Definition JDAQClock.hh:162
Timeslice with random data.
void recycle(const double T_ns)
Recycle time slice by randomly shuffling time intervals of data.