text:
      PROGRAM FLXINT
      USE MOD_ZA  ! HYCOM array I/O interface
      IMPLICIT NONE
C
C     FLUX ARRAYS.
C
      INTEGER, ALLOCATABLE :: MSK(:,:)
      REAL*4,  ALLOCATABLE :: PLON(:,:),PLAT(:,:),XAF(:,:),YAF(:,:)
      REAL*4,  ALLOCATABLE :: TAM(:,:),HAM(:,:),
     +                        FRM(:,:),FPM(:,:),PCM(:,:)
C
      LOGICAL   LARCTIC
      INTEGER   IWI,JWI
      REAL*4    XFIN,YFIN,DXIN,DYIN
      REAL*4    PLATIJ,YFMIN,YFMAX,XAMAX,XAMIN,YAMAX,YAMIN
      CHARACTER PREAMBL(5)*79
C
      REAL*8,    ALLOCATABLE :: WLAT(:),WLON(:)
      REAL*4,    ALLOCATABLE :: TAIRIN(:,:),HAIRIN(:,:),
     +                          FLXRIN(:,:),FLXPIN(:,:),
     +                          PCIPIN(:,:)
      INTEGER*2, ALLOCATABLE ::  FLDI2(:,:)
      INTEGER,   ALLOCATABLE :: MASKIN(:,:)
C
C     INTERPOLATION ARRAYS.
C
      INTEGER IBD(4)
C
      REAL*4,  ALLOCATABLE :: FXI(:),FYI(:),WK(:)
      REAL*4,  ALLOCATABLE :: TAIRI(:,:),HAIRI(:,:),
     +                        FLXRI(:,:),FLXPI(:,:),
     +                        PCIPI(:,:)
      REAL*4,  ALLOCATABLE :: WTAIR3(:,:,:)
C
C     NAMELIST.
C
      INTEGER          JPR
      COMMON/NPROCS/   JPR
      SAVE  /NPROCS/
C
      CHARACTER*79     CTITLE
      CHARACTER*6      CNAME(5)
      NAMELIST/AFTITL/ CTITLE,CNAME
      REAL*8           FSTART,WSTART,TSTART,TMAX,
     +                 BIASPC,BIASRD,PCMEAN,RDMEAN
      REAL*4           HMKS,RMKS,PMKS,SEAMSK
      NAMELIST/AFTIME/ FSTART,WSTART,TSTART,TMAX,
     +                 BIASPC,BIASRD,PCMEAN,RDMEAN,
     +                 HMKS,RMKS,PMKS,
     +                 SEAMSK
      INTEGER          IFFILE,IFTYPE,INTERP,INTMSK,NGLOBE
      NAMELIST/AFFLAG/ IFFILE,IFTYPE,INTERP,INTMSK,NGLOBE,JPR
C
C**********
C*
C 1)  FROM UNFORMATTED Ta-Ha-Qr-Qp-Pc DATA ON ITS NATIVE GRID,
C      CREATE A MODEL GRID FLUX FILE SUITABLE FOR INPUT TO
C      THE HYCOM OCEAN MODEL OVER THE GIVEN REGION.
C
C     INTERPOLATION IS EITHER PIECEWISE BILINEAR, OR PIECEWISE CUBIC
C      BESSEL (I.E. CUBIC HERMITE INTERPOLATION WITH DERIVATIVES
C      APPROXIMATED BY CENTERED DIFFERENCES), OR PIECEWISE BICUBIC
C      (WITH DERIVATIVES APPROXIMATED BY CENTERED DIFFERENCES), OR
C      CUBIC SPLINE.
C
C 2)  PARAMETERS:
C
C     NATIVE FLUX GRID SPECIFICATION, SEE (4):
C
C        IWI    = 1ST DIMENSION OF FLUX GRID
C        JWI    = 2ND DIMENSION OF FLUX GRID
C        XFIN   = LONGITUDE OF 1ST FLUX GRID POINT
C        YFIN   = LATITUDE  OF 1ST FLUX GRID POINT
C        DXIN   = FLUX LONGITUDINAL GRID SPACING
C        DYIN   = FLUX LATITUDINAL  GRID SPACING
C                  =0.0; GAUSSIAN GRID WITH JWI/2 NODES PER HEMISPHERE.
C
C 3)  NAMELIST INPUT:
C
C     /AFTITL/
C        CTITLE - ONE (79-CHARACTER) LINE TITLE.
C        CNAME  - TWO  (6-CHARACTER) LINE NAMEs (default "tauewd","taunwd").
C
C     /AFTIME/
C        FSTART - TIME OF HEAT FLUX START                  (DAYS)
C        WSTART - TIME OF WIND START, IGNORED HERE         (DAYS)
C        TMAX   - TIME OF END   OF CURRENT INTEGRATION     (DAYS)
C        TSTART - TIME OF START OF CURRENT INTEGRATION     (DAYS)
C        BIASPC - PRECIPITATION  CORRECTION TO CLOSE SALINITY (M/S)
C                  =0.0; NO CORRECTION (DEFAULT)
C        PCMEAN - PRECIPITATION  MEAN FOR MULTIPLICATIVE CORRECTION
C                  =0.0; USE ADDITIVE CORRECTION (DEFAULT)
C        BIASRD - RADIATION FLUX CORRECTION TO CLOSE HEAT (W/M**2)
C                  =0.0; NO CORRECTION (DEFAULT)
C        RDMEAN - RADIATION FLUX MEAN FOR MULTIPLICATIVE CORRECTION
C                  =0.0; USE ADDITIVE CORRECTION (DEFAULT)
C        HMKS   - SCALE FACTOR FROM INPUT HAIR UNITS TO KG/KG.
C                  =1.0; NO SCALING (DEFAULT)
C        RMKS   - SCALE FACTOR FROM INPUT FLXR UNITS TO W/M**2 INTO OCEAN.
C                  =1.0; NO SCALING (DEFAULT)
C        PMKS   - SCALE FACTOR FROM INPUT PCIP UNITS TO M/S    INTO OCEAN.
C                  =1.0; NO SCALING (DEFAULT)
C        SEAMSK - CUTOVER FROM SEA TO LAND IN MASK (0.0 to 1.0, DEFAULT 0.5)
C
C     /AFFLAG/
C        IFFILE - MODEL FLUX INPUT FILE FLAG
C                    =3; ONE MONTHLY FLUX FILE (DEFAULT)
C                    =5; ONE FLUX FILE, ACTUAL FLUX DAY
C                         FLUX RECORDS CONTAIN THE ACTUAL 'FLUX TIME'
C        IFTYPE - INPUT FILE TYPE
C                    =5; Ta-Ha-Qr-Qp-Pc (DEFAULT)
C                    =4; Ta-Ha-Qr-Qp
C                    =3; Ta-Ha
C                    =2; Qr
C                    =1; Pc
C        INTERP - INTERPOLATION FLAG.
C                    =0; PIECEWISE BI-LINEAR
C                    =1; CUBIC SPLINE (DEFAULT)
C                    =2; PIECEWISE BESSEL
C                    =3; PIECEWISE BI-CUBIC
C        INTMSK - ALLOW FOR NATIVE GRID LAND-MASK
C                    =0; NO MASK (DEFAULT)
C                    =1; MASK IS >=SEAMSK OVER OCEAN, <SEAMSK OVER LAND
C                    =2; MASK IS <=SEAMSK OVER OCEAN, >SEAMSK OVER LAND
C        NGLOBE - NEAR-GLOBAL INPUT FLAG
C                    =0; STANDARD (E.G. FULLY-GLOBAL) INPUT (DEFAULT)
C                    =1; NEAR-GLOBAL INPUT
C                    =2; REGIONAL INPUT
C
C     IF A LAND-MASK IS USED, INTERP MUST BE 0 (PIECEWISE LINEAR).
C
C     NAMELIST /AFTIME/ IS PATTERNED AFTER /XXTIME/ SO THAT THE
C      MODEL''S STANDARD AWK-BASED RUN SCRIPT CUSTOMIZER CAN ALSO
C      WORK FOR THE FLUX GENERATION SCRIPT.  IN PARTICULAR, 'WSTART'
C      IS READ IN, BUT NOT USED.
C
C     BIASPC AND BIASRD FOLLOW HYCOM SIGN CONVENTIONS (+VE MEANS
C      GAIN BY OCEAN).
C
C 4)  INPUT:
C        ON UNIT  5:    NAMELIST /AFTITL/, /AFTIME/
C        ON CDF070:     netCDF NATIVE MASK FILE, SEE (5).
C        ON CDF071-099: netCDF NATIVE FLUX FILE, SEE (5).
C     OUTPUT:
C        ON UNIT 10:    UNFORMATTED MODEL    TAIR FILE, SEE (6).
C        ON UNIT 11:    UNFORMATTED MODEL    HAIR FILE, SEE (6).
C        ON UNIT 12:    UNFORMATTED MODEL    FLXR FILE, SEE (6).
C        ON UNIT 13:    UNFORMATTED MODEL    FLXP FILE, SEE (6).
C        ON UNIT 14:    UNFORMATTED MODEL    PCIP FILE, SEE (6).
C
C 5)  THE INPUT FLUX FIELDS, FROM ENVIRONMENT VARIABLES CDF071 TO CDF099,
C      ARE ON THEIR 'NATIVE' LAT-LON GRID.  THE NETCDF FILE LAYOUT IS:
C
C --- dimensions:
C ---         MT = UNLIMITED ;
C ---         Latitude =
C ---         Longitude =
C --- variables:
C ---         double MT(MT) ;
C ---                 MT:long_name = "time" ;
C ---                 MT:units = "days since 1900-12-31 00:00:00" ;
C ---                 MT:calendar = "gregorian" ;
C ---                 MT:axis = "T" ;
C ---                 MT:next_MT =
C ---         double Latitude(Latitude) ;
C ---                 Latitude:standard_name = "latitude" ;
C ---                 Latitude:units = "degrees_north" ;
C ---         double Longitude(Longitude) ;
C ---                 Longitude:standard_name = "longitude" ;
C ---                 Longitude:units = "degrees_east" ;
C ---         short precip(MT, Latitude, Longitude) ;
C ---                 precip:time_average_bounds =
C ---                 precip:long_name = " precipitation" ;
C ---                 precip:standard_name = "lwe_precipitation_rate" ;
C ---                 precip:units = "m/s" ;
C ---                 precip:scale_factor = 7.275958e-12f ;
C ---                 precip:add_offset = 2.384186e-07f ;
C ---
C --- // global attributes:
C ---                 :Conventions = "CF-1.0" ;
C ---                 :title =
C
C     THIS EXAMPLE IS FOR CNAME(5)="precip".  IT WILL BE CONVERTED FROM
C      INTEGER*2 TO REAL*4 USING scale_factor AND add_offset.
C     time_average_bounds INDICATES ANY TIME AVERAGING OF THE FIELD, AND
C      IS IGNORED HERE.  TYPICAL VALUES, ALWAYS IN DAYS, ARE:
C         0.0,0.0:       point values, no time averaging
C        -0.0416667,0.0: hrly average, reference time at end of interval
C        -0.5,0.5:       daily running average, centered in time
C        -15.25,15.25:   monthly average, centered in time,
C                         actual interval may vary by month
C     Longitude IS ASSUMED TO BE UNIFORM SPACED.
C     Latitude  IS ASSUMED TO BE UNIFORM SPACED OR GAUSSIAN.
C     MT:calendar IS ASSUMED TO BE "gregorian" OR "standard".
C     MT IS THE NUMBER OF FLUX RECORDS, AND MT(N) IS THE DATE OF FLUX
C      "RECORD" N W.R.T. JAN 1ST 1901 (days since 1900-12-31 00:00:00).
C     next_MT IS THE EXPECTED FLUX DATE FOR THE RECORD FOLLOWING THE
C      LAST RECORD IN THE FILE, I.E. FOLLOWING MT(MT).
C     BY CONVENTION, FLUX DATES BEFORE JAN 1ST 1905 (MT(1)<1462.0)
C      INDICATE A CLIMATOLOGY.  IN SUCH CASES THE CLIMATOLOGY'S LENGTH
C      (PERIOD) IS GIVEN BY next_MT-MT(1).
C     Ta-Ha-Qr-Qp-Pc FILES ALSO WORK FOR Ta-Ha-Qr-Qp AND Ta-Ha,
C      AND Ta-Ha-Qr-Qp FILES ALSO WORK FOR Ta-Ha.
C
C 6)  THE OUTPUT HEAT FLUXES ARE AT EVERY GRID POINT OF THE MODEL'S
C     'P' GRID.  ARRAY SIZE IS 'IDM' BY 'JDM'.
C
C 7)  SEVERAL STATISTICS ARE WRITTEN OUT IN ORDER TO CHECK THE
C      INTERPOLATION BETWEEN VARIOUS MACHINES.  MIN, MAX, MEAN AND RMS
C      OF THE ENTIRE BASIN ARE OUTPUT FOR (TAIR, HAIR, FLXR) FOR EACH
C      RECORD.  NOTE HOWEVER THAT THESE VALUES MAY NOT REPRESENT THE
C      STATISTICS OF THE FLUXS AS SEEN BY THE MODEL, IF THE INPUT FLUX
C      DATA HAS NON-REALISTIC VALUES OVER LAND.  IT IS UP TO THE USER
C      TO CHECK THE LOG FILES FOR CONSISTENCY BETWEEN MACHINES.
C
C 8)  ALAN J. WALLCRAFT,  PLANNING SYSTEMS INC.,  OCTOBER 1995.
C      BASED ON EARILER VERSIONS BY SEVERAL AUTHORS.
C      UPDATED FOR NETCDF INPUT IN JULY 2010.
C*
C**********
C
      EXTERNAL LINEAR, AVERMS,MINMAX
C
      INTEGER    MAXREC
      REAL*4     ZERO,RADIAN
      PARAMETER (MAXREC=190000)
      PARAMETER (ZERO=0.0, RADIAN=57.2957795)
C
      CHARACTER*80 CLINE
      INTEGER IFREC,IUREC,IWIX,NREC
      REAL*8  WDAY8,WDAY8I
      REAL*4  WDAY(MAXREC+1),WSTRT,WEND
      REAL*4  HMINA,HMINB,HMAXA,HMAXB
C
      LOGICAL LINTERP(5)
      INTEGER I,II,J,JJ,KREC
      REAL*4  BIASRD4,BIASPC4,FLXPOLD,SCALEPC,SCALERD
      REAL*4  XFD,YFD,DXD,DYD
      REAL*4  FDY,WYR,XLIN,XFDX,XOV,
     +        XMIN,XMAX,XAVE,XRMS
C
C --- MODEL ARRAYS.
C
      CALL XCSPMD  !define idm,jdm
      ALLOCATE(  MSK(IDM,JDM),
     +          PLON(IDM,JDM),
     +          PLAT(IDM,JDM),
     +           XAF(IDM,JDM),
     +           YAF(IDM,JDM),
     +           TAM(IDM,JDM),
     +           HAM(IDM,JDM),
     +           FRM(IDM,JDM),
     +           FPM(IDM,JDM),
     +           PCM(IDM,JDM) )
