ARJO SEGERS
KNMI
Date: 3 March 2009
Provides a module:
use GribFile, only : ...
The module provides a structure type containing all what should
be known to open, read, and close a grib file:
type(TGribFile) :: gribfile
The gribfile is opened either for reading or writing with
an 'Init' routine, which is in fact an interface to 'PbOpen'.
A type TgribFile contains buffers to store one record
of a grib file (a 'grib message').
1. To open a grib file for reading, use:
call Init( gribfile, filename, 'r', status )
type(TGribFile), intent(inout) :: gribfile
character(len=*), intent(in) :: filename
character(len=*), intent(in) :: mode
integer, intent(out) :: status
To read next record in buffer, use:
call ReadRecord( gribfile, status )
Contents of the message is extracted with the 'Get' routine,
for example:
call Get( gribfile, status, pid=thepid )
Complete syntaxis and description of arguments
(see also the 'GRIB SECTIONS' part later on):
call Get( gribfile, status, &
model_id, pid, &
levtype, level, hyb_a, hyb_b, &
reftime, timerange, &
gridtype, &
ll, &
lon_first, lon_last, lon_inc, lon_n, &
lat_first, lat_last, lat_inc, lat_n, &
T, sh, &
N, gg )
integer, intent(out), optional :: model_id
! Element 3 of section 1.
integer, intent(out), optional :: pid
! Parameter identification.
! See the integer parameters 'pid_' in the code.
! Element 6 of section 1.
integer, intent(out), optional :: levtype, level
! Elements 7-9 of section 1.
! --> Predefined values for level type:
! integer, parameter :: levtype_sfc = 1 ! surface
! integer, parameter :: levtype_hyb = 109 ! hybrid level
real, intent(out), optional :: hyb_a(:), hyb_b(:)
! Vertical parameters in real section 2
integer, intent(out), optional :: reftime(5)
! Reference time: a 5 element integer array with
! year_inc_century, month, day, hour, minutes
! Elements 10-14,21 of section 1.
integer, intent(out), optional :: timerange(4)
! Time range: a 4 element integer array with
! time range unit, value 1, value 2, indicator
! Elements 15-18 of section 1.
integer, intent(out), optional :: gridtype
! Element 1 of section 2.
! --> Predefined values:
! integer, parameter :: gridtype_ll = 0 ! lat/lon
! integer, parameter :: gridtype_gg = 4 ! gaussian grid
! integer, parameter :: gridtype_sh = 50 ! spectral
integer, intent(out), optional :: lon_first, lon_last, lon_inc, lon_n
integer, intent(out), optional :: lat_first, lon_last, lat_inc, lat_n
! Elements of section 2
! NOTE: values for lat/lon in mili degrees !
real, intent(out), optional :: ll(:,:)
! Value of lat/lon grid (stored from west->east, north->south!)
complex, intent(out), optional :: sh(:)
integer, intent(out), optional :: T
! spectral coefficients;
! size should match with spectral truncation 'T'
real, intent(out), optional :: gg(:)
integer, intent(out), optional :: N
! data values on Gaussian grid;
! number of values should match grid size 'N'
Contents could be checked with the 'Check' routine,
for example:
call Check( gribfile, status, pid=154 )
Complete syntaxis:
call Check( gribfile, status, debug=1, &
model_id, pid, &
levtype, level, &
reftime, timerange, &
gridtype, &
lon_first, lon_last, lon_inc, lon_n, &
lat_first, lat_last, lat_inc, lat_n, &
T )
type(TGribFile), intent(in) :: gribfile
integer, intent(out) :: status
integer, intent(in), optional :: debug
integer, intent(in), optional :: model_id
integer, intent(in), optional :: pid
integer, intent(in), optional :: levtype, level
integer, intent(in), optional :: reftime(5)
integer, intent(in), optional :: timerange(4)
integer, intent(in), optional :: gridtype
integer, intent(in), optional :: lon_first, lon_last, lon_inc, lon_n
integer, intent(in), optional :: lat_first, lat_last, lat_inc, lat_n
integer, intent(in), optional :: T
integer, intent(in), optional :: status
integer, intent(out), optional :: status
See the 'Get' command above for a description.
In case of one or more fields not matching the grib field:
o return status is equal to number of tests failed;
o info about the failed test is printed if debug is present and >0 .
Note: hybride coefficients can not checked due to a
lack of precission ... Similar for grid values (ll, gg, and sh).
The grib file is closed with a 'Done' routine
(interface to 'PbClose'):
call Done( gribfile, status )
2. To open a grib file for writing, use:
call Init( gribfile, 'test.gb', 'w', status )
For safety, set all sections in the buffer to zero
or an other safe value:
call Clear( gribfile, status )
Fill some or all of the appropriate fields with
the 'Set' command:
call Set( gribfile, status, &
model_id, pid, &
levtype, hyb_a, hyb_b, level, &
reftime, timerange, &
gridtype, &
ll, lon_n, lat_n, &
lon_first, lon_inc, lon_last, &
lat_first, lat_inc, lat_last, &
scanning_mode, nbits )
type(TGribFile), intent(inout) :: gribfile
integer, intent(out) :: status
integer, intent(in), optional :: debug
integer, intent(in), optional :: model_id
integer, intent(in), optional :: pid
integer, intent(in), optional :: levtype
integer, intent(in), optional :: level
real, intent(in), optional :: hyb_a(:), hyb_b(:)
integer, intent(in), optional :: reftime(5) ! (/yy,mm,dd,hh,min/)
integer, intent(in), optional :: timerange(4)
integer, intent(in), optional :: gridtype
real, intent(in), optional :: ll(:,:)
integer, intent(in), optional :: lon_n, lat_n
integer, intent(in), optional :: lon_first, lon_inc, lon_last ! mili degree
integer, intent(in), optional :: lat_first, lat_inc, lat_last ! mili degree
See the 'Get' command for a description of most of the fields.
Only used as arguments to 'Set':
integer, intent(in), optional :: scanning_mode
! Something with storage order; does not work proper yet.
integer, intent(in), optional :: nbits
! Number of bits used to store a real value; default 24 .
To write the record to the file, use:
call WriteRecord( gribfile, status )
3. To check wether a file is assigned to a grib structure,
use the logical function 'Opened' :
if ( Opened(gribfile) ) ...
logical :: Opened
type(TgribFile), intent(in) :: gribfile
4. To copy the buffer from one to another sturcture, use:
call CopySections( gribIn, gribOut, status )
type(TGribFile), intent(in) :: gribIn
type(TGribFile), intent(out) :: gribOut
to copy all sections except real section 4 (the actual data);
to copy the data too, use :
call CopyAllSections( gribIn, gribOut, status )
type(TGribFile), intent(in) :: gribIn
type(TGribFile), intent(out) :: gribOut
Compile together with the GribEx or Emos library from ECMWF:
f90 -o test.exe ... file_grib.o ... -l gribex
See also GRIBEX
at the ECMWF website.
Overview of some useful entries in grib files. For code tables, see:
Section 1 - Product Definition Section.
---------------------------------------
1. Code Table 2 Version Number. 128
2. Originating centre identifier. 98
3. Model identification. 199
4. Grid definition. 255
5. Flag (Code Table 1) 10000000
6. Parameter identifier (Code Table 2).
7. Type of level (Code Table 3). 109
8. Value 1 of level (Code Table 3).
9. Value 2 of level (Code Table 3).
10. Year of reference time of data. 100
11. Month of reference time of data. 8
12. Day of reference time of data. 10
13. Hour of reference time of data. 12
14. Minute of reference time of data. 0
15. Time unit (Code Table 4). 1
16. Time range one. 0
17. Time range two. 0
18. Time range indicator (Code Table 5) 0
21. Century of reference time of data.
Section 2 - Grid Description Section.
-------------------------------------
Southern latitudes and Western longitudes are negative.
Unit is mili degrees.
1. Data represent type = lat/long (Table 6) 0
2. Number of points along a parallel. 360
3. Number of points along a meridian. 181
4. Latitude of first grid point. 90000
5. Longitude of first grid point. -180000
6. Resolution and components flag. 10000000
7. Latitude of last grid point. -90000
8. Longitude of last grid point. 179000
9. i direction (East-West) increment. 1000
10. j direction (North-South) increment. 1000
11. Scanning mode flags (Code Table 8) 00000000
1. Data represent type = spectral (Table 6) 50
2. J - Pentagonal resolution parameter. 106 <--- T
3. K - Pentagonal resolution parameter. 106
4. M - Pentagonal resolution parameter. 106
5. Representation type (Table 9) 1
6. Representation mode (Table 10). 2
1. Data represent type = gaussian lat/lon (Table 6) 4
2. Number of points along a parallel <0=reduced, >0=regular
3. Number of points along a meridian. 160
4. Latitude of first grid point. 89141
5. Longitude of first grid point. 0
6. Resolution and components flag. 00000000
7. Latitude of last grid point. -89141
8. Longitude of last grid point. 358878
9. i direction (East-West) increment Not given
10. Number of parallels between pole and equator. 80 <--- N
11. Scanning mode flags (Code Table 8) 00000000
23-22*2N. Number of points along a parallel for reduced grid
12. Number of vertical coordinate parameters. 122
(parameters stored in real part of section 2, starting form index 11)
Section 4 - Binary Data Section.
-------------------------------------
1. Number of data values coded/decoded. 11556
2. Number of bits per data value. 16
3. Type of data (0=grid pt, 128=spectral). 128
4. Type of packing (0=simple, 64=complex). 64
5. Type of data (0=float, 32=integer). 0
6. Additional flags (0=none, 16=present). 0
7. Reserved. 0
8. Number of values (0=single, 64=matrix). 0
9. Secondary bit-maps (0=none, 32=present). 0
10. Values width (0=constant, 16=variable). 0
11. Byte offset of start of packed data (N). 2214
12. Power (P * 1000). 500
13. Pentagonal resolution parameter J for subset. 20
14. Pentagonal resolution parameter K for subset. 20
15. Pentagonal resolution parameter M for subset. 20
3. Type of data (0=grid pt, 128=spectral). 0
4. Type of packing (0=simple, 64=complex). 0
5. Type of data (0=float, 32=integer). 0
6. Additional flags (0=none, 16=present). 0
7. Reserved. 0
8. Number of values (0=single, 64=matrix). 0
9. Secondary bit-maps (0=none, 32=present). 0
10. Values width (0=constant, 16=variable). 0