function [y0,w0,T,fi0,delta,S,TT,YY,YY2] = w3(y,t)
% YY  - вихідний ряд затухаючих коливань зі збільшеною кількістю проміжних точок
% YY2 - модельний ряд затухаючих коливань зі збільшеною кількістю проміжних точок

[m,n] = size(y);
[c,d] = size(t);

if n~=d disp('Matrix dimensions must agree.'), end;

% знаходимо експоненційну огигнаючу та delta
ln_y = log(y);
p = polyfit(t,ln_y,1);
pv =polyval(p,t);
delta = abs(p(1));

% збільшуємо кількість точок за допомогою сплайнів
TT=[min(t):0.001:max(t)];
YY=spline(t,y,TT);
[mm,nn]=size(YY);

% знаходимо точки перетину графіка з осями
A=sign(YY); 
B=A(1:nn-1).*A(2:nn); 
I=find(B<0);
% вибираємо перші 5 точок
I=I(1:5);
TTT=TT(I);
% визначаємо період T
T=(sum(TTT(2:5)-TTT(1:4)))/2;

% знаходимо частоту
w0 = 2*pi/T;

% визначаємо початкову швидкість
v = polyfit(TT(1:3),YY(1:3),1);
v0 = v(1);
%find y0, fi0
%w2 = sqrt(w0^2 - delta^2)

% знаходимо початкову амплітуду
first = v0 + delta*y(1);
y0 = sqrt(y(1)^2 + (first)^2)/w0;
second = y(1)*sqrt(w0);

% знаходимо початкову фазу
if ((first >= 0)&(second >= 0))|((first >= 0)&(second <= 0)), fi0 = atan(second/first);
elseif (first <= 0)&(second >= 0) fi0 = atan(second/first) + pi;
elseif (first <= 0)&(second <= 0) fi0 = atan(second/first) - pi;
end;

% розраховуємо точки модельного графіка
YY2=y0*exp(-delta*TT).*sin(w0*TT+fi0);
L=max(YY2)/max(YY);
YY2=YY2/L;

% шукаємо суму квадратів невязок S
S = (YY2-YY)*(YY2-YY)';

% save 'YY2.txt' YY2 -ASCII
% save 'YY.txt' YY -ASCII
