Python 简单的信号处理

混合正弦信号处理

import numpy as np
import matplotlib.pyplot as plt
from numpy.fft import rfft,rfftfreq

Pi2 = 2 * np.pi
framerate = 10000

def gene_fake_signal():
    #采样率
    global framerate
    #采样时长
    duration = 0.01
    #样本数
    n = round(duration*framerate)
    #采样时间点
    ts = np.arange(n)/framerate
    #频率
    freq1 = 1900
    freq2 = 2300
    #振幅
    A1 = 1
    A2 = 2
    #相位偏移
    offset1 = 10
    offset2 = 20

    ys1 = A1 * np.cos(Pi2*freq1*ts  + offset1)
    ys2 = A2 * np.cos(Pi2*freq2*ts + offset2)
    ys = ys1 + ys2
    return ts,ys

ts,ys = gene_fake_signal()

fig = plt.figure(figsize=(6,3))
ax1 = plt.subplot(121)
ax1.plot(ts,ys,"r-")
ax1.set_title("original signal")

ax2 = plt.subplot(122)

#rfft 与 rfftfreq
hs = np.abs(rfft(ys))
fs = rfftfreq(len(ys),1/framerate)

Data = np.array([[hs[i],fs[i]] for i in range(len(hs))])
temp = Data[:,0]
index = np.lexsort((temp,))
Data = Data[index]

ax2.plot(fs,hs)
plt.savefig("0.jpg")
plt.pause(0.01)

       [1.79449838e-13, 2.40000000e+03],
       [2.61232155e-13, 1.00000000e+03],
       [5.00000000e+01, 1.90000000e+03],
       [1.00000000e+02, 2.30000000e+03]])

随机噪声对正弦信号干扰

示例：单通道+噪声干扰

import numpy as np
import matplotlib.pyplot as plt
from numpy.fft import rfft,rfftfreq

Pi2 = 2 * np.pi
framerate = 10000

def gene_fake_signal():
    #采样率
    global framerate
    #采样时长
    duration = 0.01
    #样本数
    n = round(duration*framerate)
    #采样时间点
    ts = np.arange(n)/framerate
    #频率
    freq1 = 1900
    #振幅
    A1 = 1
    #相位偏移
    offset1 = 10

    rs = np.random.uniform(0,1,n)
    ys1 = A1 * np.cos(Pi2*freq1*ts  + offset1)
    ys = ys1 + rs 
    return ts,ys

ts,ys = gene_fake_signal()

fig = plt.figure(figsize=(6,3))
ax1 = plt.subplot(121)
ax1.plot(ts,ys,"r-")
ax1.set_title("original signal")

ax2 = plt.subplot(122)
hs = np.abs(rfft(ys))
fs = rfftfreq(len(ys),1/framerate)
Data = np.array([[hs[i],fs[i]] for i in range(len(hs))])
temp = Data[:,0]
index = np.lexsort((temp,))
Data = Data[index]

ax2.plot(fs,hs)
plt.savefig("0.jpg")
plt.pause(0.01)

示例：多通道+噪声干扰

import numpy as np
import matplotlib.pyplot as plt
from numpy.fft import rfft,rfftfreq

Pi2 = 2 * np.pi
framerate = 10000

def gene_fake_signal():
    #采样率
    global framerate
    #采样时长
    duration = 0.01
    #样本数
    n = round(duration*framerate)
    #采样时间点
    ts = np.arange(n)/framerate
    #频率
    freq1 = 1900
    freq2 = 2000
    freq3 = 2500
    #振幅
    A1 = 1
    A2 = 2
    A3 = 3
    #相位偏移
    offset1 = 10
    offset2 = 20
    offset3 = 30
    
    rs = np.random.uniform(0,0.1,n)
    ys1 = A1 * np.cos(Pi2*freq1*ts+offset1)
    ys2 = A2 * np.cos(Pi2*freq2*ts+offset2)
    ys3 = A3 * np.cos(Pi2*freq3*ts+offset3)
    
    ys = ys1 + ys2 + ys3 + rs 
    return ts,ys

ts,ys = gene_fake_signal()

fig = plt.figure(figsize=(6,3))
ax1 = plt.subplot(121)
ax1.plot(ts,ys,"r-")
ax1.set_title("original signal")

ax2 = plt.subplot(122)

hs = np.abs(rfft(ys))
fs = rfftfreq(len(ys),1/framerate)

Data = np.array([[hs[i],fs[i]] for i in range(len(hs))])
temp = Data[:,0]
index = np.lexsort((temp,))
Data = Data[index]

