Sipeed Tang Nano 4k z ADC - Gowin FPGA Designer

ethanak · Marzec 7, 2022

Nie o to mi chodzi. Pytanie: czy Twój program (i przetwornik) będzie pobierał dane dokładnie co 1/30 msec (czy ile tam chcesz), czy może będzie to 30k razy na sekundę, ale odstępy między kolejnymi samplami nie będą identyczne? Bo jeśli nie będą, to pewnie trzeba będzie obciąć do 8 bitów

Przy okazji, 13 bitów wystarczy, chyba że chcesz się bawić w sprawdzanie, czy facet na kotle nie pomylił się o centymetr

FlyingDutch · Marzec 7, 2022

7 minut temu, ethanak napisał:

Nie o to mi chodzi. Pytanie: czy Twój program (i przetwornik) będzie pobierał dane dokładnie co 1/30 msec (czy ile tam chcesz), czy może będzie to 30k razy na sekundę, ale odstępy między kolejnymi samplami nie będą identyczne? Bo jeśli nie będą, to pewnie trzeba będzie obciąć do 8 bitów

Przy okazji, 13 bitów wystarczy, chyba że chcesz się bawić w sprawdzanie, czy facet na kotle nie pomylił się o centymetr

@ethanak,

jak już pisałem trzeba prześledzić jak wygląda w czasie sygnał DRDY (wyjście przetwornika), który mówi (zbocze opadające), że kolejna próbka jest gotowa. Zrobię to i będę mógł powiedzieć, jak jest regularność pobierania próbek. Zajrzę też do noty katalogowej układu ADS1256, może tam też uda się znaleźć odpowiedź na to pytanie.

Pozdrawiam

ethanak · Marzec 7, 2022

Ech... tylko po co się w to bawić? Czy na pewno ten przetwornik jest przeznaczony do sygnału audio? I czy na pewno chcesz zapalać samochód na pych twierdząc, że jak go już się nauczysz tak odpalać to się zainteresujesz, co to jest ten rozrusznik?

FlyingDutch · Marzec 7, 2022

47 minut temu, ethanak napisał:

Ech... tylko po co się w to bawić? Czy na pewno ten przetwornik jest przeznaczony do sygnału audio? I czy na pewno chcesz zapalać samochód na pych twierdząc, że jak go już się nauczysz tak odpalać to się zainteresujesz, co to jest ten rozrusznik?

Hej,

poszukałem w karcie katalogowej i tam piszą, że efektywna liczba bitów dla prędkości próbkowania 30Ksps wynosi 20, oraz, że nie występuje gubienie kodów i nieliniowość jest bardzo niska. Zobaczymy jak będzie wyglądała policzona FFT dla znanego przebiegu (w funkcji czasu). Przebieg sygnału DRDY na analizatorze też wygląda bardzo regularnie.

BTW: chcę zdobyć jakieś praktyczne doświadczenie z tym typem przetwornika. Spróbuję potem z tym przetwornikiem I2S i porównam rezultaty (na początku nawet zakładałem, że pierwsze próby zrobię z przetwornikiem z I2S, ale płytka z ADS1256 dotarła pierwsza).

Pozdrawiam

Edytowano Marzec 7, 2022 przez FlyingDutch

FlyingDutch · Marzec 7, 2022

Spróbowałem policzyć tym programem FFT o długości 2048 punktów (na razie z tablicy w kodzie):

#include "FFT.h" // include the library
#include "FFT_signal.h"


void setup() {
  Serial.begin(115200); // use the serial port
  char print_buf[300];
  fft_config_t *real_fft_plan = fft_init(FFT_N, FFT_REAL, FFT_FORWARD, fft_input, fft_output);

  for (int k = 0 ; k < FFT_N ; k++)
    real_fft_plan->input[k] = (float)fft_signal[k];

  long int t1 = micros();
  // Execute transformation
  fft_execute(real_fft_plan);
  
  // Print the output
  for (int k = 1 ; k < real_fft_plan->size / 2 ; k++)
  {
    /*The real part of a magnitude at a frequency is followed by the corresponding imaginary part in the output*/
    float mag = sqrt(pow(real_fft_plan->output[2*k],2) + pow(real_fft_plan->output[2*k+1],2))/1;
    float freq = k*1.0/TOTAL_TIME;
//    sprintf(print_buf,"%f Hz : %f", freq, mag);
//    Serial.println(print_buf);
    if(mag > max_magnitude)
    {
        max_magnitude = mag;
        fundamental_freq = freq;
    }
  }
  long int t2 = micros();
  
  Serial.println();
  /*Multiply the magnitude of the DC component with (1/FFT_N) to obtain the DC component*/
  sprintf(print_buf,"DC component : %f g\n", (real_fft_plan->output[0])/10000/FFT_N);  // DC is at [0]
  Serial.println(print_buf);

  /*Multiply the magnitude at all other frequencies with (2/FFT_N) to obtain the amplitude at that frequency*/
  sprintf(print_buf,"Fundamental Freq : %f Hz\t Mag: %f g\n", fundamental_freq, (max_magnitude/10000)*2/FFT_N);
  Serial.println(print_buf);

  Serial.print("Time taken: ");Serial.print((t2-t1)*1.0/1000);Serial.println(" milliseconds!");
  
  // Clean up at the end to free the memory allocated
  fft_destroy(real_fft_plan);
  
}

void loop() {
  
}

Dla takiego sygnału:

#define FFT_N 2048 // Must be a power of 2
#define TOTAL_TIME 9.391904 //The time in which data was captured. This is equal to FFT_N/sampling_freq


float fft_input[FFT_N];
float fft_output[FFT_N];

float max_magnitude = 0;
float fundamental_freq = 0;


