Integral field spectroscopy can map astronomical objects spatially and spectroscopically. Due to instrumental and atmospheric effects, it is common for integral field instruments to yield a sampling of the sky image that is both irregular and wavelength-dependent. Most subsequent analysis procedures require a regular, wavelength independent sampling (for example a fixed rectangular grid), and thus an initial step of fundamental importance is to resample the data onto a new grid. The best possible resampling would produce a well-sampled image, with a resolution equal to that imposed by the intrinsic spatial resolution of the instrument, telescope, and atmosphere, and with no statistical correlations between neighboring pixels. A standard method in the field to produce a regular set of samples from an irregular set of samples is Shepards method, but Shepards method typically yields images with a degraded resolution and large statistical correlations between pixels. Here we introduce a new method, which improves on Shepards method in both these respects. We apply this method to data from the Mapping Nearby Galaxies at Apache Point Observatory survey, part of Sloan Digital Sky Survey IV, demonstrating a full-width half maximum close to that of the intrinsic resolution (and ~ 16% better than Shepards method) and low statistical correlations between pixels. These results nearly achieve the ideal resampling. This method can have broader applications to other integral field data sets and to other astronomical data sets (such as dithered images) with irregular sampling.