CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Protected Member Functions | Private Types | Private Member Functions | Private Attributes
DreamSD Class Reference

#include <DreamSD.h>

Inheritance diagram for DreamSD:
CaloSD SensitiveCaloDetector Observer< const BeginOfRun * > Observer< const BeginOfEvent * > Observer< const BeginOfTrack * > Observer< const EndOfTrack * > Observer< const EndOfEvent * > SensitiveDetector

Public Member Functions

 DreamSD (const std::string &, const edm::EventSetup &, const SensitiveDetectorCatalog &, edm::ParameterSet const &, const SimTrackManager *)
 
uint32_t setDetUnitId (const G4Step *) override
 
 ~DreamSD () override
 
- Public Member Functions inherited from CaloSD
 CaloSD (const std::string &aSDname, const edm::EventSetup &es, const SensitiveDetectorCatalog &clg, edm::ParameterSet const &p, const SimTrackManager *, float timeSlice=1., bool ignoreTkID=false)
 
void clear () override
 
void clearHits () override
 
void DrawAll () override
 
void EndOfEvent (G4HCofThisEvent *eventHC) override
 
void fillHits (edm::PCaloHitContainer &, const std::string &) override
 
void Initialize (G4HCofThisEvent *HCE) override
 
void PrintAll () override
 
G4bool ProcessHits (G4Step *step, G4TouchableHistory *) override
 
bool ProcessHits (G4GFlashSpot *aSpot, G4TouchableHistory *) override
 
void reset () override
 
 ~CaloSD () override
 
- Public Member Functions inherited from SensitiveCaloDetector
 SensitiveCaloDetector (const std::string &iname, const edm::EventSetup &es, const SensitiveDetectorCatalog &clg, edm::ParameterSet const &p)
 
- Public Member Functions inherited from SensitiveDetector
void EndOfEvent (G4HCofThisEvent *eventHC) override
 
const std::vector< std::string > & getNames () const
 
void Initialize (G4HCofThisEvent *eventHC) override
 
bool isCaloSD () const
 
 SensitiveDetector (const std::string &iname, const edm::EventSetup &es, const SensitiveDetectorCatalog &, edm::ParameterSet const &p, bool calo)
 
 ~SensitiveDetector () override
 
- Public Member Functions inherited from Observer< const BeginOfRun * >
 Observer ()
 
void slotForUpdate (const BeginOfRun * iT)
 
virtual ~Observer ()
 
- Public Member Functions inherited from Observer< const BeginOfEvent * >
 Observer ()
 
void slotForUpdate (const BeginOfEvent * iT)
 
virtual ~Observer ()
 
- Public Member Functions inherited from Observer< const BeginOfTrack * >
 Observer ()
 
void slotForUpdate (const BeginOfTrack * iT)
 
virtual ~Observer ()
 
- Public Member Functions inherited from Observer< const EndOfTrack * >
 Observer ()
 
void slotForUpdate (const EndOfTrack * iT)
 
virtual ~Observer ()
 
- Public Member Functions inherited from Observer< const EndOfEvent * >
 Observer ()
 
void slotForUpdate (const EndOfEvent * iT)
 
virtual ~Observer ()
 

Protected Member Functions

double getEnergyDeposit (const G4Step *) override
 
void initRun () override
 
- Protected Member Functions inherited from CaloSD
bool checkHit ()
 
CaloG4HitcreateNewHit (const G4Step *, const G4Track *)
 
virtual void endEvent ()
 
virtual bool filterHit (CaloG4Hit *, double)
 
double getAttenuation (const G4Step *aStep, double birk1, double birk2, double birk3) const
 
virtual uint16_t getDepth (const G4Step *)
 
virtual bool getFromLibrary (const G4Step *step)
 
int getNumberOfHits ()
 
double getResponseWt (const G4Track *)
 
virtual int getTrackID (const G4Track *)
 
bool hitExists (const G4Step *)
 
virtual void initEvent (const BeginOfEvent *)
 
void processHit (const G4Step *step)
 
void resetForNewPrimary (const G4Step *)
 
void setNumberCheckedHits (int val)
 
void setParameterized (bool val)
 
G4ThreeVector setToGlobal (const G4ThreeVector &, const G4VTouchable *) const
 
G4ThreeVector setToLocal (const G4ThreeVector &, const G4VTouchable *) const
 
virtual int setTrackID (const G4Step *)
 
void setUseMap (bool val)
 
void update (const BeginOfRun *) override
 This routine will be called when the appropriate signal arrives. More...
 
void update (const BeginOfEvent *) override
 This routine will be called when the appropriate signal arrives. More...
 
void update (const BeginOfTrack *trk) override
 This routine will be called when the appropriate signal arrives. More...
 
void update (const EndOfTrack *trk) override
 This routine will be called when the appropriate signal arrives. More...
 
