65.9K
CodeProject 正在变化。 阅读更多。
Home

从 WAV 到图形再到 MIDI

Andy Allinger
starIconstarIconstarIconstarIconstarIcon

5.00/5 (23投票s)

2022 年 5 月 2 日

公共领域

21分钟阅读

viewsIcon

18977

downloadIcon

498

将音频转换为图像。（编辑图像。）将图像转换为音乐。

目录

引言

将录制音乐转录成乐谱是一个难题。音符很难从泛音和回声中分辨出来。另一方面，将独奏音频转换为 MIDI 则相当简单。

一种可能的方法是通过频谱图将音频转换为图像。然后手动编辑图像以消除泛音和回声，将最细微的步骤留给人脑。编辑图像的独奏线可以转换为 MIDI。 SPECTR 和 SCRIBE 程序旨在用作 MIDI 语音控制器，类似于键盘或呼吸控制器。

程序 SPECTR 在中音区生成高分辨率频谱图，共 32 个音符，从高音谱号下方的 F 到高音谱号上方的 C。

接下来，用户编辑图像以去除泛音。由于这是一个标准图像文件，因此在此步骤中可以执行任何可以对图像进行的编辑操作。

程序 SCRIBE 将编辑后的图像转换为标准 MIDI 文件。 SCRIBE 的输出由于表情和音高弯曲的许多变化而变长，每秒 240 帧。该程序的缺点是它只转录独奏，并且手动编辑图像工作量很大。

这是一个命令提示符程序。它是免费的，并且属于公共领域。

一些子程序本身就很有趣。子程序 FRAT 实现了频率重新分配变换，它比傅里叶变换提供更高的精度。 GLOM 将序列元素合并为具有相似值的片段。 SMOOTH 是一种非线性平滑器，用于从数据序列中去除异常值。

发行版中还包含一些优秀的遗留代码

  • SOLV 和 ZERORA 基于 L F Shampine 和 H A Watts 的 ZEROIN，用于求解方程。
  • Alfred H. Morris, Jr. 的 ISHELL 对整数数组进行排序。
  • Robert Renka 的 RORDER 将排列应用于浮点数组。

第一部分：构建更好的频谱图

声音的著名图形表示是频谱图。在该视图中，时间表示在水平方向上，频率表示在垂直方向上。图像中的每条垂直线都是声音简短部分的傅里叶变换。该线内的每个像素都是傅里叶变换的系数。

有关傅里叶变换的更多信息，请参阅 CodeProject 文章

恒定 Q 变换

离散傅里叶变换接受一个序列并返回一个系数数组（参见附录）。原始信号表示为频率为 \(f_o, 2f_o, \ldots, i f_o, \ldots, \frac{N}{2} f_o \) 的正弦波和余弦波之和，其中 \(f_o\) 是基频，即采样率除以变换中包含的样本数。设置傅里叶变换的设计方程是

$MR = NF$

其中

  • M 是一个整数
  • R 是每秒循环的采样率
  • N 是变换中的样本数
  • F 是要测量的频率

傅里叶分析的局限性在于，通常，时间分辨率和频率分辨率之间存在权衡。为了使频率更紧密地间隔开，\(f_o\) 必须更小。对于给定的采样率，这要求 N 更大。您寻求事件频率的精度越高，您确定事件发生时间的精度就越低。

傅里叶变换的本质是在一组等间隔频率处找到振幅。对于音乐分析而言，这不是最佳选择。音阶的音调根据几何级数间隔开。具体来说，每个音调的频率都比前一个音调大 \(\sqrt[12]{2}\) 倍。

图中的垂直条表示音乐音调频率的指数增长。水平线表示均匀间隔的傅里叶变换箱。足以分离高频音符的变换大小将无法分离低频音符。足以分离低频音符的变换大小对于高频音符具有过高的分辨率。必须选择大的变换大小。这会以时间分辨率降低和计算量增加为代价。

