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