Hide keyboard shortcuts

Hot-keys on this page

r m x p   toggle line displays

j k   next/prev highlighted chunk

0   (zero) top of page

1   (one) first highlighted chunk

1""" 

2Discrete Fourier Transforms 

3 

4Routines in this module: 

5 

6fft(a, n=None, axis=-1) 

7ifft(a, n=None, axis=-1) 

8rfft(a, n=None, axis=-1) 

9irfft(a, n=None, axis=-1) 

10hfft(a, n=None, axis=-1) 

11ihfft(a, n=None, axis=-1) 

12fftn(a, s=None, axes=None) 

13ifftn(a, s=None, axes=None) 

14rfftn(a, s=None, axes=None) 

15irfftn(a, s=None, axes=None) 

16fft2(a, s=None, axes=(-2,-1)) 

17ifft2(a, s=None, axes=(-2, -1)) 

18rfft2(a, s=None, axes=(-2,-1)) 

19irfft2(a, s=None, axes=(-2, -1)) 

20 

21i = inverse transform 

22r = transform of purely real data 

23h = Hermite transform 

24n = n-dimensional transform 

252 = 2-dimensional transform 

26(Note: 2D routines are just nD routines with different default 

27behavior.) 

28 

29""" 

30__all__ = ['fft', 'ifft', 'rfft', 'irfft', 'hfft', 'ihfft', 'rfftn', 

31 'irfftn', 'rfft2', 'irfft2', 'fft2', 'ifft2', 'fftn', 'ifftn'] 

32 

33import functools 

34 

35from numpy.core import asarray, zeros, swapaxes, conjugate, take, sqrt 

36from . import _pocketfft_internal as pfi 

37from numpy.core.multiarray import normalize_axis_index 

38from numpy.core import overrides 

39 

40 

41array_function_dispatch = functools.partial( 

42 overrides.array_function_dispatch, module='numpy.fft') 

43 

44 

45# `inv_norm` is a float by which the result of the transform needs to be 

46# divided. This replaces the original, more intuitive 'fct` parameter to avoid 

47# divisions by zero (or alternatively additional checks) in the case of 

48# zero-length axes during its computation. 

49def _raw_fft(a, n, axis, is_real, is_forward, inv_norm): 

50 axis = normalize_axis_index(axis, a.ndim) 

51 if n is None: 

52 n = a.shape[axis] 

53 

54 if n < 1: 

55 raise ValueError("Invalid number of FFT data points (%d) specified." 

56 % n) 

57 

58 fct = 1/inv_norm 

59 

60 if a.shape[axis] != n: 

61 s = list(a.shape) 

62 index = [slice(None)]*len(s) 

63 if s[axis] > n: 

64 index[axis] = slice(0, n) 

65 a = a[tuple(index)] 

66 else: 

67 index[axis] = slice(0, s[axis]) 

68 s[axis] = n 

69 z = zeros(s, a.dtype.char) 

70 z[tuple(index)] = a 

71 a = z 

72 

73 if axis == a.ndim-1: 

74 r = pfi.execute(a, is_real, is_forward, fct) 

75 else: 

76 a = swapaxes(a, axis, -1) 

77 r = pfi.execute(a, is_real, is_forward, fct) 

78 r = swapaxes(r, axis, -1) 

79 return r 

80 

81 

82def _unitary(norm): 

83 if norm is None: 

84 return False 

85 if norm=="ortho": 

86 return True 

87 raise ValueError("Invalid norm value %s, should be None or \"ortho\"." 

88 % norm) 

89 

90 

91def _fft_dispatcher(a, n=None, axis=None, norm=None): 

92 return (a,) 

93 

94 

95@array_function_dispatch(_fft_dispatcher) 

96def fft(a, n=None, axis=-1, norm=None): 

97 """ 

98 Compute the one-dimensional discrete Fourier Transform. 

99 

100 This function computes the one-dimensional *n*-point discrete Fourier 

101 Transform (DFT) with the efficient Fast Fourier Transform (FFT) 

102 algorithm [CT]. 

103 

104 Parameters 

105 ---------- 

106 a : array_like 

107 Input array, can be complex. 

108 n : int, optional 

109 Length of the transformed axis of the output. 

110 If `n` is smaller than the length of the input, the input is cropped. 

111 If it is larger, the input is padded with zeros. If `n` is not given, 

112 the length of the input along the axis specified by `axis` is used. 

113 axis : int, optional 

114 Axis over which to compute the FFT. If not given, the last axis is 

115 used. 

116 norm : {None, "ortho"}, optional 

117 .. versionadded:: 1.10.0 

118 

119 Normalization mode (see `numpy.fft`). Default is None. 

120 

121 Returns 

122 ------- 

123 out : complex ndarray 

124 The truncated or zero-padded input, transformed along the axis 

125 indicated by `axis`, or the last one if `axis` is not specified. 

126 

127 Raises 

128 ------ 

129 IndexError 

130 if `axes` is larger than the last axis of `a`. 

131 

132 See Also 

133 -------- 

134 numpy.fft : for definition of the DFT and conventions used. 

135 ifft : The inverse of `fft`. 

136 fft2 : The two-dimensional FFT. 

137 fftn : The *n*-dimensional FFT. 

138 rfftn : The *n*-dimensional FFT of real input. 

139 fftfreq : Frequency bins for given FFT parameters. 

140 

141 Notes 

142 ----- 

143 FFT (Fast Fourier Transform) refers to a way the discrete Fourier 

144 Transform (DFT) can be calculated efficiently, by using symmetries in the 

145 calculated terms. The symmetry is highest when `n` is a power of 2, and 

146 the transform is therefore most efficient for these sizes. 

147 

148 The DFT is defined, with the conventions used in this implementation, in 

149 the documentation for the `numpy.fft` module. 

150 

151 References 

152 ---------- 

153 .. [CT] Cooley, James W., and John W. Tukey, 1965, "An algorithm for the 

154 machine calculation of complex Fourier series," *Math. Comput.* 

155 19: 297-301. 

156 

157 Examples 

158 -------- 

159 >>> np.fft.fft(np.exp(2j * np.pi * np.arange(8) / 8)) 

160 array([-2.33486982e-16+1.14423775e-17j, 8.00000000e+00-1.25557246e-15j, 

161 2.33486982e-16+2.33486982e-16j, 0.00000000e+00+1.22464680e-16j, 

162 -1.14423775e-17+2.33486982e-16j, 0.00000000e+00+5.20784380e-16j, 

163 1.14423775e-17+1.14423775e-17j, 0.00000000e+00+1.22464680e-16j]) 

164 

165 In this example, real input has an FFT which is Hermitian, i.e., symmetric 

166 in the real part and anti-symmetric in the imaginary part, as described in 

167 the `numpy.fft` documentation: 

168 

169 >>> import matplotlib.pyplot as plt 

170 >>> t = np.arange(256) 

171 >>> sp = np.fft.fft(np.sin(t)) 

172 >>> freq = np.fft.fftfreq(t.shape[-1]) 

173 >>> plt.plot(freq, sp.real, freq, sp.imag) 

174 [<matplotlib.lines.Line2D object at 0x...>, <matplotlib.lines.Line2D object at 0x...>] 

175 >>> plt.show() 

176 

177 """ 

