← Back to branch summary

~ubuntu-branches/ubuntu/maverick/commons-math/maverick

« back to all changes in this revision

Viewing changes to src/java/org/apache/commons/math/transform/FastCosineTransformer.java

  • Committer: Bazaar Package Importer
  • Author(s): Damien Raude-Morvan
  • Date: 2009-08-22 01:13:25 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20090822011325-hi4peq1ua5weguwn
Tags: 2.0-1
* New upstream release.
* Set Maintainer field to Debian Java Team
* Add myself as Uploaders
* Switch to Quilt patch system:
  - Refresh all patchs
  - Remove B-D on dpatch, Add B-D on quilt
  - Include patchsys-quilt.mk in debian/rules
* Bump Standards-Version to 3.8.3:
  - Add a README.source to describe patch system
* Maven POMs:
  - Add a Build-Depends-Indep dependency on maven-repo-helper
  - Use mh_installpom and mh_installjar to install the POM and the jar to the
    Maven repository
* Use default-jdk/jre:
  - Depends on java5-runtime-headless
  - Build-Depends on default-jdk
  - Use /usr/lib/jvm/default-java as JAVA_HOME
* Move api documentation to /usr/share/doc/libcommons-math-java/api
* Build-Depends on junit4 instead of junit

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/java/org/apache/commons/math/transform/FastCosineTransformer.java
1
 
/*
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
8
 
 *
9
 
 *      http://www.apache.org/licenses/LICENSE-2.0
10
 
 *
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.
16
 
 */
17
 
package org.apache.commons.math.transform;
18
 
 
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;
23
 
 
24
 
/**
25
 
 * Implements the <a href="http://documents.wolfram.com/v5/Add-onsLinks/
26
 
 * StandardPackages/LinearAlgebra/FourierTrig.html">Fast Cosine Transform</a>
27
 
 * for transformation of one-dimensional data sets. For reference, see
28
 
 * <b>Fast Fourier Transforms</b>, ISBN 0849371635, chapter 3.
29
 
 * <p>
30
 
 * FCT is its own inverse, up to a multiplier depending on conventions.
31
 
 * The equations are listed in the comments of the corresponding methods.</p>
32
 
 * <p>
33
 
 * Different from FFT and FST, FCT requires the length of data set to be
34
 
 * power of 2 plus one. Users should especially pay attention to the
35
 
 * function transformation on how this affects the sampling.</p>
36
 
 *
37
 
 * @version $Revision: 620312 $ $Date: 2008-02-10 12:28:59 -0700 (Sun, 10 Feb 2008) $
38
 
 * @since 1.2
39
 
 */
40
 
public class FastCosineTransformer implements Serializable {
41
 
 
42
 
    /** serializable version identifier */
43
 
    static final long serialVersionUID = -7673941545134707766L;
44
 
 
45
 
    /**
46
 
     * Construct a default transformer.
47
 
     */
48
 
    public FastCosineTransformer() {
49
 
        super();
50
 
    }
51
 
 
52
 
    /**
53
 
     * Transform the given real data set.
54
 
     * <p>
55
 
     * The formula is $ F_n = (1/2) [f_0 + (-1)^n f_N] +
56
 
     *                        \Sigma_{k=0}^{N-1} f_k \cos(\pi nk/N) $
57
 
     * </p>
58
 
     * 
59
 
     * @param f the real data array to be transformed
60
 
     * @return the real transformed array
61
 
     * @throws MathException if any math-related errors occur
62
 
     * @throws IllegalArgumentException if any parameters are invalid
63
 
     */
64
 
    public double[] transform(double f[]) throws MathException,
65
 
        IllegalArgumentException {
66
 
 
67
 
        return fct(f);
68
 
    }
69
 
 
70
 
    /**
71
 
     * Transform the given real function, sampled on the given interval.
72
 
     * <p>
73
 
     * The formula is $ F_n = (1/2) [f_0 + (-1)^n f_N] +
74
 
     *                        \Sigma_{k=0}^{N-1} f_k \cos(\pi nk/N) $
75
 
     * </p>
76
 
     * 
77
 
     * @param f the function to be sampled and transformed
78
 
     * @param min the lower bound for the interval
79
 
     * @param max the upper bound for the interval
80
 
     * @param n the number of sample points
81
 
     * @return the real transformed array
82
 
     * @throws MathException if any math-related errors occur
83
 
     * @throws IllegalArgumentException if any parameters are invalid
84
 
     */
85
 