一个直接的解决方案 [Brown, 1991] 是采用不同大小的变换来观察不同的频率。这就是恒定 Q 变换，它不是一种新的数学技术，而是傅里叶变换的有针对性的使用。它允许在整个频谱上以最佳方式选择频率-时间分辨率权衡。在程序 SPECTR 中，M 选择为 17，因为 16:17 和 17:18 的比率近似等于音阶音调之间的间隔。

离散傅里叶变换 (DFT) 远不及快速傅里叶变换 (FFT) 快，无法找到从 \(1 \ldots \frac{N}{2}\) 的所有系数。但是 DFT 可用于查找单个系数，而 FFT 总是计算 N 个系数。许多 FFT 实现强制 N 为 2 的幂或小素数的乘积。DFT 对 N 没有限制。通过适当选择 M 和 N，DFT 允许测量任何一组频率，直至奈奎斯特极限 \(\frac{R}{2}\)。

窗函数

傅里叶分析假设要分析的信号是一个完美周期函数，其周期等于变换大小。要分析的真实信号从不满足此假设。为了减少引入的误差，需要通过将样本乘以窗函数来使要分析的样本的开始和结束部分逐渐减小。窗函数定义在 \(0 \ldots 2\pi\) 上，并在其端点处等于零，在中心处具有最大值。已经定义了许多窗函数，但它们的形状非常相似。

窗函数的选择对分析结果有明确的影响。在去除噪声和振幅测量的准确性之间存在权衡。这可以在[Heinzel, Rudiger, & Schilling, 2002]中描述的方法计算的一些窗函数在频域中的图中看出。

从视觉上看，函数之间的差异可以描述为中心瓣的平坦度或锐度之间的选择。汉宁窗的中心瓣相当窄，它可以更好地检测频率。Salvatore 窗的中心瓣平坦，它可以更准确地测量振幅。

Nuttall 窗在中心 bin 上具有 < 1dB 衰减的平坦度，这意味着它对于音频应用来说足够平坦，因为这大约是人耳的刚好可察觉的差异。为了不给用户带来过多的选项，该项目选择了 Nuttall 窗并将其并入子程序 FRAT。可以通过编辑源代码中的参数值来选择不同的窗函数。

窗函数的实际效果或许可以通过示例来展示。一个非常尖锐的窗函数，如 Hamming，具有良好的分辨率

而一个平顶窗可能相当模糊

然而，比较并非完全有利于尖锐窗。它们似乎会产生分叉伪影，即单个音高（错误地？）分裂成多个频率

平顶窗受这种效应的影响要小得多

频率重新分配

傅里叶变换获得的频率分辨率在一般情况下是最好的。然而，通常在频谱的某个区间内只存在一个频率分量。在这种情况下，可以避免时间-频率权衡，从而获得更好的频率分辨率。

其思想是，由于已知分量的正弦和余弦系数，因此可知合成波的相位。频率是相位的时导数。公式 [Lazzarini, 2011] 的推导需要复杂的分析，我不敢假装理解，但结果易于使用。令 \(u_j\) 表示窗函数，\(w_j\) 表示窗函数相对于索引 _j_ 的导数，\(y_j\) 表示输入信号。要找到能量 \(e\) 和频率 \(\hat{f}\)，计算

$\begin{aligned} a & = \frac{1}{N}\sum^N_{j=1} u_j y_j \cos(\frac{2\pi}{N}Mj) \\ b & = \frac{1}{N}\sum^N_{j=1} u_j y_j \sin(\frac{2\pi}{N}Mj) \\ c & = \frac{1}{N}\sum^N_{j=1} w_j y_j \cos(\frac{2\pi}{N}Mj) \\ d & = \frac{1}{N}\sum^N_{j=1} w_j y_j \sin(\frac{2\pi}{N}Mj) \\ e & = a^2 + b^2 \\ f & = \frac{MR}{N} \\ \hat{f} & = f + \left(\frac{ad-bc}{e}\right)\left(\frac{R}{2\pi}\right) \end{aligned}$

振幅是能量的平方根。

其他注意事项

为了获得最佳结果，输入信号应通过高通滤波器以去除低于基频的频率。为了提高效率，信号进行降采样以保持 N ≤ 6000。在降采样之前，信号应进行低通滤波至降低采样率的一半。数字滤波器已在我之前的文章 可自定义的巴特沃斯数字滤波器 中讨论过。

