Automated DEM generation from SPOT imagery
Extraction of seed points
We generate seed points automatically based on some constraints. For every point (X,Y) on the earth
surface, the radius R of the search window is set. If the initial altitude H at (X,Y) can be estimated, then we
can search for the best position of its corresponding point in the left image within the range H-R to H+R.
Then according to the point’s position (X,Y) and its altitude H-R and H+R, we can calculate its position on the
left image and the right image. The cross correlation score is then calculated. So a sequence of correlation
score is achieved. We called this sequence of correlation scores a correlation curve (figure 1).
Through analyzing this correlation curve, we use several constraints to extract seed points. The constraints
are set such that only points with correlation values greater than a predetermined threshold and whose
correlation curves have one and only one peak not located at the boundary are chosen (e.g. figure 1 þ, þ).
The texture index of gray variation of the matching patches must also lies above a threshold. This constraint
is necessary to eliminate spurious matching in areas whose texture index is very small (e.g. water), but with
a high correlation score.
Altitude Estimating
Once a seed point has been extracted, its neighbor points’ altitude can be estimated. Suppose that point i is
a seed point, and its altitude is Hi. Point j is its neighbor point. If the gradient of point i to point j is K, and the
distance between i and j is S, we may estimate that the altitude range of point j is around Hi -S*K to Hi+S*K.
Thus how to choose the gradient limit is the key to estimating pixels’ altitudes. Generally a little altitude range
is expected. The less the gradient, the less is the altitude range. According to Lu et. al. [3], the gradient of
point i to point j should meet the qualification of a smooth terrain. They consider the terrain to be not smooth
such that point i and point j are not compatible if the gradient of point i to point j is greater than 1. This
assumption is theoretically correct. However, if there is some error in Hi, the estimated altitude range (Hi-S*K,
Hi+S*K) may not be correct. The risk that the estimated altitude range does not contain the correspondence
point will then increase.
Once the altitude range is wrong, no matter how accurate the image matching is, it is not possible to get the
correct matching result. A gradient limit of 1 may not be suitable and a gradient limit of 2 is used by some
researchers [11]. As the gradient limit increases, the estimated altitude range will become greater and
greater. The increasing altitude range will result in an increase in the probability of false matching which is
not desirable. In fact, no matter what gradient limit is used, we do not know whether the estimated altitude
range is correct. The only thing we can do is to minimize the probability of wrong estimation and the
probability of incorrect matching.
In order to reduce the probability of wrong estimation, we adopt a muti-point estimating schedule. For one
seed point, all its neighbor points’ altitudes are estimated. For each point to be estimated, if there are more
than 1 seed points around it, then the weighted mean value of all the seed points’ altitude is used.
Geometric Rectification to Matching Patch
Because of different sensor attitudes, the same matching patch has different shapes on the left image and
the right image. For a rectangular matching patch in the left image, its corresponding matching patch is no
longer rectangular (figure 2). The distortion mainly is caused by rotation and different image scale. In order to
reduce the matching error and increase the stability of image matching, geometric rectification is performed
on the matching patch before matching. For a pixel on the left image, we may use equation (2) to calculate
its epipolar curve on the right image [7],
where i and j are respectively the row and column indices of a pixel on the left image, i’ and j’ are row and
column indices of the corresponding pixel on the right image. Ai are coefficients based on the sensor model.
A short segment of this epipolar curve can be approximated as a straight line. The rotation angle can then be
calculated using two points on this epipolar curve. Euation (3) is used to calculate the image ratio [2].
where a is the incidence angle of left image and b is the incidence angle of right image respectively (figure 3).