178 

179 a = asarray(a) 

180 if n is None: 

181 n = a.shape[axis] 

182 inv_norm = 1 

183 if norm is not None and _unitary(norm): 

184 inv_norm = sqrt(n) 

185 output = _raw_fft(a, n, axis, False, True, inv_norm) 

186 return output 

187 

188 

189@array_function_dispatch(_fft_dispatcher) 

190def ifft(a, n=None, axis=-1, norm=None): 

191 """ 

192 Compute the one-dimensional inverse discrete Fourier Transform. 

193 

194 This function computes the inverse of the one-dimensional *n*-point 

195 discrete Fourier transform computed by `fft`. In other words, 

196 ``ifft(fft(a)) == a`` to within numerical accuracy. 

197 For a general description of the algorithm and definitions, 

198 see `numpy.fft`. 

199 

200 The input should be ordered in the same way as is returned by `fft`, 

201 i.e., 

202 

203 * ``a[0]`` should contain the zero frequency term, 

204 * ``a[1:n//2]`` should contain the positive-frequency terms, 

205 * ``a[n//2 + 1:]`` should contain the negative-frequency terms, in 

206 increasing order starting from the most negative frequency. 

207 

208 For an even number of input points, ``A[n//2]`` represents the sum of 

209 the values at the positive and negative Nyquist frequencies, as the two 

210 are aliased together. See `numpy.fft` for details. 

211 

212 Parameters 

213 ---------- 

214 a : array_like 

215 Input array, can be complex. 

216 n : int, optional 

217 Length of the transformed axis of the output. 

218 If `n` is smaller than the length of the input, the input is cropped. 

219 If it is larger, the input is padded with zeros. If `n` is not given, 

220 the length of the input along the axis specified by `axis` is used. 

221 See notes about padding issues. 

222 axis : int, optional 

223 Axis over which to compute the inverse DFT. If not given, the last 

224 axis is used. 

225 norm : {None, "ortho"}, optional 

226 .. versionadded:: 1.10.0 

227 

228 Normalization mode (see `numpy.fft`). Default is None. 

229 

230 Returns 

231 ------- 

232 out : complex ndarray 

233 The truncated or zero-padded input, transformed along the axis 

234 indicated by `axis`, or the last one if `axis` is not specified. 

235 

236 Raises 

237 ------ 

238 IndexError 

239 If `axes` is larger than the last axis of `a`. 

240 

241 See Also 

242 -------- 

243 numpy.fft : An introduction, with definitions and general explanations. 

244 fft : The one-dimensional (forward) FFT, of which `ifft` is the inverse 

245 ifft2 : The two-dimensional inverse FFT. 

246 ifftn : The n-dimensional inverse FFT. 

247 

248 Notes 

249 ----- 

250 If the input parameter `n` is larger than the size of the input, the input 

251 is padded by appending zeros at the end. Even though this is the common 

252 approach, it might lead to surprising results. If a different padding is 

253 desired, it must be performed before calling `ifft`. 

254 

255 Examples 

256 -------- 

257 >>> np.fft.ifft([0, 4, 0, 0]) 

258 array([ 1.+0.j, 0.+1.j, -1.+0.j, 0.-1.j]) # may vary 

259 

260 Create and plot a band-limited signal with random phases: 

261 

262 >>> import matplotlib.pyplot as plt 

263 >>> t = np.arange(400) 

264 >>> n = np.zeros((400,), dtype=complex) 

265 >>> n[40:60] = np.exp(1j*np.random.uniform(0, 2*np.pi, (20,))) 

266 >>> s = np.fft.ifft(n) 

267 >>> plt.plot(t, s.real, 'b-', t, s.imag, 'r--') 

268 [<matplotlib.lines.Line2D object at ...>, <matplotlib.lines.Line2D object at ...>] 

269 >>> plt.legend(('real', 'imaginary')) 

270 <matplotlib.legend.Legend object at ...> 

271 >>> plt.show() 

272 

273 """ 

274 a = asarray(a) 

275 if n is None: 

276 n = a.shape[axis] 

277 if norm is not None and _unitary(norm): 

278 inv_norm = sqrt(max(n, 1)) 

279 else: 

280 inv_norm = n 

281 output = _raw_fft(a, n, axis, False, False, inv_norm) 

282 return output 

283 

284 

285 

286@array_function_dispatch(_fft_dispatcher) 

287def rfft(a, n=None, axis=-1, norm=None): 

288 """ 

289 Compute the one-dimensional discrete Fourier Transform for real input. 

290 

291 This function computes the one-dimensional *n*-point discrete Fourier 

292 Transform (DFT) of a real-valued array by means of an efficient algorithm 

293 called the Fast Fourier Transform (FFT). 

294 

295 Parameters 

296 ---------- 

297 a : array_like 

298 Input array 

299 n : int, optional 

300 Number of points along transformation axis in the input to use. 

301 If `n` is smaller than the length of the input, the input is cropped. 

302 If it is larger, the input is padded with zeros. If `n` is not given, 

303 the length of the input along the axis specified by `axis` is used. 

304 axis : int, optional 

305 Axis over which to compute the FFT. If not given, the last axis is 

306 used. 

307 norm : {None, "ortho"}, optional 

308 .. versionadded:: 1.10.0 

309 

310 Normalization mode (see `numpy.fft`). Default is None. 

311 

312 Returns 

313 ------- 

314 out : complex ndarray 

315 The truncated or zero-padded input, transformed along the axis 

316 indicated by `axis`, or the last one if `axis` is not specified. 

317 If `n` is even, the length of the transformed axis is ``(n/2)+1``. 

318 If `n` is odd, the length is ``(n+1)/2``. 

319 

320 Raises 

321 ------ 

322 IndexError 

323 If `axis` is larger than the last axis of `a`. 

324 

325 See Also 

326 -------- 

327 numpy.fft : For definition of the DFT and conventions used. 

328 irfft : The inverse of `rfft`. 

329 fft : The one-dimensional FFT of general (complex) input. 

330 fftn : The *n*-dimensional FFT. 

331 rfftn : The *n*-dimensional FFT of real input. 

332 

333 Notes 

334 ----- 

335 When the DFT is computed for purely real input, the output is 

336 Hermitian-symmetric, i.e. the negative frequency terms are just the complex 

337 conjugates of the corresponding positive-frequency terms, and the 

338 negative-frequency terms are therefore redundant. This function does not 

339 compute the negative frequency terms, and the length of the transformed 

340 axis of the output is therefore ``n//2 + 1``. 

341 

342 When ``A = rfft(a)`` and fs is the sampling frequency, ``A[0]`` contains 

343 the zero-frequency term 0*fs, which is real due to Hermitian symmetry. 

344 

345 If `n` is even, ``A[-1]`` contains the term representing both positive 

346 and negative Nyquist frequency (+fs/2 and -fs/2), and must also be purely 

347 real. If `n` is odd, there is no term at fs/2; ``A[-1]`` contains 

348 the largest positive frequency (fs/2*(n-1)/n), and is complex in the 

349 general case. 

350 

351 If the input `a` contains an imaginary part, it is silently discarded. 

352 

353 Examples 

354 -------- 

355 >>> np.fft.fft([0, 1, 0, 0]) 

356 array([ 1.+0.j, 0.-1.j, -1.+0.j, 0.+1.j]) # may vary 

357 >>> np.fft.rfft([0, 1, 0, 0]) 

358 array([ 1.+0.j, 0.-1.j, -1.+0.j]) # may vary 

359 

360 Notice how the final element of the `fft` output is the complex conjugate 

361 of the second element, for real input. For `rfft`, this symmetry is 

362 exploited to compute only the non-negative frequency terms. 

363 

364 """ 

