CMS 3D CMS Logo

List of all members | Public Member Functions | Protected Member Functions | Private Types | Private Member Functions | Private Attributes | Static 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
 
bool ProcessHits (G4GFlashSpot *aSpot, G4TouchableHistory *) override
 
G4bool ProcessHits (G4Step *step, 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 ::EndOfEvent *) override
 
void update (const BeginOfEvent *) override
 This routine will be called when the appropriate signal arrives. More...
 
void update (const BeginOfRun *) 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 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)
 
void fillMap (const std::string &, double, double)
 
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 dd4hep_
 
bool doCherenkov_
 
G4MaterialPropertiesTable * materialPropertiesTable_
 
int nphotons_
 
bool readBothSide_
 
int side_
 
double slopeLY_
 
bool useBirk_
 
DimensionMap xtalLMap_
 

Static Private Attributes

static constexpr double k_ScaleFromDD4HepToG4 = 10.0
 
static constexpr double k_ScaleFromDDDToG4 = 1.0
 

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

◆ DimensionMap

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

Definition at line 32 of file DreamSD.h.

◆ Doubles

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

Definition at line 31 of file DreamSD.h.

Constructor & Destructor Documentation

◆ DreamSD()

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

Definition at line 35 of file DreamSD.cc.

40  : CaloSD(name, es, clg, p, manager) {
41  edm::ParameterSet m_EC = p.getParameter<edm::ParameterSet>("ECalSD");
42  useBirk_ = m_EC.getParameter<bool>("UseBirkLaw");
43  doCherenkov_ = m_EC.getParameter<bool>("doCherenkov");
44  birk1_ = m_EC.getParameter<double>("BirkC1") * (g / (MeV * cm2));
45  birk2_ = m_EC.getParameter<double>("BirkC2");
46  birk3_ = m_EC.getParameter<double>("BirkC3");
47  slopeLY_ = m_EC.getParameter<double>("SlopeLightYield");
48  readBothSide_ = m_EC.getUntrackedParameter<bool>("ReadBothSide", false);
49  dd4hep_ = m_EC.getUntrackedParameter<bool>("DD4Hep", false);
50 
51  chAngleIntegrals_.reset(nullptr);
52 
53  edm::LogVerbatim("EcalSim") << "Constructing a DreamSD with name " << GetName()
54  << "\nDreamSD:: Use of Birks law is set to " << useBirk_
55  << " with three constants kB = " << birk1_ << ", C1 = " << birk2_ << ", C2 = " << birk3_
56  << "\n Slope for Light yield is set to " << slopeLY_
57  << "\n Parameterization of Cherenkov is set to " << doCherenkov_
58  << ", readout both sides is " << readBothSide_ << " and dd4hep flag " << dd4hep_;
59 
60  initMap(name, es);
61 }

References birk1_, birk2_, birk3_, chAngleIntegrals_, dd4hep_, doCherenkov_, g, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), initMap(), MeV, Skims_PA_cff::name, AlCaHLTBitMon_ParallelJobs::p, readBothSide_, slopeLY_, and useBirk_.

◆ ~DreamSD()

DreamSD::~DreamSD ( )
inlineoverride

Definition at line 22 of file DreamSD.h.

22 {}

Member Function Documentation

◆ cherenkovDeposit_()

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

Returns the total energy due to Cherenkov radiation.

Definition at line 239 of file DreamSD.cc.

