CN104897925A - 超声波风速风向测量装置及测量方法 - Google Patents

Abstract

本发明属于超声波测量领域，涉及一种超声波风速风向测量装置及测量方法；克服了现有技术存在的超声波测风仪器成本较高、对射式仪器机械结构较大、反射式仪器精度有限的问题；所述超声波风速风向测量装置，包括超声波传感器阵列、A/D转换器组、上位机和显示器；超声波传感器阵列由六个超声波传感器组成；所述超声波风速风向测量装置的测量方法，包括以下步骤：1、上位机控制0号超声波传感器发射超声波信号，发射信号为S(t)，发送超声波频率f为200kHz；2、利用A/D转换器组对1号～5号超声波传感器接收到的信号进行采样，输入到上位机中进行算法运算；3、利用上位机控制显示模块，输出显示计算得到的实时风速风向值。

Description

超声波风速风向测量装置及测量方法

技术领域

本发明属于超声波测量领域，涉及一种超声波风速风向测量装置及测量方法。

背景技术

在风速风向测量方面，当前常用的仪器有翼状风速计、杯状风速计、热敏风速计和超声波测风仪等。翼状风速计和杯状风速计使用方便，但他们的机械摩擦阻力和惰性较大，只适用于测定风速较大的情况，对于一些需要测量较小风速的场合则无法试用；常用的机械式仪器因为存在很多的运动部件，所以均存在结构磨损、使用寿命不长、无法在恶劣环境下使用、维护成本高等缺点，此外，受到机械结构以及机械式测风仪器测量原理的限制，机械测风仪器的灵敏度也较低；热敏式风速仪由于自身测量原理的限制，受到环境温度影响较大，因此其适用场合有限；而现有的超声测风仪器大部分均采用的是基于时间差的测量方法，因此，超声测风仪器系统的测量精度完全取决于时间差的测量精度，常用的超声测风仪器采用的超声波传感器类型为对射式和反射式两种，对射式超声测风仪器机械机构较大，不利于安装使用，反射式超声测风仪器在测量精度方面会有一定的限制和误差，而且，超声测风仪器采用的超声波传感器均为收发一体的，这类传感器的市场价格较高，因此，超声测风仪器的成本也较高。

发明内容

本发明所要解决的技术问题是克服了现有技术存在的超声波测风仪器成本较高、对射式仪器机械结构较大、反射式仪器精度有限的问题，提供了一种超声波风速风向测量装置及测量方法；尤其涉及一种利用超声波传感器阵列和阵列信号处理算法来实现风速风向测量的方法。

为解决上述技术问题，本发明是采用如下技术方案实现的，结合附图说明如下：

一种超声波风速风向测量装置，包括超声波传感器阵列、A/D转换器组、上位机和显示器；

所述超声波传感器阵列由六个超声波传感器组成，六个超声波传感器中一个为发射超声波传感器，另外五个为接收超声波传感器；一个发射超声波传感器通过信号线与所述上位机相连，五个接收超声波传感器分别通过信号线与所述A/D转换器组中五个转换器相连；

所述显示器通过信号线与所述上位机相连。

技术方案中所述发射超声波传感器为0号超声波传感器，型号为MA40S4S；所述五个接收超声波传感器为1号超声波传感器～5号超声波传感器，型号为MA40B8S；1号～5号超声波传感器分别通过信号线与A/D转换器组中的1号A/D转换器～5号A/D转换器相连；0号～5号超声波传感器共地。

技术方案中所述0号超声波传感器～5号超声波传感器位于同一个水平面上，所述1号～5号超声波传感器均匀排列在一条直线上，每两个相邻超声波传感器之间间距为d；0号超声波传感器与5号超声波传感器的连线垂直于1号～5号超声波传感器排列的直线，且0号超声波传感器与5号超声波传感器的距离为4d；

所述1号A/D转换器～5号A/D转换器均相同。

