CN108174196B - 基于距离加权的成像系统调制传递函数测量方法 - Google Patents

Abstract

本发明公开了一种基于距离加权的成像系统调制传递函数测量方法，包括如下步骤：1)在拍摄得到的图像中选取合适的刃边区域，使用传感器的光电转换函数对图像数据进行线性化处理，得到待测刃边图像；2)在步骤1)得到的刃边图像中，通过逐行寻找每行边缘扩散函数的中值点得到刃边边缘位置；3)对步骤2)得到的刃边边缘位置进行最小二乘拟合，得到刃边位置函数；4)计算图像上每个像素到刃边的距离，根据距离加权得到等距采样的ESF；5)对步骤4)得到的ESF进行求导得到线扩散函数，对LSF进行傅里叶变换，得到成像系统的调制传递函数。本发明实现对含刃边区域的图像进行MTF测量，具有不受刃边角度限制、抗噪性好、计算准确、稳定性好等优点。

Description

基于距离加权的成像系统调制传递函数测量方法

技术领域

本发明属于成像系统质量评价领域，特别涉及一种针对光学数码成像系统的调制传递函数测量方法。

背景技术

在实际成像中，图像质量往往受到多方面的影响：成像系统的像差和衍射效应；传感器的分辨率、非线性响应和各类噪声；拍摄时外界的抖动和其他干扰等。对于成像系统成像质量评价，调制传递函数MTF(Modulation Transfer Function)是一个重要指标。它客观反映了不同空间频率信号经过成像系统后的衰减情况，代表成像过程中成像系统对输入信号的传递特性，是目前国际上通用的评定成像系统性能的指标之一。根据成像退化理论，如果成像系统的MTF可以精确测得，那么可以从退化图像中恢复得到真实图像。因此，对成像系统进行MTF测量具有非常重要的意义。

目前，针对数码成像系统的MTF测量，由于靶标选取条件相对宽松(人工靶标或合乎要求的刃边目标)，受噪声等因素干扰较小等因素，通常采用倾斜刃边法。ISO12233将倾斜刃边法作为电子静态图像相机分辨率测试的标准方法。

实际使用中，倾斜刃边法存在一定限制：数码成像是离散的等距采样，而刃边图像往刃边直线投影得到的边缘扩散函数的采样点间距为非等距，当使用ISO12233的合并取均值方法时，ESF采样结果会出现偏差，导致计算结果出现一定的偏差；当刃边图像存在噪声时，测量得到的ESF也必然被噪声所污染，求导得到LSF的过程会进一步放大噪声，导致测量结果失真；此外，计算误差随刃边角度变大，而实际应用中不容易控制刃边的角度。

现有技术为提高测量的准确度和稳定性，常常通过构造ESF的函数模型，对上采样的ESF数据进行非线性拟合，再用于下一步计算。这种方法可以提高计算的稳定性，但由于函数模型限制不能很好地逼近实际ESF，从而影响测量结果的准确度。另外，刃边角度对MTF测量结果的影响在现有方法中也无法得到很好的解决。

发明内容

本发明解决的技术问题是：针对倾斜刃边法测量数码成像系统MTF时受刃边角度和噪声影响，MTF测量结果不准确，提出一种基于距离加权的成像系统调制传递函数测量方法。

本发明的目的是通过以下技术方案来实现的：一种基于距离加权的成像系统调制传递函数测量方法，包括以下步骤：

1)在成像系统拍摄得到的图像中选取合适的刃边区域，使用成像系统的传感器的光电转换函数(OECF)对图像数据进行线性化处理，得到待测刃边图像；

2)在步骤1)得到的刃边图像中，通过逐行寻找每行边缘扩散函数(ESF)的中值点得到刃边边缘位置；

3)对步骤2)得到的刃边边缘位置进行最小二乘拟合，得到刃边位置函数；

4)计算图像上每个像素到刃边的距离，根据距离加权得到等距采样的ESF；

5)对步骤4)得到的ESF进行求导得到边缘扩散函数(LSF)，对LSF进行傅里叶变换，得到成像系统的调制传递函数。

进一步地，所述步骤2)中的逐行寻找每行ESF的中值点，通过以下步骤来确定：

2-1)逐行求导获得LSF，乘以中央对称的Hamming窗，积分得到加窗ESF；

2-2)计算步骤2-1)得到的各行加窗ESF的末像素的平均值，取其一半为中值，插值获取每行加窗ESF中值点所在位置；

2-3)以步骤2-2)的中值点位置为对称中心，逐行计算新的Hamming窗，并与每行的LSF相乘，积分得到加窗ESF；

