FEDRA emulsion software from the OPERA Collaboration
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
EdbEDATrackSelection Class Reference

#include <EdbEDATrackSet.h>

Inheritance diagram for EdbEDATrackSelection:
Collaboration diagram for EdbEDATrackSelection:

Public Member Functions

void DoSelection (TObjArray *tracksbase, TObjArray *tracks)
 
 EdbEDATrackSelection ()
 
bool GetImpactSearch ()
 
int GetNsegCut ()
 
double GetPHCut ()
 
double GetPHDTRMS ()
 
EdbRunTrackingGetRunTracking ()
 
int GetSideOutPlate ()
 
int GetUpstreamPlate ()
 
int ImpactSearch (EdbTrackP *t)
 
int Neighborhood (EdbSegP *s, double *dmin=NULL)
 
void Reset ()
 
void SetAngle (double tx, double ty)
 
void SetAngularCut (bool b=kTRUE)
 
void SetClearPrevious (bool b=kTRUE)
 
void SetCondBTDefault (EdbScanCond &cond)
 
void SetCondMTDefault (EdbScanCond &cond)
 
void SetCondTrackingDefault ()
 
void SetDT (double dtx, double dty=-1)
 
void SetDX (double dx)
 
void SetImpactSearch (bool b=kTRUE, EdbVertex *v=NULL)
 
void SetNeighborSearch (bool b=kTRUE, TObjArray *selected=NULL, double dzup=-1, double dzdown=-1)
 
void SetNsegCut (int nseg)
 
void SetPHCut (double phcut=0.0)
 
void SetPHDTRMS (double slope)
 
void SetRunTracking (EdbRunTracking &rt)
 
void SetSideOut (bool b=kTRUE)
 
void SetSideOutPlate (int npl)
 
void SetUpstreamPlate (int ipl)
 
int SideOut (EdbTrackP *t)
 
virtual ~EdbEDATrackSelection ()
 

Public Attributes

EdbRunTracking eRunTracking
 Run Tracking for microtrack search. More...
 
EdbVertexeVertex
 

Private Attributes

int eAngularCut
 
int eClearPrevious
 
int eImpactSearch
 
double eNeighborDzDown
 
double eNeighborDzUp
 
int eNeighborSearch
 
int eNsegCut
 
double ePHCut
 
double ePHDTRMS
 
int ePlateDown
 
int ePlateUp
 
TObjArray * eSelected
 
int eSideOut
 
int eSideOutPlate
 
double eTolDTX
 
double eTolDTY
 
double eTolDX
 
double eTX
 
double eTY
 

Constructor & Destructor Documentation

◆ EdbEDATrackSelection()

EdbEDATrackSelection::EdbEDATrackSelection ( )
inline
64{ Reset(); }
void Reset()
Definition: EdbEDATrackSet.h:67

◆ ~EdbEDATrackSelection()

virtual EdbEDATrackSelection::~EdbEDATrackSelection ( )
inlinevirtual
65{}

Member Function Documentation

◆ DoSelection()