/* Dummy data (Output of an accelerometer)
 * Frequency: 5 Hz
 * Amplitude: 0.25g
*/
double fft_signal[FFT_N] = {
11100,10600,11200,11700,12200,12900,12900,13100,11600,12900,13600,13000,12100,11700,11000,11500,10400,10800,9800,9800,9400,9700,8700,8700,8200,8200,7700,7800,7200,7200,6900,5900,5600,7400,8100,8300,8600,9400,9300,9600,10300,11000,10700,11200,11200,11200,11800,11700,12900,12700,12700,11900,11500,11100,11700,13700,12900,11700,10600,11000,10900,11100,9800,9600,9900,9400,9000,8300,8100,7600,7100,7600,7000,9700,8200,8200,7600,7000,7500,8500,8400,9400,9100,9400,10200,10400,10100,10600,11000,11700,11700,12700,12100,12200,12000,12200,14200,13700,12400,11100,11000,11500,10800,10700,10500,10400,9800,9400,9300,9300,9200,9000,8000,7600,7600,7400,7500,3600,6100,6700,9000,8500,7800,8000,8800,9100,10000,9600,10900,10600,11000,11400,11800,11700,11100,12100,12600,13400,11000,12800,12100,13100,12800,12100,11300,10900,11200,10500,10200,9500,9900,9800,9000,8600,8300,8000,7900,7300,6800,6700,9500,9300,7500,7800,7100,8000,8400,8700,8900,10000,10300,10300,9900,10300,11300,11800,12000,11200,12200,13100,13400,13100,13500,12900,13100,10900,11900,11500,11700,10400,10400,10700,9900,10200,9600,9300,8700,8500,8500,8300,8000,7700,6700,7200,4800,6300,8400,8000,7600,8100,9500,9600,9600,9600,10700,10400,11000,10800,11300,11700,11200,11900,12500,12800,12500,10900,12500,12600,11900,11800,11800,11200,11500,11000,10600,10500,9900,9800,9300,8700,9200,8600,7900,7600,7500,6600,7500,7400,8800,6200,6500,7500,8400,8700,8700,9400,9500,10000,10400,10900,11200,11500,11500,11700,11900,13000,12000,13700,12900,11600,13600,12100,11300,11500,11800,11600,10500,10900,10600,10600,9900,9000,9000,9000,9100,8000,7700,7700,8100,6700,6800,4300,7300,7400,8500,7500,8500,9100,9900,8900,9600,9700,10200,10700,10900,11500,11300,11600,11800,13300,13000,12900,11100,12500,12300,13100,12200,11800,11000,11200,11000,10700,10000,10400,9900,9400,8800,8800,8400,8300,8100,7300,6800,7400,8700,5700,5900,9500,8900,9000,8200,8800,8600,9400,10400,10500,10200,10500,12300,10700,12100,11900,13100,13300,12600,12200,11700,11800,11800,12700,12200,11800,11000,10900,10500,10000,9900,9600,10000,9200,9100,8500,8200,7800,7700,7300,6800,8600,9100,8000,8200,7400,7500,8300,8900,9100,9500,9900,10400,10600,11200,11000,11800,11700,11800,12700,12700,12800,12700,13300,13700,13700,10800,10900,11400,11600,11400,10900,10000,10300,10000,9400,9200,9100,8900,8700,8300,8000,7600,6700,7300,5400,6100,6200,8600,8500,8200,8300,8100,9900,9800,10000,9800,10300,10700,12000,12100,12100,11700,12600,12600,12800,11500,11900,11100,12000,12900,12300,11500,11200,11100,10500,10200,9700,9700,10300,9300,8900,8300,8200,7700,7500,7000,6900,9700,8600,6000,6000,7100,7400,8400,9200,9100,9700,9800,10600,10200,11300,11300,11600,12000,11400,12200,12700,12700,12800,12800,13700,13800,11800,11700,11700,11600,10900,10800,10400,10400,9400,9300,8900,9300,9000,8000,7900,8000,7500,6700,8800,4500,6700,7100,8300,7400,8600,9100,9400,9700,10000,9500,9800,11300,11900,12000,11100,11000,12100,11800,12600,12700,11300,12300,11800,13400,12700,11000,11000,11200,11100,11000,9900,9900,9800,9300,9200,8600,8600,8300,8000,7800,7100,7500,9000,9700,8600,6400,7200,7300,9100,9200,9100,9700,9400,10100,10700,11500,11800,11200,12100,11600,13100,12500,13000,12900,13100,13900,12300,11900,11200,11700,11300,10900,10300,10300,10800,10300,9500,8700,9700,9500,8300,8200,7600,6900,7500,7100,4400,6400,6100,8800,8200,8300,9000,9400,9700,10100,10500,10100,11400,11700,12000,11900,11000,11800,12300,13000,13100,11300,12100,11500,12900,12900,12800,10700,11200,11000,10800,9900,9400,9600,9800,9100,8400,8000,7500,7900,6900,7400,7100,9600,7600,6900,6400,7300,8500,8400,9000,9400,9400,9600,10300,10500,10500,11200,11400,11800,11800,12900,12500,13000,11900,13500,12800,12700,11000,12100,11400,11400,11000,11100,10300,10500,9300,9700,9500,9400,8100,8000,7700,8300,7500,6800,6800,4400,8000,7800,7900,8200,7900,9100,9200,9300,9300,10800,11200,10900,11100,10900,11900,12100,12100,12300,12500,12800,10900,13400,13500,11800,11400,11000,10800,10600,10000,9800,10500,9900,10000,9400,9000,8900,8700,8500,7800,7900,6900,6900,7800,7800,6500,6400,7500,7500,8000,8000,8300,9600,9800,10600,10200,10300,11600,11300,12100,10800,12200,12400,12400,12400,12700,11900,11800,11900,12400,12100,11400,11200,10700,10500,9700,9600,9700,9600,9000,8800,8100,8000,8500,7200,7300,7000,4000,7200,8900,8700,7800,8100,9000,8900,9300,9500,10300,10600,10800,10700,10900,11200,12100,12500,12500,12400,12700,10600,13200,13500,13300,12100,12200,11200,11500,11000,10400,10000,9800,9700,9400,9100,9000,8600,8500,7900,7400,6800,7300,9000,7100,7400,6800,7400,7600,9100,8700,9500,9700,10400,10700,10100,11000,10800,11800,11600,11300,12200,12300,12800,11300,12000,11100,11300,13300,12800,11200,10700,10900,10900,10400,10100,9300,9700,9500,9200,8300,8200,7600,8000,7000,7000,8100,8200,8300,7800,6800,7200,8500,9000,8900,9300,9500,10600,10900,11300,10500,10800,11400,10900,13100,12700,13100,11900,12500,14000,13600,11300,11000,12000,11600,11100,11800,9800,9500,9400,9800,9300,9000,8800,8200,7600,8100,7400,7400,7300,3700,7900,7800,9200,8300,8400,8900,8800,9700,9700,9600,10400,10800,11300,10900,11500,11500,12600,12400,12400,13200,11100,13200,12700,12100,12000,11800,11500,11700,11200,10700,9800,10100,9600,9700,9500,8600,8300,8400,7900,6800,7600,6600,8600,8600,6700,6800,7200,7800,8400,8900,10100,10400,10500,10700,11100,10300,10600,11800,12100,11600,12400,12500,12600,12100,13400,12700,12800,11000,12000,12100,11600,10400,10600,11000,10400,9500,9800,9100,9000,8800,7800,8000,8000,7500,7000,6900,5900,6200,8700,7600,7400,8400,8400,8800,9200,9600,9700,10800,11700,12100,11500,11500,11500,12100,11700,13100,12800,11300,12600,12800,12900,12200,11500,11200,11700,10900,10500,9800,9600,10100,9000,8800,8600,8300,8300,7600,7100,6800,7000,9800,7800,7800,6300,7100,7400,8300,8100,9200,9500,10200,10000,10500,11300,10900,12200,12300,11400,12800,12500,12900,12400,11400,11000,11500,10700,12200,12600,11400,10500,10900,10900,10300,10200,9600,9400,8900,9000,8100,8200,8000,7600,6900,7400,4900,8300,9000,8100,7700,7700,8100,9700,10300,9600,10000,9700,10400,11200,11100,11200,11200,11400,12400,12000,13100,11300,12300,12800,12700,13100,12200,11700,10900,10800,10500,10300,10100,9400,10100,9400,9000,8400,8200,7600,7600,7600,7400,9600,8800,8200,8300,7000,7800,8300,8600,8100,8500,9000,10500,10600,10800,10700,11600,11900,12100,12700,12900,12800,12300,13300,13900,13900,11000,11100,11100,11900,10900,11000,10700,10000,9800,9300,9100,9400,8900,8000,8100,8000,7500,7200,6700,5900,8000,5500,8700,8100,8600,7700,8200,9000,9700,10500,10600,10700,11000,11100,11600,11600,11700,13100,12700,12500,11100,12000,11700,12800,12800,12100,11200,10900,10900,10400,10500,9900,9200,9700,9300,9400,7900,7800,7700,7500,7300,6800,9600,7900,8400,8900,7400,8100,7700,8800,8800,9100,9500,9800,11200,11200,11300,11000,11000,11600,11500,12700,12700,12300,12400,12500,13800,11500,11600,11200,12100,11600,10800,9700,9900,9800,9500,9200,9000,8100,8600,8200,8300,7500,7500,6900,7300,5800,5700,6900,7900,8800,8400,8300,9600,9600,10900,10200,10500,10400,12000,11800,12300,12000,12200,12400,12600,11300,11900,11200,11400,13300,12600,11600,10900,11000,10900,10200,9900,9700,9100,9100,8700,8500,8100,7900,7000,7000,7100,9200,9000,8700,7000,7400,7000,8900,9000,9400,9600,9900,10200,10700,11300,11000,11600,11500,12300,13100,12700,12400,12600,13400,14000,12700,11700,11200,12200,12100,10600,10900,9800,10700,9400,9300,9100,9400,8900,8800,8200,8000,7600,7400,6900,4400,6500,7100,9100,8200,8000,8700,9100,9700,10400,9400,10200,10900,10500,12000,12100,11700,11200,12700,12400,12500,11000,11800,11800,11900,13000,11800,11400,11000,10900,11100,9900,10400,9500,10000,9100,8800,8200,8300,7500,7200,6800,7000,9700,9600,7700,7900,6500,7800,8500,8800,9200,9400,9400,10300,10300,11100,11400,11500,11900,11800,11700,12700,12800,12800,12700,13400,13500,11200,11700,11800,11300,10700,10400,10500,10700,10000,9200,8900,9100,9300,8600,7500,7900,7300,7300,7200,4600,7400,6900,8200,7300,8300,8400,9300,9700,9800,10200,10600,10900,12100,11600,11400,11600,11600,12100,12400,12800,10600,11400,12000,13500,12100,11900,10600,10600,11000,10500,10400,9600,9800,9800,9000,9200,8300,8300,7500,7200,6900,7300,10000,8400,9100,6400,6900,7300,8000,9000,9300,9500,10100,10500,10200,10700,11300,10900,11900,12000,12300,12700,12500,12900,12400,13400,12200,11100,11500,11600,11500,10600,11000,10700,10100,10100,9000,8400,9000,9300,8500,8300,8000,7300,7200,7200,5000,9000,7600,8600,7100,8000,8400,8800,9700,9400,9700,10800,10600,11300,11100,11400,11700,11100,12000,12600,13200,10800,12000,11400,13500,13700,12000,10600,10900,10800,10400,10100,9500,9800,9400,8700,8400,8800,8300,7800,7400,7000,7100,8600,6800,5800,6900,7600,8200,8500,8500,8800,10000,10300,10600,10100,10700,11300,12400,11600,12800,12400,12300,12100,11900,12700,11200,11800,12600,11500,12000,11100,10600,10800,10500,9600,9500,9500,9100,9000,8600,8400,8100,8300,7500,7000,7500,5400,8300,8600,7800,7800,7700,8700,9600,9800,8800,9900,10700,10800,11700,11000,11200,12100,11700,12200,12200,13400,11000,14000,13200,12200,10800,10600,11200,11600,10900,10000,10400,10000,9600,10100,8500,9000,8800,8900,8300,7600,7200,7300,4100,6000,6100,6700,8100,8400,8500,9200,9500,10000,9800,9600,10600,11100,11800,11200,11600,11000,12200,13000,13000,11500,12500,11400,12400,12400,12900,11300,11500,10500,10900,10600,9900,10500,9400,9100,8800,8800,8300,8200,7600,7700,7100,7600,8500,8000,9200,7100,7300,7900,8500,9100,9500,9000,10100,9800,11400,11100,11500,11300,11600,12200,12400,12300,12100,12200,13700,13700,11000,11500,11900,11600,10800,10800,10000,9800,10400,9300,9400,9000,8800,8500,8000,8300,7400,7600,7600,4700,6200,6100,8600,7900,8500,8700,9300,9700,9600,9800,10100,11900,11500,11400,11900,11300,12400,12300,12400,12000,10900,13400,11400,13200,12600,12500,11000,11300,10900,10500,10300,10100,9500,9700,9000,9000,8100,8100,7600,7500,7200,6700,9400,8800,8500,6800,7600,7700,7900,8800,9300,8300,10400,10100,11000,11000,11100,11600,11900,12500,11800,12700,13300,13000,13000,13600,13300,12100,10900,11100,12000,10800,11200,10600,10100,10000,9700,9500,9300,8900,8100,8000,8400,7200,7500,7200,5200,5700,5700,8000,7800,9100,8800,8700,9000,9300,10000,9900,11100,10900,11700,11600,11100,11600,12100,12500,13200,11100,12100,11300,11300,14000,12800,11000,10600,11600,10100,10500,9900,9500,9500,9200,9200,8400,7800,7500,7900,7700,7000,9100,9200,9100,7800,7300,7200,7900,9100,9100,9100,9600,9700,10300,11000,11200,10900,11000,11800,12300,12800,13100,13100,11700,13300,13800,12100,11100,10900,11200,11000,10800,10100,10100,9900,8700,9500,8800,8600,8400,7700,7500,7100,7800,7000,6500,5300,5300,5900,7700,9000,8800,9300,9800,9600,10000,10000,9900,10900,10900,11700,12000,11700,12100,12600,12800,11300,12100,11600,11400,13200,12800,10900,10900,10500,10500,10100,9800,9700,9600,9300,8100,8400,7900,8100,7400,7400,7000,8000,7700,8600,8400,7500,7700,8200,9200,9300,9200,9500,10000,11000,11500,11700,10900,11800,12200,12100,12500,12200,12800,12300,14200,14100,11600,11100,10900,11400,11400,10800,9800,10100,9800,9300,8900,9100,9200,8500,8200,8200,8100,7500,7100,7700,5400,5900,6600,6900,8700,8700,8000,9700,9600,10300,9800,10000,10700,11600,12000,11700,12000,12500,12500,12600,12000,12700,13500,11700,10800,11100,12400,11200,10800,10400,10100,10200,10100,9400,9300,9300,9300,8200,8100,7800,7400,7500,7500,3800,7100,6900,8600,8400,8700,8900,9200,9500,10100,10400,10400,10500,10600,12000,11700,11700,11000
};

Wynik jest jak na ekranie poniżej:

Liczenie FFT o długości 2048 punktów zajęło 6.74 milisekundy, a wyniki są zgodne z podanym w tablicy sygnałem. Nie mam żadnego programowanego generatora sygnałowego, żeby podać do przetwornika ADC jakiś ciekawszy sygnał - muszę się zastanowić jak go wygenerować. Mam jakiegoś bardzo prymitywnego DAC'a (12 bit) - muszę go podłączyć do Arduino UNO i wygenerować jakiś ciekawszy sygnał audio.

Program do liczenia FFT z tej strony WWW:

https://medium.com/swlh/how-to-perform-fft-onboard-esp32-and-get-both-frequency-and-amplitude-45ec5712d7da#9ebf

Pozdrawiam

Edytowano Marzec 7, 2022 przez FlyingDutch

ethanak · Marzec 7, 2022

Zrób coś takiego:

Znajdź w necie jakiegoś sampla trąbki (albo innwj blachy)

Zrób analizę

Zrób syntezę

Posłuchaj co wyszło

FlyingDutch · Marzec 7, 2022

Przed chwilą, ethanak napisał:

Zrób coś takiego:

Znajdź w necie jakiegoś sampla trąbki (albo innwj blachy)

Zrób analizę

Zrób syntezę

Posłuchaj co wyszło

@ethanak,

wygenerować tablicę Float z przebiegiem sygnału mogę z łatwością, ale chodzi mi o to aby liczyć FFT z sampli z przetwornika ADS1256 - więc muszę mieć sygnał analogowy.

Pozdrawiam

ethanak · Marzec 7, 2022

No to podłącz wyjście jakiejś porządnej karty dźwiękowej do przetwornika, puść tam jakiś sygnał 24 bit na 64 kHz, lepiej raczej w domowych warunkach nie znajdziesz.

FlyingDutch · Marzec 7, 2022

9 minut temu, ethanak napisał:

No to podłącz wyjście jakiejś porządnej karty dźwiękowej do przetwornika, puść tam jakiś sygnał 24 bit na 64 kHz, lepiej raczej w domowych warunkach nie znajdziesz.

Ten 12-to bitowy DAC mi wystarczy, ale to już jutro - na dzisiaj kończę

Pozdrawiam

FlyingDutch · Marzec 8, 2022

Cześć,

dzisiaj udało mi się wygenerować sygnał audio na "Analog Test Shield" (układ MCP4725 12-to bitowy DAC) podłączony do "Arduino UNO". Tutaj zdjęcie całego układu (ESP32, płytka z ADC ADS1256, Arduino UNO+"Analog Test Shield"):

Tutaj link do opisu "Analog Test Shield":

https://www.waveshare.com/wiki/Analog_Test_Shield

Za pomocą programu na Arduino UNO dla tego shield'a wygenerowałem sygnał sinusoidalny o f=100Hz zanieczyszczony sygnałem o mniejszej amplitudzie i wyższej częstotliwości - patrz poniżej (przebieg sygnału na oscyloskopie):

Te nierównomierności, które widać na przebiegu to właśnie wpływ przebiegu o większej częstotliwości nałożonego na sygnał 100Hz. Tutaj kod programu dla Arduino UNO do generacji tego przebiegu:

/*
  Analog input, PWM output for LED, DAC output for Speaker
  created 25 May. 2015
  by WaveShare
 */
#include <Wire.h>
#define MCP4725_WRITE  0x60 //7-bit slave address 

