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

#include <EdbPlateTracking.h>

Inheritance diagram for EdbPlateTracking:
Collaboration diagram for EdbPlateTracking:

Public Member Functions

 EdbPlateTracking ()
 
 EdbPlateTracking (EdbPattern &S1, EdbPattern &S2, EdbSegP &prediction)
 
 EdbPlateTracking (EdbPattern &S1, EdbPattern &S2, EdbSegP &prediction, EdbPlateP &plate)
 
int ExtrapolateCond (EdbScanCond &inputcond, int flag, EdbScanCond &outputcond)
 
int FindBestCandidate (EdbPattern &fndbt, EdbSegP &fnd, EdbPattern &cnd, float wmin, float wmindegrad, float chi2max)
 
int FindBestCandidateDS (EdbPattern &fndbt, EdbSegP &fnd, EdbPattern &cnd, float wmin, float wmindegrad, float chi2max, EdbSegP &spred, float maxdmin)
 
int FindCandidateMT (EdbPattern &fnds1, EdbPattern &fnds2, EdbSegP &fnd)
 
int FindCandidateMTDS (EdbPattern &fnds1, EdbPattern &fnds2, EdbSegP &fnd, EdbSegP &spred, float maxdmin)
 
int FindCandidates (EdbSegP &spred, EdbPattern &fndbt, EdbPattern &fnds1, EdbPattern &fnds2)
 
int FindCompliments (EdbSegP &s, EdbPattern &pat, TObjArray &found, float chi2max, TArrayF &chiarr)
 
int FindPrediction (EdbSegP &spred, EdbSegP &fndbt, EdbSegP &fnds1, EdbSegP &fnds2, EdbSegP &snewpred)
 
int FindPredictionDS (EdbSegP &spred, EdbSegP &fndbt, EdbSegP &fnds1, EdbSegP &fnds2, EdbSegP &snewpred, float maxdmin)
 
int FindTrack (EdbTrackP &pred, EdbTrackP &found, EdbPlateP &plate)
 
bool GetSBtreeEntry (int entry, TTree &tsbt)
 
void Print ()
 
void Set0 ()
 
void SetCondBT (EdbScanCond &cond)
 
void SetCondMT (EdbScanCond &cond)
 
void SetPred (const EdbSegP &pred)
 
void TransformFromPlateRS ()
 
int UpdateFlag (int flag, int status)
 
bool UpdateSBtree (TTree &tsbt, int idp[4], int idf[4])
 
virtual ~EdbPlateTracking ()
 

Static Public Member Functions

static void CloseSBtree (TTree *tree)
 
static int GetBTHoles (int flag)
 
static int GetMTHoles (int flag)
 
static TTree * InitSBtree (const char *file_name="sbt.root", const char *mode="RECREATE")
 

Public Attributes

Float_t eChi2MaxBT
 (1.5) maximum chi2 accepted between prediction and basetrack candidates More...
 
Float_t eChi2MaxMT
 (1.6) maximum chi2 accepted between prediction and microtrack candidates More...
 
EdbScanCond eCondBT
 
EdbScanCond eCondMT
 ------— processing parameters (can be default) and itermeadiate results ---------------------— More...
 
Float_t eDegradPos
 SigmaX = SigmaX(0) + degradPos * mth. More...
 
Float_t eDegradSlope
 SigmaTX = SigmaTX(0) + degradSlope * bth. More...
 
Float_t eDeltaR
 (20) More...
 
Int_t eIdf [4]
 
Int_t eIdp [4]
 to read from sbt More...
 
EdbSegP eNext
 next prediction More...
 
EdbPlateP ePlate
 plate geometry and correction parameters to be applied to prediction More...
 
EdbSegP ePred
 ------— input data ------------------------— More...
 
bool ePredictionScan
 if true use GetPatternDataForPrediction( spred.ID(), side, pat ); in FindCandidates (default is false) More...
 
Float_t ePreliminaryChi2MaxMT
 (1.6) / microtracks and basetracks selection More...
 
Float_t ePulsMinBT
 (18) More...
 
Float_t ePulsMinDegradBT
 (0) More...
 
Float_t ePulsMinDegradMT
 (0) More...
 
Float_t ePulsMinMT
 (10) mimimal number of grains accepted to select microtracks More...
 
EdbSegP eS
 ------— output: main result ------------------------— More...
 
EdbSegP eS1
 
EdbPattern eS1cnd
 microtrack candidates passed all cuts More...
 
EdbPattern eS1pre
 the result of the selection of microtracks ordered by chi square More...
 
EdbSegP eS2
 found segments More...
 
EdbPattern eS2cnd
 
EdbPattern eS2pre
 
EdbPattern eScnd
 basetracks candidates passed all cuts More...
 
EdbPattern eSide1
 side 1 microtracks in plate reference system More...
 
EdbPattern eSide2
 side 2 microtracks in plate reference system More...
 
EdbPattern eSpre
 the result of the selection of basetracks ordered by chi square More...
 
Int_t eStatus
 -1-nothing, 0-bt, 1-mt1, 2-mt2 More...
 

Detailed Description

//////////////////////////////////////////////////////////////////////// // EdbPlateTracking // // ////////////////////////////////////////////////////////////////////////

Constructor & Destructor Documentation

◆ EdbPlateTracking() [1/3]

EdbPlateTracking::EdbPlateTracking ( )
inline
77 {
78 Set0();
79 };
void Set0()
Definition: EdbPlateTracking.cxx:28

◆ EdbPlateTracking() [2/3]

EdbPlateTracking::EdbPlateTracking ( EdbPattern S1,
EdbPattern S2,
EdbSegP prediction,
EdbPlateP plate 
)
inline
82 {
83
84 for(int i=0;i<S1.N();i++)
86 for(int i=0;i<S2.N();i++)
88
89 eSide1.SetZ(S1.Z());
90 eSide2.SetZ(S2.Z());
91
92 ePred.Copy(prediction);
94 Set0();
95 };
void Copy(EdbPlateP &p)
Definition: EdbBrick.cxx:24
EdbPattern eSide1
side 1 microtracks in plate reference system
Definition: EdbPlateTracking.h:26
EdbPlateP ePlate
plate geometry and correction parameters to be applied to prediction
Definition: EdbPlateTracking.h:28
EdbSegP ePred
------— input data ------------------------—
Definition: EdbPlateTracking.h:25
EdbPattern eSide2
side 2 microtracks in plate reference system
Definition: EdbPlateTracking.h:27
void Copy(const EdbSegP &s)
Definition: EdbSegP.cxx:105
Float_t Z() const
Definition: EdbPattern.h:84
void SetZ(float z)
Definition: EdbPattern.h:41
Int_t N() const
Definition: EdbPattern.h:86
EdbSegP * GetSegment(int i) const
Definition: EdbPattern.h:66
EdbSegP * AddSegment(int i, EdbSegP &s)
Definition: EdbPattern.cxx:72
Int_t plate
Definition: merge_Energy_SytematicSources_Electron.C:1

◆ EdbPlateTracking() [3/3]

EdbPlateTracking::EdbPlateTracking ( EdbPattern S1,
EdbPattern S2,
EdbSegP prediction 
)
inline
98 {
99
100 for(int i=0;i<S1.N();i++)
101 eSide1.AddSegment(*(S1.GetSegment(i)));
102 for(int i=0;i<S2.N();i++)
103 eSide2.AddSegment(*(S2.GetSegment(i)));
104
105 eSide1.SetZ(S1.Z());
106 eSide2.SetZ(S2.Z());
107
108 ePred.Copy(prediction);
109 Set0();
110 };

◆ ~EdbPlateTracking()

