PWGJE/AliAnalysisTaskAj.cxx

   1
   2 // ******************************************
   3 // This task searches for events with large dijet imbalance
   4 // and then looks to the jet structure of the b-t-b jets.
   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
  24 #include "TChain.h"
  25 #include "TTree.h"
  26 #include "TMath.h"
  27 #include "TH1F.h"
  28 #include "TH2F.h"
  29 #include "TH3F.h"
  30 #include "THnSparse.h"
  31 #include "TMatrixD.h"
  32 #include "TMatrixDSym.h"
  33 #include "TMatrixDSymEigen.h"
  34
  35 #include "TCanvas.h"
  36
  37 #include "AliLog.h"
  38
  39 #include "AliAnalysisTask.h"
  40 #include "AliAnalysisManager.h"
  41
  42 #include "AliVEvent.h"
  43 #include "AliESDEvent.h"
  44 #include "AliESDInputHandler.h"
  45 #include "AliCentrality.h"
  46 #include "AliAnalysisHelperJetTasks.h"
  47 #include "AliInputEventHandler.h"
  48 #include "AliAODJetEventBackground.h"
  49 #include "AliAnalysisTaskFastEmbedding.h"
  50 #include "AliAODEvent.h"
  51 #include "AliAODHandler.h"
  52 #include "AliAODJet.h"
  53
  54 #include "AliAnalysisTaskAj.h"
  55
  56 using std::cout;
  57 using std::endl;
  58
  59 ClassImp(AliAnalysisTaskAj)
  60
  61 AliAnalysisTaskAj::AliAnalysisTaskAj() :
  62 AliAnalysisTaskSE(),
  63 fESD(0x0),
  64 fAODIn(0x0),
  65 fAODOut(0x0),
  66 fAODExtension(0x0),
  67 fBackgroundBranch(""),
  68 fNonStdFile(""),
  69 fIsPbPb(kTRUE),
  70 fOfflineTrgMask(AliVEvent::kAny),
  71 fMinContribVtx(1),
  72 fVtxZMin(-10.),
  73 fVtxZMax(10.),
  74 fEvtClassMin(0),
  75 fEvtClassMax(4),
  76 fFilterMask(0),
  77 fRadioFrac(0.2),
  78 fMinDist(0.1),
  79 fCentMin(0.),
  80 fCentMax(100.),
  81 fNInputTracksMin(0),
  82 fNInputTracksMax(-1),
  83 fAngStructCloseTracks(0),
  84 fCheckMethods(0),
  85 fDoEventMixing(0),
  86 fJetEtaMin(-.5),
  87 fJetEtaMax(.5),
  88 fNevents(0x0),
  89 fTindex(0x0),
  90 fJetPtMin(20.),
  91 fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
  92 fJetPtFractionMin(0.5),
  93 fNMatchJets(4),
  94 fMatchMaxDist(0.8),
  95 fKeepJets(kFALSE),
  96 fkNbranches(2),
  97 fkEvtClasses(12),
  98 fOutputList(0x0),
  99 fbEvent(kTRUE),
 100 fHistEvtSelection(0x0),
 101 fh2Pt1Pt2C10(0x0),
 102 fh2Pt1Pt2C40(0x0),
 103 fh3LocalCoordinates(0x0),
 104 fh2Sum2pt20(0x0),
 105 fh2Sum4pt20(0x0),
 106 fh2Sum6pt20(0x0),
 107 fh2Sum8pt20(0x0),
 108 fh2Sum12pt20(0x0),
 109 fh2Sum2lpt20(0x0),
 110 fh2Sum4lpt20(0x0),
 111 fh2Sum6lpt20(0x0),
 112 fh2Sum8lpt20(0x0),
 113 fh2Sum12lpt20(0x0),
 114 fh2Sum2pt40(0x0),
 115 fh2Sum4pt40(0x0),
 116 fh2Sum6pt40(0x0),
 117 fh2Sum8pt40(0x0),
 118 fh2Sum12pt40(0x0),
 119 fh2Sum2lpt40(0x0),
 120 fh2Sum4lpt40(0x0),
 121 fh2Sum6lpt40(0x0),
 122 fh2Sum8lpt40(0x0),
 123 fh2Sum12lpt40(0x0),
 124 fhnDeltaR(0x0)
 125  {
 126    // default Constructor
 127
 128
 129    fJetBranchName[0] = "";
 130    fJetBranchName[1] = "";
 131
 132    fListJets[0] = new TList;
 133    fListJets[1] = new TList;
 134 }
 135
 136 AliAnalysisTaskAj::AliAnalysisTaskAj(const char *name) :
 137 AliAnalysisTaskSE(name),
 138 fESD(0x0),
 139 fAODIn(0x0),
 140 fAODOut(0x0),
 141 fAODExtension(0x0),
 142 fBackgroundBranch(""),
 143 fNonStdFile(""),
 144 fIsPbPb(kTRUE),
 145 fOfflineTrgMask(AliVEvent::kAny),
 146 fMinContribVtx(1),
 147 fVtxZMin(-10.),
 148 fVtxZMax(10.),
 149 fEvtClassMin(0),
 150 fEvtClassMax(4),
 151 fFilterMask(0),
 152 fRadioFrac(0.2),
 153 fMinDist(0.1),
 154 fCentMin(0.),
 155 fCentMax(100.),
 156 fNInputTracksMin(0),
 157 fNInputTracksMax(-1),
 158 fAngStructCloseTracks(0),
 159 fCheckMethods(0),
 160 fDoEventMixing(0),
 161 fJetEtaMin(-.5),
 162 fJetEtaMax(.5),
 163 fNevents(0x0),
 164 fTindex(0x0),
 165 fJetPtMin(20.),
 166 fJetTriggerExcludeMask(AliAODJet::kHighTrackPtTriggered),
 167 fJetPtFractionMin(0.5),
 168 fNMatchJets(4),
 169 fMatchMaxDist(0.8),
 170 fKeepJets(kFALSE),
 171 fkNbranches(2),
 172 fkEvtClasses(12),
 173 fOutputList(0x0),
 174 fbEvent(kTRUE),
 175 fHistEvtSelection(0x0),
 176 fh2Pt1Pt2C10(0x0),
 177 fh2Pt1Pt2C40(0x0),
 178 fh3LocalCoordinates(0x0),
 179 fh2Sum2pt20(0x0),
 180 fh2Sum4pt20(0x0),
 181 fh2Sum6pt20(0x0),
 182 fh2Sum8pt20(0x0),
 183 fh2Sum12pt20(0x0),
 184 fh2Sum2lpt20(0x0),
 185 fh2Sum4lpt20(0x0),
 186 fh2Sum6lpt20(0x0),
 187 fh2Sum8lpt20(0x0),
 188 fh2Sum12lpt20(0x0),
 189 fh2Sum2pt40(0x0),
 190 fh2Sum4pt40(0x0),
 191 fh2Sum6pt40(0x0),
 192 fh2Sum8pt40(0x0),
 193 fh2Sum12pt40(0x0),
 194 fh2Sum2lpt40(0x0),
 195 fh2Sum4lpt40(0x0),
 196 fh2Sum6lpt40(0x0),
 197 fh2Sum8lpt40(0x0),
 198 fh2Sum12lpt40(0x0),
 199 fhnDeltaR(0x0)
 200  {
 201    // Constructor
 202
 203
 204
 205
 206
 207
 208
 209    fJetBranchName[0] = "";
 210    fJetBranchName[1] = "";
 211
 212    fListJets[0] = new TList;
 213    fListJets[1] = new TList;
 214
 215    DefineOutput(1, TList::Class());
 216 }
 217
 218 AliAnalysisTaskAj::~AliAnalysisTaskAj()
 219 {
 220    delete fListJets[0];
 221    delete fListJets[1];
 222 }
 223
 224 void AliAnalysisTaskAj::SetBranchNames(const TString &branch1, const TString &branch2)
 225 {
 226    fJetBranchName[0] = branch1;
 227    fJetBranchName[1] = branch2;
 228 }
 229
 230 void AliAnalysisTaskAj::Init()
 231 {
 232
 233    // check for jet branches
 234    if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
 235       AliError("Jet branch name not set.");
 236    }
 237
 238 }
 239
 240 void AliAnalysisTaskAj::UserCreateOutputObjects()
 241 {
 242    // Create histograms
 243    // Called once
 244    OpenFile(1);
 245    if(!fOutputList) fOutputList = new TList;
 246    fOutputList->SetOwner(kTRUE);
 247
 248    Bool_t oldStatus = TH1::AddDirectoryStatus();
 249    TH1::AddDirectory(kFALSE);
 250
 251
 252    fHistEvtSelection = new TH1I("fHistEvtSelection", "event selection", 6, -0.5, 5.5);
 253    fHistEvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
 254    fHistEvtSelection->GetXaxis()->SetBinLabel(2,"events IN");
 255    fHistEvtSelection->GetXaxis()->SetBinLabel(3,"event selection (rejected)");
 256    fHistEvtSelection->GetXaxis()->SetBinLabel(4,"vertex cut (rejected)");
 257    fHistEvtSelection->GetXaxis()->SetBinLabel(5,"centrality (rejected)");
 258    fHistEvtSelection->GetXaxis()->SetBinLabel(6,"multiplicity (rejected)");
 259
 260      UInt_t entries = 0; // bit coded, see GetDimParams() below
 261      entries = 1<<0 | 1<<1 | 1<<2 | 1<<3 | 1<<4 | 1<<5 | 1<<6 |1<<7|1<<8;
 262      fhnDeltaR = NewTHnSparseF("fhnDeltaR", entries);
 263
 264
 265      //change binning in centrality
 266      Double_t *xPt0 = new Double_t[4];
 267      xPt0[0] = 0.;
 268      for(int i = 1; i<=3;i++){
 269       if(xPt0[i-1]<10)xPt0[i] = xPt0[i-1] + 10; // 1 - 5
 270       else if(xPt0[i-1]<40)xPt0[i] = xPt0[i-1] + 30; // 5 - 12
 271       else xPt0[i] = xPt0[i-1] + 60.; // 18
 272      }
 273      fhnDeltaR->SetBinEdges(0,xPt0);
 274      delete [] xPt0;
 275
 276
 277
 278      //change binning in jTtrack
 279      Double_t *xPt3 = new Double_t[3];
 280      xPt3[0] = 0.;
 281      for(int i = 1; i<=2;i++){
 282       if(xPt3[i-1]<1)xPt3[i] = xPt3[i-1] + 1.; // 1 - 5
 283       else xPt3[i] = xPt3[i-1] + 149.; // 18
 284      }
 285      fhnDeltaR->SetBinEdges(3,xPt3);
 286      delete [] xPt3;
 287
 288
 289
 290      //change binning in pTtrack
 291      Double_t *xPt4 = new Double_t[5];
 292      xPt4[0] = 0.;
 293      for(int i = 1; i<=4;i++){
 294       if(xPt4[i-1]<0.4)xPt4[i] = xPt4[i-1] + 0.4; // 1 - 5
 295       else if(xPt4[i-1]<3)xPt4[i] = xPt4[i-1] + 2.6; // 5 - 12
 296       else if(xPt4[i-1]<10)xPt4[i] = xPt4[i-1] + 7.4; // 5 - 12
 297       else xPt4[i] = xPt4[i-1] + 150.; // 18
 298      }
 299     fhnDeltaR->SetBinEdges(4,xPt4);
 300     delete [] xPt4;
 301
 302
 303
 304
 305
 306
 307     //change binning in HTI
 308      Double_t *xPt5 = new Double_t[4];
 309      xPt5[0] = 0.;
 310      for(int i = 1; i<=3;i++){
 311       if(xPt5[i-1]<10)xPt5[i] = xPt5[i-1] + 5; // 10 - 12
 312       else xPt5[i] = xPt5[i-1] + 40.; // 13
 313      }
 314     fhnDeltaR->SetBinEdges(8,xPt5);
 315     delete [] xPt5;
 316
 317
 318
 319     fh2Pt1Pt2C10 = new TH2F("Dijet spectra central","",20,0.,200.,20,0.,200.);
 320     fh2Pt1Pt2C40 = new TH2F("Dijet spectra peripheral","",20,0.,200.,20,0.,200.);
 321     fh3LocalCoordinates = new TH3F("Local coordinates","",10,-2,2,10,-2,2,10,0,100);
 322     fh2Sum2pt20 = new TH2F("pL R<0.2 pt20","",10,0.,1.,100,0.,100.);
 323     fh2Sum4pt20 = new TH2F("pL R<0.4 pt20","",10,0.,1.,100,0.,100.);
 324     fh2Sum6pt20 = new TH2F("pL R<0.6 pt20","",10,0.,1.,100,0.,100.);
 325     fh2Sum8pt20 = new TH2F("pL R<0.8 pt20","",10,0.,1.,100,0.,100.);
 326     fh2Sum12pt20 = new TH2F("pL R<1.2 pt20","",10,0.,1.,100,0.,100.);
 327     fh2Sum2lpt20 = new TH2F("pL R<0.2 low pt pt20","",10,0.,1.,100,0,100);
 328     fh2Sum4lpt20 = new TH2F("pL R<0.4 low pt pt20","",10,0.,1.,100,0,100);
 329     fh2Sum6lpt20 = new TH2F("pL R<0.6 low pt pt20","",10,0.,1.,100,0,100);
 330     fh2Sum8lpt20 = new TH2F("pL R<0.8 low pt pt20","",10,0.,1.,100,0,100);
 331     fh2Sum12lpt20 = new TH2F("pL R<1.2 low pt pt20","",10,0.,1.,100,0,100);
 332     fh2Sum2pt40 = new TH2F("pL R<0.2 pt40","",10,0.,1.,100,0.,100.);
 333     fh2Sum4pt40 = new TH2F("pL R<0.4 pt40","",10,0.,1.,100,0.,100.);
 334     fh2Sum6pt40 = new TH2F("pL R<0.6 pt40","",10,0.,1.,100,0.,100.);
 335     fh2Sum8pt40 = new TH2F("pL R<0.8 pt40","",10,0.,1.,100,0.,100.);
 336     fh2Sum12pt40 = new TH2F("pL R<1.2 pt40","",10,0.,1.,100,0.,100.);
 337     fh2Sum2lpt40 = new TH2F("pL R<0.2 low pt pt40","",10,0.,1.,100,0,100);
 338     fh2Sum4lpt40 = new TH2F("pL R<0.4 low pt pt40","",10,0.,1.,100,0,100);
 339     fh2Sum6lpt40 = new TH2F("pL R<0.6 low pt pt40","",10,0.,1.,100,0,100);
 340     fh2Sum8lpt40 = new TH2F("pL R<0.8 low pt pt40","",10,0.,1.,100,0,100);
 341     fh2Sum12lpt40 = new TH2F("pL R<1.2 low pt pt40","",10,0.,1.,100,0,100);
 342
 343
 344    fOutputList->Add(fHistEvtSelection);
 345    fOutputList->Add(fh2Pt1Pt2C10);
 346    fOutputList->Add(fh2Pt1Pt2C40);
 347    fOutputList->Add(fh3LocalCoordinates);
 348
 349    fOutputList->Add(fh2Sum2pt20);
 350    fOutputList->Add(fh2Sum4pt20);
 351    fOutputList->Add(fh2Sum6pt20);
 352    fOutputList->Add(fh2Sum8pt20);
 353    fOutputList->Add(fh2Sum12pt20);
 354    fOutputList->Add(fh2Sum2lpt20);
 355    fOutputList->Add(fh2Sum4lpt20);
 356    fOutputList->Add(fh2Sum6lpt20);
 357    fOutputList->Add(fh2Sum8lpt20);
 358    fOutputList->Add(fh2Sum12lpt20);
 359
 360    fOutputList->Add(fh2Sum2pt40);
 361    fOutputList->Add(fh2Sum4pt40);
 362    fOutputList->Add(fh2Sum6pt40);
 363    fOutputList->Add(fh2Sum8pt40);
 364    fOutputList->Add(fh2Sum12pt40);
 365    fOutputList->Add(fh2Sum2lpt40);
 366    fOutputList->Add(fh2Sum4lpt40);
 367    fOutputList->Add(fh2Sum6lpt40);
 368    fOutputList->Add(fh2Sum8lpt40);
 369    fOutputList->Add(fh2Sum12lpt40);
 370
 371    fOutputList->Add(fhnDeltaR);
 372
 373
 374
 375    // fOutputList->Add(fh3specbiased);
 376
 377    // =========== Switch on Sumw2 for all histos ===========
 378    for (Int_t i=0; i<fOutputList->GetEntries(); ++i) {
 379       TH1 *h1 = dynamic_cast<TH1*>(fOutputList->At(i));
 380       if (h1){
 381          h1->Sumw2();
 382          continue;
 383       }
 384     THnSparse *hn = dynamic_cast<THnSparse*>(fOutputList->At(i));
 385       if (hn){
 386          hn->Sumw2();
 387       }
 388    }
 389    TH1::AddDirectory(oldStatus);
 390
 391    PostData(1, fOutputList);
 392 }
 393
 394 void AliAnalysisTaskAj::UserExec(Option_t *)
 395 {
 396
 397
 398    if(!strlen(fJetBranchName[0].Data()) || !strlen(fJetBranchName[1].Data())){
 399       AliError("Jet branch name not set.");
 400       return;
 401    }
 402
 403    fESD=dynamic_cast<AliESDEvent*>(InputEvent());
 404    if (!fESD) {
 405       AliError("ESD not available");
 406       fAODIn = dynamic_cast<AliAODEvent*>(InputEvent());}
 407       fAODOut = dynamic_cast<AliAODEvent*>(AODEvent());
 408
 409       static AliAODEvent* aod = 0;
 410        // take all other information from the aod we take the tracks from
 411        if(!aod){
 412        if(!fESD)aod = fAODIn;
 413        else aod = fAODOut;}
 414
 415
 416
 417     if(fNonStdFile.Length()!=0){
 418     // case that we have an AOD extension we need can fetch the jets from the extended output
 419     AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
 420     fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);
 421     if(!fAODExtension){
 422       if(fDebug>1)Printf("AODExtension found for %s",fNonStdFile.Data());
 423     }}
 424
 425
 426
 427
 428
 429    // -- event selection --
 430    fHistEvtSelection->Fill(1); // number of events before event selection
 431
 432    // physics selection
 433    AliInputEventHandler* inputHandler = (AliInputEventHandler*)
 434    ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
 435    cout<<inputHandler->IsEventSelected()<<" "<<fOfflineTrgMask<<endl;
 436    if(!(inputHandler->IsEventSelected() & fOfflineTrgMask)){
 437       if(fDebug) Printf(" Trigger Selection: event REJECTED ... ");
 438       fHistEvtSelection->Fill(2);
 439       PostData(1, fOutputList);
 440       return;
 441    }
 442
 443    // vertex selection
 444    if(!aod){
 445      if(fDebug) Printf("%s:%d No AOD",(char*)__FILE__,__LINE__);
 446      fHistEvtSelection->Fill(3);
 447       PostData(1, fOutputList);
 448    }
 449    AliAODVertex* primVtx = aod->GetPrimaryVertex();
 450
 451    if(!primVtx){
 452      if(fDebug) Printf("%s:%d No primVtx",(char*)__FILE__,__LINE__);
 453      fHistEvtSelection->Fill(3);
 454      PostData(1, fOutputList);
 455      return;
 456    }
 457
 458    Int_t nTracksPrim = primVtx->GetNContributors();
 459    if ((nTracksPrim < fMinContribVtx) ||
 460          (primVtx->GetZ() < fVtxZMin) ||
 461          (primVtx->GetZ() > fVtxZMax) ){
 462       if(fDebug) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ());
 463       fHistEvtSelection->Fill(3);
 464       PostData(1, fOutputList);
 465       return;
 466    }
 467
 468    // event class selection (from jet helper task)
 469    Int_t eventClass = AliAnalysisHelperJetTasks::EventClass();
 470    if(fDebug) Printf("Event class %d", eventClass);
 471    if (eventClass < fEvtClassMin || eventClass > fEvtClassMax){
 472       fHistEvtSelection->Fill(4);
 473       PostData(1, fOutputList);
 474       return;
 475    }
 476
 477    // centrality selection
 478    AliCentrality *cent = 0x0;
 479    Double_t centValue = 0.;
 480    if(fESD) {cent = fESD->GetCentrality();
 481      if(cent) centValue = cent->GetCentralityPercentile("V0M");}
 482    else     centValue=((AliVAODHeader*)aod->GetHeader())->GetCentrality();
 483
 484    if(fDebug) printf("centrality: %f\n", centValue);
 485    //  if (centValue < fCentMin || centValue > fCentMax){
 486    //  fHistEvtSelection->Fill(4);
 487    //  PostData(1, fOutputList);
 488    //  return;
 489    // }
 490
 491
 492    fHistEvtSelection->Fill(0);
 493    // accepted events
 494    // -- end event selection --
 495
 496    // get background
 497    AliAODJetEventBackground* externalBackground = 0;
 498    if(fAODOut&&!externalBackground&&fBackgroundBranch.Length()){
 499       externalBackground =  (AliAODJetEventBackground*)(fAODOut->FindListObject(fBackgroundBranch.Data()));
 500       if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
 501    }
 502    if(fAODExtension&&!externalBackground&&fBackgroundBranch.Length()){
 503      externalBackground =  (AliAODJetEventBackground*)(fAODExtension->GetAOD()->FindListObject(fBackgroundBranch.Data()));
 504       if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
 505    }
 506
 507     if(fAODIn&&!externalBackground&&fBackgroundBranch.Length()){
 508       externalBackground =  (AliAODJetEventBackground*)(fAODIn->FindListObject(fBackgroundBranch.Data()));
 509       if(!externalBackground)Printf("%s:%d Background branch not found %s",(char*)__FILE__,__LINE__,fBackgroundBranch.Data());;
 510     }
 511   Float_t rho = 0;
 512    if(externalBackground)rho = externalBackground->GetBackground(0);
 513
 514
 515    // fetch jets
 516    TClonesArray *aodJets[2];
 517    aodJets[0]=0;
 518    if(fAODOut&&!aodJets[0]){
 519    aodJets[0] = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fJetBranchName[0].Data()));
 520    aodJets[1] = dynamic_cast<TClonesArray*>(fAODOut->FindListObject(fJetBranchName[1].Data()));  }
 521    if(fAODExtension && !aodJets[0]){
 522    aodJets[0] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[0].Data()));
 523    aodJets[1] = dynamic_cast<TClonesArray*>(fAODExtension->GetAOD()->FindListObject(fJetBranchName[1].Data()));  }
 524    if(fAODIn&&!aodJets[0]){
 525    aodJets[0] = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fJetBranchName[0].Data()));
 526    aodJets[1] = dynamic_cast<TClonesArray*>(fAODIn->FindListObject(fJetBranchName[1].Data())); }
 527
 528    //Double_t ptsub[aodJets[0]->GetEntriesFast()];
 529    //Int_t inord[aodJets[0]->GetEntriesFast()];
 530    //for(Int_t n=0;n<aodJets[0]->GetEntriesFast();n++){
 531    //  ptsub[n]=0;
 532    //  inord[n]=0;}
 533
 534    TList ParticleList;
 535    Int_t nT = GetListOfTracks(&ParticleList);
 536      for (Int_t iJetType = 0; iJetType < 2; iJetType++) {
 537       fListJets[iJetType]->Clear();
 538       if (!aodJets[iJetType]) continue;
 539
 540       if(fDebug) Printf("%s: %d jets",fJetBranchName[iJetType].Data(),aodJets[iJetType]->GetEntriesFast());
 541
 542
 543       for (Int_t iJet = 0; iJet < aodJets[iJetType]->GetEntriesFast(); iJet++) {
 544          AliAODJet *jet = dynamic_cast<AliAODJet*>((*aodJets[iJetType])[iJet]);
 545          if (jet) fListJets[iJetType]->Add(jet);
 546          // if(iJetType==0){
 547          // ptsub[iJet]=jet->Pt()-rho*jet->EffectiveAreaCharged();}
 548       }}
 549
 550
 551    Double_t eta1=0;
 552    Int_t selec=-1;
 553    Double_t ptmax=-10;
 554    Double_t areaj=0;
 555    //   Double_t phij=0;
 556    Double_t etaj=0;
 557    Double_t ptj=0;
 558    Double_t ptcorrj=0;
 559    for(Int_t i=0; i<fListJets[0]->GetEntries(); ++i){
 560            AliAODJet* jetj = (AliAODJet*)(fListJets[0]->At(i));
 561            etaj  = jetj->Eta();
 562            //           phij  = jetj->Phi();
 563            ptj = jetj->Pt();
 564            if(ptj==0) continue;
 565            areaj = jetj->EffectiveAreaCharged();
 566            ptcorrj=ptj-rho*areaj;
 567            if(ptcorrj<=0) continue;
 568            if((etaj<fJetEtaMin)||(eta1>fJetEtaMax)) continue;
 569            if(ptcorrj>ptmax){ptmax=ptcorrj;
 570                              selec=i;}}
 571    ///hardest jet selected
 572     if(selec<0){PostData(1, fOutputList);
 573                return;}
 574     AliAODJet* jet1 = (AliAODJet*)(fListJets[0]->At(selec));
 575     //What is the hardest constituent track?
 576                        AliAODTrack* leadtrack1;
 577                        Int_t ippt=0;
 578                        Double_t ppt=-10;
 579                        TRefArray *genTrackList = jet1->GetRefTracks();
 580                        Int_t nTracksGenJet = genTrackList->GetEntriesFast();
 581                        AliAODTrack* genTrack;
 582                        for(Int_t ir=0; ir<nTracksGenJet; ++ir){
 583                        genTrack = (AliAODTrack*)(genTrackList->At(ir));
 584                        if(genTrack->Pt()>ppt){ppt=genTrack->Pt();
 585                          ippt=ir;}}
 586                         leadtrack1=(AliAODTrack*)(genTrackList->At(ippt));
 587                         //If it doesn't exist or if it is greater that 100 GeV, discard.
 588                         if(!leadtrack1)  {PostData(1, fOutputList);
 589                           return;}
 590                         if(leadtrack1->Pt()>=100.){ PostData(1, fOutputList);
 591                           return;}
 592
 593            //Look to the back-to-back jet
 594            Int_t btb=-1;
 595            for(Int_t j=1;j<fListJets[0]->GetEntries();j++){
 596            if(j==selec) continue;
 597            AliAODJet* jetb = (AliAODJet*)(fListJets[0]->At(j));
 598            etaj  = jetb->Eta();
 599            //           phij  = jetb->Phi();
 600            ptj   = jetb->Pt();
 601            if(ptj<=0) continue;
 602            areaj = jetb->EffectiveAreaCharged();
 603            ptcorrj=ptj-rho*areaj;
 604            if(ptcorrj<=0) continue;
 605            if((etaj<fJetEtaMin)||(etaj>fJetEtaMax)) continue;
 606            Double_t dphij=RelativePhi(jetb->Phi(),jet1->Phi());
 607            if(TMath::Abs(dphij)>TMath::Pi()-0.2) { btb=j;
 608              break;}}
 609
 610             AliAODJet* jet2 = (AliAODJet*)(fListJets[0]->At(btb));
 611            //the back-to-back jet is also identified
 612
 613            if(btb<0){PostData(1, fOutputList);
 614                           return;}
 615
 616             Double_t ptcorr1=jet1->Pt()-rho*jet1->EffectiveAreaCharged();
 617             Double_t ptcorr2=jet2->Pt()-rho*jet2->EffectiveAreaCharged();
 618
 619             if(centValue<10.) fh2Pt1Pt2C10->Fill(ptcorr1,ptcorr2);
 620             if(centValue>40.) fh2Pt1Pt2C40->Fill(ptcorr1,ptcorr2);
 621
 622
 623            Double_t px2=jet2->Px();
 624            Double_t py2=jet2->Py();
 625            Double_t pz2=jet2->Pz();
 626            Double_t phi2=jet2->Phi();
 627            Double_t eta2=jet2->Eta();
 628            //Once we have have a dijet event,look to the structure of the back-to-back jet:
 629
 630       TVector3  ppJ1(px2, py2, pz2);
 631       TVector3  ppJ3(- px2 * pz2, - py2 * pz2, px2 * px2 + py2 * py2);
 632       ppJ3.SetMag(1.);
 633       TVector3  ppJ2(-py2, px2, 0);
 634       ppJ2.SetMag(1.);
 635       Float_t mxx    = 0.;
 636       Float_t myy    = 0.;
 637       Float_t mxy    = 0.;
 638       Int_t   nc     = 0;
 639       Float_t sump2  = 0.;
 640       Float_t ptMax  = 0.;
 641       Float_t etaMax = 0.;
 642       Float_t phiMax = 0.;
 643       //      Int_t   iMax   = -1;
 644       //1st loop over all tracks
 645          for(int it = 0;it<nT;++it){
 646          AliVParticle *track = (AliVParticle*)ParticleList.At(it);
 647          TVector3 pp(track->Px(), track->Py(), track->Pz());
 648           Float_t phi = track->Phi();
 649           Float_t eta = track->Eta();
 650           Float_t pt  = track->Pt();
 651           Float_t jT  = pp.Perp(ppJ1);
 652           Float_t deta = eta - eta2;
 653           Float_t dphi = phi - phi2;
 654           if (dphi >   TMath::Pi()) dphi =   2. * TMath::Pi() - dphi;
 655           if (dphi < - TMath::Pi()) dphi = - 2. * TMath::Pi() - dphi;
 656           Float_t r = TMath::Sqrt(dphi * dphi + deta * deta);
 657           /////////To compute the TM axis, we use particles with large jT
 658
 659               //cout<<"Jt spectrum............"<<ptbig<<" "<<jT<<endl;
 660           /////////We compute the TM with large jT particles only
 661              if((jT >1.)&&(r<1.2)){
 662             //longitudinal and perpendicular component of the track pT in the
 663             //local frame
 664               TVector3 pLong = pp.Dot(ppJ1) / ppJ1.Mag2() * ppJ1;
 665               TVector3 pPerp = pp - pLong;
 666               //projection onto the two perpendicular vectors defined above
 667               Float_t ppjX = pPerp.Dot(ppJ2);
 668               Float_t ppjY = pPerp.Dot(ppJ3);
 669               Float_t ppjT = TMath::Sqrt(ppjX * ppjX + ppjY * ppjY);
 670               //components of the sphericity matrix
 671               mxx += (ppjX * ppjX / ppjT);
 672               myy += (ppjY * ppjY / ppjT);
 673               mxy += (ppjX * ppjY / ppjT);
 674               nc++;
 675               sump2 += ppjT;
 676               // max pt
 677               if (pt > ptMax) {
 678                   ptMax  = pt;
 679                   //              iMax   = it;
 680                   etaMax = deta;
 681                   phiMax = dphi;
 682               }
 683                 } // R < 0.4
 684          } // 1st Track Loop
 685
 686   // At this point we have mxx, myy, mxy
 687       if (nc == 0) return;
 688 // Shericity Matrix
 689       const Double_t ele[4] = {mxx / sump2, mxy / sump2, mxy / sump2, myy / sump2};
 690       TMatrixDSym m0(2,ele);
 691 // Find eigenvectors
 692       TMatrixDSymEigen m(m0);
 693       TVectorD eval(2);
 694       TMatrixD evecm = m.GetEigenVectors();
 695       eval  = m.GetEigenValues();
 696 // Largest eigenvector
 697       Int_t jev = 0;
 698       if (eval[0] < eval[1]) jev = 1;
 699       TVectorD evec0(2);
 700 // Principle axis
 701       evec0 = TMatrixDColumn(evecm, jev);
 702       TVector2 evec(evec0[0], evec0[1]);
 703 // Principle axis from leading partice
 704       Float_t phiM = TMath::ATan2(phiMax, etaMax);
 705       TVector2 evecM(TMath::Cos(phiM), TMath::Sin(phiM));
 706       Float_t phistM = evecM.DeltaPhi(evec);
 707       if (TMath::Abs(phistM) > TMath::Pi()/2.) evec*=(-1.);
 708
 709       //////we have now the direction
 710       /////along which the sum of the projections of the particle
 711       ///momentum is higher.
 712       Double_t sum2lpt20=0;
 713       Double_t sum4lpt20=0;
 714       Double_t sum6lpt20=0;
 715       Double_t sum8lpt20=0;
 716       Double_t sum12lpt20=0;
 717       Double_t sum2pt20=0;
 718       Double_t sum4pt20=0;
 719       Double_t sum6pt20=0;
 720       Double_t sum8pt20=0;
 721       Double_t sum12pt20=0;
 722
 723        Double_t sum2lpt40=0;
 724       Double_t sum4lpt40=0;
 725       Double_t sum6lpt40=0;
 726       Double_t sum8lpt40=0;
 727       Double_t sum12lpt40=0;
 728       Double_t sum2pt40=0;
 729       Double_t sum4pt40=0;
 730       Double_t sum6pt40=0;
 731       Double_t sum8pt40=0;
 732       Double_t sum12pt40=0;
 733
 734           for (Int_t ip = 0; ip < nT; ip++) {
 735           AliVParticle *track = (AliVParticle*)ParticleList.At(ip);
 736           TVector3 pp(track->Px(), track->Py(), track->Pz());
 737           Float_t phi = track->Phi();
 738           Float_t eta = track->Eta();
 739           Float_t pt  = track->Pt();
 740           Float_t jT  = pp.Perp(ppJ1);
 741           Double_t deta = eta - eta2;
 742           Double_t deltaPhi = phi - phi2;
 743           Double_t r = TMath::Sqrt(deltaPhi * deltaPhi + deta * deta);
 744
 745             if(ptcorr2>20. && ptcorr2<40.){
 746             if(pt<0.4){
 747                  if(r<0.2) sum2lpt20=sum2lpt20-1.*pt*TMath::Cos(phi-jet1->Phi());
 748                  if(r<0.4) sum4lpt20=sum4lpt20-1.*pt*TMath::Cos(phi-jet1->Phi());
 749                  if(r<0.6) sum6lpt20=sum6lpt20-1.*pt*TMath::Cos(phi-jet1->Phi());
 750                  if(r<0.8) sum8lpt20=sum8lpt20-1.*pt*TMath::Cos(phi-jet1->Phi());
 751                  if(r<1.2) sum12lpt20=sum12lpt20-1.*pt*TMath::Cos(phi-jet1->Phi());}
 752
 753                  if(r<0.2) sum2pt20=sum2pt20-1.*pt*TMath::Cos(phi-jet1->Phi());
 754                  if(r<0.4) sum4pt20=sum4pt20-1.*pt*TMath::Cos(phi-jet1->Phi());
 755                  if(r<0.6) sum6pt20=sum6pt20-1.*pt*TMath::Cos(phi-jet1->Phi());
 756                  if(r<0.8) sum8pt20=sum8pt20-1.*pt*TMath::Cos(phi-jet1->Phi());
 757                  if(r<1.2) sum12pt20=sum12pt20-1.*pt*TMath::Cos(phi-jet1->Phi());}
 758
 759
 760                 if(ptcorr2>40. && ptcorr2<60.){
 761             if(pt<0.4){
 762                  if(r<0.2) sum2lpt40=sum2lpt40-1.*pt*TMath::Cos(phi-jet1->Phi());
 763                  if(r<0.4) sum4lpt40=sum4lpt40-1.*pt*TMath::Cos(phi-jet1->Phi());
 764                  if(r<0.6) sum6lpt40=sum6lpt40-1.*pt*TMath::Cos(phi-jet1->Phi());
 765                  if(r<0.8) sum8lpt40=sum8lpt40-1.*pt*TMath::Cos(phi-jet1->Phi());
 766                  if(r<1.2) sum12lpt40=sum12lpt40-1.*pt*TMath::Cos(phi-jet1->Phi());}
 767
 768                  if(r<0.2) sum2pt40=sum2pt40-1.*pt*TMath::Cos(phi-jet1->Phi());
 769                  if(r<0.4) sum4pt40=sum4pt40-1.*pt*TMath::Cos(phi-jet1->Phi());
 770                  if(r<0.6) sum6pt40=sum6pt40-1.*pt*TMath::Cos(phi-jet1->Phi());
 771                  if(r<0.8) sum8pt40=sum8pt40-1.*pt*TMath::Cos(phi-jet1->Phi());
 772                  if(r<1.2) sum12pt40=sum12pt40-1.*pt*TMath::Cos(phi-jet1->Phi());}
 773
 774
 775
 776
 777           TVector3 pLong = pp.Dot(ppJ1) / ppJ1.Mag2() * ppJ1;
 778           TVector3 pPerp = pp - pLong;
 779           Float_t ppjX = pPerp.Dot(ppJ2);
 780           Float_t ppjY = pPerp.Dot(ppJ3);
 781           TVector2 vr(ppjX, ppjY) ;
 782           //and this is the angle between the particle and the TM axis.
 783           Float_t phistr = evec.DeltaPhi(vr);
 784
 785           //Double_t phistr=vr.Phi()-evec.Phi();
 786
 787           if(centValue<10.) fh3LocalCoordinates->Fill(ppjX,ppjY,ptcorr2);
 788           Double_t deltaEta = eta2-track->Eta();
 789
 790
 791           if(deltaPhi<-0.5*TMath::Pi()) deltaPhi+=2.*TMath::Pi();
 792           if(deltaPhi>3./2.*TMath::Pi()) deltaPhi-=2.*TMath::Pi();
 793           Double_t jetEntries[9] = {centValue,ptcorr1,ptcorr2,jT,pt,deltaEta,deltaPhi,phistr,ptMax};
 794           fhnDeltaR->Fill(jetEntries);
 795           }
 796           if(centValue<10.){
 797             if(ptcorr2>20.){
 798               if(ptcorr1>80.){
 799           Double_t aj=(ptcorr1-ptcorr2)/(ptcorr1+ptcorr2);
 800
 801           fh2Sum2pt20->Fill(aj,sum2pt20);
 802           fh2Sum4pt20->Fill(aj,sum4pt20);
 803           fh2Sum6pt20->Fill(aj,sum6pt20);
 804           fh2Sum8pt20->Fill(aj,sum8pt20);
 805           fh2Sum12pt20->Fill(aj,sum12pt20);
 806           fh2Sum2lpt20->Fill(aj,sum2lpt20);
 807           fh2Sum4lpt20->Fill(aj,sum4lpt20);
 808           fh2Sum6lpt20->Fill(aj,sum6lpt20);
 809           fh2Sum8lpt20->Fill(aj,sum8lpt20);
 810           fh2Sum12lpt20->Fill(aj,sum12lpt20);
 811
 812           fh2Sum2pt40->Fill(aj,sum2pt40);
 813           fh2Sum4pt40->Fill(aj,sum4pt40);
 814           fh2Sum6pt40->Fill(aj,sum6pt40);
 815           fh2Sum8pt40->Fill(aj,sum8pt40);
 816           fh2Sum12pt40->Fill(aj,sum12pt40);
 817           fh2Sum2lpt40->Fill(aj,sum2lpt40);
 818           fh2Sum4lpt40->Fill(aj,sum4lpt40);
 819           fh2Sum6lpt40->Fill(aj,sum6lpt40);
 820           fh2Sum8lpt40->Fill(aj,sum8lpt40);
 821           fh2Sum12lpt40->Fill(aj,sum12lpt40);
 822               }}}
 823
 824
 825
 826
 827
 828
 829
 830    PostData(1, fOutputList);
 831 }
 832
 833 void AliAnalysisTaskAj::Terminate(const Option_t *)
 834 {
 835    // Draw result to the screen
 836    // Called once at the end of the query
 837
 838    if (!GetOutputData(1))
 839    return;
 840 }
 841
 842
 843
 844
 845
 846
 847
 848
 849
 850
 851
 852 Int_t  AliAnalysisTaskAj::GetListOfTracks(TList *list){
 853
 854     Int_t iCount = 0;
 855
 856      AliAODEvent *aod = 0;
 857      if(!fESD)aod = fAODIn;
 858      else aod = fAODOut;
 859      if(!aod) return iCount;
 860
 861     for(int it = 0;it < aod->GetNumberOfTracks();++it){
 862       AliAODTrack *tr = dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
 863       if(!tr) AliFatal("Not a standard AOD");
 864       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
 865       if(TMath::Abs(tr->Eta())>0.9)continue;
 866       if(tr->Pt()<0.15)continue;
 867       list->Add(tr);
 868       //cout<<fAOD->GetNumberOfTracks()<<" "<<tr->Pt()<<endl;
 869       iCount++;
 870     }
 871
 872
 873   return iCount;
 874
 875 }
 876
 877    Int_t  AliAnalysisTaskAj::GetHardestTrackBackToJet(AliAODJet *jetbig){
 878
 879       Int_t index=-1;
 880      AliAODEvent *aod = 0;
 881      if(!fESD)aod = fAODIn;
 882      else aod = fAODOut;
 883      if(!aod) return index;
 884
 885
 886     Double_t ptmax=-10;
 887     Double_t dphi=0;
 888     //    Double_t dif=0;
 889     Int_t iCount=0;
 890     for(int it = 0;it < aod->GetNumberOfTracks();++it){
 891       AliAODTrack *tr = dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
 892       if(!tr) AliFatal("Not a standard AOD");
 893       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
 894       if(TMath::Abs(tr->Eta())>0.9)continue;
 895       if(tr->Pt()<0.15)continue;
 896       iCount=iCount+1;
 897       dphi=RelativePhi(tr->Phi(),jetbig->Phi());
 898       if(TMath::Abs(dphi)<TMath::Pi()-0.2) continue;
 899       if(tr->Pt()>ptmax){ ptmax=tr->Pt();
 900       index=iCount-1;
 901       //  dif=dphi;
 902       }}
 903
 904       return index;
 905
 906    }
 907
 908
 909
 910
 911
 912
 913
 914
 915
 916  Int_t  AliAnalysisTaskAj::GetListOfTracksCloseToJet(TList *list,AliAODJet *jetbig){
 917
 918     Int_t iCount = 0;
 919      AliAODEvent *aod = 0;
 920      if(!fESD)aod = fAODIn;
 921      else aod = fAODOut;
 922      if(!aod) return iCount;
 923
 924
 925
 926     for(int it = 0;it < aod->GetNumberOfTracks();++it){
 927       AliAODTrack *tr = dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
 928       if(!tr) AliFatal("Not a standard AOD");
 929       if((fFilterMask>0)&&!(tr->TestFilterBit(fFilterMask)))continue;
 930       if(TMath::Abs(tr->Eta())>0.9)continue;
 931       if(tr->Pt()<0.15)continue;
 932       Double_t disR=jetbig->DeltaR(tr);
 933       if(disR>0.8)  continue;
 934       list->Add(tr);
 935       //cout<<fAOD->GetNumberOfTracks()<<" "<<tr->Pt()<<endl;
 936       iCount++;
 937     }
 938
 939    list->Sort();
 940    return iCount;
 941
 942 }
 943
 944
 945
 946
 947
 948
 949
 950
 951
 952
 953
 954 Int_t AliAnalysisTaskAj::GetNInputTracks()
 955 {
 956
 957    Int_t nInputTracks = 0;
 958       AliAODEvent *aod = 0;
 959      if(!fESD)aod = fAODIn;
 960      else aod = fAODOut;
 961      if(!aod) return nInputTracks;
 962
 963
 964    TString jbname(fJetBranchName[1]);
 965    //needs complete event, use jets without background subtraction
 966    for(Int_t i=1; i<=3; ++i){
 967       if(jbname.Contains(Form("B%d",i))) jbname.ReplaceAll(Form("B%d",i),"B0");
 968    }
 969    // use only HI event
 970    if(jbname.Contains("AODextraonly")) jbname.ReplaceAll("AODextraonly","AOD");
 971    if(jbname.Contains("AODextra")) jbname.ReplaceAll("AODextra","AOD");
 972
 973    if(fDebug) Printf("Multiplicity from jet branch %s", jbname.Data());
 974    TClonesArray *tmpAODjets = dynamic_cast<TClonesArray*>(aod->FindListObject(jbname.Data()));
 975    if(!tmpAODjets){
 976       Printf("Jet branch %s not found", jbname.Data());
 977       Printf("AliAnalysisTaskAj::GetNInputTracks FAILED");
 978       return -1;
 979    }
 980
 981    for (Int_t iJet=0; iJet<tmpAODjets->GetEntriesFast(); iJet++){
 982       AliAODJet *jet = dynamic_cast<AliAODJet*>((*tmpAODjets)[iJet]);
 983       if(!jet) continue;
 984       TRefArray *trackList = jet->GetRefTracks();
 985       Int_t nTracks = trackList->GetEntriesFast();
 986       nInputTracks += nTracks;
 987       if(fDebug) Printf("#jet%d: %d tracks", iJet, nTracks);
 988    }
 989    if(fDebug) Printf("---> input tracks: %d", nInputTracks);
 990
 991    return nInputTracks;
 992 }
 993
 994
 995
 996 Double_t AliAnalysisTaskAj::RelativePhi(Double_t mphi,Double_t vphi){
 997
 998   if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi());
 999   else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
