Topography Fitering

Introduction

We have access to a 1/4 degree digital elevation model that we can use to specify the model's topography. However, when working on a grid coarser than 1/4 degree, we must filter high spatial frequency modes that are not resolved by the grid, but which could give incorrect elevations through aliasing.

The ideal filter would give zero amplitude to all modes that are not resolved well by the grid while retaining full amplitude of those that are. It would also do it efficiently, for cases where the grid is evolving, and would be able to accommodate our potentially irregularly shaped grid.

What is "resolved well"? Certainly any mode with wavelength shorter than 2 x grid-spacing is not resolved. In addition, a numerical model probably does not simulate well dynamic behavior on spatial scales less than about 3 - 4 x grid-spacing, so these scalles should probably be eliminated, too, or at least have small amplitude.

Barnes scheme filtering

One approach to filtering is to obtain topographic elevation at the model's grid points using a weighted-average interpolation scheme, such as the Barnes scheme, which is evaluated here.

The evaluation will be in one dimension, though it could readily be extending to two-dimensional fields. We assume that the 1/4 degree topography field in 1-D is transformed into a fourier representation. The question to ask is then: how much does the fiter suppress the amplitude of each fourier mode? Ideally the filter will remove all short waves (setting their amplitude to zero) while leaving longer wavelengths untouched.

We apply the Barnes scheme to waves of different wavelengths, in each case interpolating to the point of maximum value. Viewing all waves as cosine waves varying in the x-direction with wavelength L, we interpolate to the point x=0. So long as our grid is aligned with the sine waves (grid points symmetric about extreme values of +1 and -1 and the zero crossings), then examining behaviors for x=0 is sufficient for evaluating the smoothing effect of the interpolation. There will also be no phase shift in that minima, maxima and zero crossings will remain in the same locations. For irregularly spaced grids, this will not be so true.

Our grid has spacing dx, which we might view as the 1/4 degree grid spacing of our topography data. Examples are given here of contributions to the summation in the numerator of the Barnes scheme's interpolation equation for two different wavelengths: the shortest possible on the grid ("2-delta" wave) and one with a wavelength 256 times the grid spacing. For our nominal 1/4 degree input data grid, this would be a wavelength of 64 degrees, or roughly 1/6 of a latitude circle.

The first figure is for the case where D=dx. Contributions drop off rapidly with distance. There should be little smoothing of the longest waves, since their interpolated value is determined only by the wave values in the immediate neighborhood and does not include wave values of opposite sign. Note that although continuous lines are drawn, the contributions to the summation come from discrete points that can be identified by the extrema of the L=2dx curve.

The second figure shows the contributions to the summation when D=16*dx. Contributions come from many more grid points. The net effect will be to smooth wave amplitudes, especially for the shortest waves. The third figure zooms in on a portion of the second.

The last figure show the response functions for different choices of D. (Note the logarithmic scale for the x-axis.) These show how much amplitudes of the shortest waves are suppressed by the interpolation's smoothing. The explicit marks on the curves are for wavelengths (in dx units) = 2, 3, 4, 5, 6, 8, 16, 32, 64, 128 and 256. Ideally, the filtering will remove waves the grid cannot resolve but leave waves it can untouched. This is not done perfectly, but nontheless, the longest waves retain nearly all of their amplitude.

One recommendation would be to construct the model's topography using D = grid spacing at each location. This might smooth too much, as wavelengths up to 5dx lose over half their amplitude and only wavelengths longer than about 16dx retain over 90% of their original amplitude. On the other hand, one might argue that the 2dx wave retains too much of its original amplitude

An alternative would be to use D = (grid spacing)/2. Viewing the figure with "2dx" as the grid spacing, the response function would be given by the D=dx curve with marked points suitably rescaled in term of their wavelength. Now wavelengths longer than about 8 x (grid spacing) retain 90% or more of their amplitude, but on the other hand, wavelengths of (1/2) x (grid spacing) may be retained and would alias onto longer modes.


Go to CDGA main page.