const int analogIn = A0;  // Analog input
const int analogOut = 6;  // Analog output (PWM)

int i,j,count;
int Value = 0;        // value read from analog input
const unsigned short Sine12bit[32] = {	 //data for sine wave - zanieczyszczona
                      2047, 2447-400, 2831, 3185-400, 3498, 3750-400, 3939, 4056-400, 4095, 4056-400,
                      3939, 3750-400, 3495, 3185-400, 2831, 2447-400, 2047, 1647-400, 1263, 909-400, 
                      599-400, 344, 155, 38, 0, 38, 155, 344, 599, 909, 1263, 1647};

void setup() {
  // initialize serial communications at 9600 bps:
  //Serial.begin(9600);
  // initialize I2C bus:
  Wire.begin();
}

void loop() {
  
  //count++;               //counting
  for(i = 0;i<32;i++)
  {
    MCPWriteDAC(Sine12bit[i]);     //DAC output for Speaker
  }
  
//  if(count%10 ==0)             
//  {
    //analog output(PWM) for LED,maps the result to a range from 0 to 255
//    analogWrite(analogOut, map(Sine12bit[j++], 0, 4096, 0, 255));
//    if(j == 32)j = 0;
//    if(count%100 == 0)
//    {
//      count = 0;
//      //read analog input
//      Value = analogRead(analogIn); 
//      //maps the result to a range from 0 to 3300mv
//      Value = map(Value, 0, 1024, 0, 5000);
//      //print analog input
//      Serial.print("The analog input is: " );
//      Serial.print(Value);
//      Serial.println(" mv" );
//    }
//  }
}

void MCPWriteDAC(unsigned short value)    //sent 12-bit value to MCP4725
{
    char data[2];
    data[0] = value >> 8;
    data[1] = value;
    Wire.beginTransmission(MCP4725_WRITE); // transmit to device #4
    Wire.write(data[0]);
    Wire.write(data[1]);
    Wire.endTransmission();    // stop transmitting
}

Wyjście przetwornika DAC (sygnał analogowy) jest podłączone do pierwszego kanału płytki z przetwornikiem ADC (ADS1256), gdzie w czasie rzeczywistym są zbierane sample z częstotliwością próbkowania 30 KHz (30ksps). Zmodyfikowałem program do liczenia FFT na ESP32 tak, że liczona jest FFT dla 256 punktów i w wyniku otrzymujemy 256 częstotliwości z zakresu od 0Hz do 15KHz. FFT jest liczone z 256 sampli zebranych przetwornikiem ADS1256. Zebranie 256 sampli z prędkością 30Ksps zajmuje: 256/30000 s co daje około 8,5333 ms. z tego wynika, że w ciągu sekundy możemy zebrać próbki i policzyć FFT (8.5ms+0,17ms) 114 razy (0,17 ms to czas liczenia 256 punktowej FFT na ESP32). Patrz zrzut ekranu z wynikami obliczenia FFT:

A tutaj wyniki liczenia FFT tekstowo (z Serial Monitora Arduino IDE):

ets Jun  8 2016 00:22:57

rst:0x1 (POWERON_RESET),boot:0x13 (SPI_FAST_FLASH_BOOT)
ets Jun  8 2016 00:22:57

rst:0x10 (RTCWDT_RTC_RESET),boot:0x13 (SPI_FAST_FLASH_BOOT)
configsip: 0, SPIWP:0xee
clk_drv:0x00,q_drv:0x00,d_drv:0x00,cs0_drv:0x00,hd_drv:0x00,wp_drv:0x00
mode:DIO, clock div:1
load:0x3fff0018,len:4
load:0x3fff001c,len:1216
ho 0 tail 12 room 4
load:0x40078000,len:10944
load:0x40080400,len:6388
entry 0x400806b4
E (131) psram: PSRAM ID read error: 0xffffffff
Starting ADC
ADC Started
Channel set to Single end ch0
58.593754 Hz : 90.054092
117.187508 Hz : 136.518585
175.781250 Hz : 31.235022
234.375015 Hz : 20.636230
292.968750 Hz : 14.748680
351.562500 Hz : 10.471104
410.156281 Hz : 10.598127
468.750031 Hz : 9.879042
527.343750 Hz : 6.423238
585.937500 Hz : 6.201964
644.531250 Hz : 6.285541
703.125000 Hz : 3.427619
761.718750 Hz : 3.034126
820.312561 Hz : 3.983618
878.906311 Hz : 2.359965
937.500061 Hz : 2.312025
996.093811 Hz : 3.658482
1054.687500 Hz : 2.859877
1113.281250 Hz : 2.632012
1171.875000 Hz : 3.362010
1230.468750 Hz : 2.717731
1289.062500 Hz : 2.295434
1347.656250 Hz : 2.136851
1406.250000 Hz : 1.470130
1464.843750 Hz : 3.274125
1523.437500 Hz : 0.713512
1582.031372 Hz : 3.794071
1640.625122 Hz : 1.176832
1699.218872 Hz : 2.227799
1757.812622 Hz : 2.330910
1816.406372 Hz : 1.945227
1875.000122 Hz : 2.092054
1933.593872 Hz : 1.949055
1992.187622 Hz : 1.802135
2050.781250 Hz : 1.859932
2109.375000 Hz : 1.611560
2167.968750 Hz : 1.708672
2226.562500 Hz : 1.832446
2285.156250 Hz : 1.604206
2343.750000 Hz : 1.681069
2402.343750 Hz : 1.754111
2460.937500 Hz : 1.527530
2519.531250 Hz : 1.500246
2578.125000 Hz : 1.516142
2636.718750 Hz : 1.375365
2695.312500 Hz : 1.355924
2753.906250 Hz : 1.351367
2812.500000 Hz : 1.323049
2871.093750 Hz : 1.320552
2929.687500 Hz : 1.313711
2988.281250 Hz : 1.319631
3046.875000 Hz : 1.256514
3105.468750 Hz : 1.306583
3164.062744 Hz : 1.837343
3222.656494 Hz : 0.966339
3281.250244 Hz : 1.098205
3339.843994 Hz : 1.068151
3398.437744 Hz : 1.068684
3457.031494 Hz : 1.113237
3515.625244 Hz : 1.046626
3574.218994 Hz : 0.982782
3632.812744 Hz : 1.029194
3691.406494 Hz : 0.985107
3750.000244 Hz : 0.932194
3808.593994 Hz : 0.981902
3867.187744 Hz : 0.958111
3925.781494 Hz : 0.931370
3984.375244 Hz : 0.953631
4042.968994 Hz : 0.907713
4101.562500 Hz : 0.908188
4160.156250 Hz : 0.934787
4218.750000 Hz : 0.860148
4277.343750 Hz : 0.928298
4335.937500 Hz : 0.988026
4394.531250 Hz : 0.887219
4453.125000 Hz : 0.963869
4511.718750 Hz : 1.044403
4570.312500 Hz : 0.838399
4628.906250 Hz : 0.997085
4687.500000 Hz : 0.891364
4746.093750 Hz : 0.822948
4804.687500 Hz : 0.790764
4863.281250 Hz : 0.838347
4921.875000 Hz : 0.766782
4980.468750 Hz : 0.774266
5039.062500 Hz : 0.836908
5097.656250 Hz : 0.795133
5156.250000 Hz : 0.799860
5214.843750 Hz : 0.832920
5273.437500 Hz : 0.797654
5332.031250 Hz : 0.794718
5390.625000 Hz : 0.804790
5449.218750 Hz : 0.780731
5507.812500 Hz : 0.786722
5566.406250 Hz : 0.788067
5625.000000 Hz : 0.755925
5683.593750 Hz : 0.818313
5742.187500 Hz : 0.788407
5800.781250 Hz : 0.764223
5859.375000 Hz : 0.772662
5917.968750 Hz : 0.770496
5976.562500 Hz : 0.775823
6035.156250 Hz : 0.843242
6093.750000 Hz : 0.640025
6152.343750 Hz : 0.710130
6210.937500 Hz : 0.719463
6269.531738 Hz : 0.708062
6328.125488 Hz : 0.724338
6386.719238 Hz : 0.713423
6445.312988 Hz : 0.724907
6503.906738 Hz : 0.724680
6562.500488 Hz : 0.715388
6621.094238 Hz : 0.724927
6679.687988 Hz : 0.714322
6738.281738 Hz : 0.705258
6796.875488 Hz : 0.715228
6855.469238 Hz : 0.695962
6914.062988 Hz : 0.695944
6972.656738 Hz : 0.714210
7031.250488 Hz : 0.693164
7089.844238 Hz : 0.707362
7148.437988 Hz : 0.729178
7207.031738 Hz : 0.670764
7265.625488 Hz : 0.752986
7324.219238 Hz : 0.646362
7382.812988 Hz : 0.715567
7441.406738 Hz : 0.675138

DC component : 0.000122 g

Fundamental Freq : 117.187508 Hz	 Mag: 0.000107 g

Time taken: 0.17 milliseconds!

Jak widać na wynikach częstotliwość prążka widma o największej amplitudzie to 117,18 Hz - właśnie ten efekt chciałem zobaczyć na wynikach FFT (wynika on z nałożenia sygnału o mniejszej amplitudzie na sinusoidę 100Hz).

Sam program do liczenia FFT (i zbierania próbek z ADS1256) składa się z trzech plików:

1) FFT.h - zawiera kod algorytmu FFT

2) FFT_signal.h - zawiera tablicę fft_signal[FFT_N] (double ) do której są zbierane sample i stałą FFT_N z liczbą sampli oraz TOTAL_TIME który mówi w jakim czasie te sample zostały zebrane

3) FFT_on_ESP32_Arduino.ino - główny program do liczenia FFT

Tutaj kod FFT.h (algorytm FFT):

/*

  ESP32 FFT
  =========

  This provides a vanilla radix-2 FFT implementation and a test example.

  Author
  ------

  This code was written by [Robin Scheibler](http://www.robinscheibler.org) during rainy days in October 2017.

  License
  -------

  Copyright (c) 2017 Robin Scheibler

  Permission is hereby granted, free of charge, to any person obtaining a copy
  of this software and associated documentation files (the "Software"), to deal
  in the Software without restriction, including without limitation the rights
  to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
  copies of the Software, and to permit persons to whom the Software is
  furnished to do so, subject to the following conditions:

  The above copyright notice and this permission notice shall be included in all
  copies or substantial portions of the Software.

  THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
  IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
  FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
  AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
  LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
  OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
  SOFTWARE.

*/
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <complex.h>

typedef enum
{
  FFT_REAL,
  FFT_COMPLEX
} fft_type_t;

typedef enum
{
  FFT_FORWARD,
  FFT_BACKWARD
} fft_direction_t;

#define FFT_OWN_INPUT_MEM 1
#define FFT_OWN_OUTPUT_MEM 2

typedef struct
{
  int size;  // FFT size
  float *input;  // pointer to input buffer
  float *output; // pointer to output buffer
  float *twiddle_factors;  // pointer to buffer holding twiddle factors
  fft_type_t type;   // real or complex
  fft_direction_t direction; // forward or backward
  unsigned int flags; // FFT flags
} fft_config_t;


fft_config_t *fft_init(int size, fft_type_t type, fft_direction_t direction, float *input, float *output);
void fft_destroy(fft_config_t *config);
void fft_execute(fft_config_t *config);
void fft(float *input, float *output, float *twiddle_factors, int n);
void ifft(float *input, float *output, float *twiddle_factors, int n);
void rfft(float *x, float *y, float *twiddle_factors, int n);
void irfft(float *x, float *y, float *twiddle_factors, int n);
void fft_primitive(float *x, float *y, int n, int stride, float *twiddle_factors, int tw_stride);
void split_radix_fft(float *x, float *y, int n, int stride, float *twiddle_factors, int tw_stride);
void ifft_primitive(float *input, float *output, int n, int stride, float *twiddle_factors, int tw_stride);
void fft8(float *input, int stride_in, float *output, int stride_out);
void fft4(float *input, int stride_in, float *output, int stride_out);


#define TWO_PI 6.28318530
#define USE_SPLIT_RADIX 1
#define LARGE_BASE_CASE 1



