1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
17 ///////////////////////////////////////////////////////////////////////////////
19 // Component for redoing the reconstruction from the clusters and tracks
21 // The new calibration data used
23 // In reality it overwrites the content of the ESD
28 gSystem->Load("libANALYSIS");
29 gSystem->Load("libTPCcalib");
31 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
32 gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
33 AliXRDPROOFtoolkit tool;
34 TChain * chainCl = tool.MakeChain("calib.txt","Clusters",0,1);
36 TChain * chainTr = tool.MakeChain("calib.txt","Tracks",0,1);
45 // marian.ivanov@cern.ch
47 #include "AliTPCcalibCalib.h"
50 #include "TTreeStream.h"
52 #include "TTimeStamp.h"
53 #include "AliESDEvent.h"
54 #include "AliESDfriend.h"
55 #include "AliESDtrack.h"
56 #include "AliTracker.h"
57 #include "AliTPCClusterParam.h"
59 #include "AliTPCcalibDB.h"
60 #include "AliTPCTransform.h"
61 #include "AliTPCRecoParam.h"
62 #include "AliTPCclusterMI.h"
63 #include "AliTPCseed.h"
64 #include "AliTPCPointCorrection.h"
66 ClassImp(AliTPCcalibCalib)
68 AliTPCcalibCalib::AliTPCcalibCalib():
70 fApplyExBCorrection(1), // apply ExB correction
71 fApplyTOFCorrection(1), // apply TOF correction
72 fApplyPositionCorrection(1), // apply position correction
73 fApplySectorAlignment(1), // apply sector alignment
74 fApplyRPhiCorrection(1), // apply R-Phi correction
75 fApplyRCorrection(1) // apply Radial correction
84 AliTPCcalibCalib::AliTPCcalibCalib(const Text_t *name, const Text_t *title)
86 fApplyExBCorrection(1), // apply ExB correction
87 fApplyTOFCorrection(1), // apply TOF correction
88 fApplyPositionCorrection(1), // apply position correction
89 fApplySectorAlignment(1), // apply sector alignment
90 fApplyRPhiCorrection(1), // apply R-Phi correction
91 fApplyRCorrection(1) // apply Radial correction
98 AliTPCcalibCalib::AliTPCcalibCalib(const AliTPCcalibCalib&calib):
99 AliTPCcalibBase(calib),
100 fApplyExBCorrection(calib.GetApplyExBCorrection()),
101 fApplyTOFCorrection(calib.GetApplyTOFCorrection()),
102 fApplyPositionCorrection(calib.GetApplyPositionCorrection()),
103 fApplySectorAlignment(calib.GetApplySectorAlignment()),
104 fApplyRPhiCorrection(calib.GetApplyRPhiCorrection()),
105 fApplyRCorrection(calib.GetApplyRCorrection())
113 AliTPCcalibCalib &AliTPCcalibCalib::operator=(const AliTPCcalibCalib&calib){
117 ((AliTPCcalibBase *)this)->operator=(calib);
122 AliTPCcalibCalib::~AliTPCcalibCalib() {
129 void AliTPCcalibCalib::Process(AliESDEvent *event){
136 AliESDfriend *ESDfriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
141 if (GetDebugLevel()>20) printf("Hallo world: Im here\n");
142 Int_t ntracks=event->GetNumberOfTracks();
143 //AliTPCcalibDB::Instance()->SetExBField(fMagF);
149 for (Int_t i=0;i<ntracks;++i) {
150 AliESDtrack *track = event->GetTrack(i);
151 AliESDfriendTrack *friendTrack = ESDfriend->GetTrack(i);
153 const AliExternalTrackParam * trackIn = track->GetInnerParam();
154 const AliExternalTrackParam * trackOut = track->GetOuterParam();
155 AliExternalTrackParam * tpcOut = (AliExternalTrackParam *)friendTrack->GetTPCOut();
156 if (!trackIn) continue;
157 if (!trackOut) continue;
158 if (!tpcOut) continue;
159 TObject *calibObject;
160 AliTPCseed *seed = 0;
161 for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) {
162 if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
165 RefitTrack(track, seed,event->GetMagneticField());
166 (*tpcOut)=*(track->GetOuterParam());
171 Bool_t AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed, Float_t magesd){
174 // if magesd==0 forget the curvature
177 // 0 - Setup transform object
179 AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform() ;
180 transform->SetCurrentRun(fRun);
181 transform->SetCurrentTimeStamp((UInt_t)fTime);
182 if(!fApplyExBCorrection) { // disable ExB correction in transform
183 if(transform->GetCurrentRecoParam())
184 transform->GetCurrentRecoParamNonConst()->SetUseExBCorrection(0);
186 if(!fApplyTOFCorrection) { // disable TOF correction in transform
187 if(transform->GetCurrentRecoParam())
188 transform->GetCurrentRecoParamNonConst()->SetUseTOFCorrection(kFALSE);
192 // First apply calibration
194 AliTPCPointCorrection * corr = AliTPCPointCorrection::Instance();
195 for (Int_t irow=0;irow<159;irow++) {
196 AliTPCclusterMI *cluster=seed->GetClusterPointer(irow);
197 if (!cluster) continue;
198 AliTPCclusterMI cl0(*cluster);
199 Double_t x[3]={cluster->GetRow(),cluster->GetPad(),cluster->GetTimeBin()};
200 Int_t i[1]={cluster->GetDetector()};
202 transform->Transform(x,i,0,1);
204 // get position correction
207 if (cluster->GetDetector()>35) ipad=1;
208 Float_t dy =AliTPCClusterParam::SPosCorrection(0,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax());
209 Float_t dz =AliTPCClusterParam::SPosCorrection(1,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax());
210 // if(fApplyPositionCorrection) {
215 // Apply sector alignment
217 Double_t dxq = AliTPCPointCorrection::SGetCorrectionSector(0,cluster->GetDetector()%36,cluster->GetX(),
218 cluster->GetY(),cluster->GetZ());
219 Double_t dyq = AliTPCPointCorrection::SGetCorrectionSector(1,cluster->GetDetector()%36,cluster->GetX(),
220 cluster->GetY(),cluster->GetZ());
221 Double_t dzq = AliTPCPointCorrection::SGetCorrectionSector(2,cluster->GetDetector()%36,cluster->GetX(),
222 cluster->GetY(),cluster->GetZ());
223 if (fApplySectorAlignment){
229 // Apply r-phi correction - To be done on track level- knowing the track angle !!!
232 corr->RPhiCOGCorrection(cluster->GetDetector(),cluster->GetRow(), cluster->GetPad(),
233 cluster->GetY(),cluster->GetY(), cluster->GetZ(), 0., cluster->GetMax(),2.5);
235 Double_t corrR = corr->CorrectionOutR0(kFALSE,kFALSE,cluster->GetX(),cluster->GetY(),cluster->GetZ(),cluster->GetDetector());
237 if (fApplyRPhiCorrection){
238 if (cluster->GetY()>0) x[1]+=corrclY; // rphi correction
239 if (cluster->GetY()<0) x[1]-=corrclY; // rphi correction
241 if (fApplyRCorrection){
242 x[0]+=corrR; // radial correction
252 TTreeSRedirector *cstream = GetDebugStreamer();
254 (*cstream)<<"Clusters"<<
255 "run="<<fRun<< // run number
256 "event="<<fEvent<< // event number
257 "time="<<fTime<< // time stamp of event
258 "trigger="<<fTrigger<< // trigger
259 "triggerClass="<<&fTriggerClass<< // trigger
260 "mag="<<fMagF<< // magnetic field
277 Int_t ncl = seed->GetNumberOfClusters();
279 for (Int_t i=0;i<15;i++) covar[i]=0;
282 covar[5]=10.*10./(64.*64.);
283 covar[9]=10.*10./(64.*64.);
285 if (TMath::Abs(magesd)<0.05) {
286 covar[14]=0.025*0.025;
291 AliExternalTrackParam * trackInOld = (AliExternalTrackParam*)track->GetInnerParam();
292 AliExternalTrackParam * trackOutOld = (AliExternalTrackParam*)track->GetOuterParam();
295 AliExternalTrackParam trackIn = *trackOutOld;
296 if (TMath::Abs(magesd)<0.05) {
297 ((Double_t&)(trackIn.GetParameter()[4]))=0.000000001;
299 trackIn.ResetCovariance(20.);
300 trackIn.AddCovariance(covar);
302 Int_t nclIn=0,nclOut=0;
307 for (Int_t irow=159; irow>0; irow--){
308 AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
310 if (cl->GetX()<80) continue;
311 Int_t sector = cl->GetDetector();
312 Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackIn.GetAlpha();
313 if (TMath::Abs(dalpha)>0.01){
314 if (!trackIn.Rotate(TMath::DegToRad()*(sector%18*20.+10.))) break;
316 Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()};
317 Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation
318 AliTPCseed::GetError(cl, &trackIn,cov[0],cov[2]);
322 Double_t bz = AliTracker::GetBz(xyz);
324 if (!trackIn.PropagateTo(r[0],bz)) continue;
325 if (RejectCluster(cl,&trackIn)) continue;
327 trackIn.Update(&r[1],cov);
330 AliExternalTrackParam trackOut = trackIn;
331 if (TMath::Abs(magesd)<0.05) {
332 ((Double_t&)(trackOut.GetParameter()[4]))=0.000000001;
334 trackOut.ResetCovariance(20.);
335 trackOut.AddCovariance(covar);
339 //Bool_t lastEdge=kFALSE;
340 for (Int_t irow=0; irow<160; irow++){
341 AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
343 if (cl->GetX()<80) continue;
344 Int_t sector = cl->GetDetector();
345 Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackOut.GetAlpha();
347 if (TMath::Abs(dalpha)>0.01){
348 if (!trackOut.Rotate(TMath::DegToRad()*(sector%18*20.+10.))) break;
350 Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()};
352 Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation
353 AliTPCseed::GetError(cl, &trackOut,cov[0],cov[2]);
356 trackOut.GetXYZ(xyz);
357 Double_t bz = AliTracker::GetBz(xyz);
358 if (!trackOut.PropagateTo(r[0],bz)) continue;
359 if (RejectCluster(cl,&trackOut)) continue;
361 trackOut.Update(&r[1],cov);
362 //if (cl->GetType()<0) lastEdge=kTRUE;
363 //if (cl->GetType()>=0) lastEdge=kFALSE;
370 trackIn.ResetCovariance(10.);
372 // Refit in one more time
374 for (Int_t irow=159; irow>0; irow--){
375 AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
377 if (cl->GetX()<80) continue;
378 Int_t sector = cl->GetDetector();
379 Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackIn.GetAlpha();
380 if (TMath::Abs(dalpha)>0.01){
381 if (!trackIn.Rotate(TMath::DegToRad()*(sector%18*20.+10.))) break;
383 Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()};
384 Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation
385 AliTPCseed::GetError(cl, &trackIn,cov[0],cov[2]);
389 Double_t bz = AliTracker::GetBz(xyz);
391 if (!trackIn.PropagateTo(r[0],bz)) continue;
392 if (RejectCluster(cl,&trackIn)) continue;
394 trackIn.Update(&r[1],cov);
398 trackIn.Rotate(trackInOld->GetAlpha());
399 trackOut.Rotate(trackOutOld->GetAlpha());
401 trackInOld->GetXYZ(xyz);
402 Double_t bz = AliTracker::GetBz(xyz);
403 trackIn.PropagateTo(trackInOld->GetX(),bz);
405 trackOutOld->GetXYZ(xyz);
406 bz = AliTracker::GetBz(xyz);
407 trackOut.PropagateTo(trackOutOld->GetX(),bz);
411 TTreeSRedirector *cstream = GetDebugStreamer();
413 (*cstream)<<"Tracks"<<
414 "run="<<fRun<< // run number
415 "event="<<fEvent<< // event number
416 "time="<<fTime<< // time stamp of event
417 "trigger="<<fTrigger<< // trigger
418 "triggerClass="<<&fTriggerClass<< // trigger
419 "mag="<<fMagF<< // magnetic field
423 "TrIn0.="<<trackInOld<<
424 "TrOut0.="<<trackOutOld<<
425 "TrIn1.="<<&trackIn<<
426 "TrOut1.="<<&trackOut<<
431 // And now rewrite ESDtrack and TPC seed
434 (*trackInOld) = trackIn;
435 (*trackOutOld) = trackOut;
436 AliExternalTrackParam *t = &trackIn;
437 track->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
438 seed->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
439 seed->SetNumberOfClusters((nclIn+nclOut)/2);
445 Bool_t AliTPCcalibCalib::RejectCluster(AliTPCclusterMI* cl, AliExternalTrackParam * param){
447 // check the acceptance of cluster
448 // Cut on edge effects
450 Float_t kEdgeCut=2.5;
453 Bool_t isReject = kFALSE;
454 Float_t edgeY = cl->GetX()*TMath::Tan(TMath::Pi()/18);
455 Float_t dist = edgeY - TMath::Abs(cl->GetY());
456 if (param) dist = TMath::Abs(edgeY - TMath::Abs(param->GetY()));
457 if (dist<kEdgeCut) isReject=kTRUE;
459 Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation
460 AliTPCseed::GetError(cl, param,cov[0],cov[2]);
461 if (param->GetSigmaY2()<0 || param->GetSigmaZ2()<0){
462 AliError("Wrong parameters");
465 Double_t py = (cl->GetY()-param->GetY())/TMath::Sqrt(cov[0]*cov[0]+param->GetSigmaY2());
466 Double_t pz = (cl->GetZ()-param->GetZ())/TMath::Sqrt(cov[2]*cov[2]+param->GetSigmaZ2());
468 if ((py*py+pz*pz)>kSigmaCut*kSigmaCut) isReject=kTRUE;