FEDRA emulsion software from the OPERA Collaboration
EdbPatCouple Class Reference

#include <EdbPVRec.h>

Inheritance diagram for EdbPatCouple:
Collaboration diagram for EdbPatCouple:

Public Member Functions

EdbSegCoupleAddSegCouple (int id1, int id2)
 
int Align (int alignFlag)
 
void CalculateAffXY (int alignFlag)
 
void CalculateAffXYZ (float z, int alignFlag)
 
int CheckSegmentsDuplication (EdbPattern *pat)
 
float Chi2A (EdbSegCouple *scp, int iprob=1)
 
float Chi2A (EdbSegP *s1, EdbSegP *s2, int iprob=1)
 
float Chi2KF (EdbSegCouple *scp)
 
int CHI2mode () const
 
float Chi2Pz0 (EdbSegCouple *scp)
 
void ClearSegCouples ()
 
EdbScanCondCond ()
 
int CutCHI2P (float chimax)
 
int DiffPat (EdbPattern *pat1, EdbPattern *pat2, Long_t vdiff[4])
 
int DiffPatCell (TIndexCell *cel1, TIndexCell *cel2, Long_t vdiff[4])
 
 EdbPatCouple ()
 
void FillCell_XYaXaY (EdbPattern *pat, EdbScanCond *cond, float dz, float stepx, float stepy)
 
void FillCell_XYaXaY (EdbScanCond *cond, float zlink, int id=0)
 
int FillCHI2 ()
 
int FillCHI2P ()
 
int FindOffset (EdbPattern *pat1, EdbPattern *pat2, Long_t vdiff[4])
 
int FindOffset0 (float xmax, float ymax)
 
int FindOffset01 (float xmax, float ymax)
 
int FindOffset1 (float xmax, float ymax)
 
EdbAffine2DGetAff ()
 
EdbSegCoupleGetSegCouple (int i) const
 
int ID1 () const
 
int ID2 () const
 
int LinkFast ()
 
int LinkSlow (float chi2max)
 
int Ncouples () const
 
float OffsetTX () const
 
float OffsetTY () const
 
float OffsetX () const
 
float OffsetY () const
 
EdbPatternPat1 ()
 
EdbPatternPat2 ()
 
void PrintCouples ()
 
void RemoveSegCouple (EdbSegCouple *sc)
 
int SelectIsolated ()
 
void SetCHI2mode (int m)
 
void SetCond (EdbScanCond *cond)
 
void SetID (int id1, int id2)
 
void SetOffset (float o1, float o2, float o3, float o4)
 
void SetOffsetsMax (float ox, float oy)
 
void SetPat1 (EdbPattern *pat1)
 
void SetPat2 (EdbPattern *pat2)
 
void SetSigma (float s1, float s2, float s3, float s4)
 
void SetZlink (float z)
 
float SigmaTX () const
 
float SigmaTY () const
 
float SigmaX () const
 
float SigmaY () const
 
int SortByCHI2P ()
 
float Zlink () const
 
 ~EdbPatCouple ()
 

Private Attributes

EdbAffine2DeAff
 affine transformation as: pat2 = pat1->Transform(aff) More...
 
float eChi2Max
 
int eCHI2mode
 algorithm used for chi2 calculation More...
 
EdbScanCondeCond
 
int eCoupleType
 1 - up/down link; 0 - plate-to-plate More...
 
Int_t eID [2]
 id-s of patterns in volume More...
 
Float_t eOffset [4]
 
EdbPatternePat1
 
EdbPatternePat2
 
TObjArray * eSegCouples
 
Float_t eSigma [4]
 expected pat<->pat sigma after alignment More...
 
float eXoffsetMax
 
float eYoffsetMax
 
float eZlink
 at this z patterns should be linked (aligned) More...
 

Constructor & Destructor Documentation

◆ EdbPatCouple()

EdbPatCouple::EdbPatCouple ( )

◆ ~EdbPatCouple()

EdbPatCouple::~EdbPatCouple ( )

111{
112 if(eSegCouples) {
113 eSegCouples->Delete();
114 SafeDelete(eSegCouples);
115 }
116}
TObjArray * eSegCouples
Definition: EdbPVRec.h:44

Member Function Documentation

◆ AddSegCouple()

EdbSegCouple * EdbPatCouple::AddSegCouple ( int  id1,
int  id2 
)

120{
121 if(!eSegCouples) eSegCouples = new TObjArray(1000);
122 EdbSegCouple *sc = new EdbSegCouple(id1,id2);
123 eSegCouples->Add(sc);
124 return sc;
125}
Definition: EdbSegCouple.h:17

◆ Align()

int EdbPatCouple::Align ( int  alignFlag)

