Skip to content

Commit

Permalink
Clean up Hcal rec constants and use Hcal parameters
Browse files Browse the repository at this point in the history
  • Loading branch information
@ianna
ianna committed Jul 15, 2015
1 parent ef71279 commit ec2c9c1
Show file tree
Hide file tree
Showing 2 changed files with 35 additions and 110 deletions.
25 changes: 7 additions & 18 deletions Geometry/HcalCommonData/interface/HcalDDDRecConstants.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -56,23 +56,23 @@ class HcalDDDRecConstants {
std::vector<HcalEtaBin> getEtaBins(const int itype) const;
std::pair<double,double> getEtaPhi(int subdet, int ieta, int iphi) const;
std::pair<int,int> getEtaRange(const int i) const
{return std::pair<int,int>(iEtaMin[i],iEtaMax[i]);}
const std::vector<double> & getEtaTable() const {return etaTable;}
{return std::pair<int,int>(hpar->etaMin[i],hpar->etaMax[i]);}
const std::vector<double> & getEtaTable() const {return hpar->etaTable;}
const std::vector<double> & getEtaTableHF() const {return hpar->etaTableHF;}
std::pair<double,double> getEtaLimit(const int i) const
{return std::pair<double,double>(etaTable[i],etaTable[i+1]);}
{return std::pair<double,double>(hpar->etaTable[i],hpar->etaTable[i+1]);}
HcalID getHCID(int subdet, int ieta, int iphi, int lay,
int idepth) const;
int getMaxDepth(const int type) const {return maxDepth[type];}
int getNEta() const {return nEta;}
int getMaxDepth(const int type) const {return hpar->maxDepth[type];}
int getNEta() const {return hpar->etagroup.size();}
double getPhiBin(const int i) const {return phibin[i];}
double getPhiOff(const int i) const {return hpar->phioff[i];}
const std::vector<double> & getPhiOffs() const {return hpar->phioff;}
const std::vector<double> & getPhiTable() const {return phibin;}
const std::vector<double> & getPhiTableHF() const {return phibinHF;}
const std::vector<double> & getPhiTableHF() const {return hpar->phitable;}
double getRZ(int subdet, int ieta, int depth) const;
std::vector<HcalActiveLength> getThickActive(const int type) const;
int getTopoMode() const {return modeTopo_;}
int getTopoMode() const {return hpar->topologyMode;}
std::vector<HcalCellType> HcalCellTypes(HcalSubdetector) const;
unsigned int numberOfCells(HcalSubdetector) const;
unsigned int nCells(HcalSubdetector) const;
Expand All @@ -83,26 +83,15 @@ class HcalDDDRecConstants {
unsigned int layerGroupSize( unsigned int eta ) const;
unsigned int layerGroup( unsigned int eta, unsigned int i ) const;

bool tobeInitialized;
const HcalParameters *hpar;
const HcalDDDSimConstants &hcons;
int modeTopo_; // Mode for topology
std::vector<int> etaGroup; // Eta Grouping
std::vector<std::pair<int,int> > etaSimValu; // eta ranges at Sim stage
std::vector<double> etaTable; // Eta table (HB+HE)
std::vector<int> ietaMap; // Map Sim level ieta to Rec level ieta
std::vector<int> iEtaMin, iEtaMax; // Minimum and maximum eta
std::vector<int> maxDepth; // Maximum depth in HB/HE/HF/HO
int nEta; // Number of bins in eta for HB and HE
std::vector<int> phiGroup; // Eta Grouping
std::vector<double> phibin; // Phi step for all eta bins (HB, HE, HO)
std::vector<double> phibinHF; // Phi step for all eta bins (HF)
std::vector<int> phiUnitS; // Phi unit at SIM stage
std::vector<int> nOff; // Speical eta bin #'s in barrel and endcap
std::vector<std::pair<double,double> > gconsHB; // Geometry constatnts HB
std::vector<std::pair<double,double> > gconsHE; // Geometry constatnts HE
int nModule[2], nHalves[2]; // Modules, Halves for HB/HE
enum { kHOSizePreLS1 = 2160, kHFSizePreLS1 = 1728 } ;
};

#endif
120 changes: 28 additions & 92 deletions Geometry/HcalCommonData/src/HcalDDDRecConstants.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,7 @@
#include "CLHEP/Units/GlobalSystemOfUnits.h"

//#define DebugLog
enum { kHOSizePreLS1 = 2160, kHFSizePreLS1 = 1728 } ;

HcalDDDRecConstants::HcalDDDRecConstants(const HcalParameters* hp,
const HcalDDDSimConstants& hc) :
Expand All @@ -30,14 +31,14 @@ HcalDDDRecConstants::getEtaBins(const int itype) const {
std::vector<HcalDDDRecConstants::HcalEtaBin> bins;
unsigned int type = (itype == 0) ? 0 : 1;
unsigned int lymax = (type == 0) ? 17 : 19;
for (int ieta = iEtaMin[type]; ieta <= iEtaMax[type]; ++ieta) {
for (int ieta = hpar->etaMin[type]; ieta <= hpar->etaMax[type]; ++ieta) {
int nfi = (int)((20.001*nModule[itype]*CLHEP::deg)/phibin[ieta-1]);
HcalDDDRecConstants::HcalEtaBin etabin = HcalDDDRecConstants::HcalEtaBin(ieta,etaTable[ieta-1], etaTable[ieta], nfi, hpar->phioff[type], phibin[ieta-1]);
HcalDDDRecConstants::HcalEtaBin etabin = HcalDDDRecConstants::HcalEtaBin(ieta,hpar->etaTable[ieta-1], hpar->etaTable[ieta], nfi, hpar->phioff[type], phibin[ieta-1]);
int dstart = -1;
if (layerGroupSize( ieta-1 ) > 0) {
int lmin(0), lmax(0);
int dep = layerGroup( ieta-1, 0 );
if (type == 1 && ieta == iEtaMin[type]) dep = 3;
if (type == 1 && ieta == hpar->etaMin[type]) dep = 3;
unsigned lymx0 = (layerGroupSize( ieta-1 ) > lymax) ? lymax : layerGroupSize( ieta-1 );
for (unsigned int l=0; l<lymx0; ++l) {
if ((int)layerGroup( ieta-1, l ) == dep) {
Expand All @@ -50,12 +51,12 @@ HcalDDDRecConstants::getEtaBins(const int itype) const {
lmax = l;
dep = layerGroup( ieta-1, l );
}
if (type == 0 && ieta == iEtaMax[type] && dep > 2) break;
if (type == 0 && ieta == hpar->etaMax[type] && dep > 2) break;
}
if (lmax >= lmin) {
if (ieta+1 == hpar->noff[1]) {
} else if (ieta == hpar->noff[1]) {
HcalDDDRecConstants::HcalEtaBin etabin0 = HcalDDDRecConstants::HcalEtaBin(ieta-1,etaTable[ieta-2], etaTable[ieta], nfi, hpar->phioff[type], phibin[ieta-1]);
HcalDDDRecConstants::HcalEtaBin etabin0 = HcalDDDRecConstants::HcalEtaBin(ieta-1,hpar->etaTable[ieta-2], hpar->etaTable[ieta], nfi, hpar->phioff[type], phibin[ieta-1]);
etabin0.depthStart = dep;
etabin0.layer.push_back(std::pair<int,int>(lmin,lmax));
bins.push_back(etabin0);
Expand Down Expand Up @@ -96,16 +97,16 @@ HcalDDDRecConstants::getEtaPhi(int subdet, int ieta, int iphi) const {
(subdet == static_cast<int>(HcalOuter))) { // Use Eta Table
int unit = (int)(phibin[ietaAbs-1]/fiveDegInRad+0.5);
int kphi = (unit == 2) ? ((iphi-1)/2 + 1) : iphi;
double foff = (ietaAbs <= iEtaMax[0]) ? hpar->phioff[0] : hpar->phioff[1];
eta = 0.5*(etaTable[ietaAbs-1]+etaTable[ietaAbs]);
double foff = (ietaAbs <= hpar->etaMax[0]) ? hpar->phioff[0] : hpar->phioff[1];
eta = 0.5*(hpar->etaTable[ietaAbs-1]+hpar->etaTable[ietaAbs]);
phi = foff + (kphi-0.5)*phibin[ietaAbs-1];
} else {
ietaAbs -= iEtaMin[3];
int unit = (int)(phibinHF[ietaAbs-1]/fiveDegInRad+0.5);
ietaAbs -= hpar->etaMin[3];
int unit = (int)(hpar->phitable[ietaAbs-1]/fiveDegInRad+0.5);
int kphi = (unit == 4) ? ((iphi-3)/4 + 1) : ((iphi-1)/2 + 1);
double foff = (unit > 2) ? hpar->phioff[4] : hpar->phioff[2];
eta = 0.5*(hpar->etaTableHF[ietaAbs-1]+hpar->etaTableHF[ietaAbs]);
phi = foff + (kphi-0.5)*phibinHF[ietaAbs-1];
phi = foff + (kphi-0.5)*hpar->phitable[ietaAbs-1];
}
if (ieta < 0) eta = -eta;
if (phi > M_PI) phi -= (2*M_PI);
Expand Down Expand Up @@ -141,7 +142,7 @@ HcalDDDRecConstants::getHCID(int subdet, int ieta, int iphi, int lay,
unit = hcons.unitPhi(phibin[eta-1]);
phi = hcons.phiNumber(phi0,unit);
depth = layerGroup( eta-1, lay-1 );
if (eta == iEtaMin[1]) {
if (eta == hpar->etaMin[1]) {
if (subdet == static_cast<int>(HcalBarrel)) {
if (depth > 2) depth = 2;
} else {
Expand Down Expand Up @@ -171,7 +172,7 @@ double HcalDDDRecConstants::getRZ(int subdet, int ieta, int depth) const {
#ifdef DebugLog
int lay(0);
#endif
if (ietaAbs < nEta) {
if (ietaAbs < hpar->etaMax[1]) {
for (unsigned int k=0; k< layerGroupSize( ietaAbs-1 ); ++k) {
if (depth == (int)layerGroup( ietaAbs-1, k )) {
rz = ((subdet == static_cast<int>(HcalBarrel)) ? (gconsHB[k].first) :
Expand Down Expand Up @@ -307,7 +308,6 @@ unsigned int HcalDDDRecConstants::numberOfCells(HcalSubdetector subdet) const {
} else {
return hcons.numberOfCells(subdet);
}

}

unsigned int HcalDDDRecConstants::nCells(HcalSubdetector subdet) const {
Expand All @@ -334,46 +334,11 @@ unsigned int HcalDDDRecConstants::nCells() const {
}

void HcalDDDRecConstants::initialize(void) {

//Topology Mode
modeTopo_ = hpar->topologyMode;

//Eta grouping
nEta = (int)(hpar->etagroup.size());
if (nEta != (int)(hpar->phigroup.size())) {
edm::LogError("HCalGeom") << "HcalDDDRecConstants: sizes of the vectors "
<< " etaGroup (" << nEta << ") and phiGroup ("
<< hpar->phigroup.size() << ") do not match";
throw cms::Exception("DDException") << "HcalDDDRecConstants: inconsistent array sizes" << nEta << ":" << hpar->phigroup.size();
}

//Layer grouping
for (int i=0; i<nEta; ++i) {
// unsigned int k = hcons.findLayer(i+1, hpar->layerGroupEtaRec);
// if (k < hpar->layerGroupEtaRec.size()) {
// layerGroup[i] = hpar->layerGroupEtaRec[k].layerGroup;
// } else {
// layerGroup[i] = layerGroup[i-1];
// }
#ifdef DebugLog
std::cout << "HcalDDDRecConstants:Read LayerGroup" << i << ":";
for (unsigned int k=0; k<layerGroupSize( i ); k++)
std::cout << " [" << k << "] = " << layerGroup( i, k );
std::cout << std::endl;
#endif
}

iEtaMin = hpar->etaMin;
iEtaMax = hpar->etaMax;
iEtaMin[0] = 1;
iEtaMax[1] = nEta-1;
iEtaMax[2] = iEtaMin[2]+(int)(hpar->rTable.size())-2;

// First eta table
etaTable.clear(); ietaMap.clear(); etaSimValu.clear();
ietaMap.clear(); etaSimValu.clear();
int ieta(0), ietaHB(0), ietaHE(0);
etaTable.push_back(hpar->etaTable[ieta]);
for (int i=0; i<nEta; ++i) {
for( unsigned int i=0; i<hpar->etagroup.size(); ++i ) {
int ef = ieta+1;
ieta += (hpar->etagroup[i]);
if (ieta >= (int)(hpar->etaTable.size())) {
Expand All @@ -385,21 +350,19 @@ void HcalDDDRecConstants::initialize(void) {
<< " at index " << i
<< " of etaTable from SimConstant";
} else {
etaTable.push_back(hpar->etaTable[ieta]);
etaSimValu.push_back(std::pair<int,int>(ef,ieta));
}
for (int k=0; k<(hpar->etagroup[i]); ++k) ietaMap.push_back(i+1);
if (ieta <= iEtaMax[0]) ietaHB = i+1;
if (ieta <= iEtaMin[1]) ietaHE = i+1;
iEtaMax[1] = i+1;
for (int k=0; k<hpar->etagroup[i]; ++k) ietaMap.push_back(i+1);
if (ieta <= hpar->etaMax[0]) ietaHB = i+1;
if (ieta <= hpar->etaMin[1]) ietaHE = i+1;
}
iEtaMin[1] = ietaHE;
iEtaMax[0] = ietaHB;
assert( hpar->etaMin[1] == ietaHE );
assert( hpar->etaMax[0] == ietaHB );

// Then Phi bins
ieta = 0;
phibin.clear(); phiUnitS.clear();
for (int i=0; i<nEta; ++i) {
for (unsigned int i=0; i<hpar->etagroup.size(); ++i) {
double dphi = (hpar->phigroup[i])*(hpar->phibin[ieta]);
phibin.push_back(dphi);
ieta += (hpar->etagroup[i]);
Expand All @@ -408,11 +371,11 @@ void HcalDDDRecConstants::initialize(void) {
int unit = hcons.unitPhi(hpar->phibin[i-1]);
phiUnitS.push_back(unit);
}
phibinHF = hpar->phitable;

#ifdef DebugLog
std::cout << "Modified eta/deltaphi table for " << nEta << " bins" << std::endl;
for (int i=0; i<nEta; ++i)
std::cout << "Eta[" << i << "] = " << etaTable[i] << ":" << etaTable[i+1]
std::cout << "Modified eta/deltaphi table for " << hpar->etagroup.size() << " bins" << std::endl;
for (unsigned int i=0; i<hpar->etagroup.size(); ++i)
std::cout << "Eta[" << i << "] = " << hpar->etaTable[i] << ":" << hpar->etaTable[i+1]
<< ":" << etaSimValu[i].first << ":" << etaSimValu[i].second
<< " PhiBin[" << i << "] = " << phibin[i]/CLHEP::deg <<std::endl;
std::cout << "PhiUnitS";
Expand All @@ -424,38 +387,11 @@ void HcalDDDRecConstants::initialize(void) {
std::cout << " [" << i << "] = " << hpar->etaTableHF[i];
std::cout << std::endl;
std::cout << "PhiBinHF";
for (unsigned int i=0; i<phibinHF.size(); ++i)
std::cout << " [" << i << "] = " << phibinHF[i];
for (unsigned int i=0; i<hpar->phitable.size(); ++i)
std::cout << " [" << i << "] = " << hpar->phitable[i];
std::cout << std::endl;
#endif

//Now the depths
maxDepth = hpar->maxDepth;
maxDepth[0] = maxDepth[1] = 0;
for (int i=0; i<nEta; ++i) {
unsigned int imx = layerGroupSize( i );
int laymax = (imx > 0) ? layerGroup( i, imx-1 ) : 0;
if (i < iEtaMax[0]) {
int laymax0 = (imx > 16) ? layerGroup( i, 16 ) : laymax;
if (i+1 == iEtaMax[0] && laymax0 > 2) laymax0 = 2;
#ifdef DebugLog
std::cout << "HB " << i << " " << imx << " " << laymax << " " << laymax0 << std::endl;
#endif
if (maxDepth[0] < laymax0) maxDepth[0] = laymax0;
}
if (i >= iEtaMin[1]-1 && i < iEtaMax[1]) {
#ifdef DebugLog
std::cout << "HE " << i << " " << imx << " " << laymax << std::endl;
#endif
if (maxDepth[1] < laymax) maxDepth[1] = laymax;
}
}
#ifdef DebugLog
for (int i=0; i<4; ++i)
std::cout << "Detector Type[" << i << "] iEta " << iEtaMin[i] << ":"
<< iEtaMax[i] << " MaxDepth " << maxDepth[i] << std::endl;
#endif

//Now the geometry constants
nModule[0] = hpar->modHB[0];
nHalves[0] = hpar->modHB[1];
Expand All @@ -470,7 +406,7 @@ void HcalDDDRecConstants::initialize(void) {
std::cout << "rHB[" << i << "] = " << gconsHB[i].first << " +- "
<< gconsHB[i].second << std::endl;
#endif
nModule[1] = hpar->modHE[0];;
nModule[1] = hpar->modHE[0];
nHalves[1] = hpar->modHE[1];
for (unsigned int i=0; i<hpar->zHE.size(); ++i) {
gconsHE.push_back(std::pair<double,double>(hpar->zHE[i]/CLHEP::cm,
Expand Down

0 comments on commit ec2c9c1

Please sign in to comment.