PHOTOGRAMMETRIC PROCESSING OF PLANETARY LINEAR PUSHBROOM IMAGES BASED ON APPROXIMATE ORTHOPHOTOS

It is still a great challenging task to efficiently produce planetary mapping products from orbital remote sensing images. There are many disadvantages in photogrammetric processing of planetary stereo images, such as lacking ground control information and informative features. Among which, image matching is the most difficult job in planetary photogrammetry. This paper designs a photogrammetric processing framework for planetary remote sensing images based on approximate orthophotos. Both tie points extraction for bundle adjustment and dense image matching for generating digital terrain model (DTM) are performed on approximate orthophotos. Since most of planetary remote sensing images are acquired by linear scanner cameras, we mainly deal with linear pushbroom images. In order to improve the computational efficiency of orthophotos generation and coordinates transformation, a fast back-projection algorithm of linear pushbroom images is introduced. Moreover, an iteratively refined DTM and orthophotos scheme was adopted in the DTM generation process, which is helpful to reduce search space of image matching and improve matching accuracy of conjugate points. With the advantages of approximate orthophotos, the matching results of planetary remote sensing images can be greatly improved. We tested the proposed approach with Mars Express (MEX) High Resolution Stereo Camera (HRSC) and Lunar Reconnaissance Orbiter (LRO) Narrow Angle Camera (NAC) images. The preliminary experimental results demonstrate the feasibility of the proposed approach.


INTRODUCTION
Planetary photogrammetry is the main method to produce digital orthophoto map (DOM) and digital terrain model (DTM) of celestial bodies.Planetary mapping products are foundation of landing sites selection and mission path planning (Di et al., 2008;Kirk et al., 2008;Naß et al., 2017), which are also essential to planetary scientific research such as geology and geomorphology (Wang et al., 2017;Wardlaw et al., 2018).Currently, there are millions of planetary remote sensing images acquired by orbital imaging instruments.However, only a small portion of these images are fully processed, due to lacking professional people and powerful software tools.The well-known open source software in planetary mapping community, namely Integrated Software for Imagers and Spectrometers (ISIS) developed by United States Geological Survey (USGS), supports many planetary missions (Edmundson et al., 2013).As a matter of fact, ISIS is the de facto standard for preprocessing of planetary images.The Ames Research Centre of National Aeronautics and Space Administration (NASA) also developed a suite of open source 3D surface reconstruction software to facilitate planetary mapping (Shean et al., 2016).On the other hand, commercial digital photogrammetric workstation (DPW) such as SoftCopy Exploitation Tool SET (SOCET SET) provided by BAE Systems is also used to process planetary images with necessary format conversion (Scholten et al., 2005;Kirk et al., 2008;Gwinner et al., 2009).Though commercial DPWs are more robust and efficient than open source software, they are mainly developed for earth observation remote sensing images and only support limited planetary images such as High Resolution Imaging Science Experiment (HiRISE) aboard Mars * Corresponding author Reconnaissance Orbiter (MRO) and Narrow Angle Cameras (NACs) aboard Lunar Reconnaissance Orbiter (LRO).Orbital imaging instruments used for planetary mapping always adopt linear pushbroom cameras, for example NACs on LRO, HiRISE on MRO, Mars Orbiter Camera (MOC) on Mars Global Surveyor (MGS) and High Resolution Stereo Camera (HRSC) on Mars Express (MEX).It is widely acknowledged that the photogrammetric processing of linear pushbroom images is more complicated than conventional frame images.Moreover, the imaging instruments used for planetary mapping differ greatly, which makes the processing of planetary images more complex.Hence, the developed processing method for specific imaging instrument maybe not suitable for other imaging instruments.In a word, photogrammetric processing of planetary remote sensing images is still a very challenging task.Though there are available open source and commercial software tools, it is still meaningful to develop processing techniques of planetary images.
In this paper, we propose a general photogrammetric processing framework for planetary linear pushbroom images.The proposed approach is characterized by image matching based on approximate orthophotos.Here, "approximate orthophotos" means that orbital photographs are rectified with a rough DTM such as Mars Orbiter Laser Altimeter (MOLA) DTM or a preliminary auto-generated DTM rather than final high resolution DTM products.It is noted that the HRSC team also adopted the concept of approximate orthophotos to derive Mars mapping products (Albertz et al., 2005;Scholten et al., 2005).The different points of our method are: (1) a fast back-projection algorithm of linear pushbroom images is introduced to improve the computational efficiency of orthophotos generation and coordinates transformation of conjugate points; (2) an iteratively refined DTM and DOM strategy is developed to accurately determine the approximate values of conjugate points.

