Artifacts appear when using an external high res (10m) DEM with Range Doppler Terrain Correction

Hello there,

I am trying to to terrain correct S1 images. No matter whether they are SLC or GRD, I have the same problem:

If I choose an autodownload DEM (SRTM 3s or SRTM 1s) I get a perfectly terrain corrected image.
However, if I try to use and external 10m DEM, a series of narrow horizontal and vertical stripes appear. They can be better seen at the local incidence angle band. But if I choose a lower res 25m external DEM, such artifacts disappear.

If there any limitation to the resolution of the external DEMs to use? Why does it happen? Has anyone encountered a similiar issue?

The attached images correspond to a very hilly area scene after applying the Range Doppler terrain correction

INCIDENCE ANGLE FOR 10m External DEM, 25m External DEM, SRTM 1sec:

BETA0 FOR 10m External DEM, 25m External DEM, SRTM 1sec

However, such artifacts do not propagate to the elevation band:

Changing the interpolation method of the external DEM could help. Have you tried if it makes a difference?

Hello there,

Yes, I did try several methods, but the artifacts will still appear. With nearest neighbourgh interpolation their impact is even higher, and with the rest of methods available the artifacts will be similar…

At some point I check whether it has something to do with the ratio between DEM resolution and choosen resulting image resolution, but seems not to be related. With the 25m external DEM, choosing 5m, 10m or 25m will not result in artifacts. With the 10m external DEM, artifacts will appear…

Is the projection of the DEM WGS84?

You can resample it to the target resolution outside SNAP and try if it changes.

Yes, projection is WGS84 , with geographic coordinates and heigh in meters. The thing is the 25m DEM and the 10m DEM are identical, exept for the resolution. However, the 25m works fine, whereas the 10m one shows the artifact…

After TC with both 10m and 25m DEM, the generated Elevation bands are, for both DEM, flawless…So somehow the processing is indeed able to read the heigh information, but maybe it fails when computing local normals and/or local incidence angles for the 10m DEM…

GRADUAL INTERPOLATION FOR THE SAME DEM 10m, 15m, 25 m:

The artifacts gradually dissappear…

So it seems somehow a 2D harmonic-like artifact is appearing , and its impact is bigger the smaller the external DEM resolution is.
Maybe there is some kind of spatial filtering at some point, and Nyqist condition is not being fullfiled for the 10m DEM?

Browsing through the source code, I see tha the TC uses methods from the SARGeocoding class:

There are two functions for computing the local incidence angle, named computeLocalIncidenceAngle, each of them with their own set of invoking parameters.

The 2nd one uses both localDEM and DEM as input parameters. I was wondering if the artifact could somehow originate here, when the altitude alt variable takes values alternatively from DEM or from localDEM:

public static void computeLocalIncidenceAngle(
final LocalGeometry lg, final double demNoDataValue, final boolean saveLocalIncidenceAngle,
final boolean saveProjectedLocalIncidenceAngle, final boolean saveSigmaNought, final int x0,
final int y0, final int x, final int y, final double[][] localDEM, final double[] localIncidenceAngles,
final TileGeoreferencing tileGeoRef, ElevationModel dem) throws Exception {

    // Note: For algorithm and notation of the following implementation, please see Andrea's email dated
    //       May 29, 2009 and Marcus' email dated June 3, 2009, or see Eq (14.10) and Eq (14.11) on page
    //       321 and 323 in "SAR Geocoding - Data and Systems".
    //       The Cartesian coordinate (x, y, z) is represented here by a length-3 array with element[0]
    //       representing x, element[1] representing y and element[2] representing z.
    final int yy = y - y0;
    final int xx = x - x0;
    final int maxX = localDEM[0].length - 1;
    final int maxY = localDEM.length - 1;
    final int numN = 3;
    final GeoPos geo = new GeoPos();
    Double alt;
    double rightPointHeight = 0, leftPointHeight = 0, upPointHeight = 0, downPointHeight = 0;
    int cnt = 0;
    for (int n = 0; n < numN; ++n) {
        if (xx + n > maxX) {
            tileGeoRef.getGeoPos(xx + n, yy, geo);
            alt = dem.getElevation(geo);
        } else {
            alt = localDEM[yy][xx + n];
        }
        if (!alt.equals(demNoDataValue)) {
            rightPointHeight += alt;
            ++cnt;
        }
    }

I am not familiar with how this code works, so it is just a wild guess…