void EdbEDATrackSelection::DoSelection ( TObjArray *  tracksbase,
TObjArray *  tracks 
)
126 {
127 if(tracksbase==NULL||tracks==NULL) return;
128 if(eClearPrevious) tracks->Clear();
129
132 if(eVertex==NULL) {
134 printf("please select a vertex.\n");
135 }
136 }
137
138 printf("applying cuts, Nseg %d, PlateUp = %d \n",eNsegCut, ePlateUp);
139 if( eSelected==NULL ) eNeighborSearch = 0;
140 if( eNeighborSearch ){
141 printf( "Neighbor search On. tolerance : dTX = %.3lf rad, dTY = %.3lf rad, dmin = %.0lf micron\n",eTolDTX,eTolDTY, eTolDX);
142 }
143 if( eAngularCut ){
144 printf( "Angular cut On. tolerance : dTX = %.3lf rad, dTY = %.3lf rad\n", eTolDTX, eTolDTY);
145 }
146 if( eImpactSearch ){
147 printf( "Impact cut On. tolerance : IP < %.0lf micron\n", eTolDX);
148 }
149
150 if( eSideOut ){
151 printf( "Penetration and Side-Out rejection : require %d scanning area upstream\n", eSideOutPlate);
152 }
153 int nentr = tracksbase->GetEntries();
154 for(int i=0;i<nentr;i++){
155 if(i%1000==0) printf("%d / %d. %d tracks\r", i, tracksbase->GetEntriesFast(), tracks->GetEntriesFast());
156 EdbTrackP *t = (EdbTrackP *) tracksbase->At(i);
157 if(NULL==t) continue;
158 // Selections
159
160 // Nseg cut
161 if( t->N() < eNsegCut) continue;
162
163 // Max plate cut
164 int ipl_first = t->GetSegmentFirst()->Plate();
165
166 if( ipl_first < ePlateUp) continue;
167
168 // Reject Side-out tracks
169 if( eSideOut ){
170 if( SideOut(t) ) continue;
171 }
172
173 // Neighbor Search. use eTolDX, eTolDT
174 if( eNeighborSearch) {
175 if( Neighborhood(t) == 0 ) continue;
176 }
177
178 // Impact Search. use eTolDX
179 if( eImpactSearch ){
180 if( ImpactSearch (t) == 0 ) continue;
181 }
182
183 // Angular Cut
184 if( eAngularCut ){
185 if( fabs(t->TX() - eTX) > eTolDTX) continue;
186 if( fabs(t->TY() - eTY) > eTolDTY) continue;
187 }
188
189 // PH Cut
190 double ph = t->Wgrains()/t->N();
191 if(ePHCut > ph ) continue;
192
193 // PH - Angle deviation cut
194 double dtrms = DTRMS(t);
195 if(dtrms> ePHDTRMS*(ph-ePHCut)) continue;
196
197 if(tracks->FindObject(t)==NULL) tracks->AddLast(t);
198 }
199 printf("%d tracks out of %d tracks\n", tracks->GetEntries(), tracksbase->GetEntries());
200}
EdbEDA * gEDA
Definition: EdbEDA.C:3
Expr< UnaryOp< Fabs< T >, Expr< A, T, D >, T >, T, D > fabs(const Expr< A, T, D > &rhs)
Definition: UnaryOperators.hh:96
EdbVertex * GetSelectedVertex(void)
Definition: EdbEDA.h:439
double eTY
Definition: EdbEDATrackSet.h:50
int SideOut(EdbTrackP *t)
Definition: EdbEDATrackSet.C:48
int eNeighborSearch
Definition: EdbEDATrackSet.h:41
TObjArray * eSelected
Definition: EdbEDATrackSet.h:57
double ePHCut
Definition: EdbEDATrackSet.h:53
int ImpactSearch(EdbTrackP *t)
double eTolDX
Definition: EdbEDATrackSet.h:46
int eImpactSearch
Definition: EdbEDATrackSet.h:42
int ePlateUp
Definition: EdbEDATrackSet.h:39
double eTolDTY
Definition: EdbEDATrackSet.h:48
int eAngularCut
Definition: EdbEDATrackSet.h:43
int eClearPrevious
Definition: EdbEDATrackSet.h:52
int eNsegCut
Definition: EdbEDATrackSet.h:38
int eSideOut
Definition: EdbEDATrackSet.h:44
double ePHDTRMS
Definition: EdbEDATrackSet.h:54
EdbVertex * eVertex
Definition: EdbEDATrackSet.h:61
double eTolDTX
Definition: EdbEDATrackSet.h:47
int eSideOutPlate
Definition: EdbEDATrackSet.h:45
int Neighborhood(EdbSegP *s, double *dmin=NULL)
Definition: EdbEDATrackSet.C:20
double eTX
Definition: EdbEDATrackSet.h:49
Definition: EdbPattern.h:113
TTree * t
Definition: check_shower.C:4
TTree * tracks
Definition: check_tr.C:19
double DTRMS(EdbTrackP *t)
Definition: EdbEDAUtil.C:431
#define NULL
Definition: nidaqmx.h:84

◆ GetImpactSearch()

bool EdbEDATrackSelection::GetImpactSearch ( )
inline
109{ return eImpactSearch;}

◆ GetNsegCut()

