I am new to this area, but i am wondering if someone else has as well encountered the error message by plotting the scan_strategy-function?
Plotting the scan strategy with: ax = wradlib.vis.plot_scan_strategy(ranges, elevs, sitecoords, beamwidth=vol.how.beamwH) works fine. But if I want to use:
ax = wradlib.vis.plot_scan_strategy(ranges, elevs, sitecoords, cg=True, terrain=True) → I get: TypeError: download_srtm() takes 2 positional arguments but 3 were given.
Now, ax = wradlib.vis.plot_scan_strategy(ranges, elevs, sitecoords, beamwidth=vol.how.beamwH, cg=True) runs fine BUT the altitude is shifted by -1000 as MSL is now there.
everything works “fine”, but the MSL is still leveled at -1km. WHY? It is the same in the last figure of the wradlib-link you added above, and I don’t understand why it should be like this?
Furthermore, if i change az to something different than az=180, it gives me the following error message:
not sure, it would be great, if you can make a little MCVE from your example. So at least give the input data you have used. This will help to reproduce issue and find a way around it.
This might be a bug using units='km'. Does it work if you omit that kwarg? Another possible cause might be broken site_coords.
@lou Code looks good. The RuntimeError: No usable source images. might it have to do with the coverage of the SRTM data? IIRC it’s only covered to 60deg north.
Do you have a local DEM covering whole of Finland? We might just skip the automated SRTM download and prepare the terrain ourselves from local DEM data. If you need help with that you might have a glimpse at Beam Blockage Calculation using a DEM — wradlib how to interpolate a given DEM
onto the radar polar layout.
The other thing, maxalt would still have to be given in meters (even if you specify ‘km’).
The altitude-problem is probably a bug in the curvelinear grid code, as the issue also shows up in the docs. Good catch! Would you mind opening an issue on the wradlib github issue tracker Issues · wradlib/wradlib · GitHub?
Confirmed, it’s a bug, we get this offset of 1000m because we use two different earth radius in calculation and plotting (6370000 and 6371000). If you feel capable enough you can change the values to 6371000 as a quick fix:
Indeed, I forgot about that the SRTM data covers only the area up to 60deg north. I don’t have a local DEM covering Finland - I’ll have a look if I can find one and implement it.
Good to know about maxalt!
I think I don’t feel capable of fixing this in the code as I don’t know what the right way would be to do so.
I downloaded some Geo-Tiff files that I need. But I could not figure out how I should now process these data in order to implement it into the scan_strategy-function? Would you mind to help me out?
And, @kmuehlbauer, is there a way to add a star or a dot on the MSL (e.g. a location that is 60km away from the radar station) with some text?
ax = wradlib.vis.plot_scan_strategy(ranges,
elevs,
sitecoords,
maxalt=15000.0,
cg=True,
az=10,
cmap="viridis",
units='km',
terrain=False)
# get polar parasite axis
paax = ax.parasites[1]
# use this as earth radius
er = 6371000
# the distance in m
rng = 60000
# calculate bin distance
ade = wradlib.georef.bin_distance(rng, 0, sitecoords[2], ke=1.0)
# get angle
thetap = -np.degrees(ade / er) + 90.0
# plot the stuff at the wanted location
paax.plot(thetap, er, "r*", linewidth=3, label="RAD", zorder=3)
paax.text(thetap, er, "Radar", zorder=3)
Takes a while, need to find a bit of time doing that. You might explore GDAL merge gdal_merge.py — GDAL documentation on how to merge the tiles. After that just follow-up with the preprocessing steps Beam Blockage Calculation using a DEM — wradlib how to interpolate the DEM onto the radar polar grid. Then just feed the wanted azimuth to terrain (ax = wrl.vis.plot_scan_strategy(ranges, elevs, site, units="km", terrain=polarvalues[wanted_azi]).