536{
537 int npat =0;
538
539 SetZlink( (Pat1()->Z()+Pat2()->Z())/2. ); // link at central z
540
541 Log(2,"/nEdbPatCouple::Align","patterns: %d (%d) and %d (%d) at Zlink = %f",
542 Pat1()->ID(), Pat1()->N(), Pat2()->ID(),Pat2()->N(), Zlink() );
543
545
546 Long_t vdiff[4]={0,0,0,0};
547 int nitr=5;
548
549 FillCell_XYaXaY(Cond(), Zlink(), 2 ); Pat2()->Cell()->DropCouples(4);
550
551 for(int i=0; i<nitr; i++) {
552 FillCell_XYaXaY(Cond(), Zlink(), 1 ); Pat1()->Cell()->DropCouples(4);
553 npat = DiffPat(Pat1(),Pat2(),vdiff);
554
555 CalculateAffXYZ(Zlink(), alignFlag);
556
557 Pat1()->Transform(GetAff());
558 if( Pat1()->DiffAff(GetAff()) < Cond()->SigmaX(0)/20. ) break; // stop iterations
559 }
560
561 vdiff[0]=vdiff[1]=vdiff[2]=vdiff[3]=1;
562
563 npat = DiffPat( Pat1(), Pat2(), vdiff );
564 FillCHI2P();
565 SortByCHI2P();
566 CutCHI2P(1.5);
568 npat = Ncouples(); // the selected couples for the final transformation
569 CalculateAffXYZ(Zlink(),alignFlag);
570 Pat1()->Transform(GetAff());
571 Pat1()->SetNAff(npat);
572 EdbAffine2D affA;
573 affA.Set( GetAff()->A11(), GetAff()->A12(), GetAff()->A21(), GetAff()->A22(),0,0 );
574
575 //psel1->CalculateAXAY(psel2,affA);
576 Pat1()->TransformA(&affA);
577 GetAff()->Reset();
578
579
580 //psel2->CalculateAXAY(psel1,affA);
581 //pat2->TransformA(affA);
582
583 return npat;
584}
bool Log(int level, const char *location, const char *fmt,...)
Definition: EdbLog.cxx:75
Int_t npat
Definition: Xi2HatStartScript.C:33
Definition: EdbAffine.h:17
void Reset()
Definition: EdbAffine.cxx:72
void Set(EdbAffine2D &a)
Definition: EdbAffine.h:36
float eXoffsetMax
Definition: EdbPVRec.h:51
EdbPattern * Pat1()
Definition: EdbPVRec.h:89
EdbPattern * Pat2()
Definition: EdbPVRec.h:90
EdbAffine2D * GetAff()
Definition: EdbPVRec.h:76
int FillCHI2P()
Definition: EdbPVRec.cxx:266
EdbScanCond * Cond()
Definition: EdbPVRec.h:77
int SelectIsolated()
Definition: EdbPVRec.cxx:472
int FindOffset01(float xmax, float ymax)
Definition: EdbPVRec.cxx:191
int DiffPat(EdbPattern *pat1, EdbPattern *pat2, Long_t vdiff[4])
Definition: EdbPVRec.cxx:698
void SetZlink(float z)
Definition: EdbPVRec.h:69
float Zlink() const
Definition: EdbPVRec.h:79
float SigmaX() const
Definition: EdbPVRec.h:139
void FillCell_XYaXaY(EdbScanCond *cond, float zlink, int id=0)
Definition: EdbPVRec.cxx:807
int SortByCHI2P()
Definition: EdbPVRec.cxx:491
int CutCHI2P(float chimax)
Definition: EdbPVRec.cxx:452
void CalculateAffXYZ(float z, int alignFlag)
Definition: EdbPVRec.cxx:139
float eYoffsetMax
Definition: EdbPVRec.h:52
int Ncouples() const
Definition: EdbPVRec.h:80
void SetNAff(int n)
Definition: EdbPattern.h:311
TIndexCell * Cell() const
Definition: EdbPattern.h:321
virtual void Transform(const EdbAffine2D *a)
Definition: EdbVirtual.cxx:155
void TransformA(const EdbAffine2D *affA)
Definition: EdbPattern.cxx:324
int DropCouples(int level)
Definition: TIndexCell.cpp:206
Double_t Z
Definition: tlg2couples.C:104

◆ CalculateAffXY()

void EdbPatCouple::CalculateAffXY ( int  alignFlag)

calculate affine transformation for SELECTED couples at Z in the center

174{
176 float z = (Pat2()->Z()-Pat1()->Z())/2.;
177 CalculateAffXYZ(z,flag);
178 Pat1()->Transform(GetAff());
179 GetAff()->Reset();
180}
Float_t Z() const
Definition: EdbPattern.h:84

◆ CalculateAffXYZ()

void EdbPatCouple::CalculateAffXYZ ( float  z,
int  alignFlag 
)

calculate affine transformation for SELECTED couples at given z if flag==0 (default) - permit patterns deformation

140{
143
144 if(!eAff) eAff = new EdbAffine2D();
145 else eAff->Reset();
146
147 EdbSegCouple *scp=0;
148 EdbSegP *s1, *s2;
149 int ncp=Ncouples();
150
151 TArrayF x1(ncp);
152 TArrayF y1(ncp);
153 TArrayF x2(ncp);
154 TArrayF y2(ncp);
155
156 Float_t dz1,dz2;
157
158 for( int i=0; i<ncp; i++ ) {
159 scp=GetSegCouple(i);
160 s1 = Pat1()->GetSegment(scp->ID1());
161 s2 = Pat2()->GetSegment(scp->ID2());
162 dz1 = z - s1->Z();
163 dz2 = z - s2->Z();
164 x1[i] = s1->X() + dz1*s1->TX();
165 y1[i] = s1->Y() + dz1*s1->TY();
166 x2[i] = s2->X() + dz2*s2->TX();
167 y2[i] = s2->Y() + dz2*s2->TY();
168 }
169 eAff->Calculate( ncp,x1.fArray,y1.fArray,x2.fArray,y2.fArray, flag);
170}
Int_t Calculate(EdbPointsBox2D *b1, EdbPointsBox2D *b2)
Definition: EdbAffine.cxx:231
EdbSegCouple * GetSegCouple(int i) const
Definition: EdbPVRec.h:85
EdbAffine2D * eAff
affine transformation as: pat2 = pat1->Transform(aff)
Definition: EdbPVRec.h:42
int ID2() const
Definition: EdbSegCouple.h:56
int ID1() const
Definition: EdbSegCouple.h:55
Definition: EdbSegP.h:21
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 * GetSegment(int i) const
Definition: EdbPattern.h:66
EdbSegP * s1
Definition: tlg2couples.C:29
EdbSegP * s2
Definition: tlg2couples.C:30

◆ CheckSegmentsDuplication()

int EdbPatCouple::CheckSegmentsDuplication ( EdbPattern pat)

