File indexing completed on 2021-07-20 22:46:48
0001 #include "SimG4CMS/ShowerLibraryProducer/interface/FiberSD.h"
0002 #include "SimDataFormats/CaloHit/interface/HFShowerPhoton.h"
0003 #include "DataFormats/Math/interface/Point3D.h"
0004 #include "Geometry/HcalCommonData/interface/HcalDDDSimConstants.h"
0005 #include "Geometry/HcalCommonData/interface/HcalSimulationConstants.h"
0006 #include "Geometry/Records/interface/HcalSimNumberingRecord.h"
0007
0008 #include "G4VPhysicalVolume.hh"
0009 #include "G4PVPlacement.hh"
0010 #include "G4HCofThisEvent.hh"
0011 #include "G4TouchableHistory.hh"
0012 #include "G4Track.hh"
0013 #include "G4Step.hh"
0014 #include "G4VSolid.hh"
0015 #include "G4DynamicParticle.hh"
0016 #include "G4ParticleDefinition.hh"
0017 #include "G4SDManager.hh"
0018 #include "G4ios.hh"
0019
0020
0021
0022 FiberSD::FiberSD(const std::string& iname,
0023 const HcalSimulationConstants* hsps,
0024 const HcalDDDSimConstants* hdc,
0025 const SensitiveDetectorCatalog& clg,
0026 edm::ParameterSet const& p,
0027 const SimTrackManager* manager)
0028 : SensitiveCaloDetector(iname, clg), theShower(nullptr), theHCID(-1), theHC(nullptr) {
0029 edm::LogVerbatim("FiberSim") << "FiberSD : Instantiating for " << iname;
0030
0031 theShower = new HFShower(iname, hdc, hsps->hcalsimpar(), p, 1);
0032 }
0033
0034 FiberSD::~FiberSD() {
0035 delete theShower;
0036 delete theHC;
0037 }
0038
0039 void FiberSD::Initialize(G4HCofThisEvent* HCE) {
0040 edm::LogVerbatim("FiberSim") << "FiberSD : Initialize called for " << GetName() << " in collection "
0041 << HCE;
0042 theHC = new FiberG4HitsCollection(GetName(), collectionName[0]);
0043 if (theHCID < 0)
0044 theHCID = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]);
0045 HCE->AddHitsCollection(theHCID, theHC);
0046 edm::LogVerbatim("FiberSim") << "FiberSD : Add hit collectrion for " << collectionName[0] << ":"
0047 << theHCID << ":" << theHC;
0048 }
0049
0050 G4bool FiberSD::ProcessHits(G4Step* aStep, G4TouchableHistory*) {
0051
0052 double zoffset = 1000;
0053 std::vector<HFShower::Hit> hits = theShower->getHits(aStep, true, zoffset);
0054
0055 if (!hits.empty()) {
0056 std::vector<HFShowerPhoton> thePE;
0057 for (unsigned int i = 0; i < hits.size(); i++) {
0058
0059 HFShowerPhoton pe = HFShowerPhoton(
0060 hits[i].position.x(), hits[i].position.y(), hits[i].position.z(), hits[i].wavelength, hits[i].time);
0061 thePE.push_back(pe);
0062 }
0063 int trackID = aStep->GetTrack()->GetTrackID();
0064 G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
0065 const G4VTouchable* touch = preStepPoint->GetTouchable();
0066 G4LogicalVolume* lv = touch->GetVolume(0)->GetLogicalVolume();
0067 int depth = (touch->GetReplicaNumber(0)) % 10;
0068 int detID = setDetUnitId(aStep);
0069 math::XYZPoint theHitPos(
0070 preStepPoint->GetPosition().x(), preStepPoint->GetPosition().y(), preStepPoint->GetPosition().z());
0071
0072
0073 FiberG4Hit* aHit = new FiberG4Hit(lv, detID, depth, trackID);
0074 #ifdef EDM_ML_DEBUG
0075 edm::LogVerbatim("FiberSim") << "FiberSD :hit size " << hits.size() << " npe" << aHit->npe();
0076 edm::LogVerbatim("FiberSim") << "FiberSD :pre hit position " << aHit->hitPos();
0077 #endif
0078 aHit->setNpe(hits.size());
0079 aHit->setPos(theHitPos);
0080 aHit->setTime(preStepPoint->GetGlobalTime());
0081 aHit->setPhoton(thePE);
0082 #ifdef EDM_ML_DEBUG
0083 edm::LogVerbatim("FiberSim") << "FiberSD :ShowerPhoton position " << thePE[0].x() << " "
0084 << thePE[0].y() << " " << thePE[0].z();
0085
0086 edm::LogVerbatim("FiberSim") << "FiberSD: Hit created at " << lv->GetName()
0087 << " DetID: " << aHit->towerId() << " Depth: " << aHit->depth()
0088 << " Track ID: " << aHit->trackId() << " Nb. of Cerenkov Photons: " << aHit->npe()
0089 << " Time: " << aHit->time() << " at " << aHit->hitPos();
0090 for (unsigned int i = 0; i < thePE.size(); i++)
0091 edm::LogVerbatim("FiberSim") << "FiberSD: PE[" << i << "] " << thePE[i];
0092 #endif
0093 theHC->insert(aHit);
0094 }
0095 return true;
0096 }
0097
0098 void FiberSD::EndOfEvent(G4HCofThisEvent* HCE) {
0099 edm::LogVerbatim("FiberSim") << "FiberSD: finds " << theHC->entries() << " hits";
0100 clear();
0101 edm::LogVerbatim("FiberSim") << "theHC entries = " << theHC->entries();
0102 }
0103
0104 void FiberSD::clear() {}
0105
0106 void FiberSD::DrawAll() {}
0107
0108 void FiberSD::PrintAll() {}
0109
0110 void FiberSD::update(const BeginOfJob* job) {}
0111
0112 void FiberSD::update(const BeginOfRun*) {}
0113
0114 void FiberSD::update(const BeginOfEvent*) {}
0115
0116 void FiberSD::update(const ::EndOfEvent*) {}
0117
0118 void FiberSD::clearHits() {}
0119
0120 uint32_t FiberSD::setDetUnitId(const G4Step* aStep) {
0121 const G4VTouchable* touch = aStep->GetPreStepPoint()->GetTouchable();
0122 int fibre = (touch->GetReplicaNumber(1)) % 10;
0123 int cell = (touch->GetReplicaNumber(2));
0124 int tower = (touch->GetReplicaNumber(3));
0125 return ((tower * 1000 + cell) * 10 + fibre);
0126 }
0127
0128 void FiberSD::fillHits(edm::PCaloHitContainer&, const std::string&) {}