CMS 3D CMS Logo
List of all members | Classes | Public Member Functions | Private Attributes
ZdcShowerLibrary Class Reference

#include <ZdcShowerLibrary.h>

Classes

struct  Hit
 

Public Member Functions

int encodePartID (G4int parCode)
 
int getEnergyFromLibrary (const G4ThreeVector &posHit, const G4ThreeVector &momDir, double energy, G4int parCode, HcalZDCDetId::Section section, bool side, int channel)
 
std::vector< Hit > & getHits (const G4Step *aStep, bool &ok)
 
void initRun (G4ParticleTable *theParticleTable)
 
int photonFluctuation (double eav, double esig, double edis)
 
 ZdcShowerLibrary (const std::string &name, edm::ParameterSet const &p)
 
 ~ZdcShowerLibrary ()=default
 

Private Attributes

G4int anuePDG
 
G4int anumuPDG
 
G4int anutauPDG
 
G4int emPDG
 
G4int epPDG
 
G4int etaPDG
 
G4int gammaPDG
 
G4int geantinoPDG
 
std::vector< ZdcShowerLibrary::Hit > hits
 
int npe
 
G4int nuePDG
 
G4int numuPDG
 
G4int nutauPDG
 
G4int pi0PDG
 
bool verbose
 

Detailed Description

Definition at line 22 of file ZdcShowerLibrary.h.

Constructor & Destructor Documentation

◆ ZdcShowerLibrary()

ZdcShowerLibrary::ZdcShowerLibrary ( const std::string &  name,
edm::ParameterSet const &  p 
)

Definition at line 19 of file ZdcShowerLibrary.cc.

References edm::ParameterSet::getUntrackedParameter(), hits, npe, AlCaHLTBitMon_ParallelJobs::p, and verbose.

19  {
20  edm::ParameterSet m_HS = p.getParameter<edm::ParameterSet>("ZdcShowerLibrary");
21  verbose = m_HS.getUntrackedParameter<int>("Verbosity", 0);
22 
23  npe = 9; // number of channels or fibers where the energy will be deposited
24  hits.reserve(npe);
25 }
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
ZdcShowerLibrary::hits
std::vector< ZdcShowerLibrary::Hit > hits
Definition: ZdcShowerLibrary.h:57
AlCaHLTBitMon_ParallelJobs.p
def p
Definition: AlCaHLTBitMon_ParallelJobs.py:153

◆ ~ZdcShowerLibrary()

ZdcShowerLibrary::~ZdcShowerLibrary ( )
default

Member Function Documentation

◆ encodePartID()

int ZdcShowerLibrary::encodePartID ( G4int  parCode)

Definition at line 272 of file ZdcShowerLibrary.cc.

References emPDG, epPDG, and gammaPDG.

Referenced by getEnergyFromLibrary(), and getHits().

272  {
273  G4int iparCode = 1;
274  if (parCode == emPDG || parCode == epPDG || parCode == gammaPDG) {
275  iparCode = 0;
276  } else {
277  return iparCode;
278  }
279  return iparCode;
280 }

◆ getEnergyFromLibrary()

int ZdcShowerLibrary::getEnergyFromLibrary ( const G4ThreeVector &  posHit,
const G4ThreeVector &  momDir,
double  energy,
G4int  parCode,
HcalZDCDetId::Section  section,
bool  side,
int  channel 
)

Definition at line 148 of file ZdcShowerLibrary.cc.

References encodePartID(), hcalRecHitTable_cff::energy, fact, HMesonGammaMonitor_cff::nphotons, photonFluctuation(), funct::pow(), hgcalPlots::section, theXChannelBoundaries, and ecalPiZeroTask_cfi::X0.

Referenced by getHits().

