A standardized approach for quantitative characterization of impact crater topography

Icarus 241 (2014) 114–129

Contents lists available at ScienceDirect

Icarus journal homepage: www.elsevier.com/locate/icarus

A standardized approach for quantitative characterization of impact crater topography P. Mahanti a,⇑, M.S. Robinson a, D.C. Humm b, J.D. Stopar a a b

School of Earth and Space Exploration, Arizona State University, Tempe, AZ 85287, USA Space Instrument Calibration Consulting, Annapolis, MD 21401, USA

a r t i c l e

i n f o

Article history: Received 1 December 2013 Revised 20 June 2014 Accepted 22 June 2014 Available online 3 July 2014 Keywords: Cratering Data reduction techniques Impact processes Moon Image processing

a b s t r a c t Historically, topographic profiles provided a quantitative means to investigate the morphology and formation processes for impact craters, although no generic mathematical framework was developed to reduce profiles to morphology descriptors. Only need-specific polynomial expressions were utilized in previous studies, thus no standardized automated comparison of craters exists. We employ a Chebyshev polynomial function approximation to describe crater forms in a quantitative and repeatable manner. We show that the Chebyshev polynomials return coefficients that are relatable to crater morphologic characteristics, thus providing a standardized mathematical means for describing crater forms. Ó 2014 Elsevier Inc. All rights reserved.

1. Introduction The geometry and surface texture of craters record details of cratering mechanisms, the properties of the materials involved and the gradual modification of crater morphology with time. Quantification of crater morphology is an essential approach towards understanding impact phenomena. Meaningful data mining of crater topography and derived features requires a standardized numerical methodology. However, despite extensive previous work classifying crater forms from crater elevation profiles, the topic of generic numerical representation of crater elevation profiles has not received adequate attention such that a standard analytical method has not arisen. Previous crater morphometric studies relied on measures of crater features such as apparent depth, maximum depth, central depth, crater rim-to-rim diameter, mean diameter at the ambient plane, wall slopes and circularity index (Pommerol et al., 2012; Sori and Zuber, 2012, 2010; Talpe et al., 2011; Blaser et al., 2011; Mazarico et al., 2010; Simpson et al., 2008; Bray et al., 2008; Elkins-Tanton et al., 2004; Craddock and Howard, 2000; Garvin et al., 1999, 1998; Cintala and Grieve, 2010; Garvin and Frawley, 1998; Pike, 1977a, 1974a; Roddy, 1976; Ross, 1968; Baldwin, 1963), along with qualitative accounts of other morphological features such as the characteristics of the crater rim, presence of ⇑ Corresponding author. E-mail address: [email protected] (P. Mahanti). http://dx.doi.org/10.1016/j.icarus.2014.06.023 0019-1035/Ó 2014 Elsevier Inc. All rights reserved.

central peaks, floor textures, local terrain slope and deviations from symmetry (Pommerol et al., 2012; Talpe et al., 2011; Bray et al., 2008; Cook et al., 2000; Roddy, 1976; Schultz, 1976; Pike, 1974a; Hartmann, 1972). In many of these studies, crater forms represented by cross-sectional elevation profiles were quantitatively analyzed with measurements of the features mentioned above. The profiles were obtained either directly from 2-dimensional elevation data (elevation vs sample interval) or were sampled for a desired surface track from 3-dimensional digital elevation models (DEMs) (elevation vs sample position in a 2-dimensional reference frame). The majority of these previous studies addressed the group-behavior of crater morphology, often relating crater size descriptors (depth and diameter, for example) with simple geometric relations (power law is one example). Such relationships for crater populations are useful but cumbersome while relating one morphological parameter to another for individual craters. Crater topography shapes evident from the elevation profiles were also used for qualitative classification of age and crater type (Robbins and Hynek, 2012a,b; Sori and Zuber, 2012; Stopar et al., 2010; Talpe et al., 2011; Bray et al., 2008; Cintala and Grieve, 2010; Pike, 1977a; Smith and Hartnell, 1978; Schultz, 1976; Hartmann, 1972; Ross, 1968). In simulation studies, elevation profile shapes aided the investigation of formation and geological degradation mechanisms (Bray et al., 2008; Elkins-Tanton et al., 2004; Craddock and Howard, 2000; Roddy, 1976; Hartmann, 1972). Fitting crater profiles with polynomials is not uncommon and, in the past, arbitrary polynomials (rational and irrational) were

115

P. Mahanti et al. / Icarus 241 (2014) 114–129

used to satisfy need-specific crater profile analysis (for example Craddock and Howard, 2000, Simon and de Bruyn, 2007). In addition, polynomial relationships between crater feature heights (for example depth and rim height) and spatial extent in the sampling axis (for example radial distance from crater center, diameter (D)) were used implicitly in parabolic fits (Talpe et al., 2011) and power laws (Baldwin, 1963; Pike, 1977a) to describe specific relations (for example depth-to-diameter) for crater populations. However, in these polynomial fitting methods there is no consensus on the polynomial (degree or type), the data fitting methods, or the methods for extraction of crater profile features from the fitted polynomial. While a need-specific polynomial selection approach works for specific studies, the lack of standardization across different studies defeats the purpose of having a mathematical representation for the crater profile and makes quantitative comparisons and classification inefficient, if at all possible. Ordinary polynomials and other non-standardized mathematical representations of craters are not particularly desirable for describing craters by translating the shape information from the profile (in observation space) to a set of morphological descriptors (in feature space). This is because often the translated feature space representations are not compact (have a large number of subcomponents), the morphological descriptors obtained are inter-dependent (cannot be obtained selectively and (or) their magnitudes are correlated) and are inefficient as crater profile shape descriptors (physical relationship between the coefficients and shape of the crater profile is complicated). Fast, efficient, algorithm-based determination of morphological parameters from crater profiles is desirable for standardized description of crater forms with the advent of high resolution topographic data. A compact feature vector obtained from a standardized mathematical representation of crater profiles enables a quantitative comparison of crater forms within a global population, and from body-to-body. Obtaining these coefficients should be quick and procedural and the coefficient magnitudes must correlate strongly with the profile shape. Moreover, the physical relevance, particularly which coefficient (or coefficients) affects which part of the shape, should be intuitive. We show that Chebyshev function approximation of crater profiles meets all desired requirements for standardized representation of crater profiles for morphological feature extraction and comparison. Crater profiles obtained from raster topographic data (Tran et al., 2010; Scholten et al., 2012) from the Lunar Reconnaissance Orbiter Camera (LROC) (Robinson et al., 2010) validates the proposed method in this work. However, the concept and methods described here are equally applicable to other planetary topography datasets for the Moon (Smith et al., 1997, 2010), Mars (McEwen et al., 2007; Zuber et al., 1992) and Mercury (Cavanaugh et al., 2007).

2. Methods 2.1. Standardized extraction of crater profiles Local DEMs derived from the LROC Narrow Angle Camera (NAC) (Tran et al., 2010) and global lunar DEM derived from LROC Wide Angle Camera (WAC) (Scholten et al., 2012) images form the basis of this study. Elevation profiles of large craters (D > 2 km) were extracted from the WAC global DEM (sampled at 200 m/pixel) while for smaller craters (D52 km), the NAC DEMs (sampled at 5 m/pixel) were used. In total, 765 craters (of diameters 100 m to 145 km) were analyzed in this work. A standardized method (enabling ease of comparison across studies) was used to obtain the crater profiles (Fig. 1). The elevation profile was selected along the maximum gradient of the

surface outside the crater and through the crater center. The crater profile length was fixed to twice the crater rim-to-rim diameter. Once the profiles were obtained, Chebyshev function approximations of the crater profiles and the Chebyshev approximation coefficients were computed by the use of inner products (Eq. 6, Section 2.2.3). Prior to obtaining the coefficients, the horizontal axes of the crater profiles were rescaled between 1. This standardized extraction of crater profiles is simple and can either be followed manually or implemented in a semi-automated or fullyautomated software tool. 2.2. Standardized representation of crater profiles with Chebyshev polynomials Two distinct approaches are possible for quantitatively characterizing an impact crater profile and associated morphologic features. One approach is to model parts of the crater profile or to model relationships between aspects of the profile (such as crater diameter, depth, rim elevation). Assumption of a model (such as crater scaling laws) is necessary in this approach that has been historically taken and usually applied to a group of similar craters. The other approach (taken in this work) is to use function approximation to describe a crater elevation profile. In this approach, no model is assumed a priori and the desired accuracy in the representation of the crater profile can be tuned by increasing the order of approximation, as described further below. In general, polynomial approximations provide computationally efficient representations of complex functions, such as crater profiles and representing a complex function with fewer number of polynomial terms is desirable (simpler representation and computational advantage). However, during the approximation process, if a function is approximated by a higher order polynomial and then only the first few terms are used (to get fewer terms), a truncation error occurs that is dependent on the properties of the polynomial used for approximation. Chebyshev polynomials produce minimum truncation error and are commonly used for function approximations. We use Chebyshev polynomials of the ‘first kind’ (or type I), commonly denoted as T n ðxÞ (Fig. 2) for this work because the approximation coefficients are simpler to use than those of other Chebyshev polynomials. A more detailed rationale is provided in the following sections. 2.2.1. Chebyshev polynomials Chebyshev polynomials are a series of orthogonal polynomials that are commonly used to retrieve least-square-error function approximations. Orthogonal polynomial series members have a unique shape (not correlated with any other in the series) and satisfy a property that their inner-product (similar to dot product for a vector, for details see 2.3.3) is zero. The unique shapes of the member polynomials can be used as fundamental building blocks for a complex function (Fig. 1a or Fig. 3 left side). Chebyshev polynomials can be defined recursively where the expression of one polynomial leads to the next higher order polynomial. The polynomial for degree ‘n’ is denoted T n ðxÞ (x is the sample or reference axis), and is given by the recurrence relation:

T nþ1 ðxÞ ¼ 2xT n ðxÞ  T n1 ðxÞ;

j xj 6 1

ð1Þ