2-4)计算步骤2-3)得到的各行加窗ESF的末像素的平均值，取其一半为中值，插值获取每行加窗ESF中值点的位置。

进一步地，所述步骤4)中的根据距离加权得到等距采样的ESF，通过以下步骤来确定：

4-1)确定过采样倍数f，根据与刃边距离大小，将图像像素分到间隔为dx的容器(bin)中；间隔dx计算公式为：

其中n为刃边图像列数；容器中央位置x_u为：

4-2)每个容器内所有像素按与容器中央的距离计算权重，计算公式为：

其中x_u为容器中央位置，i为容器内像素序号，x(x_u，i)为容器中像素与容器中央的距离，x(x_u，i)为每个像素的权重；

4-3)加权得到以容器中央位置为采样点的等距采样ESF，计算公式为：

esf(x_u)是以容器中央位置为采样点的等距采样的ESF。

进一步地，所述步骤2-4)中，中值v_mid的计算公式如下：

其中j为行序号，m为刃边图像行数，n为刃边图像列数，ESF(j，n)为第j行第i列的加窗ESF。

进一步地，所述步骤2-4)中，插值获取每行加窗ESF中值点的位置的计算公式如下：

其中l是ESF中小于v_mid且最接近v_mid的采样点位置，满足以下关系：

ESF(l)＜v_mid≤ESF(l+1)

本发明的有益效果是：本发明实现对含刃边区域的图像进行MTF测量，通过逐行寻找ESF中值点求解出精确的边缘位置，利用多像素加权获得等距采样的ESF，从而减轻刃边角度和图像噪声的影响，使得计算更为准确，计算结果稳定性更好。

附图说明

图1为本发明方法流程示意图；

图2为中央对称的Hamming窗函数的示意图；

图3为以中值点位置为对称中心的Hamming窗函数的示意图；

图4为图像所有像素构成的非等距采样ESF的示意图；

图5为容器内各个像素的权重分布示意图；

图6为加权得到等距ESF的示意图；

图7为无噪声的20°倾角刃边图像的MTF测量结果对比图；

图8为40dB高斯白噪声下5°倾角刃边图像的MTF测量相对误差分析图。

具体实施方式

以下结合附图和实施例详细说明本发明技术方案。参见图1，实施例的流程可以分为五个步骤：

步骤1：获取待测刃边图像，包括以下子步骤：

1.1通过待测成像系统进行成像，获取待测图像；

1.2在待测图像中手动或者通过算法自动选取大小为m×n像素的刃边区域；

1.3对选取的刃边图像，使用传感器的OECF进行线性化处理，得到大小为m×n像素的待测刃边图像；

步骤2：确定刃边边缘位置，包括以下子步骤：

2.1逐行求导获得序列长度为n的LSF，与中央对称的序列长度为n的Hamming窗函数点乘，得到序列长度为n的加窗LSF。其中Hamming窗函数计算公式为：

其中j为行序号，i为列序号。图2为中央对称的Hamming窗函数的示意图。

2.2对步骤2.1得到的加窗LSF进行积分，得到加窗ESF，计算公式为：

2.3计算步骤2.2得到的各行加窗ESF的末像素的平均值，取其一半得到中值，计算公式为：

2.4通过插值获取每行中值点所在位置：

其中l是ESF中小于v_mid且最接近v_mid的采样点位置，满足以下关系：

ESF(l)＜v_mid≤ESF(l+1)

2.5以步骤2.4得到的中值点位置为对称中心，逐行计算新的Hamming窗与每行的LSF相乘，得到序列长度为n的加窗LSF。其中Hamming窗函数计算公式为：

其中L_mid是对称中心的位置。图3为以中值点位置为对称中心的Hamming窗函数的示意图。

2.6对步骤2.5得到的加窗LSF进行积分，得到加窗ESF；

2.7计算步骤2.6得到的每行加窗ESF的末像素的平均值，取其一半得到中值；

2.8通过插值获取每行中值点所在位置；

步骤3：拟合得到刃边位置函数。通过步骤2计算得到每行的准确刃边边缘位置，对刃边位置进行最小二乘拟合，得到刃边位置函数：

A*x+B*y+C＝0

其中x和y分别为行序号和列序号。

步骤4：计算图像各像素到刃边的距离，加权得到等距采样的ESF。包括以下子步骤：

4.1计算图像各像素到刃边的距离，计算公式为：

其中x(j，i)和y(j，i)分别为图像像素的行序号和列序号。图4为图像所有像素构成的非等距采样ESF的示意图，横坐标为图像各像素到刃边的距离，纵坐标为图像像素值。

