Fix #956: ausgab calls for photonuclear with xcse

Add ausgab calls before and after a photonuclear interaction, in order
to handle cross section enhancement (xcse) correctly. Without these
calls, there is no compensation for the enhanced cross section, hence
the photonuclear cross section is effectively increased by the xcse
factor.

Also add the photonuclear event trapping in Fano calculations, for
consistency. However, the Fano calculation fails when photonuclear
attenuation is turned on anyway, because the photon is discarded
without depositing its energy.

Also update a couple loops that print out the ausgab flags (IAUSFL)
values, to include the photonuclear attenuation option value.

HEN_HOUSE/omega/beamnrc/beamnrc.mortran

@@ -4279,10 +4279,13 @@ IF( use_cs_enhance ) [
  IF( iausfl(21) = 0 ) cse_return(21) = .true.;
  IF( iausfl(24) = 0 ) cse_return(24) = .true.;
  IF( iausfl(25) = 0 ) cse_return(25) = .true.;
+ IF( iausfl(30) = 0 ) cse_return(30) = .true.;
+ IF( iausfl(31) = 0 ) cse_return(31) = .true.;
  IAUSFL(16) = 1; IAUSFL(17) = 1; " before/after PAIR "
  IAUSFL(18) = 1; IAUSFL(19) = 1; " before/after COMPTON "
  IAUSFL(20) = 1; IAUSFL(21) = 1; " before/after PHOTO "
  IAUSFL(24) = 1; IAUSFL(25) = 1; " before/after RAYLEIGH "
+ IAUSFL(30) = 1; IAUSFL(31) = 1; " before/after PHOTONUC "
 ]
 
 "Ali:BCSE"
