C$$$  MAIN PROGRAM DOCUMENTATION BLOCK
C                .      .    .                                       .
C MAIN PROGRAM:  RAOBDUMP    Dump & tutorial program for MADIS raob
C   PRGMMR: M. BARTH         ORG: FSL        DATE: 01-07-31
C
C ABSTRACT: Dump MADIS raob station info, observations, and quality
C   control (QC) information for a single time, based on the user's 
C   selections specified in a text parameter file.  The program is also 
C   meant as a tutorial to demonstrate how a user can use the MADIS 
C   library calls to write his or her own programs.
C
C PROGRAM HISTORY LOG:
C   01-07-30  M. BARTH  Original version -- MADIS Version 2.0
C   02-04-30  M. BARTH  V2.1 - Changed comment re available datasets 
C   03-07-02  M. BARTH  V2.3 - Changed comment re available datasets 
C   03-10-02  M. BARTH  V2.4 - Changed comment re available datasets 
C   04-06-17  M. BARTH  V2.5 - Changed comment re available datasets 
C   04-10-15  L. BENJAM V2.7 - Changed comment re available datasets 
C   05-09-16  L. BENJAMIN V2.8 - Changed RAOB Integrated levels 
C                                   from 256 to 285 
C                              - Changed RAOB Significant levels 
C                                   from 150 to 175
C                              - Changed RAOB Wind levels 
C                                   from 76 to 80
C   05-10-27  L. BENJAMIN V2.9 - Changed comment re available datasets 
C   05-11-18  L. BENJAMIN V3.0 - Changed comment re available datasets 
C   06-01-25  M. BARTH    V3.1 - Changed RAOB Wind levels from 80 to
C                                200, and integrated levels from 285
C                                to 405
C   08-07-25  L. BENJAMIN V3.7 - Changed comment re available datasets 
C
C USAGE:
C   INPUT FILES: (for the entire program)
C
C     raobdump.par     Parameter file (in current directory) with
C                      user's selections.
C
C     MADIS initialization files in MADIS_STATIC directory   
C
C     MADIS observation files under MADIS_DATA directories
C
C   OUTPUT FILES:
C
C     raobdump.txt     Text output of station info, obs & QC
C
C     standard output  Log of what was done, errors
C
C SUBPROGRAMS CALLED:
C     LIBRARY: - MADISLIB:  MINIT, MSETWIN, MSETDOM, MSETQC, MSETCOR,
C                           MTRNTIM, MRAOBSTA, MGETIJ, MGETVCN, MGENVCN,
C                           MGETRAOB, MGETRAOBC, M_MKFNAM
C
C REMARKS:
C
C   1) This program can be used to exercise all options included in the
C      MADIS library for raob data.  It is also meant to be used as
C      a tutorial to demonstrate how a user can write his or her own
C      programs.  
C
C   2) See the INSTALL.* files in the MADIS doc directory for 
C      instructions on how to install and build the MADIS library and 
C      programs (INSTALL.unix for Unix/Linux, INSTALL.windows for
C      Windows 9x/Me/2000/NT).
C
C   3) See the raobdump_usage.txt file in the MADIS doc directory for
C      instructions on how to set up a parameter file and run the
C      program.
C    
C ATTRIBUTES:
C   LANGUAGE: FORTRAN 77 + EXTENSIONS
C   MACHINE: UNIX
C
C$$$
C-----------------------------------------------------------------------
C
C     Declarations notes:
C
C     1) Key constants:
C
C       - 0              successful status return
C       - 1              internal MADIS raob dataset code number
C       - 999999.0       missing data flag
C       - 1000           maximum number of radiosonde stations
C
C     2) In addition to dimensioning arrays for the maximum number of
C        1000 stations, they should be dimensioned for the maximum 
C        number of levels in the desired sounding type:
C
C        Sounding Type          Maximum Number of Levels
C        -------------          ------------------------
C        Integrated                      405
C        Mandatory Levels                 22
C        Significant Levels              175
C        Wind Levels                     200
C        Maximum Wind Levels               4
C        Tropopause Levels                 4

      implicit none

      integer success_p,raob_ds_p
      real miss_p
      parameter (success_p = 0)
      parameter (raob_ds_p = 1)
      parameter (miss_p = 999999.0)

      integer istatus,i,j,k,l,line
      integer nsta,nvar,nobs,idb
      character*10 stanam(1000)
      character*9 atime
      character*6 db
      character*9 timeob(1000),reltime(1000)
      integer wmoid(1000)
      real lat(1000),lon(1000),elev(1000)
      integer nlevels(1000)
      integer vcn,sounding_type(4)
      character*7 vcode(4)
      character*19 vcstr(4)
      character*1 cdummy
      character*13 tstr1,tstr2
      character*70 vline
      integer minbck,minfwd,recwin,qctyp,idmtyp,tform,ctime
      character*10 singsta
      integer corr_type
      real dummy,ri(1000),rj(1000)
      integer idummy,time,ising
      real rpar(7)
      integer ipar(2)
      character*63 twstr(0:4)
      integer twstrlen(0:4)
      character*16 qcstr(0:2)

      real ltype(405,1000),stype(1000)

      real rdata(405,1000,4)
      integer qca(405,1000,4),qcr(405,1000,4)
      character*1 qcd(405,1000,4)

      real rdata_man(22,1000,4)
      integer qca_man(22,1000,4),qcr_man(22,1000,4)
      character*1 qcd_man(22,1000,4)

      real rdata_sigt(175,1000,4)
      integer qca_sigt(175,1000,4),qcr_sigt(175,1000,4)
      character*1 qcd_sigt(175,1000,4)

      real rdata_sigw(200,1000,4)
      integer qca_sigw(200,1000,4),qcr_sigw(200,1000,4)
      character*1 qcd_sigw(200,1000,4)

      real rdata_mxwnd(4,1000,4)
      integer qca_mxwnd(4,1000,4),qcr_mxwnd(4,1000,4)
      character*1 qcd_mxwnd(4,1000,4)

      real rdata_trop(4,1000,4)
      integer qca_trop(4,1000,4),qcr_trop(4,1000,4)
      character*1 qcd_trop(4,1000,4)

      equivalence (rdata,rdata_man,rdata_sigt,rdata_sigw,
     1             rdata_mxwnd,rdata_trop)
      equivalence (qca,qca_man,qca_sigt,qca_sigw,
     1             qca_mxwnd,qca_trop)
      equivalence (qcr,qcr_man,qcr_sigt,qcr_sigw,
     1             qcr_mxwnd,qcr_trop)
      equivalence (qcd,qcd_man,qcd_sigt,qcd_sigw,
     1             qcd_mxwnd,qcd_trop)

      character*7 level_type_vcode,release_time_vcode,sonde_type_vcode
      character*1 ltstr(6)

      data level_type_vcode/'LEVTYPE'/
      data release_time_vcode/'RELTIME'/
      data sonde_type_vcode/'SONDET'/
      data ltstr/'s','m','t','w','x','r'/
      data twstr/
     1 'Return all records within the file containing nominal time',  
     2 'Return one record per fixed station, closest to nominal time',
     3'Return one record per fixed station, closest to start of window',
     4 'Return one record per fixed station, closest to end of window',
     5 'Return all records in window'/
      data twstrlen/58,60,63,61,28/
      data qcstr/'No Filter','Coarse','Screened'/

