//  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 [eps0,eps0p,eps0pp] = CL_mod_meanObliquity(jj_tai,sec_tai, flag)
// Mean obliquity (IERS 1996) - DEPRECATED
//
// Calling Sequence
// [eps0,eps0p,eps0pp] = CL_mod_meanObliquity(jj_tai,sec_tai [,flag])
//
// Description
// <itemizedlist><listitem>
// <p>This function is deprecated. </p>
// <p></p></listitem>
// <listitem>
// <p>Computes the mean obliquity and its time derivatives (see IERS standards 1996).</p> 
// </listitem>
// </itemizedlist>
//
// Parameters
// jj_tai : (integer) Modified (1950.0) Julian day (1xN)
// sec_tai : Seconds  (TAI time scale) (1xN)
// flag : (string, optional) "p" -> eps0 and eps0p; "s" -> eps0, eps0p and eps0pp (default depends on number of outputs)
// eps0 : Mean obliquity [rad] (1xN)
// eps0p : (optional) First time derivative of [rad/s] (1xN)
// eps0pp : (optional) Second time derivative of eps0 [rad/s^2] (1xN)
//
// Authors
// CNES - DCT/SB
//
// Bibliography
// 1) IERS Conventions (1996), Dennis D. McCarthy
// 2) Explanatory Supplement to the Astronomical Almanac, Seidelman (1992)
//
// See also
// CL_mod_nutationAngles
// CL_mod_nutationArg
// CL_mod_precessionAngles
//
// Examples
// jj_tai = [19500:1:20500];
// ss_tai = 43200*ones(jj_tai);
// [eps0,eps0p,eps0pp]=CL_mod_meanObliquity(jj_tai,ss_tai,"s");
// plot(jj_tai'+ss_tai'/86400,eps0'*180/%pi);
//

// Declarations:


// Code:

CL__warnDeprecated(); // deprecated function

[lhd,rhs]=argn();

if (~exists('flag','local'))
  select lhs
    case 1
      flag='n';
    case 2
      flag='p';
    case 3
      flag='s';
    else
      CL__error("check number of output arguments");
  end
end

if (rhs>3)|(rhs==1) CL__error("check number of input arguments"); end,

if (flag~="s" & flag~="p" & flag~="n") then
  CL__error("value of parameter ''flag'' = "+flag+" unknown");
end

J2000 = 36525. / 2.;
TT_MOINS_TAI = 32.184;
t = (jj_tai - J2000 + (TT_MOINS_TAI+sec_tai)/86400) / 36525;
t = t(:)';

eps0 = 84381.448d0 - 46.8150d0*t - 0.00059d0*t.^2 + 0.001813d0*t.^3;
SEC2RAD = (%pi / 180) / 3600;
eps0 = eps0*SEC2RAD;
eps0p=[]; eps0pp=[];
if (flag=="p" | flag=="s") then
  eps0p = -46.8150d0 - 2*0.00059d0*t + 3*0.001813d0*t.^2;
  eps0p = eps0p*SEC2RAD/(86400*36525);
end
if flag=="s" then
  eps0pp = -2*0.00059d0 + 6*0.001813d0*t;
  eps0pp = eps0pp*SEC2RAD/((86400*36525)^2);
end

endfunction