C
C     NAMELIST INPUT.
C
      CALL ZHOPEN(6, 'FORMATTED', 'UNKNOWN', 0)
C
      CTITLE = ' '
      CNAME(1) = 'airtmp'
      CNAME(2) = 'vapmix'
      CNAME(3) = 'radflx'
      CNAME(4) = 'shwflx'
      CNAME(5) = 'precip'
      WRITE(6,*) 'READING /AFTITL/'
      CALL ZHFLSH(6)
      READ( 5,AFTITL)
      WRITE(6,AFTITL)
C
      FSTART = 0.0
      WSTART = 0.0
      TSTART = 0.0
      TMAX   = 0.0
      BIASPC = 0.0
      BIASRD = 0.0
      PCMEAN = 0.0
      RDMEAN = 0.0
      HMKS   = 1.0
      RMKS   = 1.0
      PMKS   = 1.0
      SEAMSK = 0.5
      WRITE(6,*) 'READING /AFTIME/'
      CALL ZHFLSH(6)
      READ( 5,AFTIME)
      WRITE(6,AFTIME)
      WRITE(6,*)
      CALL ZHFLSH(6)
C
      IFFILE = 3
      IFTYPE = 5
      INTERP = 1
      INTMSK = 0
      NGLOBE = 0
      JPR    = 8
      WRITE(6,*) 'READING /AFFLAG/'
      CALL ZHFLSH(6)
      READ( 5,AFFLAG)
      WRITE(6,AFFLAG)
      WRITE(6,*)
      CALL ZHFLSH(6)
C
      IF     (INTMSK.NE.0 .AND. INTERP.NE.0) THEN
        WRITE(6,'(/ a /)')
     &    'error - INTERP must be xero if INTMSK is non-zero'
        CALL ZHFLSH(6)
        STOP
      ENDIF
C
C     GRID INPUT.
C
      CALL ZAIOST
C
      CALL ZHOPNC(21, 'regional.grid.b', 'FORMATTED', 'OLD', 0)
      CALL ZAIOPF('regional.grid.a', 'OLD', 21)
C
      READ(21,*) ! skip idm
      READ(21,*) ! skip jdm
      READ(21,*) ! skip mapflg
      READ(21,'(A)') CLINE
      I = INDEX(CLINE,'=')
      READ (CLINE(I+1:),*)   HMINB,HMAXB
      CALL ZAIORD(PLON,MSK,.FALSE., HMINA,HMAXA, 21)
      IF     (ABS(HMINA-HMINB).GT.ABS(HMINB)*1.E-4 .OR.
     &        ABS(HMAXA-HMAXB).GT.ABS(HMAXB)*1.E-4     ) THEN
        WRITE(6,'(/ a / a,1p3e14.6 / a,1p3e14.6 /)')
     &    'error - .a and .b grid files not consistent (plon):',
     &    '.a,.b min = ',HMINA,HMINB,HMINA-HMINB,
     &    '.a,.b max = ',HMAXA,HMAXB,HMAXA-HMAXB
        CALL ZHFLSH(6)
        STOP
      ENDIF
C
      READ(21,'(A)') CLINE
      I = INDEX(CLINE,'=')
      READ (CLINE(I+1:),*)   HMINB,HMAXB
      CALL ZAIORD(PLAT,MSK,.FALSE., HMINA,HMAXA, 21)
      IF     (ABS(HMINA-HMINB).GT.ABS(HMINB)*1.E-4 .OR.
     &        ABS(HMAXA-HMAXB).GT.ABS(HMAXB)*1.E-4     ) THEN
        WRITE(6,'(/ a / a,1p3e14.6 / a,1p3e14.6 /)')
     &    'error - .a and .b grid files not consistent (plat):',
     &    '.a,.b min = ',HMINA,HMINB,HMINA-HMINB,
     &    '.a,.b max = ',HMAXA,HMAXB,HMAXA-HMAXB
        CALL ZHFLSH(6)
        STOP
      ENDIF
C
      CLOSE(UNIT=21)
      CALL ZAIOCL(21)
C
C     CHECK FOR AN ARCTIC BIPOLAR PATCH
C
      LARCTIC = PLON(3*IDM/4+1,JDM-1).EQ.PLON(IDM/4,JDM) .AND.
     &          PLAT(3*IDM/4+1,JDM-1).EQ.PLAT(IDM/4,JDM)
C
      WRITE(6,*)
      WRITE(6,*) 'LARCTIC = ',LARCTIC
      WRITE(6,*)
      CALL ZHFLSH(6)
C
C     INITIALIZE OUTPUT.
C
      WRITE(6,6000) 'OUTPUT:',CTITLE
      CALL ZHFLSH(6)
C
      CALL ZAIOPN('NEW', 10)
      CALL ZAIOPN('NEW', 11)
      CALL ZAIOPN('NEW', 12)
      CALL ZAIOPN('NEW', 13)
      CALL ZAIOPN('NEW', 14)
C
      CALL ZHOPEN(10, 'FORMATTED', 'NEW', 0)
      CALL ZHOPEN(11, 'FORMATTED', 'NEW', 0)
      CALL ZHOPEN(12, 'FORMATTED', 'NEW', 0)
      CALL ZHOPEN(13, 'FORMATTED', 'NEW', 0)
      CALL ZHOPEN(14, 'FORMATTED', 'NEW', 0)
C
      PREAMBL(1) = CTITLE
      PREAMBL(2) = ' '
      PREAMBL(3) = ' '
      PREAMBL(4) = ' '
      WRITE(PREAMBL(5),'(A,2I5)')
     +        'i/jdm =',
     +       IDM,JDM
      WRITE(10,4101) PREAMBL
      WRITE(11,4101) PREAMBL
C
      IF     (IFTYPE.NE.3) THEN
      IF     (BIASRD.EQ.0.0) THEN
        PREAMBL(2) = 'No radiation flux correction'
      ELSEIF (RDMEAN.EQ.0.0) THEN
        WRITE(PREAMBL(2),'(A,F7.1,A)')
     +        'Shortwave correction is',
     +         BIASRD,' w/m**2 into the ocean'
      ELSE
        WRITE(PREAMBL(2),'(A,F7.1,A,F7.1,A)')
     +        'Shortwave correction is',
     +         BIASRD,' w/m**2 into the ocean (w.r.t.',
     +         RDMEAN,' w/m**2)'
      ENDIF
      WRITE(12,4101) PREAMBL
      WRITE(13,4101) PREAMBL
C
      IF     (IFTYPE.EQ.4 .OR. IFTYPE.EQ.2) THEN
        PREAMBL(2) = 'Zero precipitation'
      ELSEIF (BIASPC.EQ.0.0) THEN
        PREAMBL(2) = 'No precipitation correction'
      ELSEIF (PCMEAN.EQ.0.0) THEN
        WRITE(PREAMBL(2),'(A,F12.8,A)')
     +        'Precipitation correction  is',1000.0*BIASPC,
     +        ' mm/s into the ocean'
      ELSE
        WRITE(PREAMBL(2),'(A,F12.8,A,F12.8,A)')
     +    'Precipitation correction  is',
     +     1000.0*BIASPC,' mm/s into the ocean (w.r.t.',
     +     1000.0*PCMEAN,' mm/s)'
      ENDIF
      WRITE(14,4101) PREAMBL
C
      IF     (BIASRD.EQ.0.0) THEN
        PREAMBL(2) = 'No radiation flux correction'
      ELSEIF (RDMEAN.EQ.0.0) THEN
        WRITE(PREAMBL(2),'(A,F7.1,A)')
     +        'Shortwave correction is',
     +         BIASRD,' w/m**2 into the ocean'
      ELSE
        WRITE(PREAMBL(2),'(A,F7.1,A,F7.1,A)')
     +        'Shortwave correction is',
     +         BIASRD,' w/m**2 into the ocean (w.r.t.',
     +         RDMEAN,' w/m**2)'
      ENDIF
      IF     (IFTYPE.EQ.4 .OR. IFTYPE.EQ.2) THEN
        PREAMBL(3) = 'Zero precipitation'
      ELSEIF (BIASPC.EQ.0.0) THEN
        PREAMBL(3) = 'No precipitation correction'
      ELSEIF (PCMEAN.EQ.0.0) THEN
        WRITE(PREAMBL(3),'(A,F12.8,A)')
     +        'Precipitation correction  is',1000.0*BIASPC,
     +        ' mm/s into the ocean'
      ELSE
        WRITE(PREAMBL(3),'(A,F12.8,A,F12.8,A)')
     +    'Precipitation correction  is',
     +     1000.0*BIASPC,' mm/s into the ocean (w.r.t.',
     +     1000.0*PCMEAN,' mm/s)'
      ENDIF
      WRITE(6,*)
      WRITE(6, 4101) PREAMBL
      WRITE(6,*)
      ENDIF !iftype/=3
C
C     NATIVE GRID ARRAYS.
C
      CALL FREADI(IWI,JWI)
C
      ALLOCATE(   WLON(IWI+1)   )
      ALLOCATE(   WLAT(    JWI) )
      ALLOCATE(  FLDI2(IWI,JWI) )
      ALLOCATE( TAIRIN(IWI,JWI),
     +          HAIRIN(IWI,JWI),
     +          FLXRIN(IWI,JWI),
     +          FLXPIN(IWI,JWI),
     +          PCIPIN(IWI,JWI) )
      IF     (INTMSK.NE.0) THEN
        ALLOCATE( MASKIN(IWI,JWI) )
      ENDIF
C
      ALLOCATE( TAIRI(IWI+5,JWI+5),
     +          HAIRI(IWI+5,JWI+5),
     +          FLXRI(IWI+5,JWI+5),
     +          FLXPI(IWI+5,JWI+5),
     +          PCIPI(IWI+5,JWI+5) )  !+5 needed for bicubc
      IF     (INTERP.EQ.1) THEN
        ALLOCATE( FXI(IWI+4),
     +            FYI(JWI+4),
     +            WK(3*(IWI+JWI+8)+1) )
        ALLOCATE( WTAIR3(IWI+4,JWI+4,3) )
      ENDIF
C
C     INITIALIZE FOR NATIVE FLUXS.
C
      IF     (IFFILE.EQ.3) THEN
        WSTRT = -1.0
        WEND  = -1.0
      ELSE
        WSTRT = FSTART
        WEND  = FSTART + (TMAX - TSTART)
      ENDIF
      CALL FREAD0(NREC,IUREC,IFREC,WDAY,MAXREC, WSTRT,WEND,
     +            IWI,JWI,WLON,WLAT)
C
      IF     (IFFILE.EQ.3 .AND. NREC.GT.12) THEN
        WRITE(6,9200)
        CALL ZHFLSH(6)
        STOP
      ENDIF
C
C     DEFINE THE GRID COORDINATES.
C
      WLON(IWI+1) = WLON(IWI) + (WLON(2)-WLON(1))
      IF     (WLON(IWI+1)-WLON(1).GE.359.9) THEN
        IF     (ABS(WLON(IWI+1)-WLON(1)-360.0) .GT. 0.01) THEN
          WRITE(6,9050)
          CALL ZHFLSH(6)
          STOP
        ENDIF
        IWIX = IWI + 1
        IBD(1) = 3
        IBD(2) = 3
        IBD(3) = 2
        IBD(4) = 2
      ELSE
        IWIX = IWI
        IBD(1) = 2
        IBD(2) = 2
        IBD(3) = 2
        IBD(4) = 2
      ENDIF
C
C     CONVERT HYCOM LON,LAT TO CLIMATOLOGY ARRAY COORDS.
C     ASSUMES A UNIFORM LONGITUDINAL GRID SPACING.
C
      DXIN  = WLON(2) - WLON(1)
      XFIN  = WLON(1)
      XLIN  = WLON(IWIX)
      XAMIN = 2*IWI
      XAMAX = 0
      DO J= 1,JDM
        DO I= 1,IDM
          XOV = PLON(I,J)
          IF     (XOV.LT.XFIN) THEN
            XOV = XOV + 360.0
          ELSEIF (XOV.GE.XLIN) THEN
            XOV = XOV - 360.0
          ENDIF
C
          XAF(I,J) = 3.0 + (XOV - XFIN)/DXIN
C
*         IF     (MOD(J,100).EQ.1 .OR. J.EQ.JDM) THEN
*           IF     (MOD(I,MAX(10,IDM/20)).EQ.1 .OR. I.EQ.IDM) THEN
*             WRITE(6,'("I,J,LONV,XAF =",2I5,2F10.3)') I,J,XOV,XAF(I,J)
*             CALL ZHFLSH(6)
*           ENDIF
*         ENDIF
          XAMIN  = MIN( XAMIN, XAF(I,J) )
          XAMAX  = MAX( XAMAX, XAF(I,J) )
        ENDDO
      ENDDO
C
      IF     (NGLOBE.NE.2) THEN
        YFMIN = WLAT(1)
        YFMAX = WLAT(JWI)
      ELSE  ! non-global native grid, inactivate YFMIN,YFMAX
        YFMIN = -90.0
        YFMAX =  90.0
      ENDIF
      YAMIN = 2*JWI
      YAMAX = 0
      DO I= 1,IDM
        JJ = 1
        DO J= 1,JDM
          IF     (I.GT.1 .AND. ABS(PLAT(I,J)-PLAT(1,J)).LT.0.005) THEN
            YAF(I,J) = YAF(1,J)
            CYCLE
          ENDIF
C
          PLATIJ = MIN(YFMAX,MAX(YFMIN,PLAT(I,J)))
          IF     (PLATIJ.LE.YFMIN) THEN
            YAF(I,J) = 3.0
          ELSEIF (PLATIJ.GE.YFMAX) THEN
            YAF(I,J) = 1.999 + JWI
          ELSE
            IF     (WLAT(JJ).LE.PLATIJ  .AND.
     &              PLATIJ.LT.WLAT(JJ+1)     ) THEN  !try jj from last loop
              YAF(I,J) = 2.0 + JJ + (PLATIJ    -WLAT(JJ))/
     &                              (WLAT(JJ+1)-WLAT(JJ))
            ELSE
              YAF(I,J) = 0.0
              DO JJ= 1,JWI-1
                IF     (WLAT(JJ).LE.PLATIJ  .AND.
     &                  PLATIJ.LT.WLAT(JJ+1)     ) THEN
                  YAF(I,J) = 2.0 + JJ + (PLATIJ    -WLAT(JJ))/
     &                                  (WLAT(JJ+1)-WLAT(JJ))
                  EXIT
                ENDIF
              ENDDO !jj
            ENDIF !try old jj:else
          ENDIF
C
*           IF     (MOD(I,100).EQ.1 .OR. I.EQ.IDM) THEN
            IF     (I.EQ.1) THEN
              IF     (MOD(J,MAX(10,JDM/20)).EQ.1 .OR. J.EQ.JDM) THEN
                WRITE(6,'("I,J,LATU,YAF =",2I5,2F10.3)')
     +                     I,J,PLAT(I,J),YAF(I,J)
                CALL ZHFLSH(6)
              ENDIF
            ENDIF
