Production of Reduced Data Records for the Phoenix Atmospheric Structure Experiment Paul Withers (a) and D. C. Catling (b) (a) Center for Space Physics, Boston University, 725 Commonwealth Avenue, Boston, MA 02215, USA. (b) Department of Earth and Space Sciences, University of Washington, Box 351310, Seattle, WA 98195-1310, USA. Abstract The purpose of this report is to describe the methodology used to produce Reduced Data Records (RDRs) for the Phoenix Atmospheric Structure Experiment (ASE) from its Experimental Data Records (EDRs). These RDRs include vertical profiles of atmospheric density, pressure, and temperature. During its brief flight through the atmosphere of Mars, Phoenix recorded acceleration and angular velocity data using accelerometers and gyroscopes within an inertial measurement unit. These time series data, which constitute the EDRs, are archived at the NASA Planetary Data System (PDS). The archived EDRs require correction because of errors associated with the use of a preliminary version of a transformation matrix between two reference frames. The EDRs, which were acquired at a high rate of 200 Hz, are noisy, but this can be improved by smoothing at the cost of worsened vertical resolution. Smoothing of the axial acceleration requires care. We have described how the running arithmetic mean of the axial acceleration is a biased estimate of the true acceleration and explained how this bias can be corrected. The accelerations and angular velocities (EDRs) have been used to reconstruct the trajectory of Phoenix through the atmosphere and the associated atmospheric structure (RDRs). The position and velocity of Phoenix were reconstructed by forward integration of the equations of motion and its orientation history was determined directly using onboard angular velocity measurements. The density profile along the trajectory was reconstructed from about 130 km to parachute deployment, about 13.5 km above the surface, using the drag equation. The corresponding pressure profile was obtained from the equation of hydrostatic equilibrium and the corresponding temperature profile was obtained from the ideal gas law. The angle of attack of Phoenix during atmospheric flight was determined directly using angular velocity measurements and indirectly using the ratio of normal to axial accelerations and the spacecraft's aerodynamic database. The former was typically 1-2 degrees greater than the latter, a discrepancy for which we do not have a good explanation. 1. Introduction Data from the entry, descent, and landing (EDL) of the Phoenix spacecraft, have been used to obtain a profile of martian atmospheric density, pressure, and temperature from about 130 km to about 13.5 km above the surface. The thermal structure of the martian atmosphere is sensitive to radiative forcing from suspended dust and to diabatic heating associated with atmospheric dynamics [ZUREKETAL1992, LEOVY2001]. It is also perturbed by a wide variety of waves and tides [LEOVY&ZUREK1979, BANFIELDETAL2000, WITHERSETAL2003A]. This is the first such profile of martian atmospheric structure from the polar regions. The atmospheric processes that can be observed in such profiles were discussed by [MAGALHAESETAL1999], who also compared the advantages and disadvantages of this measurement technique to those of other techniques. NASA's Mars Scout program selected the Phoenix mission for flight in 2003. Its "science mission focuses on providing the ground truth for the 2002 Odyssey discovery of massive ice deposits hidden under surface soils in the circumpolar regions" [SMITHETAL2008]. The objectives of this mission were: (1) to study the history of the ground-ice and its emplacement mechanisms, (2) to address the effect that subsurface ice has on the local surface geomorphology, (3) to characterize the climate and local weather of the landing site, and (4) to address the habitability of the icy soil [SMITHETAL2008, ARVIDSONETAL2009, SMITHETAL2009]. Much of the design of the Phoenix spacecraft, including those aspects relevant for cruise and EDL, derived from the mothballed "Mars Surveyor 2001 Lander" [GUINNETAL2008, DESAIETAL2008, GROVERETAL2008]. Phoenix launched on 4 August 2007 and successfully landed on Mars in the Vastitas Borealis or northern plains region on 25 May 2008. The time and position of landing are reported in Table 1. [TABLE 1] Phoenix carried two Honeywell inertial measurement units (IMUs, model number YG9666BC) [TAYLORETAL2008]. IMU-A recorded data during EDL; IMU-B, flown as a backup, was never turned on after launch and is not discussed further. Each IMU consisted of three identical single-axis Allied Signal Q-Flex accelerometers and three identical Honeywell GG1320 RLG ring-laser gyroscopes. These were mounted orthogonally. Accelerometers measure acceleration and gyroscopes measure angular velocity, but the IMU was configured to output time-integrated quantities, specifically velocity change and angle change, every 5 ms. These measurements, which constitute the Experimental Data Records (EDRs) of the Atmospheric Structure Experiment (ASE), were archived by the NASA Planetary Data System (PDS) with associated documentation [CATLINGETAL2008]. [DESAIETAL2008] and [GROVERETAL2008] also provide valuable background information on the IMU, spacecraft and EDL sequence. Analyses of the IMU data for engineering purposes have been undertaken [DESAIETAL2008, BLANCHARD2009]. The purpose of this report is to describe and discuss a reconstruction of the Phoenix EDL trajectory and the associated atmospheric structure that has produced a set of Reduced Data Records (RDRs) for the Phoenix ASE investigation. Section 2 introduces the principles of trajectory reconstruction, Section 3 describes how the Phoenix EDR data were processed, Section 4 describes how the attitude of Phoenix was determined, Section 5 describes the state of Phoenix at atmospheric entry, Section 6 describes the trajectory reconstruction for Phoenix, Section 7 describes the atmospheric structure reconstruction for Phoenix, Section 8 describes some consistency checks on this reconstruction, Section 9 describes possible future work, and Section 10 summarizes the main conclusions of this work. Appendix A describes how errors in the EDR data were corrected and Appendix B describes how noise in the axial accelerations was reduced by smoothing. 2. Introduction to Trajectory Reconstruction The trajectory reconstruction process used for Phoenix was similar to those applied to earlier missions [WITHERS&SMITH2006A, WITHERSETAL2003B] and references therein. First, clean time series of angular velocities and accelerations at the centre of mass in a suitable frame were produced. Second, spacecraft orientation, position and velocity at the time of atmospheric entry were obtained. Third, the equations of motion were integrated forward in time from entry using the initial conditions, angular velocities, accelerations and spatially-varying gravitational field to determine spacecraft orientation, position and velocity as functions of time. The existence of reliable angular velocity data from the Phoenix gyroscopes enabled direct determination of spacecraft orientation, unlike the indirect methods used on some earlier spacecraft [SPENCERETAL1999]. 3. Preparation of Data Characteristics of the IMU and its measurements are summarized in file CATALOG/INST.CAT of [CATLINGETAL2008]. The data used in this report come from the archived EDRs; specifically, the measured velocity increments and angle increments reported in file DATA/IMU_A_EDR_M.TAB of [CATLINGETAL2008]. Zeroes were reported for all IMU data at samples 12031 to 12040, a duration of 50 ms, due to a downlink transmission problem. This data gap was filled by interpolation in this report. These data were converted into accelerations and angular velocities by dividing by the relevant time interval, 5 ms. Many different coordinate frames were associated with Phoenix, as discussed in file DOCUMENT/EDRSIS.PDF of [CATLINGETAL2008], [BLANCHARD2009], and the Phoenix frame kernel distributed by NAIF (JPL's Navigation and Ancillary Information Facility) in SPICE format (ftp://naif.jpl.nasa.gov/pub/naif/PHOENIX/kernels/fk/phx_v06.tf). The archived data were transformed into the desired frame, the "C frame", as described in Appendix A. This frame, also called the "Phoenix cruise frame", is defined in Appendix A. Appendix A also discusses a small, but significant, correction to the EDR accelerations and angular velocities. The symmetry axis of the spacecraft, along which the largest accelerations were experienced, is the +/- X-axis in the "C frame". The 5 ms interval between successive raw data points corresponds at high elevations to a change in elevation of 6 m, far smaller than any relevant atmospheric lengthscale. Accordingly, noise in the data was reduced by averaging. Accelerations along the Y and Z-axes of the "C frame" were smoothed with a 256-point running mean. These axes are normal to the spacecraft's symmetry axis. The same uncertainty, 0.0025 m s-2, was assigned to each of these smoothed data points. This value corresponds to the standard deviation of the smoothed normal accelerations at high elevation. All uses of uncertainties described in this report were based on Monte Carlo techniques assuming normal distributions. Accelerations along the X-axis of the "C frame", the spacecraft's symmetry axis, were smoothed using a more complicated algorithm. The guiding principles behind this algorithm were (A) the need for balance between reduction in noise and acceptable vertical resolution and (B) consideration of the appropriateness of an unweighted arithmetic mean. In theory, the axial accelerations at high elevations should increase exponentially with time (e.g. Figure 1 of [MAGALHAESETAL1999]) due to the descent at constant speed into an exponentially-varying atmosphere that occurs prior to peak deceleration. In such circumstances, the arithmetic mean of a subset of successive accelerations is not the best possible estimate of the true acceleration at the centre of this interval. The axial accelerations were smoothed as described in Appendix B. The standard deviation of axial accelerations recorded at high elevation, prior to the onset of detectable aerodynamic drag, is a measure of the uncertainty or error in the acceleration data. This standard deviation depends on the size of the smoothing window. Analysis of these data (not shown) indicated that the dependence of the standard deviation, sigma_aa, on the number of points in the smoothing window, N, satisfies: Equation 1 sigma_aa = a0 / N where a0 was found to equal 0.3 m s-2. Although this value of a0 is not derived directly from the instrumental digitization, it is similar to the ratio of the digitization limit of 2.7 mm s-1 to the 5 ms interval between successive data points (0.54 m s-2). Note that the 1/N dependence in Equation 1 differs from the 1/(the square root of N) trend expected for values sampled from a normal distribution. The uncertainties in the smoothed axial accelerations were assumed to equal the relevant standard deviation, so the uncertainty in any given measurement is determined by the smoothing window used. Variations in the width of the window used to smooth the axial accelerations, which are outlined in Appendix B, are shown in Figure 2. [FIGURE 1] [FIGURE 2] At high elevations, the unsmoothed angular velocities had a standard deviation of 10^(-3) rad s-1. This was adopted as the uncertainty of the unsmoothed angular velocities. Angular velocities about the X, Y and Z-axes of the "C frame" were then smoothed with a 101-point running mean. Measured accelerations are valid at the location of the IMU, whereas the reconstruction process requires accelerations at the location of the centre of mass (CM). Angular corrections were investigated using the positional information reported in the file DOCUMENT/EDRSIS.PDF of [CATLINGETAL2008], but neglected because their effects were negligible below 100 km. The appropriate correction is vector omega-dot x vector r + vector omega x (vector omega x vector r) where vector omega is the angular velocity of the spacecraft as measured by the IMU gyroscopes, vector r is the separation of the IMU from the centre of mass, and a dot indicates a derivative with respect to time [TOLSONETAL1999]. The typical magnitude of omega was 10^(-2) rad s-1, the typical magnitude of omega-dot was 10^(-2) rad s-2, and the magnitude of r was about 1 m. Hence the two angular correction terms were 10^(-2) m s-2 and 10^(-4) m s-2, respectively. Comparison with Figure 1 shows that these corrections are small. Finally, data points prior to the time of entry were discarded. This completed the preparation of the IMU data. 4. Attitude Determination The aerodynamic accelerations were measured in a frame defined relative to the Phoenix spacecraft, but the trajectory reconstruction process must use a frame that is defined relative to Mars. For example, knowledge of the aerodynamic accelerations in a spacecraft-fixed frame does not provide knowledge as to whether the aerodynamic acceleration is directed up or down relative to the surface of Mars. The attitude of Phoenix must be known in a frame that can be related to Mars itself. There are a multitude of possible ways to describe the orientation of a rigid body. In this report, we use quaternions to describe the relationship of the "C frame" to the "J2000 frame", which is defined by the Earth's mean equator and equinox at the J2000.0 epoch. We adopted the convention that the quaternion [q0, q1, q2, q3] represents a rotation about the axis parallel to the vector qv = [q1, q2, q3] by angle 2 cos^(-1) q0. The time dependence of the quaternion that describes the orientation of the "C frame" relative to the "J2000" frame satisfies: Equation 2 dq0/dt = 0.5 q0 vector omega - vector qv x vector omega Equation 3 dvector qv / dt = -0.5 vector qv . vector omega where omega is the angular velocity of the spacecraft as measured by the IMU gyroscopes. The vector aerodynamic acceleration in the "J2000" frame, vector a_J2000, can be derived from the vector aerodynamic acceleration in the "C frame", vector a_C, as follows: Equation 4 vector a_J2000 = [1-2(q2^(2)+q3^(2)) 2(q2 q1-q0 q3) 2(q3 q1+q0 q3) ] vector a_C [2(q1 q2+q0 q3) 1-2(q1^(2)+q3^(2)) 2(q3 q2-q0 q1) ] [2(q1 q3-q0 q2) 2(q2 q3+q0 q1) 1-2(q1^(2)+q2^(2))] 5. Entry State Analysis of radio communications between Earth and Phoenix prior to entry were used to determine its trajectory prior to atmospheric entry. Since atmospheres do not have sharp edges, "atmospheric entry" is a slightly arbitrary concept. The Phoenix project adopted a distance of 3522.2 km from the centre of mass of Mars as its entry interface. The time, position and velocity at the entry interface are listed in Table 2 with their uncertainties. For computational convenience, the uncertainty on entry time, formally 0.002 s, was set to zero and the uncertainty on entry radius, formally zero, was set to the radial distance travelled by Phoenix at the entry radius during 0.002 s, or about 2.5 m. [TABLE 2] The orientation of Phoenix during its cruise to Mars was monitored and reconstructed by the spacecraft team. This was archived at the NAIF node of the PDS in the form of SPICE kernels (ftp://naif.jpl.nasa.gov/pub/naif/PHOENIX/kernels/ck/). The quaternion at the entry interface was obtained from the relevant SPICE kernels and is listed in Table 3. [TABLE 3] 6. Reconstructed Trajectory The equations of motion were integrated forward in time using methods described in previous papers and associated archives [WITHERSETAL2003B, WITHERS&SMITH2006A] to obtain time series of position and velocity. The only significant difference was that the transformation of acceleration measurements from a spacecraft-fixed frame into a planet-fixed frame used a direct method enabled by the derived quaternions (Section 4), rather than an indirect method based upon acceleration ratios [WITHERSETAL2003B]. The reconstructed trajectory is shown in Figures 3-6. All "altitudes" used in this report are radial distances above 3376.3 km (Table 1). Uncertainties in all the results reported in this work were obtained using Monte Carlo techniques (1000 reconstructions) and normally distributed uncertainties. The duration of the smoothing window used for the axial accelerations, which is shown in Figure 2, can be expressed as a vertical resolution if multiplied by the rate of change of altitude, as shown in Figure 7. [FIGURE 3] [FIGURE 4] [FIGURE 5] [FIGURE 6] [FIGURE 7] The angle of attack, namely the angle between the spacecraft symmetry axis and the vector velocity of the atmosphere relative to the spacecraft, was also determined using the assumption that the atmosphere rotates with the solid body of the planet. This directly determined angle of attack is labelled alpha_D. Results are shown in Figure 8. This angle can be determined directly if gyroscopic data are available, as they are for Phoenix, but must be determined indirectly if they are not [SPENCERETAL1999, WITHERSETAL2003B, WITHERS&SMITH2006A]. [FIGURE 8] The moment of parachute deployment is clearly visible in the time series of accelerations. The time and reconstructed position of parachute deployment are listed in Table 4. The moment of impact is also clearly visible in the time series of accelerations. The time and reconstructed position of the first ground contact are listed in Table 5. Results are consistent with [DESAIETAL2008] and [BLANCHARD2009], as well as the coordinates reported in Table 1. In particular, the speed at first ground contact differs from zero by only 0.1% of the entry speed, which indicates that the inferred orientation history of Phoenix is reliable. We do not address the near-surface orientation in this work, but it should be noted that a similar reconstruction by [BLANCHARD2009] found that the reconstructed orientation appeared unreliable near the surface where it varied rapidly, whereas the safe landing of Phoenix on its three legs implies that its orientation must have been relatively stable at this time. One possible explanation for this discrepancy is that near-surface winds, whose magnitude and direction are not known from observations, were not included in the reconstruction. [TABLE 4] [TABLE 5] 7. Reconstructed Atmospheric Structure Atmospheric density, pressure and temperature were determined from the reconstructed trajectory using the techniques established by previous missions [WITHERSETAL2003B, WITHERS&SMITH2006A]. The atmospheric reconstruction began at approximately 130 km altitude, some distance below the entry interface. Data close to the entry interface cannot be usefully smoothed because they are too close to the start of the data series, and so are not useful for the atmospheric reconstruction. According to the file CATALOG/INST.CAT of [CATLINGETAL2008], the nominal entry mass and reference area of Phoenix were 572.743 kg and pi D*D / 4, respectively, where D = 2.65 m. Uncertainties on mass and reference area were neglected. An atmospheric mean molar mass of 43.49 g mol-1 was assumed at all altitudes [MAGALHAESETAL1999]. This is an over-estimate at high altitudes due to diffusive separation above the homopause (about 125 km). Readers wishing to adopt a different mass profile can straight-forwardly multiply the temperatures reported here by the ratio of their desired mean molar mass to 43.49 g mol-1. The aerodynamic characteristics of Phoenix were modelled using a range of theoretical and experimental data sources [EDQUISTETAL2008]. The resultant database was used to relate unknown density to known axial aerodynamic acceleration. An upper boundary condition is needed to determine pressure from density using the equation of hydrostatic equilibrium. First, a reconstruction was performed without consideration of uncertainties. An exponential fit to data between the onset of atmospheric reconstruction and 10 km below this level found a scale height equivalent to a temperature of 150 K for an isothermal atmosphere. Since the product of pressure and mean molecular mass equals the product of density, Boltzmann's constant and temperature for an ideal gas, this fitted temperature can be used in conjunction with density at the upper boundary to determine a pressure at the upper boundary. This pressure was used as the upper boundary condition in the integration. The uncertainty in this pressure was assumed to be equivalent to an uncertainty in the temperature at the upper boundary of 50 K. The effects of the assumed upper boundary condition and its uncertainty are negligible more than a few scale heights below the onset of atmospheric reconstruction. The reconstructed atmospheric profiles are shown in Figures 9-12. [FIGURE 9] [FIGURE 10] [FIGURE 11] [FIGURE 12] 8. Consistency Checks The aerodynamic characteristics of entry vehicles such as Phoenix depend on their angle of attack. Since Phoenix's angle of attack was determined directly during the trajectory reconstruction using gyroscope data, it was not necessary to indirectly determine the angle of attack prior to finding the appropriate aerodynamic coefficient for the atmospheric structure reconstruction. However, the angle of attack can still be indirectly determined from the ratio of normal to axial accelerations [SPENCERETAL1999, WITHERSETAL2003B, WITHERS&SMITH2006A]. This ratio, which can be obtained from observations, equals the ratio of the normal force coefficent to the axial force coefficient. The latter ratio increases monotonically with increasing angle of attack and can be found from the modelled aerodynamic characteristics. The inferred value of the angle of attack is therefore the value necessary to make the modelled ratio of the normal force coefficent to the axial force coefficient equal the observed ratio of normal to axial accelerations. The resultant indirectly determined angles of attack, alpha_I, are shown in Figure 13 alongside the directly determined angles of attack, alpha_D. Values of alpha_I are typically 1-2 degrees smaller than alpha_D. As a specific example, alpha_D = 4.0 degrees and alpha_I = 2.7 degrees at t=1980 sec. One possible explanation involves winds. Values of alpha_D are calculated under the assumption that wind speed is zero everywhere, which is clearly unrealistic, whereas inferred values of alpha_I reflect the actual atmospheric winds during entry. However, this explanation is not sufficient. If a uniform zonal wind is assumed, then a wind speed of 300 m s-1 is required to adjust alpha_D at t=1980 sec from 4.0 degrees to 2.7 degrees. This is unphysically fast. Also, the angle of attack at parachute deployment is found to be 46 degrees, much larger than any values consistent with the actual safe landing of Phoenix. Other possible explanations (possessing varying degrees of plausibility) are that the single axis accelerometers have some cross-axis sensitivity, that the modelled aerodynamic characteristics are incorrect, and that the transformation between the "IMU frame" and "C frame" is incorrect. Resolution of this discrepancy is beyond the scope of this report. Note that this issue is different from discrepancies between predicted and reconstructed angles of attack addressed by [DESAIETAL2008]. [FIGURE 13] Small oscillations in the angle of attack can be seen in Figure 8. They have amplitudes on the order of tenths of a degree and periods of about two seconds. The period of these oscillations depends on the atmospheric density according to [SCHOENENBERGERETAL2005]: Equation 5 Omega * Omega = - rho V*V A D Cmalpha / (2 I) where Omega is the angular frequency of the oscillations, rho is atmospheric density, V is atmosphere-relative speed, A is a reference area, D is a reference diameter, Cmalpha is the derivative of the pitching moment coefficient with respect to angle of attack, and I is a moment of inertia. If the observed period of oscillations, reconstructed density and speed, previously stated reference area and diameter, and Cmalpha is about 0.1 rad s-1 [EDQUISTETAL2008] are used with Equation 5 to calculate I, then the resultant I is within about 10% of 200 kg m2 over much of the trajectory. The actual moment of inertia is close to 200 kg m2 [PRINCEETAL2008]. 9. Future Work Scientific analysis of these results is in progress (Withers and Catling, Scientific results of the Phoenix Atmospheric Structure Experiment, in preparation). A direct radio communications link between Phoenix and Earth was maintained during EDL. The time series of frequencies recorded on the ground could be used to perform an independent trajectory and atmospheric structure reconstruction. These data, which determine the line of sight speed of Phoenix relative to Earth, are sufficient for this task if it is assumed that aerodynamic acceleration is parallel to the atmosphere-relative velocity vector. This technique may also offer the possibility of near-real-time trajectory and atmospheric structure reconstruction for future missions. A radar onboard Phoenix recorded the range to the surface for altitudes of 0-2 km. Incorporation of this altitude and descent speed data into the trajectory reconstruction would correct the current non-zero altitude at first ground contact and might also be useful in determining atmospheric properties in the dynamic environment encountered after parachute deployment. The discrepancy between alpha_D and alpha_I is puzzling. Its resolution is more likely to improve understanding of characteristics of the spacecraft than understanding of the martian environment. 10. Conclusions During its flight through the atmosphere of Mars, Phoenix recorded acceleration and angular velocity data using accelerometers and gyroscopes within its inertial measurement unit. The archived EDR data require correction due to the use of a preliminary version of a transformation matrix between two reference frames. We have described a correction method using an updated version of this transformation matrix. The EDR data, which were acquired at a high rate of 200 Hz, are noisy, but this can be improved by smoothing at the cost of worsened vertical resolution. Smoothing of the axial acceleration requires care. We have described how the running arithmetic mean of the axial acceleration is a biased estimate of the true acceleration and explained how this bias can be corrected. The processed accelerations and angular velocities have been used to reconstruct the trajectory of Phoenix through the atmosphere and the associated atmospheric structure. The position and velocity of Phoenix were reconstructed by forward integration of the equations of motion and its orientation history was determined directly using onboard angular velocity measurements. The density profile along the trajectory was reconstructured from about 130 km to parachute deployment, about 13.5 km above the surface, using the drag equation. The corresponding pressure profile was obtained from the equation of hydrostatic equilibrium and the corresponding temperature profile was obtained from the ideal gas law. The angle of attack of Phoenix during atmospheric flight was determined directly using angular velocity measurements and indirectly using the ratio of normal to axial accelerations and the spacecraft's aerodynamic database. The former was typically 1-2 degrees greater than the latter, a discrepancy for which we do not have a good explanation. Acknowledgments PW acknowledges support from NASA (NNX09AG16G). DC acknowledges past support from the UK Science and Technology Facilities Council (STFC) awarded to the University of Bristol for Phoenix lander data reduction. PW and DC also acknowledge help and assistance from many people associated with the Phoenix project. Appendix A. Reference frames The IMU data were recorded in the "local IMU coordinate frame" for IMU-A, but were transformed into the "Phoenix mechanical frame" before being archived as EDRs [CATLINGETAL2008]. These are abbreviated here as the IMU frame and M frame, respectively. The M frame is also called the PHX_LANDER frame and is assigned the NAIF identifier -84001 (ftp://naif.jpl.nasa.gov/pub/naif/PHOENIX/kernels/fk/phx_v06.tf). The symmetry axis of the spacecraft, along which the largest accelerations occurred, is the +/- Z-axis in the M frame. During the course of this project, we discovered that a preliminary version of the appropriate transformation matrix, M-subscript-IMU-superscript-M,Prelim, had been used to transform the data from the IMU frame to the M frame for archiving. The preliminary matrix corresponds to the spacecraft as designed, not as built. The elements of this matrix, which were reported in the file DOCUMENT/EDRSIS.PDF of [CATLINGETAL2008], are: Equation A1 M-subscript-IMU-superscript-M,Prelim = [ 0 cos(delta) -sin(delta)] [-1 0 0] [ 0 sin(delta) cos(delta)] where delta = 25 degrees. For the analysis described in this report, the Phoenix cruise frame was used because its orientation, but not the orientation of the M frame, was readily available at entry. This is abbreviated here as the C frame. The C frame is also called the PHX_LANDER_CRUISE frame and is assigned the NAIF identifier -84000. The transformation matrix from the M frame to the C frame, M-subscript-M-superscript-C, was reported in the file DOCUMENT/EDRSIS.PDF of [CATLINGETAL2008]. Its elements are: Equation A2 M-subscript-M-superscript-C = [ 0 0 1] [cos(epsilon) -sin(epsilon) 0] [sin(epsilon) cos(epsilon) 0] where epsilon = 120 degrees. The transformation matrix from the IMU frame to the C frame, M-subscript-IMU-superscript-C, was reported in Lockheed-Martin Interoffice Memo "PHX-SE-07-0196: Phoenix Mars Lander Alignment Results" by A. Stoltz. Its elements are: Equation A3 M-subscript-IMU-superscript-C = [0.001658000000000 0.434532000000000 0.900655000000000] [0.865583671733631 -0.451637796435926 0.216304845776395] [0.500761102355330 0.779233610045474 -0.37687298280806 ] This matrix was determined from pre-launch measurements, not design specifications. Its elements differ from those found by combining M-subscript-IMU-superscript-M,Prelim and M-subscript-M-superscript-C only at the third decimal place, but these small differences have a surprisingly large effect on the reconstructed trajectory and atmospheric structure. If the preliminary matrix was used, then angles of attack reached ten degrees at parachute deployment, first ground contact and parachute deployment occurred several kilometres higher than in our final reconstruction, and temperatures below 30 km were 5K warmer than in our final reconstruction. We recommend that designers of future atmospheric structure investigations carefully evaluate how uncertainties in the orientation of accelerometer axes in a spacecraft-fixed frame affect their final scientific data products. The updated transformation matrix from the IMU frame to the M frame, M-subscript-IMU-superscript-M,Updated, corresponds to the spacecraft as built. It satisfies: Equation A4 M-subscript-IMU-superscript-M,Updated = Transpose(M-subscript-M-superscript-C) * M-subscript-IMU-superscript-C We transformed the originally archived accelerations and angular rates (EDRs) from the preliminary version of the M frame into the IMU frame using M-subscript-IMU-superscript-M,Prelim, then from the IMU frame into the C frame using M-subscript-IMU-superscript-C. Users of the originally archived EDR data should modify the archived values using the matrix M-subscript-M,Prelim-superscript-M,Updated, which satisfies: Equation A5 M-subscript-M,Prelim-superscript-M,Updated = M-subscript-IMU-superscript-M,Updated * Transpose(M-subscript-IMU-superscript-M,Prelim) = Transpose(M-subscript-M-superscript-C) * M-subscript-IMU-superscript-C * Transpose(M-subscript-IMU-superscript-M,Prelim) The numerical values of the elements of the matrix M-subscript-M,Prelim-superscript-M,Updated are as follows: Equation A6 M-subscript-M,Prelim-superscript-M,Updated = [ 0.99991267 -0.00087993926 -0.013188274 ] [ 0.00090164837 0.99999798 0.0016463352 ] [ 0.013186473 -0.0016580001 0.99991177 ] To provide a mechanism by which data users can verify that they have modified the EDR data correctly, we now state the values of some modified datapoints. The original (from file DATA/IMU_A_EDR_M.TAB of [CATLINGETAL2008]) and modified velocity changes in m s-1 along the X, Y and Z axes of the "M frame" at the first time step (RELATIVE TIME=0.000 seconds) are: Equation A7 [ 0.00013646 ] [ 0.00013799366 ] [-0.00271028 ] --> [ -0.0027100467 ] [ 0.00006364 ] [ 6.9927455E-5 ] The corresponding angle changes (radians) are: Equation A8 [ -0.00000786 ] [ -0.0015487050 ] [ 0.00000600 ] --> [ 0.0011955575 ] [ -0.00000918 ] [ -0.0018585567 ] The original and modified velocity changes (m s-1) at this timestep are: Equation A9 [ -0.00430083 ] [ 0.0012597930 ] [ -0.00150571 ] --> [ -0.0022035232 ] [ -0.42150492 ] [ -0.42152195 ] The corresponding angle changes (radians) are: Equation A10 [ -0.00001348 ] [ -1.3027308E-5 ] [ 0.00001100 ] --> [ 1.0930251E-5 ] [ -0.00003497 ] [ -3.5162906E-5 ] Appendix B. Reducing noise in axial acceleration measurements The noise in the axial acceleration measurements must be reduced by averaging in order for uncertainties in the reconstructed atmospheric properties to be small enough that the results have scientific value. However, the normal averaging process provides a biased estimate of the true axial acceleration at the centre of the time series. We now explain the source of this bias and how the bias was corrected. The axial acceleration measured before substantial deceleration occurs increases exponentially with time [MAGALHAESETAL1999]. This can also be seen in Figure 1. This happens because the axial acceleration, aa, satisfies: Equation B1 m aa = -CA rho V*V / 2 where m is spacecraft mass, CA is the axial force coefficient, rho is atmospheric density, V is atmosphere-relative speed and A is a reference area. Values of m and A are fixed, V is effectively constant until substantial deceleration occurs and CA changes only slightly during flight through the rarefied upper atmosphere. Hence aa is proportional to rho, which depends exponentially on altitude, z, through rho is proportional to exp(-z/H) where H is the density scale height. Since the rate of change of altitude is also approximately constant, the magnitude of aa increases exponentially with time. Dropping the subscript a for convenience, we have: Equation B2 a = a0 exp(t/tau) where t is time, a0 is the acceleration at time t=0, and tau is the timescale, which equals the ratio of the atmosphere's density scale height to the rate of change of altitude. The mean acceleration between t=-tX and t=tX, amean, satisfies: Equation B3 amean = 1 /(2 tX) x integral from t=-tX to t=tX of a0 exp(t/tau) dt Equation B4 amean = (a0 tau)/(2 tX) (exp(tX/tau) - exp(-tX/tau)) Equation B5 amean = a0 (tau/tX) sinh(tX/tau) The mean acceleration is not the same as the desired acceleration at the centre of the time series. Although the arithmetic mean of the logarithm of the acceleration does, in principle, give an unbiased estimate of the acceleration at the centre of the time series, this is not useful in practice because the noisy measured accelerations can be positive or negative. A process that requires finding the logarithm of a negative number is not helpful. According to Equation B5, the mean acceleration is significantly greater than the acceleration at the centre of the time series unless tX is much less than tau. However, if the timescale tau can be determined, then the value of a0 can be found from amean and tX. If a "long" average, aL, and a "short" average, aS, are calculated over the ranges t=-2tS to t=2tS and t=-tS to t=tS, respectively, then these averages satisfy: Equation B6 aL = a0 tau/(2 tS) sinh(2 tS/tau) Equation B7 aS = a0 tau/(tS) sinh(tS/tau) The ratio of the averages satisfies: Equation B8 aL/aS = 0.5 sinh(2 tS/tau) / sinh(tS/tau) Equation B9 aL/aS = cosh(tS/tau) Equation B9 can be rearranged using the trigonometric identity cosh^(-1) x = ln ( x + (x^(2) - 1)^(1/2) ) to yield: Equation B10 tS/tau = ln( aL/aS + sqrt( (aL/aS)*(aL/aS) - 1) ) Thus the timescale tau can be determined from the two related means, aL and aS, and then used in either Equation B6 or Equation B7 to find the true acceleration at the centre of the series of data points. Figure 14 illustrates differences between aL and aS for the Phoenix data. Figure 15 shows how these differences are drastically reduced after application of the correction procedure outlined above. [FIGURE 14] [FIGURE 15] For practical application of these theoretical results, two considerations are critical. Namely, whether to apply this correction and what duration to choose for the averaging window. If tX is sufficiently small in comparison to tau, then the uncorrected average is so close to the true axial acceleration that correcting it is more likely to worsen it than improve it. The duration of the averaging window should be chosen such that the combination of uncertainty and vertical resolution is optimal for the desired application. Creation of the composite series of smoothed axial accelerations for Phoenix required several steps. To begin, several series of smoothed, uncorrected axial accelerations were calculated using a 2^(i) point running mean. These were labelled ai, where i ranged from 6 to 13. Next, several series of smoothed, corrected axial accelerations were calculated by taking the ai series, then correcting them using the ratio to the a(i+1) series. These were labelled aiC, where i ranged from 10 to 12. Finally, the composite series was constructed as follows. First, unaveraged values were used to fill the first 4096 data points. Second, a12C values were used to fill the next 1081 data points. Third, a11C values were used to fill the next 2595 data points. Fourth, a10C values were used to fill the next 1882 data points. The use of corrected values cannot continue indefinitely, since the axial acceleration eventually ceases to vary exponentially with time. Also, it becomes unnecessary when the averaging window is small by comparison to the timescale for changes in the axial acceleration. Accordingly, uncorrected averages were used after this point. Fifth, a9 values were used to fill the next 114 data points. Sixth, a8 values were used to fill the next 127 data points. Seventh, a7 values were used to fill the next 337 data points. Eighth, a6 values were used to fill all remaining data points prior to parachute deployment. The duration of the averaging window used to calculate the a6 values, 0.32 sec, corresponds to averaging over a vertical distance of less than 0.5 km. The use of a shorter averaging window to improve vertical resolution further is not justified. Ninth, unaveraged values were used to fill all data points after parachute deployment. The atmospheric structure reconstruction terminates at parachute deployment, so there is little motivation to reduce noise by averaging after this point. In all cases, the locations of the transitions between one average and another average were fixed based on consideration of the scatter in the average with the shorter averaging window. These variations in the width of the window used to smooth the axial accelerations are shown in Figures 2 and 7. The impact of these variations on the uncertainty of the smoothed axial acceleration is shown in Figure 1. Table 1. Time and position of Phoenix landing. Time(1) (UTC) 2008-05-25T23:38:24 Aerocentric latitude(2) (deg N) 68.21878 +/- 0.00006 Longitude(2) (deg E) 234.24845 +/- 0.000096 Radius(2) (km) 3376.2915 +/- 0.0014 Elevation(3) (m) -4131 Footnote (1): [SMITHETAL2009]. Footnote (2): Martin-Mur (personal communication, 28 May 2008). The landed latitude, longitude and radius, which are shown with 1 sigma uncertainties, were determined from Doppler tracking. Footnote (3): Elevation is with respect to the areoid defined by the Mars Orbiter Laser Altimeter (MOLA) investigation, specifically 16 pixels per degree gridded MOLA data acquired from http://geo.pds.nasa.gov/missions/mgs/megdr.html [SMITHETAL2003]. Table 2. Phoenix entry state(1) Quantity Value 1 sigma Uncertainty Time (UTC) 2008-05-25T23:30:57.733 0.002(7) t - tref(2) (s) 1857.733 0.002(7) Radius(3) (km) 3522.2 0(8) Areocentric latitude (deg N) 69.3660 0.04 Longitude (deg E) 197.7160 0.00031 Speed(4) (m/s) 5600.34 0.02 Flight path angle(5) (deg) 13.010 0.00015 Heading angle(6) (deg) 77.6720 0.000055 Footnote (1): Information from the file CATALOG/INST.CAT of [CATLINGETAL2008]. Footnote (2): tref, a UTC reference time, is 2008-05-25T23:00:00.000. Footnote (3): 145.9 km above the landing site at 3376.3 km. Footnote (4): Relative to a Mars-centred inertial frame [SPENCERETAL1999, WITHERSETAL2003B]. Footnote (5): Angle below horizontal of velocity vector in inertial frame. Footnote (6): Angle east of north of velocity vector in inertial frame. Footnote (7): Set to zero as discussed in the text. Footnote (8): Set to 2.5 m as discussed in the text. Table 3. Quaternion representing the rotation between the "J2000 frame" and the "C frame" at the entry interface. q0 0.050125700 q1 0.52484401 q2 0.66515558 q3 0.52876671 Table 4. Reconstructed conditions and their 1 sigma uncertainties at parachute deployment. Altitude (km) 13.54 0.38 Areocentric latitude (deg N) 68.265 0.030 Longitude (deg E) 234.034 0.058 Atmosphere-relative speed (m/s) 390.9 1.5 Angle of attack (deg) 4.90 0.71 Mach number 1.703 0.015 Table 5. Reconstructed conditions and their 1 sigma uncertainties at first ground contact. Altitude (km) 1.10 1.49 Areocentric latitude (deg N) 68.237 0.029 Longitude (deg E) 234.311 0.054 Atmosphere-relative speed (m/s) 6.1 3.6 Figure 1. Time series of smoothed axial acceleration measurements (dots). The associated 1 sigma uncertainties are shown by the dashed line. Figure 2. Time-varying width of window used to smooth axial acceleration measurements. Figure 3. Reconstructed altitude versus time for Phoenix with 1 sigma uncertainties shown by the grey envelope. Figure 4. Reconstructed altitude versus reconstructed latitude for Phoenix with 1 sigma uncertainties shown by the grey envelope. Figure 5. Reconstructed altitude versus reconstructed longitude for Phoenix with 1 sigma uncertainties shown by the grey envelope. Figure 6. Reconstructed altitude versus reconstructed atmosphere-relative speed for Phoenix with 1 sigma uncertainties shown by the grey envelope. Figure 7. Reconstructed altitude versus vertical extent of window used to smooth axial acceleration measurements. Figure 8. Time series of directly determined angle of attack, alphaD, with 1 sigma uncertainties shown by the grey envelope. Results after parachute deployment are not shown. Figure 9. Reconstructed altitude versus reconstructed density with 1 sigma uncertainties shown by the grey envelope. Figure 10. Reconstructed altitude versus reconstructed pressure with 1 sigma uncertainties shown by the grey envelope. Figure 11. Reconstructed altitude versus reconstructed temperature with 1 sigma uncertainties shown by the grey envelope. 501-point running mean values of temperatures and temperature uncertainties are shown to reduce high frequency noise. Figure 12. Reconstructed pressure versus reconstructed temperature with 1 sigma uncertainties shown by the grey envelope. 501-point running mean values of temperatures and temperature uncertainties are shown to reduce high frequency noise. Figure 13. Time series of directly determined angle of attack, alphaD, (black line) and indirectly determined angle of attack, alphaI (grey line). Results after parachute deployment are not shown. Figure 14. Two time series of axial acceleration measurements before bias correction. Grey dots indicate smoothing with a 1024 point running mean and black dots indicate smoothing with a 2048 point running mean. Figure 15. Two time series of axial acceleration measurements after bias correction. Grey dots indicate smoothing with a 1024 point running mean followed by correction using ratio to a 2048 point running mean. Black dots indicate smoothing with a 2048 point running mean followed by correction using ratio to a 4096 point running mean.