2
* Licensed to the Apache Software Foundation (ASF) under one or more
3
* contributor license agreements. See the NOTICE file distributed with
4
* this work for additional information regarding copyright ownership.
5
* The ASF licenses this file to You under the Apache License, Version 2.0
6
* (the "License"); you may not use this file except in compliance with
7
* the License. You may obtain a copy of the License at
9
* http://www.apache.org/licenses/LICENSE-2.0
11
* Unless required by applicable law or agreed to in writing, software
12
* distributed under the License is distributed on an "AS IS" BASIS,
13
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14
* See the License for the specific language governing permissions and
15
* limitations under the License.
17
package org.apache.commons.math.transform;
19
import java.io.Serializable;
20
import org.apache.commons.math.analysis.*;
21
import org.apache.commons.math.complex.*;
22
import org.apache.commons.math.MathException;
25
* Implements the <a href="http://documents.wolfram.com/v5/Add-onsLinks/
26
* StandardPackages/LinearAlgebra/FourierTrig.html">Fast Sine Transform</a>
27
* for transformation of one-dimensional data sets. For reference, see
28
* <b>Fast Fourier Transforms</b>, ISBN 0849371635, chapter 3.
30
* FST is its own inverse, up to a multiplier depending on conventions.
31
* The equations are listed in the comments of the corresponding methods.</p>
33
* Similar to FFT, we also require the length of data set to be power of 2.
34
* In addition, the first element must be 0 and it's enforced in function
35
* transformation after sampling.</p>
37
* @version $Revision: 620312 $ $Date: 2008-02-10 12:28:59 -0700 (Sun, 10 Feb 2008) $
40
public class FastSineTransformer implements Serializable {
42
/** serializable version identifier */
43
static final long serialVersionUID = -478002039949390854L;
46
* Construct a default transformer.
48
public FastSineTransformer() {
53
* Transform the given real data set.
55
* The formula is $ F_n = \Sigma_{k=0}^{N-1} f_k \sin(\pi nk/N) $
58
* @param f the real data array to be transformed
59
* @return the real transformed array
60
* @throws MathException if any math-related errors occur
61
* @throws IllegalArgumentException if any parameters are invalid
63
public double[] transform(double f[]) throws MathException,
64
IllegalArgumentException {
70
* Transform the given real function, sampled on the given interval.
72
* The formula is $ F_n = \Sigma_{k=0}^{N-1} f_k \sin(\pi nk/N) $
75
* @param f the function to be sampled and transformed
76
* @param min the lower bound for the interval
77
* @param max the upper bound for the interval
78
* @param n the number of sample points
79
* @return the real transformed array
80
* @throws MathException if any math-related errors occur
81
* @throws IllegalArgumentException if any parameters are invalid
83
public double[] transform(
84
UnivariateRealFunction f, double min, double max, int n)
85
throws MathException, IllegalArgumentException {
87
double data[] = FastFourierTransformer.sample(f, min, max, n);
93
* Transform the given real data set.
95
* The formula is $ F_n = \sqrt{2/N} \Sigma_{k=0}^{N-1} f_k \sin(\pi nk/N) $
98
* @param f the real data array to be transformed
99
* @return the real transformed array
100
* @throws MathException if any math-related errors occur
101
* @throws IllegalArgumentException if any parameters are invalid
103
public double[] transform2(double f[]) throws MathException,
104
IllegalArgumentException {
106
double scaling_coefficient = Math.sqrt(2.0 / f.length);
107
return FastFourierTransformer.scaleArray(fst(f), scaling_coefficient);
111
* Transform the given real function, sampled on the given interval.
113
* The formula is $ F_n = \sqrt{2/N} \Sigma_{k=0}^{N-1} f_k \sin(\pi nk/N) $
116
* @param f the function to be sampled and transformed
117
* @param min the lower bound for the interval
118
* @param max the upper bound for the interval
119
* @param n the number of sample points
120
* @return the real transformed array
121
* @throws MathException if any math-related errors occur
122
* @throws IllegalArgumentException if any parameters are invalid
124
public double[] transform2(
125
UnivariateRealFunction f, double min, double max, int n)
126
throws MathException, IllegalArgumentException {
128
double data[] = FastFourierTransformer.sample(f, min, max, n);
130
double scaling_coefficient = Math.sqrt(2.0 / n);
131
return FastFourierTransformer.scaleArray(fst(data), scaling_coefficient);
135
* Inversely transform the given real data set.
137
* The formula is $ f_k = (2/N) \Sigma_{n=0}^{N-1} F_n \sin(\pi nk/N) $
140
* @param f the real data array to be inversely transformed
141
* @return the real inversely transformed array
142
* @throws MathException if any math-related errors occur
143
* @throws IllegalArgumentException if any parameters are invalid
145
public double[] inversetransform(double f[]) throws MathException,
146
IllegalArgumentException {
148
double scaling_coefficient = 2.0 / f.length;
149
return FastFourierTransformer.scaleArray(fst(f), scaling_coefficient);
153
* Inversely transform the given real function, sampled on the given interval.
155
* The formula is $ f_k = (2/N) \Sigma_{n=0}^{N-1} F_n \sin(\pi nk/N) $
158
* @param f the function to be sampled and inversely transformed
159
* @param min the lower bound for the interval
160
* @param max the upper bound for the interval
161
* @param n the number of sample points
162
* @return the real inversely transformed array
163
* @throws MathException if any math-related errors occur
164
* @throws IllegalArgumentException if any parameters are invalid
166
public double[] inversetransform(
167
UnivariateRealFunction f, double min, double max, int n)
168
throws MathException, IllegalArgumentException {
170
double data[] = FastFourierTransformer.sample(f, min, max, n);
172
double scaling_coefficient = 2.0 / n;
173
return FastFourierTransformer.scaleArray(fst(data), scaling_coefficient);
177
* Inversely transform the given real data set.
179
* The formula is $ f_k = \sqrt{2/N} \Sigma_{n=0}^{N-1} F_n \sin(\pi nk/N) $
182
* @param f the real data array to be inversely transformed
183
* @return the real inversely transformed array
184
* @throws MathException if any math-related errors occur
185
* @throws IllegalArgumentException if any parameters are invalid
187
public double[] inversetransform2(double f[]) throws MathException,
188
IllegalArgumentException {
190
return transform2(f);
194
* Inversely transform the given real function, sampled on the given interval.
196
* The formula is $ f_k = \sqrt{2/N} \Sigma_{n=0}^{N-1} F_n \sin(\pi nk/N) $
199
* @param f the function to be sampled and inversely transformed
200
* @param min the lower bound for the interval
201
* @param max the upper bound for the interval
202
* @param n the number of sample points
203
* @return the real inversely transformed array
204
* @throws MathException if any math-related errors occur
205
* @throws IllegalArgumentException if any parameters are invalid
207
public double[] inversetransform2(
208
UnivariateRealFunction f, double min, double max, int n)
209
throws MathException, IllegalArgumentException {
211
return transform2(f, min, max, n);
215
* Perform the FST algorithm (including inverse).
217
* @param f the real data array to be transformed
218
* @return the real transformed array
219
* @throws MathException if any math-related errors occur
220
* @throws IllegalArgumentException if any parameters are invalid
222
protected double[] fst(double f[]) throws MathException,
223
IllegalArgumentException {
225
double A, B, x[], F[] = new double[f.length];
227
FastFourierTransformer.verifyDataSet(f);
229
throw new IllegalArgumentException
230
("The first element is not zero: " + f[0]);
233
if (N == 1) { // trivial case
238
// construct a new array and perform FFT on it
241
x[N >> 1] = 2.0 * f[N >> 1];
242
for (int i = 1; i < (N >> 1); i++) {
243
A = Math.sin(i * Math.PI / N) * (f[i] + f[N-i]);
244
B = 0.5 * (f[i] - f[N-i]);
248
FastFourierTransformer transformer = new FastFourierTransformer();
249
Complex y[] = transformer.transform(x);
251
// reconstruct the FST result for the original array
253
F[1] = 0.5 * y[0].getReal();
254
for (int i = 1; i < (N >> 1); i++) {
255
F[2*i] = -y[i].getImaginary();
256
F[2*i+1] = y[i].getReal() + F[2*i-1];