PWG4/JetTasks/AliAnalysisTaskUE.cxx

   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:$ */
  17
  18 //class AliAODInputHandler;
  19 //class AliAODHandler;
  20 //class AliLog;
  21 //class TROOT;
  22
  23 #include <TROOT.h>
  24 #include <TChain.h>
  25 #include <TFile.h>
  26 #include <TList.h>
  27 #include <TMath.h>
  28 #include <TTree.h>
  29 #include <TVector3.h>
  30 #include <TH1F.h>
  31
  32 #include "AliAnalyseUE.h"
  33 #include "AliAnalysisTaskUE.h"
  34 #include "AliHistogramsUE.h"
  35
  36 #include "AliAnalysisManager.h"
  37 #include "AliMCEventHandler.h"
  38
  39 #include "AliAnalysisHelperJetTasks.h"
  40 #include "AliAODHandler.h"
  41 #include "AliAODInputHandler.h"
  42 #include "AliGenPythiaEventHeader.h"
  43 #include "AliLog.h"
  44 #include "AliInputEventHandler.h"
  45
  46 ////////////////////////////////////////////////////////////////////////
  47 //
  48 // Analysis class for Underlying Event studies
  49 //
  50 // Look for correlations on the tranverse regions to
  51 // the leading charged jet
  52 //
  53 // This class needs as input AOD with track and Jets.
  54 // The output is a list of histograms
  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 //
  66 ////////////////////////////////////////////////////////////////////////
  67
  68
  69 ClassImp( AliAnalysisTaskUE)
  70
  71 // Define global pointer
  72 AliAnalysisTaskUE* AliAnalysisTaskUE::fgTaskUE=NULL;
  73
  74 //____________________________________________________________________
  75 AliAnalysisTaskUE:: AliAnalysisTaskUE(const char* name):
  76 AliAnalysisTask(name,""),
  77 fAnaUE(0x0),
  78 fAOD(0x0),
  79 fAODBranch("jets"),
  80 fDebug(0),
  81 fListOfHistos(0x0),
  82 fBinsPtInHist(30),
  83 fIsNorm2Area(kTRUE),
  84 fMaxJetPtInHist(300.),
  85 fMinJetPtInHist(0.),
  86 fConstrainDistance(kTRUE),
  87 fMinDistance(0.2),
  88 fSimulateChJetPt(kFALSE),
  89 fUseAliStack(kTRUE),
  90 fUseMCParticleBranch(kFALSE),
  91 fnTracksVertex(3),  // QA tracks pointing to principal vertex (= 3 default)
  92 fZVertex(5.),
  93 fAnaType(1),
  94 fConePosition(1),
  95 fConeRadius(0.7),
  96 fFilterBit(0xFF),
  97 fJetsOnFly(kFALSE),
  98 fRegionType(1),
  99 fUseChargeHadrons(kFALSE),
 100 fUseChPartJet(kFALSE),
 101 fUsePositiveCharge(kTRUE),
 102 fUseSingleCharge(kFALSE),
 103 fOrdering(1),
 104 fChJetPtMin(5.0),
 105 fJet1EtaCut(0.2),
 106 fJet2DeltaPhiCut(2.616),    // 150 degrees
 107 fJet2RatioPtCut(0.8),
 108 fJet3PtCut(15.),
 109 fTrackEtaCut(0.9),
 110 fTrackPtCut(0.),
 111 fAvgTrials(1)
 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());
 119
 120 }
 121
 122 //____________________________________________________________________
 123 AliAnalysisTaskUE:: AliAnalysisTaskUE(const AliAnalysisTaskUE & original):
 124 AliAnalysisTask(),
 125 fAnaUE(original.fAnaUE),
 126 fAOD(original.fAOD),
 127 fAODBranch(original.fAODBranch),
 128 fDebug(original.fDebug),
 129 fListOfHistos(original.fListOfHistos),
 130 fBinsPtInHist(original.fBinsPtInHist),
 131 fIsNorm2Area(original.fIsNorm2Area),
 132 fMaxJetPtInHist(original.fMaxJetPtInHist),
 133 fMinJetPtInHist(original.fMinJetPtInHist),
 134 fConstrainDistance(original.fConstrainDistance),
 135 fMinDistance(original.fMinDistance),
 136 fSimulateChJetPt(original.fSimulateChJetPt),
 137 fUseAliStack(original.fUseAliStack),
 138 fUseMCParticleBranch(original.fUseMCParticleBranch),
 139 fnTracksVertex(original.fnTracksVertex),  // QA tracks pointing to principal vertex (= 3 default)
 140 fZVertex(original.fZVertex),
 141 fAnaType(original.fAnaType),
 142 fConePosition(original.fConePosition),
 143 fConeRadius(original.fConeRadius),
 144 fFilterBit(original.fFilterBit),
 145 fJetsOnFly(original.fJetsOnFly),
 146 fRegionType(original.fRegionType),
 147 fUseChargeHadrons(original.fUseChargeHadrons),
 148 fUseChPartJet(original.fUseChPartJet),
 149 fUsePositiveCharge(original.fUsePositiveCharge),
 150 fUseSingleCharge(original.fUseSingleCharge),
 151 fOrdering(original.fOrdering),
 152 fChJetPtMin(original.fChJetPtMin),
 153 fJet1EtaCut(original.fJet1EtaCut),
 154 fJet2DeltaPhiCut(original.fJet2DeltaPhiCut),    // 150 degrees
 155 fJet2RatioPtCut(original.fJet2RatioPtCut),
 156 fJet3PtCut(original.fJet3PtCut),
 157 fTrackEtaCut(original.fTrackEtaCut),
 158 fTrackPtCut(original.fTrackPtCut),
 159 fAvgTrials(original.fAvgTrials)
 160 {
 161   // Copy constructor
 162 }
 163
 164
 165 //______________________________________________________________
 166 AliAnalysisTaskUE & AliAnalysisTaskUE::operator = (const AliAnalysisTaskUE & /*source*/)
 167 {
 168   // assignment operator
 169   return *this;
 170 }
 171
 172 //______________________________________________________________
 173 Bool_t AliAnalysisTaskUE::Notify()
 174 {
 175   //
 176   // Implemented Notify() to read the cross sections
 177   // and number of trials from pyxsec.root
 178   // Copy from AliAnalysisTaskJFSystematics
 179   fAvgTrials = 1;
 180   TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
 181   Float_t xsection = 0;
 182   Float_t trials  = 1;
 183   if(tree){
 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);
 193
 194         // construct average trials
 195         Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
 196         if(trials>=nEntries && nEntries>0.)fAvgTrials = trials/nEntries;
 197         }
 198
 199   return kTRUE;
 200 }
 201
 202 //____________________________________________________________________
 203 void AliAnalysisTaskUE::ConnectInputData(Option_t* /*option*/)
 204 {
 205   // Connect the input data
 206
 207   // We need AOD with tracks and jets.
 208   // Since AOD can be either connected to the InputEventHandler (input chain fron AOD files)
 209   // or to the OutputEventHandler ( AOD is create by a previus task in the train )
 210   // we need to check where it is and get the pointer to AODEvent in the right way
 211
 212   // Delta AODs are also accepted
 213
 214
 215   if (fDebug > 1) AliInfo("ConnectInputData() ");
 216
 217   TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
 218
 219   if( handler && handler->InheritsFrom("AliAODInputHandler") ) { // input AOD
 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 );
 242 }
 243
 244 //____________________________________________________________________
 245 void  AliAnalysisTaskUE::CreateOutputObjects()
 246 {
 247   // Create the output container
 248
 249   if (fDebug > 1) AliInfo("CreateOutPutData()");
 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);
 258   PostData(0,fListOfHistos);
 259
 260 }
 261
 262 //____________________________________________________________________
 263 void  AliAnalysisTaskUE::Exec(Option_t */*option*/)
 264 {
 265   // Trigger selection ************************************************
 266   AliInputEventHandler* inputHandler = (AliInputEventHandler*)
 267          ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
 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
 282   if ( fDebug > 3 ) AliInfo( " Processing event..." );
 283
 284   // fetch the pythia header info and get the trials
 285   AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
 286   Float_t nTrials = 1;
 287   if (mcHandler) {
 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
 299   AnalyseUE();
 300   fListOfHistos = (TList*)fAnaUE->GetHistograms()->Clone();
 301   PostData(0,fListOfHistos);
 302 }
 303
 304 //____________________________________________________________________
 305 void  AliAnalysisTaskUE::AnalyseUE()
 306 {
 307
 308
 309    // Get jets and order by pT
 310    TVector3 jetVect[3];
 311    *jetVect = fAnaUE->GetOrderedClusters(fAODBranch, fUseChPartJet, fChJetPtMin);
 312
 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                 }
 319
 320    // Select events according to analysis type
 321    if ( ! (fAnaUE->AnaTypeSelection(jetVect))) return;
 322
 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());
 331
 332         fAnaUE->FindMaxMinRegions( jetVect,  fConePosition );
 333
 334   }else {
 335
 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 )
 341
 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     }
 357
 358     fAnaUE->AnalyseMC(jetVect,mcEvent,pythiaGenHeader, fConePosition, fUseAliStack, fConstrainDistance, fMinDistance);
 359
 360   }
 361
 362   fAnaUE->FillRegions(fIsNorm2Area, jetVect);
 363
 364 }
 365
 366 //____________________________________________________________________
 367 AliAnalysisTaskUE* AliAnalysisTaskUE::Instance()
 368 {
 369
 370   //Create instance
 371   if (fgTaskUE) {
 372         return fgTaskUE;
 373   } else {
 374         fgTaskUE = new AliAnalysisTaskUE();
 375         return fgTaskUE;
 376         }
 377 }
 378
 379 //____________________________________________________________________
 380 void  AliAnalysisTaskUE::Terminate(Option_t */*option*/)
 381 {
 382
 383   // Terminate analysis
 384
 385   fListOfHistos = dynamic_cast<TList*> (GetOutputData(0));
 386   if (!fListOfHistos){
 387         AliError("Histogram List is not available");
 388         return;
 389   }
 390   if (!gROOT->IsBatch()){
 391   //call class AliHistogramsUE
 392   AliHistogramsUE *histos=new AliHistogramsUE(fListOfHistos);
 393   histos->DrawUE(0);
 394   } else {
 395         AliInfo(" Batch mode, not histograms will be shown...");
 396   }
 397
 398   if( fDebug > 1 )
 399     AliInfo("End analysis");
 400
 401 }
 402
 403 void  AliAnalysisTaskUE::WriteSettings()
 404 {
 405
 406   //Print analysis settings on screen
 407   if (fDebug > 5){
 408     AliInfo(" All Analysis Settings in Saved Tree");
 409     fAnaUE->WriteSettings();
 410   }
 411 }