@@ -38,6 +38,13 @@ namespace astro
38
38
namespace trajectory
39
39
{
40
40
41
+ // //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
42
+
43
+ using library::physics::units::Length ;
44
+ using library::physics::units::Derived ;
45
+
46
+ static const Derived::Unit GravitationalParameterSIUnit = Derived::Unit::GravitationalParameter(Length::Unit::Meter, library::physics::units::Time::Unit::Second) ;
47
+
41
48
// //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
42
49
43
50
Orbit::Orbit ( const orbit::Model& aModel,
@@ -140,7 +147,7 @@ Pass Orbit::getPassAt (
140
147
}
141
148
142
149
return this ->getPassWithRevolutionNumber (this ->getRevolutionNumberAt (anInstant)) ;
143
-
150
+
144
151
}
145
152
146
153
Pass Orbit::getPassWithRevolutionNumber ( const Integer& aRevolutionNumber ) const
@@ -181,7 +188,7 @@ Pass Orbit::getPassWithRevolutionNumber (
181
188
182
189
const auto lowerBoundMapIt = passMap_.lower_bound (aRevolutionNumber) ;
183
190
184
- if (lowerBoundMapIt != passMap_.end ())
191
+ if (lowerBoundMapIt != passMap_.end ())
185
192
{
186
193
// std::cout << "lowerBoundMap = " << lowerBoundMapIt->second << std::endl ;
187
194
@@ -202,7 +209,7 @@ Pass Orbit::getPassWithRevolutionNumber (
202
209
{
203
210
closestPassPtr = &(lowerBoundMapIt->second ) ;
204
211
}
205
-
212
+
206
213
}
207
214
208
215
}
@@ -284,7 +291,7 @@ Pass Orbit::getPassWithRevolutionNumber (
284
291
{
285
292
stepDuration = Duration::Seconds (-currentStateCoordinates_ECI_z * stepDuration.inSeconds () / (currentStateCoordinates_ECI_z - previousStateCoordinates_ECI_z)) ;
286
293
}
287
-
294
+
288
295
if (stepDuration.isZero ())
289
296
{
290
297
stepDuration = (currentStateCoordinates_ECI_z < 0.0 ) ? Duration::Nanoseconds (+1.0 ) : Duration::Nanoseconds (-1.0 ) ;
@@ -356,7 +363,7 @@ Pass Orbit::getPassWithRevolutionNumber (
356
363
}
357
364
358
365
return Pass::Undefined () ;
359
-
366
+
360
367
}
361
368
362
369
Shared<const Frame> Orbit::getOrbitalFrame ( const Orbit::FrameType& aFrameType ) const
@@ -421,14 +428,14 @@ Shared<const Frame> Orbit::getOrbitalFrame (
421
428
break ;
422
429
423
430
}
424
-
431
+
425
432
case Orbit::FrameType::LVLH:
426
433
{
427
434
428
435
// X axis along position vector
429
436
// Z axis along orbital momentum
430
437
// Y axis toward velocity vector
431
-
438
+
432
439
const Shared<const DynamicProvider> dynamicProviderSPtr = std::make_shared<const DynamicProvider>
433
440
(
434
441
[this ] (const Instant& anInstant) -> Transform // [TBI] Use shared_from_this instead
@@ -456,7 +463,7 @@ Shared<const Frame> Orbit::getOrbitalFrame (
456
463
orbitalFrameSPtr = Frame::Construct (frameName, false , Frame::GCRF (), dynamicProviderSPtr) ;
457
464
458
465
break ;
459
-
466
+
460
467
}
461
468
462
469
case Orbit::FrameType::VVLH:
@@ -465,7 +472,7 @@ Shared<const Frame> Orbit::getOrbitalFrame (
465
472
// Z axis along negative position vector
466
473
// Y axis along negative orbital momentum
467
474
// X axis toward velocity vector
468
-
475
+
469
476
const Shared<const DynamicProvider> dynamicProviderSPtr = std::make_shared<const DynamicProvider>
470
477
(
471
478
[this ] (const Instant& anInstant) -> Transform // [TBI] Use shared_from_this instead
@@ -493,7 +500,7 @@ Shared<const Frame> Orbit::getOrbitalFrame (
493
500
orbitalFrameSPtr = Frame::Construct (frameName, false , Frame::GCRF (), dynamicProviderSPtr) ;
494
501
495
502
break ;
496
-
503
+
497
504
}
498
505
499
506
case Orbit::FrameType::QSW:
@@ -571,7 +578,7 @@ Shared<const Frame> Orbit::getOrbitalFrame (
571
578
case Orbit::FrameType::VNC:
572
579
{
573
580
574
- // X axis along velocity vector
581
+ // X axis along velocity vector
575
582
// Y axis along orbital momentum
576
583
577
584
const Shared<const DynamicProvider> dynamicProviderSPtr = std::make_shared<const DynamicProvider>
@@ -764,9 +771,9 @@ Orbit Orbit::SunSynchronous (
764
771
{
765
772
766
773
using library::core::types::Uint8 ;
767
-
774
+
768
775
using library::math::obj::Vector3d ;
769
-
776
+
770
777
using library::physics::units::Mass ;
771
778
using library::physics::units::Derived ;
772
779
using library::physics::time::Scale ;
@@ -802,7 +809,7 @@ Orbit Orbit::SunSynchronous (
802
809
803
810
const Real a = aSemiMajorAxis.inMeters () ;
804
811
const Real R = aCelestialObjectSPtr->getEquatorialRadius ().inMeters () ;
805
- const Real mu = aCelestialObjectSPtr->getGravitationalParameter ().in ({ Length::Unit::Meter, Derived::Order ( 3 ), Mass::Unit::Undefined, Derived::Order::Zero (), library::physics::units::Time::Unit::Second, Derived::Order (- 2 ), Angle::Unit::Undefined, Derived::Order::Zero () } ) ;
812
+ const Real mu = aCelestialObjectSPtr->getGravitationalParameter ().in (GravitationalParameterSIUnit ) ;
806
813
const Real j2 = aCelestialObjectSPtr->getJ2 () ;
807
814
808
815
const Real T_sid = 31558149.504 ; // [s] Sidereal year
@@ -838,7 +845,7 @@ Orbit Orbit::SunSynchronous (
838
845
const Real sunEclipticLatitude_rad = Angle::Degrees (std::fmod (sunMeanLongitude_deg + 1.914666471 * std::sin (sunMeanAnomaly_rad) + 0.019994643 * std::sin (2.0 * sunMeanAnomaly_rad), 360.0 )).inRadians () ;
839
846
840
847
// Compute the equation of time
841
-
848
+
842
849
const Real equationOfTime_deg = - 1.914666471 * std::sin (sunMeanAnomaly_rad)
843
850
- 0.019994643 * std::sin (2.0 * sunMeanAnomaly_rad)
844
851
+ 2.466 * std::sin (2.0 * sunEclipticLatitude_rad)
@@ -865,30 +872,30 @@ Orbit Orbit::SunSynchronous (
865
872
environment.setInstant (anEpoch) ;
866
873
867
874
// Sun direction in GCRF
868
-
875
+
869
876
const Vector3d sunDirection_GCRF = environment.accessCelestialObjectWithName (" Sun" )->getPositionIn (Frame::GCRF ()).getCoordinates ().normalized () ;
870
877
871
878
// Desired angle between the Sun and the ascending node
872
-
879
+
873
880
const Angle alpha = Angle::Degrees ((localTime - 12.0 ) / 12.0 * 180.0 ) ;
874
881
875
882
// Sun Apparent Local Time (right ascension of the Sun in GCRF)
876
883
// https://en.wikipedia.org/wiki/Solar_time#Apparent_solar_time
877
-
884
+
878
885
const Angle apparentSolarTime = Angle::Radians (std::atan2 (sunDirection_GCRF.y (), sunDirection_GCRF.x ())) ;
879
886
880
887
// Equation of Time
881
888
// https://en.wikipedia.org/wiki/Equation_of_time
882
-
889
+
883
890
const Angle equationOfTime = calculateEquationOfTime (anEpoch) ;
884
891
885
892
// Sun Mean Local Time
886
893
// https://en.wikipedia.org/wiki/Solar_time#Mean_solar_time
887
894
888
895
const Angle meanSolarTime = apparentSolarTime + equationOfTime ;
889
-
896
+
890
897
// Right Ascension of the Ascending Node
891
-
898
+
892
899
const Angle raan = Angle::Radians (std::fmod (meanSolarTime.inRadians () + alpha.inRadians (), Real::TwoPi ())) ;
893
900
894
901
return raan ;
@@ -1077,4 +1084,4 @@ String Orbit::StringFromFrameType (
1077
1084
}
1078
1085
}
1079
1086
1080
- // //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
1087
+ // //////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
0 commit comments