本文实例讲述了java实现的傅里叶变化算法。分享给大家供大家参考,具体如下:
用java实现傅里叶变化 结果为复数形式 a+bi
废话不多说,实现代码如下,共两个class
fft.class 傅里叶变化功能实现代码
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
|
package fft.test; /************************************************************************* * compilation: javac fft.java execution: java fft n dependencies: complex.java * * compute the fft and inverse fft of a length n complex sequence. bare bones * implementation that runs in o(n log n) time. our goal is to optimize the * clarity of the code, rather than performance. * * limitations ----------- - assumes n is a power of 2 * * - not the most memory efficient algorithm (because it uses an object type for * representing complex numbers and because it re-allocates memory for the * subarray, instead of doing in-place or reusing a single temporary array) * *************************************************************************/ public class fft { // compute the fft of x[], assuming its length is a power of 2 public static complex[] fft(complex[] x) { int n = x.length; // base case if (n == 1 ) return new complex[] { x[ 0 ] }; // radix 2 cooley-tukey fft if (n % 2 != 0 ) { throw new runtimeexception( "n is not a power of 2" ); } // fft of even terms complex[] even = new complex[n / 2 ]; for ( int k = 0 ; k < n / 2 ; k++) { even[k] = x[ 2 * k]; } complex[] q = fft(even); // fft of odd terms complex[] odd = even; // reuse the array for ( int k = 0 ; k < n / 2 ; k++) { odd[k] = x[ 2 * k + 1 ]; } complex[] r = fft(odd); // combine complex[] y = new complex[n]; for ( int k = 0 ; k < n / 2 ; k++) { double kth = - 2 * k * math.pi / n; complex wk = new complex(math.cos(kth), math.sin(kth)); y[k] = q[k].plus(wk.times(r[k])); y[k + n / 2 ] = q[k].minus(wk.times(r[k])); } return y; } // compute the inverse fft of x[], assuming its length is a power of 2 public static complex[] ifft(complex[] x) { int n = x.length; complex[] y = new complex[n]; // take conjugate for ( int i = 0 ; i < n; i++) { y[i] = x[i].conjugate(); } // compute forward fft y = fft(y); // take conjugate again for ( int i = 0 ; i < n; i++) { y[i] = y[i].conjugate(); } // divide by n for ( int i = 0 ; i < n; i++) { y[i] = y[i].scale( 1.0 / n); } return y; } // compute the circular convolution of x and y public static complex[] cconvolve(complex[] x, complex[] y) { // should probably pad x and y with 0s so that they have same length // and are powers of 2 if (x.length != y.length) { throw new runtimeexception( "dimensions don't agree" ); } int n = x.length; // compute fft of each sequence,求值 complex[] a = fft(x); complex[] b = fft(y); // point-wise multiply,点值乘法 complex[] c = new complex[n]; for ( int i = 0 ; i < n; i++) { c[i] = a[i].times(b[i]); } // compute inverse fft,插值 return ifft(c); } // compute the linear convolution of x and y public static complex[] convolve(complex[] x, complex[] y) { complex zero = new complex( 0 , 0 ); complex[] a = new complex[ 2 * x.length]; // 2n次数界,高阶系数为0. for ( int i = 0 ; i < x.length; i++) a[i] = x[i]; for ( int i = x.length; i < 2 * x.length; i++) a[i] = zero; complex[] b = new complex[ 2 * y.length]; for ( int i = 0 ; i < y.length; i++) b[i] = y[i]; for ( int i = y.length; i < 2 * y.length; i++) b[i] = zero; return cconvolve(a, b); } // display an array of complex numbers to standard output public static void show(complex[] x, string title) { system.out.println(title); system.out.println( "-------------------" ); int complexlength = x.length; for ( int i = 0 ; i < complexlength; i++) { // 输出复数 // system.out.println(x[i]); // 输出幅值需要 * 2 / length system.out.println(x[i].abs() * 2 / complexlength); } system.out.println(); } /** * 将数组数据重组成2的幂次方输出 * * @param data * @return */ public static double [] pow2doublearr( double [] data) { // 创建新数组 double [] newdata = null ; int datalength = data.length; int sumnum = 2 ; while (sumnum < datalength) { sumnum = sumnum * 2 ; } int addlength = sumnum - datalength; if (addlength != 0 ) { newdata = new double [sumnum]; system.arraycopy(data, 0 , newdata, 0 , datalength); for ( int i = datalength; i < sumnum; i++) { newdata[i] = 0d; } } else { newdata = data; } return newdata; } /** * 去偏移量 * * @param originalarr * 原数组 * @return 目标数组 */ public static double [] deskew( double [] originalarr) { // 过滤不正确的参数 if (originalarr == null || originalarr.length <= 0 ) { return null ; } // 定义目标数组 double [] resarr = new double [originalarr.length]; // 求数组总和 double sum = 0d; for ( int i = 0 ; i < originalarr.length; i++) { sum += originalarr[i]; } // 求数组平均值 double aver = sum / originalarr.length; // 去除偏移值 for ( int i = 0 ; i < originalarr.length; i++) { resarr[i] = originalarr[i] - aver; } return resarr; } public static void main(string[] args) { // int n = integer.parseint(args[0]); double [] data = { - 0.35668879080953375 , - 0.6118094913035987 , 0.8534269560320435 , - 0.6699697478438837 , 0.35425500561437717 , 0.8910250650549392 , - 0.025718699518642918 , 0.07649691490732002 }; // 去除偏移量 data = deskew(data); // 个数为2的幂次方 data = pow2doublearr(data); int n = data.length; system.out.println(n + "数组n中数量...." ); complex[] x = new complex[n]; // original data for ( int i = 0 ; i < n; i++) { // x[i] = new complex(i, 0); // x[i] = new complex(-2 * math.random() + 1, 0); x[i] = new complex(data[i], 0 ); } show(x, "x" ); // fft of original data complex[] y = fft(x); show(y, "y = fft(x)" ); // take inverse fft complex[] z = ifft(y); show(z, "z = ifft(y)" ); // circular convolution of x with itself complex[] c = cconvolve(x, x); show(c, "c = cconvolve(x, x)" ); // linear convolution of x with itself complex[] d = convolve(x, x); show(d, "d = convolve(x, x)" ); } } /********************************************************************* * % java fft 8 x ------------------- -0.35668879080953375 -0.6118094913035987 * 0.8534269560320435 -0.6699697478438837 0.35425500561437717 0.8910250650549392 * -0.025718699518642918 0.07649691490732002 * * y = fft(x) ------------------- 0.5110172121330208 -1.245776663065442 + * 0.7113504894129803i -0.8301420417085572 - 0.8726884066879042i * -0.17611092978238008 + 2.4696418005143532i 1.1395317305034673 * -0.17611092978237974 - 2.4696418005143532i -0.8301420417085572 + * 0.8726884066879042i -1.2457766630654419 - 0.7113504894129803i * * z = ifft(y) ------------------- -0.35668879080953375 -0.6118094913035987 + * 4.2151962932466006e-17i 0.8534269560320435 - 2.691607282636124e-17i * -0.6699697478438837 + 4.1114763914420734e-17i 0.35425500561437717 * 0.8910250650549392 - 6.887033953004965e-17i -0.025718699518642918 + * 2.691607282636124e-17i 0.07649691490732002 - 1.4396387316837096e-17i * * c = cconvolve(x, x) ------------------- -1.0786973139009466 - * 2.636779683484747e-16i 1.2327819138980782 + 2.2180047699856214e-17i * 0.4386976685553382 - 1.3815636262919812e-17i -0.5579612069781844 + * 1.9986455722517509e-16i 1.432390480003344 + 2.636779683484747e-16i * -2.2165857430333684 + 2.2180047699856214e-17i -0.01255525669751989 + * 1.3815636262919812e-17i 1.0230680492494633 - 2.4422465262488753e-16i * * d = convolve(x, x) ------------------- 0.12722689348916738 + * 3.469446951953614e-17i 0.43645117531775324 - 2.78776395788635e-18i * -0.2345048043334932 - 6.907818131459906e-18i -0.5663280251946803 + * 5.829891518914417e-17i 1.2954076913348198 + 1.518836016779236e-16i * -2.212650940696159 + 1.1090023849928107e-17i -0.018407034687857718 - * 1.1306778366296569e-17i 1.023068049249463 - 9.435675069681485e-17i * -1.205924207390114 - 2.983724378680108e-16i 0.796330738580325 + * 2.4967811657742562e-17i 0.6732024728888314 - 6.907818131459906e-18i * 0.00836681821649593 + 1.4156564203603091e-16i 0.1369827886685242 + * 1.1179436667055108e-16i -0.00393480233720922 + 1.1090023849928107e-17i * 0.005851777990337828 + 2.512241462921638e-17i 1.1102230246251565e-16 - * 1.4986790192807268e-16i *********************************************************************/ |
complex.class 复数类
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
|
package fft.test; /****************************************************************************** * compilation: javac complex.java * execution: java complex * * data type for complex numbers. * * the data type is "immutable" so once you create and initialize * a complex object, you cannot change it. the "final" keyword * when declaring re and im enforces this rule, making it a * compile-time error to change the .re or .im instance variables after * they've been initialized. * * % java complex * a = 5.0 + 6.0i * b = -3.0 + 4.0i * re(a) = 5.0 * im(a) = 6.0 * b + a = 2.0 + 10.0i * a - b = 8.0 + 2.0i * a * b = -39.0 + 2.0i * b * a = -39.0 + 2.0i * a / b = 0.36 - 1.52i * (a / b) * b = 5.0 + 6.0i * conj(a) = 5.0 - 6.0i * |a| = 7.810249675906654 * tan(a) = -6.685231390246571e-6 + 1.0000103108981198i * ******************************************************************************/ import java.util.objects; public class complex { private final double re; // the real part private final double im; // the imaginary part // create a new object with the given real and imaginary parts public complex( double real, double imag) { re = real; im = imag; } // return a string representation of the invoking complex object public string tostring() { if (im == 0 ) return re + "" ; if (re == 0 ) return im + "i" ; if (im < 0 ) return re + " - " + (-im) + "i" ; return re + " + " + im + "i" ; } // return abs/modulus/magnitude public double abs() { return math.hypot(re, im); } // return angle/phase/argument, normalized to be between -pi and pi public double phase() { return math.atan2(im, re); } // return a new complex object whose value is (this + b) public complex plus(complex b) { complex a = this ; // invoking object double real = a.re + b.re; double imag = a.im + b.im; return new complex(real, imag); } // return a new complex object whose value is (this - b) public complex minus(complex b) { complex a = this ; double real = a.re - b.re; double imag = a.im - b.im; return new complex(real, imag); } // return a new complex object whose value is (this * b) public complex times(complex b) { complex a = this ; double real = a.re * b.re - a.im * b.im; double imag = a.re * b.im + a.im * b.re; return new complex(real, imag); } // return a new object whose value is (this * alpha) public complex scale( double alpha) { return new complex(alpha * re, alpha * im); } // return a new complex object whose value is the conjugate of this public complex conjugate() { return new complex(re, -im); } // return a new complex object whose value is the reciprocal of this public complex reciprocal() { double scale = re * re + im * im; return new complex(re / scale, -im / scale); } // return the real or imaginary part public double re() { return re; } public double im() { return im; } // return a / b public complex divides(complex b) { complex a = this ; return a.times(b.reciprocal()); } // return a new complex object whose value is the complex exponential of // this public complex exp() { return new complex(math.exp(re) * math.cos(im), math.exp(re) * math.sin(im)); } // return a new complex object whose value is the complex sine of this public complex sin() { return new complex(math.sin(re) * math.cosh(im), math.cos(re) * math.sinh(im)); } // return a new complex object whose value is the complex cosine of this public complex cos() { return new complex(math.cos(re) * math.cosh(im), -math.sin(re) * math.sinh(im)); } // return a new complex object whose value is the complex tangent of this public complex tan() { return sin().divides(cos()); } // a static version of plus public static complex plus(complex a, complex b) { double real = a.re + b.re; double imag = a.im + b.im; complex sum = new complex(real, imag); return sum; } // see section 3.3. public boolean equals(object x) { if (x == null ) return false ; if ( this .getclass() != x.getclass()) return false ; complex that = (complex) x; return ( this .re == that.re) && ( this .im == that.im); } // see section 3.3. public int hashcode() { return objects.hash(re, im); } // sample client for testing public static void main(string[] args) { complex a = new complex( 3.0 , 4.0 ); complex b = new complex(- 3.0 , 4.0 ); system.out.println( "a = " + a); system.out.println( "b = " + b); system.out.println( "re(a) = " + a.re()); system.out.println( "im(a) = " + a.im()); system.out.println( "b + a = " + b.plus(a)); system.out.println( "a - b = " + a.minus(b)); system.out.println( "a * b = " + a.times(b)); system.out.println( "b * a = " + b.times(a)); system.out.println( "a / b = " + a.divides(b)); system.out.println( "(a / b) * b = " + a.divides(b).times(b)); system.out.println( "conj(a) = " + a.conjugate()); system.out.println( "|a| = " + a.abs()); system.out.println( "tan(a) = " + a.tan()); } } |
希望本文所述对大家java程序设计有所帮助。
原文链接:https://blog.csdn.net/ffj0721/article/details/78521821