Numpy.fft.fft2() gives different result than tf.fft2d()

[Sushobhan04]

as mentioned in the issue #6401, the tf.fft2d() gives different result compared to np.fft.fft2(). Is there a reason for this ?

Note : numpy gives proper fourier transform after np.fft.fftshift(), and I have taken care of that in my code.
The differences are not visible here, but the mean squared error is significant.

Second image is the fft using tensorflow, and third one is using numpy. You can see the difference in the corners. The fourth image is the difference between the two images times 10.

(tf has some features while numpy does not)

I am working on an application which uses fft in backpropagation and thus it is of absolute importance that the fft in numpy are same as fft by tf.

My question is - Why is there a difference and how can I get the same fft as numpy ?

[gongbudaizhe]

I am not sure if this is related to this issue #10735

[aselle]

Can you give any indication about the relative error between the two versions? How much did you magnify the difference image?

[Sushobhan04]

@aselle Hey, sorry for late response. I magnified the difference by a factor of 10.

I realized that numpy uses float64 whereas tf uses float32. I could not find a way to use float64 in tf (creates a complex128, which is not supported by tf.fft2d). Could this be a possible problem ?

The difference kind of seems huge to consider floating point accuracy as the reason.

[aselle]

We use cuFFT to compute fft's, so it is a property of the cufft. You could look at its documentation to see if it has any insight (or google for cufft vs fftw). Also, you could see whether numpy.fft when run with floats matches better.

[absudabsu]

This is still an issue. I am not sure how much it will matter in the long run if we are training a network around it anyway..., but it is quite shocking that the error is so high. This also makes debugging hard, because TF's FFT is undoubtedly not doing what we think it should be. Perhaps cuFFT's algorithm is not numerically stable? The error is there even for 1D, but it is much much smaller. This makes me think that the 2D implementation is not optimal, with a lot of fill in or non-optimal number of operations. If someone can point me to the code, maybe we can debug / suggest a fix without changing too much.

Some python code to verify error using complex64 for both numpy and tf (i.e. based on float32):

## DEFINE FFTs
fft1D_np = np.fft.fft
fft1D_tf = tf.fft
ifft1D_np = np.fft.ifft
ifft1D_tf = tf.ifft
fft2D_np = np.fft.fft2
fft2D_tf = tf.fft2d
ifft2D_np = np.fft.ifft2
ifft2D_tf = tf.ifft2d

  ## TEST 2D FFT
   np.random.seed(13)
   M = np.random.random([1,2,2]).astype(np.complex64)
   M_tf = tf.placeholder(dtype=M.dtype, shape=[None,*M.shape[1:]], name='M_tf')
   eval_dict = {M_tf: M}
   sess, gsaver = nn_utils.init_network()
   Mhat = fft2D_np(M) 
   Mhat_tf = fft2D_tf(M_tf)
   Mhat_tfeval = Mhat_tf.eval(feed_dict=eval_dict)
   e = Mhat - Mhat_tfeval
   err = np.linalg.norm(e.reshape([-1,1]))
   assert( np.isclose(err, 0.0) )

Doing some prints in my session:

>>> M
array([[[ 0.77770239+0.j,  0.23754121+0.j],
        [ 0.82427853+0.j,  0.96574920+0.j]]], dtype=complex64)
>>> Mhat
array([[[ 2.80527139+0.j,  0.39869052+0.j],
        [-0.77478409+0.j,  0.68163186+0.j]]], dtype=complex64)
>>> Mhat_tfeval
array([[[ 2.80527139+0.j,  0.39869046+0.j],
        [-0.77478415+0.j,  0.68163186+0.j]]], dtype=complex64)
>>> e
array([[[  0.00000000e+00+0.j,   5.96046448e-08+0.j],
        [  5.96046448e-08+0.j,   0.00000000e+00+0.j]]], dtype=complex64)
>>> err
8.4293696e-08
>>> np.isclose(err, 0.0)
False

[rryan]

@sirgogo, can you please post the source to fft2D_np and fft2D_tf in your example?

Are your ops running on GPU or CPU?
(have you verified that with tf.ConfigProto(log_device_placement=True)?)

FWIW, NumPy is the oracle used for all of the FFT unit tests. The tolerance on those checks is 1e-4 at the moment, which is pretty high, but it's because we're using a one-size-fits all bound against Eigen's TensorFFT (CPU), and cuFFT (GPU). Please also note that the FFT unit tests include gradient tests, that verify the numerical gradient matches the symbolic gradient.

Due to differences in the floating point hardware across your CPU and GPU, the results between NumPy and cuFFT will differ by some amount for an identical sequence of floating point operations. Additionally, NumPy uses fftpack for its FFTs -- so whether your ops are running on CPU or GPU, it's a different FFT implementation from TensorFlow's. That's not to say there isn't a bug in TensorFlow's invocation of these FFT implementations or the implementations themselves!

Here are code pointers:

• Eigen TensorFFT
• TensorFlow FFT kernels. Used to invoke either Eigen TensorFFT or cuFFT.
• cuFFT documentation
• StreamExecutor invocation of cuFFT
• TensorFlow FFT gradient definitions

What is your expectation here? e.g. would you be happy if np.isclose returned True in your example? That would imply you need something like rtol=1e-5, atol=1e-8 ?

[absudabsu]

@rryan, thank you for the info and links! I will take a look shortly. Sorry for the delay, I was out of town, server went down, etc etc.