ax2.plot(fs,hs)
plt.savefig("0.jpg")
plt.pause(0.01)

       [5.65482521e-01, 3.60000000e+03],
       [5.75495030e-01, 3.90000000e+03],
       [5.48305771e+00, 0.00000000e+00],
       [4.97797723e+01, 1.90000000e+03],
       [1.00382508e+02, 2.00000000e+03],
       [1.50113033e+02, 2.50000000e+03]])

信号转换

• 现有一方程，其数值解如下：

• 给出源代码

import numpy as np
import matplotlib.pyplot as plt
from numpy.fft import rfft,rfftfreq

from math import *

def runge_kutta4(df,a,b,h,y0):
    num = len(y0)
    t = np.arange(a,b+h,h)
    w = np.zeros((t.size,num))
    w[0,:] = y0
 
    for i in range(t.size - 1):
        s0 = df(t[i], w[i,:])
        s1 = df(t[i] + h/2.,w[i,:] + h * s0 / 2.)
        s2 = df(t[i] + h/2.,w[i,:] + h * s1 / 2.)
        s3 = df(t[i+1], w[i,:] + h * s2)
        w[i+1,:] = w[i,:] + h * (s0 + 2*(s1+s2) + s3) / 6.
 
    return t,w
 
def df(t, variables):
    x,y,z,x1,y1,z1 = variables
    
    A = np.zeros((3,3))
    b = np.zeros(3)
    
    A[0,0] = 1
    A[0,1] = 0
    A[0,2] = 0
    
    A[1,0] = 0
    A[1,1] = 1
    A[1,2] = 0
 
    A[2,0] = 0
    A[2,1] = 0
    A[2,2] = 1
 
 
    r1 = 2 * 10**-3
    r = 14.5 * 10**-3
    m = 88.25 * 10**-3
    G = 53.4
    E = 268153
    k = 0.5
    u = 0.35
    A1 = r1**2
    L = 10 * 10**-2
    I = 0.4 * m * r**2
    g = 9.8
    
    b[0] = -u*g*sin(atan(x1/x/y1))-E*A1/m/L*2*(x-L/2)+(x+r)*x1**2
    b[1] = (-u*g*cos(atan(x1/x/y1))-2*x1*y1)/x 
    b[2] = -G*3.14* r1**3 * r * z/32/I/x + r*u*m*g*cos(atan(x1/x/y1))/I - k*m*g/I
    
    x2,y2,z2 = np.linalg.solve(A, b)
 
    return np.array([x1,y1,z1,x2,y2,z2])



def gene_fake_signal():
    a,b = 0.0,15.0
    h = 0.01
     
    x = 1/2 * 10**-2
    y = 0
    z = 100*2*3.14
    x1 = 0.1
    y1 = 0.000001
    z1 = 0.000001
     
    y0 = np.array([x, y, z, x1, y1, z1])
    t,w = runge_kutta4(df, a, b, h, y0)
    t = t
    x = w[:, 0]
    y = w[:, 1]
    z = w[:, 2]
    x1 = w[:, 3]
    y1 = w[:, 4]
    z1 = w[:, 5]
     
    fig,axes = plt.subplots(3,2)
    axes[0][0].plot(t, x)
    axes[0][1].plot(t, y)
    axes[1][0].plot(t, z)
    axes[1][1].plot(t, x1)
    axes[2][0].plot(t, y1)
    axes[2][1].plot(t, z1)
    plt.pause(0.01)

    return t, x



def RFFTFREQ(ts, ys):

    framerate = 100
    
    fig = plt.figure(figsize=(6, 3))
    ax1 = plt.subplot(121)
    ax1.plot(ts, ys, "r-")
    ax1.set_title("original signal")

    ax2 = plt.subplot(122)

    #rfft 与 rfftfreq
    hs = np.abs(rfft(ys))
    fs = rfftfreq(len(ys), 1/framerate)

    Data = np.array([[hs[i], fs[i]] for i in range(len(hs))])
    temp = Data[:,0]
    index = np.lexsort((temp, ))
    Data = Data[index]
    sorted_tuples = sorted(zip(fs, hs), key=lambda x: x[1], reverse=True)
    sorted_fs = [t[0] for t in sorted_tuples]
    print(sorted_fs[:5])
    ax2.plot(fs, hs)
    plt.savefig("0.jpg")
    plt.pause(0.01)


ts,ys = gene_fake_signal()
RFFTFREQ(ts, ys)

 

