IDL代码转 python问题记录

张霄柏

于 2023-08-29 16:51:31 发布

阅读量686

点赞数 17

文章标签： python scipy numpy

本文链接：https://blog.csdn.net/qq_40499108/article/details/132565130

版权

IDL 转 python 问题记录

fft(傅里叶变换)

在做某IDL代码转换为python 时发现了以下问题：
1、IDL中的fft方法结果与python中结果并不相同，不管是scipy.fft还是np.fft.fft 都存在差异
主要是结果精度存在差异

x = [4.9000001, 5.0999999, 4.7000003, 3.7000000, 3.8000000, 3.9000001, 3.7000000, 4.4000001,
     4.2000003, 3.8000000, 4.0000000, 2.8000000, 3.9000001, 4.9000001, 5.0999999, 4.7000003,
     3.7000000, 3.8000000, 3.9000001, 3.7000000, 4.4000001, 4.2000003, 3.8000000, 4.0000000,
     2.8000000, 3.9000001, 4.9000001, 5.0999999, 4.7000003, 3.7000000, 3.8000000, 3.9000001,
     3.7000000, 4.4000001, 4.2000003, 3.8000000, 4.0000000, 2.8000000, 3.9000001, 4.9000001,
     5.0999999, 4.7000003, 3.7000000, 3.8000000, 3.9000001, 3.7000000, 4.4000001, 4.2000003,
     3.8000000, 4.0000000, 2.8000000, 3.9000001, 4.9000001, 5.0999999, 4.7000003, 3.7000000,
     3.8000000, 3.9000001, 3.7000000, 4.4000001, 4.2000003, 3.8000000, 4.0000000, 2.8000000,
     3.9000001, 4.9000001, 5.0999999, 4.7000003, 3.7000000, 3.8000000, 3.9000001, 3.7000000,
     4.4000001, 4.2000003, 3.8000000, 4.0000000, 2.8000000, 3.9000001, 4.9000001, 5.0999999,
     4.7000003, 3.7000000, 3.8000000, 3.9000001, 3.7000000, 4.4000001, 4.2000003, 3.8000000,
     4.0000000, 2.8000000, 3.9000001, 4.9000001, 5.0999999, 4.7000003, 3.7000000, 3.8000000,
     3.9000001, 3.7000000, 4.4000001, 4.2000003, 3.8000000, 4.0000000, 2.8000000, 3.9000001,
     4.9000001, 5.0999999, 4.7000003, 3.7000000, 3.8000000, 3.9000001, 3.7000000, 4.4000001,
     4.2000003, 3.8000000, 4.0000000, 2.8000000, 3.9000001]

IDL 和 python numpy 分别对以上x进行fft

