Jpp
Classes | Public Member Functions | Static Public Member Functions | Private Types | Private Attributes | List of all members
JFIT::JORCAShowerFit Class Reference

#include <JORCAShowerFit.hh>

Classes

struct  JShowerHypothesis
 Container of shower hypothesis: During the fit the algorithm creates a grid in direction and energy space, for each of these points a shower hypothesis is saved: direction, energy and LogLikelihood value. More...
 

Public Member Functions

 JORCAShowerFit ()
 Default constructor. More...
 
 JORCAShowerFit (const JLANG::JSharedPointer< const JDETECTOR::JModuleRouter > &moduleRouter, const JFIT::JShowerFitParameters_t &parameters, const std::string pdfFile)
 Parameterized constructor. More...
 
 ~JORCAShowerFit ()
 
void getJEvt (const KM3NETDAQ::JDAQTimeslice &timeSliceBuildL0, JFIT::JEvt &InPreFits, JFIT::JEvt &OutFits) const
 Declaration of the member function that actually performs the reconstruction. More...
 

Static Public Member Functions

static double estimateInitialEnergy (int nHits)
 Function to estimate the initial em shower energy. More...
 
static bool sortLogLik (JShowerHypothesis A, JShowerHypothesis B)
 Function to sort different shower hypotheses by their likelihood. More...
 

Private Types

typedef JFIT::JRegressor< JFIT::JShower3EZ, JFIT::JSimplexJRegressor_t
 

Private Attributes

JRegressor_t fit_
 
JLANG::JSharedPointer< const JDETECTOR::JModuleRoutermoduleRouter_
 
JFIT::JShowerFitParameters_t parameters_
 

Detailed Description

Definition at line 63 of file JORCAShowerFit.hh.

Member Typedef Documentation

◆ JRegressor_t

Definition at line 67 of file JORCAShowerFit.hh.

Constructor & Destructor Documentation

◆ JORCAShowerFit() [1/2]

JFIT::JORCAShowerFit::JORCAShowerFit ( )
inline

Default constructor.

Definition at line 78 of file JORCAShowerFit.hh.

79  {}

◆ JORCAShowerFit() [2/2]

JFIT::JORCAShowerFit::JORCAShowerFit ( const JLANG::JSharedPointer< const JDETECTOR::JModuleRouter > &  moduleRouter,
const JFIT::JShowerFitParameters_t parameters,
const std::string  pdfFile 
)
inline

Parameterized constructor.

Parameters
moduleRouterJSharedPointer of JModuleRouter, this is built via detector file.
parametersstruct that holds default-optimized parameters for the reconstruction, which can be found in $JPP_DATA.
pdfFilePDF file

Definition at line 89 of file JORCAShowerFit.hh.

91  :
92  moduleRouter_(moduleRouter),
93  parameters_(parameters)
94  {
96  JRegressor_t::T_ns.setRange(-parameters.Tmax_ns,parameters.Tmax_ns);
97  JRegressor_t::Vmax_npe = 20.0;
98  JRegressor_t::MAXIMUM_ITERATIONS = 1000;
99 
100  fit_ = {pdfFile};
101  fit_.estimator.reset(new JMEstimatorNull());
102  }

◆ ~JORCAShowerFit()

JFIT::JORCAShowerFit::~JORCAShowerFit ( )
inline

Definition at line 104 of file JORCAShowerFit.hh.

104 {}

Member Function Documentation

◆ estimateInitialEnergy()

static double JFIT::JORCAShowerFit::estimateInitialEnergy ( int  nHits)
inlinestatic

Function to estimate the initial em shower energy.

Definition at line 124 of file JORCAShowerFit.hh.

124  {
125 
126  double Energy = exp((nHits + 26)/26);
127  if(Energy > 100) Energy = 100;
128 
129  return Energy;
130  }

◆ sortLogLik()

static bool JFIT::JORCAShowerFit::sortLogLik ( JShowerHypothesis  A,
JShowerHypothesis  B 
)
inlinestatic

Function to sort different shower hypotheses by their likelihood.

Definition at line 137 of file JORCAShowerFit.hh.

138  {
139  return(A.LogLik < B.LogLik);
140  };

◆ getJEvt()

