Jpp
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Static Public Member Functions | Private Types | Private Attributes | List of all members
JFIT::JORCAShowerPrefit Class Reference

class to handle first step of the shower reconstruction, mainly dedicated for ORCA More...

#include <JORCAShowerPrefit.hh>

Public Member Functions

 JORCAShowerPrefit ()
 Default constructor. More...
 
 JORCAShowerPrefit (const JLANG::JSharedPointer< const JDETECTOR::JModuleRouter > &moduleRouter, const JFIT::JShowerPrefitParameters_t &parameters)
 Parameterized constructor. More...
 
 ~JORCAShowerPrefit ()
 
void getJEvt (const KM3NETDAQ::JDAQTimeslice &timeSliceBuildL0, const KM3NETDAQ::JDAQTimeslice &timeSliceBuildL2, JFIT::JEvt &OutFits) const
 Declaration of Member function that actually performs the reconstruction. More...
 

Static Public Member Functions

static bool totW (const JTRIGGER::JHitR1 &first, const JTRIGGER::JHitR1 &second)
 Member function used in the predicate of std::sort, to sort hits according to ToT and W. More...
 

Private Types

typedef JFIT::JRegressor
< JPoint4D, JSimplex
JRegressor_t
 

Private Attributes

std::shared_ptr< JRegressor_tfit2_
 
JLANG::JSharedPointer< const
JDETECTOR::JModuleRouter
moduleRouter_
 
JFIT::JShowerPrefitParameters_t parameters_
 

Detailed Description

class to handle first step of the shower reconstruction, mainly dedicated for ORCA

Definition at line 55 of file JORCAShowerPrefit.hh.

Member Typedef Documentation

Definition at line 58 of file JORCAShowerPrefit.hh.

Constructor & Destructor Documentation

JFIT::JORCAShowerPrefit::JORCAShowerPrefit ( )
inline

Default constructor.

Definition at line 69 of file JORCAShowerPrefit.hh.

70  {}
JFIT::JORCAShowerPrefit::JORCAShowerPrefit ( const JLANG::JSharedPointer< const JDETECTOR::JModuleRouter > &  moduleRouter,
const JFIT::JShowerPrefitParameters_t parameters 
)
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.

Definition at line 79 of file JORCAShowerPrefit.hh.

80  :
81  moduleRouter_(moduleRouter),
82  parameters_(parameters)
83  {
86 
87  fit2_ = std::make_shared<JRegressor_t>(true);
88  }
static int debug
debug level (default is off).
Definition: JMessage.hh:43
JFIT::JShowerPrefitParameters_t parameters_
JLANG::JSharedPointer< const JDETECTOR::JModuleRouter > moduleRouter_
std::shared_ptr< JRegressor_t > fit2_
static int MAXIMUM_ITERATIONS
maximal number of iterations
Definition: JSimplex.hh:224
JFIT::JORCAShowerPrefit::~JORCAShowerPrefit ( )
inline

Definition at line 90 of file JORCAShowerPrefit.hh.

90 {}

Member Function Documentation

void JFIT::JORCAShowerPrefit::getJEvt ( const KM3NETDAQ::JDAQTimeslice timeSliceBuildL0,
const KM3NETDAQ::JDAQTimeslice timeSliceBuildL2,
JFIT::JEvt OutFits 
) const

Declaration of Member function that actually performs the reconstruction.

member function definition

Parameters
timeSliceBuildL0which is contructed via a JDAQEvent
timeSliceBuildL2which is contructed via a JDAQEvent
OutFitsnon const (output) JEvt.

Definition at line 129 of file JORCAShowerPrefit.hh.

