001/* 002 * ***** BEGIN LICENSE BLOCK ***** Version: MPL 1.1/GPL 2.0/LGPL 2.1 003 * 004 * The contents of this file are subject to the Mozilla Public License Version 005 * 1.1 (the "License"); you may not use this file except in compliance with the 006 * License. You may obtain a copy of the License at http://www.mozilla.org/MPL/ 007 * 008 * Software distributed under the License is distributed on an "AS IS" basis, 009 * WITHOUT WARRANTY OF ANY KIND, either express or implied. See the License for 010 * the specific language governing rights and limitations under the License. 011 * 012 * The Original Code is JTransforms. 013 * 014 * The Initial Developer of the Original Code is Piotr Wendykier, Emory 015 * University. Portions created by the Initial Developer are Copyright (C) 016 * 2007-2009 the Initial Developer. All Rights Reserved. 017 * 018 * Alternatively, the contents of this file may be used under the terms of 019 * either the GNU General Public License Version 2 or later (the "GPL"), or the 020 * GNU Lesser General Public License Version 2.1 or later (the "LGPL"), in which 021 * case the provisions of the GPL or the LGPL are applicable instead of those 022 * above. If you wish to allow use of your version of this file only under the 023 * terms of either the GPL or the LGPL, and not to allow others to use your 024 * version of this file under the terms of the MPL, indicate your decision by 025 * deleting the provisions above and replace them with the notice and other 026 * provisions required by the GPL or the LGPL. If you do not delete the 027 * provisions above, a recipient may use your version of this file under the 028 * terms of any one of the MPL, the GPL or the LGPL. 029 * 030 * ***** END LICENSE BLOCK ***** 031 */ 032 033package edu.emory.mathcs.jtransforms.fft; 034 035// @formatter:off 036/** 037 * <p> 038 * This is a set of utility methods for R/W access to data resulting from a call 039 * to the Fourier transform of <em>real</em> data. Memory optimized methods, 040 * namely 041 * <ul> 042 * <li>{@link DoubleFFT_2D#realForward(double[])}</li> 043 * <li>{@link DoubleFFT_2D#realForward(double[][])}</li> 044 * <li>{@link FloatFFT_2D#realForward(float[])}</li> 045 * <li>{@link FloatFFT_2D#realForward(float[][])}</li> 046 * </ul> 047 * are implemented to handle this case specifically. However, packing of the 048 * data in the data array is somewhat obscure. This class provides methods for 049 * direct access to the data, without the burden of all necessary tests. 050 * </p> 051 * <h3>Example for Fourier Transform of real, double precision 1d data</h3> 052 * <p> 053 * <pre> 054 * DoubleFFT_2D fft = new DoubleFFT_2D(rows, columns); 055 * double[] data = new double[2 * rows * columns]; 056 * ... 057 * fft.realForwardFull(data); 058 * data[r1 * 2 * columns + c1] = val1; 059 * val2 = data[r2 * 2 * columns + c2]; 060 * </pre> 061 * is equivalent to 062 * <pre> 063 * DoubleFFT_2D fft = new DoubleFFT_2D(rows, columns); 064 * RealFFTUtils_2D unpacker = new RealFFTUtils_2D(rows, columns); 065 * double[] data = new double[rows * columns]; 066 * ... 067 * fft.realForward(data); 068 * unpacker.pack(val1, r1, c1, data); 069 * val2 = unpacker.unpack(r2, c2, data, 0); 070 * </pre> 071 * Even (resp. odd) values of <code>c</code> correspond to the real (resp. 072 * imaginary) part of the Fourier mode. 073 * </p> 074 * <h3>Example for Fourier Transform of real, double precision 2d data</h3> 075 * <p> 076 * <pre> 077 * DoubleFFT_2D fft = new DoubleFFT_2D(rows, columns); 078 * double[][] data = new double[rows][2 * columns]; 079 * ... 080 * fft.realForwardFull(data); 081 * data[r1][c1] = val1; 082 * val2 = data[r2][c2]; 083 * </pre> 084 * is equivalent to 085 * <pre> 086 * DoubleFFT_2D fft = new DoubleFFT_2D(rows, columns); 087 * RealFFTUtils_2D unpacker = new RealFFTUtils_2D(rows, columns); 088 * double[][] data = new double[rows][columns]; 089 * ... 090 * fft.realForward(data); 091 * unpacker.pack(val1, r1, c1, data); 092 * val2 = unpacker.unpack(r2, c2, data, 0); 093 * </pre> 094 * Even (resp. odd) values of <code>c</code> correspond to the real (resp. 095 * imaginary) part of the Fourier mode. 096 * </p> 097 * 098 * @author Sébastien Brisard 099 * 100 */ 101// @formatter:on 102public class RealFFTUtils_2D { 103 /** The constant <code>int</code> value of 1. */ 104 private static final int ONE = 1; 105 106 /** The constant <code>int</code> value of 2. */ 107 private static final int TWO = 2; 108 109 /** The constant <code>int</code> value of 0. */ 110 private static final int ZERO = 0; 111 112 /** The size of the data in the second direction. */ 113 private final int columns; 114 115 /** The size of the data in the first direction. */ 116 private final int rows; 117 118 /** 119 * Creates a new instance of this class. The size of the underlying 120 * {@link DoubleFFT_2D} or {@link FloatFFT_2D} must be specified. 121 * 122 * @param rows 123 * number of rows 124 * @param columns 125 * number of columns 126 */ 127 public RealFFTUtils_2D(final int rows, final int columns) { 128 this.columns = columns; 129 this.rows = rows; 130 } 131 132 /** 133 * <p> 134 * Returns the 1d index of the specified 2d Fourier mode. In other words, if 135 * <code>packed</code> contains the transformed data following a call to 136 * {@link DoubleFFT_2D#realForward(double[])} or 137 * {@link FloatFFT_2D#realForward(float[])}, then the returned value 138 * <code>index</code> gives access to the <code>[r][c]</code> Fourier mode 139 * <ul> 140 * <li>if <code>index == {@link Integer#MIN_VALUE}</code>, then the Fourier 141 * mode is zero,</li> 142 * <li>if <code>index >= 0</code>, then the Fourier mode is 143 * <code>packed[index]</code>,</li> 144 * <li>if <code>index < 0</code>, then the Fourier mode is 145 * <code>-packed[-index]</code>,</li> 146 * </ul> 147 * </p> 148 * 149 * @param r 150 * the row index 151 * @param c 152 * the column index 153 * @return the value of <code>index</code> 154 */ 155 public int getIndex(final int r, final int c) { 156 final int cmod2 = c & ONE; 157 final int rmul2 = r << ONE; 158 if (r != ZERO) { 159 if (c <= ONE) { 160 if (rmul2 == rows) { 161 if (cmod2 == ONE) { 162 return Integer.MIN_VALUE; 163 } 164 return ((rows * columns) >> ONE); 165 } else if (rmul2 < rows) { 166 return columns * r + cmod2; 167 } else { 168 if (cmod2 == ZERO) { 169 return columns * (rows - r); 170 } else { 171 return -(columns * (rows - r) + ONE); 172 } 173 } 174 } else if ((c == columns) || (c == columns + ONE)) { 175 if (rmul2 == rows) { 176 if (cmod2 == ONE) { 177 return Integer.MIN_VALUE; 178 } 179 return ((rows * columns) >> ONE) + ONE; 180 } else if (rmul2 < rows) { 181 if (cmod2 == ZERO) { 182 return columns * (rows - r) + ONE; 183 } else { 184 return -(columns * (rows - r)); 185 } 186 } else { 187 return columns * r + ONE - cmod2; 188 } 189 } else if (c < columns) { 190 return columns * r + c; 191 } else { 192 if (cmod2 == ZERO) { 193 return columns * (rows + TWO - r) - c; 194 } else { 195 return -(columns * (rows + TWO - r) - c + TWO); 196 } 197 } 198 } else { 199 if ((c == ONE) || (c == columns + ONE)) { 200 return Integer.MIN_VALUE; 201 } else if (c == columns) { 202 return ONE; 203 } else if (c < columns) { 204 return c; 205 } else { 206 if (cmod2 == ZERO) { 207 return (columns << ONE) - c; 208 } else { 209 return -((columns << ONE) - c + TWO); 210 } 211 } 212 } 213 } 214 215 /** 216 * Sets the specified Fourier mode of the transformed data. The data array 217 * results from a call to {@link DoubleFFT_2D#realForward(double[])}. 218 * 219 * @param val 220 * the new value of the <code>[r][c]</code> Fourier mode 221 * @param r 222 * the row index 223 * @param c 224 * the column index 225 * @param packed 226 * the transformed data 227 * @param pos 228 * index of the first element in array <code>packed</code> 229 */ 230 public void pack(final double val, final int r, final int c, 231 final double[] packed, final int pos) { 232 final int index = getIndex(r, c); 233 if (index >= 0) { 234 packed[pos + index] = val; 235 } else if (index > Integer.MIN_VALUE) { 236 packed[pos - index] = -val; 237 } else { 238 throw new IllegalArgumentException( 239 String.format( 240 "[%d][%d] component cannot be modified (always zero)", 241 r, c)); 242 } 243 } 244 245 /** 246 * Sets the specified Fourier mode of the transformed data. The data array 247 * results from a call to {@link DoubleFFT_2D#realForward(double[][])}. 248 * 249 * @param val 250 * the new value of the <code>[r][c]</code> Fourier mode 251 * @param r 252 * the row index 253 * @param c 254 * the column index 255 * @param packed 256 * the transformed data 257 */ 258 public void pack(final double val, final int r, final int c, 259 final double[][] packed) { 260 final int index = getIndex(r, c); 261 if (index >= 0) { 262 packed[index / columns][index % columns] = val; 263 } else if (index > Integer.MIN_VALUE) { 264 packed[(-index) / columns][(-index) % columns] = -val; 265 } else { 266 throw new IllegalArgumentException( 267 String.format( 268 "[%d][%d] component cannot be modified (always zero)", 269 r, c)); 270 } 271 } 272 273 /** 274 * Sets the specified Fourier mode of the transformed data. The data array 275 * results from a call to {@link FloatFFT_2D#realForward(float[])}. 276 * 277 * @param val 278 * the new value of the <code>[r][c]</code> Fourier mode 279 * @param r 280 * the row index 281 * @param c 282 * the column index 283 * @param packed 284 * the transformed data 285 * @param pos 286 * index of the first element in array <code>packed</code> 287 */ 288 public void pack(final float val, final int r, final int c, 289 final float[] packed, final int pos) { 290 final int index = getIndex(r, c); 291 if (index >= 0) { 292 packed[pos + index] = val; 293 } else if (index > Integer.MIN_VALUE) { 294 packed[pos - index] = -val; 295 } else { 296 throw new IllegalArgumentException( 297 String.format( 298 "[%d][%d] component cannot be modified (always zero)", 299 r, c)); 300 } 301 } 302 303 /** 304 * Sets the specified Fourier mode of the transformed data. The data array 305 * results from a call to {@link FloatFFT_2D#realForward(float[][])}. 306 * 307 * @param val 308 * the new value of the <code>[r][c]</code> Fourier mode 309 * @param r 310 * the row index 311 * @param c 312 * the column index 313 * @param packed 314 * the transformed data 315 */ 316 public void pack(final float val, final int r, final int c, 317 final float[][] packed) { 318 final int index = getIndex(r, c); 319 if (index >= 0) { 320 packed[index / columns][index % columns] = val; 321 } else if (index > Integer.MIN_VALUE) { 322 packed[(-index) / columns][(-index) % columns] = -val; 323 } else { 324 throw new IllegalArgumentException( 325 String.format( 326 "[%d][%d] component cannot be modified (always zero)", 327 r, c)); 328 } 329 } 330 331 /** 332 * Returns the specified Fourier mode of the transformed data. The data 333 * array results from a call to {@link DoubleFFT_2D#realForward(double[])}. 334 * 335 * @param r 336 * the row index 337 * @param c 338 * the column index 339 * @param packed 340 * the transformed data 341 * @param pos 342 * index of the first element in array <code>packed</code> 343 * @return the value of the <code>[r][c]</code> Fourier mode 344 */ 345 public double unpack(final int r, final int c, final double[] packed, 346 final int pos) { 347 final int index = getIndex(r, c); 348 if (index >= 0) { 349 return packed[pos + index]; 350 } else if (index > Integer.MIN_VALUE) { 351 return -packed[pos - index]; 352 } else { 353 return ZERO; 354 } 355 } 356 357 /** 358 * Returns the specified Fourier mode of the transformed data. The data 359 * array results from a call to {@link DoubleFFT_2D#realForward(double[][])} 360 * . 361 * 362 * @param r 363 * the row index 364 * @param c 365 * the column index 366 * @param packed 367 * the transformed data 368 * @return the value of the <code>[r][c]</code> Fourier mode 369 */ 370 public double unpack(final int r, final int c, final double[][] packed) { 371 final int index = getIndex(r, c); 372 if (index >= 0) { 373 return packed[index / columns][index % columns]; 374 } else if (index > Integer.MIN_VALUE) { 375 return -packed[(-index) / columns][(-index) % columns]; 376 } else { 377 return ZERO; 378 } 379 } 380 381 /** 382 * Returns the specified Fourier mode of the transformed data. The data 383 * array results from a call to {@link FloatFFT_2D#realForward(float[])} 384 * . 385 * 386 * @param r 387 * the row index 388 * @param c 389 * the column index 390 * @param packed 391 * the transformed data 392 * @param pos 393 * index of the first element in array <code>packed</code> 394 * @return the value of the <code>[r][c]</code> Fourier mode 395 */ 396 public float unpack(final int r, final int c, final float[] packed, 397 final int pos) { 398 final int index = getIndex(r, c); 399 if (index >= 0) { 400 return packed[pos + index]; 401 } else if (index > Integer.MIN_VALUE) { 402 return -packed[pos - index]; 403 } else { 404 return ZERO; 405 } 406 } 407 408 /** 409 * Returns the specified Fourier mode of the transformed data. The data 410 * array results from a call to {@link FloatFFT_2D#realForward(float[][])} . 411 * 412 * @param r 413 * the row index 414 * @param c 415 * the column index 416 * @param packed 417 * the transformed data 418 * @return the value of the <code>[r][c]</code> Fourier mode 419 */ 420 public float unpack(final int r, final int c, final float[][] packed) { 421 final int index = getIndex(r, c); 422 if (index >= 0) { 423 return packed[index / columns][index % columns]; 424 } else if (index > Integer.MIN_VALUE) { 425 return -packed[(-index) / columns][(-index) % columns]; 426 } else { 427 return ZERO; 428 } 429 } 430}