remoe duplicate QA initialisation and do ESD QA for same detectors as RecPoint QA
[u/mrichter/AliRoot.git] / MUON / AliMUONQADataMakerRec.cxx
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
17 // --- ROOT system ---
18 #include <TClonesArray.h>
19 #include <TFile.h> 
20 #include <TH1F.h> 
21 #include <TH1I.h> 
22 #include <TH2F.h>
23 #include <TH3F.h> 
24 #include <TLorentzVector.h>
25
26 // --- AliRoot header files ---
27 #include "AliESDEvent.h"
28 #include "AliMUONConstants.h"  
29 #include "AliLog.h"
30 #include "AliRawReader.h"
31 #include "AliQAChecker.h"
32 #include "AliMpBusPatch.h"
33 #include "AliMUONCluster.h"  
34 #include "AliMUONRawStreamTracker.h"
35 #include "AliMUONRawStreamTrigger.h"
36
37 #include "AliMUONVClusterStore.h"
38 #include "AliMUONVCluster.h"
39 #include "AliESDMuonTrack.h"
40 #include "AliESDMuonCluster.h"
41
42 #include "AliMUONDigitMaker.h"
43 #include "AliMUONVDigitStore.h"
44 #include "AliMUONVTriggerStore.h"
45 #include "AliMUONVDigit.h"
46 #include "AliMUONLocalTrigger.h"
47
48 #include "AliMUONDDLTrigger.h"
49 #include "AliMUONRegHeader.h"
50 #include "AliMUONDarcHeader.h"
51 #include "AliMpTriggerCrate.h"
52 #include "AliMpDDLStore.h"
53 #include "AliMUONLocalStruct.h"
54 #include "AliMpLocalBoard.h"
55 #include "AliMpCDB.h"
56
57 #include "AliMpPad.h"
58 #include "AliMUONGeometryTransformer.h"
59 #include "AliMpVSegmentation.h"
60 #include "AliMpSegmentation.h"
61 #include "AliMpConstants.h"
62
63 #include "AliMUONQADataMakerRec.h"
64
65 //-----------------------------------------------------------------------------
66 /// \class AliMUONQADataMakerRec
67 ///
68 /// MUON base class for quality assurance data (histo) maker
69 ///
70 /// \author C. Finck
71
72 /// \cond CLASSIMP
73 ClassImp(AliMUONQADataMakerRec)
74 /// \endcond
75            
76 //____________________________________________________________________________ 
77 AliMUONQADataMakerRec::AliMUONQADataMakerRec() : 
78     AliQADataMakerRec(AliQA::GetDetName(AliQA::kMUON), "MUON Quality Assurance Data Maker"),
79     fDigitStore(0x0),
80     fTriggerStore(0x0),
81     fDigitMaker(0x0)
82 {
83     /// ctor
84   fDigitStore = AliMUONVDigitStore::Create("AliMUONDigitStoreV1");
85   fDigitMaker = new AliMUONDigitMaker(kTRUE,kFALSE);
86
87 }
88
89 //____________________________________________________________________________ 
90 AliMUONQADataMakerRec::AliMUONQADataMakerRec(const AliMUONQADataMakerRec& qadm) :
91     AliQADataMakerRec(qadm),
92     fDigitStore(0x0),
93     fTriggerStore(0x0),
94     fDigitMaker(0x0)
95 {
96     ///copy ctor 
97     SetName((const char*)qadm.GetName()) ; 
98     SetTitle((const char*)qadm.GetTitle()); 
99
100     // Do not copy the digit store and digit maker, but create its own ones
101     fDigitStore = AliMUONVDigitStore::Create("AliMUONDigitStoreV1");
102     fDigitMaker = new AliMUONDigitMaker(kTRUE,kFALSE);
103 }
104
105 //__________________________________________________________________
106 AliMUONQADataMakerRec& AliMUONQADataMakerRec::operator = (const AliMUONQADataMakerRec& qadm )
107 {
108   /// Assignment operator
109
110   // check assignment to self
111   if (this == &qadm) return *this;
112
113   this->~AliMUONQADataMakerRec();
114   new(this) AliMUONQADataMakerRec(qadm);
115   return *this;
116 }
117
118 //__________________________________________________________________
119 AliMUONQADataMakerRec::~AliMUONQADataMakerRec()
120 {
121     /// dtor
122     delete fDigitStore;
123     delete fTriggerStore;
124     delete fDigitMaker;
125 }
126
127 //____________________________________________________________________________ 
128 void AliMUONQADataMakerRec::EndOfDetectorCycle(AliQA::TASKINDEX_t task, TObjArray* list)
129 {
130     ///Detector specific actions at end of cycle
131
132     // Display trigger histos in a more user friendly way
133     DisplayTriggerInfo(task);
134
135     // do the QA checking
136     AliQAChecker::Instance()->Run(AliQA::kMUON, task, list) ;
137 }
138
139 //____________________________________________________________________________ 
140 void AliMUONQADataMakerRec::InitRaws()
141 {
142     /// create Raws histograms in Raws subdir
143     TH1I* h0 = new TH1I("hRawBusPatch", "buspatch distribution",  1932, 1, 1932); 
144     Add2RawsList(h0, kRawBusPatch);
145
146     TH1I* h1 = new TH1I("hRawCharge", "Charge distribution in rawdata", 4096, 0, 4095); 
147     Add2RawsList(h1, kRawCharge);
148                 
149     for (Int_t iDDL = 0; iDDL < 20; ++iDDL) 
150     {
151       TH1F* h2 = new TH1F(Form("%s%d", "hRawBusPatchDDL", iDDL), Form("%s %d","RAW Buspatch distribution for DDL", iDDL), 50, 0, 49); 
152       Add2RawsList(h2, kRawBuspatchDDL+iDDL);
153     }
154
155     TH3F* h3 = new TH3F("hTriggerScalersBendPlane", "Trigger scalers in bending plane",
156                         4, 10.5, 14.5,
157                         234, 0.5, 234.5,
158                         16, -0.5, 15.5);
159     h3->GetXaxis()->SetTitle("Chamber");
160     h3->GetYaxis()->SetTitle("Board");
161     h3->GetZaxis()->SetTitle("Strip");
162     Add2RawsList(h3, kTriggerScalersBP);
163
164     TH3F* h4 = new TH3F("hTriggerScalersNonBendPlane", "Trigger scalers in non-bending plane",
165                         4, 10.5, 14.5,
166                         234, 0.5, 234.5,
167                         16, -0.5, 15.5);
168     h4->GetXaxis()->SetTitle("Chamber");
169     h4->GetYaxis()->SetTitle("Board");
170     h4->GetZaxis()->SetTitle("Strip");
171     Add2RawsList(h4, kTriggerScalersNBP);
172
173     InitDisplayHistos(AliQA::kRAWS);
174 }
175
176 //____________________________________________________________________________ 
177 void AliMUONQADataMakerRec::InitRecPoints()
178 {
179     /// create Reconstructed Points histograms in RecPoints subdir
180     TH3F* h0 = new TH3F("hTriggerDigitsBendPlane", "Trigger digits in bending plane",
181                         4, 10.5, 14.5,
182                         234, 0.5, 234.5,
183                         16, -0.5, 15.5);
184     h0->GetXaxis()->SetTitle("Chamber");
185     h0->GetYaxis()->SetTitle("Board");
186     h0->GetZaxis()->SetTitle("Strip");
187     Add2RecPointsList(h0, kTriggerDigitsBendPlane);
188
189     TH3F* h1 = new TH3F("hTriggerDigitsNonBendPlane", "Trigger digits in non-bending plane",
190                         4, 10.5, 14.5,
191                         234, 0.5, 234.5,
192                         16, -0.5, 15.5);
193     h1->GetXaxis()->SetTitle("Chamber");
194     h1->GetYaxis()->SetTitle("Board");
195     h1->GetZaxis()->SetTitle("Strip");
196     Add2RecPointsList(h1, kTriggerDigitsNonBendPlane);
197
198     TH1F* h2 = new TH1F("hTriggeredBoards", "Triggered boards", 234, 0.5, 234.5);
199     Add2RecPointsList(h2, kTriggeredBoards);
200
201     InitDisplayHistos(AliQA::kRECPOINTS);
202 }
203
204
205 //____________________________________________________________________________ 
206 void AliMUONQADataMakerRec::InitESDs()
207 {
208     ///create ESDs histograms in ESDs subdir
209   TH1F* h0 = new TH1F("hESDnTracks", "ESDs track number distribution", 30, 0., 30.);  
210   Add2ESDsList(h0, kESDnTracks);
211
212   TH1F* h1 = new TH1F("hESDMomentum", "ESDs P distribution", 300, 0., 300) ; 
213   Add2ESDsList(h1, kESDMomentum);
214
215   TH1F* h2 = new TH1F("hESDPt", "ESDs Pt distribution", 200, 0., 50) ; 
216   Add2ESDsList(h2, kESDPt);
217
218   TH1F* h3 = new TH1F("hESDRapidity", "ESDs rapidity distribution", 200, -4.5,-2.) ; 
219   Add2ESDsList(h3, kESDRapidity);
220
221   for (Int_t i = 0; i < AliMUONConstants::NTrackingCh(); ++i) 
222   {
223     TH2F* h4 = new TH2F(Form("%s%d", "hESDClusterHitMap", i), 
224                      Form("%s %d", "ESD Clusters hit distribution for chamber", i),
225                      100, -1*AliMUONConstants::Rmax(i/2), AliMUONConstants::Rmax(i/2),
226                      100, -1*AliMUONConstants::Rmax(i/2), AliMUONConstants::Rmax(i/2)); 
227     Add2ESDsList(h4, kESDClusterHitMap+i);
228   }
229 }
230
231 //____________________________________________________________________________
232 void AliMUONQADataMakerRec::MakeRaws(AliRawReader* rawReader)
233 {
234     /// make QA for rawdata
235     Int_t busPatchId;
236     UShort_t manuId;  
237     UChar_t channelId;
238     UShort_t charge;
239
240     rawReader->Reset();
241     AliMUONRawStreamTracker rawStreamTrack(rawReader);
242     rawStreamTrack.First();
243     while( rawStreamTrack.Next(busPatchId, manuId, channelId, charge) ) {
244   
245       GetRawsData(kRawBusPatch)->Fill(busPatchId);
246       GetRawsData(kRawCharge)->Fill(charge);
247       Int_t iDDL = rawStreamTrack.GetCurentDDL();
248       GetRawsData(kRawBuspatchDDL + iDDL)->Fill(AliMpBusPatch::GetLocalBusID(busPatchId, iDDL));
249                 
250                   
251     } // Next digit
252
253
254     // Get trigger scalers
255
256     Int_t loCircuit=0;
257     AliMpCDB::LoadDDLStore();
258
259     AliMUONRawStreamTrigger rawStreamTrig(rawReader);
260     while (rawStreamTrig.NextDDL()) 
261     {
262       // If not a scaler event, do nothing
263       Bool_t scalerEvent =  rawReader->GetDataHeader()->GetL1TriggerMessage() & 0x1;
264       if(!scalerEvent) break;
265
266       AliMUONDDLTrigger* ddlTrigger = rawStreamTrig.GetDDLTrigger();
267       AliMUONDarcHeader* darcHeader = ddlTrigger->GetDarcHeader();
268
269       Int_t nReg = darcHeader->GetRegHeaderEntries();
270     
271       for(Int_t iReg = 0; iReg < nReg ;iReg++)
272       {   //reg loop
273
274         // crate info  
275         AliMpTriggerCrate* crate = AliMpDDLStore::Instance()->
276           GetTriggerCrate(rawStreamTrig.GetDDL(), iReg);
277
278         AliMUONRegHeader* regHeader =  darcHeader->GetRegHeaderEntry(iReg);
279
280         // loop over local structures
281         Int_t nLocal = regHeader->GetLocalEntries();
282         for(Int_t iLocal = 0; iLocal < nLocal; iLocal++) 
283         {
284           AliMUONLocalStruct* localStruct = regHeader->GetLocalEntry(iLocal);
285         
286           // if card exist
287           if (!localStruct) continue;
288           
289           loCircuit = crate->GetLocalBoardId(localStruct->GetId());
290
291           if ( !loCircuit ) continue; // empty slot
292
293           AliMpLocalBoard* localBoard = AliMpDDLStore::Instance()->GetLocalBoard(loCircuit, false);
294
295           // skip copy cards
296           if( !localBoard->IsNotified()) 
297             continue;
298
299           Int_t cathode = localStruct->GetComptXY()%2;
300
301           ERaw hindex = (cathode==0) ? kTriggerScalersBP : kTriggerScalersNBP;
302
303           // loop over strips
304           for (Int_t ibitxy = 0; ibitxy < 16; ++ibitxy) {
305             if(localStruct->GetXY1(ibitxy) > 0)
306               ((TH3F*)GetRawsData(hindex))->Fill(11+0, loCircuit, ibitxy, 2*localStruct->GetXY1(ibitxy));
307             if(localStruct->GetXY2(ibitxy) > 0)
308               ((TH3F*)GetRawsData(hindex))->Fill(11+1, loCircuit, ibitxy, 2*localStruct->GetXY2(ibitxy));
309             if(localStruct->GetXY3(ibitxy) > 0)
310               ((TH3F*)GetRawsData(hindex))->Fill(11+2, loCircuit, ibitxy, 2*localStruct->GetXY3(ibitxy));
311             if(localStruct->GetXY4(ibitxy) > 0)
312               ((TH3F*)GetRawsData(hindex))->Fill(11+3, loCircuit, ibitxy, 2*localStruct->GetXY4(ibitxy));
313           } // loop on strips
314         } // iLocal
315       } // iReg
316     } // NextDDL
317 }
318
319 //____________________________________________________________________________
320 void AliMUONQADataMakerRec::MakeRecPoints(TTree* clustersTree)
321 {
322   
323     /// makes data from trigger response
324       
325     // Fired pads info
326     fDigitStore->Clear();
327
328     if (!fTriggerStore) fTriggerStore = AliMUONVTriggerStore::Create(*clustersTree);
329     fTriggerStore->Clear();
330     fTriggerStore->Connect(*clustersTree, false);
331     clustersTree->GetEvent(0);
332
333     AliMUONLocalTrigger* locTrg;
334     TIter nextLoc(fTriggerStore->CreateLocalIterator());
335
336     while ( ( locTrg = static_cast<AliMUONLocalTrigger*>(nextLoc()) ) ) 
337     {
338       if (locTrg->IsNull()) continue;
339    
340       TArrayS xyPattern[2];
341       locTrg->GetXPattern(xyPattern[0]);
342       locTrg->GetYPattern(xyPattern[1]);
343
344       Int_t nBoard = locTrg->LoCircuit();
345
346       Bool_t xTrig=locTrg->IsTrigX();
347       Bool_t yTrig=locTrg->IsTrigY();
348     
349       if (xTrig && yTrig)
350         ((TH1F*)GetRecPointsData(kTriggeredBoards))->Fill(nBoard);
351
352       fDigitMaker->TriggerDigits(nBoard, xyPattern, *fDigitStore);
353     }
354
355     TIter nextDigit(fDigitStore->CreateIterator());
356     AliMUONVDigit* mDigit;
357     while ( ( mDigit = static_cast<AliMUONVDigit*>(nextDigit()) ) )
358     {
359       Int_t detElemId = mDigit->DetElemId();
360       Int_t ch = detElemId/100;
361       Int_t localBoard = mDigit->ManuId();
362       Int_t channel = mDigit->ManuChannel();
363       Int_t cathode = mDigit->Cathode();
364       ERecPoints hindex 
365         = ( cathode == 0 ) ? kTriggerDigitsBendPlane : kTriggerDigitsNonBendPlane;
366       
367       ((TH3F*)GetRecPointsData(hindex))->Fill(ch, localBoard, channel);
368     }
369 }
370
371 //____________________________________________________________________________
372 void AliMUONQADataMakerRec::MakeESDs(AliESDEvent* esd)
373 {
374     /// make QA data from ESDs
375
376     TLorentzVector v1;
377
378     Int_t nTracks = (Int_t)esd->GetNumberOfMuonTracks() ; 
379     GetESDsData(0)->Fill(nTracks);
380
381     for (Int_t iTrack = 0; iTrack < nTracks; ++iTrack) {
382
383       AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
384       muonTrack->LorentzP(v1);
385
386       GetESDsData(1)->Fill(v1.P());
387       GetESDsData(2)->Fill(v1.Pt());
388       GetESDsData(3)->Fill(v1.Rapidity());
389
390       TClonesArray clusters =  muonTrack->GetClusters();
391
392       for (Int_t iCluster = 0; iCluster <clusters.GetEntriesFast(); ++iCluster) {
393         AliESDMuonCluster* cluster = (AliESDMuonCluster*)clusters.At(iCluster);
394         GetESDsData(kESDClusterHitMap+cluster->GetChamberId())
395           ->Fill(cluster->GetX(), cluster->GetY());
396       }
397     }
398 }
399
400 //____________________________________________________________________________ 
401 void AliMUONQADataMakerRec::StartOfDetectorCycle()
402 {
403     /// Detector specific actions at start of cycle
404   
405 }
406
407 //____________________________________________________________________________ 
408 void AliMUONQADataMakerRec::DisplayTriggerInfo(AliQA::TASKINDEX_t task)
409 {
410   //
411   /// Display trigger information in a user-friendly way:
412   /// from local board and strip numbers to their position on chambers
413   //
414   if(task!=AliQA::kRECPOINTS && task!=AliQA::kRAWS) return;
415
416   AliMpCDB::LoadDDLStore();
417
418   Int_t detElemId;
419   Float_t xg1, yg1, zg1, xg2, yg2, zg2, binContent;
420   Float_t xLocal1=0., yLocal1=0., xLocal2=0., yLocal2=0., xWidth=0., yWidth=0.;
421   Float_t x1Label, y1Label, x2Label, y2Label;
422   Int_t x1Int, x2Int, y1Int, y2Int;
423   Int_t binCh, binBoard, binStrip;
424
425   AliMUONGeometryTransformer transform;
426   transform.LoadGeometryData();
427
428   TH3F* histoStrips=0x0;
429   TH2F* histoDisplayStrips=0x0;
430   TH2F* histoDisplayBoards=0x0;
431   TH1F* histoBoards=0x0;
432
433   const Float_t kShift = 0.15;
434
435   for (Int_t iCath = 0; iCath < 2; ++iCath)
436   {
437
438     if(task==AliQA::kRECPOINTS){
439       ERecPoints hindex 
440         = ( iCath == 0 ) ? kTriggerDigitsBendPlane : kTriggerDigitsNonBendPlane;
441       histoStrips = (TH3F*)GetRecPointsData(hindex);
442       histoBoards = (TH1F*)GetRecPointsData(kTriggeredBoards);
443     }
444     else if(task==AliQA::kRAWS){
445       ERaw hindex 
446         = ( iCath == 0 ) ? kTriggerScalersBP : kTriggerScalersNBP;
447       histoStrips = (TH3F*)GetRawsData(hindex);
448       if(histoStrips->GetEntries()==0) return; // No scalers found
449     }
450     
451
452     for (Int_t iChamber = 0; iChamber < 4; ++iChamber)
453     {
454       Int_t iCh = iChamber + AliMpConstants::NofTrackingChambers();
455
456       if(task==AliQA::kRECPOINTS){
457         histoDisplayStrips = (TH2F*)GetRecPointsData(kTriggerDigitsDisplay + 4*iCath + iChamber);
458         histoDisplayBoards = (TH2F*)GetRecPointsData(kTriggerBoardsDisplay);
459       }
460       else if(task==AliQA::kRAWS){
461         histoDisplayStrips = (TH2F*)GetRawsData(kTriggerScalersDisplay + 4*iCath + iChamber);
462       }
463
464       for(Int_t iBoard=1; iBoard<=234; iBoard++){
465         // get detElemId
466         detElemId = AliMpDDLStore::Instance()->GetDEfromLocalBoard(iBoard, iCh);
467
468         if (!detElemId) continue;
469
470         AliMpLocalBoard* localBoard = AliMpDDLStore::Instance()->GetLocalBoard(iBoard, false);
471
472         // skip copy cards
473         if( !localBoard->IsNotified()) 
474           continue;
475
476         const AliMpVSegmentation* seg 
477           = AliMpSegmentation::Instance()
478           ->GetMpSegmentation(detElemId, AliMp::GetCathodType(iCath));  
479         
480         // loop over strips
481         for (Int_t ibitxy = 0; ibitxy < 16; ++ibitxy) 
482         {
483           AliMpPad pad = seg->PadByLocation(AliMpIntPair(iBoard,ibitxy),kFALSE);
484                         
485           if (!pad.IsValid())
486             continue;
487
488           if(iCath==0){ // Geometry info from bending plane only
489             if(ibitxy==0) {
490               xLocal1 = pad.Position().X();
491               yLocal1 = pad.Position().Y();
492               xWidth = pad.Dimensions().X();
493               yWidth = pad.Dimensions().Y();
494             }
495             xLocal2 = pad.Position().X();
496             yLocal2 = pad.Position().Y();
497           }
498
499           Float_t dimX = pad.Dimensions().X();
500           Float_t dimY = pad.Dimensions().Y();
501
502           Float_t stripX = pad.Position().X();
503           Float_t stripY = pad.Position().Y();
504                    
505           transform.Local2Global(detElemId, stripX, stripY, 0, xg1, yg1, zg1);
506
507           x1Int = histoDisplayStrips->GetXaxis()->FindBin(xg1 - dimX + kShift);
508           y1Int = histoDisplayStrips->GetYaxis()->FindBin(yg1 - dimY + kShift);
509           x2Int = histoDisplayStrips->GetXaxis()->FindBin(xg1 + dimX - kShift);
510           y2Int = histoDisplayStrips->GetYaxis()->FindBin(yg1 + dimY - kShift);
511
512           binCh = histoStrips->GetXaxis()->FindBin(iCh+1);
513           binBoard = histoStrips->GetYaxis()->FindBin(iBoard);
514           binStrip = histoStrips->GetZaxis()->FindBin(ibitxy);
515           binContent = histoStrips->GetBinContent(binCh, binBoard, binStrip);
516
517           if(binContent==0) continue;
518
519           for(Int_t binX=x1Int; binX<=x2Int; binX++){
520             for(Int_t binY=y1Int; binY<=y2Int; binY++){       
521               histoDisplayStrips->SetBinContent(binX, binY, binContent);
522             }
523           }
524         }// ibitxy
525
526         if(task != AliQA::kRECPOINTS) continue;
527         if(iCath>0 || iChamber>0) continue; // Geometry info from bending plane only
528         transform.Local2Global(detElemId, xLocal1, yLocal1, 0, xg1, yg1, zg1);
529         transform.Local2Global(detElemId, xLocal2, yLocal2, 0, xg2, yg2, zg2);
530
531         // Per board
532         x1Label = TMath::Min(xg1,xg2) - xWidth + kShift;
533         y1Label = TMath::Min(yg1,yg2) - yWidth + kShift;
534         x2Label = TMath::Max(xg1,xg2) + xWidth - kShift;
535         y2Label = TMath::Max(yg1,yg2) + yWidth - kShift;
536
537         x1Int = histoDisplayBoards->GetXaxis()->FindBin(x1Label);
538         y1Int = histoDisplayBoards->GetYaxis()->FindBin(y1Label);
539         x2Int = histoDisplayBoards->GetXaxis()->FindBin(x2Label);
540         y2Int = histoDisplayBoards->GetYaxis()->FindBin(y2Label);
541
542         binBoard = histoBoards->GetXaxis()->FindBin(iBoard);
543         binContent = histoBoards->GetBinContent(binBoard);
544
545         if(binContent==0) continue;
546
547         for(Int_t binX=x1Int; binX<=x2Int; binX++){
548           for(Int_t binY=y1Int; binY<=y2Int; binY++){
549             histoDisplayBoards->SetBinContent(binX, binY, binContent);
550           }
551         }
552       } // iBoard
553     } // iChamber
554   } // iCath
555 }
556
557
558 //____________________________________________________________________________ 
559 void AliMUONQADataMakerRec::InitDisplayHistos(AliQA::TASKINDEX_t task)
560 {
561   //
562   /// Initialize trigger information display histograms,
563   /// using mapping and geometry
564   //
565   AliMpCDB::LoadDDLStore();
566
567   const Int_t kNumOfBoards = AliMpConstants::NofLocalBoards();
568
569   AliMUONGeometryTransformer transform;
570   transform.LoadGeometryData();
571
572   TString cathCode[2] = {"BendPlane", "NonBendPlane"};
573
574   TArrayF xLocal1[4], yLocal1[4], xLocal2[4], yLocal2[4], xWidth[4], yWidth[4];
575
576   TArrayF xAxisStrip[8];
577   TArrayF yAxisStrip[8];
578   TArrayF xAxisBoard[8];
579   TArrayF yAxisBoard[8];
580
581   TH2F* histoStrips=0x0;
582   TH2F* histoBoards=0x0;
583
584   const Float_t kResetValue=1234567.;
585    
586   for(Int_t ch=0; ch<4; ch++){
587     xLocal1[ch].Set(kNumOfBoards);
588     yLocal1[ch].Set(kNumOfBoards);
589     xLocal2[ch].Set(kNumOfBoards);
590     yLocal2[ch].Set(kNumOfBoards);
591     xWidth[ch].Set(kNumOfBoards);
592     yWidth[ch].Set(kNumOfBoards);
593   }
594
595   for(Int_t cath=0; cath<2; cath++){
596     for(Int_t ch=0; ch<4; ch++){
597       Int_t chCath = 4*cath + ch;
598       xAxisBoard[chCath].Set(60);
599       xAxisBoard[chCath].Reset(kResetValue);
600       yAxisBoard[chCath].Set(60);
601       yAxisBoard[chCath].Reset(kResetValue);
602
603       xAxisStrip[chCath].Set(700);
604       xAxisStrip[chCath].Reset(kResetValue);
605       yAxisStrip[chCath].Set(700);
606       yAxisStrip[chCath].Reset(kResetValue);
607     }
608   }
609    
610   Float_t xg1, yg1, zg1, xg2, yg2, zg2;
611
612   TString histoName, histoTitle;
613
614   const Float_t kShift = 0.;
615
616   for(Int_t iCath=0; iCath<2; iCath++){
617     for (Int_t iChamber = 0; iChamber < 4; ++iChamber) {
618       Int_t iCh = iChamber + AliMpConstants::NofTrackingChambers();
619       for(Int_t iLoc = 0; iLoc < 234; iLoc++) {  
620         Int_t iBoard = iLoc+1;
621         Int_t detElemId = AliMpDDLStore::Instance()->GetDEfromLocalBoard(iBoard, iCh);
622
623         if (!detElemId) continue;
624
625         AliMpLocalBoard* localBoard = AliMpDDLStore::Instance()->GetLocalBoard(iBoard, kFALSE);
626
627         // skip copy cards
628         if( !localBoard->IsNotified()) 
629           continue;
630
631         // get segmentation
632         const AliMpVSegmentation* seg = AliMpSegmentation::Instance()
633           ->GetMpSegmentation(detElemId,
634                               AliMp::GetCathodType(iCath));
635
636         Int_t chCath = 4*iCath + iChamber;
637         // loop over strips
638         for (Int_t ibitxy = 0; ibitxy < 16; ++ibitxy) {
639           // get pad from electronics
640           AliMpPad pad = seg->PadByLocation(AliMpIntPair(iBoard,ibitxy),kFALSE);
641
642           if (!pad.IsValid()) continue;
643
644           if(iCath==0){
645             if(ibitxy==0) {
646               xLocal1[iChamber][iLoc] = pad.Position().X();
647               yLocal1[iChamber][iLoc] = pad.Position().Y();
648               xWidth[iChamber][iLoc] = pad.Dimensions().X();
649               yWidth[iChamber][iLoc] = pad.Dimensions().Y();
650             }
651             xLocal2[iChamber][iLoc] = pad.Position().X();
652             yLocal2[iChamber][iLoc] = pad.Position().Y();
653           }
654
655           // Check fired pads
656           Float_t dimX = pad.Dimensions().X();
657           Float_t dimY = pad.Dimensions().Y();
658      
659           Float_t stripX = pad.Position().X();
660           Float_t stripY = pad.Position().Y();
661
662           transform.Local2Global(detElemId, stripX, stripY, 0, xg1, yg1, zg1);
663
664           AddSortedPoint(xg1 - dimX + kShift, xAxisStrip[chCath], kResetValue);
665           AddSortedPoint(xg1 + dimX - kShift, xAxisStrip[chCath], kResetValue);
666
667           AddSortedPoint(yg1 - dimY + kShift, yAxisStrip[chCath], kResetValue);
668           AddSortedPoint(yg1 + dimY - kShift, yAxisStrip[chCath], kResetValue);
669         } // loop on strips  
670
671         transform.Local2Global(detElemId, xLocal1[iChamber][iLoc], yLocal1[iChamber][iLoc], 0, xg1, yg1, zg1);
672         transform.Local2Global(detElemId, xLocal2[iChamber][iLoc], yLocal2[iChamber][iLoc], 0, xg2, yg2, zg2);
673
674         // Per board
675         AddSortedPoint(TMath::Min(xg1,xg2) - xWidth[iChamber][iLoc] + kShift, xAxisBoard[chCath], kResetValue);
676         AddSortedPoint(TMath::Max(xg1,xg2) + xWidth[iChamber][iLoc] - kShift, xAxisBoard[chCath], kResetValue);
677
678         AddSortedPoint(TMath::Min(yg1,yg2) - yWidth[iChamber][iLoc] + kShift, yAxisBoard[chCath], kResetValue);
679         AddSortedPoint(TMath::Max(yg1,yg2) + yWidth[iChamber][iLoc] - kShift, yAxisBoard[chCath], kResetValue);
680       } // loop on local boards 
681     } // loop on chambers
682   } // loop on cathodes
683
684   const Float_t kMinDiff = 0.1;
685
686   // Book histos
687   for(Int_t iCath=0; iCath<2; iCath++){
688     for(Int_t iChamber=0; iChamber<4; iChamber++){
689       Int_t chCath = 4*iCath + iChamber;
690       Int_t ipoint=0;
691       while(TMath::Abs(xAxisStrip[chCath][ipoint]-kResetValue)>kMinDiff) { ipoint++; }
692       xAxisStrip[chCath].Set(ipoint);
693
694       ipoint = 0;
695       while(TMath::Abs(yAxisStrip[chCath][ipoint]-kResetValue)>kMinDiff) { ipoint++; }
696       yAxisStrip[chCath].Set(ipoint);
697
698       ipoint = 0;
699       while(TMath::Abs(xAxisBoard[chCath][ipoint]-kResetValue)>kMinDiff) { ipoint++; }
700       xAxisBoard[chCath].Set(ipoint);
701
702       ipoint = 0;
703       while(TMath::Abs(yAxisBoard[chCath][ipoint]-kResetValue)>kMinDiff) { ipoint++; }
704       yAxisBoard[chCath].Set(ipoint);
705
706       if(task==AliQA::kRECPOINTS){
707         histoName = Form("hTriggerDigits%sChamber%i", cathCode[iCath].Data(), 11+iChamber);
708         histoTitle = Form("Chamber %i: Fired pads %s", 11+iChamber, cathCode[iCath].Data());
709       }
710       else if(task==AliQA::kRAWS){
711         histoName = Form("hScalers%sChamber%i", cathCode[iCath].Data(), 11+iChamber);
712         histoTitle = Form("Chamber %i: Scalers %s", 11+iChamber, cathCode[iCath].Data());
713       }
714       histoStrips = new TH2F(histoName.Data(), histoTitle.Data(), 
715                              xAxisStrip[chCath].GetSize()-1, xAxisStrip[chCath].GetArray(),
716                              yAxisStrip[chCath].GetSize()-1, yAxisStrip[chCath].GetArray());
717       histoStrips->SetXTitle("X (cm)");
718       histoStrips->SetYTitle("Y (cm)");
719
720       if(task==AliQA::kRECPOINTS){
721         Add2RecPointsList(histoStrips, kTriggerDigitsDisplay + 4*iCath + iChamber);
722       }
723       else if(task==AliQA::kRAWS){
724         Add2RawsList(histoStrips, kTriggerScalersDisplay + 4*iCath + iChamber);
725       }
726
727       if(task != AliQA::kRECPOINTS) continue;
728       if(iCath>0 || iChamber>0) continue;
729
730       histoName = "hFiredBoardsDisplay";
731       histoTitle = "Fired boards";
732       histoBoards = new TH2F(histoName.Data(), histoTitle.Data(), 
733                              xAxisBoard[chCath].GetSize()-1, xAxisBoard[chCath].GetArray(),
734                              yAxisBoard[chCath].GetSize()-1, yAxisBoard[chCath].GetArray());
735       histoBoards->SetXTitle("X (cm)");
736       histoBoards->SetYTitle("Y (cm)");
737
738       Add2RecPointsList(histoBoards, kTriggerBoardsDisplay + 4*iCath + iChamber);
739     } // loop on chamber
740   } // loop on cathode
741 }
742
743
744 //____________________________________________________________________________ 
745 Bool_t AliMUONQADataMakerRec::AddSortedPoint(Float_t currVal, TArrayF& position, const Float_t kResetValue)
746 {
747   //
748   /// Add sorted point in array according to an increasing order.
749   /// Used to build display histograms axis.
750   //
751   Int_t nEntries = position.GetSize();
752   Float_t tmp1, tmp2;
753   const Float_t kMinDiff = 0.1;
754   for(Int_t i=0; i<nEntries; i++){
755     if(TMath::Abs(position[i]-currVal)<kMinDiff) return kFALSE;
756     if(TMath::Abs(position[i]-kResetValue)<kMinDiff) {
757       position[i] = currVal;
758       return kTRUE;
759     }
760     if(currVal>position[i]) continue;
761     tmp1 = position[i];
762     position[i] = currVal;
763     for(Int_t j=i+1; j<nEntries; j++){
764       tmp2 = position[j];
765       position[j] = tmp1;
766       tmp1 = tmp2;
767       if(tmp1==kResetValue) break;
768     }
769     return kTRUE;
770   }
771   return kFALSE;
772 }
773