2 回答

• 7

  您可以使用Luis Mendo的帖子立即确定滤镜的频谱 . 对于1D滤波器使用freqz，对于2D滤波器使用freqz2 . 请注意，使用 freqz 和 freqz2 的图表是在 [-1,1] 之间的标准化频率方面 .

  但是，我和'd like to write this post to compliment Luis Mendo'显示了频谱的来源 .

  ---

  让我们从您的第一个过滤器开始：

  h1 = [-1,3,-1];
  

  如果你还记得，1D Discrete-Time Fourier Transform (DTFT)被定义为：

  

  ---

  旁白：为什么DTFT而不是离散傅里叶变换（DFT）？

  请特别注意这是 different 来自Discrete Fourier Transform (DFT) . 它们之间的区别在于DTFT是应用于离散信号的传统傅立叶变换 . 从连续的角度来看，我们应用常规傅立叶变换，其中频谱在频率方面是连续的 . 对于离散信号，这里基本上是相同的，但是变换应用于离散信号 . 输出在频率上也是连续的，附加约束是频谱是周期性的 . 对于DFT，它本质上是DTFT in frequency domain 的 sampled 版本 . 因此，它们之间的核心差异在于频域中，一个是连续的，而另一个是连续对应的采样版本 . 我们需要对DTFT进行采样的原因是，您可以将光谱实际存储在计算机上，让程序处理光谱（即MATLAB） .

  您的问题要求确定频谱是什么，从理论角度来看，我将向您展示DTFT . 当我们实际上是光谱时，无论如何我们实际上是DTFT，所以当我们看到光谱时，你会看到DFT（或多或少！） .

  可以在DSP StackExchange上找到关于它们之间差异的非常好的解释：https://dsp.stackexchange.com/questions/16586/difference-between-discrete-time-fourier-transform-and-discrete-fourier-transfor

  ---

  回到你的问题

  对于DTFT的定义， h[k] 是1D信号， omega 是以弧度定义的角频率 . 因此，您可以将滤波器视为此1D信号，当您在空间域中进行滤波时，它与获取此信号，将其转换为频域并与频域中的其他输入信号进行乘法相同 .

  因此，如果我们认为您的过滤器是对称的，则将值3定义为 centre 点 . 因此，您可以将 h[k] 视为：

  h = [-1    3  -1]
        ^    ^   ^
  k =   -1   0   1
  

  因此，频域表示仅仅是具有复指数项的滤波器系数的加权和 . 将 h[k] 替换为傅里叶变换公式，我们得到：

  

  如果你还记得Euler's formula，我们有：

  

  同理：

  

  如果我们将两个方程加在一起，我们可以重新排列并求解 cos(x) ：

  

  回到我们的1D滤波器的变换，我们可以做一些因子：

  

  请注意，域在 [-pi,pi] 之间，因为1D傅里叶变换是对称的 2*pi . 因此，如果要显示光谱，只需使用上面的域绘制并使用我刚创建的方程，绘制 absolute 值，因为光谱的幅度始终为正：

  w = linspace(-pi,pi);
  h = abs(3 - 2*cos(w));
  plot(w,h);
  title('Magnitude of 1D spectrum');
  axis([-pi, pi, 0, 5]);
  

  linspace从 -pi 到 pi 生成线性增加的数组，范围之间有100个点 . 您可以通过指定手动告诉 linspace 要生成多少点的第三个参数来覆盖此行为，但如果省略此参数，则默认值为100 .

  请注意，我已经使 y -axis从0开始扩展，因此您可以看到曲线的起始位置，并且它会上升到5，因为这是 h 的最大值 . 这就是我们得到的：

  

  实际上，上述频谱看起来像是一个高通滤波器，并且由于直流频率（ w = 0 ）的幅度为1，您实际上是在高通滤波结果的基础上添加原始信号，因此您的信号就是这样 .

  ---

  但是，您可以使用2D案例执行相同的过程它将涉及更多一点 . 在2D情况下，离散时间傅立叶变换被定义为：

  

  我们有两个独立的变量需要考虑，我们将有2D空间频率 . w1 和 k1 将沿着2D信号的 rows 运行， w2 和 k2 将沿着2D信号的 columns 运行 .

  鉴于你的面具：

  h2 = [-1 -1 -1; -1 8 -1; -1 -1 -1]
  

  由于 h2 的形状是对称的，因此值8将是位置 (w1,w2) = (0,0) . 因此，当我们用上面的等式计算光谱时，我们得到：

  

  我会省去简化的麻烦，但做一些重新排列并使用欧拉的公式，我们得到：

  

  注意：

  

  我使用上面的属性来简化 . 对于2D情况，我们将定义meshgrid坐标，两个维度的范围都在 [-pi,pi] 之间，我们将在surf的曲面图上绘制它 . 请记住使用过滤器的 absolute 值来显示幅度：

  [w1,w2] = meshgrid(linspace(-pi,pi));
  h = abs(8 - 2*cos(w1+w2) - 2*cos(w1) - 2*cos(w2) - 2*cos(w1-w2));
  surf(w1,w2,h);
  title('2D spectrum of filter');
  

  w1 和 w2 提供在阵列的每个相应空间位置水平和垂直定义的频率的唯一组合 . w1 和 w2 实际上是二维的 . 对于 w1 和 w2 的每个独特组合，我们确定光谱的大小 . 一旦我们计算出这些，我们就可以将所有这些信息放在一个漂亮的三维表面图中 .

  我们得到：

  

  请注意，这两个维度都跨越 [-pi,pi] . 如果检查频谱，您将看到DC分量被取消为0，这实际上是一个高通滤波器，而 not 是一个锐化滤波器 . 请参阅下面的说明 .

  ---

  小调

  顺便说一句，您的2D滤镜定义是 not 2D锐化滤镜 . 它只是一个2D拉普拉斯算子，它是一个边缘检测，因此是一个高通滤波器 . 它只检测边缘 . 如果你想正确地锐化图像，请确保在内核的 centre 中添加1，所以你真正想要的是：

  h2 = [-1 -1 -1; -1 9 -1; -1 -1 -1];
  

  偏移量1将确保您保持原始信号不变，但也会在原始信号的顶部添加高通结果，从而锐化您的输入 .

  回复于 2024-05-01T13:05:35+08:00
• 3

  要绘制1D FIR滤波器的频率响应，请使用freqz：

  h = [-1, 3, -1]; %// FIR impulse response
  N = 128; %// number of frequency samples
  freqz(h,1,N); %// "1" because the filter is not recursive
  

  

  像往常一样，频率轴标准化为采样频率的一半 .

  有关更多选项，请参阅documentation of freqz .

  ---

  要绘制2D FIR滤波器的频率响应，请使用freqz2：

  h = [-1, -1, -1; -1, 8, -1; -1, -1, -1]; %// FIR impulse response (filter mask)
  N = 128; %// number of frequency samples
  freqz2(h,N,N)
  

  

  频率轴标准化为采样频率的一半，如1D情况 .

  有关更多选项，请参阅documentation of freqz2 .

  回复于 2024-05-01T13:05:35+08:00