C-----------------------------------------------------------------------
C PROGRAM SETUP
C -------------
C
C     Open an output file for this program.

      open(unit=1,file='raobdump.txt',status='UNKNOWN')

C     Open the parameter file that will specify the user's selections.

      open(unit=2,file='raobdump.par',status='OLD')

C-----------------------------------------------------------------------
C MANDATORY MADIS INITIALIZATION
C -----------------------------
C
C     Mandatory initialization which must be done before
C     any other MADIS library calls.
C
C     Initialize the MADIS radiosonde subsystem
C     -----------------------------------------
C
C     The user selects the "subsystem" to be initialized by
C     specifying the following arguments to the MINIT call:
C
C     1) Name of the subsystem:
C
C        'SFC'        Surface
C        'RAOB'       Radiosonde              (used here)
C        'NPN'        NOAA Profiler Network
C        'ACARS'      Aircraft data (ACARS, AMDARS, etc.)
C        'ACARSP'     Aircraft data profiles at airports
C        'MAP'        Multi-Agency Profilers
C        'HYDRO'      Hydro
C        'RDMTR'      Radiometer
C        'SATWND'     Satellite wind
C        'SATSND'     Satellite soundings
C        'SATRAD'     Satellite radiance
C        'SNOW'       Snow
C        'WISDOM'     Wisdom balloon wind
C
C     2) Name of the database:
C
C        'FSL'        Data obtained from FSL
C        'AWIPS'      Data obtained from AWIPS
C
C        Note that the list of variables with QC found in the MADIS
C        documentation applies to the FSL database.  For the AWIPS
C        database, the only QC is on sfc mesonet data.  Otherwise, the
C        QC information returned by API calls will indicate that no
C        QC has been applied.
C
C     3) Error handling (.true./.false.):
C
C        If ".true.", then error messages will be written to standard
C        output.  If ".false.", no error messages will be output (but
C        ISTATUS will still contain an error code).  The error codes
C        can be found in errors.txt in the doc directory.

C     PAR FILE INPUT:  Database ('FSL' or 'AWIPS')

C     The necessary inputs are found in the parameter file, along with
C     a number of lines of explanation, section delimiters, etc. meant
C     to aid readability.  Throughout this program, we'll be skipping
C     non-input lines, reading necessary inputs, and checking for I/O
C     errors.  Otherwise, the validity of most of user's choices will
C     be tested by the affected MADIS library routines.  Invalid choices
C     will be reported to standard output, and in most cases, the
C     program will then stop.

      line = 0
      do i = 1, 3
         line = line + 1
         read(2,*,err=825,end=850)
      enddo
      line = line + 1
      read(2,145,err=825,end=850)db
 145  format(a)

      write(*,500)'MADIS RAOBDUMP Activity Log'
 500  format(1x,a)
      write(*,500)'---------------------------'

      write(*,501)'Database'
 501  format(/,1x,a)
      write(*,401)db
 401  format(3x,a)

      CALL MINIT('RAOB',db,.true.,istatus)
      if(istatus.ne.success_p)stop

      if(db(1:3).eq.'FSL')then
         idb = 0
      else if(db(1:5).eq.'AWIPS')then
         idb = 1
      endif

