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 (G4Step *aStep, bool &ok)
 
void initRun (G4ParticleTable *theParticleTable)
 
int photonFluctuation (double eav, double esig, double edis)
 
 ZdcShowerLibrary (std::string &name, const DDCompactView &cpv, edm::ParameterSet const &p)
 
 ~ZdcShowerLibrary ()
 

Private Attributes

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

Detailed Description

Definition at line 23 of file ZdcShowerLibrary.h.

Constructor & Destructor Documentation

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

Definition at line 19 of file ZdcShowerLibrary.cc.

References edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), hits, and npe.

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

Definition at line 29 of file ZdcShowerLibrary.cc.

29  {
30 }

Member Function Documentation

int ZdcShowerLibrary::encodePartID ( G4int  parCode)

Definition at line 284 of file ZdcShowerLibrary.cc.

References emPDG, epPDG, and gammaPDG.

Referenced by getEnergyFromLibrary(), and getHits().

284  {
285  G4int iparCode = 1;
286  if (parCode == emPDG ||
287  parCode == epPDG ||
288  parCode == gammaPDG) {
289  iparCode = 0;
290  } else { return iparCode; }
291  return iparCode;
292 }
int ZdcShowerLibrary::getEnergyFromLibrary ( const G4ThreeVector &  posHit,
const G4ThreeVector &  momDir,
double  energy,
G4int  parCode,
HcalZDCDetId::Section  section,
bool  side,
int  channel 
)

Definition at line 166 of file ZdcShowerLibrary.cc.

References encodePartID(), fact, GeV, createfilelist::int, LogDebug, photonFluctuation(), funct::pow(), theXChannelBoundaries, and X0.

Referenced by getHits().

167  {
168  int nphotons = -1;
169 
170  energy =energy/GeV;
171 
172  LogDebug("ZdcShower")
173  //std::cout
174  <<"\n ZdcShowerLibrary::getEnergyFromLibrary input/output variables:"
175  <<" phi: "<<59.2956*momDir.phi()
176  <<" theta: "<<59.2956*momDir.theta()
177  <<" xin : "<<hitPoint.x()
178  <<" yin : "<<hitPoint.y()
179  <<" zin : "<<hitPoint.z()
180  <<" track en: " <<energy<< "(GeV)"
181  <<" section: "<<section
182  <<" side: "<<side
183  <<" channel: "<< channel
184  <<" partID: "<<parCode;
185  //std::cout<<std::endl;
186 
187  // these varables are not used for now
188  //float phi = 59.2956*momDir.phi();
189  //float theta = 59.2956*momDir.theta();
190  //float zin = hitPoint.z();
191  //int isection = int(section);
192  //int iside = (side)? 1 : 2;
193 
194  int iparCode = encodePartID(parCode);
195 
196  double eav = 0.;
197  double esig = 0.;
198  double edis = 0.;
199 
200  float xin = hitPoint.x();
201  float yin = hitPoint.y();
202  float fact = 0.;
203 
204  if(section == 1 && iparCode !=0){
205  if(channel < 5 )
206  if(((theXChannelBoundaries[channel-1])< (xin + X0)) && ((xin + X0)<= theXChannelBoundaries[channel]))fact= 0.18;
207  if(channel ==5 )
208  if(theXChannelBoundaries[channel-1]< xin + X0)fact = 0.18;
209  }
210 
211  if(section == 2 && iparCode !=0){
212  if(channel == 1)fact = 0.34;
213  if(channel == 2)fact = 0.24;
214  if(channel == 3)fact = 0.17;
215  if(channel == 4)fact = 0.07;
216  }
217  if(section == 1 && iparCode ==0){
218  if(channel < 5 )
219  if(((theXChannelBoundaries[channel-1])< (xin + X0)) && ((xin + X0)<= theXChannelBoundaries[channel]))fact= 1.;
220  if(channel ==5 )
221  if(theXChannelBoundaries[channel-1]< xin + X0)fact = 1.0;
222  }
223 
224  //change to cm for parametrization
225  yin = yin/cm;
226  xin = xin/cm;
227 
228  if (iparCode==0){
229  eav = ((((((-0.0002*xin-2.0e-13)*xin+0.0022)*xin+1.0e-11)*xin-0.0217)*xin-3.0e-10)*xin+1.0028)*
230  (((0.0001*yin + 0.0056)*yin + 0.0508)*yin + 1.0)*300.0*pow((energy/300.0),0.99); // EM
231  esig = ((((((0.0005*xin - 1.0e-12)*xin - 0.0052)*xin + 5.0e-11)*xin + 0.032)*xin -
232  2.0e-10)*xin + 1.0)*(((0.0006*yin + 0.0071)*yin - 0.031)*yin + 1.0)*30.0*pow((energy/300.0),0.54);// EM
233  edis = 1.0;
234  }else{
235  eav = ((((((-0.0002*xin-2.0e-13)*xin+0.0022)*xin+1.0e-11)*xin-0.0217)*xin-3.0e-10)*xin+1.0028)*
236  (((0.0001*yin + 0.0056)*yin + 0.0508)*yin + 1.0)*300.0*pow((energy/300.0),1.12);// HD
237  esig = ((((((0.0005*xin - 1.0e-12)*xin - 0.0052)*xin + 5.0e-11)*xin + 0.032)*xin -
238  2.0e-10)*xin + 1.0)*(((0.0006*yin + 0.0071)*yin - 0.031)*yin + 1.0)*54.0*pow((energy/300.0),0.93); //HD
239  edis = 3.0;
240  }
241 
242  if(eav <0. || edis <0.){
243  LogDebug("ZdcShower")
244  <<" Negative everage energy from parametrization \n"
245  <<" xin: "<<xin<< "(cm)"
246  <<" yin: "<<yin<< "(cm)"
247  <<" track en: " <<energy<< "(GeV)"
248  <<" eaverage: "<<eav/GeV << " (GeV)"
249  <<" esigma: "<<esig/GeV << " (GeV)"
250  <<" edist: "<<edis << " (GeV)"
251  <<" dE hit: "<<nphotons/GeV<< " (GeV)";
252  return 0;
253  }
254 
255  // Convert from GeV to MeV for the code
256  eav = eav*GeV;
257  esig= esig*GeV;
258 
259  while(nphotons == -1 || nphotons > int(eav + 5.*esig))
260  nphotons = (int)(fact*photonFluctuation(eav, esig, edis));
261 
262  LogDebug("ZdcShower")
263  //std::cout
264  <<" track en: " <<energy<< "(GeV)"
265  <<" eaverage: "<<eav/GeV << " (GeV)"
266  <<" esigma: "<<esig/GeV << " (GeV)"
267  <<" edist: "<<edis << " (GeV)"
268  <<" dE hit: "<<nphotons/GeV<< " (GeV)";
269  //std::cout<<std::endl;
270 
271  return nphotons;
272 }
#define LogDebug(id)
const double GeV
Definition: MathUtil.h:16
static const double theXChannelBoundaries[]
int photonFluctuation(double eav, double esig, double edis)
const double fact
int encodePartID(G4int parCode)
static const double X0
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
std::vector< ZdcShowerLibrary::Hit > & ZdcShowerLibrary::getHits ( G4Step *  aStep,
bool &  ok 
)