365 a = asarray(a) 

366 inv_norm = 1 

367 if norm is not None and _unitary(norm): 

368 if n is None: 

369 n = a.shape[axis] 

370 inv_norm = sqrt(n) 

371 output = _raw_fft(a, n, axis, True, True, inv_norm) 

372 return output 

373 

374 

375@array_function_dispatch(_fft_dispatcher) 

376def irfft(a, n=None, axis=-1, norm=None): 

377 """ 

378 Compute the inverse of the n-point DFT for real input. 

379 

380 This function computes the inverse of the one-dimensional *n*-point 

381 discrete Fourier Transform of real input computed by `rfft`. 

382 In other words, ``irfft(rfft(a), len(a)) == a`` to within numerical 

383 accuracy. (See Notes below for why ``len(a)`` is necessary here.) 

384 

385 The input is expected to be in the form returned by `rfft`, i.e. the 

386 real zero-frequency term followed by the complex positive frequency terms 

387 in order of increasing frequency. Since the discrete Fourier Transform of 

388 real input is Hermitian-symmetric, the negative frequency terms are taken 

389 to be the complex conjugates of the corresponding positive frequency terms. 

390 

391 Parameters 

392 ---------- 

393 a : array_like 

394 The input array. 

395 n : int, optional 

396 Length of the transformed axis of the output. 

397 For `n` output points, ``n//2+1`` input points are necessary. If the 

398 input is longer than this, it is cropped. If it is shorter than this, 

399 it is padded with zeros. If `n` is not given, it is taken to be 

400 ``2*(m-1)`` where ``m`` is the length of the input along the axis 

401 specified by `axis`. 

402 axis : int, optional 

403 Axis over which to compute the inverse FFT. If not given, the last 

404 axis is used. 

405 norm : {None, "ortho"}, optional 

406 .. versionadded:: 1.10.0 

407 

408 Normalization mode (see `numpy.fft`). Default is None. 

409 

410 Returns 

411 ------- 

412 out : ndarray 

413 The truncated or zero-padded input, transformed along the axis 

414 indicated by `axis`, or the last one if `axis` is not specified. 

415 The length of the transformed axis is `n`, or, if `n` is not given, 

416 ``2*(m-1)`` where ``m`` is the length of the transformed axis of the 

417 input. To get an odd number of output points, `n` must be specified. 

418 

419 Raises 

420 ------ 

421 IndexError 

422 If `axis` is larger than the last axis of `a`. 

423 

424 See Also 

425 -------- 

426 numpy.fft : For definition of the DFT and conventions used. 

427 rfft : The one-dimensional FFT of real input, of which `irfft` is inverse. 

428 fft : The one-dimensional FFT. 

429 irfft2 : The inverse of the two-dimensional FFT of real input. 

430 irfftn : The inverse of the *n*-dimensional FFT of real input. 

431 

432 Notes 

433 ----- 

434 Returns the real valued `n`-point inverse discrete Fourier transform 

435 of `a`, where `a` contains the non-negative frequency terms of a 

436 Hermitian-symmetric sequence. `n` is the length of the result, not the 

437 input. 

438 

439 If you specify an `n` such that `a` must be zero-padded or truncated, the 

440 extra/removed values will be added/removed at high frequencies. One can 

441 thus resample a series to `m` points via Fourier interpolation by: 

442 ``a_resamp = irfft(rfft(a), m)``. 

443 

444 The correct interpretation of the hermitian input depends on the length of 

445 the original data, as given by `n`. This is because each input shape could 

446 correspond to either an odd or even length signal. By default, `irfft` 

447 assumes an even output length which puts the last entry at the Nyquist 

448 frequency; aliasing with its symmetric counterpart. By Hermitian symmetry, 

449 the value is thus treated as purely real. To avoid losing information, the 

450 correct length of the real input **must** be given. 

451 

452 Examples 

453 -------- 

454 >>> np.fft.ifft([1, -1j, -1, 1j]) 

455 array([0.+0.j, 1.+0.j, 0.+0.j, 0.+0.j]) # may vary 

456 >>> np.fft.irfft([1, -1j, -1]) 

457 array([0., 1., 0., 0.]) 

458 

459 Notice how the last term in the input to the ordinary `ifft` is the 

460 complex conjugate of the second term, and the output has zero imaginary 

461 part everywhere. When calling `irfft`, the negative frequencies are not 

462 specified, and the output array is purely real. 

463 

464 """ 

465 a = asarray(a) 

466 if n is None: 

467 n = (a.shape[axis] - 1) * 2 

468 inv_norm = n 

469 if norm is not None and _unitary(norm): 

470 inv_norm = sqrt(n) 

471 output = _raw_fft(a, n, axis, True, False, inv_norm) 

472 return output 

473 

474 

475@array_function_dispatch(_fft_dispatcher) 

476def hfft(a, n=None, axis=-1, norm=None): 

