;-------------------------------------------------------------
;+
; NAME:
;        READ_DAO
;
; PURPOSE:
;        Reads DAO I-6,(instantaneous 6h), A-6 (average 6h),
;        or A-3 (average 3-h) met fields from disk.
;
; CATEGORY:
;        CTM Tools
;
; CALLING SEQUENCE:
;        READ_DAO, FileName, Data [, Keywords ]
;
; INPUTS:
;        FILENAME  -> The name of the input file (or a file mask).
;             FILENAME is passed to OPEN_FILE.  If FILENAME is a null 
;             string or a file mask, then OPEN_FILE will open a
;             pickfile dialog box.
;
; KEYWORD PARAMETERS:
;        THISFILENAME -> Returns to the calling program the name of
;             the file found by OPEN_FILE.
;
;        SELECTION -> name of input field to be read (default UWND)
;
;        TIME -> time to be read (default 0000 UTC)
;
;        /RES2 -> Set array dimensions for 2 x 2.5 resolution
;             fields.  Default is for 4 x 5 resolution fields.
;
;        VERBOSE -> Will print information about which met fields
;             are being read in.
;
;        NLAYERS -> Specifies the number of vertical layers for the
;             DATA array.  Default is 20 layers (for GEOS-1).
;
; OUTPUTS:
;        DATA -> Will return data array which will be 2D or 3D
;             depending on the field selected.
;
; SUBROUTINES:
;        OPEN_FILE
;        GET_FREELUN (function)
;
; REQUIREMENTS:
;        None
;
; NOTES:
;   Major changes due to different purpose): 
;   - pass filename instead of logical unit
;          (on terra: /data/ctm/GEOS_4x5/GEOS_1/yymmdd.i6.4x5 [.gz] )
;
;   - return only selected information (1 array)
;
;=============================================================================
;   Excerpt from original FORTRAN programs:
;
;   I-3 Fields for GEOS-1 and GEOS-STRAT data:
;   -------------------------------------------------------------------------
;   (1) PS     : DAO field for surface pressure            (      mb      )
;   (2) ALBD   : DAO field for surface albedo              (  0.0 - 1.0   )
;   (3) LWI    : DAO field for surface land/water indices  (  1, 2, 3, 4  )
;   (4) UWND   : DAO field for 3-D U-wind field            (     m/s      )
;   (5) VWND   : DAO field for 3-D V-wind field            (     m/s      )
;   (6) TMPU   : DAO field for 3-D temperature field       (      K       )
;   (7) SPHU   : DAO field for 3-D specific humidity field ( g H20/kg air )
;   (8) PHIS   : DAO field for surface geopotential heights (  meters * g ) 
;
;   A-6 Fields for GEOS-1 and GEOS-STRAT data:
;   -------------------------------------------------------------------------
;   (1) MOISTQ  : DAO field for change in specific humidity    ( kg H20
;                   due to moist processes                       kg air / s )
;   (2) CLMOLW  : DAO field for maximum overlap cloud fraction (    0 - 1   )
;   (3) CLROLW  : DAO field for random  overlap cloud fraction (    0 - 1   )
;   (4) CLDMAS  : DAO field for cloud mass flux                (   kg/m2/s  )
;   (5) DTRAIN  : DAO field for cloud detrainment              (   kg/m2/s  )
;
;   A-3 fields for GEOS-1 oNLy:
;   ---------------------------
;   (1 ) HFLUX  : DAO field for sensible heat flux from surface     (  W/m2  )
;   (2 ) RADSWG : DAO field for incident solar radiation @ ground   (  W/m2  )
;   (3 ) PREACC : DAO field for total accum. precip. @ ground       ( mm/day )
;   (4 ) PRECON : DAO field for total conv.  precip. @ ground       ( mm/day )
;   (5 ) TS     : DAO field for surface air temperature             (   K    )
;   (6 ) RADSWT : DAO field for incident solar radiation @ atm. top (  W/m2  )
;   (7 ) USTAR  : DAO field for friction velocity                   (  m/s   )
;   (8 ) Z0     : DAO field for surface roughness height            (   m    )
;   (9 ) PBL    : DAO field for planetary boundary layer depth      (   mb   )
;   (10) CLDFRC : DAO field for 2-D cloud fraction                  ( 0 - 1  )
;   (11) U10M   : DAO field for U-wind speed at 10 meters altitude  (  m/s   )
;   (12) V10M   : DAO field for V-wind speed at 10 meters altitude  (  m/s   )
;
;   A-6 fields for GEOS-STRAT oNLy:
;   ---------------------------
;   (1 ) HFLUX  : DAO field for sensible heat flux from surface     (  W/m2  )
;   (2 ) RADSWG : DAO field for incident solar radiation @ ground   (  W/m2  )
;   (3 ) PREACC : DAO field for total accum. precip. @ ground       ( mm/day )
;   (4 ) PRECON : DAO field for total conv.  precip. @ ground       ( mm/day )
;   (5 ) TS     : DAO field for surface air temperature             (   K    )
;   (6 ) RADSWT : DAO field for incident solar radiation @ atm. top (  W/m2  )
;   (7 ) USTAR  : DAO field for friction velocity                   (  m/s   )
;   (8 ) Z0     : DAO field for surface roughness height            (   m    )
;   (9 ) PBL    : DAO field for planetary boundary layer depth      (  mb    )
;
;  ADDITIONAL NOTES:
;  (1) The I-6 fields are either INSTANTANEOUS 2-D or INSTANTANEOUS 3-D 
;      fields.  They are saved on times 00, 06, 12, and 18h GMT.  In other 
;      words, PS at 06 is the instantaneous value of PS at 6h, etc, etc.
;
;  (2) All of the A-6 fields are AVERAGE 3-D fields.  They are centered 
;      on times 00, 06, 12, and 18h GMT.  In other words, MOISTQ at 06h 
;      contains the average value of MOISTQ from 03h to 09h, etc, etc. 
;
;  (3) All of the A-3 fields are AVERAGE 2-D surface fields.  They are ended 
;      on times 00, 03, 06, 09, 12, 15, 18, and 21h GMT.  In other words, 
;      HFLUX at 3h contains the average HFLUX values from 00-03h GMT.
;
;  (4) GEOS-STRAT has no A-3 fields.  Several of the surface fields
;      that are A-3 in GEOS-1 are A-6 in GEOS-STRAT.
;=============================================================================
;
; EXAMPLES:
;        (1) Read 20-Layer GEOS-1 CLDMAS data 
;            for date 1994/09/01 and 06:00:00 hours UTC
;
;        FileName = '~/terra/CTM4/run_code/TEMP/940901.a6.4x5'
;        Read_DAO, FileName, Data, Selection='CLDMAS', Time=060000
;
;
;        (2) Read 46-layer GEOS-STRAT HFLUX data 
;            for date 1996/09/01 and 00:00:00 hours UTC
;
;        FileName = '~/terra/CTM4/run_code/TEMP/960901.a6.4x5'
;        Read_DAO, FileName, Data, Selection='HFLUX', NLayers=46
;        
;
;        (3) Read 26-layer GEOS-STRAT TMPU data 
;            for 1996/09/01 and 18:00:00 hours UTC
;
;        FileName = '~/terra/CTM4/run_code/TEMP/960901.a6.4x5'
;        Read_DAO, FileName, Data, Selection='TMPU', NLayers=26, Time=180000
;
;
; MODIFICATION HISTORY:
;        mgs, 06 Feb 1998: VERSION 1.00
;                          - translated into IDL by mgs, 06 FEB 1998
;        bmy, 16 Oct 1998: VERSION 2.00
;                          - Added NLAYERS keywords so that 46- or
;                            26-layer GEOS-STRAT data can be read in.
;        bmy, 20 Oct 1998: - Added call to GET_FREELUN
;        bmy, 21 Oct 1998: - TIME can now be passed as 00, 03, 09, 12, 
;                            or as 000000, 060000, 090000, 120000, etc.
;        bmy, 25 Nov 1998: - now returns if the user selects CANCEL
;                            in a dialog pickfile box
;                          - added THISFILENAME keyword
;                          - updated comments
;-
; Copyright (C) 1998, Martin Schultz and Bob Yantosca, 
; Harvard University
; This software is provided as is without any warranty
; whatsoever. It may be freely used, copied or distributed
; for non-commercial purposes. This copyright notice must be
; kept with any copy of this software. If this software shall
; be used commercially or sold as part of a larger package,
; please contact the author to arrange payment.
; Bugs and comments should be directed to mgs@io.harvard.edu
; or bmy@io.harvard.eduwith subject "IDL routine read_a6"
;-------------------------------------------------------------


