X-Git-Url: http://git.uio.no/git/?p=u%2Fmrichter%2FAliRoot.git;a=blobdiff_plain;f=TPC%2FAliTPCcalibAlign.cxx;h=ef84bc9b14c63e1b987f5011b9d676a8773002a4;hp=85dcbd22e52d74fadc05df505e50fd2296666906;hb=9f02b62b3f5da006ee8f2d7e1de20f0754fd9765;hpb=3326b32373bdf653ab2b6f410337d9838b637798 diff --git a/TPC/AliTPCcalibAlign.cxx b/TPC/AliTPCcalibAlign.cxx index 85dcbd22e52..ef84bc9b14c 100644 --- a/TPC/AliTPCcalibAlign.cxx +++ b/TPC/AliTPCcalibAlign.cxx @@ -133,7 +133,7 @@ #include "TFile.h" #include "TProfile.h" #include "TCanvas.h" - +#include "TDatabasePDG.h" #include "TTreeStream.h" #include "Riostream.h" @@ -181,7 +181,6 @@ AliTPCcalibAlign::AliTPCcalibAlign() fMatrixArray9(72*72), fMatrixArray6(72*72), fCombinedMatrixArray6(72), - fCompTracklet(0), // tracklet comparison fNoField(kFALSE), fXIO(0), fXmiddle(0), @@ -207,12 +206,14 @@ AliTPCcalibAlign::AliTPCcalibAlign() fXquadrant = roc->GetPadRowRadii(36,53); fXmiddle = ( roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5; fXIO = ( roc->GetPadRowRadii(0,roc->GetNRows(0)-1)+roc->GetPadRowRadii(36,0))*0.5; - fClusterDelta[0]=0; // cluster residuals - fClusterDelta[1]=0; // cluster residuals - fClusterDelta[2]=0; // cluster residuals - vertex constrained - fClusterDelta[3]=0; // cluster residuals - fClusterDelta[4]=0; // cluster residuals - ITS constrained - fClusterDelta[5]=0; // cluster residuals + fClusterDelta[0]=0; // cluster residuals - Y + fClusterDelta[1]=0; // cluster residuals - Z + + + fTrackletDelta[0]=0; // tracklet residuals + fTrackletDelta[1]=0; // tracklet residuals + fTrackletDelta[2]=0; // tracklet residuals + fTrackletDelta[3]=0; // tracklet residuals } AliTPCcalibAlign::AliTPCcalibAlign(const Text_t *name, const Text_t *title) @@ -234,7 +235,6 @@ AliTPCcalibAlign::AliTPCcalibAlign(const Text_t *name, const Text_t *title) fMatrixArray9(72*72), fMatrixArray6(72*72), fCombinedMatrixArray6(72), - fCompTracklet(0), // tracklet comparison fNoField(kFALSE), fXIO(0), fXmiddle(0), @@ -265,13 +265,11 @@ AliTPCcalibAlign::AliTPCcalibAlign(const Text_t *name, const Text_t *title) fXIO = ( roc->GetPadRowRadii(0,roc->GetNRows(0)-1)+roc->GetPadRowRadii(36,0))*0.5; fClusterDelta[0]=0; // cluster residuals fClusterDelta[1]=0; // cluster residuals - fClusterDelta[2]=0; // cluster residuals - vertex constrained - fClusterDelta[3]=0; // cluster residuals - fClusterDelta[4]=0; // cluster residuals - ITS constrained - fClusterDelta[5]=0; // cluster residuals - - + fTrackletDelta[0]=0; // tracklet residuals + fTrackletDelta[1]=0; // tracklet residuals + fTrackletDelta[2]=0; // tracklet residuals + fTrackletDelta[3]=0; // tracklet residuals } @@ -296,7 +294,6 @@ AliTPCcalibAlign::AliTPCcalibAlign(const AliTPCcalibAlign &align) fMatrixArray9(align.fMatrixArray9), fMatrixArray6(align.fMatrixArray6), fCombinedMatrixArray6(align.fCombinedMatrixArray6), - fCompTracklet(align.fCompTracklet), // tracklet comparison fNoField(align.fNoField), fXIO(align.fXIO), fXmiddle(align.fXmiddle), @@ -378,11 +375,11 @@ AliTPCcalibAlign::AliTPCcalibAlign(const AliTPCcalibAlign &align) } fClusterDelta[0]=0; // cluster residuals fClusterDelta[1]=0; // cluster residuals - fClusterDelta[2]=0; // cluster residuals - vertex constrained - fClusterDelta[3]=0; // cluster residuals - fClusterDelta[4]=0; // cluster residuals - ITS constrained - fClusterDelta[5]=0; // cluster residuals + fTrackletDelta[0]=0; // tracklet residuals + fTrackletDelta[1]=0; // tracklet residuals + fTrackletDelta[2]=0; // tracklet residuals + fTrackletDelta[3]=0; // tracklet residuals } @@ -432,15 +429,20 @@ AliTPCcalibAlign::~AliTPCcalibAlign() { fMatrixArray9.Delete(); // array of transnformtation matrix fMatrixArray6.Delete(); // array of transnformtation matrix - if (fCompTracklet) delete fCompTracklet; fArraySectorIntParam.SetOwner(kTRUE); // array of sector alignment parameters fArraySectorIntCovar.SetOwner(kTRUE); // array of sector alignment covariances fArraySectorIntParam.Delete(); // array of sector alignment parameters fArraySectorIntCovar.Delete(); // array of sector alignment covariances - for (Int_t i=0; i<6; i++){ + for (Int_t i=0; i<2; i++){ delete fClusterDelta[i]; // cluster residuals } + + for (Int_t i=0; i<4; i++){ + delete fTrackletDelta[i]; // tracklet residuals + } + + } void AliTPCcalibAlign::Process(AliESDEvent *event) { @@ -448,6 +450,7 @@ void AliTPCcalibAlign::Process(AliESDEvent *event) { // Process pairs of cosmic tracks // if (!fClusterDelta[0]) MakeResidualHistos(); + if (!fTrackletDelta[0]) MakeResidualHistosTracklet(); // fCurrentEvent=event; ExportTrackPoints(event); // export track points for external calibration @@ -469,11 +472,13 @@ void AliTPCcalibAlign::Process(AliESDEvent *event) { AliTPCseed *seed0 = 0; // friendTrack = (AliESDfriendTrack *)eESDfriend->GetTrack(i0);; + if (!friendTrack) continue; for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) { if ((seed0=dynamic_cast(calibObject))) break; } if (!seed0) continue; fCurrentTrack=track0; + fCurrentFriendTrack=friendTrack; fCurrentSeed=seed0; fCurrentEvent=event; ProcessSeed(seed0); @@ -505,10 +510,12 @@ void AliTPCcalibAlign::Process(AliESDEvent *event) { AliTPCseed *seed0 = 0,*seed1=0; // friendTrack = (AliESDfriendTrack *)eESDfriend->GetTrack(i0);; + if (!friendTrack) continue; for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) { if ((seed0=dynamic_cast(calibObject))) break; } friendTrack = (AliESDfriendTrack *)eESDfriend->GetTrack(i1);; + if (!friendTrack) continue; for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) { if ((seed1=dynamic_cast(calibObject))) break; } @@ -664,11 +671,13 @@ void AliTPCcalibAlign::ExportTrackPoints(AliESDEvent *event){ AliTPCseed *seed0 = 0,*seed1=0; // friendTrack = (AliESDfriendTrack *)eESDfriend->GetTrack(index0);; + if (!friendTrack) continue; for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) { if ((seed0=dynamic_cast(calibObject))) break; } if (index1>0){ friendTrack = (AliESDfriendTrack *)eESDfriend->GetTrack(index1);; + if (!friendTrack) continue; for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) { if ((seed1=dynamic_cast(calibObject))) break; } @@ -778,6 +787,7 @@ void AliTPCcalibAlign::ProcessSeed(AliTPCseed *seed) { // // make a kalman tracklets out of seed // + UpdateClusterDeltaField(seed); TObjArray tracklets= AliTPCTracklet::CreateTracklets(seed,AliTPCTracklet::kKalman, kFALSE,20,4); @@ -793,7 +803,7 @@ void AliTPCcalibAlign::ProcessSeed(AliTPCseed *seed) { AliTPCTracklet *t2=static_cast(tracklets[i2]); AliExternalTrackParam *common1=0,*common2=0; if (AliTPCTracklet::PropagateToMeanX(*t1,*t2,common1,common2)){ - ProcessTracklets(*common1,*common2,seed, t1->GetSector(),t2->GetSector()); + ProcessTracklets(*common1,*common2,seed, t1->GetSector(),t2->GetSector()); UpdateAlignSector(seed,t1->GetSector()); } delete common1; @@ -901,7 +911,8 @@ void AliTPCcalibAlign::ProcessTracklets(const AliExternalTrackParam &tp1, // Int_t accept = AcceptTracklet(tp1,tp2); Int_t acceptLinear = AcceptTracklet(parLine1,parLine2); - + + if (fStreamLevel>1 && seed){ TTreeSRedirector *cstream = GetDebugStreamer(); if (cstream){ @@ -947,7 +958,7 @@ void AliTPCcalibAlign::ProcessTracklets(const AliExternalTrackParam &tp1, if (TMath::Abs(parLine1[2])<0.8 &&TMath::Abs(parLine1[2])<0.8 ){ //angular cut FillHisto(parLine1,parLine2,s1,s2); ProcessAlign(parLine1,parLine2,s1,s2); - UpdateKalman(s1,s2,par1, cov1, par2, cov2); + //UpdateKalman(s1,s2,par1, cov1, par2, cov2); - OBSOLETE to be removed - 50 % of time here } } } @@ -958,7 +969,11 @@ void AliTPCcalibAlign::ProcessTracklets(const AliExternalTrackParam &tp1, // // use Kalman if mag field // - if (seed) ProcessDiff(tp1,tp2, seed,s1,s2); + if (seed) { + ProcessDiff(tp1,tp2, seed,s1,s2); + FillHisto((AliExternalTrackParam*)&tp1,(AliExternalTrackParam*)&tp2,s1,s2); + FillHisto((AliExternalTrackParam*)&tp2,(AliExternalTrackParam*)&tp1,s2,s1); + } FillHisto(t1,t2,s1,s2); ProcessAlign(t1,t2,s1,s2); } @@ -970,8 +985,8 @@ void AliTPCcalibAlign::ProcessAlign(Double_t * t1, // // Do intersector alignment // - Process12(t1,t2,GetOrMakeFitter12(s1,s2)); - Process9(t1,t2,GetOrMakeFitter9(s1,s2)); + //Process12(t1,t2,GetOrMakeFitter12(s1,s2)); + //Process9(t1,t2,GetOrMakeFitter9(s1,s2)); Process6(t1,t2,GetOrMakeFitter6(s1,s2)); ++fPoints[GetIndex(s1,s2)]; } @@ -1540,7 +1555,7 @@ TLinearFitter* AliTPCcalibAlign::GetOrMakeFitter12(Int_t s1,Int_t s2) { fitter->StoreData(kFALSE); fFitterArray12.AddAt(fitter,GetIndex(s1,s2)); counter12++; - if (GetDebugLevel()>0) cerr<<"Creating fitter12 "<0) cerr<<"Creating fitter12 "<StoreData(kFALSE); fFitterArray9.AddAt(fitter,GetIndex(s1,s2)); counter9++; - if (GetDebugLevel()>0) cerr<<"Creating fitter12 "<0) cerr<<"Creating fitter12 "<StoreData(kFALSE); fFitterArray6.AddAt(fitter,GetIndex(s1,s2)); counter6++; - if (GetDebugLevel()>0) cerr<<"Creating fitter6 "<0) cerr<<"Creating fitter6 "<GetAxis(ivar2)->SetName(axisName[ivar2].Data()); fClusterDelta[ivar]->GetAxis(ivar2)->SetTitle(axisName[ivar2].Data()); } @@ -1688,6 +1697,67 @@ void AliTPCcalibAlign::MakeResidualHistos(){ } + +void AliTPCcalibAlign::MakeResidualHistosTracklet(){ + // + // Make tracklet residual histograms + // + Double_t xminTrack[9], xmaxTrack[9]; + Int_t binsTrack[9]; + TString axisName[9],axisTitle[9]; + // + // 0 - delta of interest + // 1 - global phi in sector number as float + // 2 - local x + // 3 - local ky + // 4 - local kz + // 5 - sector 1 + // 5 - sector 0 + + axisName[0]="delta"; axisTitle[0]="#Delta (cm)"; + binsTrack[0]=60; xminTrack[0]=-0.6; xmaxTrack[0]=0.6; + // + axisName[1]="phi"; axisTitle[1]="#phi"; + binsTrack[1]=180; xminTrack[1]=-TMath::Pi(); xmaxTrack[1]=TMath::Pi(); + // + axisName[2]="localX"; axisTitle[2]="x (cm)"; + binsTrack[2]=10; xminTrack[2]=120.; xmaxTrack[2]=200.; + // + axisName[3]="kY"; axisTitle[3]="dy/dx"; + binsTrack[3]=10; xminTrack[3]=-0.5; xmaxTrack[3]=0.5; + // + axisName[4]="kZ"; axisTitle[4]="dz/dx"; + binsTrack[4]=22; xminTrack[4]=-1.1; xmaxTrack[4]=1.1; + // + axisName[5]="is1"; axisTitle[5]="is1"; + binsTrack[5]=72; xminTrack[5]=0; xmaxTrack[5]=72; + // + axisName[6]="is0"; axisTitle[6]="is0"; + binsTrack[6]=72; xminTrack[6]=0; xmaxTrack[6]=72; + + // + xminTrack[0]=-0.3; xmaxTrack[0]=0.3; + fTrackletDelta[0] = new THnSparseF("#Delta_{Y} (cm)","#Delta_{Y} (cm)", 7, binsTrack,xminTrack, xmaxTrack); + xminTrack[0]=-0.5; xmaxTrack[0]=0.5; + fTrackletDelta[1] = new THnSparseF("#Delta_{Z} (cm)","#Delta_{Z} (cm)", 7, binsTrack,xminTrack, xmaxTrack); + xminTrack[0]=-0.005; xmaxTrack[0]=0.005; + fTrackletDelta[2] = new THnSparseF("#Delta_{kY}","#Delta_{kY}", 7, binsTrack,xminTrack, xmaxTrack); + xminTrack[0]=-0.005; xmaxTrack[0]=0.005; + fTrackletDelta[3] = new THnSparseF("#Delta_{kZ}","#Delta_{kZ}", 7, binsTrack,xminTrack, xmaxTrack); + // + // + // + for (Int_t ivar=0;ivar<4;ivar++){ + for (Int_t ivar2=0;ivar2<7;ivar2++){ + fTrackletDelta[ivar]->GetAxis(ivar2)->SetName(axisName[ivar2].Data()); + fTrackletDelta[ivar]->GetAxis(ivar2)->SetTitle(axisName[ivar2].Data()); + } + } + +} + + + void AliTPCcalibAlign::FillHisto(const Double_t *t1, const Double_t *t2, Int_t s1,Int_t s2) { @@ -1720,6 +1790,55 @@ void AliTPCcalibAlign::FillHisto(const Double_t *t1, } +void AliTPCcalibAlign::FillHisto(AliExternalTrackParam *tp1, + AliExternalTrackParam *tp2, + Int_t s1,Int_t s2) { + // + // Fill residual histograms + // Track2-Track1 + const Double_t kEpsilon=0.001; + Double_t x[8]={0,0,0,0,0,0,0,0}; + AliExternalTrackParam p1(*tp1); + AliExternalTrackParam p2(*tp2); + if (s1%18==s2%18) { + // inner outer - match at the IROC-OROC boundary + if (!p1.PropagateTo(fXIO, AliTrackerBase::GetBz())) return; + } + if (!p2.Rotate(p1.GetAlpha())) return; + if (!p2.PropagateTo(p1.GetX(),AliTrackerBase::GetBz())) return; + if (TMath::Abs(p1.GetX()-p2.GetX())>kEpsilon) return; + Double_t xyz[3]; + p1.GetXYZ(xyz); + x[1]=TMath::ATan2(xyz[1],xyz[0]); + x[2]=p1.GetX(); + x[3]=0.5*(p1.GetSnp()+p2.GetSnp()); // mean snp + x[4]=0.5*(p1.GetTgl()+p2.GetTgl()); // mean tgl + x[5]=s2; + x[6]=s1; + + x[0]=p2.GetY()-p1.GetY(); + fTrackletDelta[0]->Fill(x); + x[0]=p2.GetZ()-p1.GetZ(); + fTrackletDelta[1]->Fill(x); + x[0]=p2.GetSnp()-p1.GetSnp(); + fTrackletDelta[2]->Fill(x); + x[0]=p2.GetTgl()-p1.GetTgl(); + fTrackletDelta[3]->Fill(x); + TTreeSRedirector *cstream = GetDebugStreamer(); + if (cstream){ + (*cstream)<<"trackletMatch"<< + "tp1.="<=18){ - kpar = fSectorParamC; - kcov = fSectorCovarC; - } - for (Int_t ipar=0;ipar<6;ipar++){ - Int_t isec1 = s1%18; - Int_t isec2 = s2%18; - if (s1>35) isec1+=18; - if (s2>35) isec2+=18; - param6s1(ipar)=(*kpar)(6*isec1+ipar,0); - param6s2(ipar)=(*kpar)(6*isec2+ipar,0); + if (fSectorParamA){ + TMatrixD * kpar = fSectorParamA; + TMatrixD * kcov = fSectorCovarA; + if (s1%36>=18){ + kpar = fSectorParamC; + kcov = fSectorCovarC; + } + for (Int_t ipar=0;ipar<6;ipar++){ + Int_t isec1 = s1%18; + Int_t isec2 = s2%18; + if (s1>35) isec1+=18; + if (s2>35) isec2+=18; + param6s1(ipar)=(*kpar)(6*isec1+ipar,0); + param6s2(ipar)=(*kpar)(6*isec2+ipar,0); + } } - - + Double_t dy=0, dz=0, dphi=0,dtheta=0; Double_t sy=0, sz=0, sphi=0,stheta=0; Double_t ny=0, nz=0, nphi=0,ntheta=0; Double_t chi2v12=0, chi2v9=0, chi2v6=0; - Int_t npoints=0; - TLinearFitter * fitter = 0; + // Int_t npoints=0; + // TLinearFitter * fitter = 0; if (fPoints[GetIndex(s1,s2)]>minPoints){ // // // - fitter = GetFitter12(s1,s2); - npoints = fitter->GetNpoints(); - chi2v12 = TMath::Sqrt(fitter->GetChisquare()/npoints); +// fitter = GetFitter12(s1,s2); +// npoints = fitter->GetNpoints(); +// chi2v12 = TMath::Sqrt(fitter->GetChisquare()/npoints); - // - fitter = GetFitter9(s1,s2); - npoints = fitter->GetNpoints(); - chi2v9 = TMath::Sqrt(fitter->GetChisquare()/npoints); - // - fitter = GetFitter6(s1,s2); - npoints = fitter->GetNpoints(); - chi2v6 = TMath::Sqrt(fitter->GetChisquare()/npoints); - fitter->GetParameters(param6Diff); - // - GetTransformation6(s1,s2,m6); - GetTransformation9(s1,s2,m9); - GetTransformation12(s1,s2,m12); - // - fitter = GetFitter6(s1,s2); - //fitter->FixParameter(3,0); - //fitter->Eval(); - GetTransformation6(s1,s2,m6FX); +// // +// fitter = GetFitter9(s1,s2); +// npoints = fitter->GetNpoints(); +// chi2v9 = TMath::Sqrt(fitter->GetChisquare()/npoints); +// // +// fitter = GetFitter6(s1,s2); +// npoints = fitter->GetNpoints(); +// chi2v6 = TMath::Sqrt(fitter->GetChisquare()/npoints); +// fitter->GetParameters(param6Diff); +// // +// GetTransformation6(s1,s2,m6); +// GetTransformation9(s1,s2,m9); +// GetTransformation12(s1,s2,m12); +// // +// fitter = GetFitter6(s1,s2); +// //fitter->FixParameter(3,0); +// //fitter->Eval(); +// GetTransformation6(s1,s2,m6FX); // TH1 * his=0; his = GetHisto(kY,s1,s2); @@ -1942,28 +2062,13 @@ void AliTPCcalibAlign::MakeTree(const char *fname, Int_t minPoints){ // } + AliMagF* magF= (AliMagF*)TGeoGlobalMagField::Instance()->GetField(); + if (!magF) AliError("Magneticd field - not initialized"); + Double_t bz = magF->SolenoidField()/10.; //field in T - // x2 = a00*x1 + a01*y1 + a02*z1 + a03 - // y2 = a10*x1 + a11*y1 + a12*z1 + a13 - // z2 = a20*x1 + a21*y1 + a22*z1 + a23 - // dydx2 = (a10 + a11*dydx1 + a12*dzdx1)/(a00 + a01*dydx1 + a02*dzdx1) - // dzdx2 = (a20 + a21*dydx1 + a22*dzdx1)/(a00 + a01*dydx1 + a02*dzdx1) - // - // a00 a01 a02 a03 p[0] p[1] p[2] p[9] - // a10 a11 a12 a13 ==> p[3] p[4] p[5] p[10] - // a20 a21 a22 a23 p[6] p[7] p[8] p[11] - - // - // - // dy:-(fXIO*m6.fElements[4]+m6.fElements[7]) - // - // dphi:-(m6.fElements[4]) - // - // dz:fXIO*m6.fElements[8]+m6.fElements[11] - // - // dtheta:m6.fElements[8] - // cstream<<"Align"<< + "run="<GetImpactParameters(vertexXY,vertexZ); + if (TMath::Abs(vertexXY)>kVertexCut) return; + if (TMath::Abs(vertexZ)>kVertexCut) return; + if (TMath::Abs(seed->Pt())GetNumberOfClusters()GetSnp())>kSnpCut) return; + if (!fClusterDelta[0]) MakeResidualHistos(); + + Int_t detector=-1; + // + // + AliExternalTrackParam trackIn = *(fCurrentTrack->GetInnerParam()); + AliExternalTrackParam trackOut = *(fCurrentFriendTrack->GetTPCOut()); + static Double_t mass = TDatabasePDG::Instance()->GetParticle("pi+")->Mass(); + // + Int_t ncl=0; + for (Int_t irow=0; irow<160; irow++){ + AliTPCclusterMI *cl=seed->GetClusterPointer(irow); + if (!cl) continue; + if (cl->GetX()<80) continue; + if (detector<0) detector=cl->GetDetector()%36; + if (detector!=cl->GetDetector()%36) return; // cluster from different sectors + // skip such tracks + ncl++; + } + if (nclGetClusterPointer(irow); + if (!cl) continue; + if (cl->GetX()<80) continue; + if (detector<0) detector=cl->GetDetector()%36; + Int_t sector = cl->GetDetector(); + Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackOut.GetAlpha(); + if (TMath::Abs(dalpha)>0.01){ + if (!trackOut.Rotate(TMath::DegToRad()*(sector%18*20.+10.))) break; + } + Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()}; + Double_t cov[3]={0.01,0.,0.01}; + AliTPCseed::GetError(cl, &trackOut,cov[0],cov[2]); + Double_t dedge = cl->GetX()*TMath::Tan(TMath::Pi()/18.)-TMath::Abs(trackOut.GetY()); + cov[0]+=1./(irow+1.); // bigger error at boundary + cov[0]+=1./(160.-irow); // bigger error at boundary + cov[2]+=1./(irow+1.); // bigger error at boundary + cov[2]+=1./(160.-irow); // bigger error at boundary + cov[0]+=0.5/dedge; // bigger error close to the boundary + cov[2]+=0.5/dedge; // bigger error close to the boundary + cov[0]*=cov[0]; + cov[2]*=cov[2]; + if (!AliTracker::PropagateTrackToBxByBz(&trackOut, r[0],mass,1.,kFALSE)) continue; + + if (TMath::Abs(dedge)GetY()-trackOut.GetY())GetDetector()%36!=detector) continue; + // + // fill residual histogram + // + Double_t resVector[5]; + trackOut.GetXYZ(xyz); + resVector[1]= 9.*TMath::ATan2(xyz[1],xyz[0])/TMath::Pi(); + if (resVector[1]<0) resVector[1]+=18; + resVector[2]= TMath::Sqrt(cl->GetX()*cl->GetX()+cl->GetY()*cl->GetY()); + resVector[3]= cl->GetZ()/resVector[2]; + // + resVector[0]= cl->GetY()-trackOut.GetY(); + fClusterDelta[0]->Fill(resVector); + resVector[0]= cl->GetZ()-trackOut.GetZ(); + fClusterDelta[1]->Fill(resVector); + } + // + // Refit in - store residual maps + // + for (Int_t irow=159; irow>=0; irow--){ + AliTPCclusterMI *cl=seed->GetClusterPointer(irow); + if (!cl) continue; + if (cl->GetX()<80) continue; + if (detector<0) detector=cl->GetDetector()%36; + Int_t sector = cl->GetDetector(); + Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackIn.GetAlpha(); + if (TMath::Abs(dalpha)>0.01){ + if (!trackIn.Rotate(TMath::DegToRad()*(sector%18*20.+10.))) break; + } + Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()}; + Double_t cov[3]={0.01,0.,0.01}; + AliTPCseed::GetError(cl, &trackIn,cov[0],cov[2]); + Double_t dedge = cl->GetX()*TMath::Tan(TMath::Pi()/18.)-TMath::Abs(trackIn.GetY()); + cov[0]+=1./(irow+1.); // bigger error at boundary + cov[0]+=1./(160.-irow); // bigger error at boundary + cov[2]+=1./(irow+1.); // bigger error at boundary + cov[2]+=1./(160.-irow); // bigger error at boundary + cov[0]+=0.5/dedge; // bigger error close to the boundary +- + cov[2]+=0.5/dedge; // bigger error close to the boundary +- + cov[0]*=cov[0]; + cov[2]*=cov[2]; + if (!AliTracker::PropagateTrackToBxByBz(&trackIn, r[0],mass,1.,kFALSE)) continue; + if (TMath::Abs(dedge)GetY()-trackIn.GetY())GetDetector()%36!=detector) continue; + // + // fill residual histogram + // + Double_t resVector[5]; + trackIn.GetXYZ(xyz); + resVector[1]= 9.*TMath::ATan2(xyz[1],xyz[0])/TMath::Pi(); + if (resVector[1]<0) resVector[1]+=18; + resVector[2]= TMath::Sqrt(cl->GetX()*cl->GetX()+cl->GetY()*cl->GetY()); + resVector[3]= cl->GetZ()/resVector[2]; + // + resVector[0]= cl->GetY()-trackIn.GetY(); + fClusterDelta[0]->Fill(resVector); + resVector[0]= cl->GetZ()-trackIn.GetZ(); + fClusterDelta[1]->Fill(resVector); + } + + +} + + void AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){ // - // Update Kalman filter of Alignment + // Update Kalman filter of Alignment - only setup without filed // IROC - OROC quadrants // + if (TMath::Abs(AliTracker::GetBz())>0.5) return; if (!fClusterDelta[0]) MakeResidualHistos(); const Int_t kMinClusterF=40; + const Int_t kMinClusterFit=10; const Int_t kMinClusterQ=10; // const Int_t kdrow1Fit =5; // rows to skip from fit at the end @@ -2537,8 +2809,8 @@ void AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){ // // full track fit parameters // - TLinearFitter fyf(2,"pol1"); - TLinearFitter fzf(2,"pol1"); + static TLinearFitter fyf(2,"pol1"); // make it static - too much time for comiling of formula + static TLinearFitter fzf(2,"pol1"); // calgrind recomendation TVectorD pyf(2), peyf(2),pzf(2), pezf(2); TVectorD pyfc(2),pzfc(2); Int_t nf=0; @@ -2547,6 +2819,7 @@ void AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){ // for (Int_t iter=0; iter<2; iter++){ fyf.ClearPoints(); + fzf.ClearPoints(); for (Int_t irow=kdrow0Fit;irow<159-kdrow1Fit;irow++) { AliTPCclusterMI *c=track->GetClusterPointer(irow); if (!c) continue; @@ -2567,11 +2840,14 @@ void AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){ c->GetY(),yfit, c->GetZ(), pyf[1], c->GetMax(),2.5); if (TMath::Abs(corrtrY)>kMaxCorrY) continue; } - fyf.AddPoint(x,c->GetY(),0.1); - fzf.AddPoint(x,c->GetZ(),0.1); + if (TMath::Abs(x[0])<10){ + fyf.AddPoint(x,c->GetY(),0.1); //use only middle rows+-10cm + } + fzf.AddPoint(x,c->GetZ(),0.1); } nf = fyf.GetNpoints(); - if (nfGetTOFsignal(); // tof signal + // Double_t tofSignal= fCurrentTrack->GetTOFsignal(); // tof signal TVectorF tpcPosG(3); // global position of track at the middle of fXmiddle Double_t lphi = TMath::ATan2(pyf[0],fXmiddle); // expected local phi angle - if vertex at 0 Double_t gphi = 2.*TMath::Pi()*(isec%18+0.5)/18.+lphi; // expected global phi if vertex at 0 @@ -2613,8 +2889,8 @@ void AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){ // get constrained parameters // Double_t xvertex=vPosL[0]-fXmiddle; - fyf.AddPoint(&xvertex,vPosL[1], 0.1+TMath::Abs(vImpact[0])); - fzf.AddPoint(&xvertex,vPosL[2], 0.1+TMath::Abs(vImpact[1])); + fyf.AddPoint(&xvertex,vPosL[1], 0.00001); + fzf.AddPoint(&xvertex,vPosL[2], 2.); fyf.Eval(); fyf.GetParameters(pyfc); fzf.Eval(); @@ -2625,8 +2901,15 @@ void AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){ // 1-4 OROC quadrants // 0 IROC // - TLinearFitter *fittersY[5]; - TLinearFitter *fittersZ[5]; + static TLinearFitter *fittersY[5]={0,0,0,0,0}; // calgrind recomendation - fater to clear points + static TLinearFitter *fittersZ[5]={0,0,0,0,0}; // than create the fitter + if (fittersY[0]==0){ + for (Int_t i=0;i<5;i++) { + fittersY[i] = new TLinearFitter(2,"pol1"); + fittersZ[i] = new TLinearFitter(2,"pol1"); + } + } + // Int_t npoints[5]; TVectorD paramsY[5]; TVectorD errorsY[5]; @@ -2638,14 +2921,14 @@ void AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){ Double_t chi2CZ[5]; for (Int_t i=0;i<5;i++) { npoints[i]=0; - fittersY[i] = new TLinearFitter(2,"pol1"); paramsY[i].ResizeTo(2); errorsY[i].ResizeTo(2); covY[i].ResizeTo(2,2); - fittersZ[i] = new TLinearFitter(2,"pol1"); paramsZ[i].ResizeTo(2); errorsZ[i].ResizeTo(2); covZ[i].ResizeTo(2,2); + fittersY[i]->ClearPoints(); + fittersZ[i]->ClearPoints(); } // // Update fitters @@ -2677,23 +2960,17 @@ void AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){ // // Fill THnSparse cluster residuals // use only primary candidates with ITS signal - if (nf>100&&fCurrentTrack->IsOn(0x4)&&TMath::Abs(vImpact[0])<1&&TMath::Abs(vImpact[1])<1){ + if (fCurrentTrack->IsOn(0x4)&&TMath::Abs(vImpact[0])<1&&TMath::Abs(vImpact[1])<1){ Double_t resVector[5]; resVector[1]= 9.*gphi/TMath::Pi(); - resVector[2]= c->GetX(); - resVector[3]= c->GetY()/c->GetX(); - resVector[4]= c->GetZ()/c->GetX(); + resVector[2]= TMath::Sqrt(c->GetX()*c->GetX()+c->GetY()*c->GetY()); + resVector[3]= c->GetZ()/resVector[2]; // - resVector[0]= c->GetY()-yfit; - fClusterDelta[0]->Fill(resVector); - resVector[0]= c->GetZ()-zfit; - fClusterDelta[1]->Fill(resVector); // resVector[0]= c->GetY()-yfitC; - fClusterDelta[2]->Fill(resVector); + fClusterDelta[0]->Fill(resVector); resVector[0]= c->GetZ()-zfitC; - fClusterDelta[3]->Fill(resVector); - + fClusterDelta[1]->Fill(resVector); } if (c->GetRow()GetRow()>159-kdrow1Fit) continue; @@ -2754,11 +3031,6 @@ void AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){ covZ[i](1,1)*=chi2FacZ*chi2FacZ; } } - for (Int_t i=0;i<5;i++){ - delete fittersY[i]; - delete fittersZ[i]; - } - // // void UpdateSectorKalman // @@ -2828,7 +3100,7 @@ void AliTPCcalibAlign::UpdateAlignSector(const AliTPCseed * track,Int_t isec){ "vPosG.="<<&vPosG<< // vertex position in global system "vPosL.="<<&vPosL<< // vertex position in local system "vImpact.="<<&vImpact<< // track impact parameter - "tofSignal="<