Как использовать вывод обратного преобразования Фурье?

Question

Я пытаюсь усилить ультразвуковой сигнал путем спектрального вычитания. Сигнал находится во временной области и содержит шум. Я разделил сигнал на окна Хэмминга длительностью 2 мкс и рассчитал преобразования Фурье этих кадров. Затем я выбрал 3 последовательных кадра, которые я интерпретировал как шум. Я усреднил спектры величин этих трех кадров и вычел это среднее из спектра величин каждого отдельного кадра. Затем я определил все спектры отрицательной величины как ноль и восстановил усиленное преобразование Фурье, объединив новые спектры амплитуды с фазовыми. Это дает мне серию комплексных чисел на кадр. Теперь я хотел бы преобразовать этот ряд обратно во временную область, используя обратное преобразование Фурье. Однако эта операция предоставляет мне комплексные числа, которые я не знаю, как использовать.

Я прочитал в нескольких постах, что это нормально, чтобы получить сложный результат обратного преобразования Фурье. Однако использование комплексных чисел разделено. Некоторые говорят, что пренебрегать мнимой частью, потому что она должна быть очень маленькой (15), но в моем случае она не пренебрежимо мала (0,01-0,5). Честно говоря, я просто не знаю, что теперь делать с числами, потому что ожидал, что обратное преобразование Фурье даст мне только действительное число. И надеюсь на очень маленькие мнимые части, но, к сожалению ..

# General parameters
#
total_samples = length(time_or) # Total numbers of samples in the current series
max_time = max(time_or) # Length of the measurement in microseconds
sampling_freq = 1/(max_time/1000000)*total_samples # Sampling frequency
frame_length_t = 2 # In microseconds (time)
frame_length_s = round(frame_length_t/1000000*sampling_freq) # In samples per frame
overlap = frame_length_s/2 # Overlap in number of frames, set to 50% overlap
#
# Transform the frame to frequency domain
#
fft_frames = specgram(amp, n=frame_length_s, Fs=125, window=hamming(frame_length_s), overlap=overlap)

mag_spec=abs(fft_frames[["S"]])
phase_spec=atan(Im(fft_frames[["S"]])/Re(fft_frames[["S"]]))

#
# Determine the arrival time of noise
#
cutoff= 10 #determine the percentage of the signal that has to be cut off
dnr=us_data[(length(us_data[,1])*(cutoff/100)):length(us_data[,1]), ]
noise_arr=(length(us_data[,1])-length(dnr[,1])+min(which(dnr[,2]>0.01)))*0.008
#
# Select the frames for noise spectrum estimation
#
noise_spec=0
noise_spec=mag_spec[,noise_arr]
noise_spec=noise_spec+mag_spec[, (noise_arr+1)]
noise_spec=noise_spec+mag_spec[, (noise_arr+2)]
noise_spec_check=noise_spec/3
#
# Subtract the estimated noise spectrum from every frame
#
est_mag_spec=mag_spec-noise_spec_check
est_mag_spec[est_mag_spec < 0] = 0
#
# Transform back to frequency spectrum
#
j=complex(real=0, imaginary=1)
enh_spec = est_mag_spec*exp(j*phase_spec)
#
# Transform back to time domain
#
install.packages("pracma")
library("pracma")
enh_time=fft(enh_spec[,2], inverse=TRUE)

Я надеюсь, что есть кто-нибудь с идеей о том, как обрабатывать эти сложные числа. Возможно, я ранее допустил ошибку в методе обработки, но я проверял ее несколько раз, и она кажется мне довольно серьезной. Это (последний, но последний) последний шаг процесса, и он действительно надеется получить хороший сигнал во временной области после обратного преобразования Фурье.

Answer 1

Это ваша проблема:

phase_spec=atan(Im(fft_frames[["S"]])/Re(fft_frames[["S"]]))

Здесь вы вычисляете угол в полукруге, сопоставляя вторую половину с первой. То есть вы теряете информацию.

Многие языки имеют функцию для получения фазы комплексного значения, например, в MATLAB это angle, а в Python numpy.angle.

В качестве альтернативы используйте функцию atan2, которая существует на каждом отдельном языке, который я когда-либо использовал, за исключением NumPy, который они решили назвать ее arctan2. Он вычисляет арктангенс четырех квадрантов, принимая два компонента как отдельные значения. То есть atan(y/x) совпадает с atan2(y,x), если результат находится в первых двух квадрантах.

Полагаю, вы можете сделать

phase_spec=atan2(Im(fft_frames[["S"]]), Re(fft_frames[["S"]]))

Answer 2

Важным средством устранения неполадок при преобразовании данных с использованием преобразования Фурье является представление о том, что вы можете сделать fft, затем взять эти данные и сделать обратное fft и вернуть ваши исходные данные ... Я предлагаю вам освоиться с игрушечным вводом данные во временной области ... допустим, звуковая волна 1 кГц, которая является вашими данными во временной области ... отправить ее в вызов fft, который вернет обратно массив его представления в частотной области ... ничего не делая с этими данными, отправьте их в обратный fft (ifft) ... возвращенные данные будут вашей исходной звуковой волной 1 кГц ... сделайте это сейчас, чтобы оценить ее мощность и используйте этот трюк в своем проекте, чтобы подтвердить, что вы находитесь в парке событий. .. В качестве альтернативы, если вы начинаете с данных в частотной области, вы также можете сделать это ...

freq domain data -> ifft -> time domain data -> fft -> same freq domain data

или

time domain -> fft -> freq domain -> ifft -> same time domain data

подробности см. Здесь Получить частоту с наибольшей амплитудой из FFT