virtual EdbPlateTracking::~EdbPlateTracking ( )
virtual

Member Function Documentation

◆ CloseSBtree()

void EdbPlateTracking::CloseSBtree ( TTree *  tree)
static
744{
745 tree->AutoSave();
746 TFile *f=0;
747 f = tree->GetCurrentFile();
748 if(f) {
749 f->Purge();
750 f->Close();
751 }
752 tree=0;
753}
FILE * f
Definition: RecDispMC.C:150

◆ ExtrapolateCond()

int EdbPlateTracking::ExtrapolateCond ( EdbScanCond inputcond,
int  flag,
EdbScanCond outputcond 
)

tTODO: tuning the dependency of sigma NOTED by Artem: when we do jumping we have no "holes" stored in the flag, so this function do not extrapolate errors correctly so one had to encrease the acceptance manually, possible solution can be to add total plates number before the last found track into the flag and use it only for extrapolation

109{
115
116 int bth = GetBTHoles(flag);
117 int mth = GetMTHoles(flag);
118 outputcond = inputcond;
119 outputcond.SetSigma0( inputcond.SigmaX(0) + eDegradPos * mth,
120 inputcond.SigmaY(0) + eDegradPos * mth,
121 inputcond.SigmaTX(0) + eDegradSlope * bth,
122 inputcond.SigmaTY(0) + eDegradSlope * bth );
123 return(0);
124}
static int GetBTHoles(int flag)
Definition: EdbPlateTracking.h:121
Float_t eDegradSlope
SigmaTX = SigmaTX(0) + degradSlope * bth.
Definition: EdbPlateTracking.h:67
Float_t eDegradPos
SigmaX = SigmaX(0) + degradPos * mth.
Definition: EdbPlateTracking.h:66
static int GetMTHoles(int flag)
Definition: EdbPlateTracking.h:122
float SigmaTX(float ax) const
Definition: EdbScanCond.h:106
float SigmaTY(float ay) const
Definition: EdbScanCond.h:107
float SigmaX(float ax) const
Definition: EdbScanCond.h:102
void SetSigma0(float x, float y, float tx, float ty)
Definition: EdbScanCond.h:62
float SigmaY(float ay) const
Definition: EdbScanCond.h:103

◆ FindBestCandidate()

int EdbPlateTracking::FindBestCandidate ( EdbPattern fndbt,
EdbSegP fnd,
EdbPattern cnd,
float  wmin,
float  wmindegrad,
float  chi2max 
)
323{
324 int n=0;
325 fnd.Set0();
326 fnd.SetChi2(10000.+chi2max);
327 for (int i=0; i<cand.N(); i++) {
328 EdbSegP *s = cand.GetSegment(i);
329 if ( s->W()<wmin+s->Chi2()*wmindegrad ) continue;
330 if ( s->Chi2()>chi2max ) continue;
331 n++;
332 passed.AddSegment(*s);
333 if (s->Chi2()<fnd.Chi2()) fnd.Copy(*s);
334 }
335 return n;
336}
Definition: EdbSegP.h:21
Float_t Chi2() const
Definition: EdbSegP.h:157
void SetChi2(float chi2)
Definition: EdbSegP.h:135
void Set0()
Definition: EdbSegP.cxx:30
s
Definition: check_shower.C:55

◆ FindBestCandidateDS()

int EdbPlateTracking::FindBestCandidateDS ( EdbPattern fndbt,
EdbSegP fnd,
EdbPattern cnd,
float  wmin,
float  wmindegrad,
float  chi2max,
EdbSegP spred,
float  maxdmin = 1. 
)
475{
476 int n=0;
477
478 float eX1=spred.X(),eY1=spred.Y(),eZ1=spred.Z(),eTX1=spred.TX(),eTY1=spred.TY(),dminz;
479 float s1,s2,s1bunsi,s1bunbo,s2bunsi,s2bunbo;
480 float p1x,p1y,p1z,p2x,p2y,p2z,p1p2;
481
482 fnd.Set0();
483 fnd.SetChi2(10000.+chi2max);
484
485 for (int i=0; i<cand.N(); i++) {
486 EdbSegP *s = cand.GetSegment(i);
487 float eX2=s->X(),eY2=s->Y(),eZ2=s->Z(),eTX2=s->TX(),eTY2=s->TY();
488
489 s1bunsi=(eTX2*eTX2+eTY2*eTY2+1)*(eTX1*(eX2-eX1)+eTY1*(eY2-eY1)+eZ2-eZ1) - (eTX1*eTX2+eTY1*eTY2+1)*(eTX2*(eX2-eX1)+eTY2*(eY2-eY1)+eZ2-eZ1);
490 s1bunbo=(eTX1*eTX1+eTY1*eTY1+1)*(eTX2*eTX2+eTY2*eTY2+1) - (eTX1*eTX2+eTY1*eTY2+1)*(eTX1*eTX2+eTY1*eTY2+1);
491 s2bunsi=(eTX1*eTX2+eTY1*eTY2+1)*(eTX1*(eX2-eX1)+eTY1*(eY2-eY1)+eZ2-eZ1) - (eTX1*eTX1+eTY1*eTY1+1)*(eTX2*(eX2-eX1)+eTY2*(eY2-eY1)+eZ2-eZ1);
492 s2bunbo=(eTX1*eTX1+eTY1*eTY1+1)*(eTX2*eTX2+eTY2*eTY2+1) - (eTX1*eTX2+eTY1*eTY2+1)*(eTX1*eTX2+eTY1*eTY2+1);
493 s1=s1bunsi/s1bunbo;
494 s2=s2bunsi/s2bunbo;
495 p1x=eX1+s1*eTX1;
496 p1y=eY1+s1*eTY1;
497 p1z=eZ1+s1*1;
498 p2x=eX2+s2*eTX2;
499 p2y=eY2+s2*eTY2;
500 p2z=eZ2+s2*1;
501 p1p2=sqrt( (p1x-p2x)*(p1x-p2x)+(p1y-p2y)*(p1y-p2y)+(p1z-p2z)*(p1z-p2z) );
502 dminz = eZ1-p1z;
503
504 if (p1p2>maxdmin) continue;
505
506 if ( s->W()<wmin+s->Chi2()*wmindegrad ) continue;
507 if ( s->Chi2()>chi2max ) continue;
508 n++;
509 passed.AddSegment(*s);
510 if (s->Chi2()<fnd.Chi2()) fnd.Copy(*s);
511 }
512
513 return n;
514}
Float_t TX() const
tangens = deltaX/deltaZ
Definition: EdbSegP.h:175
Float_t X() const
Definition: EdbSegP.h:173
Float_t Z() const
Definition: EdbSegP.h:153
Float_t Y() const
Definition: EdbSegP.h:174
Float_t TY() const
tangens = deltaY/deltaZ
Definition: EdbSegP.h:176
EdbSegP * s1
Definition: tlg2couples.C:29
EdbSegP * s2
Definition: tlg2couples.C:30

◆ FindCandidateMT()