*           IF     (MOD(I,1000).EQ.1 .OR. I.EQ.IDM) THEN
*             IF     (YAF(I,J).GT.JWI-1) THEN
*               WRITE(6,'("I,J,LON,LAT,X,YAF =",2I5,4F10.3)')
*    +                     I,J,PLON(I,J),PLAT(I,J),XAF(I,J),YAF(I,J)
*               CALL ZHFLSH(6)
*             ENDIF
*           ENDIF
C
            YAMIN  = MIN( YAMIN, YAF(I,J) )
            YAMAX  = MAX( YAMAX, YAF(I,J) )
          ENDDO !j
        ENDDO !i
C
      WRITE(6,6200) XAMIN,XAMAX,YAMIN,YAMAX
      CALL ZHFLSH(6)
C
C     CHECK THAT THE INTERPOLATION IS 'SAFE',
C
      IF     (INT(XAMIN).LT.3 .OR. INT(XAMAX).GT.IWI+2 .OR.
     +        INT(YAMIN).LT.3 .OR. INT(YAMAX).GT.JWI+2     ) THEN
        WRITE(6,9150)
        CALL ZHFLSH(6)
        STOP
      ENDIF
C
C     INPUT THE MASK.
C
      IF     (INTMSK.NE.0) THEN
        CALL FREADM(PCIPIN,FLDI2,IWI,JWI,70)
        IF     (INTMSK.EQ.1) THEN
          DO J= 1,JWI
            DO I= 1,IWI
              IF     (PCIPIN(I,J).GE.SEAMSK) THEN
                MASKIN(I,J) = 1  !ocean
              ELSE
                MASKIN(I,J) = 0  !land
              ENDIF
            ENDDO !i
          ENDDO !j
        ELSE !intmsk=2
          DO J= 1,JWI
            DO I= 1,IWI
              IF     (PCIPIN(I,J).LE.SEAMSK) THEN
                MASKIN(I,J) = 1  !ocean
              ELSE
                MASKIN(I,J) = 0  !land
              ENDIF
            ENDDO !i
          ENDDO !j
        ENDIF !intmsk
      ENDIF !MASK
C
C     PROCESS ALL THE FLUX RECORDS.
C
      LINTERP(1) = IFTYPE.NE.1 .AND. IFTYPE.NE.2  ! TAIR
      LINTERP(2) = IFTYPE.NE.1 .AND. IFTYPE.NE.2  ! HAIR
      LINTERP(3) = IFTYPE.NE.1 .AND. IFTYPE.NE.3  ! FLXR
      LINTERP(4) = IFTYPE.NE.1 .AND. IFTYPE.NE.2
     &                         .AND. IFTYPE.NE.3  ! FLXP
      LINTERP(5) = IFTYPE.NE.4 .AND. IFTYPE.NE.2
     &                         .AND. IFTYPE.NE.3  ! PCIP
      TAM = 0.0
      HAM = 0.0
      FRM = 0.0
      FPM = 0.0
      PCM = 0.0
C
      DO 810 KREC= 1,NREC
C
C       READ THE INPUT FLUXES.
C
        CALL FREAD1(TAIRIN,HAIRIN,FLXRIN,FLXPIN,PCIPIN,IWI,JWI,
     +              FLDI2, CNAME,
     +              IUREC,IFREC, KREC,WDAY(KREC),IFTYPE)
        IF     (INTMSK.NE.0) THEN
          IF     (IFTYPE.EQ.1) THEN
            CALL LANDFILL1(PCIPIN,MASKIN,IWI,JWI, 99,
     +                     IWIX.GT.IWI,.FALSE., KREC.EQ.1)
          ELSEIF (IFTYPE.EQ.2) THEN
            CALL LANDFILL1(FLXRIN,MASKIN,IWI,JWI, 99,
     +                     IWIX.GT.IWI,.FALSE., KREC.EQ.1)
          ELSEIF (IFTYPE.EQ.3) THEN
            CALL LANDFILL2(TAIRIN,HAIRIN,MASKIN,IWI,JWI, 99,
     +                     IWIX.GT.IWI,.FALSE., KREC.EQ.1)
          ELSE
            CALL LANDFILL5(TAIRIN,HAIRIN,FLXRIN,FLXPIN,
     +                     PCIPIN,MASKIN,IWI,JWI, 99,
     +                     IWIX.GT.IWI,.FALSE., KREC.EQ.1)
          ENDIF
        ENDIF
C
C       COPY INTO THE (LARGER) INTERPOLATION ARRAYS.
C       PMKS*PCIPIN AND RMKS*FLXPIN (AND FLXPI) MUST BE NON-NEGATIVE.
C
        IF     (PCMEAN.EQ.0.0) THEN
          BIASPC4 = BIASPC
        ELSE
          SCALEPC = 1.0 + BIASPC/PCMEAN
        ENDIF
        IF     (RDMEAN.EQ.0.0) THEN
          BIASRD4 = BIASRD
        ELSE
          SCALERD = 1.0 + BIASRD/RDMEAN
        ENDIF
        DO J= 1,JWI
          DO I= 1,IWI
            TAIRI(I+2,J+2) =          TAIRIN(I,J)
            HAIRI(I+2,J+2) =     HMKS*HAIRIN(I,J)
            IF     (RDMEAN.EQ.0.0) THEN
              FLXPI(I+2,J+2) = MAX(RMKS*FLXPIN(I,J) + BIASRD4, ZERO)
              FLXRI(I+2,J+2) =     RMKS*FLXRIN(I,J) + BIASRD4
            ELSE
              FLXPOLD        = MAX(RMKS*FLXPIN(I,J), ZERO)
              FLXPI(I+2,J+2) =     FLXPOLD*SCALERD
              FLXRI(I+2,J+2) =     RMKS*FLXRIN(I,J) +
     +                             FLXPOLD*SCALERD  -
     +                             FLXPOLD
            ENDIF
            IF     (PCMEAN.EQ.0.0) THEN
              PCIPI(I+2,J+2) = MAX(PMKS*PCIPIN(I,J) + BIASPC4, BIASPC4)
            ELSE
              PCIPI(I+2,J+2) = MAX(PMKS*PCIPIN(I,J)*SCALEPC,   ZERO)
            ENDIF
          ENDDO
        ENDDO
C
C       FILL IN THE PADDING AREA AS NECESSARY.
C
        IF     (INT(XAMAX).GE.IWI+1) THEN  !may need iwi+3 and perhaps iwi+4/5
          IF     (IWIX.GT.IWI) THEN
            DO 320 J= 3,JWI+2
              TAIRI(IWI+3,J) = TAIRI(3,J)
              TAIRI(IWI+4,J) = TAIRI(4,J)
              TAIRI(IWI+5,J) = TAIRI(5,J)
              HAIRI(IWI+3,J) = HAIRI(3,J)
              HAIRI(IWI+4,J) = HAIRI(4,J)
              HAIRI(IWI+5,J) = HAIRI(5,J)
              FLXRI(IWI+3,J) = FLXRI(3,J)
              FLXRI(IWI+4,J) = FLXRI(4,J)
              FLXRI(IWI+5,J) = FLXRI(5,J)
              FLXPI(IWI+3,J) = FLXPI(3,J)
              FLXPI(IWI+4,J) = FLXPI(4,J)
              FLXPI(IWI+5,J) = FLXPI(5,J)
              PCIPI(IWI+3,J) = PCIPI(3,J)
              PCIPI(IWI+4,J) = PCIPI(4,J)
              PCIPI(IWI+5,J) = PCIPI(5,J)
  320       CONTINUE
          ELSE
            DO 325 J= 3,JWI+2
              TAIRI(IWI+3,J) = 2.0*TAIRI(IWI+2,J) - TAIRI(IWI+1,J)
              TAIRI(IWI+4,J) = 2.0*TAIRI(IWI+3,J) - TAIRI(IWI+2,J)
              TAIRI(IWI+5,J) = 2.0*TAIRI(IWI+4,J) - TAIRI(IWI+3,J)
              HAIRI(IWI+3,J) = 2.0*HAIRI(IWI+2,J) - HAIRI(IWI+1,J)
              HAIRI(IWI+4,J) = 2.0*HAIRI(IWI+3,J) - HAIRI(IWI+2,J)
              HAIRI(IWI+5,J) = 2.0*HAIRI(IWI+4,J) - HAIRI(IWI+3,J)
              FLXRI(IWI+3,J) = 2.0*FLXRI(IWI+2,J) - FLXRI(IWI+1,J)
              FLXRI(IWI+4,J) = 2.0*FLXRI(IWI+3,J) - FLXRI(IWI+2,J)
              FLXRI(IWI+5,J) = 2.0*FLXRI(IWI+4,J) - FLXRI(IWI+3,J)
              FLXPI(IWI+3,J) = 2.0*FLXPI(IWI+2,J) - FLXPI(IWI+1,J)
              FLXPI(IWI+4,J) = 2.0*FLXPI(IWI+3,J) - FLXPI(IWI+2,J)
              FLXPI(IWI+5,J) = 2.0*FLXPI(IWI+4,J) - FLXPI(IWI+3,J)
              PCIPI(IWI+3,J) = 2.0*PCIPI(IWI+2,J) - PCIPI(IWI+1,J)
              PCIPI(IWI+4,J) = 2.0*PCIPI(IWI+3,J) - PCIPI(IWI+2,J)
              PCIPI(IWI+5,J) = 2.0*PCIPI(IWI+4,J) - PCIPI(IWI+3,J)
  325       CONTINUE
          ENDIF
        ENDIF
        IF     (INT(XAMIN).LE.4) THEN  !may need 2 and perhaps 1
          IF     (IWIX.GT.IWI) THEN
            DO 330 J= 3,JWI+2
              TAIRI(1,J) = TAIRI(IWI+1,J)
              TAIRI(2,J) = TAIRI(IWI+2,J)
              HAIRI(1,J) = HAIRI(IWI+1,J)
              HAIRI(2,J) = HAIRI(IWI+2,J)
              FLXRI(1,J) = FLXRI(IWI+1,J)
              FLXRI(2,J) = FLXRI(IWI+2,J)
              FLXPI(1,J) = FLXPI(IWI+1,J)
              FLXPI(2,J) = FLXPI(IWI+2,J)
              PCIPI(1,J) = PCIPI(IWI+1,J)
              PCIPI(2,J) = PCIPI(IWI+2,J)
  330       CONTINUE
          ELSE
            DO 335 J= 3,JWI+2
              TAIRI(1,J) = 3.0*TAIRI(3,J) - 2.0*TAIRI(4,J)
              TAIRI(2,J) = 2.0*TAIRI(3,J) -     TAIRI(4,J)
              HAIRI(1,J) = 3.0*HAIRI(3,J) - 2.0*HAIRI(4,J)
              HAIRI(2,J) = 2.0*HAIRI(3,J) -     HAIRI(4,J)
              FLXRI(1,J) = 3.0*FLXRI(3,J) - 2.0*FLXRI(4,J)
              FLXRI(2,J) = 2.0*FLXRI(3,J) -     FLXRI(4,J)
              FLXPI(1,J) = 3.0*FLXPI(3,J) - 2.0*FLXPI(4,J)
              FLXPI(2,J) = 2.0*FLXPI(3,J) -     FLXPI(4,J)
              PCIPI(1,J) = 3.0*PCIPI(3,J) - 2.0*PCIPI(4,J)
              PCIPI(2,J) = 2.0*PCIPI(3,J) -     PCIPI(4,J)
  335       CONTINUE
          ENDIF
        ENDIF
        IF     (INT(YAMAX).GE.JWI+1) THEN !may need jwi+3 and perhaps jwi+4/5
          IF     (IWIX.GT.IWI) THEN  !global grid
            IF     (DYIN.EQ.0.0 .OR.
     +              ABS(YFIN+JWI*DYIN-90.0).LE.0.1*DYIN) THEN
C ---         JWI+3 = 90N
              DO I= 1,IWI+5
                II = MOD(I-3+IWI/2+IWI,IWI)+3
                TAIRI(I,JWI+5) =      TAIRI(II,JWI+1)
                TAIRI(I,JWI+4) =      TAIRI(II,JWI+2)
                TAIRI(I,JWI+3) = 0.5*(TAIRI(I, JWI+2)+TAIRI(I,JWI+4))
                HAIRI(I,JWI+5) =      HAIRI(II,JWI+1)
                HAIRI(I,JWI+4) =      HAIRI(II,JWI+2)
                HAIRI(I,JWI+3) = 0.5*(HAIRI(I, JWI+2)+HAIRI(I,JWI+4))
                FLXRI(I,JWI+5) =      FLXRI(II,JWI+1)
                FLXRI(I,JWI+4) =      FLXRI(II,JWI+2)
                FLXRI(I,JWI+3) = 0.5*(FLXRI(I, JWI+2)+FLXRI(I,JWI+4))
                FLXPI(I,JWI+5) =      FLXPI(II,JWI+1)
                FLXPI(I,JWI+4) =      FLXPI(II,JWI+2)
                FLXPI(I,JWI+3) = 0.5*(FLXPI(I, JWI+2)+FLXPI(I,JWI+4))
                PCIPI(I,JWI+5) =      PCIPI(II,JWI+1)
                PCIPI(I,JWI+4) =      PCIPI(II,JWI+2)
                PCIPI(I,JWI+3) = 0.5*(PCIPI(I, JWI+2)+PCIPI(I,JWI+4))
*                WRITE(6,'(A,2I5,4F10.3)')
*    +            'I,II,FLXR = ',I,II,
*    +            FLXRI(I,JWI+1),FLXRI(I,JWI+2),
*    +            FLXRI(I,JWI+3),FLXRI(I,JWI+4)
              ENDDO !i
            ELSEIF (ABS(YFIN+(JWI-1)*DYIN-90.0).LE.0.1*DYIN) THEN
C ---         JWI+2 = 90N
              DO I= 1,IWI+5
                II = MOD(I-3+IWI/2+IWI,IWI)+3
                TAIRI(I,JWI+3) = TAIRI(II,JWI+1)
                TAIRI(I,JWI+4) = TAIRI(II,JWI  )
                TAIRI(I,JWI+5) = TAIRI(II,JWI-1)
                HAIRI(I,JWI+3) = HAIRI(II,JWI+1)
                HAIRI(I,JWI+4) = HAIRI(II,JWI  )
                HAIRI(I,JWI+5) = HAIRI(II,JWI-1)
                FLXRI(I,JWI+3) = FLXRI(II,JWI+1)
                FLXRI(I,JWI+4) = FLXRI(II,JWI  )
                FLXRI(I,JWI+5) = FLXRI(II,JWI-1)
                FLXPI(I,JWI+3) = FLXPI(II,JWI+1)
                FLXPI(I,JWI+4) = FLXPI(II,JWI  )
                FLXPI(I,JWI+5) = FLXPI(II,JWI-1)
                PCIPI(I,JWI+3) = PCIPI(II,JWI+1)
                PCIPI(I,JWI+4) = PCIPI(II,JWI  )
                PCIPI(I,JWI+5) = PCIPI(II,JWI-1)