一种超声波风速风向测量装置的测量方法，包括以下步骤：

步骤1：上位机控制0号超声波传感器发射超声波信号，发射信号为S(t)，发送超声波频率f为200kHz；

步骤2：利用A/D转换器组对1号～5号超声波传感器接收到的信号进行采样，输入到上位机中进行算法运算；

步骤3：利用上位机控制显示模块，输出显示计算得到的实时风速风向值；

技术方案中所述步骤2中利用A/D转换器组对1号～5号超声波传感器接收到的信号进行采样，输入到上位机中进行算法运算的运算过程如下：

第一步：

得到接收数据X(t)的协方差矩阵R：

R＝E[X(t)X^H(t)]   式(1)

其中：H表示共轭转置；

X(t)为阵列接收信号矩阵，X(t)的形式如下：

X(t)＝[x₁(t),x₂(t)...x_N(t)]^T   式(2)

其中：T表示矩阵转置；

X(t)为x_i(t)的矢量形式；

x_i(t)为第i个阵元接收到的信号，i＝1～5，x_i(t)的形式如下：

$x_i\left(t\right)=s\left(t\right)*e^{-j2{\mathrm{\πf\τ}}_i}+n_i\left(t\right)$    式(3)

其中：s(t)为发射信号；

f为超声波频率；

n_i(t)为第i个传感器上的噪声；

τ_i为第i个阵元与基准阵元之间接收信号的时延差，τ_i的形式如下：

${\mathrm{\τ}}_i=\frac{L_i}{c+V_i}-\frac{L_1}{c+V_1}$    式(4)

其中：L_i为信号源与各传感器之间的距离，选1号传感器作为基准阵元，L₁为信号源与1号传感器之间的距离；

c为声速；

V_i为风速V在信号源与第i个传感器连线上的分量，V₁为风速V在信号源与1号传感器连线上的分量；

运用解三角形的方法得出：

$L_1=4\sqrt{2}d;$ L₂＝5d； $L_3=2\sqrt{5}d;$ $L_4=\sqrt{17}d;$ L₅＝4d；

且 $V_1=(\frac{\sqrt{2}}{2}cos\mathrm{\θ}+\frac{\sqrt{2}}{2}sin\mathrm{\θ})V;V_2=(\frac{4}{5}cos\mathrm{\θ}+\frac{3}{5}sin\mathrm{\θ})V;V_3=(\frac{2\sqrt{5}}{5}cos\mathrm{\θ}+\frac{\sqrt{5}}{5}sin\mathrm{\θ})V;$ $V_4=(\frac{4\sqrt{17}}{17}cos\mathrm{\θ}+\frac{\sqrt{17}}{17}sin\mathrm{\θ})V;$ V₅＝V cosθ；

第二步：

对由式(1)得到的协方差矩阵R进行特征值分解：

R＝UΣU^H   式(5)，

其中：U为特征矢量矩阵；

Σ为由特征值组成的对角阵；

按照特征值的大小顺序，将最大特征值和对应的特征矢量看作信号子空间，将剩下的特征值对应的特征矢量看作噪声子空间，将式(5)整理成：

$R=U_S{\mathrm{\Σ}}_S{U}_S^H+U_N{\mathrm{\Σ}}_N{U}_N^H$    式(6)

其中：U_S为大特征值对应的信号子空间；

Σ_S为大特征值组成的对角阵；

U_N为小特征值对应的噪声子空间；

Σ_N为小特征值组成的对角阵；

第三步：

写出谱计算公式，谱计算公式P_MUSIC的形式如下：

$P_{MUSIC}(\mathrm{\θ},V)=\frac{1}{a^H(\mathrm{\θ},V)U_N{U}_N^{H}a(\mathrm{\θ},V)}$    式(7)

其中：U_N为小特征值对应的噪声子空间；

a(θ,V)为信号向量；

a(θ,V)由如下过程得到：

将式(2)写成矢量形式，有：

