simpleOrbit.cpp

Go to the documentation of this file.

/************************************************************************************************/
/*! \file simpleOrbit.cpp
*  \brief Nonlinear Orbit EOM
*  \author $Author: jayhawk_hokie $
*  \version $Revision: 1.3 $
*  \date    $Date: 2006/08/25 19:14:17 $
************************************************************************************************/
/*!
*
************************************************************************************************/

// Standard System includes
#include <iostream.h>
#include <iomanip.h>
#include <fstream.h>

// Open-SESSAME Includes
#include "attitude/AttitudeModels/QuaternionAngVelDynamics.h" // Dynamics model
#include "utils/RungeKuttaIntegrator.h" // Integrator function
#include "datahandling/History.h" // attitude and orbit history

#include "orbit/orbitstaterep/Keplerian.h"


// XML Parser Include
#include <dsacssinterface.h>

using namespace O_SESSAME;

#define PI M_PI

// Specify Xml File to Parse
const char *fileName = "DSACSSConfig.xml";
// Define xml config path names
const char *spacecraft = "SPACECRAFT_1";
const char *physicalProperties = "PHYSICAL_PROPERTIES";
const char *attitudeProperties = "ATTITUDE_PROPERTIES";
const char *hardwareProperties = "HARDWARE_PROPERTIES";
const char *initialConditions = "INITIAL_CONDITIONS";
const char *magnetometer = "MAGNETOMETER";
const char *attribute = "value";


Vector NullFunctor(const ssfTime& _pSSFTime, const OrbitState& _pOrbitState,
const AttitudeState& _pAttitudeState)
{
        return Vector(3);
}

