Skip to content

File Mars.cpp

File List > astrea > astro > astro > systems > planetary_bodies > Mars > Mars.cpp

Go to the documentation of this file

/*
 * The GNU Lesser General Public License (LGPL)
 *
 * Copyright (c) 2025 Jay Iuliano
 *
 * This file is part of Astrea.
 * Astrea is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License
 * as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
 * Astrea is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty
 * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details. You should
 * have received a copy of the GNU General Public License along with Astrea. If not, see <https://www.gnu.org/licenses/>.
 */

#include <astro/systems/planetary_bodies/Mars/Mars.hpp>

#include <map>

#ifdef ASTREA_BUILD_MARS_EPHEMERIS
#include <astro/state/orbital_elements/OrbitalElements.hpp>
#include <ephemerides/Mars/MarsEphemerisTable.hpp>
#endif // ASTREA_BUILD_MARS_EPHEMERIS

namespace astrea {
namespace astro {
namespace planetary_bodies {

using mp_units::non_si::day;
using mp_units::si::unit_symbols::kg;
using mp_units::si::unit_symbols::km;
using mp_units::si::unit_symbols::m;

// Altitude Conditions(TABLE 7-4, Vallado)
static const std::map<Altitude, Density> martianAtmosphere = { // km, kg/m^3
    { 2.0 * km, 1.19e-1 * kg / (pow<3>(m)) },  { 4.0 * km, 1.10e-1 * kg / (pow<3>(m)) },
    { 6.0 * km, 1.02e-1 * kg / (pow<3>(m)) },  { 8.0 * km, 9.39e-2 * kg / (pow<3>(m)) },
    { 10.0 * km, 8.64e-2 * kg / (pow<3>(m)) }, { 12.0 * km, 7.93e-2 * kg / (pow<3>(m)) },
    { 14.0 * km, 7.25e-2 * kg / (pow<3>(m)) }, { 16.0 * km, 6.61e-2 * kg / (pow<3>(m)) },
    { 18.0 * km, 6.00e-2 * kg / (pow<3>(m)) }, { 20.0 * km, 5.43e-2 * kg / (pow<3>(m)) },
    { 22.0 * km, 4.89e-2 * kg / (pow<3>(m)) }, { 24.0 * km, 3.91e-2 * kg / (pow<3>(m)) },
    { 26.0 * km, 3.32e-2 * kg / (pow<3>(m)) }, { 28.0 * km, 2.82e-2 * kg / (pow<3>(m)) },
    { 30.0 * km, 2.40e-2 * kg / (pow<3>(m)) }, { 32.0 * km, 2.04e-2 * kg / (pow<3>(m)) },
    { 34.0 * km, 1.73e-2 * kg / (pow<3>(m)) }, { 36.0 * km, 1.47e-2 * kg / (pow<3>(m)) },
    { 38.0 * km, 1.25e-2 * kg / (pow<3>(m)) }, { 40.0 * km, 1.06e-2 * kg / (pow<3>(m)) },
    { 45.0 * km, 7.03e-3 * kg / (pow<3>(m)) }, { 50.0 * km, 4.67e-3 * kg / (pow<3>(m)) },
    { 55.0 * km, 3.10e-3 * kg / (pow<3>(m)) }, { 60.0 * km, 2.06e-3 * kg / (pow<3>(m)) },
    { 65.0 * km, 1.36e-3 * kg / (pow<3>(m)) }, { 70.0 * km, 9.11e-4 * kg / (pow<3>(m)) },
    { 75.0 * km, 6.05e-4 * kg / (pow<3>(m)) }, { 80.0 * km, 4.02e-4 * kg / (pow<3>(m)) }
};

Density Mars::find_atmospheric_density(const Date& date, const Distance& altitude) const
{
    // The values up to 80 km are almost definitely wrong.I can't find any
    // sources that contradict them though.Please fix them(and the
    // associated crash radius of Mars) if you can find better numbers.
    Unitless altitudeValue = altitude / astrea::detail::distance_unit;
    if (altitude <= 80.0 * km) {
        const auto iter = martianAtmosphere.upper_bound(altitude);
        return (iter != martianAtmosphere.end()) ? iter->second : 0.0 * kg / (m * m * m);
    }
    else if (altitude < 200.0 * km) {
        return exp(-2.55314e-10 * pow<5>(altitudeValue) + 2.31927e-7 * pow<4>(altitudeValue) -
                   8.33206e-5 * pow<3>(altitudeValue) + 0.0151947 * pow<2>(altitudeValue) - 1.52799 * altitudeValue + 48.69659) *
               kg / (m * m * m);
    }
    else if (altitude < 300.0 * km) {
        return exp(2.65472e-11 * pow<5>(altitudeValue) - 2.45558e-8 * pow<4>(altitudeValue) + 6.31410e-6 * pow<3>(altitudeValue) +
                   4.73359e-4 * pow<2>(altitudeValue) - 0.443712 * altitudeValue + 23.79408) *
               kg / (m * m * m);
    }
    else {
        return 0.0 * kg / (m * m * m);
    }
}

#ifdef ASTREA_BUILD_MARS_EPHEMERIS

RadiusVector<frames::solar_system_barycenter::icrf> get_position_at(const Date& date) const
{
    const auto positionJbFromSsb = get_position_at_impl<MarsEphemerisTable, frames::solar_system_barycenter::icrf>(date);
    return positionJbFromSsb; // TODO: Add correction for Mars' position from Mars barycenter
}

#endif // ASTREA_BUILD_MARS_EPHEMERIS

} // namespace planetary_bodies
} // namespace astro
} // namespace astrea