Как реализовать алгоритм 2D RFFT на основе 1D RFFT?

Question

Я пытаюсь реализовать NumPy s rfft2(), функцию RFFT , которая поддерживает массивы с 2 измерениями, выполняя 1D RFFT для каждой строки, а затем повторное выполнение 1D RFFT для каждого столбца предыдущего результата.

Этот подход хорошо работает для реализации функции 2D FFT , как обсуждалось ранее в этом посте , но, похоже, он не работает для 2D RFFT .

Вот сценарий, который реализует пользовательскую функцию 2D FFT, которая следует этой идее, используя 1D версию FFT NumPy в качестве основы, а затем сравнивает свой результат с фактической 2D версией из NumPy:

import cmath
import numpy as np
import math

def my_fft2d(matrix):
    fft_rows = [np.fft.fft(row) for row in matrix]
    return np.transpose([np.fft.fft(row) for row in np.transpose(fft_rows)])


# initialize test data
img = np.array([[0,0,0,0], [0,1,0,0], [0,0,0,0], [0,0,0,0]])
print('img shape=', img.shape)

# perform custom FFT2D and print result
custom_result = my_fft2d(img)
print('\ncustom_result shape=', custom_result.shape)
for row in custom_result:
   print(', '.join(['%.3f + %.3fi' % (x.real, x.imag) for x in row]))

# perform numpy FFT2D and print result
numpy_result = np.fft.fft2(img)
print('\nnumpy_result shape=', numpy_result.shape)
for row in numpy_result:
   print(', '.join(['%.3f + %.3fi' % (x.real, x.imag) for x in row]))

# compare results
print('\nAre the results equivalent to NumPy?', np.allclose(custom_result, custom_result))
print('ASSERT(assert_array_almost_equal):', np.testing.assert_array_almost_equal(custom_result, custom_result))

Вывод :

img shape= (4, 4)

custom_result shape= (4, 4)
1.000 + 0.000i, 0.000 + -1.000i, -1.000 + 0.000i, 0.000 + 1.000i
0.000 + -1.000i, -1.000 + 0.000i, 0.000 + 1.000i, 1.000 + 0.000i
-1.000 + 0.000i, 0.000 + 1.000i, 1.000 + 0.000i, 0.000 + -1.000i
0.000 + 1.000i, 1.000 + 0.000i, 0.000 + -1.000i, -1.000 + 0.000i

numpy_result shape= (4, 4)
1.000 + 0.000i, 0.000 + -1.000i, -1.000 + 0.000i, 0.000 + 1.000i
0.000 + -1.000i, -1.000 + 0.000i, 0.000 + 1.000i, 1.000 + 0.000i
-1.000 + 0.000i, 0.000 + 1.000i, 1.000 + 0.000i, 0.000 + -1.000i
0.000 + 1.000i, 1.000 + 0.000i, 0.000 + -1.000i, -1.000 + 0.000i

Are the results equivalent to NumPy? True
ASSERT(assert_array_almost_equal): None

Вывод скрипта показывает, что реализация my_fft2d() совместима с np.fft.fft2().

Однако, когда тот же logi c применяется для реализации RFFT-версии преобразования, результирующий массив имеет другую форму, как демонстрирует сценарий ниже:

def my_rfft2d(matrix):
    fft_rows = [np.fft.rfft(row) for row in matrix]
    return np.transpose([np.fft.rfft(row) for row in np.transpose(fft_rows)])


# initialize test data
img = np.array([[0,0,0,0], [0,1,0,0], [0,0,0,0], [0,0,0,0]])
print('img shape=', img.shape)

# perform custom FFT2D and print result
custom_result = my_rfft2d(img)
print('\ncustom_result shape=', custom_result.shape)
for row in custom_result:
   print(', '.join(['%.3f + %.3fi' % (x.real, x.imag) for x in row]))