477 """ 

478 Compute the FFT of a signal that has Hermitian symmetry, i.e., a real 

479 spectrum. 

480 

481 Parameters 

482 ---------- 

483 a : array_like 

484 The input array. 

485 n : int, optional 

486 Length of the transformed axis of the output. For `n` output 

487 points, ``n//2 + 1`` input points are necessary. If the input is 

488 longer than this, it is cropped. If it is shorter than this, it is 

489 padded with zeros. If `n` is not given, it is taken to be ``2*(m-1)`` 

490 where ``m`` is the length of the input along the axis specified by 

491 `axis`. 

492 axis : int, optional 

493 Axis over which to compute the FFT. If not given, the last 

494 axis is used. 

495 norm : {None, "ortho"}, optional 

496 Normalization mode (see `numpy.fft`). Default is None. 

497 

498 .. versionadded:: 1.10.0 

499 

500 Returns 

501 ------- 

502 out : ndarray 

503 The truncated or zero-padded input, transformed along the axis 

504 indicated by `axis`, or the last one if `axis` is not specified. 

505 The length of the transformed axis is `n`, or, if `n` is not given, 

506 ``2*m - 2`` where ``m`` is the length of the transformed axis of 

507 the input. To get an odd number of output points, `n` must be 

508 specified, for instance as ``2*m - 1`` in the typical case, 

509 

510 Raises 

511 ------ 

512 IndexError 

513 If `axis` is larger than the last axis of `a`. 

514 

515 See also 

516 -------- 

517 rfft : Compute the one-dimensional FFT for real input. 

518 ihfft : The inverse of `hfft`. 

519 

520 Notes 

521 ----- 

522 `hfft`/`ihfft` are a pair analogous to `rfft`/`irfft`, but for the 

523 opposite case: here the signal has Hermitian symmetry in the time 

524 domain and is real in the frequency domain. So here it's `hfft` for 

525 which you must supply the length of the result if it is to be odd. 

526 

527 * even: ``ihfft(hfft(a, 2*len(a) - 2)) == a``, within roundoff error, 

528 * odd: ``ihfft(hfft(a, 2*len(a) - 1)) == a``, within roundoff error. 

529 

530 The correct interpretation of the hermitian input depends on the length of 

531 the original data, as given by `n`. This is because each input shape could 

532 correspond to either an odd or even length signal. By default, `hfft` 

533 assumes an even output length which puts the last entry at the Nyquist 

534 frequency; aliasing with its symmetric counterpart. By Hermitian symmetry, 

535 the value is thus treated as purely real. To avoid losing information, the 

536 shape of the full signal **must** be given. 

537 

538 Examples 

539 -------- 

540 >>> signal = np.array([1, 2, 3, 4, 3, 2]) 

541 >>> np.fft.fft(signal) 

542 array([15.+0.j, -4.+0.j, 0.+0.j, -1.-0.j, 0.+0.j, -4.+0.j]) # may vary 

543 >>> np.fft.hfft(signal[:4]) # Input first half of signal 

544 array([15., -4., 0., -1., 0., -4.]) 

545 >>> np.fft.hfft(signal, 6) # Input entire signal and truncate 

546 array([15., -4., 0., -1., 0., -4.]) 

547 

548 

549 >>> signal = np.array([[1, 1.j], [-1.j, 2]]) 

550 >>> np.conj(signal.T) - signal # check Hermitian symmetry 

551 array([[ 0.-0.j, -0.+0.j], # may vary 

552 [ 0.+0.j, 0.-0.j]]) 

553 >>> freq_spectrum = np.fft.hfft(signal) 

554 >>> freq_spectrum 

555 array([[ 1., 1.], 

556 [ 2., -2.]]) 

557 

558 """ 

559 a = asarray(a) 

560 if n is None: 

561 n = (a.shape[axis] - 1) * 2 

562 unitary = _unitary(norm) 

563 return irfft(conjugate(a), n, axis) * (sqrt(n) if unitary else n) 

564 

565 

566@array_function_dispatch(_fft_dispatcher) 

567def ihfft(a, n=None, axis=-1, norm=None): 

568 """ 

569 Compute the inverse FFT of a signal that has Hermitian symmetry. 

570 

571 Parameters 

572 ---------- 

573 a : array_like 

574 Input array. 

575 n : int, optional 

576 Length of the inverse FFT, the number of points along 

577 transformation axis in the input to use. If `n` is smaller than 

578 the length of the input, the input is cropped. If it is larger, 

579 the input is padded with zeros. If `n` is not given, the length of 

580 the input along the axis specified by `axis` is used. 

581 axis : int, optional 

582 Axis over which to compute the inverse FFT. If not given, the last 

583 axis is used. 

584 norm : {None, "ortho"}, optional 

585 Normalization mode (see `numpy.fft`). Default is None. 

586 

587 .. versionadded:: 1.10.0 

588 

589 Returns 

590 ------- 

591 out : complex ndarray 

592 The truncated or zero-padded input, transformed along the axis 

593 indicated by `axis`, or the last one if `axis` is not specified. 

594 The length of the transformed axis is ``n//2 + 1``. 

595 

596 See also 

597 -------- 

598 hfft, irfft 

599 

600 Notes 

601 ----- 

602 `hfft`/`ihfft` are a pair analogous to `rfft`/`irfft`, but for the 

603 opposite case: here the signal has Hermitian symmetry in the time 

604 domain and is real in the frequency domain. So here it's `hfft` for 

605 which you must supply the length of the result if it is to be odd: 

606 

607 * even: ``ihfft(hfft(a, 2*len(a) - 2)) == a``, within roundoff error, 

608 * odd: ``ihfft(hfft(a, 2*len(a) - 1)) == a``, within roundoff error. 

609 

610 Examples 

611 -------- 

612 >>> spectrum = np.array([ 15, -4, 0, -1, 0, -4]) 

613 >>> np.fft.ifft(spectrum) 

614 array([1.+0.j, 2.+0.j, 3.+0.j, 4.+0.j, 3.+0.j, 2.+0.j]) # may vary 

615 >>> np.fft.ihfft(spectrum) 

616 array([ 1.-0.j, 2.-0.j, 3.-0.j, 4.-0.j]) # may vary 

617 

618 """ 

619 a = asarray(a) 

620 if n is None: 

621 n = a.shape[axis] 

622 unitary = _unitary(norm) 

623 output = conjugate(rfft(a, n, axis)) 

624 return output * (1 / (sqrt(n) if unitary else n)) 

625 

626 

627def _cook_nd_args(a, s=None, axes=None, invreal=0): 

628 if s is None: 

629 shapeless = 1 

630 if axes is None: 

631 s = list(a.shape) 

632 else: 

633 s = take(a.shape, axes) 

634 else: 

635 shapeless = 0 

636 s = list(s) 

637 if axes is None: 

638 axes = list(range(-len(s), 0)) 

639 if len(s) != len(axes): 

640 raise ValueError("Shape and axes have different lengths.") 

641 if invreal and shapeless: 

642 s[-1] = (a.shape[axes[-1]] - 1) * 2 

643 return s, axes 

644 

645 

646def _raw_fftnd(a, s=None, axes=None, function=fft, norm=None): 

647 a = asarray(a) 

648 s, axes = _cook_nd_args(a, s, axes) 

649 itl = list(range(len(axes))) 

650 itl.reverse() 

651 for ii in itl: 

652 a = function(a, n=s[ii], axis=axes[ii], norm=norm) 

653 return a 

654 

655 

656def _fftn_dispatcher(a, s=None, axes=None, norm=None): 

657 return (a,) 

658 

659 

660@array_function_dispatch(_fftn_dispatcher) 

661def fftn(a, s=None, axes=None, norm=None): 

