Belle II Software  release-05-01-25
KLMElectronicsMapImporter.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2016 Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Kirill Chilikin, Giacomo De Pietro *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 /* Own header. */
12 #include <klm/calibration/KLMElectronicsMapImporter.h>
13 
14 /* KLM headers. */
15 #include <klm/dataobjects/KLMChannelIndex.h>
16 #include <klm/dbobjects/KLMElectronicsMap.h>
17 
18 /* Belle 2 headers. */
19 #include <framework/database/DBImportObjPtr.h>
20 #include <framework/database/IntervalOfValidity.h>
21 #include <framework/logging/Logger.h>
22 #include <rawdata/dataobjects/RawCOPPERFormat.h>
23 
24 using namespace Belle2;
25 
27  m_ElementNumbers(&(KLMElementNumbers::Instance()))
28 {
29 }
30 
32 {
33 }
34 
35 void KLMElectronicsMapImporter::setIOV(int experimentLow, int runLow,
36  int experimentHigh, int runHigh)
37 {
38  m_ExperimentLow = experimentLow;
39  m_RunLow = runLow;
40  m_ExperimentHigh = experimentHigh;
41  m_RunHigh = runHigh;
42 }
43 
45 {
46  m_ChannelMap.clear();
47 }
48 
50 {
51  const int minimalVersion = 1;
52  const int maximalVersion = 2;
53  if (version < minimalVersion || version > maximalVersion) {
54  B2FATAL("Incorrect version (" << version << ") of BKLM electronics map. "
55  "It must be from " << minimalVersion << " to " << maximalVersion);
56  }
57  int copperId = 0;
58  int slotId = 0;
59  int laneId;
60  int axisId = 0;
62  for (KLMChannelIndex bklmPlane = bklmPlanes.beginBKLM();
63  bklmPlane != bklmPlanes.endBKLM(); ++bklmPlane) {
64  int section = bklmPlane.getSection();
65  int sector = bklmPlane.getSector();
66  int layer = bklmPlane.getLayer();
67  int plane = bklmPlane.getPlane();
68 
69  if (section == BKLMElementNumbers::c_ForwardSection) {
70  if (sector == 3 || sector == 4 || sector == 5 || sector == 6)
71  copperId = 1 + BKLM_ID;
72  if (sector == 1 || sector == 2 || sector == 7 || sector == 8)
73  copperId = 2 + BKLM_ID;
74  }
76  if (sector == 3 || sector == 4 || sector == 5 || sector == 6)
77  copperId = 3 + BKLM_ID;
78  if (sector == 1 || sector == 2 || sector == 7 || sector == 8)
79  copperId = 4 + BKLM_ID;
80  }
81 
82  if (sector == 3 || sector == 4 || sector == 5 || sector == 6)
83  slotId = sector - 2;
84  if (sector == 1 || sector == 2)
85  slotId = sector + 2;
86  if (sector == 7 || sector == 8)
87  slotId = sector - 6;
88 
90  laneId = layer + 5;
91  axisId = plane;
92  } else {
93  laneId = layer;
94  if (plane == BKLMElementNumbers::c_ZPlane)
95  axisId = 1;
96  if (plane == BKLMElementNumbers::c_PhiPlane)
97  axisId = 0;
98  }
99 
100  int MaxiChannel = BKLMElementNumbers::getNStrips(
101  section, sector, layer, plane);
102 
103  bool dontFlip = false;
104  if (section == BKLMElementNumbers::c_ForwardSection &&
105  (sector == 7 || sector == 8 || sector == 1 || sector == 2))
106  dontFlip = true;
107  if (section == BKLMElementNumbers::c_BackwardSection &&
108  (sector == 4 || sector == 5 || sector == 6 || sector == 7))
109  dontFlip = true;
110 
111  for (int iStrip = 1; iStrip <= MaxiChannel; iStrip++) {
112  int channelId = iStrip;
113 
114  if (!(dontFlip && layer >= BKLMElementNumbers::c_FirstRPCLayer && plane == BKLMElementNumbers::c_PhiPlane))
115  channelId = MaxiChannel - iStrip + 1;
116 
117  if (plane == BKLMElementNumbers::c_PhiPlane) {
118  if (layer == 1)
119  channelId += 4;
120  if (layer == 2)
121  channelId += 2;
122  }
123 
124  if (plane == BKLMElementNumbers::c_ZPlane) {
126  int channelCheck = channelId;
129  if (layer == 1) {
130  if (version == 1) {
131  if (channelCheck > 0 && channelCheck < 9)
132  channelId = 9 - channelId;
133  if (channelCheck > 8 && channelCheck < 24)
134  channelId = 54 - channelId;
135  if (channelCheck > 23 && channelCheck < 39)
136  channelId = 54 - channelId;
137  } else {
138  if (channelCheck > 0 && channelCheck < 9)
139  channelId = 9 - channelId; // 8 : 1
140  if (channelCheck > 8 && channelCheck < 24)
141  channelId = 39 - channelId; // 30 : 16
142  if (channelCheck > 23 && channelCheck < 39)
143  channelId = 69 - channelId; // 45 : 31
144  }
145  }
146  if (layer == 2) {
147  if (channelCheck > 0 && channelCheck < 10)
148  channelId = 10 - channelId;
149  if (channelCheck > 9 && channelCheck < 24)
150  channelId = 40 - channelId;
151  if (channelCheck > 23 && channelCheck < 39)
152  channelId = 69 - channelId;
153  }
154  } else { // All the sectors except the chimney one
155  if (channelCheck > 0 && channelCheck < 10)
156  channelId = 10 - channelId;
157  if (channelCheck > 9 && channelCheck < 25)
158  channelId = 40 - channelId;
159  if (channelCheck > 24 && channelCheck < 40)
160  channelId = 70 - channelId;
161  if (channelCheck > 39 && channelCheck < 55)
162  channelId = 100 - channelId;
163  }
164  }
165  }
166 
167  uint16_t detectorChannel = m_ElementNumbers->channelNumberBKLM(
168  section, sector, layer, plane, iStrip);
169  m_ChannelMap.insert(
170  std::pair<uint16_t, KLMElectronicsChannel>(
171  detectorChannel,
172  KLMElectronicsChannel(copperId, slotId, laneId, axisId, channelId)));
173  }
174  }
175 }
176 
178 {
179  int segment, strip;
180  /* Segment is 0-based here. */
181  segment = (stripSoftware - 1) / EKLMElementNumbers::getNStripsSegment();
182  /* Order of segment readout boards in the firmware is opposite. */
183  segment = 4 - segment;
184  strip = segment * EKLMElementNumbers::getNStripsSegment() +
185  (stripSoftware - 1) % EKLMElementNumbers::getNStripsSegment() + 1;
186  return strip;
187 }
188 
190  int section, int sector, int layer, int copper, int slot, int lane)
191 {
192  for (int plane = 1; plane <= EKLMElementNumbers::getMaximalPlaneNumber(); ++plane) {
193  for (int strip = 1; strip <= EKLMElementNumbers::getMaximalStripNumber(); ++strip) {
194  int axis = plane - 1;
195  int channel = getEKLMStripFirmwareBySoftware(strip);
196  uint16_t detectorChannel = m_ElementNumbers->channelNumberEKLM(
197  section, sector, layer, plane, strip);
198  m_ChannelMap.insert(
199  std::pair<uint16_t, KLMElectronicsChannel>(
200  detectorChannel,
201  KLMElectronicsChannel(EKLM_ID + copper, slot, lane, axis, channel)));
202  }
203  }
204 }
205 
207  int section, int sector, int layer, int plane, int segment,
208  int firmwareSegment)
209 {
210  std::map<uint16_t, KLMElectronicsChannel>:: iterator it;
211  for (int strip = 1; strip <= EKLMElementNumbers::getNStripsSegment(); ++strip) {
212  int stripPlane = strip + EKLMElementNumbers::getNStripsSegment() * (segment - 1);
213  uint16_t detectorChannel = m_ElementNumbers->channelNumberEKLM(
214  section, sector, layer, plane, stripPlane);
215  it = m_ChannelMap.find(detectorChannel);
216  if (it == m_ChannelMap.end())
217  B2FATAL("The KLM electronics map is not loaded or incomplete.");
218  int channel = (firmwareSegment - 1) * EKLMElementNumbers::getNStripsSegment() +
219  strip;
220  it->second.setChannel(channel);
221  }
222 }
223 
225 {
226  const int minimalVersion = 1;
227  const int maximalVersion = 2;
228  if (version < minimalVersion || version > maximalVersion) {
229  B2FATAL("Incorrect version (" << version << ") of EKLM electronics map. "
230  "It must be from " << minimalVersion << " to " << maximalVersion);
231  }
232  /* Backward section. */
233  /* Sector 1. */
243  /* Wrong connection was fixed between phase 2 and phase 3. */
244  if (mc || (version >= 2)) {
247  } else {
250  }
252  /* Sector 2. */
253  if (version == 1) {
266  } else {
279  }
280  /* Sector 3. */
281  if (version == 1) {
294  } else {
307  }
308  /* Sector 4. */
321  /* Forward section. */
322  /* Sector 1. */
337  /* Sector 2. */
352  /* Sector 3. */
367  /* Sector 4. */
382  /* Switch of internal cables. */
383  if (!mc) {
384  setChannelsEKLMSegment(1, 1, 5, 1, 1, 4);
385  setChannelsEKLMSegment(1, 1, 5, 1, 2, 5);
386  setChannelsEKLMSegment(1, 1, 5, 2, 1, 4);
387  setChannelsEKLMSegment(1, 1, 5, 2, 2, 5);
388  setChannelsEKLMSegment(2, 1, 10, 2, 3, 2);
389  setChannelsEKLMSegment(2, 1, 10, 2, 4, 3);
390  }
391 }
392 
394  int subdetector, int section, int sector, int layer, int lane)
395 {
396  std::map<uint16_t, KLMElectronicsChannel>::iterator it;
397  int minimalPlane = m_ElementNumbers->getMinimalPlaneNumber(subdetector);
398  KLMChannelIndex klmChannel(subdetector, section, sector, layer, minimalPlane, 1);
399  KLMChannelIndex klmModule(klmChannel);
401  KLMChannelIndex klmNextModule(klmModule);
402  ++klmNextModule;
403  for (; klmChannel != klmNextModule; ++klmChannel) {
404  uint16_t channel = klmChannel.getKLMChannelNumber();
405  it = m_ChannelMap.find(channel);
406  if (it == m_ChannelMap.end())
407  B2FATAL("The KLM electronics map is not loaded or incomplete.");
408  it->second.setLane(lane);
409  }
410 }
411 
413  int subdetector, int section, int sector, int layer, int plane, int lane)
414 {
415  std::map<uint16_t, KLMElectronicsChannel>::iterator it;
416  KLMChannelIndex klmChannel(subdetector, section, sector, layer, plane, 1);
417  KLMChannelIndex klmPlane(klmChannel);
419  KLMChannelIndex klmNextPlane(klmPlane);
420  ++klmNextPlane;
421  for (; klmChannel != klmNextPlane; ++klmChannel) {
422  uint16_t channel = klmChannel.getKLMChannelNumber();
423  it = m_ChannelMap.find(channel);
424  if (it == m_ChannelMap.end())
425  B2FATAL("The KLM electronics map is not loaded or incomplete.");
426  it->second.setLane(lane);
427  }
428 }
429 
431 {
432  DBImportObjPtr<KLMElectronicsMap> electronicsMap;
433  electronicsMap.construct();
434  std::map<uint16_t, KLMElectronicsChannel>::iterator it;
435  for (it = m_ChannelMap.begin(); it != m_ChannelMap.end(); ++it) {
436  electronicsMap->addChannel(
437  it->first,
438  it->second.getCopper(),
439  it->second.getSlot(),
440  it->second.getLane(),
441  it->second.getAxis(),
442  it->second.getChannel());
443  }
446  electronicsMap.import(iov);
447 }
Belle2::KLMElectronicsMapImporter::clearElectronicsMap
void clearElectronicsMap()
Clear electronics map (to be able to import its multiple versions).
Definition: KLMElectronicsMapImporter.cc:44
Belle2::IntervalOfValidity
A class that describes the interval of experiments/runs for which an object in the database is valid.
Definition: IntervalOfValidity.h:35
Belle2::KLMElectronicsMapImporter::addEKLMLane
void addEKLMLane(int section, int sector, int layer, int copper, int slot, int lane)
Add EKLM electronics map lane.
Definition: KLMElectronicsMapImporter.cc:189
Belle2::KLMElectronicsMapImporter::getEKLMStripFirmwareBySoftware
int getEKLMStripFirmwareBySoftware(int stripSoftware) const
Get EKLM firmware strip number by software strip number.
Definition: KLMElectronicsMapImporter.cc:177
Belle2::KLMChannelIndex::endBKLM
KLMChannelIndex & endBKLM()
Last channel for BKLM.
Definition: KLMChannelIndex.cc:195
Belle2::KLMElectronicsMapImporter::m_RunHigh
int m_RunHigh
High run.
Definition: KLMElectronicsMapImporter.h:163
Belle2::KLMChannelIndex::c_IndexLevelLayer
@ c_IndexLevelLayer
Layer.
Definition: KLMChannelIndex.h:52
Belle2::KLMElementNumbers::getMinimalPlaneNumber
int getMinimalPlaneNumber(int subdetector) const
Get minimal plane number.
Definition: KLMElementNumbers.cc:244
Belle2::KLMElectronicsMapImporter::KLMElectronicsMapImporter
KLMElectronicsMapImporter()
Constructor.
Definition: KLMElectronicsMapImporter.cc:26
Belle2::KLMElectronicsMapImporter::m_ElementNumbers
const KLMElementNumbers * m_ElementNumbers
Element numbers.
Definition: KLMElectronicsMapImporter.h:148
Belle2::BKLMElementNumbers::c_ZPlane
@ c_ZPlane
Z.
Definition: BKLMElementNumbers.h:81
Belle2::BKLMElementNumbers::c_ForwardSection
@ c_ForwardSection
Forward.
Definition: BKLMElementNumbers.h:47
Belle2::DBImportObjPtr::construct
void construct(Args &&... params)
Construct an object of type T in this DBImportObjPtr using the provided constructor arguments.
Definition: DBImportObjPtr.h:57
Belle2::BKLMElementNumbers::c_ChimneySector
@ c_ChimneySector
Chimney sector: BB3 in 1-based notation; BB2 in 0-based notation.
Definition: BKLMElementNumbers.h:61
Belle2::DBImportBase::import
bool import(const IntervalOfValidity &iov)
Import the object to database.
Definition: DBImportBase.cc:38
Belle2::KLMElectronicsMapImporter::setChannelsEKLMSegment
void setChannelsEKLMSegment(int section, int sector, int layer, int plane, int segment, int firmwareSegment)
Set channels for EKLM segment.
Definition: KLMElectronicsMapImporter.cc:206
Belle2::KLMElementNumbers::channelNumberBKLM
uint16_t channelNumberBKLM(int section, int sector, int layer, int plane, int strip) const
Get channel number for BKLM.
Definition: KLMElementNumbers.cc:50
Belle2::EKLMElementNumbers::c_BackwardSection
@ c_BackwardSection
Backward.
Definition: EKLMElementNumbers.h:44
Belle2::KLMElectronicsMapImporter::loadEKLMElectronicsMap
void loadEKLMElectronicsMap(int version, bool mc)
Load EKLM electronics map.
Definition: KLMElectronicsMapImporter.cc:224
Belle2::KLMChannelIndex::getKLMChannelNumber
uint16_t getKLMChannelNumber() const
Get KLM channel number.
Definition: KLMChannelIndex.cc:145
Belle2::KLMChannelIndex::beginBKLM
KLMChannelIndex beginBKLM()
First channel for BKLM.
Definition: KLMChannelIndex.cc:189
Belle2::KLMElectronicsMapImporter::m_ChannelMap
std::map< uint16_t, KLMElectronicsChannel > m_ChannelMap
Data for creation of the electronics map.
Definition: KLMElectronicsMapImporter.h:151
Belle2::KLMChannelIndex::c_IndexLevelPlane
@ c_IndexLevelPlane
Plane.
Definition: KLMChannelIndex.h:55
Belle2::KLMElectronicsMapImporter::m_RunLow
int m_RunLow
Low run.
Definition: KLMElectronicsMapImporter.h:157
Belle2::EKLMElementNumbers::getMaximalStripNumber
static constexpr int getMaximalStripNumber()
Get maximal strip number.
Definition: EKLMElementNumbers.h:329
Belle2::BKLMElementNumbers::c_BackwardSection
@ c_BackwardSection
Backward.
Definition: BKLMElementNumbers.h:44
Belle2::KLMElectronicsMap::addChannel
void addChannel(uint16_t detectorChannel, int copper, int slot, int lane, int axis, int channel)
Add channel.
Definition: KLMElectronicsMap.cc:47
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::KLMElectronicsMapImporter::loadBKLMElectronicsMap
void loadBKLMElectronicsMap(int version)
Load BKLM electronics map.
Definition: KLMElectronicsMapImporter.cc:49
Belle2::DBImportObjPtr
Class for importing a single object to the database.
Definition: DBImportObjPtr.h:33
Belle2::KLMElectronicsMapImporter::importElectronicsMap
void importElectronicsMap()
Import electronics map.
Definition: KLMElectronicsMapImporter.cc:430
Belle2::KLMElectronicsMapImporter::m_ExperimentHigh
int m_ExperimentHigh
High experiment.
Definition: KLMElectronicsMapImporter.h:160
Belle2::EKLMElementNumbers::getMaximalPlaneNumber
static constexpr int getMaximalPlaneNumber()
Get maximal plane number.
Definition: EKLMElementNumbers.h:313
Belle2::KLMElectronicsMapImporter::~KLMElectronicsMapImporter
~KLMElectronicsMapImporter()
Destructor.
Definition: KLMElectronicsMapImporter.cc:31
Belle2::EKLMElementNumbers::getNStripsSegment
static constexpr int getNStripsSegment()
Get number of strips in a segment.
Definition: EKLMElementNumbers.h:401
Belle2::KLMElementNumbers::channelNumberEKLM
uint16_t channelNumberEKLM(int section, int sector, int layer, int plane, int strip) const
Get channel number for EKLM.
Definition: KLMElementNumbers.cc:64
Belle2::KLMChannelIndex
KLM channel index.
Definition: KLMChannelIndex.h:33
Belle2::KLMElectronicsChannel
BKLM electronics channel.
Definition: KLMElectronicsChannel.h:33
Belle2::KLMElectronicsMapImporter::setLane
void setLane(int subdetector, int section, int sector, int layer, int lane)
Set non-default lane for all channels in a module.
Definition: KLMElectronicsMapImporter.cc:393
Belle2::KLMElementNumbers
KLM element numbers.
Definition: KLMElementNumbers.h:37
Belle2::BKLMElementNumbers::c_FirstRPCLayer
@ c_FirstRPCLayer
First RPC layer.
Definition: BKLMElementNumbers.h:71
Belle2::EKLMElementNumbers::c_ForwardSection
@ c_ForwardSection
Forward.
Definition: EKLMElementNumbers.h:47
Belle2::KLMElectronicsMapImporter::setIOV
void setIOV(int experimentLow, int runLow, int experimentHigh, int runHigh)
Set interval of validity.
Definition: KLMElectronicsMapImporter.cc:35
Belle2::KLMChannelIndex::setIndexLevel
void setIndexLevel(enum IndexLevel indexLevel)
Set index level.
Definition: KLMChannelIndex.cc:67
Belle2::KLMElectronicsMapImporter::m_ExperimentLow
int m_ExperimentLow
Low experiment.
Definition: KLMElectronicsMapImporter.h:154
Belle2::BKLMElementNumbers::getNStrips
static int getNStrips(int section, int sector, int layer, int plane)
Get number of strips.
Definition: BKLMElementNumbers.cc:108
Belle2::BKLMElementNumbers::c_PhiPlane
@ c_PhiPlane
Phi.
Definition: BKLMElementNumbers.h:84