Comparison of MTF measurements using edge method:

towards reference data set

Françoise Viallefont-Robinet, Dennis Helder, Renaud Fraisse, Amy Newbury,
Frans van den Bergh, Donghan Lee, Sébastien Saunier

To cite this version:

Françoise Viallefont-Robinet, Dennis Helder, Renaud Fraisse, Amy Newbury, Frans van den Bergh,
et al.. Comparison of MTF measurements using edge method: towards reference data set. Optics
Express, 2018, 26 (26), pp.33625-33648. �hal-02055611�

HAL Id: hal-02055611

https://hal.science/hal-02055611v1
Submitted on 4 Mar 2019

HAL is a multi-disciplinary open access L’archive ouverte pluridisciplinaire HAL, est

archive for the deposit and dissemination of sci- destinée au dépôt et à la diffusion de documents
entific research documents, whether they are pub- scientifiques de niveau recherche, publiés ou non,
lished or not. The documents may come from émanant des établissements d’enseignement et de
teaching and research institutions in France or recherche français ou étrangers, des laboratoires
abroad, or from public or private research centers. publics ou privés.
Comparison of MTF measurements using edge
method: towards reference data set
FRANÇOISE VIALLEFONT-ROBINET,1,* DENNIS HELDER,2 RENAUD FRAISSE,3
AMY NEWBURY,4 FRANS VAN DEN BERGH,5 DONGHAN LEE,6 and
SÉBASTIEN SAUNIER7
1
ONERA DOTA Université de Toulouse F-31055, Toulouse, France
2
SDSU, Brookings, SD 57007, United States
3
AIRBUS Defence and Space, 31402 Toulouse Cedex 4, France
4
DIGITAL GLOBE, 1300 West 120th Ave, Westminster, CO 80234, United States
5
Remote Sensing Research Unit, CSIR Meraka Institute, Pretoria, South Africa
6
Korea Aerospace Research Institute (KARI), 169-84 Gwahak-ro Yuseong-gu Daejeon 34133, South
Korea
7
TELESPAZIO VEGA France, 26, avenue Jean-François Champollion - BP 52309 31023 Toulouse
Cedex 1, France
*Francoise.Viallefont@onera.fr

