Python量化噪声卷积信号和傅里叶时频分析

2024-04-10 23:12

本文主要是介绍Python量化噪声卷积信号和傅里叶时频分析,希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!

🎯要点

  1. 信号生成、绘制和推理图形结果:

    1. 连续时域数字信号:🎯正弦信号:🖊通用表达式、🖊不同相位式、🖊三相式。🎯指数信号:📌实数指数信号情景:如:油价指数增长函数、放射性衰减指数衰减函数,📌复数指数信号情景:傅里叶变换、正弦波形和余弦波形、线性时不变系统的特征函数,🖊实数指数、🖊PN结二极管的正向特性、🖊放射性衰减、🖊复数指数函数、🖊从复数指数提取正弦和余弦函数、🖊正弦信号与增长和衰减的实指数信号相乘。🎯非平稳信号、平稳信号、线性调频信号、方波、三角波、锯齿波、sinc 、脉冲信号和高斯信号。

    2. 离散时域信号:🎯脉冲函数列、🎯梳状函数的频域、🎯时域混叠、🎯频域混叠、🎯听觉混叠、🎯中平量化器传输特性、🎯中平量化器量化输入锯齿信号、🎯中平量化器的传输特性、🎯中升量化器量化输入锯齿信号、🎯语音信号量化、🎯均匀中平量化图像、🎯 μ \mu μ​-法则压扩编码输入信号、🎯输入信号和量化信号间的误差、🎯量化误差概率密度函数、🎯正弦信号零阶保持插值、🎯一阶保持插值、🎯正弦信号的理想或 Sinc 插值、🎯零阶保持与 Sinc 插值的比较。🎯单位样本序列生成方式:🖊使用逻辑运算、🖊使用Python库中内置函数、🖊单元步序列、🖊单元斜坡信号、🖊实数和复数指数信号,🎯数学运算:🖊幅度缩放、🖊幅度漂移、🖊信号乘积、🖊信号叠加、🖊上/下取样、🖊时移运算、:🖊使用Python内置/无内置函数的时间反转、🖊单位样本卷积、🖊信号与移位单位样本信号的卷积、🖊卷积特征演示:交换、结合和分布、🖊低通滤波系数的方波卷积、🖊高通滤波系数的方波卷积、🎯信号相似性:🖊正弦波和余弦波的自相关和互相关、🖊正弦波与其自身和噪声信号的自相关、🖊使用自相关的延迟估计。

    3. 离散时域系统:🎯解差分方程、🎯确定状态空间转换关系、🎯计算系统的脉冲和阶跃响应、🖊计算阶跃响应、🖊从阶跃响应获取脉冲响应、🎯绘制:🖊非递归系统的零点、🖊非递归系统的幅度和相位响应

    4. 数控信号处理应用:🎯语音识别、🎯电能传输质量、🎯潘-汤普金斯算法心电图检测、🎯混合模式信号链中噪声分析。

🍇Python傅立叶分析加速度计信号去噪

傅里叶分析的思想是任何时间序列都可以分解为不同频率的谐波的积分和。因此,理论上,我们可以利用多个谐波来生成任何信号。在区间 ( − T / 2 < t < T / 2 - T / 2< t < T / 2 T/2<t<T/2) 上定义的任意时间函数 f ( t ) 的傅立叶级数
f ( t ) = a 0 + ∑ n = 1 ∞ a n c o s ( 2 n π t T ) + ∑ n = 1 ∞ b n s i n ( 2 n π t T ) f(t) = a_0 + \sum_{n=1}^{\infty} a_n cos(\frac{2n\pi t}{T}) + \sum_{n=1}^{\infty} b_n sin(\frac{2n\pi t}{T}) f(t)=a0+n=1ancos(T2t)+n=1bnsin(T2t)
在上式中,我们可以看到 s i n ( 2 n π t T ) sin(\frac{2n\pi t}{T}) sin(T2t) c o s ( 2 n π t T ) cos(\frac{2n\pi t}{T}) cos(T2t) 是周期性的,周期为 T。

我们将使用傅立叶分析进行滤波,假设噪声在时域中与信号重叠,但在频域中不那么重叠。

import pandas as pd
import os, sys
import numpy as np
import matplotlib.pyplot as plt
plt.rcParams['figure.figsize'] = [10,6]
plt.rcParams.update({'font.size': 18})
plt.style.use('seaborn')dt = 0.001
t = np.arange(0, 1, dt)
signal = np.sin(2*np.pi*50*t) + np.sin(2*np.pi*120*t) #composite signal
signal_clean = signal #copy for later comparison
signal = signal + 2.5 * np.random.randn(len(t))
minsignal, maxsignal = signal.min(), signal.max()

