[发明专利]一种冲击噪声环境下的量子万有引力搜索动态DOA估计方法

权利要求书

1.一种冲击噪声环境下的量子万有引力搜索动态DOA估计方法，其特征在于： 

步骤一，对空间D个信源信号进行快拍采样和加权处理，由M个全向阵元构成的无模糊线阵，D个中心频率相同的远场窄带信号从方向{θ₁,θ₂,…,θ_D}入射到阵列上，阵列接收第k次快拍采样数据为x(k)＝As(k)+n(k)，式中x(k)＝[x₁(k),x₂(k),…,x_M(k)]^T为阵列的接收数据矢量；A(θ)＝[a(θ₁) a(θ₂) … a(θ_D)]是信号导向矢量矩阵，θ＝(θ₁,θ₂,…,θ_D)是角度矢量，θ_d（d＝1,2,…,D）是第d个信源的来波方向；s(t)＝[s₁(t),s₂(t),…,s_D(t)]^T为信号矢量； 

n(k)＝[n₁(k),n₂(k),…,n_M(k)]^T是独立同分布的满足SαS分布的冲击噪声向量；第d个导向矢量为a(θ_d)＝[exp(-jω₀τ_1d),exp(-jω₀τ_2d),…,exp(-jω₀τ_Md)]^Τ，d＝1,2,…,D，式中， c为光速，f为入射信号的频率，λ为入射信号的波长，(l＝1,2,…,M)表示第d个信号到达第l个阵元时相对于参考阵元的时延，E_l为第l个阵元的位置，接收数据的加权信号为h(k)＝x(k)/max{|x₁(k)|^0.7,|x₂(k)|^0.7,…,|x_M(k)|^0.7}，max{}为取最大值函数，定义阵元接收数据之间的加权信号协方差为C(k)，第一次采样数据的加权信号协方差为C(1)＝h(1)h^H(1)，根据经典信息论准则估计信源个数D，t代表量子万有引力搜索机制的迭代次数，初始时设t＝0； 

步骤二，量子万有引力搜索机制的所有物体初始状态确定，由N个物体组成的系统，每个物体在D维搜索空间中运动，第i个物体的量子位置定义如下： 

(i＝1,2,…,N)，为有效量子位置，为辅助量子位置，0≤y_id(t)≤1(d＝1,2,…,D)表示第i个物体在第d维的有效量子位，z_id为第i个物体在第d维的辅助量子位，[y_id(t)]²+[z_id(t)]²＝1，把第i个物体的量子位置第一行即有效量子位置影射到定义区间，就是物体的当前位置q_i(t)＝[q_i1(t),q_i2(t),…,q_iD(t)]，对于所有的i＝1,2,…,N，q_id(t)＝y_id(t)[u_d(k)-l_d(k)]+l_d(k)，d＝1,2,…,D，第i个物体对应的速度为v_i(t)＝[v_i1(t),v_i2(t),…,v_iD(t)]，-0.2≤v_id(t)≤0.2(d＝1,2,…,D)表示第i个物体的第d维速度，t代表迭代次数，量子位置在定义的量子域[0，1]间随机初始化，在[-0.2，0.2]随机初始 化速度的每一维，确定初始搜索空间为D个角度的搜索区间向量，u_d(k)和l_d(k)分别为第k次采样数据处理第d维的角度搜索区间的上限和下限值，最大迭代次数设置为，其中为取整函数，zz为取整倍数； 

步骤三，第i个物体位置q_i(t)＝[q_i1(t),q_i2(t),…,q_iD(t)]的加权信号协方差极大似然方程的目标函数值为，正交投影矩阵 ，tr()表示矩阵的求迹运算，根据目标函数 评价物体位置的优劣，代表有效量子位置 映射到搜索范围中的位置计算的目标函数值，即把影射到参数范围后，使用fit[q_i(t)]计算物体的目标函数值，第i个物体到现在为止所经历的最优量子位置的第一行定义为该物体的局部最优位置，记作p_i(t)＝[p_i1(t),p_i2(t),…,p_iD(t)]，p_id(t)为第i个物体到第t次迭代为止所经历最优有效第d维量子位，d＝1,2,…,D，所有物体到现在所经历最优量子位置的第一行记作全局最优量子位置，也就是目标函数值最大的有效量子位置，记作b(t)＝[b₁(t),b₂(t),…,b_D(t)]，b_d(t)为所有物体到第t次迭代为止所经历最优有效第d维量子位，d＝1,2,…,D，确定最好目标函数值best(t)和最差目标函数值worst(t)， 

