Jpp  19.0.0
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JShowerDirectionPrefit.hh
Go to the documentation of this file.
1 #ifndef JSHOWERDIRECTIONPREFIT_INCLUDE
2 #define JSHOWERDIRECTIONPREFIT_INCLUDE
3 
4 #include <string>
5 #include <iostream>
6 #include <set>
7 #include <vector>
8 #include <algorithm>
9 #include <memory>
10 
13 
14 #include "JTrigger/JHit.hh"
15 #include "JTrigger/JTimeslice.hh"
16 #include "JTrigger/JHitL0.hh"
17 #include "JTrigger/JHitL1.hh"
18 #include "JTrigger/JHitR1.hh"
19 #include "JTrigger/JBuildL0.hh"
20 #include "JTrigger/JBuildL2.hh"
21 #include "JTrigger/JAlgorithm.hh"
22 #include "JTrigger/JMatch3G.hh"
23 #include "JTrigger/JBind2nd.hh"
24 
26 
28 #include "JFit/JFitToolkit.hh"
29 #include "JFit/JPoint4D.hh"
30 #include "JFit/JModel.hh"
31 
32 #include "JAAnet/JAAnetToolkit.hh"
33 
34 #include "JReconstruction/JEvt.hh"
37 
38 #include "JDetector/JDetector.hh"
40 
41 #include "JPhysics/KM3NeT.hh"
42 
43 #include "JGeometry3D/JVector3D.hh"
45 #include "JGeometry3D/JOmega3D.hh"
46 #include "JGeometry3D/JTrack3D.hh"
47 #include "JGeometry3D/JTrack3E.hh"
48 #include "JGeometry3D/JShower3E.hh"
49 
50 #include "JLang/JSharedPointer.hh"
51 
52 /**
53  * \author adomi, vcarretero
54  */
55 namespace JRECONSTRUCTION {}
56 namespace JPP { using namespace JRECONSTRUCTION; }
57 
58 namespace JRECONSTRUCTION {
59 
62  using JFIT::JRegressor;
63  using JFIT::JEnergy;
64  using JFIT::JShower3EZ;
66 
67  /**
68  * class to handle the direction fit of the shower reconstruction, mainly dedicated for ORCA
69  */
72  public JRegressor<JShower3EZ, JAbstractMinimiser>
73  {
74 
76  using JRegressor_t::operator();
79 
80  public:
81 
82  /**
83  * Parameterized constructor
84  *
85  * \param parameters struct that holds default-optimized parameters for the reconstruction, available in $JPP_DATA.
86  * \param router module router, this is built via detector file.
87  * \param summary summary router
88  * \param pdfFile PDF file
89  * \param debug debug
90  */
92  const JModuleRouter& router,
93  const JSummaryRouter& summary,
94  const std::string pdfFile,
95  const int debug = 0):
97  JRegressor_t(pdfFile),
98  omega(parameters.scanAngle_deg * JMATH::PI / 180.0),
99  router(router),
100  summary(summary)
101  {
102  using namespace JPP;
103 
105  JRegressor_t::T_ns.setRange(parameters.TMin_ns, parameters.TMax_ns);
107 
108  if (Emin_GeV > Emax_GeV || En <= 1) {
109  THROW(JException, "Invalid energy input " << Emin_GeV << ' ' << Emax_GeV << ' ' << En);
110  }
111 
112  const double base = std::pow((Emax_GeV / Emin_GeV), 1.0 / (En - 1));
113 
114  for (int i = 0; i != En; ++i) {
115  Ev.push_back(Emin_GeV * std::pow(base, i));
116  }
117  }
118 
119  /**
120  * Declaration of the member function that actually performs the reconstruction
121  *
122  * \param event DAQ event
123  * \param in input fits
124  * \return result fits
125  */
127  {
128  using namespace std;
129  using namespace JPP;
130 
131  typedef vector<JHitL0> JDataL0_t;
132  JBuildL0<JHitL0> buildL0;
133 
134  JEvt out;
135 
137 
138  JDataL0_t dataL0;
139 
140  buildL0(JDAQTimeslice(event, true), router, back_inserter(dataL0));
141 
142  for (JEvt::const_iterator shower = in.begin(); shower != in.end(); ++shower) {
143 
144  const JVertex3D Vertex(getPosition(*shower),shower->getT());
146 
147  const JFIT::JModel<JPoint4D> match(Vertex, DMax_m, JRegressor_t::T_ns);
148 
149  for (JDataL0_t::const_iterator i = dataL0.begin(); i != dataL0.end(); ++i) {
150  if (match(*i)) {
151  top.insert(i->getPMTIdentifier());
152  }
153  }
154 
155  const JDetectorSubset_t subdetector(detector, Vertex.getPosition(), DMax_m);
156 
158  for (JDetectorSubset_t::const_iterator module = subdetector.begin(); module != subdetector.end(); ++module) {
159  const JDAQSummaryFrame& frame = summary.getSummaryFrame(module->getID());
160 
161  JModule dom(*module);
162 
163  for (size_t i = 0; i != dom.size(); ++i) {
164 
165  if (getDAQStatus(frame, *module, i) &&
166  getPMTStatus(frame, *module, i) &&
167  frame[i].is_valid() &&
168  !module->getPMT(i).has(PMT_DISABLE)) {
169 
170  const JDAQPMTIdentifier id(module->getID(), i);
171 
172  const double rate_Hz = summary.getRate(id);
173  const size_t count = top.count(id);
174 
175  data.push_back(JPMTW0(dom.getPMT(i), rate_Hz, count));
176  }
177  }
178  }
179 
180  const JShower3EZ sh(JVertex3D(JVector3D(0,0,0), sh.getT()), JVersor3Z(), 1);
181 
182  for (JOmega3D_t::const_iterator dir = omega.begin(); dir != omega.end(); ++dir) {
183  const JRotation3D R(*dir);
184  vector<double> chi2v(En, 0.0);
185  for (vector<JPMTW0>::const_iterator i = data.begin(); i != data.end(); ++i) {
186  JPMTW0 pmt = *i;
187  pmt.rotate(R);
188  JNPE_t::result_type H1 = (*this).getH1(sh, pmt);
189  JNPE_t::result_type H0 = (*this).getH0(pmt.getR()); // getH0 = Get background hypothesis value for time integrated PDF.
190  const bool hit = pmt.getN() != 0;
191  for(size_t j=0;j!=chi2v.size();++j){
192  chi2v[j]+= getChi2(get_value(H0+Ev[j]*H1), hit); // H1 is linear with E
193  }
194  }
195  //Store only the lowest chi2 for a given direction
196  auto minChi2 = std::min_element(chi2v.begin(), chi2v.end());
197  JShower3E sh_fit(Vertex, JDirection3D(*dir), Ev[minChi2-chi2v.begin()]);
198  out.push_back(getFit(JHistory(shower->getHistory()).add(JSHOWERDIRECTIONPREFIT), sh_fit, getQuality(*minChi2), data.size(), sh_fit.getE()));
199  }
200  }
201 
202  return out;
203 
204  }
207  };
208 }
209 
210 #endif
Auxiliary methods to evaluate Poisson probabilities and chi2.
static int debug
debug level (default is off).
Definition: JMessage.hh:45
Template definition of a data regressor of given model.
Definition: JRegressor.hh:70
double TMax_ns
maximum time for local coincidences [ns]
General exception.
Definition: JException.hh:24
JEvt operator()(const KM3NETDAQ::JDAQEvent &event, const JFIT::JEvt &in)
Declaration of the member function that actually performs the reconstruction.
JRegressor< JShower3EZ, JAbstractMinimiser > JRegressor_t
double getR() const
Get rate.
Definition: JPMTW0.hh:56
Data structure for direction in three dimensions.
Definition: JDirection3D.hh:33
class to handle the direction fit of the shower reconstruction, mainly dedicated for ORCA ...
Data structure for a composite optical module.
Definition: JModule.hh:67
double getQuality(const double chi2, const int NDF)
Get quality of fit.
Algorithms for hit clustering and sorting.
Regressor function object for JShower3EZ fit using Abstract minimiser, that just computes the chi2 wi...
Template specialisation of L0 builder for JHitL0 data type.
Definition: JBuildL0.hh:102
Auxiliary class for handling PMT geometry, rate and response.
Definition: JPMTW0.hh:22
Detector data structure.
Definition: JDetector.hh:89
Router for direct addressing of module data in detector data structure.
Rotation matrix.
Definition: JRotation3D.hh:111
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Definition: JException.hh:712
double getRate() const
Get default rate.
then usage $script< input file >[option] nPossible options count
Definition: JVolume1D.sh:31
*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
Direction set covering (part of) solid angle.
Definition: JOmega3D.hh:66
const JDAQSummaryFrame & getSummaryFrame() const
Get default summary frame.
Basic data structure for time and time over threshold information of hit.
3D track with energy.
Definition: JTrack3E.hh:30
Data structure for detector geometry and calibration.
static JTimeRange T_ns
Time window with respect to Cherenkov hypothesis [ns].
Basic data structure for L0 hit.
Auxiliary class to extract a subset of optical modules from a detector.
Properties of KM3NeT PMT and deep-sea water.
Data structure for fit of straight line in positive z-direction with energy.
Definition: JShower3EZ.hh:27
JAxis3D & rotate(const JRotation3D &R)
Rotate axis.
Definition: JAxis3D.hh:225
Detector file.
Definition: JHead.hh:226
Data structure for vector in three dimensions.
Definition: JVector3D.hh:34
Data storage class for rate measurements of all PMTs in one module.
static const int JSHOWERDIRECTIONPREFIT
static const int PMT_DISABLE
KM3NeT Data Definitions v3.4.0-8-ge14cb17 https://git.km3net.de/common/km3net-dataformat.
Definition: pmt_status.hh:12
JShowerDirectionPrefit(const JShowerDirectionPrefitParameters_t &parameters, const JModuleRouter &router, const JSummaryRouter &summary, const std::string pdfFile, const int debug=0)
Parameterized constructor.
JPosition3D getPosition(const Vec &pos)
Get position.
Definition of hit and track types and auxiliary methods for handling Monte Carlo data.
T pow(const T &x, const double y)
Power .
Definition: JMath.hh:97
Data time slice.
static const double PI
Mathematical constants.
Router for fast addressing of summary data in KM3NETDAQ::JDAQSummaryslice data structure as a functio...
const JPMT & getPMT(const int index) const
Get PMT.
Definition: JModule.hh:172
Detector subset without binary search functionality.
static double Vmax_npe
Maximal integral of PDF [npe].
Reduced data structure for L1 hit.
then JCookie sh JDataQuality D $DETECTOR_ID R
Definition: JDataQuality.sh:41
bool getPMTStatus(const JStatus &status)
Test status of PMT.
JFIT::JHistory JHistory
Definition: JHistory.hh:354
const JClass_t & getReference() const
Get reference to object.
Definition: JReference.hh:38
Data structure for set of track fit results.
JFit getFit(const int id, const JMODEL::JString &string)
Get fit parameters of string.
int getN() const
Get number of hits.
Definition: JPMTW0.hh:67
then fatal The output file must have the wildcard in the e g root fi eval JPrintDetector a $DETECTOR O IDENTIFIER eval JPrintDetector a $DETECTOR O SUMMARY JAcoustics sh $DETECTOR_ID source JAcousticsToolkit sh CHECK_EXIT_CODE typeset A EMITTERS get_tripods $WORKDIR tripod txt EMITTERS get_transmitters $WORKDIR transmitter txt EMITTERS for EMITTER in
Definition: JCanberra.sh:48
Data regression method for JFIT::JShower3EZ.
Data structure for fit of energy.
Definition: JEnergy.hh:28
int j
Definition: JPolint.hh:792
bool getDAQStatus(const JDAQFrameStatus &frame, const JStatus &status)
Test status of DAQ.
Template specialisation of class JModel to match hit with bright point.
Definition: JFit/JModel.hh:121
double getChi2(const double P)
Get chi2 corresponding to given probability.
Definition: JFitToolkit.hh:56
Data structure for normalised vector in positive z-direction.
Definition: JVersor3Z.hh:39
JResultEvaluator< JResult_t >::result_type get_value(const JResult_t &value)
Helper method to recursively evaluate a to function value.
Definition: JResult.hh:998
Abstract minimiser.
Definition: JRegressor.hh:25
double TMin_ns
minimum time for local coincidences [ns]
Match operator for Cherenkov light from shower in any direction.
Basic data structure for L1 hit.
int debug
debug level