void update (const ::EndOfEvent *) override
 
void updateHit (CaloG4Hit *)
 
- Protected Member Functions inherited from SensitiveDetector
TrackInformationcmsTrackInformation (const G4Track *aTrack)
 
Local3DPoint ConvertToLocal3DPoint (const G4ThreeVector &point) const
 
Local3DPoint FinalStepPosition (const G4Step *step, coordinates) const
 
Local3DPoint InitialStepPosition (const G4Step *step, coordinates) const
 
Local3DPoint LocalPostStepPosition (const G4Step *step) const
 
Local3DPoint LocalPreStepPosition (const G4Step *step) const
 
void NaNTrap (const G4Step *step) const
 
void setNames (const std::vector< std::string > &)
 
- Protected Member Functions inherited from Observer< const EndOfEvent * >
virtual void update (const EndOfEvent *)=0
 This routine will be called when the appropriate signal arrives. More...
 

Private Types

typedef std::map< G4LogicalVolume *, DoublesDimensionMap
 
typedef std::pair< double, double > Doubles
 

Private Member Functions

double cherenkovDeposit_ (const G4Step *aStep)
 Returns the total energy due to Cherenkov radiation. More...
 
double crystalLength (G4LogicalVolume *) const
 
double crystalWidth (G4LogicalVolume *) const
 
double curve_LY (const G4Step *, int)
 
double getAverageNumberOfPhotons_ (const double charge, const double beta, const G4Material *aMaterial, const G4MaterialPropertyVector *rIndex)
 Returns average number of photons created by track. More...
 
double getPhotonEnergyDeposit_ (const G4ParticleMomentum &p, const G4ThreeVector &x, const G4Step *aStep)
 Returns energy deposit for a given photon. More...
 
void initMap (const std::string &, const edm::EventSetup &)
 
bool setPbWO2MaterialProperties_ (G4Material *aMaterial)
 Sets material properties at run-time... More...
 

Private Attributes

double birk1
 
double birk2
 
double birk3
 
std::unique_ptr< G4PhysicsOrderedFreeVector > chAngleIntegrals_
 Table of Cherenkov angle integrals vs photon momentum. More...
 
bool doCherenkov_
 
G4MaterialPropertiesTable * materialPropertiesTable
 
int nphotons_
 
bool readBothSide_
 
int side
 
double slopeLY
 
bool useBirk
 
DimensionMap xtalLMap
 

Additional Inherited Members

- Protected Types inherited from SensitiveDetector
enum  coordinates { WorldCoordinates, LocalCoordinates }
 
- Protected Attributes inherited from CaloSD
CaloG4HitcurrentHit
 
CaloHitID currentID
 
float edepositEM
 
float edepositHAD
 
double eminHit
 
double energyCut
 
G4ThreeVector entranceLocal
 
G4ThreeVector entrancePoint
 
bool forceSave
 
float incidentEnergy
 
double kmaxIon
 
double kmaxNeutron
 
double kmaxProton
 
G4ThreeVector posGlobal
 
CaloHitID previousID
 
bool suppressHeavy
 
double tmaxHit
 

Detailed Description

Definition at line 15 of file DreamSD.h.

Member Typedef Documentation

typedef std::map<G4LogicalVolume *, Doubles> DreamSD::DimensionMap
private

Definition at line 32 of file DreamSD.h.

typedef std::pair<double, double> DreamSD::Doubles
private

Definition at line 31 of file DreamSD.h.

Constructor & Destructor Documentation

DreamSD::DreamSD ( const std::string &  name,
const edm::EventSetup es,
const SensitiveDetectorCatalog clg,
edm::ParameterSet const &  p,
const SimTrackManager manager 
)

Definition at line 31 of file DreamSD.cc.

References birk1, birk2, birk3, chAngleIntegrals_, doCherenkov_, g, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), initMap(), MeV, readBothSide_, slopeLY, and useBirk.

