//  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 [oe2] = CL_stela_convert(cmd, type_oe, cjd, oe, params, frame, ut1_tref, tt_tref)
// Mean/osculating conversion using STELA
//
// Calling Sequence
// oe2 = CL_stela_convert(cmd, type_oe, cjd, oe, params [, frame, ut1_tref, tt_tref])
//
// Description
// <itemizedlist>
// <listitem>
// <p>Convert orbital elements from mean to osculating or vice versa. </p>
// <p></p>
// <p>The input orbital elements are internally converted to the adequate type and frame before 
// the call to STELA. The results are then converted back to 
// the initial type and frame. </p>
// <p></p>
// </listitem>
// <listitem>
// <p><b>Notes</b>:</p>
// <p>- The integration time step parameter (in the "params" structure) is not used but must exist. </p>
// <p></p>
// </listitem>
// <listitem>
// <p>See the <link linkend="STELA">STELA</link> page for more details.</p>
// <p></p>
// </listitem>
// </itemizedlist>
// 
// Parameters
// cmd: (string) Wanted operation: "mean2osc" or "osc2mean" 
// type_oe: (string) Type of orbital elements: "kep", "cir", "cireq", "equin", "pv". 
// cjd: epochs (CJD, time scale: TREF) (1xN)
// oe: Mean or osculating orbital elements (6xN) 
// params: (structure) Propagation model parameters. 
// frame: (string, optional) Input/output frame. Default is "ECI" 
// ut1_tref: (optional) UT1-TREF [seconds]. Default is %CL_UT1_TREF (1x1)
// tt_tref: (optional) TT-TREF [seconds]. Default is %CL_TT_TREF (1x1)
// oe2: Osculating or Mean orbital elements (6xN)
//
// Authors
// CNES - DCT/SB/MS
//
// Examples
//
// // Date (cjd, TREF) and Keplerian mean orbital elements
// cjd = 20000; 
// mean_kep = [7.e6; 1.e-3; 1.7; 1; 2; 3]; 
//
// // STELA model parameters (default values)
// params = CL_stela_params(); 
//
// // Convert to osculating elements
// osc_kep = CL_stela_convert("mean2osc", "kep", cjd, mean_kep, params);
//
// // Check: convert back to mean elements
// mean_kep2 = CL_stela_convert("osc2mean", "kep", cjd, osc_kep, params); 
// mean_kep2 - mean_kep   // => 0


// Declarations:

// Code
if (~exists("frame", "local")); frame="ECI"; 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 

// check extension module
CL__checkExtensionModule(); 

if (typeof(frame) <> "string" | size(frame, "*") <> 1)
  CL__error("Invalid argument type or size for frame"); 
end

if (size(cjd,1) <> 1 | size(oe,6) <> 1 | size(cjd,2) <> size(oe,2))
  CL__error("Invalid arguments sizes"); 
end

// converts orbital elements to "cireq" type
// NB: uses "mu" from STELA
// converts frame to MOD if necessary
// osculating are supposed valid in same frame as assumed by stela_extrap (for consistency)
//

// STELA constants (mu)
CL__stela_manageConst(); 
cst = CLx_stela_getConst(1); 


// if not "natural" frames => compute frame conversion for all dates
convfr = %f; 
if (~(frame == "CIRS" | (frame == "ECI" & CL_configGet("ECI_FRAME") == "CIRS")))
  [M, omega] = CL_fr_convertMat(frame, "CIRS", cjd, cjd(1), ut1_tref=ut1_tref,  tt_tref=tt_tref); 
  convfr = %t; 
end

if (convfr)
  pv = CL_oe_convert(type_oe, "pv", oe, cst.mu); 
  [p, v] = CL_rot_pvConvert(pv(1:3,:), pv(4:6,:), M, omega); 
  x = CL_oe_convert("pv", "cireq", [p; v], cst.mu); 

else
  // transformation through "pv" since direct conversion may not exist
  pv = CL_oe_convert(type_oe, "pv", oe, cst.mu); 
  x = CL_oe_convert("pv", "cireq", pv, cst.mu); 
end

[x2, info] = CLx_stela_convert(cmd, cjd, x, params, ut1_tref, tt_tref); 

// check status
if (info.status < 0)
  CL__error(info.err_msg); 
end

if (convfr)
  // change direction
  [M, omega] = CL_rot_compose(M, omega, -1); 

  pv = CL_oe_convert("cireq", "pv", x2, cst.mu); 
  [p, v] = CL_rot_pvConvert(pv(1:3,:), pv(4:6,:), M, omega); 
  oe2 = CL_oe_convert("pv", type_oe, [p;v], cst.mu); 

else
   // cireq type -> initial type
  oe2 = CL_oe_convert("cireq", type_oe, x2, cst.mu); 
  
end

      
endfunction