//  Copyright (c) CNES  2008
//
//  This software is part of CelestLab, a CNES toolbox for Scilab
//
//  This software is governed by the CeCILL  license under French law and
//  abiding by the rules of distribution of free software.  You can  use,
//  modify and/ or redistribute the software under the terms of the CeCILL
//  license as circulated by CEA, CNRS and INRIA at the following URL
//  'http://www.cecill.info'.

function [dcir] = CL_op_shortPeriodsJ2Circ(sma,inc,pso, er,mu,j2)
// Short period terms due to J2 + central force
//
// Calling Sequence
// dcir = CL_op_shortPeriodsJ2Circ(sma,inc,pso [,er,mu,j2])
//
// Description
// <itemizedlist>
// <listitem>
// <p>Computes the short term effects (periods less than one orbit) due to J2 and the central force 
// on orbital elements for nearly circular orbits. </p>
// <p></p></listitem>
// <listitem> 
// <p>The formulas used are the following:</p>  
// <p><inlinemediaobject><imageobject><imagedata fileref="short_periods_j2.gif"/></imageobject></inlinemediaobject></p>
// <p></p></listitem>
// <listitem>
// <p>The orbital elements used are the following:</p>  
// <p><inlinemediaobject><imageobject><imagedata fileref="cir_par.gif"/></imageobject></inlinemediaobject></p>
// <p></p></listitem>
// <listitem> 
// <p>Notes:</p> 
// <p>- The osculating orbital elements are obtained by adding the mean orbital elements and <b>dcir</b>.</p>
// <p>- The short term effects due to J2 squared are not included in the formulas.</p> 
// <p></p></listitem>
// </itemizedlist>
//
// Parameters
// sma: Mean semi-major axis [m] (1xN or 1x1) 
// inc: Mean inclination [rad] (1xN or 1x1) 
// pso: Mean argument of latitude (1xN or 1x1) 
// er: (optional) Equatorial radius [m]. Default is %CL_eqRad
// mu: (optional) Gravitational constant [m^3/s^2]. NOT USED.
// j2: (optional) Second zonal harmonic term). Default is %CL_j1jn(2)
// dcir: J2 periodic (period less than 1 orbit) effects on orbital elements (a, ex, ey, inc, raan, arg of latitude) (6xN)
//
// Authors
// CNES - DCT/SB 
//
// Bibliography
// 1) Mecanique Spatiale, Cnes - Cepadues Editions, Tome I, section 13.4.1.1.1 (p 804) - with typo corrected. 
//
// See also
// CL_ex_eckHech
//
// Examples
// sma = 7.e6;
// inc = CL_deg2rad(98);
// pso = linspace(0,2*%pi,13); 
// CL_op_shortPeriodsJ2Circ(sma, inc, pso)


// Declarations:

// Code:
if (~exists("er", "local")); er = CL__dataGetEnv("eqRad"); end
if (~exists("j2", "local")); j2 = CL__dataGetEnv("j1jn", 2); end

[sma,inc,pso] = CL__checkInputs(sma,1,inc,1,pso,1); 
 
cosi  = cos(inc);
sini  = sin(inc);
sini2 = sini.^2;
q     = 1.5 * j2 * (er ./ sma).^2; 


// K : only term in J2 (not J2^2)
// initial formula : K     = q ./ (1 + q .* (3 - 4 * sini2));
K = q; 

C     = cos([pso; 2*pso; 3*pso]); 
S     = sin([pso; 2*pso; 3*pso]); 

// Effect on orbital elements : 
// sma, ex, ey, inc, raan, arg of latitude

dcir = [ ...
    K .* ( sma .* sini2 .* C(2,:) ); ...                                 
    K .* ( (1 - 1.25 * sini2) .* C(1,:) + (7/12) * sini2 .* C(3,:) ); ...
    K .* ( (1 - 1.75 * sini2) .* S(1,:) + (7/12) * sini2 .* S(3,:) ); ...
    K .* ( 0.5 * sini .* cosi .* C(2,:) ); ...                       
    K .* ( 0.5 * cosi .* S(2,:) ); ...                                 
    K .* ( (1.25 * sini2 - 0.5) .* S(2,:) ) ...
]; 

endfunction