CN107886555A - 一种射电干涉阵列分布式图像重建方法 - Google Patents

Abstract

本发明公开一种射电干涉阵列分布式图像重建方法，将射电干涉阵列图像重建的逆问题转换为等式约束条件下的凸优化问题，并使用近似交替方向乘子算法进行求解，最终重建的图像的质量比现有的算法高；通过将数据分割成多份，能够精确快速的分割数据；并使用多计算节点进行计算，具有较高的并行性能，可拓展性强。本发明不需要任何的数据后处理步骤，直接得到重建图像，最终得到的图像的动态范围更大，对延展结构的天体显示更加清晰，相比传统的方法更加高效。

Description

一种射电干涉阵列分布式图像重建方法

技术领域

本发明涉及射电天文成像技术领域，具体涉及一种射电干涉阵列分布式图像重建方法。

背景技术

射电干涉仪的诞生使得射电观测的灵敏度和分辨率有了大幅的提高，该技术有力的推动了天文学、宇宙学和天体物理学诸多方向的研究与发展。在宏伟的科学目标的驱动下，进入21世纪射电干涉阵列的规模不断扩大，例如荷兰的低频阵列(Low FrequencyArray，LOFAR)，核心阵有24个单元(阵元)，每个单元48个天线组成，海外阵有14个单元，每个单元96个天线组成；中国新疆地区建设的21cm阵共有81个单元，每个单元由127个天线组成；美国-澳大利亚联合建成的默奇森广角阵列(Murchison Widefield Array,MWA)共有128个单元，每个单元由16个天线组成；以及预计2019年开始建设的平方公里阵列(SKA)将共有512个单元，每个单元由约256个天线组成。这些干涉阵列共同的特点是产生的原始数据是海量的，且最终得到天文图像的像素尺寸巨大。尤其是SKA，它预期将产生千兆像素的图像并达到六或七个数量级动态范围，这将是从分辨率和灵敏度方面改善目前仪器约两个数量级。获得的数据量如此巨大，用来解决图像重建的方法必须是快速并具有良好的伸缩性。

射电干涉阵列图像的重建是利用射电干涉阵列观测得到的原始数据重建真实天空亮度分布，这个逆问题是一个病态逆问题。现阶段大多是采用洁化算法CLEAN及其延伸算法来实现图像的重建，主要的延伸算法有Cotton&Schwab提出的CLEAN算法简称CS-CLEAN和多尺度CLEAN算法MS-CLEAN。CLEAN和CS-CLEAN算法对点源的重建能够有很好效果，对延展结构的射电源图像重建的效果较差；MS-CLEAN虽然对点源和延展源均有良好的重建效果，但是它的尺度函数是固定的，难以处理大图像。总得来说，这些算法都是依赖本地贪婪迭代过程，不适合大规模并行计算和分布式计算。

发明内容

本发明所要解决的是现有的射电干涉阵列图像重建方法存在运行效率低、可拓展性差、以及依赖于复杂的参数设置的问题，提供一种射电干涉阵列分布式图像重建方法，其在天文大数据下能够快速的重建出精确的天文图像，并且具有良好的可拓展性能(或可伸缩性能)。

为解决上述问题，本发明是通过以下技术方案实现的：

一种射电干涉阵列分布式图像重建方法，包括步骤如下：

步骤1、中央节点将初始化参数广播并发送给计算节点；

步骤2、中央节点根据计算节点的个数n_d，将射电干涉仪的测量数据分割为n_d份；并将这n_d份数据块分别发送给n_d个计算节点；

步骤3、每个计算节点j根据所得到的数据块顺序计算当前迭代t下的松弛变量拉格朗日乘子和中间解并将所计算出的中间解发送给中央节点；

步骤4、中央节点根据各个计算节点返回的中间解计算当前迭代t下的梯度下降解

步骤5、中央节点根据所计算出的梯度下降解使用双重前向-后向方法计算当前迭代t下的最终解x^(t)；

