DSFU data visualization on Google Map using MIKEIO

[saiyara11]

Hi,
I am trying to create a map with dfsu result file and want to have google street map as background using MIKE IO. I am facing problem with the projection. My dfsu result file always generates in BTM projection system. Is there any way to change the projection of dfsu file to WGS84 programmatically to keep it similar to the map?

[ecomodeller]

@clemenscremer do you have any suggestion on how to get a correctly aligned background map with Dfsu data in projected coordinates?

[saiyara11]

Thank you so much for your response. I have followed your code but getting an error CRSError: The WKT could not be parsed. OGR Error code 6
Can you help me with this?

[clemenscremer]

You are very welcome. This looks like a problem with rasterio but is hard to diagnose without any further information.
I have just run the code in another one of my environments and it works flawlessly using rasterio v. 1.3.8 and contextily v. 1.3.0.

I would suggest to try upgrading rasterio and see if this fixes your error.

[saiyara11]

I have upgraded rasterio to 1.3.8 but not still not working. I have added the code and an attachment with the screenshot of the problem.

` import mikeio
import contextily as ctx
import matplotlib.pyplot as plt
import rasterio
def find_epsg(dfs):
"""finds epsg code in mike dataset"""
import re
proj_string = dfs.geometry.projection
match = re.search(r'EPSG",(\d+)', proj_string)
if match:
epsg = int(match.group(1))
else:
epsg = None

return epsg

fn_dfsu = r"C:\Users\FFWS\Desktop\Boundary_Generation_For_2D\V1054_SWRM_SW_2D PS_LR_Dike_Updt_statBaseFlood
Map_Day05.dfsu"
zoom = 12

dfs = mikeio.read(fn_dfsu, items="Maximum water depth")

fig, ax = plt.subplots(1,1, figsize=(14,8))
dfs[0].plot(cmap="hsv", ax=ax)
source = ctx.providers.Esri.WorldImagery # ctx.providers.OpenStreetMap.Mapnik
epsg = find_epsg(dfs) # find epsg code in mike dataset
ctx.add_basemap(ax, crs=f"epsg:{epsg}", source=source, zoom=zoom)`

[clemenscremer]

Might be a problem with the regex or the epsg then (just have added it for convenience and it works with my dfsu files but it might not generalize well). I suggest you specify epsg directly instead of using find_epsg(dfs): epsg = 3857.

If that doesn't work you might need to read on reprojections https://pygis.io/docs/d_raster_crs_intro.html

[ecomodeller]

https://spatialreference.org/ is a good resource.

[saiyara11]

Thank you all for your support. I have successfully added the background map. I have just added epsg=9678 to the code and it worked fine. thanks again.

[clemenscremer]

hey @saiyara11 : great to hear! congrats and thanks for letting us know about that fix :)