Add equality constraints for predicted data

[Esteban82]

Note: edited description after discussion

It would be helpful if we could add equality constraints for the predicted data when fitting a model. For example, when gridding bathymetry we may want to constrain the values along a coastline to be exactly zero. This can be general and used to fix predicted values at any set of points.

This could probably be done with Lagrange multipliers to set the constraints:

The B matrix would be the spline Jacobian for the constrained points (the coastline) and v is a vector with the desired values at those points (zero for the bathymetry case).

---

Original question

When interpolating bathymetrical data (like from single-beam surveys) it is very useful to use the coastline as well. So my idea is to have a tool or option that automatically download the near coastline (maybe the GMT database could be use) and add it as data (with a very high weight?). This could be done for the sea and also for lakes/rivers.

I did it to create a bathymetrical map of the Lago Fagnano and it improve the results very well (i.e. the results were much logical from a geological point of view).

And maybe, just maybe, with a lot of caution. the topographic data from a DEM also could be used as an input to interpolate bathymetrical data. It is very logical that the shape of the topography goes on under the water (specially for lakes in mountains areas like the Fagnano).

If you are interested, I could share some of my GMT scripts (which are very basic) which I for the Fagnano.

Are you willing to help implement and maintain this feature? Yes

[leouieda]

Hi @Esteban82 I'm not entirely sure what you mean by using the coastlines for the interpolation. Do you mean cutting the grid based on the coastline data (like our convexhull_mask)? If so, there is already this package for doing that: https://regionmask.readthedocs.io/en/stable/defined_landmask.html So you could generate the grid with Verde and cut out the lake.

[Esteban82]

Hi @leouieda. My idea is to use the coastline as a border condition to grid the data. For example, if I make a bathymetric grid of the coast of California the coastline could be used as an input and the grid should be forced to have values with 0 over the coastline. Am I clear?

[leouieda]

Ah, right! That sounds like a really good idea. More generally, this would require adding equality constraints to our gridding methods. Currently they don't support this and I'd have to think about what that would look like in the API (do we pass the constraints to fit or when we create the Spline?). The maths might be a bit more complicated as well but I'll have a look.

I'll rename this issue to reflect this. It could be generalized to any type of equality constraints, not just coastlines.

[Esteban82]

Ok, great. And also as I said before maybe also the data onshore (from a DEM) could be include. The idea is that the relief onshore should go on under the water. It is not very logical to have two different type of landform separated by the coastline.

In other words, if I have a very abrupt relief outside a lake for example, it could be expected (in general) to have abrupt relief inside the lake. And if the relief outside is very smooth, it should goes on inside the lake.

[leouieda]

In other words, if I have a very abrupt relief outside a lake for example, it could be expected (in general) to have abrupt relief inside the lake. And if the relief outside is very smooth, it should goes on inside the lake.

That is much easier to do, actually. You could add DEM points from outside to the bathymetry data before gridding. That would make the near shore predictions a lot better since they would have data from both sides.

[Esteban82]

Yes, it is much easier. Although it should be done with caution. There can be sedimentary process inside the water that change the landform.