Jpp master_rocky-44-g75b7c4f75
the software that should make you happy
Loading...
Searching...
No Matches
JRECONSTRUCTION::JMuonPrefit Struct Reference

Wrapper class to make pre-fit of muon trajectory. More...

#include <JMuonPrefit.hh>

Inheritance diagram for JRECONSTRUCTION::JMuonPrefit:
JRECONSTRUCTION::JMuonPrefitParameters_t JFIT::JEstimator< JModel_t > TObject

Classes

struct  cmz
 Auxiliary data structure for sorting of hits. More...
 

Public Types

typedef JEstimator< JLine1ZJEstimator_t
 
typedef JTRIGGER::JHitR1 hit_type
 
typedef std::vector< hit_typebuffer_type
 

Public Member Functions

 JMuonPrefit (const JMuonPrefitParameters_t &parameters, const JModuleRouter &router, const int debug=0)
 Constructor.
 
 JMuonPrefit (const JMuonPrefitParameters_t &parameters, const JModuleRouter &router, const JOmega3D &omega, const int debug=0)
 Constructor.
 
JEvt operator() (const KM3NETDAQ::JDAQEvent &event)
 Fit function.
 
JEvt operator() (const buffer_type &dataL0, const buffer_type &dataL1)
 Fit function.
 
void reset ()
 Reset fit parameters.
 
bool equals (const JMuonPrefitParameters_t &parameters) const
 Equality.
 
 ClassDef (JMuonPrefitParameters_t, 1)
 

Public Attributes

const JModuleRouterrouter
 
JOmega3D omega
 
int debug
 
int factoryLimit
 factory limit for combinatorics
 
int NMaxHits
 maximal number of hits
 
double sigma_ns
 time resolution [ns]
 
double gridAngle_deg
 grid angle for directions [deg]
 
bool useL0
 option for L0 hit use
 
int numberOfOutliers
 maximum number of outliers
 
size_t numberOfPrefits
 number of prefits
 
double DZMax
 maximal slope for downward pointing solutions
 
size_t numberOfDZMax
 additional number of downward pointing solutions
 
double TMaxLocal_ns
 time window for local coincidences [ns]
 
double ctMin
 minimal cosine space angle between PMT axes
 
double roadWidth_m
 road width [m]
 
double Qmin
 minimal quality step
 

Static Public Attributes

static const struct JRECONSTRUCTION::JMuonPrefit::cmz cmz
 

Private Member Functions

void configure ()
 Configure internal buffer(s).
 

Private Attributes

buffer_type data
 
JMatrixNZ V
 
JVectorNZ Y
 

Detailed Description

Wrapper class to make pre-fit of muon trajectory.

The JMuonPrefit fit is used to generate start values for subsequent fits (usually JMuonSimplex and JMuonGandalf).
To this end, a scan of directions is made and the time and transverse positions of the track are fitted for each direction (JFIT::JEstimator<JLine1Z>).
The directions are spaced by the parameters JMuonPrefitParameters_t::gridAngle_deg.
This angle corresponds to the envisaged angular accuracy of the result.
The probability that one of the results is less than this angle away from the correct value, multiple start values should be considered (JMuonPrefitParameters_t::numberOfPrefits).
Note that the CPU time scales with the inverse of the square of this angle.
The chi-squared is based on the time residuals.

Definition at line 75 of file JMuonPrefit.hh.

Member Typedef Documentation

◆ JEstimator_t

◆ hit_type

◆ buffer_type

Constructor & Destructor Documentation

◆ JMuonPrefit() [1/2]

JRECONSTRUCTION::JMuonPrefit::JMuonPrefit ( const JMuonPrefitParameters_t & parameters,
const JModuleRouter & router,
const int debug = 0 )
inline

Constructor.

Parameters
parametersparameters
routermodule router
debugdebug

Definition at line 92 of file JMuonPrefit.hh.

