CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
HCalSD.cc
Go to the documentation of this file.
1 // File: HCalSD.cc
3 // Description: Sensitive Detector class for calorimeters
5 
16 
19 
20 #include "G4LogicalVolumeStore.hh"
21 #include "G4LogicalVolume.hh"
22 #include "G4Step.hh"
23 #include "G4Track.hh"
24 #include "G4ParticleTable.hh"
25 
26 #include <iostream>
27 #include <fstream>
28 #include <iomanip>
29 
30 //#define DebugLog
31 
32 HCalSD::HCalSD(G4String name, const DDCompactView & cpv,
34  edm::ParameterSet const & p, const SimTrackManager* manager) :
35  CaloSD(name, cpv, clg, p, manager,
36  p.getParameter<edm::ParameterSet>("HCalSD").getParameter<int>("TimeSliceUnit"),
37  p.getParameter<edm::ParameterSet>("HCalSD").getParameter<bool>("IgnoreTrackID")),
38  numberingFromDDD(0), numberingScheme(0), showerLibrary(0), hfshower(0),
39  showerParam(0), showerPMT(0), showerBundle(0){
40 
41  //static SimpleConfigurable<bool> on1(false, "HCalSD:UseBirkLaw");
42  //static SimpleConfigurable<double> bk1(0.013, "HCalSD:BirkC1");
43  //static SimpleConfigurable<double> bk2(0.0568,"HCalSD:BirkC2");
44  //static SimpleConfigurable<double> bk3(1.75, "HCalSD:BirkC3");
45  // Values from NIM 80 (1970) 239-244: as implemented in Geant3
46  //static SimpleConfigurable<bool> on2(true,"HCalSD:UseShowerLibrary");
47 
49  useBirk = m_HC.getParameter<bool>("UseBirkLaw");
50  birk1 = m_HC.getParameter<double>("BirkC1")*(g/(MeV*cm2));
51  birk2 = m_HC.getParameter<double>("BirkC2");
52  birk3 = m_HC.getParameter<double>("BirkC3");
53  useShowerLibrary = m_HC.getParameter<bool>("UseShowerLibrary");
54  useParam = m_HC.getParameter<bool>("UseParametrize");
55  bool testNumber = m_HC.getParameter<bool>("TestNumberingScheme");
56  usePMTHit = m_HC.getParameter<bool>("UsePMTHits");
57  betaThr = m_HC.getParameter<double>("BetaThreshold");
58  eminHitHB = m_HC.getParameter<double>("EminHitHB")*MeV;
59  eminHitHE = m_HC.getParameter<double>("EminHitHE")*MeV;
60  eminHitHO = m_HC.getParameter<double>("EminHitHO")*MeV;
61  eminHitHF = m_HC.getParameter<double>("EminHitHF")*MeV;
62  useFibreBundle = m_HC.getParameter<bool>("UseFibreBundleHits");
63  useHF = m_HC.getUntrackedParameter<bool>("UseHF",true);
64  bool forTBH2 = m_HC.getUntrackedParameter<bool>("ForTBH2",false);
65  useLayerWt = m_HC.getUntrackedParameter<bool>("UseLayerWt",false);
66  std::string file = m_HC.getUntrackedParameter<std::string>("WtFile","None");
67 
68 #ifdef DebugLog
69  LogDebug("HcalSim") << "***************************************************"
70  << "\n"
71  << "* *"
72  << "\n"
73  << "* Constructing a HCalSD with name " << name << "\n"
74  << "* *"
75  << "\n"
76  << "***************************************************";
77 #endif
78  edm::LogInfo("HcalSim") << "HCalSD:: Use of HF code is set to " << useHF
79  << "\nUse of shower parametrization set to "
80  << useParam << "\nUse of shower library is set to "
81  << useShowerLibrary << "\nUse PMT Hit is set to "
82  << usePMTHit << " with beta Threshold "<< betaThr
83  << "\nUSe of FibreBundle Hit set to "<<useFibreBundle
84  << "\n Use of Birks law is set to "
85  << useBirk << " with three constants kB = "
86  << birk1 << ", C1 = " << birk2 << ", C2 = " << birk3;
87  edm::LogInfo("HcalSim") << "HCalSD:: Suppression Flag " << suppressHeavy
88  << " protons below " << kmaxProton << " MeV,"
89  << " neutrons below " << kmaxNeutron << " MeV and"
90  << " ions below " << kmaxIon << " MeV\n"
91  << " Threshold for storing hits in HB: "
92  << eminHitHB << " HE: " << eminHitHE << " HO: "
93  << eminHitHO << " HF: " << eminHitHF;
94 
95  numberingFromDDD = new HcalNumberingFromDDD(name, cpv);
96  HcalNumberingScheme* scheme;
97  if (testNumber || forTBH2)
98  scheme = dynamic_cast<HcalNumberingScheme*>(new HcalTestNumberingScheme(forTBH2));
99  else
100  scheme = new HcalNumberingScheme();
101  setNumberingScheme(scheme);
102 
103  const G4LogicalVolumeStore * lvs = G4LogicalVolumeStore::GetInstance();
104  std::vector<G4LogicalVolume *>::const_iterator lvcite;
105  G4LogicalVolume* lv;
106  std::string attribute, value;
107  if (useHF) {
108  if (useParam) showerParam = new HFShowerParam(name, cpv, p);
109  else {
110  if (useShowerLibrary) showerLibrary = new HFShowerLibrary(name, cpv, p);
111  hfshower = new HFShower(name, cpv, p, 0);
112  }
113 
114  // HF volume names
115  attribute = "Volume";
116  value = "HF";
117  DDSpecificsFilter filter0;
118  DDValue ddv0(attribute,value,0);
120  DDFilteredView fv0(cpv);
121  fv0.addFilter(filter0);
122  hfNames = getNames(fv0);
123  fv0.firstChild();
124  DDsvalues_type sv0(fv0.mergedSpecifics());
125  std::vector<double> temp = getDDDArray("Levels",sv0);
126  edm::LogInfo("HcalSim") << "HCalSD: Names to be tested for " << attribute
127  << " = " << value << " has " << hfNames.size()
128  << " elements";
129  for (unsigned int i=0; i<hfNames.size(); i++) {
130  G4String namv = hfNames[i];
131  lv = 0;
132  for (lvcite = lvs->begin(); lvcite != lvs->end(); lvcite++)
133  if ((*lvcite)->GetName() == namv) {
134  lv = (*lvcite);
135  break;
136  }
137  hfLV.push_back(lv);
138  int level = static_cast<int>(temp[i]);
139  hfLevels.push_back(level);
140  edm::LogInfo("HcalSim") << "HCalSD: HF[" << i << "] = " << hfNames[i]
141  << " LV " << hfLV[i] << " at level "
142  << hfLevels[i];
143  }
144 
145  // HF Fibre volume names
146  value = "HFFibre";
147  DDSpecificsFilter filter1;
148  DDValue ddv1(attribute,value,0);
150  DDFilteredView fv1(cpv);
151  fv1.addFilter(filter1);
152  fibreNames = getNames(fv1);
153  edm::LogInfo("HcalSim") << "HCalSD: Names to be tested for " << attribute
154  << " = " << value << ":";
155  for (unsigned int i=0; i<fibreNames.size(); i++) {
156  G4String namv = fibreNames[i];
157  lv = 0;
158  for (lvcite = lvs->begin(); lvcite != lvs->end(); lvcite++)
159  if ((*lvcite)->GetName() == namv) {
160  lv = (*lvcite);
161  break;
162  }
163  fibreLV.push_back(lv);
164  edm::LogInfo("HcalSim") << "HCalSD: (" << i << ") " << fibreNames[i]
165  << " LV " << fibreLV[i];
166  }
167 
168  // HF PMT volume names
169  value = "HFPMT";
170  DDSpecificsFilter filter3;
171  DDValue ddv3(attribute,value,0);
173  DDFilteredView fv3(cpv);
174  fv3.addFilter(filter3);
175  std::vector<G4String> pmtNames = getNames(fv3);
176  edm::LogInfo("HcalSim") << "HCalSD: Names to be tested for " << attribute
177  << " = " << value << " have " << pmtNames.size()
178  << " entries";
179  for (unsigned int i=0; i<pmtNames.size(); i++) {
180  G4String namv = pmtNames[i];
181  lv = 0;
182  for (lvcite = lvs->begin(); lvcite != lvs->end(); lvcite++)
183  if ((*lvcite)->GetName() == namv) {
184  lv = (*lvcite);
185  break;
186  }
187  pmtLV.push_back(lv);
188  edm::LogInfo("HcalSim") << "HCalSD: (" << i << ") " << pmtNames[i]
189  << " LV " << pmtLV[i];
190  }
191  if (pmtNames.size() > 0) showerPMT = new HFShowerPMT (name, cpv, p);
192 
193  // HF Fibre bundles
194  value = "HFFibreBundleStraight";
195  DDSpecificsFilter filter4;
196  DDValue ddv4(attribute,value,0);
198  DDFilteredView fv4(cpv);
199  fv4.addFilter(filter4);
200  std::vector<G4String> fibreNames = getNames(fv4);
201  edm::LogInfo("HcalSim") << "HCalSD: Names to be tested for " << attribute
202  << " = " << value << " have " << fibreNames.size()
203  << " entries";
204  for (unsigned int i=0; i<fibreNames.size(); i++) {
205  G4String namv = fibreNames[i];
206  lv = 0;
207  for (lvcite = lvs->begin(); lvcite != lvs->end(); lvcite++)
208  if ((*lvcite)->GetName() == namv) {
209  lv = (*lvcite);
210  break;
211  }
212  fibre1LV.push_back(lv);
213  edm::LogInfo("HcalSim") << "HCalSD: (" << i << ") " << fibreNames[i]
214  << " LV " << fibre1LV[i];
215  }
216 
217  value = "HFFibreBundleConical";
218  DDSpecificsFilter filter5;
219  DDValue ddv5(attribute,value,0);
221  DDFilteredView fv5(cpv);
222  fv5.addFilter(filter5);
223  fibreNames = getNames(fv5);
224  edm::LogInfo("HcalSim") << "HCalSD: Names to be tested for " << attribute
225  << " = " << value << " have " << fibreNames.size()
226  << " entries";
227  for (unsigned int i=0; i<fibreNames.size(); i++) {
228  G4String namv = fibreNames[i];
229  lv = 0;
230  for (lvcite = lvs->begin(); lvcite != lvs->end(); lvcite++)
231  if ((*lvcite)->GetName() == namv) {
232  lv = (*lvcite);
233  break;
234  }
235  fibre2LV.push_back(lv);
236  edm::LogInfo("HcalSim") << "HCalSD: (" << i << ") " << fibreNames[i]
237  << " LV " << fibre2LV[i];
238  }
239  if (fibre1LV.size() > 0 || fibre2LV.size() > 0)
240  showerBundle = new HFShowerFibreBundle (name, cpv, p);
241  }
242 
243  //Material list for HB/HE/HO sensitive detectors
244  attribute = "ReadOutName";
245  DDSpecificsFilter filter2;
246  DDValue ddv2(attribute,name,0);
248  DDFilteredView fv2(cpv);
249  fv2.addFilter(filter2);
250  bool dodet = fv2.firstChild();
251 
253  //Layer0 Weight
254  layer0wt = getDDDArray("Layer0Wt",sv);
255  edm::LogInfo("HcalSim") << "HCalSD: " << layer0wt.size() << " Layer0Wt";
256  for (unsigned int it=0; it<layer0wt.size(); it++)
257  edm::LogInfo("HcalSim") << "HCalSD: [" << it << "] " << layer0wt[it];
258 
259  const G4MaterialTable * matTab = G4Material::GetMaterialTable();
260  std::vector<G4Material*>::const_iterator matite;
261  while (dodet) {
262  const DDLogicalPart & log = fv2.logicalPart();
263  G4String namx = log.name().name();
264  if (!isItHF(namx) && !isItFibre(namx)) {
265  bool notIn = true;
266  for (unsigned int i=0; i<matNames.size(); i++) {
267  if (!strcmp(matNames[i].c_str(),log.material().name().name().c_str())){
268  notIn = false;
269  break;
270  }
271  }
272  if (notIn) {
273  namx = log.material().name().name();
274  matNames.push_back(namx);
275  G4Material* mat;
276  for (matite = matTab->begin(); matite != matTab->end(); matite++) {
277  if ((*matite)->GetName() == namx) {
278  mat = (*matite);
279  break;
280  }
281  }
282  materials.push_back(mat);
283  }
284  }
285  dodet = fv2.next();
286  }
287 
288  edm::LogInfo("HcalSim") << "HCalSD: Material names for " << attribute
289  << " = " << name << ":";
290  for (unsigned int i=0; i<matNames.size(); i++)
291  edm::LogInfo("HcalSim") << "HCalSD: (" << i << ") " << matNames[i]
292  << " pointer " << materials[i];
293 
294  mumPDG = mupPDG = 0;
295 
296  if (useLayerWt) readWeightFromFile(file);
297 
298  for (int i=0; i<9; i++) {
299  hit_[i] = time_[i]= dist_[i] = 0;
300  }
301 
302 #ifdef DebugLog
304 
305  if ( tfile.isAvailable() ) {
306  static std::string labels[9] = {"HB", "HE", "HO", "HF Absorber", "HF PMT", "HF Absorber Long", "HF Absorber Short", "HF PMT Long", "HF PMT Short"};
307  char name[20], title[60];
308  for (int i=0; i<9; i++) {
309  sprintf (title, "Hit energy in %s", labels[i].c_str());
310  sprintf (name, "HCalSDHit%d", i);
311  hit_[i] = tfile->make<TH1F>(name, title, 2000, 0., 2000.);
312  sprintf (title, "Energy (MeV)");
313  hit_[i]->GetXaxis()->SetTitle(title);
314  hit_[i]->GetYaxis()->SetTitle("Hits");
315  sprintf (title, "Time of the hit in %s", labels[i].c_str());
316  sprintf (name, "HCalSDTime%d", i);
317  time_[i] = tfile->make<TH1F>(name, title, 2000, 0., 2000.);
318  sprintf (title, "Time (ns)");
319  time_[i]->GetXaxis()->SetTitle(title);
320  time_[i]->GetYaxis()->SetTitle("Hits");
321  sprintf (title, "Longitudinal profile in %s", labels[i].c_str());
322  sprintf (name, "HCalSDDist%d", i);
323  dist_[i] = tfile->make<TH1F>(name, title, 2000, 0., 2000.);
324  sprintf (title, "Distance (mm)");
325  dist_[i]->GetXaxis()->SetTitle(title);
326  dist_[i]->GetYaxis()->SetTitle("Hits");
327  }
328  }
329 #endif
330 }
331 
334  if (numberingScheme) delete numberingScheme;
335  if (showerLibrary) delete showerLibrary;
336  if (hfshower) delete hfshower;
337  if (showerParam) delete showerParam;
338  if (showerPMT) delete showerPMT;
339  if (showerBundle) delete showerBundle;
340 }
341 
342 bool HCalSD::ProcessHits(G4Step * aStep, G4TouchableHistory * ) {
343  // TimeMe t1( theHitTimer, false);
344 
345  NaNTrap( aStep ) ;
346 
347  if (aStep == NULL) {
348  return true;
349  } else {
350  G4LogicalVolume* lv =
351  aStep->GetPreStepPoint()->GetPhysicalVolume()->GetLogicalVolume();
352  G4String nameVolume = lv->GetName();
353  if (isItHF(aStep)) {
354  G4int parCode =aStep->GetTrack()->GetDefinition()->GetPDGEncoding();
355  if (useParam) {
356 #ifdef DebugLog
357  LogDebug("HcalSim") << "HCalSD: Hit from parametrization in "
358  << nameVolume << " for Track "
359  << aStep->GetTrack()->GetTrackID()
360  <<" (" << aStep->GetTrack()->GetDefinition()->GetParticleName() << ")";
361 #endif
362  getFromParam(aStep);
363  } else {
364  bool notaMuon = true;
365  if (parCode == mupPDG || parCode == mumPDG ) notaMuon = false;
366  if (useShowerLibrary && notaMuon) {
367 #ifdef DebugLog
368  LogDebug("HcalSim") << "HCalSD: Starts shower library from "
369  << nameVolume << " for Track "
370  << aStep->GetTrack()->GetTrackID()
371  <<" (" << aStep->GetTrack()->GetDefinition()->GetParticleName() << ")";
372 #endif
373  getFromLibrary(aStep);
374  } else if (isItFibre(lv)) {
375 #ifdef DebugLog
376  LogDebug("HcalSim") << "HCalSD: Hit at Fibre in " << nameVolume
377  << " for Track "
378  << aStep->GetTrack()->GetTrackID()
379  <<" (" << aStep->GetTrack()->GetDefinition()->GetParticleName() << ")";
380 #endif
381  hitForFibre(aStep);
382  }
383  }
384  } else if (isItPMT(lv)) {
385 #ifdef DebugLog
386  LogDebug("HcalSim") << "HCalSD: Hit from PMT parametrization from "
387  << nameVolume << " for Track "
388  << aStep->GetTrack()->GetTrackID() << " ("
389  << aStep->GetTrack()->GetDefinition()->GetParticleName() << ")";
390 #endif
391  if (usePMTHit && showerPMT) getHitPMT(aStep);
392  } else if (isItStraightBundle(lv) || isItConicalBundle(lv)) {
393 #ifdef DebugLog
394  LogDebug("HcalSim") << "HCalSD: Hit from FibreBundle from "
395  << nameVolume << " for Track "
396  << aStep->GetTrack()->GetTrackID() << " ("
397  << aStep->GetTrack()->GetDefinition()->GetParticleName() << ")";
398 #endif
401  } else {
402 #ifdef DebugLog
403  LogDebug("HcalSim") << "HCalSD: Hit from standard path from "
404  << nameVolume << " for Track "
405  << aStep->GetTrack()->GetTrackID() << " ("
406  << aStep->GetTrack()->GetDefinition()->GetParticleName() << ")";
407 #endif
408  if (getStepInfo(aStep)) {
409 #ifdef DebugLog
410  if (edepositEM+edepositHAD > 0)
411  plotProfile(aStep, aStep->GetPreStepPoint()->GetPosition(),
413  aStep->GetPostStepPoint()->GetGlobalTime(), 0);
414 #endif
415  if (hitExists() == false && edepositEM+edepositHAD>0.)
417  }
418  }
419  return true;
420  }
421 }
422 
423 double HCalSD::getEnergyDeposit(G4Step* aStep) {
424 
425  double destep = aStep->GetTotalEnergyDeposit();
426  double weight = 1;
427  const G4VTouchable* touch = aStep->GetPreStepPoint()->GetTouchable();
428  int depth = (touch->GetReplicaNumber(0))%10;
429  int det = ((touch->GetReplicaNumber(1))/1000)-3;
430  if (depth==0 && (det==0 || det==1)) weight = layer0wt[det];
431  if (useLayerWt) {
432  int lay = (touch->GetReplicaNumber(0)/10)%100 + 1;
433  G4ThreeVector hitPoint = aStep->GetPreStepPoint()->GetPosition();
434  weight = layerWeight(det+3, hitPoint, depth, lay);
435  }
436 
437  if (suppressHeavy) {
438  G4Track* theTrack = aStep->GetTrack();
439  TrackInformation * trkInfo = (TrackInformation *)(theTrack->GetUserInformation());
440  if (trkInfo) {
441  int pdg = theTrack->GetDefinition()->GetPDGEncoding();
442  if (!(trkInfo->isPrimary())) { // Only secondary particles
443  double ke = theTrack->GetKineticEnergy()/MeV;
444  if (((pdg/1000000000 == 1 && ((pdg/10000)%100) > 0 &&
445  ((pdg/10)%100) > 0)) && (ke<kmaxIon)) weight = 0;
446  if ((pdg == 2212) && (ke < kmaxProton)) weight = 0;
447  if ((pdg == 2112) && (ke < kmaxNeutron)) weight = 0;
448 #ifdef DebugLog
449  if (weight == 0)
450  edm::LogInfo("HcalSim") << "HCalSD:Ignore Track "
451  << theTrack->GetTrackID() << " Type "
452  << theTrack->GetDefinition()->GetParticleName()
453  << " Kinetic Energy " << ke << " MeV";
454 #endif
455  }
456  }
457  }
458 #ifdef DebugLog
459  double weight0 = weight;
460 #endif
461  if (useBirk) {
462  G4Material* mat = aStep->GetPreStepPoint()->GetMaterial();
463  if (isItScintillator(mat))
464  weight *= getAttenuation(aStep, birk1, birk2, birk3);
465  }
466  double wt1 = getResponseWt(theTrack);
467 #ifdef DebugLog
468  edm::LogInfo("HcalSim") << "HCalSD: Detector " << det+3 << " Depth " << depth
469  << " weight " << weight0 << " " << weight << " "
470  << wt1;
471 #endif
472  return weight*wt1*destep;
473 }
474 
475 uint32_t HCalSD::setDetUnitId(G4Step * aStep) {
476 
477  G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
478  const G4VTouchable* touch = preStepPoint->GetTouchable();
479  G4ThreeVector hitPoint = preStepPoint->GetPosition();
480 
481  int depth = (touch->GetReplicaNumber(0))%10 + 1;
482  int lay = (touch->GetReplicaNumber(0)/10)%100 + 1;
483  int det = (touch->GetReplicaNumber(1))/1000;
484 
485  return setDetUnitId (det, hitPoint, depth, lay);
486 }
487 
489  if (scheme != 0) {
490  edm::LogInfo("HcalSim") << "HCalSD: updates numbering scheme for "
491  << GetName() << "\n";
492  if (numberingScheme) delete numberingScheme;
493  numberingScheme = scheme;
494  }
495 }
496 
498 
499  G4ParticleTable * theParticleTable = G4ParticleTable::GetParticleTable();
500  G4String particleName;
501  mumPDG = theParticleTable->FindParticle(particleName="mu-")->GetPDGEncoding();
502  mupPDG = theParticleTable->FindParticle(particleName="mu+")->GetPDGEncoding();
503 #ifdef DebugLog
504  edm::LogInfo("HcalSim") << "HCalSD: Particle code for mu- = " << mumPDG
505  << " for mu+ = " << mupPDG;
506 #endif
507  if (showerLibrary) showerLibrary->initRun(theParticleTable);
508  if (showerParam) showerParam->initRun(theParticleTable);
509 }
510 
511 bool HCalSD::filterHit(CaloG4Hit* aHit, double time) {
512 
513  double threshold=0;
514  DetId theId(aHit->getUnitID());
515  switch (theId.subdetId()) {
516  case HcalBarrel:
517  threshold = eminHitHB; break;
518  case HcalEndcap:
519  threshold = eminHitHE; break;
520  case HcalOuter:
521  threshold = eminHitHO; break;
522  case HcalForward:
523  threshold = eminHitHF; break;
524  default:
525  break;
526  }
527  return ((time <= tmaxHit) && (aHit->getEnergyDeposit() > threshold));
528 }
529 
530 
531 uint32_t HCalSD::setDetUnitId (int det, G4ThreeVector pos, int depth,
532  int lay=1) {
533 
534  uint32_t id = 0;
535  if (numberingFromDDD) {
536  //get the ID's as eta, phi, depth, ... indices
538  depth, lay);
539  //get the ID
540  if (numberingScheme) id = numberingScheme->getUnitID(tmp);
541  }
542  return id;
543 }
544 
545 std::vector<double> HCalSD::getDDDArray(const std::string & str,
546  const DDsvalues_type & sv) {
547 
548 #ifdef DebugLog
549  LogDebug("HcalSim") << "HCalSD:getDDDArray called for " << str;
550 #endif
551  DDValue value(str);
552  if (DDfetch(&sv,value)) {
553 #ifdef DebugLog
554  LogDebug("HcalSim") << value;
555 #endif
556  const std::vector<double> & fvec = value.doubles();
557  int nval = fvec.size();
558  if (nval < 1) {
559  edm::LogError("HcalSim") << "HCalSD : # of " << str << " bins " << nval
560  << " < 2 ==> illegal";
561  throw cms::Exception("Unknown", "HCalSD")
562  << "nval < 2 for array " << str << "\n";
563  }
564 
565  return fvec;
566  } else {
567  edm::LogError("HcalSim") << "HCalSD : cannot get array " << str;
568  throw cms::Exception("Unknown", "HCalSD")
569  << "cannot get array " << str << "\n";
570  }
571 }
572 
573 std::vector<G4String> HCalSD::getNames(DDFilteredView& fv) {
574 
575  std::vector<G4String> tmp;
576  bool dodet = fv.firstChild();
577  while (dodet) {
578  const DDLogicalPart & log = fv.logicalPart();
579  bool ok = true;
580  for (unsigned int i=0; i<tmp.size(); i++) {
581  if (!strcmp(tmp[i].c_str(),log.name().name().c_str())) {
582  ok = false;
583  break;
584  }
585  }
586  if (ok) tmp.push_back(log.name().name());
587  dodet = fv.next();
588  }
589  return tmp;
590 }
591 
592 bool HCalSD::isItHF (G4Step * aStep) {
593 
594  const G4VTouchable* touch = aStep->GetPreStepPoint()->GetTouchable();
595  int levels = ((touch->GetHistoryDepth())+1);
596  for (unsigned int it=0; it < hfNames.size(); it++) {
597  if (levels >= hfLevels[it]) {
598  G4LogicalVolume* lv = touch->GetVolume(levels-hfLevels[it])->GetLogicalVolume();
599  if (lv == hfLV[it]) return true;
600  }
601  }
602  return false;
603 }
604 
605 bool HCalSD::isItHF (G4String name) {
606 
607  std::vector<G4String>::const_iterator it = hfNames.begin();
608  for (; it != hfNames.end(); it++)
609  if (name == *it) return true;
610  return false;
611 }
612 
613 bool HCalSD::isItFibre (G4LogicalVolume* lv) {
614 
615  std::vector<G4LogicalVolume*>::const_iterator ite = fibreLV.begin();
616  for (; ite != fibreLV.end(); ite++)
617  if (lv == *ite) return true;
618  return false;
619 }
620 
621 bool HCalSD::isItFibre (G4String name) {
622 
623  std::vector<G4String>::const_iterator it = fibreNames.begin();
624  for (; it != fibreNames.end(); it++)
625  if (name == *it) return true;
626  return false;
627 }
628 
629 bool HCalSD::isItPMT (G4LogicalVolume* lv) {
630 
631  std::vector<G4LogicalVolume*>::const_iterator ite = pmtLV.begin();
632  for (; ite != pmtLV.end(); ite++)
633  if (lv == *ite) return true;
634  return false;
635 }
636 
637 bool HCalSD::isItStraightBundle (G4LogicalVolume* lv) {
638 
639  std::vector<G4LogicalVolume*>::const_iterator ite = fibre1LV.begin();
640  for (; ite != fibre1LV.end(); ite++)
641  if (lv == *ite) return true;
642  return false;
643 }
644 
645 bool HCalSD::isItConicalBundle (G4LogicalVolume* lv) {
646 
647  std::vector<G4LogicalVolume*>::const_iterator ite = fibre2LV.begin();
648  for (; ite != fibre2LV.end(); ite++)
649  if (lv == *ite) return true;
650  return false;
651 }
652 
653 bool HCalSD::isItScintillator (G4Material* mat) {
654 
655  std::vector<G4Material*>::const_iterator ite = materials.begin();
656  for (; ite != materials.end(); ite++)
657  if (mat == *ite) return true;
658  return false;
659 }
660 
661 void HCalSD::getFromLibrary (G4Step* aStep) {
662 
663  preStepPoint = aStep->GetPreStepPoint();
664  theTrack = aStep->GetTrack();
665  int det = 5;
666  bool ok;
667 
668  std::vector<HFShowerLibrary::Hit> hits = showerLibrary->getHits(aStep, ok);
669 
670  double etrack = preStepPoint->GetKineticEnergy();
671  int primaryID = setTrackID(aStep);
672  /*
673  int primaryID = 0;
674  if (etrack >= energyCut) {
675  primaryID = theTrack->GetTrackID();
676  } else {
677  primaryID = theTrack->GetParentID();
678  if (primaryID == 0) primaryID = theTrack->GetTrackID();
679  }
680  */
681 
682  // Reset entry point for new primary
683  posGlobal = preStepPoint->GetPosition();
684  resetForNewPrimary(posGlobal, etrack);
685 
686  G4int particleCode = theTrack->GetDefinition()->GetPDGEncoding();
687  if (particleCode==emPDG || particleCode==epPDG || particleCode==gammaPDG) {
688  edepositEM = 1.*GeV; edepositHAD = 0.;
689  } else {
690  edepositEM = 0.; edepositHAD = 1.*GeV;
691  }
692 
693 #ifdef DebugLog
694  edm::LogInfo("HcalSim") << "HCalSD::getFromLibrary " <<hits.size()
695  << " hits for " << GetName() << " of " << primaryID
696  << " with " << theTrack->GetDefinition()->GetParticleName()
697  << " of " << preStepPoint->GetKineticEnergy()/GeV
698  << " GeV";
699 #endif
700  for (unsigned int i=0; i<hits.size(); i++) {
701  G4ThreeVector hitPoint = hits[i].position;
702  int depth = hits[i].depth;
703  double time = hits[i].time;
704  unsigned int unitID = setDetUnitId(det, hitPoint, depth);
705  currentID.setID(unitID, time, primaryID, 0);
706 #ifdef DebugLog
707  plotProfile(aStep, hitPoint, 1.0*GeV, time, depth);
708 #endif
709 
710  // check if it is in the same unit and timeslice as the previous one
711  if (currentID == previousID) {
713  } else {
714  if (!checkHit()) currentHit = createNewHit();
715  }
716  }
717 
718  //Now kill the current track
719  if (ok) {
720  theTrack->SetTrackStatus(fStopAndKill);
721  G4TrackVector tv = *(aStep->GetSecondary());
722  for (unsigned int kk=0; kk<tv.size(); kk++) {
723  if (tv[kk]->GetVolume() == preStepPoint->GetPhysicalVolume())
724  tv[kk]->SetTrackStatus(fStopAndKill);
725  }
726  }
727 }
728 
729 void HCalSD::hitForFibre (G4Step* aStep) {
730 
731  preStepPoint = aStep->GetPreStepPoint();
732  theTrack = aStep->GetTrack();
733  int primaryID = setTrackID(aStep);
734 
735  int det = 5;
736  std::vector<HFShower::Hit> hits = hfshower->getHits(aStep);
737 
738  G4int particleCode = theTrack->GetDefinition()->GetPDGEncoding();
739  if (particleCode==emPDG || particleCode==epPDG || particleCode==gammaPDG) {
740  edepositEM = 1.*GeV; edepositHAD = 0.;
741  } else {
742  edepositEM = 0.; edepositHAD = 1.*GeV;
743  }
744 
745 #ifdef DebugLog
746  edm::LogInfo("HcalSim") << "HCalSD::hitForFibre " << hits.size()
747  << " hits for " << GetName() << " of " << primaryID
748  << " with " << theTrack->GetDefinition()->GetParticleName()
749  << " of " << preStepPoint->GetKineticEnergy()/GeV
750  << " GeV in detector type " << det;
751 #endif
752  if (hits.size() > 0) {
753 
754  G4ThreeVector hitPoint = preStepPoint->GetPosition();
755  int depth =
756  (preStepPoint->GetTouchable()->GetReplicaNumber(0))%10;
757  unsigned int unitID = setDetUnitId(det, hitPoint, depth);
758 
759  for (unsigned int i=0; i<hits.size(); i++) {
760  double time = hits[i].time;
761  currentID.setID(unitID, time, primaryID, 0);
762 #ifdef DebugLog
763  plotProfile(aStep, hitPoint, 1.0*GeV, time, depth);
764 #endif
765  // check if it is in the same unit and timeslice as the previous one
766  if (currentID == previousID) {
768  } else {
769  posGlobal = preStepPoint->GetPosition();
770  if (!checkHit()) currentHit = createNewHit();
771  }
772  }
773  }
774 
775 }
776 
777 void HCalSD::getFromParam (G4Step* aStep) {
778 
779  std::vector<HFShowerParam::Hit> hits = showerParam->getHits(aStep);
780  int nHit = static_cast<int>(hits.size());
781 
782  if (nHit > 0) {
783  preStepPoint = aStep->GetPreStepPoint();
784  int primaryID = setTrackID(aStep);
785 
786  int det = 5;
787 #ifdef DebugLog
788  edm::LogInfo("HcalSim") << "HCalSD::getFromParam " << nHit << " hits for "
789  << GetName() << " of " << primaryID << " with "
790  << aStep->GetTrack()->GetDefinition()->GetParticleName()
791  << " of " << preStepPoint->GetKineticEnergy()/GeV
792  << " GeV in detector type " << det;
793 #endif
794  for (int i = 0; i<nHit; i++) {
795  G4ThreeVector hitPoint = hits[i].position;
796  int depth = hits[i].depth;
797  double time = hits[i].time;
798  unsigned int unitID = setDetUnitId(det, hitPoint, depth);
799  currentID.setID(unitID, time, primaryID, 0);
800  edepositEM = hits[i].edep*GeV;
801  edepositHAD = 0;
802 #ifdef DebugLog
803  plotProfile(aStep, hitPoint, edepositEM, time, depth);
804 #endif
805 
806  // check if it is in the same unit and timeslice as the previous one
807  if (currentID == previousID) {
809  } else {
810  posGlobal = preStepPoint->GetPosition();
811  if (!checkHit()) currentHit = createNewHit();
812  }
813  }
814  }
815 }
816 
817 void HCalSD::getHitPMT (G4Step* aStep) {
818 
819  preStepPoint = aStep->GetPreStepPoint();
820  theTrack = aStep->GetTrack();
821  double edep = showerPMT->getHits(aStep);
822 
823  if (edep >= 0) {
824  double etrack = preStepPoint->GetKineticEnergy();
825  int primaryID = 0;
826  if (etrack >= energyCut) {
827  primaryID = theTrack->GetTrackID();
828  } else {
829  primaryID = theTrack->GetParentID();
830  if (primaryID == 0) primaryID = theTrack->GetTrackID();
831  }
832  // Reset entry point for new primary
833  posGlobal = preStepPoint->GetPosition();
834  resetForNewPrimary(posGlobal, etrack);
835 
836  //
837  int det = static_cast<int>(HcalForward);
838  G4ThreeVector hitPoint = preStepPoint->GetPosition();
839  double rr = (hitPoint.x()*hitPoint.x() + hitPoint.y()*hitPoint.y());
840  double phi = (rr == 0. ? 0. :atan2(hitPoint.y(),hitPoint.x()));
841  double etaR = showerPMT->getRadius();
842  int depth = 1;
843  if (etaR < 0) {
844  depth = 2;
845  etaR =-etaR;
846  }
847  if (hitPoint.z() < 0) etaR =-etaR;
848 #ifdef DebugLog
849  edm::LogInfo("HcalSim") << "HCalSD::Hit for Detector " << det << " etaR "
850  << etaR << " phi " << phi/deg << " depth " <<depth;
851 #endif
852  double time = (aStep->GetPostStepPoint()->GetGlobalTime());
853  uint32_t unitID = 0;
854  if (numberingFromDDD) {
856  depth,1);
857  if (numberingScheme) unitID = numberingScheme->getUnitID(tmp);
858  }
859  currentID.setID(unitID, time, primaryID, 1);
860 
861  edepositHAD = aStep->GetTotalEnergyDeposit();
862  edepositEM =-edepositHAD + (edep*GeV);
863 #ifdef DebugLog
864  plotProfile(aStep, hitPoint, edep*GeV, time, depth);
865  double beta = preStepPoint->GetBeta();
866  LogDebug("HcalSim") << "HCalSD::getHitPMT 1 hit for " << GetName()
867  << " of " << primaryID << " with "
868  << theTrack->GetDefinition()->GetParticleName()
869  << " of " << preStepPoint->GetKineticEnergy()/GeV
870  << " GeV with velocity " << beta << " UnitID "
871  << std::hex << unitID << std::dec;
872 #endif
873  // check if it is in the same unit and timeslice as the previous one
874  if (currentID == previousID) {
876  } else {
877  if (!checkHit()) currentHit = createNewHit();
878  }
879  }
880 }
881 
882 void HCalSD::getHitFibreBundle (G4Step* aStep, bool type) {
883 
884  preStepPoint = aStep->GetPreStepPoint();
885  theTrack = aStep->GetTrack();
886  double edep = showerBundle->getHits(aStep, type);
887 
888  if (edep >= 0) {
889  double etrack = preStepPoint->GetKineticEnergy();
890  int primaryID = 0;
891  if (etrack >= energyCut) {
892  primaryID = theTrack->GetTrackID();
893  } else {
894  primaryID = theTrack->GetParentID();
895  if (primaryID == 0) primaryID = theTrack->GetTrackID();
896  }
897  // Reset entry point for new primary
898  posGlobal = preStepPoint->GetPosition();
899  resetForNewPrimary(posGlobal, etrack);
900 
901  //
902  int det = static_cast<int>(HcalForward);
903  G4ThreeVector hitPoint = preStepPoint->GetPosition();
904  double rr = (hitPoint.x()*hitPoint.x() + hitPoint.y()*hitPoint.y());
905  double phi = (rr == 0. ? 0. :atan2(hitPoint.y(),hitPoint.x()));
906  double etaR = showerBundle->getRadius();
907  int depth = 1;
908  if (etaR < 0) {
909  depth = 2;
910  etaR =-etaR;
911  }
912  if (hitPoint.z() < 0) etaR =-etaR;
913 #ifdef DebugLog
914  edm::LogInfo("HcalSim") << "HCalSD::Hit for Detector " << det << " etaR "
915  << etaR << " phi " << phi/deg << " depth " <<depth;
916 #endif
917  double time = (aStep->GetPostStepPoint()->GetGlobalTime());
918  uint32_t unitID = 0;
919  if (numberingFromDDD) {
921  depth,1);
922  if (numberingScheme) unitID = numberingScheme->getUnitID(tmp);
923  }
924  if (type) currentID.setID(unitID, time, primaryID, 3);
925  else currentID.setID(unitID, time, primaryID, 2);
926 
927  edepositHAD = aStep->GetTotalEnergyDeposit();
928  edepositEM =-edepositHAD + (edep*GeV);
929 #ifdef DebugLog
930  plotProfile(aStep, hitPoint, edep*GeV, time, depth);
931  double beta = preStepPoint->GetBeta();
932  LogDebug("HcalSim") << "HCalSD::getHitFibreBundle 1 hit for " << GetName()
933  << " of " << primaryID << " with "
934  << theTrack->GetDefinition()->GetParticleName()
935  << " of " << preStepPoint->GetKineticEnergy()/GeV
936  << " GeV with velocity " << beta << " UnitID "
937  << std::hex << unitID << std::dec;
938 #endif
939  // check if it is in the same unit and timeslice as the previous one
940  if (currentID == previousID) {
942  } else {
943  if (!checkHit()) currentHit = createNewHit();
944  }
945  } // non-zero energy deposit
946 }
947 
948 int HCalSD::setTrackID (G4Step* aStep) {
949 
950  theTrack = aStep->GetTrack();
951 
952  double etrack = preStepPoint->GetKineticEnergy();
953  TrackInformation * trkInfo = (TrackInformation *)(theTrack->GetUserInformation());
954  int primaryID = trkInfo->getIDonCaloSurface();
955  if (primaryID == 0) {
956 #ifdef DebugLog
957  edm::LogInfo("HcalSim") << "HCalSD: Problem with primaryID **** set by "
958  << "force to TkID **** " <<theTrack->GetTrackID();
959 #endif
960  primaryID = theTrack->GetTrackID();
961  }
962 
963  if (primaryID != previousID.trackID())
964  resetForNewPrimary(preStepPoint->GetPosition(), etrack);
965 
966  return primaryID;
967 }
968 
969 void HCalSD::readWeightFromFile(std::string fName) {
970 
971  std::ifstream infile;
972  int entry=0;
973  infile.open(fName.c_str(), std::ios::in);
974  if (infile) {
975  int det, zside, etaR, phi, lay;
976  double wt;
977  while (infile >> det >> zside >> etaR >> phi >> lay >> wt) {
978  uint32_t id = HcalTestNumbering::packHcalIndex(det,zside,1,etaR,phi,lay);
979  layerWeights.insert(std::pair<uint32_t,double>(id,wt));
980  entry++;
981 #ifdef DebugLog
982  edm::LogInfo("HcalSim") << "HCalSD::readWeightFromFile:Entry " << entry
983  << " ID " << std::hex << id << std::dec << " ("
984  << det << "/" << zside << "/1/" << etaR << "/"
985  << phi << "/" << lay << ") Weight " << wt;
986 #endif
987  }
988  infile.close();
989  }
990  edm::LogInfo("HcalSim") << "HCalSD::readWeightFromFile: reads " << entry
991  << " weights from " << fName;
992  if (entry <= 0) useLayerWt = false;
993 }
994 
995 double HCalSD::layerWeight(int det, G4ThreeVector pos, int depth, int lay) {
996 
997  double wt = 1.;
998  if (numberingFromDDD) {
999  //get the ID's as eta, phi, depth, ... indices
1001  depth, lay);
1002  uint32_t id = HcalTestNumbering::packHcalIndex(tmp.subdet, tmp.zside,
1003  1, tmp.etaR, tmp.phis,
1004  tmp.lay);
1005  std::map<uint32_t,double>::const_iterator ite = layerWeights.find(id);
1006  if (ite != layerWeights.end()) wt = ite->second;
1007 #ifdef DebugLog
1008  edm::LogInfo("HcalSim") << "HCalSD::layerWeight: ID " << std::hex << id
1009  << std::dec << " (" << tmp.subdet << "/"
1010  << tmp.zside << "/1/" << tmp.etaR << "/"
1011  << tmp.phis << "/" << tmp.lay << ") Weight " <<wt;
1012 #endif
1013  }
1014  return wt;
1015 }
1016 
1017 void HCalSD::plotProfile(G4Step* aStep, G4ThreeVector global, double edep,
1018  double time, int id) {
1019 
1020  const G4VTouchable* touch = aStep->GetPreStepPoint()->GetTouchable();
1021  static G4String modName[8] = {"HEModule", "HVQF", "HBModule", "MBAT", "MBBT", "MBBTC", "MBBT_R1P", "MBBT_R1M"};
1022  G4ThreeVector local;
1023  bool found=false;
1024  double depth=-2000;
1025  int idx = 4;
1026  for (int n=0; n<touch->GetHistoryDepth(); n++) {
1027  G4String name = touch->GetVolume(n)->GetName();
1028 #ifdef DebugLog
1029  LogDebug("HcalSim") << "plotProfile Depth " << n << " Name " << name;
1030 #endif
1031  for (unsigned int ii=0; ii<8; ii++) {
1032  if (name == modName[ii]) {
1033  found = true;
1034  int dn = touch->GetHistoryDepth() - n;
1035  local = touch->GetHistory()->GetTransform(dn).TransformPoint(global);
1036  if (ii == 0) {depth = local.z() - 4006.5; idx = 1;}
1037  else if (ii == 1) {depth = local.z() + 825.0; idx = 3;}
1038  else if (ii == 2) {depth = local.x() - 1775.; idx = 0;}
1039  else {depth = local.y() + 15.; idx = 2;}
1040  break;
1041  }
1042  }
1043  if (found) break;
1044  }
1045  if (!found) depth = std::abs(global.z()) - 11500;
1046 #ifdef DebugLog
1047  LogDebug("HcalSim") << "plotProfile Found " << found << " Global " << global << " Local " << local << " depth " << depth << " ID " << id << " EDEP " << edep << " Time " << time;
1048 #endif
1049  if (hit_[idx] != 0) hit_[idx]->Fill(edep);
1050  if (time_[idx] != 0) time_[idx]->Fill(time,edep);
1051  if (dist_[idx] != 0) dist_[idx]->Fill(depth,edep);
1052  int jd = 2*idx + id - 7;
1053  if (jd >= 0 && jd < 4) {
1054  jd += 5;
1055  if (hit_[jd] != 0) hit_[jd]->Fill(edep);
1056  if (time_[jd] != 0) time_[jd]->Fill(time,edep);
1057  if (dist_[jd] != 0) dist_[jd]->Fill(depth,edep);
1058  }
1059 }
#define LogDebug(id)
float edepositEM
Definition: CaloSD.h:119
const double beta
double energyCut
Definition: CaloSD.h:121
type
Definition: HCALResponse.h:22
T getParameter(std::string const &) const
std::vector< double > layer0wt
Definition: HCalSD.h:85
T getUntrackedParameter(std::string const &, T const &) const
int i
Definition: DBlmapReader.cc:9
bool useParam
Definition: HCalSD.h:82
double eminHitHE
Definition: HCalSD.h:83
HFShowerParam * showerParam
Definition: HCalSD.h:77
TH1F * time_[9]
Definition: HCalSD.h:95
G4int emPDG
Definition: CaloSD.h:134
double getRadius()
Definition: HFShowerPMT.cc:150
void getFromLibrary(G4Step *step)
Definition: HCalSD.cc:661
const DDLogicalPart & logicalPart() const
The logical-part of the current node in the filtered-view.
const std::vector< double > & doubles() const
a reference to the double-valued values stored in the given instance of DDValue
Definition: DDValue.cc:151
std::vector< Hit > getHits(G4Step *aStep, bool &ok, bool onlyLong=false)
double kmaxNeutron
Definition: CaloSD.h:132
virtual uint32_t getUnitID(const HcalNumberingFromDDD::HcalID id)
void addFilter(const DDFilter &, log_op op=AND)
void updateHit(CaloG4Hit *)
Definition: CaloSD.cc:453
const N & name() const
Definition: DDBase.h:88
std::vector< Hit > getHits(G4Step *aStep)
int getIDonCaloSurface() const
bool useLayerWt
Definition: HCalSD.h:80
Definition: CaloSD.h:42
bool isItStraightBundle(G4LogicalVolume *)
Definition: HCalSD.cc:637
virtual bool ProcessHits(G4Step *, G4TouchableHistory *)
Definition: HCalSD.cc:342
bool useFibreBundle
Definition: HCalSD.h:80
double betaThr
Definition: HCalSD.h:81
HFShowerFibreBundle * showerBundle
Definition: HCalSD.h:79
std::vector< G4LogicalVolume * > fibre1LV
Definition: HCalSD.h:93
std::vector< G4LogicalVolume * > fibreLV
Definition: HCalSD.h:89
double eminHitHB
Definition: HCalSD.h:83
bool useShowerLibrary
Definition: HCalSD.h:82
std::vector< double > getDDDArray(const std::string &, const DDsvalues_type &)
Definition: HCalSD.cc:545
#define abs(x)
Definition: mlp_lapack.h:159
double getHits(G4Step *aStep, bool type)
#define NULL
Definition: scimark2.h:8
virtual uint32_t setDetUnitId(G4Step *step)
Definition: HCalSD.cc:475
double birk2
Definition: HCalSD.h:81
int setTrackID(G4Step *step)
Definition: HCalSD.cc:948
bool usePMTHit
Definition: HCalSD.h:82
G4ThreeVector posGlobal
Definition: CaloSD.h:111
void setNumberingScheme(HcalNumberingScheme *)
Definition: HCalSD.cc:488
type of data representation of DDCompactView
Definition: DDCompactView.h:81
double birk1
Definition: HCalSD.h:81
double kmaxProton
Definition: CaloSD.h:132
bool DDfetch(const DDsvalues_type *, DDValue &)
helper for retrieving DDValues from DDsvalues_type *.
Definition: DDsvalues.cc:102
const DDMaterial & material() const
Returns a reference object of the material this LogicalPart is made of.
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
bool isItHF(G4Step *)
Definition: HCalSD.cc:592
G4bool checkHit()
Definition: CaloSD.cc:343
static TrackerG4SimHitNumberingScheme & numberingScheme(const DDCompactView &cpv, const GeometricDet &det)
double eminHitHF
Definition: HCalSD.h:83
std::vector< Hit > getHits(G4Step *aStep)
Definition: HFShower.cc:69
bool isItConicalBundle(G4LogicalVolume *)
Definition: HCalSD.cc:645
static uint32_t packHcalIndex(int det, int z, int depth, int eta, int phi, int lay)
TH1F * hit_[9]
Definition: HCalSD.h:95
virtual void initRun()
Definition: HCalSD.cc:497
std::vector< G4LogicalVolume * > pmtLV
Definition: HCalSD.h:93
std::vector< G4LogicalVolume * > hfLV
Definition: HCalSD.h:86
double kmaxIon
Definition: CaloSD.h:132
bool suppressHeavy
Definition: CaloSD.h:130
bool isItScintillator(G4Material *)
Definition: HCalSD.cc:653
bool next()
set current node to the next node in the filtered tree
HFShower * hfshower
Definition: HCalSD.h:76
bool useBirk
Definition: HCalSD.h:80
float edepositHAD
Definition: CaloSD.h:119
double birk3
Definition: HCalSD.h:81
int trackID() const
Definition: CaloHitID.h:25
virtual double getEnergyDeposit(G4Step *)
Definition: HCalSD.cc:423
std::pair< std::string, MonitorElement * > entry
Definition: ME_MAP.h:8
G4int epPDG
Definition: CaloSD.h:134
std::vector< std::pair< unsigned int, DDValue > > DDsvalues_type
std::maps an index to a DDValue. The index corresponds to the index assigned to the name of the std::...
Definition: DDsvalues.h:19
G4int gammaPDG
Definition: CaloSD.h:134
bool isAvailable() const
Definition: Service.h:47
TH1F * dist_[9]
Definition: HCalSD.h:95
virtual ~HCalSD()
Definition: HCalSD.cc:332
A DDLogicalPart aggregates information concerning material, solid and sensitveness ...
Definition: DDLogicalPart.h:95
HcalNumberingScheme * numberingScheme
Definition: HCalSD.h:74
void getHitPMT(G4Step *step)
Definition: HCalSD.cc:817
CaloHitID previousID
Definition: CaloSD.h:115
double getAttenuation(G4Step *aStep, double birk1, double birk2, double birk3)
Definition: CaloSD.cc:480
CaloG4Hit * currentHit
Definition: CaloSD.h:126
std::vector< int > hfLevels
Definition: HCalSD.h:88
virtual std::vector< std::string > getNames()
void NaNTrap(G4Step *step)
G4Track * theTrack
Definition: CaloSD.h:116
virtual G4bool getStepInfo(G4Step *aStep)
Definition: CaloSD.cc:267
bool isItPMT(G4LogicalVolume *)
Definition: HCalSD.cc:629
G4int mupPDG
Definition: HCalSD.h:84
double eminHitHO
Definition: HCalSD.h:83
double tmaxHit
Definition: CaloSD.h:121
void setID(uint32_t unitID, double timeSlice, int trackID, uint16_t depth=0)
Definition: CaloHitID.cc:44
G4StepPoint * preStepPoint
Definition: CaloSD.h:118
tuple labels
Definition: L1TDQM_cfg.py:62
Definition: DetId.h:20
static const G4LogicalVolume * GetVolume(const std::string &name)
Log< T >::type log(const T &t)
Definition: Log.h:22
CaloHitID currentID
Definition: CaloSD.h:115
void readWeightFromFile(std::string)
Definition: HCalSD.cc:969
void plotProfile(G4Step *step, G4ThreeVector pos, double edep, double time, int id)
Definition: HCalSD.cc:1017
void getHitFibreBundle(G4Step *step, bool type)
Definition: HCalSD.cc:882
bool isPrimary() const
DDsvalues_type mergedSpecifics() const
void hitForFibre(G4Step *step)
Definition: HCalSD.cc:729
double getHits(G4Step *aStep)
Definition: HFShowerPMT.cc:104
void initRun(G4ParticleTable *theParticleTable)
HcalNumberingFromDDD * numberingFromDDD
Definition: HCalSD.h:73
std::vector< std::vector< double > > tmp
Definition: MVATrainer.cc:100
double getResponseWt(G4Track *)
Definition: CaloSD.cc:629
void initRun(G4ParticleTable *)
std::vector< G4String > fibreNames
Definition: HCalSD.h:90
list infile
Definition: EdgesToViz.py:90
G4bool hitExists()
Definition: CaloSD.cc:319
bool firstChild()
set the current node to the first child ...
std::vector< G4Material * > materials
Definition: HCalSD.h:92
virtual bool filterHit(CaloG4Hit *, double)
Definition: HCalSD.cc:511
T * make() const
make new ROOT object
void setCriteria(const DDValue &nameVal, comp_op, log_op l=AND, bool asString=true, bool merged=true)
Definition: DDFilter.cc:285
double layerWeight(int, G4ThreeVector, int, int)
Definition: HCalSD.cc:995
std::vector< G4String > hfNames
Definition: HCalSD.h:87
G4int mumPDG
Definition: HCalSD.h:84
void resetForNewPrimary(G4ThreeVector, double)
Definition: CaloSD.cc:468
std::map< uint32_t, double > layerWeights
Definition: HCalSD.h:94
tuple level
Definition: testEve_cfg.py:34
HcalID unitID(int det, CLHEP::Hep3Vector pos, int depth, int lay=-1) const
HFShowerPMT * showerPMT
Definition: HCalSD.h:78
uint32_t getUnitID() const
Definition: CaloG4Hit.h:69
HCalSD(G4String, const DDCompactView &, SensitiveDetectorCatalog &, edm::ParameterSet const &, const SimTrackManager *)
Definition: HCalSD.cc:32
HFShowerLibrary * showerLibrary
Definition: HCalSD.h:75
std::vector< G4LogicalVolume * > fibre2LV
Definition: HCalSD.h:93
const std::string & name() const
Returns the name.
Definition: DDName.cc:87
void getFromParam(G4Step *step)
Definition: HCalSD.cc:777
CaloG4Hit * createNewHit()
Definition: CaloSD.cc:375
std::vector< G4String > matNames
Definition: HCalSD.h:91
bool isItFibre(G4LogicalVolume *)
Definition: HCalSD.cc:613
double getEnergyDeposit() const
Definition: CaloG4Hit.h:81
bool useHF
Definition: HCalSD.h:82
Definition: DDAxes.h:10
The DDGenericFilter is a runtime-parametrized Filter looking on DDSpecifcs.
Definition: DDFilter.h:37