CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
CSCGeometryBuilder.cc
Go to the documentation of this file.
1 #include "CSCGeometryBuilder.h"
2 
4 
8 
10 
11 #include <vector>
12 
13 CSCGeometryBuilder::CSCGeometryBuilder() : myName("CSCGeometryBuilder"){}
14 
15 
17 
18 
19 void CSCGeometryBuilder::build( const std::shared_ptr<CSCGeometry>& theGeometry
20  , const RecoIdealGeometry& rig
21  , const CSCRecoDigiParameters& cscpars ) {
22 
23  std::vector<float> fpar;
24  std::vector<float> gtran;
25  std::vector<float> grmat;
26  std::vector<float> fupar;
27  std::vector<double>::const_iterator it, endIt;
28  const std::vector<DetId>& detids(rig.detIds());
29 
30  for ( size_t idt = 0; idt < detids.size(); ++idt) {
31  CSCDetId detid = CSCDetId( detids[idt] );
32  int jstation = detid.station();
33  int jring = detid.ring();
34 
35  endIt = rig.shapeEnd(idt);
36  fpar.clear();
37  for ( it = rig.shapeStart(idt); it != endIt; ++it) {
38  fpar.emplace_back( (float)(*it) );
39  }
40 
41  gtran.clear();
42  endIt = rig.tranEnd(idt);
43  for ( it = rig.tranStart(idt); it != endIt; ++it ) {
44  gtran.emplace_back((float)(*it));
45  }
46  grmat.clear();
47  endIt = rig.rotEnd(idt);
48  for ( it = rig.rotStart(idt) ; it != endIt; ++it ) {
49  grmat.emplace_back((float)(*it));
50  }
51 
52  // get the chamber type from existing info
53  int chamberType = CSCChamberSpecs::whatChamberType( jstation, jring );
54  size_t cs = 0;
55  assert ( !cscpars.pChamberType.empty() );
56  while (cs < cscpars.pChamberType.size() && chamberType != cscpars.pChamberType[cs]) {
57  ++cs;
58  }
59  assert ( cs != cscpars.pChamberType.size() );
60 
61  // check the existence of the specs for this type WHY? Remove it...
62  size_t fu, numfuPars;
64  fu = cscpars.pUserParOffset[cs];
65  numfuPars = fu + 1 + size_t(cscpars.pfupars[fu]);
66 
67  // upars from db are now uparvals + wg info so we need to unwrap only part here first...
68  LogTrace(myName) << myName << ": I think I have " << cscpars.pUserParSize[cs] << " values in pfupars (uparvals)." << std::endl;
69  LogTrace(myName) << myName << ": For fupar I will start at " << cscpars.pUserParOffset[cs] + 1
70  << " in pfupars and go to " << numfuPars << "." << std::endl;
71  for ( ++fu; fu < numfuPars; ++fu ) {
72  LogTrace(myName) << myName << ": pfupars[" << fu << "]=" << cscpars.pfupars[fu] << std::endl;
73  fupar.emplace_back(cscpars.pfupars[fu]);
74  }
75  // now, we need to start from "here" at fu to go on and build wg...
76  wg.wireSpacing = cscpars.pfupars[fu++];
77  wg.alignmentPinToFirstWire = cscpars.pfupars[fu++];
78  wg.numberOfGroups = int(cscpars.pfupars[fu++]);
79  wg.narrowWidthOfWirePlane = cscpars.pfupars[fu++];
80  wg.wideWidthOfWirePlane = cscpars.pfupars[fu++];
81  wg.lengthOfWirePlane = cscpars.pfupars[fu++];
82  size_t numgrp = static_cast<size_t>(cscpars.pfupars[fu]);
83  size_t maxFu = fu + 1 + numgrp;
84  fu++;
85  for ( ;fu < maxFu; ++fu ) {
86  wg.wiresInEachGroup.emplace_back(int(cscpars.pfupars[fu]));
87  }
88  maxFu = fu + numgrp;
89  for ( ;fu < maxFu; ++fu ) {
90  wg.consecutiveGroups.emplace_back(int(cscpars.pfupars[fu]));
91  }
92 
93  if ( wg.numberOfGroups != 0 ) {
94  LogTrace(myName) << myName << ": TotNumWireGroups = " << wg.numberOfGroups ;
95  LogTrace(myName) << myName << ": WireSpacing = " << wg.wireSpacing ;
96  LogTrace(myName) << myName << ": AlignmentPinToFirstWire = " << wg.alignmentPinToFirstWire ;
97  LogTrace(myName) << myName << ": Narrow width of wire plane = " << wg.narrowWidthOfWirePlane ;
98  LogTrace(myName) << myName << ": Wide width of wire plane = " << wg.wideWidthOfWirePlane ;
99  LogTrace(myName) << myName << ": Length in y of wire plane = " << wg.lengthOfWirePlane ;
100  LogTrace(myName) << myName << ": wg.consecutiveGroups.size() = " << wg.consecutiveGroups.size() ;
101  LogTrace(myName) << myName << ": wg.wiresInEachGroup.size() = " << wg.wiresInEachGroup.size() ;
102  LogTrace(myName) << myName << ": \tNumGroups\tWiresInGroup" ;
103  for (size_t i = 0; i < wg.consecutiveGroups.size(); i++) {
104  LogTrace(myName) << myName << " \t" << wg.consecutiveGroups[i] << "\t\t" << wg.wiresInEachGroup[i] ;
105  }
106  } else {
107  LogTrace(myName) << myName << ": DDD is MISSING SpecPars for wire groups" ;
108  }
109  LogTrace(myName) << myName << ": end of wire group info. " ;
110 
111  // Are we going to apply centre-to-intersection offsets, even if values exist in the specs file?
112  if ( !theGeometry->centreTIOffsets() ) fupar[30] = 0.; // reset to zero if flagged 'off'
113 
114  buildChamber (theGeometry, detid, fpar, fupar, gtran, grmat, wg ); //, cscpars.pWGPs[cs] );
115  fupar.clear();
116  }
117 }
118 
120  const std::shared_ptr<CSCGeometry>& theGeometry // the geometry container
121  , CSCDetId chamberId // the DetId for this chamber
122  , const std::vector<float>& fpar // volume parameters hB, hT. hD, hH
123  , const std::vector<float>& fupar // user parameters
124  , const std::vector<float>& gtran // translation vector
125  , const std::vector<float>& grmat // rotation matrix
126  , const CSCWireGroupPackage& wg // wire group info
127  ) {
128 
129  LogTrace(myName) << myName << ": entering buildChamber" ;
130 
131  int jend = chamberId.endcap();
132  int jstat = chamberId.station();
133  int jring = chamberId.ring();
134  int jch = chamberId.chamber();
135  int jlay = chamberId.layer();
136 
137  if (jlay != 0 ) edm::LogWarning(myName) << "Error! CSCGeometryBuilderFromDDD was fed layer id = " << jlay << "\n";
138 
139  const size_t kNpar = 4;
140  if ( fpar.size() != kNpar )
141  edm::LogError(myName) << "Error, expected npar="
142  << kNpar << ", found npar=" << fpar.size() << std::endl;
143 
144  LogTrace(myName) << myName << ": E" << jend << " S" << jstat << " R" << jring <<
145  " C" << jch << " L" << jlay ;
146  LogTrace(myName) << myName << ": npar=" << fpar.size() << " hB=" << fpar[0]
147  << " hT=" << fpar[1] << " hD=" << fpar[2] << " hH=" << fpar[3] ;
148  LogTrace(myName) << myName << ": gtran[0,1,2]=" << gtran[0] << " " << gtran[1] << " " << gtran[2] ;
149  LogTrace(myName) << myName << ": grmat[0-8]=" << grmat[0] << " " << grmat[1] << " " << grmat[2] << " "
150  << grmat[3] << " " << grmat[4] << " " << grmat[5] << " "
151  << grmat[6] << " " << grmat[7] << " " << grmat[8] ;
152  LogTrace(myName) << myName << ": nupar=" << fupar.size() << " upar[0]=" << fupar[0]
153  << " upar[" << fupar.size()-1 << "]=" << fupar[fupar.size()-1];
154 
155 
156  const CSCChamber* chamber = theGeometry->chamber( chamberId );
157  if ( chamber ){
158  }
159  else { // this chamber not yet built/stored
160 
161  LogTrace(myName) << myName <<": CSCChamberSpecs::build requested for ME" << jstat << jring ;
162  int chamberType = CSCChamberSpecs::whatChamberType( jstat, jring );
163  const CSCChamberSpecs* aSpecs = theGeometry->findSpecs( chamberType );
164  if ( !fupar.empty() && aSpecs == nullptr ) {
165  // make new one:
166  aSpecs = theGeometry->buildSpecs (chamberType, fpar, fupar, wg);
167  } else if ( fupar.empty() && aSpecs == nullptr ) {
168  edm::LogError(myName) << "SHOULD BE THROW? Error, wg and/or fupar size are 0 BUT this Chamber Spec has not been built!";
169  }
170 
171  // Build a Transformation out of GEANT gtran and grmat...
172  // These are used to transform a point in the local reference frame
173  // of a subdetector to the global frame of CMS by
174  // (grmat)^(-1)*local + (gtran)
175  // The corresponding transformation from global to local is
176  // (grmat)*(global - gtran)
177 
178  Surface::RotationType aRot( grmat[0], grmat[1], grmat[2],
179  grmat[3], grmat[4], grmat[5],
180  grmat[6], grmat[7], grmat[8] );
181 
182  // This rotation from GEANT considers the detector face as the x-z plane.
183  // We want this to be the local x-y plane.
184  // Furthermore, the -z_global endcap has LH local coordinates, since it is built
185  // in GEANT as a *reflection* of the +z_global endcap.
186  // So we need to rotate, and in -z flip local x.
187 
188  // aRot.rotateAxes will transform aRot in place so that it becomes
189  // applicable to the new local coordinates: detector face in x-y plane
190  // looking out along z, in either endcap.
191 
192  // The interface for rotateAxes specifies 'new' X,Y,Z but the
193  // implementation deals with them as the 'old'.
194 
195  Basic3DVector<float> oldX( 1., 0., 0. );
196  Basic3DVector<float> oldY( 0., 0., -1. );
197  Basic3DVector<float> oldZ( 0., 1., 0. );
198 
199  if ( gtran[2]<0. ) oldX *= -1;
200 
201  aRot.rotateAxes( oldX, oldY, oldZ );
202 
203  // Need to know z of layers w.r.t to z of centre of chamber.
204 
205  float frameThickness = fupar[31]/10.; // mm -> cm
206  float gapThickness = fupar[32]/10.; // mm -> cm
207  float panelThickness = fupar[33]/10.; // mm -> cm
208  float zAverageAGVtoAF = fupar[34]/10.; // mm -> cm
209 
210  float layerThickness = gapThickness; // consider the layer to be the gas gap
211  float layerSeparation = gapThickness + panelThickness; // centre-to-centre of neighbouring layers
212 
213  float chamberThickness = 7.*panelThickness + 6.*gapThickness + 2.*frameThickness ; // chamber frame thickness
214  float hChamberThickness = chamberThickness/2.; // @@ should match value returned from DDD directly
215 
216  // distAverageAGVtoAF is offset between centre of chamber (AF) and (L1+L6)/2 (average AGVs)
217  // where AF = AluminumFrame and AGV=ActiveGasVolume (volume names in DDD).
218  // It is signed based on global z values: zc - (zl1+zl6)/2
219 
220  // Local z values w.r.t. AF...
221  // z of wires in layer 1 = z_w1 = +/- zAverageAGVtoAF + 2.5*layerSeparation; // layer 1 is at most +ve local z
222  // The sign in '+/-' depends on relative directions of local and global z.
223  // It is '-' if they are the same direction, and '+' if opposite directions.
224  // z of wires in layer N = z_wN = z_w1 - (N-1)*layerSeparation;
225  // z of strips in layer N = z_sN = z_wN + gapThickness/2.; @@ BEWARE: need to check if it should be '-gapThickness/2' !
226 
227  // Set dimensions of trapezoidal chamber volume
228  // N.B. apothem is 4th in fpar but 3rd in ctor
229 
230  // hChamberThickness and fpar[2] should be the same - but using the above value at least shows
231  // how chamber structure works
232 
233  TrapezoidalPlaneBounds* bounds = new TrapezoidalPlaneBounds( fpar[0], fpar[1], fpar[3], hChamberThickness );
234 
235  // Centre of chamber in z is specified in DDD
236  Surface::PositionType aVec( gtran[0], gtran[1], gtran[2] );
237 
238  Plane::PlanePointer plane = Plane::build(aVec, aRot, bounds);
239 
240  CSCChamber* chamber = new CSCChamber( plane, chamberId, aSpecs );
241  theGeometry->addChamber( chamber );
242 
243  LogTrace(myName) << myName << ": Create chamber E" << jend << " S" << jstat
244  << " R" << jring << " C" << jch
245  << " z=" << gtran[2]
246  << " t/2=" << fpar[2] << " (DDD) or " << hChamberThickness
247  << " (specs) adr=" << chamber ;
248 
249  // Create the component layers of this chamber
250  // We're taking the z as the z of the wire plane within the layer (middle of gas gap)
251 
252  // Specify global z of layer by offsetting from centre of chamber: since layer 1
253  // is nearest to IP in stations 1/2 but layer 6 is nearest in stations 3/4,
254  // we need to adjust sign of offset appropriately...
255  int localZwrtGlobalZ = +1;
256  if ( (jend==1 && jstat<3 ) || ( jend==2 && jstat>2 ) ) localZwrtGlobalZ = -1;
257  int globalZ = +1;
258  if ( jend == 2 ) globalZ = -1;
259 
260  LogTrace(myName) << myName << ": layerSeparation=" << layerSeparation
261  << ", zAF-zAverageAGV=" << zAverageAGVtoAF
262  << ", localZwrtGlobalZ=" << localZwrtGlobalZ
263  << ", gtran[2]=" << gtran[2] ;
264 
265  for ( short j = 1; j <= 6; ++j ) {
266 
267  CSCDetId layerId = CSCDetId( jend, jstat, jring, jch, j );
268 
269  // extra-careful check that we haven't already built this layer
270  const CSCLayer* cLayer = dynamic_cast<const CSCLayer*> (theGeometry->idToDet( layerId ) );
271 
272  if ( cLayer == nullptr ) {
273 
274  // build the layer - need the chamber's specs and an appropriate layer-geometry
275  const CSCChamberSpecs* aSpecs = chamber->specs();
276  const CSCLayerGeometry* geom =
277  (j%2 != 0) ? aSpecs->oddLayerGeometry( jend ) :
278  aSpecs->evenLayerGeometry( jend );
279 
280  // Build appropriate BoundPlane, based on parent chamber, with gas gap as thickness
281 
282  // centre of chamber is at global z = gtran[2]
283  float zlayer = gtran[2] - globalZ*zAverageAGVtoAF + localZwrtGlobalZ*(3.5-j)*layerSeparation;
284 
285  Surface::RotationType chamberRotation = chamber->surface().rotation();
286  Surface::PositionType layerPosition( gtran[0], gtran[1], zlayer );
287  std::array<const float, 4> const & dims = geom->parameters(); // returns hb, ht, d, a
288  // dims[2] = layerThickness/2.; // half-thickness required and note it is 3rd value in vector
289  TrapezoidalPlaneBounds* bounds = new TrapezoidalPlaneBounds( dims[0], dims[1], dims[3], layerThickness/2. );
290  Plane::PlanePointer plane = Plane::build(layerPosition, chamberRotation, bounds);
291 
292  CSCLayer* layer = new CSCLayer( plane, layerId, chamber, geom );
293 
294  LogTrace(myName) << myName << ": Create layer E" << jend << " S" << jstat
295  << " R" << jring << " C" << jch << " L" << j
296  << " z=" << zlayer
297  << " t=" << layerThickness << " or " << layer->surface().bounds().thickness()
298  << " adr=" << layer << " layerGeom adr=" << geom ;
299 
300  chamber->addComponent(j, layer);
301  theGeometry->addLayer( layer );
302  }
303  else {
304  edm::LogError(myName) << ": ERROR, layer " << j <<
305  " for chamber = " << ( chamber->id() ) <<
306  " already exists: layer address=" << cLayer <<
307  " chamber address=" << chamber << "\n";
308  }
309 
310  } // layer construction within chamber
311  } // chamber construction
312 }
int chamber() const
Definition: CSCDetId.h:68
void buildChamber(const std::shared_ptr< CSCGeometry > &theGeometry, CSCDetId chamberId, const std::vector< float > &fpar, const std::vector< float > &fupar, const std::vector< float > &gtran, const std::vector< float > &grmat, const CSCWireGroupPackage &wg)
Build one CSC chamber, and its component layers, and add them to the geometry.
auto_ptr< ClusterSequence > cs
std::vector< double >::const_iterator rotEnd(size_t ind) const
virtual const std::array< const float, 4 > parameters() const
static int whatChamberType(int istation, int iring)
std::vector< double >::const_iterator rotStart(size_t ind) const
const Bounds & bounds() const
Definition: Surface.h:120
#define nullptr
int layer() const
Definition: CSCDetId.h:61
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:42
TkRotation & rotateAxes(const Basic3DVector< T > &newX, const Basic3DVector< T > &newY, const Basic3DVector< T > &newZ)
std::vector< float > pfupars
int endcap() const
Definition: CSCDetId.h:93
std::vector< double >::const_iterator tranEnd(size_t ind) const
std::vector< int > pChamberType
static PlanePointer build(Args &&...args)
Definition: Plane.h:33
std::vector< double >::const_iterator tranStart(size_t ind) const
const CSCLayerGeometry * evenLayerGeometry(int iendcap) const
std::vector< int > pUserParOffset
void build(const std::shared_ptr< CSCGeometry > &theGeometry, const RecoIdealGeometry &rig, const CSCRecoDigiParameters &cscpars)
Build the geometry.
const std::string myName
#define LogTrace(id)
const std::vector< DetId > & detIds() const
int ring() const
Definition: CSCDetId.h:75
virtual float thickness() const =0
std::vector< int > pUserParSize
CSCGeometryBuilder()
Constructor.
const CSCLayerGeometry * oddLayerGeometry(int iendcap) const
Accessors for LayerGeometry&#39;s.
std::vector< double >::const_iterator shapeEnd(size_t ind) const
std::vector< double >::const_iterator shapeStart(size_t ind) const
virtual ~CSCGeometryBuilder()
Destructor.
int station() const
Definition: CSCDetId.h:86