Page History

EFAS data's  x and  y coordinates are not projected coordinates but matrix indexes (i, j), It is necessary to download the upstream area static file that contains the projected coordinates and replace it in EFAS.

retrieve,
class=ce,
date=2020-10-14,
expver=1,
hdate= 20070304/20070307,#20070304/20070307/20070311/20070314/20070318/20070321/20070325/20070328/20070401/20070404/20070408/20070411/20070415/20070418/20070422/20070425/20070429/20070502/20070506/20070509/20070513/20070516/20070520/20070523/20070527/20070530,
levtype=sfc,
model=lisflood,
number=0/1/2/3/4/5/6/7/8/9/10,
origin=ecmf,
param=240023,
step=0/6/12/18,
stream=efrf,
time=00:00:00,
type=pf,
target="mars_efas_reforecast.grib"

import cdsapi

c = cdsapi.Client()

c.retrieve(
    'efas-reforecast',
    {
        'format': 'grib',
        'product_type': 'ensemble_perturbed_reforecasts',
        'variable': 'river_discharge_in_the_last_6_hours',
        'model_levels': 'surface_level',
        'hyear': '2007',
        'hmonth': 'march',
        'hday': [
            '04', '07',
        ],
        'leadtime_hour': [
            '0', '12', '18',
            '6',
        ],
    },
    'efas_reforecast.grib')

import xarray as xr
import pandas as pd
from pyproj import Transformer,CRS

parameter = "dis06"

ds = xr.open_dataset("mars_efas_reforecast.grib", engine="cfgrib")
df = pd.read_csv("GRDC.csv")
uparea = xr.open_dataset("ec_uparea4.0.nc") # the upstream area

# replace x, y
ds["x"] = uparea["x"]
ds["y"] = uparea["y"]

# define reprojection parameters
laea_proj = CRS.from_proj4("+proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs")
transformer = Transformer.from_crs('epsg:4326', laea_proj, always_xy=True)

total = len(df)

rows = []
count = 0
for lon, lat, id in zip(df.long, df.lat, df.grdc_no):
    x1, y1 = transformer.transform(lon, lat)
    extracted = ds.sel(x=x1, y=y1, number = 1, method="nearest")[parameter]
    df_temp = extracted.drop_vars(["surface", "number"]).to_dataframe().reset_index()
    df_temp["grdc"] = str(id)
    df_temp = df_temp.drop("step", axis=1)
    df_temp = df_temp.set_index(["grdc","time"])
    rows.append(df_temp)
    count += 1
    print(f"progress: {count/total*100} %")

out = pd.concat(rows)
out.to_csv("extracted.csv", index="grdc")