MAPPING OF PLANETARY SURFACE AGE BASED ON CRATER STATISTICS OBTAINED BY AN AUTOMATIC DETECTION ALGORITHM

The analysis of the impact crater size-frequency distribution (CSFD) is a well-established approach to the determination of the age of planetary surfaces. Classically, estimation of the CSFD is achieved by manual crater counting and size determination in spacecraft images, which, however, becomes very time-consuming for large surface areas and/or high image resolution. With increasing availability of high-resolution (nearly) global image mosaics of planetary surfaces, a variety of automated methods for the detection of craters based on image data and/or topographic data have been developed. In this contribution a template-based crater detection algorithm is used which analyses image data acquired under known illumination conditions. Its results are used to establish the CSFD for the examined area, which is then used to estimate the absolute model age of the surface. The detection threshold of the automatic crater detection algorithm is calibrated based on a region with available manually determined CSFD such that the age inferred from the manual crater counts corresponds to the age inferred from the automatic crater detection results. With this detection threshold, the automatic crater detection algorithm can be applied to a much larger surface region around the calibration area. The proposed age estimation method is demonstrated for a Kaguya Terrain Camera image mosaic of 7.4 m per pixel resolution of the floor region of the lunar crater Tsiolkovsky, which consists of dark and flat mare basalt and has an area of nearly 10,000 km. The region used for calibration, for which manual crater counts are available, has an area of 100 km. In order to obtain a spatially resolved age map, CSFDs and surface ages are computed for overlapping quadratic regions of about 4.4 x 4.4 km2 size offset by a step width of 74 m. Our constructed surface age map of the floor of Tsiolkovsky shows age values of typically 3.2-3.3 Ga, while for small regions lower (down to 2.9 Ga) and higher (up to 3.6 Ga) age values can be observed. It is known that CSFD-derived absolute model ages can exhibit variations although the surface has a constant age. However, for four 10-20 km sized regions in the eastern part of the crater floor our map shows age values differing by several hundred Ma from the typical age of the crater floor, where the same regions are also discernible in Clementine UV/VIS color ratio image data probably due to compositional variations, such that the age differences of these four regions may be real.


INTRODUCTION
Impact crater counting is one of the most important tools for estimating the geologic age of a planetary surface, i.e. the period of time that has passed since the last resurfacing event.
Recent planetary spacecraft missions have provided enormous amounts of high-resolution image data of global coverage for various solar system bodies, allowing for the construction of crater catalogues comprising large numbers of small craters (Stepinski et al., 2012;Robbins, 2009Robbins, , 2016)).
Impacts of small bodies on planetary surfaces lead to the formation of craters (Hörz et al., 1991).Impact craters are especially abundant on Mercury, the Moon and Mars (Grieve et al., 2007), where they persist over long periods of time due to the lack of plate tectonics, an atmosphere, and life.For planetary bodies without an atmosphere, the areal impact crater density increases with increasing surface age, i.e. densely cratered surface parts are usually older than surface parts with a low abundance of craters.The impact crater size-frequency distribution (CSFD), denoting the diameter-dependent areal density of impact craters for a given surface part, allows for an estimation of the so-called absolute model age (AMA) of the surface (Hartmann, 1999;Hartmann and Neukum, 2001;Michael and Neukum, 2010;Michael et al., 2012;Michael, 2013Michael, , 2014Michael, , 2015;;Hiesinger et al., 2011).
In this paper an automatic crater detection algorithm (CDA) is applied to the mare-like floor region of the lunar farside crater Tsiolkovsky, in order to obtain crater counts of similar accuracy as manual crater counts.A surface age map is created based on the automatically obtained crater counts, resulting in a spatially resolved map of the AMA of the examined surface area.

