PWG/CaloTrackCorrBase/AliAnalysisTaskCounter.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 /* $Id: $ */
  16
  17 //_________________________________________________________________________
  18 // Count events with different selections
  19 //
  20 // It produces a histogram with the number of events with 9 bins:
  21 // 0: all events (that passed the physics selection if it was on)
  22 // 1: same but cross check that event pointer did exist
  23 // 2: passes vertex cut
  24 // 3: passes track number cut, tracks for eta < 0.8
  25 // 4: 3 && 2
  26 // 5: pass VAND
  27 // 6: 5 && 2
  28 // 7: 5 && 3
  29 // 8: 5 && 3 && 2
  30 // 9: not pileup from SPD
  31 // 10: Good vertex
  32 // 11: 10 && 5
  33 // 12: 10 && 3
  34 // 13: 10 && 2
  35 // 14: 10 && 2 && 3 && 5
  36 // 15: 10 && 9
  37 // 16: 9  && 5
  38 //
  39 // Author: Gustavo Conesa Balbastre (LPSC)
  40 //
  41 //_________________________________________________________________________
  42
  43 #include "TH2F.h"
  44 #include "AliAODHeader.h"
  45 #include "AliTriggerAnalysis.h"
  46 #include "AliESDEvent.h"
  47 #include "AliAODEvent.h"
  48 #include "AliESDtrackCuts.h"
  49 #include "AliAnalysisManager.h"
  50 #include "AliInputEventHandler.h"
  51
  52 #include "AliAnalysisTaskCounter.h"
  53 ClassImp(AliAnalysisTaskCounter)
  54
  55 //________________________________________________________________________
  56 AliAnalysisTaskCounter::AliAnalysisTaskCounter(const char *name)
  57 : AliAnalysisTaskSE(name),
  58   fAcceptFastCluster(kTRUE),
  59   fZVertexCut(10.),
  60   fTrackMultEtaCut(0.8),
  61   fCaloFilterPatch(kFALSE),
  62   fOutputContainer(0x0),
  63   fESDtrackCuts(AliESDtrackCuts::GetStandardITSTPCTrackCuts2010()),
  64   fTriggerAnalysis (new AliTriggerAnalysis),
  65   fhNEvents(0),
  66   fhXVertex(0),    fhYVertex(0),    fhZVertex(0),
  67   fhXGoodVertex(0),fhYGoodVertex(0),fhZGoodVertex(0)
  68 {
  69   //ctor
  70   DefineOutput(1, TList::Class());
  71 }
  72
  73 //________________________________________________________________________
  74 AliAnalysisTaskCounter::AliAnalysisTaskCounter()
  75   : AliAnalysisTaskSE("DefaultAnalysis_AliAnalysisTaskCounter"),
  76     fAcceptFastCluster(kTRUE),
  77     fZVertexCut(10.),
  78     fTrackMultEtaCut(0.8),
  79     fCaloFilterPatch(kFALSE),
  80     fOutputContainer(0x0),
  81     fESDtrackCuts(AliESDtrackCuts::GetStandardITSTPCTrackCuts2010()),
  82     fTriggerAnalysis (new AliTriggerAnalysis),
  83     fhNEvents(0),
  84     fhXVertex(0),    fhYVertex(0),    fhZVertex(0),
  85     fhXGoodVertex(0),fhYGoodVertex(0),fhZGoodVertex(0)
  86 {
  87   // ctor
  88   DefineOutput(1, TList::Class());
  89 }
  90
  91 //__________________________________________________
  92 AliAnalysisTaskCounter::~AliAnalysisTaskCounter()
  93 {
  94   //Destructor
  95   if(fOutputContainer){
  96     fOutputContainer->Delete() ;
  97     delete fOutputContainer ;
  98   }
  99
 100   if(fESDtrackCuts)    delete fESDtrackCuts;
 101   if(fTriggerAnalysis) delete fTriggerAnalysis;
 102
 103 }
 104
 105
 106 //-------------------------------------------------------------------
 107 void AliAnalysisTaskCounter::UserCreateOutputObjects()
 108 {
 109   // Init histograms
 110
 111   fOutputContainer = new TList();
 112
 113   fhZVertex     = new TH1F("hZVertex", " Z vertex distribution"   , 200 , -50 , 50  ) ;
 114   fhZVertex->SetXTitle("v_{z} (cm)");
 115   fOutputContainer->Add(fhZVertex);
 116
 117   fhZGoodVertex     = new TH1F("hZGoodVertex", " Good Z vertex distribution"   , 200 , -50 , 50  ) ;
 118   fhZGoodVertex->SetXTitle("v_{z} (cm)");
 119   fOutputContainer->Add(fhZGoodVertex);
 120
 121   fhXVertex     = new TH1F("hXVertex", " X vertex distribution"   , 200 , -2 , 2  ) ;
 122   fhXVertex->SetXTitle("v_{x} (cm)");
 123   fOutputContainer->Add(fhXVertex);
 124
 125   fhXGoodVertex     = new TH1F("hXGoodVertex", " Good X vertex distribution"   , 200 , -2 , 2  ) ;
 126   fhXGoodVertex->SetXTitle("v_{x} (cm)");
 127   fOutputContainer->Add(fhXGoodVertex);
 128
 129   fhYVertex     = new TH1F("hYVertex", " Y vertex distribution"   , 200 , -2 , 2  ) ;
 130   fhYVertex->SetXTitle("v_{y} (cm)");
 131   fOutputContainer->Add(fhYVertex);
 132
 133   fhYGoodVertex     = new TH1F("hYGoodVertex", " Good Y vertex distribution"   , 200 , -2 , 2  ) ;
 134   fhYGoodVertex->SetXTitle("v_{y} (cm)");
 135   fOutputContainer->Add(fhYGoodVertex);
 136
 137
 138   fhNEvents = new TH1I("hNEvents", "Number of analyzed events", 21, 0, 21) ;
 139   fhNEvents->SetXTitle("Selection");
 140   fhNEvents->SetYTitle("# events");
 141   fhNEvents->GetXaxis()->SetBinLabel(1 ,"1  = PS");
 142   fhNEvents->GetXaxis()->SetBinLabel(2 ,"2  = 1  & ESD");
 143   fhNEvents->GetXaxis()->SetBinLabel(3 ,"3  = 2  & |Z|<10");
 144   fhNEvents->GetXaxis()->SetBinLabel(4 ,"4  = 2  & |eta|<0.8");
 145   fhNEvents->GetXaxis()->SetBinLabel(5 ,"5  = 3  & 4");
 146   fhNEvents->GetXaxis()->SetBinLabel(6 ,"6  = 2  & V0AND");
 147   fhNEvents->GetXaxis()->SetBinLabel(7 ,"7  = 3  & 6");
 148   fhNEvents->GetXaxis()->SetBinLabel(8 ,"8  = 4  & 6");
 149   fhNEvents->GetXaxis()->SetBinLabel(9 ,"9  = 5  & 6");
 150   fhNEvents->GetXaxis()->SetBinLabel(10,"10 = 2  & not pileup");
 151   fhNEvents->GetXaxis()->SetBinLabel(11,"11 = 2  & good vertex");
 152   fhNEvents->GetXaxis()->SetBinLabel(12,"12 = 3  & 11");
 153   fhNEvents->GetXaxis()->SetBinLabel(13,"13 = 4  & 11");
 154   fhNEvents->GetXaxis()->SetBinLabel(14,"14 = 6  & 11");
 155   fhNEvents->GetXaxis()->SetBinLabel(15,"15 = 9  & 11");
 156   fhNEvents->GetXaxis()->SetBinLabel(16,"16 = 10 & 11");
 157   fhNEvents->GetXaxis()->SetBinLabel(17,"17 = 6  & 10");
 158   fhNEvents->GetXaxis()->SetBinLabel(17,"17 = 1  & |Z|<50");
 159   fhNEvents->GetXaxis()->SetBinLabel(18,"18 = Reject EMCAL 1");
 160   fhNEvents->GetXaxis()->SetBinLabel(19,"19 = 18 & 2");
 161   fhNEvents->GetXaxis()->SetBinLabel(20,"20 = Reject EMCAL 2");
 162   fhNEvents->GetXaxis()->SetBinLabel(21,"20 = 20 & 2");
 163
 164   fOutputContainer->Add(fhNEvents);
 165
 166   fOutputContainer->SetOwner(kTRUE);
 167
 168   PostData(1,fOutputContainer);
 169
 170 }
 171
 172 //________________________________________________________________________
 173 void AliAnalysisTaskCounter::UserExec(Option_t *)
 174 {
 175   // Main loop
 176   // Called for each event
 177
 178   //printf("___ Event __ %d __\n",(Int_t)Entry());
 179
 180   fhNEvents->Fill(0.5);
 181
 182   AliVEvent * event = InputEvent();
 183   if (!event) {
 184     printf("AliAnalysisTaskCounter::UserExec() - ERROR: event not available \n");
 185     return;
 186   }
 187
 188   AliESDEvent * esdevent = dynamic_cast<AliESDEvent*> (event);
 189   AliAODEvent * aodevent = dynamic_cast<AliAODEvent*> (event);
 190
 191   TString triggerclasses = "";
 192   if(esdevent) triggerclasses = esdevent->GetFiredTriggerClasses();
 193   if(aodevent) triggerclasses = aodevent->GetFiredTriggerClasses();
 194
 195   if (triggerclasses.Contains("FAST") && !triggerclasses.Contains("ALL") && !fAcceptFastCluster) {
 196     //printf("Do not count events from fast cluster, trigger name %s\n",triggerclasses.Data());
 197     return;
 198   }
 199
 200   fhNEvents->Fill(1.5);
 201
 202   //Initialize bools
 203   Bool_t bSelectVZ    = kFALSE;
 204   Bool_t bV0AND       = kFALSE;
 205   Bool_t bPileup      = kFALSE;
 206   Bool_t bGoodV       = kFALSE;
 207   Bool_t bSelectTrack = kFALSE;
 208   Int_t  trackMult    = 0;
 209
 210   //---------------------------------
 211   //Get the primary vertex, cut on Z
 212   //---------------------------------
 213   Double_t v[3];
 214   event->GetPrimaryVertex()->GetXYZ(v) ;
 215   fhXVertex->Fill(v[0]);
 216   fhYVertex->Fill(v[1]);
 217   fhZVertex->Fill(v[2]);
 218
 219   if(TMath::Abs(v[2]) < fZVertexCut) {
 220     bSelectVZ=kTRUE;
 221     fhNEvents->Fill(2.5);
 222   }
 223   //else printf("Vertex out %f \n",v[2]);
 224
 225
 226   //--------------------------------------------------
 227   //Tweak for calorimeter only productions
 228   //--------------------------------------------------
 229   if(fCaloFilterPatch && !esdevent){
 230     if(event->GetNumberOfCaloClusters() > 0) {
 231       AliVCluster * calo = event->GetCaloCluster(0);
 232       if(calo->GetNLabels() == 4){
 233         Int_t * selection = calo->GetLabels();
 234         bPileup   = selection[0];
 235         bGoodV    = selection[1];
 236         bV0AND    = selection[2];
 237         trackMult = selection[3];
 238         //if(selection[0] || selection[1] || selection[2])
 239         //printf(" pu %d, gv %d, v0 %d, track mult %d\n ", selection[0], selection[1], selection[2], selection[3]);
 240         if(trackMult > 0 )
 241           bSelectTrack = kFALSE;
 242       } else {
 243         //First filtered AODs, track multiplicity stored there.
 244         trackMult = (Int_t) ((AliAODHeader*)fInputEvent->GetHeader())->GetCentrality();
 245       }
 246     }else{//at least one cluster
 247         //printf("AliAnalysisTaskCounter::UserExec() - No clusters in event\n");
 248         //Remove events with  vertex (0,0,0), bad vertex reconstruction
 249         if(TMath::Abs(v[0]) < 1.e-6 && TMath::Abs(v[1]) < 1.e-6 && TMath::Abs(v[2]) < 1.e-6) bGoodV = kFALSE;
 250
 251         //First filtered AODs, track multiplicity stored there.
 252         trackMult = (Int_t) ((AliAODHeader*)fInputEvent->GetHeader())->GetCentrality();
 253     }
 254   }
 255   else {
 256     //--------------------------------------------------
 257     //Count tracks, cut on number of tracks in eta < 0.8
 258     //--------------------------------------------------
 259     Int_t nTracks   = event->GetNumberOfTracks() ;
 260     for (Int_t itrack =  0; itrack <  nTracks; itrack++) {////////////// track loop
 261       AliVTrack * track = (AliVTrack*)event->GetTrack(itrack) ; // retrieve track from esd
 262
 263       //Only for ESDs
 264       if(esdevent && !fESDtrackCuts->AcceptTrack((AliESDtrack*)track)) continue;
 265
 266       //Do not count tracks out of acceptance cut
 267       if(TMath::Abs(track->Eta())< fTrackMultEtaCut) trackMult++;
 268     }
 269   }
 270
 271   //printf("AliAnalysisTaskCounter::UserExec() - Track Mult %d \n",trackMult);
 272
 273   //--------------------------------------------------
 274   // At least one track
 275   //--------------------------------------------------
 276   if (trackMult > 0) {
 277     bSelectTrack = kTRUE;
 278                   fhNEvents->Fill(3.5);
 279     if(bSelectVZ) fhNEvents->Fill(4.5);
 280   }
 281
 282   //---------------------------------
 283   // V0AND
 284   //---------------------------------
 285   if(esdevent && !fCaloFilterPatch) bV0AND = fTriggerAnalysis->IsOfflineTriggerFired(esdevent, AliTriggerAnalysis::kV0AND);
 286   //else if(aodevent  && !fCaloFilterPatch) bV0AND = //FIXME FOR AODs
 287
 288   if(bV0AND)
 289   {
 290                                     fhNEvents->Fill(5.5);
 291     if (bSelectVZ)                  fhNEvents->Fill(6.5);
 292     if (bSelectTrack)               fhNEvents->Fill(7.5);
 293     if (bSelectVZ &&  bSelectTrack) fhNEvents->Fill(8.5);
 294   }
 295
 296   //---------------------------------
 297   // Pileup
 298   //---------------------------------
 299   if(!fCaloFilterPatch)
 300     bPileup = event->IsPileupFromSPD(3, 0.8, 3., 2., 5.); //Default values, if not it does not compile
 301   //bPileup = event->IsPileupFromSPD();
 302   if (!bPileup){
 303                 fhNEvents->Fill(9.5);
 304     if(bV0AND)  fhNEvents->Fill(16.5);
 305   }
 306
 307   //---------------------------------
 308   // Good vertex
 309   //---------------------------------
 310   if(esdevent && !fCaloFilterPatch) bGoodV = CheckForPrimaryVertex();
 311   if(bGoodV)
 312   {
 313     fhXGoodVertex->Fill(v[0]);
 314     fhYGoodVertex->Fill(v[1]);
 315     fhZGoodVertex->Fill(v[2]);
 316
 317                      fhNEvents->Fill(10.5);
 318     if(bSelectVZ)    fhNEvents->Fill(11.5);
 319     if(bSelectTrack) fhNEvents->Fill(12.5);
 320     if(bV0AND)       fhNEvents->Fill(13.5);
 321     if(bSelectVZ && bSelectTrack && bV0AND)
 322                      fhNEvents->Fill(14.5);
 323     if(!bPileup)     fhNEvents->Fill(15.5);
 324   }
 325
 326   //printf("AliAnalysisTaskCounter::UserExec() : z vertex %d, good vertex %d, v0and %d, pile up %d, track mult %d\n ", bSelectVZ, bGoodV, bV0AND, bPileup, trackMult);
 327
 328   // Events that could be rejected in EMCAL
 329   // LHC11a, SM4 and some SM3 events cut with this
 330   Bool_t bEMCALRejected = kFALSE;
 331   for (Int_t i = 0; i < InputEvent()->GetNumberOfCaloClusters(); i++)
 332   {
 333     AliVCluster *clus = InputEvent()->GetCaloCluster(i);
 334     if(clus->IsEMCAL()){
 335       if ((clus->E() > 500 && clus->GetNCells() > 200 ) || clus->GetNCells() > 200)  {
 336
 337         //printf("Counter: Reject event with cluster: E %f, ncells %d\n",clus->E(),clus->GetNCells());
 338
 339                          fhNEvents->Fill(17.5);
 340         if(bSelectVZ)    fhNEvents->Fill(18.5);
 341         bEMCALRejected = kTRUE;
 342         break;
 343       }
 344     }
 345   }
 346
 347   //LHC11a, 3 last runs, cut with this
 348   if(!bEMCALRejected){
 349     // Count number of cells in SM3 with energy larger than 0.1, cut on this number
 350     Int_t ncellsSM3 = 0;
 351     Int_t ncellsSM4 = 0;
 352     for(Int_t icell = 0; icell < event->GetEMCALCells()->GetNumberOfCells(); icell++){
 353       if(event->GetEMCALCells()->GetAmplitude(icell) > 0.1 && event->GetEMCALCells()->GetCellNumber(icell)/(24*48)==3) ncellsSM3++;
 354       if(event->GetEMCALCells()->GetAmplitude(icell) > 0.1 && event->GetEMCALCells()->GetCellNumber(icell)/(24*48)==4) ncellsSM4++;
 355     }
 356
 357     Int_t ncellcut = 21;
 358     if(triggerclasses.Contains("EMC")) ncellcut = 35;
 359
 360     if( ncellsSM3 >= ncellcut || ncellsSM4 >= 100 ){
 361
 362       //printf("Counter: reject event with ncells in SM3: ncells %d\n",ncells);
 363
 364                        fhNEvents->Fill(19.5);
 365       if(bSelectVZ)    fhNEvents->Fill(20.5);
 366     }
 367
 368   }
 369
 370   PostData(1,fOutputContainer);
 371
 372 }
 373
 374 //____________________________________________________________________________
 375 Bool_t AliAnalysisTaskCounter::CheckForPrimaryVertex(){
 376   //Check if the vertex was well reconstructed, copy from V0Reader of conversion group
 377   //It only works for ESDs
 378
 379   AliESDEvent * event = dynamic_cast<AliESDEvent*> (InputEvent());
 380   if(!event) return 0;
 381
 382   if(event->GetPrimaryVertexTracks()->GetNContributors() > 0) {
 383     return 1;
 384   }
 385
 386   if(event->GetPrimaryVertexTracks()->GetNContributors() < 1) {
 387     // SPD vertex
 388     if(event->GetPrimaryVertexSPD()->GetNContributors() > 0) {
 389       //cout<<"spd vertex type::"<< fESDEvent->GetPrimaryVertex()->GetName() << endl;
 390       return 1;
 391
 392     }
 393     if(event->GetPrimaryVertexSPD()->GetNContributors() < 1) {
 394       //      cout<<"bad vertex type::"<< fESDEvent->GetPrimaryVertex()->GetName() << endl;
 395       return 0;
 396     }
 397   }
 398   return 0;
 399   //return fInputEvent->GetPrimaryVertex()->GetNContributors()>0;
 400 }
 401
 402
 403
 404 //_____________________________________________________
 405 void AliAnalysisTaskCounter::FinishTaskOutput()
 406 {
 407   // Put in the output some event summary histograms
 408
 409   AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
 410   AliInputEventHandler *inputH = dynamic_cast<AliInputEventHandler*>(am->GetInputEventHandler());
 411   if (!inputH) return;
 412   TH2F *histStat = dynamic_cast<TH2F*>(inputH->GetStatistics());
 413   TH2F *histBin0 = dynamic_cast<TH2F*>(inputH->GetStatistics("BIN0"));
 414
 415   if(histStat)
 416     fOutputContainer->Add(histStat);
 417   else if(DebugLevel() > 1)
 418     printf("AliAnalysisTaskCounter::FinishTaskOutput() - Stat histogram not available check, \n if ESDs, that AliPhysicsSelection was on, \n if AODs, if EventStat_temp.root exists \n");
 419
 420   if(histBin0)
 421     fOutputContainer->Add(histBin0);
 422
 423 }