]> git.uio.no Git - u/mrichter/AliRoot.git/blame - MUON/AliMUONCheck.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
- Moving class from base to sim
[u/mrichter/AliRoot.git] / MUON / AliMUONCheck.cxx
CommitLineData
70b4a8d6 1/**************************************************************************
2* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3* *
4* Author: The ALICE Off-line Project. *
5* Contributors are mentioned in the code where appropriate. *
6* *
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**************************************************************************/
15
16// $Id$
17
9780bd43 18/// \class AliMUONCheck
19///
20/// This class check the ESD tree, providing the matching with the trigger
21/// response and designing useful plots (Pt, Y, ITS vertex, multiplicity).
22/// Note that there is a special flag to turn on for pdc06 production.
23/// It also checks the TR tree giving hit densities on the two first and
24/// last planes of the spectrometer as well as the time of flight on these planes.
25/// MUONkine() provides event stack and check if the event are generated with
26/// at least one muon and two muons (for PDC06).
27/// DumpDigit() as a replacement of the function from MUONCheck.C macro.
28///
29/// \author Frederic Yermia, INFN Torino
30
70b4a8d6 31#include "AliMUONCheck.h"
9780bd43 32#include "AliMUONData.h"
33#include "AliMUONDigit.h"
34#include "AliMUONConstants.h"
35#include "AliMUONTrack.h"
36#include "AliMUONTrackParam.h"
37#include "AliMUONTrackExtrap.h"
38
39#include "AliMpSegmentation.h"
40#include "AliMpVSegmentation.h"
41#include "AliMpDEManager.h"
70b4a8d6 42
9780bd43 43#include "AliRunLoader.h"
70b4a8d6 44#include "AliLoader.h"
9780bd43 45#include "AliStack.h"
46#include "AliTrackReference.h"
47#include "AliTracker.h"
48#include "AliESD.h"
49#include "AliESDMuonTrack.h"
50#include "AliMagFMaps.h"
70b4a8d6 51#include "AliLog.h"
70b4a8d6 52
9780bd43 53#include <TSystem.h>
54#include <TCanvas.h>
55#include <TLorentzVector.h>
56#include <TFile.h>
57#include <TH1.h>
58#include <TParticle.h>
70b4a8d6 59
5398f946 60/// \cond CLASSIMP
70b4a8d6 61ClassImp(AliMUONCheck)
5398f946 62/// \endcond
5a560b6a 63AliMUONCheck::AliMUONCheck(const char* galiceFile, const char* esdFile,Int_t firstEvent, Int_t lastEvent,const char* outDir)
64: TObject(),
65 fFileName(galiceFile),
66 fFileNameSim(galiceFile),
67 fesdFileName(esdFile),
68 foutDir(outDir),
69 fFirstEvent(firstEvent),
70 fLastEvent(lastEvent),
71 fRunLoader(0x0),
72 fRunLoaderSim(0x0),
73 fData(0x0),
74 fDataSim(0x0),
75 fTree(0)
76{
77 /// ctor
78
79 fRunLoader = AliRunLoader::Open(fFileName.Data(),"MUONFolder","READ");
80 if (!fRunLoader)
81 {
82 AliError(Form("Error opening %s file \n",fFileName.Data()));
83 }
84 else
85 {
86 fLoader = fRunLoader->GetLoader("MUONLoader");
87 if ( fLoader )
88 {
89 fData = new AliMUONData(fLoader,"MUON","MUON");
90 }
91 else
92 {
93 AliError(Form("Could get MUONLoader"));
94 }
95 }
96
97 fRunLoaderSim = AliRunLoader::Open(fFileNameSim.Data(),"MUONFolderSim","READ");
98 if (!fRunLoaderSim)
99 {
100 AliError(Form("Error opening %s file \n",fFileNameSim.Data()));
101 }
102 else
103 {
104 fLoaderSim = fRunLoaderSim->GetLoader("MUONLoader");
105 if ( fLoaderSim )
106 {
107 fDataSim = new AliMUONData(fLoaderSim,"MUON","MUON");
108 }
109 else
110 {
111 AliError(Form("Could get MUONLoader"));
112 }
113 }
114
115 char command[120];
116 sprintf(command,"rm -rf %s", foutDir);
117 gSystem->Exec(command);
118 gSystem->mkdir(foutDir);
70b4a8d6 119
5a560b6a 120}
70b4a8d6 121//_____________________________________________________________________________
5a560b6a 122AliMUONCheck::AliMUONCheck(const char* galiceFile, const char* galiceFileSim, const char* esdFile,Int_t firstEvent, Int_t lastEvent,const char* outDir)
70b4a8d6 123: TObject(),
9780bd43 124 fFileName(galiceFile),
5a560b6a 125 fFileNameSim(galiceFileSim),
9780bd43 126 fesdFileName(esdFile),
127 foutDir(outDir),
128 fFirstEvent(firstEvent),
129 fLastEvent(lastEvent),
130 fRunLoader(0x0),
5a560b6a 131 fRunLoaderSim(0x0),
9780bd43 132 fData(0x0),
5a560b6a 133 fDataSim(0x0),
9780bd43 134 fTree(0)
70b4a8d6 135{
71a2d3aa 136 /// ctor
9780bd43 137
70b4a8d6 138 fRunLoader = AliRunLoader::Open(fFileName.Data(),"MUONFolder","READ");
139 if (!fRunLoader)
140 {
141 AliError(Form("Error opening %s file \n",fFileName.Data()));
142 }
143 else
144 {
9780bd43 145 fLoader = fRunLoader->GetLoader("MUONLoader");
146 if ( fLoader )
70b4a8d6 147 {
9780bd43 148 fData = new AliMUONData(fLoader,"MUON","MUON");
70b4a8d6 149 }
150 else
151 {
152 AliError(Form("Could get MUONLoader"));
153 }
154 }
9780bd43 155
5a560b6a 156 fRunLoaderSim = AliRunLoader::Open(fFileNameSim.Data(),"MUONFolderSim","READ");
157 if (!fRunLoaderSim)
158 {
159 AliError(Form("Error opening %s file \n",fFileNameSim.Data()));
160 }
161 else
162 {
163 fLoaderSim = fRunLoaderSim->GetLoader("MUONLoader");
164 if ( fLoaderSim )
165 {
166 fDataSim = new AliMUONData(fLoaderSim,"MUON","MUON");
167 }
168 else
169 {
170 AliError(Form("Could get MUONLoader"));
171 }
172 }
173
9780bd43 174 char command[120];
175 sprintf(command,"rm -rf %s", foutDir);
176 gSystem->Exec(command);
177 gSystem->mkdir(foutDir);
178
179}
180
70b4a8d6 181//_____________________________________________________________________________
182AliMUONCheck::~AliMUONCheck()
183{
5398f946 184/// Destructor
9780bd43 185 fRunLoader->UnloadAll();
5a560b6a 186 fRunLoaderSim->UnloadAll();
9780bd43 187 delete fRunLoader;
70b4a8d6 188 delete fData;
5a560b6a 189 delete fRunLoaderSim;
190 delete fDataSim;
70b4a8d6 191}
192
193//_____________________________________________________________________________
194void
9780bd43 195AliMUONCheck::CheckESD(Bool_t pdc06TriggerResponse)
70b4a8d6 196{
71a2d3aa 197 /// Check ESD files
5398f946 198
70b4a8d6 199 if ( !IsValid() ) return;
200
9780bd43 201 // Histograms
202 TH1F * fhMUONVertex ; //!
203 TH1F * fhMUONMult ; //!
204
205 // create histograms
206 fhMUONVertex = new TH1F("hMUONVertex","ITS Vertex" ,100, -25., 25.);
207 fhMUONMult = new TH1F("hMUONMult" ,"Multiplicity of ESD tracks",10, -0.5, 9.5);
70b4a8d6 208
9780bd43 209 TH1F *hY = new TH1F("hY","Rapidity",100,-5.,-1.);
210 TH1F *hPt = new TH1F("hPt","Pt",100, 0.,20.);
70b4a8d6 211
9780bd43 212 // ------------->open the ESD file
213 TFile* esdFile = TFile::Open(fesdFileName.Data());
214
215 if (!esdFile || !esdFile->IsOpen()) {
216 AliError(Form("Error opening %s file \n",fesdFileName.Data()));
217 return ;}
218
219 Int_t fSPLowpt=0 ; //!
220 Int_t fSPHighpt=0 ; //!
221 Int_t fSPAllpt=0 ; //!
222 Int_t fSMLowpt=0 ; //!
223 Int_t fSMHighpt =0 ; //!
224 Int_t fSMAllpt=0 ; //!
225 Int_t fSULowpt=0 ; //!
226 Int_t fSUHighpt=0 ; //!
227 Int_t fSUAllpt=0 ; //!
228 Int_t fUSLowpt=0 ; //!
229 Int_t fUSHighpt=0 ; //!
230 Int_t fUSAllpt=0 ; //!
231 Int_t fLSLowpt=0 ; //!
232 Int_t fLSHighpt=0 ; //!
233 Int_t fLSAllpt=0 ; //!
234
235 Int_t fSLowpt=0 ; //!
236 Int_t fSHighpt=0 ; //!
237
238 Int_t fnTrackTrig=0 ; //!
239 Int_t ftracktot=0 ; //!
240 Int_t effMatch=0 ; //!
241
242 TLorentzVector fV1;
243 Float_t muonMass = 0.105658389;
244 Double_t thetaX, thetaY, pYZ;
245 Double_t fPxRec1, fPyRec1, fPzRec1, fE1;
246 Int_t fZ1;
247
248 AliESD* fESD = new AliESD();
249 TTree* fTree = (TTree*) esdFile->Get("esdTree");
250 if (!fTree) {
251 Error("CheckESD", "no ESD tree found");
252 AliError(Form("CheckESD", "no ESD tree found"));
253 return ;
254 }
255 fTree->SetBranchAddress("ESD", &fESD);
256
257 Int_t fnevents = fRunLoader->GetNumberOfEvents();
258 Int_t endOfLoop = fLastEvent+1;
259
260 if ( fLastEvent == -1 ) endOfLoop = fnevents;
261 Int_t ievent=0;
262 Int_t nev=0;
263
264 for (ievent = fFirstEvent; ievent < endOfLoop; ++ievent )
70b4a8d6 265 {
266 fRunLoader->GetEvent(ievent);
9780bd43 267 nev++;
268
269 fTree->GetEvent(ievent);
270 if (!fESD) {
271 Error("CheckESD", "no ESD object found for event %d", ievent);
272 return ;
273 }
274 AliESDVertex* vertex = (AliESDVertex*) fESD->AliESD::GetVertex();
275
276 Double_t zVertex = 0. ;
277 if (vertex) zVertex = vertex->GetZv();
278
279 Int_t nTracks = (Int_t)fESD->GetNumberOfMuonTracks() ;
280 ULong64_t trigword=fESD->GetTriggerMask();
281
282 if(pdc06TriggerResponse){
283 if (trigword & 0x01) {
284 fSPLowpt++;
285 }
286
287 if (trigword & 0x02){
288 fSPHighpt++;
289 }
290 if (trigword & 0x04){
291 fSPAllpt++;
292 }
293 if (trigword & 0x08){
294 fSMLowpt++;
295 }
296 if (trigword & 0x010){
297 fSMHighpt++;
298 }
299 if (trigword & 0x020){
300 fSMAllpt++;
301 }
302 if (trigword & 0x040){
303 fSULowpt++;
304 }
305 if (trigword & 0x080){
306 fSUHighpt++;
307 }
308 if (trigword & 0x100){
309 fSUAllpt++;
310 }
311 if (trigword & 0x200){
312 fUSLowpt++;
313 }
314
315 if (trigword & 0x400){
316 fUSHighpt++;
317 }
318 if (trigword & 0x800){
319 fUSAllpt++;
320 }
321 if (trigword & 0x1000){
322 fLSLowpt++;
323 }
324
325 if (trigword & 0x2000){
326 fLSHighpt++;
327 }
328
329 if (trigword & 0x4000){
330 fLSAllpt++;
331 }
332 }// if pdc06TriggerResponse
333 else {
334 if (trigword & 0x01) {
335 fSLowpt++;
336 }
337
338 if (trigword & 0x02){
339 fSHighpt++;
340 }
341 if (trigword & 0x04){
342 fLSLowpt++;
343 }
344 if (trigword & 0x08){
345 fLSHighpt++;
346 }
347 if (trigword & 0x010){
348 fUSLowpt++;
349 }
350 if (trigword & 0x020){
351 fUSHighpt++;
352 }
353 }
354
355 Int_t tracktrig=0;
356 Int_t iTrack1 ;
357 for (iTrack1 = 0; iTrack1<nTracks; iTrack1++) { //1st loop
358 AliESDMuonTrack* muonTrack = fESD->GetMuonTrack(iTrack1);
359 ftracktot++;
360
361 thetaX = muonTrack->GetThetaX();
362 thetaY = muonTrack->GetThetaY();
363 pYZ = 1./TMath::Abs(muonTrack->GetInverseBendingMomentum());
70b4a8d6 364
9780bd43 365 fPzRec1 = - pYZ / TMath::Sqrt(1.0 + TMath::Tan(thetaY)*TMath::Tan(thetaY));
366 fPxRec1 = fPzRec1 * TMath::Tan(thetaX);
367 fPyRec1 = fPzRec1 * TMath::Tan(thetaY);
368 fZ1 = Int_t(TMath::Sign(1.,muonTrack->GetInverseBendingMomentum()));
369 fE1 = TMath::Sqrt(muonMass * muonMass + fPxRec1 * fPxRec1 + fPyRec1 * fPyRec1 + fPzRec1 * fPzRec1);
370 fV1.SetPxPyPzE(fPxRec1, fPyRec1, fPzRec1, fE1);
371 // -----------> transverse momentum
372 Float_t pt1 = fV1.Pt();
373 // ----------->Rapidity
374 Float_t y1 = fV1.Rapidity();
70b4a8d6 375
9780bd43 376 if(muonTrack->GetMatchTrigger()) {
377 fnTrackTrig++;
378 tracktrig++;
379 }
380 hY->Fill(y1);
381 hPt->Fill(pt1);
382
383 }// loop on track
384
385 fhMUONVertex->Fill(zVertex) ;
386 fhMUONMult->Fill(Float_t(nTracks)) ;
387
388 } // loop over events
389
390 AliInfo(Form("Terminate %s:", GetName())) ;
391
392 effMatch=100*fnTrackTrig/ftracktot;
393
394 if(pdc06TriggerResponse){
395 printf("=================================================================\n") ;
396 printf("================ %s ESD SUMMARY ==============\n", GetName()) ;
397 printf(" \n") ;
398 printf(" Total number of processed events %d \n", nev) ;
399 printf("\n") ;
400 printf("\n") ;
401 printf("Table 1: \n") ;
402 printf(" Global Trigger output Low pt High pt All\n") ;
403 printf(" number of Single Plus :\t");
404 printf("%i\t%i\t%i\t", fSPLowpt, fSPHighpt, fSPAllpt) ;
405 printf("\n");
406 printf(" number of Single Minus :\t");
407 printf("%i\t%i\t%i\t", fSMLowpt, fSMHighpt, fSMAllpt) ;
408 printf("\n");
409 printf(" number of Single Undefined :\t");
410 printf("%i\t%i\t%i\t", fSULowpt, fSUHighpt, fSUAllpt) ;
411 printf("\n");
412 printf(" number of UnlikeSign pair :\t");
413 printf("%i\t%i\t%i\t", fUSLowpt, fUSHighpt, fUSAllpt) ;
414 printf("\n");
415 printf(" number of LikeSign pair :\t");
416 printf("%i\t%i\t%i\t", fLSLowpt, fLSHighpt, fLSAllpt) ;
417 printf("\n");
418 printf("===================================================\n") ;
419 printf("\n") ;
420 printf("matching efficiency with the trigger for single tracks = %2d %% \n", effMatch);
421 printf("================================================================\n") ; printf("\n") ;
422
423 }//if(pdc06TriggerResponse)
424
425 gSystem->cd(foutDir);
426
427 FILE *outtxt=fopen("output.txt","a");
428 freopen("output.txt","a",outtxt);
429
430 if(pdc06TriggerResponse){
431 fprintf(outtxt," \n");
432 fprintf(outtxt,"===================================================\n");
433 fprintf(outtxt,"================ ESD SUMMARY ==============\n");
434 fprintf(outtxt," \n");
435 fprintf(outtxt," Total number of processed events %d \n", nev);
436 fprintf(outtxt,"\n");
437 fprintf(outtxt,"\n");
438 fprintf(outtxt,"Table 1: \n");
439 fprintf(outtxt," Global Trigger output Low pt High pt All\n");
440 fprintf(outtxt," number of Single Plus :\t");
441 fprintf(outtxt,"%i\t%i\t%i\t",fSPLowpt,fSPHighpt,fSPAllpt);
442 fprintf(outtxt,"\n");
443 fprintf(outtxt," number of Single Minus :\t");
444 fprintf(outtxt,"%i\t%i\t%i\t",fSMLowpt,fSMHighpt,fSMAllpt);
445 fprintf(outtxt,"\n");
446 fprintf(outtxt," number of Single Undefined :\t");
447 fprintf(outtxt,"%i\t%i\t%i\t",fSULowpt,fSUHighpt,fSUAllpt);
448 fprintf(outtxt,"\n");
449 fprintf(outtxt," number of UnlikeSign pair :\t");
450 fprintf(outtxt,"%i\t%i\t%i\t",fUSLowpt,fUSHighpt,fUSAllpt);
451 fprintf(outtxt,"\n");
452 fprintf(outtxt," number of LikeSign pair :\t");
453 fprintf(outtxt,"%i\t%i\t%i\t",fLSLowpt,fLSHighpt, fLSAllpt);
454 fprintf(outtxt,"\n");
455 fprintf(outtxt,"===================================================\n");
456 fprintf(outtxt,"\n");
457 fprintf(outtxt,"matching efficiency with the trigger for single tracks = %2d %% \n", effMatch);
458 }//if(pdc06TriggerResponse)
459
460 else {
461
462 fprintf(outtxt," \n");
463 fprintf(outtxt,"===================================================\n");
464 fprintf(outtxt,"================ ESD SUMMARY ==============\n");
465 fprintf(outtxt," \n");
466 fprintf(outtxt," Total number of processed events %d \n", nev);
467 fprintf(outtxt,"\n");
468 fprintf(outtxt,"\n");
469 fprintf(outtxt,"Table 1: \n");
470 fprintf(outtxt," Global Trigger output Low pt High pt \n");
471 fprintf(outtxt," number of Single :\t");
472 fprintf(outtxt,"%i\t%i\t",fSLowpt,fSHighpt);
473 fprintf(outtxt,"\n");
474 fprintf(outtxt," number of UnlikeSign pair :\t");
475 fprintf(outtxt,"%i\t%i\t",fUSLowpt,fUSHighpt);
476 fprintf(outtxt,"\n");
477 fprintf(outtxt," number of LikeSign pair :\t");
478 fprintf(outtxt,"%i\t%i\t",fLSLowpt,fLSHighpt);
479 fprintf(outtxt,"\n");
480 fprintf(outtxt,"===================================================\n");
481 fprintf(outtxt,"\n");
482 fprintf(outtxt,"matching efficiency with the trigger for single tracks = %2d %% \n", effMatch);
483 }//else
484 fclose(outtxt);
485
486 TCanvas * c1 = new TCanvas("c1", "ESD", 400, 10, 600, 700) ;
487 c1->Divide(1,2) ;
488 c1->cd(1) ;
489 fhMUONVertex->Draw() ;
490 c1->cd(2) ;
491 fhMUONMult->Draw() ;
492 c1->Print("VertexAndMul.eps") ;
493 TCanvas *c2 = new TCanvas("c2","ESD",400,10,600,700);
494 c2->Divide(1,2);
495 c2->cd(1);
496 hY->Draw();
497 c2->cd(2);
498 hPt->Draw();
499 c2->Print("YandPt.eps") ;
70b4a8d6 500}
9780bd43 501
502//_____________________________________________________________________________
503void
504AliMUONCheck::CheckKine()
505{
71a2d3aa 506 /// Check Stack
9780bd43 507 if ( !IsValid() ) return;
508
509 // Stack of particle for each event
510 AliStack* stack;
511
5a560b6a 512 Int_t fnevents = fRunLoaderSim->GetNumberOfEvents();
513 fRunLoaderSim->LoadKinematics("READ");
9780bd43 514
515 Int_t endOfLoop = fLastEvent+1;
516
517 if ( fLastEvent == -1 ) endOfLoop = fnevents;
518
519 Int_t ievent=0;
520 Int_t nev=0;
521 Int_t nmu=0;
522 Int_t nonemu=0;
523 Int_t ndimu=0;
524 Int_t npa=0;
525 Int_t npb=0;
526
527 for (ievent = fFirstEvent; ievent < endOfLoop; ++ievent ) {
528
529 Int_t nmu2=0;
530 nev++;
531
5a560b6a 532 fRunLoaderSim->GetEvent(ievent);
533 stack = fRunLoaderSim->Stack();
9780bd43 534 npa = stack->GetNprimary();
535 npb = stack->GetNtrack();
536 printf("Primary particles %i \n",npa);
537 printf("Sec particles %i \n",npb);
538 printf("=================================================================\n") ;
539 printf("Primary particles listing: \n");
540 printf("=================================================================\n") ;
541 for(Int_t i=0; i<npa; i++) {
542
543 TParticle *p = stack->Particle(i);
544 p->Print("");
545 Int_t pdg=p->GetPdgCode();
546
547
548 if (abs(pdg) == 13) {
549 nmu2++;
550 }
551 }
552 printf("=================================================================\n") ;
553 printf("=================================================================\n") ;
554
555 printf("Secondaries particles listing: \n");
556 printf("=================================================================\n") ;
557 for(Int_t i=npa; i<npb; i++) {
558
559 stack->Particle(i)->Print("");
560 }
561
562 printf("=================================================================\n") ;
563 printf(">>> Event %d, Number of primary particles is %d \n",ievent, npa);
564 printf(">>> Event %d, Number of secondary articles is %d \n",ievent, npb-npa);
565 printf("=================================================================\n");
566 if(nmu2>0){
567 printf(">>> Okay!!! Event %d with at least one muon on primary stack! \n",ievent);
568 nonemu++;
569 }
570
571 if(nmu2==0){
572 printf(">>> Warning!!! Event %d without muon on primary stack! \n",ievent);
573 nmu++;
574 }
575
576 if(nmu2>1){
577 printf(">>> Okay!!! Event %d with at least two muons on primary stack! \n",ievent);
578 ndimu++;
579 }
580 printf("=================================================================\n");
581 printf(" \n");
582 printf(" \n") ;
583 }//ievent
584
5a560b6a 585 fRunLoaderSim->UnloadKinematics();
9780bd43 586
587 printf("=================================================================\n") ;
588 printf(" Total number of processed events %d \n", nev) ;
589 printf(" \n") ;
590
591 if(nmu>0){
592 printf("---> WARNING!!! <---\n");
593 printf(" %i events without muon on primary stack \n",nmu);
594 }
595
596 if(nmu==0){
597 printf("---> OKAY!!! <---\n");
598 printf(" %i events generated with at least one muon on primary stack \n",nonemu);
599 }
600 if(ndimu>0){
601 printf("---> OKAY!!! <---\n");
602 printf(" %i events generated with at least two muons on primary stack \n",ndimu);
603 }
604
605 printf(" \n") ;
606 printf("*** Leaving MuonKine() *** \n");
607 printf("**************************************************************** \n");
608
609 gSystem->cd(foutDir);
610 FILE *outtxt=fopen("output.txt","a");
611 freopen("output.txt","a",outtxt);
612 fprintf(outtxt," \n");
613 fprintf(outtxt,"=================================================================\n");
614 fprintf(outtxt,"================ MUONkine SUMMARY ================\n");
615 fprintf(outtxt,"\n");
616 fprintf(outtxt,"=================================================================\n");
617 fprintf(outtxt," Total number of processed events %d \n", nev) ;
618 fprintf(outtxt," \n");
619
620 if(nmu>0){
621 fprintf(outtxt," ---> WARNING!!! <--- \n");
622 fprintf(outtxt," %i events without muon on primary stack \n",nmu);
623 }
624
625 if(nmu==0){
626 fprintf(outtxt," ---> OKAY!!! <--- \n");
627 fprintf(outtxt," %i events generated with at least one muon on primary stack \n",nonemu);
628 }
629
630 if(ndimu>0){
631 fprintf(outtxt," ---> OKAY!!! <--- \n");
632 fprintf(outtxt," %i events generated with at least two muons on primary stack \n",ndimu);
633 }
634
635 fprintf(outtxt," \n") ;
636 fprintf(outtxt,"*** Leaving MuonKine() *** \n");
637 fprintf(outtxt,"**************************************************************** \n");
638 fclose(outtxt);
639
5a560b6a 640 fRunLoaderSim->UnloadKinematics();
9780bd43 641}
642
643//_____________________________________________________________________________
644void
645AliMUONCheck::CheckTrackRef()
646{
71a2d3aa 647 /// Check TrackRef files
9780bd43 648
649 if ( !IsValid() ) return;
650 Int_t flag11=0,flag12=0,flag13=0,flag14=0;
651
652 TH1F *tof01= new TH1F("tof01","TOF for first tracking plane",100,0.,100);
653 tof01->SetXTitle("tof (ns)");
654 TH1F *tof14= new TH1F("tof14","TOF for MT22",100,0.,100);
655 tof14->SetXTitle("tof (ns)");
656
657 TH1F *hitDensity[4];
658 hitDensity[0] = new TH1F("TR_dhits01","",30,0,300);
659 hitDensity[0]->SetFillColor(3);
660 hitDensity[0]->SetXTitle("R (cm)");
661 hitDensity[1] = new TH1F("TR_dhits10","",30,0,300);
662 hitDensity[1]->SetFillColor(3);
663 hitDensity[1]->SetXTitle("R (cm)");
664 hitDensity[2] = new TH1F("TR_dhits11","",30,0,300);
665 hitDensity[2]->SetFillColor(3);
666 hitDensity[2]->SetXTitle("R (cm)");
667 hitDensity[3] = new TH1F("TR_dhits14","",30,0,300);
668 hitDensity[3]->SetFillColor(3);
669 hitDensity[3]->SetXTitle("R (cm)");
670 Int_t fnevents = fRunLoader->GetNumberOfEvents();
671
5a560b6a 672 fRunLoaderSim->LoadTrackRefs();
9780bd43 673 Int_t endOfLoop = fLastEvent+1;
674
675 if ( fLastEvent == -1 ) endOfLoop = fnevents;
676
677 Int_t ievent=0;
678 Int_t nev=0;
679 Int_t ntot=fLastEvent+1-fFirstEvent;
680 for (ievent = fFirstEvent; ievent < endOfLoop; ++ievent ) {
5a560b6a 681 fRunLoaderSim->GetEvent(ievent);
9780bd43 682 Int_t save=-99;
683 nev++;
5a560b6a 684 TTree *tTR = fRunLoaderSim->TreeTR();
9780bd43 685 Int_t nentries = (Int_t)tTR->GetEntries();
686 TClonesArray *fRefArray = new TClonesArray("AliTrackReference");
687 TBranch *branch = tTR->GetBranch("MUON");
688 branch->SetAddress(&fRefArray);
689
690 for(Int_t l=0; l<nentries; l++)
691 {
692 if(!branch->GetEvent(l)) continue;
693 Int_t nnn = fRefArray->GetEntriesFast();
694
695 for(Int_t k=0; k<nnn; k++)
696 {
697 AliTrackReference *tref = (AliTrackReference*)fRefArray->UncheckedAt(k);
698 Int_t label = tref->GetTrack();
699 Float_t x = tref->X(); // x-pos of hit
700 Float_t y = tref->Y(); // y-pos
701 Float_t z = tref->Z();
702
703 Float_t r=TMath::Sqrt(x*x+y*y);
704 Float_t time = tref->GetTime();
705
706 Float_t wgt=1/(2*10*TMath::Pi()*r)/(ntot);
707
708 if (save!=label){
709 save=label;
710 flag11=0;
711 flag12=0;
712 flag13=0;
713 flag14=0;
714 }
715
716 if (save==label){
717
718 //Ch 1, z=-526.16
719 if (z<=-521&& z>=-531&&flag11==0){
720 flag11=1;
721 hitDensity[0]->Fill(r,wgt);
722 tof01->Fill(1000000000*time,1);
723 };
724
725 //Ch 10, z=-1437.6
726 if (z<=-1432&&z>=-1442&&flag12==0){
727 flag12=1;
728 hitDensity[1]->Fill(r,wgt);
729 }
730
731 //Ch 11, z=-1603.5
732 if (z<=-1598&& z>=-1608&&flag13==0){
733 flag13=1;
734 hitDensity[2]->Fill(r,wgt);
735 };
736
737 //ch 14 z=-1720.5
738 if(z<=-1715&&z>=-1725&&flag14==0){
739 flag14=1;
740 hitDensity[3]->Fill(r,wgt);
741 tof14->Fill(1000000000*time,1);
742 };
743
744 }//if save==label
745
746 }//hits de tTR
747
748 }//entree de tTR
749
750 fRefArray->Delete();
751 delete fRefArray;
752 }//evt loop
753
5a560b6a 754 fRunLoaderSim->UnloadTrackRefs();
9780bd43 755 gSystem->cd(foutDir);
756 TCanvas *c6 = new TCanvas("c6","TOF",400,10,600,700);
757 c6->Divide(1,2);
758 c6->cd(1);
759
760 tof01->Draw();
761 c6->cd(2);
762 tof14->Draw();
763 c6->Print("tof_on_trigger.ps");
764
765 TCanvas *c5 = new TCanvas("c5","TRef:Hits Density",400,10,600,700);
766 c5->Divide(2,2);
767 c5->cd(1);
768 hitDensity[0]->Draw();
769 c5->cd(2);
770 hitDensity[1]->Draw();
771 c5->cd(3);
772 hitDensity[2]->Draw();
773 c5->cd(4);
774 hitDensity[3]->Draw();
775 c5->Print("TR_Hit_densities.ps");
776 printf("=================================================================\n") ;
777 printf("================ %s Tref SUMMARY ==============\n", GetName()) ;
778 printf(" \n") ;
779 printf(" Total number of processed events %d \n", nev) ;
780 printf("*** Leaving TRef() *** \n");
781 printf("*************************************************** \n");
782
5a560b6a 783 fRunLoaderSim->UnloadTrackRefs();
9780bd43 784}
785
786//_____________________________________________________________________________
787void
788AliMUONCheck::CheckOccupancy(Bool_t perDetEle) const
789{
790/// Check occupancy for the first event selected
791
792 // Loading MUON subsystem
5a560b6a 793 fLoaderSim->LoadDigits("READ");
794 AliMUONData muondata(fLoaderSim,"MUON","MUON");
9780bd43 795
796 AliMUONDigit * mDigit =0x0;
797 const AliMpVSegmentation * segbend = 0x0;
798 const AliMpVSegmentation * segnonbend = 0x0;
799 AliMpIntPair pad(0,0);
800
801 Int_t dEoccupancyBending[14][26];
802 Int_t dEoccupancyNonBending[14][26];
803 Int_t cHoccupancyBending[14];
804 Int_t cHoccupancyNonBending[14];
805 Int_t totaloccupancyBending =0;
806 Int_t totaloccupancyNonBending =0;
807
808 Int_t dEchannelsBending[14][26];
809 Int_t dEchannelsNonBending[14][26];
810 Int_t cHchannelsBending[14];
811 Int_t cHchannelsNonBending[14];
812 Int_t totalchannelsBending =0;
813 Int_t totalchannelsNonBending =0;
814
815 Int_t nchambers = AliMUONConstants::NCh(); ;
816 for (Int_t ichamber=0; ichamber<nchambers; ichamber++) {
817 cHchannelsBending[ichamber]=0;
818 cHchannelsNonBending[ichamber]=0;
819
820 for (Int_t idetele=0; idetele<26; idetele++) {
821 Int_t detele = 100*(ichamber +1)+idetele;
822 dEchannelsBending[ichamber][idetele]=0;
823 dEchannelsNonBending[ichamber][idetele]=0;
824 dEoccupancyBending[ichamber][idetele]=0;
825 dEoccupancyNonBending[ichamber][idetele]=0;
826 if ( AliMpDEManager::IsValidDetElemId(detele) ) {
827 segbend = AliMpSegmentation::Instance()
828 ->GetMpSegmentation(detele, AliMp::kCath0);
829 segnonbend = AliMpSegmentation::Instance()
830 ->GetMpSegmentation(detele, AliMp::kCath1);
831 if (AliMpDEManager::GetPlaneType(detele, AliMp::kCath0) != AliMp::kBendingPlane ) {
832 const AliMpVSegmentation* tmp = segbend;
833 segbend = segnonbend;
834 segnonbend = tmp;
835 }
836
837 for (Int_t ix=0; ix<=segbend->MaxPadIndexX(); ix++) {
838 for (Int_t iy=0; iy<=segbend->MaxPadIndexY(); iy++) {
839 pad.SetFirst(ix);
840 pad.SetSecond(iy);
841 if( segbend->HasPad(pad) ) {
842 dEchannelsBending[ichamber][idetele]++;
843 cHchannelsBending[ichamber]++;
844 totalchannelsBending++;
845 }
846 }
847 }
848 for (Int_t ix=0; ix<=segnonbend->MaxPadIndexX(); ix++) {
849 for (Int_t iy=0; iy<=segnonbend->MaxPadIndexY(); iy++) {
850 pad.SetFirst(ix);
851 pad.SetSecond(iy);
852 if(segnonbend->HasPad(pad)) {
853 dEchannelsNonBending[ichamber][idetele]++;
854 cHchannelsNonBending[ichamber]++;
855 totalchannelsNonBending++;
856 }
857 }
858 }
859 if (perDetEle) printf(">>> Detection element %4d has %5d channels in bending and %5d channels in nonbending \n",
860 detele, dEchannelsBending[ichamber][idetele], dEchannelsNonBending[ichamber][idetele] );
861 }
862 }
863 printf(">>> Chamber %2d has %6d channels in bending and %6d channels in nonbending \n",
864 ichamber+1, cHchannelsBending[ichamber], cHchannelsNonBending[ichamber]);
865 }
866 printf(">>Spectrometer has %7d channels in bending and %7d channels in nonbending \n",
867 totalchannelsBending, totalchannelsNonBending);
868
869 // Get event
870 printf(">>> Event %d \n", fFirstEvent);
5a560b6a 871 fRunLoaderSim->GetEvent(fFirstEvent);
9780bd43 872 muondata.SetTreeAddress("D");
873 muondata.GetDigits();
874
875 // Loop on chambers
876 for (Int_t ichamber=0; ichamber<nchambers; ichamber++) {
877 cHoccupancyBending[ichamber] = 0;
878 cHoccupancyNonBending[ichamber]= 0;
879
880 // Loop on digits
881 Int_t ndigits = (Int_t) muondata.Digits(ichamber)->GetEntriesFast();
882 for (Int_t idigit=0; idigit<ndigits; idigit++) {
883 mDigit = static_cast<AliMUONDigit*>(muondata.Digits(ichamber)->At(idigit));
884 Int_t detele = mDigit->DetElemId();
885 Int_t idetele = detele-(ichamber+1)*100;
886 if ( mDigit->Cathode() == 0 ) {
887 cHoccupancyBending[ichamber]++;
888 dEoccupancyBending[ichamber][idetele]++;
889 totaloccupancyBending++;
890 }
891 else {
892 cHoccupancyNonBending[ichamber]++;
893 dEoccupancyNonBending[ichamber][idetele]++;
894 totaloccupancyNonBending++;
895 }
896 }
897
898 printf(">>> Chamber %2d nChannels Bending %5d nChannels NonBending %5d \n",
899 ichamber+1,
900 cHoccupancyBending[ichamber],
901 cHoccupancyNonBending[ichamber]);
902 printf(">>> Chamber %2d Occupancy Bending %5.2f %% Occupancy NonBending %5.2f %% \n",
903 ichamber+1,
904 100.*((Float_t) cHoccupancyBending[ichamber])/((Float_t) cHchannelsBending[ichamber]),
905 100.*((Float_t) cHoccupancyNonBending[ichamber])/((Float_t) cHchannelsBending[ichamber]) );
906
907
908 for(Int_t idetele=0; idetele<26; idetele++) {
909 Int_t detele = idetele + 100*(ichamber+1);
910 if ( AliMpDEManager::IsValidDetElemId(detele) ) {
911 if (perDetEle) {
912 printf(">>> DetEle %4d nChannels Bending %5d nChannels NonBending %5d \n",
913 idetele+100*(ichamber+1),
914 dEoccupancyBending[ichamber][idetele],
915 dEoccupancyNonBending[ichamber][idetele]);
916 printf(">>> DetEle %4d Occupancy Bending %5.2f %% Occupancy NonBending %5.2f %% \n",
917 idetele+100*(ichamber+1),
918 100.*((Float_t) dEoccupancyBending[ichamber][idetele])/((Float_t) dEchannelsBending[ichamber][idetele]),
919 100.*((Float_t) dEoccupancyNonBending[ichamber][idetele])/((Float_t) dEchannelsBending[ichamber][idetele]));
920 }
921 }
922 }
923 } // end chamber loop
924 printf(">>> Muon Spectrometer Occupancy Bending %5.2f %% Occupancy NonBending %5.2f %% \n",
925 100.*((Float_t) totaloccupancyBending)/((Float_t) totalchannelsBending),
926 100.*((Float_t) totaloccupancyNonBending)/((Float_t) totalchannelsNonBending) );
927 muondata.ResetDigits();
928 // } // end cathode loop
5a560b6a 929 fLoaderSim->UnloadDigits();
9780bd43 930}
931
932//_____________________________________________________________________________
933void
934AliMUONCheck::CheckRecTracks () const
935{
936/// Reads and dumps rec tracks objects
937
938 // waiting for mag field in CDB
939 AliInfoStream() << "Loading field map...\n";
940 if (!AliTracker::GetFieldMap()) {
941 AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 1, 1., 10., AliMagFMaps::k5kG);
942 AliTracker::SetFieldMap(field, kFALSE);
943 }
944
945 // Loading data
946 fLoader->LoadTracks("READ");
947
948 Int_t endOfLoop = fLastEvent+1;
949 if ( fLastEvent == -1 ) endOfLoop = fRunLoader->GetNumberOfEvents();
950
951 for (Int_t ievent=fFirstEvent; ievent<endOfLoop; ievent++) {
952 fRunLoader->GetEvent(ievent);
953
954 fData->SetTreeAddress("RT");
955 fData->GetRecTracks();
956 TClonesArray* recTracks = fData->RecTracks();
957
958 Int_t nrectracks = (Int_t) recTracks->GetEntriesFast(); //
959 printf(">>> Event %d, Number of Recconstructed tracks %d \n",ievent, nrectracks);
960
961 // Set the magnetic field for track extrapolations
962 AliMUONTrackExtrap::SetField(AliTracker::GetFieldMap());
963
964 // Loop over tracks
965 for (Int_t iRecTracks = 0; iRecTracks < nrectracks; iRecTracks++) {
966 AliMUONTrack* recTrack = (AliMUONTrack*) recTracks->At(iRecTracks);
967 AliMUONTrackParam* trackParam = (AliMUONTrackParam*) (recTrack->GetTrackParamAtHit())->First();
968 AliMUONTrackExtrap::ExtrapToZ(trackParam,0.);
969 recTrack->Print("full");
970 }
971 fData->ResetRecTracks();
972 }
973 fLoader->UnloadTracks();
974}
975
976//_____________________________________________________________________________
977void AliMUONCheck::SetEventsToCheck(Int_t firstEvent, Int_t lastEvent)
978{
979/// Set first and last event number to check
980
981 fFirstEvent = firstEvent;
982 fLastEvent = lastEvent;
983}