一种用于遥感影像处理的时空邻域高斯加权中值滤波方法技术

本发明专利技术公开了一种用于遥感影像处理的时空邻域高斯加权中值滤波方法，根据不同时间序列SAR影像的观测值同时与其时间邻域和空间邻域的观测值高度相关，且时间间隔越大、空间距离越大相关性越小，提出根据其时间和空间距离进行高斯加权的中值滤波算法，并根据SAR影像的统计分布特征，预估中值所在的区间，缩小中值搜索范围，加快中值求解速度，本发明专利技术能较好保持影像上地物细节特征，抑制斑点噪声，同时较大程度去除时序影像中车、船等临时地物的影响，解决了时序影像斑点噪声抑制问题。

【技术实现步骤摘要】

本专利技术属于遥感影像处理
,是一种合成孔径雷达(Synthetic Aperture Radar,SAR)影像相干斑噪声抑制方法，涉及一种对重轨时间序列SAR影像中的相干斑噪声进行滤除的新方法，具体涉及一种应用于时序SAR影像斑点噪声抑制的时空邻域高斯加权中值滤波快速算法。
技术介绍
合成孔径雷达不受天气地理和时间等因素的限制，能够对地面进行高分辨率成像，并且具有一定的穿透力，因而被广泛应用于军事侦察、资源探测、环境监测、测绘制图、地理国情监测等对地遥感应用中，随着多颗雷达卫星的成功发射，时序雷达影像成为当前遥感应用的一大重要方面。由于SAR基于相干成像机理进行工作，导致SAR影像中存在严重的斑点噪声，对SAR图像的解译判读和信息提取带来非常大的影响。早期的SAR影像斑点噪声抑制是通过空间多视平均处理来实现的，该方法的主要缺点是牺牲了影像的空间分辨率，不适用于当前的高分辨率SAR系统。目前主要采用空域滤波方法来实现斑点噪声抑制，即一般是利用一个滑动窗口，对窗口内的像素进行邻域处理得到窗口中心点的像素值。而中值滤波是一种非线性滤波算法，能较好地抑制斑点噪声、椒盐噪声等非高斯噪声，保持影像细节信息，在SAR影像中具有广泛的应用。现有技术中，中值滤波仅运用在空间邻域，部分算法引入了不同形状和大小的邻域窗口，并采用高斯等方法进行了加权，但此类方法的不足是：没有充分利用时间序列观测值之间的相关性，且计算效率较低。其它一些考虑时间域相关性的滤波方法，多对各时相同等对待，未充分顾及时间间隔越大其相关性越低的问题，且不能克服临时性车、船等着某一时相偶然出现地物的影响，因此滤波结果对后续分类解译效果的提升依然有限。这些改进算法均未能彻底解决时序SAR影像中斑点噪声滤除的问题。本专利技术即针对现有技术的不足而研究提出。
技术实现思路
本专利技术一种时序SAR时空邻域高斯加权中值滤波斑噪抑制快速算法，充分利用SAR时间序列影像中时间邻域和空间邻域观测值的相关性，发挥中值滤波算法对斑点噪声等非高斯噪声抑制和细节信息保持的良好性能，在抑制斑点噪声的同时，抑制时序影像中临时地物的影响，有效解决现有方法未能很好解决时序SAR影像中斑点噪声快速滤波的问题。为解决上述技术问题，本专利技术一种时序SAR时空邻域高斯加权中值滤波斑噪抑制快速算法，时序影像的表示形式为强度影像，宽和高分别为nw和nh个像素，共有nt个时相，其中第k个时相的观测时刻记为Tk，通过滑动窗口遍历每个像素(i,j,k)，其中(i,j)表示影像的行列号，记滤波前像素(i,j,k)的强度影像为P(i,j,k)，对其时空邻域进行处理即得到像素(i,j,k)的滤波结果包括以下具体步骤：步骤1：设置时空邻域窗口参数：即将空间邻域窗口的大小设置为Ws，时间邻域窗口的大小为Wt，由此得到通过滑动窗口遍历每个像素(i,j,k)时的时空邻域大小为Ws×Ws×Wt,满足：51≥Ws≥3，且Ws为奇数；且Wt至少为最大相邻时相观测时间间隔的3倍。步骤2：设置时空邻域高斯加权函数的参数：即空间邻域加权参数σs>0，且2σs<Ws/2；时间邻域加权参数σt>0，且2σt<Wt/2。步骤3：初始化当前处理影像的时相，k＝1。步骤4：计算当前时相滑动窗口的高斯权矩阵。步骤4.1：计算参与当前时相滤波的相邻时相，按照时间先后顺序遍历所有时相的观测时间Tt，当时相t与当前时相k的时间间隔Tt-Tk≤Wt/2时，将时相t加入时相k时间邻域集合，并记时相k的时间邻域集合为Snt,k，Snt,k的大小即为时间邻域，时相个数为Wtk。步骤:4.2：计算当前处理时相k的时空邻域滑动窗口内的三维权矩阵各元素的值w(dx,dy,dt),满足： w ( d x , d y , d t ) = exp ( - d x 2 + d y 2 2 σ s 2 ) exp ( - d t 2 2 σ t 2 ) ; ]]>其中dx,dy分别是方形模板上(x,y)处到模板中心的行方向和列方向的距离，其中x,y＝1,2,...,Ws；dt为时相k的邻近时相t的观测时刻到时相k的时间间隔dt＝Tt-Tk，且Tt∈Snt,k，t＝1,2,...,Wtk。由此得到大小为Ws×Ws×Wtk的权矩阵W0。步骤4.3：计算三维权矩阵元素的总和，即总权值wall满足：其中w(i,j,k)是三维权矩阵W0第i行第j列第k个时相的权重。步骤4.4：计算归一化权矩阵W,满足：步骤5：初始化当前处理像素(ic,jc)为第一个像素，即(ic,jc)＝(0,0)，并将滑动窗口的中心与(ic,jc)对齐，开始处理第一行第一列像素。步骤6：根据时空邻域内最大值和最小值确定，并缩小的中值搜索范围。步骤6.1：记当前处理像素时空邻域Ws×Ws×Wtk个影像强度值组成的集合为SN，计算SN的最大值Pmax,c和最小值Pmin,c。步骤6.2：计算当前处理像素时空邻域内的影像幅度值的最大值和最小值步骤6.3：由于SAR幅度影像的视数越大，越接近正态分布，其中值与均值越接近，当nL>100，幅度影像近似服从正态分布，其中值在均值附近的概率极大；随着视数降低，其中值逐渐向均值的左边移动。且通过实验数据统计表明，缩小后的中值搜索范围可根据不同的视数nL进行计算： P 本文档来自技高网...

