CN115278496A - 用于麦克风阵列测量的稀疏声源辨识方法、系统 - Google Patents

用于麦克风阵列测量的稀疏声源辨识方法、系统 Download PDF

Info

Publication number
CN115278496A
CN115278496A CN202210223139.7A CN202210223139A CN115278496A CN 115278496 A CN115278496 A CN 115278496A CN 202210223139 A CN202210223139 A CN 202210223139A CN 115278496 A CN115278496 A CN 115278496A
Authority
CN
China
Prior art keywords
sound source
sound
matrix
vector
intensity
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202210223139.7A
Other languages
English (en)
Inventor
孙仕林
王天杨
褚福磊
谭建鑫
井延伟
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Hebei Jiantou New Energy Co ltd
Tsinghua University
Original Assignee
Hebei Jiantou New Energy Co ltd
Tsinghua University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Hebei Jiantou New Energy Co ltd, Tsinghua University filed Critical Hebei Jiantou New Energy Co ltd
Priority to CN202210223139.7A priority Critical patent/CN115278496A/zh
Publication of CN115278496A publication Critical patent/CN115278496A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • H—ELECTRICITY
    • H04—ELECTRIC COMMUNICATION TECHNIQUE
    • H04R—LOUDSPEAKERS, MICROPHONES, GRAMOPHONE PICK-UPS OR LIKE ACOUSTIC ELECTROMECHANICAL TRANSDUCERS; DEAF-AID SETS; PUBLIC ADDRESS SYSTEMS
    • H04R29/00—Monitoring arrangements; Testing arrangements
    • H04R29/004—Monitoring arrangements; Testing arrangements for microphones
    • H04R29/005—Microphone arrays
    • G—PHYSICS
    • G10—MUSICAL INSTRUMENTS; ACOUSTICS
    • G10L—SPEECH ANALYSIS TECHNIQUES OR SPEECH SYNTHESIS; SPEECH RECOGNITION; SPEECH OR VOICE PROCESSING TECHNIQUES; SPEECH OR AUDIO CODING OR DECODING
    • G10L25/00—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00
    • G10L25/27—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00 characterised by the analysis technique
    • G—PHYSICS
    • G10—MUSICAL INSTRUMENTS; ACOUSTICS
    • G10L—SPEECH ANALYSIS TECHNIQUES OR SPEECH SYNTHESIS; SPEECH RECOGNITION; SPEECH OR VOICE PROCESSING TECHNIQUES; SPEECH OR AUDIO CODING OR DECODING
    • G10L25/00—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00
    • G10L25/48—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00 specially adapted for particular use
    • G10L25/51—Speech or voice analysis techniques not restricted to a single one of groups G10L15/00 - G10L21/00 specially adapted for particular use for comparison or discrimination

Landscapes

  • Engineering & Computer Science (AREA)
  • Signal Processing (AREA)
  • Health & Medical Sciences (AREA)
  • Physics & Mathematics (AREA)
  • Acoustics & Sound (AREA)
  • Computational Linguistics (AREA)
  • Audiology, Speech & Language Pathology (AREA)
  • Human Computer Interaction (AREA)
  • Multimedia (AREA)
  • General Health & Medical Sciences (AREA)
  • Otolaryngology (AREA)
  • Circuit For Audible Band Transducer (AREA)

Abstract

本发明提供一种用于麦克风阵列测量的稀疏声源辨识方法，首先利用离散傅里叶变换将信号转换到频域，构成频域声压矩阵，并在频域声压矩阵中选择信号成分，以构成声压向量，再根据预设规则选择声源重构区域，将声源重构区域划分离散网格，根据麦克风阵列中阵元的位置坐标和声压向量建立声传播模型，而后根据声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对声学逆问题进行变形以获取变形表达式，采用交替分裂的方式对变形表达式进行求解，以获取待求的目标声源强度，并根据目标声源强度形成声像，完成声源辨识；如此，便于对声学逆问题进行求解，提高计算效率；提高辨识精度。

Description

