我在做科研项目时,遇到过使用相位信息,但是相位被折叠(wrapped)的情况,相位计算公式如下:
Θ
=
a
r
c
t
a
n
(
b
a
)
\Theta = arctan(\frac{b}{a})
Θ=arctan(ab)
其中b、a分别代表虚部和实部,由于反正切函数的区间是
(
−
π
2
,
π
2
)
(-\frac{\pi}{2},\frac{\pi}{2})
(−2π,2π),故实际计算相位时会将结果折叠在这个区间内。如下图,这张图使用的是matlab
的atan
函数计算并绘制的,可见,相位范围被折叠在
(
−
π
2
,
π
2
)
(-\frac{\pi}{2},\frac{\pi}{2})
(−2π,2π)区间内,且在某些点出现了跳跃。
仔细观察跳跃点,不难发现,其跳跃要么是
+
π
+\pi
+π,要么是
−
π
-\pi
−π,通过前后点的差值可以判断是否发生跳跃。现在为了说明怎么相位展开,先作以下符号说明:
x
w
(
n
)
x_{w}(n)
xw(n)表示被折叠的相位信号,
x
u
(
n
)
x_{u}(n)
xu(n)表示展开的相位信号,则相位展开步骤如下:
- 令 x u ( n ) x_{u}(n) xu(n) = x w ( n ) x_{w}(n) xw(n)
- 从 x w ( n ) x_{w}(n) xw(n)的第二个点开始,计算当前点与前一个点差
- 若差大于 π 2 \frac{\pi}{2} 2π,则 x u ( n ) x_{u}(n) xu(n)当前点及后续所有点减去 π \pi π
- 若差小于 − π 2 -\frac{\pi}{2} −2π,则 x u ( n ) x_{u}(n) xu(n)当前点及后续所有点加上 π \pi π
- 重复上述步骤,直至遍历
x
u
(
n
)
x_{u}(n)
xu(n)
对应的matlab
代码如下:
xu = xw;
for i = 2:length(xw)
diff = xw(i) - xw(i-1);
if diff > pi/2
xu(i:end) = xu(i:end) - pi;
elseif diff < -pi/2
xu(i:end) = xu(i:end) + pi;
end
end
展开后相位如下图所示,可见展开后的相位是连续的。
另:若使用atan2
计算相位,上述代码需做一些改变:所有的
π
\pi
π–>
2
π
2\pi
2π,所有的
π
2
\frac{\pi}{2}
2π–>
π
\pi
π,具体如下代码所示。
xu = xw;
for i = 2:length(xw)
diff = xw(i) - xw(i-1);
if diff > pi
xu(i:end) = xu(i:end) - 2*pi;
elseif diff < -pi
xu(i:end) = xu(i:end) + 2*pi;
end
end