From Bogdan: new files / new implementation of the MUON visualization.
[u/mrichter/AliRoot.git] / EVE / Alieve / MUONData.cxx
CommitLineData
3626c858 1//
2// Sources:
3//
4// GetTrackerMapping = AliMUONDigitMaker::GetMapping
5// GetTriggerMapping = AliMUONDigitMaker::TriggerDigits
6// GetTriggerChamber = AliMUONDigitMaker::GetTriggerChamber
7// LoadRawTracker = MUONRawStreamTracker.C
8// LoadRawTrigger = MUONRawStreamTrigger.C
9//
10
11#include "MUONData.h"
12
13#include <Alieve/MUONChamberData.h>
14#include <Alieve/EventAlieve.h>
15
16#include <AliRawReader.h>
17#include <AliRawReaderFile.h>
18#include <AliRawReaderDate.h>
19#include <AliRawReaderRoot.h>
20
21#include <AliTracker.h>
22#include <AliMagFMaps.h>
23
24#include <AliMUONTrack.h>
25#include <AliMUONTrackParam.h>
26#include <AliMUONDigit.h>
27#include <AliMUONRawStreamTracker.h>
28#include <AliMUONRawStreamTrigger.h>
29#include <AliMUONDDLTracker.h>
30#include <AliMUONBlockHeader.h>
31#include <AliMUONDspHeader.h>
32#include <AliMUONBusStruct.h>
33#include <AliMUONDDLTrigger.h>
34#include <AliMUONDarcHeader.h>
35#include <AliMUONRegHeader.h>
36#include <AliMUONLocalStruct.h>
37#include <AliMUONTriggerCrateStore.h>
38#include <AliMUONTriggerCrate.h>
39#include <AliMUONLocalTriggerBoard.h>
40#include <AliMUONTriggerCircuit.h>
41#include <mapping/AliMpBusPatch.h>
42#include <mapping/AliMpVSegmentation.h>
43#include <mapping/AliMpSegmentation.h>
44#include <mapping/AliMpPad.h>
45#include <mapping/AliMpDEManager.h>
46
47#include "TTree.h"
48#include "TString.h"
49#include "TMatrixD.h"
50#include "TClonesArray.h"
51
52using namespace Reve;
53using namespace Alieve;
54
55//______________________________________________________________________
56// MUONData
57//
58
59ClassImp(MUONData)
60
61AliRawReader* MUONData::fgRawReader = 0;
62AliMUONRawStreamTracker* MUONData::fgRawStreamTracker = 0;
63AliMUONRawStreamTrigger* MUONData::fgRawStreamTrigger = 0;
64AliMpBusPatch* MUONData::fgBusPatchManager = 0;
65
66//______________________________________________________________________
67MUONData::MUONData() :
68 fChambers(14),
69 fNTracks(0),
70 fTrackPoints(0),
71 fNPoints(0)
72{
73 //
74 // Constructor
75 //
76
77 CreateAllChambers();
78
79}
80
81//______________________________________________________________________
82MUONData::~MUONData()
83{
84 //
85 // Destructor
86 //
87
88 DeleteAllChambers();
89
90 delete [] fTrackPoints;
91
92 fTrackPoints = 0;
93
94}
95
96//______________________________________________________________________
97MUONData::MUONData(const MUONData &mdata) :
98 TObject(mdata)
99{
100 //
101 // Copy constructor
102 //
103
104}
105
106//______________________________________________________________________
107MUONData& MUONData::operator=(const MUONData &mdata)
108{
109 //
110 // Assignment operator
111 //
112
113 if (this != &mdata) {
114
115 }
116
117 return *this;
118
119}
120
121//______________________________________________________________________
122void MUONData::CreateChamber(Int_t chamber)
123{
124 //
125 // create data for the chamber with id=chamber (0 to 13)
126 //
127
128 if (fChambers[chamber] == 0)
129 fChambers[chamber] = new MUONChamberData(chamber);
130
131}
132
133//______________________________________________________________________
134void MUONData::CreateAllChambers()
135{
136 //
137 // create all 14 chambers data
138 //
139
140 for (Int_t c = 0; c < 14; ++c)
141 CreateChamber(c);
142
143}
144
145//______________________________________________________________________
146void MUONData::DropAllChambers()
147{
148 //
149 // release data from all chambers
150 //
151
152 for (Int_t c = 0; c < 14; ++c) {
153
154 if (fChambers[c] != 0)
155 fChambers[c]->DropData();
156
157 }
158
159}
160
161//______________________________________________________________________
162void MUONData::DeleteAllChambers()
163{
164 //
165 // delete all chambers data
166 //
167
168 for (Int_t c = 0; c < 14; ++c) {
169
170 delete fChambers[c];
171 fChambers[c] = 0;
172
173 }
174
175}
176
177//______________________________________________________________________
178void MUONData::LoadTracks(TTree* tree)
179{
180 //
181 // load tracks from the TreeT and creates the track points array
182 // the structure of fTrackPoints:
183 // 0, 0, 0, - new track
184 // px, py, pz - from track param at vertex
185 // x, y, z - track param at vertex
186 // x, y, z - track param at hits
187 // ..........
188 //
189 // 0, 0, 0 - new track
190 // ..........
191 //
192
193 Int_t maxTrackHits = 1+2*10+4;
194
195 TClonesArray *tracks = 0;
196 tree->SetBranchAddress("MUONTrack",&tracks);
197 tree->GetEntry(0);
198
199 Int_t ntracks = tracks->GetEntriesFast();
200 printf("Found %d tracks. \n",ntracks);
201
202 fNTracks = ntracks;
203
204 Int_t maxTrackPoints = (3+3+3*maxTrackHits)*ntracks;
205 fTrackPoints = new Float_t [maxTrackPoints];
206
207 TMatrixD smatrix(2,2);
208 TMatrixD sums(2,1);
209 TMatrixD res(2,1);
210
211 Float_t xRec, xRec0;
212 Float_t yRec, yRec0;
213 Float_t zRec, zRec0;
214 Float_t px0, py0, pz0;
215
216 Float_t zg[4] = { -1603.5, -1620.5, -1703.5, -1720.5 };
217
218 AliMUONTrack *mt;
219 Int_t count = 0;
220 for (Int_t n = 0; n < ntracks; n++) {
221
222 if (count >= maxTrackPoints) continue;
223 fTrackPoints[3*count ] = 0.0;
224 fTrackPoints[3*count+1] = 0.0;
225 fTrackPoints[3*count+2] = 0.0;
226 count++;
227
228 mt = (AliMUONTrack*) tracks->At(n);
229
230 printf("Match trigger %d \n",mt->GetMatchTrigger());
231
232 AliMUONTrackParam *trackParam = mt->GetTrackParamAtVertex();
233 xRec0 = trackParam->GetNonBendingCoor();
234 yRec0 = trackParam->GetBendingCoor();
235 zRec0 = trackParam->GetZ();
236
237 px0 = trackParam->Px();
238 py0 = trackParam->Py();
239 pz0 = trackParam->Pz();
240
241 if (count >= maxTrackPoints) continue;
242 fTrackPoints[3*count ] = px0;
243 fTrackPoints[3*count+1] = py0;
244 fTrackPoints[3*count+2] = pz0;
245 count++;
246
247 if (count >= maxTrackPoints) continue;
248 fTrackPoints[3*count ] = xRec0;
249 fTrackPoints[3*count+1] = yRec0;
250 fTrackPoints[3*count+2] = zRec0;
251 count++;
252
253 Float_t xr[20], yr[20], zr[20];
254 for (Int_t i = 0; i < 10; i++) xr[i]=yr[i]=zr[i]=0.0;
255
256 Int_t nTrackHits = mt->GetNTrackHits();
257 printf("Nhits = %d \n",nTrackHits);
258 TClonesArray* trackParamAtHit;
259 for (Int_t iHit = 0; iHit < nTrackHits; iHit++){
260 trackParamAtHit = mt->GetTrackParamAtHit();
261 trackParam = (AliMUONTrackParam*) trackParamAtHit->At(iHit);
262 xRec = trackParam->GetNonBendingCoor();
263 yRec = trackParam->GetBendingCoor();
264 zRec = trackParam->GetZ();
265
266 //printf("Hit %d x %f y %f z %f \n",iHit,xRec,yRec,zRec);
267
268 xr[iHit] = xRec;
269 yr[iHit] = yRec;
270 zr[iHit] = zRec;
271
272 if (count >= maxTrackPoints) continue;
273 fTrackPoints[3*count ] = xRec;
274 fTrackPoints[3*count+1] = yRec;
275 fTrackPoints[3*count+2] = zRec;
276 count++;
277
278 }
279
280 Float_t xrc[20], yrc[20], zrc[20];
281 Int_t nrc = 0;
282 if (mt->GetMatchTrigger() && 1) {
283
284 for (Int_t i = 0; i < nTrackHits; i++) {
285 if (TMath::Abs(zr[i]) > 1000.0) {
286 //printf("Hit %d x %f y %f z %f \n",iHit,xr[i],yr[i],zr[i]);
287 xrc[nrc] = xr[i];
288 yrc[nrc] = yr[i];
289 zrc[nrc] = zr[i];
290 nrc++;
291 }
292 }
293
294 if (nrc < 2) continue;
295
296 Double_t xv, yv;
297 Float_t ax, bx, ay, by;
298
299 // fit x-z
300 smatrix.Zero();
301 sums.Zero();
302 for (Int_t i = 0; i < nrc; i++) {
303 xv = (Double_t)zrc[i];
304 yv = (Double_t)xrc[i];
305 //printf("x-z: xv %f yv %f \n",xv,yv);
306 smatrix(0,0) += 1.0;
307 smatrix(1,1) += xv*xv;
308 smatrix(0,1) += xv;
309 smatrix(1,0) += xv;
310 sums(0,0) += yv;
311 sums(1,0) += xv*yv;
312 }
313 res = smatrix.Invert() * sums;
314 ax = res(0,0);
315 bx = res(1,0);
316
317 // fit y-z
318 smatrix.Zero();
319 sums.Zero();
320 for (Int_t i = 0; i < nrc; i++) {
321 xv = (Double_t)zrc[i];
322 yv = (Double_t)yrc[i];
323 //printf("y-z: xv %f yv %f \n",xv,yv);
324 smatrix(0,0) += 1.0;
325 smatrix(1,1) += xv*xv;
326 smatrix(0,1) += xv;
327 smatrix(1,0) += xv;
328 sums(0,0) += yv;
329 sums(1,0) += xv*yv;
330 }
331 res = smatrix.Invert() * sums;
332 ay = res(0,0);
333 by = res(1,0);
334
335 Float_t xtc, ytc, ztc;
336 for (Int_t ii = 0; ii < 4; ii++) {
337
338 ztc = zg[ii];
339 ytc = ay+by*zg[ii];
340 xtc = ax+bx*zg[ii];
341
342 //printf("tc: x %f y %f z %f \n",xtc,ytc,ztc);
343
344 if (count >= maxTrackPoints) continue;
345 fTrackPoints[3*count ] = xtc;
346 fTrackPoints[3*count+1] = ytc;
347 fTrackPoints[3*count+2] = ztc;
348 count++;
349
350 }
351
352 } // end match trigger
353
354 }
355
356 fNPoints = 3*count;
357
358 printf("MUONData found %d track points. \n",fNPoints);
359
360}
361
362//______________________________________________________________________
363void MUONData::LoadDigits(TTree* tree)
364{
365 //
366 // load digits from the TreeD
367 //
368
369 Char_t branchname[30];
370 TClonesArray *digits = 0;
371 Int_t ndigits;
372 AliMUONDigit *mdig;
373 Int_t cathode, detElemId, ix, iy, charge;
374
375 for (Int_t c = 0; c < 14; ++c) {
376
377 if (fChambers[c] == 0) continue;
378 sprintf(branchname,"MUONDigits%d",c+1);
379 tree->SetBranchAddress(branchname,&digits);
380 tree->GetEntry(0);
381
382 ndigits = digits->GetEntriesFast();
383
384 for (Int_t id = 0; id < ndigits; id++) {
385 mdig = (AliMUONDigit*)digits->UncheckedAt(id);
386
387 cathode = mdig->Cathode();
388 ix = mdig->PadX();
389 iy = mdig->PadY();
390 detElemId = mdig->DetElemId();
391 charge = (Int_t)mdig->Signal();
392
393 if (c > 9) {
394 //printf("cha %d deid %d cath %1d ix %d iy %d q %d \n",c,detElemId,cathode,ix,iy,charge);
395 }
396
397 fChambers[c]->RegisterDigit(detElemId,cathode,ix,iy,charge);
398
399 } // end digits loop
400
401 }
402
403}
404
405//______________________________________________________________________
406void MUONData::LoadRaw(TString fileName)
407{
408 //
409 // load raw data from fileName; tracker and trigger data
410 //
411
412 if (fgRawReader == 0) {
413 // check extention to choose the rawdata file format
414 if (fileName.EndsWith("/")) {
415 fgRawReader = new AliRawReaderFile(fileName); // DDL files
416 } else if (fileName.EndsWith(".root")) {
417 fgRawReader = new AliRawReaderRoot(fileName); // ROOT file
418 } else if (!fileName.IsNull()) {
419 fgRawReader = new AliRawReaderDate(fileName); // DATE file
420 }
421 fgRawStreamTracker = new AliMUONRawStreamTracker(fgRawReader);
422 fgRawStreamTrigger = new AliMUONRawStreamTrigger(fgRawReader);
423 fgBusPatchManager = new AliMpBusPatch();
424 fgBusPatchManager->ReadBusPatchFile();
425 }
426
427 LoadRawTracker();
428 LoadRawTrigger();
429
430}
431
432//______________________________________________________________________
433void MUONData::LoadRawTracker()
434{
435 //
436 // load raw data for the tracking chambers
437 //
438
439 fgRawReader->RewindEvents();
440
441 AliMUONDigit* digit = new AliMUONDigit();
442
443 Int_t maxEvent = 1000;
444 Int_t minDDL = 0, maxDDL = 19;
445 Int_t cathode, detElemId, ix, iy, iChamber;
446
447 AliMUONDDLTracker* ddlTracker = 0x0;
448 AliMUONBlockHeader* blkHeader = 0x0;
449 AliMUONDspHeader* dspHeader = 0x0;
450 AliMUONBusStruct* busStruct = 0x0;
451
452 Int_t iEvent = 0;
453 Int_t dataSize, buspatchId;
454
455 Event* aevent = Alieve::gEvent;
456
457 while (fgRawReader->NextEvent()) {
458
459 if (iEvent != aevent->GetEventId()) {
460 iEvent++;
461 continue;
462 }
463
464 if (iEvent == maxEvent)
465 break;
466
467 // read DDL while < 20 DDL
468 while(fgRawStreamTracker->NextDDL()) {
469
470 if (fgRawStreamTracker->GetDDL() < minDDL ||
471 fgRawStreamTracker->GetDDL() > maxDDL)
472 continue;
473
474 //printf("\niDDL %d\n", fgRawStreamTracker->GetDDL());
475
476 ddlTracker = fgRawStreamTracker->GetDDLTracker();
477
478 // loop over block structure
479 Int_t nBlock = ddlTracker->GetBlkHeaderEntries();
480 for(Int_t iBlock = 0; iBlock < nBlock ;iBlock++){
481
482 blkHeader = ddlTracker->GetBlkHeaderEntry(iBlock);
483 //printf("Block Total length %d\n",blkHeader->GetTotalLength());
484
485 // loop over DSP structure
486 Int_t nDsp = blkHeader->GetDspHeaderEntries();
487 for(Int_t iDsp = 0; iDsp < nDsp ;iDsp++){ //DSP loop
488
489 dspHeader = blkHeader->GetDspHeaderEntry(iDsp);
490 // printf("Dsp length %d even word %d\n",dspHeader->GetTotalLength(), dspHeader->GetEventWord());
491
492 // loop over BusPatch structure
493 Int_t nBusPatch = dspHeader->GetBusPatchEntries();
494 for(Int_t iBusPatch = 0; iBusPatch < nBusPatch; iBusPatch++) {
495
496 busStruct = dspHeader->GetBusPatchEntry(iBusPatch);
497
498 //printf("busPatchId %d", busStruct->GetBusPatchId());
499 //printf(" BlockId %d", busStruct->GetBlockId());
500 //printf(" DspId %d\n", busStruct->GetDspId());
501
502 // loop over data
503 dataSize = busStruct->GetLength();
504 buspatchId = busStruct->GetBusPatchId();
505 for (Int_t iData = 0; iData < dataSize; iData++) {
506
507 Int_t manuId = busStruct->GetManuId(iData);
508 Int_t channelId = busStruct->GetChannelId(iData);
509 Int_t charge = busStruct->GetCharge(iData);
510 //printf("manuId: %d, channelId: %d charge: %d\n", manuId, channelId, charge);
511 // set digit charge
512 digit->SetSignal(charge);
513 digit->SetPhysicsSignal(charge);
514 digit->SetADC(charge);
515 // Get Back the hits at pads
516 Int_t error;
517 error = GetTrackerMapping(buspatchId,manuId,channelId,digit);
518 if (error) {
519 printf("Mapping Error\n");
520 continue;
521 }
522
523 cathode = digit->Cathode();
524 ix = digit->PadX();
525 iy = digit->PadY();
526 detElemId = digit->DetElemId();
527 charge = (Int_t)digit->Signal();
528 iChamber = detElemId/100 - 1;
529
530 fChambers[iChamber]->RegisterDigit(detElemId,cathode,ix,iy,charge);
531
532 } // iData
533 } // iBusPatch
534 } // iDsp
535 } // iBlock
536 } // NextDDL
537
538 break;
539
540 } // end event loop
541
542 delete digit;
543
544}
545
546//______________________________________________________________________
547void MUONData::LoadRawTrigger()
548{
549 //
550 // load raw data for the trigger chambers
551 //
552
553 fgRawReader->RewindEvents();
554
555 Int_t maxEvent = 1000;
556 Int_t minDDL = 0, maxDDL = 1;
557 Int_t detElemId, iChamber, cathode, charge, ix, iy;
558
559 AliMUONDDLTrigger* ddlTrigger = 0x0;
560 AliMUONDarcHeader* darcHeader = 0x0;
561 AliMUONRegHeader* regHeader = 0x0;
562 AliMUONLocalStruct* localStruct = 0x0;
563
564 // crate manager
565 AliMUONTriggerCrateStore* crateManager = new AliMUONTriggerCrateStore();
566 crateManager->ReadFromFile();
567
568 // Loop over events
569 Int_t iEvent = 0;
570 TList digitList;
571
572 Event* aevent = Alieve::gEvent;
573
574 while (fgRawReader->NextEvent()) {
575
576 if (iEvent != aevent->GetEventId()) {
577 iEvent++;
578 continue;
579 }
580
581 if (iEvent == maxEvent)
582 break;
583
584 // read DDL while < 2 DDL
585 while(fgRawStreamTrigger->NextDDL()) {
586
587 if (fgRawStreamTrigger->GetDDL() < minDDL ||
588 fgRawStreamTrigger->GetDDL() > maxDDL)
589 continue;
590
591 //printf("\niDDL %d\n", fgRawStreamTrigger->GetDDL());
592
593 ddlTrigger = fgRawStreamTrigger->GetDDLTrigger();
594 darcHeader = ddlTrigger->GetDarcHeader();
595
596 //printf("Global output %x\n", (Int_t)darcHeader->GetGlobalOutput());
597
598 // loop over regional structures
599 Int_t nReg = darcHeader->GetRegHeaderEntries();
600 for(Int_t iReg = 0; iReg < nReg ;iReg++){ //REG loop
601
602 //printf("RegionalId %d\n", iReg);
603
604 regHeader = darcHeader->GetRegHeaderEntry(iReg);
605 // printf("Reg length %d\n",regHeader->GetHeaderLength());
606
607 // crate info
608 AliMUONTriggerCrate* crate = crateManager->Crate(fgRawStreamTrigger->GetDDL(), iReg);
609 TObjArray *boards = crate->Boards();
610
611 // loop over local structures
612 Int_t nLocal = regHeader->GetLocalEntries();
613 for(Int_t iLocal = 0; iLocal < nLocal; iLocal++) {
614
615 localStruct = regHeader->GetLocalEntry(iLocal);
616
617 // check if trigger
618 if (localStruct->GetTriggerY() == 0) { // no empty data
619
620 // local trigger circuit number
621 AliMUONLocalTriggerBoard* localBoard = (AliMUONLocalTriggerBoard*)boards->At(iLocal+1);
622
623 //printf("LocalId %d\n", localStruct->GetId());
624 /*
625 Int_t iLocCard = localBoard->GetNumber();
626 Int_t loStripX = (Int_t)localStruct->GetXPos();
627 Int_t loStripY = (Int_t)localStruct->GetYPos();
628 Int_t loDev = (Int_t)localStruct->GetXDev();
629 */
630 //printf("iLocCard: %d, XPos: %d, YPos: %d Dev: %d\n", iLocCard, loStripX, loStripY, loDev);
631
632 digitList.Clear();
633 if ( GetTriggerMapping(localBoard, localStruct, digitList) ) {
634 for (Int_t iEntry = 0; iEntry < digitList.GetEntries(); iEntry++) {
635
636 AliMUONDigit* digit = (AliMUONDigit*)digitList.At(iEntry);
637 cathode = digit->Cathode();
638 ix = digit->PadX();
639 iy = digit->PadY();
640 detElemId = digit->DetElemId();
641 charge = (Int_t)digit->Signal();
642 iChamber = detElemId/100 - 1;
643
644 //printf("cha %d deid %d cath %1d ix %d iy %d q %d \n",iChamber,detElemId,cathode,ix,iy,charge);
645
646 fChambers[iChamber]->RegisterDigit(detElemId,cathode,ix,iy,charge);
647
648 }
649
650 }
651
652 }
653 } // iLocal
654 } // iReg
655 } // NextDDL
656
657 break;
658
659 } // end event loop
660
661 delete crateManager;
662
663}
664
665//______________________________________________________________________
666Int_t MUONData::GetTrackerMapping(Int_t buspatchId, UShort_t manuId, UChar_t channelId, AliMUONDigit* digit)
667{
668 //
669 // decode digits mapping for the tracking chambers
670 //
671
672 // getting DE from buspatch
673 Int_t detElemId = fgBusPatchManager->GetDEfromBus(buspatchId);
674 //AliDebug(3,Form("detElemId: %d busPatchId %d\n", detElemId, buspatchId));
675
676 const AliMpVSegmentation* seg = AliMpSegmentation::Instance()->GetMpSegmentationByElectronics(detElemId, manuId);
677 AliMpPad pad = seg->PadByLocation(AliMpIntPair(manuId,channelId),kTRUE);
678
679 if (!pad.IsValid())
680 {
681 printf("No pad for detElemId: %d, busPatchId %d, manuId: %d, channelId: %d\n",detElemId, buspatchId, manuId, channelId);
682
683 return 1;
684 } // return error
685
686 // Getting padX, padY and cathode number.
687 Int_t padX = pad.GetIndices().GetFirst();
688 Int_t padY = pad.GetIndices().GetSecond();
689 Int_t iCath = AliMpDEManager::GetCathod(detElemId,seg->PlaneType());
690
691 // storing into digits
692 digit->SetPadX(padX);
693 digit->SetPadY(padY);
694 digit->SetCathode(iCath);
695 digit->SetDetElemId(detElemId);
696 digit->SetElectronics(manuId,channelId);
697
698 //printf("detElemId: %d, busPatchId %d, manuId: %d, channelId: %d, padx: %d pady %d\n",detElemId, buspatchId, manuId, channelId, padX, padY);
699
700 return 0;
701
702}
703
704//______________________________________________________________________
705Int_t MUONData::GetTriggerMapping(AliMUONLocalTriggerBoard* localBoard,
706 AliMUONLocalStruct* localStruct,
707 TList& digitList)
708{
709 //
710 // decode digits mapping for the trigger chambers
711 //
712
713 Int_t detElemId;
714 Int_t nBoard;
715 Int_t iCath = -1;
716 Int_t iChamber = 0;
717 Int_t xyPattern = 0;
718
719 // loop over x1-4 and y1-4
720 for (Int_t icase = 0; icase < 8; icase++) {
721
722 // get chamber, cathode and associated trigger response pattern
723 GetTriggerChamber(localStruct, xyPattern, iChamber, iCath, icase);
724
725 if (!xyPattern) continue;
726
727 // get detElemId
728 AliMUONTriggerCircuit triggerCircuit;
729 detElemId = triggerCircuit.DetElemId(iChamber, localBoard->GetName());
730 nBoard = localBoard->GetNumber();
731
732 const AliMpVSegmentation* seg = AliMpSegmentation::Instance()->GetMpSegmentation(detElemId, iCath);
733
734 // loop over the 16 bits of pattern
735 for (Int_t ibitxy = 0; ibitxy < 16; ibitxy++) {
736
737 if ((xyPattern >> ibitxy) & 0x1) {
738
739 // not quite sure about this
740 Int_t offset = 0;
741 if (iCath && localBoard->GetSwitch(6)) offset = -8;
742
743 AliMpPad pad = seg->PadByLocation(AliMpIntPair(nBoard,ibitxy+offset),kTRUE);
744
745 AliMUONDigit* digit = new AliMUONDigit();
746 if (!pad.IsValid()) {
747 AliWarning(Form("No pad for detElemId: %d, nboard %d, ibitxy: %d\n",
748 detElemId, nBoard, ibitxy));
749 continue;
750 } //
751
752 Int_t padX = pad.GetIndices().GetFirst();
753 Int_t padY = pad.GetIndices().GetSecond();
754
755 // file digit
756 digit->SetSignal(1);
757 digit->SetPadX(padX);
758 digit->SetPadY(padY);
759 digit->SetCathode(iCath);
760 digit->SetDetElemId(detElemId);
761 digit->SetElectronics(nBoard, ibitxy);
762 digitList.Add(digit);
763
764 }// xyPattern
765 }// ibitxy
766 }// case
767
768 return 1;
769
770}
771
772//____________________________________________________________________
773void MUONData::GetTriggerChamber(AliMUONLocalStruct* localStruct, Int_t& xyPattern, Int_t& iChamber, Int_t& iCath, Int_t icase)
774{
775 //
776 // extract digits pattern
777 //
778
779 // get chamber & cathode number, (chamber starts at 0 !)
780 switch(icase) {
781 case 0:
782 xyPattern = localStruct->GetX1();
783 iCath = 0;
784 iChamber = 10;
785 break;
786 case 1:
787 xyPattern = localStruct->GetX2();
788 iCath = 0;
789 iChamber = 11;
790 break;
791 case 2:
792 xyPattern = localStruct->GetX3();
793 iCath = 0;
794 iChamber = 12;
795 break;
796 case 3:
797 xyPattern = localStruct->GetX4();
798 iCath = 0;
799 iChamber = 13;
800 break;
801 case 4:
802 xyPattern = localStruct->GetY1();
803 iCath = 1;
804 iChamber = 10;
805 break;
806 case 5:
807 xyPattern = localStruct->GetY2();
808 iCath = 1;
809 iChamber = 11;
810 break;
811 case 6:
812 xyPattern = localStruct->GetY3();
813 iCath = 1;
814 iChamber = 12;
815 break;
816 case 7:
817 xyPattern = localStruct->GetY4();
818 iCath = 1;
819 iChamber = 13;
820 break;
821 }
822
823}
824
825//______________________________________________________________________
826MUONChamberData* MUONData::GetChamberData(Int_t chamber)
827{
828 //
829 // return chamber data
830 //
831
832 if (chamber < 0 || chamber > 13) return 0;
833
834 //if (fChambers[chamber] == 0) CreateChamber(chamber);
835
836 return fChambers[chamber];
837
838}