62 int main(
int argc,
char **argv)
69 JLimit_t& numberOfEvents = inputFile.getLimit();
74 int numberOfTimeslices;
85 JParser<> zap(
"Example program to search for correlations between triggered events and timeslice data.");
90 zap[
'n'] =
make_field(numberOfEvents) = JLimit::max();
93 zap[
'N'] =
make_field(numberOfTimeslices) = 100;
96 zap[
'C'] =
make_field(selector) =
"", getROOTClassSelection<JDAQTimesliceTypes_t>();
102 zap.
read(argc, argv);
104 catch(
const exception& error) {
105 FATAL(error.what() << endl);
130 const Double_t
xmin = -(numberOfTimeslices + 1) *
getFrameTime() - 0.5*binWidth_ns;
131 const Double_t
xmax = +(numberOfTimeslices + 1) *
getFrameTime() + 0.5*binWidth_ns;
132 const Int_t nx = (Int_t) ((
xmax -
xmin) / binWidth_ns);
134 JManager_t manager(
new TH1D(
"M_%", NULL, nx,
xmin,
xmax));
142 if (selector ==
"") {
150 FATAL(
"No timeslice data." << endl);
153 NOTICE(
"Selected class " << ps->getClassName() << endl);
159 ps->configure(inputFile);
180 ps->setLimit(inputFile.getLimit());
195 for (JDetector::const_iterator module =
detector.begin(); module !=
detector.end(); ++module) {
196 oosDomsId[module->getID()] =
false;
202 STATUS(
"event: " << setw(10) << in.getCounter() <<
'\r');
DEBUG(endl);
209 bool test_OOS =
false;
214 if(oosDomsId[hit->getModuleID()]) {
220 if(test_OOS)
continue;
225 buffer[
event->getFrameIndex()].push_back(t0);
232 while (ps->hasNext()) {
234 STATUS(
"event: " << setw(10) << ps->getCounter() <<
'\r');
DEBUG(endl);
238 map_type::const_iterator p = buffer.lower_bound(timeslice->
getFrameIndex() - numberOfTimeslices);
239 map_type::const_iterator q = buffer.upper_bound(timeslice->
getFrameIndex() + numberOfTimeslices);
241 int number_of_events = 0;
243 for (map_type::const_iterator i = p; i != q; ++i) {
244 number_of_events += i->second.size();
247 if (number_of_events == 0) {
251 for (JDAQTimeslice::const_iterator frame = timeslice->begin(); frame != timeslice->end(); ++frame) {
255 buildL1(*frame, router.
getModule(frame->getModuleID()), back_inserter(
data));
257 TH1D* h1 = manager[frame->getModuleID()];
261 const double t1 = *hit + frame->getFrameIndex() *
getFrameTime();
263 for (map_type::const_iterator i = p; i != q; ++i) {
264 for (map_type::mapped_type::const_iterator
j = i->second.begin();
j != i->second.end(); ++
j) {
266 const double t0 = *
j;
279 Double_t most_OOS_peak = 0;
280 Int_t most_OOS_ID = 0;
282 TF1
f1(
"f1",
"[0]*exp(-0.5*(x-[1])*(x-[1])/([2]*[2])) + [3]");
290 for (JDetector::const_iterator module =
detector.begin(); module !=
detector.end(); ++module) {
294 status2[module->getID()] = DEFAULT;
296 JManager_t::iterator p = manager.find(module->getID());
298 if (p == manager.end() || p->second->GetEntries() == 0) {
300 status2[module->getID()] = NODATA;
305 TH1D* h1 = p->second;
311 Double_t
sigma = 250.0;
314 for (
int i = 1; i != h1->GetNbinsX(); ++i) {
316 const Double_t
x = h1->GetBinCenter (i);
317 const Double_t
y = h1->GetBinContent(i);
327 ymin /= h1->GetNbinsX();
329 f1.SetParameter(0, ymax);
330 f1.SetParameter(1, x0);
331 if (binWidth_ns <
sigma)
334 f1.FixParameter(2, binWidth_ns/sqrt(12.0));
335 f1.SetParameter(3, ymin);
337 for (Int_t i = 0; i !=
f1.GetNpar(); ++i) {
338 f1.SetParError(i, 0.0);
343 h1->Fit(&
f1, option.c_str(),
"same", x0 - 5 *
sigma, x0 + 5 *
sigma);
346 f1.GetParameter(0) >=
f1.GetParameter(3));
348 if (status) status2[module->getID()] = IN_SYNC;
355 for (Int_t i = 1,
ns = -(numberOfTimeslices + 1); i <= h1->GetNbinsX(); ++i) {
357 const Double_t
x = h1->GetBinCenter (i);
358 const Double_t
y = h1->GetBinContent(i);
368 DEBUG(
"Module " << setw(8) << module->getID() << endl);
372 const Double_t
y = bg.getTotal() * i->second.getCount() / bg.getCount();
373 const Double_t P = TMath::PoissonI(i->second.getTotal(),
y);
375 DEBUG(
"Peak at " << setw(4) << i->first <<
" [frame time] "
377 <<
FIXED(7,1) << i->second.getTotal() <<
" / "
381 if (P < Pmin && i->second.getTotal() >
y)
387 if (!(sn.size() == 1 &&
388 sn.begin()->first == 0)) {
390 WARNING(
"Module " << setw(8) << module->getID() <<
" number of peaks " << sn.size() << endl);
394 const Double_t noise = bg.getTotal() * i->second.getCount() / bg.getCount();
396 WARNING(
"Peak at " << setw(4) << i->first <<
" [frame time] "
398 <<
FIXED(7,1) << i->second.getTotal() <<
" / "
399 <<
FIXED(7,1) << noise << endl);
401 peaksPerDoms[module->getID()].push_back(sn.size());
402 peaksPerDoms[module->getID()].push_back(i->first);
403 peaksPerDoms[module->getID()].push_back(i->second.getTotal());
404 peaksPerDoms[module->getID()].push_back(noise);
407 status2[module->getID()] = (sn.size() == 1 ? OUT_SYNC :
ERROR);
414 if ((
f1.GetParameter(0) > most_OOS_peak) && !(oosDomsId[module->getID()])){
416 most_OOS_peak =
f1.GetParameter(0);
417 most_OOS_ID = module->getID();
422 if (most_OOS_ID != 0) {
423 oosDomsId[most_OOS_ID] =
true;
426 for (JDetector::const_iterator module =
detector.begin(); module !=
detector.end(); ++module) {
427 manager[module->getID()]->Reset(
"ICESM");
429 peaksPerDoms.clear();
437 NOTICE(
"Managing non-working and OOS DOMs." << endl);
443 parameters.
load(pmtFile.c_str());
449 if (i->second != IN_SYNC) {
451 NOTICE(
"Module " << setw(8) << i->first <<
" set QEs of all PMTs to zero." << endl);
459 ofstream out(pmtFile.c_str());
463 out << parameters << endl;
474 ofstream stream(peaks);
476 stream <<
"#DOM NB_PEAKS TIMESLICE_SHIFT SIGNAL NOISE\n";
478 for (map_type::const_iterator dom = peaksPerDoms.begin(); dom != peaksPerDoms.end(); ++dom) {
479 stream << dom->first <<
" ";
481 for (map_type::mapped_type::const_iterator
data = dom->second.begin();
data != dom->second.end(); ++
data) {
483 ((i%4 == 0) ? stream <<
"\n" : stream <<
" ");
495 if (i->second == IN_SYNC) {
498 if (i->second != IN_SYNC &&
499 i->second != NODATA) {
504 NOTICE(
"Number of modules out/in sync " << nout <<
'/' << nin << endl);
506 QAQC(nin <<
' ' << nout << endl);
Direct access to PMT data in detector data structure for DAQ hits.
KM3NeT DAQ constants, bit handling, etc.
Data structure for detector geometry and calibration.
Dynamic ROOT object management.
General purpose messaging.
#define DEBUG(A)
Message macros.
#define QAQC(A)
QA/QC output macro.
int qaqc
QA/QC file descriptor.
Utility class to parse command line options.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
I/O formatting auxiliaries.
Scanning of objects from a single file according a format that follows from the extension of each fil...
ROOT TTree parameter settings of various packages.
int main(int argc, char **argv)
Simple wrapper around JModuleRouter class for direct addressing of PMT data in detector data structur...
const JModule & getModule(const JDAQKeyHit &hit) const
Get module parameters.
const JPMT & getPMT(const JDAQKeyHit &hit) const
Get PMT parameters.
Auxiliary class for map of PMT parameters.
Auxialiary class to assert type conversion.
The template JSinglePointer class can be used to hold a pointer to an object.
Utility class to parse command line options.
int read(const int argc, const char *const argv[])
Parse the program's command line options.
Auxiliary class to manage set of compatible ROOT objects (e.g. histograms) using unique keys.
Template definition for direct access of elements in ROOT TChain.
int getFrameIndex() const
Get frame index.
JTriggerCounter_t next()
Increment trigger counter.
const JPolynome f1(1.0, 2.0, 3.0)
Function.
double getTime(const Hit &hit)
Get true time of hit.
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
KM3NeT DAQ data structures and auxiliaries.
double getFrameTime()
Get frame time duration.
double getMaximalTime(const double R_Hz)
Get maximal time for given rate.
static const int NUMBER_OF_PMTS
Total number of PMTs in module.
std::map< int, range_type > map_type
Auxiliary data structure for floating point format specification.
Auxiliary data structure to unify weights of acoustics data according to the number of pings per emit...
void load(const char *file_name)
Load from input file.
Auxiliary class for a type holder.
Auxiliary class to select ROOT class based on class name.
Auxiliary class to select JTreeScanner based on ROOT class name.
Auxiliary class for defining the range of iterations of objects.
Auxiliary data structure for L1 build parameters.
Auxiliary data structure for floating point format specification.