more general ordering of input vertical grids in fast_interp3D

[NoeLahaye]

Hi I noticed something that I was not expecting: it seems like interp2z requires both vertical grids to be ordered the same way. It looks like this behaviour is inherited from fast_interp3D.interp routine. Would you mind check this?
https://nbviewer.jupyter.org/github/NoeLahaye/bacasable/blob/master/debug_interpz_order.ipynb
@apatlpo , @slgentil

What is strange is that in my examples, sometimes I get Nan (I would get wrong values if I use boundary extrapolation), while I never have any target value that is outside the range of initial values (for the z grid), so I hope there is not something worse behind this

All the best,

Noé

[apatlpo]

The code that performs the interpolation is fairly simple:

croco/crocosi/fast_interp3D.c

Lines 28 to 45 in 3200403

	for (k=0;k<NZT;k++)
	{
	if (ztk[k] <= zsk[0]){
	if (b_extrap == 1) vtk[k] = vsk[0];
	else if (b_extrap == 2) vtk[k] = lin_int(zsk[0],vsk[0],zsk[1],vsk[1],ztk[k]);
	else vtk[k] = NPY_NAN;
	}
	else if (ztk[k] >= zsk[NZ-1]){
	if (t_extrap == 1) vtk[k] = vsk[NZ-1];
	else if (t_extrap == 2) vtk[k] = lin_int(zsk[NZ-2],vsk[NZ-2],zsk[NZ-1],vsk[NZ-1],ztk[k]);
	else vtk[k] = NPY_NAN;
	}
	else {
	while (ztk[k] > zsk[k1]) k1++;
	k0 = k1-1;
	vtk[k] = lin_int(zsk[k0],vsk[k0],zsk[k1],vsk[k1],ztk[k]);
	}
	}

and consistent with the notebook.
The question that follows is then: what do we do about it?

Allowing for grids that are reversed is certainly feasible.
But do we want to also allow unsorted vertical grid points?

[apatlpo]

My preference would be for the former.

[NoeLahaye]

Ah OK! I could have read more carefully fast_insterp3D.c. Maybe the least we could do, if we don't want to modify the routine, is to add a mention in the docstrings or even add a test in interp2z?

[apatlpo]

Agreed, could you please push a commit updating the doc?

Could you also modify the title of this issue with something like: "more general orders of input vertical grids for fast_interp3D " ?

We'll list this issue as a new feature request to remember a PR should be proposed at some point in order to allow interpolation with vertical grids with opposite orders.

[NoeLahaye]

done

[slgentil]

There's a new method in xgcm to interpolate (see xgcm tutorial in crocosi : xgcmgrid.transform). The target coordinate must be monotonically increasing or decreasing. Performance seems quite good. May be it could help...
@NoeLahaye , @apatlpo