//  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 [s, sdot] = CL__iers_s2006(jd,x,y,xdot,ydot, cder)
// CIO locator (s) as function of CIP's X,Y coordinates (IAU 2006)
//
// Calling Sequence
// s = CL__iers_s2006(jd, x, y)
// [s, sdot] = CL__iers_s2006(jd, x , y, xdot, ydot [,cder])
//
// Description
// <itemizedlist><listitem>
// <p>The CIO locator s, positioning the Celestial Intermediate Origin on
// the equator of the Celestial Intermediate Pole, given the CIP's X,Y
// coordinates. Compatible with IAU2006 precession and IAU2000AR06 nutation models.</p>
// <p></p></listitem>
// </itemizedlist>
//
// Parameters
// jd: Two-part julian day (Time scale: TT) (2xN)
// x: X coordinate of CIP (1xN)
// y: Y coordinate of CIP (1xN)
// xdot: Time derivative of X coordinate of CIP (1xN)
// ydot: Time derivative of Y coordinate of CIP (1xN)
// cder: (boolean, optional) Option to compute derivatives. If cder is %f, the derivatives will be set to []. Default is %t. (1x1)
// s: CIO locator (1xN)
// sdot: Time derivative of CIO locator (1xN)
//
// Authors
// CNES - DCT/SB
//
// Bibliography
// 1) Technical Note 36, IERS, 2010


// Declarations:

// Code:

if (~exists("cder","local")); cder = %t; end;
if (argn(1) <= 1); cder = %f; end;

// xdot and ydot need to be provided in order to compute sdot
if (cder & argn(2) <= 3)
  CL__error("Invalid number of input arguments");
end

// Conversion factors
ARCSEC_TO_RAD = %pi / (180*3600);
CENT_TO_SEC = 86400 * 36525; 

// Julian centuries since J2000 (12h00)
t = ((jd(1,:) - 2451545) + jd(2,:)) / 36525;

// Fundamental arguments (units: rad and rad/s)
[Fls, Flsdot] = CL__iers_fundArgsLS(jd, cder);
[Fpl, Fpldot] = CL__iers_fundArgsPL(jd, cder);


// Loads series
fname = fullfile(CL_home(), "data", "frame", "iers_series_s_2006_2000AR06.dat"); 
if (~isfile(fname))
  CL__error("File: " + fname + ": not found"); 
end
load(fname, "data"); 

// s and sdot: initially arcsec and arcsec / century, 
// converted to rad and rad/s
// Note: fund args derivatives converted to rad/century
// because unit of t = century

[s, sdot] = CL__iers_seriesEval(data, t, [Fls; Fpl], [Flsdot; Fpldot]*CENT_TO_SEC, cder); 

s = s * ARCSEC_TO_RAD - x .* y / 2; 

if (cder)
  sdot = sdot * ARCSEC_TO_RAD / CENT_TO_SEC -(x.*ydot + xdot.*y)/2;; 
end

endfunction