X(t)＝a(θ,V)S(t)+N(t)   式(8)

其中：N(t)＝[n₁(t),n₂(t),...,n_N(t)]^T为噪声矩阵，a(θ,V)即为信号向量，其形式为：

$a(\mathrm{\θ},V)={\[1,e^{-j2{\mathrm{\πf\τ}}_2},...,e^{-j2{\mathrm{\πf\τ}}_5}\]}^T$    式(9)

其中：T表示矩阵转置；τ_i由式(4)得到，所以有：

$a(\mathrm{\θ},V)=\left(\begin{array}{c}1\\ \exp \[-j2\mathrm{\π}f(\frac{25d}{5c+(4\cos \mathrm{\θ}+3\sin \mathrm{\θ})V}-\frac{8d}{\sqrt{2}c+(\cos \mathrm{\θ}+\sin \mathrm{\θ})V})\]\\ \exp \[-j2\mathrm{\π}f(\frac{10d}{\sqrt{5}c+(2\cos \mathrm{\θ}+\sin \mathrm{\θ})V}-\frac{8d}{\sqrt{2}c+(\cos \mathrm{\θ}+\sin \mathrm{\θ})V})\]\\ \exp \[-j2\mathrm{\π}f(\frac{17d}{\sqrt{17}c+(4\cos \mathrm{\θ}+\sin \mathrm{\θ})V}-\frac{8d}{\sqrt{2}c+(\cos \mathrm{\θ}+\sin \mathrm{\θ})V})\]\\ \exp \[-j2\mathrm{\π}f(\frac{4d}{c+V\cos \mathrm{\θ}}-\frac{8d}{\sqrt{2}c+(\cos \mathrm{\θ}+\sin \mathrm{\θ})V})\]\end{array}\right)$    式(10)

将式(10)代入式(7)即可得到谱计算公式；

第四步：

由第三步得知，谱计算公式P_MUSIC是风速V和风向角θ的二维函数；在可测风速风向范围内，使V和θ不断地变化，按照谱计算公式来计算谱函数，谱函数P_MUSIC峰值对应的风速V和风向θ即为实时的风速和风向值。

与现有技术相比本发明的有益效果是：

1、本方法通过对超声波传感器阵列接收到的信号进行处理，增强所需要的有用信号，抑制无用的干扰和噪声，并提取有用的信号特征及信号所包含的信息来测量风速和风向，通过算法的仿真实验结果可以看出，运用该方法进行风速风向的测量，可以极大地提高超声波测风系统的抗干扰能力和测量精度。

2、本发明可扩展性强，在已有系统的基础之上增加超声波传感器的数量，构建不同种类的超声波传感器阵列，可以实现更大范围内的高精度、同步多点风速风向的测量。

3、本发明采用单发、单收超声波传感器，在一定程度上降低了系统成本。

附图说明

下面结合附图对本发明作进一步的说明：

图1是本发明所述的超声波风速风向测量装置结构原理示意框图；

图2是本发明所述的超声波风速风向测量装置中传感器阵列结构示意图；

图3是超声波传感器阵列测量风速风向系统总体流程图；

图4是多重信号分类算法示意图；

图5是风速估计的均方根误差；

图6是风向估计的均方根误差。

具体实施方式

下面结合附图对本发明作详细的描述：

参阅图1，为本发明的系统整体结构框图。本系统包括超声波传感器阵列、A/D转换器组、上位机与显示器。

超声波传感器阵列由6个超声波传感器组成，0号超声波传感器为发射超声波传感器，型号为MA40S4S，0号超声波传感器通过信号线与上位机相连，上位机通过信号线控制0号超声波传感器发送超声波；1号超声波传感器～5号超声波传感器为接收超声波传感器，型号为MA40B8S，1号～5号超声波传感器分别通过信号线与A/D转换器组中的1号A/D转换器～5号A/D转换器相连；0号～5号超声波传感器共地；0号～5号超声波传感器位置摆放参阅图2，本发明采用的超声波传感器阵列结构为均匀线阵形式，0号超声波传感器～5号超声波传感器位于同一个水平面上，0号超声波传感器位于5号传感器正上方，距离为4d；1号～5号超声波传感器均匀排列在一条直线上，每两个相邻超声波传感器之间间距为d。

