2f:null 30:["$","title",null,{"children":"Mars ephemerides computation from osculating orbital elements (2000-2050) | Christophe Yamahata | Observable"}] 35:T7d0,{ let dotplot = Plot.plot({ height: 470, width: 850, grid: true, marginLeft: 50, marginRight: 40, marginTop: 50, marginBottom: 75, insetTop: 0, insetBottom: 0, style: { fontSize: "17px", background: "#fff", color: "#000" }, x: { label: " ← Right ascension (RA), α ", tickPadding: 6, labelOffset: 64, transform: f => f/6, domain: [360/6,0], tickFormat: x => `${Math.round(6*x.toFixed(3))}°`, ticks: 10, tickSize: 9, tickRotate: 0 }, y: { label: "↑ Declination, δ", labelOffset: 35, tickFormat: x => x<=0 ?`${x}°`: `+${x}°`, ticks: 10, tickSize: 6, domain: [-30, 30] }, caption: html`
Figure 1: Theoretical astronomical coordinates of Mars for the period around the selected date.
⬤ Position of Mars at the selected date:
${date}.
┄ After selected date / ─ Before selected date
`, marks: [ Plot.frame({stroke: "black"}), Plot.ruleY([0]), Plot.line(ephemeris_mars_360_before, {x: "RA_deg", y: "Dec_deg", strokeWidth: 1.8, stroke: "crimson"} ), Plot.line(ephemeris_mars_360_after, {x: "RA_deg", y: "Dec_deg", strokeWidth: 2, strokeDasharray: "4,3", stroke: "crimson"} ), Plot.dot(ephemeris_mars_360_before, {x: coordinates_mars_selected_date.RA_deg, y: coordinates_mars_selected_date.Dec_deg, fill: "brown", r: 4.5}), Plot.text(ephemeris_mars_360_before, {x: coordinates_mars_selected_date.RA_deg , y: coordinates_mars_selected_date.Dec_deg, text: d => d.label, dx: 2, dy: -28, fill: "brown", fontSize: "28px"}), ] }) // Label in bold let select_labels = d3.select(dotplot) .selectAll("svg > g > text") .style("font-weight", "bold"); return dotplot; }36:T7f9,{ let dotplot = Plot.plot({ height: 470, width: 850, grid: true, marginLeft: 50, marginRight: 40, marginTop: 50, marginBottom: 75, insetTop: 0, insetBottom: 0, style: { fontSize: "17px", background: "#fff", color: "#000" }, x: { label: " Year → ", tickPadding: 6, labelOffset: 64, tickFormat: x => Math.round(x) == x? `Jan. \n ${x}` : `${x}`, tickSize: 9, ticks: 4, domain: [Math.floor(ephemeris_mars[0].decimal_year) , Math.floor(ephemeris_mars[0].decimal_year)+3], tickRotate: -35 }, y: { label: "↑ Distance Mars–Earth (AU)", labelOffset: 35, tickFormat: "", ticks: 8, tickSize: 6, domain: [0, 3] }, caption: html`
Figure 2: Theoretical distance of Mars from Earth for the period around the selected date.
⬤ Distance of Mars at the selected date:
${date}.
┄ After selected date / ─ Before selected date
`, marks: [ Plot.frame({stroke: "black"}), Plot.ruleY([0]), Plot.line(ephemeris_mars_360_before, {x: "decimal_year", y: "distance", strokeWidth: 1.8, stroke: "crimson"} ), Plot.line(ephemeris_mars_360_after, {x: "decimal_year", y: "distance", strokeWidth: 2, strokeDasharray: "4,3", stroke: "crimson"} ), Plot.dot(ephemeris_mars_360_before, {x: coordinates_mars_selected_date.decimal_year, y: coordinates_mars_selected_date.distance, fill: "brown", r: 4.5}), Plot.text(ephemeris_mars_360_before, {x: coordinates_mars_selected_date.decimal_year , y: coordinates_mars_selected_date.distance, text: d => d.label, dx: 2, dy: -28, fill: "brown", fontSize: "28px"}), ] }) // Label in bold let select_labels = d3.select(dotplot) .selectAll("svg > g > text") .style("font-weight", "bold"); return dotplot; }37:T578,

Heliocentric coordinates

Having found the true anomaly, ${tex`v`}, and the radius vector, ${tex`r`}, of the planet, we can now calculate the position of the planet with respect to the plane of the ecliptic (coordinates ${tex`X`}, ${tex`Y`}, ${tex`Z`}):
${tex`X = r [\cos(o) \times \cos(v + p - o) - \sin(o) \times \sin(v + p - o) \times \cos(i)]`}
${tex`Y = r [\sin(o) \times \cos(v + p - o) + \cos(o) \times \sin(v + p - o) \times \cos(i)]`}
${tex`Z = r [\sin(v + p - o) \times \sin(i)]`}

For the selected date, the distances are:
${tex`X`}_Earth = ${X(0,d).toFixed(4)} AU ${tex`X`}_Mars = ${X(1,d).toFixed(4)} AU
${tex`Y`}_Earth = ${Y(0,d).toFixed(4)} AU ${tex`Y`}_Mars = ${Y(1,d).toFixed(4)} AU
${tex`Z`}_Earth = ${Z(0,d).toFixed(4)} AU ${tex`Z`}_Mars = ${Z(1,d).toFixed(4)} AU 38:T5ba,// Osculating orbital elements for Earth-Moon barycenter and Mars barycenter // Date: 2024-06-01 // source: de440.bsp (JPL ephemeris DE440) osculating_orbital_elements = { var array = []; var data = { 0: { ID: "earth", name: "Earth-Moon barycenter", JD: 2460462.5, a: 1.00000, // Semi-major axis (astronomical unit, AU) e: 0.01669, // Eccentricity i: 0.00323, // Inclination (degrees) n: 0.98561, // Mean motion (degrees per day) o: 174.32115, // Longitude of the ascending node (degrees) w: 288.73527, // Argument of periapsis (degrees) MA: 146.57663, // Mean anomaly (degrees) p: 103.05641, // Longitude of periapsis (degrees) L: 249.63305 // Mean longitude (degrees) }, 1: { ID: "mars", name: "Mars barycenter", JD: 2460462.5, a: 1.52363, // Semi-major axis (astronomical unit, AU) e: 0.09328, // Eccentricity i: 1.84786, // Inclination (degrees) n: 0.52407, // Mean motion (degrees per day) o: 49.48904, // Longitude of the ascending node (degrees) w: 286.68524, // Argument of periapsis (degrees) MA: 12.33995, // Mean anomaly (degrees) p: 336.17428, // Longitude of periapsis (degrees) L: 348.51423 // Mean longitude (degrees) } }; Object.keys(data).forEach(key => { array.push(data[key]); }); return array; }39:Tca6,python_script = ` """ The following Python script was used to obtain the osculating orbital elements listed here: https://observablehq.com/@christophe-yamahata/mars-ephemerides-from-osculating-orbital-elements Author: C. Yamahata, May 2024 """ from skyfield.api import load from skyfield.elementslib import osculating_elements_of from skyfield.data.spice import inertial_frames def get_planet_elements(planets, planet, utc_time): """Get the osculating orbital elements of the specified planet. Skyfield project: https://pypi.org/project/skyfield/ Args: planets: Skyfield planets object. planet (str): Name of the planet. utc_time: Time object representing the Universal Time Coordinated. Returns: Osculating elements object. """ # Osculating orbital elements calculated with respect to the ecliptic plane ecliptic = inertial_frames['ECLIPJ2000'] Sun = planets['sun'] Planet = planets[planet] position = (Planet - Sun).at(utc_time) elements = osculating_elements_of(position, ecliptic) return elements def print_planet_elements(elements, planet, utc_time): """Print the osculating orbital elements of the specified planet.""" i = elements.inclination.degrees e = elements.eccentricity a = elements.semi_major_axis.au n = elements.mean_motion_per_day.degrees o = elements.longitude_of_ascending_node.degrees w = elements.argument_of_periapsis.degrees p = elements.longitude_of_periapsis.degrees M = elements.mean_anomaly.degrees L = elements.mean_longitude.degrees print('\\n-------------------------------------------') print(f'Planet: \\t {planet.capitalize()}') print(f'Reference date: {utc_time.utc_datetime().strftime("%Y-%m-%d")}') print(f'TDB Julian date: {utc_time.tdb:.2f}') print('\\nVariable \\t Value') print('-------------------------------------------') print(f'Semi-major axis: a = {a:.5f} au') print(f'Eccentricity: e = {e:.5f}') print(f'Inclination: i = {i:.5f}°') print(f'Mean motion: n = {n:.5f}°/d') print(f'Longitude of \\n ascending node: o = {o:.5f}°') print(f'Argument of \\n periapsis: w = {w:.5f}°') print(f'Mean anomaly: M = {M:.5f}°') print(f'Longitude of \\n periapsis: p = o + w = {p:.5f}°') print(f'Mean longitude: L = p + M = {L:.5f}°') print('-------------------------------------------\\n') # Reference date used: June 1, 2024 (2024-06-01) ref_year = 2024 ref_month = 6 ref_day = 1 timescale = load.timescale() utc_time = timescale.utc(ref_year, ref_month, ref_day, 0, 0) # Load the JPL ephemeris DE440 planets = load('de440.bsp') planet_names = ['earth barycenter', 'mars barycenter'] planet = planet_names[1] # 0: Earth, 1: Mars elements = get_planet_elements(planets, planet, utc_time) print_planet_elements(elements, planet, utc_time) # Calculate Earth obliquity using its osculating elements earth_position = (planets['earth barycenter'] - planets['sun']).at(utc_time) earth_elements = osculating_elements_of(earth_position) earth_obliquity = earth_elements.inclination.degrees print(f'Earth obliquity: {earth_obliquity:.5f}°') `3a:[{"id":"3fb700912c68df2326731ddcb9493ec1261c39856eaecc559af18bdbd796e5ddaf5fc2f16efd2d1d381aee261635cc619c7849b7452bd1ab84e4b61653653292","href":"https://static.observableusercontent.com/files/3fb700912c68df2326731ddcb9493ec1261c39856eaecc559af18bdbd796e5ddaf5fc2f16efd2d1d381aee261635cc619c7849b7452bd1ab84e4b61653653292","name":"Ephemerides_Mars_osculating_orbital_elements.PNG","type":"image/png","lastModified":1716594657598.615,"size":181006},{"id":"9784e920ae48e6b0ab6212332b0308d803a0a75d676f3f03534c99ed7ded56156dc457f11f851e3e299f1852bc26866f3e9c72c22e57aeaa7df6c9e52269d962","href":"https://static.observableusercontent.com/files/9784e920ae48e6b0ab6212332b0308d803a0a75d676f3f03534c99ed7ded56156dc457f11f851e3e299f1852bc26866f3e9c72c22e57aeaa7df6c9e52269d962","name":"Ephemerides_osculating_orbital_elements_Mars.xlsx","type":"application/vnd.openxmlformats-officedocument.spreadsheetml.sheet","lastModified":1716594746893.203,"size":1106064}] 33:{"id":"c982c52817935346","slug":"mars-ephemerides-from-osculating-orbital-elements","title":"Mars ephemerides computation from osculating orbital elements (2000-2050)","thumbnail":"8c3504c40f6322941f9cbb858f932227f9486dd61e4f2ada69c8310a5b9fc406","access_level":"public","create_time":"$D2024-05-23T13:22:12.950Z","update_time":"$D2026-04-20T20:31:42.725Z","view_time":null,"first_public_time":"$D2024-05-24T00:53:02.431Z","has_unpublished_edits":null,"fork":false,"stars":0,"starred":false,"owner":{"type":"team","id":"dd67891d326bcfb1","login":"christophe-yamahata","name":"Christophe Yamahata","avatar_url":"https://avatars-next.observableusercontent.com/8256ac716aa8395babba412b2835913d186874fa4190df57a4a28f2e60129876","tier":"starter_2024"},"version":1353,"publish_version":null,"role":null,"authored":false,"authors":[{"type":"individual","id":"06f083eba60f5a76","login":"christophe-yamahata","name":"Christophe Yamahata","avatar_url":"https://avatars-next.observableusercontent.com/8256ac716aa8395babba412b2835913d186874fa4190df57a4a28f2e60129876","tier":null}],"members":[{"type":"individual","id":"06f083eba60f5a76","login":"christophe-yamahata","name":"Christophe Yamahata","avatar_url":"https://avatars-next.observableusercontent.com/8256ac716aa8395babba412b2835913d186874fa4190df57a4a28f2e60129876","tier":null,"role":"editor"},{"type":"team","id":"dd67891d326bcfb1","login":"christophe-yamahata","name":"Christophe Yamahata","avatar_url":"https://avatars-next.observableusercontent.com/8256ac716aa8395babba412b2835913d186874fa4190df57a4a28f2e60129876","tier":"starter_2024","role":"viewer"}],"comments":[],"custom_thumbnail":true,"trashed":false,"lapsed":null,"subscribed":null,"forks":1,"collection":{"id":"6f8011976c15eeed","slug":"astronomy","title":"Astronomy","owner":{"id":"dd67891d326bcfb1","login":"christophe-yamahata","name":"Christophe Yamahata"}},"type":"notebook","body":{"title":"Mars ephemerides computation from osculating orbital elements (2000-2050)","cells":[{"id":0,"value":"# Mars ephemerides computation from osculating orbital elements (${Year_min}-${Year_max})","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":943,"value":"

\n Right Ascension and Declination are essential astronomical coordinates that provide a consistent way to locate objects in the sky. For planets and stars, they are independent of the observer's location on Earth. \n

This notebook computes the Right Ascension (RA) and Declination (Dec) for the planet Mars using osculating orbital elements retrieved for ${format_date(UT_from_JD(JD_ref))}. The calculations are made around the date selected in the calendar displayed below (i.e., ${format_date(UT_from_JD(JD_selected_date))}). You may compare the results obtained with this method with those shown in this notebook.\n

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":294,"value":"
\n

Osculating orbital elements

\n The osculating orbital elements can be extracted from the NASA horizons database.\n Hereafter, we provide the values extracted from the JPL ephemeris DE440 for Earth and Mars for the following date:\n
${format_date(UT_from_JD(JD_ref))} \n
(Julian Date, JD = ${JD_ref})\n

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":69,"value":"{ \n\n // This formatting is required to display the values properly in the table.\n const keys = ['a', 'e', 'i', 'n', 'o', 'w', 'MA', 'p', 'L'];\n const format = Object.fromEntries(keys.map(key => [key, x => x.toFixed(5)]));\n \n return Inputs.table(osculating_orbital_elements, {\n columns: [\n \"name\",\n \"a\", \n \"e\",\n \"i\",\n \"n\",\n \"o\",\n \"w\",\n \"MA\",\n \"p\",\n \"L\",\n ],\n \n header: {\n name: \"Planet\",\n ...variables_short,\n },\n \n align: {\n name: \"left\",\n a: \"right\",\n e: \"right\",\n i: \"right\",\n n: \"right\",\n o: \"right\",\n w: \"right\",\n MA: \"right\",\n p: \"right\",\n L: \"right\",\n },\n\n format: format,\n \n width: {\n name: 135,\n },\n\n maxWidth: 950,\n maxHeight: 100,\n height: 80\n})\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":247,"value":"{\n var data = '';\n let joinedVariables = {};\n\n // Join keys of variables_short with those of variable_long\n for (const key of Object.keys(variables_short)) {\n joinedVariables[key] = {\n short: variables_short[key],\n long: variables_long[key]\n };\n data += `${joinedVariables[key].short}: \n ${joinedVariables[key].long}
`;\n }\n return html`

${data}

`;\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":373,"value":"
\n

Mars position

The method used to calculate the planets positions (Earth and Mars) from the osculated orbital elements is described here.\n

\nBy default, this notebook calculates the position of Mars for the current date. You can also manually enter a date or select one from the calendar (year range: ${Year_min}-${Year_max}).\n

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":379,"value":"viewof date = Inputs.date({label: htl.html`Input date:`, value: Date.now(), required: true, min: Year_min+\"-01-01\", max: Year_max+\"-12-31\"})","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":450,"value":"The calculations details below are for the following date:\n
${format_date(date)} \n
(Julian Date, JD = ${JD_selected_date})","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":833,"value":"

\nFor this date, the calculated Right Ascension (RA) and the Declination (Dec) of Mars are: \n
RA = ${RA(d)} \n
Dec = ${DEC(d)} \n

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":841,"value":"

\n The calculated distance Mars-Earth is: \n
${tex`D`}_Mars-Earth = ${distance_Mars_Earth(d).toFixed(4)} AU\n

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":931,"value":"
\n

Right Ascension and Declination of Mars (♂)

\n The following graph plots the theoertical path of Mars around the selected date, as calculated from the osculating orbital elements.\n Notice the occasional loops, which occur because Mars undergoes apparent retrograde motion approximately once every 26 months (~2.1 years).\n

\n

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":893,"value":"$35","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":1317,"value":"
","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":1313,"value":"### Mars retrograde periods (2020 – 2040)\n\n| Year (start) | Period | \n|--------------|----------------------------------|\n| 2020 | Sept. 9 – Nov. 13 | \n| 2022 | Oct. 30, 2022 – Jan. 12, 2023 |\n| 2024 | Dec. 6, 2024 – Feb. 23, 2025 |\n| 2027 | Jan. 10 – April 1 |\n| 2029 | Feb. 14 – March 5 |\n| 2031 | March 28 – June 13 |\n| 2033 | May 26 – Aug. 1 |\n| 2035 | Aug. 15 – Oct. 15 |\n| 2037 | Oct. 12 – Dec. 23 |\n| 2039 | Nov. 23, 2039 – Feb. 9, 2040 |","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":974,"value":"
\n

Distance of Mars (♂) from Earth

\n The following graph plots the theoretical distance Mars–Earth around the selected date, as calculated from the osculating orbital elements.\n Compare this plot with the previous one: Mars is a superior planet that exhibits apparent retrograde motion around the time of opposition, when Earth is between Mars and the Sun (i.e., during their closest approach).\n

\n

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":981,"value":"$36","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":417,"value":"

Computation from the osculating orbital elements

Mean Anomaly

\n The Mean Anomaly is calculated using the following formula:\n
${tex`M = n \\times d + L - p`}\n
where ${tex`d`} is the number of days from the refence date.\n
Here, ${tex`d = ${JD_selected_date} - ${JD_ref} = ${d}`} day${d>1 || d<-1? html`s`: html``}.\n

\n For the selected date, the mean anomalies are: \n
${tex`M`}_Earth = ${MA(0,d).toFixed(2)}°\n
${tex`M`}_Mars = ${MA(1,d).toFixed(2)}°\n","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":472,"value":"

True Anomaly

\n While the mean Anomaly is calculated considering the orbit to be circular, the True Anomaly gives the actual position of the planet as it considers orbit to be elliptical. We convert Mean Anomaly to True Anomaly, ${tex`v`}, using the following formula:\n
${tex`v = M + \\frac{180}{\\pi} \\times [ (2e - \\frac{1}{4}e^3) \\times \\sin(M) \n + \\frac{5}{4}e^2 \\times \\sin(2M) + \\frac{13}{12} e^3 \\times \\sin(3M) ]`}\n

\n For the selected date, the true anomalies are: \n
${tex`v`}_Earth = ${TA(0,d).toFixed(2)}°\n
${tex`v`}_Mars = ${TA(1,d).toFixed(2)}°\n","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":528,"value":"

Distance to the Sun

\n The Sun occupies one of the foci of the ellipse. The distance, ${tex`r`}, from the planet to the Sun is given by a simple formula based on the geometry of the ellipse:\n
${tex`r = \\frac{a (1 - e^2)}{ 1 + e \\cos(v)}`}\n

\n For the selected date, the distances are: \n
${tex`r`}_Earth = ${r_AU(0,d).toFixed(4)} AU\n
${tex`r`}_Mars = ${r_AU(1,d).toFixed(4)} AU\n","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":546,"value":"$37","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":592,"value":"

Geocentric ecliptic coordinates

\n Changing the origin of the coordinate system from the Sun to the Earth is straightforward. We just need to subtract the Earth's coordinates from those of the planet Mars.\n
The geocentric ecliptic coordinates of Mars are (coordinates ${tex`X'`}, ${tex`Y'`}, ${tex`Z'`}):\n
${tex`X' = X_{Mars} - X_{Earth}`}\n
${tex`Y' = Y_{Mars} - Y_{Earth}`}\n
${tex`Z' = Z_{Mars} - Z_{Earth}`}\n

\n For the selected date, the geocentric ecliptic cordinates of Mars are: \n
${tex`X'`} = ${X_prime(d).toFixed(4)} AU\n
${tex`Y'`} = ${Y_prime(d).toFixed(4)} AU\n
${tex`Z'`} = ${Z_prime(d).toFixed(4)} AU\n","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":625,"value":"

Geocentric equatorial coordinates

\n To change the coordinate system from geocentric ecliptic to geocentric equatorial is just a matter of a rotation around the X-axis by an angle equal to the obliquity of the ecliptic, ${tex`\\epsilon`} = ${obliquity_ecliptic_Earth.toFixed(5)}°.\n
The geocentric equatorial coordinates of Mars are (coordinates ${tex`X_{q}`}, ${tex`Y_{q}`}, ${tex`Z_{q}`}):\n
${tex`X_{q} = X'`}\n
${tex`Y_{q} = Y' \\cos(\\epsilon) - Z' \\sin(\\epsilon)`}\n
${tex`Z_{q} = Y' \\sin(\\epsilon) + Z' \\cos(\\epsilon)`}\n

\n For the selected date, the geocentric ecliptic cordinates of Mars are: \n
${tex`X_{q}`} = ${X_q(d).toFixed(4)} AU\n
${tex`Y_{q}`} = ${Y_q(d).toFixed(4)} AU\n
${tex`Z_{q}`} = ${Z_q(d).toFixed(4)} AU\n","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":664,"value":"

Distance of Mars from Earth

\nThe distance of Mars from Earth is calculated as follows:\n
${tex`D`}_Mars-Earth = ${tex`(X_{q}^2+Y_{q}^2+Z_{q}^2)^{1/2}`}\n

\n For the selected date, the distance Mars-Earth is: \n
${tex`D`}_Mars-Earth = ${distance_Mars_Earth(d).toFixed(4)} AU\n

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":724,"value":"

Right Ascension

\n The Right Ascension of Mars (RA, ${tex`\\alpha`}) is calculated as follows:\n
${tex`\\alpha`} = ${tex`\\arctan(\\frac{Y_q}{X_q})`}, with ${tex`0° \\le \\alpha < 360°`}.\n
Note that RA is usually expressed in hours:\n
0h ≤ RA < 24h (i.e., angle in degrees divided by 15).\n

\n For the selected date, the Right Ascension of Mars is: \n
${tex`\\alpha`} = ${alpha(d).toFixed(4)}° ${tex` \\quad \\Leftrightarrow \\quad`} RA = ${(alpha(d)/15).toFixed(4)} h = ${RA(d)} \n

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":793,"value":"

Declination

\n The Declination of Mars (Dec, ${tex`\\delta`}) is calculated as follows:\n
${tex`\\delta = \\arctan(\\frac{Z_q}{(X_q^2 + Y_q^2)^{1/2}})`}, with ${tex`-90° \\le \\delta < 90°`}\n

\n For the selected date, the Declination of Mars is: \n
${tex`\\delta`} = ${delta(d).toFixed(4)}° ${tex` \\quad \\Leftrightarrow \\quad`} Dec = ${DEC(d)} \n

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":1140,"value":"
\n

Useful links

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":1142,"value":"

\n Planet positions using elliptical orbits [URL visited in May 2024] \n
→ This webpage details the steps required to calculate ephemerides from osculating orbital elements. \n
NASA HORIZONS system. \n
→ Access this database to obtain the osculating orbital elements\n
\n Skyfield Python package / Skyfield documentation [URL visited in May 2024] \n
→ Skyfield can be used as an alternative to the NASA HORIZONS database for calculating osculating orbital elements relative to the ecliptic plane. The data in this notebook were obtained using Skyfield and the JPL ephemeris DE440 (see the Python script below).

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":1238,"value":"
\n

Excel file

📈For your convenience, I have included an Excel file (see data pane) that replicates the calculations performed in this notebook.

\n","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":1277,"value":"

\n
The Excel file replicates the calculations performed in this notebook.","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":1179,"value":"
\n

Python script (Skyfield)

\nThe following Python script was used to obtain the osculating orbital elements for Mars and Earth. By default, the parameters are calculated with respect to the equatorial plane. In this Python script, the osculating orbital elements are calculated relative to the ecliptic plane.\n

","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":1217,"value":"Inputs.button(\"Copy script to clipboard\", {value: null, reduce: () => navigator.clipboard.writeText(python_script)})","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":1152,"value":"html`

\n${python_script}\n

`","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":320,"value":"

","pinned":false,"hidden":false,"mode":"html","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":5,"value":"## Appendices","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined","output":"","database":"$undefined"},{"id":408,"value":"// Mean Anomaly\nMA = function(planet,d) { \n var n = osculating_orbital_elements[planet].n;\n var L = osculating_orbital_elements[planet].L;\n var p = osculating_orbital_elements[planet].p;\n var MA = n * d + L - p;\n MA = MA - 360*Math.floor(MA/360) // 0° ≤ MA < 360°\n return MA; \n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":501,"value":"// True Anomaly\nTA = function(planet,d) { \n var M = MA(planet,d);\n var e = osculating_orbital_elements[planet].e\n var TA = M + (180/Math.PI)*[ (2*e-Math.pow(e,3)/4)*Math.sin(M*(Math.PI/180)) \n + (5/4)*Math.pow(e,2)*Math.sin(2*M*(Math.PI/180))\n + (13/12)*Math.pow(e,3)*Math.sin(3*M*(Math.PI/180)) ];\n TA = TA - 360*Math.floor(TA/360) // 0° ≤ TA < 360°\n return TA; \n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":534,"value":"r_AU = function(planet,d) { \n var v = TA(planet,d);\n var e = osculating_orbital_elements[planet].e;\n var a = osculating_orbital_elements[planet].a;\n var r = a * (1- Math.pow(e,2)) / (1 + e*Math.cos(v*Math.PI/180));\n return r; \n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":572,"value":"X = function(planet,d) { \n var r = r_AU(planet,d);\n var o = osculating_orbital_elements[planet].o; \n var v = TA(planet,d);\n var p = osculating_orbital_elements[planet].p;\n var i = osculating_orbital_elements[planet].i;\n var X = r * [ Math.cos(o*Math.PI/180)*Math.cos((v+p-o)*Math.PI/180) \n - Math.sin(o*Math.PI/180)*Math.sin((v+p-o)*Math.PI/180)*Math.cos(i*Math.PI/180) ];\n return X; \n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":584,"value":"Y = function(planet,d) { \n var r = r_AU(planet,d);\n var o = osculating_orbital_elements[planet].o; \n var v = TA(planet,d);\n var p = osculating_orbital_elements[planet].p;\n var i = osculating_orbital_elements[planet].i;\n var Y = r * [ Math.sin(o*Math.PI/180)*Math.cos((v+p-o)*Math.PI/180) \n + Math.cos(o*Math.PI/180)*Math.sin((v+p-o)*Math.PI/180)*Math.cos(i*Math.PI/180) ];\n return Y; \n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":587,"value":"Z = function(planet,d) { \n var r = r_AU(planet,d);\n var o = osculating_orbital_elements[planet].o; \n var v = TA(planet,d);\n var p = osculating_orbital_elements[planet].p;\n var i = osculating_orbital_elements[planet].i;\n var Z = r * [ Math.sin((v+p-o)*Math.PI/180)*Math.sin(i*Math.PI/180) ];\n return Z; \n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":613,"value":"X_prime = function(d) { \n return X(1,d)-X(0,d);\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":619,"value":"Y_prime = function(d) {\n return Y(1,d)-Y(0,d);\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":621,"value":"Z_prime = function(d) {\n return Z(1,d)-Z(0,d);\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":648,"value":"X_q = function(d) {\n var X_q = X_prime(d);\n return X_q;\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":652,"value":"Y_q = function(d) {\n var Y_q = Y_prime(d)*Math.cos(obliquity_ecliptic_Earth*Math.PI/180)-Z_prime(d)*Math.sin(obliquity_ecliptic_Earth*Math.PI/180);\n return Y_q;\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":656,"value":"Z_q = function(d) {\n var Z_q = Y_prime(d)*Math.sin(obliquity_ecliptic_Earth*Math.PI/180)+Z_prime(d)*Math.cos(obliquity_ecliptic_Earth*Math.PI/180);\n return Z_q;\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":667,"value":"distance_Mars_Earth = function(d) {\n var distance = Math.sqrt( Math.pow(X_q(d),2) + Math.pow(Y_q(d),2) + Math.pow(Z_q(d),2)) \n return distance;\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":748,"value":"alpha = function(d) {\n var alpha = Math.atan(Y_q(d)/X_q(d))*180/Math.PI;\n alpha = X_q(d) <0? alpha + 180: alpha; // If X_q is negative then add 180 degrees to alpha \n alpha = X_q(d) >0 && Y_q(d) <0 ? alpha + 360: alpha; // If X_q is positive and Yq is negative then add 360 degrees to alpha\n return alpha;\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":786,"value":"RA = function(d) {\n var alpha_deg = alpha(d);\n \n function convertToHMS(decimalHours) {\n const totalSeconds = decimalHours * 3600;\n \n const hours = Math.floor(totalSeconds / 3600);\n const minutes = Math.floor((totalSeconds % 3600) / 60);\n const seconds = (totalSeconds % 60).toFixed(2); // Rounds to 2 decimal places\n \n // Add leading zeros\n const paddedHours = String(hours).padStart(2, '0');\n const paddedMinutes = String(minutes).padStart(2, '0');\n const paddedSeconds = String(seconds).padStart(5, '0'); // Pad to 5 characters including the decimal point\n\n return `${paddedHours}h ${paddedMinutes}m ${paddedSeconds}s`;\n }\n \n return convertToHMS(alpha_deg/15);\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":808,"value":"delta = function(d) {\n var delta = Math.atan(Z_q(d)/Math.sqrt(Math.pow(X_q(d),2)+Math.pow(Y_q(d),2)))*180/Math.PI;\n return delta;\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":819,"value":"DEC = function(d) {\n\n var delta_deg = delta(d);\n \n function convertToDMS(decimalDegrees) {\n const isNegative = decimalDegrees < 0;\n decimalDegrees = Math.abs(decimalDegrees);\n \n const degrees = Math.floor(decimalDegrees);\n const minutesDecimal = (decimalDegrees - degrees) * 60;\n const minutes = Math.floor(minutesDecimal);\n const seconds = ((minutesDecimal - minutes) * 60).toFixed(2); // Rounds to 2 decimal places\n \n // Add leading zeros\n const paddedMinutes = String(minutes).padStart(2, '0');\n const paddedSeconds = String(seconds).padStart(5, '0'); // Pad to 5 characters including the decimal point\n \n const sign = isNegative ? '-' : '+';\n return `${sign}${degrees}° ${paddedMinutes}' ${paddedSeconds}\"`;\n }\n\n return convertToDMS(delta_deg);\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":300,"value":"// Unix timestamp (UT) calculated from the Julian Date (JD)\nUT_from_JD = JD => (JD - 2440587.5) * 86400000","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":399,"value":"UT_from_date = date => Math.floor(new Date(date).getTime())","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":951,"value":"format_date = d3.utcFormat(\"%B %e, %Y\")","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":1049,"value":"function JDToDecimalYear(JD) {\n // Define the reference Julian Date for January 1, 2000 (J2000.0)\n const J2000 = 2451545.0;\n\n // Calculate the number of days since J2000\n const daysSinceJ2000 = JD - J2000;\n\n // Calculate the corresponding year (approximately)\n const approximateYear = 2000 + daysSinceJ2000 / 365.25;\n\n // Determine the integer year\n const year = Math.floor(approximateYear);\n\n // Calculate the Julian Date for January 1st of the determined year\n const JD_Jan1 = new Date(Date.UTC(year, 0, 1)).getTime() / 86400000 + 2440587.5;\n\n // Calculate the total number of days in the determined year\n const isLeapYear = (year % 4 === 0 && year % 100 !== 0) || (year % 400 === 0);\n const totalDays = isLeapYear ? 366 : 365;\n\n // Calculate the day of the year\n const dayOfYear = JD - JD_Jan1 + 1;\n\n // Calculate the decimal year\n const decimalYear = year + (dayOfYear - 1) / totalDays;\n\n return decimalYear;\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":875,"value":"dates = function(i) {\n var UT = Date.parse(date);\n var date_i = new Date(UT+86400000*i)\n return date_i\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":364,"value":"// Julian Date calculated from calculated from the Unix timestamp (UT)\nJD = UT => Math.round(UT / 86400 / 1000,0) + 2440587.5","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":316,"value":"// Julian Date (JD) used for the osculating orbital elements listed in the table\nJD_ref = osculating_orbital_elements[0].JD","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":393,"value":"JD_selected_date = JD(UT_from_date(date))","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":419,"value":"d = JD_selected_date - JD_ref","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":868,"value":"ephemeris_mars = {\n\nvar data = []\n var i_min = -Math.round(span_days*0.4);\n var i_max = Math.round(span_days*0.6);\n for(var i = i_min; i < i_max; i++) {\n data[i-i_min] = {\n date: dates(i),\n decimal_year: parseFloat(JDToDecimalYear(JD_selected_date+i).toFixed(3)),\n RA_deg: parseFloat(alpha(d+i).toFixed(3)),\n Dec_deg: parseFloat(delta(d+i).toFixed(3)),\n distance: parseFloat(distance_Mars_Earth(d+i).toFixed(3)),\n }\n }\n return data\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":898,"value":"ephemeris_mars_360_before = ephemeris_mars.filter(d => d.date <= date).map(d => d.RA_deg <= 2 || d.RA_deg >= 358 ? NaN : d);","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":905,"value":"ephemeris_mars_360_after = ephemeris_mars.filter(d => d.date >= date).map(d => d.RA_deg <= 2 || d.RA_deg >= 358 ? NaN : d);","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":908,"value":"coordinates_mars_selected_date = {\n let data = ephemeris_mars.filter(d => d.date <= date && d.date >= date)[0];\n data.label = \"♂\";\n return data;\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":890,"value":"span_days = 600","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":382,"value":"Year_min = 2000","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":385,"value":"Year_max = 2050","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":234,"value":"variables_short = {\n return {\n a: \"a (AU)\",\n e: \"e\",\n i: \"i (°)\",\n n: \"n (°/day)\",\n o: \"o (°)\",\n w: \"w (°)\",\n MA: \"M (°)\",\n p: \"p (°)\",\n L: \"L (°)\",\n }\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":242,"value":"variables_long = {\n return {\n a: \"Semi-major axis\",\n e: \"Eccentricity\",\n i: \"Inclination\",\n n: \"Mean motion\",\n o: \"Longitude of the ascending node\",\n w: \"Argument of periapsis\",\n MA: \"Mean anomaly\",\n p: \"Longitude of periapsis\",\n L: \"Mean longitude\",\n }\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":632,"value":"// Obliquity of the ecliptic: Earth's axial tilt\n// Date: 2024-06-01\n// source: de440.bsp (JPL ephemeris DE440)\nobliquity_ecliptic_Earth = 23.43607","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":8,"value":"$38","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"},{"id":1212,"value":"$39","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined","output":"$undefined","database":"$undefined"}],"files":"$W3a","stdlib":"1","titling":"auto"}} 34:"$33:body"