662 """ 

663 Compute the N-dimensional discrete Fourier Transform. 

664 

665 This function computes the *N*-dimensional discrete Fourier Transform over 

666 any number of axes in an *M*-dimensional array by means of the Fast Fourier 

667 Transform (FFT). 

668 

669 Parameters 

670 ---------- 

671 a : array_like 

672 Input array, can be complex. 

673 s : sequence of ints, optional 

674 Shape (length of each transformed axis) of the output 

675 (``s[0]`` refers to axis 0, ``s[1]`` to axis 1, etc.). 

676 This corresponds to ``n`` for ``fft(x, n)``. 

677 Along any axis, if the given shape is smaller than that of the input, 

678 the input is cropped. If it is larger, the input is padded with zeros. 

679 if `s` is not given, the shape of the input along the axes specified 

680 by `axes` is used. 

681 axes : sequence of ints, optional 

682 Axes over which to compute the FFT. If not given, the last ``len(s)`` 

683 axes are used, or all axes if `s` is also not specified. 

684 Repeated indices in `axes` means that the transform over that axis is 

685 performed multiple times. 

686 norm : {None, "ortho"}, optional 

687 .. versionadded:: 1.10.0 

688 

689 Normalization mode (see `numpy.fft`). Default is None. 

690 

691 Returns 

692 ------- 

693 out : complex ndarray 

694 The truncated or zero-padded input, transformed along the axes 

695 indicated by `axes`, or by a combination of `s` and `a`, 

696 as explained in the parameters section above. 

697 

698 Raises 

699 ------ 

700 ValueError 

701 If `s` and `axes` have different length. 

702 IndexError 

703 If an element of `axes` is larger than than the number of axes of `a`. 

704 

705 See Also 

706 -------- 

707 numpy.fft : Overall view of discrete Fourier transforms, with definitions 

708 and conventions used. 

709 ifftn : The inverse of `fftn`, the inverse *n*-dimensional FFT. 

710 fft : The one-dimensional FFT, with definitions and conventions used. 

711 rfftn : The *n*-dimensional FFT of real input. 

712 fft2 : The two-dimensional FFT. 

713 fftshift : Shifts zero-frequency terms to centre of array 

714 

715 Notes 

716 ----- 

717 The output, analogously to `fft`, contains the term for zero frequency in 

718 the low-order corner of all axes, the positive frequency terms in the 

719 first half of all axes, the term for the Nyquist frequency in the middle 

720 of all axes and the negative frequency terms in the second half of all 

721 axes, in order of decreasingly negative frequency. 

722 

723 See `numpy.fft` for details, definitions and conventions used. 

724 

725 Examples 

726 -------- 

727 >>> a = np.mgrid[:3, :3, :3][0] 

728 >>> np.fft.fftn(a, axes=(1, 2)) 

729 array([[[ 0.+0.j, 0.+0.j, 0.+0.j], # may vary 

730 [ 0.+0.j, 0.+0.j, 0.+0.j], 

731 [ 0.+0.j, 0.+0.j, 0.+0.j]], 

732 [[ 9.+0.j, 0.+0.j, 0.+0.j], 

733 [ 0.+0.j, 0.+0.j, 0.+0.j], 

734 [ 0.+0.j, 0.+0.j, 0.+0.j]], 

735 [[18.+0.j, 0.+0.j, 0.+0.j], 

736 [ 0.+0.j, 0.+0.j, 0.+0.j], 

737 [ 0.+0.j, 0.+0.j, 0.+0.j]]]) 

738 >>> np.fft.fftn(a, (2, 2), axes=(0, 1)) 

739 array([[[ 2.+0.j, 2.+0.j, 2.+0.j], # may vary 

740 [ 0.+0.j, 0.+0.j, 0.+0.j]], 

741 [[-2.+0.j, -2.+0.j, -2.+0.j], 

742 [ 0.+0.j, 0.+0.j, 0.+0.j]]]) 

743 

744 >>> import matplotlib.pyplot as plt 

745 >>> [X, Y] = np.meshgrid(2 * np.pi * np.arange(200) / 12, 

746 ... 2 * np.pi * np.arange(200) / 34) 

747 >>> S = np.sin(X) + np.cos(Y) + np.random.uniform(0, 1, X.shape) 

748 >>> FS = np.fft.fftn(S) 

749 >>> plt.imshow(np.log(np.abs(np.fft.fftshift(FS))**2)) 

750 <matplotlib.image.AxesImage object at 0x...> 

751 >>> plt.show() 

752 

753 """ 

754 

755 return _raw_fftnd(a, s, axes, fft, norm) 

756 

757 

758@array_function_dispatch(_fftn_dispatcher) 

759def ifftn(a, s=None, axes=None, norm=None): 

760 """ 

761 Compute the N-dimensional inverse discrete Fourier Transform. 

762 

763 This function computes the inverse of the N-dimensional discrete 

764 Fourier Transform over any number of axes in an M-dimensional array by 

765 means of the Fast Fourier Transform (FFT). In other words, 

766 ``ifftn(fftn(a)) == a`` to within numerical accuracy. 

767 For a description of the definitions and conventions used, see `numpy.fft`. 

768 

769 The input, analogously to `ifft`, should be ordered in the same way as is 

770 returned by `fftn`, i.e. it should have the term for zero frequency 

771 in all axes in the low-order corner, the positive frequency terms in the 

772 first half of all axes, the term for the Nyquist frequency in the middle 

773 of all axes and the negative frequency terms in the second half of all 

774 axes, in order of decreasingly negative frequency. 

775 

776 Parameters 

777 ---------- 

778 a : array_like 

779 Input array, can be complex. 

780 s : sequence of ints, optional 

781 Shape (length of each transformed axis) of the output 

782 (``s[0]`` refers to axis 0, ``s[1]`` to axis 1, etc.). 

783 This corresponds to ``n`` for ``ifft(x, n)``. 

784 Along any axis, if the given shape is smaller than that of the input, 

785 the input is cropped. If it is larger, the input is padded with zeros. 

786 if `s` is not given, the shape of the input along the axes specified 

787 by `axes` is used. See notes for issue on `ifft` zero padding. 

788 axes : sequence of ints, optional 

789 Axes over which to compute the IFFT. If not given, the last ``len(s)`` 

790 axes are used, or all axes if `s` is also not specified. 

791 Repeated indices in `axes` means that the inverse transform over that 

792 axis is performed multiple times. 

793 norm : {None, "ortho"}, optional 

794 .. versionadded:: 1.10.0 

795 

796 Normalization mode (see `numpy.fft`). Default is None. 

797 

798 Returns 

799 ------- 

800 out : complex ndarray 

801 The truncated or zero-padded input, transformed along the axes 

802 indicated by `axes`, or by a combination of `s` or `a`, 

803 as explained in the parameters section above. 

804 

805 Raises 

806 ------ 

807 ValueError 

808 If `s` and `axes` have different length. 

809 IndexError 

810 If an element of `axes` is larger than than the number of axes of `a`. 

811 

812 See Also 

813 -------- 

814 numpy.fft : Overall view of discrete Fourier transforms, with definitions 

815 and conventions used. 

816 fftn : The forward *n*-dimensional FFT, of which `ifftn` is the inverse. 

817 ifft : The one-dimensional inverse FFT. 

818 ifft2 : The two-dimensional inverse FFT. 

819 ifftshift : Undoes `fftshift`, shifts zero-frequency terms to beginning 

820 of array. 

821 

822 Notes 

823 ----- 

824 See `numpy.fft` for definitions and conventions used. 

825 

826 Zero-padding, analogously with `ifft`, is performed by appending zeros to 

827 the input along the specified dimension. Although this is the common 

828 approach, it might lead to surprising results. If another form of zero 

829 padding is desired, it must be performed before `ifftn` is called. 

830 

831 Examples 

832 -------- 

833 >>> a = np.eye(4) 

834 >>> np.fft.ifftn(np.fft.fftn(a, axes=(0,)), axes=(1,)) 

835 array([[1.+0.j, 0.+0.j, 0.+0.j, 0.+0.j], # may vary 

836 [0.+0.j, 1.+0.j, 0.+0.j, 0.+0.j], 

837 [0.+0.j, 0.+0.j, 1.+0.j, 0.+0.j], 

838 [0.+0.j, 0.+0.j, 0.+0.j, 1.+0.j]]) 

839 

840 

841 Create and plot an image with band-limited frequency content: 

842 

843 >>> import matplotlib.pyplot as plt 

844 >>> n = np.zeros((200,200), dtype=complex) 

845 >>> n[60:80, 20:40] = np.exp(1j*np.random.uniform(0, 2*np.pi, (20, 20))) 

846 >>> im = np.fft.ifftn(n).real 

847 >>> plt.imshow(im) 

848 <matplotlib.image.AxesImage object at 0x...> 

849 >>> plt.show() 

850 

851 """ 