You can run the example with the following simple definitions: I have updated that post with this

fft1D_np = np.fft.fft
fft1D_tf = tf.fft
ifft1D_np = np.fft.ifft
ifft1D_tf = tf.ifft
fft2D_np = np.fft.fft2
fft2D_tf = tf.fft2d
ifft2D_np = np.fft.ifft2
ifft2D_tf = tf.ifft2d

I am running with tensorflow-GPU. This is the only way I was able to get any of the tf.fft-like functions to work. It complains if you are using the CPU-only version (at least for 2D). Is this sufficient to assume the ops are running on the GPU? I can try the log_device_placement setting.

Yes, I would expect that the maximum relative error needs to be lower. That would probably make me / most folks happy. Thanks for your insight into this. Just to clarify, you are saying this is an issue with Eigen's TensorFFT / cuFFT, and not something that would be fixed by compiling from source?

[Rassibassi]

Hi,
using tensorflow to do some signal processing and the fft problem has appeared. I did some tests and hope they are helpful. I also created a repository with the example signal where the fft really goes wrong using the cpu, error after fft is >1e+14. With a gpu it looks good.

Code:
(https://github.com/rassibassi/fftTensorflow)

import tensorflow as tf
import scipy.io
import numpy as np

print('tf:',tf.__version__)
print('scipy:',scipy.__version__)
print('np:',np.__version__)

# numpy fft
signal = scipy.io.loadmat('signal.mat')['signal']
signal_f = np.fft.fft(signal)

# tf fft
tf_signal = tf.placeholder(shape=signal.shape,dtype=tf.complex64)
tf_signal_f = tf.fft(tf_signal)
session = tf.InteractiveSession()
out_signal_f, out_signal = session.run([tf_signal_f,tf_signal], feed_dict={tf_signal:signal})

errorBefore = np.mean(np.power(np.abs(out_signal-signal),2))
errorAfter = np.mean(np.power(np.abs(out_signal_f-signal_f),2))
print('Error before fft:',errorBefore)
print('Error after fft:',errorAfter)

unix ubuntu cpu:

tf: 1.4.0
scipy: 0.19.1
np: 1.13.1
Error before fft: 6.43755807539e-22
Error after fft: 4.21685308781e+14

scientific linux hpc-cluster cpu:

tf: 1.3.0
scipy: 0.17.1
np: 1.11.0
Error before fft: 6.43755807539e-22
Error after fft: 4.21685308781e+14

scientific linux hpc-cluster gpu (Tesla K80) with cudnn/v5.1:

('tf:', '1.0.1')
('scipy:', '0.18.1')
('np:', '1.11.2')
('Error before fft:', 6.4375580753914162e-22)
('Error after fft:', 2.8628744462657876e-15)

[rryan]

Thanks very much for the helpful repro @Rassibassi! I can also reproduce the issue.

This suggests an issue with Eigen TensorFFT or the TensorFlow plumbing that uses it. Interestingly, if I trim your 66528 tensor to length 65536 (power of 2) the results match numpy to 1e-8, but I get high error at length 65535. I get similar results when I trim the input tensor to length 32768/32767, but not when I trim to 16384/16383. Maybe suggests something is overflowing? We'll take a look.

cc: @rmlarsen

('tf:', '1.4.1')
('scipy:', '1.0.0')
('np:', '1.13.3')
('signal=', (1, 65536))
('Error before fft:', 6.4749546908751922e-22)
('Error after fft:', 1.0376122868548086e-08)

('tf:', '1.4.1')
('scipy:', '1.0.0')
('np:', '1.13.3')
('signal=', (1, 65535))
('Error before fft:', 6.4750239427498006e-22)
('Error after fft:', 55656499607300.164)

('tf:', '1.4.1')
('scipy:', '1.0.0')
('np:', '1.13.3')
('signal=', (1, 32768))
('Error before fft:', 6.1796690528036584e-22)
('Error after fft:', 7.996054658436124e-10)

('tf:', '1.4.1')
('scipy:', '1.0.0')
('np:', '1.13.3')
('signal=', (1, 32767))
('Error before fft:', 6.1793787363606864e-22)
('Error after fft:', 1754207811720195.0)

('tf:', '1.4.1')
('scipy:', '1.0.0')
('np:', '1.13.3')
('signal=', (1, 16384))
('Error before fft:', 6.0547451539629847e-22)
('Error after fft:', 1.1566757798862581e-10)

('tf:', '1.4.1')
('scipy:', '1.0.0')
('np:', '1.13.3')
('signal=', (1, 16383))
('Error before fft:', 6.055074302285619e-22)
('Error after fft:', 8.9115448530747228e-06)

[edumotya]

Hi,
I've gotten high errors also running one dimensional tf.rfft() on CPU over large signals (>280000 samples). On the other hand, after trimming the tensor length to the next power of two the results are quite better (in terms of error), as pointed out by @rryan .

[rryan]

We think we tracked the issue (or part of the issue) down to this code in Eigen that computes twiddle factors for non-power-of-2 FFTs via Bluestein's algorithm.

https://bitbucket.org/eigen/eigen/src/034b6c3e101792a3cc3ccabd9bfaddcabe85bb58/unsupported/Eigen/CXX11/src/Tensor/TensorFFT.h?at=default&fileviewer=file-view-default#TensorFFT.h-230:249

We'll submit a fix to Eigen, then update the version of Eigen that TensorFlow uses.