CN110954303B - 一种基于高分辨率遥感图像基准的mtf自动测算方法 - Google Patents

一种基于高分辨率遥感图像基准的mtf自动测算方法 Download PDF

Info

Publication number
CN110954303B
CN110954303B CN201911142854.2A CN201911142854A CN110954303B CN 110954303 B CN110954303 B CN 110954303B CN 201911142854 A CN201911142854 A CN 201911142854A CN 110954303 B CN110954303 B CN 110954303B
Authority
CN
China
Prior art keywords
image
mtf
resolution
remote sensing
calculating
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.)
Active
Application number
CN201911142854.2A
Other languages
English (en)
Other versions
CN110954303A (zh
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.)
Southeast University
Original Assignee
Southeast 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 Southeast University filed Critical Southeast University
Priority to CN201911142854.2A priority Critical patent/CN110954303B/zh
Publication of CN110954303A publication Critical patent/CN110954303A/zh
Application granted granted Critical
Publication of CN110954303B publication Critical patent/CN110954303B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01MTESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
    • G01M11/00Testing of optical apparatus; Testing structures by optical methods not otherwise provided for
    • G01M11/02Testing optical properties
    • G01M11/0292Testing optical properties of objectives by measuring the optical modulation transfer function
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00Image enhancement or restoration
    • G06T5/10Image enhancement or restoration using non-spatial domain filtering
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10032Satellite or aerial image; Remote sensing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20024Filtering details
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20048Transform domain processing
    • G06T2207/20056Discrete and fast Fourier transform, [DFT, FFT]

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Chemical & Material Sciences (AREA)
  • Analytical Chemistry (AREA)
  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Image Analysis (AREA)

Abstract

本发明提供一种基于高分辨率遥感图像基准的MTF自动测算方法,主要解决传统方法在遥感图像在轨动态MTF检测中,图像中必须包含特定冲激响应的问题,以及在此过程中人工参与造成的效率和精度问题。本发明的核心机制是,采用与待检测图像同一成像区域的高分辨率遥感图像,依据简化的光学成像模型对高分辨率图像进行退化,通过构造一种健壮的退化图像和低分辨率图像之间分辨率相似度函数为评价函数,在参数空间中搜索图像退化的最优参数,从而确定待检测图像的MTF值。本发明的精确度、健壮度高,不受图像间平移偏差影响,受波段间光谱差异影响较小,运算速度快,可达到实用水平。适用于遥感及医学影像的MTF检测及空间尺度匹配。

Description

