Como Utilizar a FFT com o PIC32MX - DSP Library

Neste post você encontra informações de como utilizar a Fast Fourier Transformer (FFT) em seu projeto com PIC32MX.

A família de microcontroladores PIC32MX não possui instruções de DSP. Operações de DSP estão disponíveis na família PIC32MZ. Aqui iremos apresentar a utilização da função mips_fft16 ( ) disponível na biblioteca DSP Library para o PIC32MX. Utilizaremos a linguagem C com o compilador MPLAB XC32 no ambiente de desenvolvimento MPLAB X.

Não entraremos em detalhes da teoria da Transformada Rápida de Fourier - FFT. Partiremos direto para a utilização da função mips_fft16 ( ). Para saber mais sobre a FFT você pode acessar estes posts: "Como Utilizar a FFT com o dsPIC - DSP Library" e "FFT-PIC - Um Analisador de Espectro para a Rede de Energia Elétrica ".

No início do nosso projeto devemos incluir os headers "dsplib_dsp.h" e "fftc.h" da biblioteca DSP, conforme mostra o Código 1. O arquivo "fftc.h" contém os coeficientes utilizados nos cálculos da FFT.

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#include <p32xxxx.h>
#include <xc.h>         // comentar esta linha para usar o MPLAB C32
#include <plib.h>

#include "dsplib_dsp.h"
#include "fftc.h"

Código 1 - Incluindo headers da DSP Library ao projeto.

Para o nosso projeto temos que definir o número de "pontos" da FFT, conforme o Código 2.

/*         FFT com 128 pontos       */

#define N    128        // N = 2^7 = 128

#define fftc fft16c128  // from fftc.h, for N = 128

Código 2 - Define 128 pontos da FFT.

128 amostras no domínio do tempo (lidas pelo conversor analógico-digital, por exemplo), se transformarão em 64 pontos no domínio da frequência.

Dados da FFT deste exemplo:
- Frequência de amostragem (Fs): 7,68 ksps;
- Pontos da FFT: 128;
- log2 (N) = 7;
- Passo entre as frequências (Fo): 60 Hz.

A frequência de amostragem (Fs) do sinal de entrada irá determinar o espaçamento entre as frequências de cada ponto de saída (Fo).

A frequência de amostragem (Fs) do nosso sinal gerado (a onda quadrada) é de 7680 Hz, portanto Fo é 60 Hz, ou seja, 7680 sps / 128 pontos da FFT!

Isto quer dizer que temos na saída o valor do espectro das frequências múltiplas de 60 Hz, de 0 à 3840 Hz (60 Hz * 64 pontos no domínio da frequência).

Para o sinal de entrada que será aplicado a FFT iremos utilizar um vetor com dados de uma onda quadrada de 780 Hz e frequência de amostragem 7680 Hz. Esta onda foi gerada com o auxílio do Software dsPICworks, da Microchip. Para mais informações de como gerar este vetor de dados acesse "Como Utilizar a FFT com o dsPIC - DSP Library". O vetor sinal_Entrada [ ] do tipo int16, com 128 amostragens é visto no Código 3.

int16 sinal_Entrada [128]= {                   //vetor que ira armazenar os dados de entrada
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,        //valores iniciais declarados para simulacao
0x8001, 0x8001, 0x8001, 0x8001, 0x8001,        //e teste do programa
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,        //square signal 780 Hz, 7680 sps, -1.0 a +1.0
0x8001, 0x8001, 0x8001, 0x8001, 0x8001,
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,
0x8001, 0x8001, 0x8001, 0x8001, 0x8001,
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,
0x8001, 0x8001, 0x8001, 0x8001, 0x8001,
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,
0x8001, 0x8001, 0x8001, 0x8001, 0x8001,
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,
0x8001, 0x8001, 0x8001, 0x8001, 0x8001,
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,
0x8001, 0x8001, 0x8001, 0x8001, 0x8001,
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,
0x8001, 0x8001, 0x8001, 0x8001, 0x8001,
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,
0x8001, 0x8001, 0x8001, 0x8001, 0x8001,
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,
0x8001, 0x8001, 0x8001, 0x8001, 0x8001,
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,
0x8001, 0x8001, 0x8001, 0x8001, 0x8001,
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,
0x8001, 0x8001, 0x8001, 0x8001, 0x8001,
0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF, 0x7FFF,
0x8001, 0x8001, 0x8001};

 