int EdbPlateTracking::FindCandidateMT ( EdbPattern fnds1,
EdbPattern fnds2,
EdbSegP fnd 
)
302{
303 EdbSegP s1,s2;
306 // Log(2,"FindCandidateMT","Found %d+%d microtrack candidates after cuts",n1,n2);
307
308 if( n1==0&&n2==0 ) return 0;
309
310 if( s1.Chi2() <= s2.Chi2() ) {
311 fnd.Copy(s1);
312 return 1;
313 }
314 else {
315 fnd.Copy(s2);
316 return 2;
317 }
318 return 0;
319}
int FindBestCandidate(EdbPattern &fndbt, EdbSegP &fnd, EdbPattern &cnd, float wmin, float wmindegrad, float chi2max)
Definition: EdbPlateTracking.cxx:322
Float_t eChi2MaxMT
(1.6) maximum chi2 accepted between prediction and microtrack candidates
Definition: EdbPlateTracking.h:61
EdbPattern eS1cnd
microtrack candidates passed all cuts
Definition: EdbPlateTracking.h:63
Float_t ePulsMinMT
(10) mimimal number of grains accepted to select microtracks
Definition: EdbPlateTracking.h:59
EdbPattern eS2cnd
Definition: EdbPlateTracking.h:64
Float_t ePulsMinDegradMT
(0)
Definition: EdbPlateTracking.h:54

◆ FindCandidateMTDS()

int EdbPlateTracking::FindCandidateMTDS ( EdbPattern fnds1,
EdbPattern fnds2,
EdbSegP fnd,
EdbSegP spred,
float  maxdmin = 1. 
)
452{
453 EdbSegP s1,s2;
454
455 int n1=FindBestCandidateDS( fnds1, s1, eS1cnd, ePulsMinMT, ePulsMinDegradMT, eChi2MaxMT, spred, maxdmin);
456 int n2=FindBestCandidateDS( fnds2, s2, eS2cnd, ePulsMinMT, ePulsMinDegradMT, eChi2MaxMT, spred, maxdmin);
457
458 // Log(2,"FindCandidateMTDS","Found %d+%d microtrack candidates after cuts",n1,n2);
459
460 if( n1==0&&n2==0 ) return 0;
461
462 if( s1.Chi2() <= s2.Chi2() ) {
463 fnd.Copy(s1);
464 return 1;
465 }
466 else {
467 fnd.Copy(s2);
468 return 2;
469 }
470 return 0;
471}
int FindBestCandidateDS(EdbPattern &fndbt, EdbSegP &fnd, EdbPattern &cnd, float wmin, float wmindegrad, float chi2max, EdbSegP &spred, float maxdmin)
Definition: EdbPlateTracking.cxx:474

◆ FindCandidates()

int EdbPlateTracking::FindCandidates ( EdbSegP spred,
EdbPattern fndbt,
EdbPattern fnds1,
EdbPattern fnds2 
)

Find microtracks and basetracks for the prediction segment "spred" Selection criteria: 1) select all microtracks with (puls >= ePreliminaryPulsMinMT) (6) and (chi2 < ePreliminaryChi2MaxMT) for both side 2) select all basetracks by using the selected microtracks, with (chi2_bt < eChi2MaxBT) Microtracks and basetracks in output are sorted according to chi2 (ascending order = the first one is the best one)

Input: spred - track prediction Output: fnds1 - microtracks having puls>=ePreliminaryPulsMinMT and chi2<ePreliminaryChi2MaxMT in the top side of the emulsion plate fnds2 - microtracks having puls>=ePreliminaryPulsMinMT and chi2<ePreliminaryChi2MaxMT in the bottom side of the emulsion plate fndbt - basetracks built from the microtracks and (chi2_bt < eChi2MaxBT)

173{
188
189
190 EdbScanCond condBT;
191 ExtrapolateCond(eCondBT,spred.Flag(),condBT);
192 if (gEDBDEBUGLEVEL>=3) condBT.Print();
193
194 //spred.SetZ(107.); // TODO! sonoqui
195 spred.SetErrors();
196 condBT.FillErrorsCov( spred.TX(), spred.TY(), spred.COV() );
197
198 //EdbPVRec aview; //with 2 patterns of preselected microtracks
199 //aview.AddPattern( new EdbPattern(0,0,214)); // TODO! sequence??
200 //aview.AddPattern( new EdbPattern(0,0,0) );
201
202
203 for(int side=1; side<=2; side++) {
204
205
206 if(side==1)
207 {
208 for(int i=0; i<eSide1.N(); i++) {
210 s->SetErrors();
211 eCondMT.FillErrorsCov( spred.TX(), spred.TY(), s->COV() );
212 }
213 eSide1.FillCell(10,10,0.01,0.01); //divide view on this cells
214 TArrayF chi2arr(10000); //TODO!
215 TObjArray found;
216 FindCompliments( spred, eSide1, found, ePreliminaryChi2MaxMT, chi2arr );
217 for(int j=0; j<found.GetEntries(); j++) {
218 EdbSegP *s = (EdbSegP *)(found.At(j));
219 s->SetChi2(chi2arr[j]);
220 if (side==1) fnds1.AddSegment(*s);
221 else if(side==2) fnds2.AddSegment(*s);
222 }
223 }
224 else
225 {
226 for(int i=0; i<eSide2.N(); i++) {
228 s->SetErrors();
229 eCondMT.FillErrorsCov( spred.TX(), spred.TY(), s->COV() );
230 }
231 eSide2.FillCell(10,10,0.01,0.01); //divide view on this cells
232 TArrayF chi2arr(10000); //TODO!
233 TObjArray found;
234 FindCompliments( spred, eSide2, found, ePreliminaryChi2MaxMT, chi2arr );
235 for(int j=0; j<found.GetEntries(); j++) {
236 EdbSegP *s = (EdbSegP *)(found.At(j));
237 s->SetChi2(chi2arr[j]);
238 if (side==1) fnds1.AddSegment(*s);
239 else if(side==2) fnds2.AddSegment(*s);
240 }
241 }
242
243
244 }//end loop on two emulsion layers
245
246 // filling fndbt
247
248 EdbPattern bt;
249
250 for(int is1=0; is1<fnds1.N(); is1++) {
251 for(int is2=0; is2<fnds2.N(); is2++) {
252
253 EdbSegP *s1 = fnds1.GetSegment(is1);
254 EdbSegP *s2 = fnds2.GetSegment(is2);
255
256 float dx1=s1->X()-(spred.X()+spred.TX()*(s1->Z()-spred.Z()));
257 float dy1=s1->Y()-(spred.Y()+spred.TY()*(s1->Z()-spred.Z()));
258 float dx2=s2->X()-(spred.X()+spred.TX()*(s2->Z()-spred.Z()));
259 float dy2=s2->Y()-(spred.Y()+spred.TY()*(s2->Z()-spred.Z()));
260 float r = Sqrt( (dx1-dx2)*(dx1-dx2) + (dy1-dy2)*(dy1-dy2) );
261
262 if(r<eDeltaR) { // has good BT
263 EdbSegP s3;
264 s3.Copy(spred);
265 s3.SetX( 0.5*(s1->X() + s2->X()) );
266 s3.SetY( 0.5*(s1->Y() + s2->Y()) );
267 s3.SetZ( 0.5*(s1->Z() + s2->Z()) );
268 s3.SetTX( (s2->X() - s1->X()) / (s2->Z() - s1->Z()) );
269 s3.SetTY( (s2->Y() - s1->Y()) / (s2->Z() - s1->Z()) );
270 s3.SetW(s1->W()+s2->W());
271 s3.SetVid(s1->ID(),s2->ID());//ale, antonia
272
273 s3.SetFlag(is2*10000+is1);
274
275 EdbSegP s4(spred);
276 float chi = EdbTrackFitter::Chi2Seg(&s4, &s3); // depends on the sequence (best defined should be first)
277 //float chi = EdbTrackFitter::Chi2SegM(spred, s3, s4, eCondBT, eCondBT);
278 if(chi<eChi2MaxBT) {
279 s3.SetChi2(chi);
280 bt.AddSegment(s3);
281 }
282 }
283 }
284 }
285
286 TArrayF chi2arr(bt.N());
287 for (int i=0;i<bt.N();i++) chi2arr[i]=bt.GetSegment(i)->Chi2();
288 TArrayI ind(bt.N());
289 TMath::Sort(bt.N(),chi2arr.GetArray(),ind.GetArray(),false);
290 for(int i=0; i<bt.N(); i++) {
291 EdbSegP *s = bt.GetSegment(ind[i]);
292 fndbt.AddSegment(*s);
293 }
294
295 // Log(2,"FindCandidates","Found %d basetrack candidate in %d+%d preselected microtracks",fndbt.N(),fnds1.N(),fnds2.N());
296
297 return 1; //TODO!
298}
Definition: EdbPattern.h:273
void FillCell(float stepx, float stepy, float steptx, float stepty)
Definition: EdbPattern.cxx:1323
Float_t eChi2MaxBT
(1.5) maximum chi2 accepted between prediction and basetrack candidates
Definition: EdbPlateTracking.h:55
Float_t eDeltaR
(20)
Definition: EdbPlateTracking.h:49
Float_t ePreliminaryChi2MaxMT
(1.6) / microtracks and basetracks selection
Definition: EdbPlateTracking.h:44
EdbScanCond eCondBT
Definition: EdbPlateTracking.h:39
EdbScanCond eCondMT
------— processing parameters (can be default) and itermeadiate results ---------------------—
Definition: EdbPlateTracking.h:38
int ExtrapolateCond(EdbScanCond &inputcond, int flag, EdbScanCond &outputcond)
Definition: EdbPlateTracking.cxx:108
int FindCompliments(EdbSegP &s, EdbPattern &pat, TObjArray &found, float chi2max, TArrayF &chiarr)
Definition: EdbPlateTracking.cxx:127
Definition: EdbScanCond.h:10
void Print() const
Definition: EdbScanCond.cxx:50
void FillErrorsCov(float tx, float ty, TMatrixD &cov)
Definition: EdbScanCond.cxx:161
void SetY(Float_t y)
Definition: EdbSegP.h:178
void SetErrors()
Definition: EdbSegP.h:90
void SetTX(Float_t tx)
Definition: EdbSegP.h:179
Int_t ID() const
Definition: EdbSegP.h:147
void SetX(Float_t x)
Definition: EdbSegP.h:177
TMatrixD & COV() const
Definition: EdbSegP.h:123
void SetZ(float z)
Definition: EdbSegP.h:125
void SetTY(Float_t ty)
other functions
Definition: EdbSegP.h:180
void SetW(float w)
Definition: EdbSegP.h:132
Float_t W() const
Definition: EdbSegP.h:151
void SetVid(int vid, int sid)
Definition: EdbSegP.h:137
void SetFlag(int flag)
Definition: EdbSegP.h:130
Int_t Flag() const
Definition: EdbSegP.h:149
static float Chi2Seg(EdbSegP *s1, EdbSegP *s2)
Definition: EdbTrackFitter.cxx:62
gEDBDEBUGLEVEL
Definition: energy.C:7
void r(int rid=2)
Definition: test.C:201

