PWGUD/dNdPt/AliPtResolAnalysis.cxx

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