Optimum Interpolation of Imaging Microwave Radiometer Data

.4h.<lract-A technique to interpolate imaging microwave radiometer data is presented as an application of the Backus-Gilhert theory. The interpolation process is optimal in the sense that it attempts to preserve the spatial resolution of the antenna gain function associated with the sampled radiometer data. The technique is applied to the Special Sen­ sor Microwavellmager (SSM/I) satellite data, and is found to enhance the high-resolution features of the imagery. The technique is expected to be a useful tool in creating images for detailed analysis.


T. INTRODUCTION
I N many situations, proper imaging of microwave radi ometer data requires an antenna scan rate and associated data sampling rate consistcnt with collecting all thc spa tial frequency information sensed by the antenna.If the antenna is treated as having a band-limited spatial fre quency response W km -I [1], then the scan and sampling rates should be selected to be greater than or equal to the minimum sampling rate 2 W samples per km often referred to as the Nyquist rate [2].Due to sampling limitations, the amount of data collected may be inadequate to con duct detailed image analysis.In addition, the spatially sampled grid may not be compatible with the desired im age projection coordinates.Thus, interpolating.resam piing, or smoothing of the radiometric data may bccome important parts of the image processing function.
It is, of course, desirable that the interpolation process preserve the spatial resolution of the sampled radiometric data as well as avoid introducing data processing artifacts that could alter or destroy the high spatial frequency con tent of the image.In addition, the interpolation should avoid degrading the radiometric image with excessive amplification of instrument noise appearing in the sam pled data.
The interpolation filter, as it is commonly referred to, may be viewed as "optimal" if the rcsulting interpolated sampling geometry.The intent of this paper is to show that the Backus-Gilbert theory may be applied to closely approximate an optimal linear interpolation filter for im aging microwave radiometer data.
The paper draws upon the theoretical results presented by Stogryn L3 J who applied the Backus-Gilbert method ology to attempt beam shaping of the antenna gain func tion at the sampled data points as well as remove cross polarization coupling.Although the emphasis here is on interpolation, the mathematics of the Backus-Gilbert ap proach is quite general and may be applied to a number of image processing problems.For example, Stogryn [3] notes the application to smoothing a high spatial resolu tion image to match the resulting resolution to that of a specified lower resolution image, and [4] shows the im provement in retrieving sea surface temperature from multifrequency radiometer data when commensurate spa tial resolutions are maintained in the retrievals.This paper does not address the subject of increasing the spatial resolution of the interpolated radiometer data beyond that associated with the antenna 3 dB beam di mneter.Attempts to extract frequency information from an image beyond the standard diffraction limits of the an tenna aperture and obtain higher resolution (i.e., "super resolution") image features has received considerable at tention in the literature [5]- [7].Unfortunately, most methods that offer higher spatial resolution incur rela tively large amplification of sensor noise and disturbing artifacts in the processed image.

