#include "w3macros.h"
!/ ------------------------------------------------------------------- /
      PROGRAM W3OUTP
!/
!/                  +-----------------------------------+
!/                  | WAVEWATCH III           NOAA/NCEP |
!/                  |           H. L. Tolman            |
!/                  |            J.H. Alves             |
!/                  |             A. Chawla             |
!/                  |            F. Ardhuin             |
!/                  |             E. Rogers             |
!/                  |            T. Campbell            |
!/                  |                        FORTRAN 90 |
!/                  | Last update :         27-Aug-2015 |
!/                  +-----------------------------------+
!/
!/    14-Jan-1999 : Final FORTRAN 77                    ( version 1.18 )
!/    21-Jan-2000 : Upgrade to FORTRAN 90               ( version 2.00 )
!/    14-Feb-2000 : Exact nonlinear interactions        ( version 2.01 )
!/    09-Jan-2001 : U* bug fix in tabular output        ( version 2.05 )
!/    25-Jan-2001 : Flat grid version                   ( version 2.06 )
!/    02-Feb-2001 : Xnl version 3.0                     ( version 2.07 )
!/    11-Jun-2001 : Clean up                            ( version 2.11 )
!/    11-Oct-2001 : Clean up, X*, Y* in tables          ( version 2.14 )
!/    13-Nov-2002 : Add stress vector                   ( version 3.00 )
!/    27-Nov-2002 : First version of VDIA and MDIA      ( version 3.01 )
!/    24-Dec-2004 : Multiple grid version.              ( version 3.06 )
!/    17-Apr-2006 : Filter for directional spread.      ( version 3.09 )
!/    23-Jun-2006 : Linear input added.                 ( version 3.09 )
!/    28-Jun-2006 : Adding file name preamble.          ( version 3.09 )
!/    03-Jul-2006 : Separate flux modules.              ( version 3.09 )
!/    28-Oct-2006 : Add partitioning option.            ( version 3.10 )
!/    24-Mar-2007 : Add pars for entire spectrum.       ( version 3.11 )
!/    25-Apr-2007 : Battjes-Janssen Sdb added.          ( version 3.11 )
!/                  (J. H. Alves)
!/    08-Aug-2007 : Creation of buoy log file added     ( version 3.12 )
!/                  (switch O14 -- A. Chawla)
!/    09-Oct-2007 : WAM 4+ Sin and Sds added.           ( version 3.13 )
!/                  (F. Ardhuin)
!/    09-Oct-2007 : Experimental Sbs (BS1) added.       ( version 3.13 )
!/                  (F. Ardhuin)
!/    09-Apr-2008 : Adding an additional output for     ( version 3.12 )
!/                  WMO standard (A. Chawla)
!/    29-Apr-2008 : Adjust format partition output.     ( version 3.14 )
!/    29-May-2009 : Preparing distribution version.     ( version 3.14 )
!/    30-Oct-2009 : Implement run-time grid selection.  ( version 3.14 )
!/                  (W. E. Rogers & T. J. Campbell, NRL)
!/    04-Mar-2010 : Added partitions bulletin output.   ( version 3.14 )
!/                  (J. H. Alves)
!/    20-Apr-2010 : Fix initialization of USTAR.      ( version 3.14.1 )
!/    16-Jul-2012 : Move GMD (SNL3) and nonlinear filter (SNLS)
!/                  from 3.15 (HLT).                    ( version 4.08 )
!/    23-Aug-2012 : Adding movable bed friction BT4     ( version 4.08 )
!/    26-Dec-2012 : Modified obsolete declarations.     ( version 4.11 )
!/    10-Sep-2013 : Implement second order correction   ( version 4.12 )
!/                  (F. Ardhuin)
!/    06-Feb-2014 : Fix header format in part. files.   ( version 4.18 )
!/    27-Aug-2015 : Sice add as additional output       ( version 5.10 )
!/                  (in source terms)
!/    27-Jun-2017 : Expanding WMO table to 2 digits JHA ( version 6.02 )
!/    18-Aug-2018 : S_{ice} IC5 (Q. Liu)                ( version 6.06 )
!/
!/    Copyright 2009-2014 National Weather Service (NWS),
!/       National Oceanic and Atmospheric Administration.  All rights
!/       reserved.  WAVEWATCH III is a trademark of the NWS. 
!/       No unauthorized use without permission.
!/
!  1. Purpose :
!
!     Post-processing of point output.
!
!  2. Method :
!
!     Data is read from the grid output file out_pnt.ww3 (raw data)
!     and from the file ww3_outp.inp ( NDSI, output requests ).
!     Model definition and raw data files are read using WAVEWATCH III
!     subroutines.
!
!     Output types ITYPE :            Sub-type OTYPE :
!     --------------------           -----------------
!       0 : Check file.
!       1 : Spectra.
!                                      1 : Print plots.
!                                      2 : Table of 1-D spectra
!                                      3 : Transfer file
!       2 : Table of mean wave parameters
!                                      1 : Depth, current,  wind
!                                      2 : Mean wave pars.
!                                      3 : Nondimensional pars. (U*)
!                                      4 : Nondimensional pars. (U10)
!                                      5 : Validation table
!                                      6 : WMO standard output
!       3 : Source terms
!                                      1 : Print plots.
!                                      2 : Table of 1-D S(f).
!                                      3 : Table of 1-D time scales.
!                                      4 : Transfer file.
!
!       4 : Partitioning and bulletins
!                                      1 : Spectral partitions table
!                                      2 : Bulletins ASCII format
!                                      3 : Bulletins CSV format
!                                      4 : Bulletins CSV & ASCII format
!  3. Parameters :
!
!  4. Subroutines used :
!
!      Name      Type  Module   Description
!     ----------------------------------------------------------------
!      W3NMOD    Subr. W3GDATMD Set number of model.
!      W3SETG    Subr.   Id.    Point to selected model.
!      W3NDAT    Subr. W3WDATMD Set number of model for wave data.
!      W3SETW    Subr.   Id.    Point to selected model for wave data.
!      W3NAUX    Subr. W3ADATMD Set number of model for aux data.
!      W3SETA    Subr.   Id.    Point to selected model for aux data.
!      W3NOUT    Subr. W3ODATMD Set number of model for output.
!      W3SETO    Subr.   Id.    Point to selected model for output.
!      ITRACE    Subr. W3SERVMD Subroutine tracing initialization.
!      STRACE    Subr.   Id.    Subroutine tracing.
!      NEXTLN    Subr.   Id.    Get next line from input filw
!      EXTCDE    Subr.   Id.    Abort program as graceful as possible.
!      STME21    Subr. W3TIMEMD Convert time to string.
!      TICK21    Subr.   Id.    Advance time.
!      DSEC21    Func.   Id.    Difference between times.
!      W3IOGR    Subr. W3IOGRMD Reading/writing model definition file.
!      W3IOPO    Subr. W3IOPOMD Reading/writing raw point output file.
!      W3EXPO    Subr. Internal Execute point output.
!      W3BULL    Subr. W3BULLMD Generate buletins from spectral part.
!     ----------------------------------------------------------------
!
!  5. Called by :
!
!     None, stand-alone program.
!
!  6. Error messages :
!
!     Checks on input, checks in W3IOxx.
!
!  7. Remarks :
!
!     - Tables written to file 'tabNN.ww3', where NN is the
!       unit umber (NDSTAB).
!     - Transfder file written to ww3.yymmddhh.spc with multiple
!       spectra and times in file. yymmddhh relates to first
!       output (NDSTAB).
!     - !/IC1 !/IC2 !/IC3 !/IC4 !/IC5 are not included in dissipation term
!       FIXME: ICE is a dummy variable at the moment
!              Include ice parameters in point output file out_pnt.ww3
!              Ice coupling to SIN, SDS and SIC similar to w3srcemd.ftn
!
!  8. Structure :
!
!     See source code.
!
!  9. Switches :
!
!       !/S    Enable subroutine tracing.
!
!       !/NCO  NCEP NCO modifications for operational implementation.
!
!       !/O14  Buoy log file generation.
!
! 10. Source code :
!
!/ ------------------------------------------------------------------- /
      USE CONSTANTS
!/
!     USE W3GDATMD, ONLY: W3NMOD, W3SETG
      USE W3WDATMD, ONLY: W3SETW, W3NDAT
!/NL1      USE W3ADATMD, ONLY: W3SETA, W3NAUX
      USE W3ODATMD, ONLY: W3SETO, W3NOUT
      USE W3IOGRMD, ONLY: W3IOGR
      USE W3IOPOMD, ONLY: W3IOPO
      USE W3SERVMD, ONLY : ITRACE, NEXTLN, EXTCDE
!/S      USE W3SERVMD, ONLY : STRACE
      USE W3TIMEMD, ONLY: STME21, TICK21, DSEC21
!/
      USE W3GDATMD
      USE W3WDATMD, ONLY: TIME
      USE W3ODATMD, ONLY: NDSE, NDST, NDSO, NOPTS, PTLOC, PTNME,     &
                          DPO, WAO, WDO, ASO, CAO, CDO, SPCO, FNMPRE,&
                          ICEO, ICEHO, ICEFO, DIMP
      USE W3BULLMD, ONLY: NPTAB, NFLD, NPMAX, BHSMIN, BHSDROP, IYY,  &
                          HST, TPT, DMT, ASCBLINE, CSVBLINE
!/NCO       USE W3BULLMD, ONLY: CASCBLINE
!/O14      USE W3ODATMD, ONLY: GRDID
!/IG1       USE W3GIG1MD, ONLY: W3ADDIG
!/IG1       USE W3CANOMD, ONLY: W3ADD2NDORDER
!
      IMPLICIT NONE
!/
!/ ------------------------------------------------------------------- /
!/ Local parameters
!/
      INTEGER                 :: NDSI, NDSM, NDSOP,  NDSTRC, NTRACE,  &
                                 IERR, I, TOUT(2), NOUT, TDUM(2),     &
                                 NREQ, IPOINT, ITYPE, OTYPE, NDSTAB,  &
                                 IOTEST, IK, ITH, IOUT, J, DIMXP,     &
                                 NDSBUL, NDSCSV, ICSV, IJ
!/NCO       INTEGER                 :: NDSCBUL
      INTEGER                 :: ISCALE = 0
      INTEGER                 :: TIMEV(2)
!/O14      INTEGER                 :: NDBO
!/S      INTEGER, SAVE           :: IENT   = 0
      REAL                    :: DTREQ, SCALE1, SCALE2, DTEST
      REAL                    :: M2KM
      REAL, ALLOCATABLE       :: XPART(:,:)
      LOGICAL                 :: FLFORM, FLSRCE(7)
      LOGICAL, ALLOCATABLE    :: FLREQ(:)
      CHARACTER               :: COMSTR*1, IDTIME*23, IDDDAY*11,      &
                                 TABNME*9, TFNAME*64
      CHARACTER(LEN=25)       :: IDSRCE(7)
      CHARACTER               :: HSTR*6, HTYPE*3
      CHARACTER(LEN=32)       :: prefix
      integer                 :: nargs
      LOGICAL                 :: PROCESS_POINT_ONLY
      INTEGER                 :: ACTIVE_POINT, J_START, J_END
!/BT2      REAL                    :: SBTC2
!/
!/ ------------------------------------------------------------------- /
!/
      DATA IDSRCE / 'Spectrum                 ' ,                     &
                    'Wind-wave interactions   ' ,                     &
                    'Nonlinear interactions   ' ,                     &
                    'Dissipation              ' ,                     &
                    'Wave-bottom interactions ' ,                     &
                    'Wave-ice interactions    ' ,                     &
                    'Sum of selected sources  ' /
      FLSRCE = .FALSE.

      nargs=command_argument_count()
      if(nargs==0) then
        prefix="gfswave."
      elseif(nargs==1) then
        call get_command_argument(1,prefix)
        prefix=trim(prefix)//'.'
      else
        WRITE(NDSO, '(A)') 'ww3_outp: Argument count not recognized.'
        WRITE(NDSO, '(A)') 'ww3_outp: Assuming first arg is experiment prefix.'
        call get_command_argument(1,prefix)
        prefix=trim(prefix)//'.'
      endif
!
!/NCO/!     CALL W3TAGB('WAVESPEC',1998,0007,0050,'NP21   ')
!
!--- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! 1.  IO set-up.
!
      CALL W3NMOD ( 1, 6, 6 )
      CALL W3SETG ( 1, 6, 6 )
      CALL W3NDAT (    6, 6 )
      CALL W3SETW ( 1, 6, 6 )
!/NL1      CALL W3NAUX (    6, 6 )
!/NL1      CALL W3SETA ( 1, 6, 6 )
      CALL W3NOUT (    6, 6 )
      CALL W3SETO ( 1, 6, 6 )
!
      NDSI   = 10
      NDSM   = 20
      NDSOP  = 20
      NDSBUL = 0
!/NCO       NDSCBUL = 0
!
      NDSTRC =  6
      NTRACE = 10
      CALL ITRACE ( NDSTRC, NTRACE )

!
!/S      CALL STRACE (IENT, 'W3OUTP')
!
!/NCO/!
!/NCO/! Redo according to NCO
!/NCO/!
!/NCO      NDSI   = 11
!/NCO      NDSO   =  6
!/NCO      NDSE   = NDSO
!/NCO      NDST   = NDSO
!/NCO      NDSM   = 12
!/NCO      NDSOP  = 13
!/O14      NDBO   = 14
!/NCO      NDSTRC = NDSO
!
      WRITE (NDSO,900)
