%  AQUI ANALIZAMOS QUE PINTA TIENE EL SONIDO DE VARIOS DIAPASONES EN EL
%  DOMIINIO TEMPORAL Y EN EL DE FOURIER.
%
% NOTA: Todo lo de abajo asume que se han guardado los ficheritos en un
% determinado directorio que naturalmente en tu ordenador no existe (a no
% ser que estes escribiendo en el mio)
%

% Leemos los ficheritos que hemos bajado de http://www.xtec.es/centres/a8019411/caixa/ondas.htm

[Y27,FS,NBITS]=wavread('C:\Documents and Settings\jesusjesus\Mis documentos\jesus\clases\percepcion\27.wav');
[Y100,FS,NBITS]=wavread('C:\Documents and Settings\jesusjesus\Mis documentos\jesus\clases\percepcion\100.wav');
[Y200,FS,NBITS]=wavread('C:\Documents and Settings\jesusjesus\Mis documentos\jesus\clases\percepcion\200.wav');
[Y440,FS,NBITS]=wavread('C:\Documents and Settings\jesusjesus\Mis documentos\jesus\clases\percepcion\440.wav');
% Nunca fiarse de lo que se trae uno por internet, que son unos burros!!!
Y440=Y440(272:end-271);
[Y1000,FS,NBITS]=wavread('C:\Documents and Settings\jesusjesus\Mis documentos\jesus\clases\percepcion\1000.wav');
[Y3000,FS,NBITS]=wavread('C:\Documents and Settings\jesusjesus\Mis documentos\jesus\clases\percepcion\3000.wav');

Y27=Y27(1:21644);
Y100=Y100(1:21644);
Y200=Y200(1:21644);
Y1000=Y1000(1:21644);
Y3000=Y3000(1:21644);

% Normalizamos para igualar energias

Y27=Y27/max(Y27);
Y100=Y100/max(Y100);
Y200=Y200/max(Y200);
Y440=Y440/max(Y440);
Y1000=Y1000/max(Y1000);
Y3000=Y3000/max(Y3000);

% Vale, veamos las ondas...
duracion=length(Y27)/FS;

figure(1),plot(linspace(0,length(Y27)/FS,length(Y27)),Y27),xlabel('t (seg)'),ylabel('Presion (Bares)'),title('Sinusoide de 27 Hz (Dom. Temporal)'),axis([0 duracion/20 -1.3 1.3])
figure(2),plot(linspace(0,length(Y27)/FS,length(Y27)),Y100),xlabel('t (seg)'),ylabel('Presion (Bares)'),title('Sinusoide de 100 Hz (Dom. Temporal)'),axis([0 duracion/20 -1.3 1.3])
figure(3),plot(linspace(0,length(Y27)/FS,length(Y27)),Y200),xlabel('t (seg)'),ylabel('Presion (Bares)'),title('Sinusoide de 200 Hz (Dom. Temporal)'),axis([0 duracion/20 -1.3 1.3])
figure(4),plot(linspace(0,length(Y27)/FS,length(Y27)),Y440),xlabel('t (seg)'),ylabel('Presion (Bares)'),title('Sinusoide de 440 Hz (Dom. Temporal)'),axis([0 duracion/20 -1.3 1.3])
figure(5),plot(linspace(0,length(Y27)/FS,length(Y27)),Y1000),xlabel('t (seg)'),ylabel('Presion (Bares)'),title('Sinusoide de 1000 Hz (Dom. Temporal)'),axis([0 duracion/20 -1.3 1.3])
figure(6),plot(linspace(0,length(Y27)/FS,length(Y27)),Y3000),xlabel('t (seg)'),ylabel('Presion (Bares)'),title('Sinusoide de 3000 Hz (Dom. Temporal)'),axis([0 duracion/20 -1.3 1.3])

% Bueno, despues de los manejos, las ondas duran 1 seg. aprox.
%
% Por tanto si queremos amuermar al personal con musicas de 15 seg. las
% repetiremos 15 veces! (si cambias el 15 por el numero que te salga de los
% huevos, las haras tan largas como quieras).


