Range-Doppler Mapping of Space-Based Targets Using the JRO 50 MHz Radar

Abstract

The Jicamarca Radio Observatory (JRO) two-dimensional square array radar system operating at ~ 6-m wavelength was used to study the Moon and low Earth orbit satellites using the Range-Doppler inverse synthetic aperture radar technique also known as Delay-Doppler imaging. The radar data was collected on Oct 21, 2015. A circularly polarized coded pulse was transmitted from a quarter-array antenna segment during lunar transit over JRO. Dual-linear polarization receive systems were employed on two quarter-array segments and on two 1/64th array modules giving the longest possible baselines across the transit path. A Range-Doppler mapping technique that uses the rotational motion of the targets and an autofocusing motion and ionospheric delay compensation technique has been implemented to generate the two-dimensional maps of the point-target (Satellite) and range-spread target (Moon). A review of our technique and the maps obtained from these observations is presented herein. Range-Doppler maps of the Moon and satellites are instructive with regards to possible further improvement of the technique, especially regarding ionospheric compensation.

1 Introduction

The Jicamarca Radio Observatory (JRO) main radar antenna system is a fixed-pointing (near zenith), square-array used mainly as an Incoherent (as well as Coherent) Scatter Radar (ISR). It is located near Lima, Peru at the geomagnetic equator and it is used mainly to study the equatorial upper atmosphere and ionosphere. Usage examples include radar imaging of ionospheric irregularities and interferometric coherent radar imaging of meteors (Chau et al. 2008; Zhu et al. 2016; Malhotra et al. 2007). Delay-Doppler, also known as Range-Doppler, mapping of possible solar echoes was also attempted at JRO without much success as described in Coles et al. (2006). Here we report on an investigation into extensions to the imaging capabilities of this high power large aperture (HPLA) radar system. In particular, we report a new approach to imaging the Moon as augmented using ionospheric information gleaned from simultaneous LEO (Low Earth Orbit) and MEO (Medium Earth Orbit) satellite observations. In this task, we employed a modern signal processing technique—Range-Doppler Inverse Synthetic Aperture Radar (ISAR) mapping—to compensate for the lack of tracking capability and to partially correct for magnetoionic propagation effects encountered as both the Moon and the satellites transited through the wavevector-perpendicular-B (geomagnetic field) viewing region.

Prior to this work, JRO radar was used to determine the angular scattering law and the albedo of the Moon at the 6-m wavelength (Klemperer 1965; Hagfors et al. 1969). In Klemperer (1965) the antenna array configuration is described for detecting the scattered power from the subradar point to the lunar-limb in order to estimate the angular scattering power variations across the lunar surface. Whereas Hagfors et al. (1969) utilized the radiation and cosmic noise temperature fluctuations during the lunar transit over JRO to determine the dielectric properties of the lunar surface. In this paper, in addition to the scattered power variation as a function of “lunar depth”, we utilize the coherent phase information of the recorded, dual-orthogonal polarization voltages and develop an autofocusing-imaging technique to identify, map and detect the features on the lunar surface. In this we also employ a “zero-sidelobe” (truncated) IIR (Infinite Impulse Response) “un-matched” decoding filter (Lehtinen et al. 2004) to produce higher resolution images.

Additionally, while observing the Moon, LEO Satellite transits through the fixed JRO beam have been observed as described in Gao and Mathews (2015a). These ad hoc satellite observations were used in calibrating the phase and antenna transmission pattern. In this paper, we extend this calibration procedure in using the coherent phase information to image the satellite transits at meter wavelengths. Limitation on the radar receiver bandwidth (BW = 1 MHz) and hence on the imaging range resolution (150 m) makes these satellite traces an example for point-target imaging for JRO radar. Further, as the satellites are in readily characterized orbits and the radar observations are coherent across the entire transit, the actual orbit can, in principle, be resolved to about a wavelength. In classical terms, the single pulse resolution of a point target is set by the sampling rate given matching bandwidth whereas the point-target orbital resolution is established via the phase coherence of the radar pulse returns during the entire few seconds of the satellite overpass. Both ionospheric propagation effects and the signal-to-noise ratio of the satellite return signal limit the aforementioned coherence. To within the coherence issues just outlined, only one gravitational orbit is consistent with the observed satellite trajectory and can, in principle, be determined as part of the optimized inversion process. Thus, LEO satellites act as far-field specular targets that are useful in estimating and removing the strong ionospheric effects that would otherwise degrade the imaging. Imaging of a Mir space station satellite using a wideband (200–400 MHz) satellite-tracking radar system at VHF frequency (300 MHz) is described in Schmidt (2000). They reported tracking Mir using the published satellite ephemeris and employed it further to compensate for the dispersive effects of the ionosphere. In this paper, we demonstrate a similar sort of imaging algorithm that automatically compensates ionospheric delay as well as the polarization effects using the JRO narrow-band transit radar.

