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