0

I want to compile a meteorological dataset out of DWD weather data (https://opendata.dwd.de/climate_environment/REA/COSMO_REA6/daily/2D/) in .grb format with one dimensional coordinates (x and y) into a spatial subset in the netCDF fileformat with twodimensional lat lon coordinates. And I need to do it in python and on windows. So tools like "ecCodes" or "CDO" are not viable. I have tried to do this with xarray since it is able to read .grb files, but the very different structure of netCDF and GRIB, seems to make it difficult for me.

The error message I get:

KeyError: "no index found for coordinate 'latitude'"

But when I look at the dataset there should be longitude and latitude:

Out[24]: 
<xarray.Dataset>
Dimensions:     (time: 31, y: 824, x: 848)
Coordinates:
  * time        (time) datetime64[ns] 2008-01-01T01:00:00 ... 2008-01-31T01:0...
    step        timedelta64[ns] 1 days
    surface     float64 0.0
    latitude    (y, x) float64 21.95 21.97 21.99 22.01 ... 66.76 66.73 66.7
    longitude   (y, x) float64 ...
    valid_time  (time) datetime64[ns] 2008-01-02T01:00:00 ... 2008-02-01T01:0...
Dimensions without coordinates: y, x
Data variables:
    ASWDIFD_S   (time, y, x) float32 ...
Attributes:
    GRIB_edition:            1
    GRIB_centre:             edzw
    GRIB_centreDescription:  Offenbach
    GRIB_subCentre:          255
    Conventions:             CF-1.7
    institution:             Offenbach
    history:                 2022-12-09T13:41 GRIB to CDM+CF via cfgrib-0.9.1...

My current code:

from tqdm import tqdm
import xarray as xr
from tqdm import tqdm
import os 

folders = [x for x in os.listdir() ]
for var in tqdm(folders,desc=" outer", position=0):
    variables = [x for x in os.listdir(var) if not ".idx" in x]
    for i, file in tqdm(enumerate(variables), desc="inner", position = 1):
        if i == 0:
            try:
                ds = xr.open_dataset(var+"/"+file)
                ds = ds.sel(latitude = slice(53.6287, 54.95), longitude = slice(9.35, 14.685))
            except:
                ds = xr.open_dataset(var+"/"+file,engine="cfgrib")
                ds = ds.sel(latitude = slice(53.6287, 54.95), longitude = slice(9.35, 14.685))
            i += 1
        else:
            try:
                ds   = xr.concat([ds, xr.open_dataset(var+"/"+file)], dim = "time")
            except:
                ds   = xr.concat([ds, xr.open_dataset(var+"/"+file,engine="cfgrib")], dim = "time")
    
    ds.to_netcdf("COSMO_REA6_"+var+"_2008-2019.nc")
b_robran
  • 1
  • 1
  • Welcome to SO. Can you possibly simplify your code? The for loops etc. don't help anyone solve the problem. Can you rewrite it for a single file? It looks to me like more than half of your lines are irrelevant to the problem. – Robert Wilson Dec 09 '22 at 15:13
  • Unfortunatly I did not find a solution. My data provider reprocessed everything for me in netcdf format with one dimensional lat and lon coordinates. – b_robran Jul 12 '23 at 08:19

0 Answers0