239  {
240  double cherenkovEnergy = 0;
242  return cherenkovEnergy;
243  G4Material *material = aStep->GetTrack()->GetMaterial();
244 
245  // Retrieve refractive index
246  G4MaterialPropertyVector *Rindex = materialPropertiesTable_->GetProperty("RINDEX");
247  if (Rindex == nullptr) {
248  edm::LogWarning("EcalSim") << "Couldn't retrieve refractive index";
249  return cherenkovEnergy;
250  }
251 
252  // V.Ivanchenko - temporary close log output for 9.5
253  // Material refraction properties
254  int Rlength = Rindex->GetVectorLength() - 1;
255  double Pmin = Rindex->Energy(0);
256  double Pmax = Rindex->Energy(Rlength);
257 #ifdef EDM_ML_DEBUG
258  edm::LogVerbatim("EcalSim") << "Material properties: \n Pmin = " << Pmin << " Pmax = " << Pmax;
259 #endif
260  // Get particle properties
261  const G4StepPoint *pPreStepPoint = aStep->GetPreStepPoint();
262  const G4StepPoint *pPostStepPoint = aStep->GetPostStepPoint();
263  const G4ThreeVector &x0 = pPreStepPoint->GetPosition();
264  G4ThreeVector p0 = aStep->GetDeltaPosition().unit();
265  const G4DynamicParticle *aParticle = aStep->GetTrack()->GetDynamicParticle();
266  const double charge = aParticle->GetDefinition()->GetPDGCharge();
267  // beta is averaged over step
268  double beta = 0.5 * (pPreStepPoint->GetBeta() + pPostStepPoint->GetBeta());
269  double BetaInverse = 1.0 / beta;
270 
271 #ifdef EDM_ML_DEBUG
272  edm::LogVerbatim("EcalSim") << "Particle properties: \n charge = " << charge << " beta = " << beta;
273 #endif
274 
275  // Now get number of photons generated in this step
276  double meanNumberOfPhotons = getAverageNumberOfPhotons_(charge, beta, material, Rindex);
277  if (meanNumberOfPhotons <= 0.0) { // Don't do anything
278 #ifdef EDM_ML_DEBUG
279  edm::LogVerbatim("EcalSim") << "Mean number of photons is zero: " << meanNumberOfPhotons << ", stopping here";
280 #endif
281  return cherenkovEnergy;
282  }
283 
284  // number of photons is in unit of Geant4...
285  meanNumberOfPhotons *= aStep->GetStepLength();
286 
287  // Now get a poisson distribution
288  int numPhotons = static_cast<int>(G4Poisson(meanNumberOfPhotons));
289  // edm::LogVerbatim("EcalSim") << "Number of photons = " << numPhotons;
290  if (numPhotons <= 0) {
291 #ifdef EDM_ML_DEBUG
292  edm::LogVerbatim("EcalSim") << "Poission number of photons is zero: " << numPhotons << ", stopping here";
293 #endif
294  return cherenkovEnergy;
295  }
296 
297  // Material refraction properties
298  double dp = Pmax - Pmin;
299  double maxCos = BetaInverse / (*Rindex)[Rlength];
300  double maxSin2 = (1.0 - maxCos) * (1.0 + maxCos);
301 
302  // Finally: get contribution of each photon
303  for (int iPhoton = 0; iPhoton < numPhotons; ++iPhoton) {
304  // Determine photon momentum
305  double randomNumber;
306  double sampledMomentum, sampledRI;
307  double cosTheta, sin2Theta;
308 
309  // sample a momentum (not sure why this is needed!)
310  do {
311  randomNumber = G4UniformRand();
312  sampledMomentum = Pmin + randomNumber * dp;
313  sampledRI = Rindex->Value(sampledMomentum);
314  cosTheta = BetaInverse / sampledRI;
315 
316  sin2Theta = (1.0 - cosTheta) * (1.0 + cosTheta);
317  randomNumber = G4UniformRand();
318 
319  } while (randomNumber * maxSin2 > sin2Theta);
320 
321  // Generate random position of photon on cone surface
322  // defined by Theta
323  randomNumber = G4UniformRand();
324 
325  double phi = twopi * randomNumber;
326  double sinPhi = sin(phi);
327  double cosPhi = cos(phi);
328 
329  // Create photon momentum direction vector
330  // The momentum direction is still w.r.t. the coordinate system where the
331  // primary particle direction is aligned with the z axis
332  double sinTheta = sqrt(sin2Theta);
333  double px = sinTheta * cosPhi;
334  double py = sinTheta * sinPhi;
335  double pz = cosTheta;
336  G4ThreeVector photonDirection(px, py, pz);
337 
338  // Rotate momentum direction back to global (crystal) reference system
339  photonDirection.rotateUz(p0);
340 
341  // Create photon position and momentum
342  randomNumber = G4UniformRand();
343  G4ThreeVector photonPosition = x0 + randomNumber * aStep->GetDeltaPosition();
344  G4ThreeVector photonMomentum = sampledMomentum * photonDirection;
345 
346  // Collect energy on APD
347  cherenkovEnergy += getPhotonEnergyDeposit_(photonMomentum, photonPosition, aStep);
348  }
349  return cherenkovEnergy;
350 }

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

Referenced by getEnergyDeposit().

◆ crystalLength()

double DreamSD::crystalLength ( G4LogicalVolume *  lv) const
private

Definition at line 219 of file DreamSD.cc.

219  {
220  double length = -1.;
221  DimensionMap::const_iterator ite = xtalLMap_.find(lv);
222  if (ite != xtalLMap_.end())
223  length = ite->second.first;
224  return length;
225 }

References xtalLMap_.

Referenced by curve_LY(), and getPhotonEnergyDeposit_().

◆ crystalWidth()

