一种粒子优化加权信号子空间拟合的阵列直接自定位方法

未命名 08-15 阅读:153 评论:0


1.本发明属于无源定位技术领域,具体地,涉及一种粒子优化加权信号子空间拟合的阵列直接自定位方法。


背景技术:

2.定位技术是信号处理领域的一个重要研究方向,在各种工业领域有着广泛的应用。在现代社会的发展中,目标位置的确定是智能感知技术的基础,利用无线传感器节点进行精确的定位成为很多技术的核心,所以定位技术一直是学术界和工业界的研究热点方向。经过多年的发展,定位技术算法在越来越多的场景中得到大规模应用。但是,科技发展带来了新应用场景的涌现,这些变化都对定位技术提出了更高的要求。
3.自定位技术是近几年人们广泛关注的话题,传统的自定位技术使用rss、tdoa、aoa这些两步定位技术进行辅助定位,通过最小二乘定位算法实现对阵列的实时定位,但是从信息论的角度,两步定位方法只利用了数据中的角度信息,因而两步定位方法的结果存在信息的损失。
4.粒子滤波算法是一种适用性广的误差修正方法,既能够很好地运用于线性系统,也能够改善非线性系统的性能,所以被广泛地应用于传感器误差的修正之中。在直接定位中,网格对空间的量化也是一种感知误差,所以使用粒子滤波算法可以有效减少网格量化的误差。因此,研究能够使用粒子优化加权信号子空间拟合的阵列直接定位算法具有重要的意义。


技术实现要素:

5.本发明针对现有技术中的不足,提供一种粒子优化加权信号子空间拟合的阵列直接自定位方法,基于加权信号子空间拟合方法,利用粒子滤波对定位结果进行优化,实现对信源位置的精准估计。
6.为实现上述目的,本发明采用以下技术方案:一种粒子优化加权信号子空间拟合的阵列直接自定位方法,具体包括以下步骤:
7.步骤1:在定位区域分布l个信源以及用于接收信源信号的阵列;所述阵列为由m个阵元组成的圆阵;
8.步骤2:通过圆阵采集数据的模型对多信源信号进行数据采集;
9.步骤3:对定位区域进行网格点的均匀选取,并在每个网格点根据圆阵采集数据的模型构造所有信源的特征矩阵;
10.步骤4:对采集的数据计算协方差矩阵,并进行特征值分解,得到信号子空间、噪声子空间和对角矩阵,并计算加权信号子空间拟合所需要的权值矩阵;
11.步骤5:通过每个网格点构造的权值矩阵、信号子空间以及所有信源的特征矩阵进行加权信号子空间函数值拟合,对每个网格点的加权信号子空间函数值进行比较,选取函数值最大的点作为圆阵所在位置的初始估计;
12.步骤6:在初始估计位置附近随机选取粒子,粒子分布满足均值为0、方差为σ2的正态分布关系;
13.步骤7:对各个粒子位置通过粒子加权信号子空间函数进行权重评估;
14.步骤8:归一化权重评估值,计算权重对粒子分布位置的加权平均和,得到圆阵位置估计结果。
15.进一步地,步骤2中利用圆阵对多信源信号进行数据采集的过程为:
16.x=as+n
17.其中,x为圆阵采集的信号,s为信源信号矢量,n为噪声信号矢量,a为圆阵的阵列流型,a=[a1,

,a
l
,

,a
l
],a
l
为第l个信源入射到阵列方位角形成的导向矢量,λ为信源发射信号波长,r为圆阵的半径,β
l
为圆阵测量的第l个信源的俯仰角,γ
m-1
为第m个阵元与参考阵元相对于圆阵中心位置的夹角,γ
m-1
=2π(m-1)/m,j为虚数单位,θ
l
表示圆阵测得的第l个辐射源的入射方位角。
[0018]
进一步地,步骤3中每个网格点构造的所有信源的特征矩阵为:
[0019]
φi=[φ
1,i
,


l,i
,


l,i
]
[0020][0021]
其中,φi为第i个网格点构造的所有信源的特征矩阵,φ
l,i
为第i个网格点中构造的第l个信源的特征矩阵,λ为信源发射信号波长,r为圆阵的半径,β
l
为圆阵测量的第l个信源的俯仰角,η
l,i
为在第i个网格点位置第l个信源和正北方向的夹角,为圆阵中参考阵元与正北方向的夹角,γ
m-1
为第m个阵元与参考阵元相对于圆阵中心位置的夹角,γ
m-1
=2π(m-1)/m,j为虚数单位,η
l,i
为在第i个网格点位置第l个信源和正北方向的夹角,(xi,yi)为第i个网格
点的位置坐标,(p
lx
,p
ly
)为第l个信源的位置坐标。
[0022]
进一步地,步骤4包括如下子步骤:
[0023]
步骤4.1:对采集的数据计算协方差矩阵,对协方差矩阵进行特征值分解,得到信号子空间、噪声子空间和对角矩阵:
[0024][0025]
其中,r为采集数据的协方差矩阵,e()为数学期望,ds为l
×
l维的对角矩阵,ds的对角元素由特征值分解所得的最大的l个特征值构成;dn为由m-l个最小的特征值构成的对角矩阵;us为信号子空间,信号子空间由l个最大的特征值对应的特征向量构成;un为噪声子空间,噪声子空间由m-l个最小的特征值对应的特征向量构成;
[0026]
步骤4.2:通过由最大的l个特征值构成的对角矩阵ds以及由m-l个最小的特征值构成的对角矩阵计算加权信号子空间拟合所需要的权值矩阵其中,i为单位矩阵,σ2为噪声功率,由dn对角线元素的平均值构成。
[0027]
进一步地,步骤5中加权信号子空间函数值的拟合过程为:
[0028][0029]
其中,f(pi)为第i个网格点的加权信号子空间函数值,tr()为矩阵的迹运算,i是单位矩阵,us为信号子空间,w为权值矩阵,φi为第i个网格点构造的所有信源的特征矩阵。
[0030]
进一步地,步骤6中随机选取的j个粒子的位置为paj=[paxj,payj]
t
,并且粒子的位置分布满足如下关系:
[0031][0032][0033]
其中,||||2为二范数运算,n(0,σ2)表示服从均值为0、方差为σ2的正态分布关系,为圆阵所在位置的初始估计。
[0034]
进一步地,步骤7中粒子加权信号子空间函数为:
[0035][0036]
其中,c(paj)为粒子paj对应的权重评估数值,us为信号子空间,w为权值矩阵,φ
′j为第j个粒子构造的所有信源的特征矩阵,φ
′j=[φ

1,j
,



l,j
,



l,j
],φ