Código 3 - Vetor com dados de entrada para a simulação e cálculo da FFT.

 A seguir precisamos declarar os vetores e variáveis utilizadas para o cálculo da função mips_fft16 ( ), conforme o Código 4.

    int log2N = 7;      // log2(128) = 7        N = 128 = (1<<log2N)     N = 2^7 = 128
    int16c din[N];
    int16c dout[N];
    int16c scratch[N];
    int16c *p_fftc = &fftc;

Código 4 - Declaração das variáveis utilizadas no cálculo da FFT.

Os vetores utilizados pela função mips_fft16 ( ) são do tipo int16c, ou seja, inteiro 16-bit complexo, com partes real e imaginária.

No início da função main ( ) iremos transferir os coeficientes para o cálculo da FFT, conhecidos como Twiddle Factors, da memória de programa para a memória RAM, para execução mais rápida do cálculo. Veja o Código 5.

    p_fftc = &fftc;                      //inicializa ponteiro

    for (i = 0; i < N / 2; i++)          //move twiddleFactors da flash para a RAM
    {
        twiddleFactors [i] = *p_fftc;
        *p_fftc++;                       //incrementa o ponteiro
    }

Código 5 - Movendo os coeficientes para a memória RAM.

Agora podemos partir para um loop while onde continuamente o PIC32 irá mover os dados contidos no vetor sinal_Entrada [ ] para o vetor din [ ] e limitar a escala entre -0.5 e +0.5. Note que sinal_Entrada [ ] pode conter dados vindos de um conversor analógico-digital ADC, por exemplo, e din [ ] é o vetor com dados complexos utilizados nos cálculos da FFT.

   while (TRUE) {

        //move buffer do sinal de entrada (real)para vetor da FFT (complexo)
        for (i = 0; i < N; i++)
        {
            din[i].re = sinal_Entrada [i];  //parte real
            din[i].im = 0;                  //parte imaginaria
            din[i].re = din[i].re >> 1;     //divide pela metade, escalona para -0.5 a +0.5
        }
          
        ...

  

Código 6 - Preparando o vetor de dados complexos de entrada para os cálculos da FFT.

A chamada para a função mips_fft16 ( ), que aplica a Transformada Rápida de Fourier ao sinal, é mostrada no Código 7. Detalhes sobre os argumentos desta função você encontra em "PIC32 DSP Library".

        // Executa o Calculo da FFT
        mips_fft16 (dout, din, twiddleFactors, scratch, log2N);

Código 7 - Cálculo da FFT.


O resultado da FFT está contido no vetor complexo dout [ ].

Para exibir o resultado gráfico da FFT comumente se computa o single sided FFT, a magnitude quadrática do espectro de frequência, ou o módulo do número complexo. O Código 8 apresenta as fórmulas que devem ser aplicadas a cada posição do vetor de dados.

        /*     Computa o single sided FFT    */
        // Calcula a soma do quadrado das partes real e imaginária transformando vetor complexo em um vetor real.
  
        square_Magnitude = parte_Real^2 + parte_Imaginaria^2
        

        /*          Calcula o Modulo         */
        // Extrai a raiz quadrada da soma dos quadrados da parte real e imaginaria
        
        modulo = sqrtf(parte_Real^2 + parte_Imaginaria^2)  

        
  

Código 8 - Cálculo da magnitude quadrática do resultado da FFT.

Por fim, podemos calcular o espectro da frequência de maior magnitude do sinal amostrado, conforme exemplificado no Código 9.

  
       // Encontra a frequency Bin ( indice da frequencia) que possui a maior energia
        picoFrequenciaBin = vectorMaxBin (PONTOS_FFT/2,squareMagnitude);

        // Calcula a frequencia em Hz com maior espectro
        picoFrequencia = picoFrequenciaBin*(SAMPLING_FREQ/PONTOS_FFT);

Código 9 - Cálculo do espectro da frequência de maior magnitude.

O programa agora pode ser testado.


Para confirmar se o seu projeto está correto, verifique se a variável picoFrequencia apresenta o valor de 780 [Hz] após os cálculos. Esta é a frequência da onda quadrada gerada no dsPICwork para a simulação dos cálculos da FFT.

A posição do vetor com o espectro de frequência de maior valor é SquareMagnitude [13].

Veja os cálculos finais:

Fs               /    PONTOS_FFT              =      Fo
7680 sps    /     128 pontos da FFT         =     60 Hz

60 Hz * [13]       =        780 Hz!


O Vídeo 1 mostra este código rodando no MPLAB X em modo simulação passo-a-passo.

Vídeo 1 - Simulação dos cálculos da FFT com o PIC32MX.

O Vídeo 2 mostra uma aplicação da FFT com o PIC32MX para a identificação da frequência de um sinal. O display gráfico mostra os dados do vetor squareMagnitude [ ] como resultado da FFT, o espectro das frequências.

Vídeo 2 - Aplicação da FFT para identificar a frequência de um sinal.

Para quem gosta de música e áudio, o Vídeo 3 mostra uma aplicação bem legal para a FFT: um Analisador de Espectro de Áudio para a animação de um MP3 Player com o PIC32MX!

 Vídeo 3 - PIC32 MP3 Player & Audio Spectrum Analyser com a FFT.

Foi utilizada a função mips_fft16 para o cálculo da FFT com 128 pontos.

A FFT é testada com o arquivo de áudio Sine3kHz.MP3 contendo uma onda senoidal variando de 20 Hz à 3 kHz.

Referências:

HAYKIN, Simon S.; VAN VEEN, Barry. Sinais e sistemas. Porto Alegre: Bookman, 2001. 668p.

MICROCHIP. Code Examples. CE018 - Using the Fast Fourier Transform (FFT) for Frequency Detection. Disponível em: http://ww1.microchip.com/downloads/en/DeviceDoc/CE018_FFT_DSPlib_101007.zip

MICROCHIP. Software Libraries. Disponível em: http://www.microchip.com/SoftwareLib.aspx

MICROCHIP. PIC32 DSP Library. Disponível em: http://www.microchip.com/stellent/idcplg?IdcService=SS_GET_PAGE&nodeId=2680&dDocName=en552827

FFT-PIC - Um Analisador de Espectro para a Rede de Energia Elétrica

 Uma aplicação da Transformada Rápida de Fourier, dsPIC e o Processamento Digital de Sinal

 Muitos estudantes de Engenharia se deparam no decorrer de seu curso com as fórmulas de Euller, as Séries e as benditas Transformadas de Fourier. Muitos se perguntam pra que aprender tantos cálculos? Vou utilizar isso algum dia? A ideia básica deste projeto foi mostrar de forma prática e didática uma aplicação da Transformada Rápida de Fourier.

Antes de descrever o funcionamento do circuito, vamos relembrar um pouco da Transforma de Fourier. Basicamente um sinal pode ser representado de duas formas: uma no domínio do tempo e outra no domínio da frequência. O sinal no domínio do tempo pode ser obtido naturalmente como por exemplo um sinal sonoro, uma música, etc. Em algumas ocasiões a análise de um sinal se torna melhor compreendida se analisada no domínio da frequência, como as harmônicas presentes em uma rede elétrica.

            A Transformada de Fourier (FT) é utilizada para converter uma função ou sinal do domínio do tempo para o domínio da frequência e vice-versa.

A Transformada Discreta de Fourier (DFT) é aplicada para sinais discretos no tempo (dados coletados por um conversor analógico digital por exemplo).

A Transformada Rápida de Fourier (FFT) foi desenvolvida para uma execução mais rápida dos cálculos da FT possibilitando o seu uso na computação. Dentre as suas principais aplicações estão a filtragem digital, reconhecimento de padrões, bargraph para sinais de áudio, equalizadores, eliminações de ruído e interferência de imagens.

A FFT representa a soma de uma série de ondas senoidais de diferentes frequências, fases e amplitudes:

 Onde:

 

N= 2^n e n pode ser escrito como N = 2M, sendo M um inteiro positivo. 

                Para entender a diferença do sinal analisado no domínio do tempo e no domínio da frequência vamos a um exemplo prático utilizando o software dsPICworks, da Microchip. Com ele geramos alguns sinais. 

                Na Figura 1a temos um sinal senoidal de 60 Hz (domínio do tempo, eixo x do gráfico).

 Figura 1a - Sinal Senoidal de 60 Hz no tempo

              Na Figura 1b um sinal senoidal de 600 Hz (domínio do tempo).  

  Figura 1b - Sinal Senoidal de 600 Hz no tempo

               Efetuando a soma dos sinais senoidais A e B obtemos o sinal C mostrado na Figura 1c (domínio do tempo).

 Figura 1c - Soma dos Sinais Senoidais de 60 Hz e 600 Hz no tempo

              O resultado da Transformada Rápida de Fourier aplicada ao sinal C é exibido na Figura 1D. O sinal agora é visto no domínio da frequência (eixo x do gráfico). Observe 2 picos no gráfico, que representam a energia do espectro de 60 Hz e de 600 Hz.

 Figura 1d - Resultado da FFT aplicado ao sinal

             

          A análise de um sinal no domínio da frequência se torna muito interessante e/ou mais amigável conforme a aplicação, por exemplo, em um VU Meter, equalizador gráfico de áudio ou na análise das frequências harmônicas da rede de energia elétrica (múltiplas de 60 Hz), que aliás, é o foco deste nosso projeto de exemplo para a aplicação da FFT.

O desenvolvimento da linha dsPIC30F da Microchip possibilitou a utilização de algoritmos complexos para os projetos devido à alta capacidade de processamento de dados (30 MIPS – Mega instruções por segundo) aliada às instruções de DSP que otimizam o processamento digital de sinais e facilitam a execução de programas em tempo real.

O circuito proposto faz uso de um microcontrolador & DSP (Controlador Digital de Sinal -DSC) dsPIC30F4013 que coleta um sinal na faixa de 0 à 3780 Hz, aplica a FFT, calcula a frequência de maior amplitude e mostra um gráfico com o espectro de frequência do sinal, com as 63 múltiplas harmônicas de 60 Hz.

O programa faz uso da biblioteca “DSP Library” da Microchip que possui os algoritmos necessários para a execução da FFT.

A figura 2 mostra o diagrama do circuito. São utilizados poucos componentes externos além do dsPIC e do display LCD. O circuito funciona independente do display e da comunicação serial, ou seja, o leitor pode montar somente a parte do display sem a comunicação serial (MAX232), ou vice-versa, conforme os recursos que tiver.

 

  Figura 2 – Diagrama do Analisador de Espectro FFT-PIC

Funcionamento do circuito

O sinal de entrada é amostrado pela entrada analógica AN2 do dsPIC em uma frequência de amostragem de 7,68 Khz (fs), lembrando que o máximo que a entrada suporta é 5 volts. O capacitor C1 faz o acoplamento do sinal e os resistores R3 e R4 centralizam o sinal em 2.5 Vcc. O ADC de 12 bits é configurado para mostrar o sinal no formato fracionário com sinal, isto é, 0 - 5 V é representado numa escala de -1,0 à 1,0.

È executado o algoritmo da FFT sobre 128 amostras do sinal no tempo, obtendo a amplitude de 64 frequência do sinal. A frequência de amostragem e o número de ‘pontos’ amostrados foi calculado para se obter as frequências múltiplas de 60 Hz (fo), a qual é de interesse em muitas áreas da elétrica, principalmente para a análise de interferências harmônicas nas redes de energia elétrica.

A frequência de maior espectro é mostrada em um display LCD 16X2 juntamente com o seu módulo. As amplitudes das frequências são enviadas para um computador através da saída UART2, o MAX232 converte os níveis TTL para RS232 para serem enviados para a porta serial do PC. O programa Analisador de Espectro FFT-PIC desenvolvido em Visual Basic recebe os dados do dsPIC e mostra-os na forma de gráfico de barras para a visualização do espectro de frequências do sinal amostrado.

O programa do dsPIC

O código que roda no dsPIC foi desenvolvido em linguagem C no ambiente de desenvolvimento MPLAB IDE juntamente com o compilador MPLAB C30 da Microchip. O programa inicia na função Main () com a configuração dos PORTS de I/O e inicialização das variáveis. Os periféricos são configurados e iniciados com as funções:

-  uart_init (), inicia a comunicação UART2/RS232 a 9600 bps, 8 bits de dados, sem paridade e 1 bit de parada.

 -init_adc12 (), inicia o conversor analógico digital.

- init_timer1 (), configura a interrupção do Timer1 para uma amostragem de 7,68 ksps com o ADC.

- lcd_init (), configura o display LCD para comunicação com 4 bits de dados.

São enviadas as mensagens iniciais para o display com a função envia_mensagem (TXPtr), o programa passa então para um loop infinito aguardando o fim das 128 amostragens do sinal de entrada. Estes dados são escalonados para -0,5 a 0,5 e submetidos à Transformada Rápida de Fourier pela função FFTComplexIP (). A ordem dos dados de saída, já no domínio da freqüência são organizados pela função BitReverseComplex (), o sinal complexo é convertido para real através de SquareMagnitudeCplx (). A função VetorMax () retorna o índice da frequência de maior espectro. Estas últimas 4 funções pertencem à biblioteca DSP Library. É calculada a frequência de maior espectro e enviada para o LCD assim como o seu módulo.

A função transmite_dados () envia os dados das 64 frequências para o computador, são enviados primeiro a parte inteira, depois a parte decimal e 3 bytes para controle da comunicação.

O programa faz uso do temporizador WDT configurado para 4096 ms. Um LED, ligado no pino RF0, pisca intermitente quando o programa está sendo executado normalmente.

A figura 3 mostra o protótipo do circuito amostrando um sinal de 1020 Hz.

 

 Figura 3 - Circuito do FFT-PIC

O programa do Analisador de Espectro FFT-PIC

O programa do Analisador de Espectro FFT-PIC foi desenvolvido no Visual Basic 6.0 e testado no Windows XP SP2. Com ele é possível visualizar o espectro das frequências do sinal que foi amostrado pelo dsPIC.

A comunicação com o DSC é feita através do controle MSComm que acessa a porta COM1 do computador. A comunicação serial está configurada para 9600 bps, 8, N, 1.

O formulário possui menus e botões para conectar e desconectar a porta serial, com indicação de status da comunicação e um indicador de erros nos dados.

O pacote de dados recebido do dsPIC é mostrado em um gráfico de barras. Para gerar o gráfico foi utilizado o objeto ChartSpace, pertencente ao Office Web Components. 2 combobox permitem ajustar os eixos X e Y do gráfico.

O programa inicia com a rotina Load_Form () que carrega o formulário e inicia o gráfico. Cada byte recebido pela serial chama o evento MSComm1_OnComm (). Quando o pacote de dados for recebido a rotina carregavalorFFT () atualiza o gráfico com os novos valores recebidos dos espectros das frequências, numa escala de amplitude de 0 à 100, e a frequência de maior módulo. Uma caixa de texto permite visualizar os dados que são recebidos.

Na figura 4 pode ser visto o programa com um sinal amostrado de 1020 Hz.

 Figura 4 - Programa Analisador de Espectro FFT-PIC desenvolvido no VB

O circuito pode ser testado utilizando um programa gerador de sinais para PC, que pode ser encontrado na internet. Para teste dos cálculos da FFT da DSP Library foi gerado um sinal senoidal variando de 0 à 3780 Hz em 15 segundos na saída de som para fone de ouvido do computador, lembrando que a entrada analógica do dsPIC suporta de 0 à 5 Vcc.A seguir você pode assistir ao vídeo mostrando o circuito sendo testado.
 

 

 Vídeo Teste do Circuito do FFT-PIC com um Gerador de Sinal

As linhas dos programas do dsPIC e do VB foram comentadas para facilitar na compreensão do código. Todos os cálculos referentes às configurações dos registradores da comunicação UART, do conversor A/D, da interrupção do temporizador Timer1 e do temporizador do WDT estão escritos junto às suas funções no programa principal ‘main.c’ . O arquivo ‘fuses.c’ contém a programação dos bits de configuração do dispositivo como a ativação do oscilador interno, a multiplicação do clock (PLL) por 16 para a operação a 30 MIPS e atribuição do tempo do WDT.

Este software foi desenvolvido para mostrar um exemplo das diversas aplicações para a Transformada de Fourier e os recursos oferecidos pelo dsPIC, e pode ser utilizado para detectar a frequência de um sinal (em passos de 60 Hz) e também para visualizar o seu espectro de frequência, sendo muito útil para analisar componentes de interferência e harmônicas em um sinal.

   
ATENÇÂO: A entrada analógica do circuito NÂO deve ser conectada diretamente à rede de energia elétrica ou tomada pois não suporta 127 V / 220 V. Esta tensão deve ser abaixada por um transformador para uma tensão compatível com a entrada do circuito: 0 à 5 Vcc, conforme identificado no diagrama elétrico.

 
Referência Bibliográfica

www.microchip.com

MACORATTI. J. C. Desenvolvimento de programas em VB. Disponível em:  http://www.macoratti.net

SIMON S. HAYKIN, BARRY VAN VEEN. Sinais e Sistemas.

MICROCHIP. Code Examples. CE018 - Using the Fast Fourier Transform (FFT) for Frequency Detection. Disponível em: http://ww1.microchip.com/downloads/en/DeviceDoc/CE018_FFT_DSPlib_101007.zip

 * Projeto publicado na revista Saber Eletrônica 448.