【技术保护点】
时序SAR时空邻域高斯加权中值滤波斑噪抑制快速算法,时序影像的表示形式为强度影像，宽和高分别为nw和nh个像素，共有nt个时相，其中第k个时相的观测时刻记为Tk，通过滑动窗口遍历每个像素(i,j,k)，其中(i,j)表示影像的行列号，记滤波前像素(i,j,k)的强度影像为P(i,j,k)，对其时空邻域进行处理即得到像素(i,j,k)的滤波结果其特征在于包括以下步骤：步骤1：设置时空邻域窗口参数：即将空间邻域窗口的大小设置为Ws，时间邻域窗口的大小为Wt，由此得到通过滑动窗口遍历每个像素(i,j,k)时的时空邻域大小为Ws×Ws×Wt,满足：51≥Ws≥3，且Ws为奇数；且Wt至少为最大相邻时相观测时间间隔的3倍；步骤2：设置时空邻域高斯加权函数的参数：即空间邻域加权参数σs>0，且2σs<Ws/2；时间邻域加权参数σt>0，且2σt<Wt/2；步骤3：初始化当前处理影像的时相，k＝1；步骤4：计算当前时相滑动窗口的高斯权矩阵；步骤5：初始化当前处理像素(ic,jc)为第一个像素，即(ic,jc)＝(0,0)，并将滑动窗口的中心与(ic,jc)对齐，开始处理第一行第一列像素；步骤6：根据时空邻域内最大值和最小值确定，并缩小的中值搜索范围；步骤7：判断中值是否位于缩小的搜索范围步骤8：在步骤7确定的搜索区间J＝(Jmin,Jmax)内进行二分法递归快速搜索高斯时空加权的中值；步骤9：保存结果步骤10：当前时相像素全部处理完成，进入步骤11；否则移动滑动窗口中心到下一个像素，更新当前处理像素(ic,jc)，进入步骤6；步骤11：当时相k影像的全部像素都处理完时，更新当前处理时相k＝k+1；如果k≤nt进入步骤4，开始对下一时相进行处理，否则处理结束。...

