1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 // Cut on the Event at reconstructed level: for the moment
16 // the requirements on the number of charged tracks and on
17 // the vertex position and resolution are implemented
18 // The argument of IsSelected member function (passed object) is cast into
19 // an AliESDEvent. In the future may be modified to use AliVEvent interface
20 // and include more cut variables.
21 // The class derives from AliCFCutBase
22 // Author:S.Arcelli Silvia.Arcelli@cern.ch
29 #include "AliAnalysisUtils.h"
30 #include "AliESDEvent.h"
31 #include "AliESDVertex.h"
32 #include "AliAODEvent.h"
33 #include "AliAODVertex.h"
34 #include "AliHFEextraEventCuts.h"
35 ClassImp(AliHFEextraEventCuts)
36 //____________________________________________________________________
37 AliHFEextraEventCuts::AliHFEextraEventCuts() :
40 fRequireVtxCuts(kFALSE),
46 fCheckCorrelationSPDVtx(0),
47 fVtxResolution(kFALSE),
56 //____________________________________________________________________
57 AliHFEextraEventCuts::AliHFEextraEventCuts(Char_t* name, Char_t* title) :
58 AliCFCutBase(name,title),
60 fRequireVtxCuts(kFALSE),
66 fCheckCorrelationSPDVtx(0),
67 fVtxResolution(kFALSE),
76 fAnalysisUtils = new AliAnalysisUtils;
79 //____________________________________________________________________
80 AliHFEextraEventCuts::AliHFEextraEventCuts(const AliHFEextraEventCuts& c) :
83 fRequireVtxCuts(c.fRequireVtxCuts),
86 fVtxNCtrbMin(c.fVtxNCtrbMin),
87 fVtxMixed(c.fVtxMixed),
89 fCheckCorrelationSPDVtx(c.fCheckCorrelationSPDVtx),
90 fVtxResolution(c.fVtxResolution),
91 fPApileupCut(c.fPApileupCut),
97 for (Int_t i=0; i<c.kNCuts; i++){
98 for (Int_t j=0; j<c.kNStepQA; j++){
99 if(c.fhQA[i][j]) fhQA[i][j] = (TH1F*)c.fhQA[i][j]->Clone();
102 fAnalysisUtils = new AliAnalysisUtils;
105 //____________________________________________________________________
106 AliHFEextraEventCuts::~AliHFEextraEventCuts() {
111 for (Int_t i=0; i<kNCuts; i++){
112 for (Int_t j=0; j<kNStepQA; j++){
113 if(fhQA[i][j]) delete fhQA[i][j];
117 if(fBitMap) delete fBitMap;
118 if(fAnalysisUtils) delete fAnalysisUtils;
120 //_____________________________________________________________________________
121 void AliHFEextraEventCuts::Initialise()
129 // sets pointers to histos to zero
132 for(Int_t i=0; i<kNCuts; i++){
133 for(Int_t j =0; j<kNStepQA; j++){
139 //____________________________________________________________________
140 AliHFEextraEventCuts& AliHFEextraEventCuts::operator=(const AliHFEextraEventCuts& c)
143 // Assignment operator
146 AliCFCutBase::operator=(c) ;
147 fRequireVtxCuts=c.fRequireVtxCuts;
150 fVtxNCtrbMin=c.fVtxNCtrbMin;
151 fVtxMixed=c.fVtxMixed;
153 fCheckCorrelationSPDVtx=c.fCheckCorrelationSPDVtx;
154 fVtxResolution = c.fVtxResolution;
156 if(!fAnalysisUtils) fAnalysisUtils = new AliAnalysisUtils;
157 fPApileupCut = c.fPApileupCut;
160 for (Int_t i=0; i<c.kNCuts; i++){
161 for (Int_t j=0; j<c.kNStepQA; j++){
162 if(c.fhQA[i][j]) fhQA[i][j] = (TH1F*)c.fhQA[i][j]->Clone();
170 //____________________________________________________________________
171 Bool_t AliHFEextraEventCuts::IsSelected(TObject* obj) {
173 //Check if the requested cuts are passed
177 SelectionBitMap(obj);
179 if (fIsQAOn) FillHistograms(obj,0);
180 Bool_t isSelected = kTRUE;
182 for (UInt_t icut=0; icut<fBitMap->GetNbits();icut++)
183 if(!fBitMap->TestBitNumber(icut)) isSelected = kFALSE;
185 if (!isSelected) return kFALSE ;
186 if (fIsQAOn) FillHistograms(obj,1);
190 //____________________________________________________________________
191 void AliHFEextraEventCuts::SelectionBitMap(TObject* obj) {
193 //cut on the number of charged tracks and on the event vertex.
194 //so far specific to AliESDEvents
197 //Check if the requested cuts are passed and return a bitmap
199 AliVEvent *inputEvent = dynamic_cast<AliVEvent *>(obj);
201 AliDebug(1, "Not a virtual event");
202 for(Int_t j=0;j<kNCuts;j++)fBitMap->SetBitNumber(j,kFALSE);
205 const AliVVertex *vtxTracks = GetPrimaryVertexTracks(inputEvent),
206 *vtxSPD = GetPrimaryVertexSPD(inputEvent),
209 //first assume the event will be accepted:
210 for(Int_t j=0;j<kNCuts;j++)fBitMap->SetBitNumber(j,kTRUE);
213 // Use mixed vertex: Prefer vertex with tracks, in case not available use SPD vertex
214 if(vtxTracks && vtxTracks->GetNContributors() > 0) vtxPrim = vtxTracks;
215 else if(vtxSPD && vtxSPD->GetNContributors() > 0) vtxPrim = vtxSPD;
217 if(vtxSPD && vtxSPD->GetNContributors () > 0) vtxPrim = vtxSPD;
219 if(vtxTracks && vtxTracks->GetNContributors() > 0) vtxPrim = vtxTracks;
222 // No primary vertex: Reject event
223 for(Int_t j=1;j<kNCuts;j++)fBitMap->SetBitNumber(j,kFALSE);
224 AliWarning("Cannot get vertex, skipping event");
227 // Check standard vertex cuts using the primary vertex
228 if(fVtxZMin > - 1000. && fVtxZMax < 1000.){
229 // Primary vertex z cut required
230 // Pick up the position and uncertainties
232 vtxPos[0] = vtxPrim->GetX();
233 vtxPos[1] = vtxPrim->GetY();
234 vtxPos[2] = vtxPrim->GetZ();
235 if (vtxPos[2]>fVtxZMax || vtxPos[2]<fVtxZMin)
236 fBitMap->SetBitNumber(kVtxPosZ,kFALSE);
238 if(fVtxNCtrbMin > 0){
239 // cut required if the cut value is set to something larger than 0
240 // same effect as setting the min. cut value to 0
241 Int_t nCtrb = vtxPrim->GetNContributors();
242 if (nCtrb<fVtxNCtrbMin)
243 fBitMap->SetBitNumber(kVtxNCtrb,kFALSE);
246 // check vertex correlation cut
247 if(fCheckCorrelationSPDVtx){
248 if(vtxTracks && vtxTracks->GetNContributors() && vtxSPD && vtxSPD->GetNContributors()){
249 // Both vertices available, check correlation
250 if(TMath::Abs(vtxTracks->GetZ() - vtxSPD->GetZ()) >= 0.5) fBitMap->SetBitNumber(kCorrelation, kFALSE);
252 // No correlation available: set the cut to false
253 fBitMap->SetBitNumber(kCorrelation, kFALSE);
256 // check cut on the SPD vertex resolution
259 TString vtxTyp = vtxSPD->GetTitle();
261 vtxSPD->GetCovarianceMatrix(cov);
262 Double_t zRes = TMath::Sqrt(cov[5]);
263 if (vtxTyp.Contains("vertexer:Z") && (zRes>0.25)) fBitMap->SetBitNumber(kResolution, kFALSE);
265 fBitMap->SetBitNumber(kResolution, kFALSE);
268 // check pA pileup cut
270 if(fAnalysisUtils->IsPileUpEvent(inputEvent)) fBitMap->SetBitNumber(kpApileup,kFALSE);
274 //_____________________________________________________________________________
275 void AliHFEextraEventCuts::FillHistograms(TObject* obj, Bool_t b)
278 // fill the QA histograms
282 // cast TObject into VEvent
283 AliVEvent* inputEvent = dynamic_cast<AliESDEvent *>(obj);
284 if (!inputEvent) return;
286 // index = 0: fill histograms before cuts
287 // index = 1: fill histograms after cuts
289 index = ((b) ? 1 : 0);
292 const AliVVertex *vtxTracks = GetPrimaryVertexTracks(inputEvent),
293 *vtxSPD = GetPrimaryVertexSPD(inputEvent),
296 //look at vertex parameters:
298 if(vtxTracks && vtxTracks->GetNContributors() > 0) vtxPrim = vtxTracks;
299 else if(vtxSPD && vtxSPD->GetNContributors()) vtxPrim = vtxSPD;
305 // vertex position and uncertainties
306 fhQA[kVtxPosZ] [index]->Fill(vtxPrim->GetZ());
307 fhQA[kVtxNCtrb][index]->Fill(vtxPrim->GetNContributors());
308 // SPD Vertex Position Correlation
309 if(vtxTracks && vtxSPD){
310 fhQA[kCorrelation][index]->Fill(vtxTracks->GetZ()-vtxSPD->GetZ());
314 vtxSPD->GetCovarianceMatrix(cov);
315 fhQA[kResolution][index]->Fill(TMath::Sqrt(cov[5]));
318 if(fAnalysisUtils->IsPileUpEvent(inputEvent)) isPileup = 1,
319 fhQA[kpApileup][index]->Fill(isPileup);
322 //____________________________________________________________________
323 void AliHFEextraEventCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t *bins)
326 //setting x-axis bin limits of QA histogram fhQA[index]
329 for(Int_t i=0;i<kNStepQA;i++){
330 if(!fhQA[index][i]){AliWarning("non-existing histogram!");
333 fhQA[index][i]->GetXaxis()->Set(nbins,bins);
336 //____________________________________________________________________
337 void AliHFEextraEventCuts::SetHistogramBins(Int_t index, Int_t nbins, Double_t xmin, Double_t xmax)
340 //setting x-axis bins and range of QA histogram fhQA[index]
343 for(Int_t i=0;i<kNStepQA;i++){
344 if(!fhQA[index][i]){AliWarning("non-existing histogram!");
347 fhQA[index][i]->GetXaxis()->Set(nbins,xmin,xmax);
351 //_____________________________________________________________________________
352 void AliHFEextraEventCuts::DefineHistograms() {
354 // histograms for cut variables
359 AliInfo(Form("No QA histos requested, Please first set the QA flag on!"));
363 // book QA histograms
365 for (Int_t i=0; i<kNStepQA; i++) {
366 if (i==0) snprintf(str,5," ");
367 else snprintf(str,5,"_cut");
369 fhQA[kVtxPosZ][i] = new TH1F(Form("%s_Vtx_Pos_Z%s",GetName(),str), "",200,-50.,50.);
370 fhQA[kVtxNCtrb][i] = new TH1F(Form("%s_Vtx_N_Ctrb%s",GetName(),str), "",1000,0.,1000.);
371 fhQA[kCorrelation][i] = new TH1F(Form("%s_SPDCorrelation_%s",GetName(),str), "",200, -10., 10);
372 fhQA[kResolution][i] = new TH1F(Form("%s_SPDResolution_%s",GetName(), str), "", 100, 0., 1.);
373 fhQA[kpApileup][i] = new TH1F(Form("%s_pApileup_%s",GetName(),str), "", 2, 0., 2.);
375 fhQA[kVtxPosZ][i] ->SetXTitle("Vertex Position Z (cm)");
376 fhQA[kVtxNCtrb][i] ->SetXTitle("Number of contributors");
379 for(Int_t i=0; i<kNCuts; i++) fhQA[i][1]->SetLineColor(color);
383 //_____________________________________________________________________________
384 void AliHFEextraEventCuts::AddQAHistograms(TList *qaList) {
386 // saves the histograms in a TList
391 for (Int_t j=0; j<kNStepQA; j++) {
392 for(Int_t i=0; i<kNCuts; i++)
393 qaList->Add(fhQA[i][j]);
397 //_____________________________________________________________________________
398 const AliVVertex *AliHFEextraEventCuts::GetPrimaryVertexSPD(const AliVEvent * const inputEvent){
400 // Get SPD vertex from event
402 const AliVVertex *spdvtx(NULL);
403 const AliESDEvent *esd(NULL);
404 const AliAODEvent *aod(NULL);
405 if((esd = dynamic_cast<const AliESDEvent *>(inputEvent))){
406 spdvtx = esd->GetPrimaryVertexSPD();
407 } else if((aod = dynamic_cast<const AliAODEvent *>(inputEvent))){
408 spdvtx = aod->GetPrimaryVertexSPD();
413 //_____________________________________________________________________________
414 const AliVVertex *AliHFEextraEventCuts::GetPrimaryVertexTracks(const AliVEvent *const inputEvent){
416 // Get Primary Vertex from tracks
418 const AliVVertex *trkvtx(NULL);
419 const AliESDEvent *esd(NULL);
420 const AliAODEvent *aod(NULL);
421 if((esd = dynamic_cast<const AliESDEvent *>(inputEvent))){
422 trkvtx = esd->GetPrimaryVertexTracks();
423 } else if((aod = dynamic_cast<const AliAODEvent *>(inputEvent))){
424 const AliVVertex *vtxTmp = aod->GetPrimaryVertex();
425 // check whether the primary vertex is the vertex from tracks
426 TString vtxTtl = vtxTmp->GetTitle();
427 if(vtxTtl.Contains("VertexerTracks")){