Jpp  debug
the software that should make you happy
JAcousticsMonitor_short.cc
Go to the documentation of this file.
1 #include <iostream>
2 #include <iomanip>
3 #include <limits>
4 #include <map>
5 
6 #include "TROOT.h"
7 #include "TFile.h"
8 #include "TH1D.h"
9 #include "TH2D.h"
10 #include "TGraph.h"
11 #include "TGraphErrors.h"
12 
14 #include "JTools/JHashMap.hh"
15 #include "JTools/JRange.hh"
16 #include "JTools/JQuantile.hh"
17 
18 #include "JDetector/JDetector.hh"
20 
22 #include "JSupport/JTreeScanner.hh"
23 
24 #include "JROOT/JRootToolkit.hh"
25 #include "JROOT/JManager.hh"
26 
27 #include "JLang/JPredicate.hh"
28 
29 #include "JAcoustics/JEvent.hh"
30 #include "JAcoustics/JSupport.hh"
31 
32 #include "Jeep/JContainer.hh"
33 #include "Jeep/JPrint.hh"
34 #include "Jeep/JParser.hh"
35 #include "Jeep/JMessage.hh"
36 
37 
38 /**
39  * \file
40  *
41  * Example program to monitor acoustic events rate per emitter and receiver and run.
42  * Output: histogram of the acoustic event rate per emitter and run.
43  * \author mdejong, cgatius
44  */
45 int main(int argc, char **argv)
46 {
47  using namespace std;
48  using namespace JPP;
49 
50  JMultipleFileScanner<> inputFile;
51  JLimit_t& numberOfEvents = inputFile.getLimit();
52  string outputFile;
53  string detectorFile;
54  int debug;
55  double lifetime_s;
56 
57  try {
58 
59  JParser<> zap("Example program to monitor acoustic events.");
60 
61  zap['f'] = make_field(inputFile, "output of JAcousticEventBuilder[.sh]");
62  zap['n'] = make_field(numberOfEvents) = JLimit::max();
63  zap['o'] = make_field(outputFile) = "monitor.root";
64  zap['a'] = make_field(detectorFile);
65  zap['d'] = make_field(debug) = 2;
66  zap['l'] = make_field(lifetime_s, "run lifetime in seconds");
67  zap(argc, argv);
68  }
69  catch(const exception &error) {
70  FATAL(error.what() << endl);
71  }
72 
74 
75  try {
76  load(detectorFile, detector);
77  }
78  catch(const JException& error) {
79  FATAL(error);
80  }
81 
82  JHashMap<int, JLocation> receivers;
83 
84  for (JDetector::const_iterator i = detector.begin(); i != detector.end(); ++i) {
85 
86  receivers[i->getID()] = i->getLocation();
87  }
88 
89  const JHashCollection<int> string(make_array(detector.begin(), detector.end(), &JModule::getString));
90  const JRange <int> floor (make_array(detector.begin(), detector.end(), &JModule::getFloor));
91 
92  JManager<int, TH2D> H2(new TH2D("H[%].event", NULL,
93  string.size(), -0.5, string.size() - 0.5,
94  floor.getUpperLimit() + 1, - 0.5, floor.getUpperLimit() + 0.5));
95 
96  for (Int_t i = 1; i <= H2->GetXaxis()->GetNbins(); ++i) {
97  H2->GetXaxis()->SetBinLabel(i, MAKE_CSTRING(string.at(i-1)));
98  }
99 
100  for (Int_t i = 1; i <= H2->GetYaxis()->GetNbins(); ++i) {
101  H2->GetYaxis()->SetBinLabel(i, MAKE_CSTRING(i-1));
102  }
103 
105 
106  JTreeScanner_t in(inputFile);
107 
108  for (JTreeScanner_t::iterator event = in.begin(); event != in.end(); ++event) {
109 
110  JQuantile QD("", true), QQ;
111 
112  TH2D* h2 = H2[event->getID()];
113 
114  map<int, set<double> > buffer;
115 
116  for (JEvent::const_iterator i = event->begin(); i != event->end(); ++i) {
117  buffer[i->getID()].insert(i->getQ());
118  }
119 
120  for (map<int, set<double> >::const_iterator i = buffer.begin(); i != buffer.end(); ++i) {
121 
122  if (receivers.has(i->first)) {
123 
124  const JLocation& location = receivers[i->first];
125 
126  h2->Fill((double) string.getIndex(location.getString()), (double) location.getFloor(), 1.0/lifetime_s);
127  }
128  }
129  }
130  STATUS(endl);
131 
132  TFile out(outputFile.c_str(), "recreate");
133 
134  out << H2;
135 
136  out.Write();
137  out.Close();
138 }
int main(int argc, char **argv)
Acoustic event.
ROOT TTree parameter settings.
Container I/O.
string outputFile
Data structure for detector geometry and calibration.
General purpose class for a hash collection of unique elements.
General purpose class for hash map of unique elements.
Dynamic ROOT object management.
General purpose messaging.
#define STATUS(A)
Definition: JMessage.hh:63
#define FATAL(A)
Definition: JMessage.hh:67
int debug
debug level
Definition: JSirene.cc:69
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:2158
I/O formatting auxiliaries.
#define MAKE_CSTRING(A)
Make C-string.
Definition: JPrint.hh:136
Auxiliary class to define a range between two values.
Detector data structure.
Definition: JDetector.hh:96
Logical location of module.
Definition: JLocation.hh:39
int getFloor() const
Get floor number.
Definition: JLocation.hh:145
int getString() const
Get string number.
Definition: JLocation.hh:134
General exception.
Definition: JException.hh:24
Utility class to parse command line options.
Definition: JParser.hh:1714
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
Base class for JTreeScanner.
Definition: JTreeScanner.hh:56
Template definition for direct access of elements in ROOT TChain.
bool has(const T &value) const
Test whether given value is present.
T getUpperLimit() const
Get upper limit.
Definition: JRange.hh:213
JMODEL::JString getString(const JFit &fit)
Get model parameters of string.
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
const array_type< JValue_t > & make_array(const JValue_t(&array)[N])
Method to create array of values.
Definition: JVectorize.hh:54
int getIndex()
Get index for user I/O manipulation.
Definition: JManip.hh:26
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Definition: JSTDTypes.hh:14
Detector file.
Definition: JHead.hh:227
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:45
General purpose class for hash map of unique keys.
Definition: JHashMap.hh:75
Auxiliary data structure for running average, standard deviation and quantiles.
Definition: JQuantile.hh:46