text:
C$Procedure PCKW20 ( PCK, write segment, type 20 )

      SUBROUTINE PCKW20 (  HANDLE,  CLSSID,  FRAME,   FIRST,
     .                     LAST,    SEGID,   INTLEN,  N,
     .                     POLYDG,  CDATA,   ASCALE,  TSCALE,
     .                     INITJD,  INITFR                   )

C$ Abstract
C
C     Write a type 20 segment to a PCK file.
C
C$ Disclaimer
C
C     THIS SOFTWARE AND ANY RELATED MATERIALS WERE CREATED BY THE
C     CALIFORNIA INSTITUTE OF TECHNOLOGY (CALTECH) UNDER A U.S.
C     GOVERNMENT CONTRACT WITH THE NATIONAL AERONAUTICS AND SPACE
C     ADMINISTRATION (NASA). THE SOFTWARE IS TECHNOLOGY AND SOFTWARE
C     PUBLICLY AVAILABLE UNDER U.S. EXPORT LAWS AND IS PROVIDED "AS-IS"
C     TO THE RECIPIENT WITHOUT WARRANTY OF ANY KIND, INCLUDING ANY
C     WARRANTIES OF PERFORMANCE OR MERCHANTABILITY OR FITNESS FOR A
C     PARTICULAR USE OR PURPOSE (AS SET FORTH IN UNITED STATES UCC
C     SECTIONS 2312-2313) OR FOR ANY PURPOSE WHATSOEVER, FOR THE
C     SOFTWARE AND RELATED MATERIALS, HOWEVER USED.
C
C     IN NO EVENT SHALL CALTECH, ITS JET PROPULSION LABORATORY, OR NASA
C     BE LIABLE FOR ANY DAMAGES AND/OR COSTS, INCLUDING, BUT NOT
C     LIMITED TO, INCIDENTAL OR CONSEQUENTIAL DAMAGES OF ANY KIND,
C     INCLUDING ECONOMIC DAMAGE OR INJURY TO PROPERTY AND LOST PROFITS,
C     REGARDLESS OF WHETHER CALTECH, JPL, OR NASA BE ADVISED, HAVE
C     REASON TO KNOW, OR, IN FACT, SHALL KNOW OF THE POSSIBILITY.
C
C     RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF
C     THE SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY
C     CALTECH AND NASA FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE
C     ACTIONS OF RECIPIENT IN THE USE OF THE SOFTWARE.
C
C$ Required_Reading
C
C     DAF
C     NAIF_IDS
C     TIME
C     PCK
C     ROTATION
C
C$ Keywords
C
C     ORIENTATION
C
C$ Declarations

      IMPLICIT NONE
      INCLUDE 'spk20.inc'

      INTEGER               HANDLE
      INTEGER               CLSSID
      CHARACTER*(*)         FRAME
      DOUBLE PRECISION      FIRST
      DOUBLE PRECISION      LAST
      CHARACTER*(*)         SEGID
      DOUBLE PRECISION      INTLEN
      INTEGER               N
      INTEGER               POLYDG
      DOUBLE PRECISION      CDATA  ( * )
      DOUBLE PRECISION      ASCALE
      DOUBLE PRECISION      TSCALE
      DOUBLE PRECISION      INITJD
      DOUBLE PRECISION      INITFR

C$ Brief_I/O
C
C   Variable  I/O  Description
C     --------  ---  --------------------------------------------------
C     HANDLE     I   Handle of PCK file open for writing.
C     CLSSID     I   NAIF PCK frame class ID.
C     FRAME      I   Reference frame name.
C     FIRST      I   Start time of interval covered by segment.
C     LAST       I   End time of interval covered by segment.
C     SEGID      I   Segment identifier.
C     INTLEN     I   Length of time covered by logical record (days).
C     N          I   Number of logical records in segment.
C     POLYDG     I   Chebyshev polynomial degree.
C     CDATA      I   Array of Chebyshev coefficients and angles.
C     ASCALE     I   Angular scale of data.
C     TSCALE     I   Time scale of data.
C     INITJD     I   Integer part of begin time (TDB Julian date) of
C                    first record.
C     INITFR     I   Fractional part of begin time (TDB Julian date) of
C                    first record.
C     MAXDEG     P   Maximum allowed degree of Chebyshev expansions.
C     TOLSCL     P   Tolerance scale factor for coverage bound checking.
C
C$ Detailed_Input
C
C     HANDLE         is the DAF handle of a PCK file to which a type 20
C                    segment is to be added. The PCK file must be open
C                    for writing.
C
C     CLSSID         is the integer NAIF PCK frame class ID code of the
C                    reference frame whose orientation relative to its
C                    base frame is described by the segment to be
C                    created. See the Frames Required Reading for
C                    details.
C
C     FRAME          is the NAIF name for a reference frame relative to
C                    which the orientation information for CLSSID is
C                    specified. This frame is called the "base frame."
C
C     FIRST,
C     LAST           are the start and stop times of the time interval
C                    over which the segment defines the orientation of
C                    the reference frame identified by CLSSID.
C
C     SEGID          is a segment identifier. A PCK segment identifier
C                    may contain up to 40 characters.
C
C     INTLEN         is the length of time, in TDB Julian days, covered
C                    by each set of Chebyshev polynomial coefficients
C                    (each logical record).
C
C     N              is the number of logical records to be stored in
C                    the segment. There is one logical record for each
C                    time period. Each logical record contains three
C                    sets of Chebyshev coefficients---one for each
C                    coordinate---and three position vector components.
C
C     POLYDG         is the degree of each set of Chebyshev
C                    polynomials, i.e. the number of Chebyshev
C                    coefficients per angle minus one. POLYDG must be
C                    less than or equal to the parameter MAXDEG.
C
C     CDATA          is an array containing sets of Chebyshev
C                    polynomial coefficients and angles to be placed in
C                    the new segment of the PCK file. The Chebyshev
C                    coefficients represent Euler angle rates; the
C                    angles are values of the Euler angles at each
C                    interval midpoint. The angular and time units of
C                    the data are defined by the inputs ASCALE and
C                    TSCALE, which are described below.
C
C                    The Euler angles represent the orientation of the
C                    reference frame designated by CLSSID relative to
C                    its base frame. The angles, which are numbered
C                    according to their ordinal position in the logical
C                    records, define a transformation matrix R as
C                    follows:
C
C                       R = [ A*ANGLE_3 ]  [ A*ANGLE_2 ]  [ A*ANGLE_1 ]
C                                        3              1              3
C
C                    where A is the angular scale ASCALE. Here the
C                    notation
C
C                       [ THETA ]
C                                i
C
C                    denotes a reference frame rotation of THETA
C                    radians in the right-hand sense about the ith
C                    coordinate axis. See the Rotation Required Reading
C                    for further discussion of this notation.
C
C                    The matrix R transforms vectors expressed in the
C                    base frame to vectors expressed in the frame
C                    associated with CLSSID by left multiplication:
C
C                       V       = R * V
C                        CLSSID        FRAME
C
C                    In cases where the frame designated by CLSSID
C                    (which we'll abbreviate as "the CLSSID frame") is
C                    a body-fixed, right-handed frame with its +Z axis
C                    aligned with a body's north pole, the orientation
C                    angles are related to right ascension (RA) and
C                    declination (DEC) of the CLSSID frame's north
C                    pole, and prime meridian orientation (W), by the
C                    equations
C
C                       ANGLE_1 * ASCALE = RA   + pi/2 radians
C                       ANGLE_2 * ASCALE = pi/2 - DEC  radians
C                       ANGLE_3 * ASCALE = W           radians
C
C                    The coefficients and angles are stored in CDATA in
C                    order as follows:
C
C                       the (POLYDG + 1) coefficients for the rate of
C                       the first angle of the first logical record,
C                       followed by the value of the first angle at the
C                       first interval midpoint.
C
C                       the coefficients for the rate of the second
C                       angle of the first logical record, followed by
C                       the value of the second angle at the first
C                       interval midpoint.
C
C                       the coefficients for the rate of the third
C                       angle of the first logical record, followed by
C                       the value of the third angle at the first
C                       interval midpoint.
C
C                       the (degree + 1) coefficients for the rate of
C                       the first angle of the second logical record,
C                       followed by the value of the first angle at the
C                       second interval midpoint.
C
C                       and so on.
C
C                    The logical data records are stored contiguously:
C
C                       +----------+
C                       | Record 1 |
C                       +----------+
C                       | Record 2 |
C                       +----------+
C                           ...
C                       +----------+
C                       | Record N |
C                       +----------+
C
C                    The contents of an individual record are:
C
C                       +--------------------------------------+
C                       | Coeff set for ANGLE_1 rate           |
C                       +--------------------------------------+
C                       | ANGLE_1                              |
C                       +--------------------------------------+
C                       | Coeff set for ANGLE_2 rate           |
C                       +--------------------------------------+
C                       | ANGLE_2                              |
C                       +--------------------------------------+
C                       | Coeff set for ANGLE_3 rate           |
C                       +--------------------------------------+
C                       | ANGLE_3                              |
C                       +--------------------------------------+
C
C                   Each coefficient set has the structure:
C
C                       +--------------------------------------+
C                       | Coefficient of T_0                   |
C                       +--------------------------------------+
C                       | Coefficient of T_1                   |
C                       +--------------------------------------+
C                                         ...
C                       +--------------------------------------+
C                       | Coefficient of T_POLYDG              |
C                       +--------------------------------------+
C
C                    Where T_n represents the Chebyshev polynomial
C                    of the first kind of degree n.
C
C     ASCALE,
C     TSCALE         are, respectively, the angular scale of the input
C                    angle and angular rate data in radians, and the
C                    time scale of the input rate data in TDB
C                    seconds.
C
C                    For example, if the input angular data have units
C                    of degrees, ASCALE should be set to the number of
C                    radians in one degree. If the input rate data have
C                    time units of Julian days, then TSCALE should be
C                    set to the number of seconds per Julian day
C                    (86400).
C
C
C     INITJD         is the integer part of the Julian ephemeris date
C                    of initial epoch of the first record. INITJD may
C                    be less than, equal to, or greater than the
C                    initial epoch.
C
C     INITFR         is the fractional part of the Julian ephemeris
C                    date of initial epoch of the first record. INITFR
C                    has units of Julian days. INITFR has magnitude
C                    strictly less than 1 day. The sum
C
C                       INITJD + INITFR
C
C                    equals the Julian ephemeris date of the initial
C                    epoch of the first record.
C
C
C$ Detailed_Output
C
C     None. This routine writes data to a PCK file.
C
C$ Parameters
C
C     MAXDEG         is the maximum allowed degree of the input
C                    Chebyshev expansions. MAXDEG is declared in the
C                    Fortran INCLUDE file pck20.inc.
C
C     TOLSCL         is a tolerance scale factor (also called a
C                    "relative tolerance") used for time coverage
C                    bound checking. TOLSCL is unitless. TOLSCL
C                    produces a tolerance value via the formula
C
C                       TOL = TOLSCL * MAX( ABS(FIRST), ABS(LAST) )
C
C                    where FIRST and LAST are the coverage time bounds
C                    of a type 20 segment, expressed as seconds past
C                    J2000 TDB.
C
C                    The resulting parameter TOL is used as a tolerance
C                    for comparing the input segment descriptor time
C                    bounds to the first and last epoch covered by the
C                    sequence of time intervals defined by the inputs
C                    to PCKW20:
C
C                       INITJD
C                       INITFR
C                       INTLEN
C                       N
C
C                    TOLSCL is declared in the Fortran INCLUDE file
C                    pck20.inc.
C
C                    See the Exceptions section below for a description
C                    of the error check using this tolerance.
C
C$ Exceptions
C
C     1)  If the number of sets of coefficients is not positive
C         SPICE(INVALIDCOUNT) is signaled.
C
C     2)  If the interval length is not positive, SPICE(INTLENNOTPOS)
C         is signaled.
C
C     3)  If the name of the reference frame is not recognized,
C         SPICE(INVALIDREFFRAME) is signaled.
C
C     4)  If segment stop time is not greater then the begin time,
C         SPICE(BADDESCRTIMES) is signaled.
C
C     5)  If the start time of the first record exceeds the descriptor
C         begin time by more than a computed tolerance, or if the end
C         time of the last record precedes the descriptor end time by
C         more than a computed tolerance, the error SPICE(COVERAGEGAP)
C         is signaled. See the Parameters section above for a
C         description of the tolerance.
C
C     6)  If the input degree POLYDG is less than 0 or greater than
C         MAXDEG, the error SPICE(INVALIDDEGREE) is signaled.
C
C     7)  If the last non-blank character of SEGID occurs past index
C         40, or if SEGID contains any nonprintable characters, the
C         error will be diagnosed by a routine in the call tree of this
C         routine.
C
C     8)  If either the angle or time scale is non-positive, the
C         error SPICE(NONPOSITIVESCALE) will be signaled.
C
C$ Files
C
C     A new type 20 PCK segment is written to the PCK file attached
C     to HANDLE.
C
C$ Particulars
C
C     This routine writes a PCK type 20 data segment to the designated
C     PCK file, according to the format described in the PCK Required
C     Reading.
C
C     Each segment can contain data for only one reference frame
C     and base frame. The Chebyshev polynomial degree and length
C     of time covered by each logical record are also fixed. However,
C     an arbitrary number of logical records of Chebyshev polynomial
C     coefficients can be written in each segment.  Minimizing the
C     number of segments in a PCK file will help optimize how the
C     SPICE system accesses the file.
C
C$ Examples
C
C     Suppose that you have in an array CDATA sets of Chebyshev
C     polynomial coefficients and angles representing the orientation
C     of the moon, relative to the J2000 reference frame, and you want
C     to put these into a type 20 segment in an existing PCK file. The
C     following code could be used to add one new type 20 segment. To
C     add multiple segments, put the call to PCKW20 in a loop.
C
C     C
C     C      First open the PCK file and get a handle for it.
C     C
C            CALL DAFOPW ( PCKNAM, HANDLE )
C
C     C
C     C      Create a segment identifier.
C     C
C            SEGID = 'MY_SAMPLE_PCK_TYPE_20_SEGMENT'
C
C     C
C     C      Note that the interval length INTLEN has units
C     C      of Julian days. The start time of the first record
C     C      is expressed using two inputs: integer and fractional
C     C      portions of the Julian ephemeris date of the start
C     C      time.
C     C
C     C      The PCK frame class ID code is stored in the
C     C      variable CLSSID. This ID must be associated in
C     C      with a PCK frame; usually such an association is
C     C      made via a frame kernel.
C     C
C     C      Write the segment.
C     C
C            CALL PCKW20 ( HANDLE, CLSSID, 'J2000', FIRST,
C          .               LAST,   SEGID,  INTLEN,  N,
C          .               POLYDG, CDATA,  ASCALE,  TSCALE
C          .               INITJD, INITFR                  )
C
C     C
C     C      Close the file.
C     C
C            CALL DAFCLS ( HANDLE )
C
C$ Restrictions
C
C     None.
C
C$ Literature_References
C
C     None.
C
C$ Author_and_Institution
C
C     N.J. Bachman (JPL)
C     K.S. Zukor   (JPL)
C
C$ Version
C
C-    SPICELIB Version 1.0.0, 17-JAN-2014 (NJB) (KSZ)
C
C-&