Abstract: A sensor’s spatial resolution has traditionally been a difficult concept to define, but
all would agree that it is inextricably linked to the Ground Sampling Distance (GSD) and
Instantaneous Field of View (IFOV) of an imaging sensor system. As a measure of the
geospatial quality of imagery, the Modulation Transfer Function (MTF) of the system is often
used along with the signal-to-noise ratio (SNR). However, their calculation is not fully
standardized. Further, consistent measurements and comparisons are often hard to obtain.
Therefore, in the Infrared and Visible Optical Sensors (IVOS) subgroup of the Working
Group on Calibration Validation (WGCV) of the Committee for Earth Observation Satellites
(CEOS), a team from various countries and professional entities who are involved in MTF
measurement was established to address the issue of on-orbit MTF measurements and
comparisons. As a first step, a blind comparison of MTF measurements based on the slanted
edge approach has been undertaken. A set of both artificial and actual satellite edge images
was developed and a first comparison of processing results was generated. In all, seven
organizations contributed to the experiment and several significant results were generated in
2016. No single participant produced the best results for all test images as measured by either
the closest to the mean result, or closest to the truth for the synthetic test images. In addition,
close estimates of the MTF value at Nyquist did not ensure the accuracy of other MTF values
at other spatial frequencies. Some algorithm results showed that the accuracy of their
estimates depended upon the type of MTF curve that was being analyzed. After the initial
analysis, participants were allowed to modify their methodology and reprocess the test images
since, in several cases, the results contained errors. Results from the second iteration, in 2017,
verified that the anomalies in the experiment’s first iteration were due to errors in either
coding or methodology, or both. One organization implemented a third trial to fix software
errors. This emphasizes the importance of fully understanding both methodology and
implementation, in order to ensure accurate and repeatable results. To extend this comparison
study, a reference data set, which is composed of edge images and corresponding MTF
curves, will be built. A broader audience will be able to access the edge images through the
CEOS CalVal Portal (http://calvalportal.ceos.org/. This paper, which is associated with the
reference data set, can serve as a new tool to either implement or check, or both, the MTF
measurement that relies on the slanted edge method.

©2018 Optical Society of America under the terms of the OSA Open Access Publishing Agreement
1. Introduction
The geo-spatial quality of a sensor and its imagery often revolves, at least in part, around the
concept of the spatial resolution of a sensor which is often reduced to the Ground Sampling
Distance (GSD) associated with the Instantaneous Field Of View (IFOV) defined by the pixel
size. However, spatial resolution, and hence geo-spatial quality, is more complex than this.
Most agree that the effective spatial resolution is due to three (or four) features of the sensor:
the IFOV (and the GSD if different), the Modulation Transfer Function (MTF) and the signal
to noise ratio (SNR). The MTF is often used as a measure of image sharpness [1,2]. This
important parameter for image quality has to be checked on orbit in order to be sure that
launch vibrations, transition from air to vacuum, or thermal state have not degraded the
sharpness of the images. In some cases, it can lead to a refocusing decision.
This paper deals with one of the methods used for on-orbit MTF assessment, called the
edge method, the knife-edge method, or the slanted-edge method. This method is widely used
for laboratory measurements and may be implemented in various manners. For on-orbit MTF
assessment, it requires a slanted edge as explained in section 2. It has been used for numerous
space sensors such as Landsat TM [3], MOS-1 MESSR [4], IKONOS [5], SPOT5 [6], and
more recently Sentinel2 MSI [7].
In the framework of the Infrared and Visible Optical Sensors (IVOS) subgroup of the
Working Group on Calibration Validation (WGCV) of the Committee for Earth Observation
Satellites (CEOS), a team of people from various countries and professional entities, who are
involved in MTF measurement, has been created to address a variety of issues regarding the
geo-spatial quality of optical satellite imagery. One of the first efforts of this group has been
to compare processing methodologies for the edge method of MTF estimation. For this
comparison experiment, the team was composed of Frans van den Bergh from CSIR, Renaud
Fraisse from Airbus DS, Dennis Helder from SDSU, Dong Han Lee from KARI, Amy
Newbury and Robert Kudola from Digital Globe, Sébastien Saunier from Telespazio and
Françoise Viallefont-Robinet from ONERA
This paper presents the method and its various implementations followed by the
comparison experiment. The first results, obtained with a blind test approach, were analyzed.
This exercise was an opportunity to correct or improve the software of each participant. Thus,
a second run was performed by most of the participants in order to improve the results,
leading to a second comparison. For two test cases, a third and final run was performed by
one of the outliers. All comparisons are presented and commented on.
2. Edge method
2.1 Theory
Considering the sensor as a linear system without spatial variation of its response (shift
invariant), the relation between the radiances (or top of atmosphere reflectances) of the
landscape and the image is simply:
i ( x , y ) = l ( x , y ) ⊗ h( x , y ) (1)

where i(x,y) stands for the image,

l(x,y) stands for the landscape,
h(x,y) is the Point Spread Function of the sensor,
⊗ is the convolution integral.
For sensors using a CCD in the image plane, the system is no longer strictly shift
invariant; but it can nevertheless be described using the usual theory for regions of the
imaging array without loss of generality. The sampling done by the CCD can be written as a
multiplication by a Dirac comb.
Completing Eq. (1), it becomes:
 i( x, y ) = [l ( x, y ) ⊗ h( x, y)] ⋅ comb( x / px , y / p y ) (2)

where px is the size of the IFOV and py the GSD in the case of a pushbroom imaging system.
A classical way to deal with a convolution product is to apply a Fourier Transform, which
leads to:

I ( f x , f y ) =  L( f x , f y ) ⋅ H ( f x , f y )  ⊗ comb( f x / f sx , f y / f sy ) (3)

where I(fx,fy) stands for the Fourier Transform of the image,

L(fx,fy) stands for the Fourier Transform of the landscape,
H(fx,fy) is the Transfer Function of the sensor,
fsx = 1/px is the sampling frequency for the fx axis,
fsy = 1/py is the sampling frequency for the fy axis.
The sensor behavior is known to be similar to a low-pass filter without phase shift [8].
This is why the Optical Transfer Function is usually reduced to the Modulation Transfer
Function defined as the modulus of the Optical Transfer Function normalized by the zero
frequency component.
For the edge method, the landscape is close to a Heaviside function:
l ( x, y) = a ⋅ hea( x) ⋅ unit ( y ) + b ⋅ unit ( x) ⋅ unit ( y ) (4)

hea(x) being the Heaviside function centered on x = 0,

unit(x) = 1 for all values of x.
In this case, Eq. (2) becomes:
i( x, y ) = [ a ⋅ hea( x) ⋅ unit ( y ) ⊗ h( x, y) + b ⋅ unit ( x) ⋅ unit ( y)] ⋅ comb( x / px , y / p y ) (5)

As b ⋅ unit ( x) ⋅ unit ( y ) ⊗ h( x, y ) = b = b ⋅ unit ( x) ⋅ unit ( y )

Equation (5) can be rewritten as follows:
i( x, y ) = i( x) = [ a ⋅ hea( x) ⊗ LSF ( x) + b ⋅ unit ( x)] ⋅ comb( x / px ) (6)

The convolution by the comb produces aliasing. One way to overcome this problem is to use
an edge with a slight inclination relative to the row or column direction [9]. This is used to
build an oversampled 1-D edge image as illustrated in Fig. 1.
So, the 1-D edge corresponds to:
i ( z ) = a ⋅ hea( z ) ⊗ LSF ( z ) + b ⋅ unit ( z ) (7)

where z is measured in the direction perpendicular to the edge, which nearly coincides with
the x axis when the edge is vertically oriented as modelled here.
 Fig. 1. Principlee of 1-D oversamp
pled edge responsee (L represents the radiance).

2.2 Impleme
entation strate
egies
Implementatio
ons of the slaanted-edge meethod usually have three iddentifiable steeps: edge
modeling, Edge Spread Fun onstruction, andd the calculatioon of an MTF using the
nction (ESF) co
sampled ESF. The implementation of each of these steps can vary significantly, as discussed
in greater detail below. Some common traits of the participant implementations will be
discussed in this section.
Without loss of generality it will be assumed that the Region Of Interest (ROI) containing
the edge transition will be processed on a row-by-row basis with the slanted edge oriented
nearly vertically. Let R denote the ROI such that i(x, y) represents the intensity of the image
at coordinates (x, y) for all (x, y) ∊ R.
2.2.1 Edge modeling
The construction of an oversampled ESF, as illustrated in Fig. 1, requires that the location of
the edge in each image row is known with sub-pixel accuracy. This can be as simple as
calculating the centroid of the discrete derivative of the intensity of each row, as suggested in
the ISO 12233 standard [10, Appendix D], but such a method will be sensitive to image noise
and target non-uniformity. A more robust method is to fit a parametric- or spline function to
each row, using the inflection point of the fitted function as the sub-pixel edge location
estimate.
Once the location of the edge has been estimated in each row, a linear function is typically
fitted across all the per-row results to obtain a more accurate model describing the sub-pixel
location of the edge, as illustrated in section 3.1 in Fig. 3. If the physical target edge is curved,
the linear edge model can be replaced with a low-order polynomial to accommodate the
curvature [11].
2.2.2 ESF construction
A simplified example of the construction of an ESF is illustrated in Fig. 1, where the
oversampled ESF is obtained by interleaving the intensity values of each row of the ROI.
Constructing the ESF involves projecting the 2-D image intensity values i(x,y) onto a 1-D
representation i(z). The magnitude of z represents the shortest distance from a pixel at
coordinates (x,y) to the slanted edge. The oversampled ESF can be constructed using an
extension of the method described in the ISO 12233 standard [10, Appendix D], or using the
alternative method described by Kohm [12].
In the ISO 12233-based approach the ROI is a rectangle aligned with the rows and
columns of the image, as shown in Fig. 1. For each pixel the value z is calculated as:
z = [ x − e( y )] ⋅ cos(θ ) (8)

where e(y) denotes the location of the edge in row y, as predicted by the edge model, and θ
denotes the relative edge angle. The cosine factor transforms a distance measured along a row
into a distance measured perpendicularly to the edge.
With Kohm’s method, the ROI is a rectangle that is aligned with the edge itself, as shown

in Fig. 2. A unit vector perpendicular to the edge, n , is used to calculate z such that:
  
z = n ⋅ ( p x , y − p x0 , y0 ) (9)

where (x0,y0) are the coordinates of an arbitrary point on the edge.

 Fig. 2.
2 An alternative method
m of construccting the oversamppled ESF. The edgge (1) is identifiedd
during
g the edge modeliing step. A vector perpendicular too the edge (2) is constructed. Eachh
point (3) in the ROI is projected
p onto the perpendicular
p vecttor (4).

Both thesee approaches yield

y an oversaampled ESF inn the form of a set of tupless {[z, i(x,
x, y) ∊ R. In geeneral the z vaalues are not unniformly spaceed, which makes the set
y)]} for all (x
{[z, i(x, y)]} unsuitable fo or the next steep in the slannted-edge algorrithm involvinng a Fast
Fourier Transsform (FFT).
The ESF canc be resamplled to a uniform m spacing by ffitting a modell to the set {[z, i(x, y)]}
followed by sampling
s the model
m at the dessired spacing ((1/4th pixel, foor example). Thhis model
can be piecew wise, like the LOESS
L method d employed by Kohm [12], oor global like a spline or
another suitabble parametric function. A diifferent approaach, inspired byy the method ssuggested
in the ISO 12 2233 standard, is to bin the seet {[z, i(x, y)]}} into uniformlly spaced bins,, but with
additional low w-pass filterin
ng to further reduce
r the imppact of image noise and tarrget non-
uniformity.
2.2.3 MTF ca
alculation and
d correction
Once the ESF F has been resaampled to a uniform spacing,, the final step is to compute the MTF
either using th
he derivative method,
m or the spectral
s methood [9]. The deriivative methodd involves
computing th he finite-difference derivativ ve approximatiion of the uniiformly-spacedd ESF to
produce the Line
L Spread Fu unction (LSF). The MTF is thhen calculated as the magnituude of the
FFT of the LS SF, normalized d by the zero-ffrequency compponent. Care m must be taken tto correct
the MTF derrived this way y to compensaate for the finiite-difference approximationn used to
obtain the LSF F [13].
The specttral approach bypasses LSF F estimation ccompletely, instead of remooving the
contribution of
o the ideal steep edge (the ex xpected targett scene) from tthe observed E ESF, both
are transformed to the frequ uency domain, to obtain the ssystem MTF. T The spectral appproach is
described in detail
d in section
n 3.1.
3. Implementations of the edge method
3.1 ONERA
The ONERA implementation follows the spectral approach [9], which is less widely known
than the derivative method [4,5,14].
In this case, the finite number of samples has to be taken into account. This can be done as
follows:
i( z ) = [ a ⋅ hea( z ) ⊗ LSF ( z ) + b ⋅ unit ( z )] ⋅ w( z ) (10)

noting w(z) as the window corresponding to the finite interval.

In the Fourier domain, the relation becomes:
I ( f z ) = [ a ⋅ Hea( f z ) ⋅ H ( f z ,0) + b ⋅ δ ( f z )] ⊗ W ( f z ) (11)

The window has to be chosen so that Hea(fz)⊗W (fz) ≠ 0 for all frequencies and not far from
Hea(fz).
After removal of the background b, the following ratio gives the transfer function:
I ( f z ) / [ a ⋅ Hea( f z ) ⊗ W ( f z )] ≈ H ( f z , 0) (12)
ONERA tool follows the 3 steps general implementation.
For edge modeling, or in other words the 2-D to 1-D transformation, each row is
interpolated using a spline function and the inflection point is computed. A straight line is
fitted on the set of inflection points and may be used (depending on the user choice) to replace
the positions found as shown in Fig. 3. The inverse of the slope of the straight line provides
the oversampling rate for building the ESF. At this stage, the 2-D edge image is split into a
list of rows where each ESF location is indicated by the position of the edge e(y).
The second step aims at mixing the rows to construct the oversampled ESF. For the
inclination, in the ideal case of Fig. 1, the sequence spans a whole number (4 in this case) of
rows. However, it is not possible to have a whole number of rows corresponding to the
sequence for any edge. Thus, the user chooses the oversampling rate Nr, not too far from Ns
(Ns being the whole number closest to the actual number of rows of the sequence) and then
each row of the sequence is properly phased in the Nr grid. Some cases of Nr ≠ Ns, may lead
to incomplete sampling.
A user can choose between linear interpolation or model fitting options to construct an
ESF with regular sampling. In the model fitting case, all the edge samples are replaced by the
samples deduced from a parametric transfer function model as explained in [15]. The model
fitting implementation eliminates noise due to non-uniformity of the dark or light areas of the
ROI, which is an added advantage even when the sampling is regular.
Once the regularly sampled ESF is built, the parameters a and b, as well as the position z
= z0 of the Heaviside function, are assessed. Several approaches are available to the user for
selection of a, b, and z0:
• the first value of the edge for a and the last value for b,
• the values at a distance of -d (relative to z0) for a and + d for b,
• the value corresponding to (a + b)/2 for z0,
• the value corresponding to the inflection point for z0
Once the parameters a, b and z0 are found, the distribution a.Hea(z) + b is built and drawn
over the ESF. The background b is then removed.
ESF may be artificially extended, as shown Fig. 4 in order to reach the desired sampling
rate for the MTF curve and to limit windowing effects.
 A Hann window is applied to the ESF and to the distribution. An FFT is applied to the
windowed ESF and to the windowed distribution. The ratio of the modulus of each spectrum
is computed and normalized to obtain the MTF according to Eq. (12). The fz frequencies are
converted to fx frequencies. Once the MTF is computed, there is a possibility to fit the model
described in [15] to the values obtained and to resample the MTF in order to obtain the
desired frequency step.

Fig. 3. Location of the edge and angle assessment.

Fig. 4. Measured edge after truncation and extension and corresponding MTF curve.

3.2 SDSU
Based on the theory presented in Section 2, the South Dakota State University
implementation follows a series of steps designed to be largely target and model independent
so the algorithm will work with a variety of edge target types. The steps are described as
follows.
First, it is important that, to the degree possible, an edge target is used that has a proper
orientation with respect to the sampling grid of the satellite sensor. Some combinations of
edge angles and sampling grid orientation will result in data that are not reasonably uniformly
distributed as shown in Fig. 1. When this happens, the result is that the steep part of the edge
is represented by data points that are clustered together. This leads to very poor LSF/MTF
estimates. Extensive modeling has indicated that relative angles of 6-8 degrees are optimal.
This has the added advantage of orienting the edge nearly orthogonal to the sampling grid,
and minimizes the correction necessary to obtain LSF/MTF estimates in the typical along-
track and cross-track directions in which sensor specifications are often given.
A second critical aspect of target development is that the length of the edge should span a
sufficient number of rows (or columns) in the image so that the oversampling process
produces enough samples for edge reconstruction to be accurate. Empirical analysis has
indicated that a minimum of 20 cross sections of the edge should be obtained for accurate
results.
Signal-to-noise ratio (SNR) has a significant bearing on PSF/MTF estimation accuracy.
For purposes of LSF/MTF estimation, SNR can be defined as the ratio of the edge height to
the average of the standard deviations of the region on either side of the edge. Modeling has
indicated that SNR > 50 produces accurate and consistent results.
The first step to developing an oversampled ESF is to estimate the edge location from
each slice of data across the edge. A simple, but accurate, approach is to fit a Fermi function
to the data of the form of Eq. (13):
b−d
f ( x) = d + (13)
1 + e − s ( x − e ( y ))
where x represents the pixel locations for row y, d is the bias level, b is the magnitude of the
bright side of the edge, e(y) is the edge location and s represents the steepness of the edge.
This approach does not model ESF which have ringing in them, but will still fit the steep part
of the ESF well and give good estimates of the precise edge location which is the goal of this
step. To determine the parameter values, a common optimization algorithm, such as the
Levenberg-Marquardt algorithm, is employed. The output of this step is the parameter, e(y),
which provides, for each row, subpixel estimate of the true edge location using the integer-
based grid of the edge image input data.
An oversampled but irregularly spaced ESF is constructed using Eq. (8). Only the data
within a distance of five pixels from the edge are retained, based on typical LSF width of two
pixels.
The truncated, oversampled ESF is simultaneously filtered (to reduce high frequency
noise) and resampled (to obtain uniformly spaced samples) using a non-linear modified
Savitzky-Golay filter. In this approach, a window of data is selected (typically set to two
pixels in length) and a low-order polynomial is fitted to the data using a linear least squares
approach. Through extensive modeling, it was determined that a fourth order polynomial is
optimal. The output of the filter is the value of the polynomial at the center of the window.
The window is then shifted by the amount of the desired output sampling interval, and the
process is repeated. It is recommended to oversample by a factor of 10 or more. This step
produces as an output a uniformly oversampled ESF that has been smoothed for high
frequency noise but, at the same time, has not been modified significantly at frequencies
below the Nyquist frequency.
Because smoothing has already been done to the ESF, a simple first order differencing
approach is employed to obtain the LSF. The final step is to obtain the MTF by applying the
Fourier transform (via the Fast Fourier Transform) to the LSF and normalizing by the zeroth
value.
The various steps in the SDSU MTF estimation process are illustrated Fig. 5. In the upper
left corner is an actual satellite image of a slant edge that was obtained from deployed tarps.
The upper right chart shows the oversampled ESF that has been produced after application of
the modified Savitzky-Golay filter. Note the uniform spacing of the data, even in the steep
region of the edge response. At this point an estimation of the SNR for this edge response can
be calculated. The lower left plot shows the line spread function after simple first order
differentiation
n. Here the datta values have been suppresssed and a simpple curve conneecting the
data points is shown for claarity. Note the noise present on either side of the LSF, deespite the
observation thhat SNR is greeater than 90 for
fo this particullar example. TThe correspondding MTF
function is shhown in the lowwer right whicch indicates forr this examplee an MTF at N
Nyquist of
0.1117.

Fig. 5.
5 Illustration of SD
DSU LSF/MTF esstimation process.

3.3 Airbus DS
D
Airbus DS has developed and used its edge method implementatioon to evaluatee all high
resolution satellites launched over the last 20 years. It iss used operationnally for all thhe in-orbit
follow-on monitoring
m of Airbus DS satellites, inccluding SPOT T and Sentinel2. The
implementatio on follows thee familiar stepss: edge modellling, ESF moddelling, and fiinal MTF
calculation ussing the LSF.
The edge model suppo orts both lineaar target edges and curved edges, allowiing more
accurate modeling of naturaal edges that arre not necessarrily straight. Thhe sub-pixel loocation of
the edge in eaach image row is estimated by y low order po lynomial fittinng.
After the edge
e model paarameters have been obtainedd, the oversamppled ESF is coonstructed
using Eq. (8)), after aligninng the individu ual per-row eddge profiles uusing the edgee location
predicted by the (potentiallly curved) edg ge model. A uuniformly spacced ESF is obttained by
fitting a suitab
ble model to th
he irregularly spaced
s ESF, foollowed by sam mpling the ESF model at
the desired saampling interv vals. Differentt models can bbe used, somee being theoreetical like
sigmoid, oth her being link ked to the physics
p (introoducing know wledge on thee optical
combination).. Optional outlier rejection caan be activatedd during the ESSF model fittinng stage if
needed.
The LSF isi obtained by numerical diffferentiation of the uniformly spaced ESF. T The MTF
a the normalizzed magnitude of the FFT of tthe LSF.
is computed as
The robusstness of the im
mplementation is improved byy careful filteriing and/or moddel fitting
at the variouss stages of the algorithm to reduce
r the imppact of measurrement noise aand target
non-uniformitty. The Airbus DS implemen ntation steps aree illustrated Figg. 6.