A/D转换器组中包括5个完全相同的A/D转换器，1号A/D转换器～5号A/D转换器分别通过信号线与1号～5号超声波传感器相连，用来采集超声波接收信号并经A/D转换后，通过信号线输入到上位机中。

上位机通过信号线控制0号传感器发送超声波，A/D转换器组将1号～5号超声波传感器接收到超声波信号转换为数字信号后，输入到上位机中进行算法运算，并将计算结果通过显示器显示出来。

显示器通过信号线与上位机相连，用来显示风速风向的计算结果，对实时风速风向信息进行监控。

利用超声波传感器阵列结合多重信号分类算法计算风速风向的步骤。

参阅图3，为本发明计算风速风向方法流程图。现归纳利用超声波传感器均匀线阵结合多重信号分类算法计算风速风向的过程如下：

第一步：上位机控制0号超声波传感器发射超声波信号，发射信号为S(t)，发送超声波频率f为200kHz。

第二步：利用A/D转换器组对超声波传感器1～5接收到的信号进行采样，输入到上位机中进行算法运算，该算法计算过程如下：

首先设计阵列结构的近场声传播模型。

参阅图4，为本发明采用的多重信号分类算法示意图。当信号位于远场时，信号到达阵列时的波阵面可以近似看成平面波，这时等间距线阵中相邻阵元接收到的信号时间延迟大小只与信号的方向角有关。而对于本系统来说，信号源处于近场情况，到达波不能看成平面波，而应看成球面波，等间距均匀线阵中相邻阵元接收到信号时间延迟的大小，不仅与信号的方向角有关，还与信号源到阵列的距离有关，又因为发射源位置已知，因此阵列接收矢量仅与风速风向有关。

本系统采用的均匀线阵由5个阵元组成，阵元间距为d，选1号超声波传感器作为基准阵元，信号源与各阵元的距离为L_i(i＝1～5)，c为声速，设风速为V，风向角为θ，V_i为风速在信号源与第i个阵元连线上的分量。所以，第i个阵元与基准阵元接收信号的时延差：

${\mathrm{\τ}}_i=\frac{L_i}{c+V_i}-\frac{L_1}{c+V_1}$    式(1)

由图3的几何关系，运用解三角形的方法很容易得出：L₂＝5d； $L_4=\sqrt{17}d;$ L₅＝4d；且 $V_1=(\frac{\sqrt{2}}{2}cos\mathrm{\θ}+\frac{\sqrt{2}}{2}sin\mathrm{\θ})V;V_2=(\frac{4}{5}cos\mathrm{\θ}+\frac{3}{5}sin\mathrm{\θ})V;$ $V_3=(\frac{2\sqrt{5}}{5}cos\mathrm{\θ}+\frac{\sqrt{5}}{5}sin\mathrm{\θ})V;V_4=(\frac{4\sqrt{17}}{17}cos\mathrm{\θ}+\frac{\sqrt{17}}{17}sin\mathrm{\θ})V;$ V₅＝V cosθ。

那么有：

τ₁＝0   式(2)

${\mathrm{\τ}}_2=\frac{25d}{5c+(4cos\mathrm{\θ}+3sin\mathrm{\θ})V}-\frac{8d}{\sqrt{2}c+(\cos \mathrm{\θ}+sin\mathrm{\θ})V}$    式(3)

${\mathrm{\τ}}_3=\frac{10d}{\sqrt{5}c+(2\cos \mathrm{\θ}+sin\mathrm{\θ})V}-\frac{8d}{\sqrt{2}c+(\cos \mathrm{\θ}+sin\mathrm{\θ})V}$    式(4)