Basic Principle
Since planetary remote sensing images have little texture features, image matching is the most difficult task in planetary photogrammetry.Conventional image matching methods including both area-based matching (ABM) and feature-based matching (FBM) usually deliver unsatisfactory results for planetary images.Thus, many researchers proposed optimized image matching algorithms for planetary remote sensing images.Semi-global matching (SGM) algorithm won great success in earth observation fields (Hirschmüller, 2008), which was also used to process HRSC images.(Sidiropoulos et al., 2017) proposed a ring matching method to perform coregistration of multi-instrument planetary images.They are employing this method to coregister Martian images using HRSC as a baseline.(Geng et al., 2017) developed a pixel-level image matching method for HRSC images based on approximate orthophotos.(Wu et al., 2018) proposed an illumination-invariant SIFT matching method for high-resolution planetary images and acquired 40%-60% more matches than the classical SIFT method.(Tao et al., 2018) introduced an automated multi-resolution DTM processing chain based on ASP, called the Co-registration Ames Stereo Pipeline (ASP) Gotcha Optimised (CASP-GO).CASP-GO can be used to produce planet-wide MRO Context Camera (CTX) and HiRISE DTMs.Additionally, Shape from Shading (SfS) and Shape and Albedo from Shading (SAfS) are also employed to derive or refine high-resolution DTM (Albertz et al., 2005;Wu et al., 2017).However, the image matching problem of planetary photogrammetry is not entirely solved.
In order to efficiently derive reliable and accurate conjugate points of planetary images, we apply approximate orthophotos in the whole photogrammetric processing procedures.This means that both sparse tie points extraction for bundle adjustment (BA) and dense image matching for generating DTM are carried out on approximate orthophotos.Compared with original images, approximate orthophotos have many advantages in image matching: (1) geometric distortions caused by perspective projection and relief are removed greatly; (2) original images with different resolution are scaled to identical ground sample distance (GSD); (3) very good approximate values of conjugate points can be provided.These features of approximate orthophotos ease the image matching greatly.Additionally, without man-made structures, the terrain relief of planetary surface show better continuity than earth.This characteristic makes the approximate orthophotos more suitable for image matching of planetary images.

Coordinates Transformation from Approximate Orthophotos to Original Images
Though image matching is performed on approximate orthophotos, both bundle adjustment and surface reconstruction require original image point coordinates.Thus, the coordinates transformation from approximate orthophotos to original images must be solved.Since the reference DTM used in the orthophotos generation process is available, an image point on approximate orthophoto corresponds to a ground point.For linear pushbroom images, one scan line is acquired at an instant of time and each scan line has its own six exterior orientation (EO) parameters.
Consequently, given a ground point we should firstly determine the exact scan line before applying collinearity equation.This requires complicated iterative computation.Moreover, due to millions of conjugate points need such coordinates transformation, computational efficiency need to be taken into account.Here, the coordinates transformation from 3D ground point to 2D image point is referred to as back-projection, which indicates the inverse direction of an image ray.To improve performance, a fast back-projection algorithm for linear pushbroom images based on object space geometric constraints is adopted (Wang et al., 2009).Since orthophotos are always generated with indirect method, the fast back-projection algorithm can also improve the computational efficiency of geometric rectification process.(Geng et al., 2013) presented that the computational efficiency of back-projection for planetary linear pushbroom images can reach one million points per second with a normal hardware environment.This satisfies the photogrammetric processing requirements very well.

