Autoformer¶

# 获取幅度谱 - 这里乘以2是因为能量在正负频率上分布
# 对于rfft，我们只看到正频率部分
magnitude_spectrum = torch.abs(signal_rfft) * 2 / N  # 归一化
# 对于第0个和第N/2个频率分量不需要乘以2
if N % 2 == 0:  # 如果N是偶数
    magnitude_spectrum[0] /= 2
    magnitude_spectrum[-1] /= 2
magnitude_spectrum

# 获取频率轴
frequencies = torch.fft.rfftfreq(N, 1.0/sample_rate)
frequencies

# 详细打印结果
print(f"采样率: {sample_rate}Hz")
print(f"信号长度: {N}点")
print(f"频率分辨率: {sample_rate/N}Hz")
print("\n前25个频率点及其幅度:")
for i in range(min(25, len(frequencies))):
    print(f"频率 {frequencies[i]:.1f}Hz: 幅度 {magnitude_spectrum[i]:.4f}")

# 可视化结果    
plt.figure(figsize=(15, 10))

# 原始信号
plt.subplot(3, 1, 1)
plt.plot(t.numpy(), signal.numpy())
plt.title('原始信号 - 10Hz (幅度1.0) + 20Hz (幅度0.5)的正弦波 1*sin(2π*10*t) + 0.5*sin(2π*20*t)')
plt.xlabel('时间 (秒)')
plt.ylabel('幅度')
plt.grid(True)

# 频谱幅度 - 线性刻度
plt.subplot(3, 1, 2)
plt.plot(frequencies.numpy(), magnitude_spectrum.numpy())
plt.title('频谱幅度 (线性刻度)')
plt.xlabel('频率 (Hz)')
plt.ylabel('幅度')
plt.grid(True)
# 标注10Hz和20Hz的峰值
plt.axvline(x=10, color='r', linestyle='--', alpha=0.7)
plt.axvline(x=20, color='g', linestyle='--', alpha=0.7)
plt.text(10, 0.9, '10Hz', color='r')
plt.text(20, 0.45, '20Hz', color='g')

# 频谱幅度 - 对数刻度，便于观察小幅值
plt.subplot(3, 1, 3)
plt.semilogy(frequencies.numpy(), magnitude_spectrum.numpy())
plt.title('频谱幅度 (对数刻度)')
plt.xlabel('频率 (Hz)')
plt.ylabel('幅度 (log)')
plt.grid(True)
plt.axvline(x=10, color='r', linestyle='--', alpha=0.7)
plt.axvline(x=20, color='g', linestyle='--', alpha=0.7)

plt.tight_layout()
plt.show()

# 理论验证 - 重建信号
# 创建理想频谱
ideal_spectrum = torch.zeros_like(signal_rfft)
ideal_spectrum[10] = amplitude1 * N / 2 * np.exp(-0.5j * np.pi)  # 10Hz, 相位为-π/2
ideal_spectrum[20] = amplitude2 * N / 2 * np.exp(-0.5j * np.pi)  # 20Hz, 相位为-π/2

# 逆变换重建
reconstructed = torch.fft.irfft(ideal_spectrum, n=N)

# 比较原始信号和重建信号
plt.figure(figsize=(12, 6))
plt.plot(t.numpy(), signal.numpy(), 'b-', label='原始信号')
plt.plot(t.numpy(), reconstructed.numpy(), 'r--', label='从频谱重建的信号')
plt.title('原始信号与重建信号比较')
plt.legend()
plt.grid(True)
plt.show()

# 计算重建误差
error = torch.max(torch.abs(signal - reconstructed))
print(f"\n重建误差(最大绝对差): {error:.10f}")

from scipy.fft import fft
fft([4,3,2,1])

