File indexing completed on 2026-03-30 07:50:48
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038
0039
0040
0041
0042
0043 #include "Run.hh"
0044
0045 #include "DetectorConstruction.hh"
0046 #include "PrimaryGeneratorAction.hh"
0047
0048 #include "G4AnalysisManager.hh"
0049 #include "G4EmCalculator.hh"
0050 #include "G4H2O.hh"
0051 #include "G4Molecule.hh"
0052 #include "G4MoleculeCounter.hh"
0053 #include "G4MoleculeGun.hh"
0054 #include "G4MoleculeTable.hh"
0055 #include "G4SystemOfUnits.hh"
0056 #include "G4UnitsTable.hh"
0057 #include "G4ios.hh"
0058
0059 #include <G4Scheduler.hh>
0060 #include <iomanip>
0061
0062
0063
0064 Run::Run(DetectorConstruction* det) : G4Run(), fDetector(det), fParticle(0)
0065 {
0066 fSoma3DEdep = new G4double[fDetector->GetnbSomacomp()];
0067 for (G4int i = 0; i < fDetector->GetnbSomacomp(); i++) {
0068 fSoma3DEdep[i] = 0.;
0069 }
0070 fDend3DEdep = new G4double[fDetector->GetnbDendritecomp()];
0071 for (G4int i = 0; i < fDetector->GetnbDendritecomp(); i++) {
0072 fDend3DEdep[i] = 0.;
0073 }
0074 fAxon3DEdep = new G4double[fDetector->GetnbAxoncomp()];
0075 for (G4int i = 0; i < fDetector->GetnbAxoncomp(); i++) {
0076 fAxon3DEdep[i] = 0.;
0077 }
0078
0079 fEdepAll = fEdepAll_err = fEdepMedium = fEdepMedium_err = fEdepSlice = fEdepSlice_err =
0080 fEdepSoma = fEdepSoma_err = 0.;
0081 fEdepDend = fEdepDend_err = fEdepAxon = fEdepAxon_err = fEdepNeuron = fEdepNeuron_err = 0.;
0082 }
0083
0084
0085
0086 Run::~Run()
0087 {
0088 delete[] fSoma3DEdep;
0089 delete[] fDend3DEdep;
0090 delete[] fAxon3DEdep;
0091 }
0092
0093
0094
0095 void Run::SetPrimary(G4ParticleDefinition* particle, G4double energy)
0096 {
0097 fParticle = particle;
0098 fEkin = energy;
0099 }
0100
0101
0102
0103 void Run::AddPrimaryLET(G4double let)
0104 {
0105 fLET += let;
0106 fLET2 += let * let;
0107 }
0108
0109
0110
0111 void Run::SetTrackLength(G4double t)
0112 {
0113 ftrackLength = t;
0114 fTrackLen += t;
0115 fTrackLen2 += t * t;
0116 }
0117
0118
0119
0120 void Run::CountProcesses(const G4VProcess* process)
0121 {
0122 G4String procName = process->GetProcessName();
0123 std::map<G4String, G4int>::iterator it = fProcCounter.find(procName);
0124 if (it == fProcCounter.end()) {
0125 fProcCounter[procName] = 1;
0126 }
0127 else {
0128 fProcCounter[procName]++;
0129 }
0130 }
0131
0132
0133
0134 void Run::ParticleCount(G4String name, G4double Ekin)
0135 {
0136 std::map<G4String, ParticleData>::iterator it = fParticleDataMap1.find(name);
0137 if (it == fParticleDataMap1.end()) {
0138 fParticleDataMap1[name] = ParticleData(1, Ekin, Ekin, Ekin);
0139 }
0140 else {
0141 ParticleData& data = it->second;
0142 data.fCount++;
0143 data.fEmean += Ekin;
0144
0145 G4double emin = data.fEmin;
0146 if (Ekin < emin) data.fEmin = Ekin;
0147 G4double emax = data.fEmax;
0148 if (Ekin > emax) data.fEmax = Ekin;
0149 }
0150 }
0151
0152
0153
0154 void Run::ParticleCountNeuron(G4String name, G4double Ekin)
0155 {
0156 std::map<G4String, ParticleData>::iterator it = fParticleDataMap2.find(name);
0157 if (it == fParticleDataMap2.end()) {
0158 fParticleDataMap2[name] = ParticleData(1, Ekin, Ekin, Ekin);
0159 }
0160 else {
0161 ParticleData& data = it->second;
0162 data.fCount++;
0163 data.fEmean += Ekin;
0164
0165 G4double emin = data.fEmin;
0166 if (Ekin < emin) data.fEmin = Ekin;
0167 G4double emax = data.fEmax;
0168 if (Ekin > emax) data.fEmax = Ekin;
0169 }
0170 }
0171
0172
0173
0174
0175
0176
0177
0178
0179
0180
0181
0182 void Run::MoleculeCountNeuron(G4Molecule* molecule)
0183 {
0184 G4String moleculename = molecule->GetName();
0185 std::map<G4String, G4int>::iterator it = fMoleculeNumber.find(moleculename);
0186 if (it == fMoleculeNumber.end()) {
0187 fMoleculeNumber[moleculename] = 1;
0188 }
0189 else {
0190 fMoleculeNumber[moleculename]++;
0191 }
0192 }
0193
0194
0195
0196 void Run::AddEflow(G4double eflow)
0197 {
0198 fEnergyFlow += eflow;
0199 fEnergyFlow2 += eflow * eflow;
0200 }
0201
0202
0203
0204 void Run::Merge(const G4Run* run)
0205 {
0206 const Run* localRun = static_cast<const Run*>(run);
0207
0208
0209
0210 fParticle = localRun->fParticle;
0211 fEkin = localRun->fEkin;
0212 ftrackLength = localRun->ftrackLength;
0213 fTrackLen += localRun->fTrackLen;
0214 fTrackLen2 += localRun->fTrackLen2;
0215 fLET += localRun->fLET;
0216 fLET2 += localRun->fLET2;
0217
0218
0219 std::map<G4String, G4int>::const_iterator itp;
0220 for (itp = localRun->fProcCounter.begin(); itp != localRun->fProcCounter.end(); ++itp) {
0221 G4String procName = itp->first;
0222 G4int localCount = itp->second;
0223 if (fProcCounter.find(procName) == fProcCounter.end()) {
0224 fProcCounter[procName] = localCount;
0225 }
0226 else {
0227 fProcCounter[procName] += localCount;
0228 }
0229 }
0230
0231
0232 std::map<G4String, ParticleData>::const_iterator itc;
0233 for (itc = localRun->fParticleDataMap1.begin(); itc != localRun->fParticleDataMap1.end(); ++itc) {
0234 G4String name = itc->first;
0235 const ParticleData& localData = itc->second;
0236 if (fParticleDataMap1.find(name) == fParticleDataMap1.end()) {
0237 fParticleDataMap1[name] =
0238 ParticleData(localData.fCount, localData.fEmean, localData.fEmin, localData.fEmax);
0239 }
0240 else {
0241 ParticleData& data = fParticleDataMap1[name];
0242 data.fCount += localData.fCount;
0243 data.fEmean += localData.fEmean;
0244 G4double emin = localData.fEmin;
0245 if (emin < data.fEmin) data.fEmin = emin;
0246 G4double emax = localData.fEmax;
0247 if (emax > data.fEmax) data.fEmax = emax;
0248 }
0249 }
0250
0251
0252 std::map<G4String, ParticleData>::const_iterator itn;
0253 for (itn = localRun->fParticleDataMap2.begin(); itn != localRun->fParticleDataMap2.end(); ++itn) {
0254 G4String name = itn->first;
0255 const ParticleData& localData = itn->second;
0256 if (fParticleDataMap2.find(name) == fParticleDataMap2.end()) {
0257 fParticleDataMap2[name] =
0258 ParticleData(localData.fCount, localData.fEmean, localData.fEmin, localData.fEmax);
0259 }
0260 else {
0261 ParticleData& data = fParticleDataMap2[name];
0262 data.fCount += localData.fCount;
0263 data.fEmean += localData.fEmean;
0264 G4double emin = localData.fEmin;
0265 if (emin < data.fEmin) data.fEmin = emin;
0266 G4double emax = localData.fEmax;
0267 if (emax > data.fEmax) data.fEmax = emax;
0268 }
0269 }
0270
0271 std::map<G4String, G4int>::const_iterator itm;
0272 for (itm = localRun->fMoleculeNumber.begin(); itm != localRun->fMoleculeNumber.end(); ++itm) {
0273 G4String moleculeName = itm->first;
0274 G4int localCount = itm->second;
0275 if (fMoleculeNumber.find(moleculeName) == fMoleculeNumber.end()) {
0276 fMoleculeNumber[moleculeName] = localCount;
0277 }
0278 else {
0279 fMoleculeNumber[moleculeName] += localCount;
0280 }
0281 }
0282
0283
0284
0285 for (G4int i = 0; i < fDetector->GetnbSomacomp(); i++) {
0286 fSoma3DEdep[i] += localRun->fSoma3DEdep[i];
0287 }
0288 for (G4int i = 0; i < fDetector->GetnbDendritecomp(); i++) {
0289 fDend3DEdep[i] += localRun->fDend3DEdep[i];
0290 }
0291 for (G4int i = 0; i < fDetector->GetnbAxoncomp(); i++) {
0292 fAxon3DEdep[i] += localRun->fAxon3DEdep[i];
0293 }
0294
0295
0296
0297 fEdepAll += localRun->fEdepAll;
0298 fEdepAll_err += localRun->fEdepAll_err;
0299 fEdepMedium += localRun->fEdepMedium;
0300 fEdepMedium_err += localRun->fEdepMedium_err;
0301 fEdepSlice += localRun->fEdepSlice;
0302 fEdepSlice_err += localRun->fEdepSlice_err;
0303 fEdepSoma += localRun->fEdepSoma;
0304 fEdepSoma_err += localRun->fEdepSoma_err;
0305 fEdepDend += localRun->fEdepDend;
0306 fEdepDend_err += localRun->fEdepDend_err;
0307 fEdepAxon += localRun->fEdepAxon;
0308 fEdepAxon_err += localRun->fEdepAxon_err;
0309 fEdepNeuron += localRun->fEdepNeuron;
0310 fEdepNeuron_err += localRun->fEdepNeuron_err;
0311
0312 fEnergyFlow += localRun->fEnergyFlow;
0313 fEnergyFlow2 += localRun->fEnergyFlow2;
0314
0315 G4Run::Merge(run);
0316 }
0317
0318
0319
0320 void Run::EndOfRun()
0321 {
0322 G4int prec = 5, wid = prec + 2;
0323 G4int dfprec = G4cout.precision(prec);
0324
0325
0326 G4String Particle = fParticle->GetParticleName();
0327 G4double GunArea;
0328
0329 G4Material* material = fDetector->GetTargetMaterial();
0330
0331
0332
0333 fTrackLen /= numberOfEvent;
0334 fTrackLen2 /= numberOfEvent;
0335 G4double rms = fTrackLen2 - fTrackLen * fTrackLen;
0336 if (rms > 0.)
0337 rms = std::sqrt(rms);
0338 else
0339 rms = 0.;
0340
0341 G4int TotNbofEvents = numberOfEvent;
0342 G4double EdepTotal = fEdepAll;
0343 G4double EdepTotal2 = fEdepAll_err;
0344 EdepTotal /= TotNbofEvents;
0345 EdepTotal2 /= TotNbofEvents;
0346 G4double rmst = EdepTotal2 - EdepTotal * EdepTotal;
0347 if (rmst > 0.)
0348 rmst = std::sqrt(rmst);
0349 else
0350 rmst = 0.;
0351
0352
0353 G4EmCalculator emCalculator;
0354 G4double dEdxFull = 0.;
0355 if (fParticle->GetPDGCharge() != 0.) {
0356 dEdxFull = emCalculator.ComputeTotalDEDX(fEkin, fParticle, material);
0357 }
0358
0359
0360 G4double meandEdx = (EdepTotal / fTrackLen) / (keV / um);
0361 G4double meandEdxerr = (rmst / fTrackLen) / (keV / um);
0362 G4int ncols = 0;
0363 G4float tmp, En;
0364 G4int Ntrav = 0;
0365 FILE* fp = fopen("OutputPerEvent.out", "r");
0366 while (1) {
0367 ncols = fscanf(fp, " %f %f %f %f %f %f %f %f", &tmp, &tmp, &tmp, &tmp, &En, &tmp, &tmp, &tmp);
0368 if (ncols < 0) break;
0369 if (En > 0) Ntrav++;
0370 }
0371 fclose(fp);
0372
0373 GunArea = fDetector->GetTotSurfMedium();
0374
0375 G4double FluenceDoseperBeam = 0.160218 * (dEdxFull) / (GunArea * std::pow(10, 18));
0376
0377 G4cout << "\n ======= The summary of simulation results 'neuron' ========\n";
0378 G4cout << "\n Primary particle = " << Particle
0379 << "\n Kinetic energy of beam = " << fEkin / MeV << " A*MeV"
0380 << "\n Particle traversals the neuron = " << Ntrav << " of " << numberOfEvent
0381 << "\n Full LET of beam as formulas = " << dEdxFull / (keV / um) << " keV/um"
0382 << "\n Mean LET of beam as simulation = " << meandEdx << " +- " << meandEdxerr
0383 << " keV/um"
0384 << "\n Mean track length of beam = " << fTrackLen / um << " +- " << rms << " um"
0385 << "\n Particle fluence = " << numberOfEvent / (GunArea / (cm * cm))
0386 << " particles/cm^2"
0387 << "\n Fluence dose (full) = "
0388 << numberOfEvent * FluenceDoseperBeam / (joule / kg) << " Gy"
0389 << "\n Fluence dose ber beam = " << FluenceDoseperBeam / (joule / kg) << " Gy"
0390 << G4endl;
0391
0392 if (numberOfEvent == 0) {
0393 G4cout.precision(dfprec);
0394 return;
0395 }
0396
0397
0398
0399 G4cout << "\n List of generated physical process:" << G4endl;
0400
0401 G4int index = 0;
0402 std::map<G4String, G4int>::iterator it;
0403 for (it = fProcCounter.begin(); it != fProcCounter.end(); it++) {
0404 G4String procName = it->first;
0405 G4int count = it->second;
0406 G4String space = " ";
0407 if (++index % 1 == 0) space = "\n";
0408 G4cout << " " << std::setw(20) << procName << "=" << std::setw(7) << count << space;
0409 }
0410 G4cout << G4endl;
0411
0412
0413
0414 G4cout << "\n List of generated particles outside neuron structure:" << G4endl;
0415
0416 std::map<G4String, ParticleData>::iterator itc;
0417 for (itc = fParticleDataMap1.begin(); itc != fParticleDataMap1.end(); itc++) {
0418 G4String name = itc->first;
0419 ParticleData data = itc->second;
0420 G4int count = data.fCount;
0421 G4double eMean = data.fEmean / count;
0422 G4double eMin = data.fEmin;
0423 G4double eMax = data.fEmax;
0424
0425 G4double Eflow = data.fEmean / TotNbofEvents;
0426
0427 G4cout << " " << std::setw(13) << name << " : " << std::setw(7) << count
0428 << " Emean = " << std::setw(wid) << G4BestUnit(eMean, "Energy") << "\t( "
0429 << G4BestUnit(eMin, "Energy") << " --> " << G4BestUnit(eMax, "Energy")
0430 << ") \tEflow/event = " << G4BestUnit(Eflow, "Energy") << G4endl;
0431 }
0432
0433
0434
0435 G4cout << "\n Number of secondary particles inside neuron structure:" << G4endl;
0436
0437 std::map<G4String, ParticleData>::iterator itn;
0438 for (itn = fParticleDataMap2.begin(); itn != fParticleDataMap2.end(); itn++) {
0439 G4String name = itn->first;
0440 ParticleData data = itn->second;
0441 G4int count = data.fCount;
0442
0443 G4cout << " " << std::setw(13) << name << " : " << std::setw(7) << count << G4endl;
0444 }
0445
0446
0447
0448 G4cout << "\n Number of molecular products inside neuron structure:"
0449 << "\n time: 1 ps - 10 ps " << G4endl;
0450
0451 G4int ind = 0;
0452 std::map<G4String, G4int>::iterator itm;
0453 for (itm = fMoleculeNumber.begin(); itm != fMoleculeNumber.end(); itm++) {
0454 G4String moleculeName = itm->first;
0455 G4int count = itm->second;
0456 G4String space = " ";
0457 if (++ind % 3 == 0) space = "\n";
0458
0459 G4cout << " " << std::setw(13) << moleculeName << " : " << std::setw(7) << count << G4endl;
0460 }
0461
0462
0463 G4cout << "\n Total energy (MeV) deposition :" << G4endl;
0464
0465 G4cout << " " << std::setw(13) << "All volume: " << std::setw(7) << fEdepAll / MeV << "\n "
0466 << " " << std::setw(13) << "Bounding slice: " << std::setw(7)
0467 << (fEdepSlice + fEdepNeuron) / MeV << "\n "
0468 << " " << std::setw(13) << "Neuron: " << std::setw(7) << fEdepNeuron / MeV << "\n "
0469 << " " << std::setw(13) << "Soma: " << std::setw(7) << fEdepSoma / MeV << "\n "
0470 << " " << std::setw(13) << "Dendrites: " << std::setw(7) << fEdepDend / MeV << "\n "
0471 << " " << std::setw(13) << "Axon: " << std::setw(7) << fEdepAxon / MeV << G4endl;
0472
0473
0474
0475
0476 G4int somaCompHit = 0;
0477 G4double somaCompEdep = 0.;
0478 G4double somaCompDose = 0.;
0479 G4double somaCompEdep2 = 0.;
0480 G4double somaCompDose2 = 0.;
0481
0482 remove("Soma3DEdep.out");
0483 for (G4int i = 0; i < fDetector->GetnbSomacomp(); i++) {
0484 if (fSoma3DEdep[i] > 0.0) {
0485 somaCompHit++;
0486 somaCompEdep += fSoma3DEdep[i];
0487 somaCompEdep2 += fSoma3DEdep[i] * fSoma3DEdep[i];
0488
0489 std::ofstream WriteDataInSoma("Soma3DEdep.out", std::ios::app);
0490
0491 WriteDataInSoma << fDetector->GetPosSomacomp(i).x() << '\t' << " "
0492 << fDetector->GetPosSomacomp(i).y() << '\t' << " "
0493 << fDetector->GetPosSomacomp(i).z() << '\t'
0494 << " "
0495
0496 << fSoma3DEdep[i] / keV << '\t' << " " << G4endl;
0497 }
0498 }
0499
0500
0501 G4double rmsEdepS = 0.;
0502 G4double rmsDoseS = 0.;
0503 if (somaCompHit > 0) {
0504 somaCompEdep /= somaCompHit;
0505 somaCompEdep2 /= somaCompHit;
0506 rmsEdepS = somaCompEdep2 - somaCompEdep * somaCompEdep;
0507 if (rmsEdepS > 0.)
0508 rmsEdepS = std::sqrt(rmsEdepS);
0509 else
0510 rmsEdepS = 0.;
0511 somaCompDose /= somaCompHit;
0512 somaCompDose2 /= somaCompHit;
0513 rmsDoseS = somaCompDose2 - somaCompDose * somaCompDose;
0514 if (rmsDoseS > 0.)
0515 rmsDoseS = std::sqrt(rmsDoseS);
0516 else
0517 rmsDoseS = 0.;
0518 }
0519
0520 G4int DendCompHit = 0;
0521 G4double DendCompEdep = 0.;
0522 G4double DendCompDose = 0.;
0523 G4double DendCompEdep2 = 0.;
0524 G4double DendCompDose2 = 0.;
0525 remove("Dend3DEdep.out");
0526 std::ofstream WriteDataInDend("Dend3DEdep.out", std::ios::app);
0527 for (G4int i = 0; i < fDetector->GetnbDendritecomp(); i++) {
0528 if (fDend3DEdep[i] > 0.0) {
0529 DendCompHit++;
0530 DendCompEdep += fDend3DEdep[i];
0531 DendCompEdep2 += fDend3DEdep[i] * fDend3DEdep[i];
0532
0533 WriteDataInDend
0534
0535 << fDetector->GetPosDendcomp(i).x() << '\t' << " " << fDetector->GetPosDendcomp(i).y()
0536 << '\t' << " " << fDetector->GetPosDendcomp(i).z() << '\t' << " "
0537 << fDetector->GetDistADendSoma(i) << '\t' << " " << fDetector->GetDistBDendSoma(i) << '\t'
0538 << " " << fDend3DEdep[i] / keV << '\t' << " " << G4endl;
0539 }
0540 }
0541
0542 G4double rmsEdepD = 0.;
0543 G4double rmsDoseD = 0.;
0544 if (DendCompHit > 0) {
0545 DendCompEdep /= DendCompHit;
0546 DendCompEdep2 /= DendCompHit;
0547 rmsEdepD = DendCompEdep2 - DendCompEdep * DendCompEdep;
0548 if (rmsEdepD > 0.)
0549 rmsEdepD = std::sqrt(rmsEdepD);
0550 else
0551 rmsEdepD = 0.;
0552 DendCompDose /= DendCompHit;
0553 DendCompDose2 /= DendCompHit;
0554 rmsDoseD = DendCompDose2 - DendCompDose * DendCompDose;
0555 if (rmsDoseD > 0.)
0556 rmsDoseD = std::sqrt(rmsDoseD);
0557 else
0558 rmsDoseD = 0.;
0559 }
0560
0561 G4int AxonCompHit = 0;
0562 G4double AxonCompEdep = 0.;
0563 G4double AxonCompDose = 0.;
0564 G4double AxonCompEdep2 = 0.;
0565 G4double AxonCompDose2 = 0.;
0566 remove("Axon3DEdep.out");
0567 std::ofstream WriteDataInAxon("Axon3DEdep.out", std::ios::app);
0568 for (G4int i = 0; i < fDetector->GetnbAxoncomp(); i++) {
0569 if (fAxon3DEdep[i] > 0.0) {
0570 AxonCompHit++;
0571 AxonCompEdep += fAxon3DEdep[i];
0572 AxonCompEdep2 += fAxon3DEdep[i] * fAxon3DEdep[i];
0573
0574 WriteDataInAxon
0575
0576 << fDetector->GetPosAxoncomp(i).x() << '\t' << " " << fDetector->GetPosAxoncomp(i).y()
0577 << '\t' << " " << fDetector->GetPosAxoncomp(i).z() << '\t' << " "
0578 << fDetector->GetDistAxonsoma(i) << '\t' << " " << fAxon3DEdep[i] / keV << '\t' << " "
0579 << G4endl;
0580 }
0581 }
0582
0583 G4double rmsEdepA = 0.;
0584 G4double rmsDoseA = 0.;
0585 if (AxonCompHit > 0) {
0586 AxonCompEdep /= AxonCompHit;
0587 AxonCompEdep2 /= AxonCompHit;
0588 rmsEdepA = AxonCompEdep2 - AxonCompEdep * AxonCompEdep;
0589 if (rmsEdepA > 0.)
0590 rmsEdepA = std::sqrt(rmsEdepA);
0591 else
0592 rmsEdepA = 0.;
0593 AxonCompDose /= AxonCompHit;
0594 AxonCompDose2 /= AxonCompHit;
0595 rmsDoseA = AxonCompDose2 - AxonCompDose * AxonCompDose;
0596 if (rmsDoseA > 0.)
0597 rmsDoseA = std::sqrt(rmsDoseA);
0598 else
0599 rmsDoseA = 0.;
0600 }
0601
0602 G4cout << "\n Number of compartments traversed by particle tracks :" << G4endl;
0603 G4cout << " " << std::setw(13) << "Soma: " << std::setw(7) << somaCompHit
0604 << " of total: " << fDetector->GetnbSomacomp() << "\n "
0605 << " " << std::setw(13) << "Dendrites: " << std::setw(7) << DendCompHit
0606 << " of total: " << fDetector->GetnbDendritecomp() << "\n "
0607 << " " << std::setw(13) << "Axon: " << std::setw(7) << AxonCompHit
0608 << " of total: " << fDetector->GetnbAxoncomp() << "\n " << G4endl;
0609 G4cout << "\n Dendritic (or Axon) compartmental energy deposits \n"
0610 << " at the distance from Soma have been written into *.out files:"
0611 << "\n Dend3DEdep.out, Axon3DEdep.out, Soma3DEdep.out"
0612 << "\n Outputs of energy deposition per event written in data file:"
0613 << "\n OutputPerEvent.out" << G4endl;
0614
0615
0616 fProcCounter.clear();
0617 fParticleDataMap2.clear();
0618
0619
0620 G4cout.precision(dfprec);
0621 }
0622
0623