Definition at line 56 of file ZdcShowerLibrary.cc.

References ZdcShowerLibrary::Hit::DeEM, ZdcShowerLibrary::Hit::DeHad, particleFlowClusterECALTimeSelected_cfi::depth, ZdcShowerLibrary::Hit::depth, ZdcShowerLibrary::Hit::detID, HcalZDCDetId::EM, encodePartID(), ZdcShowerLibrary::Hit::entryLocal, geantinoPDG, getEnergyFromLibrary(), HcalZDCDetId::HAD, hits, mps_fire::i, LogDebug, npe, ZdcShowerLibrary::Hit::position, vertexPlots::section, theXChannelBoundaries, theZHadChannelBoundaries, theZSectionBoundaries, ZdcShowerLibrary::Hit::time, HiIsolationCommonParameters_cff::track, X0, geometryCSVtoXML::xx, Y0, geometryCSVtoXML::yy, Z0, and geometryCSVtoXML::zz.

Referenced by ZdcSD::getFromLibrary().

56  {
57 
58  G4StepPoint * preStepPoint = aStep->GetPreStepPoint();
59  G4StepPoint * postStepPoint = aStep->GetPostStepPoint();
60  G4Track * track = aStep->GetTrack();
61 
62  const G4DynamicParticle *aParticle = track->GetDynamicParticle();
63  const G4ThreeVector& momDir = aParticle->GetMomentumDirection();
64  double energy = preStepPoint->GetKineticEnergy();
65  G4ThreeVector hitPoint = preStepPoint->GetPosition();
66  const G4ThreeVector& hitPointOrig = preStepPoint->GetPosition();
67  G4int parCode = track->GetDefinition()->GetPDGEncoding();
68 
69  hits.clear();
70 
71  ok = false;
72  if (parCode == geantinoPDG)
73  return hits;
74  ok = true;
75 
76  G4ThreeVector pos;
77  G4ThreeVector posLocal;
78  double tSlice = (postStepPoint->GetGlobalTime())/nanosecond;
79 
80  int nHit = 0;
82  bool side = false;
83  int channel = 0;
84  double xx,yy,zz;
85  double xxlocal, yylocal, zzlocal;
86 
87  ZdcShowerLibrary::Hit oneHit;
88  side = (hitPointOrig.z() > 0.) ? true : false;
89 
90  float xWidthEM = fabs(theXChannelBoundaries[0] - theXChannelBoundaries[1]);
91  float zWidthEM = fabs(theZSectionBoundaries[0] - theZSectionBoundaries[1]);
92  float zWidthHAD = fabs(theZHadChannelBoundaries[0] -theZHadChannelBoundaries[1]);
93 
94  for (int i = 0; i < npe; i++) {
95  if(i < 5){
96  section = HcalZDCDetId::EM;
97  channel = i+1;
98  xxlocal = theXChannelBoundaries[i]+(xWidthEM/2.);
99  xx = xxlocal + X0;
100  yy = 0.0;
101  yylocal = yy + Y0;
102  zzlocal = theZSectionBoundaries[0]+(zWidthEM/2.);
103  zz = (hitPointOrig.z() > 0.) ? zzlocal + Z0 : zzlocal - Z0;
104  pos = G4ThreeVector(xx,yy,zz);
105  posLocal = G4ThreeVector(xxlocal,yylocal,zzlocal);
106  }
107  if(i > 4){
108  section = HcalZDCDetId::HAD;
109  channel = i - 4;
110  xxlocal = 0.0;
111  xx = xxlocal + X0;
112  yylocal = 0;
113  yy = yylocal + Y0;
114  zzlocal = (hitPointOrig.z() > 0.) ?
115  theZHadChannelBoundaries[i-5] + (zWidthHAD/2.) : theZHadChannelBoundaries[i-5] - (zWidthHAD/2.);
116  zz = (hitPointOrig.z() > 0.) ? zzlocal + Z0 : zzlocal - Z0;
117  pos = G4ThreeVector(xx,yy,zz);
118  posLocal = G4ThreeVector(xxlocal,yylocal,zzlocal);
119  }
120 
121  oneHit.position = pos;
122  oneHit.entryLocal = posLocal;
123  oneHit.depth = channel;
124  oneHit.time = tSlice;
125  oneHit.detID = HcalZDCDetId(section,side,channel);
126 
127 
128  // Note: coodinates of hit are relative to center of detector (X0,Y0,Z0)
129  hitPoint.setX(hitPointOrig.x()-X0);
130  hitPoint.setY(hitPointOrig.y()-Y0);
131  double setZ= (hitPointOrig.z()> 0.) ? hitPointOrig.z()- Z0 : fabs(hitPointOrig.z()) - Z0;
132  hitPoint.setZ(setZ);
133 
134  int dE = getEnergyFromLibrary(hitPoint,momDir,energy,parCode,section,side,channel);
135 
136  int iparCode = encodePartID(parCode);
137  if (iparCode == 0 ) {
138  oneHit.DeEM = dE;
139  oneHit.DeHad = 0.;
140  } else {
141  oneHit.DeEM = 0;
142  oneHit.DeHad = dE;
143  }
144 
145  hits.push_back(oneHit);
146 
147  LogDebug("ZdcShower")
148  //std::cout
149  << "\nZdcShowerLibrary:Generated Hit " << nHit
150  <<" orig hit pos " << hitPointOrig
151  <<" orig hit pos local coord" << hitPoint
152  <<" new position " << (hits[nHit].position)
153  <<" Channel " <<(hits[nHit].depth)
154  <<" side "<< side
155  <<" Time " <<(hits[nHit].time)
156  <<" DetectorID " << (hits[nHit].detID)
157  <<" Had Energy " << (hits[nHit].DeHad)
158  <<" EM Energy " << (hits[nHit].DeEM)
159  <<"\n";
160  nHit++;
161  }
162  return hits;
163 }
#define LogDebug(id)
G4ThreeVector entryLocal
static const double Z0
int getEnergyFromLibrary(const G4ThreeVector &posHit, const G4ThreeVector &momDir, double energy, G4int parCode, HcalZDCDetId::Section section, bool side, int channel)
static const double theXChannelBoundaries[]
std::vector< ZdcShowerLibrary::Hit > hits
static const double theZHadChannelBoundaries[]
static const double theZSectionBoundaries[]
static const double Y0
int encodePartID(G4int parCode)
static const double X0
void ZdcShowerLibrary::initRun ( G4ParticleTable *  theParticleTable)