fft_config_t *fft_init(int size, fft_type_t type, fft_direction_t direction, float *input, float *output)
{
  /*
   * Prepare an FFT of correct size and types.
   *
   * If no input or output buffers are provided, they will be allocated.
   */
  int k,m;

  fft_config_t *config = (fft_config_t *)malloc(sizeof(fft_config_t));

  // Check if the size is a power of two
  if ((size & (size-1)) != 0)  // tests if size is a power of two
    return NULL;

  // start configuration
  config->flags = 0;
  config->type = type;
  config->direction = direction;
  config->size = size;

  // Allocate and precompute twiddle factors
  config->twiddle_factors = (float *)malloc(2 * config->size * sizeof(float));

  float two_pi_by_n = TWO_PI / config->size;

  for (k = 0, m = 0 ; k < config->size ; k++, m+=2)
  {
    config->twiddle_factors[m] = cosf(two_pi_by_n * k);    // real
    config->twiddle_factors[m+1] = sinf(two_pi_by_n * k);  // imag
  }

  // Allocate input buffer
  if (input != NULL)
    config->input = input;
  else 
  {
    if (config->type == FFT_REAL)
      config->input = (float *)malloc(config->size * sizeof(float));
    else if (config->type == FFT_COMPLEX)
      config->input = (float *)malloc(2 * config->size * sizeof(float));

    config->flags |= FFT_OWN_INPUT_MEM;
  }

  if (config->input == NULL)
    return NULL;

  // Allocate output buffer
  if (output != NULL)
    config->output = output;
  else
  {
    if (config->type == FFT_REAL)
      config->output = (float *)malloc(config->size * sizeof(float));
    else if (config->type == FFT_COMPLEX)
      config->output = (float *)malloc(2 * config->size * sizeof(float));

    config->flags |= FFT_OWN_OUTPUT_MEM;
  }

  if (config->output == NULL)
    return NULL;

  return config;
}

void fft_destroy(fft_config_t *config)
{
  if (config->flags & FFT_OWN_INPUT_MEM)
    free(config->input);

  if (config->flags & FFT_OWN_OUTPUT_MEM)
    free(config->output);

  free(config->twiddle_factors);
  free(config);
}

void fft_execute(fft_config_t *config)
{
  if (config->type == FFT_REAL && config->direction == FFT_FORWARD)
    rfft(config->input, config->output, config->twiddle_factors, config->size);
  else if (config->type == FFT_REAL && config->direction == FFT_BACKWARD)
    irfft(config->input, config->output, config->twiddle_factors, config->size);
  else if (config->type == FFT_COMPLEX && config->direction == FFT_FORWARD)
    fft(config->input, config->output, config->twiddle_factors, config->size);
  else if (config->type == FFT_COMPLEX && config->direction == FFT_BACKWARD)
    ifft(config->input, config->output, config->twiddle_factors, config->size);
}

void fft(float *input, float *output, float *twiddle_factors, int n)
{
  /*
   * Forward fast Fourier transform
   * DIT, radix-2, out-of-place implementation
   *
   * Parameters
   * ----------
   *  input (float *)
   *    The input array containing the complex samples with
   *    real/imaginary parts interleaved [Re(x0), Im(x0), ..., Re(x_n-1), Im(x_n-1)]
   *  output (float *)
   *    The output array containing the complex samples with
   *    real/imaginary parts interleaved [Re(x0), Im(x0), ..., Re(x_n-1), Im(x_n-1)]
   *  n (int)
   *    The FFT size, should be a power of 2
   */

#if USE_SPLIT_RADIX
  split_radix_fft(input, output, n, 2, twiddle_factors, 2);
#else
  fft_primitive(input, output, n, 2, twiddle_factors, 2);
#endif
}

void ifft(float *input, float *output, float *twiddle_factors, int n)
{
  /*
   * Inverse fast Fourier transform
   * DIT, radix-2, out-of-place implementation
   *
   * Parameters
   * ----------
   *  input (float *)
   *    The input array containing the complex samples with
   *    real/imaginary parts interleaved [Re(x0), Im(x0), ..., Re(x_n-1), Im(x_n-1)]
   *  output (float *)
   *    The output array containing the complex samples with
   *    real/imaginary parts interleaved [Re(x0), Im(x0), ..., Re(x_n-1), Im(x_n-1)]
   *  n (int)
   *    The FFT size, should be a power of 2
   */
  ifft_primitive(input, output, n, 2, twiddle_factors, 2);
}

void rfft(float *x, float *y, float *twiddle_factors, int n)
{

  // This code uses the two-for-the-price-of-one strategy
#if USE_SPLIT_RADIX
  split_radix_fft(x, y, n / 2, 2, twiddle_factors, 4);
#else
  fft_primitive(x, y, n / 2, 2, twiddle_factors, 4);
#endif

  // Now apply post processing to recover positive
  // frequencies of the real FFT
  float t = y[0];
  y[0] = t + y[1];  // DC coefficient
  y[1] = t - y[1];  // Center coefficient

  // Apply post processing to quarter element
  // this boils down to taking complex conjugate
  y[n/2+1] = -y[n/2+1];

  // Now process all the other frequencies
  int k;
  for (k = 2 ; k < n / 2 ; k += 2)
  {
    float xer, xei, x0r, xoi, c, s, tr, ti;

    c = twiddle_factors[k];
    s = twiddle_factors[k+1];
    
    // even half coefficient
    xer = 0.5 * (y[k] + y[n-k]);
    xei = 0.5 * (y[k+1] - y[n-k+1]);

    // odd half coefficient
    x0r = 0.5 * (y[k+1] + y[n-k+1]);
    xoi = - 0.5 * (y[k] - y[n-k]);

    tr =  c * x0r + s * xoi;
    ti = -s * x0r + c * xoi;

    y[k]   = xer + tr;
    y[k+1] = xei + ti;

    y[n-k]   =   xer - tr;
    y[n-k+1] = -(xei - ti);
  }
}

void irfft(float *x, float *y, float *twiddle_factors, int n)
{
  /*
   * Destroys content of input vector
   */
  int k;

  // Here we need to apply a pre-processing first
  float t = x[0];
  x[0] = 0.5 * (t + x[1]);
  x[1] = 0.5 * (t - x[1]);

  x[n/2+1] = -x[n/2+1];

  for (k = 2 ; k < n / 2 ; k += 2)
  {
    float xer, xei, x0r, xoi, c, s, tr, ti;

    c = twiddle_factors[k];
    s = twiddle_factors[k+1];

    xer = 0.5 * (x[k] + x[n-k]);
    tr  = 0.5 * (x[k] - x[n-k]);

    xei = 0.5 * (x[k+1] - x[n-k+1]);
    ti  = 0.5 * (x[k+1] + x[n-k+1]);

    x0r = c * tr - s * ti;
    xoi = s * tr + c * ti;

    x[k]   = xer - xoi;
    x[k+1] = x0r + xei;

    x[n-k]   = xer + xoi;
    x[n-k+1] = x0r - xei;
  }

  ifft_primitive(x, y, n / 2, 2, twiddle_factors, 4);
}

void fft_primitive(float *x, float *y, int n, int stride, float *twiddle_factors, int tw_stride)
{
  /*
   * This code will compute the FFT of the input vector x
   *
   * The input data is assumed to be real/imag interleaved
   *
   * The size n should be a power of two
   *
   * y is an output buffer of size 2n to accomodate for complex numbers
   *
   * Forward fast Fourier transform
   * DIT, radix-2, out-of-place implementation
   *
   * For a complex FFT, call first stage as:
   * fft(x, y, n, 2, 2);
   *
   * Parameters
   * ----------
   *  x (float *)
   *    The input array containing the complex samples with
   *    real/imaginary parts interleaved [Re(x0), Im(x0), ..., Re(x_n-1), Im(x_n-1)]
   *  y (float *)
   *    The output array containing the complex samples with
   *    real/imaginary parts interleaved [Re(x0), Im(x0), ..., Re(x_n-1), Im(x_n-1)]
   *  n (int)
   *    The FFT size, should be a power of 2
   *  stride (int)
   *    The number of elements to skip between two successive samples
   *  tw_stride (int)
   *    The number of elements to skip between two successive twiddle factors
   */
  int k;
  float t;

#if LARGE_BASE_CASE
  // End condition, stop at n=8 to avoid one trivial recursion
  if (n == 8)
  {
    fft8(x, stride, y, 2);
    return;
  }
#else
  // End condition, stop at n=2 to avoid one trivial recursion
  if (n == 2)
  {
    y[0] = x[0] + x[stride];
    y[1] = x[1] + x[stride + 1];
    y[2] = x[0] - x[stride];
    y[3] = x[1] - x[stride + 1];
    return;
  }
#endif

  // Recursion -- Decimation In Time algorithm
  fft_primitive(x, y, n / 2, 2 * stride, twiddle_factors, 2 * tw_stride);             // even half
  fft_primitive(x + stride, y+n, n / 2, 2 * stride, twiddle_factors, 2 * tw_stride);  // odd half

  // Stitch back together

  // We can a few multiplications in the first step
  t = y[0];
  y[0] = t + y[n];
  y[n] = t - y[n];

  t = y[1];
  y[1] = t + y[n+1];
  y[n+1] = t - y[n+1];

  for (k = 1 ; k < n / 2 ; k++)
  {
    float x1r, x1i, x2r, x2i, c, s;
    c = twiddle_factors[k * tw_stride];
    s = twiddle_factors[k * tw_stride + 1];

    x1r = y[2 * k];
    x1i = y[2 * k + 1];
    x2r =  c * y[n + 2 * k] + s * y[n + 2 * k + 1];
    x2i = -s * y[n + 2 * k] + c * y[n + 2 * k + 1];

    y[2 * k] = x1r + x2r;
    y[2 * k + 1] = x1i + x2i;

    y[n + 2 * k] = x1r - x2r;
    y[n + 2 * k + 1] = x1i - x2i;
  }

}