n. THEORY
Consider a scanning microwave radiometer sensing the earth's upwelling brightness temperature field shown in Fig. 1.The sampled radiometer output at time t expressed as an antenna temperature TA provides a measurement of the brightness temperature TB incident on the antenna, and it may be expressed as where '�A (t) is the effective pointing direction of the an tenna associated with the radiometer sample at time t, h is the impulse response of the radiometer receiver presam which describes the angular weighting of radiation inci dent on the antenna in direction of unit vector -s with the antenna in direction of unit vectors o It is assumed that h is normalized so that J':'rn dt h = I, and that G is suffi ciently directive so that the angular integration may be restricted to the earth's surface and JE dO G = 1.It is also assumed that G has negligible cross-polararization cou pling or that a correction has been made to remove its effect on TA (see [3]) so that (1) applies to the principal polarization of interest.In general, the earth's brightness temperature field TB depends on the direction s, the posi tion vector p, as well as time.The effective pointing di rection ,I: A(t) may be defined as '�A ( t) = [rn dr' h (t -t') s) t ' ) , (2 To simplify (\). the angular integration is transformed to an area integration on the surface of the earth, and the time and spatial integrations are interchanged to yield where unit vector p points toward area element dA shown in Fig. 1.Expression (3) assumes an idealized spherical model for the earth.The quantity within the brackets in (3) may be identified as the effective areal weighting given by the antenna to the brightness temperature field associ ated with area dA at position p in direction s( t') at time t'.
To proceed much further, additional infonnation is needed concerning the directional and time dependence of TB• Following the development in [3], it will be assumed that over the time period of interest, i.e., the measure ment interval associated with the collection of radiometer temperature upwelling from area dA at p is essentially independent of direction - § over the angular regions of interest.and the time variation of TB may be neglected: As noted in [31.this approximation may be reasonably good for many cases; it neveltheless constitutes a funda mental limitation in any analysis which docs not account for the directional dependence and potential rapid time dependencc of TB.
With (4), (3) simplifies to 7:4(P A ) = r dA G(PA' P) TB (P) (5) JE where position vector P A is associated with direction SA, and G is given by (6) and describes the effective antenna gain function associ ated with the antenna temperature sample TAo It should be noted that G takes into account the motion of the sensor platform, the characteristics of the antenna scan geome try, the effects of the radiometer low-pass filter, and. of course, the actual antenna gain function G.
To preserve the imaging properties of the antenna, the radiometer low-pass filter h must be selected to "match" the spatial frequency response of the antenna G in (6).
This may be accomplished by noting that the low-pass filter operates on the gain function G primarily in the along scan direction.Thus, knowledge of the maximum spatial frequency sensed by G in this direction and the antenna scan rate permits a determination of the maximum fre quency of the wavefoml input to the filter.Of course, the filter should also have as narrow a passband as possible to minimize the effects of sensor noise.For an ideal low pass filtcr with cutoff frequency equal to the product of the maximum spatial frequency sensed by the antenna (km -I ) and the scan rate (km / s ), the effective gain G (6) reduces to (7) which is simply the instantaneous antenna gain function at sample time t.A number of performance parameters must be addressed in the process of selecting a realizable low-pass lilter that adequately approximates (7).These parameters include the distortion introduced by the filter due to nonuniform attenuation and group delay within the desired passband as well as the noise equivalent band width of the filter.It should be noted that it is desirable to select h so that (7) is maintained during the mcasure ments.The resulting sampled data may be smoothed with a narrower filter during image processing if desired.Now suppose that a set of antenna temperature samples TA ( Pi) at positions Pi. i = I .... • N exist and that it is desired to estimate the antenna temperature at position P d which lies within the region containing Pi' The Backus-Gilbert theory may be applied to estimate the an tenna temperature fA ( p d ) at P d using a weighted linear combination of the neighboring antenna temperature sam ples, where {a,} denotes the set of coefficients of the interpo lation fi ltcr.Substituting ( 5) into ( 8) and interchanging the sum and has undesirable sidelobe levels.then a functional form may be selected for J which emphasizes the region of the sidelobes.Depending on the form selected, the reduction of the sidelobe level may result in some widening or slight distortion of the main beam.For the results presented herein, J is set to unity.
The solution for the vector a with components {ai} may be written as where g -I is the inverse of matrix g defi ned by compo nents ( 15) integral operations yields and it and v are vectors with components where Gj may be identified as the interpolated antenna gain function N b UiG (P;, pl.i= 1 Now if T4 (Pd) represents the "true" antenna temper atures associated with P d, then the error in estimate (9) may be written as In the evcnt it is possible to select the set of {a,} so that Gr closely approximates G, then e will be minimized.If this is to occur, it is clear that Gj must represent a prop crly normalized antenna gain function so that As noted in [3], many criteria may be posed to mini mize e, with perhaps the most useful ones presented 1TI terms of the integral Qd = LdAl G (Pd. til -,� a ; G (p;, p)r J(Pd' p) (13)   where J is a specified function, and the set of weighting coefficients {ai} are selected to minimize Q d subject to the normalization constraint (12).J may be selected to emphasize regions of importance within the integration.For example, in the event the interpolated gain function Thc rcadcr may rccognize that criterion (13) with con straint (12) is equivalent to an unbiased minimum mean squared estimator of TA for the situation in which TB is treated as a homogeneous random process [8] having a uniform power spectral density.That is, expressing TB, (18) where < > denotes the expectation operator and OTB is the zero mean random component of TB, then under the above assumption for TB, the mean squared error of e (11) may be written as The first term on the right-hand side of ( 19) is the basis of the estimator and it vanishes if (12) holds.The min imization of the second term is equivale'ht to minimizing (13) with Finally. it should be noted that the instrument noise ap pearing in the sampled antenna temperatures propagates to the interpolated temperature (9).As presented in [3], a minimization criteria may be established to include a lin-car weighted sum of Q d and the variance of the instrument noise in the estimate (9) (see [13, eq. ( 13)]).The ap proach taken here is to select a; which minimize Q d sub ject to (12) without consideration of instrument noise.However, it is important to examine the resulting variance of the interpolated antenna temperatures e2 for each set of coefficients: (20 where E is the covariance matrix of the instrument noise process appearing in the sampled antenna temperatures, and is a function of system noise temperature, predetee tion bandwidth, and receiver low-pass filter.When the noise in each sample is uncorrelated and equal, then e2 is proportional to the sum of the squares of a; times the variance of the instrument noisc.

III. EXAMPLES
As an application of the above theory, consider the problem of interpolating the satellite microwave radiome ter data collected by the Special Sensor Microwave/Ima ger (SSM/I).The SSMII is a seven-channel four-fre quency linearly polarized conically scanning microwave radiometer on board the joint United States Air Forcel Navy Defense Meteorological Satellite Program (DMSP) spacecraft.The satellite operates in a circular sun-syn chronous near polar orbit at a nominal altitude of 833 km and an inclination of 98.8° and orbit period of 102 min.A description of the SSM!I instrument and the environ mental products retrievablc from its data is given in the SSMII User's Guide [9] and the SSMf} ealIVal Final Re port [l0].Suffice it to say that the antenna system, which consists of a single offset parabolic reflector illuminated by a single wide-band multichannel feedhom, produces seven nearly circular coaxial beams with 3 dB beam widths of 1.9°, 1. 06°, and 0.4° with dual linear polari zations at channel frequencies of 19.35, 37.0, and 85 GHz, respectively, and 1.6° with single polarization at 22.235 GHz.During the conical scan, the antenna beams trace out a circular arc on the earth's surface at a constant nominal 53.10 incidence angle across a swath width of 1400 km.The antenna scan rate of 30 rpm and channel data sampling rates result in a 12.5 km sampling grid for the highest resolution 85 GHz channels along the circular arc, and a 12.5 km sampling increment between succes sive scans along the direction of the sub satellite ground track.For the lower resolution channels.the spatial sam pling grid is 25 km along the circular arc and 25 km be tween successive scans.
The conical scan geometry is needed to retrieve envi ronmental products.At the same time thi, geometry pro vides a uniform antenna gain function for each sample antenna temperature.However, the conical scan intro duces a dcgree of complexity in the relative geometry of the sampled temperatures when viewed in a coordinate 8m system on the earth's surface.For example.Fig. 2 illus trates the relative positions of the 85 GHz samples for regions near the start of scan, middle of scan, and end of scan in a latitudeflongitude coordinate system.The cir cles identify the 3 dB contour of the antenna gain function projected on the earth's surface and shows the variation of the relative sample geometry and beam overlapping across the swath.
It should be noted that the SSMfI employs an integrate and dump filter for the radiometer low-pass filter which, in accordance with (6), increases the effective antenna beam diameter primarily in the along scan direction.An integration time of 4.22 ms is used at 85 GHz and trans lates to an azimuthal integration region of 0.8 0. For the remaining channels, thc integration times are 8.44 ms with an 1.6° azimuthal region.The 85 GHz data are sampled at 4.22-ms interval and the remaining channels at 8.44 ms.The increase in beam diamctcr at 85 GHz is just enough to cause the projected 3 dB contour to be nearly circular.The 3 dB contours for the remaining channels are considerably more elliptical with corresponding more beam overlapping of 3 dB contours than shown in Fig. 2.
Since thc integrate and dump filter is not matched to the antenna gain function, the full resolution of the gain func tion ( 7) is not realized in thc along scan direction.The loss in resolution is signifi cant at both 37 and 85 GHz and was implemented in the SSMfI to meet the spacecraft data rate of 3.6 kh / s.The results presented below apply to the effective gain factor ( 6) with the integrated and dump fil ter.
To illustrate the efficacy of the interpolation process presented above, consider the problem of interpolating the 85 GHz antenna temperature at a point located midway between successive scans and midway between along scan samples.The point marked by "X" in Fig. 2 identifies the interpolation point of interest for each of the three scan regions.Since the antenna gains are reasonably directive, the set of neighboring sampled antenna tcmperatures used in the interpolation may be restricted to 3-4 samples in each direction about point X.For definiteness, consider a four X four array.From Fig. 2 it is clear that the neigh boring samples (lying within the dashed line) form a square array near the subsatellite track and parallelogram arrays near the edges of the swath.
Although the "true" antenna temperature at point X is not available, the appropriateness of the interpolated tem perature may be judged on the basis of how closely the interpolated antenna gain function Of (9) approximates the "desired" antenna gain function O.In the results pre sented below, the computation of 0 (6), and, hence 0" include simulations of the satellite motion over a spherical earth, the SSMfI scan geometry, the radiometry integrate and dump low-pass filter, and actual antenna gain func tions as determined from measurements on an antcnna range.A situation similar to that presented above for the 85 GHz data exists when interpolating the 37 GHz antenna temperatures.The problem stems again from a spatial sampling grid of 25 km that is not consistent with the spa tial frequency response of the antenna at 37 GHz.The grid is approximately twice as large as needed to satisfy the Nyquist criterion.However, the interpolation of the 19 GHz data reveal that the interpolated gain functions exhibit reasonably good agreement with the desired gain.
For example, Fig. 10 presents results obtained using a 4    This arises from the fact that the antenna gain function is considerably narrower in the along scan direction, and consequently the 25 km sampling increment is too large to meet the Nyquist criteria in this direction.As noted in the discussion following (13), the level of the sidelobes may be reduced, if desired, by selecting a form for J which emphasizes the region containing the sidelobes.Table I presents the interpolation coefficients for case of Figs.For some situations.this may be a highly desirable fea ture.For example.imagery similar to that of Fig.
data correspond to what would have been measured had the radiometer actually made the measurements.Cer tainly the resulting interpolated data would possess the spatial resolution of the orginally sampled data, avoid am plification of sensor noise.and, in addition, retain any peculiar features associated with the antenna spatial fre quency response (e.g., beam assymetry) or with the data Manuscript received July 17, Iljglj: revised September 21.1989.The author was with the Space Sensing Branch .Naval Center for Space Technology, Naval Research Laboratory, Washington, DC 20375-5000.He is now with Acrojet ElectroSystcms, 1100 West Hollyvale Street, P.O.Box 296, Azusa, CA 91702.IEEE Log Number 9037460.
Fig. I. Geometry and integration elements for antenna temperature (equa tion (I)).

e
= I TA(Pd) -fA(Pd)I I LdA[ G (Pd' p) -G j (Pd, p)]7s(p)l• (1\) u; = L dA G (Pi, ti) Vi = I. dA G(Pi, j:i)G(P".j:i) Jf and the superscript T denotes transpose.(16) (17) Note that in the event the interpolation point P d coin cides with one of the sampled points P b the interpolated temperature reduces to the sampled temperature at Pk since. in this case, Q d vanishes with ak = 1 and a; = 0 (i *-k).

x 4
array of samples to interpolate the 19 GHz (vertical polarization) at a point having the same relative position midway between along scan samples and midway be-

10
except for the sidelobes appearing in the H-cut pattern.

Fig. 14
Fig. 14 clearly shows a considerable improvement over Fig. IS in image quality.This is especially noticeable in the identification of fine details along the Atlantic coast line.the shoreline regions around the Great Lakes.and along the cloud boundaries in the lower part of the image.

Fig. 14
Fig.14presents an image that attempts to preserve the spatial resolution of the original sampled radiometric data.
14 has been employed with an overlay of a world coastline da tabase to analyze and correct geolocation errors appearing in the SSM!! data [12].In addition.the interpolation pro cess described herein has been used to create imagery of tropical cyc loncs.ocean surface wind fields [13 J. and sea ice boundaries [IOJ for detailed analyses, as well as for comparisons of polar sea ice concentration maps with vi sual and infrared imagery from the Advanced Very High Resolution Radiometer (AVHRR) [14J.
809 hg. 15. 85 GHz image (horizontal polarization) with simple pixel repli cation lilter: orbit 4263.April 17. 1988.I V. COKC'I.lJSIOKSAn optimum interpolation technique has been presented to interpolate.resample.and smooth imaging microwave radiometer data.The technique was applied to the satellite Special Sensor Microwave/Imager data.and results were presented which show the efficacy of the methodology in attempting to preserve the spatial resolution of the origi nally sampled data and.hence. in enhancing the quality of the imagery.It is expected that proper tiltering.sam pling.and interpolation of satellite data will continue to receive increasing attention as analysts extract informa tion from high-resolution microwave radiometric im agery .