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 fAOD(0x0),
  61 fAODExtension(0x0),
  62 fBackgroundBranch(""),
  63 fNonStdFile(""),
  64 fIsPbPb(kTRUE),
  65 fOfflineTrgMask(AliVEvent::kAny),
  66 fMinContribVtx(1),
  67 fVtxZMin(-10.),
  68 fVtxZMax(10.),
  69 fEvtClassMin(0),
  70 fEvtClassMax(4),
  71 fFilterMask(0),
  72 fRadioFrac(0.2),
  73 fMinDist(0.1),
  74 fCentMin(0.),
  75 fCentMax(100.),
  76 fNInputTracksMin(0),
  77 fNInputTracksMax(-1),
  78 fAngStructCloseTracks(0),
  79 fCheckMethods(0),
  80 fDoEventMixing(0),
  81 fJetEtaMin(-.5),
  82 fJetEtaMax(.5),
  83 fNevents(0x0),
  84 fTindex(0x0),
  85 fTrigBufferIndex(0x0),
  86 fJetPtMin(20.),
  87 fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
  88 fJetPtFractionMin(0.5),
  89 fNMatchJets(4),
  90 fMatchMaxDist(0.8),
  91 fKeepJets(kFALSE),
  92 fkNbranches(2),
  93 fkEvtClasses(12),
  94 fOutputList(0x0),
  95 fbEvent(kTRUE),
  96 fHistEvtSelection(0x0),
  97 fhnDeltaR(0x0),
  98 fhnMixedEvents(0x0),
  99 fh2JetCoreMethod1C10(0x0),
 100 fh2JetCoreMethod2C10(0x0),
 101 fh2JetCoreMethod1C20(0x0),
 102 fh2JetCoreMethod2C20(0x0),
 103 fh2JetCoreMethod1C30(0x0),
 104 fh2JetCoreMethod2C30(0x0),
 105 fh2JetCoreMethod1C60(0x0),
 106 fh2JetCoreMethod2C60(0x0),
 107 fh2AngStructpt1C10(0x0),
 108 fh2AngStructpt2C10(0x0),
 109 fh2AngStructpt3C10(0x0),
 110 fh2AngStructpt4C10(0x0),
 111 fh2AngStructpt1C20(0x0),
 112 fh2AngStructpt2C20(0x0),
 113 fh2AngStructpt3C20(0x0),
 114 fh2AngStructpt4C20(0x0),
 115 fh2AngStructpt1C30(0x0),
 116 fh2AngStructpt2C30(0x0),
 117 fh2AngStructpt3C30(0x0),
 118 fh2AngStructpt4C30(0x0),
 119 fh2AngStructpt1C60(0x0),
 120 fh2AngStructpt2C60(0x0),
 121 fh2AngStructpt3C60(0x0),
 122 fh2AngStructpt4C60(0x0),
 123 fh3spectriggered(0x0),
 124 fh3specbiased(0x0),
 125 fh3spectot(0x0)
 126
 127
 128 {
 129    // default Constructor
 130
 131
 132  // Trigger buffer.
 133    for(Int_t i=0; i<10; i++) {
 134                 for(Int_t j=0; j<7; j++) {
 135                         fTrigBuffer[i][j]=0;
 136                 }
 137         }
 138
 139
 140
 141
 142
 143    fJetBranchName[0] = "";
 144    fJetBranchName[1] = "";
 145
 146    fListJets[0] = new TList;
 147    fListJets[1] = new TList;
 148 }
 149
 150 AliAnalysisTaskJetCore::AliAnalysisTaskJetCore(const char *name) :
 151 AliAnalysisTaskSE(name),
 152 fESD(0x0),
 153 fAOD(0x0),
 154 fAODExtension(0x0),
 155 fBackgroundBranch(""),
 156 fNonStdFile(""),
 157 fIsPbPb(kTRUE),
 158 fOfflineTrgMask(AliVEvent::kAny),
 159 fMinContribVtx(1),
 160 fVtxZMin(-10.),
 161 fVtxZMax(10.),
 162 fEvtClassMin(0),
 163 fEvtClassMax(4),
 164 fFilterMask(0),
 165 fRadioFrac(0.2),
 166 fMinDist(0.1),
 167 fCentMin(0.),
 168 fCentMax(100.),
 169 fNInputTracksMin(0),
 170 fNInputTracksMax(-1),
 171 fAngStructCloseTracks(0),
 172 fCheckMethods(0),
 173 fDoEventMixing(0),
 174 fJetEtaMin(-.5),
 175 fJetEtaMax(.5),
 176 fNevents(0x0),
 177 fTindex(0x0),
 178 fTrigBufferIndex(0x0),
 179 fJetPtMin(20.),
 180 fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
 181 fJetPtFractionMin(0.5),
 182 fNMatchJets(4),
 183 fMatchMaxDist(0.8),
 184 fKeepJets(kFALSE),
 185 fkNbranches(2),
 186 fkEvtClasses(12),
 187 fOutputList(0x0),
 188 fbEvent(kTRUE),
 189 fHistEvtSelection(0x0),
 190 fhnDeltaR(0x0),
 191 fhnMixedEvents(0x0),
 192 fh2JetCoreMethod1C10(0x0),
 193 fh2JetCoreMethod2C10(0x0),
 194 fh2JetCoreMethod1C20(0x0),
 195 fh2JetCoreMethod2C20(0x0),
 196 fh2JetCoreMethod1C30(0x0),
 197 fh2JetCoreMethod2C30(0x0),
 198 fh2JetCoreMethod1C60(0x0),
 199 fh2JetCoreMethod2C60(0x0),
 200 fh2AngStructpt1C10(0x0),
 201 fh2AngStructpt2C10(0x0),
 202 fh2AngStructpt3C10(0x0),
 203 fh2AngStructpt4C10(0x0),
 204 fh2AngStructpt1C20(0x0),
 205 fh2AngStructpt2C20(0x0),
 206 fh2AngStructpt3C20(0x0),
 207 fh2AngStructpt4C20(0x0),
 208 fh2AngStructpt1C30(0x0),
 209 fh2AngStructpt2C30(0x0),
 210 fh2AngStructpt3C30(0x0),
 211 fh2AngStructpt4C30(0x0),
 212 fh2AngStructpt1C60(0x0),
 213 fh2AngStructpt2C60(0x0),
 214 fh2AngStructpt3C60(0x0),
 215 fh2AngStructpt4C60(0x0),
 216 fh3spectriggered(0x0),
 217 fh3specbiased(0x0),
 218 fh3spectot(0x0)
 219
 220  {
 221    // Constructor
 222
 223
 224     for(Int_t i=0; i<10; i++) {
 225        for(Int_t j=0; j<7; j++) {
 226             fTrigBuffer[i][j]=0;
 227                 }
 228     }
 229
 230
 231
 232    fJetBranchName[0] = "";
 233    fJetBranchName[1] = "";
 234
 235    fListJets[0] = new TList;
 236    fListJets[1] = new TList;
 237
 238    DefineOutput(1, TList::Class());
 239 }
 240
 241 AliAnalysisTaskJetCore::~AliAnalysisTaskJetCore()
 242 {
 243    delete fListJets[0];
 244    delete fListJets[1];
 245 }
 246
 247 void AliAnalysisTaskJetCore::SetBranchNames(const TString &branch1, const TString &branch2)
 248 {
 249    fJetBranchName[0] = branch1;
 250    fJetBranchName[1] = branch2;
 251 }
 252
 253 void AliAnalysisTaskJetCore::Init()
 254 {
 255
 256    // check for jet branches
 257    if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
 258       AliError("Jet branch name not set.");
 259    }
 260
 261 }
 262
 263 void AliAnalysisTaskJetCore::UserCreateOutputObjects()
 264 {
 265    // Create histograms
 266    // Called once
 267    OpenFile(1);
 268    if(!fOutputList) fOutputList = new TList;
 269    fOutputList->SetOwner(kTRUE);
 270
 271    Bool_t oldStatus = TH1::AddDirectoryStatus();
 272    TH1::AddDirectory(kFALSE);
 273
 274
 275    fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 6, -0.5, 5.5);
 276    fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
 277    fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN");
 278    fHistEvtSelection->GetXaxis()->SetBinLabel(3,"event selection (rejected)");
 279    fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)");
 280    fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)");
 281    fHistEvtSelection->GetXaxis()->SetBinLabel(6,"multiplicity (rejected)");
 282
 283      UInt_t entries = 0; // bit coded, see GetDimParams() below
 284      entries = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6 |1<<7;
 285      fhnDeltaR = NewTHnSparseF("fhnDeltaR", entries);
 286
 287
 288      UInt_t cifras = 0; // bit coded, see GetDimParams() below
 289      cifras = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6 |1<<7;
 290      fhnMixedEvents = NewTHnSparseF("fhnMixedEvents", cifras);
 291
 292     if(fCheckMethods){
 293
 294
 295     fh2JetCoreMethod1C10 = new TH2F("JetCoreMethod1C10","",150, 0., 150.,100, 0., 1.5);
 296     fh2JetCoreMethod2C10 = new TH2F("JetCoreMethod2C10","",150, 0., 150.,100, 0., 1.5);
 297     fh2JetCoreMethod1C20 = new TH2F("JetCoreMethod1C20","",150, 0., 150.,100, 0., 1.5);
 298     fh2JetCoreMethod2C20 = new TH2F("JetCoreMethod2C20","",150, 0., 150.,100, 0., 1.5);
 299     fh2JetCoreMethod1C30 = new TH2F("JetCoreMethod1C30","",150, 0., 150.,100, 0., 1.5);
 300     fh2JetCoreMethod2C30 = new TH2F("JetCoreMethod2C30","",150, 0., 150.,100, 0., 1.5);
 301     fh2JetCoreMethod1C60 = new TH2F("JetCoreMethod1C60","",150, 0., 150.,100, 0., 1.5);
 302     fh2JetCoreMethod2C60 = new TH2F("JetCoreMethod2C60","",150, 0., 150.,100, 0., 1.5);}
 303
 304
 305     if(fAngStructCloseTracks>0){
 306     fh2AngStructpt1C10 = new TH2F("Ang struct pt1 C10","",15,0.,1.5,150,0.,10.);
 307     fh2AngStructpt2C10 = new TH2F("Ang struct pt2 C10","",15,0.,1.5,150,0.,10.);
 308     fh2AngStructpt3C10 = new TH2F("Ang struct pt3 C10","",15,0.,1.5,150,0.,10.);
 309     fh2AngStructpt4C10 = new TH2F("Ang struct pt4 C10","",15,0.,1.5,150,0.,10.);
 310     fh2AngStructpt1C20 = new TH2F("Ang struct pt1 C20","",15,0.,1.5,150,0.,10.);
 311     fh2AngStructpt2C20 = new TH2F("Ang struct pt2 C20","",15,0.,1.5,150,0.,10.);
 312     fh2AngStructpt3C20 = new TH2F("Ang struct pt3 C20","",15,0.,1.5,150,0.,10.);
 313     fh2AngStructpt4C20 = new TH2F("Ang struct pt4 C20","",15,0.,1.5,150,0.,10.);
 314     fh2AngStructpt1C30 = new TH2F("Ang struct pt1 C30","",15,0.,1.5,150,0.,10.);
 315     fh2AngStructpt2C30 = new TH2F("Ang struct pt2 C30","",15,0.,1.5,150,0.,10.);
 316     fh2AngStructpt3C30 = new TH2F("Ang struct pt3 C30","",15,0.,1.5,150,0.,10.);
 317     fh2AngStructpt4C30 = new TH2F("Ang struct pt4 C30","",15,0.,1.5,150,0.,10.);
 318     fh2AngStructpt1C60 = new TH2F("Ang struct pt1 C60","",15,0.,1.5,150,0.,10.);
 319     fh2AngStructpt2C60 = new TH2F("Ang struct pt2 C60","",15,0.,1.5,150,0.,10.);
 320     fh2AngStructpt3C60 = new TH2F("Ang struct pt3 C60","",15,0.,1.5,150,0.,10.);
 321     fh2AngStructpt4C60 = new TH2F("Ang struct pt4 C60","",15,0.,1.5,150,0.,10.); }
 322     fh3spectriggered = new TH3F("Triggered spectrum","",10,0,100,50,0.,200,50,0.,50.);
 323     fh3specbiased = new TH3F("Biased spectrum","",10,0,100,50,0.,200.,50,0.,50.);
 324     fh3spectot = new TH3F("Total spectrum","",50,0.,200.,50,0.,50.,50,0.,50.);
 325
 326    fOutputList->Add(fHistEvtSelection);
 327
 328    fOutputList->Add(fhnDeltaR);
 329
 330    fOutputList->Add(fhnMixedEvents);
 331
 332
 333
 334       if(fCheckMethods){
 335
 336       fOutputList->Add(fh2JetCoreMethod1C10);
 337       fOutputList->Add(fh2JetCoreMethod2C10);
 338       fOutputList->Add(fh2JetCoreMethod1C20);
 339       fOutputList->Add(fh2JetCoreMethod2C20);
 340       fOutputList->Add(fh2JetCoreMethod1C30);
 341       fOutputList->Add(fh2JetCoreMethod2C30);
 342       fOutputList->Add(fh2JetCoreMethod1C60);
 343       fOutputList->Add(fh2JetCoreMethod2C60);}
 344
 345
 346
 347
 348
 349         if(fAngStructCloseTracks>0){
 350        fOutputList->Add(fh2AngStructpt1C10);
 351        fOutputList->Add(fh2AngStructpt2C10);
 352        fOutputList->Add(fh2AngStructpt3C10);
 353        fOutputList->Add(fh2AngStructpt4C10);
 354        fOutputList->Add(fh2AngStructpt1C20);
 355        fOutputList->Add(fh2AngStructpt2C20);
 356        fOutputList->Add(fh2AngStructpt3C20);
 357        fOutputList->Add(fh2AngStructpt4C20);
 358        fOutputList->Add(fh2AngStructpt1C30);
 359        fOutputList->Add(fh2AngStructpt2C30);
 360        fOutputList->Add(fh2AngStructpt3C30);
 361        fOutputList->Add(fh2AngStructpt4C30);
 362        fOutputList->Add(fh2AngStructpt1C60);
 363        fOutputList->Add(fh2AngStructpt2C60);
 364        fOutputList->Add(fh2AngStructpt3C60);
 365        fOutputList->Add(fh2AngStructpt4C60);}
 366        fOutputList->Add(fh3spectriggered);
 367        fOutputList->Add(fh3specbiased);
 368        fOutputList->Add(fh3spectot);
 369
 370    // =========== Switch on Sumw2 for all histos ===========
 371    for (Int_t i=0; i<fOutputList->GetEntries(); ++i) {
 372       TH1 *h1 = dynamic_cast<TH1*>(fOutputList->At(i));
 373       if (h1){
 374          h1->Sumw2();
 375          continue;
 376       }
 377     THnSparse *hn = dynamic_cast<THnSparse*>(fOutputList->At(i));
 378       if (hn){
 379          hn->Sumw2();
 380       }
 381    }
 382    TH1::AddDirectory(oldStatus);
 383
 384    PostData(1, fOutputList);
 385 }
 386
 387 void AliAnalysisTaskJetCore::UserExec(Option_t *)
 388 {
 389
 390
 391    if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
 392       AliError("Jet branch name not set.");
 393       return;
 394    }
 395
 396    fESD=dynamic_cast<AliESDEvent*>(InputEvent());
 397    if (!fESD) {
 398       AliError("ESD not available");
 399       fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
 400    } else {
 401       fAOD = dynamic_cast<AliAODEvent*>(AODEvent());
 402    }
 403
 404     if(fNonStdFile.Length()!=0){
 405     // case that we have an AOD extension we need can fetch the jets from the extended output
 406     AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
 407     fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
 408     if(!fAODExtension){
 409       if(fDebug>1)Printf("AODExtension found for %s",fNonStdFile.Data());
 410     }}
 411
 412
 413
 414
 415
 416    // -- event selection --
 417    fHistEvtSelection->Fill(1); // number of events before event selection
 418
 419    // physics selection
 420    AliInputEventHandler* inputHandler = (AliInputEventHandler*)
 421    ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
 422    cout<<inputHandler->IsEventSelected()<<" "<<fOfflineTrgMask<<endl;
 423    if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){
 424       if(fDebug) Printf(" Trigger Selection: event REJECTED ... ");
 425       fHistEvtSelection->Fill(2);
 426       PostData(1, fOutputList);
 427       return;
 428    }
 429
 430    // vertex selection
 431    if(!fAOD){
 432      if(fDebug) Printf("%s:%d No AOD",(char*)__FILE__,__LINE__);
 433      fHistEvtSelection->Fill(3);
 434       PostData(1, fOutputList);
 435    }
 436    AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
 437
 438    if(!primVtx){
 439      if(fDebug) Printf("%s:%d No primVtx",(char*)__FILE__,__LINE__);
 440      fHistEvtSelection->Fill(3);
 441       PostData(1, fOutputList);
 442    }
 443
 444    Int_t nTracksPrim = primVtx->GetNContributors();
 445    if ((nTracksPrim < fMinContribVtx) ||
 446          (primVtx->GetZ() < fVtxZMin) ||
 447          (primVtx->GetZ() > fVtxZMax) ){
 448       if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
 449       fHistEvtSelection->Fill(3);
 450       PostData(1, fOutputList);
 451       return;
 452    }
 453
 454    // event class selection (from jet helper task)
 455    Int_t eventClass = AliAnalysisHelperJetTasks::EventClass();
 456    if(fDebug) Printf("Event class %d", eventClass);
 457    if (eventClass < fEvtClassMin || eventClass > fEvtClassMax){
 458       fHistEvtSelection->Fill(4);
 459       PostData(1, fOutputList);
 460       return;
 461    }
 462
 463    // centrality selection
 464    AliCentrality *cent = 0x0;
 465    Double_t centValue = 0.;
 466    if(fESD) cent = fESD->GetCentrality();
 467    if(cent) centValue = cent->GetCentralityPercentile("V0M");
 468    if(fDebug) printf("centrality: %f\n", centValue);
 469    if (centValue < fCentMin || centValue > fCentMax){
 470       fHistEvtSelection->Fill(4);
 471       PostData(1, fOutputList);
 472       return;
 473    }
 474
 475
 476    fHistEvtSelection->Fill(0);
 477    // accepted events
 478    // -- end event selection --
 479
 480    // get background
 481    AliAODJetEventBackground* externalBackground = 0;
 482    if(fAOD&&!externalBackground&&fBackgroundBranch.Length()){
 483       externalBackground =  (AliAODJetEventBackground*)(fAOD->FindListObject(fBackgroundBranch.Data()));
 484       if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
 485    }
 486    if(fAODExtension&&!externalBackground&&fBackgroundBranch.Length()){
 487      externalBackground =  (AliAODJetEventBackground*)(fAODExtension->GetAOD()->FindListObject(fBackgroundBranch.Data()));
 488       if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
 489    }
 490
 491    Float_t rho = 0;
 492    if(externalBackground)rho = externalBackground->GetBackground(0);
 493
 494
 495    // fetch jets
 496    TClonesArray *aodJets[2];
 497    aodJets[0]=0;
 498    if(fAOD&&!aodJets[0]){
 499    aodJets[0] = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fJetBranchName[0].Data()));
 500    aodJets[1] = dynamic_cast<TClonesArray*>(fAOD->FindListObject(fJetBranchName[1].Data()));  }
 501    if(fAODExtension && !aodJets[0]){
 502    aodJets[0] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[0].Data()));
 503    aodJets[1] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[1].Data()));  }
 504
 505    //Double_t ptsub[aodJets[0]->GetEntriesFast()];
 506    //Int_t inord[aodJets[0]->GetEntriesFast()];
 507    //for(Int_t n=0;n<aodJets[0]->GetEntriesFast();n++){
 508    //  ptsub[n]=0;
 509    //  inord[n]=0;}
 510
 511    TList ParticleList;
 512    Int_t nT = GetListOfTracks(&ParticleList);
 513      for (Int_t iJetType = 0; iJetType < 2; iJetType++) {
 514       fListJets[iJetType]->Clear();
 515       if (!aodJets[iJetType]) continue;
 516
 517       if(fDebug) Printf("%s: %d jets",fJetBranchName[iJetType].Data(),aodJets[iJetType]->GetEntriesFast());
 518
 519
 520       for (Int_t iJet = 0; iJet < aodJets[iJetType]->GetEntriesFast(); iJet++) {
 521          AliAODJet *jet = dynamic_cast<AliAODJet*>((*aodJets[iJetType])[iJet]);
 522          if (jet) fListJets[iJetType]->Add(jet);
 523          // if(iJetType==0){
 524          // ptsub[iJet]=jet->Pt()-rho*jet->EffectiveAreaCharged();}
 525       }}
 526
 527    Double_t etabig=0;
 528    Double_t ptbig=0;
 529    Double_t areabig=0;
 530    Double_t phibig=0.;
 531    Double_t etasmall=0;
 532    Double_t ptsmall=0;
 533    Double_t areasmall=0;
 534    //Double_t distr=0.;
 535    Double_t phismall=0.;
 536
 537
 538    // Double_t up1[15]={0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
 539    // Double_t up2[15]={0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
 540    // Double_t up3[15]={0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
 541    // Double_t up4[15]={0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
 542    // Double_t down1[15]={0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
 543    // Double_t down2[15]={0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
 544    // Double_t down3[15]={0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
 545    // Double_t down4[15]={0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
 546    Int_t iCount=0;
 547    Int_t trigJet=-1;
 548    Int_t trigBBTrack=-1;
 549    Int_t trigInTrack=-1;
 550
 551    for(Int_t i=0; i<fListJets[0]->GetEntries(); ++i){
 552            AliAODJet* jetbig = (AliAODJet*)(fListJets[0]->At(i));
 553            etabig  = jetbig->Eta();
 554            phibig  = jetbig->Phi();
 555            ptbig   = jetbig->Pt();
 556            if(ptbig==0) continue;
 557            areabig = jetbig->EffectiveAreaCharged();
 558            Double_t ptcorr=ptbig-rho*areabig;
 559            if(ptcorr<=0) continue;
 560            if((etabig<fJetEtaMin)||(etabig>fJetEtaMax)) continue;
 561                    Double_t dismin=100.;
 562                    Double_t ptmax=-10.;
 563                    Int_t index1=-1;
 564                    Int_t index2=-1;
 565
 566                Int_t point=GetHardestTrackBackToJet(jetbig);
 567                AliVParticle *partback = (AliVParticle*)ParticleList.At(point);
 568                if(!partback) continue;
 569                    fh3spectriggered->Fill(centValue,ptcorr,partback->Pt());
 570                        AliAODTrack* leadtrack;
 571                        Int_t ippt=0;
 572                        Double_t ppt=-10;
 573                        TRefArray *genTrackList = jetbig->GetRefTracks();
 574                        Int_t nTracksGenJet = genTrackList->GetEntriesFast();
 575                        AliAODTrack* genTrack;
 576                        for(Int_t ir=0; ir<nTracksGenJet; ++ir){
 577                        genTrack = (AliAODTrack*)(genTrackList->At(ir));
 578                        if(genTrack->Pt()>ppt){ppt=genTrack->Pt();
 579                        ippt=ir;}}
 580                         leadtrack=(AliAODTrack*)(genTrackList->At(ippt));
 581                         if(!leadtrack) continue;
 582                         fh3specbiased->Fill(centValue,ptcorr,leadtrack->Pt());
 583                         if(centValue<10)fh3spectot->Fill(ptcorr,leadtrack->Pt(),partback->Pt());
 584                         //store one trigger info
 585                         if((partback->Pt()>10.)&&(iCount==0)){
 586                           trigJet=i;
 587                           trigBBTrack=point;
 588                           trigInTrack=ippt;
 589                           iCount=iCount+1;}
 590
 591
 592                   if(fCheckMethods){
 593                   for(Int_t j=0; j<fListJets[1]->GetEntries(); ++j){
 594                   AliAODJet* jetsmall = (AliAODJet*)(fListJets[1]->At(j));
 595                   etasmall  = jetsmall->Eta();
 596                   phismall = jetsmall->Phi();
 597                   ptsmall   = jetsmall->Pt();
 598                   areasmall = jetsmall->EffectiveAreaCharged();
 599                   Double_t tmpDeltaR=(phismall-phibig)*(phismall-phibig)+(etasmall-etabig)*(etasmall-etabig);
 600                   tmpDeltaR=TMath::Sqrt(tmpDeltaR);
 601                      //Fraction in the jet core
 602                     if((ptsmall>ptmax)&&(tmpDeltaR<=fRadioFrac)){ptmax=ptsmall;
 603                     index2=j;}
 604                     if(tmpDeltaR<=dismin){ dismin=tmpDeltaR;
 605                       index1=j;}} //en of loop over R=0.2 jets
 606                 //method1:most concentric jet=core
 607                 if(dismin<fMinDist){ AliAODJet* jetmethod1 = (AliAODJet*)(fListJets[1]->At(index1));
 608                   if(centValue<10) fh2JetCoreMethod1C10->Fill(ptcorr,jetmethod1->Pt()/ptbig);
 609                   if((centValue>20)&&(centValue<40)) fh2JetCoreMethod1C20->Fill(ptcorr,jetmethod1->Pt()/ptbig);
 610                   if((centValue>30)&&(centValue<60)) fh2JetCoreMethod1C30->Fill(ptcorr,jetmethod1->Pt()/ptbig);
 611                   if(centValue>60) fh2JetCoreMethod1C60->Fill(ptcorr,jetmethod1->Pt()/ptbig); }
 612                 //method2:hardest contained jet=core
 613                 if(index2!=-1){
 614                   AliAODJet* jetmethod2 = (AliAODJet*)(fListJets[1]->At(index2));
 615                   if(centValue<10) fh2JetCoreMethod2C10->Fill(ptcorr,jetmethod2->Pt()/ptbig);
 616                   if((centValue>20)&&(centValue<40)) fh2JetCoreMethod2C20->Fill(ptcorr,jetmethod2->Pt()/ptbig);
 617                   if((centValue>30)&&(centValue<60)) fh2JetCoreMethod2C30->Fill(ptcorr,jetmethod2->Pt()/ptbig);
 618                   if(centValue>60) fh2JetCoreMethod2C60->Fill(ptcorr,jetmethod2->Pt()/ptbig); }}
 619
 620           for(int it = 0;it<nT;++it){
 621           AliVParticle *part = (AliVParticle*)ParticleList.At(it);
 622           Double_t deltaR = jetbig->DeltaR(part);
 623           Double_t deltaEta = etabig-part->Eta();
 624           Double_t deltaPhi=phibig-part->Phi();
 625           if(deltaPhi<-0.5*TMath::Pi()) deltaPhi+=2.*TMath::Pi();
 626           if(deltaPhi>3./2.*TMath::Pi()) deltaPhi-=2.*TMath::Pi();
 627           Double_t jetEntries[8] = {centValue,ptcorr,part->Pt(),deltaR,deltaEta,deltaPhi,leadtrack->Pt(),partback->Pt()};                     fhnDeltaR->Fill(jetEntries);
 628           }
 629           //end of track loop
 630           }
 631           //end of jet loop
 632
 633
 634
 635
 636           if(fDoEventMixing){
 637             //check before if the trigger exists
 638             // fTrigBuffer[i][0] = zvtx
 639             // fTrigBuffer[i][1] = phi
 640             // fTrigBuffer[i][2] = eta
 641             // fTrigBuffer[i][3] = pt_jet
 642             // fTrigBuffer[i][4] = pt_trig
 643             // fTrigBuffer[i][5]= pt_track_in
 644             // fTrigBuffer[i][6]= centrality
 645             if(fTindex==11) fTindex=0;
 646             if(fTrigBuffer[fTindex][3]>0){
 647             if (TMath::Abs(fTrigBuffer[fTindex][0]-primVtx->GetZ()<2.)){
 648             if (TMath::Abs(fTrigBuffer[fTindex][6]-centValue<10)){
 649
 650                         for(int it = 0;it<nT;++it){
 651                         AliVParticle *part = (AliVParticle*)ParticleList.At(it);                        Double_t DPhi = fTrigBuffer[fTindex][1] - part->Phi();
 652                         Double_t DEta = fTrigBuffer[fTindex][2] - part->Eta();
 653                         Double_t DR=TMath::Sqrt(DPhi*DPhi+DEta*DEta);
 654                         if(DPhi<-0.5*TMath::Pi()) DPhi+=2.*TMath::Pi();
 655                         if(DPhi>3./2.*TMath::Pi()) DPhi-=2.*TMath::Pi();
 656                         Double_t triggerEntries[8] = {centValue,fTrigBuffer[fTindex][3],part->Pt(),DR,DEta,DPhi,fTrigBuffer[fTindex][4],fTrigBuffer[fTindex][5]};
 657                         fhnMixedEvents->Fill(triggerEntries);
 658                         }
 659                         fNevents=fNevents+1;
 660                         if(fNevents==9) {fTindex=fTindex+1;
 661                         fNevents=0;}
 662             }}}
 663
 664
 665                // Copy the triggers from the current event into the buffer.
 666                //again, only if the trigger exists:
 667                 if(trigJet>-1){
 668                 AliAODJet* jetT = (AliAODJet*)(fListJets[0]->At(trigJet));
 669                 AliVParticle *partL = (AliVParticle*)ParticleList.At(trigInTrack);
 670                 AliVParticle *partT = (AliVParticle*)ParticleList.At(trigBBTrack);
 671                 fTrigBuffer[fTrigBufferIndex][0] = primVtx->GetZ();
 672                 fTrigBuffer[fTrigBufferIndex][1] = jetT->Phi();
 673                 fTrigBuffer[fTrigBufferIndex][2] = jetT->Eta();
 674                 fTrigBuffer[fTrigBufferIndex][3] = jetT->Pt()-rho*jetT->EffectiveAreaCharged();
 675                 fTrigBuffer[fTrigBufferIndex][4] = partT->Pt();
 676                 fTrigBuffer[fTrigBufferIndex][5] = partL->Pt();
 677                 fTrigBuffer[fTrigBufferIndex][6] = centValue;
 678                 fTrigBufferIndex++;
 679                 if(fTrigBufferIndex==9) fTrigBufferIndex=0;
 680                 }
 681           }
 682
 683
 684
 685
 686
 687      //////////////////ANGULAR STRUCTURE//////////////////////////////////////
 688
 689      //tracks up to R=0.8 distant from the jet axis
 690    //   if(fAngStructCloseTracks==1){
 691    //    TList CloseTrackList;
 692    //    Int_t nn=GetListOfTracksCloseToJet(&CloseTrackList,jetbig);
 693    //    Double_t difR=0.04;
 694    //    for(Int_t l=0;l<15;l++){
 695    //    Double_t rr=l*0.1+0.1;
 696    //     for(int it = 0;it<nn;++it){
 697    //         AliVParticle *part1 = (AliVParticle*)CloseTrackList.At(it);
 698    //         for(int itu=it+1;itu<CloseTrackList.GetEntries();itu++){
 699    //         AliVParticle *part2 = (AliVParticle*)CloseTrackList.At(itu);
 700    //         Double_t ptm=part1->Pt();
 701    //         Double_t ptn=part2->Pt();
 702    //         Double_t Rnm = (part1->Eta()-part2->Eta())*(part1->Eta()-part2->Eta())+(part1->Phi()-part2->Phi())*(part1->Phi()-part2->Phi());
 703    //                    Rnm=TMath::Sqrt(Rnm);
 704    //                    Double_t deltag=(1./(TMath::Sqrt(2*TMath::Pi())*difR))*TMath::Exp(-1.*(rr-Rnm)*(rr-Rnm)/(2.*difR*difR));
 705    //                    Double_t stepf=0.5*(1.+TMath::Erf((rr-Rnm)/(TMath::Sqrt(2.)*difR)));
 706    //                    if((ptcorr<85.) && (ptcorr>=70.)){up1[l]=up1[l]+ptm*ptn*Rnm*Rnm*deltag;
 707    //                                                   down1[l]=down1[l]+ptm*ptn*Rnm*Rnm*stepf;}
 708    //                    if((ptcorr<100.) && (ptcorr>=85.)){up2[l]=up2[l]+ptm*ptn*Rnm*Rnm*deltag;
 709    //                                                   down2[l]=down2[l]+ptm*ptn*Rnm*Rnm*stepf;}
 710    //                    if((ptcorr<120.) && (ptcorr>=100.)){up3[l]=up3[l]+ptm*ptn*Rnm*Rnm*deltag;
 711    //                                                    down3[l]=down3[l]+ptm*ptn*Rnm*Rnm*stepf;}
 712    //                    if((ptcorr<140.) && (ptcorr>=120.)){up4[l]=up4[l]+ptm*ptn*Rnm*Rnm*deltag;
 713    //                   down4[l]=down4[l]+ptm*ptn*Rnm*Rnm*stepf;}}}}
 714    //   }
 715
 716    //   //only jet constituents
 717    //    if(fAngStructCloseTracks==2){
 718
 719    //    Double_t difR=0.04;
 720    //    for(Int_t l=0;l<15;l++){
 721    //    Double_t rr=l*0.1+0.1;
 722
 723
 724    //    AliAODTrack* part1;
 725    //    AliAODTrack* part2;
 726
 727    //        TRefArray *genTrackListb = jetbig->GetRefTracks();
 728    //        Int_t nTracksGenJetb = genTrackListb->GetEntriesFast();
 729
 730
 731
 732    //        for(Int_t it=0; it<nTracksGenJetb; ++it){
 733    //           part1 = (AliAODTrack*)(genTrackListb->At(it));
 734    //         for(Int_t itu=0; itu<nTracksGenJetb; ++itu){
 735    //           part2 = (AliAODTrack*)(genTrackListb->At(itu));
 736    //         Double_t ptm=part1->Pt();
 737    //         Double_t ptn=part2->Pt();
 738    //         Double_t Rnm = (part1->Eta()-part2->Eta())*(part1->Eta()-part2->Eta())+(part1->Phi()-part2->Phi())*(part1->Phi()-part2->Phi());
 739    //                    Rnm=TMath::Sqrt(Rnm);
 740    //                    Double_t deltag=(1./(TMath::Sqrt(2*TMath::Pi())*difR))*TMath::Exp(-1.*(rr-Rnm)*(rr-Rnm)/(2.*difR*difR));
 741    //                    Double_t stepf=0.5*(1.+TMath::Erf((rr-Rnm)/(TMath::Sqrt(2.)*difR)));
 742    //                    if((ptcorr<85.) && (ptcorr>=70.)){up1[l]=up1[l]+ptm*ptn*Rnm*Rnm*deltag;
 743    //                                                   down1[l]=down1[l]+ptm*ptn*Rnm*Rnm*stepf;}
 744    //                    if((ptcorr<100.) && (ptcorr>=85.)){up2[l]=up2[l]+ptm*ptn*Rnm*Rnm*deltag;
 745    //                                                   down2[l]=down2[l]+ptm*ptn*Rnm*Rnm*stepf;}
 746    //                    if((ptcorr<120.) && (ptcorr>=100.)){up3[l]=up3[l]+ptm*ptn*Rnm*Rnm*deltag;
 747    //                                                    down3[l]=down3[l]+ptm*ptn*Rnm*Rnm*stepf;}
 748    //                    if((ptcorr<140.) && (ptcorr>=120.)){up4[l]=up4[l]+ptm*ptn*Rnm*Rnm*deltag;
 749    //                   down4[l]=down4[l]+ptm*ptn*Rnm*Rnm*stepf;}}}}}
 750    // }
 751      // //end loop over R=0.4 jets
 752      // if(fAngStructCloseTracks>0){
 753      // for(Int_t l=0;l<15;l++){
 754      // Double_t rr=l*0.1+0.1;
 755      //    if(down1[l]!=0){
 756      //         if(centValue<10.)fh2AngStructpt1C10->Fill(rr,rr*up1[l]/down1[l]);
 757      //    if(centValue>20. && centValue<40.) fh2AngStructpt1C20->Fill(rr,rr*up1[l]/down1[l]);
 758      //    if(centValue>30. && centValue<60.) fh2AngStructpt1C30->Fill(rr,rr*up1[l]/down1[l]);
 759      //    if(centValue>60.) fh2AngStructpt1C60->Fill(rr,rr*up1[l]/down1[l]);}
 760      //    if(down2[l]!=0){
 761      //         if(centValue<10.) fh2AngStructpt2C10->Fill(rr,rr*up2[l]/down2[l]);
 762      //    if(centValue>20. && centValue<40.) fh2AngStructpt2C20->Fill(rr,rr*up2[l]/down2[l]);
 763      //    if(centValue>30. && centValue<60.) fh2AngStructpt2C30->Fill(rr,rr*up2[l]/down2[l]);
 764      //    if(centValue>60.) fh2AngStructpt2C60->Fill(rr,rr*up2[l]/down2[l]);}
 765      //    if(down3[l]!=0){
 766      //         if(centValue<10.) fh2AngStructpt3C10->Fill(rr,rr*up3[l]/down3[l]);
 767      //    if(centValue>20. && centValue<40.) fh2AngStructpt3C20->Fill(rr,rr*up3[l]/down3[l]);
 768      //    if(centValue>30. && centValue<60.) fh2AngStructpt3C30->Fill(rr,rr*up3[l]/down3[l]);
 769      //    if(centValue>60.) fh2AngStructpt3C60->Fill(rr,rr*up3[l]/down3[l]);}
 770      //    if(down4[l]!=0){
 771      //         if(centValue<10.) fh2AngStructpt4C10->Fill(rr,rr*up4[l]/down4[l]);
 772      //    if(centValue>20. && centValue<40.) fh2AngStructpt4C20->Fill(rr,rr*up4[l]/down4[l]);
 773      //    if(centValue>30. && centValue<60.) fh2AngStructpt4C30->Fill(rr,rr*up4[l]/down4[l]);
 774      //    if(centValue>60.) fh2AngStructpt4C60->Fill(rr,rr*up4[l]/down4[l]);}}}
 775
 776
 777
 778
 779
 780
 781
 782    PostData(1, fOutputList);
 783 }
 784
 785 void AliAnalysisTaskJetCore::Terminate(const Option_t *)
 786 {
 787    // Draw result to the screen
 788    // Called once at the end of the query
 789
 790    if (!GetOutputData(1))
 791    return;
 792 }
 793
 794
 795
 796
 797
 798
 799
 800
 801
 802
 803
 804 Int_t  AliAnalysisTaskJetCore::GetListOfTracks(TList *list){
 805
 806     Int_t iCount = 0;
 807
 808
 809     for(int it = 0;it < fAOD->GetNumberOfTracks();++it){
 810       AliAODTrack *tr = fAOD->GetTrack(it);
 811       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
 812       if(TMath::Abs(tr->Eta())>0.9)continue;
 813       if(tr->Pt()<0.15)continue;
 814       list->Add(tr);
 815       //cout<<fAOD->GetNumberOfTracks()<<" "<<tr->Pt()<<endl;
 816       iCount++;
 817     }
 818
 819
 820   return iCount;
 821
 822 }
 823
 824    Int_t  AliAnalysisTaskJetCore::GetHardestTrackBackToJet(AliAODJet *jetbig){
 825
 826
 827     Int_t index=-1;
 828     Double_t ptmax=-10;
 829     Double_t dphi=0;
 830     Double_t dif=0;
 831     Int_t iCount=0;
 832     for(int it = 0;it < fAOD->GetNumberOfTracks();++it){
 833       AliAODTrack *tr = fAOD->GetTrack(it);
 834       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
 835       if(TMath::Abs(tr->Eta())>0.9)continue;
 836       if(tr->Pt()<0.15)continue;
 837       iCount=iCount+1;
 838       dphi=RelativePhi(tr->Phi(),jetbig->Phi());
 839       if(TMath::Abs(dphi)<TMath::Pi()-0.2) continue;
 840       if(tr->Pt()>ptmax){ ptmax=tr->Pt();
 841       index=iCount-1;
 842       dif=dphi;  }}
 843
 844       return index;
 845
 846    }
 847
 848
 849
 850
 851
 852
 853
 854
 855
 856  Int_t  AliAnalysisTaskJetCore::GetListOfTracksCloseToJet(TList *list,AliAODJet *jetbig){
 857
 858     Int_t iCount = 0;
 859
 860
 861     for(int it = 0;it < fAOD->GetNumberOfTracks();++it){
 862       AliAODTrack *tr = fAOD->GetTrack(it);
 863       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
 864       if(TMath::Abs(tr->Eta())>0.9)continue;
 865       if(tr->Pt()<0.15)continue;
 866       Double_t disR=jetbig->DeltaR(tr);
 867       if(disR>0.8)  continue;
 868       list->Add(tr);
 869       //cout<<fAOD->GetNumberOfTracks()<<" "<<tr->Pt()<<endl;
 870       iCount++;
 871     }
 872
 873    list->Sort();
 874    return iCount;
 875
 876 }
 877
 878
 879
 880
 881
 882
 883
 884
 885
 886
 887
 888 Int_t AliAnalysisTaskJetCore::GetNInputTracks()
 889 {
 890
 891    Int_t nInputTracks = 0;
 892
 893    TString jbname(fJetBranchName[1]);
 894    //needs complete event, use jets without background subtraction
 895    for(Int_t i=1; i<=3; ++i){
 896       if(jbname.Contains(Form("B%d",i))) jbname.ReplaceAll(Form("B%d",i),"B0");
 897    }
 898    // use only HI event
 899    if(jbname.Contains("AODextraonly")) jbname.ReplaceAll("AODextraonly","AOD");
 900    if(jbname.Contains("AODextra")) jbname.ReplaceAll("AODextra","AOD");
 901
 902    if(fDebug) Printf("Multiplicity from jet branch %s", jbname.Data());
 903    TClonesArray *tmpAODjets = dynamic_cast<TClonesArray*>(fAOD->FindListObject(jbname.Data()));
 904    if(!tmpAODjets){
 905       Printf("Jet branch %s not found", jbname.Data());
 906       Printf("AliAnalysisTaskJetCore::GetNInputTracks FAILED");
 907       return -1;
 908    }
 909
 910    for (Int_t iJet=0; iJet<tmpAODjets->GetEntriesFast(); iJet++){
 911       AliAODJet *jet = dynamic_cast<AliAODJet*>((*tmpAODjets)[iJet]);
 912       if(!jet) continue;
 913       TRefArray *trackList = jet->GetRefTracks();
 914       Int_t nTracks = trackList->GetEntriesFast();
 915       nInputTracks += nTracks;
 916       if(fDebug) Printf("#jet%d: %d tracks", iJet, nTracks);
 917    }
 918    if(fDebug) Printf("---> input tracks: %d", nInputTracks);
 919
 920    return nInputTracks;
 921 }
 922
 923
 924
 925 Double_t AliAnalysisTaskJetCore::RelativePhi(Double_t mphi,Double_t vphi){
 926
 927   if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi());
 928   else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
 929   if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi());
 930   else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
 931   double dphi = mphi-vphi;
 932   if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
 933   else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
 934   return dphi;//dphi in [-Pi, Pi]
 935 }
 936
 937
 938
 939 THnSparse* AliAnalysisTaskJetCore::NewTHnSparseF(const char* name, UInt_t entries)
 940 {
 941    // generate new THnSparseF, axes are defined in GetDimParams()
 942
 943    Int_t count = 0;
 944    UInt_t tmp = entries;
 945    while(tmp!=0){
 946       count++;
 947       tmp = tmp &~ -tmp;  // clear lowest bit
 948    }
 949
 950    TString hnTitle(name);
 951    const Int_t dim = count;
 952    Int_t nbins[dim];
 953    Double_t xmin[dim];
 954    Double_t xmax[dim];
 955
 956    Int_t i=0;
 957    Int_t c=0;
 958    while(c<dim && i<32){
 959       if(entries&(1<<i)){
 960
 961          TString label("");
 962          GetDimParams(i, label, nbins[c], xmin[c], xmax[c]);
 963          hnTitle += Form(";%s",label.Data());
 964          c++;
 965       }
 966
 967       i++;
 968    }
 969    hnTitle += ";";
 970
 971    return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
 972 }
 973
 974 void AliAnalysisTaskJetCore::GetDimParams(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
 975 {
 976    // stores label and binning of axis for THnSparse
 977
 978    const Double_t pi = TMath::Pi();
 979
 980    switch(iEntry){
 981
 982    case 0:
 983       label = "V0 centrality (%)";
 984
 985          nbins = 10;
 986          xmin = 0.;
 987          xmax = 100.;
 988          break;
 989
 990
 991    case 1:
 992       label = "corrected jet pt";
 993          nbins = 50;
 994          xmin = 0.;
 995          xmax = 200.;
 996           break;
 997
 998
 999    case 2:
1000       label = "track pT";
1001
1002          nbins = 1000;
1003          xmin = 0.;
1004          xmax = 50;
1005          break;
1006
1007
1008     case 3:
1009       label = "deltaR";
1010       nbins = 15;
1011       xmin = 0.;
1012       xmax = 1.5;
1013       break;
1014
1015
1016
1017    case 4:
1018       label = "deltaEta";
1019       nbins = 30;
1020       xmin = -1.5;
1021       xmax = 1.5;
1022       break;
1023
1024
1025   case 5:
1026       label = "deltaPhi";
1027       nbins = 90;
1028       xmin = -0.5*pi;
1029       xmax = 1.5*pi;
1030       break;
1031
1032
1033
1034     case 6:
1035       label = "leading track";
1036       nbins = 50;
1037       xmin = 0;
1038       xmax = 50;
1039       break;
1040
1041      case 7:
1042
1043       label = "trigger track";
1044       nbins =50;
1045       xmin = 0;
1046       xmax = 50;
1047       break;
1048
1049
1050
1051
1052
1053
1054
1055
1056    }
1057
1058 }
1059