#include <G4VeLowEnergyLoss.hh>
Inheritance diagram for G4VeLowEnergyLoss:
Definition at line 77 of file G4VeLowEnergyLoss.hh.
G4VeLowEnergyLoss::G4VeLowEnergyLoss | ( | const G4String & | , | |
G4ProcessType | aType = fNotDefined | |||
) |
Definition at line 87 of file G4VeLowEnergyLoss.cc.
References G4VeLowEnergyLoss().
Referenced by G4VeLowEnergyLoss().
00088 : G4VContinuousDiscreteProcess(aName, aType), 00089 lastMaterial(0), 00090 imat(-1), 00091 f1Fluct(0),f2Fluct(0),e1Fluct(0),e2Fluct(0),rateFluct(0), 00092 ipotFluct(0),e1LogFluct(-1),e2LogFluct(-1),ipotLogFluct(-1), 00093 nmaxCont1(4), 00094 nmaxCont2(16) 00095 {;}
G4VeLowEnergyLoss::G4VeLowEnergyLoss | ( | G4VeLowEnergyLoss & | ) |
Definition at line 105 of file G4VeLowEnergyLoss.cc.
References G4VeLowEnergyLoss().
00106 : G4VContinuousDiscreteProcess(right), 00107 lastMaterial(0), 00108 imat(-1), 00109 f1Fluct(0),f2Fluct(0),e1Fluct(0),e2Fluct(0),rateFluct(0), 00110 ipotFluct(0),e1LogFluct(-1),e2LogFluct(-1),ipotLogFluct(-1), 00111 nmaxCont1(4), 00112 nmaxCont2(16) 00113 {;}
G4VeLowEnergyLoss::~G4VeLowEnergyLoss | ( | ) | [virtual] |
virtual G4VParticleChange* G4VeLowEnergyLoss::AlongStepDoIt | ( | const G4Track & | track, | |
const G4Step & | Step | |||
) | [pure virtual] |
Reimplemented from G4VContinuousDiscreteProcess.
G4PhysicsTable * G4VeLowEnergyLoss::BuildInverseRangeTable | ( | G4PhysicsTable * | theRangeTable, | |
G4PhysicsTable * | theRangeCoeffATable, | |||
G4PhysicsTable * | theRangeCoeffBTable, | |||
G4PhysicsTable * | theRangeCoeffCTable, | |||
G4PhysicsTable * | theInverseRangeTable, | |||
G4double | Tmin, | |||
G4double | Tmax, | |||
G4int | nbin | |||
) | [static, protected] |
Definition at line 524 of file G4VeLowEnergyLoss.cc.
References G4PhysicsTable::clearAndDestroy(), G4PhysicsVector::GetLowEdgeEnergy(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), G4PhysicsTable::insert(), and G4PhysicsVector::PutValue().
00533 { 00534 G4bool b; 00535 00536 G4int numOfCouples = G4ProductionCutsTable::GetProductionCutsTable()->GetTableSize(); 00537 00538 if(theInverseRangeTable) 00539 { theInverseRangeTable->clearAndDestroy(); 00540 delete theInverseRangeTable; } 00541 theInverseRangeTable = new G4PhysicsTable(numOfCouples); 00542 00543 // loop for materials 00544 for (G4int i=0; i<numOfCouples; i++) 00545 { 00546 00547 G4PhysicsVector* pv = (*theRangeTable)[i]; 00548 size_t nbins = pv->GetVectorLength(); 00549 G4double elow = pv->GetLowEdgeEnergy(0); 00550 G4double ehigh = pv->GetLowEdgeEnergy(nbins-1); 00551 G4double rlow = pv->GetValue(elow, b); 00552 G4double rhigh = pv->GetValue(ehigh, b); 00553 00554 //rhigh *= std::exp(std::log(rhigh/rlow)/((G4double)(nbins-1))); 00555 00556 G4PhysicsLogVector* v = new G4PhysicsLogVector(rlow, rhigh, nbins-1); 00557 00558 v->PutValue(0,elow); 00559 G4double energy1 = elow; 00560 G4double range1 = rlow; 00561 G4double energy2 = elow; 00562 G4double range2 = rlow; 00563 size_t ilow = 0; 00564 size_t ihigh; 00565 00566 for (size_t j=1; j<nbins; j++) { 00567 00568 G4double range = v->GetLowEdgeEnergy(j); 00569 00570 for (ihigh=ilow+1; ihigh<nbins; ihigh++) { 00571 energy2 = pv->GetLowEdgeEnergy(ihigh); 00572 range2 = pv->GetValue(energy2, b); 00573 if(range2 >= range || ihigh == nbins-1) { 00574 ilow = ihigh - 1; 00575 energy1 = pv->GetLowEdgeEnergy(ilow); 00576 range1 = pv->GetValue(energy1, b); 00577 break; 00578 } 00579 } 00580 00581 G4double e = std::log(energy1) + std::log(energy2/energy1)*std::log(range/range1)/std::log(range2/range1); 00582 00583 v->PutValue(j,std::exp(e)); 00584 } 00585 theInverseRangeTable->insert(v); 00586 00587 } 00588 return theInverseRangeTable ; 00589 }
G4PhysicsTable * G4VeLowEnergyLoss::BuildLabTimeTable | ( | G4PhysicsTable * | theDEDXTable, | |
G4PhysicsTable * | theLabTimeTable, | |||
G4double | Tmin, | |||
G4double | Tmax, | |||
G4int | nbin | |||
) | [static, protected] |
Definition at line 268 of file G4VeLowEnergyLoss.cc.
References G4PhysicsTable::clearAndDestroy(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), and G4PhysicsTable::insert().
00273 { 00274 00275 G4int numOfCouples = G4ProductionCutsTable::GetProductionCutsTable()->GetTableSize(); 00276 00277 if(theLabTimeTable) 00278 { theLabTimeTable->clearAndDestroy(); 00279 delete theLabTimeTable; } 00280 theLabTimeTable = new G4PhysicsTable(numOfCouples); 00281 00282 00283 for (G4int J=0; J<numOfCouples; J++) 00284 { 00285 G4PhysicsLogVector* aVector; 00286 00287 aVector = new G4PhysicsLogVector(lowestKineticEnergy, 00288 highestKineticEnergy,TotBin); 00289 00290 BuildLabTimeVector(theDEDXTable, 00291 lowestKineticEnergy,highestKineticEnergy,TotBin,J,aVector); 00292 theLabTimeTable->insert(aVector); 00293 00294 00295 } 00296 return theLabTimeTable ; 00297 }
G4PhysicsTable * G4VeLowEnergyLoss::BuildProperTimeTable | ( | G4PhysicsTable * | theDEDXTable, | |
G4PhysicsTable * | ProperTimeTable, | |||
G4double | Tmin, | |||
G4double | Tmax, | |||
G4int | nbin | |||
) | [static, protected] |
Definition at line 301 of file G4VeLowEnergyLoss.cc.
References G4PhysicsTable::clearAndDestroy(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), and G4PhysicsTable::insert().
00306 { 00307 00308 G4int numOfCouples = G4ProductionCutsTable::GetProductionCutsTable()->GetTableSize(); 00309 00310 if(theProperTimeTable) 00311 { theProperTimeTable->clearAndDestroy(); 00312 delete theProperTimeTable; } 00313 theProperTimeTable = new G4PhysicsTable(numOfCouples); 00314 00315 00316 for (G4int J=0; J<numOfCouples; J++) 00317 { 00318 G4PhysicsLogVector* aVector; 00319 00320 aVector = new G4PhysicsLogVector(lowestKineticEnergy, 00321 highestKineticEnergy,TotBin); 00322 00323 BuildProperTimeVector(theDEDXTable, 00324 lowestKineticEnergy,highestKineticEnergy,TotBin,J,aVector); 00325 theProperTimeTable->insert(aVector); 00326 00327 00328 } 00329 return theProperTimeTable ; 00330 }
G4PhysicsTable * G4VeLowEnergyLoss::BuildRangeCoeffATable | ( | G4PhysicsTable * | theRangeTable, | |
G4PhysicsTable * | theCoeffATable, | |||
G4double | Tmin, | |||
G4double | Tmax, | |||
G4int | nbin | |||
) | [static, protected] |
Definition at line 649 of file G4VeLowEnergyLoss.cc.
References G4PhysicsTable::clearAndDestroy(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), G4PhysicsVector::GetValue(), G4PhysicsTable::insert(), and G4PhysicsVector::PutValue().
00655 { 00656 00657 G4int numOfCouples = G4ProductionCutsTable::GetProductionCutsTable()->GetTableSize(); 00658 00659 if(theRangeCoeffATable) 00660 { theRangeCoeffATable->clearAndDestroy(); 00661 delete theRangeCoeffATable; } 00662 theRangeCoeffATable = new G4PhysicsTable(numOfCouples); 00663 00664 G4double RTable = std::exp(std::log(highestKineticEnergy/lowestKineticEnergy)/TotBin) ; 00665 G4double R2 = RTable*RTable ; 00666 G4double R1 = RTable+1.; 00667 G4double w = R1*(RTable-1.)*(RTable-1.); 00668 G4double w1 = RTable/w , w2 = -RTable*R1/w , w3 = R2/w ; 00669 G4double Ti , Tim , Tip , Ri , Rim , Rip , Value ; 00670 G4bool isOut; 00671 00672 // loop for materials 00673 for (G4int J=0; J<numOfCouples; J++) 00674 { 00675 G4int binmax=TotBin ; 00676 G4PhysicsLinearVector* aVector = 00677 new G4PhysicsLinearVector(0.,binmax, TotBin); 00678 Ti = lowestKineticEnergy ; 00679 G4PhysicsVector* rangeVector= (*theRangeTable)[J]; 00680 00681 for ( G4int i=0; i<=TotBin; i++) 00682 { 00683 Ri = rangeVector->GetValue(Ti,isOut) ; 00684 if ( i==0 ) 00685 Rim = 0. ; 00686 else 00687 { 00688 Tim = Ti/RTable ; 00689 Rim = rangeVector->GetValue(Tim,isOut); 00690 } 00691 if ( i==TotBin) 00692 Rip = Ri ; 00693 else 00694 { 00695 Tip = Ti*RTable ; 00696 Rip = rangeVector->GetValue(Tip,isOut); 00697 } 00698 Value = (w1*Rip + w2*Ri + w3*Rim)/(Ti*Ti) ; 00699 00700 aVector->PutValue(i,Value); 00701 Ti = RTable*Ti ; 00702 } 00703 00704 theRangeCoeffATable->insert(aVector); 00705 } 00706 return theRangeCoeffATable ; 00707 }
G4PhysicsTable * G4VeLowEnergyLoss::BuildRangeCoeffBTable | ( | G4PhysicsTable * | theRangeTable, | |
G4PhysicsTable * | theCoeffBTable, | |||
G4double | Tmin, | |||
G4double | Tmax, | |||
G4int | nbin | |||
) | [static, protected] |
Definition at line 711 of file G4VeLowEnergyLoss.cc.
References G4PhysicsTable::clearAndDestroy(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), G4PhysicsVector::GetValue(), G4PhysicsTable::insert(), and G4PhysicsVector::PutValue().
00717 { 00718 00719 G4int numOfCouples = G4ProductionCutsTable::GetProductionCutsTable()->GetTableSize(); 00720 00721 if(theRangeCoeffBTable) 00722 { theRangeCoeffBTable->clearAndDestroy(); 00723 delete theRangeCoeffBTable; } 00724 theRangeCoeffBTable = new G4PhysicsTable(numOfCouples); 00725 00726 G4double RTable = std::exp(std::log(highestKineticEnergy/lowestKineticEnergy)/TotBin) ; 00727 G4double R2 = RTable*RTable ; 00728 G4double R1 = RTable+1.; 00729 G4double w = R1*(RTable-1.)*(RTable-1.); 00730 G4double w1 = -R1/w , w2 = R1*(R2+1.)/w , w3 = -R2*R1/w ; 00731 G4double Ti , Tim , Tip , Ri , Rim , Rip , Value ; 00732 G4bool isOut; 00733 00734 // loop for materials 00735 for (G4int J=0; J<numOfCouples; J++) 00736 { 00737 G4int binmax=TotBin ; 00738 G4PhysicsLinearVector* aVector = 00739 new G4PhysicsLinearVector(0.,binmax, TotBin); 00740 Ti = lowestKineticEnergy ; 00741 G4PhysicsVector* rangeVector= (*theRangeTable)[J]; 00742 00743 for ( G4int i=0; i<=TotBin; i++) 00744 { 00745 Ri = rangeVector->GetValue(Ti,isOut) ; 00746 if ( i==0 ) 00747 Rim = 0. ; 00748 else 00749 { 00750 Tim = Ti/RTable ; 00751 Rim = rangeVector->GetValue(Tim,isOut); 00752 } 00753 if ( i==TotBin) 00754 Rip = Ri ; 00755 else 00756 { 00757 Tip = Ti*RTable ; 00758 Rip = rangeVector->GetValue(Tip,isOut); 00759 } 00760 Value = (w1*Rip + w2*Ri + w3*Rim)/Ti; 00761 00762 aVector->PutValue(i,Value); 00763 Ti = RTable*Ti ; 00764 } 00765 theRangeCoeffBTable->insert(aVector); 00766 } 00767 return theRangeCoeffBTable ; 00768 }
G4PhysicsTable * G4VeLowEnergyLoss::BuildRangeCoeffCTable | ( | G4PhysicsTable * | theRangeTable, | |
G4PhysicsTable * | theCoeffCTable, | |||
G4double | Tmin, | |||
G4double | Tmax, | |||
G4int | nbin | |||
) | [static, protected] |
Definition at line 772 of file G4VeLowEnergyLoss.cc.
References G4PhysicsTable::clearAndDestroy(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), G4PhysicsVector::GetValue(), G4PhysicsTable::insert(), and G4PhysicsVector::PutValue().
00778 { 00779 00780 G4int numOfCouples = G4ProductionCutsTable::GetProductionCutsTable()->GetTableSize(); 00781 00782 if(theRangeCoeffCTable) 00783 { theRangeCoeffCTable->clearAndDestroy(); 00784 delete theRangeCoeffCTable; } 00785 theRangeCoeffCTable = new G4PhysicsTable(numOfCouples); 00786 00787 G4double RTable = std::exp(std::log(highestKineticEnergy/lowestKineticEnergy)/TotBin) ; 00788 G4double R2 = RTable*RTable ; 00789 G4double R1 = RTable+1.; 00790 G4double w = R1*(RTable-1.)*(RTable-1.); 00791 G4double w1 = 1./w , w2 = -RTable*R1/w , w3 = RTable*R2/w ; 00792 G4double Ti , Tim , Tip , Ri , Rim , Rip , Value ; 00793 G4bool isOut; 00794 00795 // loop for materials 00796 for (G4int J=0; J<numOfCouples; J++) 00797 { 00798 G4int binmax=TotBin ; 00799 G4PhysicsLinearVector* aVector = 00800 new G4PhysicsLinearVector(0.,binmax, TotBin); 00801 Ti = lowestKineticEnergy ; 00802 G4PhysicsVector* rangeVector= (*theRangeTable)[J]; 00803 00804 for ( G4int i=0; i<=TotBin; i++) 00805 { 00806 Ri = rangeVector->GetValue(Ti,isOut) ; 00807 if ( i==0 ) 00808 Rim = 0. ; 00809 else 00810 { 00811 Tim = Ti/RTable ; 00812 Rim = rangeVector->GetValue(Tim,isOut); 00813 } 00814 if ( i==TotBin) 00815 Rip = Ri ; 00816 else 00817 { 00818 Tip = Ti*RTable ; 00819 Rip = rangeVector->GetValue(Tip,isOut); 00820 } 00821 Value = w1*Rip + w2*Ri + w3*Rim ; 00822 00823 aVector->PutValue(i,Value); 00824 Ti = RTable*Ti ; 00825 } 00826 theRangeCoeffCTable->insert(aVector); 00827 } 00828 return theRangeCoeffCTable ; 00829 }
G4PhysicsTable * G4VeLowEnergyLoss::BuildRangeTable | ( | G4PhysicsTable * | theDEDXTable, | |
G4PhysicsTable * | theRangeTable, | |||
G4double | Tmin, | |||
G4double | Tmax, | |||
G4int | nbin | |||
) | [static, protected] |
Definition at line 134 of file G4VeLowEnergyLoss.cc.
References G4PhysicsTable::clearAndDestroy(), G4PhysicsTable::insert(), and G4PhysicsTable::length().
00140 { 00141 00142 G4int numOfCouples = theDEDXTable->length(); 00143 00144 if(theRangeTable) 00145 { theRangeTable->clearAndDestroy(); 00146 delete theRangeTable; } 00147 theRangeTable = new G4PhysicsTable(numOfCouples); 00148 00149 // loop for materials 00150 00151 for (G4int J=0; J<numOfCouples; J++) 00152 { 00153 G4PhysicsLogVector* aVector; 00154 aVector = new G4PhysicsLogVector(lowestKineticEnergy, 00155 highestKineticEnergy,TotBin); 00156 BuildRangeVector(theDEDXTable,lowestKineticEnergy,highestKineticEnergy, 00157 TotBin,J,aVector); 00158 theRangeTable->insert(aVector); 00159 } 00160 return theRangeTable ; 00161 }
virtual G4double G4VeLowEnergyLoss::GetContinuousStepLimit | ( | const G4Track & | track, | |
G4double | previousStepSize, | |||
G4double | currentMinimumStep, | |||
G4double & | currentSafety | |||
) | [pure virtual] |
Implements G4VContinuousDiscreteProcess.
G4double G4VeLowEnergyLoss::GetLossWithFluct | ( | const G4DynamicParticle * | aParticle, | |
const G4MaterialCutsCouple * | couple, | |||
G4double | MeanLoss, | |||
G4double | step | |||
) | [protected] |
Definition at line 833 of file G4VeLowEnergyLoss.cc.
References e1Fluct, e1LogFluct, e2Fluct, e2LogFluct, f1Fluct, f2Fluct, G4Poisson(), G4UniformRand, G4Material::GetElectronDensity(), G4IonisParamMat::GetEnergy0fluct(), G4IonisParamMat::GetEnergy1fluct(), G4IonisParamMat::GetEnergy2fluct(), G4IonisParamMat::GetF1fluct(), G4IonisParamMat::GetF2fluct(), G4MaterialCutsCouple::GetIndex(), G4Material::GetIonisation(), G4DynamicParticle::GetKineticEnergy(), G4IonisParamMat::GetLogEnergy1fluct(), G4IonisParamMat::GetLogEnergy2fluct(), G4IonisParamMat::GetLogMeanExcEnergy(), G4MaterialCutsCouple::GetMaterial(), G4IonisParamMat::GetMeanExcitationEnergy(), G4ProductionCutsTable::GetProductionCutsTable(), G4IonisParamMat::GetRateionexcfluct(), imat, ipotFluct, ipotLogFluct, lastMaterial, nmaxCont1, nmaxCont2, ParticleMass, and rateFluct.
00839 { 00840 static const G4double minLoss = 1.*eV ; 00841 static const G4double probLim = 0.01 ; 00842 static const G4double sumaLim = -std::log(probLim) ; 00843 static const G4double alim=10.; 00844 static const G4double kappa = 10. ; 00845 static const G4double factor = twopi_mc2_rcl2 ; 00846 const G4Material* aMaterial = couple->GetMaterial(); 00847 00848 // check if the material has changed ( cache mechanism) 00849 00850 if (aMaterial != lastMaterial) 00851 { 00852 lastMaterial = aMaterial; 00853 imat = couple->GetIndex(); 00854 f1Fluct = aMaterial->GetIonisation()->GetF1fluct(); 00855 f2Fluct = aMaterial->GetIonisation()->GetF2fluct(); 00856 e1Fluct = aMaterial->GetIonisation()->GetEnergy1fluct(); 00857 e2Fluct = aMaterial->GetIonisation()->GetEnergy2fluct(); 00858 e1LogFluct = aMaterial->GetIonisation()->GetLogEnergy1fluct(); 00859 e2LogFluct = aMaterial->GetIonisation()->GetLogEnergy2fluct(); 00860 rateFluct = aMaterial->GetIonisation()->GetRateionexcfluct(); 00861 ipotFluct = aMaterial->GetIonisation()->GetMeanExcitationEnergy(); 00862 ipotLogFluct = aMaterial->GetIonisation()->GetLogMeanExcEnergy(); 00863 } 00864 G4double threshold,w1,w2,C, 00865 beta2,suma,e0,loss,lossc,w; 00866 G4double a1,a2,a3; 00867 G4int p1,p2,p3; 00868 G4int nb; 00869 G4double Corrfac, na,alfa,rfac,namean,sa,alfa1,ea,sea; 00870 // G4double dp1; 00871 G4double dp3; 00872 G4double siga ; 00873 00874 // shortcut for very very small loss 00875 if(MeanLoss < minLoss) return MeanLoss ; 00876 00877 // get particle data 00878 G4double Tkin = aParticle->GetKineticEnergy(); 00879 00880 // G4cout << "MGP -- Fluc Tkin " << Tkin/keV << " keV " << " MeanLoss = " << MeanLoss/keV << G4endl; 00881 00882 threshold = (*((G4ProductionCutsTable::GetProductionCutsTable()) 00883 ->GetEnergyCutsVector(1)))[imat]; 00884 G4double rmass = electron_mass_c2/ParticleMass; 00885 G4double tau = Tkin/ParticleMass, tau1 = tau+1., tau2 = tau*(tau+2.); 00886 G4double Tm = 2.*electron_mass_c2*tau2/(1.+2.*tau1*rmass+rmass*rmass); 00887 00888 // G4cout << "MGP Particle mass " << ParticleMass/MeV << " Tm " << Tm << G4endl; 00889 00890 if(Tm > threshold) Tm = threshold; 00891 beta2 = tau2/(tau1*tau1); 00892 00893 // Gaussian fluctuation ? 00894 if(MeanLoss >= kappa*Tm || MeanLoss <= kappa*ipotFluct) 00895 { 00896 G4double electronDensity = aMaterial->GetElectronDensity() ; 00897 siga = std::sqrt(Tm*(1.0-0.5*beta2)*step* 00898 factor*electronDensity/beta2) ; 00899 do { 00900 loss = G4RandGauss::shoot(MeanLoss,siga) ; 00901 } while (loss < 0. || loss > 2.0*MeanLoss); 00902 return loss ; 00903 } 00904 00905 w1 = Tm/ipotFluct; 00906 w2 = std::log(2.*electron_mass_c2*tau2); 00907 00908 C = MeanLoss*(1.-rateFluct)/(w2-ipotLogFluct-beta2); 00909 00910 a1 = C*f1Fluct*(w2-e1LogFluct-beta2)/e1Fluct; 00911 a2 = C*f2Fluct*(w2-e2LogFluct-beta2)/e2Fluct; 00912 a3 = rateFluct*MeanLoss*(Tm-ipotFluct)/(ipotFluct*Tm*std::log(w1)); 00913 00914 suma = a1+a2+a3; 00915 00916 loss = 0. ; 00917 00918 if(suma < sumaLim) // very small Step 00919 { 00920 e0 = aMaterial->GetIonisation()->GetEnergy0fluct(); 00921 // G4cout << "MGP e0 = " << e0/keV << G4endl; 00922 00923 if(Tm == ipotFluct) 00924 { 00925 a3 = MeanLoss/e0; 00926 00927 if(a3>alim) 00928 { 00929 siga=std::sqrt(a3) ; 00930 p3 = std::max(0,G4int(G4RandGauss::shoot(a3,siga)+0.5)); 00931 } 00932 else p3 = G4Poisson(a3); 00933 00934 loss = p3*e0 ; 00935 00936 if(p3 > 0) loss += (1.-2.*G4UniformRand())*e0 ; 00937 // G4cout << "MGP very small step " << loss/keV << G4endl; 00938 } 00939 else 00940 { 00941 // G4cout << "MGP old Tm = " << Tm << " " << ipotFluct << " " << e0 << G4endl; 00942 Tm = Tm-ipotFluct+e0 ; 00943 00944 // MGP ---- workaround to avoid log argument<0, TO BE CHECKED 00945 if (Tm <= 0.) 00946 { 00947 loss = MeanLoss; 00948 p3 = 0; 00949 // G4cout << "MGP correction loss = MeanLoss " << loss/keV << G4endl; 00950 } 00951 else 00952 { 00953 a3 = MeanLoss*(Tm-e0)/(Tm*e0*std::log(Tm/e0)); 00954 00955 // G4cout << "MGP new Tm = " << Tm << " " << ipotFluct << " " << e0 << " a3= " << a3 << G4endl; 00956 00957 if(a3>alim) 00958 { 00959 siga=std::sqrt(a3) ; 00960 p3 = std::max(0,G4int(G4RandGauss::shoot(a3,siga)+0.5)); 00961 } 00962 else 00963 p3 = G4Poisson(a3); 00964 //G4cout << "MGP p3 " << p3 << G4endl; 00965 00966 } 00967 00968 if(p3 > 0) 00969 { 00970 w = (Tm-e0)/Tm ; 00971 if(p3 > nmaxCont2) 00972 { 00973 // G4cout << "MGP dp3 " << dp3 << " p3 " << p3 << " " << nmaxCont2 << G4endl; 00974 dp3 = G4double(p3) ; 00975 Corrfac = dp3/G4double(nmaxCont2) ; 00976 p3 = nmaxCont2 ; 00977 } 00978 else 00979 Corrfac = 1. ; 00980 00981 for(G4int i=0; i<p3; i++) loss += 1./(1.-w*G4UniformRand()) ; 00982 loss *= e0*Corrfac ; 00983 // G4cout << "MGP Corrfac = " << Corrfac << " e0 = " << e0/keV << " loss = " << loss/keV << G4endl; 00984 } 00985 } 00986 } 00987 00988 else // not so small Step 00989 { 00990 // excitation type 1 00991 if(a1>alim) 00992 { 00993 siga=std::sqrt(a1) ; 00994 p1 = std::max(0,int(G4RandGauss::shoot(a1,siga)+0.5)); 00995 } 00996 else 00997 p1 = G4Poisson(a1); 00998 00999 // excitation type 2 01000 if(a2>alim) 01001 { 01002 siga=std::sqrt(a2) ; 01003 p2 = std::max(0,int(G4RandGauss::shoot(a2,siga)+0.5)); 01004 } 01005 else 01006 p2 = G4Poisson(a2); 01007 01008 loss = p1*e1Fluct+p2*e2Fluct; 01009 01010 // smearing to avoid unphysical peaks 01011 if(p2 > 0) 01012 loss += (1.-2.*G4UniformRand())*e2Fluct; 01013 else if (loss>0.) 01014 loss += (1.-2.*G4UniformRand())*e1Fluct; 01015 01016 // ionisation ....................................... 01017 if(a3 > 0.) 01018 { 01019 if(a3>alim) 01020 { 01021 siga=std::sqrt(a3) ; 01022 p3 = std::max(0,int(G4RandGauss::shoot(a3,siga)+0.5)); 01023 } 01024 else 01025 p3 = G4Poisson(a3); 01026 01027 lossc = 0.; 01028 if(p3 > 0) 01029 { 01030 na = 0.; 01031 alfa = 1.; 01032 if (p3 > nmaxCont2) 01033 { 01034 dp3 = G4double(p3); 01035 rfac = dp3/(G4double(nmaxCont2)+dp3); 01036 namean = G4double(p3)*rfac; 01037 sa = G4double(nmaxCont1)*rfac; 01038 na = G4RandGauss::shoot(namean,sa); 01039 if (na > 0.) 01040 { 01041 alfa = w1*G4double(nmaxCont2+p3)/(w1*G4double(nmaxCont2)+G4double(p3)); 01042 alfa1 = alfa*std::log(alfa)/(alfa-1.); 01043 ea = na*ipotFluct*alfa1; 01044 sea = ipotFluct*std::sqrt(na*(alfa-alfa1*alfa1)); 01045 lossc += G4RandGauss::shoot(ea,sea); 01046 } 01047 } 01048 01049 nb = G4int(G4double(p3)-na); 01050 if (nb > 0) 01051 { 01052 w2 = alfa*ipotFluct; 01053 w = (Tm-w2)/Tm; 01054 for (G4int k=0; k<nb; k++) lossc += w2/(1.-w*G4UniformRand()); 01055 } 01056 } 01057 01058 loss += lossc; 01059 } 01060 } 01061 01062 return loss ; 01063 }
virtual G4double G4VeLowEnergyLoss::GetMeanFreePath | ( | const G4Track & | track, | |
G4double | previousStepSize, | |||
G4ForceCondition * | condition | |||
) | [pure virtual] |
Implements G4VContinuousDiscreteProcess.
virtual G4VParticleChange* G4VeLowEnergyLoss::PostStepDoIt | ( | const G4Track & | track, | |
const G4Step & | Step | |||
) | [pure virtual] |
Reimplemented from G4VContinuousDiscreteProcess.
void G4VeLowEnergyLoss::SetEnlossFluc | ( | G4bool | value | ) | [static] |
Definition at line 120 of file G4VeLowEnergyLoss.cc.
References EnlossFlucFlag.
00121 { 00122 EnlossFlucFlag = value; 00123 }
void G4VeLowEnergyLoss::SetRndmStep | ( | G4bool | value | ) | [static] |
Definition at line 115 of file G4VeLowEnergyLoss.cc.
References rndmStepFlag.
00116 { 00117 rndmStepFlag = value; 00118 }
Definition at line 125 of file G4VeLowEnergyLoss.cc.
References c1lim, c2lim, c3lim, dRoverRange, and finalRange.
00126 { 00127 dRoverRange = c1; 00128 finalRange = c2; 00129 c1lim=dRoverRange; 00130 c2lim=2.*(1-dRoverRange)*finalRange; 00131 c3lim=-(1.-dRoverRange)*finalRange*finalRange; 00132 }
G4double G4VeLowEnergyLoss::c1lim = dRoverRange [static, protected] |
G4double G4VeLowEnergyLoss::c2lim = 2.*(1.-dRoverRange)*finalRange [static, protected] |
G4double G4VeLowEnergyLoss::c3lim = -(1.-dRoverRange)*finalRange*finalRange [static, protected] |
G4double G4VeLowEnergyLoss::dRoverRange = 20*perCent [static, protected] |
G4double G4VeLowEnergyLoss::e1Fluct [protected] |
G4double G4VeLowEnergyLoss::e1LogFluct [protected] |
G4double G4VeLowEnergyLoss::e2Fluct [protected] |
G4double G4VeLowEnergyLoss::e2LogFluct [protected] |
G4bool G4VeLowEnergyLoss::EnlossFlucFlag = true [static, protected] |
G4double G4VeLowEnergyLoss::f1Fluct [protected] |
G4double G4VeLowEnergyLoss::f2Fluct [protected] |
G4double G4VeLowEnergyLoss::finalRange = 200*micrometer [static, protected] |
G4int G4VeLowEnergyLoss::imat [protected] |
G4double G4VeLowEnergyLoss::ipotFluct [protected] |
G4double G4VeLowEnergyLoss::ipotLogFluct [protected] |
const G4Material* G4VeLowEnergyLoss::lastMaterial [protected] |
G4double G4VeLowEnergyLoss::ltauhigh [static, protected] |
Definition at line 231 of file G4VeLowEnergyLoss.hh.
G4double G4VeLowEnergyLoss::ltaulow [static, protected] |
Definition at line 231 of file G4VeLowEnergyLoss.hh.
const G4int G4VeLowEnergyLoss::nmaxCont1 [protected] |
const G4int G4VeLowEnergyLoss::nmaxCont2 [protected] |
G4double G4VeLowEnergyLoss::ParticleMass [static, protected] |
G4double G4VeLowEnergyLoss::rateFluct [protected] |
G4bool G4VeLowEnergyLoss::rndmStepFlag = false [static, protected] |
G4double G4VeLowEnergyLoss::tauhigh [static, protected] |
Definition at line 231 of file G4VeLowEnergyLoss.hh.
G4double G4VeLowEnergyLoss::taulow [static, protected] |
Definition at line 231 of file G4VeLowEnergyLoss.hh.