由于频段略有重叠，频率可能会乱序。插入排序对于几乎已排序的数据是高效的。子例程 ISORTR 被动地按频率排序。然后调用 Robert Renka 的子例程 RORDER 将排序应用于每个分析帧的振幅和频率。

插值

FRAT 返回的频率和振幅必须用于设置像素强度。这需要插值。有几种方法可以做到这一点，但 Lanczos 似乎被认为是最好的。该方程是

$ L(x) = \left\{ \begin{array}{cl} 1 & x = 0 \\ \frac{\sigma\sin(\pi x)\sin(\frac{\pi x}{\sigma})}{\pi^2 x^2} & -\sigma < x < +\sigma, x \neq 0 \\ 0 & \mbox{otherwise} \\ \end{array}\right\}$

L 是一个因子，它表示根据像素位置和估计频率（两者都以像素测量）之间的距离 x 来缩放振幅的程度。参数 σ 为 2 或 3。在 SPECTR 中选择 σ = 2 以避免 σ = 3 时出现的光晕效应。

最后，取振幅的对数，并将结果缩放到 0 到 255 的范围，即 PGM 文件中像素强度的范围。

插曲：编辑图形

学习阅读频谱图需要一些练习。水平线是单个频率的纯音

吉他等弦乐器产生稳定的频率。人声则不那么稳定。正如 Johnny Cash 的“Delia's Gone”的上述片段所示，声线并不那么直。垂直线是脉冲。它们听起来像咔哒声或砰砰声。*脉冲是所有频率同时进行的活动。*以下是使用 Sonic Visualiser 制作的手掌拍击的频谱图

模糊区域是静态噪声

许多乐器在基频之上产生泛音。这些表现为平行线

有时会有回声。这将在原始声音的右侧显示为较淡的线

识别旋律线，并消除泛音、回声、静电和咔嗒声。它不需要精细的触感。我推荐许多编辑程序中提供的黑色矩形工具。用纯黑色填充除了旋律之外的所有内容。如果您遗漏了某个点，程序 SCRIBE 可以处理一些零散的像素。

第二部分：从旋律线到 MIDI

撤销 Lanczos 插值

读取编辑后的 _ .PGM_ 文件后，必须将像素位置和强度转换回频率和振幅。在时间点只有两个活动像素的特殊情况下，可以精确地完成此操作。令 (x~o~,y~o~) 表示原始的未知频率和振幅。令 (x~1~,y~1~) 和 (x~2~,y~2~) 表示图像中的像素索引和强度，以线性光测量。在每个点应用 Lanczos 插值公式得到

$ y_1 = \frac{\sigma \sin(\pi(x_1-x_o))\sin(\frac{\pi}{\sigma}(x_1-x_o))}{\pi^2 (x_1-x_o)^2}y_o \\ y_2 = \frac{\sigma \sin(\pi(x_2-x_o))\sin(\frac{\pi}{\sigma}(x_2-x_o))}{\pi^2 (x_2-x_o)^2}y_o $

解出 _y~o~_

$ y_o = \frac{\pi^2 (x_1-x_o)^2}{\sigma\sin(\pi(x_1-x_o))\sin(\frac{\pi}{\sigma}(x_1-x_o))}y_1 = \frac{\pi^2(x_2-x_o)^2}{\sigma\sin(\pi(x_2-x_o))\sin(\frac{\pi}{\sigma}(x_2-x_o))}y_2 $

因此，问题已简化为求解一个未知数的方程

$ \frac{(x_1-x_o)^2}{\sin(\pi(x_1-x_o))\sin(\frac{\pi}{\sigma}(x_1-x_o))}y_1 - \frac{(x_2-x_o)^2}{\sin(\pi(x_2-x_o))\sin(\frac{\pi}{\sigma}(x_2-x_o))}y_2 = 0 $

子程序 solv 基于 Shampine 和 Watts 的 ZEROIN。在此处调用它以数值求解 _x~o~_。现在代入 _y~o~_ 的表达式

