Initial commit from Polytechnique Montreal

1 files changed, 776 insertions, 0 deletions

diff --git a/Dragon/src/LIBWIM.f b/Dragon/src/LIBWIM.f
new file mode 100644
index 0000000..b70c13d
--- /dev/null
+++ b/Dragon/src/LIBWIM.f
@@ -0,0 +1,776 @@
+*DECK LIBWIM
+      SUBROUTINE LIBWIM(IPLIB,IPRINT,NAMFIL,NGROUP,NBISO,NL,ISONAM,
+     >                  ISONRF,IPISO,ISHINA,TN,SN,SB,MASKI,NGF,NGFR)
+*
+*-----------------------------------------------------------------------
+*
+*Purpose:
+* Transcription of the interpolated microscopic xs read from a
+* microscopic xs library in WIMS-AECL format to LCM data structures.
+*
+*Copyright:
+* Copyright (C) 2002 Ecole Polytechnique de Montreal
+* This library is free software; you can redistribute it and/or
+* modify it under the terms of the GNU Lesser General Public
+* License as published by the Free Software Foundation; either
+* version 2.1 of the License, or (at your option) any later version
+*
+*Author(s): G. Marleau
+*
+*Parameters: input
+* IPLIB   pointer to the internal library.
+* IPRINT  print flag.
+* NAMFIL  WIMS-EACL file name.
+* NGROUP  number of groups.
+* NBISO   number of isotopes.
+* NL      number of Legendre scattering order:
+*         =1 isotropic;
+*         =2 linearly anisotropic.
+* ISONAM  local isotope names.
+* ISONRF  library isotope names.
+* IPISO   pointer array towards microlib isotopes.
+* ISHINA  self-shielding isotope names.
+* TN      isotope tempterature.
+* SN      dilution xs.
+* SB      Livolant-Jeanpierre dilution xs.
+* MASKI   logical mask for processing isotope.
+*
+*Parameters: output
+* NGF     number of fast groups without self-shielding.
+* NGFR    number of fast and resonance groups.
+*
+*-----------------------------------------------------------------------
+*
+      USE GANLIB
+      IMPLICIT NONE
+*----
+*  SUBROUTINE ARGUMENTS
+*----
+      INTEGER    NDPROC
+      PARAMETER (NDPROC=10)
+      TYPE(C_PTR) IPLIB,IPISO(NBISO)
+      INTEGER    IPRINT,NGROUP,NBISO,NL,ISONAM(3,NBISO),ISONRF(3,NBISO),
+     >           ISHINA(3,NBISO),NGF,NGFR
+      CHARACTER  NAMFIL*8
+      LOGICAL    MASKI(NBISO)
+      REAL       TN(NBISO),SN(NGROUP,NBISO),SB(NGROUP,NBISO)
+*----
+* FUNCTIONS
+*----
+      DOUBLE PRECISION XDRCST
+*----
+*  INTERNAL PARAMETERS
+*----
+      INTEGER    IOUT,ITLIB,MAXRES,MAXTEM,MAXDIL,NOTX
+      REAL       CONVM
+      PARAMETER (IOUT=6,ITLIB=1,MAXRES=50,MAXTEM=20,MAXDIL=20,NOTX=-1)
+*----
+*  LOCAL VARIABLES
+*----
+      CHARACTER  NAMDXS(NDPROC)*6,HNAMIS*12,HNISOR*12,HSHIR*8,
+     >           README*96,FMT*6
+      INTEGER    IHGAR(24),IP1,NPROC,IUNIT,KDROPN,II,NEL,NGR,NGTHER,
+     >           MXSCT,NGX,IG,ILOCX,ILOCY,ILOCS,NRDT,JSO,ITC,IDRES,IEL,
+     >           IRISO,IENDF,NF,NSCT,NTMP,IREC,JJJ,IACT,ITMP,ITXS,NTYP,
+     >           LSUBTR,LSUBZ,LRESND,IGRF,IGR,NRES,IGF,JRES,KRES,NTMPR,
+     >           NDILR,NTD,ITT,IRRICS,ILL,IGRL,IG1,IP0
+      REAL       TMPT(MAXTEM),DILT(MAXTEM),RS1(3*MAXRES),XSCOR(4),
+     >           AWJSO,RIND,XRS1,ASIGPL
+      DOUBLE PRECISION TERP(MAXTEM)
+      TYPE(C_PTR) KPLIB
+*----
+*  WIMS-AECL LIBRARY PARAMETERS
+*   IUTYPE : TYPE OF FILE = 4 (DA)
+*   LRIND  : LENGHT RECORD ON DA FILE = 256
+*   IACTO  : OPEN ACTION = 2 (READ ONLY)
+*   IACTC  : CLOSE ACTION = 2 (KEEP)
+*   MAXISO : MAX. NB. OF ISO = 246
+*   NCT    : NUMBER OF C*8 IN TITLE = 10
+*   LPZ    : LENGTH OF WIMS PARAMETER ARRAY = 9
+*   LMASTB : LENGTH OF MST TAB = MAXISO+9
+*   LMASIN : LENGTH OF MST IDX = LMASTB-4
+*   LGENTB : LENGTH OF GEN TAB = 6
+*   LGENIN : LENGTH OF GEN IDX = LGENTB
+*   LSUBTB : LENGTH OF SUB TAB = 6*MAXTEM+21-5+12
+*   LSUBIN : LENGTH OF SUB IDX = LSUBTB-12
+*   LRESTB : LENGTH OF RES TAB = 5*MAXRES
+*   LRESIN : LENGTH OF RES IDX = LRESIN
+*   MASTER : MASTER INDEX ARRAY
+*   GENINX : GENERAL INDEX ARRAY
+*   SUBINX : SUB INDEX ARRAY GENERAL
+*   SUBINR : SUB INDEX ARRAY RESONANCE
+*   RESINX : RESONANCE INDEX ARRAY
+*   IWISO  : ID OF ISOTOPE
+*   CWISO  : ISOTOPE NAMES
+*   MASTER : MASTER INDEX ARRAY
+*   GENINX : GENERAL INDEX ARRAY
+*   SUBINX : SUB INDEX ARRAY
+*----
+      INTEGER    IUTYPE,LRIND,IACTO,IACTC,MAXISO,NCT,LPZ,LMASTB,
+     >           LMASIN,LGENTB,LGENIN,LSUBTB,LSUBIN,LRESTB,
+     >           LRESIN,ILONG,ITYLCM
+      PARAMETER (IUTYPE=4,LRIND=256,IACTO=2,IACTC=1,
+     >           MAXISO=246,NCT=10,LPZ=9,LMASTB=MAXISO+9,
+     >           LMASIN=LMASTB-4,LGENTB=6,LGENIN=LGENTB,
+     >           LSUBTB=6*MAXTEM+28,LSUBIN=LSUBTB-12,
+     >           LRESTB=MAXRES*5,LRESIN=LRESTB)
+      CHARACTER  CWISO(MAXISO)*8,CTITLE(NCT)*8
+      INTEGER    MASTER(LMASTB),GENINX(LGENTB),SUBINX(LSUBTB),
+     >           SUBINR(LSUBTB),RESINX(LRESTB),NXS(MAXTEM),
+     >           ITITLE(2*NCT),NPZ(LPZ),IWISO(2*MAXISO)
+      REAL       AWR
+      INTEGER      IPRLOC
+      EQUIVALENCE (SUBINX(LSUBIN+3),AWR)
+*----
+*  ALLOCATABLE ARRAYS
+*----
+      INTEGER, ALLOCATABLE, DIMENSION(:) :: ITYPRO
+      REAL, ALLOCATABLE, DIMENSION(:) :: DELTA,XSSCMP,ENER,TMPXS,TMPSC,
+     > RRI,RIT,DSIGPL
+      REAL, ALLOCATABLE, DIMENSION(:,:) :: XSREC,XSOUT,GAR
+      REAL, ALLOCATABLE, DIMENSION(:,:,:) :: SCAT
+*----
+*  DATA
+*----
+      SAVE    NAMDXS
+      DATA    NAMDXS  /'NTOT0 ','TRANC ','NUSIGF','NFTOT ','CHI   ',
+     >                 'NU    ','NG    ','N2N   ','NGOLD ','NWT0  '/
+*----
+*  SCRATCH STORAGE ALLOCATION
+*   ITYPRO  cross section processed
+*   DELTA   lethergy
+*   XSREC   general xs vector
+*   SCAT    complete scattering matrix SCAT(JG,IG) (from IG to JG)
+*   XSSCMP  compress scattering for transfer
+*   XSOUT   self shielding parameter
+*   GAR     intermediate xs vector:
+*           GAR(I,1): library fission spectrum;
+*           GAR(I,2): potential scattering xs
+*----
+      ALLOCATE(ITYPRO(NL))
+      ALLOCATE(DELTA(NGROUP),XSREC(NGROUP,NDPROC+NL),
+     >         SCAT(NGROUP,NGROUP,NL),XSSCMP(NGROUP*(NGROUP+2)),
+     >         XSOUT(NGROUP,7),GAR(NGROUP,2))
+*----
+*  OPEN WIMSLIB AND READ TITLE
+*  READ GENERAL DIMENSIONING
+*----
+      IPRLOC=IPRINT
+      IF(IPRINT .LT. 20) IPRLOC=0
+      CONVM=REAL(XDRCST('Neutron mass','amu'))
+      IP0=NDPROC+1
+      IP1=NDPROC+2
+      NPROC=NDPROC+NL
+      IUNIT=KDROPN(NAMFIL,IACTO,IUTYPE,LRIND)
+      IF(IUNIT.LE.0) CALL XABORT('LIBWIM: WIMS-AECL LIBRARY '//
+     >    NAMFIL//' CANNOT BE OPENED FOR MIXS')
+      CALL OPNIND(IUNIT,MASTER,LMASTB)
+      CALL REDIND(IUNIT,MASTER,LMASIN,GENINX,LGENTB,1)
+      CALL REDIND(IUNIT,MASTER,LMASIN,ITITLE,2*NCT,2)
+      CALL UPCKIC(ITITLE(1),CTITLE(1),NCT)
+      WRITE(README(9:96),'(6H FROM ,10A8,A2)')
+     >  (CTITLE(II),II=1,NCT),'  '
+      IF(IPRINT.GE.5) THEN
+        WRITE(IOUT,6000) NAMFIL
+        WRITE(IOUT,'(1X,10A8)') (CTITLE(II),II=1,NCT)
+      ENDIF
+      CALL REDIND(IUNIT,GENINX,LGENIN,NPZ,LPZ,1)
+      IF(NPZ(2).NE.NGROUP) THEN
+        WRITE(IOUT,9001) NGROUP,NPZ(2)
+        CALL XABORT('LIBWIM: INVALID NUMBER OF GROUPS')
+      ENDIF
+      NEL=NPZ(1)
+      NGF=NPZ(4)
+      NGR=NPZ(5)
+      NGTHER=NPZ(6)
+      NGFR=NGF+NGR
+      MXSCT=NGROUP*(NGROUP+2)
+      IF(NGFR+NGTHER.NE.NGROUP) THEN
+        WRITE(IOUT,9001) NGROUP,NGFR+NGTHER
+        CALL XABORT('LIBWIM: INVALID NUMBER OF GROUPS')
+      ENDIF
+      IF(NEL.GT.MAXISO) THEN
+        WRITE(IOUT,9003) MAXISO,NEL
+        CALL XABORT('LIBWIM: INVALID NUMBER OF ISOTOPES')
+      ENDIF
+      ALLOCATE(DSIGPL(NGR))
+*----
+*  READ ISOTOPES NAMES
+*----
+      CALL REDIND(IUNIT,GENINX,LGENIN,IWISO,2*NEL,3)
+      CALL UPCKIC(IWISO(1),CWISO(1),NEL)
+      CALL REDIND(IUNIT,GENINX,LGENIN,IWISO,NEL,2)
+*----
+*  READ GROUP STRUCTURE
+*----
+      ALLOCATE(ENER(NGROUP+1))
+      CALL REDIND(IUNIT,GENINX,LGENIN,ENER,NGROUP+1,4)
+      IF(ENER(NGROUP+1).EQ.0.0) ENER(NGROUP+1)=1.0E-5
+      CALL LCMPUT(IPLIB,'ENERGY',NGROUP+1,2,ENER)
+      NGX=0
+      DO 100 IG=1,NGROUP
+        IF(NGX.EQ.0.AND.ENER(IG+1).LT.4.0) NGX=IG-1
+        DELTA(IG)=LOG(ENER(IG)/ENER(IG+1))
+ 100  CONTINUE
+      CALL LCMPUT(IPLIB,'DELTAU',NGROUP,2,DELTA)
+      DEALLOCATE(ENER)
+*----
+*  INITIALIZE ALL XSREC
+*  READ FISSION SPECTRUM
+*----
+      GAR(:NGROUP,1)=0.0
+      CALL REDIND(IUNIT,GENINX,LGENIN,GAR(:,1),NPZ(3),5)
+*----
+*  ALLOCATE MEMORY FOR TEMPERATURE DEPENDENT XS
+*  AND FOR RESONANCE CALCULATION
+*----
+      ALLOCATE(TMPXS(5*NGROUP),TMPSC(NGROUP*NGROUP),
+     > RRI(MAXDIL*MAXTEM*2),RIT(MAXDIL))
+      ILOCX=0
+      ILOCY=NGFR
+      ILOCS=0
+      NRDT=NGTHER-1
+*----
+*  READ THROUGH DRAGON FILE AND ACCUMULATE CROSS SECTIONS FOR
+*  CROSS SECTION ARE SAVED ONLY IF ISOTOPE IS USED
+*----
+      DO 110 JSO=1,NBISO
+        IF(.NOT.MASKI(JSO)) GO TO 115
+*----
+*  LOCATE ISOTOPE
+*----
+        WRITE(HNAMIS,'(3A4)') (ISONAM(ITC,JSO),ITC=1,3)
+        WRITE(HNISOR,'(3A4)') (ISONRF(ITC,JSO),ITC=1,3)
+        WRITE(HSHIR,'(2A4)') (ISHINA(ITC,JSO),ITC=1,2)
+        IDRES=INDEX(HSHIR,'.')
+        IF(IDRES.GT.0) THEN
+          WRITE(FMT,'(2H(F,I1,3H.1))') IDRES+1
+          READ(HSHIR,FMT) RIND
+        ENDIF
+        IRISO=0
+        DO 120 IEL=1,NEL
+          IF(CWISO(IEL).EQ.HNISOR(1:8)) THEN
+            IRISO=IEL
+            IF(IDRES.EQ.0) THEN
+              RIND=FLOAT(IWISO(IRISO))
+            ENDIF
+            GO TO 125
+          ENDIF
+ 120    CONTINUE
+        WRITE(IOUT,9002) HNISOR,NAMFIL
+        CALL XABORT('LIBWIM: ISOTOPE NOT FOUND ON LIBRARY')
+ 125    CONTINUE
+        IF(IPRINT.GE.5) WRITE(IOUT,6001) HNAMIS
+        XSREC(:NGROUP,:NPROC)=0.0
+        SCAT(:NGROUP,:NGROUP,:NL)=0.0
+*----
+*  READ SUB INDEX ASSOCIATED WITH ISOTOPE
+*----
+        CALL REDIND(IUNIT,MASTER,LMASIN,SUBINX,LSUBTB,IRISO+4)
+*----
+*  FOR ENDF/B-VI LIBRARY : IENDF = 2
+*  FOR ENDF/B-V  LIBRARY : IENDF = 1
+*  FOR WINFRITH  LIBRARY : IENDF = 0
+*----
+        IENDF=SUBINX(LSUBIN+12)
+        AWJSO=AWR/CONVM
+*----
+*  FAST AND/OR RESONANCE XS
+*----
+        CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(NGF+1:,9),NGR,9)
+        DSIGPL(:NGR)=0.0
+        IF(IENDF.EQ.0) THEN
+          CALL REDIND(IUNIT,SUBINX,LSUBIN,GAR(NGF+1:,2),NGR,2)
+          DO 130 IG=NGF+1,NGFR
+            DSIGPL(IG-NGF)=GAR(IG,2)*XSREC(IG,9)
+ 130      CONTINUE
+        ENDIF
+        CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(:,8),NGF,13)
+        NF=SUBINX(LSUBIN+5)
+        IF(NF.GT.1) THEN
+          CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(:,3),NGFR,10)
+          CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(:,4),NGFR,12)
+          DO 135 IG=1,NGROUP
+            XSREC(IG,5)=GAR(IG,1)
+ 135      CONTINUE
+        ENDIF
+        NSCT=SUBINX(LSUBIN+8)
+        IF(NSCT.GT.MXSCT) THEN
+          WRITE(IOUT,9004) NSCT,MXSCT
+          CALL XABORT('LIBWIM: INVALID SCATTERING MATRIX -E/0')
+        ENDIF
+*----
+*  READ AND DECOMPRESS P0 SCATTERING CROSS SECTIONS
+*  COMPUTE TOTAL P0 SCATTERING OUT OF GROUP
+*----
+        CALL REDIND(IUNIT,SUBINX,LSUBIN,XSSCMP,NSCT,14)
+        CALL LIBWSC(NGROUP,1,NGFR,NSCT,XSSCMP,SCAT(1,1,1),XSREC(1,IP0))
+*----
+*  FOR IENDF=2   READ XS FOR NG AND TOTAL
+*  FOR IENDF=0,1 READ XS FOR NG AND TRANSPORT
+*----
+        IF(IENDF.GE.2) THEN
+*----
+*  READ TOTAL XS FOR IENDF=2
+*----
+          CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(:,2),NGFR,5)
+          CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(:,7),NGFR,5)
+        ELSE
+*----
+*  COMPUTE TOTAL XS FOR IENDF=0,1
+*----
+          CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(:,2),NGFR,4)
+          CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(:,7),NGFR,6)
+        ENDIF
+*----
+*  THERMAL XS
+*----
+        NTMP=SUBINX(LSUBIN+6)
+        IF(NTMP.GT.MAXTEM) THEN
+          CALL XABORT('LIBWIM: INVALID MAXTEM FOR P0.')
+        ENDIF
+        CALL REDIND(IUNIT,SUBINX,LSUBIN,NXS,NTMP,3)
+        CALL REDIND(IUNIT,SUBINX,LSUBIN,TMPT,NTMP,15)
+        IREC=16
+        IF(NTMP.EQ.1) THEN
+          IF(IPRINT.GE.100) THEN
+            WRITE(IOUT,6200) TN(JSO)
+          ENDIF
+          IREC=IREC+2
+          IF(NF.GT.1) THEN
+            CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(NGFR+1:,3),
+     >        NGTHER,IREC)
+            CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(NGFR+1:,4),
+     >        NGTHER,IREC+1)
+          ENDIF
+          IREC=IREC+2
+          IF(NXS(1).GT.MXSCT) THEN
+            WRITE(IOUT,9004) NXS(1),MXSCT
+            CALL XABORT('LIBWIM: INVALID SCATTERING MATRIX -T/0')
+          ENDIF
+          CALL REDIND(IUNIT,SUBINX,LSUBIN,XSSCMP,NXS(1),IREC)
+          IREC=IREC+1
+*----
+*  DECOMPRESS P0 SCATTERING CROSS SECTIONS AND COMPUTE
+*  P0 SCATTERING OUT OF GROUP
+*----
+          CALL LIBWSC(NGROUP,NGFR+1,NGROUP,NXS(1),XSSCMP,
+     >                SCAT(1,1,1),XSREC(1,IP0))
+          IF(IENDF.GE.2) THEN
+            CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(NGFR+1:,2),
+     >                  NGTHER,IREC-4)
+          ELSE
+            CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(NGFR+1:,2),
+     >                  NGTHER,IREC-5)
+          ENDIF
+          CALL REDIND(IUNIT,SUBINX,LSUBIN,XSREC(NGFR+1:,7),
+     >                NGTHER,IREC-4)
+        ELSE IF(NTMP.GT.1) THEN
+*----
+*  AVALUATE LAGRANGIAN INTERPOLATION FACTOR FOR
+*  AVAILABLE TEMPERATURES (ORDER NOTX) AND INTERPOLATE.
+*----
+          CALL LIBLEX(NTMP,TN(JSO),TMPT,NOTX,TERP)
+          IF(IPRINT.GE.100) THEN
+            WRITE(IOUT,6201) TN(JSO)
+            WRITE(IOUT,6202) (TMPT(JJJ),JJJ=1,NTMP)
+            WRITE(IOUT,6203) (TERP(JJJ),JJJ=1,NTMP)
+          ENDIF
+          NRDT=NGTHER-1
+          IACT=1
+          DO 140 ITMP=1,NTMP
+            IF(TERP(ITMP).EQ.0.0D0) THEN
+              IREC=IREC+5
+            ELSE
+              IREC=IREC+2
+              IF(NF.GT.1) THEN
+               CALL REDIND(IUNIT,SUBINX,LSUBIN,TMPXS(ILOCY+2*NGROUP+1:),
+     >                     NGTHER,IREC)
+               CALL REDIND(IUNIT,SUBINX,LSUBIN,TMPXS(ILOCY+3*NGROUP+1:),
+     >                     NGTHER,IREC+1)
+              ELSE
+                TMPXS(ILOCY+2*NGROUP+1:ILOCY+2*NGROUP+NGTHER)=0.0
+                TMPXS(ILOCY+3*NGROUP+1:ILOCY+3*NGROUP+NGTHER)=0.0
+              ENDIF
+              IREC=IREC+2
+              IF(NXS(ITMP).GT.MXSCT) THEN
+                WRITE(IOUT,9004) NXS(ITMP),MXSCT
+                CALL XABORT('LIBWIM: INVALID SCATTERING MATRIX -T/0')
+              ENDIF
+              CALL REDIND(IUNIT,SUBINX,LSUBIN,XSSCMP,NXS(ITMP),IREC)
+              IREC=IREC+1
+*----
+*  DECOMPRESS P0 SCATTERING CROSS SECTIONS AND COMPUTE
+*  P0 SCATTERING OUT OF GROUP
+*----
+              CALL LIBWSC(NGROUP,NGFR+1,NGROUP,NXS(ITMP),XSSCMP,
+     >                    TMPSC(ILOCS+1),TMPXS(ILOCX+4*NGROUP+1))
+              IF(IENDF.GE.2) THEN
+                CALL REDIND(IUNIT,SUBINX,LSUBIN,TMPXS(ILOCY+1:),
+     >                      NGTHER,IREC-4)
+              ELSE
+                CALL REDIND(IUNIT,SUBINX,LSUBIN,TMPXS(ILOCY+1:),
+     >                      NGTHER,IREC-5)
+              ENDIF
+              CALL REDIND(IUNIT,SUBINX,LSUBIN,TMPXS(ILOCY+NGROUP+1:),
+     >                    NGTHER,IREC-4)
+              ITXS=1
+              CALL LIBWTE(IACT,ITXS,NGROUP,NGTHER,1,NF,TERP(ITMP),
+     >                    SCAT(1,1,1),XSREC(1,IP0),XSREC(1,7),
+     >                    XSREC(1,3),XSREC(1,4),XSREC(1,2),TMPXS,TMPSC)
+              IACT=2
+            ENDIF
+ 140      CONTINUE
+        ENDIF
+*----
+*  BUILT CROSS SECTION FROM INFORMATION IN NG WHICH IS
+*  CURRENTLY ABSORPTION AND SCATTERING OUT OF GROUP
+*  COMPUTE REAL NG CROSS SECTION WHICH IS
+*  CURRENT NG (ABSORPTION)-FISSION+N2N
+*  SINCE ABSORPTION IS DEFINED AS
+*  TOTAL-SIGS WHERE SIGS CONTAINE 2*N2N SINCE A N2N CONTRIBUTION
+*  PRODUCES AN EQUIVALENT OF 2 NEUTRON BY DIFFUSION
+*----
+        DO 150 IG=1,NGROUP
+          XSREC(IG,1)=XSREC(IG,7)+XSREC(IG,IP0)
+          IF(NF.GT.1) THEN
+            XSREC(IG,7)=XSREC(IG,7)+XSREC(IG,8)-XSREC(IG,4)
+          ELSE
+            XSREC(IG,7)=XSREC(IG,7)+XSREC(IG,8)
+          ENDIF
+          IF(XSREC(IG,4).NE.0) THEN
+            XSREC(IG,6)=XSREC(IG,3)/XSREC(IG,4)
+          ELSE
+            XSREC(IG,6)=0
+          ENDIF
+ 150    CONTINUE
+        IF(IENDF.LT.2) THEN
+*----
+*  COMPUTE TRANSPORT CORRECTION AND STORE IN TRAN
+*----
+          DO 151 IG=1,NGROUP
+            XSREC(IG,2)=XSREC(IG,1)-XSREC(IG,2)
+ 151      CONTINUE
+        ENDIF
+*----
+*  SELF SHIELDING DATA
+*----
+        NTYP=1
+        XSCOR(1)=0.0
+        IF(SUBINX(LSUBIN+5).EQ.3) THEN
+          NTYP=2
+          XSCOR(2)=0.0
+        ENDIF
+*----
+*  MODIFIED SUB IDX LENGTH FOR RESONANCE
+*----
+        LSUBTR=NGR+7
+        LSUBZ=NGR+1
+        CALL REDIND(IUNIT,MASTER,LMASIN,SUBINR,LSUBTR,NEL+5)
+*----
+*  MODIFIED RES IDX LENGTH FOR RESONANCE
+*----
+        LRESND=SUBINR(NGR+6)
+        IF(NTYP.EQ.2.AND.SUBINR(NGR+7).EQ.1) THEN
+          NTYP=3
+          XSCOR(3)=0.0
+        ENDIF
+        XSCOR(4)=0.0
+        IGRF=NGF
+        KRES=0
+        DO 300 IGR=1,NGR
+          IGRF=IGRF+1
+          CALL REDIND(IUNIT,SUBINR,LSUBZ,RESINX,LRESND+1,IGR)
+          NRES=RESINX(LRESND+1)
+          IF(NRES.GT.MAXRES) THEN
+            WRITE(IOUT,9005) NRES,MAXRES
+            CALL XABORT('LIBWIM: INVALID NUMBER OF RESONANCE')
+          ENDIF
+          IF(IGR.EQ.1) THEN
+            CALL REDIND(IUNIT,RESINX,LRESND,RS1,3*NRES,1)
+            DO 314 IGF=1,NGFR
+              XSOUT(IGF,3)=XSREC(IGF,IP0)
+              XSOUT(IGF,4)=1.0
+              XSOUT(IGF,5)=1.0
+ 314        CONTINUE
+*----
+*  IDENTIFY SELF SHIELDING RESONNANT ISOTOPE
+*----
+            DO 310 JRES=1,NRES
+              IF(IDRES.EQ.0) THEN
+                XRS1=FLOAT(INT((RS1(3*(JRES-1)+1)+0.01)*10.)
+     >              -INT(RS1(3*(JRES-1)+1)+0.01)*10)/10.+0.02
+                XRS1=ABS(RS1(3*(JRES-1)+1)-XRS1-RIND)
+              ELSE
+                XRS1=ABS(RS1(3*(JRES-1)+1)-RIND)
+              ENDIF
+              IF(XRS1.LE.0.01) THEN
+                KRES=JRES
+                NTMPR=INT(RS1(3*(KRES-1)+2)+0.1)
+                NDILR=INT(RS1(3*(KRES-1)+3)+0.1)
+                IF(NTMPR.GT.MAXTEM) THEN
+                  WRITE(IOUT,9006) NTMPR,MAXTEM
+                  CALL XABORT('LIBWIM: INVALID NUMBER OF RES TEMP')
+                ELSE IF(NDILR.GT.MAXTEM) THEN
+                  WRITE(IOUT,9007) NDILR,MAXTEM
+                  CALL XABORT('LIBWIM: INVALID NUMBER OF RES DIL')
+                ENDIF
+                NTD=NDILR*NTMPR
+                IF(IPRINT.GE.5) THEN
+                  WRITE(IOUT,6002) RS1(3*(JRES-1)+1)
+                ENDIF
+                CALL REDIND(IUNIT,RESINX,LRESND,TMPT,NTMPR,2+5*(KRES-1))
+                CALL REDIND(IUNIT,RESINX,LRESND,DILT,NDILR,3+5*(KRES-1))
+                IF(IPRINT.GE.100) THEN
+                  WRITE(IOUT,6003) TN(JSO)
+                  WRITE(IOUT,6008) (TMPT(ITT),ITT=1,NTMPR)
+                  WRITE(IOUT,6004) SN(IGRF,JSO),DSIGPL(IGR)
+                  WRITE(IOUT,6008) (DILT(ITT),ITT=1,NDILR)
+                ENDIF
+                DO 312 II=1,NTMPR
+                  TMPT(II)=SQRT(TMPT(II))
+ 312            CONTINUE
+                DO 313 II=1,NDILR
+                  IF(DILT(II)-DSIGPL(IGR).GT.0.0) THEN
+                    DILT(II)=SQRT(DILT(II)-DSIGPL(IGR))
+                  ELSE
+                    DILT(II)=0.0
+                  ENDIF
+ 313            CONTINUE
+                GO TO 311
+              ENDIF
+ 310        CONTINUE
+*----
+*  NO SELF SHIELDING DATA FOR THIS ISOTOPE EXIT TO 301
+*----
+            XSREC(:NGROUP,10)=0.0
+            GO TO 301
+          ENDIF
+*----
+*  READ SELF SHIELDING DATA FOR THIS ISOTOPE
+*----
+ 311      CONTINUE
+*----
+*  READ FLUX FOR THIS RESONANCE INTEGRAL
+*----
+          IF(IENDF.GE.2) THEN
+*----
+*  READ TOTAL RR AND FLUX
+*----
+            CALL REDIND(IUNIT,RESINX,LRESND,RRI,2*NTD,4+5*(KRES-1))
+            CALL LIBWRI(NTMPR,NDILR,TN(JSO),SN(IGRF,JSO),TMPT,
+     >                  DILT,RRI(1),RIT,XSOUT(IGRF,1),XSCOR(1))
+            CALL LIBWRI(NTMPR,NDILR,TN(JSO),SN(IGRF,JSO),TMPT,DILT,
+     >                  RRI(NTD+1),RIT,XSOUT(IGRF,4),XSCOR(4))
+          ELSE
+*----
+*  READ TOTAL RR
+*----
+            CALL REDIND(IUNIT,RESINX,LRESND,RRI,NTD,4+5*(KRES-1))
+            CALL LIBWRI(NTMPR,NDILR,TN(JSO),SN(IGRF,JSO),TMPT,
+     >                  DILT,RRI(1),RIT,XSOUT(IGRF,1),XSCOR(1))
+            IF(IPRINT.GE.100) THEN
+              WRITE(IOUT,6005) XSOUT(IGRF,1)
+              WRITE(IOUT,6008) (RRI(ITT),ITT=1,NTD)
+            ENDIF
+          ENDIF
+          IF(NTYP.GE.2) THEN
+*----
+*  READ FISSION RR
+*----
+            CALL REDIND(IUNIT,RESINX,LRESND,RRI,NTD,5+5*(KRES-1))
+            CALL LIBWRI(NTMPR,NDILR,TN(JSO),SN(IGRF,JSO),TMPT,DILT,
+     >                  RRI(1),RIT,XSOUT(IGRF,2),XSCOR(2))
+            IF(IPRINT.GE.100) THEN
+              WRITE(IOUT,6006) XSOUT(IGRF,2)
+              WRITE(IOUT,6008) (RRI(ITT),ITT=1,NTD)
+            ENDIF
+            IF(NTYP.GE.3) THEN
+*----
+*  READ SCATTERING RR
+*----
+              CALL REDIND(IUNIT,RESINX,LRESND,RRI,NTD,
+     >                    6+5*(KRES-1))
+              IRRICS=0
+              ASIGPL=0.0
+              DO 340 ILL=1,NTD
+                ASIGPL=ASIGPL+RRI(IRRICS+1)
+                IRRICS=IRRICS+1
+ 340          CONTINUE
+              IF(ASIGPL.GT.0.0) THEN
+                CALL LIBWRI(NTMPR,NDILR,TN(JSO),SN(IGRF,JSO),TMPT,
+     >                      DILT,RRI(1),RIT,XSOUT(IGRF,3),XSCOR(3))
+               IF(IPRINT.GE.100) THEN
+                  WRITE(IOUT,6007) XSOUT(IGRF,3)
+                  WRITE(IOUT,6008) (RRI(ITT),ITT=1,NTD)
+                ENDIF
+              ENDIF
+             ENDIF
+          ENDIF
+ 300    CONTINUE
+*----
+*  CORRECT CROSS SECTIONS FOR CURRENT GROUP
+*----
+        IGRL=IGRF
+        IGRF=NGF+1
+        CALL LIBWRE(NTYP,IPRINT,ITLIB,NGROUP,1,IGRF,IGRL,NGR,
+     >              SCAT,XSREC(1,IP0),XSREC(1,1),XSREC(1,7),
+     >              XSREC(1,3),XSREC(1,4),XSREC(1,6),
+     >              DELTA,SN(1,JSO),SB(1,JSO),XSOUT,XSCOR,
+     >              DSIGPL)
+*----
+*  PRINT CROSS SECTIONS IF REQUIRED
+*----
+        IF(IPRINT.GE.5) THEN
+          WRITE(IOUT,6100)
+          DO 400 IG1=NGF+1,NGFR
+            WRITE(IOUT,6101) IG1,SN(IG1,JSO),SB(IG1,JSO),
+     >                       XSOUT(IG1,4),XSREC(IG1,1),
+     >                       XSREC(IG1,IP0),XSREC(IG1,3),
+     >                       XSREC(IG1,9)
+ 400      CONTINUE
+        ENDIF
+*----
+*  SET NWT0 THE RESONANCE FLUX WEIGHTING
+*----
+        XSREC(:NGROUP,10)=1.0
+        DO 401 IG1=NGF+1,NGFR
+          XSREC(IG1,10)=XSOUT(IG1,4)
+ 401    CONTINUE
+ 301    CONTINUE
+*----
+*  P1 SCATTERING
+*----
+        IF(NL.EQ.2) THEN
+          IREC=16+NTMP*5
+          NTMP=SUBINX(LSUBIN+10)
+          IF(NTMP+1.GT.MAXTEM) THEN
+            CALL XABORT('LIBWIM: INVALID MAXTEM FOR P1.')
+          ELSE IF(NTMP.GT.0) THEN
+            CALL REDIND(IUNIT,SUBINX,LSUBIN,NXS,NTMP+1,7)
+            CALL REDIND(IUNIT,SUBINX,LSUBIN,TMPT,NTMP,5)
+            NSCT=NXS(NTMP+1)
+            IF(NSCT.GT.MXSCT) THEN
+              WRITE(IOUT,9004) NSCT,MXSCT
+              CALL XABORT('LIBWIM: INVALID SCATTERING MATRIX -E/1')
+            ENDIF
+            CALL REDIND(IUNIT,SUBINX,LSUBIN,XSSCMP,NSCT,IREC)
+            IREC=IREC+1
+*----
+*  DECOMPRESS P1 SCATTERING CROSS SECTIONS AND COMPUTE
+*  P1 SCATTERING OUT OF GROUP
+*----
+            CALL LIBWSC(NGROUP,1,NGFR,NSCT,XSSCMP,
+     >                  SCAT(1,1,2),XSREC(1,IP1))
+          ENDIF
+          IF(NTMP.EQ.1) THEN
+            IF(NXS(1).GT.MXSCT) THEN
+              WRITE(IOUT,9004) NXS(1),MXSCT
+              CALL XABORT('LIBWIM: INVALID SCATTERING MATRIX -T/1')
+            ENDIF
+            CALL REDIND(IUNIT,SUBINX,LSUBIN,XSSCMP,NXS(1),IREC)
+            IREC=IREC+1
+*----
+*  DECOMPRESS P1 SCATTERING CROSS SECTIONS AND COMPUTE
+*  P1 SCATTERING OUT OF GROUP
+*----
+            CALL LIBWSC(NGROUP,NGFR+1,NGROUP,NXS(1),XSSCMP,
+     >                  SCAT(1,1,2),XSREC(1,IP1))
+          ELSE IF(NTMP.GT.1) THEN
+*----
+*  AVALUATE LAGRANGIAN INTERPOLATION FACTOR FOR
+*  AVAILABLE TMPTERATURES (ORDER NOTX) AND INTERPOLATE.
+*----
+            CALL LIBLEX(NTMP,TN(JSO),TMPT,NOTX,TERP)
+            NRDT=NGTHER-1
+            IACT=1
+            DO 170 ITMP=1,NTMP
+              IF(TERP(ITMP).EQ.0.0D0) THEN
+                IREC=IREC+1
+              ELSE
+                IF(NXS(ITMP).GT.MXSCT) THEN
+                  WRITE(IOUT,9004) NXS(ITMP),MXSCT
+                  CALL XABORT('LIBWIM: INVALID SCATTERING MATRIX -T/1')
+                ENDIF
+                CALL REDIND(IUNIT,SUBINX,LSUBIN,XSSCMP,NXS(ITMP),IREC)
+                IREC=IREC+1
+*----
+*  DECOMPRESS P1 SCATTERING CROSS SECTIONS AND COMPUTE
+*  P1 SCATTERING OUT OF GROUP
+*----
+                CALL LIBWSC(NGROUP,NGFR+1,NGROUP,NXS(ITMP),XSSCMP,
+     >                      TMPSC(ILOCS+1),TMPXS(ILOCX+4*NGROUP+1))
+                ITXS=2
+                CALL LIBWTE(IACT,ITXS,NGROUP,NGTHER,1,NF,TERP(ITMP),
+     >                      SCAT(1,1,2),XSREC(1,IP1),XSREC(1,7),
+     >                      XSREC(1,3),XSREC(1,4),XSREC(1,2),
+     >                      TMPXS,TMPSC)
+                IACT=2
+              ENDIF
+ 170        CONTINUE
+          ENDIF
+        ENDIF
+*----
+*  SAVE MAIN CROSS SECTIONS ON LCM
+*----
+        KPLIB=IPISO(JSO) ! set JSO-th isotope
+        CALL LCMPTC(KPLIB,'ALIAS',12,HNAMIS)
+        CALL LCMPUT(KPLIB,'AWR',1,2,AWJSO)
+        CALL XDRLGS(KPLIB,1,IPRLOC,0,NL-1,1,NGROUP,XSREC(1,NDPROC+1),
+     >              SCAT,ITYPRO)
+        CALL XDRLXS(KPLIB,1,IPRLOC,NDPROC,NAMDXS,1,NGROUP,XSREC)
+        CALL LCMLEN(KPLIB,'NTOT0',ILONG,ITYLCM)
+        IF(ILONG.EQ.0) CALL LCMPUT(KPLIB,'NTOT0',NGROUP,2,XSREC(1,1))
+        WRITE(README(:8),'(A8)') HNAMIS(1:8)
+        READ(README,'(24A4)') (IHGAR(II),II=1,24)
+        CALL LCMPUT(KPLIB,'README',24,3,IHGAR)
+        IF(IPRINT.GE.100) CALL LCMLIB(KPLIB)
+ 115  CONTINUE
+ 110  CONTINUE
+      DEALLOCATE(RIT,RRI,TMPSC,TMPXS)
+      CALL CLSIND(IUNIT)
+*----
+*  SCRATCH STORAGE DEALLOCATION
+*----
+      DEALLOCATE(DSIGPL,GAR,XSOUT,XSSCMP,SCAT,XSREC,DELTA)
+      DEALLOCATE(ITYPRO)
+*----
+*  RETURN
+*----
+      RETURN
+*----
+*  FORMAT
+*----
+ 9001 FORMAT(/'  NUMBER OF GROUPS SPECIFIED :',I10/
+     >        ' NUMBER OF GROUPS IN LIBRARY :',I10)
+ 9002 FORMAT(/' LIBWIM: MATERIAL/ISOTOPE ',A12,' IS MISSING ON WIMS',
+     >        ' FILE NAME ',A8)
+ 9003 FORMAT(/' MAXIMUM NUMBER OF ISOTOPE SPECIFIED :',I10/
+     >        '        NUMBER OF ISOTOPE IN LIBRARY :',I10)
+ 9004 FORMAT(/'        DIMENSION OF SCATTERING MATRIX  :',I10/
+     >        ' MAXIMUM DIMENSION OF SCATTERING MATRIX :',I10)
+ 9005 FORMAT(/'         NUMBER OF RESONANT ISOTOPES :',I10/
+     >        ' MAXIMUM NUMBER OF RESONANT ISOTOPES :',I10)
+ 9006 FORMAT(/'         NUMBER OF RESONANT TEMPERATURE :',I10/
+     >        ' MAXIMUM NUMBER OF RESONANT TEMPERATURE :',I10)
+ 9007 FORMAT(/'         NUMBER OF RESONANT DILUTION :',I10/
+     >        ' MAXIMUM NUMBER OF RESONANT DILUTION :',I10)
+ 6000 FORMAT(' READING WIMS-AECL LIBRARY NAME ',A8)
+ 6001 FORMAT('   PROCESSING ISOTOPE/MATERIAL = ',A12)
+ 6002 FORMAT('   SELF SHIELDING ISOTOPE      = ',F9.3)
+ 6003 FORMAT('   RESONANCE TEMPERATURE       = ',1P,E15.7)
+ 6004 FORMAT('   RESONANCE DILUTIONS         = ',1P,2E15.7)
+ 6005 FORMAT('   ABSORPTION RATE             = ',1P,E15.7)
+ 6006 FORMAT('   FISSION RATE                = ',1P,E15.7)
+ 6007 FORMAT('   SCATTERING RATE             = ',1P,E15.7)
+ 6008 FORMAT(1P,5E15.7)
+ 6100 FORMAT(/5X,'GROUP',10X,'DILUT',13X,'SB',11X,'NPHI',10X,'NTOT0',
+     > 11X,'SIGS',9X,'NUSIGF',10X,'NGOLD')
+ 6101 FORMAT(5X,I5,1P,8E15.5)
+ 6200 FORMAT(' TEMPERATURE = ',F10.5,10X,
+     >  ' CROSS SECTION TABULATED AT A SINGLE TEMPERATURE')
+ 6201 FORMAT(' TEMPERATURE = ',F10.5,10X,
+     >  ' CROSS SECTION TABULATED AT MULTIPLE TEMPERATURES')
+ 6202 FORMAT(' TABULATION TEMPERATURES= ',/(5F15.5))
+ 6203 FORMAT(' INTERPOLATION FACTORS  = ',1P,/(5E15.5))
+      END