我不明白评论。
你当然可以这样做。只需了解 DFT 的含义、如何计算 DFT bin 值以及如何将这些 bin 值解释为连续傅立叶级数系数。
首先,您正在查看的平面是复平面。您的点是一组$N$离散样本。每个样本都是一个复杂的点。因此,您所拥有的是重复复杂信号的一个周期的表示。图中的间距不是那么重要。
- 任何$N$个点的序列都可以用N个系数精确地在样本点处表示
问题是:“你的形状是否允许你忽略系数,所以你的系数要少得多?”
答案是:“取决于形状。” 因此,开始丢弃最小的幅度系数,看看精度会受到多大影响。
当您构建傅里叶级数时,您需要将未归一化的 DFT 系数除以$N$。您还希望将 DFT 的上半部分算作负频率,因此$N-1$对应于$-1$等。
所以基本上你是在对离散序列进行 DFT,然后使用系数重建插值。
希望这可以帮助。
赛德
我把这个放在我的答案中,因为我不想触发将此对话转移到聊天空间(我不同意 BTW 的政策)。
手头的问题不仅仅是一个封闭的图形是否可以参数化,这个问题被牢固地设置为本轮求和的应用(你知道,在哥白尼改变参考系之前,行星运动是如何被建模的)。是的,还有其他方法可以参数化圆周运动而不是正弦和余弦,但它们很笨拙。
还有其他方法可以退后一步并将图形作为一个整体参数化,并且不需要它是周期性的。我想到了勒让德基础。DFT 方法就是这样发生的,它本质上是周期性的。
在我看来,OP 认为可以绘制任意图形(在限制内)很酷(我也是),并试图理解本轮的概念与 DFT 的关系。
让我们做一些数学运算以使其更清楚。使用传统的归一化和符号,DFT 是:
$$ X[k] = \sum_{n=0}^{N-1} x[n] e^{-ink \frac{2\pi}{N} } $$
由于$x[n]$是已知的,$X[k]$现在也是已知的。现在,让我们看一下反面:
$$ x[n] = \frac{1}{N} \sum_{k=0}^{N-1} X[k] e^{ink \frac{2\pi}{N} } $$
如果我们简单地允许$n$为实值并将逆 DFT 定义视为一个连续方程,我们会遇到超过奈奎斯特频率的问题。在离散情况下没有区别,因为它们将在样本点处匹配。在这两者之间,确实如此。因此,总和必须移动到以 DC bin 为中心。(假设$N$是偶数,否则可以类似地计算出来)
$$ x[n] = \frac{1}{N} \sum_{k=-N/2+1}^{N/2} X[k] e^{ink \frac{2\pi}{N } } $$
该方程也可以分为实部和虚部:
$$ \Re(x(n)) = \sum_{k=-N/2+1}^{N/2} X[k] \frac{1}{N} \cos( nk \frac{2\ pi}{N} ) $$
$$ \Im(x(n)) = \sum_{k=-N/2+1}^{N/2} X[k] \frac{1}{N} \罪(NK \frac{2\pi}{N})$$
我会使用$x$和$y$,但$x$已经被占用了。这些方程显然是傅里叶级数的形式,系数为$ X[k] / N $。
我不是想在这里教育 r bj,我知道他对这些东西了如指掌。我只是说引入替代参数化或替代坐标系会分散手头的核心问题。
是的,刚刚接受 r bj 的教育。也感谢您的编辑。
事实上,奈奎斯特项应该分成两半,结果是这两个本轮将相互抵消虚部并加倍实部。由于图中没有很多曲折,我预计这个系数的大小会很低。
只是为了好玩,我写了一个小 Gambas 程序来演示数学。你可以在这里找到它:
https://forum.gambas.one/viewtopic.php?f=4&t=725
我还做了一点手绘的 Pi 符号。当然,它看起来有点醉了,但它仍然证明了这一点。
根据要求,这里有一点角落处理。角落的工作比我预期的要好。我认为这个例子真正体现了我之前所说的真正有趣的问题是沿着图形找到最接近的点放置。
在这个讨论中我没有明确说的是$ e^{i\theta} $的复数值沿着复单位圆移动,因此是一个本轮模型,所以循环内的每个产品代表如果要制作本轮动画,则在该时间点其各自本轮(即线段)的半径位置。半径的长度是系数的大小,因为$ e^{i\theta} $的大小始终为 1。
Complex.Polar(r, theta) = r * e ^ { i theta }
您可能会发现我的这篇文章有助于理解这些材料:
我不是 MATLAB 的粉丝(主要是因为对基于数组的使用极其短视),所以我不会评论您的伪代码。相反,这是我的代码,它实际计算给定“n”处的插值。
你可以点击链接自己下载(我只是放了新版本,允许在同一张图中有多个图形)。如果你有 Linux,你可以安装 Gambas (PPA:gambas-team/gambas3) 来运行它。
[注:代码中的n只是一个迭代器,t是实际的n,我懒得编辑代码。]
.
.
.
对于 n = 0 到 myPoints.Count * 100 - 1
t = n / 100
p = 计算(t, w)
Paint.Arc(p.Real, p.Imag, 1)
Paint.Fill()
下一个
.
.
.
'================================================== =======================
Public Sub Calculate(ArgT As Float, ArgDFT As Vector) As Complex
Dim k, N 作为整数
暗淡如复杂
将 a、b 调暗为浮点数
N = ArgDFT.Count
b = ArgT * Pi(2) / N
如果偶数(N)然后
GoSub EvenCase
别的
GoSub OddCase
万一
返回 p
'---------------------------------------------------- ----------------------
偶数情况:
p = ArgDFT[0] + ArgDFT[N / 2] * Cos(ArgT * Pi)
对于 k = 1 到 N / 2 - 1
a = b * k
p += ArgDFT[k] * Complex.Polar(1.0, a)
p += ArgDFT[N - k] * Complex.Polar(1.0, -a)
下一个
返回
'---------------------------------------------------- ----------------------
奇案:
p = ArgDFT[0]
对于 k = 1 到 (N - 1) / 2
a = b * k
p += ArgDFT[k] * Complex.Polar(1.0, a)
p += ArgDFT[N - k] * Complex.Polar(1.0, -a)
下一个
返回
结尾
'================================================== =======================
对 r bj 的反驳:
罗伯特,我非常不同意你的一些主张。
1) 插值点(以及它们通过 LineTo 调用形成的路径)将遵循您提供点的任何顺序
2)与我的答案相比,“t”使用 0 到$2\pi$范围会混淆问题,其中“t”的范围从 0 到 N,即与离散尺度相同的尺度,仅包括介于两者之间的实际值整数。这反映了您的参考框架是连续案例。[不再相关,我已经放弃了 t]
3) 将 (x,y) 视为一个向量,而不是单个复数值 x + iy,将参数化分为两个独立的问题,这两个问题不必通过相同的方法进行参数化。只有在复杂的价值解释中,作为这个问题核心的本轮概念才有意义。
4) 将拐角处的点捆绑在一起,中间没有足够数量的点,会导致拐角处的超限。这就是我在上一张图中添加第四个数字的原因。
5)您对$a_k$和$b_k$的定义没有意义,因为没有给出连续函数,只有一组样本点。因此,傅立叶系数应该使用离散定义来计算,即求和而不是积分。你本末倒置。使用不同的参数化,例如勒让德,您不会有超出范围的重复模式,不一定会在点之间匹配,但您会在所有样本点处匹配。
标题问题是:“如何获得傅立叶系数以使用 DFT 绘制任何形状?”
答案是:“归一化的 DFT bin 值是傅立叶系数。”
换句话说,简单地将中心逆 DFT 中的整数“n”替换为连续实值变量将产生插值结果。没有比这更优雅的了。我的代码就是对此的一种表达。我假设 OP 将在 MATLAB 中实现它(进行必要的索引调整)。
你让这种方式比它需要的更复杂。
这是偶数情况循环的等效编码,以阐明“k”的含义及其范围。
对于 k = -N / 2 + 1 到 -1
a = b * k
p += ArgDFT[k+N] * Complex.Polar(1.0, a)
下一个
对于 k = 1 到 N / 2 - 1
a = b * k
p += ArgDFT[k] * Complex.Polar(1.0, a)
下一个
这个是给 Olli 的,使用 N = 9。如果图形确实是一个三角形,你可以看到,通过适当的点放置,可以找到一个更好的拟合,它也能击中所有的点。当然,可以添加更多点(更多本轮)来获得更紧密的拟合。
这里的辅助问题(前面提到过),我认为这将是你要解决的问题,Olli,是如何将样本点放置在底层连续图形上以最小化“涟漪”或“溢出”。
当您将 DFT 的上半部分视为正频率而不是负频率时,就会发生这种情况。您可以清楚地看到所有点都被击中,但在结果之间并不是我们想要的。也许有一些新颖的应用程序会有所帮助。
我这样做是为了回应 Olli 提出的是否只能用正频率来完成的挑战。也许,如果将实部和虚部分开并且将余弦级数用于两个参数化,但我认为这会阻碍问题的意图,并且它不再是本轮实现。
我最初的倾向是拒绝。我认为这个问题相当于“你能从顺时针开瓶器的总和中构造一个逆时针开瓶器吗?” 也许有一个无穷大的数字,我在数学中看到了太多关于接近无穷大的奇怪事情来排除它,但在这里我什至无法想象一个近似的序列。
克里斯,
我不反对变量的大写字母。事实上,我喜欢使用$S_n$作为我的信号值。我不喜欢对信号使用小写“x”,对 DFT 使用大写“X”。对我来说,这还不够,因为它们描述的是两个完全不同的领域。此外,“X”是小写和大写最相似的字母之一,这使得它们在手写数学中更难区分。
在这种情况下,我们基本上有三个尺度(或功能域):
1) n 在样本点的整数上从 0 变为 N-1(对于输入点和逆 DFT 的输出)
2) k 在逆 DFT 定义中的整数上从 0 变为 N-1,然后移动半帧以消除上半部分的锯齿
3) t 从 0 变为$2\pi$是系列解决方案的域(你和罗伯特)和(0 到 N - 我的代码中的 1/100)
所以,是的,你在 T 域中使用 K 是在误导。
在我的代码中,ArgDFT 是 1/N 归一化 DFT,而 ArgT 是我原来的“t”参数,它与“n”的比例相同,但是是连续的。我在代码中的“b”和你的“t”一样。
流程总结:
当您对样本序列进行 1/N 归一化 DFT 时,您会同时找到将通过所有点的连续傅里叶级数的系数。(为什么应该使用 1/N 归一化的一个强有力的论据)。
级数解的域可以通过变量替换重新缩放:
$$ n = t \cdot \frac{N}{2\pi} $$
转化为被解释为连续函数的逆 DFT。
$$ x(n) = \frac{1}{N} \sum_{k} X[k] e^{ink \frac{2\pi}{N} } $$
$$ z(t) = x(t \cdot \frac{N}{2\pi}) = \frac{1}{N} \sum_{k} X[k] e^{it \cdot \frac{ N}{2\pi} k \frac{2\pi}{N} } $$
$$ z(t) = \sum_{k} \frac{X[k]}{N} e^{ikt } $$
这就是连续插值路径的级数解。它只是 t 的函数。如果需要,您可以将其与 t 区分开来计算您的“笔速度”。
很明显,您现在已经理解了我所说的“点放置问题”的含义,而且看起来 Olli 也对它产生了兴趣。
如果您还没有,我建议您重新阅读此线程中的所有内容。有了更好的理解基础,已经说过的话应该更有意义。
结语:关于情况的不同观点,这里的许多人都熟悉。但是,它不会产生傅立叶系数。
本来就是一堆烂摊子。
$$ z(t) = \sum_{k} \frac{X[k]}{N} e^{ikt } $$
$$ X[k] = \sum_{n=0}^{N-1} x[n] e^{-ink \frac{2\pi}{N} } $$
$$ z(t) = \frac{1}{N} \sum_{k} \sum_{n=0}^{N-1} x[n] e^{-ink \frac{2\pi}{ N} } e^{ikt } $$
$$ z(t) = \frac{1}{N} \sum_{n=0}^{N-1} x[n] \sum_{k} e^{ik ( t - \frac{n}{ N}2\pi )} $$
$$ t_n = t - \frac{n}{N}2\pi $$
$$ D(t_n) = \sum_{k} e^{ik t_n } $$
$$ z(t) = \frac{1}{N} \sum_{n=0}^{N-1} x[n] D(t_n) $$
奇数情况:$k = -(N-1)/2 \to (N-1)/2$
令$l = k + (N-1)/2$变为$0 \到 N-1$
$$ k = l - (N-1)/2 $$
$$ \begin{aligned} D(t_n) &= \sum_{l=0}^{N-1} e^{i ( l - (N-1)/2 ) t_n } \\ &= \sum_{ l=0}^{N-1} e^{il t_n } e^{-i \frac{N-1}{2} t_n } \\ &= e^{-i \frac{N-1}{ 2} t_n} \sum_{l=0}^{N-1} (e^{i t_n })^l \\ &= e^{-i \frac{N-1}{2} t_n} \frac {1 - (e^{i t_n })^N }{ 1 - e^{i t_n } } \\ &= e^{-i \frac{N-1}{2} t_n} \left[ \frac {e^{i t_n N / 2 } } { e^{i t_n / 2 } } \cdot \frac{ e^{-i t_n N / 2 } - e^{i t_n N/2 } }{ e^ {-i t_n / 2 } - e^{i t_n / 2 } } \right] \\ &= \frac{e^{i t_n N / 2 } - e^{-i t_n N / 2 }} { e ^{i t_n / 2 } - e^{-i t_n / 2 } } \\ &= \frac{ 2i \cdot \sin( N t_n / 2 ) } { 2i \cdot \sin( t_n / 2 ) } \ \ &= \frac{ \sin( N t_n / 2 ) } { \sin( t_n / 2 ) } \end{对齐} $$
$$ z(t) = \frac{1}{N} \sum_{n=0}^{N-1} x[n] \frac{ \sin( N t_n / 2 ) } { \sin( t_n / 2) } $$
$$ z(t) = \sum_{n=0}^{N-1} x[n] \frac{ \sin( N (t - \frac{n}{N}2\pi) / 2 ) } { N \sin( (t - \frac{n}{N}2\pi) / 2 ) } $$
请注意,商是实数值,因此可以将其视为权重,然后求和是样本点集的时变加权平均值。
结语二
在对其他问题进行了大量讨论之后,很明显 Nyquist bin 应该在负频率和正频率之间平均分配。
偶数情况:$k = 1/2 ( N/2 \text{ and } -N/2 ), -N/2 + 1 \to N/2 - 1 $
让$l = k + N/2 - 1 $去$0 \到 N-2$
$$ k = l - N/2 + 1 $$
$$ \begin{aligned} D(t_n) &= \frac{1}{2} \left[ e^{i ( N/2 ) t_n } + e^{i (-N/2 ) t_n } \right ] + \sum_{l=0}^{N-2} e^{i ( l - N/2 + 1 ) t_n } \\ &= \cos \left( \frac{N}{2} t_n \right ) + \sum_{l=0}^{N-2} e^{il t_n } e^{i (- N/2 + 1 ) t_n } \\ &= \cos \left( \frac{N}{ 2} t_n \right) + e^{i (- N/2 + 1 ) t_n } \sum_{l=0}^{N-2} (e^{i t_n })^l \\ &= \cos \left( \frac{N}{2} t_n \right) + e^{i (- N/2 + 1 ) t_n } \frac{1 - (e^{i t_n })^{N-1} } { 1 - e^{i t_n } } \\ &= \cos \left( \frac{N}{2} t_n \right) + e^{i (- N/2 + 1 ) t_n } \left[ \ frac{e^{i t_n ( N - 1 ) / 2 } } { e^{i t_n / 2 } } \cdot \frac{ e^{-i t_n ( N - 1 ) / 2 } - e^{i t_n ( N - 1 ) / 2 } }{ e^{-i t_n / 2 } - e^{i t_n / 2 } } \right] \\ &= \cos \left( \frac{N}{2} t_n \right) + \frac{e^{i t_n ( N - 1 ) / 2 } - e^{-i t_n ( N - 1 ) / 2 } } { e^{i t_n / 2 } - e^{-i t_n / 2 } } \\ &= \cos \left( \frac{N}{2} t_n \right) + \frac{ 2i \cdot \ sin( t_n ( N - 1 ) / 2 ) } { 2i \cdot \sin( t_n / 2 ) } \\ &= \cos \left( \frac{N}{2} t_n \right) + \frac{ \ sin( t_n N /2 ) \cos( t_n / 2 ) - \cos( t_n N /2 ) \sin( t_n / 2 ) } { \sin( t_n / 2 ) } \\ &= \cos \left( \ frac{N}{2} t_n \right) + \frac{ \sin( t_n N /2 ) } { \sin( t_n / 2 ) } \cos( t_n / 2 ) - \cos( t_n N /2 ) \ \ &= \frac{ \sin( N t_n/2 ) }{ \sin( t_n / 2 ) } \cos( t_n / 2 ) \end{对齐} $$= \cos \left( \frac{N}{2} t_n \right) + \frac{ \sin( t_n N /2 ) \cos( t_n / 2 ) - \cos( t_n N /2 ) \sin( t_n / 2 ) } { \sin( t_n / 2 ) } \\ &= \cos \left( \frac{N}{2} t_n \right) + \frac{ \sin( t_n N /2 ) } { \sin ( t_n / 2 ) } \cos( t_n / 2 ) - \cos( t_n N /2 ) \\ &= \frac{ \sin( N t_n/2 ) }{ \sin( t_n / 2 ) } \cos( t_n / 2 ) \end{对齐} $$= \cos \left( \frac{N}{2} t_n \right) + \frac{ \sin( t_n N /2 ) \cos( t_n / 2 ) - \cos( t_n N /2 ) \sin( t_n / 2 ) } { \sin( t_n / 2 ) } \\ &= \cos \left( \frac{N}{2} t_n \right) + \frac{ \sin( t_n N /2 ) } { \sin ( t_n / 2 ) } \cos( t_n / 2 ) - \cos( t_n N /2 ) \\ &= \frac{ \sin( N t_n/2 ) }{ \sin( t_n / 2 ) } \cos( t_n / 2 ) \end{对齐} $$
对于正负奈奎斯特项,上述推导可以使用 1/2 和 1/2 以外的系数来完成,但最后发生的简化不会发生,表达式会更复杂。同样非常清楚的是,如果$x[n]$的集合是真实的,则插值不一定是真实的。对于 1/2 和 1/2,插值都是实数。
则连续插值函数为:
$$ \begin{aligned} z(t) &= \frac{1}{N} \sum_{n=0}^{N-1} x[n] \left[ \frac{ \sin( N t_n / 2 ) }{ \sin( t_n / 2 ) } \right] \cos( t_n / 2 ) \\ &= \sum_{n=0}^{N-1} x[n] \left[ \frac{ \ sin( N (t - \frac{n}{N}2\pi) / 2 ) } { N \sin( (t - \frac{n}{N}2\pi) / 2 ) } \right] \ cos( (t - \frac{n}{N}2\pi) / 2 ) \\ &= \sum_{n=0}^{N-1} x[n] \frac{ \sin( N (t - \frac{n}{N}2\pi) / 2 ) } { N \tan( (t - \frac{n}{N}2\pi) / 2 ) } \end{aligned} $$
值得注意的是,该公式将奇数情况版本与前两行中应用的简单“窗口函数”相匹配。最后一个匹配 R BJ 的给定公式,格式更简洁。
查看 N = 2 的情况
$$ \begin{aligned} z(t) &= x[0] \left[ \cos^2( t / 2 ) \right] + x[1] \left[ \cos^2( (t - \pi ) / 2 ) \right] \\ &= x[0] \left[ \frac{ \cos( t ) + 1 }{2} \right] + x[1] \left[ \frac{ \cos( t - \pi ) + 1 }{2} \right] \\ &= \frac{1}{2} ( x[0] + x[1] ) + \frac{1}{2} ( x[0] - x[1] ) \cos( t ) \end{对齐} $$
这意味着交替序列 1、-1、1、-1 被插值为:
$$ z(t) = \cos(t) $$
这是奈奎斯特分裂假设的直接结果。