Jpp  16.0.0-rc.1
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JMuonVeto.cc
Go to the documentation of this file.
1 
2 #include <string>
3 #include <iostream>
4 #include <iomanip>
5 #include <vector>
6 #include <set>
7 
11 #include "JDAQ/JDAQEventIO.hh"
12 #include "JDAQ/JDAQTimesliceIO.hh"
14 
15 #include "JDetector/JDetector.hh"
18 
19 #include "JTrigger/JHit.hh"
20 #include "JTrigger/JFrame.hh"
21 #include "JTrigger/JTimeslice.hh"
22 #include "JTrigger/JHitL0.hh"
23 #include "JTrigger/JHitL1.hh"
24 #include "JTrigger/JHitR1.hh"
25 #include "JTrigger/JBuildL0.hh"
27 
31 #include "JSupport/JSupport.hh"
32 #include "JSupport/JMeta.hh"
33 
34 #include "JFit/JLine1Z.hh"
35 #include "JFit/JFitToolkit.hh"
36 #include "JFit/JTimeRange.hh"
37 
38 #include "JReconstruction/JEvt.hh"
40 
41 #include "Jeep/JParser.hh"
42 #include "Jeep/JMessage.hh"
43 
44 
45 /**
46  * \file
47  *
48  * Program to determine veto of muon trajectory.
49  * \author mdejong
50  */
51 int main(int argc, char **argv)
52 {
53  using namespace std;
54  using namespace JPP;
55  using namespace KM3NETDAQ;
56 
57  typedef vector<JTimeRange> JVeto_t;
58  typedef JParallelFileScanner< JTypeList<JDAQEvent, JEvt> > JParallelFileScanner_t;
59  typedef JParallelFileScanner_t::multi_pointer_type multi_pointer_type;
61 
62  JParallelFileScanner_t inputFile;
64  JLimit_t& numberOfEvents = inputFile.getLimit();
65  string detectorFile;
66  double roadWidth_m;
67  double R_Hz;
68  size_t numberOfPrefits;
69  JVeto_t veto;
70  bool reprocess;
71  int debug;
72 
73  try {
74 
75  JParser<> zap("Program to determine veto of muon trajectory.");
76 
77  zap['f'] = make_field(inputFile);
78  zap['o'] = make_field(outputFile) = "veto.root";
79  zap['a'] = make_field(detectorFile);
80  zap['n'] = make_field(numberOfEvents) = JLimit::max();
81  zap['R'] = make_field(roadWidth_m) = numeric_limits<double>::max();
82  zap['B'] = make_field(R_Hz) = 6.0e3;
83  zap['N'] = make_field(numberOfPrefits) = 1;
84  zap['V'] = make_field(veto);
85  zap['r'] = make_field(reprocess);
86  zap['d'] = make_field(debug) = 1;
87 
88  zap(argc, argv);
89  }
90  catch(const exception& error) {
91  FATAL(error.what() << endl);
92  }
93 
94 
95  cout.tie(&cerr);
96 
98 
99  try {
100  load(detectorFile, detector);
101  }
102  catch(const JException& error) {
103  FATAL(error);
104  }
105 
106  const JModuleRouter router(detector);
107 
108  typedef vector<JHitL0> JDataL0_t;
109  const JBuildL0<JHitL0> buildL0;
110 
111  double VS = 0.0; // veto [s]
112 
113  for (vector<JTimeRange>::const_iterator vs = veto.begin(); vs != veto.end(); ++vs) {
114  VS += vs->getLength() * 1.0e-9;
115  }
116 
117 
118  outputFile.open();
119 
120  outputFile.put(JMeta(argc, argv));
121 
122  while (inputFile.hasNext()) {
123 
124  STATUS("event: " << setw(10) << inputFile.getCounter() << '\r'); DEBUG(endl);
125 
126  multi_pointer_type ps = inputFile.next();
127 
128  const JDAQEvent* tev = ps;
129  const JEvt* in = ps;
130 
131  JEvt cp = *in;
132  JEvt out;
133 
134  cp.select(numberOfPrefits, qualitySorter);
135 
136 
137  JDataL0_t dataL0;
138 
139  buildL0(*tev, router, true, back_inserter(dataL0));
140 
141 
142  for (JEvt::iterator track = cp.begin(); track != cp.end(); ++track) {
143 
144  const JRotation3D R (getDirection(*track));
145  const JLine1Z tz(getPosition (*track).rotate(R), track->getT());
146 
147  // hit selection based on start value
148 
150 
151  for (JDataL0_t::const_iterator i = dataL0.begin(); i != dataL0.end(); ++i) {
152 
153  JHitL0 hit(*i);
154 
155  hit.rotate(R);
156 
157  bool is_veto = false;
158 
159  if (tz.getDistance(hit) <= roadWidth_m) {
160 
161  const double t1 = hit.getT() - tz.getT(hit);
162 
163  for (vector<JTimeRange>::const_iterator vs = veto.begin(); vs != veto.end() && !is_veto; ++vs) {
164  is_veto = (*vs)(t1);
165  }
166  }
167 
168  if (is_veto) {
169  top.insert(hit.getPMTIdentifier());
170  }
171  }
172 
173 
174  double npe = 0.0;
175  int count = 0;
176 
177  for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) {
178 
179  JPosition3D pos(module->getPosition());
180 
181  pos.rotate(R);
182 
183  if (tz.getDistance(pos) <= roadWidth_m) {
184 
185  for (unsigned int i = 0; i != module->size(); ++i) {
186 
187  const JDAQPMTIdentifier id(module->getID(), i);
188 
189  npe += R_Hz * VS;
190  count += top.count(id);
191  }
192  }
193  }
194 
195  out.push_back(JFit(*track).add(JMUONVETO));
196 
197  out.rbegin()->setW(track->getW());
198  out.rbegin()->setW(JVETO_NPE, npe);
199  out.rbegin()->setW(JVETO_NUMBER_OF_HITS, count);
200  }
201 
202  // apply default sorter
203 
204  sort(out.begin(), out.end(), qualitySorter);
205 
206  outputFile.put(out);
207  }
208  STATUS(endl);
209 
211 
212  io >> outputFile;
213 
214  outputFile.close();
215 }
Auxiliary methods to evaluate Poisson probabilities and chi2.
Auxiliary class for ROOT I/O of application specific meta data.
Definition: JMeta.hh:70
double getT() const
Get calibrated time of hit.
Object writing to file.
Utility class to parse command line options.
Definition: JParser.hh:1500
General exception.
Definition: JException.hh:23
int main(int argc, char *argv[])
Definition: Main.cc:15
ROOT TTree parameter settings of various packages.
static const int JVETO_NPE
number of photo-electrons from JVeto.cc
Template specialisation of L0 builder for JHitL0 data type.
Definition: JBuildL0.hh:102
const JDAQPMTIdentifier & getPMTIdentifier() const
Get PMT identifier.
#define STATUS(A)
Definition: JMessage.hh:63
Detector data structure.
Definition: JDetector.hh:89
Recording of objects on file according a format that follows from the file name extension.
Router for direct addressing of module data in detector data structure.
Rotation matrix.
Definition: JRotation3D.hh:111
General purpose class for parallel reading of objects from a single file or multiple files...
static const int JMUONVETO
Basic data structure for time and time over threshold information of hit.
JFit & add(const int type)
Add event to history.
string outputFile
Data structure for detector geometry and calibration.
bool qualitySorter(const JFit &first, const JFit &second)
Comparison of fit results.
Basic data structure for L0 hit.
Type list.
Definition: JTypeList.hh:22
Scanning of objects from a single file according a format that follows from the extension of each fil...
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:41
JAxis3D & rotate(const JRotation3D &R)
Rotate axis.
Definition: JAxis3D.hh:225
Detector file.
Definition: JHead.hh:224
JDirection3D getDirection(const Vec &dir)
Get direction.
Acoustic event fit.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1961
Parallel scanning of objects from a single file or multiple files according a format that follows fro...
ROOT I/O of application specific meta data.
JPosition3D getPosition(const Vec &pos)
Get position.
static const int JVETO_NUMBER_OF_HITS
number of hits from JVeto.cc
int debug
debug level
Definition: JSirene.cc:63
then usage $script[distance] fi case set_variable R
Definition: JDrawLED.sh:43
General purpose messaging.
#define FATAL(A)
Definition: JMessage.hh:67
Direct access to module in detector data structure.
Reduced data structure for L1 hit.
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
std::vector< int > count
Definition: JAlgorithm.hh:180
Utility class to parse command line options.
Data structure for L0 hit.
Definition: JHitL0.hh:27
then usage $script< input_file >< detector_file > fi set_variable OUTPUT_DIR set_variable SELECTOR JDAQTimesliceL1 set_variable DEBUG case set_variable DEBUG
Data structure for fit of straight line paralel to z-axis.
Definition: JLine1Z.hh:27
then cp
Object reading from a list of files.
Data structure for position in three dimensions.
Definition: JPosition3D.hh:36
const JLimit & getLimit() const
Get limit.
Definition: JLimit.hh:73
do set_variable DETECTOR_TXT $WORKDIR detector
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
JPosition3D & rotate(const JRotation3D &R)
Rotate.
Definition: JPosition3D.hh:186
Basic data structure for L1 hit.