You can not select more than 25 topics Topics must start with a chinese character,a letter or number, can include dashes ('-') and can be up to 35 characters long.

149 lines
3.9 kB

  1. __author__ = 'Jo Schlemper'
  2. import numpy as np
  3. sqrt = np.sqrt
  4. from numpy.fft import fft, fft2, ifft2, ifft, ifftshift, fftshift
  5. def fftc(x, axis=-1, norm='ortho'):
  6. ''' expect x as m*n matrix '''
  7. return fftshift(fft(ifftshift(x, axes=axis), axis=axis, norm=norm), axes=axis)
  8. def ifftc(x, axis=-1, norm='ortho'):
  9. ''' expect x as m*n matrix '''
  10. return fftshift(ifft(ifftshift(x, axes=axis), axis=axis, norm=norm), axes=axis)
  11. def fft2c(x):
  12. '''
  13. Centered fft
  14. Note: fft2 applies fft to last 2 axes by default
  15. :param x: 2D onwards. e.g: if its 3d, x.shape = (n,row,col). 4d:x.shape = (n,slice,row,col)
  16. :return:
  17. '''
  18. # axes = (len(x.shape)-2, len(x.shape)-1) # get last 2 axes
  19. axes = (-2, -1) # get last 2 axes
  20. res = fftshift(fft2(ifftshift(x, axes=axes), norm='ortho'), axes=axes)
  21. return res
  22. def ifft2c(x):
  23. '''
  24. Centered ifft
  25. Note: fft2 applies fft to last 2 axes by default
  26. :param x: 2D onwards. e.g: if its 3d, x.shape = (n,row,col). 4d:x.shape = (n,slice,row,col)
  27. :return:
  28. '''
  29. axes = (-2, -1) # get last 2 axes
  30. res = fftshift(ifft2(ifftshift(x, axes=axes), norm='ortho'), axes=axes)
  31. return res
  32. def fourier_matrix(rows, cols):
  33. '''
  34. parameters:
  35. rows: number or rows
  36. cols: number of columns
  37. return unitary (rows x cols) fourier matrix
  38. '''
  39. # from scipy.linalg import dft
  40. # return dft(rows,scale='sqrtn')
  41. col_range = np.arange(cols)
  42. row_range = np.arange(rows)
  43. scale = 1 / np.sqrt(cols)
  44. coeffs = np.outer(row_range, col_range)
  45. fourier_matrix = np.exp(coeffs * (-2. * np.pi * 1j / cols)) * scale
  46. return fourier_matrix
  47. def inverse_fourier_matrix(rows, cols):
  48. return np.array(np.matrix(fourier_matrix(rows, cols)).getH())
  49. def flip(m, axis):
  50. """
  51. ==== > Only in numpy 1.12 < =====
  52. Reverse the order of elements in an array along the given axis.
  53. The shape of the array is preserved, but the elements are reordered.
  54. .. versionadded:: 1.12.0
  55. Parameters
  56. ----------
  57. m : array_like
  58. Input array.
  59. axis : integer
  60. Axis in array, which entries are reversed.
  61. Returns
  62. -------
  63. out : array_like
  64. A view of `m` with the entries of axis reversed. Since a view is
  65. returned, this operation is done in constant time.
  66. See Also
  67. --------
  68. flipud : Flip an array vertically (axis=0).
  69. fliplr : Flip an array horizontally (axis=1).
  70. Notes
  71. -----
  72. flip(m, 0) is equivalent to flipud(m).
  73. flip(m, 1) is equivalent to fliplr(m).
  74. flip(m, n) corresponds to ``m[...,::-1,...]`` with ``::-1`` at position n.
  75. Examples
  76. --------
  77. >>> A = np.arange(8).reshape((2,2,2))
  78. >>> A
  79. array([[[0, 1],
  80. [2, 3]],
  81. [[4, 5],
  82. [6, 7]]])
  83. >>> flip(A, 0)
  84. array([[[4, 5],
  85. [6, 7]],
  86. [[0, 1],
  87. [2, 3]]])
  88. >>> flip(A, 1)
  89. array([[[2, 3],
  90. [0, 1]],
  91. [[6, 7],
  92. [4, 5]]])
  93. >>> A = np.random.randn(3,4,5)
  94. >>> np.all(flip(A,2) == A[:,:,::-1,...])
  95. True
  96. """
  97. if not hasattr(m, 'ndim'):
  98. m = np.asarray(m)
  99. indexer = [slice(None)] * m.ndim
  100. try:
  101. indexer[axis] = slice(None, None, -1)
  102. except IndexError:
  103. raise ValueError("axis=%i is invalid for the %i-dimensional input array"
  104. % (axis, m.ndim))
  105. return m[tuple(indexer)]
  106. def rot90_nd(x, axes=(-2, -1), k=1):
  107. """Rotates selected axes"""
  108. def flipud(x):
  109. return flip(x, axes[0])
  110. def fliplr(x):
  111. return flip(x, axes[1])
  112. x = np.asanyarray(x)
  113. if x.ndim < 2:
  114. raise ValueError("Input must >= 2-d.")
  115. k = k % 4
  116. if k == 0:
  117. return x
  118. elif k == 1:
  119. return fliplr(x).swapaxes(*axes)
  120. elif k == 2:
  121. return fliplr(flipud(x))
  122. else:
  123. # k == 3
  124. return fliplr(x.swapaxes(*axes))

简介

No Description

Python

贡献者 (1)