Benchmarking interpolation functions.

James William Fletcher
4 min readSep 3, 2021

--

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.

--

--