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()
137 fhMUONVertex = new TH1F("hMUONVertex","ITS Vertex" ,100, -25., 25.);
138 fhMUONMult = new TH1F("hMUONMult" ,"Multiplicity of ESD tracks",10, -0.5, 9.5);
139 fhPt = new TH1F("hPt","Pt",100, 0.,20.);
140 fhY = new TH1F("hY","Rapidity",100,-5.,-1.);
141 fheffMatchT = new TH1F("heff_matchT","Trigger Matching Efficiency",100, 0.,100.);
142 fhSLowpt = new TH1F("hSLowpt","Single Low Pt Response (%)",101, 0.,101.);
143 fhUSLowpt = new TH1F("hUSLowpt","Unlike Sign Low Pt Response (%)",101, 0.,101.);
144 fhUSHighpt = new TH1F("hUSHighpt","Unlike Sign High Pt Response (%)",101, 0.,101.);
145 fhLSLowpt = new TH1F("hLSLowpt","Like Sign Low Pt Response (%)",101, 0.,101.);
146 fhLSHighpt = new TH1F("hLSHighpt","Like Sign High Pt Response (%)",101, 0.,101.);
147 fhChi2 = new TH1F("hChi2","Chi2 by d.o.f.",100, 0.,20.);
148 fhChi2match = new TH1F("hChi2match","Chi2 of trig/track matching",100, 0.,20.);
149 // create output container
151 fOutputContainer = new TObjArray(12) ;
152 fOutputContainer->SetName(GetName()) ;
153 fOutputContainer->AddAt(fhMUONVertex, 0) ;
154 fOutputContainer->AddAt(fhMUONMult, 1) ;
155 fOutputContainer->AddAt(fhPt, 2) ;
156 fOutputContainer->AddAt(fhY, 3) ;
157 fOutputContainer->AddAt(fheffMatchT, 4) ;
158 fOutputContainer->AddAt(fhSLowpt, 5) ;
159 fOutputContainer->AddAt(fhUSLowpt, 6) ;
160 fOutputContainer->AddAt(fhUSHighpt, 7) ;
161 fOutputContainer->AddAt(fhLSLowpt, 8) ;
162 fOutputContainer->AddAt(fhLSHighpt, 9) ;
163 fOutputContainer->AddAt(fhChi2, 10) ;
164 fOutputContainer->AddAt( fhChi2match, 11) ;
167 //______________________________________________________________________________
168 void AliMUONQATask::Exec(Option_t *)
170 // Processing of one event
174 Long64_t entry = fChain->GetReadEntry() ;
177 AliError("fESD is not connected to the input!") ;
181 if ( !((entry-1)%100) )
182 AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ;
184 // ************************ MUON *************************************
186 const AliESDVertex* vertex = dynamic_cast<const AliESDVertex*>(fESD->GetVertex()) ;
188 Double_t zVertex = 0. ;
190 zVertex = vertex->GetZv() ;
192 Int_t nTracks = fESD->GetNumberOfMuonTracks() ;
194 ULong64_t trigWord = fESD->GetTriggerMask() ;
196 if (trigWord & 0x80) {
199 if (trigWord & 0x100){
202 if (trigWord & 0x200){
205 if (trigWord & 0x400){
208 if (trigWord & 0x800){
212 Int_t tracktrig = 0 ;
215 for (iTrack1 = 0 ; iTrack1 < nTracks ; iTrack1++) { //1st loop
216 AliESDMuonTrack* muonTrack = fESD->GetMuonTrack(iTrack1) ;
218 fthetaX = muonTrack->GetThetaX();
219 fthetaY = muonTrack->GetThetaY();
220 fpYZ = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
221 fPzRec1 = - fpYZ / TMath::Sqrt(1.0 + TMath::Tan(fthetaY)*TMath::Tan(fthetaY));
222 fPxRec1 = fPzRec1 * TMath::Tan(fthetaX);
223 fPyRec1 = fPzRec1 * TMath::Tan(fthetaY);
224 fZ1 = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
225 fE1 = TMath::Sqrt(fmuonMass * fmuonMass + fPxRec1 * fPxRec1 + fPyRec1 * fPyRec1 + fPzRec1 * fPzRec1);
226 fV1.SetPxPyPzE(fPxRec1, fPyRec1, fPzRec1, fE1);
228 // -----------> transverse momentum
229 Float_t pt1 = fV1.Pt();
230 // ----------->Rapidity
231 Float_t y1 = fV1.Rapidity();
233 if(muonTrack->GetMatchTrigger()) {
236 Float_t Chi2match = muonTrack->GetChi2MatchTrigger();
237 fhChi2match->Fill(Chi2match);
240 Float_t fitfmin = muonTrack->GetChi2();
241 Int_t ntrackhits = muonTrack->GetNHit();
242 Float_t Chi2= fitfmin / (2.0 * ntrackhits - 5);
249 fhMUONVertex->Fill(zVertex) ;
250 fhMUONMult->Fill(Float_t(nTracks)) ;
252 PostData(0, fOutputContainer);
255 //______________________________________________________________________________
256 void AliMUONQATask::Terminate(Option_t *)
258 // Processing when the event loop is ended
259 Int_t fSLowPt = fSLowpt;
261 fSLowPt = 100 * fSLowpt / fnevents ;
262 fhSLowpt->Fill(fSLowPt); }
263 Int_t fUSLowPt = fUSLowpt;
265 fUSLowPt = 100 * fUSLowpt / fnevents ;
266 fhUSLowpt->Fill(fUSLowPt); }
267 Int_t fUSHighPt = fUSHighpt;
269 fUSHighPt = 100 * fUSHighpt / fnevents ;
270 fhUSHighpt->Fill(fUSHighPt); }
271 Int_t fLSLowPt = fLSLowpt;
273 fLSLowPt = 100 * fLSLowpt / fnevents ;
274 fhLSLowpt->Fill(fLSLowPt); }
275 Int_t fLSHighPt = fLSHighpt;
277 fLSHighPt = 100 * fLSHighpt / fnevents ;
278 fhLSHighpt->Fill(fLSHighPt); }
280 Int_t effMatch = -1 ;
282 effMatch = 100 * fnTrackTrig / ftracktot ;
283 fheffMatchT->Fill(effMatch);}
285 Bool_t problem = kFALSE ;
286 AliInfo(Form(" *** %s Report:", GetName())) ;
288 fOutputContainer = (TObjArray*)GetOutputData(0);
289 fhMUONVertex = (TH1F*)fOutputContainer->At(0);
290 fhMUONMult = (TH1F*)fOutputContainer->At(1);
291 fhPt = (TH1F*)fOutputContainer->At(2);
292 fhY = (TH1F*)fOutputContainer->At(3);
293 fheffMatchT=(TH1F*)fOutputContainer->At(4);
294 fhSLowpt=(TH1F*)fOutputContainer->At(5);
295 fhUSLowpt=(TH1F*)fOutputContainer->At(6);
296 fhUSHighpt=(TH1F*)fOutputContainer->At(7);
297 fhLSLowpt=(TH1F*)fOutputContainer->At(8);
298 fhLSHighpt=(TH1F*)fOutputContainer->At(9);
300 printf(" Total number of processed events %d \n", fnevents) ;
303 printf(" Table 1: \n") ;
304 printf(" ===================================================\n") ;
305 printf(" Global Trigger output Low pt High pt \n") ;
306 printf(" number of Single :\t");
307 printf(" %i\t", fSLowpt) ;
309 printf(" number of UnlikeSign pair :\t");
310 printf(" %i\t%i\t", fUSLowpt, fUSHighpt) ;
312 printf(" number of LikeSign pair :\t");
313 printf(" %i\t%i\t", fLSLowpt, fLSHighpt) ;
315 printf(" matching efficiency with the trigger for single tracks = %2d %% \n", effMatch);
318 TCanvas * cMUON1 = new TCanvas("cMUON1", "MUON ESD Vert & Mult", 400, 10, 600, 700) ;
319 cMUON1->Divide(1,2) ;
321 fhMUONVertex->SetXTitle("Vz (cm)");
322 fhMUONVertex->Draw() ;
324 fhMUONMult->SetXTitle(" Track Multiplicity");
326 cMUON1->Print("MUON1.eps") ;
328 TCanvas * cMUON2 = new TCanvas("cMUON2", "MUON ESD Pt & Y", 400, 10, 600, 700) ;
329 cMUON2->Divide(1,2) ;
331 fhPt->SetXTitle("Pt (GeV)");
336 cMUON2->Print("MUON2.eps") ;
338 TCanvas * cMUON3 = new TCanvas("cMUON3", "Track Chi2 by dof and Chi2 of trigger/track matching ", 400, 10, 600, 700) ;
339 cMUON3->Divide(1,2) ;
341 fhChi2->SetXTitle("Chi2 by d.o.f.");
344 fhChi2match->SetXTitle("Chi2 of trig/track matching");
346 cMUON3->Print("MUON3.eps") ;
348 TCanvas * cMUON4 = new TCanvas("cMUON4", "Trigger Matching and Trigger Response (%)", 400, 10, 600, 700) ;
349 cMUON4->Divide(2,3) ;
351 fheffMatchT->SetXTitle("%");
352 fheffMatchT->Draw() ;
354 fhSLowpt->SetXTitle("%");
357 fhUSLowpt->SetXTitle("%");
360 fhUSHighpt->SetXTitle("%");
363 fhLSLowpt->SetXTitle("%");
366 fhLSHighpt->SetXTitle("%");
368 cMUON4->Print("MUON4.eps") ;
371 sprintf(line, ".!tar -zcf %s.tar.gz *.eps", GetName()) ;
372 gROOT->ProcessLine(line);
373 sprintf(line, ".!rm -fR *.eps");
374 gROOT->ProcessLine(line);
376 AliInfo(Form("!!! All the eps files are in %s.tar.gz !!!", GetName())) ;
380 report="Problems found, please check!!!";
384 AliInfo(Form("*** %s Summary Report: %s \n",GetName(), report)) ;