mirror/OpenSpace
1
0
Fork 0
mirror of https://github.com/OpenSpace/OpenSpace.git synced 2026-05-04 18:11:01 -05:00
Code Issues Packages Projects Releases Wiki Activity

Add habitable zone representation for exoplanet systems

Browse Source
This commit is contained in:
Emma Broman
2020-12-11 13:48:56 +01:00
parent 51869784fd
commit 05651d1aa3
4 changed files with 134 additions and 25 deletions
Download Patch File Download Diff File Expand all files Collapse all files
modules/exoplanets/exoplanetshelper.cpp
+35
View File
@@ -136,6 +136,41 @@ glm::dmat3 computeSystemRotation(glm::dvec3 starPosition) {
);
}
std::optional<glm::vec2> computeHabitableZone(float teff, float luminosity) {
// The formula only consider stars with teff in range [2600, 7200] K
if (teff > 7200.f || teff < 2600.f) {
return std::nullopt;
}
struct Coefficients {
float seffSun;
float a, b, c, d;
};
// Coefficients for planets of 1 Earth mass. Received from:
// https://depts.washington.edu/naivpl/sites/default/files/HZ_coefficients.dat
constexpr Coefficients coefficients[] = {
// Inner boundary - Runaway greenhouse
{1.10700E+00, 1.33200E-04, 1.58000E-08, -8.30800E-12, -1.93100E-15},
// Outer boundary - Maximum greenhouse
{3.56000E-01, 6.17100E-05, 1.69800E-09, -3.19800E-12, -5.57500E-16}
};
const float tstar = teff - 5780.f;
const float tstar2 = tstar * tstar;
glm::vec2 distances;
for (int i = 0; i < 2; ++i) {
const Coefficients& coeffs = coefficients[i];
float seff = coeffs.seffSun + (coeffs.a * tstar) + (coeffs.b * tstar2) +
(coeffs.c * tstar * tstar2) + (coeffs.d * tstar2 * tstar2);
distances[i] = std::pow(luminosity / seff, 0.5f);
}
return std::optional<glm::vec2>(std::move(distances));
}
std::string createIdentifier(std::string name) {
std::replace(name.begin(), name.end(), ' ', '_');
std::replace(name.begin(), name.end(), '.', '-');
modules/exoplanets/exoplanetshelper.h
+20 -7
View File
@@ -26,6 +26,7 @@
#define __OPENSPACE_MODULE_EXOPLANETS___EXOPLANETSHELPER___H__
#include <ghoul/glm.h>
#include <optional>
#include <string>
#include <vector>
@@ -60,9 +61,9 @@ struct ExoplanetDataEntry {
float rStar; // Estimated radius of the star in solar radii
float rStarUpper; // Upper uncertainty of estimated star radius
float rStarLower; // Lower uncertainty of estimated star radius
float luminosity; // Star luminosity, in units of solar luminosities [log(Solar)]
float luminosityUpper; // Upper uncertainty of star luminosity [log(Solar)]
float luminosityLower; // Lower uncertainty of star luminosity [log(Solar)]
float luminosity; // Star luminosity, in units of solar luminosities
float luminosityUpper; // Upper uncertainty of star luminosity
float luminosityLower; // Lower uncertainty of star luminosity
float teff; // Star's effective temperature in Kelvin
float teffUpper; // Upper uncertainty of effective temperature
float teffLower; // Lower uncertainty of effective temperature
@@ -79,10 +80,10 @@ struct ExoplanetDataEntry {
struct StarData {
glm::vec3 position = glm::vec3(std::numeric_limits<float>::quiet_NaN()); // In parsec
float radius = std::numeric_limits<float>::quiet_NaN(); // In solar radii
float radius = std::numeric_limits<float>::quiet_NaN(); // In solar radii
float bv = std::numeric_limits<float>::quiet_NaN();
float teff = std::numeric_limits<float>::quiet_NaN(); // In Kelvin
float luminosity = std::numeric_limits<float>::quiet_NaN(); // In log(Solar)
float teff = std::numeric_limits<float>::quiet_NaN(); // In Kelvin
float luminosity = std::numeric_limits<float>::quiet_NaN(); // In solar luminosities
};
struct ExoplanetSystem {
@@ -100,12 +101,24 @@ bool hasSufficientData(const ExoplanetDataEntry& p);
// Compute star color in RGB from b-v color index
glm::vec3 starColor(float bv);
glm::dmat4 computeOrbitPlaneRotationMatrix(float i, float bigom, float omega);
glm::dmat4 computeOrbitPlaneRotationMatrix(float i, float bigom = 180.f,
float omega = 90.f);
// Rotate the original coordinate system (where x is pointing to First Point of Aries)
// so that x is pointing from star to the sun.
glm::dmat3 computeSystemRotation(glm::dvec3 starPosition);
/**
* Compute the inner and outer boundary of the habitable zone of a star, accordring to
* formula and coefficients by Kopparapu et al. (2015) https://arxiv.org/abs/1404.5292
*
* \param teff The effective temperature of the star, in Kelvin
* \param luminosity The luminosity of the star, in solar luminosities
* \return A vec2 with the lower and upper boundary in atronomical units, if a habitable
zone could be computed. Otherwise an std::nullopt
*/
std::optional<glm::vec2> computeHabitableZone(float teff, float luminosity);
// Create an identifier without whitespaces
std::string createIdentifier(std::string name);
modules/exoplanets/exoplanetsmodule_lua.inl
+73 -15
View File
@@ -55,6 +55,12 @@ constexpr const char* NoDataTextureFile =
"${SYNC}/http/exoplanets_textures/1/grid-32.png";
constexpr const char* DiscTextureFile =
"${SYNC}/http/exoplanets_textures/1/disc_texture.png";
constexpr const char* HabitableZoneTextureFile =
"${SYNC}/http/exoplanets_textures/1/hz_disc_texture.png";
const float AstronomicalUnit = static_cast<float>(distanceconstants::AstronomicalUnit);
const float SolarRadius = static_cast<float>(distanceconstants::SolarRadius);
const float JupiterRadius = static_cast<float>(distanceconstants::JupiterRadius);
ExoplanetSystem findExoplanetSystemInData(std::string_view starName) {
ExoplanetSystem system;
@@ -164,10 +170,9 @@ void createExoplanetSystem(const std::string& starName) {
const glm::dvec3 starPos =
static_cast<glm::dvec3>(starPosInParsec) * distanceconstants::Parsec;
const glm::dmat3 exoplanetSystemRotation = computeSystemRotation(starPos);
const float solarRadius = static_cast<float>(distanceconstants::SolarRadius);
// Star
float radiusInMeter = solarRadius;
float radiusInMeter = SolarRadius;
if (!std::isnan(system.starData.radius)) {
radiusInMeter *= system.starData.radius;
}
@@ -242,7 +247,7 @@ void createExoplanetSystem(const std::string& starName) {
// Planets
for (size_t i = 0; i < system.planetNames.size(); i++) {
ExoplanetDataEntry planet = system.planetsData[i];
ExoplanetDataEntry& planet = system.planetsData[i];
const std::string planetName = system.planetNames[i];
if (std::isnan(planet.ecc)) {
@@ -271,28 +276,24 @@ void createExoplanetSystem(const std::string& starName) {
sEpoch = "2009-05-19T07:11:34.080";
}
const float astronomicalUnit =
static_cast<float>(distanceconstants::AstronomicalUnit);
const float jupiterRadius = static_cast<float>(distanceconstants::JupiterRadius);
float planetRadius;
std::string enabled;
if (std::isnan(planet.r)) {
if (std::isnan(planet.rStar)) {
planetRadius = planet.a * 0.001f * astronomicalUnit;
planetRadius = planet.a * 0.001f * AstronomicalUnit;
}
else {
planetRadius = planet.rStar * 0.1f * solarRadius;
planetRadius = planet.rStar * 0.1f * SolarRadius;
}
enabled = "false";
}
else {
planetRadius = static_cast<float>(planet.r) * jupiterRadius;
planetRadius = static_cast<float>(planet.r) * JupiterRadius;
enabled = "true";
}
const float periodInSeconds = static_cast<float>(planet.per * SecondsPerDay);
const float semiMajorAxisInMeter = planet.a * astronomicalUnit;
const float semiMajorAxisInMeter = planet.a * AstronomicalUnit;
const float semiMajorAxisInKm = semiMajorAxisInMeter * 0.001f;
const std::string planetIdentifier = createIdentifier(planetName);
@@ -365,13 +366,12 @@ void createExoplanetSystem(const std::string& starName) {
bool hasLowerAUncertainty = !std::isnan(planet.aLower);
if (hasUpperAUncertainty && hasLowerAUncertainty) {
// Get the orbit plane of the planet trail orbit from the KeplerTranslation
const glm::dmat4 orbitPlaneRotationMatrix = computeOrbitPlaneRotationMatrix(
const glm::dmat4 rotation = computeOrbitPlaneRotationMatrix(
planet.i,
planet.bigOmega,
planet.omega
);
const glm::dmat3 rotation = orbitPlaneRotationMatrix;
const glm::dmat3 rotationMat3 = static_cast<glm::dmat3>(rotation);
const std::string discNode = "{"
"Identifier = '" + planetIdentifier + "_Disc',"
@@ -391,7 +391,7 @@ void createExoplanetSystem(const std::string& starName) {
"Transform = {"
"Rotation = {"
"Type = 'StaticRotation',"
"Rotation = " + ghoul::to_string(rotation) + ""
"Rotation = " + ghoul::to_string(rotationMat3) + ""
"}"
"},"
"GUI = {"
@@ -406,6 +406,64 @@ void createExoplanetSystem(const std::string& starName) {
);
}
}
// Habitable Zone
bool hasTeff = !std::isnan(system.starData.teff);
bool hasLuminosity = !std::isnan(system.starData.luminosity);
std::optional<glm::vec2> zone = std::nullopt;
if (hasTeff && hasLuminosity) {
zone = computeHabitableZone(system.starData.teff, system.starData.luminosity);
}
if (zone.has_value()) {
float meanInclination = 0.f;
for (const ExoplanetDataEntry& p : system.planetsData) {
meanInclination += p.i;
}
meanInclination /= static_cast<float>(system.planetsData.size());
const glm::dmat4 rotation = computeOrbitPlaneRotationMatrix(meanInclination);
const glm::dmat3 rotationMat3 = static_cast<glm::dmat3>(rotation);
glm::vec2 limits = zone.value();
float half = 0.5f * (limits[1] - limits[0]);
float center = limits[0] + half;
const std::string zoneDiscNode = "{"
"Identifier = '" + starIdentifier + "_HZ_Disc',"
"Parent = '" + starIdentifier + "',"
"Enabled = true,"
"Renderable = {"
"Type = 'RenderableOrbitDisc',"
"Rotation = {"
"Type = 'StaticRotation',"
"Rotation = " + ghoul::to_string(rotationMat3) + ""
"},"
"Texture = openspace.absPath('" + HabitableZoneTextureFile + "'),"
"Size = " + std::to_string(center * AstronomicalUnit) + ","
"Eccentricity = 0,"
"Offset = { " +
std::to_string(half) + ", " +
std::to_string(half) +
"}," //min / max extend
"Opacity = 0.05"
"},"
"Transform = {"
"Rotation = {"
"Type = 'StaticRotation',"
"Rotation = " + ghoul::to_string(rotationMat3) + ""
"}"
"},"
"GUI = {"
"Name = '" + starName + " Habitable Zone',"
"Path = '" + guiPath + "'"
"}"
"}";
openspace::global::scriptEngine->queueScript(
"openspace.addSceneGraphNode(" + zoneDiscNode + ");",
scripting::ScriptEngine::RemoteScripting::Yes
);
}
}
int addExoplanetSystem(lua_State* L) {
modules/exoplanets/tasks/exoplanetsdatapreparationtask.cpp
+6 -3
View File
@@ -298,13 +298,16 @@ void ExoplanetsDataPreparationTask::perform(
}
// Star luminosity
else if (column == "st_lum") {
p.luminosity = readFloatData(data);
float dataInLogSolar = readFloatData(data);
p.luminosity = std::pow(10, dataInLogSolar);
}
else if (column == "st_lumerr1") {
p.luminosityUpper = readFloatData(data);
float dataInLogSolar = readFloatData(data);
p.luminosityUpper = std::pow(10, dataInLogSolar);
}
else if (column == "st_lumerr2") {
p.luminosityLower = -readFloatData(data);
float dataInLogSolar = readFloatData(data);
p.luminosityLower = -std::pow(10, dataInLogSolar);
}
// Is the planet orbiting a binary system?
else if (column == "cb_flag") {