COCKOS
CONFEDERATED FORUMS
Cockos : REAPER : NINJAM : Forums
Forum Home : Register : FAQ : Members List : Search :
Old 02-03-2020, 08:22 AM   #1
vasyan
Human being with feelings
 
vasyan's Avatar
 
Join Date: Oct 2019
Location: Sayan Mountains, Siberia
Posts: 22
Default How to use the WDL FFT library?

Unfortunately, I did not find any documentation and user manuals for the WDL FFT library. I tried the simplest conversion of the signal to the spectrum and vice versa, having previously applied the window function, but in the resulting array there were unreal numbers (NAN).
I suspect that such a transformation requires some feature that I do not know about. However, at the very beginning of working with the library, it called the library initialization function WDL_fft_init();
Can I see somewhere a simple working example of using this FFT library (namely, the implementation of a band-pass filter is of interest)?
Or be able to view user manuals or documentation for this library.
Thanks!
vasyan is offline   Reply With Quote
Old 02-04-2020, 04:52 AM   #2
vasyan
Human being with feelings
 
vasyan's Avatar
 
Join Date: Oct 2019
Location: Sayan Mountains, Siberia
Posts: 22
Default

No one knows?
Some libraries have this order of data output after FFT:
f - pointer on the destination array (frequencies),
where f [0...length(x)/2] = real values,
and f [length(x)/2+1...length(x)-1] = imaginary values.
However, nothing is said about this in a brief annotation to this library.
Other thoughts?
vasyan is offline   Reply With Quote
Old 02-04-2020, 05:18 AM   #3
Tale
Human being with feelings
 
Tale's Avatar
 
Join Date: Jul 2008
Location: The Netherlands
Posts: 3,645
Default

Here is some code that takes the "real" FFT, and then prints magnitude and phase for each bin. I hope this helps.

Code:
#include <math.h>
#include <stdio.h>

#include "WDL/fft.h"

int main()
{
	WDL_fft_init();

	const double srate = 44100.0;
	static const int len = 1024;

	WDL_FFT_REAL buf[len];

	// Fill buf[] with test signal.
	for (int i = 0; i < len; ++i)
	{
		const double t = (double)i / len;
		const double saw = 2*t - 1;
		buf[i] = (WDL_FFT_REAL)saw;
	}

	// Forward FFT
	WDL_real_fft(buf, len, 0);

	// Scale after forward FFT.
	for (int i = 0; i < len; ++i)
	{
		buf[i] /= (WDL_FFT_REAL)len;
	}

	// Bins are stored out of order, so use permutation table to fix.
	const int* order = WDL_fft_permute_tab(len / 2);

	printf("Bin\tFrequency\tMagnitude\tPhase\n");
	for (int i = 0; i <= len / 2; ++i)
	{
		WDL_FFT_COMPLEX* bin = (WDL_FFT_COMPLEX*)buf + order[i];

		double re, im;
		if (i == 0)
		{
			// DC (0 Hz)
			re = bin->re;
			im = 0.0;
		}
		else if (i == len / 2)
		{
			 // Nyquist frequency
			re = buf[1]; // i.e. DC bin->im
			im = 0.0;
		}
		else
		{
			re = bin->re;
			im = bin->im;
		}

		const double mag = sqrt(re*re + im*im);
		const double phase = atan2(im, re);

		const double freq = (double)i / len * srate;
		printf("%d\t%f\t%f\t%f\n", i, freq, mag, phase);
	}

	// Inverse FFT
	WDL_real_fft(buf, len, 1);

	// Scale after inverse FFT.
	for (int i = 0; i < len; ++i)
	{
		buf[i] *= (WDL_FFT_REAL)0.5;
	}

	// Now we're back to original test signal.
}

Last edited by Tale; 02-04-2020 at 05:31 AM. Reason: Typo in comment
Tale is offline   Reply With Quote
Old 02-05-2020, 07:54 AM   #4
vasyan
Human being with feelings
 
vasyan's Avatar
 
Join Date: Oct 2019
Location: Sayan Mountains, Siberia
Posts: 22
Default

Quote:
Originally Posted by Tale View Post
Here is some code that takes the "real" FFT, and then prints magnitude and phase for each bin. I hope this helps.
Thanks!
That way everything works. However, if I want to change the type of data as double precision, then I declare at the very beginning:
Code:
#define WDL_FFT_REALSIZE 8
The result is a significant discrepancy between the results, which are completely unacceptable.
Or is this algorithm not intended for double precision numbers?
vasyan is offline   Reply With Quote
Old 02-05-2020, 08:33 AM   #5
Tale
Human being with feelings
 
Tale's Avatar
 
Join Date: Jul 2008
Location: The Netherlands
Posts: 3,645
Default

Quote:
Originally Posted by vasyan View Post
However, if I want to change the type of data as double precision, then I declare at the very beginning:
Code:
#define WDL_FFT_REALSIZE 8
The result is a significant discrepancy between the results, which are completely unacceptable.
Or is this algorithm not intended for double precision numbers?
I think the problem is that you're now compiling your source file with WDL_FFT_REALSIZE == 8, but you're still compiling fft.c with WDL_FFT_REALSIZE == 4. It would therefore be better to define WDL_FFT_REALSIZE at the project level, so all source files use the same settings.
Tale is offline   Reply With Quote
Old 02-06-2020, 05:00 AM   #6
vasyan
Human being with feelings
 
vasyan's Avatar
 
Join Date: Oct 2019
Location: Sayan Mountains, Siberia
Posts: 22
Default

Quote:
Originally Posted by Tale View Post
I think the problem is that you're now compiling your source file with WDL_FFT_REALSIZE == 8, but you're still compiling fft.c with WDL_FFT_REALSIZE == 4. It would therefore be better to define WDL_FFT_REALSIZE at the project level, so all source files use the same settings.
Yes, you're right, I had to explicitly set the type in the fft.h file
Code:
 typedef double WDL_FFT_REAL;
then everything began to work perfectly.
Thank you for the consultation, you simply provided invaluable assistance in solving the problem!
vasyan is offline   Reply With Quote
Old 02-06-2020, 07:11 AM   #7
Tale
Human being with feelings
 
Tale's Avatar
 
Join Date: Jul 2008
Location: The Netherlands
Posts: 3,645
Default

Quote:
Originally Posted by vasyan View Post
Yes, you're right, I had to explicitly set the type in the fft.h file
Code:
 typedef double WDL_FFT_REAL;
then everything began to work perfectly.
Yeah, that should work. Do note that it might be more convenient add "WDL_FFT_REALSIZE=8" to your project settings or makefile. That way you don't have to edit fft.h.

Quote:
Originally Posted by vasyan View Post
Thank you for the consultation, you simply provided invaluable assistance in solving the problem!
You're welcome.
Tale is offline   Reply With Quote
Reply

Thread Tools
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off

Forum Jump


All times are GMT -7. The time now is 05:56 AM.


Powered by vBulletin® Version 3.8.11
Copyright ©2000 - 2024, vBulletin Solutions Inc.