Как извлечь первый компонент БПФ из 4D Image - PullRequest
1 голос
/ 05 марта 2020
  • У меня 4D (2D + срезы по оси z + временные рамки) серого изображения для сердца, бьющегося в разные моменты.

  • Мне нравится брать Преобразование Фурье по оси времени (для каждого среза отдельно) и анализировать фундаментальную гармонику c (также называемую компонентой H1, где H обозначает Гильбертово пространство), чтобы я мог определяет области пикселей, соответствующие ROI, которые показывают наиболее сильный отклик на кардию c частота .

  • Я использую python для этой цели, и я пытался сделайте это с помощью следующего кода, но я не уверен, что это правильный способ сделать это, потому что я не знаю, как определить частоту среза, чтобы сохранить только основную гармонию c.

Эта ссылка на изображение, с которым я имею дело

import nibabel as nib
import numpy as np
import matplotlib.pyplot as plt

img = nib.load('patient057_4d.nii.gz')
f = np.fft.fft2(img)
#  Move the DC component of the FFT output to the center of the spectrum
fshift = np.fft.fftshift(f)
fshift_orig = fshift.copy()
# logarithmic transformation
 magnitude_spectrum = 20*np.log(np.abs(fshift))
# Create mask
rows, cols = img.shape
crow, ccol = int(rows/2), int(cols/2)
# Use mask to remove low frequency components
dist1 = 20
dist2 = 10
fshift[crow-dist1:crow+dist1, ccol-dist1:ccol+dist1] = 0
#fshift[crow-dist2:crow+dist2, ccol-dist2:ccol+dist2] = fshift_orig[crow-dist2:crow+dist2, ccol-dist2:ccol+dist2] 

# logarithmic transformation
magnitude_spectrum1 = 20*np.log(np.abs(fshift)) 
f_ishift = np.fft.ifftshift(fshift)
# inverse Fourier transform
img_back = np.fft.ifft2(f_ishift)
# get rid of imaginary part by abs
img_back = np.abs(img_back)
plt.figure(num = 'Im_Back')
plt.imshow(abs(fshift[:,:,2,2]).astype('uint8'),cmap='gray')
plt.show()

1 Ответ

1 голос
/ 07 марта 2020
  • Решение состояло в том, чтобы взять трехмерное преобразование Фурье для каждого среза отдельно, а затем выбрать только 2-й компонент преобразования, чтобы преобразовать его обратно в пространственное пространство, и все.
  • Преимущество это обнаружение того, что что-то движется вдоль третьей оси (время в моем случае).
for sl in range(img.shape[2]):
   #-----Fourier--H1-----------------------------------------
   # ff1[:, :, 1] H1 compnent 1, if 0 then DC
   ff1 = FFT.fftn(img[:,:,sl,:])
   fh = np.absolute(FFT.ifftn(ff1[:, :, 1])) 

   #-----Fourier--H1-----------------------------------------
...