1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 // --- MUON header files ---
19 #include "AliMUONTriggerDisplay.h"
21 #include "AliMpDDLStore.h"
22 #include "AliMpVSegmentation.h"
23 #include "AliMpSegmentation.h"
24 #include "AliMpConstants.h"
26 #include "AliMpLocalBoard.h"
29 // --- AliRoot header files ---
32 // --- ROOT system ---
36 //-----------------------------------------------------------------------------
37 /// \class AliMUONTriggerDisplay
39 /// MUON base class for converting histos as a function of strip/board/slat number
40 /// into display histos showing the detection element position in the
43 /// Input histos can be given as:
44 /// - TH2 with x -> board # [1-234] and y -> strip # in board [0-15]. Option: kDisplayStrips
45 /// - TH1 with x -> board # [1-234] Option: kDisplayBoards
46 /// - TH1 with x -> slat # [0-17] Option: kDisplaySlats
51 ClassImp(AliMUONTriggerDisplay)
54 //____________________________________________________________________________
55 AliMUONTriggerDisplay::AliMUONTriggerDisplay() :
61 //__________________________________________________________________
62 AliMUONTriggerDisplay::~AliMUONTriggerDisplay()
68 //____________________________________________________________________________
69 TH2* AliMUONTriggerDisplay::GetEmptyDisplayHisto(TString displayHistoName, EDisplayType displayType,
70 Int_t cathode, Int_t chamber,
71 TString displayHistoTitle)
74 /// Return the display histogram with optimized binning
75 /// but do not fill it.
77 TH2F* displayHisto = new TH2F();
79 InitOrDisplayTriggerInfo(0x0, displayHisto, displayType,
80 cathode, chamber, displayHistoName, displayHistoTitle);
86 //____________________________________________________________________________
87 TH2* AliMUONTriggerDisplay::GetBoardNumberHisto(TString displayHistoName,
89 TString displayHistoTitle)
92 /// Return the display histogram with optimized binning
93 /// and fill it with the board number
96 const Int_t kNboards = AliMpConstants::NofLocalBoards();
97 TH1F* inputHisto = new TH1F("boardNumbers","Board Numbers",kNboards,0.5,(Float_t)kNboards + 0.5);
98 for(Int_t ibin=1; ibin<=kNboards; ibin++){
99 inputHisto->Fill(ibin,ibin);
102 TH2F* displayHisto = (TH2F*)GetEmptyDisplayHisto(displayHistoName, kDisplayBoards, 0, chamber, displayHistoTitle);
103 FillDisplayHistogram(inputHisto,displayHisto,kDisplayBoards,0,chamber,kNumbered);
107 displayHisto->SetStats(kFALSE);
112 //____________________________________________________________________________
113 TH2* AliMUONTriggerDisplay::GetDisplayHistogram(TH1* inputHisto, TString displayHistoName,
114 EDisplayType displayType, Int_t cathode,
115 Int_t chamber, TString displayHistoTitle,
116 EDisplayOption displayOpt)
119 /// Get histogram displaying the information contained in the input histogram
121 TH2* displayHisto = GetEmptyDisplayHisto(displayHistoName, displayType,
122 cathode, chamber, displayHistoTitle);
124 FillDisplayHistogram(inputHisto, displayHisto, displayType, cathode, chamber, displayOpt);
129 //____________________________________________________________________________
130 Bool_t AliMUONTriggerDisplay::FillDisplayHistogram(TH1* inputHisto, TH2* displayHisto,
131 EDisplayType displayType, Int_t cathode,
132 Int_t chamber, EDisplayOption displayOpt)
135 /// Fill a previously initialized display histogram
136 /// with the information contained in inputHisto.
137 /// To get initialized display, please use GetEmptyDisplayHisto method
139 return InitOrDisplayTriggerInfo(inputHisto, displayHisto, displayType,
140 cathode, chamber, "", "",displayOpt);
143 //____________________________________________________________________________
144 Bool_t AliMUONTriggerDisplay::InitOrDisplayTriggerInfo(TH1* inputHisto, TH2* displayHisto,
145 EDisplayType displayType,
146 Int_t cathode, Int_t chamber,
147 TString displayHistoName, TString displayHistoTitle,
148 EDisplayOption displayOpt)
151 /// Initialize trigger information display histograms using mapping
152 /// Trigger information is displayed in a user-friendly way:
153 /// from local board and strip numbers to their position on chambers.
157 if ( ! AliMpSegmentation::Instance(kFALSE) ) {
159 if ( ! AliMpCDB::LoadDDLStore() ) {
160 AliError("Could not access mapping from OCDB !");
161 AliError("Histograms are not initialized !");
166 Int_t iCh = (chamber > AliMpConstants::NofTrackingChambers()) ? chamber - 11 : chamber;
167 Int_t iCath = cathode;
174 Float_t yOffsetLine, xOffsetLine = 0.;
176 const Float_t kResetValue=1234567.;
180 xAxisBoard.Reset(kResetValue);
182 yAxisBoard.Reset(kResetValue);
185 xAxisStrip.Reset(kResetValue);
187 yAxisStrip.Reset(kResetValue);
189 else if(!displayHisto){
190 AliWarning("Display histogram not initialized. Please initialize it first!");
194 Float_t xWidth, yWidth, yWidthSlat=0., xWidthCol=0.;
196 Float_t x1b=0., x2b=0., y1b=0., y2b=0.;
197 Float_t xcPad, ycPad;
201 const Float_t kShiftB = 0.5;
202 const Float_t kShiftS = 0.1;
204 const Float_t kShiftX = (iCath==0) ? kShiftB : kShiftS;
205 const Float_t kShiftY = (iCath==0) ? kShiftS : kShiftB;
206 const Float_t kShiftEl = (displayType==kDisplaySlats) ? 0.01 : kShiftB;
208 Int_t iChamber = iCh + AliMpConstants::NofTrackingChambers();
209 for(Int_t iLoc = 0; iLoc < AliMpConstants::NofLocalBoards(); iLoc++) {
210 Int_t iBoard = iLoc+1;
211 Int_t detElemId = AliMpDDLStore::Instance()->GetDEfromLocalBoard(iBoard, iChamber);
213 if (!detElemId) continue;
215 AliMpLocalBoard* localBoard = AliMpDDLStore::Instance()->GetLocalBoard(iBoard, kFALSE);
218 if( !localBoard->IsNotified())
222 const AliMpVSegmentation* seg[2] = {
223 AliMpSegmentation::Instance()->GetMpSegmentation(detElemId, AliMp::GetCathodType(0)),
224 AliMpSegmentation::Instance()->GetMpSegmentation(detElemId, AliMp::GetCathodType(1))};
227 AliMpPad pad1 = seg[1]->PadByLocation(AliMpIntPair(iBoard,0),kFALSE);
228 yWidthSlat = pad1.Dimensions().Y();
229 AliMpPad pad0 = seg[0]->PadByLocation(AliMpIntPair(iBoard,0),kFALSE);
230 xOffsetLine = TMath::Abs(pad0.Position().X()) + pad0.Dimensions().X();
231 xWidthCol = 2.* pad0.Dimensions().X();
234 // Get ideal global position of DetElemId center
235 slat = detElemId%100;
236 line = (4 + slat)%18;
242 yOffsetLine = (Float_t)(line - 4) * 2. * yWidthSlat;
244 for(Int_t cath=0; cath<AliMpConstants::NofCathodes(); cath++){
246 // necessary because strip info is read for each cathode
247 // board info is read only from cathode 0
250 for (Int_t iStrip = 0; iStrip < 16; ++iStrip) {
251 // get pad from electronics
254 if (cath && localBoard->GetSwitch(6)) offset = -8;
256 AliMpPad pad = seg[cath]->PadByLocation(AliMpIntPair(iBoard,iStrip+offset),kFALSE);
258 if (!pad.IsValid()) continue;
260 xWidth = pad.Dimensions().X();
261 yWidth = pad.Dimensions().Y();
262 xcPad = sign * (pad.Position().X() + xOffsetLine);
263 if(line==4) xcPad += 0.75 * sign * xWidthCol;
264 ycPad = pad.Position().Y() + yOffsetLine;
267 x1 = xcPad - xWidth + kShiftX;
268 y1 = ycPad - yWidth + kShiftY;
269 x2 = xcPad + xWidth - kShiftX;
270 y2 = ycPad + yWidth - kShiftY;
273 AddSortedPoint(x1, xAxisStrip, kResetValue);
274 AddSortedPoint(x2, xAxisStrip, kResetValue);
276 AddSortedPoint(y1, yAxisStrip, kResetValue);
277 AddSortedPoint(y2, yAxisStrip, kResetValue);
279 else if(displayType==kDisplayStrips)
280 FillBins(inputHisto, displayHisto, iBoard, iStrip, x1, x2, y1, y2, kShiftX, kShiftY, displayOpt);
284 if(iStrip+offset==0) {
285 x1b = xcPad - xWidth + kShiftEl;
286 y1b = ycPad - yWidth + kShiftEl;
288 x2b = xcPad + xWidth - kShiftEl;
289 y2b = ycPad + yWidth - kShiftEl;
293 // if iCath==0 strip info and board info are both filled -> break!
294 // if iCath==1 board info is filled at cath==0. Strip info to be filled at cath==1
296 } // loop on cathodes
300 AddSortedPoint(x1b, xAxisBoard, kResetValue);
301 AddSortedPoint(x2b, xAxisBoard, kResetValue);
303 AddSortedPoint(y1b, yAxisBoard, kResetValue);
304 AddSortedPoint(y2b, yAxisBoard, kResetValue);
306 else if(displayType==kDisplayBoards)
307 FillBins(inputHisto, displayHisto, iBoard, -1, x1b, x2b, y1b, y2b, kShiftEl, kShiftEl, displayOpt);
308 else if(displayType==kDisplaySlats)
309 FillBins(inputHisto, displayHisto, slat, -1, x1b, x2b, y1b, y2b, kShiftEl, kShiftEl, displayOpt);
310 } // loop on local boards
312 if(inputHisto) return kTRUE;
314 displayHisto->Reset();
317 const Float_t kMinDiff = 0.1;
319 TArrayD* currArray[4] = {&xAxisStrip, &yAxisStrip,
320 &xAxisBoard, &yAxisBoard};
321 for(Int_t iaxis=0; iaxis<4; iaxis++){
323 while(TMath::Abs((*currArray[iaxis])[ipoint]-kResetValue)>kMinDiff) {
326 if(ipoint>currArray[iaxis]->GetSize()-2)
327 AliWarning(Form("Array size (%i) lower than the number of points!", currArray[iaxis]->GetSize()));
328 currArray[iaxis]->Set(ipoint);
333 displayHisto->SetBins(xAxisStrip.GetSize()-1, xAxisStrip.GetArray(),
334 yAxisStrip.GetSize()-1, yAxisStrip.GetArray());
338 displayHisto->SetBins(xAxisBoard.GetSize()-1, xAxisBoard.GetArray(),
339 yAxisBoard.GetSize()-1, yAxisBoard.GetArray());
343 displayHisto->SetName(displayHistoName.Data());
344 displayHisto->SetTitle(displayHistoTitle.Data());
345 displayHisto->SetXTitle("X (cm)");
346 displayHisto->SetYTitle("Y (cm)");
347 //displayHisto->SetStats(kFALSE);
353 //____________________________________________________________________________
354 Bool_t AliMUONTriggerDisplay::AddSortedPoint(Float_t currVal, TArrayD& position, const Float_t kResetValue)
357 /// Add sorted point in array according to an increasing order.
358 /// Used to build display histograms axis.
360 Int_t nEntries = position.GetSize()-1;
362 const Float_t kMinDiff = 0.1;
363 for(Int_t i=0; i<nEntries; i++){
364 if(TMath::Abs(position[i]-currVal)<kMinDiff) return kFALSE;
365 if(TMath::Abs(position[i]-kResetValue)<kMinDiff) {
366 position[i] = currVal;
369 if(currVal>position[i]) continue;
371 position[i] = currVal;
372 for(Int_t j=i+1; j<nEntries; j++){
376 if(tmp1==kResetValue) break;
384 //____________________________________________________________________________
385 void AliMUONTriggerDisplay::FillBins(TH1* inputHisto, TH2* displayHisto,
386 Int_t iElement1, Int_t iElement2,
387 Float_t x1, Float_t x2, Float_t y1, Float_t y2,
388 const Float_t kShiftX, const Float_t kShiftY,
389 EDisplayOption displayOpt)
392 /// Given the bin in inputHisto, search the corresponding bins
393 /// in display histo and fill it.
395 TString className = inputHisto->ClassName();
397 Float_t binContent=0;
398 Int_t binX = inputHisto->GetXaxis()->FindBin(iElement1);
399 if(className.Contains("2")) {
400 binY = inputHisto->GetYaxis()->FindBin(iElement2);
401 binContent = inputHisto->GetBinContent(binX, binY);
403 else binContent = inputHisto->GetBinContent(binX);
406 if(displayOpt==kShowZeroes) binContent = 1e-5;
410 Int_t binX1 = displayHisto->GetXaxis()->FindBin(x1 + 0.01*kShiftX);
411 Int_t binX2 = displayHisto->GetXaxis()->FindBin(x2 - 0.01*kShiftX);
412 Int_t binY1 = displayHisto->GetYaxis()->FindBin(y1 + 0.01*kShiftY);
413 Int_t binY2 = displayHisto->GetYaxis()->FindBin(y2 - 0.01*kShiftY);
415 if(displayOpt==kNumbered) {
416 Int_t meanBin = (binX1+binX2)/2;
420 meanBin = (binY1+binY2)/2;
425 for(Int_t ibinx=binX1; ibinx<=binX2; ibinx++){
426 for(Int_t ibiny=binY1; ibiny<=binY2; ibiny++){
427 displayHisto->SetBinContent(ibinx,ibiny,binContent);