Writing a Dfs2 file from python does not retain the x and y coordinates

[mlinds]

Describe the bug Hi mikeio community"

I was given a batch of Dfs2 files that contain an a projection strong, a correct origin point in the local CRS, and an orientation angle. However, They only have arbitrary coordinates as seen below:

I want to be able to interpolate variables based on local coordinates. So I went about trying to create a script that takes a dfs2 file as input, and writes a new one with the correct x and y values as output. Easy enough.

To Reproduce Steps to reproduce the behavior:

I am able to create new geometry with correct coordinates.

x0,y0 = dfs2file.geometry.origin
nx = dfs2file.geometry.nx
ny = dfs2file.geometry.ny
dx = dfs2file.geometry.dx
dy = dfs2file.geometry.dy
orient = dfs2file.geometry.orientation
new_geom = mikeio.spatial.Grid2D(x0=x0,y0=y0,nx=nx,ny=ny,dy=dy,dx=dx,orientation=orient,origin=dfs2file.geometry.origin,projection=dfs2file.geometry.projection)

After that I verify that it works:

Cool! Now, its time to set that geometry as the geometry of the Dfs2 file.

dfs2file.geometry = new_geom
# get the underlying dataset object
ds = dfs2file.read()

ds.to_dfs('my_corrected_file.dfs2')

Unfortunately, if I then load the new dfs and look at the x and y coordinates, they come out wrong.

Expected behavior I would expect that the newly-written dfs2 file has the same geometry object as the python object that was used to write it. Maybe I misunderstand the API design, but to me this seems counterintuitive to the point of being a bug.

System information:

• Python version: 3.11.6
• MIKE IO version: 1.6.3

[ecomodeller]

The read method reads the geometry from the file. You have to set the geometry on the dataset, before writing to a new file.

Give it a try!

[mlinds]

@ecomodeller Thanks for taking a look.

For me, the geometry of the dataset ds is correct when I check it. Also, I am not sure how to set the geometry of a dataset. If I try to set it by assignment it raises an attribute error. Is there a method I can use to set the geometry?

[ecomodeller]

@ecomodeller Thanks for taking a look.

For me, the geometry of the dataset ds is correct when I check it. Also, I am not sure how to set the geometry of a dataset. If I try to set it by assignment it raises an attribute error. Is there a method I can use to set the geometry?

True, the geometry has to be set on the DataArray. (Dataset is a collection of DataArrays)

Take a look at this material from the recent MIKE IO course https://dhi.github.io/getting-started-with-mikeio/dfs2.html That might help you understand the API.

[mlinds]

Thanks. I've looked at the individual DataArray and I'm running into what is seemingly the same problem. I set the geometry at the DataArray level in a loop. Then I create a new Dataset from scratch.

ds = dfs2file.read()

dalist = []
for da in ds:
    da.geometry = new_geom
    dalist.append(da)
new_ds = mikeio.Dataset(dalist)
# I can verify that all dataarrays in this new dataset have the correct coordinates
new_ds.to_dfs('fixed_dfs.dfs2')

I have checked all dataarrays in this dataset and they all have the correct coordinates. When I save the Dfs2 and then reopen it, the coordinates are one again arbitrary

Any ideas?

[ecomodeller]

Can you create minimal reproducible example with a 2x3 grid?

Or run the code in a debug and step through the lines of this function https://github.com/DHI/mikeio/blob/645213bb31940c3f94ede894057f489371e1141d/mikeio/dfs/_dfs2.py#L35

and see if you can spot the problem.

[ecomodeller]

@mlinds I am not sure if you still have a problem or not?