*                WRITE(6,'(A,2I5,4F10.3)')
*    +            'I,II,FLXR = ',I,II,
*    +            FLXRI(I,JWI+1),FLXRI(I,JWI+2),
*    +            FLXRI(I,JWI+3),FLXRI(I,JWI+4)
              ENDDO !i
            ELSE
              DO I= 1,IWI+5
                TAIRI(I,JWI+3) = 2.0*TAIRI(I,JWI+2) - TAIRI(I,JWI+1)
                TAIRI(I,JWI+4) = 2.0*TAIRI(I,JWI+3) - TAIRI(I,JWI+2)
                TAIRI(I,JWI+5) = 2.0*TAIRI(I,JWI+4) - TAIRI(I,JWI+3)
                HAIRI(I,JWI+3) = 2.0*HAIRI(I,JWI+2) - HAIRI(I,JWI+1)
                HAIRI(I,JWI+4) = 2.0*HAIRI(I,JWI+3) - HAIRI(I,JWI+2)
                HAIRI(I,JWI+5) = 2.0*HAIRI(I,JWI+4) - HAIRI(I,JWI+3)
                FLXRI(I,JWI+3) = 2.0*FLXRI(I,JWI+2) - FLXRI(I,JWI+1)
                FLXRI(I,JWI+4) = 2.0*FLXRI(I,JWI+3) - FLXRI(I,JWI+2)
                FLXRI(I,JWI+5) = 2.0*FLXRI(I,JWI+4) - FLXRI(I,JWI+3)
                FLXPI(I,JWI+3) = 2.0*FLXPI(I,JWI+2) - FLXPI(I,JWI+1)
                FLXPI(I,JWI+4) = 2.0*FLXPI(I,JWI+3) - FLXPI(I,JWI+2)
                FLXPI(I,JWI+5) = 2.0*FLXPI(I,JWI+4) - FLXPI(I,JWI+3)
                PCIPI(I,JWI+3) = 2.0*PCIPI(I,JWI+2) - PCIPI(I,JWI+1)
                PCIPI(I,JWI+4) = 2.0*PCIPI(I,JWI+3) - PCIPI(I,JWI+2)
                PCIPI(I,JWI+5) = 2.0*PCIPI(I,JWI+4) - PCIPI(I,JWI+3)
*                WRITE(6,'(A,2I5,4F10.3)')
*    +            'I,II,FLXR = ',I,II,
*    +            FLXRI(I,JWI+1),FLXRI(I,JWI+2),
*    +            FLXRI(I,JWI+3),FLXRI(I,JWI+4)
              ENDDO !i
            ENDIF
          ELSE  !non-global grid
            DO 345 I= 1,IWI+5
              TAIRI(I,JWI+3) = 2.0*TAIRI(I,JWI+2) - TAIRI(I,JWI+1)
              TAIRI(I,JWI+4) = 2.0*TAIRI(I,JWI+3) - TAIRI(I,JWI+2)
              TAIRI(I,JWI+5) = 2.0*TAIRI(I,JWI+4) - TAIRI(I,JWI+3)
              HAIRI(I,JWI+3) = 2.0*HAIRI(I,JWI+2) - HAIRI(I,JWI+1)
              HAIRI(I,JWI+4) = 2.0*HAIRI(I,JWI+3) - HAIRI(I,JWI+2)
              HAIRI(I,JWI+5) = 2.0*HAIRI(I,JWI+4) - HAIRI(I,JWI+3)
              FLXRI(I,JWI+3) = 2.0*FLXRI(I,JWI+2) - FLXRI(I,JWI+1)
              FLXRI(I,JWI+4) = 2.0*FLXRI(I,JWI+3) - FLXRI(I,JWI+2)
              FLXRI(I,JWI+5) = 2.0*FLXRI(I,JWI+4) - FLXRI(I,JWI+3)
              FLXPI(I,JWI+3) = 2.0*FLXPI(I,JWI+2) - FLXPI(I,JWI+1)
              FLXPI(I,JWI+4) = 2.0*FLXPI(I,JWI+3) - FLXPI(I,JWI+2)
              FLXPI(I,JWI+5) = 2.0*FLXPI(I,JWI+4) - FLXPI(I,JWI+3)
              PCIPI(I,JWI+3) = 2.0*PCIPI(I,JWI+2) - PCIPI(I,JWI+1)
              PCIPI(I,JWI+4) = 2.0*PCIPI(I,JWI+3) - PCIPI(I,JWI+2)
              PCIPI(I,JWI+5) = 2.0*PCIPI(I,JWI+4) - PCIPI(I,JWI+3)
  345       CONTINUE
          ENDIF
        ENDIF
        IF     (INT(YAMIN).LE.4) THEN  !may need 2 and perhaps 1
          IF     (IWIX.GT.IWI) THEN  !global grid
            IF     (DYIN.EQ.0.0 .OR.
     +              ABS(YFIN-DYIN+90.0).LE.0.1*DYIN) THEN
C ---         2 = 90S
              DO I= 1,IWI+5
                II = MOD(I-3+IWI/2+IWI,IWI)+3
                TAIRI(I,1) =      TAIRI(II,3)
                TAIRI(I,2) = 0.5*(TAIRI(I, 1)+TAIRI(I,3))
                HAIRI(I,1) =      HAIRI(II,3)
                HAIRI(I,2) = 0.5*(HAIRI(I, 1)+HAIRI(I,3))
                FLXRI(I,1) =      FLXRI(II,3)
                FLXRI(I,2) = 0.5*(FLXRI(I, 1)+FLXRI(I,3))
                FLXPI(I,1) =      FLXPI(II,3)
                FLXPI(I,2) = 0.5*(FLXPI(I, 1)+FLXPI(I,3))
                PCIPI(I,1) =      PCIPI(II,3)
                PCIPI(I,2) = 0.5*(PCIPI(I, 1)+PCIPI(I,3))
              ENDDO !i
            ELSEIF (ABS(YFIN+90.0).LE.0.1*DYIN) THEN
C ---         3 = 90S
              DO I= 1,IWI+5
                II = MOD(I-3+IWI/2+IWI,IWI)+3
                TAIRI(I,1) = TAIRI(II,5)
                TAIRI(I,2) = TAIRI(II,4)
                HAIRI(I,1) = HAIRI(II,5)
                HAIRI(I,2) = HAIRI(II,4)
                FLXRI(I,1) = FLXRI(II,5)
                FLXRI(I,2) = FLXRI(II,4)
                FLXPI(I,1) = FLXPI(II,5)
                FLXPI(I,2) = FLXPI(II,4)
                PCIPI(I,1) = PCIPI(II,5)
                PCIPI(I,2) = PCIPI(II,4)
              ENDDO !i
            ELSE
              DO I= 1,IWI+5
                TAIRI(I,1) = 3.0*TAIRI(I,3) - 2.0*TAIRI(I,4)
                TAIRI(I,2) = 2.0*TAIRI(I,3) -     TAIRI(I,4)
                HAIRI(I,1) = 3.0*HAIRI(I,3) - 2.0*HAIRI(I,4)
                HAIRI(I,2) = 2.0*HAIRI(I,3) -     HAIRI(I,4)
                FLXRI(I,1) = 3.0*FLXRI(I,3) - 2.0*FLXRI(I,4)
                FLXRI(I,2) = 2.0*FLXRI(I,3) -     FLXRI(I,4)
                FLXPI(I,1) = 3.0*FLXPI(I,3) - 2.0*FLXPI(I,4)
                FLXPI(I,2) = 2.0*FLXPI(I,3) -     FLXPI(I,4)
                PCIPI(I,1) = 3.0*PCIPI(I,3) - 2.0*PCIPI(I,4)
                PCIPI(I,2) = 2.0*PCIPI(I,3) -     PCIPI(I,4)
              ENDDO !i
            ENDIF
          ELSE  !non-global grid
            DO 355 I= 1,IWI+5
              TAIRI(I,1) = 3.0*TAIRI(I,3) - 2.0*TAIRI(I,4)
              TAIRI(I,2) = 2.0*TAIRI(I,3) -     TAIRI(I,4)
              HAIRI(I,1) = 3.0*HAIRI(I,3) - 2.0*HAIRI(I,4)
              HAIRI(I,2) = 2.0*HAIRI(I,3) -     HAIRI(I,4)
              FLXRI(I,1) = 3.0*FLXRI(I,3) - 2.0*FLXRI(I,4)
              FLXRI(I,2) = 2.0*FLXRI(I,3) -     FLXRI(I,4)
              FLXPI(I,1) = 3.0*FLXPI(I,3) - 2.0*FLXPI(I,4)
              FLXPI(I,2) = 2.0*FLXPI(I,3) -     FLXPI(I,4)
              PCIPI(I,1) = 3.0*PCIPI(I,3) - 2.0*PCIPI(I,4)
              PCIPI(I,2) = 2.0*PCIPI(I,3) -     PCIPI(I,4)
  355       CONTINUE
          ENDIF
        ENDIF
C
C       INTERPOLATE FROM NATIVE TO MODEL FLUX GRIDS.
C
        IF     (INTERP.EQ.0) THEN
          IF     (LINTERP(1)) THEN
            CALL LINEAR(TAM,XAF,YAF,IDM,IDM,JDM,
     +                  TAIRI,IWI+5,IWI+4,JWI+4)
          ENDIF
          IF     (LINTERP(2)) THEN
            CALL LINEAR(HAM,XAF,YAF,IDM,IDM,JDM,
     +                  HAIRI,IWI+5,IWI+4,JWI+4)
          ENDIF
          IF     (LINTERP(3)) THEN
            CALL LINEAR(FRM,XAF,YAF,IDM,IDM,JDM,
     +                  FLXRI,IWI+5,IWI+4,JWI+4)
          ENDIF
          IF     (LINTERP(4)) THEN
            CALL LINEAR(FPM,XAF,YAF,IDM,IDM,JDM,
     +                  FLXPI,IWI+5,IWI+4,JWI+4)
          ENDIF
          IF     (LINTERP(5)) THEN
            CALL LINEAR(PCM,XAF,YAF,IDM,IDM,JDM,
     +                  PCIPI,IWI+5,IWI+4,JWI+4)
          ENDIF
        ELSEIF (INTERP.EQ.2) THEN
          IF     (LINTERP(1)) THEN
            CALL BESSEL(TAM,XAF,YAF,IDM,IDM,JDM,
     +                  TAIRI,IWI+5,IWI+4,JWI+4)
          ENDIF
          IF     (LINTERP(2)) THEN
            CALL BESSEL(HAM,XAF,YAF,IDM,IDM,JDM,
     +                  HAIRI,IWI+5,IWI+4,JWI+4)
          ENDIF
          IF     (LINTERP(3)) THEN
            CALL BESSEL(FRM,XAF,YAF,IDM,IDM,JDM,
     +                  FLXRI,IWI+5,IWI+4,JWI+4)
          ENDIF
          IF     (LINTERP(4)) THEN
            CALL BESSEL(FPM,XAF,YAF,IDM,IDM,JDM,
     +                  FLXPI,IWI+5,IWI+4,JWI+4)
          ENDIF
          IF     (LINTERP(5)) THEN
            CALL BESSEL(PCM,XAF,YAF,IDM,IDM,JDM,
     +                  PCIPI,IWI+5,IWI+4,JWI+4)
          ENDIF
        ELSEIF (INTERP.EQ.3) THEN
          IF     (LINTERP(1)) THEN
            CALL BICUBC(TAM,XAF,YAF,IDM,IDM,JDM,
     +                  TAIRI,IWI+5,IWI+5,JWI+5)
          ENDIF
          IF     (LINTERP(2)) THEN
            CALL BICUBC(HAM,XAF,YAF,IDM,IDM,JDM,
     +                  HAIRI,IWI+5,IWI+5,JWI+5)
          ENDIF
          IF     (LINTERP(3)) THEN
            CALL BICUBC(FRM,XAF,YAF,IDM,IDM,JDM,
     +                  FLXRI,IWI+5,IWI+5,JWI+5)
          ENDIF
          IF     (LINTERP(4)) THEN
            CALL BICUBC(FPM,XAF,YAF,IDM,IDM,JDM,
     +                  FLXPI,IWI+5,IWI+5,JWI+5)
          ENDIF
          IF     (LINTERP(5)) THEN
            CALL BICUBC(PCM,XAF,YAF,IDM,IDM,JDM,
     +                  PCIPI,IWI+5,IWI+5,JWI+5)
          ENDIF
        ELSE
          IF     (LINTERP(1)) THEN
            CALL CUBSPL(TAM,XAF,YAF,IDM,IDM,JDM,
     +                  TAIRI,IWI+5,IWI+4,JWI+4, IBD, FXI,FYI,WTAIR3,WK)
          ENDIF
          IF     (LINTERP(2)) THEN
            CALL CUBSPL(HAM,XAF,YAF,IDM,IDM,JDM,
     +                  HAIRI,IWI+5,IWI+4,JWI+4, IBD, FXI,FYI,WTAIR3,WK)
          ENDIF
          IF     (LINTERP(3)) THEN
            CALL CUBSPL(FRM,XAF,YAF,IDM,IDM,JDM,
     +                  FLXRI,IWI+5,IWI+4,JWI+4, IBD, FXI,FYI,WTAIR3,WK)
          ENDIF
          IF     (LINTERP(4)) THEN
            CALL CUBSPL(FPM,XAF,YAF,IDM,IDM,JDM,
     +                  FLXPI,IWI+5,IWI+4,JWI+4, IBD, FXI,FYI,WTAIR3,WK)
C
C           CUBIC SPLINE CAN CREATE NEW (NEGATIVE) MINIMA.
C
            DO I= 1,IDM
              DO J= 1,JDM
                FPM(I,J) = MAX( FPM(I,J), ZERO )
              ENDDO
            ENDDO
          ENDIF
          IF     (LINTERP(5)) THEN
            CALL CUBSPL(PCM,XAF,YAF,IDM,IDM,JDM,
     +                  PCIPI,IWI+5,IWI+4,JWI+4, IBD, FXI,FYI,WTAIR3,WK)
C
C           CUBIC SPLINE CAN CREATE NEW MINIMA.
C
            DO I= 1,IDM
              DO J= 1,JDM
                PCM(I,J) = MAX( PCM(I,J), BIASPC4 )
              ENDDO
            ENDDO
          ENDIF
        ENDIF
C
C       WRITE OUT STATISTICS.
C
        IF     (LINTERP(1)) THEN
        CALL ARCUPD(TAM,IDM,JDM, LARCTIC,.TRUE.)
        CALL MINMAX(TAM,IDM,JDM, XMIN,XMAX)
        CALL AVERMS(TAM,IDM,JDM, XAVE,XRMS)
        WRITE(6,8100) 'TAIR', XMIN,XMAX,XAVE,XRMS
        ENDIF
