diff --git a/HEN_HOUSE/egs++/ausgab_objects/beam_dose_scoring/Makefile b/HEN_HOUSE/egs++/ausgab_objects/beam_dose_scoring/Makefile
index da5df1d20..bca2d49c2 100644
--- a/HEN_HOUSE/egs++/ausgab_objects/beam_dose_scoring/Makefile
+++ b/HEN_HOUSE/egs++/ausgab_objects/beam_dose_scoring/Makefile
@@ -40,7 +40,7 @@ my_deps = $(common_ausgab_deps)
extra_dep = $(addprefix $(DSOLIBS), $(my_deps))
#redefined from egspp.spec
-INC2 = -I$(IEGS2) -I..$(DSEP).. -I$(HEN_HOUSE)/user_codes/beampp -I$(HEN_HOUSE)interface
+INC2 = -I$(IEGS2) -I..$(DSEP).. -I$(HEN_HOUSE)interface
include $(SPEC_DIR)egspp_libs.spec
diff --git a/HEN_HOUSE/egs++/egs_advanced_application.cpp b/HEN_HOUSE/egs++/egs_advanced_application.cpp
index 678f46d42..f4001ee04 100644
--- a/HEN_HOUSE/egs++/egs_advanced_application.cpp
+++ b/HEN_HOUSE/egs++/egs_advanced_application.cpp
@@ -649,7 +649,7 @@ int EGS_AdvancedApplication::helpInit(EGS_Input *transportp, bool do_hatch) {
egsInformation("%3d %-24s AE=%7.4f AP=%7.4f %d\n",j,
geometry->getMediumName(j),the_thresh->ae[imed],
the_thresh->ap[imed],imed);
- if (Emax > the_thresh->ue[imed]-0.511 ||
+ if (Emax > the_thresh->ue[imed]-the_useful->rm ||
Emax > the_thresh->up[imed]) {
egsInformation(" upper limits (UE=%g UP=%g) not enough for "
"maximum source energy (%g)\n",the_thresh->ue[imed],
diff --git a/HEN_HOUSE/egs++/egs_advanced_application.h b/HEN_HOUSE/egs++/egs_advanced_application.h
index 8d4d7f0b3..2efd90b24 100644
--- a/HEN_HOUSE/egs++/egs_advanced_application.h
+++ b/HEN_HOUSE/egs++/egs_advanced_application.h
@@ -25,6 +25,7 @@
#
# Contributors: Ernesto Mainegra-Hing
# Frederic Tessier
+# Reid Townson
#
###############################################################################
*/
diff --git a/HEN_HOUSE/egs++/sources/egs_beam_source/Makefile b/HEN_HOUSE/egs++/sources/egs_beam_source/Makefile
index bad8fc4d1..8760e171b 100644
--- a/HEN_HOUSE/egs++/sources/egs_beam_source/Makefile
+++ b/HEN_HOUSE/egs++/sources/egs_beam_source/Makefile
@@ -39,6 +39,10 @@ lib_files = egs_beam_source
my_deps = $(common_source_deps)
extra_dep = $(addprefix $(DSOLIBS), $(my_deps))
+#redefined from egspp.spec
+INC2 = -I$(IEGS2) -I..$(DSEP).. -I$(HEN_HOUSE)interface -I.$(DSEP)
+
include $(SPEC_DIR)egspp_libs.spec
$(make_depend)
+
diff --git a/HEN_HOUSE/egs++/sources/egs_beam_source/array_sizes.h b/HEN_HOUSE/egs++/sources/egs_beam_source/array_sizes.h
new file mode 100644
index 000000000..0ccd372d5
--- /dev/null
+++ b/HEN_HOUSE/egs++/sources/egs_beam_source/array_sizes.h
@@ -0,0 +1,47 @@
+/*
+###############################################################################
+#
+# EGSnrc egs++ egs_chamber application array sizes headers
+# Copyright (C) 2015 National Research Council Canada
+#
+# This file is part of EGSnrc.
+#
+# EGSnrc is free software: you can redistribute it and/or modify it under
+# the terms of the GNU Affero General Public License as published by the
+# Free Software Foundation, either version 3 of the License, or (at your
+# option) any later version.
+#
+# EGSnrc is distributed in the hope that it will be useful, but WITHOUT ANY
+# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+# FOR A PARTICULAR PURPOSE. See the GNU Affero General Public License for
+# more details.
+#
+# You should have received a copy of the GNU Affero General Public License
+# along with EGSnrc. If not, see .
+#
+###############################################################################
+#
+# Author: Reid Townson, 2016
+#
+# Contributors:
+#
+###############################################################################
+#
+# Defines he maximum number of media (MXMED) and the maximum number of
+# particles on the stack (MXSTACK). This file gets included by the egsnrc
+# fortran subroutines (egsnrc_$my_machine.F), the base application
+# (egs_simple_application.cpp or egs_advanced_application.cpp in
+# $HEN_HOUSE/egs++), and possibly the user code, if it uses the particle
+# stack or one of the structures that depends on the maximum number of media.
+#
+###############################################################################
+*/
+
+
+#ifndef ARRAY_SIZES_
+#define ARRAY_SIZES_
+
+#define MXMED 0
+#define MXSTACK 0
+
+#endif
diff --git a/HEN_HOUSE/egs++/sources/egs_beam_source/egs_beam_source.cpp b/HEN_HOUSE/egs++/sources/egs_beam_source/egs_beam_source.cpp
index f00cdd99c..8af8bd913 100644
--- a/HEN_HOUSE/egs++/sources/egs_beam_source/egs_beam_source.cpp
+++ b/HEN_HOUSE/egs++/sources/egs_beam_source/egs_beam_source.cpp
@@ -25,6 +25,7 @@
#
# Contributors: Blake Walters
# Frederic Tessier
+# Reid Townson
#
###############################################################################
*/
@@ -40,6 +41,7 @@
#include "egs_functions.h"
#include "egs_library.h"
#include "egs_application.h"
+#include "egs_interface2.h"
#include
#include
@@ -201,7 +203,7 @@ EGS_I64 EGS_BeamSource::getNextParticle(EGS_RandomGenerator *, int &q,
// " x=(%g,%g,%g) latch=%d count=%lld\n",te,tq,twt,tx,ty,tz,
// tlatch,count);
if (tq) {
- te -= 0.5110034;
+ te -= the_useful->rm;
}
ok = true;
if (te > Emax) {
diff --git a/HEN_HOUSE/egs++/sources/egs_phsp_source/Makefile b/HEN_HOUSE/egs++/sources/egs_phsp_source/Makefile
index 5d62234d4..4bdc03393 100644
--- a/HEN_HOUSE/egs++/sources/egs_phsp_source/Makefile
+++ b/HEN_HOUSE/egs++/sources/egs_phsp_source/Makefile
@@ -39,6 +39,9 @@ lib_files = egs_phsp_source
my_deps = $(common_source_deps)
extra_dep = $(addprefix $(DSOLIBS), $(my_deps))
+#redefined from egspp.spec
+INC2 = -I$(IEGS2) -I..$(DSEP).. -I$(HEN_HOUSE)interface -I.$(DSEP)
+
include $(SPEC_DIR)egspp_libs.spec
$(make_depend)
diff --git a/HEN_HOUSE/egs++/sources/egs_phsp_source/array_sizes.h b/HEN_HOUSE/egs++/sources/egs_phsp_source/array_sizes.h
new file mode 100644
index 000000000..0ccd372d5
--- /dev/null
+++ b/HEN_HOUSE/egs++/sources/egs_phsp_source/array_sizes.h
@@ -0,0 +1,47 @@
+/*
+###############################################################################
+#
+# EGSnrc egs++ egs_chamber application array sizes headers
+# Copyright (C) 2015 National Research Council Canada
+#
+# This file is part of EGSnrc.
+#
+# EGSnrc is free software: you can redistribute it and/or modify it under
+# the terms of the GNU Affero General Public License as published by the
+# Free Software Foundation, either version 3 of the License, or (at your
+# option) any later version.
+#
+# EGSnrc is distributed in the hope that it will be useful, but WITHOUT ANY
+# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+# FOR A PARTICULAR PURPOSE. See the GNU Affero General Public License for
+# more details.
+#
+# You should have received a copy of the GNU Affero General Public License
+# along with EGSnrc. If not, see .
+#
+###############################################################################
+#
+# Author: Reid Townson, 2016
+#
+# Contributors:
+#
+###############################################################################
+#
+# Defines he maximum number of media (MXMED) and the maximum number of
+# particles on the stack (MXSTACK). This file gets included by the egsnrc
+# fortran subroutines (egsnrc_$my_machine.F), the base application
+# (egs_simple_application.cpp or egs_advanced_application.cpp in
+# $HEN_HOUSE/egs++), and possibly the user code, if it uses the particle
+# stack or one of the structures that depends on the maximum number of media.
+#
+###############################################################################
+*/
+
+
+#ifndef ARRAY_SIZES_
+#define ARRAY_SIZES_
+
+#define MXMED 0
+#define MXSTACK 0
+
+#endif
diff --git a/HEN_HOUSE/egs++/sources/egs_phsp_source/egs_phsp_source.cpp b/HEN_HOUSE/egs++/sources/egs_phsp_source/egs_phsp_source.cpp
index 5fcb35354..7d8da9fa7 100644
--- a/HEN_HOUSE/egs++/sources/egs_phsp_source/egs_phsp_source.cpp
+++ b/HEN_HOUSE/egs++/sources/egs_phsp_source/egs_phsp_source.cpp
@@ -25,6 +25,7 @@
#
# Contributors: Ernesto Mainegra-Hing
# Frederic Tessier
+# Reid Townson
#
###############################################################################
*/
@@ -38,6 +39,7 @@
#include "egs_phsp_source.h"
#include "egs_input.h"
#include "egs_functions.h"
+#include "egs_interface2.h"
EGS_PhspSource::EGS_PhspSource(const string &phsp_file,
const string &Name, EGS_ObjectFactory *f) : EGS_BaseSource(Name,f) {
@@ -396,7 +398,7 @@ void EGS_PhspSource::readParticle() {
}
first = false;
if (p.q) {
- p.E -= 0.5110034;
+ p.E -= the_useful->rm;
Nreuse = Nreuse_e;
}
else {
diff --git a/HEN_HOUSE/gui/egs_gui/pegs_page.cpp b/HEN_HOUSE/gui/egs_gui/pegs_page.cpp
index 0553c0f70..d9b69cc9d 100644
--- a/HEN_HOUSE/gui/egs_gui/pegs_page.cpp
+++ b/HEN_HOUSE/gui/egs_gui/pegs_page.cpp
@@ -23,7 +23,7 @@
#
# Author: Ernesto Mainegra-Hing, 2015
#
-# Contributors:
+# Contributors: Reid Townson
#
###############################################################################
*/
@@ -56,7 +56,7 @@ QStringList *elements = 0;
TableEventHandler *tfilter = 0;
-const double rm=0.5110034;
+const double rm=0.5109989461;
Element element_data[] = {
{1,"H",1.0079,19.2,8.3748e-05},
diff --git a/HEN_HOUSE/interface/egs_c_interface2.macros b/HEN_HOUSE/interface/egs_c_interface2.macros
index ad4115c32..08b063415 100644
--- a/HEN_HOUSE/interface/egs_c_interface2.macros
+++ b/HEN_HOUSE/interface/egs_c_interface2.macros
@@ -25,6 +25,7 @@
" "
" Contributors: Ernesto Mainegra-Hing "
" Blake Walters "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -444,7 +445,7 @@ REPLACE {$COMIN-SET-DEFAULTS;} WITH {;
};
REPLACE {$COMIN-INIT-COMPT;} WITH {;
;integer max_med; parameter (max_med = MXMED);
- ;COMIN/COMPTON-DATA,BREMPR,EDGE,MEDIA,MISC,EGS-IO/;
+ ;COMIN/COMPTON-DATA,BREMPR,EDGE,MEDIA,MISC,USEFUL,EGS-IO/;
};
REPLACE {$COMIN-COMPT;} WITH {;
;integer max_med; parameter (max_med = MXMED);
@@ -459,7 +460,7 @@ REPLACE {$COMIN-ELECTR;} WITH {;
};
REPLACE {$COMIN-MSCATI;} WITH {;
;integer max_med; parameter (max_med = MXMED);
- ;COMIN/BOUNDS,ELECIN,MEDIA,MISC,RANDOM,ET-Control,EGS-IO/;
+ ;COMIN/BOUNDS,ELECIN,MEDIA,MISC,RANDOM,ET-Control,USEFUL,EGS-IO/;
};
REPLACE {$COMIN-MOLLER;} WITH {;
;integer max_med; parameter (max_med = MXMED);
diff --git a/HEN_HOUSE/omega/beamnrc/beamnrc.mortran b/HEN_HOUSE/omega/beamnrc/beamnrc.mortran
index 1b4b251c8..b26c20899 100644
--- a/HEN_HOUSE/omega/beamnrc/beamnrc.mortran
+++ b/HEN_HOUSE/omega/beamnrc/beamnrc.mortran
@@ -36,6 +36,7 @@
" Iwan Kawrakow "
" Ernesto Mainegra-Hing "
" Frederic Tessier "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -5716,8 +5717,8 @@ IF (IARG = 5) ["Particle step just taken"
IF(W(NP)>0)["electron going forward, split it"
count_esplit=count_esplit+1;
IF(IWATCH=1|IWATCH=2)[
- OUTPUT NP,E(NP)-0.511,IQ(NP),IR(NP),X(NP),Y(NP),Z(NP),U(NP),V(NP),
- W(NP),LATCH(NP),WT(NP);
+ OUTPUT NP,E(NP)-PRM,IQ(NP),IR(NP),X(NP),Y(NP),Z(NP),
+ U(NP),V(NP),W(NP),LATCH(NP),WT(NP);
(' Electron',T36,':',I5,F9.3,2I4,3F8.3,3F7.3,I10,1PE10.3);
OUTPUT NBRSPL;(' about to be split ',I10, ' times (DBS).');
IF(IRAD_DBS=1)[
@@ -5731,8 +5732,8 @@ IF (IARG = 5) ["Particle step just taken"
ANG_DBS=2.*3.1415926/FLOAT(NBRSPL);
]
IF(IWATCH=1|IWATCH=2)[
- OUTPUT NP,E(NP)-0.511,IQ(NP),IR(NP),X(NP),Y(NP),Z(NP),U(NP),V(NP),
- W(NP),LATCH(NP),WT(NP);
+ OUTPUT NP,E(NP)-PRM,IQ(NP),IR(NP),X(NP),Y(NP),Z(NP),
+ U(NP),V(NP),W(NP),LATCH(NP),WT(NP);
(T36,':',I5,F9.3,2I4,3F8.3,3F7.3,I10,1PE10.3);
]
DO I=1,NBRSPL-1[
@@ -5750,8 +5751,8 @@ IF (IARG = 5) ["Particle step just taken"
U(NP)=U(NP-1);V(NP)=V(NP-1);
]
IF(IWATCH=1|IWATCH=2)[
- OUTPUT NP,E(NP)-0.511,IQ(NP),IR(NP),X(NP),Y(NP),Z(NP),U(NP),V(NP),
- W(NP),LATCH(NP),WT(NP);
+ OUTPUT NP,E(NP)-PRM,IQ(NP),IR(NP),X(NP),Y(NP),Z(NP),
+ U(NP),V(NP),W(NP),LATCH(NP),WT(NP);
(T36,':',I5,F9.3,2I4,3F8.3,3F7.3,I10,1PE10.3);
]
]
@@ -5874,7 +5875,7 @@ IF (IARG = 5) ["Particle step just taken"
"Energy fluence:"
SCFLU_TMP(IQL,ISCORE,ISZ,IPAR+3)=SCFLU_TMP(IQL,ISCORE,ISZ,IPAR+3)+
- WEIGHT*(E(NP)-0.511*ABS(IQ(NP)))/
+ WEIGHT*(E(NP)-PRM*ABS(IQ(NP)))/
MAX(0.08716,ABS(WL));
"Angle:"
"change to degrees on output"
@@ -5894,7 +5895,7 @@ IF (IARG = 5) ["Particle step just taken"
SCFLU_COV(IQL,ISCORE,ISZ,IPAR/2+1)=SCFLU_COV(IQL,ISCORE,ISZ,IPAR/2+1)+
SCFLU_TMP(IQL,ISCORE,ISZ,IPAR+3)*
SCFLU_TMP(IQL,ISCORE,ISZ,IPAR+2);
- SCFLU_TMP(IQL,ISCORE,ISZ,IPAR+3)=WEIGHT*(E(NP)-0.511*ABS(IQ(NP)))/
+ SCFLU_TMP(IQL,ISCORE,ISZ,IPAR+3)=WEIGHT*(E(NP)-PRM*ABS(IQ(NP)))/
MAX(0.08716,ABS(WL));
"Angle: also score angle*no. (weighted) now for covariance"
@@ -5958,12 +5959,12 @@ IF (IARG = 5) ["Particle step just taken"
IF (IO_OPT = 0)|((IO_OPT = 3)&(IHSTRY <= 100000))|(IO_OPT=4)[
"phase-space output"
NPPHSP(ISCORE)=NPPHSP(ISCORE)+1;"increment # of particles in file"
- IF((E(NP)-0.511*ABS(IQ(NP)))>EKMAXPHSP(ISCORE))[
- EKMAXPHSP(ISCORE)=E(NP)-0.511*ABS(IQ(NP));
+ IF((E(NP)-PRM*ABS(IQ(NP)))>EKMAXPHSP(ISCORE))[
+ EKMAXPHSP(ISCORE)=E(NP)-PRM*ABS(IQ(NP));
]
IF(IQ(NP)=-1)[
- IF((E(NP)-0.511)****************************
diff --git a/HEN_HOUSE/omega/progs/beamdp/beamdp.mortran b/HEN_HOUSE/omega/progs/beamdp/beamdp.mortran
index 7cf378343..af3a8f83d 100644
--- a/HEN_HOUSE/omega/progs/beamdp/beamdp.mortran
+++ b/HEN_HOUSE/omega/progs/beamdp/beamdp.mortran
@@ -30,6 +30,7 @@
" Zdenko Sego "
" Ernesto Mainegra-Hing "
" Frederic Tessier "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -2125,7 +2126,7 @@ integer lnblnk1;
]
IF(IQ ~= 0)[ "We need kinetic energy only"
- EI=EI-0.5110034;
+ EI=EI-0.5109989461;
]
"processing the ph-sp data"
@@ -4835,14 +4836,14 @@ LOOP["read phase-space data from the data file"
]
IF(IQ = -1)[ "We need kinetic energy only"
- EI=EI-0.5110034;
+ EI=EI-0.5109989461;
N_ph_sp_e = N_ph_sp_e + 1; E_ph_sp_e = E_ph_sp_e+EI;
]
ELSEIF(IQ = 0)[
N_ph_sp_g = N_ph_sp_g + 1; E_ph_sp_g = E_ph_sp_g+EI;
]
ELSEIF(IQ = +1)[
- EI=EI-0.5110034;
+ EI=EI-0.5109989461;
N_ph_sp_p = N_ph_sp_p + 1; E_ph_sp_p = E_ph_sp_p+EI;
]
ELSE[
diff --git a/HEN_HOUSE/omega/progs/beamdp/beammodel_macros.mortran b/HEN_HOUSE/omega/progs/beamdp/beammodel_macros.mortran
index 099073307..3bf57d0cb 100644
--- a/HEN_HOUSE/omega/progs/beamdp/beammodel_macros.mortran
+++ b/HEN_HOUSE/omega/progs/beamdp/beammodel_macros.mortran
@@ -27,6 +27,7 @@
" Iwan Kawrakow "
" Zdenko Sego "
" Dave Rogers "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -314,10 +315,10 @@ REPLACE {;$BEAMMODEL-SOURCE4-OUTPUT} WITH {
}
;
REPLACE{$BEAMMODEL-SOURCE4-ENERGY} WITH {
- ;IF(iqin ~= 0)[etotin = einsrc +0.511;]
+ ;IF(iqin ~= 0)[etotin = einsrc +0.5109989461;]
ELSE[etotin = einsrc;]
esrc = esrc + einsrc;
- IF(iqin =1)[esrc = esrc + 0.5110034;]
+ IF(iqin =1)[esrc = esrc + 0.5109989461;]
}
;
REPLACE {$BEAMMODEL_SOURCE_OUTPUT} WITH
diff --git a/HEN_HOUSE/omega/progs/readphsp/readphsp.mortran b/HEN_HOUSE/omega/progs/readphsp/readphsp.mortran
index cd1998b8a..fe5d15ee0 100644
--- a/HEN_HOUSE/omega/progs/readphsp/readphsp.mortran
+++ b/HEN_HOUSE/omega/progs/readphsp/readphsp.mortran
@@ -29,6 +29,7 @@
" Daryoush Sheikh-Bagheri "
" Blake Walters "
" Iwan Kawrakow "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -633,7 +634,7 @@ LOOP[ "read and write one particle at a time
ELSE[" PAW FORMAT"
"For this output only, we need to recover the original"
"phase space values, not the compressed values"
- IF( abs(IQ) = 1) [ ESHORT=ESHORT-0.511; ]
+ IF( abs(IQ) = 1) [ ESHORT=ESHORT-0.5109989461; ]
"until Feb 1995, this was if(iq = -1) => e+ had total energy"
"RWN XTUPLE(1) = ESHORT; XTUPLE(2) = X; XTUPLE(3) = Y; XTUPLE(4) = U;
"RWN XTUPLE(5) = V;
diff --git a/HEN_HOUSE/src/egs_utilities.mortran b/HEN_HOUSE/src/egs_utilities.mortran
index 6a62c05fa..3da17a149 100644
--- a/HEN_HOUSE/src/egs_utilities.mortran
+++ b/HEN_HOUSE/src/egs_utilities.mortran
@@ -27,6 +27,7 @@
" Ernesto Mainegra-Hing "
" Blake Walters "
" Frederic Tessier "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -938,7 +939,12 @@ nmed = $default_nmed;
kmpi = 12; kmpo = 8; dunit = 1;
rng_seed = 999999;
latchi = 0;
-rm = 0.5110034; rmt2 = 2*rm; rmsq = rm*rm;
+
+" The rest mass value is as recommended by CODATA 2014
+" http://physics.nist.gov/cgi-bin/cuu/Value?mec2mev
+rm = 0.5109989461;
+rmt2 = 2*rm; rmsq = rm*rm;
+
pi = 4*datan(1d0); twopi = 2*pi; pi5d2 = 2.5*pi;
nbr_split = 1;
i_play_RR = 0; i_survived_RR = 0; prob_RR = -1; n_RR_warning = 0;
diff --git a/HEN_HOUSE/src/egsnrc.macros b/HEN_HOUSE/src/egsnrc.macros
index e713dc830..f284a4241 100644
--- a/HEN_HOUSE/src/egsnrc.macros
+++ b/HEN_HOUSE/src/egsnrc.macros
@@ -27,6 +27,7 @@
" Contributors: Blake Walters "
" Ernesto Mainegra-Hing "
" Frederic Tessier "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -1128,9 +1129,9 @@ REPLACE {$COMIN-SET-DEFAULTS;} WITH {
EGS-VARIANCE-REDUCTION,EGS-IO,Spin-Data,EII-DATA,rayleigh_inputs,
EMF-INPUTS/;};
REPLACE {$COMIN-INIT-COMPT;} WITH {
- ;COMIN/COMPTON-DATA,BREMPR,EDGE,MEDIA,MISC,EGS-IO/;};
+ ;COMIN/COMPTON-DATA,BREMPR,EDGE,MEDIA,MISC,USEFUL,EGS-IO/;};
REPLACE {$COMIN-MSCATI;} WITH {
- ;COMIN/BOUNDS,ELECIN,MEDIA,MISC,RANDOM,ET-Control,EGS-IO/;};
+ ;COMIN/BOUNDS,ELECIN,MEDIA,MISC,RANDOM,ET-Control,USEFUL,EGS-IO/;};
REPLACE {$COMIN-INIT-TRIPLET;} WITH {
;COMIN/BREMPR,EGS-IO,MEDIA,TRIPLET-DATA,USEFUL/;};
REPLACE {$COMIN-GET-TRANSPORTP;} WITH {
@@ -2140,7 +2141,7 @@ REPLACE {$RE-EVALUATE-DEDX(#);} WITH
lelktmp = max(1,lelktmp);
IF(lelec < 0)[$EVALUATE dedxmid USING ededx(elktmp);]
ELSE [$EVALUATE dedxmid USING pdedx(elktmp);]
- dedx = rhof*dedxmid*(1+0.17408298*({P1}/(1-{P1})/(eke+0.5110034))**2);
+ dedx = rhof*dedxmid*(1+0.17408298*({P1}/(1-{P1})/(eke+PRM))**2);
"0.17408298 is 2/3*m**2"
{P1} = 2*{P1};
}
@@ -3443,7 +3444,7 @@ REPLACE {$RADC_WARNING;} WITH {;
};
REPLACE {$RADC_HATCH;} WITH {$RADC_WARNING;}
REPLACE {$COMIN-RADC-INIT;} WITH {
- ;COMIN/RAD_COMPTON,EGS-IO,COMPTON-DATA,MEDIA,PHOTIN/;
+ ;COMIN/RAD_COMPTON,EGS-IO,COMPTON-DATA,MEDIA,PHOTIN,USEFUL/;
};
REPLACE {$COMIN-RADC-SAMPLE;} WITH {
;COMIN/RAD_COMPTON,RANDOM,STACK,USEFUL,EGS-IO/;
diff --git a/HEN_HOUSE/src/egsnrc.mortran b/HEN_HOUSE/src/egsnrc.mortran
index e7a8dfee4..b5ca3738e 100644
--- a/HEN_HOUSE/src/egsnrc.mortran
+++ b/HEN_HOUSE/src/egsnrc.mortran
@@ -1840,6 +1840,11 @@ $INTEGER lnblnk1;" in house lnblnk function becuase not all compilers"
DATA MDLABL/$S' MEDIUM='/,LMDL/8/,LMDN/24/,DUNITO/1./;
DATA I1ST/1/,NSINSS/37/,MXSINC/$MXSINC/,ISTEST/0/,NRNA/1000/;
+" SET UP ENERGY PRECISION VARIABLES"
+PRM=RM; "PRECISE REST MASS"
+PRMT2=2.D0*PRM; "TWICE THE PRECISION REST MASS"
+PZERO=0.0D0; "PRECISE ZERO"
+
$INIT-PEGS4-VARIABLES;
" FORMAT STATEMENTS USED MULTIPLE TIMES IN SETUP"
@@ -1852,11 +1857,6 @@ IF (I1ST.NE.0)[ I1ST=0;"RESET FIRST TIME FLAG"
$HATCH-USER-INPUT-INIT;
-" SET UP ENERGY PRECISION VARIABLES"
-PRM=RM; "PRECISE REST MASS"
-PRMT2=2.D0*PRM; "TWICE THE PRECISION REST MASS"
-PZERO=0.0D0; "PRECISE ZERO"
-
" NOW CONSTRUCT PIECEWISE LINEAR FIT TO SINE FUNCTION OVER THE"
" INTERVAL (0,5*PI/2). DIVIDE THIS INTERVAL INTO MXSINC SUB-"
" INTERVALS. EACH OF THESE SUBINTERVALS IS THEN SUBDIVIDED INTO"
@@ -2413,9 +2413,8 @@ implicit none;
$COMIN-INIT-COMPT;
-$INTEGER i,j,iz,medium,nsh,j_l,j_h;
-$REAL aux,pztot,rm,atav;
-parameter (rm=0.5110034);
+$INTEGER i,j,iz,nsh,j_l,j_h;
+$REAL aux,pztot,atav;
$REAL aux_erf,erf1; "using erf1, provided with EGSnrc, because some"
"compiler don't have an intrinsic error function"
$LOGICAL getd;
@@ -2436,7 +2435,7 @@ DO j=1,$MXTOTSH [
read($INCOHUNIT,*) iz_array(j),shn_array(j),ne_array(j),
Jo_array(j),be_array(j);
Jo_array(j) = Jo_array(j)*137.;
- be_array(j) = be_array(j)*1e-6/rm;
+ be_array(j) = be_array(j)*1e-6/PRM;
aux_erf = 0.70710678119*(1+0.3*Jo_array(j));
erfJo_array(j) = 0.82436063535*(erf1(aux_erf)-1);
"0.82436063535 is exp(0.5)/2"
@@ -2472,7 +2471,7 @@ DO medium = 1,nmed [
eno_array(i,medium) = eno_array(i,medium)/pztot;
$egs_info('(i3,i4,i3,f9.5,e10.3,f10.3)',
i,j,shn_array(j),eno_array(i,medium),
- Jo_array(j),be_array(j)*rm*1000.);
+ Jo_array(j),be_array(j)*PRM*1000.);
eno_array(i,medium) = -eno_array(i,medium);
eno_atbin_array(i,medium) = i;
]
@@ -2682,9 +2681,9 @@ subroutine mscati;
implicit none;
$REAL ededx,ei,eil,eip1,eip1l,si,sip1,eke,elke,aux,ecutmn,tstbm,tstbmn;
-$REAL p2,beta2,rm,dedx0,ekef,elkef,estepx,ektmp,elktmp,chi_a2;
+$REAL p2,beta2,dedx0,ekef,elkef,estepx,ektmp,elktmp,chi_a2;
$INTEGER
- i,leil,leip1l,medium,neke,lelke,lelkef,lelktmp;
+ i,leil,leip1l,neke,lelke,lelkef,lelktmp;
$LOGICAL ise_monoton, isp_monoton;
$declare_write_buffer;
@@ -2719,7 +2718,7 @@ IF( skindepth_for_bca <= 1e-4 ) [
$egs_info(*,' minimum ECUT found: ',ecutmn);
tstbmn = 1e30;
DO medium = 1,nmed [
- tstbm = (ecutmn-0.5110034)*(ecutmn+0.5110034)/ecutmn**2;
+ tstbm = (ecutmn-prm)*(ecutmn+prm)/ecutmn**2;
tstbm = blcc(medium)*tstbm*(ecutmn/xcc(medium))**2;
aux = Log(tstbm);
IF( aux > 300 ) $egs_info(*,'aux > 300 ? ',aux);
@@ -2799,7 +2798,6 @@ $egs_info(*,' Initializing tmxs for estepe = ',estepe,' and ximax = ',ximax);
$egs_info(*,' ');
"Determine upper limit in step size for multiple scattering
-rm = 0.5110034;
DO medium = 1,nmed
[
" Calculate range array first "
@@ -3242,7 +3240,7 @@ subroutine init_spin;
implicit none;
$declare_max_medium;
-COMIN/Spin-Data,ELECIN,MEDIA,BREMPR,THRESH,EGS-IO/;
+COMIN/Spin-Data,ELECIN,MEDIA,BREMPR,THRESH,USEFUL,EGS-IO/;
"BREMPR is needed for the elemental composition"
$REAL eta_array(0:1,$0-MAXE_SPI1),
@@ -3260,7 +3258,7 @@ real*4 dum1,dum2,dum3,aux_o; " These must be 32 bit floats!!!!!"
real*4 fmax_array(0:$MAXQ_SPIN);
integer*2 i2_array(512),ii2;
-$INTEGER medium,iq,i,j,k,i_ele,iii,iZ,iiZ,n_ener,n_q,n_point,je,neke,
+$INTEGER iq,i,j,k,i_ele,iii,iZ,iiZ,n_ener,n_q,n_point,je,neke,
ndata,leil,length,ii4,irec;
character spin_file*256;
@@ -3325,7 +3323,7 @@ dbeta2 = (b2spin_max - b2spin_min)/n_ener;
beta2 = b2spin_min; earray(n_ener+1) = espin_max;
DO i=n_ener+2,2*n_ener+1 [
beta2 = beta2 + dbeta2;
- IF( beta2 < 0.999 ) [ earray(i) = 511.0034*(1/sqrt(1-beta2)-1); ]
+ IF( beta2 < 0.999 ) [ earray(i) = prm*1000.0*(1/sqrt(1-beta2)-1); ]
ELSE [ earray(i) = 50585.1; ]
$FOOL-INTEL-OPTIMIZER(25) i,earray(i);
]
@@ -3384,7 +3382,7 @@ DO medium = 1,NMED [
tmp_4 = aux_o; call egs_swap_4(c_4); aux_o = tmp_4;
]
eta_array(iq,i)=eta_array(iq,i)+tmp*Log(Z23*aux_o);
- tau = earray(i)/511.0034; "energy in the file is in keV"
+ tau = earray(i)/prm*0.001; "energy in the file is in keV"
beta2 = tau*(tau+2)/(tau+1)**2;
eta = Z23/(fine*TF_constant)**2*aux_o/4/tau/(tau+2);
c_array(iq,i)=c_array(iq,i)+
@@ -3427,7 +3425,7 @@ DO medium = 1,NMED [
" Process eta_array, c_array and g_array to their final form "
$FOOL-INTEL-OPTIMIZER(25) 'Spin corrections as read in from file';
DO i=0,$MAXE_SPI1 [
- tau = earray(i)/511.0034; beta2 = tau*(tau+2)/(tau+1)**2;
+ tau = earray(i)/prm*0.001; beta2 = tau*(tau+2)/(tau+1)**2;
DO iq=0,1 [
aux_o = Exp(eta_array(iq,i)/sum_Z2)/(fine*TF_constant)**2;
eta_array(iq,i) = 0.26112447*aux_o*blcc(medium)/xcc(medium);
@@ -3452,7 +3450,7 @@ DO medium = 1,NMED [
aae = (eil-espml)*dleneri; je = aae; aae = aae - je;
]
ELSE [
- tau = e/0.5110034; beta2 = tau*(tau+2)/(tau+1)**2;
+ tau = e/prm; beta2 = tau*(tau+2)/(tau+1)**2;
aae = (beta2 - b2spin_min)*dbeta2i;
je = aae; aae = aae - je; je = je + $MAXE_SPIN + 1;
]
@@ -3471,7 +3469,7 @@ DO medium = 1,NMED [
aae = (eil-espml)*dleneri; je = aae; aae = aae - je;
]
ELSE [
- tau = e/0.5110034; beta2 = tau*(tau+2)/(tau+1)**2;
+ tau = e/prm; beta2 = tau*(tau+2)/(tau+1)**2;
aae = (beta2 - b2spin_min)*dbeta2i;
je = aae; aae = aae - je; je = je + $MAXE_SPIN + 1;
]
@@ -3583,11 +3581,11 @@ DO medium = 1,NMED [
"Now substract scattering power that is already taken into account in"
"discrete Moller/Bhabha events"
- tauc = te(medium)/0.5110034;
+ tauc = te(medium)/prm;
si1e = 1;
DO i=1,neke-1 [
eil = (i+1 - eke0(medium))/eke1(medium);
- e = Exp(eil); leil=i+1; tau=e/0.5110034;
+ e = Exp(eil); leil=i+1; tau=e/prm;
IF( tau > 2*tauc ) [
$EVALUATE sig USING esig(eil);
$EVALUATE dedx USING ededx(eil);
@@ -3646,7 +3644,7 @@ subroutine init_spin_old;
implicit none;
$declare_max_medium;
-COMIN/Spin-Data,ELECIN,MEDIA,BREMPR,THRESH,EGS-IO/;
+COMIN/Spin-Data,ELECIN,MEDIA,BREMPR,THRESH,USEFUL,EGS-IO/;
"BREMPR is needed for the elemental composition"
$REAL eta_array(0:1,$0-MAXE_SPI1),
@@ -3659,7 +3657,7 @@ $REAL eta_array(0:1,$0-MAXE_SPI1),
af($0-MAXE_SPI1),bf($0-MAXE_SPI1),cf($0-MAXE_SPI1),
df($0-MAXE_SPI1),spline;
-$INTEGER medium,iq,i,j,k,i_ele,iii,iZ,iiZ,n_ener,n_q,n_point,je,neke,
+$INTEGER iq,i,j,k,i_ele,iii,iZ,iiZ,n_ener,n_q,n_point,je,neke,
ndata,leil,length,want_spin_unit,spin_unit,egs_get_unit;
character spin_file*256;
@@ -3729,7 +3727,7 @@ DO medium = 1,NMED [
read(spin_unit,'(a,g14.6)') string,earray(i);
read(spin_unit,*) dum1,dum2,dum3,aux_o;
eta_array(iq,i)=eta_array(iq,i)+tmp*Log(Z23*aux_o);
- tau = earray(i)/511.0034; "energy in the file is in keV"
+ tau = earray(i)/prm*0.001; "energy in the file is in keV"
beta2 = tau*(tau+2)/(tau+1)**2;
eta = Z23/(fine*TF_constant)**2*aux_o/4/tau/(tau+2);
c_array(iq,i)=c_array(iq,i)+
@@ -3766,7 +3764,7 @@ DO medium = 1,NMED [
" Process eta_array, c_array and g_array to their final form "
DO i=0,$MAXE_SPI1 [
- tau = earray(i)/511.0034; beta2 = tau*(tau+2)/(tau+1)**2;
+ tau = earray(i)/prm*0.001; beta2 = tau*(tau+2)/(tau+1)**2;
DO iq=0,1 [
aux_o = Exp(eta_array(iq,i)/sum_Z2)/(fine*TF_constant)**2;
eta_array(iq,i) = 0.26112447*aux_o*blcc(medium)/xcc(medium);
@@ -3796,7 +3794,7 @@ DO medium = 1,NMED [
aae = (eil-espml)*dleneri; je = aae; aae = aae - je;
]
ELSE [
- tau = e/0.5110034; beta2 = tau*(tau+2)/(tau+1)**2;
+ tau = e/prm; beta2 = tau*(tau+2)/(tau+1)**2;
aae = (beta2 - b2spin_min)*dbeta2i;
je = aae; aae = aae - je; je = je + $MAXE_SPIN + 1;
]
@@ -3814,7 +3812,7 @@ DO medium = 1,NMED [
aae = (eil-espml)*dleneri; je = aae; aae = aae - je;
]
ELSE [
- tau = e/0.5110034; beta2 = tau*(tau+2)/(tau+1)**2;
+ tau = e/prm; beta2 = tau*(tau+2)/(tau+1)**2;
aae = (beta2 - b2spin_min)*dbeta2i;
je = aae; aae = aae - je; je = je + $MAXE_SPIN + 1;
]
@@ -3907,11 +3905,11 @@ DO medium = 1,NMED [
"Now substract scattering power that is already taken into account in"
"discrete Moller/Bhabha events"
- tauc = te(medium)/0.5110034;
+ tauc = te(medium)/prm;
si1e = 1;
DO i=1,neke-1 [
eil = (i+1 - eke0(medium))/eke1(medium);
- e = Exp(eil); leil=i+1; tau=e/0.5110034;
+ e = Exp(eil); leil=i+1; tau=e/prm;
IF( tau > 2*tauc ) [
$EVALUATE sig USING esig(eil);
$EVALUATE dedx USING ededx(eil);
@@ -3972,7 +3970,7 @@ end;
" "
subroutine msdist_pII
(
- e0,eloss,tustep,rhof,medium,qel,spin_effects,u0,v0,w0,x0,y0,z0, "Inputs
+ e0,eloss,tustep,rhof,med,qel,spin_effects,u0,v0,w0,x0,y0,z0, "Inputs
us,vs,ws,xf,yf,zf,ustep "Outputs
);
@@ -3993,7 +3991,7 @@ $REAL
z0 "initial z-position
;
$INTEGER
- medium,"medium number
+ med,"medium number
qel "=0 for e-, =1 for e+, needed for spin effects
;
$LOGICAL
@@ -4079,7 +4077,9 @@ $INTEGER
lelke
;
$declare_max_medium;
-;COMIN/ELECIN,THRESH,UPHIOT,RANDOM,CH-Steps/;
+;COMIN/ELECIN,THRESH,UPHIOT,RANDOM,CH-Steps,USEFUL/;
+
+medium = med;
count_pII_steps = count_pII_steps + 1;
blccc = blcc(medium);
@@ -4087,7 +4087,7 @@ xcccc = xcc(medium);
"Commonly used factors
e = e0 - 0.5*eloss;
-tau = e/0.5110034;
+tau = e/prm;
tau2 = tau*tau;
epsilon = eloss/e0;
epsilonp= eloss/e;
@@ -8133,6 +8133,7 @@ real*8 sum0,a,b,x1,x2,pow_x1,pow_x2,dle,e,xmax,
anorm,anorm1,anorm2,w,dw,xold,t,aux;
$INTEGER i,j,k,ibin;
+;COMIN/USEFUL/;
write(*,'(a$)') ' preparing data for Rayleigh sampling ... ';
@@ -8162,7 +8163,7 @@ DO i=1,ndat-1 [
**************************************************************/
dle = log(emax/emin)/(ne-1); i = 1;
DO j=1,ne [
- e = emin*exp(dle*(j-1)); xmax = 20.607544d0*2*e/0.5110034d0;
+ e = emin*exp(dle*(j-1)); xmax = 20.607544d0*2*e/prm;
DO k=i,ndat-1 [
IF( xmax >= x(k) & xmax < x(k+1) ) EXIT;
]
@@ -8237,11 +8238,12 @@ $REAL function egs_KN_sigma0(e);
"=========================================================================="
implicit none;
$REAL e;
-$REAL rm,con,ko,c1,c2,c3,eps1,eps2;
-data rm/0.5110034/,con/0.1274783851/;
-ko = e/rm;
+$REAL con,ko,c1,c2,c3,eps1,eps2;
+data con/0.1274783851/;
+;COMIN/USEFUL/;
+ko = e/prm;
IF( ko < 0.01 ) [
- egs_KN_sigma0 = 8.*con/3.*(1-ko*(2-ko*(5.2-13.3*ko)))/rm;
+ egs_KN_sigma0 = 8.*con/3.*(1-ko*(2-ko*(5.2-13.3*ko)))/prm;
return;
]
c1 = 1./(ko*ko); c2 = 1. - 2*(1+ko)*c1; c3 = (1+2*ko)*c1;
@@ -8255,9 +8257,10 @@ $REAL function egs_KN_sigma1(e);
"=========================================================================="
implicit none;
$REAL e;
-$REAL rm,con,ko,c1,c2,c3,eps1,eps2;
-data rm/0.5110034/,con/0.1274783851/;
-ko = e/rm; c1 = 1./(ko*ko); c2 = 1. - 2*(1+ko)*c1; c3 = (1+2*ko)*c1;
+$REAL con,ko,c1,c2,c3,eps1,eps2;
+data con/0.1274783851/;
+;COMIN/USEFUL/;
+ko = e/prm; c1 = 1./(ko*ko); c2 = 1. - 2*(1+ko)*c1; c3 = (1+2*ko)*c1;
eps2 = 1; eps1 = 1./(1+2*ko);
egs_KN_sigma1 = c1*(1./eps1-1./eps2);
egs_KN_sigma1 = egs_KN_sigma1 + log(eps2/eps1)*(c2 - c1) - c2*(eps2-eps1);
@@ -8272,13 +8275,14 @@ subroutine egsi_get_data(flag,iunit,n,ne,zsorted,pz_sorted,ge1,ge0,data);
"=========================================================================="
implicit none;
COMIN/EGS-IO/;
-$REAL rm,eth;
+$REAL eth;
$INTEGER flag,iunit,n,ne;
$REAL ge1,ge0,zsorted(*),pz_sorted(*),data(*);
$REAL etmp($MXINPUT),ftmp($MXINPUT);
$REAL gle,sig,p,e;
$INTEGER i,j,k,kk,iz,iz_old,ndat,iiz;
-data rm/0.5110034/;
+
+;COMIN/USEFUL/;
"Ali:photonuc. The whole routine is commented out and re-written
"to accommodate reading photonuclear cross sections. A copy of the
diff --git a/HEN_HOUSE/src/emf_macros.mortran b/HEN_HOUSE/src/emf_macros.mortran
index 1fe18e7fa..948905b2c 100644
--- a/HEN_HOUSE/src/emf_macros.mortran
+++ b/HEN_HOUSE/src/emf_macros.mortran
@@ -27,6 +27,7 @@
" Amir Keyvanloo "
" Ernesto Mainegra-Hing "
" Blake Walters "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -386,7 +387,7 @@ ELSE[
REPLACE {$GET-EM-FIELD(#,#,#,#,#);} WITH {;
{P1}X0=0.0;{P1}Y0=0.0;{P1}Z0=0.0;
-{P2}X0=0.0;{P2}Y0=0.0;{P2}Z0=3.0*3.0/0.511;
+{P2}X0=0.0;{P2}Y0=0.0;{P2}Z0=3.0*3.0/PRM;
}
;
@@ -399,16 +400,16 @@ REPLACE {$GET-POTENTIAL(#,#,#,#);} WITH {;
{P1}=-GAMMA0*BETA20*RM*(EX0*{P2}+EY0*{P3}+EZ0*{P4});
}
-"3.0*3.0/0.511 explanation "
+"3.0*3.0/PRM explanation "
"nT = 3.0 (this example) = Magnetic field strength in Tesla V.s/m**2 "
"3.0 = Speed of light, in m/s *10**8 "
-"0.511 = rest mass energy of the electron in units eV*10**6 "
+"PRM = rest mass energy of the electron in units eV*10**6 "
" "
-" (nT)[V.s/m**2](3*10**8[m/s])*(e)/0.511*10**6[e.V] "
+" (nT)[V.s/m**2](3*10**8[m/s])*(e)/PRM*10**6[e.V] "
" = "
-" (nT)*3*100/0.511[1/m] "
+" (nT)*3*100/PRM[1/m] "
" = "
-" (nT)*3.0/0.511[1/cm] "
+" (nT)*3.0/PRM[1/cm] "
"**************************************************************************"
"End of example is for a constant B-field of 3 Tesla parallel to the z-axis"
diff --git a/HEN_HOUSE/src/get_media_inputs.mortran b/HEN_HOUSE/src/get_media_inputs.mortran
index 6b8e55a8f..2728dc747 100644
--- a/HEN_HOUSE/src/get_media_inputs.mortran
+++ b/HEN_HOUSE/src/get_media_inputs.mortran
@@ -23,7 +23,7 @@
" "
" Authors: Blake Walters, 2013 "
" "
-" Contributors: "
+" Contributors: Reid Townson "
" "
"#############################################################################"
@@ -663,7 +663,7 @@ $INTEGER ounit;
$declare_max_medium;
$COMIN-GET-TRANSPORTP;
-COMIN/MEDINP,ELECIN,THRESH,ELEMTB/;
+COMIN/MEDINP,ELECIN,THRESH,ELEMTB,USEFUL/;
$INTEGER ival,ival_media,ival_medfile,i,j,k,ival_ae,ival_ue,ival_ap,ival_up,
ival_rho,ival_elements,ival_rhoz,ival_iunrst,ival_iaprim,ival_gasp,
@@ -692,7 +692,7 @@ $REAL4 ZTBL;
$REAL EKE,ELKE,TMXSO,DEDXE,DEDXP,EFRACT,SIGE,SIGP,BREME,BREMP,ETAB(16),
EIE,PLOTE($MAXPOINTS),PLOTEM($MAXPOINTS),PLOTEEN($MAXPOINTS),
PLOTEMP($MAXPOINTS), PLOTEMS($MAXPOINTS);
-$INTEGER IPLOTE,IFLAG1,IFLAG2,MEDIUM,LELKE;
+$INTEGER IPLOTE,IFLAG1,IFLAG2,LELKE;
CHARACTER*60 GRAPHTITLE,XAXIS,YAXISPcom,YAXISPmfp,YAXISE,YAXISEmfp,
SUBTITLE,SERIES;
DATA ETAB/1.,1.25,1.5,1.75,2.,2.5,3.,3.5,4.,4.5,5.,5.5,6.,7.,8.,9./;
@@ -836,7 +836,7 @@ nvalue(ival) = 1; "1 input"
type(ival) = 1;
value_min(ival) = 0;
value_max(ival) = 999.;
-default(ival) = 50.511;
+default(ival) = 50 + prm;
ival = ival + 1;
ival_up = ival;
@@ -867,8 +867,9 @@ IF(error_flags(ival_ue)=0)[
ue_tmp=value(ival_ue,1);
]
ELSE[
- $WRITE_MEDERR(' Warning: UE for media not supplied. Will use 50.511 MeV');
- ue_tmp=50.511;
+ $WRITE_MEDERR(' Warning: UE for media not supplied. Will use
+ 50.5109989461 MeV');
+ ue_tmp=50 + prm;
]
IF(error_flags(ival_up)=0)[
up_tmp=value(ival_up,1);
@@ -1648,14 +1649,18 @@ DO i=1,nmed[
DO j=1,8[
DO k=1,16[
EKE=ETAB(k)*10.**(j-4);
- IF(EKE <= AE(1)-0.511) [
- IF(IFLAG1 = 0) [ IFLAG1=1; EKE=AE(1)-0.511; ] ELSE [ EKE=0.0; ]
+ IF(EKE <= AE(1)-PRM) [
+ IF(IFLAG1 = 0) [
+ IFLAG1=1; EKE=AE(1)-PRM;
+ ] ELSE [ EKE=0.0; ]
]
- IF(EKE > UE(1)-0.511) [
- IF(IFLAG2 = 0) [ IFLAG2=1; EKE=UE(1)-0.511; ] ELSE [ EKE=1.E30; ]
+ IF(EKE > UE(1)-PRM) [
+ IF(IFLAG2 = 0) [
+ IFLAG2=1; EKE=UE(1)-PRM;
+ ] ELSE [ EKE=1.E30; ]
]
- EIE=EKE+0.511; TMXSO=0.0; DEDXE=0.0;
+ EIE=EKE+PRM; TMXSO=0.0; DEDXE=0.0;
DEDXP=0.0; EFRACT=0.0;
IF(EIE >= AE(1)-0.0001 & EIE <= UE(1)+0.001) [
ELKE=LOG(EKE);
diff --git a/HEN_HOUSE/src/nrcaux.mortran b/HEN_HOUSE/src/nrcaux.mortran
index ec3849b5f..4bf7703e8 100644
--- a/HEN_HOUSE/src/nrcaux.mortran
+++ b/HEN_HOUSE/src/nrcaux.mortran
@@ -25,6 +25,7 @@
" Dave Rogers, 2000 "
" "
" Contributors: Ernesto Mainegra-Hing "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -156,7 +157,7 @@ IF((IWATCH = 4) & (IARG >= 0) & (IARG ~= 5)) [ "GRAPHICS OUTPUT"
IF(IARG = 5 | IARG < 0) RETURN;
IF(IWATCH = 4) RETURN; "NONE OF THE REST NEEDED FOR GRAPHICS OUTPUT"
-KE=E(NP);IF(IQ(NP).NE.0)[KE=E(NP)-0.511;]
+KE=E(NP);IF(IQ(NP).NE.0)[KE=E(NP)-PRM;]
IF(IARG = 0 & IWATCH = 2)[
$CNTOUT(NP);(T11,'STEP ABOUT TO OCCUR', T36,':');
diff --git a/HEN_HOUSE/src/pegs4_macros.mortran b/HEN_HOUSE/src/pegs4_macros.mortran
index 0ee72fcb7..8ed1e39b7 100644
--- a/HEN_HOUSE/src/pegs4_macros.mortran
+++ b/HEN_HOUSE/src/pegs4_macros.mortran
@@ -23,7 +23,7 @@
" "
" Authors: Blake Walters, 2013 "
" "
-" Contributors: "
+" Contributors: Reid Townson "
" "
"#############################################################################"
" "
@@ -400,7 +400,6 @@ REPLACE {$INIT-PEGS4-VARIABLES;} WITH {
" AN AVOGADRO'S NUMBER "
PIP=3.1415926536;
C=2.997925E+10;
-RME=9.1091E-28;
HBAR=1.05450E-27;
ECGS=4.80298E-10;
EMKS=1.60210E-19;
@@ -415,9 +414,10 @@ RADDEG=180./PIP;
FSC = ECGS**2/(HBAR*C);
" 1.E+7 IS THE NUMBER OF ERGS PER JOULE "
ERGMEV = (1.E+6)*(EMKS*1.E+7);
+RME = PRM/C**2*ERGMEV;
+"We are using RMP instead of PRM because single precision is required"
+RMP = PRM;
R0 = (ECGS**2)/(RME*C**2);
-RMP = RME*C**2/ERGMEV;
-RMPT2 = RMP*2.0;
RMPSQ = RMP*RMP;
A22P9 = RADDEG*SQRT(4.*PIP*AN)*ECGS**2/ERGMEV;
A6680 = 4.0*PIP*AN*(HBAR/(RME*C))**2*(0.885**2/(1.167*1.13));
diff --git a/HEN_HOUSE/src/rad_compton1.mortran b/HEN_HOUSE/src/rad_compton1.mortran
index 259b0b98e..1d8e5d8e4 100644
--- a/HEN_HOUSE/src/rad_compton1.mortran
+++ b/HEN_HOUSE/src/rad_compton1.mortran
@@ -23,7 +23,7 @@
" "
" Author: Iwan Kawrakow, 2005 "
" "
-" Contributors: "
+" Contributors: Reid Townson "
" "
"#############################################################################"
" "
@@ -229,7 +229,7 @@ DO i=0,ne [
close(radc_unit);
$egs_info(*,' OK');
-radc_le1 = -log(radc_emin*0.5110034)*radc_dlei;
+radc_le1 = -log(radc_emin*prm)*radc_dlei;
"Change the cross section and branching ratios"
$egs_info('(a,$)',' Modifying cross sections and branching ratios ...');
@@ -240,7 +240,7 @@ DO medium=1,nmed [
$EVALUATE gmfp USING gmfp(gle);
$EVALUATE gbr1 USING gbr1(gle);
$EVALUATE gbr2 USING gbr2(gle);
- IF( exp(gle) > radc_emin*0.5110034 ) [
+ IF( exp(gle) > radc_emin*prm ) [
acheck = radc_dlei*gle + radc_le1;
icheck = acheck; acheck = acheck - icheck; acheck1 = 1 - acheck;
frad = (radc_sigs(icheck)+radc_sigd(icheck))*acheck1 +
diff --git a/HEN_HOUSE/user_codes/dosrznrc/dosrznrc.mortran b/HEN_HOUSE/user_codes/dosrznrc/dosrznrc.mortran
index dc9d8818e..d9698d560 100644
--- a/HEN_HOUSE/user_codes/dosrznrc/dosrznrc.mortran
+++ b/HEN_HOUSE/user_codes/dosrznrc/dosrznrc.mortran
@@ -4946,8 +4946,8 @@ IF(IFULL=2)["for pulse height distribution"
"Set up energy bins for peak efficiency calculation"
"This only makes sense if beam is monoenergetic"
- DFEN(1,2)=EPHTOP-DELTAE;DFEN(2,2)=DFEN(1,2)-0.511;
- DFEN(3,2)=DFEN(1,2)-1.022;DFEN(4,2)=0.511-DELTAE;
+ DFEN(1,2)=EPHTOP-DELTAE;DFEN(2,2)=DFEN(1,2)-PRM;
+ DFEN(3,2)=DFEN(1,2)-1.022;DFEN(4,2)=PRM-DELTAE;
"I.E. WE SET LOWER ENERGIES FOR FULL ENERGY, SINGLE ESCAPE, DOUBLE"
"ESCAPE AND THE 511 KEV LINE"
DO IPK = 1,4 [
@@ -5010,10 +5010,10 @@ VALUES_SOUGHT(IVAL)='ESAVEIN';
NVALUE(IVAL)=1;
TYPE(IVAL)=1;
VALUE_MIN(IVAL)=0.;
-VALUE_MAX(IVAL)=EIN+0.511; "EIN is max k.e. of particles set in srcrz for"
- "ISOURC 21,22 and ine ensrc for the others"
-DEFAULT(IVAL) =EIN+0.511; "This means that if a too high number is input,"
- "it is reduced to the max needed."
+VALUE_MAX(IVAL)=EIN+PRM; "EIN is max k.e. of particles set in "
+ "srcrz for ISOURC 21,22 and ine ensrc for the others"
+DEFAULT(IVAL) =EIN+PRM; "This means that if a too high number "
+ "is input, it is reduced to the max needed."
"********************************************************************"
DELIMETER='VARIANCE REDUCTION';
diff --git a/HEN_HOUSE/user_codes/dosxyznrc/dosxyznrc.mortran b/HEN_HOUSE/user_codes/dosxyznrc/dosxyznrc.mortran
index c2170eb5a..7a0e774f1 100644
--- a/HEN_HOUSE/user_codes/dosxyznrc/dosxyznrc.mortran
+++ b/HEN_HOUSE/user_codes/dosxyznrc/dosxyznrc.mortran
@@ -34,6 +34,7 @@
" Iwan Kawrakow "
" Frederic Tessier "
" Ernesto Mainegra-Hing "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -1910,7 +1911,7 @@ ELSE[
"set initial energy - used only for enflag=0 (mono-energetic source)"
" reset for other sources elsewhere"
-IF(iqin = 0)[etotin = ein;] ELSE [etotin = ein + 0.511;]
+IF(iqin = 0)[etotin = ein;] ELSE [etotin = ein + prm;]
OUTPUT61 NCASE,IWATCH,TIMMAX,IXXIN,JXXIN,BEAM_SIZE,ISMOOTH,IRESTART,IDAT,
IREJECT,ESAVE_GLOBAL,NRCYCL,IPARALLEL,PARNUM,n_split,ihowfarless,i_phsp_out;
@@ -2656,17 +2657,17 @@ DO ibatch = 1,$NBATCH["Break into $NBATCH batches"
"electrons and photons in the total incident E"
IF(iqin =1)[ esrc = esrc + weight*1.0220068;
"we include 2 0.511 MeV photons from e+ - e- annihilation"]
- IF(iqin ~= 0)[ etotin = esrc_sp + 0.5110034;"particle total E"]
+ IF(iqin ~= 0)[ etotin = esrc_sp + prm;"particle total E"]
ELSE[etotin = esrc_sp;]
]
"-----------------------------------------------------------------"
ELSEIF(enflag = 2 | enflag = 3)["phase-space input or BEAM sim."
"-----------------------------------------------------------------"
etotin = einsrc; "einsrc is the total energy of the particle"
- IF(iqin = -1) [esrc = esrc + weight*(etotin - 0.5110034);
+ IF(iqin = -1) [esrc = esrc + weight*(etotin - prm);
"remove e- rest mass"]
ELSEIF(iqin = 0)[esrc = esrc + weight*etotin; ]
- ELSEIF(iqin = 1)[esrc = esrc + weight*(etotin + 0.5110034);
+ ELSEIF(iqin = 1)[esrc = esrc + weight*(etotin + prm);
"add 0.511 for the positron but was - till June 97"
]
]
@@ -2686,7 +2687,7 @@ DO ibatch = 1,$NBATCH["Break into $NBATCH batches"
IF(xin*xin+yin*yin<0.31831)[
IF(iqin ~= 0)[
planarfe=planarfe+weight;
- planarefe=planarefe+(etotin-0.5110034)*weight;
+ planarefe=planarefe+(etotin-prm)*weight;
]
ELSE["photons"
planarfp=planarfp+weight;
@@ -4630,12 +4631,11 @@ end;
subroutine modify_tmxs(mindel);
implicit none;
-;COMIN/BOUNDS,ELECIN,MEDIA,MISC,RANDOM,ET-Control/;
-$REAL rm,eil,ei,p2,beta2,chi_a2,elke,eke,dedx0,estepx,si,sip1,tmxs_orig,
+;COMIN/BOUNDS,ELECIN,MEDIA,MISC,RANDOM,ET-Control,USEFUL/;
+$REAL eil,ei,p2,beta2,chi_a2,elke,eke,dedx0,estepx,si,sip1,tmxs_orig,
mindel;
-$INTEGER medium,i,neke,leil,lelke;
+$INTEGER i,neke,leil,lelke;
-rm = 0.5110034;
DO medium = 1,nmed [
neke = meke(medium);
diff --git a/HEN_HOUSE/user_codes/egs_fac/egs_range_rejection.cpp b/HEN_HOUSE/user_codes/egs_fac/egs_range_rejection.cpp
index 97fd6ef12..1111bc935 100644
--- a/HEN_HOUSE/user_codes/egs_fac/egs_range_rejection.cpp
+++ b/HEN_HOUSE/user_codes/egs_fac/egs_range_rejection.cpp
@@ -23,7 +23,7 @@
#
# Author: Iwan Kawrakow, 2008
#
-# Contributors:
+# Contributors: Reid Townson, 2017
#
###############################################################################
*/
@@ -67,7 +67,7 @@ EGS_RangeRejection* EGS_RangeRejection::getRangeRejection(EGS_Input *input,
if( !result->cgeom ) egsWarning("\n\n********** no geometry named"
" %s exists => using region-by-region rejection only\n");
}
- if( result->Esave <= 0.511 && result->type == RangeDiscard ) {
+ if( result->Esave <= the_useful->prm && result->type == RangeDiscard ) {
egsWarning("\n\n********* rr_flag = 1 but Esave = 0 =>"
" not using range rejection\n\n");
delete result; result = 0;
diff --git a/HEN_HOUSE/user_codes/examin/examin.mortran b/HEN_HOUSE/user_codes/examin/examin.mortran
index 9fa9b176f..4837165a6 100644
--- a/HEN_HOUSE/user_codes/examin/examin.mortran
+++ b/HEN_HOUSE/user_codes/examin/examin.mortran
@@ -30,6 +30,7 @@
" Iwan Kawrakow "
" Ernesto Mainegra-Hing "
" Frederic Tessier "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -318,7 +319,7 @@
call egs_fdate(nunit);
write(nunit,:LAB1a:);
:LAB1a:format(/' ',79('='));
- WRITE(NUNIT,:LABL2:) AE(1)-0.511,UE(1)-0.511,AP(1),UP(1);
+ WRITE(NUNIT,:LABL2:) AE(1)-PRM,UE(1)-PRM,AP(1),UP(1);
:LABL2:FORMAT(/' Electron kinetic energy range:',T35,2F12.3,' MeV'/
' Photon energy range: ',T35,2F12.3,' MeV');
WRITE(NUNIT,:LABL3:) RLC(MEDIUM),RHO(MEDIUM);
@@ -489,7 +490,8 @@
call egs_fdate(nunit);
write(nunit,:lab2a:);
:lab2a:format(/' ',79('='));
- WRITE(NUNIT,:LABL13:) AE(1)-0.511,UE(1)-0.511,AP(1),UP(1);
+ WRITE(NUNIT,:LABL13:)
+ AE(1)-PRM,UE(1)-PRM,AP(1),UP(1);
:LABL13:FORMAT(/' Electron kinetic energy range:',T35,2F12.3,' MeV'/
' Photon energy range:' ,T35,2F12.3,' MeV');
WRITE(NUNIT,:LABL14:) RLC(MEDIUM),RHO(MEDIUM),200.*TEFF0(1);
@@ -552,23 +554,27 @@
OUTPUT;(' Electron kinetic energy( MeV): ',$);
READ(5,'(F12.0)',ERR=:label2:,END= :ELECTRON:) EKE;
- IF (EKE > UE(1)-0.511003 | EKE < AE(1)-0.511003) [
- OUTPUT AE(1)-0.511,UE(1)-0.511;
+ IF (EKE > UE(1)-PRM | EKE < AE(1)-PRM) [
+ OUTPUT AE(1)-PRM,UE(1)-PRM;
(T20,' must be in range',2(1Pe10.3));
]
IF (EKE = 0.0) [STOP; ]
]
ELSE [ EKE=ETAB(I2)*10.**(I1-4); ]
- IF(EKE <= AE(1)-0.511) [
- IF(IFLAG1 = 0) [ IFLAG1=1; EKE=AE(1)-0.511; ] ELSE [ EKE=0.0; ]
+ IF(EKE <= AE(1)-PRM) [
+ IF(IFLAG1 = 0) [
+ IFLAG1=1; EKE=AE(1)-PRM;
+ ] ELSE [ EKE=0.0; ]
]
- IF(EKE > UE(1)-0.511) [
- IF(IFLAG2 = 0) [ IFLAG2=1; EKE=UE(1)-0.511; ] ELSE [ EKE=1.E30; ]
+ IF(EKE > UE(1)-PRM) [
+ IF(IFLAG2 = 0) [
+ IFLAG2=1; EKE=UE(1)-PRM;
+ ] ELSE [ EKE=1.E30; ]
]
- EIE=EKE+0.511; TMXS=0.0; DEDXE=0.0;
+ EIE=EKE+PRM; TMXS=0.0; DEDXE=0.0;
DEDXP=0.0; EFRACT=0.0;
IF(EIE >= AE(1)-0.0001 & EIE <= UE(1)+0.001) [
diff --git a/HEN_HOUSE/user_codes/flurznrc/flurznrc.mortran b/HEN_HOUSE/user_codes/flurznrc/flurznrc.mortran
index 0e165ca70..0f76ecb4c 100644
--- a/HEN_HOUSE/user_codes/flurznrc/flurznrc.mortran
+++ b/HEN_HOUSE/user_codes/flurznrc/flurznrc.mortran
@@ -31,6 +31,7 @@
" Blake Walters "
" Iwan Kawrakow "
" Ernesto Mainegra-Hing "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -2794,7 +2795,7 @@ IF(IFLDIF = 1)[
]
IF(SLOTE > 0.0)["EQUAL BIN WIDTHS"
- /EBINW(1,1),EBINW(1,3)/ = SLOTE - (ECUTMX-0.5110034);
+ /EBINW(1,1),EBINW(1,3)/ = SLOTE - (ECUTMX-PRM);
IF(EBINW(1,1) <= 0.0)[/EBINW(1,1),EBINW(1,3)/=1.E10;]
"above in case large ECUT means electron bins don't exist"
EBINW(1,2) = SLOTE - PCUTMX;
@@ -2805,7 +2806,7 @@ IF(SLOTE > 0.0)["EQUAL BIN WIDTHS"
]
ELSE ["SLOTE=0.0 and SLOTE<0 CASE"
- /EBINW(1,1),EBINW(1,3)/ = BINTOP(1) - (ECUTMX-0.5110034);
+ /EBINW(1,1),EBINW(1,3)/ = BINTOP(1) - (ECUTMX-PRM);
EBINW(1,2) = BINTOP(1) - PCUTMX;
DO IB = 2,MAXIB[DO IQL=1,3[
EBINW(IB,IQL) = BINTOP(IB)-BINTOP(IB-1) ;
@@ -3178,7 +3179,9 @@ IF(IARG = 0)["transport step about to happen"
"IBDOWN is the actual bin being considered, IBTOP to 1 "
IBDOWN = IBTOP - IB +1;
- IF(IBDOWN=1)[BINBOT=ECUT(IRL)-0.5110034;] ELSE[BINBOT=BINTOP(IBDOWN-1);]
+ IF(IBDOWN=1)[
+ BINBOT=ECUT(IRL)-PRM;
+ ] ELSE[BINBOT=BINTOP(IBDOWN-1);]
IF(IB=1)["bin with initial energy" BINFR(IB)=(EKIN-BINBOT)/EDEP;]
ELSE [ BINFR(IB) = (BINTOP(IBDOWN) - BINBOT)/EDEP;]
@@ -3189,7 +3192,7 @@ IF(IARG = 0)["transport step about to happen"
" bin, but it will be reduced to the correct value below "
IF( (BINFR(IB) >= 1.0 - FRSUM) | IBDOWN=1 |
- (BINBOT<=ECUT(IRL)-0.5110034) )[
+ (BINBOT<=ECUT(IRL)-PRM) )[
"this must be the last bin - "
"first condition is true if the final energy in the "
"step is greater than BINBOT, second condition means there"
@@ -3198,7 +3201,7 @@ IF(IARG = 0)["transport step about to happen"
"in first and last cases, we must correct BINFR "
IF((BINFR(IB) >= 1.0-FRSUM) |
- (BINBOT<=ECUT(IRL)-0.5110034 & IBDOWN>1))[
+ (BINBOT<=ECUT(IRL)-PRM & IBDOWN>1))[
"final energy is in this bin"
BINFR(IB) = 1.0 - FRSUM; ]
@@ -3967,14 +3970,14 @@ ELSEIF (SLOTE = -999) ["set up linear/log bins for SPR calcns"
"10% of bins cover top 10% of spectrum"
"rest are equal logarithmic"
LINBIN = int(float($EBIN)/10. + 0.5); "no of bins to cover 10%"
- DELLIN = (EIN - (ECUTIN-0.5110034))*0.10/dble(LINBIN);
+ DELLIN = (EIN - (ECUTIN-PRM))*0.10/dble(LINBIN);
OUTPUT LINBIN,DELLIN;(/' ***SLOTE=-999***', I4,' BINS OF',F10.4,' MEV'/);
MAXIB = $EBIN; BINTOP(MAXIB) = EIN;
DO I =1,linbin[
BINTOP(MAXIB-I) = EIN-(dble(I)*DELLIN);]
"now set up logarithmic bins"
- RATIO=(BINTOP(MAXIB-LINBIN)/(ECUTIN-0.5110034))**(1./dble(MAXIB-LINBIN));
- BINTOP(1) = (ECUTIN - 0.5110034)*RATIO;
+ RATIO=(BINTOP(MAXIB-LINBIN)/(ECUTIN-PRM))**(1./dble(MAXIB-LINBIN));
+ BINTOP(1) = (ECUTIN - PRM)*RATIO;
DO IB = 2,MAXIB-LINBIN[
BINTOP(IB) = BINTOP(IB-1)*RATIO;]
OUTPUT MAXIB,ECUTIN,(BINTOP(IB),IB=1,MAXIB);(/I4,' SPR ENERGY ',
@@ -3986,8 +3989,8 @@ ELSEIF (SLOTE < 0.0) [ "SET UP -SLOTE EQUAL LOGARITHMIC BINS"
' ASKED FOR REDUCED TO MAX OF ', I4,' ****');]
" maximum energy in problem is EIN, and we use electron lowest"
" energy to set the bins equally from ECUTIN to EIN"
- RATIO = (EIN /(ECUTIN-0.5110034))**(1./dble(MAXIB));
- BINTOP(1) = (ECUTIN-0.5110034)*RATIO;
+ RATIO = (EIN /(ECUTIN-PRM))**(1./dble(MAXIB));
+ BINTOP(1) = (ECUTIN-PRM)*RATIO;
DO IB=2,MAXIB-1 [
BINTOP(IB) = RATIO * BINTOP(IB-1);]
BINTOP(MAXIB) = EIN; "TO AVOID ROUND-OFF ERRORS"
@@ -4063,10 +4066,10 @@ VALUES_SOUGHT(IVAL)='ESAVEIN';
NVALUE(IVAL)=1;
TYPE(IVAL)=1;
VALUE_MIN(IVAL)=0.;
-VALUE_MAX(IVAL)=EIN+0.511; "EIN is max k.e. of particles set in srcrz for"
- "ISOURC 21,22 and ine ensrc for the others"
-DEFAULT(IVAL)=EIN+0.511; "defaults to max energy in case too high a value"
- " is input, => still do range rej"
+VALUE_MAX(IVAL)=EIN+PRM; "EIN is max k.e. of particles set in srcrz"
+ "for ISOURC 21,22 and ine ensrc for the others"
+DEFAULT(IVAL)=EIN+PRM; "defaults to max energy in case too high "
+ " a value is input, => still do range rej"
DELIMETER='VARIANCE REDUCTION';
$GET_INPUTS(NUM_BREMSPLIT,NUM_ESAVEIN);
@@ -4872,7 +4875,7 @@ CHARACTER*1 a(3);
character backslash;
LOGICAL NEGVAL;
-;COMIN/GEOM,BOUNDS,IODAT1,IODAT2,PRINTC,PLOTC,SCORE/;
+;COMIN/GEOM,BOUNDS,IODAT1,IODAT2,PRINTC,PLOTC,SCORE,USEFUL/;
"March 9, 2012 DR added BOUNDS to correct first energy in point plot"
@@ -5197,7 +5200,7 @@ DO IQ=1,4 [ "Q=4 is for electrons and positrons"
XCOORD(IB) = (BINTOP(IB) + PCUT(2))/2.
]
ELSE ["a charged particle"
- XCOORD(IB) = (BINTOP(IB) + (ECUT(2)-0.5110034))/2.
+ XCOORD(IB) = (BINTOP(IB) + (ECUT(2)-PRM))/2.
]
]"end of IB=1 block"
ELSE [
@@ -5255,7 +5258,7 @@ DO IQ=1,4 [ "Q=4 is for electrons and positrons"
"We now assume that the cutoffs in region 2 are the same"
"everywhere"
"Had to add comin BOUNDS to PLOTSN to get ECUT/PCUT"
- IF(IQ=2)[BIN_CUT=PCUT(2);] ELSE [BIN_CUT=ECUT(2)-0.5110034;]
+ IF(IQ=2)[BIN_CUT=PCUT(2);] ELSE [BIN_CUT=ECUT(2)-PRM;]
"Convert regional number(IRL) to a character string"
IRLR=IRL; $CONVERT_INT(IRLR)_TO_CHAR(CH_IRL);
diff --git a/HEN_HOUSE/user_codes/g/g.mortran b/HEN_HOUSE/user_codes/g/g.mortran
index d714456af..6986a2319 100644
--- a/HEN_HOUSE/user_codes/g/g.mortran
+++ b/HEN_HOUSE/user_codes/g/g.mortran
@@ -234,7 +234,7 @@ CALL HATCH;
;OUTPUT;('\f Start g value calculation: Version 1.14 '//);
OUTPUT (media(j,1),j=1,24);(/' MEDIUM is: ',24A1/);
-;OUTPUT AE(1)-0.511, AP(1);
+;OUTPUT AE(1)-PRM, AP(1);
(' knock-on electrons can be created and any electron followed down to'
/T40,F8.3,' MeV kinetic energy'/
' brem photons can be created and any photon followed down to ',
diff --git a/HEN_HOUSE/user_codes/sprrznrc/sprrznrc.mortran b/HEN_HOUSE/user_codes/sprrznrc/sprrznrc.mortran
index 44c2f3827..2de22f9e8 100644
--- a/HEN_HOUSE/user_codes/sprrznrc/sprrznrc.mortran
+++ b/HEN_HOUSE/user_codes/sprrznrc/sprrznrc.mortran
@@ -29,6 +29,7 @@
" Aaron Merovitz "
" Iwan Kawrakow "
" Ernesto Mainegra-Hing "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -4014,10 +4015,10 @@ VALUES_SOUGHT(IVAL)='ESAVEIN';
NVALUE(IVAL)=1;
TYPE(IVAL)=1;
VALUE_MIN(IVAL)=0.;
-VALUE_MAX(IVAL)=EIN+0.511; "EIN is max k.e. of particles set in srcrz for"
- "ISOURC 21,22 and ine ensrc for the others"
-DEFAULT(IVAL) =EIN+0.511; "This means that if a too high number is input,"
- "it is reduced to the max needed."
+VALUE_MAX(IVAL)=EIN+PRM; "EIN is max k.e. of particles set in "
+ "srcrz for ISOURC 21,22 and ine ensrc for the others"
+DEFAULT(IVAL) =EIN+PRM; "This means that if a too high number is "
+ "input, it is reduced to the max needed."
DELIMETER='VARIANCE REDUCTION';
$GET_INPUTS(NUM_FORCE,NUM_ESAVEIN);
diff --git a/HEN_HOUSE/user_codes/tutor1/tutor1.mortran b/HEN_HOUSE/user_codes/tutor1/tutor1.mortran
index 65cc819bb..d2cbeb0c1 100644
--- a/HEN_HOUSE/user_codes/tutor1/tutor1.mortran
+++ b/HEN_HOUSE/user_codes/tutor1/tutor1.mortran
@@ -25,6 +25,7 @@
" "
" Contributors: Iwan Kawrakow "
" Frederic Tessier "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -62,13 +63,14 @@ REPLACE {$CALL-HOWNEAR(#);} WITH {
}
-;COMIN/BOUNDS,GEOM,MEDIA,MISC,THRESH/;"Note ; before COMIN"
+;COMIN/BOUNDS,GEOM,MEDIA,MISC,THRESH,USEFUL/;"Note ; before COMIN"
" The above expands into a series of COMMON statements"
" BOUNDS contains ECUT and PCUT"
" GEOM passes info to our HOWFAR routine"
" MEDIA contains the array MEDIA"
" MISC contains MED"
" THRESH contains AE and AP"
+" USEFUL contains PRM"
"---------------------------------------------------------------------"
"STEP 2 PRE-HATCH-CALL-INITIALIZATION "
@@ -103,7 +105,7 @@ PCUT(2)=0.1;" terminate photon histories at 0.1 MeV in the plate"
CALL HATCH;" pick up cross section data for TA"
" data file must be assigned to unit 12"
-;OUTPUT AE(1)-0.511, AP(1);
+;OUTPUT AE(1)-PRM, AP(1);
(/' knock-on electrons can be created and any electron followed down to'
/T40,F8.3,' MeV kinetic energy'/
' brem photons can be created and any photon followed down to ',
@@ -129,7 +131,7 @@ ZBOUND=0.1;" plate is 1 mm thick"
"perpendicular to the slab"
IQIN=-1;" incident charge - electrons"
-EIN=20.511;" 20 MeV kinetic energy"
+EIN=20.5109989461;" 20 MeV kinetic energy"
/XIN,YIN,ZIN/=0.0;" incident at origin"
/UIN,VIN/=0.0;WIN=1.0;" moving along Z axis"
IRIN=2;" starts in region 2, could be 1"
@@ -178,7 +180,7 @@ IF(IARG = 3)[
ANGLE=ACOS(W(NP))*180./3.14159;"angle w.r.t. Z axis in degrees"
IF(IQ(NP) = 0)[EKINE=E(NP);]
-ELSE [EKINE=E(NP)-0.511;]"get kinetic energy"
+ELSE [EKINE=E(NP)-0.510998946131;]"get kinetic energy"
OUTPUT EKINE,IQ(NP),ANGLE;(T21,F10.3,T33,I10,T49,F10.1);]
RETURN;END;"END OF AUSGAB"
diff --git a/HEN_HOUSE/user_codes/tutor2/tutor2.mortran b/HEN_HOUSE/user_codes/tutor2/tutor2.mortran
index d3d1a7821..8d2e5ba4f 100644
--- a/HEN_HOUSE/user_codes/tutor2/tutor2.mortran
+++ b/HEN_HOUSE/user_codes/tutor2/tutor2.mortran
@@ -24,6 +24,7 @@
" Author: Dave Roger, 1985 "
" "
" Contributors: Iwan Kawrakow "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -58,7 +59,7 @@ REPLACE {$CALL-HOWNEAR(#);} WITH {
"Define a COMMON for scoring in AUSGAB"
REPLACE {;COMIN/SCORE/;} WITH {;COMMON/SCORE/ESCORE(3);REAL*8 ESCORE;}
-;COMIN/BOUNDS,GEOM,MEDIA,MISC,SCORE, THRESH/;"NOTE ; BEFORE COMIN"
+;COMIN/BOUNDS,GEOM,MEDIA,MISC,SCORE,THRESH,USEFUL/;"NOTE ; BEFORE COMIN"
" The above expands into a series of COMMON statements"
" BOUNDS contains ECUT and PCUT"
" GEOM passes info to our HOWFAR routine"
@@ -67,6 +68,7 @@ REPLACE {;COMIN/SCORE/;} WITH {;COMMON/SCORE/ESCORE(3);REAL*8 ESCORE;}
" SCORE contains the scoring array ESCORE"
" THRESH contains various thresholds"
" THRESH CONTAINS AE AND AP"
+" USEFUL contains PRM"
"---------------------------------------------------------------------"
"STEP 2 PRE-HATCH-CALL-INITIALIZATION "
@@ -100,7 +102,7 @@ PCUT(2)=0.1;" terminate photon histories at 0.1 MeV in the plate"
CALL HATCH;" pick up cross section data for TA"
" data file must be assigned to unit 12"
-;OUTPUT AE(1)-0.511, AP(1);
+;OUTPUT AE(1)-PRM, AP(1);
(/' knock-on electrons can be created and any electron followed down to'
/T40,F8.3,' MeV kinetic energy'/
' brem photons can be created and any photon followed down to ',
@@ -124,7 +126,7 @@ DO I=1,3[ESCORE(I)=0.0;]"zero scoring array before starting"
"perpendicular to the slab"
IQIN=-1;" incident charge - electrons"
-EIN=20.511;" 20 MeV kinetic energy"
+EIN=20 + PRM;" 20 MeV kinetic energy"
/XIN,YIN,ZIN/=0.0;" incident at origin"
/UIN,VIN/=0.0;WIN=1.0;" moving along Z axis"
IRIN=2;" starts in region 2, could be 1"
@@ -141,7 +143,7 @@ DO I=1,NCASE [
"---------------------------------------------------------------------"
"STEP 8 OUTPUT-OF-RESULTS "
"---------------------------------------------------------------------"
-ANORM = 100./((EIN+FLOAT(IQIN)*0.511)*FLOAT(NCASE));
+ANORM = 100./((EIN+FLOAT(IQIN)*PRM)*FLOAT(NCASE));
"normalize to % of total input energy"
TOTAL=0.0; DO I=1,3[TOTAL=TOTAL+ESCORE(I);]
OUTPUT NCASE, (ESCORE(I)*ANORM, I=1,3),TOTAL*ANORM;
diff --git a/HEN_HOUSE/user_codes/tutor3/tutor3.mortran b/HEN_HOUSE/user_codes/tutor3/tutor3.mortran
index e32e11ffb..5fd90e014 100644
--- a/HEN_HOUSE/user_codes/tutor3/tutor3.mortran
+++ b/HEN_HOUSE/user_codes/tutor3/tutor3.mortran
@@ -24,6 +24,7 @@
" Author: Dave Roger, 1985 "
" "
" Contributors: Iwan Kawrakow "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -63,7 +64,7 @@ REPLACE {;COMIN/SCORE/;} WITH {
;COMMON/SCORE/EHIST,EBIN($EBIN);
$REAL EHIST,EBIN;}
-;COMIN/BOUNDS,GEOM,MEDIA,MISC,SCORE, THRESH/;"Note ; before COMIN"
+;COMIN/BOUNDS,GEOM,MEDIA,MISC,SCORE,THRESH,USEFUL/;"Note ; before COMIN"
" The above expands into a series of COMMON statements"
" BOUNDS contains ECUT and PCUT"
" GEOM passes info to our HOWFAR routine"
@@ -71,6 +72,7 @@ REPLACE {;COMIN/SCORE/;} WITH {
" MISC contains MED"
" SCORE contains the scoring variables EHIST and EBIN"
" THRESH contains various thresholds"
+" USEFUL contains PRM"
"---------------------------------------------------------------------"
"STEP 2 PRE-HATCH-CALL-INITIALIZATION "
@@ -103,7 +105,7 @@ PCUT(2)=0.1;" terminate photon histories at 0.1 MeV in the detector"
CALL HATCH;" pick up cross section data for TA"
" data file must be assigned to unit 12"
-;OUTPUT AE(1)-0.511, AP(1);
+;OUTPUT AE(1)-PRM, AP(1);
(/' knock-on electrons can be created and any electron followed down to'
/T40,F8.3,' MeV kinetic energy'/
' brem photons can be created and any photon followed down to ',
diff --git a/HEN_HOUSE/user_codes/tutor4/tutor4.mortran b/HEN_HOUSE/user_codes/tutor4/tutor4.mortran
index e8156f724..ac9d718c8 100644
--- a/HEN_HOUSE/user_codes/tutor4/tutor4.mortran
+++ b/HEN_HOUSE/user_codes/tutor4/tutor4.mortran
@@ -24,6 +24,7 @@
" Author: Dave Roger, 2000 "
" "
" Contributors: Iwan Kawrakow "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -65,7 +66,7 @@ REPLACE {;COMIN/SCORE/;} WITH {
;COMMON/SCORE/ESCORE(3),IWATCH; $INTEGER IWATCH; REAL*8 ESCORE;
}
-;COMIN/BOUNDS,GEOM,MEDIA,MISC,SCORE,STACK,THRESH/;"NOTE ; BEFORE COMIN"
+;COMIN/BOUNDS,GEOM,MEDIA,MISC,SCORE,STACK,THRESH,USEFUL/;"NOTE ; BEFORE COMIN"
" The above expands into a series of COMMON statements"
" BOUNDS contains ECUT and PCUT"
" GEOM passes info to our HOWFAR routine"
@@ -74,6 +75,7 @@ REPLACE {;COMIN/SCORE/;} WITH {
" SCORE contains the scoring array ESCORE and IWATCH"
" STACK to provide LATCHI"
" THRESH contains various thresholds"
+" USEFUL contains PRM"
"---------------------------------------------------------------------"
"STEP 2 PRE-HATCH-CALL-INITIALIZATION "
@@ -106,7 +108,7 @@ PCUT(2)=0.1;" terminate photon histories at 0.1 MeV in the plate"
CALL HATCH;" pick up cross section data for TA"
" data file must be assigned to unit 12"
-;OUTPUT AE(1)-0.511, AP(1);
+;OUTPUT AE(1)-PRM, AP(1);
(/' knock-on electrons can be created and any electron followed down to'
/T40,F8.3,' MeV kinetic energy'/
' brem photons can be created and any photon followed down to ',
@@ -139,7 +141,7 @@ CALL WATCH(-99,IWATCH); "Initializes calls to AUSGAB for WATCH"
"perpendicular to the slab"
IQIN=-1;" incident charge - electrons"
-EIN=20.511;EI=20.0;" 20 MeV kinetic energy"
+EIN=20 + PRM;EI=20.0;" 20 MeV kinetic energy"
/XIN,YIN,ZIN/=0.0;" incident at origin"
/UIN,VIN/=0.0;WIN=1.0;" moving along Z axis"
IRIN=2;" starts in region 2, could be 1"
@@ -163,7 +165,7 @@ DO I=1,NCASE [
"---------------------------------------------------------------------"
"STEP 8 OUTPUT-OF-RESULTS "
"---------------------------------------------------------------------"
-ANORM = 100./((EIN+FLOAT(IQIN)*0.511)*FLOAT(NCASE));
+ANORM = 100./((EIN+FLOAT(IQIN)*PRM)*FLOAT(NCASE));
"normalize to % of total input energy"
TOTAL=0.0; DO I=1,3[TOTAL=TOTAL+ESCORE(I);]
OUTPUT (ESCORE(I)*ANORM, I=1,3),TOTAL*ANORM;
diff --git a/HEN_HOUSE/user_codes/tutor5/tutor5.mortran b/HEN_HOUSE/user_codes/tutor5/tutor5.mortran
index 8dedcec7a..191d7f845 100644
--- a/HEN_HOUSE/user_codes/tutor5/tutor5.mortran
+++ b/HEN_HOUSE/user_codes/tutor5/tutor5.mortran
@@ -24,6 +24,7 @@
" Author: Dave Roger, 1985 "
" "
" Contributors: Iwan Kawrakow "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -65,7 +66,7 @@ REPLACE {;COMIN/SCORE/;} WITH {
REPLACE {$CALL-HOWNEAR(#);} WITH {
;CALL HOWNEAR({P1},X(NP),Y(NP),Z(NP),IRL);}
-;COMIN/BOUNDS,EPCONT,GEOM,MEDIA,MISC,SCORE,STACK,THRESH/;
+;COMIN/BOUNDS,EPCONT,GEOM,MEDIA,MISC,SCORE,STACK,THRESH,USEFUL/;
" the above expands into a series of COMMON statements"
%E "tutor5.mortran"
@@ -98,7 +99,7 @@ IRAYLR(2)=1;" turn on Rayleigh scattering in the slab "
CALL HATCH;" pick up cross section data for H2O"
" data file must be assigned to unit 12"
-;OUTPUT AE(1)-0.511, AP(1);
+;OUTPUT AE(1)-PRM, AP(1);
(/' Knock-on electrons can be created and any electron followed down to'
/T40,F8.3,' MeV kinetic energy'/
' Brem photons can be created and any photon followed down to '
diff --git a/HEN_HOUSE/user_codes/tutor6/tutor6.mortran b/HEN_HOUSE/user_codes/tutor6/tutor6.mortran
index b5bd093ad..0c1a05e3b 100644
--- a/HEN_HOUSE/user_codes/tutor6/tutor6.mortran
+++ b/HEN_HOUSE/user_codes/tutor6/tutor6.mortran
@@ -25,6 +25,7 @@
" "
" Contributors: Blake Walters "
" Iwan Kawrakow "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -133,7 +134,7 @@ call inputs;
;OUTPUT;('\f Start tutor6'//' CALL HATCH to get cross-section data'/);
CALL HATCH;
-;OUTPUT AE(1)-0.511, AP(1);
+;OUTPUT AE(1)-PRM, AP(1);
(/' knock-on electrons can be created and any electron followed down to'
/T40,F8.3,' MeV kinetic energy'/
' brem photons can be created and any photon followed down to ',
diff --git a/HEN_HOUSE/user_codes/tutor7/tutor7.mortran b/HEN_HOUSE/user_codes/tutor7/tutor7.mortran
index 740f530a5..254e58809 100644
--- a/HEN_HOUSE/user_codes/tutor7/tutor7.mortran
+++ b/HEN_HOUSE/user_codes/tutor7/tutor7.mortran
@@ -25,6 +25,7 @@
" "
" Contributors: Dave Rogers "
" Blake Walters "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -156,7 +157,7 @@ $RNG-INITIALIZATION;
;OUTPUT;('\f Start tutor7'//' CALL HATCH to get cross-section data'/);
CALL HATCH;
-;OUTPUT AE(1)-0.511, AP(1);
+;OUTPUT AE(1)-PRM, AP(1);
(/' knock-on electrons can be created and any electron followed down to'
/T40,F8.3,' MeV kinetic energy'/
' brem photons can be created and any photon followed down to ',
diff --git a/HEN_HOUSE/utils/iaea_phsp_macros.mortran b/HEN_HOUSE/utils/iaea_phsp_macros.mortran
index e94b78611..9eb000e56 100644
--- a/HEN_HOUSE/utils/iaea_phsp_macros.mortran
+++ b/HEN_HOUSE/utils/iaea_phsp_macros.mortran
@@ -25,6 +25,7 @@
" "
" Contributors: Iwan Kawrakow "
" Ernesto Mainegra-Hing "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -161,7 +162,7 @@ IF({P2}=0)["do not output if NPASS=1"
IHSTRY_PHSP({P3})={P4}; "reset ihstry"
"JWU: IAEA phsp uses kinetic energy!"
IF( ({P6}=1) | ({P6}=-1) )[
- ESHORT = ESHORT - 0.5110034;
+ ESHORT = ESHORT - 0.5109989461;
]
$IAEA_PARSE_FOR_WRITE({P8});
IF(i_iaea_open_for_write=1)[
@@ -334,7 +335,7 @@ ELSEIF(iaea_n_stat>=0)[
{P5}=iaea_typ_q(iaea_q_index);
"JWU: IAEA phsp uses kinetic energy!"
IF( ({P5}=1) | ({P5}=-1) )[
- ESHORT = ESHORT + 0.5110034;
+ ESHORT = ESHORT + 0.5109989461;
]
{P6}=ESHORT;
IF(iaea_i_zlast=-99)[
diff --git a/HEN_HOUSE/utils/nrcaux.mortran b/HEN_HOUSE/utils/nrcaux.mortran
index 7ada13193..7b09b5600 100644
--- a/HEN_HOUSE/utils/nrcaux.mortran
+++ b/HEN_HOUSE/utils/nrcaux.mortran
@@ -25,6 +25,7 @@
" Dave Rogers, 2000 "
" "
" Contributors: Blake Walters "
+" Reid Townson "
" "
"#############################################################################"
" "
@@ -165,7 +166,7 @@ IF((IWATCH = 4) & (IARG >= 0) & (IARG ~= 5)) [ "GRAPHICS OUTPUT"
IF(IARG = 5 | IARG < 0) RETURN;
IF(IWATCH = 4) RETURN; "NONE OF THE REST NEEDED FOR GRAPHICS OUTPUT"
-KE=E(NP);IF(IQ(NP).NE.0)[KE=E(NP)-0.511;]
+KE=E(NP);IF(IQ(NP).NE.0)[KE=E(NP)-PRM;]
IF(IARG = 0 & IWATCH = 2)[
$CNTOUT(NP);(T11,'STEP ABOUT TO OCCUR', T36,':');