# perform numpy FFT2D and print results
numpy_result = np.fft.rfft2(img)
print('\nnumpy_result shape=', numpy_result.shape)
for row in numpy_result:
   print(', '.join(['%.3f + %.3fi' % (x.real, x.imag) for x in row]))

Вывод :

img shape= (4, 4)
C:\Users\username\AppData\Roaming\Python\Python37\site-packages\numpy\fft\_pocketfft.py:77: ComplexWarning: Casting complex values to real discards the imaginary part
  r = pfi.execute(a, is_real, is_forward, fct)

custom_result shape= (3, 3)
1.000 + 0.000i, 0.000 + 0.000i, -1.000 + 0.000i
0.000 + -1.000i, 0.000 + 0.000i, 0.000 + 1.000i
-1.000 + 0.000i, 0.000 + 0.000i, 1.000 + 0.000i

numpy_result shape= (4, 3)
1.000 + 0.000i, 0.000 + -1.000i, -1.000 + 0.000i
0.000 + -1.000i, -1.000 + 0.000i, 0.000 + 1.000i
-1.000 + 0.000i, 0.000 + 1.000i, 1.000 + 0.000i
0.000 + 1.000i, 1.000 + 0.000i, 0.000 + -1.000i

Как видите, в выводе есть две проблемы:

• предупреждение от numpy жалуется на что-то Я не совсем уверен, как исправить;
• пользовательская реализация 2D RFFT возвращает результат с меньшим количеством строк, чем тот, который возвращается np.fft.rfft2();

Как я могу исправьте эту проблему и сделайте my_rfft2d() совместимым с np.fft.rfft2()?

Answer 1

Как сказал комментатор. Тебе нужно сдать второй раз. Это связано с тем, что вывод rfft строк сложен. Это решает комплексную и реальную ошибку, а также проблему формы.

import numpy as np

def my_rfft2d(matrix):
    fft_rows = [np.fft.rfft(row) for row in matrix]
    return np.transpose([np.fft.fft(row) for row in np.transpose(fft_rows)])


# initialize test data
img = np.array([[0,0,0,0], [0,1,0,0], [0,0,0,0], [0,0,0,0]])
print('img shape=', img.shape)

# perform custom FFT2D and print result
custom_result = my_rfft2d(img)
print('\ncustom_result shape=', custom_result.shape)
for row in custom_result:
   print(', '.join(['%.3f + %.3fi' % (x.real, x.imag) for x in row]))

# perform numpy FFT2D and print results
numpy_result = np.fft.rfft2(img)
print('\nnumpy_result shape=', numpy_result.shape)
for row in numpy_result:
   print(', '.join(['%.3f + %.3fi' % (x.real, x.imag) for x in row]))

Вывод:

custom_result shape= (4, 3)
1.000 + 0.000i, 0.000 + -1.000i, -1.000 + 0.000i
0.000 + -1.000i, -1.000 + 0.000i, 0.000 + 1.000i
-1.000 + 0.000i, 0.000 + 1.000i, 1.000 + 0.000i
0.000 + 1.000i, 1.000 + 0.000i, 0.000 + -1.000i

numpy_result shape= (4, 3)
1.000 + 0.000i, 0.000 + -1.000i, -1.000 + 0.000i
0.000 + -1.000i, -1.000 + 0.000i, 0.000 + 1.000i
-1.000 + 0.000i, 0.000 + 1.000i, 1.000 + 0.000i
0.000 + 1.000i, 1.000 + 0.000i, 0.000 + -1.000i

Answer 2

Как я сказал в своем комментарии, после взятия rfft строк вы должны взять fft вместо rfft, потому что результат rfft в целом сложен.

Я не Не знаю, почему вы пытаетесь go реально, но если вы действительно хотите go все реально, вы должны использовать DCT (дискретное косинусное преобразование) вместо FFT, потому что выход DCT является реальным. Вы можете использовать тот же подход, что и при вычислении 2D БПФ выше, потому что вы можете разложить 2D DCT аналогичным образом.