CRATER DETECTION ALGORITHMS
In the domain of planetary science, crater counts are typically performed by manual inspection of orbital images (e.g., Michael and Neukum, 2010).On the other hand, an important research topic in remote sensing is the development of crater detection algorithms (CDAs) which automatically determine the locations and sizes of craters in images (Salamanićcar andLončarić, 2008, 2010).Numerous CDAs have been developed as a consequence of the high importance of knowledge about the impact crater distribution for remote geologic studies of planetary surfaces (Stepinski et al., 2012).However, many works about CDAs are limited to the demonstration of the accuracy of a particular algorithm on a small set of test images displaying quite simple and clearly pronounced craters, while meaningful studies in the field of planetary science require CDAs achieving a high performance also for less obvious, e.g., degraded, impact craters (Stepinski et al., 2012).
Since the start of the last decade, there was a further increase in research works regarding the automatic detection of craters (Salamunićcar and Lončarić, 2010).Existing CDAs can be divided into two different groups.CDAs of the first group detect craters in images, while CDAs of the second group detect craters based on digital elevation models (Salamunićcar and Lončarić, 2010;Di et al., 2014).
According to Stepinski et al. (2012), some image-based CDAs consist of techniques which process the image to detect the characteristic circular or elliptical shapes of crater rims, followed by a matching stage such as the Hough transform (Hart, 2009), while other CDAs employ machine learning approaches for crater detection.These systems are trained using image data sets labelled by human experts, and the detector performance is tested on new image data not used during the training stage (Stepinski et al., 2012, see also Chung et al., 2014).Examples for classifier architectures are boosting techniques allowing for a combination of various feature detectors (Wetzler et al., 2005;Salamunićcar and Lončarić, 2010;Di et al., 2014;Stepinski et al., 2012) and support vector machines (Wetzler et al., 2005;Chung et al., 2014).
In this work, the template matching based CDA introduced by Grumpe and Wöhler (2013) with the ability of automatic detection of small craters (<10 image pixels diameter) in lunar images will be used for crater detection.This template-based algorithm will be applied to orbital image data of the mare-like floor region of lunar crater Tsiolkovsky.The main step in the estimation of the absolute surface age is the statistical analysis of the CSFD constructed using the image-based CDA for determining the locations and diameter values for the derivation of the AMA.Accordingly, a spatially resolved AMA map of the whole floor of Tsiolkovsky will be constructed.

TEMPLATE-BASED AUTOMATIC CRATER DETECTION
The algorithm by Grumpe and Wöhler (2013) used for the analyses in this paper is based on six generic 3D crater models.These six templates were created based on track data of the Lunar Orbiter Laser Altimeter (LOLA) instrument (Smith et al., 2010) of three craters of about 8 km diameter in the vicinity of lunar crater Plato, making use of the known illumination conditions and the known reflectance behavior of the lunar surface, which allow for generating artificial crater template images.
The craters are characterized by three basic shapes: the simple bowl shape, a transient shape with a flat floor, and a complex shape with a central mound.By dividing each of these three profiles into two halves by vertical splitting and rotating each "half-crater" around the vertical axis, a set of six 3D models was generated by Grumpe and Wöhler (2013).Figure 1 shows cross-sections of the resulting crater models.
Artificially illuminating the 3D models using the Hapke model (Hapke, 1984(Hapke, , 2002) ) for the known illumination conditions then resulted in a set of six crater template images.Six of these templates are shown in Figure 2. All available templates were scaled to sizes between 5 and 200 image pixels in order to allow for the detection of craters across a broad range of diameters based on the normalized cross-correlation coefficient between the templates and the image parts being analyzed.
The template matching stage is followed by a fusion process that replaces multiple detections of the same crater at slightly different positions or with slightly different diameters by the average of these positions or diameters.

THE REGION OF INTEREST: TSIOLKOVSKY
The crater Tsiolkovsky (Figure 3) is one of a small number of lava-filled lunar farside craters, thus having a very dark and smooth crater floor.It is located at (20° S, 129° E) and has a diameter of approximately 200 km.The central peak located in the northern part of the mare-filled region has a "W" shape (Tyrie, 1988a).
Previous studies regarding Tsiolkovsky crater suggest ages around 3.8 Ga (Walker and El-Baz, 1982) and 3.5±0.1 Ga (Tyrie, 1988b) for the floor region.The image used in this work was a mosaic of images taken by the Terrain Camera (TC) of the lunar orbiter spacecraft Kaguya (Haruyama et al., 2014) with a resolution of 7.4 m/pixel.The dark mare fill inside Tsiolkovsky crater is our area of interest, as it appears to be geologically homogeneous.
The estimation of the absolute model age (AMA) mainly depends on the number of craters per diameter interval and per area, the CSFD.We chose a surface part of around 100 km 2 size inside the crater Tsiolkovsky as a reference dataset for the crater detection, which has been evaluated manually by Pasckert et al. (2015).The crater diameter interval considered for estimating the AMA ranges from 128 m to 1000 m.Smaller craters were not taken into account because the detection rate of the CDA decreases with decreasing diameter for diameters less than about 100 m, leading to an artificial "roll-off" of the CSFD.
Like virtually all other object detection algorithms, our template matching based CDA relies on a threshold value which controls the sensitivity of the detector.We adjusted this threshold value such that for a given reference area with available manual crater counts, it yields an AMA that comes as close as possible to the manually determined AMA.
To validate this approach, the 100 km 2 reference area has been divided into three parts of equal size, for each of which the absolute difference between CDA-based and manually determined AMA was minimized by a quasi-Newton method.The resulting threshold is then tested on the remaining two areas.The number of craters in each sub-area detected by the CDA and by manual counting, respectively, is summarized in Table 1.The derived CDA-based AMA values between 3.00 Ga and 3.44 Ga are always consistent with the manually determined AMA values (Table 2).
The "optimal" threshold value was chosen to be the arithmetic mean of the three obtained thresholds.Consequently, the CDA was applied to the total 100 km² reference area with a crosscorrelation threshold value of 0.6568, resulting in an AMA of 3.21±0.13Ga, being very similar to the AMA of Ga obtained by Pasckert et al. (2015) based on manual counting.