1000   if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi());
1001   else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
1002   double dphi = mphi-vphi;
1003   if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
1004   else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
1005   return dphi;//dphi in [-Pi, Pi]
1006 }
1007
1008
1009
1010 THnSparse* AliAnalysisTaskAj::NewTHnSparseF(const char* name, UInt_t entries)
1011 {
1012    // generate new THnSparseF, axes are defined in GetDimParams()
1013
1014    Int_t count = 0;
1015    UInt_t tmp = entries;
1016    while(tmp!=0){
1017       count++;
1018       tmp = tmp &~ -tmp;  // clear lowest bit
1019    }
1020
1021    TString hnTitle(name);
1022    const Int_t dim = count;
1023    Int_t nbins[dim];
1024    Double_t xmin[dim];
1025    Double_t xmax[dim];
1026
1027    Int_t i=0;
1028    Int_t c=0;
1029    while(c<dim && i<32){
1030       if(entries&(1<<i)){
1031
1032          TString label("");
1033          GetDimParams(i, label, nbins[c], xmin[c], xmax[c]);
1034          hnTitle += Form(";%s",label.Data());
1035          c++;
1036       }
1037
1038       i++;
1039    }
1040    hnTitle += ";";
1041
1042    return new THnSparseF(name, hnTitle.Data(), dim, nbins, xmin, xmax);
1043 }
1044
1045 void AliAnalysisTaskAj::GetDimParams(Int_t iEntry, TString &label, Int_t &nbins, Double_t &xmin, Double_t &xmax)
1046 {
1047    // stores label and binning of axis for THnSparse
1048
1049    const Double_t pi = TMath::Pi();
1050
1051    switch(iEntry){
1052
1053    case 0:
1054       label = "V0 centrality (%)";
1055
1056          nbins = 3;
1057          xmin = 0.;
1058          xmax = 100.;
1059          break;
1060
1061
1062    case 1:
1063       label = "corrected jet pt1";
1064          nbins = 10;
1065          xmin = 0.;
1066          xmax = 200.;
1067           break;
1068
1069       case 2:
1070       label = "corrected jet pt2";
1071          nbins = 10;
1072          xmin = 0.;
1073          xmax = 200.;
1074           break;
1075
1076
1077
1078    case 3:
1079       label = "track jT";
1080
1081          nbins = 2;
1082          xmin = 0.;
1083          xmax = 150.;
1084          break;
1085
1086    case 4:
1087       label = "track pT";
1088
1089          nbins = 4;
1090          xmin = 0.;
1091          xmax = 150.;
1092          break;
1093
1094
1095       case 5:
1096       label = "deltaEta";
1097       nbins = 8;
1098       xmin = -1.6;
1099       xmax = 1.6;
1100       break;
1101
1102
1103       case 6:
1104       label = "deltaPhi";
1105       nbins = 60;
1106       xmin = -0.5*pi;
1107       xmax = 1.5*pi;
1108       break;
1109
1110       case 7:
1111       label = "deltaPhiTM";
1112       nbins = 60;
1113       xmin = 0.;
1114       xmax = 1.3*pi;
1115       break;
1116
1117
1118
1119       case 8:
1120       label = "leading track";
1121       nbins = 3;
1122       xmin = 0;
1123       xmax = 50;
1124       break;
1125
1126
1127
1128
1129
1130
1131
1132
1133
1134
1135    }
1136
1137 }
1138