621{
622 TIndexCell *cell=pat->Cell();
623 int n1 = cell->GetEntries();
624 TIndexCell *c1,*c2,*c3,*c4;
625 EdbSegP *s1=0,*s2=0;
626 for( int i1=0; i1<n1; i1++) {
627 c1 = cell->At(i1);
628 int n2 = c1->GetEntries();
629 for( int i2=0; i2<n2; i2++) {
630 c2 = c1->At(i2);
631 int n3 = c2->GetEntries();
632 for( int i3=0; i3<n3; i3++) {
633 c3 = c2->At(i3);
634 int n4 = c3->GetEntries();
635 for( int i4=0; i4<n4; i4++) {
636 c4 = c3->At(i4);
637 int N=c4->N();
638 if(N>1) {
639 for(int j1=0; j1<N-1; j1++) {
640 s1=pat->GetSegment(c4->At(j1)->Value());
641 for(int j2=j1+1; j2<N; j2++) {
642 s2=pat->GetSegment(c4->At(j2)->Value());
643
644 if(gDIFF) gDIFF->Fill( s1->X(),s1->Y(),s1->TX(),s1->TY(),s1->W(),
645 s2->X(),s2->Y(),s2->TX(),s2->TY(),s2->W(),
646 s1->Z(),
647 s1->Aid(0), s1->Aid(1), s2->Aid(0), s2->Aid(1)
648 );
649 if( s1->Flag()<0 ) continue;
650 if( s2->Flag()<0 ) continue;
651 if( s1->Aid(0) == s2->Aid(0)&& s1->Aid(1) == s2->Aid(1) && s1->Side()==s2->Side() ) continue;
652 if( TMath::Abs(s2->X()-s1->X()) > 3.) continue;
653 if( TMath::Abs(s2->Y()-s1->Y()) > 3.) continue;
654 s2->W()>s1->W() ? s1->SetFlag(-1):s2->SetFlag(-1);
655 }
656 }
657 }
658 }
659 }
660 }
661 }
662 return 0;
663}
TNtuple * gDIFF
Definition: EdbLog.cxx:26
Int_t Side() const
mandatory virtual functions:
Definition: EdbSegP.h:170
Float_t W() const
Definition: EdbSegP.h:151
Int_t Aid(int i) const
Definition: EdbSegP.h:169
void SetFlag(int flag)
Definition: EdbSegP.h:130
Int_t Flag() const
Definition: EdbSegP.h:149
sort collection with attributes
Definition: TIndexCell.h:19
Long_t Value() const
Definition: TIndexCell.h:79
Int_t GetEntries() const
Definition: TIndexCell.h:82
TIndexCell const * At(Int_t narg, Int_t vind[]) const
Definition: TIndexCell.cpp:519
Int_t N() const
Definition: TIndexCell.cpp:344
TCanvas * c1
Definition: energy.C:13
TCanvas * c2
Definition: energy.C:26

◆ Chi2A() [1/2]

float EdbPatCouple::Chi2A ( EdbSegCouple scp,
int  iprob = 1 
)

332{
333 EdbSegP *s1 = Pat1()->GetSegment(scp->ID1());
334 EdbSegP *s2 = Pat2()->GetSegment(scp->ID2());
335 float chi2 = Chi2A(s1,s2, iprob);
336
337 SafeDelete(scp->eS);
338 EdbSegP *s = scp->eS=new EdbSegP();
339 s->Set( 0, // id will be assigned on writing into the tree
340 (s1->X()+s2->X())/2.,
341 (s1->Y()+s2->Y())/2.,
342 (s1->X()-s2->X())/(s1->Z()-s2->Z()),
343 (s1->Y()-s2->Y())/(s1->Z()-s2->Z()),
344 s1->W()+s2->W(),0
345 );
346 s->SetZ( (s2->Z()+s1->Z())/2 );
347 s->SetChi2( chi2 );
348 s->SetDZ( s2->Z()-s1->Z() );
349 s->SetVolume( s1->Volume()+s2->Volume() );
350 s->SetMC(s1->MCEvt(),s1->MCTrack());
351
352 if(s1->EMULDigitArray())
353 for(int i=0; i<s1->EMULDigitArray()->GetEntries(); i++) s->addEMULDigit(s1->EMULDigitArray()->At(i));
354 if(s2->EMULDigitArray())
355 for(int i=0; i<s2->EMULDigitArray()->GetEntries(); i++) s->addEMULDigit(s2->EMULDigitArray()->At(i));
356
357 return chi2;
358}
cout<< tr-> GetEntries()<< endl
float Chi2A(EdbSegCouple *scp, int iprob=1)
Definition: EdbPVRec.cxx:331
EdbSegP * eS
the result of the fit
Definition: EdbSegCouple.h:27
Float_t Volume() const
Definition: EdbSegP.h:158
TRefArray * EMULDigitArray() const
Definition: EdbSegP.h:82
Int_t MCTrack() const
Definition: EdbSegP.h:146
Int_t MCEvt() const
Definition: EdbSegP.h:145
s
Definition: check_shower.C:55
Float_t chi2
Definition: testBGReduction_By_ANN.C:14

◆ Chi2A() [2/2]

float EdbPatCouple::Chi2A ( EdbSegP s1,
EdbSegP s2,
int  iprob = 1 
)

fast estimation of chi2 in the special case when the position errors of segments are negligible in respect to angular errors: sigmaXY/dz << sigmaTXY application: up/down linking, alignment (offset search)

All calculation are done in the track plane which remove the dependency of the polar angle (phi)

362{
370
371 TVector3 v1,v2,v;
372 v1.SetXYZ( s1->TX(), s1->TY() , -1. );
373 v2.SetXYZ( s2->TX(), s2->TY() , -1. );
374 v.SetXYZ( -( s2->X() - s1->X() ),
375 -( s2->Y() - s1->Y() ),
376 -( s2->Z() - s1->Z() ) );
377
378 float phi = v.Phi();
379 v.RotateZ( -phi );
380 v1.RotateZ( -phi );
381 v2.RotateZ( -phi );
382
383 float dz = v.Z();
384 float tx = v.X()/dz;
385 float ty = v.Y()/dz;
386 float stx = eCond->SigmaTX(tx);
387 float sty = eCond->SigmaTY(ty);
388
389 float prob1=1., prob2=1.;
390 if(iprob) {
391 prob1 = eCond->ProbSeg(tx,ty,s1->W());
392 prob2 = eCond->ProbSeg(tx,ty,s2->W());
393 }
394 float dtx1 = (v1.X()-tx)*(v1.X()-tx)/stx/stx/prob1;
395 float dty1 = (v1.Y()-ty)*(v1.Y()-ty)/sty/sty/prob1;
396 float dtx2 = (v2.X()-tx)*(v2.X()-tx)/stx/stx/prob2;
397 float dty2 = (v2.Y()-ty)*(v2.Y()-ty)/sty/sty/prob2;
398
399 float chi2a=TMath::Sqrt(dtx1+dty1+dtx2+dty2)/2.;
400 return chi2a;
401}
brick dz
Definition: RecDispMC.C:107
EdbScanCond * eCond
Definition: EdbPVRec.h:47
float SigmaTX(float ax) const
Definition: EdbScanCond.h:106
float SigmaTY(float ay) const
Definition: EdbScanCond.h:107
float ProbSeg(float tx, float ty, float puls) const
Definition: EdbScanCond.cxx:119