步骤6、判断是否满足收敛条件，即当前迭代t达到最大迭代次数或者当前迭代t下的最终解x^(t)使得每个数据块均满足||Φ_jx^(t)-y_j||₂≤∈_j；

如果未满足收敛条件，则迭代次数t+1，并返回步骤3；

如果满足收敛条件，则当前迭代t下的最终解x^(t)即为重建的图像；

其中，Φ_j表示第j个数据块的测量矩阵，y_j表示第j个数据块的可见度数据，∈_j表示第j个数据块的噪声上限，j∈{1,2,...,n_d}，n_d表示计算节点的个数。

上述步骤2中，对测量数据进行分割的具体步骤如下：

步骤2.1、求出计算节点个数n_d的所有约数，并将这些约数按从小到大的顺序排列为一维数组n_div；

步骤2.2、找出一维数组n_div中最先大于或等于的数v_no，并令u_no＝n_d/v_no；

步骤2.3、对测量数据的uw投影的秩统计量R_uw进行归一化并量化后，得到uw投影的量化数据s_u：

s_u＝ceil(u_no×R_uw/l_uw)

其中，ceil表示向右取整函数，R_uw表示测量数据的uw投影的秩统计量，l_uw表示测量数据的uw投影的长度；

步骤2.4、将uw投影的量化数据s_u分别等于1,2,...,u_no的下标记录下来，并根据这些下标将测量数据的投影数据uw和vw，以及可见度数据y均分割成u_no份；

步骤2.5、对测量数据的vw投影的秩统计量R_vw进行归一化并量化后，得到vw投影的量化数据s_v：

s_v＝ceil(v_no×R_vw/l_vw)

其中，ceil表示向右取整函数，R_vw表示测量数据的vw投影的秩统计量，l_vw表示测量数据的vw投影的长度；

步骤2.6、将vw投影的量化数据s_v分别等于1,2,...,v_no的下标记录下来，并根据这些下标将步骤2.4所得的被分割成u_no份的测量数据的投影数据uw和vw，以及可见度数据y均再进一步分割成v_no份，最终得到n_d份数据块。

上述步骤3中，第j个计算节点在当前迭代t下的松弛变量为：

其中，表示第j个计算节点求解最小化问题的近似算子，G_j表示第j个计算节点的紧支撑核，M_j表示第j个计算节点的数据选取矩阵，F表示傅里叶变换矩阵，Z表示补零矩阵，x^(t-1)表示上次迭代t-1下的最终解，y_j表示第j个数据块的可见度数据，表示第j个计算节点在上次迭代t-1下的拉格朗日乘子。

上述步骤3中，第j个计算节点在当前迭代t下的拉格朗日乘子为：

其中，表示第j个计算节点在上次迭代t-1下的拉格朗日乘子，β表示常数，G_j表示第j个计算节点的紧支撑核，M_j表示第j个计算节点的数据选取矩阵，F表示傅里叶变换矩阵，Z表示补零矩阵，x^(t-1)表示上次迭代t-1下的最终解，y_j表示第j个数据块的可见度数据，表示第j个计算节点在当前迭代t下的松弛变量。

上述步骤3中，第j个计算节点在当前迭代t下的中间解为：

其中，G_j表示第j个计算节点的紧支撑核，M_j表示第j个计算节点的数据选取矩阵，F表示傅里叶变换矩阵，Z表示补零矩阵，x^(t-1)表示上次迭代t-1下的最终解，y_j表示第j个数据块的可见度数据，表示第j个计算节点在当前迭代t下的松弛变量，表示第j个计算节点在当前迭代t下的拉格朗日乘子，上标^T表示转置。

上述步骤4中，当前迭代t下的梯度下降解为：

