//  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 [vmin,vmax,nb1,nb2] = CL_graduate(xmin,xmax,ns1,ns2)
// Simple graduation calculation
//
// Calling Sequence
// [vmin,vmax,nb1,nb2] = CL_graduate(xmin,xmax,ns1,ns2)
//
// Description
// <itemizedlist>
// <listitem>
// <p>Computes uniformly distributed graduation values. </p>
// <p>The steps are computed such that the corresponding values (vmin+k*(vmax-vmin)/nb1, k=0:nb1) 
// are of the form: k*10^n (i: integer), with k in [1,2,3,5,6]. The substeps are then chosen 
// accordingly. For instance, if k=1, possible values for nb2 are 1, 2, 5, 10, and so on. The number of 
// substeps is limited to 12. </p>
// </listitem>
// </itemizedlist>
//
// Parameters
// xmin: User data minimum value (1x1).
// xmax: User data maximum value (1x1).
// ns1: (integer, optional) Approximate number of steps wanted. Default is 7. 
// ns2: (integer, optional) Approximate number of sub-steps wanted. Default is 5. 
// vmin: Interval minimum value (vmin <= xmin).
// vmax: Interval maximum value (vmax >= xmax).
// nb1: (integer) Number of steps.
// nb2: (integer) Number of sub-steps.
//
// Authors
// CNES - DCT/SB
//
// See also
// CL_autoLevels
//
// Examples
// [vmin, vmax, nb1, nb2] = CL_graduate(0, 6)
// [vmin, vmax, nb1, nb2] = CL_graduate(1.14, 2.27, 5, 5)

// Declarations:


// Code:

if (~exists('ns1', 'local')); ns1=7; end
if (~exists('ns2', 'local')); ns2=5; end

// allows for values around ns1
ns1min = max(ns1-2, 1);
ns1max = max(ns1+2, 1);

if (xmax < xmin)
  x = xmin; 
  xmin = xmax; 
  xmax = x; 
end

ampl = xmax - xmin;

// If amplitude == 0 (xmin == xmax) 
// => find appropriate solution (security)
if (ampl <= 0)
  xmin = xmin - 0.1*xmin;
  xmax = xmax + 0.1*xmax;
  ampl = xmax - xmin;
  if (ampl <= 0) // if xmax == xmin == 0
    xmin = -1;
    xmax = 1;
    ampl = xmax - xmin;
  end
end

// sort the steps according to the difference with ns1 
ns1_all = (ns1min:ns1max);  
[X, I] = gsort(abs(ns1_all-ns1), "c", "i"); 
pas = ampl ./ ns1_all(I); 

fact = (10 .^ floor(log10(pas)));
pas = pas ./ fact;

// NB: value 10 for security - should not be necessary
valref = [1;2;3;5;6;10]; // step = valref(k) * 10^n. 

// A : contains all combinations
A = abs( ones(valref) * pas - valref * ones(pas) ); 

// NB: minimum on a matrix:
// - 2nd output argument returned by min is such that ij(1) = row index and ij(2) = column index
// - min returns only one minimum (the first encountered)
[delta, ij] = min(A);
pas = valref(ij(1)) * fact(ij(2));

vmin = floor(xmin/pas) * pas;
vmax = ceil(xmax/pas) * pas;
nb1 = round((vmax-vmin)/pas);

// substeps
Asub = [ 1,2,5,10;  // values for steps multiple of 1
         1,2,4,10;  // values for steps multiple of 2
         1,3,6,12;  // values for steps multiple of 3
         1,1,5,10;  // values for steps multiple of 5
         1,3,6,12;  // values for steps multiple of 6
         1,2,5,10 ] // values for steps multiple of 10

i = ij(1);
// NB: minimum on a vector ==> k is (1x1)
[v, k] = min(abs(Asub(i,:) - ns2));
nb2 = Asub(i,k);

endfunction