2022-8-30到9-02学习笔记

1、查找python函数源码方法

函数名.file,返回函数在包中的位置，然后可以一级一级的追踪；
例子：

import matplotlib.pyplot as plt
print(plt.__file__)返回：'C:\\Users\\Alien\\AppData\\Local\\Programs\\Python\\Python37\\lib\\site-packages\\matplotlib\\pyplot.py'

2、判断某些元素是否在某些元素中

借用numpy的cbook._check_in_list（），但是需要研究怎么用~

补充_check_in_list的用法

# TODO 补充cbook._check_in_list怎么用
import matplotlib.cbook as cbook

3、判断是否为函数对象

callable()函数可以判断输入对象是否是函数
例子：

callable(123)
返回 Falsedef fun_help():return 1
callable(fun_help)
返回 True

4、plt.psd()的计算方法

（1）将数据按照NFFT的长度分块
（2）对每块数据进行滤波（窗函数默认hanning窗，窗系数大小为NFFT长度）
（3）每块数据*滤波器系数？？？？？？？？？
（3）对每块数据进行FFT处理
（4）用共轭的方法求功率
（5）功率/FS求出功率谱密度
（6）功率谱密度/窗函数功率
（7）用np.roll对结果进行翻转
（9）某一个频点的功率：用所有blk在该频点的平均值来替代

代码计算实例：

