HighLAND
ConstituentsUtils.cxx
1 #include "ConstituentsUtils.hxx"
2 #include <TVector3.h>
3 #include <stdio.h>
4 #include <math.h>
5 #include <iostream>
6 #include <typeinfo>
7 #include "BaseDataClasses.hxx"
8 
9 #include "SubDetId.hxx"
10 
11 #include <bitset>
12 
13 namespace anaUtils {}
14 
15 
16 //**************************************************
17 Int_t anaUtils::GetOneSegmentPerTPC(AnaTPCParticleB* in[], Int_t nseg_in, AnaTPCParticleB* out[]) {
18 //**************************************************
19 
20  // This method takes as input all TPC segments in a track and returns an array with only one segment per TPC, the one with more nodes
21 
22  Int_t nnodes_max[3]={0,0,0};
23  Int_t itrack_nnodes_max[3]={-1,-1,-1};
24  for(Int_t iseg=0; iseg< nseg_in; iseg++){
25  Int_t tpc = SubDetId::GetTPC(in[iseg]->Detector)-1;
26  if (in[iseg]->NNodes > nnodes_max[tpc]){
27  nnodes_max[tpc] = in[iseg]->NNodes;
28  itrack_nnodes_max[tpc]=iseg;
29  }
30  }
31 
32  int nseg_out = 0;
33  for(Int_t i=0;i< 3;i++){
34  if (itrack_nnodes_max[i]!=-1)
35  out[nseg_out++]=in[itrack_nnodes_max[i]];
36  }
37 
38  return nseg_out;
39 }
40 
41 //**************************************************
43 //**************************************************
44 
45  // returns the TPC closest to the track start point
46  // simply use the linear distance for now
47 
48  /*
49  SubDetId::SubDetEnum tpc_closest = SubDetId::kInvalid;
50 
51  Float_t dist = 9999999.;
52 
53  for(int i = 0; i < track.nTPCSegments; ++i){
54  AnaTPCParticleB* tpc_track = track.TPCSegments[i];
55  Float_t dist_tmp = GetSeparationSquared(track.PositionStart, tpc_track->PositionStart);
56  if(dist_tmp < dist){
57  dist = dist_tmp;
58  tpc_closest = SubDetId::GetSubdetectorEnum(tpc_track->Detector);
59  }
60  }
61  return tpc_closest;
62  */
63 
65 
66  if (TPCSegment) return SubDetId::GetSubdetectorEnum(TPCSegment->Detector);
67  else return SubDetId::kInvalid;
68 }
69 
70 //**************************************************
72 //**************************************************
73  if (det == SubDetId::kInvalid) {
74  return 0;
75  }
76  if(!SubDetId::GetDetectorUsed(track.Detector, det)){
77  return 0;
78  }
79 
80  int count = 0;
81 
82  // Return segments for complete detector subsystems (all TPC etc.) first
83  switch(det){
84  case SubDetId::kTPC :
85  std::copy(&track.TPCSegments[0], &track.TPCSegments[track.nTPCSegments], segments);
86  return track.nTPCSegments;
87  break;
88  case SubDetId::kFGD :
89  std::copy(&track.FGDSegments[0], &track.FGDSegments[track.nFGDSegments], segments);
90  return track.nFGDSegments;
91  break;
92 
93  case SubDetId::kECAL :
94  std::copy(&track.ECALSegments[0], &track.ECALSegments[track.nECALSegments], segments);
95  return track.nECALSegments;
96  break;
97 
98  case SubDetId::kSMRD :
99  std::copy(&track.SMRDSegments[0], &track.SMRDSegments[track.nSMRDSegments], segments);
100  return track.nSMRDSegments;
101  break;
102 
103  case SubDetId::kP0D :
104  std::copy(&track.P0DSegments[0], &track.P0DSegments[track.nP0DSegments], segments);
105  return track.nP0DSegments;
106  break;
107  default:
108 
109  if(SubDetId::IsTPCDetector(det)){
110  for(int i = 0; i < track.nTPCSegments; ++i){
111  AnaTPCParticleB* tpc_track = track.TPCSegments[i];
112  if(SubDetId::GetDetectorUsed(tpc_track->Detector, det)){
113  segments[count] = tpc_track;
114  ++count;
115  }
116  }
117  }
118  else if(SubDetId::IsFGDDetector(det)){
119  for(int i = 0; i < track.nFGDSegments; ++i){
120  AnaFGDParticleB* fgd_track = track.FGDSegments[i];
121  if(SubDetId::GetDetectorUsed(fgd_track->Detector, det)){
122  segments[count] = fgd_track;
123  ++count;
124  }
125  }
126  }
127  else if(SubDetId::IsECALDetector(det)){
128  for(int i = 0; i < track.nECALSegments; ++i){
129  AnaECALParticleB* ecal_track = track.ECALSegments[i];
130  if(SubDetId::GetDetectorUsed(ecal_track->Detector, det)){
131  segments[count] = ecal_track;
132  ++count;
133  }
134  }
135  }
136  else if(SubDetId::IsSMRDDetector(det)){
137  for(int i = 0; i < track.nSMRDSegments; ++i){
138  AnaSMRDParticleB* smrd_track = track.SMRDSegments[i];
139  if(SubDetId::GetDetectorUsed(smrd_track->Detector, det)){
140  segments[count] = smrd_track;
141  ++count;
142  }
143  }
144  }
145  return count;
146  }
147  return count;
148 }
149 
150 //**************************************************
152 //**************************************************
153 
155 }
156 
157 //**************************************************
159 //**************************************************
160 
161  int tpc_closest = SubDetId::kInvalid;
162  int tpc = SubDetId::kInvalid;
163 
164  AnaParticleB* subtrack[3] = {NULL, NULL, NULL};
165 
166  Float_t dist = 9999999.;
167  int nNodes[3] = {0,0,0};
168 
169  for(int i = 0; i < track.nTPCSegments; ++i){
170  AnaTPCParticleB* tpc_track = track.TPCSegments[i];
171  Float_t dist_tmp = std::min(
172  GetSeparationSquared(pos, tpc_track->PositionStart),
173  GetSeparationSquared(pos, tpc_track->PositionEnd)
174  );
175  tpc = SubDetId::GetTPC(tpc_track->Detector);
176 
177  if(dist_tmp < dist){
178  dist = dist_tmp;
179  tpc_closest = tpc;
180  }
181  // TPC number is not zero ordered
182  if(tpc_track->NNodes > nNodes[tpc-1]){
183  nNodes[tpc-1] = tpc_track->NNodes;
184  subtrack[tpc-1] = tpc_track;
185  }
186  }
187 
188  if(tpc_closest != (int)SubDetId::kInvalid) return subtrack[tpc_closest-1];
189 
190  return NULL;
191 }
192 
193 //**************************************************
195 //**************************************************
196 
197  if (det == SubDetId::kInvalid) {
198  return NULL;
199  }
200  if(!SubDetId::GetDetectorUsed(track.Detector, det)){
201  return NULL;
202  }
203 
204  int nNodes = 0;
205  AnaParticleB* subtrack = NULL;
206 
207  switch(det){
208  case SubDetId::kTPC :
209  for(int i = 0; i < track.nTPCSegments; ++i){
210  AnaTPCParticleB* tpc_track = track.TPCSegments[i];
211  if(tpc_track->NNodes > nNodes){
212  nNodes = tpc_track->NNodes;
213  subtrack = tpc_track;
214  }
215  }
216  return subtrack;
217  break;
218  case SubDetId::kFGD :
219  for(int i = 0; i < track.nFGDSegments; ++i){
220  AnaFGDParticleB* fgd_track = track.FGDSegments[i];
221  if(fgd_track->NNodes > nNodes){
222  nNodes = fgd_track->NNodes;
223  subtrack = fgd_track;
224  }
225  }
226  return subtrack;
227  break;
228  case SubDetId::kECAL :
229  for(int i = 0; i < track.nECALSegments; ++i){
230  AnaECALParticleB* tpc_track = track.ECALSegments[i];
231  if(tpc_track->NNodes > nNodes){
232  nNodes = tpc_track->NNodes;
233  subtrack = tpc_track;
234  }
235  }
236  return subtrack;
237  break;
238  case SubDetId::kSMRD :
239  for(int i = 0; i < track.nSMRDSegments; ++i){
240  AnaSMRDParticleB* smrd_track = track.SMRDSegments[i];
241  if(smrd_track->NNodes > nNodes){
242  nNodes = smrd_track->NNodes;
243  subtrack = smrd_track;
244  }
245  }
246  return subtrack;
247  break;
248  case SubDetId::kP0D :
249  for(int i = 0; i < track.nP0DSegments; ++i){
250  AnaP0DParticleB* p0d_track = track.P0DSegments[i];
251  if(p0d_track->NNodes > nNodes){
252  nNodes = p0d_track->NNodes;
253  subtrack = p0d_track;
254  }
255  }
256  return subtrack;
257  break;
258  default:
259  if(SubDetId::IsTPCDetector(det)){
260  for(int i = 0; i < track.nTPCSegments; ++i){
261  AnaTPCParticleB* tpc_track = track.TPCSegments[i];
262  if(SubDetId::GetDetectorUsed(tpc_track->Detector, det)){
263  if(tpc_track->NNodes > nNodes){
264  nNodes = tpc_track->NNodes;
265  subtrack = tpc_track;
266  }
267  }
268  }
269  return subtrack;
270  }
271  else if(SubDetId::IsFGDDetector(det)){
272  for(int i = 0; i < track.nFGDSegments; ++i){
273  AnaFGDParticleB* fgd_track = track.FGDSegments[i];
274  if(SubDetId::GetDetectorUsed(fgd_track->Detector, det)){
275  if(fgd_track->NNodes > nNodes){
276  nNodes = fgd_track->NNodes;
277  subtrack = fgd_track;
278  }
279  }
280  }
281  return subtrack;
282  }
283 
284  else if(SubDetId::IsP0DDetector(det)){
285  for(int i = 0; i < track.nP0DSegments; ++i){
286  AnaP0DParticleB* p0d_track = track.P0DSegments[i];
287  if(SubDetId::GetDetectorUsed(p0d_track->Detector, det)){
288  if(p0d_track->NNodes > nNodes){
289  nNodes = p0d_track->NNodes;
290  subtrack = p0d_track;
291  }
292  }
293  }
294  return subtrack;
295  }
296  else if(SubDetId::IsECALDetector(det)){
297  for(int i = 0; i < track.nECALSegments; ++i){
298  AnaECALParticleB* ecal_track = track.ECALSegments[i];
299  if(SubDetId::GetDetectorUsed(ecal_track->Detector, det)){
300  if(ecal_track->NNodes > nNodes){
301  nNodes = ecal_track->NNodes;
302  subtrack = ecal_track;
303  }
304  }
305  }
306  return subtrack;
307  }
308  else if(SubDetId::IsSMRDDetector(det)){
309  for(int i = 0; i < track.nSMRDSegments; ++i){
310  AnaSMRDParticleB* smrd_track = track.SMRDSegments[i];
311  if(SubDetId::GetDetectorUsed(smrd_track->Detector, det)){
312  if(smrd_track->NNodes > nNodes){
313  nNodes = smrd_track->NNodes;
314  subtrack = smrd_track;
315  }
316  }
317  }
318  return subtrack;
319  }
320  return NULL;
321  }
322 }
323 
324 //**************************************************
326 //**************************************************
327 
328  if(SubDetId::IsTPCDetector(det)){
329  for(int i = 0; i < track.nTPCSegments; ++i){
330  AnaTPCParticleB* tpc_track = track.TPCSegments[i];
331  if(SubDetId::GetDetectorUsed(tpc_track->Detector, det)){
332  return tpc_track;
333  }
334  }
335  }
336 
337  if(SubDetId::IsFGDDetector(det)){
338  for(int i = 0; i < track.nFGDSegments; ++i){
339  AnaFGDParticleB* fgd_track = track.FGDSegments[i];
340  if(SubDetId::GetDetectorUsed(fgd_track->Detector, det)){
341  return fgd_track;
342  }
343  }
344  }
345 
346  if(SubDetId::IsP0DDetector(det)){
347  for(int i = 0; i < track.nP0DSegments; ++i){
348  AnaP0DParticleB* p0d_track = track.P0DSegments[i];
349  if(SubDetId::GetDetectorUsed(p0d_track->Detector, det)){
350  return p0d_track;
351  }
352  }
353  }
354 
355  if(SubDetId::IsECALDetector(det)){
356  for(int i = 0; i < track.nECALSegments; ++i){
357  AnaECALParticleB* ecal_track = track.ECALSegments[i];
358  if(SubDetId::GetDetectorUsed(ecal_track->Detector, det)){
359  return ecal_track;
360  }
361  }
362  }
363 
364  if(SubDetId::IsSMRDDetector(det)){
365  for(int i = 0; i < track.nSMRDSegments; ++i){
366  AnaSMRDParticleB* smrd_track = track.SMRDSegments[i];
367  if(SubDetId::GetDetectorUsed(smrd_track->Detector, det)){
368  return smrd_track;
369  }
370  }
371  }
372  return NULL;
373 }
374 
375 //**************************************************
377 //**************************************************
378  for(Int_t it = 0; it != static_cast<Int_t>(SubDetId::kInvalidSubdetector); it++) {
379  SubDetId::SubDetEnum det = static_cast<SubDetId::SubDetEnum>(it);
380  if (anaUtils::InDetVolume(det, pos))
381  return det;
382  }
383  return SubDetId::kInvalid;
384 }
385 
386 //**************************************************
387 bool anaUtils::InDetVolume(SubDetId::SubDetEnum det, const Float_t* pos){
388 //**************************************************
389 
390  Float_t null[3] = {0.,0.,0.};
391 
392  //account for a case when a "general" volume is provided
393  switch(det){
394  case SubDetId::kFGD:
395  return (InFiducialVolume(SubDetId::kFGD1, pos, null, null) || InFiducialVolume(SubDetId::kFGD2, pos, null, null));
396  break;
397  case SubDetId::kFGD1:
398  return anaUtils::InFiducialVolume(SubDetId::kFGD1, pos, null, null);
399  break;
400  case SubDetId::kFGD2:
401  return anaUtils::InFiducialVolume(SubDetId::kFGD2, pos, null, null);
402  break;
403  case SubDetId::kTPC:
404  return (
405  InFiducialVolume(SubDetId::kTPC1, pos, null, null) ||
406  InFiducialVolume(SubDetId::kTPC2, pos, null, null) ||
407  InFiducialVolume(SubDetId::kTPC3, pos, null, null)
408  );
409  break;
410  case SubDetId::kTPC1:
411  return anaUtils::InFiducialVolume(SubDetId::kTPC1, pos, null, null);
412  break;
413  case SubDetId::kTPC2:
414  return anaUtils::InFiducialVolume(SubDetId::kTPC2, pos, null, null);
415  break;
416  case SubDetId::kTPC3:
417  return anaUtils::InFiducialVolume(SubDetId::kTPC3, pos, null, null);
418  break;
419  case SubDetId::kP0D:
420  return anaUtils::InFiducialVolume(SubDetId::kP0D, pos, null, null);
421  break;
422  case SubDetId::kECAL:
423  return (
424  anaUtils::InFiducialVolume(SubDetId::kTopTECAL, pos, null, null) ||
425  anaUtils::InFiducialVolume(SubDetId::kBottomTECAL, pos, null, null) ||
426  anaUtils::InFiducialVolume(SubDetId::kLeftTECAL, pos, null, null) ||
427  anaUtils::InFiducialVolume(SubDetId::kRightTECAL, pos, null, null) ||
428  anaUtils::InFiducialVolume(SubDetId::kTopPECAL, pos, null, null) ||
429  anaUtils::InFiducialVolume(SubDetId::kBottomPECAL, pos, null, null) ||
430  anaUtils::InFiducialVolume(SubDetId::kLeftPECAL, pos, null, null) ||
431  anaUtils::InFiducialVolume(SubDetId::kRightPECAL, pos, null, null) ||
432  anaUtils::InFiducialVolume(SubDetId::kDSECAL, pos, null, null)
433  );
434  break;
435  case SubDetId::kTECAL:
436  return (
437  anaUtils::InFiducialVolume(SubDetId::kTopTECAL, pos, null, null) ||
438  anaUtils::InFiducialVolume(SubDetId::kBottomTECAL, pos, null, null) ||
439  anaUtils::InFiducialVolume(SubDetId::kLeftTECAL, pos, null, null) ||
440  anaUtils::InFiducialVolume(SubDetId::kRightTECAL, pos, null, null)
441  );
442  break;
443  case SubDetId::kTopTECAL:
444  return anaUtils::InFiducialVolume(SubDetId::kTopTECAL, pos, null, null);
445  break;
446  case SubDetId::kBottomTECAL:
447  return anaUtils::InFiducialVolume(SubDetId::kBottomTECAL, pos, null, null);
448  break;
449  case SubDetId::kLeftTECAL:
450  return anaUtils::InFiducialVolume(SubDetId::kLeftTECAL, pos, null, null);
451  break;
452  case SubDetId::kRightTECAL:
453  return anaUtils::InFiducialVolume(SubDetId::kRightTECAL, pos, null, null);
454  break;
455  case SubDetId::kPECAL:
456  return (
457  anaUtils::InFiducialVolume(SubDetId::kTopPECAL, pos, null, null) ||
458  anaUtils::InFiducialVolume(SubDetId::kBottomPECAL, pos, null, null) ||
459  anaUtils::InFiducialVolume(SubDetId::kLeftPECAL, pos, null, null) ||
460  anaUtils::InFiducialVolume(SubDetId::kRightPECAL, pos, null, null)
461  );
462  break;
463  case SubDetId::kTopPECAL:
464  return anaUtils::InFiducialVolume(SubDetId::kTopPECAL, pos, null, null);
465  break;
466  case SubDetId::kBottomPECAL:
467  return anaUtils::InFiducialVolume(SubDetId::kBottomPECAL, pos, null, null);
468  break;
469  case SubDetId::kLeftPECAL:
470  return anaUtils::InFiducialVolume(SubDetId::kLeftPECAL, pos, null, null);
471  break;
472  case SubDetId::kRightPECAL:
473  return anaUtils::InFiducialVolume(SubDetId::kRightPECAL, pos, null, null);
474  break;
475  case SubDetId::kSMRD:
476  return (
477  anaUtils::InFiducialVolume(SubDetId::kTopSMRD, pos, null, null) ||
478  anaUtils::InFiducialVolume(SubDetId::kBottomSMRD, pos, null, null) ||
479  anaUtils::InFiducialVolume(SubDetId::kLeftSMRD, pos, null, null) ||
480  anaUtils::InFiducialVolume(SubDetId::kRightSMRD, pos, null, null)
481  );
482  break;
483  case SubDetId::kTopSMRD:
484  return anaUtils::InFiducialVolume(SubDetId::kTopSMRD, pos, null, null);
485  break;
486  case SubDetId::kBottomSMRD:
487  return anaUtils::InFiducialVolume(SubDetId::kBottomSMRD, pos, null, null);
488  break;
489  case SubDetId::kLeftSMRD:
490  return anaUtils::InFiducialVolume(SubDetId::kLeftSMRD, pos, null, null);
491  break;
492  case SubDetId::kRightSMRD:
493  return anaUtils::InFiducialVolume(SubDetId::kRightSMRD, pos, null, null);
494  break;
495  case SubDetId::kDSECAL:
496  return anaUtils::InFiducialVolume(SubDetId::kDSECAL, pos, null, null);
497  break;
498  default:
499  std::cout << "Warning: anaUtils::InDetVolume() No Volume set for " << det << std::endl;
500  return false;
501  break;
502  }
503 
504 }
505 
506 //**************************************************
507 bool anaUtils::InFiducialVolume(SubDetId::SubDetEnum det, const Float_t* pos){
508 //**************************************************
509 
510 
511 
512 
513 
514  Float_t null[3] = {0.,0.,0.};
515  switch(det){
516  case SubDetId::kFGD:
517  return (InFiducialVolume(SubDetId::kFGD1,pos) || InFiducialVolume(SubDetId::kFGD2,pos));
518  break;
519  case SubDetId::kFGD1:
520  return anaUtils::InFiducialVolume(SubDetId::kFGD1, pos, FVDef::FVdefminFGD1, FVDef::FVdefmaxFGD1);
521  break;
522  case SubDetId::kFGD2:
523  return anaUtils::InFiducialVolume(SubDetId::kFGD2, pos, FVDef::FVdefminFGD2, FVDef::FVdefmaxFGD2);
524  break;
525  case SubDetId::kTPC:
526  return (
527  InFiducialVolume(SubDetId::kTPC1, pos, null, null) ||
528  InFiducialVolume(SubDetId::kTPC2, pos, null, null) ||
529  InFiducialVolume(SubDetId::kTPC3, pos, null, null)
530  );
531  break;
532  case SubDetId::kTPC1:
533  return anaUtils::InFiducialVolume(SubDetId::kTPC1, pos, null, null);
534  break;
535  case SubDetId::kTPC2:
536  return anaUtils::InFiducialVolume(SubDetId::kTPC2, pos, null, null);
537  break;
538  case SubDetId::kTPC3:
539  return anaUtils::InFiducialVolume(SubDetId::kTPC3, pos, null, null);
540  break;
541  case SubDetId::kP0D:
542  return anaUtils::InFiducialVolume(SubDetId::kP0D, pos, FVDef::FVdefminP0D, FVDef::FVdefmaxP0D);
543  break;
544  case SubDetId::kECAL:
545  return (
546  anaUtils::InFiducialVolume(SubDetId::kTopTECAL, pos, null, null) ||
547  anaUtils::InFiducialVolume(SubDetId::kBottomTECAL, pos, null, null) ||
548  anaUtils::InFiducialVolume(SubDetId::kLeftTECAL, pos, null, null) ||
549  anaUtils::InFiducialVolume(SubDetId::kRightTECAL, pos, null, null) ||
550  anaUtils::InFiducialVolume(SubDetId::kTopPECAL, pos, null, null) ||
551  anaUtils::InFiducialVolume(SubDetId::kBottomPECAL, pos, null, null) ||
552  anaUtils::InFiducialVolume(SubDetId::kLeftPECAL, pos, null, null) ||
553  anaUtils::InFiducialVolume(SubDetId::kRightPECAL, pos, null, null) ||
554  anaUtils::InFiducialVolume(SubDetId::kDSECAL, pos, null, null)
555  );
556  break;
557  case SubDetId::kTECAL:
558  return (
559  anaUtils::InFiducialVolume(SubDetId::kTopTECAL, pos, null, null) ||
560  anaUtils::InFiducialVolume(SubDetId::kBottomTECAL, pos, null, null) ||
561  anaUtils::InFiducialVolume(SubDetId::kLeftTECAL, pos, null, null) ||
562  anaUtils::InFiducialVolume(SubDetId::kRightTECAL, pos, null, null)
563  );
564  break;
565  case SubDetId::kTopTECAL:
566  return anaUtils::InFiducialVolume(SubDetId::kTopTECAL, pos, null, null);
567  break;
568  case SubDetId::kBottomTECAL:
569  return anaUtils::InFiducialVolume(SubDetId::kBottomTECAL, pos, null, null);
570  break;
571  case SubDetId::kLeftTECAL:
572  return anaUtils::InFiducialVolume(SubDetId::kLeftTECAL, pos, null, null);
573  break;
574  case SubDetId::kRightTECAL:
575  return anaUtils::InFiducialVolume(SubDetId::kRightTECAL, pos, null, null);
576  break;
577  case SubDetId::kPECAL:
578  return (
579  anaUtils::InFiducialVolume(SubDetId::kTopPECAL, pos, null, null) ||
580  anaUtils::InFiducialVolume(SubDetId::kBottomPECAL, pos, null, null) ||
581  anaUtils::InFiducialVolume(SubDetId::kLeftPECAL, pos, null, null) ||
582  anaUtils::InFiducialVolume(SubDetId::kRightPECAL, pos, null, null)
583  );
584  break;
585  case SubDetId::kTopPECAL:
586  return anaUtils::InFiducialVolume(SubDetId::kTopPECAL, pos, null, null);
587  break;
588  case SubDetId::kBottomPECAL:
589  return anaUtils::InFiducialVolume(SubDetId::kBottomPECAL, pos, null, null);
590  break;
591  case SubDetId::kLeftPECAL:
592  return anaUtils::InFiducialVolume(SubDetId::kLeftPECAL, pos, null, null);
593  break;
594  case SubDetId::kRightPECAL:
595  return anaUtils::InFiducialVolume(SubDetId::kRightPECAL, pos, null, null);
596  break;
597  case SubDetId::kSMRD:
598  return (
599  anaUtils::InFiducialVolume(SubDetId::kTopSMRD, pos, null, null) ||
600  anaUtils::InFiducialVolume(SubDetId::kBottomSMRD, pos, null, null) ||
601  anaUtils::InFiducialVolume(SubDetId::kLeftSMRD, pos, null, null) ||
602  anaUtils::InFiducialVolume(SubDetId::kRightSMRD, pos, null, null)
603  );
604  break;
605  case SubDetId::kTopSMRD:
606  return anaUtils::InFiducialVolume(SubDetId::kTopSMRD, pos, null, null);
607  break;
608  case SubDetId::kBottomSMRD:
609  return anaUtils::InFiducialVolume(SubDetId::kBottomSMRD, pos, null, null);
610  break;
611  case SubDetId::kLeftSMRD:
612  return anaUtils::InFiducialVolume(SubDetId::kLeftSMRD, pos, null, null);
613  break;
614  case SubDetId::kRightSMRD:
615  return anaUtils::InFiducialVolume(SubDetId::kRightSMRD, pos, null, null);
616  break;
617  case SubDetId::kDSECAL:
618  return anaUtils::InFiducialVolume(SubDetId::kDSECAL, pos, null, null);
619  break;
620  default:
621  std::cout << "Warning: anaUtils::InFiducialVolume() No Fiducial Volume set for " << det << std::endl;
622  return false;
623  break;
624  }
625 }
626 
627 //**************************************************
628 bool anaUtils::InFiducialVolume(SubDetId::SubDetEnum det, const Float_t* pos, const Float_t* FVdefmin, const Float_t* FVdefmax){
629 //**************************************************
630 
631  switch(det){
632  case SubDetId::kFGD1:
633  if (pos[0] > DetDef::fgd1min[0]+FVdefmin[0] &&
634  pos[0] < DetDef::fgd1max[0]-FVdefmax[0] &&
635  pos[1] > DetDef::fgd1min[1]+FVdefmin[1] &&
636  pos[1] < DetDef::fgd1max[1]-FVdefmax[1] &&
637  pos[2] > DetDef::fgd1min[2]+FVdefmin[2] &&
638  pos[2] < DetDef::fgd1max[2]-FVdefmax[2])
639  return true;
640  break;
641  case SubDetId::kFGD2:
642  if (pos[0] > DetDef::fgd2min[0]+FVdefmin[0] &&
643  pos[0] < DetDef::fgd2max[0]-FVdefmax[0] &&
644  pos[1] > DetDef::fgd2min[1]+FVdefmin[1] &&
645  pos[1] < DetDef::fgd2max[1]-FVdefmax[1] &&
646  pos[2] > DetDef::fgd2min[2]+FVdefmin[2] &&
647  pos[2] < DetDef::fgd2max[2]-FVdefmax[2])
648  return true;
649  break;
650  case SubDetId::kTPC1:
651  if (pos[0] > DetDef::tpc1min[0]+FVdefmin[0] &&
652  pos[0] < DetDef::tpc1max[0]-FVdefmax[0] &&
653  pos[1] > DetDef::tpc1min[1]+FVdefmin[1] &&
654  pos[1] < DetDef::tpc1max[1]-FVdefmax[1] &&
655  pos[2] > DetDef::tpc1min[2]+FVdefmin[2] &&
656  pos[2] < DetDef::tpc1max[2]-FVdefmax[2])
657  return true;
658  break;
659  case SubDetId::kTPC2:
660  if (pos[0] > DetDef::tpc2min[0]+FVdefmin[0] &&
661  pos[0] < DetDef::tpc2max[0]-FVdefmax[0] &&
662  pos[1] > DetDef::tpc2min[1]+FVdefmin[1] &&
663  pos[1] < DetDef::tpc2max[1]-FVdefmax[1] &&
664  pos[2] > DetDef::tpc2min[2]+FVdefmin[2] &&
665  pos[2] < DetDef::tpc2max[2]-FVdefmax[2])
666  return true;
667  break;
668  case SubDetId::kTPC3:
669  if (pos[0] > DetDef::tpc3min[0]+FVdefmin[0] &&
670  pos[0] < DetDef::tpc3max[0]-FVdefmax[0] &&
671  pos[1] > DetDef::tpc3min[1]+FVdefmin[1] &&
672  pos[1] < DetDef::tpc3max[1]-FVdefmax[1] &&
673  pos[2] > DetDef::tpc3min[2]+FVdefmin[2] &&
674  pos[2] < DetDef::tpc3max[2]-FVdefmax[2])
675  return true;
676  break;
677  case SubDetId::kP0D:
678  if (pos[0] > DetDef::p0dmin[0]+FVdefmin[0] &&
679  pos[0] < DetDef::p0dmax[0]-FVdefmax[0] &&
680  pos[1] > DetDef::p0dmin[1]+FVdefmin[1] &&
681  pos[1] < DetDef::p0dmax[1]-FVdefmax[1] &&
682  pos[2] > DetDef::p0dmin[2]+FVdefmin[2] &&
683  pos[2] < DetDef::p0dmax[2]-FVdefmax[2])
684  return true;
685  break;
686  //DsECal (Ecal)
687  case SubDetId::kDSECAL:
688  if (pos[0] > DetDef::dsecalmin[0]+FVdefmin[0] &&
689  pos[0] < DetDef::dsecalmax[0]-FVdefmax[0] &&
690  pos[1] > DetDef::dsecalmin[1]+FVdefmin[1] &&
691  pos[1] < DetDef::dsecalmax[1]-FVdefmax[1] &&
692  pos[2] > DetDef::dsecalmin[2]+FVdefmin[2] &&
693  pos[2] < DetDef::dsecalmax[2]-FVdefmax[2])
694  return true;
695  break;
696  //BarrelECal (TEcal)
697  case SubDetId::kLeftTECAL:
698  if (pos[0] > DetDef::tecalLmin[0]+FVdefmin[0] &&
699  pos[0] < DetDef::tecalLmax[0]-FVdefmax[0] &&
700  pos[1] > DetDef::tecalLmin[1]+FVdefmin[1] &&
701  pos[1] < DetDef::tecalLmax[1]-FVdefmax[1] &&
702  pos[2] > DetDef::tecalLmin[2]+FVdefmin[2] &&
703  pos[2] < DetDef::tecalLmax[2]-FVdefmax[2])
704  return true;
705  break;
706  case SubDetId::kRightTECAL:
707  if (pos[0] > DetDef::tecalRmin[0]+FVdefmin[0] &&
708  pos[0] < DetDef::tecalRmax[0]-FVdefmax[0] &&
709  pos[1] > DetDef::tecalRmin[1]+FVdefmin[1] &&
710  pos[1] < DetDef::tecalRmax[1]-FVdefmax[1] &&
711  pos[2] > DetDef::tecalRmin[2]+FVdefmin[2] &&
712  pos[2] < DetDef::tecalRmax[2]-FVdefmax[2])
713  return true;
714  break;
715  case SubDetId::kTopTECAL:
716  if (pos[0] > DetDef::tecalTLmin[0]+FVdefmin[0] &&
717  pos[0] < DetDef::tecalTLmax[0]-FVdefmax[0] &&
718  pos[1] > DetDef::tecalTLmin[1]+FVdefmin[1] &&
719  pos[1] < DetDef::tecalTLmax[1]-FVdefmax[1] &&
720  pos[2] > DetDef::tecalTLmin[2]+FVdefmin[2] &&
721  pos[2] < DetDef::tecalTLmax[2]-FVdefmax[2])
722  return true;
723  //two parts symmetric w.r.t to Z axis
724  if (pos[0] > DetDef::tecalTRmin[0]+FVdefmax[0] &&
725  pos[0] < DetDef::tecalTRmax[0]-FVdefmin[0] &&
726  pos[1] > DetDef::tecalTRmin[1]+FVdefmin[1] &&
727  pos[1] < DetDef::tecalTRmax[1]-FVdefmax[1] &&
728  pos[2] > DetDef::tecalTRmin[2]+FVdefmin[2] &&
729  pos[2] < DetDef::tecalTRmax[2]-FVdefmax[2])
730  return true;
731  break;
732  case SubDetId::kBottomTECAL:
733  if (pos[0] > DetDef::tecalBLmin[0]+FVdefmin[0] &&
734  pos[0] < DetDef::tecalBLmax[0]-FVdefmax[0] &&
735  pos[1] > DetDef::tecalBLmin[1]+FVdefmin[1] &&
736  pos[1] < DetDef::tecalBLmax[1]-FVdefmax[1] &&
737  pos[2] > DetDef::tecalBLmin[2]+FVdefmin[2] &&
738  pos[2] < DetDef::tecalBLmax[2]-FVdefmax[2])
739  return true;
740  //two parts symmetric w.r.t to Z axis
741  if (pos[0] > DetDef::tecalBRmin[0]+FVdefmax[0] &&
742  pos[0] < DetDef::tecalBRmax[0]-FVdefmin[0] &&
743  pos[1] > DetDef::tecalBRmin[1]+FVdefmin[1] &&
744  pos[1] < DetDef::tecalBRmax[1]-FVdefmax[1] &&
745  pos[2] > DetDef::tecalBRmin[2]+FVdefmin[2] &&
746  pos[2] < DetDef::tecalBRmax[2]-FVdefmax[2])
747  return true;
748  break;
749  //P0DECal (PEcal)
750  case SubDetId::kLeftPECAL:
751  if (pos[0] > DetDef::pecalLmin[0]+FVdefmin[0] &&
752  pos[0] < DetDef::pecalLmax[0]-FVdefmax[0] &&
753  pos[1] > DetDef::pecalLmin[1]+FVdefmin[1] &&
754  pos[1] < DetDef::pecalLmax[1]-FVdefmax[1] &&
755  pos[2] > DetDef::pecalLmin[2]+FVdefmin[2] &&
756  pos[2] < DetDef::pecalLmax[2]-FVdefmax[2])
757  return true;
758  break;
759  case SubDetId::kRightPECAL:
760  if (pos[0] > DetDef::pecalRmin[0]+FVdefmin[0] &&
761  pos[0] < DetDef::pecalRmax[0]-FVdefmax[0] &&
762  pos[1] > DetDef::pecalRmin[1]+FVdefmin[1] &&
763  pos[1] < DetDef::pecalRmax[1]-FVdefmax[1] &&
764  pos[2] > DetDef::pecalRmin[2]+FVdefmin[2] &&
765  pos[2] < DetDef::pecalRmax[2]-FVdefmax[2])
766  return true;
767  break;
768  case SubDetId::kTopPECAL:
769  if (pos[0] > DetDef::pecalTLmin[0]+FVdefmin[0] &&
770  pos[0] < DetDef::pecalTLmax[0]-FVdefmax[0] &&
771  pos[1] > DetDef::pecalTLmin[1]+FVdefmin[1] &&
772  pos[1] < DetDef::pecalTLmax[1]-FVdefmax[1] &&
773  pos[2] > DetDef::pecalTLmin[2]+FVdefmin[2] &&
774  pos[2] < DetDef::pecalTLmax[2]-FVdefmax[2])
775  return true;
776  //two parts symmetric w.r.t to Z axis
777  if (pos[0] > DetDef::pecalTRmin[0]+FVdefmax[0] &&
778  pos[0] < DetDef::pecalTRmax[0]-FVdefmin[0] &&
779  pos[1] > DetDef::pecalTRmin[1]+FVdefmin[1] &&
780  pos[1] < DetDef::pecalTRmax[1]-FVdefmax[1] &&
781  pos[2] > DetDef::pecalTRmin[2]+FVdefmin[2] &&
782  pos[2] < DetDef::pecalTRmax[2]-FVdefmax[2])
783  return true;
784  break;
785 
786  case SubDetId::kBottomPECAL:
787  if (pos[0] > DetDef::pecalBLmin[0]+FVdefmin[0] &&
788  pos[0] < DetDef::pecalBLmax[0]-FVdefmax[0] &&
789  pos[1] > DetDef::pecalBLmin[1]+FVdefmin[1] &&
790  pos[1] < DetDef::pecalBLmax[1]-FVdefmax[1] &&
791  pos[2] > DetDef::pecalBLmin[2]+FVdefmin[2] &&
792  pos[2] < DetDef::pecalBLmax[2]-FVdefmax[2])
793  return true;
794  //two parts symmetric w.r.t to Z axis
795  if (pos[0] > DetDef::pecalBRmin[0]+FVdefmax[0] &&
796  pos[0] < DetDef::pecalBRmax[0]-FVdefmin[0] &&
797  pos[1] > DetDef::pecalBRmin[1]+FVdefmin[1] &&
798  pos[1] < DetDef::pecalBRmax[1]-FVdefmax[1] &&
799  pos[2] > DetDef::pecalBRmin[2]+FVdefmin[2] &&
800  pos[2] < DetDef::pecalBRmax[2]-FVdefmax[2])
801  return true;
802  break;
803  //SMRD
804  case SubDetId::kLeftSMRD:
805  if (pos[1] > DetDef::smrd15Lmin[1]+FVdefmin[1] &&
806  pos[1] < DetDef::smrd15Lmax[1]-FVdefmax[1] &&
807  pos[2] > DetDef::smrd15Lmin[2]+FVdefmin[2] &&
808  pos[2] < DetDef::smrd78Lmax[2]-FVdefmax[2]){
809 
810  if (pos[2] > DetDef::smrd15Lmin[2] &&
811  pos[2] < DetDef::smrd6Lmin[2] &&
812  pos[0] > DetDef::smrd15Lmin[0]+FVdefmin[0] &&
813  pos[0] < DetDef::smrd15Lmax[0]-FVdefmax[0])
814  return true;
815 
816  if (pos[2] > DetDef::smrd6Lmin[2] &&
817  pos[2] < DetDef::smrd78Lmin[2] &&
818  pos[0] > DetDef::smrd6Lmin[0]+FVdefmin[0] &&
819  pos[0] < DetDef::smrd6Lmax[0]-FVdefmax[0])
820  return true;
821 
822  if (pos[2] > DetDef::smrd78Lmin[2] &&
823  pos[2] < DetDef::smrd78Lmax[2] &&
824  pos[0] > DetDef::smrd78Lmin[0]+FVdefmin[0] &&
825  pos[0] < DetDef::smrd78Lmax[0]-FVdefmax[0])
826  return true;
827  }
828  break;
829 
830  case SubDetId::kRightSMRD:
831  if (pos[1] > DetDef::smrd15Rmin[1]+FVdefmin[1] &&
832  pos[1] < DetDef::smrd15Rmax[1]-FVdefmax[1] &&
833  pos[2] > DetDef::smrd15Rmin[2]+FVdefmin[2] &&
834  pos[2] < DetDef::smrd78Rmax[2]-FVdefmax[2]){
835 
836  if (pos[2] > DetDef::smrd15Rmin[2] &&
837  pos[2] < DetDef::smrd6Rmin[2] &&
838  pos[0] > DetDef::smrd15Rmin[0]+FVdefmin[0] &&
839  pos[0] < DetDef::smrd15Rmax[0]-FVdefmax[0])
840  return true;
841 
842  if (pos[2] > DetDef::smrd6Rmin[2] &&
843  pos[2] < DetDef::smrd78Rmin[2] &&
844  pos[0] > DetDef::smrd6Rmin[0]+FVdefmin[0] &&
845  pos[0] < DetDef::smrd6Rmax[0]-FVdefmax[0])
846  return true;
847 
848  if (pos[2] > DetDef::smrd78Rmin[2] &&
849  pos[2] < DetDef::smrd78Rmax[2] &&
850  pos[0] > DetDef::smrd78Rmin[0]+FVdefmin[0] &&
851  pos[0] < DetDef::smrd78Rmax[0]-FVdefmax[0])
852  return true;
853  }
854  break;
855  case SubDetId::kTopSMRD:
856  if (pos[0] > DetDef::smrdTLmin[0]+FVdefmin[0] &&
857  pos[0] < DetDef::smrdTLmax[0]-FVdefmax[0] &&
858  pos[1] > DetDef::smrdTLmin[1]+FVdefmin[1] &&
859  pos[1] < DetDef::smrdTLmax[1]-FVdefmax[1] &&
860  pos[2] > DetDef::smrdTLmin[2]+FVdefmin[2] &&
861  pos[2] < DetDef::smrdTLmax[2]-FVdefmax[2])
862  return true;
863  //two parts symmetric w.r.t to Z axis
864  if (pos[0] > DetDef::smrdTRmin[0]+FVdefmax[0] &&
865  pos[0] < DetDef::smrdTRmax[0]-FVdefmin[0] &&
866  pos[1] > DetDef::smrdTRmin[1]+FVdefmin[1] &&
867  pos[1] < DetDef::smrdTRmax[1]-FVdefmax[1] &&
868  pos[2] > DetDef::smrdTRmin[2]+FVdefmin[2] &&
869  pos[2] < DetDef::smrdTRmax[2]-FVdefmax[2])
870  return true;
871  break;
872  case SubDetId::kBottomSMRD:
873  if (pos[0] > DetDef::smrdBLmin[0]+FVdefmin[0] &&
874  pos[0] < DetDef::smrdBLmax[0]-FVdefmax[0] &&
875  pos[1] > DetDef::smrdBLmin[1]+FVdefmin[1] &&
876  pos[1] < DetDef::smrdBLmax[1]-FVdefmax[1] &&
877  pos[2] > DetDef::smrdBLmin[2]+FVdefmin[2] &&
878  pos[2] < DetDef::smrdBLmax[2]-FVdefmax[2])
879  return true;
880  //two parts symmetric w.r.t to Z axis
881  if (pos[0] > DetDef::smrdBRmin[0]+FVdefmax[0] &&
882  pos[0] < DetDef::smrdBRmax[0]-FVdefmin[0] &&
883  pos[1] > DetDef::smrdBRmin[1]+FVdefmin[1] &&
884  pos[1] < DetDef::smrdBRmax[1]-FVdefmax[1] &&
885  pos[2] > DetDef::smrdBRmin[2]+FVdefmin[2] &&
886  pos[2] < DetDef::smrdBRmax[2]-FVdefmax[2])
887  return true;
888  break;
889  default:
890  std::cout << "Error: anaUtils::InFiducialVolume() given an unknown subdetector enumeration: " << det << std::endl;
891 
892  }
893 
894  return false;
895 
896 }
897 
898 
899 //**************************************************
900 int anaUtils::GetAllChargedTrajInTPCInBunch(const AnaEventB& event, AnaTrueParticleB* chargedtrajInBunch[]) {
901 //**************************************************
902  int count = 0;
903  for(Int_t i=0;i< event.nTrueParticles;i++){
904  if(!event.TrueParticles[i]->TrueVertex) continue;
905  if(event.TrueParticles[i]->TrueVertex->Bunch!=event.Bunch) continue;
906  if(event.TrueParticles[i]->Charge==0)continue;
907  Float_t dist=-9999999;
908  for(Int_t idet=0;idet<event.TrueParticles[i]->nDetCrossings;idet++){
909  //i.e crossing the active part of the tpc
910  if(SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet]->Detector, SubDetId::kTPC) && event.TrueParticles[i]->DetCrossings[idet]->InActive) {
912  if(sep>dist) dist=sep;
913  }
914  }
915  // 30* 30 originally
916  if((dist)>900 && event.TrueParticles[i]->Momentum>5){//bigger than 3 TPC hits (30*30 is faster that sqrt(dist)), and momentum > 5 MeV
917  chargedtrajInBunch[count] = event.TrueParticles[i];
918  ++count;
919  }
920  }
921  return count;
922 }
923 
924 //**************************************************
926 //**************************************************
927  /*
928  * we need here to select in-FGD tracks that potentially should have been reconstruced
929  * by FGD iso recon (the function name is confusing);
930  * this involves putting some min requirements for the true tracks:
931  * since FGD iso recon requires a track to extend for at least 4 Z layers (i.e. having hits in five consequitive layers)
932  * in order to be reconstruced this requirement should be applied for the true tracks as well.
933  * In principle one can use the geometry info to retrieve layers that true entrance and exit point correspond to
934  * but it can be time taking, so we use an approximation: a true trajectory should have a length in Z at least of the one of 4 FGD layers:
935  * so 4 cm
936 
937  */
938 
939  int count = 0;
940  for (Int_t i = 0; i < event.nTrueParticles; i++) {
941  if(!event.TrueParticles[i]->TrueVertex) continue;
942  if(event.TrueParticles[i]->TrueVertex->Bunch!=event.Bunch) continue;
943  if(event.TrueParticles[i]->Charge==0)continue;
944  Float_t dist = -9999999;
945  for (Int_t idet = 0; idet < event.TrueParticles[i]->nDetCrossings; idet++) {
946  // i.e crossing the active part of the FGD
947  if (SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet]->Detector, SubDetId::kFGD)){
948  if (SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet]->Detector, det) && event.TrueParticles[i]->DetCrossings[idet]->InActive) {
949  //the separation should be done using the z position, since the fgd is separated by layer in z,
950  //making the z position the reconstructed quantity to be cut on
951  Float_t sep = fabs(event.TrueParticles[i]->DetCrossings[idet]->EntrancePosition[2] - event.TrueParticles[i]->DetCrossings[idet]->ExitPosition[2]);
952  if(sep>dist) dist=sep;
953  }
954  }
955  }
956  // apply the cut (this cut is only valid for FGD!)
957  if (dist > 40){
958  chargedtrajInBunch[count] = event.TrueParticles[i];
959  ++count;
960  }
961  }
962 
963  return count;
964 }
965 
966 //**************************************************
968 //**************************************************
969  AnaTrueParticleB* trajInBunchInFgdx[NMAXTRUEPARTICLES];
970  Int_t ntrajInBunchInFgdx = anaUtils::GetAllChargedTrajInFGDInBunch(event, trajInBunchInFgdx,det);
971 
972  Int_t count = 0;
973  for (Int_t i = 0; i < ntrajInBunchInFgdx; i++) {
974  Float_t dist=-999999.;
975  for(Int_t idet=0;idet<trajInBunchInFgdx[i]->nDetCrossings;idet++){
976  //i.e crossing the active part of the tpc
977  if(SubDetId::GetDetectorUsed(trajInBunchInFgdx[i]->DetCrossings[idet]->Detector, SubDetId::kTPC) && trajInBunchInFgdx[i]->DetCrossings[idet]->InActive) {
978  Float_t sep = GetSeparationSquared(trajInBunchInFgdx[i]->DetCrossings[idet]->EntrancePosition, trajInBunchInFgdx[i]->DetCrossings[idet]->ExitPosition);
979 
980  if(sep>dist) dist=sep;
981  }
982  }
983 
984  bool cross_tpc = false;
985  // 250*250 originally
986  if(dist>62500)//bigger than the ~1/4 of the width of the TPC
987  cross_tpc = true;
988 
989  if (!cross_tpc){
990  chargedtrajInBunch[count] = trajInBunchInFgdx[i];
991  ++count;
992  }
993  }
994 
995  return count;
996 }
997 
998 //**************************************************
1000 //**************************************************
1001 
1002  int count = 0;
1003 
1004  for(Int_t i=0;i<event.nTrueParticles;i++){
1005  if(!event.TrueParticles[i]->TrueVertex) continue;
1006  if(event.TrueParticles[i]->TrueVertex->Bunch!=event.Bunch) continue;
1007  if(event.TrueParticles[i]->Charge==0)continue;
1008 
1009  Float_t dist=-9999999;
1010  for(Int_t idet=0;idet<event.TrueParticles[i]->nDetCrossings;idet++){
1011  //i.e crossing the active part of one of the FGDs
1012  if(SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet]->Detector, SubDetId::kFGD)){
1013  for(Int_t idet1=0;idet1<event.TrueParticles[i]->nDetCrossings;idet1++){
1014  //look for TPC1-FGD1, FGD1-TPC2, TPC2-FGD2, FGD2-TPC3 trajectories
1015  if((SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet]->Detector, SubDetId::kFGD1) &&
1016  (SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet1]->Detector, SubDetId::kTPC1) || SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet1]->Detector, SubDetId::kTPC2))) ||
1017  (SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet]->Detector, SubDetId::kFGD2) &&
1018  (SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet1]->Detector, SubDetId::kTPC2) || SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet1]->Detector, SubDetId::kTPC3))))
1019  {
1020  Float_t sep = GetSeparationSquared(event.TrueParticles[i]->DetCrossings[idet1]->EntrancePosition, event.TrueParticles[i]->DetCrossings[idet1]->ExitPosition);
1021  if(sep>dist) dist=sep;
1022  }
1023  }
1024  }
1025  }
1026 
1027  // 10*10 originally, now 100
1028  if(dist>100){
1029  chargedtrajInBunch[count] = event.TrueParticles[i];
1030  ++count;
1031  }
1032  }
1033 
1034  return count;
1035 
1036 }
1037 
1038 //**************************************************
1040 //**************************************************
1041 
1042  int count = 0;
1043 
1044  for(Int_t i=0;i< event.nTrueParticles;i++){
1045  if(!event.TrueParticles[i]->TrueVertex) continue;
1046  if(event.TrueParticles[i]->TrueVertex->Bunch!=event.Bunch) continue;
1047  if(event.TrueParticles[i]->Charge==0)continue;
1048 
1049  Float_t dist=0;
1050  for(Int_t idet=0;idet<event.TrueParticles[i]->nDetCrossings;idet++){
1051  //i.e crossing the active part of the tpc
1052  if(SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet]->Detector, SubDetId::kTPC) && event.TrueParticles[i]->DetCrossings[idet]->InActive) {
1054  if(sep>dist) dist=sep;
1055  }
1056  }
1057  //250*250 originally
1058  if(dist>62500){//bigger than the ~1/4 of the width of the TPC
1059  chargedtrajInBunch[count] = event.TrueParticles[i];
1060  ++count;
1061  }
1062 
1063  }
1064  return count;
1065 }
1066 
1067 int anaUtils::GetAllChargedTrajInP0DInBunch(const AnaEventB& event, AnaTrueParticleB* chargedTrajInBunch[])
1068 {
1069  int count = 0;
1070  for (Int_t i = 0; i<event.nTrueParticles; i++)
1071  {
1072  if (!event.TrueParticles[i]->TrueVertex) continue;
1073  if (event.TrueParticles[i]->TrueVertex->Bunch!=event.Bunch) continue;
1074  if (event.TrueParticles[i]->Charge==0) continue;
1075  Float_t dist = -9999999;
1076  for (Int_t idet = 0; idet<event.TrueParticles[i]->nDetCrossings;idet++)
1077  {
1078  if (SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet]->Detector,SubDetId::kP0D) && event.TrueParticles[i]->DetCrossings[idet]->InActive){
1080  if (sep>dist) dist = sep;
1081  }
1082  }
1083  if ( (dist)>900 && event.TrueParticles[i]->Momentum>5){
1084  chargedTrajInBunch[count] = event.TrueParticles[i];
1085  count++;
1086  }
1087  if (count >= (int) NMAXTRUEPARTICLES) break;
1088  }
1089 
1090  return count;
1091 
1092 }
1093 
1094 int anaUtils::GetAllChargedTrajInP0DAndTPCInBunch(const AnaEventB& event, AnaTrueParticleB* chargedTrajInBunch[])
1095 {
1096  int count = 0;
1097  for (Int_t i = 0; i<event.nTrueParticles; i++)
1098  {
1099  if (!event.TrueParticles[i]->TrueVertex) continue;
1100  if (event.TrueParticles[i]->TrueVertex->Bunch!=event.Bunch) continue;
1101  if (event.TrueParticles[i]->Charge==0) continue;
1102  Float_t dist = -9999999;
1103  for (Int_t idet = 0; idet<event.TrueParticles[i]->nDetCrossings;idet++)
1104  {
1105  if (SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet]->Detector,SubDetId::kP0D) && event.TrueParticles[i]->DetCrossings[idet]->InActive){
1106  for (Int_t idet1=0; idet1<event.TrueParticles[i]->nDetCrossings; idet1++)
1107  {
1108  if (SubDetId::GetDetectorUsed(event.TrueParticles[i]->DetCrossings[idet1]->Detector,SubDetId::kTPC1) && event.TrueParticles[i]->DetCrossings[idet1]->InActive){
1110  if (sep>dist) dist = sep;
1111  }//in TPC1
1112  }//Loop over detector crossings
1113  }//in P0D
1114  }//Loop over detector crossings
1115  if ( (dist)>100 && event.TrueParticles[i]->Momentum>5){
1116  chargedTrajInBunch[count] = event.TrueParticles[i];
1117  count++;
1118  }
1119  if (count >= (int) NMAXTRUEPARTICLES) break;
1120  }
1121 
1122  return count;
1123 
1124 
1125 }
1126 
1127 int anaUtils::GetAllChargedTrajInP0DAndNoTPCInBunch(const AnaEventB& event, AnaTrueParticleB* chargedTrajInBunch[])
1128 {
1129  AnaTrueParticleB* trajInP0D[NMAXTRUEPARTICLES];
1130  Int_t ntrajInP0D = GetAllChargedTrajInP0DInBunch(event,trajInP0D);
1131  Int_t count = 0;
1132  for (Int_t i = 0; i < ntrajInP0D; i++){
1133  Float_t dist=-999999.;
1134  for (Int_t idet=0; idet<trajInP0D[i]->nDetCrossings; idet++)
1135  {
1136  if (SubDetId::GetDetectorUsed(trajInP0D[i]->DetCrossings[idet]->Detector,SubDetId::kTPC1) && trajInP0D[i]->DetCrossings[idet]->InActive){
1137  Float_t sep = anaUtils::GetSeparationSquared(trajInP0D[i]->DetCrossings[idet]->EntrancePosition,trajInP0D[i]->DetCrossings[idet]->ExitPosition);
1138  if (sep>dist) dist=sep;
1139 
1140 
1141  }
1142  }
1143  bool cross_tpc = false;
1144  if (dist>900)//originally 62500=250**2: 1/4 width of TPC
1145  cross_tpc = true;
1146  if (!cross_tpc){
1147  chargedTrajInBunch[count] = trajInP0D[i];
1148  count++;
1149  if (count>=(int)NMAXTRUEPARTICLES) break;
1150  }
1151 
1152  }
1153  return count;
1154 
1155 
1156 }
1157 
1158 
1159 //**************************************************
1161 //**************************************************
1162 
1163  int count = 0;
1164  for (int it = 0; it < event.nParticles; it++) {
1165  AnaTrackB* track = static_cast<AnaTrackB*>(event.Particles[it]);
1166  if (!SubDetId::GetDetectorUsed(track->Detector, SubDetId::kTPC) && SubDetId::GetDetectorUsed(track->Detector, fgddet)) {
1167  selTracks[count] = track;
1168  ++count;
1169  }
1170  }
1171 
1172  // Sort by decreasing number of hits
1173  std::sort(&selTracks[0], &selTracks[count], AnaParticleB::CompareNHits);
1174 
1175  return count;
1176 }
1177 
1178 
1179 //**************************************************
1181 //**************************************************
1182 
1183  int count = 0;
1184  for (int it = 0; it < event.nParticles; it++) {
1185  AnaTrackB* track = static_cast<AnaTrackB*>(event.Particles[it]);
1186  if (SubDetId::GetDetectorUsed(track->Detector, SubDetId::kTPC) || SubDetId::GetDetectorUsed(track->Detector, SubDetId::kFGD)) {
1187  selTracks[count] = track;
1188  ++count;
1189  }
1190  }
1191 
1192  // Sort by decreasing number of hits
1193  std::sort(&selTracks[0], &selTracks[count], AnaParticleB::CompareNHits);
1194 
1195  return count;
1196 }
1197 
1198 //**************************************************
1199 bool anaUtils::HasTrackUsingTPC(const AnaEventB& event) {
1200 //**************************************************
1201  return anaUtils::HasTrackUsingDet(event, SubDetId::kTPC);
1202 }
1203 
1204 //**************************************************
1205 int anaUtils::GetAllTracksUsingP0D(const AnaEventB& event, AnaTrackB* selTracks[]) {
1206 //**************************************************
1207  return GetAllTracksUsingDet(event, SubDetId::kP0D, selTracks);
1208 }
1209 
1210 //**************************************************
1211 int anaUtils::GetAllTracksUsingFGD(const AnaEventB& event, AnaTrackB* selTracks[]) {
1212 //**************************************************
1213  return GetAllTracksUsingDet(event, SubDetId::kFGD, selTracks);
1214 }
1215 
1216 //**************************************************
1217 int anaUtils::GetAllTracksUsingTPC(const AnaEventB& event, AnaTrackB* selTracks[]) {
1218 //**************************************************
1219  return GetAllTracksUsingDet(event, SubDetId::kTPC, selTracks);
1220 }
1221 
1222 //**************************************************
1223 int anaUtils::GetAllTracksUsingECAL(const AnaEventB& event, AnaTrackB* selTracks[]) {
1224 //**************************************************
1225  return GetAllTracksUsingDet(event, SubDetId::kECAL, selTracks);
1226 }
1227 
1228 //**************************************************
1229 int anaUtils::GetAllTracksUsingSMRD(const AnaEventB& event, AnaTrackB* selTracks[]) {
1230 //**************************************************
1231  return GetAllTracksUsingDet(event, SubDetId::kSMRD, selTracks);
1232 }
1233 
1234 //**************************************************
1236 //**************************************************
1237 
1238  int count = 0;
1239 
1240  SubDetId::SubDetEnum dets[2];
1241  dets[0] = SubDetId::kTPC;
1242  dets[1] = det;
1243 
1244  // Loop over all tracks
1245  for (int it = 0; it < event.nParticles; it++) {
1246  AnaTrackB* track = static_cast<AnaTrackB*>(event.Particles[it]);
1247  if (anaUtils::TrackUsesDets(*track, dets, 2)){
1248  count ++;
1249  }
1250  }
1251 
1252  return count;
1253 }
1254 
1255 //********************************************************************
1257 //********************************************************************
1258  return (SubDetId::GetDetectorUsed(track->Detector, SubDetId::kTECAL) || SubDetId::GetDetectorUsed(track->Detector, SubDetId::kDSECAL));
1259 }
1260 
1261 //********************************************************************
1263 //********************************************************************
1264 
1265  int count = 0;
1266  for (Int_t i = 0; i < track->nECALSegments; i++) {
1267  AnaECALParticleB* ecal = track->ECALSegments[i];
1268  if (SubDetId::GetDetectorUsed(ecal->Detector, SubDetId::kTECAL) || SubDetId::GetDetectorUsed(ecal->Detector, SubDetId::kDSECAL)) {
1269  ecals[count] = ecal;
1270  count ++;
1271  }
1272  }
1273 
1274  return count;
1275 }
1276 
1277 //********************************************************************
1279 //********************************************************************
1280  if (track->nTPCSegments == 0) {
1281  return NULL;
1282  }
1283 
1284  return track->TPCSegments[track->nTPCSegments-1];
1285 }
1286 
1287 //********************************************************************
1289 //********************************************************************
1290  AnaECALParticleB* upstream = NULL;
1291 
1292  Float_t pos_z = 1e4;
1293 
1294  if (track) {
1295  for (Int_t i = 0; i < track->nECALSegments; i++) {
1296  AnaECALParticleB* ecal = track->ECALSegments[i];
1297  if(!ecal) continue;
1298  if (ecal->PositionStart[2] < pos_z) {
1299  upstream = ecal;
1300  pos_z = ecal->PositionStart[2];
1301  }
1302  if (ecal->PositionEnd[2] < pos_z) {
1303  upstream = ecal;
1304  pos_z = ecal->PositionEnd[2];
1305  }
1306  }
1307  }
1308 
1309  return upstream;
1310 }
1311 
1312 //********************************************************************
1313 //bool anaUtils::IsEcalContained(AnaEcalTrackEcalPID* ecal) {
1314 //********************************************************************
1315 /*
1316  TVector3 pos = ecal->ShowerPosition;
1317 
1318  if (ecal->Detector == subDetId::kDSECAL) {
1319  return (fabs(pos.X()) < 900 && fabs(pos.Y()) < 900);
1320  } else if (ecal->Detector == subDetId::kLeftTECAL || ecal->Detector == subDetId::kRightTECAL) {
1321  return (fabs(pos.Z()) < 2900 && fabs(pos.Y()) < 900);
1322  } else {
1323  return (fabs(pos.Z()) < 2900 && fabs(pos.X()) < 900);
1324  }
1325  }
1326 */
1327 
1328 
1329 //********************************************************************
1330 Int_t anaUtils::GetFgdModuleType(bool IsMC, const Float_t* pos, SubDetId::SubDetEnum det, bool includeGlueSkin){
1331 //********************************************************************
1332 
1333  if (det == SubDetId::kFGD1) std::cout << "FGD1 not supported" << std::endl;
1334  if ( ! anaUtils::InFiducialVolume(det,pos)) return 7; // here CATOUTFV is not known
1335  else return GetFgdModuleTypeNoFV(IsMC, pos, includeGlueSkin);
1336 }
1337 
1338 
1339 //********************************************************************
1340 Int_t anaUtils::GetFgdModuleTypeNoFV(bool IsMC, const Float_t* pos, bool includeGlueSkin){
1341 //********************************************************************
1342 
1343  Int_t material = GetFgdMaterialNoFV(IsMC, pos);
1344 
1345  if (includeGlueSkin) {
1346  if (material >= 10 && material <= 14) return 0; // 1st half XY module
1347  if (material >= 15 && material <= 19) return 1; // 2nd half XY module
1348  if (material >= 20 && material <= 26) return 2; // water module
1349  }
1350  else {
1351  // OLD WAY, until nd280Highland2 v1r3
1352  if (material == 10) return 0; // X layer
1353  if (material == 15) return 1; // Y layer
1354  if (material == 20) return 2; // water/PC(polycarbonate) panel
1355  }
1356 
1357  return 3; // gaps, should be empty if includeGlueSkin=true (default)
1358 }
1359 
1360 
1361 //********************************************************************
1362 Int_t anaUtils::GetFgdMaterial(bool IsMC, const Float_t* pos, SubDetId::SubDetEnum det){
1363 //********************************************************************
1364 
1365  if (det == SubDetId::kFGD1) std::cout << "FGD1 not supported" << std::endl;
1366  if ( ! anaUtils::InFiducialVolume(det,pos)) return 7; // here CATOUTFV is not known
1367  else return GetFgdMaterialNoFV(IsMC, pos);
1368 }
1369 
1370 
1371 //********************************************************************
1372 Int_t anaUtils::GetFgdMaterialNoFV(bool IsMC, const Float_t* pos){
1373 //********************************************************************
1374 
1375  // to do: do it also for FGD1
1376  Float_t barCenter[14];
1377  Float_t waterCenter[6];
1378  if (IsMC) {
1379  for (Int_t i=0; i<14; i++) barCenter[i] = DetDef::fgd2BarCenterMC[i];
1380  for (Int_t i=0; i<6; i++) waterCenter[i] = DetDef::fgd2WaterCenterMC[i];
1381  } else {
1382  for (Int_t i=0; i<14; i++) barCenter[i] = DetDef::fgd2BarCenterDATA[i];
1383  for (Int_t i=0; i<6; i++) waterCenter[i] = DetDef::fgd2WaterCenterDATA[i];
1384  }
1385 
1386  // this has to be done like this because Float_t appoximation leaves gaps
1387  Float_t midpos_preX[7];
1388  Float_t midpos_postY[7];
1389  Float_t midpos_XY[7];
1390  Int_t ixymod=0;
1391 
1392  // look in scintillator layers
1393  Float_t initZ, finalZ;
1394  for (Int_t ibar=0; ibar<14; ibar++) {
1395 
1396  initZ = barCenter[ibar] - (DetDef::fgdXYBarWidth/2);
1397  finalZ = barCenter[ibar] + (DetDef::fgdXYBarWidth/2);
1398 
1399  if (ibar%2 == 0) {
1400  if (pos[2] >= initZ && pos[2] <= finalZ) return 10; // X layer
1401  initZ = initZ - DetDef::fgdXYGlueWidth;
1402  if (pos[2] >= initZ && pos[2] <= finalZ) return 11; // glue before X layer
1403  initZ = initZ - DetDef::fgdXYSkinWidth;
1404  if (pos[2] >= initZ && pos[2] <= finalZ) return 12; // G10 skin before X layer
1405  initZ = initZ - (DetDef::fgdWaterAirWidth/2);
1406  if (pos[2] >= initZ && pos[2] <= finalZ) return 13; // half of air before X layer
1407  finalZ = finalZ + (DetDef::fgdXYMiddleGlueWidth/2);
1408  if (pos[2] >= initZ && pos[2] <= finalZ) return 14; // half of glue between X and Y
1409 
1410  midpos_preX[ixymod] = initZ;
1411  midpos_XY[ixymod] = finalZ;
1412 
1413  } else if (ibar%2 == 1) {
1414  if (pos[2] >= initZ && pos[2] <= finalZ) return 15; // Y layer
1415  finalZ = finalZ + DetDef::fgdXYGlueWidth;
1416  if (pos[2] >= initZ && pos[2] <= finalZ) return 16; // glue after Y layer
1417  finalZ = finalZ + DetDef::fgdXYSkinWidth;
1418  if (pos[2] >= initZ && pos[2] <= finalZ) return 17; // G10 skin after Y layer
1419  finalZ = finalZ + (DetDef::fgdXYAirWidth/2);
1420  if (pos[2] >= initZ && pos[2] <= finalZ) return 18; // half of air after Y layer
1421 // initZ = initZ - (DetDef::fgdXYMiddleGlueWidth/2);
1422 // if (pos[2] >= initZ && pos[2] <= finalZ) return 19; // half of glue between X and Y
1423 
1424  if (pos[2] >= midpos_XY[ixymod] && pos[2] <= finalZ) return 19; // half of glue between X and Y
1425 
1426  midpos_postY[ixymod] = finalZ;
1427  ixymod++;
1428  }
1429  }
1430 
1431  // look in water modules
1432  for (Int_t ibar=0; ibar<6; ibar++) {
1433  initZ = waterCenter[ibar] - (DetDef::fgdWaterPCWidth/2);
1434  finalZ = waterCenter[ibar] + (DetDef::fgdWaterPCWidth/2);
1435  if (pos[2] >= initZ && pos[2] <= finalZ) return 20; // water/PC(polycarbonate) panel
1436  initZ = initZ - DetDef::fgdWaterGlueWidth;
1437  if (pos[2] >= initZ && pos[2] <= finalZ) return 21; // glue before water/PC
1438  initZ = initZ - DetDef::fgdWaterSkinWidth;
1439  if (pos[2] >= initZ && pos[2] <= finalZ) return 22; // polypropylene skin before water/PC
1440 // initZ = initZ - (DetDef::fgdXYAirWidth/2);
1441 // if (pos[2] >= initZ && pos[2] <= finalZ) return 23; // half of air before water/PC
1442 
1443  if (pos[2] >= midpos_postY[ibar] && pos[2] <= finalZ) return 23; // half of air before water/PC
1444 
1445  finalZ = finalZ + DetDef::fgdWaterGlueWidth;
1446  if (pos[2] >= initZ && pos[2] <= finalZ) return 24; // glue after water/PC
1447  finalZ = finalZ + DetDef::fgdWaterSkinWidth;
1448  if (pos[2] >= initZ && pos[2] <= finalZ) return 25; // polypropilene skin after water/PC
1449 // finalZ = finalZ + (DetDef::fgdWaterAirWidth/2);
1450 // if (pos[2] >= initZ && pos[2] <= finalZ) return 26; // half of air after water/PC
1451 
1452  if (pos[2] >= initZ && pos[2] <= midpos_preX[ibar+1]) return 26; // half of air after water/PC
1453  }
1454 
1455  std::cout << "Error in GetFgdMaterialNoFV: material not found for this position." << std::endl;
1456  return -1;
1457 }
1458 
1459 
1460 //********************************************************************
1461 anaUtils::massComponentEnum anaUtils::GetMassComponent(bool IsMC, const Float_t* pos){
1462 //********************************************************************
1463 
1464  if (anaUtils::InDetVolume(SubDetId::kFGD1,pos))
1465  return anaUtils::kFGD1;
1466  else if (anaUtils::InDetVolume(SubDetId::kFGD2, pos)) {
1467  if (anaUtils::GetFgdModuleTypeNoFV(IsMC,pos,true) == 2)
1468  return anaUtils::kFGD2watermodules;
1469  else
1470  return anaUtils::kFGD2xymodules;
1471  }
1472  /*
1473  else if (anaUtils::InDetVolume(SubDetId::kP0D, pos)){
1474  if (anaUtils::GetP0DModuleTypeNoFV(event,pos,true) == 2)
1475  return anaUtils::kP0DWater;
1476  else
1477  return anaUtils::kP0DAir;
1478  }
1479  */
1480  else return kInvalid;
1481 }
int GetAllTracksUsingECAL(const AnaEventB &event, AnaTrackB *selTracks[])
bool InActive
If the particle passes through an active part of the subdetector.
AnaTrueVertexB * TrueVertex
Pointer to the AnaTrueVertexB of the interaction that created this AnaTrueParticleB.
int nP0DSegments
How many P0D tracks are associated with this track.
Float_t PositionStart[4]
The reconstructed start position of the particle.
AnaP0DParticleB * P0DSegments[NMAXP0DS]
The P0D segments that contributed to this global track.
unsigned long Detector
int nDetCrossings
The number of DetCrossing objects.
unsigned long Detector
bool UsesTrackerDsEcal(AnaTrackB *track)
Whether this track has an ecal segment in the TrackerECal or DsECal.
Float_t ExitPosition[4]
for each subdetector tell the exit position
static bool IsFGDDetector(SubDetId::SubDetEnum det)
Check if a detector enumeration refers to a FGD or not.
Definition: SubDetId.cxx:126
AnaECALParticleB * ECALSegments[NMAXECALS]
The ECAL segments that contributed to this global track.
int nECALSegments
How many ECAL tracks are associated with this track.
int GetAllTracksUsingTPC(const AnaEventB &event, AnaTrackB *selTracks[])
AnaSMRDParticleB * SMRDSegments[NMAXSMRDS]
The SMRD segments that contributed to this global track.
static bool IsP0DDetector(SubDetId::SubDetEnum det)
Check if a detector enumeration refers to a SMRDP0D or not.
Definition: SubDetId.cxx:146
Int_t NNodes
The number of nodes in the reconstructed object.
SubDetId::SubDetEnum GetClosestTPC(const AnaTrackB &track)
For tracks that start in the FGD, get the closest TPC in the direction of the track.
Representation of an ECAL segment of a global track.
AnaTPCParticleB * TPCSegments[NMAXTPCS]
The TPC segments that contributed to this global track.
int GetAllChargedTrajInFGDAndNoTPCInBunch(const AnaEventB &event, AnaTrueParticleB *chargedtrajInBunch[], SubDetId::SubDetEnum det)
static bool GetDetectorUsed(unsigned long BitField, SubDetId::SubDetEnum det)
Method to see if a certain subdetector or subdetector system is used.
Definition: SubDetId.cxx:40
massComponentEnum GetMassComponent(bool IsMC, const Float_t *pos)
Get the detector component in which the position is.
Float_t Momentum
The initial momentum of the true particle.
SubDetId::SubDetEnum GetDetector(const Float_t *pos)
Return the detector in which the position is.
int GetAllTracksUsingFGD(const AnaEventB &event, AnaTrackB *selTracks[])
AnaParticleB * GetSegmentInDet(const AnaTrackB &track, SubDetId::SubDetEnum det)
float GetSeparationSquared(const Float_t *pos1, const Float_t *pos2)
Calculate the distance between two points.
int GetAllChargedTrajInTPCFGDInBunch(const AnaEventB &event, AnaTrueParticleB *chargedtrajInBunch[])
int GetAllChargedTrajInTPCInBunch(const AnaEventB &event, AnaTrueParticleB *traj[])
static bool CompareNHits(const AnaParticleB *t1, const AnaParticleB *t2)
Compare the NHits of two particles. Return whether t1 is higher NHts than t2.
int GetAllTracksUsingFGDorTPC(const AnaBunchB &bunch, AnaTrackB *selTracks[])
Int_t GetFgdModuleType(bool IsMC, const Float_t *pos, SubDetId::SubDetEnum det, bool includeGlueSkin=true)
AnaECALParticleB * GetMostUpstreamECalSegment(AnaTrackB *track)
int nSMRDSegments
How many SMRD tracks are associated with this track.
int GetAllTracksUsingP0D(const AnaEventB &event, AnaTrackB *selTracks[])
int nTPCSegments
How many TPC tracks are associated with this track.
Representation of a true Monte Carlo trajectory/particle.
Representation of an SMRD segment of a global track.
int GetNTracksUsingTPCAndDet(const AnaEventB &event, SubDetId::SubDetEnum det)
Returns the number of tracks using both the TPC and the subdetector &#39;det&#39;.
int GetSegmentsInDet(const AnaTrackB &track, SubDetId::SubDetEnum det, AnaParticleB *selTracks[])
Float_t EntrancePosition[4]
for each subdetector tell the entrance position
static bool IsTPCDetector(SubDetId::SubDetEnum det)
Check if a detector enumeration refers to a TPC or not.
Definition: SubDetId.cxx:122
AnaDetCrossingB ** DetCrossings
static bool IsECALDetector(SubDetId::SubDetEnum det)
Check if a detector enumeration refers to a ECAL or not.
Definition: SubDetId.cxx:130
SubDetEnum
Enumeration of all detector systems and subdetectors.
Definition: SubDetId.hxx:25
int GetAllTracksUsingFGDAndNoTPC(const AnaBunchB &bunch, AnaTrackB *selTracks[])
Representation of a global track.
Int_t Bunch
The index of this bunch (0-7).
int GetAllBigEnoughChargedTrajInTPCInBunch(const AnaEventB &event, AnaTrueParticleB *chargedtrajInBunch[])
static SubDetId::SubDetEnum GetSubdetectorEnum(unsigned long BitField)
Get the single subdetector that this track is from.
Definition: SubDetId.cxx:165
Representation of a TPC segment of a global track.
bool InDetVolume(SubDetId::SubDetEnum det, const Float_t *pos)
Representation of a FGD segment of a global track.
int nFGDSegments
How many FGD tracks are associated with this track.
AnaTrueParticleB ** TrueParticles
The true MC particles used in this spill.
int GetAllTracksUsingSMRD(const AnaEventB &event, AnaTrackB *selTracks[])
AnaTPCParticleB * GetTPCBackSegment(const AnaTrackB *track)
Get the most dowstream TPC segment of the track.
bool TrackUsesDets(const AnaTrackB &track, SubDetId::SubDetEnum dets[], int nDets)
static int GetTPC(unsigned long BitField)
Definition: SubDetId.cxx:74
AnaParticleB * GetSegmentWithMostNodesInClosestTpc(const AnaTrackB &track)
Combined function to address NuMu selection needs as efficiently as possible - gets the TPC segment w...
Representation of a P0D segment of a global track.
int GetAllTracksUsingDet(const AnaBunchB &bunch, SubDetId::SubDetEnum det, AnaTrackB *selTracks[])
Float_t Charge
The true charge of the particle.
int GetAllChargedTrajInFGDInBunch(const AnaEventB &event, AnaTrueParticleB *traj[], SubDetId::SubDetEnum det)
This namespace contains useful functions for analyses related to kinematics.
bool InFiducialVolume(SubDetId::SubDetEnum det, const Float_t *pos, const Float_t *FVdefmin, const Float_t *FVdefmax)
AnaParticleB * GetSegmentWithMostNodesInDet(const AnaTrackB &track, SubDetId::SubDetEnum det)
Method to get the subtrack with most nodes in a given detector.
Representation of a reconstructed particle (track or shower).
AnaFGDParticleB * FGDSegments[NMAXFGDS]
The FGD segments that contributed to this global track.
int GetTrackerDsEcals(AnaTrackB *track, AnaECALParticleB *selTracks[])
static bool IsSMRDDetector(SubDetId::SubDetEnum det)
Check if a detector enumeration refers to a SMRD or not.
Definition: SubDetId.cxx:142
Float_t PositionEnd[4]
The reconstructed end position of the particle.