Jpp  18.2.0-rc.1
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JKexing.cc
Go to the documentation of this file.
1 #include <string>
2 #include <iostream>
3 #include <fstream>
4 #include <iomanip>
5 #include <vector>
6 
7 #include "TROOT.h"
8 #include "TFile.h"
9 #include "TH1D.h"
10 #include "TF1.h"
11 
12 #include "JSupernova.hh"
13 
15 #include "JDAQ/JDAQEvaluator.hh"
16 
17 #include "JDetector/JDetector.hh"
19 
20 #include "JROOT/JManager.hh"
21 
23 #include "JTools/JQuantile.hh"
24 
26 #include "JSupport/JTreeScanner.hh"
28 #include "JSupport/JSupport.hh"
29 #include "JSupport/JMeta.hh"
30 
32 #include "JMath/JMathToolkit.hh"
33 
34 #include "Jeep/JPrint.hh"
35 #include "Jeep/JParser.hh"
36 #include "Jeep/JMessage.hh"
37 
39 
40 using namespace std;
41 using namespace JPP;
42 using namespace KM3NETDAQ;
43 using namespace JSUPERNOVA;
44 
45 /**
46  * \author mlincett
47  * \file
48  * Example application to test supernova trigger on run files.
49  */
50 int main(int argc, char **argv)
51 {
52 
53  typedef JRange<int> JRange_t;
54 
55  JMultipleFileScanner<> inputFile;
56  JLimit_t& numberOfEvents = inputFile.getLimit();
57  string outputFile;
58  string detectorFile;
59  double TMax_ns;
60  double TVeto_ns;
61  JROOTClassSelector selector;
62  JDAQHitToTSelector totSelector_ns;
63  int debug;
64  int preTriggerThreshold;
65  JRange_t M;
66 
67 
68  // Check input parameters
69 
70  try {
71 
72  JParser<> zap("Example program test supernova triggers.");
73 
74  zap['f'] = make_field(inputFile);
75  zap['o'] = make_field(outputFile) = "kexing.root";
76  zap['n'] = make_field(numberOfEvents) = JLimit::max();
77  zap['a'] = make_field(detectorFile);
78  zap['T'] = make_field(TMax_ns) = 10.0;
79  zap['V'] = make_field(TVeto_ns) = 1000.0;
80  zap['C'] = make_field(selector) = getROOTClassSelection<JDAQTimesliceTypes_t>();
81  zap['d'] = make_field(debug) = 1;
82  zap['M'] = make_field(M) = JRange<int>(6,10);
83  zap['t'] = make_field(preTriggerThreshold) = 4;
84  zap['t'] = make_field(totSelector_ns) = JDAQHitToTSelector(4, 255);
85 
86  zap(argc, argv);
87  }
88  catch(const exception &error) {
89  FATAL(error.what() << endl);
90  }
91 
92 
94 
95  try {
96  load(detectorFile, detector);
97  }
98  catch(const JException& error) {
99  FATAL(error);
100  }
101 
102  // Configure input streams
103 
105 
107 
108  JTreeScannerInterface<JDAQTimeslice>* pts = zmap[selector];
109 
110  pts->configure(inputFile);
111 
112  // Configure routers
113 
114  const JModuleRouter moduleRouter(detector);
115 
116  const JDAQHitRouter hitRouter(detector);
117 
118  // -----------------------------------
119  // STEP 1: building vetoes from events
120  // -----------------------------------
121 
122  TH1D* h_vtr = new TH1D("VetoTimeRange","VetoTimeRange", 10000, 0, 10000);
123 
125 
126  map<int, JVetoSet> triggeredEvents;
127 
128  for (; evIn.hasNext(); ) {
129 
130  STATUS("event: " << setw(10) << evIn.getCounter() << '\r'); DEBUG(endl);
131 
132  JDAQEvent* event = evIn.next();
133 
134  JVeto vp(*event, hitRouter);
135 
136  triggeredEvents[event->getFrameIndex()].push_back(vp);
137 
138  h_vtr->Fill(vp.getLength());
139 
140  }
141 
142  STATUS(triggeredEvents.size() << " JDAQEvents loaded in veto buffer." << endl);
143 
144  //--------------------------------
145  // STEP 2: timeslice processing
146  // -------------------------------
147 
148  // 0 -> no filter
149  // 1 -> count once per track
150  // 2 -> suppress track
151  // 3 -> suppress based on trigger veto
152 
153  // Output data structures
154 
155  typedef JManager<int, TH1D> JManager_t;
156 
157  JManager_t SNT(new TH1D("SNT_F%", NULL, 100, 0.0, 100));
158 
159  JManager_t MUL(new TH1D("MUL_F%", NULL, 1 + 31, -0.5, 31 + 1 - 0.5));
160 
161  const int nStages = 5;
162  vector<vector<int> > trgHistory(nStages);
163 
164  // Loop
165 
166  counter_type counter = 0;
167 
168  for ( ; pts->hasNext() && counter != inputFile.getLimit(); ++counter) {
169 
170  STATUS("timeslice: " << setw(10) << counter << '\r'); DEBUG(endl);
171 
172  const JDAQTimeslice* timeslice = pts->next();
173 
174  int fIndex = timeslice->getFrameIndex();
175 
176  JDataSN preTrigger(TMax_ns, preTriggerThreshold);
177 
178  preTrigger(timeslice, moduleRouter, totSelector_ns);
179 
180  JTriggerSN trigger(TVeto_ns);
181 
182  trigger(preTrigger);
183 
184  // generate and configure event-based veto
185 
186  JVetoSet veto;
187 
188  if (triggeredEvents.count(fIndex)) {
189  veto = triggeredEvents.at(fIndex);
190  }
191 
192  // trigger.setVeto(veto);
193 
194  // count trigger
195 
196  JSNFilterM trgF0(M, 0);
197  JSNFilterM trgF1(M, 1);
198 
199  JSNFilterMV trgFV(M, veto);
200 
201  JRange<int> A = JRange<int>(2,31);
202  JSNFilterMV trgAV(A, veto);
203 
204  int rawCount = count_if(preTrigger.begin(), preTrigger.end(), trgF0);
205 
206  int trgCountF0 = count_if(trigger.begin(), trigger.end(), trgF0);
207  int trgCountF1 = count_if(trigger.begin(), trigger.end(), trgF1);
208  int trgCountFV = count_if(trigger.begin(), trigger.end(), trgFV);
209  int domCountF1 = trigger.getModules(trgF1).size();
210 
211  trgHistory[0].push_back(rawCount);
212  trgHistory[1].push_back(trgCountF0);
213  trgHistory[2].push_back(trgCountF1);
214  trgHistory[3].push_back(trgCountFV);
215  trgHistory[4].push_back(domCountF1);
216 
217  trigger.fill(MUL[1], trgF0);
218  trigger.fill(MUL[2], trgF1);
219  trigger.fill(MUL[3], trgFV);
220  trigger.fill(MUL[4], trgAV);
221 
222  }
223 
224  //-------------------------------------
225  // STEP 3: generating trigger summaries
226  // ------------------------------------
227 
228  for (int i = 0; i < nStages; i++) {
229  for (unsigned j = 0; j < trgHistory[i].size(); j++) {
230  SNT[i]->Fill(trgHistory[i][j]);
231  }
232  }
233 
234  SNT[0]->SetTitle("M[6,10] count before clustering");
235  SNT[1]->SetTitle("M[6,10] count after clustering");
236  SNT[2]->SetTitle("M[6,10] count after track self-veto");
237  SNT[3]->SetTitle("M[6,10] count after track trigger-veto");
238  SNT[4]->SetTitle("M[6,10] count after track self-veto, unique modules");
239 
240 
241  if (outputFile != "") {
242 
243  TFile out(outputFile.c_str(), "RECREATE");
244 
245  h_vtr->Write();
246  SNT.Write(out);
247  MUL.Write(out);
248 
249  putObject(out, JMeta(argc,argv));
250 
251  out.Close();
252  }
253 }
254 
255 
double getLength()
Get length of veto time range.
Definition: JSupernova.hh:115
Auxiliary class for ROOT I/O of application specific meta data.
Definition: JMeta.hh:70
Utility class to parse command line options.
Definition: JParser.hh:1514
General exception.
Definition: JException.hh:24
int main(int argc, char *argv[])
Definition: Main.cc:15
ROOT TTree parameter settings of various packages.
Auxiliary methods for geometrical methods.
JModuleSet getModules(JRange< int > A=JRange< int >(6, 10))
Get triggered modules after track veto.
Definition: JSupernova.hh:454
Auxiliary class to define a veto time window on a set of optical modules.
Definition: JSupernova.hh:80
#define STATUS(A)
Definition: JMessage.hh:63
Detector data structure.
Definition: JDetector.hh:89
Auxiliary class to select ROOT class based on class name.
SN filter based on veto window.
Definition: JSupernova.hh:364
Router for direct addressing of module data in detector data structure.
Auxiliary interface for direct access of elements in ROOT TChain.
Long64_t counter_type
Type definition for counter.
Dynamic ROOT object management.
Auxiliary class for a type holder.
Definition: JType.hh:19
string outputFile
Template definition for direct access of elements in ROOT TChain.
Data structure for detector geometry and calibration.
Simple wrapper around JModuleRouter class for direct addressing of PMT data in detector data structur...
int getFrameIndex() const
Get frame index.
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:41
Type definition of range.
Definition: JHead.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.
void fill(TH1D *out, const JSNFilter &F)
Fill histogram with multiplicity spectrum resulting from given filter.
Definition: JSupernova.hh:494
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1989
Auxiliary class to apply the supernova trigger to SN data.
Definition: JSupernova.hh:391
ROOT I/O of application specific meta data.
Data time slice.
SN filter based on multiplicity selection optional suppression of multi-module coincidences WARNING: ...
Definition: JSupernova.hh:329
bool putObject(TDirectory &dir, const TObject &object)
Write object to ROOT directory.
General purpose messaging.
Auxiliary class to select JTreeScanner based on ROOT class name.
#define FATAL(A)
Definition: JMessage.hh:67
Scanning of objects from multiple files according a format that follows from the extension of each fi...
Auxiliary class to build the supernova trigger dataset.
Definition: JSupernova.hh:134
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
Utility class to parse command line options.
Auxiliary class to select DAQ hits based on time-over-treshold value.
int j
Definition: JPolint.hh:703
do set_variable DETECTOR_TXT $WORKDIR detector
source $JPP_DIR setenv csh $JPP_DIR &dev null eval JShellParser o a A
int debug
debug level
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62
Auxiliary class to manage a set of vetoes.
Definition: JSupernova.hh:257