The JRO 50 MHz radar system is unique in that the very large, crossed-dipole array enables polarization measurements with relatively low cross-polarization coupling. In describing received signal polarization, we refer to the polarization direction of the echoes that matches the transmitted circular polarization wave as “depolarized” (same-sense or non-specular) and the echoes in the opposite sense of the transmitted wave polarization as “polarized”. The polarized echo is that expected for a specular or mirror-like reflection (terminology as defined in Campbell et al. 2007). Interestingly, most of our observations occur in the near k⊥B (radar wavevector direction perpendicular to the geomagnetic field) direction which complicates the interpretation of the received echo polarization due to the ionospheric manifestation of the Cotton–Mouton effect (CM; Yeh et al. 1999a). The measurements reported herein are, we believe, the first to report the CM effect on lunar and satellite echoes. We extend the past work at JRO to include the depolarized angular scattering law measurements that yield unique insight into the polarization and scattering characteristics of the lunar (sub) surface at ~ 6-m wavelength.

In the next section, we describe the experimental set-up that was used to observe the lunar and satellite transits during the same observation period. We then describe the imaging methodology employed for both the scenarios. Lastly, we present the resultant maps produced along with discussion of the ionospheric propagation effects and where to go next.

2 Experiment Set-Up

The JRO main antenna consists of a 288×288 m square array composed of 18,432 cross-polarized (COCO; coaxial-colinear) dipoles. The array is subdivided into four Quarters (North, South, West, East), each quarter consisting of 4 × 4 modules as shown in Fig. 1. Each module has “up” and “down” polarization receivers. The “up” polarization is approximately parallel to the northeast direction and is further denoted as Array Horizontal (AH). The “down” polarization is parallel to the southeast direction and is further denoted as Array Vertical (AV).

Fig. 1