$ y_o = \frac{\pi^2(x_1-x_o)^2}{\sigma\sin(\pi(x_1-x_o))\sin(\frac{\pi}{\sigma}(x_1-x_o))} $

如果一个帧中有两个以上的像素处于活动状态，则调用子例程 QWM 作为备用。它将频率估计为加权平均值，并将振幅选择为最大强度。

分离声音和静音

接下来，必须为每个帧决定是否有音符发声。设置固定阈值不会产生好的结果。微小的亮点或暗点会变成许多短暂、不愉快的音符。声音和静音之间的转换只应在持续合理时间的情况下发生。动态规划方法为这类问题提供了最佳解决方案。[Temperley, 2001]

维特比算法是动态规划最简单的形式。它在序列的每一步，对于每个可能的事件，询问前一步最可能发生哪个事件？每个步骤的答案都被保留下来，在序列的末尾，可以追溯到最可能的事件链。子程序 DYRB 实现了维特比算法，用于二元状态（发声或静音）和实值观察序列（强度）的特殊情况。

我不知道将音符概率推断为像素强度函数的正确统计模型。正态分布似乎不适用于此应用。在正态模型下，以零强度和最大强度为极值，一半强度是声音或静音之间的无差别。这不是我们想要的；无差别值应该更接近零。相反，选择的惩罚函数是对强度除以参考值的对数，为声音或静音添加或减去。

通过选取一个任意的小时间单位来选择过渡惩罚。选择 0.035 秒的值以与程序 MELISMA [Sleator & Temperley, 2000] 中的值匹配，他们称之为“pip”单位。过渡惩罚等于最大观察惩罚的 pip 持续时间。

因此，动态规划的最优性声明必须附带强烈的免责声明。它是“最优的”，鉴于假定的观察惩罚和转换惩罚。

LULU 平滑器

接下来，使用非线性平滑器修剪异常值数据。传统上，移动平均线用于此目的。然而，它们是低通滤波器。如上所示，脉冲包含所有频率。因此，移动平均平滑器会模糊突然的变化。非线性平滑器可以在不模糊的情况下移除异常点。[Jankowitz, 2007] 子例程 SMOOTH 实现了 Jankowitz “B”滤波器。它完全基于最小和最大操作，并且不进行任何算术运算。

聚合音符

采用简单贪婪的方法合并音符，这将导致每个像素频率与片段平均值之间的平方差总和增加最少。这不是一个最优算法，但它高效且已被证明在实践中效果良好。[Terzi & Tsaparas, 2006] 随着每个片段的增长，其平均频率会发生变化，因此与相邻片段合并可能导致的误差变化也会发生变化。由于聚合是数据聚类的 1 维情况，因此 Ward 公式适用。[Theodoridis & Koutroumbas, 2009]

$ e_{ij} = \frac{n_i n_j}{n_i + n_j} d_{ij}$

其中 e 是误差，d 是片段平均值之间的差异，n 是每个片段中的对象数量。这些更新可以通过将潜在的误差变化放入堆中并将最小的变化筛选到顶部来高效处理。

子例程 GLOM 接受一个实值序列和三个参数作为输入。 ADPARM 是相邻片段平均值之间的阈值差异。 MDPARM 是片段内最大值和最小值之间范围的阈值。 K 是所需的片段数量。片段会一直合并，直到满足停止条件。即，每对相邻片段之间的差异都超过 ADPARM，或者在不超过 MDPARM 的情况下无法合并任何片段，或者已达到片段数量 K。

表达式公式

在程序 SPECTR 中，取傅里叶系数的对数，并将结果缩放至像素强度范围。根据 MIDI 标准，音符力度应具有指数解释。因此，力度可以直接与音符中的最大像素强度成正比，\( V = \frac{127}{255}Y \)，其中 Y 是 0 到 255 的强度，V 是 0 到 127 的力度。

MIDI 表达式的情况更困难。表达式默认最大，并设置为在每个音符内创建渐弱。请参阅源代码以获取结果的推导