◆ FindCompliments()

int EdbPlateTracking::FindCompliments ( EdbSegP s,
EdbPattern pat,
TObjArray &  found,
float  chi2max,
TArrayF &  chiarr 
)

return found sorted by increasing chi2

128{
130
131 int nfound=0;
132 int maxcand=chiarr.GetSize();
133 TArrayF chi2arr(maxcand);
134 TObjArray arr(maxcand);
135 TArrayI ind(maxcand);
136
137 int nseg = pat.FindCompliments(s,arr,30,200); // acceptance (prelim): s.SX()*30; s.STX*200
138 // printf("\nnseg = %d\n",nseg);
139 if(nseg>maxcand) {
140 printf("Warning!: Too many segments %d, accept only the first %d \n", nseg, maxcand);
141 nseg = maxcand;
142 }
143 if(nseg<=0) return 0;
144
145 EdbSegP *s2=0;
146 for(int j=0; j<nseg; j++) {
147 s2 = (EdbSegP *)arr.At(j);
148 EdbSegP s3;
149 s3.Copy(s);
150 chi2arr[j] = EdbTrackFitter::Chi2Seg(&s3, s2);
151 //chi2arr[j] = EdbTrackFitter::Chi2SegM(s, *s2,s3,eCondBT,eCondMT);
152 }
153 TMath::Sort(nseg,chi2arr.GetArray(),ind.GetArray(),0);
154 // printf("pred = %f %f %f %f\n",s.X(),s.Y(),s.TX(),s.TY());
155 for(int j=0; j<nseg; j++) {
156 s2 = (EdbSegP *)arr.At(ind[j]);
157 // printf("j = %d, ind = %d, chi2 = %f %f\n",j,ind[j],chi2arr[ind[j]],chi2max);
158 // printf("Tx = %f %f %f %f\n",s2->X(),s2->Y(),s2->TX(),s2->TY());
159 if(chi2arr[ind[j]] > chi2max ) break;
160 chiarr[j] = chi2arr[ind[j]];
161 s2->SetMC(s.MCEvt(),s.MCTrack());
162 found.Add(s2);
163 nfound++;
164 }
165
166 // printf("nfound = %d\n",nfound);
167 return nfound;
168}
int FindCompliments(EdbSegP &s, TObjArray &arr, float nsig, float nsigt)
Definition: EdbPattern.cxx:1354
void SetMC(int mEvt, int mTrack)
Definition: EdbSegP.h:141

◆ FindPrediction()

int EdbPlateTracking::FindPrediction ( EdbSegP spred,
EdbSegP fndbt,
EdbSegP fnds1,
EdbSegP fnds2,
EdbSegP snewpred 
)

Select the best (micro or base) track matching with the prediction and prepare for a new search.

Selection criteria: 1) Call FindCandidates having the list of basetrack and microtrack candidates 2) Call FindCandidateBT which looks for the best basetrack, if any 3) If no basetrack is found, call FindCandidateMT which looks fot the best microtrack, if any. Microtracks accepted shold satisfy the following cut: (puls >= ePulsMinMT) (10) and (chi2 < eChi2MaxMT)

Input: spred - track prediction Output:

  • if a basetrack is found (status 0): fndbt - basetrack found fnds1 - microtrack top contained in the found basetrack fnds2 - microtrack bottom contained in the found basetrack snewpred - fndbt with flag equal to 0
  • if a microtrack top is found (status 1): fndbt - dummy fnds1 - microtrack top found fnds2 - dummy snewpred - a track with slopes from the prediction spred and positions from an extrapolation of fnds1. The flag is updated by UpdateFlag
  • if a microtrack bottom is found (status 2): fndbt - dummy fnds1 - dummy fnds2 - microtrack bottom found snewpred - a track with slopes from the prediction spred and positions from an extrapolation of fnds2. The flag is updated by UpdateFlag
  • if nothing is found (status -1): fndbt - dummy fnds1 - dummy fnds2 - dummy snewpred - the prediction spred. The flag is updated by UpdateFlag

Return: -1: no track found 0: basetrack found 1: microtrack top found 2: microtrack bottom found