C-----------------------------------------------------------------------
C OPTIONAL MADIS SETUP CALLS
C -------------------------
C
C     Optional calls to limit which obs MADIS will return.
C     These can be selected at any time in the code,
C     and will then take effect upon the next call to
C     MRAOBSTA (see below).  The default behavior (what
C     happens if you don't call the routine) is described
C     below.
C
C     Select time window
C     ------------------
C
C     The user selects the nominal time of the data to be returned
C     when using MADIS routines to read obs and station information.
C     Time "windowing" is used so that the MADIS routines will
C     know how close in time each station's observation time needs
C     to be to the nominal time in order to be included, and whether
C     or not duplicate records from the same station that are within
C     the time window should be included.  For radiosonde obs, the
C     "observation time" is sometimes referred to as the "synoptic
C     time", as contrasted with the "release time" of the balloon.
C     For the twice-daily North American radiosonde launches, the
C     observation times are 0000 and 1200 UTC (and the release
C     times are typically +- 70 minutes of the observation time).
C
C     The time window is specified with the MSETWIN subroutine, which
C     allows the user to select the number of minutes before and
C     after the nominal time to be used for the window, and how
C     duplicates should be handled.  For more details, see msetwin.txt
C     in the doc directory.
C
C     If MSETWIN is not called, the default will be to return all
C     records (including multiple reports from the same station)
C     where the observation time is during the hour specified in the
C     file's time.  This will return all records in a standard FSL
C     hourly RAOB file, e.g., with a nominal time of 1200 UTC
C     data time-stamped 1200 through 1259 will be returned.
C
C     The following sample MSETWIN call uses values that will return
C     obs within +- 60 minutes of the nominal time, and if multiple
C     reports from the same station are found, only the latest one
C     will be returned.  If your application runs hourly, you might
C     want to specify that when multiple reports are found, the
C     one closest to the nominal time be returned.  For other
C     applications it's suggested that using the latest record
C     within the window will provide better results in cases where
C     a second radiosonde launch has been necessary as a result of
C     equipment problems, or thunderstorms in the area. 
C
C     Because many radiosondes tend to be available with observation
C     times up to an hour earlier than the 0000 and 1200 UTC periods,
C     it's recommended that you don't use the default window (which
C     will only give you the 0000 or 1200 reports), but instead use
C     the +- 60 minute window shown in the sample call.
C
C     Sample MSETWIN call:
C
C     minbck = -60
C     minfwd = 60
C     recwin = 3
C     CALL MSETWIN(minbck,minfwd,recwin,istatus)

C     PAR FILE INPUT:  Time window parameters (minbck, minfwd, recwin)
C                      [set all to 0 for default]
      do i = 1, 3
         line = line + 1
         read(2,*,err=825,end=850)
      enddo
      line = line + 1
      read(2,*,err=825,end=850)minbck
      line = line + 1
      read(2,*,err=825,end=850)minfwd
      line = line + 1
      read(2,*,err=825,end=850)recwin

      write(*,501)'Time Window'

      if(recwin.ge.0.and.recwin.le.4)then

         write(*,401)twstr(recwin)(1:twstrlen(recwin))

         if(minbck.ne.0.or.minfwd.ne.0.or.recwin.ne.0)then

            write(*,502)'Start minute relative to nominal time ',minbck
 502        format(3x,a,i4)
            write(*,502)'End minute relative to nominal time   ',minfwd

            CALL MSETWIN(minbck,minfwd,recwin,istatus)
            if(istatus.ne.success_p)stop

         endif

      else

         write(*,500)'RAOBDUMP: INVALID WINDOW RECORD TYPE'
         stop

      endif

C     Select a horizontal domain
C     --------------------------
C
C     The horizontal domain can be limited to match one of two map
C     projections (polar stereographic, Lambert conformal conic) or it
C     can be limited to a latitude/longitude rectangle.  Only stations
C     within the domain will be returned.  If one of the map projections
C     is selected, any winds (of any form) returned by subsequent
C     MGETRAOB calls will be rotated to match the projection.  See
C     msetdom.txt in the doc directory for an explanation of the 
C     parameters needed to define the domains.
C
C     If MSETDOM is not called, all available stations will be
C     returned, regardless of their location.
C
C     The example given here is for the CONUS 212 Lambert conformal 
C     conic projection.  
C
C     CALL MSETDOM(3,40635.25,12.19,-133.459,-95.0,25.0,25.0,dummy,
C    1             185,129,istatus)

C     PAR FILE INPUT:  Domain Filter (0 - don't filter, 
C                                     1 - lat/lon corners,
C                                     2 - Polar Stereographic grid,
C                                     3 - Lambert Conformal Conic grid)
      do i = 1, 7
         line = line + 1
         read(2,*,err=825,end=850)
      enddo

      line = line + 1
      read(2,*,err=825,end=850)idmtyp

C     We'll read in the parameters appropriate for the selected
C     domain filter, and skip over all of the other lines in the
C     rest of this section of the parameter file.

      write(*,501)'Domain Filter'

      if(idmtyp.eq.1)then

         do i = 1, 6
            line = line + 1
            read(2,*,err=825,end=850)
         enddo
         do i = 1, 4
            line = line + 1
            read(2,*,err=825,end=850)rpar(i)
         enddo
         do i = 1, 23
            line = line + 1
            read(2,*,err=825,end=850)
         enddo

         write(*,401)'Latitude/Longitude Corners'
         write(*,402)'SW corner latitude (north)                 ',
     1               rpar(1)
 402     format(3x,a,f13.6)
         write(*,402)'SW corner longitude (east)                 ',
     1               rpar(2)
         write(*,402)'NE corner latitude (north)                 ',
     1               rpar(3)
         write(*,402)'NE corner longitude (east)                 ',
     1               rpar(4)

         CALL MSETDOM(idmtyp,rpar(1),rpar(2),rpar(3),rpar(4),dummy,
     1                dummy,dummy,idummy,idummy,istatus)

      else if(idmtyp.eq.2)then

         do i = 1, 13
            line = line + 1
            read(2,*,err=825,end=850)
         enddo
         do i = 1, 7
            line = line + 1
            read(2,*,err=825,end=850)rpar(i)
         enddo
         do i = 1, 2
            line = line + 1
            read(2,*,err=825,end=850)ipar(i)
         enddo
         do i = 1, 11
            line = line + 1
            read(2,*,err=825,end=850)
         enddo

         write(*,401)'Polar Stereographic Grid'
         write(*,402)'Grid box size (meters)                     ',
     1               rpar(1)
         write(*,402)'Latitude (north) of 1st grid point         ',
     1               rpar(2)
         write(*,402)'Longitude (east) of 1st grid point         ',
     1               rpar(3)
         write(*,402)'Orientation longitude (east)               ',
     1               rpar(4)
         write(*,402)'I-coordinate of pole                       ',
     1               rpar(5)
         write(*,402)'J-coordinate of pole                       ',
     1               rpar(6)
         write(*,402)'Latitude (north) at which X-Y scale is true',
     1               rpar(7)
         write(*,403)'Number of grid points (X)                  ',
     1               ipar(1)
         write(*,403)'Number of grid points (Y)                  ',
     1               ipar(2)
 403     format(3x,a,i6)   

         CALL MSETDOM(idmtyp,rpar(1),rpar(2),rpar(3),rpar(4),rpar(5),
     1                rpar(6),rpar(7),ipar(1),ipar(2),istatus)

      else if(idmtyp.eq.3)then

         do i = 1, 25
            line = line + 1
            read(2,*,err=825,end=850)
         enddo
         do i = 1, 6
            line = line + 1
            read(2,*,err=825,end=850)rpar(i)
         enddo
         do i = 1, 2
            line = line + 1
            read(2,*,err=825,end=850)ipar(i)
         enddo

         write(*,401)'Lambert Conformal Conic Grid'
         write(*,402)'Grid box size (meters)                     ',
     1               rpar(1)
         write(*,402)'Latitude (north) of 1st grid point         ',
     1               rpar(2)
         write(*,402)'Longitude (east) of 1st grid point         ',
     1               rpar(3)
         write(*,402)'Orientation longitude (east)               ',
     1               rpar(4)
         write(*,402)'Latitude (north) of 1st tangent            ',
     1               rpar(5)
         write(*,402)'Latitude (north) of 2nd tangent            ',
     1               rpar(6)
         write(*,403)'Number of grid points (X)                  ',
     1               ipar(1)
         write(*,403)'Number of grid points (Y)                  ',
     1               ipar(2)

         CALL MSETDOM(idmtyp,rpar(1),rpar(2),rpar(3),rpar(4),rpar(5),
     1                rpar(6),dummy,ipar(1),ipar(2),istatus)

      else if(idmtyp.eq.0)then

         write(*,401)'No Filter'

         do i = 1, 33
            line = line + 1
            read(2,*,err=825,end=850)
         enddo
         
      else

         write(*,500)'RAOBDUMP: INVALID DOMAIN TYPE'
         stop

      endif
      if(istatus.ne.success_p)stop

C     Select QC level of data to be returned
C     --------------------------------------
C
C     The primary variables in all MADIS datasets have been screened
C     by automated quality control (QC) checks.  Different variables
C     have different levels of QC, or "stages", applied to them.
C     For example, temperature may have level 1 checks (validity) and
C     level 2 checks (hydrostatic, superadiabatic lapse rate).  Note
C     that each level of QC incorporates the previous levels, i.e., an
C     ob that has passed level 2 has also passed any applicable level
C     1 checks.  See raob_qc_variable_list.txt in the doc directory
C     for more information.
C
C     Each call to MGETRAOB will return obs and associated QC
C     information.  If the user doesn't want to interpret the QC info,
C     he or she can limit the returned obs to selectable levels of
C     quality.  In that case, obs that don't meet the selected QC level
C     will be replaced with missing data flags, and the QC information
C     itself doesn't need to be examined by the calling program.
C
C     If MSETQC is not called, all obs will be returned (along with QC
C     info).
C
C     There are three ways that the user can request and process the QC 
C     information:
C
C     1) Don't call MSETQC at all, and examine the QC info before
C        deciding whether or not to use each observation.  This
C        allows the calling program the maximum amount of control
C        over what's available.
C
C     2) If you know ahead of time that the same level of QC is
C        available for all variables you'll be using, call MSETQC
C        one time to select the desired QC level.  Note that the
C        level selected by MSETQC will stay in effect until the
C        next time MSETQC is called.
C
C     3) To require that all variables pass the highest level
C        of QC that are specified for them, call MSETQC one time
C        and select QC level 99.
C
C     User's of MADIS data who want the data filtered by QC results,
C     but don't want to interpret the QC information themselves, should
C     call MSETQC one time and select QC level 99.