$ E = 127 \sqrt{\epsilon}^{\scriptstyle (\frac{V}{127} - \frac{Y}{255})} $

其中 \(\epsilon = \frac{1}{1024}\)，是程序 SPECTR 中的分辨率参数。

未涵盖的主题

读写 _ .WAV_、_ .PGM_ 和 _ .MID_ 文件格式并不是很有趣。有关它们的信息，请参阅参考文献。

构建说明

使用 gcc 编译 C 代码，执行

gcc -o spectr -O2 rorder.c safdiv.c isortr.c search.c filter.c frat.c spectr.c -lm
gcc -o scribe -O2 solv.c ishell.c rorder.c safdiv.c isortr.c search.c sift.c glom.c dyrb.c smooth.c scribe.c -lm

使用 gfortran 编译 FORTRAN，执行

gfortran -o spectr -O2 -fno-automatic wrapper.f both.f rorder.f safdiv.f isortr.f search.f filter.f frat.f spectr.f
gfortran -o scribe -O2 -fno-automatic wrapper.f both.f zerora.f ishell.f safdiv.f search.f sift.f glom.f dyrb.f smooth.f scribe.f

这有点长，但没有依赖项，没有 makefile，也没有头文件。

Using the Code

用法是

     spectr something.wav something.pgm
     scribe something.pgm something.mid

其中

  • _ .wav_ 是标准 Microsoft 波形音频
  • _ .pgm_ 是便携式灰度图图像文件
  • _ .mid_ 是您的输出

您可以对图像进行任何类型的变换，但其高度必须保持 640 像素。请记住，SCRIBE 仅转录独奏。

这两个程序都有一个 --tune 选项来移调。所以你可以做一些事情，比如

     spectr --tune +24.5 deep.wav alto.pgm
     scribe --tune -24.5 alto.pgm deep.mid

这将男中音向上移两个八度加一个四分音符到频谱图的范围，然后在写入 MIDI 时将其移调回来。

程序 SCRIBE 具有 --patch、--text、--copyright 和 --seqname 选项，因此

     scribe --patch 68 --seqname "Hacker's Blues" --copyright "(C)2021 Joe" 
            --text "nobody reads the comments" blues.pgm blues.mid

将音色设置为双簧管，为音轨命名，添加版权声明和注释。

程序 SPECTR 有一个 --swab 选项，如果您的 _ .WAV_ 是一个 32 位浮点数，并且字节顺序与计算机的默认值相反，您可能需要它。

示例问题

人声

请看这个由 Katy 演唱的样本，由 digifishmusic 提供。[Katy, 2007]

普通的频谱图（用 sox 制作）相当模糊。

Spectrogram made with sox

使用 SPECTR 制作的高分辨率频谱图显示了颤音的细节。

接下来，使用 Windows 画图或任何类似的程序编辑图像以去除泛音。

擦除除一条平行线之外的所有线条

因为这是一个标准图像文件，所以在此步骤中可以执行任何可以对图像进行的编辑操作。现在使用 SCRIBE 处理图像以生成 MIDI 文件。

作为比较，ofoct.com 上的免费文件转换器提供的结果可用。此服务免费，无需在您的计算机上安装任何东西。两个 MIDI 都准确捕捉了演奏。竞争对手的 MIDI 非常简洁。另一方面，它包含 167 个音符。 SCRIBE 的输出由于表达方式的多次变化而变得很长，但只有 28 个音符。从这个意义上说，它是一个更准确的转录。

器乐

长笛演奏 [kerri, 2007] 由 SCAT 系统处理。频谱图的一个细节

编辑图像后

结果是相当不错的转录，尽管一些装饰音已经丢失。

发布历史

  • 3 2023 年 8 月 12 日：聚合寻求最小平方误差。
  • 2 2022 年 4 月 30 日：非线性平滑，逆插值，无 _ .exe_ 文件
  • 1 2021 年 10 月 10 日：添加 readme、FORTRAN 源代码和 _ .exe_ 文件
  • 0 2021 年 10 月 6 日：首次发布到 MidKar 小组

