Resampling using interpolations other than bilinear

[andrevitorelli]

Describe the feature and the current behavior/state.
Addition of interpolators beyond bilinear in resampler_ops.

As detailed in #2535, gradients of the bilinear interpolation fail when applied to an integer-pixel warp. This is expected, as in bilinear interpolation, pixel positions will correspond to non-differentiable corners in the interpolated function. This is an issue when using tensorflow addons for our(@andrevitorelli, @EiffL, @b-remy, @aguinot, @Dr-Zero) projects, at CosmoStat/autometacal and DifferentiableUniverseInitiative/GalFlow to implement the same techniques as esheldon/ngmix with TensorFlow.

Relevant information

• Are you willing to contribute to it (yes/no): Yes. We have code running and we'd like to open the PR, from andrevitorelli/addons/new_kernels
• Are you willing to maintain it going forward? (yes/no): Yes
• Is there a relevant academic paper? (if so, where): Yes. The most relevant paper is the GalSim paper, the framework on top of which ngmix works. Our most used interpolation kernel (quintic) comes from Bernstein & Gruen 2014
• Does the relevant academic paper exceed 50 citations? (yes/no): Yes.
• Is there already an implementation in another framework? (if so, where): Yes, in galsim-developers/GalSim, without automatic differentiation or GPU support.
• Was it part of tf.contrib? (if so, where): No

Which API type would this fall under (layer, metric, optimizer, etc.)
custom_ops
Who will benefit with this feature?
Anyone using resampling or affine transformations of images through resampling.
Any other info.
In our fork, we currently have implementations of cubic and quintic interpolations, though we would like to have more in the future (such as Lanczos).