C     PAR FILE INPUT:  QC level (0-none, 1-coarse, 2-screened, 
C                                99-highest possible)
      do i = 1, 3
         line = line + 1
         read(2,*,err=825,end=850)
      enddo
      line = line + 1
      read(2,*,err=825,end=850)qctyp

      write(*,501)'QC Filter'
      if(db(1:3).eq.'FSL')then
         if(qctyp.ge.0.and.qctyp.le.3)then
            write(*,401)qcstr(qctyp)
         else if(qctyp.eq.99)then
            write(*,401)'Highest Possible'
         else
            write(*,500)'RAOBDUMP: INVALID QC TYPE'
            stop
         endif
      else if(db(1:5).eq.'AWIPS')then
         write(*,401)'AWIPS database - QC not available'
         qctyp = 0
      endif

C     Select the QC level to be returned.

      CALL MSETQC(qctyp,istatus)
      if(istatus.ne.success_p)stop

C     Select Correction Options
C     -------------------------
C
C     Mandatory level heights and temperatures have a hydrostatic
C     check applied to them that, given certain circumstances, can
C     have sufficient confidence to recommend a corrected value
C     be used instead of the original value, which has failed the
C     check.  The actual data value, however, isn't changed in the
C     MADIS database, and the user has control over whether or not
C     the correction is made before the data are returned to the
C     user.  Note that this choice also applies to variables that
C     are derived in MGETRAOB (e.g., interpolated pressure and
C     height, dew point temperature and relative humidity).
C
C     If MSETCOR is not called, the default behavior is to apply
C     any recommended hydrostatic check corrections.

      write(*,501)'Hydrostatic Correction'