double DreamSD::crystalWidth ( G4LogicalVolume *  lv) const
private

Definition at line 228 of file DreamSD.cc.

228  {
229  double width = -1.;
230  DimensionMap::const_iterator ite = xtalLMap_.find(lv);
231  if (ite != xtalLMap_.end())
232  width = ite->second.second;
233  return width;
234 }

References ApeEstimator_cff::width, and xtalLMap_.

Referenced by getPhotonEnergyDeposit_().

◆ curve_LY()

double DreamSD::curve_LY ( const G4Step *  aStep,
int  flag 
)
private

Definition at line 192 of file DreamSD.cc.

192  {
193  auto const stepPoint = aStep->GetPreStepPoint();
194  auto const lv = stepPoint->GetTouchable()->GetVolume(0)->GetLogicalVolume();
195  G4String nameVolume = lv->GetName();
196 
197  double weight = 1.;
198  G4ThreeVector localPoint = setToLocal(stepPoint->GetPosition(), stepPoint->GetTouchable());
199  double crlength = crystalLength(lv);
200  double localz = localPoint.x();
201  double dapd = 0.5 * crlength - flag * localz; // Distance from closest APD
202  if (dapd >= -0.1 || dapd <= crlength + 0.1) {
203  if (dapd <= 100.)
204  weight = 1.0 + slopeLY_ - dapd * 0.01 * slopeLY_;
205  } else {
206  edm::LogWarning("EcalSim") << "DreamSD: light coll curve : wrong distance "
207  << "to APD " << dapd << " crlength = " << crlength << " crystal name = " << nameVolume
208  << " z of localPoint = " << localz << " take weight = " << weight;
209  }
210 #ifdef EDM_ML_DEBUG
211  edm::LogVerbatim("EcalSim") << "DreamSD, light coll curve : " << dapd << " crlength = " << crlength
212  << " crystal name = " << nameVolume << " z of localPoint = " << localz
213  << " take weight = " << weight;
214 #endif
215  return weight;
216 }

References crystalLength(), RemoveAddSevLevel::flag, CaloSD::setToLocal(), slopeLY_, and mps_merge::weight.

Referenced by getEnergyDeposit().

◆ fillMap()

void DreamSD::fillMap ( const std::string &  name,
double  length,
double  width 
)
private

Definition at line 174 of file DreamSD.cc.

174  {
175  const G4LogicalVolumeStore *lvs = G4LogicalVolumeStore::GetInstance();
176  std::vector<G4LogicalVolume *>::const_iterator lvcite;
177  G4LogicalVolume *lv = nullptr;
178  for (lvcite = lvs->begin(); lvcite != lvs->end(); lvcite++) {
179  if ((*lvcite)->GetName() == static_cast<G4String>(name)) {
180  lv = (*lvcite);
181  break;
182  }
183  }
184 #ifdef EDM_ML_DEBUG
185  edm::LogVerbatim("EcalSim") << "DreamSD::fillMap (for " << name << " Logical Volume " << lv << " Length " << length
186  << " Width " << width;
187 #endif
188  xtalLMap_.insert(std::pair<G4LogicalVolume *, Doubles>(lv, Doubles(length, width)));
189 }

References Skims_PA_cff::name, ApeEstimator_cff::width, and xtalLMap_.

Referenced by initMap().

◆ getAverageNumberOfPhotons_()

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 355 of file DreamSD.cc.

