%--------------------------------------------------------------------------

% I N S T I T U T O P O L I T É C N I C O N A C I O N A L
% Unidad Profesional Interdisciplinaria de Ingeniería Campus Zacatecas
%--------------------------------------------------------------------------
% Unidad de aprendizaje: Procesamiento Digital de Señales
% Carrera: Ingeniería en Sistemas Computacionales
% Profesor: Dr. Teodoro Ibarra Pérez
% Alumno: Alonso Domínguez López
%==========================================================================
%------ Amplificacion de señales y ruido ---------
%==========================================================================
% Objetivo de la práctica: Aprender a generar señales con ruido
%--------------------------------------------------------------------------
% Nombre del archivo: e13_Domínguez_Alonso_5CM1.m
% Fecha de creación: 09 de octubre de 2023
% Software: Matlab R2022a
% Versión: 1.0
% Notas: comentarios adicionales que considere importantes
%--------------------------------------------------------------------------

%Konnichiwa

%De cada 6 comparaciones que haga, 4 las pueda identificar al 100 porciento, un
error
%de 33.33 por ciento
%clear
Fs = 100; %Frecuencia de muestreo
Ts = 1/Fs; %Periodo de muestreo
A1 = 5; %Amplitud
A2 = 10;
A3 = 15;
f1 = 1;
f2 = 2;
f3 = 3;
f4 = 4;
n = 0:Fs-1; %Tamaño del vector

%Señal 1
%[y1,Fs] = audioread('palabra.wav');
y11 = (A1+0.2)*sin(2*pi*f1*n*Ts)+(A1-0.3)*sin(2*pi*f2*n*Ts); %Componente de 1Hz y
de 2Hz
y12 = (A1-0.2)*sin(2*pi*f3*n*Ts)+(A1+0.1)*sin(2*pi*f4*n*Ts); %Componente de 3Hz y
de 4z
y1 = y11+y12;
%Señal 2
y21 = (A2+0.1)*sin(2*pi*(f1+8)*n*Ts)+(A2-0.3)*sin(2*pi*(f2+8)*n*Ts); %Componente de
9Hz y de 10Hz
y22 = (A2+0.3)*sin(2*pi*(f3+7)*n*Ts)+(A2+0.4)*sin(2*pi*(f4+8)*n*Ts); %Componente de
10Hz y de 12z
y2 = y21+y22;
%Señal 3
y31 = (A3-0.2)*sin(2*pi*(f1+16)*n*Ts)+(A3-0.2)*sin(2*pi*(f2+16)*n*Ts); %Componente
de 17Hz y de 18Hz
y32 = (A3+0.2)*sin(2*pi*(f3+16)*n*Ts)+(A3+0.4)*sin(2*pi*(f4+16)*n*Ts); %Componente
de 19Hz y de 20hz
y3=y31+y32;
%Transformada rápida de Fourier
a11 = fft(y11);%Primer coeficiente;
a12 = fft(y12);
a21 = fft(y21);
a22 = fft(y22);
a31 = fft(y31);
a32 = fft(y32);
%Se obtienen las apmitudes de las componentes
va11 = abs(a11);
va12 = abs(a12);
va21 = abs(a21);
va22 = abs(a22);
va31 = abs(a31);
va32 = abs(a32);
%Se normaliza el vector de amplitud a su valor original
norm_a11 = va11/(Fs/2);
norm_a12 = va12/(Fs/2);
norm_a21 = va21/(Fs/2);
norm_a22 = va22/(Fs/2);
norm_a31 = va31/(Fs/2);
norm_a32 = va32/(Fs/2);
%Se reduce a la mitad y se recorre a la izquierda
n = (0:length(n)/2-1);
%Se reduce a la mitad
norm_a11 = norm_a11(1:length(norm_a11)/2);
norm_a12 = norm_a12(1:length(norm_a12)/2);
norm_a21 = norm_a21(1:length(norm_a21)/2);
norm_a22 = norm_a22(1:length(norm_a22)/2);
norm_a31 = norm_a31(1:length(norm_a31)/2);
norm_a32 = norm_a32(1:length(norm_a32)/2);
%Filtrado de componentes
%Diseñar un filtro, hacemos una funcion que recibe un parámetro
b = 0.001;
f_a11= filtrar(norm_a11,b);
f_a12= filtrar(norm_a12,b);
f_a21= filtrar(norm_a21,b);
f_a22= filtrar(norm_a22,b);
f_a31= filtrar(norm_a31,b);
f_a32= filtrar(norm_a32,b);
%Dibujar la gráfica
plot(n,norm_a11+norm_a12,'ro',n,norm_a21+norm_a22,'bo',n,norm_a31+norm_a32,'go');
P = [norm_a11(:) norm_a12(:) norm_a21(:) norm_a22(:) norm_a31(:) norm_a32(:)];