94 :
95 JMuonPrefitParameters_t(parameters),
97 omega (parameters.gridAngle_deg * JMATH::PI/180.0),
98 debug (debug)
99 {
100 configure();
101 }
static const double PI
Mathematical constants.
double gridAngle_deg
grid angle for directions [deg]
const JModuleRouter & router
void configure()
Configure internal buffer(s).

◆ JMuonPrefit() [2/2]

JRECONSTRUCTION::JMuonPrefit::JMuonPrefit ( const JMuonPrefitParameters_t & parameters,
const JModuleRouter & router,
const JOmega3D & omega,
const int debug = 0 )
inline

Constructor.

Parameters
routermodule router
parametersparameters
omegadirections
debugdebug

Definition at line 112 of file JMuonPrefit.hh.

115 :
116 JMuonPrefitParameters_t(parameters),
117 router(router),
118 omega (omega),
119 debug (debug)
120 {
121 configure();
122 }

Member Function Documentation

◆ operator()() [1/2]

JEvt JRECONSTRUCTION::JMuonPrefit::operator() ( const KM3NETDAQ::JDAQEvent & event)
inline

Fit function.

Parameters
eventevent
Returns
fit results

Definition at line 131 of file JMuonPrefit.hh.

132 {
133 using namespace std;
134 using namespace JPP;
135
136 const JBuildL0<hit_type> buildL0;
138
139 buffer_type dataL0;
140 buffer_type dataL1;
141
142 buildL2(event, router, !useL0, back_inserter(dataL1));
143
144 // 3D cluster of unique optical modules
145
147
148 sort(dataL1.begin(), dataL1.end(), hit_type::compare);
149
150 buffer_type::iterator __end = dataL1.end();
151
152 __end = unique(dataL1.begin(), __end, equal_to<JDAQModuleIdentifier>());
153
154 __end = clusterizeWeight(dataL1.begin(), __end, match3B);
155
156 dataL1.erase(__end, dataL1.end());
157
158
159 if (useL0) {
160
161 buildL0(event, router, true, back_inserter(dataL0));
162
163 __end = dataL0.end();
164
165 for (buffer_type::iterator i = dataL0.begin(); i != __end; ) {
166
167 if (match3B.count(*i, dataL1.begin(), dataL1.end()) != 0)
168 ++i;
169 else
170 swap(*i, *--__end);
171 }
172
173 dataL0.erase(__end, dataL0.end());
174 }
175
176 return (*this)(dataL0, dataL1);
177 }
Template L0 hit builder.
Definition JBuildL0.hh:38
Template L2 builder.
Definition JBuildL2.hh:49
static const struct JTRIGGER::JHitR1::compare compare
3D match criterion with road width.
Definition JMatch3B.hh:36
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
static const struct JTRIGGER::clusterizeWeight clusterizeWeight
double ctMin
minimal cosine space angle between PMT axes
double TMaxLocal_ns
time window for local coincidences [ns]
Data structure for L2 parameters.

◆ operator()() [2/2]

JEvt JRECONSTRUCTION::JMuonPrefit::operator() ( const buffer_type & dataL0,
const buffer_type & dataL1 )
inline

Fit function.

Parameters
dataL0L0 hit data
dataL1L1 hit data
Returns
fit results

Definition at line 187 of file JMuonPrefit.hh.