154  {
155  int nphotons = -1;
156 
157  energy = energy / GeV;
158 
159  edm::LogVerbatim("ZdcShower") << "\n ZdcShowerLibrary::getEnergyFromLibrary input/output variables:"
160  << " phi: " << momDir.phi() / CLHEP::deg << " theta: " << momDir.theta() / CLHEP::deg
161  << " xin : " << hitPoint.x() << " yin : " << hitPoint.y() << " zin : " << hitPoint.z()
162  << " track en: " << energy << "(GeV)"
163  << " section: " << section << " side: " << side << " channel: " << channel
164  << " partID: " << parCode;
165 
166  // these varables are not used for now
167  //float phi = momDir.phi() / CLHEP::deg;
168  //float theta = momDir.theta() / CLHEP::deg;
169  //float zin = hitPoint.z();
170  //int isection = int(section);
171  //int iside = (side)? 1 : 2;
172 
173  int iparCode = encodePartID(parCode);
174 
175  double eav = 0.;
176  double esig = 0.;
177  double edis = 0.;
178 
179  float xin = hitPoint.x();
180  float yin = hitPoint.y();
181  float fact = 0.;
182 
183  if (section == 1 && iparCode != 0) {
184  if (channel < 5)
185  if (((theXChannelBoundaries[channel - 1]) < (xin + X0)) && ((xin + X0) <= theXChannelBoundaries[channel]))
186  fact = 0.18;
187  if (channel == 5)
188  if (theXChannelBoundaries[channel - 1] < xin + X0)
189  fact = 0.18;
190  }
191 
192  if (section == 2 && iparCode != 0) {
193  if (channel == 1)
194  fact = 0.34;
195  if (channel == 2)
196  fact = 0.24;
197  if (channel == 3)
198  fact = 0.17;
199  if (channel == 4)
200  fact = 0.07;
201  }
202  if (section == 1 && iparCode == 0) {
203  if (channel < 5)
204  if (((theXChannelBoundaries[channel - 1]) < (xin + X0)) && ((xin + X0) <= theXChannelBoundaries[channel]))
205  fact = 1.;
206  if (channel == 5)
207  if (theXChannelBoundaries[channel - 1] < xin + X0)
208  fact = 1.0;
209  }
210 
211  //change to cm for parametrization
212  yin = yin / cm;
213  xin = xin / cm;
214 
215  if (iparCode == 0) {
216  eav = ((((((-0.0002 * xin - 2.0e-13) * xin + 0.0022) * xin + 1.0e-11) * xin - 0.0217) * xin - 3.0e-10) * xin +
217  1.0028) *
218  (((0.0001 * yin + 0.0056) * yin + 0.0508) * yin + 1.0) * 300.0 * pow((energy / 300.0), 0.99); // EM
219  esig = ((((((0.0005 * xin - 1.0e-12) * xin - 0.0052) * xin + 5.0e-11) * xin + 0.032) * xin - 2.0e-10) * xin + 1.0) *
220  (((0.0006 * yin + 0.0071) * yin - 0.031) * yin + 1.0) * 30.0 * pow((energy / 300.0), 0.54); // EM
221  edis = 1.0;
222  } else {
223  eav = ((((((-0.0002 * xin - 2.0e-13) * xin + 0.0022) * xin + 1.0e-11) * xin - 0.0217) * xin - 3.0e-10) * xin +
224  1.0028) *
225  (((0.0001 * yin + 0.0056) * yin + 0.0508) * yin + 1.0) * 300.0 * pow((energy / 300.0), 1.12); // HD
226  esig = ((((((0.0005 * xin - 1.0e-12) * xin - 0.0052) * xin + 5.0e-11) * xin + 0.032) * xin - 2.0e-10) * xin + 1.0) *
227  (((0.0006 * yin + 0.0071) * yin - 0.031) * yin + 1.0) * 54.0 * pow((energy / 300.0), 0.93); //HD
228  edis = 3.0;
229  }
230 
231  if (eav < 0. || edis < 0.) {
232  edm::LogVerbatim("ZdcShower") << " Negative everage energy from parametrization \n"
233  << " xin: " << xin << "(cm)"
234  << " yin: " << yin << "(cm)"
235  << " track en: " << energy << "(GeV)"
236  << " eaverage: " << eav / CLHEP::GeV << " (GeV)"
237  << " esigma: " << esig / CLHEP::GeV << " (GeV)"
238  << " edist: " << edis << " (GeV)"
239  << " dE hit: " << nphotons / CLHEP::GeV << " (GeV)";
240  return 0;
241  }
242 
243  // Convert from GeV to MeV for the code
244  eav = eav * CLHEP::GeV;
245  esig = esig * CLHEP::GeV;
246 
247  while (nphotons == -1 || nphotons > int(eav + 5. * esig))
248  nphotons = static_cast<int>(fact * photonFluctuation(eav, esig, edis));
249 
250  edm::LogVerbatim("ZdcShower") << " track en: " << energy << "(GeV)"
251  << " eaverage: " << eav / CLHEP::GeV << " (GeV)"
252  << " esigma: " << esig / CLHEP::GeV << " (GeV)"
253  << " edist: " << edis << " (GeV)"
254  << " dE hit: " << nphotons / CLHEP::GeV << " (GeV)";
255 
256  return nphotons;
257 }
edm::LogVerbatim
Log< level::Info, true > LogVerbatim
Definition: MessageLogger.h:134
ecalPiZeroTask_cfi.X0
X0
Definition: ecalPiZeroTask_cfi.py:25
hcalRecHitTable_cff.energy
energy
Definition: hcalRecHitTable_cff.py:13
theXChannelBoundaries
static const double theXChannelBoundaries[]
Definition: ZdcHardcodeGeometryData.h:35
HMesonGammaMonitor_cff.nphotons
nphotons
Definition: HMesonGammaMonitor_cff.py:10
ZdcShowerLibrary::photonFluctuation
int photonFluctuation(double eav, double esig, double edis)
Definition: ZdcShowerLibrary.cc:259
fact
const double fact
Definition: NuclearInteractionFTFSimulator.cc:74
ZdcShowerLibrary::encodePartID
int encodePartID(G4int parCode)
Definition: ZdcShowerLibrary.cc:272
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29