l,j
为第l个信源相对于第j个粒子的构造矩阵:
[0037][0038]
其中,为圆阵中参考阵元与正北方向的夹角,λ为信源发射信号波长,r为圆阵的半径,β
l
为圆阵测量的第l个信源的俯仰角,j为虚数单位,γ
m-1
为第m个阵元与参考阵元相对于圆阵中心位置的夹角,γ
m-1
=2π(m-1)/m,η

l,j
为在第j个粒子位置第l个信源和正北方向的夹角,(p
lx
,p
ly
)为第l个信源的位置坐标。
[0039]
进一步地,步骤8中权重的归一化过程为:
[0040][0041]
其中,c'(paj)为第j个粒子归一化权重评估值。
[0042]
进一步地,所述圆阵位置估计结果为:
[0043][0044]
与现有技术相比,本发明的有益效果是:本发明通过对阵列接收信号协方差矩阵进行特征值分解,得到信号子空间、噪声子空间及相应的特征值对角矩阵,使用权值矩阵弥补信号子空间理论值和实际值的差距;并通过信源的特征矩阵对信号子空间的位置信息进行拟合,相比于传统方法的多次参数估计计算,本发明方法减少了参数多次计算带来的累计误差,能够有效提高自定位的精度;同时,本发明使用粒子滤波方法进行自定位结果的优化,通过对初始估计位置随机选取周围粒子点作为样本,计算粒子的代价函数作为概率密度函数的近似,使用均值代替积分运算,使得定位估计值的误差方差最小,从而减少了网格量化误差。
附图说明
[0045]
图1是本发明粒子优化加权信号子空间拟合的阵列直接自定位方法的流程图;
[0046]
图2是本发明粒子优化加权信号子空间拟合的阵列直接自定位方法的场景图;
[0047]
图3是本发明粒子优化加权信号子空间拟合的阵列直接自定位方法的仿真效果与wssf方法的对比图;
[0048]
图4是本发明粒子优化加权信号子空间拟合的直接自定位方法在不同信噪比下的阵列位置估计性能比较图。
具体实施方式
[0049]
下面结合附图对本发明的技术方案作进一步详细地说明。
[0050]
如图1为本发明粒子优化加权信号子空间拟合的阵列直接自定位方法的流程图,该阵列直接自定位方法具体包括如下步骤:
[0051]
步骤1:在定位区域分布l个信源以及用于接收信源信号的阵列;本发明中阵列为由m个阵元组成的圆阵,相比于其它阵型,圆阵的能够得到目标信源的二维测向数据,并且使用阵元数量相对而言最少;
[0052]
步骤2:通过圆阵采集数据的模型对多信源信号进行数据采集;
[0053]
本发明中利用圆阵对多信源信号进行数据采集的过程为:
[0054]
x=as+n
[0055]
其中,x为圆阵采集的信号,s为信源信号矢量,n为噪声信号矢量,a为圆阵的阵列流型,a=[a1,

,a
l
,

,a
l
],a
l
为第l个信源入射到阵列方位角形成的导向矢量,λ为信源发射信号波长,r为圆阵的半径,β
l
为圆阵测量的第l个信源的俯仰角,γ
m-1
为第m个阵元与参考阵元相对于圆阵中心位置的夹角,γ
m-1
=2π(m-1)/m,j为虚数单位,θ
l
表示圆阵测得的第l个辐射源的入射方位角。
[0056]
步骤3:对定位区域进行网格点的均匀选取,并在每个网格点根据圆阵采集数据的模型构造所有信源的特征矩阵,使用该特征矩阵能够很好地对网格点相对于信源位置约束关系进行描述;
[0057]
本发明中每个网格点构造的所有信源的特征矩阵为:
[0058]
φi=[φ
1,i
,


l,i
,


l,i
]
[0059]
[0060]
其中,φi为第i个网格点构造的所有信源的特征矩阵,φ
l,i
为第i个网格点中构造的第l个信源的特征矩阵,λ为信源发射信号波长,r为圆阵的半径,β
l
为圆阵测量的第l个信源的俯仰角,在二维测向的时候为β
l
=90
°
;为圆阵中参考阵元与正北方向的夹角,γ
m-1
为第m个阵元和参考阵元相对于圆阵中心位置的夹角,γ
m-1
=2π(m-1)/m,j为虚数单位,η
l,i
为在第i个网格点位置第l个信源和正北方向的夹角,在第i个网格点位置第l个信源和正北方向的夹角,(xi,yi)为第i个网格点的位置坐标,(p
lx
,p
ly
)为第l个信源的位置坐标。
[0061]
步骤4:对采集的数据计算协方差矩阵,并进行特征值分解,得到信号子空间、噪声子空间和对角矩阵,并计算加权信号子空间拟合所需要的权值矩阵,由于采集信号的数据有限性、实际环境中存在的噪声等因素的影响,根据阵列采集的数据计算得到的信号子空间与理论值是存在差异的,通过权值矩阵可以弥补这一差异;具体包括如下子步骤:
[0062]
步骤4.1:对采集的数据计算协方差矩阵,对协方差矩阵进行特征值分解,得到信号子空间、噪声子空间和对角矩阵:
[0063][0064]
其中,r为采集数据的协方差矩阵,e()为数学期望,ds为l
×
l维的对角矩阵,ds的对角元素由特征值分解所得的最大的l个特征值构成;dn为由m-l个最小的特征值构成的对角矩阵;us为信号子空间,信号子空间由l个最大的特征值对应的特征向量构成;un为噪声子空间,噪声子空间由m-l个最小的特征值对应的特征向量构成;
[0065]
步骤4.2:通过由最大的l个特征值构成的对角矩阵ds以及由m-l个最小的特征值构成的对角矩阵计算加权信号子空间拟合所需要的权值矩阵其中,i为单位矩阵,σ2为噪声功率,由dn对角线元素的平均值构成。
[0066]
步骤5:通过每个网格点构造的权值矩阵、信号子空间以及所有信源的特征矩阵进行加权信号子空间函数值拟合,对每个网格点的加权信号子空间函数值进行比较,函数值表示的是网格点与阵列真实位置的匹配程度,函数值越大,匹配程度也就越大,所以选取函数值最大的点作为圆阵所在位置的初始估计;
[0067]
本发明中加权信号子空间函数值的拟合过程为:
[0068]
根据阵列信号处理理论,信号子空间张成的空间与阵列流型张成的空间是同一个空间,可知:
[0069]
span{us}=span{a}
[0070]
此时也就存在一个满秩矩阵t,使得
[0071]us
=at
[0072]
由于实际环境中采样次数的有限性、噪声环境的干扰,都会造成信号子空间不能够完全和理论重合,所以可以使用权值矩阵w来弥补信息的缺失,从而可以得到下式的关系:
[0073]usw1/2
=at
[0074]
根据上式两边的关系,将求解阵列位置的问题转换成最小二乘拟合问题:
[0075][0076]
其中,表示阵列初始位置的估计值,坐标表示为表示阵列初始位置的估计值,坐标表示为表示满秩矩阵t的估计值;
[0077]
先固信源的特征矩阵φi,对未知参数进行如下的估计:
[0078][0079]
结合上面两式,可以得到:
[0080][0081]
因此,可以根据上式对网格点进行遍历,在网格点pi处,得到如下的加权信号子空间拟合函数:
[0082][0083]
其中,tr()为矩阵的迹运算,i是单位矩阵。
[0084]
步骤6:在初始估计位置附近随机选取粒子,粒子分布满足均值为0、方差为σ2的正态分布关系,相比于随机均匀分布,这样的正态分布能够保证粒子覆盖的位置信息丰富的同时,也能够约束粒子始终分布在初始估计位置附近;具体地,随机选取的j个粒子的位置为paj=[paxj,payj]
t
,并且粒子的位置分布满足如下关系:
[0085][0086][0087]
其中,||||2为二范数运算,n(0,σ2)表示服从均值为0、方差为σ2的正态分布关系,为圆阵所在位置的初始估计。
[0088]
步骤7:对各个粒子位置通过粒子加权信号子空间函数进行权重评估,对“粒子是阵列位置”的事件进行概率密度函数计算;
[0089]
本发明中粒子加权信号子空间函数为:
[0090][0091]
其中,c(paj)为粒子paj对应的权重评估数值,us为信号子空间,w为权值矩阵,φ
′j为第j个粒子构造的所有信源的特征矩阵,φ
′j=[φ

1,j
,



l,j
,



l,j
],φ

