Sound2Light auf RGB(W)-LED-Leisten
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

fft.c 3.6KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151
  1. /*
  2. * vim: sw=2 ts=2 expandtab
  3. *
  4. * THE PIZZA-WARE LICENSE" (derived from "THE BEER-WARE LICENCE"):
  5. * <cfr34k@tkolb.de> wrote this file. As long as you retain this notice you can
  6. * do whatever you want with this stuff. If we meet some day, and you think
  7. * this stuff is worth it, you can buy me a pizza in return. - Thomas Kolb
  8. */
  9. #include <stdio.h>
  10. #include <math.h>
  11. #include <time.h>
  12. #include <stdint.h>
  13. #include "config.h"
  14. #include "fft.h"
  15. double hanning_buffer[BLOCK_LEN];
  16. int lookup_table[BLOCK_LEN];
  17. void init_fft(void) {
  18. int i = 0;
  19. int ri, b;
  20. for(i = 0; i < BLOCK_LEN; i++) {
  21. hanning_buffer[i] = 0.5 * (1 - cos(2 * M_PI * i / BLOCK_LEN));
  22. }
  23. for(i = 0; i < BLOCK_LEN; i++) {
  24. ri = 0;
  25. for(b = 0; b < FFT_EXPONENT; b++)
  26. ri |= ((i >> b) & 1) << (FFT_EXPONENT - b - 1);
  27. lookup_table[i] = ri;
  28. }
  29. }
  30. void complex_to_absolute(double *re, double *im, double *result) {
  31. int i;
  32. for(i = 0; i < DATALEN; i++)
  33. {
  34. result[i] = sqrt( re[i]*re[i] + im[i]*im[i] );
  35. }
  36. }
  37. void apply_hanning(sample *dftinput) {
  38. int i;
  39. for(i = 0; i < BLOCK_LEN; i++)
  40. dftinput[i] *= hanning_buffer[i];
  41. }
  42. void fft_transform(sample *samples, double *resultRe, double *resultIm) {
  43. int i;
  44. int layer, part, element;
  45. int num_parts, num_elements;
  46. int left, right;
  47. double x_left_re, x_left_im, x_right_re, x_right_im;
  48. double param;
  49. double sinval, cosval;
  50. // re-arrange the input array according to the lookup table
  51. // and store it into the real output array (as the input is obviously real).
  52. // zero the imaginary output
  53. for(i = 0; i < BLOCK_LEN; i++)
  54. {
  55. resultRe[lookup_table[i]] = samples[i];
  56. resultIm[i] = 0;
  57. }
  58. // walk layers
  59. for(layer = 0; layer < FFT_EXPONENT; layer++)
  60. {
  61. // number of parts in current layer
  62. num_parts = 1 << (FFT_EXPONENT - layer - 1);
  63. // walk parts of layer
  64. for(part = 0; part < num_parts; part++)
  65. {
  66. // number of elements in current part
  67. num_elements = (1 << layer);
  68. // walk elements in part
  69. for(element = 0; element < num_elements; element++)
  70. {
  71. // calculate index of element in left and right half of part
  72. left = (1 << (layer + 1)) * part + element;
  73. right = left + (1 << layer);
  74. // buffer the elements for the calculation
  75. x_left_re = resultRe[left];
  76. x_left_im = resultIm[left];
  77. x_right_re = resultRe[right];
  78. x_right_im = resultIm[right];
  79. // precalculate the parameter for sin and cos
  80. param = -M_PI * element / (1 << layer);
  81. // precalculate sinus and cosinus values for param
  82. sinval = sin(param);
  83. cosval = cos(param);
  84. // combine the values according to a butterfly diagram
  85. resultRe[left] = x_right_re + x_left_re * cosval - x_left_im * sinval;
  86. resultIm[left] = x_right_im + x_left_im * cosval + x_left_re * sinval;
  87. resultRe[right] = x_right_re - x_left_re * cosval + x_left_im * sinval;
  88. resultIm[right] = x_right_im - x_left_im * cosval - x_left_re * sinval;
  89. }
  90. }
  91. }
  92. }
  93. uint32_t find_loudest_frequency(double *absFFT) {
  94. int maxPos = 0;
  95. double maxVal = 0;
  96. int i;
  97. for(i = 0; i < BLOCK_LEN; i++) {
  98. if(absFFT[i] > maxVal) {
  99. maxPos = i;
  100. maxVal = absFFT[i];
  101. }
  102. }
  103. return (double)maxPos * SAMPLE_RATE / BLOCK_LEN;
  104. }
  105. double get_energy_in_band(double *fft, uint32_t minFreq, uint32_t maxFreq) {
  106. int firstBlock = minFreq * BLOCK_LEN / SAMPLE_RATE;
  107. int lastBlock = maxFreq * BLOCK_LEN / SAMPLE_RATE;
  108. int i;
  109. double energy = 0;
  110. for(i = firstBlock; i < lastBlock; i++) {
  111. energy += fft[i];
  112. }
  113. return energy;
  114. }