with T 0 ðxÞ ¼ 1 and T 1 ðxÞ ¼ x being the first two member polynomials (higher order polynomials shown in Table 1). There are two main types of Chebyshev polynomials. Throughout this work we use the term ‘Chebyshev polynomial’ or T n ðxÞ to refer to type I Chebyshev polynomials. Type II Chebyshev polynomials have similar recurrence relations (as the first kind) given by U nþ1 ðxÞ ¼ 2xU n ðxÞ  U n1 ðxÞ with U 0 ðxÞ ¼ 1 and U 1 ðxÞ ¼ 2x. Type III and Type IV Chebyshev polynomials also exist but further discussion of the Chebyshev polynomials other than the first kind and

116

P. Mahanti et al. / Icarus 241 (2014) 114–129

Fig. 1. (A) First five Chebyshev basis functions and a crater profile having a local gradient and central peak. Weighted Chebyshev basis functions (dotted line) contribute to the total crater profile in accordance to their shape resemblance to different portions of the crater profile. (B) The crater profile (gray arrow, length is twice crater rim-to-rim diameter) is obtained along the maximum gradient direction of the pre-impact surface. Physical relationship between Chebyshev coefficients and some crater features: Depth (C 0  C 2 ), average profile elevation (C 0 ), and pre-impact plane height (C 0 þ C 2 ). (C) Tycho crater its elevation profile; The profile is extracted from the DEM (black line) and the profile sampling axis is rescaled for the Chebyshev approximation. The approximated profile (blue) is obtained from the true profile (black) using only 8 Chebyshev coefficients (shown in gray). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Fig. 3. Rationale for selecting type I Chebyshev polynomials for crater elevation profile approximation.

Fig. 2. The shapes of Chebyshev polynomials and their use to construct any arbitrary crater elevation profile shape. Shapes for only the first 9 Chebyshev polynomials are shown (and the 33rd). For each of the polynomials, the maximum and minimum values are +1 and 1.

other orthogonal polynomials is beyond the scope of this work and can be found elsewhere (see Gautschi, 2004, p. 28).

Chebyshev polynomials are a series of (function) shapes and by expanding a series of Chebyshev polynomials, a crater profile can be approximately reproduced (synthesized). While synthesizing an arbitrary complex function, such as a crater profile, each of the component polynomials is scaled by a multiplier (coefficient) and summed (Fig. 2). The magnitude of the multiplier determines the contribution of an individual polynomial shape. For example, if a

P. Mahanti et al. / Icarus 241 (2014) 114–129 Table 1 First 9 Chebyshev polynomials. N

T N ðxÞ

T N ð1Þ

T N ð0Þ

T N ð1Þ

0 1 2 3 4 5 6 7 8

1 x 2x2  1 4x3  3x 8x4  8x2 þ 1 16x5  20x3 þ 5x 32x6  48x4 þ 18x2  1 64x7  112x5 þ 56x3  7x 128x8  256x6 þ 160x4  32x2 þ 1

1 1 1 1 1 1 1 1 1

1 0 1 0 1 0 1 0 1

1 1 1 1 1 1 1 1 1

C n T n ðxÞ

ð2Þ

n¼0

where C represents coefficients and TðxÞ represents basis functions. Basis functions (in this case) are lower degree polynomial functions of the sample position. A regular (also called as basic or ordinary) polynomial expressed in the above form is a similar linear sum

pðxÞ ¼

M X n¼0

an xn ¼ a0 þ a1 x þ a2 x2 þ    þ aM xM

M X C n T n ðxÞ

ð4Þ

n¼0

2.2.2. Rationale for using Chebyshev polynomials for this work In the introduction (Section 1), the historical use of polynomials to interpret crater elevation profiles was discussed and earlier in this section we have discussed our intent of using function approximation for characterizing impact crater elevation profiles (summarized in Fig. 3). It turns out that, apart from the computational ease of using polynomials, their use to approximate an arbitrary function defined on an interval is guaranteed to be efficient. The Weierstass theorem (Rivlin, 2003; Hewitt, 1947) provides confidence that an efficient approximation (convergence of the polynomial representation to the true function) is mathematically guaranteed in the context of the current work (where the crater elevation profile is the true function and is being represented by polynomials). Note that the scaling of the sample axis between 1 that follows the crater elevation profile selection confines the continuous elevation function over a fixed interval of sample positions. A typical polynomial pðxÞ of one variable (with positive integer exponents only) and M þ 1 terms can be expressed as a linear sum M X

where an are the coefficients or scalar multipliers and the basis functions are powers of x (1; x; x2 ; . . .). Also, from our earlier discussion, an arbitrary elevation profile function f ðxÞ defined over an interval can be approximated by a polynomial. Combining these two aspects we get the expression for the Chebyshev approximation of an arbitrary crater elevation profile as:

f ðxÞ ffi PðxÞ ¼

component shape is not present, the corresponding multiplier is zero. Such a collection of fundamental component functions (shapes) that can be used to synthesize a complex function is called a ‘basis set’. The Chebyshev polynomials constitute a basis set, and it can be mathematically shown that any arbitrary continuous function within the jxj 6 1 interval can be synthesized using this basis set (Gautschi, 2004; Nguyen et al., 2011). The analysis of any arbitrary function shape (Fig. 2) can be thought of as a reverse process to synthesis, where the function is broken down into fundamental, unique component shapes. If the component shapes are pre-decided, that is, if the basis set is known, then only the multipliers (coefficients) need to be computed. This procedure is called function approximation, and when the Chebyshev basis set is used, this is called the Chebyshev approximation process. As more and more component shapes (and coefficients) are selected to approximate an arbitrary function, the magnitude of the difference (residual) between the target function shape (a crater profile) and its approximation (linear combination of component shapes) decreases until the desired order of approximation (or residual magnitude) is reached. One of the main goals of this work is to show that the coefficients of the basis set are strongly correlated to morphological descriptors and can be used in classifying craters based on their elevation profiles.

pðxÞ ¼

117

ð3Þ

It may be noted that an approximation sign was used in the above expression and not an equality. This is because the sum on the right-hand-side has finite number of terms (terms after M þ 1 were truncated) and this is essential for computational purposes. However, the truncation implies an approximation error whose magnitude is equal to the combined higher order terms not included. When an arbitrary elevation profile is expressed using a regular polynomial (as in some historical approaches) then the basis functions are powers of x (T n ðxÞ ¼ xn ). The problem with the basis functions of regular polynomials is that they are highly correlated (for example, if x increases, x2 increases) unless the argument (x) has a wide variation. So, even if one additional polynomial term is added, the coefficients of all the polynomial terms can change substantially even if the additional term adds an insignificant contribution. In order to have useful standardized representation of the crater profile by polynomials, this is not desirable as the definition of the coefficients would keep changing with the order of polynomial. Minimal correlation between the basis functions and minimal dependence of the coefficient values with the polynomial order is required so that irrespective of the polynomial order used, users can compare corresponding numbered coefficients. The above condition is possible with polynomials having an orthogonal basis set such that there is zero correlation between different functions from the basis set (inner product is zero). Additionally, the choice of orthogonal polynomials ensures minimum mean squared error during the approximation. There are many systems of orthogonal polynomial functions and of those Chebyshev polynomials are particularly well-suited for the proposed work for approximating other functions (Beckmann, 1973; Gautschi, 2004). This is because the smallest maximum deviation from the true function (in this case the crater elevation profile) at every point on the reference axis is obtained only for the Chebyshev polynomials (Mason and Handscomb, 2010; Powell, 1967). Note that mean squared error is an error measure over the entire profile shape and says nothing about local conformity to the target function. Additionally, an upper bound on the maximum approximation error can be immediately computed when Chebyshev polynomials are used in an approximation. Further, any of the Chebyshev polynomial forms (first or second kind) can be chosen for this work (since one form can be represented in terms of the other). However the Chebyshev polynomials of the first kind is chosen as all the polynomials have a maximum value of +1 and a minimum value of 1 which simplifies computations and comparisons (scaling factor on the elevation axis is fixed) across component shapes. 2.2.3. Relevant properties of Chebyshev polynomials and Chebyshev approximation In this section some of the properties of Chebyshev polynomial which make them special for the current work are discussed. Note that unless otherwise stated: (a) n refers to the degree of a particular Chebyshev basis, so T n ðxÞ ¼ 2x2  1 implies n ¼ 2 because the degree of the polynomial basis is 2. (b) M is the order of approximation and the upper limit of the approximation sum. It implies M þ 1 terms in the approximation and M ¼ 0 corresponds to the first basis function. An order M approximation also means that

118

P. Mahanti et al. / Icarus 241 (2014) 114–129

Table 2 Tycho crater elevation profile approximation using ordinary and Chebyshev polynomials. N

Chebyshev polynomial approximation Coefficients

0 1 2 3 4 5 6 7 8

Regular polynomial approximation Total RMSE

Coefficients

M=4

M=8

M=4

M=8

M=4

M=8

M=4

4.0 0.1 11.8 0.6 6.1 0.2 3.9 0.8 0.3

4.0 0.1 11.8 0.6 6.1

14.9 14.9 9.1 9.1 6.8 6.7 5.5 5.4 5.4

14.9 14.9 9.1 9.1 6.8

34.3 6.9 185.0 63.8 250.1 126.7 41.8 69.6 67.7

29.3 2.3 102.0 3.2 68.6

33.8 34.1 83.1 88.8 59.5 59.8 44.4 30.7 5.4

29.4 29.4 36.7 36.7 6.7

