TPC/AliTPCcalibTime.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 Comments to be written here:
  18
  19 1. What do we calibrate.
  20
  21   Time dependence of gain and drift velocity in order to account for changes in: temperature, pressure, gas composition.
  22
  23   AliTPCcalibTime *calibTime = new AliTPCcalibTime("cosmicTime","cosmicTime",0, 1213.9e+06, 1213.96e+06, 0.04e+04, 0.04e+04);
  24
  25 2. How to interpret results
  26
  27 3. Simple example
  28
  29   a) determine the required time range:
  30
  31   AliXRDPROOFtoolkit tool;
  32   TChain * chain = tool.MakeChain("pass2.txt","esdTree",0,6000);
  33   chain->Draw("GetTimeStamp()")
  34
  35   b) analyse calibration object on Proof in calibration train
  36
  37   AliTPCcalibTime *calibTime = new AliTPCcalibTime("cosmicTime","cosmicTime", StartTimeStamp, EndTimeStamp, IntegrationTimeVdrift);
  38
  39   c) plot results
  40   .x ~/NimStyle.C
  41   gSystem->Load("libANALYSIS");
  42   gSystem->Load("libTPCcalib");
  43
  44   TFile f("CalibObjects.root");
  45   AliTPCcalibTime *cal = (AliTPCcalibTime *)f->Get("TPCCalib")->FindObject("calibTime");
  46   cal->GetHistoDrift("all")->Projection(2,0)->Draw()
  47   cal->GetFitDrift("all")->Draw("lp")
  48
  49 4. Analysis using debug streamers.
  50
  51   gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
  52   gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
  53   AliXRDPROOFtoolkit tool;
  54   TChain * chainTime = tool.MakeChain("time.txt","timeInfo",0,10200);
  55   chainTime->Lookup();
  56 */
  57
  58 #include "Riostream.h"
  59 #include "TChain.h"
  60 #include "TTree.h"
  61 #include "TH1F.h"
  62 #include "TH2F.h"
  63 #include "TH3F.h"
  64 #include "THnSparse.h"
  65 #include "TList.h"
  66 #include "TMath.h"
  67 #include "TCanvas.h"
  68 #include "TFile.h"
  69 #include "TF1.h"
  70 #include "TVectorD.h"
  71 #include "TProfile.h"
  72 #include "TGraphErrors.h"
  73 #include "TCanvas.h"
  74
  75 #include "AliTPCclusterMI.h"
  76 #include "AliTPCseed.h"
  77 #include "AliESDVertex.h"
  78 #include "AliESDEvent.h"
  79 #include "AliESDfriend.h"
  80 #include "AliESDInputHandler.h"
  81 #include "AliAnalysisManager.h"
  82
  83 #include "AliTracker.h"
  84 #include "AliMagF.h"
  85 #include "AliTPCCalROC.h"
  86
  87 #include "AliLog.h"
  88
  89 #include "AliTPCcalibTime.h"
  90
  91 #include "TTreeStream.h"
  92 #include "AliTPCTracklet.h"
  93 #include "TTimeStamp.h"
  94 #include "AliTPCcalibDB.h"
  95 #include "AliTPCcalibLaser.h"
  96 #include "AliDCSSensorArray.h"
  97 #include "AliDCSSensor.h"
  98
  99 ClassImp(AliTPCcalibTime)
 100
 101
 102 AliTPCcalibTime::AliTPCcalibTime()
 103   :AliTPCcalibBase(),
 104    fLaser(0),       // pointer to laser calibration
 105    fDz(0),          // current delta z
 106    fCutMaxD(2),        // maximal distance in rfi ditection
 107    fCutMaxDz(20),      // maximal distance in rfi ditection
 108    fCutTheta(0.03),    // maximal distan theta
 109    fCutMinDir(-0.99),  // direction vector products
 110    fCutTracks(10),
 111    fMapDz(0),          //NEW! Tmap of V drifts for different triggers
 112    fTimeBins(0),
 113    fTimeStart(0),
 114    fTimeEnd(0),
 115    fPtBins(0),
 116    fPtStart(0),
 117    fPtEnd(0),
 118    fVdriftBins(0),
 119    fVdriftStart(0),
 120    fVdriftEnd(0),
 121    fRunBins(0),
 122    fRunStart(0),
 123    fRunEnd(0)
 124 //   fBinsVdrift(fTimeBins,fPtBins,fVdriftBins),
 125 //   fXminVdrift(fTimeStart,fPtStart,fVdriftStart),
 126 //   fXmaxVdrift(fTimeEnd,fPtEnd,fVdriftEnd)
 127
 128 {
 129   AliInfo("Default Constructor");
 130   for (Int_t i=0;i<3;i++) {
 131     fHistVdriftLaserA[i]=0;
 132     fHistVdriftLaserC[i]=0;
 133   }
 134   for (Int_t i=0;i<5;i++) {
 135     fCosmiMatchingHisto[i]=0;
 136   }
 137 }
 138
 139
 140 AliTPCcalibTime::AliTPCcalibTime(const Text_t *name, const Text_t *title, UInt_t StartTime, UInt_t EndTime, Int_t deltaIntegrationTimeVdrift)
 141   :AliTPCcalibBase(),
 142    fLaser(0),       // pointer to laser calibration
 143    fDz(0),          // current delta z
 144    fCutMaxD(2),        // maximal distance in rfi ditection
 145    fCutMaxDz(20),      // maximal distance in rfi ditection
 146    fCutTheta(0.03),    // maximal distan theta
 147    fCutMinDir(-0.99),  // direction vector products
 148    fCutTracks(10),
 149    fMapDz(0),          //NEW! Tmap of V drifts for different triggers
 150    fTimeBins(0),
 151    fTimeStart(0),
 152    fTimeEnd(0),
 153    fPtBins(0),
 154    fPtStart(0),
 155    fPtEnd(0),
 156    fVdriftBins(0),
 157    fVdriftStart(0),
 158    fVdriftEnd(0),
 159    fRunBins(0),
 160    fRunStart(0),
 161    fRunEnd(0)
 162 {
 163
 164   SetName(name);
 165   SetTitle(title);
 166   for (Int_t i=0;i<3;i++) {
 167     fHistVdriftLaserA[i]=0;
 168     fHistVdriftLaserC[i]=0;
 169   }
 170
 171   AliInfo("Non Default Constructor");
 172
 173   fTimeBins   =(EndTime-StartTime)/deltaIntegrationTimeVdrift;
 174   fTimeStart  =StartTime; //(((TObjString*)(mapGRP->GetValue("fAliceStartTime")))->GetString()).Atoi();
 175   fTimeEnd    =EndTime;   //(((TObjString*)(mapGRP->GetValue("fAliceStopTime")))->GetString()).Atoi();
 176   fPtBins     = 200;
 177   fPtStart    = -0.04;
 178   fPtEnd      =  0.04;
 179   fVdriftBins = 200;
 180   fVdriftStart= -20.0/500.0;
 181   fVdriftEnd  =  20.0/500.0;
 182   fRunBins    = 100000;
 183   fRunStart   = -0.5;
 184   fRunEnd     = 0.5;
 185
 186   Int_t    binsVdriftLaser[4] = {fTimeBins , fPtBins , fVdriftBins*20, fRunBins };
 187   Double_t xminVdriftLaser[4] = {fTimeStart, fPtStart, fVdriftStart  , fRunStart};
 188   Double_t xmaxVdriftLaser[4] = {fTimeEnd  , fPtEnd  , fVdriftEnd    , fRunEnd  };
 189
 190   for (Int_t i=0;i<3;i++) {
 191     fHistVdriftLaserA[i] = new THnSparseF("HistVdriftLaser","HistVdriftLaser;time;p/T ratio;Vdrift;run",4,binsVdriftLaser,xminVdriftLaser,xmaxVdriftLaser);
 192     fHistVdriftLaserC[i] = new THnSparseF("HistVdriftLaser","HistVdriftLaser;time;p/T ratio;Vdrift;run",4,binsVdriftLaser,xminVdriftLaser,xmaxVdriftLaser);
 193   }
 194   fBinsVdrift[0] = fTimeBins;
 195   fBinsVdrift[1] = fPtBins;
 196   fBinsVdrift[2] = fVdriftBins;
 197   fBinsVdrift[3] = fRunBins;
 198   fXminVdrift[0] = fTimeStart;
 199   fXminVdrift[1] = fPtStart;
 200   fXminVdrift[2] = fVdriftStart;
 201   fXminVdrift[3] = fRunStart;
 202   fXmaxVdrift[0] = fTimeEnd;
 203   fXmaxVdrift[1] = fPtEnd;
 204   fXmaxVdrift[2] = fVdriftEnd;
 205   fXmaxVdrift[3] = fRunEnd;
 206
 207   fMapDz=new TMap();
 208
 209   for (Int_t i=0;i<5;i++) {
 210     fCosmiMatchingHisto[i]=new TH1F("Cosmics matching","Cosmics matching",100,0,0);
 211   }
 212 }
 213
 214 AliTPCcalibTime::~AliTPCcalibTime(){
 215   //
 216   // Destructor
 217   //
 218   for (Int_t i=0;i<3;i++){
 219     delete fHistVdriftLaserA[i];
 220     delete fHistVdriftLaserC[i];
 221   }
 222   fMapDz->SetOwner();
 223   fMapDz->Delete();
 224   delete fMapDz;
 225   for (Int_t i=0;i<5;i++) {
 226     delete fCosmiMatchingHisto[i];
 227   }
 228 }
 229
 230 void AliTPCcalibTime::ResetCurrent(){
 231   //
 232   // reset current values
 233   //
 234   fDz=0;          // current delta z
 235 }
 236
 237 Bool_t AliTPCcalibTime::IsLaser(AliESDEvent *event){
 238   return ((event->GetFiredTriggerClasses()).Contains("0LSR")==1);
 239 }
 240
 241 void AliTPCcalibTime::Process(AliESDEvent *event){
 242   if(!event) return;
 243   if (event->GetNumberOfTracks()<2) return;
 244   ResetCurrent();
 245
 246 //  if(IsLaser(event))
 247   ProcessLaser (event);
 248 //  else
 249   ProcessCosmic(event);
 250 }
 251
 252 void AliTPCcalibTime::ProcessLaser(AliESDEvent *event){
 253   //
 254   // Fit drift velocity using laser
 255   //
 256   // 0. cuts
 257   const Int_t    kMinTracks     = 20;    // minimal number of laser tracks
 258   const Int_t    kMinTracksSide = 10;    // minimal number of tracks per side
 259   const Float_t  kMaxRMS        = 0.1;   // maximal RMS of tracks
 260   const Float_t  kMaxDeltaZ     = 3.;    // maximal deltaZ A-C side
 261   const Float_t  kMaxDeltaV     = 0.01;  // maximal deltaV A-C side
 262   const Float_t  kMaxDeltaY     = 2.;    // maximal deltaY A-C side
 263   /*
 264     TCut cutRMS("sqrt(laserA.fElements[4])<0.1&&sqrt(laserC.fElements[4])<0.1");
 265     TCut cutZ("abs(laserA.fElements[0]-laserC.fElements[0])<3");
 266     TCut cutV("abs(laserA.fElements[1]-laserC.fElements[1])<0.01");
 267     TCut cutY("abs(laserA.fElements[2]-laserC.fElements[2])<2");
 268     TCut cutAll = cutRMS+cutZ+cutV+cutY;
 269   */
 270   if (event->GetNumberOfTracks()<kMinTracks) return;
 271   //
 272   if(!fLaser) fLaser = new AliTPCcalibLaser("laserTPC","laserTPC",kFALSE);
 273   fLaser->Process(event);
 274   if (fLaser->GetNtracks()<kMinTracks) return;   // small amount of tracks cut
 275   if (fLaser->fFitAside->GetNrows()==0) return;  // no fit A side
 276   if (fLaser->fFitCside->GetNrows()==0) return;  // no fit C side
 277   //
 278   // debug streamer  - activate stream level
 279   // Use it for tuning of the cuts
 280   //
 281   if (fStreamLevel>0){
 282     printf("Trigger: %s\n",event->GetFiredTriggerClasses().Data());
 283
 284     TTreeSRedirector *cstream = GetDebugStreamer();
 285     if (cstream){
 286       TTimeStamp tstamp(fTime);
 287       Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
 288       Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
 289       Double_t ptrelative0   = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
 290       Double_t ptrelative1   = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
 291       Double_t temp0         = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
 292       Double_t temp1         = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
 293       TVectorD vecGoofie(20);
 294       AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
 295       if (goofieArray){
 296         for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
 297           AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
 298           if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
 299         }
 300       }
 301       (*cstream)<<"laserInfo"<<
 302         "run="<<fRun<<              //  run number
 303         "event="<<fEvent<<          //  event number
 304         "time="<<fTime<<            //  time stamp of event
 305         "trigger="<<fTrigger<<      //  trigger
 306         "mag="<<fMagF<<             //  magnetic field
 307         // Environment values
 308         "press0="<<valuePressure0<<
 309         "press1="<<valuePressure1<<
 310         "pt0="<<ptrelative0<<
 311         "pt1="<<ptrelative1<<
 312         "temp0="<<temp0<<
 313         "temp1="<<temp1<<
 314         "vecGoofie.=<<"<<&vecGoofie<<
 315         //laser
 316         "laserA.="<<fLaser->fFitAside<<
 317         "laserC.="<<fLaser->fFitCside<<
 318         "laserAC.="<<fLaser->fFitACside<<
 319         "trigger="<<event->GetFiredTriggerClasses()<<
 320         "\n";
 321     }
 322   }
 323   //
 324   // Apply custs
 325   //
 326   if ((*fLaser->fFitAside)[3] <kMinTracksSide) return; // enough tracks A side
 327   if ((*fLaser->fFitCside)[3]<kMinTracksSide) return;  // enough tracks C side
 328   //
 329   if (TMath::Abs((*fLaser->fFitAside)[0]-(*fLaser->fFitCside)[0])>kMaxDeltaZ) return;
 330   if (TMath::Abs((*fLaser->fFitAside)[2]-(*fLaser->fFitCside)[2])>kMaxDeltaY) return;
 331   if (TMath::Abs((*fLaser->fFitAside)[1]-(*fLaser->fFitCside)[1])>kMaxDeltaV) return;
 332   if (TMath::Sqrt(TMath::Abs((*fLaser->fFitAside)[4]))>kMaxRMS) return;
 333   if (TMath::Sqrt(TMath::Abs((*fLaser->fFitCside)[4]))>kMaxRMS) return;
 334   //
 335   // fill histos
 336   //
 337   TVectorD vdriftA(5), vdriftC(5),vdriftAC(5);
 338   vdriftA=*(fLaser->fFitAside);
 339   vdriftC=*(fLaser->fFitCside);
 340   vdriftAC=*(fLaser->fFitACside);
 341   Int_t npointsA=0, npointsC=0;
 342   Float_t chi2A=0, chi2C=0;
 343   npointsA= TMath::Nint(vdriftA[3]);
 344   chi2A= vdriftA[4];
 345   npointsC= TMath::Nint(vdriftC[3]);
 346   chi2C= vdriftC[4];
 347
 348   TTimeStamp tstamp(fTime);
 349   Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
 350   Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
 351
 352   Double_t vecVdriftLaserA[4]={fTime,(ptrelative0+ptrelative1)/2.0,1./((*(fLaser->fFitAside))[1])-1,event->GetRunNumber()};
 353   Double_t vecVdriftLaserC[4]={fTime,(ptrelative0+ptrelative1)/2.0,1./((*(fLaser->fFitCside))[1])-1,event->GetRunNumber()};
 354
 355   for (Int_t i=0;i<3;i++){
 356     if (i==0){ //z0 shift
 357       vecVdriftLaserA[3]=(*(fLaser->fFitAside))[0]/250.;
 358       vecVdriftLaserA[3]=(*(fLaser->fFitCside))[0]/250.;
 359     }
 360     if (i==1){ //vdrel shift
 361       vecVdriftLaserA[3]=1./(*(fLaser->fFitAside))[1]-1.;
 362       vecVdriftLaserA[3]=1./(*(fLaser->fFitCside))[1]-1.;
 363     }
 364     if (i==2){ //gy shift - full gy - full drift
 365       vecVdriftLaserA[3]=(*(fLaser->fFitAside))[2]/250.;
 366       vecVdriftLaserA[3]=(*(fLaser->fFitCside))[2]/250.;
 367     }
 368     fHistVdriftLaserA[i]->Fill(vecVdriftLaserA);
 369     fHistVdriftLaserC[i]->Fill(vecVdriftLaserC);
 370   }
 371 }
 372
 373 void AliTPCcalibTime::ProcessCosmic(AliESDEvent *event){
 374   if (!event) {
 375     Printf("ERROR: ESD not available");
 376     return;
 377   }
 378   if (event->GetTimeStamp() == 0 ) {
 379     Printf("no time stamp!");
 380     return;
 381   }
 382
 383   //fd
 384   // Find cosmic pairs
 385   //
 386   // Track0 is choosen in upper TPC part
 387   // Track1 is choosen in lower TPC part
 388   //
 389   Int_t ntracks=event->GetNumberOfTracks();
 390   if (ntracks==0) return;
 391   if (ntracks > fCutTracks) return;
 392
 393   if (GetDebugLevel()>1) printf("Hallo world: Im here\n");
 394   AliESDfriend *ESDfriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
 395
 396   TObjArray  tpcSeeds(ntracks);
 397   Double_t vtxx[3]={0,0,0};
 398   Double_t svtxx[3]={0.000001,0.000001,100.};
 399   AliESDVertex vtx(vtxx,svtxx);
 400   //
 401   // track loop
 402   //
 403   for (Int_t i=0;i<ntracks;++i) {
 404     AliESDtrack *track = event->GetTrack(i);
 405
 406     const AliExternalTrackParam * trackIn = track->GetInnerParam();
 407     const AliExternalTrackParam * trackOut = track->GetOuterParam();
 408     if (!trackIn) continue;
 409     if (!trackOut) continue;
 410
 411     AliESDfriendTrack *friendTrack = ESDfriend->GetTrack(i);
 412     TObject *calibObject;
 413     AliTPCseed *seed = 0;
 414     for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
 415     if (seed) tpcSeeds.AddAt(seed,i);
 416   }
 417   if (ntracks<2) return;
 418   //
 419   // Find pairs
 420   //
 421   for (Int_t i=0;i<ntracks;++i) {
 422     AliESDtrack *track0 = event->GetTrack(i);
 423     // track0 - choosen upper part
 424     if (!track0) continue;
 425     if (!track0->GetOuterParam()) continue;
 426     if (track0->GetOuterParam()->GetAlpha()<0) continue;
 427     Double_t d1[3];
 428     track0->GetDirection(d1);
 429     for (Int_t j=0;j<ntracks;++j) {
 430       if (i==j) continue;
 431       AliESDtrack *track1 = event->GetTrack(j);
 432       //track 1 lower part
 433       if (!track1) continue;
 434       if (!track1->GetOuterParam()) continue;
 435       if (track1->GetOuterParam()->GetAlpha()>0) continue;
 436       //
 437       Double_t d2[3];
 438       track1->GetDirection(d2);
 439
 440       AliTPCseed * seed0 = (AliTPCseed*) tpcSeeds.At(i);
 441       AliTPCseed * seed1 = (AliTPCseed*) tpcSeeds.At(j);
 442       if (! seed0) continue;
 443       if (! seed1) continue;
 444       Float_t dir = (d1[0]*d2[0] + d1[1]*d2[1] + d1[2]*d2[2]);
 445       Float_t dist0  = track0->GetLinearD(0,0);
 446       Float_t dist1  = track1->GetLinearD(0,0);
 447       //
 448       // conservative cuts - convergence to be guarantied
 449       // applying before track propagation
 450       if (TMath::Abs(dist0+dist1)>fCutMaxD) continue;   // distance to the 0,0
 451       if (dir>fCutMinDir) continue;               // direction vector product
 452       Float_t bz = AliTracker::GetBz();
 453       Float_t dvertex0[2];   //distance to 0,0
 454       Float_t dvertex1[2];   //distance to 0,0
 455       track0->GetDZ(0,0,0,bz,dvertex0);
 456       track1->GetDZ(0,0,0,bz,dvertex1);
 457       if (TMath::Abs(dvertex0[1])>250) continue;
 458       if (TMath::Abs(dvertex1[1])>250) continue;
 459       //
 460       //
 461       //
 462       Float_t dmax = TMath::Max(TMath::Abs(dist0),TMath::Abs(dist1));
 463       AliExternalTrackParam param0(*track0);
 464       AliExternalTrackParam param1(*track1);
 465       //
 466       // Propagate using Magnetic field and correct fo material budget
 467       //
 468       AliTracker::PropagateTrackTo(&param0,dmax+1,0.0005,3,kTRUE);
 469       AliTracker::PropagateTrackTo(&param1,dmax+1,0.0005,3,kTRUE);
 470       //
 471       // Propagate rest to the 0,0 DCA - z should be ignored
 472       //
 473       //Bool_t b0 = ;
 474       param0.PropagateToDCA(&vtx,bz,1000);
 475       //Bool_t b1 =
 476       param1.PropagateToDCA(&vtx,bz,1000);
 477       //
 478       param0.GetDZ(0,0,0,bz,dvertex0);
 479       param1.GetDZ(0,0,0,bz,dvertex1);
 480       //
 481       Double_t xyz0[3];//,pxyz0[3];
 482       Double_t xyz1[3];//,pxyz1[3];
 483       param0.GetXYZ(xyz0);
 484       param1.GetXYZ(xyz1);
 485       Bool_t isPair = IsPair(&param0,&param1);
 486
 487       Double_t z0 = track0->GetOuterParam()->GetZ();
 488       Double_t z1 = track1->GetOuterParam()->GetZ();
 489
 490       Double_t z0inner = track0->GetInnerParam()->GetZ();
 491       Double_t z1inner = track1->GetInnerParam()->GetZ();
 492
 493       if (isPair && z0>0 && z0inner>0 && z1<0 && z1inner<0) {
 494         if (track0->GetTPCNcls() > 80) {
 495           fDz = param0.GetZ() - param1.GetZ();
 496           TTimeStamp tstamp(fTime);
 497           Double_t ptrelative0 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
 498           Double_t ptrelative1 = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
 499           Double_t vecVdrift[4]={fTime,(ptrelative0+ptrelative1)/2.0,fDz/500.0,event->GetRunNumber()};
 500           THnSparse* curHist=0;
 501
 502           curHist=(THnSparseF*)(fMapDz->GetValue(event->GetFiredTriggerClasses()));
 503           if(!curHist){
 504             curHist=new THnSparseF(event->GetFiredTriggerClasses(),"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
 505             fMapDz->Add(new TObjString(event->GetFiredTriggerClasses()),curHist);
 506           }
 507           curHist->Fill(vecVdrift);
 508
 509           curHist=(THnSparseF*)(fMapDz->GetValue("all"));
 510           if(!curHist){
 511             curHist=new THnSparseF("all","HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
 512             fMapDz->Add(new TObjString("all"),curHist);
 513           }
 514           curHist->Fill(vecVdrift);
 515         }
 516       }
 517     } // end 2nd order loop
 518   } // end 1st order loop
 519
 520   if (fStreamLevel>0){
 521     TTreeSRedirector *cstream = GetDebugStreamer();
 522     if (cstream){
 523       TTimeStamp tstamp(fTime);
 524       Float_t valuePressure0 = AliTPCcalibDB::GetPressure(tstamp,fRun,0);
 525       Float_t valuePressure1 = AliTPCcalibDB::GetPressure(tstamp,fRun,1);
 526       Double_t ptrelative0   = AliTPCcalibDB::GetPTRelative(tstamp,fRun,0);
 527       Double_t ptrelative1   = AliTPCcalibDB::GetPTRelative(tstamp,fRun,1);
 528       Double_t temp0         = AliTPCcalibDB::GetTemperature(tstamp,fRun,0);
 529       Double_t temp1         = AliTPCcalibDB::GetTemperature(tstamp,fRun,1);
 530       TVectorD vecGoofie(20);
 531       AliDCSSensorArray* goofieArray = AliTPCcalibDB::Instance()->GetGoofieSensors(fRun);
 532       if (goofieArray){
 533         for (Int_t isensor=0; isensor<goofieArray->NumSensors();isensor++){
 534           AliDCSSensor *gsensor = goofieArray->GetSensor(isensor);
 535           if (gsensor) vecGoofie[isensor]=gsensor->GetValue(tstamp);
 536         }
 537       }
 538       (*cstream)<<"timeInfo"<<
 539         "run="<<fRun<<              //  run number
 540         "event="<<fEvent<<          //  event number
 541         "time="<<fTime<<            //  time stamp of event
 542         "trigger="<<fTrigger<<      //  trigger
 543         "mag="<<fMagF<<             //  magnetic field
 544         // Environment values
 545         "press0="<<valuePressure0<<
 546         "press1="<<valuePressure1<<
 547         "pt0="<<ptrelative0<<
 548         "pt1="<<ptrelative1<<
 549         "temp0="<<temp0<<
 550         "temp1="<<temp1<<
 551         "vecGoofie.=<<"<<&vecGoofie<<
 552         //
 553         // accumulated values
 554         //
 555         "fDz="<<fDz<<          //! current delta z
 556         "trigger="<<event->GetFiredTriggerClasses()<<
 557         "\n";
 558     }
 559   }
 560   printf("Trigger: %s\n",event->GetFiredTriggerClasses().Data());
 561 }
 562
 563 void AliTPCcalibTime::Analyze(){}
 564
 565 THnSparse* AliTPCcalibTime::GetHistoDrift(TObjString* name){
 566   return (THnSparseF*)(fMapDz->GetValue(name));
 567 }
 568
 569 THnSparse* AliTPCcalibTime::GetHistoDrift(const char* name){
 570   TObjString* objName=new TObjString(name);
 571   THnSparse* histoDrift=0;
 572   if(objName){
 573     histoDrift=GetHistoDrift(objName);
 574     delete objName;
 575     objName=0;
 576   }
 577   return histoDrift;
 578 }
 579
 580 TMap* AliTPCcalibTime::GetHistoDrift(){
 581   return fMapDz;
 582 }
 583
 584 TGraphErrors* AliTPCcalibTime::GetGraphDrift(TObjString* name){
 585   THnSparse* histoDrift=GetHistoDrift(name);
 586   TGraphErrors*    graphDrift=0;
 587   if(histoDrift) graphDrift=FitSlices(histoDrift,2,0,400,100,0.05,0.95, kTRUE);
 588   return graphDrift;
 589 }
 590
 591 TGraphErrors* AliTPCcalibTime::GetGraphDrift(const char* name){
 592   TObjString* objName=new TObjString(name);
 593   TGraphErrors* graphDrift=0;
 594   if(objName){
 595     graphDrift=GetGraphDrift(objName);
 596     delete objName;
 597     objName=0;
 598   }
 599   return graphDrift;
 600 }
 601
 602 TMap* AliTPCcalibTime::GetGraphDrift(){
 603   TMap* mapGraphDrift=new TMap();
 604   TIterator* iterator=fMapDz->MakeIterator();
 605   iterator->Reset();
 606   TPair* addPair=0;
 607   while((addPair=(TPair*)(fMapDz->FindObject(iterator->Next())))) mapGraphDrift->Add((TObjString*)addPair->Key(), GetGraphDrift((TObjString*)addPair->Key()));
 608   return mapGraphDrift;
 609 }
 610
 611 TGraph* AliTPCcalibTime::GetFitDrift(TObjString* name){
 612   TGraphErrors* graphDrift=GetGraphDrift(name);
 613   TGraph* fitDrift=0;
 614   if(graphDrift && graphDrift->GetN()){
 615     AliSplineFit fit;
 616     fit.SetGraph(graphDrift);
 617     fit.SetMinPoints(graphDrift->GetN()+1);
 618     fit.InitKnots(graphDrift,2,0,0.001);
 619     fit.SplineFit(0);
 620     fitDrift=fit.MakeGraph(graphDrift->GetX()[0],graphDrift->GetX()[graphDrift->GetN()-1],50000,0);
 621     delete graphDrift;
 622     graphDrift=0;
 623   }
 624   return fitDrift;
 625 }
 626
 627 TGraph* AliTPCcalibTime::GetFitDrift(const char* name){
 628   TObjString* objName=new TObjString(name);
 629   TGraph* fitDrift=0;
 630   if(objName){
 631     fitDrift=GetFitDrift(objName);
 632     delete objName;
 633     objName=0;
 634   }
 635   return fitDrift;
 636 }
 637
 638 TMap* AliTPCcalibTime::GetFitDrift(){
 639   TMap* mapFitDrift=new TMap();
 640   TIterator* iterator = fMapDz->MakeIterator();
 641   iterator->Reset();
 642   TPair* addPair=0;
 643   while((addPair=(TPair*)(fMapDz->FindObject(iterator->Next())))) mapFitDrift->Add((TObjString*)addPair->Key(), GetFitDrift((TObjString*)addPair->Key()));
 644   return mapFitDrift;
 645 }
 646
 647 Long64_t AliTPCcalibTime::Merge(TCollection *li) {
 648
 649   TIterator* iter = li->MakeIterator();
 650   AliTPCcalibTime* cal = 0;
 651
 652   while ((cal = (AliTPCcalibTime*)iter->Next())) {
 653     if (!cal->InheritsFrom(AliTPCcalibTime::Class())) {
 654       Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
 655       return -1;
 656     }
 657     for (Int_t imeas=0; imeas<3; imeas++){
 658       if (cal->GetHistVdriftLaserA(imeas) && cal->GetHistVdriftLaserA(imeas)){
 659         fHistVdriftLaserA[imeas]->Add(cal->GetHistVdriftLaserA(imeas));
 660         fHistVdriftLaserC[imeas]->Add(cal->GetHistVdriftLaserC(imeas));
 661       }
 662     }
 663     TMap * addMap=cal->GetHistoDrift();
 664     if(!addMap) return 0;
 665     TIterator* iterator = addMap->MakeIterator();
 666     iterator->Reset();
 667     TPair* addPair=0;
 668     while((addPair=(TPair *)(addMap->FindObject(iterator->Next())))){
 669       THnSparse* addHist=dynamic_cast<THnSparseF*>(addPair->Value());
 670       if (!addHist) continue;
 671       addHist->Print();
 672       THnSparse* localHist=dynamic_cast<THnSparseF*>(fMapDz->GetValue(addHist->GetName()));
 673       if(!localHist){
 674         localHist=new THnSparseF(addHist->GetName(),"HistVdrift;time;p/T ratio;Vdrift;run",4,fBinsVdrift,fXminVdrift,fXmaxVdrift);
 675         fMapDz->Add(new TObjString(addHist->GetName()),localHist);
 676       }
 677       localHist->Add(addHist);
 678     }
 679     if (cal->GetCosmiMatchingHisto(0)) fCosmiMatchingHisto[0]->Add(cal->GetCosmiMatchingHisto(0));
 680     if (cal->GetCosmiMatchingHisto(1)) fCosmiMatchingHisto[1]->Add(cal->GetCosmiMatchingHisto(1));
 681     if (cal->GetCosmiMatchingHisto(2)) fCosmiMatchingHisto[2]->Add(cal->GetCosmiMatchingHisto(2));
 682     if (cal->GetCosmiMatchingHisto(3)) fCosmiMatchingHisto[3]->Add(cal->GetCosmiMatchingHisto(3));
 683     if (cal->GetCosmiMatchingHisto(4)) fCosmiMatchingHisto[4]->Add(cal->GetCosmiMatchingHisto(4));
 684   }
 685   return 0;
 686 }
 687
 688 Bool_t  AliTPCcalibTime::IsPair(AliExternalTrackParam *tr0, AliExternalTrackParam *tr1){
 689   /*
 690   // 0. Same direction - OPOSITE  - cutDir +cutT
 691   TCut cutDir("cutDir","dir<-0.99")
 692   // 1.
 693   TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.03")
 694   //
 695   // 2. The same rphi
 696   TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<5")
 697   //
 698   TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<1&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10");
 699   // 1/Pt diff cut
 700   */
 701   const Double_t *p0 = tr0->GetParameter();
 702   const Double_t *p1 = tr1->GetParameter();
 703   fCosmiMatchingHisto[0]->Fill(p0[0]+p1[0]);
 704   fCosmiMatchingHisto[1]->Fill(p0[1]-p1[1]);
 705   fCosmiMatchingHisto[2]->Fill(tr0->GetAlpha()+tr1->GetAlpha());
 706   fCosmiMatchingHisto[3]->Fill(p0[3]+p1[3]);
 707   fCosmiMatchingHisto[4]->Fill(p0[4]+p1[4]);
 708
 709   if (TMath::Abs(p0[3]+p1[3])>fCutTheta) return kFALSE;
 710   if (TMath::Abs(p0[0]+p1[0])>fCutMaxD)  return kFALSE;
 711   if (TMath::Abs(p0[1]-p1[1])>fCutMaxDz)  return kFALSE;
 712   Double_t d0[3], d1[3];
 713   tr0->GetDirection(d0);
 714   tr1->GetDirection(d1);
 715   if (d0[0]*d1[0] + d0[1]*d1[1] + d0[2]*d1[2] >fCutMinDir) return kFALSE;
 716
 717   return kTRUE;
 718 }
 719