PIC32: ¿Cómo usar la biblioteca DSP?

Estoy usando una placa Microstick II con el PIC32MX250F128B incluido y me gustaría realizar un procesamiento de audio en tiempo real. Uso una entrada analógica para obtener el sonido de entrada y una salida PWM, en combinación con un filtro de paso bajo seguido de un amplificador, como mi salida de audio (conectada a auriculares).

Lo que me gustaría hacer es realizar un procesamiento de señal más avanzado, como filtrado, etc. El sitio web de Microchip dice que hay una biblioteca DSP dedicada a eso. El manual de lib se puede encontrar aquí , pero no contiene mucha información, solo algunos prototipos de funciones...

En una primera vez, me gustaría realizar un filtrado de paso bajo y paso de banda. Hasta ahora, logré realizar fft usando el siguiente código:

// include header files
#include <plib.h>
#include <p32xxxx.h>
#include <dsplib_dsp.h>
#include <fftc.h>

// Config Bits
#pragma config FNOSC = FRCPLL       // Internal Fast RC oscillator (8 MHz) w/ PLL
#pragma config FPLLIDIV = DIV_2     // Divide FRC before PLL (now 4 MHz)
#pragma config FPLLMUL = MUL_20     // PLL Multiply (now 80 MHz)
#pragma config FPLLODIV = DIV_2     // Divide After PLL (now 40 MHz)
#pragma config FWDTEN = OFF         // Watchdog Timer Disabled
#pragma config ICESEL = ICS_PGx1    // ICE/ICD Comm Channel Select (pins 4,5)
#pragma config JTAGEN = OFF         // Disable JTAG
#pragma config FSOSCEN = OFF        // Disable Secondary Oscillator

// Defines

#define fftc fft16c64 // from fftc.h, for N = 64

#define SYSCLK (40000000L)

#define SAMPLES 64          

#define PWM_FREQ   48000            // Output PWM frequency
#define DUTY_CYCLE  1               

int tab[SAMPLES];
int i = 0;

int analogRead(char analogPIN)
{
    AD1CHS = analogPIN << 16;       // AD1CHS<16:19> controls which analog pin goes to the ADC

    AD1CON1bits.SAMP = 1;           // Sampling
    while(AD1CON1bits.SAMP);        // wait until acquisition is done
    while(!AD1CON1bits.DONE);       // wait until conversion done

    return ADC1BUF0;                
}

void adcConfigureManual()
{
    AD1CON1CLR = 0x8000;    // disable ADC before configuration

    AD1CON1 = 0x00E0;       // internal counter ends sampling and starts conversion (manual sample)
    AD1CON2 = 0;            // AD1CON2<15:13> set voltage reference to pins AVSS/AVDD

    // Found on the web (todo: check the datasheet)
    AD1CON3 = 0x0f01;       // TAD = 4*TPB, acquisition time = 15*TAD
} 

int main( void)
{
    SYSTEMConfigPerformance(SYSCLK);

    // Set OC1 to pin 2 with peripheral pin select
    RPA0Rbits.RPA0R = 0x0005;

    // Configure standard PWM mode for output compare module 1 
    OC1CON = 0x0006;

    for(i = 0; i<SAMPLES; i++)
        tab[i] = 0;

    // From datasheet: 
    // PR = [FPB / (PWM Frequency * TMR Prescale Value)] – 1
    PR2 = (SYSCLK / PWM_FREQ) - 1;

    // Initial duty cycle value
    OC1RS = (PR2 + 1) * ((float)DUTY_CYCLE / 100);

    T2CONSET = 0x8000;      // Enable Timer2, prescaler 1:1
    OC1CONSET = 0x8000;     // Enable Output Compare Module 1

        // Configure pins as analog inputs
        ANSELBbits.ANSB3 = 1;   // set RB3 (AN5) to analog
        TRISBbits.TRISB3 = 1;   // set RB3 as an input
        TRISBbits.TRISB5 = 0;   // set RB5 as an output (note RB5 is a digital only pin)

        adcConfigureManual();   // Configure ADC
        AD1CON1SET = 0x8000;    // Enable ADC

        int pos=0, dat = 0;

        int log2N = 6; // log2(64) = 6
        int N = 1 << log2N; // N = 2^6 = 64
        int din[N];
        int dout[N];
        int dout2[N];
        int scratch[N];

    while(1)
    {
            //foo = analogRead 5); // note that we call pin AN5 (RB3) by it's analog number
            dat = analogRead(5);

            //dat += din[pos]*3/10;
            //if(dat > 1023) dat -= 512;

            din[pos] =  dat;

            mips_fft16(dout, din, fftc, scratch, 1);
            mips_fft16(dout2, dout, fftc, scratch, 1);


            if(++pos >= SAMPLES) pos = 0;

            OC1RS = (PR2 + 1) * ( ((float)dout2[pos])/1023); // Write new duty cycle
    }

    return 0;
}

