Jpp  master_rocky-37-gf0c5bc59d
the software that should make you happy
Functions
JSpark.cc File Reference

Example program to check for sparks. More...

#include <string>
#include <iostream>
#include <iomanip>
#include <set>
#include <vector>
#include <cmath>
#include "TROOT.h"
#include "TFile.h"
#include "TH1D.h"
#include "TH2D.h"
#include "km3net-dataformat/offline/Head.hh"
#include "km3net-dataformat/offline/Evt.hh"
#include "km3net-dataformat/online/JDAQTriggerMask.hh"
#include "JAAnet/JAAnetToolkit.hh"
#include "JAAnet/JHead.hh"
#include "JAAnet/JHeadToolkit.hh"
#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JDetector/JStringRouter.hh"
#include "JDynamics/JDynamics.hh"
#include "JROOT/JManager.hh"
#include "JSupport/JMultipleFileScanner.hh"
#include "JSupport/JFileRecorder.hh"
#include "JSupport/JSupport.hh"
#include "JSupport/JMeta.hh"
#include "JSupport/JMonteCarloFileSupportkit.hh"
#include "Jeep/JProperties.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"

Go to the source code of this file.

Functions

int main (int argc, char **argv)
 

Detailed Description

Example program to check for sparks.

Author
mdejong

Definition in file JSpark.cc.

Function Documentation

◆ main()

int main ( int  argc,
char **  argv 
)

Definition at line 47 of file JSpark.cc.

