unexpected behavior for add_contour method

[simonweppe]

Hi @SorooshMani-NOAA ,

Im running into a strange issue when using the add_contour method.

I have defined a shape file that I want to use as my domain, and have tif file with bathy information from which I want to define my size function.
See code below:

hfun_rasters = [Raster(p) for p in dem_paths] #  if memory issues chunk_size=500, overlap=50
size_mesh_min = 100
size_mesh_max = 1000
exp_rate = 0.005

hfun = Hfun(
    deepcopy(hfun_rasters),
    base_shape=base_gs.unary_union,
    base_shape_crs=base_gs.crs,
    hmin=size_mesh_min, 
    hmax=size_mesh_max,
    nprocs = 1, 
    #method='fast'
    )

hfun.add_constant_value(value=size_mesh_max) 
hfun.add_constant_value(value=size_mesh_min, lower_bound=-3, upper_bound=3) # constant element size for size function between lower_bound and upper bound
hfun.add_contour(level=-3, expansion_rate= exp_rate , target_size=size_mesh_min) # element size starts at 100m at depth = -3 below sea level, then expands

Im getting a mesh that is looking good and as expected, expect near the outer ocean boundary where it's refined to size_mesh_min. This seems to be due to the add_contour step, which specifies a size of 100 along the outer edge of my mesh, while it should relax to much larger element size given it's deep through there. As if the computed -3 contour was going around the outer edge of my grid (it's not the case in the tiff).
I've tried a few things and can get more relaxed resolution there using other approaches (e.g. add_topo_bound_constraint, add_subtidal_flow_limiter), but as soon as I use add_contour, I get that odd high res zone.

Am I missing something obvious ? See some snapshots below.

[SorooshMani-NOAA]

Hi, I think something has changed in matplotlib recently. Can you try using matplotlib version 3.5.2 and see if this resolves the issue or not?

[SorooshMani-NOAA]

@felicio93 do you think this is related to #144? Have you had the chance to explore that other ticket?

[simonweppe]

Thanks for the fast reply @SorooshMani-NOAA ..yes I saw this one other ticket but Im using matplotlib 3.5.1 so shouldnt be the issue.

[SorooshMani-NOAA]

Another possibility is that when the base shape "cuts" the DEM it creates nan regions. I remember seeing something similar to what you see a couple of years ago when I was implementing this. I don't remember if I ever had a good solution, something like this might work:

Try using the baseshape to create your geometry, but buffer the base shape before passing it to create hfun object. This way, the invalid contours due to cutting will fall outside the domain. Note that you need to play with the buffer size to get a good results with this approach.

[felicio93]

@SorooshMani-NOAA and @simonweppe,

I agree with I don't think this is related to the matlab problem, at least it is not the same problem.

I tried recreating 2 meshes using the same script I was using before (both cases worked). This is what I get now:

case1 - CA Mesh:

case2 - AK Mesh:

@SorooshMani-NOAA, do you know if something has changed with our dependencies?

[felicio93]

FYI, the problem here is not the same of #144. I just tested and the ocsmesh "get_contour" raster function is working as expected.

[SorooshMani-NOAA]

@felicio93 thanks for testing, I think the problem might then be what I described in #154 (comment)

[simonweppe]

Thanks @SorooshMani-NOAA @felicio93 I will try your workaround

[felicio93]

Yes, you are right!

I tested rerunning with the buffer as:

hfun = Hfun(
    hfun_rast_list,
    base_shape=domain.buffer(0.1).unary_union,
    base_shape_crs=geom.crs,
    hmin=100, hmax=2000,
    method='fast')

For the sake of documentation, this is how this looks like:

[simonweppe]

Great - thanks for the quick replies and support on this.

[SorooshMani-NOAA]

Can we close this ticket now?

[simonweppe]

Yes you can, sorry I should have closed it back then - thanks for the help on this!

…

On Mon, 22 Jul 2024, 16:14 Felicio Cassalho, ***@***.***> wrote: Closed #154 <#154> as completed. — Reply to this email directly, view it on GitHub <#154 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ABXQTUZYQQT7BVSH3SRUXYDZNUHV5AVCNFSM6AAAAABJGFZLSCVHI2DSMVQWIX3LMV45UABCJFZXG5LFIV3GK3TUJZXXI2LGNFRWC5DJN5XDWMJTGYYDEMBVHAYDMNY> . You are receiving this because you were mentioned.Message ID: ***@***.***>