C$ Index_Entries
C
C     write pck type_20 data segment
C
C-&

C
C     SPICELIB functions
C
      DOUBLE PRECISION      J2000
      DOUBLE PRECISION      SPD

      LOGICAL               FAILED
      LOGICAL               RETURN

C
C     Local Parameters
C

C
C     DTYPE is the PCK data type.
C
      INTEGER               DTYPE
      PARAMETER           ( DTYPE   =   20 )

C
C     ND is the number of double precision components in a PCK
C     segment descriptor. PCK uses ND = 2.
C
      INTEGER               ND
      PARAMETER           ( ND      =   2 )

C
C     NI is the number of integer components in a PCK segment
C     descriptor. PCK uses NI = 5.
C
      INTEGER               NI
      PARAMETER           ( NI      =   5 )

C
C     NS is the size of a packed PCK segment descriptor.
C
      INTEGER               NS
      PARAMETER           ( NS      =   5 )

C
C     SIDLEN is the maximum number of characters allowed in an
C     PCK segment identifier.
C
      INTEGER               SIDLEN
      PARAMETER           ( SIDLEN  =  40 )

C
C     Local variables
C
      CHARACTER*(SIDLEN)    ETSTR
      CHARACTER*(SIDLEN)    NETSTR

      DOUBLE PRECISION      BTIME
      DOUBLE PRECISION      DCD   ( ND )
      DOUBLE PRECISION      DESCR ( NS )
      DOUBLE PRECISION      LTIME
      DOUBLE PRECISION      NUMREC
      DOUBLE PRECISION      TOL

      INTEGER               ICD   ( NI )
      INTEGER               NINREC
      INTEGER               REFCOD

C
C
C     Standard SPICE error handling.
C
      IF ( RETURN () ) THEN
         RETURN
      END IF

      CALL CHKIN ( 'PCKW20' )

