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,""),
58 fmuonMass(0.105658389),
81 // Input slot #0 works with an Ntuple
82 DefineInput(0, TChain::Class());
83 // Output slot #0 writes into a TH1 container
84 DefineOutput(0, TObjArray::Class()) ;
89 //______________________________________________________________________________
90 AliMUONQATask::~AliMUONQATask()
93 fOutputContainer->Clear() ;
94 delete fOutputContainer ;
110 //______________________________________________________________________________
111 void AliMUONQATask::ConnectInputData(const Option_t*)
113 // Initialisation of branch container and histograms
115 AliInfo(Form("*** Initialization of %s", GetName())) ;
118 fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
120 AliError(Form("Input 0 for %s not found\n", GetName()));
124 // One should first check if the branch address was taken by some other task
125 char ** address = (char **)GetBranchAddress(0, "ESD");
127 fESD = (AliESD*)(*address);
130 SetBranchAddress(0, "ESD", &fESD);
134 //________________________________________________________________________
135 void AliMUONQATask::CreateOutputObjects()
141 fhMUONVertex = new TH1F("hMUONVertex","ITS Vertex" ,100, -25., 25.);
142 fhMUONMult = new TH1F("hMUONMult" ,"Multiplicity of ESD tracks",10, -0.5, 9.5);
143 fhPt = new TH1F("hPt","Pt",100, 0.,20.);
144 fhY = new TH1F("hY","Rapidity",100,-5.,-1.);
145 fheffMatchT = new TH1F("heff_matchT","Trigger Matching Efficiency",100, 0.,100.);
146 fhSLowpt = new TH1F("hSLowpt","Single Low Pt Response (%)",101, 0.,101.);
147 fhUSLowpt = new TH1F("hUSLowpt","Unlike Sign Low Pt Response (%)",101, 0.,101.);
148 fhUSHighpt = new TH1F("hUSHighpt","Unlike Sign High Pt Response (%)",101, 0.,101.);
149 fhLSLowpt = new TH1F("hLSLowpt","Like Sign Low Pt Response (%)",101, 0.,101.);
150 fhLSHighpt = new TH1F("hLSHighpt","Like Sign High Pt Response (%)",101, 0.,101.);
151 fhChi2 = new TH1F("hChi2","Chi2 by d.o.f.",100, 0.,20.);
152 fhChi2match = new TH1F("hChi2match","Chi2 of trig/track matching",100, 0.,20.);
153 // create output container
155 fOutputContainer = new TObjArray(12) ;
156 fOutputContainer->SetName(GetName()) ;
157 fOutputContainer->AddAt(fhMUONVertex, 0) ;
158 fOutputContainer->AddAt(fhMUONMult, 1) ;
159 fOutputContainer->AddAt(fhPt, 2) ;
160 fOutputContainer->AddAt(fhY, 3) ;
161 fOutputContainer->AddAt(fheffMatchT, 4) ;
162 fOutputContainer->AddAt(fhSLowpt, 5) ;
163 fOutputContainer->AddAt(fhUSLowpt, 6) ;
164 fOutputContainer->AddAt(fhUSHighpt, 7) ;
165 fOutputContainer->AddAt(fhLSLowpt, 8) ;
166 fOutputContainer->AddAt(fhLSHighpt, 9) ;
167 fOutputContainer->AddAt(fhChi2, 10) ;
168 fOutputContainer->AddAt( fhChi2match, 11) ;
171 //______________________________________________________________________________
172 void AliMUONQATask::Exec(Option_t *)
174 // Processing of one event
178 Long64_t entry = fChain->GetReadEntry() ;
181 AliError("fESD is not connected to the input!") ;
185 if ( !((entry-1)%100) )
186 AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ;
188 // ************************ MUON *************************************
190 const AliESDVertex* vertex = dynamic_cast<const AliESDVertex*>(fESD->GetVertex()) ;
192 Double_t zVertex = 0. ;
194 zVertex = vertex->GetZv() ;
196 Int_t nTracks = fESD->GetNumberOfMuonTracks() ;
198 ULong64_t trigWord = fESD->GetTriggerMask() ;
200 if (trigWord & 0x80) {
203 if (trigWord & 0x100){
206 if (trigWord & 0x200){
209 if (trigWord & 0x400){
212 if (trigWord & 0x800){
216 Int_t tracktrig = 0 ;
219 for (iTrack1 = 0 ; iTrack1 < nTracks ; iTrack1++) { //1st loop
220 AliESDMuonTrack* muonTrack = fESD->GetMuonTrack(iTrack1) ;
222 fthetaX = muonTrack->GetThetaX();
223 fthetaY = muonTrack->GetThetaY();
224 fpYZ = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
225 fPzRec1 = - fpYZ / TMath::Sqrt(1.0 + TMath::Tan(fthetaY)*TMath::Tan(fthetaY));
226 fPxRec1 = fPzRec1 * TMath::Tan(fthetaX);
227 fPyRec1 = fPzRec1 * TMath::Tan(fthetaY);
228 fZ1 = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
229 fE1 = TMath::Sqrt(fmuonMass * fmuonMass + fPxRec1 * fPxRec1 + fPyRec1 * fPyRec1 + fPzRec1 * fPzRec1);
230 fV1.SetPxPyPzE(fPxRec1, fPyRec1, fPzRec1, fE1);
232 // -----------> transverse momentum
233 Float_t pt1 = fV1.Pt();
234 // ----------->Rapidity
235 Float_t y1 = fV1.Rapidity();
237 if(muonTrack->GetMatchTrigger()) {
240 Float_t Chi2match = muonTrack->GetChi2MatchTrigger();
241 fhChi2match->Fill(Chi2match);
244 Float_t fitfmin = muonTrack->GetChi2();
245 Int_t ntrackhits = muonTrack->GetNHit();
246 Float_t Chi2= fitfmin / (2.0 * ntrackhits - 5);
253 fhMUONVertex->Fill(zVertex) ;
254 fhMUONMult->Fill(Float_t(nTracks)) ;
256 PostData(0, fOutputContainer);
259 //______________________________________________________________________________
260 void AliMUONQATask::Terminate(Option_t *)
262 // Processing when the event loop is ended
263 Int_t fSLowPt = fSLowpt;
265 fSLowPt = 100 * fSLowpt / fnevents ;
266 fhSLowpt->Fill(fSLowPt); }
267 Int_t fUSLowPt = fUSLowpt;
269 fUSLowPt = 100 * fUSLowpt / fnevents ;
270 fhUSLowpt->Fill(fUSLowPt); }
271 Int_t fUSHighPt = fUSHighpt;
273 fUSHighPt = 100 * fUSHighpt / fnevents ;
274 fhUSHighpt->Fill(fUSHighPt); }
275 Int_t fLSLowPt = fLSLowpt;
277 fLSLowPt = 100 * fLSLowpt / fnevents ;
278 fhLSLowpt->Fill(fLSLowPt); }
279 Int_t fLSHighPt = fLSHighpt;
281 fLSHighPt = 100 * fLSHighpt / fnevents ;
282 fhLSHighpt->Fill(fLSHighPt); }
284 Int_t effMatch = -1 ;
286 effMatch = 100 * fnTrackTrig / ftracktot ;
287 fheffMatchT->Fill(effMatch);}
289 Bool_t problem = kFALSE ;
290 AliInfo(Form(" *** %s Report:", GetName())) ;
292 fOutputContainer = (TObjArray*)GetOutputData(0);
293 fhMUONVertex = (TH1F*)fOutputContainer->At(0);
294 fhMUONMult = (TH1F*)fOutputContainer->At(1);
295 fhPt = (TH1F*)fOutputContainer->At(2);
296 fhY = (TH1F*)fOutputContainer->At(3);
297 fheffMatchT=(TH1F*)fOutputContainer->At(4);
298 fhSLowpt=(TH1F*)fOutputContainer->At(5);
299 fhUSLowpt=(TH1F*)fOutputContainer->At(6);
300 fhUSHighpt=(TH1F*)fOutputContainer->At(7);
301 fhLSLowpt=(TH1F*)fOutputContainer->At(8);
302 fhLSHighpt=(TH1F*)fOutputContainer->At(9);
304 printf(" Total number of processed events %d \n", fnevents) ;
307 printf(" Table 1: \n") ;
308 printf(" ===================================================\n") ;
309 printf(" Global Trigger output Low pt High pt \n") ;
310 printf(" number of Single :\t");
311 printf(" %i\t", fSLowpt) ;
313 printf(" number of UnlikeSign pair :\t");
314 printf(" %i\t%i\t", fUSLowpt, fUSHighpt) ;
316 printf(" number of LikeSign pair :\t");
317 printf(" %i\t%i\t", fLSLowpt, fLSHighpt) ;
319 printf(" matching efficiency with the trigger for single tracks = %2d %% \n", effMatch);
322 TCanvas * cMUON1 = new TCanvas("cMUON1", "MUON ESD Vert & Mult", 400, 10, 600, 700) ;
323 cMUON1->Divide(1,2) ;
325 fhMUONVertex->SetXTitle("Vz (cm)");
326 fhMUONVertex->Draw() ;
328 fhMUONMult->SetXTitle(" Track Multiplicity");
330 cMUON1->Print("MUON1.eps") ;
332 TCanvas * cMUON2 = new TCanvas("cMUON2", "MUON ESD Pt & Y", 400, 10, 600, 700) ;
333 cMUON2->Divide(1,2) ;
335 fhPt->SetXTitle("Pt (GeV)");
340 cMUON2->Print("MUON2.eps") ;
342 TCanvas * cMUON3 = new TCanvas("cMUON3", "Track Chi2 by dof and Chi2 of trigger/track matching ", 400, 10, 600, 700) ;
343 cMUON3->Divide(1,2) ;
345 fhChi2->SetXTitle("Chi2 by d.o.f.");
348 fhChi2match->SetXTitle("Chi2 of trig/track matching");
350 cMUON3->Print("MUON3.eps") ;
352 TCanvas * cMUON4 = new TCanvas("cMUON4", "Trigger Matching and Trigger Response (%)", 400, 10, 600, 700) ;
353 cMUON4->Divide(2,3) ;
355 fheffMatchT->SetXTitle("%");
356 fheffMatchT->Draw() ;
358 fhSLowpt->SetXTitle("%");
361 fhUSLowpt->SetXTitle("%");
364 fhUSHighpt->SetXTitle("%");
367 fhLSLowpt->SetXTitle("%");
370 fhLSHighpt->SetXTitle("%");
372 cMUON4->Print("MUON4.eps") ;
375 sprintf(line, ".!tar -zcf %s.tar.gz *.eps", GetName()) ;
376 gROOT->ProcessLine(line);
377 sprintf(line, ".!rm -fR *.eps");
378 gROOT->ProcessLine(line);
380 AliInfo(Form("!!! All the eps files are in %s.tar.gz !!!", GetName())) ;
384 report="Problems found, please check!!!";
388 AliInfo(Form("*** %s Summary Report: %s \n",GetName(), report.Data())) ;