852 

853 return _raw_fftnd(a, s, axes, ifft, norm) 

854 

855 

856@array_function_dispatch(_fftn_dispatcher) 

857def fft2(a, s=None, axes=(-2, -1), norm=None): 

858 """ 

859 Compute the 2-dimensional discrete Fourier Transform 

860 

861 This function computes the *n*-dimensional discrete Fourier Transform 

862 over any axes in an *M*-dimensional array by means of the 

863 Fast Fourier Transform (FFT). By default, the transform is computed over 

864 the last two axes of the input array, i.e., a 2-dimensional FFT. 

865 

866 Parameters 

867 ---------- 

868 a : array_like 

869 Input array, can be complex 

870 s : sequence of ints, optional 

871 Shape (length of each transformed axis) of the output 

872 (``s[0]`` refers to axis 0, ``s[1]`` to axis 1, etc.). 

873 This corresponds to ``n`` for ``fft(x, n)``. 

874 Along each axis, if the given shape is smaller than that of the input, 

875 the input is cropped. If it is larger, the input is padded with zeros. 

876 if `s` is not given, the shape of the input along the axes specified 

877 by `axes` is used. 

878 axes : sequence of ints, optional 

879 Axes over which to compute the FFT. If not given, the last two 

880 axes are used. A repeated index in `axes` means the transform over 

881 that axis is performed multiple times. A one-element sequence means 

882 that a one-dimensional FFT is performed. 

883 norm : {None, "ortho"}, optional 

884 .. versionadded:: 1.10.0 

885 

886 Normalization mode (see `numpy.fft`). Default is None. 

887 

888 Returns 

889 ------- 

890 out : complex ndarray 

891 The truncated or zero-padded input, transformed along the axes 

892 indicated by `axes`, or the last two axes if `axes` is not given. 

893 

894 Raises 

895 ------ 

896 ValueError 

897 If `s` and `axes` have different length, or `axes` not given and 

898 ``len(s) != 2``. 

899 IndexError 

900 If an element of `axes` is larger than than the number of axes of `a`. 

901 

902 See Also 

903 -------- 

904 numpy.fft : Overall view of discrete Fourier transforms, with definitions 

905 and conventions used. 

906 ifft2 : The inverse two-dimensional FFT. 

907 fft : The one-dimensional FFT. 

908 fftn : The *n*-dimensional FFT. 

909 fftshift : Shifts zero-frequency terms to the center of the array. 

910 For two-dimensional input, swaps first and third quadrants, and second 

911 and fourth quadrants. 

912 

913 Notes 

914 ----- 

915 `fft2` is just `fftn` with a different default for `axes`. 

916 

917 The output, analogously to `fft`, contains the term for zero frequency in 

918 the low-order corner of the transformed axes, the positive frequency terms 

919 in the first half of these axes, the term for the Nyquist frequency in the 

920 middle of the axes and the negative frequency terms in the second half of 

921 the axes, in order of decreasingly negative frequency. 

922 

923 See `fftn` for details and a plotting example, and `numpy.fft` for 

924 definitions and conventions used. 

925 

926 

927 Examples 

928 -------- 

929 >>> a = np.mgrid[:5, :5][0] 

930 >>> np.fft.fft2(a) 

931 array([[ 50. +0.j , 0. +0.j , 0. +0.j , # may vary 

932 0. +0.j , 0. +0.j ], 

933 [-12.5+17.20477401j, 0. +0.j , 0. +0.j , 

934 0. +0.j , 0. +0.j ], 

935 [-12.5 +4.0614962j , 0. +0.j , 0. +0.j , 

936 0. +0.j , 0. +0.j ], 

937 [-12.5 -4.0614962j , 0. +0.j , 0. +0.j , 

938 0. +0.j , 0. +0.j ], 

939 [-12.5-17.20477401j, 0. +0.j , 0. +0.j , 

940 0. +0.j , 0. +0.j ]]) 

941 

942 """ 

943 

944 return _raw_fftnd(a, s, axes, fft, norm) 

945 

946 

947@array_function_dispatch(_fftn_dispatcher) 

948def ifft2(a, s=None, axes=(-2, -1), norm=None): 