4.2确定过采样倍数f，根据与刃边距离大小，将图像像素分到间隔为dx的容器中。间隔dx计算公式为：

其中n为刃边图像列数。容器中央位置x_u为：

4.3每个容器内所有像素按与容器中央的距离计算权重，计算公式为：

其中x_u为容器中央位置，i为容器内像素序号，x(x_u，i)为容器中像素与容器中央的距离，计算公式为：

x(x_u，i)＝d(x_u，i)-x_u

图5为容器内各个像素的权重分布示意图。像素离容器中央位置越近，像素值与中央位置真实值越接近，因此权重越大。

4.4加权得到以容器中央位置为采样点的等距采样ESF，计算公式为：

图6为加权得到等距ESF的示意图。可以看到，受噪声影响，像素投影得到的非等距ESF出现锯齿状波动，进一步求取LSF时将放大噪声。通过使用加权方法，可以有效去除噪声影响，获得准确的ESF。

步骤5：对步骤4得到的ESF进行求导得到LSF，对LSF进行傅里叶变换，得到成像系统的调制传递函数。图7是对无噪声的20°倾角刃边图像进行MTF测量的结果，可以看到ISO12233方法出现较大偏差，而本方法与理论真值更为接近，测量结果更为准确。

图8为40dB高斯白噪声水平下本方法与ISO12233方法对5°倾角刃边图像MTF计算结果的相对误差分析。考虑噪声随机性，对100幅噪声图像进行计算，显示结果为相对误差平均值。其中相对误差RE的计算公式为：

其中MTF_cal(u)为各频率的MTF实测值，MTF_true(u)为各频率的MTF理论真值。可以看到在不同频率处，本发明方法都更有效地抑制噪声影响，计算更加稳定，计算结果更加精确。

Claims (4)

1.一种基于距离加权的成像系统调制传递函数测量方法，其特征在于，该方法包括以下步骤：

1)在成像系统拍摄得到的图像中选取合适的刃边区域，使用成像系统的传感器的光电转换函数(OECF)对图像数据进行线性化处理，得到待测刃边图像；

2)在步骤1)得到的刃边图像中，通过逐行寻找每行边缘扩散函数(ESF)的中值点得到刃边边缘位置；

3)对步骤2)得到的刃边边缘位置进行最小二乘拟合，得到刃边位置函数；

4)计算图像上每个像素到刃边的距离，根据距离加权得到等距采样的ESF，通过以下步骤来确定：

4-1)确定过采样倍数f，根据与刃边距离大小，将图像像素分到间隔为dx的容器(bin)中；间隔dx计算公式为：

其中n为刃边图像列数；容器中央位置x_u为：

4-2)每个容器内所有像素按与容器中央的距离计算权重，计算公式为：

其中x_u为容器中央位置，i为容器内像素序号，x(x_u,i)为容器中像素与容器中央的距离，w(x_u,i)为每个像素的权重；

4-3)加权得到以容器中央位置为采样点的等距采样ESF，计算公式为：

esf(x_u)是以容器中央位置为采样点的等距采样的ESF；

5)对步骤4)得到的ESF进行求导得到线扩散函数(LSF)，对LSF进行傅里叶变换，得到成像系统的调制传递函数。

2.如权利要求1所述的基于距离加权的成像系统调制传递函数测量方法，其特征在于：所述步骤2)中的逐行寻找每行ESF的中值点，通过以下步骤来确定：

2-1)逐行求导获得LSF，乘以中央对称的Hamming窗，积分得到加窗ESF；

2-2)计算步骤2-1)得到的各行加窗ESF的末像素的平均值，取其一半为中值，插值获取每行加窗ESF中值点所在位置；

2-3)以步骤2-2)的中值点位置为对称中心，逐行计算新的Hamming窗，并与每行的LSF相乘，积分得到加窗ESF；

2-4)计算步骤2-3)得到的各行加窗ESF的末像素的平均值，取其一半为中值，插值获取每行加窗ESF中值点的位置。

3.如权利要求2所述的基于距离加权的成像系统调制传递函数测量方法，其特征在于：所述步骤2-4)中，中值v_mid的计算公式如下：

其中j为行序号，m为刃边图像行数，n为刃边图像列数，ESF(j,n)为第j行第n列的加窗ESF。

4.如权利要求3所述的基于距离加权的成像系统调制传递函数测量方法，其特征在于：所述步骤2-4)中，插值获取每行加窗ESF中值点的位置的计算公式如下：

其中l是ESF中小于v_mid且最接近v_mid的采样点位置，满足以下关系：

ESF(l)＜v_mid≤ESF(l+1)。