Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

GFS: Missing levels for sigmaLayer and depthBelowLandLayer #195

Closed
maxhollmann opened this issue Jan 23, 2021 · 5 comments
Closed

GFS: Missing levels for sigmaLayer and depthBelowLandLayer #195

maxhollmann opened this issue Jan 23, 2021 · 5 comments
Labels
enhancement New feature or request

Comments

@maxhollmann
Copy link

maxhollmann commented Jan 23, 2021

cfgrib seems to miss some levels for the level types sigmaLayer and depthBelowLandLayer.

I'm working with this grib file: https://noaa-gfs-bdp-pds.s3.amazonaws.com/gfs.20210115/00/gfs.t00z.pgrb2.0p25.f003

wgrib2 -lev ~/Downloads/gfs.t00z.pgrb2.0p25.f003 shows four levels for each of these level types:

425:224793214:0-0.1 m below ground
427:225650288:0.1-0.4 m below ground
429:226508191:0.4-1 m below ground
431:227363072:1-2 m below ground

558:309354159:0.33-1 sigma layer
559:310143600:0.44-1 sigma layer
560:310922234:0.72-0.94 sigma layer
561:311736866:0.44-0.72 sigma layer

Reading the same file with cfgrib results in only two:

import cfgrib

ds = cfgrib.open_datasets("gfs.t00z.pgrb2.0p25.f003")

d = next(d for d in ds if 'depthBelowLandLayer' in d.variables)
d.depthBelowLandLayer
# <xarray.IndexVariable 'depthBelowLandLayer' (depthBelowLandLayer: 2)>
# array([0, 1])
# Attributes:
#     long_name:      soil depth
#     units:          m
#     positive:       down
#     standard_name:  depth

d = next(d for d in ds if 'sigmaLayer' in d.variables)
d.sigmaLayer
# <xarray.IndexVariable 'sigmaLayer' (sigmaLayer: 2)>
# array([0, 1])
# Attributes:
#     long_name:  original GRIB coordinate for key: level(sigmaLayer)
#     units:      1

d.r.shape
# (2, 721, 1440)
@shahramn
Copy link
Collaborator

shahramn commented Jan 23, 2021

Several of the messages in this file use non-standard i.e. local values for "typeOfFirstFixedSurface"
and therefore the ecCodes key "typeOfLevel" returns "unknown".

% grib_ls -p count,typeOfFirstFixedSurface,typeOfSecondFixedSurface,typeOfLevel gfs.t00z.pgrb2.0p25.f003 |\
               grep unknown
8                         220                       255                       unknown
9                         220                       255                       unknown
10                        220                       255                       unknown
476                       200                       255                       unknown
477                       200                       255                       unknown
478                       200                       255                       unknown
479                       200                       255                       unknown
480                       214                       255                       unknown
548                       204                       255                       unknown

The key typeOfFirstFixedSurface and typeOfSecondFixedSurface use Code Table 4.5 in which we have:

...
184 184 Grid tile glacier ice and inland ice fraction as a model surface
# 185-191 Reserved 
# 192-254 Reserved for local use 
255 255 Missing 

So the values used in this file are in the range 192-254 i.e. they are local

@maxhollmann
Copy link
Author

Right, but that's not what I mean here. The level types sigmaLayer and depthBelowLandLayer both get recognized. The values for these two levels are different when parsed by eccodes than by wgrib2:

grib_ls -p count,typeOfLevel,level gfs.t00z.pgrb2.0p25.f003 | grep -E "depthBelowLandLayer|sigmaLayer":

424          depthBelowLandLayer  0           
425          depthBelowLandLayer  0           
426          depthBelowLandLayer  0           
427          depthBelowLandLayer  0           
428          depthBelowLandLayer  0           
429          depthBelowLandLayer  0           
430          depthBelowLandLayer  1           
431          depthBelowLandLayer  1           
558          sigmaLayer   0           
559          sigmaLayer   0           
560          sigmaLayer   1           
561          sigmaLayer   0           

wgrib2 -lev gfs.t00z.pgrb2.0p25.f003 | grep -E "sigma layer|below ground":

424:224283344:0-0.1 m below ground
425:224793214:0-0.1 m below ground
426:225146662:0.1-0.4 m below ground
427:225650288:0.1-0.4 m below ground
428:226004347:0.4-1 m below ground
429:226508191:0.4-1 m below ground
430:226864070:1-2 m below ground
431:227363072:1-2 m below ground
558:309354159:0.33-1 sigma layer
559:310143600:0.44-1 sigma layer
560:310922234:0.72-0.94 sigma layer
561:311736866:0.44-0.72 sigma layer

@shahramn
Copy link
Collaborator

OK. Thanks for the clarification.
You have to decode the "level" keys as DOUBLEs rather than INTEGERs
Try:

% grib_ls -p count,typeOfLevel,topLevel:d,bottomLevel:d gfs.t00z.pgrb2.0p25.f003 |\
    grep -E "depthBelowLandLayer|sigmaLayer"
424          depthBelowLandLayer  0            0.1
425          depthBelowLandLayer  0            0.1
426          depthBelowLandLayer  0.1          0.4
427          depthBelowLandLayer  0.1          0.4
428          depthBelowLandLayer  0.4          1
429          depthBelowLandLayer  0.4          1
430          depthBelowLandLayer  1            2
431          depthBelowLandLayer  1            2
558          sigmaLayer   0.33         1
559          sigmaLayer   0.44         1
560          sigmaLayer   0.72         0.94
561          sigmaLayer   0.44         0.72

Here I am using the keys "topLevel" and "bottomLevel" decoded as doubles (using the ":d" qualifier)

@alexamici
Copy link
Contributor

@shahramn then this looks like an actual bug (or missing feature) in cfgrib as it only requests the level key and it accepts the default data type, that is int.

Thanks for the clarification.

@alexamici alexamici added the enhancement New feature or request label Jan 28, 2021
point9repeating added a commit to point9repeating/cfgrib that referenced this issue Mar 4, 2021
…attributes

and if the level is actually an integer cast it back to integer

fixes ecmwf#195
@alexamici
Copy link
Contributor

Levels will be read as floats in the next release.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
enhancement New feature or request
Projects
None yet
Development

Successfully merging a pull request may close this issue.

3 participants