132 {
133  const double Tmax_ns = parameters_.Tmax_ns; // [ns]
134  const double ctMin = parameters_.ctMin; // optimal value
135  const double Dmax_m = parameters_.Dmax_m; // [m] optimal value
137  const JTOOLS::JRange<double> time_grid = parameters_.time_grid;
138  const double pos_step = parameters_.pos_step; // [m] density of the grid
139  const double time_step = parameters_.time_step; // [ns] density of the grid
140  const double sigma_ns = parameters_.sigma_ns;
141  const int numberOfOutliers = parameters_.numberOfOutliers;
142  const bool useL0 = parameters_.useL0;
143 
144  const int FACTORY_LIMIT = 40;
145  const double STANDARD_DEVIATIONS = 1.0;
146  const int NUMBER_OF_L1HITS = (useL0 ? 1 : 4);
147 
148  typedef std::vector<JTRIGGER::JHitL0> JDataL0_t;
149  typedef std::vector<JTRIGGER::JHitL1> JDataL1_t;
150 
152  const JTRIGGER::JBuildL2<JTRIGGER::JHitL1> buildL2(2, Tmax_ns, ctMin);
153  const JTRIGGER::JMatch3G<JTRIGGER::JHitL1> match3G(Dmax_m, Tmax_ns); // CAUSALITY FOR SHOWERS
154  const JFIT::JHitL1Comparator compare;
155 
156  //define empty vectors to fill with JHitL0 and JHitL1
157  JDataL0_t dataL0;
158  JDataL0_t dataL0_selected;
159  JDataL1_t data;
160 
161  buildL0(timeSliceBuildL0, *moduleRouter_, back_inserter(dataL0));
162  buildL2(timeSliceBuildL1, *moduleRouter_, back_inserter(data));
163 
164  JDataL1_t::iterator __end = clusterizeWeight(data.begin(), data.end(), match3G);
165 
166  if(std::distance(data.begin(), __end) >= NUMBER_OF_L1HITS){
167 
168  JDataL1_t dataL1;
169  std::copy(data.begin(), data.end(), std::back_inserter(dataL1));
170  JDataL1_t::iterator __end1 = dataL1.end();
171 
172  __end1 = clusterizeWeight(dataL1.begin(), __end1, match3G);
173 
174  if (useL0) {
175 
176  dataL1.erase(__end1, dataL1.end());
177  dataL1.reserve(dataL1.size() + dataL0.size());
178  __end1 = dataL1.end();
179 
180  for (JDataL0_t::const_iterator i = dataL0.begin(); i != dataL0.end(); ++i) {
181 
182  const JTRIGGER::JHitL1 hit(*i);
183 
184  if (std::count_if(dataL1.begin(), __end1, JTRIGGER::JBind2nd(match3G,hit)) >= NUMBER_OF_L1HITS){
185  dataL1.push_back(hit);
186  }
187  }
188 
189  std::sort(__end1, dataL1.end(), compare);
190  __end1 = clusterize(__end1, dataL1.end(), match3G);
191 
192  }
193 
194  // 4D fit
196  JFIT::JPoint4D vx;
197  double chi2 = std::numeric_limits<double>::max();
198  int NDF = 0;
199  int N = 0;
200 
201  if(std::distance(dataL1.begin(), __end1) <= FACTORY_LIMIT){
202 
203  double ymin = std::numeric_limits<double>::max();
204 
205  JDataL1_t::iterator __end2 = __end1;
206 
207  for (int n = 0; n <= numberOfOutliers && std::distance(dataL1.begin(), __end2) >
209 
210  std::sort(dataL1.begin(), __end1, compare);
211 
212  do {
213  try {
214 
215  fit(dataL1.begin(), __end2);
216  double y = JFIT::getChi2(fit, dataL1.begin(), __end2, sigma_ns);
217 
218  if (y <= 0.0) {
219  std::cout << "chi2(1) " << y << std::endl; // <- it was used WARNING, which gives a bug.
220  } else if (y < ymin) {
221  ymin = y;
222  vx = fit;
223  NDF = std::distance(dataL1.begin(), __end2) - JFIT::JEstimator<JFIT::JPoint4D>::NUMBER_OF_PARAMETERS;
224  N = JFIT::getCount(dataL1.begin(), __end2);
225  chi2 = y;
226  }
227  }
228  catch(JLANG::JException& error) {}
229 
230  } while (JTOOLS::next_permutation(dataL1.begin(), __end2, __end1, compare));
231 
232  ymin -= STANDARD_DEVIATIONS * STANDARD_DEVIATIONS;
233  }
234 
235  } else {
236 
237  try{
238  fit(dataL1.begin(), __end1);
239  NDF = std::distance(dataL1.begin(), __end1) - JEstimator<JFIT::JPoint4D>::NUMBER_OF_PARAMETERS;
240  N = JFIT::getCount(dataL1.begin(), __end1);
241  vx = fit;
242  chi2 = JFIT::getChi2(fit, dataL1.begin(), __end1, sigma_ns);
243  }
244  catch(const JLANG::JException& error){}
245  }
246 
247  if (NDF >= 0) {
248 
249  for(double x = pos_grid.getLowerLimit(); x <= pos_grid.getUpperLimit(); x += pos_step){
250  for(double y = pos_grid.getLowerLimit(); y <= pos_grid.getUpperLimit(); y += pos_step){
251  for(double z = pos_grid.getLowerLimit(); z <= pos_grid.getUpperLimit(); z += pos_step){
252  for(double t = time_grid.getLowerLimit(); t <= time_grid.getUpperLimit(); t += time_step){
253 
254  const JGEOMETRY3D::JVector3D pos_shifted(vx.getX() + x, vx.getY() + y, vx.getZ() + z);
255  const JFIT::JPoint4D point_shifted(pos_shifted, vx.getT() + t);
256 
257  dataL0_selected.clear();
258 
259  //select a subcluster of L0s
260  const JFIT::JModel<JFIT::JPoint4D> match(vx, 80, JTOOLS::JTimeRange(-100, +100));
261 
262  for (JDataL0_t::const_iterator i = dataL0.begin(); i != dataL0.end(); ++i) {
263 
264  const JHitL0 hit(*i);
265 
266  if(match(hit)){
267  dataL0_selected.push_back(hit);
268  }
269  }
270 
271  double simplex_step = 2; // [m]
272  fit2_->step.resize(4);
273  fit2_->step[0] = JFIT::JPoint4D(JVector3D(simplex_step, 0.0, 0.0), 0.0);
274  fit2_->step[1] = JFIT::JPoint4D(JVector3D(0.0, simplex_step, 0.0), 0.0);
275  fit2_->step[2] = JFIT::JPoint4D(JVector3D(0.0, 0.0, simplex_step), 0.0);
276  fit2_->step[3] = JFIT::JPoint4D(JVector3D(0.0, 0.0, 0.0), 20);
277 
278  chi2 = (*fit2_)(point_shifted, dataL0_selected.begin(), dataL0_selected.end());
279 
280  const JGEOMETRY3D::JVector3D fitPos(fit2_->value);
281  const JGEOMETRY3D::JTrack3D fitResult(fitPos, JGEOMETRY3D::JVersor3Z(), fit2_->value.getT());
282 
283  const double energy(0);
285  fitResult, getQuality(chi2, (useL0 ? N : NDF)),
286  NDF, energy,JFIT::JFitStatus_t::OKAY );
287 
288  OutFits.push_back(outFit);
289  OutFits.rbegin()->setW(13, chi2);
290  OutFits.rbegin()->setW(14, NDF);
291 
292  }
293  }
294  }
295  }
296  }
297  std::sort(OutFits.begin(), OutFits.end(), qualitySorter);
298  }else{
299  std::cout<<"Too few hits " << std::endl;// <- it was used DEBUG, which gives a bug.
300  }
301 
302 }
General exception.
Definition: JException.hh:40
3G match criterion.
Definition: JMatch3G.hh:29
T getLowerLimit() const
Get lower limit.
Definition: JRange.hh:180
Data structure for L1 hit.
Definition: JHitL1.hh:34
JBinder2nd< JHit_t > JBind2nd(const JMatch< JHit_t > &match, const JHit_t &second)
Auxiliary method to create JBinder2nd object.
Definition: JBind2nd.hh:70
JFIT::JShowerPrefitParameters_t parameters_
Data structure for vertex fit.
Definition: JPoint4D.hh:22
int getCount(const JHitL0 &hit)
Get hit count.
Definition: JEvtToolkit.hh:728
JShowerPrefit.cc.
Container for historical events.
Definition: JHistory.hh:95
Template definition of linear fit.
Definition: JEstimator.hh:25
JHitIterator_t clusterizeWeight(JHitIterator_t __begin, JHitIterator_t __end, const JMatch< JHit_t > &match)
Partition data according given binary match operator.
Definition: JAlgorithm.hh:310
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
Data structure for track fit results.
Definition: JEvt.hh:32
JLANG::JSharedPointer< const JDETECTOR::JModuleRouter > moduleRouter_
Auxiliary class for permutations of L1 hits.
Definition: JEvtToolkit.hh:501
JHitIterator_t clusterize(JHitIterator_t __begin, JHitIterator_t __end, const JMatch< JHit_t > &match, const int Nmin=1)
Partition data according given binary match operator.
Definition: JAlgorithm.hh:45
Template L2 builder.
Definition: JBuildL2.hh:47
Data structure for vector in three dimensions.
Definition: JVector3D.hh:32
T getUpperLimit() const
Get upper limit.
Definition: JRange.hh:191
std::shared_ptr< JRegressor_t > fit2_
Auxiliary class to match data points with given model.
Definition: JModel.hh:24
double getQuality(const double chi2, const int NDF)
Get quality of fit.
Definition: JEvtToolkit.hh:195
void copy(const Head &from, JHead &to)
Copy header from from to to.
Definition: JHead.cc:40
bool qualitySorter(const JFit &first, const JFit &second)
Comparison of fit results.
Definition: JEvtToolkit.hh:223
bool next_permutation(T __begin, T __last, T __end, JComparator_t compare, std::bidirectional_iterator_tag)
Implementation of method next_permutation for bidirectional iterators.
Definition: JPermutation.hh:20
Template L0 hit builder.
Definition: JBuildL0.hh:35
double getChi2(const double P)
Get chi2 corresponding to given probability.
Definition: JFitToolkit.hh:80
Data structure for normalised vector in positive z-direction.
Definition: JVersor3Z.hh:36
static bool JFIT::JORCAShowerPrefit::totW ( const JTRIGGER::JHitR1 first,
const JTRIGGER::JHitR1 second 
)
inlinestatic

Member function used in the predicate of std::sort, to sort hits according to ToT and W.

Parameters
firstlhs hit
secondrhs hit

Definition at line 113 of file JORCAShowerPrefit.hh.

115  {
116  return (first.getToT() + first.getW() > second.getToT() + second.getW());
117  }
double getToT() const
Get calibrated time over threshold of hit.
Definition: JHit.hh:157
double getW() const
Get weight.
Definition: JHitR1.hh:194

Member Data Documentation

std::shared_ptr<JRegressor_t> JFIT::JORCAShowerPrefit::fit2_
private

Definition at line 59 of file JORCAShowerPrefit.hh.

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

Definition at line 61 of file JORCAShowerPrefit.hh.

JFIT::JShowerPrefitParameters_t JFIT::JORCAShowerPrefit::parameters_
private

Definition at line 62 of file JORCAShowerPrefit.hh.


The documentation for this class was generated from the following file: