Амплитуда и фазовый спектр. Сдвиг фазы без изменения амплитуды

У меня есть данные с эквивалентными интервалами и соответствующими измерениями в соответствующих точках. В качестве примера приведу отрывок из имеющихся у меня данных:

y =[2.118, 2.1289, 2.1374, 2.1458, 2.1542, 2.1615, 2.1627, 2.165 2.1687...]

интервал между точками 0,1

Итак, что мне нужно получить из данных, это спектр амплитуды (амплитуда против частоты), а также фазовый спектр (угол фазы против частоты). Вдобавок я должен сдвинуть фазу данных на отрицательные 90 градусов (-pi / 2).

После сдвига фазы и оставления амплитуды нетронутой мне нужно сделать обратный fft и получить новый сигнал. Я хочу сделать это на Python.

Не могли бы вы привести мне пример того, как это сделать.

Код, который я использовал, был взят из другого вопроса SO, но я сделал некоторые изменения

## Perform FFT WITH SCIPY
signalFFT = np.fft.fft(y)

## Get Power Spectral Density
signalPSD = np.abs(signalFFT) ** 2
signalPhase = np.angle(signalFFT)

## Shift the phase py +90 degrees
new_signalPhase =(180/np.pi)*np.angle(signalFFT)+90 


## Get frequencies corresponding to signal 
fftFreq = np.fft.fftfreq(len(signalPSD), 0.1)

## Get positive half of frequencies
i = fftFreq>0

##
plt.figurefigsize=(8,4)
#plt.plot(fftFreq[i], 10*np.log10(signalPSD[i]));

plt.plot(fftFreq[i], new_signalPhase[i]);
plt.ylim(-200, 200);
plt.xlabel('Frequency Hz');
plt.ylabel('Phase Angle')
plt.grid()
plt.show()

Проблема в том, что я хочу восстановить сигнал с теми же амплитудами, но со сдвигом фазы. Я знаю, что ответ связан с чем-то, связанным с ошибкой, но как мне подготовить для этого данные? Не могли бы вы посоветовать мне дальнейшие шаги.

выход

Привет, Кямран, и добро пожаловать в ТАК! Задавая здесь вопрос, важно показать нам, что вы уже пытались решить для своих проблем, и очень конкретно указать, с чем вы застряли. Прочтите как мне задать хороший вопрос?. Это поможет нам помочь вам. Сейчас нам трудно сказать, не понимаете ли вы, как работает БПФ, или просто не знаете, как его кодировать.

Floris 05.09.2018 11:16

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

Scott Stensland 05.09.2018 19:31

Кямран, пожалуйста, посмотрите мой ответ и код ниже. Сигнал с фазовым сдвигом создается тремя строками кода, преобразованием Фурье, умножением на exp (фаза i) и обратным преобразованием.

DrM 12.09.2018 17:26

@DrM большое спасибо большое. Поправьте меня, если я ошибаюсь, умножение на exp (i-фаза) - это просто сложение в степени exp (i w t), я имею в виду в выражении DFT. Я правильно понимаю?

Kamran Abbasov 12.09.2018 17:29

exp (i w t) x exp (i phase) равно exp [i (w t + phi)], то есть сдвигает фазу. Другими словами, после умножения сигнала, преобразованного Фурье, на exp (i-фаза), у вас есть FT сигнала, сдвинутого на ту же самую фазу. Это работает для любого входного сигнала.

DrM 12.09.2018 17:40

Попробуйте, и вы убедитесь, что это работает. И, пожалуйста, не забудьте отметить это правильно и / или проголосовать за него.

DrM 12.09.2018 17:45

@DrM есть ли способ не нормализовать newSignal? Кажется, что новостной сигнал помещается около 0 и не соответствует исходным значениям. Я приложил вывод в описании проблемы выше.

Kamran Abbasov 12.09.2018 20:36

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