342{
386
387 //EdbPattern vfndbt,vfnds1,vfnds2;
388
389 eStatus = -1;
390 SetPred(spred);
391 FindCandidates( spred, eSpre, eS1pre, eS2pre ); //eSpre , eS1pre , eS2pre riempiti
392
393 EdbSegP fnd;
395 if ( nbt > 0 ) {
396 eS.Copy(fnd);
397 eS1.Copy(*(eS1pre.GetSegment(fnd.Flag()%10000)));
398 eS2.Copy(*(eS2pre.GetSegment(fnd.Flag()/10000)));
399 eS.SetFlag(0);
400 eNext.Copy(fnd);
401 eNext.SetFlag(UpdateFlag(spred.Flag(),0)); // if bt found : bth=0, mth=0, tb=0
402 eStatus = 0;
403 }
404
405 int if_mt = FindCandidateMT(eS1pre,eS2pre,fnd);
406 if(eStatus!=-1) goto RESUME;
407
408 switch(if_mt) {
409 case 0: // find nothing
410 eNext.Copy(spred);
411 eNext.SetFlag(UpdateFlag(spred.Flag(),-1)); // hole: if not found: bth++, mth++, tb= keep last value
412 eNext.SetW(0);
413 eStatus = -1; goto RESUME;
414 case 1: // best microtrack is on the 1-st side
415 eS1.Copy(fnd);
416 eNext.Copy(spred);
417 eNext.SetX( fnd.X() + spred.TX()*(spred.Z()-fnd.Z()) );
418 eNext.SetY( fnd.Y() + spred.TY()*(spred.Z()-fnd.Z()) );
419 eNext.SetZ(spred.Z());
420 eNext.SetFlag(UpdateFlag(spred.Flag(),1)); // if mt found : bth++, mth=0, tb=1
421 eNext.SetW(fnd.W());
422 eNext.SetID(fnd.ID());//ale,antonia
423 eStatus = 1; goto RESUME;
424 case 2: // best microtrack is on the 2-d side
425 eS2.Copy(fnd);
426 eNext.Copy(spred);
427 eNext.SetX( fnd.X() + spred.TX()*(spred.Z()-fnd.Z()) );
428 eNext.SetY( fnd.Y() + spred.TY()*(spred.Z()-fnd.Z()) );
429 eNext.SetZ(spred.Z());
430 eNext.SetFlag(UpdateFlag(spred.Flag(),2)); // if mt found : bth++, mth=0, tb=2
431 eNext.SetW(fnd.W());
432 eNext.SetID(fnd.ID());//ale,antonia
433 eStatus = 2; goto RESUME;
434 }
435
436 RESUME:
437
438 /* Log(2,"FindPrediction","status = %d, good candidates [s:s1:s2] %d:%d:%d ; preliminary [s:s1:s2] %d:%d:%d",
439 eStatus,
440 eScnd.N(),eS1cnd.N(),eS2cnd.N(),
441 eSpre.N(),eS1pre.N(),eS2pre.N()
442 );*/
443 snewpred.Copy(eNext);
444 fndbt.Copy(eS);
445 fnds1.Copy(eS1);
446 fnds2.Copy(eS2);
447 return eStatus;
448}
EdbPattern eSpre
the result of the selection of basetracks ordered by chi square
Definition: EdbPlateTracking.h:51
Float_t ePulsMinDegradBT
(0)
Definition: EdbPlateTracking.h:60
int FindCandidateMT(EdbPattern &fnds1, EdbPattern &fnds2, EdbSegP &fnd)
Definition: EdbPlateTracking.cxx:301
Int_t eStatus
-1-nothing, 0-bt, 1-mt1, 2-mt2
Definition: EdbPlateTracking.h:34
EdbSegP eNext
next prediction
Definition: EdbPlateTracking.h:33
EdbSegP eS2
found segments
Definition: EdbPlateTracking.h:32
EdbSegP eS1
Definition: EdbPlateTracking.h:32
int UpdateFlag(int flag, int status)
Definition: EdbPlateTracking.cxx:89
EdbPattern eS2pre
Definition: EdbPlateTracking.h:47
EdbSegP eS
------— output: main result ------------------------—
Definition: EdbPlateTracking.h:32
void SetPred(const EdbSegP &pred)
Definition: EdbPlateTracking.cxx:53
EdbPattern eS1pre
the result of the selection of microtracks ordered by chi square
Definition: EdbPlateTracking.h:46
Float_t ePulsMinBT
(18)
Definition: EdbPlateTracking.h:53
int FindCandidates(EdbSegP &spred, EdbPattern &fndbt, EdbPattern &fnds1, EdbPattern &fnds2)
Definition: EdbPlateTracking.cxx:172
EdbPattern eScnd
basetracks candidates passed all cuts
Definition: EdbPlateTracking.h:57
void SetID(int id)
Definition: EdbSegP.h:128

◆ FindPredictionDS()

int EdbPlateTracking::FindPredictionDS ( EdbSegP spred,
EdbSegP fndbt,
EdbSegP fnds1,
EdbSegP fnds2,
EdbSegP snewpred,
float  maxdmin = 1. 
)

Select the best (micro or base) track matching with the prediction and prepare for a new search.

Selection criteria: 1) Call FindCandidates having the list of basetrack and microtrack candidates 2) Call FindCandidateBT which looks for the best basetrack, if any 3) If no basetrack is found, call FindCandidateMT which looks fot the best microtrack, if any. Microtracks accepted shold satisfy the following cut: (puls >= ePulsMinMT) (10) and (chi2 < eChi2MaxMT)

Input: spred - track prediction Output:

  • if a basetrack is found (status 0): fndbt - basetrack found fnds1 - microtrack top contained in the found basetrack fnds2 - microtrack bottom contained in the found basetrack snewpred - fndbt with flag equal to 0
  • if a microtrack top is found (status 1): fndbt - dummy fnds1 - microtrack top found fnds2 - dummy snewpred - a track with slopes from the prediction spred and positions from an extrapolation of fnds1. The flag is updated by UpdateFlag
  • if a microtrack bottom is found (status 2): fndbt - dummy fnds1 - dummy fnds2 - microtrack bottom found snewpred - a track with slopes from the prediction spred and positions from an extrapolation of fnds2. The flag is updated by UpdateFlag
  • if nothing is found (status -1): fndbt - dummy fnds1 - dummy fnds2 - dummy snewpred - the prediction spred. The flag is updated by UpdateFlag

Return: -1: no track found 0: basetrack found 1: microtrack top found 2: microtrack bottom found

