Matlab psf2otf与fft2函数的关系

落日映苍穹つ 2022-08-05 07:54 277阅读 0赞

**psf2otf**

--------------------

> Convert point-spread function to optical transfer function

**Syntax**

> OTF = psf2otf(PSF)  
> OTF = psf2otf(PSF,OUTSIZE)

**Description**

> OTF = psf2otf(PSF) computes the fast Fourier transform (FFT) of the  
> point-spread function (PSF) array and creates the optical transfer  
> function array, OTF, that is not influenced by the PSF off-centering.  
> By default, the OTF array is the same size as the PSF array.
> 
> OTF = psf2otf(PSF,OUTSIZE) converts the PSF array into an OTF array,  
> where OUTSIZE specifies the size of the OTF array. OUTSIZE cannot be  
> smaller than the PSF array size in any dimension.
> 
> To ensure that the OTF is not altered because of PSF off-centering,  
> psf2otf postpads the PSF array (down or to the right) with 0’s to  
> match dimensions specified in OUTSIZE, then circularly shifts the  
> values of the PSF array up (or to the left) until the central pixel  
> reaches (1,1) position.
> 
> Note that this function is used in image convolution/deconvolution  
> when the operations involve the FFT.

**fft2**

--------------------

> 2-D fast Fourier transform

**Syntax**

> Y = fft2(X)  
> Y = fft2(X,m,n)

两者之间的官方解释也说了，psf2otf，其作用是将一个空间点扩散函数转换为频谱面的光学传递函数，执行的也是对PSF的FFT变换。真的只是这样吗？

--------------------

**大胆假设 小心求证**  
我们用matlab验证一下，首先我们生成PSF函数（其实就是空间卷积嘛，囧），随便来一个拉普拉斯模板吧，So easy。

l = [0 -1 0;-1 4 -1;0 -1 0]

*l*=⎡⎣⎢0−10−14−10−10⎤⎦⎥

接下来，对  *l* 进行傅里叶变换，这里为了显示方便，我简单选择m，n都是3。

F = fft2(l,3,3)

*F*=⎡⎣⎢0−1.5000−2.5981*i*−1.5000\+2.5981*i*−1.5000−2.5981*i*−3.0000\+5.1962*i*6.0000−1.5000\+2.5981*i*6.0000−3.0000−5.1962*i*⎤⎦⎥

如果使用psf2otf函数呢

P = psf2otf(l,[3,3])

*P*=⎡⎣⎢033366666⎤⎦⎥

这里可以看出两者结果是不同的，看下面解释吧，我就不翻译了。  
官方解释如下：

> **To ensure that the OTF is not altered because of PSF off-centering,**  
> **psf2otf postpads the PSF array (down or to the right) with 0’s to**  
> **match dimensions specified in OUTSIZE, then circularly shifts the**  
> **values of the PSF array up (or to the left) until the central pixel**  
> **reaches (1,1) position.**

我们试一下，将*l*的中心元素移动到矩阵的左上角。

S = circshift(l,-floor(size(l)/2))

*S*=⎡⎣⎢4−1−1−100−100⎤⎦⎥

然后再进行  *f**f**t* 变换

FS = fft2(S)

*F**S*=⎡⎣⎢033366666⎤⎦⎥

这样两者就一样了。  
实际上psf2otf中最后还舍弃由于浮点数运算带来的微小误差。  
Matlab中psf2otf源码如下：

[psf, psfSize, outSize] = ParseInputs(varargin{:});
    
    if  ~all(psf(:)==0),
    
       % Pad the PSF to outSize
       padSize = outSize - psfSize;
       psf     = padarray(psf, padSize, 'post');
    
       % Circularly shift otf so that the "center" of the PSF is at the
       % (1,1) element of the array.
       psf    = circshift(psf,-floor(psfSize/2));
    
       % Compute the OTF
       otf = fftn(psf);
    
       % Estimate the rough number of operations involved in the 
       % computation of the FFT.
       nElem = prod(psfSize);
       nOps  = 0;
       for k=1:ndims(psf)
          nffts = nElem/psfSize(k);
          nOps  = nOps + psfSize(k)*log2(psfSize(k))*nffts; 
       end
    
       % Discard the imaginary part of the psf if it's within roundoff error.
       if max(abs(imag(otf(:))))/max(abs(otf(:))) <= nOps*eps
          otf = real(otf);
       end
    else
       otf = zeros(outSize);
    end

## 更多阅读 ##

[http://www.imageprocessingplace.com/forumphpBB3/phpBB3/viewtopic.php?f=1&t=9][http_www.imageprocessingplace.com_forumphpBB3_phpBB3_viewtopic.php_f_1_t_9]

[http_www.imageprocessingplace.com_forumphpBB3_phpBB3_viewtopic.php_f_1_t_9]: http://www.imageprocessingplace.com/forumphpBB3/phpBB3/viewtopic.php?f=1&t=9