Jpp  16.0.0-rc.2
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JHydrophone.cc
Go to the documentation of this file.
1 #include <iostream>
2 #include <iomanip>
3 
4 #include "TROOT.h"
5 #include "TFile.h"
6 #include "TH1D.h"
7 
8 #include "JDB/JToAshort.hh"
9 #include "JDB/JSupport.hh"
10 
11 #include "JDetector/JDetector.hh"
13 #include "JDetector/JTripod.hh"
14 #include "JDetector/JHydrophone.hh"
15 
17 #include "JSupport/JTreeScanner.hh"
18 
19 #include "JTools/JHashMap.hh"
20 #include "JTools/JRange.hh"
21 #include "JLang/JComparator.hh"
22 #include "JLang/JComparison.hh"
23 
24 #include "JROOT/JManager.hh"
25 #include "JROOT/JRootToolkit.hh"
26 
27 #include "JAcoustics/JEmitterID.hh"
29 #include "JAcoustics/JEmitter.hh"
30 #include "JAcoustics/JReceiver.hh"
34 #include "JAcoustics/JHit.hh"
35 #include "JAcoustics/JEvent.hh"
36 #include "JAcoustics/JSupport.hh"
37 
38 #include "Jeep/JContainer.hh"
39 #include "Jeep/JPrint.hh"
40 #include "Jeep/JParser.hh"
41 #include "Jeep/JMessage.hh"
42 
43 namespace {
44 
45  /**
46  * Auxiliary data structure for organistation of histograms.
47  */
48  struct key_type :
49  public std::pair<int, int>
50  {
51  /**
52  * Constructor.
53  *
54  * \param first first value
55  * \param second second value
56  */
57  key_type(const int first,
58  const int second) :
59  std::pair<int, int>(first, second)
60  {}
61 
62 
63  /**
64  * Read key from input.
65  *
66  * \param in input stream
67  * \param key key
68  * \return input stream
69  */
70  friend inline std::istream& operator>>(std::istream& in, key_type& key)
71  {
72  in >> key.first;
73  in >> key.second;
74 
75  return in;
76  }
77 
78 
79  /**
80  * Write key to output.
81  *
82  * \param out output stream
83  * \param key key
84  * \return output stream
85  */
86  friend inline std::ostream& operator<<(std::ostream& out, const key_type& key)
87  {
88  out << key.first;
89  out << '.';
90  out << key.second;
91 
92  return out;
93  }
94 
95  /**
96  * Less=than operator.
97  *
98  * \param first first key
99  * \param second second key
100  * \return true if first key less than second key; else false
101  */
102  friend inline bool operator<(const key_type& first, const key_type& second)
103  {
104  if (first.first == second.first)
105  return first.second < second.second;
106  else
107  return first.first < second.first;
108  }
109  };
110 }
111 
112 
113 /**
114  * \file
115  *
116  * Example program to plot hydrophone data.
117  * \author mdejong
118  */
119 int main(int argc, char **argv)
120 {
121  using namespace std;
122  using namespace JPP;
123 
124  typedef JContainer< vector<JTripod> > tripods_container;
125  typedef JContainer< vector<JHydrophone> > hydrophones_container;
126 
128  string detectorFile;
129  JLimit_t& numberOfEvents = inputFile.getLimit();
130  string outputFile;
131  JSoundVelocity V = getSoundVelocity; // default sound velocity
132  tripods_container tripods; // tripods
133  hydrophones_container hydrophones; // hydrophones
134  double Q;
135  int debug;
136 
137  try {
138 
139  JParser<> zap("Example program to plot hydrophone data.");
140 
141  zap['f'] = make_field(inputFile);
142  zap['n'] = make_field(numberOfEvents) = JLimit::max();
143  zap['o'] = make_field(outputFile) = "hydrophone.root";
144  zap['a'] = make_field(detectorFile);
145  zap['V'] = make_field(V, "sound velocity") = JPARSER::initialised();
146  zap['T'] = make_field(tripods, "tripod data");
147  zap['H'] = make_field(hydrophones, "hydrophone data") = JPARSER::initialised();
148  zap['W'] = make_field(getEmitterID, "waveform identification data") = JPARSER::initialised();
149  zap['Q'] = make_field(Q, "quality") = 0.0;
150  zap['d'] = make_field(debug) = 2;
151 
152  zap(argc, argv);
153  }
154  catch(const exception &error) {
155  FATAL(error.what() << endl);
156  }
157 
158 
160 
161  try {
162  load(detectorFile, detector);
163  }
164  catch(const JException& error) {
165  FATAL(error);
166  }
167 
168  V.set(detector.getUTMZ());
169 
170  JHashMap<int, JReceiver> receivers;
171  JHashMap<int, JEmitter> emitters;
172 
173  for (hydrophones_container::const_iterator i = hydrophones.begin(); i != hydrophones.end(); ++i) {
174 
175  for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) {
176 
177  if (i->getLocation() == module->getLocation()) {
178 
179  receivers[module->getID()] = JReceiver(module->getID(),
180  module->getPosition() + i->getPosition(),
181  module->getT0() * 1.0e-9);
182  break;
183  }
184  }
185  }
186 
187  for (tripods_container::const_iterator i = tripods.begin(); i != tripods.end(); ++i) {
188  emitters[i->getID()] = JEmitter(i->getID(),
189  i->getUTMPosition() - detector.getUTMPosition());
190  }
191 
192 
193  const JRange<double> T0(-0.05, +0.05);
194  const JRange<double> T1(-0.01, +0.01);
195 
196  TH1D h0(MAKE_CSTRING("Q0 " << T0), NULL, 100, 0.0, 8.0);
197  TH1D h1(MAKE_CSTRING("Q1 " << T1), NULL, 100, 0.0, 8.0);
198 
199  JManager<key_type, TH1D> H1(new TH1D("[%]", NULL, 50000, -5.0e-2, +5.0e-2));
200 
201 
202  while (inputFile.hasNext()) {
203 
204  STATUS("event: " << setw(10) << inputFile.getCounter() << '\r'); DEBUG(endl);
205 
206  const JEvent* evt = inputFile.next();
207 
208  if (!evt->empty() && emitters.has(evt->getID())) {
209 
211 
212  for (JEvent::const_iterator hit = evt->begin(); hit != evt->end(); ++hit) {
213  if (receivers.has(hit->getID())) {
214  buffer[hit->getID()].push_back(*hit);
215  }
216  }
217 
218  for (map<int, vector<JTransmission> >::iterator i = buffer.begin(); i != buffer.end(); ++i) {
219 
220  sort(i->second.begin(), i->second.end(), make_comparator(&JTransmission::getQ, JComparison::gt()));
221 
222  vector<JTransmission>::const_iterator hit = i->second.begin();
223 
224  if (hit->getQ() >= Q) {
225 
226  const JPosition3D& p1 = emitters [evt->getID()].getPosition();
227  const JPosition3D& p2 = receivers[hit->getID()].getPosition();
228 
229  const double D = p2.getDistance(p1);
230  const double Vi = V.getInverseVelocity(D, p1.getZ(), p2.getZ());
231 
232  const double t0 = evt->begin()->getToE() + D * Vi;
233  const double t1 = hit->getToA() - t0;
234 
235  H1[key_type(hit->getID(), evt->getID())]->Fill(t1);
236 
237  const double Q = log10(hit->getQ());
238 
239  if (T0(t1)) { h0.Fill(Q); }
240  if (T1(t1)) { h1.Fill(Q); }
241  }
242  }
243  }
244  }
245  STATUS(endl);
246 
247  TFile out(outputFile.c_str(), "recreate");
248 
249  out << h0 << h1 << H1;
250 
251  out.Write();
252  out.Close();
253 }
Utility class to parse command line options.
Definition: JParser.hh:1500
Acoustic hit.
General exception.
Definition: JException.hh:23
Acoustic receiver.
Definition: JReceiver.hh:27
Q(UTCMax_s-UTCMin_s)-livetime_s
int main(int argc, char *argv[])
Definition: Main.cc:15
TPaveText * p1
Greater than.
Definition: JComparison.hh:73
JComparator< JResult_t T::*, JComparison::lt > make_comparator(JResult_t T::*member)
Helper method to create comparator between values of data member.
Definition: JComparator.hh:185
Sound velocity.
General purpose class for hash map of unique keys.
Definition: JHashMap.hh:72
bool operator<(const Head &first, const Head &second)
Less than operator.
Definition: JHead.hh:1678
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
#define MAKE_CSTRING(A)
Make C-string.
Definition: JPrint.hh:151
then for HISTOGRAM in h0 h1
Definition: JMatrixNZ.sh:71
Acoustic event.
Empty structure for specification of parser element that is initialised (i.e. does not require input)...
Definition: JParser.hh:66
Dynamic ROOT object management.
V(JDAQEvent-JTriggerReprocessor)*1.0/(JDAQEvent+1.0e-10)
string outputFile
Acoustic emitter.
double getDistance(const JVector3D &pos) const
Get distance to point.
Definition: JVector3D.hh:270
Data structure for detector geometry and calibration.
then echo The file $DIR KM3NeT_00000001_00000000 root already please rename or remove it first
Data structure for hydrophone.
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:41
static const JSoundVelocity getSoundVelocity(1541.0,-17.0e-3,-2000.0)
Function object for velocity of sound.
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
Acoustics support kit.
Acoustic emitter.
Definition: JEmitter.hh:27
Acoustics toolkit.
Auxiliary wrapper for I/O of container with optional comment (see JComment).
Definition: JContainer.hh:39
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1961
set_variable E_E log10(E_{fit}/E_{#mu})"
void Write(TDirectory &out, const bool wm=false)
Write objects to file.
Definition: JManager.hh:295
Emitter identification.
ROOT TTree parameter settings.
int debug
debug level
Definition: JSirene.cc:63
double getQ(const double D_m, const double f_kHz, const double d_m)
Get relative quality for given frequency at given distance.
General purpose messaging.
Implementation for depth dependend velocity of sound.
#define FATAL(A)
Definition: JMessage.hh:67
Scanning of objects from multiple files according a format that follows from the extension of each fi...
Acoustic receiver.
p2
Definition: module-Z:fit.sh:74
std::istream & operator>>(std::istream &in, JAANET::JHead &header)
Read header from input.
Definition: JHead.hh:1693
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
Auxiliary class to define a range between two values.
General purpose class for object reading from a list of file names.
Utility class to parse command line options.
Acoustic transceiver.
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.
std::ostream & operator<<(std::ostream &stream, const CLBCommonHeader &header)
Data structure for position in three dimensions.
Definition: JPosition3D.hh:36
const JLimit & getLimit() const
Get limit.
Definition: JLimit.hh:73
static JEmitterID getEmitterID
Function object for emitter identification.
Definition: JEmitterID.hh:119
do set_variable DETECTOR_TXT $WORKDIR detector
int getID() const
Get identifier.
Data structure for tripod.
then fatal Not enough tripods
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
do echo Generating $dir eval D
Definition: JDrawLED.sh:53
double getZ() const
Get z position.
Definition: JVector3D.hh:115
Container I/O.