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 "AliFMDDebug.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(Bool_t forceReInit, UInt_t what)
118 // Initialize the parameters manager. We need to get stuff from the
120 if (forceReInit) fIsInit = kFALSE;
122 if (what & kPulseGain) InitPulseGain();
123 if (what & kPedestal) InitPedestal();
124 if (what & kDeadMap) InitDeadMap();
125 if (what & kSampleRate) InitSampleRate();
126 if (what & kZeroSuppression) InitZeroSuppression();
127 if (what & kAltroMap) InitAltroMap();
132 //__________________________________________________________________
133 #define DET2IDX(det,ring,sec,str) \
134 (det * 1000 + (ring == 'I' ? 0 : 512) + str)
136 //__________________________________________________________________
138 AliFMDParameters::Draw(Option_t* option)
142 kPulseGain, // Path to PulseGain calib object
143 kThreshold, // Path to PulseGain calib object
144 kPedestal, // Path to Pedestal calib object
145 kPedestalWidth, // Path to Pedestal calib object
146 kDead, // Path to Dead calib object
147 kSampleRate, // Path to SampleRate calib object
148 kAltroMap, // Path to AltroMap calib object
149 kZeroSuppression, // Path to ZeroSuppression cal object
150 kMinStripRange, // Path to strip range cal object
151 kMaxStripRange // Path to strip range cal object
155 if (opt.Contains("dead", TString::kIgnoreCase))
157 else if (opt.Contains("threshold",TString::kIgnoreCase))
159 else if (opt.Contains("gain",TString::kIgnoreCase))
161 else if (opt.Contains("pedestal",TString::kIgnoreCase))
163 else if (opt.Contains("noise",TString::kIgnoreCase))
164 what = kPedestalWidth;
165 else if (opt.Contains("zero",TString::kIgnoreCase))
166 what = kZeroSuppression;
167 else if (opt.Contains("rate",TString::kIgnoreCase))
169 else if (opt.Contains("min",TString::kIgnoreCase))
170 what = kMinStripRange;
171 else if (opt.Contains("max",TString::kIgnoreCase))
172 what = kMaxStripRange;
173 else if (opt.Contains("map",TString::kIgnoreCase))
176 Warning("Draw", "unknown parameter: %s\n\tShould be one of\n\t"
177 "dead, threshold, gain, pedestal, noise, zero, rate, "
183 TArrayD xbins(3 * 512 + 2 * 256 + 5);
186 for (UShort_t det = 1; det <= 3; det++) {
187 UShort_t nRings = (det == 1 ? 1 : 2);
188 for (UShort_t iring = 0; iring < nRings; iring++) {
189 UShort_t nStrip = (iring == 0 ? 512 : 256);
190 Char_t ring = (iring == 0 ? 'I' : 'O');
191 for (UShort_t str = 0; str < nStrip; str++) {
192 // UShort_t nSec = (iring == 0 ? 20 : 40);
193 // Char_t ring = (iring == 0 ? 'I' : 'O');
194 // for (UShort_t sec = 0; sec < nSec; sec++) {
195 Int_t idx = DET2IDX(det, ring, 0, str);
196 // Int_t idx = DET2IDX(det, ring, sec, 0);
198 xbins[i-1] = idx - .5;
209 for (Int_t i = 0; i < ybins.fN; i++) ybins[i] = Float_t(i - .5);
210 TH2D* hist = new TH2D("calib", Form("Calibration %s", option),
211 xbins.fN-1, xbins.fArray,
212 ybins.fN-1, ybins.fArray);
213 hist->GetXaxis()->SetTitle("1000 #times detector + 512 #times ring + strip");
214 hist->GetYaxis()->SetTitle("sector");
216 // hist->Draw("Lego");
219 for (UShort_t det = 1; det <= 3; det++) {
220 UShort_t nRings = (det == 1 ? 1 : 2);
221 for (UShort_t iring = 0; iring < nRings; iring++) {
222 UShort_t nSector = (iring == 0 ? 20 : 40);
223 UShort_t nStrip = (iring == 0 ? 512 : 256);
224 Char_t ring = (iring == 0 ? 'I' : 'O');
225 for (UShort_t sec = 0; sec < nSector; sec++) {
226 for (UShort_t str = 0; str < nStrip; str++) {
227 Int_t idx = DET2IDX(det, ring, sec, str);
231 case kPulseGain: // Path to PulseGain calib object
232 val = GetPulseGain(det,ring,sec,str); break;
233 case kThreshold: // Path to PulseGain calib object
234 val = GetThreshold(); break;
235 case kPedestal: // Path to Pedestal calib object
236 val = GetPedestal(det,ring,sec,str); break;
237 case kPedestalWidth: // Path to Pedestal calib object
238 val = GetPedestalWidth(det,ring,sec,str); break;
239 case kDead: // Path to Dead calib object
240 val = IsDead(det,ring,sec,str); break;
241 case kSampleRate: // Path to SampleRate calib object
242 val = GetSampleRate(det,ring,sec,str); break;
243 case kAltroMap: // Path to AltroMap calib object
244 Detector2Hardware(det,ring,sec,str, ddl, addr);
246 case kZeroSuppression: // Path to ZeroSuppression cal object
247 val = GetZeroSuppression(det,ring,sec,str); break;
248 case kMinStripRange: // Path to strip range cal object
249 val = GetMinStrip(det,ring,sec,str); break;
250 case kMaxStripRange: // Path to strip range cal object
251 val = GetMaxStrip(det,ring,sec,str); break;
253 hist->Fill(idx,sec,val);
254 // hist->Fill(idx,str,val);
262 //__________________________________________________________________
264 AliFMDParameters::Print(Option_t* option) const
266 // Print information.
267 // If option contains an 'A' then everything is printed.
268 // If the option contains the string "FMD" the function will search
269 // for detector, ring, sector, and strip numbers to print, in the
272 // FMD<detector><ring>[<sector>,<string>]
274 // The wild card '*' means all of <detector>, <ring>, <sector>, or
277 Bool_t showStrips = opt.Contains("a", TString::kIgnoreCase);
278 UShort_t ds[] = { 1, 2, 3, 0 };
279 Char_t rs[] = { 'I', 'O', '\0' };
280 UShort_t minStrip = 0;
281 UShort_t maxStrip = 512;
282 UShort_t minSector = 0;
283 UShort_t maxSector = 40;
286 if (opt.Contains("fmd",TString::kIgnoreCase)) {
288 size_t i = opt.Index("fmd",TString::kIgnoreCase);
289 size_t j = opt.Index("]",TString::kIgnoreCase);
300 std::stringstream s(opt(i+4, j-i-3).Data());
301 while (state != kEnd) {
302 Char_t tmp = s.peek();
303 if (tmp == ' ' || tmp == '\t') {
308 case kReadDet: { // First, try to kRead the detector
309 if (tmp == '*') s.get();
318 state = (s.bad() ? kEnd : kReadRing);
320 case kReadRing: { // Then try to read the ring;
323 if (ring != '*' && !s.bad()) {
327 state = (s.bad() ? kEnd : kReadLbrack);
329 case kReadLbrack: { // Try to read a left bracket
332 state = (s.bad() ? kEnd : kReadSector);
334 case kReadSector: { // Try to read a sector
335 if (tmp == '*') s.get();
344 state = (s.bad() ? kEnd : kReadComma);
346 case kReadComma: { // Try to read a left bracket
349 state = (s.bad() ? kEnd : kReadStrip);
351 case kReadStrip: { // Try to read a strip
352 if (tmp == '*') s.get();
361 state = (s.bad() ? kEnd : kReadRbrack);
363 case kReadRbrack: { // Try to read a left bracket
375 while ((det = *(dp++))) {
379 while ((ring = *(rp++))) {
380 if (det == 1 && ring == 'O') continue;
381 UShort_t min = GetMinStrip(det, ring, 0, 0);
382 UShort_t max = GetMaxStrip(det, ring, 0, 0);
383 UShort_t rate = GetSampleRate(det, ring, 0, 0);
384 std::cout << "FMD" << det << ring
386 << std::setw(3) << min << ","
387 << std::setw(3) << max << " Rate: "
388 << std::setw(2) << rate << std::endl;
390 if (!showStrips) continue;
391 UShort_t nSec = ( ring == 'I' ? 20 : 40 );
392 UShort_t nStr = ( ring == 'I' ? 512 : 256 );
393 for (UShort_t sec = minSector; sec < maxSector && sec < nSec; sec++) {
395 << " Strip | Pedestal | Gain | ZS thr. | Address\n"
396 << "--------+-------------------+------------+---------+---------"
398 for (UShort_t str = minStrip; str < nStr && str < maxStrip; str++) {
399 if (str == minStrip) std::cout << std::setw(3) << sec << ",";
400 else std::cout << " ";
401 std::cout << std::setw(3) << str << " | ";
402 if (IsDead(det, ring, sec, str)) {
403 std::cout << "dead" << std::endl;
407 Detector2Hardware(det, ring, sec, str, ddl, addr);
408 std::cout << std::setw(7) << GetPedestal(det, ring, sec, str)
409 << "+/-" << std::setw(7)
410 << GetPedestalWidth(det, ring, sec, str)
411 << " | " << std::setw(10)
412 << GetPulseGain(det, ring, sec, str)
413 << " | " << std::setw(7)
414 << GetZeroSuppression(det, ring, sec, str)
415 << " | 0x" << std::hex << std::setw(4)
416 << std::setfill('0') << ddl << ",0x" << std::setw(3)
417 << addr << std::dec << std::setfill(' ') << std::endl;
421 << "============================================================="
428 //__________________________________________________________________
430 AliFMDParameters::SetStripRange(UShort_t min, UShort_t max)
432 // Set fixed strip range
433 fFixedMinStrip = min;
434 fFixedMaxStrip = max;
437 //__________________________________________________________________
439 AliFMDParameters::InitPulseGain()
441 // Get pulse gain from CDB or used fixed
442 AliCDBManager* cdb = AliCDBManager::Instance();
443 AliCDBEntry* gain = cdb->Get(fgkPulseGain);
445 AliFatal(Form("No %s found in CDB, perhaps you need to "
446 "use AliFMDCalibFaker?", fgkPulseGain));
450 AliFMDDebug(1, ("Got gain from CDB"));
451 fPulseGain = dynamic_cast<AliFMDCalibGain*>(gain->GetObject());
452 if (!fPulseGain) AliFatal("Invalid pulser gain object from CDB");
454 //__________________________________________________________________
456 AliFMDParameters::InitPedestal()
458 // Initialize the pedestals from CDB
459 AliCDBManager* cdb = AliCDBManager::Instance();
460 AliCDBEntry* pedestal = cdb->Get(fgkPedestal);
462 AliFatal(Form("No %s found in CDB, perhaps you need to "
463 "use AliFMDCalibFaker?", fgkPedestal));
466 AliFMDDebug(1, ("Got pedestal from CDB"));
467 fPedestal = dynamic_cast<AliFMDCalibPedestal*>(pedestal->GetObject());
468 if (!fPedestal) AliFatal("Invalid pedestal object from CDB");
471 //__________________________________________________________________
473 AliFMDParameters::InitDeadMap()
475 // Get Dead-channel-map from CDB
476 AliCDBManager* cdb = AliCDBManager::Instance();
477 AliCDBEntry* deadMap = cdb->Get(fgkDead);
479 AliFatal(Form("No %s found in CDB, perhaps you need to "
480 "use AliFMDCalibFaker?", fgkDead));
483 AliFMDDebug(1, ("Got dead map from CDB"));
484 fDeadMap = dynamic_cast<AliFMDCalibDeadMap*>(deadMap->GetObject());
485 if (!fDeadMap) AliFatal("Invalid dead map object from CDB");
488 //__________________________________________________________________
490 AliFMDParameters::InitZeroSuppression()
492 // Get 0-suppression from CDB
493 AliCDBManager* cdb = AliCDBManager::Instance();
494 AliCDBEntry* zeroSup = cdb->Get(fgkZeroSuppression);
496 AliFatal(Form("No %s found in CDB, perhaps you need to "
497 "use AliFMDCalibFaker?", fgkZeroSuppression));
500 AliFMDDebug(1, ("Got zero suppression from CDB"));
502 dynamic_cast<AliFMDCalibZeroSuppression*>(zeroSup->GetObject());
503 if (!fZeroSuppression)AliFatal("Invalid zero suppression object from CDB");
506 //__________________________________________________________________
508 AliFMDParameters::InitSampleRate()
510 // get Sample rate from CDB
511 AliCDBManager* cdb = AliCDBManager::Instance();
512 AliCDBEntry* sampRat = cdb->Get(fgkSampleRate);
514 AliFatal(Form("No %s found in CDB, perhaps you need to "
515 "use AliFMDCalibFaker?", fgkSampleRate));
518 AliFMDDebug(1, ("Got zero suppression from CDB"));
519 fSampleRate = dynamic_cast<AliFMDCalibSampleRate*>(sampRat->GetObject());
520 if (!fSampleRate) AliFatal("Invalid zero suppression object from CDB");
523 //__________________________________________________________________
525 AliFMDParameters::InitAltroMap()
527 // Get hardware mapping from CDB
528 AliCDBManager* cdb = AliCDBManager::Instance();
529 AliCDBEntry* hwMap = cdb->Get(fgkAltroMap);
531 AliFatal(Form("No %s found in CDB, perhaps you need to "
532 "use AliFMDCalibFaker?", fgkAltroMap));
533 fAltroMap = new AliFMDAltroMapping;
536 AliFMDDebug(1, ("Got ALTRO map from CDB"));
537 fAltroMap = dynamic_cast<AliFMDAltroMapping*>(hwMap->GetObject());
539 AliFatal("Invalid ALTRO map object from CDB");
540 fAltroMap = new AliFMDAltroMapping;
544 //__________________________________________________________________
546 AliFMDParameters::InitStripRange()
548 // Get strips read-out from CDB
549 AliCDBManager* cdb = AliCDBManager::Instance();
550 AliCDBEntry* range = cdb->Get(fgkStripRange);
552 AliFatal(Form("No %s found in CDB, perhaps you need to "
553 "use AliFMDCalibFaker?", fgkStripRange));
556 AliFMDDebug(1, ("Got strip range from CDB"));
557 fStripRange = dynamic_cast<AliFMDCalibStripRange*>(range->GetObject());
558 if (!fStripRange) AliFatal("Invalid strip range object from CDB");
562 //__________________________________________________________________
564 AliFMDParameters::GetThreshold() const
566 // Get threshold from CDB
567 if (!fPulseGain) return fFixedThreshold;
568 return fPulseGain->Threshold();
571 //__________________________________________________________________
573 AliFMDParameters::GetPulseGain(UShort_t detector, Char_t ring,
574 UShort_t sector, UShort_t strip) const
576 // Returns the pulser calibrated gain for strip # strip in sector #
577 // sector or ring id ring of detector # detector.
579 // For simulation, this is normally set to
582 // ------------------ * MIP_Energy_Loss
583 // ALTRO_channel_size
586 if (fFixedPulseGain <= 0)
587 fFixedPulseGain = fVA1MipRange * GetEdepMip() / fAltroChannelSize;
588 return fFixedPulseGain;
590 AliFMDDebug(50, ("pulse gain for FMD%d%c[%2d,%3d]=%f",
591 detector, ring, sector, strip,
592 fPulseGain->Value(detector, ring, sector, strip)));
593 return fPulseGain->Value(detector, ring, sector, strip);
596 //__________________________________________________________________
598 AliFMDParameters::IsDead(UShort_t detector, Char_t ring,
599 UShort_t sector, UShort_t strip) const
601 // Check if the channel is dead
602 if (!fDeadMap) return kFALSE;
603 AliFMDDebug(50, ("Dead for FMD%d%c[%2d,%3d]=%s",
604 detector, ring, sector, strip,
605 fDeadMap->operator()(detector, ring, sector, strip) ?
607 return fDeadMap->operator()(detector, ring, sector, strip);
610 //__________________________________________________________________
612 AliFMDParameters::GetZeroSuppression(UShort_t detector, Char_t ring,
613 UShort_t sector, UShort_t strip) const
615 // Get zero suppression threshold
616 if (!fZeroSuppression) return fFixedZeroSuppression;
617 // Need to map strip to ALTRO chip.
618 AliFMDDebug(50, ("zero sup. for FMD%d%c[%2d,%3d]=%f",
619 detector, ring, sector, strip,
620 fZeroSuppression->operator()(detector, ring,
622 return fZeroSuppression->operator()(detector, ring, sector, strip/128);
625 //__________________________________________________________________
627 AliFMDParameters::GetSampleRate(UShort_t det, Char_t ring, UShort_t sector,
631 if (!fSampleRate) return fFixedSampleRate;
632 // Need to map sector to digitizier card.
633 UInt_t ret = fSampleRate->Rate(det, ring, sector, str);
634 AliFMDDebug(50, ("Sample rate for FMD%d%c[%2d,%3d]=%d",
635 det, ring, sector, str, ret));
639 //__________________________________________________________________
641 AliFMDParameters::GetMinStrip(UShort_t det, Char_t ring, UShort_t sector,
644 // Get strip range read out
645 if (!fStripRange) return fFixedMinStrip;
646 // Need to map sector to digitizier card.
647 UInt_t ret = fStripRange->Min(det, ring, sector, str);
648 AliFMDDebug(50, ("Min strip # for FMD%d%c[%2d,%3d]=%d",
649 det, ring, sector, str, ret));
653 //__________________________________________________________________
655 AliFMDParameters::GetMaxStrip(UShort_t det, Char_t ring, UShort_t sector,
658 // Get strip range read out
659 if (!fStripRange) return fFixedMaxStrip;
660 // Need to map sector to digitizier card.
661 UInt_t ret = fStripRange->Max(det, ring, sector, str);
662 AliFMDDebug(50, ("Max strip # for FMD%d%c[%2d,%3d]=%d",
663 det, ring, sector, str, ret));
667 //__________________________________________________________________
669 AliFMDParameters::GetPedestal(UShort_t detector, Char_t ring,
670 UShort_t sector, UShort_t strip) const
673 if (!fPedestal) return fFixedPedestal;
674 AliFMDDebug(50, ("pedestal for FMD%d%c[%2d,%3d]=%f",
675 detector, ring, sector, strip,
676 fPedestal->Value(detector, ring, sector, strip)));
677 return fPedestal->Value(detector, ring, sector, strip);
680 //__________________________________________________________________
682 AliFMDParameters::GetPedestalWidth(UShort_t detector, Char_t ring,
683 UShort_t sector, UShort_t strip) const
686 if (!fPedestal) return fFixedPedestalWidth;
687 AliFMDDebug(50, ("pedetal width for FMD%d%c[%2d,%3d]=%f",
688 detector, ring, sector, strip,
689 fPedestal->Width(detector, ring, sector, strip)));
690 return fPedestal->Width(detector, ring, sector, strip);
693 //__________________________________________________________________
695 AliFMDParameters::GetAltroMap() const
697 // Get the hardware address to detector index map
702 //__________________________________________________________________
704 AliFMDParameters::Hardware2Detector(UInt_t ddl, UInt_t board,
705 UInt_t chip, UInt_t chan,
706 UShort_t& det, Char_t& ring,
707 UShort_t& sec, UShort_t& str) const
709 // Map hardware address to detector index
710 if (!fAltroMap) return kFALSE;
711 return fAltroMap->Hardware2Detector(ddl,board,chip,chan, det,ring,sec,str);
713 //__________________________________________________________________
715 AliFMDParameters::Hardware2Detector(UInt_t ddl, UInt_t addr,
716 UShort_t& det, Char_t& ring,
717 UShort_t& sec, UShort_t& str) const
719 // Map hardware address to detector index
720 if (!fAltroMap) return kFALSE;
721 return fAltroMap->Hardware2Detector(ddl, addr, det, ring, sec, str);
724 //__________________________________________________________________
726 AliFMDParameters::Detector2Hardware(UShort_t det, Char_t ring,
727 UShort_t sec, UShort_t str,
728 UInt_t& ddl, UInt_t& board,
729 UInt_t& chip, UInt_t& chan) const
731 // Map detector index to hardware address
732 if (!fAltroMap) return kFALSE;
733 return fAltroMap->Detector2Hardware(det,ring,sec,str, ddl,board,chip,chan);
736 //__________________________________________________________________
738 AliFMDParameters::Detector2Hardware(UShort_t det, Char_t ring,
739 UShort_t sec, UShort_t str,
740 UInt_t& ddl, UInt_t& addr) const
742 // Map detector index to hardware address
743 if (!fAltroMap) return kFALSE;
744 return fAltroMap->Detector2Hardware(det, ring, sec, str, ddl, addr);
748 //__________________________________________________________________
750 AliFMDParameters::GetEdepMip() const
752 // Get energy deposited by a MIP in the silicon sensors
754 AliFMDGeometry* fmd = AliFMDGeometry::Instance();
755 fEdepMip = (fkSiDeDxMip
756 * fmd->GetRing('I')->GetSiThickness()
757 * fmd->GetSiDensity());
766 //____________________________________________________________________