◆ Chi2KF()

float EdbPatCouple::Chi2KF ( EdbSegCouple scp)

exact estimation of chi2 with full fitting procedure Applications: up/down linking, alignment (offset search)

All calculation are done with full corellation matrix for both segments (dependency on the polar angle (phi) is taken into account)

305{
311
312 EdbSegP *s1 = Pat1()->GetSegment(scp->ID1());
313 EdbSegP *s2 = Pat2()->GetSegment(scp->ID2());
314 s1->SetErrors();
315 s2->SetErrors();
316 eCond->FillErrorsCov( s1->TX(), s1->TY(), s1->COV() );
317 eCond->FillErrorsCov( s2->TX(), s2->TY(), s2->COV() );
318
319 SafeDelete(scp->eS);
320 scp->eS=new EdbSegP(*s2);
321 float chi2 = EdbTrackFitter::Chi2Seg(scp->eS, s1);
322 scp->eS->PropagateToCOV( 0.5*(s1->Z()+s2->Z()) );
323
324 if(s1->EMULDigitArray())
325 for(int i=0; i<s1->EMULDigitArray()->GetEntries(); i++) scp->eS->addEMULDigit(s1->EMULDigitArray()->At(i));
326
327 return chi2;
328}
void FillErrorsCov(float tx, float ty, TMatrixD &cov)
Definition: EdbScanCond.cxx:161
void SetErrors()
Definition: EdbSegP.h:90
void addEMULDigit(TObject *a)
Definition: EdbSegP.h:77
TMatrixD & COV() const
Definition: EdbSegP.h:123
void PropagateToCOV(float z)
Definition: EdbSegP.cxx:253
static float Chi2Seg(EdbSegP *s1, EdbSegP *s2)
Definition: EdbTrackFitter.cxx:62

◆ CHI2mode()

int EdbPatCouple::CHI2mode ( ) const
inline
104{ return eCHI2mode; }
int eCHI2mode
algorithm used for chi2 calculation
Definition: EdbPVRec.h:56

◆ Chi2Pz0()

float EdbPatCouple::Chi2Pz0 ( EdbSegCouple scp)

405{
406 float sx=.5, sy=.5, sz=3., dz=44.;
407 float sa;
408 float a1,a2;
409 float tx,ty;
410
411 TVector3 v1,v2,v;
412
413 EdbSegP *s1 = Pat1()->GetSegment(scp->ID1());
414 EdbSegP *s2 = Pat2()->GetSegment(scp->ID2());
415
416 tx = (s1->TX()+s2->TX())/2.;
417 ty = (s1->TY()+s2->TY())/2.;
418
419 v1.SetXYZ( s1->TX()/sx, s1->TY()/sy , 1./sz );
420 v2.SetXYZ( s2->TX()/sx, s2->TY()/sy , 1./sz );
421
422 v = v1+v2;
423 v *= .5;
424
425 a1 = v.Angle(v1);
426 a2 = v.Angle(v2);
427
428 sa = sz/dz*TMath::Cos(v.Theta());
429
430 float chi2 = TMath::Sqrt( (a1*a1 + a2*a2)/2. ) /
431 eCond->ProbSeg( tx, ty, (s1->W()+s2->W())/2. ) / sa;
432
433 SafeDelete(scp->eS);
434 EdbSegP *s = scp->eS=new EdbSegP();
435 s->Set( 0, // id will be assigned on writing into the tree
436 (s1->X()+s2->X())/2.,
437 (s1->Y()+s2->Y())/2.,
438 tx,
439 ty,
440 s1->W()+s2->W(),0
441 );
442 s->SetZ( (s2->Z()+s1->Z())/2 );
443 s->SetDZ( (s2->DZ()+s1->DZ())/2. );
444 s->SetChi2( chi2 );
445 s->SetVolume( s1->Volume()+s2->Volume() );
446 s->SetMC(s1->MCEvt(),s1->MCTrack());
447
448 return chi2;
449}
Float_t DZ() const
Definition: EdbSegP.h:154
static Float_t Angle(const EdbSegP &s1, const EdbSegP &s2)
Definition: EdbSegP.cxx:444

◆ ClearSegCouples()

void EdbPatCouple::ClearSegCouples ( )
inline
82{ if(eSegCouples) eSegCouples->Clear(); }

◆ Cond()

EdbScanCond * EdbPatCouple::Cond ( )
inline
77{return eCond;}

◆ CutCHI2P()

int EdbPatCouple::CutCHI2P ( float  chi2max)

453{
454 TObjArray *sCouples = new TObjArray();
455 EdbSegCouple *sc = 0;
456 int ncp = Ncouples();
457
458 if(gEDBDEBUGLEVEL>1) printf("CutCHI2P (%4.1f): %d -> ", chi2max,ncp );
459
460 for( int i=ncp-1; i>=0; i-- ) {
461 sc = GetSegCouple(i);
462 if(sc->CHI2P()<=chi2max) sCouples->Add(sc);
463 else SafeDelete(sc);
464 }
465 SafeDelete(eSegCouples);
466 eSegCouples=sCouples;
467 if(gEDBDEBUGLEVEL>1) printf("%d\n", Ncouples() );
468 return Ncouples();
469}
float CHI2P() const
Definition: EdbSegCouple.h:62
gEDBDEBUGLEVEL
Definition: energy.C:7