36  : CaloSD(name, es, clg, p, manager) {
37  edm::ParameterSet m_EC = p.getParameter<edm::ParameterSet>("ECalSD");
38  useBirk = m_EC.getParameter<bool>("UseBirkLaw");
39  doCherenkov_ = m_EC.getParameter<bool>("doCherenkov");
40  birk1 = m_EC.getParameter<double>("BirkC1") * (g / (MeV * cm2));
41  birk2 = m_EC.getParameter<double>("BirkC2");
42  birk3 = m_EC.getParameter<double>("BirkC3");
43  slopeLY = m_EC.getParameter<double>("SlopeLightYield");
44  readBothSide_ = m_EC.getUntrackedParameter<bool>("ReadBothSide", false);
45 
46  chAngleIntegrals_.reset(nullptr);
47 
48  edm::LogInfo("EcalSim") << "Constructing a DreamSD with name " << GetName() << "\n"
49  << "DreamSD:: Use of Birks law is set to " << useBirk
50  << " with three constants kB = " << birk1 << ", C1 = " << birk2 << ", C2 = " << birk3 << "\n"
51  << " Slope for Light yield is set to " << slopeLY << "\n"
52  << " Parameterization of Cherenkov is set to " << doCherenkov_
53  << " and readout both sides is " << readBothSide_;
54 
55  initMap(name, es);
56 }
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
bool readBothSide_
Definition: DreamSD.h:51
double birk1
Definition: DreamSD.h:52
void initMap(const std::string &, const edm::EventSetup &)
Definition: DreamSD.cc:108
double birk3
Definition: DreamSD.h:52
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e g
Definition: Activities.doc:4
const double MeV
CaloSD(const std::string &aSDname, const edm::EventSetup &es, const SensitiveDetectorCatalog &clg, edm::ParameterSet const &p, const SimTrackManager *, float timeSlice=1., bool ignoreTkID=false)
Definition: CaloSD.cc:27
double birk2
Definition: DreamSD.h:52
std::unique_ptr< G4PhysicsOrderedFreeVector > chAngleIntegrals_
Table of Cherenkov angle integrals vs photon momentum.
Definition: DreamSD.h:59
bool doCherenkov_
Definition: DreamSD.h:51
double slopeLY
Definition: DreamSD.h:53
bool useBirk
Definition: DreamSD.h:51
DreamSD::~DreamSD ( )
inlineoverride

Definition at line 22 of file DreamSD.h.

References getEnergyDeposit(), initRun(), and setDetUnitId().

22 {}

Member Function Documentation

double DreamSD::cherenkovDeposit_ ( const G4Step *  aStep)
private

Returns the total energy due to Cherenkov radiation.

Definition at line 198 of file DreamSD.cc.

References zMuMuMuonUserData::beta, ALCARECOTkAlJpsiMuMu_cff::charge, funct::cos(), Calorimetry_cff::dp, getAverageNumberOfPhotons_(), getPhotonEnergyDeposit_(), LogDebug, materialPropertiesTable, phi, dataAnalyzerFineBiningParameters_cff::Pmax, multPhiCorr_741_25nsDY_cfi::px, multPhiCorr_741_25nsDY_cfi::py, funct::sin(), and mathSSE::sqrt().

Referenced by getEnergyDeposit().

