| GISdevelopment.net ---> AARS ---> ACRS 1994 ---> Mapping from Space |
DEM Generation using JERSI
OPS Stereo Data
Mitsuharu Tokunaga, Seiichi
Hara
Central Computer Services Co., Ltd.
6-41-10, Kameido, Koto-Ku, Tokyo 136, Japan
TEL + 81-3-5626-7739 FAX + 81-5626-7792
E-Mail : PXF06057@niftyserve.or.jp
Shigetoshi Ishikawa
Earth Remote Sensing Data Analysis Centre
3-12-1, Kachidoki, Chuo-ku, Tokyo 104, Japan
Shunji Murai
Asian Institute of Technology
P.O. Box 2754, Bangkok 10501, Thailand
Central Computer Services Co., Ltd.
6-41-10, Kameido, Koto-Ku, Tokyo 136, Japan
TEL + 81-3-5626-7739 FAX + 81-5626-7792
E-Mail : PXF06057@niftyserve.or.jp
Shigetoshi Ishikawa
Earth Remote Sensing Data Analysis Centre
3-12-1, Kachidoki, Chuo-ku, Tokyo 104, Japan
Shunji Murai
Asian Institute of Technology
P.O. Box 2754, Bangkok 10501, Thailand
Abstract
The Method for DEM generation using JERS1 OPS stereo data was developed. JERS1 satellite was launched by NASDA in 1992. The ground resolution of OPS data is about 18m on the cross track and about 24m on the along track. Also the OPS is able to perform stereoscopic observation from the same orbit with 0.3 B/H ratio.
DEM were generated the following method. (1) The orientation elements were determined without GCP. (2) DEM were by automatically stereo matching. (3) DEM were adjusted using a few GCP.
The Accuracy of DEM was 32.1m on the test site. Therefore, the generated DEM has the potentiality of drawing 1:200,000 scale topographic map.
1. Introduction
The Japanese Earth Resources Satellite-1 (JERS1) that was development by the National Space Development Agency of Japan (NASDA) and Ministry of International Trade and Industry (MITI) was carried on a Synthetic Aperture Radar (SAR) and an Optical Sensor (OPS).
The OPS is passive linear array sensor that measures in the visible, near infrared, and short wave infrared. Each spectral band has a 4096 elements linear CCD. The band wave length shows in Table 1.
Table 1 Observation band
The Visible and Near Infrared Radiometer (VNIR) can produce stereoscopic images using band3 and band4 which looks 15.3 degrees forward in the orbit plane. It means that Base and height ratio (B/H) is 0.3
The swath width is about 75km and the ground resolution is about 18.3m on the cross track and about 24.2m on the along track.
This stereoscopic system has the advantage of the stereo matching as compared with the side looking observation system. Because the stereoscopic images is observed almost same time. Therefore it is expected that mismatching points are decreased.
The purpose of this study is to develop the method of DEM generation without orientation using the GCPs. Because it is difficult to get many GCPs in our target area where isn't prepared to the topographic map and the DEM generation also can't be automatically calculated by using GCPs. Automatic generation was given priority over accuracy of DEM.
2. Method of DEM generation
2.1 Orientation Elements
The Satellite's orbit and orientation elements were calculated from ephemeris data to be recorded in CCT without the GCPs.
The satellite's orbit are recorded to the following approximate of degree three.
X = a0 + a1 t + a2 t2
+ a3 t3
Y = b0 + b1 t + b2 t2 + b3 t3
Z = c0 + c1 t + c2 t2 + c3 t3
where :
X,Y,Z: the location of satellite in the earth fixed coordinate system
Xaxis is from centre of the earth to the equator through meridian zero.
Zaxis is along the earth's spin axis.
Y axis is the orthogonal axis on XZ plane by a right hand coordinate system.
an , bn , cn : coefficient
t : progressing time from a standard time
and the satellite's orbit are recorded to the following approximate Fourier series.
where :
q : yaw, roll, pitch angle.
an , bn : coefficient
T : approximate width on the time (32 sec)
t : progressing time from a standard time
When the setellite's orientation on the image coordinate is calculated, t is defined as the function of image line number.
t = F(L)
where :
L : line number on the image
t : the progressing time on the L line.
The geometry of observation by a linear array sensor is shown in Fig. 1. There kind of coordinate system were defined as the earth fixed coordinate system, the satellite orbit coordinate system and the satellite fixed coordinate system.
Figure 1 Geometry of Observation
The earth fixed coordinate system is defined that the origin is the center of the earth, Xe axis is from center of the earth to the equator through meridian zero, Ze axis is along the earth's spin axis and Ye axis is the orthogonal axis on XeZe plane by a right hand coordinate system.
The satellite orbit coordinate system is defined that the origin is the center of the satellite, Xs axis is the flying direction, Zs axis is the nadir direction and Ys axis is the orthogonal axis on XsZs plane a right hand coordinate system.
This coordinate is calculated from the location and velocity of satellite. The satellite fixed coordinate system is defined as the local attitude on the satellite orbit coordinate system. This coordinate is recorded as the ancillary data in the CCT.
The relationship of the each coordinate system are shown the follows.
Where :
(XE,YE, ZE) : the target point on the ground in the earth fixed coordinate.
(X0,Y0, Z0) : the satellite location in the earth fixed coordinate.
(X1,Y1, Z1) : the target point on the ground in the satellite fixed coordinate.
(X2,Y2, Z2) : the target point on the ground in the satellite fixed coordinate with view angle
R (K0, j0, w0) : The rotation matrix from the earth fixed coordinate to the satellite orbit coordinate.
R (K1, j1, w1) : The rotation matrix from the satellite orbit coordinate to the satellite fixed coordinate.
q : the observation angle, q=0 (nadir view) q=15.33 degree (forward view)
As the center of projection, the target point in the image and the target point on the ground are alignment, the following equation are satisfied from satisfied from collinear condition.
Where :
f : the sensor focus ( = 2153mm )
Further; the relationship of the sensor coordinate and the image coordinate shows the following equation.
Where :
x y : the location of the target point in the sensor coordinate x is the constant zero, because the sensor is a linear array sensor.
P, L : the number of the target point in the image coordinate
N : the number of the CCD sensor ( = 4096 )
D : an interval of a CCD sensor ( = 7 mm)
2.2 Stereo Matching
The stereo matching were carried out by using the area correlation method on the two phases.
The approximate matching was carried out using the low frequencies images in the first phase. The low frequencies images were the two value images transformed by the difference from the original image to the local average.
The detail matching was carried out using the high frequencies images in the second phase. The high frequencies images were the enhanced images by added the difference from the original image to the local average. The next equation shows the calculation of the average.
Where :
p(x, y) : original image value
m(x, y) : local average
N : window size for local average
the first phase : 25 - 35 the second phase 19: - 15
The neighborhood of matching point detected in the first phase were researched using the correlation window.
This method has the advantage that the size of correlation window is more smaller than other method, so that the computation time is reduced.
3 Generation of DEM
3.1 Test site
The test site is shown in table 2
Table 2 Test Site
3.2 Evaluation for Generated DEM
The 28 ground control points selected in the test site for the conparison the generated DEM with the ground control points.
These points were got the absolute elevation using 1/50,000 topographical map.
3.3 Results
The results of the generated DEM are shown in Table 3 and Fig. 2 the generated DEM had the offset with vias and the gain with trend by compared with GCPs.
Table 3 Calculation Accuracy of Elevation
Figure 2 Contour Lines Generated from JERS - 1 OPS Stereo Data (left) and Those of 1/200000 Map (right)
It's for this reson that the DEM were generated by orientation without GCP. Therefore the relative elevation error of DEM as well as the absolute error were evaluated.
The bias with offset was defined as the difference from the average of DEM to absolute elevation data. So the trend with gain was defined as the incline of plane, because it is assumed that DEM is inclined along a plane.
That plane is defined as follows.
Dhxy = aX + bY + c
Where :
Dhxy : the elevation error in the target point ( X, Y, )
a,b,c : coefficients are calculated by the method of least squares.
4. Conclusion
The absolute elevation error weren't got good accuracy. It's for this reson that the mis matching points were existed in the clout area and the orientation recorded in CCT hasn't so high accuracy. The orientation errors weren't still investigated. Therefore, the elevation error effected by the orientation error weren't also investigated.
However, the relative elevation accuracy are as same as the theoetical value that calculated by the ground resolution devided (B/H ratio)*2. [the theoretical value is 30.5 m - 40.3 m ] if some GCPs are got in the target area, the high accuracy DEM are generated b corrected the bias and the trens. But it isn't nessesary to go many GCPs.
The accuracy of generated DEM has the potentiality of drawing 1:200,000 scale topgraphic map.
Also it isn't nessesary to be operated by a man, Therefore DEM is generated automatically. It's useful to generated the DEM for global data.
Reference
The Method for DEM generation using JERS1 OPS stereo data was developed. JERS1 satellite was launched by NASDA in 1992. The ground resolution of OPS data is about 18m on the cross track and about 24m on the along track. Also the OPS is able to perform stereoscopic observation from the same orbit with 0.3 B/H ratio.
DEM were generated the following method. (1) The orientation elements were determined without GCP. (2) DEM were by automatically stereo matching. (3) DEM were adjusted using a few GCP.
The Accuracy of DEM was 32.1m on the test site. Therefore, the generated DEM has the potentiality of drawing 1:200,000 scale topographic map.
1. Introduction
The Japanese Earth Resources Satellite-1 (JERS1) that was development by the National Space Development Agency of Japan (NASDA) and Ministry of International Trade and Industry (MITI) was carried on a Synthetic Aperture Radar (SAR) and an Optical Sensor (OPS).
The OPS is passive linear array sensor that measures in the visible, near infrared, and short wave infrared. Each spectral band has a 4096 elements linear CCD. The band wave length shows in Table 1.
| Band | Wave Length (mm) |
| 1 | 0.52 - 0.60 |
| 2 | 0.63 - 0.69 |
| 3 | 0.76 - 0.86 |
| 4 | 0.76 - 0.86 |
| 5 | 1.60 - 1.71 |
| 6 | 2.01 - 2.12 |
| 7 | 2.13 - 2.25 |
| 8 | 2.27 - 2.40 |
The Visible and Near Infrared Radiometer (VNIR) can produce stereoscopic images using band3 and band4 which looks 15.3 degrees forward in the orbit plane. It means that Base and height ratio (B/H) is 0.3
The swath width is about 75km and the ground resolution is about 18.3m on the cross track and about 24.2m on the along track.
This stereoscopic system has the advantage of the stereo matching as compared with the side looking observation system. Because the stereoscopic images is observed almost same time. Therefore it is expected that mismatching points are decreased.
The purpose of this study is to develop the method of DEM generation without orientation using the GCPs. Because it is difficult to get many GCPs in our target area where isn't prepared to the topographic map and the DEM generation also can't be automatically calculated by using GCPs. Automatic generation was given priority over accuracy of DEM.
2. Method of DEM generation
2.1 Orientation Elements
The Satellite's orbit and orientation elements were calculated from ephemeris data to be recorded in CCT without the GCPs.
The satellite's orbit are recorded to the following approximate of degree three.
Y = b0 + b1 t + b2 t2 + b3 t3
Z = c0 + c1 t + c2 t2 + c3 t3
where :
X,Y,Z: the location of satellite in the earth fixed coordinate system
Xaxis is from centre of the earth to the equator through meridian zero.
Zaxis is along the earth's spin axis.
Y axis is the orthogonal axis on XZ plane by a right hand coordinate system.
an , bn , cn : coefficient
t : progressing time from a standard time
and the satellite's orbit are recorded to the following approximate Fourier series.
where :
q : yaw, roll, pitch angle.
an , bn : coefficient
T : approximate width on the time (32 sec)
t : progressing time from a standard time
When the setellite's orientation on the image coordinate is calculated, t is defined as the function of image line number.
t = F(L)
where :
L : line number on the image
t : the progressing time on the L line.
The geometry of observation by a linear array sensor is shown in Fig. 1. There kind of coordinate system were defined as the earth fixed coordinate system, the satellite orbit coordinate system and the satellite fixed coordinate system.
Figure 1 Geometry of Observation
The earth fixed coordinate system is defined that the origin is the center of the earth, Xe axis is from center of the earth to the equator through meridian zero, Ze axis is along the earth's spin axis and Ye axis is the orthogonal axis on XeZe plane by a right hand coordinate system.
The satellite orbit coordinate system is defined that the origin is the center of the satellite, Xs axis is the flying direction, Zs axis is the nadir direction and Ys axis is the orthogonal axis on XsZs plane a right hand coordinate system.
This coordinate is calculated from the location and velocity of satellite. The satellite fixed coordinate system is defined as the local attitude on the satellite orbit coordinate system. This coordinate is recorded as the ancillary data in the CCT.
The relationship of the each coordinate system are shown the follows.
Where :
(XE,YE, ZE) : the target point on the ground in the earth fixed coordinate.
(X0,Y0, Z0) : the satellite location in the earth fixed coordinate.
(X1,Y1, Z1) : the target point on the ground in the satellite fixed coordinate.
(X2,Y2, Z2) : the target point on the ground in the satellite fixed coordinate with view angle
R (K0, j0, w0) : The rotation matrix from the earth fixed coordinate to the satellite orbit coordinate.
R (K1, j1, w1) : The rotation matrix from the satellite orbit coordinate to the satellite fixed coordinate.
q : the observation angle, q=0 (nadir view) q=15.33 degree (forward view)
As the center of projection, the target point in the image and the target point on the ground are alignment, the following equation are satisfied from satisfied from collinear condition.
Where :
f : the sensor focus ( = 2153mm )
Further; the relationship of the sensor coordinate and the image coordinate shows the following equation.
Where :
x y : the location of the target point in the sensor coordinate x is the constant zero, because the sensor is a linear array sensor.
P, L : the number of the target point in the image coordinate
N : the number of the CCD sensor ( = 4096 )
D : an interval of a CCD sensor ( = 7 mm)
2.2 Stereo Matching
The stereo matching were carried out by using the area correlation method on the two phases.
The approximate matching was carried out using the low frequencies images in the first phase. The low frequencies images were the two value images transformed by the difference from the original image to the local average.
The detail matching was carried out using the high frequencies images in the second phase. The high frequencies images were the enhanced images by added the difference from the original image to the local average. The next equation shows the calculation of the average.
Where :
p(x, y) : original image value
m(x, y) : local average
N : window size for local average
the first phase : 25 - 35 the second phase 19: - 15
The neighborhood of matching point detected in the first phase were researched using the correlation window.
This method has the advantage that the size of correlation window is more smaller than other method, so that the computation time is reduced.
3 Generation of DEM
3.1 Test site
The test site is shown in table 2
| Band3 ( nadir ) | Band 4 ( forward ) | |
| GRS | 78 - 248 | 78 - 246 |
| Observation Date | March 18, 1993 | March 18, 1993 |
| Processing Data Level | Level 10 | Level 10 |
| Gain Mode | High Gain | High Gain |
3.2 Evaluation for Generated DEM
The 28 ground control points selected in the test site for the conparison the generated DEM with the ground control points.
These points were got the absolute elevation using 1/50,000 topographical map.
3.3 Results
The results of the generated DEM are shown in Table 3 and Fig. 2 the generated DEM had the offset with vias and the gain with trend by compared with GCPs.
| Positive Maximum Error (m) | Negative Maximum Error (m) | Ave. of Absolute Difference | R.E.M.S. (m) | Bias (m) | Trend X Direction (m/pixel) | Trend Y direction (m/line) | |
| no corre- ction |
1159.2 | 0.0 | 1084.1 | 1084.8 | 1083.1 | 0.0303 | 0.0028 |
| bias | 75.1 | -58.6 | 31.6 | 36.8 | 0.0 | 0.0303 | 0.0028 |
| bias and trend | 82.7 | -53.5 | 22.8 | 30.0 | 0.0 | 0.0 | 0.0 |
| 1 pixel = 18.3m, 1 line = 24.2m | |||||||
Figure 2 Contour Lines Generated from JERS - 1 OPS Stereo Data (left) and Those of 1/200000 Map (right)
It's for this reson that the DEM were generated by orientation without GCP. Therefore the relative elevation error of DEM as well as the absolute error were evaluated.
The bias with offset was defined as the difference from the average of DEM to absolute elevation data. So the trend with gain was defined as the incline of plane, because it is assumed that DEM is inclined along a plane.
That plane is defined as follows.
Where :
Dhxy : the elevation error in the target point ( X, Y, )
a,b,c : coefficients are calculated by the method of least squares.
4. Conclusion
The absolute elevation error weren't got good accuracy. It's for this reson that the mis matching points were existed in the clout area and the orientation recorded in CCT hasn't so high accuracy. The orientation errors weren't still investigated. Therefore, the elevation error effected by the orientation error weren't also investigated.
However, the relative elevation accuracy are as same as the theoetical value that calculated by the ground resolution devided (B/H ratio)*2. [the theoretical value is 30.5 m - 40.3 m ] if some GCPs are got in the target area, the high accuracy DEM are generated b corrected the bias and the trens. But it isn't nessesary to go many GCPs.
The accuracy of generated DEM has the potentiality of drawing 1:200,000 scale topgraphic map.
Also it isn't nessesary to be operated by a man, Therefore DEM is generated automatically. It's useful to generated the DEM for global data.
Reference
- Y. Mukai, T. Sugimura, K. Arai, Generation of DEM by the method of Correcting the Elevations of Ground Control Points using System Corrected SPOT Data, Vol.12, No.3, Journal of Remote Sensing Society of Japan, 1992
- Y. Mukia. T. Sugimura. H. Kimura. Early Evaluation of JERS-1 OPS Data and their Distribution, Vol.12, No.3, Journal of Remote Sensing Society of Japan, 1992.
- H. Hino, H. Ono, Development Status of JERS-1 Mission Equipment Initiaal Checkout of OPS, MDT and MDR, Vol.12, No.3, Journal of Remote Sensing Society of Japan, 1992
- S. Murai and Others, Analytical Photogrammetry, JSPRS, 1988
- S. Murai and R. Shibasaki, Automated Generation of Digital Terrain Model (D.T.M.) using Linear Array Sensor Data, Vol.23, No.3, Photogrammetry and Remote Sensing.