[u/mrichter/AliRoot.git] / ANALYSIS / AliAnalysisTaskPIDqa.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* $Id: AliAnalysisTaskPIDqa.cxx 43811 2010-09-23 14:13:31Z wiechula $ */
#include <TList.h>
#include <TVectorD.h>
#include <TObjArray.h>
#include <TH2.h>
#include <TFile.h>
#include <TPRegexp.h>
#include <TChain.h>
#include <TF1.h>
#include <TSpline.h>

#include <AliAnalysisManager.h>
#include <AliInputEventHandler.h>
#include <AliVEventHandler.h>
#include <AliVEvent.h>
#include <AliVParticle.h>
#include <AliVTrack.h>
#include <AliLog.h>
#include <AliPID.h>
#include <AliPIDResponse.h>
#include <AliITSPIDResponse.h>
#include <AliTPCPIDResponse.h>
#include <AliTRDPIDResponse.h>
#include <AliTOFPIDResponse.h>

#include "AliAnalysisTaskPIDqa.h"


ClassImp(AliAnalysisTaskPIDqa)

//______________________________________________________________________________
AliAnalysisTaskPIDqa::AliAnalysisTaskPIDqa():
AliAnalysisTaskSE(),
fPIDResponse(0x0),
fListQA(0x0),
fListQAits(0x0),
fListQAtpc(0x0),
fListQAtrd(0x0),
fListQAtof(0x0)
{
  //
  // Dummy constructor
  //
}

//______________________________________________________________________________
AliAnalysisTaskPIDqa::AliAnalysisTaskPIDqa(const char* name):
AliAnalysisTaskSE(name),
fPIDResponse(0x0),
fListQA(0x0),
fListQAits(0x0),
fListQAtpc(0x0),
fListQAtrd(0x0),
fListQAtof(0x0)
{
  //
  // Default constructor
  //
  DefineInput(0,TChain::Class());
  DefineOutput(1,TList::Class());
}

//______________________________________________________________________________
AliAnalysisTaskPIDqa::~AliAnalysisTaskPIDqa()
{
  //
  // Destructor
  //

}

//______________________________________________________________________________
void AliAnalysisTaskPIDqa::UserCreateOutputObjects()
{
  //
  // Create the output QA objects
  //

  AliLog::SetClassDebugLevel("AliAnalysisTaskPIDqa",10);

  //input hander
  AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
  AliInputEventHandler *inputHandler=dynamic_cast<AliInputEventHandler*>(man->GetInputEventHandler());
  if (!inputHandler) AliFatal("Input handler needed");

  //pid response object
  fPIDResponse=inputHandler->GetPIDResponse();
  if (!fPIDResponse) AliFatal("PIDResponse object was not created");
  
  //
  fListQA=new TList;
  fListQA->SetOwner();
  
  fListQAits=new TList;
  fListQAits->SetOwner();
  fListQAits->SetName("ITS");
  
  fListQAtpc=new TList;
  fListQAtpc->SetOwner();
  fListQAtpc->SetName("TPC");
  
  fListQAtrd=new TList;
  fListQAtrd->SetOwner();
  fListQAtrd->SetName("TRD");
  
  fListQAtof=new TList;
  fListQAtof->SetOwner();
  fListQAtof->SetName("TOF");
  
  fListQA->Add(fListQAits);
  fListQA->Add(fListQAtpc);
  fListQA->Add(fListQAtrd);
  fListQA->Add(fListQAtof);

  SetupITSqa();
  SetupTPCqa();
  SetupTRDqa();
  SetupTOFqa();

  PostData(1,fListQA);
}


//______________________________________________________________________________
void AliAnalysisTaskPIDqa::UserExec(Option_t */*option*/)
{
  //
  // Setup the PID response functions and fill the QA histograms
  //

  AliVEvent *event=InputEvent();
  if (!event||!fPIDResponse) return;

  
  FillITSqa();
  FillTPCqa();
  FillTOFqa();

  PostData(1,fListQA);
}