C
C     The number of sets of coefficients must be positive.
C
      IF ( N .LE. 0 ) THEN

         CALL SETMSG ( 'The number of sets of coordinate '
     .   //            'coefficients is not positive. N = #.' )
         CALL ERRINT ( '#', N                                 )
         CALL SIGERR ( 'SPICE(INVALIDCOUNT)'                  )
         CALL CHKOUT ( 'PCKW20'                               )
         RETURN

      END IF

C
C     Make sure that the degree of the interpolating polynomials is
C     in range.
C
      IF (  ( POLYDG .LT. 0 ) .OR. ( POLYDG .GT. MAXDEG )  )  THEN

         CALL SETMSG ( 'The interpolating polynomials have degree #; '
     .   //            'the valid degree range is [0, #].'             )
         CALL ERRINT ( '#', POLYDG                                     )
         CALL ERRINT ( '#', MAXDEG                                     )
         CALL SIGERR ( 'SPICE(INVALIDDEGREE)'                          )
         CALL CHKOUT ( 'PCKW20'                                        )
         RETURN

      END IF

C
C     The interval length must be positive.
C
      IF ( INTLEN .LE. 0 ) THEN

         CALL SETMSG ( 'The interval length is not positive.'
     .   //            'N = #'                                )
         CALL ERRDP  ( '#', INTLEN                            )
         CALL SIGERR ( 'SPICE(INTLENNOTPOS)'                  )
         CALL CHKOUT ( 'PCKW20'                               )
         RETURN

      END IF

C
C     Get the NAIF integer code for the reference frame.
C
      CALL NAMFRM ( FRAME, REFCOD )

      IF ( REFCOD .EQ. 0 ) THEN

         CALL SETMSG ( 'The reference frame # is not supported.' )
         CALL ERRCH  ( '#', FRAME                                )
         CALL SIGERR ( 'SPICE(INVALIDREFFRAME)'                  )
         CALL CHKOUT ( 'PCKW20'                                  )
         RETURN

      END IF

C
C     The segment stop time must be greater than the begin time.
C
      IF ( FIRST .GE. LAST ) THEN

         CALL SETMSG ( 'The segment start time: # (# TDB) is '
     .   //            'not less than the segment end time: '
     .   //            '(# TDB).'                               )
         CALL ETCAL  ( FIRST, ETSTR                             )
         CALL ERRCH  ( '#',   ETSTR                             )
         CALL ERRDP  ( '#',   FIRST                             )
         CALL ETCAL  ( LAST,  NETSTR                            )
         CALL ERRCH  ( '#',   NETSTR                            )
         CALL ERRDP  ( '#',   LAST                              )
         CALL SIGERR ( 'SPICE(BADDESCRTIMES)'                   )
         CALL CHKOUT ( 'PCKW20'                                 )
         RETURN

      END IF

C
C     The angle and time scales must be positive.
C
      IF ( ASCALE .LE. 0.D0 ) THEN

         CALL SETMSG ( 'The angle scale is not positive.'  //
     .                 'ASCALE = #'                           )
         CALL ERRDP  ( '#', ASCALE                            )
         CALL SIGERR ( 'SPICE(NONPOSITIVESCALE)'              )
         CALL CHKOUT ( 'PCKW20'                               )
         RETURN

      END IF

      IF ( TSCALE .LE. 0.D0 ) THEN

         CALL SETMSG ( 'The time scale is not positive.'     //
     .                 'TSCALE = #'                           )
         CALL ERRDP  ( '#', TSCALE                            )
         CALL SIGERR ( 'SPICE(NONPOSITIVESCALE)'              )
         CALL CHKOUT ( 'PCKW20'                               )
         RETURN

      END IF

C
C     The begin time of the first record must be less than or equal
C     to the begin time of the segment. Convert the two-part input
C     epoch to seconds past J2000 for the purpose of this check.
C
      BTIME = SPD()  *  (  ( INITJD - J2000() ) + INITFR  )

      LTIME = BTIME  +  ( N * INTLEN * SPD() )

C
C     Compute the tolerance to use for descriptor time bound checks.
C
      TOL = TOLSCL * MAX( ABS(BTIME), ABS(LTIME) )

      IF ( FIRST .LT. BTIME - TOL ) THEN

         CALL SETMSG ( 'The segment descriptor start time # is too '
     .   //            'much less than the beginning time of the '
     .   //            'segment data # (in seconds past J2000: #). '
     .   //            'The difference is # seconds; the '
     .   //            'tolerance is # seconds.'                     )
         CALL ETCAL  ( FIRST, ETSTR                                  )
         CALL ERRCH  ( '#',   ETSTR                                  )
         CALL ETCAL  ( BTIME, ETSTR                                  )
         CALL ERRCH  ( '#',   ETSTR                                  )
         CALL ERRDP  ( '#',   FIRST                                  )
         CALL ERRDP  ( '#',   BTIME-FIRST                            )
         CALL ERRDP  ( '#',   TOL                                    )
         CALL SIGERR ( 'SPICE(COVERAGEGAP)'                          )
         CALL CHKOUT ( 'PCKW20'                                      )
         RETURN

      END IF

C
C     The end time of the final record must be greater than or
C     equal to the end time of the segment.
C
      IF ( LAST .GT. LTIME + TOL ) THEN

         CALL SETMSG ( 'The segment descriptor end time # is too '
     .   //            'much greater than the end time of the segment '
     .   //            'data # (in seconds past J2000: #). The '
     .   //            'difference is # seconds; the tolerance is # '
     .   //            'seconds.'                                    )
         CALL ETCAL  ( LAST,  ETSTR                                  )
         CALL ERRCH  ( '#',   ETSTR                                  )
         CALL ETCAL  ( LTIME, ETSTR                                  )
         CALL ERRCH  ( '#',   ETSTR                                  )
         CALL ERRDP  ( '#',   LAST                                   )
         CALL ERRDP  ( '#',   LAST-LTIME                             )
         CALL ERRDP  ( '#',   TOL                                    )
         CALL SIGERR ( 'SPICE(COVERAGEGAP)'                          )
         CALL CHKOUT ( 'PCKW20'                                      )
         RETURN

      END IF

C
C     Now check the validity of the segment identifier.
C
      CALL CHCKID ( 'PCK segment identifier', SIDLEN, SEGID )

      IF ( FAILED() ) THEN
         CALL CHKOUT ( 'PCKW20' )
         RETURN
      END IF

C
C     Store the start and end times to be associated
C     with this segment.
C
      DCD(1) = FIRST
      DCD(2) = LAST

C
C     Create the integer portion of the descriptor.
C
      ICD(1) = CLSSID
      ICD(2) = REFCOD
      ICD(3) = DTYPE

C
C     Pack the segment descriptor.
C
      CALL DAFPS ( ND, NI, DCD, ICD, DESCR )

C
C     Begin a new segment of PCK type 20 form:
C
C        Record 1
C        Record 2
C        ...
C        Record N
C        ASCALE     ( angular scale in radians )
C        TSCALE     ( time scale in seconds )
C        INITJD     ( integer part of initial epoch of first record,
C                     expressed as a TDB Julian date )
C        INITFR     ( fractional part of initial epoch, in units of
C                     TDB Julian days )
C        INTLEN     ( length of interval covered by each record, in
C                     units of TDB Julian days )
C        RSIZE      ( number of data elements in each record )
C        N          ( number of records in segment )
C
C     Each record will have the form:
C
C        ANGLE_1 coefficients
C        ANGLE_1 angle at interval midpoint
C        ANGLE_2 coefficients
C        ANGLE_2 angle at interval midpoint
C        ANGLE_3 coefficients
C        ANGLE_3 angle at interval midpoint
C
C
      CALL DAFBNA ( HANDLE, DESCR, SEGID )

C
C     Calculate the number of entries in a record.
C
      NINREC = ( POLYDG + 2 ) * 3

C
C     Fill segment with N records of data.
C
      CALL DAFADA ( CDATA, N * NINREC )

C
C     Store the angle and time scales.
C
      CALL DAFADA ( ASCALE, 1 )
      CALL DAFADA ( TSCALE, 1 )

C
C     Store the integer and fractional parts of the initial epoch of
C     the first record.
C
      CALL DAFADA ( INITJD, 1 )
      CALL DAFADA ( INITFR, 1 )

C
C     Store the length of interval covered by each record.
C
      CALL DAFADA ( INTLEN, 1 )

C
C     Store the size of each record (total number of array elements).
C     Note that this size is smaller by 2 than the size of a type 2
C     record of the same degree, since the record coverage midpoint
C     and radius are not stored.
C
      CALL DAFADA ( DBLE( NINREC ), 1 )

C
C     Store the number of records contained in the segment.
C
      NUMREC = N
      CALL DAFADA ( NUMREC, 1 )

C
C     End this segment.
C
      CALL DAFENA

      CALL CHKOUT ( 'PCKW20' )
      RETURN
      END

max_stars_repo_name: agforero/FTFramework

max_issues_repo_name: agforero/fortran-testing-framework

max_forks_repo_name: agforero/fortran-testing-framework


text:
c     Subroutine to find basic stats for month values
c     AJ_Kettle, Nov30/2017

      SUBROUTINE calc_basic_stats_mon(f_ndflag,l_monrec_common,
     +  f_monrec_airt_k,f_monrec_wspd_ms,
     +  f_monrec_slpres_hpa,f_monrec_ppt_mm,

     +  f_monavg_airt_ngood,f_monavg_airt_nbad,
     +  f_monavg_airt_min_k,f_monavg_airt_max_k,
     +  f_monavg_airt_avg_k,
     +  f_monavg_wspd_ngood,f_monavg_wspd_nbad,
     +  f_monavg_wspd_min_ms,f_monavg_wspd_max_ms,
     +  f_monavg_wspd_avg_ms,
     +  f_monavg_slpres_ngood,f_monavg_slpres_nbad,
     +  f_monavg_slpres_min_hpa,f_monavg_slpres_max_hpa,
     +  f_monavg_slpres_avg_hpa,
     +  f_monavg_ppt_ngood,f_monavg_ppt_nbad,
     +  f_monavg_ppt_min_mm,f_monavg_ppt_max_mm,
     +  f_monavg_ppt_avg_mm)

      IMPLICIT NONE
c************************************************************************
      REAL                :: f_ndflag
      INTEGER             :: l_monrec_common

      REAL                :: f_monrec_airt_k(5000)
      REAL                :: f_monrec_wspd_ms(5000)
      REAL                :: f_monrec_slpres_hpa(5000)
      REAL                :: f_monrec_ppt_mm(5000)

      REAL                :: f_monavg_airt_ngood
      REAL                :: f_monavg_airt_nbad
      REAL                :: f_monavg_airt_min_k
      REAL                :: f_monavg_airt_max_k
      REAL                :: f_monavg_airt_avg_k

      REAL                :: f_monavg_wspd_ngood
      REAL                :: f_monavg_wspd_nbad
      REAL                :: f_monavg_wspd_min_ms
      REAL                :: f_monavg_wspd_max_ms
      REAL                :: f_monavg_wspd_avg_ms

      REAL                :: f_monavg_slpres_ngood
      REAL                :: f_monavg_slpres_nbad
      REAL                :: f_monavg_slpres_min_hpa
      REAL                :: f_monavg_slpres_max_hpa
      REAL                :: f_monavg_slpres_avg_hpa

      REAL                :: f_monavg_ppt_ngood
      REAL                :: f_monavg_ppt_nbad
      REAL                :: f_monavg_ppt_min_mm
      REAL                :: f_monavg_ppt_max_mm
      REAL                :: f_monavg_ppt_avg_mm

      INTEGER             :: ii,jj,kk,i,j,k
c************************************************************************
c     Find stats dayavg_airt
      CALL dwd_mon_stat_package(l_monrec_common,f_ndflag,
     + f_monrec_airt_k,
     + f_monavg_airt_min_k,f_monavg_airt_max_k,
     + f_monavg_airt_avg_k,
     + f_monavg_airt_ngood,f_monavg_airt_nbad)

c     Find stats dayavg_wspd
      CALL dwd_mon_stat_package(l_monrec_common,f_ndflag,
     + f_monrec_wspd_ms,
     + f_monavg_wspd_min_ms,f_monavg_wspd_max_ms,
     + f_monavg_wspd_avg_ms,
     + f_monavg_wspd_ngood,f_monavg_wspd_nbad)

c     Find stats dayavg_slpres
      CALL dwd_mon_stat_package(l_monrec_common,f_ndflag,
     + f_monrec_slpres_hpa,
     + f_monavg_slpres_min_hpa,f_monavg_slpres_max_hpa,
     + f_monavg_slpres_avg_hpa,
     + f_monavg_slpres_ngood,f_monavg_slpres_nbad)

c     Find stats dayavg_ppt
      CALL dwd_mon_stat_package(l_monrec_common,f_ndflag,
     + f_monrec_ppt_mm,
     + f_monavg_ppt_min_mm,f_monavg_ppt_max_mm,
     + f_monavg_ppt_avg_mm,
     + f_monavg_ppt_ngood,f_monavg_ppt_nbad)
c************************************************************************
c      print*,'airt_k=',f_monavg_airt_min_k,f_monavg_airt_max_k,
c     + f_monavg_airt_avg_k,
c     + f_monavg_airt_ngood,f_monavg_airt_nbad
c      print*,'wspd_ms=',f_monavg_wspd_min_ms,f_monavg_wspd_max_ms,
c     + f_monavg_wspd_avg_ms,
c     + f_monavg_wspd_ngood,f_monavg_wspd_nbad
c      print*,'slpres_hpa=',
c     + f_monavg_slpres_min_hpa,f_monavg_slpres_max_hpa,
c     + f_monavg_slpres_avg_hpa,
c     + f_monavg_slpres_ngood,f_monavg_slpres_nbad
c      print*,'ppt_mm=',f_monavg_ppt_min_mm,f_monavg_ppt_max_mm,
c     + f_monavg_ppt_avg_mm,
c     + f_monavg_ppt_ngood,f_monavg_ppt_nbad

      RETURN
      END

ext: f

max_stars_repo_name: ajkettle/glamod-nuim

max_issues_repo_name: ajkettle/glamod-nuim

max_forks_repo_name: ajkettle/glamod-nuim


text:
#redirect Habitat for Humanity

Welcome to the wiki Its obvious you are with Habitat for Humanity and the page for that already exists!
Please consider the importance of using your RealName   Daubert

max_stars_repo_name: voflo/Search

max_issues_repo_name: voflo/Search

max_forks_repo_name: voflo/Search


text:
C$Procedure      EKNSEG ( EK, number of segments in file )

      INTEGER FUNCTION EKNSEG ( HANDLE )

C$ Abstract
C
C     Return the number of segments in a specified EK.
C
C$ Disclaimer
C
C     THIS SOFTWARE AND ANY RELATED MATERIALS WERE CREATED BY THE
C     CALIFORNIA INSTITUTE OF TECHNOLOGY (CALTECH) UNDER A U.S.
C     GOVERNMENT CONTRACT WITH THE NATIONAL AERONAUTICS AND SPACE
C     ADMINISTRATION (NASA). THE SOFTWARE IS TECHNOLOGY AND SOFTWARE
C     PUBLICLY AVAILABLE UNDER U.S. EXPORT LAWS AND IS PROVIDED "AS-IS"
C     TO THE RECIPIENT WITHOUT WARRANTY OF ANY KIND, INCLUDING ANY
C     WARRANTIES OF PERFORMANCE OR MERCHANTABILITY OR FITNESS FOR A
C     PARTICULAR USE OR PURPOSE (AS SET FORTH IN UNITED STATES UCC
C     SECTIONS 2312-2313) OR FOR ANY PURPOSE WHATSOEVER, FOR THE
C     SOFTWARE AND RELATED MATERIALS, HOWEVER USED.
C
C     IN NO EVENT SHALL CALTECH, ITS JET PROPULSION LABORATORY, OR NASA
C     BE LIABLE FOR ANY DAMAGES AND/OR COSTS, INCLUDING, BUT NOT
C     LIMITED TO, INCIDENTAL OR CONSEQUENTIAL DAMAGES OF ANY KIND,
C     INCLUDING ECONOMIC DAMAGE OR INJURY TO PROPERTY AND LOST PROFITS,
C     REGARDLESS OF WHETHER CALTECH, JPL, OR NASA BE ADVISED, HAVE
C     REASON TO KNOW, OR, IN FACT, SHALL KNOW OF THE POSSIBILITY.
C
C     RECIPIENT BEARS ALL RISK RELATING TO QUALITY AND PERFORMANCE OF
C     THE SOFTWARE AND ANY RELATED MATERIALS, AND AGREES TO INDEMNIFY
C     CALTECH AND NASA FOR ALL THIRD-PARTY CLAIMS RESULTING FROM THE
C     ACTIONS OF RECIPIENT IN THE USE OF THE SOFTWARE.
C
C$ Required_Reading
C
C     EK
C
C$ Keywords
C
C     EK
C     UTILITY
C
C$ Declarations

      INCLUDE 'ekfilpar.inc'

      INTEGER               HANDLE

C$ Brief_I/O
C
C     Variable  I/O  Description
C     --------  ---  --------------------------------------------------
C     HANDLE     I   EK file handle.
C
C     The function returns the number of segments in the specified
C     E-kernel.
C
C$ Detailed_Input
C
C     HANDLE         is the handle of an EK file opened for read
C                    access.
C
C$ Detailed_Output
C
C     The function returns the number of segments in the specified
C     E-kernel.
C
C$ Parameters
C
C     None.
C
C$ Exceptions
C
C     1)  If HANDLE is invalid, the error will be diagnosed by routines
C         called by this routine.  EKNSEG will return the value zero.
C
C     2)  If an I/O error occurs while trying to read the EK, the error
C         will be diagnosed by routines called by this routine.
C         EKNSEG will return the value zero.
C
C$ Files
C
C     See the description of HANDLE in $Detailed_Input.
C
C$ Particulars
C
C     This routine is used to support the function of summarizing an
C     EK file.  Given the number of segments in the file, a program
C     can use EKSSUM in a loop to summarize each of them.
C
C$ Examples
C
C     1)  Open an EK file and count the segments in it.
C
C            CALL EKOPR ( EKNAME, HANDLE )
C            N = EKNSEG  ( HANDLE )
C
C$ Restrictions
C
C     None.
C
C$ Literature_References
C
C     None.
C
C$ Author_and_Institution
C
C     N.J. Bachman   (JPL)
C
C$ Version
C
C-    Beta Version 1.0.0, 26-SEP-1995 (NJB)
C
C-&

C$ Index_Entries
C
C     return number of segments in an E-kernel
C
C-&

C
C     SPICELIB functions
C
      LOGICAL               FAILED
      LOGICAL               RETURN

C
C     Non-SPICELIB functions
C
      INTEGER               ZZEKTRBS
      INTEGER               ZZEKTRSZ

C
C     Local variables
C
      INTEGER               BASE
      INTEGER               TREE

C
C     Set a default value for EKNSEG.
C
      EKNSEG  =  0

C
C     Standard SPICE error handling.
C
      IF ( RETURN () ) THEN
         RETURN
      ELSE
         CALL CHKIN ( 'EKNSEG' )
      END IF

C
C     Make sure this is a paged DAS EK.
C
      CALL ZZEKPGCH ( HANDLE, 'READ' )

      IF ( FAILED() ) THEN
         CALL CHKOUT ( 'EKNSEG' )
         RETURN
      END IF

C
C     Obtain the base address of the first integer page.
C
      BASE = ZZEKTRBS ( 1 )

C
C     Look up the head node of the segment tree.
C
      CALL DASRDI ( HANDLE, BASE+SGTIDX, BASE+SGTIDX, TREE )

C
C     Get the entry count for the segment tree.
C
      EKNSEG  =  ZZEKTRSZ ( HANDLE, TREE )

      CALL CHKOUT ( 'EKNSEG' )
      RETURN
      END


text:
! This automatically generated Fortran wrapper file allows codes
! written in Fortran to be called directly from C and translates all
! C-style arguments into expected Fortran-style arguments (with
! assumed size, local type declarations, etc.).

MODULE C_SWAP
USE ISO_FORTRAN_ENV , ONLY : REAL64 , INT64
  IMPLICIT NONE

CONTAINS

  SUBROUTINE C_SWAP_I64(V1, V2) BIND(C)
    USE SWAP, ONLY: SWAP_I64
    IMPLICIT NONE
    INTEGER(KIND=INT64), INTENT(INOUT) :: V1
    INTEGER(KIND=INT64), INTENT(INOUT) :: V2

    CALL SWAP_I64(V1, V2)
  END SUBROUTINE C_SWAP_I64

  SUBROUTINE C_SWAP_R64(V1, V2) BIND(C)
    USE SWAP, ONLY: SWAP_R64
    IMPLICIT NONE
    REAL(KIND=REAL64), INTENT(INOUT) :: V1
    REAL(KIND=REAL64), INTENT(INOUT) :: V2

    CALL SWAP_R64(V1, V2)
  END SUBROUTINE C_SWAP_R64

END MODULE C_SWAP

MODULE C_FAST_SELECT
USE ISO_FORTRAN_ENV , ONLY : REAL64 , INT64
USE SWAP , ONLY : SWAP_I64 , SWAP_R64
  IMPLICIT NONE

CONTAINS

  SUBROUTINE C_ARGSELECT_R64(VALUES_DIM_1, VALUES, INDICES_DIM_1, INDICES, K, DIVISOR_PRESENT, DIVISOR, MAX_SIZE_PRESENT, MAX_SIZE)&
& BIND(C)
    USE FAST_SELECT, ONLY: ARGSELECT_R64
    IMPLICIT NONE
    INTEGER, INTENT(IN) :: VALUES_DIM_1
    REAL(KIND=REAL64), INTENT(INOUT), DIMENSION(VALUES_DIM_1) :: VALUES
    INTEGER, INTENT(IN) :: INDICES_DIM_1
    INTEGER(KIND=INT64), INTENT(INOUT), DIMENSION(INDICES_DIM_1) :: INDICES
    INTEGER(KIND=INT64), INTENT(IN) :: K
    LOGICAL, INTENT(IN) :: DIVISOR_PRESENT
    INTEGER(KIND=INT64), INTENT(IN) :: DIVISOR
    LOGICAL, INTENT(IN) :: MAX_SIZE_PRESENT
    INTEGER(KIND=INT64), INTENT(IN) :: MAX_SIZE

    IF (DIVISOR_PRESENT) THEN
      IF (MAX_SIZE_PRESENT) THEN
        CALL ARGSELECT_R64(VALUES=VALUES, INDICES=INDICES, K=K, DIVISOR=DIVISOR, MAX_SIZE=MAX_SIZE)
      ELSE
        CALL ARGSELECT_R64(VALUES=VALUES, INDICES=INDICES, K=K, DIVISOR=DIVISOR)
      END IF
    ELSE
      IF (MAX_SIZE_PRESENT) THEN
        CALL ARGSELECT_R64(VALUES=VALUES, INDICES=INDICES, K=K, MAX_SIZE=MAX_SIZE)
      ELSE
        CALL ARGSELECT_R64(VALUES=VALUES, INDICES=INDICES, K=K)
      END IF
    END IF
  END SUBROUTINE C_ARGSELECT_R64

  SUBROUTINE C_ARGSELECT_I64(VALUES_DIM_1, VALUES, INDICES_DIM_1, INDICES, K, DIVISOR_PRESENT, DIVISOR, MAX_SIZE_PRESENT, MAX_SIZE)&
& BIND(C)
    USE FAST_SELECT, ONLY: ARGSELECT_I64
    IMPLICIT NONE
    INTEGER, INTENT(IN) :: VALUES_DIM_1
    INTEGER(KIND=INT64), INTENT(INOUT), DIMENSION(VALUES_DIM_1) :: VALUES
    INTEGER, INTENT(IN) :: INDICES_DIM_1
    INTEGER(KIND=INT64), INTENT(INOUT), DIMENSION(INDICES_DIM_1) :: INDICES
    INTEGER(KIND=INT64), INTENT(IN) :: K
    LOGICAL, INTENT(IN) :: DIVISOR_PRESENT
    INTEGER(KIND=INT64), INTENT(IN) :: DIVISOR
    LOGICAL, INTENT(IN) :: MAX_SIZE_PRESENT
    INTEGER(KIND=INT64), INTENT(IN) :: MAX_SIZE

    IF (DIVISOR_PRESENT) THEN
      IF (MAX_SIZE_PRESENT) THEN
        CALL ARGSELECT_I64(VALUES=VALUES, INDICES=INDICES, K=K, DIVISOR=DIVISOR, MAX_SIZE=MAX_SIZE)
      ELSE
        CALL ARGSELECT_I64(VALUES=VALUES, INDICES=INDICES, K=K, DIVISOR=DIVISOR)
      END IF
    ELSE
      IF (MAX_SIZE_PRESENT) THEN
        CALL ARGSELECT_I64(VALUES=VALUES, INDICES=INDICES, K=K, MAX_SIZE=MAX_SIZE)
      ELSE
        CALL ARGSELECT_I64(VALUES=VALUES, INDICES=INDICES, K=K)
      END IF
    END IF
  END SUBROUTINE C_ARGSELECT_I64

END MODULE C_FAST_SELECT

MODULE C_FAST_SORT
USE ISO_FORTRAN_ENV , ONLY : REAL64 , INT64
USE SWAP , ONLY : SWAP_I64 , SWAP_R64
  IMPLICIT NONE

CONTAINS

  SUBROUTINE C_ARGSORT_R64(VALUES_DIM_1, VALUES, INDICES_DIM_1, INDICES, MIN_SIZE_PRESENT, MIN_SIZE) BIND(C)
    USE FAST_SORT, ONLY: ARGSORT_R64
    IMPLICIT NONE
    INTEGER, INTENT(IN) :: VALUES_DIM_1
    REAL(KIND=REAL64), INTENT(INOUT), DIMENSION(VALUES_DIM_1) :: VALUES
    INTEGER, INTENT(IN) :: INDICES_DIM_1
    INTEGER(KIND=INT64), INTENT(INOUT), DIMENSION(INDICES_DIM_1) :: INDICES
    LOGICAL, INTENT(IN) :: MIN_SIZE_PRESENT
    INTEGER(KIND=INT64), INTENT(IN) :: MIN_SIZE

    IF (MIN_SIZE_PRESENT) THEN
      CALL ARGSORT_R64(VALUES=VALUES, INDICES=INDICES, MIN_SIZE=MIN_SIZE)
    ELSE
      CALL ARGSORT_R64(VALUES=VALUES, INDICES=INDICES)
    END IF
  END SUBROUTINE C_ARGSORT_R64

  SUBROUTINE C_INSERTION_ARGSORT_R64(VALUES_DIM_1, VALUES, INDICES_DIM_1, INDICES) BIND(C)
    USE FAST_SORT, ONLY: INSERTION_ARGSORT_R64
    IMPLICIT NONE
    INTEGER, INTENT(IN) :: VALUES_DIM_1
    REAL(KIND=REAL64), INTENT(INOUT), DIMENSION(VALUES_DIM_1) :: VALUES
    INTEGER, INTENT(IN) :: INDICES_DIM_1
    INTEGER(KIND=INT64), INTENT(INOUT), DIMENSION(INDICES_DIM_1) :: INDICES

    CALL INSERTION_ARGSORT_R64(VALUES, INDICES)
  END SUBROUTINE C_INSERTION_ARGSORT_R64

  SUBROUTINE C_ARGPARTITION_R64(VALUES_DIM_1, VALUES, INDICES_DIM_1, INDICES, LEFT) BIND(C)
    USE FAST_SORT, ONLY: ARGPARTITION_R64
    IMPLICIT NONE
    INTEGER, INTENT(IN) :: VALUES_DIM_1
    REAL(KIND=REAL64), INTENT(INOUT), DIMENSION(VALUES_DIM_1) :: VALUES
    INTEGER, INTENT(IN) :: INDICES_DIM_1
    INTEGER(KIND=INT64), INTENT(INOUT), DIMENSION(INDICES_DIM_1) :: INDICES
    INTEGER(KIND=INT64) :: LEFT

    LEFT = ARGPARTITION_R64(VALUES, INDICES)
  END SUBROUTINE C_ARGPARTITION_R64

END MODULE C_FAST_SORT

max_stars_repo_name: thchang/SortVectors

max_issues_repo_name: thchang/SortVectors

max_forks_repo_name: thchang/SortVectors


text:
program where_test
    use, intrinsic :: iso_fortran_env, only : sp => REAL32, error_unit
    implicit none
    integer, parameter :: N = 5
    real, parameter :: eps = 1.0e-2_sp
    real(kind=sp), dimension(N, N) :: A
    integer, dimension(N, N) :: Pos, Neg, Zero
    real(kind=sp) :: r
    integer :: i, j
    character(len=40) :: r_fmt_str
    character(len=40) :: i_fmt_str

    write (r_fmt_str, "('(', I2, 'F7.2)')") N
    write (i_fmt_str, "('(', I2, 'I7)')") N

    do i = 1, N
        do j = 1, N
            call random_number(r)
            A(i, j) = 2.0_sp*(r - 0.5_sp)
        end do
    end do

    where (A > eps)
        Pos = 1
        Neg = 0
        Zero = 0
    elsewhere (A < eps)
        Pos = 0
        Neg = 1
        Zero = 0
    elsewhere
        Pos = 0
        Neg = 0
        Zero = 1
    end where

    print '(A)', 'matrix:'
    do i = 1, N
        print r_fmt_str, A(i, :)
    end do
    print '(2A)', '---------', 'positive'
    do i = 1, N
        print i_fmt_str, Pos(i, :)
    end do
    print '(2A)', '---------', 'negative'
    do i = 1, N
        print i_fmt_str, Neg(i, :)
    end do
    print '(2A)', '---------', 'zero'
    do i = 1, N
        print i_fmt_str, Zero(i, :)
    end do

    if (sum(Pos + Neg + Zero) /= N*N) then
        write (unit=error_unit, fmt='(A)') '# error: not all covered'
        stop 1
    end if

end program where_test

max_stars_repo_name: Gjacquenot/training-material

max_issues_repo_name: Gjacquenot/training-material

max_forks_repo_name: Gjacquenot/training-material


text:
! wind_realtime_lib.for - routines for accessing wind/waves realtime tm stream

!------------------------------------------------------------------------------
! Sets up realtime-specific internal data processing elements.
!------------------------------------------------------------------------------
	integer*4	function	w_setup_rt_stream(ch) !,f)
	implicit	none
	include		'wind_tm_user_def.for'
	include		'wind_realtime_def.for'
	integer*4	ch
!	character*(*)	f		! argument not used
	integer*4	ok
	integer*4	i
	integer*4	ios
	external	w_rt_close_stream
	external	w_rt_get_rec
	external	w_rt_frame_rate
	external	w_rt_scet_of_mfmf
	external	w_rt_scet_dbms_of_mfmf
	external	w_rt_get_packet
	external	w_rt_get_xtra_event_info
	external	w_rt_get_rec_by_time
	external	w_rt_get_word
	external	w_rt_get_mode
	external	w_rt_get_hk
	external	w_rt_get_dpu_mf
	external	w_rt_get_mfmf
	integer*4	wind_open_realtime
	logical*4	wind_suppress_messages

	w_setup_rt_stream = 0

	ok = wind_open_realtime(ch)
	if (ok .ne. 1) goto 10

	! these values must be non-zero: use dummy functions if necessary
	user(ch).f_get_word            = %loc(w_rt_get_word)
	user(ch).f_get_tm_mode         = %loc(w_rt_get_mode)
	user(ch).f_get_intrnl_rec      = %loc(w_rt_get_rec)
	user(ch).f_get_plain_packet    = %loc(w_rt_get_packet)
	user(ch).f_goto_time           = %loc(w_rt_get_rec_by_time)
	user(ch).f_get_mf_scet         = %loc(w_rt_scet_of_mfmf)
	user(ch).f_get_mf_scet_dbms    = %loc(w_rt_scet_dbms_of_mfmf)
	user(ch).f_get_frame_rate      = %loc(w_rt_frame_rate)
	user(ch).f_get_xtra_event_info = %loc(w_rt_get_xtra_event_info)
	user(ch).f_close_ch            = %loc(w_rt_close_stream)
	user(ch).f_get_hk              = %loc(w_rt_get_hk)
	user(ch).f_get_dpu_mf          = %loc(w_rt_get_dpu_mf)
	user(ch).f_get_rec_idx         = %loc(w_rt_get_mfmf)

	rltm(ch).wait_for_new_records	= .true.
	rltm(ch).ch = ch

	do i=0,sizeof_minor_frame-1
	   user(ch).is_a_valid_ptr(i) = rltm(ch).is_a_valid_ptr(i)
	end do

	w_setup_rt_stream = 1

	return
   1	format(1x,'W_SETUP_RT_STREAM: ', a)
  10	continue
	w_setup_rt_stream = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) ' cannot setup realtime TM stream.'
	return
	end

!------------------------------------------------------------------------------
! Handles the special cases of updating the various realtime stream positioning
! capabilities provided by wind_lib (e.g.: earliest_mfmf, latest_mfmf).
!------------------------------------------------------------------------------
	integer*4	function	w_rt_get_mfmf
	1				(ch,major,minor,position_type)
	implicit	none
	include		'parm_def.for'
	include		'glbsec_def.for'
	include		'wind_realtime_def.for'
	include		'low_byte_def.for'
	include		'wind_wnd_def.for'
	integer*4	ch
	integer*4	major
	integer*4	minor
	integer*4	position_type

	integer*4	ok
	integer*4	ios
	integer*4	i
	integer*4	w_wnd_close_stream
	integer*4	w_wnd_open_file
	logical*4	wind_suppress_messages

	w_rt_get_mfmf = 0

	if (position_type .eq. w_tk_stream_mfmf) then
	   major = glsec.last_major_validated
	   minor = glsec.last_minor_validated
	else if (position_type .eq. w_tk_latest_mfmf) then
	   major = glsec.last_major_validated
	   minor = glsec.last_minor_validated
	else if (position_type .eq. w_tk_earliest_mfmf) then
	   ! is there a recording file?
	   if (glsec.filename(1:1) .gt. ' ') then
	      if (rltm(ch).file .eq. glsec.filename) then
	         ! the offline file has already been opened and is current
	         major = wnd(ch).first_major
	         minor = wnd(ch).first_minor.i4val
	      else if (rltm(ch).file(1:1) .gt. ' ') then
	         ! the offline file is not the current recording file
	         ! so close it and open the new one
	         ok = w_wnd_close_stream(ch)
	         rltm(ch).file = ' '
	         rltm(ch).file = glsec.filename
	         ok = w_wnd_open_file(ch,rltm(ch).file)
	         if (ok .ne. 1) goto 20
	         major = wnd(ch).first_major
	         minor = wnd(ch).first_minor.i4val
	      else
	         ! the current recording file has not been opened yet,
                 ! so open it now
	         rltm(ch).file = glsec.filename
	         ok = w_wnd_open_file(ch,rltm(ch).file)
	         if (ok .ne. 1) goto 20
	         major = wnd(ch).first_major
	         minor = wnd(ch).first_minor.i4val
	      end if
	   else
	      if (rltm(ch).file(1:1) .ge. ' ') then
	         ! there was a recording file, but now there is not
	         ok = w_wnd_close_stream(ch)
	         rltm(ch).file = ' '
	      end if
	      if (glsec.number_of_interrupts .gt. 491) then
	         ! when purely online and record buffer is full
	         ! the oldest record in the ring buffer is at glsec.pointer+1,
	         ! we allow a margin of ten records because the very oldest
	         ! record is very likely to get overwritten before the user
	         ! can get it
	         i = glsec.pointer + 10
	         if (i.gt.500) i = i - 500
	         major = glsec.recs(i).major_frame
	         minor = zext(glsec.recs(i).minor_frame)
	      else if (glsec.number_of_interrupts .le. 491) then
	         ! when purely online and record buffer is partially
	         ! full we are in a startup state
	         major = glsec.recs(1).major_frame
	         minor = zext(glsec.recs(1).minor_frame)
	      end if
	   end if
	else if (position_type .eq. w_tk_current_mfmf) then
	   major = rltm(ch).rec.major_frame
	   minor = zext(rltm(ch).rec.minor_frame)
	else
	   goto 30
	end if

	w_rt_get_mfmf = 1

	return
  1	format(1x,'W_RT_GET_MFMF: ', a)
 20	continue
	w_rt_get_mfmf = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) 'cannot open current recording file.'
	write(6,'(1x,a,a)',iostat=ios) '  File: ', rltm(ch).file(:60)
	return
 30	continue
	w_rt_get_mfmf = 0
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) 'unknown MF.mf positioning token'
	return
	end