◆ getHits()

std::vector< ZdcShowerLibrary::Hit > & ZdcShowerLibrary::getHits ( const G4Step *  aStep,
bool &  ok 
)

Definition at line 48 of file ZdcShowerLibrary.cc.

References ZdcShowerLibrary::Hit::DeEM, ZdcShowerLibrary::Hit::DeHad, hcalRecHitTable_cff::depth, ZdcShowerLibrary::Hit::depth, ZdcShowerLibrary::Hit::detID, HcalZDCDetId::EM, encodePartID(), hcalRecHitTable_cff::energy, ZdcShowerLibrary::Hit::entryLocal, geantinoPDG, getEnergyFromLibrary(), HcalZDCDetId::HAD, hits, mps_fire::i, npe, convertSQLiteXML::ok, ZdcShowerLibrary::Hit::position, hgcalPlots::section, theXChannelBoundaries, theZHadChannelBoundaries, theZSectionBoundaries, ZdcShowerLibrary::Hit::time, HLT_2023v12_cff::track, funct::true, ecalPiZeroTask_cfi::X0, geometryCSVtoXML::xx, Y0, geometryCSVtoXML::yy, Z0, and geometryCSVtoXML::zz.

48  {
49  G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
50  G4StepPoint* postStepPoint = aStep->GetPostStepPoint();
51  G4Track* track = aStep->GetTrack();
52 
53  const G4DynamicParticle* aParticle = track->GetDynamicParticle();
54  const G4ThreeVector& momDir = aParticle->GetMomentumDirection();
55  double energy = preStepPoint->GetKineticEnergy();
56  G4ThreeVector hitPoint = preStepPoint->GetPosition();
57  const G4ThreeVector& hitPointOrig = preStepPoint->GetPosition();
58  G4int parCode = track->GetDefinition()->GetPDGEncoding();
59 
60  hits.clear();
61 
62  ok = false;
63  if (parCode == geantinoPDG)
64  return hits;
65  ok = true;
66 
67  G4ThreeVector pos;
68  G4ThreeVector posLocal;
69  double tSlice = (postStepPoint->GetGlobalTime()) / nanosecond;
70 
71  int nHit = 0;
72  HcalZDCDetId::Section section;
73  bool side = false;
74  int channel = 0;
75  double xx, yy, zz;
76  double xxlocal, yylocal, zzlocal;
77 
78  ZdcShowerLibrary::Hit oneHit;
79  side = (hitPointOrig.z() > 0.) ? true : false;
80 
81  float xWidthEM = fabs(theXChannelBoundaries[0] - theXChannelBoundaries[1]);
82  float zWidthEM = fabs(theZSectionBoundaries[0] - theZSectionBoundaries[1]);
83  float zWidthHAD = fabs(theZHadChannelBoundaries[0] - theZHadChannelBoundaries[1]);
84 
85  for (int i = 0; i < npe; i++) {
86  if (i < 5) {
87  section = HcalZDCDetId::EM;
88  channel = i + 1;
89  xxlocal = theXChannelBoundaries[i] + (xWidthEM / 2.);
90  xx = xxlocal + X0;
91  yy = 0.0;
92  yylocal = yy + Y0;
93  zzlocal = theZSectionBoundaries[0] + (zWidthEM / 2.);
94  zz = (hitPointOrig.z() > 0.) ? zzlocal + Z0 : zzlocal - Z0;
95  pos = G4ThreeVector(xx, yy, zz);
96  posLocal = G4ThreeVector(xxlocal, yylocal, zzlocal);
97  }
98  if (i > 4) {
99  section = HcalZDCDetId::HAD;
100  channel = i - 4;
101  xxlocal = 0.0;
102  xx = xxlocal + X0;
103  yylocal = 0;
104  yy = yylocal + Y0;
105  zzlocal = (hitPointOrig.z() > 0.) ? theZHadChannelBoundaries[i - 5] + (zWidthHAD / 2.)
106  : theZHadChannelBoundaries[i - 5] - (zWidthHAD / 2.);
107  zz = (hitPointOrig.z() > 0.) ? zzlocal + Z0 : zzlocal - Z0;
108  pos = G4ThreeVector(xx, yy, zz);
109  posLocal = G4ThreeVector(xxlocal, yylocal, zzlocal);
110  }
111 
112  oneHit.position = pos;
113  oneHit.entryLocal = posLocal;
114  oneHit.depth = channel;
115  oneHit.time = tSlice;
116  oneHit.detID = HcalZDCDetId(section, side, channel);
117 
118  // Note: coodinates of hit are relative to center of detector (X0,Y0,Z0)
119  hitPoint.setX(hitPointOrig.x() - X0);
120  hitPoint.setY(hitPointOrig.y() - Y0);
121  double setZ = (hitPointOrig.z() > 0.) ? hitPointOrig.z() - Z0 : fabs(hitPointOrig.z()) - Z0;
122  hitPoint.setZ(setZ);
123 
124  int dE = getEnergyFromLibrary(hitPoint, momDir, energy, parCode, section, side, channel);
125 
126  int iparCode = encodePartID(parCode);
127  if (iparCode == 0) {
128  oneHit.DeEM = dE;
129  oneHit.DeHad = 0.;
130  } else {
131  oneHit.DeEM = 0;
132  oneHit.DeHad = dE;
133  }
134 
135  hits.push_back(oneHit);
136 
137  edm::LogVerbatim("ZdcShower") << "\nZdcShowerLibrary:Generated Hit " << nHit << " orig hit pos " << hitPointOrig
138  << " orig hit pos local coord" << hitPoint << " new position "
139  << (hits[nHit].position) << " Channel " << (hits[nHit].depth) << " side " << side
140  << " Time " << (hits[nHit].time) << " DetectorID " << (hits[nHit].detID)
141  << " Had Energy " << (hits[nHit].DeHad) << " EM Energy " << (hits[nHit].DeEM)
142  << "\n";
143  nHit++;
144  }
145  return hits;
146 }
edm::LogVerbatim
Log< level::Info, true > LogVerbatim
Definition: MessageLogger.h:134
mps_fire.i
i
Definition: mps_fire.py:429
ZdcShowerLibrary::Hit::entryLocal
G4ThreeVector entryLocal
Definition: ZdcShowerLibrary.h:30
ecalPiZeroTask_cfi.X0
X0
Definition: ecalPiZeroTask_cfi.py:25
hcalRecHitTable_cff.energy
energy
Definition: hcalRecHitTable_cff.py:13
theXChannelBoundaries
static const double theXChannelBoundaries[]
Definition: ZdcHardcodeGeometryData.h:35
ZdcShowerLibrary::getEnergyFromLibrary
int getEnergyFromLibrary(const G4ThreeVector &posHit, const G4ThreeVector &momDir, double energy, G4int parCode, HcalZDCDetId::Section section, bool side, int channel)
Definition: ZdcShowerLibrary.cc:148
theZHadChannelBoundaries
static const double theZHadChannelBoundaries[]
Definition: ZdcHardcodeGeometryData.h:38
ZdcShowerLibrary::hits
std::vector< ZdcShowerLibrary::Hit > hits
Definition: ZdcShowerLibrary.h:57
Z0
static const double Z0
Definition: ZdcHardcodeGeometryData.h:25
hcalRecHitTable_cff.depth
depth
Definition: hcalRecHitTable_cff.py:17
ZdcShowerLibrary::encodePartID
int encodePartID(G4int parCode)
Definition: ZdcShowerLibrary.cc:272
Y0
static const double Y0
Definition: ZdcHardcodeGeometryData.h:26
theZSectionBoundaries
static const double theZSectionBoundaries[]
Definition: ZdcHardcodeGeometryData.h:34
HLT_2023v12_cff.track
track
Definition: HLT_2023v12_cff.py:9624