189 {
190 using namespace std;
191 using namespace JPP;
192
193 const double STANDARD_DEVIATIONS = 3.0; // [unit]
194 const double HIT_OFF = 1.0e3 * sigma_ns * sigma_ns; // [ns^2]
195
197
198 data.reserve(dataL0.size() +
199 dataL1.size());
200
201 JEvt out;
202
203 for (JOmega3D_t::const_iterator dir = omega.begin(); dir != omega.end(); ++dir) {
204
205 const JRotation3D R(*dir);
206
207
208 buffer_type::iterator __end = copy(dataL1.begin(), dataL1.end(), data.begin());
209
210 for (buffer_type::iterator i = data.begin(); i != __end; ++i) {
211 i->rotate(R);
212 }
213
214
215 // reduce data
216
217 if (distance(data.begin(), __end) > NMaxHits) {
218
219 buffer_type::iterator __p = data.begin();
220
221 advance(__p, NMaxHits);
222
223 partial_sort(data.begin(), __p, __end, cmz);
224
225 __end = __p;
226 }
227
228
229 // 1D cluster
230
231 __end = clusterizeWeight(data.begin(), __end, match1D);
232
233 if (useL0) {
234
235 buffer_type::iterator p = __end; // begin L0 data
236 buffer_type::iterator q = copy(dataL0.begin(), dataL0.end(), p); // end L0 data
237
238 for (buffer_type::iterator i = p; i != q; ++i) {
239
240 if (find_if(data.begin(), __end, make_predicate(&hit_type::getModuleID, i->getModuleID())) == __end) {
241
242 i->rotate(R);
243
244 if (match1D.count(*i, data.begin(), __end) != 0) {
245 *p = *i;
246 ++p;
247 }
248 }
249 }
250
251 __end = clusterize(__end, p, match1D);
252 }
253
254
255 if (distance(data.begin(), __end) <= NUMBER_OF_PARAMETERS) {
256 continue;
257 }
258
259
260 // 1D fit
261
262 JLine1Z tz;
263 double chi2 = numeric_limits<double>::max();
264 int NDF = distance(data.begin(), __end) - NUMBER_OF_PARAMETERS;
265 int N = getCount(data.begin(), __end);
266
267
268 if (distance(data.begin(), __end) <= factoryLimit) {
269
270 int number_of_outliers = numberOfOutliers;
271
272 if (number_of_outliers > NDF - 1) {
273 number_of_outliers = NDF - 1;
274 }
275
276 double ymin = numeric_limits<double>::max();
277
278 buffer_type::iterator __end1 = __end;
279
280 for (int n = 0; n <= number_of_outliers; ++n, --__end1) {
281
282 sort(data.begin(), __end, hit_type::compare);
283
284 do {
285 /*
286 if (getNumberOfStrings(router, data.begin(), __end1) < 2) {
287 continue;
288 }
289 */
290 try {
291
292 (*this)(data.begin(), __end1);
293
294 V.set(*this, data.begin(), __end1, gridAngle_deg, sigma_ns);
295 Y.set(*this, data.begin(), __end1);
296
297 V.invert();
298
299 double y = getChi2(Y, V);
300
301 if (y <= -(STANDARD_DEVIATIONS * STANDARD_DEVIATIONS)) {
302
303 WARNING(endl << "chi2(1) " << y << endl);
304
305 } else {
306
307 if (y < 0.0) {
308 y = 0.0;
309 }
310
311 if (y < ymin) {
312 ymin = y;
313 tz = *this;
314 chi2 = ymin;
315 NDF = distance(data.begin(), __end1) - NUMBER_OF_PARAMETERS;
316 N = getCount(data.begin(), __end1);
317 }
318 }
319 }
320 catch(const exception& error) {}
321
322 } while (next_permutation(data.begin(), __end1, __end, hit_type::compare));
323
324 ymin -= STANDARD_DEVIATIONS * STANDARD_DEVIATIONS;
325 }
326
327 } else {
328
329 const int number_of_outliers = NDF - 1;
330
331 try {
332
333 (*this)(data.begin(), __end);
334
335 V.set(*this, data.begin(), __end, gridAngle_deg, sigma_ns);
336 Y.set(*this, data.begin(), __end);
337
338 V.invert();
339
340 for (int n = 0; n <= number_of_outliers; ++n) {
341
342 double ymax = 0.0;
343 int k = -1;
344
345 for (size_t i = 0; i != Y.size(); ++i) {
346
347 double y = getChi2(Y, V, i);
348
349 if (y > ymax) {
350 ymax = y;
351 k = i;
352 }
353 }
354
355 if (ymax < STANDARD_DEVIATIONS * STANDARD_DEVIATIONS) {
356 break;
357 }
358
359 try {
360
361 V.update(k, +HIT_OFF);
362
363 this->update(data.begin(), __end, V);
364
365 Y.set(*this, data.begin(), __end);
366
367 tz = *this;
368 NDF -= 1;
369 N -= getCount(data[k]);
370 }
371 catch(const exception& error) {
372
373 V.update(k, -HIT_OFF);
374
375 static_cast<JLine1Z&>(*this) = tz;
376
377 Y.set(*this, data.begin(), __end);
378
379 break;
380 }
381 }
382
383 chi2 = getChi2(Y, V);
384 tz = *this;
385 }
386 catch(const exception& error) {}
387 }
388
389 if (chi2 != numeric_limits<double>::max()) {
390
391 tz.rotate_back(R);
392
393 out.push_back(getFit(JHistory(JMUONPREFIT), tz, *dir, getQuality(chi2, N, NDF), NDF));
394 }
395 }
396
397
398 if (numberOfPrefits > 0) {
399
400 JEvt::iterator __end = out.end();
401
402 if (Qmin > 0) { // sort distinct maxima
403
404 __end = out.begin(); // output
405
406 for (JEvt::iterator p1 = out.begin(); p1 != out.end() && (size_t) distance(out.begin(), __end) < numberOfPrefits; ) {
407
408 JEvt::iterator p2 = p1;
409
410 for (JEvt::iterator i = p1; i != out.end(); ++i) {
411 if (i->getQ() > p2->getQ()) {
412 p2 = i;
413 }
414 }
415
416 swap(*p2, *p1);
417
418 p2 = p1++;
419
420 sort(p1, out.end(), JPointing(*p2));
421
422 for (double Q = p2->getQ(); p1 != out.end() && (p1->getQ() >= p2->getQ() - Qmin || p1->getQ() <= Q); Q = (p1++)->getQ()) {}
423
424 swap(*(__end++), *p2);
425 }
426
427 } else if (numberOfPrefits < out.size()) { // sort subset
428
429 advance(__end = out.begin(), numberOfPrefits);
430
431 partial_sort(out.begin(), __end, out.end(), qualitySorter);
432
433 } else { // sort all
434
435 sort(out.begin(), __end, qualitySorter);
436 }
437
438 // add downward pointing solutions if available but not yet sufficient
439
440 int nz = numberOfDZMax - count_if(out.begin(), __end, make_predicate(&JFit::getDZ, DZMax, JComparison::le()));
441
442 if (nz > 0) {
443
444 JEvt::iterator __p = __end;
445 JEvt::iterator __q = __end = partition(__p, out.end(), make_predicate(&JFit::getDZ, DZMax, JComparison::le()));
446
447 if (nz < distance(__p, __q)) {
448
449 advance(__end = __p, nz);
450
451 partial_sort(__p, __end, __q, qualitySorter);
452
453 } else {
454
455 sort(__p, __end, qualitySorter);
456 }
457 }
458
459 out.erase(__end, out.end());
460
461 } else {
462
463 sort(out.begin(), out.end(), qualitySorter);
464 }
465
466 return out;
467 }
TPaveText * p1
#define WARNING(A)
Definition JMessage.hh:65
std::vector< T >::difference_type distance(typename std::vector< T >::const_iterator first, typename PhysicsEvent::const_iterator< T > second)
Specialisation of STL distance.
double getDZ() const
Get Z-slope.
Data structure for fit of straight line paralel to z-axis.
Definition JLine1Z.hh:29
void set(const JVector3D &pos, T __begin, T __end, const double alpha, const double sigma)
Set co-variance matrix.
Definition JMatrixNZ.hh:85
void set(const JLine1Z &track, T __begin, T __end)
Set time residual vector.
Definition JVectorNZ.hh:68
JPosition3D & rotate_back(const JRotation3D &R)
Rotate back.
Auxiliary class to compare fit results with respect to a reference direction (e.g....
1D match criterion.
Definition JMatch1D.hh:33
int getModuleID() const
Get module identifier.
static const int JMUONPREFIT
void copy(const Head &from, JHead &to)
Copy header from from to to.
Definition JHead.cc:162
double getChi2(const double P)
Get chi2 corresponding to given probability.
size_t getCount(const array_type< T > &buffer, const JCompare_t &compare)
Count number of unique values.
double getQuality(const double chi2, const int N, const int NDF)
Get quality of fit.
JFIT::JHistory JHistory
Definition JHistory.hh:354
bool qualitySorter(const JFit &first, const JFit &second)
Comparison of fit results.
JFit getFit(const JHistory &history, const JTrack3D &track, const double Q, const int NDF, const double energy=0.0, const int status=SINGLE_STAGE)
Get fit.
counter_type advance(counter_type &counter, const counter_type value, const counter_type limit=std::numeric_limits< counter_type >::max())
Advance counter.
const int n
Definition JPolint.hh:786
bool next_permutation(T __begin, T __last, T __end, JComparator_t compare, std::bidirectional_iterator_tag)
Implementation of method next_permutation for bidirectional iterators.
static const struct JTRIGGER::clusterize clusterize
Acoustic event fit.
void update(const size_t k, const double value)
Update inverted matrix at given diagonal element.
Definition JMatrixNS.hh:446
void invert()
Invert matrix according LDU decomposition.
Definition JMatrixNS.hh:75
double DZMax
maximal slope for downward pointing solutions
int factoryLimit
factory limit for combinatorics
size_t numberOfDZMax
additional number of downward pointing solutions
static const struct JRECONSTRUCTION::JMuonPrefit::cmz cmz

◆ configure()

void JRECONSTRUCTION::JMuonPrefit::configure ( )
inlineprivate

Configure internal buffer(s).

Definition at line 500 of file JMuonPrefit.hh.

501 {
502 using namespace JPP;
503
506 }
const JClass_t & getReference() const
Get reference to object.
Definition JReference.hh:38
int getNumberOfPMTs(const JModule &module)
Get number of PMTs.
static const JModuleCounter getNumberOfModules
Function object to count unique modules.

◆ reset()

void JRECONSTRUCTION::JMuonPrefitParameters_t::reset ( )
inlineinherited

Reset fit parameters.

Definition at line 39 of file JMuonPrefitParameters_t.hh.

40 {
41 factoryLimit = 8;
42 NMaxHits = 50;
43 sigma_ns = 5;
44 gridAngle_deg = 1;
45 useL0 = false;
47 numberOfPrefits = 12;
48 DZMax = 0.0;
49 numberOfDZMax = 1;
50 TMaxLocal_ns = 18.0;
51 ctMin = 0;
52 roadWidth_m = 200;
53 Qmin = 0.0;
54 }

◆ equals()

bool JRECONSTRUCTION::JMuonPrefitParameters_t::equals ( const JMuonPrefitParameters_t & parameters) const
inlineinherited

Equality.

Parameters
parametersfit parameters
Returns
true if equals; else false

Definition at line 62 of file JMuonPrefitParameters_t.hh.

63 {
64 return (this->factoryLimit == parameters.factoryLimit &&
65 this->NMaxHits == parameters.NMaxHits &&
66 this->sigma_ns == parameters.sigma_ns &&
67 this->gridAngle_deg == parameters.gridAngle_deg &&
68 this->useL0 == parameters.useL0 &&
69 this->numberOfOutliers == parameters.numberOfOutliers &&
70 this->numberOfPrefits == parameters.numberOfPrefits &&
71 this->DZMax == parameters.DZMax &&
72 this->numberOfDZMax == parameters.numberOfDZMax &&
73 this->TMaxLocal_ns == parameters.TMaxLocal_ns &&
74 this->roadWidth_m == parameters.roadWidth_m &&
75 this->Qmin == parameters.Qmin);
76 }

◆ ClassDef()

JRECONSTRUCTION::JMuonPrefitParameters_t::ClassDef ( JMuonPrefitParameters_t ,
1  )
inherited

Member Data Documentation

◆ cmz

const struct JRECONSTRUCTION::JMuonPrefit::cmz JRECONSTRUCTION::JMuonPrefit::cmz
static

◆ router

const JModuleRouter& JRECONSTRUCTION::JMuonPrefit::router

Definition at line 492 of file JMuonPrefit.hh.

◆ omega

JOmega3D JRECONSTRUCTION::JMuonPrefit::omega

Definition at line 493 of file JMuonPrefit.hh.

◆ debug

int JRECONSTRUCTION::JMuonPrefit::debug

Definition at line 494 of file JMuonPrefit.hh.

◆ data

buffer_type JRECONSTRUCTION::JMuonPrefit::data
private

Definition at line 508 of file JMuonPrefit.hh.

◆ V

JMatrixNZ JRECONSTRUCTION::JMuonPrefit::V
private

Definition at line 509 of file JMuonPrefit.hh.

◆ Y

JVectorNZ JRECONSTRUCTION::JMuonPrefit::Y
private

Definition at line 510 of file JMuonPrefit.hh.

◆ factoryLimit

int JRECONSTRUCTION::JMuonPrefitParameters_t::factoryLimit
inherited

factory limit for combinatorics

Definition at line 80 of file JMuonPrefitParameters_t.hh.

◆ NMaxHits

int JRECONSTRUCTION::JMuonPrefitParameters_t::NMaxHits
inherited

maximal number of hits

Definition at line 81 of file JMuonPrefitParameters_t.hh.

◆ sigma_ns

double JRECONSTRUCTION::JMuonPrefitParameters_t::sigma_ns
inherited

time resolution [ns]

Definition at line 82 of file JMuonPrefitParameters_t.hh.

◆ gridAngle_deg

double JRECONSTRUCTION::JMuonPrefitParameters_t::gridAngle_deg
inherited

grid angle for directions [deg]

Definition at line 83 of file JMuonPrefitParameters_t.hh.

◆ useL0

bool JRECONSTRUCTION::JMuonPrefitParameters_t::useL0
inherited

option for L0 hit use

Definition at line 84 of file JMuonPrefitParameters_t.hh.

◆ numberOfOutliers

int JRECONSTRUCTION::JMuonPrefitParameters_t::numberOfOutliers
inherited

maximum number of outliers

Definition at line 85 of file JMuonPrefitParameters_t.hh.

◆ numberOfPrefits

size_t JRECONSTRUCTION::JMuonPrefitParameters_t::numberOfPrefits
inherited

number of prefits

Definition at line 86 of file JMuonPrefitParameters_t.hh.

◆ DZMax

double JRECONSTRUCTION::JMuonPrefitParameters_t::DZMax
inherited

maximal slope for downward pointing solutions

Definition at line 87 of file JMuonPrefitParameters_t.hh.

◆ numberOfDZMax

size_t JRECONSTRUCTION::JMuonPrefitParameters_t::numberOfDZMax
inherited

additional number of downward pointing solutions

Definition at line 88 of file JMuonPrefitParameters_t.hh.

◆ TMaxLocal_ns

double JRECONSTRUCTION::JMuonPrefitParameters_t::TMaxLocal_ns
inherited

time window for local coincidences [ns]

Definition at line 89 of file JMuonPrefitParameters_t.hh.

◆ ctMin

double JRECONSTRUCTION::JMuonPrefitParameters_t::ctMin
inherited

minimal cosine space angle between PMT axes

Definition at line 90 of file JMuonPrefitParameters_t.hh.

◆ roadWidth_m

double JRECONSTRUCTION::JMuonPrefitParameters_t::roadWidth_m
inherited

road width [m]

Definition at line 91 of file JMuonPrefitParameters_t.hh.

◆ Qmin

double JRECONSTRUCTION::JMuonPrefitParameters_t::Qmin
inherited

minimal quality step

Definition at line 92 of file JMuonPrefitParameters_t.hh.


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