int EdbEDATrackSelection::GetNsegCut ( )
inline
114{ return eNsegCut;}

◆ GetPHCut()

double EdbEDATrackSelection::GetPHCut ( )
inline
115{ return ePHCut;}

◆ GetPHDTRMS()

double EdbEDATrackSelection::GetPHDTRMS ( )
inline
116{ return ePHDTRMS;}

◆ GetRunTracking()

EdbRunTracking & EdbEDATrackSelection::GetRunTracking ( )
inline
122{ return eRunTracking;}
EdbRunTracking eRunTracking
Run Tracking for microtrack search.
Definition: EdbEDATrackSet.h:62

◆ GetSideOutPlate()

int EdbEDATrackSelection::GetSideOutPlate ( )
inline
117{ return eSideOutPlate;}

◆ GetUpstreamPlate()

int EdbEDATrackSelection::GetUpstreamPlate ( )
inline
113{ return ePlateUp;}

◆ ImpactSearch()

int EdbEDATrackSelection::ImpactSearch ( EdbTrackP t)

◆ Neighborhood()

int EdbEDATrackSelection::Neighborhood ( EdbSegP s,
double *  dmin = NULL 
)
20 {
21 int flag=0;
22 double dminz;
23 if(dmin) *dmin=1e5;
24 for(int i=0;i<eSelected->GetEntries();i++){
25 EdbSegP *ss = (EdbSegP *)eSelected->At(i);
26 if(ss==s) return 1;
27
28 if( fabs(s->TX()-ss->TX()) > eTolDTX) continue; // angular cut
29 if( fabs(s->TY()-ss->TY()) > eTolDTY) continue;
30
31 double d = CalcDmin( ss, s, &dminz);
32 if(dmin) if(*dmin>d) *dmin=d;
33// if(d<10) printf("%lf %lf %lf %lf\n", d, dminz, eNeighborDzDown, eNeighborDzUp);
34
35 if( d > eTolDX) continue;
36
37 if( -eNeighborDzDown<dminz&&dminz<eNeighborDzUp) flag++;
38 else if(dminz<-eNeighborDzDown){
39 if( CalcDistance(ss, s, ss->Z()+eNeighborDzDown)<eTolDX) flag++;
40 }
41 else if(eNeighborDzUp<dminz){
42 if( CalcDistance(ss, s, ss->Z()-eNeighborDzUp)<eTolDX) flag++;
43 }
44 }
45 return flag;
46}
void d()
Definition: RecDispEX.C:381
double eNeighborDzDown
Definition: EdbEDATrackSet.h:56
double eNeighborDzUp
Definition: EdbEDATrackSet.h:55
Definition: EdbSegP.h:21
s
Definition: check_shower.C:55
ss
Definition: energy.C:62
double CalcDistance(EdbSegP *s1, EdbSegP *s2, double z)
Definition: EdbEDAUtil.C:422
double CalcDmin(EdbSegP *seg1, EdbSegP *seg2, double *dminz=NULL)
Definition: EdbEDAUtil.C:239

◆ Reset()