    public double[] transform(
86
 
        UnivariateRealFunction f, double min, double max, int n)
87
 
        throws MathException, IllegalArgumentException {
88
 
 
89
 
        double data[] = FastFourierTransformer.sample(f, min, max, n);
90
 
        return fct(data);
91
 
    }
92
 
 
93
 
    /**
94
 
     * Transform the given real data set.
95
 
     * <p>
96
 
     * The formula is $ F_n = \sqrt{1/2N} [f_0 + (-1)^n f_N] +
97
 
     *                        \sqrt{2/N} \Sigma_{k=0}^{N-1} f_k \cos(\pi nk/N) $
98
 
     * </p>
99
 
     * 
100
 
     * @param f the real data array to be transformed
101
 
     * @return the real transformed array
102
 
     * @throws MathException if any math-related errors occur
103
 
     * @throws IllegalArgumentException if any parameters are invalid
104
 
     */
105
 
    public double[] transform2(double f[]) throws MathException,
106
 
        IllegalArgumentException {
107
 
 
108
 
        double scaling_coefficient = Math.sqrt(2.0 / (f.length-1));
109
 
        return FastFourierTransformer.scaleArray(fct(f), scaling_coefficient);
110
 
    }
111
 
 
112
 
    /**
113
 
     * Transform the given real function, sampled on the given interval.
114
 
     * <p>
115
 
     * The formula is $ F_n = \sqrt{1/2N} [f_0 + (-1)^n f_N] +
116
 
     *                        \sqrt{2/N} \Sigma_{k=0}^{N-1} f_k \cos(\pi nk/N) $
117
 
     *
118
 
     * </p>
119
 
     * 
120
 
     * @param f the function to be sampled and transformed
121
 
     * @param min the lower bound for the interval
122
 
     * @param max the upper bound for the interval
123
 
     * @param n the number of sample points
124
 
     * @return the real transformed array
125
 
     * @throws MathException if any math-related errors occur
126
 
     * @throws IllegalArgumentException if any parameters are invalid
127
 
     */
128
 
    public double[] transform2(
129
 
        UnivariateRealFunction f, double min, double max, int n)
130
 
        throws MathException, IllegalArgumentException {
131
 
 
132
 
        double data[] = FastFourierTransformer.sample(f, min, max, n);
133
 
        double scaling_coefficient = Math.sqrt(2.0 / (n-1));
134
 
        return FastFourierTransformer.scaleArray(fct(data), scaling_coefficient);
135
 
    }
136
 
 
137
 
    /**
138
 
     * Inversely transform the given real data set.
139
 
     * <p>
140
 
     * The formula is $ f_k = (1/N) [F_0 + (-1)^k F_N] +
141
 
     *                        (2/N) \Sigma_{n=0}^{N-1} F_n \cos(\pi nk/N) $
142
 
     * </p>
143
 
     * 
144
 
     * @param f the real data array to be inversely transformed
145
 
     * @return the real inversely transformed array
146
 
     * @throws MathException if any math-related errors occur
147
 
     * @throws IllegalArgumentException if any parameters are invalid
148
 
     */
149
 
    public double[] inversetransform(double f[]) throws MathException,
150
 
        IllegalArgumentException {
151
 
 
152
 
        double scaling_coefficient = 2.0 / (f.length - 1);
153
 
        return FastFourierTransformer.scaleArray(fct(f), scaling_coefficient);
154
 
    }
155
 
 
156
 
    /**
157
 
     * Inversely transform the given real function, sampled on the given interval.
158
 
     * <p>
159
 
     * The formula is $ f_k = (1/N) [F_0 + (-1)^k F_N] +
160
 
     *                        (2/N) \Sigma_{n=0}^{N-1} F_n \cos(\pi nk/N) $
161
 
     * </p>
162
 
     * 
163
 
     * @param f the function to be sampled and inversely transformed
164
 
     * @param min the lower bound for the interval
165
 
     * @param max the upper bound for the interval
166
 
     * @param n the number of sample points
167
 
     * @return the real inversely transformed array
168
 
     * @throws MathException if any math-related errors occur
169
 
     * @throws IllegalArgumentException if any parameters are invalid
170
 
     */
171
 
    public double[] inversetransform(
172
 
        UnivariateRealFunction f, double min, double max, int n)
173
 
        throws MathException, IllegalArgumentException {
174
 
 
175
 
        double data[] = FastFourierTransformer.sample(f, min, max, n);
176
 
        double scaling_coefficient = 2.0 / (n - 1);
177
 
        return FastFourierTransformer.scaleArray(fct(data), scaling_coefficient);
178
 
    }
179
 
 
180
 
