Jpp  15.0.3
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JAcousticsMonitor.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 "JAcoustics/JEvent.hh"
28 #include "JAcoustics/JSupport.hh"
30 
31 #include "Jeep/JContainer.hh"
32 #include "Jeep/JPrint.hh"
33 #include "Jeep/JParser.hh"
34 #include "Jeep/JMessage.hh"
35 
36 
37 /**
38  * \file
39  *
40  * Example program to monitor acoustic events.
41  * \author mdejong
42  */
43 int main(int argc, char **argv)
44 {
45  using namespace std;
46  using namespace JPP;
47 
48  JMultipleFileScanner<> inputFile;
49  JLimit_t& numberOfEvents = inputFile.getLimit();
50  string outputFile;
51  string detectorFile;
53  int debug;
54 
55  try {
56 
57  JParser<> zap("Example program to monitor acoustic events.");
58 
59  zap['f'] = make_field(inputFile, "output of JAcousticEventBuilder[.sh]");
60  zap['n'] = make_field(numberOfEvents) = JLimit::max();
61  zap['o'] = make_field(outputFile) = "monitor.root";
62  zap['a'] = make_field(detectorFile);
63  zap['@'] = make_field(parameters, "trigger parameters") = JPARSER::initialised();
64  zap['d'] = make_field(debug) = 2;
65 
66  zap(argc, argv);
67  }
68  catch(const exception &error) {
69  FATAL(error.what() << endl);
70  }
71 
73 
74  try {
75  load(detectorFile, detector);
76  }
77  catch(const JException& error) {
78  FATAL(error);
79  }
80 
81  JHashMap<int, JLocation> receivers;
82 
83  for (JDetector::const_iterator i = detector.begin(); i != detector.end(); ++i) {
84  receivers[i->getID()] = i->getLocation();
85  }
86 
87 
88  const JHashCollection<int> string(make_array(detector.begin(), detector.end(), &JModule::getString));
89  const JRange <int> floor (make_array(detector.begin(), detector.end(), &JModule::getFloor));
90 
91 
92  JManager<int, TH1D> HA(new TH1D("H[%].size", NULL, detector.size() + 1, -0.5, detector.size() + 0.5));
93  JManager<int, TH1D> HB(new TH1D("H[%].overlays", NULL, 101, -0.5, 100.5));
94  JManager<int, TH1D> HC(new TH1D("H[%].time", NULL, 800, 0.0, 4.0));
95  JManager<int, TH1D> HD(new TH1D("H[%].rms", NULL, 200, 0.0, 5.0e-3));
96  JManager<int, TH1D> HE(new TH1D("H[%].deviation", NULL, 200, 0.0, 50.0e-3));
97  JManager<int, TH1D> HQ(new TH1D("H[%].quality", NULL, 200, 0.0, 8.0));
98  JManager<int, TH1D> H1(new TH1D("H[%].multiplicity", NULL, 101, -0.5, 100.5));
99  JManager<int, TH2D> H2(new TH2D("H[%].event", NULL,
100  string.size(), -0.5, string.size() - 0.5,
101  floor.getUpperLimit() + 1, - 0.5, floor.getUpperLimit() + 0.5));
102  JManager<int, TH2D> H3((TH2D*) H2->Clone("H[%].doubles"));
103 
104  JManager<int, TGraph> GA(new TGraph(), "G[%].size");
105  JManager<int, TGraph> GD(new TGraph(), "G[%].rms");
106  JManager<int, TGraph> GE(new TGraph(), "G[%].deviation");
107  JManager<int, TGraph> GQ(new TGraph(), "G[%].quality");
108 
109  for (Int_t i = 1; i <= H2->GetXaxis()->GetNbins(); ++i) {
110  H2->GetXaxis()->SetBinLabel(i, MAKE_CSTRING(string.at(i-1)));
111  H3->GetXaxis()->SetBinLabel(i, MAKE_CSTRING(string.at(i-1)));
112  }
113 
114  for (Int_t i = 1; i <= H2->GetYaxis()->GetNbins(); ++i) {
115  H2->GetYaxis()->SetBinLabel(i, MAKE_CSTRING(i-1));
116  H3->GetYaxis()->SetBinLabel(i, MAKE_CSTRING(i-1));
117  }
118 
120 
121  JTreeScanner_t in(inputFile);
122 
124 
125  for (JTreeScanner_t::iterator event = in.begin(); event != in.end(); ++event) {
126 
127  STATUS("event " << setw(8) << event->getCounter() << '\r'); DEBUG(endl);
128 
129  JEvent evt = *event;
130 
131  evt.remove(parameters.numberOfOutliers, parameters.stdev);
132 
133  const double t1 = evt.begin()->getToE();
134 
135  JQuantile QD, QQ;
136 
137  for (JEvent::const_iterator i = evt.begin(); i != evt.end(); ++i) {
138  QD.put(i->getToE() - t1);
139  QQ.put(log10(i->getQ()));
140  }
141 
142  HA[evt.getID()]->Fill((double) evt.size());
143  HB[evt.getID()]->Fill((double) evt.getOverlays());
144 
145  if (toe.has(evt.getID())) {
146  HC[evt.getID()]->Fill(log10(evt.begin()->getToE() - toe[evt.getID()]));
147  }
148 
149  HD[evt.getID()]->Fill(QD.getSTDev());
150  HE[evt.getID()]->Fill(QD.getDeviation());
151 
152  AddPoint(GA[evt.getID()], t1, evt.size());
153  AddPoint(GD[evt.getID()], t1, QD.getSTDev());
154  AddPoint(GE[evt.getID()], t1, QD.getDeviation());
155  AddPoint(GQ[evt.getID()], t1, QQ.getMean());
156 
157  toe[evt.getID()] = evt.begin()->getToE();
158 
159  TH1D* h1 = H1[evt.getID()];
160  TH1D* hq = HQ[evt.getID()];
161  TH2D* h2 = H2[evt.getID()];
162  TH2D* h3 = H3[evt.getID()];
163 
164  map<int, set<double> > buffer;
165 
166  for (JEvent::const_iterator i = evt.begin(); i != evt.end(); ++i) {
167  buffer[i->getID()].insert(i->getQ());
168  }
169 
170  for (map<int, set<double> >::const_iterator i = buffer.begin(); i != buffer.end(); ++i) {
171 
172  h1->Fill((double) i->second.size());
173 
174  hq->Fill(log10(*(i->second.rbegin())));
175 
176  if (receivers.has(i->first)) {
177 
178  const JLocation& location = receivers[i->first];
179 
180  h2->Fill((double) string.getIndex(location.getString()), (double) location.getFloor(), 1.0);
181 
182  if (i->second.size() > 1u) {
183  h3->Fill((double) string.getIndex(location.getString()), (double) location.getFloor(), 1.0);
184  }
185  }
186  }
187  }
188  STATUS(endl);
189 
190 
191  TFile out(outputFile.c_str(), "recreate");
192 
193  out << HA << HB << HC << HD << HE << HQ << H1 << H2 << H3;
194  out << GA << GD << GE << GQ;
195 
196  out.Write();
197  out.Close();
198 }
Utility class to parse command line options.
Definition: JParser.hh:1500
General exception.
Definition: JException.hh:23
int getOverlays() const
Get overlays.
int main(int argc, char *argv[])
Definition: Main.cc:15
void AddPoint(TGraph *g1, const Double_t x, const Double_t y)
Add point to TGraph.
int getFloor() const
Get floor number.
Definition: JLocation.hh:145
General purpose class for hash map of unique keys.
Definition: JHashMap.hh:72
double getSTDev() const
Get standard deviation.
Definition: JQuantile.hh:266
Auxiliary data structure for running average, standard deviation and quantiles.
Definition: JQuantile.hh:43
General purpose class for hash map of unique elements.
#define STATUS(A)
Definition: JMessage.hh:63
ROOT TTree parameter settings.
Detector data structure.
Definition: JDetector.hh:89
*fatal Wrong number of arguments esac JCookie sh typeset Z DETECTOR typeset Z SOURCE_RUN typeset Z TARGET_RUN set_variable PARAMETERS_FILE $WORKDIR parameters
Definition: diff-Tuna.sh:38
#define MAKE_CSTRING(A)
Make C-string.
Definition: JPrint.hh:151
then for HISTOGRAM in h0 h1
Definition: JMatrixNZ.sh:71
Acoustic event.
General purpose class for a hash collection of unique elements.
Empty structure for specification of parser element that is initialised (i.e. does not require input)...
Definition: JParser.hh:66
Dynamic ROOT object management.
string outputFile
Template definition for direct access of elements in ROOT TChain.
Definition: JTreeScanner.hh:91
Data structure for detector geometry and calibration.
double getMean() const
Get mean value.
Definition: JQuantile.hh:252
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:41
Auxiliary class to manage set of compatible ROOT objects (e.g. histograms) using unique keys...
Definition: JManager.hh:43
I/O formatting auxiliaries.
Detector file.
Definition: JHead.hh:224
Logical location of module.
Definition: JLocation.hh:37
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1961
const array_type< JValue_t > & make_array(const JValue_t(&array)[N])
Method to create array of values.
Definition: JVectorize.hh:54
T getUpperLimit() const
Get upper limit.
Definition: JRange.hh:213
void put(const double x, const double w=1.0)
Put value.
Definition: JQuantile.hh:133
void Write(TDirectory &out, const bool wm=false)
Write objects to file.
Definition: JManager.hh:295
int getIndex()
Get index for user I/O manipulation.
Definition: JManip.hh:26
int debug
debug level
Definition: JSirene.cc:63
General purpose messaging.
#define FATAL(A)
Definition: JMessage.hh:67
Scanning of objects from multiple files according a format that follows from the extension of each fi...
Base class for JTreeScanner.
Definition: JTreeScanner.hh:52
double getDeviation(const bool relative=true) const
Get maximal deviation from average.
Definition: JQuantile.hh:281
int getString() const
Get string number.
Definition: JLocation.hh:134
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
Acoustic trigger parameters.
void remove(const size_t ns, const double stdev)
Remove outliers.
Auxiliary class to define a range between two values.
Utility class to parse command line options.
then usage $script< input_file >< detector_file > fi set_variable OUTPUT_DIR set_variable SELECTOR JDAQTimesliceL1 set_variable DEBUG case set_variable DEBUG
Acoustic event.
const JLimit & getLimit() const
Get limit.
Definition: JLimit.hh:73
do set_variable DETECTOR_TXT $WORKDIR detector
int getID() const
Get identifier.
double u[N+1]
Definition: JPolint.hh:739
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable INPUT_FILE $argv[2] eval JPrintDetector a $DETECTOR O IDENTIFIER eval JPrintDetector a $DETECTOR O SUMMARY source JAcoustics sh $DETECTOR_ID CHECK_EXIT_CODE typeset A TRIPODS get_tripods $WORKDIR tripod txt TRIPODS for EMITTER in
Definition: JCanberra.sh:42
Container I/O.