358  {
359  const double rFact = 369.81 / (eV * cm);
360 
361  if (beta <= 0.0)
362  return 0.0;
363 
364  double BetaInverse = 1. / beta;
365 
366  // Vectors used in computation of Cerenkov Angle Integral:
367  // - Refraction Indices for the current material
368  // - new G4PhysicsOrderedFreeVector allocated to hold CAI's
369 
370  // Min and Max photon momenta
371  int Rlength = Rindex->GetVectorLength() - 1;
372  double Pmin = Rindex->Energy(0);
373  double Pmax = Rindex->Energy(Rlength);
374 
375  // Min and Max Refraction Indices
376  double nMin = (*Rindex)[0];
377  double nMax = (*Rindex)[Rlength];
378 
379  // Max Cerenkov Angle Integral
380  double CAImax = chAngleIntegrals_.get()->GetMaxValue();
381 
382  double dp = 0., ge = 0., CAImin = 0.;
383 
384  // If n(Pmax) < 1/Beta -- no photons generated
385  if (nMax < BetaInverse) {
386  }
387 
388  // otherwise if n(Pmin) >= 1/Beta -- photons generated
389  else if (nMin > BetaInverse) {
390  dp = Pmax - Pmin;
391  ge = CAImax;
392  }
393  // If n(Pmin) < 1/Beta, and n(Pmax) >= 1/Beta, then
394  // we need to find a P such that the value of n(P) == 1/Beta.
395  // Interpolation is performed by the GetPhotonEnergy() and
396  // GetProperty() methods of the G4MaterialPropertiesTable and
397  // the GetValue() method of G4PhysicsVector.
398  else {
399  Pmin = Rindex->Value(BetaInverse);
400  dp = Pmax - Pmin;
401  // need boolean for current implementation of G4PhysicsVector
402  // ==> being phased out
403  double CAImin = chAngleIntegrals_->Value(Pmin);
404  ge = CAImax - CAImin;
405  }
406 
407  // Calculate number of photons
408  double numPhotons = rFact * charge / eplus * charge / eplus * (dp - ge * BetaInverse * BetaInverse);
409 
410 #ifdef EDM_ML_DEBUG
411  edm::LogVerbatim("EcalSim") << "@SUB=getAverageNumberOfPhotons\nCAImin = " << CAImin << "\nCAImax = " << CAImax
412  << "\ndp = " << dp << ", ge = " << ge << "\nnumPhotons = " << numPhotons;
413 #endif
414  return numPhotons;
415 }

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

Referenced by cherenkovDeposit_().

◆ getEnergyDeposit()

double DreamSD::getEnergyDeposit ( const G4Step *  aStep)
overrideprotectedvirtual

Reimplemented from CaloSD.

Definition at line 64 of file DreamSD.cc.

64  {
65  // take into account light collection curve for crystals
66  double weight = curve_LY(aStep, side_);
67  if (useBirk_)
69  double edep = aStep->GetTotalEnergyDeposit() * weight;
70 
71  // Get Cerenkov contribution
72  if (doCherenkov_) {
73  edep += cherenkovDeposit_(aStep);
74  }
75 #ifdef EDM_ML_DEBUG
76  edm::LogVerbatim("EcalSim") << "DreamSD:: " << aStep->GetPreStepPoint()->GetPhysicalVolume()->GetName() << " Side "
77  << side_ << " Light Collection Efficiency " << weight << " Weighted Energy Deposit "
78  << edep / CLHEP::MeV << " MeV";
79 #endif
80  return edep;
81 }

References birk1_, birk2_, birk3_, cherenkovDeposit_(), curve_LY(), doCherenkov_, CaloSD::getAttenuation(), MeV, side_, useBirk_, and mps_merge::weight.

◆ getPhotonEnergyDeposit_()

double DreamSD::getPhotonEnergyDeposit_ ( const G4ParticleMomentum &  p,
const G4ThreeVector &  x,
const G4Step *  aStep 
)
private

Returns energy deposit for a given photon.

Definition at line 501 of file DreamSD.cc.

501  {
502  double energy = 0;
503 
504  // Crystal dimensions
505 
506  // edm::LogVerbatim("EcalSim") << p << x;
507 
508  // 1. Check if this photon goes straight to the APD:
509  // - assume that APD is at x=xtalLength/2.0
510  // - extrapolate from x=x0 to x=xtalLength/2.0 using momentum in x-y
511 
512  G4StepPoint *stepPoint = aStep->GetPreStepPoint();
513  G4LogicalVolume *lv = stepPoint->GetTouchable()->GetVolume(0)->GetLogicalVolume();
514  G4String nameVolume = lv->GetName();
515 
516  double crlength = crystalLength(lv);
517  double crwidth = crystalWidth(lv);
518  double dapd = 0.5 * crlength - x.x(); // Distance from closest APD
519  double y = p.y() / p.x() * dapd;
520 
521 #ifdef EDM_ML_DEBUG
522  edm::LogVerbatim("EcalSim") << "Distance to APD: " << dapd << " - y at APD: " << y;
523 #endif
524  // Not straight: compute probability
525  if (std::abs(y) > crwidth * 0.5) {
526  }
527 
528  // 2. Retrieve efficiency for this wavelength (in nm, from MeV)
529  double waveLength = p.mag() * 1.239e8;
530 
531  energy = p.mag() * PMTResponse::getEfficiency(waveLength);
532 #ifdef EDM_ML_DEBUG
533  edm::LogVerbatim("EcalSim") << "Wavelength: " << waveLength << " - Energy: " << energy;
534 #endif
535  return energy;
536 }

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

Referenced by cherenkovDeposit_().

◆ initMap()

void DreamSD::initMap ( const std::string &  sd,
const edm::EventSetup es 
)
private

Definition at line 115 of file DreamSD.cc.

