TPC/AliTPCcalibCalib.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   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  **************************************************************************/
  15
  16
  17 ///////////////////////////////////////////////////////////////////////////////
  18 //                                                                           //
  19 //   Component for redoing the reconstruction from the clusters and tracks
  20 //
  21 //   The new calibration data used
  22 //
  23 //   In reality it overwrites the content of the ESD
  24 //
  25
  26 /*
  27
  28   gSystem->Load("libANALYSIS");
  29   gSystem->Load("libTPCcalib");
  30   //
  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);
  35   chainCl->Lookup();
  36   TChain * chainTr = tool.MakeChain("calib.txt","Tracks",0,1);
  37   chainTr->Lookup();
  38
  39
  40
  41 */
  42
  43
  44
  45 //  marian.ivanov@cern.ch
  46 //
  47 #include "AliTPCcalibCalib.h"
  48 #include "TSystem.h"
  49 #include "TFile.h"
  50 #include "TTreeStream.h"
  51 #include "AliLog.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"
  58
  59 #include "AliTPCcalibDB.h"
  60 #include "AliTPCTransform.h"
  61 #include "AliTPCclusterMI.h"
  62 #include "AliTPCseed.h"
  63
  64 ClassImp(AliTPCcalibCalib)
  65
  66 AliTPCcalibCalib::AliTPCcalibCalib():
  67     AliTPCcalibBase()
  68 {
  69   //
  70   // Constructor
  71   //
  72 }
  73
  74
  75 AliTPCcalibCalib::AliTPCcalibCalib(const Text_t *name, const Text_t *title)
  76   :AliTPCcalibBase()
  77 {
  78   SetName(name);
  79   SetTitle(title);
  80 }
  81
  82
  83 AliTPCcalibCalib::AliTPCcalibCalib(const AliTPCcalibCalib&calib):
  84   AliTPCcalibBase(calib)
  85 {
  86   //
  87   // copy constructor
  88   //
  89 }
  90
  91 AliTPCcalibCalib &AliTPCcalibCalib::operator=(const AliTPCcalibCalib&calib){
  92   //
  93   //
  94   //
  95   ((AliTPCcalibBase *)this)->operator=(calib);
  96   return *this;
  97 }
  98
  99
 100 AliTPCcalibCalib::~AliTPCcalibCalib() {
 101   //
 102   // destructor
 103   //
 104 }
 105
 106
 107 void     AliTPCcalibCalib::Process(AliESDEvent *event){
 108   //
 109   //
 110   //
 111   if (!event) {
 112     return;
 113   }
 114   AliESDfriend *ESDfriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
 115   if (!ESDfriend) {
 116    return;
 117   }
 118
 119   if (GetDebugLevel()>20) printf("Hallo world: Im here\n");
 120   Int_t ntracks=event->GetNumberOfTracks();
 121   //AliTPCcalibDB::Instance()->SetExBField(fMagF);
 122
 123   //
 124   //
 125   //
 126
 127   for (Int_t i=0;i<ntracks;++i) {
 128     AliESDtrack *track = event->GetTrack(i);
 129     const AliExternalTrackParam * trackIn = track->GetInnerParam();
 130     const AliExternalTrackParam * trackOut = track->GetOuterParam();
 131     if (!trackIn) continue;
 132     if (!trackOut) continue;
 133
 134     AliESDfriendTrack *friendTrack = ESDfriend->GetTrack(i);
 135     TObject *calibObject;
 136     AliTPCseed *seed = 0;
 137     for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) {
 138       if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
 139     }
 140     if (!seed) continue;
 141     RefitTrack(track, seed);
 142   }
 143   return;
 144 }
 145
 146 Bool_t  AliTPCcalibCalib::RefitTrack(AliESDtrack * track, AliTPCseed *seed){
 147   //
 148   // Refit track
 149   //
 150
 151   //
 152   // First apply calibration
 153   //
 154
 155   AliTPCTransform *transform = AliTPCcalibDB::Instance()->GetTransform();
 156   for (Int_t irow=0;irow<159;irow++) {
 157     AliTPCclusterMI *cluster=seed->GetClusterPointer(irow);
 158     if (!cluster) continue;
 159     AliTPCclusterMI cl0(*cluster);
 160     Double_t x[3]={cluster->GetRow(),cluster->GetPad(),cluster->GetTimeBin()};
 161     Int_t i[1]={cluster->GetDetector()};
 162     transform->Transform(x,i,0,1);
 163     //
 164     // get position correction
 165     //
 166     Int_t ipad=0;
 167     if (cluster->GetDetector()>35) ipad=1;
 168     Float_t dy =AliTPCClusterParam::SPosCorrection(0,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax());
 169     Float_t dz =AliTPCClusterParam::SPosCorrection(1,ipad,cluster->GetPad(),cluster->GetTimeBin(),cluster->GetZ(),cluster->GetSigmaY2(),cluster->GetSigmaZ2(),cluster->GetMax());
 170     //x[1]-=dy;
 171     //x[2]-=dz;
 172
 173     //
 174     cluster->SetX(x[0]);
 175     cluster->SetY(x[1]);
 176     cluster->SetZ(x[2]);
 177     if (fStreamLevel>2){
 178       TTreeSRedirector *cstream = GetDebugStreamer();
 179       if (cstream){
 180         (*cstream)<<"Clusters"<<
 181           "run="<<fRun<<              //  run number
 182           "event="<<fEvent<<          //  event number
 183           "time="<<fTime<<            //  time stamp of event
 184           "trigger="<<fTrigger<<      //  trigger
 185           "mag="<<fMagF<<             //  magnetic field
 186           "cl0.="<<&cl0<<
 187           "cl.="<<cluster<<
 188           "cy="<<dy<<
 189           "cz="<<dz<<
 190           "\n";
 191       }
 192     }
 193   }
 194   Int_t ncl = seed->GetNumberOfClusters();
 195   Double_t covar[15];
 196   for (Int_t i=0;i<15;i++) covar[i]=0;
 197   covar[0]=10.*10.;
 198   covar[2]=10.*10.;
 199   covar[5]=10.*10./(64.*64.);
 200   covar[9]=10.*10./(64.*64.);
 201   covar[14]=1*1;
 202
 203   //
 204   // And now do refit
 205   //
 206   AliExternalTrackParam * trackInOld  = (AliExternalTrackParam*)track->GetInnerParam();
 207   AliExternalTrackParam * trackOutOld = (AliExternalTrackParam*)track->GetOuterParam();
 208
 209   AliExternalTrackParam trackIn  = *trackOutOld;
 210   AliExternalTrackParam trackOut = *trackInOld;
 211   trackIn.ResetCovariance(200.);
 212   trackOut.ResetCovariance(200.);
 213   trackIn.AddCovariance(covar);
 214   trackOut.AddCovariance(covar);
 215   Double_t xyz[3];
 216   Int_t nclIn=0,nclOut=0;
 217   //
 218   // Refit out
 219   //
 220   for (Int_t irow=0; irow<160; irow++){
 221     AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
 222     if (!cl) continue;
 223     if (cl->GetX()<80) continue;
 224     Int_t sector = cl->GetDetector();
 225     Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackOut.GetAlpha();
 226
 227     if (TMath::Abs(dalpha)>0.01)
 228       trackOut.Rotate(TMath::DegToRad()*(sector%18*20.+10.));
 229     Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()};
 230
 231     Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation
 232     AliTPCseed::GetError(cl, &trackOut,cov[0],cov[2]);
 233     cov[0]*=cov[0];
 234     cov[2]*=cov[2];
 235     trackOut.GetXYZ(xyz);
 236     Double_t bz = AliTracker::GetBz(xyz);
 237     if (trackOut.PropagateTo(r[0],bz)) nclOut++;
 238     if (RejectCluster(cl,&trackOut)) continue;
 239     trackOut.Update(&r[1],cov);
 240   }
 241   //
 242   // Refit in
 243   //
 244
 245   for (Int_t irow=159; irow>0; irow--){
 246     AliTPCclusterMI *cl=seed->GetClusterPointer(irow);
 247     if (!cl) continue;
 248     if (cl->GetX()<80) continue;
 249     Int_t sector = cl->GetDetector();
 250     Float_t dalpha = TMath::DegToRad()*(sector%18*20.+10.)-trackIn.GetAlpha();
 251     if (TMath::Abs(dalpha)>0.01)
 252       trackIn.Rotate(TMath::DegToRad()*(sector%18*20.+10.));
 253     Double_t r[3]={cl->GetX(),cl->GetY(),cl->GetZ()};
 254     Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation
 255     AliTPCseed::GetError(cl, &trackIn,cov[0],cov[2]);
 256     cov[0]*=cov[0];
 257     cov[2]*=cov[2];
 258     trackIn.GetXYZ(xyz);
 259     Double_t bz = AliTracker::GetBz(xyz);
 260
 261     if (trackIn.PropagateTo(r[0],bz)) nclIn++;
 262     if (RejectCluster(cl,&trackIn)) continue;
 263     trackIn.Update(&r[1],cov);
 264   }
 265   trackIn.Rotate(trackInOld->GetAlpha());
 266   trackOut.Rotate(trackOutOld->GetAlpha());
 267   //
 268   trackInOld->GetXYZ(xyz);
 269   Double_t bz = AliTracker::GetBz(xyz);
 270   trackIn.PropagateTo(trackInOld->GetX(),bz);
 271   //
 272   trackOutOld->GetXYZ(xyz);
 273   bz = AliTracker::GetBz(xyz);
 274   trackOut.PropagateTo(trackOutOld->GetX(),bz);
 275
 276
 277   if (fStreamLevel>0){
 278     TTreeSRedirector *cstream = GetDebugStreamer();
 279     if (cstream){
 280       (*cstream)<<"Tracks"<<
 281         "run="<<fRun<<              //  run number
 282         "event="<<fEvent<<          //  event number
 283         "time="<<fTime<<            //  time stamp of event
 284         "trigger="<<fTrigger<<      //  trigger
 285         "mag="<<fMagF<<             //  magnetic field
 286         "nclIn="<<nclIn<<
 287         "nclOut="<<nclOut<<
 288         "ncl="<<ncl<<
 289         "TrIn0.="<<trackInOld<<
 290         "TrOut0.="<<trackOutOld<<
 291         "TrIn1.="<<&trackIn<<
 292         "TrOut1.="<<&trackOut<<
 293         "\n";
 294     }
 295   }
 296   //
 297   // And now rewrite ESDtrack
 298   //
 299
 300   (*trackInOld)  = trackIn;
 301   (*trackOutOld) = trackOut;
 302   AliExternalTrackParam *t = &trackIn;
 303   track->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
 304   return kTRUE;
 305 }
 306
 307
 308
 309 Bool_t AliTPCcalibCalib::RejectCluster(AliTPCclusterMI* cl, AliExternalTrackParam * param){
 310   //
 311   // check the acceptance of cluster
 312   // Cut on edge effects
 313   //
 314   Bool_t isReject = kFALSE;
 315   Float_t edgeY = cl->GetX()*TMath::Tan(TMath::Pi()/18);
 316   Float_t dist  = edgeY - TMath::Abs(cl->GetY());
 317   if (param)  dist  = TMath::Abs(edgeY - TMath::Abs(param->GetY()));
 318   if (dist<3) isReject=kTRUE;
 319   if (cl->GetType()<0) isReject=kTRUE;
 320   return isReject;
 321 }
 322
 323
 324