现在有两个长度相同的Python 列表fs hs ，fs 和 hs 一一对应。fs中元素表示对应位置上hs列表中的元素的编号。现在将hs 中的元素从大到小排序，给出顺序改变后的fs列表

fs = [1, 2, 3, 4, 5, 6, 7, 8, 9, 10]
hs = [423, 324, 27, 43, 32, 2, 74, 7, 53, 4]

sorted_tuples = sorted(zip(fs, hs), key=lambda x: x[1], reverse=True)
sorted_fs = [t[0] for t in sorted_tuples]

>>> sorted_fs
[1, 2, 7, 9, 4, 5, 3, 8, 10, 6]
>>> 
 

with water

import numpy as np
import matplotlib.pyplot as plt
from numpy.fft import rfft,rfftfreq

from math import *

def runge_kutta4(df,a,b,h,y0):
    num = len(y0)
    t = np.arange(a,b+h,h)
    w = np.zeros((t.size,num))
    w[0,:] = y0
 
    for i in range(t.size - 1):
        s0 = df(t[i], w[i,:])
        s1 = df(t[i] + h/2.,w[i,:] + h * s0 / 2.)
        s2 = df(t[i] + h/2.,w[i,:] + h * s1 / 2.)
        s3 = df(t[i+1], w[i,:] + h * s2)
        w[i+1,:] = w[i,:] + h * (s0 + 2*(s1+s2) + s3) / 6.
 
    return t,w
 
def df(t, variables):
    x,y,z,x1,y1,z1 = variables
    
    A = np.zeros((3,3))
    b = np.zeros(3)
    
    A[0,0] = 1
    A[0,1] = 0
    A[0,2] = 0
    
    A[1,0] = 0
    A[1,1] = 1
    A[1,2] = 0
 
    A[2,0] = 0
    A[2,1] = 0
    A[2,2] = 1
 
 
    r1 = 2 * 10**-3
    r = 14.5 * 10**-3
    m = 88.25 * 10**-3
    G = 53.4
    E = 268153
    k = 0.5
    u = 0.35
    A1 = r1**2
    L = 10 * 10**-2
    I = 0.4 * m * r**2
    g = 9.8
    
    b[0] = -u*g*sin(atan(x1/x/y1))-E*A1/m/L*2*(x-L/2)+(x+r)*x1**2
    b[1] = (-u*g*cos(atan(x1/x/y1))-2*x1*y1)/x 
    b[2] = -G*3.14* r1**3 * r * z/32/I/x + r*u*m*g*cos(atan(x1/x/y1))/I - k*m*g/I
    
    x2,y2,z2 = np.linalg.solve(A, b)
 
    return np.array([x1,y1,z1,x2,y2,z2])

def gene_fake_signal():
    a,b = 0.0,15.0
    h = 0.01
     
    x = 1/2 * 10**-2
    y = 0
    z = 100*2*3.14
    x1 = 0.1
    y1 = 0.000001
    z1 = 0.000001
     
    y0 = np.array([x, y, z, x1, y1, z1])
    t,w = runge_kutta4(df, a, b, h, y0)
    t = t
    x = w[:, 0]
    y = w[:, 1]
    z = w[:, 2]
    x1 = w[:, 3]
    y1 = w[:, 4]
    z1 = w[:, 5]
     
    fig,axes = plt.subplots(3,2)
    axes[0][0].plot(t, x)
    axes[0][1].plot(t, y)
    axes[1][0].plot(t, z)
    axes[1][1].plot(t, x1)
    axes[2][0].plot(t, y1)
    axes[2][1].plot(t, z1)
    plt.pause(0.01)

    return t, x

def RFFTFREQ(ts, ys):

    framerate = 100
    
    fig = plt.figure(figsize=(6, 3))
    ax1 = plt.subplot(121)
    ax1.plot(ts, ys, "r-")
    ax1.set_title("original signal")

    ax2 = plt.subplot(122)

    #rfft 与 rfftfreq
    hs = np.abs(rfft(ys))
    fs = rfftfreq(len(ys), 1/framerate)

    Data = np.array([[hs[i], fs[i]] for i in range(len(hs))])
    temp = Data[:,0]
    index = np.lexsort((temp, ))
    Data = Data[index]
    sorted_tuples = sorted(zip(fs, hs), key=lambda x: x[1], reverse=True)
    sorted_fs = [t[0] for t in sorted_tuples]
    print(sorted_fs[:5])
    ax2.plot(fs, hs)
    plt.savefig("0.jpg")
    plt.pause(0.01)

ts,ys = gene_fake_signal()
RFFTFREQ(ts, ys)