PWG0/multPbPb/AliAnalysisTaskTriggerStudy.cxx

   1 // AliAnalysisTaskTriggerStudy
   2
   3 // Author: Michele Floris, CERN
   4 // TODO:
   5 // - Add chi2/cluster plot for primary, secondaries and fakes
   6
   7
   8 #include "AliAnalysisTaskTriggerStudy.h"
   9 #include "AliESDInputHandler.h"
  10 #include "AliHistoListWrapper.h"
  11 #include "AliAnalysisManager.h"
  12 #include "AliMCEvent.h"
  13 #include "AliStack.h"
  14 #include "TH1I.h"
  15 #include "TH3D.h"
  16 #include "AliMCParticle.h"
  17 #include "AliGenEventHeader.h"
  18 #include "AliESDCentrality.h"
  19
  20 #include <iostream>
  21 #include "AliTriggerAnalysis.h"
  22 #include "AliMultiplicity.h"
  23 #include "TFile.h"
  24 #include "AliLog.h"
  25 #include "AliESDtrackCuts.h"
  26 #include "AliESDVZERO.h"
  27 #include "TH2F.h"
  28
  29 using namespace std;
  30
  31 ClassImp(AliAnalysisTaskTriggerStudy)
  32
  33 //const char * AliAnalysisTaskTriggerStudy::kVDNames[] = {"C0SM1","C0SM2","C0VBA","C0VBC","C0OM2"};
  34 const char * AliAnalysisTaskTriggerStudy::kVDNames[] = {"C0SM1","C0SM2","C0VBA","C0VBC"};//,"C0OM2"};
  35
  36 AliAnalysisTaskTriggerStudy::AliAnalysisTaskTriggerStudy()
  37 : AliAnalysisTaskSE("TaskTriggerStudy"),
  38   fESD(0),fHistoList(0),fIsMC(0),fTriggerAnalysis(0),fHistoSuffix(""),fNTrackletsCut(1000000),fNTrackletsCutKine(100),fRejectBGWithV0(0)
  39 {
  40   // constructor
  41
  42   DefineOutput(1, AliHistoListWrapper::Class());
  43
  44 }
  45 AliAnalysisTaskTriggerStudy::AliAnalysisTaskTriggerStudy(const char * name)
  46   : AliAnalysisTaskSE(name),
  47     fESD(0),fHistoList(0),fIsMC(0),fTriggerAnalysis(0),fHistoSuffix(""),fNTrackletsCut(1000000),fNTrackletsCutKine(100),fRejectBGWithV0(0)
  48 {
  49   //
  50   // Standard constructur which should be used
  51   //
  52
  53   DefineOutput(1, AliHistoListWrapper::Class());
  54
  55 }
  56
  57 AliAnalysisTaskTriggerStudy::AliAnalysisTaskTriggerStudy(const AliAnalysisTaskTriggerStudy& obj) :
  58   AliAnalysisTaskSE(obj) ,fESD (0), fIsMC(0), fTriggerAnalysis(0),fHistoSuffix(""),fNTrackletsCut(1000000),fNTrackletsCutKine(100),fRejectBGWithV0(0)
  59 {
  60   //copy ctor
  61   fESD = obj.fESD ;
  62   fHistoList = obj.fHistoList;
  63   fTriggerAnalysis = obj.fTriggerAnalysis;
  64   fHistoSuffix = obj.fHistoSuffix;
  65   fNTrackletsCut = obj.fNTrackletsCut;
  66   fNTrackletsCutKine = obj.fNTrackletsCutKine;
  67   fRejectBGWithV0 = obj.fRejectBGWithV0;
  68 }
  69
  70 AliAnalysisTaskTriggerStudy::~AliAnalysisTaskTriggerStudy(){
  71   // destructor
  72
  73   if(!AliAnalysisManager::GetAnalysisManager()->IsProofMode()) {
  74     if(fHistoList) {
  75       delete fHistoList;
  76       fHistoList = 0;
  77     }
  78     if(fTriggerAnalysis) {
  79       delete fTriggerAnalysis;
  80       fHistoList = 0;
  81     }
  82   }
  83   // Histo list should not be destroyed: fListWrapper is owner!
  84
  85 }
  86 void AliAnalysisTaskTriggerStudy::UserCreateOutputObjects()
  87 {
  88   // Called once
  89   fHistoList = new AliHistoListWrapper("histoList","histogram list for trigger studies");
  90   fTriggerAnalysis = new AliTriggerAnalysis();
  91 }
  92
  93
  94 void AliAnalysisTaskTriggerStudy::UserExec(Option_t *)
  95 {
  96   // User code
  97
  98   // FIXME: make sure you have the right cuts here
  99
 100   /* PostData(0) is taken care of by AliAnalysisTaskSE */
 101   PostData(1,fHistoList);
 102
 103   fESD = dynamic_cast<AliESDEvent*>(fInputEvent);
 104   if (strcmp(fESD->ClassName(),"AliESDEvent")) {
 105     AliFatal("Not processing ESDs");
 106   }
 107
 108   // get the multiplicity object
 109   const AliMultiplicity* mult = fESD->GetMultiplicity();
 110   Int_t ntracklets = mult->GetNumberOfTracklets();
 111   // Get Number of tracks
 112   Int_t ntracks    = AliESDtrackCuts::GetReferenceMultiplicity(fESD,kTRUE); // tpc only
 113
 114   // Get V0 Multiplicity
 115   AliESDVZERO* esdV0 = fESD->GetVZEROData();
 116   Float_t multV0A=esdV0->GetMTotV0A();
 117   Float_t multV0C=esdV0->GetMTotV0C();
 118   Float_t multV0 = multV0A+multV0C;
 119
 120   // Get number of clusters in layer 1
 121   Float_t outerLayerSPD = mult->GetNumberOfITSClusters(1);
 122   Float_t innerLayerSPD = mult->GetNumberOfITSClusters(0);
 123   Float_t totalClusSPD = outerLayerSPD+innerLayerSPD;
 124
 125   GetHistoSPD1  ("All", "All events before any selection")->Fill(outerLayerSPD);
 126
 127
 128   // Physics selection
 129   Bool_t isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB);
 130   if(!isSelected) return;
 131   // FIXME trigger classes
 132   // if ( !(fESD->IsTriggerClassFired("CMBS2A-B-NOPF-ALL")|| fESD->IsTriggerClassFired("CMBS2C-B-NOPF-ALL") || fESD->IsTriggerClassFired("CMBAC-B-NOPF-ALL")) ) return;
 133
 134   GetHistoSPD1  ("AllPSNoPrino", "All events after physsel, before Francesco")->Fill(outerLayerSPD);
 135
 136
 137   // Francesco's cuts
 138   const AliESDVertex * vtxESDTPC= fESD->GetPrimaryVertexTPC();
 139   if(vtxESDTPC->GetNContributors()<1) return;
 140   if (vtxESDTPC->GetNContributors()<(-10.+0.25*fESD->GetMultiplicity()->GetNumberOfITSClusters(0)))     return;
 141   const AliESDVertex * vtxESDSPD= fESD->GetPrimaryVertexSPD();
 142   Float_t tpcContr=vtxESDTPC->GetNContributors();
 143
 144
 145
 146   // GetT0 Stuff
 147   const Double32_t *meanT0   = fESD->GetT0TOF();
 148   const Double32_t  meanT0A  = 0.001* meanT0[1];
 149   const Double32_t  meanT0C  = 0.001* meanT0[2];
 150   const Double32_t  meanT0AC = 0.001* meanT0[0];
 151   Double32_t  T0Vertex = fESD->GetT0zVertex();
 152   //  Double32_t  *ampT0 =Esdevent ->GetT0amplitude();
 153
 154 //   cut1yesF = ( (meanC < 95. && meanA < 95.) && (meanC < -2.) ) && francescoscut
 155 // cut1notF = ( (meanC < 95. && meanA < 95.) && (meanC < -2.) ) && ! francescoscut
 156 // cut2 = ( (meanC < 95. && meanA < 95.) && ( (meanC-meanA) <=-0.7) && meanC > -2) )
 157 // cut3 = ( (meanC < 95. && meanA < 95.) && ( (meanC-meanA) < 0.7 && (meanC-meanA) > -0.7 ) )
 158 //   cut4 = ( (meanC < 95. && meanA < 95.) && (meanA < -2.)
 159
 160   Bool_t cut1T0 =  ( (meanT0C < 95. && meanT0A < 95.) && (meanT0C < -2.) );
 161   Bool_t cut2T0 = ( (meanT0C < 95. && meanT0A < 95.) && ( (meanT0C-meanT0A) <=-0.7) && meanT0C > -2) ;
 162
 163   if(ntracklets > fNTrackletsCut) return;
 164
 165   // Reset histo suffix and fill reference histograms without any suffix
 166   fHistoSuffix = "";
 167
 168   // Fast or in the outer layer
 169   Int_t nFastOrOnline  = fTriggerAnalysis->SPDFiredChips(fESD, 1, 0, 2); // offline
 170   Int_t nFastOrOffline = fTriggerAnalysis->SPDFiredChips(fESD, 0, 0, 2); // online
 171
 172   Bool_t c0sm1 = nFastOrOffline >= 1;
 173   Bool_t c0sm2 = nFastOrOffline >= 2;
 174   Bool_t c0sm3 = nFastOrOffline >= 3;
 175   Bool_t c0sm4 = nFastOrOffline >= 4;
 176   Bool_t c0sm5 = nFastOrOffline >= 5;
 177
 178   // V0 triggers
 179   Bool_t c0v0A       = fTriggerAnalysis->IsOfflineTriggerFired(fESD, AliTriggerAnalysis::kV0A);
 180   Bool_t c0v0C       = fTriggerAnalysis->IsOfflineTriggerFired(fESD, AliTriggerAnalysis::kV0C);
 181   Bool_t v0AHW     = (fTriggerAnalysis->V0Trigger(fESD, AliTriggerAnalysis::kASide, kTRUE) == AliTriggerAnalysis::kV0BB);// should replay hw trigger
 182   Bool_t v0CHW     = (fTriggerAnalysis->V0Trigger(fESD, AliTriggerAnalysis::kCSide, kTRUE) == AliTriggerAnalysis::kV0BB);// should replay hw trigger
 183
 184   Bool_t v0ABG = fTriggerAnalysis->IsOfflineTriggerFired(fESD, AliTriggerAnalysis::kV0ABG);
 185   Bool_t v0CBG = fTriggerAnalysis->IsOfflineTriggerFired(fESD, AliTriggerAnalysis::kV0CBG);
 186   Bool_t v0BG  = v0ABG || v0CBG;
 187
 188
 189   // TOF triggers
 190   // FIXME: move to triggeranalysis?
 191   AliESDHeader*h = fESD->GetHeader(); // taken the header from AliESDEvent
 192   Bool_t c0OM2 = h->IsTriggerInputFired("0OM2"); // thr >= 2 (input 19)
 193   Bool_t c0OM3 = h->IsTriggerInputFired("0OM3"); // thr >= 3 (input 20)
 194
 195   // ZDC triggers
 196   Bool_t zdcA   = kFALSE;
 197   Bool_t zdcC   = kFALSE;
 198   Bool_t zdcBar = kFALSE;
 199
 200   if (!fIsMC) {
 201     // If it's data, we use the TDCs
 202     zdcA   = fTriggerAnalysis->ZDCTDCTrigger(fESD, AliTriggerAnalysis::kASide, kTRUE, kFALSE) ;
 203     zdcC   = fTriggerAnalysis->ZDCTDCTrigger(fESD, AliTriggerAnalysis::kCSide, kTRUE, kFALSE) ;
 204     zdcBar = fTriggerAnalysis->ZDCTDCTrigger(fESD, AliTriggerAnalysis::kCentralBarrel) ;
 205   } else {
 206     // If it's MC, we use the energy
 207     Double_t minEnergy = 0;
 208     AliESDZDC *esdZDC = fESD->GetESDZDC();
 209     Double_t  zNCEnergy = esdZDC->GetZDCN1Energy();
 210     Double_t  zPCEnergy = esdZDC->GetZDCP1Energy();
 211     Double_t  zNAEnergy = esdZDC->GetZDCN2Energy();
 212     Double_t  zPAEnergy = esdZDC->GetZDCP2Energy();
 213     zdcA = (zNAEnergy>minEnergy || zPAEnergy>minEnergy);
 214     zdcC = (zNCEnergy>minEnergy || zPCEnergy>minEnergy);
 215   }
 216
 217
 218   // Some macros for the online triggers
 219   Bool_t cMBS2A = c0sm2 && c0v0A;
 220   Bool_t cMBS2C = c0sm2 && c0v0C;
 221   Bool_t cMBAC  = c0v0A && c0v0C;
 222
 223
 224   Bool_t vdArray[kNVDEntries];
 225   vdArray[kVDC0MBS1] = c0sm1;
 226   vdArray[kVDC0MBS2] = c0sm2;
 227   vdArray[kVDC0VBA]  = c0v0A;
 228   vdArray[kVDC0VBC]  = c0v0C;
 229   //vdArray[kVDC0OM2]  = c0OM2;
 230
 231   // Plots requested by Jurgen on 18/11/2010
 232   // FIXME: will skip everything else
 233
 234   GetHistoSPD1  ("PhysSel", "All events after physics selection and Francesco's cut")->Fill(outerLayerSPD);
 235   GetHistoTracks("PhysSel", "All events after physics selection and Francesco's cut")->Fill(ntracks);
 236   GetHistoV0M   ("PhysSel", "All events after physics selection and Francesco's cut")->Fill(multV0);
 237   if(c0v0A && c0v0C) {
 238     GetHistoSPD1  ("V0AND", "V0A & V0C")->Fill(outerLayerSPD);
 239     GetHistoTracks("V0AND", "V0A & V0C")->Fill(ntracks);
 240     GetHistoV0M   ("V0AND", "V0A & V0C")->Fill(multV0);
 241   }
 242   if(zdcA && zdcC) {
 243     GetHistoSPD1  ("ZDCAND", "ZDCA & ZDCC")->Fill(outerLayerSPD);
 244     GetHistoTracks("ZDCAND", "ZDCA & ZDCC")->Fill(ntracks);
 245     GetHistoV0M   ("ZDCAND", "ZDCA & ZDCC")->Fill(multV0);
 246   }
 247   if((c0v0A && c0v0C) && !(zdcA && zdcC)) {
 248     GetHistoSPD1  ("V0ANDNOTZDCAND", "(V0A & V0C) & !(ZDCA & ZDCC)")->Fill(outerLayerSPD);
 249     GetHistoTracks("V0ANDNOTZDCAND", "(V0A & V0C) & !(ZDCA & ZDCC)")->Fill(ntracks);
 250     GetHistoV0M   ("V0ANDNOTZDCAND", "(V0A & V0C) & !(ZDCA & ZDCC)")->Fill(multV0);
 251   }
 252   if((c0v0A && c0v0C) && !zdcA && !zdcC) {
 253     GetHistoSPD1  ("V0ANDNOTZDC", "(V0A & V0C) & !ZDCA & !ZDCC)")->Fill(outerLayerSPD);
 254     GetHistoTracks("V0ANDNOTZDC", "(V0A & V0C) & !ZDCA & !ZDCC)")->Fill(ntracks);
 255     GetHistoV0M   ("V0ANDNOTZDC", "(V0A & V0C) & !ZDCA & !ZDCC)")->Fill(multV0);
 256 }
 257   if((c0v0A && c0v0C) && ((!zdcA && zdcC) || (zdcA && !zdcC))) {
 258     GetHistoSPD1  ("V0ANDONEZDC", "(V0A & V0C) &  ((!ZDCA && ZDCC) || (ZDCA && !ZDCC)))")->Fill(outerLayerSPD);
 259     GetHistoTracks("V0ANDONEZDC", "(V0A & V0C) &  ((!ZDCA && ZDCC) || (ZDCA && !ZDCC)))")->Fill(ntracks);
 260     GetHistoV0M   ("V0ANDONEZDC", "(V0A & V0C) &  ((!ZDCA && ZDCC) || (ZDCA && !ZDCC)))")->Fill(multV0);
 261   }
 262
 263   if(((c0v0A && !c0v0C) || (!c0v0A && c0v0C)) && (zdcA && zdcC)) {
 264     GetHistoSPD1  ("V0ONEZDCBOTH", "((V0A && !V0C) || (!V0A && V0C)) && (ZDCA && ZDCC)")->Fill(outerLayerSPD);
 265     GetHistoTracks("V0ONEZDCBOTH", "((V0A && !V0C) || (!V0A && V0C)) && (ZDCA && ZDCC)")->Fill(ntracks);
 266     GetHistoV0M   ("V0ONEZDCBOTH", "((V0A && !V0C) || (!V0A && V0C)) && (ZDCA && ZDCC)")->Fill(multV0);
 267   }
 268
 269   if((c0v0A && c0v0C) && (zdcA && zdcC)) {
 270     GetHistoSPD1  ("V0ANDZDCAND", "(V0A & V0C) & (ZDCA & ZDCC)")->Fill(outerLayerSPD);
 271     GetHistoTracks("V0ANDZDCAND", "(V0A & V0C) & (ZDCA & ZDCC)")->Fill(ntracks);
 272     GetHistoV0M   ("V0ANDZDCAND", "(V0A & V0C) & (ZDCA & ZDCC)")->Fill(multV0);
 273   }
 274
 275   if((c0v0A && zdcA && !zdcC && !c0v0C) || (c0v0C && zdcC && !zdcA && !c0v0A)) {
 276     GetHistoSPD1  ("OneSided", "(V0A & ZDCA & !ZDCC & !V0C) || (V0C & ZDCC & !ZDCA & !V0A)")->Fill(outerLayerSPD);
 277     GetHistoTracks("OneSided", "(V0A & ZDCA & !ZDCC & !V0C) || (V0C & ZDCC & !ZDCA & !V0A)")->Fill(ntracks);
 278     GetHistoV0M   ("OneSided", "(V0A & ZDCA & !ZDCC & !V0C) || (V0C & ZDCC & !ZDCA & !V0A)")->Fill(multV0);
 279   }
 280
 281   // GetHistoCorrelationSPDTPCVz("All", "After physics selection and Francesco's cut")->Fill(vtxESDSPD->GetZ(),vtxESDTPC->GetZ());
 282   // if(cut1T0)   GetHistoCorrelationSPDTPCVz("Cut1T0", "T0 Cut 1")->Fill(vtxESDSPD->GetZ(),vtxESDTPC->GetZ());
 283   // if(cut2T0)   GetHistoCorrelationSPDTPCVz("Cut2T0", "T0 Cut 2")->Fill(vtxESDSPD->GetZ(),vtxESDTPC->GetZ());
 284
 285   // GetHistoCorrelationContrTPCSPDCls("All", "After physics selection and Francesco's cut")->Fill(totalClusSPD,tpcContr);
 286   // if(cut1T0)   GetHistoCorrelationContrTPCSPDCls("Cut1T0", "T0 Cut 1")->Fill(totalClusSPD,tpcContr);
 287   // if(cut2T0)   GetHistoCorrelationContrTPCSPDCls("Cut2T0", "T0 Cut 2")->Fill(totalClusSPD,tpcContr);
 288
 289   // GetHistoCorrelationTrackletsSPDCls("All", "After physics selection and Francesco's cut")->Fill(totalClusSPD,ntracklets);
 290   // if(cut1T0)   GetHistoCorrelationTrackletsSPDCls("Cut1T0", "T0 Cut 1")->Fill(totalClusSPD,ntracklets);
 291   // if(cut2T0)   GetHistoCorrelationTrackletsSPDCls("Cut2T0", "T0 Cut 2")->Fill(totalClusSPD,ntracklets);
 292
 293
 294   return; // FIXME
 295
 296
 297   // Reject background
 298   if (v0BG && fRejectBGWithV0) {
 299     cout << "Rejection BG" << endl;
 300
 301     return;
 302   }
 303   // Fill global histos
 304   GetHistoTracklets("all","All events")->Fill(ntracklets);
 305   FillTriggerOverlaps("All", "All Events",vdArray);
 306
 307   // Fill some combination of trigger classes
 308   Bool_t cmbs1aOnline = fESD->IsTriggerClassFired("CMBS1A-B-NOPF-ALL");
 309   Bool_t cmbs1cOnline = fESD->IsTriggerClassFired("CMBS1C-B-NOPF-ALL");
 310   Bool_t cmbacOnline  = fESD->IsTriggerClassFired("CMBAC-B-NOPF-ALL");
 311
 312   Bool_t twoOutOfThree = kFALSE;
 313   if (cmbs1aOnline || cmbs1cOnline ||cmbacOnline) twoOutOfThree = kTRUE;
 314   if (twoOutOfThree) GetHistoTracklets("All_TwoOutOfThree" ,"Events 2-out-of-3 online" )->Fill(ntracklets);
 315
 316
 317   // loop over trigger classes in the event
 318   TObjArray * tokens = 0;
 319   if(fIsMC) {
 320     // in case of montecarlo I override the trigger class
 321     tokens = new TObjArray;
 322     tokens->SetOwner();
 323     //    tokens->Add(new TObjString("CINT1B-ABCE-NOPF-ALL"));
 324     tokens->Add(new TObjString("MC"));
 325   }
 326   else {
 327     TString trgClasses = fESD->GetFiredTriggerClasses();
 328     tokens = trgClasses.Tokenize("  ");
 329   }
 330   TIterator  * iter = (TIterator*) tokens->MakeIterator();
 331
 332   TString classes = fESD->GetFiredTriggerClasses();
 333
 334   // if (classes.Contains("SMH")) {
 335   //    tokens->Print();
 336   //    cout << classes.Data() << endl;
 337   // }
 338   //  iter->Reset();
 339   //Int_t itoken = 0;
 340   TObjString * tok=0;
 341   while((tok = (TObjString*) iter->Next())){
 342     // clean up trigger name
 343     TString trg = tok->GetString();
 344     trg.Strip(TString::kTrailing, ' ');
 345     trg.Strip(TString::kLeading, ' ');
 346
 347     fHistoSuffix = "_";
 348     fHistoSuffix += trg;
 349
 350     //    cout << itoken++ << " " << trg.Data() << endl;
 351     // continue;
 352     //    if (trg.Contains("SMH")) cout << itoken++ << " " << trg.Data() << endl;
 353
 354     // Fill tracklets
 355     GetHistoTracklets("All" ,"Events no offline trigger" )->Fill(ntracklets);
 356
 357
 358     // Fill histograms mismatchs
 359     // TODO: check mismatch trigger class
 360     if(nFastOrOffline != nFastOrOnline) {
 361       GetHistoTracklets("mismatchingFastOr", "Events where fast or offline differs from fast-or online")->Fill(ntracklets);
 362     }
 363
 364     if (c0v0A != v0AHW){
 365       GetHistoTracklets("mismatchingV0A", "Events where V0A offline differs from V0A online")->Fill(ntracklets);
 366     }
 367
 368     if (c0v0C != v0CHW){
 369       GetHistoTracklets("mismatchingV0C", "Events where V0C offline differs from V0C online")->Fill(ntracklets);
 370     }
 371
 372     // Fill a tracklet histo for each trigger type
 373     if(c0sm1)  GetHistoTracklets("c0sm1" ,"Events were trigger c0sm1 fired" )->Fill(ntracklets);
 374     if(c0sm2)  GetHistoTracklets("c0sm2" ,"Events were trigger c0sm2 fired" )->Fill(ntracklets);
 375     if(c0sm3)  GetHistoTracklets("c0sm3" ,"Events were trigger c0sm3 fired" )->Fill(ntracklets);
 376     if(c0sm4)  GetHistoTracklets("c0sm4" ,"Events were trigger c0sm4 fired" )->Fill(ntracklets);
 377     if(c0sm5)  GetHistoTracklets("c0sm5" ,"Events were trigger c0sm5 fired" )->Fill(ntracklets);
 378     if(c0OM2)  GetHistoTracklets("c0OM2" ,"Events were trigger c0OM2 fired" )->Fill(ntracklets);
 379     if(c0OM3)  GetHistoTracklets("c0OM3" ,"Events were trigger c0OM3 fired" )->Fill(ntracklets);
 380     if(c0v0A)  GetHistoTracklets("c0v0A" ,"Events were trigger c0v0A fired" )->Fill(ntracklets);
 381     if(c0v0C)  GetHistoTracklets("c0v0C" ,"Events were trigger c0v0C fired" )->Fill(ntracklets);
 382     if(cMBS2A) GetHistoTracklets("cMBS2A","Events were trigger cMBS2A fired")->Fill(ntracklets);
 383     if(cMBS2C) GetHistoTracklets("cMBS2C","Events were trigger cMBS2C fired")->Fill(ntracklets);
 384     if(cMBAC ) GetHistoTracklets("cMBAC", "Events were trigger cMBAC  fired")->Fill(ntracklets);
 385     if(zdcA )  GetHistoTracklets("cZDCA", "Events were trigger cZDCA  fired")->Fill(ntracklets);
 386     if(zdcC )  GetHistoTracklets("cZDCC", "Events were trigger cZDCC  fired")->Fill(ntracklets);
 387     if(!zdcA && !zdcC )  GetHistoTracklets("NoZDC", "Events were zdc trigger don't  fired")->Fill(ntracklets);
 388     if( (zdcA && !zdcC) || (!zdcA && zdcC) )  GetHistoTracklets("OneSideZDC", "Events were only 1 ZDC  fired")->Fill(ntracklets);
 389     if( (zdcA && zdcC) )  GetHistoTracklets("TwoSideZDC", "Events were both  ZDC  fired")->Fill(ntracklets);
 390
 391     if(twoOutOfThree) {
 392       if(!zdcA && !zdcC )  GetHistoTracklets("NoZDC_TwoOutOfThree", "Events were zdc trigger don't  fired")->Fill(ntracklets);
 393       if( (zdcA && !zdcC) || (!zdcA && zdcC) )  GetHistoTracklets("OneSideZDC_TwoOutOfThree", "Events were only 1 ZDC  fired")->Fill(ntracklets);
 394       if( (zdcA && zdcC) )  GetHistoTracklets("TwoSideZDC_TwoOutOfThree", "Events were both  ZDC  fired")->Fill(ntracklets);
 395     }
 396     if(zdcBar) GetHistoTracklets("cZDCBar","Events were trigger cZDCB  fired")->Fill(ntracklets);
 397     //  if() GetHistoTracklets("","Events were trigger  fired");
 398
 399     // Fill trigger overlaps
 400     FillTriggerOverlaps("All", "All Events in trigger class",vdArray);
 401
 402
 403     // Fill kinematic variables for peripheral events
 404     static AliESDtrackCuts * cuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(); // FIXME: change this ?
 405     if (ntracklets < fNTrackletsCutKine) {
 406       // 1. Loop on tracklets
 407       for(Int_t itracklet = 0; itracklet < ntracklets; itracklet++){
 408         GetHistoEta("All", Form("Tracklet #eta distribution, ntracklets < %d",fNTrackletsCutKine))->Fill(mult->GetEta(itracklet));
 409       }
 410       // 2. loop on tracks
 411       TObjArray * goodTracks = cuts->GetAcceptedTracks(fESD);
 412       TIterator * iterTracks = goodTracks->MakeIterator();
 413       AliESDtrack * esdTrack = 0;
 414       while ((esdTrack = (AliESDtrack*) iterTracks->Next())) {
 415         GetHistoPt("All", Form("Tracklet p_{T} distribution, ntracklets < %d",fNTrackletsCutKine))->Fill(esdTrack->Pt());
 416       }
 417       delete goodTracks;
 418     }
 419
 420
 421   }
 422   delete tokens;
 423
 424     // if (fIsMC) {
 425
 426
 427   //   if (!fMCEvent) {
 428   //     AliError("No MC info found");
 429   //   } else {
 430
 431   //     //loop on the MC event
 432   //     //      Int_t nMCTracks = fMCEvent->GetNumberOfTracks();
 433   //     Int_t offset    = fMCEvent->GetPrimaryOffset();
 434   //     Int_t nMCTracks = fMCEvent->GetNumberOfPrimaries()+offset;
 435   //     for (Int_t ipart=offset; ipart<nMCTracks; ipart++) {
 436
 437   //    AliMCParticle *mcPart  = (AliMCParticle*)fMCEvent->GetTrack(ipart);
 438
 439   //    // We don't care about neutrals and non-physical primaries
 440   //    if(mcPart->Charge() == 0) continue;
 441
 442   // PHYSICAL PRIMARY
 443   //    // Get MC vertex
 444     //  TArrayF   vertex;
 445   //    fMCEvent->GenEventHeader()->PrimaryVertex(vertex);
 446   //    Float_t zv = vertex[2];
 447   //    //      Float_t zv = vtxESD->GetZ();
 448   //    // Fill generated histo
 449   //    hTracks[AliAnalysisMultPbTrackHistoManager::kHistoGen]->Fill(mcPart->Pt(),mcPart->Eta(),zv);
 450
 451   //     }
 452   //   }
 453   // }
 454
 455
 456
 457
 458 }
 459
 460 void   AliAnalysisTaskTriggerStudy::Terminate(Option_t *){
 461   // terminate
 462   // Save output in a more friendly format
 463   fHistoList = dynamic_cast<AliHistoListWrapper*> (GetOutputData(1));
 464   if (!fHistoList){
 465     Printf("ERROR: fHistoList not available");
 466     return;
 467   }
 468   TFile * f = new TFile("trigger_study.root", "recreate");
 469   fHistoList->GetList()->Write();
 470   f->Close();
 471
 472 }
 473
 474 TH1 *   AliAnalysisTaskTriggerStudy::GetHistoTracklets(const char * name, const char * title){
 475   // Book histo of events vs ntracklets, if needed
 476
 477   TString hname = "hTracklets_";
 478   hname+=name;
 479   hname+=fHistoSuffix;
 480   TH1 * h = (TH1*) fHistoList->GetList()->FindObject(hname.Data());
 481
 482   if(!h) {
 483     AliInfo(Form("Booking histo %s",hname.Data()));
 484     Bool_t oldStatus = TH1::AddDirectoryStatus();
 485     TH1::AddDirectory(kFALSE);
 486     // h = new TH1F (hname.Data(), title, 1000, -0.5, 10000-0.5);
 487     h = new TH1F (hname.Data(), title, 4000, -0.5, 4000-0.5);
 488     h->Sumw2();
 489     h->SetXTitle("ntracklets");
 490     fHistoList->GetList()->Add(h);
 491     TH1::AddDirectory(oldStatus);
 492   }
 493   return h;
 494 }
 495
 496
 497 TH1 *   AliAnalysisTaskTriggerStudy::GetHistoSPD1(const char * name, const char * title){
 498   // Book histo of events vs ntracklets, if needed
 499
 500   TString hname = "hSPD1_";
 501   hname+=name;
 502   hname+=fHistoSuffix;
 503   TH1 * h = (TH1*) fHistoList->GetList()->FindObject(hname.Data());
 504
 505   // const Int_t nbin = 118;
 506   // Float_t bins[nbin] = {0};
 507   // bins[0] = 0;
 508   // for(Int_t ibin = 1; ibin <= nbin; ibin++){
 509   //   if (ibin < 100)
 510   //     bins[ibin] = bins [ibin-1]+1;
 511   //   else if (ibin < 1000)
 512   //     bins[ibin] = bins [ibin-1]+100;
 513   //   else if (ibin < 10000)
 514   //     bins[ibin] = bins [ibin-1]+1000;
 515   // }
 516
 517
 518   if(!h) {
 519     AliInfo(Form("Booking histo %s",hname.Data()));
 520     Bool_t oldStatus = TH1::AddDirectoryStatus();
 521     TH1::AddDirectory(kFALSE);
 522     //    h = new TH1F (hname.Data(), title, nbin, bins);
 523     //h = new TH1F (hname.Data(), title, 100, -0.5, 100-0.5);
 524     h = new TH1F (hname.Data(), title, 10000, -0.5, 10000-0.5);
 525     h->Sumw2();
 526     h->SetXTitle("SPD1");
 527     fHistoList->GetList()->Add(h);
 528     TH1::AddDirectory(oldStatus);
 529   }
 530   return h;
 531 }
 532
 533 TH1 *   AliAnalysisTaskTriggerStudy::GetHistoV0M(const char * name, const char * title){
 534   // Book histo of events vs ntracklets, if needed
 535
 536   TString hname = "hV0M_";
 537   hname+=name;
 538   hname+=fHistoSuffix;
 539   TH1 * h = (TH1*) fHistoList->GetList()->FindObject(hname.Data());
 540
 541   if(!h) {
 542     AliInfo(Form("Booking histo %s",hname.Data()));
 543     Bool_t oldStatus = TH1::AddDirectoryStatus();
 544     TH1::AddDirectory(kFALSE);
 545     // h = new TH1F (hname.Data(), title, 1000, -0.5, 10000-0.5);
 546     h = new TH1F (hname.Data(), title, 300, -0.5, 300-0.5);
 547     h->Sumw2();
 548     h->SetXTitle("V0M");
 549     fHistoList->GetList()->Add(h);
 550     TH1::AddDirectory(oldStatus);
 551   }
 552   return h;
 553 }
 554
 555
 556 TH1 *   AliAnalysisTaskTriggerStudy::GetHistoTracks(const char * name, const char * title){
 557   // Book histo of events vs ntracklets, if needed
 558
 559   TString hname = "hTPCTracks_";
 560   hname+=name;
 561   hname+=fHistoSuffix;
 562   TH1 * h = (TH1*) fHistoList->GetList()->FindObject(hname.Data());
 563
 564   if(!h) {
 565     AliInfo(Form("Booking histo %s",hname.Data()));
 566     Bool_t oldStatus = TH1::AddDirectoryStatus();
 567     TH1::AddDirectory(kFALSE);
 568     // h = new TH1F (hname.Data(), title, 1000, -0.5, 10000-0.5);
 569     h = new TH1F (hname.Data(), title, 50, -0.5, 50-0.5);
 570     h->Sumw2();
 571     h->SetXTitle("ntracks");
 572     fHistoList->GetList()->Add(h);
 573     TH1::AddDirectory(oldStatus);
 574   }
 575   return h;
 576 }
 577
 578
 579 TH1 *   AliAnalysisTaskTriggerStudy::GetHistoEta(const char * name, const char * title){
 580   // Book histo of events vs ntracklets, if needed
 581
 582   TString hname = "hEta_";
 583   hname+=name;
 584   hname+=fHistoSuffix;
 585   TH1 * h = (TH1*) fHistoList->GetList()->FindObject(hname.Data());
 586
 587   if(!h) {
 588     AliInfo(Form("Booking histo %s",hname.Data()));
 589     Bool_t oldStatus = TH1::AddDirectoryStatus();
 590     TH1::AddDirectory(kFALSE);
 591     // h = new TH1F (hname.Data(), title, 1000, -0.5, 10000-0.5);
 592     h = new TH1F (hname.Data(), title, 20, -2, 2);
 593     h->Sumw2();
 594     h->SetXTitle("#eta");
 595     fHistoList->GetList()->Add(h);
 596     TH1::AddDirectory(oldStatus);
 597   }
 598   return h;
 599 }
 600
 601 TH1 *   AliAnalysisTaskTriggerStudy::GetHistoPt(const char * name, const char * title){
 602   // Book histo of pt distribution, if needed
 603
 604   TString hname = "hPt_";
 605   hname+=name;
 606   hname+=fHistoSuffix;
 607   TH1 * h = (TH1*) fHistoList->GetList()->FindObject(hname.Data());
 608
 609   if(!h) {
 610     AliInfo(Form("Booking histo %s",hname.Data()));
 611     Bool_t oldStatus = TH1::AddDirectoryStatus();
 612     TH1::AddDirectory(kFALSE);
 613     // h = new TH1F (hname.Data(), title, 1000, -0.5, 10000-0.5);
 614     h = new TH1F (hname.Data(), title, 100, -0.5, 1-0.5);
 615     h->Sumw2();
 616     h->SetXTitle("p_{T}");
 617     fHistoList->GetList()->Add(h);
 618     TH1::AddDirectory(oldStatus);
 619   }
 620   return h;
 621 }
 622
 623 TH1 *   AliAnalysisTaskTriggerStudy::GetHistoCorrelationSPDTPCVz(const char * name, const char * title){
 624   // Book histo of events vs ntracklets, if needed
 625
 626   TString hname = "hTPCvsSPD_";
 627   hname+=name;
 628   hname+=fHistoSuffix;
 629   TH1 * h = (TH1*) fHistoList->GetList()->FindObject(hname.Data());
 630
 631   if(!h) {
 632     AliInfo(Form("Booking histo %s",hname.Data()));
 633     Bool_t oldStatus = TH1::AddDirectoryStatus();
 634     TH1::AddDirectory(kFALSE);
 635     // h = new TH1F (hname.Data(), title, 1000, -0.5, 10000-0.5);
 636     h = new TH2F (hname.Data(), title, 80, -20, 20,  80, -20, 20);
 637     h->Sumw2();
 638     h->SetXTitle("SPD Vz");
 639     h->SetYTitle("TPC Vz");
 640     fHistoList->GetList()->Add(h);
 641     TH1::AddDirectory(oldStatus);
 642   }
 643   return h;
 644 }
 645
 646 TH1 *   AliAnalysisTaskTriggerStudy::GetHistoCorrelationContrTPCSPDCls(const char * name, const char * title){
 647   // Book histo of events vs ntracklets, if needed
 648
 649   TString hname = "hContrTPCvsSPDCls_";
 650   hname+=name;
 651   hname+=fHistoSuffix;
 652   TH1 * h = (TH1*) fHistoList->GetList()->FindObject(hname.Data());
 653
 654   if(!h) {
 655     AliInfo(Form("Booking histo %s",hname.Data()));
 656     Bool_t oldStatus = TH1::AddDirectoryStatus();
 657     TH1::AddDirectory(kFALSE);
 658     // h = new TH1F (hname.Data(), title, 1000, -0.5, 10000-0.5);
 659     h = new TH2F (hname.Data(), title, 1000, 0, 9000,  1000, 0, 5000);
 660     h->Sumw2();
 661     h->SetXTitle("SPD Clusters");
 662     h->SetYTitle("TPC V Contributors");
 663     fHistoList->GetList()->Add(h);
 664     TH1::AddDirectory(oldStatus);
 665   }
 666   return h;
 667 }
 668 TH1 *   AliAnalysisTaskTriggerStudy::GetHistoCorrelationTrackletsSPDCls(const char * name, const char * title){
 669   // Book histo of events vs ntracklets, if needed
 670
 671   TString hname = "hTrackletsvsSPDCls_";
 672   hname+=name;
 673   hname+=fHistoSuffix;
 674   TH1 * h = (TH1*) fHistoList->GetList()->FindObject(hname.Data());
 675
 676   if(!h) {
 677     AliInfo(Form("Booking histo %s",hname.Data()));
 678     Bool_t oldStatus = TH1::AddDirectoryStatus();
 679     TH1::AddDirectory(kFALSE);
 680     // h = new TH1F (hname.Data(), title, 1000, -0.5, 10000-0.5);
 681     h = new TH2F (hname.Data(), title, 1000, 0, 9000,  1000, 0, 5000);
 682     h->Sumw2();
 683     h->SetXTitle("SPD Clusters");
 684     h->SetYTitle("N tracklets");
 685     fHistoList->GetList()->Add(h);
 686     TH1::AddDirectory(oldStatus);
 687   }
 688   return h;
 689 }
 690
 691
 692 void AliAnalysisTaskTriggerStudy::FillTriggerOverlaps (const char * name, const char * title, Bool_t * vdArray){
 693   //Fills a histo with the different trigger statistics in a venn like diagramm. Books it if needed.
 694
 695   // Get or book histo
 696   TString hname = "hTrigStat_";
 697   hname+=name;
 698   hname+=fHistoSuffix;
 699   TH1 * h = (TH1*) fHistoList->GetList()->FindObject(hname.Data());
 700
 701   if(!h) {
 702     AliInfo(Form("Booking histo %s",hname.Data()));
 703     Bool_t oldStatus = TH1::AddDirectoryStatus();
 704     TH1::AddDirectory(kFALSE);
 705     Int_t nbins = 0;
 706     for(Int_t ientry = 0; ientry < kNVDEntries; ientry++){
 707       nbins = nbins | (1<<ientry);
 708     }
 709     //    cout << "NBINS " << nbins << endl;
 710     nbins = nbins+1;
 711     h = new TH1I (hname, title, nbins, -0.5, nbins-0.5);
 712     fHistoList->GetList()->Add(h);
 713     TH1::AddDirectory(oldStatus);
 714
 715     // we look at the combinations of n triggers
 716     // We set a bit for each trigger to fill the diagram
 717     // This is much simpler and faster than any recursive function
 718     h->GetXaxis()->SetBinLabel(1,"NONE");
 719     for(Int_t ibin = 1; ibin < nbins; ibin++){
 720       TString binname = "";
 721       Bool_t first = kTRUE;
 722       for(Int_t ivdentry = 0; ivdentry < kNVDEntries; ivdentry++){
 723         if (ibin & (1<<ivdentry)) {
 724           if(!first) binname += " & ";
 725           binname += kVDNames[ivdentry];
 726           first=kFALSE;
 727         }
 728       }
 729       h->GetXaxis()->SetBinLabel(ibin+1,binname.Data());
 730     }
 731
 732   }
 733
 734   UInt_t mask = 0;
 735   for(Int_t ivdentry = 0; ivdentry < kNVDEntries; ivdentry++){
 736     if(vdArray[ivdentry]) {
 737       mask  = mask | (1<<ivdentry);
 738       //      cout << " 1 "   ;
 739     } //else cout << " 0 ";
 740   }
 741   //  cout << hex << " = " << mask << endl;
 742
 743   h->Fill(mask);
 744
 745 }