步骤四，计算在t次迭代的每个物体质量：M_ai＝M_pi＝M_ii＝M_i,i＝1,2,…,N, M_pi表示被作用物体i的质量，M_ai表示作用到物体i的质量，M_ii为物体i的引力质量，M_i(t)为第i个物体的惯性质量，计算N个物体的万有引力的合力，对于第d维，第j个物体作用在第i个物体上的引力大小为 引力系数T为最大迭代次数，R_ij(t)是第j个物体和第i个物体的欧氏距离，G₀和γ为常系数，R_ij(t)＝||p_i(t),p_j(t)||₂，对于 第d维，作用在第i个物体上的总作用力i＝1,2,…,N，d＝1,2,…,D，rand_j是[0，1]间的随机数，更新所有物体的加速度、速度和量子位置，第i个物体的第d维加速度为前N/2个物体速度根据引力方程进行更新，则第i个物体的第d维速度更新为 i＝1,2,…,N/2，d＝1,2,…,D；其他N/2个物体速度根据群行为进行更新，则第i个物体的第d维速度更新为v_id(t+1)＝c₄r₄[b_d(t)-y_id(t)]+c₅r₅[p_id(t)-y_id(t)]+w_tv_id(t)，w_t从初次迭代的0.9线形递减到最大迭代的0.5，i＝N/2+1,N/2+2,…,N，d＝1,2,…,D，r₁、r₂、r₃、r₄和r₅都是[0，1]之间的均匀随机数，加权常数c₁＝c₂＝c₄＝c₅＝0.78，c₃＝0.1；对于v_id(t+1)，若超出边界值，将其限制在边界，即v_id(t+1)＞0.2，v_id(t+1)＝0.2，若v_id(t+1)<-0.2，v_id(t+1)＝-0.2，量子位置为i＝1,2,…,N，其中 d＝1,2,…,D，abs()为量子位取绝对值函数； 

步骤五，量子位置y_i(t+1)影射到定义区间q_i(t+1)，计算目标函数值， 

i＝1,2,…,N，d＝1,2,…,D，若 

fit[q_i(t+1)]＞fitq[p_i(t)]，则令局部最优量子位置否则, 

p_i(t+1)＝p_i(t)，把最优的局部最优量；子位置设置成全局最优量子位置为 

b(t+1)＝[b₁(t+1),b₂(t+1),…,b_D(t+1)]，确定当前最好目标函数值best(t+1)和最差目标函数值worst(t+1)，

步骤六，判断是否达到最大迭代次数，若是，记录最优位置，执行步骤七；否则，t＝t+1，返回步骤四； 

步骤七，快拍采样新数据x(k+1)＝[x₁(k+1),x₂(k+1),…,x_M(k+1)]^T，对采样数据进行加权处 理：h(k+1)＝x(k+1)/max{|x₁(k+1)|^0.7,|x₂(k+1)|^0.7,…,|x_M(k+1)|^0.7}，计算当前采样的加权信号协方差：其中是新增加的第k+1个采样数据的加权信号协方差的增量，更新搜索空间为Z(k+1)， d＝1,2,…,D，β为收敛因子，决定了搜索空间的收敛速度；常数r为搜索空间在锁定状态下的搜索半径；为第d个方向在第k次采样时的估计值；为第d个方向在第k次采样时搜索空间的中心值，其更新公式为其中δ为遗传因子； 

步骤八，如果达到最大跟踪次数，执行步骤九；否则，设k＝k+1，t＝0，返回步骤二继续估计动态目标下一个时刻的方向； 

步骤九，得到所有快拍采样下的全局最优位置就是检测的动态目标方向值，输出动态跟踪结果。