11 #include "CLHEP/Units/GlobalPhysicalConstants.h"
21 firstHERing_(999), lastHERing_(0),
22 firstHFRing_(29), lastHFRing_(41),
23 firstHORing_(1), lastHORing_(15),
24 firstHEDoublePhiRing_(999), firstHEQuadPhiRing_(999),
25 firstHFQuadPhiRing_(40), firstHETripleDepthRing_(999),
101 const double fiveDegInRad = 2*
M_PI/72;
115 std::cout <<
"Set segmentation for ring " <<
ring <<
" with "
116 << segmentation.size() <<
" elements:";
117 for (
unsigned int k=0;
k<segmentation.size(); ++
k)
138 mode_(mode), triggerMode_(tmode),
139 firstHBRing_(1), lastHBRing_(16),
140 firstHERing_(16), lastHERing_(29),
141 firstHFRing_(29), lastHFRing_(41),
142 firstHORing_(1), lastHORing_(15),
143 firstHEDoublePhiRing_((mode==HcalTopologyMode::
H2 || mode==HcalTopologyMode::
H2HE)?(22):(21)),
144 firstHEQuadPhiRing_(999), firstHFQuadPhiRing_(40),
145 firstHETripleDepthRing_((mode==HcalTopologyMode::
H2 || mode==HcalTopologyMode::
H2HE)?(24):(27)),
146 singlePhiBins_(
IPHI_MAX), doublePhiBins_(36),
147 maxDepthHB_(maxDepthHB), maxDepthHE_(maxDepthHE), maxDepthHF_(2),
148 etaHE2HF_(30), etaHF2HE_(29), maxPhiHE_(
IPHI_MAX),
149 HBSize_(kHBSizePreLS1),
150 HESize_(kHESizePreLS1),
153 HTSize_(kHTSizePreLS1),
154 numberOfShapes_(( mode==HcalTopologyMode::
SLHC ) ? 500 : 87 ) {
177 edm::LogWarning(
"CaloTopology") <<
"This is an incomplete constructor of HcalTopology - be warned that many functionalities will not be there - revert from this - get from EventSetup";
202 if (
id.iphi()<1 ||
id.iphi()>
IPHI_MAX ||
id.ieta()==0)
return false;
203 if (
id.
depth() != 0)
return false;
206 (
id.ietaAbs()>32))
return false;
208 if (
id.ietaAbs()>ietaMax && ((
id.iphi()%4)!=1))
return false;
211 if (
id.ietaAbs()<30 ||
id.ietaAbs()>41)
return false;
212 if (
id.ietaAbs()>29 && ((
id.iphi()%2)==0))
return false;
213 if (
id.ietaAbs()>39 && ((
id.iphi()%4)!=3))
return false;
221 switch (
id.subdet()) {
254 std::vector<DetId> vNeighborsDetId;
257 if (neighbors[
i].oldFormat()) neighbors[
i].changeForm();
258 vNeighborsDetId.push_back(
DetId(neighbors[
i].rawId()));
260 return vNeighborsDetId;
264 std::vector<DetId> vNeighborsDetId;
267 if (neighbors[
i].oldFormat()) neighbors[
i].changeForm();
268 vNeighborsDetId.push_back(
DetId(neighbors[
i].rawId()));
270 return vNeighborsDetId;
274 std::vector<DetId> vNeighborsDetId;
278 vNeighborsDetId.push_back(
DetId(neighbor.
rawId()));
280 return vNeighborsDetId;
284 std::vector<DetId> vNeighborsDetId;
288 vNeighborsDetId.push_back(
DetId(neighbor.
rawId()));
290 return vNeighborsDetId;
295 std::vector<DetId> vNeighborsDetId;
298 vNeighborsDetId.push_back(neighbor);
300 return vNeighborsDetId;
305 std::vector<DetId> vNeighborsDetId;
308 vNeighborsDetId.push_back(neighbor);
310 return vNeighborsDetId;
330 int depth_l=
std::min(depth1,depth2);
331 int depth_h=
std::max(depth1,depth2);
334 for (
int ieta=ieta_l; ieta<=ieta_h; ieta++)
335 for (
int iphi=iphi_l; iphi<=iphi_h; iphi++)
374 const int ie (
id.ietaAbs());
375 const int ip (
id.iphi());
378 return ( ( ip >= 1 ) &&
385 ( ( ( ie == 15 ) || ( ie == 16 ) ) &&
405 ( ip%2 == 1 ) ) ) ) ||
416 ( ip%4 == 3 ) ) ) ) ) ) ;
423 int aieta=
id.ietaAbs();
424 int depth=
id.depth();
428 if ((ieta==0 || iphi<=0 || iphi>maxPhi) || aieta>
maxEta_) ok =
false;
434 (depth<hcons_->getMinDepth(0,aieta))) ok=
false;
436 if (aieta>
lastHBRing() || depth>2 || (aieta<=14 && depth>1)) ok=
false;
442 (depth<hcons_->getMinDepth(1,aieta))) {
450 if (depth < 1) ok =
false;
484 switch (
id.subdet()) {
488 else neighbor=
HcalDetId(
id.subdet(),
id.ieta(),
id.iphi()+1,
id.
depth());
493 else neighbor=
HcalDetId(
id.subdet(),
id.ieta(),
id.iphi()+4,
id.
depth());
496 else neighbor=
HcalDetId(
id.subdet(),
id.ieta(),
id.iphi()+2,
id.
depth());
499 else neighbor=
HcalDetId(
id.subdet(),
id.ieta(),
id.iphi()+1,
id.
depth());
505 else neighbor=
HcalDetId(
id.subdet(),
id.ieta(),
id.iphi()+4,
id.
depth());
508 else neighbor=
HcalDetId(
id.subdet(),
id.ieta(),
id.iphi()+2,
id.
depth());
510 if (!
validRaw(neighbor)) ok =
false;
522 switch (
id.subdet()) {
526 else neighbor=
HcalDetId(
id.subdet(),
id.ieta(),
id.iphi()-1,
id.
depth());
531 else neighbor=
HcalDetId(
id.subdet(),
id.ieta(),
id.iphi()-4,
id.
depth());
534 else neighbor=
HcalDetId(
id.subdet(),
id.ieta(),
id.iphi()-2,
id.
depth());
537 else neighbor=
HcalDetId(
id.subdet(),
id.ieta(),
id.iphi()-1,
id.
depth());
543 else neighbor=
HcalDetId(
id.subdet(),
id.ieta(),
id.iphi()-4,
id.
depth());
546 else neighbor=
HcalDetId(
id.subdet(),
id.ieta(),
id.iphi()-2,
id.
depth());
548 if (!
validRaw(neighbor)) ok =
false;
569 int aieta=
id.ietaAbs();
574 neighbors[0]=
HcalDetId(
id.subdet(),(aieta+1)*
id.
zside(),((
id.iphi()==1)?(71):(
id.iphi()-2)),
id.
depth());
576 neighbors[0]=
HcalDetId(
id.subdet(),(aieta+1)*
id.
zside(),((
id.iphi()==1)?(71):(
id.iphi()-2)),
id.
depth());
584 if (!
valid(neighbors[0])) n=0;
591 int aieta=
id.ietaAbs();
607 }
else if (aieta==1) {
616 if (!
valid(neighbors[0]) && n==2) {
617 if (!
valid(neighbors[1])) n=0;
620 neighbors[0]=neighbors[1];
623 if (n==2 && !
valid(neighbors[1])) n=1;
624 if (n==1 && !
valid(neighbors[0])) n=0;
631 int & nDepthBins,
int & startingBin)
const {
662 }
else if (etaRing==17) {
680 std::cerr <<
"Bad HCAL subdetector " << subdet << std::endl;
687 int ieta = detId.
ieta();
690 int nDepthBins, startingBin;
695 if (depth > nDepthBins) {
708 (ieta > 0) ? ++ieta : --ieta;
713 (ieta > 0) ? --ieta : ++ieta;
726 int ieta = detId.
ieta();
729 int nDepthBins, startingBin;
745 (ieta > 0) ? --ieta : ++ieta;
746 }
else if (depth <= 0) {
751 ieta = (ieta > 0) ? etaRing : -etaRing;
782 static const double twopi =
M_PI+
M_PI;
822 static const double twopi =
M_PI+
M_PI;
838 if (phi<0.0) phi += twopi;
839 if (phi>twopi) phi -= twopi;
840 int phibin(1),
unit(1);
853 if (unit == 2) iphi = (phibin-1)*2+1;
854 else if (unit == 4) iphi = (phibin-1)*4+3;
863 throw cms::Exception(
"HcalTopology") <<
"No depth segmentation found for ring" <<
ring;
867 readoutDepths = pos->second;
875 std::vector<int> readoutDepths;
877 int d1 = std::lower_bound(readoutDepths.begin(), readoutDepths.end(),
depth) - readoutDepths.begin();
878 int d2 = std::upper_bound(readoutDepths.begin(), readoutDepths.end(),
depth) - readoutDepths.begin();
879 return std::pair<int, int>(d1, d2);
903 int ieta = (keta > 0) ? keta : -keta;
912 if ((ietal < (
int)(
etaTable.size())) && (ieta > 0))
913 return std::pair<double,double>(
etaTable[ieta-1],etaTable[ietal]);
915 return std::pair<double,double>(0,0);
922 const int ip (hid.
iphi() ) ;
923 const int ie (hid.
ietaAbs() ) ;
925 const int zn (hid.
zside() < 0 ? 1 : 0 ) ;
927 ( ip - 1 )*18 +
dp - 1 + ie - ( ie<16 ? 1 : 0 ) + zn*
kHBhalf :
929 2*
kHBhalf + ( ip - 1 )*8 + ( ip/2 )*20 +
930 ( ( ie==16 || ie==17 ) ? ie - 16 :
931 ( ( ie>=18 && ie<=20 ) ? 2 + 2*( ie - 18 ) +
dp - 1 :
932 ( ( ie>=21 && ie<=26 ) ? 8 + 2*( ie - 21 ) +
dp - 1 :
933 ( ( ie>=27 && ie<=28 ) ? 20 + 3*( ie - 27 ) +
dp - 1 :
934 26 + 2*( ie - 29 ) +
dp - 1 ) ) ) ) + zn*
kHEhalf :
939 ( ( ip - 1 )/4 )*4 + ( ( ip - 1 )/2 )*22 +
940 2*( ie - 29 ) + (
dp - 1 ) + zn*
kHFhalf : 0xFFFFFFFFu ) ) ) ) ;
947 const int ip (hid.
iphi() ) ;
948 const int ie (hid.
ietaAbs() ) ;
950 const int zn (hid.
zside() < 0 ? 1 : 0 ) ;
951 unsigned int retval = 0xFFFFFFFFu;
953 retval=( ip - 1 )*18 +
dp - 1 + ie - ( ie<16 ? 1 : 0 ) + zn*
kHBhalf;
964 const int ip (hid.
iphi() ) ;
965 const int ie (hid.
ietaAbs() ) ;
967 const int zn (hid.
zside() < 0 ? 1 : 0 ) ;
968 unsigned int retval = 0xFFFFFFFFu;
970 retval=( ip - 1 )*8 + ( ip/2 )*20 +
971 ( ( ie==16 || ie==17 ) ? ie - 16 :
972 ( ( ie>=18 && ie<=20 ) ? 2 + 2*( ie - 18 ) +
dp - 1 :
973 ( ( ie>=21 && ie<=26 ) ? 8 + 2*( ie - 21 ) +
dp - 1 :
974 ( ( ie>=27 && ie<=28 ) ? 20 + 3*( ie - 27 ) +
dp - 1 :
975 26 + 2*( ie - 29 ) +
dp - 1 ) ) ) ) + zn*
kHEhalf;
986 const int ip (hid.
iphi() ) ;
987 const int ie (hid.
ietaAbs() ) ;
988 const int zn (hid.
zside() < 0 ? 1 : 0 ) ;
990 unsigned int retval = 0xFFFFFFFFu;
992 retval=( ip - 1 )*15 + ( ie - 1 ) + zn*
kHOhalf;
1002 const int ip (hid.
iphi() ) ;
1003 const int ie (hid.
ietaAbs() ) ;
1005 const int zn (hid.
zside() < 0 ? 1 : 0 ) ;
1007 unsigned int retval = 0xFFFFFFFFu;
1009 retval = ( ( ip - 1 )/4 )*4 + ( ( ip - 1 )/2 )*22 +
1010 2*( ie - 29 ) + (
dp - 1 ) + zn*
kHFhalf;
1012 retval=
dp-1+2*(ip-1);
1022 unsigned int ietaAbs = tid.
ietaAbs();
1023 unsigned int iphi = tid.
iphi();
1024 unsigned int ivers = tid.
version();
1028 if ((iphi-1)%4==0) index = (iphi-1)*32 + (ietaAbs-1) - (12*((iphi-1)/4));
1029 else index = (iphi-1)*28 + (ietaAbs-1) + (4*(((iphi-1)/4)+1));
1030 if (zside == -1) index +=
kHThalf;
1034 index += (36 * (ietaAbs - 30) + ((iphi - 1)/2));
1043 int ieta = tid.
ieta();
1044 int iphi = tid.
iphi();
1046 unsigned int index=0xFFFFFFFFu;
1056 index = ((iphi+1)/4-1) + 18*channel + 27*(ieta+1);
1061 if (channel>2) channel-=1;
1062 index = ((iphi+1)/4-1) + 18*channel + 54*(ieta+1) + 108;
1066 if (channel==8) channel = 2;
1068 index = (iphi-1)/18 + 4*channel + 6*(ieta+1) + 324;
1073 index = (ieta+2) + 420;
1077 if (ieta<0) index = ((iphi+1)/12-1) + 36*channel + 6*(ieta+2) + 348;
1078 else if (ieta>0) index = ((iphi+1)/12-1) + 36*channel + 6*(ieta+2) + 6 + 348;
1079 else index = ((iphi+1)/6-1) + 36*channel + 6*(ieta+2) + 348;
1090 if (
abs(ieta)==4) index = ((iphi-1)%36) + (((zside+1)*36)/2) + 72 + 425;
1091 else index = (iphi-1) + (36*(zside+1)*2) + 425;
1099 unsigned int retval(0);
1129 << std::hex << retval <<
std::dec << std::endl;
1140 int in ( denseid ) ;
1147 iz = ( in<
kHFhalf ? 1 : -1 ) ;
1151 ip += 1 + ( in>21 ? 2 : 0 ) ;
1152 if( 3 == ip%4 ) in -= 22 ;
1158 iz = ( in<
kHOhalf ? 1 : -1 ) ;
1162 ie = 1 + ( in - 15*( ip - 1 ) ) ;
1163 }
else if ( in > 2*
kHBhalf - 1 ) {
1166 iz = ( in<
kHEhalf ? 1 : -1 ) ;
1171 if( 0 == ip%2 ) in %= 28 ;
1172 ie = 15 + ( in<2 ? 1 + in : 2 +
1173 ( in<20 ? 1 + ( in - 2 )/2 : 9 +
1174 ( in<26 ? 1 + ( in - 20 )/3 : 3 ) ) ) ;
1177 ( in<20 ? 1 + ( in - 2 )%2 :
1178 ( in<26 ? 1 + ( in - 20 )%3 :
1179 ( 1 + ( in - 26 )%2 ) ) ) ) ) ;
1181 iz = ( in<
kHBhalf ? 1 : -1 ) ;
1196 if (denseid <
ncells()) {
1204 if (ie > 12) {ie = 54 -ie; iz = -1;}
1205 else {ie += 29; iz = 1;}
1212 if (ie > 14) {ie = 30 -ie; iz = -1;}
1213 else {ie += 1; iz = 1;}
1214 }
else if (denseid >= (
HBSize_)) {
1237 <<
" : " << hid.
rawId() <<
std::dec <<
" | " << hid << std::endl;
int getNPhi(const int type) const
std::vector< int > unitPhiHF
int zside() const
get the z-side of the tower (1/-1)
int decIEta(const HcalDetId &id, HcalDetId neighbors[2]) const
unsigned int numberOfShapes_
unsigned int detId2denseIdPreLS1(const DetId &id) const
virtual unsigned int detId2denseId(const DetId &id) const
return a linear packed id
int getDepthEta16(int i) const
std::vector< double > dPhiTableHF
void excludeSubdetector(HcalSubdetector subdet)
HcalSubdetector subdet() const
get the subdetector
std::vector< HcalDDDRecConstants::HcalEtaBin > etaBinsHB_
int phiBin(HcalSubdetector subdet, int etaRing, double phi) const
std::vector< int > unitPhi
CalibDetType calibFlavor() const
get the flavor of this calibration detid
virtual std::vector< DetId > down(const DetId &id) const
int zside() const
get the z-side of the cell (1/-1)
bool validDetId(HcalSubdetector subdet, int ieta, int iphi, int depth) const
virtual std::vector< DetId > south(const DetId &id) const
unsigned int detId2denseIdHT(const DetId &id) const
return a linear packed id from HT
unsigned int detId2denseIdHF(const DetId &id) const
return a linear packed id from HF
void exclude(const HcalDetId &id)
const std::vector< double > & getEtaTableHF() const
HcalTopologyMode::TriggerMode triggerMode_
int nPhiBins(int etaRing) const
how many phi segments in this ring
int incIEta(const HcalDetId &id, HcalDetId neighbors[2]) const
bool decrementDepth(HcalDetId &id) const
HcalTopology(const HcalDDDRecConstants *hcons)
unsigned int detId2denseIdHB(const DetId &id) const
return a linear packed id from HB
bool decIPhi(const HcalDetId &id, HcalDetId &neighbor) const
bool validHcal(const HcalDetId &id) const
bool neighbor(int endcap, int sector, int SectIndex, int id, int sub, int station)
std::vector< HcalDDDRecConstants::HcalEtaBin > etaBinsHE_
int ieta() const
get the rbx name (if relevant)
int getMinDepth(const int itype, const int ieta) const
const std::vector< double > & getPhiOffs() const
virtual int topoVersion() const
return a version which identifies the given topology
int getNoff(const int i) const
int getTriggerMode() const
uint32_t rawId() const
get the raw id
int getMaxDepth(const int type) const
bool isExcluded(const HcalDetId &id) const
int depth() const
get the tower depth
bool incIPhi(const HcalDetId &id, HcalDetId &neighbor) const
bool validRaw(const HcalDetId &id) const
int decAIEta(const HcalDetId &id, HcalDetId neighbors[2]) const
int maxDepth(HcalSubdetector subdet) const
HcalTopologyMode::Mode mode_
int etaRing(HcalSubdetector subdet, double eta) const
eta and phi index from eta, phi values
void setDepthSegmentation(unsigned ring, const std::vector< int > &readoutDepths)
int firstHETripleDepthRing() const
int ieta() const
get the cell ieta
Abs< T >::type abs(const T &t)
virtual std::vector< DetId > west(const DetId &id) const
void depthBinInformation(HcalSubdetector subdet, int etaRing, int &nDepthBins, int &startingBin) const
finds the number of depth bins and which is the number to start with
virtual DetId denseId2detId(unsigned int) const
return a linear packed id
std::vector< double > etaTableHF
int firstHEDoublePhiRing() const
SegmentationMap depthSegmentation_
std::vector< double > etaTable
int ietaAbs() const
get the absolute value of the cell ieta
int iphi() const
get the low-edge iphi (if relevant)
bool validDetIdPreLS1(const HcalDetId &id) const
int iphi() const
get the cell iphi
bool incrementDepth(HcalDetId &id) const
std::pair< int, int > segmentBoundaries(unsigned ring, unsigned depth) const
int zside() const
get the sign of ieta (+/-1)
std::vector< double > phioff
const std::vector< double > & getPhiTable() const
void getDepthSegmentation(unsigned ring, std::vector< int > &readoutDepths) const
virtual std::vector< DetId > north(const DetId &id) const
int version() const
get the version code for the trigger tower
virtual bool valid(const DetId &id) const
const std::vector< double > & getEtaTable() const
std::pair< int, int > getEtaRange(const int i) const
int firstHFQuadPhiRing() const
Geom::Phi< T > phi() const
virtual std::vector< DetId > up(const DetId &id) const
const std::vector< int > & getDepth(const unsigned int i) const
std::vector< double > dPhiTable
int cboxChannel() const
get the calibration box channel (if relevant)
std::pair< double, double > etaRange(HcalSubdetector subdet, int ieta) const
TString units(TString variable, Char_t axis)
const std::vector< double > & getPhiTableHF() const
int firstHETripleDepthRing_
unsigned int detId2denseIdHO(const DetId &id) const
return a linear packed id from HO
std::vector< HcalDetId > exclusionList_
unsigned int detId2denseIdHE(const DetId &id) const
return a linear packed id from HE
int firstHEQuadPhiRing() const
HcalSubdetector hcalSubdet() const
get the HcalSubdetector (if relevant)
virtual std::vector< DetId > east(const DetId &id) const
volatile std::atomic< bool > shutdown_flag false
int maxHFDepth(int ieta, int iphi) const
int incAIEta(const HcalDetId &id, HcalDetId neighbors[2]) const
unsigned int detId2denseIdCALIB(const DetId &id) const
return a linear packed id from CALIB
int ietaAbs() const
get the absolute value of the tower ieta
int iphi() const
get the tower iphi
virtual unsigned int ncells() const
return a count of valid cells (for dense indexing use)
double etaMax(HcalSubdetector subdet) const
int firstHEDoublePhiRing_
const HcalDDDRecConstants * hcons_
static const int IPHI_MAX
std::vector< HcalEtaBin > getEtaBins(const int itype) const
bool validHT(const HcalTrigTowerDetId &id) const