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 **************************************************************************/
16 #include "AliESDMuonTrackCuts.h"
18 #include <AliESDMuonTrack.h>
20 #include <AliESDEvent.h>
26 #include <TDirectory.h>
28 //____________________________________________________________________
29 ClassImp(AliESDMuonTrackCuts)
32 const Char_t* AliESDMuonTrackCuts::fgkCutNames[kNCuts] = {
42 //____________________________________________________________________
43 AliESDMuonTrackCuts::AliESDMuonTrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title),
67 //##############################################################################
68 // setting default cuts
80 //_____________________________________________________________________________
81 AliESDMuonTrackCuts::AliESDMuonTrackCuts(const AliESDMuonTrackCuts &c) : AliAnalysisCuts(c),
103 ((AliESDMuonTrackCuts &) c).Copy(*this);
106 AliESDMuonTrackCuts::~AliESDMuonTrackCuts()
111 for (Int_t i=0; i<2; i++) {
119 delete fhCutStatistics;
120 if (fhCutCorrelation)
121 delete fhCutCorrelation;
124 void AliESDMuonTrackCuts::Init()
127 // sets everything to zero
144 fHistogramsOn = kFALSE;
146 for (Int_t i=0; i<2; ++i)
152 fhCutCorrelation = 0;
155 //_____________________________________________________________________________
156 AliESDMuonTrackCuts &AliESDMuonTrackCuts::operator=(const AliESDMuonTrackCuts &c)
159 // Assignment operator
161 if (this != &c) ((AliESDMuonTrackCuts &) c).Copy(*this);
165 //_____________________________________________________________________________
166 void AliESDMuonTrackCuts::Copy(TObject &c) const
172 AliESDMuonTrackCuts& target = (AliESDMuonTrackCuts &) c;
176 target.fPMin = fPMin;
177 target.fPMax = fPMax;
178 target.fPtMin = fPtMin;
179 target.fPtMax = fPtMax;
180 target.fPxMin = fPxMin;
181 target.fPxMax = fPxMax;
182 target.fPyMin = fPyMin;
183 target.fPyMax = fPyMax;
184 target.fPzMin = fPzMin;
185 target.fPzMax = fPzMax;
186 target.fEtaMin = fEtaMin;
187 target.fEtaMax = fEtaMax;
188 target.fRapMin = fRapMin;
189 target.fRapMax = fRapMax;
191 target.fHistogramsOn = fHistogramsOn;
193 for (Int_t i=0; i<2; ++i)
195 if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone();
196 if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone();
199 if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
200 if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
205 //_____________________________________________________________________________
206 Long64_t AliESDMuonTrackCuts::Merge(TCollection* list) {
207 // Merge a list of AliESDMuonTrackCuts objects with this (needed for PROOF)
208 // Returns the number of merged objects (including this)
219 TIterator* iter = list->MakeIterator();
222 // collection of measured and generated histograms
224 while ((obj = iter->Next())) {
226 AliESDMuonTrackCuts* entry = dynamic_cast<AliESDMuonTrackCuts*>(obj);
230 if (!entry->fHistogramsOn)
233 for (Int_t i=0; i<2; i++) {
234 fhPt[i]->Add(entry->fhPt[i]);
235 fhEta[i]->Add(entry->fhEta[i]);
238 fhCutStatistics->Add(entry->fhCutStatistics);
239 fhCutCorrelation ->Add(entry->fhCutCorrelation);
247 void AliESDMuonTrackCuts::EnableNeededBranches(TTree* tree)
249 // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack
251 tree->SetBranchStatus("fTracks.fFlags", 1);
252 tree->SetBranchStatus("fTracks.fP*", 1);
253 tree->SetBranchStatus("fTracks.fR*", 1); //detector response probability
256 //____________________________________________________________________
258 AliESDMuonTrackCuts::AcceptTrack(AliESDMuonTrack* esdMuTrack) {
260 // figure out if the tracks survives all the track cuts defined
262 // the different kinematic values are first
263 // retrieved from the track. then it is found out what cuts the
264 // track did not survive and finally the cuts are imposed.
266 // getting the kinematic variables of the track
267 // (assuming the mass is known)
269 esdMuTrack->PxPyPz(p);
270 Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
271 Float_t pt = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
272 Float_t energy = TMath::Sqrt(TMath::Power(esdMuTrack->M(),2) + TMath::Power(momentum,2));
275 //y-eta related calculations
278 if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
279 eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
280 if((energy != TMath::Abs(p[2]))&&(momentum != 0))
281 y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
283 //########################################################################
287 for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
289 // track kinematics cut
290 if((momentum < fPMin) || (momentum > fPMax))
292 if((pt < fPtMin) || (pt > fPtMax))
294 if((p[0] < fPxMin) || (p[0] > fPxMax))
296 if((p[1] < fPyMin) || (p[1] > fPyMax))
298 if((p[2] < fPzMin) || (p[2] > fPzMax))
300 if((eta < fEtaMin) || (eta > fEtaMax))
302 if((y < fRapMin) || (y > fRapMax))
306 for (Int_t i=0; i<kNCuts; i++)
307 if (cuts[i]) cut = kTRUE;
309 //########################################################################
310 // filling histograms
312 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
315 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
317 for (Int_t i=0; i<kNCuts; i++) {
319 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
321 for (Int_t j=i; j<kNCuts; j++) {
322 if (cuts[i] && cuts[j]) {
323 Float_t xC = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
324 Float_t yC = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
325 fhCutCorrelation->Fill(xC, yC);
335 //########################################################################
337 if (cut) return kFALSE;
339 //########################################################################
340 // filling histograms after cut
351 //____________________________________________________________________
352 TObjArray* AliESDMuonTrackCuts::GetAcceptedTracks(AliESD* esd)
355 // returns an array of all tracks that pass the cuts
358 TObjArray* acceptedTracks = new TObjArray();
360 // loop over esd tracks
361 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
362 AliESDMuonTrack* track = esd->GetMuonTrack(iTrack);
364 if (AcceptTrack(track))
365 acceptedTracks->Add(track);
368 return acceptedTracks;
371 //____________________________________________________________________
372 Int_t AliESDMuonTrackCuts::CountAcceptedTracks(AliESD* esd)
375 // returns an the number of tracks that pass the cuts
380 // loop over esd tracks
381 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
382 AliESDMuonTrack* track = esd->GetMuonTrack(iTrack);
384 if (AcceptTrack(track))
391 //____________________________________________________________________
392 TObjArray* AliESDMuonTrackCuts::GetAcceptedTracks(AliESDEvent* esd)
395 // returns an array of all tracks that pass the cuts
398 TObjArray* acceptedTracks = new TObjArray();
400 // loop over esd tracks
401 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
402 AliESDMuonTrack* track = esd->GetMuonTrack(iTrack);
404 if (AcceptTrack(track))
405 acceptedTracks->Add(track);
408 return acceptedTracks;
411 //____________________________________________________________________
412 Int_t AliESDMuonTrackCuts::CountAcceptedTracks(AliESDEvent* esd)
415 // returns an the number of tracks that pass the cuts
420 // loop over esd tracks
421 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
422 AliESDMuonTrack* track = esd->GetMuonTrack(iTrack);
424 if (AcceptTrack(track))
431 //____________________________________________________________________
432 void AliESDMuonTrackCuts::DefineHistograms(Int_t color) {
434 // diagnostics histograms are defined
439 Bool_t oldStatus = TH1::AddDirectoryStatus();
440 TH1::AddDirectory(kFALSE);
442 //###################################################################################
443 // defining histograms
445 fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
447 fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
448 fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
450 fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
452 for (Int_t i=0; i<kNCuts; i++) {
453 fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
454 fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
455 fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
458 fhCutStatistics ->SetLineColor(color);
459 fhCutCorrelation ->SetLineColor(color);
460 fhCutStatistics ->SetLineWidth(2);
461 fhCutCorrelation ->SetLineWidth(2);
464 for (Int_t i=0; i<2; i++) {
465 if (i==0) sprintf(str," ");
466 else sprintf(str,"_cut");
468 fhPt[i] = new TH1F(Form("pt%s",str) ,"p_{T} distribution;p_{T} (GeV/c)",500,0.0,100.0);
469 fhEta[i] = new TH1F(Form("eta%s",str) ,"#eta distribution;#eta",40,-2.0,2.0);
472 TH1::AddDirectory(oldStatus);
475 //____________________________________________________________________
476 Bool_t AliESDMuonTrackCuts::LoadHistograms(const Char_t* dir)
479 // loads the histograms from a file
480 // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
486 if (!gDirectory->cd(dir))
489 fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
490 fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
493 for (Int_t i=0; i<2; i++) {
496 gDirectory->cd("before_cuts");
501 gDirectory->cd("after_cuts");
505 fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("pt%s",str)));
506 fhEta[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("eta%s",str)));
508 gDirectory->cd("../");
511 gDirectory->cd("..");
516 //____________________________________________________________________
517 void AliESDMuonTrackCuts::SaveHistograms(const Char_t* dir) {
519 // saves the histograms in a directory (dir)
522 if (!fHistogramsOn) {
523 AliDebug(0, "Histograms not on - cannot save histograms!!!");
530 gDirectory->mkdir(dir);
533 gDirectory->mkdir("before_cuts");
534 gDirectory->mkdir("after_cuts");
536 fhCutStatistics->Write();
537 fhCutCorrelation->Write();
539 for (Int_t i=0; i<2; i++) {
541 gDirectory->cd("before_cuts");
543 gDirectory->cd("after_cuts");
548 gDirectory->cd("../");
551 gDirectory->cd("../");
554 //____________________________________________________________________
555 void AliESDMuonTrackCuts::DrawHistograms()
557 gStyle->SetPalette(1);
558 gStyle->SetFrameFillColor(10);
559 gStyle->SetCanvasColor(10);
561 TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Cut Results", 800, 500);
562 canvas1->Divide(2, 1);
565 fhCutStatistics->SetStats(kFALSE);
566 fhCutStatistics->LabelsOption("v");
567 gPad->SetBottomMargin(0.3);
568 fhCutStatistics->Draw();
571 fhCutCorrelation->SetStats(kFALSE);
572 fhCutCorrelation->LabelsOption("v");
573 gPad->SetBottomMargin(0.3);
574 gPad->SetLeftMargin(0.3);
575 fhCutCorrelation->Draw("COLZ");
577 canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));