Fig. 6. Illu
ustration of Airbus DS operational M
MTF measurement tool.

3.4 Digital Globe

G
Digital Globee uses the Imagge Quality Meetric Toolkit too measure MTF F using a grouund-based
tilted MTF target made up ofo black and white
w square arreas that are appproximately 220 m on a
side. The toolkit was develloped by Harriis Corporation and is a plugg-in to the Exeelis ENVI
image processsing softwaree. The toolkiit takes the measured edgge response aand uses
enhancementss to the method dology describ
bed by Taitian in his 1965 paaper [16] to coonvert the
measured edg ge response to MTF
M as a funcction of samplinng frequency. On-orbit meassurements
of Nyquist MTF
M using this tool have been consistent too within aboutt 0.006 (1σ) foor periods
well over 1 yeear.
3.5 CSIR
The CSIR’s im mplementation n of the slanted
d-edge method,, as used in thiss study, was innspired by
the method of o Kohm [12]]. The only additional
a userr input to thee algorithm iss a mask
specifying thee region of intterest (ROI) foor each edge too be analyzed. The ROI is tyypically a
rectangular reegion that is aliigned with the edge as in Fig.. 2.
For the edge
e model estimation step p, the edge loocation and oorientation is eestimated
iteratively ussing weighted d Principal Co omponent Annalysis (PCA) applied to tthe (x,y)
coordinates of
o the samples within the RO OI, where the weight of eacch sample is thhe image
gradient magn nitude raised to
t the fourth power. Once ann initial estimaate of the edgee centroid
and edge normal vector is obtained, the image gradient magnitude values are projected onto
the normal and grouped into coarse bins by their signed distance from the edge. Outliers
within each bin are flagged using a variant of Tukey’s quartile test. The PCA and outlier
rejection steps are iterated until the outliers stabilize. In this context, the PCA provides a total-
least-squares estimate of the edge parameters.
For the ESF/LSF construction step, following the method of Kohm outlined in section 2.2,
the coordinates of the samples within the ROI are projected onto the edge normal to yield the
set{[z, i(x, y)]} using Eq. (9). Outliers identified in the edge parameter estimation stage are
excluded. An ESF with a regular sample spacing of 1/8th pixel is constructed by weighted
binning of the set {[z, i(x, y)]}. To avoid poor sampling caused by certain edge slopes (e.g.,
1/2 or 1/4), the value of each ESF bin k is calculated according to Eq. (14):

 f [ z ( x, y ) − m ] i ( x, y )