!------------------------------------------------------------------------------
! Performs channel closing operations specific to the realtime TM stream.
!------------------------------------------------------------------------------
	integer*4	function	w_rt_close_stream(ch)
	implicit	none
	include		'parm_def.for'
	include		'glbsec_def.for'
	include		'wind_realtime_def.for'
	integer*4	ch
	integer*4	ok
	integer*4	w_wnd_close_stream
	integer*4	ios

	if (rltm(ch).file(1:1) .gt. ' ') then
	   ok = w_wnd_close_stream(ch)
	   if (ok .ne. 1) then
	      write(6,1,iostat=ios) 'cannot close associated offline file.'
	      write(6,'(1x,a,a)',iostat=ios) '  File: ', rltm(ch).file(:60)
	   end if
	end if

	rltm(ch).file = ' '
	rltm(ch).ch = 0

	w_rt_close_stream = 1
	return
  1	format(1x,'W_RT_CLOSE_STREAM: ', a, :, a)
	end

!------------------------------------------------------------------------------
! Returns one word of data from the realtime stream.
!------------------------------------------------------------------------------
	integer*4	function	w_rt_get_word
	1				(ch,major,minor,word,buf)
	implicit	none
	include		'parm_def.for'
	include		'wind_record_def.for'
	include		'wind_realtime_def.for'
	include		'wind_return_code_def.for'
	integer*4	ch				! TM channel number
	integer*4	major				! major frame number
	integer*4	minor				! minor frame number
	integer*4	word				! word in mnr fr to get
	byte		buf				! byte size buffer

	integer*4	ok
	integer*4	ios
	logical*4	wind_suppress_messages
	integer*4	w_rt_set_rec

	w_rt_get_word = 0

	ok = w_rt_set_rec(ch,major,minor)
	if (.not. ok) goto 10

	buf = rltm(ch).rec.data(word)

	w_rt_get_word = 1
	return
   1	format(1x,'W_RT_GET_WORD: ', a)
 10	w_rt_get_word = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) ' cannot get realtime record.'
	return
	end

