True anomaly works

This commit is contained in:
Nick 2022-04-10 22:54:55 -04:00
parent 71971d0522
commit 967ea5e31e
2 changed files with 73 additions and 26 deletions

View File

@ -17,24 +17,19 @@ namespace NewHorizons.Builder.Orbital
{
public static InitialMotion Make(GameObject body, AstroObject primaryBody, AstroObject secondaryBody, OWRigidbody OWRB, OrbitModule orbit)
{
body.SetActive(false);
InitialMotion initialMotion = body.AddComponent<InitialMotion>();
return Update(initialMotion, body, primaryBody, secondaryBody, OWRB, orbit);
}
public static float SiderealPeriodToAngularSpeed(float siderealPeriod)
{
return siderealPeriod == 0 ? 0f : 2f * Mathf.PI / (siderealPeriod * 60f);
}
// This bit makes the initial motion not try to calculate the orbit velocity itself for reasons
initialMotion._orbitImpulseScalar = 0f;
public static InitialMotion Update(InitialMotion initialMotion, GameObject body, AstroObject primaryBody, AstroObject secondaryBody, OWRigidbody OWRB, OrbitModule orbit)
{
// Rotation
initialMotion._initAngularSpeed = SiderealPeriodToAngularSpeed(orbit.SiderealPeriod);
initialMotion._initAngularSpeed = orbit.SiderealPeriod == 0 ? 0f : 2f * Mathf.PI / (orbit.SiderealPeriod * 60f);
//initialMotion._primaryBody = primaryBody?.GetAttachedOWRigidbody();
var rotationAxis = Quaternion.AngleAxis(orbit.AxialTilt + 90f, Vector3.right) * Vector3.up;
body.transform.rotation = Quaternion.FromToRotation(Vector3.up, rotationAxis);
initialMotion._rotationAxis = rotationAxis;
initialMotion._orbitImpulseScalar = 0f;
if (!orbit.IsStatic && primaryBody != null)
{
initialMotion.SetPrimaryBody(primaryBody.GetAttachedOWRigidbody());
@ -45,10 +40,18 @@ namespace NewHorizons.Builder.Orbital
var secondaryGravity = new Gravity(secondaryBody.GetGravityVolume());
var velocity = orbit.GetOrbitalParameters(primaryGravity, secondaryGravity).InitialVelocity;
initialMotion._initLinearDirection = velocity.normalized;
initialMotion._initLinearSpeed = velocity.magnitude;
initialMotion._cachedInitVelocity = velocity + primaryBody?.GetComponent<InitialMotion>()?.GetInitVelocity() ?? Vector3.zero;
initialMotion._isInitVelocityDirty = false;
}
}
else
{
initialMotion._initLinearDirection = Vector3.forward;
initialMotion._initLinearSpeed = 0f;
}
body.SetActive(true);
return initialMotion;
}

View File

@ -4,6 +4,7 @@ using System.Linq;
using System.Text;
using System.Threading.Tasks;
using UnityEngine;
using Logger = NewHorizons.Utility.Logger;
namespace NewHorizons.Components.Orbital
{
@ -41,35 +42,78 @@ namespace NewHorizons.Components.Orbital
var power = primaryGravity.Power;
var period = (float) (GravityVolume.GRAVITATIONAL_CONSTANT * (primaryMass + secondaryMass) / (4 * Math.PI * Math.PI * Math.Pow(semiMajorAxis, power)));
var f = Mathf.Deg2Rad * trueAnomaly; // True anomaly in radians
var eccentricAnomaly = Mathf.Rad2Deg * Mathf.Atan2(Mathf.Sqrt(1 - eccentricity * eccentricity) * Mathf.Sin(f), 1 + eccentricity * Mathf.Cos(f));
var meanAnomaly = eccentricAnomaly - eccentricity * Mathf.Sin(Mathf.Deg2Rad * eccentricAnomaly);
// All in radians
var f = Mathf.Deg2Rad * trueAnomaly;
var eccentricAnomaly = Mathf.Atan2(Mathf.Sqrt(1 - (eccentricity * eccentricity)) * Mathf.Sin(f), eccentricity + Mathf.Cos(f));
var meanAnomaly = eccentricAnomaly - eccentricity * Mathf.Sin(eccentricAnomaly);
orbitalParameters.Period = period;
orbitalParameters.MeanAnomaly = meanAnomaly;
orbitalParameters.EccentricAnomaly = eccentricAnomaly;
orbitalParameters.MeanAnomaly = meanAnomaly * Mathf.Rad2Deg;
orbitalParameters.EccentricAnomaly = eccentricAnomaly * Mathf.Rad2Deg;
// Position
var ap = Mathf.Deg2Rad * argumentOfPeriapsis; // Argument of periapsis in radians
var la = Mathf.Deg2Rad * longitudeOfAscendingNode; // Longitude of ascending node in radians
var i = Mathf.Deg2Rad * inclination; // Inclination in radians
var p = semiMajorAxis * (1 - eccentricity * eccentricity); // Semi-latus rectum
var p = semiMajorAxis * (1 - (eccentricity * eccentricity)); // Semi-latus rectum
var r = p / (1 + eccentricity * Mathf.Cos(f));
var mu = GravityVolume.GRAVITATIONAL_CONSTANT * primaryMass;
var h = Mathf.Sqrt(mu * p);
var v = Mathf.Sqrt(mu * ((2 / r) - (1 / semiMajorAxis)));
var semiMinorAxis = semiMajorAxis * Mathf.Sqrt(1 - (eccentricity * eccentricity));
var focusDistance = Mathf.Sqrt((semiMajorAxis * semiMajorAxis) - (semiMinorAxis * semiMinorAxis));
var x = semiMajorAxis * Mathf.Cos(eccentricAnomaly) - focusDistance;
var y = semiMinorAxis * Mathf.Sin(eccentricAnomaly);
var R1 = Quaternion.AngleAxis(longitudeOfAscendingNode, Vector3.up);
var R2 = Quaternion.AngleAxis(inclination, Vector3.forward);
var R3 = Quaternion.AngleAxis(argumentOfPeriapsis, Vector3.up);
var n_p = new Vector2(x, y).normalized;
var n_v = new Vector2(-y, x).normalized;
// Idk why y and z get swapped just go with it
var pos = new Vector3(r * n_p.x, r * n_p.y, 0f);
var vel = new Vector3(v * n_v.x, 0f, -v * n_v.y);
orbitalParameters.InitialPosition = R1 * R2 * R3 * pos;
orbitalParameters.InitialVelocity = R1 * R2 * R3 * vel;
Logger.Log($"POSITION: {orbitalParameters.InitialPosition}, {orbitalParameters.InitialVelocity}, {n_p}, {n_v}");
/*
var x = r * (Mathf.Cos(la) * Mathf.Cos(ap + f) - Mathf.Sin(la) * Mathf.Sin(ap + f) * Mathf.Cos(i));
var y = r * (Mathf.Sin(la) * Mathf.Cos(ap + f) - Mathf.Cos(la) * Mathf.Sin(ap + f) * Mathf.Cos(i));
var y = r * (Mathf.Sin(la) * Mathf.Cos(ap + f) + Mathf.Cos(la) * Mathf.Sin(ap + f) * Mathf.Cos(i));
var z = r * (Mathf.Sin(i) * Mathf.Sin(ap + f));
orbitalParameters.InitialPosition = new Vector3(x, y, z);
// Velocity
var coefficient = h * eccentricity * Mathf.Sin(f) / (r * p);
var vx = x * coefficient - (h / r) * (Mathf.Cos(la) * Mathf.Sin(ap + f) + Mathf.Sin(la) * Mathf.Cos(ap + f) * Mathf.Cos(i));
var vy = y * coefficient - (h / r) * (Mathf.Sin(la) * Mathf.Sin(ap + f) - Mathf.Cos(la) * Mathf.Cos(ap + f) * Mathf.Cos(i));
var vz = z * coefficient + (h / r) * (Mathf.Sin(i) * Mathf.Cos(ap + f));
orbitalParameters.InitialVelocity = new Vector3(vx, vy, vz);
*/
/*
// Where x points towards the periapsis
var ox = semiMajorAxis * Mathf.Cos(f);
var oy = semiMinorAxis * Mathf.Sin(f);
var o = new Vector3(ox, 0, oy);
var mu = GravityVolume.GRAVITATIONAL_CONSTANT * primaryMass;
var r = o.magnitude;
var o_dot_coeff = Mathf.Sqrt(mu * ((2/r) - (1/semiMajorAxis)));
var ox_dot = o_dot_coeff * -Mathf.Sin(eccentricAnomaly);
var oy_dot = o_dot_coeff * Mathf.Cos(eccentricAnomaly);
var o_dot = new Vector3(ox_dot, 0, oy_dot);
// Do some rotations
var R1 = Quaternion.AngleAxis(Mathf.Rad2Deg * -longitudeOfAscendingNode, Vector3.up);
var R2 = Quaternion.AngleAxis(Mathf.Rad2Deg * -inclination, Vector3.forward);
var R3 = Quaternion.AngleAxis(Mathf.Rad2Deg * -argumentOfPeriapsis, Vector3.up);
orbitalParameters.InitialPosition = R1 * R2 * R3 * o;
orbitalParameters.InitialVelocity = R1 * R2 * R3 * o_dot;
*/
return orbitalParameters;
}