void EdbEDATrackSelection::Reset ( )
inline
67 {
68 eNsegCut = 3;
69 ePlateUp = 0;
70 ePlateDown = 0;
72 eImpactSearch = 0;
73 eAngularCut = 0;
74 eSideOut = 0;
75 eSideOutPlate = 3;
76 eTolDX = 200.;
77 eTolDTX = 2.0;
78 eTolDTY = 2.0;
79 eTX = 0.0;
80 eTY = 0.0;
82 ePHCut = 0.0;
83 ePHDTRMS = 0.05;
84 eNeighborDzUp = 5000.;
85 eNeighborDzDown = 5000.;
86 eVertex = NULL;
90 }
int ePlateDown
Definition: EdbEDATrackSet.h:40
void SetCondTrackingDefault()
Definition: EdbEDATrackSet.h:124
bool ePredictionScan
if true use GetPatternDataForPrediction( spred.ID(), side, pat ); in FindCandidates (default is false...
Definition: EdbRunTracking.h:62

◆ SetAngle()

void EdbEDATrackSelection::SetAngle ( double  tx,
double  ty 
)
inline
105{ eTX=tx; eTY=ty;}

◆ SetAngularCut()

void EdbEDATrackSelection::SetAngularCut ( bool  b = kTRUE)
inline
102{ eAngularCut = b;}

◆ SetClearPrevious()

void EdbEDATrackSelection::SetClearPrevious ( bool  b = kTRUE)
inline
106{ eClearPrevious = b;}

◆ SetCondBTDefault()

void EdbEDATrackSelection::SetCondBTDefault ( EdbScanCond cond)
inline
159 {
160 cond.SetSigma0( 10., 10., 0.007, 0.007 ); // sigma0 "x, y, tx, ty" at zero angle
161 cond.SetDegrad( 2. ); // sigma(tx) = sigma0*(1+degrad*tx)
162 cond.SetBins(0, 0, 0, 0); // bins in [sigma] for checks
163 cond.SetPulsRamp0( 5., 5. ); // in range (Pmin:Pmax) Signal/All is nearly linear
164 cond.SetPulsRamp04( 5., 5. );
165 cond.SetChi2Max( 6.5 );
166 cond.SetChi2PMax( 6.5 );
167 cond.SetRadX0( 5810. );
168 cond.SetName("OPERA_basetrack");
169 }
void SetPulsRamp0(float p1, float p2)
Definition: EdbScanCond.h:74
void SetChi2Max(float chi2)
Definition: EdbScanCond.h:83
void SetDegrad(float d)
Definition: EdbScanCond.h:71
void SetSigma0(float x, float y, float tx, float ty)
Definition: EdbScanCond.h:62
void SetBins(float bx, float by, float btx, float bty)
Definition: EdbScanCond.h:65
void SetPulsRamp04(float p1, float p2)
Definition: EdbScanCond.h:75
void SetRadX0(float x0)
Definition: EdbScanCond.h:57
void SetChi2PMax(float chi2)
Definition: EdbScanCond.h:84

◆ SetCondMTDefault()

void EdbEDATrackSelection::SetCondMTDefault ( EdbScanCond cond)
inline
147 {
148 cond.SetSigma0( 1., 1., 0.010, 0.010 ); // sigma0 "x, y, tx, ty" at zero angle
149 cond.SetDegrad( 5. ); // sigma(tx) = sigma0*(1+degrad*tx)
150 cond.SetBins(0, 0, 0, 0); //??? // bins in [sigma] for checks
151 cond.SetPulsRamp0( 5., 5. ); // in range (Pmin:Pmax) Signal/All is nearly linear
152 cond.SetPulsRamp04( 5., 5. );
153 cond.SetChi2Max( 6.5 );
154 cond.SetChi2PMax( 6.5 );
155 cond.SetRadX0( 5810. );
156 cond.SetName("OPERA_microtrack");
157 }

◆ SetCondTrackingDefault()

void EdbEDATrackSelection::SetCondTrackingDefault ( )
inline
124 {
127 eRunTracking.eDeltaR = 20.;
128
131
135
139
142
145 }
void SetCondBTDefault(EdbScanCond &cond)
Definition: EdbEDATrackSet.h:159
void SetCondMTDefault(EdbScanCond &cond)
Definition: EdbEDATrackSet.h:147
Float_t ePulsMinDegradBT
(0)
Definition: EdbRunTracking.h:45
Float_t ePulsMinBT
(18)
Definition: EdbRunTracking.h:38
Float_t eDeltaR
(20)
Definition: EdbRunTracking.h:34
EdbScanCond eCondMT
conditions for microtracks
Definition: EdbRunTracking.h:23
Float_t ePreliminaryPulsMinMT
(6) _ preliminary cuts to microtracks candidates for
Definition: EdbRunTracking.h:28
Float_t eChi2MaxMT
(1.6) maximum chi2 accepted between prediction and microtrack candidates
Definition: EdbRunTracking.h:46
Float_t eDeltaRview
(400)
Definition: EdbRunTracking.h:26
Float_t ePreliminaryChi2MaxMT
(1.6) / microtracks and basetracks selection
Definition: EdbRunTracking.h:29
Float_t ePulsMinMT
(10) mimimal number of grains accepted to select microtracks
Definition: EdbRunTracking.h:44
Float_t eDegradSlope
SigmaTX = SigmaTX(0) + degradSlope * bth.
Definition: EdbRunTracking.h:57
EdbScanCond eCondBT
conditions for basetracks
Definition: EdbRunTracking.h:24
Float_t eChi2MaxBT
(1.5) maximum chi2 accepted between prediction and basetrack candidates
Definition: EdbRunTracking.h:40
Float_t eDegradPos
SigmaX = SigmaX(0) + degradPos * mth.
Definition: EdbRunTracking.h:56
Float_t ePulsMinDegradMT
(0)
Definition: EdbRunTracking.h:39
Float_t eDeltaTheta
(0.15) slope acceptance
Definition: EdbRunTracking.h:27

◆ SetDT()

void EdbEDATrackSelection::SetDT ( double  dtx,
double  dty = -1 
)
inline
104{ eTolDTX=dtx; eTolDTY=dty; if(eTolDTY<0) eTolDTY=eTolDTX;}

◆ SetDX()

void EdbEDATrackSelection::SetDX ( double  dx)
inline
103{ eTolDX=dx;}

◆ SetImpactSearch()

void EdbEDATrackSelection::SetImpactSearch ( bool  b = kTRUE,
EdbVertex v = NULL 
)
inline
94{ eImpactSearch = b; eVertex = v;}

◆ SetNeighborSearch()

void EdbEDATrackSelection::SetNeighborSearch ( bool  b = kTRUE,
TObjArray *  selected = NULL,
double  dzup = -1,
double  dzdown = -1 
)
inline
95 {
96 eNeighborSearch = b; eSelected=selected;
97 if(dzup>=0) eNeighborDzUp = dzup;
98 if(dzdown>=0) eNeighborDzDown = dzdown;
99 }

◆ SetNsegCut()

void EdbEDATrackSelection::SetNsegCut ( int  nseg)
inline
91{ eNsegCut = nseg;}

◆ SetPHCut()

void EdbEDATrackSelection::SetPHCut ( double  phcut = 0.0)
inline
107{ ePHCut=phcut;}

◆ SetPHDTRMS()

void EdbEDATrackSelection::SetPHDTRMS ( double  slope)
inline
108{ePHDTRMS=slope;}

◆ SetRunTracking()

void EdbEDATrackSelection::SetRunTracking ( EdbRunTracking rt)
inline
121{eRunTracking=rt;}

◆ SetSideOut()

void EdbEDATrackSelection::SetSideOut ( bool  b = kTRUE)
inline
100{ eSideOut = b;}

◆ SetSideOutPlate()

void EdbEDATrackSelection::SetSideOutPlate ( int  npl)
inline
101{ eSideOutPlate = npl;}

◆ SetUpstreamPlate()

void EdbEDATrackSelection::SetUpstreamPlate ( int  ipl)
inline
92{ ePlateUp = ipl;}

◆ SideOut()

int EdbEDATrackSelection::SideOut ( EdbTrackP t)

Return 1, if the track goes side-out.

  • if require 3 scanning area upstream from the 1st segment.
    EdbPVRec *pvr = gEDA->GetPVR();
    if(pvr==NULL) return 0;
48 {
53
54 int sideout=0;
55
56 // Upstream Side-out check.
57 int npl = 0;
58 EdbSegP *s1 = t->GetSegmentFirst();
59 EdbSegP s=*s1;
60
61 // original. but very slow. EdbPointsBox2D::Xmax().....
62 /*
63 for(int pid=s1->PID()-1; pid>=0; pid--){
64 EdbPattern *p = pvr->GetPattern(pid);
65 if(p->N()==0) continue;
66 s.PropagateTo(p->Z());
67 if(s.X()<p->Xmin()) continue;
68 if(s.X()>p->Xmax()) continue;
69 if(s.Y()<p->Ymin()) continue;
70 if(s.Y()>p->Ymax()) continue;
71
72 npl++;
73 }
74 */
75 /*
76 for(int pid=s1->PID()-1; pid>=0; pid--){
77 EdbEDAArea *a = gEDA->GetAreaSet()->GetArea(pid);
78 s.PropagateTo(a->Z());
79 if(s.X()<a->Xmin()) continue;
80 if(s.X()>a->Xmax()) continue;
81 if(s.Y()<a->Ymin()) continue;
82 if(s.Y()>a->Ymax()) continue;
83
84 npl++;
85 }
86 */
87 for(int ipl=s1->Plate()-1; ipl>0; ipl--){
89 if(a==NULL) continue;
90 s.PropagateTo(a->Z());
91 if(s.X()<a->Xmin()) continue;
92 if(s.X()>a->Xmax()) continue;
93 if(s.Y()<a->Ymin()) continue;
94 if(s.Y()>a->Ymax()) continue;
95
96 npl++;
97 }
98
99
100 if(npl<eSideOutPlate) sideout++;
101 /*
102 // Downstream Side-out check.
103 // - if 0 or only 1 scanning area downstream from the last segment. if the propagated position is in the scanning area of the last plate, eceptionally doesn't cut the track.
104
105 npl = 0;
106 s1 = t->GetSegmentLast();
107 s=*s1;
108 if(s1->PID()==pvr->Npatterns()-1) npl=2;
109 for(int pid=s1->PID()+1; pid<pvr->Npatterns(); pid++){
110 EdbPattern *p = pvr->GetPattern(pid);
111 if(p->N()==0) continue;
112 s.PropagateTo(p->Z());
113 if(s.X()<p->Xmin()) continue;
114 if(s.X()>p->Xmax()) continue;
115 if(s.Y()<p->Ymin()) continue;
116 if(s.Y()>p->Ymax()) continue;
117
118 npl++;
119 if(pid==pvr->Npatterns()-1) npl++; // if the propagated position is in the scanning area of last plate, eceptionally doesn't cut the track.
120 }
121 if(npl<=1) sideout++;
122 */
123 return sideout;
124}
void a()
Definition: check_aligned.C:59
EdbEDAArea * GetAreaIPL(int ipl)
Definition: EdbEDASets.h:135
Definition: EdbEDASets.h:7
EdbEDAAreaSet * GetAreaSet()
Definition: EdbEDATrackSet.h:602
EdbEDATrackSet * GetTrackSet(int i)
Definition: EdbEDA.h:617
Int_t Plate() const
Definition: EdbSegP.h:159
EdbSegP * s1
Definition: tlg2couples.C:29

Member Data Documentation

◆ eAngularCut

int EdbEDATrackSelection::eAngularCut
private

◆ eClearPrevious

int EdbEDATrackSelection::eClearPrevious
private

◆ eImpactSearch

int EdbEDATrackSelection::eImpactSearch
private

◆ eNeighborDzDown

double EdbEDATrackSelection::eNeighborDzDown
private

◆ eNeighborDzUp

double EdbEDATrackSelection::eNeighborDzUp
private

◆ eNeighborSearch

int EdbEDATrackSelection::eNeighborSearch
private

◆ eNsegCut

int EdbEDATrackSelection::eNsegCut
private

◆ ePHCut

double EdbEDATrackSelection::ePHCut
private

◆ ePHDTRMS

double EdbEDATrackSelection::ePHDTRMS
private

◆ ePlateDown

int EdbEDATrackSelection::ePlateDown
private

◆ ePlateUp

int EdbEDATrackSelection::ePlateUp
private

◆ eRunTracking

EdbRunTracking EdbEDATrackSelection::eRunTracking

Run Tracking for microtrack search.

◆ eSelected

TObjArray* EdbEDATrackSelection::eSelected
private

◆ eSideOut

int EdbEDATrackSelection::eSideOut
private

◆ eSideOutPlate

int EdbEDATrackSelection::eSideOutPlate
private

◆ eTolDTX

double EdbEDATrackSelection::eTolDTX
private

◆ eTolDTY

double EdbEDATrackSelection::eTolDTY
private

◆ eTolDX

double EdbEDATrackSelection::eTolDX
private

◆ eTX

double EdbEDATrackSelection::eTX
private

◆ eTY

double EdbEDATrackSelection::eTY
private

◆ eVertex

EdbVertex* EdbEDATrackSelection::eVertex

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