pro Read_DAO, FileName, Data,                        $
              ThisFileName=ThisFileName,             $
              Selection=Selection,  Time=Time,       $
              Res2=Res2,            Verbose=Verbose, $
              NLayers=NLayers,      NYMD=NYMD,       $
              NHMS=NHMS,            _EXTRA=e

   ; Pass external functions
   FORWARD_FUNCTION Get_FreeLun

   ; Default value for NLAYERS: 20 for GEOS-1, 46 for GEOS-STRAT
   if ( N_Elements( NLayers ) eq 0 ) then NLayers = 20
    
   ; Define necessary parameters
   Name = BytArr(8)             ; FORTRAN char *8
   NI = 72                      ; number of longitude boxes   (2x2.5 = 144)
   NJ = 46                      ; number of latitude boxes    (2x2.5 =  91)
   NL = NLayers                 ; vertical levels
   if ( Keyword_Set( Res2 ) ) then begin
      NI = 144 & NJ = 91 
   end
 
   XYMD = 0.                    ; value for year month day
   XHMS = 0.                    ; value for hour minute second
   Found = 0                    ; flag to indicate that selected field 
                                ;  was Found
 
   ; Set default selection to MOISTQ
   if ( N_Elements( Selection ) le 0 ) then Selection = 'PS'        
   selection = strupcase(selection)
 
   ; Set default time to 000000 UTC
   if ( N_Elements( Time ) le 0 ) then Time = 0. 
 
   ; Times are keyed by 000000, 030000, 060000, etc.
   ; so adjust for this if necessary
   XTime = Time
   if ( XTime lt 10000. ) then XTime = XTime * 10000.

   ; File Found, now open and read
   Open_File, FileName, Get_FreeLun( Ilun ),            $
      /F77_Unformatted, FileName=ThisFileName,          $
      Title='Select a met fields file (CANCEL quits!)', _EXTRA=e

   ; Return if THISFILENAME is a null string
   if ( ThisFileName eq '' ) then return

   ; Read through all of the fields in the file
   while ( not EOF( Ilun ) and not Found ) do begin
      Readu, Ilun, Name
      StrName = StrCompress( Name, /Remove_all )
 
      if ( Keyword_Set( Verbose ) ) then print,'read in label : ',strname
      
      ; STRNAME will be in uppercase, as read in from the file...
      case ( StrName ) of

         ; I-6 fields for GEOS-1 and GEOS-STRAT
         'PS' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end

         'ALBD' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end

         'PHIS' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end

         'LWI' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end

         'UWND' : begin
            Data = FltArr(NI,NJ,NL)
            ReadU,Ilun,XYMD,XHMS,Data
         end

         'VWND' : begin
            Data = FltArr(NI,NJ,NL)
            ReadU,Ilun,XYMD,XHMS,Data
         end

         'TMPU' : begin
            Data = FltArr(NI,NJ,NL)
            ReadU,Ilun,XYMD,XHMS,Data
         end

         'SPHU' : begin
            Data = FltArr(NI,NJ,NL)
            ReadU,Ilun,XYMD,XHMS,Data
         end
         
         ; A-6 Fields for both GEOS-1 and GEOS-STRAT...
         'MOISTQ' : begin
            Data = FltArr(NI,NJ,NL)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         'CLDMAS' : begin
            Data = FltArr(NI,NJ,NL)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         'DTRAIN' : begin
            Data = FltArr(NI,NJ,NL)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         'CLMOLW' : begin
            Data = FltArr(NI,NJ,NL)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         'CLROLW' : begin
            Data = FltArr(NI,NJ,NL)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         ; A-3 fields for GEOS-1
         ; but A-6 fields for GEOS-STRAT
         'HFLUX' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         'RADSWG' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         'PREACC' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         'PRECON' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         'TS' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         'RADSWT' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         'USTAR' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         'Z0' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         'PBL' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end
 
         ; A-3 fields for GEOS-1 only
         'CLDFRC' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end

         'U10M' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end

         'V10M' : begin
            Data = FltArr(NI,NJ)
            ReadU,Ilun,XYMD,XHMS,Data
         end

         else : begin
            Message, 'Invalid selection!', /Continue
            stop
         end
      endcase
 
 
      ; see if information is what we wanted
      ; that means we are done
      if ( StrName eq Selection AND XHMS eq Xtime ) then Found = 1
      
   endwhile
 
   Print,'read (selection, date and time) : ', StrName, XYMD, XHMS
 
   ; return NYMD and NHMS as read from the file as longwords
   NYMD = Long( XYMD )
   NHMS = Long( XHMS )

   Free_LUN, Ilun
 
   return
end