PWGJE/AliAnalysisTaskJetCore.cxx

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