Skip to content

Commit

Permalink
Tracklet fitting tuning from Alex: ref: https://indico.cern.ch/event/…
Browse files Browse the repository at this point in the history
  • Loading branch information
shahor02 committed May 23, 2014
1 parent d0e35a0 commit 6107975
Show file tree
Hide file tree
Showing 5 changed files with 145 additions and 33 deletions.
55 changes: 55 additions & 0 deletions TRD/AliTRDrecoParam.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -63,6 +63,7 @@ AliTRDrecoParam::AliTRDrecoParam()
,fNumberOfConfigs(3)
,fFlags(0)
,fRawStreamVersion("DEFAULT")
,fdzdxXcrossFactor(0.)
,fMinMaxCutSigma(4.)
,fMinLeftRightCutSigma(8.)
,fClusMaxThresh(4.5)
Expand Down Expand Up @@ -100,6 +101,7 @@ AliTRDrecoParam::AliTRDrecoParam()
SetImproveTracklets();
SetLUT();
SetTailCancelation();
SetTrackletParams();
}

//______________________________________________________________
Expand Down Expand Up @@ -133,6 +135,7 @@ AliTRDrecoParam::AliTRDrecoParam(const AliTRDrecoParam &ref)
,fNumberOfConfigs(ref.fNumberOfConfigs)
,fFlags(ref.fFlags)
,fRawStreamVersion(ref.fRawStreamVersion)
,fdzdxXcrossFactor(ref.fdzdxXcrossFactor)
,fMinMaxCutSigma(ref.fMinMaxCutSigma)
,fMinLeftRightCutSigma(ref.fMinLeftRightCutSigma)
,fClusMaxThresh(ref.fClusMaxThresh)
Expand All @@ -149,6 +152,12 @@ AliTRDrecoParam::AliTRDrecoParam(const AliTRDrecoParam &ref)
memcpy(fTCParams, ref.fTCParams, 8*sizeof(Double_t));
memcpy(fPIDThreshold, ref.fPIDThreshold, AliTRDCalPID::kNMom*sizeof(Double_t));
memcpy(fStreamLevel, ref.fStreamLevel, kTRDreconstructionTasks * sizeof(Int_t));

// tracklet params
memcpy(fdzdxCorrFactor, ref.fdzdxCorrFactor, 2*sizeof(Double_t));
memcpy(fdzdxCorrRCbias, ref.fdzdxCorrRCbias, 2*sizeof(Double_t));
memcpy(fYcorrTailCancel, ref.fdzdxCorrRCbias, 6*sizeof(Double_t));
memcpy(fS2Ycorr, ref.fS2Ycorr, 2*sizeof(Double_t));
}

//______________________________________________________________
Expand Down Expand Up @@ -188,6 +197,7 @@ AliTRDrecoParam& AliTRDrecoParam::operator=(const AliTRDrecoParam &ref)
fNumberOfConfigs = ref.fNumberOfConfigs;
fFlags = ref.fFlags;
fRawStreamVersion = ref.fRawStreamVersion;
fdzdxXcrossFactor = ref.fdzdxXcrossFactor;
fMinMaxCutSigma = ref.fMinMaxCutSigma;
fMinLeftRightCutSigma = ref.fMinLeftRightCutSigma;
fClusMaxThresh = ref.fClusMaxThresh;
Expand All @@ -201,6 +211,12 @@ AliTRDrecoParam& AliTRDrecoParam::operator=(const AliTRDrecoParam &ref)
memcpy(fTCParams, ref.fTCParams, 8*sizeof(Double_t));
memcpy(fPIDThreshold, ref.fPIDThreshold, AliTRDCalPID::kNMom*sizeof(Double_t));
memcpy(fStreamLevel, ref.fStreamLevel, kTRDreconstructionTasks * sizeof(Int_t));

// tracklet params
memcpy(fdzdxCorrFactor, ref.fdzdxCorrFactor, 2*sizeof(Double_t));
memcpy(fdzdxCorrRCbias, ref.fdzdxCorrRCbias, 2*sizeof(Double_t));
memcpy(fYcorrTailCancel, ref.fdzdxCorrRCbias, 6*sizeof(Double_t));
memcpy(fS2Ycorr, ref.fS2Ycorr, 2*sizeof(Double_t));
return *this;
}

