HIGH RESOLUTION SEAMLESS DOM GENERATION OVER CHANG ' E-5 LANDING AREA USING LROC NAC IMAGES

Chang’e-5, China’s first sample return lunar mission, will be launched in 2019, and the planned landing area is near Mons Rümker in Oceanus Procellarum. High-resolution and high-precision mapping of the landing area is of great importance for supporting scientific analysis and safe landing. This paper proposes a systematic method for large area seamless digital orthophoto map (DOM) generation, and presents the mapping result of Chang’e-5 landing area using over 700 LROC NAC images. The developed method mainly consists of two stages of data processing: stage 1 includes subarea block adjustment with rational function model (RFM) and seamless subarea DOM generation; stage 2 includes whole area adjustment through registration of the subarea DOMs with thin plate spline model and seamless DOM mosaicking. The resultant seamless DOM coves a large area (20 longitude  4 latitude) and is tied to the widely used reference DEM SLDEM2015. As a result, the RMS errors of the tie points are all around half pixel in image space, indicating a high internal precision; the RMS errors of the control points are about one grid cell size of SLDEM2015, indicating that the resultant DOM is tied to SLDEM2015 well. * Corresponding author


INTRODUCTION
The Chinese Lunar Exploration Program started in 2004 and includes three phases, named as orbital missions, soft landers/rovers and sample return missions (Ouyang et al., 2010).The first two phases have been successfully achieved by orbital missions of Chang'e-1 (CE-1) and Chang'e-2 (CE-2), lander and rover mission Chang'e-3.The last phase is being implemented by Chang'e-5 (CE-5) mission, which will be launched in 2019 and return at least 2 kilograms of lunar soil and rock samples back to earth.The planned landing area of CE-5 is near Mons Rümker in Oceanus Procellarum, a large area of lunar mare in the northwest region of the Moon (Xinhuanet, 2017;Gbtimes, 2017;Zhao et al., 2017).
High resolution landing area mapping is of great importance for landed missions to support landing trajectory planning and detailed analysis of potential sampling sites.At present, the available lunar global image mosaic products include Clementine global mosaic (100 m/pixel) (Robinson et al, 1999), LROC WAC global mosaic (100 m/pixel) (ASU, 2011), CE-1 CCD camera global mosaic (120 m/pixel) (Li et al, 2010), and CE-2 global mosaic (50 m/pixel released) (Li et al, 2015), and so on.These maps are useful for general analysis of the landing area, but may not be sufficient for detailed analysis and mission planning of CE-5 for their low resolutions.
Lunar Reconnaissance Orbiter Camera (LROC) Narrow Angle Camera (NAC) provides lunar orbiter images with highest resolutions of 0.5 ~ 2 m, and has covered nearly the whole lunar surface.Some featured NAC mosaics (Klem et al, 2014) and RDR products (Tran et al., 2010) have been generated and released by the LROC team, and most of them are of limited coverage area.It is valuable to produce a high-resolution seamless digital orthophoto map (DOM) over the planned CE-5 landing area using LROC NAC images to support detailed morphological and geological analysis, and mission planning.Such a large area seamless mapping is challenging due to the widespread geometric and radiometric inconsistencies among the images.
In this research, we have developed a systematic method and procedure to produce high-resolution seamless DOM of CE-5 landing area using over 700 LROC NAC images.The resultant DOM coves a large area (20 longitude  4 latitude) and is tied to the widely used reference DEM -SLDEM2015.

