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 /// An analysis task to check the MUON data in simulated data
19 /// This class checks out the ESD tree, providing the matching with
20 /// the trigger,trigger responses for Low and High Pt cuts
21 /// (in Single, Unlike Sign and Like Sign) and gives Pt, Y, ITS vertex
22 /// and multiplicity distributions. All results are in histogram form.
23 /// The output is a root file and eps files in MUON.tar.gz.
25 //*-- Frederic Yermia, yermia@to.infn.it
26 //////////////////////////////////////////////////////////////////////////////
27 //////////////////////////////////////////////////////////////////////////////
34 #include <TLorentzVector.h>
35 #include "AliMUONQATask.h"
38 #include "AliESDVertex.h"
39 #include "AliESDMuonTrack.h"
40 #include <TLorentzVector.h>
41 //______________________________________________________________________________
42 AliMUONQATask::AliMUONQATask(const char *name) :
43 AliAnalysisTask(name,""),
54 fmuonMass(0.105658389),
77 // Input slot #0 works with an Ntuple
78 DefineInput(0, TChain::Class());
79 // Output slot #0 writes into a TH1 container
80 DefineOutput(0, TObjArray::Class()) ;
85 //______________________________________________________________________________
86 AliMUONQATask::~AliMUONQATask()
89 fOutputContainer->Clear() ;
90 delete fOutputContainer ;
106 //______________________________________________________________________________
107 void AliMUONQATask::ConnectInputData(const Option_t*)
109 // Initialisation of branch container and histograms
111 AliInfo(Form("*** Initialization of %s", GetName())) ;
114 fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
116 AliError(Form("Input 0 for %s not found\n", GetName()));
120 // One should first check if the branch address was taken by some other task
121 char ** address = (char **)GetBranchAddress(0, "ESD");
123 fESD = (AliESD*)(*address);
126 SetBranchAddress(0, "ESD", &fESD);
130 //________________________________________________________________________
131 void AliMUONQATask::CreateOutputObjects()
134 fhMUONVertex = new TH1F("hMUONVertex","ITS Vertex" ,100, -25., 25.);
135 fhMUONMult = new TH1F("hMUONMult" ,"Multiplicity of ESD tracks",10, -0.5, 9.5);
136 fhPt = new TH1F("hPt","Pt",100, 0.,20.);
137 fhY = new TH1F("hY","Rapidity",100,-5.,-1.);
138 fheffMatchT = new TH1F("heff_matchT","Trigger Matching Efficiency",100, 0.,100.);
139 fhSLowpt = new TH1F("hSLowpt","Single Low Pt Response (%)",101, 0.,101.);
140 fhUSLowpt = new TH1F("hUSLowpt","Unlike Sign Low Pt Response (%)",101, 0.,101.);
141 fhUSHighpt = new TH1F("hUSHighpt","Unlike Sign High Pt Response (%)",101, 0.,101.);
142 fhLSLowpt = new TH1F("hLSLowpt","Like Sign Low Pt Response (%)",101, 0.,101.);
143 fhLSHighpt = new TH1F("hLSHighpt","Like Sign High Pt Response (%)",101, 0.,101.);
144 fhChi2 = new TH1F("hChi2","Chi2 by d.o.f.",100, 0.,20.);
145 fhChi2match = new TH1F("hChi2match","Chi2 of trig/track matching",100, 0.,20.);
146 // create output container
148 fOutputContainer = new TObjArray(12) ;
149 fOutputContainer->SetName(GetName()) ;
150 fOutputContainer->AddAt(fhMUONVertex, 0) ;
151 fOutputContainer->AddAt(fhMUONMult, 1) ;
152 fOutputContainer->AddAt(fhPt, 2) ;
153 fOutputContainer->AddAt(fhY, 3) ;
154 fOutputContainer->AddAt(fheffMatchT, 4) ;
155 fOutputContainer->AddAt(fhSLowpt, 5) ;
156 fOutputContainer->AddAt(fhUSLowpt, 6) ;
157 fOutputContainer->AddAt(fhUSHighpt, 7) ;
158 fOutputContainer->AddAt(fhLSLowpt, 8) ;
159 fOutputContainer->AddAt(fhLSHighpt, 9) ;
160 fOutputContainer->AddAt(fhChi2, 10) ;
161 fOutputContainer->AddAt( fhChi2match, 11) ;
164 //______________________________________________________________________________
165 void AliMUONQATask::Exec(Option_t *)
167 // Processing of one event
171 Long64_t entry = fChain->GetReadEntry() ;
174 AliError("fESD is not connected to the input!") ;
178 if ( !((entry-1)%100) )
179 AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ;
181 // ************************ MUON *************************************
183 const AliESDVertex* vertex = dynamic_cast<const AliESDVertex*>(fESD->GetVertex()) ;
185 Double_t zVertex = 0. ;
187 zVertex = vertex->GetZv() ;
189 Int_t nTracks = fESD->GetNumberOfMuonTracks() ;
191 ULong64_t trigWord = fESD->GetTriggerMask() ;
193 if (trigWord & 0x80) {
196 if (trigWord & 0x100){
199 if (trigWord & 0x200){
202 if (trigWord & 0x400){
205 if (trigWord & 0x800){
209 Int_t tracktrig = 0 ;
212 for (iTrack1 = 0 ; iTrack1 < nTracks ; iTrack1++) { //1st loop
213 AliESDMuonTrack* muonTrack = fESD->GetMuonTrack(iTrack1) ;
215 fthetaX = muonTrack->GetThetaX();
216 fthetaY = muonTrack->GetThetaY();
217 fpYZ = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
218 fPzRec1 = - fpYZ / TMath::Sqrt(1.0 + TMath::Tan(fthetaY)*TMath::Tan(fthetaY));
219 fPxRec1 = fPzRec1 * TMath::Tan(fthetaX);
220 fPyRec1 = fPzRec1 * TMath::Tan(fthetaY);
221 fZ1 = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
222 fE1 = TMath::Sqrt(fmuonMass * fmuonMass + fPxRec1 * fPxRec1 + fPyRec1 * fPyRec1 + fPzRec1 * fPzRec1);
223 fV1.SetPxPyPzE(fPxRec1, fPyRec1, fPzRec1, fE1);
225 // -----------> transverse momentum
226 Float_t pt1 = fV1.Pt();
227 // ----------->Rapidity
228 Float_t y1 = fV1.Rapidity();
230 if(muonTrack->GetMatchTrigger()) {
233 Float_t Chi2match = muonTrack->GetChi2MatchTrigger();
234 fhChi2match->Fill(Chi2match);
237 Float_t fitfmin = muonTrack->GetChi2();
238 Int_t ntrackhits = muonTrack->GetNHit();
239 Float_t Chi2= fitfmin / (2.0 * ntrackhits - 5);
246 fhMUONVertex->Fill(zVertex) ;
247 fhMUONMult->Fill(Float_t(nTracks)) ;
249 PostData(0, fOutputContainer);
252 //______________________________________________________________________________
253 void AliMUONQATask::Terminate(Option_t *)
255 // Processing when the event loop is ended
256 Int_t fSLowPt = fSLowpt;
258 fSLowPt = 100 * fSLowpt / fnevents ;
259 fhSLowpt->Fill(fSLowPt); }
260 Int_t fUSLowPt = fUSLowpt;
262 fUSLowPt = 100 * fUSLowpt / fnevents ;
263 fhUSLowpt->Fill(fUSLowPt); }
264 Int_t fUSHighPt = fUSHighpt;
266 fUSHighPt = 100 * fUSHighpt / fnevents ;
267 fhUSHighpt->Fill(fUSHighPt); }
268 Int_t fLSLowPt = fLSLowpt;
270 fLSLowPt = 100 * fLSLowpt / fnevents ;
271 fhLSLowpt->Fill(fLSLowPt); }
272 Int_t fLSHighPt = fLSHighpt;
274 fLSHighPt = 100 * fLSHighpt / fnevents ;
275 fhLSHighpt->Fill(fLSHighPt); }
277 Int_t effMatch = -1 ;
279 effMatch = 100 * fnTrackTrig / ftracktot ;
280 fheffMatchT->Fill(effMatch);}
282 AliInfo(Form(" *** %s Report:", GetName())) ;
284 fOutputContainer = (TObjArray*)GetOutputData(0);
285 fhMUONVertex = (TH1F*)fOutputContainer->At(0);
286 fhMUONMult = (TH1F*)fOutputContainer->At(1);
287 fhPt = (TH1F*)fOutputContainer->At(2);
288 fhY = (TH1F*)fOutputContainer->At(3);
289 fheffMatchT=(TH1F*)fOutputContainer->At(4);
290 fhSLowpt=(TH1F*)fOutputContainer->At(5);
291 fhUSLowpt=(TH1F*)fOutputContainer->At(6);
292 fhUSHighpt=(TH1F*)fOutputContainer->At(7);
293 fhLSLowpt=(TH1F*)fOutputContainer->At(8);
294 fhLSHighpt=(TH1F*)fOutputContainer->At(9);
296 printf(" Total number of processed events %d \n", fnevents) ;
299 printf(" Table 1: \n") ;
300 printf(" ===================================================\n") ;
301 printf(" Global Trigger output Low pt High pt \n") ;
302 printf(" number of Single :\t");
303 printf(" %i\t", fSLowpt) ;
305 printf(" number of UnlikeSign pair :\t");
306 printf(" %i\t%i\t", fUSLowpt, fUSHighpt) ;
308 printf(" number of LikeSign pair :\t");
309 printf(" %i\t%i\t", fLSLowpt, fLSHighpt) ;
311 printf(" matching efficiency with the trigger for single tracks = %2d %% \n", effMatch);
314 TCanvas * cMUON1 = new TCanvas("cMUON1", "MUON ESD Vert & Mult", 400, 10, 600, 700) ;
315 cMUON1->Divide(1,2) ;
317 fhMUONVertex->SetXTitle("Vz (cm)");
318 fhMUONVertex->Draw() ;
320 fhMUONMult->SetXTitle(" Track Multiplicity");
322 cMUON1->Print("MUON1.eps") ;
324 TCanvas * cMUON2 = new TCanvas("cMUON2", "MUON ESD Pt & Y", 400, 10, 600, 700) ;
325 cMUON2->Divide(1,2) ;
327 fhPt->SetXTitle("Pt (GeV)");
332 cMUON2->Print("MUON2.eps") ;
334 TCanvas * cMUON3 = new TCanvas("cMUON3", "Track Chi2 by dof and Chi2 of trigger/track matching ", 400, 10, 600, 700) ;
335 cMUON3->Divide(1,2) ;
337 fhChi2->SetXTitle("Chi2 by d.o.f.");
340 fhChi2match->SetXTitle("Chi2 of trig/track matching");
342 cMUON3->Print("MUON3.eps") ;
344 TCanvas * cMUON4 = new TCanvas("cMUON4", "Trigger Matching and Trigger Response (%)", 400, 10, 600, 700) ;
345 cMUON4->Divide(2,3) ;
347 fheffMatchT->SetXTitle("%");
348 fheffMatchT->Draw() ;
350 fhSLowpt->SetXTitle("%");
353 fhUSLowpt->SetXTitle("%");
356 fhUSHighpt->SetXTitle("%");
359 fhLSLowpt->SetXTitle("%");
362 fhLSHighpt->SetXTitle("%");
364 cMUON4->Print("MUON4.eps") ;
367 sprintf(line, ".!tar -zcf %s.tar.gz *.eps", GetName()) ;
368 gROOT->ProcessLine(line);
369 sprintf(line, ".!rm -fR *.eps");
370 gROOT->ProcessLine(line);
372 AliInfo(Form("!!! All the eps files are in %s.tar.gz !!! \n", GetName())) ;