ROMS and wetting/drying masks

[kthyng]

Hi! I need to override the constant values set in consts.py for "sea_floor_depth_below_sea_level" because in ROMS the depths are negative instead of positive and the default range defined is 0 to 12000. I have tried a few approaches for this, most recently:

Initial setup:

reader = reader_ROMS_native.Reader(loc)
ocean_drift_sim = OceanDrift(loglevel=20)
ocean_drift_sim.add_reader([reader])

then:

key = "sea_floor_depth_below_sea_level"
c = {}
c[f'environment:fallback:{key}'] = ocean_drift_sim.get_configspec()[f'environment:constant:{key}']
c[f'environment:fallback:{key}'].update({"min": -12000, "max": 0})
ocean_drift_sim._add_config(c)
ocean_drift_sim.set_config(f'environment:fallback:{key}', -10000)

Also, separately, tried this:

ocean_drift_sim.get_configspec()['environment:constant:sea_floor_depth_below_sea_level']["min"] = -12000
ocean_drift_sim.get_configspec()['environment:constant:sea_floor_depth_below_sea_level']["max"] = 0

But when I actually then seeded the simulation and ran it, I still got log messages like:

11:22:54 WARNING opendrift.readers.basereader.variables:635: Invalid values (-6.475353 to -0.3697905) found for sea_floor_depth_below_sea_level, replacing with NaN
11:22:54 WARNING opendrift.readers.basereader.variables:638: (allowed range: [0, 12000])

or in the case that I changed the fallback value:

16:26:17 WARNING opendrift.readers.basereader.variables:635: Invalid values (-10000.0 to -10000.0) found for sea_floor_depth_below_sea_level, replacing with NaN
16:26:17 WARNING opendrift.readers.basereader.variables:638: (allowed range: [0, 12000])

So, neither of these approaches seemed to "stick" and it seems that what is hard-wired into "consts.py" isn't meant to be overridden. Any suggestions for how to proceed? Thanks.

[knutfrode]

I believe it should be better to update the reader to provide positive depths, as this is what is expected by various drift modules (although this is just a convention). In all ROMS native files I have seen so far, depth is positive. But if it is negative in your case, we could update the reader to reverse sign. Maybe there are some attributes to the variable to identify that it is negative?

If you could share a small sample file, I could have a look at it. Sharing output from ncdump could also be useful.

[kthyng]

As you probably know, ROMS has s_rho and s_w stretching vertical coordinates which have to be mapped to depth based on some of the grid details and the sea surface height. I guess it is a choice at that mapping stage as to whether to map them to positive or negative values, but my impression is they are usually mapped to negative values since s_rho and s_w have the "positive: up" attribute and themselves have negative values. You can look at a representation of the xarray Dataset summary for an example dataset here: http://xpublish-nwgoa.srv.axds.co/datasets/nwgoa_all/.

And I guess this must be happening in the ROMS reader itself, but I didn't change anything in there to change it to be positive or negative in particular. I need to look a little more closely at that code.

[kthyng]

The ROMS reader does seem to think the depths should be negative, for example:

https://github.com/OpenDrift/opendrift/blob/15546beb529f0d94629dc971a00b6edcb80dbf2b/opendrift/readers/reader_ROMS_native.py#L74-L78

unless I am over-interpreting those lines as being more meaningful than they are.

[knutfrode]

Yes, the vertical coordinate z is always negative (by convention in OpenDrift), and when seeding, you should use negative z values.

However, the variable sea_floor_depth_below_sea_level is the absolute value of the depth below sea level. Although the CF standard name table does not say explicitly, I believe this always have to be a positive value. At least I have not seen any cases or model outputs where it is negative.

[knutfrode]

z has to be negative inn the ocean, as OpenDrift also supports drift in air. Thus z=-10 is in ocean, and z=10 means that particle will follow wind.

[kthyng]

@knutfrode Ah I realize now I conflated the depths of the drifter with the overall seafloor depth. Thank you for that clarification! I dug in and see the reason for some seafloor depths being negative: this model was run with wetting and drying on. Negative grid cells are where wetting and drying is allowed to occur, and then also in permanently masked areas (those all have a fixed negative value). Has opendrift been used in models/areas with wetting and drying?

[kthyng]