Expand Down Expand Up @@ -322,3 +338,42 @@ void AliTRDrecoParam::SetPIDLQslices(Int_t s)
}
}

//___________________________________________________
void AliTRDrecoParam::SetTrackletParams(Double_t *par)
{
// Load tracklet reconstruction parameters. If none are set use defaults
if(par){
// correct dzdx for the bias in z
fdzdxCorrFactor[0] = par[0]; // !RC
fdzdxCorrFactor[1] = par[1]; // RC
// correct dzdx in RC tracklets for the bias in cluster attachment
fdzdxCorrRCbias[0] = par[2]; // dz/dx > 0
fdzdxCorrRCbias[1] = par[3]; // dz/dx < 0
/// correct x_cross for the bias in dzdx
fdzdxXcrossFactor = par[4];
// y linear q/pt correction due to wrong tail cancellation.
fYcorrTailCancel[0][0] = par[5]; fYcorrTailCancel[0][1] = par[6]; // opposite sign !RC
fYcorrTailCancel[1][0] = par[7]; fYcorrTailCancel[1][1] = par[8]; // same sign !RC
fYcorrTailCancel[2][0] = par[9]; fYcorrTailCancel[2][1] = par[10]; // RC
// inflation factor of error parameterization in r-phi due to wrong estimation of residuals.
fS2Ycorr[0] = par[11]; // opposite sign,
fS2Ycorr[1] = par[12]; // same sign

} else {
// correct dzdx for the bias in z
fdzdxCorrFactor[0] = 1.09; // !RC
fdzdxCorrFactor[1] = 1.05; // RC
// correct dzdx in RC tracklets for the bias in cluster attachment
fdzdxCorrRCbias[0] = 0.; // dz/dx > 0
fdzdxCorrRCbias[1] = -0.012; // dz/dx < 0
/// correct x_cross for the bias in dzdx
fdzdxXcrossFactor = 0.14;
// y linear q/pt correction due to wrong tail cancellation.
fYcorrTailCancel[0][0] = 0.; fYcorrTailCancel[0][1] = 0.027; // opposite sign !RC
fYcorrTailCancel[1][0] = 0.04; fYcorrTailCancel[1][1] = 0.027; // same sign !RC
fYcorrTailCancel[2][0] = 0.013;fYcorrTailCancel[2][1] = 0.018; // RC
// inflation factor of error parameterization in r-phi due to wrong estimation of residuals.
fS2Ycorr[0] = 5.; // opposite sign,
fS2Ycorr[1] = 10; // same sign
}
}
45 changes: 34 additions & 11 deletions TRD/AliTRDrecoParam.h
Original file line number Diff line number Diff line change
Expand Up @@ -63,8 +63,8 @@ class AliTRDrecoParam : public AliDetectorRecoParam
Double_t GetNMeanClusters() const { return fkNMeanClusters; }
Double_t GetNSigmaClusters() const { return fkNSigmaClusters; }
Double_t GetFindableClusters() const { return fkFindable; }
inline Int_t GetPIDLQslices() const;
inline AliTRDPIDResponse::ETRDPIDMethod GetPIDmethod() const;
Int_t GetPIDLQslices() const;
AliTRDPIDResponse::ETRDPIDMethod GetPIDmethod() const;
Double_t GetMaxTheta() const { return fkMaxTheta; }
Double_t GetMaxPhi() const { return fkMaxPhi; }
Double_t GetPlaneQualityThreshold() const { return fkPlaneQualityThreshold; }
Expand All @@ -77,9 +77,9 @@ class AliTRDrecoParam : public AliDetectorRecoParam
Double_t GetRoad2z() const { return fkRoad2z; }
Double_t GetRoadzMultiplicator() const { return fkRoadzMultiplicator; }
Double_t GetTrackLikelihood() const { return fkTrackLikelihood; }
inline void GetSysCovMatrix(Double_t *sys) const;
inline void GetTCParams(Double_t *par) const;
inline Int_t GetStreamLevel(ETRDReconstructionTask task) const;
void GetSysCovMatrix(Double_t *sys) const;
void GetTCParams(Double_t *par) const;
Int_t GetStreamLevel(ETRDReconstructionTask task) const;
const TString *GetRawStreamVersion() const{ return &fRawStreamVersion; };
Double_t GetMinMaxCutSigma() const { return fMinMaxCutSigma; };
Double_t GetMinLeftRightCutSigma() const { return fMinLeftRightCutSigma; };
Expand All @@ -90,6 +90,13 @@ class AliTRDrecoParam : public AliDetectorRecoParam
Int_t GetNumberOfPostsamples() const { return fNumberOfPostsamples;}
Int_t GetNumberOfSeedConfigs() const { return fNumberOfConfigs;}
Int_t GetRecEveryNTB() const { return fRecEveryNTB; }
// Tracklet parameters
Double_t GetCorrDZDXbiasRC(Bool_t dzdx) const { return fdzdxCorrRCbias[dzdx];}
Double_t GetCorrDZDX(Bool_t rc) const { return fdzdxCorrFactor[rc];}
Double_t GetCorrDZDXxcross() const { return fdzdxXcrossFactor;}
void GetYcorrTailCancel(Int_t it, Double_t par[2]) const;
Double_t GetS2Ycorr(Bool_t sgn) const { return fS2Ycorr[sgn];}