参考文献

  1. Jef Poskanzer 撰写的“PGM 格式规范”，1989 年。在线
  2. Judith C. Brown 撰写的“恒定 Q 频谱变换的计算”，《美国声学学会杂志》第 89 卷第 1 期，1991 年 1 月。可在 researchgate 上获取。
  3. MIDI 制造商协会撰写的“标准 MIDI 文件规范”，1996 年。在线
  4. Craig Stuart Sapp 撰写的“Microsoft WAVE 音频文件格式”，1997 年。在线
  5. 基本音乐结构的认知，David Temperley，MIT Press，2001 年。出版商提供第一章。
  6. G. Heinzel、A. Rudiger 和 R. Schilling 撰写的“通过离散傅里叶变换 (DFT) 进行频谱和频谱密度估计，包括窗函数的综合列表和一些新的平顶窗”，Max-Planck-Institut fur Gravitationsphysik，2002 年 2 月 15 日。可在 researchgate 上获取。
  7. Robert Bristow-Johnson 撰写的“音频 EQ 双二阶滤波器系数的食谱公式”，[2005？]。可在 musicsdp 上获取。
  8. John Stensby 撰写的“巴特沃斯低通滤波器”，2005 年 10 月 19 日。可从阿拉巴马大学获取。
  9. Evimaria Terzi 和 Panayiotis Tsaparas 撰写的“序列分割的有效算法”，《第六届 SIAM 国际数据挖掘会议论文集》，2006 年 4 月。可在 researchgate 上获取。
  10. LULU 平滑器的一些统计方面，Maria Dorothea Jankowitz，论文，斯泰伦博斯大学，2007 年 12 月。可从 斯泰伦博斯大学 获取。
  11. kerri 撰写的“闹鬼峡谷长笛”，2007 年。在 freesound。
  12. Katy 演唱的“Katy Sings LaaOooAaa”，由 digifishmusic 录制，2007 年。在 freesound。
  13. 模式识别 Sergios Theodoridis 和 Konstantinos Koutroumbas，Academic Press / Elsevier，2009 年。
  14. Bear 文件转换器, 2009.
  15. Victor Lazzarini 撰写的“相位声码器编程”，收录于 Richard Boulanger 和 Victor Lazzarini 编辑的《音频编程书》。麻省理工学院出版社，2011 年。

附录：离散傅里叶变换的简短推导

假设在某个持续时间 T 内，以 N 个等间隔的间隔 \((y_0, y_1, y_2, \ldots, y_j, \ldots, y_n) \) 测量离散信号。假设可以将此信号表示为正弦波和余弦波之和

$ y_j = \sum^{N/2}_{i=1} (a_i \cos( 2 \pi f_i t_j) + b_i \sin(2 \pi f_i t_j)) $

对于 \(i = 1 \ldots \frac{N}{2} \)，代入 \(t_j = \frac{j}{N}T\) 和 \(f_i = i f_o = \frac{i}{T}\)，因此

$ y_j = \sum^{N/2}_{i=1} (a_i \cos(\frac{2\pi}{N}ij) + b_i \sin(\frac{2\pi}{N}ij)) $

写出表格

