new UE Task from Sara
[u/mrichter/AliRoot.git] / PWG4 / JetTasks / AliAnalysisTaskUE.cxx
CommitLineData
f3050824 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id:$ */
6ef5bfa4 17
00a1afbd 18//class AliAODInputHandler;
19//class AliAODHandler;
20//class AliLog;
21//class TROOT;
22
f3050824 23#include <TROOT.h>
f3050824 24#include <TChain.h>
25#include <TFile.h>
26#include <TList.h>
f3050824 27#include <TMath.h>
6ef5bfa4 28#include <TTree.h>
00a1afbd 29#include <TVector3.h>
30#include <TH1F.h>
f3050824 31
00a1afbd 32#include "AliAnalyseUE.h"
f3050824 33#include "AliAnalysisTaskUE.h"
00a1afbd 34#include "AliHistogramsUE.h"
35
f3050824 36#include "AliAnalysisManager.h"
37#include "AliMCEventHandler.h"
f3050824 38
623ed0a0 39#include "AliAnalysisHelperJetTasks.h"
00a1afbd 40#include "AliAODHandler.h"
41#include "AliAODInputHandler.h"
42#include "AliGenPythiaEventHeader.h"
f3050824 43#include "AliLog.h"
00a1afbd 44#include "AliInputEventHandler.h"
f3050824 45
00a1afbd 46////////////////////////////////////////////////////////////////////////
47//
f3050824 48// Analysis class for Underlying Event studies
49//
50// Look for correlations on the tranverse regions to
51// the leading charged jet
52//
02480db2 53// This class needs as input AOD with track and Jets.
54// The output is a list of histograms
f3050824 55//
56// AOD can be either connected to the InputEventHandler
57// for a chain of AOD files
58// or
59// to the OutputEventHandler
60// for a chain of ESD files, so this case class should be
61// in the train after the Jet finder
62//
63// Arian.Abrahantes.Quintana@cern.ch
64// Ernesto.Lopez.Torres@cern.ch
65//
00a1afbd 66////////////////////////////////////////////////////////////////////////
f3050824 67
68
69ClassImp( AliAnalysisTaskUE)
70
00a1afbd 71// Define global pointer
72AliAnalysisTaskUE* AliAnalysisTaskUE::fgTaskUE=NULL;
f3050824 73
74//____________________________________________________________________
75AliAnalysisTaskUE:: AliAnalysisTaskUE(const char* name):
94ee88d4 76AliAnalysisTask(name,""),
00a1afbd 77fAnaUE(0x0),
6ef5bfa4 78fAOD(0x0),
00a1afbd 79fAODBranch("jets"),
80fDebug(0),
6ef5bfa4 81fListOfHistos(0x0),
82fBinsPtInHist(30),
baf36668 83fIsNorm2Area(kTRUE),
00a1afbd 84fMaxJetPtInHist(300.),
85fMinJetPtInHist(0.),
623ed0a0 86fConstrainDistance(kTRUE),
87fMinDistance(0.2),
88fSimulateChJetPt(kFALSE),
89fUseAliStack(kTRUE),
00a1afbd 90fUseMCParticleBranch(kFALSE),
94ee88d4 91fnTracksVertex(3), // QA tracks pointing to principal vertex (= 3 default)
92fZVertex(5.),
6ef5bfa4 93fAnaType(1),
623ed0a0 94fConePosition(1),
00a1afbd 95fConeRadius(0.7),
96fFilterBit(0xFF),
97fJetsOnFly(kFALSE),
98fRegionType(1),
9d9a3c85 99fUseChargeHadrons(kFALSE),
00a1afbd 100fUseChPartJet(kFALSE),
6ef5bfa4 101fUsePositiveCharge(kTRUE),
00a1afbd 102fUseSingleCharge(kFALSE),
6ef5bfa4 103fOrdering(1),
623ed0a0 104fChJetPtMin(5.0),
6ef5bfa4 105fJet1EtaCut(0.2),
106fJet2DeltaPhiCut(2.616), // 150 degrees
107fJet2RatioPtCut(0.8),
108fJet3PtCut(15.),
6ef5bfa4 109fTrackEtaCut(0.9),
00a1afbd 110fTrackPtCut(0.),
111fAvgTrials(1)
f3050824 112{
113 // Default constructor
114 // Define input and output slots here
115 // Input slot #0 works with a TChain
116 DefineInput(0, TChain::Class());
117 // Output slot #0 writes into a TList container
118 DefineOutput(0, TList::Class());
94ee88d4 119
f3050824 120}
121
00a1afbd 122//____________________________________________________________________
123AliAnalysisTaskUE:: AliAnalysisTaskUE(const AliAnalysisTaskUE & original):
124AliAnalysisTask(),
125fAnaUE(original.fAnaUE),
126fAOD(original.fAOD),
127fAODBranch(original.fAODBranch),
128fDebug(original.fDebug),
129fListOfHistos(original.fListOfHistos),
130fBinsPtInHist(original.fBinsPtInHist),
131fIsNorm2Area(original.fIsNorm2Area),
132fMaxJetPtInHist(original.fMaxJetPtInHist),
133fMinJetPtInHist(original.fMinJetPtInHist),
134fConstrainDistance(original.fConstrainDistance),
135fMinDistance(original.fMinDistance),
136fSimulateChJetPt(original.fSimulateChJetPt),
137fUseAliStack(original.fUseAliStack),
138fUseMCParticleBranch(original.fUseMCParticleBranch),
139fnTracksVertex(original.fnTracksVertex), // QA tracks pointing to principal vertex (= 3 default)
140fZVertex(original.fZVertex),
141fAnaType(original.fAnaType),
142fConePosition(original.fConePosition),
143fConeRadius(original.fConeRadius),
144fFilterBit(original.fFilterBit),
145fJetsOnFly(original.fJetsOnFly),
146fRegionType(original.fRegionType),
147fUseChargeHadrons(original.fUseChargeHadrons),
148fUseChPartJet(original.fUseChPartJet),
149fUsePositiveCharge(original.fUsePositiveCharge),
150fUseSingleCharge(original.fUseSingleCharge),
151fOrdering(original.fOrdering),
152fChJetPtMin(original.fChJetPtMin),
153fJet1EtaCut(original.fJet1EtaCut),
154fJet2DeltaPhiCut(original.fJet2DeltaPhiCut), // 150 degrees
155fJet2RatioPtCut(original.fJet2RatioPtCut),
156fJet3PtCut(original.fJet3PtCut),
157fTrackEtaCut(original.fTrackEtaCut),
158fTrackPtCut(original.fTrackPtCut),
159fAvgTrials(original.fAvgTrials)
160{
161 // Copy constructor
162}
163
164
165//______________________________________________________________
166AliAnalysisTaskUE & AliAnalysisTaskUE::operator = (const AliAnalysisTaskUE & /*source*/)
167{
168 // assignment operator
169 return *this;
170}
171
623ed0a0 172//______________________________________________________________
173Bool_t AliAnalysisTaskUE::Notify()
174{
175 //
176 // Implemented Notify() to read the cross sections
177 // and number of trials from pyxsec.root
519378fb 178 // Copy from AliAnalysisTaskJFSystematics
856f9829 179 fAvgTrials = 1;
623ed0a0 180 TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
519378fb 181 Float_t xsection = 0;
00a1afbd 182 Float_t trials = 1;
623ed0a0 183 if(tree){
00a1afbd 184 TFile *curfile = tree->GetCurrentFile();
185 if (!curfile) {
186 Error("Notify","No current file");
187 return kFALSE;
188 }
189
190 AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,trials);
191
192 fAnaUE->FillXsec("<#sigma>",xsection);
9d9a3c85 193
00a1afbd 194 // construct average trials
195 Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
196 if(trials>=nEntries && nEntries>0.)fAvgTrials = trials/nEntries;
197 }
c4396fd4 198
623ed0a0 199 return kTRUE;
200}
201
f3050824 202//____________________________________________________________________
203void AliAnalysisTaskUE::ConnectInputData(Option_t* /*option*/)
204{
6ef5bfa4 205 // Connect the input data
f3050824 206
6ef5bfa4 207 // We need AOD with tracks and jets.
208 // Since AOD can be either connected to the InputEventHandler (input chain fron AOD files)
1e996f55 209 // or to the OutputEventHandler ( AOD is create by a previus task in the train )
6ef5bfa4 210 // we need to check where it is and get the pointer to AODEvent in the right way
f3050824 211
00a1afbd 212 // Delta AODs are also accepted
c4396fd4 213
214
215 if (fDebug > 1) AliInfo("ConnectInputData() ");
f3050824 216
217 TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
218
c4396fd4 219 if( handler && handler->InheritsFrom("AliAODInputHandler") ) { // input AOD
00a1afbd 220 fAOD = ((AliAODInputHandler*)handler)->GetEvent();
221 if(!fJetsOnFly){
222 if (fDebug > 1) AliInfo(" ==== Tracks and Jets from AliAODInputHandler");
223 }else{
224 AliInfo(" ==== Tracks from AliAODInputHandler / Jets on-the-fly");
225 }
226 } else { //output AOD
227 handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
228 if( handler && handler->InheritsFrom("AliAODHandler") ) {
229 fAOD = ((AliAODHandler*)handler)->GetAOD();
230 if (!fJetsOnFly){
231 if (fDebug > 1) AliInfo(" ==== Tracks and Jets from AliAODHandler");
232 } else {
233 if (fDebug > 1) AliInfo(" ==== Tracks from AliAODHandler / Jets on-the-fly");
234 }
235 }else {
236 AliFatal("I can't get any AOD Event Handler");
237 return;
238 }
239 }
240
241 fAnaUE->Initialize( *this );
f3050824 242}
243
244//____________________________________________________________________
245void AliAnalysisTaskUE::CreateOutputObjects()
246{
6ef5bfa4 247 // Create the output container
00a1afbd 248
f3050824 249 if (fDebug > 1) AliInfo("CreateOutPutData()");
00a1afbd 250
251 //Initialize AliAnalysisUE, a class implementing the main analysis algorithms
252 fAnaUE = new AliAnalyseUE();
253
254 //Initialize output histograms
255 fAnaUE->CreateHistograms(fBinsPtInHist, fMinJetPtInHist, fMaxJetPtInHist);
256 fListOfHistos = (TList*)fAnaUE->GetHistograms()->Clone();
257 fListOfHistos->SetOwner(kTRUE);
94ee88d4 258 PostData(0,fListOfHistos);
1e996f55 259
f3050824 260}
261
f3050824 262//____________________________________________________________________
263void AliAnalysisTaskUE::Exec(Option_t */*option*/)
264{
00a1afbd 265 // Trigger selection ************************************************
266 AliInputEventHandler* inputHandler = (AliInputEventHandler*)
216601f0 267 ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
00a1afbd 268 if( !inputHandler->InheritsFrom("AliAODInputHandler") ) { // input AOD
269 if (inputHandler->IsEventSelected()) {
270 if (fDebug > 1) AliInfo(" Trigger Selection: event ACCEPTED ... ");
271 } else {
272 if (fDebug > 1) AliInfo(" Trigger Selection: event REJECTED ... ");
273 return;
274 }
275 }
276 // Event selection (vertex) *****************************************
277
278 if(!fAnaUE->VertexSelection(fAOD,fnTracksVertex,fZVertex)) return;
279
280 // Execute analysis for current event ******************************
281
6ef5bfa4 282 if ( fDebug > 3 ) AliInfo( " Processing event..." );
00a1afbd 283
519378fb 284 // fetch the pythia header info and get the trials
285 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
286 Float_t nTrials = 1;
c4396fd4 287 if (mcHandler) {
00a1afbd 288 AliMCEvent* mcEvent = mcHandler->MCEvent();
289 if (mcEvent) {
290 AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
291 if(pythiaGenHeader){
292 nTrials = pythiaGenHeader->Trials();
293 }
294 }
295 }
296 fAnaUE->FillTrials("#sum{ntrials}",fAvgTrials);
297
298 //analyse the event
02480db2 299 AnalyseUE();
00a1afbd 300 fListOfHistos = (TList*)fAnaUE->GetHistograms()->Clone();
301 PostData(0,fListOfHistos);
f3050824 302}
f3050824 303
304//____________________________________________________________________
305void AliAnalysisTaskUE::AnalyseUE()
306{
00a1afbd 307
6ef5bfa4 308
00a1afbd 309 // Get jets and order by pT
310 TVector3 jetVect[3];
311 *jetVect = fAnaUE->GetOrderedClusters(fAODBranch, fUseChPartJet, fChJetPtMin);
c4396fd4 312
00a1afbd 313 if( jetVect[0].Pt() < 0. ) {
314 if( fDebug > 1 ) AliInfo("\n Skipping Event, not jet found...");
315 return;
316 } else {
317 if (fDebug >1 ) AliInfo(Form("\n Pt Leading Jet = %6.1f eta=%5.3f ", jetVect[0].Pt(), jetVect[0].Eta() ));
318 }
c4396fd4 319
00a1afbd 320 // Select events according to analysis type
321 if ( ! (fAnaUE->AnaTypeSelection(jetVect))) return;
02480db2 322
00a1afbd 323 // Find max and min regions with real tracks
324 if (!fUseMCParticleBranch){
325 // Primary vertex distribution
326 AliAODVertex* vertex = (AliAODVertex*)fAOD->GetPrimaryVertex();
327 fAnaUE->FillVertex(vertex->GetNContributors());
328
329 // Fill leading "jet" histogram
330 fAnaUE->FillLeadingJet(jetVect[0].Pt());
02480db2 331
00a1afbd 332 fAnaUE->FindMaxMinRegions( jetVect, fConePosition );
baf36668 333
00a1afbd 334 }else {
6ef5bfa4 335
623ed0a0 336 // this is the part we only use when we have MC information
337 // More than a test for values of it also resumes the reconstruction efficiency of jets
338 // As commented bellow if available for the data, we try to pair reconstructed jets with simulated ones
339 // afterwards we kept angular variables of MC jet to perform UE analysis over MC particles
340 // TODO: Handle Multiple jet environment. 06/2009 just suited for inclusive jet condition ( fAnaType = 1 )
6ef5bfa4 341
623ed0a0 342 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
343 if (!mcHandler) {
344 Printf("ERROR: Could not retrieve MC event handler");
345 return;
346 }
347
348 AliMCEvent* mcEvent = mcHandler->MCEvent();
349 if (!mcEvent) {
350 Printf("ERROR: Could not retrieve MC event");
351 return;
352 }
353 AliGenPythiaEventHeader* pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
354 if(!pythiaGenHeader){
355 return;
356 }
f3050824 357
00a1afbd 358 fAnaUE->AnalyseMC(jetVect,mcEvent,pythiaGenHeader, fConePosition, fUseAliStack, fConstrainDistance, fMinDistance);
623ed0a0 359
f3050824 360 }
f3050824 361
00a1afbd 362 fAnaUE->FillRegions(fIsNorm2Area, jetVect);
6ef5bfa4 363
f3050824 364}
365
366//____________________________________________________________________
00a1afbd 367AliAnalysisTaskUE* AliAnalysisTaskUE::Instance()
368{
f3050824 369
00a1afbd 370 //Create instance
371 if (fgTaskUE) {
372 return fgTaskUE;
373 } else {
374 fgTaskUE = new AliAnalysisTaskUE();
375 return fgTaskUE;
376 }
02480db2 377}
6ef5bfa4 378
379//____________________________________________________________________
00a1afbd 380void AliAnalysisTaskUE::Terminate(Option_t */*option*/)
6ef5bfa4 381{
94ee88d4 382
00a1afbd 383 // Terminate analysis
6ef5bfa4 384
00a1afbd 385 fListOfHistos = dynamic_cast<TList*> (GetOutputData(0));
386 if (!fListOfHistos){
387 AliError("Histogram List is not available");
388 return;
6ef5bfa4 389 }
00a1afbd 390 if (!gROOT->IsBatch()){
391 //call class AliHistogramsUE
392 AliHistogramsUE *histos=new AliHistogramsUE(fListOfHistos);
393 histos->DrawUE(0);
6ef5bfa4 394 } else {
00a1afbd 395 AliInfo(" Batch mode, not histograms will be shown...");
6ef5bfa4 396 }
00a1afbd 397
6ef5bfa4 398 if( fDebug > 1 )
399 AliInfo("End analysis");
00a1afbd 400
6ef5bfa4 401}
f3050824 402
6ef5bfa4 403void AliAnalysisTaskUE::WriteSettings()
1e996f55 404{
00a1afbd 405
406 //Print analysis settings on screen
407 if (fDebug > 5){
6ef5bfa4 408 AliInfo(" All Analysis Settings in Saved Tree");
00a1afbd 409 fAnaUE->WriteSettings();
6ef5bfa4 410 }
411}