Develop additional binning kernels

[eguil]

First implementation of binDensity uses a linear interpolation (numpy.interp). This has two limitations:

• high CPU (loop on lat/lon required)
• less control on binning of extensive variables (e.g. transport).

The goal here is to add other binning kernels and limits the use of loops.

[durack1]

@eguil while we're also thinking about performance, the actual vertical/cast interpolation should also be considered, there are now quite a number of available methods - see https://github.com/eguil/Density_bining/issues/42#issuecomment-327609810

[eguil]

Maybe we can use ESMF for 4D interpolation: https://github.com/raphaeldussin/ESMPy3/blob/master/examples/ungridded_dimension_regrid.py Any thoughts Paul ?

[durack1]

@eguil the dev folks are currently working on porting CDMS2 to python3, and as part of this testing of the ESMF 7.x branches are underway. We should be able to engage the ESMF dev folks directly on this question, it's an expansion of functionality that would benefit all cdms users..

@dnadeau4 @doutriaux1 pinging you both here

[doutriaux1]

@durack1 what about Python3? I do not see the connection.

[durack1]

@doutriaux1 the ESMF 6 -> 7 update was what I was attempting to point out.. This is all being done as an important part of the cdms move from python2 -> python3. And as @dnadeau4 is interacting with the ESMF folks on that work, it makes sense that we could also engage with them on this 4D interpolation question

[doutriaux1]

ah ok. ESMF/ESMPY 7.0 is already in. Since 2.10

[eguil]

Exchange with one ESMF developer:

Dear Eric,

I haven't used ESMF remapping for vertical interpolation, although it seems to be able to handle 3D grids.

Unfortunately the documentation is lacking some precious information. One way of doing it faster and at relatively low coding time would be to code or copy your favorite cubic spline or higher order scheme in fortran and use f2py to call it from python. I had tremendous speedups using this trick, especially when loops are involved.

Best, Raphael Dussin (raphael.dussin@rutgers.edu)

---

Dear Raphael,

A group of us are trying to optimize a code which bins 4D ocean data on density (going from [z,y,z,t] to [x,y,Sigma,t]). At the moment we are doing a horribly expensive loop on [x,y] using numpy.interp which eats 80% of the elapse time of the program (see line 675-679 in https://github.com/eguil/Density_bining/blob/master/binDensity.py).

We already use ESMFpy for horizontal interpolation and are exploring to see how we could also use ESMF (and its parallel capabilities) for the vertical interpolation. It seems this program would take us in the right direction: https://github.com/raphaeldussin/ESMPy3/blob/master/examples/ungridded_dimension_regrid.py but before pursuing we want to ensure we are on the right track.

Any help would be much welcomed ! Thanks in advance and best wishes, Eric Guilyardi

[doutriaux1]

@nadeau4 i think you have done 3D regrid right?

[dnadeau4]

yes.

[durack1]

@dnadeau4 is it fast?

[dnadeau4]

As fast as ESMF can be.

[durack1]

@dnadeau4 ok I shall pose my question better. How is the performance of the 4D (or 3D) interpolation in your user experience, is it fast and efficient, or relatively slow and resource hungry? Are there other options that you can contrast ESMF against?

[durack1]

@eguil this comment was interesting mom-ocean/MOM5#202 - looks like there are other folks interested in this

[eguil]

added function binIndices in libDensity.py and testing on simple cases using test_binnoloop.py locally