同步两个相位不同的嘈杂正弦信号

Question

同步两个相位不同的嘈杂正弦信号

3

我目前有两个不同相位的嘈杂正弦波，我正在寻找一个将它们设置为同相位的解决方案。首先，我试图通过查找上升趋势的零交叉值来调整正弦波，然而，由于信号中存在一些噪声，这种方法不可行。我制作了一个小的模拟示例:

time = np.arange(0, 60, 0.005)
amp = 2
F = 0.25
signal1 = amp*np.sin(2*np.pi*F*time + np.pi)
signal2 = amp*np.sin(2*np.pi*F*time + np.pi/3)

mu, sigma = 0.2, 0.1 # mean and standard deviation
noise1 = np.random.normal(mu, sigma, len(signal1))
noise2 = np.random.normal(mu, sigma, len(signal2))


noisy_signal1 = signal1 + noise1
noisy_signal2 = signal2 + noise2

plt.figure()
plt.plot(time, noisy_signal1)
plt.plot(time, noisy_signal2)
plt.show()

当然实际上，我不知道相移是哪个，但我预先知道每个正弦波的振幅和频率。我想要实现的是，从运行代码时可以看到的图形，转换成这样的图形：

编辑1：使用相关函数可能是一个适合的解决方案。然而，由于信号是嘈杂的，最大相关可能会在其中一个信号的几个周期之前找到，因此在进行移位时会丢失很多数据点。我特别寻求一种识别相位并只移动一个正弦波的解决方案。

- Marc Schwambach

我更新了答案。如果已知或计算出频率，则只需要删除与最大位置模数和一个周期中的样本数相对应的样本数。 - Deepstop

更新后，使用用于计算频率的FFT结果来提供相位差，这将快几个数量级。 - Deepstop

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Deepstop · Accepted Answer

感谢您提出这个问题。我正在研究交叉相关的应用，这是我第一次尝试编写一些代码。

交叉相关重复卷积两个信号并提供（长）数组作为输出，每个数组表示在任何给定时间差异下两个信号的相关程度。对于像您的正弦信号这样的信号，几个周期可能就足够了，而不是将整个数组进行相关，降采样也可能是合理的，以提高执行时间。此外，我使用的是numpy correlate函数而不是scipy函数，后者使用FFT来更快地处理大型数组。这应该是一个迷人的探索领域。

要移动第二个信号的相位，首先将信号进行交叉相关，然后找到具有最大值的点。第二个数组通过该最大值的索引旋转，以使两个数组同步，从而达到最大的相关性。旋转假定完整的周期，因此移位和填充零可能更通用。

这只是3行代码（numpy真是太棒了），插入在绘图之前。可能还有我没有考虑到的边缘情况，但这正是您想要的方式。

corr = np.correlate(noisy_signal1, noisy_signal2, 'full')
corr_max_idx = np.argmax(corr)
noisy_signal2 = np.roll(noisy_signal2, corr_max_idx)

要获得一个非常相似的结果且损失不到半个周期，您只需在除以样本频率时找到 argmax 计算的余数即可。我们可以从您的示例程序中了解到频率是多少，但我还是使用FFT进行计算。虽然 FFT 不是计算单个频率的推荐方法，但我发现信号太嘈杂了，无法通过零交叉点可靠地计算每个周期的数量。

以下是我的计算方法。我稍微改变了示例代码，以便将采样率保存供后续使用，尽管我最终没有用到它。

import numpy as np
from scipy import fft
from matplotlib import pyplot as plt

sr = 200
time = np.arange(0, 60, 1/sr)

amp = 2
F = 0.25
omega = 2 * np.pi * F

signal1 = amp*np.sin(omega * time + np.pi)
signal2 = amp*np.sin(omega * time + np.pi/3)

mu, sigma = 0.2, 0.1 # mean and standard deviation
noise1 = np.random.normal(mu, sigma, len(signal1))
noise2 = np.random.normal(mu, sigma, len(signal2))


noisy_signal1 = signal1 + noise1
noisy_signal2 = signal2 + noise2

# Get frequency in samples
fft1 = fft.fft(noisy_signal1)
fft2 = fft.fft(noisy_signal2)
max_f = np.argmax(np.abs(fft1))
f = len(noisy_signal1) / max_f

# Perform correlation to determine best fit for phase difference (in samples)
corr = np.correlate(noisy_signal1, noisy_signal2, 'full')
delta_p = int(round(np.argmax(corr) % f))

# Chop phase difference (maximum 1/2 cycle)
fi = int(round(f,0))
if delta_p < fi/ 2:
    noisy_signal1 = noisy_signal1[delta_p:]
    noisy_signal2 = noisy_signal2[:-delta_p]
    time = time[delta_p:]
else:
    noisy_signal1 = noisy_signal1[:delta_p - fi]
    noisy_signal2 = noisy_signal2[fi - delta_p:]
    time = time[:delta_p - fi]

plt.figure()
plt.plot(time, noisy_signal1)
plt.plot(time, noisy_signal2)
plt.show()

为了消除计算成本高昂的相关性，我们也可以使用用于计算频率的FFT，将两行相关性代码更改为：

# Do the same calculation with the FFT
p1 = np.angle(fft1[max_f])
p2 = np.angle(fft2[max_f])
delta_p = int((f * ((p2 - p1) / (2 * np.pi))) % f)

如果已知频率，不需要相关性或FFT的另一种可能性是使用单个bin FFT，然后测量相位角之间的差异以获取相位差，例如：

mix = np.exp(-1j * omega * time)
p1 = np.angle(noisy_signal1 * mix)
p2 = np.angle(noisy_signal2 * mix)
delta_p = int((f * (np.mean(p2 - p1) / (2 * np.pi))) % f)