198  {
199  double cherenkovEnergy = 0;
201  return cherenkovEnergy;
202  G4Material *material = aStep->GetTrack()->GetMaterial();
203 
204  // Retrieve refractive index
205  G4MaterialPropertyVector *Rindex = materialPropertiesTable->GetProperty("RINDEX");
206  if (Rindex == nullptr) {
207  edm::LogWarning("EcalSim") << "Couldn't retrieve refractive index";
208  return cherenkovEnergy;
209  }
210 
211  // V.Ivanchenko - temporary close log output for 9.5
212  // Material refraction properties
213  int Rlength = Rindex->GetVectorLength() - 1;
214  double Pmin = Rindex->Energy(0);
215  double Pmax = Rindex->Energy(Rlength);
216  LogDebug("EcalSim") << "Material properties: "
217  << "\n"
218  << " Pmin = " << Pmin << " Pmax = " << Pmax;
219 
220  // Get particle properties
221  const G4StepPoint *pPreStepPoint = aStep->GetPreStepPoint();
222  const G4StepPoint *pPostStepPoint = aStep->GetPostStepPoint();
223  const G4ThreeVector &x0 = pPreStepPoint->GetPosition();
224  G4ThreeVector p0 = aStep->GetDeltaPosition().unit();
225  const G4DynamicParticle *aParticle = aStep->GetTrack()->GetDynamicParticle();
226  const double charge = aParticle->GetDefinition()->GetPDGCharge();
227  // beta is averaged over step
228  double beta = 0.5 * (pPreStepPoint->GetBeta() + pPostStepPoint->GetBeta());
229  double BetaInverse = 1.0 / beta;
230 
231  LogDebug("EcalSim") << "Particle properties: "
232  << "\n"
233  << " charge = " << charge << " beta = " << beta;
234 
235  // Now get number of photons generated in this step
236  double meanNumberOfPhotons = getAverageNumberOfPhotons_(charge, beta, material, Rindex);
237  if (meanNumberOfPhotons <= 0.0) { // Don't do anything
238  LogDebug("EcalSim") << "Mean number of photons is zero: " << meanNumberOfPhotons << ", stopping here";
239  return cherenkovEnergy;
240  }
241 
242  // number of photons is in unit of Geant4...
243  meanNumberOfPhotons *= aStep->GetStepLength();
244 
245  // Now get a poisson distribution
246  int numPhotons = static_cast<int>(G4Poisson(meanNumberOfPhotons));
247  // edm::LogVerbatim("EcalSim") << "Number of photons = " << numPhotons;
248  if (numPhotons <= 0) {
249  LogDebug("EcalSim") << "Poission number of photons is zero: " << numPhotons << ", stopping here";
250  return cherenkovEnergy;
251  }
252 
253  // Material refraction properties
254  double dp = Pmax - Pmin;
255  double maxCos = BetaInverse / (*Rindex)[Rlength];
256  double maxSin2 = (1.0 - maxCos) * (1.0 + maxCos);
257 
258  // Finally: get contribution of each photon
259  for (int iPhoton = 0; iPhoton < numPhotons; ++iPhoton) {
260  // Determine photon momentum
261  double randomNumber;
262  double sampledMomentum, sampledRI;
263  double cosTheta, sin2Theta;
264 
265  // sample a momentum (not sure why this is needed!)
266  do {
267  randomNumber = G4UniformRand();
268  sampledMomentum = Pmin + randomNumber * dp;
269  sampledRI = Rindex->Value(sampledMomentum);
270  cosTheta = BetaInverse / sampledRI;
271 
272  sin2Theta = (1.0 - cosTheta) * (1.0 + cosTheta);
273  randomNumber = G4UniformRand();
274 
275  } while (randomNumber * maxSin2 > sin2Theta);
276 
277  // Generate random position of photon on cone surface
278  // defined by Theta
279  randomNumber = G4UniformRand();
280 
281  double phi = twopi * randomNumber;
282  double sinPhi = sin(phi);
283  double cosPhi = cos(phi);
284 
285  // Create photon momentum direction vector
286  // The momentum direction is still w.r.t. the coordinate system where the
287  // primary particle direction is aligned with the z axis
288  double sinTheta = sqrt(sin2Theta);
289  double px = sinTheta * cosPhi;
290  double py = sinTheta * sinPhi;
291  double pz = cosTheta;
292  G4ThreeVector photonDirection(px, py, pz);
293 
294  // Rotate momentum direction back to global (crystal) reference system
295  photonDirection.rotateUz(p0);
296 
297  // Create photon position and momentum
298  randomNumber = G4UniformRand();
299  G4ThreeVector photonPosition = x0 + randomNumber * aStep->GetDeltaPosition();
300  G4ThreeVector photonMomentum = sampledMomentum * photonDirection;
301 
302  // Collect energy on APD
303  cherenkovEnergy += getPhotonEnergyDeposit_(photonMomentum, photonPosition, aStep);
304  }
305  return cherenkovEnergy;
306 }
#define LogDebug(id)
double getPhotonEnergyDeposit_(const G4ParticleMomentum &p, const G4ThreeVector &x, const G4Step *aStep)
Returns energy deposit for a given photon.
Definition: DreamSD.cc:458
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
T sqrt(T t)
Definition: SSEVec.h:19
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
double getAverageNumberOfPhotons_(const double charge, const double beta, const G4Material *aMaterial, const G4MaterialPropertyVector *rIndex)
Returns average number of photons created by track.
Definition: DreamSD.cc:311
G4MaterialPropertiesTable * materialPropertiesTable
Definition: DreamSD.h:60
double DreamSD::crystalLength ( G4LogicalVolume *  lv) const
private

Definition at line 178 of file DreamSD.cc.

References xtalLMap.

Referenced by curve_LY(), and getPhotonEnergyDeposit_().

178  {
179  double length = -1.;
180  DimensionMap::const_iterator ite = xtalLMap.find(lv);
181  if (ite != xtalLMap.end())
182  length = ite->second.first;
183  return length;
184 }
DimensionMap xtalLMap
Definition: DreamSD.h:54
double DreamSD::crystalWidth ( G4LogicalVolume *  lv) const
private

Definition at line 187 of file DreamSD.cc.

References ApeEstimator_cff::width, and xtalLMap.

Referenced by getPhotonEnergyDeposit_().

187  {
188  double width = -1.;
189  DimensionMap::const_iterator ite = xtalLMap.find(lv);
190  if (ite != xtalLMap.end())
191  width = ite->second.second;
192  return width;
193 }
DimensionMap xtalLMap
Definition: DreamSD.h:54
double DreamSD::curve_LY ( const G4Step *  aStep,
int  flag 
)
private

Definition at line 153 of file DreamSD.cc.

References crystalLength(), LogDebug, CaloSD::setToLocal(), slopeLY, and mps_merge::weight.

Referenced by getEnergyDeposit().