int main()
{

        /******************************************************************
        *  Load Configuration Parameters from xml file 
        *******************************************************************/
        // Declare TiXmlDocument
        TiXmlDocument config( fileName );
        // Load Xml File
        bool loadOkay = config.LoadFile();
        checkLoadFile(loadOkay, fileName, config);
        // Declare TiXmlHandle
        TiXmlHandle docHandle( &config );
        int response;

        /* Get mass from xml config file */
        double mass;
        response = docHandle.FirstChild( spacecraft ).FirstChild( physicalProperties ).
        Child( "MASS", 0 ).Element() -> QueryDoubleAttribute(attribute, &mass);
        checkResponse(response);

        /* Get inertia matrix from xml config file */
        Matrix inertia = simulatorInertia(docHandle, spacecraft);

        /* Get initial angular rate from xml config file */
        Vector angularRateBody(3); // wb_b
        for(int iterator=0; iterator<3; iterator++)
        {
                response = docHandle.FirstChild( spacecraft ).FirstChild( attitudeProperties ).
                FirstChild( initialConditions ).Child( "ANGULAR_RATE_BODY", iterator ).Element() -> 
                QueryDoubleAttribute( "value", &angularRateBody(VectorIndexBase + iterator) );
                checkResponse(response);
        }
        
        /* Get initial attitude from xml config file */
        Vector attitudeVector(4); // quaternion
        for(int iterator=0; iterator<4; iterator++)
        {
                response = docHandle.FirstChild( spacecraft ).FirstChild( attitudeProperties ).
                FirstChild( initialConditions ).Child( "QUATERNION", iterator ).Element() -> 
                QueryDoubleAttribute( "value", &attitudeVector(VectorIndexBase + iterator) );
                checkResponse(response);
        }
        Quaternion quaternion(attitudeVector);

        /* Get initial orbit information from xml config file */
        Vector orbitVector(6);
        Keplerian orbit;
        string str = docHandle.FirstChild( spacecraft ).FirstChild( "ORBIT_PROPERTIES" ).
        FirstChild( initialConditions ).Child( "KEPLERIAN", 0 ).Element() -> 
        Attribute( "selection" );
        if ("ON" == str) // Keplarian orbit parameters (a, e, i, OMEGA, omega, f)
        {
                for(int iterator=0; iterator<6; iterator++)
                {
                        response = docHandle.FirstChild( spacecraft ).FirstChild( "ORBIT_PROPERTIES" ).
                        FirstChild( initialConditions ).FirstChild( "KEPLERIAN" ).Child( "VECTOR", iterator ).
                        Element() -> QueryDoubleAttribute( "value", &orbitVector(VectorIndexBase + iterator) );
                        checkResponse(response);
                }
                // Convert from deg to rad
                for(int iterator=3; iterator<7; iterator++)
                {
                        orbitVector(iterator) = Deg2Rad(orbitVector(iterator));
                }
                orbit.SetKeplerianRepresentationTrueAnamoly(orbitVector); // store orbitVector in orbit

                cout << "\nSemimajor Axis (km): " << orbit.GetSemimajorAxis()
                << " Eccentricity (~): " << orbit.GetEccentricity()
                << " Inclination (rad): " << orbit.GetInclination()
                << " Long Ascending Node (rad): " << orbit.GetLongAscNode()
                << " Argument of Perigee (rad): " << orbit.GetArgPerigee()
                << " True Anomaly (rad): " << orbit.GetTrueAnomaly() 
                << endl << endl;
        }

        /* Get numerical integration step size (sec) */
        double integrationStepSize;
        response = docHandle.FirstChild( spacecraft ).FirstChild( "TIME" ).
        Child( "INTEGRATION_STEP_SIZE", 0 ).Element() -> QueryDoubleAttribute(attribute, &integrationStepSize);
        checkResponse(response);
        
        /* Get number of orbits before simulation ends. */
        double numberOrbits;
        response = docHandle.FirstChild( spacecraft ).FirstChild( "TIME" ).
        Child( "NUMBER_OF_ORBITS", 0 ).Element() -> QueryDoubleAttribute(attribute, &numberOrbits);
        checkResponse(response);


        /******************************************************************
        *  Numerical Simulation 
        *******************************************************************/

        /* Start of execution time */   
        tick(); // (sec)

        /* Create attitude object */
        AttitudeState spacecraft_1;

        /* Create */
        SpecificFunctor AttitudeForcesFunctor(&NullFunctor);

        /* Create integrator object */
        RungeKuttaIntegrator attitudeIntegrator;

        /* Set number of integration steps */
        double orbitPeriod = 2*PI/( orbit.GetMeanMotion() ); // (sec)
        int integrationSteps = orbitPeriod * numberOrbits / integrationStepSize;
        attitudeIntegrator.SetNumSteps( integrationSteps );

        /* Create a container class to hold integration times */
        vector<ssfTime> integrationTimes;
        /* Simulation Time (sec) */
        ssfTime begin(0);
        ssfTime end( begin + (orbitPeriod * numberOrbits) ); 
        /* set integration times */
        integrationTimes.push_back(begin);
        integrationTimes.push_back(end);

        DirectionCosineMatrix Rbo; // initialize Rbo
        Rbo.Set( quaternion ); // create DCM from quaternion

        Vector angularRateInertial(3);
        angularRateInertial = angularRateBody; 
        //angularRateInertial = angularRateBody  - ( orbit.GetMeanMotion() ) * Rbo( _(i+0,i+2),_(i+1,i+1) ); // wbi_bo

        /* Set attitude and angular rates */
        spacecraft_1.SetRotation(Rbo);
        spacecraft_1.SetAngularVelocity(angularRateInertial);

        Matrix history = attitudeIntegrator.Integrate(integrationTimes, &AttituteDynamics_QuaternionAngVel, spacecraft_1.GetState(),
        NULL, NULL, inertia, AttitudeForcesFunctor);

        /* Save output */
        ofstream outputFile;
        outputFile.open("DataFiles/GravGradNonLin.dat");
        outputFile << history;
        outputFile.close();

        /* Output execution time to screen */
        cout << "The numerical simulation excution time was " << tock() << " seconds. \n \n" << endl;

        return(0);
}


// Do not change the comments below - they will be added automatically by CVS
/*****************************************************************************
*       $Log: simpleOrbit.cpp,v $
*       Revision 1.3  2006/08/25 19:14:17  jayhawk_hokie
*       Modified keplerian representation.
*
*       Revision 1.2  2006/08/11 19:50:28  jayhawk_hokie
*       Added config parameters.
*
*       Revision 1.1  2006/08/08 22:11:43  jayhawk_hokie
*       Initial submission.
*
*
*
*
******************************************************************************/

---