//
//////////////////////////////////////////////////////////////////////////////
-#include "TMatrixDSym.h"
-#include "TMatrixD.h"
-#include "TArrayI.h"
-#include "AliTrackFitterRieman.h"
+#include <TArrayI.h>
+#include <TLinearFitter.h>
+#include <TMath.h>
+#include <TMatrixD.h>
+#include <TMatrixDSym.h>
+#include <TRandom.h>
+#include <TTreeStream.h>
+
#include "AliLog.h"
-#include "TTreeStream.h"
-#include "AliRieman.h"
#include "AliLog.h"
+#include "AliRieman.h"
+#include "AliTrackFitterRieman.h"
ClassImp(AliTrackFitterRieman)
AliTrackFitterRieman::AliTrackFitterRieman():
AliTrackFitter(),
+ fBCorrection(kFALSE),
fAlpha(0.),
fNUsed(0),
fConv(kFALSE),
AliTrackFitterRieman::AliTrackFitterRieman(AliTrackPointArray *array, Bool_t owner):
AliTrackFitter(array,owner),
+ fBCorrection(kFALSE),
fAlpha(0.),
fNUsed(0),
fConv(kFALSE),
AliTrackFitterRieman::AliTrackFitterRieman(const AliTrackFitterRieman &rieman):
AliTrackFitter(rieman),
+ fBCorrection(rieman.fBCorrection),
fAlpha(rieman.fAlpha),
fNUsed(rieman.fNUsed),
fConv(rieman.fConv),
if(this==&rieman) return *this;
((AliTrackFitter *)this)->operator=(rieman);
+ fBCorrection = rieman.fBCorrection;
fAlpha = rieman.fAlpha;
fNUsed = rieman.fNUsed;
fConv = rieman.fConv;
// at the space-point reference plane. The reference plane is
// found using the covariance matrix of the point
// (assuming sigma(x)=0 at the reference coordinate system.
+ Int_t debugLevel = AliLog::GetDebugLevel("","AliTrackFitterRieman");
+ Float_t debugRatio = 1./(1.+debugLevel);
- Reset();
+ const Int_t kMinPoints =1;
Int_t npoints = fPoints->GetNPoints();
- if (fPoints && AliLog::GetDebugLevel("","AliTrackFitterRieman")>1){
+ if ( npoints<fMinNPoints) return kFALSE;
+ //
+ // fast count points
+ Int_t countFit = 0;
+ for (Int_t ifit=0; ifit<volIdsFit->GetSize(); ifit++){
+ Int_t volIdFit = volIdsFit->At(ifit);
+ for (Int_t ipoint = 0; ipoint < npoints; ipoint++)
+ if (volIdFit==fPoints->GetVolumeID()[ipoint]) countFit++;
+ }
+ if (countFit<fMinNPoints) return kFALSE;
+ //
+ Int_t countPoint = 0;
+ for (Int_t jpoint=0; jpoint<volIds->GetSize(); jpoint++){
+ Int_t volIdPoint = volIds->At(jpoint);
+ for (Int_t ipoint = 0; ipoint < npoints; ipoint++)
+ if (volIdPoint==fPoints->GetVolumeID()[ipoint]) countPoint++;
+ }
+ if (countPoint<kMinPoints) return kFALSE;
+ //
+ //
+
+ Reset();
+
+ if (fPoints && AliLog::GetDebugLevel("","AliTrackFitterRieman")>1&& gRandom->Rndm()<debugRatio){
Int_t nVol = volIds->GetSize();
Int_t nVolFit = volIdsFit->GetSize();
Int_t volId = volIds->At(0);
"fPoints.="<<fPoints<< // input points
"\n";
}
- if (npoints < 3) return kFALSE;
Bool_t isAlphaCalc = kFALSE;
AliTrackPoint p,plocal;
}
}
- if (AliLog::GetDebugLevel("","AliTrackFitterRieman")>0){
+ if (AliLog::GetDebugLevel("","AliTrackFitterRieman")>0 && gRandom->Rndm()<debugRatio){
//
// Debug Info
//
AliTrackPointArray *lPVolId = new AliTrackPointArray(npVolId);
AliTrackPointArray *lPTrack = new AliTrackPointArray(npVolId);
+ AliTrackPointArray *lPTrackE = new AliTrackPointArray(npVolId);
AliRieman * residual = fRieman->MakeResiduals();
- AliTrackPoint p2local;
for (Int_t ipoint = 0; ipoint < npVolId; ipoint++){
- Int_t index = pindex[ipoint];
- fPoints->GetPoint(p,index);
- if (GetPCA(p,p2)) {
- plocal = p.Rotate(fAlpha);
- // plocal.Rotate(fAlpha);
- Float_t xyz[3],xyz2[3];
- p2local = plocal;
- plocal.GetXYZ(xyz);
- xyz2[0] = xyz[0];
- xyz2[1] = GetYat(xyz[0]);
- xyz2[2] = GetZat(xyz[0]);
- p2local.SetXYZ(xyz2);
- lPVolId->AddPoint(ipoint,&plocal);
- lPTrack->AddPoint(ipoint,&p2local);
- }
+ AliTrackPoint p0, p0local;
+ AliTrackPoint pFit, pFitlocal, pFitLocalE;
+ fPVolId->GetPoint(p0,ipoint);
+ Float_t lAngle = p0.GetAngle();
+ p0local= p0.MasterToLocal();
+ fPTrack->GetPoint(pFit,ipoint);
+ pFitlocal= pFit.Rotate(lAngle);
+ //
+ Float_t xyz[3], cov[6];
+ xyz[0] = pFitlocal.GetX();
+ xyz[1] = pFitlocal.GetY();
+ xyz[2] = pFitlocal.GetZ();
+ for (Int_t icov=0; icov<6; icov++) cov[icov]=0;
+ cov[3] = GetErrY2at(xyz[0]);
+ cov[5] = GetErrZ2at(xyz[0]);
+ pFitLocalE.SetXYZ(xyz,cov);
+ //
+ lPVolId->AddPoint(ipoint,&p0local);
+ lPTrack->AddPoint(ipoint,&pFitlocal);
+ lPTrackE->AddPoint(ipoint,&pFitLocalE);
}
//
// debug info
Int_t volId = volIds->At(0);
Int_t modId =0;
Int_t layer = AliAlignObj::VolUIDToLayer(volId,modId);
+ Int_t volIdFit = volIdsFit->At(0);
+ Int_t modIdFit =0;
+ Int_t layerFit = AliAlignObj::VolUIDToLayer(volIdFit,modIdFit);
(*fDebugStream)<<"Fit"<<
"VolId="<<volId<< // volume ID
"Layer="<<layer<< // layer ID
"Module="<<modId<< // module ID
+ "LayerFit="<<layerFit<< // layer ID fit
+ "ModuleFit="<<modIdFit<< // module ID fit
"NVol="<<nVol<< // number of volumes
"NvolFit="<<nVolFit<< // number of volumes to fit
"Points0.="<<fPVolId<< // original points
"Points1.="<<fPTrack<< // fitted points
"LPoints0.="<<lPVolId<< // original points - local frame
"LPoints1.="<<lPTrack<< // fitted points - local frame
+ "LPointsE.="<<lPTrackE<< // fitted points with ext error - local frame
"Rieman.="<<this<< // original rieman fit
"Res.="<<residual<< // residuals of rieman fit
"\n";
fNUsed = fRieman->GetN();
fConv = kTRUE;
}
+ //
+ //
+ //
+ TLinearFitter fitY(3,"pol2");
+ TLinearFitter fitZ(3,"pol2");
+ for (Int_t ip=0; ip<fRieman->GetN();ip++){
+ Double_t x = fRieman->GetX()[ip];
+ fitY.AddPoint(&x,fRieman->GetY()[ip]-fRieman->GetYat(x),1);
+ fitZ.AddPoint(&x,fRieman->GetZ()[ip]-fRieman->GetZat(x),1);
+ }
+ fitY.Eval();
+ fitZ.Eval();
+ for (Int_t iparam=0; iparam<3; iparam++){
+ fCorrY[iparam]=fitY.GetParameter(iparam);
+ fCorrZ[iparam]=fitZ.GetParameter(iparam);
+ }
+ fCorrY[3]=fitY.GetChisquare()/Float_t(fRieman->GetN()-3);
+ fCorrZ[3]=fitZ.GetChisquare()/Float_t(fRieman->GetN()-3);
}
+
+
//_____________________________________________________________________________
Bool_t AliTrackFitterRieman::GetPCA(const AliTrackPoint &p, AliTrackPoint &p2) const
{
cov[3] = sigmax*dydx*dydx; cov[4] = sigmax*dydx*dzdx;
cov[5] = sigmax*dzdx*dzdx;
+ Double_t sigmay2 = GetErrY2at(x2);
+ Double_t sigmaz2 = GetErrZ2at(x2);
+ cov[3] += sigmay2;
+ cov[5] += sigmaz2;
+
+
Float_t newcov[6];
newcov[0] = cov[0]*cos*cos-
2*cov[1]*sin*cos+
newcov[5] = cov[5];
p2.SetXYZ(xsecond,ysecond,zsecond,newcov);
- if (AliLog::GetDebugLevel("","AliTrackFitterRieman")>1){
+ Int_t debugLevel = AliLog::GetDebugLevel("","AliTrackFitterRieman");
+ Float_t debugRatio = 1./(1.+debugLevel);
+ if (AliLog::GetDebugLevel("","AliTrackFitterRieman")>0 && gRandom->Rndm()<debugRatio){
AliTrackPoint lp0(p);
AliTrackPoint lp2(p2);
- if (0)(*fDebugStream)<<"PCA"<<
+ (*fDebugStream)<<"PCA"<<
"P0.="<<&lp0<<
"P2.="<<&lp2<<
"\n";
}
return kTRUE;
}
+
+Double_t AliTrackFitterRieman::GetYat(Double_t x) const {
+ //
+ // get y position at given point
+ //
+ Double_t correction=0;
+ if (fBCorrection){ // systematic effect correction
+ correction = fCorrY[0]+fCorrY[1]*x +fCorrY[2]*x*x;
+ }
+ return fRieman->GetYat(x)+correction;
+}
+
+Double_t AliTrackFitterRieman::GetZat(Double_t x) const {
+ //
+ // get z position at given point
+ //
+ Double_t correction=0;
+ if (fBCorrection){ // systematic effect correction
+ correction = fCorrZ[0]+fCorrZ[1]*x +fCorrZ[2]*x*x;
+ }
+ return fRieman->GetZat(x)+correction;
+}
+
+Double_t AliTrackFitterRieman::GetErrY2at(Double_t x) const {
+ //
+ // get estimate of extrapolation error
+ //
+ Double_t error = fRieman->GetErrY(x);
+ Double_t correction=0;
+ if (fBCorrection){ // everestimate error due systematic effect
+ error *=fCorrY[3];
+ correction = fCorrY[0]+fCorrY[1]*x +fCorrY[2]*x*x;
+ correction *=correction;
+ }
+ return TMath::Sqrt(error+correction);
+}
+
+Double_t AliTrackFitterRieman::GetErrZ2at(Double_t x) const {
+ //
+ // get estimate of extrapolation error
+ //
+ Double_t error = fRieman->GetErrZ(x)*fCorrZ[3];
+ Double_t correction=0;
+ if (fBCorrection){
+ error *= fCorrZ[3];
+ correction = fCorrZ[0]+fCorrZ[1]*x +fCorrZ[2]*x*x;
+ correction*= correction;
+ }
+ return TMath::Sqrt(error+correction);
+}
+