  Numerical simulations of merging black hole binaries produce the most
accurate gravitational waveforms. The availability of hundreds of these
numerical relativity (NR) waveforms, often containing many higher spherical
harmonic modes, allows one to study many aspects of gravitational waves.
Amongst these are the response of data analysis pipelines, the calibration of
semi-analytical models, the building of reduced-order surrogates, the
estimation of the parameters of detected gravitational waves, and the
composition of public catalogs of NR waveform data. The large number of
generated NR waveforms consequently requires efficient data storage and
handling, especially since many more waveforms will be generated at an
increased rate in the forthcoming years. In addition, gravitational wave data
analyses often require the NR waveforms to be interpolated and uniformly
resampled at high sampling rates. Previously, this resulted in very large data
files (up to $\sim$ several GB) in memory-intensive operations, which is
unfeasible when confronted with hundreds of multi-modal NR waveforms. To handle
these challenges, we present a simple and efficient method to significantly
\emph{compress} the original waveform data sets while accurately reproducing
the original data via spline interpolation. The method is generically
applicable to relatively smooth, one-dimensional datasets and uses a greedy
algorithm to determine the most relevant subset of the full data such that a
spline interpolant of a specified polynomial degree will represent the original
data to within a requested point-wise tolerance. We find significant
compression of the original NR data sets presented here. These compressed data
sets can then be evaluated fast and efficiently and resampled as desired.