115  {
116  if (dd4hep_) {
118  es.get<IdealGeometryRecord>().get(cpv);
119  const cms::DDFilter filter("ReadOutName", sd);
120  cms::DDFilteredView fv((*cpv), filter);
121  while (fv.firstChild()) {
122  std::string name = static_cast<std::string>(cms::dd::noNamespace(fv.name()));
123  std::vector<double> paras(fv.parameters());
124 #ifdef EDM_ML_DEBUG
125  edm::LogVerbatim("EcalSim") << "DreamSD::initMap (for " << sd << "): Solid " << name << " Shape "
126  << cms::dd::name(cms::DDSolidShapeMap, fv.shape()) << " Parameter 0 = " << paras[0];
127 #endif
128  // Set length to be the largest size, width the smallest
129  std::sort(paras.begin(), paras.end());
130  double length = 2.0 * k_ScaleFromDD4HepToG4 * paras.back();
131  double width = 2.0 * k_ScaleFromDD4HepToG4 * paras.front();
132  fillMap(name, length, width);
133  }
134  } else {
136  es.get<IdealGeometryRecord>().get(cpv);
137  DDSpecificsMatchesValueFilter filter{DDValue("ReadOutName", sd, 0)};
138  DDFilteredView fv((*cpv), filter);
139 
140  bool dodet = fv.firstChild();
141  while (dodet) {
142  const DDSolid &sol = fv.logicalPart().solid();
143  std::vector<double> paras(sol.parameters());
144  std::string name = sol.name().name();
145 #ifdef EDM_ML_DEBUG
146  edm::LogVerbatim("EcalSim") << "DreamSD::initMap (for " << sd << "): Solid " << name << " Shape " << sol.shape()
147  << " Parameter 0 = " << paras[0];
148 #endif
149  // Set length to be the largest size, width the smallest
150  std::sort(paras.begin(), paras.end());
151  double length = 2.0 * k_ScaleFromDDDToG4 * paras.back();
152  double width = 2.0 * k_ScaleFromDDDToG4 * paras.front();
153  fillMap(name, length, width);
154  dodet = fv.next();
155  }
156  }
157 #ifdef EDM_ML_DEBUG
158  edm::LogVerbatim("EcalSim") << "DreamSD: Length Table for ReadOutName = " << sd << ":";
159 #endif
160  DimensionMap::const_iterator ite = xtalLMap_.begin();
161  int i = 0;
162  for (; ite != xtalLMap_.end(); ite++, i++) {
163  G4String name = "Unknown";
164  if (ite->first != nullptr)
165  name = (ite->first)->GetName();
166 #ifdef EDM_ML_DEBUG
167  edm::LogVerbatim("EcalSim") << " " << i << " " << ite->first << " " << name << " L = " << ite->second.first
168  << " W = " << ite->second.second;
169 #endif
170  }
171 }

References dd4hep_, cms::DDSolidShapeMap, fillMap(), ALCARECOTkAlBeamHalo_cff::filter, DDFilteredView::firstChild(), cms::DDFilteredView::firstChild(), edm::EventSetup::get(), get, mps_fire::i, k_ScaleFromDD4HepToG4, k_ScaleFromDDDToG4, DDFilteredView::logicalPart(), Skims_PA_cff::name, cms::dd::name(), DDName::name(), DDBase< N, C >::name(), cms::DDFilteredView::name(), DDFilteredView::next(), cms::dd::noNamespace(), DDSolid::parameters(), cms::DDFilteredView::parameters(), sd, DDSolid::shape(), cms::DDFilteredView::shape(), DDLogicalPart::solid(), AlCaHLTBitMon_QueryRunRegistry::string, ApeEstimator_cff::width, and xtalLMap_.

Referenced by DreamSD().

◆ initRun()

void DreamSD::initRun ( )
overrideprotectedvirtual

Reimplemented from CaloSD.

Definition at line 84 of file DreamSD.cc.

84  {
85  // Get the material and set properties if needed
86  DimensionMap::const_iterator ite = xtalLMap_.begin();
87  const G4LogicalVolume *lv = (ite->first);
88  G4Material *material = lv->GetMaterial();
89  edm::LogVerbatim("EcalSim") << "DreamSD::initRun: Initializes for material " << material->GetName() << " in "
90  << lv->GetName();
91  materialPropertiesTable_ = material->GetMaterialPropertiesTable();
93  if (!setPbWO2MaterialProperties_(material)) {
94  edm::LogWarning("EcalSim") << "Couldn't retrieve material properties table\n Material = " << material->GetName();
95  }
96  materialPropertiesTable_ = material->GetMaterialPropertiesTable();
97  }
98 }