//______________________________________________________________________________
void AliAnalysisTaskPIDqa::FillITSqa()
{
  AliVEvent *event=InputEvent();
  
  Int_t ntracks=event->GetNumberOfTracks();
  for(Int_t itrack = 0; itrack < ntracks; itrack++){
    AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
    ULong_t status=track->GetStatus();
    // not that nice. status bits not in virtual interface
    // ITS refit + ITS pid
    if (!( ( (status&0x0004)==0x0004 ) && ( (status&0x0008)==0x0008 ) )) return;
    Double_t mom=track->P();
    
    for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
      TH2 *h=(TH2*)fListQAits->At(ispecie);
      if (!h) continue;
      Double_t nSigma=fPIDResponse->NumberOfSigmasITS(track, (AliPID::EParticleType)ispecie);
      h->Fill(mom,nSigma);
    }
    
    TH2 *h=(TH2*)fListQAits->At(AliPID::kSPECIES);
    if (h) {
      Double_t sig=track->GetITSsignal();
      h->Fill(mom,sig);
    }
    
  }
}

//______________________________________________________________________________
void AliAnalysisTaskPIDqa::FillTPCqa()
{
  AliVEvent *event=InputEvent();
  
  Int_t ntracks=event->GetNumberOfTracks();
  for(Int_t itrack = 0; itrack < ntracks; itrack++){
    AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
    
    //
    //basic track cuts
    //
    ULong_t status=track->GetStatus();
    // not that nice. status bits not in virtual interface
    // TPC refit + ITS refit + TPC pid
    if (!( (status&0x0040)==0x0040) || !( (status&0x0004)==0x0004) ||  !((status&0x0080)==0x0080) ) return;
    
    Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
    Float_t  ratioCrossedRowsOverFindableClustersTPC = 1.0;
    if (track->GetTPCNclsF()>0) {
      ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
    }
    
    if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
    
    Double_t mom=track->GetTPCmomentum();
    
    for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
      TH2 *h=(TH2*)fListQAtpc->At(ispecie);
      if (!h) continue;
      Double_t nSigma=fPIDResponse->NumberOfSigmasTPC(track, (AliPID::EParticleType)ispecie);
      h->Fill(mom,nSigma);
    }
    
    TH2 *h=(TH2*)fListQAtpc->At(AliPID::kSPECIES);
    if (h) {
      Double_t sig=track->GetTPCsignal();
      h->Fill(mom,sig);
    }
  }
}

//______________________________________________________________________________
void AliAnalysisTaskPIDqa::FillTOFqa()
{
  AliVEvent *event=InputEvent();
  
  Int_t ntracks=event->GetNumberOfTracks();
  for(Int_t itrack = 0; itrack < ntracks; itrack++){
    AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
    
    //
    //basic track cuts
    //
    ULong_t status=track->GetStatus();
    // not that nice. status bits not in virtual interface
    // TPC refit + ITS refit +
    // TOF out + TOFpid +
    // kTIME
    if (!((status&0x0040)==0x0040) || !((status&0x0004)==0x0004)
        || !((status&0x2000)==0x2000) || !((status&0x8000)==0x8000)
        || !((status&0x80000000)==0x80000000) ) return;
    
    Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
    Float_t  ratioCrossedRowsOverFindableClustersTPC = 1.0;
    if (track->GetTPCNclsF()>0) {
      ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
    }
    
    if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
    
    
    Double_t mom=track->P();
    
    for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
      TH2 *h=(TH2*)fListQAtof->At(ispecie);
      if (!h) continue;
      Double_t nSigma=fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)ispecie);
      h->Fill(mom,nSigma);
    }
    
    TH2 *h=(TH2*)fListQAtof->At(AliPID::kSPECIES);
    if (h) {
      Double_t sig=track->GetTOFsignal();
      h->Fill(mom,sig);
    }
    
  }
}

