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