一种基于高分辨率遥感图像基准的MTF自动测算方法
技术领域
本发明涉及属于航天光学遥感成像及地面处理系统和数字图像处理技术领域,特别是涉及一种基于高分辨率遥感图像基准的MTF自动测算方法。
背景技术
卫星遥感技术可以对地球表面进行大规模全覆盖的观测,在测绘、气象、海洋、农业、自然资源调查、灾害监测、国防安全等诸多方面均发挥巨大作用。由于光学遥感成像系统在遥感图像成像过程中,会由于光学系统特性、探测器积分同步、传感器平台运动、震颤、大气状况、杂光等因素造成图像模糊,这些因素合并起来对冲击输入的空间响应称为点扩散函数(Point Spread Function,PSF),其在频率域的幅度称为调制传递函数(ModulationTransfer Function,MTF).从影响因素来看,光学系统特性、传感器探测器积分等因素构成了静态MTF,合并其他因素则构成了卫星在轨动态MTF.MTF决定了遥感图像的实际空间分辨率(effective spatial resolution),直接影响遥感图像增强、融合与同化、多尺度分析、尺度转化、目标提取等任务的精度。
静态MTF在卫星发射前进行测量,但卫星在太空运行后,参数可能会发生改变;随着服役时间增加,各项性能也可能发生衰退。震颤、大气状况等影响则无法预测,因此除了建设使用地面标准测试场或对特定区域的观测,在轨动态MTF测量的主流方法是假定在光学图像中存在一些对容易建模的“尖锐”信号的响应,通过对这部分信号的处理来提取MTF,如刀刃法对图像中分布接近顺轨和跨轨方向的强烈边缘进行提取,在与边缘垂直方向上进行微分获取线扩散函数(Line Spread Function,LSF),对LSF进行微分操作获取PSF,然后经过傅里叶变换提取MTF.这类方法对图像内容有一定的要求,同时需要人工参与进行特征选取,具有一定的主观性,影响了MTF测量的自动化程度和有效性。
卫星数据提供商提供的数据产品,一般会对全色数据进行MTF补偿,但绝大多数遥感成像系统的MTF补偿方法都没有公开。因此,在处理多光谱/高光谱数据时,利用较待处理图像具有更高空间分辨率的遥感图像作为参考,特别是同平台经过了MTF补偿的全色数据,对待处理图像的MTF进行测量,是一条值得探索的思路。如何在MTF测量中避免冲击法、刃边法等对图像内容的要求,同时避免人工干预,是一个需要解决的问题。
发明内容
为了现有方法测量在轨动态MTF时,对图像内容的限制问题,以及如何避免人工干预进行全自动MTF测量的问题,本发明提供一种基于高分辨率遥感图像基准的MTF自动测算方法,包括下述步骤,其特征在于:
步骤1:数据预处理,包括高分辨率图像Ih和待检测图像Il对齐、裁剪共同覆盖区域;
步骤2:参数初始化,设置高、低分辨率图像的名义分辨率比值r,初始搜索步长S,终止阈值Sth,退化参数m;
步骤3:对待检测图像Tl施行FFT,获取其幅度谱
Figure GDA0002945669270000021
步骤4:采用高斯模型根据m构造相应的低通滤波器G,并利用G对Ih进行低通滤波;
步骤5:将步骤4获得的图像使用盒采样方法向下r倍重采样至Il尺寸,获取Ih的退化图像Id
步骤6:计算Il相对于Id的MTF值,包括以下步骤,
步骤6.1.初始化相对光学MTF值mr
步骤6.2.对Id施行FFT,获取其幅度谱
Figure GDA0002945669270000022
步骤6.3.采用高斯模型根据mr构造相应的频率域滤波器H;
步骤6.4.利用梯度下降法计算令分辨率相似度SSR(Il,HId,1-H)达到最大的mr值,迭代终止条件为S<Sth
步骤6.5.根据mr和m更新m值为m←m·mr
步骤7:调整m值,m更新如下m←m+0.05,重复步骤4-6一次后转入下一步骤,将频谱混叠对测量结果的影响降至可接受的程度,从而得到对Il光学MTF更准确的估计;
步骤8:根据m计算待检测图像Il的光学MTF和总MTF并输出结果。
作为本发明进一步改进,步骤1对高分辨率图像Ih和待检测图像Il进行对齐,主要是方向对齐,使二者没有角度偏转;共同覆盖区域无需精确重叠,裁剪图像的尺寸比例与分辨率比例相同即可;需尽量避开图像边缘空值区域。
作为本发明进一步改进,步骤2中采用高、低分辨率图像的名义分辨率比值r,如按图像尺寸比值计算,需确保所处理图像未经过图像缩放处理;退化参数m∈[0.7,1)用于搜索Il的光学MTF,推荐初始值为0.8。
作为本发明进一步改进,步骤3和步骤6.2中对待检测图像Il和退化图像Id分别施行FFT,仅提取其幅度谱参与后续运算,不保留相位谱,以提高鲁棒性,并且减少运算量。
作为本发明进一步改进,步骤4对Ih进行低通滤波以减弱降采样中频谱混叠对测量的影响,滤波器G定义为
Figure GDA0002945669270000023
σs为空间域高斯滤波器的尺度参数。
Figure GDA0002945669270000031
滤波过程为
Figure GDA0002945669270000032
作为本发明进一步改进,步骤5采用盒采样方法对步骤4获得的图像向下r倍重采样至Il尺寸,获取Ih的退化图像Id,其效果为合并模拟探测器MTF退化操作和降采样操作,降低算法复杂度,提高效率。
作为本发明进一步改进,步骤6采用了一种计算图像间分辨率相似度指标,该指标定义如下:对于相同尺寸N的图像Ix和Iy,幅度谱分别为X和Y,在给定的频率域滤波器H下,二者的分辨率相似度SSR定义为
Figure GDA0002945669270000033
其中Xk,Yk,Hk分别为X,Y,H的第k个元素。
作为本发明进一步改进,步骤6在低分辨率图像尺度上进行MTF测量,可以避免上述采样方法对测量结果的影响,并减少运算量,提高运算速度。
作为本发明进一步改进,步骤6.4计算令分辨率相似度SSR(Il,HId,1-H)达到最大的mr值,作为对Il和Id的光学MTF值比值的估计。
作为本发明进一步改进,将MTF分解为光学MTF和探测器MTF两部分,截断频率处光学MTF部分Mopt估计为
Mopt=m
探测器MTF部分估计为
Figure GDA0002945669270000034
则截断频率处的总MTF为二者乘积
MTF=MoptMdet
本发明与现有技术相比,改进主要体现在:首先,消除了MTF测量方法对检测图像内容的依赖,可以应用于任意图像,提高了MTF测量的时效性;其次,本方法无需人工参与,避免了主观性和不确定性;再次,本方法提出了一种在频率域中计算的分辨率相似度指数,和现有的基于图像内容相关系数相比,不受配准平移偏差的影响,受波段间光谱差异影响小,测量结果的稳定性好,同时运算速度快;最后,本方法将MTF分解为包含大气等影响的光学MTF部分和探测器MTF部分,对二者分别建模,提高了MTF测量的精度。
具体实施方式
下面结合具体实施方式对本发明作进一步详细描述:
本发明提供一种基于高分辨率遥感图像基准的MTF自动测算方法,本发明采用高分辨率遥感图像作为基准,将动态MTF简化为包含了大气影响的光学MTF和探测器MTF两部分,通过不同参数设置进行图像退化;提出一种分辨率相似度度量函数,用于评价不同参数下的退化图像与待检测图像之间的分辨率相似度,利用梯度下降法搜索相似度达到最大的参数设置,从而计算对应的MTF值。
作为本发明具体实施例,本发明采用真实WorldView 2星载遥感多光谱和全色光图像,利用全色波段图像(PAN)作为基准,对多光谱波段2图像进行MTF测量。实施本发明包括以下步骤:
步骤一:数据预处理:由于采用同平台同时相的PAN波段图像作为高分辨率基准图像,因此无需对齐,直接按分辨率比例(r=4)进行裁边即可,裁边宽度100像元;
步骤二:参数初始化:设置高、低分辨率图像的名义分辨率比值r,初始搜索步长S,终止阈值S_th,退化参数m,MATLAB代码如下:
r=4;
S=0.05;
S_th=0.001;
m=0.8
步骤三:对MS波段图像进行FFT,获取功率谱。MATLAB代码如下:
Fa_ms=abs(fft2(I_MS);
步骤四:构造低通滤波器对PAN波段图像进行低通滤波,并采用盒采样方法向下4倍重采样至MS图像尺寸,获取退化图像I_d。MATLAB代码示意如下:
N=41;
alpha=sqrt((N*(0.5/r))^2/(-2*log(m)));
H=fspecial('gaussian',N,alpha);
Hd=H./max(H(:));
h=fwind1(Hd,kaiser(N));
I_d=double(imfilter(I_PAN,real(h),'circular'));
I_d=imresize(I_d,1/r,'box');
步骤五:计算I_MS相对于I_d的MTF比值m_r,初始估计值为m_r=0.5.
步骤六:提取I_d的幅度谱Fa_d如下:
Fa_d=abs(fft2(I_d);
步骤七:根据m_r构造频率域滤波器H.MATLAB代码示意如下:
M,N=size(I_d);
sigma_x=N*(0.5/r)./sqrt(-2*log(m_r));
sigma_y=M*(0.5/r)./sqrt(-2*log(m_r));
H_x=fspecial('gaussian',[1,N],sigma_x);
H_y=fspecial('gaussian',[M,1],sigma_y);
H=H_y*H_x;
步骤八:计算SSR值,并通过梯度下降法求得最大值对应的m_r值。SSR计算MATLAB代码如下:
H=ifftshift(H);
fx=Fa_d(:).*H(:).*(1-H(:));
fy=Fa_ms(:).*(1-H(:));
fx(1)=0;
fy(1)=0;
nsum_xy=real(fx*fy');
nsum_xx=real(fx*fx');
nsum_yy=real(fy*fy');
ssr=(nsum_xy)/(sqrt(nsum_xx)*sqrt(nsum_yy));
步骤九:更新m_r。更新方法如下:
m_r=m_r*m+0.05;
步骤十:重复步骤七和步骤八一次;
步骤十一:计算传感器MTF并与m_r合并,输出总MTF。MATLAB代码示意如下:
MTF_opt=m_r;
MTF_det=sinc(0.5);
MTF=MTF_opt*MTF_det;
disp(MTF)。
以上所述,仅是本发明的较佳实施例而已,并非是对本发明作任何其他形式的限制,而依据本发明的技术实质所作的任何修改或等同变化,仍属于本发明所要求保护的范围。

Claims (10)

1.一种基于高分辨率遥感图像基准的MTF自动测算方法,包括下述步骤,其特征在于:
步骤1:数据预处理,包括高分辨率图像Ih和待检测图像Il对齐、裁剪共同覆盖区域;
步骤2:参数初始化,设置高、低分辨率图像的名义分辨率比值r,初始搜索步长S,终止阈值Sth,退化参数m;
步骤3:对待检测图像Il施行FFT,获取其幅度谱
Figure FDA0002945669260000011
步骤4:采用高斯模型根据m构造相应的低通滤波器G,并利用G对Ih进行低通滤波;
步骤5:将步骤4获得的图像使用盒采样方法向下r倍重采样至Il尺寸,获取Ih的退化图像Id
步骤6:计算Il相对于Id的MTF值,包括以下步骤,
步骤6.1.初始化相对光学MTF值mr
步骤6.2.对Id施行FFT,获取其幅度谱
Figure FDA0002945669260000012
步骤6.3.采用高斯模型根据mr构造相应的频率域滤波器H;
步骤6.4.利用梯度下降法计算令分辨率相似度SSR(Il,HId,1-H)达到最大的mr值,迭代终止条件为S<Sth
步骤6.5.根据mr和m更新m值为m←m·mr
步骤7:调整m值,m更新如下m←m+0.05,重复步骤4-6一次后转入下一步骤,将频谱混叠对测量结果的影响降至可接受的程度,从而得到对Il光学MTF更准确的估计;
步骤8:根据m计算待检测图像Il的光学MTF和总MTF并输出结果。
2.根据权利要求1所述的一种基于高分辨率遥感图像基准的MTF自动测算方法,其特征在于:步骤1对高分辨率图像Ih和待检测图像Il进行对齐,主要是方向对齐,使二者没有角度偏转;共同覆盖区域无需精确重叠,裁剪图像的尺寸比例与分辨率比例相同即可;需尽量避开图像边缘空值区域。
3.根据权利要求1所述的一种基于高分辨率遥感图像基准的MTF自动测算方法,其特征在于:步骤2中采用高、低分辨率图像的名义分辨率比值r,如按图像尺寸比值计算,需确保所处理图像未经过图像缩放处理;退化参数m∈[0.7,1)用于搜索Il的光学MTF,推荐初始值为0.8。
4.根据权利要求1所述的一种基于高分辨率遥感图像基准的MTF自动测算方法,其特征在于:步骤3和步骤6.2中对待检测图像Il和退化图像Id分别施行FFT,仅提取其幅度谱参与后续运算,不保留相位谱,以提高鲁棒性,并且减少运算量。
5.根据权利要求1所述的一种基于高分辨率遥感图像基准的MTF自动测算方法,其特征在于:步骤4对Ih进行低通滤波以减弱降采样中频谱混叠对测量的影响,低通滤波器G定义为
Figure FDA0002945669260000021
σs为空间域高斯滤波器的尺度参数;
Figure FDA0002945669260000022
滤波过程为
Figure FDA0002945669260000023
6.根据权利要求1所述的一种基于高分辨率遥感图像基准的MTF自动测算方法,其特征在于:步骤5采用盒采样方法对步骤4获得的图像向下r倍重采样至Il尺寸,获取Ih的退化图像Id,其效果为合并模拟探测器MTF退化操作和降采样操作,降低算法复杂度,提高效率。
7.根据权利要求1所述的一种基于高分辨率遥感图像基准的MTF自动测算方法,其特征在于:步骤6采用了一种计算图像间分辨率相似度指标,该指标定义如下:对于相同尺寸N的图像Ix和Iy,幅度谱分别为X和Y,在给定的频率域滤波器H下,二者的分辨率相似度SSR定义为
Figure FDA0002945669260000024
其中Xk,Yk,Hk分别为X,Y,H的第k个元素。
8.根据权利要求1所述的一种基于高分辨率遥感图像基准的MTF自动测算方法,其特征在于:步骤6在低分辨率图像尺度上进行MTF测量,可以避免上述采样方法对测量结果的影响,并减少运算量,提高运算速度。
9.根据权利要求1所述的一种基于高分辨率遥感图像基准的MTF自动测算方法,其特征在于:步骤6.4计算令分辨率相似度SSR(Il,HId,1-H)达到最大的mr值,作为对Il和Id的光学MTF值比值的估计。
10.根据权利要求1所述的一种基于高分辨率遥感图像基准的MTF自动测算方法,其特征在于:将MTF分解为光学MTF和探测器MTF两部分,截断频率处光学MTF部分Mopt估计为Mopt=m
探测器MTF部分估计为
Figure FDA0002945669260000025
则截断频率处的总MTF为二者乘积
MTF=MoptMdet
CN201911142854.2A 2019-11-20 2019-11-20 一种基于高分辨率遥感图像基准的mtf自动测算方法 Active CN110954303B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911142854.2A CN110954303B (zh) 2019-11-20 2019-11-20 一种基于高分辨率遥感图像基准的mtf自动测算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911142854.2A CN110954303B (zh) 2019-11-20 2019-11-20 一种基于高分辨率遥感图像基准的mtf自动测算方法

Publications (2)

Publication Number Publication Date
CN110954303A CN110954303A (zh) 2020-04-03
CN110954303B true CN110954303B (zh) 2021-05-18

Family

ID=69978036

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911142854.2A Active CN110954303B (zh) 2019-11-20 2019-11-20 一种基于高分辨率遥感图像基准的mtf自动测算方法

Country Status (1)

Country Link
CN (1) CN110954303B (zh)

Family Cites Families (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101441764B (zh) * 2008-12-31 2011-01-26 中国资源卫星应用中心 一种mtfc遥感图像复原方法
CN101793599A (zh) * 2010-03-29 2010-08-04 中国科学院对地观测与数字地球科学中心 非理想标靶情况下mtf参数测试方法
CN104298844B (zh) * 2014-05-23 2017-04-12 中国科学院光电研究院 获取点阵法测量光学遥感载荷在轨mtf测量精度的方法
JP6370646B2 (ja) * 2014-09-05 2018-08-08 日本放送協会 Mtf測定装置
CN104318526B (zh) * 2014-10-20 2017-02-15 南京理工大学 基于mtf的遥感器在轨自动优化在轨参数的方法
CN104680138B (zh) * 2015-02-09 2017-11-03 北京空间飞行器总体设计部 基于特性参数相关性的卫星图像自动判读系统和方法
CN104820980B (zh) * 2015-04-15 2018-07-24 北京空间机电研究所 自适应高精度mtf测量方法
CN105976317B (zh) * 2016-04-28 2019-03-19 中国科学院遥感与数字地球研究所 一种图像空间退化模拟方法及系统
CN110020993B (zh) * 2018-11-28 2021-01-26 北京理工大学 基于一种靶标的空间重采样gf-4卫星的mtf计算方法

Also Published As

Publication number Publication date
CN110954303A (zh) 2020-04-03

Similar Documents

Publication Publication Date Title
CN108921885B (zh) 一种综合三类数据源联合反演森林地上生物量的方法
Baltsavias et al. Digital surface modelling by airborne laser scanning and digital photogrammetry for glacier monitoring
CN111781146B (zh) 利用高分辨率卫星光学影像的波浪参数反演方法
CN107145891B (zh) 一种基于遥感影像的水体提取方法及系统
CN110456352B (zh) 一种基于相干系数阈值的冰川识别方法
CN104217426A (zh) 一种基于ENVISAT ASAR与Landsat TM遥感数据面向对象提取水体的方法
CN112529788B (zh) 一种基于薄云厚度图估计的多光谱遥感图像薄云去除方法
CN110109118B (zh) 一种森林冠层生物量的预测方法
CN110174673B (zh) 一种利用时序接力干涉图叠加高效减弱大气相位影响的方法
CN109388887A (zh) 一种地面沉降影响因素定量分析方法及系统
CN110703244A (zh) 基于遥感数据识别城区内水体的方法和装置
CN114091274A (zh) 一种滑坡易发性评价方法及系统
CN114202535A (zh) 作物种植面积提取方法及装置
Zhang et al. Reconstruction of GF-1 soil moisture observation based on satellite and in situ sensor collaboration under full cloud contamination
CN111144350B (zh) 一种基于参考底图的遥感影像定位精度评价方法
CN110954303B (zh) 一种基于高分辨率遥感图像基准的mtf自动测算方法
Belfiore et al. Orthorectification and pan-sharpening of worldview-2 satellite imagery to produce high resolution coloured ortho-photos
CN112946647A (zh) 大气误差改正InSAR干涉图堆叠地质灾害普查方法和装置
CN117058522A (zh) 一种融合光谱神经网络和叶面积指数的遥感地表生态变化检测方法
CN117078033A (zh) 一种台风灾害影响评估方法、装置及计算机
CN116773516A (zh) 一种基于遥感数据的土壤碳含量分析系统
CN114239379A (zh) 一种基于形变检测的输电线路地质灾害分析方法及系统
CN114545410A (zh) 基于合成孔径雷达双极化数据相干性的作物倒伏监测方法
Cheng et al. Generation of pixel-level SAR image time series using a locally adaptive matching technique
Lombardi et al. Accuracy of high resolution CSK interferometric Digital Elevation Models

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
GR01 Patent grant
GR01 Patent grant