◆ initRun()

void ZdcShowerLibrary::initRun ( G4ParticleTable *  theParticleTable)

Definition at line 27 of file ZdcShowerLibrary.cc.

References anuePDG, anumuPDG, anutauPDG, emPDG, epPDG, etaPDG, gammaPDG, geantinoPDG, nuePDG, numuPDG, nutauPDG, and pi0PDG.

27  {
28  G4String parName;
29  emPDG = theParticleTable->FindParticle(parName = "e-")->GetPDGEncoding();
30  epPDG = theParticleTable->FindParticle(parName = "e+")->GetPDGEncoding();
31  gammaPDG = theParticleTable->FindParticle(parName = "gamma")->GetPDGEncoding();
32  pi0PDG = theParticleTable->FindParticle(parName = "pi0")->GetPDGEncoding();
33  etaPDG = theParticleTable->FindParticle(parName = "eta")->GetPDGEncoding();
34  nuePDG = theParticleTable->FindParticle(parName = "nu_e")->GetPDGEncoding();
35  numuPDG = theParticleTable->FindParticle(parName = "nu_mu")->GetPDGEncoding();
36  nutauPDG = theParticleTable->FindParticle(parName = "nu_tau")->GetPDGEncoding();
37  anuePDG = theParticleTable->FindParticle(parName = "anti_nu_e")->GetPDGEncoding();
38  anumuPDG = theParticleTable->FindParticle(parName = "anti_nu_mu")->GetPDGEncoding();
39  anutauPDG = theParticleTable->FindParticle(parName = "anti_nu_tau")->GetPDGEncoding();
40  geantinoPDG = theParticleTable->FindParticle(parName = "geantino")->GetPDGEncoding();
41  edm::LogVerbatim("ZdcShower") << "ZdcShowerLibrary: Particle codes for e- = " << emPDG << ", e+ = " << epPDG
42  << ", gamma = " << gammaPDG << ", pi0 = " << pi0PDG << ", eta = " << etaPDG
43  << ", geantino = " << geantinoPDG << "\n nu_e = " << nuePDG
44  << ", nu_mu = " << numuPDG << ", nu_tau = " << nutauPDG << ", anti_nu_e = " << anuePDG
45  << ", anti_nu_mu = " << anumuPDG << ", anti_nu_tau = " << anutauPDG;
46 }
edm::LogVerbatim
Log< level::Info, true > LogVerbatim
Definition: MessageLogger.h:134

