Geant4 Cross Reference |
1 // 2 // ******************************************************************** 3 // * License and Disclaimer * 4 // * * 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. * 10 // * * 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitation of liability. * 17 // * * 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************************************** 25 // 26 // ******************************************************************** 27 // 28 // CaTS (Calorimetry and Tracking Simulation) 29 // 30 // Authors : Hans Wenzel 31 // Soon Yung Jun 32 // (Fermi National Accelerator Laboratory) 33 // 34 // History 35 // October 18th, 2021 : first implementation 36 // 37 // ******************************************************************** 38 // 39 /// \file PhotonSD.cc 40 /// \brief Implementation of the CaTS::PhotonSD class 41 42 // Geant4 headers 43 #include "G4VProcess.hh" 44 #include "G4OpticalPhoton.hh" 45 #include "G4HCofThisEvent.hh" 46 #include "G4Step.hh" 47 #include "G4SDManager.hh" 48 #include "ConfigurationManager.hh" 49 // project headers 50 #include "PhotonSD.hh" 51 #ifdef WITH_G4OPTICKS 52 # include "G4Opticks.hh" 53 # include "TrackInfo.hh" 54 # include "OpticksGenstep.h" 55 # include "OpticksFlags.hh" 56 # include "G4OpticksHit.hh" 57 #endif 58 59 PhotonSD::PhotonSD(G4String name) 60 : G4VSensitiveDetector(name) 61 { 62 G4String HCname = name + "_HC"; 63 collectionName.insert(HCname); 64 G4cout << collectionName.size() << " PhotonSD name: " << name 65 << " collection Name: " << HCname << G4endl; 66 fHCID = -1; 67 verbose = ConfigurationManager::getInstance()->isEnable_verbose(); 68 } 69 70 void PhotonSD::Initialize(G4HCofThisEvent* hce) 71 { 72 fPhotonHitsCollection = 73 new PhotonHitsCollection(SensitiveDetectorName, collectionName[0]); 74 if(fHCID < 0) 75 { 76 if(verbose) 77 G4cout << "PhotonSD::Initialize: " << SensitiveDetectorName << " " 78 << collectionName[0] << G4endl; 79 fHCID = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]); 80 } 81 hce->AddHitsCollection(fHCID, fPhotonHitsCollection); 82 } 83 84 G4bool PhotonSD::ProcessHits(G4Step* aStep, G4TouchableHistory*) 85 { 86 G4Track* theTrack = aStep->GetTrack(); 87 // we only deal with optical Photons: 88 if(theTrack->GetDefinition() != G4OpticalPhoton::OpticalPhotonDefinition()) 89 { 90 return false; 91 } 92 G4double theEdep = theTrack->GetTotalEnergy() / CLHEP::eV; 93 const G4VProcess* thisProcess = theTrack->GetCreatorProcess(); 94 G4String processname; 95 if(thisProcess != NULL) 96 processname = thisProcess->GetProcessName(); 97 else 98 processname = "No Process"; 99 unsigned theCreationProcessid; 100 if(processname == "Cerenkov") 101 { 102 theCreationProcessid = 0; 103 } 104 else if(processname == "Scintillation") 105 { 106 theCreationProcessid = 1; 107 } 108 else 109 { 110 theCreationProcessid = -1; 111 } 112 PhotonHit* newHit = new PhotonHit( 113 0, theCreationProcessid, etolambda(theEdep), theTrack->GetGlobalTime(), 114 aStep->GetPostStepPoint()->GetPosition(), 115 aStep->GetPostStepPoint()->GetMomentumDirection(), 116 aStep->GetPostStepPoint()->GetPolarization()); 117 fPhotonHitsCollection->insert(newHit); 118 theTrack->SetTrackStatus(fStopAndKill); 119 return true; 120 } 121 122 void PhotonSD::EndOfEvent(G4HCofThisEvent*) 123 { 124 if(verbose) 125 { 126 G4int NbHits = fPhotonHitsCollection->entries(); 127 G4cout << " PhotonSD::EndOfEvent Number of PhotonHits: " << NbHits 128 << G4endl; 129 } 130 } 131 #ifdef WITH_G4OPTICKS 132 void PhotonSD::AddOpticksHits() 133 { 134 G4Opticks* g4ok = G4Opticks::Get(); 135 bool way_enabled = g4ok->isWayEnabled(); 136 unsigned num_hits = g4ok->getNumHit(); 137 if(verbose) 138 G4cout << "PhotonSD::AddOpticksHits PhotonHits: " << num_hits << G4endl; 139 G4OpticksHit hit; 140 G4OpticksHitExtra hit_extra; 141 G4OpticksHitExtra* hit_extra_ptr = way_enabled ? &hit_extra : NULL; 142 for(unsigned i = 0; i < num_hits; i++) 143 { 144 g4ok->getHit(i, &hit, hit_extra_ptr); 145 PhotonHit* newHit = 146 new PhotonHit(i, 0, hit.wavelength, hit.time, hit.global_position, 147 hit.global_direction, hit.global_polarization); 148 fPhotonHitsCollection->insert(newHit); 149 } 150 if(verbose) 151 G4cout << "AddOpticksHits size: " << fPhotonHitsCollection->entries() 152 << G4endl; 153 } 154 #endif 155