# 模拟Autoformer的时间序列分解函数
def series_decomp(data, kernel_size=25):
    # 使用卷积实现移动平均
    weights = torch.ones(kernel_size) / kernel_size
    weights = weights.view(1, 1, -1).to(torch.float32)
    data_pad = F.pad(torch.from_numpy(data).float().view(1, 1, -1), 
                    (kernel_size//2, kernel_size//2), mode='replicate')
    trend = F.conv1d(data_pad, weights).view(-1).numpy()
    seasonal = data - trend
    return seasonal, trend

# 使用FFT分析季节性成分中的谐波
def analyze_harmonics(seasonal, sampling_rate=1):
    # 执行FFT
    fft_result = np.fft.rfft(seasonal)
    # 计算频率
    freqs = np.fft.rfftfreq(len(seasonal), d=1/sampling_rate)
    # 计算幅度
    magnitudes = np.abs(fft_result) * 2 / len(seasonal)
    # 计算相位 (转换为度)
    phases = np.angle(fft_result) * 180 / np.pi
    return freqs, magnitudes, phases

# 分解时间序列
seasonal, trend = series_decomp(power_demand)

# 分析谐波
freqs, magnitudes, phases = analyze_harmonics(seasonal)

# 找出主要谐波 (幅度最大的几个)
n_harmonics = 5
top_indices = np.argsort(magnitudes)[-n_harmonics:][::-1]
main_freqs = freqs[top_indices]
main_periods = 1/main_freqs[main_freqs > 0]  # 避免除以零
main_magnitudes = magnitudes[top_indices]
main_phases = phases[top_indices]

# 重构信号，仅使用主要谐波
reconstructed_seasonal = np.zeros_like(seasonal)
t = np.arange(len(seasonal))
for i, idx in enumerate(top_indices):
    if freqs[idx] > 0:  # 跳过直流分量
        reconstructed_seasonal += main_magnitudes[i] * np.cos(2*np.pi*freqs[idx]*t + np.radians(main_phases[i]))

# 可视化
plt.figure(figsize=(18, 12))

# 原始时间序列
plt.subplot(3, 2, 1)
plt.plot(hours, power_demand)
plt.title('原始电力需求时间序列')
plt.xlabel('小时')
plt.ylabel('功率 (MW)')
plt.grid(True)

# 趋势和季节性分解
plt.subplot(3, 2, 2)
plt.plot(hours, trend, label='趋势')
plt.plot(hours, seasonal, label='季节性')
plt.title('时间序列分解')
plt.xlabel('小时')
plt.ylabel('功率 (MW)')
plt.legend()
plt.grid(True)

# 频谱分析
plt.subplot(3, 2, 3)
plt.stem(freqs[:50], magnitudes[:50])  # 只显示前50个频率
plt.title('季节性成分的频谱')
plt.xlabel('频率 (周期/小时)')
plt.ylabel('幅度')
plt.grid(True)

# 相位分析
plt.subplot(3, 2, 4)
plt.stem(freqs[:50], phases[:50])
plt.title('季节性成分的相位谱')
plt.xlabel('频率 (周期/小时)')
plt.ylabel('相位 (度)')
plt.grid(True)

# 主要谐波
plt.subplot(3, 2, 5)
for i, idx in enumerate(top_indices):
    if freqs[idx] > 0:  # 跳过直流分量
        period = 1/freqs[idx]
        plt.plot(hours, main_magnitudes[i] * np.cos(2*np.pi*freqs[idx]*t + np.radians(main_phases[i])),
                 label=f'{period:.1f}小时周期')
plt.title('主要谐波')
plt.xlabel('小时')
plt.ylabel('幅度')
plt.legend()
plt.grid(True)

# 原始季节性与重构比较
plt.subplot(3, 2, 6)
plt.plot(hours, seasonal, label='原始季节性')
plt.plot(hours, reconstructed_seasonal, label='谐波重构', linestyle='--')
plt.title('季节性成分与谐波重构比较')
plt.xlabel('小时')
plt.ylabel('幅度')
plt.legend()
plt.grid(True)

plt.tight_layout()
plt.show()

# 打印主要谐波信息
print("主要谐波:")
for i, idx in enumerate(top_indices):
    if freqs[idx] > 0:
        print(f"- 周期: {1/freqs[idx]:.1f}小时, 幅度: {main_magnitudes[i]:.2f}, 相位: {main_phases[i]:.1f}°")

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.font_manager import FontProperties
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.colors as colors

# 设置中文显示
font_path = '/System/Library/Fonts/Supplemental/Songti.ttc'  # macOS 中的宋体字体
chinese_font = FontProperties(fname=font_path)
plt.rcParams['axes.unicode_minus'] = False  # 解决负号显示问题

# 参数设置
T = 6  # 周期
omega = (2 * np.pi) / T  # 角频率
n = 2  # 周期个数
t = np.arange(0, n * T + 0.1, 0.1)  # 时间序列，从0到n*T，步长为0.1

# 计算复指数信号
f = np.exp(1j * omega * t)
real_part = np.real(f)  # 实部
imag_part = np.imag(f)  # 虚部

# 创建三维图像
fig = plt.figure(figsize=(20, 16))
ax = fig.add_subplot(111, projection='3d')

# 绘制三维曲线
spiral, = ax.plot(t, imag_part, real_part, linewidth=3, color='blue', label='复指数信号')

# 创建投影平面
# xy平面 (z=0)
max_t = n * T
xy_x, xy_y = np.meshgrid(np.linspace(0, max_t, 20), np.linspace(-1.5, 1.5, 20))
xy_z = np.zeros_like(xy_x)
xy_plane = ax.plot_surface(xy_x, xy_y, xy_z, alpha=0.2, color='lightblue', label='xy平面')

# xz平面 (y=0)
xz_x, xz_z = np.meshgrid(np.linspace(0, max_t, 20), np.linspace(-1.5, 1.5, 20))
xz_y = np.zeros_like(xz_x)
xz_plane = ax.plot_surface(xz_x, xz_y, xz_z, alpha=0.2, color='lightgreen', label='xz平面')

# yz平面 (x=0)
yz_y, yz_z = np.meshgrid(np.linspace(-1.5, 1.5, 20), np.linspace(-1.5, 1.5, 20))
yz_x = np.zeros_like(yz_y)
yz_plane = ax.plot_surface(yz_x, yz_y, yz_z, alpha=0.2, color='lightyellow', label='yz平面')

# 绘制投影曲线
xy_proj, = ax.plot(t, imag_part, np.zeros_like(t), 'g--', linewidth=2, label='xy平面投影')
xz_proj, = ax.plot(t, np.zeros_like(t), real_part, 'r--', linewidth=2, label='xz平面投影')
yz_proj, = ax.plot(np.zeros_like(t), imag_part, real_part, 'm--', linewidth=2, label='yz平面投影')

# 绘制单位圆
theta = np.linspace(0, 2*np.pi, 100)
y_circle = np.sin(theta)
z_circle = np.cos(theta)
t_circle = np.zeros_like(theta)
unit_circle, = ax.plot(t_circle, y_circle, z_circle, 'k-', linewidth=2, label='单位圆')

# 添加阴影效果
# 选择几个关键点来画阴影
shadow_indices = np.arange(0, len(t), 10)  # 每10个点取一个
for i in shadow_indices:
    # 从复平面到时间轴的竖线
    ax.plot([0, t[i]], [0, 0], [0, 0], 'k:', alpha=0.3)
    
    # 从时间轴到xy投影的线
    ax.plot([t[i], t[i]], [0, imag_part[i]], [0, 0], 'g:', alpha=0.6)
    
    # 从时间轴到xz投影的线
    ax.plot([t[i], t[i]], [0, 0], [0, real_part[i]], 'r:', alpha=0.6)
    
    # 从投影点到空间曲线的连线
    ax.plot([t[i], t[i]], [imag_part[i], imag_part[i]], [0, real_part[i]], 'b-', alpha=0.4)
    ax.plot([t[i], t[i]], [0, imag_part[i]], [real_part[i], real_part[i]], 'b-', alpha=0.4)
    
    # 在0点处画一个到当前点的扇形阴影
    if i > 0:
        # 从原点到复平面上的点的扇形
        angle = np.linspace(0, omega * t[i], 20)
        x_fan = np.zeros_like(angle)
        y_fan = np.sin(angle)
        z_fan = np.cos(angle)
        ax.plot(x_fan, y_fan, z_fan, 'r-', alpha=0.3)
        
        # 填充扇形区域
        ax.plot_surface(
            np.outer(x_fan, np.ones(2)),
            np.outer(y_fan, np.linspace(0, 1, 2)),
            np.outer(z_fan, np.linspace(0, 1, 2)),
            color='yellow', alpha=0.1
        )

# 标记特殊点
special_times = [0, T/4, T/2, 3*T/4, T]
for st in special_times:
    idx = np.argmin(np.abs(t - st))
    ax.scatter([t[idx]], [imag_part[idx]], [real_part[idx]], color='red', s=100, label=f't={st}' if st == 0 else "")
    ax.text(t[idx]+0.2, imag_part[idx], real_part[idx], f't={st}', fontsize=12)

# 设置坐标轴标签
ax.set_xlabel('x轴：时间 t', fontproperties=chinese_font, fontsize=14)
ax.set_ylabel('y 轴：虚部 Im(f) sin(1*omega*t)', fontproperties=chinese_font, fontsize=14)
ax.set_zlabel('z 轴：实部 Re(f) cos(1*omega*t)', fontproperties=chinese_font, fontsize=14)

# 设置坐标轴范围
ax.set_xlim(0, max_t)
ax.set_ylim(-1.5, 1.5)
ax.set_zlim(-1.5, 1.5)

# 添加网格
ax.grid(True)

# 设置标题
ax.set_title(f'复指数信号 $f=e^{{i\\omega t}}$ 的三维可视化，周期T={T}，n={n}个周期',
            fontproperties=chinese_font, fontsize=16)

# 添加0点标记
ax.scatter([0], [0], [0], color='black', s=100, marker='o')
ax.text(0.2, 0, 0, '原点', fontsize=12, fontproperties=chinese_font)

# 添加图例
# ax.legend(loc='upper right', prop=chinese_font)  # 图例太多可能会挡住图像

# 添加数学公式说明
formula_text = (
    f"复指数信号: $f = e^{{i\\omega t}}$，其中 $\\omega = \\frac{{2\\pi}}{{T}} = \\frac{{2\\pi}}{{{T}}}$\n"
    f"实部: $\\cos(\\omega t)$，虚部: $\\sin(\\omega t)$\n"
    f"周期: T = {T}，总时间: {n*T}"
)
fig.text(0.5, 0.02, formula_text, ha='center', fontproperties=chinese_font, fontsize=14)

# 设置更好的视角
ax.view_init(elev=25, azim=40)

# 调整图形布局
plt.tight_layout(rect=[0, 0.05, 1, 0.95])
plt.show()

# 绘制平面投影的详细视图
plt.figure(figsize=(18, 12))

# 1. 时间域信号
plt.subplot(2, 2, 1)
plt.plot(t, real_part, 'b-', label=r'$\cos(\omega t)$ (实部)')
plt.plot(t, imag_part, 'r-', label=r'$\sin(\omega t)$ (虚部)')
plt.grid(True)
plt.legend(prop=chinese_font)
plt.title('复指数信号的实部和虚部', fontproperties=chinese_font)
plt.xlabel('时间 t', fontproperties=chinese_font)
plt.ylabel('幅度', fontproperties=chinese_font)

# 2. 复平面上的轨迹
plt.subplot(2, 2, 2)
plt.plot(real_part, imag_part, 'g-')
plt.grid(True)
plt.title('复平面轨迹 (单位圆)', fontproperties=chinese_font)
plt.xlabel('实部', fontproperties=chinese_font)
plt.ylabel('虚部', fontproperties=chinese_font)
plt.axis('equal')  # 保持坐标轴比例一致

# 添加阴影扇区
for st in special_times:
    idx = np.argmin(np.abs(t - st))
    angle = np.linspace(0, omega * t[idx], 100)
    x_fan = np.cos(angle)
    y_fan = np.sin(angle)
    plt.fill(x_fan, y_fan, 'yellow', alpha=0.1)
    plt.plot(x_fan, y_fan, 'k--', alpha=0.3)
    
# 标记特殊点
for st in special_times:
    idx = np.argmin(np.abs(t - st))
    plt.plot(real_part[idx], imag_part[idx], 'ro')
    plt.annotate(f't = {st}', 
                xy=(real_part[idx], imag_part[idx]),
                xytext=(real_part[idx]*1.2, imag_part[idx]*1.2),
                fontproperties=chinese_font,
                fontsize=10)

# 3. 三维螺旋的俯视图
plt.subplot(2, 2, 3)
plt.plot(t, imag_part, 'g-')
plt.grid(True)
plt.title('三维螺旋的俯视图 (t-虚部 平面)', fontproperties=chinese_font)
plt.xlabel('时间 t', fontproperties=chinese_font)
plt.ylabel('虚部', fontproperties=chinese_font)

# 添加从时间轴到曲线的垂线阴影
shadow_indices = np.arange(0, len(t), 20)
for i in shadow_indices:
    plt.plot([t[i], t[i]], [0, imag_part[i]], 'k:', alpha=0.5)
    
# 标记特殊点
for st in special_times:
    idx = np.argmin(np.abs(t - st))
    plt.plot(t[idx], imag_part[idx], 'ro')
    plt.annotate(f't = {st}', 
                xy=(t[idx], imag_part[idx]),
                xytext=(t[idx]+0.2, imag_part[idx]+0.1),
                fontproperties=chinese_font,
                fontsize=10)

# 4. 三维螺旋的侧视图
plt.subplot(2, 2, 4)
plt.plot(t, real_part, 'r-')
plt.grid(True)
plt.title('三维螺旋的侧视图 (t-实部 平面)', fontproperties=chinese_font)
plt.xlabel('时间 t', fontproperties=chinese_font)
plt.ylabel('实部', fontproperties=chinese_font)

# 添加从时间轴到曲线的垂线阴影
for i in shadow_indices:
    plt.plot([t[i], t[i]], [0, real_part[i]], 'k:', alpha=0.5)
    
# 标记特殊点
for st in special_times:
    idx = np.argmin(np.abs(t - st))
    plt.plot(t[idx], real_part[idx], 'ro')
    plt.annotate(f't = {st}', 
                xy=(t[idx], real_part[idx]),
                xytext=(t[idx]+0.2, real_part[idx]+0.1),
                fontproperties=chinese_font,
                fontsize=10)

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.font_manager import FontProperties
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.colors as colors
from matplotlib import cm

# 设置中文显示
font_path = '/System/Library/Fonts/Supplemental/Songti.ttc'  # macOS 中的宋体字体
chinese_font = FontProperties(fname=font_path)
plt.rcParams['axes.unicode_minus'] = False  # 解决负号显示问题

# 参数设置
T = 6  # 周期
omega = (2 * np.pi) / T  # 角频率
n = 2  # 周期个数
t = np.arange(0, n * T + 0.1, 0.1)  # 时间序列，从0到n*T，步长为0.1

# 计算复指数信号
f = np.exp(2j * omega * t)  # 注意这里使用2j
real_part = np.real(f)  # 实部
imag_part = np.imag(f)  # 虚部

# 创建三维图像
fig = plt.figure(figsize=(20, 16))
ax = fig.add_subplot(111, projection='3d')

# 绘制三维曲线
spiral, = ax.plot(t, imag_part, real_part, linewidth=3, color='blue', label='复指数信号')

# 参数设置
max_t = n * T
t_grid = np.linspace(0, max_t, 100)
xy_range = np.linspace(-1.5, 1.5, 100)
zero_grid = np.zeros_like(t_grid)

# 创建投影平面网格
t_mesh, im_mesh = np.meshgrid(t_grid, xy_range)
re_mesh = np.zeros_like(t_mesh)

# XY平面（地面投影）- 生成实际投影曲面
Z_xy = np.zeros_like(t_mesh)  # z=0平面
# 创建投影数据
t_proj = t.reshape(-1, 1)
imag_proj = imag_part.reshape(-1, 1)
points_xy = np.hstack((t_proj, imag_proj))

# 绘制XY平面基础
xy_plane = ax.plot_surface(t_mesh, im_mesh, re_mesh, 
                          color='lightblue', alpha=0.2)

# 在XY平面上绘制投影曲线
xy_proj, = ax.plot(t, imag_part, np.zeros_like(t), 'g-', linewidth=3)

# 绘制从曲线到XY平面的垂直线（"阴影线"）
for i in range(0, len(t), 10):
    ax.plot([t[i], t[i]], [imag_part[i], imag_part[i]], [real_part[i], 0], 'k--', alpha=0.3)

# XZ平面（后墙投影）
t_mesh, re_mesh = np.meshgrid(t_grid, xy_range)
im_mesh = np.zeros_like(t_mesh)  # y=0平面

# 绘制XZ平面
xz_plane = ax.plot_surface(t_mesh, im_mesh, re_mesh,
                          color='lightgreen', alpha=0.2)

# 在XZ平面上绘制投影曲线 
xz_proj, = ax.plot(t, np.zeros_like(t), real_part, 'r-', linewidth=3)

# 绘制从曲线到XZ平面的垂直线
for i in range(0, len(t), 10):
    ax.plot([t[i], t[i]], [0, imag_part[i]], [real_part[i], real_part[i]], 'k--', alpha=0.3)

# YZ平面（侧墙投影）- 复平面
im_mesh, re_mesh = np.meshgrid(xy_range, xy_range)
t_mesh = np.zeros_like(im_mesh)  # x=0平面

# 绘制YZ平面
yz_plane = ax.plot_surface(t_mesh, im_mesh, re_mesh,
                          color='lightyellow', alpha=0.2)

# 在YZ平面上绘制单位圆 
theta = np.linspace(0, 2*np.pi, 100)
y_circle = np.sin(2*theta)  # 这里使用2*theta匹配f=exp(2j*omega*t)
z_circle = np.cos(2*theta)
t_circle = np.zeros_like(theta)
unit_circle, = ax.plot(t_circle, y_circle, z_circle, 'k-', linewidth=2)

# 在复平面(YZ平面)上创建时间刻度的投影点
# 这些点表示复指数在特定时间点的值在复平面上的投影
special_times = np.linspace(0, T, 13)  # 更多的点以显示轨迹
for st in special_times:
    idx = np.argmin(np.abs(t - st))
    ax.scatter([0], [imag_part[idx]], [real_part[idx]], color='purple', s=30)
    # 连接原点到复平面上的点
    ax.plot([0, 0], [0, imag_part[idx]], [0, real_part[idx]], 'purple', alpha=0.3)

# 填充扇区阴影 - 在复平面(YZ平面)上
# 选择几个代表性时间点
for st in [T/4, T/2, 3*T/4, T]:
    idx = np.argmin(np.abs(t - st))
    # 计算扇区范围内的角度
    angle = np.linspace(0, 2*omega*t[idx], 50)  # 注意这里使用2*omega匹配f=exp(2j*omega*t)
    # 创建扇区点
    x_fan = np.zeros_like(angle)
    y_fan = np.sin(angle)
    z_fan = np.cos(angle)
    # 填充扇区
    fan_x = np.vstack((x_fan, x_fan))
    fan_y = np.vstack((np.zeros_like(angle), y_fan))
    fan_z = np.vstack((np.zeros_like(angle), z_fan))
    ax.plot_surface(fan_x, fan_y, fan_z, color='yellow', alpha=0.1)
    # 绘制扇区边界
    ax.plot(x_fan, y_fan, z_fan, 'k:', alpha=0.5)

# 标记特殊点
highlight_times = [0, T/4, T/2, 3*T/4, T]
for st in highlight_times:
    idx = np.argmin(np.abs(t - st))
    # 在螺旋上标记点
    ax.scatter([t[idx]], [imag_part[idx]], [real_part[idx]], color='red', s=100)
    ax.text(t[idx]+0.2, imag_part[idx], real_part[idx], f't={st}', fontsize=12)
    
    # 在xy平面投影上标记点
    ax.scatter([t[idx]], [imag_part[idx]], [0], color='green', s=80)
    
    # 在xz平面投影上标记点
    ax.scatter([t[idx]], [0], [real_part[idx]], color='red', s=80)
    
    # 在yz平面复平面上标记点
    ax.scatter([0], [imag_part[idx]], [real_part[idx]], color='purple', s=80)

# 添加时间轴上的原点标记
ax.scatter([0], [0], [0], color='black', s=100)
ax.text(0.2, 0, 0, '原点', fontsize=12, fontproperties=chinese_font)

# 设置坐标轴标签
ax.set_xlabel('x轴：时间 t', fontproperties=chinese_font, fontsize=14)
ax.set_ylabel('y 轴：虚部 Im(f) sin(2*omega*t)', fontproperties=chinese_font, fontsize=14)
ax.set_zlabel('z 轴：实部 Re(f) cos(2*omega*t)', fontproperties=chinese_font, fontsize=14)

# 设置坐标轴范围
ax.set_xlim(0, max_t)
ax.set_ylim(-1.5, 1.5)
ax.set_zlim(-1.5, 1.5)

# 添加网格
ax.grid(True)

# 设置标题 - 注意这里是2j而不是j
ax.set_title(f'复指数信号 $f=e^{{2i\\omega t}}$ 的三维可视化，周期T={T/2}，n={n}个周期',
            fontproperties=chinese_font, fontsize=16)

# 添加数学公式说明
formula_text = (
    f"复指数信号: $f = e^{{2i\\omega t}}$，其中 $\\omega = \\frac{{2\\pi}}{{T}} = \\frac{{2\\pi}}{{{T}}}$\n"
    f"实部: $\\cos(2\\omega t)$，虚部: $\\sin(2\\omega t)$\n"
    f"实际周期: T/{2} = {T/2}，总时间: {n*T}"
)
fig.text(0.5, 0.02, formula_text, ha='center', fontproperties=chinese_font, fontsize=14)

# 设置更好的视角
ax.view_init(elev=25, azim=40)

# 调整图形布局
plt.tight_layout(rect=[0, 0.05, 1, 0.95])
plt.show()

# 绘制平面投影的详细视图
plt.figure(figsize=(18, 12))

# 1. 时间域信号
plt.subplot(2, 2, 1)
plt.plot(t, real_part, 'b-', label=r'$\cos(2\omega t)$ (实部)')
plt.plot(t, imag_part, 'r-', label=r'$\sin(2\omega t)$ (虚部)')
plt.grid(True)
plt.legend(prop=chinese_font)
plt.title('复指数信号的实部和虚部', fontproperties=chinese_font)
plt.xlabel('时间 t', fontproperties=chinese_font)
plt.ylabel('幅度', fontproperties=chinese_font)

# 2. 复平面上的轨迹
plt.subplot(2, 2, 2)
plt.plot(real_part, imag_part, 'g-')
plt.grid(True)
plt.title('复平面轨迹 (单位圆)', fontproperties=chinese_font)
plt.xlabel('实部', fontproperties=chinese_font)
plt.ylabel('虚部', fontproperties=chinese_font)
plt.axis('equal')  # 保持坐标轴比例一致

# 添加阴影扇区
for st in highlight_times:
    idx = np.argmin(np.abs(t - st))
    angle = np.linspace(0, 2*omega * t[idx], 100)  # 注意这里用2*omega
    x_fan = np.cos(angle)
    y_fan = np.sin(angle)
    plt.fill(x_fan, y_fan, 'yellow', alpha=0.1)
    plt.plot(x_fan, y_fan, 'k--', alpha=0.3)
    
# 标记特殊点
for st in highlight_times:
    idx = np.argmin(np.abs(t - st))
    plt.plot(real_part[idx], imag_part[idx], 'ro')
    plt.annotate(f't = {st}', 
                xy=(real_part[idx], imag_part[idx]),
                xytext=(real_part[idx]*1.2, imag_part[idx]*1.2),
                fontproperties=chinese_font,
                fontsize=10)

# 3. 时间-虚部平面
plt.subplot(2, 2, 3)
plt.plot(t, imag_part, 'g-')
plt.grid(True)
plt.title('三维螺旋的俯视图 (t-虚部 平面)', fontproperties=chinese_font)
plt.xlabel('时间 t', fontproperties=chinese_font)
plt.ylabel('虚部', fontproperties=chinese_font)

# 添加从时间轴到曲线的垂线阴影
shadow_indices = np.arange(0, len(t), 20)
for i in shadow_indices:
    plt.plot([t[i], t[i]], [0, imag_part[i]], 'k:', alpha=0.5)
    
# 标记特殊点
for st in highlight_times:
    idx = np.argmin(np.abs(t - st))
    plt.plot(t[idx], imag_part[idx], 'ro')
    plt.annotate(f't = {st}', 
                xy=(t[idx], imag_part[idx]),
                xytext=(t[idx]+0.2, imag_part[idx]+0.1),
                fontproperties=chinese_font,
                fontsize=10)

# 4. 时间-实部平面
plt.subplot(2, 2, 4)
plt.plot(t, real_part, 'r-')
plt.grid(True)
plt.title('三维螺旋的侧视图 (t-实部 平面)', fontproperties=chinese_font)
plt.xlabel('时间 t', fontproperties=chinese_font)
plt.ylabel('实部', fontproperties=chinese_font)

# 添加从时间轴到曲线的垂线阴影
for i in shadow_indices:
    plt.plot([t[i], t[i]], [0, real_part[i]], 'k:', alpha=0.5)
    
# 标记特殊点
for st in highlight_times:
    idx = np.argmin(np.abs(t - st))
    plt.plot(t[idx], real_part[idx], 'ro')
    plt.annotate(f't = {st}', 
                xy=(t[idx], real_part[idx]),
                xytext=(t[idx]+0.2, real_part[idx]+0.1),
                fontproperties=chinese_font,
                fontsize=10)

plt.tight_layout()
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.font_manager import FontProperties
import matplotlib as mpl

# 设置中文显示
font_path = '/System/Library/Fonts/Supplemental/Songti.ttc'  # macOS 中的宋体字体
chinese_font = FontProperties(fname=font_path)
plt.rcParams['axes.unicode_minus'] = False  # 解决负号显示问题

# 创建数据
theta = np.linspace(0, 2*np.pi, 1000)
x = np.cos(theta)
y = np.sin(theta)

# 复平面上的特殊点
special_points = {
    0: [1, 0, r'$e^{i \cdot 0} = 1$', 'red'],
    np.pi/2: [0, 1, r'$e^{i \cdot \pi/2} = i$', 'blue'],
    np.pi: [-1, 0, r'$e^{i \cdot \pi} = -1$', 'green'],
    3*np.pi/2: [0, -1, r'$e^{i \cdot 3\pi/2} = -i$', 'purple'],
    2*np.pi: [1, 0, r'$e^{i \cdot 2\pi} = 1$', 'orange']
}

# 创建图形
plt.figure(figsize=(10, 10))
ax = plt.subplot(111)

# 绘制单位圆
plt.plot(x, y, 'k-', lw=1.5, alpha=0.7)

# 绘制坐标轴
plt.axhline(y=0, color='gray', linestyle='-', alpha=0.7)
plt.axvline(x=0, color='gray', linestyle='-', alpha=0.7)

# 绘制单位圆上的特殊点
for angle, (px, py, text, color) in special_points.items():
    plt.plot(px, py, 'o', markersize=10, color=color)
    plt.annotate(text, xy=(px, py), xytext=(px*1.5, py*1.5),
                fontsize=12, color=color,
                arrowprops=dict(facecolor=color, shrink=0.05, width=2, alpha=0.7))

# 为某个任意点绘制复指数可视化
t_example = np.pi/4  # 可以更改这个值看不同角度
x_example = np.cos(t_example)
y_example = np.sin(t_example)

# 绘制到任意点的半径
plt.plot([0, x_example], [0, y_example], 'r--', lw=2)

# 标记角度
arc = np.linspace(0, t_example, 100)
arc_x = 0.2 * np.cos(arc)
arc_y = 0.2 * np.sin(arc)
plt.plot(arc_x, arc_y, 'b-', lw=1.5)
plt.text(0.25 * np.cos(t_example/2), 0.25 * np.sin(t_example/2), 
         r'$\theta$', fontsize=14, ha='center', va='center')

# 添加坐标标签
plt.text(1.05, 0, r'$Re$', fontsize=14)
plt.text(0, 1.05, r'$Im$', fontsize=14)

# 设置绘图范围和标题
plt.grid(True, linestyle='--', alpha=0.7)
plt.axis('equal')
plt.xlim(-2, 2)
plt.ylim(-2, 2)
plt.title('欧拉公式 $e^{i\\theta} = \\cos\\theta + i\\sin\\theta$', fontproperties=chinese_font, fontsize=16)

# 添加公式注释
formula_text = (
    r"欧拉公式: $e^{i\theta} = \cos\theta + i\sin\theta$"
    "\n"
    r"当 $\theta = \pi$ 时: $e^{i\pi} + 1 = 0$"
)
plt.figtext(0.5, 0.01, formula_text, ha='center', fontproperties=chinese_font, fontsize=14)

plt.tight_layout(rect=[0, 0.05, 1, 0.95])
plt.show()

# 另外绘制欧拉公式分量图
t = np.linspace(0, 2*np.pi, 1000)
cos_t = np.cos(t)
sin_t = np.sin(t)

plt.figure(figsize=(12, 8))

# 实部和虚部作为函数图像
plt.subplot(211)
plt.plot(t, cos_t, 'b-', label=r'$\cos(\theta)$ (实部)')
plt.plot(t, sin_t, 'r-', label=r'$\sin(\theta)$ (虚部)')
plt.grid(True)
plt.legend(prop=chinese_font)
plt.title('欧拉公式的实部和虚部', fontproperties=chinese_font)
plt.xlabel(r'$\theta$')
plt.ylabel(r'$\cos(\theta), \sin(\theta)$')

# 参数方程表示
plt.subplot(212)
plt.plot(cos_t, sin_t, 'g-')
plt.grid(True)
plt.title('欧拉公式在复平面上的轨迹(单位圆)', fontproperties=chinese_font)
plt.xlabel(r'$\cos(\theta)$ (实部)')
plt.ylabel(r'$\sin(\theta)$ (虚部)')
plt.axis('equal')

plt.tight_layout()
plt.show()