${\mathrm{\τ}}_4=\frac{17d}{\sqrt{17}c+(4\cos \mathrm{\θ}+sin\mathrm{\θ})V}-\frac{8d}{\sqrt{2}c+(\cos \mathrm{\θ}+sin\mathrm{\θ})V}$    式(5)

${\mathrm{\τ}}_5=\frac{4d}{c+Vcos\mathrm{\θ}}-\frac{8d}{\sqrt{2}c+(\cos \mathrm{\θ}+sin\mathrm{\θ})V}$    式(6)

第i个阵元接收信号为：

$x_i\left(t\right)=s\left(t\right)*e^{-j2{\mathrm{\πf\τ}}_i}+n_i\left(t\right)$    式(7)

x_i(t)为第i个阵元接收到的信号；f为超声波频率，n_i(t)为第i个传感器上的噪声。将式(7)写成矢量形式有，

X(t)＝a(θ,V)S(t)＋N(t)   式(8)

式中，X(t)＝[x₁(t),x₂(t)...x_N(t)]^T为阵列接收信号矩阵，N(t)＝[n₁(t),n₂(t),...,n_N(t)]^T为噪声矩阵，信号向量为：

$a(\mathrm{\θ},V)={\[1,e^{-j2{\mathrm{\πf\τ}}_2},...,e^{-j2{\mathrm{\πf\τ}}_5}\]}^T$    式(9)

式(9)中的τ₂～τ₅可由式(2)～(6)得到，所以：

$a(\mathrm{\θ},V)=\left(\begin{array}{c}1\\ \exp \[-j2\mathrm{\π}f(\frac{25d}{5c+(4\cos \mathrm{\θ}+3\sin \mathrm{\θ})V}-\frac{8d}{\sqrt{2}c+(\cos \mathrm{\θ}+\sin \mathrm{\θ})V})\]\\ \exp \[-j2\mathrm{\π}f(\frac{10d}{\sqrt{5}c+(2\cos \mathrm{\θ}+\sin \mathrm{\θ})V}-\frac{8d}{\sqrt{2}c+(\cos \mathrm{\θ}+\sin \mathrm{\θ})V})\]\\ \exp \[-j2\mathrm{\π}f(\frac{17d}{\sqrt{17}c+(4\cos \mathrm{\θ}+\sin \mathrm{\θ})V}-\frac{8d}{\sqrt{2}c+(\cos \mathrm{\θ}+\sin \mathrm{\θ})V})\]\\ \exp \[-j2\mathrm{\π}f(\frac{4d}{c+V\cos \mathrm{\θ}}-\frac{8d}{\sqrt{2}c+(\cos \mathrm{\θ}+\sin \mathrm{\θ})V})\]\end{array}\right)$    式(10)

由式(10)可以看出，方向矢量是风速V和风向θ的二维函数。

多重信号分类(MUSIC)算法的基本思想是将任意阵列输出数据的协方差矩阵进行特征分解，从而得到与信号分量相对应的信号子空间和与信号相正交的噪声子空间，利用这两个空间的正交性来估计信号的参数。

对于本系统，信号源个数为1，因此，A(θ,V)＝a(θ,V)。假定各阵元的噪声是均值为零的平稳白噪声过程，方差为σ²，并且噪声之间不相关，且与信号不相关。因此接收数据X(t)的协方差矩阵为：

R＝E[X(t)X^H(t)]   式(11)

其中上标H表示共轭转置。

R＝AR_sA^H+σ²I   式(12)

由于信号与噪声相互独立，数据协方差矩阵可以分解为与信号、噪声相关的两部分，其中R_s是信号的协方差矩阵。

阵列协方差矩阵R可以划分为两个空间：

$R=U_S{\mathrm{\Σ}}_S{U}_S^H+U_N{\mathrm{\Σ}}_N{U}_N^H$    式(13)