!------------------------------------------------------------------------------
! Gets the sc TM mode and converts it to the equivalent cdhf value
!------------------------------------------------------------------------------
	integer*4	function	w_rt_get_mode
	1			(ch,major,minor,mode)
	implicit	none
	include		'parm_def.for'
	include		'wind_record_def.for'
	include		'wind_realtime_def.for'
	include		'wind_return_code_def.for'
	include		'rt_mode_def.for'
	integer*4	ch			! index into user's structure
	integer*4	major			! major frame number
	integer*4	minor			! minor frame number
	integer*4	mode			! tm mode
	integer*4	imode
	integer*4	mode_byte
	integer*4	speed
	logical*4	wind_suppress_messages
	integer*4	ok
	integer*4	w_rt_set_rec
	integer*4	ios
	integer*4	mnr
	integer*4	i
	integer*4	w_rt_tm_speed
	logical*4	caller_wants_speed

	w_rt_get_mode = 0

	mnr = minor
	i = mod(mnr,5)
	if (i .ne. 0) then
	   mnr = mnr - i
	end if

	ok = w_rt_set_rec(ch,major,mnr)
	if (ok.ne.1) goto 10

	imode = rltm(ch).rec.data(4)
	imode = imode .and. '00f0'x
	caller_wants_speed = .false.
	goto 1000

	!----------------------------------------------------------------------
	entry	w_rt_tm_speed(ch,major,mode_byte, speed)
	imode = mode_byte .and. '00f0'x
	caller_wants_speed = .true.

 1000	continue

	if (imode .eq. s1) then
	   if (caller_wants_speed) then
	      speed = 1
	   else
	      mode = science_1x
	   end if
	else if (imode .eq. s2) then
	   if (caller_wants_speed) then
	      speed = 2
	   else
	      mode = science_2x
	   end if
	else if (imode .eq. m1) then
	   if (caller_wants_speed) then
	      speed = 1
	   else
	      mode = maneuver_1x
	   end if
	else if (imode .eq. m2) then
	   if (caller_wants_speed) then
	      speed = 2
	   else
	      mode = maneuver_2x
	   end if
	else if (imode .eq. c1) then
	   if (caller_wants_speed) then
	      speed = 1
	   else
	      mode = contingency_1x
	   end if
	else if (imode .eq. c2) then
	   if (caller_wants_speed) then
	      speed = 2
	   else
	      mode = contingency_2x
	   end if
	else
	   ! just pass on the funny value
	   if (.not. wind_suppress_messages(ch)) then
	      write(6,1,iostat=ios) 'unrecognized mode: ', imode
	   end if
	   if (caller_wants_speed) then
	      speed = 9
	   else
	      mode = imode
	   end if
	end if

	w_rt_get_mode = 1

	return
  1	format(1x,'W_RT_GET_MODE: ', a, i)
 10	continue
	w_rt_get_mode = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) 'error getting record, ok= ', ok
	return
	end

!------------------------------------------------------------------------------
! Sets the current realtime record to specified frame numbers.  May call
! w_rt_get_rec to copy specified record to internal buffer.
!------------------------------------------------------------------------------
	integer*4	function	w_rt_set_rec(ch,major,minor)
	implicit	none
	include		'parm_def.for'
	include		'wind_record_def.for'
	include		'wind_realtime_def.for'
	include		'low_byte_def.for'
	include		'wind_return_code_def.for'
	integer*4	ch
	integer*4	major
	record		/low_byte/	minor
	integer*4	ios
	integer*4	ok
	integer*4	w_rt_get_rec
	logical*4	wind_suppress_messages

	w_rt_set_rec = 0

	if (rltm(ch).rec.major_frame .eq. major .and.
	1   rltm(ch).rec.minor_frame .eq. minor.b) then
	   w_rt_set_rec = 1
	   return
	end if

	ok = w_rt_get_rec(ch, major, minor, rltm(ch).rec)
	if (ok .ne. 1) goto 10

	w_rt_set_rec = 1

	return
  1	format(1x,'W_RT_SET_REC: ', a)
 10	continue
	w_rt_set_rec = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) ' cannot set to specified frame.'
	return
	end

!------------------------------------------------------------------------------
! Copies the indicated record into the record buffer of the specified
! channel.  If the desired record is in the TM gatherer's global section
! the record is simply copied.  If the desired record has already been written
! to disk, then the record is read from disk.  If the desired record is yet
! to be received, the "wait_for_new_records" flag for the channel is checked.
! If wait_for_new_records is true, we use the current bit rate to calculate
! a time span to wait and reiterate.  When wait_for_new_records is false,
! a "not_waiting" error is returned.
!------------------------------------------------------------------------------
	integer*4	function	w_rt_get_rec
	1				(ch,major,minor,wr)
	implicit	none
	include		'parm_def.for'
	include		'glbsec_def.for'
	include		'wind_realtime_def.for'
	include		'low_byte_def.for'
	include		'wind_return_code_def.for'
	integer*4	ch
	integer*4	major
	record		/low_byte/	minor
	record 		/wind_record/	wr	! main wind lib internal buffer
	integer*4	ok
	integer*4	diff
	integer*4	wind_tm_delta_mfmf
	integer*4	ptr, vmajor, vminor
	real*4		wait_time
	integer*4	o
	logical*4	wind_suppress_messages
	integer*4	w_rt_get_offline_rec

	w_rt_get_rec = 0

	if (rltm(ch).rec.major_frame .eq. major .and.
	1   rltm(ch).rec.minor_frame .eq. minor.b) then
	   wr = rltm(ch).rec
	   w_rt_get_rec = 1
	   return
	end if

 100	ptr = glsec.pointer
	if (ptr.lt. 1) goto 10
	vmajor = glsec.recs(ptr).major_frame
	vminor = zext(glsec.recs(ptr).minor_frame)

	! routine wind_tm_delta_mfmf checks ranges of first 4 args
	ok = wind_tm_delta_mfmf(major,minor,vmajor,vminor,diff)
	if (.not. ok) goto 30

	if (diff .eq. 0) then
	   ! current record is the desired record
	   wr = glsec.recs(ptr)
	   w_rt_get_rec = 1
	else if (diff .lt. 0) then
	   ! desired record is in the future
	   if (rltm(ch).wait_for_new_records) then
	      ! estimate the wait time
	      wait_time = glsec.frame_rate * abs(diff)
	      wait_time = amax1(wait_time, 0.1)
	      wait_time = amin1(wait_time, 1.0)
	      call lib$wait(wait_time)
	      goto 100
	   else
	      goto 40
	   end if
	else if (diff .lt. 490) then
	   ! desired record is in the global section
	   ptr = ptr - diff
	   if (ptr .lt. 1) ptr = 500 + ptr
	   wr = glsec.recs(ptr)
	   if (wr.minor_frame .ne. minor.b) goto 50
	   w_rt_get_rec = 1
	else
	   ! the desired record is on disk if recording, update
	   ok = w_rt_get_offline_rec(ch,major,minor,wr)
	   if (ok .ne. 1) goto 60
	   w_rt_get_rec = 1
	end if

	return
  1	format(1x,'W_RT_GET_REC: ', a, a)

 10	w_rt_get_rec = w_no_realtime
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=o)
	1 'Global section pointer is zero, is GSE up? Telemetry?'
	return

 30	w_rt_get_rec = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=o) 'Invalid frame number(s).'
	return

 4	format(1x,'W_RT_GET_REC: ',a,' major=',i,',minor=',i3,'.')
 40	w_rt_get_rec = w_no_record_waiting
	if (wind_suppress_messages(ch)) return
	write(6,4,iostat=o) 'NOT waiting for record',
	1	major, minor.i4val
	return

 5	format(1x,a,i3,a,i3,a,i3,'.')
 50	w_rt_get_rec = w_unexpected_record
	if (wind_suppress_messages(ch)) return
	write(6,5,iostat=o) 'W_RT_GET_REC: Unexpected record, wanted mf ',
	1	 minor.i4val, ' but got ',
	1	 zext(rltm(ch).rec.minor_frame), ', ptr=', ptr
	return

 60	w_rt_get_rec = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=o) 'cannot get specified offline disk record.'
	return
	end

!------------------------------------------------------------------------------
! This routine attempts to get the specified record from the offline .wnd file.
! There are several possibliities:
!	1) the specified record is found in the current offline file
!	2) disk recording is not enabled, the record is ungettable
!	3) disk recording is enabled, but recording has stopped and restarted
!	   since the user first opened the realtime stream, leaving the
!	   the desired record in an (as yet) unopened file.
!------------------------------------------------------------------------------
	integer*4	function	w_rt_get_offline_rec(ch,major,minor,wr)
	implicit	none
	include		'parm_def.for'
	include		'glbsec_def.for'
	include		'wind_realtime_def.for'
	include		'low_byte_def.for'
	include		'wind_wnd_def.for'
	include		'wind_return_code_def.for'
	integer*4	ch
	integer*4	major
	integer*4	minor
	record /wind_record/	wr		! main wind lib internal buffer

	integer*4	w_wnd_get_rec_xc
	integer*4	w_wnd_close_stream
	integer*4	w_wnd_open_file
	logical*4	wind_suppress_messages
	integer*4	ok
	integer*4	ios

	w_rt_get_offline_rec = 0

	! case of not currently recording
	if (glsec.filename(1:1) .le. ' ') then
	   ! not recording
	   if (rltm(ch).file(1:1) .gt. ' ') then
	      ! there was an offline recording file, close it now
	      ok = w_wnd_close_stream(ch)
	      rltm(ch).file = ' '
	   end if
	   goto 10
	end if

	! case of recording, but this is the first disk access
	if (rltm(ch).file(1:1) .le. ' ') then
	    rltm(ch).file = glsec.filename
	    ok = w_wnd_open_file(ch,rltm(ch).file)
	    if (ok .ne. 1) goto 20
	end if

	! update the end-of-file frame numbers
	if (wnd(ch).last_major .lt. glsec.last_major_written) then
	   wnd(ch).last_major = glsec.last_major_written
	   wnd(ch).last_minor.b = glsec.last_minor_written_b
	end if

	ok = w_wnd_get_rec_xc(ch, major, minor, wr)

	if (ok .ne. 1) then
	   ! has recording restarted?
	   if (glsec.filename .ne. rltm(ch).file) then
	      ! recording has restarted
	      ok = w_wnd_close_stream(ch)
	      rltm(ch).file = glsec.filename
	      ok = w_wnd_open_file(ch,rltm(ch).file)
	      if (ok .ne. 1) goto 20
	      ok = w_wnd_get_rec_xc(ch, major, minor, wr)
	      if (ok .ne. 1) goto 30
	   end if
	end if

	w_rt_get_offline_rec = 1

	return
  1	format(1x,'W_RT_GET_OFFLINE_REC: ', a, a)
  2	format(1x,'W_RT_GET_OFFLINE_REC: ', a, /,
	1	'User Ch:',i2,1x,i12,'.',i3.3)

  10	w_rt_get_offline_rec = w_not_recording
	if (wind_suppress_messages(ch)) return
	type 2, 'Record not in buffer and disk recording not enabled.',
	1	ch, major, minor
	return

  20	w_rt_get_offline_rec = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) 'cannot open offline file named:'
	write(6,'(1x,a)',iostat=ios) rltm(ch).file(:60)
	rltm(ch).file = ' '
	return

  30	w_rt_get_offline_rec = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) 'cannot get offline record:'
	return
	end

!------------------------------------------------------------------------------
! This routine enables the user to control wether or not to wait for
! new telemetry to come in when the user requests a minor frame that is
! yet to be received.
!------------------------------------------------------------------------------
	integer*4	function	wind_tm_set_wait(ch)
	implicit	none
	integer*4	ch
	include		'parm_def.for'
	include		'wind_record_def.for'
	include		'wind_realtime_def.for'
	include		'wind_return_code_def.for/nolist'
	logical*4	wind_suppress_messages

	wind_tm_set_wait = 0
	if (ch .lt. 1 .or. ch .gt. max_channels) goto 8

	rltm(ch).wait_for_new_records = .true.

	wind_tm_set_wait = 1
	return
  1	format(1x,'WIND_TM_SET_WAIT: ', a)
 08	wind_tm_set_wait = w_bad_channel
	if (wind_suppress_messages(ch)) return
	type 1, 'invalid channel number.'
	return
	end

