半波还原整波——封装函数

import matplotlib.pyplot as plt
import numpy as np
from scipy.fftpack import fft, ifft #FFT库
import xlrd     #读取excel库

# 打开采样数据文件
file_name = r'C:\Users\Chaver\Desktop\task\sampled.xlsx'
wb = xlrd.open_workbook(file_name)
# 根据 sheet 名称获取内容
sh1 = wb.sheet_by_name('3V-R-65Hz-data1')

# 读取采集A、B、C 线电压值 （A:Uab，B:Ubc,C:Uca ）
phase_A = sh1.col_values(1,1,257) # A相电压
phase_B = sh1.col_values(2,1,257) # B相电压


# 上零点函数
def find_up(data):
    up = []
    for x in range(2,len(data)-1):
        if data[x-2] < 50 and data[x-1] < 50 and data[x] < 50 and data[x+1] - data[x] > 50:
            up.append(x)
    return up


# 下零点函数
def find_dw(data):
    dw = []
    for x in range(1,len(data)-2):
        if data[x+2] < 50 and data[x+1] < 50 and data[x] < 50 and data[x] - data[x-1] < -40:
            dw.append(x)
    return dw


# 获取相位差
def angle_error(phase_A,phase_B):
    # 加窗 Hanning
    phase_A_window = np.hanning(len(phase_A))*phase_A
    phase_B_window = np.hanning(len(phase_B))*phase_B
    # fft计算
    F_phase_A = fft(phase_A_window)
    F_phase_B = fft(phase_B_window)
    # 角度计算
    angle_A = np.angle(F_phase_A,deg=True)
    angle_B = np.angle(F_phase_B,deg=True)
    # 相位差
    Delta_BA = np.subtract(angle_B,angle_A)
    Delta_BA = [x+360 if x<-180 else 360-x if x>180 else x for x in Delta_BA]
    return Delta_BA


# 获取还原前角度差
angle_original = angle_error(phase_A, phase_B)
print("还原前角度差：", angle_original[1:4])

# 获取A,B零点位置
up_A = find_up(phase_A)
dw_A = find_dw(phase_A)
up_B = find_up(phase_B)
dw_B = find_dw(phase_B)


# 还原全波
# 原始数据近零数值置零
def pre_data(data):
    for x in range(len(data)):
        if data[x] < 50:
            data[x] = 0
    return data


# 获取完整的第一个半波
def zero(D, up, dw):
    for x in range(len(D)):
        if x < up or x > dw:
            D[x] = 0
    return D


# 还原全波主函数
def wave_recovery(data_1,up,dw):
    if up[0] < dw[0]:
        half_length = dw[0] - up[0]
        # 生成新列表data_2用来存放平移反转后的列表data_1
        data_2 = [0] * dw[0] + [None] * (len(data_1) - dw[0])
        for x in range(dw[0], len(data_1)):
            data_2[x] = data_1[x - half_length]
            data_2[x] = -data_2[x]
        # 上零点左侧曲线复原，data_3
        data_3 = [0] * up[0]
        for x in range(up[0], 2 * up[0]):
            data_3[x - up[0]] = -data_1[x]
        # 将元素逆序呈现
        data_3 = data_3[::-1]
        # 补齐列表长度，添加零元素
        data_3 = data_3 + [0] * (len(data_1) - len(data_3))
        # 合并data_1与data_2
        sum_1 = np.add(data_1, data_2)
        phase_A = np.add(sum_1, data_3)
    else:
        half_length = up[0] - dw[0]
        data_2 = [0] * dw[1] + [None] * (len(data_1) - dw[1])
        for x in range(dw[1], len(data_1)):
            data_2[x] = data_1[x - half_length]
            data_2[x] = -data_2[x]

        # 合并data_1和data_2，不知为什么一定要放在data_3的for循环前？
        sum_data_1_2 = np.add(data_1, data_2)
        # 获取第一个完整半波
        data_3 = zero(data_1, up[0], dw[1])
        # 将第一个半波向左平移half_length反转,得到data_4
        data_4 = [0] * dw[0] + [0] * (len(data_3) - dw[0])
        for x in range(dw[0], up[0]):
            data_4[x] = data_3[x + half_length]
            data_4[x] = -data_4[x]
        # 合并data_1与data_2,data_4
        phase_A = np.add(sum_data_1_2, data_4)
    return phase_A


# 获取data1
phase_A_pre = pre_data(phase_A)
phase_B_pre = pre_data(phase_B)

# 获取还原后数据
phase_A_recovery=wave_recovery(phase_A_pre,up_A,dw_A)
phase_B_recovery=wave_recovery(phase_B_pre,up_B,dw_B)

plt.subplot(2,1,1)
plt.plot(phase_A_recovery,marker='o')
plt.subplot(2,1,2)
plt.plot(phase_B_recovery,marker='o')

# 比较还原前与还原后相位差
angle_recovery = angle_error(phase_A_recovery,phase_B_recovery)
print("还原后角度差: ", angle_recovery[1:4])

plt.show()