Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer << 3 // * DISCLAIMER * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th << 5 // * The following disclaimer summarizes all the specific disclaimers * 6 // * the Geant4 Collaboration. It is provided << 6 // * of contributors to this software. The specific disclaimers,which * 7 // * conditions of the Geant4 Software License << 7 // * govern, are listed with their locations in: * 8 // * LICENSE and available at http://cern.ch/ << 8 // * http://cern.ch/geant4/license * 9 // * include a list of copyright holders. << 10 // * 9 // * * 11 // * Neither the authors of this software syst 10 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 11 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 12 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 13 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file << 14 // * use. * 16 // * for the full disclaimer and the limitatio << 17 // * 15 // * * 18 // * This code implementation is the result << 16 // * This code implementation is the intellectual property of the * 19 // * technical work of the GEANT4 collaboratio << 17 // * GEANT4 collaboration. * 20 // * By using, copying, modifying or distri << 18 // * By copying, distributing or modifying the Program (or any work * 21 // * any work based on the software) you ag << 19 // * based on the Program) you indicate your acceptance of this * 22 // * use in resulting scientific publicati << 20 // * statement, and all its terms. * 23 // * acceptance of all terms of the Geant4 Sof << 24 // ******************************************* 21 // ******************************************************************** 25 // 22 // 26 // G4ParticleChangeForLoss class implementatio << 27 // 23 // 28 // Author: Hisaya Kurashige, 23 March 1998 << 24 // $Id: G4ParticleChangeForLoss.cc,v 1.13 2004/06/15 08:17:38 vnivanch Exp $ 29 // Revision: Vladimir Ivantchenko, 16 January << 25 // GEANT4 tag $Name: geant4-07-00-cand-01 $ 30 // ------------------------------------------- << 26 // 31 << 27 // >> 28 // -------------------------------------------------------------- >> 29 // GEANT 4 class implementation file >> 30 // >> 31 // ------------------------------------------------------------ >> 32 // Implemented for the new scheme 23 Mar. 1998 H.Kurahige >> 33 // -------------------------------------------------------------- >> 34 // >> 35 // Modified: >> 36 // 16.01.04 V.Ivanchenko update for model variant of energy loss >> 37 // >> 38 // ------------------------------------------------------------ >> 39 // 32 #include "G4ParticleChangeForLoss.hh" 40 #include "G4ParticleChangeForLoss.hh" 33 #include "G4SystemOfUnits.hh" << 41 #include "G4Track.hh" >> 42 #include "G4Step.hh" >> 43 #include "G4DynamicParticle.hh" 34 #include "G4ExceptionSeverity.hh" 44 #include "G4ExceptionSeverity.hh" 35 45 36 // ------------------------------------------- << 46 G4ParticleChangeForLoss::G4ParticleChangeForLoss():G4VParticleChange() 37 G4ParticleChangeForLoss::G4ParticleChangeForLo << 38 { 47 { 39 // Disable flag that is enabled in G4VPartic << 48 theSteppingControlFlag = NormalCondition; 40 debugFlag = false; 49 debugFlag = false; 41 SetNumberOfSecondaries(1); << 50 #ifdef G4VERBOSE >> 51 if (verboseLevel>2) { >> 52 G4cout << "G4ParticleChangeForLoss::G4ParticleChangeForLoss() " << G4endl; >> 53 } >> 54 #endif 42 } 55 } 43 56 44 // ------------------------------------------- << 57 G4ParticleChangeForLoss::~G4ParticleChangeForLoss() 45 void G4ParticleChangeForLoss::DumpInfo() const << 46 { 58 { 47 // use base-class DumpInfo << 59 #ifdef G4VERBOSE 48 G4VParticleChange::DumpInfo(); << 60 if (verboseLevel>2) { >> 61 G4cout << "G4ParticleChangeForLoss::~G4ParticleChangeForLoss() " << G4endl; >> 62 } >> 63 #endif >> 64 } >> 65 >> 66 G4ParticleChangeForLoss::G4ParticleChangeForLoss( >> 67 const G4ParticleChangeForLoss &right): G4VParticleChange(right) >> 68 { >> 69 if (verboseLevel>1) { >> 70 G4cout << "G4ParticleChangeForLoss:: copy constructor is called " << G4endl; >> 71 } >> 72 currentTrack = right.currentTrack; >> 73 proposedKinEnergy = right.proposedKinEnergy; >> 74 currentCharge = right.currentCharge; >> 75 //theProposedWeight = right.theProposedWeight; >> 76 proposedMomentumDirection = right.proposedMomentumDirection; >> 77 } 49 78 50 G4long oldprc = G4cout.precision(8); << 79 // assignment operator 51 G4cout << " --------------------------- << 80 G4ParticleChangeForLoss & G4ParticleChangeForLoss::operator=( 52 G4cout << " G4ParticleChangeForLoss p << 81 const G4ParticleChangeForLoss &right) 53 G4cout << " Charge (eplus) : " << s << 82 { 54 << currentCharge / eplus << G4endl; << 83 if (verboseLevel>1) { 55 G4cout << " Kinetic Energy (MeV): " < << 84 G4cout << "G4ParticleChangeForLoss:: assignment operator is called " << G4endl; 56 << proposedKinEnergy / MeV << G4endl; << 85 } 57 G4cout << " Momentum Direct - x : " < << 86 if (this != &right) 58 << proposedMomentumDirection.x() << G << 87 { 59 G4cout << " Momentum Direct - y : " < << 88 theListOfSecondaries = right.theListOfSecondaries; 60 << proposedMomentumDirection.y() << G << 89 theSizeOftheListOfSecondaries = right.theSizeOftheListOfSecondaries; 61 G4cout << " Momentum Direct - z : " < << 90 theNumberOfSecondaries = right.theNumberOfSecondaries; 62 << proposedMomentumDirection.z() << G << 91 theStatusChange = right.theStatusChange; 63 G4cout.precision(oldprc); << 92 theLocalEnergyDeposit = right.theLocalEnergyDeposit; >> 93 theSteppingControlFlag = right.theSteppingControlFlag; >> 94 >> 95 currentTrack = right.currentTrack; >> 96 proposedKinEnergy = right.proposedKinEnergy; >> 97 currentCharge = right.currentCharge; >> 98 //theProposedWeight = right.theProposedWeight; >> 99 proposedMomentumDirection = right.proposedMomentumDirection; >> 100 } >> 101 return *this; 64 } 102 } 65 103 66 // ------------------------------------------- << 104 //---------------------------------------------------------------- >> 105 // methods for updating G4Step >> 106 // >> 107 67 G4Step* G4ParticleChangeForLoss::UpdateStepFor 108 G4Step* G4ParticleChangeForLoss::UpdateStepForAlongStep(G4Step* pStep) 68 { 109 { 69 const G4StepPoint* pPreStepPoint = pStep->Ge << 110 // A physics process always calculates the final state of the >> 111 // particle relative to the initial state at the beginning >> 112 // of the Step, i.e., based on information of G4Track (or >> 113 // equivalently the PreStepPoint). >> 114 // So, the differences (delta) between these two states have to be >> 115 // calculated and be accumulated in PostStepPoint. >> 116 >> 117 G4StepPoint* pPreStepPoint = pStep->GetPreStepPoint(); 70 G4StepPoint* pPostStepPoint = pStep->GetPost 118 G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint(); 71 119 72 // accumulate change of the kinetic energy << 120 // calculate new kinetic energy 73 G4double preKinEnergy = pPreStepPoint->GetKi << 121 G4double kinEnergy = pPostStepPoint->GetKineticEnergy() 74 G4double kinEnergy = << 122 + (proposedKinEnergy - pPreStepPoint->GetKineticEnergy()); 75 pPostStepPoint->GetKineticEnergy() + (prop << 123 76 << 124 // update kinetic energy and momentum direction 77 pPostStepPoint->SetCharge(currentCharge); << 125 if (kinEnergy < 0.0) { 78 << 126 theLocalEnergyDeposit += kinEnergy; 79 // calculate velocity << 127 kinEnergy = 0.0; 80 if(kinEnergy > 0.0) << 81 { << 82 pPostStepPoint->SetKineticEnergy(kinEnergy << 83 << 84 // assuming that mass>0, zero mass particl << 85 pPostStepPoint->SetVelocity(CLHEP::c_light << 86 } << 87 else << 88 { << 89 pPostStepPoint->SetKineticEnergy(0.0); << 90 pPostStepPoint->SetVelocity(0.0); << 91 } << 92 << 93 if(isParentWeightProposed) << 94 { << 95 pPostStepPoint->SetWeight(theParentWeight) << 96 } 128 } 97 129 98 pStep->AddTotalEnergyDeposit(theLocalEnergyD << 130 pPostStepPoint->SetKineticEnergy( kinEnergy ); 99 pStep->AddNonIonizingEnergyDeposit(theNonIon << 131 pPostStepPoint->SetCharge( currentCharge ); >> 132 >> 133 // update weight >> 134 // this feature is commented out, it should be overwritten in case >> 135 // if energy loss processes will use biasing >> 136 // G4double newWeight = theProposedWeight/(pPreStepPoint->GetWeight())*(pPostStepPoint->GetWeight()); >> 137 // pPostStepPoint->SetWeight( newWeight ); >> 138 >> 139 >> 140 // Not necessary to check now >> 141 //#ifdef G4VERBOSE >> 142 // if (debugFlag) CheckIt(*aTrack); >> 143 //#endif >> 144 >> 145 // Update the G4Step specific attributes >> 146 pStep->AddTotalEnergyDeposit( theLocalEnergyDeposit ); 100 return pStep; 147 return pStep; 101 } 148 } 102 149 103 // ------------------------------------------- << 104 G4Step* G4ParticleChangeForLoss::UpdateStepFor 150 G4Step* G4ParticleChangeForLoss::UpdateStepForPostStep(G4Step* pStep) 105 { 151 { >> 152 // A physics process always calculates the final state of the >> 153 // particle relative to the initial state at the beginning >> 154 // of the Step, i.e., based on information of G4Track (or >> 155 // equivalently the PreStepPoint). >> 156 // So, the differences (delta) between these two states have to be >> 157 // calculated and be accumulated in PostStepPoint. >> 158 106 G4StepPoint* pPostStepPoint = pStep->GetPost 159 G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint(); 107 160 108 pPostStepPoint->SetCharge(currentCharge); << 161 pPostStepPoint->SetKineticEnergy( proposedKinEnergy ); 109 pPostStepPoint->SetMomentumDirection(propose << 162 pPostStepPoint->SetCharge( currentCharge ); 110 if(proposedKinEnergy > 0.0) << 163 pPostStepPoint->SetMomentumDirection( proposedMomentumDirection ); 111 { << 164 112 pPostStepPoint->SetKineticEnergy(proposedK << 165 // update weight >> 166 // this feature is commented out, it should be overwritten in case >> 167 // if energy loss processes will use biasing >> 168 // pPostStepPoint->SetWeight( theProposedWeight ); >> 169 >> 170 // Not necessary to check now >> 171 //#ifdef G4VERBOSE >> 172 // if (debugFlag) CheckIt(*aTrack); >> 173 //#endif >> 174 >> 175 // Update the G4Step specific attributes >> 176 pStep->AddTotalEnergyDeposit( theLocalEnergyDeposit ); >> 177 return pStep; >> 178 } >> 179 >> 180 //---------------------------------------------------------------- >> 181 // methods for printing messages >> 182 // >> 183 >> 184 void G4ParticleChangeForLoss::DumpInfo() const >> 185 { >> 186 // use base-class DumpInfo >> 187 G4VParticleChange::DumpInfo(); >> 188 >> 189 G4cout.precision(3); >> 190 G4cout << " Charge (eplus) : " >> 191 << std::setw(20) << currentCharge/eplus >> 192 << G4endl; >> 193 G4cout << " Kinetic Energy (MeV): " >> 194 << std::setw(20) << proposedKinEnergy/MeV >> 195 << G4endl; >> 196 G4cout << " Momentum Direct - x : " >> 197 << std::setw(20) << proposedMomentumDirection.x() >> 198 << G4endl; >> 199 G4cout << " Momentum Direct - y : " >> 200 << std::setw(20) << proposedMomentumDirection.y() >> 201 << G4endl; >> 202 G4cout << " Momentum Direct - z : " >> 203 << std::setw(20) << proposedMomentumDirection.z() >> 204 << G4endl; >> 205 } >> 206 >> 207 G4bool G4ParticleChangeForLoss::CheckIt(const G4Track& aTrack) >> 208 { >> 209 G4bool itsOK = true; >> 210 G4bool exitWithError = false; >> 211 >> 212 G4double accuracy; 113 213 114 // assuming that mass>0, zero mass particl << 214 // Energy should not be lager than initial value 115 pPostStepPoint->SetVelocity(CLHEP::c_light << 215 accuracy = ( proposedKinEnergy - aTrack.GetKineticEnergy())/MeV; >> 216 if (accuracy > accuracyForWarning) { >> 217 #ifdef G4VERBOSE >> 218 G4cout << "G4ParticleChangeForLoss::CheckIt: "; >> 219 G4cout << "KinEnergy become larger than the initial value!" << G4endl; >> 220 G4cout << " Difference: " << accuracy << "[MeV] " <<G4endl; >> 221 #endif >> 222 itsOK = false; >> 223 if (accuracy > accuracyForException) exitWithError = true; 116 } 224 } 117 else << 225 118 { << 226 // dump out information of this particle change 119 pPostStepPoint->SetKineticEnergy(0.0); << 227 #ifdef G4VERBOSE 120 pPostStepPoint->SetVelocity(0.0); << 228 if (!itsOK) { >> 229 G4cout << "G4ParticleChangeForLoss::CheckIt " << G4endl; >> 230 DumpInfo(); 121 } 231 } 122 pPostStepPoint->SetPolarization(proposedPola << 232 #endif 123 233 124 if(isParentWeightProposed) << 234 // Exit with error 125 { << 235 if (exitWithError) { 126 pPostStepPoint->SetWeight(theParentWeight) << 236 G4Exception("G4ParticleChangeForLoss::CheckIt", >> 237 "400", >> 238 EventMustBeAborted, >> 239 "energy was illegal"); 127 } 240 } 128 241 129 pStep->AddTotalEnergyDeposit(theLocalEnergyD << 242 //correction 130 pStep->AddNonIonizingEnergyDeposit(theNonIon << 243 if (!itsOK) { 131 return pStep; << 244 proposedKinEnergy = aTrack.GetKineticEnergy(); >> 245 } >> 246 >> 247 itsOK = (itsOK) && G4VParticleChange::CheckIt(aTrack); >> 248 return itsOK; 132 } 249 } 133 250