Iteratively Refined DTM and DOM
In theory, the pixel coordinates of conjugate points on two overlapping orthophotos should be the same.As a matter of fact, since the geometric accuracy of the initial EO and rough DTM are limited, pixel coordinates of conjugate points on approximate orthophotos always show displacements in the order of several pixels.Note that the initial EO data can be updated through bundle adjustment and the derived DTM can be refined through hierarchical image matching.Thus, new approximate orthophotos can be generated using the updated EO data and refined DTM.This will result in smaller pixel coordinates displacements of conjugate points on approximate orthophotos.Hence, through iteratively refined DTM and DOM, more accurate approximate values of conjugate points are available.Obviously, the more accurate the approximate values of conjugate points, the faster the conjugate points be determined.

Photogrammetric Processing Framework
A general rigorous geometric model (RGM) for planetary linear pushbroom images is introduced into the framework.Thus, different linear pushbroom images such as MEX HRSC and LRO NAC images can be processed with a generic method.This reduces the software development burden and increase reusability.The normalized cross-correlation method is used to perform image matching.Then, least-squares matching (LSM) is applied to achieve sub-pixel accuracy.It is widely accepted that the initial EO parameters show inconsistencies of several pixels.Thus, BA is essential to improve the geometric accuracy.In the BA process, a piecewise polynomial model (PPM) is adopted (Gruen et al., 2003) as trajectory model of linear pushbroom images.The BA program is implemented based on an open source photogrammetric software DGAP (Stallman, 2018).The widely used hierarchical image matching strategy is adopted in the DTM generation process.Moreover, DTMs and approximate orthophotos are refined at each pyramid level to reduce the search space of image matching.The flowchart of the photogrammetric processing framework is presented in Figure 1.The detailed processing procedures are described as follows.

PDS
(1) In the data input process, planetary stereo images of regions of interest are selected firstly.Then, the corresponding Spacecraft Planet Instrument C-matrix Events (SPICE) kernels are determined using ISIS.A rough DTM such as MOLA DTM or LOLA DTM is prepared for later use, for example generating approximate orthophotos and collecting Ground Control Points (GCPs).All photogrammetric calculations refer to Cartesian coordinates system, and the coordinates system of the reference DTM should conform with that used in the photogrammetric processing procedures.
(2) The RGM for planetary linear pushbroom images are established.We adopt the file formats used for airborne linear camera ADS40 to implement the RGM (Tempelmann et al., 2000).More specifically, for each linear array, a camera file (.cam) containing each pixel's focal plane coordinates and an orientation data file (.odf) containing each scan line's EO parameters are derived from SPICE kernels (Acton et al., 2016).Meanwhile, the original Planetary Data System (PDS) format images are converted to common Tag Image File Format (TIFF) images using the famous open source image processing software tool, namely Geospatial Data Abstraction Library (GDAL).
(3) Approximate orthophotos are generated using a rough reference DTM.This reference DTM are stored for subsequent coordinates transformation from approximate orthophotos to original images.
(4) The preparations for BA include GCP collection, tie points generation and weight setting.Ground control information are manually selected at flat area using the reference DTM data.Tie points are matched firstly on approximate orthophotos and then converted to original images.The weight values should properly reflect the accuracy of the observations of different type.Thus, they should be assigned carefully according to empirical value.The observation equations are constructed using extended collinearity equation.Then a Least Square (LS) method is utilized to refine EO parameters.
(5) In the DTM generation process, dense image matching is performed on approximate orthophotos at each pyramid level.Then, the forward intersection is carried out to derive 3D object points.Moreover, both the DTM and the approximate orthophotos are refined iteratively.More specially, the approximate orthophotos at pyramid level i are generated using the DTM derived at pyramid level i+1.At full resolution level, a few manual editing of the auto-generated DTM is usually required.At last, the final mapping products including DTM and DOM are derived.

EXPERIMENTAL RESULTS AND DISCUSSION
The software development was carried out using Visual Studio 2013 and Qt 5.4.2 on the Windows 7 platform.The BA module was implemented based on the open source software DGAP.The hardware configurations of the experimental environment were Intel Core i5 CPU with 8 GB RAM capacity.The test datasets include MEX HRSC and LRO NAC images.LRO NACs have two cameras, namely NAC-Left and NAC-Right, whereas only the NAC-Left images were used to carry out experiments.The radius of the reference sphere for Moon and Mars is 1737.4kmand 3396.19kmrespectively.

