BASS(同步信号块处理)方法:利用相邻的数据块的不同相位差值除以时间得到更精确的载波多普勒变化。和用pll最大不同是,不用反馈。
1.对每段数据,我们可以用dft的定义来得到频域的值:
X
(
k
)
=
∑
n
=
0
N
−
1
x
(
n
)
e
−
j
2
π
n
k
/
N
X(k)=\sum_{n=0}^{N-1} x(n) e^{-j 2 \pi n k / N}
X(k)=n=0∑N−1x(n)e−j2πnk/N
现在,假设x(n)是复数信号
e
j
(
2
π
f
n
+
π
/
4
)
{{\rm{e}}^{j(2\pi fn+\pi/4)}}
ej(2πfn+π/4),初相为
π
/
4
\pi/4
π/4,那么最大幅度的|X(ki)|就对应这个正弦信号的频率。同时我们还可以容易得到余弦信号的初始相位:
θ
=
tan
−
1
(
Im
[
X
(
k
i
)
]
Re
[
X
(
k
i
)
]
)
\theta=\tan ^{-1}\left(\frac{\operatorname{Im}\left[X\left(k_{i}\right)\right]}{\operatorname{Re}\left[X\left(k_{i}\right)\right]}\right)
θ=tan−1(Re[X(ki)]Im[X(ki)])
证明:我们假设f=0.25,采样频率fs=1,采样点数N=100,作DFT:
X
(
k
)
=
∑
n
=
0
N
−
1
x
(
n
)
e
−
j
2
π
n
k
/
N
=
∑
n
=
0
N
−
1
e
j
(
2
π
0.25
n
+
π
/
4
)
e
−
j
2
π
n
k
/
100
\begin{array}{l} X(k) = \sum\limits_{n = 0}^{N - 1} x (n){e^{ - j2\pi nk/N}}\\ = \sum\limits_{n = 0}^{N - 1} {{{\rm{e}}^{j(2\pi 0.25n+\pi/4)}}} {e^{ - j2\pi nk/100}} \end{array}
X(k)=n=0∑N−1x(n)e−j2πnk/N=n=0∑N−1ej(2π0.25n+π/4)e−j2πnk/100
我们知道X(100)对应于fs=1处(这里也可以看出来matlabFFT其实没有求到fs处,只求到X(99)处还差一个点。从另一方面也可以看出来,FFT分辨率是fs/N ,第一个点X(0)频率对应0,那么第100个点X(99)对应99fs/N啊,不是fs。即范围为半开半闭区间[0,2
π
\pi
π),那么X(25)就对应于f=0.25处,所以我们可以直接求X(25)
X
(
25
)
=
∑
n
=
0
100
−
1
x
(
n
)
e
−
j
2
π
n
25
/
N
=
∑
n
=
0
100
−
1
e
j
(
2
π
0.25
n
+
π
/
4
)
e
−
j
2
π
n
25
/
100
=
∑
n
=
0
100
−
1
e
j
[
2
π
(
0.25
−
25
/
100
)
n
+
π
/
4
]
=
e
j
π
/
4
∑
n
=
0
100
−
1
e
j
[
2
π
(
0.25
−
25
/
100
)
n
]
=
e
j
π
/
4
\begin{array}{l} X(25) = \sum\limits_{n = 0}^{100 - 1} x (n){e^{ - j2\pi n25/N}}\\ = \sum\limits_{n = 0}^{100 - 1} {{{\rm{e}}^{j(2\pi 0.25n + \pi /4)}}} {e^{ - j2\pi n25/100}}\\ = \sum\limits_{n = 0}^{100 - 1} {{{\rm{e}}^{j\left[ {2\pi (0.25 - 25/100)n + \pi /4} \right]}}} \\ = {e^{j\pi /4}}\sum\limits_{n = 0}^{100 - 1} {{{\rm{e}}^{j\left[ {2\pi (0.25 - 25/100)n} \right]}}} \\ = {e^{j\pi /4}} \end{array}
X(25)=n=0∑100−1x(n)e−j2πn25/N=n=0∑100−1ej(2π0.25n+π/4)e−j2πn25/100=n=0∑100−1ej[2π(0.25−25/100)n+π/4]=ejπ/4n=0∑100−1ej[2π(0.25−25/100)n]=ejπ/4
最后结果是
e
j
π
/
4
{e^{j\pi /4}}
ejπ/4,其幅值确实是归一化最大值,相位确实是原信号初相.
postscript:原来我以为matlab求FFT的结果范围是[0,
2
π
2\pi
2π] (即[0,fs]),其实不是,是[0,fs)。
通常情况下,如果对输入信号完全不知道,那么就需要计算所有的k值(0~N-1),然后比较求最大幅度的值,这就必须借助FFT减少计算次数;如果输入信号的频率知道一些,即:有1khz或100hz左右的分辨率也就是误差啦,这时候,我们就没必要计算所有的k值,完全可以根据目前的频率值(有误差)/频率分辨率=k值,然后直接用dft计算一个X(k),X(k)值即可找到最大的幅度值的位置。但值得注意的是,我们这里不关心幅度,而是关心相位,也就是说,连续的时间段内,我们用上面的dft的方法都可以得到一个X(k),虽然不同时间段计算出来的|X(k)|相同,但是θ(k)肯定不同,这个θ(k)的差值就可以得到更精确的多普勒变化,即:
f
=
δ
θ
2
π
m
=
θ
m
+
n
−
θ
n
2
π
m
f=\frac{\delta \theta}{2\pi m}=\frac{\theta_{m+n}-\theta_{n}}{2\pi m}
f=2πmδθ=2πmθm+n−θn
这里m是时间,单位是秒。
我们知道这么一个对应关系,如果
2
π
2\pi
2π对应于信号的一个周期T=1/f,那么
δ
θ
\delta \theta
δθ就对应
δ
t
\delta t
δt,这里的
δ
t
\delta t
δt就是上面的m,即:
2
π
T
=
2
π
1
f
=
δ
θ
m
=
>
f
=
δ
θ
2
π
m
\begin{array}{l} \frac{{2\pi }}{T} = \frac{{2\pi }}{{\frac{1}{f}}} = \frac{{\delta \theta }}{m}\\ = > f = \frac{{\delta \theta }}{{2\pi m}} \end{array}
T2π=f12π=mδθ=>f=2πmδθ
BASS 捕获:
第一步:伪码并行算法,1ms数据获得分辨率为1k多普勒的捕获结果和伪码捕获结果
第二步:剥离伪码,用10ms数据FFT得到分辨率为100Hz的多普勒结果
第三步:BASS法,用每毫秒数据的初始相位差得到更精确的频率估计值。
对每毫秒数据作DFT
X
(
k
)
=
∑
n
=
0
N
−
1
x
(
n
)
e
−
j
2
π
n
k
/
N
X(k)=\sum_{n=0}^{N-1} x(n) e^{-j 2 \pi n k / N}
X(k)=n=0∑N−1x(n)e−j2πnk/N
然后
θ
=
tan
−
1
(
Im
[
X
(
k
i
)
]
Re
[
X
(
k
i
)
]
)
\theta=\tan ^{-1}\left(\frac{\operatorname{Im}\left[X\left(k_{i}\right)\right]}{\operatorname{Re}\left[X\left(k_{i}\right)\right]}\right)
θ=tan−1(Re[X(ki)]Im[X(ki)])
得到每毫秒的初相,然后利用BASS公式得到频率值
有一些详细问题比如有导航信息信息位翻转问题等在此不详表。
参考:1.《BASS软件GPS数据处理方法》白伟华
2.https://blog.csdn.net/xinqrs01/article/details/55798482