◆ DiffPat()

int EdbPatCouple::DiffPat ( EdbPattern pat1,
EdbPattern pat2,
Long_t  vdiff[4] 
)
700{
701 if(!pat1) return 0;
702 if(!pat2) return 0;
703 TIndexCell *c1=pat1->Cell();
704 TIndexCell *c2=pat2->Cell();
705 int ncouples = 0;
706
707 ncouples = DiffPatCell(c1,c2,vdiff);
708
709 return ncouples;
710}
int DiffPatCell(TIndexCell *cel1, TIndexCell *cel2, Long_t vdiff[4])
Definition: EdbPVRec.cxx:713

◆ DiffPatCell()

int EdbPatCouple::DiffPatCell ( TIndexCell cel1,
TIndexCell cel2,
Long_t  vdiff[4] 
)

Input: 2 cells as "x:y:ax:ay:entry" where entry is local id vdiff - vector of differences

Output: 1 patterns as: "entry1:entry2"

Action: all entries from both patterns satisfied to current vdiff for each cell are selected. Entries could be taken more than once

tTODO: move this algorithm to TIndexCell (n-dim)

715{
723
725
726 int ncouples = 0;
727 int npat = 0;
729
730 TIndexCell *c1[4],*c2[4];
731
732 Long_t v[2]={0,0};
733 Long_t val0[6]={0,0,0,0,0,0};
734 Long_t val[5]={0,0,0,0,0};
735 int vind[5]={-1,-1,-1,-1,-1};
736
737 //EdbSegP *s1,*s2;
738
739 int ncel1, nc10, nc11, nc12, nc13, nc23;
740
741 ncel1 = cel1->GetEntries();
742 for( vind[0]=0; vind[0]<ncel1; vind[0]++) { //x1
743 c1[0] = cel1->At(vind[0]);
744 val0[0] = c1[0]->Value();
745 for( val[0] = val0[0]-vdiff[0]; val[0]<=val0[0]+vdiff[0]; val[0]++ ) { //x2
746 c2[0] = cel2->Find(val[0]);
747 if(!c2[0]) continue;
748
749 nc10 = c1[0]->GetEntries();
750 for( vind[1]=0; vind[1]<nc10; vind[1]++) { //y1
751 c1[1] = c1[0]->At(vind[1]);
752 val0[1] = c1[1]->Value();
753 for( val[1] = val0[1]-vdiff[1]; val[1]<=val0[1]+vdiff[1]; val[1]++ ) { //y2
754 c2[1] = c2[0]->Find(val[1]);
755 if(!c2[1]) continue;
756
757 nc11 = c1[1]->GetEntries();
758 for( vind[2]=0; vind[2]<nc11; vind[2]++) { //ax1
759 c1[2] = c1[1]->At(vind[2]);
760 val0[2] = c1[2]->Value();
761 for( val[2] = val0[2]-vdiff[2]; val[2]<=val0[2]+vdiff[2]; val[2]++ ) { //ax2
762 c2[2] = c2[1]->Find(val[2]);
763 if(!c2[2]) continue;
764
765 nc12 = c1[2]->GetEntries();
766 for( vind[3]=0; vind[3]<nc12; vind[3]++) { //ay1
767 c1[3] = c1[2]->At(vind[3]);
768 val0[3] = c1[3]->Value();
769 for( val[3] = val0[3]-vdiff[3]; val[3]<=val0[3]+vdiff[3]; val[3]++ ) { //ay2
770 c2[3] = c2[2]->Find(val[3]);
771 if(!c2[3]) continue;
772
773 npat++;
774
775 nc13 = c1[3]->GetEntries();
776 for(int ie1=0; ie1<nc13; ie1++) {
777 nc23 = c2[3]->GetEntries();
778 for(int ie2=0; ie2<nc23; ie2++) {
779
780 ncouples++;
781
782 v[0] = c1[3]->At(ie1)->Value();
783 v[1] = c2[3]->At(ie2)->Value();
784
785// s1 = Pat1()->GetSegment(v[0]);
786// s2 = Pat2()->GetSegment(v[1]);
787// if(IsCompatible(s1,s2,eCond) //TODO
788
789 AddSegCouple((int)v[0],(int)v[1]);
790
791 }
792 }
793
794 }
795 }
796 }
797 }
798 }
799 }
800 }
801 }
802
803 return npat;
804}
EdbSegCouple * AddSegCouple(int id1, int id2)
Definition: EdbPVRec.cxx:119
void ClearSegCouples()
Definition: EdbPVRec.h:82
TIndexCell * Find(Int_t narg, Long_t varg[]) const
Definition: TIndexCell.cpp:613

◆ FillCell_XYaXaY() [1/2]

void EdbPatCouple::FillCell_XYaXaY ( EdbPattern pat,
EdbScanCond cond,
float  dz,
float  stepx,
float  stepy 
)

fill cells according to Sigma(angle) functions here cells size do not depends on dz!

824{
827
828 TIndexCell *cell = pat->Cell();
829 if(cell) cell->Drop();
830 pat->SetStep(stepx,stepy,0,0);
831
832 float x,y,tx,ty;
833 Long_t val[5]; // x,y,ax,ay,i
834 EdbSegP *p;
835 int npat = pat->N();
836 Log(3,"EdbPatCouple::FillCell_XYaXaY","npat = %d",npat);
837 for(int i=0; i<npat; i++ ) {
838 p = pat->GetSegment(i);
839 if( p->Track() >-1 ) continue; //to check side effects!!!
840 tx = p->TX();
841 ty = p->TY();
842 x = p->X() + tx*dz;
843 y = p->Y() + ty*dz;
844 val[0]= (Long_t)(x / stepx );
845 val[1]= (Long_t)(y / stepy );
846 val[2]= (Long_t)(tx/ cond->StepTX(tx) );
847 val[3]= (Long_t)(ty/ cond->StepTY(ty) );
848 val[4]= (Long_t)(i);
849 cell->Add(5,val);
850 }
851 cell->Sort();
852 //cell->PrintStat();
853 //cell->SetName("x:y:tx:ty:entry");
854}
void SetStep(float stepx, float stepy, float steptx, float stepty)
Definition: EdbPattern.h:299
float StepTX(float tx) const
Definition: EdbScanCond.h:93
float StepTY(float ty) const
Definition: EdbScanCond.h:94
Int_t N() const
Definition: EdbPattern.h:86
void Drop()
Definition: TIndexCell.cpp:237
void Sort(Int_t upto=kMaxInt)
Definition: TIndexCell.cpp:539
Int_t Add(Int_t narg, Long_t varg[])
Definition: TIndexCell.cpp:602
p
Definition: testBGReduction_AllMethods.C:8