(       4.0692306,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)( -1.0188829e-009,  2.0377657e-009)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(  3.0566487e-009,      0.00000000)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)(      0.11039685,     -0.11208207)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
( -2.0377657e-009,  2.0377657e-009)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(  1.0188829e-009,  1.0188829e-009)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.17974262,     -0.25503349)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)( -4.0755315e-009, -8.1510629e-009)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
(      0.00000000, -1.0188829e-009)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(     0.095397085,    -0.093490288)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)( -5.0944142e-009,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
(     0.043500885,    -0.049952786)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(  2.0377657e-009,      0.00000000)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)(     0.027866472,      0.14249153)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,  4.0755315e-009)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(    -0.041519184,     0.029256972)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
(  8.1510629e-009,  8.1510629e-009)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(    -0.041519184,    -0.029256970)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000, -2.0377657e-009)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)(  4.0755315e-009,  4.0755315e-009)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
(     0.027866475,     -0.14249153)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,  8.8237845e-010)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)(     0.043500885,     0.049952783)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
( -3.0566487e-009,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(     0.095397085,     0.093490295)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)( -2.0377657e-009,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.17974260,      0.25503352)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(  2.0377657e-009, -1.0188829e-009)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)(      0.00000000, -2.0377657e-009)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)
(      0.11039685,      0.11208207)(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(  3.0566487e-009,      0.00000000)
(      0.00000000,      0.00000000)(      0.00000000,      0.00000000)(      0.00000000, -8.8237845e-010)(      0.00000000,      0.00000000)
(      0.00000000,      0.00000000)

#scipy.fft结果如下：
[ 4.06923084e+00-0.00000000e+00j  1.53638023e-16-3.92973826e-17j
  8.08033680e-17-1.81927335e-17j  4.93432455e-17-4.17519770e-17j
  6.57404929e-17-6.28012122e-17j  1.89157497e-17-2.36632406e-17j
  5.31388798e-17-6.07301484e-17j  3.16782483e-17-6.32308867e-17j
  4.35835329e-18-3.81146049e-17j  1.10396863e-01-1.12082056e-01j
 -3.13103043e-17-2.64775867e-17j  1.95684825e-17-3.18205311e-17j
 -3.03650742e-17+4.55476113e-17j -6.01708749e-17+3.11564011e-17j
  2.98116109e-17+5.39934398e-17j  1.89781714e-18+3.79563427e-18j
 -5.33681292e-18+1.59327504e-17j  4.35979066e-17-3.52710449e-17j
  1.79742572e-01-2.55033542e-01j  2.37276415e-17+4.85720902e-18j
  2.71165942e-17-5.25152192e-18j -2.83309999e-17-1.51825371e-17j
 -1.38105390e-17-4.41470666e-17j  1.11557712e-17-2.31618207e-17j
  3.79563427e-17-2.84672570e-17j  2.78165467e-17-1.60769532e-17j
  1.94464013e-16-5.85548803e-17j  9.53970774e-02-9.34903004e-02j
  2.45818818e-17-4.20305559e-17j -3.27494439e-17-3.53013764e-17j
  2.52798885e-17-1.31484628e-17j  6.57387881e-19+2.37851367e-17j
  5.92149728e-18+1.47884340e-17j  9.48908568e-18+0.00000000e+00j
  3.38403945e-17-1.49438944e-17j  4.49022424e-18-1.95693911e-17j
  4.35008342e-02-4.99528109e-02j -2.50236702e-18+9.71781184e-18j
  2.52673699e-17+3.71730724e-17j -4.85841187e-16+2.10375405e-16j
  6.02195660e-17+2.52293924e-17j  1.56409473e-17+2.07937310e-17j
  1.89781714e-17+1.89781714e-17j  3.78415115e-18-1.37026930e-17j
 -1.85954358e-17-2.66755659e-18j  2.78664449e-02+1.42491514e-01j
 -3.33684426e-18+2.73640327e-17j -1.64026440e-17+1.63521304e-17j
 -4.04624256e-17-1.31484628e-17j -9.16495971e-18-4.86027966e-18j
  1.84776012e-17-1.02100248e-17j -3.79563427e-17+6.07301484e-17j
 -1.34293138e-16-8.97112814e-17j -6.38709495e-18+8.36527950e-17j
 -4.15191607e-02+2.92569972e-02j  4.77778421e-18+1.61100016e-17j
  3.85340493e-17+4.80584924e-17j -2.03407428e-18+1.51825371e-17j
  3.14135988e-18-3.82950490e-18j  3.14135988e-18+3.82950490e-18j
 -2.03407428e-18-1.51825371e-17j  3.85340493e-17-4.80584924e-17j
  4.77778421e-18-1.61100016e-17j -4.15191607e-02-2.92569972e-02j
 -6.38709495e-18-8.36527950e-17j -1.34293138e-16+8.97112814e-17j
 -3.79563427e-17-6.07301484e-17j  1.84776012e-17+1.02100248e-17j
 -9.16495971e-18+4.86027966e-18j -4.04624256e-17+1.31484628e-17j
 -1.64026440e-17-1.63521304e-17j -3.33684426e-18-2.73640327e-17j
  2.78664449e-02-1.42491514e-01j -1.85954358e-17+2.66755659e-18j
  3.78415115e-18+1.37026930e-17j  1.89781714e-17-1.89781714e-17j
  1.56409473e-17-2.07937310e-17j  6.02195660e-17-2.52293924e-17j
 -4.85841187e-16-2.10375405e-16j  2.52673699e-17-3.71730724e-17j
 -2.50236702e-18-9.71781184e-18j  4.35008342e-02+4.99528109e-02j
  4.49022424e-18+1.95693911e-17j  3.38403945e-17+1.49438944e-17j
  9.48908568e-18-0.00000000e+00j  5.92149728e-18-1.47884340e-17j
  6.57387881e-19-2.37851367e-17j  2.52798885e-17+1.31484628e-17j
 -3.27494439e-17+3.53013764e-17j  2.45818818e-17+4.20305559e-17j
  9.53970774e-02+9.34903004e-02j  1.94464013e-16+5.85548803e-17j
  2.78165467e-17+1.60769532e-17j  3.79563427e-17+2.84672570e-17j
  1.11557712e-17+2.31618207e-17j -1.38105390e-17+4.41470666e-17j
 -2.83309999e-17+1.51825371e-17j  2.71165942e-17+5.25152192e-18j
  2.37276415e-17-4.85720902e-18j  1.79742572e-01+2.55033542e-01j
  4.35979066e-17+3.52710449e-17j -5.33681292e-18-1.59327504e-17j
  1.89781714e-18-3.79563427e-18j  2.98116109e-17-5.39934398e-17j
 -6.01708749e-17-3.11564011e-17j -3.03650742e-17-4.55476113e-17j
  1.95684825e-17+3.18205311e-17j -3.13103043e-17+2.64775867e-17j
  1.10396863e-01+1.12082056e-01j  4.35835329e-18+3.81146049e-17j
  3.16782483e-17+6.32308867e-17j  5.31388798e-17+6.07301484e-17j
  1.89157497e-17+2.36632406e-17j  6.57404929e-17+6.28012122e-17j
  4.93432455e-17+4.17519770e-17j  8.08033680e-17+1.81927335e-17j
  1.53638023e-16+3.92973826e-17j]


  #np.fft.fft结果如下
[ 4.06923084+0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.11039686-0.11208206j  0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.17974257-0.25503354j  0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.09539708-0.0934903j   0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.04350083-0.04995281j  0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.02786644+0.14249151j  0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
 -0.04151916+0.029257j    0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
 -0.04151916-0.029257j    0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.02786644-0.14249151j  0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.04350083+0.04995281j  0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.09539708+0.0934903j   0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.17974257+0.25503354j  0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.11039686+0.11208206j  0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j
  0.        +0.j          0.        +0.j          0.        +0.j  ]

通过多次试验可知，IDL中fft会自动进行归一化，也就是结果值/采样数据长度，所以python中的fft方法必须使用np.fft.fft（x,norm=‘forward’）,设置关键字 norm=‘forward’ 可以保证结果大部分一致，但是IDL中计算结果保留精度明显高于python，还是有部分结果不能完全一致，也尝试去设置数据为双精度浮点型，但是不能成功，这个问题未能解决，有知道的大佬可以指导一下。

滑动平均

IDL中smooth和python中多种数据平滑方式都不一致，故在此记录一下。

IDL

IDL 中，可以使用 SMOOTH 函数对一组数据进行平滑处理。SMOOTH 函数的基本语法如下：

css
Copy
result = SMOOTH(data, width [, /EDGE_TRIM])
其中，data 表示需要平滑处理的数据，可以是一维或多维数组；width 表示平滑窗口的宽度，必须是正整数；/EDGE_TRIM 表示是否对边缘进行修整，如果指定了该选项，则会在计算平均值时排除边缘数据。

以下是使用 SMOOTH 函数对一组数据进行平滑处理的示例代码：


; 创建一组数据
data = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10]