Illustration of JRO main array receiver configuration. Channel A, B, E, F are module-array receivers with A, E receiving AH polarization and B, F receiving AV polarization. C, D, G, H are quarter-array receivers with C, G receiving AH polarization and B, F receiving AV polarization. Right-hand circular polarization was transmitted and the both linear AH and AV polarizations were received. Magnetic declination, the angle difference between the true north (N) and magnetic North, is 1.65°W (https://www.ngdc.noaa.gov/geomag-web/#declination)

The approximate orbital period of the Moon around the Earth is 27.3 days. Given the lunar orbital inclination to the ecliptic 5.14° and the Earth’s 23.44° spin axis tilt to the ecliptic, the lunar transit zenith angle at the equator varies from 28.58° North to 18.3° South and from 40.53°N to 6.35°S at JRO (antenna center at geographical 11.9516°S, 76.8743°W). Thus, lunar transits near beam center of the nearly zenith pointing JRO antenna can be observed approximately twice per month. The lunar transit over the JRO radar is along the east–west direction. Figure 2 shows the transit path of the Moon on 21-10-2015 overlaid on the quarter-antenna transmission radiation pattern. Figure 2 also indicates the locus of radar pointing perpendicular to the geomagnetic field (indicated elsewhere as k⊥B) at different altitudes. Our initial timing and “tracking” data was obtained using the elevation data of the Moon over JRO main radar antenna center from the NASA Horizons ephemeris data (http://ssd.jpl.nasa.gov/horizons.cgi). At the transit time of 19:04 h on 21-10-2015 the elevation value is ~ 88° and the expected duration of the lunar echoes is estimated to be ~ 25 min by taking the beam-edge elevation angle to be about between 85° and 88° range (approximate module antenna array radar beamwidth) from the ephemeris data.

Fig. 2

The lunar transit path (Orange line) through the JRO quarter-array radiation (transmit) pattern on the observation date of 21 October 2015 (Ref: http://jro.igp.gob.pe/programs/overjro/overjro.py). The entire ~ 20 min, east-to-west transit occurred near 1900 h. The thick black line represents the locus where the radar-pointing direction is perpendicular to the geomagnetic field direction at three altitudes. This set-up shows that lunar transit occurs near k⊥B direction resulting in Cotton–Mouton effect on the received lunar echoes. See Gao and Mathews (2015b) for details on JRO module, quarter-array, and full-array beam patterns

Both polarizations (AH- and AV-linear polarizations combined at − 90° to right-hand circular polarization; 600 kW peak power per polarization) of the East quarter-array antenna were used for transmission. This approach was used instead of the whole radar array as the quarter-array has a wider beamwidth (~ 2°) than the full antenna (~ 1.1°), hence yielding a longer observing window as transit occurred (see Table 1). We transmitted using circular polarization in order to mitigate ionospheric Faraday rotation effects (Davies 1990). This only partially succeeded as we describe later. We employed a 169-baud, nested Barker code at 10 μs per baud to achieve both a higher Signal-to-Noise Ratio (SNR) as well as an optimum sidelobe level. This code was generated using two-nested (Kronecker product) Barker 13-baud codes. Additional details about the nested Barker coding are given in Levanon 2005; Richards 2005; Vierinen 2012. Full radar parameters are given in Table 1. The InterPulse Period (IPP; alternatively, the pulse repetition interval) was adjusted to 39 ms and phased such that the received lunar echo of ~ 11.2 ms duration falls in between two transmission pulses. The lunar echo phasing was also adjusted to occur at a different apparent range from typical LEO satellite ranges. The maximum receiver bandwidth available at JRO is 1 MHz and is complex sampled at 1 MHz to achieve the best intrinsic range resolution. The receiver configuration of the JRO main radar is shown in Fig. 1. Quarter-array receiver beamwidth (~ 2.2°) is larger and comparable to the angular beamwidth of Moon i.e. ~ 0.5° and hence suitable to image the whole visible lunar surface with high SNR however with north–south hemisphere ambiguity (Thompson 1970). Therefore, the far-most North and South corner array-modules are also considered such that ambiguity can be removed using interferometry. The lunar interferometry results will be discussed in the future paper. The ~ 7º beamwidth module array elements aid in detecting high-inclination LEO satellites over wider transit paths.

Table 1 Radar parameters for the observations reported herein

3 Lunar Crossing Satellites at JRO

Lunar radar echoes were observed for approximately 20 min during the transit over JRO (~ 19:04 LT—Local Time). The measured lunar range is ambiguous to an integer number of IPPs, 39 ms. Since we know from the a priori ephemeris that the round-trip light time to the lunar sub-radar point is ~ 2.43 s, our measured apparent range of 19 ms (~ 2850 km) is in the 62nd IPP as shown in Fig. 3. A satellite transit of nearly 7 s (~ 19:06:31 LT) was also observed at ~ 515 km range in all the four receivers. After an online database search (www.calsky.com), this LEO satellite was identified as the Cosmos 1626. The AH-polarization range-time-intensity (RTI) plot from the North module receiver of the satellite echoes is shown in the bottom panel of Fig. 3. We also observed two more satellite echoes at approximately 18:51 LT and 18:53 LT and identified them as a DMSP and an Iridium satellite, respectively, by matching characteristics with the satellite information available at Calsky.com. The methodology described below is applied individually to both the Moon and Cosmos 1626, using the respective signals from the North- and South-module channels.

Fig. 3

The top panel is the range time intensity (RTI) plot of lunar echoes over the entire transit. The total JRO-lunar round-trip delay is ~ 2.43 s. As the InterPulse Period (IPP) is 39 ms, the lunar echo that reaches the receiver after the delay is seen in 62nd IPP from the transmit time and positioned at ~ 3000 km (~ 19 ms delay) range. The bottom panel is the RTI plot of the lunar crossing Satellite at JRO. Matched filtering is used to decode the data shown in these plots with the resultant decoding side lobes clearly visible. This decoding approach is later contrasted with the (nearly) zero-sidelobe alternative decoding approach

4 Methods

The processing steps include preprocessing polarization decomposition, parameter estimation for autofocusing followed by the Range-Doppler processing steps that include range alignment, phase adjustment and azimuth FFT (Fast Fourier Transform) to generate the Range-Doppler image. The Range-Doppler algorithm steps are similar to the standard motion compensation path for generating Range-Doppler maps of targets [see Fig. 3.6 in Chen and Martorella (2014)].

4.1 Preprocessing

A necessary signal calibration procedure, the preprocessing step, is performed on the lunar return signal to obtain the depolarized target vector (further defined below). In standard lunar radar mapping—e.g., datasets collected while transmitting from Arecibo and receiving at Green Bank observatories at 70 cm wavelength—the orthogonal receiver polarizations were used directly to get the polarized and depolarized maps (Campbell et al. 2007). As defined in Sect. 1, in the usual case, the depolarized vector refers to the polarization sense of radar scattering from the lunar surface (range-spread target) that is in the same sense as that of the transmitted signal. Whereas, in our case, even though we collect the data from orthogonal linear polarizations, the presence of geomagnetic field (around 1° dip angle) at the JRO location, causes the polarization ellipse of the signal to Faraday rotate as the wave propagates through the ionosphere. Faraday rotation as such does not affect circular polarization however it is especially important at VHF in the presence of a “turbulent” ionosphere and much more so when the wavevector is very close to perpendicular to B due to the Cotton–Mouton effect. This variation in polarization is highly sensitive to the aspect angle between the direction of wave propagation vector (k) and the geomagnetic field vector (Segre 1999; Yeh et al. 1999a, b).

In the case of the JRO lunar observations reported here, the wave vector is nearly perpendicular to the geomagnetic field throughout the transit as shown in Fig. 2. This results in the conversion of the transmitted circular polarization to near linear polarization in a propagation process known as the Cotton–Mouton effect which is applicable for small aspect angles (< 0.5°) (Milla 2010). Basically, the Cotton–Mouton effect stems from the (decomposed) parallel-B linear polarization wave E-field component not “seeing” the B-field while the perpendicular-B E-field component is Faraday rotated until it is parallel-B. This magnetoionic process of course occurs on both the transmitted and scattered waves. Thus, the orientation angle and shape of the polarization ellipse are dependent on magnetic field strength, radar wavelength, and the total ionospheric electron content on the ray path as well as on the original transmitted polarization and the polarization signature of the scattering process. For detailed mathematical formulation that shows the dependence of the parameters, see Eqs. (6.6)–(6.8) given in Segre (1999) and theoretical explanation of polarization change at various aspect angles in Yeh et al. (1999b).

Additionally, as the ionosphere at the geomagnetic equator is electrodynamically very “turbulent”, ionospheric effects vary throughout the observation period. In addressing this issue as well as the Cotton–Mouton effect, we developed an optimization/calibration procedure that was used to derive the best possible depolarized vector from the received data itself. In this procedure, the received polarized vectors (AH, AV) are decomposed to estimate the axes of the received signal polarization ellipse that best represent depolarized and polarized component vectors. This decomposition is performed such that the difference of the power values of the matched filtered decoded leading-edge of the target is maximized. The axis that has minimum power value is defined as the depolarized vector and the maximum power value as the polarized vector.

The satellite we observed transited from north-to-south across the beam (main lobe) within a short time (~ 7 s) compared to the lunar transit (~ 20 min). The projected satellite track across the antenna array is shown in Fig. 4. This transit path is approximately determined using the projected latitude/longitude “trajectory” of the satellite derived from the two-line orbit elements corresponding to the observing time as available at the Calsky.com website. The propagation path aspect angle with respect to the geomagnetic field varies for each IPP as the target passes through the beam as shown in Fig. 4. In this case, polarization shifts from circular to linear as the propagating wave direction nears the perpendicular-to-B direction and then back to circular (Milla 2010). This can be observed as the low-to-high-to-low variation in axial ratio values of the polarization ellipse as shown in Fig. 5. The change in signal polarization throughout the transit also adds additional phase variation to each IPP in the received orthogonal components resulting in less coherent phase information. To compensate for the effect a common polarization vector with respect to all IPPs received was created using the same optimization technique described above. As the transit orbits of the other two satellite events we observed did not provide suitable aspect angle variation with respect to the geomagnetic field, they are not discussed here. The derived depolarized vector is used to estimate the autofocus parameters in the next step.

Fig. 4

Cosmos 1626 satellite transit path and time over JRO radar as derived from the two-line orbit elements (TLE) available from the Calsky.com online website source. The timing accuracy is ~ 230 ms. This is obtained by comparing the timing of the minimum range in the radar data (19 h 06 m 31.77 s) with that derived from the TLE data (19 h 06 m 32 s). The aspect angle variation (angle between receiver position located at axes origin and B-field at the given latitude–longitude) during the transit path is shown in the top label. Axes origin is the JRO main radar center coordinates. The B-field estimates are obtained from: http://jro.igp.gob.pe/programs/overjro/overjro.py

Fig. 5

The top panel gives the Left- and Right- circular polarization response (LCP, RCP) of the leading-edge echo (sidelobe-free decoding) of the satellite event observed in the north- and south-module receivers. As the transmitted signal was RCP, the (expected) specular return is LCP. The bottom panel is the axial ratio (AR = LCP/RCP) observed in both the north- and south-module receivers. The maximum axial ratio spike at ~ 13 s is observed when the aspect angle is nearly 90° (shown in Fig. 4). This represents the Cotton–Mouton effect region which results in largely linear polarization as discussed in the text. Note that a small AR represents circular polarization. In the ideal case of linear polarization, AR is infinite

4.2 Autofocus Parameter Estimation Procedure

In the Range-Doppler mapping technique, the net translational motion of the target results in the blurring if phase (Doppler) distortion is not corrected (Chen and Martorella 2014). It is the rotational motion of the target that is exploited to get the Doppler dimension of the map. In conventional processing of lunar radar data with mapping as the goal, the relative motion between the radar and the target is known from ephemeris information and is used to remove the smearing/blurring effect to develop the focused maps. Further, Doppler offset has in the past been hardware compensated in real-time (Pettengill et al. 1974; Campbell et al. 2007; Margot et al. 2000). As we observe the Moon only as it transits the radar, the range-walk and range-rate phase (Doppler) corrections are estimated from the received data itself and compensated (aligned and adjusted) accordingly. This has been referred to as the autofocusing technique (Chen and Martorella 2014). Note that this approach automatically compensates for the delay due to the ionosphere that is significant and potentially quite variable especially at the JRO VHF observing wavelength (~ 6 m) where magnetoionic effects are very apparent.

Different translational motion compensation techniques such as Prominent Point Processing (PPP), Phase Gradient Algorithm (PGA), Entropy Minimization Algorithm, and Maximum Contrast based compensation techniques are described in the Range-Doppler Inverse Synthetic Aperture Radar literature (ISAR) (e.g., Chen and Martorella 2014; Chen and Andrews 1980; Berizzi et al. 2004; Küçükkiliç 2006). The signal processing technique based on the Prominent Point Detection technique shown in Fig. 6a is implemented to determine the range/delay and range-rate/Doppler of the target as a function of (transmitted pulse) time. An issue is that even a small outlier result in the range detection of a prominent point—in this case the lunar subradar point—results in a small change in Doppler offset which highly blurs the phase-sensitive resultant image (ill-conditioned system). This is because lunar limb-to-limb Doppler bandwidth is very small—1.25 Hz at 6-m wavelength. Therefore, selected IPP’s are chosen to derive the Doppler fit parameters that result in best resultant map based on visual inspection. As more than one receiver is available, data from North- and South- module receivers are combined to determine the range fit and Doppler shift variation and to establish common focusing parameters for both maps.

Fig. 6

The top panel is the flow chart for the autofocus processing used to estimate the sub-radar “point” Range and Doppler parameters. The bottom panel is the Range-Doppler image formation signal-processing path

4.3 Range-Doppler Processing

The Range-Doppler processing flow chart is shown in Fig. 6b. Initially, the received raw data from both linear polarizations are rotated by an angle estimated using the optimized technique described in Sect. 4.1. In the next step, range compression is performed by using the principle of sidelobe free decoding filter described in Lehtinen et al. (2004). It is possible to find the unique sidelobe-free (IIR: Infinite Impulse Response) decoding filter provided the modulation function of the code has no zeros in the frequency domain. However, the impulse response of such a side lobe-free filter is necessarily truncated. In our case, a filter length of 28.5 ms that yields the maximum SNR value with minimal sidelobes was chosen. The transmitted code, decoding filter and the impulse response of this filter are shown in Fig. 7.

Fig. 7

The 169-bit code (two nested, Barker 13 codes) with 10 μs baud length transmitted code (top). The sidelobe-free IIR compression filter is shown in the middle panel. The IIR truncated filter length is 28.5 ms. Bottom left is the impulse response with the matched filter and to the right with the sidelobe free compression filter. The principle is discussed in Lehtinen et al. (2004). An annotated Mathematica notebook giving these signal processing algorithms are available on request

The autofocus range parameter corresponding to the leading edge (subradar point) of the target are interpolated at an order of 3 and are used to align the ranges and perform the coarse range alignment step such that the leading edge from each IPP occurs at the same range gate. The phase adjustment step described in Kesaraju et al. (2016) is implemented next using the estimated autofocus Doppler parameter values interpolated at an order of 3 (in standard processing, Doppler parameters are obtained from ephemeris data). Once the entire decoded dataset is compensated for the translational motion, the Doppler processing stage involves Fourier transforming the entire voltage sequence at each range gate yielding an estimate of the Doppler power spectrum at each range. This step is also known as the azimuth (Fast) Fourier transform. Lastly for the range spread target, the consecutive- IPP, power spectrum values in all range gates are incoherently averaged to get the desired resolution Range-Doppler map. Details surrounding these steps are discussed below.

5 Results and Discussion

5.1 Autofocus Parameter Estimation

The leading edge of lunar echoes observed for approximately 500 s (~ 12,800 individual radar pulse returns) was used to estimate the autofocus parameters. Figure 8 shows the quadratic fit to the range of the leading-edge point derived from both the North and South module depolarized channel data. Selected IPP’s from both receiver channels are combined such that common focusing parameters—as defined in Eq. (1)—is established. As per Eq. (1), the linear term is the sub-microsecond alignment needed for the proper range and phase alignment of all radar “pulses” reflected from the Moon. Subsequently, the differentiated linear equation of the quadratic fit in Eq. (1) results in the estimated Doppler shift of the lunar leading-edge point as given by Eq. (2). Similarly, Fig. 9 shows the range fit and Eqs. (3), (4) give range and Doppler fit for the satellite event. Note, that for an IPP of 39 ms, the maximum unambiguous velocity is 38 m/s that is of order of the projected (near zenith) orbital velocity (~ 7 km/s) of satellites in the LEO orbit.

Fig. 8

A quadratic fit to the range values of the detected depolarized (non-specular) lunar leading edge point. This is estimated using the signal processing path outlined in Fig. 6a. The range fit equation is derived using the combined data from both module receivers as given by Eq. (1). The extreme outlier points are not included in the least-squares quadratic fit to this data. The delay range outlier points seem to be due to ionospheric effects

Fig. 9

A quadratic fit to the range values of the depolarized (non-specular) satellite leading edge point. This is estimated using the signal processing mentioned in Fig. 6a from North and South module receivers. The range fit equation derived by combining the data from both receivers shown in Eq. (3)

Lunar parameters:

$$\begin{array}{*{20}c} {r_{m} \left( {\text{megameters}} \right) = 365.439 - 69.039 \times 10^{ - 6} t + 1.57946 \times 10^{ - 8} t^{2} } \\ \end{array}$$

(1)

$$\begin{array}{*{20}c} {f_{dm} \left( {\text{Hz}} \right) = - 22.9921 + 0.0105202 t} \\ \end{array}$$

(2)

Satellite parameters:

$$\begin{array}{*{20}c} {r_{s} \left( {\text{megameters}} \right) = 0.523214 - 0.0012111 t + 0.0000498118 t^{2} } \\ \end{array}$$

(3)

$$\begin{array}{*{20}c} {f_{ds} \left( {\text{Hz}} \right) = - 413.224 + 33.9879 t} \\ \end{array}$$

(4)

As seen in Fig. 8, possible reasons for the outliers in the range estimation of the leading edge could be due to rapidly changing specularity because of the lunar libration and/or due to the influence of the ionospheric path (group delay) scintillation. The estimated range values of the lunar leading-edge point (Fig. 8) are compared with the range values generated using the ephemeris data provided by the NASA Horizons system service and is available online at http://ssd.jpl.nasa.gov/horizons.cgi. Figure 10 shows the offset between both the values throughout the observation period. The main reason for such an offset is the significant two-way group delay of the VHF signals during the ionospheric propagation. The approximate (note that this ignores Faraday rotation effects) one-way group delay in meters ($\Delta \tau$) due to the presence of ionosphere is estimated to be (Giffard 1999; Schmidt 2000)

$$\begin{array}{*{20}c} {\Delta \tau = 40.3\frac{{N_{T} }}{{f^{2} }} } \\ \end{array}$$

(5)

where $N_{T}$ (TEC units of 10¹⁶ electrons-m⁻²) is the total free electron content in the ionosphere and f is the frequency of the propagated signal. From Fig. 10, root-mean-square (RMS) error of the ephemeris delay and estimated delay from the observation, i.e., the one-way group delay ($\Delta \tau$), is calculated to be 38.9 μs yielding an apparent 1.85 TECU ionosphere electron content. Here, the range derived from the ephemeris has ~ 1 m accuracy (ftp://ssd.jpl.nasa.gov/pub/ssd/Horizons_doc.pdf) that is equivalent to sub-microsecond accuracy level for the estimated TEC content. However, the mean radar range to the lunar subradar point is less accurate and difficult to quantify per Fig. 8.

Fig. 10

Estimated lunar sub-radar point range (megameters) derived from both the north and south module depolarized vector using autofocus processing is represented as black fit. This is compared with the sub-radar point range provided by the NASA horizons ephemeris system (blue line). The range offset is due to the VHF signal propagation through the ionosphere. The resultant difference corresponds to total electron content of 1.85 TEC units (10¹⁶ m⁻²). The delay outliers among the individual data points are also assumed to be due to propagation effects from ionospheric “blobs” and scintillations

Further, the estimated one-way group delay/range from the lunar radar echoes is subtracted from the observed range of the satellite to estimate the true range. In our case, to determine the satellite orbital speed from the radar data, the angular trajectory of the satellite within the radar beam is estimated using the technique described in Gao and Mathews (2015b) considering the two perpendicular baselines available in this receiver configuration—baseline-1-North–South module pair and the baseline 2-East–West Quarter. This results in an orbital speed of 7.84 km/s. Whereas, gravitational orbital speed equation [see Eq. (3) in Gao and Mathews (2015a)] using the true range of the satellite results in ~ 7.9 km/s orbital speed. These results confirm the calibration of the radar, as well as the orbit of the satellite, even in the presence of ionospheric propagation effects that distort the phase and delay terms of the received radar echoes.

Additionally, vertical TEC data from the worldwide GPS receiver network is maintained online in the MIT Haystack Madrigal distribution database system (Rideout and Coster 2006). Even though the exact TEC values from this network were not available at the radar location during the observing period, we can infer that the average TEC from the GPS receiver network relatively nearby JRO during the whole observation period to be 1.4 TECU. The small-scale ionospheric irregularities affect the ionospheric delay and thus blur the resultant Range-Doppler map although the autofocusing step described above serves to mitigate this issue.

5.2 Range-Doppler Maps

Satellite events, due to the intrinsic resolution of 1.5 km (10 μs baud-length), are seen as the point targets in the range dimension. As pointed out earlier, the combined two-body orbit and coherent radar transit observations yield a net orbital accuracy of the order of a wavelength in range with sufficient observing bandwidth and power. Doppler resolution (ΔD) and cross-range resolution (Δr) is dependent on the net observing time/coherence time (T_obs), the rotation rate of the target (ω) and wavelength (λ) as shown in Eq. (6) (Margot et al. 2000)

$$\begin{array}{*{20}c} {\Delta D = \frac{1}{{T_{obs} }} = \frac{2 \omega \Delta r}{\lambda } } \\ \end{array}$$

(6)

From the Calsky online website source, the angular velocity (rate of axial rotation) of the observed satellite is known to be 0.86°/s which gives the cross-range resolution at the transmitted VHF frequency to be approximately 20 m for the observing period of 10 s. As the dimensions of the satellite are ~ 6 m × 2 m the resultant map will be single pixel image as shown in Fig. 11. The left panel of Fig. 11 shows the depolarized Range-Doppler map when matched filtering is used for range compression and the right panel of Fig. 11 shows the result obtained when sidelobe free filter (~ IIR) is used to map the point target. Therefore, the lunar-crossing satellite radar data validates the removal of ringing effect in the Range-Doppler maps along the range dimension when the sidelobe-free filter is used. Per Fig. 5b, the satellite echoes also confirm the Cotton–Mouton effect thus validating our approach to the much more complex lunar polarization effects.

Fig. 11

Example Range-Doppler map of the Cosmos 1626 lunar crossing satellite. The left and right panels show the same satellite echoes decoded using matched and IIR filters, respectively. The zero-range-sidelobe IIR filter is clearly necessary for both the satellite and lunar mapping observations

For the lunar maps, a net observing period of 6000 IPP’s equivalent to T _obs  = 234 s are considered. The observed signal coherence time needed to achieve the Doppler resolution required is apparently limited by the frequency stability of the clock used by the JRO radar with ionospheric variations likely a contributing issue given the Fig. 8 results. The Trimble clock (Thunderbolt E GPS disciplined double ovenized clock) used at JRO has frequency stability (Allan deviation) of order 10⁻¹¹ (Trimble 2011). The average coherency period is estimated to be 250 s at the 50 MHz frequency and is calculated using the Eq. (9) of Vierinen and Lehtinen (2009). This net integration time (T _obs ) of 234 s yields a frequency resolution of ~ 0.004 Hz which corresponds to cross range resolution of ~ 11 km. This considers the limb-to-limb (2R _m ) lunar Doppler bandwidth to be ~ 1.25 Hz (at 50 MHz) which then corresponds to Doppler spatial resolution ~ 0.36 mHz/km where R _m is 1738 km lunar radius and Doppler bandwidth given by 4ω R _m /λ. Therefore, accuracy of Doppler fit parameters (Eq. 2) should be ~ 3 mHz to achieve 11 km cross range resolution without smearing effect.

The North module receiver Range-Doppler-Intensity (RDI) map—where Range is the range beyond the subradar point—of depolarized lunar surface scattered power is given in Fig. 12. This map is the result of averaging in the range dimension to arrive at resolutions of ~ 10 × 11 km along the range and cross range directions, respectively. Features labeled in Fig. 12 were found by comparison with the lunar optical map from the Lunar Reconnaissance Orbiter Camera (LROC) that is available online at http://webmap.lroc.asu.edu/lunaserv.html.

Fig. 12

Lunar Range-Doppler-Intensity (RDI) map (includes north–south ambiguity) of normalized (to subradar point power) depolarized (non-specular) lunar surface scattered power. This map is derived from 6000 radar pulse returns (over ~ 234 s). The signal processing technique is described in the text. The labeled features are identified to be Pythagoras (A), Carpenter (B), Philolais (C) and Anaxagoras (D), Copernicus (E) craters and Mare Nectaris (F) on the lunar surface

Even though the scattering properties of the individual lunar features cannot yet be used for the further analysis from the north–south ambiguous map, the normalized angular scattering law observation can be derived by averaging the power values at each range/delay from all the observed data. Figure 13 shows the normalized echo power (dB) versus the lunar delay depth from both the orthogonal polarizations—i.e., depolarized and polarized obtained after the preprocessing step. The maximum return power observed at the initial delay depth is due to the specular (mirror-like) reflections from the sub-radar/leading-edge region. As the 6-m wavelength radar wave penetrates deep into the regolith, depolarized (non-specular) scattering from the sub-radar region is relatively small compared with polarized (specular) reflections as seen in Fig. 13. Nonetheless, polarized scattering dominates depolarized scattering at all lunar delay depths. While keeping in mind that the JRO array polarization separation is likely no more than ~ 25 dB, this result points to the scattering complexity of the lunar surface/regolith structure. Figure 4 of Pettengill and Henry (1962) shows a similar angular (lunar depth) scattering variation of the lunar echo observed at 440 MHz. The polarized (opposite-sense) signal at 440 MHz also dominates the depolarized signal at all the delay depths.

Fig. 13

Normalized lunar power echo versus delay obtained from all the received echoes for two optimally estimated orthogonal polarizations at 50 MHz. The depolarized echo displays lower net power than the polarized echo with respect to delay similar to the results obtained at 440 MHz described in Pettengill and Henry (1962)

To further validate that the sub-radar region acts as a specular (mirror-like) reflector, the polarization ratio (Eq. 9, Pettengill and Henry 1962)

$$\begin{array}{*{20}c} {P = \frac{{I_{P} - I_{DP} }}{{I_{P} + I_{DP} }}} \\ \end{array} .$$

(7)

for both the North and South modules is given in Fig. 14. Here, I _P represents the polarized and I _DP represents the depolarized vector. It can be inferred that the scattering surface tends to have diffuse scattering characteristics at greater delays. This is likely the result of the net scattering from many scattering centers that include meter-scale boulders in the regolith. Similar conclusions were inferred from a radar polarization study at 440 MHz (see Fig. 5 in Pettengill and Henry 1962). However, they find that the minimum polarization value at the largest delay (~ 11 ms) is greater than 0.4 which is larger than the value (~ 0.34) observed in our data. A similar set of measurements and polarization ratio response of lunar power from the orthogonal linear polarization for circular polarized transmitted signal at 23 cm and 3.8 cm wavelength is discussed in Hagfors et al. (1965) and Hagfors (1967). These measurements were obtained at 1.6 ms delay resolution and only up to one tenth of the diameter of the lunar disk (~ 1 ms delay depth). Therefore, these values were not used for comparison analysis here. We can attribute increase in depolarization at 6 m wavelength to the subsurface multiple scattering features observed due to the greater regolith penetration depth at 50 MHz. However, it might be premature to conclude that there are fundamental differences in scattering properties at this point. This last conclusion is cautious as ionospheric (magneto-ionic) effects are quite important and even dominant in the Cotton–Mouton region as discussed in Sect. 4.1. Additionally, fading of the received echoes—i.e., variations due to lunar librational motion (Evans et al. 1959) as well as ionospheric “turbulence” is in the range of a few seconds or faster—certainly effects these results. Therefore, future observations will be directed towards better understanding of ionospheric propagation effects to achieve true orthogonal depolarized and polarized components and validate the results presented here.

Fig. 14

Polarization ratio versus lunar delay obtained from all the received echoes at 50 MHz. Polarization ratio is defined as per Eq. (7) and lies in the range of 0–1. A polarization ratio near 1 is expected for mirror-like returns

6 Conclusion

In this paper, we report on a comprehensive analysis of JRO 50 MHz radar echoes observed from the lunar surface as well as from a lunar-crossing satellite. We demonstrate suitable observational and processing techniques for mapping LEO satellites. We especially consider ionospheric propagation effects on signal polarization as the magnetoionic Cotton–Mouton effect is a dominant feature when the radar wave vector (k) is very nearly perpendicular to B. We further present an extension of this technique to range-spread targets in the Delay (Range)-Doppler domain. Previously, as described in Davis and Rohlfs (1964) and Davis et al. (1965), many attempts were made to obtain lunar and planetary radar echoes at long wavelength (in the meter range) resulting, however, in a challenge to mitigate ionosphere effects. Range-Doppler maps at longer wavelength also require longer coherence time in obtaining the Doppler dimension of the maps. In our case, in addition to the challenges imposed by the longer wavelength, as the JRO radar does not have tracking capability, we had to develop and apply the offline focus methodologies unlike the conventional processing techniques (when tracking the lunar sub-radar point or other features using ephemeris data) used to obtain the lunar maps from Arecibo, Haystack, GreenBank (bi-static, receive-only) radar observatories. We also derived and removed the lunar range-rate Doppler (full offset), a function that was accomplished in hardware or using ephemeris data previously. Secondly, the location of JRO near the equator results in smaller aspect angles with the Earth’s magnetic field that effects the polarization of the received radar echoes. However, JRO is the largest high power large aperture radar (HPLA) radar in the world (as of present) operating at or near 50 MHz frequency that can provide the required SNR to detect and map the lunar surface. Thus, offering observational insights not available from other radars. A parallel calibration technique was developed using LEO satellites detected during the lunar observations.

The feasibility of obtaining scientifically relevant lunar radar (scattering) maps at JRO even with a small lunar limb-to-limb Doppler Bandwidth of ~ 1.25 Hz, no transit tracking capability and the additional magneto-ionospheric propagation effects is discussed with the context of modern Range-Doppler autofocus signal processing techniques and the observational setup introduced in this paper. The satellite result presented in this paper (refer to Sect. 4.1 and Fig. 5) provides the first scientific Cotton–Mouton effect results obtained via radar observations of satellites. Moreover, variation of TEC results derived from the lunar echoes is shown. Lastly, the lunar radar map we derive show that the depolarization (which can be attributed to diffuse scattering) is stronger at greater delays than compared to that observed at 440 MHz by Pettengill and Henry (1962). We tentatively attribute this to the greater penetration of the 6-m wavelength radar signal into the regolith with sub-surface specular scattering away from the radar resulting in a net increase in relative depolarization of the received signal.

Future efforts will concentrate on module receiver interferometry to remove the north–south ambiguity while maintaining maximum polarization separation. Modern super-resolution spectral estimation techniques like Capon’s spectral algorithm (Capon 1983) and amplitude and phase estimation (APES) approach (Li and Stoica 1996) will also be explored. Perhaps most importantly, further lunar observations at different lunar phases and that avoid the CM region will allow for better estimates of the circular polarization ratio for individual features of the Range-Doppler map. Lastly, the planned addition of a modest tracking capability to the JRO radar will allow for substantially better lunar imaging.