l,j
为第l个信源相对于第j个粒子的构造矩阵:
[0092][0093]
其中,为圆阵中参考阵元与正北方向的夹角,λ为信源发射信号波长,r为圆阵的半径,β
l
为圆阵测量的第l个信源的俯仰角,j为虚数单位,γ
m-1
为第m个阵元与参考阵元相对于圆阵中心位置的夹角,γ
m-1
=2π(m-1)/m,η

l,j
为在第j个粒子位置第l个信源和正北方向的夹角,(p
lx
,p
ly
)为第l个信源的位置坐标。
[0094]
步骤8:归一化权重评估值计算权重对粒子分布位置的加权平均和,得到圆阵位置估计结果从而实现对信源位置的精准估计。
[0095]
为验证本发明粒子优化加权信号子空间拟合的阵列直接自定位方法的有效性,下面通过matlab仿真分析进行证明,性能估计指标为均方根误差(root mean square error,rmse),定义为:
[0096][0097]
式中,k表示蒙特卡洛仿真次数,qk表示第k次蒙特卡洛实验的估计值和真实值。
[0098]
如图2为仿真实验中构建的本发明的场景图,空间中有l个信源,定位区域有一个圆阵接收来自信源的信号数据。
[0099]
图3为本发明粒子优化加权信号子空间拟合的阵列直接自定位方法的仿真效果与加权信号子空间拟合(weighted signal subspace fitting,wssf)方法的对比图,在本次仿真实验中选取的圆阵阵元m=7,圆阵半径r=0.3m;网格范围为长100m、宽100m,网格间距为5m,有三个信源,信源位置为(25m,25m),(50m,75m),(75m,25m),快拍数设置为2048,信噪比条件设置为15db。可以看出本发明的直接自定位方法预测的圆阵位置更接近于真实的圆阵位置,表明通过粒子优化的方法能够更加准确地对圆阵位置进行估计。
[0100]
图4为本发明方法在不同信噪比下的信源位置估计性能比较图,仿真中选取的圆阵阵元m=7,圆阵半径r=0.3m;网格范围为长300m、宽300m,网格间距为2m。在仿真中信源位置随机、圆阵位置随机,快拍数设置为1024,蒙特卡罗实验次数为500。从仿真结果可以看
出,相比于music(multiple signal classification)、isf(initial signal fitting)、wssf方法,在不同信噪比情况下,本发明方法都可以更加准确地估计出信源位置,并且随着信噪比条件的改善,优化性能也会得到更大的提高。
[0101]
以上仅是本发明的优选实施方式,本发明的保护范围并不仅局限于上述实施方式,凡属于本发明思路下的技术方案均属于本发明的保护范围。应当指出,对于本技术领域的普通技术人员来说,在不脱离本发明原理前提下的若干改进和润饰,应视为本发明的保护范围。