Definition at line 32 of file ZdcShowerLibrary.cc.

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

Referenced by ZdcSD::initRun().

32  {
33  G4String parName;
34  emPDG = theParticleTable->FindParticle(parName="e-")->GetPDGEncoding();
35  epPDG = theParticleTable->FindParticle(parName="e+")->GetPDGEncoding();
36  gammaPDG = theParticleTable->FindParticle(parName="gamma")->GetPDGEncoding();
37  pi0PDG = theParticleTable->FindParticle(parName="pi0")->GetPDGEncoding();
38  etaPDG = theParticleTable->FindParticle(parName="eta")->GetPDGEncoding();
39  nuePDG = theParticleTable->FindParticle(parName="nu_e")->GetPDGEncoding();
40  numuPDG = theParticleTable->FindParticle(parName="nu_mu")->GetPDGEncoding();
41  nutauPDG= theParticleTable->FindParticle(parName="nu_tau")->GetPDGEncoding();
42  anuePDG = theParticleTable->FindParticle(parName="anti_nu_e")->GetPDGEncoding();
43  anumuPDG= theParticleTable->FindParticle(parName="anti_nu_mu")->GetPDGEncoding();
44  anutauPDG= theParticleTable->FindParticle(parName="anti_nu_tau")->GetPDGEncoding();
45  geantinoPDG= theParticleTable->FindParticle(parName="geantino")->GetPDGEncoding();
46  LogDebug("ZdcShower")<< "ZdcShowerLibrary: Particle codes for e- = " << emPDG
47  << ", e+ = " << epPDG << ", gamma = " << gammaPDG
48  << ", pi0 = " << pi0PDG << ", eta = " << etaPDG
49  << ", geantino = " << geantinoPDG << "\n nu_e = "
50  << nuePDG << ", nu_mu = " << numuPDG << ", nu_tau = "
51  << nutauPDG << ", anti_nu_e = " << anuePDG
52  << ", anti_nu_mu = " << anumuPDG << ", anti_nu_tau = "
53  << anutauPDG;
54 }
#define LogDebug(id)
int ZdcShowerLibrary::photonFluctuation ( double  eav,
double  esig,
double  edis 
)

