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 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)
 
int photonFluctuation (double eav, double esig, double edis)
 
 ZdcShowerLibrary (const std::string &name, const DDCompactView &cpv, edm::ParameterSet const &p)
 
 ~ZdcShowerLibrary ()
 

Private Attributes

std::vector< ZdcShowerLibrary::Hit > hits
 
int npe
 
bool verbose
 

Detailed Description

Definition at line 23 of file ZdcShowerLibrary.h.

Constructor & Destructor Documentation

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

Definition at line 20 of file ZdcShowerLibrary.cc.

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

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

Definition at line 29 of file ZdcShowerLibrary.cc.

29  {
30 }

Member Function Documentation

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

Definition at line 142 of file ZdcShowerLibrary.cc.

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

Referenced by getHits().

143  {
144  int nphotons = -1;
145 
146  energy =energy/GeV;
147 
148  LogDebug("ZdcShower")
149  <<"\n ZdcShowerLibrary::getEnergyFromLibrary input/output variables:"
150  <<" phi: "<<59.2956*momDir.phi()
151  <<" theta: "<<59.2956*momDir.theta()
152  <<" xin : "<<hitPoint.x()
153  <<" yin : "<<hitPoint.y()
154  <<" zin : "<<hitPoint.z()
155  <<" track en: " <<energy<< "(GeV)"
156  <<" section: "<<section
157  <<" side: "<<side
158  <<" channel: "<< channel
159  <<" partID: "<<parCode;
160 
161  double eav = 0.;
162  double esig = 0.;
163  double edis = 0.;
164 
165  float xin = hitPoint.x();
166  float yin = hitPoint.y();
167  float fact = 0.;
168 
169  bool isEM = G4TrackToParticleID::isGammaElectronPositron(parCode);
170 
171  if(section == 1 && !isEM){
172  if(channel < 5 )
173  if(((theXChannelBoundaries[channel-1])< (xin + X0)) && ((xin + X0)<= theXChannelBoundaries[channel]))fact= 0.18;
174  if(channel ==5 )
175  if(theXChannelBoundaries[channel-1]< xin + X0)fact = 0.18;
176  }
177 
178  if(section == 2 && !isEM){
179  if(channel == 1)fact = 0.34;
180  if(channel == 2)fact = 0.24;
181  if(channel == 3)fact = 0.17;
182  if(channel == 4)fact = 0.07;
183  }
184  if(section == 1 && isEM){
185  if(channel < 5 )
186  if(((theXChannelBoundaries[channel-1])< (xin + X0)) && ((xin + X0)<= theXChannelBoundaries[channel]))fact= 1.;
187  if(channel ==5 )
188  if(theXChannelBoundaries[channel-1]< xin + X0)fact = 1.0;
189  }
190 
191  //change to cm for parametrization
192  yin = yin/cm;
193  xin = xin/cm;
194 
195  if (isEM){
196  eav = ((((((-0.0002*xin-2.0e-13)*xin+0.0022)*xin+1.0e-11)*xin-0.0217)*xin-3.0e-10)*xin+1.0028)*
197  (((0.0001*yin + 0.0056)*yin + 0.0508)*yin + 1.0)*300.0*pow((energy/300.0),0.99); // EM
198  esig = ((((((0.0005*xin - 1.0e-12)*xin - 0.0052)*xin + 5.0e-11)*xin + 0.032)*xin -
199  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
200  edis = 1.0;
201  }else{
202  eav = ((((((-0.0002*xin-2.0e-13)*xin+0.0022)*xin+1.0e-11)*xin-0.0217)*xin-3.0e-10)*xin+1.0028)*
203  (((0.0001*yin + 0.0056)*yin + 0.0508)*yin + 1.0)*300.0*pow((energy/300.0),1.12);// HD
204  esig = ((((((0.0005*xin - 1.0e-12)*xin - 0.0052)*xin + 5.0e-11)*xin + 0.032)*xin -
205  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
206  edis = 3.0;
207  }
208 
209  if(eav <0. || edis <0.){
210  LogDebug("ZdcShower")
211  <<" Negative everage energy from parametrization \n"
212  <<" xin: "<<xin<< "(cm)"
213  <<" yin: "<<yin<< "(cm)"
214  <<" track en: " <<energy<< "(GeV)"
215  <<" eaverage: "<<eav/GeV << " (GeV)"
216  <<" esigma: "<<esig/GeV << " (GeV)"
217  <<" edist: "<<edis << " (GeV)"
218  <<" dE hit: "<<nphotons/GeV<< " (GeV)";
219  return 0;
220  }
221 
222  // Convert from GeV to MeV for the code
223  eav = eav*GeV;
224  esig= esig*GeV;
225 
226  while(nphotons == -1 || nphotons > int(eav + 5.*esig))
227  nphotons = (int)(fact*photonFluctuation(eav, esig, edis));
228 
229  LogDebug("ZdcShower")
230  //std::cout
231  <<" track en: " <<energy<< "(GeV)"
232  <<" eaverage: "<<eav/GeV << " (GeV)"
233  <<" esigma: "<<esig/GeV << " (GeV)"
234  <<" edist: "<<edis << " (GeV)"
235  <<" dE hit: "<<nphotons/GeV<< " (GeV)";
236 
237  return nphotons;
238 }
LogDebug
#define LogDebug(id)
Definition: MessageLogger.h:670
GeV
const double GeV
Definition: MathUtil.h:16
theXChannelBoundaries
static const double theXChannelBoundaries[]
Definition: CastorGeometryData.h:41
ZdcShowerLibrary::photonFluctuation
int photonFluctuation(double eav, double esig, double edis)
Definition: ZdcShowerLibrary.cc:240
fact
const double fact
Definition: NuclearInteractionFTFSimulator.cc:74
G4TrackToParticleID::isGammaElectronPositron
static bool isGammaElectronPositron(int pdgCode)
Definition: G4TrackToParticleID.cc:17
X0
static const double X0
Definition: CastorGeometryData.h:31
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
std::vector< ZdcShowerLibrary::Hit > & ZdcShowerLibrary::getHits ( const G4Step *  aStep,
bool &  ok 
)

Definition at line 32 of file ZdcShowerLibrary.cc.

References funct::abs(), ZdcShowerLibrary::Hit::DeEM, ZdcShowerLibrary::Hit::DeHad, egammaForCoreTracking_cff::depth, ZdcShowerLibrary::Hit::depth, ZdcShowerLibrary::Hit::detID, HcalZDCDetId::EM, ZdcShowerLibrary::Hit::entryLocal, getEnergyFromLibrary(), HcalZDCDetId::HAD, hits, mps_fire::i, G4TrackToParticleID::isGammaElectronPositron(), G4TrackToParticleID::isStableHadronIon(), LogDebug, npe, ZdcShowerLibrary::Hit::position, trackingPlots::section, theXChannelBoundaries, theZHadChannelBoundaries, theZSectionBoundaries, ZdcShowerLibrary::Hit::time, HiIsolationCommonParameters_cff::track, X0, geometryCSVtoXML::xx, Y0, geometryCSVtoXML::yy, Z0, and geometryCSVtoXML::zz.

32  {
33 
34  const G4StepPoint * preStepPoint = aStep->GetPreStepPoint();
35  const G4StepPoint * postStepPoint = aStep->GetPostStepPoint();
36  const G4Track * track = aStep->GetTrack();
37 
38  const G4DynamicParticle *aParticle = track->GetDynamicParticle();
39  const G4ThreeVector& momDir = aParticle->GetMomentumDirection();
40  double energy = preStepPoint->GetKineticEnergy();
41  G4ThreeVector hitPoint = preStepPoint->GetPosition();
42  const G4ThreeVector& hitPointOrig = preStepPoint->GetPosition();
43  G4int parCode = track->GetDefinition()->GetPDGEncoding();
44 
45  hits.clear();
46 
47  ok = false;
48  bool isEM = G4TrackToParticleID::isGammaElectronPositron(parCode);
49  bool isHad = G4TrackToParticleID::isStableHadronIon(track);
50  if (!isEM && !isHad)
51  return hits;
52  ok = true;
53 
54  G4ThreeVector pos;
55  G4ThreeVector posLocal;
56  double tSlice = (postStepPoint->GetGlobalTime())/nanosecond;
57 
58  int nHit = 0;
59  HcalZDCDetId::Section section;
60  bool side = false;
61  int channel = 0;
62  double xx,yy,zz;
63  double xxlocal, yylocal, zzlocal;
64 
65  ZdcShowerLibrary::Hit oneHit;
66  side = (hitPointOrig.z() > 0.) ? true : false;
67 
68  float xWidthEM = std::abs(theXChannelBoundaries[0] - theXChannelBoundaries[1]);
69  float zWidthEM = std::abs(theZSectionBoundaries[0] - theZSectionBoundaries[1]);
70  float zWidthHAD = std::abs(theZHadChannelBoundaries[0] -theZHadChannelBoundaries[1]);
71 
72  for (int i = 0; i < npe; i++) {
73  if(i < 5){
74  section = HcalZDCDetId::EM;
75  channel = i+1;
76  xxlocal = theXChannelBoundaries[i]+(xWidthEM/2.);
77  xx = xxlocal + X0;
78  yy = 0.0;
79  yylocal = yy + Y0;
80  zzlocal = theZSectionBoundaries[0]+(zWidthEM/2.);
81  zz = (hitPointOrig.z() > 0.) ? zzlocal + Z0 : zzlocal - Z0;
82  pos = G4ThreeVector(xx,yy,zz);
83  posLocal = G4ThreeVector(xxlocal,yylocal,zzlocal);
84  }
85  if(i > 4){
86  section = HcalZDCDetId::HAD;
87  channel = i - 4;
88  xxlocal = 0.0;
89  xx = xxlocal + X0;
90  yylocal = 0;
91  yy = yylocal + Y0;
92  zzlocal = (hitPointOrig.z() > 0.) ?
93  theZHadChannelBoundaries[i-5] + (zWidthHAD/2.) : theZHadChannelBoundaries[i-5] - (zWidthHAD/2.);
94  zz = (hitPointOrig.z() > 0.) ? zzlocal + Z0 : zzlocal - Z0;
95  pos = G4ThreeVector(xx,yy,zz);
96  posLocal = G4ThreeVector(xxlocal,yylocal,zzlocal);
97  }
98 
99  oneHit.position = pos;
100  oneHit.entryLocal = posLocal;
101  oneHit.depth = channel;
102  oneHit.time = tSlice;
103  oneHit.detID = HcalZDCDetId(section,side,channel);
104 
105 
106  // Note: coodinates of hit are relative to center of detector (X0,Y0,Z0)
107  hitPoint.setX(hitPointOrig.x()-X0);
108  hitPoint.setY(hitPointOrig.y()-Y0);
109  double setZ= (hitPointOrig.z()> 0.) ? hitPointOrig.z()- Z0 : fabs(hitPointOrig.z()) - Z0;
110  hitPoint.setZ(setZ);
111 
112  int dE = getEnergyFromLibrary(hitPoint,momDir,energy,parCode,section,side,channel);
113 
114  if (isEM) {
115  oneHit.DeEM = dE;
116  oneHit.DeHad = 0.;
117  } else {
118  oneHit.DeEM = 0;
119  oneHit.DeHad = dE;
120  }
121 
122  hits.push_back(oneHit);
123 
124  LogDebug("ZdcShower")
125  << "\nZdcShowerLibrary:Generated Hit " << nHit
126  <<" orig hit pos " << hitPointOrig
127  <<" orig hit pos local coord" << hitPoint
128  <<" new position " << (hits[nHit].position)
129  <<" Channel " <<(hits[nHit].depth)
130  <<" side "<< side
131  <<" Time " <<(hits[nHit].time)
132  <<" DetectorID " << (hits[nHit].detID)
133  <<" Had Energy " << (hits[nHit].DeHad)
134  <<" EM Energy " << (hits[nHit].DeEM)
135  <<"\n";
136  nHit++;
137  }
138  return hits;
139 }
LogDebug
#define LogDebug(id)
Definition: MessageLogger.h:670
mps_fire.i
i
Definition: mps_fire.py:338
ZdcShowerLibrary::Hit::entryLocal
G4ThreeVector entryLocal
Definition: ZdcShowerLibrary.h:35
Z0
static const double Z0
Definition: CastorGeometryData.h:32
ZdcShowerLibrary::getEnergyFromLibrary
int getEnergyFromLibrary(const G4ThreeVector &posHit, const G4ThreeVector &momDir, double energy, G4int parCode, HcalZDCDetId::Section section, bool side, int channel)
Definition: ZdcShowerLibrary.cc:142
theXChannelBoundaries
static const double theXChannelBoundaries[]
Definition: CastorGeometryData.h:41
ZdcShowerLibrary::hits
std::vector< ZdcShowerLibrary::Hit > hits
Definition: ZdcShowerLibrary.h:54
theZHadChannelBoundaries
static const double theZHadChannelBoundaries[]
Definition: CastorGeometryData.h:42
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
theZSectionBoundaries
static const double theZSectionBoundaries[]
Definition: CastorGeometryData.h:40
G4TrackToParticleID::isStableHadronIon
static bool isStableHadronIon(const G4Track *)
Definition: G4TrackToParticleID.cc:40
HiIsolationCommonParameters_cff.track
track
Definition: HiIsolationCommonParameters_cff.py:9
G4TrackToParticleID::isGammaElectronPositron
static bool isGammaElectronPositron(int pdgCode)
Definition: G4TrackToParticleID.cc:17
Y0
static const double Y0
Definition: CastorGeometryData.h:33
X0
static const double X0
Definition: CastorGeometryData.h:31
egammaForCoreTracking_cff.depth
depth
Definition: egammaForCoreTracking_cff.py:29
int ZdcShowerLibrary::photonFluctuation ( double  eav,
double  esig,
double  edis 
)

Definition at line 240 of file ZdcShowerLibrary.cc.

References createfilelist::int.

Referenced by getEnergyFromLibrary().

240  {
241  int nphot=0;
242  double efluct = 0.;
243  if(edis == 1.0) efluct = eav+esig*CLHEP::RandGaussQ::shoot();
244  if(edis == 3.0) efluct = eav+esig*CLHEP::RandLandau::shoot();
245  nphot = int(efluct);
246  if(nphot<0)nphot = 0;
247  return nphot;
248 }

Member Data Documentation

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

Definition at line 54 of file ZdcShowerLibrary.h.

Referenced by getHits(), and ZdcShowerLibrary().

int ZdcShowerLibrary::npe
private

Definition at line 53 of file ZdcShowerLibrary.h.

Referenced by getHits(), and ZdcShowerLibrary().

bool ZdcShowerLibrary::verbose
private

Definition at line 51 of file ZdcShowerLibrary.h.

Generated for CMSSW Reference Manual by doxygen 1.8.11