技术特征:
1.一种粒子优化加权信号子空间拟合的阵列直接自定位方法,其特征在于,具体包括以下步骤:步骤1:在定位区域分布l个信源以及用于接收信源信号的阵列;所述阵列为由m个阵元组成的圆阵;步骤2:通过圆阵采集数据的模型对多信源信号进行数据采集;步骤3:对定位区域进行网格点的均匀选取,并在每个网格点根据圆阵采集数据的模型构造所有信源的特征矩阵;步骤4:对采集的数据计算协方差矩阵,并进行特征值分解,得到信号子空间、噪声子空间和对角矩阵,并计算加权信号子空间拟合所需要的权值矩阵;步骤5:通过每个网格点构造的权值矩阵、信号子空间以及所有信源的特征矩阵进行加权信号子空间函数值拟合,对每个网格点的加权信号子空间函数值进行比较,选取函数值最大的点作为圆阵所在位置的初始估计;步骤6:在初始估计位置附近随机选取粒子,粒子分布满足均值为0、方差为σ2的正态分布关系;步骤7:对各个粒子位置通过粒子加权信号子空间函数进行权重评估;步骤8:归一化权重评估值,计算权重对粒子分布位置的加权平均和,得到圆阵位置估计结果。2.根据权利要求1所述的一种粒子优化加权信号子空间拟合的阵列直接自定位方法,其特征在于:步骤2中利用圆阵对多信源信号进行数据采集的过程为:x=as+n其中,x为圆阵采集的信号,s为信源信号矢量,n为噪声信号矢量,a为圆阵的阵列流型,a=[a1,

,a
l
,

,a
l
],a
l
为第l个信源入射到阵列方位角形成的导向矢量,λ为信源发射信号波长,r为圆阵的半径,β
l
为圆阵测量的第l个信源的俯仰角,γ
m-1
为第m个阵元与参考阵元相对于圆阵中心位置的夹角,γ
m-1
=2π(m-1)/m,j为虚数单位,θ
l
表示圆阵测得的第l个辐射源的入射方位角。3.根据权利要求1所述的一种优化加权信号子空间拟合的阵列直接自定位方法,其特征在于,步骤3中每个网格点构造的所有信源的特征矩阵为:φ
i
=[φ
1,i
,


l,i
,


l,i
]
其中,φ
i
为第i个网格点构造的所有信源的特征矩阵,φ
l,i
为第i个网格点中构造的第l个信源的特征矩阵,λ为信源发射信号波长,r为圆阵的半径,β
l
为圆阵测量的第l个信源的俯仰角,η
l,i
为在第i个网格点位置第l个信源和正北方向的夹角,为圆阵中参考阵元与正北方向的夹角,γ
m-1
为第m个阵元与参考阵元相对于圆阵中心位置的夹角,γ
m-1
=2π(m-1)/m,j为虚数单位,η
l,i
为在第i个网格点位置第l个信源和正北方向的夹角,(x
i
,y
i
)为第i个网格点的位置坐标,(p
lx
,p
ly
)为第l个信源的位置坐标。4.根据权利要求1所述的一种粒子优化加权信号子空间拟合的阵列直接自定位方法,其特征在于:步骤4包括如下子步骤:步骤4.1:对采集的数据计算协方差矩阵,对协方差矩阵进行特征值分解,得到信号子空间、噪声子空间和对角矩阵:其中,r为采集数据的协方差矩阵,e()为数学期望,d
s
为l
×
l维的对角矩阵,d
s
的对角元素由特征值分解所得的最大的l个特征值构成;d
n
为由m-l个最小的特征值构成的对角矩阵;u
s
为信号子空间,信号子空间由l个最大的特征值对应的特征向量构成;u
n
为噪声子空间,噪声子空间由m-l个最小的特征值对应的特征向量构成;步骤4.2:通过由最大的l个特征值构成的对角矩阵d
s
以及由m-l个最小的特征值构成的对角矩阵计算加权信号子空间拟合所需要的权值矩阵其中,i为单位矩阵,σ2为噪声功率,由d
n
对角线元素的平均值构成。5.根据权利要求1所述的一种粒子优化加权信号子空间拟合的直接自定位方法,其特征在于:步骤5中加权信号子空间函数值的拟合过程为:其中,f(p
i
)为第i个网格点的加权信号子空间函数值,tr()为矩阵的迹运算,i是单位矩阵,u
s
为信号子空间,w为权值矩阵,φ
i
为第i个网格点构造的所有信源的特征矩阵。6.根据权利要求1所述的一种粒子优化加权信号子空间拟合的阵列直接自定位方法,其特征在于:步骤6中随机选取的j个粒子的位置为pa
j
=[pax
j
,pay
j
]
t
,并且粒子的位置分
布满足如下关系:布满足如下关系:其中,||||2为二范数运算,n(0,σ2)表示服从均值为0、方差为σ2的正态分布关系,为圆阵所在位置的初始估计。7.根据权利要求6所述的一种粒子优化加权信号子空间拟合的阵列直接自定位方法,其特征在于:步骤7中粒子加权信号子空间函数为:其中,c(pa
j
)为粒子pa
j
对应的权重评估数值,u
s
为信号子空间,w为权值矩阵,φ