References materialPropertiesTable_, setPbWO2MaterialProperties_(), and xtalLMap_.

◆ setDetUnitId()

uint32_t DreamSD::setDetUnitId ( const G4Step *  aStep)
overridevirtual

Implements CaloSD.

Definition at line 101 of file DreamSD.cc.

101  {
102  const G4VTouchable *touch = aStep->GetPreStepPoint()->GetTouchable();
103  uint32_t id = (touch->GetReplicaNumber(1)) * 10 + (touch->GetReplicaNumber(0));
104  side_ = readBothSide_ ? -1 : 1;
105  if (side_ < 0) {
106  ++id;
107  }
108 #ifdef EDM_ML_DEBUG
109  edm::LogVerbatim("EcalSim") << "DreamSD:: ID " << id;
110 #endif
111  return id;
112 }

References triggerObjects_cff::id, readBothSide_, and side_.

◆ setPbWO2MaterialProperties_()

bool DreamSD::setPbWO2MaterialProperties_ ( G4Material *  aMaterial)
private

Sets material properties at run-time...

Definition at line 420 of file DreamSD.cc.

420  {
421  std::string pbWO2Name("E_PbWO4");
422  if (pbWO2Name != aMaterial->GetName()) { // Wrong material!
423  edm::LogWarning("EcalSim") << "This is not the right material: "
424  << "expecting " << pbWO2Name << ", got " << aMaterial->GetName();
425  return false;
426  }
427 
428  G4MaterialPropertiesTable *table = new G4MaterialPropertiesTable();
429 
430  // Refractive index as a function of photon momentum
431  // FIXME: Should somehow put that in the configuration
432  const int nEntries = 14;
433  double PhotonEnergy[nEntries] = {1.7712 * eV,
434  1.8368 * eV,
435  1.90745 * eV,
436  1.98375 * eV,
437  2.0664 * eV,
438  2.15625 * eV,
439  2.25426 * eV,
440  2.3616 * eV,
441  2.47968 * eV,
442  2.61019 * eV,
443  2.75521 * eV,
444  2.91728 * eV,
445  3.09961 * eV,
446  3.30625 * eV};
447  double RefractiveIndex[nEntries] = {2.17728,
448  2.18025,
449  2.18357,
450  2.18753,
451  2.19285,
452  2.19813,
453  2.20441,
454  2.21337,
455  2.22328,
456  2.23619,
457  2.25203,
458  2.27381,
459  2.30282,
460  2.34666};
461 
462  table->AddProperty("RINDEX", PhotonEnergy, RefractiveIndex, nEntries);
463  aMaterial->SetMaterialPropertiesTable(table); // FIXME: could this leak? What does G4 do?
464 
465  // Calculate Cherenkov angle integrals:
466  // This is an ad-hoc solution (we hold it in the class, not in the material)
467  chAngleIntegrals_.reset(new G4PhysicsOrderedFreeVector());
468 
469  int index = 0;
470  double currentRI = RefractiveIndex[index];
471  double currentPM = PhotonEnergy[index];
472  double currentCAI = 0.0;
473  chAngleIntegrals_.get()->InsertValues(currentPM, currentCAI);
474  double prevPM = currentPM;
475  double prevCAI = currentCAI;
476  double prevRI = currentRI;
477  while (++index < nEntries) {
478  currentRI = RefractiveIndex[index];
479  currentPM = PhotonEnergy[index];
480  currentCAI = 0.5 * (1.0 / (prevRI * prevRI) + 1.0 / (currentRI * currentRI));
481  currentCAI = prevCAI + (currentPM - prevPM) * currentCAI;
482 
483  chAngleIntegrals_.get()->InsertValues(currentPM, currentCAI);
484 
485  prevPM = currentPM;
486  prevCAI = currentCAI;
487  prevRI = currentRI;
488  }
489 
490 #ifdef EDM_ML_DEBUG
491  edm::LogVerbatim("EcalSim") << "Material properties set for " << aMaterial->GetName();
492 #endif
493  return true;
494 }

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

Referenced by initRun().

Member Data Documentation

◆ birk1_

double DreamSD::birk1_
private

Definition at line 56 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

◆ birk2_

double DreamSD::birk2_
private

Definition at line 56 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

◆ birk3_

double DreamSD::birk3_
private

Definition at line 56 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

◆ chAngleIntegrals_

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

Table of Cherenkov angle integrals vs photon momentum.

Definition at line 63 of file DreamSD.h.

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

◆ dd4hep_

bool DreamSD::dd4hep_
private