◆ FillCell_XYaXaY() [2/2]

void EdbPatCouple::FillCell_XYaXaY ( EdbScanCond cond,
float  zlink,
int  id = 0 
)

fill cells according to Sigma(angle) functions at z=zlink

808{
810
811 float dz1 = zlink - Pat1()->Z();
812 float dz2 = zlink - Pat2()->Z();
813 float dz = TMath::Max( TMath::Abs(dz1), TMath::Abs(dz2) );
814 float stepx = cond->StepX(dz);
815 float stepy = cond->StepY(dz);
816
817 if(id==0||id==1) FillCell_XYaXaY( Pat1(), cond, dz1, stepx, stepy );
818 if(id==0||id==2) FillCell_XYaXaY( Pat2(), cond, dz2, stepx, stepy );
819}
float StepY(float dz) const
Definition: EdbScanCond.cxx:97
float StepX(float dz) const
Definition: EdbScanCond.cxx:90

◆ FillCHI2()

int EdbPatCouple::FillCHI2 ( )

final chi2 calculation based on the linked track

289{
291
292 EdbSegCouple *scp=0;
293 float chi2;
294 int ncp = Ncouples();
295 for( int i=0; i<ncp; i++ ) {
296 scp=GetSegCouple(i);
297 chi2 = Chi2A(scp,0);
298 scp->SetCHI2(chi2);
299 }
300 return Ncouples();
301}
void SetCHI2(float chi2)
Definition: EdbSegCouple.h:47

◆ FillCHI2P()

int EdbPatCouple::FillCHI2P ( )

fast chi2 calculation used for couples selection

267{
269
270 EdbSegCouple *scp=0;
271 float chi2;
272 int ncp = Ncouples();
273
274 Log(3,"EdbPatCouple::FillCHI2P","eCHI2mode = %d\n",eCHI2mode);
275 for( int i=0; i<ncp; i++ ) {
276 scp=GetSegCouple(i);
277
278 if(eCHI2mode==2) chi2 = Chi2Pz0(scp);
279 else if(eCHI2mode==3) chi2 = Chi2KF(scp);
280 else chi2 = Chi2A(scp);
281
282 scp->SetCHI2P(chi2);
283 }
284 return Ncouples();
285}
float Chi2Pz0(EdbSegCouple *scp)
Definition: EdbPVRec.cxx:404
float Chi2KF(EdbSegCouple *scp)
Definition: EdbPVRec.cxx:304
void SetCHI2P(float chi2)
Definition: EdbSegCouple.h:48

◆ FindOffset()

int EdbPatCouple::FindOffset ( EdbPattern pat1,
EdbPattern pat2,
Long_t  vdiff[4] 
)
208{
209 float dz=pat2->Z()-pat1->Z();
210
211 float voff[2]={0,0};
212 float stepx = Cond()->StepX(dz/2);
213 float stepy = Cond()->StepY(dz/2);
214 int nx = (int)( eXoffsetMax/stepx );
215 int ny = (int)( eYoffsetMax/stepy );
216
217 if(nx==0&&ny==0) return 2;
218 Log(2,"EdbPatCouple::FindOffset","( %d x %d ) with steps %8.3f %8.3f \t%f %f",
219 2*nx+1,2*ny+1,stepx,stepy,Cond()->StepTX(0),Cond()->StepTY(0) );
220
221 FillCell_XYaXaY(pat1,Cond(),dz/2.,stepx,stepy);
222 Log(3,"EdbPatCouple::FindOffset"," drop couples: %d",pat1->Cell()->DropCouples(4) );
223 FillCell_XYaXaY(pat2,Cond(),-dz/2.,stepx,stepy);
224 Log(3,"EdbPatCouple::FindOffset"," drop couples: %d",pat2->Cell()->DropCouples(4) );
225
226 Long_t vshift[4] = {0,0,0,0};
227 TIndexCell *c1=pat1->Cell();
228 TIndexCell *c2=pat2->Cell();
229
230 int npat0 = 0;
231 int npat = 0;
232
233 for(int iy=-ny; iy<=ny; iy++ ) {
234 for(int ix=-nx; ix<=nx; ix++ ) {
235 vshift[0] = ix;
236 vshift[1] = iy;
237 c1->Shift(2,vshift);
238
239 npat = c1->ComparePatterns(4,vdiff,c2);
240 if(gEDBDEBUGLEVEL>1) printf("%5d",npat);
241
242 if(npat>npat0) {
243 npat0 = npat;
244 voff[0] = ix*stepx;
245 voff[1] = iy*stepy;
246 }
247
248 vshift[0] = -ix;
249 vshift[1] = -iy;
250 c1->Shift(2,vshift);
251 }
252 if(gEDBDEBUGLEVEL>1) printf("\n");
253 }
254
255 Log(2,"EdbPatCouple::FindOffset","Offset: %f %f npat = %d", voff[0], voff[1], npat0 );
256
257 EdbAffine2D aff;
258 aff.ShiftX(voff[0]);
259 aff.ShiftY(voff[1]);
260 pat1->Transform(&aff);
261
262 return npat0;
263}
void ShiftX(float d)
Definition: EdbAffine.h:64
void ShiftY(float d)
Definition: EdbAffine.h:65

◆ FindOffset0()

