+void AliTPCcalibAlign::UpdateClusterDeltaField(const AliTPCseed * seed){
+ //
+ // Update the cluster residula histograms for setup with field
+ // Kalman track fitting is used
+ // Only high momenta primary tracks used
+ //
+ // 1. Apply selection
+ // 2. Refit the track - in-out
+ // 3. Refit the track - out-in
+ // 4. Combine In and Out track - - fil cluster residuals
+ //
+ const Double_t kPtCut=1.0; // pt
+ const Double_t kSnpCut=0.2; // snp cut
+ const Double_t kNclCut=120; //
+ const Double_t kVertexCut=1;
+ const Double_t kMaxDist=0.5; // max distance between tracks and cluster
+ const Double_t kEdgeCut = 2.5;
+ const Double_t kDelta2=0.2*0.2; // initial increase in covar matrix
+ const Double_t kSigma=0.3; // error increase towards edges of TPC
+ const Double_t kSkipBoundary=7.5; // skip track updates in the boundary IFC,OFC, IO
+ //
+ if (!fCurrentTrack) return;
+ if (!fCurrentFriendTrack) return;
+ Float_t vertexXY=0,vertexZ=0;
+ fCurrentTrack->GetImpactParameters(vertexXY,vertexZ);
+ if (TMath::Abs(vertexXY)>kVertexCut) return;
+ if (TMath::Abs(vertexZ)>kVertexCut) return;
+ if (TMath::Abs(seed->Pt())<kPtCut) return;
+ if (seed->GetNumberOfClusters()<kNclCut) return;
+ if (TMath::Abs(seed->GetSnp())>kSnpCut) return;
+ if (!fClusterDelta[0]) MakeResidualHistos();
+ //
+ AliExternalTrackParam fitIn[160];
+ AliExternalTrackParam fitOut[160];
+ AliTPCROC * roc = AliTPCROC::Instance();
+ Double_t xmiddle = ( roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5;
+ Double_t xDiff = ( -roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5;
+ Double_t xIFC = ( roc->GetPadRowRadii(0,0));
+ Double_t xOFC = ( roc->GetPadRowRadii(36,roc->GetNRows(36)-1));
+ //
+ Int_t detector=-1;
+ //
+ //
+ AliExternalTrackParam trackIn = *(fCurrentTrack->GetInnerParam());
+ AliExternalTrackParam trackOut = *(fCurrentFriendTrack->GetTPCOut());
+ trackIn.ResetCovariance(10);
+ trackOut.ResetCovariance(10);
+ Double_t *covarIn = (Double_t*)trackIn.GetCovariance();
+ Double_t *covarOut = (Double_t*)trackOut.GetCovariance();
+ covarIn[0]+=kDelta2; covarIn[2]+=kDelta2;
+ covarIn[5]+=kDelta2/(100.*100.); covarIn[9]=kDelta2/(100.*100.);
+ covarIn[14]+=kDelta2/(5.*5.);
+ covarOut[0]+=kDelta2; covarOut[2]+=kDelta2;
+ covarOut[5]+=kDelta2/(100.*100.); covarOut[9]=kDelta2/(100.*100.);
+ covarOut[14]+=kDelta2/(5.*5.);
+ //
+ 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 (ncl<kNclCut) return;
+ Int_t nclIn=0,nclOut=0;
+ Double_t xyz[3];
+ //
+ // Refit out - store residual maps
+ //
+ 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;
+ Int_t sector = cl->GetDetector();
+ Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackOut.GetAlpha();
+ if (cl->GetDetector()%36!=detector) continue;
+ 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.1,0.,0.1};
+ Double_t dedge = cl->GetX()*TMath::Tan(TMath::Pi()/18.)-TMath::Abs(trackOut.GetY());
+ Double_t dmiddle = TMath::Abs(cl->GetX()-xmiddle)/xDiff;
+ dmiddle*=dmiddle;
+ //
+ cov[0]+=kSigma*dmiddle; // bigger error at boundary
+ cov[0]+=kSigma*dmiddle; // bigger error at boundary
+ cov[2]+=kSigma*dmiddle; // bigger error at boundary
+ cov[2]+=kSigma*dmiddle; // bigger error at boundary
+ cov[0]+=kSigma/dedge; // bigger error close to the boundary
+ cov[2]+=kSigma/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)<kEdgeCut) continue;
+ //
+ Bool_t doUpdate=kTRUE;
+ if (TMath::Abs(cl->GetX()-xIFC)<kSkipBoundary) doUpdate=kFALSE;
+ if (TMath::Abs(cl->GetX()-xOFC)<kSkipBoundary) doUpdate=kFALSE;
+ if (TMath::Abs(cl->GetX()-fXIO)<kSkipBoundary) doUpdate=kFALSE;
+ //
+ if (TMath::Abs(cl->GetY()-trackOut.GetY())<kMaxDist){
+ nclOut++;
+ if (doUpdate) trackOut.Update(&r[1],cov);
+ }
+ fitOut[irow]=trackOut;
+ }
+
+ //
+ // 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 (cl->GetDetector()%36!=detector) continue;
+ 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.1,0.,0.1};
+ Double_t dedge = cl->GetX()*TMath::Tan(TMath::Pi()/18.)-TMath::Abs(trackIn.GetY());
+ Double_t dmiddle = TMath::Abs(cl->GetX()-xmiddle)/xDiff;
+ dmiddle*=dmiddle;
+ //
+ cov[0]+=kSigma*dmiddle; // bigger error at boundary
+ cov[0]+=kSigma*dmiddle; // bigger error at boundary
+ cov[2]+=kSigma*dmiddle; // bigger error at boundary
+ cov[2]+=kSigma*dmiddle; // bigger error at boundary
+ cov[0]+=kSigma/dedge; // bigger error close to the boundary
+ cov[2]+=kSigma/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)<kEdgeCut) continue;
+ Bool_t doUpdate=kTRUE;
+ if (TMath::Abs(cl->GetX()-xIFC)<kSkipBoundary) doUpdate=kFALSE;
+ if (TMath::Abs(cl->GetX()-xOFC)<kSkipBoundary) doUpdate=kFALSE;
+ if (TMath::Abs(cl->GetX()-fXIO)<kSkipBoundary) doUpdate=kFALSE;
+ if (TMath::Abs(cl->GetY()-trackIn.GetY())<kMaxDist){
+ nclIn++;
+ if (doUpdate) trackIn.Update(&r[1],cov);
+ }
+ fitIn[irow]=trackIn;
+ }
+ //
+ //
+ for (Int_t irow=159; irow>=0; irow--){
+ //
+ // Update kalman - +- direction
+ // Store cluster residuals
+ AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
+ if (!cl) continue;
+ if (cl->GetX()<80) continue;
+ if (detector<0) detector=cl->GetDetector()%36;
+ if (cl->GetDetector()%36!=detector) continue;
+ AliExternalTrackParam trackSmooth = fitIn[irow];
+ AliTrackerBase::UpdateTrack(trackSmooth, fitOut[irow]);
+ //
+ Double_t resVector[5];
+ trackSmooth.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()-trackSmooth.GetY();
+ fClusterDelta[0]->Fill(resVector);
+ resVector[0]= cl->GetZ()-trackSmooth.GetZ();
+ fClusterDelta[1]->Fill(resVector);
+ }
+
+}
+
+