x, y
k

e= (14)
 f [ z ( x, y ) − m ]
x, y
k

where mk denotes the midpoint value (distance from edge) of bin k, and f(d) a low-pass kernel
function. In practice, the function f(z) = exp(−13|z|) works well. The effect of this low-pass
filter must be removed from the final MTF by dividing the measured MTF by the Fourier
transform of f(z), i.e., 132/(132 + fz2) in its normalized form. The proposed kernel f(z) is wide
enough so that even if the edge angle is 26.565 degrees (a slope of 1/2), none of the ESF bins
will have a zero denominator.
The notion of filtering the ESF during construction is taken one step further by switching
to a different low-pass kernel for the tails of the ESF, similar to the method proposed by
Williams and Burns [14]. In particular, using f(z) = rect(k·z) with k decreasing with distance
from the edge reduces the impact of noise on the eventual MTF measurement. The starting
locations of the ESF tails are defined relative to the 10% and 90% quantiles of a heavily
smoothed temporary ESF. No correction of the MTF is applied for these ESF-tail low-pass
kernels.
The LSF is constructed by a finite-difference approximation of the ESF derivative. The
ESF-tail smoothing described previously obviates the need for windowing of the LSF before
applying the FFT, since the LSF tails naturally taper to zero with sufficient smoothing.
For the final MTF calculation step, the FFT of the LSF is computed, followed by
normalization. The appropriate sinc(c·fz) correction is applied to the final MTF to compensate
for the finite-difference approximation. To reduce the variance of the MTF estimates further,
a variable-width Savitzky-Golay filter is applied. The width of the filter increases gradually at
higher normalized frequency values.
3.6 KARI
The KARI slanted-edge implementation has been developed to measure the spatial quality of
the KOMPSAT image data starting from the SDSU algorithm [5]. The same implementation
was also used for ground testing before launch.
The selection of the inputs for the current KARI algorithm are the ones leading to is the
largest values of Relative Edge Response (RER), Function Width at half maximum (FWHM)
for LSF and MTF.
The edge modeling involves fitting an edge model, starting with an initial estimate of the
edge location within each row of the ROI that is obtained by finding the pair of adjacent
pixels with the largest difference. This estimate is refined by computing the inflection point of
a cubic polynomial fitted through the four values surrounding the initial estimate. The overall
edge model is obtained by fitting a linear function across all the rows through the refined edge
locations.
 In the ESSF construction n step, the irrregularly spacced oversampleed ESF is obtained as
described in Section
S 2.2, usiing Eq. (8). It involves trimm
ming the ESF tto contain onlyy the edge
transition regiion in order too limit the imp pact of noise annd target non-uuniformity as shown in
Fig. 7(3). Thee edge transitiion region of the t ESF is inittially estimated in three stepps: first a
cubic smooth hing spline is fitted to the oversampled
o E
ESF, then the L LSF is obtaineed as the
derivative of the smoothed d spline, and laastly the kneee points of thee LSF are founnd as the
extrema of the second derivative
d of the LSF. Thhese LSF knnee points coorrespond
approximately y to the pointss of maximum curvature in tthe LSF, and ttherefore correespond to
the knee pointts of the ESF, illustrated
i withh the short dashhed blue lines in Fig. 7(4:a). The knee
points are exttended outward ds from the ed dge by one pixxel to define thhe final edge transition
region, illustrrated as the innterval between n the dotted b lue lines of Fig. 7(4:b,c). F From here
onwards, only y this central part of the ES SF is considerred. Another ccubic smoothinng spline
(MATLAB’s csaps routine with p = 0.98)) is fitted to thee trimmed ESF F, which is sammpled at a
uniform spaciing determined d by the desired oversamplinng factor. The ffirst-order diffference of
the uniformlyy sampled trimm med ESF yield ds the LSF.
In the lastt step the FFT of the LSF iss computed, annd the resultingg magnitude vvalues are
normalized too obtain the MT TF.
The follow wing metrics area calculated in addition too the MTF: R Relative Edge R Response
(RER), LSF width
w at 25%, 50%
5 (FWHM),, and 80% of m maximum as illlustrated in Figg. 8.

Fig. 7. Edge Spread Function

F from Edgee data by KARI.

KARI’s im
mplementation adds five consstraints designeed to achieve m
more accurate rresults:
(1) Edgee step magnitu
ude: greater th
han 2000 DN
N (KOMPSAT
T-3 is a 14 biit sensor,
0~166383 DN).
(2) Edge angle: greater than or equal to
t 5°.
(3) Numb
ber of rows in ROI:
R at least 15 lines.
(4) RMSE ge model fit: at most 3 pixels..
E of linear edg
(5) MTF must be strictlly non-increasiing up to the N
Nyquist frequeency. A violatioon of this
rule is indicated by
y the red arrow in Fig. 8.
 Fig. 8. Constrain
nt in KARI implem
mentation.

3.7 Telespaz
zio
The edge MT TF measuremen nt method pressented here haas been developped and operatted in the
context of ESSA contribution n to the ALOS S PRISM calibbration campaiign [17]. The iimproved
method and th he basis of the algorithm usedd in this compaarison exercisee are described below.
Step 1 - Edge identification and orrientation: the main objectiive of this staage is to
determine thee location of thhe edge and itss angle with thhe Along-track (AL) and Acrross-track
(AC)) directioons, as shown in Fig. 9. A seecond objectivve is to check eeach image linne/column
in terms of no
oise and contraast affecting thee edge profile. Beside the usee of an edge m
model, this
preparatory taask requires maanual operationns and, in partiicular, visual innspection.
 Fig. 9.
9 As the initial sttep, the operator identifies suitablee key marks whicch will be used too
estimaate parameters of the two edge linees (depicted in yelllow). Key marks are select at eachh
side of
o the image. The location of a given n edge line is estim
mated at two distincct points (depictedd
in bluue) and then linee parameters com mputed. The mainn objective of thi s procedure is too
compu ute the angle betwween Along-track (AL),
( Across-trackk (AC) directions, respectively, andd
the veertical/horizontal directions.
d

Step 2 - Edge
E spread fuunction constru uction: oversam
mpling of the iinput target is therefore
accomplished d using the meethod proposed d by Kohm [112], as describbed in section 2.2. The
algorithm pro ojects edge pro ofiles (in AL, AC directionss) onto the peerpendicular line to the
edge, using Eq.
E (9) as illusttrated in Fig. 3.3 The perpenddicular directioon is deduced from the
angle discusseed above. All edge
e profiles are
a projected too construct a noon-equally spaaced ESF.
The purpose of o the next step
p is to provide an ESF with unniform samplinng.
Step 3 - ESFE Modeling g: different meethods exist fofor the ESF m modeling, and these are
divided into two main caategories: mod deling based on a non-parrametric approoach and
modeling baseed on a parameetric approach.. These approaaches were testeed for our test case and,
finally, the paarametric curvee method proposed in [18] w was selected forr its robustnesss to noise
and flexibilitty as discusseed above. Parrameters of thhe following ffunction are eestimated
iteratively.
 x − a3 
f ( x) = a1 ⋅ erf   + a2 (15)
σ 2 
where a3 is thhe parameter designating
d the inflection poinnt.
Some sho ortcomings hav ve been observ ved with this model: asymm metric ESFs ccannot be
modelled, and d the model fits poorly in thee corners (kneee points) of the ESF. For thesse reasons
a second curvve fitting, appliied to locally deform
d the curvve in the transiition area, is peerformed.
Finally, the ESF parametric model is summ marized with thhe following foormula:
3
e ( z ) = f ( z ) + w( z − a3 ) c j ( z − a3 ) 2 j −1
esf (16)
j =1

where w( z ) ) is the Hann window

w and thee c j coefficiennts are estimateed by using leaast square
adjustment method. Results are shown in Fig.
F 10.
 Fig. 10.
1 ESF parametric function (depictted in green) withh parameters estim
mated by using thee
non-equally spaced ESF
F data points (depiicted in blue).

This stagee requires seveeral manual operations to set up the transitiion area extentt in order
to control thee deformation applied. Quallity informatioon such as Reelative Edge R Response,
Signal to No oise Ratio and d the L2 Norm m residual aree used to asseess the validitty of the
procedure. Viisual inspectionn remains an im
mportant aspecct.
Step 4 - Modulation
M Traansfer Functionn Calculation: Numerical diffferentiation off the ESF
discussed abo ove leads to th he LSF as sho own in Fig. 111(a). Special attention is paaid when
trimming the LSF to obtain n sufficient poiints in the MT
TF. An example of the MTF,, with 0.1
frequency stepp, is shown in Fig. 11(b) belo
ow.

Fig. 11. (a) LSF producced with the first derivative

d of ESF ccurve, (b) MTF noormalized moduluss
produuced from the LSF..

4. Comparis
son experime
ent
Airbus DS, Digital
D Globe annd CSIR have provided
p sampple image data of slanted-edgge targets.
The availablee data set is large: 20 edg ges. A subset of six edges was selectedd for this
experiment ass shown Fig. 12 m are simulatedd and the remaaining two are extracted
2. Four of them
from actual im
mages.
 Fig. 12. Ov
verview of the datta set.

Concernin ng simulations,, two images have

h been geneerated thanks tto analytic Poinnt Spread
Functions, wiithout noise, acccording to [19 9]. The angle iis small: less tthan 5°. The tw
wo others
are as represeentative as po ossible of real acquisitions iin term of sennsor performannces. The
pattern is a sllanted edge of angle close to o 15° and size 50 meters. Thhe simulated M MTF takes
into account typical
t optical diffraction and aberrations, the detector M MTF and some dynamic
contribution. Then representtative noise is added (photonn noise plus darrk noise) with a SNR of
150 for the bright
b area. Twwo values of GSD
G have beeen provided: 1m with high M MTF and
30cm with low w MTF.
The data were put on the t Cal/Val Po ortal (calvalpoortal.ceos.org) with, for the moment,
access limited d to the MTF project team. Each particippant was askedd to obtain thhe sample
images, to pro ocess them and d to send results to ONERA A. ONERA perrformed the coomparison
using a blind experimental design: each participant
p wass assigned a leetter from the alphabet,
which was sh hared with on nly that particcipant, allowinng participantss to identify ttheir own
results, but no
ot those of the others.
o
The first results,
r obtaineed with this bllind test approoach were colleected in 2016. To draw
further benefiit from the experiment, a seco ond run was prroposed in ordeer to correct orr improve
the processess and/or the inputs
i used. This
T second rrun was perfoormed by mosst of the
participants leeading to a seccond comparison in 2017. Foor the StdSysteem cases, an aadditional
third run was performed by participant A.
Table 1. suummarizes thee history of resu ult submission s.
Table 1. Result submission h
history
Data namee Meassured by
StdSystem_
_1 All (77), reprocessed by A, B, C, D and F
StdSystem_
_30 All (77), reprocessed by A, B, C, D and F
Apnn 5, repprocessed by B, D and F
Cgpnn 5, repprocessed by B, D and F
14oct…_P3 6, repprocessed by B, C, D and F
15aug…_P
P3 5, repprocessed by B, C, D and F

