#include "AliESDtrack.h"
#include "AliTracker.h"
#include "AliTPCClusterParam.h"
+#include "AliTPCParam.h"
#include "AliTPCcalibDB.h"
#include "AliTPCTransform.h"
+#include "AliTPCRecoParam.h"
#include "AliTPCclusterMI.h"
#include "AliTPCseed.h"
-
+#include "AliTPCPointCorrection.h"
+#include <TGeoManager.h>
+#include <TGeoPhysicalNode.h>
+#include "TDatabasePDG.h"
ClassImp(AliTPCcalibCalib)
AliTPCcalibCalib::AliTPCcalibCalib():
- AliTPCcalibBase()
+AliTPCcalibBase(),
+ fApplyExBCorrection(1), // apply ExB correction
+ fApplyTOFCorrection(1), // apply TOF correction
+ fApplyPositionCorrection(1), // apply position correction
+ fApplySectorAlignment(1), // apply sector alignment
+ fApplyRPhiCorrection(1), // apply R-Phi correction
+ fApplyRCorrection(1) // apply Radial correction
+
{
//
// Constructor
AliTPCcalibCalib::AliTPCcalibCalib(const Text_t *name, const Text_t *title)
- :AliTPCcalibBase()
+ :AliTPCcalibBase(),
+ fApplyExBCorrection(1), // apply ExB correction
+ fApplyTOFCorrection(1), // apply TOF correction
+ fApplyPositionCorrection(1), // apply position correction
+ fApplySectorAlignment(1), // apply sector alignment
+ fApplyRPhiCorrection(1), // apply R-Phi correction
+ fApplyRCorrection(1) // apply Radial correction
{
SetName(name);
SetTitle(title);
AliTPCcalibCalib::AliTPCcalibCalib(const AliTPCcalibCalib&calib):
- AliTPCcalibBase(calib)
+ AliTPCcalibBase(calib),
+ fApplyExBCorrection(calib.GetApplyExBCorrection()),
+ fApplyTOFCorrection(calib.GetApplyTOFCorrection()),
+ fApplyPositionCorrection(calib.GetApplyPositionCorrection()),
+ fApplySectorAlignment(calib.GetApplySectorAlignment()),
+ fApplyRPhiCorrection(calib.GetApplyRPhiCorrection()),
+ fApplyRCorrection(calib.GetApplyRCorrection())
+
{
//
// copy constructor
if (GetDebugLevel()>20) printf("Hallo world: Im here\n");
Int_t ntracks=event->GetNumberOfTracks();
- AliTPCcalibDB::Instance()->SetExBField(fMagF);
+ //AliTPCcalibDB::Instance()->SetExBField(fMagF);
//
//
//
for (Int_t i=0;i<ntracks;++i) {
- AliESDtrack *track = event->GetTrack(i);
- const AliExternalTrackParam * trackIn = track->GetInnerParam();
+ AliESDtrack *track = event->GetTrack(i);
+ AliESDfriendTrack *friendTrack = (AliESDfriendTrack*) ESDfriend->GetTrack(i);
+ if (!friendTrack) continue;
+ //track->SetFriendTrack(friendTrack);
+ fCurrentFriendTrack=friendTrack;
+ const AliExternalTrackParam * trackIn = track->GetInnerParam();
const AliExternalTrackParam * trackOut = track->GetOuterParam();
+ AliExternalTrackParam * tpcOut = (AliExternalTrackParam *)friendTrack->GetTPCOut();
if (!trackIn) continue;
if (!trackOut) continue;
-
- AliESDfriendTrack *friendTrack = ESDfriend->GetTrack(i);
+ if (!tpcOut) continue;
TObject *calibObject;
AliTPCseed *seed = 0;
for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) {
if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
}
if (!seed) continue;
- RefitTrack(track, seed);
+ RefitTrack(track, seed, event->GetMagneticField());
+ (*tpcOut)=*(track->GetOuterParam());
}
return;
}
-Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){
+Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed, Float_t magesd){
//
// Refit track
+ // if magesd==0 forget the curvature
+
+ //
+ // 0 - Setup transform object
//
+ const Double_t kxIFC = 83.; // position of IFC
+ const Double_t kxOFC = 250.; // position of OFC
+ const Double_t kaFC = 1.; // amplitude
+ const Double_t ktFC = 5.0; // slope of error
+ //cov[0]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC));
+ //cov[2]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC));
+
+ static Int_t streamCounter=0;
+ streamCounter++;
+ AliESDfriendTrack *friendTrack = fCurrentFriendTrack;
+
+ AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;
+ AliTPCParam *param = AliTPCcalibDB::Instance()->GetParameters();
+ transform->SetCurrentRun(fRun);
+ transform->SetCurrentTimeStamp((UInt_t)fTime);
+ if(!fApplyExBCorrection) { // disable ExB correction in transform
+ if(transform->GetCurrentRecoParam())
+ transform->GetCurrentRecoParamNonConst()->SetUseExBCorrection(0);
+ }
+ if(!fApplyTOFCorrection) { // disable TOF correction in transform
+ if(transform->GetCurrentRecoParam())
+ transform->GetCurrentRecoParamNonConst()->SetUseTOFCorrection(kFALSE);
+ }
//
// First apply calibration
//
-
- AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform();
+ // AliTPCPointCorrection * corr = AliTPCPointCorrection::Instance();
for (Int_t irow=0;irow<159;irow++) {
AliTPCclusterMI *cluster=seed->GetClusterPointer(irow);
if (!cluster) continue;
AliTPCclusterMI cl0(*cluster);
Double_t x[3]={cluster->GetRow(),cluster->GetPad(),cluster->GetTimeBin()};
Int_t i[1]={cluster->GetDetector()};
+
transform->Transform(x,i,0,1);
//
// get position correction
//
Int_t ipad=0;
if (cluster->GetDetector()>35) ipad=1;
- Float_t dy =AliTPCClusterParam::SPosCorrection(0,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax());
- Float_t dz =AliTPCClusterParam::SPosCorrection(1,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax());
- //x[1]-=dy;
- //x[2]-=dz;
-
+ Float_t dy =0;//AliTPCClusterParam::SPosCorrection(0,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax());
+ Float_t dz =0;//AliTPCClusterParam::SPosCorrection(1,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax());
//
cluster->SetX(x[0]);
cluster->SetY(x[1]);
cluster->SetZ(x[2]);
- if (fStreamLevel>2){
+
+ //
+ // Apply alignemnt
+ //
+ if (transform->GetCurrentRecoParam()->GetUseSectorAlignment()){
+ if (!param->IsGeoRead()) param->ReadGeoMatrices();
+ TGeoHMatrix *mat = param->GetClusterMatrix(cluster->GetDetector());
+ //TGeoHMatrix mat;
+ Double_t pos[3]= {cluster->GetX(),cluster->GetY(),cluster->GetZ()};
+ Double_t posC[3]={cluster->GetX(),cluster->GetY(),cluster->GetZ()};
+ if (mat) mat->LocalToMaster(pos,posC);
+ else{
+ // chack Loading of Geo matrices from GeoManager - TEMPORARY FIX
+ }
+ cluster->SetX(posC[0]);
+ cluster->SetY(posC[1]);
+ cluster->SetZ(posC[2]);
+ }
+
+
+
+ if (fStreamLevel>2 && streamCounter<20*fStreamLevel ){
+ // dump debug info if required
TTreeSRedirector *cstream = GetDebugStreamer();
if (cstream){
(*cstream)<<"Clusters"<<
"event="<<fEvent<< // event number
"time="<<fTime<< // time stamp of event
"trigger="<<fTrigger<< // trigger
+ "triggerClass="<<&fTriggerClass<< // trigger
"mag="<<fMagF<< // magnetic field
"cl0.="<<&cl0<<
"cl.="<<cluster<<
}
}
}
+ //
+ //
+ //
Int_t ncl = seed->GetNumberOfClusters();
+ const Double_t kResetCov=4.;
+ const Double_t kSigma=5.;
Double_t covar[15];
for (Int_t i=0;i<15;i++) covar[i]=0;
- covar[0]=10.*10.;
- covar[2]=10.*10.;
- covar[5]=10.*10./(64.*64.);
- covar[9]=10.*10./(64.*64.);
- covar[14]=1*1;
-
+ covar[0]=kSigma*kSigma;
+ covar[2]=kSigma*kSigma;
+ covar[5]=kSigma*kSigma/Float_t(ncl*ncl);
+ covar[9]=kSigma*kSigma/Float_t(ncl*ncl);
+ covar[14]=0.2*0.2;
+ if (TMath::Abs(magesd)<0.05) {
+ covar[14]=0.025*0.025;
+ }
//
// And now do refit
//
AliExternalTrackParam * trackInOld = (AliExternalTrackParam*)track->GetInnerParam();
- AliExternalTrackParam * trackOutOld = (AliExternalTrackParam*)track->GetOuterParam();
+ AliExternalTrackParam * trackOuter = (AliExternalTrackParam*)track->GetOuterParam();
+ AliExternalTrackParam * trackOutOld = (AliExternalTrackParam *)friendTrack->GetTPCOut();
+ Double_t mass = TDatabasePDG::Instance()->GetParticle("pi+")->Mass();
+
+
AliExternalTrackParam trackIn = *trackOutOld;
- AliExternalTrackParam trackOut = *trackInOld;
- trackIn.ResetCovariance(200.);
- trackOut.ResetCovariance(200.);
+ trackIn.ResetCovariance(kResetCov);
trackIn.AddCovariance(covar);
- trackOut.AddCovariance(covar);
+ if (TMath::Abs(magesd)<0.05) {
+ ((Double_t&)(trackIn.GetParameter()[4]))=0.000000001;
+ ((Double_t&)(trackIn.GetCovariance()[14]))=covar[14]; // fix the line
+ }
+
Double_t xyz[3];
Int_t nclIn=0,nclOut=0;
+ //
+ // Refit in
+ //
+ for (Int_t irow=159; irow>0; irow--){
+ AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
+ if (!cl) continue;
+ if (cl->GetX()<80) continue;
+ 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}; //TODO: correct error parametrisation
+ AliTPCseed::GetError(cl, &trackIn,cov[0],cov[2]);
+ cov[0]*=cov[0];
+ cov[2]*=cov[2];
+ cov[0]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC));
+ cov[2]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC));
+ trackIn.GetXYZ(xyz);
+ // Double_t bz = AliTracker::GetBz(xyz);
+
+ // if (!trackIn.PropagateTo(r[0],bz)) continue;
+ if (!AliTracker::PropagateTrackToBxByBz(&trackIn, r[0],mass,1.,kFALSE)) continue;
+
+ if (RejectCluster(cl,&trackIn)) continue;
+ nclIn++;
+ trackIn.Update(&r[1],cov);
+ }
+ //
+ AliExternalTrackParam trackOut = trackIn;
+ trackOut.ResetCovariance(kResetCov);
+ trackOut.AddCovariance(covar);
+ if (TMath::Abs(magesd)<0.05) {
+ ((Double_t&)(trackOut.GetParameter()[4]))=0.000000001;
+ ((Double_t&)(trackOut.GetCovariance()[14]))=covar[14]; // fix the line
+ }
+
//
// Refit out
//
+ //Bool_t lastEdge=kFALSE;
for (Int_t irow=0; irow<160; irow++){
AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
if (!cl) continue;
Int_t sector = cl->GetDetector();
Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackOut.GetAlpha();
- if (TMath::Abs(dalpha)>0.01)
- trackOut.Rotate(TMath::DegToRad()*(sector%18*20.+10.));
+ 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}; //TODO: correct error parametrisation
+
+ Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation
+ AliTPCseed::GetError(cl, &trackOut,cov[0],cov[2]);
+ cov[0]*=cov[0];
+ cov[2]*=cov[2];
+ cov[0]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC));
+ cov[2]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC));
trackOut.GetXYZ(xyz);
- Double_t bz = AliTracker::GetBz(xyz);
- if (trackOut.PropagateTo(r[0],bz)) nclOut++;
+ //Double_t bz = AliTracker::GetBz(xyz);
+ // if (!trackOut.PropagateTo(r[0],bz)) continue;
+ if (!AliTracker::PropagateTrackToBxByBz(&trackOut, r[0],mass,1.,kFALSE)) continue;
+
if (RejectCluster(cl,&trackOut)) continue;
- trackOut.Update(&r[1],cov);
+ nclOut++;
+ trackOut.Update(&r[1],cov);
+ //if (cl->GetType()<0) lastEdge=kTRUE;
+ //if (cl->GetType()>=0) lastEdge=kFALSE;
}
//
- // Refit in
//
-
+ //
+ nclIn=0;
+ trackIn = trackOut;
+ trackIn.ResetCovariance(kResetCov);
+ if (TMath::Abs(magesd)<0.05) {
+ ((Double_t&)(trackIn.GetParameter()[4]))=0.000000001;
+ ((Double_t&)(trackIn.GetCovariance()[14]))=covar[14]; // fix the line
+ }
+ //
+ // Refit in one more time
+ //
for (Int_t irow=159; irow>0; irow--){
AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
if (!cl) continue;
if (cl->GetX()<80) continue;
Int_t sector = cl->GetDetector();
Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackIn.GetAlpha();
- if (TMath::Abs(dalpha)>0.01)
- trackIn.Rotate(TMath::DegToRad()*(sector%18*20.+10.));
+ 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}; //TODO: correct error parametrisation
- trackOut.GetXYZ(xyz);
- Double_t bz = AliTracker::GetBz(xyz);
+ AliTPCseed::GetError(cl, &trackIn,cov[0],cov[2]);
+ cov[0]*=cov[0];
+ cov[2]*=cov[2];
+ cov[0]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC));
+ cov[2]+= kaFC*(TMath::Exp(-TMath::Abs(cl->GetX()-kxIFC)/ktFC)+TMath::Exp(-TMath::Abs(cl->GetX()-kxOFC)/ktFC));
+
+ trackIn.GetXYZ(xyz);
+ //Double_t bz = AliTracker::GetBz(xyz);
+
+ // if (!trackIn.PropagateTo(r[0],bz)) continue;
+ if (!AliTracker::PropagateTrackToBxByBz(&trackIn, r[0],mass,1,kFALSE)) continue;
- if (trackIn.PropagateTo(r[0],bz)) nclIn++;
if (RejectCluster(cl,&trackIn)) continue;
+ nclIn++;
trackIn.Update(&r[1],cov);
}
+
+
trackIn.Rotate(trackInOld->GetAlpha());
trackOut.Rotate(trackOutOld->GetAlpha());
//
trackOut.PropagateTo(trackOutOld->GetX(),bz);
- if (fStreamLevel>0){
+ if (fStreamLevel>0 && streamCounter<100*fStreamLevel){
TTreeSRedirector *cstream = GetDebugStreamer();
if (cstream){
(*cstream)<<"Tracks"<<
"event="<<fEvent<< // event number
"time="<<fTime<< // time stamp of event
"trigger="<<fTrigger<< // trigger
+ "triggerClass="<<&fTriggerClass<< // trigger
"mag="<<fMagF<< // magnetic field
"nclIn="<<nclIn<<
"nclOut="<<nclOut<<
}
}
//
- // And now rewrite ESDtrack
+ // And now rewrite ESDtrack and TPC seed
//
(*trackInOld) = trackIn;
(*trackOutOld) = trackOut;
+ (*trackOuter) = trackOut;
AliExternalTrackParam *t = &trackIn;
- track->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
+ //track->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
+ seed->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
+ seed->SetNumberOfClusters((nclIn+nclOut)/2);
return kTRUE;
}
// check the acceptance of cluster
// Cut on edge effects
//
+ if (!param) return kTRUE;
+ Float_t kEdgeCut=2.5;
+ Float_t kSigmaCut=6;
+
Bool_t isReject = kFALSE;
Float_t edgeY = cl->GetX()*TMath::Tan(TMath::Pi()/18);
Float_t dist = edgeY - TMath::Abs(cl->GetY());
- if (param) dist = TMath::Abs(edgeY - TMath::Abs(param->GetY()));
- if (dist<3) isReject=kTRUE;
- if (cl->GetType()<0) isReject=kTRUE;
+ dist = TMath::Abs(edgeY - TMath::Abs(param->GetY()));
+ if (dist<kEdgeCut) isReject=kTRUE;
+
+ Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation
+ AliTPCseed::GetError(cl, param,cov[0],cov[2]);
+ if (param->GetSigmaY2()<0 || param->GetSigmaZ2()<0){
+ AliError("Wrong parameters");
+ return kFALSE;
+ }
+ Double_t py = (cl->GetY()-param->GetY())/TMath::Sqrt(cov[0]*cov[0]+param->GetSigmaY2());
+ Double_t pz = (cl->GetZ()-param->GetZ())/TMath::Sqrt(cov[2]*cov[2]+param->GetSigmaZ2());
+ //
+ if ((py*py+pz*pz)>kSigmaCut*kSigmaCut) isReject=kTRUE;
+
return isReject;
}
+