File indexing completed on 2023-03-17 13:02:42
0001 #include "CSCGeometryBuilder.h"
0002
0003 #include <Geometry/CSCGeometry/interface/CSCGeometry.h>
0004
0005 #include <DataFormats/DetId/interface/DetId.h>
0006 #include <DataFormats/MuonDetId/interface/CSCDetId.h>
0007 #include <Geometry/CSCGeometry/interface/CSCWireGroupPackage.h>
0008
0009 #include <FWCore/MessageLogger/interface/MessageLogger.h>
0010
0011 #include <vector>
0012
0013 CSCGeometryBuilder::CSCGeometryBuilder() : myName("CSCGeometryBuilder") {}
0014
0015 CSCGeometryBuilder::~CSCGeometryBuilder() {}
0016
0017 void CSCGeometryBuilder::build(CSCGeometry& theGeometry,
0018 const RecoIdealGeometry& rig,
0019 const CSCRecoDigiParameters& cscpars) {
0020 std::vector<float> fpar;
0021 std::vector<float> gtran;
0022 std::vector<float> grmat;
0023 std::vector<float> fupar;
0024 std::vector<double>::const_iterator it, endIt;
0025 const std::vector<DetId>& detids(rig.detIds());
0026
0027 for (size_t idt = 0; idt < detids.size(); ++idt) {
0028 CSCDetId detid = CSCDetId(detids[idt]);
0029 int jstation = detid.station();
0030 int jring = detid.ring();
0031
0032 endIt = rig.shapeEnd(idt);
0033 fpar.clear();
0034 for (it = rig.shapeStart(idt); it != endIt; ++it) {
0035 fpar.emplace_back((float)(*it));
0036 }
0037
0038 gtran.clear();
0039 endIt = rig.tranEnd(idt);
0040 for (it = rig.tranStart(idt); it != endIt; ++it) {
0041 gtran.emplace_back((float)(*it));
0042 }
0043 grmat.clear();
0044 endIt = rig.rotEnd(idt);
0045 for (it = rig.rotStart(idt); it != endIt; ++it) {
0046 grmat.emplace_back((float)(*it));
0047 }
0048
0049
0050 int chamberType = CSCChamberSpecs::whatChamberType(jstation, jring);
0051 size_t cs = 0;
0052 assert(!cscpars.pChamberType.empty());
0053 while (cs < cscpars.pChamberType.size() && chamberType != cscpars.pChamberType[cs]) {
0054 ++cs;
0055 }
0056 assert(cs != cscpars.pChamberType.size());
0057
0058
0059 size_t fu, numfuPars;
0060 CSCWireGroupPackage wg;
0061 fu = cscpars.pUserParOffset[cs];
0062 numfuPars = fu + 1 + size_t(cscpars.pfupars[fu]);
0063
0064
0065 LogTrace(myName) << myName << ": I think I have " << cscpars.pUserParSize[cs] << " values in pfupars (uparvals)."
0066 << std::endl;
0067 LogTrace(myName) << myName << ": For fupar I will start at " << cscpars.pUserParOffset[cs] + 1
0068 << " in pfupars and go to " << numfuPars << "." << std::endl;
0069 for (++fu; fu < numfuPars; ++fu) {
0070 LogTrace(myName) << myName << ": pfupars[" << fu << "]=" << cscpars.pfupars[fu] << std::endl;
0071 fupar.emplace_back(cscpars.pfupars[fu]);
0072 }
0073
0074 wg.wireSpacing = cscpars.pfupars[fu++];
0075 wg.alignmentPinToFirstWire = cscpars.pfupars[fu++];
0076 wg.numberOfGroups = int(cscpars.pfupars[fu++]);
0077 wg.narrowWidthOfWirePlane = cscpars.pfupars[fu++];
0078 wg.wideWidthOfWirePlane = cscpars.pfupars[fu++];
0079 wg.lengthOfWirePlane = cscpars.pfupars[fu++];
0080 size_t numgrp = static_cast<size_t>(cscpars.pfupars[fu]);
0081 size_t maxFu = fu + 1 + numgrp;
0082 fu++;
0083 for (; fu < maxFu; ++fu) {
0084 wg.wiresInEachGroup.emplace_back(int(cscpars.pfupars[fu]));
0085 }
0086 maxFu = fu + numgrp;
0087 for (; fu < maxFu; ++fu) {
0088 wg.consecutiveGroups.emplace_back(int(cscpars.pfupars[fu]));
0089 }
0090
0091 if (wg.numberOfGroups != 0) {
0092 LogTrace(myName) << myName << ": TotNumWireGroups = " << wg.numberOfGroups;
0093 LogTrace(myName) << myName << ": WireSpacing = " << wg.wireSpacing;
0094 LogTrace(myName) << myName << ": AlignmentPinToFirstWire = " << wg.alignmentPinToFirstWire;
0095 LogTrace(myName) << myName << ": Narrow width of wire plane = " << wg.narrowWidthOfWirePlane;
0096 LogTrace(myName) << myName << ": Wide width of wire plane = " << wg.wideWidthOfWirePlane;
0097 LogTrace(myName) << myName << ": Length in y of wire plane = " << wg.lengthOfWirePlane;
0098 LogTrace(myName) << myName << ": wg.consecutiveGroups.size() = " << wg.consecutiveGroups.size();
0099 LogTrace(myName) << myName << ": wg.wiresInEachGroup.size() = " << wg.wiresInEachGroup.size();
0100 LogTrace(myName) << myName << ": \tNumGroups\tWiresInGroup";
0101 for (size_t i = 0; i < wg.consecutiveGroups.size(); i++) {
0102 LogTrace(myName) << myName << " \t" << wg.consecutiveGroups[i] << "\t\t" << wg.wiresInEachGroup[i];
0103 }
0104 } else {
0105 LogTrace(myName) << myName << ": DDD is MISSING SpecPars for wire groups";
0106 }
0107 LogTrace(myName) << myName << ": end of wire group info. ";
0108
0109
0110 if (!theGeometry.centreTIOffsets())
0111 fupar[30] = 0.;
0112
0113 buildChamber(theGeometry, detid, fpar, fupar, gtran, grmat, wg);
0114 fupar.clear();
0115 }
0116 }
0117
0118 void CSCGeometryBuilder::buildChamber(CSCGeometry& theGeometry,
0119 CSCDetId chamberId,
0120 const std::vector<float>& fpar,
0121 const std::vector<float>& fupar,
0122 const std::vector<float>& gtran,
0123 const std::vector<float>& grmat,
0124 const CSCWireGroupPackage& wg
0125 ) {
0126 LogTrace(myName) << myName << ": entering buildChamber";
0127
0128 int jend = chamberId.endcap();
0129 int jstat = chamberId.station();
0130 int jring = chamberId.ring();
0131 int jch = chamberId.chamber();
0132 int jlay = chamberId.layer();
0133
0134 if (jlay != 0)
0135 edm::LogWarning(myName) << "Error! CSCGeometryBuilderFromDDD was fed layer id = " << jlay << "\n";
0136
0137 const size_t kNpar = 4;
0138 if (fpar.size() != kNpar)
0139 edm::LogError(myName) << "Error, expected npar=" << kNpar << ", found npar=" << fpar.size() << std::endl;
0140
0141 LogTrace(myName) << myName << ": E" << jend << " S" << jstat << " R" << jring << " C" << jch << " L" << jlay;
0142 LogTrace(myName) << myName << ": npar=" << fpar.size() << " hB=" << fpar[0] << " hT=" << fpar[1] << " hD=" << fpar[2]
0143 << " hH=" << fpar[3];
0144 LogTrace(myName) << myName << ": gtran[0,1,2]=" << gtran[0] << " " << gtran[1] << " " << gtran[2];
0145 LogTrace(myName) << myName << ": grmat[0-8]=" << grmat[0] << " " << grmat[1] << " " << grmat[2] << " " << grmat[3]
0146 << " " << grmat[4] << " " << grmat[5] << " " << grmat[6] << " " << grmat[7] << " " << grmat[8];
0147 LogTrace(myName) << myName << ": nupar=" << fupar.size() << " upar[0]=" << fupar[0] << " upar[" << fupar.size() - 1
0148 << "]=" << fupar[fupar.size() - 1];
0149
0150 const CSCChamber* chamber = theGeometry.chamber(chamberId);
0151 if (chamber) {
0152 } else {
0153
0154 LogTrace(myName) << myName << ": CSCChamberSpecs::build requested for ME" << jstat << jring;
0155 int chamberType = CSCChamberSpecs::whatChamberType(jstat, jring);
0156 const CSCChamberSpecs* aSpecs = theGeometry.findSpecs(chamberType);
0157 if (!fupar.empty() && aSpecs == nullptr) {
0158
0159 aSpecs = theGeometry.buildSpecs(chamberType, fpar, fupar, wg);
0160 } else if (fupar.empty() && aSpecs == nullptr) {
0161 edm::LogError(myName)
0162 << "SHOULD BE THROW? Error, wg and/or fupar size are 0 BUT this Chamber Spec has not been built!";
0163 }
0164
0165
0166
0167
0168
0169
0170
0171
0172 Surface::RotationType aRot(
0173 grmat[0], grmat[1], grmat[2], grmat[3], grmat[4], grmat[5], grmat[6], grmat[7], grmat[8]);
0174
0175
0176
0177
0178
0179
0180
0181
0182
0183
0184
0185
0186
0187
0188 Basic3DVector<float> oldX(1., 0., 0.);
0189 Basic3DVector<float> oldY(0., 0., -1.);
0190 Basic3DVector<float> oldZ(0., 1., 0.);
0191
0192 if (gtran[2] < 0.)
0193 oldX *= -1;
0194
0195 aRot.rotateAxes(oldX, oldY, oldZ);
0196
0197
0198
0199 float frameThickness = fupar[31] / 10.;
0200 float gapThickness = fupar[32] / 10.;
0201 float panelThickness = fupar[33] / 10.;
0202 float zAverageAGVtoAF = fupar[34] / 10.;
0203
0204 float layerThickness = gapThickness;
0205 float layerSeparation = gapThickness + panelThickness;
0206
0207 float chamberThickness = 7. * panelThickness + 6. * gapThickness + 2. * frameThickness;
0208 float hChamberThickness = chamberThickness / 2.;
0209
0210
0211
0212
0213
0214
0215
0216
0217
0218
0219
0220
0221
0222
0223
0224
0225
0226
0227 TrapezoidalPlaneBounds* bounds = new TrapezoidalPlaneBounds(fpar[0], fpar[1], fpar[3], hChamberThickness);
0228
0229
0230 Surface::PositionType aVec(gtran[0], gtran[1], gtran[2]);
0231
0232 Plane::PlanePointer plane = Plane::build(aVec, aRot, bounds);
0233
0234 CSCChamber* chamber = new CSCChamber(plane, chamberId, aSpecs);
0235 theGeometry.addChamber(chamber);
0236
0237 LogTrace(myName) << myName << ": Create chamber E" << jend << " S" << jstat << " R" << jring << " C" << jch
0238 << " z=" << gtran[2] << " t/2=" << fpar[2] << " (DDD) or " << hChamberThickness
0239 << " (specs) adr=" << chamber;
0240
0241
0242
0243
0244
0245
0246
0247 int localZwrtGlobalZ = +1;
0248 if ((jend == 1 && jstat < 3) || (jend == 2 && jstat > 2))
0249 localZwrtGlobalZ = -1;
0250 int globalZ = +1;
0251 if (jend == 2)
0252 globalZ = -1;
0253
0254 LogTrace(myName) << myName << ": layerSeparation=" << layerSeparation << ", zAF-zAverageAGV=" << zAverageAGVtoAF
0255 << ", localZwrtGlobalZ=" << localZwrtGlobalZ << ", gtran[2]=" << gtran[2];
0256
0257 for (short j = 1; j <= 6; ++j) {
0258 CSCDetId layerId = CSCDetId(jend, jstat, jring, jch, j);
0259
0260
0261 const CSCLayer* cLayer = dynamic_cast<const CSCLayer*>(theGeometry.idToDet(layerId));
0262
0263 if (cLayer == nullptr) {
0264
0265 const CSCChamberSpecs* aSpecs = chamber->specs();
0266 const CSCLayerGeometry* geom = (j % 2 != 0) ? aSpecs->oddLayerGeometry(jend) : aSpecs->evenLayerGeometry(jend);
0267
0268
0269
0270
0271 float zlayer = gtran[2] - globalZ * zAverageAGVtoAF + localZwrtGlobalZ * (3.5 - j) * layerSeparation;
0272
0273 Surface::RotationType chamberRotation = chamber->surface().rotation();
0274 Surface::PositionType layerPosition(gtran[0], gtran[1], zlayer);
0275 std::array<const float, 4> const& dims = geom->parameters();
0276
0277 TrapezoidalPlaneBounds* bounds = new TrapezoidalPlaneBounds(dims[0], dims[1], dims[3], layerThickness / 2.);
0278 Plane::PlanePointer plane = Plane::build(layerPosition, chamberRotation, bounds);
0279
0280 CSCLayer* layer = new CSCLayer(plane, layerId, chamber, geom);
0281
0282 LogTrace(myName) << myName << ": Create layer E" << jend << " S" << jstat << " R" << jring << " C" << jch
0283 << " L" << j << " z=" << zlayer << " t=" << layerThickness << " or "
0284 << layer->surface().bounds().thickness() << " adr=" << layer << " layerGeom adr=" << geom;
0285
0286 chamber->addComponent(j, layer);
0287 theGeometry.addLayer(layer);
0288 } else {
0289 edm::LogError(myName) << ": ERROR, layer " << j << " for chamber = " << (chamber->id())
0290 << " already exists: layer address=" << cLayer << " chamber address=" << chamber << "\n";
0291 }
0292
0293 }
0294 }
0295 }