I can be able to add Wradlib environment, Data that I am using is from DWD which refers to Germany. data type is Radolan binary data. I can able to read the data and able to plot it using Matplotlib.
I am attaching the graph and also shapefile of my region just for an overview. How can I extract data for specific region and also can I be able to look into the values with in the grids?
the code that I used to plot this is,
import glob
import os
import wradlib as wrl
import warnings
warnings.filterwarnings(“ignore”)
import matplotlib.pyplot as pl
import numpy as np
import xarray as xr
from mpl_toolkits import mplot3d
import time
import os
from os import listdir
from os.path import isfile, join
mypath = os.getcwd()
mypath = mypath + ‘\combine_test’
print(mypath)
os.environ[“WRADLIB_DATA”] = mypath
onlyfiles = [f for f in listdir(mypath) if isfile(join(mypath, f))]
onlyfiles.pop(0)
onlyfiles.pop(0)
def plot_radolan(data, attrs, grid, clabel=None):
fig = pl.figure(figsize=(10,8))
ax = fig.add_subplot(111, aspect=‘equal’)
x = grid[:,:,0]
y = grid[:,:,1]
pm = ax.pcolormesh(x, y, data, cmap=‘viridis’)
cb = fig.colorbar(pm, shrink=0.75)
cb.set_label(clabel)
pl.xlabel(“x [km]”)
pl.ylabel(“y [km]”)
pl.title(‘{0} Product\n{1}’.format(attrs[‘producttype’],
attrs[‘datetime’].isoformat()))
pl.xlim((x[0,0],x[-1,-1]))
pl.ylim((y[0,0],y[-1,-1]))
pl.grid(color=‘r’)
radolan_grid_xy = wrl.georef.get_radolan_grid(900,900)
radolan_egrid_xy = wrl.georef.get_radolan_grid(1500,1400)
radolan_wgrid_xy = wrl.georef.get_radolan_grid(1100, 900)
x = radolan_grid_xy[:,:,0]
y = radolan_grid_xy[:,:,1]
xe = radolan_egrid_xy[:,:,0]
ye = radolan_egrid_xy[:,:,1]
xw = radolan_wgrid_xy[:,:,0]
yw = radolan_wgrid_xy[:,:,1]
fpath = “raa01-rw_10000-0701010050-dwd—bin.gz”
f = wrl.util.get_wradlib_data_file(fpath)
img,meta = wrl.io.read_radolan_composite(f)
#print(onlyfiles)
flist = [wrl.util.get_wradlib_data_file(f) for f in onlyfiles]
comp3 = wrl.io.open_radolan_mfdataset(flist)
comp3.RW.plot(col=‘time’)
latitude = 54.307857
Replace with the desired latitude value
longitude = 10.117256 # Replace with the desired longitude value
radolan_location = wrl.georef.get_radolan_location()
x, y = wrl.georef.reproject(latitude, longitude,
projection_source=wrl.georef.epsg_to_osr(4326), # EPSG 4326 is WGS84
projection_target=wrl.georef.create_osr(“dwd-radolan”))
Find the nearest grid point in the Radolan data
x, y = wrl.util.find_nearest(polar_coords[…, 0], polar_coords[…, 1], radolan_coords[…, 0], radolan_coords[…, 1])
x = int(np.round(x))
y = int(np.round(y))
x = max(0, min(x, img.shape[1] - 1))
y = max(0, min(y, img.shape[0] - 1))
data = np.ma.masked_equal(img, -9999)
Extract the Radolan data value at the specific location
radolan_value = img[y, x]
plot_radolan(data, meta, radolan_grid_xy, clabel=‘mm * h-1’)