C
        IF     (LINTERP(2)) THEN
        CALL ARCUPD(HAM,IDM,JDM, LARCTIC,.TRUE.)
        CALL MINMAX(HAM,IDM,JDM, XMIN,XMAX)
        CALL AVERMS(HAM,IDM,JDM, XAVE,XRMS)
        WRITE(6,8100) 'HAIR', XMIN,XMAX,XAVE,XRMS
        ENDIF
C
        IF     (LINTERP(3)) THEN
        CALL ARCUPD(FRM,IDM,JDM, LARCTIC,.TRUE.)
        CALL MINMAX(FRM,IDM,JDM, XMIN,XMAX)
        CALL AVERMS(FRM,IDM,JDM, XAVE,XRMS)
        WRITE(6,8100) 'FLXR', XMIN,XMAX,XAVE,XRMS
        ENDIF
C
        IF     (LINTERP(4)) THEN
        CALL ARCUPD(FPM,IDM,JDM, LARCTIC,.TRUE.)
        CALL MINMAX(FPM,IDM,JDM, XMIN,XMAX)
        CALL AVERMS(FPM,IDM,JDM, XAVE,XRMS)
        WRITE(6,8100) 'FLXP', XMIN,XMAX,XAVE,XRMS
        ENDIF
C
        IF     (LINTERP(5)) THEN
        CALL ARCUPD(PCM,IDM,JDM, LARCTIC,.TRUE.)
        CALL MINMAX(PCM,IDM,JDM, XMIN,XMAX)
        CALL AVERMS(PCM,IDM,JDM, XAVE,XRMS)
        WRITE(6,8100) 'PCIP', 1000.0*XMIN,1000.0*XMAX,
     +                        1000.0*XAVE,1000.0*XRMS
        ENDIF
C
C       WRITE OUT HYCOM FLUXS.
C
        WDAY8  = WDAY(KREC)
        WDAY8I = WDAY(KREC+1)
        IF     (WDAY8I-WDAY8.GT.0.03D0) THEN
C         CORRECT WIND DAYS TO NEAREST HOUR
          WDAY8  = NINT(WDAY8 *24.D0)/24.D0
          WDAY8I = NINT(WDAY8I*24.D0)/24.D0
        ENDIF
        WDAY8I = WDAY8I-WDAY8
C
        IF     (LINTERP(1)) THEN
        CALL ZAIOWR(TAM,MSK,.FALSE., XMIN,XMAX, 10, .FALSE.)
        IF     (IFFILE.EQ.3) THEN
          WRITE(10,4102) '  airtmp',KREC,XMIN,XMAX
        ELSE
          WRITE(10,4112) '  airtmp',WDAY8,WDAY8I,XMIN,XMAX
        ENDIF
        CALL ZHFLSH(10)
        ENDIF
C
        IF     (LINTERP(2)) THEN
        CALL ZAIOWR(HAM,MSK,.FALSE., XMIN,XMAX, 11, .FALSE.)
        IF     (IFFILE.EQ.3) THEN
          WRITE(11,4102) '  vapmix',KREC,XMIN,XMAX
        ELSE
          WRITE(11,4112) '  vapmix',WDAY8,WDAY8I,XMIN,XMAX
        ENDIF
        CALL ZHFLSH(11)
        ENDIF
C
        IF     (LINTERP(3)) THEN
        CALL ZAIOWR(FRM,MSK,.FALSE., XMIN,XMAX, 12, .FALSE.)
        IF     (IFFILE.EQ.3) THEN
          WRITE(12,4102) '  radflx',KREC,XMIN,XMAX
        ELSE
          WRITE(12,4112) '  radflx',WDAY8,WDAY8I,XMIN,XMAX
        ENDIF
        CALL ZHFLSH(12)
        ENDIF
C
        IF     (LINTERP(4)) THEN
        CALL ZAIOWR(FPM,MSK,.FALSE., XMIN,XMAX, 13, .FALSE.)
        IF     (IFFILE.EQ.3) THEN
          WRITE(13,4102) '  shwflx',KREC,XMIN,XMAX
        ELSE
          WRITE(13,4112) '  shwflx',WDAY8,WDAY8I,XMIN,XMAX
        ENDIF
        CALL ZHFLSH(13)
        ENDIF
C
        IF     (LINTERP(5)) THEN
        CALL ZAIOWR(PCM,MSK,.FALSE., XMIN,XMAX, 14, .FALSE.)
        IF     (IFFILE.EQ.3) THEN
          WRITE(14,4102) '  precip',KREC,XMIN,XMAX
        ELSE
          WRITE(14,4112) '  precip',WDAY8,WDAY8I,XMIN,XMAX
        ENDIF
        CALL ZHFLSH(14)
        ENDIF
C
        IF     (IFFILE.EQ.3) THEN
          WRITE(6,6300) KREC,WDAY8
        ELSE
          CALL WNDAY(WDAY8, WYR,FDY)
          WRITE(6,6350) KREC,WDAY8,FDY,NINT(WYR)
        ENDIF
        CALL ZHFLSH(6)
  810 CONTINUE
C
      CALL ZAIOCL(10)
      CLOSE( UNIT=10)
      CALL ZAIOCL(11)
      CLOSE( UNIT=11)
      CALL ZAIOCL(12)
      CLOSE( UNIT=12)
      CALL ZAIOCL(13)
      CLOSE( UNIT=13)
      CALL ZAIOCL(14)
      CLOSE( UNIT=14)
C
C     SUMMARY.
C
      WDAY8  = WDAY(1)
      IF     (WDAY8I.GT.0.03D0) THEN
C       CORRECT WIND DAY TO NEAREST HOUR
        WDAY8  = NINT(WDAY8*24.D0)/24.D0
      ENDIF
      CALL WNDAY(WDAY8, WYR,FDY)
      IF     (WYR.LT.1904.5) THEN
        IF     (WDAY8I.GT.0.03D0) THEN
          WDAY8I = WDAY(NREC+1)
          WDAY8I = NINT(WDAY8I*24.D0)/24.D0
        ELSE
          WDAY8I = WDAY(NREC+1)
        ENDIF
        WRITE(6,6400) NREC,FDY,NINT(WYR),WDAY8I-WDAY8
      ELSE
        IF     (WDAY8I.GT.0.03D0) THEN
          WDAY8I = WDAY(NREC)
          WDAY8I = NINT(WDAY8I*24.D0)/24.D0
        ELSE
          WDAY8I = WDAY(NREC)
        ENDIF
        WRITE(6,6450) NREC,FDY,NINT(WYR),WDAY8I-WDAY8
      ENDIF
      CALL ZHFLSH(6)
      STOP
C
 4101 FORMAT(A79)
 4102 FORMAT(A,': month,range = ',I2.2,1P2E16.7)
 4112 FORMAT(A,': day,span,range =',F12.5,F10.6,1P2E16.7)
 6000 FORMAT(1X,A,2X,A //)
 6200 FORMAT(/ 1X,'MIN,MAX I COORDS = ',F8.2,',',F8.2
     +       / 1X,'MIN,MAX J COORDS = ',F8.2,',',F8.2 /)
 6300 FORMAT(10X,'WRITING FLUX RECORD',I3,'    FDAY =',F10.3 /)
 6350 FORMAT(10X,'WRITING FLUX RECORD',I5,
     +           '    FDAY =',F10.3,
     +            '  FDATE =',F8.3,'/',I4 /)
 6400 FORMAT(I5,' RECORD CLIMATOLOGY STARTING ON',F7.2,'/',I4,
     +   ' COVERING',F9.2,' DAYS')
 6450 FORMAT(I5,' FLUX RECORDS STARTING ON',F7.2,'/',I4,
     +   ' COVERING',F9.2,' DAYS')
 8100 FORMAT(1X,A,': MIN=',F13.8,' MAX=',F13.8,
     +             ' AVE=',F13.8,' RMS=',F13.8)
 9050 FORMAT(// 20X,'**********  ERROR - ',
     +   'INPUT IS GLOBAL AND IWI*DXIN IS NOT 360 DEGREES  ****' //)
 9150 FORMAT(// 20X,'**********  ERROR - ',
     +   'THE OUTPUT GRID IS NOT INSIDE THE INPUT GRID  **********' //)
 9200 FORMAT(// 20X,'**********  ERROR - ',
     +   'IFFILE=3 BUT RECORDS ARE NOT MONTHLY  **********' //)
C     END OF PROGRAM FLXINT.
      END
      SUBROUTINE WNDAY(WDAY, YEAR,DAY)
      IMPLICIT NONE
      REAL*8 WDAY
      REAL*4 YEAR,DAY
C
C**********
C*
C  1) CONVERT 'FLUX DAY' INTO JULIAN DAY AND YEAR.
C
C  2) THE 'FLUX DAY' IS THE NUMBER OF DAYS SINCE 001/1901 (WHICH IS
C      FLUX DAY 1.0).
C     FOR EXAMPLE:
C      A) YEAR=1901.0 AND DAY=1.0, REPRESENTS 0000Z HRS ON 001/1901
C         SO WDAY WOULD BE 1.0.
C      B) YEAR=1901.0 AND DAY=2.5, REPRESENTS 1200Z HRS ON 002/1901
C         SO WDAY WOULD BE 2.5.
C     YEAR MUST BE NO LESS THAN 1901.0, AND NO GREATER THAN 2099.0.
C     NOTE THAT YEAR 2000 IS A LEAP YEAR (BUT 1900 AND 2100 ARE NOT).
C
C  3) ALAN J. WALLCRAFT, PLANNING SYSTEMS INC., FEBRUARY 1993.
C*
C**********
C
      INTEGER IYR,NLEAP
      REAL*8  WDAY1
C
C     FIND THE RIGHT YEAR.
C
      IYR   = (WDAY-1.0)/365.25
      NLEAP = IYR/4
      WDAY1 = 365.0*IYR + NLEAP + 1.0
      DAY   = WDAY - WDAY1 + 1.0
      IF     (WDAY1.GT.WDAY) THEN
        IYR   = IYR - 1
      ELSEIF (DAY.GE.367.0) THEN
        IYR   = IYR + 1
      ELSEIF (DAY.GE.366.0 .AND. MOD(IYR,4).NE.3) THEN
        IYR   = IYR + 1
      ENDIF
      NLEAP = IYR/4
      WDAY1 = 365.0*IYR + NLEAP + 1.0
C
C     RETURN YEAR AND JULIAN DAY.
C
      YEAR = 1901 + IYR
      DAY  = WDAY - WDAY1 + 1.0
      RETURN
C     END OF WNDAY.
      END
      SUBROUTINE FREADI(IWI,JWI)
      USE netcdf   ! NetCDF fortran 90 interface
      IMPLICIT NONE
C
      INTEGER IWI,JWI
C
C**********
C*
C  1)  INITIALIZE ARRAY SIZES FOR READING NATIVE FLUXS.
C
C      SEE 'FREAD0' FOR HEADER INITIALIZATION.
C      SEE 'FREAD1' FOR READING ACTUAL FLUX RECORDS.
C
C  2) ON EXIT:
C      IWI = LONGIUDINAL ARRAY DIMENSION
C      JWI = LATITUDINAL ARRAY DIMENSION
C
C  3) ALAN J. WALLCRAFT,  JULY 2010.
C*
C*********
C
      CHARACTER*240 CFILE
      INTEGER       ncFID,ncDID,ncVID
C
C     FIRST INPUT FLUX FILE.
C
      CFILE = ' '
      CALL GETENV('CDF071',CFILE)
      IF     (CFILE.EQ.' ') THEN
        WRITE(0,*) 'kp_nc: no CDF071 environment variable'
        CALL EXIT(1)
        STOP
      ENDIF
      call nchek('nf90_open',
     &            nf90_open(trim(CFILE), nf90_nowrite, ncFID))
C
      call nchek('nf90_inq_dimid-Latitude',
     &            nf90_inq_dimid(        ncFID, 'Latitude',ncDID))
      call nchek('nf90_inquire_dimension-Latitude',
     &            nf90_inquire_dimension(ncFID,            ncDID,
     &                                                 len=JWI))
C
      call nchek('nf90_inq_dimid-Longitude',
     &            nf90_inq_dimid(        ncFID, 'Longitude',ncDID))
      call nchek('nf90_inquire_dimension-Longitude',
     &            nf90_inquire_dimension(ncFID,             ncDID,
     &                                                  len=IWI))
      call nchek("nf90_close",
     &            nf90_close(ncFID))
      RETURN
C     END OF FREADI.
      END
      SUBROUTINE FREAD0(NREC,IUREC,IFREC,WDAY,MAXREC, WSTART,WEND,
     +                  IWI,JWI,WLON,WLAT)
      USE netcdf   ! NetCDF fortran 90 interface
      IMPLICIT NONE
C
      INTEGER NREC,IUREC,IFREC,MAXREC
      INTEGER IWI,JWI
      REAL*8  WLON(IWI),WLAT(JWI)
      REAL*4  WSTART,WEND
      REAL*4  WDAY(MAXREC+1)
C
C**********
C*
C  1)  INITIALIZE FOR READING NATIVE FLUXS.
C
C      SEE 'FREAD1' FOR READING ACTUAL FLUX RECORDS.
C
C  2) ON EXIT:
C      NREC  = NUMBER OF FLUX RECORDS REQUIRED
C      IUREC = UNIT NUMBER   OF FIRST INPUT RECORD (71...99)
C      IFREC = RECORD NUMBER OF FIRST INPUT RECORD
C      WDAY  = FLUX DAYS FOR RECORDS 1...NREC
C
C      IF TSTART=-1.0, THE INPUT IS A SINGLE CLIMATOLOGY FILE AND
C       THE ENTIRE FILE IS OUTPUT.
C      OTHERWISE, ONLY THE RECORDS THAT SPAN WSTART TO WEND ARE
C       OUTPUT (AND LISTED IN WDAY) BUT THIS CAN INVOLVE MULTIPLE
C       INPUT FILES AND THE INITIAL FLUX DAY NEED NOT BE IN THE FIRST
C       (UNIT 71) FLUX FILE.
C
C  3) ALAN J. WALLCRAFT,
C      UPDATED FOR NETCDF INPUT IN JULY 2010.
C*
C*********
C
      CHARACTER*6       CENV
      CHARACTER*240     CFILE
      INTEGER           ncFID,ncDID,ncVID
C
      CHARACTER*40     CWNAME
      INTEGER          IR,IFILE,IUNIT,NFREC
      REAL*4           WFDAY(9000)
      DOUBLE PRECISION TIME(9000),TIME_NEXT
*     write(6,*) 'FREAD0 - MAXREC,WSTART,WEND = ',
*    +           MAXREC,WSTART,WEND
*     call zhflsh(6)
C
C     FIRST INPUT FLUX FILE.
C
      CFILE = ' '
      CALL GETENV('CDF071',CFILE)
      IF     (CFILE.EQ.' ') THEN
        WRITE(0,*) 'kp_nc: no CDF071 environment variable'
        CALL EXIT(1)
        STOP
      ENDIF
      call nchek('nf90_open',
     &            nf90_open(trim(CFILE), nf90_nowrite, ncFID))
