Importing Radolan data into arcGIS using wradlib

[lord rufusus]

Hey Guys,

I am still a student and new to the topic. I want to use RADOLAN Data (Index of /weather/radar/composit/rv/ (dwd.de)) in ArcGIS Pro. The DWD told me that I have to parse the radolan data format first using wradlib. Is that possible?

In the examples in the documentation the data has always bin.bz format but my data is in .tar.bz2 format. I am honestly not even sure if I have to unzip it first.  Never worked with a similar Datatype before.

Hope you can help me with some ideas or example code. Thanks in Advance!

Paul

[Kai Muehlbauer]

Hi Paul,

yes, parsing the file using wradlib is possible. Unfortunately wradlib
only supports gzip compressed files out of the box.

Also your archive files (tar.bz2) seems to contain multiple data files.
My suggestion would be to decompress your archive and then use one of
the available wradlib radolan readers (as in the examples).

HTH,
Kai

Am 13.09.22 um 11:49 schrieb lord rufusus:

> Hey Guys,
>
> I am still a student and new to the topic. I want to use RADOLAN Data (Index
> of /weather/radar/composit/rv/ (dwd.de)

> <https://opendata.dwd.de/weather/radar/composit/rv/>) in ArcGIS Pro. The

[lord rufusus]

Thx Kai,

I can parse the radolan data into GeoTiff Format that I can visualize in ArcGis Pro, using the following code:

filename = wrl.util.get_wradlib_data_file(file_name)
data_raw, meta = wrl.io.read_radolan_composite(file_name)

# This is the RADOLAN projection
proj_osr = wrl.georef.create_osr("dwd-radolan")

# Get projected RADOLAN coordinates for corner definition
xy_raw = wrl.georef.get_radolan_grid(900, 900)

data, xy = wrl.georef.set_raster_origin(data_raw, xy_raw, "upper")

# create 3 bands
data = np.stack((data))
#data = np.stack((data, data + 100, data + 1000))
ds = wrl.georef.create_raster_dataset(data, xy, projection=proj_osr)
wrl.io.write_raster_dataset( file_name  , ds, "GTiff")

Unfortunately something goes wrong with the projection. Do you have any idea, how to fix that?

[Kai Muehlbauer]

Hi,

you have the wrong shape applied and thus get the wrong gridsize. That
leads to the shifted location.

Here is a minimal example to produce georeferenced geotiffs from your data:

filename = 'DE1200_RV2209200825.tar.bz2'
with tarfile.open(filename, "r:bz2") as tar:
for fname in tar.getnames():
f = tar.extractfile(fname)
data_raw, meta = wrl.io.read_radolan_composite(f)
xy_raw = wrl.georef.get_radolan_grid(*data_raw.shape)

data, xy = wrl.georef.set_raster_origin(data_raw, xy_raw, "upper")

proj_osr = wrl.georef.create_osr("dwd-radolan")

ds = wrl.georef.create_raster_dataset(data, xy,
projection=proj_osr)

wrl.io.write_raster_dataset(f"{fname}.tif", ds, "GTiff")

HTH,
Kai

[lord rufusus]

Hey Kai,

vielen Dank schonmal für die schnellen Rückmeldungen!

Aber keine Ahnung was ich falsch mache. Wenn ich deinen Code so 1:1 kopiere, passen zwar die Proportionen besser, aber es ist noch weiter nördlich als vorher. Habe auch schon probiert, die Daten schon mit wradlib in 4326 zu transformieren, aber da liegt das auch falsch.

Sowohl Map als auch das Raster sind in der Radolan Projektion.

Du hast nicht zufällig noch eine Idee woran das liegen kann? Tut mir auch leid dass ich da so nerve, aber ich brauch das eigentlich für meine Bachelorarbeit ;) 

LG Paul

[Kai Muehlbauer]

Hi Paul,

I'll stay in English, maybe others have the same problem, too.

I've used QGIS to inspect the create geotiff, but at my place it looks
good. I've also checked after transforming to UTM32N (EPSG 25832) and it
looks similar, no shifts.

filename = 'DE1200_RV2209200825.tar.bz2'
with tarfile.open(filename, "r:bz2") as tar:
for fname in tar.getnames():
f = tar.extractfile(fname)
data_raw, meta = wrl.io.read_radolan_composite(f)
xy_raw = wrl.georef.get_radolan_grid(*data_raw.shape)
data, xy = wrl.georef.set_raster_origin(data_raw, xy_raw, "upper")
proj_osr = wrl.georef.create_osr("dwd-radolan")

proj_25832 = wrl.georef.epsg_to_osr(25832)
xy = wrl.georef.reproject(xy, projection_source=proj_osr,
projection_target=proj_25832)
ds = wrl.georef.create_raster_dataset(data, xy,
projection=proj_25832)

wrl.io.write_raster_dataset(f"{fname}.tif", ds, "GTiff")

Maybe your ArcGIS can't work with the polar stereographic transform. Can
you check with a non-proprietary solution like QGIS?

HTH,
Kai


Am 22.09.22 um 15:45 schrieb lord rufusus:

[lord rufusus]

Hey Kai,

thx a lot. Using your code with the projection in ETRS it also works in ArcGis. :) No idea what I have done without you!!!

Falls du je in Oldenburg stranden solltest, sag Bescheid, Bier geht auf mich ;)

Lg Paul