//  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 [M, omega] = CL__fr_gcrs2cirs(jd, args, comega)
// GCRS to CIRS transformation matrix and angular velocity vector
//
// Calling Sequence
// [M, omega] = CL__fr_gcrs2cirs(jd, args [,comega])
//
// Description
// <itemizedlist><listitem>
// <p>Computes the frame transformation matrix <b>M</b> from GCRS to CIRS</p>
// <p>By convention, multiplying <b>M</b> by coordinates relative to GCRS yields coordinates relative to CIRS.</p> 
// <p>Optionaly computes the angular velocity vector <b>omega</b> of CIRS relative to GCRS, with coordinates relative to GCRS.
// See <link linkend="Data types">Data types</link> for more details on the definition of angular velocity vectors and frame transformation matrix.</p>
// <p></p></listitem>
// <listitem>
// <p>args is a structure containing the following fields:</p> 
// <p>- "model": (string) Available values are :
// <p>-          - "classic" for IAU2006 precession and IAU2000AR06 nutation, using classical angles.</p>
// <p>-          - "series" for IAU2006 precession and IAU2000AR06 nutation, using series.</p>
// <p>-          - "interp" for interpolation of IAU2006 precession and IAU2000AR06 nutation 
//                  for dates inside [2000-2100]. "classic" values are used for dates outside this range</p>
// <p>- "dx06" : Correction to the X coordinate of CIP. (radians)</p>
// <p>- "dy06" : Correction to the Y coordinate of CIP. (radians)</p>
// <p></p></listitem>
// <listitem>
// <p>See <link linkend="Reference frames">Reference frames</link> for more details on the definition of reference frames.</p> 
// </listitem>
// </itemizedlist>
//
// Parameters
// jd: Two-part julian day (Time scale: TT) (1xN)
// comega: (boolean, optional) Option to compute omega. If comega is %f, omega will be set to []. Default is %t. (1x1)
// args: Structure of arguments. See description for more details.
// M: Transformation matrix (3x3xN)
// omega: (optional) Angular velocity vector [rad/s] (3xN)
//
// Authors
// CNES - DCT/SB
//
// Bibliography
// 1) Technical Note 36, IERS, 2010, p48

// Declarations:

// Code:

if ~isstruct(args); CL__error("args should be a structure"); end;

if (~isfield(args,"model")); CL__error("Field(s) missing in args"); end;
if (~isfield(args,"dx06") | ~isfield(args,"dy06"))
  CL__error("Field(s) missing in args"); 
end

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

if (args.model == "interp")
  [x, y, s, xdot, ydot, sdot] = CL__iers_xys_interp(jd, comega);
elseif (args.model == "classic" | args.model == "series")
  [x, y, s, xdot, ydot, sdot] = CL__iers_xys_block(jd, args.model, comega);    
else
  CL__error("Unknown model");
end

// CIP corrections 
x = x + args.dx06;
y = y + args.dy06;

// GCRS to CIRS matrix + angular velocity
[M,omega] = CL__iers_xys2Mat(x, y, s, xdot, ydot, sdot, comega);

endfunction