518{
562
563 //EdbPattern vfndbt,vfnds1,vfnds2;
564
565 eStatus = -1;
566 SetPred(spred);
567 FindCandidates( spred, eSpre, eS1pre, eS2pre ); //eSpre , eS1pre , eS2pre riempiti
568
569
570 EdbSegP fnd;
571 int nbt = FindBestCandidateDS(eSpre,fnd, eScnd, ePulsMinBT, ePulsMinDegradBT, eChi2MaxBT, spred, maxdmin);
572 if ( nbt > 0 ) {
573 eS.Copy(fnd);
574 eS1.Copy(*(eS1pre.GetSegment(fnd.Flag()%10000)));
575 eS2.Copy(*(eS2pre.GetSegment(fnd.Flag()/10000)));
576 eS.SetFlag(0);
577 eNext.Copy(fnd);
578 eNext.SetFlag(UpdateFlag(spred.Flag(),0)); // if bt found : bth=0, mth=0, tb=0
579 eStatus = 0;
580 }
581
582 int if_mt = FindCandidateMTDS(eS1pre,eS2pre,fnd,spred,maxdmin);
583 if(eStatus!=-1) goto RESUME;
584
585 switch(if_mt) {
586 case 0: // find nothing
587 eNext.Copy(spred);
588 eNext.SetFlag(UpdateFlag(spred.Flag(),-1)); // hole: if not found: bth++, mth++, tb= keep last value
589 eNext.SetW(0);
590 eStatus = -1; goto RESUME;
591 case 1: // best microtrack is on the 1-st side
592 eS1.Copy(fnd);
593 eNext.Copy(spred);
594 eNext.SetX( fnd.X() + spred.TX()*(spred.Z()-fnd.Z()) );
595 eNext.SetY( fnd.Y() + spred.TY()*(spred.Z()-fnd.Z()) );
596 eNext.SetZ(spred.Z());
597 eNext.SetFlag(UpdateFlag(spred.Flag(),1)); // if mt found : bth++, mth=0, tb=1
598 eNext.SetW(fnd.W());
599 eNext.SetID(fnd.ID());//ale,antonia
600 eStatus = 1; goto RESUME;
601 case 2: // best microtrack is on the 2-d side
602 eS2.Copy(fnd);
603 eNext.Copy(spred);
604 eNext.SetX( fnd.X() + spred.TX()*(spred.Z()-fnd.Z()) );
605 eNext.SetY( fnd.Y() + spred.TY()*(spred.Z()-fnd.Z()) );
606 eNext.SetZ(spred.Z());
607 eNext.SetFlag(UpdateFlag(spred.Flag(),2)); // if mt found : bth++, mth=0, tb=2
608 eNext.SetW(fnd.W());
609 eNext.SetID(fnd.ID());//ale,antonia
610 eStatus = 2; goto RESUME;
611 }
612
613 RESUME:
614
615 /* Log(2,"FindPredictionDS","status = %d, good candidates [s:s1:s2] %d:%d:%d ; preliminary [s:s1:s2] %d:%d:%d",
616 eStatus,
617 eScnd.N(),eS1cnd.N(),eS2cnd.N(),
618 eSpre.N(),eS1pre.N(),eS2pre.N()
619 );*/
620 snewpred.Copy(eNext);
621 fndbt.Copy(eS);
622 fnds1.Copy(eS1);
623 fnds2.Copy(eS2);
624 return eStatus;
625}
int FindCandidateMTDS(EdbPattern &fnds1, EdbPattern &fnds2, EdbSegP &fnd, EdbSegP &spred, float maxdmin)
Definition: EdbPlateTracking.cxx:451

◆ FindTrack()

int EdbPlateTracking::FindTrack ( EdbTrackP pred,
EdbTrackP found,
EdbPlateP plate 
)

look for tracks in this plate track - input track in brick RS - will be updated on output plate - all plate parameters including affine transformation plate-to-brick

771{
775
776 int status=-100;
777 float DZmax = 1350*100;
778
779 EdbAffine2D p2b(*(plate.GetAffineXY())); // from plate to brick
780 EdbAffine2D b2p(p2b); b2p.Invert(); // from brick to plate
781
782 EdbSegP ps;
783 float dz = pred.MakePredictionTo(plate.Z(), ps);
784
785 if(Abs(dz)>DZmax) return status;
786 if(GetBTHoles(pred.Flag())>5) return status;
787 if(GetMTHoles(pred.Flag())>3) return status;
788
789 ps.SetFlag(pred.Flag());
790 ps.Transform(&b2p); // plate.Transoform(seg) ???
791
792 EdbSegP fndbt, fnds1, fnds2, snewpred;
793 status = FindPrediction( ps, fndbt, fnds1, fnds2, snewpred ); // -1: not found; 0-bt, 1-bot, 2-top
794 found.SetFlag(snewpred.Flag());
795
796 TransformFromPlateRS(); // transform all components into brick RS
797
798 if(status>=0) {
799 found.AddSegment( new EdbSegP(eNext) );
800 found.AddSegmentF( new EdbSegP(ePred) ); // add prediction as "fitted segment" because it is an extrapolation
801 found.SetSegmentsTrack();
802 }
803
804 // Log(2,"EdbPlateTracking::FindTracks","status = %d",status);
805 return status;
806}
brick dz
Definition: RecDispMC.C:107
Definition: EdbAffine.h:17
int FindPrediction(EdbSegP &spred, EdbSegP &fndbt, EdbSegP &fnds1, EdbSegP &fnds2, EdbSegP &snewpred)
Definition: EdbPlateTracking.cxx:341
void TransformFromPlateRS()
Definition: EdbPlateTracking.cxx:809
virtual void Transform(const EdbAffine2D *a)
void AddSegmentF(EdbSegP *s)
Definition: EdbPattern.h:233
void AddSegment(EdbSegP *s)
Definition: EdbPattern.h:214
float MakePredictionTo(Float_t z, EdbSegP &ss)
Definition: EdbPattern.cxx:1050
int SetSegmentsTrack(int id)
Definition: EdbPattern.h:246
float DZmax
Definition: check_vertex.C:28

◆ GetBTHoles()

static int EdbPlateTracking::GetBTHoles ( int  flag)
inlinestatic
121{ return(flag/10000); }

◆ GetMTHoles()

static int EdbPlateTracking::GetMTHoles ( int  flag)
inlinestatic
122{ return((flag/100)%100); }

◆ GetSBtreeEntry()

bool EdbPlateTracking::GetSBtreeEntry ( int  entry,
TTree &  tsbt 
)
712{
713 EdbSegP *s_pred = &ePred, *s_bt = &eS, *s_mt1=&eS1, *s_mt2=&eS2, *s_next=&eNext;
714
715 tsbt.SetBranchAddress("idpred",eIdp);
716 tsbt.SetBranchAddress("idfound",eIdf);
717 tsbt.SetBranchAddress("stat",&eStatus);
718 tsbt.SetBranchAddress("pred.", &s_pred);
719 tsbt.SetBranchAddress("s.", &s_bt);
720 tsbt.SetBranchAddress("s1.", &s_mt1);
721 tsbt.SetBranchAddress("s2.", &s_mt2);
722 tsbt.SetBranchAddress("next.", &s_next);
723
724 TClonesArray *s_cnd = eScnd.GetSegments();
725 tsbt.SetBranchAddress("scnd", &s_cnd);
726 TClonesArray *s1_cnd = eS1cnd.GetSegments();
727 tsbt.SetBranchAddress("s1cnd", &s1_cnd);
728 TClonesArray *s2_cnd = eS2cnd.GetSegments();
729 tsbt.SetBranchAddress("s2cnd", &s2_cnd);
730
731 TClonesArray *s_pre = eSpre.GetSegments();
732 tsbt.SetBranchAddress("spre", &s_pre);
733 TClonesArray *s1_pre = eS1pre.GetSegments();
734 tsbt.SetBranchAddress("s1pre", &s1_pre);
735 TClonesArray *s2_pre = eS2pre.GetSegments();
736 tsbt.SetBranchAddress("s2pre", &s2_pre);
737
738 tsbt.GetEntry(entry);
739 return true;
740}
TLegendEntry * entry
Definition: Canv_SYSTEMATICS_ALLCOMBINED__RMSEnergy__vs__Energy__ELECTRON.C:130
Int_t eIdp[4]
to read from sbt
Definition: EdbPlateTracking.h:69
Int_t eIdf[4]
Definition: EdbPlateTracking.h:70
TClonesArray * GetSegments() const
Definition: EdbPattern.h:69

◆ InitSBtree()