48 {
49  using namespace std;
50  using namespace JPP;
51  using namespace KM3NETDAQ;
52 
54 
55  JMultipleFileScanner<Evt> inputFile;
58  JLimit_t& numberOfEvents = inputFile.getLimit();
59  string detectorFile;
60  double Tmax_s;
61  bool use_weights;
62  JDAQTriggerMask trigger_mask = TRIGGER_MASK_ON;
63  set<double> Dmax_m = { 0.3, 2.0 };
64  int debug;
65 
66  try {
67 
68  JProperties properties;
69 
70  properties.insert(gmake_property(trigger_mask));
71  properties.insert(gmake_property(Dmax_m));
72 
73  JParser<> zap("Example program to check for sparks.");
74 
75  zap['f'] = make_field(inputFile);
76  zap['n'] = make_field(numberOfEvents) = JLimit::max();
77  zap['a'] = make_field(detectorFile);
78  zap['@'] = make_field(properties) = JPARSER::initialised();
79  zap['+'] = make_field(calibrationFile) = JPARSER::initialised();
80  zap['T'] = make_field(Tmax_s) = 100.0;
81  zap['o'] = make_field(outputFile);
82  zap['W'] = make_field(use_weights);
83  zap['d'] = make_field(debug) = 2;
84 
85  zap(argc, argv);
86  }
87  catch(const exception& error) {
88  FATAL(error.what() << endl);
89  }
90 
91 
93 
94  try {
95  load(detectorFile, detector);
96  }
97  catch(const JException& error) {
98  FATAL(error);
99  }
100 
101  sort(detector.begin(), detector.end(), make_comparator(&JModule::getLocation));
102 
103  unique_ptr<JDynamics> dynamics;
104 
105  try {
106 
107  dynamics.reset(new JDynamics(detector, Tmax_s));
108 
109  dynamics->load(calibrationFile);
110  }
111  catch(const exception& error) {
112  if (!calibrationFile.empty()) {
113  FATAL(error.what());
114  }
115  }
116 
117 
118  const floor_range range = getRangeOfFloors(detector);
119  const JStringRouter router(detector);
120 
121  TH1D h0("h0", NULL, 50, -1.0, 2.0);
122 
123  JManager<int, TH2D> H2(new TH2D("h2[%]", NULL, router.size(), -0.5, router.size() - 0.5, range.getUpperLimit(), 1 - 0.5, range.getUpperLimit() + 0.5));
124 
125  for (Int_t i = 1; i <= H2->GetXaxis()->GetNbins(); ++i) {
126  H2->GetXaxis()->SetBinLabel(i, MAKE_CSTRING(router.at(i-1)));
127  }
128  for (Int_t i = 1; i <= H2->GetYaxis()->GetNbins(); ++i) {
129  H2->GetYaxis()->SetBinLabel(i, MAKE_CSTRING(i));
130  }
131 
132 
133  outputFile.open();
134  outputFile.put(JMeta(argc, argv));
135 
136  double W = 1.0;
137 
138  for (string file_name = ""; inputFile.hasNext(); ) {
139 
140  STATUS("event: " << setw(10) << inputFile.getCounter() << '\r'); DEBUG(endl);
141 
142  const Evt* evt = inputFile.next();
143 
144  if (JDAQTriggerMask(evt->trigger_mask).hasTriggerMask(trigger_mask)) {
145 
146  if (dynamics) {
147  dynamics->update(evt->t.GetSec());
148  }
149 
150  if (use_weights) {
151 
152  if (!evt->w.empty()) {
154  }
155 
156  } else {
157 
158  if (file_name != inputFile.getFilename()) {
159 
160  file_name = inputFile.getFilename();
161 
162  try {
163 
164  const JAANET::JHead head(getHeader(file_name));
165 
166  if (is_mupage(head)) { W = head.DAQ.livetime_s / head.livetime.numberOfSeconds; }
167  if (is_pure_noise(head)) { W = head.DAQ.livetime_s / head.K40.livetime_s; }
168  }
169  catch(const exception& error) {}
170  }
171  }
172 
173  double Dmin = numeric_limits<double>::max();
174  JLocation location(-1, -1);
175 
176  if (has_reconstructed_jppshower(*evt)) {
177 
178  const Trk trk = get_best_reconstructed_jppshower(*evt);
179  const JTrack3D ts = getTrack(trk);
180 
181  if (trk.status != TRK_ST_UNDEFINED) {
182 
183  for (const auto& module : (dynamics ? dynamics->getDetector() : detector)) {
184 
185  if (module.getFloor() != 0) {
186 
187  const double D = getDistance(module.getPosition(), ts.getPosition());
188 
189  if (D < Dmin) {
190  Dmin = D;
191  location = module.getLocation();
192  }
193  }
194  }
195  }
196 
197  if (Dmin != numeric_limits<double>::max()) {
198 
199  h0.Fill(log10(Dmin), W);
200 
201  const int i = distance(Dmax_m.begin(), Dmax_m.lower_bound(Dmin));
202 
203  H2 ->Fill((double) router.getIndex(location.getString()), (double) location.getFloor(), W);
204  H2[i]->Fill((double) router.getIndex(location.getString()), (double) location.getFloor(), W);
205  }
206  }
207  }
208  //outputFile.put(*evt);
209  }
210  STATUS(endl);
211 
212 
213  outputFile.put(h0);
214  outputFile.put(*H2);
215 
216  for (auto& h2 : H2) {
217  outputFile.put(*h2.second);
218  }
219 
220  outputFile.close();
221 }
string outputFile
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62
#define STATUS(A)
Definition: JMessage.hh:63
#define FATAL(A)
Definition: JMessage.hh:67
int debug
debug level
Definition: JSirene.cc:69
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:2142
#define MAKE_CSTRING(A)
Make C-string.
Definition: JPrint.hh:72
#define gmake_property(A)
macros to convert (template) parameter to JPropertiesElement object
std::vector< T >::difference_type distance(typename std::vector< T >::const_iterator first, typename PhysicsEvent::const_iterator< T > second)
Specialisation of STL distance.
Monte Carlo run header.
Definition: JHead.hh:1236
Detector data structure.
Definition: JDetector.hh:96
Logical location of module.
Definition: JLocation.hh:40
Utility class to parse parameter values.
Definition: JProperties.hh:501
const JPosition3D & getPosition() const
Get position.
Definition: JPosition3D.hh:130
General exception.
Definition: JException.hh:24
Utility class to parse command line options.
Definition: JParser.hh:1698
Auxiliary class to manage set of compatible ROOT objects (e.g. histograms) using unique keys.
Definition: JManager.hh:47
Object writing to file.
virtual bool hasNext() override
Check availability of next element.
counter_type getCounter() const
Get counter.
virtual const pointer_type & next() override
Get next element.
const std::string & getFilename() const
Get current file name.
T getUpperLimit() const
Get upper limit.
Definition: JRange.hh:213
Auxiliary class for trigger mask.
bool hasTriggerMask(const JDAQTriggerMask &mask) const
Has trigger bit pattern.
JTrack3E getTrack(const Trk &track)
Get track.
bool is_pure_noise(const JHead &header)
Check for generator.
bool is_mupage(const JHead &header)
Check for generator.
Definition: JHeadToolkit.hh:85
floor_range getRangeOfFloors(const JDetector &detector)
Get range of floors.
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
JComparator< JResult_t T::*, JComparison::lt > make_comparator(JResult_t T::*member)
Helper method to create comparator between values of data member.
double getDistance(const JFirst_t &first, const JSecond_t &second)
Get distance between objects.
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Head getHeader(const JMultipleFileScanner_t &file_list)
Get Monte Carlo header.
KM3NeT DAQ data structures and auxiliaries.
Definition: DataQueue.cc:39
static const JDAQTriggerMask TRIGGER_MASK_ON
Trigger mask on;.
Definition: JSTDTypes.hh:14
The Evt class respresent a Monte Carlo (MC) event as well as an offline event.
Definition: Evt.hh:21
std::vector< double > w
MC: Weights w[0]=w1, w[1]=w2, w[2]=w3 (see e.g. Tag list or km3net-dataformat/definitions)
Definition: Evt.hh:42
ULong64_t trigger_mask
trigger mask from raw data (i.e. the trigger bits)
Definition: Evt.hh:30
TTimeStamp t
UTC time of the timeslice, or the event_time for MC. (default: 01 Jan 1970 00:00:00)
Definition: Evt.hh:33
Detector file.
Definition: JHead.hh:227
Router for mapping of string identifier to index.
Dynamic detector calibration.
Definition: JDynamics.hh:84
Type list.
Definition: JTypeList.hh:23
Empty structure for specification of parser element that is initialised (i.e. does not require input)...
Definition: JParser.hh:68
Auxiliary class for defining the range of iterations of objects.
Definition: JLimit.hh:45
Auxiliary class for ROOT I/O of application specific meta data.
Definition: JMeta.hh:72
The Trk class represents a Monte Carlo (MC) particle as well as a reconstructed track/shower.
Definition: Trk.hh:15
int status
MC status code, see km3net-dataformat/definitions/trkmembers.csv for values.
Definition: Trk.hh:28
const Trk & get_best_reconstructed_jppshower(const Evt &evt)
Get best reconstructed shower.
bool has_reconstructed_jppshower(const Evt &evt)
Test whether given event has a track with shower reconstruction.
static const int TRK_ST_UNDEFINED
status was not defined for this MC track (all reco tracks have this value)
Definition: trkmembers.hh:14
static const int WEIGHTLIST_RUN_BY_RUN_WEIGHT
w[1]*DAQ_livetime / MC_evts_summary::n_gen (gseagen; [GeV m2 sr s]) or DAQ_livetime/MC_evts_summary::...
Definition: weightlist.hh:18