PWGJE/PWGJE/AliAnalysisTaskJetCore.cxx

   1 // ******************************************
   2 // This task computes several jet observables like
   3 // the fraction of energy in inner and outer coronnas,
   4 // jet-track correlations,triggered jet shapes and
   5 // correlation strength distribution of particles inside jets.
   6 // Author: lcunquei@cern.ch
   7 // *******************************************
   8
   9
  10 /**************************************************************************
  11  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  12  *                                                                        *
  13  * Author: The ALICE Off-line Project.                                    *
  14  * Contributors are mentioned in the code where appropriate.              *
  15  *                                                                        *
  16  * Permission to use, copy, modify and distribute this software and its   *
  17  * documentation strictly for non-commercial purposes is hereby granted   *
  18  * without fee, provided that the above copyright notice appears in all   *
  19  * copies and that both the copyright notice and this permission notice   *
  20  * appear in the supporting documentation. The authors make no claims     *
  21  * about the suitability of this software for any purpose. It is          *
  22  * provided "as is" without express or implied warranty.                  *
  23  **************************************************************************/
  24
  25
  26 #include "TChain.h"
  27 #include "TTree.h"
  28 #include "TMath.h"
  29 #include "TH1F.h"
  30 #include "TH2F.h"
  31 #include "TH3F.h"
  32 #include "THnSparse.h"
  33 #include "TCanvas.h"
  34 #include "TRandom3.h"
  35 #include "AliLog.h"
  36
  37 #include "AliAnalysisTask.h"
  38 #include "AliAnalysisManager.h"
  39
  40 #include "AliVEvent.h"
  41 #include "AliESDEvent.h"
  42 #include "AliESDInputHandler.h"
  43 #include "AliCentrality.h"
  44 #include "AliAnalysisHelperJetTasks.h"
  45 #include "AliInputEventHandler.h"
  46 #include "AliAODJetEventBackground.h"
  47 #include "AliAODMCParticle.h"
  48 #include "AliAnalysisTaskFastEmbedding.h"
  49 #include "AliAODEvent.h"
  50 #include "AliAODHandler.h"
  51 #include "AliAODJet.h"
  52
  53 #include "AliAnalysisTaskJetCore.h"
  54
  55 using std::cout;
  56 using std::endl;
  57
  58 ClassImp(AliAnalysisTaskJetCore)
  59
  60 AliAnalysisTaskJetCore::AliAnalysisTaskJetCore() :
  61 AliAnalysisTaskSE(),
  62 fESD(0x0),
  63 fAODIn(0x0),
  64 fAODOut(0x0),
  65 fAODExtension(0x0),
  66 fBackgroundBranch(""),
  67 fNonStdFile(""),
  68 fIsPbPb(kTRUE),
  69 fOfflineTrgMask(AliVEvent::kAny),
  70 fMinContribVtx(1),
  71 fVtxZMin(-10.),
  72 fVtxZMax(10.),
  73 fEvtClassMin(0),
  74 fEvtClassMax(4),
  75 fFilterMask(0),
  76 fFilterMaskBestPt(0),
  77 fFilterType(0),
  78 fRadioFrac(0.2),
  79 fMinDist(0.1),
  80 fCentMin(0.),
  81 fCentMax(100.),
  82 fNInputTracksMin(0),
  83 fNInputTracksMax(-1),
  84 fRequireITSRefit(0),
  85 fApplySharedClusterCut(0),
  86 fAngStructCloseTracks(0),
  87 fCheckMethods(0),
  88 fDoEventMixing(0),
  89 fFlagPhiBkg(0),
  90 fFlagEtaBkg(0),
  91 fFlagJetHadron(0),
  92 fDodiHadron(0),
  93 fFrac(0.8),
  94 fTTLowRef(11),
  95 fTTUpRef(13),
  96 fTTLowSig(15),
  97 fTTUpSig(19),
  98 fHardest(0),
  99 fFlagRandom(0),
 100 fFlagOnlyRecoil(0),
 101 fFlagOnlyHardest(1),
 102 fTrackTypeRec(kTrackUndef),
 103 fRPAngle(0),
 104 fNRPBins(50),
 105 fSemigoodCorrect(0),
 106 fHolePos(4.71),
 107 fHoleWidth(0.2),
 108 fJetEtaMin(-.5),
 109 fJetEtaMax(.5),
 110 fNevents(0),
 111 fTindex(0),
 112 fTrigBufferIndex(0),
 113 fCountAgain(0),
 114 fJetPtMin(20.),
 115 fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
 116 fJetPtFractionMin(0.5),
 117 fNMatchJets(4),
 118 fMatchMaxDist(0.8),
 119 fKeepJets(kFALSE),
 120 fRunAnaAzimuthalCorrelation(kFALSE),
 121 fkNbranches(2),
 122 fkEvtClasses(12),
 123 fOutputList(0x0),
 124 fbEvent(kTRUE),
 125 fHistEvtSelection(0x0),
 126 fhnDeltaR(0x0),
 127 fhnMixedEvents(0x0),
 128 fh2JetCoreMethod1C10(0x0),
 129 fh2JetCoreMethod2C10(0x0),
 130 fh2JetCoreMethod1C20(0x0),
 131 fh2JetCoreMethod2C20(0x0),
 132 fh2JetCoreMethod1C30(0x0),
 133 fh2JetCoreMethod2C30(0x0),
 134 fh2JetCoreMethod1C60(0x0),
 135 fh2JetCoreMethod2C60(0x0),
 136 fh3JetTrackC3060(0x0),
 137 fh3JetTrackC20(0x0),
 138 fh2AngStructpt1C10(0x0),
 139 fh2AngStructpt2C10(0x0),
 140 fh2AngStructpt3C10(0x0),
 141 fh2AngStructpt4C10(0x0),
 142 fh2AngStructpt1C20(0x0),
 143 fh2AngStructpt2C20(0x0),
 144 fh2AngStructpt3C20(0x0),
 145 fh2AngStructpt4C20(0x0),
 146 fh2AngStructpt1C30(0x0),
 147 fh2AngStructpt2C30(0x0),
 148 fh2AngStructpt3C30(0x0),
 149 fh2AngStructpt4C30(0x0),
 150 fh2AngStructpt1C60(0x0),
 151 fh2AngStructpt2C60(0x0),
 152 fh2AngStructpt3C60(0x0),
 153 fh2AngStructpt4C60(0x0),
 154 fh1TrigRef(0x0),
 155 fh1TrigSig(0x0),
 156 fh1TrackPhiDistance(0x0),
 157 fh1TrackRDistance(0x0),
 158 fh2Ntriggers(0x0),
 159 fh2Ntriggers2C10(0x0),
 160 fh2Ntriggers2C20(0x0),
 161 fh3JetDensity(0x0),
 162 fh3JetDensityA4(0x0),
 163 fh2RPJetsC10(0x0),
 164 fh2RPJetsC20(0x0),
 165 fh2RPTC10(0x0),
 166 fh2RPTC20(0x0),
 167 fHJetSpec(0x0),
 168 fhTTPt(0x0),
 169 fHJetPhiCorr(0x0),
 170 fRandom(0x0)
 171
 172 {
 173    // default Constructor
 174
 175
 176  // Trigger buffer.
 177    for(Int_t i=0; i<10; i++) {
 178                 for(Int_t j=0; j<6; j++) {
 179                         fTrigBuffer[i][j]=0;
 180                 }
 181         }
 182
 183
 184
 185
 186
 187    fJetBranchName[0] = "";
 188    fJetBranchName[1] = "";
 189
 190    fListJets[0] = new TList;
 191    fListJets[1] = new TList;
 192 }
 193
 194 AliAnalysisTaskJetCore::AliAnalysisTaskJetCore(const char *name) :
 195 AliAnalysisTaskSE(name),
 196 fESD(0x0),
 197 fAODIn(0x0),
 198 fAODOut(0x0),
 199 fAODExtension(0x0),
 200 fBackgroundBranch(""),
 201 fNonStdFile(""),
 202 fIsPbPb(kTRUE),
 203 fOfflineTrgMask(AliVEvent::kAny),
 204 fMinContribVtx(1),
 205 fVtxZMin(-10.),
 206 fVtxZMax(10.),
 207 fEvtClassMin(0),
 208 fEvtClassMax(4),
 209 fFilterMask(0),
 210 fFilterMaskBestPt(0),
 211 fFilterType(0),
 212 fRadioFrac(0.2),
 213 fMinDist(0.1),
 214 fCentMin(0.),
 215 fCentMax(100.),
 216 fNInputTracksMin(0),
 217 fNInputTracksMax(-1),
 218 fRequireITSRefit(0),
 219 fApplySharedClusterCut(0),
 220 fAngStructCloseTracks(0),
 221 fCheckMethods(0),
 222 fDoEventMixing(0),
 223 fFlagPhiBkg(0),
 224 fFlagEtaBkg(0),
 225 fFlagJetHadron(0),
 226 fDodiHadron(0),
 227 fFrac(0.8),
 228 fTTLowRef(11),
 229 fTTUpRef(13),
 230 fTTLowSig(15),
 231 fTTUpSig(19),
 232 fHardest(0),
 233 fFlagRandom(0),
 234 fFlagOnlyRecoil(0),
 235 fFlagOnlyHardest(1),
 236 fTrackTypeRec(kTrackUndef),
 237 fRPAngle(0),
 238 fNRPBins(50),
 239 fSemigoodCorrect(0),
 240 fHolePos(4.71),
 241 fHoleWidth(0.2),
 242 fJetEtaMin(-.5),
 243 fJetEtaMax(.5),
 244 fNevents(0),
 245 fTindex(0),
 246 fTrigBufferIndex(0),
 247 fCountAgain(0),
 248 fJetPtMin(20.),
 249 fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
 250 fJetPtFractionMin(0.5),
 251 fNMatchJets(4),
 252 fMatchMaxDist(0.8),
 253 fKeepJets(kFALSE),
 254 fRunAnaAzimuthalCorrelation(kFALSE),
 255 fkNbranches(2),
 256 fkEvtClasses(12),
 257 fOutputList(0x0),
 258 fbEvent(kTRUE),
 259 fHistEvtSelection(0x0),
 260 fhnDeltaR(0x0),
 261 fhnMixedEvents(0x0),
 262 fh2JetCoreMethod1C10(0x0),
 263 fh2JetCoreMethod2C10(0x0),
 264 fh2JetCoreMethod1C20(0x0),
 265 fh2JetCoreMethod2C20(0x0),
 266 fh2JetCoreMethod1C30(0x0),
 267 fh2JetCoreMethod2C30(0x0),
 268 fh2JetCoreMethod1C60(0x0),
 269 fh2JetCoreMethod2C60(0x0),
 270 fh3JetTrackC3060(0x0),
 271 fh3JetTrackC20(0x0),
 272 fh2AngStructpt1C10(0x0),
 273 fh2AngStructpt2C10(0x0),
 274 fh2AngStructpt3C10(0x0),
 275 fh2AngStructpt4C10(0x0),
 276 fh2AngStructpt1C20(0x0),
 277 fh2AngStructpt2C20(0x0),
 278 fh2AngStructpt3C20(0x0),
 279 fh2AngStructpt4C20(0x0),
 280 fh2AngStructpt1C30(0x0),
 281 fh2AngStructpt2C30(0x0),
 282 fh2AngStructpt3C30(0x0),
 283 fh2AngStructpt4C30(0x0),
 284 fh2AngStructpt1C60(0x0),
 285 fh2AngStructpt2C60(0x0),
 286 fh2AngStructpt3C60(0x0),
 287 fh2AngStructpt4C60(0x0),
 288 fh1TrigRef(0x0),
 289 fh1TrigSig(0x0),
 290 fh1TrackPhiDistance(0x0),
 291 fh1TrackRDistance(0x0),
 292 fh2Ntriggers(0x0),
 293 fh2Ntriggers2C10(0x0),
 294 fh2Ntriggers2C20(0x0),
 295 fh3JetDensity(0x0),
 296 fh3JetDensityA4(0x0),
 297 fh2RPJetsC10(0x0),
 298 fh2RPJetsC20(0x0),
 299 fh2RPTC10(0x0),
 300 fh2RPTC20(0x0),
 301 fHJetSpec(0x0),
 302 fhTTPt(0x0),
 303 fHJetPhiCorr(0x0),
 304 fRandom(0x0)
 305
 306  {
 307    // Constructor
 308
 309
 310     for(Int_t i=0; i<10; i++) {
 311        for(Int_t j=0; j<6; j++) {
 312             fTrigBuffer[i][j]=0;
 313                 }
 314     }
 315
 316
 317
 318    fJetBranchName[0] = "";
 319    fJetBranchName[1] = "";
 320
 321    fListJets[0] = new TList;
 322    fListJets[1] = new TList;
 323
 324    DefineOutput(1, TList::Class());
 325 }
 326
 327 AliAnalysisTaskJetCore::~AliAnalysisTaskJetCore()
 328 {
 329    delete fListJets[0];
 330    delete fListJets[1];
 331 }
 332
 333 void AliAnalysisTaskJetCore::SetBranchNames(const TString &branch1, const TString &branch2)
 334 {
 335    fJetBranchName[0] = branch1;
 336    fJetBranchName[1] = branch2;
 337 }
 338
 339 void AliAnalysisTaskJetCore::Init()
 340 {
 341
 342    // check for jet branches
 343    if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
 344       AliError("Jet branch name not set.");
 345    }
 346
 347 }
 348
 349 void AliAnalysisTaskJetCore::UserCreateOutputObjects()
 350 {
 351    // Create histograms
 352    // Called once
 353    OpenFile(1);
 354    if(!fOutputList) fOutputList = new TList;
 355    fOutputList->SetOwner(kTRUE);
 356
 357    Bool_t oldStatus = TH1::AddDirectoryStatus();
 358    TH1::AddDirectory(kFALSE);
 359
 360
 361         // set seed
 362    fRandom = new TRandom3(0);
 363
 364
 365
 366
 367    fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 6, -0.5, 5.5);
 368    fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
 369    fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN");
 370    fHistEvtSelection->GetXaxis()->SetBinLabel(3,"event selection (rejected)");
 371    fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)");
 372    fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)");
 373    fHistEvtSelection->GetXaxis()->SetBinLabel(6,"multiplicity (rejected)");
 374
 375      UInt_t entries = 0; // bit coded, see GetDimParams() below
 376      entries = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6 |1<<7;
 377      fhnDeltaR = NewTHnSparseF("fhnDeltaR", entries);
 378
 379      //change binning in pTtrack
 380      Double_t *xPt3 = new Double_t[10];
 381      xPt3[0] = 0.;
 382      for(int i = 1; i<=9;i++){
 383       if(xPt3[i-1]<2)xPt3[i] = xPt3[i-1] + 0.4; // 1 - 5
 384       else if(xPt3[i-1]<11)xPt3[i] = xPt3[i-1] + 3; // 5 - 12
 385       else xPt3[i] = xPt3[i-1] + 150.; // 18
 386      }
 387     fhnDeltaR->SetBinEdges(2,xPt3);
 388     delete [] xPt3;
 389
 390     //change binning in HTI
 391      Double_t *xPt4 = new Double_t[14];
 392      xPt4[0] = 0.;
 393      for(int i = 1; i<=13;i++){
 394       if(xPt4[i-1]<10)xPt4[i] = xPt4[i-1] + 1; // 1 - 10
 395       else if(xPt4[i-1]<20)xPt4[i] = xPt4[i-1] + 5; // 10 - 12
 396       else xPt4[i] = xPt4[i-1] + 30.; // 13
 397      }
 398     fhnDeltaR->SetBinEdges(6,xPt4);
 399     delete [] xPt4;
 400
 401
 402
 403
 404
 405
 406
 407      if(fDoEventMixing){
 408      UInt_t cifras = 0; // bit coded, see GetDimParams() below
 409      cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<7;
 410      fhnMixedEvents = NewTHnSparseF("fhnMixedEvents", cifras);}
 411
 412     if(fCheckMethods){
 413     fh2JetCoreMethod1C10 = new TH2F("JetCoreMethod1C10","",150, 0., 150.,100, 0., 1.5);
 414     fh2JetCoreMethod2C10 = new TH2F("JetCoreMethod2C10","",150, 0., 150.,100, 0., 1.5);
 415     fh2JetCoreMethod1C20 = new TH2F("JetCoreMethod1C20","",150, 0., 150.,100, 0., 1.5);
 416     fh2JetCoreMethod2C20 = new TH2F("JetCoreMethod2C20","",150, 0., 150.,100, 0., 1.5);
 417     fh2JetCoreMethod1C30 = new TH2F("JetCoreMethod1C30","",150, 0., 150.,100, 0., 1.5);
 418     fh2JetCoreMethod2C30 = new TH2F("JetCoreMethod2C30","",150, 0., 150.,100, 0., 1.5);
 419     fh2JetCoreMethod1C60 = new TH2F("JetCoreMethod1C60","",150, 0., 150.,100, 0., 1.5);
 420     fh2JetCoreMethod2C60 = new TH2F("JetCoreMethod2C60","",150, 0., 150.,100, 0., 1.5);}
 421
 422     fh3JetTrackC3060=new TH3F("JetTrackC3060","",50,0,50,150,0.,150.,35,0.,3.5);
 423     fh3JetTrackC20=new TH3F("JetTrackC20","",50,0,50,150,0.,150.,35,0.,3.5);
 424     if(fAngStructCloseTracks>0){
 425     fh2AngStructpt1C10 = new TH2F("Ang struct pt1 C10","",15,0.,1.5,150,0.,10.);
 426     fh2AngStructpt2C10 = new TH2F("Ang struct pt2 C10","",15,0.,1.5,150,0.,10.);
 427     fh2AngStructpt3C10 = new TH2F("Ang struct pt3 C10","",15,0.,1.5,150,0.,10.);
 428     fh2AngStructpt4C10 = new TH2F("Ang struct pt4 C10","",15,0.,1.5,150,0.,10.);
 429     fh2AngStructpt1C20 = new TH2F("Ang struct pt1 C20","",15,0.,1.5,150,0.,10.);
 430     fh2AngStructpt2C20 = new TH2F("Ang struct pt2 C20","",15,0.,1.5,150,0.,10.);
 431     fh2AngStructpt3C20 = new TH2F("Ang struct pt3 C20","",15,0.,1.5,150,0.,10.);
 432     fh2AngStructpt4C20 = new TH2F("Ang struct pt4 C20","",15,0.,1.5,150,0.,10.);
 433     fh2AngStructpt1C30 = new TH2F("Ang struct pt1 C30","",15,0.,1.5,150,0.,10.);
 434     fh2AngStructpt2C30 = new TH2F("Ang struct pt2 C30","",15,0.,1.5,150,0.,10.);
 435     fh2AngStructpt3C30 = new TH2F("Ang struct pt3 C30","",15,0.,1.5,150,0.,10.);
 436     fh2AngStructpt4C30 = new TH2F("Ang struct pt4 C30","",15,0.,1.5,150,0.,10.);
 437     fh2AngStructpt1C60 = new TH2F("Ang struct pt1 C60","",15,0.,1.5,150,0.,10.);
 438     fh2AngStructpt2C60 = new TH2F("Ang struct pt2 C60","",15,0.,1.5,150,0.,10.);
 439     fh2AngStructpt3C60 = new TH2F("Ang struct pt3 C60","",15,0.,1.5,150,0.,10.);
 440     fh2AngStructpt4C60 = new TH2F("Ang struct pt4 C60","",15,0.,1.5,150,0.,10.); }
 441
 442
 443     fh1TrigRef=new TH1D("Trig Ref","",10,0.,10);
 444     fh1TrigSig=new TH1D("Trig Sig","",10,0.,10);
 445     fh1TrackPhiDistance=new TH1D("PhiDistance","",35,0.,3.5);
 446     fh1TrackRDistance=new TH1D("RDistance","",10,0.,1);
 447     fh2Ntriggers=new TH2F("# of triggers","",100,0.,100.,50,0.,50.);
 448     fh2Ntriggers2C10=new TH2F("# of triggers2C10","",50,0.,50.,50,0.,50.);
 449     fh2Ntriggers2C20=new TH2F("# of triggers2C20","",50,0.,50.,50,0.,50.);
 450     fh3JetDensity=new TH3F("Jet density vs mutliplicity A>0.07","",100,0.,4000.,100,0.,5.,25,0.,50.);
 451     fh3JetDensityA4=new TH3F("Jet density vs multiplicity A>0.4","",100,0.,4000.,100,0.,5.,25,0.,50.);
 452     fh2RPJetsC10=new TH2F("RPJetC10","",35,0.,3.5,100,0.,100.);
 453     fh2RPJetsC20=new TH2F("RPJetC20","",35,0.,3.5,100,0.,100.);
 454     fh2RPTC10=new TH2F("RPTriggerC10","",35,0.,3.5,50,0.,50.);
 455     fh2RPTC20=new TH2F("RPTriggerC20","",35,0.,3.5,50,0.,50.);
 456
 457
 458
 459
 460    fOutputList->Add(fHistEvtSelection);
 461
 462    fOutputList->Add(fhnDeltaR);
 463
 464    fOutputList->Add(fhnMixedEvents);
 465
 466
 467
 468       if(fCheckMethods){
 469
 470       fOutputList->Add(fh2JetCoreMethod1C10);
 471       fOutputList->Add(fh2JetCoreMethod2C10);
 472       fOutputList->Add(fh2JetCoreMethod1C20);
 473       fOutputList->Add(fh2JetCoreMethod2C20);
 474       fOutputList->Add(fh2JetCoreMethod1C30);
 475       fOutputList->Add(fh2JetCoreMethod2C30);
 476       fOutputList->Add(fh2JetCoreMethod1C60);
 477       fOutputList->Add(fh2JetCoreMethod2C60);}
 478
 479       fOutputList->Add(fh3JetTrackC3060);
 480       fOutputList->Add(fh3JetTrackC20);
 481
 482
 483
 484
 485
 486         if(fAngStructCloseTracks>0){
 487        fOutputList->Add(fh2AngStructpt1C10);
 488        fOutputList->Add(fh2AngStructpt2C10);
 489        fOutputList->Add(fh2AngStructpt3C10);
 490        fOutputList->Add(fh2AngStructpt4C10);
 491        fOutputList->Add(fh2AngStructpt1C20);
 492        fOutputList->Add(fh2AngStructpt2C20);
 493        fOutputList->Add(fh2AngStructpt3C20);
 494        fOutputList->Add(fh2AngStructpt4C20);
 495        fOutputList->Add(fh2AngStructpt1C30);
 496        fOutputList->Add(fh2AngStructpt2C30);
 497        fOutputList->Add(fh2AngStructpt3C30);
 498        fOutputList->Add(fh2AngStructpt4C30);
 499        fOutputList->Add(fh2AngStructpt1C60);
 500        fOutputList->Add(fh2AngStructpt2C60);
 501        fOutputList->Add(fh2AngStructpt3C60);
 502        fOutputList->Add(fh2AngStructpt4C60);}
 503
 504
 505
 506         fOutputList->Add(fh1TrigRef);
 507         fOutputList->Add(fh1TrigSig);
 508         fOutputList->Add(fh1TrackPhiDistance);
 509         fOutputList->Add(fh1TrackRDistance);
 510         fOutputList->Add(fh2Ntriggers);
 511         fOutputList->Add(fh2Ntriggers2C10);
 512         fOutputList->Add(fh2Ntriggers2C20);
 513         fOutputList->Add(fh3JetDensity);
 514         fOutputList->Add(fh3JetDensityA4);
 515         fOutputList->Add(fh2RPJetsC10);
 516         fOutputList->Add(fh2RPJetsC20);
 517          fOutputList->Add(fh2RPTC10);
 518         fOutputList->Add(fh2RPTC20);
 519
 520         const Int_t dimSpec = 5;
 521         const Int_t nBinsSpec[dimSpec]     = {100,6, 140, 50, fNRPBins};
 522         const Double_t lowBinSpec[dimSpec] = {0,0,-80, 0, 0};
 523         const Double_t hiBinSpec[dimSpec]  = {100,1, 200, 50,  static_cast<Double_t>(fNRPBins)};
 524         fHJetSpec = new THnSparseF("fHJetSpec","Recoil jet spectrum",dimSpec,nBinsSpec,lowBinSpec,hiBinSpec);
 525
 526              //change binning in jet area
 527      Double_t *xPt6 = new Double_t[7];
 528      xPt6[0] = 0.;
 529      xPt6[1]=0.07;
 530      xPt6[2]=0.2;
 531      xPt6[3]=0.4;
 532      xPt6[4]=0.6;
 533      xPt6[5]=0.8;
 534      xPt6[6]=1;
 535     fHJetSpec->SetBinEdges(1,xPt6);
 536     delete [] xPt6;
 537
 538
 539
 540
 541
 542         fOutputList->Add(fHJetSpec);
 543
 544
 545         if(fRunAnaAzimuthalCorrelation)
 546           {
 547             fhTTPt = new TH2F("fhTTPt","Trigger track p_{T} vs centrality",10,0,100,100,0,100);
 548             fOutputList->Add(fhTTPt);
 549
 550             const Int_t dimCor = 5;
 551             const Int_t nBinsCor[dimCor]     = {50, 200, 100,              8,   10};
 552             const Double_t lowBinCor[dimCor] = {0,  -50, -0.5*TMath::Pi(), 0,   0};
 553             const Double_t hiBinCor[dimCor]  = {50, 150, 1.5*TMath::Pi(),  0.8, 100};
 554             fHJetPhiCorr = new THnSparseF("fHJetPhiCorr","TT p_{T} vs jet p_{T} vs dPhi vs area vs centrality",dimCor,nBinsCor,lowBinCor,hiBinCor);
 555             fOutputList->Add(fHJetPhiCorr);
 556           }
 557
 558
 559
 560
 561    // =========== Switch on Sumw2 for all histos ===========
 562    for (Int_t i=0; i<fOutputList->GetEntries(); ++i) {
 563       TH1 *h1 = dynamic_cast<TH1*>(fOutputList->At(i));
 564       if (h1){
 565          h1->Sumw2();
 566          continue;
 567       }
 568     THnSparse *hn = dynamic_cast<THnSparse*>(fOutputList->At(i));
 569       if (hn){
 570          hn->Sumw2();
 571       }
 572    }
 573    TH1::AddDirectory(oldStatus);
 574
 575    PostData(1, fOutputList);
 576 }
 577
 578 void AliAnalysisTaskJetCore::UserExec(Option_t *)
 579 {
 580
 581
 582    if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
 583       AliError("Jet branch name not set.");
 584       return;
 585    }
 586
 587    fESD=dynamic_cast<AliESDEvent*>(InputEvent());
 588    if (!fESD) {
 589       AliError("ESD not available");
 590       fAODIn = dynamic_cast<AliAODEvent*>(InputEvent());
 591    }
 592       fAODOut = dynamic_cast<AliAODEvent*>(AODEvent());
 593
 594        static AliAODEvent* aod = 0;
 595        // take all other information from the aod we take the tracks from
 596        if(!aod){
 597        if(!fESD)aod = fAODIn;
 598        else aod = fAODOut;}
 599
 600
 601
 602     if(fNonStdFile.Length()!=0){
 603     // case that we have an AOD extension we need can fetch the jets from the extended output
 604     AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
 605     fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
 606     if(!fAODExtension){
 607       if(fDebug>1)Printf("AODExtension found for %s",fNonStdFile.Data());
 608     }}
 609
 610
 611    // -- event selection --
 612    fHistEvtSelection->Fill(1); // number of events before event selection
 613
 614
 615        Bool_t selected=kTRUE;
 616        selected = AliAnalysisHelperJetTasks::Selected();
 617        if(!selected){
 618       // no selection by the service task, we continue
 619        PostData(1,fOutputList);
 620        return;}
 621
 622
 623
 624   // physics selection: this is now redundant, all should appear as accepted after service task selection
 625    AliInputEventHandler* inputHandler = (AliInputEventHandler*)
 626    ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
 627          std::cout<<inputHandler->IsEventSelected()<<" "<<fOfflineTrgMask<<std::endl;
 628    if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){
 629       if(fDebug) Printf(" Trigger Selection: event REJECTED ... ");
 630       fHistEvtSelection->Fill(2);
 631       PostData(1, fOutputList);
 632       return;
 633    }
 634
 635
 636
 637    // vertex selection
 638    if(!aod){
 639      if(fDebug) Printf("%s:%d No AOD",(char*)__FILE__,__LINE__);
 640      fHistEvtSelection->Fill(3);
 641       PostData(1, fOutputList);
 642    }
 643    AliAODVertex* primVtx = aod->GetPrimaryVertex();
 644
 645    if(!primVtx){
 646      if(fDebug) Printf("%s:%d No primVtx",(char*)__FILE__,__LINE__);
 647      fHistEvtSelection->Fill(3);
 648      PostData(1, fOutputList);
 649      return;
 650    }
 651
 652    Int_t nTracksPrim = primVtx->GetNContributors();
 653    if ((nTracksPrim < fMinContribVtx) ||
 654          (primVtx->GetZ() < fVtxZMin) ||
 655          (primVtx->GetZ() > fVtxZMax) ){
 656       if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
 657       fHistEvtSelection->Fill(3);
 658       PostData(1, fOutputList);
 659       return;
 660    }
 661
 662    // event class selection (from jet helper task)
 663    Int_t eventClass = AliAnalysisHelperJetTasks::EventClass();
 664    if(fDebug) Printf("Event class %d", eventClass);
 665    if (eventClass < fEvtClassMin || eventClass > fEvtClassMax){
 666       fHistEvtSelection->Fill(4);
 667       PostData(1, fOutputList);
 668       return;
 669    }
 670
 671    // centrality selection
 672    AliCentrality *cent = 0x0;
 673    Double_t centValue = 0.;
 674    if(fIsPbPb){
 675    if(fESD) {cent = fESD->GetCentrality();
 676      if(cent) centValue = cent->GetCentralityPercentile("V0M");}
 677    else     centValue=((AliVAODHeader*)aod->GetHeader())->GetCentrality();
 678
 679    if(fDebug) printf("centrality: %f\n", centValue);
 680       if (centValue < fCentMin || centValue > fCentMax){
 681       fHistEvtSelection->Fill(4);
 682       PostData(1, fOutputList);
 683       return;
 684       }}
 685
 686
 687    fHistEvtSelection->Fill(0);
 688    // accepted events
 689    // -- end event selection --
 690
 691    // get background
 692    AliAODJetEventBackground* externalBackground = 0;
 693    if(fAODOut&&!externalBackground&&fBackgroundBranch.Length()){
 694       externalBackground =  (AliAODJetEventBackground*)(fAODOut->FindListObject(fBackgroundBranch.Data()));
 695       if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
 696    }
 697    if(fAODExtension&&!externalBackground&&fBackgroundBranch.Length()){
 698      externalBackground =  (AliAODJetEventBackground*)(fAODExtension->GetAOD()->FindListObject(fBackgroundBranch.Data()));
 699       if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
 700    }
 701
 702     if(fAODIn&&!externalBackground&&fBackgroundBranch.Length()){
 703       externalBackground =  (AliAODJetEventBackground*)(fAODIn->FindListObject(fBackgroundBranch.Data()));
 704       if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
 705     }
 706
 707    Float_t rho = 0;
 708
 709    if(fFlagRandom==0){
 710      if(externalBackground)rho = externalBackground->GetBackground(0);}
 711    if(fFlagRandom==2){
 712       if(externalBackground)rho = externalBackground->GetBackground(2);}
 713    if(fFlagRandom==3){
 714       if(externalBackground)rho = externalBackground->GetBackground(3);}
 715    // fetch jets
 716    TClonesArray *aodJets[2];
 717    aodJets[0]=0;
 718    if(fAODOut&&!aodJets[0]){
 719    aodJets[0] = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fJetBranchName[0].Data()));
 720    aodJets[1] = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fJetBranchName[1].Data()));  }
 721    if(fAODExtension && !aodJets[0]){
 722    aodJets[0] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[0].Data()));
 723    aodJets[1] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[1].Data()));  }
 724      if(fAODIn&&!aodJets[0]){
 725    aodJets[0] = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fJetBranchName[0].Data()));
 726    aodJets[1] = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fJetBranchName[1].Data()));  }
 727
 728
 729
 730    Int_t nT=0;
 731    TList ParticleList;
 732    Double_t minT=0;
 733    Double_t maxT=0;
 734    Int_t number=0;
 735    Double_t dice=fRandom->Uniform(0,1);
 736    if(dice>fFrac){ minT=fTTLowRef;
 737                    maxT=fTTUpRef;}
 738    if(dice<=fFrac){minT=fTTLowSig;
 739                    maxT=fTTUpSig;}
 740
 741
 742
 743    if(fHardest==1 || fHardest==2) nT = GetListOfTracks(&ParticleList);
 744    if(fHardest==0) nT=SelectTrigger(&ParticleList,minT,maxT,number);
 745    if(nT<0){
 746    PostData(1, fOutputList);
 747    return;}
 748
 749       if(dice>fFrac) fh1TrigRef->Fill(number);
 750       if(dice<=fFrac)fh1TrigSig->Fill(number);
 751
 752      for (Int_t iJetType = 0; iJetType < 2; iJetType++) {
 753       fListJets[iJetType]->Clear();
 754       if (!aodJets[iJetType]) continue;
 755       if(fDebug) Printf("%s: %d jets",fJetBranchName[iJetType].Data(),aodJets[iJetType]->GetEntriesFast());
 756       for (Int_t iJet = 0; iJet < aodJets[iJetType]->GetEntriesFast(); iJet++) {
 757          AliAODJet *jet = dynamic_cast<AliAODJet*>((*aodJets[iJetType])[iJet]);
 758          if (jet) fListJets[iJetType]->Add(jet);
 759          // if(iJetType==0){
 760          // ptsub[iJet]=jet->Pt()-rho*jet->EffectiveAreaCharged();}
 761       }}
 762
 763    Double_t etabig=0;
 764    Double_t ptbig=0;
 765    Double_t areabig=0;
 766    Double_t phibig=0.;
 767    Double_t etasmall=0;
 768    Double_t ptsmall=0;
 769    //   Double_t areasmall=0;
 770    Double_t phismall=0.;
 771
 772
 773
 774    Int_t iCount=0;
 775    Int_t trigJet=-1;
 776    Int_t trigBBTrack=-1;
 777           // Int_t trigInTrack=-1;
 778    fRPAngle = ((AliVAODHeader*)aod->GetHeader())->GetEventplane();
 779
 780    if(fHardest==0 || fHardest==1){
 781    AliVParticle *partback = (AliVParticle*)ParticleList.At(nT);
 782    if(!partback){
 783    PostData(1, fOutputList);
 784    return;}
 785
 786     if(fSemigoodCorrect){
 787    Double_t disthole=RelativePhi(partback->Phi(),fHolePos);
 788    if(TMath::Abs(disthole)+fHoleWidth>TMath::Pi()-0.6){
 789    PostData(1, fOutputList);
 790    return;}}
 791
 792
 793     }
 794
 795
 796    for(Int_t tt=0;tt<ParticleList.GetEntries();tt++){
 797      if(fHardest==0||fHardest==1){if(tt!=nT) continue;}
 798      AliVParticle *partback = (AliVParticle*)ParticleList.At(tt);
 799      if(!partback) continue;
 800      if(partback->Pt()<8) continue;
 801
 802    Double_t accep=2.*TMath::Pi()*1.8;
 803    Int_t injet4=0;
 804    Int_t injet=0;
 805
 806    if(fSemigoodCorrect){
 807    Double_t disthole=RelativePhi(partback->Phi(),fHolePos);
 808    if(TMath::Abs(disthole)+fHoleWidth>TMath::Pi()-0.6) continue;}
 809
 810
 811
 812    fh2Ntriggers->Fill(centValue,partback->Pt());
 813    Double_t phiBinT = RelativePhi(partback->Phi(),fRPAngle);
 814    if(centValue<20.) fh2RPTC20->Fill(TMath::Abs(phiBinT),partback->Pt());
 815    if(centValue<10.) fh2RPTC10->Fill(TMath::Abs(phiBinT),partback->Pt());
 816
 817
 818
 819    for(Int_t i=0; i<fListJets[0]->GetEntries(); ++i){
 820            AliAODJet* jetbig = (AliAODJet*)(fListJets[0]->At(i));
 821            etabig  = jetbig->Eta();
 822            phibig  = jetbig->Phi();
 823            ptbig   = jetbig->Pt();
 824            if(ptbig==0) continue;
 825            Double_t phiBin = RelativePhi(phibig,fRPAngle);
 826            areabig = jetbig->EffectiveAreaCharged();
 827            Double_t ptcorr=ptbig-rho*areabig;
 828            if((etabig<fJetEtaMin)||(etabig>fJetEtaMax)) continue;
 829            if(areabig>=0.07) injet=injet+1;
 830            if(areabig>=0.4) injet4=injet4+1;
 831            Double_t dphi=RelativePhi(partback->Phi(),phibig);
 832
 833            if(fFlagEtaBkg==1){
 834            Double_t etadif= partback->Eta()-etabig;
 835            if(TMath::Abs(etadif)<=0.5){
 836
 837            if(centValue<20.) fh3JetTrackC20->Fill(partback->Pt(),ptcorr,TMath::Abs(dphi));
 838            if(centValue>30. && centValue<60.) fh3JetTrackC3060->Fill(partback->Pt(),ptcorr,TMath::Abs(dphi));}}
 839            if(fFlagEtaBkg==0){
 840            if(centValue<20.) fh3JetTrackC20->Fill(partback->Pt(),ptcorr,TMath::Abs(dphi));
 841            if(centValue>30. && centValue<60.) fh3JetTrackC3060->Fill(partback->Pt(),ptcorr,TMath::Abs(dphi));}
 842
 843
 844            if(fFlagJetHadron==0){
 845            if(fFlagPhiBkg==1) if((TMath::Abs(dphi)<TMath::Pi()/2.-0.1)||(TMath::Abs(dphi)>TMath::Pi()/2.+0.1)) continue;
 846            if(fFlagPhiBkg==0) if(TMath::Abs(dphi)<TMath::Pi()-0.6) continue;
 847            if(fFlagPhiBkg==2) if(TMath::Abs(dphi)<TMath::Pi()-0.7) continue;
 848            if(fFlagPhiBkg==3) if(TMath::Abs(dphi)<TMath::Pi()-0.5) continue;}
 849
 850            if(fFlagJetHadron!=0) if(TMath::Abs(dphi)>0.4) continue;
 851
 852
 853            if(centValue<10.) fh2RPJetsC10->Fill(TMath::Abs(phiBin), ptcorr);
 854            if(centValue<20.) fh2RPJetsC20->Fill(TMath::Abs(phiBin), ptcorr);
 855                    Double_t dismin=100.;
 856                    Double_t ptmax=-10.;
 857                    Int_t index1=-1;
 858                    Int_t index2=-1;
 859
 860                    Float_t phitt=partback->Phi();
 861                    if(phitt<0)phitt+=TMath::Pi()*2.;
 862                    Int_t phiBintt = GetPhiBin(phitt-fRPAngle);
 863
 864                    Double_t fillspec[] = {centValue,jetbig->EffectiveAreaCharged(),ptcorr,partback->Pt(), static_cast<Double_t>(phiBintt)};
 865                   fHJetSpec->Fill(fillspec);
 866
 867
 868
 869                    if(ptcorr<=0) continue;
 870
 871                        AliAODTrack* leadtrack=0;
 872                        Int_t ippt=0;
 873                        Double_t ppt=-10;
 874                        if(fFlagJetHadron==0){
 875                          TRefArray *genTrackList = jetbig->GetRefTracks();
 876                          Int_t nTracksGenJet = genTrackList->GetEntriesFast();
 877                          AliAODTrack* genTrack;
 878                          for(Int_t ir=0; ir<nTracksGenJet; ++ir){
 879                            genTrack = (AliAODTrack*)(genTrackList->At(ir));
 880                            if(genTrack->Pt()>ppt){ppt=genTrack->Pt();
 881                              ippt=ir;}}
 882                          leadtrack=(AliAODTrack*)(genTrackList->At(ippt));
 883                          if(!leadtrack) continue;
 884                        }
 885
 886
 887
 888                        AliVParticle* leadtrackb=0;
 889                        if(fFlagJetHadron!=0){
 890                          Int_t nTb = GetHardestTrackBackToJet(jetbig);
 891                          leadtrackb = (AliVParticle*)ParticleList.At(nTb);
 892                          if(!leadtrackb) continue;
 893                        }
 894
 895
 896
 897
 898
 899                        //store one trigger info
 900                        if(iCount==0){
 901                          trigJet=i;
 902                          trigBBTrack=nT;
 903                          //                      trigInTrack=ippt;
 904                          iCount=iCount+1;}
 905
 906
 907                        if(fCheckMethods){
 908                          for(Int_t j=0; j<fListJets[1]->GetEntries(); ++j){
 909                            AliAODJet* jetsmall = (AliAODJet*)(fListJets[1]->At(j));
 910                            etasmall  = jetsmall->Eta();
 911                            phismall = jetsmall->Phi();
 912                            ptsmall   = jetsmall->Pt();
 913                            //                      areasmall = jetsmall->EffectiveAreaCharged();
 914                            Double_t tmpDeltaR=(phismall-phibig)*(phismall-phibig)+(etasmall-etabig)*(etasmall-etabig);
 915                            tmpDeltaR=TMath::Sqrt(tmpDeltaR);
 916                            //Fraction in the jet core
 917                            if((ptsmall>ptmax)&&(tmpDeltaR<=fRadioFrac)){ptmax=ptsmall;
 918                              index2=j;}
 919                     if(tmpDeltaR<=dismin){ dismin=tmpDeltaR;
 920                       index1=j;}} //en of loop over R=0.2 jets
 921                 //method1:most concentric jet=core
 922                 if(dismin<fMinDist){ AliAODJet* jetmethod1 = (AliAODJet*)(fListJets[1]->At(index1));
 923                   if(centValue<10) fh2JetCoreMethod1C10->Fill(ptcorr,jetmethod1->Pt()/ptbig);
 924                   if((centValue>20)&&(centValue<40)) fh2JetCoreMethod1C20->Fill(ptcorr,jetmethod1->Pt()/ptbig);
 925                   if((centValue>30)&&(centValue<60)) fh2JetCoreMethod1C30->Fill(ptcorr,jetmethod1->Pt()/ptbig);
 926                   if(centValue>60) fh2JetCoreMethod1C60->Fill(ptcorr,jetmethod1->Pt()/ptbig); }
 927                 //method2:hardest contained jet=core
 928                 if(index2!=-1){
 929                   AliAODJet* jetmethod2 = (AliAODJet*)(fListJets[1]->At(index2));
 930                   if(centValue<10) fh2JetCoreMethod2C10->Fill(ptcorr,jetmethod2->Pt()/ptbig);
 931                   if((centValue>20)&&(centValue<40)) fh2JetCoreMethod2C20->Fill(ptcorr,jetmethod2->Pt()/ptbig);
 932                   if((centValue>30)&&(centValue<60)) fh2JetCoreMethod2C30->Fill(ptcorr,jetmethod2->Pt()/ptbig);
 933                   if(centValue>60) fh2JetCoreMethod2C60->Fill(ptcorr,jetmethod2->Pt()/ptbig); }}
 934                   if(centValue<10&&leadtrack) fh2Ntriggers2C10->Fill(leadtrack->Pt(),partback->Pt());
 935                   if(centValue<20&&leadtrack) fh2Ntriggers2C20->Fill(leadtrack->Pt(),partback->Pt());
 936          if(fDoEventMixing==0 && fFlagOnlyRecoil==0){
 937          for(int it = 0;it<ParticleList.GetEntries();++it){
 938           AliVParticle *part = (AliVParticle*)ParticleList.At(it);
 939           Double_t deltaR = jetbig->DeltaR(part);
 940           Double_t deltaEta = etabig-part->Eta();
 941
 942           Double_t deltaPhi=phibig-part->Phi();
 943           if(deltaPhi<-0.5*TMath::Pi()) deltaPhi+=2.*TMath::Pi();
 944           if(deltaPhi>3./2.*TMath::Pi()) deltaPhi-=2.*TMath::Pi();
 945           Double_t pTcont=0;
 946           if(fFlagJetHadron==0) pTcont=leadtrack->Pt();
 947           if(fFlagJetHadron!=0) pTcont=leadtrackb->Pt();
 948            Double_t jetEntries[8] = {centValue,ptcorr,part->Pt(),deltaR,deltaEta,deltaPhi,pTcont,partback->Pt()};
 949            fhnDeltaR->Fill(jetEntries);}
 950
 951
 952           }
 953
 954
 955           //end of track loop, we only do it if EM is switched off
 956
 957
 958
 959
 960
 961
 962
 963
 964
 965    }
 966    if(injet>0) fh3JetDensity->Fill(ParticleList.GetEntries(),injet/accep,partback->Pt());
 967    if(injet4>0)fh3JetDensityA4->Fill(ParticleList.GetEntries(),injet4/accep,partback->Pt());
 968           //end of jet loop
 969
 970    //}
 971
 972
 973           if(fDoEventMixing>0){
 974             //check before if the trigger exists
 975             // fTrigBuffer[i][0] = zvtx
 976             // fTrigBuffer[i][1] = phi
 977             // fTrigBuffer[i][2] = eta
 978             // fTrigBuffer[i][3] = pt_jet
 979             // fTrigBuffer[i][4] = pt_trig
 980             // fTrigBuffer[i][5]= centrality
 981             if(fTindex==10) fTindex=0;
 982             if(fTrigBuffer[fTindex][3]>0){
 983             if (TMath::Abs(fTrigBuffer[fTindex][0]-primVtx->GetZ()<2.)){
 984             if (TMath::Abs(fTrigBuffer[fTindex][5]-centValue<5)){
 985
 986                         for(int it = 0;it<nT;++it){
 987                         AliVParticle *part = (AliVParticle*)ParticleList.At(it);
 988                         Double_t DPhi = fTrigBuffer[fTindex][1] - part->Phi();
 989                         Double_t DEta = fTrigBuffer[fTindex][2] - part->Eta();
 990                         Double_t DR=TMath::Sqrt(DPhi*DPhi+DEta*DEta);
 991                         if(DPhi<-0.5*TMath::Pi()) DPhi+=2.*TMath::Pi();
 992                         if(DPhi>3./2.*TMath::Pi()) DPhi-=2.*TMath::Pi();
 993                         Double_t triggerEntries[7] = {centValue,fTrigBuffer[fTindex][3],part->Pt(),DR,DEta,DPhi,fTrigBuffer[fTindex][4]};
 994                         fhnMixedEvents->Fill(triggerEntries);
 995                         }
 996                         fNevents=fNevents+1;
 997                         if(fNevents==10) fTindex=fTindex+1;
 998             }}}
 999
1000                if(fTindex==10&&fNevents==10) fCountAgain=0;
1001
1002                // Copy the triggers from the current event into the buffer.
1003                //again, only if the trigger exists:
1004                if(fCountAgain==0){
1005                 if(trigJet>-1){
1006                 AliAODJet* jetT = (AliAODJet*)(fListJets[0]->At(trigJet));                      AliVParticle *partT = (AliVParticle*)ParticleList.At(trigBBTrack);
1007                 fTrigBuffer[fTrigBufferIndex][0] = primVtx->GetZ();
1008                 fTrigBuffer[fTrigBufferIndex][1] = jetT->Phi();
1009                 fTrigBuffer[fTrigBufferIndex][2] = jetT->Eta();
1010                 fTrigBuffer[fTrigBufferIndex][3] = jetT->Pt()-rho*jetT->EffectiveAreaCharged();
1011                 fTrigBuffer[fTrigBufferIndex][4] = partT->Pt();
1012                 fTrigBuffer[fTrigBufferIndex][5] = centValue;
1013                 fTrigBufferIndex++;
1014                 if(fTrigBufferIndex==9) {fTrigBufferIndex=0;
1015                                          fCountAgain=1;}
1016                 }
1017                }
1018
1019           }
1020
1021           /////////////////////////////////////////////////////////////////////////////
1022           ////////////////////// Rongrong's analysis //////////////////////////////////
1023           if(fRunAnaAzimuthalCorrelation)
1024             {
1025               fhTTPt->Fill(centValue,partback->Pt());
1026               for(Int_t ij=0; ij<fListJets[0]->GetEntries(); ij++)
1027                 {
1028                   AliAODJet* jet = (AliAODJet*)(fListJets[0]->At(ij));
1029                   Double_t jetPt   = jet->Pt();
1030                   Double_t jetEta  = jet->Eta();
1031                   Double_t jetPhi  = jet->Phi();
1032                   if(jetPt==0) continue;
1033                   if((jetEta<fJetEtaMin)||(jetEta>fJetEtaMax)) continue;
1034                   Double_t jetArea = jet->EffectiveAreaCharged();
1035                   Double_t jetPtCorr=jetPt-rho*jetArea;
1036                   Double_t dPhi=jetPhi-partback->Phi();
1037                   if(dPhi>2*TMath::Pi()) dPhi -= 2*TMath::Pi();
1038                   if(dPhi<-2*TMath::Pi()) dPhi += 2*TMath::Pi();
1039                   if(dPhi<-0.5*TMath::Pi()) dPhi += 2*TMath::Pi();
1040                   if(dPhi>1.5*TMath::Pi()) dPhi -= 2*TMath::Pi();
1041
1042                   Double_t fill[] = {partback->Pt(),jetPtCorr,dPhi,jetArea,centValue};
1043                   fHJetPhiCorr->Fill(fill);
1044                 }
1045             }
1046           /////////////////////////////////////////////////////////////////////////////
1047           /////////////////////////////////////////////////////////////////////////////
1048
1049
1050      //////////////////ANGULAR STRUCTURE//////////////////////////////////////
1051
1052      //tracks up to R=0.8 distant from the jet axis
1053    //   if(fAngStructCloseTracks==1){
1054    //    TList CloseTrackList;
1055    //    Int_t nn=GetListOfTracksCloseToJet(&CloseTrackList,jetbig);
1056    //    Double_t difR=0.04;
1057    //    for(Int_t l=0;l<15;l++){
1058    //    Double_t rr=l*0.1+0.1;
1059    //     for(int it = 0;it<nn;++it){
1060    //         AliVParticle *part1 = (AliVParticle*)CloseTrackList.At(it);
1061    //         for(int itu=it+1;itu<CloseTrackList.GetEntries();itu++){
1062    //         AliVParticle *part2 = (AliVParticle*)CloseTrackList.At(itu);
1063    //         Double_t ptm=part1->Pt();
1064    //         Double_t ptn=part2->Pt();
1065    //         Double_t Rnm = (part1->Eta()-part2->Eta())*(part1->Eta()-part2->Eta())+(part1->Phi()-part2->Phi())*(part1->Phi()-part2->Phi());
1066    //                    Rnm=TMath::Sqrt(Rnm);
1067    //                    Double_t deltag=(1./(TMath::Sqrt(2*TMath::Pi())*difR))*TMath::Exp(-1.*(rr-Rnm)*(rr-Rnm)/(2.*difR*difR));
1068    //                    Double_t stepf=0.5*(1.+TMath::Erf((rr-Rnm)/(TMath::Sqrt(2.)*difR)));
1069    //                    if((ptcorr<85.) && (ptcorr>=70.)){up1[l]=up1[l]+ptm*ptn*Rnm*Rnm*deltag;
1070    //                                                   down1[l]=down1[l]+ptm*ptn*Rnm*Rnm*stepf;}
1071    //                    if((ptcorr<100.) && (ptcorr>=85.)){up2[l]=up2[l]+ptm*ptn*Rnm*Rnm*deltag;
1072    //                                                   down2[l]=down2[l]+ptm*ptn*Rnm*Rnm*stepf;}
1073    //                    if((ptcorr<120.) && (ptcorr>=100.)){up3[l]=up3[l]+ptm*ptn*Rnm*Rnm*deltag;
1074    //                                                    down3[l]=down3[l]+ptm*ptn*Rnm*Rnm*stepf;}
1075    //                    if((ptcorr<140.) && (ptcorr>=120.)){up4[l]=up4[l]+ptm*ptn*Rnm*Rnm*deltag;
1076    //                   down4[l]=down4[l]+ptm*ptn*Rnm*Rnm*stepf;}}}}
1077    //   }
1078
1079    //   //only jet constituents
1080    //    if(fAngStructCloseTracks==2){
1081
1082    //    Double_t difR=0.04;
1083    //    for(Int_t l=0;l<15;l++){
1084    //    Double_t rr=l*0.1+0.1;
1085
1086
1087    //    AliAODTrack* part1;
1088    //    AliAODTrack* part2;
1089
1090    //        TRefArray *genTrackListb = jetbig->GetRefTracks();
1091    //        Int_t nTracksGenJetb = genTrackListb->GetEntriesFast();
1092
1093
1094
1095    //        for(Int_t it=0; it<nTracksGenJetb; ++it){
1096    //           part1 = (AliAODTrack*)(genTrackListb->At(it));
1097    //         for(Int_t itu=0; itu<nTracksGenJetb; ++itu){
1098    //           part2 = (AliAODTrack*)(genTrackListb->At(itu));
1099    //         Double_t ptm=part1->Pt();
1100    //         Double_t ptn=part2->Pt();
1101    //         Double_t Rnm = (part1->Eta()-part2->Eta())*(part1->Eta()-part2->Eta())+(part1->Phi()-part2->Phi())*(part1->Phi()-part2->Phi());
1102    //                    Rnm=TMath::Sqrt(Rnm);
1103    //                    Double_t deltag=(1./(TMath::Sqrt(2*TMath::Pi())*difR))*TMath::Exp(-1.*(rr-Rnm)*(rr-Rnm)/(2.*difR*difR));
1104    //                    Double_t stepf=0.5*(1.+TMath::Erf((rr-Rnm)/(TMath::Sqrt(2.)*difR)));
1105    //                    if((ptcorr<85.) && (ptcorr>=70.)){up1[l]=up1[l]+ptm*ptn*Rnm*Rnm*deltag;
1106    //                                                   down1[l]=down1[l]+ptm*ptn*Rnm*Rnm*stepf;}
1107    //                    if((ptcorr<100.) && (ptcorr>=85.)){up2[l]=up2[l]+ptm*ptn*Rnm*Rnm*deltag;
1108    //                                                   down2[l]=down2[l]+ptm*ptn*Rnm*Rnm*stepf;}
1109    //                    if((ptcorr<120.) && (ptcorr>=100.)){up3[l]=up3[l]+ptm*ptn*Rnm*Rnm*deltag;
1110    //                                                    down3[l]=down3[l]+ptm*ptn*Rnm*Rnm*stepf;}
1111    //                    if((ptcorr<140.) && (ptcorr>=120.)){up4[l]=up4[l]+ptm*ptn*Rnm*Rnm*deltag;
1112    //                   down4[l]=down4[l]+ptm*ptn*Rnm*Rnm*stepf;}}}}}
1113    // }
1114      // //end loop over R=0.4 jets
1115      // if(fAngStructCloseTracks>0){
1116      // for(Int_t l=0;l<15;l++){
1117      // Double_t rr=l*0.1+0.1;
1118      //    if(down1[l]!=0){
1119      //         if(centValue<10.)fh2AngStructpt1C10->Fill(rr,rr*up1[l]/down1[l]);
1120      //    if(centValue>20. && centValue<40.) fh2AngStructpt1C20->Fill(rr,rr*up1[l]/down1[l]);
1121      //    if(centValue>30. && centValue<60.) fh2AngStructpt1C30->Fill(rr,rr*up1[l]/down1[l]);
1122      //    if(centValue>60.) fh2AngStructpt1C60->Fill(rr,rr*up1[l]/down1[l]);}
1123      //    if(down2[l]!=0){
1124      //         if(centValue<10.) fh2AngStructpt2C10->Fill(rr,rr*up2[l]/down2[l]);
1125      //    if(centValue>20. && centValue<40.) fh2AngStructpt2C20->Fill(rr,rr*up2[l]/down2[l]);
1126      //    if(centValue>30. && centValue<60.) fh2AngStructpt2C30->Fill(rr,rr*up2[l]/down2[l]);
1127      //    if(centValue>60.) fh2AngStructpt2C60->Fill(rr,rr*up2[l]/down2[l]);}
1128      //    if(down3[l]!=0){
1129      //         if(centValue<10.) fh2AngStructpt3C10->Fill(rr,rr*up3[l]/down3[l]);
1130      //    if(centValue>20. && centValue<40.) fh2AngStructpt3C20->Fill(rr,rr*up3[l]/down3[l]);
1131      //    if(centValue>30. && centValue<60.) fh2AngStructpt3C30->Fill(rr,rr*up3[l]/down3[l]);
1132      //    if(centValue>60.) fh2AngStructpt3C60->Fill(rr,rr*up3[l]/down3[l]);}
1133      //    if(down4[l]!=0){
1134      //         if(centValue<10.) fh2AngStructpt4C10->Fill(rr,rr*up4[l]/down4[l]);
1135      //    if(centValue>20. && centValue<40.) fh2AngStructpt4C20->Fill(rr,rr*up4[l]/down4[l]);
1136      //    if(centValue>30. && centValue<60.) fh2AngStructpt4C30->Fill(rr,rr*up4[l]/down4[l]);
1137      //    if(centValue>60.) fh2AngStructpt4C60->Fill(rr,rr*up4[l]/down4[l]);}}}
1138
1139
1140
1141    }
1142
1143
1144
1145    PostData(1, fOutputList);
1146 }
1147
1148 void AliAnalysisTaskJetCore::Terminate(const Option_t *)
1149 {
1150    // Draw result to the screen
1151    // Called once at the end of the query
1152
1153    if (!GetOutputData(1))
1154    return;
1155 }
1156
1157
1158
1159
1160
1161
1162 Int_t  AliAnalysisTaskJetCore::GetListOfTracks(TList *list){
1163
1164       Int_t iCount = 0;
1165       AliAODEvent *aod = 0;
1166
1167      if(!fESD)aod = fAODIn;
1168      else aod = fAODOut;
1169
1170      if(!aod)return 0;
1171
1172      Int_t index=-1;
1173      Double_t ptmax=-10;
1174
1175
1176
1177      for(int it = 0;it < aod->GetNumberOfTracks();++it){
1178       AliAODTrack *tr = dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
1179       if(!tr) AliFatal("Not a standard AOD");
1180       Bool_t bGood = false;
1181       if(fFilterType == 0)bGood = true;
1182       else if(fFilterType == 1)bGood = tr->IsHybridTPCConstrainedGlobal();
1183       else if(fFilterType == 2)bGood = tr->IsHybridGlobalConstrainedGlobal();
1184       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
1185       if(fRequireITSRefit==1){if((tr->GetStatus()&AliESDtrack::kITSrefit)==0)continue;}
1186       if(bGood==false) continue;
1187       if (fApplySharedClusterCut) {
1188            Double_t frac = Double_t(tr->GetTPCnclsS()) /Double_t(tr->GetTPCncls());
1189            if (frac > 0.4) continue;
1190       }
1191      if(TMath::Abs(tr->Eta())>0.9)continue;
1192       if(tr->Pt()<0.15)continue;
1193       list->Add(tr);
1194       iCount++;
1195       if(fFilterType==2 && fFilterMaskBestPt>0){// only set the trigger track index for good quality tracks
1196         if(tr->TestFilterBit(fFilterMaskBestPt)){
1197           if(tr->Pt()>ptmax){
1198             ptmax=tr->Pt();
1199             index=iCount-1;
1200           }
1201         }
1202       }
1203       else{
1204         if(tr->Pt()>ptmax){
1205           ptmax=tr->Pt();
1206           index=iCount-1;
1207         }
1208       }
1209      }
1210
1211
1212     // else if (type == kTrackAODMCCharged) {
1213     // TClonesArray *tca = dynamic_cast<TClonesArray*>(aod->FindListObject(AliAODMCParticle::StdBranchName()));
1214     // if(!tca)return iCount;
1215     // for(int it = 0;it < tca->GetEntriesFast();++it){
1216     //   AliAODMCParticle *part = dynamic_cast<AliAODMCParticle*>(tca->At(it));
1217     //   if(!part)continue;
1218     //   if(part->Pt()<0.15)continue;
1219     //   if(!part->IsPhysicalPrimary())continue;
1220     //   if(part->Charge()==0)continue;
1221     //   if(TMath::Abs(part->Eta())>0.9)continue;
1222     //   list->Add(part);
1223     //   iCount++;
1224     //   if(part->Pt()>ptmax){ ptmax=part->Pt();
1225     //  index=iCount-1;}}}
1226       return index;
1227
1228 }
1229
1230
1231
1232 Int_t  AliAnalysisTaskJetCore::SelectTrigger(TList *list,Double_t minT,Double_t maxT,Int_t &number){
1233      Int_t iCount = 0;
1234      AliAODEvent *aod = 0;
1235      if(!fESD)aod = fAODIn;
1236      else aod = fAODOut;
1237      if(!aod)return 0;
1238      Int_t index=-1;
1239      Int_t triggers[100];
1240      Int_t triggers2[100];
1241      for(Int_t cr=0;cr<100;cr++){triggers[cr]=-1;
1242        triggers2[cr]=-1;}
1243      Int_t im=0;
1244      Int_t im2=0;
1245      for(int it = 0;it < aod->GetNumberOfTracks();++it){
1246       AliAODTrack *tr = dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
1247       if(!tr) {
1248         AliFatal("Not a standard AOD");
1249         return -1;
1250       }
1251       Bool_t bGood = false;
1252       if(fFilterType == 0)bGood = true;
1253       else if(fFilterType == 1)bGood = tr->IsHybridTPCConstrainedGlobal();
1254       else if(fFilterType == 2)bGood = tr->IsHybridGlobalConstrainedGlobal();
1255       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
1256        if(fRequireITSRefit==1){if((tr->GetStatus()&AliESDtrack::kITSrefit)==0)continue;}
1257       if(bGood==false) continue;
1258       if (fApplySharedClusterCut) {
1259            Double_t frac = Double_t(tr->GetTPCnclsS()) /Double_t(tr->GetTPCncls());
1260            if (frac > 0.4) continue;
1261       }
1262       if(TMath::Abs(tr->Eta())>0.9)continue;
1263       if(tr->Pt()<0.15)continue;
1264       list->Add(tr);
1265       iCount++;
1266
1267       if(tr->Pt()>=minT && tr->Pt()<maxT){
1268         if(tr->Pt()<20){
1269         triggers[im]=iCount-1;
1270         im=im+1;}
1271         if(tr->Pt()>=20.){triggers2[im2]=iCount-1;
1272           im2=im2+1;}
1273
1274       }}
1275
1276       number=im;
1277       Int_t rd=0;
1278       if(im2==0) rd=0;
1279       if(im2>0) rd=fRandom->Integer(im2);
1280       index=triggers2[rd];
1281       if(index==-1) return index;
1282       AliVParticle *tr1 = (AliVParticle*)list->At(index);
1283       if(!tr1) return index;
1284
1285
1286       for(Int_t kk=0;kk<im;kk++){
1287         //if(kk==rd) continue;
1288         if(index==triggers[kk]) continue;
1289         Int_t lab=triggers[kk];
1290          AliVParticle *tr2 = (AliVParticle*)list->At(lab);
1291          if(!tr2) continue;
1292        Double_t detat=tr1->Eta()-tr2->Eta();
1293        Double_t dphit=RelativePhi(tr1->Phi(),tr2->Phi());
1294        Double_t deltaRt=TMath::Sqrt(detat*detat+dphit*dphit);
1295        fh1TrackPhiDistance->Fill(TMath::Abs(dphit));
1296        fh1TrackRDistance->Fill(deltaRt);
1297
1298        if(fDodiHadron==1)  if(deltaRt>0.4) number=number-1;
1299        if(fDodiHadron==2) if((deltaRt>0.4) && (TMath::Abs(dphit)>TMath::Pi()-0.6)) number=number-1;}
1300
1301
1302
1303
1304
1305
1306       return index;
1307
1308 }
1309
1310
1311
1312
1313
1314
1315
1316
1317
1318
1319
1320
1321
1322
1323
1324
1325    Int_t  AliAnalysisTaskJetCore::GetHardestTrackBackToJet(AliAODJet *jetbig){
1326
1327     AliAODEvent *aod = 0;
1328     if(!fESD)aod = fAODIn;
1329     else aod = fAODOut;
1330     Int_t index=-1;
1331     Double_t ptmax=-10;
1332     Double_t dphi=0;
1333     // Double_t dif=0;
1334     Int_t iCount=0;
1335     for(int it = 0;it < aod->GetNumberOfTracks();++it){
1336       AliAODTrack *tr = dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
1337       if(!tr) AliFatal("Not a standard AOD");
1338       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
1339       if(TMath::Abs(tr->Eta())>0.9)continue;
1340       if(tr->Pt()<0.15)continue;
1341       iCount=iCount+1;
1342       dphi=RelativePhi(tr->Phi(),jetbig->Phi());
1343       if(TMath::Abs(dphi)<TMath::Pi()-0.6) continue;
1344       if(tr->Pt()>ptmax){ ptmax=tr->Pt();
1345       index=iCount-1;
1346       //  dif=dphi;
1347       }}
1348
1349       return index;
1350
1351    }
1352
1353
1354
1355
1356
1357
1358
1359
1360
1361  Int_t  AliAnalysisTaskJetCore::GetListOfTracksCloseToJet(TList *list,AliAODJet *jetbig){
1362
1363     Int_t iCount = 0;
1364       AliAODEvent *aod = 0;
1365      if(!fESD)aod = fAODIn;
1366      else aod = fAODOut;
1367
1368       for(int it = 0;it < aod->GetNumberOfTracks();++it){
1369       AliAODTrack *tr = dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
1370       if(!tr) AliFatal("Not a standard AOD");
1371       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
1372       if(TMath::Abs(tr->Eta())>0.9)continue;
1373       if(tr->Pt()<0.15)continue;
1374       Double_t disR=jetbig->DeltaR(tr);
1375       if(disR>0.8)  continue;
1376       list->Add(tr);
1377       //cout<<fAOD->GetNumberOfTracks()<<" "<<tr->Pt()<<endl;
1378       iCount++;
1379     }
1380
1381    list->Sort();
1382    return iCount;
1383
1384 }
1385
1386
1387
1388
1389
1390
1391
1392
1393
1394
1395
1396 Int_t AliAnalysisTaskJetCore::GetNInputTracks()
1397 {
1398
1399    Int_t nInputTracks = 0;
1400      AliAODEvent *aod = 0;
1401      if(!fESD)aod = fAODIn;
1402      else aod = fAODOut;
1403    TString jbname(fJetBranchName[1]);
1404    //needs complete event, use jets without background subtraction
1405    for(Int_t i=1; i<=3; ++i){
1406       if(jbname.Contains(Form("B%d",i))) jbname.ReplaceAll(Form("B%d",i),"B0");
1407    }
1408    // use only HI event
1409    if(jbname.Contains("AODextraonly")) jbname.ReplaceAll("AODextraonly","AOD");
1410    if(jbname.Contains("AODextra")) jbname.ReplaceAll("AODextra","AOD");
1411
1412    if(fDebug) Printf("Multiplicity from jet branch %s", jbname.Data());
1413    TClonesArray *tmpAODjets = dynamic_cast<TClonesArray*>(aod->FindListObject(jbname.Data()));
1414    if(!tmpAODjets){
1415       Printf("Jet branch %s not found", jbname.Data());
1416       Printf("AliAnalysisTaskJetCore::GetNInputTracks FAILED");
1417       return -1;
1418    }
1419
1420    for (Int_t iJet=0; iJet<tmpAODjets->GetEntriesFast(); iJet++){
1421       AliAODJet *jet = dynamic_cast<AliAODJet*>((*tmpAODjets)[iJet]);
1422       if(!jet) continue;
1423       TRefArray *trackList = jet->GetRefTracks();
1424       Int_t nTracks = trackList->GetEntriesFast();
1425       nInputTracks += nTracks;
1426       if(fDebug) Printf("#jet%d: %d tracks", iJet, nTracks);
1427    }
1428    if(fDebug) Printf("---> input tracks: %d", nInputTracks);
1429
1430    return nInputTracks;
1431 }
1432
1433
1434
1435 Double_t AliAnalysisTaskJetCore::RelativePhi(Double_t mphi,Double_t vphi){
1436
1437   if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi());
1438   else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
1439   if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi());
1440   else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
1441   double dphi = mphi-vphi;
1442   if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
1443   else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
1444   return dphi;//dphi in [-Pi, Pi]
1445 }
1446
1447 Int_t AliAnalysisTaskJetCore::GetPhiBin(Double_t phi)
1448 {
1449     Int_t phibin=-1;
1450     if(!(TMath::Abs(phi)<=2*TMath::Pi())){AliError("phi w.r.t. RP out of defined range");return -1;}
1451     Double_t phiwrtrp=TMath::ACos(TMath::Abs(TMath::Cos(phi)));
1452     phibin=Int_t(fNRPBins*phiwrtrp/(0.5*TMath::Pi()));
1453     if(phibin<0||phibin>=fNRPBins){AliError("Phi Bin not defined");}
1454     return phibin;
1455 }
1456
1457
1458
1459
1460 THnSparse* AliAnalysisTaskJetCore::NewTHnSparseF(const char* name, UInt_t entries)
1461 {
1462    // generate new THnSparseF, axes are defined in GetDimParams()
1463
1464    Int_t count = 0;
1465    UInt_t tmp = entries;
1466    while(tmp!=0){
1467       count++;
1468       tmp = tmp &~ -tmp;  // clear lowest bit
1469    }
1470
1471    TString hnTitle(name);
1472    const Int_t dim = count;
1473    Int_t nbins[dim];
1474    Double_t xmin[dim];
1475    Double_t xmax[dim];
1476
1477    Int_t i=0;
1478    Int_t c=0;
1479    while(c<dim && i<32){
1480       if(entries&(1<<i)){
1481
1482          TString label("");
1483          GetDimParams(i, label, nbins[c], xmin[c], xmax[c]);
1484          hnTitle += Form(";%s",label.Data());
1485          c++;
1486       }
1487
1488       i++;
1489    }
1490    hnTitle += ";";
1491
1492    return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
1493 }
1494
1495 void AliAnalysisTaskJetCore::GetDimParams(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
1496 {
1497    // stores label and binning of axis for THnSparse
1498
1499    const Double_t pi = TMath::Pi();
1500
1501    switch(iEntry){
1502
1503    case 0:
1504       label = "V0 centrality (%)";
1505
1506          nbins = 10;
1507          xmin = 0.;
1508          xmax = 100.;
1509          break;
1510
1511
1512    case 1:
1513       label = "corrected jet pt";
1514          nbins = 20;
1515          xmin = 0.;
1516          xmax = 200.;
1517           break;
1518
1519
1520    case 2:
1521       label = "track pT";
1522
1523          nbins = 9;
1524          xmin = 0.;
1525          xmax = 150;
1526          break;
1527
1528
1529     case 3:
1530       label = "deltaR";
1531       nbins = 15;
1532       xmin = 0.;
1533       xmax = 1.5;
1534       break;
1535
1536
1537
1538    case 4:
1539       label = "deltaEta";
1540       nbins = 8;
1541       xmin = -1.6;
1542       xmax = 1.6;
1543       break;
1544
1545
1546   case 5:
1547       label = "deltaPhi";
1548       nbins = 90;
1549       xmin = -0.5*pi;
1550       xmax = 1.5*pi;
1551       break;
1552
1553
1554
1555     case 6:
1556       label = "leading track";
1557       nbins = 13;
1558       xmin = 0;
1559       xmax = 50;
1560       break;
1561
1562      case 7:
1563
1564       label = "trigger track";
1565       nbins =10;
1566       xmin = 0;
1567       xmax = 50;
1568       break;
1569
1570
1571
1572
1573
1574
1575
1576
1577    }
1578
1579 }
1580