◆ photonFluctuation()

int ZdcShowerLibrary::photonFluctuation ( double  eav,
double  esig,
double  edis 
)

Definition at line 259 of file ZdcShowerLibrary.cc.

References createfilelist::int.

Referenced by getEnergyFromLibrary().

259  {
260  int nphot = 0;
261  double efluct = 0.;
262  if (edis == 1.0)
263  efluct = eav + esig * CLHEP::RandGaussQ::shoot();
264  if (edis == 3.0)
265  efluct = eav + esig * CLHEP::RandLandau::shoot();
266  nphot = int(efluct);
267  if (nphot < 0)
268  nphot = 0;
269  return nphot;
270 }

Member Data Documentation

◆ anuePDG

G4int ZdcShowerLibrary::anuePDG
private

Definition at line 55 of file ZdcShowerLibrary.h.

Referenced by initRun().

◆ anumuPDG

G4int ZdcShowerLibrary::anumuPDG
private

Definition at line 55 of file ZdcShowerLibrary.h.

Referenced by initRun().

◆ anutauPDG

G4int ZdcShowerLibrary::anutauPDG
private

Definition at line 55 of file ZdcShowerLibrary.h.

Referenced by initRun().

◆ emPDG

G4int ZdcShowerLibrary::emPDG
private

