; docformat = 'rst'
;
; NAME: 
;    shapeCreateEllipsoid
; PURPOSE:
;    Build a 3-D shape model from tri-axial diameters
;+
; :Description:
;    Build a 3-D shape model from tri-axial diameters
;
; :Categories:
;    Ephemerides, Shape 3-D
;
; :Params:
;    a: in, optional, type=float, default=1
;      Equatorial diameter
;    b: in, optional, type=float, default=a
;      Second equatorial diameter
;    c: in, optional, type=float, default=b
;      Polar diameter
;
; :Returns: A structure containing the shape model, fields are similar
;    to readVer function::
;       .VERTEX:  A structure with Information on Vertex::
;          .N   = Number of Vertex
;          .XYZ = Cartesian coordinates of the vertex
;       .FACET:  A structure with Information on Facet::
;          .N   = Number of Facets
;          .CON = Connectivity array ([N_vertex, ID_vertex....])
;
; :Keywords:
;    nbLong: in, optional, type=int, default=36
;      Number of evaluation points along longitudes
;    nbLat: in, optional, type=int, default=9
;      Number of evaluation points along latitudes
;    dump: in, optional, type=string
;      If set, export the shape model into a file at dump path.
;
; :Uses:
;    writeVer
;
; :Examples:
;    Create a 3x2x1 ellipsoid and plot it from the top::
;      IDL>  model = shapeCreateEllipsoid(3.,2.,1.)
;      IDL>  cgPlot, model.vertex.xyz[0,*], model.vertex.xyz[1,*],psym=4,/iso
;
;    Create an oblate spheroid 1.1x1.1x1 and write it on disk::
;      IDL>  model = shapeCreateEllipsoid(1.1,1.1,1, dump='oblate.ver')
;
; :Author:
;    B.Carry (OCA)
;
; :History:
;    Change History::
;      Original Version written in June 2017, B. Carry (OCA), 
;-
function shapeCreateEllipsoid, a, b, c, nbLong=nbLong, nbLat=nbLat, dump=dump
  COMPILE_OPT hidden, idl2
;-----------------------------------------------------------------------------------------------;
;-----------------------------------------------------------------------------------------------;
;---  TAG ---    I -- Initialization And Input Verification              -----------------------;
;-----------------------------------------------------------------------------------------------;
;-----------------------------------------------------------------------------------------------;
 ;--I.1-- Ellipsoid Default Diameters  ---------------------------------------------------------
 ;--I.1.1-- Set Values
  if not keyword_set(a) then a=1.
  if not keyword_set(b) then b=a
  if not keyword_set(c) then c=b

 ;--I.1.2-- Compute Square Values
  a2=float(a*a)
  b2=float(b*b)
  c2=float(c*c)
  

 ;--I.2-- Number of Evaluation Points  ---------------------------------------------------------
  if not keyword_set(nbLong) then nbLong=72
  if not keyword_set(nbLat)  then nbLat =18


 ;--I.3-- Spherical Coordinates  ---------------------------------------------------------------
 ;--I.3.1-- Define Angles
  theta=2*!PI*findgen(nbLong)/nbLong
  phi  =!PI*findgen(nbLat)/nbLat

 ;--I.3.2-- Compute Trogonometric Values
  cosTheta = cos(theta)
  sinTheta = sin(theta)
  cosPhi   = cos( phi )
  sinPhi   = sin( phi )


;-----------------------------------------------------------------------------------------------;
;-----------------------------------------------------------------------------------------------;
;---  TAG ---   II -- Create an Ellispoid Shape Model                    -----------------------;
;-----------------------------------------------------------------------------------------------;
;-----------------------------------------------------------------------------------------------;
 ;--II.1-- Create Triangular Facets  -----------------------------------------------------------
 ;--II.1.1-- Northern Pole
  xyz = [0, 0, c]

 ;--II.1.2-- Mid-Latitude Points
  for kLat=1, nbLat-2 do begin
    for kLon=0, nbLong-1 do begin
      x = a * cosTheta[kLon] * sinPhi[kLat]
      y = b * sinTheta[kLon] * sinPhi[kLat]
      z = c *                  cosPhi[kLat]
      xyz = [ [xyz], [x,y,z] ]
    endfor
  endfor

 ;--II.1.3-- Southern Pole
  xyz = [ [xyz], [0,0,-c] ]


 ;--II.2-- Create Radius Vector Array  ---------------------------------------------------------
  mesh_obj, 0, oVert, oConn, xyz
    nbFacet = n_elements(oConn)/4L
    nbVertex= n_elements(oVert)/3L


 ;--II.3-- Store 3-D Model in IDL Structure  ---------------------------------------------------
  vertexStr = { n:nbVertex, xyz:oVert }
  facetStr  = { n:nbFacet, con:oConn}
  model = {vertex:vertexStr, facet:facetStr }


 ;--II.6-- Write Model on Disk  [Optional]  ----------------------------------------------------
  if keyword_set(dump) then writeVer, dump, model


 ;--II.7-- Return Model Structure  -------------------------------------------------------------
  return, model

end