DATA
The planned CE-5 landing area covers an large area of 49 ~ 69 west longitude and 41~ 45 north latitude, or roughly 443.7km  121.4km.As of December 2017, there are 2299 NAC images within or overlapping with this area.We select 765 images for the mapping task, they are generally of good quality, similar solar azimuth angles, and similar resolutions.The EDR level NAC images are downloaded from the PDS website and are all pre-processed with Integrated System for Imagers and Spectrometers (ISIS) to attach the SPICE kernels (NAIF, 2014) with "spiceinit", radiometric corrections are performed with "lronaccal" command and echo effects are removed with "lronacecho" command.
SLDEM2015 is a lunar DEM generated by co-registration and combining SELENE TC DEM with LRO laser altimetric data.It has a horizontal resolution of 512 pixels per degree (~60 m at the equator) and a 3 to 4 m root-mean-square (RMS) elevation residuals to LOLA profiles (Barker et al., 2016).SLDEM2015 is downloaded from the PDS LOLA data node (http://imbrium.mit.edu/EXTRAS/SLDEM2015/), and is used as the reference for geopositioning and DOM generation.LROC WAC global mosaic is downloaded from LROC web site (http://wms.lroc.asu.edu/lroc/view_rdr/WAC_GLOBAL) and is used as a bridge for matching between NAC images and SLDEM2015.

METHODOLOGY
A two-stage processing method is developed in this research for large area controlled seamless DOM generation and a partitioning strategy is used to ensure both mapping precision and computational efficiency.Figure 1 shows the flowchart of the large area DOM generation process.In stage 1, within each subarea, a planar block adjustment is applied to improve the relative precision among NAC images and improve the absolute accuracy with respect to the control points acquired from SLDEM2015.The rational function model (RFM) of each image is used and refined in the block adjustment.DOM of each image is automatically generated and DOM of the subarea is mosaicked subsequently.After the block adjustment, the inconsistencies between neighboring NAC images within each subarea have been effectively removed.Due to the precision limitation of the control sources, some inconsistencies between the neighboring subareas still exist.Therefore, in stage 2, a thin plate spline (TPS) model based image registration is conducted to the subarea DOMs, based on which the seamless DOM of the whole landing area is produced.

Geometric Modeling of LROC NAC Imagery
A geometric model of the orbital imagery establishes the relationship between image-space coordinates and object-space coordinates and it is the foundation of the block adjustment and image ortho-rectification.
The rigorous sensor model (RSM) of LROC NAC imagery is established based on the collinearity equations with exterior orientation (EO) and interior orientation (IO) parameters, which are retrieved from the relevant SPICE kernels (Tran et al., 2010;Speyerer et al., 2016;Wu et al., 2017;Henriksen et al., 2016).The RFM is a mathematical fitting of the RSM, usually with a RMS residual of 1/100 pixel level.It establishes the relationship between image-space coordinates and object-space coordinates with the ratios of polynomials, and has advantages of platform independence, simple form, and high calculation speed (Di et al., 2003;Liu et al., 2016).
The RFM can be represented as ( , , ) where a 1 , a 2 … to a 20 are the coefficients of the polynomial function P i , named as the rational polynomial coefficients (RPCs).The RPCs are solved by least-squares fitting using large numbers of virtual control points, which are derived by the RSM of the image.

Subarea Block Adjustment
Due to orbit and attitude errors, the RSM, as well as the fitted RFM, may not be sufficiently accurate for mapping applications.
The geopositioning errors are reflected by the inconsistencies between the coordinates of the same surface points measured from different images.To generate high precision seamless DOM, these inconsistencies must be reduced to a sub-pixel level.This is achieved by block adjustment of the LROC NAC images within each subarea, which is the first stage processing of the whole procedure.
In traditional photogrammetric block adjustment, 3-D ground coordinates of the tie points are calculated and refined during the adjustment (Gwinner et al., 2010;Wu et al., 2014;Di et al., 2014).However, when the stereo convergence angle is too small (e.g., < 10), the iteration of 3D block adjustment will be hard to converge and the calculated ground height may be abnormal.Such situation is very common in our large area mapping using LROC NAC images.That is, either the stereo convergence angles are very small in overlapping regions of neighboring images, or there is no stereo coverage at all in most of the regions.To tackle this problem, we propose a planer block adjustment method based on RFM for each subarea, in which ground planar coordinates of the tie points are calculated and refined, while the height coordinates are interpolated from the SLDEM.Meanwhile, instead of recalculating the RPCs, we use an affine transformation model to compensate the systematic errors in image space in the planer block adjustment.
Tie points between every two overlapping NAC images are obtained automatically by image matching of feature points.Control points are acquired interactively between NAC images and the SLDEM with the help of the WAC mosaic.The control points should be evenly distributed in the subarea but it is not necessary to pick control points in every image.
The error equation of the planer block adjustment can be expressed as where, 0 1 2 ,, e e e and 0 1 2 ,, f f f are affine transformation parameters in sample and line directions respectively, sample and line are the back-projected image coordinates of a ground point with original RFM as shown in Equation (1), while x and y are corresponding measured coordinates in image space.The error equation is no-linear with respect to the ground coordinates and thus it is linearized using Taylor series expansion as shown in Equation ( 4).There are two types of unknown parameters: affine transformation parameters and ground planar coordinates (lat, lon) of the tie points.The height coordinates will be interpolated from the reference DEM iteratively.
For control points, it is only necessary to calculate affine transformation parameters and the error equation can be expressed as For each tie point and control point, error equations can be constructed and the unknown values can be solved by the leastsquares algorithm iteratively until the correction values are smaller than empirically given thresholds.
Through the planer block adjustment, the geopositioning inconsistencies among the images within each subarea are effectively reduced to sub-pixel level.Subsequently, based on the refined RFMs, the DOM for each image is generated and a seamless DOM for the subarea is produced through mosaicking .

TPS-based Registration of Subarea DOMs
Inconsistencies still exist among subarea DOMs because of the limited control precision in the subarea block adjustment.To eliminate these geometric inconsistencies, a TPS (Wahba, 1990) based registration among the subarea DOMs are proposed as the second stage processing.
The TPS model is given by in Equation ( 6) (Shen et al, 2017), where, m is the number of control points,  and  0 ,  1 ,  2 and  j in Equation ( 6) are the coefficients that need to be estimated by minimizing the weighted sum of error measure E(f) and roughness measure R(f), as shown below, where  is the smoothing parameter which is restricted to be non-negative.When =0, there is no smoothness constraint and the fitted TPS model passes exactly through all the control points.If    , the coefficients  j will be zero and TPS will reduce to the affine transformation model.In this research,  is set to be 0.
In the TPS-based registration of subarea DOMs, the DOM of the middle subarea is chosen as the reference.The DOMs of adjacent subareas are matched with the reference DOM, and the corresponding points are used as control points for computation of the TPS parameters.Generally, for the adjacent DOM on the left of the reference DOM (Left DOM in Figure 2), the matched control points are distributed in the right-side margin of the DOM; for the adjacent DOM on the of the reference DOM (Right DOM in Figure 2), the matched control points are distributed in the left-side margin of the DOM.In order to avoid over correction of the DOM and restrain the geometric error in local area, auxiliary control points are generated in the opposite-side margin.The TPS parameters for geometric rectification of the DOM are determined by the control points and auxiliary control point together.As a result, coordinates of the control points in the DOM are transformed to that of the reference DOM, while the coordinates of the auxiliary control points in the DOM remain unchanged; the coordinates in between the control points and the auxiliary control points are transformed smoothly using the TPS parameters.Then, the rectified DOMs are then considered as references and the registration process repeats for the rest of the DOMs.In the end, a seamless DOM mosaic of the whole landing area is produced.

Planar Block Adjustment Results
To maintain both the adjustment precision and computational efficiency, the whole CE-5 landing area is divided into 10 subareas in longitude direction.Each subarea covers about 3 with an overlapping of 1.The borders of the subareas are depicted in Figure 3 with red and blue polygons for easy distinguishing.The subareas from the left to the right are named as Part1, Part2, …, and Part10, respectively.The number of LROC NAC images in each subarea is listed in the second column of Table 1.Within each subarea, around 100 control points are extracted to link the NAC images to the reference DEM, and tens of thousands of tie points are extracted to link all the NAC images to form the subarea block.
Figure 3.The borders of the ten subareas The planar block adjustment is applied to each subarea separately and the RMS errors of the control points and tie points are displayed in Table 1.The unit of RMS errors is the NAC image pixel, which is set as 1.5 meters.It can be seen in Table 1 that the RMS errors of the tie points of all the 10 blocks in image space are around half pixel, indicating that the geometric inconsistencies among the NAC images within the subareas have been effectively removed or reduced.The RMS errors of the control points of all blocks are around 30 pixels, which are about one grid cell size of SLDEM2015, indicating that the subarea DOMs are tied to SLDEM2015 well.

TPS-based DOM Registration Results
Restricted by the low resolution of control sources, positional deviations still exist in the DOMs of adjacent subareas.For example, the left image in Figure 5 shows the positional differences between Part1 and Part2 DOMs.The deviations are removed by the TPS-based registration of the adjacent subarea DOMs.The registration is conducted subarea by subarea using the procedure detailed in Section 3.3.The right image in Figure 5 shows that after registration the two DOMs can be seamlessly overlaid and mosaicked.The registration results can be quantitatively evaluated by the coordinate differences of the check points in the overlapping region of two adjacent subarea DOMs.The Part1 and part2 are also taken as an example.Ten automatically matched check points are evenly distributed in the overlapping area.The coordinate differences in latitude, longitude directions and the overall planar coordinate differences are listed in  After registration and rectification of the subarea DOMs, seamless DOM mosaic of the whole area is generated through mosaicking.The final DOM of the landing area is geometrically seamless and radiometrically homogeneous.Figure 6

CONCLUSIONS
This paper presented a systematic method for large area controlled seamless DOM generation using hundreds of LROC NAC images.The developed method mainly consists of two stages of data processing.The first stage is the planar block adjustment with RFM, which can reduce the geometric inconsistencies among the NAC images to sub-pixel level and also tie the images to the reference DEM.The planer block adjustment is applied to subareas covering tens to over one hundred NAC images.The second stage is TPS-based registration of the subarea DOMs, which can reduce the positional differences among the adjacent subarea DOMs and ensure seamless DOM mosaicking in the whole area.
Based on the developed method, a high resolution seamless DOM has been generated in the planned landing area of Chang'e-5 mission using 765 LROC NAC images.As a result, the RMS errors of the tie points are all around half pixel in image space, indicating a high internal precision; the RMS errors of the control points are about one grid cell size of SLDEM2015, indicating that the resultant DOM is tied to SLDEM2015 well.The planar precision after TPS-based DOM registration is generally within one pixel of the output DOM.The final DOM of Chang'e-5 landing area has a spatial resolution of 1.5 meters and the image size is 224721 columns 44945 rows, covering 49 -69 west longitude and 41 -45 north latitude.The DOM is valuable for detailed morphological study of the landing area.The developed method can be applied to other interesting large areas of the moon surface to produce high-resolution seamless DOMs using LROC NAC images.

Figure 1 .
Figure 1.Flowchart of large area seamless DOM generation

F
are the residuals of the image coordinates.

jr
is the Euclidean distance between an arbitrary image point( , )   xy and the jth control point ( , ) jj xy as shown in Equation (7),

Figure 2 .
Figure 2. Illustration of control points (circles) and auxiliary control points (crosses) in TPS-based DOM registration

Figure 4 .
Figure 4. Comparison of the seams between the DOMs of neighboring NAC images before (left) and after (middle and right) planar block adjustment.The right figure shows result of further radiometric correction based on histogram matching.

Figure 4
Figure 4 compares the seams between two neighboring images before and after block adjustment.The left figure shows parts the DOMs generated using NAC images M1158132972L and M1160489335L with the original RPCs, while the middle one shows the DOMs generated using the same images with the

Figure 5 .
Figure 5.Comparison of the geometric deviation between the adjacent subarea DOMs (part1 and part2) before (left) and after (right) TPS-based registration.
is a zooming-out view of the resultant seamless DOM mosaic of Chang'e-5 landing area.The DOM has a spatial resolution of 1.5 meters and the image size is 224721 columns 44945 rows.

Figure 6 .
Figure 6.Zooming-out view of the resultant seamless DOM mosaic of Chang'e-5 landing area.

Table 1 .
No. of images and block adjustment precisions within the ten subareas

Table 2 .
The planar precision after the TPS-based registration is generally within one pixel (1.5m) of the output DOM with the largest error of 2.02 pixels, indicating an satisfactory registration