void split_radix_fft(float *x, float *y, int n, int stride, float *twiddle_factors, int tw_stride)
{
  /*
   * This code will compute the FFT of the input vector x
   *
   * The input data is assumed to be real/imag interleaved
   *
   * The size n should be a power of two
   *
   * y is an output buffer of size 2n to accomodate for complex numbers
   *
   * Forward fast Fourier transform
   * Split-Radix
   * DIT, radix-2, out-of-place implementation
   *
   * For a complex FFT, call first stage as:
   * fft(x, y, n, 2, 2);
   *
   * Parameters
   * ----------
   *  x (float *)
   *    The input array containing the complex samples with
   *    real/imaginary parts interleaved [Re(x0), Im(x0), ..., Re(x_n-1), Im(x_n-1)]
   *  y (float *)
   *    The output array containing the complex samples with
   *    real/imaginary parts interleaved [Re(x0), Im(x0), ..., Re(x_n-1), Im(x_n-1)]
   *  n (int)
   *    The FFT size, should be a power of 2
   *  stride (int)
   *    The number of elements to skip between two successive samples
   *  twiddle_factors (float *)
   *    The array of twiddle factors
   *  tw_stride (int)
   *    The number of elements to skip between two successive twiddle factors
   */
  int k;

#if LARGE_BASE_CASE
  // End condition, stop at n=2 to avoid one trivial recursion
  if (n == 8)
  {
    fft8(x, stride, y, 2);
    return;
  }
  else if (n == 4)
  {
    fft4(x, stride, y, 2);
    return;
  }
#else
  // End condition, stop at n=2 to avoid one trivial recursion
  if (n == 2)
  {
    y[0] = x[0] + x[stride];
    y[1] = x[1] + x[stride + 1];
    y[2] = x[0] - x[stride];
    y[3] = x[1] - x[stride + 1];
    return;
  }
  else if (n == 1)
  {
    y[0] = x[0];
    y[1] = x[1];
    return;
  }
#endif

  // Recursion -- Decimation In Time algorithm
  split_radix_fft(x, y, n / 2, 2 * stride, twiddle_factors, 2 * tw_stride);
  split_radix_fft(x + stride, y + n, n / 4, 4 * stride, twiddle_factors, 4 * tw_stride);
  split_radix_fft(x + 3 * stride, y + n + n / 2, n / 4, 4 * stride, twiddle_factors, 4 * tw_stride);

  // Stitch together the output
  float u1r, u1i, u2r, u2i, x1r, x1i, x2r, x2i;
  float t;

  // We can save a few multiplications in the first step
  u1r = y[0];
  u1i = y[1];
  u2r = y[n / 2];
  u2i = y[n / 2 + 1];

  x1r = y[n];
  x1i = y[n + 1];
  x2r = y[n / 2 + n];
  x2i = y[n / 2 + n + 1];

  t = x1r + x2r;
  y[0] = u1r + t;
  y[n]     = u1r - t;

  t = x1i + x2i;
  y[1] = u1i + t;
  y[n + 1] = u1i - t;

  t = x2i - x1i;
  y[n / 2]     = u2r - t;
  y[n + n / 2]     = u2r + t;

  t = x1r - x2r;
  y[n / 2 + 1] = u2i - t;
  y[n + n / 2 + 1] = u2i + t;

  for (k = 1 ; k < n / 4 ; k++)
  {
    float u1r, u1i, u2r, u2i, x1r, x1i, x2r, x2i, c1, s1, c2, s2;
    c1 = twiddle_factors[k * tw_stride];
    s1 = twiddle_factors[k * tw_stride + 1];
    c2 = twiddle_factors[3 * k * tw_stride];
    s2 = twiddle_factors[3 * k * tw_stride + 1];

    u1r = y[2 * k];
    u1i = y[2 * k + 1];
    u2r = y[2 * k + n / 2];
    u2i = y[2 * k + n / 2 + 1];

    x1r =  c1 * y[n + 2 * k] + s1 * y[n + 2 * k + 1];
    x1i = -s1 * y[n + 2 * k] + c1 * y[n + 2 * k + 1];
    x2r =  c2 * y[n / 2 + n + 2 * k] + s2 * y[n / 2 + n + 2 * k + 1];
    x2i = -s2 * y[n / 2 + n + 2 * k] + c2 * y[n / 2 + n + 2 * k + 1];

    t = x1r + x2r;
    y[2 * k]     = u1r + t;
    y[2 * k + n]     = u1r - t;

    t = x1i + x2i;
    y[2 * k + 1] = u1i + t;
    y[2 * k + n + 1] = u1i - t;

    t = x2i - x1i;
    y[2 * k + n / 2]     = u2r - t;
    y[2 * k + n + n / 2]     = u2r + t;

    t = x1r - x2r;
    y[2 * k + n / 2 + 1] = u2i - t;
    y[2 * k + n + n / 2 + 1] = u2i + t;
  }

}


void ifft_primitive(float *input, float *output, int n, int stride, float *twiddle_factors, int tw_stride)
{

#if USE_SPLIT_RADIX
  split_radix_fft(input, output, n, stride, twiddle_factors, tw_stride);
#else
  fft_primitive(input, output, n, stride, twiddle_factors, tw_stride);
#endif

  int ks;

  int ns = n * stride;

  // reverse all coefficients from 1 to n / 2 - 1
  for (ks = stride ; ks < ns / 2 ; ks += stride)
  {
    float t;

    t = output[ks];
    output[ks] = output[ns-ks];
    output[ns-ks] = t;

    t = output[ks+1];
    output[ks+1] = output[ns-ks+1];
    output[ns-ks+1] = t;
  }

  // Apply normalization
  float norm = 1. / n;
  for (ks = 0 ; ks < ns ; ks += stride)
  {
    output[ks]   *= norm;
    output[ks+1] *= norm;
  }

}

inline void fft8(float *input, int stride_in, float *output, int stride_out)
{
  /*
   * Unrolled implementation of FFT8 for a little more performance
   */
  float a0r, a1r, a2r, a3r, a4r, a5r, a6r, a7r;
  float a0i, a1i, a2i, a3i, a4i, a5i, a6i, a7i;
  float b0r, b1r, b2r, b3r, b4r, b5r, b6r, b7r;
  float b0i, b1i, b2i, b3i, b4i, b5i, b6i, b7i;
  float t;
  float sin_pi_4 = 0.7071067812;

  a0r = input[0];
  a0i = input[1];
  a1r = input[stride_in];
  a1i = input[stride_in+1];
  a2r = input[2*stride_in];
  a2i = input[2*stride_in+1];
  a3r = input[3*stride_in];
  a3i = input[3*stride_in+1];
  a4r = input[4*stride_in];
  a4i = input[4*stride_in+1];
  a5r = input[5*stride_in];
  a5i = input[5*stride_in+1];
  a6r = input[6*stride_in];
  a6i = input[6*stride_in+1];
  a7r = input[7*stride_in];
  a7i = input[7*stride_in+1];

  // Stage 1

  b0r = a0r + a4r;
  b0i = a0i + a4i;

  b1r = a1r + a5r;
  b1i = a1i + a5i;

  b2r = a2r + a6r;
  b2i = a2i + a6i;

  b3r = a3r + a7r;
  b3i = a3i + a7i;

  b4r = a0r - a4r;
  b4i = a0i - a4i;

  b5r = a1r - a5r;
  b5i = a1i - a5i;
  // W_8^1 = 1/sqrt(2) - j / sqrt(2)
  t = b5r + b5i;
  b5i = (b5i - b5r) * sin_pi_4;
  b5r = t * sin_pi_4;

  // W_8^2 = -j
  b6r = a2i - a6i;
  b6i = a6r - a2r;

  b7r = a3r - a7r;
  b7i = a3i - a7i;
  // W_8^3 = -1 / sqrt(2) + j / sqrt(2)
  t = sin_pi_4 * (b7i - b7r);
  b7i = - (b7r + b7i) * sin_pi_4;
  b7r = t;

  // Stage 2

  a0r = b0r + b2r;
  a0i = b0i + b2i;

  a1r = b1r + b3r;
  a1i = b1i + b3i;

  a2r = b0r - b2r;
  a2i = b0i - b2i;

  // * j
  a3r = b1i - b3i;
  a3i = b3r - b1r;

  a4r = b4r + b6r;
  a4i = b4i + b6i;

  a5r = b5r + b7r;
  a5i = b5i + b7i;

  a6r = b4r - b6r;
  a6i = b4i - b6i;

  // * j
  a7r = b5i - b7i;
  a7i = b7r - b5r;

  // Stage 3

  // X[0]
  output[0] = a0r + a1r;
  output[1] = a0i + a1i;

  // X[4]
  output[4*stride_out] = a0r - a1r;
  output[4*stride_out+1] = a0i - a1i;

  // X[2]
  output[2*stride_out] = a2r + a3r;
  output[2*stride_out+1] = a2i + a3i;

  // X[6]
  output[6*stride_out] = a2r - a3r;
  output[6*stride_out+1] = a2i - a3i;

  // X[1]
  output[stride_out] = a4r + a5r;
  output[stride_out+1] = a4i + a5i;

  // X[5]
  output[5*stride_out] = a4r - a5r;
  output[5*stride_out+1] = a4i - a5i;

  // X[3]
  output[3*stride_out] = a6r + a7r;
  output[3*stride_out+1] = a6i + a7i;

  // X[7]
  output[7*stride_out] = a6r - a7r;
  output[7*stride_out+1] = a6i - a7i;

}

inline void fft4(float *input, int stride_in, float *output, int stride_out)
{
  /*
   * Unrolled implementation of FFT4 for a little more performance
   */
  float t1, t2;

  t1 = input[0] + input[2*stride_in];
  t2 = input[stride_in] + input[3*stride_in];
  output[0] = t1 + t2;
  output[2*stride_out] = t1 - t2;

  t1 = input[1] + input[2*stride_in+1];
  t2 = input[stride_in+1] + input[3*stride_in+1];
  output[1] = t1 + t2;
  output[2*stride_out+1] = t1 - t2;

  t1 = input[0] - input[2*stride_in];
  t2 = input[stride_in+1] - input[3*stride_in+1];
  output[stride_out] = t1 + t2;
  output[3*stride_out] = t1 - t2;

  t1 = input[1] - input[2*stride_in+1];
  t2 = input[3*stride_in] - input[stride_in];
  output[stride_out+1] = t1 + t2;
  output[3*stride_out+1] = t1 - t2;
}

Tutaj kod pliku FFT_signal.h:

#define FFT_N 256 // Must be a power of 2
#define TOTAL_TIME 0.008533333 //The time in which data was captured. This is equal to FFT_N/sampling_freq
//9.391904 was time

float fft_input[FFT_N];
float fft_output[FFT_N];

float max_magnitude = 0;
float fundamental_freq = 0;


/* Dummy data (Output of an accelerometer)
 * Frequency: 5 Hz
 * Amplitude: 0.25g
*/
double fft_signal[FFT_N];

A tutaj kod głównego programu FFT_on_ESP32_Arduino.ino (C++):

#include "FFT.h" // include the library
#include "FFT_signal.h"
#include <SPI.h>
#include <ADS1256.h>

float clockMHZ = 7.68; // crystal frequency used on ADS1256
float vRef = 2.5; // voltage reference

// Construct and init ADS1256 object
ADS1256 adc(clockMHZ,vRef,false); // RESETPIN is permanently tied to 3.3v

float sensor1;
bool startSampling=false;
bool endSampling=false;
volatile int licSample=0;

void setup() { //---------------------------------------------------------------------------------------
  Serial.begin(115200); // use the serial port
  Serial.println("Starting ADC");

  // start the ADS1256 with data rate of 15 SPS and gain x1
  adc.begin(ADS1256_DRATE_30000SPS,ADS1256_GAIN_1,false); 

  Serial.println("ADC Started");
  // Set MUX Register to AINO so it start doing the ADC conversion
  Serial.println("Channel set to Single end ch0");
  adc.setChannel(0);
  licSample=0;
  startSampling=true;
}//---------------------------------------------------------------------------------------

void loop() {

  if (startSampling) {
   adc.waitDRDY(); // wait for DRDY to go low before next register read
   sensor1 = adc.readCurrentChannel(); // read as voltage according to gain and vref
  
   //Serial.println(sensor1 , 10); // Print as decimal, 10 decimal places
   fft_signal[licSample]=(double)sensor1;
   licSample++;
   if (licSample>=256) {
     licSample=0;
     startSampling=false;
     endSampling=true;
   }
  }

  if (endSampling) {
    char print_buf[300];
    fft_config_t *real_fft_plan = fft_init(FFT_N, FFT_REAL, FFT_FORWARD, fft_input, fft_output);
  
    for (int k = 0 ; k < FFT_N ; k++)
      real_fft_plan->input[k] = (float)fft_signal[k];
  
    long int t1 = micros();
    // Execute transformation
    fft_execute(real_fft_plan);

    //sprintf(print_buf,"sizeFFT= %d", real_fft_plan->size);
    //Serial.println(print_buf);

    long int t2 = micros();
    
    // Print the output
    for (int k = 1 ; k < real_fft_plan->size / 2 ; k++)
    {
      /*The real part of a magnitude at a frequency is followed by the corresponding imaginary part in the output*/
      float mag = sqrt(pow(real_fft_plan->output[2*k],2) + pow(real_fft_plan->output[2*k+1],2))/1;
      float freq = k*1.0/(2*TOTAL_TIME);
      sprintf(print_buf,"%f Hz : %f", freq, mag);
      Serial.println(print_buf);
      if(mag > max_magnitude)
      {
          max_magnitude = mag;
          fundamental_freq = freq;
      }
    }   
    
    Serial.println();
    /*Multiply the magnitude of the DC component with (1/FFT_N) to obtain the DC component*/
    sprintf(print_buf,"DC component : %f g\n", (real_fft_plan->output[0])/10000/FFT_N);  // DC is at [0]
    Serial.println(print_buf);
  
    /*Multiply the magnitude at all other frequencies with (2/FFT_N) to obtain the amplitude at that frequency*/
    sprintf(print_buf,"Fundamental Freq : %f Hz\t Mag: %f g\n", fundamental_freq, (max_magnitude/10000)*2/FFT_N);
    Serial.println(print_buf);
  
    Serial.print("Time taken: ");Serial.print((t2-t1)*1.0/1000);Serial.println(" milliseconds!");
    
    // Clean up at the end to free the memory allocated
    fft_destroy(real_fft_plan);
    endSampling=false;
    }
  
}