DrM 12.09.2018 22:41
Почему в Python есть оператор "pass"?
Почему в Python есть оператор "pass"?
Оператор pass в Python - это простая концепция, которую могут быстро освоить даже новички без опыта программирования.
Некоторые методы, о которых вы не знали, что они существуют в Python
Некоторые методы, о которых вы не знали, что они существуют в Python
Python - самый известный и самый простой в изучении язык в наши дни. Имея широкий спектр применения в области машинного обучения, Data Science,...
Основы Python Часть I
Основы Python Часть I
Вы когда-нибудь задумывались, почему в программах на Python вы видите приведенный ниже код?
LeetCode - 1579. Удаление максимального числа ребер для сохранения полной проходимости графа
LeetCode - 1579. Удаление максимального числа ребер для сохранения полной проходимости графа
Алиса и Боб имеют неориентированный граф из n узлов и трех типов ребер:
Оптимизация кода с помощью тернарного оператора Python
Оптимизация кода с помощью тернарного оператора Python
И последнее, что мы хотели бы показать вам, прежде чем двигаться дальше, это
Советы по эффективной веб-разработке с помощью Python
Советы по эффективной веб-разработке с помощью Python
Как веб-разработчик, Python может стать мощным инструментом для создания эффективных и масштабируемых веб-приложений.
3
8
7 370
1
Перейти к ответу Данный вопрос помечен как решенный

Ответы 1

Ответ принят как подходящий

Вот ваш код с некоторыми изменениями. Мы применяем преобразование Фурье, сдвигаем по фазе преобразованный сигнал, а затем выполняем обратное преобразование Фурье, чтобы получить сдвинутый по фазе сигнал во временной области.

Обратите внимание, что преобразования выполняются с помощью rfft () и irfft (), а фазовый сдвиг выполняется простым умножением преобразованных данных на cmath.rect (1., фаза). Фазовый сдвиг эквивалентен умножению комплексно преобразованного сигнала на exp (i * phase).

На графиках на левой панели мы показываем исходный и смещенный сигналы. Новый сигнал выдвинут на 90 градусов. На правой панели мы показываем спектр мощности на левой оси. В этом примере у нас есть мощность на одной частоте. Фаза изображена на правой оси. Но опять же, поскольку у нас есть мощность только на одной частоте, фазовый спектр показывает шум повсюду.

#!/usr/bin/python

import matplotlib.pyplot as plt
import numpy as np
import cmath

# Generate a model signal
t0 = 1250.0
dt = 0.152
freq = (1./dt)/128

t = np.linspace( t0, t0+1024*dt, 1024, endpoint=False )
signal = np.sin( t*(2*np.pi)*freq )

## Fourier transform of real valued signal
signalFFT = np.fft.rfft(signal)

## Get Power Spectral Density
signalPSD = np.abs(signalFFT) ** 2
signalPSD /= len(signalFFT)**2

## Get Phase
signalPhase = np.angle(signalFFT)

## Phase Shift the signal +90 degrees
newSignalFFT = signalFFT * cmath.rect( 1., np.pi/2 )

## Reverse Fourier transform
newSignal = np.fft.irfft(newSignalFFT)

## Uncomment this line to restore the original baseline
# newSignal += signalFFT[0].real/len(signal)


# And now, the graphics -------------------

## Get frequencies corresponding to signal 
fftFreq = np.fft.rfftfreq(len(signal), dt)

plt.figure( figsize=(10, 4) )

ax1 = plt.subplot( 1, 2, 1 )
ax1.plot( t, signal, label='signal')
ax1.plot( t, newSignal, label='new signal')
ax1.set_ylabel( 'Signal' )
ax1.set_xlabel( 'time' )
ax1.legend()

ax2 = plt.subplot( 1, 2, 2 )
ax2.plot( fftFreq, signalPSD )
ax2.set_ylabel( 'Power' )
ax2.set_xlabel( 'frequency' )

ax2b = ax2.twinx()
ax2b.plot( fftFreq, signalPhase, alpha=0.25, color='r' )
ax2b.set_ylabel( 'Phase', color='r' )


plt.tight_layout()

plt.show()

А вот и графический вывод.

Raw and shifted signals, power and phase spectrum

Другие вопросы по теме