1"""
2Discrete Fourier Transforms - basic.py
3"""
4# Created by Pearu Peterson, August,September 2002
5__all__ = ['fft','ifft','fftn','ifftn','rfft','irfft',
6           'fft2','ifft2']
7
8from scipy.fft import _pocketfft
9from .helper import _good_shape
10
11
12def fft(x, n=None, axis=-1, overwrite_x=False):
13    """
14    Return discrete Fourier transform of real or complex sequence.
15
16    The returned complex array contains ``y(0), y(1),..., y(n-1)``, where
17
18    ``y(j) = (x * exp(-2*pi*sqrt(-1)*j*np.arange(n)/n)).sum()``.
19
20    Parameters
21    ----------
22    x : array_like
23        Array to Fourier transform.
24    n : int, optional
25        Length of the Fourier transform. If ``n < x.shape[axis]``, `x` is
26        truncated. If ``n > x.shape[axis]``, `x` is zero-padded. The
27        default results in ``n = x.shape[axis]``.
28    axis : int, optional
29        Axis along which the fft's are computed; the default is over the
30        last axis (i.e., ``axis=-1``).
31    overwrite_x : bool, optional
32        If True, the contents of `x` can be destroyed; the default is False.
33
34    Returns
35    -------
36    z : complex ndarray
37        with the elements::
38
39            [y(0),y(1),..,y(n/2),y(1-n/2),...,y(-1)]        if n is even
40            [y(0),y(1),..,y((n-1)/2),y(-(n-1)/2),...,y(-1)]  if n is odd
41
42        where::
43
44            y(j) = sum[k=0..n-1] x[k] * exp(-sqrt(-1)*j*k* 2*pi/n), j = 0..n-1
45
46    See Also
47    --------
48    ifft : Inverse FFT
49    rfft : FFT of a real sequence
50
51    Notes
52    -----
53    The packing of the result is "standard": If ``A = fft(a, n)``, then
54    ``A[0]`` contains the zero-frequency term, ``A[1:n/2]`` contains the
55    positive-frequency terms, and ``A[n/2:]`` contains the negative-frequency
56    terms, in order of decreasingly negative frequency. So ,for an 8-point
57    transform, the frequencies of the result are [0, 1, 2, 3, -4, -3, -2, -1].
58    To rearrange the fft output so that the zero-frequency component is
59    centered, like [-4, -3, -2, -1,  0,  1,  2,  3], use `fftshift`.
60
61    Both single and double precision routines are implemented. Half precision
62    inputs will be converted to single precision. Non-floating-point inputs
63    will be converted to double precision. Long-double precision inputs are
64    not supported.
65
66    This function is most efficient when `n` is a power of two, and least
67    efficient when `n` is prime.
68
69    Note that if ``x`` is real-valued, then ``A[j] == A[n-j].conjugate()``.
70    If ``x`` is real-valued and ``n`` is even, then ``A[n/2]`` is real.
71
72    If the data type of `x` is real, a "real FFT" algorithm is automatically
73    used, which roughly halves the computation time. To increase efficiency
74    a little further, use `rfft`, which does the same calculation, but only
75    outputs half of the symmetrical spectrum. If the data is both real and
76    symmetrical, the `dct` can again double the efficiency by generating
77    half of the spectrum from half of the signal.
78
79    Examples
80    --------
81    >>> from scipy.fftpack import fft, ifft
82    >>> x = np.arange(5)
83    >>> np.allclose(fft(ifft(x)), x, atol=1e-15)  # within numerical accuracy.
84    True
85
86    """
87    return _pocketfft.fft(x, n, axis, None, overwrite_x)
88
89
90def ifft(x, n=None, axis=-1, overwrite_x=False):
91    """
92    Return discrete inverse Fourier transform of real or complex sequence.
93
94    The returned complex array contains ``y(0), y(1),..., y(n-1)``, where
95
96    ``y(j) = (x * exp(2*pi*sqrt(-1)*j*np.arange(n)/n)).mean()``.
97
98    Parameters
99    ----------
100    x : array_like
101        Transformed data to invert.
102    n : int, optional
103        Length of the inverse Fourier transform.  If ``n < x.shape[axis]``,
104        `x` is truncated. If ``n > x.shape[axis]``, `x` is zero-padded.
105        The default results in ``n = x.shape[axis]``.
106    axis : int, optional
107        Axis along which the ifft's are computed; the default is over the
108        last axis (i.e., ``axis=-1``).
109    overwrite_x : bool, optional
110        If True, the contents of `x` can be destroyed; the default is False.
111
112    Returns
113    -------
114    ifft : ndarray of floats
115        The inverse discrete Fourier transform.
116
117    See Also
118    --------
119    fft : Forward FFT
120
121    Notes
122    -----
123    Both single and double precision routines are implemented. Half precision
124    inputs will be converted to single precision. Non-floating-point inputs
125    will be converted to double precision. Long-double precision inputs are
126    not supported.
127
128    This function is most efficient when `n` is a power of two, and least
129    efficient when `n` is prime.
130
131    If the data type of `x` is real, a "real IFFT" algorithm is automatically
132    used, which roughly halves the computation time.
133
134    Examples
135    --------
136    >>> from scipy.fftpack import fft, ifft
137    >>> import numpy as np
138    >>> x = np.arange(5)
139    >>> np.allclose(ifft(fft(x)), x, atol=1e-15)  # within numerical accuracy.
140    True
141
142    """
143    return _pocketfft.ifft(x, n, axis, None, overwrite_x)
144
145
146def rfft(x, n=None, axis=-1, overwrite_x=False):
147    """
148    Discrete Fourier transform of a real sequence.
149
150    Parameters
151    ----------
152    x : array_like, real-valued
153        The data to transform.
154    n : int, optional
155        Defines the length of the Fourier transform. If `n` is not specified
156        (the default) then ``n = x.shape[axis]``. If ``n < x.shape[axis]``,
157        `x` is truncated, if ``n > x.shape[axis]``, `x` is zero-padded.
158    axis : int, optional
159        The axis along which the transform is applied. The default is the
160        last axis.
161    overwrite_x : bool, optional
162        If set to true, the contents of `x` can be overwritten. Default is
163        False.
164
165    Returns
166    -------
167    z : real ndarray
168        The returned real array contains::
169
170          [y(0),Re(y(1)),Im(y(1)),...,Re(y(n/2))]              if n is even
171          [y(0),Re(y(1)),Im(y(1)),...,Re(y(n/2)),Im(y(n/2))]   if n is odd
172
173        where::
174
175          y(j) = sum[k=0..n-1] x[k] * exp(-sqrt(-1)*j*k*2*pi/n)
176          j = 0..n-1
177
178    See Also
179    --------
180    fft, irfft, scipy.fft.rfft
181
182    Notes
183    -----
184    Within numerical accuracy, ``y == rfft(irfft(y))``.
185
186    Both single and double precision routines are implemented. Half precision
187    inputs will be converted to single precision. Non-floating-point inputs
188    will be converted to double precision. Long-double precision inputs are
189    not supported.
190
191    To get an output with a complex datatype, consider using the newer
192    function `scipy.fft.rfft`.
193
194    Examples
195    --------
196    >>> from scipy.fftpack import fft, rfft
197    >>> a = [9, -9, 1, 3]
198    >>> fft(a)
199    array([  4. +0.j,   8.+12.j,  16. +0.j,   8.-12.j])
200    >>> rfft(a)
201    array([  4.,   8.,  12.,  16.])
202
203    """
204    return _pocketfft.rfft_fftpack(x, n, axis, None, overwrite_x)
205
206
207def irfft(x, n=None, axis=-1, overwrite_x=False):
208    """
209    Return inverse discrete Fourier transform of real sequence x.
210
211    The contents of `x` are interpreted as the output of the `rfft`
212    function.
213
214    Parameters
215    ----------
216    x : array_like
217        Transformed data to invert.
218    n : int, optional
219        Length of the inverse Fourier transform.
220        If n < x.shape[axis], x is truncated.
221        If n > x.shape[axis], x is zero-padded.
222        The default results in n = x.shape[axis].
223    axis : int, optional
224        Axis along which the ifft's are computed; the default is over
225        the last axis (i.e., axis=-1).
226    overwrite_x : bool, optional
227        If True, the contents of `x` can be destroyed; the default is False.
228
229    Returns
230    -------
231    irfft : ndarray of floats
232        The inverse discrete Fourier transform.
233
234    See Also
235    --------
236    rfft, ifft, scipy.fft.irfft
237
238    Notes
239    -----
240    The returned real array contains::
241
242        [y(0),y(1),...,y(n-1)]
243
244    where for n is even::
245
246        y(j) = 1/n (sum[k=1..n/2-1] (x[2*k-1]+sqrt(-1)*x[2*k])
247                                     * exp(sqrt(-1)*j*k* 2*pi/n)
248                    + c.c. + x[0] + (-1)**(j) x[n-1])
249
250    and for n is odd::
251
252        y(j) = 1/n (sum[k=1..(n-1)/2] (x[2*k-1]+sqrt(-1)*x[2*k])
253                                     * exp(sqrt(-1)*j*k* 2*pi/n)
254                    + c.c. + x[0])
255
256    c.c. denotes complex conjugate of preceding expression.
257
258    For details on input parameters, see `rfft`.
259
260    To process (conjugate-symmetric) frequency-domain data with a complex
261    datatype, consider using the newer function `scipy.fft.irfft`.
262
263    Examples
264    --------
265    >>> from scipy.fftpack import rfft, irfft
266    >>> a = [1.0, 2.0, 3.0, 4.0, 5.0]
267    >>> irfft(a)
268    array([ 2.6       , -3.16405192,  1.24398433, -1.14955713,  1.46962473])
269    >>> irfft(rfft(a))
270    array([1., 2., 3., 4., 5.])
271
272    """
273    return _pocketfft.irfft_fftpack(x, n, axis, None, overwrite_x)
274
275
276def fftn(x, shape=None, axes=None, overwrite_x=False):
277    """
278    Return multidimensional discrete Fourier transform.
279
280    The returned array contains::
281
282      y[j_1,..,j_d] = sum[k_1=0..n_1-1, ..., k_d=0..n_d-1]
283         x[k_1,..,k_d] * prod[i=1..d] exp(-sqrt(-1)*2*pi/n_i * j_i * k_i)
284
285    where d = len(x.shape) and n = x.shape.
286
287    Parameters
288    ----------
289    x : array_like
290        The (N-D) array to transform.
291    shape : int or array_like of ints or None, optional
292        The shape of the result. If both `shape` and `axes` (see below) are
293        None, `shape` is ``x.shape``; if `shape` is None but `axes` is
294        not None, then `shape` is ``numpy.take(x.shape, axes, axis=0)``.
295        If ``shape[i] > x.shape[i]``, the ith dimension is padded with zeros.
296        If ``shape[i] < x.shape[i]``, the ith dimension is truncated to
297        length ``shape[i]``.
298        If any element of `shape` is -1, the size of the corresponding
299        dimension of `x` is used.
300    axes : int or array_like of ints or None, optional
301        The axes of `x` (`y` if `shape` is not None) along which the
302        transform is applied.
303        The default is over all axes.
304    overwrite_x : bool, optional
305        If True, the contents of `x` can be destroyed. Default is False.
306
307    Returns
308    -------
309    y : complex-valued N-D NumPy array
310        The (N-D) DFT of the input array.
311
312    See Also
313    --------
314    ifftn
315
316    Notes
317    -----
318    If ``x`` is real-valued, then
319    ``y[..., j_i, ...] == y[..., n_i-j_i, ...].conjugate()``.
320
321    Both single and double precision routines are implemented. Half precision
322    inputs will be converted to single precision. Non-floating-point inputs
323    will be converted to double precision. Long-double precision inputs are
324    not supported.
325
326    Examples
327    --------
328    >>> from scipy.fftpack import fftn, ifftn
329    >>> y = (-np.arange(16), 8 - np.arange(16), np.arange(16))
330    >>> np.allclose(y, fftn(ifftn(y)))
331    True
332
333    """
334    shape = _good_shape(x, shape, axes)
335    return _pocketfft.fftn(x, shape, axes, None, overwrite_x)
336
337
338def ifftn(x, shape=None, axes=None, overwrite_x=False):
339    """
340    Return inverse multidimensional discrete Fourier transform.
341
342    The sequence can be of an arbitrary type.
343
344    The returned array contains::
345
346      y[j_1,..,j_d] = 1/p * sum[k_1=0..n_1-1, ..., k_d=0..n_d-1]
347         x[k_1,..,k_d] * prod[i=1..d] exp(sqrt(-1)*2*pi/n_i * j_i * k_i)
348
349    where ``d = len(x.shape)``, ``n = x.shape``, and ``p = prod[i=1..d] n_i``.
350
351    For description of parameters see `fftn`.
352
353    See Also
354    --------
355    fftn : for detailed information.
356
357    Examples
358    --------
359    >>> from scipy.fftpack import fftn, ifftn
360    >>> import numpy as np
361    >>> y = (-np.arange(16), 8 - np.arange(16), np.arange(16))
362    >>> np.allclose(y, ifftn(fftn(y)))
363    True
364
365    """
366    shape = _good_shape(x, shape, axes)
367    return _pocketfft.ifftn(x, shape, axes, None, overwrite_x)
368
369
370def fft2(x, shape=None, axes=(-2,-1), overwrite_x=False):
371    """
372    2-D discrete Fourier transform.
373
374    Return the 2-D discrete Fourier transform of the 2-D argument
375    `x`.
376
377    See Also
378    --------
379    fftn : for detailed information.
380
381    Examples
382    --------
383    >>> from scipy.fftpack import fft2, ifft2
384    >>> y = np.mgrid[:5, :5][0]
385    >>> y
386    array([[0, 0, 0, 0, 0],
387           [1, 1, 1, 1, 1],
388           [2, 2, 2, 2, 2],
389           [3, 3, 3, 3, 3],
390           [4, 4, 4, 4, 4]])
391    >>> np.allclose(y, ifft2(fft2(y)))
392    True
393    """
394    return fftn(x,shape,axes,overwrite_x)
395
396
397def ifft2(x, shape=None, axes=(-2,-1), overwrite_x=False):
398    """
399    2-D discrete inverse Fourier transform of real or complex sequence.
400
401    Return inverse 2-D discrete Fourier transform of
402    arbitrary type sequence x.
403
404    See `ifft` for more information.
405
406    See also
407    --------
408    fft2, ifft
409
410    Examples
411    --------
412    >>> from scipy.fftpack import fft2, ifft2
413    >>> y = np.mgrid[:5, :5][0]
414    >>> y
415    array([[0, 0, 0, 0, 0],
416           [1, 1, 1, 1, 1],
417           [2, 2, 2, 2, 2],
418           [3, 3, 3, 3, 3],
419           [4, 4, 4, 4, 4]])
420    >>> np.allclose(y, fft2(ifft2(y)))
421    True
422
423    """
424    return ifftn(x,shape,axes,overwrite_x)
425