式(13)中，U_S代表信号子空间，U_N代表噪声子空间。在理想条件下，数据空间中的信号子空间和和噪声子空间是相互正交的，即信号子空间的导向矢量也与噪声子空间正交：

a^H(θ,V)U_N＝0   式(14)

由于噪声存在，a^H(R₁,θ)与U_N并不能完全正交，也就是说式(14)并不成立，实际求DOA是以最小优化搜索实现的，当信号源位于近场时，二维MUSIC定位算法的谱估计公式为：

$P_{MUSIC}=\frac{1}{a^H(\mathrm{\θ},V)U_N{U}_N^{H}a(\mathrm{\θ},V)}$    式(15)

对目标(θ,V)进行二维谱峰搜索，得到定位声图，谱峰对应的二维参数即为风速V以及风向θ。

经过以上讨论，将第二步中进行的多重信号分类算法运算步骤归纳如下：

(1)根据5个超声波传感器接收到的信号矢量，得到协方差矩阵的估计值：

$R=\frac{1}{5}\underset{n=1}{\overset{5}{\Σ}}x\left(n\right)x^H\left(n\right)$    式(16)

对式(16)得到的协方差矩阵进行特征值分解：R＝UΣU^H。

(2)按照特征值的大小顺序，把最大特征值和对应的特征矢量看作信号子空间，把剩下的特征值对应的特征矢量看作噪声子空间。即

(3)使(θ,V)变化，按照来计算谱函数，通过寻求峰值得到实时的风速和风向值。

第三步，利用上位机控制显示模块，输出显示计算得到的实时风速风向值。

实施例

根据本发明提供的估计风速风向的方法步骤，设计仿真实验，分析随信噪比变化时，风速风向估计的均方根误差，以此来说明该方法具有较好的抗干扰能力和较高的测量精度。

实验条件：均匀线阵超声波传感器阵列由5个超声波传感器构成，传感器间距d＝0.05m。超声波传感器0作为声源，发射超声波入射到超声波传感器阵列中，超声波传感器0位于超声波传感器5正上方，距离为4d＝0.2m，风速风向信息为：风速＝10.2m/s；风向＝141°。超声波阵元噪声为加性高斯白噪声。快拍数为100，信噪比从0dB变化到25dB，间隔为5dB。对每种情况下做20次独立的仿真实验，估计性能使用均方根误差(RMSE)来评判，信号参数估计的RMSE定义如下：

$RMSE=\sqrt{E\left\{{(\hat{a}-a_{true})}^2\right\}}$    式(17)

式中为参数的估计值，a_true为参数的真实值。

仿真结果如图5、图6所示。

参阅图5，在信噪比为10dB时，本发明提供的方法估计出的风速信息，误差低于0.2m/s，满足精度要求。

参阅图6，在信噪比为10dB时，本发明提供的方法估计出的风向信息，误差低于1°，满足精度要求。

综上所述，利用本发明提供的基于均匀线阵结合多重信号分类算法估计风向风速的方法，是可行的，且在较低信噪比的情况下，其估计精度仍是可观的。

Claims (5)

1.一种超声波风速风向测量装置，其特征在于：包括超声波传感器阵列、A/D转换器组、上位机和显示器；

所述超声波传感器阵列由六个超声波传感器组成，六个超声波传感器中一个为发射超声波传感器，另外五个为接收超声波传感器；一个发射超声波传感器通过信号线与所述上位机相连，五个接收超声波传感器分别通过信号线与所述A/D转换器组中五个转换器相连；

所述显示器通过信号线与所述上位机相连。

2.根据权利要求1所述的一种超声波风速风向测量装置，其特征在于：

所述发射超声波传感器为0号超声波传感器，型号为MA40S4S；所述五个接收超声波传感器为1号超声波传感器～5号超声波传感器，型号为MA40B8S；1号～5号超声波传感器分别通过信号线与A/D转换器组中的1号A/D转换器～5号A/D转换器相连；0号～5号超声波传感器共地。