C     PAR FILE INPUT:  Correct mandatory H, T (0-yes, 1-no)

      do i = 1, 6
         line = line + 1
         read(2,*,err=825,end=850)
      enddo
      line = line + 1
      read(2,*,err=825,end=850)corr_type
      if(corr_type.eq.0)then
         write(*,401)'Yes'
      else if(corr_type.eq.1)then
         write(*,401)'No'
      endif

      CALL MSETCOR(corr_type,istatus)
      if(istatus.ne.success_p)stop

C-----------------------------------------------------------------------
C SELECT NOMINAL TIME OF DATA TO BE RETURNED 
C ------------------------------------------
C
C     The MADIS library's native time zone and format is Universal 
C     Coordinated Time, using a Julian date, in the form of "YYJJJHHMM", 
C     where:
C
C     YY  = 2-digit year (good from 1980 - 2179)
C     JJJ = Julian date of the year (001-366)
C     HH  = Hour (00-23)
C     MM  = Minute (00-59)
C
C     MADIS also provides support for a "Month/Day" date format, using
C     "YYYYMMDD_HHMM", where:
C
C     YYYY = 4-digit year (still only good from 1980 - 2179)
C     MMDD = Month/Day 
C     HHMM = Hour/Minute 
C
C     The MADIS library includes a basic set of routines that the caller
C     can use to:
C
C     - Add & subtract time in the native format.  (MCHGTIM)
C     - Translate the native format to integers representing year, 
C       month, day, hour, and minute.  (MINTTIM)
C     - Translate integers representing year, month, day, hour and
C       minute into the native format.  (MCHRTIM)
C     - Translate between Julian and Month/Day formats (MTRNTIM)
C
C     See time_routines.txt in the doc directory for details on how 
C     to use these calls.
C
C     Also note that nominal times can be specified at any minute, not
C     just at the top of the hour (minute 00).

C     PAR FILE INPUT:  Time format to use (0 - Julian, 1 - Month/Day)
      
      write(*,501)'Nominal Time'

      do i = 1, 4
         line = line + 1
         read(2,*,err=825,end=850)
      enddo
      line = line + 1
      read(2,*,err=825,end=850)tform
      if(tform.ne.0.and.tform.ne.1)then
         write(*,500)'RAOBDUMP: INVALID TIME TRANSLATION FORMAT'
         stop
      endif

C     PAR FILE INPUT:  Nominal time (0 - specified in par file,
C                                    1 - use current time)
      line = line + 1
      read(2,*,err=825,end=850)
      line = line + 1
      read(2,*,err=825,end=850)ctime
      line = line + 1
      read(2,*,err=825,end=850)

C     PAR FILE INPUT:  Nominal time (in selected format)

      line = line + 1
      read(2,144,err=825,end=850)tstr1
 144  format(a)

C     If the user wants to use the current time, get it and translate 
C     that into the chosen format (via our internal formats).  
C     Otherwise we have the chosen format in TSTR1 right now.

      if(ctime.eq.1)then
         i = TIME()
         i = i - 1325376000
         CALL M_MKFNAM(i,atime,istatus)
         if(tform.eq.1)then
            CALL MTRNTIM(atime,2,tstr1,istatus)
         else
            tstr1 = atime
         endif
      endif

      if(tform.eq.1)then
         CALL MTRNTIM(tstr1,1,atime,istatus)
         if(istatus.ne.success_p)stop
      else
         atime = tstr1
      endif

      write(*,401)tstr1

C     For most purposes, you'll want to read data for all stations
C     that meet the criteria already specified above.  For dump program
C     purposes, however, we offer the option to only get one station.

      write(*,501)'Station Selection'

