ungrib/src/ngl/g2/getgb2r.f

   1 C-----------------------------------------------------------------------
   2       SUBROUTINE GETGB2R(LUGB,CINDEX,GFLD,IRET)
   3 C$$$  SUBPROGRAM DOCUMENTATION BLOCK
   4 C
   5 C SUBPROGRAM: GETGB2R        READS AND UNPACKS A GRIB FIELD
   6 C   PRGMMR: GILBERT          ORG: W/NP11     DATE: 02-01-15
   7 C
   8 C ABSTRACT: READ AND UNPACK SECTIONS 6 AND 7 FROM A GRIB2 MESSAGE.
   9 C
  10 C   This routine assumes that the "metadata" for this field
  11 C   already exists in derived type gribfield.  Specifically,
  12 C   it requires gfld%ibmap,gfld%ngrdpts,gfld%idrtnum,gfld%idrtmpl,
  13 C   and gfld%ndpts.
  14 C
  15 C   The decoded information for the selected GRIB field
  16 C   is returned in a derived type variable, gfld.
  17 C   Gfld is of type gribfield, which is defined
  18 C   in module grib_mod, so users of this routine will need to include
  19 C   the line "USE GRIB_MOD" in their calling routine.  Each component of the
  20 C   gribfield type is described in the OUTPUT ARGUMENT LIST section below.
  21 C
  22 C PROGRAM HISTORY LOG:
  23 C   95-10-31  IREDELL
  24 C 2002-01-11  GILBERT     MODIFIED FROM GETGB1R TO WORK WITH GRIB2
  25 C
  26 C USAGE:    CALL GETGB2R(LUGB,CINDEX,GFLD,IRET)
  27 C   INPUT ARGUMENTS:
  28 C     LUGB         INTEGER UNIT OF THE UNBLOCKED GRIB DATA FILE
  29 C     CINDEX       INDEX RECORD OF THE GRIB FIELD  ( SEE DOCBLOCK OF
  30 C                  SUBROUTINE IXGB2 FOR DESCRIPTION OF AN INDEX RECORD.)
  31 C   OUTPUT ARGUMENTS:
  32 C     gfld - derived type gribfield ( defined in module grib_mod )
  33 C            ( NOTE: See Remarks Section )
  34 C        gfld%version = GRIB edition number ( currently 2 )
  35 C        gfld%discipline = Message Discipline ( see Code Table 0.0 )
  36 C        gfld%idsect() = Contains the entries in the Identification
  37 C                        Section ( Section 1 )
  38 C                        This element is actually a pointer to an array
  39 C                        that holds the data.
  40 C            gfld%idsect(1)  = Identification of originating Centre
  41 C                                    ( see Common Code Table C-1 )
  42 C                             7 - US National Weather Service
  43 C            gfld%idsect(2)  = Identification of originating Sub-centre
  44 C            gfld%idsect(3)  = GRIB Master Tables Version Number
  45 C                                    ( see Code Table 1.0 )
  46 C                             0 - Experimental
  47 C                             1 - Initial operational version number
  48 C            gfld%idsect(4)  = GRIB Local Tables Version Number
  49 C                                    ( see Code Table 1.1 )
  50 C                             0     - Local tables not used
  51 C                             1-254 - Number of local tables version used
  52 C            gfld%idsect(5)  = Significance of Reference Time (Code Table 1.2)
  53 C                             0 - Analysis
  54 C                             1 - Start of forecast
  55 C                             2 - Verifying time of forecast
  56 C                             3 - Observation time
  57 C            gfld%idsect(6)  = Year ( 4 digits )
  58 C            gfld%idsect(7)  = Month
  59 C            gfld%idsect(8)  = Day
  60 C            gfld%idsect(9)  = Hour
  61 C            gfld%idsect(10)  = Minute
  62 C            gfld%idsect(11)  = Second
  63 C            gfld%idsect(12)  = Production status of processed data
  64 C                                    ( see Code Table 1.3 )
  65 C                              0 - Operational products
  66 C                              1 - Operational test products
  67 C                              2 - Research products
  68 C                              3 - Re-analysis products
  69 C            gfld%idsect(13)  = Type of processed data ( see Code Table 1.4 )
  70 C                              0  - Analysis products
  71 C                              1  - Forecast products
  72 C                              2  - Analysis and forecast products
  73 C                              3  - Control forecast products
  74 C                              4  - Perturbed forecast products
  75 C                              5  - Control and perturbed forecast products
  76 C                              6  - Processed satellite observations
  77 C                              7  - Processed radar observations
  78 C        gfld%idsectlen = Number of elements in gfld%idsect().
  79 C        gfld%local() = Pointer to character array containing contents
  80 C                       of Local Section 2, if included
  81 C        gfld%locallen = length of array gfld%local()
  82 C        gfld%ifldnum = field number within GRIB message
  83 C        gfld%griddef = Source of grid definition (see Code Table 3.0)
  84 C                      0 - Specified in Code table 3.1
  85 C                      1 - Predetermined grid Defined by originating centre
  86 C        gfld%ngrdpts = Number of grid points in the defined grid.
  87 C        gfld%numoct_opt = Number of octets needed for each
  88 C                          additional grid points definition.
  89 C                          Used to define number of
  90 C                          points in each row ( or column ) for
  91 C                          non-regular grids.
  92 C                          = 0, if using regular grid.
  93 C        gfld%interp_opt = Interpretation of list for optional points
  94 C                          definition.  (Code Table 3.11)
  95 C        gfld%igdtnum = Grid Definition Template Number (Code Table 3.1)
  96 C        gfld%igdtmpl() = Contains the data values for the specified Grid
  97 C                         Definition Template ( NN=gfld%igdtnum ).  Each
  98 C                         element of this integer array contains an entry (in
  99 C                         the order specified) of Grid Defintion Template 3.NN
 100 C                         This element is actually a pointer to an array
 101 C                         that holds the data.
 102 C        gfld%igdtlen = Number of elements in gfld%igdtmpl().  i.e. number of
 103 C                       entries in Grid Defintion Template 3.NN
 104 C                       ( NN=gfld%igdtnum ).
 105 C        gfld%list_opt() = (Used if gfld%numoct_opt .ne. 0)  This array
 106 C                          contains the number of grid points contained in
 107 C                          each row ( or column ).  (part of Section 3)
 108 C                          This element is actually a pointer to an array
 109 C                          that holds the data.  This pointer is nullified
 110 C                          if gfld%numoct_opt=0.
 111 C        gfld%num_opt = (Used if gfld%numoct_opt .ne. 0)  The number of entries
 112 C                       in array ideflist.  i.e. number of rows ( or columns )
 113 C                       for which optional grid points are defined.  This value
 114 C                       is set to zero, if gfld%numoct_opt=0.
 115 C        gfdl%ipdtnum = Product Definition Template Number (see Code Table 4.0)
 116 C        gfld%ipdtmpl() = Contains the data values for the specified Product
 117 C                         Definition Template ( N=gfdl%ipdtnum ).  Each element
 118 C                         of this integer array contains an entry (in the
 119 C                         order specified) of Product Defintion Template 4.N.
 120 C                         This element is actually a pointer to an array
 121 C                         that holds the data.
 122 C        gfld%ipdtlen = Number of elements in gfld%ipdtmpl().  i.e. number of
 123 C                       entries in Product Defintion Template 4.N
 124 C                       ( N=gfdl%ipdtnum ).
 125 C        gfld%coord_list() = Real array containing floating point values
 126 C                            intended to document the vertical discretisation
 127 C                            associated to model data on hybrid coordinate
 128 C                            vertical levels.  (part of Section 4)
 129 C                            This element is actually a pointer to an array
 130 C                            that holds the data.
 131 C        gfld%num_coord = number of values in array gfld%coord_list().
 132 C        gfld%ndpts = Number of data points unpacked and returned.
 133 C        gfld%idrtnum = Data Representation Template Number
 134 C                       ( see Code Table 5.0)
 135 C        gfld%idrtmpl() = Contains the data values for the specified Data
 136 C                         Representation Template ( N=gfld%idrtnum ).  Each
 137 C                         element of this integer array contains an entry
 138 C                         (in the order specified) of Product Defintion
 139 C                         Template 5.N.
 140 C                         This element is actually a pointer to an array
 141 C                         that holds the data.
 142 C        gfld%idrtlen = Number of elements in gfld%idrtmpl().  i.e. number
 143 C                       of entries in Data Representation Template 5.N
 144 C                       ( N=gfld%idrtnum ).
 145 C        gfld%unpacked = logical value indicating whether the bitmap and
 146 C                        data values were unpacked.  If false,
 147 C                        gfld%bmap and gfld%fld pointers are nullified.
 148 C        gfld%expanded = Logical value indicating whether the data field
 149 C                         was expanded to the grid in the case where a
 150 C                         bit-map is present.  If true, the data points in
 151 C                         gfld%fld match the grid points and zeros were
 152 C                         inserted at grid points where data was bit-mapped
 153 C                         out.  If false, the data values in gfld%fld were
 154 C                         not expanded to the grid and are just a consecutive
 155 C                         array of data points corresponding to each value of
 156 C                         "1" in gfld%bmap.
 157 C        gfld%ibmap = Bitmap indicator ( see Code Table 6.0 )
 158 C                     0 = bitmap applies and is included in Section 6.
 159 C                     1-253 = Predefined bitmap applies
 160 C                     254 = Previously defined bitmap applies to this field
 161 C                     255 = Bit map does not apply to this product.
 162 C        gfld%bmap() = Logical*1 array containing decoded bitmap,
 163 C                      if ibmap=0 or ibap=254.  Otherwise nullified.
 164 C                      This element is actually a pointer to an array
 165 C                      that holds the data.
 166 C        gfld%fld() = Array of gfld%ndpts unpacked data points.
 167 C                     This element is actually a pointer to an array
 168 C                     that holds the data.
 169 C     IRET         INTEGER RETURN CODE
 170 C                    0      ALL OK
 171 C                    97     ERROR READING GRIB FILE
 172 C                    OTHER  GF_GETFLD GRIB UNPACKER RETURN CODE
 173 C
 174 C SUBPROGRAMS CALLED:
 175 C   BAREAD         BYTE-ADDRESSABLE READ
 176 C   GF_UNPACK6     UNAPCKS BIT_MAP SECTION
 177 C   GF_UNPACK7     UNAPCKS DATA SECTION
 178 C
 179 C REMARKS:
 180 C   DO NOT ENGAGE THE SAME LOGICAL UNIT FROM MORE THAN ONE PROCESSOR.
 181 C   THIS SUBPROGRAM IS INTENDED FOR PRIVATE USE BY GETGB2 ROUTINES ONLY.
 182 C
 183 C   Note that derived type gribfield contains pointers to many
 184 C   arrays of data.  The memory for these arrays is allocated
 185 C   when the values in the arrays are set, to help minimize
 186 C   problems with array overloading.  Because of this, users
 187 C   are encouraged to free up this memory, when it is no longer
 188 C   needed, by an explicit call to subroutine gf_free.
 189 C   ( i.e.   CALL GF_FREE(GFLD) )
 190 C
 191 C ATTRIBUTES:
 192 C   LANGUAGE: FORTRAN 90
 193 C
 194 C$$$
 195       USE GRIB_MOD
 196
 197       INTEGER,INTENT(IN) :: LUGB
 198       CHARACTER(LEN=1),INTENT(IN) :: CINDEX(*)
 199       INTEGER,INTENT(OUT) :: IRET
 200       TYPE(GRIBFIELD) :: GFLD
 201
 202       INTEGER :: LSKIP,SKIP6,SKIP7
 203       CHARACTER(LEN=1):: CSIZE(4)
 204       CHARACTER(LEN=1),ALLOCATABLE :: CTEMP(:)
 205       real,pointer,dimension(:) :: newfld
 206
 207       interface
 208          subroutine gf_unpack6(cgrib,lcgrib,iofst,ngpts,ibmap,
 209      &                         bmap,ierr)
 210            character(len=1),intent(in) :: cgrib(lcgrib)
 211            integer,intent(in) :: lcgrib,ngpts
 212            integer,intent(inout) :: iofst
 213            integer,intent(out) :: ibmap
 214            integer,intent(out) :: ierr
 215            logical*1,pointer,dimension(:) :: bmap
 216          end subroutine gf_unpack6
 217          subroutine gf_unpack7(cgrib,lcgrib,iofst,igdsnum,igdstmpl,
 218      &                         idrsnum,idrstmpl,ndpts,fld,ierr)
 219            character(len=1),intent(in) :: cgrib(lcgrib)
 220            integer,intent(in) :: lcgrib,ndpts,idrsnum,igdsnum
 221            integer,intent(inout) :: iofst
 222            integer,pointer,dimension(:) :: idrstmpl,igdstmpl
 223            integer,intent(out) :: ierr
 224            real,pointer,dimension(:) :: fld
 225          end subroutine gf_unpack7
 226       end interface
 227 C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 228 C  GET INFO
 229       NULLIFY(gfld%bmap,gfld%fld)
 230       IRET=0
 231       CALL GBYTE(CINDEX,LSKIP,4*8,4*8)
 232       CALL GBYTE(CINDEX,SKIP6,24*8,4*8)
 233       CALL GBYTE(CINDEX,SKIP7,28*8,4*8)
 234
 235 C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 236 C  READ AND UNPACK BIT_MAP, IF PRESENT
 237       IF ( gfld%ibmap.eq.0.OR.gfld%ibmap.eq.254 ) THEN
 238          ISKIP=LSKIP+SKIP6
 239          CALL BAREAD(LUGB,ISKIP,4,LREAD,CSIZE)    ! GET LENGTH OF SECTION
 240          CALL GBYTE(CSIZE,ILEN,0,32)
 241          ALLOCATE(CTEMP(ILEN))
 242          CALL BAREAD(LUGB,ISKIP,ILEN,LREAD,CTEMP)  ! READ IN SECTION
 243          IF (ILEN.NE.LREAD) THEN
 244             IRET=97
 245             DEALLOCATE(CTEMP)
 246             RETURN
 247          ENDIF
 248          IOFST=0
 249          CALL GF_UNPACK6(CTEMP,ILEN,IOFST,gfld%ngrdpts,idum,
 250      &                   gfld%bmap,ierr)
 251          IF (IERR.NE.0) THEN
 252             IRET=98
 253             DEALLOCATE(CTEMP)
 254             RETURN
 255          ENDIF
 256          DEALLOCATE(CTEMP)
 257       ENDIF
 258 C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 259 C  READ AND UNPACK DATA FIELD
 260       ISKIP=LSKIP+SKIP7
 261       CALL BAREAD(LUGB,ISKIP,4,LREAD,CSIZE)    ! GET LENGTH OF SECTION
 262       CALL GBYTE(CSIZE,ILEN,0,32)
 263       ALLOCATE(CTEMP(ILEN))
 264       CALL BAREAD(LUGB,ISKIP,ILEN,LREAD,CTEMP)  ! READ IN SECTION
 265       IF (ILEN.NE.LREAD) THEN
 266          IRET=97
 267          DEALLOCATE(CTEMP)
 268          RETURN
 269       ENDIF
 270       IOFST=0
 271       CALL GF_UNPACK7(CTEMP,ILEN,IOFST,gfld%igdtnum,gfld%igdtmpl,
 272      &                   gfld%idrtnum,gfld%idrtmpl,gfld%ndpts,
 273      &                   gfld%fld,ierr)
 274       IF (IERR.NE.0) THEN
 275          IRET=98
 276          DEALLOCATE(CTEMP)
 277          RETURN
 278       ENDIF
 279       DEALLOCATE(CTEMP)
 280 C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 281       !  If bitmap is used with this field, expand data field
 282       !  to grid, if possible.
 283       if ( gfld%ibmap .ne. 255 .AND. associated(gfld%bmap) ) then
 284             allocate(newfld(gfld%ngrdpts))
 285             !newfld=0.0
 286             !newfld=unpack(lgfld%fld,lgfld%bmap,newfld)
 287             n=1
 288             do j=1,gfld%ngrdpts
 289                 if ( gfld%bmap(j) ) then
 290                   newfld(j)=gfld%fld(n)
 291                   n=n+1
 292                 else
 293                   newfld(j)=0.0
 294                 endif
 295             enddo
 296             deallocate(gfld%fld);
 297             gfld%fld=>newfld;
 298             gfld%expanded=.true.
 299       else
 300          gfld%expanded=.true.
 301       endif
 302 C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 303       RETURN
 304       END