Finding maximum using DFT

I'm trying to find an efficient way to compute maximum of a signal using its DFT. More formally:

$$\max\left\{ \mathcal F^{-1}\left(X_k\right)\right\}, X_k\text{ is the DFT of the signal and } \mathcal F^{-1} \text{ is the IDFT}$$

The original signal $x(n)$ is real, and it has some noise when sometimes there are wide peaks.

I'm looking for a solution that is quicker than the actual IFFT since the signal is very long (but we have its DFT).

• So, wait, are you looking for its maximum in time or in frequency domain? Can you please define what $F$ and $X$ are? Mar 27, 2017 at 12:15
• so, am I understanding this correctly, you only have the DFT of your signal, but not the original signal? If so, does the signal or its DFT have any special properties (real-valuedness, symmetry, that kind of thing)? what constitutes a maximum? Maximum magnitude? Maximum Real Part? Mar 27, 2017 at 12:24
• By the way, what is "very long"? Maybe you'd want to also tell us a bit about your computing platform, because things will be different for implementation on ASIC/FPGA/GPU/CPU with SIMD/MCU… Also, it tends to be a good idea to explain where your constraint comes from (is it latency, memory or computational throughput, power consumption?) Mar 27, 2017 at 12:38
• Would you be satisfied with an estimate of the rms level or do you really need to be concerned with the value of an occasional "spike"? Mar 27, 2017 at 23:33
• Mar 29, 2017 at 11:17

This sounds like a good candidate for the Sparse Fast Fourier Transform (sFFT). Have a look at H. Hassanieh, P. Indyk, D. Katabi, and E. Price, sFFT: Sparse Fast Fourier Transform, 2012. I don't know details of the algorithm. You have the domains swapped, which should be fine as the discrete Fourier transform (DFT) and its inverse (IDFT) are nearly identical.

sFFT has been discussed also here on dsp.stackexchange.com.

• Actually I think that as is it won't help much, but there is potential for something in the theory that could help me, thanks! Mar 27, 2017 at 12:48
• I'm sorry I ever doubted you, this actually looks very fitting! thanks a lot! Mar 27, 2017 at 13:16

There is no way to compute the (exact) maximum of the time domain signal directly from its DFT (without doing an IDFT first). Without any further knowledge, the best you can do is

$$|x[n]|=\left|\frac{1}{N}\sum_{k=0}^{N-1}X[k]e^{j2\pi nk/N}\right|\le\frac{1}{N}\sum_{k=0}^{N-1}|X[k]|$$

which gives you an upper bound. However, this upper bound is usually not very tight, especially for large values of the DFT length $N$.

• Yes, this becomes an equality for any $n$ for which $X[k]$ has for all $k$ the same phase as the complex conjugate of $e^{j2\pi n k/N}.$ Apr 18, 2017 at 13:26

Generally:

There can be no shortcut here – every single element of the (inverse or forward, doesn't matter) DFT needs every element of the input. So, the FFT is already pretty much as fast as you can go.

However, you say:

The original signal $x(n)$ is real, and it has some noise when sometimes there are wide peaks.

Aha! That means that your $X$ is (hermitian) symmetric. Which means you might be able to pick an FFT algorithm that benefits from that knowledge and can omit quite a few calculation – but I don't think that'll change the general $\mathcal O(N \log N)$ complexity of the IFFT, nor the $\mathcal O(N)$ of finding the maximum.

• I actually thought the more important part is the wide peaks part, because I thought of using something like the sparse cognitive radio, but yea this helps too, thanks! Mar 27, 2017 at 12:47
• well, sparsity is something you can exploit when you know a lot about the number of base vectors you'd need to have to represent your signal. You didn't say anything about the number of your "wide peaks" (which by the way is pretty non-descriptive) nor about how long your vectors are. Mar 27, 2017 at 13:01

The maximum of $$x[n]$$, where $$X[k] = \texttt{DFT}\{x\}$$, is

\begin{align} & \texttt{max}\{x\} = \frac{1}{N}\underset{s\in[0, N-1]}{\texttt{max}}\left\{ \sum_{l=0}^{N-1} (X[k] \cdot e^{-j2\pi s k/N})[l] \right\} \tag{1} \\ & k = [0, 1, ..., N - 1] \end{align}

or in code-like syntax,

$$\texttt{max}\{x\} = \frac{1}{N}\texttt{max}\big\{ \big[\texttt{sum}\{\texttt{FFT}\{x\} \cdot e^{-j2\pi s [0:N-1] / N}\} \ \text{for}\ s=[0:N-1]\big]\big\} \tag{2}$$

The DC bin $$= \texttt{sum}\{x\}$$. Per duality, $$x[0]$$ is the DC bin for the DFT: $$x[0] = \texttt{sum}\{X\} / N$$, for any $$x$$, where $$1/N$$ is the only difference between the domains. Hence, if we shift $$\texttt{max}\{x\}$$ to index 0, it's retrievable by summing $$X$$. We don't know where max is ahead of time, so we try all possible shifts $$s$$, and by definition the max will be the max of all sums (sum at shift $$-s$$ equals $$x[s]$$).

$$\texttt{min}\{x\}$$ is supported, and so are real and complex-valued inputs. In fact any $$\texttt{func}\{x\}$$ is supported, since we're quite literally retrieving $$x$$.

This is a very slow algorithm, and I don't know if it is reducible to a satisfactory speed. If we seek to fetch the maximum by working purely in frequency, or understanding the relationship, then it's one way. I found it useful for proving that "time aliasing" cannot cause peaks in convolution.

Code validation

import numpy as np

for func in (max, min, lambda x: max(abs(x))):
# N==1 works but requires more pure coding related stuff
for N in range(2, 129):
x = np.random.randn(N) + 1j*np.random.randn(N)
xf = fft(x)

exparg = -1j*2*np.pi*np.arange(N)/N
shifts = np.array([sum(xf * np.exp(exparg * s)) for s in range(N)])
out = func(shifts) / N
assert np.allclose(out, func(x))