]> git.uio.no Git - u/mrichter/AliRoot.git/blame - MUON/AliMUONTriggerChamberEff.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
removing old AliMultiplicityCorrectionSelector
[u/mrichter/AliRoot.git] / MUON / AliMUONTriggerChamberEff.cxx
CommitLineData
f165a2d2 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
9edefa04 16/* $Id$ */
17
f165a2d2 18/// \class AliMUONTriggerChamberEff
19/// implementation of the trigger chamber efficiency determination from
20/// data, and returns the
21/// efficiencyCells.dat with the calculated efficiencies
22/// Author: Diego Stocco (Torino)
23
9edefa04 24#include <fstream>
25
26#include <TFile.h>
27#include <TH1F.h>
28#include <TMath.h>
29#include <TString.h>
30
f165a2d2 31#include "AliMUONTriggerChamberEff.h"
32#include "AliMUONDigit.h"
33#include "AliMUONConstants.h"
34#include "AliMUONGlobalTrigger.h"
35#include "AliMUONGeometryTransformer.h"
36#include "AliMUONSegmentation.h"
37#include "AliMUON.h"
38#include "AliMUONData.h"
39#include "AliMUONTriggerTrack.h"
40#include "AliMpVSegmentation.h"
41#include "AliMpSegmentation.h"
42#include "AliMpPad.h"
43#include "AliMpDEIterator.h"
44#include "AliRunLoader.h"
45#include "AliRun.h"
46
f165a2d2 47ClassImp(AliMUONTriggerChamberEff)
48
49//_____________________________________________________________________________
50AliMUONTriggerChamberEff::AliMUONTriggerChamberEff(const char* galiceFile,
51 Int_t firstEvent, Int_t lastEvent)
52: TObject(),
53 fFirstEvent(firstEvent),
54 fLastEvent(lastEvent),
55 fFirstRun(-1),
56 fLastRun(-1),
57 fRunLoader(0x0),
58 fData(0x0),
59 fReproduceTrigResponse(kFALSE),
60 fPrintInfo(kFALSE),
61 fMUON(0x0),
62 fDebugLevel(0),
63 fGaliceDir(0x0)
64{
65/// Standard constructor
66 SetGaliceFile(galiceFile);
67 ResetArrays();
68}
69
70//_____________________________________________________________________________
71AliMUONTriggerChamberEff::AliMUONTriggerChamberEff(Int_t firstRun, Int_t lastRun,
72 const char* galiceRunDir,
73 Int_t firstEvent, Int_t lastEvent)
74: TObject(),
75 fFirstEvent(firstEvent),
76 fLastEvent(lastEvent),
77 fFirstRun(firstRun),
78 fLastRun(lastRun),
79 fRunLoader(0x0),
80 fData(0x0),
81 fReproduceTrigResponse(kFALSE),
82 fPrintInfo(kFALSE),
83 fMUON(0x0),
84 fDebugLevel(0),
85 fGaliceDir(galiceRunDir)
86{
87 ResetArrays();
88 delete gAlice;
89}
90
91//_____________________________________________________________________________
92AliMUONTriggerChamberEff::~AliMUONTriggerChamberEff()
93{
94/// Destructor
95 fRunLoader->UnloadAll();
96 delete fRunLoader;
97 delete fData;
98}
99
100//_____________________________________________________________________________
101void AliMUONTriggerChamberEff::SetGaliceFile(const char *galiceFile)
102{
103 fRunLoader = AliRunLoader::Open(galiceFile,"MUONFolder","READ");
104 if (!fRunLoader)
105 {
106 AliError(Form("Error opening %s file \n",galiceFile));
107 }
108 else
109 {
110 fRunLoader->LoadgAlice();
111 gAlice = fRunLoader->GetAliRun();
112 fMUON = (AliMUON*)gAlice->GetModule("MUON");
113
114 if(fLastEvent<=0 || fLastEvent>fRunLoader->GetNumberOfEvents())fLastEvent = fRunLoader->GetNumberOfEvents()-1;
115 if(fFirstEvent<0)fFirstEvent=0;
116
117
118 AliLoader* loader = fRunLoader->GetLoader("MUONLoader");
119 if ( loader )
120 {
121 fData = new AliMUONData(loader,"MUON","MUON");
122 loader->LoadTracks("READ");
123 loader->LoadDigits("READ");
124 }
125 else
126 {
127 AliError(Form("Could get MUONLoader"));
128 }
129 }
130}
131
132//_____________________________________________________________________________
133void AliMUONTriggerChamberEff::CleanGalice()
134{
135 fRunLoader->UnloadAll();
136 delete fRunLoader;
137}
138
139//_____________________________________________________________________________
140void AliMUONTriggerChamberEff::ResetArrays()
141{
142 for(Int_t ch=0; ch<fgkNchambers; ch++){
143 for(Int_t cath=0; cath<fgkNcathodes; cath++){
144 fTrigger34[ch][cath] = 0;
145 fTrigger44[cath] = 0;
146 for(Int_t slat=0; slat<fgkNslats; slat++){
147 fInefficientSlat[ch][cath][slat] = 0;
148 fHitPerSlat[ch][cath][slat] = 0;
149 }
150 }
151 }
152}
153
154
155//_____________________________________________________________________________
156void AliMUONTriggerChamberEff::InfoDigit()
157{
158 AliMUONDigit * mDigit=0x0;
159 Int_t firstTrigCh = AliMUONConstants::NTrackingCh();
160 // Addressing
161 Int_t nchambers = AliMUONConstants::NCh();
162 fData->SetTreeAddress("D,GLT");
163
164 fData->GetDigits();
165 // Loop on chambers
166 for(Int_t ichamber=firstTrigCh; ichamber<nchambers; ichamber++) {
167 TClonesArray* digits = fData->Digits(ichamber);
168 digits->Sort();
169 Int_t ndigits = (Int_t)digits->GetEntriesFast();
170 for(Int_t idigit=0; idigit<ndigits; idigit++) {
171 mDigit = (AliMUONDigit*)digits->At(idigit);
172 mDigit->Print();
173 } // end digit loop
174 } // end chamber loop
175 fData->ResetDigits();
176}
177
178
179//_____________________________________________________________________________
180Bool_t AliMUONTriggerChamberEff::PadMatchTrack(Float_t xPad, Float_t yPad, Float_t dpx, Float_t dpy,
181 Float_t xTrackAtPad, Float_t yTrackAtPad, Int_t chamber)
182{
183 Float_t numOfHalfWidth = 5.;
184 Bool_t match = kFALSE;
185 Float_t maxDistX = dpx;
186 if(fReproduceTrigResponse && chamber>=2) maxDistX = 3.*dpx;// Non-bending plane: check the +- 1 strip between stations
187 if(!fReproduceTrigResponse)maxDistX = numOfHalfWidth*dpx;
188 Float_t maxDistY = dpy;
189 if(fReproduceTrigResponse && chamber%2) maxDistY = 3*dpy;// bending plane: check the +- 1 strip between planes in the same station
190 if(!fReproduceTrigResponse) maxDistY = numOfHalfWidth*dpy;
191 Float_t deltaX = TMath::Abs(xPad-xTrackAtPad);
192 Float_t deltaY = TMath::Abs(yPad-yTrackAtPad);
193 if(deltaX<=maxDistX && deltaY<=maxDistY)match = kTRUE;
194 return match;
195}
196
197
198//_____________________________________________________________________________
199Bool_t AliMUONTriggerChamberEff::IsDiffLocalBoard(Int_t currDetElemId, Int_t iy, Int_t detElemIdP1, Int_t iyDigitP1)
200{
201 Bool_t isDiff = kTRUE;
202 if(detElemIdP1<0 || iyDigitP1<0)return kFALSE;
203 Int_t currSlat = currDetElemId%100;
204 Int_t slatP1 = detElemIdP1%100;
205 Int_t currLoc = iy/16;
206 Int_t locP1 = iyDigitP1/16;
207 if(currSlat==slatP1 && currLoc==locP1)isDiff = kFALSE;
208 return isDiff;
209}
210
211
212//_____________________________________________________________________________
213void AliMUONTriggerChamberEff::CalculateEfficiency(Int_t trigger44, Int_t trigger34,
214 Float_t &efficiency, Float_t &error, Bool_t failuresAsInput)
215{
216 efficiency=-9.;
217 error=0.;
218 if(trigger34>0){
219 efficiency=(Double_t)trigger44/((Double_t)trigger34);
220 if(failuresAsInput)efficiency=1.-(Double_t)trigger44/((Double_t)trigger34);
221 }
222 Double_t q = TMath::Abs(1-efficiency);
223 if(efficiency<0)error=0.0;
224 else error = TMath::Sqrt(efficiency*q/((Double_t)trigger34));
225}
226
227
228//_____________________________________________________________________________
229Int_t AliMUONTriggerChamberEff::DetElemIdFromPos(Float_t x, Float_t y, Int_t chamber, Int_t cathode)
230{
231 Int_t resultingDetElemId = -1;
232 AliMpDEIterator it;
233 const AliMUONGeometryTransformer *kGeomTransformer = fMUON->GetGeometryTransformer();
234 AliMUONSegmentation *segmentation = fMUON->GetSegmentation();
235 for ( it.First(chamber-1); ! it.IsDone(); it.Next() ){
236 Int_t detElemId = it.CurrentDE();
237
238 if ( segmentation->HasDE(detElemId) ){
239 const AliMpVSegmentation* seg =
240 AliMpSegmentation::Instance()->GetMpSegmentation(detElemId,cathode);
241 if (seg){
242 Float_t deltax = seg->Dimensions().X();
243 Float_t deltay = seg->Dimensions().Y();
244 Float_t xlocal1 = -deltax;
245 Float_t ylocal1 = -deltay;
246 Float_t xlocal2 = +deltax;
247 Float_t ylocal2 = +deltay;
248 Float_t xg01, yg01, zg1, xg02, yg02, zg2;
249 kGeomTransformer->Local2Global(detElemId, xlocal1, ylocal1, 0, xg01, yg01, zg1);
250 kGeomTransformer->Local2Global(detElemId, xlocal2, ylocal2, 0, xg02, yg02, zg2);
251
252 Float_t xg1 = xg01, xg2 = xg02, yg1 = yg01, yg2 = yg02;
253
254 if(xg01>xg02){
255 xg1 = xg02;
256 xg2 = xg01;
257 }
258 if(yg01>yg02){
259 yg1 = yg02;
260 yg2 = yg01;
261 }
262
263 if(x>=xg1 && x<=xg2 && y>=yg1 && y<=yg2){
264 resultingDetElemId = detElemId;
265 break;
266 }
267 }
268 }
269 }
270 return resultingDetElemId;
271}
272
273
274//_____________________________________________________________________________
275void AliMUONTriggerChamberEff::LocalBoardFromPos(Float_t x, Float_t y, Int_t detElemId, Int_t cathode, Int_t localBoard[4])
276{
277 for(Int_t loc=0; loc<4; loc++){
278 localBoard[loc]=-1;
279 }
280 const AliMUONGeometryTransformer *kGeomTransformer = fMUON->GetGeometryTransformer();
281 Float_t xl, yl, zl;
282 kGeomTransformer->Global2Local(detElemId, x, y, 0, xl, yl, zl);
283 TVector2 pos(xl,yl);
284 const AliMpVSegmentation* seg =
285 AliMpSegmentation::Instance()->GetMpSegmentation(detElemId,cathode);
286 if (seg){
287 AliMpPad pad = seg->PadByPosition(pos,kFALSE);
288 for (Int_t loc=0; loc<pad.GetNofLocations(); loc++){
289 AliMpIntPair location = pad.GetLocation(loc);
290 localBoard[loc] = location.GetFirst();
291 }
292 }
293}
294
295
296//_____________________________________________________________________________
297void AliMUONTriggerChamberEff::PrintTrigger(AliMUONGlobalTrigger *globalTrig)
298{
299
300 printf("===================================================\n");
301 printf(" Global Trigger output\t \tLow pt\tHigh pt\n");
302
303 printf(" number of Single:\t \t");
304 printf("%i\t",globalTrig->SingleLpt());
305 printf("%i\t",globalTrig->SingleHpt());
306 printf("\n");
307
308 printf(" number of UnlikeSign pair:\t");
309 printf("%i\t",globalTrig->PairUnlikeLpt());
310 printf("%i\t",globalTrig->PairUnlikeHpt());
311 printf("\n");
312
313 printf(" number of LikeSign pair:\t");
314 printf("%i\t",globalTrig->PairLikeLpt());
315 printf("%i\t",globalTrig->PairLikeHpt());
316 printf("\n");
317 printf("===================================================\n");
318 printf("\n");
319}
320
321void AliMUONTriggerChamberEff::PerformTriggerChamberEff(const char* outputDir)
322{
323 enum {bending, nonBending, numOfCathodes};
324 Int_t evtBeforePrint = 1000;
325 Float_t rad2deg = 180./TMath::Pi();
326
327 Int_t chOrder[] = {0,2,1,3};
328 Int_t station[] = {0,0,1,1};
329 Float_t zRealMatch[fgkNchambers] = {0.0};
330 Float_t correctFactor[fgkNcathodes] = {1.};
331
332 Bool_t match[fgkNchambers][fgkNcathodes] = {{kFALSE}};
333
334 TClonesArray *RecTrigTracks = 0x0;
335 AliMUONTriggerTrack *recTrigTrack = 0x0;
336 AliMUONDigit * mDigit = 0x0;
337
338 Float_t zMeanChamber[fgkNchambers];
339 for(Int_t ch=0; ch<fgkNchambers; ch++){
340 zMeanChamber[ch] = AliMUONConstants::DefaultChamberZ(10+ch);
341 }
342
343 TClonesArray * globalTrigger = 0x0;
344 AliMUONGlobalTrigger * gloTrg = 0x0;
345
346 Int_t partNumOfTrig[fgkNchambers][fgkNcathodes] = {{0}};
347 Int_t totNumOfTrig[fgkNchambers][fgkNcathodes] = {{0}};
348 Int_t atLeast1MuPerEv[fgkNchambers][fgkNcathodes] = {{0}};
349 Int_t digitPerTrack[fgkNcathodes] = {0};
350
351 Float_t trackIntersectCh[2][fgkNchambers]={{0.0}};
352
353 Int_t slatInPlane1[2][fgkNcathodes];
354 Int_t iyDigitInPlane1[2][fgkNcathodes];
355
356 const Int_t maxNumOfTracks = 10;
357 Int_t trigScheme[maxNumOfTracks][fgkNchambers][fgkNcathodes]={{{0}}};
358 Int_t triggeredDigits[maxNumOfTracks][fgkNchambers][fgkNcathodes] = {{{-1}}};
359 Int_t slatThatTriggered[maxNumOfTracks][fgkNchambers][fgkNcathodes]={{{-1}}};
360 Int_t boardThatTriggered[maxNumOfTracks][fgkNchambers][fgkNcathodes][4]={{{{-1}}}};
361 Int_t nboard[4]={-1};
362 Int_t ineffBoard[4]={-1};
363
364 const Int_t maxNumOfDigits = 20;
365 Int_t detElOfDigitsInData[maxNumOfDigits][fgkNchambers][fgkNcathodes] = {{{-1}}};
366
367 char filename[150];
368 FileStat_t fs;
369
370 if(fFirstRun<0)fFirstRun=fLastRun=-1;
371
372 for(Int_t iRun = fFirstRun; iRun <= fLastRun; iRun++){// Loop over runs
373 // open run loader and load gAlice
374 if(fFirstRun>=0){
375 cout<<"\n\nRun = "<<iRun<<endl;
376 sprintf(filename, "%s/run%i/galice.root", fGaliceDir.Data(), iRun);
377 if(gSystem->GetPathInfo(filename,fs)){
378 cout<<"Warning: "<<filename<<" not found. Skip to next one"<<endl;
379 continue;
380 }
381 cout<<"Opening file "<<filename<<endl;
382 SetGaliceFile(filename);
383 }
384
385 for (Int_t ievent=fFirstEvent; ievent<=fLastEvent; ievent++) { // event loop
386 Bool_t isClearEvent = kTRUE;
387
388 for(Int_t ch=0; ch<fgkNchambers; ch++){
389 for(Int_t cath=0; cath<fgkNcathodes; cath++){
390 partNumOfTrig[ch][cath]=0;
391 match[ch][cath]=kFALSE;
392 for(Int_t itrack=0; itrack<maxNumOfTracks; itrack++){
393 triggeredDigits[itrack][ch][cath]=-1;
394 slatThatTriggered[itrack][ch][cath]=-1;
395 for(Int_t loc=0; loc<4; loc++){
396 boardThatTriggered[itrack][ch][cath][loc]=-1;
397 }
398 }
399 for(Int_t idig=0; idig<maxNumOfDigits; idig++){
400 detElOfDigitsInData[idig][ch][cath]=-1;
401 }
402 }
403 }
404
405 fRunLoader->GetEvent(ievent);
406 if (ievent%evtBeforePrint==0) printf("\t Event = %d\n",ievent);
407
408 fData->SetTreeAddress("RL");
409 fData->GetRecTriggerTracks();
410 RecTrigTracks = fData->RecTriggerTracks();
411 Int_t nRecTrigTracks = (Int_t) RecTrigTracks->GetEntriesFast();
412
413 fData->SetTreeAddress("D,GLT");
414 fData->GetDigits();
415
416 const AliMUONGeometryTransformer* kGeomTransformer = fMUON->GetGeometryTransformer();
417
418 for (Int_t iRecTrigTrack=0; iRecTrigTrack<nRecTrigTracks; iRecTrigTrack++) {
419 for(Int_t cath=0; cath<fgkNcathodes; cath++){
420 digitPerTrack[cath]=0;
421 for(Int_t sta=0; sta<2; sta++){
422 slatInPlane1[sta][cath] = -9999;
423 iyDigitInPlane1[sta][cath] = -9999;
424 }
425 }
426
427 Bool_t doubleCountTrack = kFALSE;
428
429 // reading info from tracks
430 recTrigTrack = (AliMUONTriggerTrack *)RecTrigTracks->At(iRecTrigTrack);
431 Float_t x11 = recTrigTrack->GetX11();// x position (info from non-bending plane)
432 Float_t y11 = recTrigTrack->GetY11();// y position (info from bending plane)
433 Float_t thetaX = recTrigTrack->GetThetax();
434 Float_t thetaY = recTrigTrack->GetThetay();
435
436 if(fDebugLevel>=3)printf("\tEvent = %i, Track = %i\npos from track: (x,y) = (%f, %f), (thetaX, thetaY) = (%f, %f)\n",ievent,iRecTrigTrack,x11,y11,thetaX*rad2deg,thetaY*rad2deg);
437
438 for(Int_t ch=0; ch<fgkNchambers; ch++) {
439 zRealMatch[ch] = zMeanChamber[ch];
440 for(Int_t cath=0; cath<fgkNcathodes; cath++){
441 trigScheme[iRecTrigTrack][ch][cath] = 0;
442 }
443 }
444
445 for(Int_t ch=0; ch<fgkNchambers; ch++) { // chamber loop
446 Int_t currCh = chOrder[ch];
447 Int_t ichamber = 10+currCh;
448 Int_t currStation = station[currCh];
449 TClonesArray* digits = fData->Digits(ichamber);
450 digits->Sort();
451 Int_t ndigits = (Int_t)digits->GetEntriesFast();
452 if(fDebugLevel>=2){
453 if(fDebugLevel<3)printf("\tEvent = %i, Track = %i\n", ievent, iRecTrigTrack);
454 printf("DigitNum: %i digits detected\n",ndigits);
455 }
456
457 for(Int_t cath=0; cath<fgkNcathodes; cath++){
458 correctFactor[cath]=1.;
459 }
460 // calculate corrections to trigger track theta
461 if(ch>=1)correctFactor[nonBending] = zMeanChamber[0]/zRealMatch[0];// corrects x position
462 if(ch>=2)correctFactor[bending] = (zMeanChamber[2] - zMeanChamber[0]) / (zRealMatch[2] - zRealMatch[0]);// corrects y position
463
464 // searching track intersection with chambers (first approximation)
465 Float_t deltaZ = zMeanChamber[currCh] - zMeanChamber[0];
466 trackIntersectCh[0][currCh] = zMeanChamber[currCh] * TMath::Tan(thetaX) * correctFactor[nonBending];// x position (info from non-bending plane)
467 trackIntersectCh[1][currCh] = y11 + deltaZ * TMath::Tan(thetaY) * correctFactor[bending];// y position (info from bending plane)
468
469 for(Int_t idigit=0; idigit<ndigits; idigit++) { // digit loop
470 mDigit = (AliMUONDigit*)digits->At(idigit);
471 for(Int_t loc=0; loc<4; loc++){
472 nboard[loc]=-1;
473 }
474
475 // searching loaded digit global position and dimension
476 Int_t detElemId = mDigit->DetElemId();
477 Int_t cathode = mDigit->Cathode();
478 Int_t ix = mDigit->PadX();
479 Int_t iy = mDigit->PadY();
480 Float_t xpad, ypad, zpad;
481 if(detElOfDigitsInData[idigit][ch][cathode]==-1)detElOfDigitsInData[idigit][ch][cathode] = detElemId;
482
483 if(fDebugLevel>=2)printf("cathode = %i\n",cathode);
484 const AliMpVSegmentation* seg =
485 AliMpSegmentation::Instance()->GetMpSegmentation(detElemId,cathode);
486
487 AliMpPad pad = seg->PadByIndices(AliMpIntPair(ix,iy),kTRUE);
488 for (Int_t loc=0; loc<pad.GetNofLocations(); loc++){
489 AliMpIntPair location = pad.GetLocation(loc);
490 nboard[loc] = location.GetFirst();
491 }
492
493 // get the pad position and dimensions
494 Float_t xlocal1 = pad.Position().X();
495 Float_t ylocal1 = pad.Position().Y();
496 Float_t dpx = pad.Dimensions().X();
497 Float_t dpy = pad.Dimensions().Y();
498
499 kGeomTransformer->Local2Global(detElemId, xlocal1, ylocal1, 0, xpad, ypad, zpad);
500
501 if(fDebugLevel>=3)printf("ch = %i\t cath = %i\tpad = (%4.1f, %4.1f, %4.1f)\tsize = (%3.1f, %3.1f)\n",currCh,cathode,xpad,ypad,zpad,dpx,dpy);
502
503 // searching track intersection with chambers (second approximation)
504 if(ch>=2){
505 deltaZ = zpad - zRealMatch[0];
506 trackIntersectCh[1][currCh] = y11 + deltaZ * TMath::Tan(thetaY) * correctFactor[bending];// y position (info from bending plane)
507 }
508
509 // deciding if digit matches track
510 Bool_t isDiffLocBoard = kFALSE;
511 if(fReproduceTrigResponse)isDiffLocBoard = IsDiffLocalBoard(detElemId, iy, slatInPlane1[currStation][cathode], iyDigitInPlane1[currStation][cathode]);
512 Bool_t matchPad = PadMatchTrack(xpad, ypad, dpx, dpy, trackIntersectCh[0][currCh], trackIntersectCh[1][currCh], currCh);
513
514 if(matchPad && ch<2){
515 slatInPlane1[currStation][cathode] = detElemId;
516 iyDigitInPlane1[currStation][cathode] = iy;
517 if(fDebugLevel>=3)printf("slatInPlane1[%i][%i] = %i\tiyDigitInPlane1[%i][%i] = %i\n",currStation,cathode,slatInPlane1[currStation][cathode],currStation,cathode,iyDigitInPlane1[currStation][cathode]);
518 }
519
520 if(isDiffLocBoard && fDebugLevel>=1)printf("\tDifferent local board\n");
521
522 match[currCh][cathode] = (matchPad && !isDiffLocBoard);
523
524 if(match[currCh][cathode]){
525 digitPerTrack[cathode]++;
526 trigScheme[iRecTrigTrack][currCh][cathode]++;
527 triggeredDigits[iRecTrigTrack][currCh][cathode] = idigit;
528 slatThatTriggered[iRecTrigTrack][currCh][cathode] = detElemId;
529 for(Int_t loc=0; loc<4; loc++){
530 boardThatTriggered[iRecTrigTrack][currCh][cathode][loc] = nboard[loc];
531 }
532 if(digitPerTrack[cathode]>4 && !fReproduceTrigResponse)isClearEvent = kFALSE;
533 }
534
535 // in case of match, store real z position of the chamber
536 if(matchPad)zRealMatch[currCh] = zpad;
537
538 } // end digit loop
539 } // end chamber loop
540
541 for(Int_t cath=0; cath<fgkNcathodes; cath++){
542 if(digitPerTrack[cath]<3 && !fReproduceTrigResponse)isClearEvent = kFALSE;
543 }
544
545 if(!isClearEvent && !fReproduceTrigResponse){
546 fData->ResetDigits();
547 continue;
548 }
549
550 Int_t commonDigits = 0;
551 Int_t doubleTrack = -1;
552 for(Int_t itrack=0; itrack<iRecTrigTrack; itrack++){
553 for(Int_t ch=0; ch<fgkNchambers; ch++){
554 if(triggeredDigits[itrack][ch][bending]==triggeredDigits[iRecTrigTrack][ch][bending])commonDigits++;
555 }
556 if(commonDigits>=2){
557 doubleCountTrack=kTRUE;
558 doubleTrack = itrack;
559 break;
560 }
561 }
562
563 for(Int_t cath=0; cath<fgkNcathodes; cath++){
564 if(!doubleCountTrack || fReproduceTrigResponse){
565 Int_t is44 = 1;
566 Bool_t goodForSlatEff=kTRUE;
567 Bool_t goodForBoardEff=kTRUE;
568 Int_t firstSlat=slatThatTriggered[iRecTrigTrack][0][cath]%100;
569 if(firstSlat<0)firstSlat=slatThatTriggered[iRecTrigTrack][1][cath]%100;
570 Int_t firstBoard=boardThatTriggered[iRecTrigTrack][0][0][0];
571 if(firstBoard<0)firstBoard=boardThatTriggered[iRecTrigTrack][1][0][0];
572 for(Int_t ch=0; ch<fgkNchambers; ch++){
573 is44 *= trigScheme[iRecTrigTrack][ch][cath];
574 Int_t currSlat=slatThatTriggered[iRecTrigTrack][ch][cath]%100;
575 if(currSlat<0)continue;
576 if(currSlat!=firstSlat)goodForSlatEff=kFALSE;
577 Bool_t atLeastOneLoc=kFALSE;
578 for(Int_t loc=0; loc<4; loc++){
579 Int_t currBoard = boardThatTriggered[iRecTrigTrack][ch][cath][loc];
580 if(currBoard==firstBoard){
581 atLeastOneLoc=kTRUE;
582 break;
583 }
584 }
585 if(!atLeastOneLoc)goodForBoardEff=kFALSE;
586 }
587 if(fDebugLevel==1)printf("\tEvent = %i, Track = %i\n", ievent, iRecTrigTrack);
588 if(is44==1){
589 fTrigger44[cath]++;
590 if(fDebugLevel>=1)printf("Trigger44[%i] = %i\n",cath,fTrigger44[cath]);
591 if(goodForSlatEff){
592 for(Int_t ch=0; ch<fgkNchambers; ch++){
593 for(Int_t slat=0; slat<fgkNslats; slat++){
594 Int_t corrDetEl = (ch+11)*100 + slat;
595 if(corrDetEl==slatThatTriggered[iRecTrigTrack][ch][cath]){
596 fHitPerSlat[ch][cath][slat]++;
597 if(fDebugLevel>=1)printf("Slat that triggered = %i\n",corrDetEl);
598 if(goodForBoardEff && firstBoard>0){
599 fHitPerBoard[ch][cath][firstBoard-1]++;
600 if(fDebugLevel>=1)printf("Board that triggered = %i\n",firstBoard);
601 }
602 else if(fDebugLevel>=1)printf("Event = %i, Track = %i: Particle crossed different boards: rejected!\n",ievent,iRecTrigTrack);
603 }
604 }
605 }
606 }
607 else printf("Event = %i, Track = %i: Particle crossed different slats: rejected!\n",ievent,iRecTrigTrack);
608 }
609 if(digitPerTrack[cath]==3){
610 for(Int_t ch=0; ch<fgkNchambers; ch++){
611 if(match[ch][cath])partNumOfTrig[ch][cath]++;
612 if(trigScheme[iRecTrigTrack][ch][cath]==0){
613 fTrigger34[ch][cath]++;
614 if(fDebugLevel>=1)printf("Trigger34[%i][%i] = %i\n",ch,cath,fTrigger34[ch][cath]);
615 if(!goodForSlatEff){
616 printf("Event %i, Track = %i: Particle crossed different slats: rejected!\n",ievent,iRecTrigTrack);
617 continue;
618 }
619 Int_t ineffSlat = DetElemIdFromPos(trackIntersectCh[0][ch], trackIntersectCh[1][ch], 11+ch, cath);
620 if(fDebugLevel>=1)printf("Slat non efficient = %i\n",ineffSlat);
621 if(ineffSlat>0){
622 Int_t slatInCh = ineffSlat%100;
623 fInefficientSlat[ch][cath][slatInCh]++;
624 for(Int_t idig=0; idig<maxNumOfDigits; idig++){
625 if(ineffSlat==detElOfDigitsInData[idig][ch][cath])cout<<"Warning: "<<ineffSlat<<" is not inefficient!!!"<<endl;
626 }
627 LocalBoardFromPos(trackIntersectCh[0][ch], trackIntersectCh[1][ch], ineffSlat, cath, ineffBoard);
628 Int_t boardNonEff=-1;
629 for(Int_t loc=0; loc<4; loc++){
630 if(ineffBoard[loc]==firstBoard){
631 boardNonEff=ineffBoard[loc];
632 break;
633 }
634 }
635 if(fDebugLevel>=1)printf("Board non efficient = %i\n",boardNonEff);
636 if(boardNonEff>0)fInefficientBoard[ch][cath][boardNonEff-1]++;
637 else if(fDebugLevel>=1){
638 printf("Inefficient board should be %i.\tBoards found:\n", firstBoard);
639 for(Int_t loc=0; loc<4; loc++){
640 printf("%i\t",ineffBoard[loc]);
641 }
642 printf("\n");
643 }
644
645 }
646 }
647 }
648 }
649 }
650 else if(doubleCountTrack){
651 if(fDebugLevel<=1)printf("\n\tEvent = %i, Track = %i: ", ievent,iRecTrigTrack);
652 printf("Double Count Track: %i similar to %i. Track rejected!\n",iRecTrigTrack, doubleTrack);
653 }
654 }
655 }// end trigger tracks loop
656 if(nRecTrigTracks<=0){
657 fData->ResetDigits();
658 continue;
659 }
660 for(Int_t ch=0; ch<fgkNchambers; ch++){
661 for(Int_t cath=0; cath<fgkNcathodes; cath++){
662 totNumOfTrig[ch][cath] += partNumOfTrig[ch][cath];
663 if(partNumOfTrig[ch][cath]>0)atLeast1MuPerEv[ch][cath]++;
664 }
665 }
666
667 if(fPrintInfo){
668 //Global trigger
669 globalTrigger = fData->GlobalTrigger();
670 Int_t nglobals = (Int_t) globalTrigger->GetEntriesFast(); // should be 1
671
672 for (Int_t iglobal=0; iglobal<nglobals; iglobal++) { // Global Trigger
673 gloTrg = (AliMUONGlobalTrigger*)globalTrigger->At(iglobal);
674 }
675 PrintTrigger(gloTrg);
676 InfoDigit();
677 cout<<"\n"<<endl;
678 }
679
680 fData->ResetDigits();
681 }// end event loop
682 if(fFirstRun>=0)CleanGalice();
683 } //end loop over run
684
685 // Write output data
686 WriteEfficiencyMap(outputDir);
687
688 WriteOutput(outputDir, totNumOfTrig, atLeast1MuPerEv);
689}
690
691//_____________________________________________________________________________
692void AliMUONTriggerChamberEff::WriteOutput(const char* outputDir, Int_t totNumOfTrig[4][2], Int_t atLeast1MuPerEv[4][2])
693{
694 char *cathodeName[fgkNcathodes]={"Bending plane", "Non-Bending plane"};
695 char *cathCode[fgkNcathodes] = {"bendPlane", "nonBendPlane"};
696
697 char outFileName[100];
698 sprintf(outFileName, "%s/triggerChamberEff.out",outputDir);
699 FILE *outfile = fopen(outFileName, "w");
700 for(Int_t cath=0; cath<fgkNcathodes; cath++){
701 fprintf(outfile,"%s:\n",cathodeName[cath]);
702 for(Int_t ch=0; ch<fgkNchambers; ch++){
703 fprintf(outfile,"Total number of muon triggers chamber 1%i = %i\n",ch+1,totNumOfTrig[ch][cath]);
704 }
705 fprintf(outfile,"\n");
706 }
707 fprintf(outfile,"\n");
708 for(Int_t cath=0; cath<fgkNcathodes; cath++){
709 fprintf(outfile,"%s:\n",cathodeName[cath]);
710 for(Int_t ch=0; ch<fgkNchambers; ch++){
711 fprintf(outfile,"At least 1 muon triggered chamber 1%i = %i\n",ch+1, atLeast1MuPerEv[ch][cath]);
712 }
713 fprintf(outfile,"\n");
714 }
715 fprintf(outfile,"\n\n");
716 for(Int_t cath=0; cath<fgkNcathodes; cath++){
717 fprintf(outfile,"%s:\n",cathodeName[cath]);
718 fprintf(outfile,"Number of triggers where all chambers counted = %i\n",fTrigger44[cath]);
719 fprintf(outfile,"\n");
720 }
721 fprintf(outfile,"\n");
722 for(Int_t cath=0; cath<fgkNcathodes; cath++){
723 fprintf(outfile,"%s:\n",cathodeName[cath]);
724 for(Int_t ch=0; ch<fgkNchambers; ch++){
725 fprintf(outfile,"Number of triggers where chamber 1%i did not count = %i\n",ch+1,fTrigger34[ch][cath]);
726 }
727 fprintf(outfile,"\n");
728 }
729
730 fprintf(outfile,"\n");
731 for(Int_t cath=0; cath<fgkNcathodes; cath++){
732 fprintf(outfile,"%s:\n",cathodeName[cath]);
733 for(Int_t ch=0; ch<fgkNchambers; ch++){
734 Int_t sumIneff = 0, sumHits = 0;
735 fprintf(outfile,"\n Chamber %1i\n", ch+1);
736 for(Int_t slat=0; slat<fgkNslats; slat++){
737 fprintf(outfile,"Number of triggers where slat %2i - did not count = %5i - was hit (hit%sCh%iSlat%i) = %5i\n",slat,fInefficientSlat[ch][cath][slat],cathCode[cath],11+ch,slat,fHitPerSlat[ch][cath][slat]);
738 sumIneff += fInefficientSlat[ch][cath][slat];
739 sumHits += fHitPerSlat[ch][cath][slat];
740 }
741 fprintf(outfile,"Number of triggers where chamber %1i - did not count = %5i - was hit (hit%sCh%i) = %5i\n",ch+1,sumIneff,cathCode[cath],11+ch,sumHits);
742 }
743 fprintf(outfile,"\n");
744 }
745 fclose(outfile);
746
747 sprintf(outFileName, "%s/triggerChamberEff.root",outputDir);
748 TFile *outputHistoFile = new TFile(outFileName,"RECREATE");
749 TDirectory *dir = gDirectory;
750
751 enum {slatIn11, slatIn12, slatIn13, slatIn14, chamberEff};
752 char *yAxisTitle = "trigger efficiency (a.u.)";
753 char *xAxisTitle = "chamber";
754
755 TH1F *histo[fgkNcathodes][fgkNchambers+1];
756 TH1F *histoBoard[fgkNcathodes][fgkNchambers];
757
758 char histoName[30];
759 char histoTitle[90];
760
761 for(Int_t cath=0; cath<fgkNcathodes; cath++){
762 for(Int_t ch=0; ch<fgkNchambers+1; ch++){
763 if(ch==chamberEff){
764 sprintf(histoName, "%sChamberEff", cathCode[cath]);
765 sprintf(histoTitle, "Chamber efficiency %s", cathCode[cath]);
766 histo[cath][ch] = new TH1F(histoName, histoTitle, fgkNchambers, 11-0.5, 15-0.5);
767 histo[cath][ch]->SetXTitle(xAxisTitle);
768 histo[cath][ch]->SetYTitle(yAxisTitle);
769 histo[cath][ch]->GetXaxis()->SetNdivisions(fgkNchambers);
770 }
771 else {
772 sprintf(histoName, "%sSlatEffChamber%i", cathCode[cath], 11+ch);
773 sprintf(histoTitle, "Chamber %i: slat efficiency %s", 11+ch, cathCode[cath]);
774 histo[cath][ch] = new TH1F(histoName, histoTitle, fgkNslats, 0-0.5, fgkNslats-0.5);
775 histo[cath][ch]->SetXTitle("slat");
776 histo[cath][ch]->SetYTitle(yAxisTitle);
777 histo[cath][ch]->GetXaxis()->SetNdivisions(fgkNslats);
778
779 sprintf(histoName, "%sBoardEffChamber%i", cathCode[cath], 11+ch);
780 sprintf(histoTitle, "Chamber %i: board efficiency %s", 11+ch, cathCode[cath]);
781 histoBoard[cath][ch] = new TH1F(histoName, histoTitle, fgkNboards, 1-0.5, fgkNboards+1.-0.5);
782 histoBoard[cath][ch]->SetXTitle("boards");
783 histoBoard[cath][ch]->SetYTitle(yAxisTitle);
784 histoBoard[cath][ch]->GetXaxis()->SetNdivisions(fgkNboards);
785 }
786 }
787 }
788
789 Float_t efficiency, efficiencyError;
790
791 for(Int_t cath=0; cath<fgkNcathodes; cath++){
792 for(Int_t ch=0; ch<fgkNchambers; ch++){
793 for(Int_t slat=0; slat<fgkNslats; slat++){
794 CalculateEfficiency(fHitPerSlat[ch][cath][slat], fHitPerSlat[ch][cath][slat]+fInefficientSlat[ch][cath][slat], efficiency, efficiencyError, kFALSE);
795 histo[cath][ch]->SetBinContent(slat+1, efficiency);
796 histo[cath][ch]->SetBinError(slat+1, efficiencyError);
797 }
798 CalculateEfficiency(fTrigger44[cath], fTrigger34[ch][cath]+fTrigger44[cath], efficiency, efficiencyError, kFALSE);
799 histo[cath][chamberEff]->SetBinContent(ch+1, efficiency);
800 histo[cath][chamberEff]->SetBinError(ch+1, efficiencyError);
801
802 for(Int_t board=0; board<fgkNboards; board++){
803 CalculateEfficiency(fHitPerBoard[ch][cath][board], fHitPerBoard[ch][cath][board]+fInefficientBoard[ch][cath][board], efficiency, efficiencyError, kFALSE);
804 histoBoard[cath][ch]->SetBinContent(board+1, efficiency);
805 histoBoard[cath][ch]->SetBinError(board+1, efficiencyError);
806 }
807 }
808 }
809
810// write all histos
811 outputHistoFile->cd();
812 dir->GetList()->Write();
813 outputHistoFile->Close();
814}
815
816
817//_____________________________________________________________________________
818void AliMUONTriggerChamberEff::WriteEfficiencyMap(const char* outputDir)
819{
820 Int_t effOutWidth=4;
821
822 Float_t efficiency, efficiencyError;
823
824 Int_t aCapo[] = {16, 38, 60, 76, 92, 108, 117, 133, 155, 177, 193, 209, 225, 234};
825
826 char filename[70];
827 sprintf(filename, "%s/efficiencyCells.dat", outputDir);
828 ofstream outFile(filename);
829 outFile << "localBoards" << endl;
830 for(Int_t ch=0; ch<fgkNchambers; ch++){
831 //Print information
832 outFile << "\n\ndetElemId:\t" << 11+ch;
833 outFile << "00";
834 for(Int_t cath=0; cath<fgkNcathodes; cath++){
835 outFile << "\n cathode:\t" << cath << endl;
836 Int_t currLine=0;
837 for(Int_t board=0; board<fgkNboards; board++){
838
839 if(board==aCapo[currLine]){
840 outFile << endl;
841 currLine++;
842 }
843 CalculateEfficiency(fHitPerBoard[ch][cath][board], fHitPerBoard[ch][cath][board]+fInefficientBoard[ch][cath][board], efficiency, efficiencyError, kFALSE);
844 outFile << " " << setw(effOutWidth) << efficiency;
845 }// loop on boards
846 outFile << endl;
847 }// loop on cathodes
848 }// loop on chambers
849}
850