Program jest na razie bardzo prymitywny (to była pierwsza próba, aby udowodnić wykonalność zadania na wybranym sprzęcie), teraz będę go modyfikował, aby był lepszy. Za pomocą "Analog Test Shield" i Arduino UNO nie uda się wygenerować przebiegu audio zawierającego wyższe częstotliwości (ograniczenia DAC i pojemności pamięci RAM ATMegi328), więc będę musiał znaleźć inne źródło sygnału. Szkoda, że wyjścia karty dźwiękowej z PC są u mnie trudno dostępne i w dużej odległości od montowanego układu (będę musiał zamówić zewnętrzną kartę dźwiękową ze złączem USB).

Potem planuję użycie z ESP32 przetwornika audio I2S - patrz link :

https://pl.aliexpress.com/item/32830812025.html?gatewayAdapt=glo2pol&spm=a2g0o.9042311.0.0.78f25c0fxqFC0h

Prawdopodobnie będzie to lepsze i tańsze rozwiązanie (jak sugerował kolega @ethanak). Jeżeli nie napotkam jakichś większych problemów z napisanie programu do obsługi tego przetwornika.

Cieszę się, bo osiągnąłem zamierzony cel: udowodniłem, że na wybranym sprzęcie(ADS1256 i ESP32) można śmiało policzyć Transformatę Fouriera z dostateczną dokładnością i w czasie spełniającym warunki, które były zadane.

BTW: @Gieneq wychodzi na to, że twój cały projekt udało by się zrobić na ESP32 (pod warunkiem, że nie zabraknie Ci pinów I/O na obsługę matrycy i innych potrzebnych rzeczy). Byłoby to także najtańsze rozwiązanie. Jak zrobię próby z przetwornikiem ADC I2S i ESP32, to mam zamiar powrócić do FPGA. Najpierw napiszę kod do podłączenia płytki z ADS1256 z układem FPGA (w Verilogu). Jak mam działające rozwiązanie w C++ (biblioteka) to dużo łatwiej będzie napisać kod w Verilog'u do tego celu (będzie potrzebny IP Core z SPI Master).

BTW: w kodzie głównego programu zakradł się błąd - w linii programu było:

 float freq = k*1.0/(2*TOTAL_TIME);

a powinno być:

 float freq = k*1.0/(1*TOTAL_TIME);

wtedy prążki widma będą do około 15Khz max (a nie 7 Khz). Za błąd przepraszam.

Pozdrawiam

Edytowano Marzec 9, 2022 przez FlyingDutch
update

FlyingDutch · Marzec 8, 2022

Cześć @ethanak,

mógłbyś napisać jak byś widział podłączenie przetwornika audio z I2S do ESP32. Wiem, że masz większe doświadczenie niż ja z samym ESP32.Wszelkie podpowiedzi jak się do tego zabrać powitam z radością .

Pozdrawiam

ethanak · Marzec 8, 2022

@FlyingDutch Powiem szczerze: zawsze działałem w drugą stronę (tzn. ESP32 generujący sygnał I2S), nawet miałem zamiar tym się pobawić ale jakoś tak wtedy coś dziwnego stało się z cenami

Spróbuj zerknąć tutaj: https://www.hackster.io/esikora/audio-visualization-with-esp32-i2s-mic-and-rgb-led-strip-4a251c

Wydaje mi się, że temat dość podobny...

Gieneq · Marzec 9, 2022

Cześć @FlyingDutch dopinałem ostatnio film na YT, wyszedł znośny ale jednak spędziłem przy nim kilkadziesiąt godzin. Tak to jest jak bierzesz się za coś o czym nie masz większego pojęcia, ale na ten rok mam sentencję z jednej animacji "Jak zawsze będziesz robić to co umiesz nigdy nie przeskoczysz samego siebie".

19 godzin temu, FlyingDutch napisał:

Szkoda, że wyjścia karty dźwiękowej z PC są u mnie trudno dostępne i w dużej odległości od montowanego układu (będę musiał zamówić zewnętrzną kartę dźwiękową ze złączem USB).

A masz w telefonie złącze Jack? Ja swoje układy testuję z generatorem online https://www.szynalski.com/tone-generator/

19 godzin temu, FlyingDutch napisał:

wychodzi na to, że twój cały projekt udało by się zrobić na ESP32

Wiem że jest to możliwe do tej pory tak było zrobione, tylko przez używanie wbudowanego przetwornika projekt tracił, bo miałem 5kHz sampling rate ograniczony szybkością przetwornika.

W ogóle bardzo dziękuję za chęci, widzę że równie wnikliwie usiadłeś do tematu

Co do kodu programu to przeniosę go do swojego programu i dam znać jak działa. U mnie jest już rozgraniczenie na 2 rdzenie. Kod rozbiłem na dość sporo plików:

ale spróbuję to jakoś opisać. Kod też jest na githubie: https://github.com/Gieneq/RGBHSV

main.cpp - inicjalizacja, wybór efektów:

#include <Arduino.h>
#include "presets.h"
#include "effects.h"
#include "ffthsv.h"

effect_handler_t fire_effect_handlers[] = {fire_effect};

void setup() {
    Serial.begin(115200);
    FFTHSV_begin();
    FFTHSV_set_effects(fire_effect_handlers, 1);
    FFTHSV_select_effect(0);
    FFTHSV_benchmark_params(10,20000,200,10000);
    FFTHSV_benchmark(false);
}


void loop() {
    FFTHSV_update();
}

Widać funkcje z benchmark w nazwie, one mogą Cię zainteresować, dlatego rpzechodzę do pliku ffthsv.cpp zawierający m.in zadania dla obu rdzeni:

#include "ffthsv.h"

/*
 * LEDS
 */
CRGB display_leds[DISPLAY_LEDS_COUNT];
CRGB base_leds[BASE_LEDS_COUNT];

/*
 * EFFECTS SELECTOR
 */
constexpr int EFFECTS_MAX_COUNT = 8;
effect_handler_t effects_options[EFFECTS_MAX_COUNT];
int effects_count;
int current_effect_index;

/*
 * COUNTERS
 */
long last_time_1;
int loops_counter_1;

long last_time_core1;
int loops_counter_core1;
double delta_time_core1;

long last_time_draw;
int loops_counter_draw;

long last_time_update;
int loops_counter_update;

/*
 * PHYSICS & ANIMATIONS
 */
double dt;
double velocities[BANDS_COUNT];
double gravity = 100;

constexpr int ENERGY_PARAMS_COUNT = 3;
double energy[ENERGY_PARAMS_COUNT];

double bars_threshold[BANDS_COUNT];
double bars_pos[BANDS_COUNT];
double diffs[BANDS_COUNT];

int bars[BANDS_COUNT];  // rendering

/*
 * BURST
 */
constexpr int BURST_BUFFER_LEN = 16;  // eval differential
int burst_buffer[BURST_BUFFER_LEN];
int burst_end_index;
int burst_front_index = BURST_BUFFER_LEN - 1;

/*
 * EFFT EVAL
 */
int test_frequency = 1; // Hz
constexpr int READING_BUFFER_SIZE = 512 * 2;
volatile int adc_readings_counter;
volatile int readings_buffer[READING_BUFFER_SIZE];
volatile int reading_buffer_index = 0;
hw_timer_t* sampling_timer = NULL;
portMUX_TYPE sampling_timerMux = portMUX_INITIALIZER_UNLOCKED;
volatile double mags[BANDS_COUNT];

/*
 * BENCHMARK
 */
bool benchmarking = false;
long benchmark_time_ms = 10000;
long benchmark_step_interval = 10;
int benchmark_steps_count = 1000;
int benchmark_step;
double benchmark_min_freq = 10;
double benchmark_max_freq = 20000;
long last_benchmark_time;
double benchmark_exp_factor = 0.001763;

/*
 * GENERATE SAMPLES
 */
static inline int calc_test_sample(int amplitude, int sample_index, const int samples_count, const double freq) {
    return int(amplitude * sin(2.0*PI*freq*sample_index / SAMPLING_FREQUENCY)); // spacing t = 1/SF
}

/*
 * USED TO PREVENT AUTORESET
 */
static void feedTheDog() {
    // feed dog 0
    TIMERG0.wdt_wprotect = TIMG_WDT_WKEY_VALUE;  // write enable
    TIMERG0.wdt_feed = 1;                        // feed dog
    TIMERG0.wdt_wprotect = 0;                    // write protect
    // feed dog 1
    TIMERG1.wdt_wprotect = TIMG_WDT_WKEY_VALUE;  // write enable
    TIMERG1.wdt_feed = 1;                        // feed dog
    TIMERG1.wdt_wprotect = 0;                    // write protect
}

/*
 * CALC ALL BAR HEIGHTS
 */
static void calc_movement(double* thrs, double* pos, double* vel, double grav,
                   int datalen, double dt, double* diffs) {
    for (int i = 0; i < BANDS_COUNT; i++) {
        pos[i] += dt * vel[i];

        if (pos[i] < 0.0) {
            pos[i] = 0.0;
            vel[i] = 0.0;
        }

        if (pos[i] > 21.0) {
            pos[i] = 21.0;
            vel[i] = 0.0;
        }
    }

    for (int i = 0; i < BANDS_COUNT; i++) {
        thrs[i] = constrain(thrs[i], 0.0, 21.0);  // target
        diffs[i] = thrs[i] - pos[i];
        // drag up
        if (diffs[i] >= 0.0) {
            pos[i] = thrs[i];
            vel[i] = 0;
        } else
            vel[i] += -grav * dt - 0.5;  // gravity drop
    }
}

/*
 * CALC BURST OF BRIGHTNES
 */
static void calc_energy_eval(int* bar_heights, double* bar_thrs, double* energy) {
    int sum_bar_thrs = 0;
    for (int display_x = 0; display_x < 3; display_x++)
        sum_bar_thrs += bar_heights[display_x];

    ++burst_end_index %= BURST_BUFFER_LEN;
    ++burst_front_index %= BURST_BUFFER_LEN;

    energy[energy_t::ENERGY_CURRENT] = double(sum_bar_thrs) / 3;
    burst_buffer[burst_end_index] = *energy;
    energy[energy_t::ENERGY_BURST] =
        burst_buffer[burst_end_index] - burst_buffer[burst_front_index];

    energy[energy_t::ENERGY_MEAN] = 0.0;
    for (int burst_i = 0; burst_i < BURST_BUFFER_LEN; burst_i++)
        energy[energy_t::ENERGY_MEAN] += burst_buffer[burst_i];
    energy[energy_t::ENERGY_MEAN] /= BURST_BUFFER_LEN;

#if DEBUG_FPS == 1
    Serial.print("TB:");
    Serial.print(total_brightness / 10.0);
    Serial.print(", ");
    Serial.print("Energy:");
    Serial.print(energy);
    Serial.print(", ");
    Serial.print("Mean_energy:");
    Serial.print(mean_energy);
    Serial.print(", ");
    Serial.print("Burst_energy:");
    Serial.println(burst_energy);
#endif
}

/*
 * CORE 0 INTERRUPT - SAMPLING ADC
 */
static void IRAM_ATTR on_timer_sample() {
    portENTER_CRITICAL_ISR(&sampling_timerMux);
    // adc_readings_counter++; //used in printing, meh
    readings_buffer[reading_buffer_index] = test_frequency <= 0 ? analogRead(A6) : calc_test_sample(512, reading_buffer_index, SAMPLES_COUNT, test_frequency);
    reading_buffer_index++;
    if (reading_buffer_index >= READING_BUFFER_SIZE) 
        reading_buffer_index = 0;
    portEXIT_CRITICAL_ISR(&sampling_timerMux);
}

/*
 * CORE 0 - SETUP, SAMPLE, CALC FFT
 */