Test Datasets
The test datasets are listed in Table 1.For each imaging instrument, one group of planetary stereo images were tested.HRSC has the ability to acquire single-pass stereo images, whereas LRO NACs require consecutive orbits observations to constitute stereopairs.

Planetary Mapping Results
For MEX HRSC images, the generated DTM and DOM are shown in Figure 2. The image resolution of the generated HRSC DOM is 13m, which is equal to the nominal resolution of original images.The generated HRSC DTM has a grid spacing of 50m.
The generated HRSC DTM is in the elevation of -4050~ -1250m.
We compared the derived HRSC DTM with HRSC Level-4 DTM through plotting terrain profiles, and the results are presented in Figure 3.

Discussion
As can be observed in Figure 2, there is a large flat area in HRSC 4165 orbit images.Though it is very difficult to find conjugate points in flat area due to lacking features and low contrast, the auto-generated DTM can express the Martian terrain relief well.This demonstrates that image matching on approximate orthophotos is beneficial to surface reconstruction of planetary images.To further evaluate the geometric accuracy of the generated DTMs, we measured the height values of 500 evenly distributed sample points for the plotted terrain profiles.Then, the height displacements were calculated and the results are presented in Table 2.As can be observed, both HRSC DTM and LRO NAC DTM show small average errors.This demonstrates that there are no obvious systematic displacements between the generated DTMs and the reference DTMs.Considering that the grid spacing of the generated HRSC DTM is similar to that of the HRSC Level-4 DTM product (75m), the small maximum errors (61.2m/1.2pixels)and standard deviations (25.3m/0.5pixels)indicate that the generated HRSC DTM conforms with the HRSC Level-4 DTM very well.Note that LOLA DTM has a resolution of 128 pixels per degree (about 250m), whereas the generated lunar DTM has a grid spacing of 5m.Thus, the large maximum errors (100.7m/20.1pixels)and standard deviation (50.1m/10.0pixels)between the generated LRO NAC DTM and LOLA DTM can be explained by the large differences of DTM resolution.Table 2.The height displacements between the generated DTMs and reference DTMs for the plotted terrain profiles.

CONCLUSION
Approximate orthophotos are extremely suitable for image matching of planetary remote sensing images.The preliminary experimental results demonstrate the feasibility of the proposed photogrammetric processing framework.However, there is still a lot of work to be done.Firstly, the in-house developed software requires further improvement, such as adding more practical features including self-calibration, adaptive image matching parameters, automatic coregistration of reference DTM, etc.Secondly, the photogrammetric processing of multi-instrument planetary images will be enhanced.Additionally, we will improve the parallel computing ability through multi-threaded programming and GPU acceleration.

For
LRO NAC images, the generated DTM and DOM are shown in Figure4.The generated LRO NAC DOM has a resolution of 1.25m.The grid spacing of the generated LRO NAC DTM is 5m.Similarly, we compared the generated lunar DTM with LOLA DTM using terrain profiles.The results are shown in Figure5.As can be observed, LRO NAC images show a thin strip shape.Thus, the terrain profiles were plotted with a north to south direction.Since LRO NAC can acquire very high resolution remote sensing images, very detailed lunar topographic information can be derived.Hence, craters with diameters of several meters can be identified.

Figure 2 .
Figure 2. The DOM and DTM results generated from MEX HRSC images.The red solid line indicates the positions of sampled points for plotting terrain profile.

Figure 3 .
Figure 3.The terrain profiles comparison between the generated Martian DTM and HRSC Level-4 DTM.

Figure 4 .
Figure 4.The DOM and DTM results generated from LRO NAC images.The red solid line indicates the positions of sampled points for plotting terrain profile.

Figure 5 .
Figure 5.The terrain profiles comparison between the generated lunar DTM and LOLA DTM

Table 1 .
Basic information of the test datasets.GSD refers to ground sample distance.