the maximum numbered Chebyshev coefficient available from the approximation is M and corresponds to the ðM þ 1Þth basis function which is a polynomial of degree M. Approximation of the Tycho crater profile (Fig. 1C) using Chebyshev and regular (ordinary) polynomials (Table 2) is used throughout this section as example and also to illustrate the properties of Chebyshev approximation. Two different metrics for evaluating the approximation accuracy are used. The first one evaluates the local performance and is called spatial relative reconstruction error (SRRE) and the other indicates performance for the entire profile and is a normalized (weighted average) version of SRRE and also will be referred to as average approximation error in this work. Details for SRRE and normalized SRRE are discussed in the analysis of approximation error (Section 3.4). 2.2.3.1. Symmetry and extrema properties of Chebyshev polynomials. Chebyshev polynomial series members are alternately even (symmetric) and odd (anti-symmetric) polynomial functions. If the value of n is even, then the corresponding Chebyshev polynomial T n is also even and a symmetric function (symmetric about the vertical axis). Similarly, T n is an odd (anti-symmetric) function if n is odd. Only the odd Chebyshev polynomials pass through the origin (they have x ¼ 0 as a root). Type I Chebyshev polynomials have the unique property that for the 1 6 x 6 1 interval, the critical (extreme) values of the polynomials are either þ1 or 1. Moreover, these extrema occur at critical reference axis positions (x ¼ 0, x ¼ þ1 and x ¼ 1), enabling easy association of Chebyshev function magnitudes with important reference axis positions (Table 1). The functional values of the Chebyshev polynomials at critical reference axis positions are important in context to the crater elevation profile analysis. The standardized crater elevation profile extraction method assigns the crater rim center to x ¼ 0; x 2 ½1; þ1 on the sample reference axis. The crater elevation profile extraction also ensures that x ¼ 1 are spatial positions at a distance equal to the rimto-rim diameter (D) from the crater rim center. Since each Chebyshev polynomial is unique in shape (basis functions), and every crater elevation profile is approximated over the same reference axis, topographic behavior of the crater elevation profile gets ‘spatially locked’ to the Chebyshev basis function magnitudes. 2.2.3.2. Orthogonality property of the Chebyshev polynomials. The Chebyshev basis is orthogonal, meaning the inner product of two Chebyshev polynomials is zero if they are not identical (they are uncorrelated). The following inner product expression summarizes the orthogonality property:

8 > < 0; m – n T m ðxj ÞT n ðxj Þ ¼ P þ 1; m ¼ n ¼ 0 > : j¼0 0:5ðP þ 1Þ; m ¼ n – 0

P X

Total RMSE

M=8

ð5Þ

Note the summation that gives the inner product above simply computes the point-to-point sum-of-product for two Chebyshev polynomials, and the result is a constant if the polynomials are same and zero otherwise. The discrete point set xj where the Chebyshev polynomials are evaluated are zeros (roots) of the Chebyshev polynomial T Nþ1 ðxÞ. Details on the selection of this discrete point set and orthogonality can be found in Mason and Handscomb (2010). It can be shown that the orthogonality property helps in evaluating the individual Chebyshev approximation coefficients as per the following expression. The summation is again computed over the discrete point set xj defined earlier.

Ck ¼

n 2 X f ðxj ÞT k ðxj Þ n þ 1 j¼0

ð6Þ

In the above expression, f ðxÞ is the function to be approximated (the crater elevation profile) and C k and T k ðxÞ are the kth Chebyshev coefficient and kth Chebyshev polynomial, respectively. Note that since the Chebyshev basis expression is already known, obtaining the Chebyshev coefficients is straight-forward and procedural and can be efficiently implemented via inner products. 2.2.3.3. Compactness of Chebyshev approximation. The compactness property of Chebyshev coefficients means that only a small set of coefficients is sufficient to represent a crater profile. The Chebyshev approximation process is such that the first coefficient is the largest possible and the contribution of the basis function progressively decreases as higher order basis functions are included in the approximation process. It may be noted however that while for an absolutely smooth continuous (ideal case) crater profile, absolute values of Chebyshev coefficients would decrease monotonically, for more realistic scenarios, this can be more of a general or group-wise trend, rather than a strict monotonic decrease. So occasionally there may be a higher numbered coefficient(s) with a value larger than the lower numbered one preceding it, but the general trend will be a decrease in absolute coefficient magnitude as coefficient number increases. The magnitude or energy of coefficients indicate the topographic variation represented by a basis function. Since most of the energy (and hence topographic variation) is concentrated in the lower numbered coefficients, the first few basis functions represent most of the topographic variation and higher numbered basis functions represent progressively subtler topography. A mathematical proof for the compactness property can be found in Gil et al. (2007). A compact representation also implies an efficient digital storage scheme for the crater elevation profile sample data. The efficiency is quantified by compression ratio which is defined as the number of Chebyshev coefficients required to reconstruct a crater elevation profile with desired level of accuracy to the actual

P. Mahanti et al. / Icarus 241 (2014) 114–129

number of elevation profile samples. Some examples of compression ratios are discussed in Section 3.4. 2.2.3.4. Order of approximation and independent nature of Chebyshev coefficients. The order of Chebyshev approximation is not fixed in the proposed approach but depends on the desired accuracy of approximation and also on retrieving desired Chebyshev coefficients. The accuracy of approximation increases as more higher spatial frequency topography is represented by including higher numbered Chebyshev coefficients. Alternately, if a particular global comparison across craters requires specific set of coefficients, for example, the sum C 4 þ C 8 þ C 12 þ C 16 , then the order of approximation is dependent on the highest numbered coefficient (M ¼ 16 in the above example). It is emphasized here that irrespective of the inclusion or noninclusion of higher order coefficients in the Chebyshev approximation, the lower-order coefficients remain the same. For example, if there are two different studies done on different crater profiles, one with M ¼ 10 and another with M ¼ 20, the first 11 coefficients can be compared across these two studies. This is also shown by the Chebyshev coefficient values (Table 2) obtained by an order 8 (M ¼ 8) and an order 4 (M ¼ 4) approximation of the Tycho crater profile (Fig. 1C). The corresponding coefficient values are identical at each N, irrespective of the order of approximation. In comparison, for the basic (ordinary) polynomial, the coefficient values completely change with order of approximation. Thus, order of Chebyshev approximation does not change the value (and the physical meaning) of the Chebyshev coefficients. Further, we note that only desired Chebyshev coefficients can be obtained without the need of obtaining all the coefficients due to the orthogonality property discussed earlier. 2.2.3.5. Truncation error and smallest maximum deviation for Chebyshev approximation. In the earlier sections, it was stated that the value of M is selected based on the desired reconstruction accuracy. Chebyshev approximation procedure allows us also to measure the truncation error caused by not considering the terms indexed M þ 2 and beyond. A simple rule of thumb is that the residual (truncation error) is the same order-of-magnitude as the M þ 1th coefficient, assuming the magnitude of the coefficients decrease geometrically as M increases. The magnitude of residuals at a pre-fixed value of M indicates the magnitude of higher frequency topographic variations that are not captured by the first M þ 1 coefficients. Note that truncation error, approximation error and residual magnitude are equivalent terms for describing the mismatch between the true crater elevation profile and its Chebyshev approximation. To illustrate the uniqueness of Chebyshev polynomials with respect to truncation error, the profile for Tycho crater (Fig. 1C) is approximated by a regular polynomial and also by a Chebyshev polynomial up to 9 coefficients (M ¼ 8). The highest degree for each approximation will be 8 (from the term x8 ) and if all the terms are chosen (N = 8, Table 2), then both Chebyshev and regular polynomial approximations have same root mean square error (RMSE) of 5:4. Similarly if the approximation is up to 5 coefficients (M ¼ 4), the highest degree for each will be 4 and if all terms are chosen for approximation, both methods will have same RMSE. However, if the M ¼ 8 approximation is truncated to first 5 terms (N ¼ 4), then the RMSE is higher for the ordinary polynomial approximation and this is true at each N (except N ¼ M when the approximation errors are equal). The approximation of the Tycho crater elevation profile shows decreasing truncation error as the number of terms selected for approximation increases (Table 2). In the typical case of Chebyshev approximation the coefficients decrease rapidly and the resulting residual error function is an oscillatory function distributed evenly

119

over the interval (1). The smooth spreading of error is an important unique property of the Chebyshev polynomial causing it to approach the best possible approximating polynomial for the crater elevation profile) of degree M having the smallest maximum deviation (from the true elevation profile).

3. Results from using Chebyshev coefficients 3.1. Physical meaning and errors associated with the Chebyshev coefficients The coefficients of Chebyshev approximation are multiplication constants for the Chebyshev basis functions. As such the Chebyshev coefficients can be understood as ‘tuning factors’ to the shape of the basis functions which are unique building blocks of the crater elevation profile. Note that the Chebyshev basis functions are always fixed and defined over a domain and a range of 1. So the maximum difference between the highest and lowest value of a Chebyshev basis function is 2. But when multiplied by a coefficient C n , this difference becomes 2C n . This is true for all the coefficients that are used in approximation and the combination of the tuned shapes (basis functions) builds the crater elevation profile. The relative importance or the physical meaning associated with a Chebyshev coefficient is governed by the shape of the corresponding basis function. Hence, for a specific location on the sample reference axis it is important to know the sign and relative magnitude of each of the basis functions being used. For example the first Chebyshev basis function is a constant = 1. So, its effect is uniform and the first Chebyshev coefficient adds an offset (C 0 ) at each spatial position. The second Chebyshev basis function vanishes at x ¼ 0, so the second Chebyshev coefficient C 1 has no control on the shape at x ¼ 0 but linearly offsets the shape at all other positions. The second Chebyshev basis is maximum at x ¼ 1 and is minimum (but does not vanish) at x ¼ 0, accordingly the coefficient C 2 strongly correlates to the overall shape at these critical points. Moreover, from the previous discussion, C 0 and C 2 are more relevant than C 1 at x ¼ 0 (at the center of the crater profile). The relative spatial importance of a Chebyshev coefficient can be qualitatively understood by visual inspection of the function shapes and quantitatively by substituting the value of x in the basis function expression. Analysis of Chebyshev basis function shapes leads to the selection of coefficients or combinations of coefficients which correlate strongly to specific crater morphological indicators (Sections 3.2 and 3.3). The errors that occur in practice for the proposed standardized method is not as much associated with the actual value of the coefficient themselves but in their use and interpretation. This is because computations of Chebyshev coefficients is robustly performed based on the orthogonal properties of the Chebyshev basis. Accordingly, negligible procedural error can be associated with the value of the coefficient themselves and the error magnitude is only marginally dependent on the computation scheme and does not vary for different crater profiles. However, the values of the coefficients are dependent on the spatial sampling resolution of the elevation profile. The differences in computed coefficient values due to differences in sampling resolution are more evident in the higher order coefficients (which represent the finer morphological details). For specific morphological descriptors, order of reconstruction dictates the error magnitude. For example, the depth of the crater is proportional to the difference C 0  C 2 . But if the order of reconstruction is increased, other even numbered coefficients can be utilized to get an even more accurate estimator for depth. As such, more accurate topographic indicators are obtained from