949 """ 

950 Compute the 2-dimensional inverse discrete Fourier Transform. 

951 

952 This function computes the inverse of the 2-dimensional discrete Fourier 

953 Transform over any number of axes in an M-dimensional array by means of 

954 the Fast Fourier Transform (FFT). In other words, ``ifft2(fft2(a)) == a`` 

955 to within numerical accuracy. By default, the inverse transform is 

956 computed over the last two axes of the input array. 

957 

958 The input, analogously to `ifft`, should be ordered in the same way as is 

959 returned by `fft2`, i.e. it should have the term for zero frequency 

960 in the low-order corner of the two axes, the positive frequency terms in 

961 the first half of these axes, the term for the Nyquist frequency in the 

962 middle of the axes and the negative frequency terms in the second half of 

963 both axes, in order of decreasingly negative frequency. 

964 

965 Parameters 

966 ---------- 

967 a : array_like 

968 Input array, can be complex. 

969 s : sequence of ints, optional 

970 Shape (length of each axis) of the output (``s[0]`` refers to axis 0, 

971 ``s[1]`` to axis 1, etc.). This corresponds to `n` for ``ifft(x, n)``. 

972 Along each axis, if the given shape is smaller than that of the input, 

973 the input is cropped. If it is larger, the input is padded with zeros. 

974 if `s` is not given, the shape of the input along the axes specified 

975 by `axes` is used. See notes for issue on `ifft` zero padding. 

976 axes : sequence of ints, optional 

977 Axes over which to compute the FFT. If not given, the last two 

978 axes are used. A repeated index in `axes` means the transform over 

979 that axis is performed multiple times. A one-element sequence means 

980 that a one-dimensional FFT is performed. 

981 norm : {None, "ortho"}, optional 

982 .. versionadded:: 1.10.0 

983 

984 Normalization mode (see `numpy.fft`). Default is None. 

985 

986 Returns 

987 ------- 

988 out : complex ndarray 

989 The truncated or zero-padded input, transformed along the axes 

990 indicated by `axes`, or the last two axes if `axes` is not given. 

991 

992 Raises 

993 ------ 

994 ValueError 

995 If `s` and `axes` have different length, or `axes` not given and 

996 ``len(s) != 2``. 

997 IndexError 

998 If an element of `axes` is larger than than the number of axes of `a`. 

999 

1000 See Also 

1001 -------- 

1002 numpy.fft : Overall view of discrete Fourier transforms, with definitions 

1003 and conventions used. 

1004 fft2 : The forward 2-dimensional FFT, of which `ifft2` is the inverse. 

1005 ifftn : The inverse of the *n*-dimensional FFT. 

1006 fft : The one-dimensional FFT. 

1007 ifft : The one-dimensional inverse FFT. 

1008 

1009 Notes 

1010 ----- 

1011 `ifft2` is just `ifftn` with a different default for `axes`. 

1012 

1013 See `ifftn` for details and a plotting example, and `numpy.fft` for 

1014 definition and conventions used. 

1015 

1016 Zero-padding, analogously with `ifft`, is performed by appending zeros to 

1017 the input along the specified dimension. Although this is the common 

1018 approach, it might lead to surprising results. If another form of zero 

1019 padding is desired, it must be performed before `ifft2` is called. 

1020 

1021 Examples 

1022 -------- 

1023 >>> a = 4 * np.eye(4) 

1024 >>> np.fft.ifft2(a) 

1025 array([[1.+0.j, 0.+0.j, 0.+0.j, 0.+0.j], # may vary 

1026 [0.+0.j, 0.+0.j, 0.+0.j, 1.+0.j], 

1027 [0.+0.j, 0.+0.j, 1.+0.j, 0.+0.j], 

1028 [0.+0.j, 1.+0.j, 0.+0.j, 0.+0.j]]) 

1029 

1030 """ 

1031 

1032 return _raw_fftnd(a, s, axes, ifft, norm) 

1033 

1034 

1035@array_function_dispatch(_fftn_dispatcher) 

1036def rfftn(a, s=None, axes=None, norm=None): 

1037 """ 

1038 Compute the N-dimensional discrete Fourier Transform for real input. 

1039 

1040 This function computes the N-dimensional discrete Fourier Transform over 

1041 any number of axes in an M-dimensional real array by means of the Fast 

1042 Fourier Transform (FFT). By default, all axes are transformed, with the 

1043 real transform performed over the last axis, while the remaining 

1044 transforms are complex. 

1045 

1046 Parameters 

1047 ---------- 

1048 a : array_like 

1049 Input array, taken to be real. 

1050 s : sequence of ints, optional 

1051 Shape (length along each transformed axis) to use from the input. 

1052 (``s[0]`` refers to axis 0, ``s[1]`` to axis 1, etc.). 

1053 The final element of `s` corresponds to `n` for ``rfft(x, n)``, while 

1054 for the remaining axes, it corresponds to `n` for ``fft(x, n)``. 

1055 Along any axis, if the given shape is smaller than that of the input, 

1056 the input is cropped. If it is larger, the input is padded with zeros. 

1057 if `s` is not given, the shape of the input along the axes specified 

1058 by `axes` is used. 

1059 axes : sequence of ints, optional 

1060 Axes over which to compute the FFT. If not given, the last ``len(s)`` 

1061 axes are used, or all axes if `s` is also not specified. 

1062 norm : {None, "ortho"}, optional 

1063 .. versionadded:: 1.10.0 

1064 

1065 Normalization mode (see `numpy.fft`). Default is None. 

1066 

1067 Returns 

1068 ------- 

1069 out : complex ndarray 

1070 The truncated or zero-padded input, transformed along the axes 

1071 indicated by `axes`, or by a combination of `s` and `a`, 

1072 as explained in the parameters section above. 

1073 The length of the last axis transformed will be ``s[-1]//2+1``, 

1074 while the remaining transformed axes will have lengths according to 

1075 `s`, or unchanged from the input. 

1076 

1077 Raises 

1078 ------ 

1079 ValueError 

1080 If `s` and `axes` have different length. 

1081 IndexError 

1082 If an element of `axes` is larger than than the number of axes of `a`. 

1083 

1084 See Also 

1085 -------- 

1086 irfftn : The inverse of `rfftn`, i.e. the inverse of the n-dimensional FFT 

1087 of real input. 

1088 fft : The one-dimensional FFT, with definitions and conventions used. 

1089 rfft : The one-dimensional FFT of real input. 

1090 fftn : The n-dimensional FFT. 

1091 rfft2 : The two-dimensional FFT of real input. 

1092 

1093 Notes 

1094 ----- 

1095 The transform for real input is performed over the last transformation 

1096 axis, as by `rfft`, then the transform over the remaining axes is 

1097 performed as by `fftn`. The order of the output is as for `rfft` for the 

1098 final transformation axis, and as for `fftn` for the remaining 

1099 transformation axes. 

1100 

1101 See `fft` for details, definitions and conventions used. 

1102 

1103 Examples 

1104 -------- 

1105 >>> a = np.ones((2, 2, 2)) 

1106 >>> np.fft.rfftn(a) 

1107 array([[[8.+0.j, 0.+0.j], # may vary 

1108 [0.+0.j, 0.+0.j]], 

1109 [[0.+0.j, 0.+0.j], 

1110 [0.+0.j, 0.+0.j]]]) 

1111 

1112 >>> np.fft.rfftn(a, axes=(2, 0)) 

1113 array([[[4.+0.j, 0.+0.j], # may vary 

1114 [4.+0.j, 0.+0.j]], 

1115 [[0.+0.j, 0.+0.j], 

1116 [0.+0.j, 0.+0.j]]]) 

1117 

1118 """ 

1119 a = asarray(a) 

1120 s, axes = _cook_nd_args(a, s, axes) 

1121 a = rfft(a, s[-1], axes[-1], norm) 

1122 for ii in range(len(axes)-1): 

1123 a = fft(a, s[ii], axes[ii], norm) 

1124 return a 

1125 

1126 

1127@array_function_dispatch(_fftn_dispatcher) 

1128def rfft2(a, s=None, axes=(-2, -1), norm=None): 