; 计算平滑后的数据
smoothed_data = SMOOTH(data, 3)

; 输出结果
PRINT, smoothed_data
输出结果为：


   1.50000   2.00000   3.00000   4.00000   5.00000   6.00000   7.00000   8.00000   9.00000   9.50000

在上述代码中，我们首先创建了一维数组 data，表示需要平滑处理的数据。然后，使用 SMOOTH 函数计算了 data 的平滑值，并将结果保存到变量 smoothed_data 中。最后，使用 IDL 的 PRINT 命令输出了 smoothed_data 的值。

需要注意的是，在使用 SMOOTH 函数计算平滑值时，如果指定了 /EDGE_TRIM 选项，则会在计算平均值时排除边缘数据。例如，以下代码将 SMOOTH 函数的 /EDGE_TRIM 选项设置为 1，表示对边缘进行修整：


; 创建一组数据
data = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10]

; 计算平滑后的数据（对边缘进行修整）
smoothed_data = SMOOTH(data, 3, /EDGE_TRIM)

; 输出结果
PRINT, smoothed_data
输出结果为：


   2.00000   3.00000   4.00000   5.00000   6.00000   7.00000   8.00000   9.00000

在上述代码中，我们将 SMOOTH 函数的 /EDGE_TRIM 选项设置为 1，表示对边缘进行修整。由于窗口大小为 $3$ ，因此在计算第一个平滑值时，只使用了前两个数的平均值。在计算最后一个平滑值时，只使用了后两个数的平均值。

