CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
HcalGeometry.cc
Go to the documentation of this file.
4 #include <algorithm>
5 
6 #include <Math/Transform3D.h>
7 #include <Math/EulerAngles.h>
8 
13 
15  m_topology(topology), m_mergePosition(topology.getMergePositionFlag()) {
16  init();
17 }
18 
20 
23  edm::LogInfo("HcalGeometry") << "HcalGeometry::init() "
24  << " HBSize " << m_topology.getHBSize()
25  << " HESize " << m_topology.getHESize()
26  << " HOSize " << m_topology.getHOSize()
27  << " HFSize " << m_topology.getHFSize();
28 
33 }
34 
36  if(!m_emptyIds.isSet()) {
37  std::unique_ptr<std::vector<DetId>> p_hbIds{new std::vector<DetId>};
38  std::unique_ptr<std::vector<DetId>> p_heIds{new std::vector<DetId>};
39  std::unique_ptr<std::vector<DetId>> p_hoIds{new std::vector<DetId>};
40  std::unique_ptr<std::vector<DetId>> p_hfIds{new std::vector<DetId>};
41  std::unique_ptr<std::vector<DetId>> p_emptyIds{new std::vector<DetId>};
42 
43  const std::vector<DetId>& baseIds (CaloSubdetectorGeometry::getValidDetIds());
44  for (unsigned int i ( 0 ) ; i != baseIds.size() ; ++i) {
45  const DetId id ( baseIds[i] );
46  if (id.subdetId() == HcalBarrel) {
47  p_hbIds->push_back( id ) ;
48  } else if (id.subdetId() == HcalEndcap) {
49  p_heIds->push_back( id ) ;
50  } else if (id.subdetId() == HcalOuter) {
51  p_hoIds->push_back( id ) ;
52  } else if (id.subdetId() == HcalForward) {
53  p_hfIds->push_back( id ) ;
54  }
55  }
56  std::sort( p_hbIds->begin(), p_hbIds->end() ) ;
57  std::sort( p_heIds->begin(), p_heIds->end() ) ;
58  std::sort( p_hoIds->begin(), p_hoIds->end() ) ;
59  std::sort( p_hfIds->begin(), p_hfIds->end() ) ;
60  p_emptyIds->resize( 0 ) ;
61 
62  m_hbIds.set(std::move(p_hbIds));
63  m_heIds.set(std::move(p_heIds));
64  m_hoIds.set(std::move(p_hoIds));
65  m_hfIds.set(std::move(p_hfIds));
66  m_emptyIds.set(std::move(p_emptyIds));
67  }
68 }
69 
70 const std::vector<DetId>&
72  int subdet ) const {
73  if( 0 != subdet && !m_hbIds.isSet() ) fillDetIds() ;
74  return ( 0 == subdet ? CaloSubdetectorGeometry::getValidDetIds() :
75  ( HcalBarrel == subdet ? *m_hbIds.load() :
76  ( HcalEndcap == subdet ? *m_heIds.load() :
77  ( HcalOuter == subdet ? *m_hoIds.load() :
78  ( HcalForward == subdet ? *m_hfIds.load() : *m_emptyIds.load() ) ) ) ) ) ;
79 }
80 
82 
83  // Now find the closest eta_bin, eta value of a bin i is average
84  // of eta[i] and eta[i-1]
85  static const double z_long=1100.0;
86  double abseta = fabs(r.eta());
87  double absz = fabs(r.z());
88 
89  // figure out subdetector, giving preference to HE in HE/HF overlap region
91  if (abseta <= m_topology.etaMax(HcalBarrel) ) {
92  bc = HcalBarrel;
93  } else if (absz >= z_long) {
94  bc = HcalForward;
95  } else if (m_topology.etaMax(HcalEndcap) ) {
96  bc = HcalEndcap;
97  } else {
98  bc = HcalForward;
99  }
100 
101  // find eta bin
102  int etaring = etaRing(bc, abseta);
103 
104  int phibin = phiBin(bc, etaring, r.phi());
105 
106  // add a sign to the etaring
107  int etabin = (r.z() > 0) ? etaring : -etaring;
108 
109  if (bc == HcalForward) {
110  static const double z_short=1137.0;
111  // Next line is premature depth 1 and 2 can coexist for large z-extent
112  // HcalDetId bestId(bc,etabin,phibin,((fabs(r.z())>=z_short)?(2):(1)));
113  // above line is no good with finite precision
114  HcalDetId bestId(bc,etabin,phibin,((fabs(r.z()) - z_short >-0.1)?(2):(1)));
115  return correctId(bestId);
116  } else {
117 
118  //Now do depth if required
119  int zside = (r.z() > 0) ? 1 : -1;
120  int dbin = m_topology.dddConstants()->getMinDepth(((int)(bc)-1),etaring,phibin,zside);
121  double pointrz(0), drz(99999.);
122  HcalDetId currentId(bc, etabin, phibin, dbin);
123  if (bc == HcalBarrel) pointrz = r.mag();
124  else pointrz = std::abs(r.z());
125  HcalDetId bestId;
126  for ( ; currentId != HcalDetId(); m_topology.incrementDepth(currentId)) {
127  const CaloCellGeometry * cell = getGeometry(currentId);
128  if (cell == 0) {
129  assert (bestId != HcalDetId());
130  break;
131  } else {
132  double rz;
133  if (bc == HcalEndcap) rz = std::abs(cell->getPosition().z());
134  else rz = cell->getPosition().mag();
135  if (std::abs(pointrz-rz)<drz) {
136  bestId = currentId;
137  drz = std::abs(pointrz-rz);
138  }
139  }
140  }
141 
142  return correctId(bestId);
143  }
144 }
145 
147  if (!m_mergePosition) {
148  return (getGeometry(id)->getPosition());
149  } else {
150  return (getGeometry(m_topology.idFront(id))->getPosition());
151  }
152 }
153 
155  if (!m_mergePosition) {
156  return (getGeometry(id)->getBackPoint());
157  } else {
158  std::vector<HcalDetId> ids;
160  return (getGeometry(ids.back())->getBackPoint());
161  }
162 }
163 
165  if (!m_mergePosition) {
166  return (getGeometry(id)->getCorners());
167  } else {
168  std::vector<HcalDetId> ids;
173  for (unsigned int k=0; k<4; ++k) {
174  mcorners[k] = mcf[k];
175  mcorners[k+4] = mcb[k+4];
176  }
177  return mcorners;
178  }
179 }
180 
181 int HcalGeometry::etaRing(HcalSubdetector bc, double abseta) const {
182  return m_topology.etaRing(bc, abseta);
183 }
184 
185 int HcalGeometry::phiBin(HcalSubdetector bc, int etaring, double phi) const {
186  return m_topology.phiBin(bc, etaring, phi);
187 }
188 
190  double dR ) const {
191  CaloSubdetectorGeometry::DetIdSet dis; // this is the return object
192 
193  if (0.000001 < dR) {
194  if (dR > M_PI/2.) {// this version needs "small" dR
195  dis = CaloSubdetectorGeometry::getCells(r, dR); // base class version
196  } else {
197  const double dR2 ( dR*dR ) ;
198  const double reta ( r.eta() ) ;
199  const double rphi ( r.phi() ) ;
200  const double lowEta ( reta - dR ) ;
201  const double highEta ( reta + dR ) ;
202  const double lowPhi ( rphi - dR ) ;
203  const double highPhi ( rphi + dR ) ;
204 
205  const double hfEtaHi (m_topology.etaMax(HcalForward));
206 
207  if (highEta > -hfEtaHi &&
208  lowEta < hfEtaHi ) { // in hcal
210 
211  for (unsigned int is ( 0 ) ; is != 4 ; ++is ) {
212  const int sign ( reta>0 ? 1 : -1 ) ;
213  const int ieta_center ( sign*etaRing( hs[is], fabs( reta ) ) ) ;
214  const int ieta_lo ( ( 0 < lowEta*sign ? sign : -sign )*etaRing( hs[is], fabs( lowEta ) ) ) ;
215  const int ieta_hi ( ( 0 < highEta*sign ? sign : -sign )*etaRing( hs[is], fabs( highEta ) ) ) ;
216  const int iphi_lo ( phiBin( hs[is], ieta_center, lowPhi ) ) ;
217  const int iphi_hi ( phiBin( hs[is], ieta_center, highPhi ) ) ;
218  const int jphi_lo ( iphi_lo>iphi_hi ? iphi_lo - 72 : iphi_lo ) ;
219  const int jphi_hi ( iphi_hi ) ;
220 
221  const int idep_lo ( 1 == is ? 4 : 1 ) ;
222  const int idep_hi ( m_topology.maxDepth(hs[is]) );
223  for (int ieta ( ieta_lo ) ; ieta <= ieta_hi ; ++ieta) {// over eta limits
224  if (ieta != 0) {
225  for (int jphi ( jphi_lo ) ; jphi <= jphi_hi ; ++jphi) { // over phi limits
226  const int iphi ( 1 > jphi ? jphi+72 : jphi ) ;
227  for (int idep ( idep_lo ) ; idep <= idep_hi ; ++idep ) {
228  const HcalDetId did ( hs[is], ieta, iphi, idep ) ;
229  if (m_topology.valid(did)) {
230  const CaloCellGeometry* cell ( getGeometry( did ) );
231  if (0 != cell ) {
232  const GlobalPoint& p ( cell->getPosition() ) ;
233  const double eta ( p.eta() ) ;
234  const double phi ( p.phi() ) ;
235  if (reco::deltaR2(eta, phi, reta, rphi ) < dR2) dis.insert( did ) ;
236  }
237  }
238  }
239  }
240  }
241  }
242  }
243  }
244  }
245  }
246  return dis;
247 }
248 
249 unsigned int HcalGeometry::getHxSize(const int type) const {
250  unsigned int hxsize(0);
251  if (!m_hbIds.isSet()) fillDetIds();
252  if (type == 1) hxsize = (m_hbIds.isSet()) ? m_hbIds->size() : 0;
253  else if (type == 2) hxsize = (m_heIds.isSet()) ? m_heIds->size() : 0;
254  else if (type == 3) hxsize = (m_hoIds.isSet()) ? m_hoIds->size() : 0;
255  else if (type == 4) hxsize = (m_hfIds.isSet()) ? m_hfIds->size() : 0;
256  else if (type == 0) hxsize = (m_emptyIds.isSet()) ? m_emptyIds->size() : 0;
257  return hxsize;
258 }
259 
261  assert( i < numberOfBarrelAlignments() ) ;
262  const int ieta ( i < numberOfBarrelAlignments()/2 ? -1 : 1 ) ;
263  const int iphi ( 1 + (4*i)%72 ) ;
264  return HcalDetId( HcalBarrel, ieta, iphi, 1 ) ;
265 }
266 
268  assert( i < numberOfEndcapAlignments() ) ;
269  const int ieta ( i < numberOfEndcapAlignments()/2 ? -16 : 16 ) ;
270  const int iphi ( 1 + (4*i)%72 ) ;
271  return HcalDetId( HcalEndcap, ieta, iphi, 1 ) ;
272 }
273 
275  assert( i < numberOfForwardAlignments() ) ;
276  const int ieta ( i < numberOfForwardAlignments()/2 ? -29 : 29 ) ;
277  const int iphi ( 1 + (4*i)%72 ) ;
278  return HcalDetId( HcalForward, ieta, iphi, 1 ) ;
279 }
280 
282  assert( i < numberOfOuterAlignments() ) ;
283  const int ring ( i/12 ) ;
284  const int ieta ( 0 == ring ? -11 :
285  1 == ring ? -5 :
286  2 == ring ? 1 :
287  3 == ring ? 5 : 11 ) ;
288  const int iphi ( 1 + ( i - ring*12 )*6 ) ;
289  return HcalDetId( HcalOuter, ieta, iphi, 4 ) ;
290 }
291 
293  assert( i < numberOfAlignments() ) ;
294 
295  const unsigned int nB ( numberOfBarrelAlignments() ) ;
296  const unsigned int nE ( numberOfEndcapAlignments() ) ;
297  const unsigned int nF ( numberOfForwardAlignments() ) ;
298 // const unsigned int nO ( numberOfOuterAlignments() ) ;
299 
300  return ( i < nB ? detIdFromBarrelAlignmentIndex( i ) :
301  i < nB+nE ? detIdFromEndcapAlignmentIndex( i - nB ) :
302  i < nB+nE+nF ? detIdFromForwardAlignmentIndex( i - nB - nE ) :
303  detIdFromOuterAlignmentIndex( i - nB - nE - nF ) ) ;
304 }
305 
307  unsigned int nD) {
308  const HcalDetId hid ( id ) ;
309  const unsigned int iphi ( hid.iphi() ) ;
310  const int ieta ( hid.ieta() ) ;
311  const unsigned int index ( ( 0 < ieta ? nD/2 : 0 ) + ( iphi + 1 )%72/4 ) ;
312  assert( index < nD ) ;
313  return index ;
314 }
315 
318 }
319 
322 }
323 
326 }
327 
329  const HcalDetId hid ( id ) ;
330  const int ieta ( hid.ieta() ) ;
331  const int iphi ( hid.iphi() ) ;
332  const int ring ( ieta < -10 ? 0 :
333  ( ieta < -4 ? 1 :
334  ( ieta < 5 ? 2 :
335  ( ieta < 11 ? 3 : 4 ) ) ) ) ;
336 
337  const unsigned int index ( 12*ring + ( iphi - 1 )/6 ) ;
338  assert( index < numberOfOuterAlignments() ) ;
339  return index ;
340 }
341 
343  assert(id.det() == DetId::Hcal) ;
344 
345  const HcalDetId hid ( id ) ;
346  bool isHB = (hid.subdet() == HcalBarrel);
347  bool isHE = (hid.subdet() == HcalEndcap);
348  bool isHF = (hid.subdet() == HcalForward);
349  // bool isHO = (hid.subdet() == HcalOuter);
350 
351  const unsigned int nB ( numberOfBarrelAlignments() ) ;
352  const unsigned int nE ( numberOfEndcapAlignments() ) ;
353  const unsigned int nF ( numberOfForwardAlignments() ) ;
354  // const unsigned int nO ( numberOfOuterAlignments() ) ;
355 
356  const unsigned int index (isHB ? alignmentBarrelIndexLocal(id) :
357  isHE ? alignmentEndcapIndexLocal(id) + nB :
358  isHF ? alignmentForwardIndexLocal( id ) + nB + nE :
359  alignmentOuterIndexLocal(id) + nB + nE + nF
360  );
361 
362  assert( index < numberOfAlignments() ) ;
363  return index ;
364 }
365 
367  return (unsigned int)DetId::Hcal - 1 ;
368 }
369 
371  const CCGFloat* pv,
372  unsigned int i,
373  Pt3D& ref) {
375 
376  if (hid.subdet() == HcalForward ) {
377  IdealZPrism::localCorners( lc, pv, ref ) ;
378  } else {
379  IdealObliquePrism::localCorners( lc, pv, ref ) ;
380  }
381 }
382 
384  const GlobalPoint& f2 ,
385  const GlobalPoint& f3 ,
386  const CCGFloat* parm ,
387  const DetId& detId) {
388 
389  assert (detId.det()==DetId::Hcal);
390 
391  const HcalDetId hid ( detId ) ;
392  unsigned int din=m_topology.detId2denseId(detId);
393 
394  edm::LogInfo("HcalGeometry") << " newCell subdet "
395  << detId.subdetId() << ", raw ID "
396  << detId.rawId() << ", hid " << hid << ", din "
397  << din << ", index ";
398 
399  if (hid.subdet()==HcalBarrel) {
400  m_hbCellVec.at( din ) = IdealObliquePrism( f1, cornersMgr(), parm ) ;
401  } else if (hid.subdet()==HcalEndcap) {
402  const unsigned int index ( din - m_hbCellVec.size() ) ;
403  m_heCellVec.at( index ) = IdealObliquePrism( f1, cornersMgr(), parm ) ;
404  } else if (hid.subdet()==HcalOuter) {
405  const unsigned int index ( din
406  - m_hbCellVec.size()
407  - m_heCellVec.size() ) ;
408  m_hoCellVec.at( index ) = IdealObliquePrism( f1, cornersMgr(), parm ) ;
409  } else {
410  const unsigned int index ( din
411  - m_hbCellVec.size()
412  - m_heCellVec.size()
413  - m_hoCellVec.size() ) ;
414  m_hfCellVec.at( index ) = IdealZPrism( f1, cornersMgr(), parm, hid.depth()==1 ? IdealZPrism::EM : IdealZPrism::HADR ) ;
415  }
416 
417  addValidID( detId ) ;
418  m_dins.push_back( din );
419 }
420 
421 const CaloCellGeometry* HcalGeometry::cellGeomPtr( unsigned int din ) const {
422  const CaloCellGeometry* cell ( 0 ) ;
423  if( m_hbCellVec.size() > din ) {
424  cell = &m_hbCellVec[ din ] ;
425  } else {
426  if (m_hbCellVec.size() + m_heCellVec.size() > din) {
427  const unsigned int index (din - m_hbCellVec.size() ) ;
428  cell = &m_heCellVec[ index ] ;
429  } else if (m_hbCellVec.size()+m_heCellVec.size()+m_hoCellVec.size() > din) {
430  const unsigned int index (din - m_hbCellVec.size() - m_heCellVec.size());
431  cell = &m_hoCellVec[ index ] ;
432  } else if (m_hbCellVec.size() + m_heCellVec.size() + m_hoCellVec.size() +
433  m_hfCellVec.size() > din) {
434  const unsigned int index (din - m_hbCellVec.size() - m_heCellVec.size() -
435  m_hoCellVec.size() ) ;
436  cell = &m_hfCellVec[ index ] ;
437  }
438  }
439 
440  return (( 0 == cell || 0 == cell->param()) ? 0 : cell ) ;
441 }
442 
446  CaloSubdetectorGeometry::IVec& dinsVec ) const {
447 
448  tVec.reserve( m_topology.ncells()*numberOfTransformParms());
449  iVec.reserve( numberOfShapes()==1 ? 1 : m_topology.ncells());
450  dVec.reserve( numberOfShapes()*numberOfParametersPerShape());
451  dinsVec.reserve( m_topology.ncells());
452 
453  for (ParVecVec::const_iterator ivv (parVecVec().begin());
454  ivv != parVecVec().end() ; ++ivv) {
455  const ParVec& pv ( *ivv ) ;
456  for (ParVec::const_iterator iv ( pv.begin() ) ; iv != pv.end() ; ++iv) {
457  dVec.push_back( *iv ) ;
458  }
459  }
460 
461  for( auto i : m_dins ) {
462  Tr3D tr ;
463  const CaloCellGeometry* ptr ( cellGeomPtr( i ) ) ;
464 
465  if (0 != ptr) {
466  dinsVec.push_back( i );
467 
468  const CCGFloat* par ( ptr->param() ) ;
469 
470  unsigned int ishape ( 9999 ) ;
471 
472  for( unsigned int ivv ( 0 ) ; ivv != parVecVec().size() ; ++ivv ) {
473  bool ok ( true ) ;
474  const CCGFloat* pv ( &(*parVecVec()[ivv].begin() ) ) ;
475  for( unsigned int k ( 0 ) ; k != numberOfParametersPerShape() ; ++k ) {
476  ok = ok && ( fabs( par[k] - pv[k] ) < 1.e-6 ) ;
477  }
478  if( ok ) {
479  ishape = ivv;
480  break ;
481  }
482  }
483  assert( 9999 != ishape ) ;
484 
485  const unsigned int nn (( numberOfShapes()==1) ? (unsigned int)1 : m_dins.size() ) ;
486  if( iVec.size() < nn ) iVec.push_back( ishape ) ;
487 
488  ptr->getTransform( tr, ( Pt3DVec* ) 0 ) ;
489 
490  if( Tr3D() == tr ) { // for preshower there is no rotation
491  const GlobalPoint& gp ( ptr->getPosition() ) ;
492  tr = HepGeom::Translate3D( gp.x(), gp.y(), gp.z() ) ;
493  }
494 
495  const CLHEP::Hep3Vector tt ( tr.getTranslation() ) ;
496  tVec.push_back( tt.x() ) ;
497  tVec.push_back( tt.y() ) ;
498  tVec.push_back( tt.z() ) ;
499  if (6 == numberOfTransformParms()) {
500  const CLHEP::HepRotation rr ( tr.getRotation() ) ;
501  const ROOT::Math::Transform3D rtr (rr.xx(), rr.xy(), rr.xz(), tt.x(),
502  rr.yx(), rr.yy(), rr.yz(), tt.y(),
503  rr.zx(), rr.zy(), rr.zz(), tt.z());
505  rtr.GetRotation( ea ) ;
506  tVec.push_back( ea.Phi() ) ;
507  tVec.push_back( ea.Theta() ) ;
508  tVec.push_back( ea.Psi() ) ;
509  }
510  }
511  }
512 }
513 
515 
516  if (m_mergePosition) {
517  HcalDetId hid(id);
518  return ((DetId)(m_topology.mergedDepthDetId(hid)));
519  } else {
520  return id;
521  }
522 }
static unsigned int numberOfBarrelAlignments()
Definition: HcalGeometry.h:59
type
Definition: HCALResponse.h:21
unsigned int getHFSize() const
Definition: HcalTopology.h:135
virtual void newCell(const GlobalPoint &f1, const GlobalPoint &f2, const GlobalPoint &f3, const CCGFloat *parm, const DetId &detId)
int i
Definition: DBlmapReader.cc:9
virtual DetId getClosestCell(const GlobalPoint &r) const
Definition: HcalGeometry.cc:81
static unsigned int alignmentTransformIndexLocal(const DetId &id)
const HcalDDDRecConstants * dddConstants() const
Definition: HcalTopology.h:161
virtual unsigned int detId2denseId(const DetId &id) const
return a linear packed id
std::vector< CCGFloat > DimVec
unsigned int getHOSize() const
Definition: HcalTopology.h:134
CaloCellGeometry::Pt3DVec Pt3DVec
Definition: HcalGeometry.cc:11
edm::AtomicPtrCache< std::vector< DetId > > m_heIds
Definition: HcalGeometry.h:135
static unsigned int alignmentEndcapIndexLocal(const DetId &id)
HcalSubdetector subdet() const
get the subdetector
Definition: HcalDetId.h:49
int phiBin(HcalSubdetector subdet, int etaRing, double phi) const
virtual void getSummary(CaloSubdetectorGeometry::TrVec &trVector, CaloSubdetectorGeometry::IVec &iVector, CaloSubdetectorGeometry::DimVec &dimVector, CaloSubdetectorGeometry::IVec &dinsVector) const
CaloTopology const * topology(0)
unsigned int getHxSize(const int type) const
static DetId detIdFromOuterAlignmentIndex(unsigned int i)
Geom::Phi< T > phi() const
Definition: PV3DBase.h:69
const HcalTopology & m_topology
Definition: HcalGeometry.h:131
std::vector< unsigned int > IVec
edm::AtomicPtrCache< std::vector< DetId > > m_emptyIds
Definition: HcalGeometry.h:138
bool isHE(int etabin, int depth)
std::vector< CCGFloat > TrVec
const_iterator begin() const
Definition: EZArrayFL.h:63
virtual const CaloCellGeometry * getGeometry(const DetId &id) const
Get the cell geometry of a given detector id. Should return false if not found.
Definition: HcalGeometry.h:102
HepGeom::Transform3D Tr3D
std::vector< IdealObliquePrism > HECellVec
Definition: HcalGeometry.h:19
std::vector< Pt3D > Pt3DVec
bool isHB(int etabin, int depth)
CaloCellGeometry::CCGFloat CCGFloat
bool m_mergePosition
Definition: HcalGeometry.h:132
bool isSet() const
virtual unsigned int numberOfParametersPerShape() const
Definition: HcalGeometry.h:38
virtual const std::vector< DetId > & getValidDetIds(DetId::Detector det=DetId::Detector(0), int subdet=0) const
Get a list of valid detector ids (for the given subdetector)
Definition: HcalGeometry.cc:71
virtual CaloSubdetectorGeometry::DetIdSet getCells(const GlobalPoint &r, double dR) const
Get a list of all cells within a dR of the given cell.
virtual const std::vector< DetId > & getValidDetIds(DetId::Detector det=DetId::Detector(0), int subdet=0) const
Get a list of valid detector ids (for the given subdetector)
HcalDetId mergedDepthDetId(const HcalDetId &id) const
Definition: HcalTopology.h:163
CaloCellGeometry::Tr3D Tr3D
Definition: HcalGeometry.cc:12
uint32_t rawId() const
get the raw id
Definition: DetId.h:43
CaloCellGeometry::Pt3DVec Pt3DVec
Definition: HcalGeometry.h:25
virtual DetIdSet getCells(const GlobalPoint &r, double dR) const
Get a list of all cells within a dR of the given cell.
const CCGFloat * param() const
T mag() const
Definition: PV3DBase.h:67
HBCellVec m_hbCellVec
Definition: HcalGeometry.h:141
HOCellVec m_hoCellVec
Definition: HcalGeometry.h:143
static void localCorners(Pt3DVec &vec, const CCGFloat *pv, Pt3D &ref)
int getMinDepth(const int itype, const int ieta, const int iphi, const int zside) const
HFCellVec m_hfCellVec
Definition: HcalGeometry.h:144
virtual void getTransform(Tr3D &tr, Pt3DVec *lptr) const
--------— only needed by specific utility; overloaded when needed -—
std::vector< IdealObliquePrism > HOCellVec
Definition: HcalGeometry.h:20
static unsigned int numberOfOuterAlignments()
Definition: HcalGeometry.h:65
int maxDepth(HcalSubdetector subdet) const
int etaRing(HcalSubdetector subdet, double eta) const
eta and phi index from eta, phi values
DetId correctId(const DetId &id) const
T z() const
Definition: PV3DBase.h:64
int ieta() const
get the cell ieta
Definition: HcalDetId.h:56
static DetId detIdFromForwardAlignmentIndex(unsigned int i)
static unsigned int alignmentTransformIndexGlobal(const DetId &id)
static unsigned int alignmentBarrelIndexLocal(const DetId &id)
HcalSubdetector
Definition: HcalAssistant.h:31
def pv(vc)
Definition: MetAnalyzer.py:6
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
static DetId detIdFromBarrelAlignmentIndex(unsigned int i)
void unmergeDepthDetId(const HcalDetId &id, std::vector< HcalDetId > &ids) const
Definition: HcalTopology.h:165
CaloCellGeometry::CCGFloat CCGFloat
Definition: HcalGeometry.cc:9
virtual DetId denseId2detId(unsigned int) const
return a linear packed id
int phiBin(HcalSubdetector bc, int etaring, double phi) const
static DetId detIdFromLocalAlignmentIndex(unsigned int i)
int etaRing(HcalSubdetector bc, double abseta) const
helper methods for getClosestCell
CaloSubdetectorGeometry::IVec m_dins
Definition: HcalGeometry.h:139
HcalGeometry(const HcalTopology &topology)
Definition: HcalGeometry.cc:14
int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:37
GlobalPoint getPosition(const DetId &id) const
void localCorners(Pt3DVec &lc, const CCGFloat *pv, unsigned int i, Pt3D &ref)
#define M_PI
int k[5][pyjets_maxn]
static unsigned int alignmentOuterIndexLocal(const DetId &id)
std::vector< IdealZPrism > HFCellVec
Definition: HcalGeometry.h:21
int iphi() const
get the cell iphi
Definition: HcalDetId.cc:103
Definition: DetId.h:18
bool incrementDepth(HcalDetId &id) const
static unsigned int numberOfEndcapAlignments()
Definition: HcalGeometry.h:61
CaloCellGeometry::Pt3D Pt3D
Definition: HcalGeometry.cc:10
AlgebraicVector EulerAngles
Definition: Definitions.h:36
void addValidID(const DetId &id)
bool set(std::unique_ptr< T > iNewValue) const
virtual ~HcalGeometry()
The HcalGeometry will delete all its cell geometries at destruction time.
Definition: HcalGeometry.cc:19
CaloCellGeometry::CornersMgr * cornersMgr()
Detector
Definition: DetId.h:24
static void localCorners(Pt3DVec &vec, const CCGFloat *pv, Pt3D &ref)
Definition: IdealZPrism.cc:122
CaloCellGeometry::CornersVec getCorners(const DetId &id) const
HepGeom::Point3D< CCGFloat > Pt3D
bool isHF(int etabin, int depth)
virtual bool valid(const DetId &id) const
std::vector< IdealObliquePrism > HBCellVec
Definition: HcalGeometry.h:18
virtual const CaloCellGeometry * cellGeomPtr(unsigned int index) const
unsigned int getHESize() const
Definition: HcalTopology.h:133
T1 deltaR2(T1 eta1, T2 phi1, T3 eta2, T4 phi2)
Definition: deltaR.h:36
HECellVec m_heCellVec
Definition: HcalGeometry.h:142
CaloCellGeometry::Tr3D Tr3D
T eta() const
Definition: PV3DBase.h:76
MgrType::const_iterator const_iterator
Definition: EZArrayFL.h:27
virtual unsigned int numberOfShapes() const
Definition: HcalGeometry.h:37
edm::AtomicPtrCache< std::vector< DetId > > m_hfIds
Definition: HcalGeometry.h:137
#define begin
Definition: vmac.h:30
const_iterator end() const
Definition: EZArrayFL.h:64
void fillDetIds() const
Definition: HcalGeometry.cc:35
static unsigned int numberOfForwardAlignments()
Definition: HcalGeometry.h:63
static unsigned int alignmentForwardIndexLocal(const DetId &id)
GlobalPoint getBackPosition(const DetId &id) const
HcalDetId idFront(const HcalDetId &id) const
Definition: HcalTopology.h:170
static DetId detIdFromEndcapAlignmentIndex(unsigned int i)
static unsigned int alignmentBarEndForIndexLocal(const DetId &id, unsigned int nD)
virtual unsigned int numberOfTransformParms() const
edm::AtomicPtrCache< std::vector< DetId > > m_hoIds
Definition: HcalGeometry.h:136
CaloCellGeometry::Pt3D Pt3D
Definition: HcalGeometry.h:24
bool withSpecialRBXHBHE() const
Definition: HcalTopology.h:162
Detector det() const
get the detector field from this detid
Definition: DetId.h:35
unsigned int getHBSize() const
Definition: HcalTopology.h:132
const GlobalPoint & getPosition() const
Returns the position of reference for this cell.
T const * load() const
edm::AtomicPtrCache< std::vector< DetId > > m_hbIds
Definition: HcalGeometry.h:134
static unsigned int numberOfAlignments()
Definition: HcalGeometry.h:69
virtual unsigned int ncells() const
return a count of valid cells (for dense indexing use)
def move(src, dest)
Definition: eostools.py:510
double etaMax(HcalSubdetector subdet) const