HLT/QA/tasks/AliAnalysisTaskHLTCentralBarrel.cxx

   1 // $Id$
   2 //**************************************************************************
   3 //* This file is property of and copyright by the ALICE HLT Project        *
   4 //* ALICE Experiment at CERN, All rights reserved.                         *
   5 //*                                                                        *
   6 //* Primary Authors: Kalliopi Kanaki <Kalliopi.Kanaki@ift.uib.no>          *
   7 //*                  for The ALICE HLT Project.                            *
   8 //*                                                                        *
   9 //* Permission to use, copy, modify and distribute this software and its   *
  10 //* documentation strictly for non-commercial purposes is hereby granted   *
  11 //* without fee, provided that the above copyright notice appears in all   *
  12 //* copies and that both the copyright notice and this permission notice   *
  13 //* appear in the supporting documentation. The authors make no claims     *
  14 //* about the suitability of this software for any purpose. It is          *
  15 //* provided "as is" without express or implied warranty.                  *
  16 //**************************************************************************
  17
  18
  19 /** @file  AliAnalysisTaskHLTCentralBarrel.cxx
  20     @author Per Ivar Lønne, Hege Erdal, Kalliopi Kanaki
  21     @date
  22     @brief An analysis task containing
  23     loops over HLT and offline ESD trees for comparing
  24     event and track properties
  25 */
  26
  27 #include <iostream>
  28
  29 #include "AliAnalysisTaskHLTCentralBarrel.h"
  30 #include "AliESDEvent.h"
  31 #include "AliESDtrack.h"
  32 #include "AliESDInputHandler.h"
  33 #include "AliTracker.h"
  34 #include "AliESDVZERO.h"
  35 #include "AliHLTGlobalTriggerDecision.h"
  36
  37 #include "TAxis.h"
  38 #include "TString.h"
  39 #include "TText.h"
  40 #include "TTimeStamp.h"
  41
  42 ClassImp(AliAnalysisTaskHLTCentralBarrel)
  43 //_______________________________________________________________________________________________//
  44
  45 AliAnalysisTaskHLTCentralBarrel::AliAnalysisTaskHLTCentralBarrel()
  46 :AliAnalysisTaskSE()
  47   ,fUseHLTTrigger(kFALSE)
  48   ,fOutputList(0)
  49   ,fEventOFF(0)
  50   ,fEventHLT(0)
  51   ,fTrackOFF(0)
  52   ,fTrackHLT(0)
  53   ,fTextBox(0)
  54 {
  55   // Constructor
  56   // Define input and output slots here
  57   // Input slot #0 works with a TChain
  58   // DefineInput(0, TChain::Class());
  59   // Output slot #0 writes into a TH1 container
  60
  61   //DefineOutput(1, TList::Class());
  62 }
  63
  64 AliAnalysisTaskHLTCentralBarrel::AliAnalysisTaskHLTCentralBarrel(const char *name)
  65 :AliAnalysisTaskSE(name)
  66   ,fUseHLTTrigger(kFALSE)
  67   ,fOutputList(0)
  68   ,fEventOFF(0)
  69   ,fEventHLT(0)
  70   ,fTrackOFF(0)
  71   ,fTrackHLT(0)
  72   ,fTextBox(0)
  73 {
  74   // Constructor
  75   // Define input and output slots here
  76   // Input slot #0 works with a TChain
  77   //DefineInput(0, TChain::Class());
  78   // Output slot #0 writes into a TH1 container
  79
  80   DefineOutput(1, TList::Class());
  81 }
  82
  83 AliAnalysisTaskHLTCentralBarrel::~AliAnalysisTaskHLTCentralBarrel(){
  84   // destructor
  85 }
  86
  87 void AliAnalysisTaskHLTCentralBarrel::UserCreateOutputObjects(){
  88   // Create THnSpare objects
  89
  90   OpenFile(1);
  91   fOutputList = new TList();
  92   fOutputList->SetName(GetName());
  93
  94   static const int sizeEvent = 6;
  95
  96   int    binsEvent[sizeEvent] = { 200, 200, 250,  2000,  2000, 2 };
  97   double minEvent [sizeEvent] = {  -2,  -2, -30,     0,     0, 0 };
  98   double maxEvent [sizeEvent] = {   2,   2,  30, 20000, 20000, 2 };
  99
 100   fEventHLT = CreateEventTHnSparse("fEventHLT",sizeEvent,binsEvent,minEvent,maxEvent);
 101   fEventOFF = CreateEventTHnSparse("fEventOFF",sizeEvent,binsEvent,minEvent,maxEvent);
 102
 103   static const int sizeTrack = 13;
 104
 105   Int_t    binsTrack[sizeTrack] = {1500, 200, 200, 200, 200,  400,  400,    3,  400,  400, 6,  2000, 2 };
 106   Double_t minTrack [sizeTrack] = {   0,   0,  -1,  -3,  -1, -100, -100, -1.5, -100, -100, 0,     0, 0 };
 107   Double_t maxTrack [sizeTrack] = { 150, 200,   4,   3,   7,  100,  100,  1.5,  100,  100, 6, 20000, 2 };
 108
 109   fTrackHLT = CreateTrackTHnSparse("fTrackHLT",sizeTrack,binsTrack,minTrack,maxTrack);
 110   fTrackOFF = CreateTrackTHnSparse("fTrackOFF",sizeTrack,binsTrack,minTrack,maxTrack);
 111
 112   fTextBox = new TText();
 113   fOutputList->Add(fEventOFF);
 114   fOutputList->Add(fEventHLT);
 115   fOutputList->Add(fTrackOFF);
 116   fOutputList->Add(fTrackHLT);
 117   fOutputList->Add(fTextBox);
 118 }
 119
 120 void AliAnalysisTaskHLTCentralBarrel::NotifyRun(){
 121
 122   AliESDEvent *esdOFF = dynamic_cast<AliESDEvent*>(InputEvent());
 123   TTimeStamp* timestamp = new TTimeStamp(esdOFF->GetTimeStamp());
 124   fTextBox->SetName("text");
 125
 126   TString s = Form("Run %d, Date %d", esdOFF->GetRunNumber(), timestamp->GetDate());
 127   fTextBox->SetTitle(s);
 128 }
 129
 130 void AliAnalysisTaskHLTCentralBarrel::UserExec(Option_t *){
 131 // see header for documentation
 132
 133   AliESDEvent *esdOFF = dynamic_cast<AliESDEvent*>(InputEvent());
 134   if(!esdOFF){
 135      printf("ERROR: offline ESD is not available.\n");
 136      return;
 137   }
 138
 139   AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*>(fInputHandler);
 140   AliESDEvent *esdHLT = NULL;
 141
 142   if(esdH) esdHLT = esdH->GetHLTEvent();
 143   if(!esdHLT){
 144      printf("ERROR: HLT ESD is not available.\n");
 145      return;
 146   }
 147
 148   if(fUseHLTTrigger && !((AliHLTGlobalTriggerDecision*)esdHLT->GetHLTTriggerDecision())->Result()) return;
 149
 150   if(fUseCentrality){
 151     Int_t centbin = CalculateCentrality(esdHLT);
 152     Printf("Centrality bin = %d", centbin);
 153   }
 154
 155
 156
 157   //============================ OFFLINE =============================//
 158
 159   const AliESDVertex *vertOFF = esdOFF->GetPrimaryVertexTracks();
 160
 161   Double_t bfield = esdOFF->GetMagneticField();
 162   Int_t nr_tracksOFF = 0;
 163
 164   if(esdOFF->GetEventSpecie()==16) return;
 165
 166   for(Int_t i=0; i<esdOFF->GetNumberOfTracks(); i++){
 167       AliESDtrack *esdTrackOFF = esdOFF->GetTrack(i);
 168       if (!esdTrackOFF) continue;
 169       if(!(esdTrackOFF->GetStatus()&AliESDtrack::kTPCin)) continue;
 170       nr_tracksOFF++;
 171   }
 172
 173   for(Int_t i=0; i<esdOFF->GetNumberOfTracks(); i++){
 174
 175       AliESDtrack *esdTrackOFF = esdOFF->GetTrack(i);
 176       if (!esdTrackOFF) continue;
 177       if(!(esdTrackOFF->GetStatus()&AliESDtrack::kTPCin)) continue;
 178
 179       //DCA calculations(from offline)
 180       Double_t x[3];
 181       Double_t b[3];
 182       esdTrackOFF->GetXYZ(x);
 183
 184       AliTracker::GetBxByBz(x,b);
 185       Bool_t isOK = esdTrackOFF->RelateToVertexTPCBxByBz(vertOFF, b, kVeryBig);
 186       if(!isOK) continue;
 187
 188       const AliExternalTrackParam *track = esdTrackOFF->GetTPCInnerParam();
 189       if(!track) continue;
 190
 191       Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
 192       if(vertOFF->GetX()==0 && vertOFF->GetY()==0 && vertOFF->GetZ()==0 ){
 193          dca[0]=-99;
 194          dca[1]=-99;
 195       }
 196       else esdTrackOFF->GetImpactParametersTPC(dca,cov);
 197
 198       Float_t DCAr =-99, DCAz = -99.;
 199       Double_t trackOFF[] = {
 200                                 TMath::Abs(esdTrackOFF->Pt())
 201                                ,esdTrackOFF->GetTPCNcls()
 202                                ,esdTrackOFF->Theta()
 203                                ,esdTrackOFF->Eta()
 204                                ,esdTrackOFF->Phi()
 205                                ,dca[0]
 206                                ,dca[1]
 207                                ,esdTrackOFF->Charge()
 208                                ,DCAr
 209                                ,DCAz
 210                                ,esdTrackOFF->GetNcls(0)
 211                                ,nr_tracksOFF
 212                                ,vertOFF->GetStatus()
 213                             };
 214       fTrackOFF->Fill(trackOFF);
 215     }
 216
 217     Double_t eventOFF[] = { vertOFF->GetX(), vertOFF->GetY(), vertOFF->GetZ(), vertOFF->GetNContributors(), nr_tracksOFF, vertOFF->GetStatus()};
 218     fEventOFF->Fill(eventOFF);
 219
 220
 221
 222   //======================================== HLT ==========================================//
 223
 224   Int_t nr_tracksHLT = 0;
 225   if(esdHLT->GetEventSpecie()==16) return;
 226   const AliESDVertex *vertHLT = esdHLT->GetPrimaryVertexTracks();
 227
 228   for(Int_t i=0; i<esdHLT->GetNumberOfTracks(); i++){
 229       AliESDtrack *esdTrackHLT = esdHLT->GetTrack(i);
 230       if (!esdTrackHLT) continue;
 231       if(!(esdTrackHLT->GetStatus()&AliESDtrack::kTPCin)) continue;
 232       nr_tracksHLT++;
 233   }
 234
 235   Double_t eventHLT[] = { vertHLT->GetX(), vertHLT->GetY(), vertHLT->GetZ(), vertHLT->GetNContributors(), nr_tracksHLT, vertHLT->GetStatus()};
 236   fEventHLT->Fill(eventHLT);
 237
 238   for(Int_t i=0; i<esdHLT->GetNumberOfTracks(); i++){
 239
 240       AliESDtrack *esdTrackHLT = esdHLT->GetTrack(i);
 241       if(!esdTrackHLT) continue;
 242       if(!(esdTrackHLT->GetStatus()&AliESDtrack::kTPCin)) continue;
 243
 244       //DCA calculations
 245       Float_t DCAr=-99;
 246       Float_t DCAz=-99;
 247       Float_t dca[2];
 248       if(vertHLT->GetX()==0 && vertHLT->GetY()==0 && vertHLT->GetZ() ==0 ){
 249         dca[0]=-99;
 250         dca[1]=-99;
 251       }
 252       else{
 253         //Calculating DCA "old" fashion
 254         esdTrackHLT->GetDZ(esdHLT->GetPrimaryVertex()->GetXv(), esdHLT->GetPrimaryVertex()->GetYv(), esdHLT->GetPrimaryVertex()->GetZv(), bfield, dca);
 255         // plotting the DCA calculated by Sergey
 256         esdTrackHLT->GetImpactParametersTPC(DCAr,DCAz);
 257       }
 258
 259       Double_t trackHLT[] = {
 260                                TMath::Abs(esdTrackHLT->Pt())
 261                               ,esdTrackHLT->GetTPCNcls()
 262                               ,esdTrackHLT->Theta()
 263                               ,esdTrackHLT->Eta()
 264                               ,esdTrackHLT->Phi()
 265                               ,DCAr
 266                               ,DCAz
 267                               ,esdTrackHLT->Charge()
 268                               ,dca[0]
 269                               ,dca[1]
 270                               ,esdTrackHLT->GetNcls(0)
 271                               ,nr_tracksHLT
 272                               ,vertHLT->GetStatus()
 273                             };
 274       fTrackHLT->Fill(trackHLT);
 275   }
 276   // Post output data.
 277   PostData(1, fOutputList);
 278 }
 279
 280 void AliAnalysisTaskHLTCentralBarrel::Terminate(Option_t *){
 281   // see header file of AliAnalysisTask for documentation
 282 }
 283
 284 Int_t AliAnalysisTaskHLTCentralBarrel::CalculateCentrality(AliESDEvent* esd){
 285 //see header for documentation
 286
 287   Int_t centrality =-1;
 288
 289   AliESDVZERO* esdV0 = esd->GetVZEROData();
 290   //AliESDZDC* esdZDC = esd->GetZDCData();
 291   //Int_t partZDC = esdZDC->GetZDCParticipants();
 292
 293   Float_t multV0 = esdV0->GetMTotV0A() + esdV0->GetMTotV0C();
 294
 295   if (      multV0 >=    0.  && multV0 <=   124.5 ) centrality = 90;
 296   else if ( multV0 >   124.5 && multV0 <=   274.5 ) centrality = 80;
 297   else if ( multV0 >   274.5 && multV0 <=   574.5 ) centrality = 70;
 298   else if ( multV0 >   574.5 && multV0 <=  1224.5 ) centrality = 60;
 299   else if ( multV0 >  1224.5 && multV0 <=  2174.5 ) centrality = 50;
 300   else if ( multV0 >  2174.5 && multV0 <=  3624.5 ) centrality = 40;
 301   else if ( multV0 >  3624.5 && multV0 <=  5574.5 ) centrality = 30;
 302   else if ( multV0 >  5574.5 && multV0 <=  8274.5 ) centrality = 20;
 303   else if ( multV0 >  8274.5 && multV0 <= 12024.5 ) centrality = 10;
 304   else if ( multV0 > 12024.5 && multV0 <= 14674.5 ) centrality = 5;
 305   else if ( multV0 > 14674.5 && multV0 <= 19449.5 ) centrality = 0;
 306
 307   return centrality;
 308 }
 309
 310 THnSparseF* AliAnalysisTaskHLTCentralBarrel::CreateEventTHnSparse(const char* name, Int_t size, const Int_t* bins, Double_t* min, Double_t* max){
 311 //see header for documentation
 312
 313   THnSparseF *thn = new THnSparseF(name,"",size,bins,min,max);
 314   thn->GetAxis(0)->SetTitle("x (cm)");
 315   thn->GetAxis(1)->SetTitle("y (cm)");
 316   thn->GetAxis(2)->SetTitle("z (cm)");
 317   thn->GetAxis(3)->SetTitle("number of contributors");
 318   thn->GetAxis(4)->SetTitle("track multiplicity");
 319   thn->GetAxis(5)->SetTitle("vertex status");
 320   return thn;
 321 }
 322
 323 THnSparseF* AliAnalysisTaskHLTCentralBarrel::CreateTrackTHnSparse(const char* name, Int_t size, const Int_t* bins, Double_t* min, Double_t* max){
 324 //see header for documentation
 325
 326   THnSparseF *thn = new THnSparseF(name,"",size,bins,min,max);
 327   thn->GetAxis(0)->SetTitle("transverse momentum");
 328   thn->GetAxis(1)->SetTitle("TPC clusters per track");
 329   thn->GetAxis(2)->SetTitle("theta");
 330   thn->GetAxis(3)->SetTitle("eta");
 331   thn->GetAxis(4)->SetTitle("phi");
 332   thn->GetAxis(5)->SetTitle("DCAr");
 333   thn->GetAxis(6)->SetTitle("DCAz");
 334   thn->GetAxis(7)->SetTitle("charge");
 335   thn->GetAxis(8)->SetTitle("DCArSG");
 336   thn->GetAxis(9)->SetTitle("DCAzSG");
 337   thn->GetAxis(10)->SetTitle("ITS clusters per track");
 338   return thn;
 339 }
 340
 341 // void AliAnalysisTaskHLTCentralBarrel::Fill(AliESDevent *esd, THnSparseF* thn){
 342 // //see header for documentation
 343 //
 344 //   int nTracks = 0;
 345 //
 346 //   for(int i=0; i<esdHLT->GetNumberOfTracks(); i++){
 347 //
 348 //       AliESDtrack *esdTrack = esd->GetTrack(i);
 349 //       if(!esdTrack) continue;
 350 //       if(!(esdTrack->GetStatus()&AliESDtrack::kTPCin)) continue;
 351 //
 352 //       //DCA calculations
 353 //       Float_t DCAr=-99;
 354 //       Float_t DCAz=-99;
 355 //       Float_t dca[2];
 356 //       if(vertHLT->GetX()==0 && vertHLT->GetY()==0 && vertHLT->GetZ() ==0 ){
 357 //      dca[0]=-99;
 358 //      dca[1]=-99;
 359 //       }
 360 //       else{
 361 //      //Calculating DCA "old" fashion
 362 //      esdTrackHLT->GetDZ(esdHLT->GetPrimaryVertex()->GetXv(), esdHLT->GetPrimaryVertex()->GetYv(), esdHLT->GetPrimaryVertex()->GetZv(), bfield, dca);
 363 //      // plotting the DCA calculated by Sergey
 364 //      esdTrackHLT->GetImpactParametersTPC(DCAr,DCAz);
 365 //       }
 366 //
 367 //       Double_t trackHLT[] = {
 368 //                              TMath::Abs(esdTrackHLT->Pt())
 369 //                             ,esdTrackHLT->GetTPCNcls()
 370 //                             ,esdTrackHLT->Theta()
 371 //                             ,esdTrackHLT->Eta()
 372 //                             ,esdTrackHLT->Phi()
 373 //                             ,DCAr
 374 //                             ,DCAz
 375 //                             ,esdTrackHLT->Charge()
 376 //                             ,dca[0]
 377 //                             ,dca[1]
 378 //                             ,esdTrackHLT->GetStatus()
 379 //                             ,esdTrackHLT->GetNcls(0)
 380 //                           };
 381 //       fTrackHLT->Fill(trackHLT);
 382 //       nr_tracksHLT++;
 383 //   }
 384 //   Double_t eventHLT[] = { vertHLT->GetX(), vertHLT->GetY(), vertHLT->GetZ(), vertHLT->GetNContributors(), nr_tracksHLT, vertHLT->GetStatus()};
 385 //   fEventHLT->Fill(eventHLT);
 386 // }