001: /*
002: * Copyright (c) 2000 Silvere Martin-Michiellot All Rights Reserved.
003: *
004: * Silvere Martin-Michiellot grants you ("Licensee") a non-exclusive,
005: * royalty free, license to use, modify and redistribute this
006: * software in source and binary code form,
007: * provided that i) this copyright notice and license appear on all copies of
008: * the software; and ii) Licensee does not utilize the software in a manner
009: * which is disparaging to Silvere Martin-Michiellot.
010: *
011: * This software is provided "AS IS," without a warranty of any kind. ALL
012: * EXPRESS OR IMPLIED CONDITIONS, REPRESENTATIONS AND WARRANTIES, INCLUDING ANY
013: * IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE OR
014: * NON-INFRINGEMENT, ARE HEREBY EXCLUDED. Silvere Martin-Michiellot
015: * AND ITS LICENSORS SHALL NOT BE LIABLE FOR ANY DAMAGES
016: * SUFFERED BY LICENSEE AS A RESULT OF USING, MODIFYING
017: * OR DISTRIBUTING THE SOFTWARE OR ITS DERIVATIVES. IN NO EVENT WILL
018: * Silvere Martin-Michiellot OR ITS LICENSORS BE LIABLE
019: * FOR ANY LOST REVENUE, PROFIT OR DATA, OR FOR DIRECT,
020: * INDIRECT, SPECIAL, CONSEQUENTIAL, INCIDENTAL OR PUNITIVE DAMAGES, HOWEVER
021: * CAUSED AND REGARDLESS OF THE THEORY OF LIABILITY, ARISING OUT OF THE USE OF
022: * OR INABILITY TO USE SOFTWARE, EVEN IF Silvere Martin-Michiellot HAS BEEN
023: * ADVISED OF THE POSSIBILITY OF SUCH DAMAGES.
024: *
025: * This software is not designed or intended for use in on-line control of
026: * aircraft, air traffic, aircraft navigation or aircraft communications; or in
027: * the design, construction, operation or maintenance of any nuclear
028: * facility. Licensee represents and warrants that it will not use or
029: * redistribute the Software for such purposes.
030: *
031: */
032:
033: // This code is repackaged after the code from Craig A. Lindley, from Digital Audio with Java
034: // Site ftp://ftp.prenhall.com/pub/ptr/professional_computer_science.w-022/digital_audio/
035: // Email
036: package com.db.media.audio.dsp.monitors;
037:
038: /**
039: * This is a Java implementation of the fast Fourier transform
040: * written by Jef Poskanzer. The copyright appears above.
041: */
042:
043: /* libfft.c - fast Fourier transform library
044: **
045: ** Copyright (C) 1989 by Jef Poskanzer.
046: **
047: ** Permission to use, copy, modify, and distribute this software and its
048: ** documentation for any purpose and without fee is hereby granted, provided
049: ** that the above copyright notice appear in all copies and that both that
050: ** copyright notice and this permission notice appear in supporting
051: ** documentation. This software is provided "as is" without express or
052: ** implied warranty.
053: */
054:
055: public class FFT {
056:
057: private static final double TWOPI = 2.0 * Math.PI;
058:
059: // Limits on the number of bits this algorithm can utilize
060: private static final int LOG2_MAXFFTSIZE = 15;
061: private static final int MAXFFTSIZE = 1 << LOG2_MAXFFTSIZE;
062:
063: // Private class data
064: private int bits;
065: private int[] bitreverse = new int[MAXFFTSIZE];
066:
067: /**
068: * FFT class constructor
069: * Initializes code for doing a fast Fourier transform
070: *
071: * @param int bits is a power of two such that 2^b is the number
072: * of samples.
073: */
074: public FFT(int bits) {
075:
076: this .bits = bits;
077:
078: if (bits > LOG2_MAXFFTSIZE) {
079: System.out.println("" + bits + " is too big");
080: System.exit(1);
081: }
082: for (int i = (1 << bits) - 1; i >= 0; --i) {
083: int k = 0;
084: for (int j = 0; j < bits; ++j) {
085: k *= 2;
086: if ((i & (1 << j)) != 0)
087: k++;
088: }
089: bitreverse[i] = k;
090: }
091:
092: }
093:
094: /**
095: * A fast Fourier transform routine
096: *
097: * @param double [] xr is the real part of the data to be transformed
098: * @param double [] xi is the imaginary part of the data to be transformed
099: * (normally zero unless inverse transofrm is effect).
100: * @param boolean invFlag which is true if inverse transform is being
101: * applied. false for a forward transform.
102: */
103: public void doFFT(double[] xr, double[] xi, boolean invFlag) {
104:
105: int n, n2, i, k, kn2, l, p;
106: double ang, s, c, tr, ti;
107:
108: n2 = (n = (1 << bits)) / 2;
109:
110: for (l = 0; l < bits; ++l) {
111: for (k = 0; k < n; k += n2) {
112: for (i = 0; i < n2; ++i, ++k) {
113: p = bitreverse[k / n2];
114: ang = TWOPI * p / n;
115: c = Math.cos(ang);
116: s = Math.sin(ang);
117: kn2 = k + n2;
118:
119: if (invFlag)
120: s = -s;
121:
122: tr = xr[kn2] * c + xi[kn2] * s;
123: ti = xi[kn2] * c - xr[kn2] * s;
124:
125: xr[kn2] = xr[k] - tr;
126: xi[kn2] = xi[k] - ti;
127: xr[k] += tr;
128: xi[k] += ti;
129: }
130: }
131: n2 /= 2;
132: }
133:
134: for (k = 0; k < n; k++) {
135: if ((i = bitreverse[k]) <= k)
136: continue;
137:
138: tr = xr[k];
139: ti = xi[k];
140: xr[k] = xr[i];
141: xi[k] = xi[i];
142: xr[i] = tr;
143: xi[i] = ti;
144: }
145:
146: // Finally, multiply each value by 1/n, if this is the forward
147: // transform.
148: if (!invFlag) {
149: double f = 1.0 / n;
150:
151: for (i = 0; i < n; i++) {
152: xr[i] *= f;
153: xi[i] *= f;
154: }
155: }
156:
157: }
158:
159: }
|