Java Apache Commons Mathを使って、高速フーリエ変換する

今回は高速フーリエ変換です。
高速フーリエ変換には、Apache Commons Mathを使います。

Apache Commons Mathの入手

mavenを使って入手

pom.xmlは、こんな感じ。
今回はフーリエ変換なので、Apache Commonsの中でも「Apache Commons Math」が対象。

<!-- https://mvnrepository.com/artifact/org.apache.commons/commons-math3 -->
<dependency>
    <groupId>org.apache.commons</groupId>
    <artifactId>commons-math3</artifactId>
    <version>3.6.1</version>
</dependency>

<groupId>org.apache.commons</groupId>

<artifactId>commons-math3</artifactId>

</dependency>

なかなか4.0は正式版が出ないねｗ

手動で入手

手動で入手するのであれば、以下のjarを入手することになります。
「Apache Commons Math」は依存関係がないので、手動でも簡単に手に入れられますね。
https://mvnrepository.com/repos/central

commons-math3-3.6.1.jar

※2023年11月現在です。

ちなみに、Apache License 2.0で提供されています。

Apache Commons Mathを使って、高速フーリエ変換するサンプル

220Hz,振幅=2と440Hz,振幅=3の正弦波を合成。
これをフーリエ変換して、周波数領域のデータを出力。
さらに逆フーリエ変換して、時間領域に戻してデータを出力します。

import org.apache.commons.math3.complex.Complex;
import org.apache.commons.math3.transform.DftNormalization;
import org.apache.commons.math3.transform.FastFourierTransformer;
import org.apache.commons.math3.transform.TransformType;

public class SinFFTTest {

	private static final int SamplesPerSec = 1024;

	public static void main(String[] args) {

		//Commons Mathの高速フーリエ変換のクラス
		FastFourierTransformer fft = new FastFourierTransformer(DftNormalization.STANDARD);

		//データ準備（A3:220Hz+A4:440Hz）
		double[] data = createSinWave();
		print("-- 時間領域（作成） --", data);

		//フーリエ変換
		Complex[] frDomain = fft.transform(data, TransformType.FORWARD);
		printFr("-- 周波数領域(Hz,絶対値) --", frDomain);

		//逆フーリエ変換
		Complex[] tmDomain = fft.transform(frDomain, TransformType.INVERSE);
		printTm("-- 時間領域（戻し） --", tmDomain);
	}

	private static double[] createSinWave() {

		double[] data = new double[SamplesPerSec];

		//220Hz,振幅=2と440Hz,振幅=3の正弦波を合成して、１秒分のデータを作成
		for (int i = 0; i<SamplesPerSec ; i++) {
			data[i] = 2.0 * Math.sin(2.0 * Math.PI * 220.0 * i / SamplesPerSec)
					+ 3.0 * Math.sin(2.0 * Math.PI * 440.0 * i / SamplesPerSec);
		}

		return data;
	}

	private static void print(String msg, double[] data) {
		System.out.println(msg);
		for (double d : data) {
			System.out.println(d);
		}
	}

	private static void printFr(String msg, Complex[] data) {
		System.out.println(msg);
		for (int i = 0; i<data.length/2 ; i++) {
			System.out.println(i + "\t" + data[i].abs());
		}
	}

	private static void printTm(String msg, Complex[] data) {
		System.out.println(msg);
		for (Complex d : data) {
			System.out.println(d.getReal());
		}
	}
}

import org.apache.commons.math3.complex.Complex;

import org.apache.commons.math3.transform.DftNormalization;

import org.apache.commons.math3.transform.FastFourierTransformer;

import org.apache.commons.math3.transform.TransformType;

public class SinFFTTest {

private static final int SamplesPerSec = 1024;

public static void main(String[] args) {

//Commons Mathの高速フーリエ変換のクラス

FastFourierTransformer fft = new FastFourierTransformer(DftNormalization.STANDARD);

//データ準備（A3:220Hz+A4:440Hz）

double[] data = createSinWave();

print("-- 時間領域（作成） --", data);

//フーリエ変換

Complex[] frDomain = fft.transform(data, TransformType.FORWARD);

printFr("-- 周波数領域(Hz,絶対値) --", frDomain);

//逆フーリエ変換

Complex[] tmDomain = fft.transform(frDomain, TransformType.INVERSE);

printTm("-- 時間領域（戻し） --", tmDomain);

}

private static double[] createSinWave() {

double[] data = new double[SamplesPerSec];

//220Hz,振幅=2と440Hz,振幅=3の正弦波を合成して、１秒分のデータを作成

for (int i = 0; i<SamplesPerSec ; i++) {

data[i] = 2.0 * Math.sin(2.0 * Math.PI * 220.0 * i / SamplesPerSec)

+ 3.0 * Math.sin(2.0 * Math.PI * 440.0 * i / SamplesPerSec);

}

return data;

}

private static void print(String msg, double[] data) {

System.out.println(msg);

for (double d : data) {

System.out.println(d);

}

private static void printFr(String msg, Complex[] data) {

System.out.println(msg);

for (int i = 0; i<data.length/2 ; i++) {

System.out.println(i + "\t" + data[i].abs());

}

private static void printTm(String msg, Complex[] data) {

System.out.println(msg);

for (Complex d : data) {

System.out.println(d.getReal());

}

実行結果

作成したデータ、フーリエ変換した周波数領域のデータ、逆フーリエ変換した時間領域のデータが出力されます。
220Hzと440Hzでスパイクが出ています。

-- 時間領域（作成） --
0.0
3.234069540367903
-1.4639211732276465
（略）
-- 周波数領域(Hz,絶対値) --
0	2.858717572861284E-12
1	2.703783409928831E-12
2	2.6013291279036116E-12
（略）
218	2.754723339373896E-12
219	6.402309347522719E-12
220	1023.9999999999935
221	6.247030050571049E-12
222	2.541112064724628E-12
（略）
438	6.858590434315256E-12
439	1.70847996035536E-11
440	1535.9999999999832
441	1.6480719359497152E-11
442	8.608039025728898E-12
（略）
-- 時間領域（戻し） --
7.570730034197228E-18
3.2340695403678734
-1.4639211732276394
（略）

-- 時間領域（作成） --

0.0

3.234069540367903

-1.4639211732276465

（略）

-- 周波数領域(Hz,絶対値) --

0 2.858717572861284E-12

1 2.703783409928831E-12

2 2.6013291279036116E-12

（略）

218 2.754723339373896E-12

219 6.402309347522719E-12

220 1023.9999999999935

221 6.247030050571049E-12

222 2.541112064724628E-12

（略）

438 6.858590434315256E-12

439 1.70847996035536E-11

440 1535.9999999999832

441 1.6480719359497152E-11

442 8.608039025728898E-12

（略）

-- 時間領域（戻し） --

7.570730034197228E-18

3.2340695403678734

-1.4639211732276394

（略）

サンプルの解説

Commons MathのFastFourierTransformerクラスを使うと、（離散の）フーリエ変換が簡単にできます。
ただし、高速フーリエ変換ですので、扱えるのは2^nサイズのデータに限られます。

時間領域から周波数領域へのフーリエ変換（行き）は、
FastFourierTransformer#transform(double[], TransformType.FORWARD)
周波数領域から時間領域への逆フーリエ変換（帰り）は、
FastFourierTransformer#transform(Complex[], TransformType.INVERSE)
を使ってます。第２引数で行き、帰りを指定しています。
ほぼ？元に戻っていますね。

ちなみに、周波数の解析範囲は、標本化定理により、サンプリング周波数の半分（ナイキスト）までです。