其中，x^(t-1)表示上次迭代t-1下的最终解，ρ表示收敛常量，Z表示补零矩阵，F表示傅里叶矩阵，M_j表示第j个计算节点的数据选取矩阵，表示第j个计算节点在当前迭代t下的中间解，n_d表示计算节点的个数，上标^T表示转置。

与现有技术相比，本发明具有如下特点：

1、将射电干涉阵列图像重建的逆问题转换为等式约束条件下的凸优化问题，并使用近似交替方向乘子算法进行求解，每次迭代需要的总计算成本较低；

2、将数据分割成多份，能够精确快速的分割数据，并使用多计算节点进行计算，具有较高的并行性能，可拓展性强；

3、不需要任何的数据后处理步骤，直接得到重建图像，最终得到的图像的动态范围更大，对延展结构的天体显示更加清晰，相比传统的方法更加高效。

附图说明

图1为一种射电干涉阵列分布式图像重建方法的原理框图。

图2为几种方法图像重建结果对比图；其中(a)本发明；(b)CLEAN；(c)CS_CLEAN；(d)MS_CLEAN。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实例，并参照附图，对本发明进一步详细说明。

本实例使用的原始数据是由Karl G.Jansky甚大阵观测超新星遗迹3C391得到的数据再使用通用天文应用软件CASA进行校准之后的数据，数据格式为Measurement Set，文件名为3c391_ctm_mosaic_spw0.ms。

首先，查看文件3c391_ctm_mosaic_spw0.ms中的数据，选择目标射电源视场，即视场ID编号为0；使用CASA软件将目标射电源视场相关的观测数据从文件3c391_ctm_mosaic_spw0.ms中分离并创建1个子文件3c391_field0.ms。

然后，读取子文件3c391_field0.ms获取观测频率通道、uvw数据、可见度数据的实部和虚部以及可见度数据对应的权重。

接着，根据观测频率通道和偏振划分uvw数据和可见度数据及权重，分离出1个偏振的uvw数据和可见度数据及权重，然后再分离出1个频率通道的UVW数据和可见度数据及权重。将这些分离后的可见度数据转换为复数形式后，与分离后的uvw数据和可见度权重一起按列方式保存到一个文本文件中。

最后，读取上述得到的文本文件uvw数据中的uv数据，计算出uv值的最大范围bmax即的最大值，计算出uw和vw：vw＝v×π/(bmax×dl)，uw＝u×π/(bmax×dl)，其中dl是1个像素的大小；读取出文本文件中的可见度数据y_I及其权重weights，将y_I与weights相乘得到新的可见度数据y。

本发明所处理的数据是已经进行了校准的数据，因此不需要考虑方向依赖响应。

假设待重建的图像(N维实数信号)在某个基(如小波基)Ψ中可以用系数向量稀疏表示(有且仅有D＜＜N个非零系数)，即：x＝Ψα,||α||₀＝D＜＜N。

本发明通过将由可见度测量值(M维实数信号)重建图像的病态逆问题：y＝Φx+n，转换为求解下式的凸优化最小化问题的解重建图像x：

其中，Ψ^T是Ψ的转置；Φ是测量矩阵(或测量算子)，其表达式为Φ＝GFZ；这里表示内插核，用于去网格；F表示傅里叶变换算子；矩阵表示补零或零填充，用于解释过采样和缩放图像以补偿不完全的内插；l表示字典Ψ下的稀疏先验范数l＝||·||₁；函数f表示图像恢复解的真实性和正数性要求，表达式为：

h函数确保数据保真约束残差在由噪声水平∈定义的范数球内，表达式为：

这里∈是残差的范数的上限值。

矩阵(算子)表示n_b个稀疏基的集合。

为了实现分布式(并行)计算本发明将原始数据分割成n_d个数据块，有：

其中，G由紧支撑核和矩阵用于选取出数据分块j部分的离散傅里叶平面。

因此，对于每个数据分块有：

y_j＝Φ_jx+n_j

其中，n_j是关于y_j的噪声。

所以射电干涉阵列分布式图像重建方法通过求解下式的最小化问题：