Definition at line 55 of file DreamSD.h.

Referenced by DreamSD(), and initMap().

◆ doCherenkov_

bool DreamSD::doCherenkov_
private

Definition at line 55 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

◆ k_ScaleFromDD4HepToG4

constexpr double DreamSD::k_ScaleFromDD4HepToG4 = 10.0
staticconstexprprivate

Definition at line 53 of file DreamSD.h.

Referenced by initMap().

◆ k_ScaleFromDDDToG4

constexpr double DreamSD::k_ScaleFromDDDToG4 = 1.0
staticconstexprprivate

Definition at line 52 of file DreamSD.h.

Referenced by initMap().

◆ materialPropertiesTable_

G4MaterialPropertiesTable* DreamSD::materialPropertiesTable_
private

Definition at line 64 of file DreamSD.h.

Referenced by cherenkovDeposit_(), and initRun().

◆ nphotons_

int DreamSD::nphotons_
private

Definition at line 66 of file DreamSD.h.

◆ readBothSide_

bool DreamSD::readBothSide_
private

Definition at line 55 of file DreamSD.h.

Referenced by DreamSD(), and setDetUnitId().

◆ side_

int DreamSD::side_
private

Definition at line 60 of file DreamSD.h.

Referenced by getEnergyDeposit(), and setDetUnitId().

◆ slopeLY_

double DreamSD::slopeLY_
private

Definition at line 57 of file DreamSD.h.

Referenced by curve_LY(), and DreamSD().

◆ useBirk_

bool DreamSD::useBirk_
private

Definition at line 55 of file DreamSD.h.

Referenced by DreamSD(), and getEnergyDeposit().

◆ xtalLMap_

DimensionMap DreamSD::xtalLMap_
private

Definition at line 58 of file DreamSD.h.

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