//______________________________________________________________________________
void AliAnalysisTaskPIDqa::SetupITSqa()
{
  //
  // Create the ITS qa objects
  //
  
  TVectorD *vX=MakeLogBinning(200,.1,30);
  
  for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
    TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_ITS_%s",AliPID::ParticleName(ispecie)),
                              Form("ITS n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
                              vX->GetNrows()-1,vX->GetMatrixArray(),
                              200,-10,10);
    fListQAits->Add(hNsigmaP);
  }
  
  
  TH2F *hSig = new TH2F("hSigP_ITS",
                        "ITS signal vs. p;p [GeV]; ITS signal [arb. units]",
                        vX->GetNrows()-1,vX->GetMatrixArray(),
                        300,0,300);
  
  fListQAits->Add(hSig);
  
}

//______________________________________________________________________________
void AliAnalysisTaskPIDqa::SetupTPCqa()
{
  //
  // Create the TPC qa objects
  //
  
  TVectorD *vX=MakeLogBinning(200,.1,30);
  
  for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
    TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TPC_%s",AliPID::ParticleName(ispecie)),
                              Form("TPC n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
                              vX->GetNrows()-1,vX->GetMatrixArray(),
                              200,-10,10);
    fListQAtpc->Add(hNsigmaP);
  }
  
  
  TH2F *hSig = new TH2F("hSigP_TPC",
                        "TPC signal vs. p;p [GeV]; TPC signal [arb. units]",
                        vX->GetNrows()-1,vX->GetMatrixArray(),
                        300,0,300);
  
  fListQAtpc->Add(hSig);
  
}

//______________________________________________________________________________
void AliAnalysisTaskPIDqa::SetupTRDqa()
{
  //
  // Create the TRD qa objects
  //
  
}

//______________________________________________________________________________
void AliAnalysisTaskPIDqa::SetupTOFqa()
{
  //
  // Create the TOF qa objects
  //
  
  TVectorD *vX=MakeLogBinning(200,.1,30);
  
  for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
    TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TOF_%s",AliPID::ParticleName(ispecie)),
                              Form("TOF n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
                              vX->GetNrows()-1,vX->GetMatrixArray(),
                              200,-10,10);
    fListQAtof->Add(hNsigmaP);
  }
  
  
  TH2F *hSig = new TH2F("hSigP_TOF",
                        "TOF signal vs. p;p [GeV]; TOF signal [arb. units]",
                        vX->GetNrows()-1,vX->GetMatrixArray(),
                        300,0,300);
  
  fListQAtof->Add(hSig);
  
}

//______________________________________________________________________________
TVectorD* AliAnalysisTaskPIDqa::MakeLogBinning(Int_t nbinsX, Double_t xmin, Double_t xmax)
{
  //
  // Make logarithmic binning
  // the user has to delete the array afterwards!!!
  //
  
  //check limits
  if (xmin<1e-20 || xmax<1e-20){
    AliError("For Log binning xmin and xmax must be > 1e-20. Using linear binning instead!");
    return MakeLinBinning(nbinsX, xmin, xmax);
  }
  if (xmax<xmin){
    Double_t tmp=xmin;
    xmin=xmax;
    xmax=tmp;
  }
  TVectorD *binLim=new TVectorD(nbinsX+1);
  Double_t first=xmin;
  Double_t last=xmax;
  Double_t expMax=TMath::Log(last/first);
  for (Int_t i=0; i<nbinsX+1; ++i){
    (*binLim)[i]=first*TMath::Exp(expMax/nbinsX*(Double_t)i);
  }
  return binLim;
}

//______________________________________________________________________________
TVectorD* AliAnalysisTaskPIDqa::MakeLinBinning(Int_t nbinsX, Double_t xmin, Double_t xmax)
{
  //
  // Make linear binning
  // the user has to delete the array afterwards!!!
  //
  if (xmax<xmin){
    Double_t tmp=xmin;
    xmin=xmax;
    xmax=tmp;
  }
  TVectorD *binLim=new TVectorD(nbinsX+1);
  Double_t first=xmin;
  Double_t last=xmax;
  Double_t binWidth=(last-first)/nbinsX;
  for (Int_t i=0; i<nbinsX+1; ++i){
    (*binLim)[i]=first+binWidth*(Double_t)i;
  }
  return binLim;
}