我们通过对不同频率(50Hz 和 120Hz)的两个正弦函数求和来创建信号。然后我们创建了一组随机噪声并将该噪声叠加到信号上。

快速傅里叶变换

n = len(t)
fhat = np.fft.fft(signal, n) 
psd = fhat * np.conj(fhat)/n
freq = (1/(dt*n)) * np.arange(n) 
idxs_half = np.arange(1, np.floor(n/2), dtype=np.int32) 

Numpy 的 fft.fft 函数使用高效的快速傅里叶变换 (FFT) 算法返回一维离散傅里叶变换。 该函数的输出是复数,我们将其与其共轭相乘以获得噪声信号的功率谱。 我们使用采样间隔 (dt) 和样本数 (n) 创建频率数组。

在上图中,我们可以看到原始信号的两个频率很突出。现在,我们可以创建一个滤波器,可以删除幅度小于阈值的所有频率。

## Filter out noise
threshold = 100
psd_idxs = psd > threshold #array of 0 and 1
psd_clean = psd * psd_idxs #zero out all the unnecessary powers
fhat_clean = psd_idxs * fhat #used to retrieve the signalsignal_filtered = np.fft.ifft(fhat_clean) #inverse fourier transform
## Visualization
fig, ax = plt.subplots(4,1)
ax[0].plot(t, signal, color='b', lw=0.5, label='Noisy Signal')
ax[0].plot(t, signal_clean, color='r', lw=1, label='Clean Signal')
ax[0].set_ylim([minsignal, maxsignal])
ax[0].set_xlabel('t axis')
ax[0].set_ylabel('Vals')
ax[0].legend()ax[1].plot(freq[idxs_half], np.abs(psd[idxs_half]), color='b', lw=0.5, label='PSD noisy')
ax[1].set_xlabel('Frequencies in Hz')
ax[1].set_ylabel('Amplitude')
ax[1].legend()ax[2].plot(freq[idxs_half], np.abs(psd_clean[idxs_half]), color='r', lw=1, label='PSD clean')
ax[2].set_xlabel('Frequencies in Hz')
ax[2].set_ylabel('Amplitude')
ax[2].legend()ax[3].plot(t, signal_filtered, color='r', lw=1, label='Clean Signal Retrieved')
ax[3].set_ylim([minsignal, maxsignal])
ax[3].set_xlabel('t axis')
ax[3].set_ylabel('Vals')
ax[3].legend()plt.subplots_adjust(hspace=0.4)
plt.savefig('signal-analysis.png', bbox_inches='tight', dpi=300)

功率阈值对真实数据进行去噪

import numpy as np
import matplotlib.pyplot as plt
plt.rcParams['figure.figsize'] = [10,6]
plt.rcParams.update({'font.size': 18})
plt.style.use('seaborn')otime = UTCDateTime('2021-04-18T22:14:37')filenameZ = 'TRCEC7A.msd'stZ = read(filenameZ)
streams = [stZ.copy()]
traces = []
for st in streams:tr = st[0].trim(otime, otime+120)traces.append(tr)delta = stZ[0].stats.delta
starttime = np.datetime64(stZ[0].stats.starttime)
endtime = np.datetime64(stZ[0].stats.endtime)
signalZ = traces[0].data/10**6
minsignal, maxsignal = signalZ.min(), signalZ.max()t = traces[0].times("utcdatetime") n = len(t)
fhat = np.fft.fft(signalZ, n) 
psd = fhat * np.conj(fhat)/n
freq = (1/(delta*n)) * np.arange(n) 
idxs_half = np.arange(1, np.floor(n/2), dtype=np.int32) 
psd_real = np.abs(psd[idxs_half]) sort_psd = np.sort(psd_real)[::-1]
threshold = sort_psd[300]
psd_idxs = psd > threshold 
psd_clean = psd * psd_idxs 
fhat_clean = psd_idxs * fhat signal_filtered = np.fft.ifft(fhat_clean) fig, ax = plt.subplots(4,1)
ax[0].plot(t, signalZ, color='b', lw=0.5, label='Noisy Signal')
ax[0].set_xlabel('t axis')
ax[0].set_ylabel('Accn in Gal')
ax[0].legend()ax[1].plot(freq[idxs_half], np.abs(psd[idxs_half]), color='b', lw=0.5, label='PSD noisy')
ax[1].set_xlabel('Frequencies in Hz')
ax[1].set_ylabel('Amplitude')
ax[1].legend()ax[2].plot(freq[idxs_half], np.abs(psd_clean[idxs_half]), color='r', lw=1, label='PSD clean')
ax[2].set_xlabel('Frequencies in Hz')
ax[2].set_ylabel('Amplitude')
ax[2].legend()ax[3].plot(t, signal_filtered, color='r', lw=1, label='Clean Signal Retrieved')
ax[3].set_ylim([minsignal, maxsignal])
ax[3].set_xlabel('t axis')
ax[3].set_ylabel('Accn in Gal')
ax[3].legend()plt.subplots_adjust(hspace=0.6)
plt.savefig('real-signal-analysis.png', bbox_inches='tight', dpi=300)
参阅一:计算思维
参阅二:亚图跨际