$ \left[ \begin{array}{llllllllll} \cos(\frac{2\pi}{N}) & \sin(\frac{2\pi}{N}) & \cos(\frac{2\pi}{N}(2)) & \sin(\frac{2\pi}{N}(2)) & \cdots & \cos(\frac{2\pi}{N}i) & \sin(\frac{2\pi}{N}i) & \cdots & \cos(\pi) & \sin(\pi) \\ \cos(\frac{2\pi}{N}(2)) & \sin(\frac{2\pi}{N}(2)) & \cos(\frac{2\pi}{N}(2)(2)) & \sin(\frac{2\pi}{N}(2)(2))& \cdots & \cos(\frac{2\pi}{N}(2)i) & \sin(\frac{2\pi}{N}(2)i) & \cdots & \cos(2\pi) & \sin(2\pi) \\ \vdots & \vdots &\vdots &\vdots & \ddots & \vdots & \vdots & \cdots & \vdots & \vdots \\ \cos(\frac{2\pi}{N}j) & \sin(\frac{2\pi}{N}j) & \cos(\frac{2\pi}{N}(2)j) & \sin(\frac{2\pi}{N}(2)j) &\cdots & \cos(\frac{2\pi}{N}ij) & \sin(\frac{2\pi}{N}ij) & \cdots & \cos(\pi j) & \sin(\pi j) \\ \vdots & \vdots &\vdots & \vdots& & \vdots & \vdots & \ddots &\vdots & \vdots \\ \cos(2\pi) & \sin(2\pi) & \cos(2\pi(2)) & \sin(2\pi(2)) & \cdots & \cos(2\pi i) & \sin(2\pi i) & \cdots & \cos(\pi N) & \sin(\pi N) \\ \end{array} \right] \left[ \begin{array}{l} a_1 \\ b_1 \\ a_2 \\ b_2 \\ \vdots \\ a_i \\ b_i \\ \vdots \\ a_\frac{N}{2} \\ b_\frac{N}{2} \\ \end{array} \right] = \left[ \begin{array}{l} y_1 \\ y_2 \\ \vdots \\ y_j \\ \vdots \\ y_N \\ \end{array} \right] $

正弦和余弦项形成一个方阵，因为有 \(\frac{N}{2}\) 对 \((a_i, b_i)\)。它也是一个正交矩阵，所以它的逆矩阵与其转置相同。请注意，无法获得超出 \(\frac{N}{2}\) 的频率，因为那样未知数将多于方程。因此，

$ \left[ \begin{array}{l} a_1 \\ b_1 \\ a_2 \\ b_2 \\ \vdots \\ a_i \\ b_i \\ \vdots \\ a_\frac{N}{2} \\ b_\frac{N}{2} \\ \end{array} \right] = \left[ \begin{array}{llllll} \cos(\frac{2\pi}{N}) & \cos(\frac{2\pi}{N}(2)) & \cdots & \cos(\frac{2\pi}{N}j) & \cdots & \cos(2\pi) \\ \sin(\frac{2\pi}{N}) & \sin(\frac{2\pi}{N}(2)) & \cdots & \sin(\frac{2\pi}{N}j) & \cdots & \sin(2\pi) \\ \vdots & \vdots & \ddots & \vdots & & \vdots \\ \cos(\frac{2\pi}{N}i) & \cos(\frac{2\pi}{N}i(2)) & \cdots & \cos(\frac{2\pi}{N}ij) & \cdots & \cos(2\pi i) \\ \sin(\frac{2\pi}{N}i) & \sin(\frac{2\pi}{N}i(2)) & \cdots & \sin(\frac{2\pi}{N}ij) & \cdots & \sin(2\pi i) \\ \vdots & \vdots &\vdots &\vdots & \ddots & \vdots \\ \cos(\pi) & \cos(2\pi) & \cdots & \cos(\pi j) &\cdots & \cos(\pi N) \\ \sin(\pi) & \sin(2\pi) & \cdots & \sin(\pi j) & \cdots & \sin(\pi N) \\ \end{array} \right] \left[ \begin{array}{l} y_1 \\ y_2 \\ \vdots \\ y_j \\ \vdots \\ y_N \\ \end{array} \right] $

用求和符号表示

$ a_i = \sum^N_{j=1} \cos(\frac{2\pi}{N} ij) y_j\\ b_i = \sum^N_{j=1} \sin(\frac{2\pi}{N} ij)y_j $

对于 \(i = 1 \ldots \frac{N}{2} \)。由于采样率 \(R = \frac{N}{T}\)，这些系数对应于频率 \(f_i = i \frac{R}{N}\)。

为了使振幅独立于变换的长度，必须将系数除以 N，从而得出以下公式

$ a_i = \frac{1}{N}\sum^N_{j=1} \cos(\frac{2\pi}{N} ij) y_j\\ b_i = \frac{1}{N}\sum^N_{j=1} \sin(\frac{2\pi}{N} ij)y_j $

请注意，离散傅里叶变换是在不使用任何虚数的情况下获得的。

© . All rights reserved.