@@ -5816,16 +5819,16 @@ IF(IBRSPL=2 & USE_REJPLN) ["DBS WITH A REJECTION PLANE"
 IF( use_cs_enhance ) [
  iicm = IR_to_CM(irl);
  IF( cs_enhance(iicm) > 1 ) [
- IF( iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23 ) [
- "Pair/Compton/Photo/Rayleigh about to occur"
+ IF( iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23 | iarg = 29 ) [
+ "Pair/Compton/Photo/Rayleigh/Photonuclear about to occur"
  np = np+1; $CHECK-STACK(np,'ausgab(cs_enhance)');
  $TRANSFER PROPERTIES TO (np) FROM (np-1);
  W(NP)=W(NP-1); U(np)=U(np-1); V(np)=V(np-1);
  E(NP)=E(NP-1);IQ(NP)=IQ(NP-1); wt(np) = wt(np)/cs_enhance(iicm);
  return;
  ]
- ELSE IF( iarg = 16 | iarg = 18 | iarg = 20 | iarg = 24 ) [
- "Pair/Compton/Photo/Rayleigh just occured."
+ ELSE IF( iarg = 16 | iarg = 18 | iarg = 20 | iarg = 24 | iarg = 30) [
+ "Pair/Compton/Photo/Rayleigh/Photonuclear just occured."
  $RANDOMSET rnno35;
  IF( rnno35*cs_enhance(iicm) > 1 ) [
  "keep the original photon and trow away scattered photons "

HEN_HOUSE/user_codes/cavrznrc/cavrznrc.mortran

@@ -2406,6 +2406,8 @@ IF(ifano = 1)
  iausfl(21) = 1; "After photoelectric"
  iausfl(24) = 1; "Before Rayleigh"
  iausfl(25) = 1; "After Rayleigh"
+ iausfl(30) = 1; "Before photonuclear"
+ iausfl(31) = 1; "After photonuclear"
 
  "AUSGAB will be responsible for throwing away any photons resulting"
  "from a primary electron. ie. True equilibtrium requires that all"
@@ -2426,6 +2428,8 @@ IF( n_split > 1 ) [
  iausfl(20) = 0;
  iausfl(24) = 0;
  iausfl(25) = 0;
+ iausfl(30) = 0;
+ iausfl(31) = 0;
 
 ]
 
@@ -2447,14 +2451,16 @@ IF( use_enhance ) [
  iausfl(21) = 1; "After photoelectric"
  iausfl(24) = 1; "Before Rayleigh"
  iausfl(25) = 1; "After Rayleigh"
+ iausfl(30) = 1; "Before photonuclear"
+ iausfl(31) = 1; "After photonuclear"
  iausfl(8) = 1; "After bremsstrahlung"
  iausfl(14) = 1; "A positron has annihilated in-flight"
  iausfl(15) = 1; "A positron has annihilated at rest"
 ]
 
 write(6,'(/a)')
 '************************ IAUSFL ************************** ';
-DO j=1,28 [
+DO j=1,30 [
  IF( iausfl(j).ne.0 ) write(6,'(i3,$)') j;
 ]
 write(6,'(/a//)')
@@ -3422,9 +3428,9 @@ IF( use_enhance | n_split > 1 ) [
 IF( use_enhance ) [ "If we use cross section enhancement, all scoring "
  " is done here and the rest of ausgab is ignored "
 
- IF (iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23)
+ IF (iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23 | iarg = 29)
  [
- "A pair/Compton/photoelectric/Rayleigh event is about to take place"
+ "Just before a pair/Compton/photoelectric/Rayleigh/photonuclear event "
  "As we have increased the photon cross section by a factor of "
  "cs_enhance, we must split the photon into a scattering portion "
  "(1/cs_enhance) and a nor-scattering portion (1-1/cs_enhance) "
@@ -3464,8 +3470,8 @@ IF( use_enhance ) [ "If we use cross section enhancement, all scoring "
  return;
  ]
 
- IF( iarg = 18 | iarg = 20 | iarg = 24 |
- " A Compton/photo-absorption/Rayleigh event just occured"
+ IF( iarg = 18 | iarg = 20 | iarg = 24 | iarg = 30 |
+ " A Compton/photo-absorption/Rayleigh/photonuclear event just occured"
  iarg = 7 | iarg = 13 | iarg = 14 ) [
  " A bremas/annihilation event just occured"
  " All scattered photons and photons originating in brems/annihilation"
@@ -3548,9 +3554,9 @@ IF(IARG = 0)["ABOUT TO TRANSPORT A PARTICLE"
 
 IF (ifano = 1)
 [
- IF (iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23)
+ IF (iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23 | iarg = 29)
  [
- "A pair/Compton/photoelectric/Rayleigh event is about to take place"
+ "Just before pair/Compton/photoelectric/Rayleigh/photonuclear event"
  np = np + 1; "Boost the stack"
  IF(np + 1 > $MXSTACK)
  [
@@ -3581,6 +3587,7 @@ IF (ifano = 1)
  IF ( iarg = 18 " Compton has occured"
  | iarg = 20 " After photo-absorption "
  | iarg = 24 " After Rayleigh "
+ | iarg = 30 " After photonuclear "
  | iarg = 7 " After brems "
  | iarg = 13 " After annihilation "
  | iarg = 14) " After annihilation at rest "

HEN_HOUSE/user_codes/cavsphnrc/cavsphnrc.mortran

@@ -1881,6 +1881,8 @@ IF(ifano = 1)
  iausfl(21) = 1; "After photoelectric"
  iausfl(24) = 1; "Before Rayleigh"
  iausfl(25) = 1; "After Rayleigh"
+ iausfl(30) = 1; "Before photonuclear"
+ iausfl(31) = 1; "After photonuclear"
 
  "AUSGAB will be responsible for throwing away any photons resulting"
  "from a primary electron. ie. True equilibtrium requires that all"
@@ -1905,6 +1907,8 @@ IF( n_split > 1 ) [
  iausfl(20) = 0;
  iausfl(24) = 0;
  iausfl(25) = 0;
+ iausfl(30) = 0;
+ iausfl(31) = 0;
 
 ]
 
@@ -1926,6 +1930,8 @@ IF( use_enhance ) [
  iausfl(21) = 1; "After photoelectric"
  iausfl(24) = 1; "Before Rayleigh"
  iausfl(25) = 1; "After Rayleigh"
+ iausfl(30) = 1; "Before photonuclear"
+ iausfl(31) = 1; "After photonuclear"
  iausfl(8) = 1; "After bremsstrahlung"
  iausfl(14) = 1; "A positron has annihilated in-flight"
  iausfl(15) = 1; "A positron has annihilated at rest"
@@ -3255,12 +3261,12 @@ IF( use_enhance | n_split > 1 ) [
 IF( use_enhance ) [ "If we use cross section enhancement, all scoring "
  " is done here and the rest of ausgab is ignored "
 
- IF (iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23)
+ IF (iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23 | iarg = 29)
  [
- "A pair/Compton/photoelectric/Rayleigh event is about to take place"
+ "Before pair/Compton/photoelectric/Rayleigh/photonuclear event. "
  "As we have increased the photon cross section by a factor of "
  "cs_enhance, we must split the photon into a scattering portion "
- "(1/cs_enhance) and a nor-scattering portion (1-1/cs_enhance) "
+ "(1/cs_enhance) and a non-scattering portion (1-1/cs_enhance) "
  "Start with placing an identical photon on the stack "
  np = np + 1;
  IF(np + 1 > $MXSTACK) [
@@ -3298,8 +3304,8 @@ IF( use_enhance ) [ "If we use cross section enhancement, all scoring "
  ]
 
 
- IF( iarg = 18 | iarg = 20 | iarg = 24 |
- " A Compton/photo-absorption/Rayleigh event just occured"
+ IF( iarg = 18 | iarg = 20 | iarg = 24 | iarg = 30 |
+ " A Compton/photo-absorption/Rayleigh/photonuclear event just occured"
  iarg = 7 | iarg = 13 | iarg = 14 ) [
  " A bremas/annihilation event just occured"
  " All scattered photons and photons originating in brems/annihilation"
@@ -3392,9 +3398,9 @@ IF(IWATCH.GT.0) [CALL WATCH(IARG,IWATCH); "SIGNAL WATCH ROUTINE IF ACTIVE"]
 
 IF (ifano = 1)
 [
- IF (iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23)
+ IF (iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23 | iarg = 29)
  [
- "A pair/Compton/photoelectric/pair event is about to take place"
+ "Just before pair/Compton/photoelectric/Rayleigh/photonuclear event"
  np = np + 1; "Boost the stack"
  IF(np + 1 > $MXSTACK)
  [
@@ -3425,6 +3431,7 @@ IF (ifano = 1)
  IF ( iarg = 18 " Compton has occured"
  | iarg = 20 " After photo-absorption "
  | iarg = 24 " After Rayleigh "
+ | iarg = 30 " After photonuclear "
  | iarg = 7 " After brems "
  | iarg = 13 " After annihilation "
  | iarg = 14) " After annihilation at rest "

HEN_HOUSE/user_codes/dosrznrc/dosrznrc.mortran

@@ -2605,6 +2605,8 @@ write(6,*) 'flagged all photon intereaction types';
  iausfl(21) = 1; "After photoelectric"
  iausfl(24) = 1; "Before Rayleigh"
  iausfl(25) = 1; "After Rayleigh"
+ iausfl(30) = 1; "Before Photonuclear"
+ iausfl(31) = 1; "After Photonuclear"
 ]
 ELSE [ DO j=1,$MXREG [ iefl(j) = 0; ] ]
 
@@ -3634,8 +3636,8 @@ MXNP=MAX(MXNP,NP);"keep track of how deep stack is"
 
 IF (ienhance = 1) [" Option to enhance photon cross section in some region"
 "write(6,*) 'in ausgab to recreate photon';"
- IF (iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23) [
- "A pair/Compton/photoelectric/pair event is about to take place"
+ IF (iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23 | iarg = 29) [
+ "Just before pair/Compton/photoelectric/Rayleigh/photonuclear event"
  np = np + 1; "Boost the stack"
  IF(np > $MXSTACK) [
  OUTPUT $MXSTACK; (
@@ -3667,7 +3669,8 @@ IF (ienhance = 1) [" Option to enhance photon cross section in some region"
  " between NPold and NP for RR"
  IF ( ( iarg = 18 " Compton has occured"
  | iarg = 20 " After photo-absorption "
- | iarg = 24) " After Rayleigh "
+ | iarg = 24 " After Rayleigh "
+ | iarg = 30) " After photonuclear "
  & ienhance = 1 )
  [
  ienhance=0;

HEN_HOUSE/user_codes/egs_chamber/egs_chamber.cpp

@@ -1871,9 +1871,11 @@ int EGS_ChamberApplication::initScoring() {
  setAusgabCall(BeforeCompton, true);
  setAusgabCall(BeforePhoto, true);
  setAusgabCall(BeforeRayleigh, true);
+ setAusgabCall(BeforePhotoNuc, true);
  setAusgabCall(AfterCompton, true);
  setAusgabCall(AfterPhoto, true);
  setAusgabCall(AfterRayleigh, true);
+ setAusgabCall(AfterPhotoNuc, true);
  setAusgabCall(AfterPair, true);
  setAusgabCall(AfterTransport, true);
  setAusgabCall(AfterBrems,true);
@@ -1944,7 +1946,7 @@ int EGS_ChamberApplication::ausgab(int iarg) {
  if( cs_enhance[ig][ir] > 1 ){
  // devide photon in interacting and non-interacting portions
  if( iarg == BeforePair || iarg == BeforeCompton ||
- iarg == BeforePhoto || iarg == BeforeRayleigh ) {
+ iarg == BeforePhoto || iarg == BeforeRayleigh || iarg == BeforePhotoNuc) {
  // increase stack
  ++the_stack->np;
  if (the_stack->np > MXSTACK)
@@ -1956,7 +1958,7 @@ int EGS_ChamberApplication::ausgab(int iarg) {
  the_stack->wt[np+1] /= cs_enhance[ig][ir];
  return 0;
  }
- if( iarg == AfterCompton || iarg == AfterPhoto || iarg == AfterRayleigh || iarg == AfterPair ) {
+ if( iarg == AfterCompton || iarg == AfterPhoto || iarg == AfterRayleigh || iarg == AfterPair || iarg == AfterPhotoNuc) {
  if( rndm->getUniform()*cs_enhance[ig][ir] < 1 ){
  // remove non-interacting portion
  the_stack->wt[the_stack->npold-2] = 0;

HEN_HOUSE/user_codes/flurznrc/flurznrc.mortran

@@ -2269,7 +2269,7 @@ IF(IPRIM=3)["need to distinguish photon primaries from secondaries"
 
 write(6,'(/a)')
 '************************ IAUSFL values ************************** ';
-DO j=1,28 [
+DO j=1,30 [
  IF( iausfl(j).ne.0 ) write(6,'(i3,$)') j;
 ]
 write(6,'(/a//)')

HEN_HOUSE/user_codes/sprrznrc/sprrznrc.mortran

@@ -1934,6 +1934,8 @@ IF(ifano = 1) [
  iausfl(21) = 1; "After photoelectric"
  iausfl(24) = 1; "Before Rayleigh"
  iausfl(25) = 1; "After Rayleigh"
+ iausfl(30) = 1; "Before photonuclear"
+ iausfl(31) = 1; "After photonuclear"
 
  "AUSGAB will be responsible for throwing away any photons resulting"
  "from a primary electron. ie. True equilibtrium requires that all"
@@ -2803,8 +2805,8 @@ IF(IARG = 0)["about to transport a particle"
 
 IF (ifano = 1) [ "option to regenerate photons after interaction"
  "and throw away scattered particles"
- IF (iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23) [
- "A pair/Compton/photoelectric/pair event is about to take place"
+ IF (iarg = 15 | iarg = 17 | iarg = 19 | iarg = 23 | iarg = 29) [
+ "Just before pair/Compton/photoelectric/pair/photonuclear event"
  np = np + 1; "Boost the stack"
  IF(np + 1 > $MXSTACK) [ "check there is space on STACK"
  OUTPUT; ( ' Fano calculation unable to boost stack.'/
@@ -2829,6 +2831,7 @@ IF (ifano = 1) [ "option to regenerate photons after interaction"
  IF ( iarg = 18 " Compton has occured"
  | iarg = 20 " After photo-absorption "
  | iarg = 24 " After Rayleigh "
+ | iarg = 30 " After Photonuclear "
  | iarg = 7 " After brems "
  | iarg = 13 " After annihilation "
  | iarg = 14) " After annihilation at rest "
@@ -4933,4 +4936,3 @@ return; end;
 ;
 "end of sprrznrc.mortran"
 ;
-