用于麦克风阵列测量的稀疏声源辨识方法、系统
技术领域
本发明涉及阵列信号处理技术领域，更为具体地，涉及一种用于麦克风阵列测量的稀疏声源辨识方法、系统。
背景技术
麦克风阵列由多个麦克风按照特定的空间分布构成，相比于利用独立麦克风的测量方式，在阵列测量过程中多个通道同步采集信号，由于阵元位置的相异性，不同麦克风采集到的声信号存在一定的相位延迟，因此不同通道之间的差异可以反映声源的空域信息。人耳能够通过听到声音感知到声源的存在，但仅靠听觉很难对声源的确切位置和强度进行准确估计，因此声源的定量辨识需要通过其它手段的辅助来完成。通过阵列信号处理技术，麦克风阵列能够获得声源的时域、空域等多维度特征，将可听的声音转化成可见的图像，从而实现声源辨识。基于麦克风阵列的声源辨识可以重建特定区域内声源的分布情况，识别异响声源的存在和性质，在汽车工业、航空航天、机械健康监测等领域具有广阔的应用前景。
稀疏声源广泛存在于现代生产、生活场景中，如汽车异响声源大多出现在中控台、车窗边缘等局部区域，由机械故障引起的冲击噪声源仅存在于运动副附近。这些声源在空域中的分布具有稀疏性，即声源辨识结果中仅存在少量的非零元素，如果直接利用常规波束形成(Williams J R.Fast Beam-Forming Algorithm[J].The Journal of theAcoustical Society of America,1968,44(5):1454-1455.)、自适应波束形成(CN201610409905、CN201911406837)等方法开展声源辨识，所得结果的分辨率和准确性较低，不能完全满足实际应用中稀疏声源辨识的需要。在稀疏声源辨识过程中，空域稀疏性可以作为先验知识提供辨识所需的约束条件，进而提高声源辨识的准确性。目前，稀疏声源辨识的基本思路为：根据声源重构区域和阵列布局建立声传播模型，并在稀疏先验知识的约束下对基于声传播模型的逆问题进行求解，从而重构声源位置和声源强度等信息。声源辨识中的稀疏约束通过声源强度向量的0范数表示，0范数的非平滑和非凸性质为声学逆问题的求解带来了困难，因此无法对声学逆问题进行直接求解。根据声学逆问题近似求解方法的不同，稀疏声源辨识方法可以分为两类：基于正交匹配追踪的方法和基于基追踪的方法。基于正交匹配追踪的方法是一种贪婪求解方法，通过逐步迭代计算的方式对声源进行重构，在每一步迭代中以局部最优解代替全局最优解，试图通过局部最优解的累积更新逼近全局最优解，因此具有较高的计算效率，求解耗时少，但由于求解策略的局限性，局部最优解往往不能精确地逼近全局最优解，这一问题导致算法收敛时所得声源强度向量的准确性难以保证，削弱了稀疏声源辨识结果的可靠性。基于基追踪的方法将原始的声学逆问题中的0范数用近似函数代替，近似函数在声学逆问题求解过程中能够近似0范数的稀疏约束特性，对解向量中的非零元素施加稀疏惩罚，从而得到尽可能稀疏的声源辨识结果；另外，近似函数相比于0范数具有更好的求解性质，比如具有凸性、平滑性等，有助于稀疏恢复过程中最优化问题的顺利求解，因此这一类方法能够得到相比于正交匹配追踪方法更准确的结果。常见的近似函数有1范数、近似0范数等。然而，由于近似函数与0范数在基本性质上存在差异，基于基追踪的稀疏声源辨识方法在结果准确性方面仍存在不足。非零元素的0范数均为同一数值，表明不同幅值的非零元素受到0范数带来的稀疏惩罚是一致的，但在上述的近似函数正则化的声学逆问题求解过程中，解向量中幅值不同的非零元素受到的稀疏惩罚程度不同。这一差异削弱了现有的基于基追踪的方法在声源辨识中的可靠性，使得声源强度存在欠估计的问题，当求解结果的动态范围较小时，强度被严重欠估计的声源将被错误地忽略，从而导致现有方法不能有效地获取声源位置和声源强度等信息。另外，在基于基追踪的稀疏声源辨识方法中，声学逆问题的求解大多通过迭代的方式进行，需要考虑结果准确性与计算效率之间的平衡，所以对于特定的近似函数，声学逆问题求解算法的构建成为稀疏声源辨识的难点之一。
因此，亟需一种提高0范数的近似效果，增强非零元素在稀疏恢复过程中所受惩罚程度的一致性，提高稀疏声源辨识结果的准确性，提高计算效率的用于麦克风阵列测量的稀疏声源辨识方法、系统。
发明内容
鉴于上述问题，本发明的目的是提供一种用于麦克风阵列测量的稀疏声源辨识方法，以解决现有技术中声源强度存在欠估计，当求解结果的动态范围较小时，强度被严重欠估计的声源将被错误地忽略，从而导致现有方法不能有效地获取声源位置和声源强度等信息的问题。
本发明提供的一种用于麦克风阵列测量的稀疏声源辨识方法，包括：
利用包含M个通道的麦克风阵列采集长度为La的声信号pIm(m＝1,2,...,M)，以构建大小为M×La的原始声压矩阵
Figure BDA0003534509740000031
对各个通道所采集的声信号pIm进行分割处理以形成长度为L、交叠率为R的信号段，并根据所述原始声压矩阵获取包含有
Figure BDA0003534509740000032
个快拍的声压矩阵
Figure BDA0003534509740000033
其中，
Figure BDA0003534509740000034
表示不大于
Figure BDA0003534509740000035
的整数，第s(s＝1,2,...,S)个快拍对应的子矩阵Ps II＝[ps,1,ps,2,...,ps,M]T大小为M×L，
Figure BDA0003534509740000036
表示
Figure BDA0003534509740000037
中的第L(1-R)(s-1)+1到第L(1-R)(s-1)+L个元素；
利用离散傅里叶变换将第s个快拍对应的信号转换到频域，构成频域声压矩阵
Figure BDA0003534509740000038
并在所述频域声压矩阵中选择第l(l＝1,2,...,L)个频率fl所对应的信号成分，以构成声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T，其中，
Figure BDA0003534509740000039
DFT表示离散傅里叶变换；pm(fl)代表
Figure BDA00035345097400000310
中频率fl对应的信号成分；
根据预设规则选择声源重构区域，将所述声源重构区域划分为N＝N1×N2个离散网格，所述离散网格中网格点的位置坐标为un＝(unx,uny,unz)(n＝1,2,...,N)；根据所述麦克风阵列中阵元的位置坐标和所述声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T建立声传播模型p＝Gq+r，其中，q为待求解的声源强度向量，r为噪声向量，G为大小为M×N的传输矩阵，且所述传输矩阵的元素为
Figure BDA0003534509740000041
其中exp()是指数函数，c是声速，
Figure BDA0003534509740000042
是虚数单位；
根据所述声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对所述声学逆问题进行变形以获取变形表达式，采用交替分裂的方式对所述变形表达式进行求解，以获取待求的目标声源强度，并根据所述目标声源强度形成声像，完成声源辨识。
优选地，在根据所述声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对所述声学逆问题进行变形以获取变形表达式的过程中，包括：
根据所述声传播模型建立基于广义最小最大凹函数的声学逆问题
Figure BDA0003534509740000043
其中，λ＞0为需要给定的正则化参数，
Figure BDA0003534509740000044
为广义最小最大凹函数，B为辅助矩阵；v为辅助向量，且为变量；
给定平衡参数0＜γ＜1，令所述辅助矩阵
Figure BDA0003534509740000045
以对所述声学逆问题进行变形确定所述声学逆问题的变形表达式：
Figure BDA0003534509740000046
将所述变形表达式作为目标函数，并将所述目标函数对于变量q的平滑项记为f1(q)，将所述目标函数对于变量q的非平滑项记为f2(q)，将所述目标函数对于变量v的平滑项记为g1(v)，将所述目标函数对于变量v的非平滑项记为g2(v)。
优选地，采用交替分裂的方式对所述变形表达式进行求解，以获取待求的目标声源强度，并根据所述目标声源强度形成声像的过程，包括：
初始化计数变量，并随机初始化中间变量
Figure BDA0003534509740000051
和
Figure BDA0003534509740000052
随着计数变量的更新根据所述变形表达式依次迭代更新计算待求的声源强度向量的值，并计算两个相邻次数所求得的声源强度向量的值的强度误差；
判断所述强度误差是否超过预设的收敛阈值；若所述强度误差超过预设的收敛阈值，则继续迭代更新计算待求的声源强度向量的值；若所述强度误差不超过所述收敛阈值，则将最后一次计算的声源强度向量的值作为待求的目标声源强度；
根据所述目标声源强度生成声像。
优选地，所述收敛阈值为0.0001。
优选地，所述交替分裂的方式为Douglas-Rachford分裂的方式。
优选地，随着计数变量的更新根据所述变形表达式依次迭代更新计算待求的声源强度向量的值，并计算两个相邻次数所求得的声源强度向量的值的强度误差的过程，包括：
计算当前计数变量下当前的声源中间变量的值和当前的声源强度向量的值；其中，当前的声源中间变量为：
Figure BDA0003534509740000053
当前的声源强度向量为：
Figure BDA0003534509740000054
令所述当前的声源强度向量中k＝k+1，计算所述当前计数变量的下一变量的源强度向量的值，并计算所述下一变量的源强度向量的值与所述当前的声源强度向量的值的强度误差eR。
优选地，所述强度误差为
Figure BDA0003534509740000055
优选地，在计算当前计数变量下当前的声源中间变量的值和当前的声源强度向量的值的过程中，
还计算当前的辅助向量的值和当前的辅助中间变量的值；其中，
当前的辅助向量的为：
Figure BDA0003534509740000061
当前的辅助中间变量为：
Figure BDA0003534509740000062
其中，
Figure BDA0003534509740000063
为预算的近端算子。
优选地，若所述强度误差不超过所述收敛阈值，则将最后一次计算的声源强度向量的值作为待求的目标声源强度；根据所述目标声源强度生成声像的过程，包括：
将最后一次计算的声源强度向量的值作为待求的目标声源强度的值；
将所述目标声源强度转换为大小为N1×N2的声源强度矩阵
Figure BDA0003534509740000064
并对所述声源强度矩阵
Figure BDA0003534509740000065
中的各个元素进行对数变换以获取目标声源矩阵
Figure BDA0003534509740000066
其中，
Figure BDA0003534509740000067
为矩阵
Figure BDA0003534509740000068
中的元素，Qref为预设的参考声源强度；
根据所述目标声源矩阵生成声像，根据所述声像进行声源辨识。
另一方面，本发明还提供一种用于麦克风阵列测量的稀疏声源辨识系统，实现如前所述的用于麦克风阵列测量的稀疏声源辨识方法，包括：
声压矩阵获取单元，用于利用包含M个通道的麦克风阵列采集长度为La的声信号pIm(m＝1,2,...,M)，以构建大小为M×La的原始声压矩阵
Figure BDA0003534509740000069
对各个通道所采集的声信号pIm进行分割处理以形成长度为L、交叠率为R的信号段，并根据所述原始声压矩阵获取包含有
Figure BDA00035345097400000610
个快拍的声压矩阵
Figure BDA00035345097400000611
其中，
Figure BDA00035345097400000612
表示不大于
Figure BDA00035345097400000613
的整数，第s(s＝1,2,...,S)个快拍对应的子矩阵
Figure BDA00035345097400000614
大小为M×L，
Figure BDA00035345097400000615
表示
Figure BDA00035345097400000616
中的第L(1-R)(s-1)+1到第L(1-R)(s-1)+L个元素；
频域转换单元，用于利用离散傅里叶变换将第s个快拍对应的信号转换到频域，构成频域声压矩阵
Figure BDA0003534509740000071
并在所述频域声压矩阵中选择第l(l＝1,2,...,L)个频率fl所对应的信号成分，以构成声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T，其中，
Figure BDA0003534509740000072
DFT表示离散傅里叶变换；pm(fl)代表
Figure BDA0003534509740000073
中频率fl对应的信号成分；
模型建立单元，用于根据预设规则选择声源重构区域，将所述声源重构区域划分为N＝N1×N2个离散网格，所述离散网格中网格点的位置坐标为un＝(unx,uny,unz)(n＝1,2,...,N)；根据所述麦克风阵列中阵元的位置坐标和所述声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T建立声传播模型p＝Gq+r，其中，q为待求解的声源强度向量，r为噪声向量，G为大小为M×N的传输矩阵，且所述传输矩阵的元素为
Figure BDA0003534509740000074
其中exp()是指数函数，c是声速，
Figure BDA0003534509740000075
是虚数单位；
辨识求解单元，用于根据所述声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对所述声学逆问题进行变形以获取变形表达式，采用交替分裂的方式对所述变形表达式进行求解，以获取待求的目标声源强度，并根据所述目标声源强度形成声像，完成声源辨识。
从上面的技术方案可知，本发明提供的用于麦克风阵列测量的稀疏声源辨识方法，首先利用麦克风阵列采集声信号，以构建大小为M×La的原始声压矩阵；对各个通道所采集的声信号进行分割处理以形成信号段，并根据原始声压矩阵获取声压矩阵，再利用离散傅里叶变换将信号转换到频域，构成频域声压矩阵，并在频域声压矩阵中选择信号成分，以构成声压向量，再根据预设规则选择声源重构区域，将声源重构区域划分离散网格，根据麦克风阵列中阵元的位置坐标和声压向量建立声传播模型，而后根据声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对声学逆问题进行变形以获取变形表达式，采用交替分裂的方式对变形表达式进行求解，以获取待求的目标声源强度，并根据目标声源强度形成声像，完成声源辨识；如此，利用广义最小最大凹函数作为正则化函数建立稀疏声源辨识方法，采用的近似函数对0范数的近似效果更好，非零元素在稀疏恢复过程中受到的惩罚程度的一致性更强，能够解决现有技术对声源强度欠估计的问题；并且声学逆问题的目标函数具有凸性，能够保证稀疏声源辨识结果的准确性，引入辅助变量后，采用交替分裂的方式对声学逆问题进行求解，建立了显式的变量迭代更新公式，能够方便地对声学逆问题进行求解，提高计算效率；而且考虑信号交叠分割策略和不同的声源频率以适应外部环境的变化，求解过程采用的稀疏先验信息能够抑制测量噪声带来的干扰，提高辨识精度。
附图说明
通过参考以下结合附图的说明书内容，并且随着对本发明的更全面理解，本发明的其它目的及结果将更加明白及易于理解。在附图中：
图1为根据本发明实施例的用于麦克风阵列测量的稀疏声源辨识方法的流程图；
图2为根据本发明实施例的用于麦克风阵列测量的稀疏声源辨识方法所获取的声源辨识结果；
图3为本发明实施例的用于麦克风阵列测量的稀疏声源辨识方法与现有技术所得到的声源辨识结果的均方根误差；
图4为根据本发明实施例的用于麦克风阵列测量的稀疏声源辨识系统的框图。
具体实施方式
非零元素的0范数均为同一数值，表明不同幅值的非零元素受到0范数带来的稀疏惩罚是一致的，但在现有的近似函数正则化的声学逆问题求解过程中，解向量中幅值不同的非零元素受到的稀疏惩罚程度不同。这一差异削弱了现有的基于基追踪的方法在声源辨识中的可靠性，使得声源强度存在欠估计的问题，当求解结果的动态范围较小时，强度被严重欠估计的声源将被错误地忽略，从而导致现有方法不能有效地获取声源位置和声源强度等信息。另外，在基于基追踪的稀疏声源辨识方法中，声学逆问题的求解大多通过迭代的方式进行，需要考虑结果准确性与计算效率之间的平衡，所以对于特定的近似函数，声学逆问题求解算法的构建成为稀疏声源辨识的难点之一。
针对上述问题，本发明提供一种用于麦克风阵列测量的稀疏声源辨识方法，以下将结合附图对本发明的具体实施例进行详细描述。
为了说明本发明提供的用于麦克风阵列测量的稀疏声源辨识方法、系统，图1-图3对本发明实施例的用于麦克风阵列测量的稀疏声源辨识方法进行了示例性标示；图4对本发明实施例的用于麦克风阵列测量的稀疏声源辨识系统进行了示例性标示。
以下示例性实施例的描述实际上仅仅是说明性的，决不作为对本发明及其应用或使用的任何限制。对于相关领域普通技术人员已知的技术和设备可能不作详细讨论，但在适当情况下，所述技术和设备应当被视为说明书的一部分。
如图1、图2、图3共同所示，本发明提供的用于麦克风阵列测量的稀疏声源辨识方法，包括：
利用包含M个通道的麦克风阵列采集长度为La的声信号pIm(m＝1,2,...,M)，以构建大小为M×La的原始声压矩阵
Figure BDA0003534509740000091
对各个通道所采集的声信号pIm进行分割处理以形成长度为L、交叠率为R的信号段，并根据所述原始声压矩阵获取包含有
Figure BDA0003534509740000092
个快拍的声压矩阵
Figure BDA0003534509740000093
其中，
Figure BDA0003534509740000094
表示不大于
Figure BDA0003534509740000095
的整数，第s(s＝1,2,...,S)个快拍对应的子矩阵Ps II＝[ps,1,ps,2,...,ps,M]T大小为M×L，
Figure BDA0003534509740000096
表示
Figure BDA0003534509740000097
中的第L(1-R)(s-1)+1到第L(1-R)(s-1)+L个元素；
利用离散傅里叶变换将第s个快拍对应的信号转换到频域，构成频域声压矩阵
Figure BDA0003534509740000098
并在所述频域声压矩阵中选择第l(l＝1,2,...,L)个频率fl所对应的信号成分，以构成声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T，其中，
Figure BDA0003534509740000101
DFT表示离散傅里叶变换；pm(fl)代表
Figure BDA0003534509740000102
中频率fl对应的信号成分；
根据预设规则选择声源重构区域，将所述声源重构区域划分为N＝N1×N2个离散网格，所述离散网格中网格点的位置坐标为un＝(unx,uny,unz)(n＝1,2,...,N)；根据所述麦克风阵列中阵元的位置坐标和所述声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T建立声传播模型p＝Gq+r，其中，q为待求解的声源强度向量，r为噪声向量，G为大小为M×N的传输矩阵，且所述传输矩阵的元素为
Figure BDA0003534509740000103
其中exp()是指数函数，c是声速，
Figure BDA0003534509740000104
是虚数单位；
根据所述声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对所述声学逆问题进行变形以获取变形表达式，采用交替分裂的方式对所述变形表达式进行求解，以获取待求的目标声源强度，并根据所述目标声源强度形成声像，完成声源辨识。
在图1、图2、图3、共同所示的实施例中，步骤S1为利用包含M个通道的麦克风阵列采集长度为La的声信号pIm(m＝1,2,...,M)，以构建大小为M×La的原始声压矩阵
Figure BDA0003534509740000105
对各个通道所采集的声信号pIm进行分割处理以形成长度为L、交叠率为R的信号段，并根据所述原始声压矩阵获取包含有
Figure BDA0003534509740000106
个快拍的声压矩阵
Figure BDA0003534509740000107
的过程。
在一个具体实施例中，首先选择通道数量M＝50，信号采集长度La＝96000，阵列拓扑形式为平面均匀圆阵，阵列直径为1m，阵列位于三维直角坐标系的x-y平面内，阵列平面垂直于z轴的正方向，阵列中心位于坐标原点，而后将各通道的信号pIm分割为长度为L、交叠率为R的信号段，由原始声压矩阵得到包含有
Figure BDA0003534509740000111
个快拍的声压矩阵
Figure BDA0003534509740000112
其中
Figure BDA0003534509740000113
表示不大于
Figure BDA0003534509740000114
的整数，第s(s＝1,2,...,S)个快拍对应的子矩阵
Figure BDA0003534509740000115
大小为M×L，
Figure BDA0003534509740000116
表示
Figure BDA0003534509740000117
中的第L(1-R)(s-1)+1到第L(1-R)(s-1)+L个元素；并且在具体实施过程中，选择信号分割长度L＝96000，交叠率R＝0.5。
在图1、图2、图3、共同所示的实施例中，步骤S2为利用离散傅里叶变换将第s个快拍对应的信号转换到频域，构成频域声压矩阵
Figure BDA0003534509740000118
并在所述频域声压矩阵中选择第l(l＝1,2,...,L)个频率fl所对应的信号成分，以构成声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T的过程。
在本实施例中，对于第s个快拍对应的信号，利用离散傅里叶变换将信号转换到频域，构成频域声压矩阵
Figure BDA0003534509740000119
其中
Figure BDA00035345097400001110
DFT表示离散傅里叶变换；选择频域声压矩阵中的第l(l＝1,2,...,L)个频率fl对应的信号成分，构成声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T，其中pm(fl)代表
Figure BDA00035345097400001111
中频率fl对应的元素；在具体实施过程中，选择分析第5个快拍对应的信号，选择分析频率fl＝4000Hz对应的信号成分。
在图1、图2、图3、共同所示的实施例中，步骤S3为根据预设规则选择声源重构区域，将所述声源重构区域划分为N＝N1×N2个离散网格，所述离散网格中网格点的位置坐标为un＝(unx,uny,unz)(n＝1,2,...,N)；根据所述麦克风阵列中阵元的位置坐标和所述声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T建立声传播模型p＝Gq+r的过程。
在本实施例中，根据预设的规则选择声源重构区域，将声源重构区域划分为N＝N1×N2个离散网格，网格点的位置坐标为un＝(unx,uny,unz)(n＝1,2,...,N)，根据麦克风阵列中阵元的位置坐标um＝(umx,umy,umz)建立声传播模型p＝Gq+r，其中q是需要求解的声源强度向量，r是噪声向量，G是大小为M×N的传输矩阵，其中的元素为
Figure BDA0003534509740000121
exp()是指数函数，c是声速，
Figure BDA0003534509740000122
是虚数单位。
在具体实施的过程中，声源位于(-0.1m,0.15m,1m)、(0.2m,0.3m,1m)，声源强度均为100dB，声源强度定义为距离声源1m处的声压级，参考值为2×10-5Pa，重构区域为平行于x-y平面的1m×1m的正方形区域，重构区域与麦克风阵列平面之间的距离为1m，重构区域的中心位于(0m,0m,1m)，选择N1＝41，N2＝41，将声源重构区域均匀划分成离散网格，声速c＝340m/s。
在图1、图2、图3、共同所示的实施例中，步骤S4为根据声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对所述声学逆问题进行变形以获取变形表达式，采用交替分裂的方式对所述变形表达式进行求解，以获取待求的目标声源强度，并根据所述目标声源强度形成声像，完成声源辨识的过程；其中，
在根据所述声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对所述声学逆问题进行变形以获取变形表达式的过程中，包括：
S411：根据所述声传播模型建立基于广义最小最大凹函数的声学逆问题
Figure BDA0003534509740000123
其中，λ＞0为需要给定的正则化参数，
Figure BDA0003534509740000124
为广义最小最大凹函数，B为辅助矩阵；v为辅助向量，且为变量；
S412：给定平衡参数0＜γ＜1，令所述辅助矩阵
Figure BDA0003534509740000125
以对所述声学逆问题进行变形确定所述声学逆问题的变形表达式：
Figure BDA0003534509740000131
S413：将所述变形表达式作为目标函数，并将所述目标函数对于变量q的平滑项记为f1(q)，将所述目标函数对于变量q的非平滑项记为f2(q)，将所述目标函数对于变量v的平滑项记为g1(v)，将所述目标函数对于变量v的非平滑项记为g2(v)。
采用交替分裂的方式对所述变形表达式进行求解，以获取待求的目标声源强度，并根据所述目标声源强度形成声像的过程，包括：
S421：初始化计数变量，并随机初始化中间变量
Figure BDA0003534509740000132
和
Figure BDA0003534509740000133
S422：随着计数变量的更新根据所述变形表达式依次迭代更新计算待求的声源强度向量的值，并计算两个相邻次数所求得的声源强度向量的值的强度误差；
S423：判断所述强度误差是否超过预设的收敛阈值；若所述强度误差超过预设的收敛阈值，则继续迭代更新计算待求的声源强度向量的值；若所述强度误差不超过所述收敛阈值，则将最后一次计算的声源强度向量的值作为待求的目标声源强度；
S424：根据所述目标声源强度生成声像。
在一个具体实施例中，该收敛阈值为0.0001；该交替分裂的方式为Douglas-Rachford分裂的方式。
在本具体实施例中，随着计数变量的更新根据所述变形表达式依次迭代更新计算待求的声源强度向量的值，并计算两个相邻次数所求得的声源强度向量的值的强度误差的过程，包括：
S4221：计算当前计数变量下当前的声源中间变量的值和当前的声源强度向量的值；其中，当前的声源中间变量为：
Figure BDA0003534509740000134
当前的声源强度向量为：
Figure BDA0003534509740000135
S4222：令所述当前的声源强度向量中k＝k+1，计算所述当前计数变量的下一变量的源强度向量的值，并计算所述下一变量的源强度向量的值与所述当前的声源强度向量的值的强度误差eR；其中，该强度误差为
Figure BDA0003534509740000141
在计算当前计数变量下当前的声源中间变量的值和当前的声源强度向量的值的过程中，
还计算当前的辅助向量的值和当前的辅助中间变量的值；其中，
当前的辅助向量的为：
Figure BDA0003534509740000142
当前的辅助中间变量为：
Figure BDA0003534509740000143
其中，
Figure BDA0003534509740000144
为预算的近端算子。
若所述强度误差不超过所述收敛阈值，则将最后一次计算的声源强度向量的值作为待求的目标声源强度；根据所述目标声源强度生成声像的过程，包括：
S4231：将最后一次计算的声源强度向量的值作为待求的目标声源强度的值；
S4232：将所述目标声源强度转换为大小为N1×N2的声源强度矩阵
Figure BDA0003534509740000145
并对所述声源强度矩阵
Figure BDA0003534509740000146
中的各个元素进行对数变换以获取目标声源矩阵
Figure BDA0003534509740000147
其中，
Figure BDA0003534509740000148
为矩阵
Figure BDA0003534509740000149
中的元素，Qref为预设的参考声源强度；
S4233：根据所述目标声源矩阵生成声像，根据所述声像进行声源辨识。
更为具体的，在一个更为具体的实施例中，步骤S4具体步骤如下：
(1)首先根据声传播模型建立基于广义最小最大凹函数的声学逆问题
Figure BDA00035345097400001410
其中λ＞0是需要给定的正则化参数，
Figure BDA00035345097400001411
是广义最小最大凹函数，B是辅助矩阵；在具体实施过程中，选择正则化参数λ＝0.1。
(2)：给定平衡参数0＜γ＜1，令辅助矩阵
Figure BDA0003534509740000151
得到声学逆问题的进一步表达式：
Figure BDA0003534509740000152
上述目标函数对于变量q的平滑项f1(q)，记上述目标函数对于变量q的非平滑项f2(q)，记上述目标函数对于变量v的平滑项g1(v)，记上述目标函数对于变量v的非平滑项g2(v)；在具体实施过程中，选择平衡参数γ＝0.7，目标函数对于变量q的平滑项
Figure BDA0003534509740000153
目标函数对于变量q的非平滑项f2(q)＝λ||q1，目标函数对于变量v的平滑项
Figure BDA0003534509740000154
目标函数对于变量v的非平滑项g2(v)＝-λv||1。
(3)：初始化计数变量k＝0，随机初始化中间变量
Figure BDA0003534509740000155
和
Figure BDA0003534509740000156
在具体实施过程中，中间变量
Figure BDA0003534509740000157
和
Figure BDA0003534509740000158
均按照均值为0，方差为1的高斯分布进行随机初始化，即其中所有的元素随机产生。
(4)：采用Douglas-Rachford分裂的方式，计算声源强度向量的更新值
Figure BDA0003534509740000159
计算中间变量的更新值
Figure BDA00035345097400001510
计算辅助向量的更新值
Figure BDA00035345097400001511
计算中间变量的更新值
Figure BDA00035345097400001512
其中，近端算子为
Figure BDA00035345097400001513
在具体实施过程中，声源强度向量的更新值
Figure BDA00035345097400001514
中间变量的更新值
Figure BDA0003534509740000161
辅助向量的更新值
Figure BDA0003534509740000162
中间变量的更新值
Figure BDA0003534509740000163
I是与GTG大小相同的单位矩阵。
(5)：令k＝k+1，计算声源强度向量在更新前后的相对误差eR，如果eR超过预先设定的收敛阈值eT，进入步骤(6)；否则，返回步骤(4)。
在具体实施过程中，相对误差位
Figure BDA0003534509740000164
选择收敛阈值eT＝0.0001。
(6)：令q＝q(k)，将声源强度向量q变换成大小为N1×N2的矩阵
Figure BDA0003534509740000165
并将矩阵
Figure BDA0003534509740000166
中的各元素进行对数变换，形成声像Q用于声源辨识。
在具体实施过程中，对数变换
Figure BDA0003534509740000167
其中
Figure BDA0003534509740000168
是矩阵
Figure BDA0003534509740000169
中的元素，Qref是参考声源强度，选择参考声源强度Qref＝0.00002。
如此获取声源强度矩阵，从而形成声像，如此将声像用于声源辨识。
如图2、图3所示，通过本发明实施例中的用于麦克风阵列测量的稀疏声源辨识方法，距离麦克风阵列1m处的两个声源能够被辨识出来，所得声像中的声源位置与实际声源位置一致，所得声像中的声源强度与实际声源强度十分接近，这一结果说明了本发明实施例中的用于麦克风阵列测量的稀疏声源辨识方法能够实现稀疏声源辨识的目标；如图3所示，本发明实施例中的技术方案与现有技术进行了对比，计算了不同方法所得的声源强度向量与真实声源强度向量间的均方误差，相比于现有技术，本发明实施例中的技术方案所得的声源辨识结果的均方误差更小，这一结果说明了本发明实施例中的用于麦克风阵列测量的稀疏声源辨识方法能够提高稀疏声源辨识的准确性。另外，在相同的硬件条件下记录了上述方法的计算时间，本发明实施例中的用于麦克风阵列测量的稀疏声源辨识方法的技术方案计算时间为0.61秒，基于1范数的方法计算时间为0.44秒，基于近似0范数的方法计算时间为0.72秒，基于正交匹配追踪的方法计算时间为0.41秒，表明本发明实施例中的用于麦克风阵列测量的稀疏声源辨识方法能够在保持较高的计算效率的条件下取得最准确的声源辨识结果，在结果准确性和计算效率方面达到更好的平衡，因此本发明具有较强的实用性。
综上所述，本发明提供的用于麦克风阵列测量的稀疏声源辨识方法，首先利用麦克风阵列采集声信号，以构建大小为M×La的原始声压矩阵；对各个通道所采集的声信号进行分割处理以形成信号段，并根据原始声压矩阵获取声压矩阵，再利用离散傅里叶变换将信号转换到频域，构成频域声压矩阵，并在频域声压矩阵中选择信号成分，以构成声压向量，再根据预设规则选择声源重构区域，将声源重构区域划分离散网格，根据麦克风阵列中阵元的位置坐标和声压向量建立声传播模型，而后根据声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对声学逆问题进行变形以获取变形表达式，采用交替分裂的方式对变形表达式进行求解，以获取待求的目标声源强度，并根据目标声源强度形成声像，完成声源辨识；如此，利用广义最小最大凹函数作为正则化函数建立稀疏声源辨识方法，采用的近似函数对0范数的近似效果更好，非零元素在稀疏恢复过程中受到的惩罚程度的一致性更强，能够解决现有技术对声源强度欠估计的问题；并且声学逆问题的目标函数具有凸性，能够保证稀疏声源辨识结果的准确性，引入辅助变量后，采用交替分裂的方式对声学逆问题进行求解，建立了显式的变量迭代更新公式，能够方便地对声学逆问题进行求解，提高计算效率；而且考虑信号交叠分割策略和不同的声源频率以适应外部环境的变化，求解过程采用的稀疏先验信息能够抑制测量噪声带来的干扰，提高辨识精度。
如图4所示，本发明还提供一种用于麦克风阵列测量的稀疏声源辨识系统100，实现如前所述的用于麦克风阵列测量的稀疏声源辨识方法，包括：
声压矩阵获取单元101，用于利用包含M个通道的麦克风阵列采集长度为La的声信号pIm(m＝1,2,...,M)，以构建大小为M×La的原始声压矩阵
Figure BDA0003534509740000171
对各个通道所采集的声信号pIm进行分割处理以形成长度为L、交叠率为R的信号段，并根据所述原始声压矩阵获取包含有
Figure BDA0003534509740000181
个快拍的声压矩阵
Figure BDA0003534509740000182
其中，
Figure BDA0003534509740000183
表示不大于
Figure BDA0003534509740000184
的整数，第s(s＝1,2,...,S)个快拍对应的子矩阵Ps II＝[ps,1,ps,2,...,ps,M]T大小为M×L，
Figure BDA0003534509740000185
表示
Figure BDA0003534509740000186
中的第L(1-R)(s-1)1到第L(1-R)(s-1)+L个元素；
频域转换单元102，用于利用离散傅里叶变换将第s个快拍对应的信号转换到频域，构成频域声压矩阵
Figure BDA0003534509740000187
并在所述频域声压矩阵中选择第l(l＝1,2,...,L)个频率fl所对应的信号成分，以构成声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T，其中，
Figure BDA0003534509740000188
DFT表示离散傅里叶变换；pm(fl)代表
Figure BDA0003534509740000189
中频率fl对应的信号成分；
模型建立单元103，用于根据预设规则选择声源重构区域，将所述声源重构区域划分为N＝N1×N2个离散网格，所述离散网格中网格点的位置坐标为un＝(unx,uny,unz)(n＝1,2,...,N)；根据所述麦克风阵列中阵元的位置坐标和所述声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T建立声传播模型p＝Gq+r，其中，q为待求解的声源强度向量，r为噪声向量，G为大小为M×N的传输矩阵，且所述传输矩阵的元素为
Figure BDA00035345097400001810
其中exp()是指数函数，c是声速，
Figure BDA00035345097400001811
是虚数单位；
辨识求解单元104，用于根据所述声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对所述声学逆问题进行变形以获取变形表达式，采用交替分裂的方式对所述变形表达式进行求解，以获取待求的目标声源强度，并根据所述目标声源强度形成声像，完成声源辨识。
具体的实施方法不做具体限制，在此不做过多赘述，可参照用于麦克风阵列测量的稀疏声源辨识方法的具体实施例。
如上所述，本发明提供的用于麦克风阵列测量的稀疏声源辨识系统，首先通过声压矩阵获取单元101利用麦克风阵列采集声信号，以构建原始声压矩阵；对各个通道所采集的声信号进行分割处理以形成信号段，并根据原始声压矩阵获取声压矩阵，再通过频域转换单元102利用离散傅里叶变换将信号转换到频域，构成频域声压矩阵，并在频域声压矩阵中选择信号成分，以构成声压向量，再通过模型建立单元103根据预设规则选择声源重构区域，将声源重构区域划分离散网格，根据麦克风阵列中阵元的位置坐标和声压向量建立声传播模型，而后通过辨识求解单元104根据声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对声学逆问题进行变形以获取变形表达式，采用交替分裂的方式对变形表达式进行求解，以获取待求的目标声源强度，并根据目标声源强度形成声像，完成声源辨识；如此，利用广义最小最大凹函数作为正则化函数建立稀疏声源辨识方法，采用的近似函数对0范数的近似效果更好，非零元素在稀疏恢复过程中受到的惩罚程度的一致性更强，能够解决现有技术对声源强度欠估计的问题；并且声学逆问题的目标函数具有凸性，能够保证稀疏声源辨识结果的准确性，引入辅助变量后，采用交替分裂的方式对声学逆问题进行求解，建立了显式的变量迭代更新公式，能够方便地对声学逆问题进行求解，提高计算效率；而且考虑信号交叠分割策略和不同的声源频率以适应外部环境的变化，求解过程采用的稀疏先验信息能够抑制测量噪声带来的干扰，提高辨识精度。
如上参照附图以示例的方式描述了根据本发明提出的用于麦克风阵列测量的稀疏声源辨识方法、系统。但是，本领域技术人员应当理解，对于上述本发明所提出的用于麦克风阵列测量的稀疏声源辨识方法、系统，还可以在不脱离本发明内容的基础上做出各种改进。因此，本发明的保护范围应当由所附的权利要求书的内容确定。