!------------------------------------------------------------------------------
! This routine enables the user to control wether or not to wait for
! new telemetry to come in when the user requests a minor frame that is
! yet to be received.
!------------------------------------------------------------------------------
	integer*4	function	wind_tm_set_nowait(ch)
	implicit	none
	integer*4	ch
	include		'parm_def.for'
	include		'wind_record_def.for'
	include		'wind_realtime_def.for'
	include		'wind_return_code_def.for/nolist'
	logical*4	wind_suppress_messages

	wind_tm_set_nowait = 0
	if (ch .lt. 1 .or. ch .gt. max_channels) goto 8

	rltm(ch).wait_for_new_records = .false.

	wind_tm_set_nowait = 1
	return
  1	format(1x,'WIND_TM_SET_NOWAIT: ', a)
 08	wind_tm_set_nowait = w_bad_channel
	if (wind_suppress_messages(ch)) return
	type 1, 'invalid channel number.'
	return
	end

!------------------------------------------------------------------------------
! Maps to the common global section used by the Gatherer process
! to buffer telemetry data as it comes in.
!------------------------------------------------------------------------------
	integer*4	function	wind_open_realtime(ch)
	implicit	none
	integer*4	ch
	include		'parm_def.for'
	include		'glbsec_def.for/nolist'
	include		'wind_realtime_def.for/nolist'
	include		'low_byte_def.for'
	include		'wind_return_code_def.for/nolist'
	include		'($secdef)'
	integer*4	addr(2)
	integer*4	sys$mgblsc
	integer*4	flags /sec$m_sysgbl/
!	parameter	null=char(0)
	integer*4	ok
	integer*4	i
	integer*4	first_time /1/
	logical*4	wind_suppress_messages

	wind_open_realtime = 0

	! map to the TM gatherer's global section
	if (first_time) then
	   addr(1) = %loc( glsec )
	   addr(2) = addr(1) + sizeof(glsec) - 1
	   ok = sys$mgblsc(addr,,,%val(flags),
	1	wind_tm_glbsec_name,,)
	   if (.not. ok) goto 10
	   first_time = 0
	end if

	! read in the most recent TM record and make it the user's current rec
	i = glsec.pointer
	if (i .lt. 1 .or. i .gt. 500) goto 20
	rltm(ch).rec = glsec.recs(i)

	rltm(ch).num_tm_ptrs	= glsec.num_tm_keeper_words
	do i=0,sizeof_minor_frame-1
	   rltm(ch).is_a_valid_ptr(i) = .false.
	end do
	do i=0,rltm(ch).num_tm_ptrs
	   rltm(ch).tm_ptrs(i) = zext(glsec.tm_keeper_words(i))
	   rltm(ch).is_a_valid_ptr(rltm(ch).tm_ptrs(i)) = .true.
	end do

	wind_open_realtime = 1

	return
  1	format(1x,'WIND_OPEN_REALTIME: ',a,:,a)
 10	wind_open_realtime = w_no_realtime
	if (wind_suppress_messages(ch)) return
	type 1, 'cannot access TM realtime data buffers.'
	type '(1x,a,z8.8)',
	1	'Cannot map '//wind_tm_glbsec_name//
	1	' global section, status: ',
	1	 ok
	type *, 'GSE software is probably not running (gss/start/all).'
	return
 20	wind_open_realtime = w_ungettable_record
	if (wind_suppress_messages(ch)) return
	type 1, 'invalid realtime ring buffer index.'
	type *, 'Is telemetry active/available?'
	return
	end

!------------------------------------------------------------------------------
! Gathers a packet of data from a .wnd file or the realtime stream.
! If argument "type" is not a valid
! packet type (between 0 and 15) the current packet, specified by major,minor
! frame number arguments, is selected.  Otherwise the type, direction, and
! seek_1st flag arguments are used to prescribe packet search and selection
! criteria.
! The user's major and minor arguments are advanced to point to the last
! frame of the selected packet on return from this routine.
!------------------------------------------------------------------------------
	integer*4	function	w_rt_get_packet
	1			 (ch,major,minor,type,direction,seek_1st,pkt)
	implicit	none
	include		'parm_def.for'
	include		'wind_record_def.for'
	include		'wind_realtime_def.for'
	include		'wind_return_code_def.for'

	integer*4	ch
	integer*4	major
	integer*4	minor
	integer*4	type
	integer*4	direction
	logical*4	seek_1st
	byte		pkt(0:*)

	integer*4	w_rt_set_rec
	integer*4	zext_bits_to_longword
	logical*4	wind_suppress_messages
	integer*4	i,j,k,n
	include		'wind_packet_addresses_def.for'
	integer*4	ok
	integer*4	mode
	integer*4	my_type
	integer*4	word
	integer*4	check

	logical*4	good_frame
	logical*4	any_bad_frame
	logical*4	science_mode
	logical*4	maneuver_mode
	logical*4	first_frame_ok
	logical*4	by_type

	w_rt_get_packet = 0

	! adjust minor frame number to point to beginning of packet
	! --this could let major,minor indicate a nonexistent record,
	! --such as could happen at the beginning of a file
	if (mod(minor,10) .ne. 0) then
	   minor = minor - mod(minor,10)
	end if

	! This loop gets the first frame of the packet, searching forwards
	! or backwards through the data looking for a specific packet type
	! or a first packet of a specific event (packet_type=event_type, here)
	! as prescribed by the calling arguments.
	by_type = type .ge. min_packet_id .and. type .le. max_packet_id
	first_frame_ok = .false.
	do while(.not. first_frame_ok)
	   ok = w_rt_set_rec(ch,major,minor)
	   if (ok.ne.1) goto 10
	   ! determine the mode: science or maneuver
	   mode = zext(rltm(ch).rec.data(4))
	   mode = mode .and. '00f0'x
	   science_mode = mode .eq. '00'x .or. mode .eq. '40'x
	   maneuver_mode = mode .eq. '10'x .or. mode .eq. '50'x
!	   if (.not. (science_mode .or. maneuver_mode)) science_mode = 1
	   if (.not. (science_mode .or. maneuver_mode)) goto 30
	   if (by_type) then
	      k = 1
	      if (maneuver_mode) k = 432
	      word = zext(rltm(ch).rec.data( p_addrs(k) ) )
	      my_type = zext_bits_to_longword(word,4,4)
	      first_frame_ok = my_type .eq. type
	      if (first_frame_ok .and. seek_1st) then
	         !first_frame_ok = word	! low bit set is 1st packet flag on vax
	         first_frame_ok = (word .and. 1) .ne. 0 ! sun fortran
	      end if
	      if (.not. first_frame_ok) then
	         if (direction .eq. forward) then
	            call wind_tm_increment_packet(major,minor)
	         else if (direction .eq. reverse) then
	            call wind_tm_decrement_packet(major,minor)
	         else
	            goto 20
	         end if
	      end if
	   else
	      first_frame_ok = .true.
	   end if
	end do

	! now gather successive frames and packet words
	k = 1
	if (maneuver_mode) k = 432
	j = 0
	any_bad_frame = .false.
	good_frame = .true.
	do i=1,n_frames_per_packet
	   if (good_frame) then
	      ! copy the data words from the TM record
	      do n=1,n_words_frame(i)
	         pkt(j) = rltm(ch).rec.data( p_addrs(k) )
	         k = k + 1
	         j = j + 1
	      end do
	   else
	      ! set all the bits in each word in the storage area derived
	      ! from this invalid/missing frame
	      do n=1,n_words_frame(i)
	         pkt(j) = 'ff'x
	         k = k + 1
	         j = j + 1
	      end do
	   end if
	   if (i .lt. n_frames_per_packet) then
	      call wind_tm_increment_mfmf(major,minor)
	      ok = w_rt_set_rec(ch,major,minor)
	      good_frame = ok .eq. 1
	      if (ok .eq. w_end_of_file) then
	         goto 40
	      else if (ok .ne. 1) then
	         ! return error code from the get_record routine
	         if (.not. any_bad_frame) w_rt_get_packet = ok
	         any_bad_frame = 1
	         if (.not. wind_suppress_messages(ch)) then
	             type 2, i-1, 'bad frame in packet, MF.mf=', major,minor
	         end if
	      end if
	   end if
	end do

	if (.not. any_bad_frame) w_rt_get_packet = 1

	call wind_plain_packet_sum(pkt,check)
	if (check .ne. 0) goto 50

	return
  1	format(1x,'W_RT_GET_PACKET: ',a, i8,'.',i3.3)
  2	format(1x,t1,'(',i2,')',1x,a,i8,'.',i3.3)
  3	format(1x,'W_RT_GET_PACKET: ',a, i)
  4	format(1x,'W_RT_GET_PACKET: ',a, z2.2)
 10	w_rt_get_packet = ok
	if (wind_suppress_messages(ch)) return
	type 1, 'cannot get 1st frame of packet, MF.mf=', major,minor
	return
 20	w_rt_get_packet = 0
	if (wind_suppress_messages(ch)) return
	type 3, 'invalid "direction" argument: ', direction
	return
 30	w_rt_get_packet = 0
	if (wind_suppress_messages(ch)) return
	type 4, 'invalid SC mode from word #4: ',
	1 rltm(ch).rec.data(4)
	return
 40	w_rt_get_packet = w_end_of_file
	if (wind_suppress_messages(ch)) return
	type 1, 'EOF detected at ', major,minor
	return
 50	w_rt_get_packet = w_bad_checksum
	if (wind_suppress_messages(ch)) return
	type 1, 'invalid checksum in packet at (major.minor): ', major, minor
	return
	end

!------------------------------------------------------------------------------
! For realtime TM access to the ring buffer global section this routine
! is not meaningful.  This routine is therefore implemented as an interface
! to the w_wnd_get_rec_by_time routine for the offline file currently
! being recorded.
!------------------------------------------------------------------------------
	integer*4	function	w_rt_get_rec_by_time
	1		(ch,arg_major,arg_minor,user_time,by_ert,cpy_to_buf,wr)
	implicit	none
	include		'parm_def.for'
	include		'glbsec_def.for'
	include		'wind_realtime_def.for'
	include		'low_byte_def.for'
	include		'wind_wnd_def.for'
	include		'wind_return_code_def.for'
	integer*4	ch		! channel id
	integer*4	arg_major	! major frame number
	record /low_byte/ arg_minor	! minor frame number
	character*(*)	user_time	! 23/24 character ERT or SCET
	logical*4	by_ert		! flag for search by ERT or SCET
	logical*4	cpy_to_buf	!
	record /wind_record/ wr		! wind_tm_lib internal buffer

	integer*4	ok
	integer*4	ios
	integer*4	i,j
	integer*4	w_rt_get_mfmf
	integer*4	w_wnd_get_rec_by_time
	logical*4	wind_suppress_messages

	w_rt_get_rec_by_time = 0

	ok = w_rt_get_mfmf(ch,i,j, w_tk_earliest_mfmf)
	if (ok .ne. 1) goto 10

	if (rltm(ch).file(1:1) .le. ' ') goto 20

	! update the end-of-file frame numbers
	if (wnd(ch).last_major .lt. glsec.last_major_written) then
	   wnd(ch).last_major = glsec.last_major_written
	   wnd(ch).last_minor.b = glsec.last_minor_written_b
	end if

	ok = w_wnd_get_rec_by_time(
	1    ch,arg_major,arg_minor,user_time,by_ert,cpy_to_buf,wr)

	w_rt_get_rec_by_time = ok

	return
  1	format(1x,'W_RT_GET_REC_BY_TIME: ', a)
 10	continue
	w_rt_get_rec_by_time = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) 'error determining offline file position.'
	return
 20	continue
	w_rt_get_rec_by_time = 0
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) 'no current offline recording file.'
	write(6,1,iostat=ios) 'realtime access by time not supported.'
	return
	end

!------------------------------------------------------------------------------
! Returns the SCET of the specified frame.
!------------------------------------------------------------------------------
	integer*4	function	w_rt_scet_of_mfmf(ch,major,minor,scet)
	implicit	none
	include		'parm_def.for'
	include		'wind_record_def.for'
	include		'wind_realtime_def.for'
	integer*4	ch
	integer*4	major
	integer*4	minor
	record /vms_64bit_time/ scet
	integer*4	ok
	integer*4	w_rt_set_rec
	integer*4	wind_suppress_messages
	integer*4	ios
	integer*4	wnd_to_dbms			! a function
	integer*4	ymd,hms,s1000
	integer*4	w_rt_scet_dbms_of_mfmf
	logical*4	do_dbms_cnvrt

	do_dbms_cnvrt = .false.
	goto 1000

	!----------------------------------------------------------------------
	entry		w_rt_scet_dbms_of_mfmf(ch,major,minor,ymd,hms,s1000)
	do_dbms_cnvrt = .true.

 1000	continue
	w_rt_scet_of_mfmf = 0

	ok = w_rt_set_rec(ch,major,minor)
	if (ok .ne. 1) goto 10

	if (do_dbms_cnvrt) then
	   ok = wnd_to_dbms(rltm(ch).rec.scet, ymd,hms,s1000)
	   w_rt_scet_dbms_of_mfmf = ok
	   return
	else
	   scet = rltm(ch).rec.scet
	end if

	w_rt_scet_of_mfmf = 1

	return
  1	format(1x,'W_RT_SCET_OF_MFMF: ', a)
 10	continue
	w_rt_scet_of_mfmf = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) 'cannot get specified realtime[/.wnd] record.'
	return
	end

