Skip to content

Commit

Permalink
Merge pull request #2079 from eggrobin/geopotential-api
Browse files Browse the repository at this point in the history
Geopotential API
  • Loading branch information
eggrobin authored Feb 17, 2019
2 parents 60d9122 + 0db2459 commit 36fe638
Show file tree
Hide file tree
Showing 3 changed files with 195 additions and 1 deletion.
75 changes: 75 additions & 0 deletions ksp_plugin/interface_external.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@

#include <string>

#include "absl/strings/str_cat.h"
#include "base/array.hpp"
#include "base/status.hpp"
#include "base/status_or.hpp"
Expand All @@ -23,6 +24,7 @@ using ksp_plugin::Vessel;
using physics::BodyCentredNonRotatingDynamicFrame;
using physics::ComputeApsides;
using physics::DiscreteTrajectory;
using physics::OblateBody;
using physics::RigidMotion;
using physics::RigidTransformation;

Expand Down Expand Up @@ -71,6 +73,79 @@ Status principia__ExternalFlowFreefall(
"|ExternalFlowFreefall| is not yet implemented"));
}

Status principia__ExternalGeopotentialGetCoefficient(
Plugin const* const plugin,
int const body_index,
int const degree,
int const order,
XY* const coefficient) {
journal::Method<journal::ExternalGeopotentialGetCoefficient> m{
{plugin,
body_index,
degree,
order},
{coefficient}};
if (plugin == nullptr) {
return m.Return(
MakeStatus(Error::INVALID_ARGUMENT, "|plugin| must not be null"));
}
if (!plugin->HasCelestial(body_index)) {
return m.Return(MakeStatus(
Error::NOT_FOUND,
absl::StrCat("No celestial with index ", body_index)));
}
if (order < 0 || order > degree) {
return m.Return(MakeStatus(
Error::INVALID_ARGUMENT,
absl::StrCat(u8"Expected 0 ≤ order ≤ degree; got degree = ",
degree, ", order = ", order)));
}
if (degree == 0) {
*coefficient = {1, 0};
return m.Return(OK());
}
auto const& body = *plugin->GetCelestial(body_index).body();
if (!body.is_oblate()) {
*coefficient = {0, 0};
return m.Return(OK());
}
auto const& oblate_body = dynamic_cast<OblateBody<Barycentric> const&>(body);
if (degree > oblate_body.geopotential_degree()) {
*coefficient = {0, 0};
return m.Return(OK());
}
*coefficient = {oblate_body.cos()[degree][order],
oblate_body.sin()[degree][order]};
return m.Return(OK());
}

Status principia__ExternalGeopotentialGetReferenceRadius(
Plugin const* const plugin,
int const body_index,
double* const reference_radius) {
journal::Method<journal::ExternalGeopotentialGetReferenceRadius> m{
{plugin,
body_index},
{reference_radius}};
if (plugin == nullptr) {
return m.Return(
MakeStatus(Error::INVALID_ARGUMENT, "|plugin| must not be null"));
}
if (!plugin->HasCelestial(body_index)) {
return m.Return(MakeStatus(
Error::NOT_FOUND,
absl::StrCat("No celestial with index ", body_index)));
}
auto const& body = *plugin->GetCelestial(body_index).body();
if (!body.is_oblate()) {
*reference_radius = body.mean_radius() / Metre;
return m.Return(OK());
}
auto const& oblate_body = dynamic_cast<OblateBody<Barycentric> const&>(body);
*reference_radius = oblate_body.reference_radius() / Metre;
return m.Return(OK());
}

Status principia__ExternalGetNearestPlannedCoastDegreesOfFreedom(
Plugin const* const plugin,
int const central_body_index,
Expand Down
70 changes: 70 additions & 0 deletions ksp_plugin_test/interface_external_test.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,7 @@ using ksp_plugin::PartId;
using ksp_plugin::FakePlugin;
using ksp_plugin::Vessel;
using physics::SolarSystem;
using quantities::Sqrt;
using quantities::si::Centi;
using quantities::si::Hour;
using quantities::si::Kilo;
Expand Down Expand Up @@ -114,5 +115,74 @@ TEST_F(InterfaceExternalTest, GetNearestPlannedCoastDegreesOfFreedom) {
Lt(1 * Centi(Metre) / Second)))));
}

TEST_F(InterfaceExternalTest, Geopotential) {
XY coefficient;
double radius;
auto status = principia__ExternalGeopotentialGetCoefficient(
&plugin_,
SolarSystemFactory::Earth,
/*degree=*/2,
/*order=*/0,
&coefficient);
EXPECT_THAT(base::Status(static_cast<base::Error>(status.error), ""), IsOk());
EXPECT_THAT(-coefficient.x * Sqrt(5), IsNear(1.08e-3));
EXPECT_THAT(coefficient.y, Eq(0));

status = principia__ExternalGeopotentialGetCoefficient(
&plugin_,
SolarSystemFactory::Earth,
/*degree=*/3,
/*order=*/1,
&coefficient);
EXPECT_THAT(base::Status(static_cast<base::Error>(status.error), ""), IsOk());
EXPECT_THAT(coefficient.x, IsNear(2.03e-6));
EXPECT_THAT(coefficient.y, IsNear(0.248e-6));

status = principia__ExternalGeopotentialGetCoefficient(
&plugin_,
SolarSystemFactory::Earth,
/*degree=*/1729,
/*order=*/163,
&coefficient);
EXPECT_THAT(base::Status(static_cast<base::Error>(status.error), ""), IsOk());
EXPECT_THAT(coefficient.x, Eq(0));
EXPECT_THAT(coefficient.y, Eq(0));

status = principia__ExternalGeopotentialGetReferenceRadius(
&plugin_,
SolarSystemFactory::Earth,
&radius);
EXPECT_THAT(base::Status(static_cast<base::Error>(status.error), ""), IsOk());
EXPECT_THAT(radius, Eq(6'378'136.3));

status = principia__ExternalGeopotentialGetCoefficient(
&plugin_,
SolarSystemFactory::Ariel,
/*degree=*/2,
/*order=*/2,
&coefficient);
EXPECT_THAT(base::Status(static_cast<base::Error>(status.error), ""), IsOk());
EXPECT_THAT(coefficient.x, Eq(0));
EXPECT_THAT(coefficient.y, Eq(0));

status = principia__ExternalGeopotentialGetCoefficient(
&plugin_,
SolarSystemFactory::Ariel,
/*degree=*/0,
/*order=*/0,
&coefficient);
EXPECT_THAT(base::Status(static_cast<base::Error>(status.error), ""), IsOk());
EXPECT_THAT(coefficient.x, Eq(1));
EXPECT_THAT(coefficient.y, Eq(0));

status = principia__ExternalGeopotentialGetReferenceRadius(
&plugin_,
SolarSystemFactory::Ariel,
&radius);
EXPECT_THAT(base::Status(static_cast<base::Error>(status.error), ""), IsOk());
EXPECT_THAT(radius, Eq(578'900));

}

} // namespace interface
} // namespace principia
51 changes: 50 additions & 1 deletion serialization/journal.proto
Original file line number Diff line number Diff line change
Expand Up @@ -151,7 +151,7 @@ message XY {
}

message Method {
extensions 5000 to 5999; // Last used: 5155.
extensions 5000 to 5999; // Last used: 5157.
}

message AdvanceTime {
Expand Down Expand Up @@ -423,6 +423,55 @@ message ExternalFlowFreefall {
optional Return return = 3;
}

// Sets |coefficient| to the normalized geopotential coefficient of the given
// |degree| and |order| of the body with index |body_index|.
// |coefficient.x| is set to Cnm, |coefficient.y| is set to Snm.
message ExternalGeopotentialGetCoefficient {
extend Method {
optional ExternalGeopotentialGetCoefficient extension = 5156;
}
message In {
required fixed64 plugin = 1 [(pointer_to) = "Plugin const",
(is_subject) = true];
required int32 body_index = 2;
required int32 degree = 3;
required int32 order = 4;
}
message Out {
required XY coefficient = 1;
}
message Return {
required Status result = 1;
}
optional In in = 1;
optional Out out = 2;
optional Return return = 3;
}

// Sets |reference_radius| to the value in metres of the reference radius of
// the geopotential model for the body with index |body_index|.
// If the body does not have a geopotential model, i.e., if it is modeled as an
// isotropic sphere, the mean radius is used.
message ExternalGeopotentialGetReferenceRadius {
extend Method {
optional ExternalGeopotentialGetReferenceRadius extension = 5157;
}
message In {
required fixed64 plugin = 1 [(pointer_to) = "Plugin const",
(is_subject) = true];
required int32 body_index = 2;
}
message Out {
required double reference_radius = 1;
}
message Return {
required Status result = 1;
}
optional In in = 1;
optional Out out = 2;
optional Return return = 3;
}

// Returns the first point of the coast phase following the given manoeuvre
// which is locally nearest to the reference position, or the nearest endpoint
// if there is no local minimum.
Expand Down

0 comments on commit 36fe638

Please sign in to comment.