3.根据权利要求2所述的一种超声波风速风向测量装置，其特征在于：

所述0号超声波传感器～5号超声波传感器位于同一个水平面上，所述1号～5号超声波传感器均匀排列在一条直线上，每两个相邻超声波传感器之间间距为d；0号超声波传感器与5号超声波传感器的连线垂直于1号～5号超声波传感器排列的直线，且0号超声波传感器与5号超声波传感器的距离为4d；

所述1号A/D转换器～5号A/D转换器均相同。

4.如权利要求1所述的一种超声波风速风向测量装置的测量方法，其特征在于，包括以下步骤：

步骤1：上位机控制0号超声波传感器发射超声波信号，发射信号为S(t)，发送超声波频率f为200kHz；

步骤2：利用A/D转换器组对1号～5号超声波传感器接收到的信号进行采样，输入到上位机中进行算法运算；

步骤3：利用上位机控制显示模块，输出显示计算得到的实时风速风向值。

5.根据权利要求4所述的一种超声波风速风向测量装置的测量方法，其特征在于：

所述步骤2中利用A/D转换器组对1号～5号超声波传感器接收到的信号进行采样，输入到上位机中进行算法运算的运算过程如下：

第一步：

得到接收数据X(t)的协方差矩阵R：

R＝E[X(t)X^H(t)]          式(1)

其中：H表示共轭转置；

X(t)为阵列接收信号矩阵，X(t)的形式如下：

X(t)＝[x₁(t),x₂(t)…x_N(t)]^T           式(2)

其中：T表示矩阵转置；

X(t)为x_i(t)的矢量形式；

x_i(t)为第i个阵元接收到的信号，i＝1～5，x_i(t)的形式如下：

             式(3)

其中：s(t)为发射信号；

f为超声波频率；

n_i(t)为第i个传感器上的噪声；

τ_i为第i个阵元与基准阵元之间接收信号的时延差，τ_i的形式如下：

               式(4)

其中：L_i为信号源与各传感器之间的距离，选1号传感器作为基准阵元，L₁为信号源与1号传感器之间的距离；

c为声速；

V_i为风速V在信号源与第i个传感器连线上的分量，V₁为风速V在信号源与1号传感器连线上的分量；

运用解三角形的方法得出：

L₂＝5d； L₅＝4d；

且

V₅＝V cosθ；

第二步：

对由式(1)得到的协方差矩阵R进行特征值分解：

R＝U∑U^H              式(5)

其中：U为特征矢量矩阵；

∑为由特征值组成的对角阵；

按照特征值的大小顺序，将最大特征值和对应的特征矢量看作信号子空间，将剩下的特征值对应的特征矢量看作噪声子空间，将式(5)整理成：

                式(6)

其中：U_S为大特征值对应的信号子空间；

∑_S为大特征值组成的对角阵；

U_N为小特征值对应的噪声子空间；

∑_N为小特征值组成的对角阵；

第三步：

写出谱计算公式，谱计算公式P_MUSIC的形式如下：

            式(7)

其中：U_N为小特征值对应的噪声子空间；

a(θ,V)为信号向量；

a(θ,V)由如下过程得到：

将式(2)写成矢量形式，有：

X(t)＝a(θ,V)S(t)+N(t)         式(8)

其中：N(t)＝[n₁(t),n₂(t),…,n_N(t)]^T为噪声矩阵，a(θ,V)即为信号向量，其形式为：

              式(9)

其中：T表示矩阵转置；τ_i由式(4)得到，所以有：

 式(10)

将式(10)代入式(7)即可得到谱计算公式；

第四步：

由第三步得知，谱计算公式P_MUSIC是风速V和风向角θ的二维函数；在可测风速风向范围内，使V和θ不断地变化，按照谱计算公式来计算谱函数，谱函数P_MUSIC峰值对应的风速V和风向θ即为实时的风速和风向值。