%Algoritmo de Clusterización Deep Learning

%cluster1 = [f_a11; f_a12; f_a21; f_a22; f_a31; f_a32];

%Comparacion
%Señal 1
% y1p = (A1+0.2)*sin(2*pi*f1*n*Ts)+(A1-0.3)*sin(2*pi*f2*n*Ts); %Componente de 1Hz y
de 2Hz
% y2p = (A1-0.1)*sin(2*pi*f3*n*Ts)+(A1+0.1)*sin(2*pi*f4*n*Ts); %Componente de 3Hz y
de 4z
% y1 = y11+y12;
% %Señal 2
% y21 = (A2+0.2)*sin(2*pi*(f1+8)*n*Ts)+(A2-0.3)*sin(2*pi*(f2+8)*n*Ts); %Componente
de 9Hz y de 10Hz
% y22 = (A2-0.1)*sin(2*pi*(f3+7)*n*Ts)+(A2+0.1)*sin(2*pi*(f4+8)*n*Ts); %Componente
de 10Hz y de 12z
% y2 = y21+y22;
% %Señal 3
% y31 = (A3+0.2)*sin(2*pi*(f1+16)*n*Ts)+(A3-0.3)*sin(2*pi*(f2+16)*n*Ts);
%Componente de 17Hz y de 18Hz
% y32 = (A3-0.1)*sin(2*pi*(f3+16)*n*Ts)+(A3+0.1)*sin(2*pi*(f4+16)*n*Ts);
%Componente de 19Hz y de 20hz
% y3=y31+y32;
y1p = (A3+0.2)*sin(2*pi*(f1+16)*n*Ts)+(A1-0.3)*sin(2*pi*(f2+16)*n*Ts); %Componente
de 17Hz y de 18Hz
y2p = (A3-0.1)*sin(2*pi*(f3+16)*n*Ts)+(A1+0.1)*sin(2*pi*(f4+16)*n*Ts); %Componente
de 19Hz y de 20z
y = y1p+y1p;

a1 = fft(y);
va1=abs(a1);
norm_a1 = va1/(Fs/2);
norm_a1 = norm_a1(1:(length(norm_a1))/2);
norm_a1 = filtrar(norm_a1,b);

t=1;
k=1;
while(true)
[g,cluster] = kmeans(P',3);
disp(['Propuestas No. ', num2str(k),' :[',num2str(g(1)),',',num2str(g(2)), ...
',',num2str(g(3)),',',num2str(g(4)),',',num2str(g ...
(5)),',',num2str(g(6)),']']);
if g(1)==1 && g(2)==1 && g(3)==2 && g(4)==2 && g(5)==3 && g(6)==3
c1=cluster;
g1=g;
break;
end
k=k+1;
end

%comments
% function [salida] = filtrar(vector, sensibilidad, frecuenciaMuestreo)
% i = 1;
% while i <= frecuenciaMuestreo/2
% if vector(i) >= sensibilidad
% salida(i) = vector(i);
% else
% salida(i) = 0;
% end
% i = i + 1;
% end
% end
%subplot(4,1,1);
%stem(n,norm_a11);
%
% subplot(4,1,2);
% stem(n,norm_a32);
% %xlim([1 4000]);
%
% subplot(4,1,3);
% stem(n,a32);
% %xlim([1 500]);
%
% subplot(4,1,4);
% stem(n,a31);