153  {
154  auto const stepPoint = aStep->GetPreStepPoint();
155  auto const lv = stepPoint->GetTouchable()->GetVolume(0)->GetLogicalVolume();
156  G4String nameVolume = lv->GetName();
157 
158  double weight = 1.;
159  G4ThreeVector localPoint = setToLocal(stepPoint->GetPosition(), stepPoint->GetTouchable());
160  double crlength = crystalLength(lv);
161  double localz = localPoint.x();
162  double dapd = 0.5 * crlength - flag * localz; // Distance from closest APD
163  if (dapd >= -0.1 || dapd <= crlength + 0.1) {
164  if (dapd <= 100.)
165  weight = 1.0 + slopeLY - dapd * 0.01 * slopeLY;
166  } else {
167  edm::LogWarning("EcalSim") << "DreamSD: light coll curve : wrong distance "
168  << "to APD " << dapd << " crlength = " << crlength << " crystal name = " << nameVolume
169  << " z of localPoint = " << localz << " take weight = " << weight;
170  }
171  LogDebug("EcalSim") << "DreamSD, light coll curve : " << dapd << " crlength = " << crlength
172  << " crystal name = " << nameVolume << " z of localPoint = " << localz
173  << " take weight = " << weight;
174  return weight;
175 }
#define LogDebug(id)
Definition: weight.py:1
G4ThreeVector setToLocal(const G4ThreeVector &, const G4VTouchable *) const
Definition: CaloSD.cc:291
double crystalLength(G4LogicalVolume *) const
Definition: DreamSD.cc:178
double slopeLY
Definition: DreamSD.h:53
double DreamSD::getAverageNumberOfPhotons_ ( const double  charge,
const double  beta,
const G4Material *  aMaterial,
const G4MaterialPropertyVector *  rIndex 
)
private

Returns average number of photons created by track.

Definition at line 311 of file DreamSD.cc.

References zMuMuMuonUserData::beta, chAngleIntegrals_, Calorimetry_cff::dp, LogDebug, metDiagnosticParameterSet_cfi::nMax, metDiagnosticParameterSet_cfi::nMin, and dataAnalyzerFineBiningParameters_cff::Pmax.

Referenced by cherenkovDeposit_().

314  {
315  const double rFact = 369.81 / (eV * cm);
316 
317  if (beta <= 0.0)
318  return 0.0;
319 
320  double BetaInverse = 1. / beta;
321 
322  // Vectors used in computation of Cerenkov Angle Integral:
323  // - Refraction Indices for the current material
324  // - new G4PhysicsOrderedFreeVector allocated to hold CAI's
325 
326  // Min and Max photon momenta
327  int Rlength = Rindex->GetVectorLength() - 1;
328  double Pmin = Rindex->Energy(0);
329  double Pmax = Rindex->Energy(Rlength);
330 
331  // Min and Max Refraction Indices
332  double nMin = (*Rindex)[0];
333  double nMax = (*Rindex)[Rlength];
334 
335  // Max Cerenkov Angle Integral
336  double CAImax = chAngleIntegrals_.get()->GetMaxValue();
337 
338  double dp = 0., ge = 0., CAImin = 0.;
339 
340  // If n(Pmax) < 1/Beta -- no photons generated
341  if (nMax < BetaInverse) {
342  }
343 
344  // otherwise if n(Pmin) >= 1/Beta -- photons generated
345  else if (nMin > BetaInverse) {
346  dp = Pmax - Pmin;
347  ge = CAImax;
348  }
349  // If n(Pmin) < 1/Beta, and n(Pmax) >= 1/Beta, then
350  // we need to find a P such that the value of n(P) == 1/Beta.
351  // Interpolation is performed by the GetPhotonEnergy() and
352  // GetProperty() methods of the G4MaterialPropertiesTable and
353  // the GetValue() method of G4PhysicsVector.
354  else {
355  Pmin = Rindex->Value(BetaInverse);
356  dp = Pmax - Pmin;
357  // need boolean for current implementation of G4PhysicsVector
358  // ==> being phased out
359  double CAImin = chAngleIntegrals_->Value(Pmin);
360  ge = CAImax - CAImin;
361  }
362 
363  // Calculate number of photons
364  double numPhotons = rFact * charge / eplus * charge / eplus * (dp - ge * BetaInverse * BetaInverse);
365 
366  LogDebug("EcalSim") << "@SUB=getAverageNumberOfPhotons"
367  << "CAImin = " << CAImin << "\n"
368  << "CAImax = " << CAImax << "\n"
369  << "dp = " << dp << ", ge = " << ge << "\n"
370  << "numPhotons = " << numPhotons;
371 
372  return numPhotons;
373 }
#define LogDebug(id)
std::unique_ptr< G4PhysicsOrderedFreeVector > chAngleIntegrals_
Table of Cherenkov angle integrals vs photon momentum.
Definition: DreamSD.h:59
double DreamSD::getEnergyDeposit ( const G4Step *  aStep)
overrideprotectedvirtual

Reimplemented from CaloSD.