DreamSD::chAngleIntegrals_
std::unique_ptr< G4PhysicsOrderedFreeVector > chAngleIntegrals_
Table of Cherenkov angle integrals vs photon momentum.
Definition: DreamSD.h:63
DreamSD::readBothSide_
bool readBothSide_
Definition: DreamSD.h:55
DDAxes::y
ApeEstimator_cff.width
width
Definition: ApeEstimator_cff.py:24
DreamSD::Doubles
std::pair< double, double > Doubles
Definition: DreamSD.h:31
mps_fire.i
i
Definition: mps_fire.py:355
metDiagnosticParameterSet_cfi.nMin
nMin
Definition: metDiagnosticParameterSet_cfi.py:10
DreamSD::fillMap
void fillMap(const std::string &, double, double)
Definition: DreamSD.cc:174
multPhiCorr_741_25nsDY_cfi.py
py
Definition: multPhiCorr_741_25nsDY_cfi.py:12
zMuMuMuonUserData.beta
beta
Definition: zMuMuMuonUserData.py:10
mps_merge.weight
weight
Definition: mps_merge.py:88
DreamSD::k_ScaleFromDDDToG4
static constexpr double k_ScaleFromDDDToG4
Definition: DreamSD.h:52
AlCaHLTBitMon_ParallelJobs.p
p
Definition: AlCaHLTBitMon_ParallelJobs.py:153
DreamSD::cherenkovDeposit_
double cherenkovDeposit_(const G4Step *aStep)
Returns the total energy due to Cherenkov radiation.
Definition: DreamSD.cc:239
DreamSD::dd4hep_
bool dd4hep_
Definition: DreamSD.h:55
CaloSD::getAttenuation
double getAttenuation(const G4Step *aStep, double birk1, double birk2, double birk3) const
Definition: CaloSD.cc:435
DreamSD::getPhotonEnergyDeposit_
double getPhotonEnergyDeposit_(const G4ParticleMomentum &p, const G4ThreeVector &x, const G4Step *aStep)
Returns energy deposit for a given photon.
Definition: DreamSD.cc:501
MeV
const double MeV
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
DDAxes::x
DreamSD::crystalWidth
double crystalWidth(G4LogicalVolume *) const
Definition: DreamSD.cc:228
DreamSD::doCherenkov_
bool doCherenkov_
Definition: DreamSD.h:55
cms::DDFilteredView
Definition: DDFilteredView.h:65
DreamSD::materialPropertiesTable_
G4MaterialPropertiesTable * materialPropertiesTable_
Definition: DreamSD.h:64
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
DDSolid::shape
DDSolidShape shape(void) const
The type of the solid.
Definition: DDSolid.cc:119
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
edm::EventSetup::get
T get() const
Definition: EventSetup.h:73
cms::DDFilter
Definition: DDFilteredView.h:54
DreamSD::initMap
void initMap(const std::string &, const edm::EventSetup &)
Definition: DreamSD.cc:115
cms::dd::name
std::string name(Mapping a, V value)
Definition: DDSolidShapes.h:31
CaloSD::setToLocal
G4ThreeVector setToLocal(const G4ThreeVector &, const G4VTouchable *) const
Definition: CaloSD.cc:291
Calorimetry_cff.dp
dp
Definition: Calorimetry_cff.py:157
DreamSD::useBirk_
bool useBirk_
Definition: DreamSD.h:55
DreamSD::getAverageNumberOfPhotons_
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:355
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
DDBase::name
const N & name() const
Definition: DDBase.h:59
cms::dd::noNamespace
std::string_view noNamespace(std::string_view)
Definition: Filter.cc:63
DreamSD::xtalLMap_
DimensionMap xtalLMap_
Definition: DreamSD.h:58
HCALHighEnergyHPDFilter_cfi.energy
energy
Definition: HCALHighEnergyHPDFilter_cfi.py:5
ALCARECOTkAlBeamHalo_cff.filter
filter
Definition: ALCARECOTkAlBeamHalo_cff.py:27
AlCaHLTBitMon_QueryRunRegistry.string
string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:256
DreamSD::setPbWO2MaterialProperties_
bool setPbWO2MaterialProperties_(G4Material *aMaterial)
Sets material properties at run-time...
Definition: DreamSD.cc:420
edm::LogWarning
Definition: MessageLogger.h:141
ALCARECOTkAlJpsiMuMu_cff.charge
charge
Definition: ALCARECOTkAlJpsiMuMu_cff.py:47
DreamSD::birk1_
double birk1_
Definition: DreamSD.h:56
edm::ParameterSet
Definition: ParameterSet.h:36
cms::DDSolidShapeMap
const std::array< const cms::dd::NameValuePair< DDSolidShape >, 19 > DDSolidShapeMap
Definition: DDSolidShapes.h:84
PMTResponse::getEfficiency
static double getEfficiency(const double &waveLengthNm)
Return efficiency for given photon wavelength (in nm)
Definition: PMTResponse.cc:6
CaloSD::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)
Definition: CaloSD.cc:27
DDName::name
const std::string & name() const
Returns the name.
Definition: DDName.cc:40
edm::LogVerbatim
Definition: MessageLogger.h:297
get
#define get
DreamSD::birk2_
double birk2_
Definition: DreamSD.h:56
edm::ESTransientHandle
Definition: ESTransientHandle.h:41
DreamSD::k_ScaleFromDD4HepToG4
static constexpr double k_ScaleFromDD4HepToG4
Definition: DreamSD.h:53
DDAxes::phi
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
multPhiCorr_741_25nsDY_cfi.px
px
Definition: multPhiCorr_741_25nsDY_cfi.py:10
DreamSD::crystalLength
double crystalLength(G4LogicalVolume *) const
Definition: DreamSD.cc:219
dataAnalyzerFineBiningParameters_cff.Pmax
Pmax
Definition: dataAnalyzerFineBiningParameters_cff.py:10
DDValue
Definition: DDValue.h:21
triggerObjects_cff.id
id
Definition: triggerObjects_cff.py:31
DDSolid::parameters
const std::vector< double > & parameters(void) const
Give the parameters of the solid.
Definition: DDSolid.cc:121
Skims_PA_cff.name
name
Definition: Skims_PA_cff.py:17
DDSolid
A DDSolid represents the shape of a part.
Definition: DDSolid.h:39
DDSpecificsMatchesValueFilter
Definition: DDFilter.h:70
AlignmentPI::index
index
Definition: AlignmentPayloadInspectorHelper.h:46
DreamSD::curve_LY
double curve_LY(const G4Step *, int)
Definition: DreamSD.cc:192
DreamSD::side_
int side_
Definition: DreamSD.h:60
sd
double sd
Definition: CascadeWrapper.h:113
DDFilteredView
Definition: DDFilteredView.h:20
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
DreamSD::slopeLY_
double slopeLY_
Definition: DreamSD.h:57
TableParser.table
table
Definition: TableParser.py:111
metDiagnosticParameterSet_cfi.nMax
nMax
Definition: metDiagnosticParameterSet_cfi.py:11
weight
Definition: weight.py:1
IdealGeometryRecord
Definition: IdealGeometryRecord.h:27
g
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
RemoveAddSevLevel.flag
flag
Definition: RemoveAddSevLevel.py:116
DreamSD::birk3_
double birk3_
Definition: DreamSD.h:56