En primer lugar, me gustaría saber si lo que estoy haciendo es correcto. Quiero decir, no hay función ifft, ¿el segundo fft realiza un ifft? En segundo lugar, al hacer esto, puedo escuchar un sonido de alta frecuencia, que puede resultar de la latencia que retrasa el cambio del ciclo de trabajo de PWM OC1RS. Cómo puedo arreglar esto ? Vi que la gente suele usar Timer2 para cambiar el ciclo de trabajo, pero cada vez que lo intento, no obtengo nada en la salida (sin sonido). ¿Cómo puedo implementar esto en mi caso?

¿Tiene una copia de "Recetas numéricas"? Es un libro muy bueno sobre estos y otros temas. Incluyendo hacer filtrado con FFTs/IFFTs. Pero no es difícil ver cómo hacerlo con fuerza bruta. Asegúrese de comprender también las funciones de las ventanas y su propósito. "Fast Fourier Transform and Its Applications" de E Brigham en prácticamente CUALQUIER edición es un libro de primer nivel sobre el tema, en general. Y su pregunta parece bastante general, por lo que también es un buen libro para obtener y leer a fondo.
He oído hablar de él, gracias por esta buena sugerencia. Pero mi pregunta es más sobre "cómo hacerlo correctamente", quiero decir en términos de rendimiento. Sé usar FFT, las usé mucho con Matlab, pero nunca con un PIC o dsPIC. Por lo tanto, me gustaría usar la función DSP lib lo más relevante posible.
El rendimiento es al menos cuatro preguntas: ¿Qué pueden lograr las bibliotecas "enlatadas"? ¿Cómo puedo cambiar la declaración del problema para hacer un mejor uso de las bibliotecas "enlatadas" para el rendimiento? ¿Qué podría lograr si escribiera un código personalizado para mis necesidades? ¿Alguno o todos estos cumplen con mis requisitos? Dado que conoce su solicitud mejor que cualquiera de nosotros, y dado que las anteriores son realmente "preguntas de solicitud", no estoy seguro de cómo ayudaremos mucho aquí sin mucha más información sobre los detalles de la solicitud. (Sin embargo, el pase bajo es fácil. ¿Qué has intentado?)
Gracias a sus comentarios, pude realizar FFT. Pero todavía tengo algunas preguntas sobre los resultados que obtengo, y particularmente sobre el uso de Timer2, que eventualmente podría solucionar el problema de latencia.

Respuestas (1)

Encuentro que esto se aprende mejor usando la herramienta DCMI y la simulación, aunque mipsfir no simula por alguna razón.

mips_fft16/32 hace fft e ifft. Es todo lo mismo, de verdad. Transfiere la entrada de un dominio a otro.

Su implementación tiene fallas. La salida del mipsfft es un número complejo, me sorprende que aún esté compilado. Haz un estudio sobre los números complejos. No veo su definición de fftc, pero un tamaño de 1 no tiene sentido en la llamada mips_fft16. Es realmente bastante simple. Si desea crear una onda sinusoidal en un búfer de muestra de 32, coloque el valor (amplitud*32) en din[1].re. din[2].re para un seno doble, y así sucesivamente. Entonces, para una frecuencia de muestreo de 16000, 16000/32 = 500 * 1 = 500 Hz. Tenga en cuenta que su salida está en dout[x].re. Entonces, una llamada a fft/ifft se parecerá más a esto para un fft de tamaño 32. mips_fft16(dout, din, fftc16c32, cero, 6);

Además, recuerda que la salida de un fft será real e imaginaria, así que si solo quieres picos, no olvides que A al cuadrado + B al cuadrado = C al cuadrado.

Real e imaginario = fase de salida, muy simplificado.

OC1RS = (PR2 + 1) * ( ((float)dout2[pos])/1023); // Write new duty cycle

Esta no es la forma correcta de hacer esto en absoluto. La actualización de OC1RS puede estar en cualquier lugar, pero sugiero colocar un temporizador dentro de una frecuencia de muestreo elegida, digamos 16000. Luego, alimente OC1RS con un puntero de datos incrementales de algún tipo para recorrer su salida.