Abrindo e filtrando arquivos de som do tipo "wav" - em Scilab e Python

Uma das frases que eu acho interessantes no filme 2001: Uma Odisseia no Espaço. O áudio apresenta algum ruído.

---

O objetivo desta postagem é mostrar como um arquivo de som do tipo "wav" pode ser aberto e manipulado em Python e em Scilab, mas é bom deixar claro que sou bem mais fluente em Scilab. Sobre arquivo no formato "wav" ver aqui. Os arquivos de teste foram escolhidos do filme 2001: Uma Odisseia no Espaço - ver aqui.

O que fazem os códigos? Basicamente, abrem um arquivo se som, mostram o sinal no tempo e na frequência (espectro - uso da fft), filtram o sinal com um filtro do tipo passa-baixa Butterworth e mostram o resultado. No código Scilab, o trecho final cria e salva um tom senoidal com um pouco de ruído.

Principais comandos Python usados:

• convolve - realiza a operação de convolução - ver aqui. 
• fft - calcula a transformada rápida de Fourier - ver aqui. 
• wave - biblioteca para manipular arquivos wav - ver aqui.  

---

Código Scilab:

///////////////////////////////////////////////////////////////
// Ler um arquivo .wav e mostrar no     //
// tempo e na frequência, filtrar depois //
///////////////////////////////////////////////////////////////
// Prof. Francisco J. A. de Aquino        //
// ver. 1.0   -   15/04/2017,  Fortaleza  //
//////////////////////////////////////////////////////////////

clc; xdel(winsid()); // fechando as janelas.

// Frases:
frase = 'This conversation can serve no purpose anymore. Goodbye.';
//frase = "We are all – by any practical definition of the words – foolproof and incapable of error."
//frase = 'Tom de teste';
//frase = "I am a HAL 9000 computer.";
//frase = "My mind is going... I can feel it.";

// Arquivos de som:
wavfile = 'C:\Users\Silvestre\Downloads\goodbye1.wav';
//wavfile = 'C:\Users\Silvestre\Downloads\foolproof.wav';
//wavfile = 'C:\Users\Silvestre\Downloads\testeff.wav';
//wavfile="C:\Users\Silvestre\Downloads\hal_9000.wav";
//wavfile="C:\Users\Silvestre\Downloads\mind_is_going.wav";

// Lendo e mostrando o arquivo de som:
[y,Fs,bits]=wavread(wavfile);
fg0=figure();
h0=uicontrol(fg0,"style","pushbutton","string","|>",'callback','sound(y,Fs)');
plot(y);  title(frase);

// Espectro do sinal:
tx = round(max(size(y))/2 - 1);
f = 1:max(size(y)); f = f - 1; f = f/max(f); f = f*Fs;
yf = abs(fft(y));
figure; plot(f(1:tx),yf(1:tx)); title('Espectro');

// Filtro passa-baixa:
[p,z,g]=iir(4,'lp','butt',[.20 .40],[.08 .03]);
n1 = convol([1 -z(1)],[1, -z(2)]);
n2 = convol([1 -z(3)],[1, -z(4)]);
nz = convol(n1,n2);

p1 = convol([1 -p(1)],[1, -p(2)]);
p2 = convol([1 -p(3)],[1, -p(4)]);
pp = real(convol(p1,p2));

// Sinal filtrado:
yff = g*filter(nz,pp,y);
yff = yff/max(abs(yff));

fg=figure();
h=uicontrol(fg,"style","pushbutton","string","Som",'callback','sound(yff,Fs)');
plot(yff); title(frase);

yfx = abs(fft(yff));
figure; plot(f(1:tx),yfx(1:tx));

/////// Gerando um tom:
p2 = 2*%pi;
frq = Fs; dt = 1/frq;
t=0:dt:2;
msom = sin(p2*t/2);
som = sin(p2*t*440) + cos(p2*t*880);
som = som.*msom;
som = som + 0.03*rand(1,max(size(som)),'n');
som = som/max(abs(som));
figure; plot(t,som); title('Teste');
wavfiles = 'C:\Users\Silvestre\Downloads\testeff.wav';
wavwrite(som, Fs, 16, wavfiles);

---

Código Python:

###########################################

## Lendo um arquivo .wav                  #

## Vendo o sinal no tempo e na frequência #

## Filtrando por um filtro passa-baixa    #

###########################################

## Prof. Francisco J. A. de Aquino        #

## Ver. 1.0 - 15/04/2017, Fortaleza - Ce  #

###########################################
import matplotlib.pyplot as plt
import numpy as np
import wave as fw
import winsound

arqff = 'C:/Users/Silvestre/Downloads/hal_9000.wav';
#arqff = 'C:/Users/Silvestre/Downloads/goodbye1.wav';

#arqff = 'C:/Users\Silvestre/Downloads/foolproof.wav';

#arqff = 'C:/Users\Silvestre/Downloads/testeff.wav';

#arqff = 'C:/Users/Silvestre/Downloads/hal_9000.wav';

#arqff = 'C:/Users/Silvestre/Downloads/mind_is_going.wav';
arquivoWav = fw.open(arqff, 'r');

print("Número canais: ", arquivoWav.getnchannels());
print("Número bytes: ", arquivoWav.getsampwidth());
print("Taxa de amostragem: ", arquivoWav.getframerate());
print("Número de frames: ", arquivoWav.getnframes());
print("Compactação: ", arquivoWav.getcompname());

tipos = np.uint8;
Damp = 256;
if arquivoWav.getsampwidth()>1:
    tipos = np.int16;
    Damp = 32760;

#frames = arquivoWav.readframes(arquivoWav.getnframes());

frames = arquivoWav.readframes(-1);
Amplitude = np.fromstring(frames, tipos)/Damp;
fim = np.size(Amplitude);
fim2 = round(fim/2 - 1);

dt = 1.0 / arquivoWav.getframerate();
tempo = np.arange(0, arquivoWav.getnframes() * dt, dt);
ffq = np.arange(0,arquivoWav.getnframes());
ffq = arquivoWav.getframerate()*ffq/max(ffq);
winsound.PlaySound(arqff, winsound.SND_FILENAME);
arquivoWav.close();

plt.plot(tempo, Amplitude); plt.grid();
plt.title(arqff);

Ampf = abs(np.fft.fft(Amplitude));
plt.figure(); plt.plot(ffq[1:fim2],Ampf[1:fim2]);

# filtragem com filtro do tipo FIR:
hh = [1.0, 4.7820952, 9.0600752, 8.0167949, 1.952694, -2.4130316, -2.0234161, 

0.1735586, 1.0120036, 0.3765159, -0.3010215, -0.3118091, -0.0008764, 0.1450086,  

0.0661129, - 0.0376649, -0.0478969,  -0.0041389, 0.0204601];
Ampc = 0.0465829*np.convolve(hh,Amplitude,'same');
plt.figure(); plt.plot(tempo,Ampc); plt.title('Sinal filtrado'); plt.grid();

Ampfc = abs(np.fft.fft(Ampc));
plt.figure();
plt.plot(ffq[1:fim2],Ampfc[1:fim2]);
plt.title('Espectro do sinal filtrado'); plt.grid();

plt.show();

---

Exemplo de filtragem:

Sinal original com ruído e espectro; sinal filtrado e espectro.