PWGJE/AliAnalysisTaskJetProperties.cxx

   1 // **************************************************************************************
   2 // *                                                                                    *
   3 // * Task for Jet properties and jet shape analysis in PWG4 Jet Task Force Train for pp *
   4 // *                                                                                    *
   5 // **************************************************************************************
   6
   7
   8 /**************************************************************************
   9  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  10  *                                                                        *
  11  * Author: The ALICE Off-line Project.                                    *
  12  * Contributors are mentioned in the code where appropriate.              *
  13  *                                                                        *
  14  * Permission to use, copy, modify and distribute this software and its   *
  15  * documentation strictly for non-commercial purposes is hereby granted   *
  16  * without fee, provided that the above copyright notice appears in all   *
  17  * copies and that both the copyright notice and this permission notice   *
  18  * appear in the supporting documentation. The authors make no claims     *
  19  * about the suitability of this software for any purpose. It is          *
  20  * provided "as is" without express or implied warranty.                  *
  21  **************************************************************************/
  22
  23 /* $Id: */
  24
  25 #include "TList.h"
  26 #include "TH1F.h"
  27 #include "TH2F.h"
  28 #include "TH3F.h"
  29 #include "TString.h"
  30 #include "THnSparse.h"
  31 #include "TProfile.h"
  32 #include "TFile.h"
  33 #include "TKey.h"
  34 #include "TRandom3.h"
  35
  36 #include "AliAODInputHandler.h"
  37 #include "AliAODHandler.h"
  38 #include "AliESDEvent.h"
  39 #include "AliAODMCParticle.h"
  40 #include "AliAODJet.h"
  41 #include "AliGenPythiaEventHeader.h"
  42 #include "AliGenHijingEventHeader.h"
  43 #include "AliInputEventHandler.h"
  44
  45 #include "AliAnalysisHelperJetTasks.h"
  46 #include "AliAnalysisManager.h"
  47 #include "AliAnalysisTaskSE.h"
  48 #include "AliVParticle.h"
  49 #include "AliAODTrack.h"
  50 #include "AliVEvent.h"
  51
  52
  53 #include "AliAnalysisTaskJetProperties.h"
  54
  55 ClassImp(AliAnalysisTaskJetProperties)
  56
  57 //_________________________________________________________________________________//
  58
  59   AliAnalysisTaskJetProperties::AliAnalysisTaskJetProperties()
  60     : AliAnalysisTaskSE()
  61     ,fESD(0)
  62     ,fAOD(0)
  63     ,fAODJets(0)
  64     ,fAODExtension(0)
  65     ,fNonStdFile("")
  66     ,fBranchJets("jets")
  67     ,fTrackType(kTrackAOD)
  68     ,fJetRejectType(0)
  69     ,fUseAODInputJets(kTRUE)
  70     ,fFilterMask(0)
  71     ,fUsePhysicsSelection(kTRUE)
  72     ,fMaxVertexZ(10)
  73     ,fNContributors(2)
  74     ,fTrackPtCut(0)
  75     ,fTrackEtaMin(0)
  76     ,fTrackEtaMax(0)
  77     ,fJetPtCut(0)
  78     ,fJetEtaMin(0)
  79     ,fJetEtaMax(0)
  80     ,fAvgTrials(0)
  81     ,fJetRadius(0.4)
  82     ,fJetList(0)
  83     ,fTrackList(0)
  84     ,fTrackListUE(0)
  85     ,fTrackListJet(0)
  86     ,fCommonHistList(0)
  87     ,fh1EvtSelection(0)
  88     ,fh1VertexNContributors(0)
  89     ,fh1VertexZ(0)
  90     ,fh1Xsec(0)
  91     ,fh1Trials(0)
  92     ,fh1PtHard(0)
  93     ,fh1PtHardTrials(0)
  94     ,fh2EtaJet(0)
  95     ,fh2PhiJet(0)
  96     ,fh2PtJet(0)
  97     ,fh1PtJet(0)
  98     ,fh2NtracksJet(0)
  99     ,fProNtracksJet(0)
 100     ,fh2EtaTrack(0)
 101     ,fh2PhiTrack(0)
 102     ,fh2PtTrack(0)
 103     ,fh2FF(0)
 104     ,fh2Ksi(0)
 105     ,fh2DelEta(0)
 106     ,fh2DelPhi(0)
 107     ,fh2DelR(0)
 108
 109     ,fh1PtLeadingJet(0)
 110     ,fh2NtracksLeadingJet(0)
 111     ,fProNtracksLeadingJet(0)
 112     ,fh2DelR80pcNch(0)
 113     ,fProDelR80pcNch(0)
 114     ,fh2DelR80pcPt(0)
 115     ,fProDelR80pcPt(0)
 116     ,fh2AreaCh(0)
 117     ,fProAreaCh(0)
 118     ,fh3PtDelRNchSum(0)
 119     ,fh3PtDelRPtSum(0)
 120     ,fProDiffJetShape(0)
 121     ,fProIntJetShape(0)
 122
 123     ,fh1PtSumInJetConeUE(0)
 124     ,fh2NtracksLeadingJetUE(0)
 125     ,fProNtracksLeadingJetUE(0)
 126     ,fh2DelR80pcNchUE(0)
 127     ,fProDelR80pcNchUE(0)
 128     ,fh2DelR80pcPtUE(0)
 129     ,fProDelR80pcPtUE(0)
 130     ,fh2AreaChUE(0)
 131     ,fProAreaChUE(0)
 132     ,fh3PtDelRNchSumUE(0)
 133     ,fh3PtDelRPtSumUE(0)
 134     ,fProDiffJetShapeUE(0)
 135     ,fProIntJetShapeUE(0)
 136
 137     ,fh1CorrJetPt(0)
 138     ,fh2CorrPtTrack1(0)
 139     ,fh2CorrFF1(0)
 140     ,fh2CorrKsi1(0)
 141     ,fh2CorrjT1(0)
 142     ,fh1JetPtvsTrkSum(0)
 143
 144 {
 145   for(Int_t ii=0; ii<13; ii++){
 146     if(ii<6){
 147       fh2CorrPt1Pt2[ii]   = NULL;
 148       fh2CorrZ1Z2[ii]     = NULL;
 149       fh2CorrKsi1Ksi2[ii] = NULL;
 150       fh2CorrjT1jT2[ii]   = NULL;
 151     }
 152
 153     fProDelRNchSum[ii]    = NULL;
 154     fProDelRPtSum[ii]     = NULL;
 155     fProDiffJetShapeA[ii] = NULL;
 156     fProIntJetShapeA[ii]  = NULL;
 157
 158     fProDelRNchSumUE[ii]    = NULL;
 159     fProDelRPtSumUE[ii]     = NULL;
 160     fProDiffJetShapeAUE[ii] = NULL;
 161     fProIntJetShapeAUE[ii]  = NULL;
 162   }//ii loop
 163   // default constructor
 164 }
 165 //_________________________________________________________________________________//
 166
 167 AliAnalysisTaskJetProperties::AliAnalysisTaskJetProperties(const char *name)
 168   : AliAnalysisTaskSE(name)
 169   ,fESD(0)
 170   ,fAOD(0)
 171   ,fAODJets(0)
 172   ,fAODExtension(0)
 173   ,fNonStdFile("")
 174   ,fBranchJets("jets")
 175   ,fTrackType(kTrackAOD)
 176   ,fJetRejectType(0)
 177   ,fUseAODInputJets(kTRUE)
 178   ,fFilterMask(0)
 179   ,fUsePhysicsSelection(kTRUE)
 180   ,fMaxVertexZ(10)
 181   ,fNContributors(2)
 182   ,fTrackPtCut(0)
 183   ,fTrackEtaMin(0)
 184   ,fTrackEtaMax(0)
 185   ,fJetPtCut(0)
 186   ,fJetEtaMin(0)
 187   ,fJetEtaMax(0)
 188   ,fAvgTrials(0)
 189   ,fJetRadius(0.4)
 190   ,fJetList(0)
 191   ,fTrackList(0)
 192   ,fTrackListUE(0)
 193   ,fTrackListJet(0)
 194   ,fCommonHistList(0)
 195   ,fh1EvtSelection(0)
 196   ,fh1VertexNContributors(0)
 197   ,fh1VertexZ(0)
 198   ,fh1Xsec(0)
 199   ,fh1Trials(0)
 200   ,fh1PtHard(0)
 201   ,fh1PtHardTrials(0)
 202   ,fh2EtaJet(0)
 203   ,fh2PhiJet(0)
 204   ,fh2PtJet(0)
 205   ,fh1PtJet(0)
 206   ,fh2NtracksJet(0)
 207   ,fProNtracksJet(0)
 208   ,fh2EtaTrack(0)
 209   ,fh2PhiTrack(0)
 210   ,fh2PtTrack(0)
 211   ,fh2FF(0)
 212   ,fh2Ksi(0)
 213   ,fh2DelEta(0)
 214   ,fh2DelPhi(0)
 215   ,fh2DelR(0)
 216
 217   ,fh1PtLeadingJet(0)
 218   ,fh2NtracksLeadingJet(0)
 219   ,fProNtracksLeadingJet(0)
 220   ,fh2DelR80pcNch(0)
 221   ,fProDelR80pcNch(0)
 222   ,fh2DelR80pcPt(0)
 223   ,fProDelR80pcPt(0)
 224   ,fh2AreaCh(0)
 225   ,fProAreaCh(0)
 226   ,fh3PtDelRNchSum(0)
 227   ,fh3PtDelRPtSum(0)
 228   ,fProDiffJetShape(0)
 229   ,fProIntJetShape(0)
 230
 231   ,fh1PtSumInJetConeUE(0)
 232   ,fh2NtracksLeadingJetUE(0)
 233   ,fProNtracksLeadingJetUE(0)
 234   ,fh2DelR80pcNchUE(0)
 235   ,fProDelR80pcNchUE(0)
 236   ,fh2DelR80pcPtUE(0)
 237   ,fProDelR80pcPtUE(0)
 238   ,fh2AreaChUE(0)
 239   ,fProAreaChUE(0)
 240   ,fh3PtDelRNchSumUE(0)
 241   ,fh3PtDelRPtSumUE(0)
 242   ,fProDiffJetShapeUE(0)
 243   ,fProIntJetShapeUE(0)
 244
 245   ,fh1CorrJetPt(0)
 246   ,fh2CorrPtTrack1(0)
 247   ,fh2CorrFF1(0)
 248   ,fh2CorrKsi1(0)
 249   ,fh2CorrjT1(0)
 250   ,fh1JetPtvsTrkSum(0)
 251 {
 252   for(Int_t ii=0; ii<13; ii++){
 253     if(ii<6){
 254       fh2CorrPt1Pt2[ii]   = NULL;
 255       fh2CorrZ1Z2[ii]     = NULL;
 256       fh2CorrKsi1Ksi2[ii] = NULL;
 257       fh2CorrjT1jT2[ii]   = NULL;
 258     }
 259
 260     fProDelRNchSum[ii]    = NULL;
 261     fProDelRPtSum[ii]     = NULL;
 262     fProDiffJetShapeA[ii] = NULL;
 263     fProIntJetShapeA[ii]  = NULL;
 264
 265     fProDelRNchSumUE[ii]    = NULL;
 266     fProDelRPtSumUE[ii]     = NULL;
 267     fProDiffJetShapeAUE[ii] = NULL;
 268     fProIntJetShapeAUE[ii]  = NULL;
 269   }//ii loop
 270   // constructor
 271   DefineOutput(1,TList::Class());
 272 }
 273 //_________________________________________________________________________________//
 274
 275 AliAnalysisTaskJetProperties::~AliAnalysisTaskJetProperties()
 276 {
 277   // destructor
 278   if(fJetList)   delete fJetList;
 279   if(fTrackList)  delete fTrackList;
 280   if(fTrackListUE) delete fTrackListUE;
 281   if(fTrackListJet) delete fTrackListJet;
 282 }
 283 //_________________________________________________________________________________//
 284
 285 Bool_t AliAnalysisTaskJetProperties::Notify()
 286 {
 287   //
 288   // Implemented Notify() to read the cross sections
 289   // and number of trials from pyxsec.root
 290   // (taken from AliAnalysisTaskJetSpectrum2)
 291   //
 292   if(fDebug > 1) Printf("AliAnalysisTaskJetProperties::Notify()");
 293
 294   TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
 295   Float_t xsection = 0;
 296   Float_t ftrials  = 1;
 297
 298   fAvgTrials = 1;
 299   if(tree){
 300     TFile *curfile = tree->GetCurrentFile();
 301     if (!curfile) {
 302       Error("Notify","No current file");
 303       return kFALSE;
 304     }
 305     if(!fh1Xsec||!fh1Trials){
 306       Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
 307       return kFALSE;
 308     }
 309     AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
 310     fh1Xsec->Fill("<#sigma>",xsection);
 311     // construct a poor man average trials
 312     Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
 313     if(ftrials>=nEntries && nEntries>0.)fAvgTrials = ftrials/nEntries;
 314   }
 315   return kTRUE;
 316 }
 317 //_________________________________________________________________________________//
 318
 319 void AliAnalysisTaskJetProperties::UserCreateOutputObjects()
 320 {
 321   //(Here, event selection part is taken from AliAnalysisTaskFramentationFunction)
 322   // create output objects
 323   if(fDebug > 1) Printf("AliAnalysisTaskJetProperties::UserCreateOutputObjects()");
 324
 325   // create list of tracks and jets
 326   fJetList = new TList();
 327   fJetList->SetOwner(kFALSE);
 328   fTrackList = new TList();
 329   fTrackList->SetOwner(kFALSE);
 330   fTrackListUE = new TList();
 331   fTrackListUE->SetOwner(kFALSE);
 332   fTrackListJet = new TList();
 333   fTrackListJet->SetOwner(kFALSE);
 334
 335   // Create histograms / output container
 336   OpenFile(1);
 337   fCommonHistList = new TList();
 338   fCommonHistList->SetOwner();
 339
 340   Bool_t oldStatus = TH1::AddDirectoryStatus();
 341   TH1::AddDirectory(kFALSE);
 342   // Histograms/TProfiles
 343   fh1EvtSelection            = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
 344   fh1EvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
 345   fh1EvtSelection->GetXaxis()->SetBinLabel(2,"event selection: rejected");
 346   fh1EvtSelection->GetXaxis()->SetBinLabel(3,"event class: rejected");
 347   fh1EvtSelection->GetXaxis()->SetBinLabel(4,"vertex Ncontr: rejected");
 348   fh1EvtSelection->GetXaxis()->SetBinLabel(5,"vertex z: rejected");
 349   fh1EvtSelection->GetXaxis()->SetBinLabel(6,"vertex type: rejected");
 350
 351
 352   fh1VertexNContributors     = new TH1F("fh1VertexNContributors", "Vertex N contributors", 2500,-.5, 2499.5);
 353   fh1VertexZ                 = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
 354   fh1Xsec                    = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
 355   fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
 356   fh1Trials                  = new TH1F("fh1Trials","trials from pyxsec.root",1,0,1);
 357   fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
 358   fh1PtHard                  = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
 359   fh1PtHardTrials            = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
 360
 361
 362   Int_t kNbinsPtSlice=20; Float_t xMinPtSlice=0.0; Float_t xMaxPtSlice=100.0;
 363   Int_t kNbinsEta=40;     Float_t xMinEta=-2.0;    Float_t xMaxEta=2.0;
 364   Int_t kNbinsPhi=90;     Float_t xMinPhi=0.0;     Float_t xMaxPhi=TMath::TwoPi();
 365   Int_t kNbinsPt=250;     Float_t xMinPt=0.0;      Float_t xMaxPt=250.0;
 366   Int_t kNbinsNtracks=50; Float_t xMinNtracks=0.0; Float_t xMaxNtracks=50.0;
 367   Int_t kNbinsFF=200;     Float_t xMinFF=-0.05;    Float_t xMaxFF=1.95;
 368   Int_t kNbinsKsi=80;     Float_t xMinKsi=0.;      Float_t xMaxKsi = 8.0;
 369   Int_t kNbinsDelR1D=50;  Float_t xMinDelR1D=0.0;  Float_t xMaxDelR1D=1.0;
 370   Int_t kNbinsjT=100;     Float_t xMinjT=0.0;      Float_t xMaxjT=10.0;
 371
 372   fh2EtaJet      = new TH2F("EtaJet","EtaJet",
 373                             kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 374                             kNbinsEta,     xMinEta,     xMaxEta);
 375   fh2PhiJet      = new TH2F("PhiJet","PhiJet",
 376                             kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 377                             kNbinsPhi,     xMinPhi,     xMaxPhi);
 378   fh2PtJet       = new TH2F("PtJet","PtJet",
 379                             kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 380                             kNbinsPt,      xMinPt,      xMaxPt);
 381   fh1PtJet       = new TH1F("PtJet1D","PtJet1D",
 382                             kNbinsPt,      xMinPt,      xMaxPt);
 383   fh2NtracksJet  = new TH2F("NtracksJet","NtracksJet",
 384                             kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 385                             kNbinsNtracks, xMinNtracks, xMaxNtracks);
 386   fProNtracksJet = new TProfile("AvgNoOfTracksJet","AvgNoOfTracksJet",
 387                                 kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 388                                 xMinNtracks, xMaxNtracks);
 389   fh2EtaTrack    = new TH2F("EtaTrack","EtaTrack",
 390                             kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 391                             kNbinsEta,     xMinEta,     xMaxEta);
 392   fh2PhiTrack    = new TH2F("PhiTrack","PhiTrack",
 393                             kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 394                             kNbinsPhi,     xMinPhi,     xMaxPhi);
 395   fh2PtTrack     = new TH2F("PtTrack","PtTrack",
 396                             kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 397                             2*kNbinsPt,      xMinPt,      xMaxPt);
 398   fh2FF          = new TH2F("FF","FF",
 399                             kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 400                             kNbinsFF,      xMinFF,      xMaxFF);
 401   fh2Ksi         = new TH2F("Ksi","Ksi",
 402                             kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 403                             kNbinsKsi,      xMinKsi,      xMaxKsi);
 404   fh2DelEta      = new TH2F("DelEta","DelEta",
 405                             kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 406                             kNbinsEta,     xMinEta,     xMaxEta);
 407   fh2DelPhi      = new TH2F("DelPhi","DelPhi",
 408                             kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 409                             kNbinsPhi,     xMinPhi,     xMaxPhi);
 410   fh2DelR        = new TH2F("DelR","DelR",
 411                             kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 412                             kNbinsDelR1D,  xMinDelR1D,  xMaxDelR1D);
 413
 414
 415
 416   Int_t kNbinsDelR=100;      Float_t xMinDelR=0.0;      Float_t xMaxDelR=1.0;
 417   Int_t kNbinsPtSliceJS=100; Float_t xMinPtSliceJS=0.0; Float_t xMaxPtSliceJS=250.0;
 418
 419   fh1PtLeadingJet       = new TH1F("PtLeadingJet","PtLeadingJet",
 420                                    kNbinsPt,      xMinPt,      xMaxPt);
 421   fh2NtracksLeadingJet  = new TH2F("NtracksLeadingJet","NtracksLeadingJet",
 422                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 423                                    kNbinsNtracks,   xMinNtracks,   xMaxNtracks);
 424   fProNtracksLeadingJet = new TProfile("AvgNoOfTracksLeadingJet","AvgNoOfTracksLeadingJet",
 425                                        kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 426                                        xMinNtracks, xMaxNtracks);
 427   fh2DelR80pcNch        = new TH2F("delR80pcNch","delR80pcNch",
 428                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 429                                    kNbinsDelR,      xMinDelR,      xMaxDelR);
 430   fProDelR80pcNch       = new TProfile("AvgdelR80pcNch","AvgdelR80pcNch",
 431                                        kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 432                                        xMinDelR,  xMaxDelR);
 433   fh2DelR80pcPt         = new TH2F("delR80pcPt","delR80pcPt",
 434                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 435                                    kNbinsDelR,      xMinDelR,      xMaxDelR);
 436   fProDelR80pcPt        = new TProfile("AvgdelR80pcPt","AvgdelR80pcPt",
 437                                        kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 438                                        xMinDelR,  xMaxDelR);
 439   fh2AreaCh             = new TH2F("Area","Area",
 440                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 441                                    72, 0.0, TMath::Pi());
 442   fProAreaCh            = new TProfile("AvgArea","AvgArea",
 443                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 444                                    0.0, TMath::Pi());
 445   fh3PtDelRNchSum       = new TH3F("PtDelRNchSum","PtDelRNchSum",
 446                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 447                                    kNbinsDelR1D,    xMinDelR1D,    xMaxDelR1D,
 448                                    kNbinsNtracks, xMinNtracks, xMaxNtracks);
 449   fh3PtDelRPtSum        = new TH3F("PtDelRPtSum","PtDelRPtSum",
 450                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 451                                    kNbinsDelR1D,    xMinDelR1D,    xMaxDelR1D,
 452                                    kNbinsPt,        xMinPt,        xMaxPt);
 453   fProDiffJetShape      = new TProfile("DiffJetShape","DiffJetShape",
 454                                   10,0.0,1.0,0.0,250.0);
 455   fProIntJetShape       = new TProfile("IntJetShape","IntJetShape",
 456                                   10,0.0,1.0,0.0,250.0);
 457
 458
 459   fh1PtSumInJetConeUE       = new TH1F("PtSumInJetConeUE","PtSumInJetConeUE",
 460                                        500,      0.,      100.);
 461   fh2NtracksLeadingJetUE  = new TH2F("NtracksLeadingJetUE","NtracksLeadingJetUE",
 462                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 463                                    kNbinsNtracks,   xMinNtracks,   xMaxNtracks);
 464   fProNtracksLeadingJetUE = new TProfile("AvgNoOfTracksLeadingJetUE","AvgNoOfTracksLeadingJetUE",
 465                                        kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 466                                        xMinNtracks, xMaxNtracks);
 467   fh2DelR80pcNchUE        = new TH2F("delR80pcNchUE","delR80pcNchUE",
 468                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 469                                    kNbinsDelR,      xMinDelR,      xMaxDelR);
 470   fProDelR80pcNchUE       = new TProfile("AvgdelR80pcNchUE","AvgdelR80pcNchUE",
 471                                        kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 472                                        xMinDelR,  xMaxDelR);
 473   fh2DelR80pcPtUE         = new TH2F("delR80pcPtUE","delR80pcPtUE",
 474                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 475                                    kNbinsDelR,      xMinDelR,      xMaxDelR);
 476   fProDelR80pcPtUE        = new TProfile("AvgdelR80pcPtUE","AvgdelR80pcPtUE",
 477                                        kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 478                                        xMinDelR,  xMaxDelR);
 479   fh2AreaChUE             = new TH2F("AreaUE","AreaUE",
 480                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 481                                    72, 0.0, TMath::Pi());
 482   fProAreaChUE            = new TProfile("AvgAreaUE","AvgAreaUE",
 483                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 484                                    0.0, TMath::Pi());
 485   fh3PtDelRNchSumUE       = new TH3F("PtDelRNchSumUE","PtDelRNchSumUE",
 486                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 487                                    kNbinsDelR1D,    xMinDelR1D,    xMaxDelR1D,
 488                                    kNbinsNtracks, xMinNtracks, xMaxNtracks);
 489   fh3PtDelRPtSumUE        = new TH3F("PtDelRPtSumUE","PtDelRPtSumUE",
 490                                    kNbinsPtSliceJS, xMinPtSliceJS, xMaxPtSliceJS,
 491                                    kNbinsDelR1D,    xMinDelR1D,    xMaxDelR1D,
 492                                    kNbinsPt,        xMinPt,        xMaxPt);
 493   fProDiffJetShapeUE      = new TProfile("DiffJetShapeUE","DiffJetShapeUE",
 494                                   10,0.0,1.0,0.0,250.0);
 495   fProIntJetShapeUE       = new TProfile("IntJetShapeUE","IntJetShapeUE",
 496                                   10,0.0,1.0,0.0,250.0);
 497
 498
 499   fh1CorrJetPt    = new TH1F("JetPtCorr","JetPtCorr",
 500                              kNbinsPt,      xMinPt,      xMaxPt);
 501   fh2CorrPtTrack1 = new TH2F("TrackPtCorr", "TrackPtCorr",
 502                              kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 503                              2*kNbinsPt,      xMinPt,      xMaxPt);
 504   fh2CorrFF1      = new TH2F("FFCorr","FFCorr",
 505                              kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 506                              kNbinsFF,      xMinFF,      xMaxFF);
 507   fh2CorrKsi1     = new TH2F("KsiCorr","KsiCorr",
 508                              kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 509                              kNbinsKsi,      xMinKsi,      xMaxKsi);
 510   fh2CorrjT1      = new TH2F("jTCorr","jTCorr",
 511                              kNbinsPtSlice, xMinPtSlice, xMaxPtSlice,
 512                              kNbinsjT,      xMinjT,      xMaxjT);
 513
 514   fh1JetPtvsTrkSum = new TH1F("diffJetPt_TrkPtSum","diffJetPt_TrkPtSum",
 515                               500, -250., 250.);
 516   TString title;
 517   for(Int_t ii=0; ii<13; ii++){
 518     if(ii==0)title = "_JetPt20to25";
 519     if(ii==1)title = "_JetPt25to30";
 520     if(ii==2)title = "_JetPt20to30";
 521
 522     if(ii==3)title = "_JetPt30to35";
 523     if(ii==4)title = "_JetPt35to40";
 524     if(ii==5)title = "_JetPt30to40";
 525
 526     if(ii==6)title = "_JetPt40to50";
 527     if(ii==7)title = "_JetPt50to60";
 528     if(ii==8)title = "_JetPt40to60";
 529
 530     if(ii==9) title = "_JetPt60to70";
 531     if(ii==10)title = "_JetPt70to80";
 532     if(ii==11)title = "_JetPt60to80";
 533
 534     if(ii==12)title = "_JetPt80to100";
 535     fProDelRNchSum[ii] = new TProfile(Form("AvgNchSumDelR%s",title.Data()),Form("AvgNchSumDelR%s",title.Data()),
 536                                       kNbinsDelR1D,  xMinDelR1D,  xMaxDelR1D,
 537                                       xMinNtracks, xMaxNtracks);
 538
 539     fProDelRPtSum[ii]  = new TProfile(Form("AvgPtSumDelR%s",title.Data()),Form("AvgPtSumDelR%s",title.Data()),
 540                                       kNbinsDelR1D,  xMinDelR1D,  xMaxDelR1D,
 541                                       xMinPt,      xMaxPt);
 542     fProDelRNchSum[ii] ->GetXaxis()->SetTitle("R");
 543     fProDelRNchSum[ii] ->GetYaxis()->SetTitle("<NchSum>");
 544     fProDelRPtSum[ii]  ->GetXaxis()->SetTitle("R");
 545     fProDelRPtSum[ii]  ->GetYaxis()->SetTitle("<PtSum>");
 546
 547     fProDiffJetShapeA[ii] = new TProfile(Form("DiffJetShape%s",title.Data()),Form("DiffJetShape%s",title.Data()),
 548                                          10,0.0,1.0,0.0,250.0);
 549     fProIntJetShapeA[ii]  = new TProfile(Form("IntJetShape%s",title.Data()),Form("IntJetShape%s",title.Data()),
 550                                          10,0.0,1.0,0.0,250.0);
 551
 552     fProDiffJetShapeA[ii]->GetXaxis()->SetTitle("R");
 553     fProDiffJetShapeA[ii]->GetYaxis()->SetTitle("Diff jet shape");
 554     fProIntJetShapeA[ii]->GetXaxis()->SetTitle("R");
 555     fProIntJetShapeA[ii]->GetYaxis()->SetTitle("Integrated jet shape");
 556
 557     fCommonHistList->Add(fProDelRNchSum[ii]);
 558     fCommonHistList->Add(fProDelRPtSum[ii]);
 559     fCommonHistList->Add(fProDiffJetShapeA[ii]);
 560     fCommonHistList->Add(fProIntJetShapeA[ii]);
 561
 562     fProDelRNchSumUE[ii] = new TProfile(Form("AvgNchSumDelR%sUE",title.Data()),Form("AvgNchSumDelR%sUE",title.Data()),
 563                                       kNbinsDelR1D,  xMinDelR1D,  xMaxDelR1D,
 564                                       xMinNtracks, xMaxNtracks);
 565
 566     fProDelRPtSumUE[ii]  = new TProfile(Form("AvgPtSumDelR%sUE",title.Data()),Form("AvgPtSumDelR%sUE",title.Data()),
 567                                       kNbinsDelR1D,  xMinDelR1D,  xMaxDelR1D,
 568                                       xMinPt,      xMaxPt);
 569     fProDelRNchSumUE[ii] ->GetXaxis()->SetTitle("R");
 570     fProDelRNchSumUE[ii] ->GetYaxis()->SetTitle("<N_{ch}^{Sum,UE}>");
 571     fProDelRPtSumUE[ii]  ->GetXaxis()->SetTitle("R");
 572     fProDelRPtSumUE[ii]  ->GetYaxis()->SetTitle("<p_{T}^{sum, UE}>");
 573
 574     fProDiffJetShapeAUE[ii] = new TProfile(Form("DiffJetShape%sUE",title.Data()),Form("DiffJetShape%sUE",title.Data()),
 575                                          10,0.0,1.0,0.0,250.0);
 576     fProIntJetShapeAUE[ii]  = new TProfile(Form("IntJetShape%sUE",title.Data()),Form("IntJetShape%sUE",title.Data()),
 577                                          10,0.0,1.0,0.0,250.0);
 578
 579     fProDiffJetShapeAUE[ii]->GetXaxis()->SetTitle("R");
 580     fProDiffJetShapeAUE[ii]->GetYaxis()->SetTitle("Diff jet shape UE");
 581     fProIntJetShapeAUE[ii]->GetXaxis()->SetTitle("R");
 582     fProIntJetShapeAUE[ii]->GetYaxis()->SetTitle("Integrated jet shape UE");
 583
 584     fCommonHistList->Add(fProDelRNchSumUE[ii]);
 585     fCommonHistList->Add(fProDelRPtSumUE[ii]);
 586     fCommonHistList->Add(fProDiffJetShapeAUE[ii]);
 587     fCommonHistList->Add(fProIntJetShapeAUE[ii]);
 588   }//ii loop
 589
 590   fh2EtaJet     ->GetXaxis()->SetTitle("JetPt");  fh2EtaJet     ->GetYaxis()->SetTitle("JetEta");
 591   fh2PhiJet     ->GetXaxis()->SetTitle("JetPt");  fh2PhiJet     ->GetYaxis()->SetTitle("JetPhi");
 592   fh2PtJet      ->GetXaxis()->SetTitle("JetPt");  fh2PtJet      ->GetYaxis()->SetTitle("JetPt");
 593   fh1PtJet      ->GetXaxis()->SetTitle("JetPt");  fh1PtJet      ->GetYaxis()->SetTitle("#jets");
 594   fh2NtracksJet ->GetXaxis()->SetTitle("JetPt");  fh2NtracksJet ->GetYaxis()->SetTitle("#tracks");
 595   fProNtracksJet->GetXaxis()->SetTitle("JetPt");  fProNtracksJet->GetYaxis()->SetTitle("AgvNtracks");
 596   fh2EtaTrack   ->GetXaxis()->SetTitle("JetPt");  fh2EtaTrack   ->GetYaxis()->SetTitle("TrackEta");
 597   fh2PhiTrack   ->GetXaxis()->SetTitle("JetPt");  fh2PhiTrack   ->GetYaxis()->SetTitle("TrackPhi");
 598   fh2PtTrack    ->GetXaxis()->SetTitle("JetPt");  fh2PtTrack    ->GetYaxis()->SetTitle("TrackPt");
 599   fh2FF         ->GetXaxis()->SetTitle("JetPt");  fh2FF         ->GetYaxis()->SetTitle("FF");
 600   fh2Ksi        ->GetXaxis()->SetTitle("JetPt");  fh2Ksi        ->GetYaxis()->SetTitle("Ksi");
 601   fh2DelEta     ->GetXaxis()->SetTitle("JetPt");  fh2DelEta     ->GetYaxis()->SetTitle("DelEta");
 602   fh2DelPhi     ->GetXaxis()->SetTitle("JetPt");  fh2DelPhi     ->GetYaxis()->SetTitle("DelPhi");
 603   fh2DelR       ->GetXaxis()->SetTitle("JetPt");  fh2DelR       ->GetYaxis()->SetTitle("DelR");
 604
 605   fh1PtLeadingJet       ->GetXaxis()->SetTitle("JetPt");  fh1PtLeadingJet       ->GetYaxis()->SetTitle("#leading jets");
 606   fh2NtracksLeadingJet  ->GetXaxis()->SetTitle("JetPt");  fh2NtracksLeadingJet  ->GetYaxis()->SetTitle("#tracks leading jet");
 607   fProNtracksLeadingJet ->GetXaxis()->SetTitle("JetPt");  fProNtracksLeadingJet ->GetYaxis()->SetTitle("AvgNtracks leading jet");
 608   fh2DelR80pcNch        ->GetXaxis()->SetTitle("JetPt");  fh2DelR80pcNch        ->GetYaxis()->SetTitle("R containing 80% of tracks");
 609   fProDelR80pcNch       ->GetXaxis()->SetTitle("JetPt");  fProDelR80pcNch       ->GetYaxis()->SetTitle("<R> containing 80% of tracks");
 610   fh2DelR80pcPt         ->GetXaxis()->SetTitle("JetPt");  fh2DelR80pcPt         ->GetYaxis()->SetTitle("R containing 80% of pT");
 611   fProDelR80pcPt        ->GetXaxis()->SetTitle("JetPt");  fProDelR80pcPt        ->GetYaxis()->SetTitle("<R> containing 80% of pT");
 612   fh2AreaCh             ->GetXaxis()->SetTitle("JetPt");  fh2AreaCh             ->GetYaxis()->SetTitle("Jet area");
 613   fProAreaCh            ->GetXaxis()->SetTitle("JetPt");  fProAreaCh            ->GetYaxis()->SetTitle("<jet area>");
 614   fh3PtDelRNchSum       ->GetXaxis()->SetTitle("JetPt");  fh3PtDelRNchSum       ->GetYaxis()->SetTitle("R");  fh3PtDelRNchSum->GetZaxis()->SetTitle("NchSum");
 615   fh3PtDelRPtSum        ->GetXaxis()->SetTitle("JetPt");  fh3PtDelRPtSum        ->GetYaxis()->SetTitle("R");  fh3PtDelRPtSum ->GetZaxis()->SetTitle("PtSum");
 616   fProDiffJetShape->GetXaxis()->SetTitle("R");
 617   fProDiffJetShape->GetYaxis()->SetTitle("Diff jet shape");
 618   fProIntJetShape->GetXaxis()->SetTitle("R");
 619   fProIntJetShape->GetYaxis()->SetTitle("Integrated jet shape");
 620
 621   fh1PtSumInJetConeUE     ->GetXaxis()->SetTitle("p_{T}^{sum, UE}(in cone R)");  fh1PtSumInJetConeUE     ->GetYaxis()->SetTitle("#leading jets");
 622   fh2NtracksLeadingJetUE  ->GetXaxis()->SetTitle("JetPt");  fh2NtracksLeadingJetUE  ->GetYaxis()->SetTitle("#tracks UE");
 623   fProNtracksLeadingJetUE ->GetXaxis()->SetTitle("JetPt");  fProNtracksLeadingJetUE ->GetYaxis()->SetTitle("AvgNtracks UE");
 624   fh2DelR80pcNchUE        ->GetXaxis()->SetTitle("JetPt");  fh2DelR80pcNchUE        ->GetYaxis()->SetTitle("R containing 80% of tracks");
 625   fProDelR80pcNchUE       ->GetXaxis()->SetTitle("JetPt");  fProDelR80pcNchUE       ->GetYaxis()->SetTitle("<R> containing 80% of tracks");
 626   fh2DelR80pcPtUE         ->GetXaxis()->SetTitle("JetPt");  fh2DelR80pcPtUE         ->GetYaxis()->SetTitle("R containing 80% of pT");
 627   fProDelR80pcPtUE        ->GetXaxis()->SetTitle("JetPt");  fProDelR80pcPtUE        ->GetYaxis()->SetTitle("<R> containing 80% of pT");
 628   fh2AreaChUE             ->GetXaxis()->SetTitle("JetPt");  fh2AreaChUE             ->GetYaxis()->SetTitle("UE area");
 629   fProAreaChUE            ->GetXaxis()->SetTitle("JetPt");  fProAreaChUE            ->GetYaxis()->SetTitle("<UE area>");
 630   fh3PtDelRNchSumUE       ->GetXaxis()->SetTitle("JetPt");  fh3PtDelRNchSumUE       ->GetYaxis()->SetTitle("R");  fh3PtDelRNchSumUE->GetZaxis()->SetTitle("NchSumUE");
 631   fh3PtDelRPtSumUE        ->GetXaxis()->SetTitle("JetPt");  fh3PtDelRPtSumUE        ->GetYaxis()->SetTitle("R");  fh3PtDelRPtSumUE ->GetZaxis()->SetTitle("PtSumUE");
 632   fProDiffJetShapeUE->GetXaxis()->SetTitle("R");
 633   fProDiffJetShapeUE->GetYaxis()->SetTitle("Diff jet shape UE");
 634   fProIntJetShapeUE->GetXaxis()->SetTitle("R");
 635   fProIntJetShapeUE->GetYaxis()->SetTitle("Integrated jet shape UE");
 636
 637   fh1CorrJetPt    ->GetXaxis()->SetTitle("JetPt");  fh1CorrJetPt    ->GetYaxis()->SetTitle("#jets");
 638   fh2CorrPtTrack1 ->GetXaxis()->SetTitle("JetPt");  fh2CorrPtTrack1 ->GetYaxis()->SetTitle("pt_track");
 639   fh2CorrFF1      ->GetXaxis()->SetTitle("JetPt");  fh2CorrFF1      ->GetYaxis()->SetTitle("z_track");
 640   fh2CorrKsi1     ->GetXaxis()->SetTitle("JetPt");  fh2CorrKsi1     ->GetYaxis()->SetTitle("ksi_track");
 641   fh2CorrjT1      ->GetXaxis()->SetTitle("JetPt");  fh2CorrjT1      ->GetYaxis()->SetTitle("jt_track");
 642   fh1JetPtvsTrkSum->GetXaxis()->SetTitle("JetPt-TrkPtSum");
 643
 644   fCommonHistList->Add(fh1EvtSelection);
 645   fCommonHistList->Add(fh1VertexNContributors);
 646   fCommonHistList->Add(fh1VertexZ);
 647   fCommonHistList->Add(fh1Xsec);
 648   fCommonHistList->Add(fh1Trials);
 649   fCommonHistList->Add(fh1PtHard);
 650   fCommonHistList->Add(fh1PtHardTrials);
 651   fCommonHistList->Add(fh2EtaJet);
 652   fCommonHistList->Add(fh2PhiJet);
 653   fCommonHistList->Add(fh2PtJet);
 654   fCommonHistList->Add(fh1PtJet);
 655   fCommonHistList->Add(fh2NtracksJet);
 656   fCommonHistList->Add(fProNtracksJet);
 657   fCommonHistList->Add(fh2EtaTrack);
 658   fCommonHistList->Add(fh2PhiTrack);
 659   fCommonHistList->Add(fh2PtTrack);
 660   fCommonHistList->Add(fh2FF);
 661   fCommonHistList->Add(fh2Ksi);
 662   fCommonHistList->Add(fh2DelEta);
 663   fCommonHistList->Add(fh2DelPhi);
 664   fCommonHistList->Add(fh2DelR);
 665
 666   fCommonHistList->Add(fh1PtLeadingJet);
 667   fCommonHistList->Add(fh2NtracksLeadingJet);
 668   fCommonHistList->Add(fProNtracksLeadingJet);
 669   fCommonHistList->Add(fh2DelR80pcNch);
 670   fCommonHistList->Add(fProDelR80pcNch);
 671   fCommonHistList->Add(fh2DelR80pcPt);
 672   fCommonHistList->Add(fProDelR80pcPt);
 673   fCommonHistList->Add(fh2AreaCh);
 674   fCommonHistList->Add(fProAreaCh);
 675   fCommonHistList->Add(fh3PtDelRNchSum);
 676   fCommonHistList->Add(fh3PtDelRPtSum);
 677   fCommonHistList->Add(fProDiffJetShape);
 678   fCommonHistList->Add(fProIntJetShape);
 679
 680   fCommonHistList->Add(fh1PtSumInJetConeUE);
 681   fCommonHistList->Add(fh2NtracksLeadingJetUE);
 682   fCommonHistList->Add(fProNtracksLeadingJetUE);
 683   fCommonHistList->Add(fh2DelR80pcNchUE);
 684   fCommonHistList->Add(fProDelR80pcNchUE);
 685   fCommonHistList->Add(fh2DelR80pcPtUE);
 686   fCommonHistList->Add(fProDelR80pcPtUE);
 687   fCommonHistList->Add(fh2AreaChUE);
 688   fCommonHistList->Add(fProAreaChUE);
 689   fCommonHistList->Add(fh3PtDelRNchSumUE);
 690   fCommonHistList->Add(fh3PtDelRPtSumUE);
 691   fCommonHistList->Add(fProDiffJetShapeUE);
 692   fCommonHistList->Add(fProIntJetShapeUE);
 693
 694   fCommonHistList->Add(fh1CorrJetPt);
 695   fCommonHistList->Add(fh2CorrPtTrack1);
 696   fCommonHistList->Add(fh2CorrFF1);
 697   fCommonHistList->Add(fh2CorrKsi1);
 698   fCommonHistList->Add(fh2CorrjT1);
 699   fCommonHistList->Add(fh1JetPtvsTrkSum);
 700
 701
 702
 703   TString titleCorr;
 704   for(Int_t jj=0; jj<6; jj++){
 705     if(jj == 0)titleCorr = "_JetPt10to20";
 706     if(jj == 1)titleCorr = "_JetPt20to30";
 707     if(jj == 2)titleCorr = "_JetPt30to40";
 708     if(jj == 3)titleCorr = "_JetPt40to60";
 709     if(jj == 4)titleCorr = "_JetPt60to80";
 710     if(jj == 5)titleCorr = "_JetPt80to100";
 711
 712     fh2CorrPt1Pt2[jj]   = new TH2F(Form("CorrPt1Pt2%s",titleCorr.Data()),Form("CorrPt1Pt2%s",titleCorr.Data()),
 713                                    2*kNbinsPt,      xMinPt,      xMaxPt,
 714                                    2*kNbinsPt,      xMinPt,      xMaxPt);
 715     fh2CorrZ1Z2[jj]     = new TH2F(Form("CorrZ1Z2%s",titleCorr.Data()),Form("CorrZ1Z2%s",titleCorr.Data()),
 716                                    kNbinsFF,      xMinFF,      xMaxFF,
 717                                    kNbinsFF,      xMinFF,      xMaxFF);
 718     fh2CorrKsi1Ksi2[jj] = new TH2F(Form("CorrKsi1Ksi2%s",titleCorr.Data()),Form("CorrKsi1Ksi2%s",titleCorr.Data()),
 719                                    kNbinsKsi,      xMinKsi,      xMaxKsi,
 720                                    kNbinsKsi,      xMinKsi,      xMaxKsi);
 721     fh2CorrjT1jT2[jj]   = new TH2F(Form("CorrjT1jT2%s",titleCorr.Data()),Form("CorrjT1jT2%s",titleCorr.Data()),
 722                                    kNbinsjT,      xMinjT,      xMaxjT,
 723                                    kNbinsjT,      xMinjT,      xMaxjT);
 724
 725
 726     fh2CorrPt1Pt2[jj]   ->GetXaxis()->SetTitle("pt_track1");
 727     fh2CorrZ1Z2[jj]     ->GetXaxis()->SetTitle("z_track1");
 728     fh2CorrKsi1Ksi2[jj] ->GetXaxis()->SetTitle("ksi_track1");
 729     fh2CorrjT1jT2[jj]   ->GetXaxis()->SetTitle("jt_track1");
 730
 731     fh2CorrPt1Pt2[jj]   ->GetYaxis()->SetTitle("pt_track2");
 732     fh2CorrZ1Z2[jj]     ->GetYaxis()->SetTitle("z_track2");
 733     fh2CorrKsi1Ksi2[jj] ->GetYaxis()->SetTitle("ksi_track2");
 734     fh2CorrjT1jT2[jj]   ->GetYaxis()->SetTitle("jt_track2");
 735
 736     fCommonHistList->Add(fh2CorrPt1Pt2[jj]);
 737     fCommonHistList->Add(fh2CorrZ1Z2[jj]);
 738     fCommonHistList->Add(fh2CorrKsi1Ksi2[jj]);
 739     fCommonHistList->Add(fh2CorrjT1jT2[jj]);
 740   }//jj loop
 741
 742
 743
 744   // =========== Switch on Sumw2 for all histos ===========
 745   for (Int_t i=0; i<fCommonHistList->GetEntries(); ++i){
 746     TH1 *h1 = dynamic_cast<TH1*>(fCommonHistList->At(i));
 747     if (h1) h1->Sumw2();
 748     else{
 749       TProfile *hPro = dynamic_cast<TProfile*>(fCommonHistList->At(i));
 750       if(hPro) hPro->Sumw2();
 751     }
 752   }
 753
 754   TH1::AddDirectory(oldStatus);
 755   PostData(1, fCommonHistList);
 756 }
 757 //_________________________________________________________________________________//
 758
 759 void AliAnalysisTaskJetProperties::Init()
 760 {
 761   // Initialization
 762   if(fDebug > 1) Printf("AliAnalysisTaskJetProperties::Init()");
 763
 764 }
 765 //_________________________________________________________________________________//
 766
 767 void AliAnalysisTaskJetProperties::UserExec(Option_t *)
 768 {
 769   // Main loop
 770   // Called for each event
 771   if(fDebug > 2) Printf("AliAnalysisTaskJetProperties::UserExec()");
 772
 773   //if(fDebug > 1) Printf("Analysis event #%5d", (Int_t) fEntry);
 774   // Trigger selection
 775
 776   AliInputEventHandler* inputHandler = (AliInputEventHandler*)
 777     ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
 778   if(!(inputHandler->IsEventSelected() & AliVEvent::kMB)){
 779     if(inputHandler->InheritsFrom("AliESDInputHandler") && fUsePhysicsSelection){ // PhysicsSelection only with ESD input
 780       fh1EvtSelection->Fill(1.);
 781       if (fDebug > 2 ) Printf(" Trigger Selection: event REJECTED ... ");
 782       PostData(1, fCommonHistList);
 783       return;
 784     }
 785   }
 786
 787   fESD = dynamic_cast<AliESDEvent*>(InputEvent());
 788   if(!fESD){
 789     if(fDebug>2) Printf("%s:%d ESDEvent not found in the input", (char*)__FILE__,__LINE__);
 790   }
 791
 792   fMCEvent = MCEvent();
 793   if(!fMCEvent){
 794     if(fDebug>2) Printf("%s:%d MCEvent not found in the input", (char*)__FILE__,__LINE__);
 795   }
 796
 797   // get AOD event from input/ouput
 798     TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
 799     if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
 800       fAOD  =  ((AliAODInputHandler*)handler)->GetEvent();
 801       if(fUseAODInputJets) fAODJets = fAOD;
 802       if (fDebug > 2)  Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
 803     }
 804     else {
 805       handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
 806       if( handler && handler->InheritsFrom("AliAODHandler") ) {
 807         fAOD = ((AliAODHandler*)handler)->GetAOD();
 808         fAODJets = fAOD;
 809         if (fDebug > 2)  Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
 810       }
 811     }
 812
 813     if(!fAODJets && !fUseAODInputJets){ // case we have AOD in input & output and want jets from output
 814     TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
 815     if( outHandler && outHandler->InheritsFrom("AliAODHandler") ) {
 816       fAODJets = ((AliAODHandler*)outHandler)->GetAOD();
 817       if (fDebug > 2)  Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
 818     }
 819     }
 820
 821     if(fNonStdFile.Length()!=0){
 822       // case we have an AOD extension - fetch the jets from the extended output
 823
 824       AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
 825       fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
 826       if(!fAODExtension){
 827         if(fDebug>2)Printf("AODExtension not found for %s",fNonStdFile.Data());
 828       }
 829     }
 830
 831     if(!fAOD){
 832       Printf("%s:%d AODEvent not found", (char*)__FILE__,__LINE__);
 833       return;
 834     }
 835     if(!fAODJets){
 836       Printf("%s:%d AODEvent with jet branch not found", (char*)__FILE__,__LINE__);
 837       return;
 838     }
 839
 840
 841     // *** vertex cut ***
 842     AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
 843     Int_t nTracksPrim = primVtx->GetNContributors();
 844     fh1VertexNContributors->Fill(nTracksPrim);
 845
 846   if (fDebug > 2) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
 847   //if(!nTracksPrim){
 848   if(nTracksPrim<fNContributors){
 849     if (fDebug > 2) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__);
 850     fh1EvtSelection->Fill(3.);
 851     PostData(1, fCommonHistList);
 852     return;
 853   }
 854
 855   fh1VertexZ->Fill(primVtx->GetZ());
 856
 857   if(TMath::Abs(primVtx->GetZ())>fMaxVertexZ){
 858     if (fDebug > 2) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
 859     fh1EvtSelection->Fill(4.);
 860     PostData(1, fCommonHistList);
 861     return;
 862   }
 863
 864   TString primVtxName(primVtx->GetName());
 865
 866   if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
 867     if (fDebug > 2) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
 868     fh1EvtSelection->Fill(5.);
 869     PostData(1, fCommonHistList);
 870     return;
 871   }
 872
 873   if (fDebug > 2) Printf("%s:%d event ACCEPTED ...",(char*)__FILE__,__LINE__);
 874   fh1EvtSelection->Fill(0.);
 875   //___ get MC information __________________________________________________________________
 876   Double_t ptHard = 0.;
 877   Double_t nTrials = 1; // trials for MC trigger weight for real data
 878
 879   if(fMCEvent){
 880     AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
 881     if(genHeader){
 882       AliGenPythiaEventHeader*  pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
 883       AliGenHijingEventHeader*  hijingGenHeader = 0x0;
 884
 885       if(pythiaGenHeader){
 886         if(fDebug>2) Printf("%s:%d pythiaGenHeader found", (char*)__FILE__,__LINE__);
 887         nTrials = pythiaGenHeader->Trials();
 888         ptHard  = pythiaGenHeader->GetPtHard();
 889
 890         fh1PtHard->Fill(ptHard);
 891         fh1PtHardTrials->Fill(ptHard,nTrials);
 892
 893
 894       } else { // no pythia, hijing?
 895
 896         if(fDebug>2) Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
 897
 898         hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
 899         if(!hijingGenHeader){
 900           Printf("%s:%d no pythiaGenHeader or hjingGenHeader found", (char*)__FILE__,__LINE__);
 901         } else {
 902           if(fDebug>2) Printf("%s:%d hijingGenHeader found", (char*)__FILE__,__LINE__);
 903         }
 904       }
 905
 906       fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
 907     }
 908   }
 909
 910   //___ fetch jets __________________________________________________________________________
 911   Int_t nJ = GetListOfJets(fJetList);
 912   Int_t nJets = 0;
 913   if(nJ>=0) nJets = fJetList->GetEntries();
 914   if(fDebug>2){
 915     Printf("%s:%d Selected jets: %d %d",(char*)__FILE__,__LINE__,nJ,nJets);
 916     if(nJ != nJets) Printf("%s:%d Mismatch Selected Jets: %d %d",(char*)__FILE__,__LINE__,nJ,nJets);
 917   }
 918   FillJetProperties(fJetList);
 919   FillJetShape(fJetList);
 920   FillJetShapeUE(fJetList);
 921   FillFFCorr(fJetList);
 922   fJetList->Clear();
 923   //Post output data.
 924   PostData(1, fCommonHistList);
 925 }//UserExec
 926 //_________________________________________________________________________________//
 927
 928 void AliAnalysisTaskJetProperties::Terminate(Option_t *)
 929 {
 930   // terminated
 931
 932   if(fDebug > 1) printf("AliAnalysisTaskJetProperties::Terminate() \n");
 933 }
 934 //_________________________________________________________________________________//
 935
 936 Int_t AliAnalysisTaskJetProperties::GetListOfJets(TList *list)
 937 {
 938   //this functionality is motivated by AliAnalysisTaskFragmentationFunction
 939   if(fDebug > 2) printf("AliAnalysisTaskJetProperties::GetListOfJets() \n");
 940   // fill list of jets selected according to type
 941   if(!list){
 942     if(fDebug>2) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
 943     return -1;
 944   }
 945
 946   if(fBranchJets.Length()==0){
 947     Printf("%s:%d no jet branch specified", (char*)__FILE__,__LINE__);
 948     if(fDebug>2)fAOD->Print();
 949     return 0;
 950   }
 951   TClonesArray *aodJets = 0;
 952   if(fBranchJets.Length())      aodJets = dynamic_cast<TClonesArray*>(fAODJets->FindListObject(fBranchJets.Data()));
 953   if(!aodJets)                  aodJets = dynamic_cast<TClonesArray*>(fAODJets->GetList()->FindObject(fBranchJets.Data()));
 954   if(fAODExtension&&!aodJets)   aodJets = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fBranchJets.Data()));
 955
 956   if(!aodJets){
 957     if(fBranchJets.Length())Printf("%s:%d no jet array with name %s in AOD", (char*)__FILE__,__LINE__,fBranchJets.Data());
 958     if(fDebug>2)fAOD->Print();
 959     return 0;
 960   }
 961
 962   Int_t nJets = 0;
 963   for(Int_t ijet=0; ijet<aodJets->GetEntries(); ijet++){
 964     AliAODJet *tmp = dynamic_cast<AliAODJet*>(aodJets->At(ijet));
 965     if(!tmp) continue;
 966     if( tmp->Pt() < fJetPtCut ) continue;
 967     if( tmp->Eta() < fJetEtaMin || tmp->Eta() > fJetEtaMax)continue;
 968     if(fJetRejectType==kReject1Track && tmp->GetRefTracks()->GetEntriesFast()==1)continue;//rejecting 1track jet if...
 969     list->Add(tmp);
 970     nJets++;
 971   }//ij loop
 972   list->Sort();
 973   return nJets;
 974 }
 975 //_________________________________________________________________________________//
 976
 977 Int_t AliAnalysisTaskJetProperties::GetListOfJetTracks(TList* list, const AliAODJet* jet)
 978 {
 979   //this functionality is motivated by AliAnalysisTaskFragmentationFunction
 980   if(fDebug > 3) printf("AliAnalysisTaskJetProperties::GetListOfJetTracks() \n");
 981
 982   // list of jet tracks from trackrefs
 983   Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
 984   Int_t NoOfTracks=0;
 985   for (Int_t itrack=0; itrack<nTracks; itrack++) {
 986     if(fTrackType==kTrackUndef){
 987       if(fDebug>3)Printf("%s:%d unknown track type %d in AOD", (char*)__FILE__,__LINE__,kTrackUndef);
 988       return 0;
 989     }//if
 990     else if(fTrackType == kTrackAOD){
 991       AliAODTrack* trackAod = dynamic_cast<AliAODTrack*>(jet->GetRefTracks()->At(itrack));
 992       if(!trackAod){
 993         AliError("expected ref track not found ");
 994         continue;
 995       }
 996       list->Add(trackAod);
 997       NoOfTracks++;
 998     }//if
 999
1000     else if(fTrackType==kTrackAODMC){
1001       AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(jet->GetRefTracks()->At(itrack));
1002       if(!trackmc){
1003         AliError("expected ref trackmc not found ");
1004         continue;
1005       }
1006       list->Add(trackmc);
1007       NoOfTracks++;
1008     }//if
1009
1010     else if (fTrackType==kTrackKine){
1011       AliVParticle* trackkine = dynamic_cast<AliVParticle*>(jet->GetRefTracks()->At(itrack));
1012       if(!trackkine){
1013         AliError("expected ref trackkine not found ");
1014         continue;
1015       }
1016       list->Add(trackkine);
1017       NoOfTracks++;
1018     }//if
1019   }//itrack loop
1020   list->Sort();
1021   return NoOfTracks;
1022 }//initial loop
1023 //_________________________________________________________________________________//
1024
1025 void AliAnalysisTaskJetProperties::FillJetProperties(TList *jetlist){
1026   //filling up the histograms jets and tracks inside jet
1027   if(fDebug > 2) printf("AliAnalysisTaskJetProperties::FillJetProperties() \n");
1028
1029   for(Int_t iJet=0; iJet < jetlist->GetEntries(); iJet++){
1030     Float_t JetPt;Float_t JetPhi;Float_t JetEta;Float_t JetE;
1031     AliAODJet *jet = dynamic_cast<AliAODJet*>(jetlist->At(iJet));
1032     if(!jet)continue;
1033     JetEta = jet->Eta();
1034     JetPhi = jet->Phi();
1035     JetPt  = jet->Pt();
1036     JetE   = jet->E();
1037     fh2EtaJet ->Fill(JetPt,JetEta);
1038     fh2PhiJet ->Fill(JetPt,JetPhi);
1039     fh2PtJet  ->Fill(JetPt,JetPt);
1040     fh1PtJet  ->Fill(JetPt);
1041     fTrackListJet->Clear();
1042     Int_t nJT = GetListOfJetTracks(fTrackListJet,jet);
1043     Int_t nJetTracks = 0;
1044     if(nJT>=0) nJetTracks = fTrackListJet->GetEntries();
1045     if(fDebug>2){
1046       Printf("%s:%d Jet tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
1047       if(nJT != nJetTracks) Printf("%s:%d Mismatch Jet Tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
1048     }
1049
1050     fh2NtracksJet  ->Fill(JetPt,fTrackListJet->GetEntries());
1051     fProNtracksJet ->Fill(JetPt,fTrackListJet->GetEntries());
1052
1053     for (Int_t j =0; j< fTrackListJet->GetEntries(); j++){
1054       if(fTrackType==kTrackUndef)continue;
1055         Float_t TrackEta=-99.0; Float_t TrackPt=-99.0; Float_t TrackPhi=-99.0;
1056         Float_t FF=-99.0;       Float_t DelEta=-99.0;  Float_t DelPhi=-99.0;
1057         Float_t DelR=-99.0;     Float_t AreaJ=-99.0;   Float_t Ksi=-99.0;
1058         if(fTrackType==kTrackAOD){
1059           AliAODTrack *trackaod = dynamic_cast<AliAODTrack*>(fTrackListJet->At(j));
1060           if(!trackaod)continue;
1061           TrackEta = trackaod->Eta();
1062           TrackPhi = trackaod->Phi();
1063           TrackPt  = trackaod->Pt();
1064         }//if kTrackAOD
1065         else if(fTrackType==kTrackAODMC){
1066           AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(fTrackListJet->At(j));
1067           if(!trackmc)continue;
1068           TrackEta = trackmc->Eta();
1069           TrackPhi = trackmc->Phi();
1070           TrackPt  = trackmc->Pt();
1071         }//if kTrackAODMC
1072         else if(fTrackType==kTrackKine){
1073           AliVParticle* trackkine = dynamic_cast<AliVParticle*>(fTrackListJet->At(j));
1074           if(!trackkine)continue;
1075           TrackEta = trackkine->Eta();
1076           TrackPhi = trackkine->Phi();
1077           TrackPt  = trackkine->Pt();
1078         }//kTrackKine
1079         if(JetPt)FF = TrackPt/JetPt;
1080         if(FF)Ksi = TMath::Log(1./FF);
1081         DelEta      = TMath::Abs(JetEta - TrackEta);
1082         DelPhi      = TMath::Abs(JetPhi - TrackPhi);
1083         if(DelPhi>TMath::Pi())DelPhi = TMath::Abs(DelPhi-TMath::TwoPi());
1084         DelR        = TMath::Sqrt(DelEta*DelEta + DelPhi*DelPhi);
1085         AreaJ       = TMath::Pi()*DelR*DelR;
1086         fh2EtaTrack ->Fill(JetPt,TrackEta);
1087         fh2PhiTrack ->Fill(JetPt,TrackPhi);
1088         fh2PtTrack  ->Fill(JetPt,TrackPt);
1089         fh2FF       ->Fill(JetPt,FF);
1090         fh2Ksi      ->Fill(JetPt,Ksi);
1091         fh2DelEta   ->Fill(JetPt,DelEta);
1092         fh2DelPhi   ->Fill(JetPt,DelPhi);
1093         fh2DelR     ->Fill(JetPt,DelR);
1094     }//track loop
1095     fTrackListJet->Clear();
1096   }//iJet loop
1097 }//FillJetProperties
1098 //_________________________________________________________________________________//
1099 void AliAnalysisTaskJetProperties::FillFFCorr(TList *jetlist){
1100   //filling up the histograms jets and tracks inside jet
1101   if(fDebug > 2) printf("AliAnalysisTaskJetProperties::FillFFCorr() \n");
1102
1103   for(Int_t iJet=0; iJet < jetlist->GetEntries(); iJet++){
1104     Float_t JetPt;Float_t JetTheta;
1105     AliAODJet *jet = dynamic_cast<AliAODJet*>(jetlist->At(iJet));
1106     if(!jet)continue;
1107     JetTheta = jet->Theta();
1108     JetPt  = jet->Pt();
1109     fh1CorrJetPt -> Fill(JetPt);
1110
1111     fTrackListJet->Clear();
1112     Int_t nJT = GetListOfJetTracks(fTrackListJet,jet);
1113     Int_t nJetTracks = 0;
1114     if(nJT>=0) nJetTracks = fTrackListJet->GetEntries();
1115     if(fDebug>2){
1116       Printf("%s:%d Jet tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
1117       if(nJT != nJetTracks) Printf("%s:%d Mismatch Jet Tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
1118     }//fDebug
1119
1120     Float_t TotPt = 0.;
1121     for (Int_t j =0; j< fTrackListJet->GetEntries(); j++){
1122       if(fTrackType==kTrackUndef)continue;
1123       Float_t TrackPt1=-99.0;  Float_t TrackPt2=-99.0;
1124       Float_t FF1=-99.0;       Float_t FF2=-99.0;
1125       Float_t Ksi1=-99.0;      Float_t Ksi2=-99.0;
1126       Float_t TrackTheta1=0.0; Float_t TrackTheta2=0.0;
1127       Float_t DelTheta1=0.0;   Float_t DelTheta2=0.0;
1128       Float_t jT1=-99.0;       Float_t jT2=-99.0;
1129       if(fTrackType==kTrackAOD){
1130         AliAODTrack *trackaod = dynamic_cast<AliAODTrack*>(fTrackListJet->At(j));
1131         if(!trackaod)continue;
1132         TrackTheta1 = trackaod->Theta();
1133         TrackPt1    = trackaod->Pt();
1134       }//if kTrackAOD
1135       else if(fTrackType==kTrackAODMC){
1136         AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(fTrackListJet->At(j));
1137         if(!trackmc)continue;
1138         TrackTheta1 = trackmc->Theta();
1139         TrackPt1    = trackmc->Pt();
1140       }//if kTrackAODMC
1141       else if(fTrackType==kTrackKine){
1142         AliVParticle* trackkine = dynamic_cast<AliVParticle*>(fTrackListJet->At(j));
1143         if(!trackkine)continue;
1144         TrackTheta1 = trackkine->Theta();
1145         TrackPt1    = trackkine->Pt();
1146       }//kTrackKine
1147       TotPt += TrackPt1;
1148       if(JetPt)FF1 = TrackPt1/JetPt;
1149       if(FF1)Ksi1  = TMath::Log(1./FF1);
1150       DelTheta1    = TMath::Abs(JetTheta - TrackTheta1);
1151       jT1          = TrackPt1 * TMath::Sin(DelTheta1);
1152       fh2CorrPtTrack1  ->Fill(JetPt,TrackPt1);
1153       fh2CorrFF1       ->Fill(JetPt,FF1);
1154       fh2CorrKsi1      ->Fill(JetPt,Ksi1);
1155       fh2CorrjT1       ->Fill(JetPt,jT1);
1156       for (Int_t jj =j+1; jj< fTrackListJet->GetEntries(); jj++){
1157         if(fTrackType==kTrackUndef)continue;
1158         if(fTrackType==kTrackAOD){
1159           AliAODTrack *trackaod2 = dynamic_cast<AliAODTrack*>(fTrackListJet->At(jj));
1160           if(!trackaod2)continue;
1161           TrackTheta2 = trackaod2->Theta();
1162           TrackPt2    = trackaod2->Pt();
1163         }//if kTrackAOD
1164         else if(fTrackType==kTrackAODMC){
1165           AliAODMCParticle* trackmc2 = dynamic_cast<AliAODMCParticle*>(fTrackListJet->At(jj));
1166           if(!trackmc2)continue;
1167           TrackTheta2 = trackmc2->Theta();
1168           TrackPt2    = trackmc2->Pt();
1169         }//if kTrackAODMC
1170         else if(fTrackType==kTrackKine){
1171           AliVParticle* trackkine2 = dynamic_cast<AliVParticle*>(fTrackListJet->At(jj));
1172           if(!trackkine2)continue;
1173           TrackTheta2 = trackkine2->Theta();
1174           TrackPt2    = trackkine2->Pt();
1175         }//kTrackKine
1176         if(JetPt)FF2 = TrackPt2/JetPt;
1177         if(FF2)Ksi2  = TMath::Log(1./FF2);
1178         DelTheta2    = TMath::Abs(JetTheta - TrackTheta2);
1179         jT2          = TrackPt2 * TMath::Sin(DelTheta2);
1180         Float_t ptmin=10.; Float_t ptmax=20.0;
1181         for(Int_t iBin=0; iBin<6; iBin++){
1182           if(iBin == 0)ptmin=10.; ptmax=20.;
1183           if(iBin == 1)ptmin=20.; ptmax=30.;
1184           if(iBin == 2)ptmin=30.; ptmax=40.;
1185           if(iBin == 3)ptmin=40.; ptmax=60.;
1186           if(iBin == 4)ptmin=60.; ptmax=80.;
1187           if(iBin == 5)ptmin=80.; ptmax=100.;
1188           if(JetPt>ptmin && JetPt <= ptmax){
1189             fh2CorrPt1Pt2[iBin]  ->Fill(TrackPt1,TrackPt2);
1190             fh2CorrZ1Z2[iBin]    ->Fill(FF1,FF2);
1191             fh2CorrKsi1Ksi2[iBin]->Fill(Ksi1,Ksi2);
1192             fh2CorrjT1jT2[iBin]  ->Fill(jT1,jT2);
1193           }//if loop
1194         }//iBin loop
1195       }//inside track loop
1196     }//outer track loop
1197     Float_t diff_JetPt_TrkPtSum = JetPt - TotPt;
1198     fh1JetPtvsTrkSum->Fill(diff_JetPt_TrkPtSum);
1199     TotPt = 0.;
1200     fTrackListJet->Clear();
1201   }//iJet loop
1202 }//FillFFCorr
1203 //_________________________________________________________________________________//
1204
1205 void AliAnalysisTaskJetProperties::FillJetShape(TList *jetlist){
1206   //filling up the histograms
1207   if(fDebug > 2) printf("AliAnalysisTaskJetProperties::FillJetShape() \n");
1208   Float_t JetEta; Float_t JetPhi; Float_t JetPt;
1209   AliAODJet *jet = dynamic_cast<AliAODJet*>(jetlist->At(0));//Leading jet only
1210   if(jet){
1211     JetEta = jet->Eta();
1212     JetPhi = jet->Phi();
1213     JetPt  = jet->Pt();
1214     fh1PtLeadingJet->Fill(JetPt);
1215     Float_t NchSumA[50]         = {0.};
1216     Float_t PtSumA[50]          = {0.};
1217     Float_t delRPtSum80pc       = 0;
1218     Float_t delRNtrkSum80pc     = 0;
1219     Float_t PtSumDiffShape[10]  = {0.0};
1220     Float_t PtSumIntShape[10]   = {0.0};
1221     Int_t kNbinsR               = 10;
1222     fTrackListJet->Clear();
1223     Int_t nJT = GetListOfJetTracks(fTrackListJet,jet);
1224     Int_t nJetTracks = 0;
1225     if(nJT>=0) nJetTracks = fTrackListJet->GetEntries();
1226     if(fDebug>3){
1227       Printf("%s:%d Jet tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
1228       if(nJT != nJetTracks) Printf("%s:%d Mismatch Jet Tracks: %d %d",(char*)__FILE__,__LINE__,nJT,nJetTracks);
1229     }
1230     fh2NtracksLeadingJet->Fill(JetPt,nJetTracks);
1231     fProNtracksLeadingJet->Fill(JetPt,nJetTracks);
1232     Int_t   *index     = new Int_t   [nJetTracks];//dynamic array containing index
1233     Float_t *delRA     = new Float_t [nJetTracks];//dynamic array containing delR
1234     Float_t *delEtaA   = new Float_t [nJetTracks];//dynamic array containing delEta
1235     Float_t *delPhiA   = new Float_t [nJetTracks];//dynamic array containing delPhi
1236     Float_t *trackPtA  = new Float_t [nJetTracks];//dynamic array containing pt-track
1237     Float_t *trackEtaA = new Float_t [nJetTracks];//dynamic array containing eta-track
1238     Float_t *trackPhiA = new Float_t [nJetTracks];//dynamic array containing phi-track
1239     for(Int_t ii=0; ii<nJetTracks; ii++){
1240       index[ii]     = 0;
1241       delRA[ii]     = 0.;
1242       delEtaA[ii]   = 0.;
1243       delPhiA[ii]   = 0.;
1244       trackPtA[ii]  = 0.;
1245       trackEtaA[ii] = 0.;
1246       trackPhiA[ii] = 0.;
1247     }//ii loop
1248
1249     for (Int_t j =0; j< nJetTracks; j++){
1250       if(fTrackType==kTrackUndef)continue;
1251
1252       Float_t TrackEta=-99.0; Float_t TrackPt=-99.0; Float_t TrackPhi=-99.0;
1253       Float_t DelEta=-99.0;  Float_t DelPhi=-99.0; Float_t DelR=-99.0; Float_t AreaJ=-99.0;
1254
1255       if(fTrackType==kTrackAOD){
1256         AliAODTrack *trackaod = dynamic_cast<AliAODTrack*>(fTrackListJet->At(j));
1257         if(!trackaod)continue;
1258         TrackEta = trackaod->Eta();
1259         TrackPhi = trackaod->Phi();
1260         TrackPt  = trackaod->Pt();
1261       }//if kTrackAOD
1262       else if(fTrackType==kTrackAODMC){
1263         AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(fTrackListJet->At(j));
1264         if(!trackmc)continue;
1265         TrackEta = trackmc->Eta();
1266         TrackPhi = trackmc->Phi();
1267         TrackPt  = trackmc->Pt();
1268       }//if kTrackAODMC
1269       else if(fTrackType==kTrackKine){
1270         AliVParticle* trackkine = dynamic_cast<AliVParticle*>(fTrackListJet->At(j));
1271         if(!trackkine)continue;
1272         TrackEta = trackkine->Eta();
1273         TrackPhi = trackkine->Phi();
1274         TrackPt  = trackkine->Pt();
1275       }//if kTrackKine
1276
1277       DelEta = TMath::Abs(JetEta - TrackEta);
1278       DelPhi = TMath::Abs(JetPhi - TrackPhi);
1279       if(DelPhi>TMath::Pi())DelPhi = TMath::Abs(DelPhi-TMath::TwoPi());
1280       DelR   = TMath::Sqrt(DelEta*DelEta + DelPhi*DelPhi);
1281       AreaJ  = TMath::Pi()*DelR*DelR;
1282
1283       fh2AreaCh ->Fill(JetPt,AreaJ);
1284       fProAreaCh->Fill(JetPt,AreaJ);
1285       delRA[j]     = DelR;
1286       delEtaA[j]   = DelEta;
1287       delPhiA[j]   = DelPhi;
1288       trackPtA[j]  = TrackPt;
1289       trackEtaA[j] = TrackEta;
1290       trackPhiA[j] = TrackPhi;
1291
1292       //calculating diff and integrated jet shapes
1293       Float_t kDeltaR = 0.1;
1294       Float_t RMin    = kDeltaR/2.0;
1295       Float_t RMax    = kDeltaR/2.0;
1296       Float_t tmpR    = 0.05;
1297       for(Int_t ii1=0; ii1<kNbinsR;ii1++){
1298         if((DelR > (tmpR-RMin)) && (DelR <=(tmpR+RMax)))PtSumDiffShape[ii1]+= TrackPt;
1299         if(DelR>0.0 && DelR <=(tmpR+RMax))PtSumIntShape[ii1]+= TrackPt;
1300         tmpR += 0.1;
1301       }//ii1 loop
1302
1303       for(Int_t ibin=1; ibin<=50; ibin++){
1304         Float_t xlow = 0.02*(ibin-1);
1305         Float_t xup  = 0.02*ibin;
1306         if( xlow <= DelR && DelR < xup){
1307           NchSumA[ibin-1]++;
1308           PtSumA[ibin-1]+= TrackPt;
1309         }//if loop
1310       }//for ibin loop
1311     }//track loop
1312     fTrackListJet->Clear();
1313
1314     //---------------------
1315     Float_t tmp1R = 0.05;
1316     for(Int_t jj1=0; jj1<kNbinsR;jj1++){
1317       if(JetPt>20 && JetPt<=100){
1318         fProDiffJetShape->Fill(tmp1R,PtSumDiffShape[jj1]/JetPt);
1319         fProIntJetShape ->Fill(tmp1R,PtSumIntShape[jj1]/JetPt);
1320       }
1321       Float_t jetPtMin0=20.0; Float_t jetPtMax0=30.0;
1322       for(Int_t k=0; k<13; k++){
1323         if(k==0) {jetPtMin0=20.0;jetPtMax0=25.0;}
1324         if(k==1) {jetPtMin0=25.0;jetPtMax0=30.0;}
1325         if(k==2) {jetPtMin0=20.0;jetPtMax0=30.0;}
1326         if(k==3) {jetPtMin0=30.0;jetPtMax0=35.0;}
1327         if(k==4) {jetPtMin0=35.0;jetPtMax0=40.0;}
1328         if(k==5) {jetPtMin0=30.0;jetPtMax0=40.0;}
1329         if(k==6) {jetPtMin0=40.0;jetPtMax0=50.0;}
1330         if(k==7) {jetPtMin0=50.0;jetPtMax0=60.0;}
1331         if(k==8) {jetPtMin0=40.0;jetPtMax0=60.0;}
1332         if(k==9) {jetPtMin0=60.0;jetPtMax0=70.0;}
1333         if(k==10){jetPtMin0=70.0;jetPtMax0=80.0;}
1334         if(k==11){jetPtMin0=60.0;jetPtMax0=80.0;}
1335         if(k==12){jetPtMin0=80.0;jetPtMax0=100.0;}
1336         if(JetPt>jetPtMin0 && JetPt<=jetPtMax0){
1337           fProDiffJetShapeA[k]->Fill(tmp1R,PtSumDiffShape[jj1]/JetPt);
1338           fProIntJetShapeA[k] ->Fill(tmp1R,PtSumIntShape[jj1]/JetPt);
1339         }//if
1340       }//k loop
1341       tmp1R +=0.1;
1342     }//jj1 loop
1343     //----------------------//
1344     Float_t PtSum = 0;
1345     Int_t NtrkSum = 0;
1346     Bool_t iflagPtSum   = kFALSE;
1347     Bool_t iflagNtrkSum = kFALSE;
1348     TMath::Sort(nJetTracks,delRA,index,0);
1349     for(Int_t ii=0; ii<nJetTracks; ii++){
1350       NtrkSum ++;
1351       PtSum += trackPtA[index[ii]];
1352       /*
1353         cout << index[ii] << "\t" <<
1354         delR[ii] << "\t" <<
1355         delEta[ii]<< "\t" <<
1356         delPhi[ii]<< "\t" <<
1357         trackPt[ii]<< "\t" <<
1358         trackEta[ii]<< "\t" <<
1359         trackPhi[ii]<< "\t DelR " <<
1360         delR[index[ii]] << endl;
1361       */
1362       if(!iflagNtrkSum){
1363         if((Float_t)NtrkSum/(Float_t)nJetTracks > 0.79){
1364           delRNtrkSum80pc = delRA[index[ii]];
1365           iflagNtrkSum = kTRUE;
1366         }//if NtrkSum
1367       }//if iflag
1368       if(!iflagPtSum){
1369         if(PtSum/JetPt >= 0.8000){
1370           delRPtSum80pc = delRA[index[ii]];
1371           iflagPtSum = kTRUE;
1372         }//if PtSum
1373       }//if iflag
1374     }//track loop 2nd
1375     delete [] index;
1376     delete [] delRA;
1377     delete [] delEtaA;
1378     delete [] delPhiA;
1379     delete [] trackPtA;
1380     delete [] trackEtaA;
1381     delete [] trackPhiA;
1382     fh2DelR80pcNch ->Fill(JetPt,delRNtrkSum80pc);
1383     fProDelR80pcNch->Fill(JetPt,delRNtrkSum80pc);
1384     fh2DelR80pcPt  ->Fill(JetPt,delRPtSum80pc);
1385     fProDelR80pcPt ->Fill(JetPt,delRPtSum80pc);
1386     for(Int_t ibin=0; ibin<50; ibin++){
1387       Float_t iR = 0.02*ibin + 0.01;
1388       fh3PtDelRNchSum->Fill(JetPt,iR,NchSumA[ibin]);
1389       fh3PtDelRPtSum ->Fill(JetPt,iR,PtSumA[ibin]);
1390       Float_t jetPtMin=20.0; Float_t jetPtMax=30.0;
1391       for(Int_t k=0; k<13; k++){
1392         if(k==0) {jetPtMin=20.0;jetPtMax=25.0;}
1393         if(k==1) {jetPtMin=25.0;jetPtMax=30.0;}
1394         if(k==2) {jetPtMin=20.0;jetPtMax=30.0;}
1395         if(k==3) {jetPtMin=30.0;jetPtMax=35.0;}
1396         if(k==4) {jetPtMin=35.0;jetPtMax=40.0;}
1397         if(k==5) {jetPtMin=30.0;jetPtMax=40.0;}
1398         if(k==6) {jetPtMin=40.0;jetPtMax=50.0;}
1399         if(k==7) {jetPtMin=50.0;jetPtMax=60.0;}
1400         if(k==8) {jetPtMin=40.0;jetPtMax=60.0;}
1401         if(k==9) {jetPtMin=60.0;jetPtMax=70.0;}
1402         if(k==10){jetPtMin=70.0;jetPtMax=80.0;}
1403         if(k==11){jetPtMin=60.0;jetPtMax=80.0;}
1404         if(k==12){jetPtMin=80.0;jetPtMax=100.0;}
1405         if(JetPt>jetPtMin && JetPt<jetPtMax){
1406           fProDelRNchSum[k]->Fill(iR,NchSumA[ibin]);
1407           fProDelRPtSum[k]->Fill(iR,PtSumA[ibin]);
1408         }//if
1409       }//k loop
1410     }//ibin loop
1411   }//if(jet)
1412 }//FillJetShape()
1413 //_________________________________________________________________________________//
1414
1415 void AliAnalysisTaskJetProperties::FillJetShapeUE(TList *jetlist){
1416   //filling up the histograms
1417   if(fDebug > 2) printf("AliAnalysisTaskJetProperties::FillJetShape() \n");
1418   AliAODJet *jet = dynamic_cast<AliAODJet*>(jetlist->At(0));//Leading jet only
1419   if(jet){
1420     Double_t jetMom[3];
1421     jet->PxPyPz(jetMom);
1422
1423     TVector3 jet3mom(jetMom);
1424     // Rotate phi and keep eta unchanged
1425
1426     const Double_t alpha = TMath::Pi()/2.;
1427     Double_t etaTilted = jet3mom.Eta();
1428     Double_t phiTilted = TVector2::Phi_0_2pi(jet3mom.Phi()) + alpha;
1429     if(phiTilted > 2*TMath::Pi()) phiTilted = phiTilted - 2*TMath::Pi();
1430
1431     Double_t JetPt  = jet->Pt();
1432
1433     Float_t NchSumA[50]         = {0.};
1434     Float_t PtSumA[50]          = {0.};
1435     Float_t delRPtSum80pc       = 0;
1436     Float_t delRNtrkSum80pc     = 0;
1437     Float_t PtSumDiffShape[10]  = {0.0};
1438     Float_t PtSumIntShape[10]   = {0.0};
1439     Int_t kNbinsR               = 10;
1440
1441     //Int_t nTracks = GetListOfTracks(fTrackList,fTrackType);
1442     fTrackList->Clear();
1443     GetListOfTracks(fTrackList,fTrackType);
1444     Double_t sumPtPerp = 0;
1445     //if(fDebug > 100) printf("Cone radius for bckg. = %f Track type = %d\n",fJetRadius, fTrackType);
1446     fTrackListUE->Clear();
1447     GetTracksTiltedwrpJetAxis(alpha, fTrackList, fTrackListUE,jet,fJetRadius,sumPtPerp);
1448     fTrackList->Clear();
1449     fh1PtSumInJetConeUE->Fill(sumPtPerp);
1450     Int_t nTracksUE = fTrackListUE->GetEntries();
1451     fh2NtracksLeadingJetUE->Fill(JetPt,nTracksUE);
1452     fProNtracksLeadingJetUE->Fill(JetPt,nTracksUE);
1453
1454     Int_t   *index     = new Int_t   [nTracksUE];//dynamic array containing index
1455     Float_t *delRA     = new Float_t [nTracksUE];//dynamic array containing delR
1456     Float_t *delEtaA   = new Float_t [nTracksUE];//dynamic array containing delEta
1457     Float_t *delPhiA   = new Float_t [nTracksUE];//dynamic array containing delPhi
1458     Float_t *trackPtA  = new Float_t [nTracksUE];//dynamic array containing pt-track
1459     Float_t *trackEtaA = new Float_t [nTracksUE];//dynamic array containing eta-track
1460     Float_t *trackPhiA = new Float_t [nTracksUE];//dynamic array containing phi-track
1461     for(Int_t ii=0; ii<nTracksUE; ii++){
1462       index[ii]     = 0;
1463       delRA[ii]     = 0.;
1464       delEtaA[ii]   = 0.;
1465       delPhiA[ii]   = 0.;
1466       trackPtA[ii]  = 0.;
1467       trackEtaA[ii] = 0.;
1468       trackPhiA[ii] = 0.;
1469     }//ii loop
1470
1471     for (Int_t j =0; j< nTracksUE; j++){
1472       if(fTrackType==kTrackUndef)continue;
1473
1474       Float_t TrackEta=-99.0; Float_t TrackPt=-99.0; Float_t TrackPhi=-99.0;
1475       Float_t DelEta=-99.0;  Float_t DelPhi=-99.0; Float_t DelR=-99.0; Float_t AreaJ=-99.0;
1476
1477       if(fTrackType==kTrackAOD){
1478         AliAODTrack *trackaod = dynamic_cast<AliAODTrack*>(fTrackListUE->At(j));
1479         if(!trackaod)continue;
1480         TrackEta = trackaod->Eta();
1481         TrackPhi = trackaod->Phi();
1482         TrackPt  = trackaod->Pt();
1483         //if(fDebug > 100) printf("FillJetShapeUE itrack, trackPt %d,  %f\n",j,TrackPt);
1484       }//if kTrackAOD
1485       else if(fTrackType==kTrackAODMC){
1486         AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(fTrackListUE->At(j));
1487         if(!trackmc)continue;
1488         TrackEta = trackmc->Eta();
1489         TrackPhi = trackmc->Phi();
1490         TrackPt  = trackmc->Pt();
1491         //if(fDebug > 100) printf("FillJetShapeUE itrack, trackPt %d,  %f\n",j,TrackPt);
1492       }//if kTrackAODMC
1493       else if(fTrackType==kTrackKine){
1494         AliVParticle* trackkine = dynamic_cast<AliVParticle*>(fTrackListUE->At(j));
1495         if(!trackkine)continue;
1496         TrackEta = trackkine->Eta();
1497         TrackPhi = trackkine->Phi();
1498         TrackPt  = trackkine->Pt();
1499         //if(fDebug > 100) printf("FillJetShapeUE itrack, trackPt %d,  %f\n",j,TrackPt);
1500       }//if kTrackKine
1501
1502       DelEta = TMath::Abs(etaTilted - TrackEta);
1503       DelPhi = TMath::Abs(phiTilted - TrackPhi);
1504       if(DelPhi>TMath::Pi())DelPhi = TMath::Abs(DelPhi-TMath::TwoPi());
1505       DelR   = TMath::Sqrt(DelEta*DelEta + DelPhi*DelPhi);
1506       AreaJ  = TMath::Pi()*DelR*DelR;
1507
1508       fh2AreaChUE ->Fill(JetPt,AreaJ);
1509       fProAreaChUE->Fill(JetPt,AreaJ);
1510
1511       delRA[j]     = DelR;
1512       delEtaA[j]   = DelEta;
1513       delPhiA[j]   = DelPhi;
1514       trackPtA[j]  = TrackPt;
1515       trackEtaA[j] = TrackEta;
1516       trackPhiA[j] = TrackPhi;
1517
1518       //calculating diff and integrated jet shapes
1519       Float_t kDeltaR = 0.1;
1520       Float_t RMin    = kDeltaR/2.0;
1521       Float_t RMax    = kDeltaR/2.0;
1522       Float_t tmpR    = 0.05;
1523       for(Int_t ii1=0; ii1<kNbinsR;ii1++){
1524         if((DelR > (tmpR-RMin)) && (DelR <=(tmpR+RMax)))PtSumDiffShape[ii1]+= TrackPt;
1525         if(DelR>0.0 && DelR <=(tmpR+RMax))PtSumIntShape[ii1]+= TrackPt;
1526         tmpR += 0.1;
1527       }//ii1 loop
1528
1529       for(Int_t ibin=1; ibin<=50; ibin++){
1530         Float_t xlow = 0.02*(ibin-1);
1531         Float_t xup  = 0.02*ibin;
1532         if( xlow <= DelR && DelR < xup){
1533           NchSumA[ibin-1]++;
1534           PtSumA[ibin-1]+= TrackPt;
1535         }//if loop
1536       }//for ibin loop
1537     }//track loop
1538     fTrackListUE->Clear();
1539
1540     //---------------------
1541     Float_t tmp1R = 0.05;
1542     for(Int_t jj1=0; jj1<kNbinsR;jj1++){
1543       if(JetPt>20 && JetPt<=100){
1544         fProDiffJetShapeUE->Fill(tmp1R,PtSumDiffShape[jj1]/JetPt);
1545         fProIntJetShapeUE ->Fill(tmp1R,PtSumIntShape[jj1]/JetPt);
1546       }
1547       Float_t jetPtMin0=20.0; Float_t jetPtMax0=30.0;
1548       for(Int_t k=0; k<13; k++){
1549         if(k==0) {jetPtMin0=20.0;jetPtMax0=25.0;}
1550         if(k==1) {jetPtMin0=25.0;jetPtMax0=30.0;}
1551         if(k==2) {jetPtMin0=20.0;jetPtMax0=30.0;}
1552         if(k==3) {jetPtMin0=30.0;jetPtMax0=35.0;}
1553         if(k==4) {jetPtMin0=35.0;jetPtMax0=40.0;}
1554         if(k==5) {jetPtMin0=30.0;jetPtMax0=40.0;}
1555         if(k==6) {jetPtMin0=40.0;jetPtMax0=50.0;}
1556         if(k==7) {jetPtMin0=50.0;jetPtMax0=60.0;}
1557         if(k==8) {jetPtMin0=40.0;jetPtMax0=60.0;}
1558         if(k==9) {jetPtMin0=60.0;jetPtMax0=70.0;}
1559         if(k==10){jetPtMin0=70.0;jetPtMax0=80.0;}
1560         if(k==11){jetPtMin0=60.0;jetPtMax0=80.0;}
1561         if(k==12){jetPtMin0=80.0;jetPtMax0=100.0;}
1562         if(JetPt>jetPtMin0 && JetPt<=jetPtMax0){
1563           fProDiffJetShapeAUE[k]->Fill(tmp1R,PtSumDiffShape[jj1]/JetPt);
1564           fProIntJetShapeAUE[k] ->Fill(tmp1R,PtSumIntShape[jj1]/JetPt);
1565         }//if
1566       }//k loop
1567       tmp1R +=0.1;
1568     }//jj1 loop
1569     //----------------------//
1570     Float_t PtSum = 0;
1571     Int_t NtrkSum = 0;
1572     Bool_t iflagPtSum   = kFALSE;
1573     Bool_t iflagNtrkSum = kFALSE;
1574     TMath::Sort(nTracksUE,delRA,index,0);
1575     for(Int_t ii=0; ii<nTracksUE; ii++){
1576       NtrkSum ++;
1577       PtSum += trackPtA[index[ii]];
1578       /*
1579         cout << index[ii] << "\t" <<
1580         delR[ii] << "\t" <<
1581         delEta[ii]<< "\t" <<
1582         delPhi[ii]<< "\t" <<
1583         trackPt[ii]<< "\t" <<
1584         trackEta[ii]<< "\t" <<
1585         trackPhi[ii]<< "\t DelR " <<
1586         delR[index[ii]] << endl;
1587       */
1588       if(!iflagNtrkSum){
1589         if((Float_t)NtrkSum/(Float_t)nTracksUE > 0.79){
1590           delRNtrkSum80pc = delRA[index[ii]];
1591           iflagNtrkSum = kTRUE;
1592         }//if NtrkSum
1593       }//if iflag
1594       if(!iflagPtSum){
1595         if(PtSum/JetPt >= 0.8000){
1596           delRPtSum80pc = delRA[index[ii]];
1597           iflagPtSum = kTRUE;
1598         }//if PtSum
1599       }//if iflag
1600     }//track loop 2nd
1601     delete [] index;
1602     delete [] delRA;
1603     delete [] delEtaA;
1604     delete [] delPhiA;
1605     delete [] trackPtA;
1606     delete [] trackEtaA;
1607     delete [] trackPhiA;
1608     fh2DelR80pcNchUE ->Fill(JetPt,delRNtrkSum80pc);
1609     fProDelR80pcNchUE->Fill(JetPt,delRNtrkSum80pc);
1610     fh2DelR80pcPtUE  ->Fill(JetPt,delRPtSum80pc);
1611     fProDelR80pcPtUE ->Fill(JetPt,delRPtSum80pc);
1612     for(Int_t ibin=0; ibin<50; ibin++){
1613       Float_t iR = 0.02*ibin + 0.01;
1614       fh3PtDelRNchSumUE->Fill(JetPt,iR,NchSumA[ibin]);
1615       fh3PtDelRPtSumUE ->Fill(JetPt,iR,PtSumA[ibin]);
1616       Float_t jetPtMin=20.0; Float_t jetPtMax=30.0;
1617       for(Int_t k=0; k<13; k++){
1618         if(k==0) {jetPtMin=20.0;jetPtMax=25.0;}
1619         if(k==1) {jetPtMin=25.0;jetPtMax=30.0;}
1620         if(k==2) {jetPtMin=20.0;jetPtMax=30.0;}
1621         if(k==3) {jetPtMin=30.0;jetPtMax=35.0;}
1622         if(k==4) {jetPtMin=35.0;jetPtMax=40.0;}
1623         if(k==5) {jetPtMin=30.0;jetPtMax=40.0;}
1624         if(k==6) {jetPtMin=40.0;jetPtMax=50.0;}
1625         if(k==7) {jetPtMin=50.0;jetPtMax=60.0;}
1626         if(k==8) {jetPtMin=40.0;jetPtMax=60.0;}
1627         if(k==9) {jetPtMin=60.0;jetPtMax=70.0;}
1628         if(k==10){jetPtMin=70.0;jetPtMax=80.0;}
1629         if(k==11){jetPtMin=60.0;jetPtMax=80.0;}
1630         if(k==12){jetPtMin=80.0;jetPtMax=100.0;}
1631         if(JetPt>jetPtMin && JetPt<jetPtMax){
1632           fProDelRNchSumUE[k]->Fill(iR,NchSumA[ibin]);
1633           fProDelRPtSumUE[k]->Fill(iR,PtSumA[ibin]);
1634         }//if
1635       }//k loop
1636     }//ibin loop
1637   }//if(jet)
1638 }//FillJetShapeUE()
1639 //_________________________________________________________________________________//
1640
1641 void AliAnalysisTaskJetProperties::GetTracksTiltedwrpJetAxis(Float_t alpha, TList* inputlist, TList* outputlist, const AliAODJet* jet, Double_t radius,Double_t& sumPt){
1642   //This part  is inherited from AliAnalysisTaskFragmentationFunction.cxx
1643   // List of tracks in cone perpendicular to the jet azimuthal direction
1644   Double_t jetMom[3];
1645   jet->PxPyPz(jetMom);
1646
1647   TVector3 jet3mom(jetMom);
1648   // Rotate phi and keep eta unchanged
1649   Double_t etaTilted = jet3mom.Eta();//no change in eta
1650   Double_t phiTilted = TVector2::Phi_0_2pi(jet3mom.Phi()) + alpha;//rotate phi by alpha
1651   if(phiTilted > 2*TMath::Pi()) phiTilted = phiTilted - 2*TMath::Pi();
1652
1653   for (Int_t itrack=0; itrack<inputlist->GetSize(); itrack++){
1654     if(fTrackType==kTrackUndef)continue;
1655     Double_t trackMom[3];
1656
1657     if(fTrackType==kTrackAOD){
1658       AliAODTrack *trackaod = dynamic_cast<AliAODTrack*>(inputlist->At(itrack));
1659       if(!trackaod)continue;
1660       trackaod->PxPyPz(trackMom);
1661       TVector3 track3mom(trackMom);
1662
1663       Double_t deta = track3mom.Eta() - etaTilted;
1664       Double_t dphi = TMath::Abs(track3mom.Phi() - phiTilted);
1665       if (dphi > TMath::Pi()) dphi = 2. * TMath::Pi() - dphi;
1666       Double_t dR = TMath::Sqrt(deta * deta + dphi * dphi);
1667
1668       if(dR<=radius){
1669         outputlist->Add(trackaod);
1670         sumPt += trackaod->Pt();
1671         //if(fDebug > 100) printf("GetTracksTiltewrpJetAxis itrack, trackPt %d,  %f\n",itrack,track3mom.Pt());
1672       }//if dR<=radius
1673     }//if kTrackAOD
1674     else if(fTrackType==kTrackAODMC){
1675       AliAODMCParticle* trackmc = dynamic_cast<AliAODMCParticle*>(inputlist->At(itrack));
1676       if(!trackmc)continue;
1677       trackmc->PxPyPz(trackMom);
1678       TVector3 track3mom(trackMom);
1679
1680       Double_t deta = track3mom.Eta() - etaTilted;
1681       Double_t dphi = TMath::Abs(track3mom.Phi() - phiTilted);
1682       if (dphi > TMath::Pi()) dphi = 2. * TMath::Pi() - dphi;
1683       Double_t dR = TMath::Sqrt(deta * deta + dphi * dphi);
1684
1685       if(dR<=radius){
1686         outputlist->Add(trackmc);
1687         sumPt += trackmc->Pt();
1688         //if(fDebug > 100) printf("GetTracksTiltewrpJetAxis itrack, trackPt %d,  %f\n",itrack,track3mom.Pt());
1689       }//if dR<=radius
1690     }//if kTrackAODMC
1691     else if(fTrackType==kTrackKine){
1692       AliVParticle* trackkine = dynamic_cast<AliVParticle*>(inputlist->At(itrack));
1693       if(!trackkine)continue;
1694       trackkine->PxPyPz(trackMom);
1695       TVector3 track3mom(trackMom);
1696
1697       Double_t deta = track3mom.Eta() - etaTilted;
1698       Double_t dphi = TMath::Abs(track3mom.Phi() - phiTilted);
1699       if (dphi > TMath::Pi()) dphi = 2. * TMath::Pi() - dphi;
1700       Double_t dR = TMath::Sqrt(deta * deta + dphi * dphi);
1701
1702       if(dR<=radius){
1703         outputlist->Add(trackkine);
1704         sumPt += trackkine->Pt();
1705         //if(fDebug > 100) printf("GetTracksTiltewrpJetAxis itrack, trackPt %d,  %f\n",itrack,track3mom.Pt());
1706       }//if dR<=radius
1707     }//kTrackKine
1708   }//itrack
1709 }//initial loop
1710 //-----------------------------------------------//
1711 Int_t AliAnalysisTaskJetProperties::GetListOfTracks(TList *list, Int_t type)
1712 {
1713   //This part  is inherited from AliAnalysisTaskFragmentationFunction.cxx (and modified)
1714   // fill list of tracks selected according to type
1715
1716   if(fDebug > 3) Printf("%s:%d Selecting tracks with %d", (char*)__FILE__,__LINE__,type);
1717
1718   if(!list){
1719     if(fDebug>3) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
1720     return -1;
1721   }
1722
1723   if(!fAOD) return -1;
1724
1725   if(!fAOD->GetTracks()) return 0;
1726
1727   if(type==kTrackUndef) return 0;
1728
1729   Int_t iCount = 0;
1730   if(type==kTrackAOD){
1731     // all rec. tracks, esd filter mask, within acceptance
1732     for(Int_t it=0; it<fAOD->GetNumberOfTracks(); ++it){
1733       AliAODTrack *tr = fAOD->GetTrack(it);
1734       if(!tr)continue;
1735       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;//selecting filtermask
1736       if(tr->Eta() < fTrackEtaMin || tr->Eta() > fTrackEtaMax) continue;
1737       if(tr->Pt()  < fTrackPtCut) continue;
1738       list->Add(tr);
1739       iCount++;
1740     }//track loop
1741   }//if type
1742   else if (type==kTrackKine){
1743     // kine particles, primaries, charged only within acceptance
1744     if(!fMCEvent) return iCount;
1745     for(Int_t it=0; it<fMCEvent->GetNumberOfTracks(); ++it){
1746       if(!fMCEvent->IsPhysicalPrimary(it))continue;//selecting only primaries
1747       AliMCParticle* part = (AliMCParticle*) fMCEvent->GetTrack(it);
1748       if(!part)continue;
1749       if(part->Particle()->GetPDG()->Charge()==0)continue;//removing neutrals
1750       if(part->Eta() < fTrackEtaMin || part->Eta() > fTrackEtaMax) continue;//eta cut
1751       if(part->Pt()  < fTrackPtCut)continue;//pt cut
1752       list->Add(part);
1753       iCount++;
1754     }//track loop
1755   }//if type
1756   else if (type==kTrackAODMC) {
1757     // MC particles (from AOD), physical primaries, charged only within acceptance
1758     if(!fAOD) return -1;
1759
1760     TClonesArray *tca = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
1761     if(!tca)return iCount;
1762
1763     for(int it=0; it<tca->GetEntriesFast(); ++it){
1764       AliAODMCParticle *part = dynamic_cast<AliAODMCParticle*>(tca->At(it));
1765       if(!part)continue;
1766       if(!part->IsPhysicalPrimary())continue;
1767       if(part->Charge()==0) continue;
1768       if(part->Eta() > fTrackEtaMax || part->Eta() < fTrackEtaMin)continue;
1769       if(part->Pt()  < fTrackPtCut) continue;
1770       list->Add(part);
1771       iCount++;
1772     }//track loop
1773   }//if type
1774
1775   list->Sort();
1776   return iCount;
1777 }
1778 //_______________________________________________________________________________