TTree * EdbPlateTracking::InitSBtree ( const char *  file_name = "sbt.root",
const char *  mode = "RECREATE" 
)
static
631{
632 const char *tree_name="sbt";
633 TTree *tree=0;
634 if (!tree) {
635 TFile *f = new TFile(file_name,mode);
636 if (f) tree = (TTree*)f->Get(tree_name);
637 if(!tree) {
638
639 f->cd();
640 tree = new TTree(tree_name,tree_name);
641 tree->SetMaxTreeSize(15000000000LL); //set 15 Gb file size limit)
642
643 Int_t idp[4], idf[4], stat;
644 EdbSegP *s_pred = 0, *s_bt = 0, *s_mt1=0, *s_mt2=0, *s_next=0;
645 TClonesArray *scnd = new TClonesArray("EdbSegP");
646 TClonesArray *s1cnd = new TClonesArray("EdbSegP");
647 TClonesArray *s2cnd = new TClonesArray("EdbSegP");
648 TClonesArray *spre = new TClonesArray("EdbSegP");
649 TClonesArray *s1pre = new TClonesArray("EdbSegP");
650 TClonesArray *s2pre = new TClonesArray("EdbSegP");
651
652 tree->Branch("idpred",idp,"idp[4]/I");
653 tree->Branch("idfound",idf,"idf[4]/I");
654 tree->Branch("stat",&stat,"stat/I");
655 tree->Branch("pred.", "EdbSegP", &s_pred);
656 tree->Branch("s.", "EdbSegP", &s_bt);
657 tree->Branch("s1.", "EdbSegP", &s_mt1);
658 tree->Branch("s2.", "EdbSegP", &s_mt2);
659 tree->Branch("next.", "EdbSegP", &s_next);
660 tree->Branch("scnd",&scnd);
661 tree->Branch("s1cnd",&s1cnd);
662 tree->Branch("s2cnd",&s2cnd);
663 tree->Branch("spre",&spre);
664 tree->Branch("s1pre",&s1pre);
665 tree->Branch("s2pre",&s2pre);
666 tree->Write();
667 }
668 }
669
670 if(!tree) Log(1,"InitSBtree","ERROR!!! can't initialize tree at %s as %s\n",file_name,mode);
671 return tree;
672}
bool Log(int level, const char *location, const char *fmt,...)
Definition: EdbLog.cxx:75

◆ Print()

void EdbPlateTracking::Print ( )
68{
69 printf("EdbPlateTracking selection criteria:\n");
70 printf("Preliminary mt selection:\n");
71 //printf(" eDeltaRview = %5.3f\n", eDeltaRview);
72 //printf(" eDeltaTheta = %5.3f\n", eDeltaTheta);
73 //printf(" ePreliminaryPulsMinMT = %3.1f\n", ePreliminaryPulsMinMT);
74 printf(" ePreliminaryChi2MaxMT = %5.3f\n", ePreliminaryChi2MaxMT);
75 printf("Preliminary bt selection:\n");
76 printf(" eDeltaR = %3.1f\n", eDeltaR);
77 printf(" eChi2MaxBT = %5.3f\n", eChi2MaxBT);
78 printf("Final bt selection:\n");
79 printf(" ePulsMinBT = %3.1f\n", ePulsMinBT);
80 printf(" ePulsMinDegradBT = %3.1f\n", ePulsMinDegradBT);
81 printf("Final mt selection:\n");
82 printf(" ePulsMinMT = %3.1f\n", ePulsMinMT);
83 printf(" ePulsMinDegradMT = %3.1f\n", ePulsMinDegradMT);
84 printf(" eChi2MaxMT = %5.3f\n", eChi2MaxMT);
85 printf("\n");
86}

◆ Set0()