    /**
181
 
     * Inversely transform the given real data set.
182
 
     * <p>
183
 
     * The formula is $ f_k = \sqrt{1/2N} [F_0 + (-1)^k F_N] +
184
 
     *                        \sqrt{2/N} \Sigma_{n=0}^{N-1} F_n \cos(\pi nk/N) $
185
 
     * </p>
186
 
     * 
187
 
     * @param f the real data array to be inversely transformed
188
 
     * @return the real inversely transformed array
189
 
     * @throws MathException if any math-related errors occur
190
 
     * @throws IllegalArgumentException if any parameters are invalid
191
 
     */
192
 
    public double[] inversetransform2(double f[]) throws MathException,
193
 
        IllegalArgumentException {
194
 
 
195
 
        return transform2(f);
196
 
    }
197
 
 
198
 
    /**
199
 
     * Inversely transform the given real function, sampled on the given interval.
200
 
     * <p>
201
 
     * The formula is $ f_k = \sqrt{1/2N} [F_0 + (-1)^k F_N] +
202
 
     *                        \sqrt{2/N} \Sigma_{n=0}^{N-1} F_n \cos(\pi nk/N) $
203
 
     * </p>
204
 
     * 
205
 
     * @param f the function to be sampled and inversely transformed
206
 
     * @param min the lower bound for the interval
207
 
     * @param max the upper bound for the interval
208
 
     * @param n the number of sample points
209
 
     * @return the real inversely transformed array
210
 
     * @throws MathException if any math-related errors occur
211
 
     * @throws IllegalArgumentException if any parameters are invalid
212
 
     */
213
 
    public double[] inversetransform2(
214
 
        UnivariateRealFunction f, double min, double max, int n)
215
 
        throws MathException, IllegalArgumentException {
216
 
 
217
 
        return transform2(f, min, max, n);
218
 
    }
219
 
 
220
 
    /**
221
 
     * Perform the FCT algorithm (including inverse).
222
 
     *
223
 
     * @param f the real data array to be transformed
224
 
     * @return the real transformed array
225
 
     * @throws MathException if any math-related errors occur
226
 
     * @throws IllegalArgumentException if any parameters are invalid
227
 
     */
228
 
    protected double[] fct(double f[]) throws MathException,
229
 
        IllegalArgumentException {
230
 
 
231
 
        double A, B, C, F1, x[], F[] = new double[f.length];
232
 
 
233
 
        int N = f.length - 1;
234
 
        if (!FastFourierTransformer.isPowerOf2(N)) {
235
 
            throw new IllegalArgumentException
236
 
                ("Number of samples not power of 2 plus one: " + f.length);
237
 
        }
238
 
        if (N == 1) {       // trivial case
239
 
            F[0] = 0.5 * (f[0] + f[1]);
240
 
            F[1] = 0.5 * (f[0] - f[1]);
241
 
            return F;
242
 
        }
243
 
 
244
 
        // construct a new array and perform FFT on it
245
 
        x = new double[N];
246
 
        x[0] = 0.5 * (f[0] + f[N]);
247
 
        x[N >> 1] = f[N >> 1];
248
 
        F1 = 0.5 * (f[0] - f[N]);   // temporary variable for F[1]
249
 
        for (int i = 1; i < (N >> 1); i++) {
250
 
            A = 0.5 * (f[i] + f[N-i]);
251
 
            B = Math.sin(i * Math.PI / N) * (f[i] - f[N-i]);
252
 
            C = Math.cos(i * Math.PI / N) * (f[i] - f[N-i]);
253
 
            x[i] = A - B;
254
 
            x[N-i] = A + B;
255
 
            F1 += C;
256
 
        }
257
 
        FastFourierTransformer transformer = new FastFourierTransformer();
258
 
        Complex y[] = transformer.transform(x);
259
 
 
260
 
        // reconstruct the FCT result for the original array
261
 
        F[0] = y[0].getReal();
262
 
        F[1] = F1;
263
 
        for (int i = 1; i < (N >> 1); i++) {
264
 
            F[2*i] = y[i].getReal();
265
 
            F[2*i+1] = F[2*i-1] - y[i].getImaginary();
266
 
        }
267
 
        F[N] = y[N >> 1].getReal();
268
 
 
269
 
        return F;
270
 
    }
271
 
}