I also have a concern that if the domain is built too large, there may be issues with the conversion from x, y to lat and lon, as I will use [md. mesh. lat md. mesh. long]=xy2ll (md. mesh. x, md. mesh. y, -1);
Two functions, md. geometry. bed=InterpFromGridToMesh (lat, lon, topg, md. mesh. lat, md. mesh. long, 0)