C
      call nchek("nf90_get_att-title",
     &            nf90_get_att(ncFID,nf90_global,
     &                         "title",
     &                         CWNAME))
C
      call nchek('nf90_inq_varid-Latitude',
     &            nf90_inq_varid(ncFID,'Latitude', ncVID))
      call nchek('nf90_get_var-Latitude',
     &            nf90_get_var(  ncFID,            ncVID,WLAT(:)))
c
      call nchek('nf90_inq_varid-Longitude',
     &            nf90_inq_varid(ncFID,'Longitude',ncVID))
      call nchek('nf90_get_var-Longitude',
     &            nf90_get_var(  ncFID,            ncVID,WLON(:)))
c
      call nchek('nf90_inq_dimid-MT',
     &            nf90_inq_dimid(        ncFID, 'MT',ncDID))
      call nchek('nf90_inquire_dimension-MT',
     &            nf90_inquire_dimension(ncFID,      ncDID,
     &                                               len=NFREC))
      IF     (NFREC.GE.9000) THEN
        WRITE(6,9150) NFREC
        CALL ZHFLSH(6)
        STOP
      ENDIF
      call nchek('nf90_inq_varid-MT',
     &            nf90_inq_varid(ncFID,'MT',ncVID))
      call nchek('nf90_get_var-MT',
     &            nf90_get_var(  ncFID,     ncVID,TIME(1:NFREC)))
      call nchek("nf90_get_att-next_MT",
     &            nf90_get_att(  ncFID,     ncVID,
     &                                      "next_MT",
     &                                      TIME_NEXT))
      call nchek("nf90_close",
     &            nf90_close(ncFID))
      WFDAY(1:NFREC) = TIME(1:NFREC)
      WFDAY(NFREC+1) = TIME_NEXT
C
      IF     (WSTART.LT.0.0) THEN
C
C       FLUX CLIMATOLOGY.
C
        IF     (MAXREC.LT.NFREC) THEN
          WRITE(6,9100) MAXREC,NFREC
          CALL ZHFLSH(6)
          STOP
        ENDIF
C
        IUREC = 71
        IFREC =  1
        NREC  = NFREC
        DO 110 IR= 1,NFREC+1
          WDAY(IR) = WFDAY(IR)
  110   CONTINUE
      ELSE
C
C       FLUXS FROM WSTART TO WEND.
C
        DO 210 IFILE= 71,99
*         write(6,*) 'DO 210 - IFILE,NFREC,WFDAY = ',
*    +               IFILE,NFREC,WFDAY(1),WFDAY(NFREC+1)
*         call zhflsh(6)
          IF     (WFDAY(1).GT.WSTART) THEN
            WRITE(6,9200) IFILE,WFDAY(1),WSTART
            CALL ZHFLSH(6)
            STOP
          ELSEIF (WFDAY(NFREC+1).LE.WSTART) THEN
C
C           NOT IN THIS FILE, READ NEXT HEADER.
C
            IF     (IFILE.NE.99) THEN
              IUNIT = IFILE + 1
              WRITE(CENV,'(a,i3.3)') 'CDF',IUNIT
              CFILE = ' '
              CALL GETENV(CENV,CFILE)
              IF     (CFILE.EQ.' ') THEN
                GOTO 950
              ENDIF
              call nchek('nf90_open',
     &                    nf90_open(trim(CFILE), nf90_nowrite, ncFID))
C
              call nchek('nf90_inq_dimid-MT',
     &                    nf90_inq_dimid(        ncFID, 'MT',ncDID))
              call nchek('nf90_inquire_dimension-MT',
     &                    nf90_inquire_dimension(ncFID,      ncDID,
     &                                                   len=NFREC))
              IF     (NFREC.GE.9000) THEN
                WRITE(6,9150) NFREC
                CALL ZHFLSH(6)
                STOP
              ENDIF
              call nchek('nf90_inq_varid-MT',
     &                    nf90_inq_varid(ncFID,'MT',ncVID))
              call nchek('nf90_get_var-MT',
     &                    nf90_get_var(  ncFID,     ncVID,
     &                                              TIME(1:NFREC)))
              call nchek("nf90_get_att-next_MT",
     &                    nf90_get_att(  ncFID,     ncVID,
     &                                              "next_MT",
     &                                              TIME_NEXT))
              call nchek("nf90_close",
     &                    nf90_close(ncFID))
              WFDAY(1:NFREC) = TIME(1:NFREC)
              WFDAY(NFREC+1) = TIME_NEXT
            ENDIF
          ELSE
C
C           INITIAL RECORD IS IN THIS FILE.
C
            IUREC = IFILE
            IF     (WFDAY(2).GT.WSTART) THEN
              NREC    = 1
              IFREC   = 1
              WDAY(1) = WFDAY(1)
*             write(6,*) 'NREC,WDAY = ',NREC,WDAY(NREC)
*             call zhflsh(6)
            ELSE
              NREC  = 0
              IFREC = NFREC
            ENDIF
            DO 220 IR= 2,NFREC
              IF     (WFDAY(IR-1).GE.WEND) THEN
                GOTO 1220
              ENDIF
              IF     (WFDAY(IR+1).GT.WSTART) THEN
                NREC       = NREC + 1
                IFREC      = MIN( IR, IFREC )
                WDAY(NREC) = WFDAY(IR)
*               write(6,*) 'NREC,WDAY = ',NREC,WDAY(NREC)
*               call zhflsh(6)
                IF     (MAXREC.LT.NREC) THEN
                  WRITE(6,9100) MAXREC,NREC
                  CALL ZHFLSH(6)
                  STOP
                ENDIF
              ENDIF
  220       CONTINUE
            IR = NFREC+1
 1220       CONTINUE
            WDAY(NREC+1) = WFDAY(IR)
C
            IF     (WDAY(NREC).LT.WEND) THEN
              DO 230 IUNIT= IFILE+1,99
                WRITE(CENV,'(a,i3.3)') 'CDF',IUNIT
                CFILE = ' '
                CALL GETENV(CENV,CFILE)
                IF     (CFILE.EQ.' ') THEN
                  GOTO 950
                ENDIF
                call nchek('nf90_open',
     &                      nf90_open(trim(CFILE),
     &                                nf90_nowrite,
     &                                ncFID))
C
                call nchek('nf90_inq_dimid-MT',
     &                      nf90_inq_dimid(        ncFID, 'MT',ncDID))
                call nchek('nf90_inquire_dimension-MT',
     &                      nf90_inquire_dimension(ncFID,      ncDID,
     &                                                     len=NFREC))
                IF     (NFREC.GE.9000) THEN
                  WRITE(6,9150) NFREC
                  CALL ZHFLSH(6)
                  STOP
                ENDIF
                call nchek('nf90_inq_varid-MT',
     &                      nf90_inq_varid(ncFID,'MT',ncVID))
                call nchek('nf90_get_var-MT',
     &                      nf90_get_var(  ncFID,     ncVID,
     &                                                TIME(1:NFREC)))
                call nchek("nf90_get_att-next_MT",
     &                      nf90_get_att(  ncFID,     ncVID,
     &                                                "next_MT",
     &                                                TIME_NEXT))
                call nchek("nf90_close",
     &                      nf90_close(ncFID))
                WFDAY(1:NFREC) = TIME(1:NFREC)
                WFDAY(NFREC+1) = TIME_NEXT
*               write(6,*) 'DO 230 - IUNIT,NFREC,WFDAY = ',
*    +                     IUNIT,NFREC,WFDAY(1),WFDAY(NFREC+1)
*               call zhflsh(6)
                IF     (WFDAY(1).NE.WDAY(NREC+1)) THEN
                  WRITE(6,9300) IUNIT,WFDAY(1),WDAY(NREC+1)
                  CALL ZHFLSH(6)
                  STOP
                ENDIF
                NREC = NREC + 1
*               write(6,*) 'NREC,WDAY = ',NREC,WDAY(NREC)
*               call zhflsh(6)
                DO 240 IR= 2,NFREC
                  IF     (WFDAY(IR-1).GE.WEND) THEN
                    GOTO 1240
                  ENDIF
                  NREC       = NREC + 1
                  WDAY(NREC) = WFDAY(IR)
*                 write(6,*) 'NREC,WDAY = ',NREC,WDAY(NREC)
*                 call zhflsh(6)
                  IF     (MAXREC.LT.NREC) THEN
                    WRITE(6,9100) MAXREC,NREC
                    CALL ZHFLSH(6)
                    STOP
                  ENDIF
  240           CONTINUE
                IR = NFREC+1
 1240           CONTINUE
                WDAY(NREC+1) = WFDAY(IR)
C
                IF     (WDAY(NREC).GE.WEND) THEN
                  GOTO 1230
                ENDIF
  230         CONTINUE
 1230         CONTINUE
            ENDIF
            GOTO 1210
          ENDIF
  210   CONTINUE
 1210   CONTINUE
      ENDIF
      RETURN
C
C     NO MORE FLUX FILES.
C
  950 CONTINUE
        WRITE(6,9500) IUNIT
        CALL ZHFLSH(6)
        STOP