其中，

l＝||·||₁，

∈_j表示第j个数据分块的噪声临界值。对于稀疏先验，范数是可分离的稀疏基可以表示为：

其中，i∈{1,...,n_b}。

本发明的主要原理是通过近似交替方向乘子算法求解上述最小化问题，其主要步骤如下：

(1)令因此上述最小化问题又可以改写为：其中，由于是指示函数的总和，我们可以重新定义为与

(2)将上述最小化问题转换为等式约束条件下的最小化问题：s.t.r＝Φx，其中，r是松弛变量。该等式的增广拉格朗日方程为：其中，向量s和常量μ均为拉格朗日乘子。

(3)上述寻找解的最小化问题最终转换为使用近似交替方向乘子算法迭代寻找增广拉格朗日方程的解：在每一次迭代中，近似交替方向乘子算法交替求解图像x的最小化问题：和松弛变量r的最小化：其中，拉格朗日乘子s使用步长为ρ的梯度上升法进行求解。根据函数表达式上述最小化问题又可以分为n_d个子问题：

(4)上式通过近似算子来计算得到，相当于计算在Φ_jx+s_j的近似算子：这个方法使得每个r_j接近Φ_jx，当x收敛于解时。收敛速度取决于拉格朗日乘子μ和与拉格朗日算子变量s最大化相联系的上升步长ρ。

(5)使用前向-后向方法求解x的最小化问题，即

(5.1)初始化：收敛常量其中是函数f(x)的近似算子，表达式为：P_C(z)＝max{z,0}。

(5.2)计算第m次迭代，n_b个次进程并行计算：每个进程计算得到的发送给主进程，其中i∈{1,...,n_b}；函数为软阈值算子，即阈值为的第m次迭代阈值算子，阈值算子的表达式为：

(5.3)主进程接收计算：

(5.4)重复2～3步骤，直至最后一次迭代或||Φ_jx-y_j||₂≤∈_j，返回最后一次迭代的

据此，本发明所设计的一种射电干涉阵列分布式图像重建方法，其具体包括步骤如下：

步骤1：中央节点初始化第1次迭代参数：初始图像x⁽¹⁾，第j数据块的松弛变量第j数据块的拉格朗日乘子第j数据块中间解常数β＝0.9，收敛常量其中表示Φ的谱范数，变量m；将初始化的参数广播并发送给计算节点。第j数据块对应第j计算节点。

在本实施例中，初始化参数：设迭代次数变量为t并初始化t＝1，初始化图像解x_sol为512×512的全零矩阵；增广拉格朗日乘子s、松弛变量r和q，它们与y数据类型相同，初始化值为均设为零；上升步数β＝0.9；根据uw和vw计算出栅格化矩阵G、傅里叶矩阵F和补零矩阵Z；计算收敛常量

步骤2：中央节点根据计算节点个数对数据y进行分块，并将分块后的数据y_j按照j∈{1,2,...,n_d}的顺序发送到n_d个计算节点。其中数据分割的具体方法如下：

步骤2.1、求出n_d的所有约数并按从小到大的顺序排列记为一维数组n_div。

步骤2.2、找出n_div中最先大于或等于的数并记为v_no，令u_no＝n_d/v_no。

步骤2.3、计算uw的秩统计量R_uw，并将秩统计量R_uw进行归一化并量化得到s_u。量化范围为：1～u_no，计算公式为：s_u＝ceil(y_no×R_uw/l_uw)，其中ceil是向右取整函数，l_uw为uw的长度。

步骤2.4、将s_u分别等于1,2,...,u_no的下标记录下来，uw、vw和y分别按照这些下标分割成u_no份。

步骤2.5、将步骤2.4中分割出的数据按照步骤2.3的计算方法计算出步骤2.4之后的vw的秩统计量R_vw，并将秩统计量R_vw进行归一化和量化得到s_v。量化范围为：1～v_no，计算公式为：s_v＝ceil(v_no×R_vw/l_vw)，其中ceil是向右取整函数，l_vw为vw的长度。