//_____________________________________________________________________________
TVectorD* AliAnalysisTaskPIDqa::MakeArbitraryBinning(const char* bins)
{
  //
  // Make arbitrary binning, bins separated by a ','
  //
  TString limits(bins);
  if (limits.IsNull()){
    AliError("Bin Limit string is empty, cannot add the variable");
    return 0x0;
  }
  
  TObjArray *arr=limits.Tokenize(",");
  Int_t nLimits=arr->GetEntries();
  if (nLimits<2){
    AliError("Need at leas 2 bin limits, cannot add the variable");
    delete arr;
    return 0x0;
  }
  
  TVectorD *binLimits=new TVectorD(nLimits);
  for (Int_t iLim=0; iLim<nLimits; ++iLim){
    (*binLimits)[iLim]=(static_cast<TObjString*>(arr->At(iLim)))->GetString().Atof();
  }
  
  delete arr;
  return binLimits;
}
Commit	Line	Data
1990d7b6	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: AliAnalysisTaskPIDqa.cxx 43811 2010-09-23 14:13:31Z wiechula $ */
	17	#include <TList.h>
	18	#include <TVectorD.h>
	19	#include <TObjArray.h>
	20	#include <TH2.h>
	21	#include <TFile.h>
	22	#include <TPRegexp.h>
	23	#include <TChain.h>
	24	#include <TF1.h>
	25	#include <TSpline.h>
	26
	27	#include <AliAnalysisManager.h>
	28	#include <AliInputEventHandler.h>
	29	#include <AliVEventHandler.h>
	30	#include <AliVEvent.h>
	31	#include <AliVParticle.h>
	32	#include <AliVTrack.h>
	33	#include <AliLog.h>
	34	#include <AliPID.h>
	35	#include <AliPIDResponse.h>
	36	#include <AliITSPIDResponse.h>
	37	#include <AliTPCPIDResponse.h>
	38	#include <AliTRDPIDResponse.h>
	39	#include <AliTOFPIDResponse.h>
	40
	41	#include "AliAnalysisTaskPIDqa.h"
	42
	43
	44	ClassImp(AliAnalysisTaskPIDqa)
	45
	46	//______________________________________________________________________________
	47	AliAnalysisTaskPIDqa::AliAnalysisTaskPIDqa():
	48	AliAnalysisTaskSE(),
	49	fPIDResponse(0x0),
	50	fListQA(0x0),
	51	fListQAits(0x0),
	52	fListQAtpc(0x0),
	53	fListQAtrd(0x0),
	54	fListQAtof(0x0)
	55	{
	56	//
	57	// Dummy constructor
	58	//
	59	}
	60
	61	//______________________________________________________________________________
	62	AliAnalysisTaskPIDqa::AliAnalysisTaskPIDqa(const char* name):
	63	AliAnalysisTaskSE(name),
	64	fPIDResponse(0x0),
