--------------------------------------------------------------------
| 135 | |
| 136 | // -------------------------------------------------------------------- |
| 137 | G4Step* G4ParticleChange::UpdateStepForAlongStep(G4Step* pStep) |
| 138 | { |
| 139 | // A physics process always calculates the final state of the |
| 140 | // particle relative to the initial state at the beginning |
| 141 | // of the Step, i.e., based on information of G4Track (or |
| 142 | // equivalently the PreStepPoint). |
| 143 | // So, the differences (delta) between these two states have to be |
| 144 | // calculated and be accumulated in PostStepPoint |
| 145 | |
| 146 | // Take note that the return type of GetMomentumDirectionChange() is a |
| 147 | // pointer to G4ParticleMometum. Also it is a normalized |
| 148 | // momentum vector |
| 149 | |
| 150 | const G4StepPoint* pPreStepPoint = pStep->GetPreStepPoint(); |
| 151 | G4StepPoint* pPostStepPoint = pStep->GetPostStepPoint(); |
| 152 | |
| 153 | // set Mass/Charge/MagneticMoment |
| 154 | pPostStepPoint->SetMass(theMassChange); |
| 155 | pPostStepPoint->SetCharge(theChargeChange); |
| 156 | pPostStepPoint->SetMagneticMoment(theMagneticMomentChange); |
| 157 | |
| 158 | // calculate new kinetic energy |
| 159 | G4double preEnergy = pPreStepPoint->GetKineticEnergy(); |
| 160 | G4double energy = |
| 161 | pPostStepPoint->GetKineticEnergy() + (theEnergyChange - preEnergy); |
| 162 | |
| 163 | // update kinetic energy and momentum direction |
| 164 | if(energy > 0.0) |
| 165 | { |
| 166 | // calculate new momentum |
| 167 | G4ThreeVector pMomentum = pPostStepPoint->GetMomentum() |
| 168 | + (CalcMomentum(theEnergyChange, theMomentumDirectionChange, |
| 169 | theMassChange) - pPreStepPoint->GetMomentum()); |
| 170 | G4double tMomentum2 = pMomentum.mag2(); |
| 171 | G4ThreeVector direction(1.0, 0.0, 0.0); |
| 172 | if(tMomentum2 > 0.) |
| 173 | { |
| 174 | direction = pMomentum / std::sqrt(tMomentum2); |
| 175 | } |
| 176 | pPostStepPoint->SetMomentumDirection(direction); |
| 177 | pPostStepPoint->SetKineticEnergy(energy); |
| 178 | |
| 179 | // if velocity is not set it should be recomputed |
| 180 | // |
| 181 | if(!isVelocityChanged) |
| 182 | { |
| 183 | if (theMassChange > 0.0) |
| 184 | { |
| 185 | theVelocityChange = CLHEP::c_light * |
| 186 | std::sqrt(energy*(energy + 2*theMassChange))/(energy + theMassChange); |
| 187 | } |
| 188 | else |
| 189 | { |
| 190 | // zero mass particle |
| 191 | theVelocityChange = CLHEP::c_light; |
| 192 | // optical photon case |
| 193 | if(theCurrentTrack->GetParticleDefinition()->GetPDGEncoding() == -22) |
| 194 | { |
no test coverage detected