120

P. Mahanti et al. / Icarus 241 (2014) 114–129

combinations of coefficients (Section 3.3) and single Chebyshev coefficients are simpler approximations (Section 3.2). It must be noted also that the analysis of crater profiles is 2-dimensional, so the coefficients describing a crater profile should not necessarily be interpreted to represent the 3-dimensional crater topography. 3.2. Lower numbered Chebyshev coefficients as morphological indicators Meaningful crater morphologic characterization is obtained from the Chebyshev approximation of crater profiles. Preliminary results are presented (Figs. 4–6 and Table 3) showing the correlation between the Chebyshev coefficients and true morphologic measurements. Fig. 4C shows the relationship between the 1st Chebyshev coefficient (C 0 ) and mean elevation for profiles obtained from 702 small (D < 200 m) craters. From previous discussions C 0 is the scaled average of the profile and hence a strong linear relationship (R-squared value  1) is expected. The effect of varying the magnitude of the 2nd Chebyshev coefficient (C 1 ) is illustrated in Fig. 5A, where a topographic gradient (plus offset) is added to a crater profile. Note that the magnitude of C 1 changes strongly (Fig. 5A bar chart: groups of first 6 coefficients shown) from the first profile (no slope added) to the third profile (maximum slope added) due to the addition of the local slope. The first coefficient (C 0 ) also changes slightly (inset, Fig. 5A bar chart) due to the offset. A slope corrected profile is easily generated by reconstructing the crater profile as per Eq. 4 but setting the 2nd coefficient (C 1 ) intentionally to zero. This procedure was applied to the profile of Giordano Bruno crater (Fig. 5B) where 33 coefficients were used for the reconstruction (M ¼ 32; C1 ¼ 0). C 2 , the 3rd Chebyshev coefficient is the magnitude for the second Chebyshev basis function (2x2  1). The domain of the

Fig. 5. Identification and correction of local slope. (A) Variation in local elevation and gradient proportionally alters the first and second Chebyshev coefficient (C 1 ) respectively (bar chart). (B) Profile of Girodano Bruno modified (slope removed) by setting C 1 ¼ 0 in the Chebyshev approximation.

Fig. 4. Examples of relationship of Chebyshev coefficients and combinations of Chebyshev coefficients to crater features. (A) Sum of every 4th Chebyshev coefficient in the series of Chebyshev coefficients ((C 4 þ C 8 þ C 12 þ C 16 ) used here) can be used to classify craters based on the heights of their central peaks relative to their depths. (Icarus is atypical) (B) Correlation between Chebyshev coefficient C 2 and crater depth (C) Correlation between average elevation and Chebyshev coefficient C 0 . For (B) and (C) the data is obtained from 702 small (D < 200 m) lunar craters in the highland ponds region (41°N, 167°E).

P. Mahanti et al. / Icarus 241 (2014) 114–129

121

Fig. 6. True crater profile and Chebyshev approximation (M = 16) for Linne and Larmor Q craters (A and B). The local relative error (percent error) indicates the spatial error magnitude variation. Magnitude of Chebyshev coefficients (contribution to profile shape) and associated power (normalized and expressed as percentage of total power) (C and D) decreases (dashed arrow direction) for higher numbered coefficients.

Table 3 Summary of Chebyshev coefficient based crater morphology indicators discussed in this work. Indicator

Description of indicated morphology

C0 C1 C2 C3 C4 I1 I2 I3 I4 I5

Average crater profile elevation Local topographic gradient Crater depth (first approximation) Crater profile asymmetry Central peak (first approximation) Elevation at crater center Central peak Crater profile elevation at left-hand-side extrema Crater profile elevation at right-hand-side extrema Crater depth (more accurate)

Chebyshev polynomial is between + 1 and 1 and so the minimum and maximum values of this expression (1 and 1) occur at 0 and 1 respectively. If the crater profile is approximated by the first 3 basis functions that is: f ðxÞ ffi C 0 þ C 1 x þ C 2 ð2x2  1Þ, then the difference between the minimum and maximum values (2C 2 ) is expected to be positively correlated with the depth of the crater. The correlation is stronger when the crater profile shape matches a parabolic profile, common to small fresh craters (Fig. 4B). The robust fit to the data is linear with a goodness-of-fit (R-squared) value of 0.9. Five craters with varying relative size of central peaks, from no central peak (Linne) to central peak higher than crater rim (Icarus) were compared using the sum of Chebyshev coefficients for the Chebyshev basis functions that have greater magnitude at the crater center (x ¼ 0; where the peak would occur, if present). One such basis function is T 4 and every 4th basis function following T 4 behaves similarly (see shapes of T 4 and T 8 in Fig. 2). Restricting

the Chebyshev approximation to 17 coefficients (M = 16), the sum C 4 þ C 8 þ C 12 þ C 16 is used to compare the height of crater central peaks relative to their depths (Fig. 4A). The bar chart results differentiate the craters distinctly, with a progression from the crater Linne (no peak) to crater Icarus (central peak summit above crater rim). 3.3. Combinations of Chebyshev coefficients as morphological indicators In addition to the shape relevance of simple Chebyshev basis functions to the crater topographic profile (Fig. 1A), the combinations of basis functions are also useful for characterizing crater elevation profiles (summarized in Table 3). For example the odd (N) basis functions have zero value at x ¼ 0. Hence, only even coefficients are relevant when a feature at the crater center (x ¼ 0) is investigated. The following identity gives the value of the crater elevation profile at x ¼ 0

I1 ¼

M 1 X k¼0

ðjÞ

k

½1 þ ð1Þk  C k ; k 2 Iþ 2

ð7Þ

 pffiffiffiffiffiffiffi where j is the unit imaginary number j ¼ 1 . If only 17 coefficients ðM ¼ 16Þ are used to approximate the crater elevation profile under analysis, the profile elevation value at the center of the crater can be approximated by C 0  C 2 þ C 4  C 6 þ C 8  C 10 þ C 12  C 14 þ C 16 . If the crater exhibits a central peak morphology, then this sum is directly proportional to the central peak height. Magnitudes for the first nine basis functions at the critical points (x ¼ 1; 0; 1) are obtained by evaluating each basis function at the critical points on the sample reference axis (Table 1). Note that rows corresponding to odd (N) Chebyshev basis functions are identical: their contribution is de-emphasized (value is

122

P. Mahanti et al. / Icarus 241 (2014) 114–129

negative) at x ¼ 1 and emphasized at x ¼ 1 (value is positive) thus making their spatial contribution from left-to-right asymmetric. Accordingly, a higher contribution from the odd (N) Chebyshev basis functions is expected for craters that have an asymmetric elevation profile. Row values are also seen to repeat for the N series: N ¼ ½2; 6; 10; . . . and also for N ¼ ½0; 4; 8; . . .. For the first set the contribution at the center is de-emphasized (value is negative) while for the second set the contribution at the center is emphasized (value is positive). Enhanced topographic activity specifically at the crater center is thus indicated by the Chebyshev coefficients for N ¼ ½0; 4; 8; . . . and a quantitative metric can be formulated from these coefficients. One simple example metric is the sum of every 4th Chebyshev coefficient starting at C 4 M

I2 ¼

k¼ 4 X C 4k ;

k 2 Iþ

ð8Þ

k¼1

which obtains the sum of all positive relief contributions at the crater center. If only 17 coefficients (M ¼ 16) are used the elevation value at the center of the crater can be approximated by C 0 þ C 4 þ C 8 þ C 12 þ C 16 ; a measure useful to describe central peaks, as was already shown in the earlier section. As discussed earlier, the crater elevation profile extremities occur at x ¼ 1, on the reference plane outside the crater rim. The Chebyshev basis functions attain extrema at these reference interval end-points (T N ð1Þ ¼ ð1ÞN ; T N ð1Þ ¼ 1) and the corresponding elevation values are easily obtained from a Chebyshev coefficient combination. Specifically, the profile elevation at lefthand-side and right-hand-side sample reference axis extremities (x ¼ 1 and x ¼ 1 respectively) are given by I3 and I4 (as below), respectively.

I3 ¼

M 1 X

T N ð1Þ ¼

n¼0

I4 ¼

M1 X n¼0

M 1 X

ð1Þk C k

ð9Þ

k¼0

T N ð1Þ ¼

M1 X

Ck

ð10Þ

k¼0

A crater elevation profile can be considered as a negative relief on a pre-impact reference (line) and identities I3 and I4 pertain to elevations at the endpoints of this reference line. The reference line can be horizontal or tilted. In the absence of strong random topographic variations and nearby craters near the extreme points (x ¼ 1), the quantity 12 ðI3 þ I4 Þ approximates the average elevation of the reference baseline (different from the crater elevation profile average indicated by C 0 ). Subtracting I1 (Eq. 7) from this average indicates the elevation difference between the average reference line elevation and the elevation at x ¼ 0. The resulting expression can be simplified to

I5 ¼

  X pffiffiffiffiffiffiffi 1 M1 k 1 þ ð1Þk 1  ðjÞ C k ; j ¼ 1 2 k¼0

Chebyshev coefficients (Section 3.3) up to the Nth coefficient, the value of N is depends on the order (M) of the Chebyshev approximation.

