Benchmarking interpolation functions.

Interpolation functions help turn low resolution discretely sampled signals back into their continuous counterparts. In this document we will be benchmarking the performance of a wide range of interpolators from MusicDSP.org, Paul Bourke, and Olli Niemitalo.

Wikipedia
In signal processing, sampling is the reduction of a continuous-time signal to a discrete-time signal. A common example is the conversion of a sound wave (a continuous signal) to a sequence of samples (a discrete-time signal). A sample is a value or set of values at a point in time and/or space.

Essentially a computer will quantize a continuous-time domain signal by sampling it at a certain frequency known as the sampling rate, this leaves us with an array of points in memory that represent what was once a continuous signal now as a discretely sampled signal. If we wish to sample between the samples in the array we will need to interpolate between two given samples, commonly the most basic approach to this is to use linear interpolation which is basically a ramping function which ramps linearly from one point to the other. But, we can improve upon this by using more complicated polynomial methods of interpolation.

The application of these functions can range from audio applications, general digital signal processing such as communications, image filtering, and more.

This document has taken a number of advanced interpolation functions and benchmarked them, when it comes to the quality of these interpolators this is not a subject we will be touching upon in this document as it is generally already covered in the documents by the respective authors.

First of all, we have a bunch of random interpolation functions cobbled together from MusicDSP.com, a website that lacks somewhat in academic reputability however that is not to say it is all bad, albeit cobbled together snippets of code from the late ’90s to date of varying quality it resembles a DSP specific version of Rosetta Code, there are some real gems to be found here and some interesting insights spanning back through the years of DSP.

Second, we have a lovely selection from Paul Bourke, first published here. This includes the classic Linear Interpolation, Cosine Interpolation, Cubic/Catmull Interpolation, and Hermite Interpolation. A really great selection of interpolators to begin with.

And thirdly we have a very in-depth set of 51 different interpolators from Olli Niemitalo, these are very interesting interpolators, personally, some of my favourite since I first discovered them back in 2011, the original paper Olli published can be read over at his website here.

The benchmarks were performed using a 32-second window of subsequent executions measured in microsecond accuracy. Thus each figure represents the number of executions performed in that 32-second window.

On each iteration the rndData() function is called to randomise the array of data points the interpolators are sampling from, this helps to ensure that the compiler does not optimise away any vital code that would damage the validity of the final results.

Here are the benchmarks, -Ofast enabled…

Executions in 32 seconds (higher is better) …
Fastest to slowest …
LinearInterpolate()-- 1,043,021,450*- A CLASS
hermite4()----------- 864,143,335 
CubicInterpolate()--- 861,368,667 
hermite1()----------- 852,650,628 
watte_42z()---------- 849,978,125 *
parabolic_42x()------ 845,329,016 *
parabolic_42z()------ 843,051,682 * - B CLASS
optimal32_42z()------ 842,327,674 
optimal16_42z()------ 842,173,284 
lagrange_43x()------- 824,777,726 
optimal8_42z()------- 818,573,293 
optimal4_43z()------- 817,554,978 
CatmullInterpolate()- 817,184,839 *
optimal2_42z()------- 817,017,749 
hermite3()----------- 815,787,715 
optimal8_43z()------- 814,932,498 
optimal_32z()-------- 814,597,788 
optimal2_43z()------- 814,205,750 
optimal16_43z()------ 813,596,244 
optimal32_43z()------ 813,219,182 
watte_42x()---------- 808,266,048 
optimal2_44z()------- 808,242,622 
HermiteInterpolate()- 803,204,436 * - C CLASS
lagrange_43z()------- 794,197,574 
optimal16_23z()------ 793,989,204 
hermite2()----------- 793,803,450 
optimal4_42z()------- 787,663,466 
optimal32_44z()------ 787,217,022 
bspline_43x()-------- 785,941,426 
bspline_43z()-------- 784,983,803 
optimal4_44z()------- 780,245,345 
optimal2_23z()------- 779,379,641 
optimal16_44z()------ 778,057,097 
optimal8_44z()------- 777,514,336 
hermite_43x()-------- 776,101,332 
optimal8_64z()------- 773,369,284 
optimal16_64z()------ 772,615,012 
optimal4_23z()------- 772,374,874 
optimal2_64z()------- 769,535,479 
hermite_43z()-------- 768,493,305 
optimal8_23z()------- 766,868,773 
optimal4_64z()------- 762,666,104 
optimal32_64z()------ 759,502,144 
osculating_45z()----- 754,361,072 
optimal16_65z()------ 746,048,605 
optimal32_65z()------ 743,849,077 *
optimal2_65z()------- 741,950,936 
optimal8_65z()------- 740,090,345 
order3Spline()------- 735,844,776 
optimal4_65z()------- 732,867,435 
hermite_63x()-------- 714,680,545 
hermite_63z()-------- 709,276,578 
osculating_45x()----- 706,546,540   - D CLASS
lagrange_65x()------- 699,041,053 
bspline_65x()-------- 688,782,582 
lagrange_65z()------- 678,868,756 
osculating_65z()----- 671,200,187 
hermite_65z()-------- 669,898,731 
hermite_65x()-------- 646,029,454 
osculating_65x()----- 626,678,371 *
CosineInterpolate()-- 555,128,841 *

An ODS file of the benchmark with comparison is available here.

I added class segmentation because the order tends to variate per benchmark but always tending to stay within a specific area of class rating, so for example, any interpolators in C class are likely to perform anywhere between 700,000,000 and 800,000,000 executions per 32-seconds. The class starts from the line below the class tag, so for example, osculating_45x() is part of the C class. Any interpolators with an asterisk * following them are ones I felt that held their position strongly within that class.

It goes without saying that LinearInterpolate() was always going to come out on top so to speak, but generally, many of the interpolators tended to perform around the 800,000,000 executions range which is only ~26% less than the linear interpolation function. I think that is brilliant and almost always worth the trade-off where higher quality interpolation is required! The optimal set of interpolators reached 732,867,435 executions at their lowest but I still think this is great performance, these are also my personal favourite when it comes down to a question of quality. The worst performing had to be the 65x and 65z variations which are short for 6-point, 5th-order coming in at the ~650,000,000 executions range, with the optimal variations achieving the best performance out of the 65z set coming in at 740,000,000 executions.

The code for the benchmark is available over at GitHub here.