Definition at line 274 of file ZdcShowerLibrary.cc.

References createfilelist::int.

Referenced by getEnergyFromLibrary().

274  {
275  int nphot=0;
276  double efluct = 0.;
277  if(edis == 1.0) efluct = eav+esig*CLHEP::RandGaussQ::shoot();
278  if(edis == 3.0) efluct = eav+esig*CLHEP::RandLandau::shoot();
279  nphot = int(efluct);
280  if(nphot<0)nphot = 0;
281  return nphot;
282 }

Member Data Documentation

G4int ZdcShowerLibrary::anuePDG
private

Definition at line 59 of file ZdcShowerLibrary.h.

Referenced by initRun().

G4int ZdcShowerLibrary::anumuPDG
private

Definition at line 59 of file ZdcShowerLibrary.h.

Referenced by initRun().

G4int ZdcShowerLibrary::anutauPDG
private

Definition at line 59 of file ZdcShowerLibrary.h.

Referenced by initRun().

G4int ZdcShowerLibrary::emPDG
private

Definition at line 57 of file ZdcShowerLibrary.h.

Referenced by encodePartID(), and initRun().

G4int ZdcShowerLibrary::epPDG
private

Definition at line 57 of file ZdcShowerLibrary.h.

Referenced by encodePartID(), and initRun().

G4int ZdcShowerLibrary::etaPDG
private

Definition at line 58 of file ZdcShowerLibrary.h.

Referenced by initRun().

G4int ZdcShowerLibrary::gammaPDG
private

Definition at line 57 of file ZdcShowerLibrary.h.

Referenced by encodePartID(), and initRun().

G4int ZdcShowerLibrary::geantinoPDG
private

Definition at line 59 of file ZdcShowerLibrary.h.

Referenced by getHits(), and initRun().

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

Definition at line 61 of file ZdcShowerLibrary.h.

Referenced by getHits(), and ZdcShowerLibrary().

G4int ZdcShowerLibrary::nuePDG
private

Definition at line 58 of file ZdcShowerLibrary.h.

Referenced by initRun().

G4int ZdcShowerLibrary::numuPDG
private

Definition at line 58 of file ZdcShowerLibrary.h.

Referenced by initRun().

G4int ZdcShowerLibrary::nutauPDG
private

Definition at line 58 of file ZdcShowerLibrary.h.

Referenced by initRun().

G4int ZdcShowerLibrary::pi0PDG
private

Definition at line 58 of file ZdcShowerLibrary.h.

Referenced by initRun().

bool ZdcShowerLibrary::verbose
private

Definition at line 56 of file ZdcShowerLibrary.h.