% -------------------------------------------------------------------------
% Least Squares_Method
% МЕТОД НАИМЕНЬШИХ КВАДРАТОВ ДЛЯ СВОБОДНО ОПРЕДЕЛЯЕМЫХ БАЗИСНЫХ ФУНКЦИЙ
% -------------------------------------------------------------------------
% Синтаксис вызова:
%                       [A]=LSM01Free(F,V)
% 
%  F - Исходная матрица табличная функция вида:
%
%    x1 x2 x3 x4 .... xM
%    y1 y2 y3 y4 .... yM
%
%  V - Литеральный вектор столбец описывающий базисные функции :
%
%    v1
%    v2
%   ....................
%    vN
% 
%  Например (подготовка тестовых данных): 
% 
%      x=[0:pi/100:2*pi];
%      y=eval('x+sin(x)-log(x+1)');
%      F=[x;y];
%
%      V = ['x^0     '; 'x       '; 'sin(x)  '; 'log(x+1)']; 
%      ВНИМАНИЕ! Длина строк всех литералов должна быть строго одинаковой и
%      при необходимости расширяться пробелами справа.
%
%      Другой вариант построения вектора V:
%      V = char('x^0', 'x', 'sin(x)', 'log(x+1)');
%
%  Вызов метода:
%  [A]=lsm01Free(F,V) 
%
%  Результат:
%  Максимальная абсолютная ошибка : 3.5527e-014
%  A =
%    0.0000
%    1.0000
%    1.0000
%   -1.0000
% -------------------------------------------------------------------------
%
function [A]=lsm01Free(F,V) 
if isempty(F) || isempty(V)
   disp('Одна из входных матриц не задана'); 
   A=NaN;
   return
end;
% Минимальный входной контроль
if (size(F,1) ~= 2)
   disp('Размеры матрицы [F] заданы неправильно'); 
   A=NaN;
   return
end;
nV=size(V,1);
if (size(F,2) < nV)
   disp('Недостаточно точек в [F] для определения коэффициентов для [V]'); 
   A=NaN;
   return
end;
% Подготовка рабочих матриц
C=zeros(nV,nV);
B=zeros(nV,1);
% Построим систему линейных уравнений для аппроксимации
% Для всех строк матрицы 'C'
for row=1:nV
    % Для всех 'x' матрицы 'F'
    B(row,1)=0;
    for k=1:size(F,2)
        x=F(1,k);
        B(row,1)=B(row,1)+eval(V(row,:))*F(2,k);  
    end;
    % Для всех колонок матрицы 'C'
    for col=1:nV
        % Для всех 'x' матрицы 'F'
        C(row,col)=0;
        for k=1:size(F,2)
            x=F(1,k);
            C(row,col)=C(row,col)+eval(V(row,:))*eval(V(col,:));  
        end;
    end;
end;
% Вычислим коэффициенты для аппроксимирующего выражения
ab=cat(2,C,B);
A = fgauss01(ab);

% Вычислим значения функции по аппроксимирующему выражению
% Для всех 'x' матрицы 'F'
for k=1:size(F,2)
    x=F(1,k);
    y(1,k)=0;
    % Для всех строк матрицы 'V'
    for row=1:nV
        y(1,k)=y(1,k)+A(row,1)*eval(V(row,:)); 
    end;   
end;

clf;
hold on;

disp('График исходной табличной функции'); 
% Построим график исходной табличной функции
plot(F(1,:),F(2,:),':b');
grid on;

disp('Далее...'); 
reply=input('Построить график по результатам аппроксимации?  y/n [y]: ','s');
if isempty(reply)
    reply = 'y';
end;
if reply == 'y'
   % Построим график по результатам аппроксимации
   plot(F(1,:),y,'-r');
end;   
grid on;
hold off;

disp('Далее...');
delta=F(2,:)-y;
maxdelta=max(abs(delta));
disp(strcat('Максимальная абсолютная ошибка : ', num2str(maxdelta)));

disp('Далее...'); 
reply=input('Построим график абсолютной ошибки аппроксимации?  y/n [y]: ','s');
if isempty(reply)
    reply = 'y';
end;
if reply == 'y'
   % Построим график абсолютной ошибки
   if max(abs(delta)) > 1e-16
      plot(F(1,:),delta);
      grid on;
   end;
end;

disp('Далее...');
reply=input('Для завершения нажмите Enter...');
close all hidden;