66int main(
int argc,
char **argv)
70 JLimit_t& numberOfEvents = inputFile.getLimit();
78 bool globalOutputOnly;
86 JParser<> zap(
"Example program to examine rates as a function of time on ms-level timescales.");
93 zap[
'T'] =
make_field(TMax_ns,
"Time window for local coincidences (if 0 run in L0 mode)") = 0.0;
94 zap[
'B'] =
make_field(binWidth_ms,
"Bin width (experimental)") = 1;
95 zap[
'V'] =
make_field(backVeto,
"Apply retroactive veto.");
96 zap[
'D'] =
make_field(mode2D,
"L1 mode: create 2D histogram with time differences of coincidences (heavy memory usage, ignored if TMax_ns = 0)");
97 zap[
'O'] =
make_field(globalOutputOnly,
"Write only aggregate histograms");
98 zap[
'M'] =
make_field(multiplicityRange,
"L1 mode: multiplicity range (ignored if TMax_ns = 0)") =
JRange<int>(2,31);
105 catch(
const exception &error) {
106 FATAL(error.what() << endl);
111 FATAL(
"Frame time must be an integer multiple of bin width");
131 pts->configure(inputFile);
133 int fEnd = pts->rbegin()->getFrameIndex();
134 int fStart = pts->begin( )->getFrameIndex();
136 if (fEnd > inputFile.getUpperLimit()) {
137 fEnd = fStart + inputFile.getUpperLimit();
141 double tStart_ms = (fStart - 1) *
getFrameTime() / 1.0e6;
143 int runNumber = pts->begin()->getRunNumber();
145 TString runTag = Form(
"%d" , runNumber);
147 double tRun_ms = tEnd_ms - tStart_ms;
149 NOTICE(
"START/END/DURATION [s]: " << tStart_ms / 1000 <<
" " << tEnd_ms / 1000 <<
" " << tRun_ms / 1000 << endl);
168 int nx = (tEnd_ms - tStart_ms) / binWidth_ms;
170 JManager_t RT_DOM(
new TH1F(
"RT_%", NULL, nx , tStart_ms, tEnd_ms));
171 JManager_t NC_DOM(
new TH1F(
"NC_%", NULL, nx / 100, tStart_ms, tEnd_ms));
173 TString rt_tag = runTag + TString(
".RT_DET_%");
174 TString nc_tag = runTag + TString(
".NC_DET_%");
176 JManager_t RT_DET(
new TH1F(rt_tag, NULL, nx, tStart_ms, tEnd_ms));
177 JManager_t NC_DET(
new TH1F(nc_tag, NULL, nx / 100, tStart_ms, tEnd_ms));
188 if (nx >= max_size) {
189 FATAL(
"2D histogram size not supported by ROOT file output; limit input size (-n) below " << floor(max_size / 100.0) << endl);
192 TString rt2d_tag = runTag + TString(
".RT2D_DET");
194 hT =
new h2d_t(rt2d_tag, NULL, nx, tStart_ms, tEnd_ms, ny, -TMax_ns, +TMax_ns);
209 if (pts->hasNext()) {
210 curr = *(pts->next());
213 FATAL(
"Input file is too short.");
218 for ( ; pts->hasNext() && counter != inputFile.getLimit(); ++counter) {
220 STATUS(
"timeslice: " << setw(10) << counter <<
'\r');
DEBUG(endl);
234 if (backVeto && ((in - ic) == 1)) {
236 summaryRouter.
update(nextSummary);
241 for (JDAQTimeslice::const_iterator frame = curr.begin(); frame != curr.end(); ++frame) {
243 const int moduleID = frame->getModuleID();
245 const string moduleLabel =
getLabel(module);
247 TH1F* RD = RT_DOM[moduleLabel];
253 if (frame->testHighRateVeto() || frame->testFIFOStatus()) {
255 veto[pmt] = ( frame->testHighRateVeto(pmt) || frame->testFIFOStatus(pmt) );
261 if (nextSummary != NULL) {
275 NC_DOM[moduleLabel]->Fill(tTimeslice_ms, count(veto.begin(), veto.end(),
false));
279 JSuperFrame2D_t& buffer2D = JSuperFrame2D_t::demultiplex(*frame, module, totSelector_ns);
281 for (JSuperFrame2D_t::iterator i = buffer2D.begin(); i != buffer2D.end(); ++i) {
282 if (veto[i->getPMTAddress()]) {
289 JSuperFrame1D_t& data = JSuperFrame1D_t::multiplex(buffer2D);
294 if (data.size() > 1) {
299 sort(data.begin(), data.end());
304 for (vector<hit_type>::const_iterator p = data.begin(); p != data.end(); p++) {
306 const double tHit_ms = tTimeslice_ms + (p->getT() / 1.0e6);
308 if (TMax_ns == 0.0) {
314 vector<hit_type>::const_iterator q = p;
316 while (++q != data.end() && q->getT() - p->getT() <= TMax_ns) {}
320 if (multiplicityRange(M)) {
328 for (vector<hit_type>::const_iterator __p = p; __p != q; ++__p) {
330 for (vector<hit_type>::const_iterator __q = __p; ++__q != q; ) {
334 hT->Fill(tHit_ms, dt, 1.0 / W);
352 if (nextSummary != NULL) {
delete nextSummary; }
365 NOTICE(
"Processing histograms." << endl);
367 for (JDetector::const_iterator module =
detector.begin(); module !=
detector.end(); ++module) {
369 string moduleLabel =
getLabel(*module);
371 if (RT_DOM.count(moduleLabel) && NC_DOM.count(moduleLabel)) {
373 TH1F* rt = RT_DOM.at(moduleLabel);
374 TH1F* nc = NC_DOM.at(moduleLabel);
376 for (
int b = 1; b <= rt->GetXaxis()->GetNbins(); b++) {
378 double r = rt->GetBinContent(b);
379 double t = rt->GetBinCenter( b);
381 RT_DET[
"SUM"]->Fill(t, r);
386 for (
int b = 1; b <= nc->GetXaxis()->GetNbins(); b++) {
388 double n = nc->GetBinContent(b);
389 double t = nc->GetBinCenter( b);
391 NC_DET[
"SUM"]->Fill(t,
n);
397 DEBUG(moduleLabel <<
" not active." << endl);
403 NOTICE(
"Writing output file" << endl);
409 NOTICE(
"Writing 1D histograms" << endl);
414 NOTICE(
"Writing 2D histogram" << endl);
420 if (!globalOutputOnly) {
422 TString dir_tag = runTag + TString(
".Modules");
424 NOTICE(
"Writing individual modules histograms" << endl);
426 TDirectory* dir = out.mkdir(dir_tag);
431 NOTICE(
"Closing file" << endl);