1129 """ 

1130 Compute the 2-dimensional FFT of a real array. 

1131 

1132 Parameters 

1133 ---------- 

1134 a : array 

1135 Input array, taken to be real. 

1136 s : sequence of ints, optional 

1137 Shape of the FFT. 

1138 axes : sequence of ints, optional 

1139 Axes over which to compute the FFT. 

1140 norm : {None, "ortho"}, optional 

1141 .. versionadded:: 1.10.0 

1142 

1143 Normalization mode (see `numpy.fft`). Default is None. 

1144 

1145 Returns 

1146 ------- 

1147 out : ndarray 

1148 The result of the real 2-D FFT. 

1149 

1150 See Also 

1151 -------- 

1152 rfftn : Compute the N-dimensional discrete Fourier Transform for real 

1153 input. 

1154 

1155 Notes 

1156 ----- 

1157 This is really just `rfftn` with different default behavior. 

1158 For more details see `rfftn`. 

1159 

1160 """ 

1161 

1162 return rfftn(a, s, axes, norm) 

1163 

1164 

1165@array_function_dispatch(_fftn_dispatcher) 

1166def irfftn(a, s=None, axes=None, norm=None): 

1167 """ 

1168 Compute the inverse of the N-dimensional FFT of real input. 

1169 

1170 This function computes the inverse of the N-dimensional discrete 

1171 Fourier Transform for real input over any number of axes in an 

1172 M-dimensional array by means of the Fast Fourier Transform (FFT). In 

1173 other words, ``irfftn(rfftn(a), a.shape) == a`` to within numerical 

1174 accuracy. (The ``a.shape`` is necessary like ``len(a)`` is for `irfft`, 

1175 and for the same reason.) 

1176 

1177 The input should be ordered in the same way as is returned by `rfftn`, 

1178 i.e. as for `irfft` for the final transformation axis, and as for `ifftn` 

1179 along all the other axes. 

1180 

1181 Parameters 

1182 ---------- 

1183 a : array_like 

1184 Input array. 

1185 s : sequence of ints, optional 

1186 Shape (length of each transformed axis) of the output 

1187 (``s[0]`` refers to axis 0, ``s[1]`` to axis 1, etc.). `s` is also the 

1188 number of input points used along this axis, except for the last axis, 

1189 where ``s[-1]//2+1`` points of the input are used. 

1190 Along any axis, if the shape indicated by `s` is smaller than that of 

1191 the input, the input is cropped. If it is larger, the input is padded 

1192 with zeros. If `s` is not given, the shape of the input along the axes 

1193 specified by axes is used. Except for the last axis which is taken to be 

1194 ``2*(m-1)`` where ``m`` is the length of the input along that axis. 

1195 axes : sequence of ints, optional 

1196 Axes over which to compute the inverse FFT. If not given, the last 

1197 `len(s)` axes are used, or all axes if `s` is also not specified. 

1198 Repeated indices in `axes` means that the inverse transform over that 

1199 axis is performed multiple times. 

1200 norm : {None, "ortho"}, optional 

1201 .. versionadded:: 1.10.0 

1202 

1203 Normalization mode (see `numpy.fft`). Default is None. 

1204 

1205 Returns 

1206 ------- 

1207 out : ndarray 

1208 The truncated or zero-padded input, transformed along the axes 

1209 indicated by `axes`, or by a combination of `s` or `a`, 

1210 as explained in the parameters section above. 

1211 The length of each transformed axis is as given by the corresponding 

1212 element of `s`, or the length of the input in every axis except for the 

1213 last one if `s` is not given. In the final transformed axis the length 

1214 of the output when `s` is not given is ``2*(m-1)`` where ``m`` is the 

1215 length of the final transformed axis of the input. To get an odd 

1216 number of output points in the final axis, `s` must be specified. 

1217 

1218 Raises 

1219 ------ 

1220 ValueError 

1221 If `s` and `axes` have different length. 

1222 IndexError 

1223 If an element of `axes` is larger than than the number of axes of `a`. 

1224 

1225 See Also 

1226 -------- 

1227 rfftn : The forward n-dimensional FFT of real input, 

1228 of which `ifftn` is the inverse. 

1229 fft : The one-dimensional FFT, with definitions and conventions used. 

1230 irfft : The inverse of the one-dimensional FFT of real input. 

1231 irfft2 : The inverse of the two-dimensional FFT of real input. 

1232 

1233 Notes 

1234 ----- 

1235 See `fft` for definitions and conventions used. 

1236 

1237 See `rfft` for definitions and conventions used for real input. 

1238 

1239 The correct interpretation of the hermitian input depends on the shape of 

1240 the original data, as given by `s`. This is because each input shape could 

1241 correspond to either an odd or even length signal. By default, `irfftn` 

1242 assumes an even output length which puts the last entry at the Nyquist 

1243 frequency; aliasing with its symmetric counterpart. When performing the 

1244 final complex to real transform, the last value is thus treated as purely 

1245 real. To avoid losing information, the correct shape of the real input 

1246 **must** be given. 

1247 

1248 Examples 

1249 -------- 

1250 >>> a = np.zeros((3, 2, 2)) 

1251 >>> a[0, 0, 0] = 3 * 2 * 2 

1252 >>> np.fft.irfftn(a) 

1253 array([[[1., 1.], 

1254 [1., 1.]], 

1255 [[1., 1.], 

1256 [1., 1.]], 

1257 [[1., 1.], 

1258 [1., 1.]]]) 

1259 

1260 """ 

1261 a = asarray(a) 

1262 s, axes = _cook_nd_args(a, s, axes, invreal=1) 

1263 for ii in range(len(axes)-1): 

1264 a = ifft(a, s[ii], axes[ii], norm) 

1265 a = irfft(a, s[-1], axes[-1], norm) 

1266 return a 

1267 

1268 

1269@array_function_dispatch(_fftn_dispatcher) 

1270def irfft2(a, s=None, axes=(-2, -1), norm=None): 

1271 """ 

1272 Compute the 2-dimensional inverse FFT of a real array. 

1273 

1274 Parameters 

1275 ---------- 

1276 a : array_like 

1277 The input array 

1278 s : sequence of ints, optional 

1279 Shape of the real output to the inverse FFT. 

1280 axes : sequence of ints, optional 

1281 The axes over which to compute the inverse fft. 

1282 Default is the last two axes. 

1283 norm : {None, "ortho"}, optional 

1284 .. versionadded:: 1.10.0 

1285 

1286 Normalization mode (see `numpy.fft`). Default is None. 

1287 

1288 Returns 

1289 ------- 

1290 out : ndarray 

1291 The result of the inverse real 2-D FFT. 

1292 

1293 See Also 

1294 -------- 

1295 irfftn : Compute the inverse of the N-dimensional FFT of real input. 

1296 

1297 Notes 

1298 ----- 

1299 This is really `irfftn` with different defaults. 

1300 For more details see `irfftn`. 

1301 

1302 """ 

1303 

1304 return irfftn(a, s, axes, norm)