一种感应器列阵来波方向的估计方法与流程

本发明属于列阵信号处理
技术领域：
，具体涉及一种感应器列阵来波方向的估计方法。
背景技术：
：在语音学、地震学、声纳和雷达系统领域中，doa(directionofarrival，doa)估计技术扮演着十分重要的角色。然而，doa估计技术中，估计感应器列阵阵元中的每一帧方向是估计感应器来波方向的基础，所以，准确的确定感应器列阵阵元中每一帧的方向对估计感应器来波方向起到重要的作用。在doa技术中目前普遍应用最广的是mvdr(minimumvariancedistortionlessresponse)空间普估计算法，该算法运算规律性强，算法鲁棒性优异，可以使用各种信号条件下的来波方向估计，但是其缺点也很明显，由于需要通过二维搜索确定来波方向，因此运算量巨大，很难集成进低功耗设备中进行实时运算，限制了该算法的进一步应用。现有超分辨率估计算法中，多重信号分类法music(multiplesignalclassification，music)根据信号特征空间和噪声特征空间正交的特性，搜索所有可能的方向候选值来确定最适合的值作为doa估计，music算法普遍适合任意列阵，但往往需要进行多维的空间搜索，之后再与基础数据进行比对，运算量相当大，实时定位困难。技术实现要素：为了解决现有技术中存在的上述问题，本发明提供了一种感应器列阵来波方向的估计方法。本发明要解决的技术问题通过以下技术方案实现：本发明实施例提供了一种感应器列阵来波方向的估计方法，包括步骤：s1.获取待测n列阵阵元，n为正整数；s2.对所述待测n列阵阵元中的每一列阵阵元接收到的信号数据进行分帧，所述每一列阵阵元均得到m帧时域数据，m为正整数；s3.对各所述每一列阵阵元的m帧时域数据进行转换，得到每一帧时域数据对应的频域数据，并设定初始方向估计值；s4.将各列阵阵元中的第一帧均作为当前帧，获取各所述当前帧的频域数据中第一频点的数据，将所述第一频点作为当前频点；s5.根据所述当前频点，计算得到所述当前频点数据的频域快拍信号；s6.通过所述频域快拍信号，计算得到所述当前频点的协方差矩阵的逆矩阵；将所述初始方向估计值作为当前方向估计值；s7.根据预设方法计算，得到所述当前频点相对于所述当前方向估计值的调整量；s8.将各所述当前帧所有频点的目标频点中除所述第一频点外的所有频点均作为当前频点，重复步骤s5-s7，计算所述当前帧中的目标频点的方向估计值调整量，并将各所述调整量求和，得到调整量累计值；s9.判断所述调整量累计值的绝对值是否大于给定的门限阈值的大小；如果大于，根据预设方法更新所述当前方向估计值，直至满足预设条件，停止更新，将当前方向估计值确定为来波方向估计值，所述来波方向估计值的方向矢量方位角与仰角确定为当前帧的信号来波方向；如果不大于，将所述当前方向估计值的方向矢量方位角与仰角确定为当前帧的信号来波方向；s10.将各列阵阵元中的每一帧均作为当前帧，重复步骤s4-s9，确定各个列阵中每一帧的来波方向。在一个具体实施方式中，所述目标频点为所述当前帧的所有频点或所述当前帧的部分频点。在一个具体实施方式中，所述待测n列阵阵元为任意几何形状的列阵阵元，其阵元位置表达式为：p＝{pi}；pi＝[xi，yi，zi]式中，p为阵元位置，i为阵元序号，i＝1-n，pi为阵元位置坐标。在一个具体实施方式中，通过快速傅氏变换运算对各所述m帧时域数据进行转换，得到每一帧时域数据对应的频域数据。在一个具体实施方式中，所述当前频点的数据的频域快拍信号的表达式为：[x(wl，1)，x(wl，2)，…，x(wl，n)]式中，x(wl，1)为当前帧第一列阵阵元的当前帧的当前频点的频域数据，x(wl，2)为当前帧第二列阵阵元的当前帧的当前频点的频域数据。x(wl，n)为当前帧第n列阵阵元的当前帧的当前频点的频域数据，wl为当前帧的当前频点。在一个具体实施方式中，所述初始方向估计值的表达式为：[θ(0)，φ(0)]＝[θinit，φinit]式中，θ为来波方向矢量方位角，ф为来波仰角。在一个具体实施方式中，所述步骤s6中，所述预设方法为：其中，ha(pij，θ，φ)＝(-xijsinθcosφ+yijcosθcosφ)he(pij，θ，φ)＝(-xijcosθsinφ+yijsinθsinφ+zijcosφ)h(φij，pij，wl)＝sin(φij+wlτij)式中pi＝[xi，yi，zi]为阵元i在直角坐标系中的位置坐标；pij＝[xij,yij,zij]为阵元i到阵元j的相对位置坐标，xij＝xi-xj；yij＝yi-yj；zij＝zi-zj；aij，фij分别为信号协方差矩阵的逆矩阵元素的幅值和相角，τi,j＝τi,0-τj,0为列阵阵元i和列阵阵元j之间的相对延时，由列阵阵元的位置坐标及信号来波方向确定，c为信号传播速度。在一个具体实施方式中，累计所述目标频点分别对应的所述调整量，得到调整量累计值的计算公式为：δθ＝∑lδθwlδφ＝∑lδφwl式中，wl为当前帧的当前频点，δθ为来波方向矢量方位角调整量累计值，δφ为仰角调整量累计值。在一个具体实施方式中，根据预设方法更新所述初始方向估计值中的所述预设方法为：θ(n+1)＝θ(n)+uδθφ(n+1)＝φ(n)+uδφ式中，u为正实数或正变量。在一个具体实施方式中，所述预设条件为：判断更新初始方向估计值的次数是否超过预设次数；如果不超过预设次数，将更新后的方向估计值作为当前方向估计值，将调整量累计值清零后，重复步骤s7-步骤s8，得到更新后的当前方向估计值对应的调整量累计值；判断更新后的方向调整量累计值是否不大于所述门限阈值；若不大于所述门限阈值，则停止更新；若大于所述门限阈值，则继续上述更新步骤；如果超过预设次数，停止更新。与现有技术相比，本发明的有益效果：本技术和已有技术相比，本发明涉及一种感应器列阵来波方向的估计方法，本发明提供的估计方法中，doa估计是直接对doa的增量进行迭代计算从而渐进收敛到目标doa值附近，在得到每一次doa估计时运算量大幅降低，其原因是其得到在计算过程中，无需对所有可能的方向进行逐一扫描计算得到空间谱后，再判断其最大值对应的方向，因此在达到相同估计精度的前提下，本发明提供的估计方法相比逐一扫描的技术从计算量上降低了2-3个数量级，运算量减少的同时，可实现实时定位；本发明提供的估计方法采用逼近算法，列阵的拓扑结构规则与否不影响本发明提供的估计方法的计算，所以，本发明提供的估计方法可以用于各种规则及不规则列阵的来波方向估计，提高了实用性。附图说明图1为本发明实施例提供的一种感应器列阵来波方向的估计方法流程图；图2为本发明实施例提供的另一种感应器列阵来波方向的估计方法流程图；图3为本发明实施例提供的预设条件流程图。具体实施方式本申请所涉及的术语解释：doa：估计技术(directionofarrival，doa)music：多重信号分类法(multiplesignalclassification，music)fft：快速傅氏变换(fastfouriertransformation，fft)，是离散傅氏变换的快速算法，它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。mvdr：自适应的空间波数谱估计算法(minimumvariancedistortionlessresponse)该算法具有较快的信干噪比意义下的收敛速度。下面结合具体实施例对本发明做进一步详细的描述，但本发明的实施方式不限于此。实施例1请参见图1，图2。本发明提供了一种感应器列阵来波方向的估计方法，包括步骤：s1.获取待测n列阵阵元，n为正整数，上述n列阵阵元的信号为机械波信号或电磁波信号，机械波信号对应使用水听器列阵，实际应用中为声纳，机械波信号对应使用传声器列阵或超声波感应器列阵，本发明提供的估计方法可适用于任意拓扑形状的感应器列阵，在方向估计结果的精度上高于传统基于空间谱多维扫描方式的精度，并且运算复杂度远远低于传统方法。s2.对所述待测n列阵阵元中的每一列阵阵元接收到的信号数据进行分帧，所述每一列阵阵元均得到m帧时域数据，m为正整数；s3.对各所述每一列阵阵元的m帧时域数据进行转换，得到每一帧时域数据对应的频域数据，并设定初始方向估计值，上述初始方向估计值为[0,0]或根据已知信息得到的粗略方向估计值；对于n阵元的列阵，将每一阵元的接收信号进行长度为n_fft个采样点分帧操作，将分帧后的一帧时域数据通过fft运算变换至频域，本发明采用fft运算，将时域数据转换为频域数据，具体的，fft运算原理为设x(n)为n项的复数序列，由dft变换，任一x(m)的计算都需要n次复数乘法和n-1次复数加法，而一次复数乘法等于四次实数乘法和两次实数加法，一次复数加法等于两次实数加法，即使把一次复数乘法和一次复数加法定义成一次“运算”(四次实数乘法和四次实数加法)，那么求出n项复数序列的x(m),即n点dft变换大约就需要n2次运算。当n＝1024点甚至更多的时候，需要n2＝1048576次运算，在fft中，利用wn的周期性和对称性，把一个n项序列(设n＝2k,k为正整数)，分为两个n/2项的子序列，每个n/2点dft变换需要(n/2)2次运算，再用n次运算把两个n/2点的dft变换组合成一个n点的dft变换。这样变换以后，总的运算次数就变成n2(n/2)2＝nn2/2。继续上面的例子，n＝1024时，总的运算次数就变成了525312次，节省了大约50％的运算量。将上述一分为二进行下去，直到分成两两一组的dft运算单元，那么n点的dft变换就只需要nlog2n次的运算，n在1024点时，运算量仅有10240次，是先前的直接算法的1％，点数越多，运算量的节约就越大，采用fft运算，提高了运算效率。s4.将各列阵阵元中的第一帧均作为当前帧，获取各所述当前帧的频域数据中第一频点的数据，将所述第一频点作为当前频点；具体的，上述第一帧为各个列阵阵元中，处于相同位置的任一帧，上述第一频点为各个列阵阵元中的第一帧中处于相同位置的频点。s5.根据所述当前频点，计算得到所述当前频点数据的频域快拍信号；具体的，首先将阵列接收数据变换到频域，取频域单快拍数据作为压缩感知的测量值，然后根据频域快拍对应的频率、搜索方位和阵列流形构造过完备的阵列流形矩阵作为压缩感知的感知矩阵，本发明提供的估计方法要求的阵元数和快拍数较少，要求的信噪比更低，恢复的目标信号更加准确，波形相关系数达到89％以上。s6.通过所述频域快拍信号，计算得到所述当前频点的协方差矩阵的逆矩阵；将所述初始方向估计值作为当前方向估计值；s7.根据预设方法计算，得到所述当前频点相对于所述当前方向估计值的调整量；在实际应用中，在给定初始方向估计值后，根据列阵几何拓扑，可以计算得到对应频点相对于当前方向估计值的调整量。s8.将各所述当前帧所有频点的目标频点中除所述第一频点外的所有频点均作为当前频点，重复步骤s5-s7，计算所述当前帧中的目标频点的方向估计值调整量，并将各所述调整量求和，得到调整量累计值；s9.判断所述调整量累计值的绝对值是否大于给定的门限阈值的大小；如果大于，根据预设方法更新所述当前方向估计值，直至满足预设条件，停止更新，将当前方向估计值确定为来波方向估计值，所述来波方向估计值的方向矢量方位角与仰角确定为当前帧的信号来波方向；如果不大于，将所述当前方向估计值的方向矢量方位角与仰角确定为当前帧的信号来波方向；s10.将各列阵阵元中的每一帧均作为当前帧，重复步骤s4-s9，确定各个列阵中每一帧的来波方向。具体的，在实际应用中，本技术和已有技术相比，本发明涉及一种感应器列阵来波方向的估计方法，本发明提供的估计方法中，doa估计是直接对doa的增量进行迭代计算从而渐进收敛到目标doa值附近，在得到每一次doa估计时运算量大幅降低，其原因是其得到在计算过程中，无需对所有可能的方向进行逐一扫描计算得到空间谱后，再判断其最大值对应的方向，因此在达到相同估计精度的前提下，本发明提供的估计方法相比逐一扫描的技术从计算量上降低了2-3个数量级，运算量减少的同时，可实现实时定位；本发明提供的估计方法采用逼近算法，列阵的拓扑结构规则与否不影响本发明提供的估计方法的计算，所以，本发明提供的估计方法可以用于各种规则及不规则列阵的来波方向估计，提高了实用性。进一步的，所述目标频点为所述当前帧的所有频点或所述当前帧的部分频点。本发明实施例中提供的目标频点可以是当前帧的所有频点，也可以是当前帧的部分频点，确定当前帧的所有频点或部分频点的方向估计值，均可以采用本发明实施例提供的方法进行估计，进一步的，所述待测n列阵阵元为任意几何形状的列阵阵元，其阵元位置表达式为：p＝{pi}；pi＝[xi，yi，zi]式中，p为阵元位置，i为阵元序号，i＝1-n，pi为阵元位置坐标。进一步的，通过快速傅氏变换运算对各所述m帧时域数据进行转换，得到每一帧时域数据对应的频域数据。进一步的，所述当前频点的数据的频域快拍信号的表达式为：[x(wl，1)，x(wl，2)，…，x(wl，n)]式中，x(wl，1)为当前帧第一列阵阵元的当前帧的当前频点的频域数据，x(wl，2)为当前帧第二列阵阵元的当前帧的当前频点的频域数据。x(wl，n)为当前帧第n列阵阵元的当前帧的当前频点的频域数据，wl为当前帧的当前频点。进一步的，所述初始方向估计值的表达式为：[θ(0)，φ(0)]＝[θinit，φinit]式中，θ为来波方向矢量方位角，ф为来波仰角。进一步的，所述预设方法为：其中，ha(pij，θ，φ)＝(-xijsinθcosφ+yijcosθcosφ)he(pij，θ，φ)＝(-xijcosθsinφ+yijsinθsinφ+zijcosφ)h(φij，pij，wl)＝sin(φij+wlτij)式中pi＝[xi，yi，zi]为阵元i在直角坐标系中的位置坐标；pij＝[xij,yij,zij]为阵元i到阵元j的相对位置坐标，xij＝xi-xj；yij＝yi-yj；zij＝zi-zj；aij，фij分别为信号协方差矩阵的逆矩阵元素的幅值和相角，τi,j＝τi,0-τj,0为列阵阵元i和列阵阵元j之间的相对延时，由列阵阵元的位置坐标及信号来波方向确定，c为信号传播速度。上述公式的等价变换均为本申请的保护范围。进一步的，累计所述目标频点分别对应的所述调整量，得到调整量累计值的计算公式为：δθ＝∑lδθwlδφ＝∑lδφwl式中，wl为当前帧的当前频点，δθ为来波方向矢量方位角的累计值，δφ为仰角累计值。进一步的，根据预设方法更新所述调整量累计值中的所述预设方法为：θ(n+1)＝θ(n)+uδθφ(n+1)＝φ(n)+uδφ式中，u为正实数或正变量。进一步的，请参见图3，所述预设条件为：判断更新初始方向估计值的次数是否超过预设次数；如果不超过预设次数，将更新后的方向估计值作为当前方向估计值，将调整量累计值清零后，重复步骤s7-步骤s8，得到更新后的当前方向估计值对应的调整量累计值；判断更新后的方向调整量累计值是否不大于所述门限阈值；若不大于所述门限阈值，则停止更新；若大于所述门限阈值，则继续上述更新步骤；如果超过预设次数，停止更新。具体的，本发明提供的方法在确定每一帧方向时采用双重约束条件，在有限次数下，对调整量累计值进行叠加，如果大于预设次数，则停止迭代计算，将当前方向估计值确定为本发明来波方向估计值，如果在有限次数内计算得到的调整量累计值小于门限阈值，同样停止更新，保证了本发明提供的估计方法正常运行。实施例2本发明实施例中，将mvdr空间谱进行逐一方向搜索来波方向估计与本发明提供的来波方向估计方法进行对比：对比条件：假设声源方向矢量方位角和仰角的范围都是[-90°+90°]，扫描间隔为1°，阵元数目n，计算范围为一个频点，本技术计算方法迭代次数为l次。对比结果如下：mvdr本发明提供的估计方法需要存储列阵流行矢量数目32761*n*32bit0bit需要存储列空间谱扫描结果的缓存32761*32bit0bit估计过程用到的复数乘法器数量32761*n*(n+1)0估计过程用到的复数加法器数量32761*n*(n+1)0估计过程用到的除法器数量327610估计过程用到的实数乘法器数量0l*n*(n+1)*5估计过程用到的实数加法器数量0l*n*(n+1)*2.5需要用到的正弦查找表0个1个通过上表我们可以看出，当本技术的迭代次数l为正常值32时，本发明提供的估计方法进行一次来波方向估计，所需乘法运算次数和加法运算次数比空间谱逐一扫描法降低了3个数量级；同时需要缓存的临时数据所占用的空间为零。本技术和已有技术相比，本发明涉及一种感应器列阵来波方向的估计方法，本发明提供的估计方法中，doa估计是直接对doa的增量进行迭代计算从而渐进收敛到目标doa值附近，在得到每一次doa估计时运算量大幅降低，其原因是其得到在计算过程中，无需对所有可能的方向进行逐一扫描计算得到空间谱后，再判断其最大值对应的方向，因此在达到相同估计精度的前提下，本发明提供的估计方法相比逐一扫描的技术从计算量上降低了2-3个数量级，运算量减少的同时，可实现实时定位；本发明提供的估计方法采用逼近算法，列阵的拓扑结构规则与否不影响本发明提供的估计方法的计算，所以，本发明提供的估计方法可以用于各种规则及不规则列阵的来波方向估计，提高了实用性。本发明提供的估计方法具有通用性。首先本技术可适用于任意拓扑结构的感应器列阵进行doa估计，感应器列阵为均匀矩形列阵、均匀环形列阵、均匀线列阵、螺旋线列阵、拓扑结构为任意无规律的列阵中的任一种。本发明提供的估计方法的通用性能够给与感应器列阵拓扑结构设计者更大的列阵设计自由度，进而使设计者可以设计出对于某种应用更加优化的列阵拓扑结构，以满足灵活和高性能的设计要求。本发明提供的估计方法对于各种类型的波也具有通用性。以机械波、声波、超声波、电磁波中的任一种为对象的doa估计，均可以采用本发明提供的估计方法进行估计。对于属于电磁波类型的雷达列阵信号，可视为窄带信号，用本发明提供的技术进行doa估计时，可将其视为单一子带上的doa估计；对于机械波类型的麦克风列阵、水听器列阵等，其信号为宽带信号的类型，本技术亦包含了将该宽带信号分解为子带信号，并进行doa估计的方法，从适用性上更加广泛。以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属
技术领域：
的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单推演或替换，都应当视为属于本发明的保护范围。当前第1页12