PWGLF/SPECTRA/ChargedHadrons/dNdPt/AliPtResolAnalysis.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 // AliPtResolAnalysis class.
  17 //
  18 // a. functionality:
  19 // - fills analysis control histograms
  20 //
  21 // b. data members:
  22 // - control histograms
  23 //
  24 // Author: J.Otwinowski 04/11/2008
  25 //------------------------------------------------------------------------------
  26
  27 #include "TH1.h"
  28 #include "TH2.h"
  29 #include "TCanvas.h"
  30 #include "THnSparse.h"
  31
  32 #include "AliHeader.h"
  33 #include "AliInputEventHandler.h"
  34 #include "AliAnalysisManager.h"
  35 #include "AliGenEventHeader.h"
  36 #include "AliStack.h"
  37 #include "AliESDEvent.h"
  38 #include "AliMCEvent.h"
  39 #include "AliESDtrackCuts.h"
  40 #include "AliLog.h"
  41 #include "AliMultiplicity.h"
  42 #include "AliTracker.h"
  43
  44 #include "AlidNdPtEventCuts.h"
  45 #include "AlidNdPtAcceptanceCuts.h"
  46 #include "AliPhysicsSelection.h"
  47 #include "AliTriggerAnalysis.h"
  48
  49 #include "AliPWG0Helper.h"
  50 #include "AlidNdPtHelper.h"
  51 #include "AliPtResolAnalysis.h"
  52
  53 using namespace std;
  54
  55 ClassImp(AliPtResolAnalysis)
  56
  57 //_____________________________________________________________________________
  58   AliPtResolAnalysis::AliPtResolAnalysis(): AlidNdPt(),
  59   fAnalysisFolder(0),
  60   fTrackParamHist(0),
  61   fTrackParamHist2(0)
  62 {
  63   // default constructor
  64   Init();
  65 }
  66
  67 //_____________________________________________________________________________
  68 AliPtResolAnalysis::AliPtResolAnalysis(Char_t* name, Char_t* title): AlidNdPt(name,title),
  69   fAnalysisFolder(0),
  70   fTrackParamHist(0),
  71   fTrackParamHist2(0)
  72 {
  73   Init();
  74 }
  75
  76 //_____________________________________________________________________________
  77 AliPtResolAnalysis::~AliPtResolAnalysis() {
  78   //
  79   // destructor
  80   //
  81   if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
  82   if(fTrackParamHist) delete fTrackParamHist; fTrackParamHist=0;
  83   if(fTrackParamHist2) delete fTrackParamHist2; fTrackParamHist2=0;
  84 }
  85
  86 //_____________________________________________________________________________
  87 void AliPtResolAnalysis::Init(){
  88   //
  89   // Generic histograms to be corrected
  90   //
  91   //1/pT:#sigma(1/pT)
  92   Int_t binsTrackParamHist[2]={400,300};
  93   Double_t minTrackParamHist[2]={0,0};
  94   Double_t maxTrackParamHist[2]={1,0.015};
  95
  96   fTrackParamHist = new THnSparseF("fTrackParamHist","1/pT:#sigma(1/pT)",2,binsTrackParamHist,minTrackParamHist,maxTrackParamHist);
  97   fTrackParamHist->GetAxis(0)->SetTitle("1/pT (GeV/c)^{-1}");
  98   fTrackParamHist->GetAxis(1)->SetTitle("#sigma(1/pT)");
  99   fTrackParamHist->Sumw2();
 100
 101   //pt:sigma(1/pT)*pT
 102   const Int_t ptNbins = 73;
 103   Double_t bins[74] = {0.0, 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0, 2.2, 2.4, 2.6, 2.8, 3.0, 3.2, 3.4, 3.6, 3.8, 4.0, 4.5, 5.0, 5.5, 6.0, 6.5, 7.0, 8.0, 9.0, 10.0, 11.0, 12.0, 13.0, 14.0, 15.0, 16.0, 18.0, 20.0, 22.0, 24.0, 26.0, 28.0, 30.0, 32.0, 34.0, 36.0, 40.0, 45.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0 };
 104
 105   Int_t binsTrackParamHist2[2]={ptNbins,200};
 106   Double_t minTrackParamHist2[2]={0,0};
 107   Double_t maxTrackParamHist2[2]={100,0.2};
 108
 109   fTrackParamHist2 = new THnSparseF("fTrackParamHist2","pT:#sigma(1/pT)*pT",2,binsTrackParamHist2,minTrackParamHist2,maxTrackParamHist2);
 110   fTrackParamHist2->SetBinEdges(0,bins);
 111   fTrackParamHist2->GetAxis(0)->SetTitle("pT (GeV/c)");
 112   fTrackParamHist2->GetAxis(1)->SetTitle("#sigma(1/pT)*pT");
 113   fTrackParamHist2->Sumw2();
 114
 115   // init folder
 116   fAnalysisFolder = CreateFolder("folderdNdPt","Analysis dNdPt Folder");
 117 }
 118
 119 //_____________________________________________________________________________
 120 void AliPtResolAnalysis::Process(AliESDEvent *const esdEvent, AliMCEvent *const mcEvent)
 121 {
 122   //
 123   // Process real and/or simulated events
 124   //
 125   if(!esdEvent) {
 126     AliDebug(AliLog::kError, "esdEvent not available");
 127     return;
 128   }
 129
 130   // get selection cuts
 131   AlidNdPtEventCuts *evtCuts = GetEventCuts();
 132   AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts();
 133   AliESDtrackCuts *esdTrackCuts = GetTrackCuts();
 134
 135   if(!evtCuts || !accCuts  || !esdTrackCuts) {
 136     AliDebug(AliLog::kError, "cuts not available");
 137     return;
 138   }
 139
 140   // trigger selection
 141   Bool_t isEventTriggered = kTRUE;
 142   AliPhysicsSelection *physicsSelection = NULL;
 143   AliTriggerAnalysis* triggerAnalysis = NULL;
 144
 145   //
 146   AliInputEventHandler* inputHandler = (AliInputEventHandler*) AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
 147   if (!inputHandler)
 148   {
 149     Printf("ERROR: Could not receive input handler");
 150     return;
 151   }
 152
 153   if(evtCuts->IsTriggerRequired())
 154   {
 155     // always MB
 156     isEventTriggered = inputHandler->IsEventSelected() & AliVEvent::kMB;
 157
 158     physicsSelection = static_cast<AliPhysicsSelection*> (inputHandler->GetEventSelection());
 159     if(!physicsSelection) return;
 160     //SetPhysicsTriggerSelection(physicsSelection);
 161
 162     if (isEventTriggered && (GetTrigger() == AliTriggerAnalysis::kV0AND)) {
 163       // set trigger (V0AND)
 164       triggerAnalysis = physicsSelection->GetTriggerAnalysis();
 165       if(!triggerAnalysis) return;
 166       isEventTriggered = triggerAnalysis->IsOfflineTriggerFired(esdEvent, GetTrigger());
 167     }
 168
 169   // get reconstructed vertex
 170   const AliESDVertex* vtxESD = 0;
 171   Bool_t isRecVertex = kFALSE;
 172   if(evtCuts->IsRecVertexRequired())
 173   {
 174     Bool_t bRedoTPCVertex = evtCuts->IsRedoTPCVertex();
 175     Bool_t bUseConstraints = evtCuts->IsUseBeamSpotConstraint();
 176     vtxESD = AlidNdPtHelper::GetVertex(esdEvent,evtCuts,accCuts,esdTrackCuts,GetAnalysisMode(),kFALSE,bRedoTPCVertex,bUseConstraints);
 177     isRecVertex = AlidNdPtHelper::TestRecVertex(vtxESD, esdEvent->GetPrimaryVertexSPD(), GetAnalysisMode(), kFALSE);
 178   }
 179
 180   Bool_t isEventOK = evtCuts->AcceptEvent(esdEvent,mcEvent,vtxESD) && isRecVertex;
 181   //printf("isEventOK %d, isEventTriggered %d \n",isEventOK, isEventTriggered);
 182   //printf("GetAnalysisMode() %d \n",GetAnalysisMode());
 183
 184   //
 185   TObjArray *allChargedTracks=0;
 186   // check event cuts
 187   if(isEventOK && isEventTriggered)
 188   {
 189     // get all charged tracks
 190     allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,GetAnalysisMode());
 191     if(!allChargedTracks) return;
 192
 193     Int_t entries = allChargedTracks->GetEntries();
 194     //printf("entries %d \n",entries);
 195
 196     // fill histograms
 197     for(Int_t i=0; i<entries;++i)
 198     {
 199       AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);
 200       if(!track) continue;
 201       if(track->Charge()==0) continue;
 202
 203       // only postive charged
 204       if(GetParticleMode() == AlidNdPtHelper::kPlus && track->Charge() < 0)
 205         continue;
 206
 207       // only negative charged
 208       if(GetParticleMode() == AlidNdPtHelper::kMinus && track->Charge() > 0)
 209         continue;
 210
 211       if(esdTrackCuts->AcceptTrack(track))
 212       {
 213         if(accCuts->AcceptTrack(track))
 214         {
 215           //Double_t x, par[5], cov[15];
 216           //track->GetExternalParameters(x, p);
 217           //track->GetExternalCovariance(cov);
 218
 219           Double_t v1[2] = {track->OneOverPt(),TMath::Sqrt(track->GetSigma1Pt2())};
 220           fTrackParamHist->Fill(v1);
 221
 222           Double_t v2[2] = {track->Pt(),track->Pt()*TMath::Sqrt(track->GetSigma1Pt2())};
 223           fTrackParamHist2->Fill(v2);
 224         }
 225       }
 226     }
 227   }
 228  }
 229 }
 230
 231 //_____________________________________________________________________________
 232 Long64_t AliPtResolAnalysis::Merge(TCollection* const list)
 233 {
 234   // Merge list of objects (needed by PROOF)
 235
 236   if (!list)
 237   return 0;
 238
 239   if (list->IsEmpty())
 240   return 1;
 241
 242   TIterator* iter = list->MakeIterator();
 243   TObject* obj = 0;
 244
 245   //
 246   //TList *collPhysSelection = new TList;
 247
 248   // collection of generated histograms
 249
 250   Int_t count=0;
 251   while((obj = iter->Next()) != 0) {
 252     AliPtResolAnalysis* entry = dynamic_cast<AliPtResolAnalysis*>(obj);
 253     if (entry == 0) continue;
 254
 255     //
 256     fTrackParamHist->Add(entry->fTrackParamHist);
 257     fTrackParamHist2->Add(entry->fTrackParamHist2);
 258   }
 259
 260 return count;
 261 }
 262
 263 //_____________________________________________________________________________
 264 void AliPtResolAnalysis::Analyse()
 265 {
 266   // Analyse histograms
 267   //
 268   TH1::AddDirectory(kFALSE);
 269   TObjArray *aFolderObj = new TObjArray;
 270   if(!aFolderObj) return;
 271
 272   //
 273   // Reconstructed event vertex
 274   //
 275
 276   // export objects to analysis folder
 277   fAnalysisFolder = ExportToFolder(aFolderObj);
 278   if(!fAnalysisFolder) {
 279     if(aFolderObj) delete aFolderObj;
 280     return;
 281   }
 282
 283   // delete only TObjArray
 284   if(aFolderObj) delete aFolderObj;
 285 }
 286
 287 //_____________________________________________________________________________
 288 TFolder* AliPtResolAnalysis::ExportToFolder(TObjArray * const array)
 289 {
 290   // recreate folder every time and export objects to new one
 291   //
 292   if(!array) return NULL;
 293
 294   AliPtResolAnalysis * comp=this;
 295   TFolder *folder = comp->GetAnalysisFolder();
 296
 297   TString name, title;
 298   TFolder *newFolder = 0;
 299   Int_t i = 0;
 300   Int_t size = array->GetSize();
 301
 302   if(folder) {
 303      // get name and title from old folder
 304      name = folder->GetName();
 305      title = folder->GetTitle();
 306
 307          // delete old one
 308      delete folder;
 309
 310          // create new one
 311      newFolder = CreateFolder(name.Data(),title.Data());
 312      newFolder->SetOwner();
 313
 314          // add objects to folder
 315      while(i < size) {
 316            newFolder->Add(array->At(i));
 317            i++;
 318          }
 319   }
 320
 321 return newFolder;
 322 }
 323
 324 //_____________________________________________________________________________
 325 TFolder* AliPtResolAnalysis::CreateFolder(TString name,TString title) {
 326 // create folder for analysed histograms
 327 //
 328 TFolder *folder = 0;
 329   folder = new TFolder(name.Data(),title.Data());
 330
 331   return folder;
 332 }