步骤2.6、将s_v分别等于1,2,...,v_no的下标记录下来，步骤2.4分割之后的每份数据即uw、vw和y按照这些下标再进一步分割成v_no份，最终得到n_d份数据。

在本实施例中，将投影坐标uw和vw和可见光数据y分割成n_d份，每份数据长度尽量相等，并发送到n_d个计算节点；计算每个数据块的噪声上限

步骤3：计算节点j顺序计算松弛变量拉格朗日乘子和中间解n_d个计算节点同时进行计算，各个计算节点得到的计算结果同时发送给中央节点。

其中，表示求解最小化问题的近似算子，G_j表示第j个计算节点的紧支撑核，F表示傅里叶变换矩阵，Z表示补零矩阵，M_j表示第j个计算节点的数据选取矩阵，y_j表示第j个数据块的可见度数据(即第j个计算节点的可见度数据)，x^(t-1)表示上次迭代t-1下的最终解，表示上次迭代t-1下的拉格朗日乘子，β表示常数，上标^T表示转置。

步骤4：中央节点计算：其中，x^(t-1)表示上次迭代t-1下的最终解，ρ表示收敛常量，Z表示补零矩阵，F表示傅里叶矩阵，M_j表示第j个计算节点的数据选取矩阵，表示当前迭代t下的中间解，n_d表示计算节点的个数，上标^T表示转置。

步骤5：中央节点使用双重前向-后向方法计算：

步骤5.1：初始化：收敛常量其中是函数f(x)的近似算子，表达式为：P_C(z)＝max{z,0}；其中Ψ在本实例中使用的是小波基。

步骤5.2：计算第m次迭代，n_b个次进程并行计算：每个进程计算得到的发送给主进程，其中i∈{1,...,n_b}；函数为软阈值算子，即阈值为的第m次迭代阈值算子，阈值算子的表达式为：

步骤5.3：主进程接收计算：

步骤5.4：重复5.2～5.3步骤，直至最后一次迭代或||Φ_jx-y_j||₂≤∈_j，返回最后一次迭代的并将赋值给x^(t)。

步骤6：重复步骤3～5，直至||Φ_jx^(t)-y_j||₂≤∈_j或最后一次迭代结束，最后得到的x^(t)即为重建的图像。

图2分别为本发明算法与CLEAN、CS-CLEAN、MS-CLEAN算法得到的图像重建结果，迭代次数均为100次。这些图像的动态范围分别是：10298.0000、65.1900、167.9851和131.2030。可以看出本发明算法的得到图像动态范围最高，图像质量最高，CLEAN算法的动态范围最小，图像质量最低，说明了CLEAN算法对延展结构的天体重建效果很差，表明本发明算法图像重建效果最好。

需要说明的是，尽管以上本发明所述的实施例是说明性的，但这并非是对本发明的限制，因此本发明并不局限于上述具体实施方式中。在不脱离本发明原理的情况下，凡是本领域技术人员在本发明的启示下获得的其它实施方式，均视为在本发明的保护之内。

Claims (6)

1.一种射电干涉阵列分布式图像重建方法，其特征是，包括步骤如下：

步骤1、中央节点将初始化参数广播并发送给计算节点；

步骤2、中央节点根据计算节点的个数n_d，将射电干涉仪的测量数据分割为n_d份；并将这n_d份数据块分别发送给n_d个计算节点；

步骤3、每个计算节点j根据所得到的数据块顺序计算当前迭代t下的松弛变量拉格朗日乘子和中间解并将所计算出的中间解发送给中央节点；

步骤4、中央节点根据各个计算节点返回的中间解计算当前迭代t下的梯度下降解

步骤5、中央节点根据所计算出的梯度下降解使用双重前向-后向方法计算当前迭代t下的最终解x^(t)；