3.4. Analysis of approximation error An accurate approximation of the crater profile enables better topographic indicators and a higher order approximation gets closer to the true crater profile. However, smooth crater profiles can be well represented with smaller Chebyshev approximation order. Chebyshev approximation of two contrasting crater profiles (Linne, 2 km diameter with symmetrical and smooth topography and Larmor Q, 19 km diameter with asymmetric and irregular topography) (Fig. 6) is performed. For reconstruction, only the first 17 Chebyshev polynomials were used (M ¼ 16) resulting in a compression ratio (coefficients to true samples) of about 99% and a spatial relative reconstruction error (SRRE) magnitude less than 20% (less than 5% for most parts). Due to the irregular shape of the Larmor Q profile, a larger average reconstruction error is observed over the entire profile (about 3% compared to 0.5% for Linne). Larger reconstruction errors occur at points where the elevation profile changes sharply such as at the crater rims. This is because, a sharp change in a function shape is associated with higher spatial frequency components that can be represented only by using higher-order Chebyshev polynomials (Fig. 2 shows the change in spatial frequency for T N as N increases). To include higher-order Chebyshev polynomials in the approximation the overall order of approximation (M) must be made sufficiently large and if the order of approximation is not large enough, reconstruction errors are observed. For analyzing the local behavior of the approximation error (at each spatial point on the profile), the error measure used (SRRE) is a ratio of the absolute value of the local error and the absolute value of the local elevation. SRRE is more sensitive at elevation values close to zero as the denominator values get smaller and so for a global measure of error (average reconstruction error across the entire profile), this metric is normalized. The average reconstruction error is computed by a weighted average of the SRRE. The weights are selected proportional to the fractional energy contributed by local elevation to the whole profile.

ð11Þ

If only 17 coefficients (M ¼ 16) are used to approximate the profile, the elevation difference between the average reference baseline elevation and the elevation at x ¼ 0 can be approximated by 2ðC 2 þ C 6 þ C 10 þ C 14 Þ and by only 2C 2 for a near parabolic profile. For small values of M; I5 is expected to show reduced correlation with elevation difference between the average reference baseline elevation and the elevation corresponding to x ¼ 0. However, for simple small craters, the contribution of higher numbered coefficients is not significant and stronger correlations are observed, as demonstrated in the earlier sections. The Chebyshev coefficient based crater morphology indicators discussed in this work (summarized in Table 3) are obtained from crater elevation profiles acquired by the standardized method. We derived the indicators I1 to I5 to represent combinations of

Fig. 7. Variation of percent approximation error (normalized SRRE) with order of approximation (M) and crater diameter (D). Value of M changes in steps of 4 from 4 (dark blue) to 32 (dark red). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

123

P. Mahanti et al. / Icarus 241 (2014) 114–129

Elevation profiles from different sized craters (total of 56 craters, 300 m < D < 140 km) were approximated using Chebyshev polynomials for 8 different approximation orders (M ¼ 4; 8; 12; 16; 20; 24; 28; 32). SRRE (expressed as percentage) was used as the comparison metric. It was found that for the craters analyzed, a choice of M ¼ 16 led to less than 10% overall approximation error for 54 of these craters (Fig. 7). Smaller craters (D < 500 m) produced lower reconstruction errors even at M ¼ 8 or M ¼ 12. For craters with diameter values between 500 m and 1 km, the SRRE values for M ¼ 4 were high. On an average, SRRE values for M ¼ 4 were even higher for craters with 1 km < D < 100 km and highest for craters with D > 100 km. At approximation order of 32, 53 of the 56 craters had SRRE < 5% (remaining 3 craters had SRRE between 6% and 8%). Similar to the decrease in SRRE values with increasing values of M, a decreasing trend is also observed at each M with decreasing crater size (D). Another observation is that the rate of decrease of SRRE values is smaller for bigger craters (in Fig. 7 compare the rate of decrease for D > 100 km and 500 m < D < 1 km). Possible reasons for this include the flattening of the crater floors which is more difficult to approximate using smoothly varying functions like the Chebyshev polynomials. Choice of an approximation order only gives a maximum bound on M beyond which the coefficients are not expected to possess sufficient energy and are not as important for morphological comparisons. It is emphasized that while M ¼ 16 provides for a good approximation, the order or approximation can vary with craters and with requirement. The dataset analyzed above only considers lunar craters and a larger number for craters are required for a more robust inference about the best choice of approximation order. Performance also needs to be evaluated against other metrics possible for reconstruction. 4. Discussion The principal advantages of the proposed standardized crater topographic profile characterization are as follows. (1) The first few Chebyshev basis functions represent different geomorphic aspects of an impact crater (Fig. 1A). Knowledge of the basis function shapes and their combinations enable intuitive understanding of the contribution of the basis functions to the overall crater shape. (2) Chebyshev approximation coefficients can be obtained independently without interfering with another coefficient and selectively such that only required topographic indicators are obtained. (3) The physical meaning of a Chebyshev approximation coefficient is independent of the total number of coefficients used for the approximation. Thus, similar numbered coefficients can always be compared across less detailed (small M) and more detailed (large M) crater representations. (4) The Chebyshev representation is compact and the coefficient absolute magnitudes decrease for higher numbered coefficients. Since most of the energy (topographic variation) is concentrated in the first few coefficients, a small set of coefficients is sufficient to represent a crater elevation profile – converting a crater elevation profile into a feature vector (with as many components or dimensions as the number of coefficients and (or) morphological indicators used). (5) The basis function inner product approach provides for a simple method of obtaining the Chebyshev coefficients and their combinations and is easily implementable as a software tool. 4.1. Uniqueness and differences of proposed method from earlier work The published methods for crater elevation profile analysis using mathematical relationships available prior to this work were briefly introduced earlier (Section 1). These methods fall into one of the following two categories. The first is where a geometric

relationship between fixed values of crater features is studied. Examples of such relations include the power law relations used by Pike (1974a,b, 1976, 1977a,b):

log ½d ¼ log ½a þ b log ½D

ð12Þ

where a; b are coefficients of the regression. Another example is the polynomial in log domain used by Baldwin (1963, 1965): 2

log ½D ¼ 0:0256ðlog ½dÞ þ log ½d þ 0:6300

ð13Þ

These relationships and similar others describe scaling relationships between fixed values of depth and diameter obtained by regression and can describe and classify groups of craters based on values of the fitted coefficients. However, these forms cannot represent the complete elevation profile information as the variables are fixed (for example d is rim-to-floor depth and D is rim-to-rim diameter). So they can represent a population of craters but not individual craters, which is possible only via an elevation profile function approximation as in the proposed method. Geometric relationships are and will be useful for describing populations of craters but not for describing the details of each crater separately. By using the proposed method, we can describe populations of craters (by similarity of coefficients or their combinations) as well as distinguish between individual craters. Further, we can reconstruct the elevation profile from stored coefficients for a more detailed investigation. The second category of mathematical relations are rational and irrational polynomial relations between crater profile elevation and radial distance from the center of crater. An example of using irrational polynomials is Simon and de Bruyn (2007)

z ¼ zc þ

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 b þ c2 r 2

ð14Þ

where z is the profile elevation and zc ; b and c are fitting constants. While the selection of such a hyperbolic parametric fit was needspecific to the study, it is found that the fit diverges at the crater rim and at the central peak (the method cannot describe the complete crater profile). Similar other mathematical expressions can be formulated for describing morphological features of a crater but their applicability is often restricted (only applicable to a particular population of craters or part of the crater) or may require a geologic interpretation during fitting (fitting model not purely based on objective rules). A more complete description is found by using piecewise rational polynomials (Craddock and Howard, 2000) or a simple second-order polynomial (quadratic fit) (Talpe et al., 2011). In the former example for a fresh crater two separate fourth-order polynomial equations are used. One for describing the profile from center to rim and another to describe the profile from rim crest to edge of ejecta. For a degraded crater, a single third-order polynomial was used. This method intuitively uses the polynomial expressions but is still a regression (describes groups of craters and not individual craters themselves). Two different morphologies are treated separately (fourth and third order expressions), and because the polynomial orders are not higher than four, fitting of craters with central peak, wall collapses or fractured floors will have large approximation errors (due to sharp elevation function changes requiring higher order terms). For the second example, the use of the quadratic term by Talpe et al. (2011) as a quantifier of the relationship between the height and width of a parabola is analogous to the use of C 2 (coefficient of a second order Chebyshev basis) in this work as a simple indicator of crater depth. Polynomials that are standard, applicable to all craters (if enough terms are used), and follow objective rules are intuitive choices for representing crater elevation profiles. Any polynomial fit of the same order will be the same, but the coefficients of standard polynomials are not as well correlated to morphological features as the Chebyshev coefficients. The Chebyshev approximation method uses objective

124

P. Mahanti et al. / Icarus 241 (2014) 114–129

rules, does not require any human intervention to calculate, and is highly extensible to a potentially huge set of craters. Coefficients (and their combinations) derived from the Chebyshev approximation for crater profiles compactly represent a crater making it possible to convert a population of craters into a table of morphological descriptors. Each row of such a table corresponds to a crater (feature vector for that crater) and each column lists either a Chebyshev coefficient or a morphological feature indicator obtained from the Chebyshev coefficients. Such a standardized tabulation enables a simplified classification scheme for craters. If only Chebyshev coefficients are used then due to compactness of representation (Section 2.2.3), the feature vector length to be used for comparison is typically small. The Chebyshev approximation for crater profiles can do everything an ordinary polynomial can, and additionally satisfies all the requirements for an effective generalized approach for crater comparison and morphological study. The difference between parametric fitting of a single complex function (fit topographic data to preconceived model) and Chebyshev function approximation (build topographic profile model as-per-data) makes Chebyshev approximation a powerful method for crater profile characterization. Useful properties of the Chebyshev polynomials leading to ease of coefficient evaluation and intuitive topographic shape relevance (distinction from trigonometric basis) make them additionally attractive. In the discussion that follows, the advantages of using the proposed representation for crater profiles is elaborated. 4.2. Physical significance of Chebyshev coefficients Coefficients obtained from the Chebyshev approximation correspond to the geomorphic characteristics of a crater profile. The higher contribution (Fig. 6C and D) for lower numbered coefficients indicate their stronger relationship with the main profile shape and the higher numbered coefficients add details to the main shape. An irregular crater profile (Fig. 6B, Larmor Q) represents a complex function (more high frequency components) shape, whereas a smooth simple crater (Fig. 6A, Linne) represents a simple function shape. Accordingly, when the same number of coefficients are used, the reconstruction error is smaller for the smooth simple crater. The relative compactness of representation is illustrated in Fig. 6C and more significant higher numbered (18 and larger) coefficients (representing high frequency components) are observed for Larmor Q. Accordingly, the total normalized contribution of the coefficients not shown in this figure (total power from coefficients 18 and larger) will be larger than that for Linne. The normalized contribution of coefficients serves as a measure of relative profile smoothness, if the profile sampling densities (meters/ pixel) are constrained to be equal. The selection M ¼ 16 was found sufficient for the characterization of all craters (via the morphologic indicators, see Table 3) considered for this work. The 1st Chebyshev coefficient (C 0 ) is positively correlated to the average topographic surface elevation for the crater rim and not just the average elevation from the crater topographic profile. This correlation increases with the radial symmetry and ‘simple’-ness of the crater and decreases as craters become complex and asymmetric. The sign of the 1st coefficient indicates whether the topography is above (+) or below () the local datum. The nature of topographic variation (around the local datum) for these two craters is indicated by the coefficient value of the 1st Chebyshev coefficient (Fig. 6C and D). The log coefficient power for C 0 is approximately 2, indicating a contribution of more than 90% of the total Chebyshev representation (signal power) for the Linne and Larmor Q crater elevation profiles. Most planetary impact craters show radial symmetry, however some craters (for example on sloping surfaces) exhibit an