Claims (10)

1.一种用于麦克风阵列测量的稀疏声源辨识方法，其特征在于，包括：
利用包含M个通道的麦克风阵列采集长度为La的声信号pIm(m＝1,2,...,M)，以构建大小为M×La的原始声压矩阵
Figure FDA0003534509730000011
对各个通道所采集的声信号pIm进行分割处理以形成长度为L、交叠率为R的信号段，并根据所述原始声压矩阵获取包含有
Figure FDA0003534509730000012
个快拍的声压矩阵
Figure FDA0003534509730000013
其中，
Figure FDA0003534509730000014
表示不大于
Figure FDA0003534509730000015
的整数，第s(s＝1,2,...,S)个快拍对应的子矩阵
Figure FDA0003534509730000016
大小为M×L，
Figure FDA0003534509730000017
表示
Figure FDA0003534509730000018
中的第L(1-R)(s-1)+1到第L(1-R)(s-1)+L个元素；
利用离散傅里叶变换将第s个快拍对应的信号转换到频域，构成频域声压矩阵
Figure FDA0003534509730000019
并在所述频域声压矩阵中选择第l(l＝1,2,...,L)个频率fl所对应的信号成分，以构成声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T，其中，
Figure FDA00035345097300000110
DFT表示离散傅里叶变换；pm(fl)代表
Figure FDA00035345097300000111
中频率fl对应的信号成分；
根据预设规则选择声源重构区域，将所述声源重构区域划分为N＝N1×N2个离散网格，所述离散网格中网格点的位置坐标为un＝(unx,uny,unz)(n＝1,2,...,N)；根据所述麦克风阵列中阵元的位置坐标和所述声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T建立声传播模型p＝Gq+r，其中，q为待求解的声源强度向量，r为噪声向量，G为大小为M×N的传输矩阵，且所述传输矩阵的元素为
Figure FDA00035345097300000112
其中exp( )是指数函数，c是声速，
Figure FDA0003534509730000021
是虚数单位；
根据所述声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对所述声学逆问题进行变形以获取变形表达式，采用交替分裂的方式对所述变形表达式进行求解，以获取待求的目标声源强度，并根据所述目标声源强度形成声像，完成声源辨识。
2.如权利要求1所述的用于麦克风阵列测量的稀疏声源辨识方法，其特征在于，在根据所述声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对所述声学逆问题进行变形以获取变形表达式的过程中，包括：
根据所述声传播模型建立基于广义最小最大凹函数的声学逆问题
Figure FDA0003534509730000022
其中，λ＞0为需要给定的正则化参数，
Figure FDA0003534509730000023
为广义最小最大凹函数，B为辅助矩阵；v为辅助向量，且为变量；
给定平衡参数0＜γ＜1，令所述辅助矩阵
Figure FDA0003534509730000024
以对所述声学逆问题进行变形确定所述声学逆问题的变形表达式：
Figure FDA0003534509730000025
将所述变形表达式作为目标函数，并将所述目标函数对于变量q的平滑项记为f1(q)，将所述目标函数对于变量q的非平滑项记为f2(q)，将所述目标函数对于变量v的平滑项记为g1(v)，将所述目标函数对于变量v的非平滑项记为g2(v)。
3.如权利要求2所述的用于麦克风阵列测量的稀疏声源辨识方法，其特征在于，采用交替分裂的方式对所述变形表达式进行求解，以获取待求的目标声源强度，并根据所述目标声源强度形成声像的过程，包括：
初始化计数变量，并随机初始化中间变量
Figure FDA0003534509730000031
和
Figure FDA0003534509730000032
随着计数变量的更新根据所述变形表达式依次迭代更新计算待求的声源强度向量的值，并计算两个相邻次数所求得的声源强度向量的值的强度误差；
判断所述强度误差是否超过预设的收敛阈值；若所述强度误差超过预设的收敛阈值，则继续迭代更新计算待求的声源强度向量的值；若所述强度误差不超过所述收敛阈值，则将最后一次计算的声源强度向量的值作为待求的目标声源强度；
根据所述目标声源强度生成声像。
4.如权利要求3所述的用于麦克风阵列测量的稀疏声源辨识方法，其特征在于，
所述收敛阈值为0.0001。
5.如权利要求3所述的用于麦克风阵列测量的稀疏声源辨识方法，其特征在于，
所述交替分裂的方式为Douglas-Rachford分裂的方式。
6.如权利要求5所述的用于麦克风阵列测量的稀疏声源辨识方法，其特征在于，随着计数变量的更新根据所述变形表达式依次迭代更新计算待求的声源强度向量的值，并计算两个相邻次数所求得的声源强度向量的值的强度误差的过程，包括：
计算当前计数变量下当前的声源中间变量的值和当前的声源强度向量的值；其中，当前的声源中间变量为：
Figure FDA0003534509730000033
当前的声源强度向量为：
Figure FDA0003534509730000034
令所述当前的声源强度向量中k＝k+1，计算所述当前计数变量的下一变量的源强度向量的值，并计算所述下一变量的源强度向量的值与所述当前的声源强度向量的值的强度误差eR。
7.如权利要求6所述的用于麦克风阵列测量的稀疏声源辨识方法，
所述强度误差为
Figure FDA0003534509730000041
8.如权利要求7所述的用于麦克风阵列测量的稀疏声源辨识方法，其特征在于，在计算当前计数变量下当前的声源中间变量的值和当前的声源强度向量的值的过程中，
还计算当前的辅助向量的值和当前的辅助中间变量的值；其中，
当前的辅助向量的为：
Figure FDA0003534509730000042
当前的辅助中间变量为：
Figure FDA0003534509730000043
其中，
Figure FDA0003534509730000044
为预算的近端算子。
9.如权利要求8所述的用于麦克风阵列测量的稀疏声源辨识方法，其特征在于，若所述强度误差不超过所述收敛阈值，则将最后一次计算的声源强度向量的值作为待求的目标声源强度；根据所述目标声源强度生成声像的过程，包括：
将最后一次计算的声源强度向量的值作为待求的目标声源强度的值；
将所述目标声源强度转换为大小为N1×N2的声源强度矩阵
Figure FDA0003534509730000045
并对所述声源强度矩阵
Figure FDA0003534509730000046
中的各个元素进行对数变换以获取目标声源矩阵
Figure FDA0003534509730000047
其中，
Figure FDA0003534509730000048
为矩阵
Figure FDA0003534509730000049
中的元素，Qref为预设的参考声源强度；
根据所述目标声源矩阵生成声像，根据所述声像进行声源辨识。
10.一种用于麦克风阵列测量的稀疏声源辨识系统，其特征在于，实现权利要求1-9任一所述的用于麦克风阵列测量的稀疏声源辨识方法，包括：
声压矩阵获取单元，用于利用包含M个通道的麦克风阵列采集长度为La的声信号pIm(m＝1,2,...,M)，以构建大小为M×La的原始声压矩阵
Figure FDA0003534509730000051
对各个通道所采集的声信号pIm进行分割处理以形成长度为L、交叠率为R的信号段，并根据所述原始声压矩阵获取包含有
Figure FDA0003534509730000052
个快拍的声压矩阵
Figure FDA0003534509730000053
其中，
Figure FDA0003534509730000054
表示不大于
Figure FDA0003534509730000055
的整数，第s(s＝1,2,...,S)个快拍对应的子矩阵
Figure FDA0003534509730000056
大小为M×L，
Figure FDA0003534509730000057
表示
Figure FDA0003534509730000058
中的第L(1-R)(s-1)+1到第L(1-R)(s-1)+L个元素；
频域转换单元，用于利用离散傅里叶变换将第s个快拍对应的信号转换到频域，构成频域声压矩阵
Figure FDA0003534509730000059
并在所述频域声压矩阵中选择第l(l＝1,2,...,L)个频率fl所对应的信号成分，以构成声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T，其中，
Figure FDA00035345097300000510
DFT表示离散傅里叶变换；pm(fl)代表
Figure FDA00035345097300000511
中频率fl对应的信号成分；
模型建立单元，用于根据预设规则选择声源重构区域，将所述声源重构区域划分为N＝N1×N2个离散网格，所述离散网格中网格点的位置坐标为un＝(unx,uny,unz)(n＝1,2,...,N)；根据所述麦克风阵列中阵元的位置坐标和所述声压向量p＝[p1(fl),p2(fl),...,pM(fl)]T建立声传播模型p＝Gq+r，其中，q为待求解的声源强度向量，r为噪声向量，G为大小为M×N的传输矩阵，且所述传输矩阵的元素为
Figure FDA00035345097300000512
其中exp( )是指数函数，c是声速，
Figure FDA00035345097300000513
是虚数单位；
辨识求解单元，用于根据所述声传播模型建立基于广义最小最大凹函数的声学逆问题，引入辅助矩阵对所述声学逆问题进行变形以获取变形表达式，采用交替分裂的方式对所述变形表达式进行求解，以获取待求的目标声源强度，并根据所述目标声源强度形成声像，完成声源辨识。
CN202210223139.7A 2022-03-07 2022-03-07 用于麦克风阵列测量的稀疏声源辨识方法、系统 Pending CN115278496A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210223139.7A CN115278496A (zh) 2022-03-07 2022-03-07 用于麦克风阵列测量的稀疏声源辨识方法、系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210223139.7A CN115278496A (zh) 2022-03-07 2022-03-07 用于麦克风阵列测量的稀疏声源辨识方法、系统