int EdbPatCouple::FindOffset0 ( float  xmax,
float  ymax 
)
184{
185 Long_t vdiff[4]={0,0,0,0};
186 SetOffsetsMax(xmax,ymax);
187 return FindOffset( Pat1(),Pat2(),vdiff);
188}
void SetOffsetsMax(float ox, float oy)
Definition: EdbPVRec.h:63
int FindOffset(EdbPattern *pat1, EdbPattern *pat2, Long_t vdiff[4])
Definition: EdbPVRec.cxx:207

◆ FindOffset01()

int EdbPatCouple::FindOffset01 ( float  xmax,
float  ymax 
)
192{
193 Long_t vdiff[4]={0,0,1,1};
194 SetOffsetsMax(xmax,ymax);
195 return FindOffset( Pat1(),Pat2(),vdiff);
196}

◆ FindOffset1()

int EdbPatCouple::FindOffset1 ( float  xmax,
float  ymax 
)
200{
201 Long_t vdiff[4]={1,1,1,1};
202 SetOffsetsMax(xmax,ymax);
203 return FindOffset( Pat1(),Pat2(),vdiff);
204}

◆ GetAff()

EdbAffine2D * EdbPatCouple::GetAff ( )
inline
76{return eAff;}

◆ GetSegCouple()

EdbSegCouple * EdbPatCouple::GetSegCouple ( int  i) const
inline
86 { return (EdbSegCouple *)(eSegCouples->At(i)); }

◆ ID1()

int EdbPatCouple::ID1 ( ) const
inline
131{ return eID[0]; }
Int_t eID[2]
id-s of patterns in volume
Definition: EdbPVRec.h:31

◆ ID2()

int EdbPatCouple::ID2 ( ) const
inline
132{ return eID[1]; }

◆ LinkFast()

int EdbPatCouple::LinkFast ( )

to link segments of already aligned patterns fast - because check couples by cells

667{
670
671 if(Pat1()->N()<1) return 0;
672 if(Pat2()->N()<1) return 0;
673 int npat =0;
674 SetZlink( (Pat2()->Z()+Pat1()->Z())/2. );
675
676 Log(2,"EdbPatCouple::LinkFast"," %3d (%7d) and %3d (%7d) at Z = %13.3f",
677 Pat1()->ID(), Pat1()->N(), Pat2()->ID(),Pat2()->N(), Zlink() );
678
680 Log(3,"EdbPatCouple::LinkFast","1");
681
682 //CheckSegmentsDuplication(Pat1());
683 //CheckSegmentsDuplication(Pat2());
684
685 Long_t vdiff[4]={1,1,1,1};
686
687 npat = DiffPat( Pat1(), Pat2(), vdiff );
688 Log(3,"EdbPatCouple::LinkFast","2");
689
690 FillCHI2P();
691 Log(3,"EdbPatCouple::LinkFast","End of linking: %3d (%7d) and %3d (%7d) at Z = %13.3f npat= %d",
692 Pat1()->ID(), Pat1()->N(), Pat2()->ID(),Pat2()->N(), Zlink(),npat );
693
694 return npat;
695}

◆ LinkSlow()

int EdbPatCouple::LinkSlow ( float  chi2max)

to link segments of already aligned patterns slow - because check all couples

588{
591
592 int npat =0;
593 EdbPattern const *pat1=Pat1(); if(!pat1) return npat;
594 EdbPattern const *pat2=Pat2(); if(!pat2) return npat;
595
597
598 Log(2,"EdbPatCouple::LinkSlow","patterns: %d (%d) and %d (%d) \n",
599 pat1->ID(), pat1->N(), pat2->ID(),pat2->N());
600
601 EdbSegCouple *c;
602 float chi2;
603
604 int n1,n2;
605 n1 = pat1->N();
606 for( int i1=0; i1<n1; i1++ ) {
607 n2 = pat2->N();
608 for( int i2=0; i2<n2; i2++ ) {
609 chi2 = Chi2A( pat1->GetSegment(i1), pat2->GetSegment(i2) );
610 if( chi2 > chi2max ) continue;
611 c=AddSegCouple(i1,i2);
612 c->SetCHI2P(chi2);
613 npat++;
614 }
615 }
616 return npat;
617}
Definition: EdbPattern.h:273
int ID() const
Definition: EdbPattern.h:319

◆ Ncouples()

int EdbPatCouple::Ncouples ( ) const
inline
80 { if(eSegCouples) return eSegCouples->GetEntriesFast();
81 else return 0; }

◆ OffsetTX()

float EdbPatCouple::OffsetTX ( ) const
inline
136{ return eOffset[2]; }
Float_t eOffset[4]
Definition: EdbPVRec.h:38

◆ OffsetTY()

float EdbPatCouple::OffsetTY ( ) const
inline
137{ return eOffset[3]; }

◆ OffsetX()

float EdbPatCouple::OffsetX ( ) const
inline
134{ return eOffset[0]; }

◆ OffsetY()

float EdbPatCouple::OffsetY ( ) const
inline
135{ return eOffset[1]; }

◆ Pat1()

EdbPattern * EdbPatCouple::Pat1 ( )
inline
89{ return ePat1; }
EdbPattern * ePat1
Definition: EdbPVRec.h:35

◆ Pat2()

EdbPattern * EdbPatCouple::Pat2 ( )
inline
90{ return ePat2; }
EdbPattern * ePat2
Definition: EdbPVRec.h:36

◆ PrintCouples()

void EdbPatCouple::PrintCouples ( )

129{
130 EdbSegCouple *sc=0;
131 int ncp = Ncouples();
132 for(int i=0; i<ncp; i++) {
133 sc = GetSegCouple(i);
134 sc->Print();
135 }
136}
void Print()
Definition: EdbSegCouple.cxx:68

◆ RemoveSegCouple()

void EdbPatCouple::RemoveSegCouple ( EdbSegCouple sc)
inline
87{ SafeDelete(sc); }

◆ SelectIsolated()

int EdbPatCouple::SelectIsolated ( )