步骤6、判断是否满足收敛条件，即当前迭代t达到最大迭代次数或者当前迭代t下的最终解x^(t)使得每个数据块均满足||Φ_jx^(t)-y_j||₂≤∈_j；

如果未满足收敛条件，则迭代次数t+1，并返回步骤3；

如果满足收敛条件，则当前迭代t下的最终解x^(t)即为重建的图像；

其中，Φ_j表示第j个数据块的测量矩阵，y_j表示第j个数据块的可见度数据，∈_j表示第j个数据块的噪声上限，j∈{1,2,...,n_d}，n_d表示计算节点的个数。

2.根据权利要求1所述一种射电干涉阵列分布式图像重建方法，其特征是，步骤2中，对测量数据进行分割的具体步骤如下：

步骤2.1、求出计算节点个数n_d的所有约数，并将这些约数按从小到大的顺序排列为一维数组n_div；

步骤2.2、找出一维数组n_div中最先大于或等于的数v_no，并令u_no＝n_d/v_no；

步骤2.3、对测量数据的uw投影的秩统计量R_uw进行归一化并量化后，得到uw投影的量化数据s_u：

s_u＝ceil(u_no×R_uw/l_uw)

其中，ceil表示向右取整函数，R_uw表示测量数据的uw投影的秩统计量，l_uw表示测量数据的uw投影的长度；

步骤2.4、将uw投影的量化数据s_u分别等于1,2,...,u_no的下标记录下来，并根据这些下标将测量数据的投影数据uw和vw，以及可见度数据y均分割成u_no份；

步骤2.5、对测量数据的vw投影的秩统计量R_vw进行归一化并量化后，得到vw投影的量化数据s_v：

s_v＝ceil(v_no×R_vw/l_vw)

其中，ceil表示向右取整函数，R_vw表示测量数据的vw投影的秩统计量，l_vw表示测量数据的vw投影的长度；

步骤2.6、将vw投影的量化数据s_v分别等于1,2,...,v_no的下标记录下来，并根据这些下标将步骤2.4所得的被分割成u_no份的测量数据的投影数据uw和vw，以及可见度数据y均再进一步分割成v_no份，最终得到n_d份数据块。

3.根据权利要求1所述一种射电干涉阵列分布式图像重建方法，其特征是，步骤3中，第第j个计算节点在当前迭代t下的松弛变量为：

<mrow> <msubsup> <mi>r</mi> <mi>j</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msub> <mi>P</mi> <msub> <mi>B</mi> <mi>j</mi> </msub> </msub> <mrow> <mo>(</mo> <msub> <mi>G</mi> <mi>j</mi> </msub> <msub> <mi>M</mi> <mi>j</mi> </msub> <msup> <mi>FZx</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>-</mo> <mn>1</mn> <mo>)</mo> </mrow> </msup> <mo>-</mo> <msub> <mi>y</mi> <mi>j</mi> </msub> <mo>+</mo> <msubsup> <mi>s</mi> <mi>j</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>-</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <mo>)</mo> </mrow> </mrow>

其中，表示第j个计算节点求解最小化问题的近似算子，G_j表示第j个计算节点的紧支撑核，M_j表示第j个计算节点的数据选取矩阵，F表示傅里叶变换矩阵，Z表示补零矩阵，x^(t-1)表示上次迭代t-1下的最终解，y_j表示第j个数据块的可见度数据，表示第j个计算节点在上次迭代t-1下的拉格朗日乘子，j∈{1,2,...,n_d}，n_d表示计算节点的个数。

4.根据权利要求1或3所述一种射电干涉阵列分布式图像重建方法，其特征是，步骤3中，第j个计算节点在当前迭代t下的拉格朗日乘子为：

