clear all; close all; N = 1000; fp = 1500; t = (0:N-1)/fp; f1 = 150; f2 = 2f1; f3 = 3f1; f4 = 4*f1;

x = sin(2pif1t) + sin(2pif2t) + sin(2pif3t) + sin(2pif4t); q = max(x); y = (1.5/q)*x + 2; max(y); subplot(321); plot(t,y); xlabel('czas [s]'); ylabel('sygnal');

Nf = N; N21 = Nf/2+1; v = fft(x,Nf); w = abs(v); f = linspace(0,fp/2,N21); subplot(322); stem(f,w(1:N21)); xlabel('czestotliwosc [Hz]'); ylabel('modul widma');

M = 101; fg = [50,150]; fgu = fg/(fp/2); h = fir1(M-1,fgu,'stop'); th = (0:M-1)/fp; % subplot(323); % plot(th,h); % xlabel('czas [s]'); % ylabel('odpowiedz impulsowa');

Nf1 = 2^nextpow2(M); N21h = Nf1/2 + 1; vx = fft(h,Nf1); wx = abs(vx); fx = linspace(0,fp/2,N21h); subplot(323); plot(fx,wx(1:N21h)); xlabel('czestotliwosc [Hz]'); ylabel('modul transmitancji');

z = filter(h,1,y);

z = z - mean(z);

subplot(324) wz = abs(fft(z,N)); stem(f,wz(1:N21)); xlabel('czestotliwosc [Hz]'); ylabel('modul widma');

figure(2); spectrogram(z);