C     PAR FILE INPUT:  Station Selection (0 - all stations, 1 - single)

      do i = 1, 3
         line = line + 1
         read(2,*,err=825,end=850)
      enddo

      line = line + 1
      read(2,*,err=825,end=850)ising
      line = line + 1
      read(2,*,err=825,end=850)
      line = line + 1
      read(2,144,err=825,end=850)singsta
      
      if(ising.eq.1)then
         write(*,401)singsta
      else
         write(*,401)'All stations'
      endif

      do i = 1, 3
         line = line + 1
         read(2,*,err=825,end=850)
      enddo

C-----------------------------------------------------------------------
C LOAD IN STATION INFO FOR SELECTED DATA
C --------------------------------------
C
C     Load in the stations covering the selected domain.  The MRAOBSTA
C     call needs to be done as the first call when getting data for a
C     new time, or whenever any of the optional setup calls are used.
C     Thereafter, MGETRAOB and MGETRAOBC can be called as many times
C     as desired in order to read variables.

      CALL MRAOBSTA(atime,nsta,stanam,wmoid,lat,lon,elev,timeob,
     1              istatus)
      if(istatus.ne.success_p)stop

C     The MGETRAOB routine is used to read all numeric variables.
C     All variables are returned as Fortran REAL type (RDATA), the
C     QC data descriptors are always CHARACTER*1 (QCD), and the
C     bit maps showing the QC checks that were applied (QCA) and
C     their results (QCR) are always INTEGER.
C
C     The returned arrays will match the records in the station info 
C     arrays returned from MRAOBSTA.  If a QC level has been selected, 
C     the obs will only be returned if they meet or exceed the selected 
C     QC level.

C     Read the level type variable.  For the FSL database, also read
C     the sonde type and release time variables.

      CALL MGETRAOB(atime,level_type_vcode,nsta,nobs,nlevels,ltype,
     1              qcd,qca,qcr,istatus)
      if(istatus.ne.success_p)stop

      if(db(1:3).eq.'FSL')then

         CALL MGETRAOB(atime,sonde_type_vcode,nsta,nobs,nlevels,stype,
     1                 qcd,qca,qcr,istatus)

C     The MGETRAOBC routine is used to read character variables.  No QC
C     is currently available for any of these variables.  The only
C     character variable is the release time, which is returned in
C     the standard MADIS 9-character time string.  Note that there are
C     no missing data flags for character data, and therefore no
C     equivalent of the NOBS argument returned from MGETRAOB; instead,
C     blanks are filled in where no genuine data exist.

         CALL MGETRAOBC(atime,release_time_vcode,nsta,reltime,istatus)

C     (We didn't check the status on the calls to get the sonde type
C     and release time because these aren't critical variables, and
C     it's conceivable that they might not exist if you're only getting 
C     data for one station.)  

      endif


C     If one of the map projections has been selected for the horizontal
C     domain, the station latitude and longitude can be converted to a 
C     floating point grid index into the 2-dimensional grid (RI,RJ).  
C     Note that this works for grids laid out this way:
C
C                        NE (I=NX,J=NY)
C
C     SW (I=1,J=1)
C
C     If we're to convert LAT/LON to RI,RJ for the user's grid, do
C     so.  Otherwise, fill in 0's so there'll be something in the
C     output line that isn't too messy.

      if(idmtyp.eq.2.or.idmtyp.eq.3)then

         do i = 1, nsta
            CALL MGETIJ(lat(i),lon(i),ri(i),rj(i),istatus)
         enddo
         
      else

         do i = 1, nsta
            ri(i) = 0
            rj(i) = 0
         enddo

      endif

C-----------------------------------------------------------------------
C LOOP THROUGH DESIRED VARIABLES:
C -------------------------------
C
C     Variables are specified by code names.  The list of code names
C     can be found in the first column of raob_variable_list.txt in the
C     doc directory.  Separate variable names are used for the different
C     sounding types, e.g., "P" will return the integrated sounding
C     pressures, and "PMAN" will return the mandatory level pressures.
C
C     For most applications, the user should use the set of variable 
C     code names that will retrieve the integrated sounding.
C
C     This program will read in up to 4 variable codes per line in the
C     parameter file, read them, output them, then move on to the next
C     set of variables until the end of the parameter file is reached.

      write(*,501)'Variables'

C     PAR FILE INPUT:  1 to 4 variable codes (free-formatted)

 100  read(2,145,end=900)vline

C     Parse the line to determine the number of variables and their
C     code names (1 to 4 code names separated by spaces).

      nvar = 0
      l = -1
      do i =  1, 4
         j = l
         CALL M_FINDWORD(vline,70,j,k,l)
         if(k.eq.0)then
            go to 101
         else
            nvar = nvar + 1
            vcode(nvar) = vline(k:l)
         endif
      enddo

 101  if(nvar.eq.0) go to 100