!------------------------------------------------------------------------------
! Determines the "frame rate" for a given minor frame by examining the mode
! bits in data word 4.
!------------------------------------------------------------------------------
	integer*4	function	w_rt_frame_rate
	1				(ch,major,minor,rate)
	implicit	none
	include		'parm_def.for'
	include		'wind_record_def.for'
	include		'wind_realtime_def.for'
	include		'wind_return_code_def.for'
	integer*4	ch			! index into user's data struc
	integer*4	major			! major frame number
	integer*4	minor			! minor frame number
	real*4		rate			! bit rate
	integer*4	ok			! for checking return statuses
	integer*4	ios
	integer*4	mnr
	integer*4	mode
	record /wind_record/ xwr
	integer*4	w_rt_set_rec
	logical*4	wind_suppress_messages
	parameter	s1='0000'x		! science 1x, 92 sec
	parameter	m1='0010'x		! maneuver 1x, 92 sec
	parameter	e ='0020'x		! engineering
	parameter	c1='0030'x		! contingency 1x
	parameter	s2='0040'x		! science 2x, 46 sec
	parameter	m2='0050'x		! maneuver 2x, 46sec
	parameter	na='0060'x		! not applicable (I think)
	parameter	c2='0070'x		! contingency 2x
	real*4		hi,lo
	parameter	(hi=46.0/250.0)
	parameter	(lo=92.0/250.0)

	w_rt_frame_rate = 0
	rate = 0.0
	mnr = minor
	xwr = rltm(ch).rec

	if (mod(mnr,5) .ne. 0) then		! frames 0,5,10,15,...240,245
	   mnr = mnr - mod(mnr,5)		! contain the mode in word 4
	end if

	! establish the current record
	ok = w_rt_set_rec(ch, major, mnr)	! position to mf with mode
	mode = zext(rltm(ch).rec.data(4))	! get the mode
	rltm(ch).rec = xwr			! restore position
	if (ok .ne. 1) goto 10			! branch on error

	mode = mode .and. '00f0'x

	if (mode .eq. s2) then
	   rate = hi
	else if (mode .eq. s1) then
	   rate = lo
	else if (mode .eq. m1) then
	   rate = lo
	else if (mode .eq. m2) then
	   rate = hi
	else if (mode .eq. c1) then
	   rate = lo
	else if (mode .eq. c2) then
	   rate = hi
	else
	   goto 20
	end if

	w_rt_frame_rate = 1

	return
  1	format(1x,'W_RT_FRAME_RATE: ', a)
 10	continue
	w_rt_frame_rate = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) 'cannot get specified realtime[/.wnd] record.'
	return
 20	continue
	w_rt_frame_rate = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) 'invalid spacecraft mode.'
	return
	end

!------------------------------------------------------------------------------
! This routine stores extra information related to the event for later
! retrieval by the user if desired.  Stored information includes:
!    1) major&minor frame numbers of the 1st frame of 1st packet of the event
!    2) Earth Receive Time of the 1st frame of 1st packet of the event in
!       DBMS format
!    3) DPU major frame number from HK enclosing 1st packet of the event
!    4) SCET of the event
!    5) ScriptPlayer TEST and STEP values
!
!------------------------------------------------------------------------------
	integer*4	function	w_rt_get_xtra_event_info
	1				(ch,major,minor)
	implicit	none
	include		'parm_def.for'
	include		'low_byte_def.for'
	include		'wind_record_def.for'
	include		'wind_extra_info_def.for'
	include		'wind_realtime_def.for'
	include		'wind_return_code_def.for'
	integer*4	ch				! TM channel number
	integer*4	major				! major frame number
	integer*4	minor				! minor frame number
	integer*4	ok				! return status var
	integer*4	w_rt_set_rec			! a function
	integer*4	wind_get_event_scet		! a function
	integer*4	wind_get_record			! a function
	integer*4	wnd_to_dbms			! a function
	logical*4	wind_suppress_messages
	record	/wind_record/	my_record
	integer*4	xmajor, xminor
	integer*4	w_get_dpu_mjr_fr_num		! a function
	integer*4	w_rt_tm_speed			! a function
	integer*4	speed
	integer*4	return_code
	integer*4	o
	parameter	w_scet_err	= '010'x
	parameter	w_scet_frctn_err= '020'x
	parameter	w_dpu_majf_err	= '040'x
	parameter	w_ert_err	= '100'x
	parameter	w_ert_frctn_err	= '200'x

	w_rt_get_xtra_event_info = 0
	return_code = 1

	! get the major and minor frame numbers of the beginning of the event
	exi(ch).major = major
	exi(ch).minor = minor - 9

	! establish beginning of event as the current record
	ok = wind_get_record(ch,major,exi(ch).minor,my_record)
	if (ok .ne. 1) goto 10

	! convert the ert time
	ok = wnd_to_dbms(my_record.gathertime,
	1	exi(ch).ert(1),
	1	exi(ch).ert(2),
	1	exi(ch).ert1000)
	if (ok .ne. 1) goto 20

	! get the script player's test and step number
	exi(ch).sp_test_number = rltm(ch).rec.sp_test_number
	exi(ch).sp_step_number = rltm(ch).rec.sp_step_number

	xmajor = major
	xminor = minor
	ok = w_get_dpu_mjr_fr_num(ch,xmajor,xminor,exi(ch).dpu_major_ert)
	if (ok .ne. 1) return_code = return_code .or. w_dpu_majf_err

	! get the bit_rate, put in event buffer
	ok = w_rt_tm_speed(ch,xmajor,zext(my_record.data(4)),speed)
	if (speed .eq. 1) then
	   exi(ch).bit_rate = 0
	else if (speed .eq. 2) then
	   exi(ch).bit_rate = 1
	else
	   exi(ch).bit_rate = 9
	end if

	! get the spacecraft event time
	xmajor = major
	xminor = minor
	ok = wind_get_event_scet(ch,xmajor,xminor)
	if (ok .ne. 1) then
	   return_code = return_code .or. w_scet_err
	   ! use the ert to grossly estimate the scet
	   exi(ch).scet(1) = exi(ch).ert(1)
	   exi(ch).scet(2) = exi(ch).ert(2)
	end if

	! restore user's stream position to end of first packet of event
	ok = w_rt_set_rec(ch,major,minor)
	if (ok .ne. 1) goto 80

	w_rt_get_xtra_event_info = return_code

	if ( (.not. wind_suppress_messages(ch)) .and. (.not.return_code)) then
	   if ((return_code .and. w_scet_err) .ne. 0) then
	      write(6,1,iostat=o) 'cannot get event SCET.'
	   end if
	   if ((return_code .and. w_ert_frctn_err) .ne. 0) then
	      write(6,1,iostat=o) 'cannot get fractional ERT.'
	   end if
	   if ((return_code .and. w_dpu_majf_err) .ne. 0) then
	      write(6,1,iostat=o) 'cannot get DPU MF.'
	   end if
	end if

	return
 1	format(1x,'W_RT_GET_XTRA_EVENT_INFO: ', a, :, a)

 10	w_rt_get_xtra_event_info = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=o) 'cannot get record for ERT.'
	return

 20	w_rt_get_xtra_event_info = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=o) 'cannot convert record ERT to DBMS format time.'
	return

 80	w_rt_get_xtra_event_info = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=o) 'cannot get record for position restoration.'
	return
	end

!------------------------------------------------------------------------------
! Returns one to N (N=w2-w1+1) bytes of WND TM HK data to caller.
! This routine returns the original (more or less) 81 HK words, 25
! packet id words, 125 bytes of minor frame quality (always zero,
! but provided for cdhf compatability in form), 1 byte counting the
! number of minor frames with fill (provided for cdhf compatablility in
! form, though the number of errors encountered gathering the first 106
! hk words), 1 byte containing the number of minor frames with synch err
! (same value as number of minor frames with fill).
!------------------------------------------------------------------------------
	integer*4	function	w_rt_get_hk(ch,mjr,w1,w2,buf)
	implicit	none
	include		'wind_os_def.for'
	include		'low_byte_def.for'
	include		'parm_def.for'
	include		'wind_record_def.for'
	include		'wind_realtime_def.for'
	include		'wind_hk_addr_def.for'
	include		'wind_return_code_def.for'
	integer*4	ch		! caller's channel number
	integer*4	mjr		! major frame number
	integer*4	w1		! first/lowest hk word number
	integer*4	w2		! second/highest hk word number
	byte		buf(*)		! caller's buffer to receive word(s)
	integer*4	ok
	integer*4	w_rt_get_word
	integer*4	w_rt_get_mode
	logical*4	wind_suppress_messages	! a function
	integer*4	ios
	integer*4	mnr
	integer*4	mode
	integer*4	i,j,k
	record /low_byte/ lb
	integer*4	n_err

	w_rt_get_hk = 0

	j = 0
	n_err = 0
	do i=w1,w2
	   j = j + 1
	   buf(j) = 0
	   if (i .ge. w_first_core_hk_word .and. 		! reg. hk
	1      i .le. w_last_core_hk_word) then
	      ok = w_rt_get_word(ch, mjr,
	1	   hk_mf_addr(i),		! minor frame number
	1	   hk_word_addr(i),		! word number in minor frame
	1	   buf(j))			! returned word
	      if (ok .ne. 1) n_err = n_err + 1
	   else if (i .ge. w_first_pktid_hk_word .and. 		! pkt id's
	1           i .le. w_last_pktid_hk_word) then
	      mnr = (i - w_first_pktid_hk_word) * 10
	      ok = w_rt_get_mode(ch,mjr,mnr,mode)
	      if (ok .eq. 1) then
	         if (mode .eq. science_1x .or.
	1            mode .eq. science_2x .or.
	1            mode .eq. contingency_1x .or.
	1            mode .eq. contingency_2x) then
	             k = 24
	         else if (mode .eq. maneuver_1x .or.
	1              mode .eq. maneuver_2x) then
	             k = 32
	         else
	             goto 35
	         end if
	         ok = w_rt_get_word(ch, mjr, mnr, k, buf(j))
	         if (ok .ne. 1) n_err = n_err + 1
	      else
	         n_err = n_err + 1
	      end if
	   else if (i .eq. w_sc_mode_hk_word) then
	      ok = w_rt_get_mode(ch,mjr,mnr,lb.i4val)
	      if (ok .ne. 1) write(6,1,iostat=ios) 'cannot get HK sc mode word.'
	      buf(j) = lb.b
	   else if (i .ge. w_first_mfq_hk_word .and.
	1           i .le. w_last_mfq_hk_word) then
	      buf(j) = 0
	   else if (i .eq. w_n_mf_w_fill_hk_word) then
	      lb.i4val = n_err
	      buf(j) = lb.b
	   else if (i .eq. w_n_mf_w_synch_err_hk_word) then
	      lb.i4val = n_err
	      buf(j) = lb.b
	   else
	      goto 13
	   end if
	end do

	if (n_err .ne. 0) then
	   write(6,4,iostat=ios) n_err, ' words with errors.'
	   if (n_err .ge. 106) return ! return err if all errors
	   w_rt_get_hk = w_partial_hk_event
	   return
	end if

	w_rt_get_hk = 1

	return
   1	format(1x,'W_RT_GET_HK: ', a, :, i)
   4	format(1x,'W_RT_GET_HK: ', i3, a)
! 10	continue
!	w_rt_get_hk = ok
!	if (wind_suppress_messages(ch)) return
!	write(6,1,iostat=ios) 'cannot goto stream position'
!	return
 13	w_rt_get_hk = w_bad_argument_value
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios) 'bad house keeping index.'
	return
! 20	w_rt_get_hk = ok
!	if (wind_suppress_messages(ch)) return
!	write(6,1,iostat=ios) 'cannot get specified core housekeeping word.'
!	return
! 30	w_rt_get_hk = ok
!	if (wind_suppress_messages(ch)) return
!	write(6,1,iostat=ios)
!	1 'cannot get sc mode for "packet id" housekeeping word.'
!	return
 35	w_rt_get_hk = ok
	if (wind_suppress_messages(ch)) return
	write(6,1,iostat=ios)
	1 'invalid sc mode for "packet id" housekeeping word.', mode
	return
! 38	w_rt_get_hk = ok
!	if (wind_suppress_messages(ch)) return
!	write(6,1,iostat=ios)
!	1 'cannot get specified "packet id" housekeeping word.'
!	return
	end

!------------------------------------------------------------------------------
	integer*4	function	w_rt_get_dpu_mf(ch,major,minor,dpu)
! Extracts the dpu major frame number from the stream at a given position.
! Attempts to estimate the value at EOF and BOF.  Restores the caller
! to the stream position at entry.
	implicit	none
	include		'parm_def.for'
	include		'low_byte_def.for'
	include		'wind_record_def.for'
	include		'wind_realtime_def.for'
	include		'wind_return_code_def.for'
	integer*4	ch				! TM channel number
	integer*4	major				! major frame number
	integer*4	minor
	integer*4	dpu				! dpu major frame #
	integer*4	ok				! return status var
	integer*4	w_rt_get_rec			! a function
	record	/low_byte/	lw
	record	/wind_record/	my_record
	integer*4	omajor,ominor
	integer*4	error_count
	integer*4	wind_get_message_state		! a function
	integer*4	wind_set_message_state		! a function

	logical*4	state
	logical*4	subtract
	logical*4	add

	lw.i4val = 0
	add      = .false.
	subtract = .false.
	omajor   = rltm(ch).rec.major_frame
	ominor   = zext(rltm(ch).rec.minor_frame)
	error_count = 0
	ok = wind_get_message_state(ch, state)
	ok = wind_set_message_state(ch, .true.)

	w_rt_get_dpu_mf = 0

 1000	continue
	if (error_count .gt. 1) goto 2000

	ok = w_rt_get_rec(ch,omajor,44,my_record)
	if (ok .eq. w_beginning_of_file) then
	   subtract = .true.
	   omajor   = omajor + 1
	   error_count = error_count + 1
	   goto 1000
	else if (ok .eq. w_end_of_file) then
	   add      = .true.
	   omajor   = omajor - 1
	   error_count = error_count + 1
	   goto 1000
	else if (ok .ne. 1) then
	   goto 114
	end if
	lw.b =  my_record.data(18)
	lw.b1 = my_record.data(17)

	ok = w_rt_get_rec(ch,omajor,114,my_record)
	if (ok .eq. w_end_of_file) then
	   add      = .true.
	   omajor   = omajor - 1
	   error_count = error_count + 1
	   goto 1000
	else if (ok .ne. 1) then
	   goto 116
	end if
	lw.b2 = my_record.data(17)

 2000	continue
	if (error_count .gt. 1) then
	   w_rt_get_dpu_mf = 0
	else if (error_count .eq. 1) then
	   if (subtract) then
	      lw.i4val = lw.i4val - 1
	   else
	      lw.i4val = lw.i4val + 1
	   end if
	   dpu = lw.i4val
	   w_rt_get_dpu_mf = 1
	else
	   dpu = lw.i4val
	   w_rt_get_dpu_mf = 1
	end if

	ok = w_rt_get_rec(ch,omajor,ominor,my_record)
	if (ok .ne. 1) goto 120

	ok = wind_set_message_state(ch, state)

	return
  1	format(1x,'W_RT_GET_DPU_MF: ',a,:,i,'.',i3.3)
 114	w_rt_get_dpu_mf = ok
	ok = wind_set_message_state(ch, state)
	if (state) return
	type 1, 'cannot get wind record for HK DPU at', omajor,44
	return