j
为第j个粒子构造的所有信源的特征矩阵,φ

j
=[φ

1,j
,...,φ

l,j
,



l,j
],φ

l,j
为第l个信源相对于第j个粒子的构造矩阵:其中,为圆阵中参考阵元与正北方向的夹角,λ为信源发射信号波长,r为圆阵的半径,β
l
为圆阵测量的第l个信源的俯仰角,j为虚数单位,γ
m-1
为第m个阵元与参考阵元相对于圆阵中心位置的夹角,γ
m-1
=2π(m-1)/m,η

l,j
为在第j个粒子位置第l个信源和正北方向的夹角,(p
lx
,p
ly
)为第l个信源的位置坐标。8.根据权利要求7所述的一种粒子优化加权信号子空间拟合的阵列直接自定位方法,其特征在于:步骤8中权重的归一化过程为:其中,c'(pa
j
)为第j个粒子归一化权重评估值。9.根据权利要求8所述的一种粒子优化加权信号子空间拟合的阵列直接自定位方法,其特征在于:所述圆阵位置估计结果为:

技术总结
本发明公开了一种粒子优化加权信号子空间拟合的阵列直接自定位方法,包括:在定位区域分布信源以及用于接收信源信号的阵列,通过阵列采集数据的模型对多信源信号进行采集;对定位区域进行网格点均匀选取,并在每个网格点构造所有信源的特征矩阵;对采集的数据计算协方差矩阵,进行特征值分解,得到信号子空间,计算加权信号子空间拟合所需要的权值矩阵,结合特征矩阵进行加权信号子空间函数值拟合,选取函数值最大的点作为圆阵所在位置的初始估计;在初始估计位置附近随机选粒子,对各个粒子位置进行权重评估,归一化权重评估值,计算权重对粒子分布位置的加权平均和,得到圆阵位置估计结果。相比于传统方法,本发明的定位估计值的误差方差最小。的误差方差最小。的误差方差最小。


技术研发人员:曹仲康 李建峰 李潘 晋本周 张小飞
受保护的技术使用者:南京航空航天大学
技术研发日:2023.05.04
技术公布日:2023/8/14
版权声明

本文仅代表作者观点,不代表航空之家立场。
本文系作者授权航家号发表,未经原创作者书面授权,任何单位或个人不得引用、复制、转载、摘编、链接或以其他任何方式复制发表。任何单位或个人在获得书面授权使用航空之家内容时,须注明作者及来源 “航空之家”。如非法使用航空之家的部分或全部内容的,航空之家将依法追究其法律责任。(航空之家官方QQ:2926969996)

飞行汽车 https://www.autovtol.com/

分享:

扫一扫在手机阅读、分享本文

相关推荐