Added task to show AliPIDCombined use
authorantoniol <antoniol@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 8 Dec 2011 19:42:30 +0000 (19:42 +0000)
committerantoniol <antoniol@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 8 Dec 2011 19:42:30 +0000 (19:42 +0000)
ANALYSIS/ANALYSISaliceLinkDef.h
ANALYSIS/AliAnalysisTaskPIDCombined.cxx [new file with mode: 0644]
ANALYSIS/AliAnalysisTaskPIDCombined.h [new file with mode: 0644]
ANALYSIS/CMakelibANALYSISalice.pkg

index 3ec2be9..c05ab3e 100644 (file)
@@ -42,6 +42,7 @@
 #pragma link C++ class AliHEPDataParser+;
 #pragma link C++ class AliAnalysisTaskPIDqa+;
 #pragma link C++ class AliAnalysisTaskBaseLine+;
+#pragma link C++ class AliAnalysisTaskPIDCombined+;
 
 #pragma link C++ class AliEventPool+;
 #pragma link C++ class AliEventPoolManager+;
diff --git a/ANALYSIS/AliAnalysisTaskPIDCombined.cxx b/ANALYSIS/AliAnalysisTaskPIDCombined.cxx
new file mode 100644 (file)
index 0000000..c1ec628
--- /dev/null
@@ -0,0 +1,321 @@
+/*************************************************************************\r
+* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *\r
+*                                                                        *\r
+* Author: The ALICE Off-line Project.                                    *\r
+* Contributors are mentioned in the code where appropriate.              *\r
+*                                                                        *\r
+* Permission to use, copy, modify and distribute this software and its   *\r
+* documentation strictly for non-commercial purposes is hereby granted   *\r
+* without fee, provided that the above copyright notice appears in all   *\r
+* copies and that both the copyright notice and this permission notice   *\r
+* appear in the supporting documentation. The authors make no claims     *\r
+* about the suitability of this software for any purpose. It is          *\r
+* provided "as is" without express or implied warranty.                  *\r
+**************************************************************************/\r
+\r
+///////////////////////////////////////////////////////////////////////////\r
+//                                                                       //\r
+//                        Basic Analysis Task                            //\r
+//                      for PID        Analysis                          //\r
+//                                                                       //\r
+///////////////////////////////////////////////////////////////////////////\r
+\r
+#include <TChain.h>\r
+#include <TH1D.h>\r
+#include <TH2D.h>\r
+\r
+#include <AliVEvent.h>\r
+#include <AliInputEventHandler.h>\r
+#include <AliAnalysisManager.h>\r
+\r
+#include <AliLog.h>\r
+#include <AliPID.h>\r
+#include <AliPIDCombined.h>\r
+#include <AliPIDResponse.h>\r
+\r
+#include "AliAnalysisTaskPIDCombined.h"\r
+\r
+const char *AliAnalysisTaskPIDCombined::fgkBinMomDesc[AliAnalysisTaskPIDCombined::kPtBins] = {\r
+  " 0 <= p < 0.5 GeV/c",\r
+  " 0.5 <= p < 0.7 GeV/c",\r
+  " 0.7 <= p < 1.0 GeV/c",\r
+  " 1.0 <= p < 1.5 GeV/c",\r
+  " 1.5 <= p < 2.0 GeV/c",\r
+  " p >= 2.0 GeV/c"\r
+};\r
+\r
+ClassImp(AliAnalysisTaskPIDCombined)\r
+\r
+//_________________________________________________________________________________\r
+AliAnalysisTaskPIDCombined::AliAnalysisTaskPIDCombined() :\r
+  AliAnalysisTaskSE(),\r
+  fHistList(),\r
+  fPIDResponse(0x0),\r
+  fPIDCombined(0x0),\r
+  fTrackCuts(0x0),\r
+  fTrackFilter(0x0)\r
+{\r
+  //\r
+  // Constructor\r
+  //\r
+}\r
+\r
+//_________________________________________________________________________________\r
+AliAnalysisTaskPIDCombined::AliAnalysisTaskPIDCombined(const char *name) :\r
+  AliAnalysisTaskSE(name),\r
+  fHistList(),\r
+  fPIDResponse(0x0),\r
+  fPIDCombined(0x0),\r
+  fTrackCuts(0x0),\r
+  fTrackFilter(0x0)\r
+{\r
+  //\r
+  // Constructor\r
+  //\r
+  DefineInput(0,TChain::Class());\r
+  DefineOutput(1, TList::Class());\r
+}\r
+\r
+//_________________________________________________________________________________\r
+void AliAnalysisTaskPIDCombined::UserCreateOutputObjects()\r
+{\r
+  //\r
+  // Initialise the framework objects\r
+  //\r
+\r
+\r
+  // ------- track cuts\r
+  fTrackCuts = new AliESDtrackCuts("fTrackCuts", "Standard");\r
+  fTrackCuts->SetAcceptKinkDaughters(kFALSE);\r
+  fTrackCuts->SetMinNClustersTPC(80);\r
+  fTrackCuts->SetMaxChi2PerClusterTPC(4);\r
+  fTrackCuts->SetMaxDCAToVertexXY(3);\r
+  fTrackCuts->SetMaxDCAToVertexZ(3);\r
+  fTrackCuts->SetRequireTPCRefit(kTRUE);\r
+  fTrackCuts->SetRequireITSRefit(kTRUE);\r
+  fTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);\r
+  fTrackFilter = new AliAnalysisFilter("trackFilter");\r
+  fTrackFilter->AddCuts(fTrackCuts);\r
+\r
+\r
+\r
+  // ------- setup PIDCombined\r
+  fPIDCombined=new AliPIDCombined;\r
+  fPIDCombined->SetDefaultTPCPriors();\r
+  fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC+AliPIDResponse::kDetTOF);\r
+\r
+  // no light nuclei\r
+  fPIDCombined->SetSelectedSpecies(AliPID::kSPECIES);\r
+\r
+\r
+  fHistList.Add(new TH1D("nEvents","Number of Evnets;Selection",2,0,2));\r
+\r
+  for (Int_t ispec=0; ispec<AliPID::kSPECIES; ++ispec){\r
+\r
\r
+    fProbTPC[ispec]=new TH2D(Form("prob%s_mom_TPC",AliPID::ParticleName(ispec)),\r
+                                   Form("%s probability vs. momentum;momentum;probability",AliPID::ParticleName(ispec)),\r
+                                   100,0.,20.,50,0.,1.);\r
+    fHistList.Add(fProbTPC[ispec]);\r
+    fProbTPCnSigma[ispec]=new TH2D(Form("prob%s_nSigma_TPC",AliPID::ParticleName(ispec)),\r
+                                   Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
+                                   20,-5.,5.,50,0.,1.);\r
+    fHistList.Add(fProbTPCnSigma[ispec]);\r
+\r
+    for (Int_t ibin=0;ibin<kPtBins;ibin++) {\r
+      fProbTPCnSigTPCMom[ibin][ispec]=new TH2D(Form("prob%s_nSigma_TPC (%s)",AliPID::ParticleName(ispec),fgkBinMomDesc[ibin]),\r
+                                              Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
+                                              20,-5.,5.,50,0.,1.);\r
+      fHistList.Add(fProbTPCnSigTPCMom[ibin][ispec]);\r
+    }\r
+\r
+\r
+\r
+    fProbTOF[ispec]=new TH2D(Form("prob%s_mom_TOF",AliPID::ParticleName(ispec)),\r
+                                   Form("%s probability vs. momentum;momentum;probability",AliPID::ParticleName(ispec)),\r
+                                   100,0.,20.,50,0.,1.);\r
+    fHistList.Add(fProbTOF[ispec]);\r
+    fProbTOFnSigma[ispec]=new TH2D(Form("prob%s_nSigma_TOF",AliPID::ParticleName(ispec)),\r
+                                   Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
+                                   20,-5.,5.,50,0.,1.);\r
+    fHistList.Add(fProbTOFnSigma[ispec]);\r
+    for (Int_t ibin=0;ibin<kPtBins;ibin++) {\r
+      fProbTOFnSigTOFMom[ibin][ispec]=new TH2D(Form("prob%s_nSigma_TOF (%s)",AliPID::ParticleName(ispec),fgkBinMomDesc[ibin]),\r
+                                              Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
+                                              20,-5.,5.,50,0.,1.);\r
+      fHistList.Add(fProbTOFnSigTOFMom[ibin][ispec]);\r
+    }\r
+\r
+\r
+\r
+    fProbTPCTOF[ispec]=new TH2D(Form("prob%s_mom_TPCTOF",AliPID::ParticleName(ispec)),\r
+                                   Form("%s probability vs. momentum;momentum;probability",AliPID::ParticleName(ispec)),\r
+                                   100,0.,20.,50,0.,1.);\r
+    fHistList.Add(fProbTPCTOF[ispec]);\r
+    fProbTPCTOFnSigmaTPC[ispec]=new TH2D(Form("probTPCTOF%s_nSigma_TPC",AliPID::ParticleName(ispec)),\r
+                                   Form("%s TPCTOF probability vs. n#sigmaTPC;n#sigma;probability",AliPID::ParticleName(ispec)),\r
+                                   20,-5.,5.,50,0.,1.);\r
+    fHistList.Add(fProbTPCTOFnSigmaTPC[ispec]);\r
+    for (Int_t ibin=0;ibin<kPtBins;ibin++) {\r
+      fProbTPCTOFnSigTPCMom[ibin][ispec]=new TH2D(Form("probTPCTOF%s_nSigma_TPC (%s)",AliPID::ParticleName(ispec),fgkBinMomDesc[ibin]),\r
+                                              Form("%s probability vs. n#sigma;n#sigma;probability",AliPID::ParticleName(ispec)),\r
+                                              20,-5.,5.,50,0.,1.);\r
+      fHistList.Add(fProbTPCTOFnSigTPCMom[ibin][ispec]);\r
+    }\r
+\r
+\r
+\r
+    // basic priors\r
+    fPriors[ispec]=new TH1F(Form("%s_priors",AliPID::ParticleName(ispec)),\r
+                           Form("%s priors vs momentum",AliPID::ParticleName(ispec)),\r
+                           100,0.,20.);\r
+    fHistList.Add(fPriors[ispec]);\r
+    switch (ispec) {\r
+    case AliPID::kElectron:\r
+      for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.02);\r
+      break;\r
+    case AliPID::kMuon:\r
+      for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.02);\r
+      break;\r
+    case AliPID::kPion:\r
+      for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.56);\r
+      break;\r
+    case AliPID::kKaon:\r
+      for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.20);\r
+      break;\r
+    case AliPID::kProton:\r
+      for (Int_t ich=1;ich<=100;ich++) fPriors[ispec]->SetBinContent(ich,0.20);\r
+      break;\r
+    default:\r
+      break;\r
+    }\r
+    fPIDCombined->SetPriorDistribution((AliPID::EParticleType)ispec,fPriors[ispec]);\r
+\r
+    // priors used\r
+    fPriorsUsed[ispec] = new TH2D(Form("%s_priorsUsed",AliPID::ParticleName(ispec)),\r
+                           Form("%s priors vs transverse momentum;p_{t} (GeV/c);priors",AliPID::ParticleName(ispec)),\r
+                                 100,0.,20.,101,0,1.01);      \r
+    fHistList.Add(fPriorsUsed[ispec]);\r
+  }\r
+\r
+\r
+\r
+  fHistList.SetOwner();\r
+  PostData(1,&fHistList);\r
+\r
+\r
+}\r
+\r
+//_________________________________________________________________________________\r
+void AliAnalysisTaskPIDCombined::UserExec(Option_t *)\r
+{\r
+  //\r
+  // Main loop. Called for every event\r
+  //\r
+  AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();\r
+  AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());\r
+  fPIDResponse=inputHandler->GetPIDResponse();\r
+  if (!fPIDResponse) AliFatal("This Task needs the PID response attached to the inputHandler");\r
+\r
+  //  Printf(" ---------------------- UserExec PID task ---------------------");\r
+  \r
+  FillHistogram("nEvents",0.);\r
+\r
+  AliVEvent *event=InputEvent();\r
+  AliVTrack *track=0x0;\r
+  Int_t ntracks=event->GetNumberOfTracks();\r
+\r
+  Double_t probTPC[AliPID::kSPECIES]={0.};\r
+  Double_t probTOF[AliPID::kSPECIES]={0.};\r
+  Double_t probTPCTOF[AliPID::kSPECIES]={0.};\r
+  \r
+  //loop over all tracks\r
+  for (Int_t itrack=0; itrack<ntracks; ++itrack){\r
+\r
+    track=(AliVTrack*)event->GetTrack(itrack);\r
+\r
+    if ( fTrackFilter->IsSelected(track) ) {\r
+\r
+      Double_t mom=track->GetTPCmomentum();\r
+      Int_t ibin=GetMomBin(mom);\r
+\r
+      fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC);\r
+      UInt_t detUsed = fPIDCombined->ComputeProbabilities(track, fPIDResponse, probTPC);\r
+     \r
+      if (detUsed != 0) {  // TPC is available\r
+       for (Int_t ispec=0; ispec<AliPID::kSPECIES; ++ispec) {\r
+         Double_t nSigmaTPC = fPIDResponse->NumberOfSigmasTPC(track,(AliPID::EParticleType)ispec);\r
+         fProbTPC[ispec]->Fill(mom,probTPC[ispec]);\r
+         fProbTPCnSigma[ispec]->Fill(nSigmaTPC,probTPC[ispec]);\r
+         fProbTPCnSigTPCMom[ibin][ispec]->Fill(nSigmaTPC,probTPC[ispec]);\r
+       }\r
+\r
+       // compute priors for TPC+TOF, even if we ask just TOF for PID\r
+       fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTOF);\r
+       detUsed = fPIDCombined->ComputeProbabilities(track, fPIDResponse, probTOF);\r
+       Double_t priors[5];     // check priors used for TOF\r
+       fPIDCombined->GetPriors(track,priors,fPIDResponse,detUsed);\r
+       for(Int_t ispec=0;ispec<5;ispec++) fPriorsUsed[ispec]->Fill(TMath::Abs(track->Pt()),priors[ispec]);\r
+\r
+       if (detUsed != 0) {  // TOF is available\r
+         for (Int_t ispec=0; ispec<AliPID::kSPECIES; ++ispec) {\r
+           Double_t nSigmaTOF = fPIDResponse->NumberOfSigmasTOF(track,(AliPID::EParticleType)ispec);\r
+           fProbTOF[ispec]->Fill(mom,probTOF[ispec]);\r
+           fProbTOFnSigma[ispec]->Fill(nSigmaTOF,probTOF[ispec]);\r
+           fProbTOFnSigTOFMom[ibin][ispec]->Fill(nSigmaTOF,probTOF[ispec]);\r
+         }\r
+       }\r
+\r
+       fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTOF|AliPIDResponse::kDetTPC);\r
+       detUsed = fPIDCombined->ComputeProbabilities(track, fPIDResponse, probTPCTOF);\r
+       if (detUsed == (UInt_t)fPIDCombined->GetDetectorMask() ) {\r
+         for (Int_t ispec=0; ispec<AliPID::kSPECIES; ++ispec) {\r
+           Double_t nSigmaTPC = fPIDResponse->NumberOfSigmasTPC(track,(AliPID::EParticleType)ispec);\r
+           fProbTPCTOF[ispec]->Fill(mom,probTPCTOF[ispec]);\r
+           fProbTPCTOFnSigmaTPC[ispec]->Fill(nSigmaTPC,probTPCTOF[ispec]);\r
+           fProbTPCTOFnSigTPCMom[ibin][ispec]->Fill(nSigmaTPC,probTPCTOF[ispec]);\r
+         }\r
+       }\r
+\r
+      }\r
+\r
+\r
+    }\r
+  }\r
+\r
+  PostData(1, &fHistList);\r
+}\r
+\r
+//_________________________________________________________________________________\r
+void AliAnalysisTaskPIDCombined::FillHistogram(const char* name, Double_t x, Double_t weight)\r
+{\r
+  //\r
+  // Fill 1D histogram by name\r
+  //\r
+  ((TH1*)fHistList.FindObject(name))->Fill(x,weight);\r
+}\r
+\r
+//_________________________________________________________________________________\r
+void AliAnalysisTaskPIDCombined::FillHistogram(const char* name, Double_t x, Double_t y, Double_t weight)\r
+{\r
+  //\r
+  // Fill 2D histogram by name\r
+  //\r
+  ((TH2*)fHistList.FindObject(name))->Fill(x,y,weight);\r
+}\r
+\r
+\r
+//_________________________________________________________________________________\r
+Int_t AliAnalysisTaskPIDCombined::GetMomBin(Float_t mom)\r
+{\r
+  //\r
+  // Given momentum return histogram to be filled\r
+  //\r
+  if (mom>0. && mom < 0.5) return 0;\r
+  if (mom>=0.5 && mom < 0.7) return 1;\r
+  if (mom>=0.7 && mom < 1.0) return 2;\r
+  if (mom>=1.0 && mom < 1.5) return 3;\r
+  if (mom>=1.5 && mom < 2.0) return 4;\r
+  return kPtBins-1;\r
+}\r
+\r
diff --git a/ANALYSIS/AliAnalysisTaskPIDCombined.h b/ANALYSIS/AliAnalysisTaskPIDCombined.h
new file mode 100644 (file)
index 0000000..e4fd958
--- /dev/null
@@ -0,0 +1,75 @@
+#ifndef ALIPIDCOMBINEDTASK_H\r
+#define ALIPIDCOMBINEDTASK_H\r
+/* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *\r
+ * See cxx source for full Copyright notice                               */\r
+\r
+//#########################################################\r
+//#                                                       # \r
+//#        Task for testing the combined PID              #\r
+//#                                                       #\r
+//#  Pietro Antonioli, INFN / Pietro.Antonioli@bo.infn.it #\r
+//#  Jens Wiechula, Uni Tübingen / Jens.Wiechula@cern.ch  #\r
+//#                                                       #\r
+//#########################################################\r
+\r
+#include <AliPID.h>\r
+#include <AliPIDResponse.h>\r
+#include <AliPIDCombined.h>\r
+\r
+#include <AliESDtrackCuts.h>\r
+#include <AliAnalysisFilter.h>\r
+\r
+#include "AliAnalysisTaskSE.h"\r
+\r
+class TH1F;\r
+class TH2D;\r
+class fPIDResponse;\r
+class fPIDCombined;\r
+\r
+class AliAnalysisTaskPIDCombined : public AliAnalysisTaskSE {\r
+  \r
+public:\r
+  static const Int_t kPtBins = 6;\r
+\r
+\r
+  AliAnalysisTaskPIDCombined();\r
+  AliAnalysisTaskPIDCombined(const char *name);\r
+  virtual ~AliAnalysisTaskPIDCombined(){;}\r
+  \r
+  virtual void  UserExec(Option_t *option);\r
+  virtual void  UserCreateOutputObjects();\r
+    \r
+\r
+private:\r
+\r
+  TList fHistList;                   //! list of histograms\r
+  TH2D *fProbTPCnSigma[AliPID::kSPECIES];    //! probabilities vs nSigma in the TPC\r
+  TH2D *fProbTOFnSigma[AliPID::kSPECIES];    //! probabilities vs nSigma  the TOF\r
+  TH2D *fProbTPCTOFnSigmaTPC[AliPID::kSPECIES];    //! comb. probabilities vs nSigma TPC\r
+  TH2D *fProbTPC[AliPID::kSPECIES];          //! probabilities vs mom in the TPC\r
+  TH2D *fProbTOF[AliPID::kSPECIES];          //! probabilities vs mom in the TOF\r
+  TH2D *fProbTPCTOF[AliPID::kSPECIES];       //! combined probabilities vs mom TPC-TOF\r
+  TH1F *fPriors[AliPID::kSPECIES];           //! priors\r
+\r
+  TH2D *fProbTPCTOFnSigTPCMom[kPtBins][AliPID::kSPECIES];  // prob. x mom. bins\r
+  TH2D *fProbTPCnSigTPCMom[kPtBins][AliPID::kSPECIES];     // prob. x mom. bins\r
+  TH2D *fProbTOFnSigTOFMom[kPtBins][AliPID::kSPECIES];     // prob. x mom. bins\r
+\r
+  TH2D *fPriorsUsed[AliPID::kSPECIES];       //! priors used\r
+\r
+  const AliPIDResponse *fPIDResponse;     //! PID response object\r
+  AliPIDCombined       *fPIDCombined;     //! combined PID object\r
+  AliESDtrackCuts      *fTrackCuts;            //! track selection\r
+  AliAnalysisFilter    *fTrackFilter;         //! track filter\r
+\r
+  AliAnalysisTaskPIDCombined(const AliAnalysisTaskPIDCombined &c);\r
+  AliAnalysisTaskPIDCombined& operator= (const AliAnalysisTaskPIDCombined &c);\r
+\r
+  void FillHistogram(const char* name, Double_t x, Double_t weight=1.);\r
+  void FillHistogram(const char* name, Double_t x, Double_t y, Double_t weight=1.);\r
+  Int_t GetMomBin(Float_t mom);\r
+  static const char* fgkBinMomDesc[kPtBins];\r
+  \r
+  ClassDef(AliAnalysisTaskPIDCombined, 2);\r
+};\r
+#endif\r
index 7a0cff2..9476d8a 100644 (file)
@@ -59,6 +59,7 @@ set ( SRCS
     AliAnalysisTaskPIDqa.cxx 
     AliAnalysisTaskBaseLine.cxx
     AliEventPoolManager.cxx
+    AliAnalysisTaskPIDCombined.cxx     
     )
 
 if( ROOTHASALIEN STREQUAL "yes")