! 115	w_rt_get_dpu_mf = ok
!	ok = wind_set_message_state(ch, state)
!	if (state) return
!	type 1, 'cannot get wind record for HK DPU at', omajor,119
!	return

 116	w_rt_get_dpu_mf = ok
	ok = wind_set_message_state(ch, state)
	if (state) return
	type 1, 'cannot get wind record for HK DPU at', omajor,114
	return

 120	w_rt_get_dpu_mf = ok
	ok = wind_set_message_state(ch, state)
	if (state) return
	type 1, 'cannot restore stream position to ', major, minor
	return
	end

ext: for


text:
      MODULE dkes_input
!     Modified indata namelist to be named dkes_input so as not to
!     conflict with VMEC input namelist.
!
!     CONTAINS INPUT PARAMETERS AND ARRAYS READ IN FROM THE DATAIN NAMELIST
!
!     mpol      : Number of Fourier harmonics used to represent the theta
!                 dependence of the distribution functions; 1 < mpol <= mpold
!     mpolb     : Number of Fourier harmonics used to represent the theta
!                 dependence of the equilibrium magnetic field; 1 < mpolb <= mpolbd
!     ntor      : Number of Fourier harmonics used to represent the zeta
!                 dependence of the distributions; 0 < ntor <= ntord
!     ntorb     : Number of Fourier harmonics used to represent the zeta
!                 dependence of the equilibrium field; 0 < ntorb <= ntorbd
!     mmnn      : Matrix containing the distribution of toroidal mode numbers
!                 (expressed in units of nzperiod) in column 1, and the
!                 number of poloidal modes associated with each toroidal
!                 mode in column 2. the poloidal mode numbers associated
!                 with each toroidal mode are stored in the same row,
!                 starting in column 3. WARNING : this spectrum must
!                 encompass the spectrum of the input magnetic field.

!     nrun      : Number of cmul = nu/v  and  efield = Erad/v values
!                 for which solutions are to be calculated; must be < krun
!     cmul      : Corresponding array of cmul values; should be non-zero.
!     efield    : Corresponding array of efield values.
!     nzperiod  : Number of toroidal field periods.
!     chip      : chip = d(chi)/d(rho), radial derivative of the poloidal flux.
!     psip      : psip = d(psi)/d(rho), radial derivative of the toroidal flux
!     btheta    : covariant poloidal component of B. In Boozer coordinates,
!                 a flux function (I, the TOROIDAL current density)
!     bzeta     : covariant toroidal component of B. In Boozer coordinates,
!                 a flux function (J, the POLOIDAL current density)
!     borbi     : Array of nonzero Fourier coefficients of |B| (or 1/|B|) defined by
!                 borbi_realspace (theta,zeta) = SUM borbi(n,m)*cos[m*theta - nzperiod*n*zeta],
!                 where m = 0 to mpolb, n = -ntorb to ntorb
!                 Note: borbi represents "|B|" or 1/|B| ("binverse"), depending on ibbi)
!                 OBSOLETE (see nvalsb array below):
!                 borbi_realspace(theta,zeta) = borbi(n,m)*cos[(m-1)*theta - nzperiod*nvalsb(n)*zeta].
!     ibbi  = 1 : borbi represents B (default, comes from boozmn file, for example)
!     ibbi  = 2 : borbi represents 1/B
!     nvalsb    : (NOW OBSOLETE) Array of ntorb toroidal mode numbers
!                 expressed in units of nzperiod. For an axisymmetric
!                 device, ntorb = 1 and nvalsb(1) = 0.

!     ipmb  = 0 : DKES computes both sine and cosine solutions (default)
!     ipmb  = 1 : DKES computes only sine            solutions.
!     ipmb  = 2 : DKES computes only cosine          solutions.

!     idisk = 0 : DKES computes all l values.
!     idisk = 1 : DKES computes only l = 0,1,2,3,4,5 (default)

!     lfout = 0 : Final distributions not written in the file DKESF. (default)
!     lfout > 0 : Final distributions     written in the file DKESF,
!                 for l+1 = 1,2,.....,lfout <= lalpha.

!     meshtz    : Determines theta and zeta meshes for Fourier transforms.
!                 See ntheta and nzeta in the MNSET subroutine. Used to
!                 ensure accuracy of spectral matrix elements.

!     lalpha    : Number of orthonormalized Legendre polynomials used to
!                 represent the dependence of the distributions on the
!                 pitch-angle variable cos(alpha); lalpha >= 6
!                 NOTE OF IMPORTANCE: the sqrt(l+1/2) factor is built in!

      USE stel_kinds
      INTEGER, PARAMETER :: mpold = 200
      INTEGER, PARAMETER :: ntord = 200
      INTEGER, PARAMETER :: mpolbd = 100
      INTEGER, PARAMETER :: ntorbd = 100
      INTEGER, PARAMETER :: bigint = 2**16
      INTEGER, PARAMETER :: krun = 100

      INTEGER :: mpol, ntor, lalpha, ipmb, meshtz, idisk,
     1   lfout, ifscl, nrun, nzperiod, mpolb, ntorb, ibbi

      INTEGER, DIMENSION(mpold + 2,ntord) :: mmnn

      REAL(rprec), DIMENSION(krun) :: cmul, efield
      REAL(rprec) :: chip, psip, btheta, bzeta
      REAL(rprec), DIMENSION(-ntorbd:ntorbd,0:mpolbd) :: borbi
      INTEGER, DIMENSION(ntorbd) :: nvalsb                              !OBSOLETE
      CHARACTER*1, PARAMETER :: dashes(131) = (/ ('-', idisk=1,131) /)
      CHARACTER*80 :: summary_file            !record file addition
      INTEGER :: itab_out                     !record file addition
      LOGICAL :: lscreen

      NAMELIST /dkes_indata/mpol, ntor, lalpha, ipmb,mmnn,meshtz,idisk,
     1   lfout, ifscl, nrun, cmul, efield, nzperiod, chip, psip, mpolb,
     2   ntorb, nvalsb, borbi, btheta, bzeta, ibbi

      CONTAINS

      SUBROUTINE read_dkes_input
      USE safe_open_mod
      USE vimatrix, ONLY: ioout, ioout_opt
      IMPLICIT NONE
C-----------------------------------------------
C   L o c a l   P a r a m e t e r s
C-----------------------------------------------
      INTEGER, PARAMETER :: idata = 7
      INTEGER, PARAMETER :: iout = 20
      INTEGER, PARAMETER :: iout_opt = 14
      REAL (rprec), PARAMETER :: zero = 0
C-----------------------------------------------
C   L o c a l   V a r i a b l e s
C-----------------------------------------------
      INTEGER :: numargs, icount, istat, index_blank, iodata
      INTEGER :: m, n, nmax, mmax
      CHARACTER :: output_file*55, opt_file*55, input_file*50
      CHARACTER*50 :: arg1(7)
      CHARACTER :: tb*1            !record file addition
      LOGICAL :: lexist
C-----------------------------------------------
C   E x t e r n a l   F u n c t i o n s
C-----------------------------------------------
      EXTERNAL getcarg, dkes_input_prepare
C-----------------------------------------------
      tb = char(9)            !record file addition
      lscreen = .true.

c  Read input data from "datain" namelist
      CALL getcarg(1, arg1(1), numargs)

      DO icount = 2, MIN(7,numargs)
         CALL getcarg(icount, arg1(icount), istat)
      END DO
      IF (numargs .ge. 4) THEN
         IF (numargs .eq. 4) arg1(5) = 'T'
         CALL dkes_input_prepare (arg1, numargs, input_file, 1)
         IF (arg1(5)(1:1).eq.'f' .or. arg1(5)(1:1).eq.'F')
     1      lscreen = .false.
      ELSE IF (numargs .eq. 1) THEN
         input_file = TRIM(arg1(1))
      ELSE
         PRINT *, ' MUST ENTER INPUT FILENAME ON COMMAND LINE'
         STOP
      END IF

      index_blank = INDEX(input_file,'input.')

      IF (index_blank .eq. 0) THEN
         INQUIRE (file=TRIM(input_file), exist=lexist)
         IF (.not.lexist) THEN
            input_file = 'input.' // TRIM(input_file)
            index_blank = 1
         ELSE
            index_blank = INDEX(input_file,'.') - 5
         END IF
      END IF

      index_blank = MAX(index_blank + 5,1)
      output_file= 'dkesout' // input_file(index_blank:)
      opt_file= 'opt_dkes' // input_file(index_blank:)             !DAS 2/21/2000
      summary_file = 'results' // input_file(index_blank:)         !record file addition

!
!     OPEN INPUT AND OUTPUT FILES FOR READING (AND WRITING OUTPUT)
!
      iodata = idata
      CALL safe_open(iodata, istat, input_file, 'old', 'formatted')
      IF (istat .ne. 0) STOP 'Error reading input file in DKES'
      ioout = iout
      CALL safe_open(ioout, istat, output_file, 'replace', 'formatted')
      IF (istat .ne. 0) STOP 'Error writing output file'
      ioout_opt = iout_opt
      CALL safe_open(ioout_opt, istat, opt_file, 'replace','formatted')
      IF (istat .ne. 0) STOP 'Error writing opt_output file'

!     Read namelist (datain) input
      nvalsb(1) = -bigint-1
      idisk = 1
      lfout = 0
      ibbi = 1
      borbi = 0
      READ (iodata, nml=dkes_indata, iostat=istat)
      IF (istat .ne. 0) STOP'Error reading dkes_indata NAMELIST in DKES'
      CLOSE (iodata)

!
!     Recompute ntorb, mpolb for new style input where
!     borbi is input with actual index value, borbi(n,m)
!
      IF (nvalsb(1) <= -bigint) THEN
         nmax = -bigint;  mmax = -bigint
         DO n = -ntorbd, ntorbd
            DO m = 0, mpolbd
               IF (borbi(n,m) .ne. zero) THEN
                  nmax = max (nmax, abs(n))
                  mmax = max (mmax, abs(m))
               END IF
            END DO
         END DO

!        User MAY input smaller values if he does not want to use entire array
         ntorb = min (abs(ntorb), nmax)
         mpolb = min (abs(mpolb-1), mmax)

      ELSE
         IF (ntorb > ntorbd) THEN
            WRITE (ioout, 45) ntorb, ntorbd
            STOP ' ntorb > ntorbd in DKES input'
         END IF
         IF (mpolb < 2) THEN
            WRITE (ioout, 20) mpolb
            STOP ' mpolb < 2 in DKES input'
         ENDIF
  20     FORMAT(' mpolb = ',i5,'  is less than 2')
  45     FORMAT(' ntorb = ',i5,'  is greater than ntorbd = ',i5)

      END IF
!
!     perform error checking on input data

      IF (nzperiod .le. 0) nzperiod = 1
      IF (ibbi<1 .or. ibbi >2) THEN
         WRITE (ioout,'(a)') ' ibbi must =1 or =2'
         STOP ' ibbi <1 or ibbi >2 in DKES input'
      END IF

      IF (mpol < 1) THEN
         WRITE (ioout, 10) mpol
  10     FORMAT(' mpol = ',i5,'  is less than 1')
         STOP ' mpol < 1 in DKES input'
      ENDIF

      IF (mpol > mpold) THEN
         WRITE (ioout, 15) mpol, mpold
  15     FORMAT(' mpol = ',i5,'  is greater than mpold = ',i5)
         STOP ' mpol > mpold in DKES input'
      ENDIF

      IF (mpolb > mpolbd) THEN
         WRITE (ioout, 25) mpolb, mpolbd
  25     FORMAT(' mpolb = ',i5,'  is greater than mpolbd = ',i5)
         STOP ' mpolb > mpolbd in DKES input'
      ENDIF

      IF (ntor < 1) THEN
         WRITE (ioout, 30) ntor
  30     FORMAT(' ntor = ',i5,'  is less than 1')
         STOP ' ntor < 1 in DKES input'
      ENDIF

      IF (ntor > ntord) THEN
         WRITE (ioout, 35) ntor, ntord
  35     FORMAT(' ntor = ',i5,'  is greater than ntord = ',i5)
         STOP ' ntor > ntord in DKES input'
      ENDIF

      IF (lalpha < 6) THEN
         WRITE (ioout, 50) lalpha
  50     FORMAT(' lalpha = ',i5,'  is less than 6')
         STOP ' lalpha < 6 in DKES input'
      ENDIF

      IF (ipmb<0 .or. ipmb>2) ipmb = 0

      meshtz = MAX(0,meshtz)

      END SUBROUTINE read_dkes_input

      END MODULE dkes_input

max_stars_repo_name: joseluisvelasco/STELLOPT

max_issues_repo_name: joseluisvelasco/STELLOPT

max_forks_repo_name: joseluisvelasco/STELLOPT