C
 9000 FORMAT(// 20X,'*****  ERROR IN FREAD0  -  ',
     +   'INPUT DOES NOT AGREE WITH PARAMETERS  *****' //
     +   1X,'IWI,JWI   = ',I5,  I10   //)
 9100 FORMAT(// 20X,'*****  ERROR IN FREAD0  -  ',
     +   'MAXREC TOO SMALL (MAXREC,NREC = ',2I6,')  *****' //)
 9150 FORMAT(// 20X,'*****  ERROR IN FREAD0  -  ',
     +   'NFREC LARGER THAN 8999 (NFREC = ',I6,')  *****' //)
 9200 FORMAT(// 20X,'*****  ERROR IN FREAD0  -  ',
     +   'INITIAL FLUX DAY ON UNIT ',I2,' IS ',F8.2,
     +   ', BUT WSTART = ',F8.2,'  *****' //)
 9300 FORMAT(// 20X,'*****  ERROR IN FREAD0  -  ',
     +   'INITIAL FLUX DAY ON UNIT ',I2,' IS ',F8.2,
     +   ', BUT SHOULD BE ',F8.2,'  *****' //)
 9500 FORMAT(// 10X,'*****  ERROR IN FREAD0  -  I/O ERROR ON UNIT',
     +   I3,'  *****' //)
C     END OF FREAD0.
      END
      SUBROUTINE FREAD1(TAIRIN,HAIRIN,FLXRIN,FLXPIN,PCIPIN,IWI,JWI,
     +                  FLDI2, CNAME,
     +                  IUREC,IFREC, KREC,WDAYK, IFTYPE)
      USE netcdf   ! NetCDF fortran 90 interface
      IMPLICIT NONE
C
      INTEGER     IWI,JWI, IUREC,IFREC,KREC,IFTYPE
      CHARACTER*6 CNAME(5)
      INTEGER*2    FLDI2(IWI,JWI)
      REAL*4      TAIRIN(IWI,JWI),HAIRIN(IWI,JWI),
     +            FLXRIN(IWI,JWI),FLXPIN(IWI,JWI),
     +            PCIPIN(IWI,JWI),WDAYK
C
C**********
C*
C  1)  READ THE 'KREC'-TH REQUIRED NATIVE FLUX RECORD.
C      THIS IS EITHER RECORD 'IFREC' ON UNIT 'IUREC', OR RECORD 1
C      ON UNIT 'IUREC'+1.
C
C      MUST BE CALLED WITH 'KREC' IN ASCENDING ORDER.
C
C      SEE 'FREAD0' FOR HEADER INITIALIZATION.
C
C  2) ON EXIT:
C      TAIRIN  = AIR TEMPERATURE FIELD
C      HAIRIN  = AIR RELATIVE HUMDITITY FIELD
C      FLXRIN  = RADIATION HEAT FLUX FIELD
C      FLXPIN  = SHORTWAVE RADIATION HEAT FLUX FIELD
C      PCIPIN  = PRECIPITATION FIELD
C       FLDI2  = WORKSPACE ARRAY
C      IUREC   = UNIT NUMBER   OF NEXT INPUT RECORD (71...99)
C      IFREC   = RECORD NUMBER OF NEXT INPUT RECORD
C
C  3) ALAN J. WALLCRAFT,  PLANNING SYSTEMS INC.,  OCTOBER 1995.
C      UPDATED FOR NETCDF INPUT IN JULY 2010.
C*
C*********
C
      INTEGER NFREC
      REAL*4  WFDAY(9000)
      SAVE    NFREC,WFDAY
C
      CHARACTER*6      CENV
      CHARACTER*240    CFILE
      INTEGER          ncFID,ncDID,ncVID
      DOUBLE PRECISION TIME(9000),TIME_NEXT
      REAL             SCALE_F,ADD_OFF
C
      CHARACTER*40 CWNAME
      INTEGER      I,IR,J
C
*     write(6,*) 'krec,ifrec,iurec = ',krec,ifrec,iurec
      IF     (KREC.GT.1 .AND. IFREC.GT.NFREC) THEN
C
C       NEW FLUX FILE.
C
        IUREC = IUREC + 1
        IFREC = 1
      ENDIF
*     write(6,*) 'krec,ifrec,iurec = ',krec,ifrec,iurec
C
      WRITE(CENV,'(a,i3.3)') 'CDF',IUREC
      CFILE = ' '
      CALL GETENV(CENV,CFILE)
      call nchek('nf90_open',
     &            nf90_open(trim(CFILE), nf90_nowrite, ncFID))
C
C     READ A HEADER IF REQUIRED.
C
      IF     (KREC.EQ.1 .OR. IFREC.EQ.1) THEN
C
C       FIRST RECORD FOR THIS FILE.
C
        call nchek('nf90_inq_dimid-MT',
     &              nf90_inq_dimid(        ncFID, 'MT',ncDID))
        call nchek('nf90_inquire_dimension-MT',
     &              nf90_inquire_dimension(ncFID,      ncDID,
     &                                             len=NFREC))
        IF     (NFREC.GE.9000) THEN
          WRITE(6,9150) NFREC
          CALL ZHFLSH(6)
          STOP
        ENDIF
        call nchek('nf90_inq_varid-MT',
     &              nf90_inq_varid(ncFID,'MT',ncVID))
        call nchek('nf90_get_var-MT',
     &              nf90_get_var(  ncFID,     ncVID,
     &                                        TIME(1:NFREC)))
        call nchek("nf90_get_att-next_MT",
     &              nf90_get_att(  ncFID,     ncVID,
     &                                        "next_MT",
     &                                        TIME_NEXT))
        WFDAY(1:NFREC) = TIME(1:NFREC)
        WFDAY(NFREC+1) = TIME_NEXT
      ENDIF
C
C     CHECK THE FLUX DAY.
C
      IF     (WDAYK.NE.WFDAY(IFREC)) THEN
        WRITE(6,9000) WFDAY(IFREC),WDAYK
        CALL ZHFLSH(6)
        STOP
      ENDIF
C
C     READ THE FLUX RECORDS.
C
      if     (IFTYPE.NE.1 .AND. IFTYPE.NE.2) THEN
C
C       TAIR
C
        call nchek('nf90_inq_varid-CNAME1',
     &              nf90_inq_varid(ncFID,CNAME(1),ncVID))
        call nchek('nf90_get_var-FLDI2',
     &              nf90_get_var(  ncFID,      ncVID,
     &                                         FLDI2(:,:),
     &                                         (/ 1,1,IFREC /) ))
        call nchek("nf90_get_att-scale_factor",
     &              nf90_get_att(  ncFID,     ncVID,
     &                                        "scale_factor",
     &                                        SCALE_F))
        call nchek("nf90_get_att-add_offset",
     &              nf90_get_att(  ncFID,     ncVID,
     &                                        "add_offset",
     &                                        ADD_OFF))
        DO I= 1,IWI
          DO J= 1,JWI
            TAIRIN(I,J) = SCALE_F*FLDI2(I,J) + ADD_OFF
          ENDDO !j
        ENDDO !i
      ELSE
        DO I= 1,IWI
          DO J= 1,JWI
            TAIRIN(I,J) = 0.0
          ENDDO !j
        ENDDO !i
      ENDIF
C
      if     (IFTYPE.NE.1 .AND. IFTYPE.NE.2) THEN
C
C       HAIR
C
        call nchek('nf90_inq_varid-CNAME2',
     &              nf90_inq_varid(ncFID,CNAME(2),ncVID))
        call nchek('nf90_get_var-FLDI2',
     &              nf90_get_var(  ncFID,      ncVID,
     &                                         FLDI2(:,:),
     &                                         (/ 1,1,IFREC /) ))
        call nchek("nf90_get_att-scale_factor",
     &              nf90_get_att(  ncFID,     ncVID,
     &                                        "scale_factor",
     &                                        SCALE_F))
        call nchek("nf90_get_att-add_offset",
     &              nf90_get_att(  ncFID,     ncVID,
     &                                        "add_offset",
     &                                        ADD_OFF))
        DO I= 1,IWI
          DO J= 1,JWI
            HAIRIN(I,J) = SCALE_F*FLDI2(I,J) + ADD_OFF
          ENDDO !j
        ENDDO !i
      ELSE
        DO I= 1,IWI
          DO J= 1,JWI
            HAIRIN(I,J) = 0.0
          ENDDO !j
        ENDDO !i
      ENDIF
C
      if     (IFTYPE.NE.1 .AND. IFTYPE.NE.3) THEN
C
C       FLXR
C
        call nchek('nf90_inq_varid-CNAME3',
     &              nf90_inq_varid(ncFID,CNAME(3),ncVID))
        call nchek('nf90_get_var-FLDI2',
     &              nf90_get_var(  ncFID,      ncVID,
     &                                         FLDI2(:,:),
     &                                         (/ 1,1,IFREC /) ))
        call nchek("nf90_get_att-scale_factor",
     &              nf90_get_att(  ncFID,     ncVID,
     &                                        "scale_factor",
     &                                        SCALE_F))
        call nchek("nf90_get_att-add_offset",
     &              nf90_get_att(  ncFID,     ncVID,
     &                                        "add_offset",
     &                                        ADD_OFF))
        DO I= 1,IWI
          DO J= 1,JWI
            FLXRIN(I,J) = SCALE_F*FLDI2(I,J) + ADD_OFF
          ENDDO !j
        ENDDO !i
      ELSE
        DO I= 1,IWI
          DO J= 1,JWI
            FLXRIN(I,J) = 0.0
          ENDDO !j
        ENDDO !i
      ENDIF
C
      if     (IFTYPE.NE.1 .AND. IFTYPE.NE.2 .AND. IFTYPE.NE.3) THEN
C
C       FLXP
C
        call nchek('nf90_inq_varid-CNAME4',
     &              nf90_inq_varid(ncFID,CNAME(4),ncVID))
        call nchek('nf90_get_var-FLDI2',
     &              nf90_get_var(  ncFID,      ncVID,
     &                                         FLDI2(:,:),
     &                                         (/ 1,1,IFREC /) ))
        call nchek("nf90_get_att-scale_factor",
     &              nf90_get_att(  ncFID,     ncVID,
     &                                        "scale_factor",
     &                                        SCALE_F))
        call nchek("nf90_get_att-add_offset",
     &              nf90_get_att(  ncFID,     ncVID,
     &                                        "add_offset",
     &                                        ADD_OFF))
        DO I= 1,IWI
          DO J= 1,JWI
            FLXPIN(I,J) = SCALE_F*FLDI2(I,J) + ADD_OFF
          ENDDO !j
        ENDDO !i
      ELSE
        DO I= 1,IWI
          DO J= 1,JWI
            FLXPIN(I,J) = 0.0
          ENDDO !j
        ENDDO !i
      ENDIF
C
      if     (IFTYPE.NE.4 .AND. IFTYPE.NE.2 .AND. IFTYPE.NE.3) THEN
C
C       PCIP
C
        call nchek('nf90_inq_varid-CNAME5',
     &              nf90_inq_varid(ncFID,CNAME(5),ncVID))
        call nchek('nf90_get_var-FLDI2',
     &              nf90_get_var(  ncFID,      ncVID,
     &                                         FLDI2(:,:),
     &                                         (/ 1,1,IFREC /) ))
        call nchek("nf90_get_att-scale_factor",
     &              nf90_get_att(  ncFID,     ncVID,
     &                                        "scale_factor",
     &                                        SCALE_F))
        call nchek("nf90_get_att-add_offset",
     &              nf90_get_att(  ncFID,     ncVID,
     &                                        "add_offset",
     &                                        ADD_OFF))
        DO I= 1,IWI
          DO J= 1,JWI
            PCIPIN(I,J) = SCALE_F*FLDI2(I,J) + ADD_OFF
          ENDDO !j
        ENDDO !i
      ELSE
        DO I= 1,IWI
          DO J= 1,JWI
            PCIPIN(I,J) = 0.0
          ENDDO !j
        ENDDO !i
      ENDIF
C
      call nchek("nf90_close",
     &            nf90_close(ncFID))
C
C     UPDATE FOR NEXT RECORD.
C
      IFREC = IFREC + 1
      RETURN
C
 9000 FORMAT(// 20X,'*****  ERROR IN FREAD1  -  ',
     +   'INPUT FLUX DAY IS ',F8.2,
     +   ', BUT SHOULD BE ',F8.2,'  *****' //)
 9150 FORMAT(// 20X,'*****  ERROR IN FREAD1  -  ',
     +   'NFREC LARGER THAN 8999 (NFREC = ',I6,')  *****' //)
C     END OF FREAD1.
      END
      SUBROUTINE FREADM(MASKIN,MASKI2,IWI,JWI,IUNIT)
      USE netcdf   ! NetCDF fortran 90 interface
      IMPLICIT NONE
C
      INTEGER   IWI,JWI, IUNIT
      INTEGER*2 MASKI2(IWI,JWI)
      REAL*4    MASKIN(IWI,JWI)
C
C**********
C*
C  1)  READ THE NATIVE LAND/SEA MASK ON UNIT 'IUNIT'.
C
C  2) ON EXIT:
C      MASKIN  = LAND/SEA MASK FIELD
C
C  3) ALAN J. WALLCRAFT,  NAVAL RESEARCH LABORATORY,  APRIL 2004.
C      UPDATED FOR NETCDF INPUT IN JULY 2010.
C*
C*********
C
      CHARACTER*6      CENV
      CHARACTER*240    CFILE
      INTEGER          ncFID,ncDID,ncVID
      REAL             SCALE_F,ADD_OFF
      INTEGER          I,IWIT,J,JWIT
C
C     NetCDF FILE FROM ENVIRONMENT VARIABLE CDFXXX
C
      WRITE(CENV,'(a,i3.3)') 'CDF',IUNIT
      CFILE = ' '
      CALL GETENV(CENV,CFILE)
      call nchek('nf90_open',
     &            nf90_open(trim(CFILE), nf90_nowrite, ncFID))
C
C     CHECK ARRAY DIMENSIONS.
C
      call nchek('nf90_inq_dimid-Latitude',
     &            nf90_inq_dimid(        ncFID, 'Latitude', ncDID))
      call nchek('nf90_inquire_dimension-Latitude',
     &            nf90_inquire_dimension(ncFID,             ncDID,
     &                                                   len=JWIT))
      call nchek('nf90_inq_dimid-Longitude',
     &            nf90_inq_dimid(        ncFID, 'Longitude',ncDID))
      call nchek('nf90_inquire_dimension-Longitude',
     &            nf90_inquire_dimension(ncFID,             ncDID,
     &                                                   len=IWIT))
      IF     (IWIT.NE.IWI .OR.
     +        JWIT.NE.JWI     ) THEN
        WRITE(6,9000) IWI,JWI
        CALL ZHFLSH(6)
        STOP
      ENDIF
C
C     NetCDF VARIABLE NAME MUST BE "lndsea"
C
      call nchek('nf90_inq_varid-CNAME',
     &            nf90_inq_varid(ncFID,'lndsea',ncVID))
      call nchek('nf90_get_var-MASKI2',
     &            nf90_get_var(  ncFID,      ncVID,
     &                                       MASKI2(:,:)))
      call nchek("nf90_get_att-scale_factor",
     &            nf90_get_att(  ncFID,     ncVID,
     &                                      "scale_factor",
     &                                      SCALE_F))
      call nchek("nf90_get_att-add_offset",
     &            nf90_get_att(  ncFID,     ncVID,
     &                                      "add_offset",
     &                                      ADD_OFF))
      call nchek("nf90_close",
     &            nf90_close(ncFID))
C
      DO I= 1,IWI
        DO J= 1,JWI
          MASKIN(I,J) = SCALE_F*MASKI2(I,J) + ADD_OFF
        ENDDO !j
      ENDDO !i
      RETURN
 9000 FORMAT(// 20X,'*****  ERROR IN FREADM  -  ',
     +   'INPUT DOES NOT AGREE WITH PARAMETERS  *****' //
     +   1X,'IWI,JWI   =',I6,  I10   //)
C     END OF FREADM.
      END
      subroutine nchek(cnf90,status)
      use netcdf   ! NetCDF fortran 90 interface
      implicit none
c
      character*(*), intent(in) :: cnf90
      integer,       intent(in) :: status
c
c     subroutine to handle NetCDF errors
c
*     if     (.TRUE. ) then !debug
      if     (.FALSE.) then !nodebug
        write(6,'(a)') trim(cnf90)
      endif

      if (status /= nf90_noerr) then
        write(6,'(/a)')   'error in profout - from NetCDF library'
        write(6,'(a/)')   trim(cnf90)
        write(6,'(a/)')   trim(nf90_strerror(status))
        stop
      end if
      end subroutine nchek

max_stars_repo_name: TillRasmussen/HYCOM-tools

max_issues_repo_name: TillRasmussen/HYCOM-tools

max_forks_repo_name: TillRasmussen/HYCOM-tools


text:
function ffun(x, y)
    implicit none
    !GCC$ ATTRIBUTES DLLEXPORT, STDCALL :: FFUN
    double precision :: x, y, ffun
    ffun = x + y * y
end function

ext: f

max_stars_repo_name: mullikine/RosettaCodeData

max_issues_repo_name: seanwallawalla-forks/RosettaCodeData

max_forks_repo_name: seanwallawalla-forks/RosettaCodeData


        !COMPILER-GENERATED INTERFACE MODULE: Fri Mar  6 15:35:15 2020
        ! This source file is for reference only and may not completely
        ! represent the generated interface used by the compiler.
        MODULE CARBT_DIAGS__genmod
          INTERFACE
            SUBROUTINE CARBT_DIAGS(ISTEP)
              INTEGER(KIND=4) :: ISTEP
            END SUBROUTINE CARBT_DIAGS
          END INTERFACE
        END MODULE CARBT_DIAGS__genmod


text:
!
!     CalculiX - A 3-dimensional finite element program
!              Copyright (C) 1998-2015 Guido Dhondt
!
!     This program is free software; you can redistribute it and/or
!     modify it under the terms of the GNU General Public License as
!     published by the Free Software Foundation(version 2);
!
!
!     This program 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 General Public License for more details.
!
!     You should have received a copy of the GNU General Public License
!     along with this program; if not, write to the Free Software
!     Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
!
      subroutine checktruecontact(ntie,tieset,tietol,elcon,itruecontact,
     &  ncmat_,ntmat_)
!
!     check whether for face-to-face penalty contact the
!     surface behavior definition is such that true contact is
!     defined and not just tied contact
!
      implicit none
!
      character*81 tieset(3,*)
!
      integer itruecontact,i,ntie,imat,ncmat_,ntmat_
!
      real*8 tietol(3,*),elcon(0:ncmat_,ntmat_,*)
!
!     if at least one tied contact is present, itruecontact
!     is set to zero and no check is performed whether tension
!     occurs in the contact areas
!
      itruecontact=1
      do i=1,ntie
         if(tieset(1,i)(81:81).eq.'C') then
            imat=int(tietol(2,i))
            if(int(elcon(3,1,imat)).eq.4) then
               itruecontact=0
               exit
            endif
         endif
      enddo
!
      return
      end

max_stars_repo_name: alleindrach/calculix-desktop

max_issues_repo_name: alleindrach/calculix-desktop

max_forks_repo_name: alleindrach/calculix-desktop


text:
      implicit none
      integer uin,uouts,uouts_cen,uouts_min,uouts_max,
     &        uoutpdf,uoutpdf_cen,uoutpdf_min,uoutpdf_max
      parameter (uin=99,uouts=98,uoutpdf=97,uouts_cen=96,
     &           uoutpdf_cen=95,uouts_min=94,uoutpdf_min=93,
     &           uouts_max=92,uoutpdf_max=91)
      integer i,j,nng,nps,lb,ipdfmin,ipdfmax,iscalevar,iPDFvar,
     &        npdf,i0,i1
      logical done
      character*10 str
      character*100 infile,outfile_s,outfile_pdf,outfilebase,var
      character*100 outfile_s_central,outfile_s_min,outfile_s_max
      character*100 outfile_pdf_central,outfile_pdf_min,outfile_pdf_max
      double precision xs(0:1000,10),ys(1000,10)
      double precision xpdf(0:1000,40),ypdf(1000,40),
     &     ysmin(1000),ysmax(1000),ypdfmin(1000),ypdfmax(1000)
      integer nc
      parameter (nc=15)
      character*(nc) scale(10),pdf(100)
      data scale/'central value  ',
     &           'muR=1.0 muF=1.0','muR=1.0 muF=2.0','muR=1.0 muF=0.5',
     &           'muR=2.0 muF=1.0','muR=2.0 muF=2.0','muR=2.0 muF=0.5',
     &           'muR=0.5 muF=1.0','muR=0.5 muF=2.0','muR=0.5 muF=0.5'/
c
      write (*,*) 'Give input .top file'
      read(*,'(a)') infile
      open(unit=uin,file=infile,status='old')
      outfilebase=infile(1:index(infile,'.top')-1)
      lb=index(outfilebase,' ')-1
c
      write (*,*) 'Type 1 for scale variations, 0 not'
      read (*,*) iscalevar
      if(iscalevar.ne.0.and.iscalevar.ne.1)then
         write(*,*)'You should type 0 or 1, not ',iscalevar
         stop
      endif
c
      write (*,*) 'Type 1 for PDF variations, 0 not'
      read (*,*) iPDFvar
      if(iPDFvar.ne.0.and.iPDFvar.ne.1)then
         write(*,*)'You should type 0 or 1, not ',iPDFvar
         stop
      endif
c
      npdf=0
      str='no PDF'
      if(iPDFvar.eq.1)then
         write (*,*) 'Give minimum and maximum uncertainty PDF sets'
         write (*,*)'(like e.g. 21101 21140)'
         read(*,*)ipdfmin,ipdfmax
         npdf=abs(ipdfmax-ipdfmin)+1
         do j=1,npdf
            write(str,'(i8)')ipdfmin-1+j
            pdf(j)='PDF='//str
         enddo
      endif
c
      if(iscalevar.eq.1)then
         outfile_s=outfilebase(1:lb)//'_scale.top'
         outfile_s_central=outfilebase(1:lb)//'_scale_central.top'
         outfile_s_min=outfilebase(1:lb)//'_scale_min.top'
         outfile_s_max=outfilebase(1:lb)//'_scale_max.top'
         open(unit=uouts,file=outfile_s,status='unknown')
         open(unit=uouts_cen,file=outfile_s_central,status='unknown')
         open(unit=uouts_min,file=outfile_s_min,status='unknown')
         open(unit=uouts_max,file=outfile_s_max,status='unknown')
         call write_top_start(uouts_cen)
         call write_top_start(uouts_min)
         call write_top_start(uouts_max)
         call write_top_start(uouts)
      endif
      if(iPDFvar.eq.1)then
         outfile_pdf=outfilebase(1:lb)//'_pdf.top'
         outfile_pdf_central=outfilebase(1:lb)//'_pdf_central.top'
         outfile_pdf_min=outfilebase(1:lb)//'_pdf_min.top'
         outfile_pdf_max=outfilebase(1:lb)//'_pdf_max.top'
         open(unit=uoutpdf,file=outfile_pdf,status='unknown')
         open(unit=uoutpdf_cen,file=outfile_pdf_central,status=
     &                                                     'unknown')
         open(unit=uoutpdf_min,file=outfile_pdf_min,status='unknown')
         open(unit=uoutpdf_max,file=outfile_pdf_max,status='unknown')
         call write_top_start(uoutpdf_cen)
         call write_top_start(uoutpdf_min)
         call write_top_start(uoutpdf_max)
         call write_top_start(uoutpdf)
      endif
c
      var='first variable to find'
      done=.false.
      do
         rewind(uin)
         call get_reference_var(uin,var,done,i0,i1)
         if(done) exit
         xs(0,1)=0d0
         call get_data(uin,xs(0,1),ys(1,1))
c Scale uncertainties
         if(iscalevar.eq.1)then
            do i=2,10
               call get_var(uin,' (  '//var(1:i1-5)//scale(i),i0,i1)
               xs(0,i)=0d0
               call get_data(uin,xs(0,i),ys(1,i))
               if (xs(0,i).ne.xs(0,1)) then
                  write (*,*) 'data not consistent... skipping plot'
                  cycle
               endif
            enddo
         endif
c PDF uncertainties
         if(iPDFvar.eq.1)then
            do i=1,npdf
               call get_var(uin,' (  '//var(1:i1-5)//pdf(i),i0,i1)
               xpdf(0,i)=0d0
               call get_data(uin,xpdf(0,i),ypdf(1,i))
               if (xpdf(0,i).ne.xs(0,1)) then
                  write (*,*) 'data not consistent... skipping plot'
                  cycle
               endif
            enddo
         endif
c Now I should have all data
c Compute max/min scale
         if(iscalevar.eq.1)then
            do i=1,nint(xs(0,1))
               ysmax(i)=-1d99
               ysmin(i)=+1d99
               do j=2,10
                  ysmax(i)=max( ysmax(i),ys(i,j) )
                  ysmin(i)=min( ysmin(i),ys(i,j) )
               enddo
            enddo
            call write_top_header(uouts_cen,var)
            call write_top_header(uouts_min,var)
            call write_top_header(uouts_max,var)
            call write_top_header(uouts,var)
            call write_top_data(uouts_cen,xs(0,1),ys(1,1))
            call write_top_data(uouts_min,xs(0,1),ysmin)
            call write_top_data(uouts_max,xs(0,1),ysmax)
            call write_top_data2(uouts,xs(0,1),ys(1,1),ysmax,ysmin)
         endif
c Compute max/min pdf according to hessian method
         if(iPDFvar.eq.1)then
            do j=1,nint(xpdf(0,1))
               ypdfmax(j)=0d0
               ypdfmin(j)=0d0
               do i=1,20
                  nps=2*i-1
                  nng=2*i
                  ypdfmin(j)=ypdfmin(j)+( max(0.d0,
     &                 ys(j,1)-ypdf(j,nps),
     &                 ys(j,1)-ypdf(j,nng)) )**2
                  ypdfmax(j)=ypdfmax(j)+( max(0.d0,
     &                 ypdf(j,nps)-ys(j,1),
     &                 ypdf(j,nng)-ys(j,1)) )**2
               enddo
               ypdfmin(j)=ys(j,1)-sqrt(ypdfmin(j))
               ypdfmax(j)=ys(j,1)+sqrt(ypdfmax(j))
            enddo
            call write_top_header(uoutpdf_cen,var)
            call write_top_header(uoutpdf_min,var)
            call write_top_header(uoutpdf_max,var)
            call write_top_header(uoutpdf,var)
            call write_top_data(uoutpdf_cen,xs(0,1),ys(1,1))
            call write_top_data(uoutpdf_min,xs(0,1),ypdfmin)
            call write_top_data(uoutpdf_max,xs(0,1),ypdfmax)
            call write_top_data2(uoutpdf,xs(0,1),ys(1,1),ypdfmax,
     &                                                   ypdfmin)
         endif
      enddo
c
      close (uin)
      if(iscalevar.eq.1)close (uouts)
      if(iPDFvar.eq.1)close (uoutpdf)
c
      return
      end

      subroutine get_reference_var(uin,var,done,i0,i1)
      implicit none
      integer uin
      character*100 buff,var
      logical found_original,done
      logical verbose
      parameter (verbose=.false.)
      integer i0,i1
      character*100 lines(200)
      integer nlines
      common/to_lines/lines, nlines
      nlines=0
      if (var.eq.'first variable to find') then
         found_original=.true.
      else
         found_original=.false.
      endif
      do
         read(uin,'(a)',end=999) buff
         i0=index(buff,'(')
         i1=index(buff,'central')
         if (found_original) then
            if (i0.ne.0.and.i1.ne.0)then
               var=buff(i0+3:i1-1)
               exit
            else
               nlines = nlines+1
               lines(nlines) = buff(1:100)
               if (index(buff,'PLOT').ne.0) nlines = 0
               if (index(buff,'NEW PLOT').ne.0) nlines = 1
            endif
         else
            if(buff(1:i1-1).eq.' (  '//var(1:i1-5))then
               found_original=.true.
            else
               nlines = nlines+1
               lines(nlines) = buff(1:100)
               if (index(buff,'PLOT').ne.0) nlines = 0
               if (index(buff,'NEW PLOT').ne.0) nlines = 1
            endif
         endif
      enddo
      read(uin,'(a)') buff
      return
 999  done=.true.
      return
      end

      subroutine get_var(uin,var,i0,i1)
      implicit none
      integer uin,i0,i1
      integer nc
      parameter(nc=15)
      logical verbose
      parameter (verbose=.false.)
      character*100 var,buff
      do
         read(uin,'(a)',end=998,err=998) buff
         if (buff(1:i1+nc-1).eq.var(1:i1+nc-1)) exit
      enddo
      read(uin,'(a)') buff
      return
 998  write (*,*) 'ERROR var not found'
      if(verbose) write (*,*) var
      stop
      end

      subroutine get_data(uin,x,y)
      implicit none
      integer uin
      double precision x(0:1000),y(1000),dy(1000)
      character*100 buff
      do
         read(uin,'(a)')buff
         if(index(buff,'HIST').eq.0) then
            x(0)=x(0)+1d0
            read(buff,*),x(nint(x(0))),y(nint(x(0))),dy(nint(x(0)))
         else
            exit
         endif
      enddo
      return
      end

      subroutine write_top_start(unit)
      implicit none
      integer unit
      write (unit,'(a)') 'SET DEVICE POSTSCRIPT ORIENT=3'
      write (unit,'(a)') 'SET FONT DUPLEX'
      return
      end

      subroutine write_top_header(unit,var)
      implicit none
      integer unit
      character*100 var
      character*100 lines(200)
      logical verbose
      parameter (verbose=.false.)
      integer nlines,i
      common/to_lines/lines, nlines
      write(unit,'(a)') 'PLOT'
      do i = 1, nlines
        if(verbose) write(*,*) lines(i)
        write(unit,'(a)') trim(lines(i))
      enddo
      write(unit,'(a)') trim(' (  '//var)
      return
      end

      subroutine write_top_data(unit,x,y)
      implicit none
      integer unit,i
      double precision x(0:1000),y(1000)
      do i=1,nint(x(0))
         write(unit,999) x(i),y(i),0d0
      enddo
      write (unit,'(a)') 'HIST SOLID'
      return
999   format(3(d16.6))
      end

      subroutine write_top_data2(unit,x,y,ymax,ymin)
      implicit none
      integer unit,i
      double precision x(0:1000),y(1000),ymax(1000),ymin(1000)
      do i=1,nint(x(0))
         write(unit,999) x(i),y(i),0d0
      enddo
      write (unit,'(a)') 'HIST SOLID RED'
      do i=1,nint(x(0))
         write(unit,999) x(i),ymax(i),0d0
      enddo
      write (unit,'(a)') 'HIST SOLID'
      do i=1,nint(x(0))
         write (unit,999) x(i),ymin(i),0d0
      enddo
      write (unit,'(a)') 'HIST SOLID'
      write (unit,'(a)') ''
      return
999   format(3(d16.6))
      end

max_stars_repo_name: JaySandesara/madminer

max_issues_repo_name: JaySandesara/madminer

max_forks_repo_name: JaySandesara/madminer


text:
Since Davis is landlocked, assume that this term refers to the illegal downloading of copyrighted material. Otherwise, be sure to see the worlds only http://www.826valencia.org/store/index.html pirate supply store next time youre in San Francisco or the Ninjas vs. Pirates at the MU.

The piracy method of choice for UC Davis is wiki:WikiPedia:BitTorrent BitTorrent (http://www.infoanarchy.org/wiki/wiki.pl?BitTorrent iA) and wiki:WikiPedia:Kazaa KaZaA (http://www.infoanarchy.org/wiki/wiki.pl?KaZaA iA).

The wiki:WikiPedia:MPAA MPAA (http://www.infoanarchy.org/wiki/wiki.pl?MPAA iA) has not yet sued any UC Davis students for movie piracy.

Aggie article on wiki:WikiPedia:RIAA RIAA (http://www.infoanarchy.org/wiki/wiki.pl?RIAA iA) http://www.californiaaggie.com/media/storage/paper981/news/2005/06/27/FrontPage/RiaaTo.Sue.Four.Ucd.Students1320868.shtml suing four UCD students.
   Anyone know the names of the RIAA Four? Are all still attending UCD?
   I heard about one of them (or at least a similar incident) a couple of months ago, though I dont have a name. An acquaintance of mine told me that a high school student had been living in the dorms and working on some kind of research project for future UC students. He was working in the same lab as my acquaintance, and he downloaded some stuff over the school network. He was caught and, after a discussion with administration and some corporate lawyers, was tentatively given a stern warning by both. According to my acquaintance, though, if the industry guys decide to press charges, the student in question will be blacklisted for life from the entire UC system. Make what you will of the story, but it sounds plausible to me.  RoyWright

max_stars_repo_name: voflo/Search

max_issues_repo_name: voflo/Search

max_forks_repo_name: voflo/Search


text:
The term FAQ or Frequently Asked Questions predates the World Wide Web.  Long before the web existed, the internet had other ways of being used.  One of the more popular was usenet, a message system that predates the web, and still exists today.  Popular at universities, most universities have a usenet server and university specific newsgroups.  UC Davis is no exception, and the ucd. series of groups have their own FAQ list.  news.ucdavis.edu is UC Davis usenet server.   Although the ucd. groups have faded, there is an historical archive of USENET at http://groupsbeta.google.com/ and you can use it to see the huge extent to which the UC Davis groups were used.

Last updated in 1996, the FAQ can be http://faqs.netesq.com/ucdavis/part1.html found here.

In it you can find thrillingly dated questions such as Why should I use the Internet?

20041221 21:06:43 nbsp Yeah, the Users/DavidPrenatt guy who wrote this FAQ went to UCD Law school and is an interesting guy.  He makes his living selling technology solutions to law firms.  I asked him about updating it once, and he said that I should by all means.  This site has filled the void left by the outofdate FAQ, however.  Even so, it would be interesting to have a pseudoFAQ where people post questions and its answered by referring people to certain links on this site.  Users/JaimeRaba

20051204 13:59:03 nbsp This page ought to be renamed to something like UCDavis Usenet FAQ? Users/AlexanderWoo


text:
Davis Pet Nanny is a local Pet Sitting / Pet Care business, serving the Davis, wiki:Woodland, wiki:Winters & Dixon areas.  Specializing in daily, weekly and outoftown care of your furry, feathered, finned and four legged loved ones.  Licensed, insured & bonded.  Animal first aid certified.  Personalized care to meet the needs of your pets.  Services include daily dog walks, exercise & socialization.  Cat care, birds, fish, rodents & exotics too!  Farm animal care.  Plant watering, mail collection and security checks while youre out of town.

20080830 18:52:23 nbsp We recently had Robin take care of our kitty while we were away for 10 days.  She did an excellent job!  She is accommodating, kind, and communicated with us through out our trip.  She also left a cute picture of him, which we thought was great touch.  And to top it off, our kitty liked her, which doesnt always happen.  We feel safe leaving for a few days now that we have a wonderful pet nanny!   Users/keriH

20120618 13:55:20 nbsp Valued customers.  As of March 2012, I am no longer offering pet sitting services in Davis.  I have relocated from the area. There are some other wonderful pet care providers in the area that can assist you, such as Elizabeths Pet Sitting, Pitter Patter Pet Care & Pet Pals Davis.  Thank you for the wonderful years of service of your pets needs!  I will miss Davis!  Robin Ray Users/rray

20120619 00:26:54 nbsp Well miss you, Robin!  Users/ScottMeehleib