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 **************************************************************************/
16 /** @file AliFMDParameters.cxx
17 @author Christian Holm Christensen <cholm@nbi.dk>
18 @date Mon Mar 27 12:44:26 2006
19 @brief Manager of FMD parameters
21 //____________________________________________________________________
23 // Forward Multiplicity Detector based on Silicon wafers.
25 // This class is a singleton that handles various parameters of
27 // The manager normally serves the parameters from the Conditions
28 // Database (CDB). These are retrivied by the member function
29 // `Init'. Optionally, the class can serve hard-coded constants, if
30 // no CDB is available.
32 #include "AliLog.h" // ALILOG_H
33 #include "AliFMDParameters.h" // ALIFMDPARAMETERS_H
34 #include "AliFMDGeometry.h" // ALIFMDGEOMETRY_H
35 #include "AliFMDRing.h" // ALIFMDRING_H
36 #include "AliFMDCalibGain.h" // ALIFMDCALIBGAIN_H
37 #include "AliFMDCalibPedestal.h" // ALIFMDCALIBPEDESTAL_H
38 #include "AliFMDCalibSampleRate.h" // ALIFMDCALIBSAMPLERATE_H
39 #include "AliFMDCalibStripRange.h" // ALIFMDCALIBSTRIPRANGE_H
40 #include "AliFMDAltroMapping.h" // ALIFMDALTROMAPPING_H
41 #include <AliCDBManager.h> // ALICDBMANAGER_H
42 #include <AliCDBEntry.h> // ALICDBMANAGER_H
43 #include <Riostream.h>
48 //====================================================================
49 ClassImp(AliFMDParameters)
51 ; // This is here to keep Emacs for indenting the next line
54 //____________________________________________________________________
55 AliFMDParameters* AliFMDParameters::fgInstance = 0;
57 //____________________________________________________________________
58 const char* AliFMDParameters::fgkPulseGain = "FMD/Calib/PulseGain";
59 const char* AliFMDParameters::fgkPedestal = "FMD/Calib/Pedestal";
60 const char* AliFMDParameters::fgkDead = "FMD/Calib/Dead";
61 const char* AliFMDParameters::fgkSampleRate = "FMD/Calib/SampleRate";
62 const char* AliFMDParameters::fgkAltroMap = "FMD/Calib/AltroMap";
63 const char* AliFMDParameters::fgkZeroSuppression = "FMD/Calib/ZeroSuppression";
64 const char* AliFMDParameters::fgkStripRange = "FMD/Calib/StripRange";
67 //____________________________________________________________________
69 AliFMDParameters::Instance()
71 // Get static instance
72 if (!fgInstance) fgInstance = new AliFMDParameters;
76 //____________________________________________________________________
77 AliFMDParameters::AliFMDParameters()
85 fFixedPedestalWidth(0),
86 fFixedZeroSuppression(0),
101 // Default constructor
103 SetAltroChannelSize();
104 SetChannelsPerAltro();
105 SetZeroSuppression();
114 //__________________________________________________________________
116 AliFMDParameters::Init()
118 // Initialize the parameters manager. We need to get stuff from the
125 InitZeroSuppression();
131 //__________________________________________________________________
132 #define DET2IDX(det,ring,sec,str) \
133 (det * 10000 + (ring == 'I' ? 0 : 1000) + str)
135 //__________________________________________________________________
137 AliFMDParameters::Draw(Option_t* option)
141 kPulseGain, // Path to PulseGain calib object
142 kThreshold, // Path to PulseGain calib object
143 kPedestal, // Path to Pedestal calib object
144 kPedestalWidth, // Path to Pedestal calib object
145 kDead, // Path to Dead calib object
146 kSampleRate, // Path to SampleRate calib object
147 kAltroMap, // Path to AltroMap calib object
148 kZeroSuppression, // Path to ZeroSuppression cal object
149 kMinStripRange, // Path to strip range cal object
150 kMaxStripRange // Path to strip range cal object
154 if (opt.Contains("dead", TString::kIgnoreCase))
156 else if (opt.Contains("threshold",TString::kIgnoreCase))
158 else if (opt.Contains("gain",TString::kIgnoreCase))
160 else if (opt.Contains("pedestal",TString::kIgnoreCase))
162 else if (opt.Contains("noise",TString::kIgnoreCase))
163 what = kPedestalWidth;
164 else if (opt.Contains("zero",TString::kIgnoreCase))
165 what = kZeroSuppression;
166 else if (opt.Contains("rate",TString::kIgnoreCase))
168 else if (opt.Contains("min",TString::kIgnoreCase))
169 what = kMinStripRange;
170 else if (opt.Contains("max",TString::kIgnoreCase))
171 what = kMaxStripRange;
172 else if (opt.Contains("map",TString::kIgnoreCase))
175 Warning("Draw", "unknown parameter: %s\n\tShould be one of\n\t"
176 "dead, threshold, gain, pedestal, noise, zero, rate, "
182 TArrayD xbins(3 * 512 + 2 * 256 + 5);
185 for (UShort_t det = 1; det <= 3; det++) {
186 UShort_t nRings = (det == 1 ? 1 : 2);
187 for (UShort_t iring = 0; iring < nRings; iring++) {
188 UShort_t nStrip = (iring == 0 ? 512 : 256);
189 Char_t ring = (iring == 0 ? 'I' : 'O');
190 for (UShort_t str = 0; str < nStrip; str++) {
191 Int_t idx = DET2IDX(det, ring, 0, str);
193 xbins[i-1] = idx - .5;
204 for (Int_t i = 0; i < 41; i++) ybins[i] = Float_t(i - .5);
205 TH2D* hist = new TH2D("calib", Form("Calibration %s", option),
206 xbins.fN-1, xbins.fArray,
207 ybins.fN-1, ybins.fArray);
209 // hist->Draw("Lego");
212 for (UShort_t det = 1; det <= 3; det++) {
213 UShort_t nRings = (det == 1 ? 1 : 2);
214 for (UShort_t iring = 0; iring < nRings; iring++) {
215 UShort_t nSector = (iring == 0 ? 20 : 40);
216 UShort_t nStrip = (iring == 0 ? 512 : 256);
217 Char_t ring = (iring == 0 ? 'I' : 'O');
218 for (UShort_t sec = 0; sec < nSector; sec++) {
219 for (UShort_t str = 0; str < nStrip; str++) {
220 Int_t idx = DET2IDX(det, ring, sec, str);
224 case kPulseGain: // Path to PulseGain calib object
225 val = GetPulseGain(det,ring,sec,str); break;
226 case kThreshold: // Path to PulseGain calib object
227 val = GetThreshold(); break;
228 case kPedestal: // Path to Pedestal calib object
229 val = GetPedestal(det,ring,sec,str); break;
230 case kPedestalWidth: // Path to Pedestal calib object
231 val = GetPedestalWidth(det,ring,sec,str); break;
232 case kDead: // Path to Dead calib object
233 val = IsDead(det,ring,sec,str); break;
234 case kSampleRate: // Path to SampleRate calib object
235 val = GetSampleRate(det,ring,sec,str); break;
236 case kAltroMap: // Path to AltroMap calib object
237 Detector2Hardware(det,ring,sec,str, ddl, addr);
239 case kZeroSuppression: // Path to ZeroSuppression cal object
240 val = GetZeroSuppression(det,ring,sec,str); break;
241 case kMinStripRange: // Path to strip range cal object
242 val = GetMinStrip(det,ring,sec,str); break;
243 case kMaxStripRange: // Path to strip range cal object
244 val = GetMaxStrip(det,ring,sec,str); break;
246 hist->Fill(idx,sec,val);
254 //__________________________________________________________________
256 AliFMDParameters::Print(Option_t* option) const
258 // Print information.
259 // If option contains an 'A' then everything is printed.
261 Bool_t showStrips = opt.Contains("a", TString::kIgnoreCase);
262 if (opt.Contains("fmd",TString::kIgnoreCase)) {
263 size_t i = opt.Index("fmd",TString::kIgnoreCase);
264 size_t j = opt.Index("]",TString::kIgnoreCase);
265 UShort_t det, sec, str;
266 Char_t ring, lbrack, rbrack, comma;
268 std::stringstream s(opt(i+4, j-i-3).Data());
269 s >> det >> ring >> lbrack >> sec >> comma >> str >> rbrack;
270 Detector2Hardware(det, ring, sec, str, ddl, addr);
272 << " Strip | Pedestal | Gain | ZS thr. | Address\n"
273 << "--------------+-------------------+------------+---------+---------"
274 << "\nFMD" << det << ring << "[" << std::setw(2) << sec << ","
275 << std::setw(3) << str << "] | "
276 << std::setw(7) << GetPedestal(det, ring, sec, str)
277 << "+/-" << std::setw(7)
278 << GetPedestalWidth(det, ring, sec, str)
279 << " | " << std::setw(10)
280 << GetPulseGain(det, ring, sec, str)
281 << " | " << std::setw(7)
282 << GetZeroSuppression(det, ring, sec, str)
283 << " | 0x" << std::hex << std::setw(4)
284 << std::setfill('0') << ddl << ",0x" << std::setw(3)
285 << addr << std::dec << std::setfill(' ') << std::endl;
288 for (UShort_t det=1 ; det <= 3; det++) {
289 std::cout << "FMD" << det << std::endl;
290 Char_t rings[] = { 'I', (det == 1 ? '\0' : 'O'), '\0' };
291 for (Char_t* ring = rings; *ring != '\0'; ring++) {
292 std::cout << " Ring " << *ring << std::endl;
293 UShort_t nSec = ( *ring == 'I' ? 20 : 40 );
294 UShort_t nStr = ( *ring == 'I' ? 512 : 256 );
295 for (UShort_t sec = 0; sec < nSec; sec++) {
296 UShort_t min = GetMinStrip(det, *ring, sec, 0);
297 UShort_t max = GetMaxStrip(det, *ring, sec, 0);
298 UShort_t rate = GetSampleRate(det, *ring, sec, 0);
299 std::cout << " Sector " << std::setw(2) << sec
300 << " Strip range: " << std::setw(3) << min << ","
301 << std::setw(3) << max << " Rate: " << std::setw(2)
302 << rate << std::endl;
303 if (!showStrips) continue;
305 << " Strip | Pedestal | Gain | ZS thr. | Address\n"
306 << "--------+-------------------+------------+---------+---------"
308 for (UShort_t str = 0; str < nStr; str++) {
309 std::cout << " " << std::setw(3) << str << " | ";
310 if (IsDead(det, *ring, sec, str)) {
311 std::cout << "dead" << std::endl;
315 Detector2Hardware(det, *ring, sec, str, ddl, addr);
316 std::cout << std::setw(7) << GetPedestal(det, *ring, sec, str)
317 << "+/-" << std::setw(7)
318 << GetPedestalWidth(det, *ring, sec, str)
319 << " | " << std::setw(10)
320 << GetPulseGain(det, *ring, sec, str)
321 << " | " << std::setw(5)
322 << GetZeroSuppression(det, *ring, sec, str)
323 << " | 0x" << std::hex << std::setw(4)
324 << std::setfill('0') << ddl << ",0x" << std::setw(3)
325 << addr << std::dec << std::setfill(' ') << std::endl;
332 //__________________________________________________________________
334 AliFMDParameters::SetStripRange(UShort_t min, UShort_t max)
336 // Set fixed strip range
337 fFixedMinStrip = min;
338 fFixedMaxStrip = max;
341 //__________________________________________________________________
343 AliFMDParameters::InitPulseGain()
345 // Get pulse gain from CDB or used fixed
346 AliCDBManager* cdb = AliCDBManager::Instance();
347 AliCDBEntry* gain = cdb->Get(fgkPulseGain);
349 AliWarning(Form("No %s found in CDB, perhaps you need to "
350 "use AliFMDCalibFaker?", fgkPulseGain));
354 AliDebug(1, Form("Got gain from CDB"));
355 fPulseGain = dynamic_cast<AliFMDCalibGain*>(gain->GetObject());
356 if (!fPulseGain) AliWarning("Invalid pulser gain object from CDB");
358 //__________________________________________________________________
360 AliFMDParameters::InitPedestal()
362 // Initialize the pedestals from CDB
363 AliCDBManager* cdb = AliCDBManager::Instance();
364 AliCDBEntry* pedestal = cdb->Get(fgkPedestal);
366 AliWarning(Form("No %s found in CDB, perhaps you need to "
367 "use AliFMDCalibFaker?", fgkPedestal));
370 AliDebug(1, Form("Got pedestal from CDB"));
371 fPedestal = dynamic_cast<AliFMDCalibPedestal*>(pedestal->GetObject());
372 if (!fPedestal) AliWarning("Invalid pedestal object from CDB");
375 //__________________________________________________________________
377 AliFMDParameters::InitDeadMap()
379 // Get Dead-channel-map from CDB
380 AliCDBManager* cdb = AliCDBManager::Instance();
381 AliCDBEntry* deadMap = cdb->Get(fgkDead);
383 AliWarning(Form("No %s found in CDB, perhaps you need to "
384 "use AliFMDCalibFaker?", fgkDead));
387 AliDebug(1, Form("Got dead map from CDB"));
388 fDeadMap = dynamic_cast<AliFMDCalibDeadMap*>(deadMap->GetObject());
389 if (!fDeadMap) AliWarning("Invalid dead map object from CDB");
392 //__________________________________________________________________
394 AliFMDParameters::InitZeroSuppression()
396 // Get 0-suppression from CDB
397 AliCDBManager* cdb = AliCDBManager::Instance();
398 AliCDBEntry* zeroSup = cdb->Get(fgkZeroSuppression);
400 AliWarning(Form("No %s found in CDB, perhaps you need to "
401 "use AliFMDCalibFaker?", fgkZeroSuppression));
404 AliDebug(1, Form("Got zero suppression from CDB"));
406 dynamic_cast<AliFMDCalibZeroSuppression*>(zeroSup->GetObject());
407 if (!fZeroSuppression)AliWarning("Invalid zero suppression object from CDB");
410 //__________________________________________________________________
412 AliFMDParameters::InitSampleRate()
414 // get Sample rate from CDB
415 AliCDBManager* cdb = AliCDBManager::Instance();
416 AliCDBEntry* sampRat = cdb->Get(fgkSampleRate);
418 AliWarning(Form("No %s found in CDB, perhaps you need to "
419 "use AliFMDCalibFaker?", fgkSampleRate));
422 AliDebug(1, Form("Got zero suppression from CDB"));
423 fSampleRate = dynamic_cast<AliFMDCalibSampleRate*>(sampRat->GetObject());
424 if (!fSampleRate) AliWarning("Invalid zero suppression object from CDB");
427 //__________________________________________________________________
429 AliFMDParameters::InitAltroMap()
431 // Get hardware mapping from CDB
432 AliCDBManager* cdb = AliCDBManager::Instance();
433 AliCDBEntry* hwMap = cdb->Get(fgkAltroMap);
435 AliWarning(Form("No %s found in CDB, perhaps you need to "
436 "use AliFMDCalibFaker?", fgkAltroMap));
437 fAltroMap = new AliFMDAltroMapping;
440 AliDebug(1, Form("Got ALTRO map from CDB"));
441 fAltroMap = dynamic_cast<AliFMDAltroMapping*>(hwMap->GetObject());
443 AliWarning("Invalid ALTRO map object from CDB");
444 fAltroMap = new AliFMDAltroMapping;
448 //__________________________________________________________________
450 AliFMDParameters::InitStripRange()
452 // Get strips read-out from CDB
453 AliCDBManager* cdb = AliCDBManager::Instance();
454 AliCDBEntry* range = cdb->Get(fgkStripRange);
456 AliWarning(Form("No %s found in CDB, perhaps you need to "
457 "use AliFMDCalibFaker?", fgkStripRange));
460 AliDebug(1, Form("Got strip range from CDB"));
461 fStripRange = dynamic_cast<AliFMDCalibStripRange*>(range->GetObject());
462 if (!fStripRange) AliWarning("Invalid strip range object from CDB");
466 //__________________________________________________________________
468 AliFMDParameters::GetThreshold() const
470 // Get threshold from CDB
471 if (!fPulseGain) return fFixedThreshold;
472 return fPulseGain->Threshold();
475 //__________________________________________________________________
477 AliFMDParameters::GetPulseGain(UShort_t detector, Char_t ring,
478 UShort_t sector, UShort_t strip) const
480 // Returns the pulser calibrated gain for strip # strip in sector #
481 // sector or ring id ring of detector # detector.
483 // For simulation, this is normally set to
486 // ------------------ * MIP_Energy_Loss
487 // ALTRO_channel_size
490 if (fFixedPulseGain <= 0)
491 fFixedPulseGain = fVA1MipRange * GetEdepMip() / fAltroChannelSize;
492 return fFixedPulseGain;
494 AliDebug(50, Form("pulse gain for FMD%d%c[%2d,%3d]=%f",
495 detector, ring, sector, strip,
496 fPulseGain->Value(detector, ring, sector, strip)));
497 return fPulseGain->Value(detector, ring, sector, strip);
500 //__________________________________________________________________
502 AliFMDParameters::IsDead(UShort_t detector, Char_t ring,
503 UShort_t sector, UShort_t strip) const
505 // Check if the channel is dead
506 if (!fDeadMap) return kFALSE;
507 AliDebug(50, Form("Dead for FMD%d%c[%2d,%3d]=%s",
508 detector, ring, sector, strip,
509 fDeadMap->operator()(detector, ring, sector, strip) ?
511 return fDeadMap->operator()(detector, ring, sector, strip);
514 //__________________________________________________________________
516 AliFMDParameters::GetZeroSuppression(UShort_t detector, Char_t ring,
517 UShort_t sector, UShort_t strip) const
519 // Get zero suppression threshold
520 if (!fZeroSuppression) return fFixedZeroSuppression;
521 // Need to map strip to ALTRO chip.
522 AliDebug(50, Form("zero sup. for FMD%d%c[%2d,%3d]=%f",
523 detector, ring, sector, strip,
524 fZeroSuppression->operator()(detector, ring,
526 return fZeroSuppression->operator()(detector, ring, sector, strip/128);
529 //__________________________________________________________________
531 AliFMDParameters::GetSampleRate(UShort_t det, Char_t ring, UShort_t sector,
535 if (!fSampleRate) return fFixedSampleRate;
536 // Need to map sector to digitizier card.
537 UInt_t ret = fSampleRate->Rate(det, ring, sector, str);
538 AliDebug(50, Form("Sample rate for FMD%d%c[%2d,%3d]=%d",
539 det, ring, sector, str, ret));
543 //__________________________________________________________________
545 AliFMDParameters::GetMinStrip(UShort_t det, Char_t ring, UShort_t sector,
548 // Get strip range read out
549 if (!fStripRange) return fFixedMinStrip;
550 // Need to map sector to digitizier card.
551 UInt_t ret = fStripRange->Min(det, ring, sector, str);
552 AliDebug(50, Form("Min strip # for FMD%d%c[%2d,%3d]=%d",
553 det, ring, sector, str, ret));
557 //__________________________________________________________________
559 AliFMDParameters::GetMaxStrip(UShort_t det, Char_t ring, UShort_t sector,
562 // Get strip range read out
563 if (!fStripRange) return fFixedMaxStrip;
564 // Need to map sector to digitizier card.
565 UInt_t ret = fStripRange->Max(det, ring, sector, str);
566 AliDebug(50, Form("Max strip # for FMD%d%c[%2d,%3d]=%d",
567 det, ring, sector, str, ret));
571 //__________________________________________________________________
573 AliFMDParameters::GetPedestal(UShort_t detector, Char_t ring,
574 UShort_t sector, UShort_t strip) const
577 if (!fPedestal) return fFixedPedestal;
578 AliDebug(50, Form("pedestal for FMD%d%c[%2d,%3d]=%f",
579 detector, ring, sector, strip,
580 fPedestal->Value(detector, ring, sector, strip)));
581 return fPedestal->Value(detector, ring, sector, strip);
584 //__________________________________________________________________
586 AliFMDParameters::GetPedestalWidth(UShort_t detector, Char_t ring,
587 UShort_t sector, UShort_t strip) const
590 if (!fPedestal) return fFixedPedestalWidth;
591 AliDebug(50, Form("pedetal width for FMD%d%c[%2d,%3d]=%f",
592 detector, ring, sector, strip,
593 fPedestal->Width(detector, ring, sector, strip)));
594 return fPedestal->Width(detector, ring, sector, strip);
597 //__________________________________________________________________
599 AliFMDParameters::GetAltroMap() const
601 // Get the hardware address to detector index map
606 //__________________________________________________________________
608 AliFMDParameters::Hardware2Detector(UInt_t ddl, UInt_t addr, UShort_t& det,
609 Char_t& ring, UShort_t& sec,
612 // Map hardware address to detector index
613 if (!fAltroMap) return kFALSE;
614 return fAltroMap->Hardware2Detector(ddl, addr, det, ring, sec, str);
617 //__________________________________________________________________
619 AliFMDParameters::Detector2Hardware(UShort_t det, Char_t ring, UShort_t sec,
620 UShort_t str, UInt_t& ddl,
623 // Map detector index to hardware address
624 if (!fAltroMap) return kFALSE;
625 return fAltroMap->Detector2Hardware(det, ring, sec, str, ddl, addr);
629 //__________________________________________________________________
631 AliFMDParameters::GetEdepMip() const
633 // Get energy deposited by a MIP in the silicon sensors
635 AliFMDGeometry* fmd = AliFMDGeometry::Instance();
636 fEdepMip = (fkSiDeDxMip
637 * fmd->GetRing('I')->GetSiThickness()
638 * fmd->GetSiDensity());
647 //____________________________________________________________________