Definition at line 53 of file ZdcShowerLibrary.h.

Referenced by encodePartID(), and initRun().

◆ epPDG

G4int ZdcShowerLibrary::epPDG
private

Definition at line 53 of file ZdcShowerLibrary.h.

Referenced by encodePartID(), and initRun().

◆ etaPDG

G4int ZdcShowerLibrary::etaPDG
private

Definition at line 54 of file ZdcShowerLibrary.h.

Referenced by initRun().

◆ gammaPDG

G4int ZdcShowerLibrary::gammaPDG
private

Definition at line 53 of file ZdcShowerLibrary.h.

Referenced by encodePartID(), and initRun().

◆ geantinoPDG

G4int ZdcShowerLibrary::geantinoPDG
private

Definition at line 55 of file ZdcShowerLibrary.h.

Referenced by getHits(), and initRun().

◆ hits

std::vector<ZdcShowerLibrary::Hit> ZdcShowerLibrary::hits
private

Definition at line 57 of file ZdcShowerLibrary.h.

Referenced by getHits(), and ZdcShowerLibrary().

◆ npe

int ZdcShowerLibrary::npe
private

Definition at line 56 of file ZdcShowerLibrary.h.

Referenced by getHits(), and ZdcShowerLibrary().

◆ nuePDG

G4int ZdcShowerLibrary::nuePDG
private

Definition at line 54 of file ZdcShowerLibrary.h.

Referenced by initRun().

◆ numuPDG

G4int ZdcShowerLibrary::numuPDG
private

Definition at line 54 of file ZdcShowerLibrary.h.

Referenced by initRun().

◆ nutauPDG

G4int ZdcShowerLibrary::nutauPDG
private

Definition at line 54 of file ZdcShowerLibrary.h.

Referenced by initRun().

◆ pi0PDG

G4int ZdcShowerLibrary::pi0PDG
private

Definition at line 54 of file ZdcShowerLibrary.h.

Referenced by initRun().

◆ verbose

bool ZdcShowerLibrary::verbose
private

Definition at line 52 of file ZdcShowerLibrary.h.

Referenced by ZdcShowerLibrary().

Generated for CMSSW Reference Manual by doxygen 1.8.14