65	fListQA(0x0),
66	fListQAits(0x0),
67	fListQAtpc(0x0),
68	fListQAtrd(0x0),
69	fListQAtof(0x0)
70	{
71	//
72	// Default constructor
73	//
74	DefineInput(0,TChain::Class());
75	DefineOutput(1,TList::Class());
76	}
77
78	//______________________________________________________________________________
79	AliAnalysisTaskPIDqa::~AliAnalysisTaskPIDqa()
80	{
81	//
82	// Destructor
83	//
84
85	}
86
87	//______________________________________________________________________________
88	void AliAnalysisTaskPIDqa::UserCreateOutputObjects()
89	{
90	//
91	// Create the output QA objects
92	//
93
94	AliLog::SetClassDebugLevel("AliAnalysisTaskPIDqa",10);
95
96	//input hander
97	AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
98	AliInputEventHandler inputHandler=dynamic_cast<AliInputEventHandler>(man->GetInputEventHandler());
99	if (!inputHandler) AliFatal("Input handler needed");
100
101	//pid response object
102	fPIDResponse=inputHandler->GetPIDResponse();
103	if (!fPIDResponse) AliFatal("PIDResponse object was not created");
104
105	//
106	fListQA=new TList;
107	fListQA->SetOwner();
108
109	fListQAits=new TList;
110	fListQAits->SetOwner();
111	fListQAits->SetName("ITS");
112
113	fListQAtpc=new TList;
114	fListQAtpc->SetOwner();
115	fListQAtpc->SetName("TPC");
116
117	fListQAtrd=new TList;
118	fListQAtrd->SetOwner();
119	fListQAtrd->SetName("TRD");
120
121	fListQAtof=new TList;
122	fListQAtof->SetOwner();
123	fListQAtof->SetName("TOF");
124
125	fListQA->Add(fListQAits);
126	fListQA->Add(fListQAtpc);
127	fListQA->Add(fListQAtrd);
128	fListQA->Add(fListQAtof);
129
130	SetupITSqa();
131	SetupTPCqa();
132	SetupTRDqa();
133	SetupTOFqa();
134
135	PostData(1,fListQA);
136	}
137
138
139	//______________________________________________________________________________
140	void AliAnalysisTaskPIDqa::UserExec(Option_t /option*/)
141	{
142	//
143	// Setup the PID response functions and fill the QA histograms
144	//
145
146	AliVEvent *event=InputEvent();
147	if (!event\|\|!fPIDResponse) return;
148
149
150	FillITSqa();
151	FillTPCqa();
152	FillTOFqa();
153
154	PostData(1,fListQA);
155	}
156
157	//______________________________________________________________________________
158	void AliAnalysisTaskPIDqa::FillITSqa()
159	{
160	AliVEvent *event=InputEvent();
161
162	Int_t ntracks=event->GetNumberOfTracks();
163	for(Int_t itrack = 0; itrack < ntracks; itrack++){
164	AliVTrack track=(AliVTrack)event->GetTrack(itrack);
165	ULong_t status=track->GetStatus();
166	// not that nice. status bits not in virtual interface
167	// ITS refit + ITS pid
168	if (!( ( (status&0x0004)==0x0004 ) && ( (status&0x0008)==0x0008 ) )) return;
169	Double_t mom=track->P();
170
171	for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
172	TH2 h=(TH2)fListQAits->At(ispecie);
173	if (!h) continue;
174	Double_t nSigma=fPIDResponse->NumberOfSigmasITS(track, (AliPID::EParticleType)ispecie);
175	h->Fill(mom,nSigma);
176	}
177
178	TH2 h=(TH2)fListQAits->At(AliPID::kSPECIES);
179	if (h) {
180	Double_t sig=track->GetITSsignal();
181	h->Fill(mom,sig);
182	}
183
184	}
185	}
186
187	//______________________________________________________________________________
188	void AliAnalysisTaskPIDqa::FillTPCqa()
189	{
190	AliVEvent *event=InputEvent();
191
192	Int_t ntracks=event->GetNumberOfTracks();
193	for(Int_t itrack = 0; itrack < ntracks; itrack++){
194	AliVTrack track=(AliVTrack)event->GetTrack(itrack);
195
196	//
197	//basic track cuts
198	//
199	ULong_t status=track->GetStatus();
200	// not that nice. status bits not in virtual interface
201	// TPC refit + ITS refit + TPC pid
202	if (!( (status&0x0040)==0x0040) \|\| !( (status&0x0004)==0x0004) \|\| !((status&0x0080)==0x0080) ) return;
203
204	Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
205	Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
206	if (track->GetTPCNclsF()>0) {
207	ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
208	}
209
210	if ( nCrossedRowsTPC<70 \|\| ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
211
212	Double_t mom=track->GetTPCmomentum();
213
214	for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
215	TH2 h=(TH2)fListQAtpc->At(ispecie);
216	if (!h) continue;
217	Double_t nSigma=fPIDResponse->NumberOfSigmasTPC(track, (AliPID::EParticleType)ispecie);
218	h->Fill(mom,nSigma);
219	}
220
221	TH2 h=(TH2)fListQAtpc->At(AliPID::kSPECIES);
222	if (h) {
223	Double_t sig=track->GetTPCsignal();
224	h->Fill(mom,sig);
225	}
226	}
227	}
228
229	//______________________________________________________________________________
230	void AliAnalysisTaskPIDqa::FillTOFqa()
231	{
232	AliVEvent *event=InputEvent();
233
234	Int_t ntracks=event->GetNumberOfTracks();
235	for(Int_t itrack = 0; itrack < ntracks; itrack++){
236	AliVTrack track=(AliVTrack)event->GetTrack(itrack);
237
238	//
239	//basic track cuts
240	//
241	ULong_t status=track->GetStatus();
242	// not that nice. status bits not in virtual interface
243	// TPC refit + ITS refit +
244	// TOF out + TOFpid +
245	// kTIME
246	if (!((status&0x0040)==0x0040) \|\| !((status&0x0004)==0x0004)
247	\|\| !((status&0x2000)==0x2000) \|\| !((status&0x8000)==0x8000)
248	\|\| !((status&0x80000000)==0x80000000) ) return;
249
250	Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
251	Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
252	if (track->GetTPCNclsF()>0) {
253	ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
254	}
255
256	if ( nCrossedRowsTPC<70 \|\| ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
257
258
259	Double_t mom=track->P();
260
261	for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
262	TH2 h=(TH2)fListQAtof->At(ispecie);
263	if (!h) continue;
264	Double_t nSigma=fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)ispecie);
265	h->Fill(mom,nSigma);
266	}
267
268	TH2 h=(TH2)fListQAtof->At(AliPID::kSPECIES);
269	if (h) {
270	Double_t sig=track->GetTOFsignal();
271	h->Fill(mom,sig);
272	}
273
274	}
275	}
276
277	//______________________________________________________________________________
278	void AliAnalysisTaskPIDqa::SetupITSqa()
279	{
280	//
281	// Create the ITS qa objects
282	//
283
284	TVectorD *vX=MakeLogBinning(200,.1,30);
285
286	for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
287	TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_ITS_%s",AliPID::ParticleName(ispecie)),
288	Form("ITS n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
289	vX->GetNrows()-1,vX->GetMatrixArray(),
290	200,-10,10);
291	fListQAits->Add(hNsigmaP);
292	}
293
294
295	TH2F *hSig = new TH2F("hSigP_ITS",
296	"ITS signal vs. p;p [GeV]; ITS signal [arb. units]",
297	vX->GetNrows()-1,vX->GetMatrixArray(),
298	300,0,300);
299
300	fListQAits->Add(hSig);
301
302	}
303
304	//______________________________________________________________________________
305	void AliAnalysisTaskPIDqa::SetupTPCqa()
306	{
307	//
308	// Create the TPC qa objects
309	//
310
311	TVectorD *vX=MakeLogBinning(200,.1,30);
312
313	for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
314	TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TPC_%s",AliPID::ParticleName(ispecie)),
315	Form("TPC n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
316	vX->GetNrows()-1,vX->GetMatrixArray(),
317	200,-10,10);
318	fListQAtpc->Add(hNsigmaP);
319	}
320
321
322	TH2F *hSig = new TH2F("hSigP_TPC",
323	"TPC signal vs. p;p [GeV]; TPC signal [arb. units]",
324	vX->GetNrows()-1,vX->GetMatrixArray(),
325	300,0,300);
326
327	fListQAtpc->Add(hSig);
328
329	}
330
331	//______________________________________________________________________________
332	void AliAnalysisTaskPIDqa::SetupTRDqa()
333	{
334	//
335	// Create the TRD qa objects
336	//
337
338	}
339
340	//______________________________________________________________________________
341	void AliAnalysisTaskPIDqa::SetupTOFqa()
342	{
343	//
344	// Create the TOF qa objects
345	//
346
347	TVectorD *vX=MakeLogBinning(200,.1,30);
348
349	for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
350	TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TOF_%s",AliPID::ParticleName(ispecie)),
351	Form("TOF n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
352	vX->GetNrows()-1,vX->GetMatrixArray(),
353	200,-10,10);
354	fListQAtof->Add(hNsigmaP);
355	}
356
357
358	TH2F *hSig = new TH2F("hSigP_TOF",
359	"TOF signal vs. p;p [GeV]; TOF signal [arb. units]",
360	vX->GetNrows()-1,vX->GetMatrixArray(),
361	300,0,300);
362
363	fListQAtof->Add(hSig);
364
365	}
366
367	//______________________________________________________________________________
368	TVectorD* AliAnalysisTaskPIDqa::MakeLogBinning(Int_t nbinsX, Double_t xmin, Double_t xmax)
369	{
370	//
371	// Make logarithmic binning
372	// the user has to delete the array afterwards!!!
373	//
374
375	//check limits
376	if (xmin<1e-20 \|\| xmax<1e-20){
377	AliError("For Log binning xmin and xmax must be > 1e-20. Using linear binning instead!");
378	return MakeLinBinning(nbinsX, xmin, xmax);
379	}
380	if (xmax<xmin){
381	Double_t tmp=xmin;
382	xmin=xmax;
383	xmax=tmp;
384	}
385	TVectorD *binLim=new TVectorD(nbinsX+1);
386	Double_t first=xmin;
387	Double_t last=xmax;
388	Double_t expMax=TMath::Log(last/first);
389	for (Int_t i=0; i<nbinsX+1; ++i){
390	(binLim)[i]=firstTMath::Exp(expMax/nbinsX*(Double_t)i);
391	}
392	return binLim;
393	}
394
395	//______________________________________________________________________________
396	TVectorD* AliAnalysisTaskPIDqa::MakeLinBinning(Int_t nbinsX, Double_t xmin, Double_t xmax)
397	{
398	//
399	// Make linear binning
400	// the user has to delete the array afterwards!!!
401	//
402	if (xmax<xmin){
403	Double_t tmp=xmin;
404	xmin=xmax;
405	xmax=tmp;
406	}
407	TVectorD *binLim=new TVectorD(nbinsX+1);
408	Double_t first=xmin;
409	Double_t last=xmax;
410	Double_t binWidth=(last-first)/nbinsX;
411	for (Int_t i=0; i<nbinsX+1; ++i){
412	(binLim)[i]=first+binWidth(Double_t)i;
413	}
414	return binLim;
415	}
416
417	//_____________________________________________________________________________
418	TVectorD* AliAnalysisTaskPIDqa::MakeArbitraryBinning(const char* bins)
419	{
420	//
421	// Make arbitrary binning, bins separated by a ','
422	//
423	TString limits(bins);
424	if (limits.IsNull()){
425	AliError("Bin Limit string is empty, cannot add the variable");
426	return 0x0;
427	}
428
429	TObjArray *arr=limits.Tokenize(",");
430	Int_t nLimits=arr->GetEntries();
431	if (nLimits<2){
432	AliError("Need at leas 2 bin limits, cannot add the variable");
433	delete arr;
434	return 0x0;
435	}
436
437	TVectorD *binLimits=new TVectorD(nLimits);
438	for (Int_t iLim=0; iLim<nLimits; ++iLim){
439	(binLimits)[iLim]=(static_cast<TObjString>(arr->At(iLim)))->GetString().Atof();
440	}
441
442	delete arr;
443	return binLimits;
444	}