Bool_t IsArgon() const { return TESTBIT(fFlags, kDriftGas); }
Bool_t IsCheckTimeConsistency() const { return kCheckTimeConsistency;}
Bool_t IsOverPtThreshold(Double_t pt) const {return Bool_t(pt>fkPtThreshold);}
Expand Down Expand Up @@ -117,7 +124,7 @@ class AliTRDrecoParam : public AliDetectorRecoParam
void SetLUT(Bool_t b=kTRUE) {if(b) SETBIT(fFlags, kLUT); else CLRBIT(fFlags, kLUT);}
void SetGAUS(Bool_t b=kTRUE) {if(b) SETBIT(fFlags, kGAUS); else CLRBIT(fFlags, kGAUS);}
void SetPIDNeuralNetwork(Bool_t b=kTRUE) {if(b) SETBIT(fFlags, kSteerPID); else CLRBIT(fFlags, kSteerPID);}
inline void SetPIDmethod(AliTRDPIDResponse::ETRDPIDMethod method);
void SetPIDmethod(AliTRDPIDResponse::ETRDPIDMethod method);
void SetPIDLQslices(Int_t s);
void SetTailCancelation(Bool_t b=kTRUE) {if(b) SETBIT(fFlags, kTailCancelation); else CLRBIT(fFlags, kTailCancelation);}
void SetXenon(Bool_t b = kTRUE) {if(b) CLRBIT(fFlags, kDriftGas); else SETBIT(fFlags, kDriftGas);}
Expand All @@ -140,12 +147,13 @@ class AliTRDrecoParam : public AliDetectorRecoParam
void SetMinLeftRightCutSigma(Float_t minLeftRightCutSigma) { fMinLeftRightCutSigma = minLeftRightCutSigma; };
void SetClusMaxThresh(Float_t thresh) { fClusMaxThresh = thresh; };
void SetClusSigThresh(Float_t thresh) { fClusSigThresh = thresh; };
inline void SetPIDThreshold(Double_t *pid);
void SetPIDThreshold(Double_t *pid);
void SetPtThreshold(Double_t pt) {fkPtThreshold = pt;}
void SetNexponential(Int_t nexp) { fTCnexp = nexp; };
inline void SetTCParams(Double_t *par);
inline void SetStreamLevel(ETRDReconstructionTask task, Int_t level);
inline void SetSysCovMatrix(Double_t *sys);
void SetTCParams(Double_t *par);
void SetTrackletParams(Double_t *par=NULL);
void SetStreamLevel(ETRDReconstructionTask task, Int_t level);
void SetSysCovMatrix(Double_t *sys);
void SetNumberOfPresamples(Int_t n) { fNumberOfPresamples = n;}
void SetNumberOfPostsamples(Int_t n) { fNumberOfPostsamples = n;}
void SetRecEveryTwoTB() { fRecEveryNTB = 2; fkNMeanClusters = 10; }
Expand Down Expand Up @@ -192,6 +200,13 @@ class AliTRDrecoParam : public AliDetectorRecoParam
// Raw Reader Params
TString fRawStreamVersion; // Raw Reader version

