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 **************************************************************************/
18 #include "AliESDMuonTrackCuts.h"
20 #include <AliESDMuonTrack.h>
22 #include <AliESDEvent.h>
28 #include <TDirectory.h>
30 //____________________________________________________________________
31 ClassImp(AliESDMuonTrackCuts)
34 const Char_t* AliESDMuonTrackCuts::fgkCutNames[kNCuts] = {
44 //____________________________________________________________________
45 AliESDMuonTrackCuts::AliESDMuonTrackCuts(const Char_t* name, const Char_t* title) : AliAnalysisCuts(name,title),
69 //##############################################################################
70 // setting default cuts
82 //_____________________________________________________________________________
83 AliESDMuonTrackCuts::AliESDMuonTrackCuts(const AliESDMuonTrackCuts &c) : AliAnalysisCuts(c),
105 ((AliESDMuonTrackCuts &) c).Copy(*this);
108 AliESDMuonTrackCuts::~AliESDMuonTrackCuts()
113 for (Int_t i=0; i<2; i++) {
121 delete fhCutStatistics;
122 if (fhCutCorrelation)
123 delete fhCutCorrelation;
126 void AliESDMuonTrackCuts::Init()
129 // sets everything to zero
146 fHistogramsOn = kFALSE;
148 for (Int_t i=0; i<2; ++i)
154 fhCutCorrelation = 0;
157 //_____________________________________________________________________________
158 AliESDMuonTrackCuts &AliESDMuonTrackCuts::operator=(const AliESDMuonTrackCuts &c)
161 // Assignment operator
163 if (this != &c) ((AliESDMuonTrackCuts &) c).Copy(*this);
167 //_____________________________________________________________________________
168 void AliESDMuonTrackCuts::Copy(TObject &c) const
174 AliESDMuonTrackCuts& target = (AliESDMuonTrackCuts &) c;
178 target.fPMin = fPMin;
179 target.fPMax = fPMax;
180 target.fPtMin = fPtMin;
181 target.fPtMax = fPtMax;
182 target.fPxMin = fPxMin;
183 target.fPxMax = fPxMax;
184 target.fPyMin = fPyMin;
185 target.fPyMax = fPyMax;
186 target.fPzMin = fPzMin;
187 target.fPzMax = fPzMax;
188 target.fEtaMin = fEtaMin;
189 target.fEtaMax = fEtaMax;
190 target.fRapMin = fRapMin;
191 target.fRapMax = fRapMax;
193 target.fHistogramsOn = fHistogramsOn;
195 for (Int_t i=0; i<2; ++i)
197 if (fhPt[i]) target.fhPt[i] = (TH1F*) fhPt[i]->Clone();
198 if (fhEta[i]) target.fhEta[i] = (TH1F*) fhEta[i]->Clone();
201 if (fhCutStatistics) target.fhCutStatistics = (TH1F*) fhCutStatistics->Clone();
202 if (fhCutCorrelation) target.fhCutCorrelation = (TH2F*) fhCutCorrelation->Clone();
207 //_____________________________________________________________________________
208 Long64_t AliESDMuonTrackCuts::Merge(TCollection* list) {
209 // Merge a list of AliESDMuonTrackCuts objects with this (needed for PROOF)
210 // Returns the number of merged objects (including this)
221 TIterator* iter = list->MakeIterator();
224 // collection of measured and generated histograms
226 while ((obj = iter->Next())) {
228 AliESDMuonTrackCuts* entry = dynamic_cast<AliESDMuonTrackCuts*>(obj);
232 if (!entry->fHistogramsOn)
235 for (Int_t i=0; i<2; i++) {
236 fhPt[i]->Add(entry->fhPt[i]);
237 fhEta[i]->Add(entry->fhEta[i]);
240 fhCutStatistics->Add(entry->fhCutStatistics);
241 fhCutCorrelation ->Add(entry->fhCutCorrelation);
249 void AliESDMuonTrackCuts::EnableNeededBranches(TTree* tree)
251 // enables the branches needed by AcceptTrack, for a list see comment of AcceptTrack
253 tree->SetBranchStatus("fTracks.fFlags", 1);
254 tree->SetBranchStatus("fTracks.fP*", 1);
255 tree->SetBranchStatus("fTracks.fR*", 1); //detector response probability
258 //____________________________________________________________________
260 AliESDMuonTrackCuts::AcceptTrack(AliESDMuonTrack* esdMuTrack) {
262 // figure out if the tracks survives all the track cuts defined
264 // the different kinematic values are first
265 // retrieved from the track. then it is found out what cuts the
266 // track did not survive and finally the cuts are imposed.
268 // getting the kinematic variables of the track
269 // (assuming the mass is known)
271 esdMuTrack->PxPyPz(p);
272 Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
273 Float_t pt = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2));
274 Float_t energy = TMath::Sqrt(TMath::Power(esdMuTrack->M(),2) + TMath::Power(momentum,2));
277 //y-eta related calculations
280 if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
281 eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
282 if((energy != TMath::Abs(p[2]))&&(momentum != 0))
283 y = 0.5*TMath::Log((energy + p[2])/(energy - p[2]));
285 //########################################################################
289 for (Int_t i=0; i<kNCuts; i++) cuts[i]=kFALSE;
291 // track kinematics cut
292 if((momentum < fPMin) || (momentum > fPMax))
294 if((pt < fPtMin) || (pt > fPtMax))
296 if((p[0] < fPxMin) || (p[0] > fPxMax))
298 if((p[1] < fPyMin) || (p[1] > fPyMax))
300 if((p[2] < fPzMin) || (p[2] > fPzMax))
302 if((eta < fEtaMin) || (eta > fEtaMax))
304 if((y < fRapMin) || (y > fRapMax))
308 for (Int_t i=0; i<kNCuts; i++)
309 if (cuts[i]) cut = kTRUE;
311 //########################################################################
312 // filling histograms
314 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n tracks")));
317 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin("n cut tracks")));
319 for (Int_t i=0; i<kNCuts; i++) {
321 fhCutStatistics->Fill(fhCutStatistics->GetBinCenter(fhCutStatistics->GetXaxis()->FindBin(fgkCutNames[i])));
323 for (Int_t j=i; j<kNCuts; j++) {
324 if (cuts[i] && cuts[j]) {
325 Float_t xC = fhCutCorrelation->GetXaxis()->GetBinCenter(fhCutCorrelation->GetXaxis()->FindBin(fgkCutNames[i]));
326 Float_t yC = fhCutCorrelation->GetYaxis()->GetBinCenter(fhCutCorrelation->GetYaxis()->FindBin(fgkCutNames[j]));
327 fhCutCorrelation->Fill(xC, yC);
337 //########################################################################
339 if (cut) return kFALSE;
341 //########################################################################
342 // filling histograms after cut
353 //____________________________________________________________________
354 TObjArray* AliESDMuonTrackCuts::GetAcceptedTracks(AliESD* esd)
357 // returns an array of all tracks that pass the cuts
360 TObjArray* acceptedTracks = new TObjArray();
362 // loop over esd tracks
363 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
364 AliESDMuonTrack* track = esd->GetMuonTrack(iTrack);
366 if (AcceptTrack(track))
367 acceptedTracks->Add(track);
370 return acceptedTracks;
373 //____________________________________________________________________
374 Int_t AliESDMuonTrackCuts::CountAcceptedTracks(AliESD* esd)
377 // returns an the number of tracks that pass the cuts
382 // loop over esd tracks
383 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
384 AliESDMuonTrack* track = esd->GetMuonTrack(iTrack);
386 if (AcceptTrack(track))
393 //____________________________________________________________________
394 TObjArray* AliESDMuonTrackCuts::GetAcceptedTracks(AliESDEvent* esd)
397 // returns an array of all tracks that pass the cuts
400 TObjArray* acceptedTracks = new TObjArray();
402 // loop over esd tracks
403 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
404 AliESDMuonTrack* track = esd->GetMuonTrack(iTrack);
406 if (AcceptTrack(track))
407 acceptedTracks->Add(track);
410 return acceptedTracks;
413 //____________________________________________________________________
414 Int_t AliESDMuonTrackCuts::CountAcceptedTracks(AliESDEvent* esd)
417 // returns an the number of tracks that pass the cuts
422 // loop over esd tracks
423 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
424 AliESDMuonTrack* track = esd->GetMuonTrack(iTrack);
426 if (AcceptTrack(track))
433 //____________________________________________________________________
434 void AliESDMuonTrackCuts::DefineHistograms(Int_t color) {
436 // diagnostics histograms are defined
441 Bool_t oldStatus = TH1::AddDirectoryStatus();
442 TH1::AddDirectory(kFALSE);
444 //###################################################################################
445 // defining histograms
447 fhCutStatistics = new TH1F("cut_statistics","cut statistics",kNCuts+4,-0.5,kNCuts+3.5);
449 fhCutStatistics->GetXaxis()->SetBinLabel(1,"n tracks");
450 fhCutStatistics->GetXaxis()->SetBinLabel(2,"n cut tracks");
452 fhCutCorrelation = new TH2F("cut_correlation","cut correlation",kNCuts,-0.5,kNCuts-0.5,kNCuts,-0.5,kNCuts-0.5);;
454 for (Int_t i=0; i<kNCuts; i++) {
455 fhCutStatistics->GetXaxis()->SetBinLabel(i+4,fgkCutNames[i]);
456 fhCutCorrelation->GetXaxis()->SetBinLabel(i+1,fgkCutNames[i]);
457 fhCutCorrelation->GetYaxis()->SetBinLabel(i+1,fgkCutNames[i]);
460 fhCutStatistics ->SetLineColor(color);
461 fhCutCorrelation ->SetLineColor(color);
462 fhCutStatistics ->SetLineWidth(2);
463 fhCutCorrelation ->SetLineWidth(2);
466 for (Int_t i=0; i<2; i++) {
467 if (i==0) snprintf(str,256," ");
468 else snprintf(str,256,"_cut");
470 fhPt[i] = new TH1F(Form("pt%s",str) ,"p_{T} distribution;p_{T} (GeV/c)",500,0.0,100.0);
471 fhEta[i] = new TH1F(Form("eta%s",str) ,"#eta distribution;#eta",40,-2.0,2.0);
474 TH1::AddDirectory(oldStatus);
477 //____________________________________________________________________
478 Bool_t AliESDMuonTrackCuts::LoadHistograms(const Char_t* dir)
481 // loads the histograms from a file
482 // if dir is empty a directory with the name of this object is taken (like in SaveHistogram)
488 if (!gDirectory->cd(dir))
491 fhCutStatistics = dynamic_cast<TH1F*> (gDirectory->Get("cut_statistics"));
492 fhCutCorrelation = dynamic_cast<TH2F*> (gDirectory->Get("cut_correlation"));
495 for (Int_t i=0; i<2; i++) {
498 gDirectory->cd("before_cuts");
503 gDirectory->cd("after_cuts");
504 snprintf(str,5,"_cut");
507 fhPt[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("pt%s",str)));
508 fhEta[i] = dynamic_cast<TH1F*> (gDirectory->Get(Form("eta%s",str)));
510 gDirectory->cd("../");
513 gDirectory->cd("..");
518 //____________________________________________________________________
519 void AliESDMuonTrackCuts::SaveHistograms(const Char_t* dir) {
521 // saves the histograms in a directory (dir)
524 if (!fHistogramsOn) {
525 AliDebug(0, "Histograms not on - cannot save histograms!!!");
532 gDirectory->mkdir(dir);
535 gDirectory->mkdir("before_cuts");
536 gDirectory->mkdir("after_cuts");
538 fhCutStatistics->Write();
539 fhCutCorrelation->Write();
541 for (Int_t i=0; i<2; i++) {
543 gDirectory->cd("before_cuts");
545 gDirectory->cd("after_cuts");
550 gDirectory->cd("../");
553 gDirectory->cd("../");
556 //____________________________________________________________________
557 void AliESDMuonTrackCuts::DrawHistograms()
559 gStyle->SetPalette(1);
560 gStyle->SetFrameFillColor(10);
561 gStyle->SetCanvasColor(10);
563 TCanvas* canvas1 = new TCanvas(Form("%s_1", GetName()), "Track Cut Results", 800, 500);
564 canvas1->Divide(2, 1);
567 fhCutStatistics->SetStats(kFALSE);
568 fhCutStatistics->LabelsOption("v");
569 gPad->SetBottomMargin(0.3);
570 fhCutStatistics->Draw();
573 fhCutCorrelation->SetStats(kFALSE);
574 fhCutCorrelation->LabelsOption("v");
575 gPad->SetBottomMargin(0.3);
576 gPad->SetLeftMargin(0.3);
577 fhCutCorrelation->Draw("COLZ");
579 canvas1->SaveAs(Form("%s_%s.gif", GetName(), canvas1->GetName()));