Geometric Conditions of Space Imagery for Mapping

GEOMETRIC CONDITIONS OF SPACE IMAGERY FOR MAPPING K. Jacobsen (*), G. Büyüksalih (**), A. Marangoz (**), U. Sefercik (**) and İ. Büyüksalih (**) *University of Hannover *[email protected] ** Zonguldak Karaelmas University ** [email protected] KEY WORDS: Satellite, Optical Image, Orientation, Geometry, Mathematical Models ABSTRACT: From the test field Zonguldak different high and very high resolution optical satellite images are available like TK350, ASTER, KOMPSAT-1, IRS-1C, SPOT 5, KVR1000, IKONOS and QuickBird. The images partially have been achieved as close to original images (level 1A) and partially projected to a plane with constant object height (level 1B). For some images, based on direct sensor information, a good image orientation is given which for some only has to be improved by a shift in X and Y, while for others only rough orientations are distributed. In addition sometimes the orientation has to be improved by additional parameters to compensate systematic geometric effects. Some orientation information of IKONOS- and QuickBird-images is available also as rational polynomial coefficients (RPCs), describing the relation between the image and the object coordinates by a ration of polynomials. The different orientation procedures are described with their advantages and disadvantages. In most cases sub-pixel accuracy was possible. The orientation of some images could be made with different procedures leading to similar results for a sufficient number of well distributed control points. But with a smaller number and also not well distributed control points quite different results have been achieved leading to the clear result that a correct mathematical model, using the available information of the image orientation should be used. This can be done with a geometric reconstruction of the image geometry or sensor oriented rational polynomial coefficients (RPCs) while for the 3D-affinity transformation more and well distributed control points are required. The DLT-method and the terrain dependent RPCs should not be used. 1. INTRODUCTION For mapping purposes only high and very high resolution space images can be used. As a rule of thumb a ground sample distance (GSD), traditionally also named as pixel size on the ground, of approximately 0.1mm in the map scale is required. That means, a map 1 : 50 000 can be generated with a GSD of 0.1mm ∗ 50 000 = 5m. Of course also maps in the scale 1 : 50 000 have been generated with SPOT 1 – 4 images having 10m GSD, but these maps do not contain the same amount of details like topographic maps based on aerial photos. By this reason only space sensors having a GSD not exceeding 15m are respected in this presentation. High resolution space images are distributed as products with different geometry, starting from close to original scenes like QuickBird Basic, OrbView Basic or SPOT, ASTER and IRS 1C level 1A, as projection to a plane with constant height (level 1B, QuickBird ORStandard, IKONOS Geo) or even as rough ortho images (QuickBird Standard). All imaging systems are equipped with a positioning system like GPS, gyroscopes and star sensors. So without control points the geolocation can be determined with accuracies depending upon the system. The information about the sensor orientation is distributed as full data set, including all elements of inner and exterior orientation, as reduced information including only few elements or as polynomial functions describing the relation between the image and the object space. For using the full accuracy potential of space images, control points are required together with the correct mathematical model. Different models are in use, based on available sensor information or just based on control points, independent upon existing orientation information. 2. IMAGING GEOMETRY All the high and very high resolution space images are from CCD-line sensors. The images are generated by the movement of the sensor or sensor view direction. Most of the sensors do not have just one CCD-line but a combination of shorter CCD-lines or small CCD-arrays. The high resolution panchromatic CCD-lines are shifted against each other and the individual colour CCD lines are shifted against this (figure 1). Figure 1: arrangement of CCD-lines in focal plane above: panchromatic below: multispectral The arrangement of the sub-images achieved by the panchromatic CCD-lines belongs to the inner orientation and the user will not see something about it. In addition, usually the matching of the corresponding subimages is in the lower sub-pixel range so that the geometry of the mosaicked image does not show any influence. This may be different for the larger offset of the colour CCD-lines. Not moving objects are fused without any problems during the pan-sharpening process. By theory only in extreme mountainous areas unimportant effects can be seen. This is different for moving objects – the time delay of the colour against the panchromatic image is causing different locations of the grey value and the colour image (figure 2). The colour is always following the grey values. This effect is unimportant for mapping because only not moving objects are used. Figure 2: pansharpened IKONOS image caused by the time delay of the colour imaging, the colour of moving objects are shifted against the grey value image The optical space sensors are located in a flying altitude corresponding to a speed of approximately 7km/sec of the image projected to the ground. So for a ground sampling distance (GSD) of 1m only 1.4msec exposure time is available. The ground sampling distance is the distance of the neighboured pixel centres on the ground, which must not be identical to the size a pixel projected to the ground because of the over- or undersampling. For the user the GSD is the visible object pixel size. 1.4msec is not a sufficient integration time for the generation of an acceptable image quality, by this reason, some of the very high resolution space sensors are equipped with time delay and integration (TDI) sensors. The TDI-sensors used in space are CCD-arrays with a small dimension in flight direction. The charge generated by the energy reflected from the ground will be shifted with the speed of the image motion to the next CCD-element and more charge can be added to the charge collected by the first CCD-element. So a larger charge can be summed up over several CCDelements. There are some limits caused by inclined view directions, so in most cases the energy is summed up over 13 CCD-elements. IKONOS, QuickBird and OrbView-3 are equipped with TDI-sensors while EROS-A and the Indian TES do not have it. They have to enlarge the integration time by a permanent rotation of the satellite during imaging (see figure 3e). Also QuickBird is using this enlargement of the integration time because the sensor originally was planned for the same flying altitude like IKONOS, but with the allowance of a smaller GSD, the flying height was reduced, resulting in a smaller pixel size on the ground and a shorter exposure time which is partially compensated by the change of the view direction during imaging, but with a quite smaller factor like for EROS-A and TES. Fig. 3a: traditional image configuration – fixed orientation in relation to orbit Fig. 3b: yaw control by SPOT 5 HRG - permanent change of view direction across orbit Fig. 3c: flexible view direction – also scan parallel to ground coordinate system Fig. 3d: flexible view direction – scan against orbit possible Fig. 3e: enlargement of integration time with factor b/a by continuous change of view direction Fig. 3f: level 1B, IKONOS Geo, QuickBird OR Standard = projection to plane with constant height The traditional CCD-line sensor satellites, like SPOT 1-4, ASTER, KOMPSAT-1, IRS-1C /1D and the HRS sensor of SPOT-5 do not change the view direction in relation to the orbit during imaging (figure 3a). SPOT 5 is using for the main imaging sensor HRG a yaw correction to compensate the effect of the earth rotation by a permanent change of the view direction across the orbit (figure 3b). The very high resolution and agile satellites like IKONOS, QuickBird, OrbView, EROS-A and TES are able to scan the earth surface in any direction by a permanent change of the satellite orientation. These satellites are equipped with high torque reaction wheels for all axes. If these reaction wheels are slowed down or accelerated, a moment will go to the satellite and it is rotating. No fuel is required for this, only electric energy coming from the solar paddles. The images are distributed as original images, geometrically just improved by the inner orientation, so it looks like from a sensor with one solid CCD-line. These original images are named level 1A in the case of SPOT, KOMPSAT, ASTER or IRS, for QuickBird it is named Basic Imagery. The next step of image product is a projection to a plane with constant height level (see figure 3f). These products are named level 1B, Carterra Geo (IKONOS) or QuickBird ORStandard. The location of the imaged points is depending upon the individual height, so the orientation process has to include also the terrain relief correction based on the individual point height and the nadir angle. For QuickBird also rough orthoimages are distributed as QuickBird Standard – they are related to the GTOPO30-DEM having just a spacing of 30 arcsec, corresponding to approximately 900m at the equator. The mathematical models have to respect the different image geometry. 3. MATHEMATICAL MODELS Different mathematical models are in use, they are ranging from exact geometric reconstruction over reconstruction of the relation image to ground included in polynomials based on correct geometric models (sensor oriented RPCs) over mathematical models not using existing sensor information to polynomial approximations just based on control points. 3.1 Original Images (Level 1A) For some original images the sensor orientation is available for a sufficient number of CCD-lines between which an interpolation is possible. The orientation information includes the projection centre and the attitudes. In addition the focal length and pixel size of the sensor is required. This orientation information is using the direct sensor orientation of the satellite, based on GPS or corresponding positional systems, gyros and star sensors. The accuracy of the satellite direct sensor information is quite different; it is ranging from +/-4m to approximately +/-800m. So control points are required for a geometric refinement and a reliability check. The orientation often is influenced by limited information about the datum of the used coordinate system – the relation of the national net coordinate system to the ITRF-frame used with WGS84. In the Hannover program BLASPO, the image geometry is reconstructed based on the given view direction, the general satellite orbit and few control points. Based on control points the attitudes and the satellite height are improved. The X- and Y-locations are fixed because they are nearly mathematical dependent upon the view direction. In addition two additional parameters for image affinity and angular affinity are required. For the 6 orientation unknowns 3 control points are necessary. More additional parameters can be introduced if geometric problems exist. 3.2 Images projected to plane with constant height (Level 1B) The image orientation of images projected to a plane with constant height (e.g. level 1B, IKONOS Geo, QuickBird ORStandard) can be based on correct mathematic models, but also on approximations. Rational Polynomial Coefficients (RPCs) from the satellite image vendors – they do describe the location of image positions as a function of the object coordinates (longitude, latitude, height) by the ration of polynomials (Grodecki 2001) – see formula 1. These sensor related RPCs are based on the direct sensor orientation of the satellite together with information about the inner orientation and do have an accuracy depending upon the quality of the direct sensor information. Third order polynomials with 20 coefficients are used, so with 80 coefficients the relation of the image coordinates to the object coordinates can be described. The RPCs have to be improved by means of control points. For IKONOS for example a simple shift is usually sufficient, for other sensors or old IKONOS images without the information of the reference height, a two-dimensional affinity transformation of the computed object coordinates to the control points is required. xij = Pi1( X , Y , Z ) j Pi 2( X , Y , Z ) j yij = Pi3( X , Y , Z ) j Pi 4( X , Y , Z ) j xij, yij =scene coordinates X,Y = geographic object coordinates Pn(X,Y,Z)j = a1 + a2∗Y + a3∗X +a4∗Z + a5∗Y∗X + a6∗Y∗Z + a7∗X∗Z + a8∗Y² + a9∗X² + a10∗Z²+ a11∗Y*X*Z + a12∗Y³ +a13∗Y∗X² + a14∗Y∗Z² + a15∗Y²∗X + a16∗X3 + a17∗X∗Z² + a18∗Y²∗Z+ a19∗X²*Z+ a20∗Z3 Formula 1: rational polynomial coefficients Reconstruction of imaging geometry: For the scene centre or the first line, the direction to the satellite is available in the image header data. This direction can be intersected with the orbit of the satellite published with its Keppler elements. Depending upon the location of an image point, the location of the corresponding projection centre on the satellite orbit together with the view direction can be computed. So the view direction from any ground point to the corresponding projection centre can be reconstructed. This method requires the same number of control points like the sensor oriented RPC-solution, that means it can be used also without control points if the direct sensor orientation is accepted as accurate enough or it requires the same additional transformation of the computed object points to the control points like the sensor oriented RPCs. Three-dimensional affinity transformation: It is not using available sensor orientation information. The 8 unknowns for the transformation of the object point coordinates to the image coordinates have to be computed based on control points located not in the same plane. At least 4 well distributed control points are required. The computed unknowns should be checked for high correlation values between the unknowns – large values are indicating numerical problems which cannot be seen at the residuals of the control points, but they may cause large geometric problems for extrapolations outside the three-dimensional area of the control points. A simple significance check of the parameters, e.g. by a Student test, is not sufficient. The 3D-affinity transformation is based on a parallel projection which is approximately available in the orbit direction but not in the direction of the CCD-line. The transformation can be improved by a correction term for the correct geometric relation of the satellite images having only a limited influence (Hanley et al 2002). xij = a1 + a2 ∗X + a3 ∗Y + a4 ∗Z yij = a5 + a6 ∗ X + a7 ∗ Y + a8 ∗ Z Formula 2: 3D-affinity transformation Direct Linear Transformation (DLT): Like the 3D-affinity transformation the DLT is not using any preinformation. The 11 unknowns for the transformation of the object point coordinates to the image coordinates have to be determined with at least 6 control points. The small field of view for high resolution satellite images together with the low object height distribution in relation to the satellite flying height is causing quite more numerical problems like the 3D-affinity transformation. The DLT is based on a perspective image geometry which is available only in the direction of the CCD-line. There is no justification for the use of this method for the orientation of satellite images having more unknowns as required for the solution. xij = L1 ∗ X + L 2 ∗ Y + L3 ∗ Z + L 4 L9 ∗ X + L10 ∗ Y + L11 ∗ Z + 1 yij = L5 ∗ X + L6 ∗ Y + L7 ∗ Z + L8 L9 ∗ X + L10 ∗ Y + L11 ∗ Z + 1 Formula 3: DLT transformation Terrain dependent RPCs: For the relation scene to object coordinates, a limited number of polynomial coefficients shown in formula 1 are calculated based on control points. The number of chosen unknowns is quite depending upon the number and distribution of the control points. Just by the residuals of the control points the effect of this method cannot be controlled. Some commercial programs including this method do not use any statistical checks for high correlations of the unknowns making the correct handling very dangerous. 4. COMPARISON OF METHODS The different mathematical models have been compared especially for IKONOS images in the Zonguldak test field in Turkey. In this area control points have been determined by GPS with a sufficient accuracy. Figure 4: IKONOS, Zonguldak 3D-affinity transformation based on 4 control points discrepancies at independent check points: RMSX=1.91m RMSY=18.53m - at the control points no discrepancies because of missing over-determination Figure 4 shows the result of the IKONOS orientation by means of the 3D-affinity transformation using 4 well distributed control points with quite different height values. Because of missing over-determination there are no discrepancies at the control points. Independent check points were leading to not acceptable results of root mean square differences of RMSX=1.91m and RMSY=18.53m. The problems have been indicated by correlation coefficients listed with r=0.999, resulting in a warning by the Hannover program TRAN3D. Most other programs do not check the numerical problems which have been caused by the fact, that the 4 control points are located nearly on a tilted plane. Also more control points located in this tilted plane would not improve the results. The orientation with a direct linear transformation resulted in similar problems which cannot be controlled just by the location and distribution of the control points. With 6 three-dimensionally well distributed control points the root mean square discrepancies at independent check points have been +/-2.4m - still too much for IKONOS (figure 5). With one more control point, the discrepancies have not been better. The geometric problems are indicated again by high correlation coefficients which have reached r=0.999. Because of this a warning was given by the used Hannover program TRAN3D. Figure 5: IKONOS, Zonguldak direct linear transformation based on 6 control points discrepancies at independent check points: RMSX/Y=2.4m With the terrain dependent RPC-solution similar problems exist like with the two previously mentioned methods. The used commercial software did not indicate any problem for the case shown in figure 6 where 8 control points in a not optimal distribution have been used. The independent check points outside the range of the control points have had discrepancies up to 500m, but also in the area located within the range of the control points extreme errors up to 50m have been present. Figure 6: IKONOS, Zonguldak terrain dependent RCP-solution with 8 control points discrepancies at independent check points With the exception of the terrain dependent RPC-solution all other methods have been tested with a different number of control points (figure 7). Caused by the number of unknowns, the 3D-affinity transformation starts at 4 control points and the DLT at 6 control points. For the sensor oriented RPC-solution and the geometric reconstruction at least one control point has been used to determine the absolute positioning including also remaining datum problems. The geometric reconstruction and the sensor oriented RPCs do show a very homogenous solution - nearly independent upon the number of control points, while the 3Daffinity transformation and the DLT must have at least an over-determination of 2 control points before showing reliable results. Even with a higher number of control points these both methods do show larger root mean square discrepancies at the independent check points. The sensor oriented RPCs are a little bellow the root mean square discrepancies of the geometric reconstruction, but both method are in the sub-pixel accuracy starting at just one control point. As a result it can be mentioned, that the direct linear transformation and the terrain dependent RPCs should not be used. The 3D-affinity transformation requires at least 3 more control points like the geometric reconstruction and the sensor oriented RPCs, in addition the unknowns of the 3D-affinity transformation have to be checked for strong correlations and the control points have to be distributed three-dimensionally. So also the 3D-affinity transformation cannot be recommended. The sensor oriented RPCs and the geometric reconstruction can be used without problems with a small number of not optimal distributed control points. Figure 7: IKONOS Zonguldak Results at independent check points for the different orientation methods as a function of the control point number. Only the case of 32 control points shows the residuals at control points. The geometric reconstruction and the sensor oriented RPCs do transform the scene coordinates to the ground coordinates using the height information for the terrain relief correction. The transformed ground coordinates are based on the accuracy level of the direct sensor orientation. The relation of the transformed coordinates to the control points can be determined by a simple shift of a two-dimensional transformation, for example a two-dimensional affinity transformation. In the Hannover programs RAPORI for the use of the RPCs and CORIKON for the geometric reconstruction, a two-dimensional affinity transformation can be used. The unknowns of the affinity transformation are checked for strong correlation and significance by a Student test. The not justified unknowns can be removed as unknowns. With these methods the required type of transformation has been checked. For the IKONOS-data in the area of Zonguldak there was no justification of an affinity transformation. With just a shift in X and Y even a better accuracy has been reached. Only based on 15 or more control points there was a negligible advantage of the two-dimensional affinity transformation against a simple shift (see figure 8). Figure 8: influence of shift and affinity transformation after terrain relief correction IKONOS, Zonguldak A similar test has been made with a QuickBird image in the same area, partially with the same control points. The QuickBird image did not show the same inner accuracy like IKONOS and a two-dimensional affinity transformation to the control points was required (table 1). After affinity transformation the geometric reconstruction showed some not negligible systematic errors which could be removed with 2 additional parameters. That means, sub-pixel accuracy is possible with QuickBird images, but at least an affinity transformation is required after the terrain relief correction needing at least 4 control points per scene. RPCs geometric reconstruction RMSX RMSY RMSX RMSY shift 1.74 0.72 1.84 0.89 2D-affinity transformation 0.40 0.59 0.81 0.66 0.48 0.46 affinity + 2 additional parameters Table 1: correction of QuickBird ORStandard after terrain relief correction – discrepancies at check points 5. RESULTS ACHIEVED WITH DIFFERENT SATELLITE IMAGES Different optical images have been analysed, including also film images from the Russian TK350, KVR1000, KFA1000 and KFA3000. The results reached with CORONA images are not included because of limited accuracy of the used control points, but a relative accuracy of these panoramic images in the range of 2-3m in X and Y and 3m in Z has been reached. The KVR1000 has a similar geometry like CORONA which can be handled in the Hannover program system BLUH. The film images do not have a pixel size, but the film resolution can be transferred to pixels with the empirical relation 1 line pair = 2 pixels. The analysed images did not have the resolution claimed by Sojuzkarta, so the resolution has been estimated and the result has been transferred to the dimension of pixels. Only the horizontal accuracy is shown in table 2. The vertical accuracy was corresponding to the achieved horizontal accuracy multiplied with the height to base relation (formula 4). SZ = h Spx Formula 4: vertical accuracy of digital images with Spx in [ground sampling distance] b RMSX / RMSY [m] RMSx‘ / RMSy‘ [ground pixel] TK 350, Zonguldak KVR 1000, level 1A, Duisburg KVR 1000, level 1B, Zonguldak KFA 1000, Hannover KFA 3000, Vienna ASTER, Zonguldak KOMPSAT-1, Zonguldak SPOT, level 1A, Hannover. SPOT 5, level 1A / 1B, Zonguldak SPOT HRS, Bavaria 8.3 3.3 10.2 6.5 2.5 10.8 8.5 4.6 5,1 (0,8) (1.6) (5.1) (1.3) (2) 0.7 1.3 0.5 1,0 6.1 0.7 / 1.1 IRS-1C, level 1A, Hannover 5.1 0.9 IRS-1C, level 1B, Zonguldak 9.1 1.6 IKONOS Geo, Zonguldak 0.7 0.7 QuickBird, ORStandard, Zonguldak 0.47 0.76 Table 2: standard deviation of ground coordinates achieved with different space images (space photos) The summary of the results shown in table 2 demonstrates, with well defined and accurate control points usually a sub-pixel accuracy of the orientation of the high and very high resolution space images is possible. In all cases with a lower quality there have been some problems with the control points. A correct mathematical model for the image orientation is required. The available information about the scene orientation should be used to lead the solution to the smallest possible number of unknowns. 6. CONCLUSION The analysis of the different data sets and mathematical solutions showed very clear, a correct mathematical model is required for the handling of the space images and the available sensor orientation should be used. All methods not using the scene orientation information do require more control points and may cause numerical problems if the solution is not checked for high correlation values indicating mathematical dependencies. The extrapolation out of the three-dimensional control point area may lead to extreme large errors for the solutions just based on control points. The direct linear transformation and the terrain dependent RPCs should not be used. Also the 3D-affinity transformation has some clear disadvantages, so the sensor oriented RPCs or the geometric reconstruction should be used for the handling of the level 1Bimages – the projection of the images to a plane with constant height. If level 1A and level 1B-images are given, the same accuracy has been reached with both. With the correct data handling and precise control points in general sub-pixel accuracy is possible. ACKNOWLEGMENTS Thanks are going to the Jülich Research Centre, Germany, and TUBITAK, Turkey, for the financial support of parts of the investigation. REFERENCES Dial, G., Grodecki, J. (2002): IKONOS Accuracy without Ground Control, Pecora 15 / Land Satellite Information IV / ISPRS Com. I, Denver 2002, on CD Grodecki, J. (2001): IKONOS Stereo Feature Extraction – RPC Approach, ASPRS annual conference St. Louis, 2001, on CD Hanley, H.B., Yamakawa. T., Fraser, C.S. (2002): Sensor Orientation for High Resolution Imagery, Pecora 15 / Land Satellite Information IV / ISPRS Com. I, Denver 2002, on CD Jacobsen, K., 1997: Geometric Aspects of High Resolution Satellite Sensors for Mapping, ASPRS Seattle 1997, on CD