Definition at line 59 of file DreamSD.cc.

References birk1, birk2, birk3, cherenkovDeposit_(), curve_LY(), doCherenkov_, CaloSD::getAttenuation(), LogDebug, MeV, side, useBirk, and mps_merge::weight.

Referenced by ~DreamSD().

59  {
60  // take into account light collection curve for crystals
61  double weight = curve_LY(aStep, side);
62  if (useBirk)
63  weight *= getAttenuation(aStep, birk1, birk2, birk3);
64  double edep = aStep->GetTotalEnergyDeposit() * weight;
65 
66  // Get Cerenkov contribution
67  if (doCherenkov_) {
68  edep += cherenkovDeposit_(aStep);
69  }
70  LogDebug("EcalSim") << "DreamSD:: " << aStep->GetPreStepPoint()->GetPhysicalVolume()->GetName() << " Side " << side
71  << " Light Collection Efficiency " << weight << " Weighted Energy Deposit " << edep / MeV
72  << " MeV";
73 
74  return edep;
75 }
#define LogDebug(id)
Definition: weight.py:1
double birk1
Definition: DreamSD.h:52
double birk3
Definition: DreamSD.h:52
const double MeV
double cherenkovDeposit_(const G4Step *aStep)
Returns the total energy due to Cherenkov radiation.
Definition: DreamSD.cc:198
int side
Definition: DreamSD.h:56
double birk2
Definition: DreamSD.h:52
double curve_LY(const G4Step *, int)
Definition: DreamSD.cc:153
double getAttenuation(const G4Step *aStep, double birk1, double birk2, double birk3) const
Definition: CaloSD.cc:434
bool doCherenkov_
Definition: DreamSD.h:51
bool useBirk
Definition: DreamSD.h:51
double DreamSD::getPhotonEnergyDeposit_ ( const G4ParticleMomentum &  p,
const G4ThreeVector &  x,
const G4Step *  aStep 
)
private

Returns energy deposit for a given photon.

Definition at line 458 of file DreamSD.cc.

References funct::abs(), crystalLength(), crystalWidth(), HCALHighEnergyHPDFilter_cfi::energy, PMTResponse::getEfficiency(), LogDebug, and y.

Referenced by cherenkovDeposit_().

458  {
459  double energy = 0;
460 
461  // Crystal dimensions
462 
463  // edm::LogVerbatim("EcalSim") << p << x;
464 
465  // 1. Check if this photon goes straight to the APD:
466  // - assume that APD is at x=xtalLength/2.0
467  // - extrapolate from x=x0 to x=xtalLength/2.0 using momentum in x-y
468 
469  G4StepPoint *stepPoint = aStep->GetPreStepPoint();
470  G4LogicalVolume *lv = stepPoint->GetTouchable()->GetVolume(0)->GetLogicalVolume();
471  G4String nameVolume = lv->GetName();
472 
473  double crlength = crystalLength(lv);
474  double crwidth = crystalWidth(lv);
475  double dapd = 0.5 * crlength - x.x(); // Distance from closest APD
476  double y = p.y() / p.x() * dapd;
477 
478  LogDebug("EcalSim") << "Distance to APD: " << dapd << " - y at APD: " << y;
479 
480  // Not straight: compute probability
481  if (std::abs(y) > crwidth * 0.5) {
482  }
483 
484  // 2. Retrieve efficiency for this wavelength (in nm, from MeV)
485  double waveLength = p.mag() * 1.239e8;
486 
487  energy = p.mag() * PMTResponse::getEfficiency(waveLength);
488 
489  LogDebug("EcalSim") << "Wavelength: " << waveLength << " - Energy: " << energy;
490 
491  return energy;
492 }
#define LogDebug(id)
static double getEfficiency(const double &waveLengthNm)
Return efficiency for given photon wavelength (in nm)
Definition: PMTResponse.cc:6
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
double crystalWidth(G4LogicalVolume *) const
Definition: DreamSD.cc:187
double crystalLength(G4LogicalVolume *) const
Definition: DreamSD.cc:178
void DreamSD::initMap ( const std::string &  sd,
const edm::EventSetup es 
)
private

Definition at line 108 of file DreamSD.cc.

References ALCARECOTkAlBeamHalo_cff::filter, DDFilteredView::firstChild(), edm::EventSetup::get(), mps_fire::i, LogDebug, DDFilteredView::logicalPart(), Skims_PA_cff::name, DDName::name(), DDBase< N, C >::name(), DDFilteredView::next(), DDSolid::parameters(), DDSolid::shape(), DDLogicalPart::solid(), ApeEstimator_cff::width, and xtalLMap.

Referenced by DreamSD().