C     Use internal MADIS routines (MGETVCN, MGENVCN) to translate the 
C     variable code name into the MADIS variable code number, then use 
C     that to determine the sounding type.  This is needed here because
C     this program is an all-purpose utility meant to handle any desired
C     sounding type.  In general, a user's MADIS program shouldn't have 
C     to determine the sounding type (as you would know what variables 
C     you're reading).
C
C        Sounding Type
C        -------------------
C        0        Integrated
C        2        Mandatory
C        3        Significant
C        4        Wind
C        5        Max wind
C        6        Tropopause
     
      do i = 1, nvar
         CALL MGETVCN(vcode(i),vcn,cdummy,istatus)
         if(istatus.ne.success_p)stop
         CALL MGENVCN(vcn,raob_ds_p,idb,dummy,sounding_type(i))
      enddo

C     Make sure the user didn't select a mismatch by mixing variables
C     with different sounding types on the same input line.  (This
C     is against the rules in this program because we're expecting
C     the same number of levels for all variables output on the same
C     line -- one output line per level will result.)

      j = sounding_type(1)
      do i = 2, nvar
         if(sounding_type(i).ne.j)then
            write(*,519)'Variable/sounding type mismatch: ',
     1                  (vcode(k),k=1,nvar)
 519        format(1x,a,4(1x,a))
            go to 100
         endif
      enddo

      do 150 i = 1, nvar

         write(*,505)vcode(i)
 505     format(3x,'Processing',5x,a)

C     Read the desired variable and QC info.  

         if(sounding_type(1).eq.0)then

            CALL MGETRAOB(atime,vcode(i),nsta,nobs,nlevels,rdata(1,1,i),
     1                    qcd(1,1,i),qca(1,1,i),qcr(1,1,i),istatus)

         else if(sounding_type(1).eq.2)then

            CALL MGETRAOB(atime,vcode(i),nsta,nobs,nlevels,
     1                    rdata_man(1,1,i),qcd_man(1,1,i),
     2                    qca_man(1,1,i),qcr_man(1,1,i),istatus)

         else if(sounding_type(1).eq.3)then

            CALL MGETRAOB(atime,vcode(i),nsta,nobs,nlevels,
     1                    rdata_sigt(1,1,i),qcd_sigt(1,1,i),
     2                    qca_sigt(1,1,i),qcr_sigt(1,1,i),istatus)


         else if(sounding_type(1).eq.4)then

            CALL MGETRAOB(atime,vcode(i),nsta,nobs,nlevels,
     1                    rdata_sigw(1,1,i),qcd_sigw(1,1,i),
     2                    qca_sigw(1,1,i),qcr_sigw(1,1,i),istatus)

         else if(sounding_type(1).eq.5)then

            CALL MGETRAOB(atime,vcode(i),nsta,nobs,nlevels,
     1                    rdata_mxwnd(1,1,i),qcd_mxwnd(1,1,i),
     2                    qca_mxwnd(1,1,i),qcr_mxwnd(1,1,i),istatus)

         else if(sounding_type(1).eq.6)then

            CALL MGETRAOB(atime,vcode(i),nsta,nobs,nlevels,
     1                    rdata_trop(1,1,i),qcd_trop(1,1,i),
     2                    qca_trop(1,1,i),qcr_trop(1,1,i),istatus)

         endif

         write(1,146)'var ',vcode(i),' total stns ',nsta,
     1               ' # non-missing obs ',nobs
 146     format(/,1x,a,a,a,i5,a,i5)

         if(istatus.ne.success_p)go to 150

C     Put the variable code name into an array that will be used as
C     a label in header lines in the output file.

         vcstr(i) = vcode(i)//' QCD QCA QCR'

 150  continue

C     Output the station info, variables and QC.

      write(1,*)


      do 200 i = 1, nsta

C     If we're getting only one station, make sure this is it.

         if(ising.eq.1.and.singsta.ne.stanam(i))go to 200

C     Use the time output format desired by the user.

         if(tform.eq.0)then

            tstr1 = '  '//timeob(i)
            tstr2 = '  '//reltime(i)

         else if(tform.eq.1)then

            CALL MTRNTIM(timeob(i),2,tstr1,istatus)
            CALL MTRNTIM(reltime(i),2,tstr2,istatus)

         endif

         if(nlevels(i).ne.0.and.db(1:3).eq.'FSL')then

            write(1,1000)stanam(i),wmoid(i),elev(i),lat(i),lon(i),ri(i),
     1                   rj(i),nlevels(i),tstr1,tstr2,nint(stype(i)),
     2                   (vcstr(j),j=1,nvar)
 1000       format(1x,'Station     WMOID  Elev(m) Lat(N) Lon(E)   ',
     1             'Grid I   Grid J  NumLev    ObTime      ',
     2             'ReleaseTime   SondeType',/,
     3             1x,a,1x,i6,1x,f7.2,1x,f6.2,1x,f7.2,1x,
     4             2(f8.3,1x),2x,i3,3x,a,1x,a,1x,i6,/,/,
     6             1x,'Lev/Type  ',4(a,5x))

         else if(nlevels(i).ne.0.and.db(1:5).eq.'AWIPS')then

            write(1,1001)stanam(i),wmoid(i),elev(i),lat(i),lon(i),ri(i),
     1                   rj(i),nlevels(i),tstr1,(vcstr(j),j=1,nvar)
 1001       format(1x,'Station     WMOID  Elev(m) Lat(N) Lon(E)   ',
     1             'Grid I   Grid J  NumLev    ObTime',/,
     3             1x,a,1x,i6,1x,f7.2,1x,f6.2,1x,f7.2,1x,
     4             2(f8.3,1x),2x,i3,3x,a,1x,/,/,
     6             1x,'Lev/Type  ',4(a,5x))

         endif

         do j = 1, nlevels(i)

            if(sounding_type(1).eq.0)then

               write(1,2000)j,ltstr(int(ltype(j,i))),(rdata(j,i,k),
     1                      qcd(j,i,k),qca(j,i,k),qcr(j,i,k),k=1,nvar)
 2000          format(1x,i3,'/',a,4(1x,f13.6,1x,a,1x,i3,1x,i3))

            else if(sounding_type(1).eq.2)then

               write(1,2000)j,ltstr(sounding_type(1)),
     1                      (rdata_man(j,i,k),qcd_man(j,i,k),
     2                       qca_man(j,i,k),qcr_man(j,i,k),k=1,nvar)

            else if(sounding_type(1).eq.3)then

               write(1,2000)j,ltstr(sounding_type(1)),
     1                      (rdata_sigt(j,i,k),qcd_sigt(j,i,k),
     2                       qca_sigt(j,i,k),qcr_sigt(j,i,k),k=1,nvar)

            else if(sounding_type(1).eq.4)then

               write(1,2000)j,ltstr(sounding_type(1)),
     1                      (rdata_sigw(j,i,k),qcd_sigw(j,i,k),
     2                       qca_sigw(j,i,k),qcr_sigw(j,i,k),k=1,nvar)

            else if(sounding_type(1).eq.5)then

               write(1,2000)j,ltstr(sounding_type(1)),
     1                      (rdata_mxwnd(j,i,k),qcd_mxwnd(j,i,k),
     2                       qca_mxwnd(j,i,k),qcr_mxwnd(j,i,k),k=1,nvar)

            else if(sounding_type(1).eq.6)then

               write(1,2000)j,ltstr(sounding_type(1)),
     1                      (rdata_trop(j,i,k),qcd_trop(j,i,k),
     2                       qca_trop(j,i,k),qcr_trop(j,i,k),k=1,nvar)
            endif

         enddo

         if(i.ne.nsta.and.nlevels(i).ne.0)write(1,*)

 200  continue

C     Process the next line of variables.

      go to 100

C-----------------------------------------------------------------------
C PROGRAM EXIT
C ------------
C
C     Error exits.

 825  write(*,503)'RAOBDUMP: READ ERROR IN PARAMETER FILE AT LINE ',line
 503  format(1x,a,i4)
      stop

 850  write(*,500)'RAOBDUMP: PREMATURE END OF PARAMETER FILE'

C     Normal exit.

 900  stop
      end
C**---------------------------------------------------------------------
C**REQUIRED STANDARD FSL DISCLAIMER (NOVEMBER 2000)
C**---------------------------------------------------------------------
C**  OPEN SOURCE LICENSE/DISCLAIMER,  FORECAST SYSTEMS LABORATORY
C**  NOAA/OAR/FSL, 325 BROADWAY BOULDER, CO 80305
C**
C**  THIS SOFTWARE IS DISTRIBUTED UNDER THE OPEN SOURCE DEFINITION, 
C**  WHICH MAY BE FOUND AT http://www.opensource.org/osd.html.
C**
C**  IN PARTICULAR, REDISTRIBUTION AND USE IN SOURCE AND BINARY FORMS,
C**  WITH OR WITHOUT MODIFICATION, ARE PERMITTED PROVIDED THAT THE 
C**  FOLLOWING CONDITIONS ARE MET:
C**
C**  - REDISTRIBUTIONS OF SOURCE CODE MUST RETAIN THIS NOTICE, THIS LIST
C**    OF CONDITIONS AND THE FOLLOWING DISCLAIMER.
C**
C**  - REDISTRIBUTIONS IN BINARY FORM MUST PROVIDE ACCESS TO THIS 
C**    NOTICE, THIS LIST OF CONDITIONS AND THE FOLLOWING DISCLAIMER, AND
C**    THE UNDERLYING SOURCE CODE.
C**
C**  - ALL MODIFICATIONS TO THIS SOFTWARE MUST BE CLEARLY DOCUMENTED, 
C**    AND ARE SOLELY THE RESPONSIBILITY OF THE AGENT MAKING THE 
C**    MODIFICATIONS.
C**
C**  - IF SIGNIFICANT MODIFICATIONS OR ENHANCEMENTS ARE MADE TO THIS
C**    SOFTWARE, THE FSL SOFTWARE POLICY MANAGER 
C**    (softwaremgr.fsl@noaa.gov) BE NOTIFIED.
C**
C**  THIS SOFTWARE AND ITS DOCUMENTATION ARE IN THE PUBLIC DOMAIN AND 
C**  ARE FURNISHED "AS IS."  THE AUTHORS, THE UNITED STATES GOVERNMENT,
C**  ITS INSTRUMENTALITIES, OFFICERS, EMPLOYEES, AND AGENTS MAKE NO 
C**  WARRANTY, EXPRESS OR IMPLIED, AS TO THE USEFULNESS OF THE SOFTWARE
C**  AND DOCUMENTATION FOR ANY PURPOSE.  THEY ASSUME NO RESPONSIBILITY 
C**  (1) FOR THE USE OF THE SOFTWARE AND DOCUMENTATION; OR (2) TO 
C**  PROVIDE TECHNICAL SUPPORT TO USERS.
C**---------------------------------------------------------------------
