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>
37 #include "AliMUONQATask.h"
40 #include "AliESDVertex.h"
41 #include "AliESDMuonTrack.h"
42 #include <TLorentzVector.h>
43 //______________________________________________________________________________
44 AliMUONQATask::AliMUONQATask(const char *name) :
45 AliAnalysisTask(name,""),
56 fmuonMass(0.105658389),
79 // Input slot #0 works with an Ntuple
80 DefineInput(0, TChain::Class());
81 // Output slot #0 writes into a TH1 container
82 DefineOutput(0, TObjArray::Class()) ;
87 //______________________________________________________________________________
88 AliMUONQATask::~AliMUONQATask()
91 fOutputContainer->Clear() ;
92 delete fOutputContainer ;
108 //______________________________________________________________________________
109 void AliMUONQATask::ConnectInputData(const Option_t*)
111 // Initialisation of branch container and histograms
113 AliInfo(Form("*** Initialization of %s", GetName())) ;
116 fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
118 AliError(Form("Input 0 for %s not found\n", GetName()));
122 // One should first check if the branch address was taken by some other task
123 char ** address = (char **)GetBranchAddress(0, "ESD");
125 fESD = (AliESD*)(*address);
128 SetBranchAddress(0, "ESD", &fESD);
132 //________________________________________________________________________
133 void AliMUONQATask::CreateOutputObjects()
139 fhMUONVertex = new TH1F("hMUONVertex","ITS Vertex" ,100, -25., 25.);
140 fhMUONMult = new TH1F("hMUONMult" ,"Multiplicity of ESD tracks",10, -0.5, 9.5);
141 fhPt = new TH1F("hPt","Pt",100, 0.,20.);
142 fhY = new TH1F("hY","Rapidity",100,-5.,-1.);
143 fheffMatchT = new TH1F("heff_matchT","Trigger Matching Efficiency",100, 0.,100.);
144 fhSLowpt = new TH1F("hSLowpt","Single Low Pt Response (%)",101, 0.,101.);
145 fhUSLowpt = new TH1F("hUSLowpt","Unlike Sign Low Pt Response (%)",101, 0.,101.);
146 fhUSHighpt = new TH1F("hUSHighpt","Unlike Sign High Pt Response (%)",101, 0.,101.);
147 fhLSLowpt = new TH1F("hLSLowpt","Like Sign Low Pt Response (%)",101, 0.,101.);
148 fhLSHighpt = new TH1F("hLSHighpt","Like Sign High Pt Response (%)",101, 0.,101.);
149 fhChi2 = new TH1F("hChi2","Chi2 by d.o.f.",100, 0.,20.);
150 fhChi2match = new TH1F("hChi2match","Chi2 of trig/track matching",100, 0.,20.);
151 // create output container
153 fOutputContainer = new TObjArray(12) ;
154 fOutputContainer->SetName(GetName()) ;
155 fOutputContainer->AddAt(fhMUONVertex, 0) ;
156 fOutputContainer->AddAt(fhMUONMult, 1) ;
157 fOutputContainer->AddAt(fhPt, 2) ;
158 fOutputContainer->AddAt(fhY, 3) ;
159 fOutputContainer->AddAt(fheffMatchT, 4) ;
160 fOutputContainer->AddAt(fhSLowpt, 5) ;
161 fOutputContainer->AddAt(fhUSLowpt, 6) ;
162 fOutputContainer->AddAt(fhUSHighpt, 7) ;
163 fOutputContainer->AddAt(fhLSLowpt, 8) ;
164 fOutputContainer->AddAt(fhLSHighpt, 9) ;
165 fOutputContainer->AddAt(fhChi2, 10) ;
166 fOutputContainer->AddAt( fhChi2match, 11) ;
169 //______________________________________________________________________________
170 void AliMUONQATask::Exec(Option_t *)
172 // Processing of one event
176 Long64_t entry = fChain->GetReadEntry() ;
179 AliError("fESD is not connected to the input!") ;
183 if ( !((entry-1)%100) )
184 AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ;
186 // ************************ MUON *************************************
188 const AliESDVertex* vertex = dynamic_cast<const AliESDVertex*>(fESD->GetVertex()) ;
190 Double_t zVertex = 0. ;
192 zVertex = vertex->GetZv() ;
194 Int_t nTracks = fESD->GetNumberOfMuonTracks() ;
196 ULong64_t trigWord = fESD->GetTriggerMask() ;
198 if (trigWord & 0x80) {
201 if (trigWord & 0x100){
204 if (trigWord & 0x200){
207 if (trigWord & 0x400){
210 if (trigWord & 0x800){
214 Int_t tracktrig = 0 ;
217 for (iTrack1 = 0 ; iTrack1 < nTracks ; iTrack1++) { //1st loop
218 AliESDMuonTrack* muonTrack = fESD->GetMuonTrack(iTrack1) ;
220 fthetaX = muonTrack->GetThetaX();
221 fthetaY = muonTrack->GetThetaY();
222 fpYZ = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
223 fPzRec1 = - fpYZ / TMath::Sqrt(1.0 + TMath::Tan(fthetaY)*TMath::Tan(fthetaY));
224 fPxRec1 = fPzRec1 * TMath::Tan(fthetaX);
225 fPyRec1 = fPzRec1 * TMath::Tan(fthetaY);
226 fZ1 = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
227 fE1 = TMath::Sqrt(fmuonMass * fmuonMass + fPxRec1 * fPxRec1 + fPyRec1 * fPyRec1 + fPzRec1 * fPzRec1);
228 fV1.SetPxPyPzE(fPxRec1, fPyRec1, fPzRec1, fE1);
230 // -----------> transverse momentum
231 Float_t pt1 = fV1.Pt();
232 // ----------->Rapidity
233 Float_t y1 = fV1.Rapidity();
235 if(muonTrack->GetMatchTrigger()) {
238 Float_t Chi2match = muonTrack->GetChi2MatchTrigger();
239 fhChi2match->Fill(Chi2match);
242 Float_t fitfmin = muonTrack->GetChi2();
243 Int_t ntrackhits = muonTrack->GetNHit();
244 Float_t Chi2= fitfmin / (2.0 * ntrackhits - 5);
251 fhMUONVertex->Fill(zVertex) ;
252 fhMUONMult->Fill(Float_t(nTracks)) ;
254 PostData(0, fOutputContainer);
257 //______________________________________________________________________________
258 void AliMUONQATask::Terminate(Option_t *)
260 // Processing when the event loop is ended
261 Int_t fSLowPt = fSLowpt;
263 fSLowPt = 100 * fSLowpt / fnevents ;
264 fhSLowpt->Fill(fSLowPt); }
265 Int_t fUSLowPt = fUSLowpt;
267 fUSLowPt = 100 * fUSLowpt / fnevents ;
268 fhUSLowpt->Fill(fUSLowPt); }
269 Int_t fUSHighPt = fUSHighpt;
271 fUSHighPt = 100 * fUSHighpt / fnevents ;
272 fhUSHighpt->Fill(fUSHighPt); }
273 Int_t fLSLowPt = fLSLowpt;
275 fLSLowPt = 100 * fLSLowpt / fnevents ;
276 fhLSLowpt->Fill(fLSLowPt); }
277 Int_t fLSHighPt = fLSHighpt;
279 fLSHighPt = 100 * fLSHighpt / fnevents ;
280 fhLSHighpt->Fill(fLSHighPt); }
282 Int_t effMatch = -1 ;
284 effMatch = 100 * fnTrackTrig / ftracktot ;
285 fheffMatchT->Fill(effMatch);}
287 Bool_t problem = kFALSE ;
288 AliInfo(Form(" *** %s Report:", GetName())) ;
290 fOutputContainer = (TObjArray*)GetOutputData(0);
291 fhMUONVertex = (TH1F*)fOutputContainer->At(0);
292 fhMUONMult = (TH1F*)fOutputContainer->At(1);
293 fhPt = (TH1F*)fOutputContainer->At(2);
294 fhY = (TH1F*)fOutputContainer->At(3);
295 fheffMatchT=(TH1F*)fOutputContainer->At(4);
296 fhSLowpt=(TH1F*)fOutputContainer->At(5);
297 fhUSLowpt=(TH1F*)fOutputContainer->At(6);
298 fhUSHighpt=(TH1F*)fOutputContainer->At(7);
299 fhLSLowpt=(TH1F*)fOutputContainer->At(8);
300 fhLSHighpt=(TH1F*)fOutputContainer->At(9);
302 printf(" Total number of processed events %d \n", fnevents) ;
305 printf(" Table 1: \n") ;
306 printf(" ===================================================\n") ;
307 printf(" Global Trigger output Low pt High pt \n") ;
308 printf(" number of Single :\t");
309 printf(" %i\t", fSLowpt) ;
311 printf(" number of UnlikeSign pair :\t");
312 printf(" %i\t%i\t", fUSLowpt, fUSHighpt) ;
314 printf(" number of LikeSign pair :\t");
315 printf(" %i\t%i\t", fLSLowpt, fLSHighpt) ;
317 printf(" matching efficiency with the trigger for single tracks = %2d %% \n", effMatch);
320 TCanvas * cMUON1 = new TCanvas("cMUON1", "MUON ESD Vert & Mult", 400, 10, 600, 700) ;
321 cMUON1->Divide(1,2) ;
323 fhMUONVertex->SetXTitle("Vz (cm)");
324 fhMUONVertex->Draw() ;
326 fhMUONMult->SetXTitle(" Track Multiplicity");
328 cMUON1->Print("MUON1.eps") ;
330 TCanvas * cMUON2 = new TCanvas("cMUON2", "MUON ESD Pt & Y", 400, 10, 600, 700) ;
331 cMUON2->Divide(1,2) ;
333 fhPt->SetXTitle("Pt (GeV)");
338 cMUON2->Print("MUON2.eps") ;
340 TCanvas * cMUON3 = new TCanvas("cMUON3", "Track Chi2 by dof and Chi2 of trigger/track matching ", 400, 10, 600, 700) ;
341 cMUON3->Divide(1,2) ;
343 fhChi2->SetXTitle("Chi2 by d.o.f.");
346 fhChi2match->SetXTitle("Chi2 of trig/track matching");
348 cMUON3->Print("MUON3.eps") ;
350 TCanvas * cMUON4 = new TCanvas("cMUON4", "Trigger Matching and Trigger Response (%)", 400, 10, 600, 700) ;
351 cMUON4->Divide(2,3) ;
353 fheffMatchT->SetXTitle("%");
354 fheffMatchT->Draw() ;
356 fhSLowpt->SetXTitle("%");
359 fhUSLowpt->SetXTitle("%");
362 fhUSHighpt->SetXTitle("%");
365 fhLSLowpt->SetXTitle("%");
368 fhLSHighpt->SetXTitle("%");
370 cMUON4->Print("MUON4.eps") ;
373 sprintf(line, ".!tar -zcf %s.tar.gz *.eps", GetName()) ;
374 gROOT->ProcessLine(line);
375 sprintf(line, ".!rm -fR *.eps");
376 gROOT->ProcessLine(line);
378 AliInfo(Form("!!! All the eps files are in %s.tar.gz !!!", GetName())) ;
382 report="Problems found, please check!!!";
386 AliInfo(Form("*** %s Summary Report: %s \n",GetName(), report.Data())) ;