I still get the wrong result like:

grdtrack [WARNING]: Some input points were outside the grid domain(s).
          0        1            2
0 -110.9536 -42.2489 -2797.394987
1 -111.3500 -42.2033 -2715.359945
2 -111.7309 -43.1706 -2711.704748
3 -112.2196 -44.3819 -2787.948153
4 -112.5505 -45.4878 -2746.303734
5 -112.6736 -45.8347 -2685.852247
6 -112.4235 -45.8561 -2785.099618
7 -112.9507 -47.2270 -2672.932596
8 -113.3124 -48.3208 -2706.345811
grdtrack [WARNING]: Some input points were outside the grid domain(s).
          0        1            2
0 -110.9536 -42.2489 -2974.656296
1 -111.3500 -42.2033 -2859.961069
2 -111.7309 -43.1706 -2711.704748
3 -112.2196 -44.3819 -2787.948153
4 -112.5505 -45.4878 -2746.303734
5 -112.6736 -45.8347 -2685.852247
6 -112.4235 -45.8561 -2785.099618
7 -112.9507 -47.2270 -2672.932596
8 -113.3124 -48.3208 -2798.693406

OK, I will have to debug that again, but running into dinner bbq. Will work on this branch later.

Oops, let's try MIN this time.

Should work now since you got -nl to work, @seisman

I have double-checked that I'm using the latest commit in this branch, but I still get the wrong result.

And I just debugged and it still says interpolation = 3 after the MIN calls with 1 in it...
RUnning out of time before dinner now, sorry.

OK, now it should work.

I'm glad that I'm no longer the only one claiming at the pad but I thought the answer was grids-must-have-a-pad.
Perhaps a better solution would be to lower down to bi-linear only when the track touches the grid borders.

The pad was introduced decades before there was a plan to have an API for externals to use. The pad made all sorts of operations requiring interpolation and sampling simple. Of course, the pad concepts collides with what externals do, and it is causing problems like the one here as well as similar situations.

The pad is especially valuable when we are working on a subset of a larger grid, since the boundary pads are actual boundary values and not derived from some mathematical condition. So when externals do a region cut and pass in a pad-less grid they have already lost that data. Oh well, GMT cannot do anything about that.

The only real solutions I see are these:

1. Expand those grids to have pads, i.e., allocate new memory and set the BCs. This annoys the memory-limited folks.
2. Rewrite those sections of GMT that actually deals with pads. I think this involves gmt_grd_project, gmt_img_project, the BCR machinery, and many functions in grdmath, for sure. There may be more; getting a complete list will be important.

Option 1 is what we currently do in general and it works fine. It allows us to preserve data BCs when they exist. Option 2 could be considered, but I think we have to do that in a clever way that does not regress to the point of not using data BCs. I will have to think about that some more but I can imagine storing (or computing) the pad as two small mini-grids of 4mx and 4my in sizes and have the bookkeeping access the right pad nodes when we exceed the interior grid during calculations. In principle that is no different than doing different things near the border since it means we will need if-tests, and bye-bye to OpenMP accelerations.

The big problem with option 2, and I have always acknowledged that, is the work that it implies. But I don't see why the bye-bye. If the loops are done from row, col = 2:end-2 and next the four edges the OMP will not know.

Because of the pad I was not able to fix some of the bugs we have open when dealing with sub-regions of images.

You are right, the inner grid loop can be OMP. But I think we agree that this would require lots of work (by whom?) vs a doubling in memory (duplicating a grid on input) for externals [which I think is required anyway in Matlab]. I just don't see the motivation. I have yet to have any tool tell me that you are out of memory on any reasonable project.
BTW, I started to look at where pads are used and the list is more exhaustive than I feared. Lots of modules operate due to pads, such as grdcut, grdpaste, and the list goes on. In addition to gmt_*.c of course.

I think option 1 via the GMT_GRID_NEEDS_PADS is the way to go, and then we can tell our grandchildren that there were actual people back then who worried about computer memory before the infiniti-crystals become common. Of I suppose it could be GMT_GRID_PAD_ONE and GMT_GRID_PAD_TWO if we want to save a tiny amount in grdgradient.

I'm fine in living it as is. And it was not me who brought the pad issue this time.

if we want to save a tiny amount in grdgradient

but again, it's not that memory that counts. It's the fact that one are obliged to duplicate the array in and out.