计算傅立叶级数以表示数据点
Computing Fourier Series to represent data points
我想计算一个函数（傅里叶级数），它通过一组给定的点。
与这里发生的事情类似 https://gofigure.impara.ai/ ，但我不希望对其进行动画处理。我只想要这个功能，这样我就可以自己画形状了。我已经阅读了很多描述它的数学资料和动画它的代码，但我正在为我的实现而苦苦挣扎。
我现在的代码如下[应该可以运行单独在python笔记本中]
import matplotlib.pyplot as plt
import scipy
import cmath
import math
import numpy as np
from matplotlib import collections as mc
import pylab as pl
midpoints = [[305.0, 244.5], [293.5, 237.0], [274.5, 367.5], [270.5, 373.5], [229.5, 391.0], [216.0, 396.0], [302.0, 269.0], [295.0, 271.0], [60.5, 146.5], [54.0, 153.0], [52.0, 167.0], [54.0, 153.0], [52.0, 167.0], [45.0, 178.0], [75.0, 76.5], [68.5, 98.5], [75.0, 76.5], [97.0, 58.5], [283.5, 357.5], [274.5, 367.5], [309.0, 255.0], [305.0, 244.5], [309.0, 255.0], [302.0, 269.0], [299.5, 291.5], [300.0, 297.0], [300.0, 297.0], [295.0, 309.5], [62.5, 105.0], [61.0, 118.5], [62.5, 105.0], [68.5, 98.5], [58.0, 139.0], [60.5, 146.5], [241.0, 111.5], [252.0, 124.5], [256.0, 132.5], [252.0, 124.5], [283.0, 356.0], [283.5, 357.5], [300.0, 290.5], [299.5, 291.5], [300.0, 290.5], [296.0, 280.5], [296.0, 280.5], [295.0, 271.0], [158.0, 387.0], [177.0, 396.5], [197.5, 402.0], [192.5, 403.0], [189.5, 400.5], [192.5, 403.0], [197.5, 402.0], [202.5, 401.0], [214.0, 395.5], [216.0, 396.0], [202.5, 401.0], [214.0, 395.5], [233.5, 375.0], [229.5, 391.0], [233.5, 375.0], [249.0, 372.5], [282.5, 340.0], [284.5, 328.0], [284.5, 328.0], [295.0, 309.5], [45.0, 178.0], [49.5, 189.0], [57.0, 198.0], [49.5, 189.0], [238.5, 108.5], [241.0, 111.5], [162.0, 57.5], [170.0, 59.0], [239.5, 204.5], [239.0, 200.0], [293.5, 237.0], [291.0, 227.5], [265.0, 229.5], [291.0, 227.5], [239.0, 189.5], [245.5, 178.0], [239.0, 189.5], [241.0, 193.5], [241.0, 193.5], [239.0, 200.0], [55.0, 119.0], [61.0, 118.5], [53.5, 134.0], [58.0, 139.0], [50.0, 129.0], [55.0, 119.0], [50.0, 129.0], [53.5, 134.0], [107.0, 46.0], [119.5, 50.5], [97.0, 54.0], [97.0, 58.5], [107.0, 46.0], [97.0, 54.0], [150.5, 377.0], [158.0, 387.0], [257.5, 367.5], [270.5, 373.5], [249.0, 372.5], [257.5, 367.5], [280.0, 349.5], [282.5, 340.0], [280.0, 349.5], [283.0, 356.0], [239.5, 90.0], [238.5, 98.0], [238.5, 108.5], [238.5, 98.0], [130.0, 49.0], [119.5, 50.5], [189.0, 65.0], [191.0, 64.5], [189.0, 65.0], [177.0, 62.5], [170.0, 59.0], [177.0, 62.5], [256.0, 132.5], [257.5, 139.5], [128.0, 361.5], [127.5, 360.0], [136.5, 382.5], [131.5, 378.5], [126.5, 370.0], [131.5, 378.5], [128.0, 361.5], [126.5, 370.0], [105.5, 343.5], [101.0, 324.5], [105.5, 343.5], [121.5, 347.5], [126.0, 353.0], [127.5, 360.0], [121.5, 347.5], [126.0, 353.0], [191.0, 64.5], [198.5, 72.0], [237.5, 83.5], [239.5, 90.0], [145.5, 49.0], [138.5, 49.0], [159.0, 57.0], [162.0, 57.5], [145.5, 49.0], [159.0, 57.0], [265.0, 229.5], [254.5, 220.0], [253.0, 216.5], [254.5, 220.0], [253.0, 216.5], [248.0, 208.5], [248.0, 208.5], [239.5, 204.5], [245.0, 173.5], [245.5, 178.0], [250.0, 158.0], [245.0, 173.5], [257.5, 139.5], [250.0, 158.0], [177.0, 396.5], [181.0, 395.5], [181.0, 395.5], [189.5, 400.5], [147.0, 377.0], [150.5, 377.0], [140.5, 381.5], [147.0, 377.0], [140.5, 381.5], [136.5, 382.5], [92.5, 313.5], [101.0, 324.5], [99.5, 290.0], [92.5, 313.5], [98.0, 271.0], [99.5, 290.0], [134.5, 47.5], [130.0, 49.0], [134.5, 47.5], [138.5, 49.0], [73.0, 222.5], [71.0, 214.5], [107.5, 246.0], [110.5, 248.0], [104.0, 266.5], [98.0, 271.0], [69.0, 209.5], [71.0, 214.5], [226.5, 87.0], [237.5, 83.5], [226.5, 87.0], [205.5, 94.5], [205.5, 94.5], [205.5, 77.5], [198.5, 72.0], [205.5, 77.5], [96.5, 244.5], [107.5, 246.0], [109.0, 265.5], [104.0, 266.5], [108.0, 263.5], [110.5, 248.0], [109.0, 265.5], [108.0, 263.5], [65.0, 205.5], [69.0, 209.5], [57.0, 198.0], [56.0, 201.5], [65.0, 205.5], [56.0, 201.5], [90.0, 240.0], [96.5, 244.5], [83.0, 224.0], [73.0, 222.5], [90.5, 231.5], [83.0, 224.0], [90.0, 240.0], [90.5, 231.5]]
x_list = [ p[0] for p in midpoints ]
y_list = [ p[1] for p in midpoints ]
complexmdpts = [ [p[0]+1j*p[1]] for p in midpoints ]
plt.scatter(x_list, y_list, s=50, marker="x", color='y')
coefs = np.fft.fftshift(scipy.fft.fft(complexmdpts))
n = len(coefs)
print("coeffs[{}]:\n{}".format(n, coefs[:5]))
#todo: sort coeffs?
# function in terms of t to trace out curve
def f(t):
ftx=0
fty=0
for i in range(-int(n/2), int(n/2)+1):
ftx+=(coefs[i]*cmath.exp(1j*2*math.pi*i*t/n)).real.tolist()[0]
fty+=(coefs[i]*cmath.exp(1j*2*math.pi*i*t/n)).imag.tolist()[0]
return [ftx/n, fty/n]
lines = [] # store computed lines segments to approximate function
pft = f(0) # compute first point
t_list = np.linspace(0, 2*math.pi, n) # compute list of dts to use when drawing
for t in t_list:
cft = f(t)
# print("f({}): {} \n".format(t, cft))
lines.append([cft,pft])
pft = cft
#draw f(t) approximation
lc = mc.LineCollection(lines)
fig, ax = pl.subplots()
ax.add_collection(lc)
ax.autoscale()
ax.margins(0.1)
这是我的输出：
我想概述的要点
我的函数逼近
我不确定快速傅里叶变换是否正确使用。从我读到的 fft 是我需要的，我移动它是因为 scipy fft returns 数组移动了，我认为我的其余代码是正确的，假设系数是正确的，这就是为什么我怀疑这些系数。
变换和系数之间是否有我遗漏的步骤？还是我的函数评估给出的系数不正确？还是我还漏掉了什么？
实际上有几个问题需要解决才能得到预期的大纲。让我们逐一讨论这些问题。
傅立叶系数计算
由于您正在计算二维数组 complexmdpts 的 FFT（每个元素都是 1 个复值的数组），scipy.fft 的默认行为是计算沿最后一个轴的 FFT .在这种情况下，这意味着您实际上正在计算 n 长度为 1 的 FFT，并且假设长度为 1 的 FFT 是恒等式，则整个计算 returns 原始数组。
一种解决方案是明确指定 axis=0：
scipy.fft(complexmdpts,axis=0)
您也可以将二维数组展平为一维数组，因为它基本上只有一维，但这会导致对基于此二维结构的其余代码进行更多更改。
系数平移
在尝试解释 FFT 系数时存在一个常见的混淆。较高指数中的系数对应于发生在奈奎斯特频率之上的环绕之后的负频率。
为了 可视化 负频率出现在其直观位置的图表上的这些系数，使用 fftshift 将较高索引中的系数与较低索引中的系数进行交换。这样负频率对应的系数出现在正频率对应的系数之前。
在您的特定情况下，f(t) 中的计算将负频率（只要 i 为负，cmath.exp 的参数也为负）与负索引处的系数相关联。由于 python 数组索引方便地使用从数组末尾倒数的元素作为负索引，因此在使用负索引进行索引时，您正确地使用了对应于负频率的系数。
所有这些都是说你不需要用 fftshift 交换系数，直接用：
获得系数
coefs = scipy.fft(complexmdpts,axis=0)
采样域
您已将 t 指定为 0 到 2*math.pi 的范围。考虑到 f(t) 的实现，域实际上在 0 到 n 之间（例如，当 i==1 时 cmath.exp 的参数从 0 到 1j*2*math.pi 为 t 从 0 到 n）。要使曲线跨越整个域，您应该相应地更新 t 的值：
t_list = np.linspace(0, n, n)
积分排序
最后，通过使用散点图绘制原始点系列，您隐藏了线图的样子：
试图获得一个平滑的函数来插入这个结果在一个图中，相似的线从轮廓的一侧跳到另一侧。如果你想捕捉这些跳跃，那么我想你已经完成了。但我怀疑您可能想要捕捉该区域的轮廓。在那种情况下，您必须对原始点进行排序，以便它们遵循大纲。
一种方法是迭代地将最近的点附加到先前选择的点。这将给出更接近轮廓的东西：
出于演示目的（即不声称这是最好的方法），您可以使用类似以下的方法进行排序：
def sort_points(points):
# pick a point
reference_point = points[0]
sorted = [reference_point]
remaining_points = range(1,len(points))
for i in range(1,len(points)):
# find the closest point to reference_point,
mindiff = np.sum(np.square(np.array(points[remaining_points[0]])-reference_point))
idx = 0
# loop over all the other remaining points
for j in range(1,len(remaining_points)):
diff = np.sum(np.square(np.array(points[remaining_points[j]])-reference_point))
if diff < mindiff:
mindiff = diff
idx = j
# found the closest: update the selected point, and add it to the list of sorted points
reference_point = points[remaining_points[idx]]
sorted.append(reference_point )
remaining_points = np.delete(remaining_points, idx)
return sorted
供参考，完整代码如下：
import matplotlib.pyplot as plt
import scipy
import cmath
import math
import numpy as np
from matplotlib import collections as mc
import pylab as pl
def sort_points(points):
# pick a point
reference_point = points[0]
sorted = [reference_point]
remaining_points = range(1,len(points))
for i in range(1,len(points)):
# find the closest point to reference_point,
mindiff = np.sum(np.square(np.array(points[remaining_points[0]])-reference_point))
idx = 0
# loop over all the other remaining points
for j in range(1,len(remaining_points)):
diff = np.sum(np.square(np.array(points[remaining_points[j]])-reference_point))
if diff < mindiff:
mindiff = diff
idx = j
# found the closest: update the selected point, and add it to the list of sorted points
reference_point = points[remaining_points[idx]]
sorted.append(reference_point )
remaining_points = np.delete(remaining_points, idx)
return sorted
midpoints = [[305.0, 244.5], [293.5, 237.0], [274.5, 367.5], [270.5, 373.5], [229.5, 391.0], [216.0, 396.0], [302.0, 269.0], [295.0, 271.0], [60.5, 146.5], [54.0, 153.0], [52.0, 167.0], [54.0, 153.0], [52.0, 167.0], [45.0, 178.0], [75.0, 76.5], [68.5, 98.5], [75.0, 76.5], [97.0, 58.5], [283.5, 357.5], [274.5, 367.5], [309.0, 255.0], [305.0, 244.5], [309.0, 255.0], [302.0, 269.0], [299.5, 291.5], [300.0, 297.0], [300.0, 297.0], [295.0, 309.5], [62.5, 105.0], [61.0, 118.5], [62.5, 105.0], [68.5, 98.5], [58.0, 139.0], [60.5, 146.5], [241.0, 111.5], [252.0, 124.5], [256.0, 132.5], [252.0, 124.5], [283.0, 356.0], [283.5, 357.5], [300.0, 290.5], [299.5, 291.5], [300.0, 290.5], [296.0, 280.5], [296.0, 280.5], [295.0, 271.0], [158.0, 387.0], [177.0, 396.5], [197.5, 402.0], [192.5, 403.0], [189.5, 400.5], [192.5, 403.0], [197.5, 402.0], [202.5, 401.0], [214.0, 395.5], [216.0, 396.0], [202.5, 401.0], [214.0, 395.5], [233.5, 375.0], [229.5, 391.0], [233.5, 375.0], [249.0, 372.5], [282.5, 340.0], [284.5, 328.0], [284.5, 328.0], [295.0, 309.5], [45.0, 178.0], [49.5, 189.0], [57.0, 198.0], [49.5, 189.0], [238.5, 108.5], [241.0, 111.5], [162.0, 57.5], [170.0, 59.0], [239.5, 204.5], [239.0, 200.0], [293.5, 237.0], [291.0, 227.5], [265.0, 229.5], [291.0, 227.5], [239.0, 189.5], [245.5, 178.0], [239.0, 189.5], [241.0, 193.5], [241.0, 193.5], [239.0, 200.0], [55.0, 119.0], [61.0, 118.5], [53.5, 134.0], [58.0, 139.0], [50.0, 129.0], [55.0, 119.0], [50.0, 129.0], [53.5, 134.0], [107.0, 46.0], [119.5, 50.5], [97.0, 54.0], [97.0, 58.5], [107.0, 46.0], [97.0, 54.0], [150.5, 377.0], [158.0, 387.0], [257.5, 367.5], [270.5, 373.5], [249.0, 372.5], [257.5, 367.5], [280.0, 349.5], [282.5, 340.0], [280.0, 349.5], [283.0, 356.0], [239.5, 90.0], [238.5, 98.0], [238.5, 108.5], [238.5, 98.0], [130.0, 49.0], [119.5, 50.5], [189.0, 65.0], [191.0, 64.5], [189.0, 65.0], [177.0, 62.5], [170.0, 59.0], [177.0, 62.5], [256.0, 132.5], [257.5, 139.5], [128.0, 361.5], [127.5, 360.0], [136.5, 382.5], [131.5, 378.5], [126.5, 370.0], [131.5, 378.5], [128.0, 361.5], [126.5, 370.0], [105.5, 343.5], [101.0, 324.5], [105.5, 343.5], [121.5, 347.5], [126.0, 353.0], [127.5, 360.0], [121.5, 347.5], [126.0, 353.0], [191.0, 64.5], [198.5, 72.0], [237.5, 83.5], [239.5, 90.0], [145.5, 49.0], [138.5, 49.0], [159.0, 57.0], [162.0, 57.5], [145.5, 49.0], [159.0, 57.0], [265.0, 229.5], [254.5, 220.0], [253.0, 216.5], [254.5, 220.0], [253.0, 216.5], [248.0, 208.5], [248.0, 208.5], [239.5, 204.5], [245.0, 173.5], [245.5, 178.0], [250.0, 158.0], [245.0, 173.5], [257.5, 139.5], [250.0, 158.0], [177.0, 396.5], [181.0, 395.5], [181.0, 395.5], [189.5, 400.5], [147.0, 377.0], [150.5, 377.0], [140.5, 381.5], [147.0, 377.0], [140.5, 381.5], [136.5, 382.5], [92.5, 313.5], [101.0, 324.5], [99.5, 290.0], [92.5, 313.5], [98.0, 271.0], [99.5, 290.0], [134.5, 47.5], [130.0, 49.0], [134.5, 47.5], [138.5, 49.0], [73.0, 222.5], [71.0, 214.5], [107.5, 246.0], [110.5, 248.0], [104.0, 266.5], [98.0, 271.0], [69.0, 209.5], [71.0, 214.5], [226.5, 87.0], [237.5, 83.5], [226.5, 87.0], [205.5, 94.5], [205.5, 94.5], [205.5, 77.5], [198.5, 72.0], [205.5, 77.5], [96.5, 244.5], [107.5, 246.0], [109.0, 265.5], [104.0, 266.5], [108.0, 263.5], [110.5, 248.0], [109.0, 265.5], [108.0, 263.5], [65.0, 205.5], [69.0, 209.5], [57.0, 198.0], [56.0, 201.5], [65.0, 205.5], [56.0, 201.5], [90.0, 240.0], [96.5, 244.5], [83.0, 224.0], [73.0, 222.5], [90.5, 231.5], [83.0, 224.0], [90.0, 240.0], [90.5, 231.5]]
midpoints = sort_points(midpoints)
x_list = [ p[0] for p in midpoints ]
y_list = [ p[1] for p in midpoints ]
complexmdpts = [ [p[0]+1j*p[1]] for p in midpoints ]
plt.scatter(x_list, y_list, s=50, marker="x", color='y')
coefs = scipy.fft(complexmdpts,axis=0)
n = len(coefs)
print("coeffs[{}]:\n{}".format(n, coefs[:5]))
# function in terms of t to trace out curve
m = n
def f(t):
ftx=0
fty=0
for i in range(-int(m/2), int(m/2)+1):
ftx+=(coefs[i]*cmath.exp(1j*2*math.pi*i*t/n)).real.tolist()[0]
fty+=(coefs[i]*cmath.exp(1j*2*math.pi*i*t/n)).imag.tolist()[0]
return [ftx/n, fty/n]
lines = [] # store computed lines segments to approximate function
pft = f(0) # compute first point
t_list = np.linspace(0, n, n) # compute list of dts to use when drawing
for t in t_list:
cft = f(t)
lines.append([cft,pft])
pft = cft
#draw f(t) approximation
lc = mc.LineCollection(lines)
fig, ax = pl.subplots()
ax.add_collection(lc)
ax.autoscale()
ax.margins(0.1)
我想计算一个函数（傅里叶级数），它通过一组给定的点。
与这里发生的事情类似 https://gofigure.impara.ai/ ，但我不希望对其进行动画处理。我只想要这个功能，这样我就可以自己画形状了。我已经阅读了很多描述它的数学资料和动画它的代码，但我正在为我的实现而苦苦挣扎。
我现在的代码如下[应该可以运行单独在python笔记本中]
import matplotlib.pyplot as plt
import scipy
import cmath
import math
import numpy as np
from matplotlib import collections as mc
import pylab as pl
midpoints = [[305.0, 244.5], [293.5, 237.0], [274.5, 367.5], [270.5, 373.5], [229.5, 391.0], [216.0, 396.0], [302.0, 269.0], [295.0, 271.0], [60.5, 146.5], [54.0, 153.0], [52.0, 167.0], [54.0, 153.0], [52.0, 167.0], [45.0, 178.0], [75.0, 76.5], [68.5, 98.5], [75.0, 76.5], [97.0, 58.5], [283.5, 357.5], [274.5, 367.5], [309.0, 255.0], [305.0, 244.5], [309.0, 255.0], [302.0, 269.0], [299.5, 291.5], [300.0, 297.0], [300.0, 297.0], [295.0, 309.5], [62.5, 105.0], [61.0, 118.5], [62.5, 105.0], [68.5, 98.5], [58.0, 139.0], [60.5, 146.5], [241.0, 111.5], [252.0, 124.5], [256.0, 132.5], [252.0, 124.5], [283.0, 356.0], [283.5, 357.5], [300.0, 290.5], [299.5, 291.5], [300.0, 290.5], [296.0, 280.5], [296.0, 280.5], [295.0, 271.0], [158.0, 387.0], [177.0, 396.5], [197.5, 402.0], [192.5, 403.0], [189.5, 400.5], [192.5, 403.0], [197.5, 402.0], [202.5, 401.0], [214.0, 395.5], [216.0, 396.0], [202.5, 401.0], [214.0, 395.5], [233.5, 375.0], [229.5, 391.0], [233.5, 375.0], [249.0, 372.5], [282.5, 340.0], [284.5, 328.0], [284.5, 328.0], [295.0, 309.5], [45.0, 178.0], [49.5, 189.0], [57.0, 198.0], [49.5, 189.0], [238.5, 108.5], [241.0, 111.5], [162.0, 57.5], [170.0, 59.0], [239.5, 204.5], [239.0, 200.0], [293.5, 237.0], [291.0, 227.5], [265.0, 229.5], [291.0, 227.5], [239.0, 189.5], [245.5, 178.0], [239.0, 189.5], [241.0, 193.5], [241.0, 193.5], [239.0, 200.0], [55.0, 119.0], [61.0, 118.5], [53.5, 134.0], [58.0, 139.0], [50.0, 129.0], [55.0, 119.0], [50.0, 129.0], [53.5, 134.0], [107.0, 46.0], [119.5, 50.5], [97.0, 54.0], [97.0, 58.5], [107.0, 46.0], [97.0, 54.0], [150.5, 377.0], [158.0, 387.0], [257.5, 367.5], [270.5, 373.5], [249.0, 372.5], [257.5, 367.5], [280.0, 349.5], [282.5, 340.0], [280.0, 349.5], [283.0, 356.0], [239.5, 90.0], [238.5, 98.0], [238.5, 108.5], [238.5, 98.0], [130.0, 49.0], [119.5, 50.5], [189.0, 65.0], [191.0, 64.5], [189.0, 65.0], [177.0, 62.5], [170.0, 59.0], [177.0, 62.5], [256.0, 132.5], [257.5, 139.5], [128.0, 361.5], [127.5, 360.0], [136.5, 382.5], [131.5, 378.5], [126.5, 370.0], [131.5, 378.5], [128.0, 361.5], [126.5, 370.0], [105.5, 343.5], [101.0, 324.5], [105.5, 343.5], [121.5, 347.5], [126.0, 353.0], [127.5, 360.0], [121.5, 347.5], [126.0, 353.0], [191.0, 64.5], [198.5, 72.0], [237.5, 83.5], [239.5, 90.0], [145.5, 49.0], [138.5, 49.0], [159.0, 57.0], [162.0, 57.5], [145.5, 49.0], [159.0, 57.0], [265.0, 229.5], [254.5, 220.0], [253.0, 216.5], [254.5, 220.0], [253.0, 216.5], [248.0, 208.5], [248.0, 208.5], [239.5, 204.5], [245.0, 173.5], [245.5, 178.0], [250.0, 158.0], [245.0, 173.5], [257.5, 139.5], [250.0, 158.0], [177.0, 396.5], [181.0, 395.5], [181.0, 395.5], [189.5, 400.5], [147.0, 377.0], [150.5, 377.0], [140.5, 381.5], [147.0, 377.0], [140.5, 381.5], [136.5, 382.5], [92.5, 313.5], [101.0, 324.5], [99.5, 290.0], [92.5, 313.5], [98.0, 271.0], [99.5, 290.0], [134.5, 47.5], [130.0, 49.0], [134.5, 47.5], [138.5, 49.0], [73.0, 222.5], [71.0, 214.5], [107.5, 246.0], [110.5, 248.0], [104.0, 266.5], [98.0, 271.0], [69.0, 209.5], [71.0, 214.5], [226.5, 87.0], [237.5, 83.5], [226.5, 87.0], [205.5, 94.5], [205.5, 94.5], [205.5, 77.5], [198.5, 72.0], [205.5, 77.5], [96.5, 244.5], [107.5, 246.0], [109.0, 265.5], [104.0, 266.5], [108.0, 263.5], [110.5, 248.0], [109.0, 265.5], [108.0, 263.5], [65.0, 205.5], [69.0, 209.5], [57.0, 198.0], [56.0, 201.5], [65.0, 205.5], [56.0, 201.5], [90.0, 240.0], [96.5, 244.5], [83.0, 224.0], [73.0, 222.5], [90.5, 231.5], [83.0, 224.0], [90.0, 240.0], [90.5, 231.5]]
x_list = [ p[0] for p in midpoints ]
y_list = [ p[1] for p in midpoints ]
complexmdpts = [ [p[0]+1j*p[1]] for p in midpoints ]
plt.scatter(x_list, y_list, s=50, marker="x", color='y')
coefs = np.fft.fftshift(scipy.fft.fft(complexmdpts))
n = len(coefs)
print("coeffs[{}]:\n{}".format(n, coefs[:5]))
#todo: sort coeffs?
# function in terms of t to trace out curve
def f(t):
ftx=0
fty=0
for i in range(-int(n/2), int(n/2)+1):
ftx+=(coefs[i]*cmath.exp(1j*2*math.pi*i*t/n)).real.tolist()[0]
fty+=(coefs[i]*cmath.exp(1j*2*math.pi*i*t/n)).imag.tolist()[0]
return [ftx/n, fty/n]
lines = [] # store computed lines segments to approximate function
pft = f(0) # compute first point
t_list = np.linspace(0, 2*math.pi, n) # compute list of dts to use when drawing
for t in t_list:
cft = f(t)
# print("f({}): {} \n".format(t, cft))
lines.append([cft,pft])
pft = cft
#draw f(t) approximation
lc = mc.LineCollection(lines)
fig, ax = pl.subplots()
ax.add_collection(lc)
ax.autoscale()
ax.margins(0.1)
这是我的输出：
我想概述的要点
我的函数逼近
我不确定快速傅里叶变换是否正确使用。从我读到的 fft 是我需要的，我移动它是因为 scipy fft returns 数组移动了，我认为我的其余代码是正确的，假设系数是正确的，这就是为什么我怀疑这些系数。
变换和系数之间是否有我遗漏的步骤？还是我的函数评估给出的系数不正确？还是我还漏掉了什么？
实际上有几个问题需要解决才能得到预期的大纲。让我们逐一讨论这些问题。
傅立叶系数计算
由于您正在计算二维数组 complexmdpts 的 FFT（每个元素都是 1 个复值的数组），scipy.fft 的默认行为是计算沿最后一个轴的 FFT .在这种情况下，这意味着您实际上正在计算 n 长度为 1 的 FFT，并且假设长度为 1 的 FFT 是恒等式，则整个计算 returns 原始数组。
一种解决方案是明确指定 axis=0：
scipy.fft(complexmdpts,axis=0)
您也可以将二维数组展平为一维数组，因为它基本上只有一维，但这会导致对基于此二维结构的其余代码进行更多更改。
系数平移
在尝试解释 FFT 系数时存在一个常见的混淆。较高指数中的系数对应于发生在奈奎斯特频率之上的环绕之后的负频率。
为了 可视化 负频率出现在其直观位置的图表上的这些系数，使用 fftshift 将较高索引中的系数与较低索引中的系数进行交换。这样负频率对应的系数出现在正频率对应的系数之前。
在您的特定情况下，f(t) 中的计算将负频率（只要 i 为负，cmath.exp 的参数也为负）与负索引处的系数相关联。由于 python 数组索引方便地使用从数组末尾倒数的元素作为负索引，因此在使用负索引进行索引时，您正确地使用了对应于负频率的系数。
所有这些都是说你不需要用 fftshift 交换系数，直接用：
coefs = scipy.fft(complexmdpts,axis=0)
采样域
您已将 t 指定为 0 到 2*math.pi 的范围。考虑到 f(t) 的实现，域实际上在 0 到 n 之间（例如，当 i==1 时 cmath.exp 的参数从 0 到 1j*2*math.pi 为 t 从 0 到 n）。要使曲线跨越整个域，您应该相应地更新 t 的值：
t_list = np.linspace(0, n, n)
积分排序
最后，通过使用散点图绘制原始点系列，您隐藏了线图的样子：
试图获得一个平滑的函数来插入这个结果在一个图中，相似的线从轮廓的一侧跳到另一侧。如果你想捕捉这些跳跃，那么我想你已经完成了。但我怀疑您可能想要捕捉该区域的轮廓。在那种情况下，您必须对原始点进行排序，以便它们遵循大纲。 一种方法是迭代地将最近的点附加到先前选择的点。这将给出更接近轮廓的东西：
出于演示目的（即不声称这是最好的方法），您可以使用类似以下的方法进行排序：
def sort_points(points):
# pick a point
reference_point = points[0]
sorted = [reference_point]
remaining_points = range(1,len(points))
for i in range(1,len(points)):
# find the closest point to reference_point,
mindiff = np.sum(np.square(np.array(points[remaining_points[0]])-reference_point))
idx = 0
# loop over all the other remaining points
for j in range(1,len(remaining_points)):
diff = np.sum(np.square(np.array(points[remaining_points[j]])-reference_point))
if diff < mindiff:
mindiff = diff
idx = j
# found the closest: update the selected point, and add it to the list of sorted points
reference_point = points[remaining_points[idx]]
sorted.append(reference_point )
remaining_points = np.delete(remaining_points, idx)
return sorted
供参考，完整代码如下：
import matplotlib.pyplot as plt
import scipy
import cmath
import math
import numpy as np
from matplotlib import collections as mc
import pylab as pl
def sort_points(points):
# pick a point
reference_point = points[0]
sorted = [reference_point]
remaining_points = range(1,len(points))
for i in range(1,len(points)):
# find the closest point to reference_point,
mindiff = np.sum(np.square(np.array(points[remaining_points[0]])-reference_point))
idx = 0
# loop over all the other remaining points
for j in range(1,len(remaining_points)):
diff = np.sum(np.square(np.array(points[remaining_points[j]])-reference_point))
if diff < mindiff:
mindiff = diff
idx = j
# found the closest: update the selected point, and add it to the list of sorted points
reference_point = points[remaining_points[idx]]
sorted.append(reference_point )
remaining_points = np.delete(remaining_points, idx)
return sorted
midpoints = [[305.0, 244.5], [293.5, 237.0], [274.5, 367.5], [270.5, 373.5], [229.5, 391.0], [216.0, 396.0], [302.0, 269.0], [295.0, 271.0], [60.5, 146.5], [54.0, 153.0], [52.0, 167.0], [54.0, 153.0], [52.0, 167.0], [45.0, 178.0], [75.0, 76.5], [68.5, 98.5], [75.0, 76.5], [97.0, 58.5], [283.5, 357.5], [274.5, 367.5], [309.0, 255.0], [305.0, 244.5], [309.0, 255.0], [302.0, 269.0], [299.5, 291.5], [300.0, 297.0], [300.0, 297.0], [295.0, 309.5], [62.5, 105.0], [61.0, 118.5], [62.5, 105.0], [68.5, 98.5], [58.0, 139.0], [60.5, 146.5], [241.0, 111.5], [252.0, 124.5], [256.0, 132.5], [252.0, 124.5], [283.0, 356.0], [283.5, 357.5], [300.0, 290.5], [299.5, 291.5], [300.0, 290.5], [296.0, 280.5], [296.0, 280.5], [295.0, 271.0], [158.0, 387.0], [177.0, 396.5], [197.5, 402.0], [192.5, 403.0], [189.5, 400.5], [192.5, 403.0], [197.5, 402.0], [202.5, 401.0], [214.0, 395.5], [216.0, 396.0], [202.5, 401.0], [214.0, 395.5], [233.5, 375.0], [229.5, 391.0], [233.5, 375.0], [249.0, 372.5], [282.5, 340.0], [284.5, 328.0], [284.5, 328.0], [295.0, 309.5], [45.0, 178.0], [49.5, 189.0], [57.0, 198.0], [49.5, 189.0], [238.5, 108.5], [241.0, 111.5], [162.0, 57.5], [170.0, 59.0], [239.5, 204.5], [239.0, 200.0], [293.5, 237.0], [291.0, 227.5], [265.0, 229.5], [291.0, 227.5], [239.0, 189.5], [245.5, 178.0], [239.0, 189.5], [241.0, 193.5], [241.0, 193.5], [239.0, 200.0], [55.0, 119.0], [61.0, 118.5], [53.5, 134.0], [58.0, 139.0], [50.0, 129.0], [55.0, 119.0], [50.0, 129.0], [53.5, 134.0], [107.0, 46.0], [119.5, 50.5], [97.0, 54.0], [97.0, 58.5], [107.0, 46.0], [97.0, 54.0], [150.5, 377.0], [158.0, 387.0], [257.5, 367.5], [270.5, 373.5], [249.0, 372.5], [257.5, 367.5], [280.0, 349.5], [282.5, 340.0], [280.0, 349.5], [283.0, 356.0], [239.5, 90.0], [238.5, 98.0], [238.5, 108.5], [238.5, 98.0], [130.0, 49.0], [119.5, 50.5], [189.0, 65.0], [191.0, 64.5], [189.0, 65.0], [177.0, 62.5], [170.0, 59.0], [177.0, 62.5], [256.0, 132.5], [257.5, 139.5], [128.0, 361.5], [127.5, 360.0], [136.5, 382.5], [131.5, 378.5], [126.5, 370.0], [131.5, 378.5], [128.0, 361.5], [126.5, 370.0], [105.5, 343.5], [101.0, 324.5], [105.5, 343.5], [121.5, 347.5], [126.0, 353.0], [127.5, 360.0], [121.5, 347.5], [126.0, 353.0], [191.0, 64.5], [198.5, 72.0], [237.5, 83.5], [239.5, 90.0], [145.5, 49.0], [138.5, 49.0], [159.0, 57.0], [162.0, 57.5], [145.5, 49.0], [159.0, 57.0], [265.0, 229.5], [254.5, 220.0], [253.0, 216.5], [254.5, 220.0], [253.0, 216.5], [248.0, 208.5], [248.0, 208.5], [239.5, 204.5], [245.0, 173.5], [245.5, 178.0], [250.0, 158.0], [245.0, 173.5], [257.5, 139.5], [250.0, 158.0], [177.0, 396.5], [181.0, 395.5], [181.0, 395.5], [189.5, 400.5], [147.0, 377.0], [150.5, 377.0], [140.5, 381.5], [147.0, 377.0], [140.5, 381.5], [136.5, 382.5], [92.5, 313.5], [101.0, 324.5], [99.5, 290.0], [92.5, 313.5], [98.0, 271.0], [99.5, 290.0], [134.5, 47.5], [130.0, 49.0], [134.5, 47.5], [138.5, 49.0], [73.0, 222.5], [71.0, 214.5], [107.5, 246.0], [110.5, 248.0], [104.0, 266.5], [98.0, 271.0], [69.0, 209.5], [71.0, 214.5], [226.5, 87.0], [237.5, 83.5], [226.5, 87.0], [205.5, 94.5], [205.5, 94.5], [205.5, 77.5], [198.5, 72.0], [205.5, 77.5], [96.5, 244.5], [107.5, 246.0], [109.0, 265.5], [104.0, 266.5], [108.0, 263.5], [110.5, 248.0], [109.0, 265.5], [108.0, 263.5], [65.0, 205.5], [69.0, 209.5], [57.0, 198.0], [56.0, 201.5], [65.0, 205.5], [56.0, 201.5], [90.0, 240.0], [96.5, 244.5], [83.0, 224.0], [73.0, 222.5], [90.5, 231.5], [83.0, 224.0], [90.0, 240.0], [90.5, 231.5]]
midpoints = sort_points(midpoints)
x_list = [ p[0] for p in midpoints ]
y_list = [ p[1] for p in midpoints ]
complexmdpts = [ [p[0]+1j*p[1]] for p in midpoints ]
plt.scatter(x_list, y_list, s=50, marker="x", color='y')
coefs = scipy.fft(complexmdpts,axis=0)
n = len(coefs)
print("coeffs[{}]:\n{}".format(n, coefs[:5]))
# function in terms of t to trace out curve
m = n
def f(t):
ftx=0
fty=0
for i in range(-int(m/2), int(m/2)+1):
ftx+=(coefs[i]*cmath.exp(1j*2*math.pi*i*t/n)).real.tolist()[0]
fty+=(coefs[i]*cmath.exp(1j*2*math.pi*i*t/n)).imag.tolist()[0]
return [ftx/n, fty/n]
lines = [] # store computed lines segments to approximate function
pft = f(0) # compute first point
t_list = np.linspace(0, n, n) # compute list of dts to use when drawing
for t in t_list:
cft = f(t)
lines.append([cft,pft])
pft = cft
#draw f(t) approximation
lc = mc.LineCollection(lines)
fig, ax = pl.subplots()
ax.add_collection(lc)
ax.autoscale()
ax.margins(0.1)