<mrow> <msubsup> <mi>s</mi> <mi>j</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>s</mi> <mi>j</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>-</mo> <mn>1</mn> <mo>)</mo> </mrow> </msubsup> <mo>+</mo> <mi>&amp;beta;</mi> <mrow> <mo>(</mo> <msub> <mi>G</mi> <mi>j</mi> </msub> <msub> <mi>M</mi> <mi>j</mi> </msub> <msup> <mi>FZx</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>-</mo> <mn>1</mn> <mo>)</mo> </mrow> </msup> <mo>-</mo> <msub> <mi>y</mi> <mi>j</mi> </msub> <mo>-</mo> <msubsup> <mi>r</mi> <mi>j</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> </msubsup> <mo>)</mo> </mrow> </mrow>

其中，表示第j个计算节点在上次迭代t-1下的拉格朗日乘子，β表示常数，G_j表示第j个计算节点的紧支撑核，M_j表示第j个计算节点的数据选取矩阵，F表示傅里叶变换矩阵，Z表示补零矩阵，x^(t-1)表示上次迭代t-1下的最终解，y_j表示第j个数据块的可见度数据，表示第j个计算节点在当前迭代t下的松弛变量，j∈{1,2,...,n_d}，n_d表示计算节点的个数。

5.根据权利要求1或4所述一种射电干涉阵列分布式图像重建方法，其特征是，步骤3中，第j个计算节点在当前迭代t下的中间解为：

<mrow> <msubsup> <mi>q</mi> <mi>j</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> </msubsup> <mo>=</mo> <msubsup> <mi>G</mi> <mi>j</mi> <mi>T</mi> </msubsup> <mrow> <mo>(</mo> <msub> <mi>G</mi> <mi>j</mi> </msub> <msub> <mi>M</mi> <mi>j</mi> </msub> <msup> <mi>FZx</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>-</mo> <mn>1</mn> <mo>)</mo> </mrow> </msup> <mo>-</mo> <msub> <mi>y</mi> <mi>j</mi> </msub> <mo>+</mo> <msubsup> <mi>r</mi> <mi>j</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> </msubsup> <mo>-</mo> <msubsup> <mi>s</mi> <mi>j</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> </msubsup> <mo>)</mo> </mrow> </mrow>

其中，G_j表示第j个计算节点的紧支撑核，M_j表示第j个计算节点的数据选取矩阵，F表示傅里叶变换矩阵，Z表示补零矩阵，x^(t-1)表示上次迭代t-1下的最终解，y_j表示第j个数据块的可见度数据，表示第j个计算节点在当前迭代t下的松弛变量，表示第j个计算节点在当前迭代t下的拉格朗日乘子，上标^T表示转置，j∈{1,2,...,n_d}，n_d表示计算节点的个数。

6.根据权利要求1所述一种射电干涉阵列分布式图像重建方法，其特征是，步骤4中，当前迭代t下的梯度下降解为：

<mrow> <msup> <mover> <mi>x</mi> <mo>~</mo> </mover> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> </msup> <mo>=</mo> <msup> <mi>x</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>-</mo> <mn>1</mn> <mo>)</mo> </mrow> </msup> <mo>-</mo> <msup> <mi>&amp;rho;Z</mi> <mi>T</mi> </msup> <msup> <mi>F</mi> <mi>T</mi> </msup> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>j</mi> <mo>=</mo> <mn>1</mn> </mrow> <msub> <mi>n</mi> <mi>d</mi> </msub> </munderover> <msubsup> <mi>M</mi> <mi>j</mi> <mi>T</mi> </msubsup> <msubsup> <mi>q</mi> <mi>j</mi> <mrow> <mo>(</mo> <mi>t</mi> <mo>)</mo> </mrow> </msubsup> </mrow>

其中，x^(t-1)表示上次迭代t-1下的最终解，ρ表示收敛常量，Z表示补零矩阵，F表示傅里叶矩阵，M_j表示第j个计算节点的数据选取矩阵，表示第j个计算节点在当前迭代t下的中间解，n_d表示计算节点的个数，上标^T表示转置。