TSIOLKOVSKY AGE MAP
To create a spatially resolved AMA map, the template matching technique has been applied to the 7.4 m/pixel resolution Kaguya image data with a 600x600 pixel (4.4x4.4 km²) wide sliding window moving over the whole image of the crater Tsiolkovsky with a step width of 10 pixels (74 m).The CDA successfully identifies a high fraction of the manually determined craters and determines reasonable diameter values, as apparent from Figure 4 (a more detailed performance analysis of the CDA is currently in progress).Slight offsets by a few pixels between the CDA-based and manual crater positions tend to occur, which do not affect the AMA estimation as it only depends on the estimated diameters and areal density of the craters.
Based on the optimal threshold value of 0.6568, CSFDs were computed for overlapping quadratic sub-regions covering the whole mare-like crater floor in order to estimate an AMA value for each sub-region, leading to a spatially resolved AMA map (Figure 5).The AMA values were obtained using our Matlab implementation of the CSFD-based method for AMA estimation (Michael and Neukum, 2010;Michael et al., 2012;Michael, 2013Michael, , 2014Michael, , 2015)).Most AMA values of the map lie between 2.9 and 3.6 Ga, with the majority of AMAs centered around 3.3 Ga.In the AMA map, the central peak, the crater walls and the parts of the crater floor not covered by mare basalt were excluded.According to stratigraphic considerations, the age of these steeply sloped regions exceeds that of the dark mare floor.
The color ratio map of the crater Tsiolkovsky (Figure 6) was derived from Clementine UV/VIS multispectral data (red channel: 750 nm/415 nm reflectance; green channel: 750 nm/950 nm reflectance; blue channel: 415 nm/750 nm reflectance, according to Pieters et al. (1994) and http://www.mapaplanet.org/explorer/help/data_set.html#moon_clementine_ratio).It can be used to identify compositionally distinct geological units on the mare-like crater floor.If the mapping of the CDA-derived AMA yields reliable results, similar units should be found within the age map.Although frame-specific calibration artefacts make it hard to unambiguously divide the color ratio image into geological units, some similarities between structures visible in Figures 5 and 6 are clearly apparent through visual comparison.Specifically, four localized AMA anomalies can be found in the eastern part of the crater floor.They are represented by the letters A, B, C and D in Figures 5 and 6, where each letter denotes a small spectrally distinct region in the color ratio image.Regions A, C and D have lower AMAs than the surrounding surface area, respectively (about 2.9-3.2Ga), while for region B a higher AMA of about 3.5 Ga has been estimated.According to Figure 6, these small regions can all be distinguished spectrally from the surrounding surface, indicating compositional variations of the mare basalt.
Figure 2. Basic crater templates.Templates 1 and 2: craters with a small central mound.Templates 3 and 4: craters with flat floor.Templates 5 and 6: simple bowl-shaped craters (Grumpe and Wöhler, 2013).Van der Bogert et al. (2015) showed for a cratered surface area of 100 km² for which the CSFD has been synthetically generated under the assumption of a uniform AMA of 4.0 Ga that estimates of the AMA for sub-regions of 4 km² size deviate by several hundred Ma from the true AMA, while accurate and reliable estimated AMA values are obtained for 10 km² subregions.Given the comparably large size of about 20 km² of the overlapping sub-regions analyzed in this work, the AMA fluctuations apparent in Figure 5 are likely to be real, which is confirmed by their correlation with spectral variations in the Clementine color ratio image.
A more detailed view of the eastern and southeastern part of the crater Tsiolkovsky is shown in Figure 7.It reveals that the two small low-AMA anomalies marked in Figures 5 and 6 by the letter D are located at the same locations as two bright details in Figure 3, corresponding to small mountains of the original crater floor that penetrate through the dark mare basalt layer.Hence, the lower estimated AMA is probably due to small craters being destroyed more easily on the steep slopes on the flanks of these mountains.The low-AMA anomalies marked by the letters A and C cannot be identified with any conspicuous surface details but correspond to dark surface parts with a visually apparent low abundance of small impact craters.In contrast, the high-AMA anomaly marked by the letter B corresponds to an area of slightly increased albedo characterized by a visually apparent exceptionally high abundance of small craters.The distribution of these craters as a cluster, however, suggests that they are secondary craters (McEwen and Bierhaus, 2006).The seemingly increased AMA of region B is thus probably not due to the real surface age but to the presence of these secondary craters (a detailed analysis of the effect of secondary craters on the estimated surface age is provided by McEwen and Bierhaus, 2006).
Hence, the AMA anomalies B and D are not due to lava flows having an age differing from the mean age of the crater floor.The negative AMA anomalies A and C might be interpreted as lava flows of a composition and age different from most of the crater floor area, which would imply the occurrence of small eruptions on the floor of Tsiolkovsky several hundred Ma after the emplacement of the majority of the mare basalts.An alternative interpretation for regions A and C is that the distribution of small craters on the floor of Tsiolkovsky may be affected in a spatially variable manner by secondary craters, which would then cause the observed AMA variations across the crater floor.

SUMMARY AND CONCLUSION
In this work, an image template matching based CDA has been applied to a Kaguya spacecraft image of the mare-like floor of the lunar crater Tsiolkovsky.The AMA of the surface has been estimated based on the CSFD inferred from the detection results of the CDA.The detection threshold was calibrated based on manual crater counts performed for a small sub-region of the crater floor.The average AMA of the crater floor is consistent with previous surface age estimates.
A spatially resolved map of the AMA has been constructed by applying the CDA to overlapping sub-regions of the crater floor.Correlations between structures in the AMA map and local variations of spectral properties in Clementine UV/VIS color ratio data could be established.Two small local negative anomalies in the AMA map could be identified with highlandlike remnants of the original crater floor penetrating the mare fill, while one positive anomaly is probably due to the presence of a cluster of secondary craters.Two other local negative AMA anomalies may correspond to lava flows of different age and composition than most of the crater floor, but they may likewise be due to the spatially varying influence of secondary craters on the CSFD of the small craters on the floor of Tsiolkovsky.
Independent of these interpretations, the proposed surface age estimation and mapping technique allows for a distinction between different geological units in terms of their AMAs.It can be applied in a straightforward manner to different surface regions, as long as manually determined AMA estimates are available for small sub-regions.To validate the results of this method in a more rigorous way, it should be tested on large mare areas on the lunar nearside whose surface ages are accurately known.Regarding the CDA, the choice of a constant detection threshold for all crater sizes as suggested in this work is the simplest possible one.Hence, since we found the sensitivity of the proposed CDA to be higher for larger than for smaller craters, we expect that an improvement of the detection results will be achieved by defining a diameter-dependent detection threshold in future extensions of the CDA.Furthermore, it will be important to extend the CDA such that it is able to distinguish secondary from primary craters based on their shape and/or spatial distribution in order to suppress the influence of secondary craters on the AMA estimation.(Speyerer et al., 2011).
Figure 1.Cross-sections of the crater models.The crater models were constructed from LOLA tracks that closely resemble crosssections of approximately 8 km diameter satellite craters of the lunar crater Plato.Each track was split at the centre of the crosssection, respectively, and the models were derived by rotating the resulting profiles around the centre of the corresponding satellite crater, respectively.(a) Model cross-sections derived from the southern half of the LOLA track.(b) Model cross-sections derived from the northern half of the LOLA track.

Figure 4 .
Figure 4. Part of the calibration area.Manually counted craters are indicated by green circles, CDA detection results by red circles.Image data: Kaguya TC (Haruyama et al., 2014).

Figure 5 .Figure 7 .
Figure 5. AMA map of the floor of the crater Tsiolkovsky.The mare-like floor region has been extracted manually.Black: no data.

Table 1 .
Results of the threshold adaptation: number of craters.

Table 2 .
Results of the threshold adaptation: estimated AMA.