108  {
110  es.get<IdealGeometryRecord>().get(cpv);
111 
112  G4String attribute = "ReadOutName";
114  DDFilteredView fv(*cpv, filter);
115  fv.firstChild();
116 
117  const G4LogicalVolumeStore *lvs = G4LogicalVolumeStore::GetInstance();
118  std::vector<G4LogicalVolume *>::const_iterator lvcite;
119  bool dodet = true;
120  while (dodet) {
121  const DDSolid &sol = fv.logicalPart().solid();
122  std::vector<double> paras(sol.parameters());
123  G4String name = sol.name().name();
124  G4LogicalVolume *lv = nullptr;
125  for (lvcite = lvs->begin(); lvcite != lvs->end(); lvcite++)
126  if ((*lvcite)->GetName() == name) {
127  lv = (*lvcite);
128  break;
129  }
130  LogDebug("EcalSim") << "DreamSD::initMap (for " << sd << "): Solid " << name << " Shape " << sol.shape()
131  << " Parameter 0 = " << paras[0] << " Logical Volume " << lv;
132  double length = 0, width = 0;
133  // Set length to be the largest size, width the smallest
134  std::sort(paras.begin(), paras.end());
135  length = 2.0 * paras.back();
136  width = 2.0 * paras.front();
137  xtalLMap.insert(std::pair<G4LogicalVolume *, Doubles>(lv, Doubles(length, width)));
138  dodet = fv.next();
139  }
140  LogDebug("EcalSim") << "DreamSD: Length Table for " << attribute << " = " << sd << ":";
141  DimensionMap::const_iterator ite = xtalLMap.begin();
142  int i = 0;
143  for (; ite != xtalLMap.end(); ite++, i++) {
144  G4String name = "Unknown";
145  if (ite->first != nullptr)
146  name = (ite->first)->GetName();
147  LogDebug("EcalSim") << " " << i << " " << ite->first << " " << name << " L = " << ite->second.first
148  << " W = " << ite->second.second;
149  }
150 }
#define LogDebug(id)
const std::vector< double > & parameters(void) const
Give the parameters of the solid.
Definition: DDSolid.cc:121
const N & name() const
Definition: DDBase.h:59
A DDSolid represents the shape of a part.
Definition: DDSolid.h:39
DDSolidShape shape(void) const
The type of the solid.
Definition: DDSolid.cc:119
std::pair< double, double > Doubles
Definition: DreamSD.h:31
DimensionMap xtalLMap
Definition: DreamSD.h:54
double sd
T get() const
Definition: EventSetup.h:73
const std::string & name() const
Returns the name.
Definition: DDName.cc:40
void DreamSD::initRun ( )
overrideprotectedvirtual

Reimplemented from CaloSD.

Definition at line 78 of file DreamSD.cc.

References materialPropertiesTable, setPbWO2MaterialProperties_(), and xtalLMap.

Referenced by ~DreamSD().

78  {
79  // Get the material and set properties if needed
80  DimensionMap::const_iterator ite = xtalLMap.begin();
81  const G4LogicalVolume *lv = (ite->first);
82  G4Material *material = lv->GetMaterial();
83  edm::LogInfo("EcalSim") << "DreamSD::initRun: Initializes for material " << material->GetName() << " in "
84  << lv->GetName();
85  materialPropertiesTable = material->GetMaterialPropertiesTable();
87  if (!setPbWO2MaterialProperties_(material)) {
88  edm::LogWarning("EcalSim") << "Couldn't retrieve material properties table\n"
89  << " Material = " << material->GetName();
90  }
91  materialPropertiesTable = material->GetMaterialPropertiesTable();
92  }
93 }
bool setPbWO2MaterialProperties_(G4Material *aMaterial)
Sets material properties at run-time...
Definition: DreamSD.cc:378
G4MaterialPropertiesTable * materialPropertiesTable
Definition: DreamSD.h:60
DimensionMap xtalLMap
Definition: DreamSD.h:54
uint32_t DreamSD::setDetUnitId ( const G4Step *  aStep)
overridevirtual

Implements CaloSD.

Definition at line 96 of file DreamSD.cc.

References triggerObjects_cff::id, LogDebug, readBothSide_, and side.

Referenced by ~DreamSD().

96  {
97  const G4VTouchable *touch = aStep->GetPreStepPoint()->GetTouchable();
98  uint32_t id = (touch->GetReplicaNumber(1)) * 10 + (touch->GetReplicaNumber(0));
99  side = readBothSide_ ? -1 : 1;
100  if (side < 0) {
101  ++id;
102  }
103  LogDebug("EcalSim") << "DreamSD:: ID " << id;
104  return id;
105 }
#define LogDebug(id)
bool readBothSide_
Definition: DreamSD.h:51
int side
Definition: DreamSD.h:56
bool DreamSD::setPbWO2MaterialProperties_ ( G4Material *  aMaterial)
private

Sets material properties at run-time...