// Tracklet parameters
Double_t fdzdxCorrFactor[2]; // correction of dzdx estimation due to z bias; [0] for !RC, [1] for RC
Double_t fdzdxCorrRCbias[2]; // correction of dzdx estimation bias for RC; [0] for dz/dx>0, [1] for dz/dx<0
Double_t fdzdxXcrossFactor; // bias in dzdx of estimated xcross [RC]
Double_t fYcorrTailCancel[3][2]; // y linear q/pt correction due to wrong tail cancellation. [0] opposite sign !RC, [1] same sign !RC, [2] RC
Double_t fS2Ycorr[2]; // inflation factor of error parameterization in r-phi due to wrong estimation of residuals. [0] opposite sign, [1] same sign

// Clusterization parameter
Double_t fMinMaxCutSigma; // Threshold sigma noise pad middle
Double_t fMinLeftRightCutSigma; // Threshold sigma noise sum pad
Expand All @@ -205,7 +220,7 @@ class AliTRDrecoParam : public AliDetectorRecoParam
Int_t fNumberOfPresamples; // number of presamples
Int_t fNumberOfPostsamples; // number of postsamples

ClassDef(AliTRDrecoParam, 12) // Reconstruction parameters for TRD detector
ClassDef(AliTRDrecoParam, 13) // Reconstruction parameters for TRD detector

};

Expand Down Expand Up @@ -257,6 +272,14 @@ inline void AliTRDrecoParam::SetTCParams(Double_t *par)
memcpy(fTCParams, par, 8*sizeof(Double_t));
}


//___________________________________________________
inline void AliTRDrecoParam::GetYcorrTailCancel(Int_t it, Double_t par[2]) const
{
if(it<0||it>2) return;
par[0] = fYcorrTailCancel[it][0]; par[1] = fYcorrTailCancel[it][1];
}