Y=Y27;
for i=1:15
    Y=[Y;Y27];
end
Y27larga=Y;
clear Y

Y=Y100;
for i=1:15
    Y=[Y;Y100];
end
Y100larga=Y;
clear Y

Y=Y200;
for i=1:15
    Y=[Y;Y200];
end
Y200larga=Y;
clear Y

Y=Y440;
for i=1:15
    Y=[Y;Y440];
end
Y440larga=Y;
clear Y

Y=Y1000;
for i=1:15
    Y=[Y;Y1000];
end
Y1000larga=Y;
clear Y

Y=Y3000;
for i=1:15
    Y=[Y;Y3000];
end
Y3000larga=Y;
clear Y

% Vale, y ahora que ya la tenemos bien larga, nos la guardamos... y con
% esto jugaremos con el windows-media-player (que es una basura).

wavwrite(Y27larga,FS,NBITS,'C:\Documents and Settings\jesusjesus\Mis documentos\jesus\clases\percepcion\27larga.wav')
wavwrite(Y100larga,FS,NBITS,'C:\Documents and Settings\jesusjesus\Mis documentos\jesus\clases\percepcion\100larga.wav')
wavwrite(Y200larga,FS,NBITS,'C:\Documents and Settings\jesusjesus\Mis documentos\jesus\clases\percepcion\200larga.wav')
wavwrite(Y440larga,FS,NBITS,'C:\Documents and Settings\jesusjesus\Mis documentos\jesus\clases\percepcion\440larga.wav')
wavwrite(Y1000larga,FS,NBITS,'C:\Documents and Settings\jesusjesus\Mis documentos\jesus\clases\percepcion\1000larga.wav')
wavwrite(Y3000larga,FS,NBITS,'C:\Documents and Settings\jesusjesus\Mis documentos\jesus\clases\percepcion\3000larga.wav')

% Para ser melenudos rockeros como DIOS manda, usaremos el Matlab, es
% decir, haremos la TF BIEN!

TFY27=abs(fftshift(fft(Y27)));
TFY100=abs(fftshift(fft(Y100)));
TFY200=abs(fftshift(fft(Y200)));
TFY440=abs(fftshift(fft(Y440)));
TFY1000=abs(fftshift(fft(Y1000)));
TFY3000=abs(fftshift(fft(Y3000)));

% En una TF discreta (de ordenador) las frecuencias van entre -FS/2 y FS/2 
% Eso es la frecuencia de NYQUIST, que era un tio de los aņos 20, y FS es
% Sampling Frequency, que en catalan se dice Frequencia de Mostreig, i es
% la frequencia a la que mostregem(?) l'ona.
% 

fre=linspace(-FS/2,FS/2,length(Y27));

figure(7),plot(fre,TFY27),xlabel('Frecuencia (Hz)'),ylabel('|TF(Presion)|'),title('Sinusoide de 27 Hz (Dom. Fourier)'),axis([-6000 6000 0 12000])
figure(8),plot(fre,TFY100),xlabel('Frecuencia (Hz)'),ylabel('|TF(Presion)|'),title('Sinusoide de 100 Hz (Dom. Fourier)'),axis([-6000 6000 0 12000])
figure(9),plot(fre,TFY200),xlabel('Frecuencia (Hz)'),ylabel('|TF(Presion)|'),title('Sinusoide de 200 Hz (Dom. Fourier)'),axis([-6000 6000 0 12000])
figure(10),plot(fre,TFY440),xlabel('Frecuencia (Hz)'),ylabel('|TF(Presion)|'),title('Sinusoide de 440 Hz (Dom. Fourier)'),axis([-6000 6000 0 12000])
figure(11),plot(fre,TFY1000),xlabel('Frecuencia (Hz)'),ylabel('|TF(Presion)|'),title('Sinusoide de 1000 Hz (Dom. Fourier)'),axis([-6000 6000 0 12000])
figure(12),plot(fre,TFY3000),xlabel('Frecuencia (Hz)'),ylabel('|TF(Presion)|'),title('Sinusoide de 3000 Hz (Dom. Fourier)'),axis([-6000 6000 0 12000])