Projeto Análise de Fourier Resolvido 2023 1

1 de 1 Aluno:_____________________________________________________ Nota: Semestre:________ Data: _____/______/________ Turma: ________ Prof. Wagner Endo Disciplina: Processamento Digital de Sinais Visto: _______________ Projeto I O problema: Aquisição e análise de Fourier de um sinal de áudio (voz-fala). O desafio do problema: Analisar e compreender as características de um sinal de áudio em módulo e fase. Objetivos da análise do projeto: • Aquisitar um sinal de áudio diretamente pelo MatLab; • Plotar e apresentar os dados da análise de Fourier do sinal: módulo e fase (Dica: selecione trechos diferentes do sinal gravado para enriquecer esta análise); • Implementar um código em MatLab que modifique somente os valores de módulo, mantendo a fase; • Implementar um código em MatLab que modifique somente os valores de fase, mantendo a módulo; • Comparar os sinais de áudio obtidos e verificar os efeitos do módulo e da fase no sinal de áudio. • Elaborar um relatório final do projeto proposto; MINISTÉRIO DA EDUCAÇÃO UNIVERSIDADE TECNOLÓGICA FEDERAL DO PARANÁ CAMPUS CORNÉLIO PROCÓPIO PROCESSAMENTO DIGITAL DE SINAIS UNIVERSIDADE TECNOLÓGICA FEDERAL DO PARANÁ PR MINISTÉRIO DA EDUCAÇÃO UNIVERSIDADE TECNOLÓGICA FEDERAL DO PARANÁ CAMPUS CORNÉLIO PROCÓPIO PROCESSAMENTO DIGITAL DE SINAIS Aluno:_____________________________________________________ Semestre:________ Data: _____/______/________ Turma: ________ Prof. Wagner Endo Disciplina: Processamento Digital de Sinais Projeto I Aquisição e análise de Fourier de um sinal de áudio (voz-fala). 14 de maio de 2023 Curitiba - PR Introdução Inúmeras ferramentas matemáticas são largamente utilizadas para ajudar a interpretar os sinais. Uma dessas ferramentas foi desenvolvida por Fourier, sendo ela a transformada de Fourier. A transformada de Fourier consiste em observar um sinal no domínio do tempo e a partir dele ter uma representação matemática deste sinal no domínio das frequências. Desta forma é possível observar o sinal de acordo com as regiões de seu espectro de frequências, podendo este ser dividido em graves, médios e agudos para sinais sonoros, mas tal estudo se expande para todo sinal que possua frequências. Sendo assim, a transformada de Fourier é utilizada para análise e tratamento de imagens, transmissão e recepção de sinais sonoros, de vídeo, analise espectral de atômica entre outros temas. A transformada de Fourier foi adaptada para o caso discreto, onde as amostras são finitas e futuramente foi desenvolvido um algoritmo que trouxe velocidade ao cálculo desta transformada por máquinas, possibilitando assim ser calculada por computadores de forma rápida e eficiente, esta é a transformada rápida de Fourier (fft). O MATLAB é uma poderosa ferramenta matemática onde é possível analisar sinais utilizando de funções matemáticas, como é o caso da fft. Neste trabalho será analisado um sinal de áudio gravado com o próprio MATLAB e realizada uma análise a partir dos seus gráficos no tempo e na frequência em magnitude e fase. Desenvolvimento Aquisição do áudio diretamente pelo MATLAB Incialmente coleta-se o áudio no MATLAB. Para isso é necessário definir alguns parâmetros, como a frequência de amostragem do sinal, o número de bits de cada amostra, o número de canais de áudio e qual será o microfone utilizado para captar o áudio. A frequência de amostragem foi definida como 20kHz, isto indica que serão capturadas 20 mil amostras do áudio a cada segundo. Pelo teorema de Nyquist, poderemos analisar o sinal entre 0 e 10kHz. O número de bits para cada amostra foi definido como 16bits, desta forma, a amostra poderá receber 216 valores diferentes. Foi utilizado apenas 1 canal, sendo este um sinal de áudio mono. Desta forma iremos analisar apensas um sinal gravado. O microfone utilizado será o padrão decido pelo próprio sistema operacional. Com essas informações é possível criar o objeto de áudio que receberá o sinal. Fs = 20000; %Frequência de amostragem do sinal nBits = 16; %Quantidade de bits para cada amostra nChannels = 1; %Número de canais sonoros ID = -1; % Dispositivo de entrada de áudio padrão recObj = audiorecorder(Fs,nBits,nChannels,ID); %Objeto onde será salvo o áudio Para realizar a gravação do sinal de áudio, basta definir o tempo de duração do áudio, informar que iniciará a gravação, gravar e informar quando por encerrado a gravação. recDuration = 4; % Duração do áudio em segundos disp("Comece a falar.") %Inicia a gravar recordblocking(recObj,recDuration); %Grava o audio no objeto disp("Fim da gravação.") %Indica que a gravação terminou Com o objeto do áudio, podemos criar o array com os valores do áudio no tempo a partir da função getaudiodata(), definir o período de cada amostra como sendo o inverso da frequência e assim gerar o sinal de tempo. y = getaudiodata(recObj); %A partir do objeto de audio cria um array com seus dados no tempo Ts=1/Fs; %O período de cada amostra é o inverso da frequencia de amostragem t=0:Ts:recDuration-Ts; %O vetor de tempo inicia em zero, com passo igual ao período de amostragem e vai até o tempo de duração do %audio menos uma amostra. Para uma análise de vários trechos do sinal, podemos repartir o sinal de tempo em 4 partes distintas: y1=y(1:recDuration/4*Fs); %Primeiro quarto do sinal y2=y(recDuration*(1/4)*Fs+1:recDuration*(2/4)*Fs); %Segundo quarto do sinal y3=y(recDuration*(2/4)*Fs+1:recDuration*(3/4)*Fs); %Terceiro quarto do sinal y4=y(recDuration*(3/4)*Fs+1:length(y)); %Ultimo quarto do sinal zeros=zeros(length(y)/4,1); %Um array só com zeros de tamanho de um quarto do sinal y1=[y1; zeros; zeros; zeros]; %Novo y1 onde é inserido zeros em seus trechos sem sinal y2=[zeros; y2; zeros; zeros]; %Novo y2 onde é inserido zeros em seus trechos sem sinal y3=[zeros; zeros; y3; zeros]; %Novo y3 onde é inserido zeros em seus trechos sem sinal y4=[zeros; zeros; zeros; y4]; %Novo y4 onde é inserido zeros em seus trechos sem sinal Podemos visualizar tais sinais a partir dos seguintes plots: figure(1) subplot(2,1,1); plot(t,y) %Plota o sinal de áudio no tempo title("Sinal de áudio gravado no tempo"); xlabel("t(s)"); ylabel("Amplitude"); grid; subplot(2,4,5) plot(t,y1) %Plota o sinal de áudio no tempo title("Sinal 1"); xlabel("t(s)"); ylabel("Amplitude"); grid; subplot(2,4,6) plot(t,y2) %Plota o sinal de áudio no tempo title("Sinal 2"); xlabel("t(s)"); ylabel("Amplitude"); grid; subplot(2,4,7) plot(t,y3) %Plota o sinal de áudio no tempo title("Sinal 3"); xlabel("t(s)"); ylabel("Amplitude"); grid; subplot(2,4,8) plot(t,y4) %Plota o sinal de áudio no tempo title("Sinal 4"); xlabel("t(s)"); ylabel("Amplitude"); grid; Plot do módulo e fase do sinal Agora com a amostra do áudio no tempo, é necessário definirmos alguns parâmetros para obtermos o espectro do sinal na frequência. Inicialmente definimos a frequência máxima a ser lida, dada pelo teorema de Nyquist e o espectro de frequências será desde a frequência máxima negativa até ela positiva. 𝐹𝑚𝑎𝑥 = 𝐹𝑠 2 −𝐹𝑚𝑎𝑥 < 𝑓 < 𝐹𝑚𝑎𝑥 O passo é dado pelo inverso do tempo de duração do sinal: 𝑃𝑎𝑠𝑠𝑜 = 1 𝑇𝑑𝑢𝑟𝑎çã𝑜 Também definimos uma tolerância máxima para representar o valor zero para garantir assim que o gráfico de amplitudes não será poluído por sinais de baixa amplitude. fMax=(Fs/2); %A frequência máxima a ser analisada é igual a frequencia de amostragem sobre 2. tol = 5; %Tolerancia do sinal para coletar sinais de ângulo. f=-fMax:1/recDuration:fMax-1/Fs;%O vetor de frequencias inicia em - Fs/2 e vai até Fs/2 %com passo igual ao inverso do tempo de duração do sinal. Podemos agora calcular a fft do sinal e assim encontrarmos o vetor de valores complexos que representa este sinal na frequência. O MATLAB encontra a partir da fft os valores das frequências entre 0 e 2𝜋 do sinal. Como a analise só é realizada entre 0 e 𝜋 ou ainda entre −𝜋 e 𝜋 considerando as frequências negativas, então temos que utilizar a função fftshift() que altera o sinal de 0 a 2𝜋 para −𝜋 a 𝜋. Após isto, podemos calcular o modulo e a fase do sinal a partir das funções abs() e angle() respectivamente. Yw=fftshift(fft(y)); %Calcula a fft de y e desloca o sinal para - pi<w<pi YwMag=abs(Yw); %Calcula o módulo do sinal Yw(abs(Yw) < tol) = 0; %Define a tolerancia maxima para o calculo da fase YwAng=angle(Yw); %Calcula a fase do sinal Yw1=fftshift(fft(y1)); Yw1Mag=abs(Yw1); Yw1(abs(Yw1) < tol) = 0; Yw1Ang=angle(Yw1); Yw2=fftshift(fft(y2)); Yw2Mag=abs(Yw2); Yw2(abs(Yw2) < tol) = 0; Yw2Ang=angle(Yw2); Yw3=fftshift(fft(y3)); Yw3Mag=abs(Yw3); Yw3(abs(Yw3) < tol) = 0; Yw3Ang=angle(Yw3); Yw4=fftshift(fft(y4)); Yw4Mag=abs(Yw4); Yw4(abs(Yw4) < tol) = 0; Yw4Ang=angle(Yw4); Com isto podemos plotar os sinais no eixo das frequências e assim realizar a análise parte a parte do sinal: figure(2) subplot(2,3,1) plot(f,YwMag) title("Amplitude na frequência"); xlabel("f(Hz)"); ylabel("Amplitude"); grid; subplot(2,3,4) plot(f,YwAng) title("Angulo na frequência"); xlabel("f(Hz)"); ylabel("\theta (rad/s)"); grid; subplot(4,3,2) plot(f,Yw1Mag) title("Amplitude Sinal 1"); xlabel("f(Hz)"); ylabel("Modulo"); grid; subplot(4,3,3) plot(f,Yw1Ang) title("Angulo Sinal 1"); xlabel("f(Hz)"); ylabel("\theta (rad/s)"); grid; subplot(4,3,5) plot(f,Yw2Mag) title("Amplitude Sinal 2"); xlabel("f(Hz)"); ylabel("Modulo"); grid; subplot(4,3,6) plot(f,Yw2Ang) title("Angulo Sinal 2"); xlabel("f(Hz)"); ylabel("\theta (rad/s)"); grid; subplot(4,3,8) plot(f,Yw3Mag) title("Amplitude Sinal 3"); xlabel("f(Hz)"); ylabel("Modulo"); grid; subplot(4,3,9) plot(f,Yw3Ang) title("Angulo Sinal 3"); xlabel("f(Hz)"); ylabel("\theta (rad/s)"); grid; subplot(4,3,11) plot(f,Yw4Mag) title("Amplitude Sinal 4"); xlabel("f(Hz)"); ylabel("Modulo"); grid; subplot(4,3,12) plot(f,Yw4Ang) title("Angulo Sinal 4"); xlabel("f(Hz)"); ylabel("\theta (rad/s)"); grid; Podemos observar que o sinal de áudio possui sinais entre 0 e 3kHz. Ampliando o sinal para observar apenas esta faixa de frequências, temos: Com este gráfico podemos notar que o sinal possui algumas frequências com maior amplitude, sendo elas os 226.5Hz e 2399.5Hz. Isto indica que no sinal de áudio gravado essas frequências estavam bem presentes. Observando parte por parte do sinal, no primeiro segundo não aparece a frequência de 2399.5Hz, apenas a de 226.5Hz e algumas frequências um pouco maiores, mas não muito, isto indica que a componente aguda do sinal ainda não estava presente, mas algumas notas de frequência um pouco maior estavam. No segundo segundo do sinal também não aparece a frequência aguda de 2399.5Hz, mas nela as frequências um pouco maiores que 226Hz e próximas a 1kHz já não aparecem mais. No terceiro segundo a frequência de 2399.5Hz aparece junto a parte do sinal anterior. Já no último segundo só temos a presença do sinal de 2399.5Hz. Isto se dá pois inicialmente ao começar a gravar o áudio a minha voz possui diferentes frequências iniciais até a faixa dos 1kHz. Mas evidenciando a tonalidade da minha voz em algo próximo a 226.5Hz. No segundo segundo, após já ter iniciado a fala, o meu tom de voz se sobressai perante os ruídos do início da fala. A partir deste momento realizei um assobio que pode ser verificado como um sinal agudo que neste caso estava na frequência de 2399.5Hz e no quarto e último segundo apenas estava presente o assobio. No gráfico das fases ainda podemos ver uma variação muito grande de fase pois um sinal de voz possui inúmeras frequências trabalhando em fases distintas e mesmo com a tolerância fica difícil notar as variações na fase. Apenas no assobio por ser um sinal muito próximo de uma senoide que verificamos os impulsos na fase. Modificação do sinal de módulo Para realizar alterações no sinal de módulo, primeiro foi criado um array com valores aleatórios entre 0 e 1 para ser multiplicado pelas frequências entre 0 e 10Khz, assim trazendo pesos distintos para cada frequência. Para isto foi utilizado a função randi com valores entre 0 e 1000 e dividindo esses valores por 1000. Com este array, multiplicamos as frequências entre -10kHz e 0 por ele e entre 0 e 10kHz por ele invertido. random=(randi(1000,1,length(YwMag)/2)/1000)'; YwMagNovo=[YwMag(1:length(YwMag)/2).*random; YwMag(length(YwMag)/2+1:length(YwMag)).*flipud(random)]; A partir desta nova função podemos verificar nosso novo gráfico: figure(3) subplot(2,1,1) plot(f,YwMag) title("Modulo do sinal original"); xlabel("f(Hz)"); ylabel("Amplitude") axis([-4000 4000 0 250]) grid; subplot(2,1,2) plot(f,YwMagNovo) title("Modulo do sinal alterado"); xlabel("f(Hz)"); ylabel("Amplitude") axis([-4000 4000 0 250]) grid; Como podemos ver o sinal sofreu redução de amplitude em diversos setores inclusive em suas frequências principais. Para realizar a análise do sinal no tempo devemos repassar este sinal novamente para o tempo com o uso da ifft(), para isto, precisamos do sinal complexo novamente. 𝑦 = 𝐴(cos(𝜃) + 𝑖 ∙ 𝑠𝑒𝑛(𝜃)) YwNovo=YwMagNovo.*complex(cos(YwAng),sin(YwAng)); Com o vetor complexo podemos calcular a ifft, para isto voltamos o sinal para o formato que o MATLAB trabalha utilizando a fftshift, calculados a fft e pegamos aparte real do sinal. Podemos escutar o sinal original e o alterado com a função sound. ynovo=real(ifft(fftshift(YwNovo))); sound(y./max(y),Fs,16) pause(recDuration); sound(ynovo./max(ynovo),Fs,16) É notável a diferença em amplitude do sinal e como partes do sinal acabaram se misturando pela normalização em amplitude. Modificação do sinal de fase Multiplicamos as frequências entre -10kHz e 0 pelo array randômico e entre 0 e 10kHz por ele invertido. Assim temos: YwAngNovo=[YwAng(1:length(YwAng)/2).*random; YwAng(length(YwAng)/2+1:length(YwAng)).*flipud(random)]; A partir desta nova função podemos verificar nosso novo gráfico: figure(4) subplot(2,1,1) plot(f,YwAng) title("Fase do sinal original"); xlabel("f(Hz)"); ylabel("Fase") axis([-4000 4000 -4 4]) grid; subplot(2,1,2) plot(f,YwAngNovo) title("Fase do sinal alterado"); xlabel("f(Hz)"); ylabel("Fase") axis([-4000 4000 -4 4]) grid; Como podemos ver o sinal sofreu redução na fase em diversos setores. Para realizar a análise do sinal no tempo devemos repassar este sinal novamente para o tempo com o uso da ifft(), para isto, precisamos do sinal complexo novamente. YwNovo2=YwMag.*complex(cos(YwAngNovo),sin(YwAngNovo)); Com o vetor complexo podemos calcular a ifft, para isto voltamos o sinal para o formato que o MATLAB trabalha utilizando a fftshift, calculados a fft e pegamos aparte real do sinal. Podemos escutar o sinal original e o alterado com a função sound. ynovo2=real(ifft(fftshift(YwNovo2))); sound(y./max(y),Fs,16) pause(recDuration); sound(ynovo2./max(ynovo2),Fs,16) É notável que o som em alguns pontos se torna nulo e em outros é mais evidenciado pois algumas componentes são canceladas pelo ajuste de fase. Também é notável que a modificação na fase do sinal embaralhou partes do sinal onde o assobio se viu presente em mais trechos do áudio. Para observar podemos recorrer aos gráficos: figure(5) subplot(3,1,1) plot(t,y) title("Sinal de áudio gravado no tempo"); xlabel("t(s)"); ylabel("Amplitude"); grid; subplot(3,1,2) plot(t,ynovo) title("Sinal alterando módulo"); xlabel("t(s)"); ylabel("Amplitude"); grid; subplot(3,1,3) plot(t,ynovo2) title("Sinal alterando fase"); xlabel("t(s)"); ylabel("Amplitude"); grid; Podemos notar que as alterações no módulo do sinal influenciaram bastante na amplitude do sinal, já na frequência ele trouxe o sinal do assobio para toda a extensão do áudio e anulou partes por cancelamento de fase.