clear;

N=1024;
% making two-sin signal
for (i=1:1:N)
    y(i)=sin(2*pi/10*(i-1))+sin(2*pi/100*(i-1));
end;

plot(y);
% ARMA process generating
ar(1)=0.5*randn(1);

for (i=2:1:N)
    ar(i)=-0.9*ar(i-1)+0.5*randn(1);
end;

plot(ar);

%making a sum
s=y+ar;

plot(s);

%output to file
%cd d:\arts\filtr\4;
fout=fopen('signal.dat','wt');
for(i=1:1:N)
 fprintf(fout,'%12.6f%12.6f\n',i,s(i)); 
end;
fclose(fout);

savefile = 'signal.mat';
save(savefile,'s');

% fast Fourier transformation
N_ft=N;

fts=fft(s,N_ft);

plot(abs(fts)/N_ft);

% frequency calculation
freq=(1:N_ft/2-1)/N_ft;

plot(freq, abs(fts(2:N_ft/2))/N_ft);
title('modul of the s fourier-transformation ')
xlabel('frequency')

% avarage substraction
abs(fts(1))/N
mean(s)
sc=s-mean(s);

% ACF calculation

for(tau=1:1:N)
 acf(tau)=0;
 for(j=1:1:N-tau)
    acf(tau)=acf(tau)+sc(j)*sc(j+tau-1);
 end;
  acf(tau)= acf(tau)/(N-tau+1);
end;
plot(acf);

plot(acf(1:100));


% Spectrum energy calculations

N_ft=N;
ftacf=fft(acf,N_ft);

plot(abs(ftacf)/N_ft);

freq=(1:N_ft/2-1)/N_ft;
plot(freq, abs(ftacf(2:N_ft/2))/N_ft);

% Window generating

N_w=254
 for(i=1:1:N_w)
      w(i)=(1-abs(i-1)/(N_w-1));
end;

plot(w)

% ACF by window multiplication
for(i=1:1:N_w)
 acfw(i)=acf(i)*w(i);
end;

plot(acfw);

% Blackman-Tukey spectrum calculation

N_ft=N_w;
ftacfw=fft(acfw,N_ft);

plot(abs(ftacfw)/N_ft);

freq=(1:N_ft/2-1)/N_ft;

plot(freq, abs(ftacfw(2:N_ft/2))/N_ft);

%wavelet-transformation of the signal

c = cwt(s,[1:128],'morl','plot');

%3d-plot
[A,B] = meshgrid([1:100],[1:128]);
for(i=1:1:128)
 for(j=1:1:100)
   Z(i,j)=c(i,j*10);
 end;
end;

surface(B,A,Z)