其他分享
首页 > 其他分享> > 如何在R中写fftshift和ifftshift?

如何在R中写fftshift和ifftshift?

作者:互联网

在numpy,我们有以下功能:

import numpy
from numpy.fft import fft2, ifft2, fftshift, ifftshift

我想在R中的R. fft中重写这些函数,就像python中的fft或fft2一样.同样对于ifft2,我们必须做fft(,inverse = T)

现在我想知道如何在R中有效地重写fftshift和ifftshift函数(对于矩阵).

解决方法:

fftshift和ifftshift背后的概念非常简单.这是我从MathWorks(MATLAB的创建者)中提取的数字:

Source: MathWorks doc page on fftshift

想象一下,您的输入2D矩阵被分割成象限.象限#1位于左上方,象限#2位于右上方,象限#3位于右下方,象限#4位于左下方.对于2D矩阵,fftshift默认交换第一和第三象限以及第二和第四象限.您可以覆盖此行为,只需在一个维度上单独执行fftshift即可.如果这样做,您将交换所谓的半空格.如果您指定沿行交换(即维度1),则矩阵的上半部分与下半部分交换.如果您指定沿列交换(即维度2),则右半部分与左半部分交换:

Source: MathWorks doc page on fftshift

默认情况下,使用fftshift链接按顺序交换尺寸1和尺寸2.如果你有一个偶数大小的矩阵,其中行和列是偶数,那么将矩阵切割成四个部分并进行交换是非常明确的.但是,如果矩阵是奇数大小,则取决于您正在查看的语言.例如,在MATLAB和Python numpy中,执行切换的位置定义为(r,c)= ceil(rows / 2),ceil(cols / 2),其中rows和cols是矩阵的行和列. r和c是交换发生的行和列.

对于ifftshift,您只需反转在fftshift上执行的操作.因此,默认操作是交换尺寸2然后尺寸1.然而,您必须重新定义切换中心对于奇数尺寸矩阵的位置.你必须使用地板,而不是ceil,因为这可以精确地确定在执行ffthift之后半空间的位置,而你现在正在撤消对原始矩阵所做的操作.因此,新的切换中心是(r,c)= floor(rows / 2),floor(cols / 2).除此之外,在单个维度之间交换的逻辑是相同的 – 只是切换中心现在已经改变.

因此,这就是我想出的.这些函数采用R中定义的矩阵以及第二个可选参数来确定要交换的维度.省略此操作会执行我刚才谈到的默认象限交换:

fftshift <- function(input_matrix, dim = -1) {

    rows <- dim(input_matrix)[1]    
    cols <- dim(input_matrix)[2]    

    swap_up_down <- function(input_matrix) {
        rows_half <- ceiling(rows/2)
        return(rbind(input_matrix[((rows_half+1):rows), (1:cols)], input_matrix[(1:rows_half), (1:cols)]))
    }

    swap_left_right <- function(input_matrix) {
        cols_half <- ceiling(cols/2)
        return(cbind(input_matrix[1:rows, ((cols_half+1):cols)], input_matrix[1:rows, 1:cols_half]))
    }

    if (dim == -1) {
        input_matrix <- swap_up_down(input_matrix)
        return(swap_left_right(input_matrix))
    }
    else if (dim == 1) {
        return(swap_up_down(input_matrix))
    }
    else if (dim == 2) {
        return(swap_left_right(input_matrix))
    }
    else {
        stop("Invalid dimension parameter")
    }
}

ifftshift <- function(input_matrix, dim = -1) {

    rows <- dim(input_matrix)[1]    
    cols <- dim(input_matrix)[2]    

    swap_up_down <- function(input_matrix) {
        rows_half <- floor(rows/2)
        return(rbind(input_matrix[((rows_half+1):rows), (1:cols)], input_matrix[(1:rows_half), (1:cols)]))
    }

    swap_left_right <- function(input_matrix) {
        cols_half <- floor(cols/2)
        return(cbind(input_matrix[1:rows, ((cols_half+1):cols)], input_matrix[1:rows, 1:cols_half]))
    }

    if (dim == -1) {
        input_matrix <- swap_left_right(input_matrix)
        return(swap_up_down(input_matrix))
    }
    else if (dim == 1) {
        return(swap_up_down(input_matrix))
    }
    else if (dim == 2) {
        return(swap_left_right(input_matrix))
    }
    else {
        stop("Invalid dimension parameter")
    }
}

