Finding whether set of coordinates are within a bounding box set by a condition

Hi,

I have two meshgrids of X and Y coordinates (UTM/polar stereographic) with a corresponding grid of elevation values.

I have a set of specific coordinates marking sites of interest and I'd like to find out if they fall within a particular elevation band.

E.g.

index1 = elevation > 3000

and I want to find if any of the sites fall in the coordinates specified by the elevation. The specific points will not correspond exactly to the points on the meshgrid.

I've tried inpolygon, filterm, etc to no avail.

Update

Added sample data: etopo_15.mat loads LAT, LON and topo_bed (elevation) matrices.

You can plot this using

figure();
pcolor(LON,LAT,topo_bed)
shading flat

% And let's say the coordinates of the sites I'm interested in are
coords = [2.43875e+02,-48.375;2.45625e+02,44.87540.125,5.375;2.47125e+02,10.625]%[LON1,LAT1;LON2,LAT2...]

I create elevation bands from

elev1 = find(topo_bed > 3000);
elev2 = find(topo_bed > 1000 & topo_bed<2000);

How do I find out if the sites in coords are in those elevations?

I think one method might be to make a polygon from the elevation criteria, and the use inpolygon? Having trouble making a polygon from the matrix though..

[darova]

Please attach the data. Can you show what are trying to do (image/picture or scheme)?

[darova]

I painted in red the region you are interested in. Is it correct?

[Elizabeth Case]

Hi darova,

The image is really small, I can't see what you've painted, could you repost it? Thanks!

[darova]

Please use right mouse button

Answer 1

[Guillaume]

Can't you just create an interpolant for your elevation and query the interpolant at your desired coordinates. It's then trivial to check if the elevation is in the desired bands:

%note that griddedInterpolant want NDGRID'ed data, not MESHGRID'ed data, so everything has to be transposed
ElevInterp = griddedInterpolant(LON.', LAT.', topo_bed.');  %using default linear interpolant here. Use whichever interpolation method you prefer

siteelevation = ElevInterp(coords)  %get elevation at site coordinates

isinband = siteelevation > 3000 | (siteelevation > 1000 & siteelevation < 2000)

[Elizabeth Case]

Hi Guillaume, this seems like a good option. I'm getting an error in gridded interopolant that might grids aren't monotonically increasing, but as far as I know/can tell, they are (see attached imshow). I'll try something along these lines and follow up here.

update: the issue is that my y (polar stereographic) values are decreasing

[Guillaume]

As long as the grid vectors are strictly monotonic, griddedInterpolant should work, it doesn't matter if it's monotonically increasing or decreasing. Certainly, I have no issue with the LAT and LON in your attached mat file.

If the coordinates are not monotonic or not gridded, you can always use scatteredInterpolant instead. The syntax is slightly different from that of griddedInterpolant but the principle is the same.