File indexing completed on 2023-03-29 01:42:50
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014 #include "Geometry/DTGeometryBuilder/src/DTGeometryBuilderFromCondDB.h"
0015
0016
0017 #include <CondFormats/GeometryObjects/interface/RecoIdealGeometry.h>
0018 #include <Geometry/DTGeometry/interface/DTGeometry.h>
0019 #include <DataFormats/MuonDetId/interface/DTChamberId.h>
0020 #include <DataFormats/MuonDetId/interface/DTSuperLayerId.h>
0021 #include <DataFormats/MuonDetId/interface/DTLayerId.h>
0022 #include "DataFormats/Math/interface/GeantUnits.h"
0023 #include <FWCore/MessageLogger/interface/MessageLogger.h>
0024
0025
0026 #include <iostream>
0027 using namespace std;
0028
0029 using namespace geant_units;
0030 using namespace geant_units::operators;
0031
0032
0033
0034
0035
0036 DTGeometryBuilderFromCondDB::DTGeometryBuilderFromCondDB() {}
0037
0038
0039 DTGeometryBuilderFromCondDB::~DTGeometryBuilderFromCondDB() {}
0040
0041
0042 void DTGeometryBuilderFromCondDB::build(const std::shared_ptr<DTGeometry>& theGeometry, const RecoIdealGeometry& rig) {
0043 #ifdef EDM_ML_DEBUG
0044 edm::LogVerbatim("DTGeometry") << "DTGeometryBuilderFromCondDB ";
0045 #endif
0046 const std::vector<DetId>& detids(rig.detIds());
0047 #ifdef EDM_ML_DEBUG
0048 edm::LogVerbatim("DTGeometry") << "size " << detids.size();
0049 #endif
0050
0051 size_t idt = 0;
0052 DTChamber* chamber(nullptr);
0053 DTSuperLayer* sl(nullptr);
0054 while (idt < detids.size()) {
0055
0056 if (int(*(rig.shapeStart(idt))) == 0) {
0057
0058
0059 if (chamber)
0060 theGeometry->add(chamber);
0061
0062 DTChamberId chid(detids[idt]);
0063 #ifdef EDM_ML_DEBUG
0064 edm::LogVerbatim("DTGeometry") << "CH: " << chid;
0065 #endif
0066 chamber = buildChamber(chid, rig, idt);
0067 } else if (int(*(rig.shapeStart(idt))) == 1) {
0068 DTSuperLayerId slid(detids[idt]);
0069 #ifdef EDM_ML_DEBUG
0070 edm::LogVerbatim("DTGeometry") << " SL: " << slid;
0071 #endif
0072 sl = buildSuperLayer(chamber, slid, rig, idt);
0073 theGeometry->add(sl);
0074 } else if (int(*(rig.shapeStart(idt))) == 2) {
0075 DTLayerId lid(detids[idt]);
0076 #ifdef EDM_ML_DEBUG
0077 edm::LogVerbatim("DTGeometry") << " LAY: " << lid;
0078 #endif
0079 DTLayer* lay = buildLayer(sl, lid, rig, idt);
0080 theGeometry->add(lay);
0081 } else {
0082 edm::LogVerbatim("DTGeometry") << "What is this?";
0083 }
0084 ++idt;
0085 }
0086 if (chamber)
0087 theGeometry->add(chamber);
0088 }
0089
0090
0091 RectangularPlaneBounds* dtGeometryBuilder::getRecPlaneBounds(const std::vector<double>::const_iterator& shapeStart) {
0092 float width = convertMmToCm(*(shapeStart));
0093 float length = convertMmToCm(*(shapeStart + 1));
0094 float thickness = convertMmToCm(*(shapeStart + 2));
0095 return new RectangularPlaneBounds(width, length, thickness);
0096 }
0097
0098 DTChamber* DTGeometryBuilderFromCondDB::buildChamber(const DetId& id, const RecoIdealGeometry& rig, size_t idt) const {
0099 DTChamberId detId(id);
0100
0101
0102
0103
0104
0105
0106
0107 RCPPlane surf(
0108 plane(rig.tranStart(idt), rig.rotStart(idt), dtGeometryBuilder::getRecPlaneBounds(++rig.shapeStart(idt))));
0109
0110 DTChamber* chamber = new DTChamber(detId, surf);
0111
0112 return chamber;
0113 }
0114
0115 DTSuperLayer* DTGeometryBuilderFromCondDB::buildSuperLayer(DTChamber* chamber,
0116 const DetId& id,
0117 const RecoIdealGeometry& rig,
0118 size_t idt) const {
0119 DTSuperLayerId slId(id);
0120
0121
0122
0123
0124
0125
0126 RCPPlane surf(
0127 plane(rig.tranStart(idt), rig.rotStart(idt), dtGeometryBuilder::getRecPlaneBounds(++rig.shapeStart(idt))));
0128
0129 DTSuperLayer* slayer = new DTSuperLayer(slId, surf, chamber);
0130
0131 #ifdef EDM_ML_DEBUG
0132 edm::LogVerbatim("DTGeometry") << "adding slayer " << slayer->id() << " to chamber " << chamber->id();
0133 #endif
0134 assert(chamber);
0135 chamber->add(slayer);
0136 return slayer;
0137 }
0138
0139 DTLayer* DTGeometryBuilderFromCondDB::buildLayer(DTSuperLayer* sl,
0140 const DetId& id,
0141 const RecoIdealGeometry& rig,
0142 size_t idt) const {
0143 DTLayerId layId(id);
0144
0145
0146
0147
0148
0149
0150 auto shapeStartPtr = rig.shapeStart(idt);
0151 RCPPlane surf(
0152 plane(rig.tranStart(idt), rig.rotStart(idt), dtGeometryBuilder::getRecPlaneBounds((shapeStartPtr + 1))));
0153
0154
0155 int firstWire = static_cast<int>(*(shapeStartPtr + 4));
0156 int WCounter = static_cast<int>(*(shapeStartPtr + 5));
0157 double sensibleLength = convertMmToCm(*(shapeStartPtr + 6));
0158 DTTopology topology(firstWire, WCounter, sensibleLength);
0159
0160 DTLayerType layerType;
0161
0162 DTLayer* layer = new DTLayer(layId, surf, topology, layerType, sl);
0163 #ifdef EDM_ML_DEBUG
0164 edm::LogVerbatim("DTGeometry") << "adding layer " << layer->id() << " to sl " << sl->id();
0165 #endif
0166 assert(sl);
0167 sl->add(layer);
0168 return layer;
0169 }
0170
0171 DTGeometryBuilderFromCondDB::RCPPlane DTGeometryBuilderFromCondDB::plane(const vector<double>::const_iterator tranStart,
0172 const vector<double>::const_iterator rotStart,
0173 Bounds* bounds) const {
0174
0175 const Surface::PositionType posResult(*(tranStart), *(tranStart + 1), *(tranStart + 2));
0176
0177 Surface::RotationType rotResult(*(rotStart + 0),
0178 *(rotStart + 1),
0179 *(rotStart + 2),
0180 *(rotStart + 3),
0181 *(rotStart + 4),
0182 *(rotStart + 5),
0183 *(rotStart + 6),
0184 *(rotStart + 7),
0185 *(rotStart + 8));
0186
0187 return RCPPlane(new Plane(posResult, rotResult, bounds));
0188 }