python

在 NumPy 中，可以使用 numpy.convolve 函数来实现一维数组的滑动平均。numpy.convolve 函数的基本语法如下：

result = numpy.convolve(data, window, mode=‘valid’) / sum(window)
其中，data 表示需要平滑处理的一维数组，window 表示平滑窗口，必须是一维数组，mode 表示卷积模式，可以是 full、same 或 valid 中的一个。在计算滑动平均时，通常使用 mode=‘valid’。

以下是使用 numpy.convolve 函数对一组数据进行滑动平均的示例代码：


import numpy as np

# 创建一组数据
data = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10]

# 创建平滑窗口
window = np.ones(3) / 3

# 计算滑动平均
smoothed_data = np.convolve(data, window, mode='valid')

# 输出结果
print(smoothed_data)
输出结果为：


[2. 3. 4. 5. 6. 7. 8. 9.]

在上述代码中，我们首先创建了一维数组 data，表示需要平滑处理的数据。然后，使用 np.ones 函数创建了一个长度为 $3$ 的一维数组 window，表示平滑窗口。接着，使用 np.convolve 函数计算了 data 的滑动平均，并将结果保存到变量 smoothed_data 中。最后，使用 Python 的 print 函数输出了 smoothed_data 的值。

需要注意的是，在使用 np.convolve 函数计算滑动平均时，由于窗口大小可能会导致边缘数据的缺失，因此通常需要在计算结果时进行修整。在上述代码中，我们使用了 mode=‘valid’ 参数来指定卷积模式为“有效卷积”，这样就可以避免边缘数据的缺失。此外，由于卷积后的结果会比原始数据短 $n - 1$ 个元素（其中 $n$ 是窗口大小），因此在计算结果时通常需要除以窗口中所有元素之和，以保证结果的总和不变。在上述代码中，我们通过除以 sum(window) 来实现了这一点。

除了 numpy.convolve 函数之外，Python 中还有其他一些常用的平滑函数，例如：

pandas.DataFrame.rolling：用于对 Pandas 数据帧进行滑动窗口计算。具体而言，可以使用 rolling 方法创建一个滑动窗口对象，然后调用其 mean、median、std 等方法计算滑动平均、中位数、标准差等指标。以下是使用 rolling 方法计算一组数据的滑动平均的示例代码：


import pandas as pd

# 创建一组数据
data = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10]

# 将数据转换为 Pandas 数据帧
df = pd.DataFrame(data)

# 计算滑动平均
smoothed_data = df.rolling(window=3).mean()

# 输出结果
print(smoothed_data)
输出结果为：


       0
0    NaN
1    NaN
2    2.0
3    3.0
4    4.0
5    5.0
6    6.0
7    7.0
8    8.0
9    9.0

scipy.signal.savgol_filter：用于进行 Savitzky-Golay 滤波，可以实现多项式拟合和平滑处理。具体而言，可以使用 savgol_filter 函数对一维或多维数组进行平滑处理，指定窗口长度、多项式次数等参数。以下是使用 savgol_filter 函数对一组数据进行平滑处理的示例代码：


from scipy.signal import savgol_filter

# 创建一组数据
data = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10]

# 计算平滑后的数据
smoothed_data = savgol_filter(data, window_length=3, polyorder=1)

# 输出结果
print(smoothed_data)
输出结果为：


[1.         2.         3.         ... [omitted]

需要注意的是，不同的平滑函数适用于不同的数据类型和处理场景，具体使用时需要根据实际情况选择合适的函数。

张霄柏

关注

17
点赞
踩
6

收藏

觉得还不错? 一键收藏
0
评论
IDL代码转 python问题记录

在做某IDL代码转换为python 时发现了以下问题：1、IDL中的fft方法结果与python中结果并不相同，不管是scipy.fft还是np.fft.fft 都存在差异主要是结果精度存在差异IDL 和 python numpy 分别对以上x进行fft。
复制链接

扫一扫