//  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_convertMat(frame1, frame2, cjd, ut1_tref, tt_tref, xp, yp, dX, dY, use_interp)
// Frame transformation matrix and angular velocity vector 
//
// Calling Sequence
// [M, omega] = CL_fr_convertMat(frame1, frame2, cjd, [ut1_tref, tt_tref, xp, yp, dX, dY, use_interp])
//
// Description
// <itemizedlist><listitem>
// <p>Computes the frame transformation from <b>frame1</b> to <b>frame2</b>, which is characterized by 
// a frame transformation matrix <b>M</b> and an angular velocity vector <b>omega</b>. </p>
// <p>By convention: </p>
// <p>- multiplying <b>M</b> by coordinates relative to <b>frame1</b> yields coordinates relative to <b>frame2</b>. </p>
// <p>- <b>omega</b> is the angular velocity vector of <b>frame2</b> wrt <b>frame1</b>, with coordinates in <b>frame1</b>. </p>
// <p></p></listitem>
// <listitem>
// <p>The available frames are: </p>
// <p><inlinemediaobject><imageobject><imagedata fileref="frame_convert.gif"/></imageobject></inlinemediaobject></p>
// <p></p>
// <p>In addition, 2 other frames (or names) are defined: ECI and ECF. </p>
// <p>See <link linkend="Reference frames">Reference frames</link> for more details on the definition of reference frames.</p> 
// <p></p></listitem>
// <listitem>
// <p><b>Notes:</b></p>
// <p>- The date 'cjd' is relative to the TREF time scale. </p>
// <p>- Which  optional arguments need to be supplied depends on the frame transformation. See 
// the frame transformation diagram for more details. </p>
// <p>- The frame names are case sensitive.</p>
// </listitem>
// </itemizedlist>
//
// Parameters
// cjd: Modified Julian date from 1950.0 (TREF time scale) (1xN or 1x1)
// frame1: (string) Initial frame (1x1)
// frame2: (string) Final frame (1x1)
// ut1_tref: (optional) UT1-TREF [seconds]. Default is %CL_UT1_TREF (1xN or 1x1)
// tt_tref: (optional) TT-TREF [seconds]. Default is %CL_TT_TREF. (1xN or 1x1)
// xp, yp: (optional) Position of the pole in ITRS [rad](1xN or 1x1)
// dX, dY: (optional) Corrections to CIP coordinates [rad] (1xN or 1x1)
// use_interp: (optional, boolean) %t to use interpolation when computing CIP coordinates. Default is %t. (1x1)
// M: Frame transformation matrix from frame1 to frame2 (3x3xN)
// omega: Angular velocity vector of frame2 wrt frame1 with coordinates in frame1 [rad/s] (3xN)
//
// Authors
// CNES - DCT/SB
//
// Examples
// // Frame transformation from Veis to EME2000 (All optional arguments to default)
// cjd_tref = 21010;
// [M,omega] = CL_fr_convertMat("Veis", "EME2000", cjd_tref);
//
// // Frame transformation from GCRS to TIRS
// cjd_tref = [21010 , 21011];
// [M,omega] = CL_fr_convertMat("GCRS", "TIRS", cjd_tref, ut1_tref=-0.2);

// Declarations:


// Code:
if (argn(2) < 3)
  CL__error('Not enough input arguments');
end

if (~exists("ut1_tref", "local")); ut1_tref = CL__dataGetEnv("UT1_TREF"); end
if (~exists("tt_tref", "local")); tt_tref = CL__dataGetEnv("TT_TREF"); end
if (~exists("xp", "local")); xp = 0; end
if (~exists("yp", "local")); yp = 0; end
if (~exists("dX", "local")); dX = 0; end
if (~exists("dY", "local")); dY = 0; end
if (~exists("use_interp", "local")); use_interp = %t; end

// Check argument sizes, and resize if necessary:
[cjd,ut1_tref,tt_tref,xp,yp,dX,dY] = CL__checkInputs(cjd,1, ut1_tref,1, tt_tref,1, xp,1, yp,1, dX,1, dY,1);

// opt = flag specifying to the internal function CL__fr_convert what needs to be computed
if (argn(1) == 1)
  opt = "mat";
else
  opt = "mat_omega";
end

// Arguments structure to be used in all sub-functions :
if (use_interp)
  model = "interp"; 
else
  model = "classic";
end

args = struct( ...
  "model", model, ...
  "precession_model", "2006", ...
  "nutation_model", "2000AR06", ...
  "ut1_tt", ut1_tref-tt_tref, ...
  "xp", xp, ...
  "yp", yp, ...
  "dx06", dX, ...
  "dy06", dY);

// Convert input date from MJD1950.0 (TREF time scale) 
// to a two part JD (TT time scale)
// ttjda + ttjdb = cjd + 2433282.5;
ttjda = floor(cjd) + 2433282;
ttjdb = cjd - floor(cjd) + 0.5 + tt_tref / 86400;
ttjd = [ttjda;ttjdb]; 
 
[M,omega] = CL__fr_convert(frame1, frame2, ttjd, [], [], opt, args);

endfunction