void EdbPlateTracking::Set0 ( )
29{
32
33 // eDeltaRview = 400;
34 // eDeltaTheta = 0.15;
35 // ePreliminaryPulsMinMT = 5;
37
38 eDeltaR = 20;
39 ePulsMinBT = 18;
41 eChi2MaxBT = 2.5;
42
43 eChi2MaxMT = 1.6;
45 ePulsMinMT = 10.;
46
47 eDegradPos = 0;
48 eDegradSlope = 0;
49 ePredictionScan = false;
50}
bool ePredictionScan
if true use GetPatternDataForPrediction( spred.ID(), side, pat ); in FindCandidates (default is false...
Definition: EdbPlateTracking.h:72
void SetDefault()
Definition: EdbScanCond.cxx:16

◆ SetCondBT()

void EdbPlateTracking::SetCondBT ( EdbScanCond cond)
inline
116{ eCondBT = cond; }

◆ SetCondMT()

void EdbPlateTracking::SetCondMT ( EdbScanCond cond)
inline
115{ eCondMT = cond; }

◆ SetPred()

void EdbPlateTracking::SetPred ( const EdbSegP pred)
54{
55 ePred.Copy(pred);
56 eNext.Set0();
57 eS.Set0(); eS1.Set0(); eS2.Set0();
58 eS1cnd.Reset();
59 eS2cnd.Reset();
60 eScnd.Reset();
61 eS1pre.Reset();
62 eS2pre.Reset();
63 eSpre.Reset();
64}
void Reset()
Definition: EdbPattern.cxx:1509

◆ TransformFromPlateRS()

void EdbPlateTracking::TransformFromPlateRS ( )
810{
811 EdbAffine2D p2b(*(ePlate.GetAffineXY())); // from plate to brick
812
813 ePred.Transform(&p2b);
814 eNext.Transform(&p2b);
815 eS.Transform(&p2b);
816 eS1.Transform(&p2b);
817 eS2.Transform(&p2b);
818 eScnd.Transform(&p2b);
819 eS1cnd.Transform(&p2b);
820 eS2cnd.Transform(&p2b);
821 eSpre.Transform(&p2b);
822 eS1pre.Transform(&p2b);
823 eS2pre.Transform(&p2b);
824
825 ePred.SetPID( ePlate.ID() );
826 eNext.SetPID( ePlate.ID() );
827 eS.SetPID( ePlate.ID() );
828 eS1.SetPID( ePlate.ID() );
829 eS2.SetPID( ePlate.ID() );
830 for(int i=0; i<eScnd.N(); i++) eScnd.GetSegment(i)->SetPID( ePlate.ID() );
831 for(int i=0; i<eS1cnd.N(); i++) eS1cnd.GetSegment(i)->SetPID( ePlate.ID() );
832 for(int i=0; i<eS2cnd.N(); i++) eS2cnd.GetSegment(i)->SetPID( ePlate.ID() );
833 for(int i=0; i<eSpre.N(); i++) eSpre.GetSegment(i)->SetPID( ePlate.ID() );
834 for(int i=0; i<eS1pre.N(); i++) eS1pre.GetSegment(i)->SetPID( ePlate.ID() );
835 for(int i=0; i<eS2pre.N(); i++) eS2pre.GetSegment(i)->SetPID( ePlate.ID() );
836
837 ePred.SetZ( ePlate.Z() );
838 eNext.SetZ( ePlate.Z() );
839 eS.SetZ( ePlate.Z() );
840 eS1.SetZ( ePlate.GetLayer(1)->Z() + ePlate.Z() );
841 eS2.SetZ( ePlate.GetLayer(2)->Z() + ePlate.Z() );
842 for(int i=0; i<eScnd.N(); i++) eScnd.GetSegment(i)->SetZ( ePlate.Z() );
843 for(int i=0; i<eS1cnd.N(); i++) eS1cnd.GetSegment(i)->SetZ( ePlate.GetLayer(1)->Z() + ePlate.Z() );
844 for(int i=0; i<eS2cnd.N(); i++) eS2cnd.GetSegment(i)->SetZ( ePlate.GetLayer(2)->Z() + ePlate.Z() );
845 for(int i=0; i<eSpre.N(); i++) eSpre.GetSegment(i)->SetZ( ePlate.Z() );
846 for(int i=0; i<eS1pre.N(); i++) eS1pre.GetSegment(i)->SetZ( ePlate.GetLayer(1)->Z() + ePlate.Z() );
847 for(int i=0; i<eS2pre.N(); i++) eS2pre.GetSegment(i)->SetZ( ePlate.GetLayer(2)->Z() + ePlate.Z() );
848
849 eNext.SetDZ( ePlate.GetLayer(0)->DZ() );
850 eS.SetDZ( ePlate.GetLayer(0)->DZ() );
851 eS1.SetDZ( ePlate.GetLayer(1)->DZ());
852 eS2.SetDZ( ePlate.GetLayer(2)->DZ());
853 for(int i=0; i<eScnd.N(); i++) eScnd.GetSegment(i)->SetDZ( ePlate.GetLayer(0)->DZ() );
854 for(int i=0; i<eS1cnd.N(); i++) eS1cnd.GetSegment(i)->SetDZ( ePlate.GetLayer(1)->DZ() );
855 for(int i=0; i<eS2cnd.N(); i++) eS2cnd.GetSegment(i)->SetDZ( ePlate.GetLayer(2)->DZ() );
856 for(int i=0; i<eSpre.N(); i++) eSpre.GetSegment(i)->SetDZ( ePlate.GetLayer(0)->DZ() );
857 for(int i=0; i<eS1pre.N(); i++) eS1pre.GetSegment(i)->SetDZ( ePlate.GetLayer(1)->DZ() );
858 for(int i=0; i<eS2pre.N(); i++) eS2pre.GetSegment(i)->SetDZ( ePlate.GetLayer(2)->DZ() );
859}
int ID() const
Definition: EdbLayer.h:73
float DZ() const
Definition: EdbLayer.h:84
EdbAffine2D * GetAffineXY()
Definition: EdbLayer.h:119
float Z() const
Definition: EdbLayer.h:77
EdbLayer * GetLayer(int i)
Definition: EdbBrick.h:29
virtual void Transform(const EdbAffine2D *a)
Definition: EdbVirtual.cxx:155
void SetPID(int pid)
Definition: EdbSegP.h:129
void SetDZ(float dz)
Definition: EdbSegP.h:126

◆ UpdateFlag()

int EdbPlateTracking::UpdateFlag ( int  flag,
int  status 
)

status: -1 -found nothing, 0-bt, 1-mt1, 2-mt2

90{
92
93 int bth = flag/10000;
94 int mth = (flag/100)%100;
95 int tb = flag%10;
96
97 switch (status) {
98 case -1: bth++; mth++; break;
99 case 0: bth=0; mth=0; tb=0; break;
100 case 1: bth++; mth=0; tb=1; break;
101 case 2: bth++; mth=0; tb=2; break;
102 }
103
104 return( bth*10000+mth*100+tb );
105}

◆ UpdateSBtree()

bool EdbPlateTracking::UpdateSBtree ( TTree &  tsbt,
int  idp[4],
int  idf[4] 
)
679{
680 EdbSegP *s_pred = &ePred, *s_bt = &eS, *s_mt1=&eS1, *s_mt2=&eS2, *s_next=&eNext;
681
682 tsbt.SetBranchAddress("idpred",idp);
683 tsbt.SetBranchAddress("idfound",idf);
684 tsbt.SetBranchAddress("stat",&eStatus);
685 tsbt.SetBranchAddress("pred.", &s_pred);
686 tsbt.SetBranchAddress("s.", &s_bt);
687 tsbt.SetBranchAddress("s1.", &s_mt1);
688 tsbt.SetBranchAddress("s2.", &s_mt2);
689 tsbt.SetBranchAddress("next.", &s_next);
690
691 TClonesArray *s_cnd = eScnd.GetSegments();
692 tsbt.SetBranchAddress("scnd", &s_cnd);
693 TClonesArray *s1_cnd = eS1cnd.GetSegments();
694 tsbt.SetBranchAddress("s1cnd", &s1_cnd);
695 TClonesArray *s2_cnd = eS2cnd.GetSegments();
696 tsbt.SetBranchAddress("s2cnd", &s2_cnd);
697
698 TClonesArray *s_pre = eSpre.GetSegments();
699 tsbt.SetBranchAddress("spre", &s_pre);
700 TClonesArray *s1_pre = eS1pre.GetSegments();
701 tsbt.SetBranchAddress("s1pre", &s1_pre);
702 TClonesArray *s2_pre = eS2pre.GetSegments();
703 tsbt.SetBranchAddress("s2pre", &s2_pre);
704
705 tsbt.Fill();
706 return true;
707}

Member Data Documentation

◆ eChi2MaxBT

Float_t EdbPlateTracking::eChi2MaxBT

(1.5) maximum chi2 accepted between prediction and basetrack candidates

◆ eChi2MaxMT

Float_t EdbPlateTracking::eChi2MaxMT

(1.6) maximum chi2 accepted between prediction and microtrack candidates

◆ eCondBT

EdbScanCond EdbPlateTracking::eCondBT

◆ eCondMT

EdbScanCond EdbPlateTracking::eCondMT

------— processing parameters (can be default) and itermeadiate results ---------------------—

◆ eDegradPos

Float_t EdbPlateTracking::eDegradPos

SigmaX = SigmaX(0) + degradPos * mth.

◆ eDegradSlope

Float_t EdbPlateTracking::eDegradSlope

SigmaTX = SigmaTX(0) + degradSlope * bth.

◆ eDeltaR

Float_t EdbPlateTracking::eDeltaR

(20)

◆ eIdf

Int_t EdbPlateTracking::eIdf[4]

◆ eIdp

Int_t EdbPlateTracking::eIdp[4]

to read from sbt

◆ eNext

EdbSegP EdbPlateTracking::eNext

next prediction

◆ ePlate

EdbPlateP EdbPlateTracking::ePlate

plate geometry and correction parameters to be applied to prediction

◆ ePred

EdbSegP EdbPlateTracking::ePred

------— input data ------------------------—

prediction to be found in this plate in Brick reference system

◆ ePredictionScan

bool EdbPlateTracking::ePredictionScan

if true use GetPatternDataForPrediction( spred.ID(), side, pat ); in FindCandidates (default is false)

◆ ePreliminaryChi2MaxMT

Float_t EdbPlateTracking::ePreliminaryChi2MaxMT

(1.6) / microtracks and basetracks selection

◆ ePulsMinBT

Float_t EdbPlateTracking::ePulsMinBT

(18)

◆ ePulsMinDegradBT

Float_t EdbPlateTracking::ePulsMinDegradBT

(0)

◆ ePulsMinDegradMT

Float_t EdbPlateTracking::ePulsMinDegradMT

(0)

◆ ePulsMinMT

Float_t EdbPlateTracking::ePulsMinMT

(10) mimimal number of grains accepted to select microtracks

◆ eS

EdbSegP EdbPlateTracking::eS

------— output: main result ------------------------—

◆ eS1

EdbSegP EdbPlateTracking::eS1

◆ eS1cnd

EdbPattern EdbPlateTracking::eS1cnd

microtrack candidates passed all cuts

◆ eS1pre

EdbPattern EdbPlateTracking::eS1pre

the result of the selection of microtracks ordered by chi square

◆ eS2

EdbSegP EdbPlateTracking::eS2

found segments

◆ eS2cnd

EdbPattern EdbPlateTracking::eS2cnd

◆ eS2pre

EdbPattern EdbPlateTracking::eS2pre

◆ eScnd

EdbPattern EdbPlateTracking::eScnd

basetracks candidates passed all cuts

◆ eSide1

EdbPattern EdbPlateTracking::eSide1

side 1 microtracks in plate reference system

◆ eSide2

EdbPattern EdbPlateTracking::eSide2

side 2 microtracks in plate reference system

◆ eSpre

EdbPattern EdbPlateTracking::eSpre

the result of the selection of basetracks ordered by chi square

◆ eStatus

Int_t EdbPlateTracking::eStatus

-1-nothing, 0-bt, 1-mt1, 2-mt2


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