# How can I interpolate between two FFTs?

I have two FFTs that represent two impulse responses. I want to interpolate between the two of them to generate a 3rd FFT in the manner shown in the diagram below. This is to be used in real time audio rendering filtering between two FFT filters.

Linear interpolation between each of the FFT's complex coefficients doesn't work as wanted.

FFT plots of two impulse responses
(x = frequency/coefficient, y = frequency magnitude)

+                                   +
|                                   |
|        +                          |
|        |                          |
|        |                          |
|        |                          |
|        |                          |
|        |                          |                   +
|        |                          |                   |
|        |                          |                   |
+--------+------------------+       +-------------------+-------+
Input A                              Input B

FFT plots of interpolation between A and B (Linear, and Wanted)
(x = frequency/coefficient, y = frequency magnitude)

+                                   +
|                                   |
|                                   |
|                                   |
|                                   |
|                                   |             +
|        +                          |             |
|        |                          |             |
|        |          +               |             |
|        |          |               |             |
+--------+----------+-------+       +-------------+-------------+
Linear Interpolation (A+B/2)            Wanted Interpolation


How do I do the interpolation as described in the Wanted Interpolation graph?

• Welcome to dsp.se! Your question is definitely interesting, but I think a bit more context is needed. Pleased edit your question to include the following: What is your application? What do you mean by interpolate between two FFT? What do your figures show (In particular, please add axis labels)? What is input to the algorithm you look for and what should the output then look like? (I think this is what you want to show with your ASCII-art) – Maximilian Matthé Feb 15 '17 at 20:45
• Done. Hopefully this has clear up your questions. – JodiTheTigger Feb 15 '17 at 21:24
• You want to sweep the frequency peaks to new frequencies over time? You'll need to do some kind of peak picking first to identify what counts as a "peak", then synthesize whatever sweep shape you want. I imagine it will sound weird. The A and B spectra always have the same number of frequencies? And it's obvious how one maps to the other? Or are you saying that you want to "expand" all frequencies simultaneously? Like a spectral stretching? – endolith Feb 15 '17 at 21:33
• Yes, I suppose I do, I want to "sweep" the peak in a at frequency fa, to the peak at b at frequency fb. Or interpolate between these two peaks, where phase, magnitude and frequency are interpolated. However this is just a simplified explination. My actual FFT has multiple peaks, and it might not be obvious which peak in A maps to which peak in B (but that's starting to get out of scope of the question). – JodiTheTigger Feb 15 '17 at 21:50
• Well, I think the matching between the peaks of A and B is the main issue here. It could be done by associating the peaks based on their sorted frequencies. However, first the peaks need to be detected. What should happen, if one FFT has more peaks than others? Once you have associations between input frequencies and output frequencies, you can create a mapping that performs a e.g. piece-wise linear mapping between input and output frequencies. – Maximilian Matthé Feb 15 '17 at 21:59

This is the solution for a linear interpolation between two frequencies as a "weighted average". You would just need to decide what your weighting actually is in terms of magnitude or dB scale for the vertical axis. X is the interpolated frequency, weighted toward the stronger signal, and Y is the linear interpolation at the new frequency location of the signal strenths a and b.

Here is a simple example demonstrating the results that would be achieved with the formulas given:

Have a look at the script below. It does not solve all your problems, but it can be a starting point:

import scipy.interpolate

# 1) Here, I define some arbitrary initial spectrum,
#    for example it consts of two paeaks
#    Here, I simplify and set the frequencies to be between 0 and 1 in 100 steps.
f = np.linspace(0,1,100)
#H_begin = lambda f: 1.5+np.sin(2*np.pi*f)
H_begin = lambda f: (abs(f-0.05) < 0.01) + (abs(f-0.5)<0.005)

# define the matching frequencies:
#  0 at the beginning remains zero after the warping
#  1 at the beginning remains one after warping
#  frequency 0.5 should be warped to 0.25
f_matches = [(0,0), (1,1), (0.5,0.25)]

# generate the according inputs and outputs of the 2d function
# that should be interpolated
f_in = np.array([x[0] for x in f_matches]*2)
t_in = np.array([0]*len(f_matches)+[1]*len(f_matches))
f_out = np.array([x[0] for x in f_matches] + [x[1] for x in f_matches])

# just printout the supporting points we have defined
for n in range(len(f_in)):
print ("X(f=%f,t=%f)=%f" % (f_in[n], t_in[n], f_out[n]))

# generate the interpolation
X = scipy.interpolate.interp2d(f_out, t_in, f_in)

# visualize the interpolation
t_range = np.linspace(0,1,100)
img = np.zeros((len(f), len(t_range)))
for n, t0 in enumerate(t_range):
img[:,n] = H_begin(X(f,t0))
plt.imshow(img)
plt.xlabel('interpolation time')
plt.ylabel('frequency');


Output:

X(f=0.000000,t=0.000000)=0.000000
X(f=1.000000,t=0.000000)=1.000000
X(f=0.500000,t=0.000000)=0.500000
X(f=0.000000,t=1.000000)=0.000000
X(f=1.000000,t=1.000000)=1.000000
X(f=0.500000,t=1.000000)=0.250000


It does the following:

Essentially, you have a 2D-interpolation problem: You have pairs of frequencies that should be warped from the first element in the pair to the second element, when the interpolation factor $t$ increases from 0 to 1. (I.e. for $t=0$ you get the initial spectrum, for $t=1$ you get the final spectrum). So, given a mapping between some supporting start-frequencies and end-frequencies, we can create a 2-dimensional interpolation. For an input frequency and interpolation time $t$, it returns to where this frequency would be be warped to.

What is not solved:

• How to find the supporting points f_match. You might look for peak-matching algorithms or come up with something that is adapted for your application
• I think the interpolation should be done in log-domain rather than in linear domain if this is about audio processing. However, the algorithm wont change, just take the log of every frequency before interpolation and convert it back to linear frequency after interpolation.
• my script assumes you have accurate information on the spectrum at all frequencies (H_begin is a function). However, in reality you only have access to discrete frequencies of the spectrum. So, if the interpolation requires you to get the value from a frequency between two bins, you'd need to interpolate the original spectrum to approximate this value.