Just wanted to follow up here — any advice on using OpenDrift with wetting and drying model output? Is there any special handling in this case? I can imagine if the cell is wet or dry then a tracking algorithm would just work as usual, but might need some special care in the case that a cell becomes dry underneath a drifter — does it become stranded in that case?

[knutfrode]

Coastline stranding is determined by the binary landmask (land_binary_mask) and not the ocean depth variable (sea_floor_depth_below_sea_level) although the two are of course related. In principle there should be nothing wrong with a time-varying landmask, but in the ROMS reader it is in fact hardcoded that landmask is assumed to be static (2D), and read only at first timestep, for performance reasons. This may however be changed if needed.

With the config setting general:coastline_action you can chose whether particles should be deactivated (stranded) or moved back to previous location when hitting coastline. https://opendrift.github.io/gallery/example_coastline_options.html

However, you'd might here instead want particles to "wait" on the dried location, and then one could (internally) instead set self.elements.moving to 0 for the elements where depth is negative, similar to what is done in the SedimentDrift module: https://opendrift.github.io/_modules/opendrift/models/sedimentdrift.html#SedimentDrift This could well be added as a general functionality in OceanDrift, perhaps with a config switch to allow deactivating this feature. I could have a look at this next week, if it is of interest. With this solution, one can keep the assumption in ROMS reader that land_binary_mask is 2D/static, as sea_floor_depth_below_sea_level is used instead.

[kthyng]

Coastline stranding is determined by the binary landmask (land_binary_mask) and not the ocean depth variable (sea_floor_depth_below_sea_level)

Ok that makes sense, though what about in the scenario that the following warning occurs which says depths have been changed to nan — does that add the nan-ed location to the landmask in effect? It looks to me from some cases I'm running and plotting that these grid cells with negative depth are in fact still being used but the warning makes it sound otherwise.

11:27:05 INFO    opendrift.models.basemodel:1989: 1999-02-01 02:00:00 - step 3 of 4 - 100 active elements (0 deactivated)
11:27:32 WARNING opendrift.readers.basereader.variables:635: Invalid values (-11.243423 to -0.09407605) found for sea_floor_depth_below_sea_level, replacing with NaN

However, you'd might here instead want particles to "wait" on the dried location, and then one could (internally) instead set self.elements.moving to 0 for the elements where depth is negative ... With this solution, one can keep the assumption in ROMS reader that land_binary_mask is 2D/static, as sea_floor_depth_below_sea_level is used instead

The variable sea_floor_depth_below_sea_level is static in time. mask_rho shows the cells that are never active as land and everything that could ever be active as active. sea_floor_depth_below_sea_level has negative depths where wetting and drying can occur, but it depends on the time-varying sea surface height on top of sea_floor_depth_below_sea_level for whether a location is wet or dry. So to try the functionality in the sediment module, which I think is a good idea, I would want to check the value of sea_floor_depth_below_sea_level + sea_surface_height to see if a cell is wet or dry.

[knutfrode]

There is a "safety-mechanism" in OpenDrift to ensure that environment variables have reasonable values, e.g. avoiding problems with corrupted file in an operational setup. sea_floor_depth_below_sea_level is expected to be between 0 and 12000: https://github.com/OpenDrift/opendrift/blob/master/opendrift/readers/basereader/consts.py#L20

One could here set the minimum value to e.g. -10 (?) to allow for drying cells.

It should be fine to check sea_floor_depth_below_sea_level + sea_surface_height, but the latter must then also be added to required variables in the actual module. It may however have 0 as a default value to be used whenever is is not available from any readers.

[kthyng]

As a follow up, I used the wet/dry mask as the land binary mask and allowed it to vary in time to solve my problem. Thanks for your guidance.

[knutfrode]

Very good!

[kthyng]

@knutfrode Well I have fixed some issues regarding the wet/dry mask. However, I am diving deeper now and I see more items. As part of this, I see that in fact I have not fixed the issue of the cells being nan'ed out. You suggest

One could here set the minimum value to e.g. -10 (?) to allow for drying cells.

which indeed is what I would like to do. Is there any way to do this without modifying the opendrift code itself?

[knutfrode]

Perhaps you can set the fallback value for seafloor depth to 0: >>> o.set_config('environment:fallback:sea_floor_depth_below_sea_level', 0) (This can also be hardcoded within your module, by setting fallback to 0 for sea_floor_depth_below_sea_level in required_variablesdictionary)