# 1、按照NFFT的大小把数据拆分成块,大小不够删掉尾部凑成整数倍
Nfft_len = 256  # psd的默认值是256
fs = 900
data_len = data_100.shape[0]
blk_num = data_len//Nfft_len
data_reshape = data_100[:blk_num*Nfft_len].reshape(blk_num,Nfft_len).T
# 2、对每块数据进行滤波（窗函数默认hanning窗，窗系数大小为NFFT长度）
window_coef = np.hanning(Nfft_len).reshape(-1,1)
# 3、每块数据*滤波器系数（因为滤波器长度和数据长度一样长，所以就直接相乘，相当于卷积）（滤波的目的是为了抗混叠，因为把数据直接拆成块，相当于加了矩形窗）
data_flt = data_reshape * window_coef
# data_flt = data_reshape * 1
# 4、对每块数据进行FFT处理
data_fft = np.fft.fft(data_flt,n=Nfft_len,axis=0)
# 5、用共轭的方法求功率
data_pow = abs(data_fft * np.conj(data_fft))   # 结果是功率，取abs是为了抛掉值为0的虚部
# 6、功率/FS求出功率谱密度
data_powfs = data_pow / fs
# 7、功率谱密度/窗函数功率
data_one_win = data_powfs / sum(abs(window_coef)**2)
# data_one_win = data_powfs / 1
# 8、用np.roll对结果进行翻转
data_fftshift = np.roll(data_one_win,-Nfft_len//2,axis=0)  # 垂直方向翻转：因为一个blk其实是一个频点周期
# 9、某一个频点的功率：用所有blk在该频点的平均值来替代
data_fft_psd = np.mean(data_fftshift,axis=1)
# 10、画图check
# fs_point = np.linspace(-fs/2,fs/2,Nfft_len)  # 错的
fs_point_temp = np.fft.fftfreq(Nfft_len, 1/fs)
fs_point = np.roll(fs_point_temp,-Nfft_len//2,axis=0)
plt.plot(fs_point,10*np.log10(data_fft_psd),'*-')  # 注意这里用的是real画图,是因为前面计算的结果是功率，所以这块虚部是0
plt.psd(data_100,Fs=fs)# 总结：用代码简单归纳fft和psd的关系
data_fft0 = data_fft[:,0]
data_pow = abs(data_fft0)**2
data_pow_fs = data_pow/fs
data_pow_coef = data_pow_fs/sum(abs(window_coef)**2)
data_shift = np.roll(data_pow_coef,128)
plt.plot(fs_point, 10*np.log10(data_shift),'*-')# 最后：psd的含义是功率谱，而fft的结果是每个频点幅值+相位，所以需要将fft的结果转换成功率然后/fs，
# 同时为了表示方便，将每个频点的功率转换成dB表示的形式

5、用noise产生信号有直流原因

根本原因是在生成noise信号时，用的是np.random.rand函数，它产生的数据是[0-1]均匀分布，所以数据的均值就是0.5，所以在0频附近会有直流，把生成的noise信号减去均值后，直流消失

# 生成底噪信号
noise_data = np.random.rand(8192)+np.random.rand(8192)*1j  # 生成的是[0-1]均匀分布,所以画出的频谱在0位置有尖
p1,f1 = plt.psd(noise_data,NFFT=1024,Fs=983.04)

data_mean = np.mean(noise_data)
data_diff_mean = noise_data - data_mean
plt.psd(data_diff_mean,NFFT=1024,Fs=983.04)

6、为什么原始带直流的noise信号用psd和scipy.singal.welch画出的结果不一样

根据原因就是两个函数对直流的默认处理不一样，体现在detrend的默认参数值，psd的默认是不处理，welch的默认是去直流

plt.subplots(3,1)
plt.subplot(311)
plt.psd(noise_data,NFFT=1024,Fs=983.04)
plt.subplot(312)
plt.psd(noise_data,NFFT=1024,Fs=983.04,detrend='constant')
plt.subplot(313)
f2,p2 = sig.welch(noise_data,983.40,window='hann',nperseg=1024,nfft=1024,return_onesided='False')  # 返回的是功率值
f2 = scipy.fft.fftshift(f2)
p2 = scipy.fft.fftshift(p2)
plt.plot(f2,10*np.log10(abs(p2)))  # 信号频响

fft与psd的关系【傅里叶变化求功率谱】相关推荐

傅里叶变化，短时傅里叶分析，小波变换
作者:咚懂咚懂咚链接:https://www.zhihu.com/question/22864189/answer/40772083 来源:知乎著作权归作者所有,转载请联系作者获得授权. 从傅里叶 ...
傅里叶变化（一）—— 复数
[参考资料] 1.万门大学:傅立叶变换.拉普拉斯变换与小波变换 [傅里叶变换系列博客] 1.傅里叶变化(一)-- 复数 2.傅里叶变换(二)-- 卷积复数的定义: z=a+ib(a,b∈R)z=a+ ...
傅里叶变化的本质：复数的实部和虚部的对应关系
之前做计算光学成像,需要用到图像的相位信息.但是设计到傅里叶变化的实部和虚部的问题的时候,发现教科书上一般来讲,只会介绍一句: 如果f(x,y)是实函数,则它的傅里叶变化就是关于原点共轭对称的: F( ...
傅里叶变化与卷积和互相关操作的转换
已知有二维信号f(x,y),g(x,y),其对应的傅里叶变化结果为F(x,y),G(x,y). 本篇文章不对公式进行推导,只是运用.有兴趣的话可以自行百度. 1.卷积与傅里叶的转换 FT{f(x,y) ...
图像处理：二维傅里叶变化的平移性_matlab实现
傅里叶变化的平移性: matlab代码验证过程实现: %% 研究傅里叶变化的平移特性 %空间域乘以exp ,频率域移动clc;clear I = imread('rice.jpg'); I = ...
matlab 范德蒙德矩阵,Matlab中fft与fwelch有什么区别？如何用fft求功率谱？
讲这个话题,就要先搞清楚频谱.功率谱的概念,可参考我的另一篇文章做信号处理的朋友应该都会fft比较熟悉,就是求傅里叶变换.我在这里也不再去讲这个函数了,但需要注意的一点:实信号的频谱关于0频对称,是 ...
fft matlab 区别,Matlab中fft与fwelch有什么区别？如何用fft求功率谱？
讲这个话题,就要先搞清楚频谱.功率谱的概念,可参考我的另一篇文章做信号处理的朋友应该都会fft比较熟悉,就是求傅里叶变换.我在这里也不再去讲这个函数了,但需要注意的一点:实信号的频谱关于0频对称,是 ...
matlab 傅里叶平移,图像处理：二维傅里叶变化的平移性_matlab实现
傅里叶变化的平移性: matlab代码验证过程实现: %% 研究傅里叶变化的平移特性 %空间域乘以exp ,频率域移动 clc;clear I = imread('rice.jpg'); I = ...
049万能图像处理小助手1.1_傅里叶变化_椒盐噪声_直方图均衡等图片批量处理
视频演示和demo仓库找049期银色子弹zg的个人空间-银色子弹zg个人主页-哔哩哔哩视频直接上效果图 049万能图像处理小助手1.1_傅里叶变化_椒盐噪声_直方图均衡等图片批量处理代码界面一 ...

fft与psd的关系【傅里叶变化求功率谱】

2022-8-30到9-02学习笔记

1、查找python函数源码方法

2、判断某些元素是否在某些元素中

3、判断是否为函数对象

4、plt.psd()的计算方法

5、用noise产生信号有直流原因

6、为什么原始带直流的noise信号用psd和scipy.singal.welch画出的结果不一样

fft与psd的关系【傅里叶变化求功率谱】相关推荐

最新文章

热门文章