Definition at line 378 of file DreamSD.cc.

References chAngleIntegrals_, LogDebug, AlCaHLTBitMon_QueryRunRegistry::string, and TableParser::table.

Referenced by initRun().

378  {
379  std::string pbWO2Name("E_PbWO4");
380  if (pbWO2Name != aMaterial->GetName()) { // Wrong material!
381  edm::LogWarning("EcalSim") << "This is not the right material: "
382  << "expecting " << pbWO2Name << ", got " << aMaterial->GetName();
383  return false;
384  }
385 
386  G4MaterialPropertiesTable *table = new G4MaterialPropertiesTable();
387 
388  // Refractive index as a function of photon momentum
389  // FIXME: Should somehow put that in the configuration
390  const int nEntries = 14;
391  double PhotonEnergy[nEntries] = {1.7712 * eV,
392  1.8368 * eV,
393  1.90745 * eV,
394  1.98375 * eV,
395  2.0664 * eV,
396  2.15625 * eV,
397  2.25426 * eV,
398  2.3616 * eV,
399  2.47968 * eV,
400  2.61019 * eV,
401  2.75521 * eV,
402  2.91728 * eV,
403  3.09961 * eV,
404  3.30625 * eV};
405  double RefractiveIndex[nEntries] = {2.17728,
406  2.18025,
407  2.18357,
408  2.18753,
409  2.19285,
410  2.19813,
411  2.20441,
412  2.21337,
413  2.22328,
414  2.23619,
415  2.25203,
416  2.27381,
417  2.30282,
418  2.34666};
419 
420  table->AddProperty("RINDEX", PhotonEnergy, RefractiveIndex, nEntries);
421  aMaterial->SetMaterialPropertiesTable(table); // FIXME: could this leak? What does G4 do?
422 
423  // Calculate Cherenkov angle integrals:
424  // This is an ad-hoc solution (we hold it in the class, not in the material)
425  chAngleIntegrals_.reset(new G4PhysicsOrderedFreeVector());
426 
427  int index = 0;
428  double currentRI = RefractiveIndex[index];
429  double currentPM = PhotonEnergy[index];
430  double currentCAI = 0.0;
431  chAngleIntegrals_.get()->InsertValues(currentPM, currentCAI);
432  double prevPM = currentPM;
433  double prevCAI = currentCAI;
434  double prevRI = currentRI;
435  while (++index < nEntries) {
436  currentRI = RefractiveIndex[index];
437  currentPM = PhotonEnergy[index];
438  currentCAI = 0.5 * (1.0 / (prevRI * prevRI) + 1.0 / (currentRI * currentRI));
439  currentCAI = prevCAI + (currentPM - prevPM) * currentCAI;
440 
441  chAngleIntegrals_.get()->InsertValues(currentPM, currentCAI);
442 
443  prevPM = currentPM;
444  prevCAI = currentCAI;
445  prevRI = currentRI;
446  }
447 
448  LogDebug("EcalSim") << "Material properties set for " << aMaterial->GetName();
449 
450  return true;
451 }
#define LogDebug(id)
std::unique_ptr< G4PhysicsOrderedFreeVector > chAngleIntegrals_
Table of Cherenkov angle integrals vs photon momentum.
Definition: DreamSD.h:59

Member Data Documentation

double DreamSD::birk1
private

Definition at line 52 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

double DreamSD::birk2
private

Definition at line 52 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

double DreamSD::birk3
private

Definition at line 52 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

std::unique_ptr<G4PhysicsOrderedFreeVector> DreamSD::chAngleIntegrals_
private

Table of Cherenkov angle integrals vs photon momentum.

Definition at line 59 of file DreamSD.h.

Referenced by DreamSD(), getAverageNumberOfPhotons_(), and setPbWO2MaterialProperties_().

bool DreamSD::doCherenkov_
private

Definition at line 51 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

G4MaterialPropertiesTable* DreamSD::materialPropertiesTable
private

Definition at line 60 of file DreamSD.h.

Referenced by cherenkovDeposit_(), and initRun().

int DreamSD::nphotons_
private

Definition at line 62 of file DreamSD.h.

bool DreamSD::readBothSide_
private

Definition at line 51 of file DreamSD.h.

Referenced by DreamSD(), and setDetUnitId().

int DreamSD::side
private

Definition at line 56 of file DreamSD.h.

Referenced by getEnergyDeposit(), and setDetUnitId().

double DreamSD::slopeLY
private

Definition at line 53 of file DreamSD.h.

Referenced by curve_LY(), and DreamSD().

bool DreamSD::useBirk
private

Definition at line 51 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

DimensionMap DreamSD::xtalLMap
private

Definition at line 54 of file DreamSD.h.

Referenced by crystalLength(), crystalWidth(), initMap(), and initRun().