!
      J      = LEN_TRIM(FNMPRE)
      OPEN (NDSI,FILE=FNMPRE(:J)//'ww3_outp.inp',STATUS='OLD',ERR=800,IOSTAT=IERR)
      READ (NDSI,'(A)',END=801,ERR=802) COMSTR
      IF (COMSTR.EQ.' ') COMSTR = '$'
      WRITE (NDSO,901) COMSTR
!
!--- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! 2.  Read model definition file.
!
      CALL W3IOGR ( 'READ', NDSM )
      WRITE (NDSO,920) GNAME
!
      IF ( FLAGLL ) THEN
          M2KM = 1.
        ELSE
          M2KM = 1.E-3
        END IF
!
      DIMXP  = ((NK+1)/2) * ((NTH-1)/2)
      ALLOCATE ( XPART(DIMP,0:DIMXP) )
      XPART  = UNDEF
!
!--- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! 3.  Read requests from input file.
!     Output times
!
      CALL NEXTLN ( COMSTR , NDSI , NDSE )
      READ (NDSI,*,END=801,ERR=802) TOUT, DTREQ, NOUT
      DTREQ  = MAX ( 0. , DTREQ )
      IF ( DTREQ.EQ.0 ) NOUT = 1
      NOUT   = MAX ( 1 , NOUT )

      CALL W3IOPO_dynm ( NDSOP, TOUT, IOTEST )
      WRITE (NDSO,930)
      DO I=1, NOPTS
        IF ( FLAGLL ) THEN
          WRITE (NDSO,931) PTNME(I), M2KM*PTLOC(1,I), M2KM*PTLOC(2,I)
        ELSE
          WRITE (NDSO,932) PTNME(I), M2KM*PTLOC(1,I), M2KM*PTLOC(2,I)
        END IF
      END DO
!
      CALL STME21 ( TOUT , IDTIME )
      WRITE (NDSO,940) IDTIME
!
      TDUM   = 0
      CALL TICK21 ( TDUM , DTREQ )
      CALL STME21 ( TDUM , IDTIME )
      IF ( DTREQ .GE. 86400. ) THEN
          WRITE (IDDDAY,'(I10,1X)') INT(DTREQ/86400.)
        ELSE
          IDDDAY = '           '
        END IF
      IDTIME = IDDDAY
      IDTIME(21:23) = '   '
      WRITE (NDSO,941) IDTIME, NOUT
!
! ... Output points
!
      ALLOCATE ( FLREQ(NOPTS) )
      FLREQ = .FALSE.
      NREQ   = 0
!
      DO I=1, NOPTS
        ! reads point index
        CALL NEXTLN ( COMSTR , NDSI , NDSE )
        READ (NDSI,*,END=801,ERR=802) IPOINT
        ! last index
        IF (IPOINT .LT. 0) THEN
          IF (I.EQ.1) THEN
            FLREQ = .TRUE.
            NREQ = NOPTS
          END IF
          EXIT
        END IF
        ! existing index in out_pnt.ww3   
        IF ( (IPOINT .GT. 0) .AND. (IPOINT .LE. NOPTS) ) THEN
          IF ( .NOT. FLREQ(IPOINT) ) THEN 
            NREQ = NREQ + 1
          END IF
          FLREQ(IPOINT) = .TRUE.
        END IF
        ! read the 'end of list' if nopts reached before it
        IF ( (IPOINT .GT. 0) .AND. (NREQ .EQ. NOPTS) ) THEN
          CALL NEXTLN ( COMSTR , NDSI , NDSE )
          READ (NDSI,*,END=801,ERR=802) IPOINT
        END IF
      END DO
      ! check if last point index is -1
      IF (IPOINT .NE. -1) THEN
        WRITE (NDSE,1007)
        CALL EXTCDE ( 47 )
      END IF

!
! ... Output type
!
      CALL NEXTLN ( COMSTR , NDSI , NDSE )
      READ (NDSI,*,END=801,ERR=802) ITYPE
!
! ... ITYPE = 0
!
      IF ( ITYPE .EQ. 0 ) THEN
!
!/O14          WRITE (NDSO,942) ITYPE, 'Generating buoy log file'
!/O14          OPEN (NDBO,FILE=FNMPRE(:J)//'buoy_log.ww3',            &
!/O14                STATUS='NEW',ERR=805,IOSTAT=IERR)
!/O14          DO I = 1,NOPTS
!/O14            WRITE(NDBO,945) I, PTNME(I), PTLOC(1,I),             &
!/O14                            PTLOC(2,I), GRDID(I)
!/O14            END DO
!/O14          CLOSE(NDBO)
!
          WRITE (NDSO,942) ITYPE, 'Checking contents of file'
          DO
            CALL STME21 ( TIME , IDTIME )
            WRITE (NDSO,948) IDTIME
            CALL W3IOPO ( 'READ', NDSOP, IOTEST )
            IF ( IOTEST .EQ. -1 ) THEN
                WRITE (NDSO,949)
                GOTO 888
              END IF
            END DO
!
! ... ITYPE = 1
!
        ELSE IF (ITYPE .EQ. 1) THEN
          WRITE (NDSO,942) ITYPE, '1-D and/or 2-D spectra'
          CALL NEXTLN ( COMSTR , NDSI , NDSE )
          READ (NDSI,*,END=801,ERR=802) OTYPE, SCALE1, SCALE2,        &
                NDSTAB, FLFORM
!/NCO          NDSTAB = 51
          IF (OTYPE .EQ. 1) THEN
              WRITE (NDSO,943) 'print plots'
              IF ( SCALE1 .LT. 0.  ) THEN
                  WRITE (NDSO,1940) '1-D'
                ELSE IF ( SCALE1 .EQ. 0.  ) THEN
                  WRITE (NDSO,1941) '1-D'
                ELSE
                  WRITE (NDSO,1942) '1-D', SCALE1
                END IF
              IF ( SCALE2 .LT. 0.  ) THEN
                  WRITE (NDSO,1940) '2-D'
                ELSE IF ( SCALE2 .EQ. 0.  ) THEN
                  WRITE (NDSO,1941) '2-D'
                ELSE
                  WRITE (NDSO,1942) '2-D', SCALE2
                END IF
            ELSE IF ( OTYPE .EQ. 2 ) THEN
              WRITE (NDSO,943) 'Table of 1-D spectral data'
              TABNME = 'tab--.ww3'
              IF ( NDSTAB.LE.0 .OR. NDSTAB.GT.99 ) NDSTAB = 51
              WRITE ( TABNME(4:5) , '(I2.2)' ) NDSTAB
              J      = LEN_TRIM(FNMPRE)
              OPEN (NDSTAB,FILE=FNMPRE(:J)//TABNME,ERR=803,IOSTAT=IERR)
              WRITE (NDSO,1947) TABNME
            ELSE IF ( OTYPE .EQ. 3 ) THEN
              WRITE (NDSO,943) 'Transfer file for each point'
              DO IJ = 1, NOPTS
                  IF (FLREQ(IJ)) THEN
                      TFNAME = TRIM(prefix)//TRIM(PTNME(IJ))//'.spec'
                      WRITE (NDSO,1943) TRIM(TFNAME), 'Transfer File'
                      J = LEN_TRIM(FNMPRE)
                      IF (FLFORM) THEN
                          OPEN(NDSTAB, FILE=TRIM(TFNAME), ERR=804, &
                               IOSTAT=IERR, FORM='UNFORMATTED')
                      WRITE (NDSTAB) 'WAVEWATCH III SPECTRA',     &
                                     NK, NTH, 1, GNAME
                      WRITE (NDSTAB) (SIG(IK)*TPIINV, IK = 1, NK)
                      WRITE (NDSTAB) (MOD(2.5*PI-TH(ITH), TPI), ITH = 1, NTH)
                  ELSE
                      OPEN(NDSTAB, FILE=TRIM(TFNAME), ERR=804, &
                           IOSTAT=IERR, FORM='FORMATTED')
                      WRITE (NDSTAB,1944) 'WAVEWATCH III SPECTRA', &
                                           NK, NTH, 1, GNAME
                      WRITE (NDSTAB,1945) (SIG(IK)*TPIINV, IK = 1, NK)
                      WRITE (NDSTAB,1946) (MOD(2.5*PI-TH(ITH), TPI), ITH= 1, NTH)
                  END IF
                  CLOSE(NDSTAB)
              END IF
              END DO 

            
            ELSE
              WRITE (NDSE,1011) OTYPE
              CALL EXTCDE ( 10 )
            END IF
!
! ... ITYPE = 2
!
        ELSE IF (ITYPE .EQ. 2) THEN
          WRITE (NDSO,942) ITYPE, 'Table of mean wave parameters'
          CALL NEXTLN ( COMSTR , NDSI , NDSE )
          READ (NDSI,*,END=801,ERR=802) OTYPE, NDSTAB
!/NCO          NDSTAB = 51
          TABNME = 'tab--.ww3'
          IF ( NDSTAB.LE.0 .OR. NDSTAB.GT.99 ) NDSTAB = 51
          WRITE ( TABNME(4:5) , '(I2.2)' ) NDSTAB
          J      = LEN_TRIM(FNMPRE)
          OPEN (NDSTAB,FILE=FNMPRE(:J)//TABNME,ERR=803,IOSTAT=IERR)
          IF ( OTYPE .EQ. 1 ) THEN
              WRITE (NDSO,2940) 'Depth, current and wind', TABNME
            ELSE IF ( OTYPE .EQ. 2 ) THEN
              WRITE (NDSO,2940) 'Mean wave parameters', TABNME
            ELSE IF ( OTYPE .EQ. 3 ) THEN
              WRITE (NDSO,2940) 'Nondimensional parameters (U*)',     &
                                 TABNME
            ELSE IF ( OTYPE .EQ. 4 ) THEN
              WRITE (NDSO,2940) 'Nondimensional parameters (U10)',    &
                                 TABNME
            ELSE IF ( OTYPE .EQ. 5 ) THEN
              WRITE (NDSO,2940) 'Validation parameters', TABNME
            ELSE IF ( OTYPE .EQ. 6 ) THEN
              WRITE (NDSO,2940) 'WMO standard mean parameters', TABNME
            ELSE
              WRITE (NDSE,1011) OTYPE
              CALL EXTCDE ( 20 )
            END IF
!
! ... ITYPE = 3
!
        ELSE IF (ITYPE .EQ. 3) THEN
          WRITE (NDSO,942) ITYPE, 'Source terms'
          CALL NEXTLN ( COMSTR , NDSI , NDSE )
          READ (NDSI,*,END=801,ERR=802) OTYPE, SCALE1, SCALE2,        &
                                        NDSTAB, FLSRCE, ISCALE, FLFORM
!/NCO          NDSTAB = 51
          ISCALE = MAX ( 0 , MIN ( 5 , ISCALE ) )
          IF ( OTYPE .EQ. 1 ) THEN
              WRITE (NDSO,943) 'Print plots'
            ELSE IF ( OTYPE .EQ. 2 ) THEN
              IF ( ISCALE .LE. 2) THEN
                  WRITE (NDSO,943) 'Tables as a function of freq.'
                ELSE
                  WRITE (NDSO,943) 'Tables as a function of f/fp.'
                END IF
              IF ( MOD(ISCALE,3) .EQ. 1 ) THEN
                  WRITE (NDSO,944) '(nondimensional based on U10)'
                ELSE IF ( MOD(ISCALE,3) .EQ. 2) THEN
                  WRITE (NDSO,944) '(nondimensional based on U*)'
                END IF
            ELSE IF ( OTYPE .EQ. 3 ) THEN
              IF ( ISCALE .LE. 2) THEN
                  WRITE (NDSO,943) 'Time scales as a function of freq.'
                ELSE
                  WRITE (NDSO,943) 'Time scales as a function of f/fp.'
                END IF
              IF ( ISCALE .EQ. 1 ) THEN
                  WRITE (NDSO,944) '(nondimensional based on U10)'
                ELSE IF ( ISCALE .EQ. 2) THEN
                  WRITE (NDSO,944) '(nondimensional based on U*)'
                END IF
            ELSE IF ( OTYPE .EQ. 4 ) THEN
              TFNAME = 'ww3.--------.src'
              WRITE (TFNAME(5:12),'(I6.6,I2.2)')                      &
                   MOD(TOUT(1),1000000), TOUT(2)/10000
              WRITE (NDSO,943) 'Transfer file'
              IF ( FLFORM ) THEN
                  WRITE (NDSO,3943) TFNAME, 'UNFORMATTED'
                  J      = LEN_TRIM(FNMPRE)
                  OPEN  (NDSTAB,FILE=FNMPRE(:J)//TFNAME,ERR=804,      &
                         IOSTAT=IERR,FORM='UNFORMATTED')
                  WRITE (NDSTAB) 'WAVEWATCH III SOURCES',             &
                                  NK, NTH, NREQ, FLSRCE
                  WRITE (NDSTAB) (SIG(IK)*TPIINV,IK=1,NK)
                  WRITE (NDSTAB) (MOD(2.5*PI-TH(ITH),TPI),ITH=1,NTH)

                ELSE
                  WRITE (NDSO,3943) TFNAME, 'FORMATTED'
                  J      = LEN_TRIM(FNMPRE)
                  OPEN  (NDSTAB,FILE=FNMPRE(:J)//TFNAME,ERR=804,      &
                         IOSTAT=IERR,FORM='FORMATTED')
                  WRITE (NDSTAB,3944) 'WAVEWATCH III SOURCES',        &
                                       NK, NTH, NREQ, FLSRCE
                  WRITE (NDSTAB,3945) (SIG(IK)*TPIINV,IK=1,NK)
                  WRITE (NDSTAB,3946)                                 &
                               (MOD(2.5*PI-TH(ITH),TPI),ITH=1,NTH)
                END IF
            ELSE
              WRITE (NDSE,1011) OTYPE
              CALL EXTCDE ( 30 )
            END IF
!
          DO I=1, 7
            IF ( FLSRCE(I) ) WRITE (NDSO,3940) IDSRCE(I)
            END DO
          WRITE (NDSO,*) ' '
!
          IF ( OTYPE .EQ. 1 ) THEN
              IF ( SCALE1 .LT. 0.  ) THEN
                  WRITE (NDSO,1940) '1-D'
                ELSE IF ( SCALE1 .EQ. 0.  ) THEN
                  WRITE (NDSO,1941) '1-D'
                ELSE
                  WRITE (NDSO,1942) '1-D', SCALE1
                END IF
              IF ( SCALE2 .LT. 0.  ) THEN
                  WRITE (NDSO,1940) '2-D'
                ELSE IF ( SCALE2 .EQ. 0.  ) THEN
                  WRITE (NDSO,1941) '2-D'
                ELSE
                  WRITE (NDSO,1942) '2-D', SCALE2
                END IF
            END IF
!
          IF ( OTYPE.EQ.2 .OR. OTYPE.EQ.3 ) THEN
              TABNME = 'tab--.ww3'
              IF ( NDSTAB.LE.0 .OR. NDSTAB.GT.99 ) NDSTAB = 51
              WRITE ( TABNME(4:5) , '(I2.2)' ) NDSTAB
              J      = LEN_TRIM(FNMPRE)
              OPEN (NDSTAB,FILE=FNMPRE(:J)//TABNME,ERR=803,IOSTAT=IERR)
              WRITE (NDSO,3941) TABNME
            END IF
!
! ... ITYPE = 4
!
        ELSE IF (ITYPE .EQ. 4) THEN
          WRITE (NDSO,942) ITYPE, 'Spectral partitions or bulletins'
          CALL NEXTLN ( COMSTR , NDSI , NDSE )
          READ (NDSI,*,END=801,ERR=802) OTYPE, NDSTAB, TIMEV, HTYPE
!/NCO          NDSTAB = 51
            IF ( OTYPE .EQ. 1 ) THEN
              WRITE (NDSO,943) 'Partitioning of spectra'
              TABNME = 'tab--.ww3'
              IF ( NDSTAB.LE.0 .OR. NDSTAB.GT.99 ) NDSTAB = 51
              WRITE ( TABNME(4:5) , '(I2.2)' ) NDSTAB
              J      = LEN_TRIM(FNMPRE)
              OPEN (NDSTAB,FILE=FNMPRE(:J)//TABNME,ERR=803,IOSTAT=IERR)
              WRITE (NDSO,1947) TABNME

            ELSEIF ( OTYPE .GE. 2 ) THEN
              IF (OTYPE .EQ. 2 .OR. OTYPE .EQ. 4 ) THEN
              WRITE (NDSO,943) 'Bulletins, ASCII format'
              J      = LEN_TRIM(FNMPRE)
              DO IJ = 1,NOPTS
                IF (FLREQ(IJ)) THEN 
                  NDSBUL = NDSTAB + (IJ - 1)
                  OPEN(NDSBUL,FILE=trim(prefix)//TRIM(PTNME(IJ))//&
                                   '.bull',ERR=803,IOSTAT=IERR)
                  WRITE (NDSO,1947) TRIM(PTNME(IJ))//'.bull'
!/NCO                   NDSCBUL = NDSTAB + (IJ - 1) + NOPTS
!/NCO                   OPEN(NDSCBUL,FILE=trim(prefix)//TRIM(PTNME(IJ))//&
!/NCO                                     '.cbull',ERR=803,IOSTAT=IERR)
!/NCO                   WRITE (NDSO,1947) TRIM(PTNME(IJ))//'.cbull'
                ENDIF
              ENDDO
              ENDIF
              IF ( OTYPE .EQ. 3 .OR. OTYPE .EQ. 4 ) THEN
              WRITE (NDSO,943) 'Bulletins, CSV format'
              J      = LEN_TRIM(FNMPRE)
              DO IJ = 1,NOPTS
                IF (FLREQ(IJ)) THEN
                  ICSV = 0
                  NDSBUL = NDSTAB + (IJ - 1)
!/NCO                   NDSCBUL = NDSTAB + (IJ - 1) + NOPTS
                  NDSCSV = NDSTAB + (IJ - 1) + 2*NOPTS
                  OPEN(NDSCSV,FILE=trim(prefix)//TRIM(PTNME(IJ))//&
                                   '.csv',ERR=803,IOSTAT=IERR)
                  WRITE (NDSO,1947) TRIM(PTNME(IJ))//'.csv'
                ENDIF
              ENDDO
              ENDIF
            ELSE
              WRITE (NDSE,1011) OTYPE
              CALL EXTCDE ( 50 )
            END IF
!
! ... ITYPE = ILLEGAL
!
        ELSE
          WRITE (NDSE,1010) ITYPE
          CALL EXTCDE ( 1 )
!
        END IF
!
! ... Output of output points
!
      WRITE (NDSO,950) NREQ
      DO I=1, NOPTS
        IF (FLREQ(I)) THEN
            IF ( FLAGLL ) THEN
                WRITE (NDSO,951) PTNME(I), M2KM*PTLOC(1,I),   &
                                           M2KM*PTLOC(2,I)
              ELSE
                WRITE (NDSO,953) PTNME(I), M2KM*PTLOC(1,I),   &
                                           M2KM*PTLOC(2,I)
              END IF
          END IF
        END DO
!
      IF ( ITYPE.EQ.3 .AND. OTYPE.EQ.4 ) WRITE (NDSO,952)
!
!--- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! 4.  Time management.
     
      PROCESS_POINT_ONLY = .FALSE.
      ACTIVE_POINT = -1

      IOUT   = 0
!
! remark: it would be better to write these warnings only if source term
!         output is requested
!/IC1      WRITE(NDSO,3960)
!/IC2      WRITE(NDSO,3960)
!/IC3      WRITE(NDSO,3960)
!/IC4      WRITE(NDSO,3960)
!/IC5      WRITE(NDSO,3960)

      DO
        DTEST  = DSEC21 ( TIME , TOUT )
        IF ( DTEST .GT. 0. ) THEN
            CALL W3IOPO_dynm ( NDSOP, TOUT, IOTEST )
            IF ( IOTEST .EQ. -1 ) THEN
                WRITE (NDSO,949)
                EXIT
            END IF
            CYCLE
        END IF
        IF ( DTEST .LT. 0. ) THEN
            CALL TICK21 ( TOUT , DTREQ )
            CYCLE
        END IF
!
        IOUT   = IOUT + 1
        CALL STME21 ( TOUT , IDTIME )
        IF ( ( ITYPE.EQ.1 .AND. OTYPE.EQ.1 ) .OR.                     &
             ( ITYPE.EQ.3 .AND. OTYPE.EQ.1 )                          &
             ) WRITE (NDSO,960) IDTIME


        IF (ITYPE .EQ. 1 .AND. OTYPE .EQ. 3) THEN
            DO IJ = 1, NOPTS
                IF (FLREQ(IJ)) THEN
                    TFNAME = TRIM(prefix)//TRIM(PTNME(IJ))//'.spec'
                    J = LEN_TRIM(FNMPRE)
                    IF (FLFORM) THEN
                        OPEN(NDSTAB, FILE=TRIM(TFNAME), STATUS='OLD', &
                             IOSTAT=IERR, FORM='UNFORMATTED', ACCESS='APPEND')
                    ELSE
                        OPEN(NDSTAB, FILE=TRIM(TFNAME), STATUS='OLD', &
                             IOSTAT=IERR, FORM='FORMATTED', ACCESS='APPEND')
                    END IF

                    PROCESS_POINT_ONLY = .TRUE.
                    ACTIVE_POINT = IJ

                    CALL W3EXPO

                    PROCESS_POINT_ONLY = .FALSE.

                    CLOSE(NDSTAB)
                END IF
            END DO
        ELSE
            CALL W3EXPO
        END IF

        CALL TICK21 ( TOUT , DTREQ )
        IF ( IOUT .GE. NOUT ) EXIT
      END DO
!
! ... ITYPE=4 & OTYPES=[2,4] requires adding lines at bottom of
!     bulletin output for compatibility with version 2.22
!
      IF (ITYPE .EQ. 4 .AND. ( OTYPE .EQ. 2 .OR. OTYPE .EQ. 4 ) ) THEN
        DO IJ = 1,NOPTS
          IF (FLREQ(IJ)) THEN
            NDSBUL = NDSTAB + (IJ - 1)
            WRITE(NDSBUL,971) 
            WRITE(NDSBUL,974) BHSDROP, BHSMIN 
!/NCO             NDSCBUL = NDSTAB + (IJ - 1) + NOPTS
!/NCO             WRITE(NDSCBUL,961)
!/NCO             WRITE(NDSCBUL,962)
            close(NDSBUL)
!/NCO            close(NDSCBUL)
            NDSCSV = NDSTAB + (IJ - 1) + 2*NOPTS
            close(NDSCSV)
          ENDIF
        ENDDO
      ENDIF
!
      GOTO 888
!
! Escape locations read errors :
!
  800 CONTINUE
      WRITE (NDSE,1000) IERR
      CALL EXTCDE ( 40 )
!
  801 CONTINUE
      WRITE (NDSE,1001)
      CALL EXTCDE ( 41 )
!
  802 CONTINUE
      WRITE (NDSE,1002) IERR
      CALL EXTCDE ( 42 )
!
  803 CONTINUE
      WRITE (NDSE,1003) IERR
      CALL EXTCDE ( 43 )
!
  804 CONTINUE
      WRITE (NDSE,1004) IERR
      CALL EXTCDE ( 44 )
!
!/O14  805 CONTINUE
!/O14      WRITE (NDSE,1005) IERR
!/O14      CALL EXTCDE ( 45 )
!
  888 CONTINUE
!
      WRITE (NDSO,999)
!
!/NCO/!     CALL W3TAGE('WAVESPEC')
!
! Formats
!
  900 FORMAT (/15X,'    *** WAVEWATCH III Point output post.***    '/ &
               15X,'==============================================='/)
  901 FORMAT ( '  Comment character is ''',A,''''/)
!
  920 FORMAT ( '  Grid name : ',A/)
!
  930 FORMAT ( '  Points in file : '/                                 &
               ' ------------------------------------')
  931 FORMAT ( '      ',A,2F10.2)
  932 FORMAT ( '      ',A,2(F8.1,'E3'))
!
  940 FORMAT (/'  Output time data : '/                               &
               ' --------------------------------------------------'/ &
               '      First time         : ',A)
  941 FORMAT ( '      Interval           : ',A/                       &
               '      Number of requests : ',I6)
  942 FORMAT (/'  Output type ',I2,' :'/                              &
               ' --------------------------------------------------'/ &
               '      ',A/)
  943 FORMAT ( '      Subtype   : ',A)
  944 FORMAT ( '                  ',A)
!/O14   945 FORMAT ( '      ',I5,3X,A,2F10.2,3X,A)
  948 FORMAT ( '      Data for ',A)
  949 FORMAT (/'      End of file reached '/)
!
  950 FORMAT (/'  Requested output for',I5,' points : '/              &
               ' --------------------------------------------------')
  951 FORMAT ( '      ',A,2F10.2)
  953 FORMAT ( '      ',A,2(F8.1,'E3'))
  952 FORMAT (/'  Output times :'/                                    &
               ' --------------------------------------------------')
!/NCO  961 FORMAT ('----------------------------------------',              &
!/NCO              '---------------------------')
!/NCO  962 FORMAT ( 'DD  = Day of Month'/                                   &
!/NCO        'HH  = Hour of Day'/                                           &
!/NCO        'HS  = Total Significant Wave Height (feet)'/                  &
!/NCO        'SS  = Significant Wave Height of separate system (feet)'/     &
!/NCO        'PP  = Peak Period of separate system (whole seconds)'/        &
!/NCO        'DDD = Mean Direction of separate system (degrees,"from")')
  971 FORMAT (' +-------+-----------+-----------------+',             &
             '-----------------+-----------------+----',              &
             '-------------+-----------------+--------',              &
             '---------+')!
  974 FORMAT ( &
       75X,'Hst : Total sigificant wave height.'/                     &
       75X,'n   : Number of fields with Hs > ',f4.2,                  &
                ' in 2-D spectrum.'/                                  &
       75X,'x   : Number of fields with Hs > ',f4.2,                  &
                ' not in table.'/                                     &
       75X,'Hs  : Significant wave height of separate wave field.'/   &
       75X,'Tp  : Peak period of separate wave field.'/               &
       75X,'dir : Mean direction of separate wave field.'/            &
       75X,'*   : Wave generation due to local wind probable.')
 
 1940 FORMAT ( '      ',A,' print plots not requested.')
 1941 FORMAT ( '      ',A,' print plots normalized.')
 1942 FORMAT ( '      Scale factor ',A,' spectrum : ',E10.3)
 1943 FORMAT ( '      File name : ',A,' (',A,')')
 1944 FORMAT ('''',A,'''',1X,3I6,1X,'''',A,'''')
 1945 FORMAT (8E10.3)
 1946 FORMAT (7E11.3)
 1947 FORMAT ( '      File name : ',A)
!
 2940 FORMAT ( '      Table output : ',A/                             &
               '      File name    : ',A)
!
 3940 FORMAT ( '                        ',A)
 3941 FORMAT ( '      File name : ',A)
 3943 FORMAT ( '      File name : ',A,' (',A,')')
 3944 FORMAT ('''',A,'''',1X,3I6,6L2)
 3945 FORMAT (8E10.3)
 3946 FORMAT (7E11.3)
!
  960 FORMAT (//'  Output for ',A/                                    &
               ' --------------------------------------------------')
!
  999 FORMAT (/'  End of program '/                                   &
               ' ========================================='/          &
               '         WAVEWATCH III Point output '/)
!
 1000 FORMAT (/' *** WAVEWATCH III ERROR IN W3OUTP : '/               &
               '     ERROR IN OPENING INPUT FILE'/                    &
               '     IOSTAT =',I5/)
!
 1001 FORMAT (/' *** WAVEWATCH III ERROR IN W3OUTP : '/               &
               '     PREMATURE END OF INPUT FILE'/)
!
 1002 FORMAT (/' *** WAVEWATCH III ERROR IN W3OUTP : '/               &
               '     ERROR IN READING FROM INPUT FILE'/               &
               '     IOSTAT =',I5/)
!
 1003 FORMAT (/' *** WAVEWATCH III ERROR IN W3OUTP : '/               &
               '     ERROR IN OPENING TABLE FILE'/                    &
               '     IOSTAT =',I5/)
!
 1004 FORMAT (/' *** WAVEWATCH III ERROR IN W3OUTP : '/               &
               '     ERROR IN OPENING IDL FILE'/                      &
               '     IOSTAT =',I5/)
!
!/O14 1005 FORMAT (/' *** WAVEWATCH III ERROR IN W3OUTP : '/          &
!/O14               '     ERROR IN OPENING BUOY LOG FILE'/            &
!/O14               '     IOSTAT =',I5/)
!
 1007 FORMAT (/' *** WAVEWATCH III ERROR IN W3OUTP : '/              &
               '     ERROR IN READING FROM INPUT FILE'/               &
               '     LAST POINT INDEX IS NOT -1'/                    &
               '     OR TOO MANY POINT INDEXES DEFINED'/)
!
 1010 FORMAT (/' *** WAVEWATCH III ERROR IN W3OUTP : '/               &
               '     ILLEGAL TYPE, ITYPE =',I4/)
!
 1011 FORMAT (/' *** WAVEWATCH III ERROR IN W3OUTP : '/               &
               '     ILLEGAL TYPE, OTYPE =',I4/)
!
 1012 FORMAT (/' *** WAVEWATCH III ERROR IN W3OUTP : '/               &
               '     MULTIPLE OUTPUT POINTS DEFINED, ITYPE =',I4,/    &
               '     ONLY SINGLE POINT ALLOWED IN THIS VERSION'/)
!/IC1 3960 FORMAT (/' *** WAVEWATCH-III WARNING IN W3OUTP :'/         &
!/IC1               '     Ice source terms !/IC1 skipped'/            &
!/IC1               '     in dissipation term.'/)
!/IC2 3960 FORMAT (/' *** WAVEWATCH-III WARNING IN W3OUTP :'/         &
!/IC2               '     Ice source terms !/IC2 skipped'/            &
!/IC2               '     in dissipation term.'/)
!/IC3 3960 FORMAT (/' *** WAVEWATCH-III WARNING IN W3OUTP :'/         &
!/IC3               '     Ice source terms !/IC3 skipped'/            &
!/IC3               '     in dissipation term.'/)
!/IC4 3960 FORMAT (/' *** WAVEWATCH-III WARNING IN W3OUTP :'/         &
!/IC4               '     Ice source terms !/IC4 skipped'/            &
!/IC4               '     in dissipation term.'/)
!/IC5 3960 FORMAT (/' *** WAVEWATCH-III WARNING IN W3OUTP :'/         &
!/IC5               '     Ice source terms !/IC5 skipped'/            &
!/IC5               '     in dissipation term.'/)
!
!/
!/ Internal subroutine W3EXPO ---------------------------------------- /
!/
      CONTAINS
!/ ------------------------------------------------------------------- /
      SUBROUTINE W3EXPO
!/
!/                  +-----------------------------------+
!/                  | WAVEWATCH III           NOAA/NCEP |
!/                  |           H. L. Tolman            |
!/                  |            J.H. Alves             |
!/                  |            F. Ardhuin             |
!/                  |             A. Chawla             |
!/                  |                        FORTRAN 90 |
!/                  | Last update :         06-Feb-2014 |
!/                  +-----------------------------------+
!/
!/    08-Jun-1999 : Final FORTRAN 77                    ( version 1.18 )
!/    21-Jan-2000 : Upgrade to FORTRAN 90               ( version 2.00 )
!/                  Massive changes to logistics
!/    09-Jan-2001 : U* bug fix in tabular output        ( version 2.05 )
!/    25-Jan-2001 : Flat grid version                   ( version 2.06 )
!/    02-Feb-2001 : Xnl version 3.0                     ( version 2.07 )
!/    11-Jun-2001 : Clean up                            ( version 2.11 )
!/    11-Oct-2001 : Clean up, X*, Y* in tables          ( version 2.14 )
!/    24-Dec-2004 : Multiple grid version.              ( version 3.06 )
!/    17-Apr-2006 : Filter for directional spread.      ( version 3.09 )
!/    23-Jun-2006 : Linear input added.                 ( version 3.09 )
!/    03-Jul-2006 : Separate flux modules.              ( version 3.09 )
!/    28-Oct-2006 : Add partitioning option.            ( version 3.10 )
!/    24-Mar-2007 : Add pars for entire spectrum.       ( version 3.11 )
!/    25-Apr-2007 : Battjes-Janssen Sdb added.          ( version 3.11 )
!/                  (J. H. Alves)
!/    09-Oct-2007 : WAM 4+ Sin and Sds added.           ( version 3.13 )
!/                  (F. Ardhuin)
!/    09-Oct-2007 : Experimental Sbs (BS1) added.       ( version 3.13 )
!/                  (F. Ardhuin)
!/    09-Apr-2008 : Adding an additional output for     ( version 3.12 )
!/                  WMO standard (A. Chawla)
!/    29-Apr-2008 : Adjust format partition output.     ( version 3.14 )
!/    01-Jul-2011 : Adding BT4                          ( version 4.01 )
!/    16-Jul-2012 : Move GMD (SNL3) and nonlinear filter (SNLS)
!/                  from 3.15 (HLT).                    ( version 4.08 )
!/    26-Dec-2012 : Modified obsolete declarations.     ( version 4.11 )
!/    06-Feb-2014 : Fix header format in part. files.   ( version 4.18 )
!/
!  1. Purpose :
!
!     Perform actual point output.
!
!  3. Parameters :
!
!  4. Subroutines used :
!
!      Name      Type  Module   Description
!     ----------------------------------------------------------------
!      W3SPRn    Subr. W3SRCnMD Mean wave parameters for use in
!                               source terms.
!      W3FLXn    Subr. W3FLXnMD Flux/stress computation.
!      W3SLNn    Subr. W3SLNnMD Linear input.
!      W3SINn    Subr. W3SRCnMD Input source term.
!      W3SDSn    Subr. W3SRCnMD Whitecapping source term
!      W3SNLn    Subr. W3SNLnMD Nonlinear interactions.
!      W3SBTn    Subr. W3SBTnMD Bottom friction source term.
!      W3SDBn    Subr. W3SBTnMD Depth induced breaking source term.
!      W3STRn    Subr. W3STRnMD Triad interaction source term.
!      W3SBSn    Subr. W3SBSnMD Bottom scattering source term.
!      W3SXXn    Subr. W3SXXnMD Unclassified source term.
!      W3PART    Sunr. W3PARTMD Spectral partitioning routine.
!      STRACE    Subr. W3SERVMD Subroutine tracing.
!      STME21    Subr. W3TIMEMD Convert time to string.
!      PRT1DS    Subr. W3ARRYMD Print plot of 1-D spectrum.
!      PRT1DM    Subr.   Id.    Print plot of several 1-D spectra.
!      PRT2DS    Subr.   Id.    Print plot of 2-D spectrum.
!      WAVNU1    Subr. W3DISPMD Solve dispersion relation.
!     ----------------------------------------------------------------
!
!  5. Called by :
!
!     Main program in which it is contained,
!
!  6. Error messages :
!
!     None.
!
!  7. Remarks :
!
!     - Spectra are relative frequency energy spectra.
!     - Note that arrays CX and CY of the main program now contain
!       the absolute current speed and direction respectively.
!
!     - BT8&9 issues : 
!
!       Q: What is the problem?
!       A: Point output of Sbot with BT8 or BT9 is not presently 
!          supported.
!
!       Q: What can a user do now?
!       A: When using BT8 or BT9 with ITYPE=3 , the
!          user should set the 5th T/F value in ww3_outp.inp for
!          ITYPE=3 to "F" like so :
!          2  1. 1. 51   T T T T F T  0  F
!        $               ^ ^ ^ ^ ^ ^ Sum of selected sources
!        $               | | | | ^ Wave-bottom interactions
!        $               | | | ^ Dissipation
!        $               | | ^ Nonlinear interactions
!        $               | ^ Wind-wave interactions
!        $               ^ Spectrum 
!          If the user really need this source function, he/she
!          needs to add test output to the mud subroutine 
!          directly
!
!       Q: Why doesn't this functionality exist?
!       A: The Sbot source function in ww3_outp was originally written 
!          with the case of BT1 in mind. BT1 uses a uniform friction 
!          factor, so it does not need any special variable for the 
!          local friction factor. BT8 and BT9 allow non-uniform mud 
!          variables (thickness, density, viscosity) and the mud 
!          subroutines are written with ww3_shel in mind, where the 
!          source function is calculated on the computational grid 
!          point IX IY. 

!       Q: How can we add this functionality?
!       A: To fix it, we would need to :
!          1) interpolate the mud variables from the computational
!             grid point IX IY to the output points (this is already
!             done now for wind, for example) (the same should probably
!             be done for the ice properties also) This would be done
!             in w3iopomd.ftn, analogous to what is done now for the
!             wind variable WAO.
!          2) manage the arrays for the new variables (mud and ice
!             properties on the output points) This would be done in
!             w3odatmd.ftn, again analogous to what is done now for the
!             wind variable WAO.
!          3) change the mud routines so that they take the local mud
!             parameters through the subroutine arguments rather than
!             taking IX IY as subroutine arguments. This would allow
!             flexibility to call the mud routine from ww3_shel or
!             ww3_outp (instead of just ww3_shel as is the case now).
!
!/---------------------------------------------------------------------/
!
!  8. Structure :
!
!     See source code.
!
!  9. Switches :
!
!       !/S      Enable subroutine tracing.
!       !/T      Enable test output.
!
!       !/FLXx   Flux/stress computation.
!       !/LNx    Linear input package
!       !/STx    Source term package
!       !/NLx    Nonlinear interaction package
!       !/BTx    Bottom friction package
!       !/ICx    S_ice source term package
!       !/DBx    Depth-induced breaking package
!       !/TRx    Triad interaction package
!       !/BSx    Bottom scattering package
!       !/XXx    Arbitrary adittional source term package
!
!       !/STAB2  Stability correction for !/ST2
!
! 10. Source code :
!
!/ ------------------------------------------------------------------- /
!/FLX1      USE W3FLX1MD
!/FLX2      USE W3FLX2MD
!/FLX3      USE W3FLX3MD
!/FLX4      USE W3FLX4MD
!/LN1      USE W3SLN1MD
!/LNX      USE W3SLNXMD
!/ST1      USE W3SRC1MD
!/ST2      USE W3SRC2MD
!/ST3      USE W3SRC3MD
!/ST4      USE W3SRC4MD, ONLY : W3SPR4, W3SIN4, W3SDS4
!/ST6      USE W3SRC6MD
!/ST6      USE W3SWLDMD, ONLY : W3SWL6
!/ST6      USE W3GDATMD, ONLY : SWL6S6
!/STX      USE W3SRCXMD
!/NL1      USE W3SNL1MD
!/NL2      USE W3SNL2MD
!/NL3      USE W3SNL3MD
!/NL4      USE W3SNL4MD
!/NLX      USE W3SNLXMD
!/NLS      USE W3SNLSMD
!/BT1      USE W3SBT1MD
!/BT2      USE W3SBT2MD
!/BT4      USE W3SBT4MD
!/BT8      USE W3SBT8MD
!/BT9      USE W3SBT9MD
!/BTX      USE W3SBTXMD
!/DB1      USE W3SDB1MD
!/DBX      USE W3SDBXMD
!/TRX      USE W3STRXMD
!/BS1      USE W3SBS1MD
!/BSX      USE W3SBSXMD
!/IS2      USE W3SIS2MD
!/IS2      USE W3GDATMD, ONLY: IICEDISP
!/XXX      USE W3SXXXMD
      USE W3PARTMD, ONLY: W3PART
      USE W3DISPMD, ONLY: WAVNU1, LIU_FORWARD_DISPERSION
!/
      USE W3ARRYMD, ONLY: PRT1DS, PRT2DS, PRT1DM
      USE W3DISPMD, ONLY: NAR1D, DFAC, N1MAX, ECG1, EWN1, DSIE
      USE W3BULLMD, ONLY: W3BULL
!/
!/ ------------------------------------------------------------------- /
!/ Local parameters
!/
      INTEGER                 :: J, I1, I2, ISP, IKM, IKL, IKH, ITH,  &
                                 IK, IH, IM, IS, IYR, IMTH, IDY, ITT, &
                                 I, NPART, IP, IX, IY, ISEA
      INTEGER, SAVE           :: IPASS  = 0
!/S      INTEGER, SAVE           :: IENT   = 0
      REAL                    :: DEPTH, SQRTH, CDIR, SIX, R1, R2,     &
                                 UDIR, UDIRR, UABS, XL, XH, XL2, XH2, &
                                 ET, EWN, ETR, ETX, ETY, EBND, EBX,   &
                                 EBY, HSIG, WLEN, TMEAN, THMEAN,      &
                                 THSPRD, EMAX, EL, EH, DENOM, FP, THP,&
                                 SPP, CD, USTAR, FACTOR, UNORM, ESTAR,&
                                 FPSTAR, FACF, FACE, FACS, HMAT, WNA, &
                                 XYZ, AGE1, AFR, AGE2, FACT, XSTAR,   &
                                 YSTAR, FHIGH, ZWND, Z0, USTD, EMEAN, &
                                 FMEAN, WNMEAN, UDIRCA, X, Y, CHARN,  &
                                 M2KM, ICEF, ICEDMAX, ICETHICK,       &
                                 ICECON
!/IS2      REAL                    :: WN_R(NK),CG_ICE(NK), ALPHA_LIU(NK)
!/ST1      REAL                    :: AMAX, FH1, FH2
!/ST2      REAL                    :: AMAX, ALPHA(NK), FPI
!/ST3      REAL                    :: AMAX, FMEANS, FMEANWS, TAUWX, TAUWY, &
!/ST3                                 TAUWNX, TAUWNY
!/ST4      REAL                    :: AMAX, FMEANS, FMEANWS, TAUWX, TAUWY, &
!/ST4                                 TAUWNX, TAUWNY, FMEAN1, WHITECAP(1:4)
!/ST6      REAL                    :: AMAX, TAUWX, TAUWY, TAUWNX, TAUWNY
!/BS1      REAL                    :: TAUSCX, TAUSCY
!/BT4      REAL                    :: D50, PSIC, BEDFORM(3), TAUBBL(2)
           REAL                    :: ICE
!/STAB2      REAL                    :: STAB0, STAB,  COR1, COR2, ASFAC,     &
!/STAB2                                 THARG1, THARG2
      REAL, SAVE              :: HSMIN  = 0.05
      REAL                    :: WN(NK), CG(NK), R(NK)
      REAL                    :: E(NK,NTH), E1(NK), APM(NK),           &
                                 THBND(NK), SPBND(NK), A(NTH,NK),      &
                                 WN2(NTH,NK)
      REAL                    :: DIA(NTH,NK), SWN(NK,NTH), SNL(NK,NTH),&
                                 SDS(NK,NTH), SBT(NK,NTH), SIS(NK,NTH),&
                                 STT(NK,NTH), DIA2(NTH,NK)
      REAL                    :: XLN(NTH,NK), XIN(NTH,NK), XNL(NTH,NK),&
                                 XTR(NTH,NK), XDS(NTH,NK), XDB(NTH,NK),&
                                 XBT(NTH,NK), XBS(NTH,NK), XXX(NTH,NK),&
                                 XIS(NTH,NK), XWL(NTH,NK)
      REAL                    :: SIN1(NK), SNL1(NK), SDS1(NK),         &
                                 SBT1(NK), STT1(NK), SIS1(NK),         &
                                 E1ALL(NK,6)
      LOGICAL                 :: LBREAK
!/ST3      LOGICAL                 :: LLWS(NSPEC)
!/ST4      LOGICAL                 :: LLWS(NSPEC)
!/ST4      REAL                    :: LAMBDA(NSPEC)
      CHARACTER               :: DTME21*23
      CHARACTER(LEN=4)         VAR1(6)
      CHARACTER(LEN=1)         IDLAT, IDLON
      CHARACTER(LEN=100)       BT8MSG
!
      DATA VAR1   / 'Sin ' , 'Snl ', 'Sds ' , 'Sbt ' , 'Sice', 'Stot' /
!/
!/ ------------------------------------------------------------------- /
!/
! 1. Initialisations
!
!/S      CALL STRACE (IENT, 'W3EXPO')
!
      IF ( FLAGLL ) THEN
          M2KM   = 1.
        ELSE
          M2KM   = 1.E-3
        END IF
!
      XL     = 1./XFR - 1.
      XH     =  XFR - 1.
      XL2    = XL**2
      XH2    = XH**2
      IPASS  = IPASS + 1
!
      IF ( ITYPE .EQ. 3 ) THEN
          XLN = 0.
          XIN = 0.
          XNL = 0.
          XTR = 0.
          XDS = 0.
          XDB = 0.
          XBT = 0.
          XBS = 0.
          XWL = 0.
          XXX = 0.
          XIS = 0.
        END IF
!
!/T      WRITE (NDST,9000) (FLREQ(J),J=1,NOPTS)
!/T      WRITE (NDST,9001) ITYPE, OTYPE, NREQ, SCALE1, SCALE2, FLSRCE
!
!     Output of time
!
      IF (  ( ITYPE.EQ.1 .AND. OTYPE.EQ.3 ) .OR.                      &
            ( ITYPE.EQ.3 .AND. OTYPE.EQ.4 ) ) THEN
          IF ( FLFORM ) THEN
              WRITE (NDSTAB) TIME
            ELSE
              WRITE (NDSTAB,900) TIME
            END IF
        END IF
!
      IF (ITYPE.EQ.2) THEN
          IF ( NREQ.EQ.1 .AND. IPASS.EQ.1 ) THEN
              IF ( OTYPE.EQ.1 ) THEN
                  WRITE (NDSTAB,1901)
                ELSE IF ( OTYPE.EQ.2 ) THEN
                  WRITE (NDSTAB,1902)
                ELSE IF ( OTYPE.EQ.3 ) THEN
                  WRITE (NDSTAB,1903)
                ELSE IF ( OTYPE.EQ.4 ) THEN
                  WRITE (NDSTAB,1904)
                ELSE IF ( OTYPE.EQ.5 ) THEN
                  WRITE (NDSTAB,1905)
                ELSE IF ( OTYPE.EQ.6 ) THEN
                  WRITE (NDSTAB,1906)
                END IF
            END IF
          IF ( NREQ.NE.1 ) THEN
              CALL STME21 ( TIME , DTME21 )
              IF ( IPASS .NE. 1 ) WRITE (NDSTAB,1910)
              IF ( OTYPE.EQ.1 ) THEN
                  IF ( FLAGLL ) THEN
                      WRITE (NDSTAB,1911) DTME21
                    ELSE
                      WRITE (NDSTAB,1711) DTME21
                    END IF
                ELSE IF ( OTYPE.EQ.2 ) THEN
                  IF ( FLAGLL ) THEN
                      WRITE (NDSTAB,1912) DTME21
                    ELSE
                      WRITE (NDSTAB,1712) DTME21
                    END IF
                ELSE IF ( OTYPE.EQ.3 ) THEN
                  WRITE (NDSTAB,1913) DTME21
                ELSE IF ( OTYPE.EQ.4 ) THEN
                  WRITE (NDSTAB,1914) DTME21
                ELSE IF ( OTYPE.EQ.5 ) THEN
                  IF ( FLAGLL ) THEN
                      WRITE (NDSTAB,1915) DTME21
                    ELSE
                      WRITE (NDSTAB,1715) DTME21
                    END IF
                ELSE IF ( OTYPE.EQ.6 ) THEN
                  IF ( FLAGLL ) THEN
                      WRITE (NDSTAB,1916) DTME21
                    ELSE
                      WRITE (NDSTAB,1716) DTME21
                    END IF
                END IF
            END IF
        END IF
!
      IF (ITYPE.EQ.3) THEN
            IF ( OTYPE .EQ. 4 ) THEN
              CALL STME21 ( TIME , DTME21 )
              WRITE (NDSO,905) DTME21
            END IF
        END IF
!
!--- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!     Loop over output points.
 
      IF (PROCESS_POINT_ONLY) THEN
          J_START = ACTIVE_POINT
          J_END = ACTIVE_POINT
      ELSE
          J_START = 1
          J_END = NOPTS
      END IF

      DO J=J_START, J_END
        IF ( FLREQ(J) ) THEN
!
!/T            WRITE (NDST,9002) PTNME(J)
!
! 2. Calculate grid parameters using and inlined version of WAVNU1.
!
            DEPTH    = MAX ( DMIN, DPO(J) )
            SQRTH    = SQRT ( DEPTH )
            UDIR     = MOD ( 270. - WDO(J)*RADE , 360. )
            UDIRCA   = WDO(J)*RADE
            UDIRR    = WDO(J)
            UABS     = MAX ( 0.001 , WAO(J) )
            CDIR     = MOD ( 270. - CDO(J)*RADE , 360. )
            ICEDMAX  = MAX ( 0., ICEFO(J))
            ICEF     = ICEDMAX
            ICETHICK = MAX (0., ICEHO(J))
            ICECON   = MAX (0., ICEO(J))
!
!/STAB2            STAB0  = ZWIND * GRAV / 273.
!/STAB2            STAB   = STAB0 * ASO(J) / MAX(5.,WAO(J))**2
!/STAB2            STAB   = MAX ( -1. , MIN ( 1. , STAB ) )
!/STAB2            THARG1 = MAX ( 0. , FFNG*(STAB-OFSTAB))
!/STAB2            THARG2 = MAX ( 0. , FFPS*(STAB-OFSTAB))
!/STAB2            COR1   = CCNG * TANH(THARG1)
!/STAB2            COR2   = CCPS * TANH(THARG2)
!/STAB2            ASFAC  = SQRT ( (1.+COR1+COR2)/SHSTAB )
!
!/T            WRITE (NDST,9010) DEPTH
            DO IK=1, NK
              SIX    = SIG(IK) * SQRTH
              I1     = INT(SIX/DSIE)
              IF (I1.LE.N1MAX) THEN
                  I2 = I1 + 1
                  R1 = SIX/DSIE - REAL(I1)
                  R2 = 1. - R1
                  WN(IK) = ( R2*EWN1(I1) + R1*EWN1(I2) ) / DEPTH
                  CG(IK) = ( R2*ECG1(I1) + R1*ECG1(I2) ) * SQRTH
                ELSE
                  WN(IK) = SIG(IK)*SIG(IK)/GRAV
                  CG(IK) = 0.5 * GRAV / SIG(IK)
                END IF
!/T              WRITE (NDST,9011) IK, TPI/SIG(IK), WN(IK), CG(IK)
!
              END DO

!
! Computes 2nd order spectrum 
!
!/IG1      IF (IGPARS(2).EQ.1) THEN 
!/IG1        IF(IGPARS(1).EQ.1) THEN 
!/IG1          CALL W3ADDIG(SPCO(:,J),DPO(J),WN,CG,0)
!/IG1        ELSE
!/IG1          CALL W3ADD2NDORDER(SPCO(:,J),DPO(J),WN,CG,0)
!/IG1          END IF
!/IG1        END IF
!
! 3.  Prepare spectra etc.
! 3.a Mean wave parameters.
!
            ET     = 0.
            EWN    = 0.
            ETR    = 0.
            ETX    = 0.
            ETY    = 0.
            DO IK=1, NK
              EBND   = 0.
              EBX    = 0.
              EBY    = 0.
              DO ITH=1, NTH
                ISP    = ITH + (IK-1)*NTH
                E(IK,ITH) = SPCO(ISP,J)
                EBND   = EBND + SPCO(ISP,J)
                EBX    = EBX  + SPCO(ISP,J)*ECOS(ITH)
                EBY    = EBY  + SPCO(ISP,J)*ESIN(ITH)
                END DO
              E1(IK) = EBND * DTH
              APM(IK)= E1(IK) / ( TPI * GRAV**2 / SIG(IK)**5  )
              IF ( E1(IK) .GT. 1.E-5) THEN
                  THBND(IK) = MOD(630.- RADE*ATAN2(EBY,EBX),360.)
                  SPBND(IK) = RADE * SQRT ( MAX ( 0. , 2.*( 1. -      &
                    SQRT( MAX(0.,(EBX**2+EBY**2)/EBND**2) ) ) ) )
                ELSE
                  THBND(IK) = -999.9
                  SPBND(IK) = -999.9
                END IF
              EBND   = E1(IK) * DSII(IK) * TPIINV
              ET     = ET  + EBND
              EWN    = EWN + EBND / WN(IK)
              ETR    = ETR + EBND / SIG(IK)
              ETX    = ETX + EBX * DSII(IK)
              ETY    = ETY + EBY * DSII(IK)
              END DO
!
! tail factors for radian action etc ...!
!
            EBND   = E1(NK) * TPIINV / ( SIG(NK) * DTH )
            ET     = ET  + FTE *EBND
            EWN    = EWN + FTWL*EBND
            ETR    = ETR + FTTR*EBND
            ETX    = DTH*ETX*TPIINV + FTE*EBX*TPIINV/SIG(NK)
            ETY    = DTH*ETY*TPIINV + FTE*EBY*TPIINV/SIG(NK)
!
            HSIG   = 4. * SQRT ( MAX(0.,ET) )
            IF ( HSIG .GT. HSMIN ) THEN
                WLEN   = EWN / ET * TPI
                TMEAN  = ETR / ET * TPI
                THMEAN = MOD ( 630. - RADE*ATAN2(ETY,ETX) , 360. )
                THSPRD = RADE * SQRT ( MAX ( 0. , 2.*( 1. - SQRT(     &
                           MAX(0.,(ETX**2+ETY**2)/ET**2) ) ) ) )
                IF ( THSPRD .LT. 0.01*RADE*DTH ) THSPRD = 0.
              ELSE
                WLEN   = 0.
                TMEAN  = 0.
                THMEAN = 0.
                THSPRD = 0.
                DO IK=1, NK
                  E1(IK) = 0.
                  DO ITH=1, NTH
                    E(IK,ITH) = 0.
                  END DO
                END DO
              END IF
!
! peak frequency
!
            EMAX   = E1(NK)
            IKM    = NK
!
            DO IK=NK-1, 1, -1
              IF ( E1(IK) .GT. EMAX ) THEN
                  EMAX   = E1(IK)
                  IKM    = IK
                END IF
              END DO
!
            IKL    = MAX (  1 , IKM-1 )
            IKH    = MIN ( NK , IKM+1 )
            EL     = E1(IKL) - E1(IKM)
            EH     = E1(IKH) - E1(IKM)
            DENOM  = XL*EH - XH*EL
!
            IF ( HSIG .GE. HSMIN ) THEN
                FP     = SIG(IKM) * ( 1. + 0.5 * ( XL2*EH - XH2*EL )  &
                            / SIGN ( MAX(ABS(DENOM),1.E-15) , DENOM ) )
                THP    = THBND(IKM)
                SPP    = SPBND(IKM)
                IF ( SPP .LT. 0.01*RADE*DTH ) SPP = 0.
              ELSE
                FP     = 0.
                THP    = 0.
                SPP    = 0.
              END IF
!
! spectral partitioning
!
            IF ( ITYPE.EQ.4 ) CALL W3PART                              &
               ( E, UABS, UDIRCA, DEPTH, WN, NPART, XPART, DIMXP )
!
! nondimensional parameters
!
            IF ( ( ITYPE.EQ.2 .AND. (OTYPE.EQ.3.OR.OTYPE.EQ.4) ) .OR. &
                 ( ITYPE.EQ.1 .AND. (OTYPE.EQ.2) ) ) THEN
!
                DO IK=1, NK
                  FACTOR = TPIINV * CG(IK) / SIG(IK)
                  DO ITH=1, NTH
                    ISP    = ITH + (IK-1)*NTH
                    A(ITH,IK)   = FACTOR * SPCO(ISP,J)
                    WN2(ITH,IK) = WN(IK)
                    END DO
                  END DO
!
!/STAB2                UABS   = UABS / ASFAC
!
!/ST0                ZWND   = 10.
!/ST1                ZWND   = 10.
!/ST2                ZWND   = ZWIND
!/ST3                ZWND   = ZZWND
!/ST3                TAUWX  = 0.
!/ST3                TAUWY  = 0.
!/ST3                LLWS(:)  = .TRUE.
!/ST4                LLWS(:)  = .TRUE.
!/ST4                ZWND   = ZZWND
!/ST4                TAUWX  = 0.
!/ST4                TAUWY  = 0.
!/ST6                ZWND   = 10.
                USTAR  = 1.
!
!/ST1                CALL W3SPR1 (A, CG, WN, EMEAN, FMEAN, WNMEAN, AMAX)
!/ST1                FP     = 0.85 * FMEAN
!/ST2                CALL W3SPR2 (A, CG, WN, DEPTH, FP , UABS, USTAR, &
!/ST2                             EMEAN, FMEAN, WNMEAN, AMAX, ALPHA, FP )
!/ST3                CALL W3SPR3 (A, CG, WN, EMEAN, FMEAN, FMEANS,       &
!/ST3                             WNMEAN, AMAX, UABS, UDIRR, USTAR, USTD,&
!/ST3                             TAUWX, TAUWY, CD, Z0, CHARN, LLWS, FMEANWS )
!/ST4                CALL W3SPR4 (A, CG, WN, EMEAN, FMEAN,  FMEAN1,        &
!/ST4                             WNMEAN, AMAX, UABS, UDIRR, USTAR, USTD,&
!/ST4                             TAUWX, TAUWY, CD, Z0, CHARN, LLWS, FMEANWS )
!/ST6                CALL W3SPR6 (A, CG, WN, EMEAN, FMEAN, WNMEAN, AMAX, FP)
!
!/FLX1                CALL W3FLX1 ( ZWND, UABS, UDIRR,                   &
!/FLX1                              USTAR, USTD, Z0, CD )
!/FLX2                CALL W3FLX2 ( ZWND, DEPTH, FP, UABS, UDIRR,        &
!/FLX2                                          USTAR, USTD, Z0, CD )
!/FLX3                CALL W3FLX3 ( ZWND, DEPTH, FP, UABS, UDIRR,        &
!/FLX3                                          USTAR, USTD, Z0, CD )
!/FLX4                CALL W3FLX4 ( ZWND, UABS, UDIRR, USTAR, USTD, Z0, CD )
!/FLXX                CALL W3FLXX
!
                DO ITT=1, 3
!/ST2                  CALL W3SIN2 (A, CG, WN2, UABS, UDIRR, CD, Z0,    &
!/ST2                                                 FPI, XIN, DIA )
!/ST2                  CALL W3SPR2 (A, CG, WN, DEPTH, FPI, UABS, USTAR, &
!/ST2                               EMEAN, FMEAN, WNMEAN, AMAX, ALPHA, FP )
!/ST3                  IX=1
!/ST3                  IY=1
!/ST3                  CALL W3SIN3 ( A, CG, WN2, UABS, USTAR, DAIR/DWAT,&
!/ST3                               ASO(J), UDIRR, Z0, CD, TAUWX, TAUWY,&
!/ST3                               TAUWNX, TAUWNY, ICE, XIN, DIA, LLWS, IX, IY )
!/ST3                  CALL W3SPR3 (A, CG, WN, EMEAN, FMEAN, FMEANS,       &
!/ST3                             WNMEAN, AMAX, UABS, UDIRR, USTAR, USTD,&
!/ST3                             TAUWX, TAUWY, CD, Z0, CHARN, LLWS, FMEANWS )
!/ST4                  IX=1
!/ST4                  IY=1
!/ST4                  CALL W3SIN4 ( A, CG, WN2, UABS, USTAR, DAIR/DWAT,&
!/ST4                               ASO(J), UDIRR, Z0, CD, TAUWX, TAUWY,&
!/ST4                               TAUWNX, TAUWNY, XIN, DIA, LLWS, IX, IY, LAMBDA )
!/ST4                  CALL W3SPR4 (A, CG, WN, EMEAN, FMEAN,  FMEAN1,      &
!/ST4                             WNMEAN, AMAX, UABS, UDIRR, USTAR, USTD,&
!/ST4                             TAUWX, TAUWY, CD, Z0, CHARN, LLWS, FMEANWS )
!/FLX2                  CALL W3FLX2 ( ZWND, DEPTH, FP, UABS, UDIRR,      &
!/FLX2                                            USTAR, USTD, Z0, CD )
!/FLX3                  CALL W3FLX3 ( ZWND, DEPTH, FP, UABS, UDIRR,      &
!/FLX3                                            USTAR, USTD, Z0, CD )
                  END DO
!
! Add alternative flux calculations here as part of !/ST2 option ....
! Also add before actual source term calculation !!!
!
!/STAB2                UABS   = UABS * ASFAC
!
              IF ( WAO(J) .LT. 0.01 ) THEN
                  UNORM  = 0.
                  ESTAR  = 0.
                  FPSTAR = 0.
                ELSE
                  IF ( OTYPE.EQ.3 ) THEN
                      UNORM  = USTAR
                    ELSE
                      UNORM  = WAO(J)
                    END IF
                  ESTAR  = ET * GRAV**2 / UNORM**4
                  FPSTAR = FP * TPIINV * UNORM / GRAV
                  XSTAR  = PTLOC(1,J) * GRAV / UNORM**2
                  YSTAR  = PTLOC(2,J) * GRAV / UNORM**2
                  IF ( FLAGLL ) THEN
                      XSTAR  = XSTAR * DERA * RADIUS &
                             * COS(PTLOC(2,J)*DERA)
                      YSTAR  = YSTAR * DERA * RADIUS
                    END IF
                END IF
!
              END IF
!
! 3.4 source terms
!
            IF ( ITYPE.EQ.3 ) THEN
!
                DO IK=1, NK
                  FACTOR = TPIINV * CG(IK) / SIG(IK)
                  DO ITH=1, NTH
                    A(ITH,IK)   = FACTOR * SPCO(ITH+(IK-1)*NTH,J)
                    WN2(ITH,IK) = WN(IK)
                    END DO
                  END DO
!
!/STAB2                UABS   = UABS / ASFAC
!
!/ST0                ZWND   = 10.
!/ST1                ZWND   = 10.
!/ST2                ZWND   = ZWIND
!/ST3                ZWND   = ZZWND
!/ST0                USTAR  = 1.
!/ST1                USTAR  = 1.
!/ST2                USTAR  = 1.
!/ST3                USTAR  = 0.
!/ST3                USTD   = 0.
!/ST3                TAUWX  = 0.
!/ST3                TAUWY  = 0.
!/ST4                ZWND   = ZZWND
!/ST4                USTAR  = 0.
!/ST4                USTD   = 0.
!/ST4                TAUWX  = 0.
!/ST4                TAUWY  = 0.
!/ST6                ZWND   = 10.
!
!/ST0                FHIGH  = SIG(NK)
!/ST1                CALL W3SPR1 (A, CG, WN, EMEAN, FMEAN, WNMEAN, AMAX)
!/ST1                FP     = 0.85 * FMEAN
!/ST1                FH1    = FXFM * FMEAN
!/ST1                FH2    = FXPM / USTAR
!/ST1                FHIGH  = MAX ( FH1 , FH2 )
!/ST2                CALL W3SPR2 (A, CG, WN, DEPTH, FP , UABS, USTAR, &
!/ST2                             EMEAN, FMEAN, WNMEAN, AMAX, ALPHA, FP )
!/ST3                CALL W3SPR3 (A, CG, WN, EMEAN, FMEAN, FMEANS,       &
!/ST3                             WNMEAN, AMAX, UABS, UDIRR, USTAR, USTD,&
!/ST3                             TAUWX, TAUWY, CD, Z0, CHARN, LLWS, FMEANWS )
!/ST4                CALL W3SPR4 (A, CG, WN, EMEAN, FMEAN,  FMEAN1,        &
!/ST4                             WNMEAN, AMAX, UABS, UDIRR, USTAR, USTD,&
!/ST4                             TAUWX, TAUWY, CD, Z0, CHARN, LLWS, FMEANWS )
!/ST6                CALL W3SPR6 (A, CG, WN, EMEAN, FMEAN, WNMEAN, AMAX, FP)
!/ST6                FHIGH  = SIG(NK)
!/STX                CALL W3SPRX
!
!/FLX1                CALL W3FLX1 ( ZWND, UABS, UDIRR,                   &
!/FLX1                              USTAR, USTD, Z0, CD )
!/FLX2                CALL W3FLX2 ( ZWND, DEPTH, FP, UABS, UDIRR,        &
!/FLX2                                          USTAR, USTD, Z0, CD )
!/FLX3                CALL W3FLX3 ( ZWND, DEPTH, FP, UABS, UDIRR,        &
!/FLX3                                          USTAR, USTD, Z0, CD )
!/FLX4                CALL W3FLX4 ( ZWND, UABS, UDIRR, USTAR, USTD, Z0, CD )
!/FLXX                CALL W3FLXX
!
                DO ITT=1, 3
!/ST2                  CALL W3SIN2 (A, CG, WN2, UABS, UDIRR, CD, Z0,    &
!/ST2                                                 FPI, XIN, DIA )
!/ST2                  CALL W3SPR2 (A, CG, WN, DEPTH, FPI, UABS, USTAR, &
!/ST2                               EMEAN, FMEAN, WNMEAN, AMAX, ALPHA, FP )
!/ST3                CALL W3SPR3 (A, CG, WN, EMEAN, FMEAN, FMEANS,       &
!/ST3                             WNMEAN, AMAX, UABS, UDIRR, USTAR, USTD,&
!/ST3                             TAUWX, TAUWY, CD, Z0, CHARN, LLWS, FMEANWS )
!/ST3                  CALL W3SIN3 ( A, CG, WN2, UABS, USTAR, DAIR/DWAT,&
!/ST3                                ASO(J), UDIRR, Z0, CD,TAUWX, TAUWY, &
!/ST3                                TAUWNX, TAUWNY, ICE, XIN, DIA, LLWS, IX, IY )
!/ST4                CALL W3SPR4 (A, CG, WN, EMEAN, FMEAN, FMEAN1,        &
!/ST4                             WNMEAN, AMAX, UABS, UDIRR, USTAR, USTD,&
!/ST4                             TAUWX, TAUWY, CD, Z0, CHARN, LLWS, FMEANWS )
!/ST4                  CALL W3SIN4 ( A, CG, WN2, UABS, USTAR, DAIR/DWAT,&
!/ST4                                ASO(J), UDIRR, Z0, CD,TAUWX, TAUWY, &
!/ST4                                TAUWNX, TAUWNY, XIN, DIA, LLWS, IX, IY, LAMBDA )
!/FLX2                  CALL W3FLX2 ( ZWND, DEPTH, FP, UABS, UDIRR,      &
!/FLX2                                            USTAR, USTD, Z0, CD )
!/FLX3                  CALL W3FLX3 ( ZWND, DEPTH, FP, UABS, UDIRR,      &
!/FLX3                                            USTAR, USTD, Z0, CD )
                  END DO
!
!/ST2                FHIGH  = XFC * FPI
!
                IF ( FLSRCE(2) ) THEN
!/LN1                    CALL W3SLN1 (WN, FHIGH, USTAR, UDIRR, XLN )
!/LNX                    CALL W3SLNX
!
!/ST1                    CALL W3SIN1 (A, WN2, USTAR, UDIRR,  XIN, DIA )
!/ST2                    CALL W3SIN2 (A, CG, WN2, UABS, UDIRR, CD, Z0,&
!/ST2                                                   FPI, XIN, DIA )
!/ST3                    CALL W3SIN3 ( A, CG, WN2, UABS, USTAR,       &
!/ST3                                  DAIR/DWAT, ASO(J), UDIRR,      &
!/ST3                            Z0, CD, TAUWX, TAUWY,TAUWNX, TAUWNY, &
!/ST3                            ICE, XIN, DIA, LLWS, IX, IY )
!/ST4                    CALL W3SIN4 ( A, CG, WN2, UABS, USTAR,       &
!/ST4                                  DAIR/DWAT, ASO(J), UDIRR,      &
!/ST4                            Z0, CD, TAUWX, TAUWY,TAUWNX, TAUWNY, &
!/ST4                            XIN, DIA, LLWS, IX, IY, LAMBDA )
!/ST6                    CALL W3SIN6 (A, CG, WN2, UABS, USTAR, UDIRR, CD,   &
!/ST6                          TAUWX, TAUWY, TAUWNX, TAUWNY, XIN, DIA )
!/STX                    CALL W3SINX
                  END IF
                IF ( FLSRCE(3) ) THEN
!/NL1                    CALL W3SNL1 ( A, CG, WNMEAN*DEPTH,  XNL, DIA )
!/NL2                    CALL W3SNL2 ( A, CG, DEPTH,         XNL, DIA )
!/NL3                    CALL W3SNL3 ( A, CG, WN, DEPTH,     XNL, DIA )
!/NL4                    CALL W3SNL4 ( A, CG, WN, DEPTH,     XNL, DIA )
!/NLX                    CALL W3SNLX
!
!/TRX                    CALL W3STRX
                  END IF
                IF ( FLSRCE(4) ) THEN
!/ST1                    CALL W3SDS1 ( A, WN2, EMEAN, FMEAN, WNMEAN,  &
!/ST1                                                        XDS, DIA )
!/ST2                    CALL W3SDS2 ( A, CG, WN, FPI, USTAR,         &
!/ST2                                  ALPHA,                XDS, DIA )
!/ST3                    CALL W3SDS3 ( A, WN, CG, EMEAN, FMEANS, WNMEAN,  &
!/ST3                                  USTAR, USTD, DEPTH, XDS, DIA, IX, IY )
!/ST4                    CALL W3SDS4 ( A, WN, CG,  USTAR, USTD, DEPTH, XDS, &
!/ST4                                  DIA, IX, IY, LAMBDA, WHITECAP )
!/ST6                    CALL W3SDS6 ( A, CG, WN,            XDS, DIA )
!/ST6                    IF (SWL6S6) CALL W3SWL6 ( A, CG, WN,  XWL, DIA )
!/STX                    CALL W3SDSX
!
!/DB1                    CALL W3SDB1 ( I, A, DEPTH, EMEAN, FMEAN,    &
!/DB1                                  WNMEAN, CG, LBREAK, XDB, DIA )
!/DBX                    CALL W3SDBX
                  END IF
                IF ( FLSRCE(5) ) THEN
!/BT1                    CALL W3SBT1 ( A, CG, WN, DEPTH,     XBT, DIA )
!/BT2                    SBTC2 = 2. * -0.067 / GRAV
!/BT2                    CALL W3SBT2 ( A, CG, WN, DEPTH, XBT, DIA, SBTC2 )
!/BT4                    IX=1    ! to be fixed later 
!/BT4                    IY=1    ! to be fixed later 
!/BT4                    ISEA=1  ! to be fixed later 
!/BT4                    D50 = SED_D50(ISEA)
!/BT4                    PSIC= SED_PSIC(ISEA)

!/BT4                    CALL W3SBT4 ( A, CG, WN, DEPTH, D50, PSIC, TAUBBL,   & 
!/BT4                                       BEDFORM, XBT, DIA, IX, IY )

                 BT8MSG='ww3_outp: ITYPE=3 with BT8 or BT9: Sbot out'//&
                        'put is not yet supported. Use "F" for the 5'//&
                        'th T/F flag.'
!/BT8                    CALL EXTCDE( 516,MSG=BT8MSG)
!/BT9                    CALL EXTCDE( 516,MSG=BT8MSG)

! For info on this issue, see : "BT8&9 issues" in "Remarks" section above.

!...broken....!/BT8                    CALL W3SBT8 ( A, DEPTH, XBT, DIA, IX, IY )
!...broken....!/BT9                    CALL W3SBT9 ( A, DEPTH, XBT, DIA, IX, IY )

!/BTX                    CALL W3SBTX

!
!/BS1                    CALL W3SBS1 ( A, CG, WN, DEPTH,              &
!/BS1                           CAO(J)*COS(CDO(J)), CAO(J)*SIN(CDO(J)), &
!/BS1                                  TAUSCX, TAUSCY,   XBS, DIA )
!/BSX                    CALL W3SBSX
!
                  END IF
!
                IF ( FLSRCE(6) ) THEN
!/IS2                   IF (IICEDISP) THEN
!/IS2                     CALL LIU_FORWARD_DISPERSION  (ICETHICK,0.,DEPTH, &
!/IS2                                                   SIG,WN_R,CG_ICE,ALPHA_LIU)
!/IS2                   ELSE
!/IS2                     WN_R=WN
!/IS2                     CG_ICE=CG
!/IS2                     END IF
!
!/IS2                  CALL W3SIS2(A, DEPTH, ICECON, ICETHICK, ICEF, ICEDMAX, &
!/IS2                              IX, IY, XIS, DIA, DIA2, WN, CG, WN_R, CG_ICE, R)    
                  END IF
!
!/XXX                CALL W3SXXX
!
!/STAB2                UABS   = UABS * ASFAC
!
                IF ( ISCALE.EQ.0 .OR. ISCALE.EQ.3 ) THEN
                    FACF   = TPIINV
                    FACE   = 1.
                    FACS   = 1.
                  ELSE IF ( ISCALE.EQ.1 .OR. ISCALE.EQ.4 ) THEN
                    FACF   = TPIINV * UABS / GRAV
                    FACE   = GRAV**3 / UABS**5
                    FACS   = GRAV**2 / UABS**4
                  ELSE IF ( ISCALE.EQ.2 .OR. ISCALE.EQ.5 ) THEN
                    FACF   = TPIINV * USTAR / GRAV
                    FACE   = GRAV**3 / USTAR**5
                    FACS   = GRAV**2 / USTAR**4
                  END IF
!
                DO IK=1, NK
                  FACTOR = TPI / CG(IK) * SIG(IK)
                  E1  (IK) = 0.
                  SIN1(IK) = 0.
                  SNL1(IK) = 0.
                  SDS1(IK) = 0.
                  SBT1(IK) = 0.
                  STT1(IK) = 0.
                  SIS1(IK) = 0.
                  DO ITH=1, NTH
                    ISP         = ITH + (IK-1)*NTH
                    E  (IK,ITH) = SPCO(ISP,J)
                    SWN(IK,ITH) = ( XLN(ITH,IK) + XIN(ITH,IK) ) * FACTOR
                    SNL(IK,ITH) = ( XNL(ITH,IK) + XTR(ITH,IK) ) * FACTOR
                    SDS(IK,ITH) = ( XDS(ITH,IK) + XDB(ITH,IK) ) * FACTOR
!/ST6               SDS(IK,ITH) =   SDS(IK,ITH) +(XWL(ITH,IK)   * FACTOR)
                    SBT(IK,ITH) = ( XBT(ITH,IK) * XBS(ITH,IK) ) * FACTOR
                    SIS(IK,ITH) = XIS(ITH,IK) * FACTOR
                    STT(IK,ITH) = SWN(IK,ITH) + SNL(IK,ITH)+SDS(IK,ITH)&
                                  + SBT(IK,ITH) + SIS(IK,ITH) &
                                  + XXX(ITH,IK) * FACTOR
                    E1  (IK) = E1  (IK) + E(IK,ITH)
                    SIN1(IK) = SIN1(IK) + SWN(IK,ITH)
                    SNL1(IK) = SNL1(IK) + SNL(IK,ITH)
                    SDS1(IK) = SDS1(IK) + SDS(IK,ITH)
                    SBT1(IK) = SBT1(IK) + SBT(IK,ITH)
                    SIS1(IK) = SIS1(IK) + SIS(IK,ITH)
                    END DO
                  E1  (IK) = E1(IK)   * DTH * FACE
                  SIN1(IK) = SIN1(IK) * DTH * FACS
                  SNL1(IK) = SNL1(IK) * DTH * FACS
                  SDS1(IK) = SDS1(IK) * DTH * FACS
                  SBT1(IK) = SBT1(IK) * DTH * FACS
                  SIS1(IK) = SIS1(IK) * DTH * FACS
                  END DO
!
                STT1       = SIN1 + SNL1 + SDS1 + SBT1 + SIS1
                E1ALL(:,1) = SIN1
                E1ALL(:,2) = SNL1
                E1ALL(:,3) = SDS1
                E1ALL(:,4) = SBT1
                E1ALL(:,5) = SIS1
                E1ALL(:,6) = STT1
!
              END IF
!
! 4.a Perform output type 1 ( print plots / tables / file )
!
            IF ( ITYPE .EQ. 1 ) THEN
!
                IF ( OTYPE .EQ. 1 ) THEN
!
                    IF ( SCALE1 .GE. 0. )                             &
                        CALL PRT1DS (NDSO, NK, E1, SIG(1:NK), 'RAD/S',&
                             17, SCALE1, 'E(f)', 'm^2s', PTNME(J) )
                    IF ( SCALE2 .GE. 0. )                             &
                        CALL PRT2DS (NDSO, NK, NK, NTH, E, SIG(1:NK), &
                             'RAD/S', 1., SCALE2, 0.0001, 'E(f,th)',  &
                             'm^2s', PTNME(J) )
                    WRITE (NDSO,910) DPO(J), UABS
                    IF ( WAO(J) .GT. 0. ) WRITE (NDSO,911) UDIR
                    WRITE (NDSO,912) ASO(J), CAO(J)
                    IF ( CAO(J) .GT. 0. ) WRITE (NDSO,913) CDIR
                    WRITE (NDSO,914) HSIG, WLEN, TMEAN, THMEAN, THSPRD
!
                  ELSE IF ( OTYPE .EQ. 2 ) THEN
!
                    CALL STME21 ( TIME , DTME21 )
                    IF ( FLAGLL ) THEN
                        WRITE (NDSTAB,920) DTME21, PTNME(J),          &
                            M2KM*PTLOC(1,J), M2KM*PTLOC(2,J),         &
                            DPO(J), USTAR, WAO(J), UDIR
                      ELSE
                        WRITE (NDSTAB,720) DTME21, PTNME(J),          &
                            M2KM*PTLOC(1,J), M2KM*PTLOC(2,J),         &
                            DPO(J), USTAR, WAO(J), UDIR
                      END IF
                    IF ( FP .EQ. 0. ) FP = SIG(NK)
                    DO IK=1, NK
                      WRITE (NDSTAB,921) TPIINV*SIG(IK), SIG(IK)/FP,  &
                          E1(IK), THBND(IK), SPBND(IK), APM(IK)
                      END DO
                    IF ( FP .EQ. SIG(NK) ) FP = 0.
                    WRITE (NDSTAB,922)
!
                  ELSE IF ( OTYPE .EQ. 3 ) THEN
!
                    IF ( FLFORM ) THEN
                        WRITE (NDSTAB) PTNME(J), PTLOC(2,J),          &
                                       PTLOC(1,J), DPO(J), WAO(J),    &
                                       UDIR, CAO(J), CDIR
                        WRITE (NDSTAB) ((E(IK,ITH),IK=1,NK),ITH=1,NTH)
                      ELSE
                        WRITE (NDSTAB,901) PTNME(J), M2KM*PTLOC(2,J), &
                                           M2KM*PTLOC(1,J), DPO(J),   &
                                           WAO(J), UDIR, CAO(J), CDIR
                        WRITE (NDSTAB,902)                            &
                              ((E(IK,ITH),IK=1,NK),ITH=1,NTH)
                      END IF
!
                  END IF
!
! 4.b Perform output type 2 ( tables )
!
              ELSE IF ( ITYPE .EQ. 2 ) THEN
!
                IF ( NREQ .EQ. 1 ) THEN
!
                    IYR    = TIME(1) / 10000
                    IMTH   = MOD ( TIME(1) , 10000 ) / 100
                    IDY    = MOD ( TIME(1) , 100 )
                    IH     = TIME(2) / 10000
                    IM     = MOD ( TIME(2) , 10000 ) / 100
                    IS     = MOD ( TIME(2) , 100 )
                    IF ( OTYPE .EQ. 1 ) THEN
                        WRITE (NDSTAB,1921) TIME(1), IH, IM, IS,      &
                               DPO(J), CAO(J), CDIR, WAO(J), UDIR
                      ELSE IF ( OTYPE .EQ. 2 ) THEN
                        WRITE (NDSTAB,1922) TIME(1), IH, IM, IS,      &
                               HSIG, WLEN, TMEAN, THMEAN, THSPRD,     &
                               FP*TPIINV, THP, SPP
                      ELSE IF ( OTYPE.EQ.3 ) THEN
                        WRITE (NDSTAB,1923) TIME(1), IH, IM, IS,      &
                           UNORM, ESTAR, FPSTAR, CD*1000., APM(NK)*100.
                      ELSE IF ( OTYPE.EQ.4 ) THEN
                        WRITE (NDSTAB,1924) TIME(1), IH, IM, IS,      &
                           UNORM, ESTAR, FPSTAR, CD*1000., APM(NK)*100.
                      ELSE IF ( OTYPE.EQ.5 ) THEN
                        HMAT   = MIN ( 100. , 3.33*GRAV*HSIG/UABS**2 )
                        IF ( HSIG .GE. HSMIN ) THEN
                            CALL WAVNU1 ( FP, DPO(J), WNA, XYZ )
                            AGE1   = MIN ( 100. , FP / WNA / UABS )
                            AFR    = TPI / TMEAN
                            CALL WAVNU1 ( AFR, DPO(J), WNA, XYZ )
                            AGE2   = MIN ( 100. , AFR / WNA / UABS )
                          ELSE
                            AGE1   = -9.99
                            AGE2   = -9.99
                          END IF
                        WRITE (NDSTAB,1925) TIME(1), IH, IM, IS,      &
                               WAO(J), UDIR, HSIG, HMAT, AGE1, AGE2,  &
                               ASO(J)
                      ELSE IF ( OTYPE.EQ.6 ) THEN
                        IF ( HSIG .GE. HSMIN ) THEN
                           WRITE (NDSTAB,1926) IYR, IMTH, IDY, IH,    &
                                 WAO(J), NINT(UDIR), HSIG, TPI / FP
                        ELSE
                           WRITE (NDSTAB,1926) IYR, IMTH, IDY, IH,    &
                                 WAO(J), NINT(UDIR), HSIG, 0.0
                        END IF
                      END IF
!
                  ELSE
!
                    IF ( OTYPE .EQ. 1 ) THEN
                        IF ( FLAGLL ) THEN
                            WRITE (NDSTAB,1931) M2KM*PTLOC(1,J),      &
                                   M2KM*PTLOC(2,J), DPO(J), CAO(J),   &
                                   CDIR, WAO(J), UDIR
                          ELSE
                            WRITE (NDSTAB,1731) M2KM*PTLOC(1,J),      &
                                   M2KM*PTLOC(2,J), DPO(J), CAO(J),   &
                                   CDIR, WAO(J), UDIR
                          END IF
                      ELSE IF ( OTYPE .EQ. 2 ) THEN
                        IF ( FLAGLL ) THEN
                            WRITE (NDSTAB,1932) M2KM*PTLOC(1,J),      &
                                   M2KM*PTLOC(2,J), HSIG, WLEN,       &
                                   TMEAN, THMEAN, THSPRD, FP*TPIINV,  &
                                   THP, SPP
                          ELSE
                            WRITE (NDSTAB,1732) M2KM*PTLOC(1,J),      &
                                   M2KM*PTLOC(2,J), HSIG, WLEN,       &
                                   TMEAN, THMEAN, THSPRD, FP*TPIINV,  &
                                   THP, SPP
                          END IF
                      ELSE IF ( OTYPE .EQ. 3 ) THEN
                        WRITE (NDSTAB,1933) 1.E-4*XSTAR,              &
                               1.E-4*YSTAR, UNORM, ESTAR, FPSTAR,     &
                               CD*1000., APM(NK)*100.
                      ELSE IF ( OTYPE .EQ. 4 ) THEN
                        WRITE (NDSTAB,1934) XSTAR, YSTAR, UNORM,      &
                               ESTAR, FPSTAR, CD*1000., APM(NK)*100.
                      ELSE IF ( OTYPE .EQ. 5 ) THEN
                        HMAT   = MIN ( 100. , 3.33*GRAV*HSIG/UABS**2 )
                        CALL WAVNU1 ( FP, DPO(J), WNA, XYZ )
                        AGE1   = MIN ( 100. , FP / WNA / UABS )
                        AFR    = TPI / TMEAN
                        CALL WAVNU1 ( AFR, DPO(J), WNA, XYZ )
                        AGE2   = MIN ( 100. , AFR / WNA / UABS )
                        IF ( FLAGLL ) THEN
                            WRITE (NDSTAB,1935) M2KM*PTLOC(1,J),      &
                                   M2KM*PTLOC(2,J), WAO(J), UDIR,     &
                                   HSIG, HMAT, AGE1, AGE2, ASO(J)
                          ELSE
                            WRITE (NDSTAB,1735) M2KM*PTLOC(1,J),      &
                                   M2KM*PTLOC(2,J), WAO(J), UDIR,     &
                                   HSIG, HMAT, AGE1, AGE2, ASO(J)
                          END IF
                      ELSE IF ( OTYPE .EQ. 6 ) THEN
                        IF ( HSIG .GE. HSMIN ) THEN
                           IF ( FLAGLL ) THEN
                               WRITE (NDSTAB,1936) M2KM*PTLOC(1,J),       &
                                      M2KM*PTLOC(2,J), WAO(J), NINT(UDIR),&
                                      HSIG, TPI / FP
                             ELSE
                               WRITE (NDSTAB,1736) M2KM*PTLOC(1,J),       &
                                      M2KM*PTLOC(2,J), WAO(J), NINT(UDIR),&
                                      HSIG, TPI / FP
                             END IF
                        ELSE
                           IF ( FLAGLL ) THEN
                               WRITE (NDSTAB,1936) M2KM*PTLOC(1,J),       &
                                      M2KM*PTLOC(2,J), WAO(J), NINT(UDIR),&
                                      HSIG, 0.0
                             ELSE
                               WRITE (NDSTAB,1736) M2KM*PTLOC(1,J),       &
                                      M2KM*PTLOC(2,J), WAO(J), NINT(UDIR),&
                                      HSIG, 0.0
                             END IF
                        END IF
                      END IF
!
                  END IF
!
! 4.c Perform output type 3 ( source terms )
!
              ELSE IF ( ITYPE .EQ. 3 ) THEN
!
                IF ( OTYPE .EQ. 1 ) THEN
!
                    IF ( SCALE1 .GE. 0. ) THEN
                        IF ( FLSRCE(1) )                              &
                            CALL PRT1DS (NDSO, NK, E1, SIG(1:NK),     &
                                 'RAD/S', 17,  0., 'E(f)', 'm^2s',    &
                                 PTNME(J) )
                        IF (FLSRCE(2) .OR. FLSRCE(3) .OR.             &
                            FLSRCE(4) .OR. FLSRCE(5) .OR.             &
                            FLSRCE(6) .OR. FLSRCE(7) )                &
                            CALL PRT1DM (NDSO, NK, 6, E1ALL, SIG(1:NK),&
                                 'RAD/S', 17, SCALE1, VAR1, 'M2',     &
                                 PTNME(J) )
                      END IF
                    IF ( SCALE2 .GE. 0. ) THEN
                        IF ( FLSRCE(1) )                              &
                            CALL PRT2DS (NDSO, NK, NK, NTH, E,        &
                                 SIG(1:NK), 'RAD/S', 1., 0., 0.0001,  &
                                 'E(f,th)', 'm^2s', PTNME(J) )
                        IF ( FLSRCE(2) )                              &
                            CALL PRT2DS (NDSO, NK, NK, NTH, SWN,      &
                                SIG(1:NK), 'RAD/S', 1., SCALE2, 0.0001,&
                                 'Sin(f,th)', 'm^2', PTNME(J) )
                        IF ( FLSRCE(3) )                              &
                            CALL PRT2DS (NDSO, NK, NK, NTH, SNL,      &
                                SIG(1:NK), 'RAD/S', 1., SCALE2, 0.0001,&
                                 'Snl(f,th)', 'm^2', PTNME(J) )
                        IF ( FLSRCE(4) )                              &
                            CALL PRT2DS (NDSO, NK, NK, NTH, SDS,      &
                                SIG(1:NK), 'RAD/S', 1., SCALE2, 0.0001,&
                                 'Sds(f,th)', 'm^2', PTNME(J) )
                        IF ( FLSRCE(5) )                              &
                            CALL PRT2DS (NDSO, NK, NK, NTH, SBT,      &
                                SIG(1:NK), 'RAD/S', 1., SCALE2, 0.0001,&
                                 'Sbt(f,th)', 'm^2', PTNME(J) )
                        IF ( FLSRCE(6) )                              &
                            CALL PRT2DS (NDSO, NK, NK, NTH, SIS,      &
                                SIG(1:NK), 'RAD/S', 1., SCALE2, 0.0001,&
                                 'Sice(f,th)', 'm^2', PTNME(J) )
                        IF ( FLSRCE(7) )                              &
                            CALL PRT2DS (NDSO, NK, NK, NTH, STT,      &
                                SIG(1:NK), 'RAD/S', 1., SCALE2, 0.0001,&
                                 'Stot(f,th)', 'm^2', PTNME(J) )
                      END IF
!
                  ELSE IF ( OTYPE .EQ. 2 ) THEN
!
                    CALL STME21 ( TIME , DTME21 )
                    IF ( FLAGLL ) THEN
                        WRITE (NDSTAB,2920) DTME21, PTNME(J),         &
                            M2KM*PTLOC(1,J), M2KM*PTLOC(2,J),         &
                            DPO(J), USTAR, WAO(J)
                      ELSE
                        WRITE (NDSTAB,2720) DTME21, PTNME(J),         &
                            M2KM*PTLOC(1,J), M2KM*PTLOC(2,J),         &
                            DPO(J), USTAR, WAO(J)
                      END IF
                    IF ( ISCALE.EQ.0 ) THEN
                        WRITE (NDSTAB,2921)
                      ELSE IF ( ISCALE.EQ.1 .OR. ISCALE.EQ.2 ) THEN
                        WRITE (NDSTAB,2922)
                      ELSE IF ( ISCALE.EQ.3 ) THEN
                        WRITE (NDSTAB,2923)
                      ELSE IF ( ISCALE.EQ.4 .OR. ISCALE.EQ.5 ) THEN
                        WRITE (NDSTAB,2924)
                      END IF
                    IF ( ISCALE.GE.3 ) FACF = 1. / FP
                    DO IK=1, NK
                      WRITE (NDSTAB,2930) FACF*SIG(IK), E1(IK),       &
                        SIN1(IK), SNL1(IK), SDS1(IK), SBT1(IK),       &
                        SIS1(IK), STT1(IK)
                     
                      END DO
                    WRITE (NDSTAB,2940)
!
                  ELSE IF ( OTYPE .EQ. 3 ) THEN
!
                    CALL STME21 ( TIME , DTME21 )
                    IF ( FLAGLL ) THEN
                        WRITE (NDSTAB,2920) DTME21, PTNME(J),         &
                            M2KM*PTLOC(1,J), M2KM*PTLOC(2,J),         &
                            DPO(J), USTAR, WAO(J)
                      ELSE
                        WRITE (NDSTAB,2720) DTME21, PTNME(J),         &
                            M2KM*PTLOC(1,J), M2KM*PTLOC(2,J),         &
                            DPO(J), USTAR, WAO(J)
                      END IF
                    IF ( ISCALE.EQ.0 ) THEN
                        WRITE (NDSTAB,2925)
                      ELSE IF ( ISCALE.EQ.1 .OR. ISCALE.EQ.2 ) THEN
                        WRITE (NDSTAB,2926)
                      ELSE IF ( ISCALE.EQ.3 ) THEN
                        WRITE (NDSTAB,2927)
                      ELSE IF ( ISCALE.EQ.4 .OR. ISCALE.EQ.5 ) THEN
                        WRITE (NDSTAB,2928)
                      END IF
!
                    IF ( ISCALE.GE.3 ) FACF = 1. / FP
                    DO IK=1, NK
                      FACT   = 1. / MAX ( 1.E-10 , E1(IK) )
                      IF ( E1(IK) .GT. 1.E-10 ) THEN
                          WRITE (NDSTAB,2931) FACF*SIG(IK), E1(IK),   &
                            FACT*SIN1(IK), FACT*SNL1(IK),             &
                            FACT*SDS1(IK), FACT*SBT1(IK),             &
                            FACT*SIS1(IK),FACT*STT1(IK)
                        ELSE
                          WRITE (NDSTAB,2931) FACF*SIG(IK), E1(IK)
                        END IF
                      END DO
                    WRITE (NDSTAB,2940)
!
                  ELSE IF ( OTYPE .EQ. 4 ) THEN
!
                    IF ( FLFORM ) THEN
                        WRITE (NDSTAB) PTNME(J), PTLOC(2,J),          &
                                       PTLOC(1,J), DPO(J), WAO(J),    &
                                       UDIR, CAO(J), CDIR
                        IF ( FLSRCE(1) ) WRITE (NDSTAB)               &
                              ((E(IK,ITH),IK=1,NK),ITH=1,NTH)
                        IF ( FLSRCE(2) ) WRITE (NDSTAB)               &
                              ((SWN(IK,ITH),IK=1,NK),ITH=1,NTH)
                        IF ( FLSRCE(3) ) WRITE (NDSTAB)               &
                              ((SNL(IK,ITH),IK=1,NK),ITH=1,NTH)
                        IF ( FLSRCE(4) ) WRITE (NDSTAB)               &
                              ((SDS(IK,ITH),IK=1,NK),ITH=1,NTH)
                        IF ( FLSRCE(5) ) WRITE (NDSTAB)               &
                              ((SBT(IK,ITH),IK=1,NK),ITH=1,NTH)
                        IF ( FLSRCE(6) ) WRITE (NDSTAB)               &
                              ((SIS(IK,ITH),IK=1,NK),ITH=1,NTH)
                        IF ( FLSRCE(7) ) WRITE (NDSTAB)               &
                              ((STT(IK,ITH),IK=1,NK),ITH=1,NTH)
                      ELSE
                        IF ( FLAGLL ) THEN
                            WRITE (NDSTAB,901) PTNME(J),              &
                                M2KM*PTLOC(2,J), M2KM*PTLOC(1,J),     &
                                DPO(J), WAO(J), UDIR, CAO(J), CDIR
                          ELSE
                            WRITE (NDSTAB,701) PTNME(J),              &
                                M2KM*PTLOC(2,J), M2KM*PTLOC(1,J),     &
                                DPO(J), WAO(J), UDIR, CAO(J), CDIR
                          END IF
                        IF ( FLSRCE(1) ) WRITE (NDSTAB,902)           &
                              ((E(IK,ITH),IK=1,NK),ITH=1,NTH)
                        IF ( FLSRCE(2) ) WRITE (NDSTAB,902)           &
                              ((SWN(IK,ITH),IK=1,NK),ITH=1,NTH)
                        IF ( FLSRCE(3) ) WRITE (NDSTAB,902)           &
                              ((SNL(IK,ITH),IK=1,NK),ITH=1,NTH)
                        IF ( FLSRCE(4) ) WRITE (NDSTAB,902)           &
                              ((SDS(IK,ITH),IK=1,NK),ITH=1,NTH)
                        IF ( FLSRCE(5) ) WRITE (NDSTAB,902)           &
                              ((SBT(IK,ITH),IK=1,NK),ITH=1,NTH)
                        IF ( FLSRCE(6) ) WRITE (NDSTAB,902)           &
                              ((SIS(IK,ITH),IK=1,NK),ITH=1, NTH) 
                        IF ( FLSRCE(7) ) WRITE (NDSTAB,902)           &
                              ((STT(IK,ITH),IK=1,NK),ITH=1,NTH)
                      END IF
!
                  END IF
!
! 4.d Perform output type 4 ( Spectral partitions and bulletins )
!
              ELSE IF ( ITYPE .EQ. 4 ) THEN
!
                  IF ( OTYPE .EQ. 1 ) THEN
!
                    IF ( FLAGLL ) THEN
                        IF ( PTLOC(1,J) .LT. 0. )                     &
                            PTLOC(1,J) = PTLOC(1,J) + 360.
                        WRITE (NDSTAB,940) TIME, M2KM*PTLOC(2,J),     &
                            M2KM*PTLOC(1,J), PTNME(J), NPART, DEPTH,  &
                            WAO(J), UDIR, CAO(J), CDIR
                      ELSE
                        WRITE (NDSTAB,943) TIME, M2KM*PTLOC(1,J),     &
                            M2KM*PTLOC(2,J), PTNME(J), NPART, DEPTH,  &
                            WAO(J), UDIR, CAO(J), CDIR
                      END IF
!                   WRITE (NDSTAB,941)
                    DO I=0, NPART
                      WRITE (NDSTAB,942) I, XPART(:,I)
                      END DO
!
                  ELSEIF ( OTYPE .GE. 2 ) THEN
                      CALL W3BULL (NPART, XPART, DIMXP, UABS,         &
                                   UDIR, J, IOUT, TIMEV )
!
                    IF ( FLAGLL ) THEN
                      X = M2KM * PTLOC(1,J)
                      Y = M2KM * PTLOC(2,J)

                      X      = MOD ( X+720. , 360. )
                      IF ( X .LE. 180. ) THEN
                        IDLON  = 'E'
                      ELSE
                        X      = 360. - X
                        IDLON  = 'W'
                      ENDIF
                      !IF ( ABS(Y) .LE. 0.0049 ) THEN
                       !IDLAT  = '-'
                      IF ( Y .GE. 0. ) THEN
                        IDLAT  = 'N'
                      ELSE
                        IDLAT  = 'S'
                        Y      = -Y
                      ENDIF
                    ELSE
                      IDLAT = ' '
                      IDLON = ' '
                    ENDIF
                      IF ( OTYPE .EQ. 2 .OR. OTYPE .EQ. 4 ) THEN
                        NDSBUL=NDSTAB + (J - 1)
!/NCO                         NDSCBUL=NDSTAB + (J - 1) + NOPTS
                        IF (IOUT .EQ. 1) THEN
                          WRITE(HSTR,'(I2,1X,A)') TIMEV(2)/10000,   &
                                   HTYPE 
                          WRITE (NDSBUL,970) PTNME(J), Y, IDLAT, X, &
                                           IDLON, GNAME, TIMEV(1),  &
                                           HSTR
                          WRITE (NDSBUL,971) 
                          WRITE (NDSBUL,972)
                          WRITE (NDSBUL,971)
!/NCO                           WRITE (NDSCBUL,960) PTNME(J), Y, IDLAT,   &
!/NCO                              X, IDLON, GNAME, TIMEV(1), HSTR
!/NCO                           WRITE (NDSCBUL,961)
                        ENDIF
                        WRITE (NDSBUL,973) ASCBLINE
!/NCO                         WRITE (NDSCBUL,963) CASCBLINE
                     ENDIF
                      IF ( OTYPE .EQ. 3 .OR. OTYPE .EQ. 4 ) THEN
                        ICSV = 0 
                        NDSCSV = NDSTAB + (J - 1) + 2*NOPTS
                        WRITE (NDSCSV,'(A664)') CSVBLINE
                      ENDIF
                  END IF
!
              END IF
! ... End of fields loop
!
          END IF
        END DO
!
      RETURN
!
! Formats
!
  900 FORMAT (I8.8,I7.6)
  901 FORMAT ('''',A10,'''',2F7.2,F10.1,2(F7.2,F6.1))
  701 FORMAT ('''',A10,'''',2(F8.1,'E3'),F10.1,2(F7.2,F6.1))
  902 FORMAT (7E11.3)
  905 FORMAT (9X,A)
  910 FORMAT (/15X,' Water depth       :',F7.1,'  (m)'/               &
               15X,' Wind speed        :',F8.2,' (m/s)')
  911 FORMAT ( 15X,' Wind direction    :',F7.1,'  (degr)')
  912 FORMAT ( 15X,' Air-sea temp. dif.:',F7.1,'  (degr)'/            &
               15X,' Current speed     :',F8.2,' (m/s)')
  913 FORMAT ( 15X,' Current direction :',F7.1,'  (degr)')
  914 FORMAT ( 15X,' Wave height       :',F8.2,' (m)'/                &
               15X,' Mean wave length  :',F6.0,'   (m)'/              &
               15X,' Mean wave period  :',F7.1,'  (s)'/               &
               15X,' Mean wave direct. :',F7.1,'  (degr)'/            &
               15X,' Direct. spread    :',F7.1,'  (degr)'/)
  920 FORMAT (' Time     : ',A/                                  &
              ' Location : ',A,'  (',2F8.2,' )'/                 &
              ' depth    : ',F7.1,'   m'/                        &
              ' U*       : ',F9.3,' m/s'/                        &
              ' U10      : ',F7.1,'   m/s'/                      &
              ' Dir U10  : ',F7.1,'   degr'//                    &
       '      f     f/fp      F(f)    theta    spr    alpha  '/  &
       '     (Hz)    (-)     (m2s)    (deg)   (deg)    (-)   '/  &
           '  --------------------------------------------------')
  720 FORMAT (' Time     : ',A/                                  &
              ' Location : ',A,'  (',2(F8.1,'E3'),' )'/          &
              ' depth    : ',F7.1,'   m'/                        &
              ' U*       : ',F9.3,' m/s'/                        &
              ' U10      : ',F7.1,'   m/s'/                      &
              ' Dir U10  : ',F7.1,'   degr'//                    &
       '      f     f/fp      F(f)    theta    spr    alpha  '/  &
       '     (Hz)    (-)     (m2s)    (deg)   (deg)    (-)   '/  &
       '  --------------------------------------------------')
  921 FORMAT (1x,F8.5,F7.3,E11.3,2F8.1,F8.4)
  922 FORMAT (' '/' ')
!
  940 FORMAT (1X,I8.8,1X,I6.6,2F8.3,2X,'''',A10,'''',            &
              1X,I3,F7.1,F5.1,f6.1,F5.2,F6.1)
  943 FORMAT (1X,I8.8,1X,I6.6,2(F8.1,'E3'),2X,'''',A10,'''',     &
              1X,I3,F7.1,F5.1,f6.1,F5.2,F6.1)
  941 FORMAT ('        hs     tp     lp       theta     sp      wf')
  942 FORMAT (I3,3F8.2,2F9.2,10F7.2)
!
! 
!/NCO  960 FORMAT ( 'Location : ',A,' (',F5.2,A,1X,F6.2,A,')'/              &
!/NCO                'Model    : ',A/                                       &
!/NCO                'Cycle    : ',I8,1X,A//                                &
!/NCO                'DDHH HS SS PP DDD SS PP DDD SS PP DDD',               &
!/NCO                       ' SS PP DDD SS PP DDD SS PP DDD')
!/NCO  961 FORMAT ('----------------------------------------',              &
!/NCO              '---------------------------')
!/NCO  963 FORMAT (A)
!
 970 FORMAT ( '  Location : ',A,' (',F5.2,A,1X,F6.2,A,')'/            &
              '  Model    : ',A/                                      &
              '  Cycle    : ',I8,1X,A)
 971 FORMAT (' +-------+-----------+-----------------+',              &
             '-----------------+-----------------+----',              &
             '-------------+-----------------+--------',              &
             '---------+')
 972 FORMAT (' | day & |  Hst  n x |    Hs   Tp  dir |',              &
                                   '    Hs   Tp  dir |',              &
                                   '    Hs   Tp  dir |',              &
                                   '    Hs   Tp  dir |',              &
                                   '    Hs   Tp  dir |',              &
                                   '    Hs   Tp  dir |'/              &
             ' |  hour |  (m)  - - |    (m)  (s) (d) |',              &
                                   '    (m)  (s) (d) |',              &
                                   '    (m)  (s) (d) |',              &
                                   '    (m)  (s) (d) |',              &
                                   '    (m)  (s) (d) |',              &
                                   '    (m)  (s) (d) |')
 973 FORMAT (1X,A)
!
 1901 FORMAT (                                                        &
       '    Date     Time        d     Uc     Dir.  U10    Dir. '/    &
       '            h  m  s     (m)   (m/s)  (d.N) (m/s)  (d.N) '/    &
       ' ---------------------------------------------------------')
 1902 FORMAT (                                                        &
       '    Date     Time        Hs     L      Tr    Dir.  Spr. ',    &
       '    fp    p_dir  p_spr'/                                      &
       '            h  m  s     (m)    (m)    (s)   (d.N)  (deg)',    &
       '   (Hz)   (d.N)  (deg)'/                                      &
       ' -------------------------------------------------------',    &
       '-----------------------')
 1903 FORMAT (                                                        &
    '    Date     Time       U*       E*        fp*       Cd   alpha'/&
    '            h  m  s   (m/s)     (-)        (-)     *1000  *100'/ &
    ' --------------------------------------------------------------')
 1904 FORMAT (                                                        &
     '    Date     Time     U10       E*        fp*       Cd   alpha'/&
     '            h  m  s  (m/s)     (-)        (-)     *1000  *100'/ &
     ' --------------------------------------------------------------')
 1905 FORMAT (                                                        &
       '   Date     Time     U10    Dir.    Hs      H*    cp/U  ',    &
           '  cm/U     Dt'/                                           &
       '                    (m/s)  (d.N)   (m)     (-)    (-)   ',    &
           '   (-)   (deg)'/                                          &
      ' --------------------------------------------------',          &
       '---------------------')
 1906 FORMAT (                                                         &
     '     Time      U10    Dir. Hs  Tp  '/         &
     '  yr mth dy h  (m/s) (d.N) (m) (s) '/         &
     ' ----------------------------------')
 1910 FORMAT ( ' '/' ' )
 1911 FORMAT (' Time : ',A//                                     &
   '    Long.    Lat.       d     Uc     Dir.  U10    Dir. '/    &
   '                       (m)   (m/s)  (d.N) (m/s)  (d.N) '/    &
   ' --------------------------------------------------------')
 1912 FORMAT (' Time : ',A//                                     &
   '    Long.    Lat.       Hs     L      Tr    Dir.  Spr. ',    &
   '    fp    p_dir  p_spr'/                                     &
   '                       (m)    (m)    (s)   (d.N)  (deg)',    &
   '   (Hz)   (d.N)  (deg)'/                                     &
   ' ------------------------------------------------------',    &
   '-----------------------')
 1711 FORMAT (' Time : ',A//                                     &
   '       X        Y          d     Uc     Dir.  U10    Dir. '/ &
   '      (m)      (m)        (m)   (m/s)  (d.N) (m/s)  (d.N) '/ &
   ' ----------------------------------------------------------')
 1712 FORMAT (' Time : ',A//                                     &
   '       X        Y          Hs     L      Tr    Dir.  Spr. ', &
   '    fp    p_dir  p_spr'/                                     &
   '      (m)      (m))       (m)    (m)    (s)   (d.N)  (deg)', &
   '   (Hz)   (d.N)  (deg)'/                                     &
   ' ------------------------------------------------------',    &
   '-------------------------')
 1913 FORMAT (' Time : ',A//                                          &
   '       X*       Y*       U*       E*        fp*       Cd   alpha'/&
   '      (-)      (-)     (m/s)     (-)        (-)     *1000  *100'/ &
   ' --------------------------------------------------------------')
 1914 FORMAT (' Time : ',A//                                          &
   '       X*       Y*     U10       E*        fp*       Cd   alpha'/ &
   '      (-)      (-)    (m/s)     (-)        (-)     *1000  *100 '/ &
   ' --------------------------------------------------------------')
 1915 FORMAT (' Time : ',A//                                     &
   '     Long.    Lat.   U10    Dir.    Hs      H*    cp/U  ',   &
       '  cm/U     Dt'/                                          &
   '                    (m/s)  (d.N)   (m)     (-)    (-)   ',   &
       '   (-)   (deg)'/                                         &
   ' -------------------------------------------------',         &
   '---------------------')
 1715 FORMAT (' Time : ',A//                                     &
   '       X        Y      U10    Dir.    Hs      H*    cp/U  ', &
       '  cm/U     Dt'/                                          &
   '      (m)      (m)    (m/s)  (d.N)   (m)     (-)    (-)   ', &
       '   (-)   (deg)'/                                         &
   ' ---------------------------------------------------',       &
   '---------------------')
 1916 FORMAT (' Time : ',A//                                     &
   '   Long.     Lat.   U10   Dir.     Hs      Tp  '/         &
   '                 (m/s)   (d.N)    (m)     (s)  '/         &
   '-----------------------------------------------')
 1716 FORMAT (' Time : ',A//                                     &
   '     X      Y      U10      Dir.      Hs       Tp  '/         &
   '    (m)    (m)    (m/s)    (d.N)     (m)      (s)  '/         &
   '---------------------------------------------------')
 1921 FORMAT ( 2X,I8.8,I3,2(1X,I2.2),F10.1,F6.2,F7.1,F6.2,F7.1)
 1922 FORMAT ( 2X,I8.8,I3,2(1X,I2.2),F9.3,F7.1,F7.2,F7.1,F7.2,        &
               F8.4,F7.1,F7.2)
 1923 FORMAT ( 2X,I8.8,I3,2(1X,I2.2),F8.4,2E11.3,2F7.3)
 1924 FORMAT ( 2X,I8.8,I3,2(1X,I2.2),F7.1,2E11.3,2F7.3)
 1925 FORMAT ( 2X,I8.8,I3,2(1X,I2.2),F7.2,F7.1,2F7.2,2F8.2,F7.1)
 1926 FORMAT ( 2X,I4,3(1X,I2),F6.2,1X,I3,2F6.2)
 1931 FORMAT ( 2X,2F8.3,F10.1,F6.2,F7.1,F6.2,F7.1)
 1932 FORMAT ( 2X,2F8.3,F9.3,F7.1,F7.2,F7.1,F7.2,                &
               F8.4,F7.1,F7.2)
 1731 FORMAT ( 2X,2(F7.1,'E3'),F10.1,F6.2,F7.1,F6.2,F7.1)
 1732 FORMAT ( 2X,2(F7.1,'E3'),F9.3,F7.1,F7.2,F7.1,F7.2,         &
               F8.4,F7.1,F7.2)
 1933 FORMAT ( 2X,2(F7.1,'E4'),F8.4,2E11.3,2F7.3)
 1934 FORMAT ( 2X,2F9.1,F7.1,2E11.3,2F7.3)
 1935 FORMAT ( 2X,2F8.3,F7.2,F7.1,2F7.2,2F8.2,F7.1)
 1735 FORMAT ( 2X,2(F7.1,'E3'),F7.2,F7.1,2F7.2,2F8.2,F7.1)
 1936 FORMAT ( 2X,2F8.3,F6.2,1X,I3,2F6.2)
 1736 FORMAT ( 2X,2(F8.2,'E3'),F6.2,1X,I3,2F6.2)
!
 2920 FORMAT (' Time     : ',A/                                  &
              ' Location : ',A,'  (',2F8.2,' )'/                 &
              ' depth    : ',F7.1,'   m'/                        &
              ' U*       : ',F9.3,' m/s'/                        &
              ' U10      : ',F7.1,'   m/s'/)
 2720 FORMAT (' Time     : ',A/                                  &
              ' Location : ',A,'  (',2(F8.1,'E3'),' )'/          &
              ' depth    : ',F7.1,'   m'/                        &
              ' U*       : ',F9.3,' m/s'/                        &
              ' U10      : ',F7.1,'   m/s'/)
 2921 FORMAT ('    f           E      ',                              &
              '    Sin        Snl        Sds        Sbt       Sice       Stot'/  &
              '   (Hz)       (m2s)    ',                              &
              '   (m2)       (m2)       (m2)       (m2)       (m2)       (m2)'/  &
              ' ------------------------------------------',               &
              '-------------------------------------------')
 2922 FORMAT ('    f*          E*    ',                               &
              '    Sin*       Snl*       Sds*       Sbt*      Sice*      Stot*'/ &
              '   (-)         (-)    ',                               &
              '    (-)        (-)        (-)        (-)        (-)        (-)'/  &
              ' ------------------------------------------',               &
              '-------------------------------------------')
 2923 FORMAT ('   f/fp         E      ',                              &
              '   Sin        Snl        Sds        Sbt       Sice      Stot'/   &
              '   (-)        (m2s)    ',                              &
              '  (m2)       (m2)       (m2)       (m2)       (m2)      (m2)'/   &
              ' ------------------------------------------',               &
              '-------------------------------------------')
 2924 FORMAT ('   f/fp         E*     ',                              &
              '   Sin*       Snl*       Sds*       Sbt*      Sice*      Stot*'/  &
              '     (-)       (-)     ',                              &
              '   (-)        (-)        (-)        (-)        (-)        (-)'/   &
              ' ------------------------------------------',               &
              '-------------------------------------------')
 2925 FORMAT ('    f         E      ',                                &
              '  Tini       Tnli       Tdsi       Tbti      Ticei      Ttoti'/   &
              '   (Hz)     (m2s)    ',                                &
              ' (1/s)      (1/s)      (1/s)      (1/s)      (1/s)      (1/s)'/   &
              ' ----------------------------------------',                 &
              '-------------------------------------------')
 2926 FORMAT ('    f*        E*     ',                                &
              ' Tini*      Tnli*      Tdsi*      Tbti*     Ticei*      Ttoti*'/   &
              '   (-)       (-)     ',                                &
              '  (-)        (-)        (-)        (-)        (-)        (-)'/    &
              ' ----------------------------------------',                 &
              '-------------------------------------------')
 2927 FORMAT ('   f/fp       E      ',                                &
              '  Tini       Tnli       Tdsi       Tbti      Ticei      Ttoti'/   &
              '   (-)       (m2s)    ',                               &
              ' (1/s)      (1/s)      (1/s)      (1/s)      (1/s)      (1/s)'/   &
              ' ----------------------------------------',                 &
              '-------------------------------------------')
 2928 FORMAT ('   f/fp       E*     ',                                &
              ' Tini*      Tnli*      Tdsi*      Tbti*     Ticei*      Ttoti*'/   &
              '   (-)       (-)     ',                                &
              '  (-)        (-)        (-)        (-)        (-)        (-)'/    &
              ' ----------------------------------------',                 &
              '-------------------------------------------')
 2930 FORMAT (1X,F6.4,2X,7E11.3)
 2931 FORMAT (1X,F6.4,7E11.3)
 2940 FORMAT ( ' '/' ' )
!
!/T 9000 FORMAT (' TEST W3EXPO : FLAGS :',40L2)
!/T 9001 FORMAT (' TEST W3EXPO : ITPYE  :',I4/                        &
!/T              '               OTPYE  :',I4/                        &
!/T              '               NREQ   :',I4/                        &
!/T              '               SCALE1 :',E10.3/                     &
!/T              '               SCALE2 :',E10.3/                     &
!/T              '               FLSRCE :',7L2)
!/T 9002 FORMAT (' TEST W3EXPO : OUTPUT POINT : ',A)
!/T 9010 FORMAT (' TEST W3EXPO : DEPTH =',F7.1,'  IK, T, K, CG :')
!/T 9011 FORMAT ('               ',I3,F8.2,F8.4,F8.2)
!/
!/ End of W3EXPO ----------------------------------------------------- /
!/
      END SUBROUTINE W3EXPO
!/
      SUBROUTINE W3IOPO_dynm ( NDSOP, TOUT, IOTST)
!/
!/                  +-----------------------------------+
!/                  | WAVEWATCH III           NOAA/NCEP |
!/                  |           H. L. Tolman            |
!/                  |                        FORTRAN 90 |
!/                  | Last update :         25-Jul-2006 |
!/                  +-----------------------------------+
!/
!/    07-Jan-1999 : Distributed FORTRAN 77 version.     ( version 1.18 )
!/    30-Dec-1999 : Upgrade to FORTRAN 90               ( version 2.00 )
!/                  Major changes to logistics.
!/    10-Nov-2004 : Multiple grid version.              ( version 3.06 )
!/    27-Jun-2006 : Adding file name preamble.          ( version 3.09 )
!/    25-Jul-2006 : Adding grid ID per point.           ( version 3.10 )
!/    27-Aug-2015 : Adding interpolation for the ice.   ( version 5.10 )
!/
!  1. Purpose :
!
!     Read/write point output.
!
!  3. Parameters :
!
!     Parameter list
!     ----------------------------------------------------------------
!       VEROPT  C*10  Private  Point output file version number.
!       IDSTR   C*32  Private  Point output file ID string.
!       INXOUT  C*(*)  I   Test string for read/write, valid are:
!                          'READ' and 'WRITE'.
!       NDSOP   Int.   I   File unit number.
!       IOTST   Int.   O   Test indictor for reading.
!                           0 : Data read.
!                          -1 : Past end of file.
!       IMOD    I(O)   I   Model number for W3GDAT etc.
!     ----------------------------------------------------------------
!
!  4. Subroutines used :
!
!     See module documentation.
!
!  5. Called by :
!
!      Name      Type  Module   Description
!     ----------------------------------------------------------------
!      W3WAVE    Subr. W3WAVEMD Actual wave model routine.
!      WW3_OUTP  Prog.   N/A    Postprocessing for point output.
!      GX_OUTP   Prog.   N/A    Grads postprocessing for point output.
!     ----------------------------------------------------------------
!
!  6. Error messages :
!
!       Tests on INXOUT, file status and on array dimensions.
!
!  7. Remarks :
!
!     - The output file has the pre-defined name 'out_pnt.FILEXT'.
!     - In MPP version of model data is supposed to be gatherd at the
!       correct processor before the routine is called.
!     - No error output filtering needed.
!
!  8. Structure :
!
!     See source code.
!
!  9. Switches :
!
!     !/SHRD  Switch for shared / distributed memory architecture.
!     !/DIST  Id.
!
!     !/S     Enable subroutine tracing.
!     !/T     Test output.
!
! 10. Source code :
!
!/ ------------------------------------------------------------------- /
      USE W3GDATMD, ONLY: W3SETG
      USE W3WDATMD, ONLY: W3SETW
      USE W3ODATMD, ONLY: W3SETO, W3DMO2
!/
      USE W3GDATMD, ONLY: NTH, NK, NSPEC, FILEXT
      USE W3WDATMD, ONLY: TIME
      USE W3ODATMD, ONLY: NDST, NDSE, IPASS => IPASS2, NOPTS, IPTINT, &
                          IL, IW, II, PTLOC, PTIFAC, DPO, WAO, WDO,   &
                          ASO, CAO, CDO, SPCO, PTNME, O2INIT, FNMPRE, &
                          GRDID, ICEO, ICEHO, ICEFO
      USE W3ODATMD, ONLY :  OFILES
!/
!/SETUP      USE W3ODATMD, ONLY: ZET_SETO
!/
      USE W3SERVMD, ONLY: EXTCDE
!/S      USE W3SERVMD, ONLY: STRACE
!
      IMPLICIT NONE
!/
!/ Private parameter statements (ID strings)
!/
      CHARACTER(LEN=10), PARAMETER :: VEROPT = 'III  2.01 '
      CHARACTER(LEN=31), PARAMETER ::                        &
                           IDSTR = 'WAVEWATCH III POINT OUTPUT FILE'
!/
!/ ------------------------------------------------------------------- /
!/ Parameter list
!/
      INTEGER, INTENT(IN)           :: NDSOP
      INTEGER, INTENT(IN)           :: TOUT(2)
      INTEGER, INTENT(OUT)          :: IOTST
!/
!/ ------------------------------------------------------------------- /
!/ local parameters
!/
      INTEGER                 :: IGRD, IERR, MK, MTH, I, J
!/S      INTEGER, SAVE           :: IENT = 0
      LOGICAL,SAVE            :: WRITE
      CHARACTER(LEN=31)       :: IDTST
      CHARACTER(LEN=10)       :: VERTST
!/
      CHARACTER(LEN=15) :: TIMETAG
      CHARACTER(LEN=120) :: FILENAME
!/
!/ ------------------------------------------------------------------- /
!/
!/S      CALL STRACE (IENT, 'W3IOPO')
      IOTST  = 0
!
! test input parameters ---------------------------------------------- *
!
      if(TOUT(2) == 0) then
        WRITE(FILENAME,'(A,"out_pnt.points." I0,A,I0,A)') &
          trim(prefix),TOUT(1),'.0',TOUT(2),'0000'
      elseif(TOUT(2) .lt. 100000) then
        WRITE(FILENAME,'(A,"out_pnt.points." I0,A,I0)') &
          trim(prefix),TOUT(1),'.0',TOUT(2)
      else
        WRITE(FILENAME,'(A,"out_pnt.points." I0,A,I0)') &
          trim(prefix),TOUT(1),'.',TOUT(2)
      endif
!
      IGRD = 1
      CALL W3SETO ( IGRD, NDSE, NDST )
      CALL W3SETG ( IGRD, NDSE, NDST )
      CALL W3SETW ( IGRD, NDSE, NDST )
!
! open file ---------------------------------------------------------- *
!
      OPEN (NDSOP,FILE=trim(FILENAME),FORM='UNFORMATTED'&
            ,ERR=800,IOSTAT=IERR,STATUS='OLD')
!
      REWIND ( NDSOP )
!
      READ (NDSOP,END=801,ERR=802,IOSTAT=IERR) IDTST, VERTST, MK, &
            MTH, NOPTS
!
      IF ( IDTST .NE. IDSTR ) THEN
        WRITE (NDSE,902) IDTST, IDSTR
        CALL EXTCDE ( 10 )
      END IF
      IF ( VERTST .NE. VEROPT ) THEN
        WRITE (NDSE,903) VERTST, VEROPT
        CALL EXTCDE ( 11 )
      END IF
      IF (NK.NE.MK .OR. NTH.NE.MTH) THEN
        WRITE (NDSE,904) MK, MTH, NK, NTH
        CALL EXTCDE ( 12 )
      END IF
      IF ( .NOT. O2INIT ) CALL W3DMO2 ( IGRD, NDSE, NDST, NOPTS )
!
      READ  (NDSOP,END=801,ERR=802,IOSTAT=IERR)               &
                    ((PTLOC(J,I),J=1,2),I=1,NOPTS), (PTNME(I),I=1,NOPTS)
!
! TIME --------------------------------------------------------------- *
!
      READ (NDSOP,END=803,ERR=802,IOSTAT=IERR) TIME
!
!/T      WRITE (NDST,9010) TIME
!
!
! Loop over spectra -------------------------------------------------- *
!
      DO I=1, NOPTS
!
        READ (NDSOP,END=801,ERR=802,IOSTAT=IERR)                 &
              IW(I), II(I), IL(I), DPO(I), WAO(I), WDO(I),      &
!/SETUP       ZET_SETO(I),                                      &
              ASO(I), CAO(I), CDO(I), ICEO(I), ICEHO(I),        &
              ICEFO(I), GRDID(I), (SPCO(J,I),J=1,NSPEC)
!
      END DO
      CLOSE (NDSOP)
!
      RETURN
!
! Escape locations read errors
!
  800 CONTINUE
      WRITE (NDSE,1000) IERR
      CALL EXTCDE ( 20 )
!
  801 CONTINUE
      WRITE (NDSE,1001)
      CALL EXTCDE ( 21 )
!
  802 CONTINUE
      WRITE (NDSE,1002) IERR
      CALL EXTCDE ( 22 )
!
  803 CONTINUE
      IOTST  = -1
!/T      WRITE (NDST,9011)
      RETURN
!
! Formats
!
  900 FORMAT (/' *** WAVEWATCH III ERROR IN W3IOPO2 :'/               &
               '     ILEGAL INXOUT VALUE: ',A/)
  901 FORMAT (/' *** WAVEWATCH III ERROR IN W3IOPO2 :'/               &
               '     MIXED READ/WRITE, LAST REQUEST: ',A/)
  902 FORMAT (/' *** WAVEWATCH III ERROR IN W3IOPO2 :'/               &
               '     ILEGAL IDSTR, READ : ',A/                        &
               '                  CHECK : ',A/)
  903 FORMAT (/' *** WAVEWATCH III ERROR IN W3IOPO2 :'/               &
               '     ILEGAL VEROPT, READ : ',A/                       &
               '                   CHECK : ',A/)
  904 FORMAT (/' *** WAVEWATCH III ERROR IN W3IOPO2 :'/               &
               '     ERROR IN SPECTRA, MK, MTH : ',2I8/               &
               '              ARRAY DIMENSIONS : ',2I8/)
!
 1000 FORMAT (/' *** WAVEWATCH III ERROR IN W3IOPO2 : '/              &
               '     ERROR IN OPENING FILE'/                          &
               '     IOSTAT =',I5/)
 1001 FORMAT (/' *** WAVEWATCH III ERROR IN W3IOPO2 : '/              &
               '     PREMATURE END OF FILE'/)
 1002 FORMAT (/' *** WAVEWATCH III ERROR IN W3IOPO2 : '/              &
               '     ERROR IN READING FROM FILE'/                     &
               '     IOSTAT =',I5/)
!
!/T 9000 FORMAT (' TEST W3IOPO : IPASS =',I4,'    INXOUT = ',A,       &
!/T              ' WRITE = ',L1,' UNIT =',I3/                         &
!/T              '               IGRD =',I3,' FEXT = ',A)

!/T 9001 FORMAT (' TEST W3IOPO : OPENING NEW FILE [',A,']')
!/T 9002 FORMAT (' TEST W3IOPO : TEST PARAMETERS:'/                   &
!/T              '       IDSTR : ',A/                                 &
!/T              '      VEROPT : ',A/                                 &
!/T              '      NK,NTH :',I5,I8/                              &
!/T              '        NOPT :',I5)
!/T 9003 FORMAT (' TEST W3IOPO : POINT LOCATION AND ID')
!/T 9004 FORMAT (3X,I4,2F10.2,2X,A)
!
!/T 9010 FORMAT (' TEST W3IOPO : TIME  :',I9.8,I7.6)
!/T 9011 FORMAT (' TEST W3IOPO : END OF FILE REACHED')
!
!/T 9020 FORMAT (' TEST W3IOPO : POINT NR.:',I5)
!/T 9021 FORMAT (' TEST W3IOPO :',2I4,2F6.3)
!/T 9022 FORMAT (' TEST W3IOPO :',4I7,2X,4I2,2X,4F5.2)
!/T 9030 FORMAT (' TEST W3IOPO :',F8.1,2(F7.2,F7.1))
!/
!/ End of W3IOPO ----------------------------------------------------- /
!/
      END SUBROUTINE W3IOPO_dynm
!/
!/ End of W3OUTP ----------------------------------------------------- /
!/
      END PROGRAM W3OUTP