From b36336c12d4bca9aedce295934e56bc0a885dce2 Mon Sep 17 00:00:00 2001 From: pranabendra Date: Sun, 7 Jun 2020 20:51:03 +0530 Subject: [PATCH 01/12] Fix DCRemoval(): for all samples in time --- src/arduinoFFT.cpp | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/arduinoFFT.cpp b/src/arduinoFFT.cpp index f552fc9..580f8ae 100644 --- a/src/arduinoFFT.cpp +++ b/src/arduinoFFT.cpp @@ -199,13 +199,13 @@ void arduinoFFT::DCRemoval() { // calculate the mean of vData double mean = 0; - for (uint16_t i = 1; i < ((this->_samples >> 1) + 1); i++) + for (uint16_t i = 0; i < this->_samples; i++) { mean += this->_vReal[i]; } mean /= this->_samples; // Subtract the mean from vData - for (uint16_t i = 1; i < ((this->_samples >> 1) + 1); i++) + for (uint16_t i = 0; i < this->_samples; i++) { this->_vReal[i] -= mean; } @@ -216,13 +216,13 @@ void arduinoFFT::DCRemoval(double *vData, uint16_t samples) // calculate the mean of vData #warning("This method is deprecated and may be removed on future revisions.") double mean = 0; - for (uint16_t i = 1; i < ((samples >> 1) + 1); i++) + for (uint16_t i = 0; i < samples; i++) { mean += vData[i]; } mean /= samples; // Subtract the mean from vData - for (uint16_t i = 1; i < ((samples >> 1) + 1); i++) + for (uint16_t i = 0; i < samples; i++) { vData[i] -= mean; } From ee0459d537fd67c202393f4d4d28152278158954 Mon Sep 17 00:00:00 2001 From: Enrique Condes Date: Tue, 6 Oct 2020 14:47:01 -0500 Subject: [PATCH 02/12] Fix for issue 56 --- Examples/FFT_01/FFT_01.ino | 1 + Examples/FFT_02/FFT_02.ino | 1 + Examples/FFT_03/FFT_03.ino | 1 + Examples/FFT_04/FFT_04.ino | 2 ++ Examples/FFT_05/FFT_05.ino | 1 + 5 files changed, 6 insertions(+) diff --git a/Examples/FFT_01/FFT_01.ino b/Examples/FFT_01/FFT_01.ino index caefe80..ffeb852 100644 --- a/Examples/FFT_01/FFT_01.ino +++ b/Examples/FFT_01/FFT_01.ino @@ -53,6 +53,7 @@ double vImag[samples]; void setup() { Serial.begin(115200); + while(!Serial); Serial.println("Ready"); } diff --git a/Examples/FFT_02/FFT_02.ino b/Examples/FFT_02/FFT_02.ino index cb97e29..b621db3 100644 --- a/Examples/FFT_02/FFT_02.ino +++ b/Examples/FFT_02/FFT_02.ino @@ -53,6 +53,7 @@ unsigned long time; void setup() { Serial.begin(115200); + while(!Serial); Serial.println("Ready"); exponent = FFT.Exponent(samples); } diff --git a/Examples/FFT_03/FFT_03.ino b/Examples/FFT_03/FFT_03.ino index 9e1640c..1e59279 100644 --- a/Examples/FFT_03/FFT_03.ino +++ b/Examples/FFT_03/FFT_03.ino @@ -47,6 +47,7 @@ void setup() { sampling_period_us = round(1000000*(1.0/samplingFrequency)); Serial.begin(115200); + while(!Serial); Serial.println("Ready"); } diff --git a/Examples/FFT_04/FFT_04.ino b/Examples/FFT_04/FFT_04.ino index 5c1fcd5..66d03b5 100644 --- a/Examples/FFT_04/FFT_04.ino +++ b/Examples/FFT_04/FFT_04.ino @@ -54,6 +54,8 @@ double vImag[samples]; void setup() { Serial.begin(115200); + while(!Serial); + Serial.println("Ready"); } void loop() diff --git a/Examples/FFT_05/FFT_05.ino b/Examples/FFT_05/FFT_05.ino index 7abe3b6..599616a 100644 --- a/Examples/FFT_05/FFT_05.ino +++ b/Examples/FFT_05/FFT_05.ino @@ -54,6 +54,7 @@ double vImag[samples]; void setup() { Serial.begin(115200); + while(!Serial); Serial.println("Ready"); } From b19bdc7b6ca400a7569a3f2af07c16e2523aab18 Mon Sep 17 00:00:00 2001 From: Enrique Condes Date: Tue, 6 Oct 2020 15:08:33 -0500 Subject: [PATCH 03/12] Fix to issue 57 --- src/arduinoFFT.cpp | 16 ++++++++++++++++ 1 file changed, 16 insertions(+) diff --git a/src/arduinoFFT.cpp b/src/arduinoFFT.cpp index 580f8ae..790da4f 100644 --- a/src/arduinoFFT.cpp +++ b/src/arduinoFFT.cpp @@ -248,7 +248,11 @@ void arduinoFFT::Windowing(uint8_t windowType, uint8_t dir) weighingFactor = 0.54 * (1.0 - cos(twoPi * ratio)); break; case FFT_WIN_TYP_TRIANGLE: // triangle (Bartlett) + #if defined(ESP8266) || defined(ESP32) + weighingFactor = 1.0 - ((2.0 * fabs(indexMinusOne - (samplesMinusOne / 2.0))) / samplesMinusOne); + #else weighingFactor = 1.0 - ((2.0 * abs(indexMinusOne - (samplesMinusOne / 2.0))) / samplesMinusOne); + #endif break; case FFT_WIN_TYP_NUTTALL: // nuttall weighingFactor = 0.355768 - (0.487396 * (cos(twoPi * ratio))) + (0.144232 * (cos(fourPi * ratio))) - (0.012604 * (cos(sixPi * ratio))); @@ -302,7 +306,11 @@ void arduinoFFT::Windowing(double *vData, uint16_t samples, uint8_t windowType, weighingFactor = 0.54 * (1.0 - cos(twoPi * ratio)); break; case FFT_WIN_TYP_TRIANGLE: // triangle (Bartlett) + #if defined(ESP8266) || defined(ESP32) + weighingFactor = 1.0 - ((2.0 * fabs(indexMinusOne - (samplesMinusOne / 2.0))) / samplesMinusOne); + #else weighingFactor = 1.0 - ((2.0 * abs(indexMinusOne - (samplesMinusOne / 2.0))) / samplesMinusOne); + #endif break; case FFT_WIN_TYP_NUTTALL: // nuttall weighingFactor = 0.355768 - (0.487396 * (cos(twoPi * ratio))) + (0.144232 * (cos(fourPi * ratio))) - (0.012604 * (cos(sixPi * ratio))); @@ -376,7 +384,11 @@ void arduinoFFT::MajorPeak(double *f, double *v) interpolatedX = ((IndexOfMaxY + delta) * this->_samplingFrequency) / (this->_samples); // returned value: interpolated frequency peak apex *f = interpolatedX; + #if defined(ESP8266) || defined(ESP32) + *v = fabs(this->_vReal[IndexOfMaxY - 1] - (2.0 * this->_vReal[IndexOfMaxY]) + this->_vReal[IndexOfMaxY + 1]); + #else *v = abs(this->_vReal[IndexOfMaxY - 1] - (2.0 * this->_vReal[IndexOfMaxY]) + this->_vReal[IndexOfMaxY + 1]); + #endif } double arduinoFFT::MajorPeak(double *vD, uint16_t samples, double samplingFrequency) @@ -424,7 +436,11 @@ void arduinoFFT::MajorPeak(double *vD, uint16_t samples, double samplingFrequenc interpolatedX = ((IndexOfMaxY + delta) * samplingFrequency) / (samples); // returned value: interpolated frequency peak apex *f = interpolatedX; + #if defined(ESP8266) || defined(ESP32) + *v = fabs(vD[IndexOfMaxY - 1] - (2.0 * vD[IndexOfMaxY]) + vD[IndexOfMaxY + 1]); + #else *v = abs(vD[IndexOfMaxY - 1] - (2.0 * vD[IndexOfMaxY]) + vD[IndexOfMaxY + 1]); + #endif } uint8_t arduinoFFT::Exponent(uint16_t value) From 566803e9ca2b6ff50895da2791fcfbf80f8fd23a Mon Sep 17 00:00:00 2001 From: Enrique Condes Date: Tue, 6 Oct 2020 15:18:15 -0500 Subject: [PATCH 04/12] Publish version 1.5.6 --- library.json | 2 +- library.properties | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/library.json b/library.json index 32c7606..19e4073 100644 --- a/library.json +++ b/library.json @@ -20,7 +20,7 @@ "email": "contact@arduinoos.com" } ], - "version": "1.5.5", + "version": "1.5.6", "frameworks": ["arduino","mbed","espidf"], "platforms": "*" } diff --git a/library.properties b/library.properties index cd86480..ea09fa2 100644 --- a/library.properties +++ b/library.properties @@ -1,5 +1,5 @@ name=arduinoFFT -version=1.5.5 +version=1.5.6 author=Enrique Condes maintainer=Enrique Condes sentence=A library for implementing floating point Fast Fourier Transform calculations on Arduino. From 54c383a62f5779cfc40559b00f97729a82a622a8 Mon Sep 17 00:00:00 2001 From: MarScaper Date: Fri, 5 Feb 2021 17:24:36 +0100 Subject: [PATCH 05/12] Improve MajorPeak with parabola I would like to submit you an improvement in frequency estimation by using a parabola estimation using three points... https://github.com/kosme/arduinoFFT/issues/41 --- README.md | 1 + keywords.txt | 3 ++- src/arduinoFFT.cpp | 47 ++++++++++++++++++++++++++++++++++++++++++++++ src/arduinoFFT.h | 3 +++ 4 files changed, 53 insertions(+), 1 deletion(-) diff --git a/README.md b/README.md index 37861ba..986a6c2 100644 --- a/README.md +++ b/README.md @@ -61,6 +61,7 @@ Calcuates the Fast Fourier Transform. Removes the DC component from the sample data. * **MajorPeak**(double *vD, uint16_t samples, double samplingFrequency); * **MajorPeak**(); +* **MajorPeakParabola**(); Looks for and returns the frequency of the biggest spike in the analyzed signal. * **Revision**(void); Returns the library revision. diff --git a/keywords.txt b/keywords.txt index 1daabde..0643678 100644 --- a/keywords.txt +++ b/keywords.txt @@ -3,7 +3,7 @@ ####################################### ####################################### -# Datatypes (KEYWORD1) +# Datatypes (KEYWORD1) ####################################### arduinoFFT KEYWORD1 @@ -19,6 +19,7 @@ Windowing KEYWORD2 Exponent KEYWORD2 Revision KEYWORD2 MajorPeak KEYWORD2 +MajorPeakParabola KEYWORD2 ####################################### # Constants (LITERAL1) diff --git a/src/arduinoFFT.cpp b/src/arduinoFFT.cpp index 790da4f..fb602a9 100644 --- a/src/arduinoFFT.cpp +++ b/src/arduinoFFT.cpp @@ -443,6 +443,53 @@ void arduinoFFT::MajorPeak(double *vD, uint16_t samples, double samplingFrequenc #endif } +double arduinoFFT::MajorPeakParabola() +{ + double maxY = 0; + uint16_t IndexOfMaxY = 0; + //If sampling_frequency = 2 * max_frequency in signal, + //value would be stored at position samples/2 + for (uint16_t i = 1; i < ((this->_samples >> 1) + 1); i++) + { + if ((this->_vReal[i-1] < this->_vReal[i]) && (this->_vReal[i] > this->_vReal[i+1])) + { + if (this->_vReal[i] > maxY) + { + maxY = this->_vReal[i]; + IndexOfMaxY = i; + } + } + } + + double freq = 0; + if( IndexOfMaxY>0 ) + { + // Assume the three points to be on a parabola + double a,b,c; + Parabola(IndexOfMaxY-1, this->_vReal[IndexOfMaxY-1], IndexOfMaxY, this->_vReal[IndexOfMaxY], IndexOfMaxY+1, this->_vReal[IndexOfMaxY+1], &a, &b, &c); + + // Peak is at the middle of the parabola + double x = -b/(2*a); + + // And magnitude is at the extrema of the parabola if you want It... + // double y = a*x*x+b*x+c; + + // Convert to frequency + freq = (x * this->_samplingFrequency) / (this->_samples); + } + + return freq; +} + +void arduinoFFT::Parabola(double x1, double y1, double x2, double y2, double x3, double y3, double *a, double *b, double *c) +{ + double reversed_denom = 1/((x1 - x2) * (x1 - x3) * (x2 - x3)); + + *a = (x3 * (y2 - y1) + x2 * (y1 - y3) + x1 * (y3 - y2)) * reversed_denom; + *b = (x3*x3 * (y1 - y2) + x2*x2 * (y3 - y1) + x1*x1 * (y2 - y3)) * reversed_denom; + *c = (x2 * x3 * (x2 - x3) * y1 + x3 * x1 * (x3 - x1) * y2 + x1 * x2 * (x1 - x2) * y3) *reversed_denom; +} + uint8_t arduinoFFT::Exponent(uint16_t value) { #warning("This method may not be accessible on future revisions.") diff --git a/src/arduinoFFT.h b/src/arduinoFFT.h index 73927dc..c3c1a0b 100644 --- a/src/arduinoFFT.h +++ b/src/arduinoFFT.h @@ -98,6 +98,8 @@ public: void MajorPeak(double *f, double *v); void Windowing(uint8_t windowType, uint8_t dir); + double MajorPeakParabola(); + private: /* Variables */ uint16_t _samples; @@ -107,6 +109,7 @@ private: uint8_t _power; /* Functions */ void Swap(double *x, double *y); + void Parabola(double x1, double y1, double x2, double y2, double x3, double y3, double *a, double *b, double *c); }; #endif From e39f2f12a8f70bc6f0c222d0e12f9511fce83958 Mon Sep 17 00:00:00 2001 From: Ivan Kravets Date: Fri, 19 Nov 2021 12:16:53 +0200 Subject: [PATCH 06/12] Declare header files for PlatformIO --- library.json | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/library.json b/library.json index 19e4073..e1a5953 100644 --- a/library.json +++ b/library.json @@ -22,5 +22,6 @@ ], "version": "1.5.6", "frameworks": ["arduino","mbed","espidf"], - "platforms": "*" + "platforms": "*", + "headers": "arduinoFFT.h" } From dec237ab147d17af2043cdbf4326c30616c7fbec Mon Sep 17 00:00:00 2001 From: Slavey Karadzhov Date: Thu, 9 Dec 2021 12:42:40 +0100 Subject: [PATCH 07/12] Fixed small spelling mistakes. --- Examples/FFT_01/FFT_01.ino | 2 +- Examples/FFT_02/FFT_02.ino | 6 +++--- Examples/FFT_03/FFT_03.ino | 2 +- Examples/FFT_04/FFT_04.ino | 4 ++-- Examples/FFT_05/FFT_05.ino | 2 +- src/arduinoFFT.cpp | 4 ++-- src/arduinoFFT.h | 2 +- 7 files changed, 11 insertions(+), 11 deletions(-) diff --git a/Examples/FFT_01/FFT_01.ino b/Examples/FFT_01/FFT_01.ino index ffeb852..39e3708 100644 --- a/Examples/FFT_01/FFT_01.ino +++ b/Examples/FFT_01/FFT_01.ino @@ -1,6 +1,6 @@ /* - Example of use of the FFT libray + Example of use of the FFT library Copyright (C) 2014 Enrique Condes This program is free software: you can redistribute it and/or modify diff --git a/Examples/FFT_02/FFT_02.ino b/Examples/FFT_02/FFT_02.ino index b621db3..e6262c9 100644 --- a/Examples/FFT_02/FFT_02.ino +++ b/Examples/FFT_02/FFT_02.ino @@ -1,8 +1,8 @@ /* - Example of use of the FFT libray to compute FFT for several signals over a range of frequencies. - The exponent is calculated once before the excecution since it is a constant. - This saves resources during the excecution of the sketch and reduces the compiled size. + Example of use of the FFT library to compute FFT for several signals over a range of frequencies. + The exponent is calculated once before the execution since it is a constant. + This saves resources during the execution of the sketch and reduces the compiled size. The sketch shows the time that the computing is taking. Copyright (C) 2014 Enrique Condes diff --git a/Examples/FFT_03/FFT_03.ino b/Examples/FFT_03/FFT_03.ino index 1e59279..b7c3fba 100644 --- a/Examples/FFT_03/FFT_03.ino +++ b/Examples/FFT_03/FFT_03.ino @@ -1,6 +1,6 @@ /* - Example of use of the FFT libray to compute FFT for a signal sampled through the ADC. + Example of use of the FFT library to compute FFT for a signal sampled through the ADC. Copyright (C) 2018 Enrique Condés and Ragnar Ranøyen Homb This program is free software: you can redistribute it and/or modify diff --git a/Examples/FFT_04/FFT_04.ino b/Examples/FFT_04/FFT_04.ino index 66d03b5..7aa5bb6 100644 --- a/Examples/FFT_04/FFT_04.ino +++ b/Examples/FFT_04/FFT_04.ino @@ -1,6 +1,6 @@ /* - Example of use of the FFT libray + Example of use of the FFT library Copyright (C) 2018 Enrique Condes This program is free software: you can redistribute it and/or modify @@ -26,7 +26,7 @@ of each of the frequencies that compose the signal are calculated. Finally, the frequency spectrum magnitudes are printed. If you use the Arduino IDE serial plotter, you will see a single spike corresponding to the 1000 Hz - frecuency. + frequency. */ #include "arduinoFFT.h" diff --git a/Examples/FFT_05/FFT_05.ino b/Examples/FFT_05/FFT_05.ino index 599616a..7050b10 100644 --- a/Examples/FFT_05/FFT_05.ino +++ b/Examples/FFT_05/FFT_05.ino @@ -1,6 +1,6 @@ /* - Example of use of the FFT libray + Example of use of the FFT library Copyright (C) 2014 Enrique Condes This program is free software: you can redistribute it and/or modify diff --git a/src/arduinoFFT.cpp b/src/arduinoFFT.cpp index fb602a9..f1e6d61 100644 --- a/src/arduinoFFT.cpp +++ b/src/arduinoFFT.cpp @@ -1,6 +1,6 @@ /* - FFT libray + FFT library Copyright (C) 2010 Didier Longueville Copyright (C) 2014 Enrique Condes @@ -230,7 +230,7 @@ void arduinoFFT::DCRemoval(double *vData, uint16_t samples) void arduinoFFT::Windowing(uint8_t windowType, uint8_t dir) {// Weighing factors are computed once before multiple use of FFT -// The weighing function is symetric; half the weighs are recorded +// The weighing function is symmetric; half the weighs are recorded double samplesMinusOne = (double(this->_samples) - 1.0); for (uint16_t i = 0; i < (this->_samples >> 1); i++) { double indexMinusOne = double(i); diff --git a/src/arduinoFFT.h b/src/arduinoFFT.h index c3c1a0b..7883145 100644 --- a/src/arduinoFFT.h +++ b/src/arduinoFFT.h @@ -1,6 +1,6 @@ /* - FFT libray + FFT library Copyright (C) 2010 Didier Longueville Copyright (C) 2014 Enrique Condes From 94453e54acb0f5213264fe812ce4964ef003f59f Mon Sep 17 00:00:00 2001 From: kosme Date: Thu, 9 Mar 2023 22:28:40 -0600 Subject: [PATCH 08/12] Bump to version 1.6 --- changeLog.txt | 3 + library.json | 2 +- library.properties | 2 +- src/arduinoFFT.cpp | 896 +++++++++++++++++++++++---------------------- src/arduinoFFT.h | 176 +++++---- 5 files changed, 579 insertions(+), 500 deletions(-) diff --git a/changeLog.txt b/changeLog.txt index e2cb0fb..33612ef 100644 --- a/changeLog.txt +++ b/changeLog.txt @@ -1,3 +1,6 @@ +03/09/23 v1.6 +Include some functionality from development branch. + 01/27/20 v1.5.5 Lookup table for constants c1 and c2 used during FFT comupting. This increases the FFT computing speed in around 5%. diff --git a/library.json b/library.json index e1a5953..eae9358 100644 --- a/library.json +++ b/library.json @@ -20,7 +20,7 @@ "email": "contact@arduinoos.com" } ], - "version": "1.5.6", + "version": "1.6", "frameworks": ["arduino","mbed","espidf"], "platforms": "*", "headers": "arduinoFFT.h" diff --git a/library.properties b/library.properties index ea09fa2..ee5f5ec 100644 --- a/library.properties +++ b/library.properties @@ -1,5 +1,5 @@ name=arduinoFFT -version=1.5.6 +version=1.6 author=Enrique Condes maintainer=Enrique Condes sentence=A library for implementing floating point Fast Fourier Transform calculations on Arduino. diff --git a/src/arduinoFFT.cpp b/src/arduinoFFT.cpp index f1e6d61..15ab81f 100644 --- a/src/arduinoFFT.cpp +++ b/src/arduinoFFT.cpp @@ -21,489 +21,529 @@ #include "arduinoFFT.h" -arduinoFFT::arduinoFFT(void) -{ // Constructor - #warning("This method is deprecated and may be removed on future revisions.") +arduinoFFT::arduinoFFT(void) { // Constructor +#warning("This method is deprecated and may be removed on future revisions.") } -arduinoFFT::arduinoFFT(double *vReal, double *vImag, uint16_t samples, double samplingFrequency) -{// Constructor - this->_vReal = vReal; - this->_vImag = vImag; - this->_samples = samples; - this->_samplingFrequency = samplingFrequency; - this->_power = Exponent(samples); +arduinoFFT::arduinoFFT(double *vReal, double *vImag, uint16_t samples, + double samplingFrequency) { // Constructor + this->_vReal = vReal; + this->_vImag = vImag; + this->_samples = samples; + this->_samplingFrequency = samplingFrequency; + this->_power = Exponent(samples); } -arduinoFFT::~arduinoFFT(void) -{ -// Destructor +arduinoFFT::~arduinoFFT(void) { + // Destructor } -uint8_t arduinoFFT::Revision(void) -{ - return(FFT_LIB_REV); +uint8_t arduinoFFT::Revision(void) { return (FFT_LIB_REV); } + +void arduinoFFT::Compute(double *vReal, double *vImag, uint16_t samples, + FFTDirection dir) { +#warning("This method is deprecated and may be removed on future revisions.") + Compute(vReal, vImag, samples, Exponent(samples), dir); } -void arduinoFFT::Compute(double *vReal, double *vImag, uint16_t samples, uint8_t dir) -{ - #warning("This method is deprecated and may be removed on future revisions.") - Compute(vReal, vImag, samples, Exponent(samples), dir); -} - -void arduinoFFT::Compute(uint8_t dir) -{// Computes in-place complex-to-complex FFT / - // Reverse bits / - uint16_t j = 0; - for (uint16_t i = 0; i < (this->_samples - 1); i++) { - if (i < j) { - Swap(&this->_vReal[i], &this->_vReal[j]); - if(dir==FFT_REVERSE) - Swap(&this->_vImag[i], &this->_vImag[j]); - } - uint16_t k = (this->_samples >> 1); - while (k <= j) { - j -= k; - k >>= 1; - } - j += k; - } - // Compute the FFT / +void arduinoFFT::Compute(FFTDirection dir) { + // Computes in-place complex-to-complex FFT / + // Reverse bits / + uint16_t j = 0; + for (uint16_t i = 0; i < (this->_samples - 1); i++) { + if (i < j) { + Swap(&this->_vReal[i], &this->_vReal[j]); + if (dir == FFT_REVERSE) + Swap(&this->_vImag[i], &this->_vImag[j]); + } + uint16_t k = (this->_samples >> 1); + while (k <= j) { + j -= k; + k >>= 1; + } + j += k; + } +// Compute the FFT #ifdef __AVR__ - uint8_t index = 0; + uint8_t index = 0; #endif - double c1 = -1.0; - double c2 = 0.0; - uint16_t l2 = 1; - for (uint8_t l = 0; (l < this->_power); l++) { - uint16_t l1 = l2; - l2 <<= 1; - double u1 = 1.0; - double u2 = 0.0; - for (j = 0; j < l1; j++) { - for (uint16_t i = j; i < this->_samples; i += l2) { - uint16_t i1 = i + l1; - double t1 = u1 * this->_vReal[i1] - u2 * this->_vImag[i1]; - double t2 = u1 * this->_vImag[i1] + u2 * this->_vReal[i1]; - this->_vReal[i1] = this->_vReal[i] - t1; - this->_vImag[i1] = this->_vImag[i] - t2; - this->_vReal[i] += t1; - this->_vImag[i] += t2; - } - double z = ((u1 * c1) - (u2 * c2)); - u2 = ((u1 * c2) + (u2 * c1)); - u1 = z; - } + double c1 = -1.0; + double c2 = 0.0; + uint16_t l2 = 1; + for (uint8_t l = 0; (l < this->_power); l++) { + uint16_t l1 = l2; + l2 <<= 1; + double u1 = 1.0; + double u2 = 0.0; + for (j = 0; j < l1; j++) { + for (uint16_t i = j; i < this->_samples; i += l2) { + uint16_t i1 = i + l1; + double t1 = u1 * this->_vReal[i1] - u2 * this->_vImag[i1]; + double t2 = u1 * this->_vImag[i1] + u2 * this->_vReal[i1]; + this->_vReal[i1] = this->_vReal[i] - t1; + this->_vImag[i1] = this->_vImag[i] - t2; + this->_vReal[i] += t1; + this->_vImag[i] += t2; + } + double z = ((u1 * c1) - (u2 * c2)); + u2 = ((u1 * c2) + (u2 * c1)); + u1 = z; + } #ifdef __AVR__ - c2 = pgm_read_float_near(&(_c2[index])); - c1 = pgm_read_float_near(&(_c1[index])); - index++; + c2 = pgm_read_float_near(&(_c2[index])); + c1 = pgm_read_float_near(&(_c1[index])); + index++; #else - c2 = sqrt((1.0 - c1) / 2.0); - c1 = sqrt((1.0 + c1) / 2.0); + c2 = sqrt((1.0 - c1) / 2.0); + c1 = sqrt((1.0 + c1) / 2.0); #endif - if (dir == FFT_FORWARD) { - c2 = -c2; - } - } - // Scaling for reverse transform / - if (dir != FFT_FORWARD) { - for (uint16_t i = 0; i < this->_samples; i++) { - this->_vReal[i] /= this->_samples; - this->_vImag[i] /= this->_samples; - } - } + if (dir == FFT_FORWARD) { + c2 = -c2; + } + } + // Scaling for reverse transform / + if (dir != FFT_FORWARD) { + for (uint16_t i = 0; i < this->_samples; i++) { + this->_vReal[i] /= this->_samples; + this->_vImag[i] /= this->_samples; + } + } } -void arduinoFFT::Compute(double *vReal, double *vImag, uint16_t samples, uint8_t power, uint8_t dir) -{ // Computes in-place complex-to-complex FFT - // Reverse bits - #warning("This method is deprecated and may be removed on future revisions.") - uint16_t j = 0; - for (uint16_t i = 0; i < (samples - 1); i++) { - if (i < j) { - Swap(&vReal[i], &vReal[j]); - if(dir==FFT_REVERSE) - Swap(&vImag[i], &vImag[j]); - } - uint16_t k = (samples >> 1); - while (k <= j) { - j -= k; - k >>= 1; - } - j += k; - } - // Compute the FFT +void arduinoFFT::Compute(double *vReal, double *vImag, uint16_t samples, + uint8_t power, FFTDirection dir) { +// Computes in-place complex-to-complex FFT +// Reverse bits +#warning("This method is deprecated and may be removed on future revisions.") + uint16_t j = 0; + for (uint16_t i = 0; i < (samples - 1); i++) { + if (i < j) { + Swap(&vReal[i], &vReal[j]); + if (dir == FFT_REVERSE) + Swap(&vImag[i], &vImag[j]); + } + uint16_t k = (samples >> 1); + while (k <= j) { + j -= k; + k >>= 1; + } + j += k; + } +// Compute the FFT #ifdef __AVR__ - uint8_t index = 0; + uint8_t index = 0; #endif - double c1 = -1.0; - double c2 = 0.0; - uint16_t l2 = 1; - for (uint8_t l = 0; (l < power); l++) { - uint16_t l1 = l2; - l2 <<= 1; - double u1 = 1.0; - double u2 = 0.0; - for (j = 0; j < l1; j++) { - for (uint16_t i = j; i < samples; i += l2) { - uint16_t i1 = i + l1; - double t1 = u1 * vReal[i1] - u2 * vImag[i1]; - double t2 = u1 * vImag[i1] + u2 * vReal[i1]; - vReal[i1] = vReal[i] - t1; - vImag[i1] = vImag[i] - t2; - vReal[i] += t1; - vImag[i] += t2; - } - double z = ((u1 * c1) - (u2 * c2)); - u2 = ((u1 * c2) + (u2 * c1)); - u1 = z; - } + double c1 = -1.0; + double c2 = 0.0; + uint16_t l2 = 1; + for (uint8_t l = 0; (l < power); l++) { + uint16_t l1 = l2; + l2 <<= 1; + double u1 = 1.0; + double u2 = 0.0; + for (j = 0; j < l1; j++) { + for (uint16_t i = j; i < samples; i += l2) { + uint16_t i1 = i + l1; + double t1 = u1 * vReal[i1] - u2 * vImag[i1]; + double t2 = u1 * vImag[i1] + u2 * vReal[i1]; + vReal[i1] = vReal[i] - t1; + vImag[i1] = vImag[i] - t2; + vReal[i] += t1; + vImag[i] += t2; + } + double z = ((u1 * c1) - (u2 * c2)); + u2 = ((u1 * c2) + (u2 * c1)); + u1 = z; + } #ifdef __AVR__ - c2 = pgm_read_float_near(&(_c2[index])); - c1 = pgm_read_float_near(&(_c1[index])); - index++; + c2 = pgm_read_float_near(&(_c2[index])); + c1 = pgm_read_float_near(&(_c1[index])); + index++; #else - c2 = sqrt((1.0 - c1) / 2.0); - c1 = sqrt((1.0 + c1) / 2.0); + c2 = sqrt((1.0 - c1) / 2.0); + c1 = sqrt((1.0 + c1) / 2.0); #endif - if (dir == FFT_FORWARD) { - c2 = -c2; - } - } - // Scaling for reverse transform - if (dir != FFT_FORWARD) { - for (uint16_t i = 0; i < samples; i++) { - vReal[i] /= samples; - vImag[i] /= samples; - } - } + if (dir == FFT_FORWARD) { + c2 = -c2; + } + } + // Scaling for reverse transform + if (dir != FFT_FORWARD) { + for (uint16_t i = 0; i < samples; i++) { + vReal[i] /= samples; + vImag[i] /= samples; + } + } } -void arduinoFFT::ComplexToMagnitude() -{ // vM is half the size of vReal and vImag - for (uint16_t i = 0; i < this->_samples; i++) { - this->_vReal[i] = sqrt(sq(this->_vReal[i]) + sq(this->_vImag[i])); - } +void arduinoFFT::ComplexToMagnitude() { + // vM is half the size of vReal and vImag + for (uint16_t i = 0; i < this->_samples; i++) { + this->_vReal[i] = sqrt(sq(this->_vReal[i]) + sq(this->_vImag[i])); + } } -void arduinoFFT::ComplexToMagnitude(double *vReal, double *vImag, uint16_t samples) -{ // vM is half the size of vReal and vImag - #warning("This method is deprecated and may be removed on future revisions.") - for (uint16_t i = 0; i < samples; i++) { - vReal[i] = sqrt(sq(vReal[i]) + sq(vImag[i])); - } +void arduinoFFT::ComplexToMagnitude(double *vReal, double *vImag, + uint16_t samples) { +// vM is half the size of vReal and vImag +#warning("This method is deprecated and may be removed on future revisions.") + for (uint16_t i = 0; i < samples; i++) { + vReal[i] = sqrt(sq(vReal[i]) + sq(vImag[i])); + } } -void arduinoFFT::DCRemoval() -{ - // calculate the mean of vData - double mean = 0; - for (uint16_t i = 0; i < this->_samples; i++) - { - mean += this->_vReal[i]; - } - mean /= this->_samples; - // Subtract the mean from vData - for (uint16_t i = 0; i < this->_samples; i++) - { - this->_vReal[i] -= mean; - } +void arduinoFFT::DCRemoval() { + // calculate the mean of vData + double mean = 0; + for (uint16_t i = 0; i < this->_samples; i++) { + mean += this->_vReal[i]; + } + mean /= this->_samples; + // Subtract the mean from vData + for (uint16_t i = 0; i < this->_samples; i++) { + this->_vReal[i] -= mean; + } } -void arduinoFFT::DCRemoval(double *vData, uint16_t samples) -{ - // calculate the mean of vData - #warning("This method is deprecated and may be removed on future revisions.") - double mean = 0; - for (uint16_t i = 0; i < samples; i++) - { - mean += vData[i]; - } - mean /= samples; - // Subtract the mean from vData - for (uint16_t i = 0; i < samples; i++) - { - vData[i] -= mean; - } +void arduinoFFT::DCRemoval(double *vData, uint16_t samples) { +// calculate the mean of vData +#warning("This method is deprecated and may be removed on future revisions.") + double mean = 0; + for (uint16_t i = 0; i < samples; i++) { + mean += vData[i]; + } + mean /= samples; + // Subtract the mean from vData + for (uint16_t i = 0; i < samples; i++) { + vData[i] -= mean; + } } -void arduinoFFT::Windowing(uint8_t windowType, uint8_t dir) -{// Weighing factors are computed once before multiple use of FFT -// The weighing function is symmetric; half the weighs are recorded - double samplesMinusOne = (double(this->_samples) - 1.0); - for (uint16_t i = 0; i < (this->_samples >> 1); i++) { - double indexMinusOne = double(i); - double ratio = (indexMinusOne / samplesMinusOne); - double weighingFactor = 1.0; - // Compute and record weighting factor - switch (windowType) { - case FFT_WIN_TYP_RECTANGLE: // rectangle (box car) - weighingFactor = 1.0; - break; - case FFT_WIN_TYP_HAMMING: // hamming - weighingFactor = 0.54 - (0.46 * cos(twoPi * ratio)); - break; - case FFT_WIN_TYP_HANN: // hann - weighingFactor = 0.54 * (1.0 - cos(twoPi * ratio)); - break; - case FFT_WIN_TYP_TRIANGLE: // triangle (Bartlett) - #if defined(ESP8266) || defined(ESP32) - weighingFactor = 1.0 - ((2.0 * fabs(indexMinusOne - (samplesMinusOne / 2.0))) / samplesMinusOne); - #else - weighingFactor = 1.0 - ((2.0 * abs(indexMinusOne - (samplesMinusOne / 2.0))) / samplesMinusOne); - #endif - break; - case FFT_WIN_TYP_NUTTALL: // nuttall - weighingFactor = 0.355768 - (0.487396 * (cos(twoPi * ratio))) + (0.144232 * (cos(fourPi * ratio))) - (0.012604 * (cos(sixPi * ratio))); - break; - case FFT_WIN_TYP_BLACKMAN: // blackman - weighingFactor = 0.42323 - (0.49755 * (cos(twoPi * ratio))) + (0.07922 * (cos(fourPi * ratio))); - break; - case FFT_WIN_TYP_BLACKMAN_NUTTALL: // blackman nuttall - weighingFactor = 0.3635819 - (0.4891775 * (cos(twoPi * ratio))) + (0.1365995 * (cos(fourPi * ratio))) - (0.0106411 * (cos(sixPi * ratio))); - break; - case FFT_WIN_TYP_BLACKMAN_HARRIS: // blackman harris - weighingFactor = 0.35875 - (0.48829 * (cos(twoPi * ratio))) + (0.14128 * (cos(fourPi * ratio))) - (0.01168 * (cos(sixPi * ratio))); - break; - case FFT_WIN_TYP_FLT_TOP: // flat top - weighingFactor = 0.2810639 - (0.5208972 * cos(twoPi * ratio)) + (0.1980399 * cos(fourPi * ratio)); - break; - case FFT_WIN_TYP_WELCH: // welch - weighingFactor = 1.0 - sq((indexMinusOne - samplesMinusOne / 2.0) / (samplesMinusOne / 2.0)); - break; - } - if (dir == FFT_FORWARD) { - this->_vReal[i] *= weighingFactor; - this->_vReal[this->_samples - (i + 1)] *= weighingFactor; - } - else { - this->_vReal[i] /= weighingFactor; - this->_vReal[this->_samples - (i + 1)] /= weighingFactor; - } - } +void arduinoFFT::Windowing(FFTWindow windowType, FFTDirection dir) { + // Weighing factors are computed once before multiple use of FFT + // The weighing function is symmetric; half the weighs are recorded + double samplesMinusOne = (double(this->_samples) - 1.0); + for (uint16_t i = 0; i < (this->_samples >> 1); i++) { + double indexMinusOne = double(i); + double ratio = (indexMinusOne / samplesMinusOne); + double weighingFactor = 1.0; + // Compute and record weighting factor + switch (windowType) { + case FFT_WIN_TYP_RECTANGLE: // rectangle (box car) + weighingFactor = 1.0; + break; + case FFT_WIN_TYP_HAMMING: // hamming + weighingFactor = 0.54 - (0.46 * cos(twoPi * ratio)); + break; + case FFT_WIN_TYP_HANN: // hann + weighingFactor = 0.54 * (1.0 - cos(twoPi * ratio)); + break; + case FFT_WIN_TYP_TRIANGLE: // triangle (Bartlett) +#if defined(ESP8266) || defined(ESP32) + weighingFactor = + 1.0 - ((2.0 * fabs(indexMinusOne - (samplesMinusOne / 2.0))) / + samplesMinusOne); +#else + weighingFactor = + 1.0 - ((2.0 * abs(indexMinusOne - (samplesMinusOne / 2.0))) / + samplesMinusOne); +#endif + break; + case FFT_WIN_TYP_NUTTALL: // nuttall + weighingFactor = 0.355768 - (0.487396 * (cos(twoPi * ratio))) + + (0.144232 * (cos(fourPi * ratio))) - + (0.012604 * (cos(sixPi * ratio))); + break; + case FFT_WIN_TYP_BLACKMAN: // blackman + weighingFactor = 0.42323 - (0.49755 * (cos(twoPi * ratio))) + + (0.07922 * (cos(fourPi * ratio))); + break; + case FFT_WIN_TYP_BLACKMAN_NUTTALL: // blackman nuttall + weighingFactor = 0.3635819 - (0.4891775 * (cos(twoPi * ratio))) + + (0.1365995 * (cos(fourPi * ratio))) - + (0.0106411 * (cos(sixPi * ratio))); + break; + case FFT_WIN_TYP_BLACKMAN_HARRIS: // blackman harris + weighingFactor = 0.35875 - (0.48829 * (cos(twoPi * ratio))) + + (0.14128 * (cos(fourPi * ratio))) - + (0.01168 * (cos(sixPi * ratio))); + break; + case FFT_WIN_TYP_FLT_TOP: // flat top + weighingFactor = 0.2810639 - (0.5208972 * cos(twoPi * ratio)) + + (0.1980399 * cos(fourPi * ratio)); + break; + case FFT_WIN_TYP_WELCH: // welch + weighingFactor = 1.0 - sq((indexMinusOne - samplesMinusOne / 2.0) / + (samplesMinusOne / 2.0)); + break; + } + if (dir == FFT_FORWARD) { + this->_vReal[i] *= weighingFactor; + this->_vReal[this->_samples - (i + 1)] *= weighingFactor; + } else { + this->_vReal[i] /= weighingFactor; + this->_vReal[this->_samples - (i + 1)] /= weighingFactor; + } + } } - -void arduinoFFT::Windowing(double *vData, uint16_t samples, uint8_t windowType, uint8_t dir) -{// Weighing factors are computed once before multiple use of FFT +void arduinoFFT::Windowing(double *vData, uint16_t samples, + FFTWindow windowType, FFTDirection dir) { +// Weighing factors are computed once before multiple use of FFT // The weighing function is symetric; half the weighs are recorded - #warning("This method is deprecated and may be removed on future revisions.") - double samplesMinusOne = (double(samples) - 1.0); - for (uint16_t i = 0; i < (samples >> 1); i++) { - double indexMinusOne = double(i); - double ratio = (indexMinusOne / samplesMinusOne); - double weighingFactor = 1.0; - // Compute and record weighting factor - switch (windowType) { - case FFT_WIN_TYP_RECTANGLE: // rectangle (box car) - weighingFactor = 1.0; - break; - case FFT_WIN_TYP_HAMMING: // hamming - weighingFactor = 0.54 - (0.46 * cos(twoPi * ratio)); - break; - case FFT_WIN_TYP_HANN: // hann - weighingFactor = 0.54 * (1.0 - cos(twoPi * ratio)); - break; - case FFT_WIN_TYP_TRIANGLE: // triangle (Bartlett) - #if defined(ESP8266) || defined(ESP32) - weighingFactor = 1.0 - ((2.0 * fabs(indexMinusOne - (samplesMinusOne / 2.0))) / samplesMinusOne); - #else - weighingFactor = 1.0 - ((2.0 * abs(indexMinusOne - (samplesMinusOne / 2.0))) / samplesMinusOne); - #endif - break; - case FFT_WIN_TYP_NUTTALL: // nuttall - weighingFactor = 0.355768 - (0.487396 * (cos(twoPi * ratio))) + (0.144232 * (cos(fourPi * ratio))) - (0.012604 * (cos(sixPi * ratio))); - break; - case FFT_WIN_TYP_BLACKMAN: // blackman - weighingFactor = 0.42323 - (0.49755 * (cos(twoPi * ratio))) + (0.07922 * (cos(fourPi * ratio))); - break; - case FFT_WIN_TYP_BLACKMAN_NUTTALL: // blackman nuttall - weighingFactor = 0.3635819 - (0.4891775 * (cos(twoPi * ratio))) + (0.1365995 * (cos(fourPi * ratio))) - (0.0106411 * (cos(sixPi * ratio))); - break; - case FFT_WIN_TYP_BLACKMAN_HARRIS: // blackman harris - weighingFactor = 0.35875 - (0.48829 * (cos(twoPi * ratio))) + (0.14128 * (cos(fourPi * ratio))) - (0.01168 * (cos(sixPi * ratio))); - break; - case FFT_WIN_TYP_FLT_TOP: // flat top - weighingFactor = 0.2810639 - (0.5208972 * cos(twoPi * ratio)) + (0.1980399 * cos(fourPi * ratio)); - break; - case FFT_WIN_TYP_WELCH: // welch - weighingFactor = 1.0 - sq((indexMinusOne - samplesMinusOne / 2.0) / (samplesMinusOne / 2.0)); - break; - } - if (dir == FFT_FORWARD) { - vData[i] *= weighingFactor; - vData[samples - (i + 1)] *= weighingFactor; - } - else { - vData[i] /= weighingFactor; - vData[samples - (i + 1)] /= weighingFactor; - } - } +#warning("This method is deprecated and may be removed on future revisions.") + double samplesMinusOne = (double(samples) - 1.0); + for (uint16_t i = 0; i < (samples >> 1); i++) { + double indexMinusOne = double(i); + double ratio = (indexMinusOne / samplesMinusOne); + double weighingFactor = 1.0; + // Compute and record weighting factor + switch (windowType) { + case FFT_WIN_TYP_RECTANGLE: // rectangle (box car) + weighingFactor = 1.0; + break; + case FFT_WIN_TYP_HAMMING: // hamming + weighingFactor = 0.54 - (0.46 * cos(twoPi * ratio)); + break; + case FFT_WIN_TYP_HANN: // hann + weighingFactor = 0.54 * (1.0 - cos(twoPi * ratio)); + break; + case FFT_WIN_TYP_TRIANGLE: // triangle (Bartlett) +#if defined(ESP8266) || defined(ESP32) + weighingFactor = + 1.0 - ((2.0 * fabs(indexMinusOne - (samplesMinusOne / 2.0))) / + samplesMinusOne); +#else + weighingFactor = + 1.0 - ((2.0 * abs(indexMinusOne - (samplesMinusOne / 2.0))) / + samplesMinusOne); +#endif + break; + case FFT_WIN_TYP_NUTTALL: // nuttall + weighingFactor = 0.355768 - (0.487396 * (cos(twoPi * ratio))) + + (0.144232 * (cos(fourPi * ratio))) - + (0.012604 * (cos(sixPi * ratio))); + break; + case FFT_WIN_TYP_BLACKMAN: // blackman + weighingFactor = 0.42323 - (0.49755 * (cos(twoPi * ratio))) + + (0.07922 * (cos(fourPi * ratio))); + break; + case FFT_WIN_TYP_BLACKMAN_NUTTALL: // blackman nuttall + weighingFactor = 0.3635819 - (0.4891775 * (cos(twoPi * ratio))) + + (0.1365995 * (cos(fourPi * ratio))) - + (0.0106411 * (cos(sixPi * ratio))); + break; + case FFT_WIN_TYP_BLACKMAN_HARRIS: // blackman harris + weighingFactor = 0.35875 - (0.48829 * (cos(twoPi * ratio))) + + (0.14128 * (cos(fourPi * ratio))) - + (0.01168 * (cos(sixPi * ratio))); + break; + case FFT_WIN_TYP_FLT_TOP: // flat top + weighingFactor = 0.2810639 - (0.5208972 * cos(twoPi * ratio)) + + (0.1980399 * cos(fourPi * ratio)); + break; + case FFT_WIN_TYP_WELCH: // welch + weighingFactor = 1.0 - sq((indexMinusOne - samplesMinusOne / 2.0) / + (samplesMinusOne / 2.0)); + break; + } + if (dir == FFT_FORWARD) { + vData[i] *= weighingFactor; + vData[samples - (i + 1)] *= weighingFactor; + } else { + vData[i] /= weighingFactor; + vData[samples - (i + 1)] /= weighingFactor; + } + } } -double arduinoFFT::MajorPeak() -{ - double maxY = 0; - uint16_t IndexOfMaxY = 0; - //If sampling_frequency = 2 * max_frequency in signal, - //value would be stored at position samples/2 - for (uint16_t i = 1; i < ((this->_samples >> 1) + 1); i++) { - if ((this->_vReal[i-1] < this->_vReal[i]) && (this->_vReal[i] > this->_vReal[i+1])) { - if (this->_vReal[i] > maxY) { - maxY = this->_vReal[i]; - IndexOfMaxY = i; - } - } - } - double delta = 0.5 * ((this->_vReal[IndexOfMaxY-1] - this->_vReal[IndexOfMaxY+1]) / (this->_vReal[IndexOfMaxY-1] - (2.0 * this->_vReal[IndexOfMaxY]) + this->_vReal[IndexOfMaxY+1])); - double interpolatedX = ((IndexOfMaxY + delta) * this->_samplingFrequency) / (this->_samples-1); - if(IndexOfMaxY==(this->_samples >> 1)) //To improve calculation on edge values - interpolatedX = ((IndexOfMaxY + delta) * this->_samplingFrequency) / (this->_samples); - // returned value: interpolated frequency peak apex - return(interpolatedX); +double arduinoFFT::MajorPeak() { + double maxY = 0; + uint16_t IndexOfMaxY = 0; + // If sampling_frequency = 2 * max_frequency in signal, + // value would be stored at position samples/2 + for (uint16_t i = 1; i < ((this->_samples >> 1) + 1); i++) { + if ((this->_vReal[i - 1] < this->_vReal[i]) && + (this->_vReal[i] > this->_vReal[i + 1])) { + if (this->_vReal[i] > maxY) { + maxY = this->_vReal[i]; + IndexOfMaxY = i; + } + } + } + double delta = + 0.5 * + ((this->_vReal[IndexOfMaxY - 1] - this->_vReal[IndexOfMaxY + 1]) / + (this->_vReal[IndexOfMaxY - 1] - (2.0 * this->_vReal[IndexOfMaxY]) + + this->_vReal[IndexOfMaxY + 1])); + double interpolatedX = + ((IndexOfMaxY + delta) * this->_samplingFrequency) / (this->_samples - 1); + if (IndexOfMaxY == + (this->_samples >> 1)) // To improve calculation on edge values + interpolatedX = + ((IndexOfMaxY + delta) * this->_samplingFrequency) / (this->_samples); + // returned value: interpolated frequency peak apex + return (interpolatedX); } -void arduinoFFT::MajorPeak(double *f, double *v) -{ - double maxY = 0; - uint16_t IndexOfMaxY = 0; - //If sampling_frequency = 2 * max_frequency in signal, - //value would be stored at position samples/2 - for (uint16_t i = 1; i < ((this->_samples >> 1) + 1); i++) { - if ((this->_vReal[i - 1] < this->_vReal[i]) && (this->_vReal[i] > this->_vReal[i + 1])) { - if (this->_vReal[i] > maxY) { - maxY = this->_vReal[i]; - IndexOfMaxY = i; - } - } - } - double delta = 0.5 * ((this->_vReal[IndexOfMaxY - 1] - this->_vReal[IndexOfMaxY + 1]) / (this->_vReal[IndexOfMaxY - 1] - (2.0 * this->_vReal[IndexOfMaxY]) + this->_vReal[IndexOfMaxY + 1])); - double interpolatedX = ((IndexOfMaxY + delta) * this->_samplingFrequency) / (this->_samples - 1); - if (IndexOfMaxY == (this->_samples >> 1)) //To improve calculation on edge values - interpolatedX = ((IndexOfMaxY + delta) * this->_samplingFrequency) / (this->_samples); - // returned value: interpolated frequency peak apex - *f = interpolatedX; - #if defined(ESP8266) || defined(ESP32) - *v = fabs(this->_vReal[IndexOfMaxY - 1] - (2.0 * this->_vReal[IndexOfMaxY]) + this->_vReal[IndexOfMaxY + 1]); - #else - *v = abs(this->_vReal[IndexOfMaxY - 1] - (2.0 * this->_vReal[IndexOfMaxY]) + this->_vReal[IndexOfMaxY + 1]); - #endif +void arduinoFFT::MajorPeak(double *f, double *v) { + double maxY = 0; + uint16_t IndexOfMaxY = 0; + // If sampling_frequency = 2 * max_frequency in signal, + // value would be stored at position samples/2 + for (uint16_t i = 1; i < ((this->_samples >> 1) + 1); i++) { + if ((this->_vReal[i - 1] < this->_vReal[i]) && + (this->_vReal[i] > this->_vReal[i + 1])) { + if (this->_vReal[i] > maxY) { + maxY = this->_vReal[i]; + IndexOfMaxY = i; + } + } + } + double delta = + 0.5 * + ((this->_vReal[IndexOfMaxY - 1] - this->_vReal[IndexOfMaxY + 1]) / + (this->_vReal[IndexOfMaxY - 1] - (2.0 * this->_vReal[IndexOfMaxY]) + + this->_vReal[IndexOfMaxY + 1])); + double interpolatedX = + ((IndexOfMaxY + delta) * this->_samplingFrequency) / (this->_samples - 1); + if (IndexOfMaxY == + (this->_samples >> 1)) // To improve calculation on edge values + interpolatedX = + ((IndexOfMaxY + delta) * this->_samplingFrequency) / (this->_samples); + // returned value: interpolated frequency peak apex + *f = interpolatedX; +#if defined(ESP8266) || defined(ESP32) + *v = fabs(this->_vReal[IndexOfMaxY - 1] - (2.0 * this->_vReal[IndexOfMaxY]) + + this->_vReal[IndexOfMaxY + 1]); +#else + *v = abs(this->_vReal[IndexOfMaxY - 1] - (2.0 * this->_vReal[IndexOfMaxY]) + + this->_vReal[IndexOfMaxY + 1]); +#endif } -double arduinoFFT::MajorPeak(double *vD, uint16_t samples, double samplingFrequency) -{ - #warning("This method is deprecated and may be removed on future revisions.") - double maxY = 0; - uint16_t IndexOfMaxY = 0; - //If sampling_frequency = 2 * max_frequency in signal, - //value would be stored at position samples/2 - for (uint16_t i = 1; i < ((samples >> 1) + 1); i++) { - if ((vD[i-1] < vD[i]) && (vD[i] > vD[i+1])) { - if (vD[i] > maxY) { - maxY = vD[i]; - IndexOfMaxY = i; - } - } - } - double delta = 0.5 * ((vD[IndexOfMaxY-1] - vD[IndexOfMaxY+1]) / (vD[IndexOfMaxY-1] - (2.0 * vD[IndexOfMaxY]) + vD[IndexOfMaxY+1])); - double interpolatedX = ((IndexOfMaxY + delta) * samplingFrequency) / (samples-1); - if(IndexOfMaxY==(samples >> 1)) //To improve calculation on edge values - interpolatedX = ((IndexOfMaxY + delta) * samplingFrequency) / (samples); - // returned value: interpolated frequency peak apex - return(interpolatedX); +double arduinoFFT::MajorPeak(double *vD, uint16_t samples, + double samplingFrequency) { +#warning("This method is deprecated and may be removed on future revisions.") + double maxY = 0; + uint16_t IndexOfMaxY = 0; + // If sampling_frequency = 2 * max_frequency in signal, + // value would be stored at position samples/2 + for (uint16_t i = 1; i < ((samples >> 1) + 1); i++) { + if ((vD[i - 1] < vD[i]) && (vD[i] > vD[i + 1])) { + if (vD[i] > maxY) { + maxY = vD[i]; + IndexOfMaxY = i; + } + } + } + double delta = + 0.5 * + ((vD[IndexOfMaxY - 1] - vD[IndexOfMaxY + 1]) / + (vD[IndexOfMaxY - 1] - (2.0 * vD[IndexOfMaxY]) + vD[IndexOfMaxY + 1])); + double interpolatedX = + ((IndexOfMaxY + delta) * samplingFrequency) / (samples - 1); + if (IndexOfMaxY == (samples >> 1)) // To improve calculation on edge values + interpolatedX = ((IndexOfMaxY + delta) * samplingFrequency) / (samples); + // returned value: interpolated frequency peak apex + return (interpolatedX); } -void arduinoFFT::MajorPeak(double *vD, uint16_t samples, double samplingFrequency, double *f, double *v) -{ - #warning("This method is deprecated and may be removed on future revisions.") - double maxY = 0; - uint16_t IndexOfMaxY = 0; - //If sampling_frequency = 2 * max_frequency in signal, - //value would be stored at position samples/2 - for (uint16_t i = 1; i < ((samples >> 1) + 1); i++) { - if ((vD[i - 1] < vD[i]) && (vD[i] > vD[i + 1])) { - if (vD[i] > maxY) { - maxY = vD[i]; - IndexOfMaxY = i; - } - } - } - double delta = 0.5 * ((vD[IndexOfMaxY - 1] - vD[IndexOfMaxY + 1]) / (vD[IndexOfMaxY - 1] - (2.0 * vD[IndexOfMaxY]) + vD[IndexOfMaxY + 1])); - double interpolatedX = ((IndexOfMaxY + delta) * samplingFrequency) / (samples - 1); - //double popo = - if (IndexOfMaxY == (samples >> 1)) //To improve calculation on edge values - interpolatedX = ((IndexOfMaxY + delta) * samplingFrequency) / (samples); - // returned value: interpolated frequency peak apex - *f = interpolatedX; - #if defined(ESP8266) || defined(ESP32) - *v = fabs(vD[IndexOfMaxY - 1] - (2.0 * vD[IndexOfMaxY]) + vD[IndexOfMaxY + 1]); - #else - *v = abs(vD[IndexOfMaxY - 1] - (2.0 * vD[IndexOfMaxY]) + vD[IndexOfMaxY + 1]); - #endif +void arduinoFFT::MajorPeak(double *vD, uint16_t samples, + double samplingFrequency, double *f, double *v) { +#warning("This method is deprecated and may be removed on future revisions.") + double maxY = 0; + uint16_t IndexOfMaxY = 0; + // If sampling_frequency = 2 * max_frequency in signal, + // value would be stored at position samples/2 + for (uint16_t i = 1; i < ((samples >> 1) + 1); i++) { + if ((vD[i - 1] < vD[i]) && (vD[i] > vD[i + 1])) { + if (vD[i] > maxY) { + maxY = vD[i]; + IndexOfMaxY = i; + } + } + } + double delta = + 0.5 * + ((vD[IndexOfMaxY - 1] - vD[IndexOfMaxY + 1]) / + (vD[IndexOfMaxY - 1] - (2.0 * vD[IndexOfMaxY]) + vD[IndexOfMaxY + 1])); + double interpolatedX = + ((IndexOfMaxY + delta) * samplingFrequency) / (samples - 1); + // double popo = + if (IndexOfMaxY == (samples >> 1)) // To improve calculation on edge values + interpolatedX = ((IndexOfMaxY + delta) * samplingFrequency) / (samples); + // returned value: interpolated frequency peak apex + *f = interpolatedX; +#if defined(ESP8266) || defined(ESP32) + *v = + fabs(vD[IndexOfMaxY - 1] - (2.0 * vD[IndexOfMaxY]) + vD[IndexOfMaxY + 1]); +#else + *v = abs(vD[IndexOfMaxY - 1] - (2.0 * vD[IndexOfMaxY]) + vD[IndexOfMaxY + 1]); +#endif } -double arduinoFFT::MajorPeakParabola() -{ - double maxY = 0; - uint16_t IndexOfMaxY = 0; - //If sampling_frequency = 2 * max_frequency in signal, - //value would be stored at position samples/2 - for (uint16_t i = 1; i < ((this->_samples >> 1) + 1); i++) - { - if ((this->_vReal[i-1] < this->_vReal[i]) && (this->_vReal[i] > this->_vReal[i+1])) - { - if (this->_vReal[i] > maxY) - { - maxY = this->_vReal[i]; - IndexOfMaxY = i; - } - } - } +double arduinoFFT::MajorPeakParabola() { + double maxY = 0; + uint16_t IndexOfMaxY = 0; + // If sampling_frequency = 2 * max_frequency in signal, + // value would be stored at position samples/2 + for (uint16_t i = 1; i < ((this->_samples >> 1) + 1); i++) { + if ((this->_vReal[i - 1] < this->_vReal[i]) && + (this->_vReal[i] > this->_vReal[i + 1])) { + if (this->_vReal[i] > maxY) { + maxY = this->_vReal[i]; + IndexOfMaxY = i; + } + } + } - double freq = 0; - if( IndexOfMaxY>0 ) - { - // Assume the three points to be on a parabola - double a,b,c; - Parabola(IndexOfMaxY-1, this->_vReal[IndexOfMaxY-1], IndexOfMaxY, this->_vReal[IndexOfMaxY], IndexOfMaxY+1, this->_vReal[IndexOfMaxY+1], &a, &b, &c); + double freq = 0; + if (IndexOfMaxY > 0) { + // Assume the three points to be on a parabola + double a, b, c; + Parabola(IndexOfMaxY - 1, this->_vReal[IndexOfMaxY - 1], IndexOfMaxY, + this->_vReal[IndexOfMaxY], IndexOfMaxY + 1, + this->_vReal[IndexOfMaxY + 1], &a, &b, &c); - // Peak is at the middle of the parabola - double x = -b/(2*a); + // Peak is at the middle of the parabola + double x = -b / (2 * a); - // And magnitude is at the extrema of the parabola if you want It... - // double y = a*x*x+b*x+c; + // And magnitude is at the extrema of the parabola if you want It... + // double y = a*x*x+b*x+c; - // Convert to frequency - freq = (x * this->_samplingFrequency) / (this->_samples); - } + // Convert to frequency + freq = (x * this->_samplingFrequency) / (this->_samples); + } - return freq; + return freq; } -void arduinoFFT::Parabola(double x1, double y1, double x2, double y2, double x3, double y3, double *a, double *b, double *c) -{ - double reversed_denom = 1/((x1 - x2) * (x1 - x3) * (x2 - x3)); +void arduinoFFT::Parabola(double x1, double y1, double x2, double y2, double x3, + double y3, double *a, double *b, double *c) { + double reversed_denom = 1 / ((x1 - x2) * (x1 - x3) * (x2 - x3)); - *a = (x3 * (y2 - y1) + x2 * (y1 - y3) + x1 * (y3 - y2)) * reversed_denom; - *b = (x3*x3 * (y1 - y2) + x2*x2 * (y3 - y1) + x1*x1 * (y2 - y3)) * reversed_denom; - *c = (x2 * x3 * (x2 - x3) * y1 + x3 * x1 * (x3 - x1) * y2 + x1 * x2 * (x1 - x2) * y3) *reversed_denom; + *a = (x3 * (y2 - y1) + x2 * (y1 - y3) + x1 * (y3 - y2)) * reversed_denom; + *b = (x3 * x3 * (y1 - y2) + x2 * x2 * (y3 - y1) + x1 * x1 * (y2 - y3)) * + reversed_denom; + *c = (x2 * x3 * (x2 - x3) * y1 + x3 * x1 * (x3 - x1) * y2 + + x1 * x2 * (x1 - x2) * y3) * + reversed_denom; } -uint8_t arduinoFFT::Exponent(uint16_t value) -{ - #warning("This method may not be accessible on future revisions.") - // Calculates the base 2 logarithm of a value - uint8_t result = 0; - while (((value >> result) & 1) != 1) result++; - return(result); +uint8_t arduinoFFT::Exponent(uint16_t value) { +#warning("This method may not be accessible on future revisions.") + // Calculates the base 2 logarithm of a value + uint8_t result = 0; + while (((value >> result) & 1) != 1) + result++; + return (result); } // Private functions -void arduinoFFT::Swap(double *x, double *y) -{ - double temp = *x; - *x = *y; - *y = temp; +void arduinoFFT::Swap(double *x, double *y) { + double temp = *x; + *x = *y; + *y = temp; } diff --git a/src/arduinoFFT.h b/src/arduinoFFT.h index 7883145..4878d01 100644 --- a/src/arduinoFFT.h +++ b/src/arduinoFFT.h @@ -19,97 +19,133 @@ */ -#ifndef arduinoFFT_h /* Prevent loading library twice */ -#define arduinoFFT_h +#ifndef ArduinoFFT_h /* Prevent loading library twice */ +#define ArduinoFFT_h #ifdef ARDUINO - #if ARDUINO >= 100 - #include "Arduino.h" - #else - #include "WProgram.h" /* This is where the standard Arduino code lies */ - #endif +#if ARDUINO >= 100 +#include "Arduino.h" #else - #include - #include - #ifdef __AVR__ - #include - #include - #endif - #include - #include "defs.h" - #include "types.h" +#include "WProgram.h" /* This is where the standard Arduino code lies */ +#endif +#else +#include +#include + +#ifdef __AVR__ +#include +#include +#endif +#include "defs.h" +#include "types.h" +#include + #endif -#define FFT_LIB_REV 0x14 +// Define this to use a low-precision square root approximation instead of the +// regular sqrt() call +// This might only work for specific use cases, but is significantly faster. +// Only works for ArduinoFFT. +// #define FFT_SQRT_APPROXIMATION + +#ifdef FFT_SQRT_APPROXIMATION +#include +#else +#define sqrt_internal sqrt +#endif + +enum class FFTDirection { Reverse, Forward }; + +enum class FFTWindow { + Rectangle, // rectangle (Box car) + Hamming, // hamming + Hann, // hann + Triangle, // triangle (Bartlett) + Nuttall, // nuttall + Blackman, // blackman + Blackman_Nuttall, // blackman nuttall + Blackman_Harris, // blackman harris + Flat_top, // flat top + Welch // welch +}; +#define FFT_LIB_REV 0x15 /* Custom constants */ -#define FFT_FORWARD 0x01 -#define FFT_REVERSE 0x00 +#define FFT_FORWARD FFTDirection::Forward +#define FFT_REVERSE FFTDirection::Reverse /* Windowing type */ -#define FFT_WIN_TYP_RECTANGLE 0x00 /* rectangle (Box car) */ -#define FFT_WIN_TYP_HAMMING 0x01 /* hamming */ -#define FFT_WIN_TYP_HANN 0x02 /* hann */ -#define FFT_WIN_TYP_TRIANGLE 0x03 /* triangle (Bartlett) */ -#define FFT_WIN_TYP_NUTTALL 0x04 /* nuttall */ -#define FFT_WIN_TYP_BLACKMAN 0x05 /* blackman */ -#define FFT_WIN_TYP_BLACKMAN_NUTTALL 0x06 /* blackman nuttall */ -#define FFT_WIN_TYP_BLACKMAN_HARRIS 0x07 /* blackman harris*/ -#define FFT_WIN_TYP_FLT_TOP 0x08 /* flat top */ -#define FFT_WIN_TYP_WELCH 0x09 /* welch */ +#define FFT_WIN_TYP_RECTANGLE FFTWindow::Rectangle /* rectangle (Box car) */ +#define FFT_WIN_TYP_HAMMING FFTWindow::Hamming /* hamming */ +#define FFT_WIN_TYP_HANN FFTWindow::Hann /* hann */ +#define FFT_WIN_TYP_TRIANGLE FFTWindow::Triangle /* triangle (Bartlett) */ +#define FFT_WIN_TYP_NUTTALL FFTWindow::Nuttall /* nuttall */ +#define FFT_WIN_TYP_BLACKMAN FFTWindow::Blackman /* blackman */ +#define FFT_WIN_TYP_BLACKMAN_NUTTALL \ + FFTWindow::Blackman_Nuttall /* blackman nuttall */ +#define FFT_WIN_TYP_BLACKMAN_HARRIS \ + FFTWindow::Blackman_Harris /* blackman harris*/ +#define FFT_WIN_TYP_FLT_TOP FFTWindow::Flat_top /* flat top */ +#define FFT_WIN_TYP_WELCH FFTWindow::Welch /* welch */ /*Mathematial constants*/ #define twoPi 6.28318531 #define fourPi 12.56637061 #define sixPi 18.84955593 #ifdef __AVR__ - static const double _c1[]PROGMEM = {0.0000000000, 0.7071067812, 0.9238795325, 0.9807852804, - 0.9951847267, 0.9987954562, 0.9996988187, 0.9999247018, - 0.9999811753, 0.9999952938, 0.9999988235, 0.9999997059, - 0.9999999265, 0.9999999816, 0.9999999954, 0.9999999989, - 0.9999999997}; - static const double _c2[]PROGMEM = {1.0000000000, 0.7071067812, 0.3826834324, 0.1950903220, - 0.0980171403, 0.0490676743, 0.0245412285, 0.0122715383, - 0.0061358846, 0.0030679568, 0.0015339802, 0.0007669903, - 0.0003834952, 0.0001917476, 0.0000958738, 0.0000479369, - 0.0000239684}; +static const double _c1[] PROGMEM = { + 0.0000000000, 0.7071067812, 0.9238795325, 0.9807852804, 0.9951847267, + 0.9987954562, 0.9996988187, 0.9999247018, 0.9999811753, 0.9999952938, + 0.9999988235, 0.9999997059, 0.9999999265, 0.9999999816, 0.9999999954, + 0.9999999989, 0.9999999997}; +static const double _c2[] PROGMEM = { + 1.0000000000, 0.7071067812, 0.3826834324, 0.1950903220, 0.0980171403, + 0.0490676743, 0.0245412285, 0.0122715383, 0.0061358846, 0.0030679568, + 0.0015339802, 0.0007669903, 0.0003834952, 0.0001917476, 0.0000958738, + 0.0000479369, 0.0000239684}; #endif class arduinoFFT { public: - /* Constructor */ - arduinoFFT(void); - arduinoFFT(double *vReal, double *vImag, uint16_t samples, double samplingFrequency); - /* Destructor */ - ~arduinoFFT(void); - /* Functions */ - uint8_t Revision(void); - uint8_t Exponent(uint16_t value); + /* Constructor */ + arduinoFFT(void); + arduinoFFT(double *vReal, double *vImag, uint16_t samples, + double samplingFrequency); + /* Destructor */ + ~arduinoFFT(void); + /* Functions */ + uint8_t Revision(void); + uint8_t Exponent(uint16_t value); - void ComplexToMagnitude(double *vReal, double *vImag, uint16_t samples); - void Compute(double *vReal, double *vImag, uint16_t samples, uint8_t dir); - void Compute(double *vReal, double *vImag, uint16_t samples, uint8_t power, uint8_t dir); - void DCRemoval(double *vData, uint16_t samples); - double MajorPeak(double *vD, uint16_t samples, double samplingFrequency); - void MajorPeak(double *vD, uint16_t samples, double samplingFrequency, double *f, double *v); - void Windowing(double *vData, uint16_t samples, uint8_t windowType, uint8_t dir); + void ComplexToMagnitude(double *vReal, double *vImag, uint16_t samples); + void Compute(double *vReal, double *vImag, uint16_t samples, + FFTDirection dir); + void Compute(double *vReal, double *vImag, uint16_t samples, uint8_t power, + FFTDirection dir); + void DCRemoval(double *vData, uint16_t samples); + double MajorPeak(double *vD, uint16_t samples, double samplingFrequency); + void MajorPeak(double *vD, uint16_t samples, double samplingFrequency, + double *f, double *v); + void Windowing(double *vData, uint16_t samples, FFTWindow windowType, + FFTDirection dir); - void ComplexToMagnitude(); - void Compute(uint8_t dir); - void DCRemoval(); - double MajorPeak(); - void MajorPeak(double *f, double *v); - void Windowing(uint8_t windowType, uint8_t dir); + void ComplexToMagnitude(); + void Compute(FFTDirection dir); + void DCRemoval(); + double MajorPeak(); + void MajorPeak(double *f, double *v); + void Windowing(FFTWindow windowType, FFTDirection dir); - double MajorPeakParabola(); + double MajorPeakParabola(); private: - /* Variables */ - uint16_t _samples; - double _samplingFrequency; - double *_vReal; - double *_vImag; - uint8_t _power; - /* Functions */ - void Swap(double *x, double *y); - void Parabola(double x1, double y1, double x2, double y2, double x3, double y3, double *a, double *b, double *c); + /* Variables */ + uint16_t _samples; + double _samplingFrequency; + double *_vReal; + double *_vImag; + uint8_t _power; + /* Functions */ + void Swap(double *x, double *y); + void Parabola(double x1, double y1, double x2, double y2, double x3, + double y3, double *a, double *b, double *c); }; #endif From fa12c7b8f02a30357049f0f2c43241f6d1fdb6e9 Mon Sep 17 00:00:00 2001 From: Bjorn <75190918+BjornTheProgrammer@users.noreply.github.com> Date: Fri, 17 Mar 2023 23:31:04 -0700 Subject: [PATCH 09/12] Updated Examples to use Non-Deprecated Methods Updated the examples FFT_01, FFT_02, FFT_03, FFT_04, and FFT_05. --- Examples/FFT_01/FFT_01.ino | 12 +++++++----- Examples/FFT_02/FFT_02.ino | 13 ++++++------- Examples/FFT_03/FFT_03.ino | 11 ++++++----- Examples/FFT_04/FFT_04.ino | 11 ++++++----- Examples/FFT_05/FFT_05.ino | 11 ++++++----- 5 files changed, 31 insertions(+), 27 deletions(-) diff --git a/Examples/FFT_01/FFT_01.ino b/Examples/FFT_01/FFT_01.ino index 39e3708..c60e054 100644 --- a/Examples/FFT_01/FFT_01.ino +++ b/Examples/FFT_01/FFT_01.ino @@ -30,7 +30,7 @@ #include "arduinoFFT.h" -arduinoFFT FFT = arduinoFFT(); /* Create FFT object */ +arduinoFFT FFT; /* These values can be changed in order to evaluate the functions */ @@ -67,21 +67,23 @@ void loop() //vReal[i] = uint8_t((amplitude * (sin((i * (twoPi * cycles)) / samples) + 1.0)) / 2.0);/* Build data displaced on the Y axis to include only positive values*/ vImag[i] = 0.0; //Imaginary part must be zeroed in case of looping to avoid wrong calculations and overflows } + + FFT = arduinoFFT(vReal, vImag, samples, samplingFrequency); /* Create FFT object */ /* Print the results of the simulated sampling according to time */ Serial.println("Data:"); PrintVector(vReal, samples, SCL_TIME); - FFT.Windowing(vReal, samples, FFT_WIN_TYP_HAMMING, FFT_FORWARD); /* Weigh data */ + FFT.Windowing(FFT_WIN_TYP_HAMMING, FFT_FORWARD); /* Weigh data */ Serial.println("Weighed data:"); PrintVector(vReal, samples, SCL_TIME); - FFT.Compute(vReal, vImag, samples, FFT_FORWARD); /* Compute FFT */ + FFT.Compute(FFT_FORWARD); /* Compute FFT */ Serial.println("Computed Real values:"); PrintVector(vReal, samples, SCL_INDEX); Serial.println("Computed Imaginary values:"); PrintVector(vImag, samples, SCL_INDEX); - FFT.ComplexToMagnitude(vReal, vImag, samples); /* Compute magnitudes */ + FFT.ComplexToMagnitude(); /* Compute magnitudes */ Serial.println("Computed magnitudes:"); PrintVector(vReal, (samples >> 1), SCL_FREQUENCY); - double x = FFT.MajorPeak(vReal, samples, samplingFrequency); + double x = FFT.MajorPeak(); Serial.println(x, 6); while(1); /* Run Once */ // delay(2000); /* Repeat after delay */ diff --git a/Examples/FFT_02/FFT_02.ino b/Examples/FFT_02/FFT_02.ino index e6262c9..367814d 100644 --- a/Examples/FFT_02/FFT_02.ino +++ b/Examples/FFT_02/FFT_02.ino @@ -23,7 +23,7 @@ #include "arduinoFFT.h" -arduinoFFT FFT = arduinoFFT(); /* Create FFT object */ +arduinoFFT FFT; /* These values can be changed in order to evaluate the functions */ @@ -31,7 +31,6 @@ These values can be changed in order to evaluate the functions const uint16_t samples = 64; const double sampling = 40; const uint8_t amplitude = 4; -uint8_t exponent; const double startFrequency = 2; const double stopFrequency = 16.4; const double step_size = 0.1; @@ -55,7 +54,6 @@ void setup() Serial.begin(115200); while(!Serial); Serial.println("Ready"); - exponent = FFT.Exponent(samples); } void loop() @@ -74,18 +72,19 @@ void loop() /*Serial.println("Data:"); PrintVector(vReal, samples, SCL_TIME);*/ time=millis(); - FFT.Windowing(vReal, samples, FFT_WIN_TYP_HAMMING, FFT_FORWARD); /* Weigh data */ + FFT = arduinoFFT(vReal, vImag, samples, sampling); /* Create FFT object */ + FFT.Windowing(FFT_WIN_TYP_HAMMING, FFT_FORWARD); /* Weigh data */ /*Serial.println("Weighed data:"); PrintVector(vReal, samples, SCL_TIME);*/ - FFT.Compute(vReal, vImag, samples, exponent, FFT_FORWARD); /* Compute FFT */ + FFT.Compute(FFT_FORWARD); /* Compute FFT */ /*Serial.println("Computed Real values:"); PrintVector(vReal, samples, SCL_INDEX); Serial.println("Computed Imaginary values:"); PrintVector(vImag, samples, SCL_INDEX);*/ - FFT.ComplexToMagnitude(vReal, vImag, samples); /* Compute magnitudes */ + FFT.ComplexToMagnitude(); /* Compute magnitudes */ /*Serial.println("Computed magnitudes:"); PrintVector(vReal, (samples >> 1), SCL_FREQUENCY);*/ - double x = FFT.MajorPeak(vReal, samples, sampling); + double x = FFT.MajorPeak(); Serial.print(frequency); Serial.print(": \t\t"); Serial.print(x, 4); diff --git a/Examples/FFT_03/FFT_03.ino b/Examples/FFT_03/FFT_03.ino index b7c3fba..ac3fa12 100644 --- a/Examples/FFT_03/FFT_03.ino +++ b/Examples/FFT_03/FFT_03.ino @@ -20,7 +20,7 @@ #include "arduinoFFT.h" -arduinoFFT FFT = arduinoFFT(); /* Create FFT object */ +arduinoFFT FFT; /* These values can be changed in order to evaluate the functions */ @@ -64,21 +64,22 @@ void loop() } microseconds += sampling_period_us; } + FFT = arduinoFFT(vReal, vImag, samples, samplingFrequency); /* Create FFT object */ /* Print the results of the sampling according to time */ Serial.println("Data:"); PrintVector(vReal, samples, SCL_TIME); - FFT.Windowing(vReal, samples, FFT_WIN_TYP_HAMMING, FFT_FORWARD); /* Weigh data */ + FFT.Windowing(FFT_WIN_TYP_HAMMING, FFT_FORWARD); /* Weigh data */ Serial.println("Weighed data:"); PrintVector(vReal, samples, SCL_TIME); - FFT.Compute(vReal, vImag, samples, FFT_FORWARD); /* Compute FFT */ + FFT.Compute(FFT_FORWARD); /* Compute FFT */ Serial.println("Computed Real values:"); PrintVector(vReal, samples, SCL_INDEX); Serial.println("Computed Imaginary values:"); PrintVector(vImag, samples, SCL_INDEX); - FFT.ComplexToMagnitude(vReal, vImag, samples); /* Compute magnitudes */ + FFT.ComplexToMagnitude(); /* Compute magnitudes */ Serial.println("Computed magnitudes:"); PrintVector(vReal, (samples >> 1), SCL_FREQUENCY); - double x = FFT.MajorPeak(vReal, samples, samplingFrequency); + double x = FFT.MajorPeak(); Serial.println(x, 6); //Print out what frequency is the most dominant. while(1); /* Run Once */ // delay(2000); /* Repeat after delay */ diff --git a/Examples/FFT_04/FFT_04.ino b/Examples/FFT_04/FFT_04.ino index 7aa5bb6..36c2ac7 100644 --- a/Examples/FFT_04/FFT_04.ino +++ b/Examples/FFT_04/FFT_04.ino @@ -31,7 +31,7 @@ #include "arduinoFFT.h" -arduinoFFT FFT = arduinoFFT(); /* Create FFT object */ +arduinoFFT FFT; /* These values can be changed in order to evaluate the functions */ @@ -68,11 +68,12 @@ void loop() //vReal[i] = uint8_t((amplitude * (sin((i * (twoPi * cycles)) / samples) + 1.0)) / 2.0);/* Build data displaced on the Y axis to include only positive values*/ vImag[i] = 0.0; //Imaginary part must be zeroed in case of looping to avoid wrong calculations and overflows } - FFT.Windowing(vReal, samples, FFT_WIN_TYP_HAMMING, FFT_FORWARD); /* Weigh data */ - FFT.Compute(vReal, vImag, samples, FFT_FORWARD); /* Compute FFT */ - FFT.ComplexToMagnitude(vReal, vImag, samples); /* Compute magnitudes */ + FFT = arduinoFFT(vReal, vImag, samples, samplingFrequency); /* Create FFT object */ + FFT.Windowing(FFT_WIN_TYP_HAMMING, FFT_FORWARD); /* Weigh data */ + FFT.Compute(FFT_FORWARD); /* Compute FFT */ + FFT.ComplexToMagnitude(); /* Compute magnitudes */ PrintVector(vReal, samples>>1, SCL_PLOT); - double x = FFT.MajorPeak(vReal, samples, samplingFrequency); + double x = FFT.MajorPeak(); while(1); /* Run Once */ // delay(2000); /* Repeat after delay */ } diff --git a/Examples/FFT_05/FFT_05.ino b/Examples/FFT_05/FFT_05.ino index 7050b10..6a87686 100644 --- a/Examples/FFT_05/FFT_05.ino +++ b/Examples/FFT_05/FFT_05.ino @@ -31,7 +31,7 @@ #include "arduinoFFT.h" -arduinoFFT FFT = arduinoFFT(); /* Create FFT object */ +arduinoFFT FFT; /* These values can be changed in order to evaluate the functions */ @@ -68,23 +68,24 @@ void loop() //vReal[i] = uint8_t((amplitude * (sin((i * (twoPi * cycles)) / samples) + 1.0)) / 2.0);/* Build data displaced on the Y axis to include only positive values*/ vImag[i] = 0.0; //Imaginary part must be zeroed in case of looping to avoid wrong calculations and overflows } + FFT = arduinoFFT(vReal, vImag, samples, samplingFrequency); /* Create FFT object */ /* Print the results of the simulated sampling according to time */ Serial.println("Data:"); PrintVector(vReal, samples, SCL_TIME); - FFT.Windowing(vReal, samples, FFT_WIN_TYP_HAMMING, FFT_FORWARD); /* Weigh data */ + FFT.Windowing(FFT_WIN_TYP_HAMMING, FFT_FORWARD); /* Weigh data */ Serial.println("Weighed data:"); PrintVector(vReal, samples, SCL_TIME); - FFT.Compute(vReal, vImag, samples, FFT_FORWARD); /* Compute FFT */ + FFT.Compute(FFT_FORWARD); /* Compute FFT */ Serial.println("Computed Real values:"); PrintVector(vReal, samples, SCL_INDEX); Serial.println("Computed Imaginary values:"); PrintVector(vImag, samples, SCL_INDEX); - FFT.ComplexToMagnitude(vReal, vImag, samples); /* Compute magnitudes */ + FFT.ComplexToMagnitude(); /* Compute magnitudes */ Serial.println("Computed magnitudes:"); PrintVector(vReal, (samples >> 1), SCL_FREQUENCY); double x; double v; - FFT.MajorPeak(vReal, samples, samplingFrequency, &x, &v); + FFT.MajorPeak(&x, &v); Serial.print(x, 6); Serial.print(", "); Serial.println(v, 6); From f9d6fa7dae2b8650a4a6c38410e78a33bbf387ee Mon Sep 17 00:00:00 2001 From: Bjorn <75190918+BjornTheProgrammer@users.noreply.github.com> Date: Fri, 17 Mar 2023 23:41:46 -0700 Subject: [PATCH 10/12] Updated Documentation to Reflect Actual Methods --- README.md | 21 +++++++++++++-------- 1 file changed, 13 insertions(+), 8 deletions(-) diff --git a/README.md b/README.md index 986a6c2..0a36864 100644 --- a/README.md +++ b/README.md @@ -19,12 +19,14 @@ To install this library, just place this entire folder as a subfolder in your Ar When installed, this library should look like: +``` Arduino\libraries\arduinoFTT (this library's folder) Arduino\libraries\arduinoFTT\arduinoFTT.cpp (the library implementation file, uses 32 bits floats vectors) Arduino\libraries\arduinoFTT\arduinoFTT.h (the library header file, uses 32 bits floats vectors) Arduino\libraries\arduinoFTT\keywords.txt (the syntax coloring file) Arduino\libraries\arduinoFTT\examples (the examples in the "open" menu) Arduino\libraries\arduinoFTT\readme.md (this file) +``` ### Building on Arduino @@ -44,28 +46,31 @@ select arduinoFTT. This will add a corresponding line to the top of your sketch * Spectrum table? ### API +The exclamation mark `!` denotes that this method is deprecated and may be removed on future revisions. -* **arduinoFFT**(void); +* **!arduinoFFT**(void); * **arduinoFFT**(double *vReal, double *vImag, uint16_t samples, double samplingFrequency); Constructor * **~arduinoFFT**(void); Destructor -* **ComplexToMagnitude**(double *vReal, double *vImag, uint16_t samples); +* **!ComplexToMagnitude**(double *vReal, double *vImag, uint16_t samples); * **ComplexToMagnitude**(); -* **Compute**(double *vReal, double *vImag, uint16_t samples, uint8_t dir); -* **Compute**(double *vReal, double *vImag, uint16_t samples, uint8_t power, uint8_t dir); +* **!Compute**(double *vReal, double *vImag, uint16_t samples, uint8_t dir); +* **!Compute**(double *vReal, double *vImag, uint16_t samples, uint8_t power, uint8_t dir); * **Compute**(uint8_t dir); Calcuates the Fast Fourier Transform. -* **DCRemoval**(double *vData, uint16_t samples); +* **!DCRemoval**(double *vData, uint16_t samples); * **DCRemoval**(); Removes the DC component from the sample data. -* **MajorPeak**(double *vD, uint16_t samples, double samplingFrequency); +* **!MajorPeak**(double *vD, uint16_t samples, double samplingFrequency); +* **!MajorPeak**(double *vD, uint16_t samples, double samplingFrequency, double *f, double *v); * **MajorPeak**(); +* **MajorPeak**(double *f, double *v); * **MajorPeakParabola**(); Looks for and returns the frequency of the biggest spike in the analyzed signal. * **Revision**(void); Returns the library revision. -* **Windowing**(double *vData, uint16_t samples, uint8_t windowType, uint8_t dir); +* **!Windowing**(double *vData, uint16_t samples, uint8_t windowType, uint8_t dir); * **Windowing**(uint8_t windowType, uint8_t dir); Performs a windowing function on the values array. The possible windowing options are: * FFT_WIN_TYP_RECTANGLE @@ -78,5 +83,5 @@ Performs a windowing function on the values array. The possible windowing option * FFT_WIN_TYP_BLACKMAN_HARRIS * FFT_WIN_TYP_FLT_TOP * FFT_WIN_TYP_WELCH -* **Exponent**(uint16_t value); +* **!Exponent**(uint16_t value); Calculates and returns the base 2 logarithm of the given value. From b0e4c69bcb0d024e2923f69b3fc29b3e03184ae6 Mon Sep 17 00:00:00 2001 From: Enrique Condes Date: Tue, 25 Jul 2023 23:39:34 +0800 Subject: [PATCH 11/12] Fix issue 84 --- Examples/FFT_02/FFT_02.ino | 2 +- library.json | 2 +- library.properties | 2 +- src/arduinoFFT.cpp | 12 +++++++----- 4 files changed, 10 insertions(+), 8 deletions(-) diff --git a/Examples/FFT_02/FFT_02.ino b/Examples/FFT_02/FFT_02.ino index 367814d..c71e588 100644 --- a/Examples/FFT_02/FFT_02.ino +++ b/Examples/FFT_02/FFT_02.ino @@ -3,7 +3,7 @@ Example of use of the FFT library to compute FFT for several signals over a range of frequencies. The exponent is calculated once before the execution since it is a constant. This saves resources during the execution of the sketch and reduces the compiled size. - The sketch shows the time that the computing is taking. + The sketch shows the time that the computation takes. Copyright (C) 2014 Enrique Condes This program is free software: you can redistribute it and/or modify diff --git a/library.json b/library.json index eae9358..23366bd 100644 --- a/library.json +++ b/library.json @@ -20,7 +20,7 @@ "email": "contact@arduinoos.com" } ], - "version": "1.6", + "version": "1.6.1", "frameworks": ["arduino","mbed","espidf"], "platforms": "*", "headers": "arduinoFFT.h" diff --git a/library.properties b/library.properties index ee5f5ec..712867e 100644 --- a/library.properties +++ b/library.properties @@ -1,5 +1,5 @@ name=arduinoFFT -version=1.6 +version=1.6.1 author=Enrique Condes maintainer=Enrique Condes sentence=A library for implementing floating point Fast Fourier Transform calculations on Arduino. diff --git a/src/arduinoFFT.cpp b/src/arduinoFFT.cpp index 15ab81f..417cd27 100644 --- a/src/arduinoFFT.cpp +++ b/src/arduinoFFT.cpp @@ -103,9 +103,10 @@ void arduinoFFT::Compute(FFTDirection dir) { } // Scaling for reverse transform / if (dir != FFT_FORWARD) { + double reciprocal = 1.0 / this->_samples; for (uint16_t i = 0; i < this->_samples; i++) { - this->_vReal[i] /= this->_samples; - this->_vImag[i] /= this->_samples; + this->_vReal[i] *= reciprocal; + this->_vImag[i] *= reciprocal; } } } @@ -169,9 +170,10 @@ void arduinoFFT::Compute(double *vReal, double *vImag, uint16_t samples, } // Scaling for reverse transform if (dir != FFT_FORWARD) { + double reciprocal = 1.0 / samples; for (uint16_t i = 0; i < samples; i++) { - vReal[i] /= samples; - vImag[i] /= samples; + vReal[i] *= reciprocal; + vImag[i] *= reciprocal; } } } @@ -535,7 +537,7 @@ uint8_t arduinoFFT::Exponent(uint16_t value) { #warning("This method may not be accessible on future revisions.") // Calculates the base 2 logarithm of a value uint8_t result = 0; - while (((value >> result) & 1) != 1) + while (value >>= 1) result++; return (result); } From c7845ab60eece7646e1a38455a0956acd72ab905 Mon Sep 17 00:00:00 2001 From: Enrique Condes Date: Fri, 29 Dec 2023 11:00:17 +0800 Subject: [PATCH 12/12] Fix bug on synthetic data generation --- Examples/FFT_01/FFT_01.ino | 6 +++--- Examples/FFT_02/FFT_02.ino | 4 ++-- Examples/FFT_04/FFT_04.ino | 6 +++--- Examples/FFT_05/FFT_05.ino | 6 +++--- README.md | 2 +- changeLog.txt | 6 ++++++ library.json | 2 +- library.properties | 2 +- src/arduinoFFT.cpp | 2 +- 9 files changed, 21 insertions(+), 15 deletions(-) diff --git a/Examples/FFT_01/FFT_01.ino b/Examples/FFT_01/FFT_01.ino index c60e054..c41d22c 100644 --- a/Examples/FFT_01/FFT_01.ino +++ b/Examples/FFT_01/FFT_01.ino @@ -60,11 +60,11 @@ void setup() void loop() { /* Build raw data */ - double cycles = (((samples-1) * signalFrequency) / samplingFrequency); //Number of signal cycles that the sampling will read + double ratio = twoPi * signalFrequency / samplingFrequency; // Fraction of a complete cycle stored at each sample (in radians) for (uint16_t i = 0; i < samples; i++) { - vReal[i] = int8_t((amplitude * (sin((i * (twoPi * cycles)) / samples))) / 2.0);/* Build data with positive and negative values*/ - //vReal[i] = uint8_t((amplitude * (sin((i * (twoPi * cycles)) / samples) + 1.0)) / 2.0);/* Build data displaced on the Y axis to include only positive values*/ + vReal[i] = int8_t(amplitude * sin(i * ratio) / 2.0);/* Build data with positive and negative values*/ + //vReal[i] = uint8_t((amplitude * (sin(i * ratio) + 1.0)) / 2.0);/* Build data displaced on the Y axis to include only positive values*/ vImag[i] = 0.0; //Imaginary part must be zeroed in case of looping to avoid wrong calculations and overflows } diff --git a/Examples/FFT_02/FFT_02.ino b/Examples/FFT_02/FFT_02.ino index c71e588..4c7ecd6 100644 --- a/Examples/FFT_02/FFT_02.ino +++ b/Examples/FFT_02/FFT_02.ino @@ -63,10 +63,10 @@ void loop() for(double frequency = startFrequency; frequency<=stopFrequency; frequency+=step_size) { /* Build raw data */ - double cycles = (((samples-1) * frequency) / sampling); + double ratio = twoPi * frequency / sampling; // Fraction of a complete cycle stored at each sample (in radians) for (uint16_t i = 0; i < samples; i++) { - vReal[i] = int8_t((amplitude * (sin((i * (twoPi * cycles)) / samples))) / 2.0); + vReal[i] = int8_t(amplitude * sin(i * ratio) / 2.0);/* Build data with positive and negative values*/ vImag[i] = 0; //Reset the imaginary values vector for each new frequency } /*Serial.println("Data:"); diff --git a/Examples/FFT_04/FFT_04.ino b/Examples/FFT_04/FFT_04.ino index 36c2ac7..be6a970 100644 --- a/Examples/FFT_04/FFT_04.ino +++ b/Examples/FFT_04/FFT_04.ino @@ -61,11 +61,11 @@ void setup() void loop() { /* Build raw data */ - double cycles = (((samples-1) * signalFrequency) / samplingFrequency); //Number of signal cycles that the sampling will read + double ratio = twoPi * signalFrequency / samplingFrequency; // Fraction of a complete cycle stored at each sample (in radians) for (uint16_t i = 0; i < samples; i++) { - vReal[i] = int8_t((amplitude * (sin((i * (twoPi * cycles)) / samples))) / 2.0);/* Build data with positive and negative values*/ - //vReal[i] = uint8_t((amplitude * (sin((i * (twoPi * cycles)) / samples) + 1.0)) / 2.0);/* Build data displaced on the Y axis to include only positive values*/ + vReal[i] = int8_t(amplitude * sin(i * ratio) / 2.0);/* Build data with positive and negative values*/ + //vReal[i] = uint8_t((amplitude * (sin(i * ratio) + 1.0)) / 2.0);/* Build data displaced on the Y axis to include only positive values*/ vImag[i] = 0.0; //Imaginary part must be zeroed in case of looping to avoid wrong calculations and overflows } FFT = arduinoFFT(vReal, vImag, samples, samplingFrequency); /* Create FFT object */ diff --git a/Examples/FFT_05/FFT_05.ino b/Examples/FFT_05/FFT_05.ino index 6a87686..a382b7f 100644 --- a/Examples/FFT_05/FFT_05.ino +++ b/Examples/FFT_05/FFT_05.ino @@ -61,11 +61,11 @@ void setup() void loop() { /* Build raw data */ - double cycles = (((samples-1) * signalFrequency) / samplingFrequency); //Number of signal cycles that the sampling will read + double ratio = twoPi * signalFrequency / samplingFrequency; // Fraction of a complete cycle stored at each sample (in radians) for (uint16_t i = 0; i < samples; i++) { - vReal[i] = int8_t((amplitude * (sin((i * (twoPi * cycles)) / samples))) / 2.0);/* Build data with positive and negative values*/ - //vReal[i] = uint8_t((amplitude * (sin((i * (twoPi * cycles)) / samples) + 1.0)) / 2.0);/* Build data displaced on the Y axis to include only positive values*/ + vReal[i] = int8_t(amplitude * sin(i * ratio) / 2.0);/* Build data with positive and negative values*/ + //vReal[i] = uint8_t((amplitude * (sin(i * ratio) + 1.0)) / 2.0);/* Build data displaced on the Y axis to include only positive values*/ vImag[i] = 0.0; //Imaginary part must be zeroed in case of looping to avoid wrong calculations and overflows } FFT = arduinoFFT(vReal, vImag, samples, samplingFrequency); /* Create FFT object */ diff --git a/README.md b/README.md index 0a36864..93cd705 100644 --- a/README.md +++ b/README.md @@ -58,7 +58,7 @@ Destructor * **!Compute**(double *vReal, double *vImag, uint16_t samples, uint8_t dir); * **!Compute**(double *vReal, double *vImag, uint16_t samples, uint8_t power, uint8_t dir); * **Compute**(uint8_t dir); -Calcuates the Fast Fourier Transform. +Calculates the Fast Fourier Transform. * **!DCRemoval**(double *vData, uint16_t samples); * **DCRemoval**(); Removes the DC component from the sample data. diff --git a/changeLog.txt b/changeLog.txt index 33612ef..e15375d 100644 --- a/changeLog.txt +++ b/changeLog.txt @@ -1,3 +1,9 @@ +28/12/23 v1.6.2 +Fix issue 52 and symplify calculation of synthetic data for examples 1, 2, 4, and 5. + +25/07/23 v1.6.1 +Code optimization for speed improvements. See issue 84 for more details. + 03/09/23 v1.6 Include some functionality from development branch. diff --git a/library.json b/library.json index 23366bd..87dadea 100644 --- a/library.json +++ b/library.json @@ -20,7 +20,7 @@ "email": "contact@arduinoos.com" } ], - "version": "1.6.1", + "version": "1.6.2", "frameworks": ["arduino","mbed","espidf"], "platforms": "*", "headers": "arduinoFFT.h" diff --git a/library.properties b/library.properties index 712867e..951628c 100644 --- a/library.properties +++ b/library.properties @@ -1,5 +1,5 @@ name=arduinoFFT -version=1.6.1 +version=1.6.2 author=Enrique Condes maintainer=Enrique Condes sentence=A library for implementing floating point Fast Fourier Transform calculations on Arduino. diff --git a/src/arduinoFFT.cpp b/src/arduinoFFT.cpp index 417cd27..02b1b67 100644 --- a/src/arduinoFFT.cpp +++ b/src/arduinoFFT.cpp @@ -292,7 +292,7 @@ void arduinoFFT::Windowing(FFTWindow windowType, FFTDirection dir) { void arduinoFFT::Windowing(double *vData, uint16_t samples, FFTWindow windowType, FFTDirection dir) { // Weighing factors are computed once before multiple use of FFT -// The weighing function is symetric; half the weighs are recorded +// The weighing function is symmetric; half the weighs are recorded #warning("This method is deprecated and may be removed on future revisions.") double samplesMinusOne = (double(samples) - 1.0); for (uint16_t i = 0; i < (samples >> 1); i++) {