%% ==================================================================
%% VENTANA SOBRE UNA SEÑAL PARA DISMINUIR EL RIPPLE EN FRECUENCIA
%%
%% Medir una señal por un tiempo determinado es equivalente a 
%% truncarla mediante una ventana rectangular. Esto genera mucho
%% ripple es frecuencia. Por los tanto, para procesar un dato es 
%% conveniente en primer término truncar la señal mediante una 
%% ventana más suave. 
%% En este ejercicio se muestra cómo se reduce el ripple de un dato
%% de entrada al truncar la señal mediante una ventana hamming en 
%% lugar de una ventana rectangular.
%% ==================================================================

clear all;clear all; close all;clf

% --------------------------
% SENAL
% --------------------------

w0 = pi/4;    % frecuencia angular digital w0 = 2*pi*f*dt
T  = 2*pi/w0; % período digital

n = (0:1:6*T); % dt==1
s = sin(w0*n);

N = length(s);

% --------------------------
%% VENTANA
% --------------------------

%% a) dando la longitud total de la ventana

%Nw = 21                   % longitud total de la ventana: 2*Mw+1
%Mw = (Nw-1)/2;            % longitud de la parte ascendente/descendente

%% b) dando la longitud de la parte ascendente/descendente
Mw = 10
Nw = 2*Mw+1;


hamm = hamming(Nw);       % ventana de hamming ;)


% --------------------------
% SENAL TRUNCADA POR VENTANA
% --------------------------

% a) usando un vector ventana de dimension igual a la senal (didactico).
win  = [ hamm(1:Mw)'  ones(1,N-2*Mw) hamm(Mw:-1:1)' ];
%swin = s .* win;

% b) usando directamente los indices (mas sencillo).

swin=s;
swin(1:Mw)     = s(1:Mw)     .* hamm(1:Mw)' ;
swin(N-Mw+1:N) = s(N-Mw+1:N) .* hamm(Mw:-1:1)' ;


% --------------------------
% ESPECTROS DE AMPLITUD
% --------------------------

% zero padding

spad=s;
swinpad=swin;

Npad          = 2**10;
spad(Npad)    = 0;
swinpad(Npad) = 0;

dw=2*pi/Npad;

wpad = (0:(Npad-1))*dw;

S    = abs(fft(spad));
SWIN = abs(fft(swinpad));

% --------------------------
% GRAFICOS
% --------------------------

figure(1); clf

% En tiempo

subplot(3,1,1)
plot(s(1:N),'Linewidth',1.5)
axis([0,N-1,-1.2,1.2])
title('Senal')

subplot(3,1,2)
plot(win(1:N),'g','Linewidth',1.5)
axis([0,N-1,-0.05,1.2])
title('Ventana (Hamming)')

subplot(3,1,3)
plot(swin(1:N),'5','Linewidth',1.5)
axis([0,N-1,-1.2,1.2])
title('Senal con Ventana')


figure(2)

subplot(2,3,1)
hold on
stem(hamm(1:Mw),'g')
plot(s(1:Mw),'b-o')
axis([1,Mw,-1.2,1.2])
hold off
subplot(2,3,2)
plot(s(Mw+1:N-(Mw-1)),'b')
axis([Mw+1,29,-1.2,1.2])
subplot(2,3,3)
plot(s(Mw+1:N-(Mw-2)),1)
hold on
stem(hamm(Mw:-1:1),'g')
plot(s(N-(Mw-1):N),'b-o')
axis([1,Mw,-1.2,1.2])
hold off

subplot(2,3,4)
plot(swin(1:Mw),'5-o')
axis([1,Mw,-1.2,1.2])
subplot(2,3,5)
plot(s(Mw+1:N-(Mw-1)),'b')
axis([Mw+1,29,-1.2,1.2])
subplot(2,3,6)
plot(swin(N-(Mw-1):N),'5-o')
axis([1,Mw,-1.2,1.2])



figure(3)

% Ticks para los gráficos
xt      = 0:pi/2:2*pi;
xtlabel = ["0";"{/Symbol p}/2";"{/Symbol p}";"3/2{/Symbol p}";"2{/Symbol p}"];

subplot(2,2,1)
plot(s(1:N),'Linewidth',1.5)
axis([0,N-1,-1.2,1.2])
title('Senal')

subplot(2,2,2)
plot(wpad,S,'r','Linewidth',1)
axis([0,2*pi,-5,max(S)])
title('Espectro de amplitud - Senal')
set(gca,"xtick",xt,"xticklabel",xtlabel);

subplot(2,2,3)
plot(swin(1:N),'5','Linewidth',1.5)
axis([0,N-1,-1.2,1.2])
title('Senal con ventana')

subplot(2,2,4)
plot(wpad,SWIN,'r','Linewidth',1)
axis([0,2*pi,-5,max(SWIN)])
title('Espectro de amplitud - Senal con ventana')
set(gca,"xtick",xt,"xticklabel",xtlabel);



% ---------------------
% SALIDA
% ---------------------

print('ventana_senal.eps','-depsc');