【技术特征摘要】
1.时序SAR时空邻域高斯加权中值滤波斑噪抑制快速算法,时序影像的表示形式为强度影像，宽和高分别为nw和nh个像素，共有nt个时相，其中第k个时相的观测时刻记为Tk，通过滑动窗口遍历每个像素(i,j,k)，其中(i,j)表示影像的行列号，记滤波前像素(i,j,k)的强度影像为P(i,j,k)，对其时空邻域进行处理即得到像素(i,j,k)的滤波结果其特征在于包括以下步骤：步骤1：设置时空邻域窗口参数：即将空间邻域窗口的大小设置为Ws，时间邻域窗口的大小为Wt，由此得到通过滑动窗口遍历每个像素(i,j,k)时的时空邻域大小为Ws×Ws×Wt,满足：51≥Ws≥3，且Ws为奇数；且Wt至少为最大相邻时相观测时间间隔的3倍；步骤2：设置时空邻域高斯加权函数的参数：即空间邻域加权参数σs>0，且2σs<Ws/2；时间邻域加权参数σt>0，且2σt<Wt/2；步骤3：初始化当前处理影像的时相，k＝1；步骤4：计算当前时相滑动窗口的高斯权矩阵；步骤5：初始化当前处理像素(ic,jc)为第一个像素，即(ic,jc)＝(0,0)，并将滑动窗口的中心与(ic,jc)对齐，开始处理第一行第一列像素；步骤6：根据时空邻域内最大值和最小值确定，并缩小的中值搜索范围；步骤7：判断中值是否位于缩小的搜索范围步骤8：在步骤7确定的搜索区间J＝(Jmin,Jmax)内进行二分法递归快速搜索高斯时空加权的中值；步骤9：保存结果步骤10：当前时相像素全部处理完成，进入步骤11；否则移动滑动窗口中心到下一个像素，更新当前处理像素(ic,jc)，进入步骤6；步骤11：当时相k影像的全部像素都处理完时，更新当前处理时相k＝k+1；如果k≤nt进入步骤4，开始对下一时相进行处理，否则处理结束。2.根据权利要求1所述时序SAR时空邻域高斯加权中值滤波斑噪抑制快速算法，其特征在于所述步骤4还包括：步骤4.1：计算参与当前时相滤波的相邻时相，按照时间先后顺序遍历所有时相的观测时间Tt，当时相t与当前时相k的时间间隔Tt-Tk≤Wt/2时，将时相t加入时相k时间邻域集合，并记时相k的时间邻域集合为Snt,k，Snt,k的大小即为时间邻域，时相个数为Wtk；步骤:4.2：计算当前处理时相k的时空邻域滑动窗口内的三维权矩阵各元素的值w(dx,dy,dt),满足： w ( d x , d y , d t ) = exp ( - d x 2 + d y 2 2 σ s 2 ) exp ( - d t 2 2 σ t 2 ) ; ]]>其中dx,dy分别是方形模板上(x,y)处到模板中心的行方向和列方向的距离，其中x,y＝1,2,...,Ws；dt为时相k的邻近时相t的观测时刻到时相k的时间间隔dt＝Tt-Tk，且Tt∈Snt,k，t＝1,2,...,Wtk，由此得到大小为Ws×Ws×Wtk的权矩阵W0；步骤4.3：计算三维权矩阵元素的总和，即总权值wall满足：其中w(i,j,k)是三维权矩阵W0第i行第j列第k个时相的权重；步骤4.4：计算归一化权矩阵W,满足：3.根据权利要求1所述时序SAR时空邻域高斯加权中值滤波斑噪抑制快速算法，其特征在于所述步骤6包括：步骤6.1：记当前处理像素时空邻域Ws×Ws×Wtk个影像强度值组成的集合为SN，计算SN的最大值Pmax,c和最小值Pmin,c；步骤6.2：计算当前处理像素时空邻域内的影像幅度值的最大值和最小值步骤6.3：由于SAR幅度影像的视数越大，越接近正态分布，其中值与均值越接近，当nL>100，幅度影像近似服从正态分布，其中值在均值附近的概率极大；随着视数降低，其中值逐渐向均值的左边移动，且通过实验数据统计表明，缩小后的中值搜索范围可根据不同的视数nL进行计算： P ^ min , c = ( A min , c + 0.15 · ( A max , c - A min , c ) ) 2 P ^ max , c = ( A min , c + 0.32 · ( A max , c - A min , c ) ) 2 , i f n L ≤ 5 , ]]> P ^ min , c = ( A min , c + 0.25 · ( A max , c - A min , c ) ) 2 P ^ max , c = ( A min , c + 0.40 · ( A max , c - A min , c ) ) 2 , i f 5 < n L ≤ 10 , ]]> P ^ min , c = ( A min , c + 0.22 · ( A max , c - A min , c ) ) 2 P ^ max , c = ( A min ...

【专利技术属性】
技术研发人员：邓少平，甘宗平，孙盛，李胜，吴泽洪，方志民，刘学林，
申请(专利权)人：甘宗平，
类型：发明
国别省市：广东;44