5. Results
For each edgee the following
g will be presen
nted:
• the 20116 and 2017 MTF M curves coorresponding reespectively to the first resullts and to
me participants ((according to T
resullts after reproceessing for som Table 1),
 • the curves of the difference relative to the mean MTF or to the known MTF,
• a table with the mean, the standard deviation and the largest difference (max-min) at
Nyquist frequency.
In all graphs and in all tables of this section, ACT stands for across track direction or axis,
and ALT stands for along track direction or axis. In 2016, the mean is computed using all the
results. In 2017, the mean is computed without the outliers.
The MTF curves will be limited to the across track direction as they provide enough
illustration of the discrepancies that were observed.
5.1 Results for StdSystem_1 edges
The 2016 and first 2017 (2017a) MTF curves and the discrepancy for the StdSystem_1 edge
are presented in Fig. 13. Another graph, Fig. 14, provides the same curves but related to last
2017 (2017b) results. Table 2 gives the 2016 and the last 2017 (2017b) MTF values obtained
at Nyquist frequency.
(a) StdSystem_1 (2016) (b) StdSystem_1 (2017a)
G 1 G
1
F F
0,8 0,8
E E
0,6 0,6
D D
ACT MTF

ACT MTF

0,4 C 0,4 C

B B
0,2 0,2
A A
0 0
Model 0 0,2 0,4 0,6 0,8 1 Model
0 0,2 0,4 0,6 0,8 1
Normalized frequency Normalized frequency

(c) Deviation to the model (d) Deviation to the model

StdSystem_1 (2016) StdSystem_1 (2017a)
G G
0,1 F 0,1 F
ACT MTF difference

ACT MTF difference

0,05 E 0,05 E

D D
0 0
0 0,2 0,4 0,6 0,8 1 C 0 0,2 0,4 0,6 0,8 1 C
-0,05 B -0,05
B

-0,1 A -0,1 A
Normalized frequency Normalized frequency

Fig. 13. Results for the StdSystem_1 edge: (a) MTF curves from 2016 runs, (b) MTF curves
from first 2017 runs, (c) deviation to the model from 2016 runs, (d) deviation to the model
from first 2017 runs.

(a) StdSystem_1 (2017b) (b) Deviation to the model

G StdSystem_1 (2017b)
G
1 F
0,1 F
0,8 E
ACT MTF difference

0,05 E
0,6 D
ACT MTF

C D
0,4 0
B 0 0,2 0,4 0,6 0,8 1 C
0,2 -0,05
A B
0 Model
-0,1 A
0 0,2 0,4 0,6 0,8 1
Normalized frequency Normalized frequency

Fig. 14. Last 2017 results for the StdSystem_1 edge: (a) MTF, (b) deviation to the model.
 Table 2. Results at Nyquist frequency for StdSystem_1 edges.
ACT 2016 ALT 2016 ACT 2017b ALT 2017b

Mean 0.29 0.29 0.29 0.30

Standard deviation 0.03 0.02 0.02 0.02
Max-min 0.08 0.07 0.04 0.04
For this first case, the 2017b results are in quite good agreement. The discrepancy
between the participants increases with the frequency but remains small up to the Nyquist
frequency. This illustrates the value of the comparison for purposes of improving algorithms
and removing errors in algorithms.
5.2 Results for StdSystem_30 edges
The 2016 and first 2017 (2017a) MTF curves and the discrepancy for the StdSystem_30 edge
are presented in Fig. 15. Another graph Fig. 16 provides the same curves but corresponding to
last 2017 (2017b) results. Table 3 shows the 2016 and last 2017 (2017b) MTF values at
Nyquist frequency.
(a) StdSystem_30 (2016) (b) StdSystem_30 (2017a)
G G
1 1
F F
0,8 0,8
E E
0,6 0,6
ACT MTF

ACT MTF

D D

0,4 C 0,4 C

B B
0,2 0,2
A A
0 0
0 0,2 0,4 0,6 0,8 1 Model 0 0,2 0,4 0,6 0,8 1 Model
Normalized frequency Normalized frequency

(c) Deviation to the model (d) Deviation to the model

StdSystem_30 (2016) StdSystem_30 (2017a)
G G

0,1 F 0,1 F
ACT MTF difference

ACT MTF difference

0,05 E 0,05 E

D D
0 0
0 0,2 0,4 0,6 0,8 1 C 0 0,2 0,4 0,6 0,8 1 C
-0,05 -0,05
B B
-0,1 A -0,1 A
Normalized frequency Normalized frequency

Fig. 15. Results for the StdSystem_30 edge: (a) MTF curves from 2016 runs, (b) MTF curves
from first 2017 (2017a) runs, (c) deviation to the model from 2016 runs, (d) deviation to the
model from first 2017 (2017a) runs.

(a) StdSystem_30 (2017b) (b) Deviation to the model

G
StdSystem_30 (2017b)
1
F G
0,1
0,8
ACT MTF difference

E F
0,05
0,6 D E
ACT MTF

C 0 D
0,4
0 0,2 0,4 0,6 0,8 1
B C
0,2 -0,05
A B
0
Model -0,1 A
0 0,2 0,4 0,6 0,8 1
Normalized frequency Normalized frequency

Fig. 16. Last 2017 (2017b) results for the StdSystem_30 edge: (a) MTF, (b) deviation to the
model.
 Table 3. Results at Nyquist frequency for StdSystem_30 edges (without D for 2017)
ACT 2016 ALT 2016 ACT 2017b ALT 2017b

Mean 0.08 0.09 0.10 0.09

Standard deviation 0.04 0.04 0.01 0.01
Max-min 0.11 0.13 0.04 0.02
Initially, for this second case, there were two singular results. The successive reprocessing
for case A clearly improved the result. Except for very low frequencies, the discrepancy
between the participants is small.
5.3 Results for apnn edge
The 2016 and 2017 MTF curves and the discrepancy for the apnn edge are shown in Fig. 17.
Table 4 gives the 2016 and 2017 MTF values at Nyquist frequency.
(a) apnn (2016) (b) apnn (2017)

1 1
G G
0,8 0,8
D F
0,6 0,6
ACT MTF

ACT MTF
D
C
0,4 0,4 C
B
0,2 0,2 B
Model
0 Model
0
0 0,2 0,4 0,6 0,8 1 0 0,2 0,4 0,6 0,8 1
Normalized frequency Normalized frequency

(c) Deviation to the model (d) Deviation to the model

apnn (2016) apnn (2017)

0,1 0,1
ACT MTF difference

G
ACT MTF difference

G
0,05 0,05 F
D
0 0 D
0 0,2 0,4 0,6 0,8 1 C 0 0,2 0,4 0,6 0,8 1
C
-0,05 -0,05
B B
-0,1 -0,1
Normalized frequency Normalized frequency

Fig. 17. Results for the apnn edge: (a) MTF curves from 2016 runs, (b) MTF curves from 2017
runs, (c) deviation to the model from 2016 runs, (d) deviation to the model from 2017 runs.
Table 4. Results at Nyquist frequency for apnn edge (without D for 2017).
ACT 2016 ALT 2016 ACT 2017 ACT 2017

Mean 0.20 - 0.21 -

Standard deviation 0.03 - 0.01 -
Max-min 0.07 - 0.02 -
For this case, the curves are in agreement except for D. It appears that this approach has
difficulties when processing this type of edge (compare to cgpnn below) which was not
solved with reprocessing in 2017.
5.4 Results for cgpnn edge
The 2016 and 2017 MTF curves and the discrepancy for the cgpnn edge are presented in Fig.
18. Table 5 shows the corresponding MTF values at Nyquist frequency.
 (a) cgpnn (2016) (b) cgpnn (2017)

1 1
G G
0,8 0,8
F
D
0,6 0,6
ACT MTF

ACT MTF
D
C
0,4 0,4 C
B
B
0,2 0,2
Model Model
0 0
0 0,2 0,4 0,6 0,8 1 0 0,2 0,4 0,6 0,8 1
Normalized frequency Normalized frequency

(c) Deviation to the model (d) Deviation to the model

cgpnn (2016) cgpnn (2017)

0,1 0,1
ACT MTF difference

ACT MTF difference

G G
0,05 0,05
F
D
0 0 D
0 0,2 0,4 0,6 0,8 1 C 0 0,2 0,4 0,6 0,8 1
C
-0,05 -0,05
B B
-0,1 -0,1
Normalized frequency Normalized frequency

Fig. 18. Results for the cgpnn edge: (a) MTF curves from 2016 runs, (b) MTF curves from
2017 runs, (c) deviation to the model from 2016 runs, (d) deviation to the model from 2017
runs.
Table 5. Results at Nyquist frequency for cgpnn edge
ACT 2016 ALT 2016 ACT 2017 ACT 2017

Mean 0.13 - 0.14 -

Standard deviation 0.01 - 0.01 -
Max-min 0.03 - 0.01 -
For this case, there is a very good agreement among all results. The only significant
deviation comes from C above the normalized frequency 0.8.
5.5 Results for 14oct_P3 edges
The 2016 and 2017 MTF curves and the discrepancy for the 14oct_P3 edge are presented in
Fig. 1920. Table 6 gives the 2016 and 2017 MTF values at Nyquist frequency.
 (a) 14oct_P3 (2016) (b) 14oct_P3 (2017)

1 1
G G
0,8 0,8 F
E
E
ACT MTF

ACT MTF
0,6 D 0,6
D
0,4 C 0,4
C
B
0,2 0,2 B
Mean Mean
0 0
0 0,2 0,4 0,6 0,8 1 0 0,2 0,4 0,6 0,8 1
Normalized frequency Normalized frequency

(c) Deviation to the mean (d) Deviation to the mean

14oct_P3 (2016) 14oct_P3 (2017)

0,1 0,1
ACT MTF difference

G G

ACT MTF difference

0,05 E
0,05 F
E
0 D 0
0 0,2 0,4 0,6 0,8 1 0 0,2 0,4 0,6 0,8 1 D
C
-0,05 -0,05 C
B B
-0,1 -0,1
Normalized frequency Normalized frequency

Fig. 19. Results for the 14oct_P3 edge: (a) MTF curves from 2016 runs, (b) MTF curves from
2017 runs, (c) deviation to the mean from 2016 runs, (d) deviation to the mean from 2017 runs.
Table 6. Results at Nyquist frequency for 14oct_P3 edges (without D for 2017)
ACT 2016 ALT 2016 ACT 2017 ALT 2017

Mean 0.12 0.14 0.12 0.13

Standard deviation 0.02 0.03 0.01 0.02
Max-min 0.05 0.07 0.04 0.04
For this first case with actual satellite image data, there is good agreement between the
results except for D at low frequencies. Once again, it looks like a problem related to the type
of edge and resulting MTF shape.
5.6 Results for 15aug_P3 edges
The 2016 and 2017 MTF curves and the discrepancy for the 15aug_P3 edge are shown in Fig.
20. Table 7 gives the corresponding MTF values at Nyquist frequency.
 (a) 15aug_P3 (2016) (b) 15aug_P3 (2017)

1 1
G G
0,8 0,8
F
D
0,6 0,6
ACT MTF

ACT MTF
D
C
0,4 0,4 C
B
0,2 0,2 B
Mean Mean
0 0
0 0,2 0,4 0,6 0,8 1 0 0,2 0,4 0,6 0,8 1
Normalized frequency Normalized frequency

(c) Deviation to the mean (d) Deviation to the mean

15aug_P3 (2016) 15aug_P3 (2017)

0,1 0,1

ACT MTF difference

G
ACT MTF difference

G
0,05 0,05 F
D
0 0 D
0 0,2 0,4 0,6 0,8 1 C 0 0,2 0,4 0,6 0,8 1
C
-0,05 -0,05
B B
-0,1 -0,1
Normalized frequency Normalized frequency

Fig. 20. Results for the 15aug_P3 edge: (a) MTF curves from 2016 runs, (b) MTF curves from
2017 runs, (c) deviation to the mean from 2016 runs, (d) deviation to the mean from 2017 runs.
Table 7. Results at Nyquist frequency for 15aug_P3 edges (without D for 2017).
ACT 2016 ALT 2016 ACT 2017 ACT 2017

Mean 0.12 0.13 0.12 0.13

Standard deviation 0.02 0.02 0.01 0.02
Max-min 0.05 0.06 0.02 0.04
For this actual case, there is a good agreement between the results except for D at low
frequencies, between 0.1 and 0.2.
9. Conclusion
A test of several algorithms derived from the widely used edge method has been performed.
For the test, a set of images of edges was created, mixing both simulations and actual images,
and was made available to the participants without any information about the edges, PSF, or
MTF. Each participant processed the edge data set to estimate the MTF curve for each edge.
Thus, the first comparison of the MTF corresponds to blind test results. A second one
occurred one year later which allowed possible improvements to the processes or the inputs.
For either the first or for the second comparison, none of the participants was able to
always produce the best estimate (the closest to the expected one for simulation or the closest
to the mean of the measurements for the actual cases). This experiment showed that, in some
cases, the error or inaccuracy may be MTF shape dependent. Thus, a validation should
include several MTF shapes. It also stressed that the results may seem to be consistent when
looking at MTF value, at Nyquist frequency, but are not always consistent for the whole
curve. Indeed, for some participants, the quality of the assessment depends strongly on the
shape of the MTF curve. All participants presented their methods and no theoretical problems
were found. The explanation of some unexpected results could possibly be a bug in the
software or some inadequate inputs. This emphasizes that a full understanding of the method
is required to obtain reliable results.
To extend this comparison study, a reference data set composed of edge image and
corresponding MTF curves will be built. It is planned to give access to edge images through
the CEOS CalVal Portal (http://calvalportal.ceos.org/) and make them available to a broader
audience. In order to promote blind testing as well as to enhance and enlarge this reference
data set, it is planned to deliver the reference MTF curves upon receipt of MTF curves from
user. Moreover, users are invited to propose new images to enlarge the data set. This paper
associated with the reference data set can be seen as a new tool to implement and/or check
MTF measurement relying on the slanted edge method.
Acknowledgments
We wish to strongly acknowledge Airbus DS, Digital Globe and CSIR for providing the edge
images and ESA for hosting the data sets on the CEOS WGCV CalValPortal.
References
1. M. R. B. Forshaw, A. Haskell, P. F. Miller, D. J. Stanley, and J. R. G. Townshend, “Spatial resolution of
remotely sensed imagery - a review paper,” Int. J. Remote Sens. 4(3), 497–520 (1983).
2. Satellite imagery. From acquisition principles to processing of optical images for observing the earth (Cepadues
Editions 2012).
3. W. H. Carnahan and G. Zhou, “Fourier transform techniques for the evaluation of the Thematic Mapper line
spread function,” Photogram. Eng. Rem. S. 52, 639–648 (1986).
4. K. Maeda, M. Kojima, and Y. Azuma, “Geometric and radiometric performance evaluation methods for marine
observation satellite-1 (MOS-1) verification program (MVP),” Acta Astronaut. 15(6-7), 297–304 (1987).
5. T. Choi, “IKONOS satellite in orbit, modulation transfer function measurement using edge and pulse methods”,
MSc Thesis, South Dakota State University (2002).
6. D. Léger, F. Viallefont, P. Déliot, and C. Valorge, “On-orbit MTF assessment of satellite cameras,” in Post-
launch calibration of satellite sensors (Morain and Budge ed., Taylor and Francis group, 2004).
7. F. Gascon, C. Bouzinac, O. Thépaut, M. Jung, B. Francesconi, J. Louis, V. Lonjou, B. Lafrance, S. Massera, A.
Gaudel-Vacaresse, F. Languille, B. Alhammoud, F. Viallefont, B. Pflug, J. Bieniarz, S. Clerc, L. Pessiot, T.
Trémas, E. Cadau, R. De Bonis, C. Isola, P. Martimort, and V. Fernandez, “Copernicus Sentinel-2A Calibration
and Products Validation Status,” Remote Sens. 9(6), 584 (2017).
8. C. L. Norton, G. C. Brooks, and R. Welch, “Optical and modulation transfer function,” Photogram. Eng. Rem. S.
43, 613–636 (1977).
9. F. Lei and H. J. Tiziani, “A comparison of methods to measure the modulation transfer function of aerial survey
lens systems from image structures,” Photogram. Eng. Rem. S. 54, 41–46 (1988).
10. ISO 12233:2017(E), “Photography – Electronic still picture imaging – Resolution and spatial frequency
responses”, International Organization for Standardization, Geneva, CH (2017).
11. D. Williams and P. D. Burns, “Evolution of slanted edge gradient SFR measurement,” IS&T/SPIE Electronic
Imaging, International Society for Optics and Photonics, San Francisco, California, USA (2014).
12. K. Kohm, “Modulation transfer function measurement method and results for the Orbview-3 high resolution
imaging satellite,” Proceedings of ISPRS, Istanbul, Turkey (2004).
13. M. Estribeau and P. Magnan, “Fast MTF measurement of CMOS imagers using ISO 12233 slanted-edge
methodology”, in Proc. SPIE 5251, Detectors and Associated Signal Processing, 243–253 (2004).
14. H. Hwang, Y.-W. Choi, S. Kwak, M. Kim, and W. Park, “MTF assessment of high resolution satellite images
using ISO 12233 slanted-edge method,” Proc. SPIE 7109, 710905 (2008).
15. F. Viallefont-Robinet and D. Léger, “Improvement of the edge method for on-orbit MTF measurement,” Opt.
Express 18(4), 3531–3545 (2010).
16. B. Tatian, “Method for Obtaining the Transfer Function from the Edge Response Function,” J. Opt. Soc. Am.
55(8), 1014 (1965).
17. S. Saunier, P. Goryl, G. Chander, R. Santer, M. Bouvet, B. Collet, A. Mambimba, and S. Kocaman Aksakal,
“Radiometric, geometric and image quality assessment of the ALOS AVNIR-2 and PRISM sensors,” TGRS
48(10), 3855–3866 (2010).
18. U. M. Leloglu and E. Tunali, “On orbit modulation transfer function estimation for Bilsat Imagers,” in Proc.
ISPRS, Paris, France, B, T04–18.1 (2006)
19. F. van den Bergh, “On the rendering of synthetic images with specific point spread functions,” 23rd Annual
Symposium of the Pattern Recognition Association of South Africa (PRASA), Pretoria, South Africa (2012).