使用 IIR 滤波器生成粉红噪声答案

【问题标题】：Pink Noise Generation with IIR filter使用 IIR 滤波器生成粉红噪声
【发布时间】：2012-07-03 15:44:40
【问题描述】：

我实现了来自 N. Jeremy Kasdin（第 825 页）pdf 的“IEEE 会议记录”中的以下代码。但我不明白这些行，因为我没有数字食谱书：

/* perform the discrete Fourier transform */  
realft (hfa,n_pts, 1); 
realft (wfa,n_pts, 1);

wfa[1]=wfa[1]*hfa[1]; 
wfa[2]=wfa[2]*hfa[2]; 

for(i=3;i<=nn;i+=2) { 
wr=wfa[i]; 
wi=wfa[i+1]; 
wfa[il=wr*hfa[i]-wi*hfa[i+1]; 
wfa[i+l]=wr*hfa[i+1]+wi*hfa[i];
}

谁能给我一些指示？

【问题讨论】：

数字食谱可在线获取nr.com
一题太贵了……
旧版本是免费的，对于这个特殊问题应该绰绰有余 - 转到nr.com/oldverswitcher.html

标签： java math signal-processing

【解决方案1】：

NR 中的realft 函数执行以下操作。你给它一个由 N 个实数组成的数组。（N 必须是 2 的幂。）它的离散傅里叶变换由 N 个遵循共轭对称关系的复数组成：F(k) 和 F(N-k) 是共轭的。特别是，F(0) 和 F(N/2) 是实数。所以realft返回N个实数，如下：F(0)，F(N/2)，F(1)的实部，F(1)的虚部，F(2)的实部，.. .，F(N/2-1)的虚部。

NR 最初都是 Fortran，并且（至少在旧版本中）使用基于 1 的索引而不是基于 0 的索引。即使在 C 语言中也是如此。这就是为什么代码从对元素 1 进行操作开始并以包含而不是排他的方式运行到 nn。

所以，您已经准备好 hfa 和 wfa 的 FT。剩下的代码只是简单地计算结果的元素乘积——前两行是简单的实数乘法，其余的是复数相乘。

我猜在那之后还有一个对realft 的调用，最后一个参数是-1（意思是进行逆运算）。所以整个事情是：hfa 和wfa 的FT；将它们逐元素相乘；逆FT。换句话说，卷积。

【讨论】：