static void core0_task(void* parameter) {
    Serial.print("FFT sampler/processor running on core ");
    Serial.println(xPortGetCoreID());
    // ok, inside infinite function
    static long last_time = millis();
    static int loops_counter = 0;
    static double delta_time;

    sampling_timer = timerBegin(0, 80, true);
    timerAttachInterrupt(sampling_timer, &on_timer_sample, true);
    timerAlarmWrite(sampling_timer, 25, true);  // 40 kHz
    timerAlarmEnable(sampling_timer);

    for (;;) {
        load_samples(reading_buffer_index, readings_buffer, READING_BUFFER_SIZE);
        calculate_bars();
        double* fft_bars = get_magnitudes();
        for (int i = 0; i < BANDS_COUNT; i++) 
            mags[i] = fft_bars[i];

#if DEBUG_FPS == 1
        // eval frequency of FFT + sampling
        delta_time = millis() - last_time;
        if (delta_time > 4000) {
            delta_time /= max(loops_counter, 1);
            Serial.print("Core ");
            Serial.print(xPortGetCoreID());
            Serial.print(" FPS=");
            Serial.print(int(1000.0 / delta_time));
            Serial.print("Hz, buffidx=");
            Serial.print(reading_buffer_index);
            Serial.print(", sr=");
            Serial.print(int(adc_readings_counter / 4.0));
            Serial.println("Hz");
            Serial.println();
            last_time = millis();
            adc_readings_counter = 0;
            loops_counter = 0;
        }
#endif

        feedTheDog();
        loops_counter++;
    }
}

/*
 * CORE 1 - SETUP
 */
void FFTHSV_begin() {
    FastLED
        .addLeds<LED_TYPE, DISPLAY_PIN, COLOR_ORDER>(display_leds,
                                                     DISPLAY_LEDS_COUNT)
        .setCorrection(TypicalSMD5050);
    FastLED.addLeds<LED_TYPE, BASE_PIN, COLOR_ORDER>(base_leds, BASE_LEDS_COUNT)
        .setCorrection(TypicalSMD5050);
    FastLED.setBrightness(MIN_BRIGHTNESS);
    FastLED.show();

    Serial.println(F("FFT Spectrum Display"));
    Serial.println(xPortGetCoreID());

    xTaskCreatePinnedToCore(core0_task, "FFT_processor", 64000, NULL, 1, NULL,
                            0);
    last_time_core1 = millis();
    last_time_draw = millis();
}

/*
 * CORE 1 - CALC EFFECTS & DRAW
 */
void FFTHSV_update() {
    //BENCHMARK
    if(benchmarking){
        if (millis() - last_benchmark_time > benchmark_step_interval) {
            int benchmark_freq = benchmark_min_freq * pow(benchmark_exp_factor, benchmark_step);
            FFTHSV_set_test_frequency(benchmark_freq);
            Serial.print("Benchmark[Hz]:");
            Serial.println(benchmark_freq);
            ++benchmark_step%=benchmark_steps_count;
            last_benchmark_time += benchmark_step_interval;
        }
    }

    // UPDATE
    if (millis() - last_time_update > UPDATE_INTERVAL_MS) {
        dt = (millis() - last_time_update) / 1000.0;  // secs

        for (int i = 0; i < BANDS_COUNT; i++) bars_threshold[i] = mags[i];

        calc_movement(bars_threshold, bars_pos, velocities, gravity,
                      BANDS_COUNT, dt, diffs);

        for (int i = 0; i < BANDS_COUNT; i++)
            bars[i] = constrain(int(bars_pos[i]), 0, 21);

        calc_energy_eval(bars, bars_threshold, energy);
        effects_options[current_effect_index](bars, bars_threshold, energy,
                                              display_leds, base_leds);

        last_time_update += UPDATE_INTERVAL_MS;
    }

    // DRAW
    if (millis() - last_time_draw > DRAW_INTERVAL_MS) {
        FastLED.show();
        last_time_draw += DRAW_INTERVAL_MS;
    }
}

void FFTHSV_set_effects(const effect_handler_t* effects, const int count) {
    effects_count = 0;
    for (int i = 0; i < count; i++) {
        effects_options[effects_count++] = effects[i];
        if (effects_count >= EFFECTS_MAX_COUNT) return;
    }
}

void FFTHSV_select_effect(int index) { current_effect_index = index; }

void FFTHSV_set_test_frequency(int freq){
    test_frequency = freq;
}

void FFTHSV_benchmark_params(double min_freq, double max_freq, int steps, long evaluation_time_ms) {
    benchmark_time_ms = evaluation_time_ms;
    benchmark_steps_count = steps;
    benchmark_step_interval = long(1.0*evaluation_time_ms/steps);
    benchmark_step = 0;
    benchmark_min_freq = min_freq;
    benchmark_max_freq = max_freq;
    benchmark_exp_factor = pow(max_freq/min_freq, 1.0/benchmark_steps_count);
}


void FFTHSV_benchmark(bool state) {
    benchmarking = state;
    if(state == true){
        last_benchmark_time = millis();
        test_frequency = benchmark_min_freq;
    } else {
        benchmark_step = -1;
        test_frequency = 0;
    }
}

// bool benchmarking = false;
// long benchmark_time_ms = 10000;
// int benchmark_steps_count = 1000;
// int benchmark_step;
// double benchmark_min_freq = 10;
// double benchmark_max_freq = 20000;
// long last_benchmark_time;
// double benchmark_exp_factor = 0.001763;

Teraz opiszę to co jest w tych 300 liniach kodu:

static inline int calc_test_sample(int amplitude, int sample_index, const int samples_count, const double freq) {
    return int(amplitude * sin(2.0*PI*freq*sample_index / SAMPLING_FREQUENCY)); // spacing t = 1/SF
}

Bardzo prosta funkcja - wybieram częstotliwość i generuję próbki. Składając wyniki funkcji można by generować zmieszane harmoniczne.

Funkcja ta występuje w zadania timera dla rdzenia 0 (ten co zazwyczaj ma obsługę WiFi)

static void IRAM_ATTR on_timer_sample() {
    portENTER_CRITICAL_ISR(&sampling_timerMux);
    // adc_readings_counter++; //used in printing, meh
    readings_buffer[reading_buffer_index] = test_frequency <= 0 ? analogRead(A6) : calc_test_sample(512, reading_buffer_index, SAMPLES_COUNT, test_frequency);
    reading_buffer_index++;
    if (reading_buffer_index >= READING_BUFFER_SIZE) 
        reading_buffer_index = 0;
    portEXIT_CRITICAL_ISR(&sampling_timerMux);
}

Jeżeli zmienna mająca testową częstotliwość zawiera bzdury to korzystam z ADC. Tu jest też mój mechanizm nieczekania na pomairy - bufor jest 2x większy, a indeks wskazuje wolną komórkę tablicy do zapisu.

W tym miejscu odpalam timer:

static void core0_task(void* parameter) {
    Serial.print("FFT sampler/processor running on core ");
    Serial.println(xPortGetCoreID());
    // ok, inside infinite function
    static long last_time = millis();
    static int loops_counter = 0;
    static double delta_time;

    sampling_timer = timerBegin(0, 80, true);
    timerAttachInterrupt(sampling_timer, &on_timer_sample, true);
    timerAlarmWrite(sampling_timer, 25, true);  // 40 kHz
    timerAlarmEnable(sampling_timer);

    for (;;) {
        load_samples(reading_buffer_index, readings_buffer, READING_BUFFER_SIZE);
        calculate_bars();
        double* fft_bars = get_magnitudes();
        for (int i = 0; i < BANDS_COUNT; i++) 
            mags[i] = fft_bars[i];

#if DEBUG_FPS == 1
        // eval frequency of FFT + sampling
        delta_time = millis() - last_time;
        if (delta_time > 4000) {
            delta_time /= max(loops_counter, 1);
            Serial.print("Core ");
            Serial.print(xPortGetCoreID());
            Serial.print(" FPS=");
            Serial.print(int(1000.0 / delta_time));
            Serial.print("Hz, buffidx=");
            Serial.print(reading_buffer_index);
            Serial.print(", sr=");
            Serial.print(int(adc_readings_counter / 4.0));
            Serial.println("Hz");
            Serial.println();
            last_time = millis();
            adc_readings_counter = 0;
            loops_counter = 0;
        }
#endif

        feedTheDog();
        loops_counter++;
    }
}

a w pętli ładuję wycinek z zapisanych próbek taką funkcją z pliku analyser.cpp :

void load_samples(int read_index, volatile int* read_buffer, const int buffersize) {
    for (int i = 0; i < SAMPLES_COUNT; i++) {
        vReal[i] = read_buffer[read_index--]; //DC removal  - 1800
        vImag[i] = 0;
        if(read_index < 0)
            read_index = buffersize-1;
    }
}

Pewnie można by to jakoś przyspieszyć - istota jest taka, że wycinek może składać się z końcówki tablicy i jej początku - więc potrzeba by 2 wskaźniki.

Wracajac, to co jest w głównej pętli core0_task pobiera próbki (nieważne jak wygenerowane) i wyznacza fft w funkcji:

void calculate_bars() {
    FFT.DCRemoval();
    // FFT.Windowing(FFT_WIN_TYP_NUTTALL, FFT_FORWARD); //nieglupie
    FFT.Windowing(FFT_WIN_TYP_HAMMING, FFT_FORWARD); //nieglupie

    FFT.Compute(FFT_FORWARD);

    FFT.ComplexToMagnitude();  // wynik w vReal
  
    for (int i = 0; i < BANDS_COUNT; i++)
      magnitude_bands[i] = 0.0;

    for (int i = 1; i < USABLE_SAMPLES_COUNT; i++) {
        if (vReal[i] > NOISE_LIMIT)
            magnitude_bands[bands512[i]] += vReal[i];
    }
    for (int i = 0; i < BANDS_COUNT; i++)
      magnitude_bands[i] *= BAR_HEIGHT_MODIFIER;
}

bazującej na bibliotece.

Teraz funkcja z rdzenia 1, która jest wykonywana domyślnie i inicjuje to co jest na rdzeniu 0:

void FFTHSV_begin() {
    FastLED
        .addLeds<LED_TYPE, DISPLAY_PIN, COLOR_ORDER>(display_leds,
                                                     DISPLAY_LEDS_COUNT)
        .setCorrection(TypicalSMD5050);
    FastLED.addLeds<LED_TYPE, BASE_PIN, COLOR_ORDER>(base_leds, BASE_LEDS_COUNT)
        .setCorrection(TypicalSMD5050);
    FastLED.setBrightness(MIN_BRIGHTNESS);
    FastLED.show();

    Serial.println(F("FFT Spectrum Display"));
    Serial.println(xPortGetCoreID());

    xTaskCreatePinnedToCore(core0_task, "FFT_processor", 64000, NULL, 1, NULL,
                            0);
    last_time_core1 = millis();
    last_time_draw = millis();
}

Za wyjatkiem składni przypięcia funkcji do rdzenia, raczej nic ciekawego.

Tu w zasadzei kończy się próbkowanie i FFT, jeszcze jest obsługa efektów. Do tego mam funkcję update która liczy różne parametryi później obsługuje funkcję przypisaną jako efekt/animacja:

void FFTHSV_update() {
    //BENCHMARK
    if(benchmarking){
        if (millis() - last_benchmark_time > benchmark_step_interval) {
            int benchmark_freq = benchmark_min_freq * pow(benchmark_exp_factor, benchmark_step);
            FFTHSV_set_test_frequency(benchmark_freq);
            Serial.print("Benchmark[Hz]:");
            Serial.println(benchmark_freq);
            ++benchmark_step%=benchmark_steps_count;
            last_benchmark_time += benchmark_step_interval;
        }
    }

    // UPDATE
    if (millis() - last_time_update > UPDATE_INTERVAL_MS) {
        dt = (millis() - last_time_update) / 1000.0;  // secs

        for (int i = 0; i < BANDS_COUNT; i++) bars_threshold[i] = mags[i];

        calc_movement(bars_threshold, bars_pos, velocities, gravity,
                      BANDS_COUNT, dt, diffs);

        for (int i = 0; i < BANDS_COUNT; i++)
            bars[i] = constrain(int(bars_pos[i]), 0, 21);

        calc_energy_eval(bars, bars_threshold, energy);
        effects_options[current_effect_index](bars, bars_threshold, energy,
                                              display_leds, base_leds);

        last_time_update += UPDATE_INTERVAL_MS;
    }

    // DRAW
    if (millis() - last_time_draw > DRAW_INTERVAL_MS) {
        FastLED.show();
        last_time_draw += DRAW_INTERVAL_MS;
    }
}