473{
474 TObjArray *sCouples = new TObjArray();
475 EdbSegCouple *sc = 0;
476 int ncp = Ncouples();
477 if(gEDBDEBUGLEVEL>1) printf("SelectIsolated: %d -> ", ncp );
478
479 for( int i=ncp-1; i>=0; i-- ) {
480 sc = GetSegCouple(i);
481 if( sc->N1tot()>1 || sc->N2tot()>1 ) {SafeDelete(sc);}
482 else sCouples->Add(sc);
483 }
484 SafeDelete(eSegCouples);
485 eSegCouples=sCouples;
486 if(gEDBDEBUGLEVEL>1) printf(" %d \n", Ncouples() );
487 return Ncouples();
488}
int N2tot() const
Definition: EdbSegCouple.h:60
int N1tot() const
Definition: EdbSegCouple.h:59

◆ SetCHI2mode()

void EdbPatCouple::SetCHI2mode ( int  m)
inline
103{ eCHI2mode=m; }

◆ SetCond()

void EdbPatCouple::SetCond ( EdbScanCond cond)
inline
74{ eCond=cond; }

◆ SetID()

void EdbPatCouple::SetID ( int  id1,
int  id2 
)
inline
62{ eID[0]=id1; eID[1]=id2; }

◆ SetOffset()

void EdbPatCouple::SetOffset ( float  o1,
float  o2,
float  o3,
float  o4 
)
inline
66 { eOffset[0]=o1; eOffset[1]=o2; eOffset[2]=o3; eOffset[3]=o4; }

◆ SetOffsetsMax()

void EdbPatCouple::SetOffsetsMax ( float  ox,
float  oy 
)
inline
64 { eXoffsetMax=ox; eYoffsetMax=oy; }

◆ SetPat1()

void EdbPatCouple::SetPat1 ( EdbPattern pat1)
inline
71{ ePat1=pat1; }

◆ SetPat2()

void EdbPatCouple::SetPat2 ( EdbPattern pat2)
inline
72{ ePat2=pat2; }

◆ SetSigma()

void EdbPatCouple::SetSigma ( float  s1,
float  s2,
float  s3,
float  s4 
)
inline
68 { eSigma[0]=s1; eSigma[1]=s2; eSigma[2]=s3; eSigma[3]=s4; }
Float_t eSigma[4]
expected pat<->pat sigma after alignment
Definition: EdbPVRec.h:40

◆ SetZlink()

void EdbPatCouple::SetZlink ( float  z)
inline
69{eZlink=z;}
float eZlink
at this z patterns should be linked (aligned)
Definition: EdbPVRec.h:49

◆ SigmaTX()

float EdbPatCouple::SigmaTX ( ) const
inline
141{ return eSigma[2]; }

◆ SigmaTY()

float EdbPatCouple::SigmaTY ( ) const
inline
142{ return eSigma[3]; }

◆ SigmaX()

float EdbPatCouple::SigmaX ( ) const
inline
139{ return eSigma[0]; }

◆ SigmaY()

float EdbPatCouple::SigmaY ( ) const
inline
140{ return eSigma[1]; }

◆ SortByCHI2P()

int EdbPatCouple::SortByCHI2P ( )

492{
493 int npat=0;
494
495 EdbSegCouple::SetSortFlag(0); // sort by CHI2P
496 eSegCouples->UnSort();
497 eSegCouples->Sort();
498
499 int np1 = ePat1->N();
500 int np2 = ePat2->N();
501 TArrayI found1(np1);
502 TArrayI found2(np2);
503 int i;
504 for(i=0; i<np1; i++) found1[i]=0;
505 for(i=0; i<np2; i++) found2[i]=0;
506
507 EdbSegCouple *sc = 0;
508
509 int ncp = Ncouples();
510 for( i=0; i<ncp; i++ ) {
511 sc = GetSegCouple(i);
512 found1[sc->ID1()]++;
513 found2[sc->ID2()]++;
514 sc->SetN1( found1[sc->ID1()] );
515 sc->SetN2( found2[sc->ID2()] );
516 npat++;
517 }
518
519 ncp = Ncouples();
520 for( i=0; i<ncp; i++ ) {
521 sc = GetSegCouple(i);
522 sc->SetN1tot( found1[sc->ID1()] );
523 sc->SetN2tot( found2[sc->ID2()] );
524 }
525
526 EdbSegCouple::SetSortFlag(1); // sort by (numbers + CHI2P)
527 eSegCouples->UnSort();
528 eSegCouples->Sort();
529
530 return npat;
531}
void SetN1(int n1)
Definition: EdbSegCouple.h:43
void SetN2(int n2)
Definition: EdbSegCouple.h:44
void SetN1tot(int n)
Definition: EdbSegCouple.h:45
void SetN2tot(int n)
Definition: EdbSegCouple.h:46
static void SetSortFlag(int s=0)
Definition: EdbSegCouple.cxx:62

◆ Zlink()

float EdbPatCouple::Zlink ( ) const
inline
79{return eZlink;}

Member Data Documentation

◆ eAff

EdbAffine2D* EdbPatCouple::eAff
private

affine transformation as: pat2 = pat1->Transform(aff)

◆ eChi2Max

float EdbPatCouple::eChi2Max
private

◆ eCHI2mode

int EdbPatCouple::eCHI2mode
private

algorithm used for chi2 calculation

◆ eCond

EdbScanCond* EdbPatCouple::eCond
private

◆ eCoupleType

int EdbPatCouple::eCoupleType
private

1 - up/down link; 0 - plate-to-plate

◆ eID

Int_t EdbPatCouple::eID[2]
private

id-s of patterns in volume

◆ eOffset

Float_t EdbPatCouple::eOffset[4]
private

maximal allowed offset before alignment "x:y:tx:ty"

◆ ePat1

EdbPattern* EdbPatCouple::ePat1
private

◆ ePat2

EdbPattern* EdbPatCouple::ePat2
private

◆ eSegCouples

TObjArray* EdbPatCouple::eSegCouples
private

◆ eSigma

Float_t EdbPatCouple::eSigma[4]
private

expected pat<->pat sigma after alignment

◆ eXoffsetMax

float EdbPatCouple::eXoffsetMax
private

◆ eYoffsetMax

float EdbPatCouple::eYoffsetMax
private

◆ eZlink

float EdbPatCouple::eZlink
private

at this z patterns should be linked (aligned)


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