asymmetric interior morphology (Plescia, 2012). If we consider a symmetric crater formed over a flat pre-impact surface the crater profile has minimal radial asymmetry and this corresponds to a small 2nd Chebyshev coefficient (C 1 ) in the Chebyshev approximation of the crater elevation profile. Conversely, for a crater on a sloping surface the 2nd Chebyshev coefficient is expected to be large. The magnitude of the 2nd Chebyshev coefficient indicates the local gradient of the pre-impact plane outside the crater rim. A local gradient can potentially impede the computation of features (Baker et al., 2012; Williams and Zuber, 1998) but the Chebyshev approximation of the crater elevation profile enables a simple correction – removal of the local gradient by setting the 2nd Chebyshev coefficient to zero (correction of Giordano Bruno crater profile in Fig. 5B and then reconstructing the crater profile (using Eq. 4). The 3rd Chebyshev coefficient (C 2 ) correlates positively with the depth of the crater. However, a more accurate indicator of the crater depth is given by the combination of Chebyshev coefficients expressed by I5 (Eq. 11, Section 3.3). Similarly, where a sum of C 2 and C 0 can be used to estimate the elevation of the pre-impact surface (Fig. 1B) a more accurate indicator of the same is provided by the magnitude 12 ðI3 þ I4 Þ, as discussed in Section 3.3. It may be noted that for radially symmetric craters the correlation between 1-dimensional crater features (average elevation, depth, pre-impact surface) obtained from profiles and 3-dimensional elevation feature measurements obtained from gridded DEMs will be stronger. Each Chebyshev basis function has a unique shape and spatially varying magnitude along the scaled sample profile axis, enabling selective analytical treatment of different segments of the crater profile. An example application was discussed earlier (Fig. 4A) where the sum of specific coefficients was utilized to sort craters based on their central peaks, showing the physical relevance between the crater profile and groups of basis functions with stronger magnitudes at the center (selection of desired location). Note that any morphological activity that modifies the central elevation feature (for example central pit or summit pit features) of the crater will affect the signs and magnitudes of the above mentioned Chebyshev coefficients (redistribution of the topographic profile energy), but such a study is currently beyond the scope and will be part of a future follow-up work. Similarly, by analysis of the basis function shapes, other coefficient groups can be effectively utilized to obtain indicators (and hence classifiers) for other crater features. 4.3. Depth-to-diameter indicator from Chebyshev coefficients So far in this paper, quantitative descriptors that are directly linked to crater morphologic features were obtained from Chebyshev coefficients, either individually or as a group. However, Chebyshev approximation also facilitates analysis of crater degradation that causes topographic variations. As an example, depthto-diameter ratio of craters   is explored and a new indicator for depth-to-diameter ratio Dd obtained from Chebyshev coefficients is introduced in this work. The definitions for depth (d) (rim-crest to floor) and diameter (D) (rim-crest diameter or rim-to-rim) values and the measurement methodology are adapted from Pike (1976). Depth-diameter relationships of lunar craters were examined at least as early as 1963, and deterministic models (with uncertainty margins) such as the power law were used in the works of Baldwin, Pike, Elachi and others (Baldwin, 1963; Baldwin, 1965; Elachi et al., 1976; Pike, 1977a). Relationships between depth and diameter were described by Baldwin (1963) as a 3rd order polynomial in the log-domain. Pike (1967), suggested a power law was a better description of the relationship between these quantities.

P. Mahanti et al. / Icarus 241 (2014) 114–129 Table 4 Selected craters, depth-to-diameter values and Chebyshev approximation based indicators of depth-to-diameter ratio. Name

d (km)

D (km)

d D

C2

eI 5

JansenE TaruntiusP Toscanelli AuwersA LickE SantosDumont Norman Rosse MostingA Hortensius ProclusA Glaisher Birt Gardner Black BinghamH Romer AbulWafa Ritz Macrobius Reichenbach Werner King Tycho Aristoteles Copernicus Theophilus Langemak Skladowska Langrenus Albategnius Neper

1.38 1.43 1.48 1.44 2.04 1.43 1.98 2.42 2.82 2.84 2.85 2.93 3.62 2.93 3.13 2.68 2.92 3.54 3.24 3.58 3.14 4.08 4.28 4.44 3.41 3.66 4.19 4.11 3.94 4.14 3.27 3.48

7.12 7.37 7.48 7.85 8.78 8.96 10.52 11.72 12.77 14.57 14.89 16.59 16.59 17.90 19.80 36.55 41.49 54.73 60.60 62.33 69.05 70.08 75.45 84.29 86.70 93.09 98.97 104.42 124.46 132.65 138.02 144.68

0.19 0.19 0.20 0.18 0.23 0.16 0.19 0.21 0.22 0.20 0.19 0.18 0.22 0.16 0.16 0.07 0.07 0.06 0.05 0.06 0.05 0.06 0.06 0.05 0.04 0.04 0.04 0.04 0.03 0.03 0.02 0.02

27.61 39.58 27.37 28.62 39.54 10.59 30.72 37.21 39.02 31.43 32.43 47.18 33.97 32.13 41.06 22.66 19.50 11.98 4.02 10.66 5.58 8.35 1.55 5.71 8.91 6.99 5.79 7.37 7.60 5.93 1.46 4.31

48.36 59.52 52.62 53.26 73.31 28.81 50.41 60.85 66.18 53.52 49.23 64.07 63.46 45.72 55.52 29.68 15.38 18.03 8.64 13.27 8.67 11.29 0.24 5.88 9.88 7.78 4.26 9.29 8.15 6.65 2.64 2.90

wide variety of craters was chosen to verify the adaptability of the Chebyshev polynomial analysis scheme. The selected crater profiles varied from V-shaped profiles to U-shaped to even more degraded examples. Some of the elevation profiles showed variations of central peak morphology and some were asymmetric with respect to the pre-impact plane. The diameters of the selected craters varied from 7 km to 160 km and their depths varied from 1 km to 4.5 km. Obtaining Chebyshev coefficients for the selected craters is exactly similar as discussed earlier, except the crater elevation profiles were first normalized by the diameter of the crater (obtained during the elevation profile extraction process). Analytically, this can be represented as M f ðxÞ X CN ffi T N ðxÞ D D N¼0

Subsequently, Apollo data were used to fit both power law (Pike, 1977b) and linear relations (Elachi et al., 1976). Historical data from an earlier study (Pike, 1976) obtained primarily from Apollo-based lunar topographic orthophotomaps (LTO) is used as the basis of reference in this work. Crater elevation profiles selected from this historical reference were extracted from LROC based DEMs and analyzed using the Chebyshev decomposition. A subset (32 craters) of the complete list (Table I in Pike (1976)) of primary impact craters was selected (Table 4) for our analysis. The complete set of craters (as well as the selected subset) comprises two sub-fields (Fig. 8), which overlap at a pronounced elbow that exists between crater diameters 10 km and 20 km. A

ð15Þ

Since the right-hand-side of the above expression is a linear sum (Chebyshev approximation process is a linear expansion), this causes the Chebyshev coefficients themselves to be normalized with respect to the crater diameter. Expanding the right-hand-side we get

~f ðxÞ ffi C 0 T 0 ðxÞ þ C 1 T 1 ðxÞ þ . . . þ C N T N ðxÞ D D D

ð16Þ

where ~f ðxÞ denotes the normalized elevation profile. A further simplified form of the above expression is

~f ðxÞ ffi

M X

e N T N ðxÞ C

ð17Þ

N¼0

where ðC~N Þ denotes the normalized Nth Chebyshev coefficient. The above expression is similar to Eq. 4 except both the elevation profile as well as the Chebyshev coefficients are now normalized, but the form of the Chebyshev approximation is unaltered. In Section 3.3, various indicators composed of multiple Chebyshev coefficients were identified. Normalization by the crater diameter prior to the Chebyshev approximation thus normalizes the indicators – an indicator of depth (I5 – Eq. 11 in Section 3.3) is now transformed to an indicator of depth-to-diameter ratio. An approximation of I5 (without normalization) for M ¼ 16 is given by the sum C 2 þ C 6 þ C 10 þ C 14 . The normalized indicator eI 5 is then given by adding the same coefficients as earlier but the effective sum is

eI 5 ffi C 2 þ C 6 þ C 10 þ C 14 D Fig. 8. Data-points (craters) from (Pike, 1976) analyzed using Chebyshev approximation. The scatter-plot shows the variation of depth-vs-diameter for all the craters on a log–log scale. Note the typical ‘elbow’ that separates the craters into two overlapping clusters (diameters below 10 km and above 20 km), each of the clusters can be fitted via a straight line. All the data-points relate to primary impact craters (Table I in Pike (1976)). The data-points marked in red are the craters that were selected (Table 4) for analysis. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

125

ð18Þ

In addition to the implicit normalization of the Chebyshev coefficients, normalization by crater diameter causes a reference axis transformation (from true spatial distance between elevation samples to a scaled ratio). Since morphological features for the same crater are proportional (a crater with large diameter is expected to have a large depth), normalized Chebyshev coefficients can be analyzed for large and small craters simultaneously over a linear scale, eliminating need for a log-based scale. A summary of the measurements obtained for the selected crater elevation profiles (Fig. 9) is presented in Table 4. The first four columns list the crater names, depth values, diameter values, and the corresponding Dd values. The fifth and sixth columns list the value of the 3rd Chebyshev coefficient (C 2 ) and eI 5 respectively. Depth-to-diameter ratio from the fourth column (Table 4) and the corresponding Dd values obtained from Pike (1976) LTO based observations are plotted against eI 5 (Fig. 10). Craters with Dd < 0:1 and those with Dd > 0:1 are observed in separate clusters along the principal diagonal of the plot. Moreover, the relationship between eI 5 and d/D is linear and thus eI 5 can be successfully utilized as an indicator for Dd . Both Dd data sets show a linear relationship

126

P. Mahanti et al. / Icarus 241 (2014) 114–129

Fig. 9. Elevation profiles of the 32 craters selected (Table 4) for Chebyshev polynomial approximation. Craters were selected from a database used to analyze depth-diameter relationships (Pike, 1976, 1977a).

Fig. 10. eI 5 (x-axis) used as an indicator of depth-to-diameter (y-axis) values. The plot shows a linear relationship between eI 5 and depth-to-diameter values. The size of the craters is indicated by color (for LRO data) – larger craters occupy the lower left quadrant of the plot. The black dots represent Dd values from earlier measurements (Pike, 1976) corresponding to the Dd values from LROC WAC global DEM based measurements. The definitions for depth (d) and diameter (D) is shown for a hypothetical crater elevation profile (inset). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

with eI 5 , but the LROC WAC global DEM based Dd values exhibit lesser spread (lower root-mean-square error (RMSE) value for the fit, also note the vertical shift of the Apollo data markers to LRO data markers) owing to the higher resolution of the topographic data. eI 5 values less than 20 indicate d values less than 0.1 and eI 5 valD ues more than 20 indicate Dd values greater than 0.1. Larger craters are correlated with smaller eI 5 and Dd values. Additionally, since Dd values are indicators of relative age of craters, eI 5 can also be used as a relative age indicator (both Dd and eI 5 magnitudes decreasing with increasing relative age). The scaled 3rd Chebyshev coefficient shows the same trend as eI 5 , but eI 5 is a more robust estimator since it includes more Chebyshev coefficients and this is illustrated in Table 4. By using the scaled C 2 value alone, the two clusters are not as well For example the scaled C 2 value for Santos  segregated.  Dumont Dd > 0:1 is less than 4 other craters in Table 4 (rows 16, 17, 18, 20) which have Dd < 0:1. Adding the other three Chebyshev coefficients, helps in separating the clusters. Among the craters analyzed, some of the crater elevation profiles were asymmetric while some showed central peak morphology. Our results indicate that such deviations from the simple symmetric bowl-shaped crater elevation profile has little effect on the correlation of Dd with eI 5 . The reason is that coefficients selected to form I5 (and hence eI 5 ) are neither effected by asymmetry (C 1 and other odd-numbered Chebyshev coefficients) nor by central peak morphology (which affects C 4 ; C 8 ; C 12 ; C 16 ). Thus, by judicious selection of coefficients, the independence of the Chebyshev coefficients (and their combinations) can be exploited effectively to independently analyze crater morphological

P. Mahanti et al. / Icarus 241 (2014) 114–129

127

Fig. 11. Approximation error (M ¼ 8) for three craters of different sizes and degradation stages. Crater A appears fresher compared to C which looks more degraded from the images.

indicators, for example, the depth of the crater can be analyzed even if the crater is on a sloped pre-impact surface and has a complex central peak.

4.4. Crater degradation and crater rim approximation As a crater degrades over time, the sharpness of the rim is reduced and the crater profile elevation function becomes smoother and can be well approximated by a smaller order polynomial (for a rough function a higher order polynomial performs better). Hence, for a fixed order of approximation, lesser approximation error (SRRE) would be observed at the crater rims for a more degraded crater. Three craters with diameters 957 m, 414 m, and 610 m were analyzed and their depth and diameters were computed from LROC NAC DEMs (Fig. 11). The change in maturity of the craters is observable from the images and also from the Dd values. The average approximation error (normalized SRRE, see Section 3.4 for description) increases with decreasing crater maturity (increasing Dd value). Selective spatial regions of the crater profile can be analyzed by the proposed Chebyshev approximation. The standardized procedure of obtaining the crater profile ensures that the crater rims are between 0.4 to 0.8 and +0.4 to +0.8 on the scaled sample reference axis. The approximation error at the rims for different order of approximation (Fig. 12) shows that fresher craters require a higher order of approximation for attaining the same average elevation error. Alternately, approximation error magnitudes at crater rims for a lower order approximation can provide indication of crater rim smoothness from which degradation and maturity may be inferred. At higher approximation orders (M > 16), average approximation errors are similar.

Fig. 12. Variation of approximation error (normalized SRRE) with change in order of approximation for craters A, B and C. More degraded craters show a quicker decrease of approximation error. For crater examples B and C the approximation error is smaller than 10% at M ¼ 4. The green shaded regions (inset, crater C approximation) indicate the spatial extent over which the average approximation error is measured. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

4.5. Uncertainties in the proposed method In the earlier sections it was shown why the proposed standardized representation of crater profile is a promising choice. In this section we discuss the limitation of the proposed method and present an analysis of three classes of uncertainties and their significance.

128

P. Mahanti et al. / Icarus 241 (2014) 114–129

The first category is the order of Chebyshev approximation and the associated approximation error. The order of approximation (M) is essentially left to the choice of the user, and so for the same extracted profile, two different users can have varying truncation (approximation) error. However, the magnitude of the error is completely characterized by the value of M and the approximation coefficients are still comparable up to the available corresponding term (or even combinations of available coefficients). For comparison across crater profiles, the order of Chebyshev approximation varies based on the smoothness of the elevation profile but in our analysis presented earlier (see Section 3.4) a value of M ¼ 16 was found to consistently generate average spatial residuals less than 5%. Similarly, by restricting the average spatial residual magnitude, a robust lower bound on M can be obtained by trial and error for a given set of crater profiles. However, in most cases the compactness property of Chebyshev approximation ensures that the value of M will not be arbitrarily high, for example, a value of M larger than 50 is usually not expected. In essence, the approximation error can be minimized by maintaining a relatively large value of M as long as the topographic data is from the same source (has same resolution and similar elevation uncertainty magnitudes). Note that here it was assumed that the standardized extraction method was followed. In our future work, planetary impact craters from various origin will be compared by the proposed method and the choice of approximation order will be discussed for a larger population of craters. The second category are the uncertainties associated with the standardized crater profile extraction method and measurement error. While the crater profile extraction method is standardized, there can be situations where more than one possible profiles can be drawn through the center of the crater and parallel to the maximum absolute gradient direction (from higher to lower). The effect increases as craters become less axisymmetric. This is a limitation but can be mitigated by deciding on one particular direction (N–S or E–W) in addition to the previous constraints. It is also possible that by the proposed standardized profile selection, a crater feature in an orthogonal direction is not accounted for in the Chebyshev approximation. Possible methods to include such a feature are averaging across multiple profiles, or considering multiple profiles to obtain multiple sets of Chebyshev coefficients. Note that the performance of the Chebyshev approximation process remains the same, irrespective of the profile selection direction. If a particular crater has a local feature in the N–S direction and not in the E–W direction, then the feature will affect more the Chebyshev coefficients for a N–S direction profile. Thus, when required, multiple topographic profiles may be used (for example 4 topographic profiles: E–W, N–S and the two diagonals), but how far this is advantageous will be analyzed in future work. Profiles drawn off-center also contribute to this second category of ambiguity. Measurement errors due to off-center profiles are expected to be enhanced for craters with central features (for example central peak) and considerably less for flat-floored craters (assuming small offsets from the center). Some uncertainties can be implicitly associated with the DEM being used and represent measurement error. In such cases, minimum elevation error for the DEM is an important consideration as it affects the magnitude of the smallest topographic feature that can be represented. The elevation uncertainties in the DEM propagate into the evaluated Chebyshev coefficients. Large coefficients (typically lower order Chebyshev coefficients) are less affected as they represent a stronger topographic variation (higher order coefficients are more affected). Variation of the topographic resolution across craters and (or) users can also lead to measurement errors for the coefficients themselves and affects the smaller magnitude coefficients the most (see Section 3.1). The effect of topographic resolution on Chebyshev coefficients will be analyzed in a future

work. As a current mitigation strategy, use of either the highest resolution topographic data available or pre-deciding on a fixed available resolution is suggested. The last category is the interpretation uncertainty associated with interpreting the values of the coefficients against the actual morphology of the crater. For example, while the value of the 3rd Chebyshev coefficient (C 2 ) reflects the crater depth (measured centrally), it is not as accurate as I5 and accuracy of I5 depends on how much the extracted topographic profile endpoints represent the elevation of the plane around the crater. The central depth is also affected by the presence of topographical features at the center of the crater, for example a central peak. While it is possible to adapt the proposed method to central features in the future, careful considerations must be made when such craters are analyzed by the proposed method and also when a large number of craters are analyzed together. In some cases, the above mentioned reasons will produce outliers as evident from the outliers in the plot of crater depth against C 2 (Fig. 4) and the plot of Dd against eI 5 (Fig. 10). The use of combinations of coefficient values for standardized comparison across craters is often more efficient than the direct use of the coefficients values. It is emphasized here that the proposed method does not extract true values of morphological features (for example crater maximum depth, central peak height) but provides an alternative to direct (and often manual) measurement by making corresponding morphological indicators available via a repeatable automatic or semi-automatic process. Note that since the entire crater elevation profile is approximated by the proposed method, precise computations of true values of morphological features is possible, at least in theory (since the information already exists) from using the Chebyshev approximation coefficients and will be investigated in future work. 5. Conclusion Unique and physically meaningful analysis of crater topography can be achieved by the use of Chebyshev polynomials. Chebyshev polynomials offer distinctive advantages of mutual independence of Chebyshev coefficients, ease of coefficient evaluation, relationship to usual crater features and compactness of representation relative to simple polynomials. Crater topographic profile analysis presented in this work shows that the first 17 Chebyshev coefficients (M ¼ 16) (or a subset of the first 17 coefficients) are sufficient to represent, characterize and even classify craters. Examples shown in this work indicate the possibility of a generic framework of automated or semi-automated crater topographic profile analysis based on Chebyshev polynomial representation. This standardized method of crater analysis shows potential for the classification, comparison and analysis of impact craters across various planets, domains and researchers - thus improving our understanding of impact phenomena from remotely sensed topography. Acknowledgments We gratefully acknowledge the work of LRO and LROC operations personnel who helped in generating the high-resolution lunar topographic data. Thanks are extended to Richard Stelling for help in data-mining. Finally, we thank the reviewers for their constructive criticism and thoughtful reviews. This work was supported by the National Aeronautics and Space Administration Lunar Reconnaissance Orbiter Project. References Baker, D., Head, J., Neumann, G., Smith, D., Zuber, M., 2012. The transition from complex craters to multi-ring basins on the Moon: Quantitative geometric

P. Mahanti et al. / Icarus 241 (2014) 114–129 properties from lunar reconnaissance orbiter lunar orbiter laser altimeter (LOLA) data. J. Geophys. Res. 117, E00H16. Baldwin, R., 1963. The Measure of the Moon. University of Chicago Press, Chicago (1). Baldwin, R., 1965. A Fundamental Survey of the Moon. McGraw-Hill, New York (1). Beckmann, P., 1973. Orthogonal Polynomials for Engineers and Physicists, vol. 153. Golem Press, Boulder, CO. Blaser, A., Holman, R., Brown, D., Willis, J., 2011. Slope diffusion models and digitally-acquired morphometric parameters yield age constraints on cinder cones, examples from the Spencer High Point and Craters of the Moon National Monument, Snake River Plain, Idaho. In: AGU (Fall Meet.), vol. 1, p. 0706 (abstracts). Bray, V., Collins, G., Morgan, J., Schenk, P., 2008. The effect of target properties on crater morphology: Comparison of central peak craters on the Moon and Ganymede. Meteorit. Planet. Sci. 43 (12), 1979–1992. Cavanaugh, J. et al., 2007. The Mercury laser altimeter instrument for the messenger mission. Space Sci. Rev. 131 (1), 451–479. Cintala, M., Grieve, R., 2010. Scaling impact melting and crater dimensions: Implications for the lunar cratering record. Meteorit. Planet. Sci. 33 (4), 889– 912. Cook, A., Watters, T., Robinson, M., Spudis, P., Bussey, D., 2000. Lunar polar topography derived from clementine stereoimages. J. Geophys. Res 105 (E5), 12–023. Craddock, R.A., Howard, A.D., 2000. Simulated degradation of lunar impact craters and a new method for age dating farside mare deposits. J. Geophys. Res.: Planets (1991–2012) 105 (E8), 20387–20401. Elachi, C., Kobrick, M., Roth, L., Tiernan, M., Brown, W., 1976. Local lunar topography from the Apollo 17 ALSE radar imagery and altimetry. Earth Moon Planets 15 (1), 119–131. Elkins-Tanton, L., Hager, B., Grove, T., 2004. Magmatic effects of the lunar late heavy bombardment. Earth Planet. Sci. Lett. 222 (1), 17–27. Garvin, J., Frawley, J., 1998. Geometric properties of martian impact craters: Preliminary. Geophys. Res. Lett. 25 (24), 4405–4408. Garvin, J., Sakimoto, S., Frawley, J., 1998. Martian impact crater ejecta topography as measured by the Mars orbiter laser altimeter (MOLA). Lunar Planet. Sci. XXIX. Garvin, J., Sakimoto, S., Schnetzler, C., Frawley, J., 1999. Global geometric properties of martian impact craters – A preliminary assessment using Mars orbiter laser altimeter (MOLA). In: The Fifth International Conference on Mars, Pasadena, CA. Gautschi, W., 2004. Orthogonal Polynomials: Computation and Approximation. Oxford University Press, USA. Gil, A., Segura, J., Temme, N.M., 2007. Numerical Methods for Special Functions. Siam. Hartmann, W., 1972. Interplanet variations in scale of crater morphology Earth, Mars, Moon. Icarus 17 (3), 707–713. Hewitt, E., 1947. Certain generalizations of the Weierstrass approximation theorem. Duke Math. J. 14 (2), 419–427. Mason, J.C., Handscomb, D.C., 2010. Chebyshev Polynomials. CRC Press. Mazarico, E. et al., 2010. Depth-diameter ratios of small craters from LOLA multibeam laser altimeter data. In: Lunar and Planetary Institute Science Conference, vol. 41, p. 2443 (abstracts). McEwen, A. et al., 2007. Mars reconnaissance orbiter’s high resolution imaging science experiment (HIRISE). J. Geophys. Res. 112 (E5), E05S02. Nguyen, N., Burken, J., Ishihara, A., 2011. Least-squares adaptive control using Chebyshev orthogonal polynomials. In: AIAA Infotech@ Aerospace Conference, AIAA-2011-1402. Pike, R., 1967. Schroeters rule and the modification of lunar crater impact morphology. J. Geophys. Res. 72 (8), 2099–2106. Pike, R., 1974a. Craters on Earth, Moon, and Mars: Multivariate classification and mode of origin. Earth Planet. Sci. Lett. 22 (3), 245–255. Pike, R., 1974b. Depth/diameter relations of fresh lunar craters: Revision from spacecraft data. Geophys. Res. Lett. 1 (7), 291–294.

129

Pike, R., 1976. Crater dimensions from Apollo data and supplemental sources. Earth Moon Planets 15 (3), 463–477. Pike, R., 1977a. Apparent depth/apparent diameter relation for lunar craters. Proc. Lunar Sci. Conf. 8, 3427–3436. Pike, R., 1977b. Size-dependence in the shape of fresh impact craters on the Moon. In: Impact and Explosion Cratering: Planetary and Terrestrial Implications, vol. 1, pp. 489–509. Plescia, J., 2012. Asymmetric lunar craters-impacts on sloping surfaces. In: 2012 GSA Annual Meeting in Charlotte. Pommerol, A., Chakraborty, S., Thomas, N., 2012. Comparative study of the surface roughness of the Moon, Mars and Mercury. Planet. Space Sci. 73 (1), 287–293. Powell, M., 1967. On the maximum errors of polynomial approximations defined by interpolation and by least squares criteria. Comput. J. 9 (4), 404–407. Rivlin, T.J., 2003. An Introduction to the Approximation of Functions. Dover Publications.com. Robbins, S.J., Hynek, B.M., 2012a. A new global database of Mars impact craters 1 km: 1. Database creation, properties, and parameters. J. Geophys. Res.: Planets (1991–2012) 117 (E5). Robbins, S.J., Hynek, B.M., 2012b. A new global database of Mars impact craters 1 km: 2. Global crater properties and regional variations of the simple-tocomplex transition diameter. J. Geophys. Res.: Planets (1991–2012) 117 (E6). Robinson, M. et al., 2010. Lunar Reconnaissance Orbiter Camera (LROC) instrument overview. Space Sci. Rev. 150 (1), 81–124. Roddy, D., 1976. High-explosive cratering analogs for bowl-shaped, central uplift, and multiring impact craters. Proc. Lunar Sci. Conf. 7, 3027–3056. Ross, H., 1968. A simplified mathematical model for lunar crater erosion. J. Geophys. Res. 73 (4), 1343–1354. Scholten, F. et al., 2012. Gld100: The near-global lunar 100 m raster DTM from LROC wac stereo image data. J. Geophys. Res. 117, E00H17. Schultz, P., 1976. Floor-fractured lunar craters. Earth Moon Planets 15 (3), 241–273. Simon, J., de Bruyn, J.R., 2007. Shape of impact craters in granular media. Phys. Rev. E 76 (4), 041306. Simpson, J., Kim, J., Muller, J., 2008. 3d crater database production on Mars by automated crater detection and data fusion. In: ISPRS Congress 2008 Proceedings. Smith, D. et al., 2010. The lunar orbiter laser altimeter investigation on the lunar reconnaissance orbiter mission. Space Sci. Rev. 150 (1), 209–241. Smith, D., Zuber, M., Neumann, G., Lemoine, F., 1997. Topography of the Moon from the clementine lidar. J. Geophys. Res. 102 (E1), 1591–1611. Smith, E.I., Hartnell, J.A., 1978. Crater size-shape profiles for the Moon and Mercury: Terrain effects and interplanetary comparisons. Moon Planets 19 (4), 479–511. Sori, M., Zuber, M., 2010. Preliminary measurement of depth-to-diameter ratios of lunar craters in the transition regime between complex craters and multiringed basins. In: Lunar and Planetary Institute Science Conference, vol. 41, p. 2202 (abstracts). Sori, M., Zuber, M., 2012. Anomalous shallowing of lunar impact craters in the south pole-aitken basin from lunar orbiter laser altimeter (LOLA) observations. In: Lunar and Planetary Institute Science Conference, vol. 43, p. 2707 (abstracts). Stopar, J., Robinson, M., Barnouin, O., Tran, T., 2010. Depths, diameters, and profiles of small lunar craters from LROC NAC stereo images. In: AGU (Fall Meet.), vol. 1, p. 1543 (abstracts). Talpe, M., Zuber, M., Clark, M., Mazarico, E., 2011. Regional cataloguing of lunar crater morphology. In: Lunar and Planetary Institute Science Conference, vol. 42, p. 2549 (abstracts). Tran, T., Rosiek, M., Howington-Kraus, E., Archinal, B., Anderson, E., Team, L., 2010. Generating digital terrain models using LROC NAC images. In: Special Joint Symposium of ISPRS Technical Commission IV and AutoCarto, pp. 15–19. Williams, K., Zuber, M., 1998. Measurement and analysis of lunar basin depths from clementine altimetry. Icarus 131 (1), 107–122. Zuber, M. et al., 1992. The Mars observer laser altimeter investigation. J. Geophys. Res. 97 (E5), 7781–7797.