Then self.environment.sea_floor_depth_below_sea_level will be set to 0 whenever an "illegal value" (i.e. <0) is encountered. And then you could use this as the condition for a drying cell in your module. You will then also get 0 if it becomes larger than maximum (12000m), but you can probably assume that this is not going to happen (unless a corrupted file).

[kthyng]

It looks like the "bad" values are always replaced with nan regardless of the fallback value:

https://github.com/OpenDrift/opendrift/blob/7aae898f383052e803e28aef6d32434327f4520d/opendrift/readers/basereader/variables.py#L642-L650

If we were to change the minimum depth value to -10 to allow for drying cells (in consts.py) would that allow the drifters to move freely in the cells that are temporarily dry? Or is the depth used for some part of the drifter movement such that having a negative depth would cause an issue?

If that would be an issue, the only thing I can think of is making the depth also vary in time (like the mask) and have it update with the water column depth (h+zeta) each time step. I might be able to make that change outside of opendrift and input the modified Dataset to opendrift (with the new version that allows for inputting Datasets to the ROMS reader) but would still need to modify the ROMS reader to allow for time-varying depth.

[knutfrode]

I believe thse NaN-values are later replaced with the fallback value? Otherwise I don't think negative ocean depth would make a problem somewhere, in contrary to positive z-value (particle depth, or rather height (positive upwards), not to be confuse with seafloor depth) for which there are some checks in some of the modules.

It would be a pity having to physically modifying the input forcing fields, as it should be possible to fix somehow with the code. Can you try to set the valid_min to -10 here: https://github.com/OpenDrift/opendrift/blob/master/opendrift/readers/basereader/consts.py#L20 If that works, I can see if there is (or make) a way to do this without needing to modify that file.

[kthyng]

@knutfrode Yes changing the consts.py code does fix the issue of the warning that my wet/dry cells have been nan'ed out.

    'sea_floor_depth_below_sea_level': {'valid_min': -10, 'valid_max': 12000,
                         'long_name': 'Depth of seafloor'},

With that change, I don't get any more warnings like

12:13:24 WARNING opendrift.readers.basereader.variables:635: Invalid values (-2.4052637 to -2.4052637) found for sea_floor_depth_below_sea_level, replacing with NaN

[knutfrode]

Good. I now made this a permanent change in OpenDrift, as it does not seem to interfer with anything else: https://github.com/OpenDrift/opendrift/commit/debdf6b53e6ef1fb6a7f7340f368de43071c2f57

[knutfrode]

Does this solve the problem, or are there more issues related to this?

[kthyng]

@knutfrode Yes thank you, this resolves this issue. I plan to submit a PR for the ROMS reader with my changes for using the wet/dry mask as well as a few other things. Still working on it for now though. Thanks for your help!

[kthyng]

I am afraid this issue of wetting/drying has again become more complex than I realized.

I was able to incorporate the wet/dry mask into the ROMS reader and you fixed the check that nan'ed out the negative depths that are grid cells that are allowed to wet and dry (actually the depth end point needs to be more like -20 instead of -10 in consts.py).

However, I hadn't noticed that when drifters go into the wet/dry cells, they are found to be under the seafloor and given the depth of the water column. I have not been able to find a satisfactory solution to this issue. I tried but failed to think of a way to differentiate between this wet/dry tidal flat situation and the situation opendrift is trying to avoid of a drifter being accidentally below the seafloor:

https://github.com/OpenDrift/opendrift/blob/668986480358eb55b211f80f1c27f2e4c17785f8/opendrift/models/basemodel/__init__.py#L697-L698

I also tried deactivating the drifters in the tidal flat when they get stranded like @knutfrode suggested above (with the thinking one of us might implement a "wait"-style command) and that worked in the sense that drifters that ran into the tidal flats (i.e. wet/dry cells) stopped there. But then I ran drifters near the seabed with the same flag (which I would have on generally if I was using it for 3D simulations in which the drifters might reach the tidal flats) and the drifters stopped if they happened to pass the seafloor (not the tidal flats) — I guess the math allows them to do so and this check is what corrects that behavior? So this seems like I cannot rely on using it for the tidal flats if it will also stop drifters at the seabed when they are near the seabed.

Do you have any corrections to what I have seen so far, or any other ideas?

Edited to add: I have realized some mistakes in my thinking since I wrote this and am working on it; I'll report back soon.

[knutfrode]

Closed by https://github.com/OpenDrift/opendrift/pull/1264