Geant4 Cross Reference |
1 // 1 2 // ******************************************* 3 // * License and Disclaimer 4 // * 5 // * The Geant4 software is copyright of th 6 // * the Geant4 Collaboration. It is provided 7 // * conditions of the Geant4 Software License 8 // * LICENSE and available at http://cern.ch/ 9 // * include a list of copyright holders. 10 // * 11 // * Neither the authors of this software syst 12 // * institutes,nor the agencies providing fin 13 // * work make any representation or warran 14 // * regarding this software system or assum 15 // * use. Please see the license in the file 16 // * for the full disclaimer and the limitatio 17 // * 18 // * This code implementation is the result 19 // * technical work of the GEANT4 collaboratio 20 // * By using, copying, modifying or distri 21 // * any work based on the software) you ag 22 // * use in resulting scientific publicati 23 // * acceptance of all terms of the Geant4 Sof 24 // ******************************************* 25 // 26 /// \file Par01/src/Par01EMShowerModel.cc 27 /// \brief Implementation of the Par01EMShower 28 // 29 // 30 // 31 #include "Par01EMShowerModel.hh" 32 33 #include "Par01EnergySpot.hh" 34 35 #include "G4Electron.hh" 36 #include "G4Gamma.hh" 37 #include "G4NistManager.hh" 38 #include "G4PhysicalConstants.hh" 39 #include "G4Positron.hh" 40 #include "G4SystemOfUnits.hh" 41 #include "G4TouchableHandle.hh" 42 #include "G4TransportationManager.hh" 43 #include "G4VSensitiveDetector.hh" 44 #include "Randomize.hh" 45 46 //....oooOO0OOooo........oooOO0OOooo........oo 47 48 Par01EMShowerModel::Par01EMShowerModel(G4Strin 49 : G4VFastSimulationModel(modelName, envelope 50 { 51 fFakeStep = new G4Step(); 52 fFakePreStepPoint = fFakeStep->GetPreStepPoi 53 fFakePostStepPoint = fFakeStep->GetPostStepP 54 fTouchableHandle = new G4TouchableHistory(); 55 fpNavigator = new G4Navigator(); 56 fNaviSetup = false; 57 fCsI = nullptr; 58 } 59 60 //....oooOO0OOooo........oooOO0OOooo........oo 61 62 Par01EMShowerModel::Par01EMShowerModel(G4Strin 63 { 64 fFakeStep = new G4Step(); 65 fFakePreStepPoint = fFakeStep->GetPreStepPoi 66 fFakePostStepPoint = fFakeStep->GetPostStepP 67 fTouchableHandle = new G4TouchableHistory(); 68 fpNavigator = new G4Navigator(); 69 fNaviSetup = false; 70 fCsI = nullptr; 71 } 72 73 //....oooOO0OOooo........oooOO0OOooo........oo 74 75 Par01EMShowerModel::~Par01EMShowerModel() 76 { 77 delete fFakeStep; 78 delete fpNavigator; 79 } 80 81 //....oooOO0OOooo........oooOO0OOooo........oo 82 83 G4bool Par01EMShowerModel::IsApplicable(const 84 { 85 return &particleType == G4Electron::Electron 86 || &particleType == G4Positron::Posit 87 || &particleType == G4Gamma::GammaDef 88 } 89 90 //....oooOO0OOooo........oooOO0OOooo........oo 91 92 G4bool Par01EMShowerModel::ModelTrigger(const 93 { 94 // Applies the parameterisation above 100 Me 95 return fastTrack.GetPrimaryTrack()->GetKinet 96 } 97 98 //....oooOO0OOooo........oooOO0OOooo........oo 99 100 void Par01EMShowerModel::DoIt(const G4FastTrac 101 { 102 // Kill the parameterised particle: 103 fastStep.KillPrimaryTrack(); 104 fastStep.ProposePrimaryTrackPathLength(0.0); 105 fastStep.ProposeTotalEnergyDeposited(fastTra 106 107 // split into "energy spots" energy accordin 108 Explode(fastTrack); 109 110 // and put those energy spots into the cryst 111 BuildDetectorResponse(); 112 } 113 114 //....oooOO0OOooo........oooOO0OOooo........oo 115 116 void Par01EMShowerModel::Explode(const G4FastT 117 { 118 //------------------------------------------ 119 // 120 //------------------------------------------ 121 122 // Reduced quantities: 123 // -- critical energy in CsI: 124 G4double Ec = 800 * MeV / (54. + 1.2); // 5 125 G4double Energy = fastTrack.GetPrimaryTrack( 126 G4double y = Energy / Ec; 127 128 // compute value of parameter "a" of longitu 129 G4double a, tmax, b(0.5), C; 130 if (fastTrack.GetPrimaryTrack()->GetDefiniti 131 C = 0.5; 132 else 133 C = -0.5; 134 tmax = 1.0 * (std::log(y) + C); 135 a = 1.0 + b * tmax; 136 137 // t : reduced quantity = z/X0: 138 G4double t, bt; 139 if (fCsI == nullptr) fCsI = G4NistManager::I 140 G4double X0 = fCsI->GetRadlen(); 141 // Moliere radius: 142 G4double Es = 21 * MeV; 143 G4double Rm = X0 * Es / Ec; 144 145 // axis of the shower, in global reference f 146 G4ThreeVector xShower, yShower, zShower; 147 zShower = fastTrack.GetPrimaryTrack()->GetMo 148 xShower = zShower.orthogonal(); 149 yShower = zShower.cross(xShower); 150 // starting point of the shower: 151 G4ThreeVector sShower = fastTrack.GetPrimary 152 153 // We shoot 100 spots of energy: 154 G4int nSpots = 100; 155 G4double deposit = Energy / double(nSpots); 156 Par01EnergySpot eSpot; 157 eSpot.SetEnergy(deposit); 158 G4ThreeVector ePoint; 159 G4double z, r, phi; 160 161 feSpotList.clear(); 162 for (int i = 0; i < nSpots; i++) { 163 // Longitudinal profile: 164 // -- shoot z according to Gamma distribut 165 bt = G4RandGamma::shoot(a, 1.0); 166 t = bt / b; 167 z = t * X0; 168 169 // transverse profile: 170 // we set 90% of energy in one Rm, 171 // the rest between 1 and 3.5 Rm: 172 G4double xr = G4UniformRand(); 173 if (xr < 0.9) 174 r = xr / 0.9 * Rm; 175 else 176 r = ((xr - 0.9) / 0.1 * 2.5 + 1.0) * Rm; 177 phi = G4UniformRand() * twopi; 178 179 // build the position: 180 ePoint = sShower + z * zShower + r * std:: 181 182 // and the energy spot: 183 eSpot.SetPosition(ePoint); 184 185 // Records the eSpot: 186 feSpotList.push_back(eSpot); 187 } 188 } 189 190 //....oooOO0OOooo........oooOO0OOooo........oo 191 192 void Par01EMShowerModel::BuildDetectorResponse 193 { 194 // Does the assignation of the energy spots 195 for (size_t i = 0; i < feSpotList.size(); i+ 196 // Draw the energy spot: 197 // feSpotList[i].Draw(); 198 // feSpotList[i].Print(); 199 200 // "converts" the energy spot into the fak 201 // G4Step to pass to sensitive detector: 202 AssignSpotAndCallHit(feSpotList[i]); 203 } 204 } 205 206 //....oooOO0OOooo........oooOO0OOooo........oo 207 208 void Par01EMShowerModel::AssignSpotAndCallHit( 209 { 210 // 211 // "converts" the energy spot into the fake 212 // G4Step to pass to sensitive detector: 213 // 214 FillFakeStep(eSpot); 215 216 // 217 // call sensitive part: taken/adapted from t 218 // Send G4Step information to Hit/Dig if the 219 // 220 G4VPhysicalVolume* pCurrentVolume = fFakeSte 221 G4VSensitiveDetector* pSensitive; 222 223 if (pCurrentVolume != nullptr) { 224 pSensitive = pCurrentVolume->GetLogicalVol 225 if (pSensitive != nullptr) { 226 pSensitive->Hit(fFakeStep); 227 } 228 } 229 } 230 231 //....oooOO0OOooo........oooOO0OOooo........oo 232 233 void Par01EMShowerModel::FillFakeStep(const Pa 234 { 235 //------------------------------------------ 236 // find in which volume the spot is. 237 //------------------------------------------ 238 if (!fNaviSetup) { 239 fpNavigator->SetWorldVolume(G4Transportati 240 ->GetNavigat 241 ->GetWorldVo 242 fpNavigator->LocateGlobalPointAndUpdateTou 243 eSpot.GetPosition(), G4ThreeVector(0., 0 244 fNaviSetup = true; 245 } 246 else { 247 fpNavigator->LocateGlobalPointAndUpdateTou 248 eSpot.GetPosition(), G4ThreeVector(0., 0 249 } 250 //-------------------------------------- 251 // Fills attribute of the G4Step needed 252 // by our sensitive detector: 253 //------------------------------------- 254 // set touchable volume at PreStepPoint: 255 fFakePreStepPoint->SetTouchableHandle(fTouch 256 // set total energy deposit: 257 fFakeStep->SetTotalEnergyDeposit(eSpot.GetEn 258 } 259