void JFIT::JORCAShowerFit::getJEvt ( const KM3NETDAQ::JDAQTimeslice timeSliceBuildL0,
JFIT::JEvt InPreFits,
JFIT::JEvt OutFits 
) const

Declaration of the member function that actually performs the reconstruction.

member function definition

Parameters
timeSliceBuildL0which is contructed via a JDAQEvent
InPreFitsinput fits, normally produced by JORCAShowerPositionFit
OutFitsoutput fits

Definition at line 164 of file JORCAShowerFit.hh.

167 {
168  using namespace std;
169 
170  const double Dmax_m = parameters_.Dmax_m;
171  const double R_Hz = parameters_.R_Hz;
172  const size_t numberOfPrefits = parameters_.numberOfPrefits;
173 
174  if ( InPreFits.empty() ) return;
175 
176  typedef vector<JTRIGGER::JHitL0> JDataL0_t;
177  typedef vector<JFIT::JPMTW0> JPMTW0_t;
179 
181 
182  JEvt::iterator __end = InPreFits.end();
183  std::copy(InPreFits.begin(), __end, back_inserter(OutFits));
184 
185  if (numberOfPrefits > 0) {
186  advance(__end = InPreFits.begin(), min(numberOfPrefits, InPreFits.size()));
187  }
188 
189  partial_sort(InPreFits.begin(), __end, InPreFits.end(), qualitySorter);
190 
191  JDataL0_t dataL0;
192  buildL0(timeSliceBuildL0, *moduleRouter_, back_inserter(dataL0));
193 
194  for (JEvt::const_iterator shower = InPreFits.begin(); shower != __end; ++shower) {
195  const JVector3D pos(getPosition(*shower));
196 
197  JDetectorSubset_t subdetector(detector, pos, Dmax_m);
198 
199  JPoint4D pt(pos, shower->getT());
200  JVector3D start_dir(0,0,0);
201 
202  // get the detector center and find the distance between the reco vertex and it.
203  const JCenter3D center(detector.begin(), detector.end());
204  const double radius = pos.getDistance(center);
205 
207 
208  double Nhits = 0;
209 
210  for (JDataL0_t::const_iterator i = dataL0.begin(); i != dataL0.end(); ++i) {
211 
212  const JHitW0 hit(*i, R_Hz);
213  const JPosition3D hit_pos(hit.getPosition());
214 
215  const double D = hit_pos.getDistance(pt.getPosition());
216  const double t_res = hit.getT() - pt.getT(hit_pos);
217 
218  const JDirection3D photonDir(hit_pos - pt.getPosition());
219  const double cosT = photonDir.getDot(hit.getDirection());
220 
221  if(D < Dmax_m && (t_res >= -25 && t_res <= 25) && (cosT >= -1 && cosT <= 0.1)){
222  top.insert(hit.getPMTIdentifier());
223  const JVector3D d(photonDir.getDX(), photonDir.getDY(), photonDir.getDZ());
224  start_dir.add(d);
225  }
226  if(D < Dmax_m && (t_res > -30 && t_res < 20)){
227  Nhits++;
228  }
229  }
230 
231  double start_E = estimateInitialEnergy(Nhits);
232 
233  // rotate subdetector wrt start_dir, subdetector is already subtracted by pos
234 
235  const JDirection3D conversion(start_dir.getX(), start_dir.getY(), start_dir.getZ());
236  const JRotation3D R(conversion);
237 
238  vector<JPMTW0> buffer;
239 
240  JFIT::JORCAShowerFit::JShowerHypothesis shower_hypothesis;
242 
243  double chi2 = 0.0;
244  double max_theta = 0, max_phi = 0, scan_step = 0;
245  double Nmin = 0;
246 
247  if(start_E > 10 && radius < 100){
248  max_theta = 30;
249  max_phi = 30;
250  scan_step = 5;
251  Nmin = 10;
252  } else {
253  max_theta = 35;
254  max_phi = 35;
255  scan_step = 5;
256  Nmin = 20;
257  }
258 
259  const JAngle3D start_dir_angles(start_dir.getX(), start_dir.getY(), start_dir.getZ());
260 
261  for(double E = -15; E <= 5; E += 7.5){ // scan in energy
262  for(double th = -max_theta; th <= max_theta; th += scan_step){
263  for(double ph = -max_phi; ph <= max_phi; ph += scan_step){
264 
265  buffer.clear();
266  chi2 = 0.0;
267 
268  //todo: use precise value of pi
269  const double theta = th * PI / 180;
270  const double phi = ph * PI / 180;
271  const double scan_E = start_E + E; // [GeV]
272 
273  const JAngle3D dir_shifted_angles(start_dir_angles.getTheta() + theta,
274  start_dir_angles.getPhi() + phi);
275  const JVector3D dir_shifted(dir_shifted_angles.getDX(),
276  dir_shifted_angles.getDY(),
277  dir_shifted_angles.getDZ());
278  const JDirection3D conversion(dir_shifted.getX(),
279  dir_shifted.getY(),
280  dir_shifted.getZ());
281  const JRotation3D R(conversion);
282  pt.rotate(R);
283 
284  for (JDetectorSubset_t::iterator module = subdetector.begin();
285  module != subdetector.end(); ++module) {
286 
287  module->rotate(R);
288  JPosition3D module_pos(module->getPosition());
289 
290  for (unsigned int i = 0; i != module->size(); ++i) {
291  const JDAQPMTIdentifier id(module->getID(), i);
292  JPMT pmt(module->getPMT(i));
293  buffer.push_back(JPMTW0(pmt, R_Hz, top.count(id)));
294  }
295  }
296 
297  for(JPMTW0_t::iterator pmt = buffer.begin(); pmt != buffer.end(); ++pmt ){
298  chi2 += fit_(JShower3EZ(pt, JVersor3Z(), JEnergy(log10(scan_E))), *pmt);
299  }
300 
301  shower_hypothesis.LogLik = chi2;
302  shower_hypothesis.Energy = scan_E;
303  shower_hypothesis.Direction = dir_shifted;
304 
305  shower_hypotheses.push_back(shower_hypothesis);
306 
307  pt.rotate_back(R);
308 
309  }
310  }
311  }
312 
313  vector<JShowerHypothesis>::iterator __end2 = shower_hypotheses.begin() + Nmin;
314  partial_sort(shower_hypotheses.begin(), __end2, shower_hypotheses.end(), sortLogLik);
315 
316  for(vector<JShowerHypothesis>::iterator hp = shower_hypotheses.begin(); hp != __end2; ++hp){
317 
318  fit_.step.resize(3);
319  const double dir_step = 0.05;
320  fit_.step[0] = JShower3EZ(JPoint4D(JVector3D(), 0.0), JVersor3Z(dir_step, 0.0), JEnergy());
321  fit_.step[1] = JShower3EZ(JPoint4D(JVector3D(), 0.0), JVersor3Z(0.0, dir_step), JEnergy());
322  fit_.step[2] = JShower3EZ(JPoint4D(JVector3D(), 0.0), JVersor3Z(), 0.05);
323 
324  const double scan_E = hp->Energy;
325 
326  buffer.clear();
327 
328  const JDirection3D startVersor(hp->Direction);
329 
330  const JDirection3D conversion(startVersor.getDX(), startVersor.getDY(), startVersor.getDZ());
331  JRotation3D R(conversion);
332  pt.rotate(R);
333 
334  for(JDetector::const_iterator module = detector.begin();
335  module != detector.end(); ++module) {
336  JPosition3D pos(module->getPosition());
337  pos.rotate(R);
338 
339  if (pt.getDistance(pos) < Dmax_m) {
340  for (unsigned int i = 0; i != module->size(); ++i) {
341  const JDAQPMTIdentifier id(module->getID(), i);
342  JPMT pmt(module->getPMT(i));
343  pmt.rotate(R);
344  buffer.push_back(JPMTW0(pmt, R_Hz, top.count(id)));
345  }
346  }
347  }
348 
349  const int NDF = buffer.size() - fit_.step.size();
350 
351  chi2 = fit_(JShower3EZ(pt, JVersor3Z(), JEnergy(log10(scan_E))), buffer.begin(), buffer.end());
352 
353  pt.rotate_back(R);
354 
355  const JVector3D resPos(fit_.value);
356  const JVersor3D resDir(fit_.value.getDX(), fit_.value.getDY(), fit_.value.getDZ());
357 
358  const JTrack3D td(resPos, resDir, shower->getT());
359  double E_reco_corrected = abs(fit_.value.getE() - 4);
360 
361  // Now that JShowerFit uses a Param PDF, the energy has to be corrected with a linear function obtained
362  // from a fit to the reconstruction error: this can be improved.
363 
364  const double p0 = -2.028;
365  const double p1 = 0.7683;
366  if(E_reco_corrected > 6){
367  E_reco_corrected = (E_reco_corrected - p0) / (1 + p1);
368  }
369 
370  JTrack3E tb(td, E_reco_corrected);
371  tb.rotate_back(R);
372 
373  const JFit &outFit = getFit(JHistory(shower->getHistory()).add(JSHOWERCOMPLETEFIT),
374  tb, -chi2, NDF, E_reco_corrected,
376 
377  OutFits.push_back(outFit);
378  }
379  }
380  sort(OutFits.begin(), OutFits.end(), qualitySorter);
381 }

Member Data Documentation

◆ fit_

JRegressor_t JFIT::JORCAShowerFit::fit_
mutableprivate

Definition at line 68 of file JORCAShowerFit.hh.

◆ moduleRouter_

JLANG::JSharedPointer<const JDETECTOR::JModuleRouter> JFIT::JORCAShowerFit::moduleRouter_
private

Definition at line 70 of file JORCAShowerFit.hh.

◆ parameters_

JFIT::JShowerFitParameters_t JFIT::JORCAShowerFit::parameters_
private

Definition at line 71 of file JORCAShowerFit.hh.


The documentation for this class was generated from the following file:
JFIT::JORCAShowerFit::JShowerHypothesis::LogLik
double LogLik
Definition: JORCAShowerFit.hh:117
JFIT::JEnergy
Data structure for fit of energy.
Definition: JEnergy.hh:24
JFIT::JHitW0
Auxiliary class for a hit with background rate value.
Definition: JHitW0.hh:21
JFIT::JORCAShowerFit::JShowerHypothesis::Direction
JVector3D Direction
Definition: JORCAShowerFit.hh:115
JGEOMETRY3D::JCenter3D
Center.
Definition: JGeometry3DToolkit.hh:29
JFIT::JORCAShowerFit::sortLogLik
static bool sortLogLik(JShowerHypothesis A, JShowerHypothesis B)
Function to sort different shower hypotheses by their likelihood.
Definition: JORCAShowerFit.hh:137
JROOT::advance
counter_type advance(counter_type &counter, const counter_type value, const counter_type limit=std::numeric_limits< counter_type >::max())
Advance counter.
Definition: JCounter.hh:35
JFIT::JORCAShowerFit::JShowerHypothesis
Container of shower hypothesis: During the fit the algorithm creates a grid in direction and energy s...
Definition: JORCAShowerFit.hh:113
JFIT::JShower3EZ
Data structure for fit of straight line in positive z-direction with energy.
Definition: JShower3EZ.hh:25
JGEOMETRY3D::JVersor3Z
Data structure for normalised vector in positive z-direction.
Definition: JVersor3Z.hh:37
JFIT::JORCAShowerFit::estimateInitialEnergy
static double estimateInitialEnergy(int nHits)
Function to estimate the initial em shower energy.
Definition: JORCAShowerFit.hh:124
JGEOMETRY3D::JVector3D::getDistance
double getDistance(const JVector3D &pos) const
Get distance to point.
Definition: JVector3D.hh:269
JFIT::JShowerFitParameters_t::numberOfPrefits
std::size_t numberOfPrefits
Definition: JShowerFitParameters_t.hh:25
JGEOMETRY3D::JDirection3D::getDot
double getDot(const JAngle3D &angle) const
Get dot product.
Definition: JDirection3D.hh:333
std::vector
Definition: JSTDTypes.hh:12
JFIT::JORCAShowerFit::fit_
JRegressor_t fit_
Definition: JORCAShowerFit.hh:68
JGEOMETRY3D::JDirection3D
Data structure for direction in three dimensions.
Definition: JDirection3D.hh:32
JSHOWERCOMPLETEFIT
static const int JSHOWERCOMPLETEFIT
Definition: reconstruction.hh:25
JGEOMETRY3D::JVersor3D
Data structure for normalised vector in three dimensions.
Definition: JVersor3D.hh:23
JFIT::JShowerFitParameters_t::Dmax_m
double Dmax_m
Definition: JShowerFitParameters_t.hh:23
JFIT::JPMTW0
Auxiliary class for handling PMT geometry, rate and response.
Definition: JPMTW0.hh:22
JAANET::copy
void copy(const Head &from, JHead &to)
Copy header from from to to.
Definition: JHead.cc:152
JGEOMETRY3D::JPosition3D::rotate
JPosition3D & rotate(const JRotation3D &R)
Rotate.
Definition: JPosition3D.hh:185
JFIT::JMEstimatorNull
Null M-estimator.
Definition: JMEstimator.hh:53
JTRIGGER::JBuildL0
Template L0 hit builder.
Definition: JBuildL0.hh:34
JGEOMETRY3D::JVector3D
Data structure for vector in three dimensions.
Definition: JVector3D.hh:33
JFIT::qualitySorter
bool qualitySorter(const JFit &first, const JFit &second)
Comparison of fit results.
Definition: JEvtToolkit.hh:223
debug
int debug
debug level
Definition: JSirene.cc:59
JLANG::JReference::getReference
const JClass_t & getReference() const
Get reference to object.
Definition: JReference.hh:38
JGEOMETRY3D::JPosition3D
Data structure for position in three dimensions.
Definition: JPosition3D.hh:35
std::multiset
Definition: JSTDTypes.hh:14
JGEOMETRY3D::JTrack3E
3D track with energy.
Definition: JTrack3E.hh:29
JGEOMETRY3D::JAxis3D::rotate
JAxis3D & rotate(const JRotation3D &R)
Rotate axis.
Definition: JAxis3D.hh:236
JDETECTOR::JPMT
Data structure for PMT geometry and calibration.
Definition: JPMT.hh:53
JFIT::JPoint4D
Data structure for vertex fit.
Definition: JPoint4D.hh:22
JGEOMETRY3D::JAngle3D
Data structure for angles in three dimensions.
Definition: JAngle3D.hh:31
JFIT::JShowerFitParameters_t::R_Hz
double R_Hz
Definition: JShowerFitParameters_t.hh:24
KM3NETDAQ::JDAQPMTIdentifier
PMT identifier.
Definition: JDAQPMTIdentifier.hh:20
JDETECTOR::JDetector
Detector data structure.
Definition: JDetector.hh:80
JAANET::detector
Detector file.
Definition: JHead.hh:130
JFIT::JORCAShowerFit::JShowerHypothesis::Energy
double Energy
Definition: JORCAShowerFit.hh:116
std
Definition: jaanetDictionary.h:36
JGEOMETRY3D::JTrack3D
3D track.
Definition: JTrack3D.hh:30
JTOOLS::PI
static const double PI
Constants.
Definition: JConstants.hh:20
JFIT::JShowerFitParameters_t::Tmax_ns
double Tmax_ns
Definition: JShowerFitParameters_t.hh:22
JFIT::getFit
JFit getFit(const JHistory &history, const JTrack3D &track, const double Q, const int NDF, const double energy=0.0, const int status=0)
Get fit.
Definition: JEvtToolkit.hh:116
p1
TPaveText * p1
Definition: JDrawModule3D.cc:35
JGEOMETRY3D::JVector3D::add
JVector3D & add(const JVector3D &vector)
Add vector.
Definition: JVector3D.hh:141
JFIT::JHistory
Container for historical events.
Definition: JHistory.hh:95
JFIT::JORCAShowerFit::moduleRouter_
JLANG::JSharedPointer< const JDETECTOR::JModuleRouter > moduleRouter_
Definition: JORCAShowerFit.hh:70
JFIT::getPosition
JPosition3D getPosition(const JFit &fit)
Get position.
Definition: JEvtToolkit.hh:63
JDETECTOR::JDetectorSubset_t
Detector subset without binary search functionality.
Definition: JDetectorSubset.hh:37
JGEOMETRY3D::JRotation3D
Rotation matrix.
Definition: JRotation3D.hh:111
JFIT::JORCAShowerFit::parameters_
JFIT::JShowerFitParameters_t parameters_
Definition: JORCAShowerFit.hh:71
JFIT::JFit
Data structure for track fit results.
Definition: JEvt.hh:31
JFIT::OKAY
Definition: JFitStatus.hh:22