effects_options jest typu:

typedef void (*effect_handler_t)(const int* bar_heights, const double* bar_thrs, double* energy, CRGB* display_leds, CRGB* base_leds);

i przypisałem jedną funkcję na starcie w pliku main.cpp:

FFTHSV_set_effects(fire_effect_handlers, 1);

funkcją:

void FFTHSV_set_effects(const effect_handler_t* effects, const int count) {
    effects_count = 0;
    for (int i = 0; i < count; i++) {
        effects_options[effects_count++] = effects[i];
        if (effects_count >= EFFECTS_MAX_COUNT) return;
    }
}

Funkcje efektów zapisałem w pliku effects, którego nagłówek wygląda tak:

#pragma once
#include "Arduino.h"
#include "FastLED.h"
#include "presets.h"

constexpr int NOISE_FLAMES_count = 128;
const uint8_t NOISE_FLAMES[] = {199, 218, 234, 245, 252, 255, 253, 248, 239, 228, 216, 204, 191, 180, 170, 163, 157, 154, 152, 152, 153, 155, 157, 159, 160, 160, 159, 157, 154, 151, 147, 143, 139, 136, 134, 133, 133, 134, 136, 138, 140, 143, 144, 146, 146, 145, 143, 140, 136, 131, 127, 122, 119, 116, 114, 114, 114, 116, 119, 123, 127, 131, 135, 137, 139, 139, 138, 136, 132, 128, 123, 119, 115, 112, 110, 110, 111, 114, 118, 124, 130, 137, 143, 148, 151, 153, 153, 151, 147, 141, 135, 127, 120, 114, 108, 105, 104, 105, 109, 115, 122, 131, 140, 148, 155, 160, 163, 162, 157, 150, 139, 126, 110, 94, 79, 64, 53, 44, 40, 40, 45, 55, 69, 87, 108, 131, 155, 178};

void fire_effect(const int* bar_heights, const double* bar_thrs, double* energy, CRGB* display_leds, CRGB* base_leds);

Ta ogromna tablica, to złączenie kilku częstotliwości wygenerowane w Pythonie - o tym już psialiśmy. Sama funkcja efektu jest dość obszerna, nie będę się w to zagłębiał:

#include "effects.h"

/*
 * TOOLS
 */
inline int get_idx(int x, int y) {
  x = DISPLAY_W - x - 1;
  #if ZIGZAGING == 1
      return (x&1 ? DISPLAY_H - y - 1 : y) + x * DISPLAY_H; //todo - tablica?
  #else
      return y + x * DISPLAY_H;
  #endif
}

void inline set_hsv(CRGB* hsv_leds, int x, int y, uint8_t h, uint8_t s, uint8_t v) {
  hsv_leds[get_idx(x, y)] = CHSV(h,s,v);
}

  /*
  * FIRE
  */
  constexpr uint8_t FIRE_TONGUE_HUE[] = {12, 16, 32, 33, 34, 35, 36, 37, 37, 37, 37, 36, 35, 34, 33, 20, 19, 18, 255, 255, 255};
  int burnt[DISPLAY_LEDS_COUNT];
  constexpr uint8_t FIRE_HUE = 22; 
  constexpr uint8_t BURNT_HSV_VALUE = 32;
  constexpr uint8_t BURNT_HUE = 254;
  constexpr long BURN_INTERVAL = 300;
  long last_burn_time;
  extern const uint8_t NOISE_FLAMES[];
  int fire_twinkle_counter;
  constexpr int FIRE_BASE_TICS = 6;
  constexpr int FIRE_RANDOM_TICS = 24;
  int total_brightness = MIN_BRIGHTNESS;
  long last_burst_time;
  constexpr int BRURST_INTERVAL = 2;

  void fire_effect(const int* bar_heights, const double* bar_thrs, double* energy, CRGB* display_leds, CRGB* base_leds){
      ++fire_twinkle_counter;

      static int reduce_burn;
      if(millis() - last_burn_time > BURN_INTERVAL){
        reduce_burn = 1;
        last_burn_time += BURN_INTERVAL;
      }

    int base_hue = constrain(int(energy[energy_t::ENERGY_MEAN]), 0, DISPLAY_H);
    for(int base_x = 0; base_x < BASE_LEDS_COUNT; base_x++)
      base_leds[base_x] = CHSV(FIRE_TONGUE_HUE[base_hue], 255,255);
    

      if(total_brightness < TOP_BRIGHTNES/2) {
        if(energy[energy_t::ENERGY_CURRENT] - energy[energy_t::ENERGY_MEAN] > 4)
          total_brightness = TOP_BRIGHTNES;
      }


      if (millis() - last_burst_time > BRURST_INTERVAL) {
        if(total_brightness > MIN_BRIGHTNESS)
          total_brightness-=5;
        if(total_brightness < MIN_BRIGHTNESS)
          total_brightness = MIN_BRIGHTNESS;

        last_burst_time += BRURST_INTERVAL;
      }

      FastLED.setBrightness(total_brightness);

      for(int display_x = 0; display_x < DISPLAY_W; display_x++) {
        for(int display_y = 0; display_y < DISPLAY_H; display_y++)
          set_hsv(display_leds, display_x, display_y, 0,0,0);

        for(int display_y = 0; display_y < bar_heights[display_x]; display_y++) {
          int factor = DISPLAY_H - (bar_heights[display_x] - display_y - 1) - 1; // dla BH=4, 17, 18, 19, 20
          uint8_t fire_hue = FIRE_TONGUE_HUE[factor];

          uint8_t fire_saturation = 255;
          if(factor < 10)
            fire_saturation -= (10-factor-1) * 28;

          uint8_t fire_value = 255;
          if(factor > 9 && factor < 21){
            fire_value = NOISE_FLAMES[(fire_twinkle_counter + display_y + display_x*6)%128];
            
            if(factor > 9 && factor < 17){
              //factor: 10...16
              double mix_factor = (factor - 10)/6.0;
              fire_value = (mix_factor * fire_value + (1.0 - mix_factor) * 255); //mix
              if(fire_value>255)
                fire_value = 255;
            }
          }

          set_hsv(display_leds, display_x, display_y, fire_hue, fire_saturation, fire_value);
          burnt[get_idx(display_x, display_y)] = FIRE_BASE_TICS + random(FIRE_RANDOM_TICS);
        } 

        for(int display_y = bar_heights[display_x]; display_y < DISPLAY_H; display_y++) {
          int idx = get_idx(display_x, display_y);
          if(burnt[idx] > 0){
            set_hsv(display_leds, display_x, display_y, BURNT_HUE, 255, BURNT_HSV_VALUE);
            if(reduce_burn == 1)
              burnt[idx]--;
          } else
            set_hsv(display_leds, display_x, display_y, 0, 0, 0);
        }
      }
  }

Jeszcze wracając do benchmarka czyli generowania sygnału. Po przełożeniu wyników FFT na kubełki, charakterystyka jest wykładnicza, zatem i benchmark też powinien być. Przy ustawianiu wartości od razu przeliczam z funkcją pow:

void FFTHSV_benchmark_params(double min_freq, double max_freq, int steps, long evaluation_time_ms) {
    benchmark_time_ms = evaluation_time_ms;
    benchmark_steps_count = steps;
    benchmark_step_interval = long(1.0*evaluation_time_ms/steps);
    benchmark_step = 0;
    benchmark_min_freq = min_freq;
    benchmark_max_freq = max_freq;
    benchmark_exp_factor = pow(max_freq/min_freq, 1.0/benchmark_steps_count);
}

Nie jestem zawodowym programistą, ale kod chyba jest dość dobry.

@FlyingDutch czyli plan mam taki, podłączam przetwornik, dodaję twój kod z jego obsługą i sprawdzam jak sobie poradzi. Ale mam jeszcze pomysł - rozumiem że ESP sobie świetnie poradzi, ale gdyby tak zostawić WiFi? WiFi zabierze rdzeń 0, a na 1 rdzeniu tego wszystkiego się nie zrobi. Jak ostatnio liczyłem, obsługa 399 diod zajmuje 12,02 ms na samo wygenerowanie przebiegu - nie wiem jak to jest realizowane, ale na rdzeniu 1 nie udało mi się wyjść ponad 80 Hz wiec chyba działa blokująco. Mam głupi pomysł, żeby dodać WiFi dokładając dodatkowy układ ESP12F, ale to jest zbrodnia. ESP8266 nie ma zabezpeiczeń, a ESP32 ma tylko, że ma też potencjał do obliczeń. Mogę teżdać 2 ESP32, a co tam

FlyingDutch · Marzec 9, 2022

1 godzinę temu, Gieneq napisał:

Cześć @FlyingDutch dopinałem ostatnio film na YT, wyszedł znośny ale jednak spędziłem przy nim kilkadziesiąt godzin. Tak to jest jak bierzesz się za coś o czym nie masz większego pojęcia, ale na ten rok mam sentencję z jednej animacji "Jak zawsze będziesz robić to co umiesz nigdy nie przeskoczysz samego siebie".

A masz w telefonie złącze Jack? Ja swoje układy testuję z generatorem online https://www.szynalski.com/tone-generator/

Wiem że jest to możliwe do tej pory tak było zrobione, tylko przez używanie wbudowanego przetwornika projekt tracił, bo miałem 5kHz sampling rate ograniczony szybkością przetwornika.

...

@FlyingDutch czyli plan mam taki, podłączam przetwornik, dodaję twój kod z jego obsługą i sprawdzam jak sobie poradzi. Ale mam jeszcze pomysł - rozumiem że ESP sobie świetnie poradzi, ale gdyby tak zostawić WiFi? WiFi zabierze rdzeń 0, a na 1 rdzeniu tego wszystkiego się nie zrobi. Jak ostatnio liczyłem, obsługa 399 diod zajmuje 12,02 ms na samo wygenerowanie przebiegu - nie wiem jak to jest realizowane, ale na rdzeniu 1 nie udało mi się wyjść ponad 80 Hz wiec chyba działa blokująco. Mam głupi pomysł, żeby dodać WiFi dokładając dodatkowy układ ESP12F, ale to jest zbrodnia. ESP8266 nie ma zabezpeiczeń, a ESP32 ma tylko, że ma też potencjał do obliczeń. Mogę teżdać 2 ESP32, a co tam

Cześć @Gieneq,

dzięki za link do generatora. Odnośnie drugiego tematu, to może warto spróbować zainstalować na ESP32 RTOS'a wtedy zadania będziesz dzielił na "wątki" a nie rdzenie, a przydziałem wątków do rdzeni zajmie się OS. Może wtedy wszystko udało by się zrobić na jednym ESP32.

Pozdrawiam

Gieneq · Marzec 9, 2022

@FlyingDutch ciekawy pomysł, akurat wczoraj YT podpowiedział mi kanał gościa który opowiada o FreeRTOS na Pico. Trochę jestem zielony w tym temacie, ale może to być to

Tylko czy komponenty ESPHome da się jakoś z tym zgrać. No trudno, będę czytał.

Sipeed Tang Nano 4k z ADC - Gowin FPGA Designer

Pomocna odpowiedź

ethanak

FlyingDutch

ethanak

FlyingDutch

FlyingDutch

ethanak

FlyingDutch

ethanak

FlyingDutch

FlyingDutch

FlyingDutch

ethanak

Gieneq

FlyingDutch

Gieneq

Bądź aktywny - zaloguj się lub utwórz konto!

Utwórz konto w ~20 sekund!

Zaloguj się

Quizy

Najnowsze posty w innych tematach

Ścieżka rozwoju dla przyszłego architekta embedded. Szukam dwóch książek.

Domowy serwer dla 6 użytkowników

Lista osób zainteresowanych płatnymi zleceniami z zakresu elektroniki i programowania

Ciekawe filmy (elektronika, mechanika, robotyka) 1 2 3 4 63

Układ ze starego radia...

Inne