/******************************************************************************
* Copyright (c) ATTX INC 2025. All Rights Reserved.
*
* This software and associated documentation (the "Software") are the 
* proprietary and confidential information of ATTX INC. The Software is 
* furnished under a license agreement between ATTX and the user organization 
* and may be used or copied only in accordance with the terms of the agreement.
* Refer to 'license/attx_license.adoc' for standard license terms.
*
* EXPORT CONTROL NOTICE: THIS SOFTWARE MAY INCLUDE CONTENT CONTROLLED UNDER THE
* INTERNATIONAL TRAFFIC IN ARMS REGULATIONS (ITAR) OR THE EXPORT ADMINISTRATION 
* REGULATIONS (EAR99). No part of the Software may be used, reproduced, or 
* transmitted in any form or by any means, for any purpose, without the express 
* written permission of ATTX INC.
******************************************************************************/
/*
WMM Magnetic field model

Author: James Tabony
*/
/*
Metadata for MS GUI:
imdata = {"exclude" : True}
*/

#ifndef MODELS_ENVIRONMENT_WORLD_MAGNETIC_FIELD_MODEL_H 
#define MODELS_ENVIRONMENT_WORLD_MAGNETIC_FIELD_MODEL_H

#include <string>

#include "simulation/Model.h"
#include "locations.h"
#include "constants/planetdefaults.h"
#include "utils/magneticFieldUtils.hpp"

namespace warptwin {

    /**
     * @brief   Spherical Harmonics Model
     * 
     * Implements the World Magnetic Model (WMM)
     * 
     * References:
     * [1] The US/UK World Magnetic Model for 2020-2025 by Arnaud Chulliat et al. (NOAA) and William Brown et al. (British Geological Survey Geomagnitsm Team)
     * [2] Pysical Geodesy by Heiskanen and Moritz 1967
     * [3] Geomagnetism and Schmidt quasi-normalization by Winch et al. 2005
     * 
     * @note This model can be expanded to a higher fidelity with the WMMHF model procuded by NOAA
     * 
     * Author: James Tabony
     * Email: james.tabony@attx.tech
     * 
    */     
    MODEL(WorldMagneticFieldModel)
    public:
        // Model params
        //         NAME                     TYPE                    DEFAULT VALUE
        START_PARAMS
            /** The filename for the file that contains the Gauss coefficients for spherical
             *  harmonic calculation. NOAA releases a new set of coefficients every 5 years on years
             *  that are divisible by 5 (e.g. 2015, 2020, 2025, ...) that are most valid from the
             *  year they are realseed until the next coefficient set release. */
            SIGNAL(filename,                std::string,            warptwinDir() + "data/WMM/WMM_2025.COF")
            /** The coefficient n which represents the variations in gravitational field in the latitudinal direction.
             *  Also the harmonic function order. This value must be chosen before executive startup. (1 <= N <= 12) */
            SIGNAL(N,                       int,                    12)
        END_PARAMS

        // Model inputs
        //         NAME                     TYPE                    DEFAULT VALUE
        START_INPUTS
            /** The body position wrt the planet in a planet-centered rotating frame (meters) */
            SIGNAL(pos_body_planet__pcr,    CartesianVector3,       CartesianVector3({0.0, 0.0, 0.0}))
        END_INPUTS

        // Model Outputs
        //         NAME                     TYPE                    DEFAULT VALUE
        START_OUTPUTS
            /** The local magnetic field vector expressed in body's local North East Down frame. (nanoTesla) */
            SIGNAL(mag_field__NED,          CartesianVector3,       CartesianVector3({0.0, 0.0, 0.0}))
            /** The local magnetic field vector expressed in ECEF frame. (nanoTesla) */
            SIGNAL(mag_field__ECEF,         CartesianVector3,       CartesianVector3({0.0, 0.0, 0.0}))
        END_OUTPUTS

        int16 activate() override;
        int16 deactivate() override;
    protected:
        int16 start() override;
        int16 execute() override;

        /// @brief Implicit returns form the read coefficients function in magneticFieldUtils.h
        /// @param _g Gauss coefficients g (nT)
        /// @param _h Gauss coefficients h (nT)
        /// @param _gdot Time derivative of Gauss coefficients g (nT/year)
        /// @param _hdot Time derivative of Gauss coefficients h (nT/year)
        /// @param _epoch Decimal year of the associated data
        double _g[NUMBER_OF_READ_COEFFICIENTS];
        double _h[NUMBER_OF_READ_COEFFICIENTS];
        double _gdot[NUMBER_OF_READ_COEFFICIENTS];
        double _hdot[NUMBER_OF_READ_COEFFICIENTS];
        double _epoch;

        // Local variables for lla coordinates
        double _lat_centric, _lat_detic;
        double _lon;

        // Dummy variable for implicit returns that can be ignored
        double _dummy;

        // Local variable for magnitude of position vector
        double _rad;

        // Magnetic field vector components in geocentric coordinates [reference 1]
        double _Xprime, _Yprime, _Zprime;

        // Temporary variable for the time difference between .COF epoch and current decimal year
        double dt;
        // Variable for warning the user if the coefficient file used is not optimal
        bool _WARN_BAD_FILE_EPOCH;

        // Temporary value for Gauss coefficients after time interpolation at a given n and m
        double _g_n_m, _h_n_m;
        // Tempory value for the index in the coefficient array that corresponds to a given n and m
        int _index;

        // Temporary variable for folding the evaluation of the Schmidt semi-normalized associated Legendre function
        double _P_n_m;

        // Variable for holding the model order on startup
        int _N;
        // Variable for warning the user if they change the model order after startup
        bool _WARN_CANT_CHANGE_N;

        // Arrays for holding the Schmidt normalized associated Legendre function outputs
        double *_P, *_dP, *_PcupS;

        // Boolean variable for tracking if at a singularity
        bool _at_singularity;

        // Temporary variables for faster computation
        double _pow_a_by_r, _cos_m_lon, _sin_m_lon;

        // Temporary DCM for coordinate transformations
        clockwerk::DCM _temp_NED_ECEF;

    };

}

#endif