在每个函数中,我定义了交换左半部分和右半部分以及上半部分和下半部分的函数.要交换左半部分和右半部分,只需确定用于执行交换的列,并使用索引通过将这两个部分连接在一起来创建新矩阵,其中前半部分是右半部分,后半部分是左半部分.在交换上半部分和下半部分但找到正确的行来执行交换时,您也会这样做.

第二个参数dim可以选择设置为1或2以交换相应的维度.请注意,fftshift和ifftshift之间的唯一区别是定义的交换中心以及默认交换两个维度时的操作顺序.其余的代码是一样的.

作为演示的一种方式,假设我已经声明了一个5 x 5数字矩阵,如下所示:

> O <- matrix(1:25, 5, 5)
> O
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    6   11   16   21
[2,]    2    7   12   17   22
[3,]    3    8   13   18   23
[4,]    4    9   14   19   24
[5,]    5   10   15   20   25

请注意,矩阵的大小在两个维度上都是奇数.执行ffthift给我们:

> P <- fftshift(O)
> P
     [,1] [,2] [,3] [,4] [,5]
[1,]   19   24    4    9   14
[2,]   20   25    5   10   15
[3,]   16   21    1    6   11
[4,]   17   22    2    7   12
[5,]   18   23    3    8   13

您可以看到相应的尺寸已被交换.用ifftshift反转这个应该给我们原来的矩阵,它确实:

> Q <- ifftshift(P)
> Q
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    6   11   16   21
[2,]    2    7   12   17   22
[3,]    3    8   13   18   23
[4,]    4    9   14   19   24
[5,]    5   10   15   20   25

当然你可以覆盖它并指定你想要交换的维度,所以假设你想只交换行:

> fftshift(O, 1)
     [,1] [,2] [,3] [,4] [,5]
[1,]    4    9   14   19   24
[2,]    5   10   15   20   25
[3,]    1    6   11   16   21
[4,]    2    7   12   17   22
[5,]    3    8   13   18   23

> ifftshift(fftshift(O, 1), 1)
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    6   11   16   21
[2,]    2    7   12   17   22
[3,]    3    8   13   18   23
[4,]    4    9   14   19   24
[5,]    5   10   15   20   25

请注意,当您在一个尺寸上交换的矩阵上执行ifftshift时,必须使用与使用fftshift时用于交换的相同尺寸,以确保返回原始矩阵.

我们也可以为第二个维度做同样的事情:

> ifftshift(O, 2)
     [,1] [,2] [,3] [,4] [,5]
[1,]   11   16   21    1    6
[2,]   12   17   22    2    7
[3,]   13   18   23    3    8
[4,]   14   19   24    4    9
[5,]   15   20   25    5   10

> ifftshift(fftshift(O, 2), 2)
     [,1] [,2] [,3] [,4] [,5]
[1,]    1    6   11   16   21
[2,]    2    7   12   17   22
[3,]    3    8   13   18   23
[4,]    4    9   14   19   24
[5,]    5   10   15   20   25

作为一个有趣的说明,我们可以验证numpy与我们上面讨论的相同.这是一个IPython会话:

In [16]: import numpy as np

In [17]: from numpy.fft import fftshift, ifftshift

In [18]: O = np.reshape(np.arange(1,26), (5,5)).T

In [19]: O
Out[19]:
array([[ 1,  6, 11, 16, 21],
       [ 2,  7, 12, 17, 22],
       [ 3,  8, 13, 18, 23],
       [ 4,  9, 14, 19, 24],
       [ 5, 10, 15, 20, 25]])

In [20]: fftshift(O)
Out[20]:
array([[19, 24,  4,  9, 14],
       [20, 25,  5, 10, 15],
       [16, 21,  1,  6, 11],
       [17, 22,  2,  7, 12],
       [18, 23,  3,  8, 13]])

In [21]: ifftshift(fftshift(O))
Out[21]:
array([[ 1,  6, 11, 16, 21],
       [ 2,  7, 12, 17, 22],
       [ 3,  8, 13, 18, 23],
       [ 4,  9, 14, 19, 24],
       [ 5, 10, 15, 20, 25]])

导入numpy以及来自numpy.fft包的fftshift和ifftshift并创建一个2D矩阵,它与你在R中定义的例子中看到的相同.然后我们调用fftshift,然后依次调整fftshift和ifftshift可以看到我们得到了与R代码中看到的相同的结果.

标签:python,r,numpy,fft,ifft
来源: https://codeday.me/bug/20190829/1757356.html