Publications (1)

Publication Number Publication Date
CN115278496A true CN115278496A (zh) 2022-11-01

Family

ID=83758602

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210223139.7A Pending CN115278496A (zh) 2022-03-07 2022-03-07 用于麦克风阵列测量的稀疏声源辨识方法、系统

Country Status (1)

Country Link
CN (1) CN115278496A (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115422693A (zh) * 2022-11-02 2022-12-02 杭州兆华电子股份有限公司 一种区域管路漏损声波检测阵列的位置设计方法及系统
CN117825898A (zh) * 2024-03-04 2024-04-05 国网浙江省电力有限公司电力科学研究院 一种gis分布式振声联合监测方法、装置及介质

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2014048399A (ja) * 2012-08-30 2014-03-17 Nippon Telegr & Teleph Corp <Ntt> 音響信号解析装置、方法、及びプログラム
CN112526451A (zh) * 2020-10-29 2021-03-19 清华大学 基于麦克风阵列成像的压缩波束形成及系统
CN113238189A (zh) * 2021-05-24 2021-08-10 清华大学 基于阵列测量和稀疏先验信息的声源辨识方法、系统

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2014048399A (ja) * 2012-08-30 2014-03-17 Nippon Telegr & Teleph Corp <Ntt> 音響信号解析装置、方法、及びプログラム
CN112526451A (zh) * 2020-10-29 2021-03-19 清华大学 基于麦克风阵列成像的压缩波束形成及系统
CN113238189A (zh) * 2021-05-24 2021-08-10 清华大学 基于阵列测量和稀疏先验信息的声源辨识方法、系统

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
MIAOFAN LI ET AL: ""Synchro-Reassigning Transform for Instantaneous Frequency Estimation and Signal Reconstructuion"", 《IEEE》, 4 August 2021 (2021-08-04) *
李禹志;王青东;郑海祥;郭各朴;马青玉;龚秀芬;: "环形稀疏声源阵列声涡旋的优化与操控", 声学学报, no. 05, 15 September 2016 (2016-09-15) *

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115422693A (zh) * 2022-11-02 2022-12-02 杭州兆华电子股份有限公司 一种区域管路漏损声波检测阵列的位置设计方法及系统
CN115422693B (zh) * 2022-11-02 2023-03-24 杭州兆华电子股份有限公司 一种区域管路漏损声波检测阵列的位置设计方法及系统
CN117825898A (zh) * 2024-03-04 2024-04-05 国网浙江省电力有限公司电力科学研究院 一种gis分布式振声联合监测方法、装置及介质
CN117825898B (zh) * 2024-03-04 2024-06-11 国网浙江省电力有限公司电力科学研究院 一种gis分布式振声联合监测方法、装置及介质

Similar Documents

Publication Publication Date Title
CN112526451B (zh) 基于麦克风阵列成像的压缩波束形成及系统
CN115278496A (zh) 用于麦克风阵列测量的稀疏声源辨识方法、系统
WO2009145310A1 (ja) 音源分離及び表示方法並びにシステム
Kujawski et al. A deep learning method for grid-free localization and quantification of sound sources
KR102316671B1 (ko) Cnn을 이용한 음향 처리방법
CN110109058A (zh) 一种平面阵列反卷积声源识别方法
CN109884592B (zh) 一种面向低频高斯噪声源的声源定位仿真方法
Wang et al. Estimation of multiple sound sources with data and model uncertainties using the EM and evidential EM algorithms
CN111798869B (zh) 一种基于双麦克风阵列的声源定位方法
Traa et al. Multichannel source separation and tracking with RANSAC and directional statistics
CN107765221A (zh) 适用于相干和非相干声源的反卷积声源成像算法
KR100921368B1 (ko) 이동형 마이크로폰 어레이를 이용한 소음원 위치 판별정밀도 개선 시스템 및 방법
JP4293986B2 (ja) 音場を表す方法及びシステム
CN109870669A (zh) 一种二维多快拍无网格压缩波束形成声源识别方法
CN109410178A (zh) 一种工件裂纹检测方法及系统
CN109923430B (zh) 用于进行相位差展开的装置及方法
CN111812581A (zh) 基于原子范数的球面阵列声源波达方向估计方法
JP7254938B2 (ja) 音響源用の結合音源定位及び分離方法
Katzberg et al. Spatial interpolation of room impulse responses using compressed sensing
Bianchi et al. Deconvolution of plenacoustic images
CN115267673A (zh) 考虑重建网格偏移的稀疏声源成像方法、系统
CN110956978B (zh) 一种基于欠定卷积混叠模型的稀疏盲分离方法
CN108594163B (zh) 用单个传声器移动结合内积运算识别声源的方法
JP4676920B2 (ja) 信号分離装置、信号分離方法、信号分離プログラム及び記録媒体
Lobato et al. Near-field acoustic holography using far-field measurements

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination