Merge pull request Tudat#50 from tudat-team/develop

Develop
DominicDirkx · Jan 24, 2022 · c28f2a3 · c28f2a3
2 parents dbbd353 + f31e986
commit c28f2a3
Show file tree

Hide file tree

Showing 19 changed files with 404 additions and 128 deletions.
diff --git a/include/tudat/astro/basic_astro/customAccelerationModel.h b/include/tudat/astro/basic_astro/customAccelerationModel.h
@@ -25,7 +25,6 @@ class CustomAccelerationModel: public basic_astrodynamics::AccelerationModel3d
             const std::function< Eigen::Vector3d( const double ) > accelerationFunction ):
         accelerationFunction_( accelerationFunction )
     {
-        updateMembers( TUDAT_NAN );
     }
 
     virtual void updateMembers( const double currentTime = TUDAT_NAN )

diff --git a/include/tudat/astro/electromagnetism/cannonBallRadiationPressureAcceleration.h b/include/tudat/astro/electromagnetism/cannonBallRadiationPressureAcceleration.h
@@ -105,7 +105,6 @@ class CannonBallRadiationPressureAcceleration: public basic_astrodynamics::Accel
           areaFunction_( areaFunction ),
           massFunction_( massFunction )
     {
-        this->updateMembers( );
     }
 
     //! Constructor taking functions pointers and constant values for parameters.
@@ -135,7 +134,6 @@ class CannonBallRadiationPressureAcceleration: public basic_astrodynamics::Accel
           areaFunction_( [ = ]( ){ return area; } ),
           massFunction_( [ = ]( ){ return mass; } )
     {
-        this->updateMembers( );
     }
 
     //! Update member variables used by the radiation pressure acceleration model.

diff --git a/include/tudat/astro/electromagnetism/solarSailAcceleration.h b/include/tudat/astro/electromagnetism/solarSailAcceleration.h
@@ -145,7 +145,6 @@ class SolarSailAcceleration: public basic_astrodynamics::AccelerationModel3d
           areaFunction_( areaFunction ),
           massFunction_( massFunction )
     {
-        this->updateMembers( );
     }
 
     //! Constructor taking functions pointers and constant values for parameters.
@@ -201,7 +200,6 @@ class SolarSailAcceleration: public basic_astrodynamics::AccelerationModel3d
           areaFunction_( [ = ]( ){ return area;} ),
           massFunction_( [ = ]( ){ return mass;} )
     {
-        this->updateMembers( );
     }
 
 
@@ -232,7 +230,6 @@ class SolarSailAcceleration: public basic_astrodynamics::AccelerationModel3d
         areaFunction_( std::bind( &RadiationPressureInterface::getArea, radiationPressureInterface ) ),
         massFunction_( massFunction )
     {
-        this->updateMembers( );
     }
 
 

diff --git a/include/tudat/astro/ephemerides/compositeEphemeris.h b/include/tudat/astro/ephemerides/compositeEphemeris.h
@@ -358,6 +358,16 @@ std::shared_ptr< Ephemeris > createReferencePointEphemeris(
         std::shared_ptr< RotationalEphemeris > bodyRotationModel,
         std::function< Eigen::Vector6d( const double& ) > referencePointRelativeStateFunction )
 {
+    if( bodyEphemeris == nullptr )
+    {
+        throw std::runtime_error( "Error when creating reference point composite ephemeris, no body ephemeris is provided" );
+    }
+
+    if( bodyRotationModel == nullptr )
+    {
+        throw std::runtime_error( "Error when creating reference point composite ephemeris, no body rotation model is provided" );
+    }
+
     typedef Eigen::Matrix< StateScalarType, 6, 1 > StateType;
 
     // Cast state fucntion of body (global) and reference point (local) into correct form.

diff --git a/include/tudat/astro/gravitation/centralGravityModel.h b/include/tudat/astro/gravitation/centralGravityModel.h
@@ -175,7 +175,6 @@ class CentralGravitationalAccelerationModel
                 positionOfBodyExertingAccelerationFunction,
                 isMutualAttractionUsed )
     {
-        this->updateMembers( );
     }
 
     //! Constructor taking position-functions for bodies, and constant gravitational parameter.
@@ -209,7 +208,6 @@ class CentralGravitationalAccelerationModel
                 positionOfBodyExertingAccelerationFunction,
                 isMutualAttractionUsed )
     {
-        this->updateMembers( );
     }
 
     //! Update members.

diff --git a/include/tudat/astro/gravitation/sphericalHarmonicsGravityModel.h b/include/tudat/astro/gravitation/sphericalHarmonicsGravityModel.h
@@ -228,7 +228,6 @@ class SphericalHarmonicsGravitationalAccelerationModel
                                      sphericalHarmonicsCache_->getMaximumDegree( ) ),
                     std::max< int >( maximumOrder_,
                                      sphericalHarmonicsCache_->getMaximumOrder( ) ) + 1 );
-        this->updateMembers( );
     }
 
     //! Constructor taking functions for position of bodies, and parameters of spherical harmonics
@@ -292,7 +291,6 @@ class SphericalHarmonicsGravitationalAccelerationModel
                                      sphericalHarmonicsCache_->getMaximumOrder( ) ) + 1 );
 
 
-        this->updateMembers( );
     }
 
     //! Get gravitational acceleration in body-fixed frame of body undergoing acceleration.

diff --git a/include/tudat/astro/orbit_determination/estimatable_parameters/estimatableParameter.h b/include/tudat/astro/orbit_determination/estimatable_parameters/estimatableParameter.h
@@ -225,6 +225,9 @@ class EstimatableParameter
         return Eigen::VectorXd::Zero( 0 );
     }
 
+    virtual void throwExceptionIfNotFullyDefined( ){ }
+
+
 
 protected:
 

diff --git a/include/tudat/astro/orbit_determination/estimatable_parameters/observationBiasParameter.h b/include/tudat/astro/orbit_determination/estimatable_parameters/observationBiasParameter.h
@@ -83,7 +83,14 @@ class ConstantObservationBiasParameter: public EstimatableParameter< Eigen::Vect
      */
     void setParameterValue( Eigen::VectorXd parameterValue )
     {
-        resetCurrentBias_( parameterValue );
+        if( resetCurrentBias_ != nullptr )
+        {
+            resetCurrentBias_( parameterValue );
+        }
+        else
+        {
+            throwExceptionIfNotFullyDefined( );
+        }
     }
 
     //! Function to retrieve the size of the parameter (equal to the size of the observable).
@@ -117,6 +124,20 @@ class ConstantObservationBiasParameter: public EstimatableParameter< Eigen::Vect
         resetCurrentBias_ = resetCurrentBias;
     }
 
+    void throwExceptionIfNotFullyDefined( )
+    {
+        if( getCurrentBias_ == nullptr || resetCurrentBias_ == nullptr )
+        {
+            throw std::runtime_error(
+                        "Error in " + getParameterTypeString( parameterName_.first ) +
+                        " of observable type " + observation_models::getObservableName(
+                            observableType_, linkEnds_.size( ) ) +
+                        " with link ends: " + observation_models::getLinkEndsString( linkEnds_ ) +
+                        " parameter not linked to bias object. Associated bias model been implemented in observation model. " +
+                        " This may be because you are resetting the parameter value before creating observation models, or because you have not defined the required bias model.");
+        }
+    }
+
     //! Function to retrieve the observation link ends for which the bias is active.
     /*!
      * Function to retrieve the observation link ends for which the bias is active.
@@ -242,12 +263,19 @@ class ArcWiseObservationBiasParameter: public EstimatableParameter< Eigen::Vecto
     {
         std::vector< Eigen::VectorXd > observationBiases;
 
-        for( int i = 0; i < numberOfArcs_; i++ )
+        if( resetBiasList_ != nullptr )
         {
-           observationBiases.push_back( parameterValue.segment( i * observableSize_, observableSize_ ) );
-        }
+            for( int i = 0; i < numberOfArcs_; i++ )
+            {
+                observationBiases.push_back( parameterValue.segment( i * observableSize_, observableSize_ ) );
+            }
 
-        resetBiasList_( observationBiases );
+            resetBiasList_( observationBiases );
+        }
+        else
+        {
+            throwExceptionIfNotFullyDefined( );
+        }
     }
 
     //! Function to retrieve the size of the parameter (equal to the size of the observable).
@@ -348,6 +376,19 @@ class ArcWiseObservationBiasParameter: public EstimatableParameter< Eigen::Vecto
         lookupScheme_ = lookupScheme;
     }
 
+    void throwExceptionIfNotFullyDefined( )
+    {
+        if( getBiasList_ == nullptr || resetBiasList_ == nullptr )
+        {
+            throw std::runtime_error(
+                        "Error in " + getParameterTypeString( parameterName_.first ) +
+                        " of observable type " + observation_models::getObservableName(
+                            observableType_, linkEnds_.size( ) ) +
+                        " with link ends: " + observation_models::getLinkEndsString( linkEnds_ ) +
+                        " parameter not linked to bias object. Has associated bias model been implemented in observation model?");
+        }
+    }
+
 protected:
 
 private:

diff --git a/include/tudat/astro/orbit_determination/podInputOutputTypes.h b/include/tudat/astro/orbit_determination/podInputOutputTypes.h
@@ -428,6 +428,18 @@ class EstimationConvergenceChecker
     unsigned int numberOfIterationsWithoutImprovement_;
 };
 
+
+inline std::shared_ptr< EstimationConvergenceChecker > estimationConvergenceChecker(
+        const unsigned int maximumNumberOfIterations = 5,
+        const double minimumResidualChange = 0.0,
+        const double minimumResidual = 1.0E-20,
+        const int numberOfIterationsWithoutImprovement = 2 )
+{
+    return std::make_shared< EstimationConvergenceChecker >(
+            maximumNumberOfIterations, minimumResidualChange, minimumResidual, numberOfIterationsWithoutImprovement );
+}
+
+
 void scaleInformationMatrixWithWeights(
         Eigen::MatrixXd& informationMatrix,
         const Eigen::VectorXd& weightsDiagonal );