这篇关于Python量化噪声卷积信号和傅里叶时频分析的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!



http://www.chinasem.cn/article/892420

相关文章

Python将博客内容html导出为Markdown格式

《Python将博客内容html导出为Markdown格式》Python将博客内容html导出为Markdown格式,通过博客url地址抓取文章,分析并提取出文章标题和内容,将内容构建成html,再转... 目录一、为什么要搞?二、准备如何搞?三、说搞咱就搞!抓取文章提取内容构建html转存markdown

Python获取中国节假日数据记录入JSON文件

《Python获取中国节假日数据记录入JSON文件》项目系统内置的日历应用为了提升用户体验,特别设置了在调休日期显示“休”的UI图标功能,那么问题是这些调休数据从哪里来呢?我尝试一种更为智能的方法:P... 目录节假日数据获取存入jsON文件节假日数据读取封装完整代码项目系统内置的日历应用为了提升用户体验,

Go标准库常见错误分析和解决办法

《Go标准库常见错误分析和解决办法》Go语言的标准库为开发者提供了丰富且高效的工具,涵盖了从网络编程到文件操作等各个方面,然而,标准库虽好,使用不当却可能适得其反,正所谓工欲善其事,必先利其器,本文将... 目录1. 使用了错误的time.Duration2. time.After导致的内存泄漏3. jsO

Python FastAPI+Celery+RabbitMQ实现分布式图片水印处理系统

《PythonFastAPI+Celery+RabbitMQ实现分布式图片水印处理系统》这篇文章主要为大家详细介绍了PythonFastAPI如何结合Celery以及RabbitMQ实现简单的分布式... 实现思路FastAPI 服务器Celery 任务队列RabbitMQ 作为消息代理定时任务处理完整

Python Websockets库的使用指南

《PythonWebsockets库的使用指南》pythonwebsockets库是一个用于创建WebSocket服务器和客户端的Python库,它提供了一种简单的方式来实现实时通信,支持异步和同步... 目录一、WebSocket 简介二、python 的 websockets 库安装三、完整代码示例1.

揭秘Python Socket网络编程的7种硬核用法

《揭秘PythonSocket网络编程的7种硬核用法》Socket不仅能做聊天室,还能干一大堆硬核操作,这篇文章就带大家看看Python网络编程的7种超实用玩法,感兴趣的小伙伴可以跟随小编一起... 目录1.端口扫描器:探测开放端口2.简易 HTTP 服务器:10 秒搭个网页3.局域网游戏:多人联机对战4.

使用Python实现快速搭建本地HTTP服务器

《使用Python实现快速搭建本地HTTP服务器》:本文主要介绍如何使用Python快速搭建本地HTTP服务器,轻松实现一键HTTP文件共享,同时结合二维码技术,让访问更简单,感兴趣的小伙伴可以了... 目录1. 概述2. 快速搭建 HTTP 文件共享服务2.1 核心思路2.2 代码实现2.3 代码解读3.

Python使用自带的base64库进行base64编码和解码

《Python使用自带的base64库进行base64编码和解码》在Python中,处理数据的编码和解码是数据传输和存储中非常普遍的需求,其中,Base64是一种常用的编码方案,本文我将详细介绍如何使... 目录引言使用python的base64库进行编码和解码编码函数解码函数Base64编码的应用场景注意

Python基于wxPython和FFmpeg开发一个视频标签工具

《Python基于wxPython和FFmpeg开发一个视频标签工具》在当今数字媒体时代,视频内容的管理和标记变得越来越重要,无论是研究人员需要对实验视频进行时间点标记,还是个人用户希望对家庭视频进行... 目录引言1. 应用概述2. 技术栈分析2.1 核心库和模块2.2 wxpython作为GUI选择的优

Python如何使用__slots__实现节省内存和性能优化

《Python如何使用__slots__实现节省内存和性能优化》你有想过,一个小小的__slots__能让你的Python类内存消耗直接减半吗,没错,今天咱们要聊的就是这个让人眼前一亮的技巧,感兴趣的... 目录背景:内存吃得满满的类__slots__:你的内存管理小助手举个大概的例子:看看效果如何?1.