//___________________________________________________
inline Int_t AliTRDrecoParam::GetPIDLQslices() const
{
Expand Down
61 changes: 42 additions & 19 deletions TRD/AliTRDseedV1.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -633,19 +633,16 @@ Double_t AliTRDseedV1::EstimatedCrossPoint(AliTRDpadPlane *pp)
fZfit[1] = (fZfit[0]+zoff)/dx;

// correct dzdx for the bias
if(!IsRowCross()) fZfit[1] *= 1.09;
else{
// correct dzdx for the bias
fZfit[1] *= 1.05;
if(fZfit[1]<0) fZfit[1] -= 0.012;

UnbiasDZDX(IsRowCross());
if(IsRowCross()){
// correct x_cross/sigma(x_cross) for the bias in dzdx
fS2Y += 0.14*TMath::Abs(fZfit[1]);
sx += 2.e-2*GetS2DZDX(fZfit[1]);
const AliTRDrecoParam* const recoParam = fkReconstructor->GetRecoParam();
if(recoParam){
fS2Y += recoParam->GetCorrDZDXxcross()*TMath::Abs(fZfit[1]);
sx += recoParam->GetCorrDZDXxcross()*recoParam->GetCorrDZDXxcross()*GetS2DZDX(fZfit[1]);
}
// correct sigma(x_cross) for the width of the crossing area
Double_t sxdz = TMath::Abs(fZfit[1])>0.05?(TMath::Exp(-1.58839-TMath::Abs(fZfit[1])*3.24116)):
(0.957043 -TMath::Abs(fZfit[1])*12.4597);
sx += sxdz;
sx += GetS2XcrossDZDX(TMath::Abs(fZfit[1]));

// estimate z and error @ anode wire
fZfit[0] += fZfit[1]*fS2Y;
Expand All @@ -654,6 +651,36 @@ Double_t AliTRDseedV1::EstimatedCrossPoint(AliTRDpadPlane *pp)
return sx;
}

//____________________________________________________________________
void AliTRDseedV1::UnbiasDZDX(Bool_t rc)
{
// correct dzdx for the bias in z according to MC
const AliTRDrecoParam* const recoParam = fkReconstructor->GetRecoParam();
if(!recoParam) return;
fZfit[1] *= recoParam->GetCorrDZDX(rc);
if(rc) fZfit[1] += recoParam->GetCorrDZDXbiasRC(fZfit[1]<0);
}

//____________________________________________________________________
Double_t AliTRDseedV1::UnbiasY(Bool_t rc, Bool_t sgn, Int_t chg)
{
// correct y coordinate for tail cancellation. This should be fixed by considering TC as a function of q/pt.
// rc : TRUE if tracklet crosses rows
// sgn : TRUE if track has same sign with magnetic field
// chg : -1 for negative particles, +1 for the rest

const AliTRDrecoParam* const recoParam = fkReconstructor->GetRecoParam();
if(!recoParam) return 0.;
Double_t par[2]={0.};
if(rc) recoParam->GetYcorrTailCancel(2, par);
else{
if(sgn && 1./fPt > 1.5) recoParam->GetYcorrTailCancel(1, par);
else if(!sgn) recoParam->GetYcorrTailCancel(0, par);
}
return par[0]+par[1]*chg/fPt;
}


//____________________________________________________________________
Float_t AliTRDseedV1::GetQperTB(Int_t tb) const
{
Expand Down Expand Up @@ -2134,7 +2161,7 @@ Bool_t AliTRDseedV1::Fit(UChar_t opt)


//____________________________________________________________________
Bool_t AliTRDseedV1::FitRobust(AliTRDpadPlane *pp, Int_t opt)
Bool_t AliTRDseedV1::FitRobust(AliTRDpadPlane *pp, Bool_t sgn, Int_t chg, Int_t opt)
{
//
// Linear fit of the clusters attached to the tracklet
Expand Down Expand Up @@ -2188,12 +2215,8 @@ Bool_t AliTRDseedV1::FitRobust(AliTRDpadPlane *pp, Int_t opt)
Int_t n(0), // clusters used in fit
row[]={-1, -1},// pad row spanned by the tracklet
col(-1); // pad column of current cluster
Double_t ycorr(0.); Int_t chg(fYref[1]>0?1:-1);
if(IsRowCross()) ycorr = 0.013 + 0.018*chg/fPt;
else{
if(chg<0 && 1./fPt > 1.5) ycorr = 0.04 + 0.027*chg/fPt;
else if(chg>0) ycorr = 0.027*chg/fPt;
}
Double_t ycorr(UnbiasY(IsRowCross(), sgn, chg)),
kS2Ycorr(recoParam->GetS2Ycorr(sgn));

AliTRDcluster *c(NULL), **jc = &fClusters[0];
for(Int_t ic=0; ic<kNtb; ic++, ++jc) {
Expand Down Expand Up @@ -2306,7 +2329,7 @@ Bool_t AliTRDseedV1::FitRobust(AliTRDpadPlane *pp, Int_t opt)
fYfit[1] = -par[1];
//printf(" yfit: %f [%f] x[%e] dydx[%f]\n", fYfit[0], par[0], fX, par[1]);
// store covariance
fCov[0] = (chg>0?5.:10.)*cov[0]; // variance of y0
fCov[0] = kS2Ycorr*cov[0]; // variance of y0
fCov[1] = kScalePulls*cov[2]; // covariance of y0, dydx
fCov[2] = kScalePulls*cov[1]; // variance of dydx
// check radial position
Expand Down
13 changes: 12 additions & 1 deletion TRD/AliTRDseedV1.h
Original file line number Diff line number Diff line change
Expand Up @@ -91,7 +91,7 @@ class AliTRDseedV1 : public AliTRDtrackletBase
void CookLabels();
Bool_t CookPID();
Bool_t Fit(UChar_t opt=0); // OBSOLETE
Bool_t FitRobust(AliTRDpadPlane *pp, Int_t opt=0);
Bool_t FitRobust(AliTRDpadPlane *pp, Bool_t sgn, Int_t chg, Int_t opt=0);
Double_t EstimatedCrossPoint(AliTRDpadPlane *pp);
Bool_t Init(const AliTRDtrackV1 *track);
void Init(const AliRieman *fit);
Expand Down Expand Up @@ -157,6 +157,7 @@ class AliTRDseedV1 : public AliTRDtrackletBase
Float_t GetS2Z() const { return fS2Z;}
Double_t GetS2DYDX(Float_t) const { return fCov[2];}
inline Double_t GetS2DZDX(Float_t) const;
inline Double_t GetS2XcrossDZDX(Double_t absdzdx) const;
Float_t GetSigmaY() const { return fS2Y > 0. ? TMath::Sqrt(fS2Y) : 0.2;}
Float_t GetSnp() const { return fYref[1]/TMath::Sqrt(1+fYref[1]*fYref[1]);}
Float_t GetTgl() const { return fZref[1]/TMath::Sqrt(1+fYref[1]*fYref[1]);}
Expand Down Expand Up @@ -215,6 +216,8 @@ class AliTRDseedV1 : public AliTRDtrackletBase

protected:
void Copy(TObject &ref) const;
void UnbiasDZDX(Bool_t rc);
Double_t UnbiasY(Bool_t rc, Bool_t sgn, Int_t chg);

private:
inline void SetN(Int_t n);
Expand Down Expand Up @@ -298,6 +301,14 @@ inline Double_t AliTRDseedV1::GetPID(Int_t is) const
return 0.;
}

//____________________________________________________________
Double_t AliTRDseedV1::GetS2XcrossDZDX(Double_t absdzdx) const
{
// correct sigma(x_cross) for the width of the crossing area
if(absdzdx>0.05) return TMath::Exp(-1.58839-absdzdx*3.24116);
else return 0.957043-absdzdx*12.4597;
}

//____________________________________________________________
Double_t AliTRDseedV1::GetS2DZDX(Float_t dzdx) const
{
Expand Down
4 changes: 2 additions & 2 deletions TRD/AliTRDtrackerV1.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -875,6 +875,7 @@ Int_t AliTRDtrackerV1::FollowBackProlongation(AliTRDtrackV1 &t)
continue;
}

Float_t prod(t.GetBz()*t.Charge());
ptrTracklet = tracklets[ily];
if(!ptrTracklet){ // BUILD TRACKLET
AliDebug(3, Form("Building tracklet det[%d]", det));
Expand Down Expand Up @@ -920,7 +921,6 @@ Int_t AliTRDtrackerV1::FollowBackProlongation(AliTRDtrackV1 &t)
// Select attachment base on track to B field sign not only track charge which is buggy
// mark kFALSE same sign tracks and kTRUE opposite sign tracks
// A.Bercuci 3.11.2011
Float_t prod(t.GetBz()*t.Charge());
if(!ptrTracklet->AttachClusters(chamber, kTRUE, prod<0.?kTRUE:kFALSE, fEventInFile)){
t.SetErrStat(AliTRDtrackV1::kNoAttach, ily);
if(debugLevel>3){
Expand Down Expand Up @@ -956,7 +956,7 @@ Int_t AliTRDtrackerV1::FollowBackProlongation(AliTRDtrackV1 &t)
// tilt correction options
// 0 : no correction
// 2 : pseudo tilt correction
if(!ptrTracklet->FitRobust(fGeom->GetPadPlane(ily, stk))){
if(!ptrTracklet->FitRobust(fGeom->GetPadPlane(ily, stk), prod>0., t.Charge())){
t.SetErrStat(AliTRDtrackV1::kNoFit, ily);
AliDebug(4, "Failed Tracklet Fit");
continue;
Expand Down

0 comments on commit 6107975

Please sign in to comment.