GRIB
to TXT converter
This page descibes the logic of the grib2txt routine, written by Hans
Gleisner (DMI), in
/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/hw_grib2ropp/hgl_script/grib2txt.
Interface:
!==============================================================================!
!
!
!
NAME
!
!
get_ecmwf_profs
!
!
!
!
DESCRIPTION
!
! Retrieve vertical profiles from ECMWF GRIB
files.
!
!
!
!
USAGE
!
! get_ecmwf_profs -gribroot
dir
&
!
!
-anfc
afstrng
&
!
!
-date dstrng -time
tstrng
&
!
!
-lat lat -lon
lon
&
!
!
-outfile
outfile
&
!
!
-datdir
datadir
&
!
!
[-ascii |
-ropp]
&
!
!
[-sil]
&
!
!
[-h]
!
!
!
! get_ecmwf_profs -gribroot
dir
&
!
!
-anfc
afstrng
&
!
!
-infile
infile
&
!
!
-datdir
datadir
&
!
!
[-ascii |
-ropp]
&
!
!
[-sil]
&
!
!
[-h]
!
!
!
!
ARGUMENTS
!
! -gribroot
<dir> Root
directory for the GRIB file
tree.
!
!
-anfc
<anfc> Analysis or
forecast: 'AN' or
'FC'
!
!
-date
<dstrng>
Date string, e.g.
'20080915'.
!
!
-time
<tstrng>
Time string, e.g. '0937' or
'122534'.
!
!
-lat
<lat> Latitude:
-90.0 to
+90.0.
!
!
-lon
<lon> Longitude:
-180.0 to
+180.0.
!
! -outfile
<outfile> Name of output
file
!
! -datdir
<datadir> Root directory
to the GSM-96 data
directory.
!
! -infile
<infile> Name of
input file (each line contains date,
time, !
!
lat, lon,
outfile).
!
!
-ascii
Output stored in a ROPP
file.
!
!
-ropp
Output stored in an ASCII
file.
!
!
-sil
Execute
silently.
!
!
-h
Causes a summary help to be
output.
!
!
!
!
DEPENDENCIES
!
! module
rom_const
!
! module
rom_utils
!
! module
rom_num
!
! module
rom_io
!
! module
grib_api
!
! module
ncdf
!
!
!
!
NOTES
!
! Already existing output files will be
overwritten.
!
!
!
!
AUTHOR
!
! Hans Gleisner
(hgl@dmi.dk)
!
!
!
!
VERSION
!
!
10-Jul-2012
!
!
!
! CODE
DESCRIPTION:
!
! Language: Fortran
90
!
!
!
!------------------------------------------------------------------------------!
Steps:
1) Read required occ dates and times.
2) Choose corresponding GRIB file.
If anfc = 'AN', this is the file with validity time within
±3h of the
occ time, in the subdirectory of gribroot called
yyyymmddhh/000.
If anfc = 'FC', this is the file with validity time within
±3h of the
occ time, in the subdirectory of gribroot called yyyymmddhh/HHH, which
contains the short range forecast from the most
recent analysis, which is valid at either 00Z or 12Z.
Thus, for example, specifying an occultation time of 0212Z on 20110508
would result in the extraction of a background from the 00Z analysis
(ie 2011050800/000) if anfc='AN', but a background from the T+12
forecast starting from the previous day's 12Z analysis (ie
2011050712/012) if
anfc='FC'.
[Explanation from
Hallgeir of DMI's ECMWF gribfile archive structure:
You may know this
already, as it is probably standard, but here comes a quick explanation
of the files in:
https://svn.grassaf.org/ropp/ropp_src/branches/dev/Share/hw_grib2ropp/grib_files
The 000 files are
analysis
files. The other files (006, 012, 018) are forecasts based on the 000
analysis in that directory. The numbers mean 006 hour-, 012 hour- and
018 hour forecasts. :)
We get forecasts at
00 and 12
hours every day, and analysis at 00, 06, 12, 18 hours every day. Our
file structure for a certain day, here January 2nd 2012, therefore
looks like this:
00 hour:
2012010200/000
2012010200/006
2012010200/012
2012010200/018
06 hour:
2012010206/000
12 hour:
2012010212/000
2012010212/006
2012010212/012
2012010212/018
18 hour:
2012010218/000
In other words, the
files you have gotten in svn is only based on one single analysis...]
All we have are:
idculv@eld037:> ls -ltR
/net/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/grib_files
/net/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/grib_files:
total 6
drwxr-xr-x 4 idculv satsense 50 Sep 17 12:15 grib2
drwxr-xr-x 4 idculv satsense 50 Sep 17 12:09 grib1
drwxr-xr-x 4 idculv satsense 50 Sep 17 12:02 grib2_double
/net/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/grib_files/grib2:
total 2
drwxr-xr-x 3 idculv satsense 106 Sep 17 12:22 2011050800
/net/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/grib_files/grib2/2011050800:
total 89390
-rw-r--r-- 1 idculv satsense 18056976 Sep 17 12:22 018
-rw-r--r-- 1 idculv satsense 18056976 Sep 17 12:22 006
-rw-r--r-- 1 idculv satsense 18056976 Sep 17 12:22 012
-rw-r--r-- 1 idculv satsense 18188228 Sep 17 12:22 000
/net/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/grib_files/grib1:
total 2
drwxr-xr-x 3 idculv satsense 106 Sep 17 12:15 2011050800
/net/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/grib_files/grib1/2011050800:
total 89294
-rw-r--r-- 1 idculv satsense 18036480 Sep 17 12:15 018
-rw-r--r-- 1 idculv satsense 18036480 Sep 17 12:15 006
-rw-r--r-- 1 idculv satsense 18036480 Sep 17 12:15 012
-rw-r--r-- 1 idculv satsense 18167620 Sep 17 12:15 000
/net/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/grib_files/grib2_double:
total 2
drwxr-xr-x 3 idculv satsense 106 Sep 17 12:09 2011070200
/net/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/grib_files/grib2_double/2011070200:
total 102086
-rw-r--r-- 1 idculv satsense 18056976 Sep 17 12:09 018
-rw-r--r-- 1 idculv satsense 28417536 Sep 17 12:09 006
-rw-r--r-- 1 idculv satsense 18056976 Sep 17 12:09 012
-rw-r--r-- 1 idculv satsense 18188228 Sep 17 12:09 000
3) Call read_ECMWF_profs. This redoes (?) the calculation of
the gribfilename, and checks for its existence. It then opens
gribfile (if available), which contains {lnsp,Φ0, T(i), q(i)},
ie 2N_lev+2 fields. Currently there is a check that it only
contains this many fields, so we'd need to switch that off for ROPP,
which should presumably be able to handle data from a standard ECMWF
GRIB file.
4) This calls read_ECMWF_3Dgrid, which redoes (?) the check that the
variable 'pv' is of size 2N_lev+2, and if so, it extracts the level
coeffs Ak and Bk from the 1st and 2nd halves of this array (in reverse
order). It then reads {lnsp, Φ0, T(i),
q(i)}, converts to {P0=exp(lnsp), Z0=Φ0/g_wmo, T, Q(g/kg)=q(kg/kg)*1000}, and currently fails if it comes across any
other fields - again, we may want to switch this off in the ROPP
implementation. It also seems to extract the lat and the lon every
time? Could be taken out of loop, since they're the same
here.
5) 3D pressure field on full levels, P, calculated according to
P(ilat,ilon,ilev) = 0.5d0*(ak(ilev) +
bk(ilev)*P0(ilat,ilon)) + 0.5d0*(ak(ilev+1) + bk(ilev+1)*P0(ilat,ilon))
6) 3D GPH field on full levels, Z, calculated by integrating R_dry
* Tvirt / g_wmo from Z0. If it didn't find Z0, then it
returns Z0=Z=mdi. There's some complicated business involving
call read_ECMWF_orography, which as far as I can tell is just looking
for Z0 again.
7) It also recalcalutes Z_grb if anfc='FC', upon return to subroutine read_ECMWF_profs.
8) It locates the target point's 4 nearest neighbours (in lat & lon) and bilinearly interpolates {Z0, P0, Z(i), P(i), T(i),
Q(i)} to it.
9) It computes orthometric (mean-sea level) heights, H, for the profiles.
a) u(iocc) = GeoidHeight(datadir, lat(iocc), lon(iocc)) =
the height of the EGM96 geoid above the WGS84 ellipsoid (ie, the
undulation). The location of the dir containing the data is under user
control. We could put this somewhere central in ROPP.
Copyright?
b) ! call curvature(lat(iocc), lon(iocc), az, Rcoc, Rc) is the
ROPP routine, with az=0. Why? - it can make a huge difference.
Presumably that's why it's not currently called. If it were, it would return Rcoc and Rc.
c) call geopot2msl_ropp(lat(iocc), u(iocc), Z(iocc,:),
H(iocc,:)). Currently this is basically the ROPP routine H =
geopotential2geometric(lat, Z), ie the geoid height isn't used. So we
convert Z wrt geoid to H wrt geoid (the MSL one), using Somigliana's
R_eff and g.
10) It then writes a profile in ascii format:
write(9,'(A)') '# Profile pressure (hPa), temperature (K), and specific humidity (g/kg).'
do ilev=1,Nlev
write(9,'(F8.3,1X,F7.2,1X,F8.5)') Pecm(i,ilev)/100.0d0, Tecm(i,ilev), Qecm(i,ilev)
enddo
(So H isn't used at present).
Example:
cd /data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/hgl_script/grib2txt
./compile_getprofs.bash
./get_ecmwf_profs.x \
-gribroot '/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/grib_files/grib2' \
-datdir '/net/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/hgl_script/ancillary' \
-anfc 'AN' \
-date '20110508' \
-time '0212' \
-lat 50.0 \
-lon -3.0 \
-ascii \
-outfile 'example1.out'
produces:
#
/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/grib_files/grib2/2011050800/000
at LAT: 50.000 LON: -3.000
# Surface geopotential (m), pressure (hPa), LS mask, and Nlev.
-1.43 1007.129 -1 91
# Profile pressure (hPa), temperature (K), and specific humidity (g/kg).
1005.936 286.18 8.83387
1003.050 286.13 8.92123
999.052 286.00 8.99149
993.804 285.83 9.07474
987.317 285.79 9.11437
979.704 285.56 9.12646
970.818 285.31 9.04204
960.400 285.02 8.94932
948.421 284.64 8.79558
934.885 284.54 8.47089
919.795 283.73 8.19792
903.156 282.82 7.83847
885.007 282.13 6.77560
865.423 280.86 6.09711
844.477 279.39 5.66150
822.233 277.96 6.10021
798.780 275.75 5.28384
774.238 274.14 5.06247
748.725 273.02 4.77243
722.344 271.41 4.45200
695.212 269.52 4.02497
667.472 267.62 3.54029
639.265 265.44 3.09754
610.708 263.43 2.81931
581.927 261.62 2.55990
553.069 259.44 2.14411
524.355 256.91 1.67860
496.066 254.23 1.12169
468.478 250.35 0.85643
441.848 247.22 0.45206
416.415 244.90 0.18884
392.282 242.20 0.07058
369.411 238.82 0.04959
347.751 235.26 0.05150
327.245 231.91 0.04244
307.838 228.77 0.03338
289.477 225.85 0.02766
272.111 223.53 0.02314
255.693 221.78 0.02003
240.176 220.33 0.01669
225.517 218.46 0.01454
211.671 216.09 0.01329
198.600 214.65 0.01097
186.264 216.03 0.00703
174.627 217.57 0.00450
163.652 220.72 0.00328
153.307 222.11 0.00311
143.559 220.27 0.00294
134.377 218.90 0.00288
125.732 218.71 0.00289
117.586 218.55 0.00290
109.887 217.88 0.00291
102.579 216.87 0.00288
95.616 216.75 0.00282
88.958 217.58 0.00281
82.572 218.56 0.00282
76.429 218.52 0.00284
70.506 217.73 0.00285
64.798 217.80 0.00286
59.315 217.86 0.00288
54.070 218.37 0.00291
49.073 219.28 0.00293
44.331 220.09 0.00294
39.848 219.71 0.00294
35.628 218.64 0.00295
31.674 218.01 0.00296
27.986 218.12 0.00299
24.565 218.94 0.00302
21.410 220.96 0.00308
18.518 223.01 0.00319
15.884 223.78 0.00326
13.504 224.73 0.00326
11.369 226.88 0.00326
9.471 231.47 0.00326
7.799 236.67 0.00330
6.342 238.83 0.00337
5.086 244.75 0.00346
4.018 249.20 0.00354
3.121 256.42 0.00362
2.380 263.96 0.00373
1.778 268.75 0.00383
1.299 271.71 0.00391
0.924 272.85 0.00396
0.640 270.36 0.00398
0.428 265.95 0.00400
0.277 260.09 0.00402
0.172 251.76 0.00401
0.101 237.05 0.00406
0.057 217.36 0.00339
0.030 201.11 0.00241
0.010 187.84 0.00156
This agrees pretty much with:
/data/nwp1/idculv/GRIB/grib_api-1.9.9_dir/bin
idculv@eld037:> grib_ls -l50,-3,1 $file |sed
-es/"2
ecmf
20110508
an
regular_ll
0
hybrid"/""/ |sort -nr
184 of 184 total grib messages in 1 files
184 of 184 grib messages in /data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/grib_files/grib2/2011050800/000
91
t
grid_simple 286.182
91
q
grid_simple 0.00883387
90
t
grid_simple 286.133
90
q
grid_simple 0.00892123
89
t
grid_simple 286.004
89
q
grid_simple 0.00899149
88
t
grid_simple 285.828
88
q
grid_simple 0.00907474
87
t
grid_simple 285.791
87
q
grid_simple 0.00911437
86
t
grid_simple 285.562
86
q
grid_simple 0.00912646
85
t
grid_simple 285.31
85
q
grid_simple 0.00904204
84
t
grid_simple 285.02
84
q
grid_simple 0.00894932
83
t
grid_simple 284.644
83
q
grid_simple 0.00879558
82
t
grid_simple 284.54
82
q
grid_simple 0.00847089
81
t
grid_simple 283.734
81
q
grid_simple 0.00819792
80
t
grid_simple 282.816
80
q
grid_simple 0.00783847
79
t
grid_simple 282.131
79
q
grid_simple 0.0067756
78
t
grid_simple 280.865
78
q
grid_simple 0.00609711
77
t
grid_simple 279.389
77
q
grid_simple 0.0056615
76
t
grid_simple 277.959
76
q
grid_simple 0.00610021
75
t
grid_simple 275.75
75
q
grid_simple 0.00528384
74
t
grid_simple 274.142
74
q
grid_simple 0.00506247
73
t
grid_simple 273.019
73
q
grid_simple 0.00477243
72
t
grid_simple 271.411
72
q
grid_simple 0.004452
71
t
grid_simple 269.523
71
q
grid_simple 0.00402497
70
t
grid_simple 267.619
70
q
grid_simple 0.00354029
69
t
grid_simple 265.44
69
q
grid_simple 0.00309754
68
t
grid_simple 263.427
68
q
grid_simple 0.00281931
67
t
grid_simple 261.623
67
q
grid_simple 0.0025599
66
t
grid_simple 259.445
66
q
grid_simple 0.00214411
65
t
grid_simple 256.911
65
q
grid_simple 0.0016786
64
t
grid_simple 254.227
64
q
grid_simple 0.00112169
63
t
grid_simple 250.347
63
q
grid_simple 0.00085643
62
t
grid_simple 247.216
62
q
grid_simple 0.000452057
61
t
grid_simple 244.898
61
q
grid_simple 0.000188837
60
t
grid_simple 242.201
60
q
grid_simple 7.05798e-05
59
t
grid_simple 238.817
59
q
grid_simple 4.95934e-05
58
t
grid_simple 235.257
58
q
grid_simple 5.14984e-05
57
t
grid_simple 231.915
57
q
grid_simple 4.24385e-05
56
t
grid_simple 228.768
56
q
grid_simple 3.33786e-05
55
t
grid_simple 225.854
55
q
grid_simple 2.766e-05
54
t
grid_simple 223.533
54
q
grid_simple 2.31432e-05
53
t
grid_simple 221.777
53
q
grid_simple 2.00323e-05
52
t
grid_simple 220.329
52
q
grid_simple 1.66893e-05
51
t
grid_simple 218.458
51
q
grid_simple 1.45435e-05
50
t
grid_simple 216.088
50
q
grid_simple 1.32918e-05
49
t
grid_simple 214.647
49
q
grid_simple 1.09673e-05
48
t
grid_simple 216.025
48
q
grid_simple 7.03335e-06
47
t
grid_simple 217.567
47
q
grid_simple 4.50015e-06
46
t
grid_simple 220.724
46
q
grid_simple 3.27826e-06
45
t
grid_simple 222.106
45
q
grid_simple 3.10689e-06
44
t
grid_simple 220.271
44
q
grid_simple 2.93553e-06
43
t
grid_simple 218.904
43
q
grid_simple 2.87592e-06
42
t
grid_simple 218.715
42
q
grid_simple 2.8871e-06
41
t
grid_simple 218.554
41
q
grid_simple 2.90386e-06
40
t
grid_simple 217.875
40
q
grid_simple 2.90945e-06
39
t
grid_simple 216.871
39
q
grid_simple 2.87965e-06
38
t
grid_simple 216.75
38
q
grid_simple 2.81602e-06
37
t
grid_simple 217.58
37
q
grid_simple 2.80837e-06
36
t
grid_simple 218.555
36
q
grid_simple 2.82151e-06
35
t
grid_simple 218.522
35
q
grid_simple 2.83583e-06
34
t
grid_simple 217.725
34
q
grid_simple 2.84932e-06
33
t
grid_simple 217.795
33
q
grid_simple 2.86236e-06
32
t
grid_simple 217.858
32
q
grid_simple 2.88346e-06
31
t
grid_simple 218.372
31
q
grid_simple 2.90988e-06
30
t
grid_simple 219.283
30
q
grid_simple 2.92827e-06
29
t
grid_simple 220.089
29
q
grid_simple 2.93849e-06
28
t
grid_simple 219.708
28
q
grid_simple 2.94274e-06
27
t
grid_simple 218.645
27
q
grid_simple 2.94806e-06
26
t
grid_simple 218.014
26
q
grid_simple 2.96126e-06
25
t
grid_simple 218.125
25
q
grid_simple 2.98719e-06
24
t
grid_simple 218.94
24
q
grid_simple 3.02117e-06
23
t
grid_simple 220.96
23
q
grid_simple 3.08349e-06
22
t
grid_simple 223.01
22
q
grid_simple 3.18947e-06
21
t
grid_simple 223.778
21
q
grid_simple 3.26119e-06
20
t
grid_simple 224.728
20
q
grid_simple 3.25947e-06
19
t
grid_simple 226.885
19
q
grid_simple 3.26296e-06
18
t
grid_simple 231.473
18
q
grid_simple 3.25986e-06
17
t
grid_simple 236.668
17
q
grid_simple 3.3009e-06
16
t
grid_simple 238.828
16
q
grid_simple 3.37098e-06
15
t
grid_simple 244.755
15
q
grid_simple 3.45707e-06
14
t
grid_simple 249.2
14
q
grid_simple 3.53864e-06
13
t
grid_simple 256.419
13
q
grid_simple 3.6194e-06
12
t
grid_simple 263.959
12
q
grid_simple 3.73099e-06
11
t
grid_simple 268.754
11
q
grid_simple 3.82798e-06
10
t
grid_simple 271.708
10
q
grid_simple 3.9072e-06
9
t
grid_simple 272.846
9
q
grid_simple 3.95614e-06
8
t
grid_simple 270.363
8
q
grid_simple 3.97863e-06
7
t
grid_simple 265.951
7
q
grid_simple 3.99851e-06
6
t
grid_simple 260.086
6
q
grid_simple 4.01725e-06
5
t
grid_simple 251.757
5
q
grid_simple 4.01412e-06
4
t
grid_simple 237.045
4
q
grid_simple 4.05898e-06
3
t
grid_simple 217.361
3
q
grid_simple 3.39441e-06
2
t
grid_simple 201.113
2
q
grid_simple 2.40852e-06
1
z
grid_simple -14.0557 = 9.80665 * -1.4332825
1
t
grid_simple 187.843
1
q
grid_simple 1.56344e-06
1
lnsp grid_simple
11.52 = ln(100710.)
edition
centre
date
dataType gridType
stepRange typeOfLevel
level
shortName packingType value
Other grid Points
Input Point: latitude=50.00 longitude=-3.00
Grid Point chosen #2 index=14577 latitude=50.00 longitude=-3.00 distance=0.00 (Km)
/data/nwp1/idculv/ROPP/ropp_src/branches/dev/Share/ic_grib2ropp/grib_files/grib2/2011050800/000
- 4 - index=14217 latitude=51.00 longitude=-3.00 distance=111.20 (Km)
- 3 - index=14218 latitude=51.00 longitude=-2.00 distance=131.79 (Km)
- 2 - index=14577 latitude=50.00 longitude=-3.00 distance=0.00 (Km)
- 1 - index=14578 latitude=50.00 longitude=-2.00 distance=71.48 (Km)
Aspects relevant to ROPP implementation:
1) Need to be able to read from general ECMWF grib2 (and grib1?) dump file.
2) Need to include GRIB API as an external library, like netCDF or
BUFR. Copyright? Currently I use grib_api_1.9.9, as 1.9.16 has
problems with DMI's example file. Under discussion with ECMWF.
3) Need to clean up the reading/checking (don't generate lat and lon
for every call: bracket this with an (IF Nmess ==1) control statement.
What's the need for call_read_ECMWF_orography? What's the
need for the Geoid Height and Curvature business? Forget about them
intially.
4) Try to write with the flexibility to read two grib files and
interpolate (linearly) in time between them (eventually). Ensure
we can handle GRIB1 and GRIB2 (and GRIB2_double?). Keep ability
to produce multiple bgr files per call. I suggest we don't bother
with Hallgeir's suggestion to allow forward modelling as part of the
script.
5) Natural place to put it: ropp_io/tools/grib2text.
Last updated 24th
September 2012
by ian.culverwell@metoffice.gov.uk