import numpy as np
import matplotlib.pyplot as plt
from matplotlib import colors
import time
from mpl_toolkits.mplot3d import Axes3D
# 解决中文显示和负号显示问题
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False
# ====== 参数设置 ======
c = 3e8 # 光速 (m/s)
f0 = 2.4e9 # 中心频率 (Hz)
N = 32 # 阵元数量
d = 0.5 * c / f0 # 阵元间距 (m)
B = 1.5e6 # 频偏带宽 (Hz)
r0 = 2500 # 目标距离 (m)
theta0 = 0 # 目标角度 (度)
r_min, r_max = 0, 5000 # 距离扫描范围
theta_min, theta_max = -90, 90 # 角度扫描范围
r_points = 200 # 距离扫描点数
theta_points = 200 # 角度扫描点数
# ====== 非线性频偏函数 ======
def frequency_offsets(offset_type, N, B):
"""计算四种非线性频偏"""
n = np.arange(N)
n_centered = n - (N - 1) / 2 # 中心对称的阵元索引
if offset_type == 'sin':
return (B / 2) * np.sin(2 * np.pi * n / N)
elif offset_type == 'log':
abs_n = np.abs(n_centered) + 1e-6 # 避免log(0)
log_norm = np.log(1 + abs_n) / np.log(1 + N / 2)
return (B / 2) * log_norm * np.sign(n_centered)
elif offset_type == 'square':
return (B / 2) * (n_centered / (N / 2)) ** 2 * np.sign(n_centered)
elif offset_type == 'reciprocal':
a = 5 # 控制曲线陡峭度
abs_n = np.abs(n_centered) + 1e-6 # 避免除零
return (B / 2) * a / (a + abs_n) * np.sign(n_centered)
else:
raise ValueError("不支持的频偏类型。请选择: 'sin', 'log', 'square', 'reciprocal'")
# ====== 计算导向矢量 (向量化优化) ======
def steering_vector(r, theta, freqs, d, c):
"""向量化计算导向矢量"""
theta_rad = np.deg2rad(theta)
n = np.arange(len(freqs))
# 向量化计算相位项
phase_range = 2 * np.pi * freqs * r / c
phase_angle = 2 * np.pi * freqs * n * d * np.sin(theta_rad) / c
phases = np.exp(-1j * (phase_range - phase_angle))
return phases
# ====== 计算波束图 (修复形状问题) ======
def calculate_beam_pattern(offset_type, r0, theta0):
print(f"计算 {offset_type} 频偏的波束图...")
start_time = time.time()
# 计算频率偏移
delta_f = frequency_offsets(offset_type, N, B)
freqs = f0 + delta_f
# 目标位置导向矢量 (权向量)
w = steering_vector(r0, theta0, freqs, d, c)
w = w.reshape(-1, 1) # 关键修复：转换为列向量 (32, 1)
# 创建扫描网格
r_scan = np.linspace(r_min, r_max, r_points)
theta_scan = np.linspace(theta_min, theta_max, theta_points)
R, Theta = np.meshgrid(r_scan, theta_scan)
# 预分配内存
beam = np.zeros_like(R, dtype=np.float64)
# 向量化计算波束响应
for i, theta in enumerate(theta_scan):
# 计算当前角度下所有距离点的导向矢量
theta_rad = np.deg2rad(theta)
n = np.arange(N)
# 距离引起的相位延迟 (形状: (N, r_points))
phase_range = 2 * np.pi * freqs[:, None] * r_scan / c
# 角度引起的相位差 (形状: (N, 1))
phase_angle = 2 * np.pi * freqs * n * d * np.sin(theta_rad) / c
phase_angle = phase_angle[:, None] # 转换为列向量
# 计算导向矩阵 (形状: (N, r_points))
A = np.exp(-1j * (phase_range - phase_angle))
# 计算波束响应 (形状: (r_points,))
# 关键修复：确保w(32,1)和A(32,200)形状兼容
beam[i, :] = np.abs(np.sum(w.conj() * A, axis=0))
print(f"计算完成，耗时: {time.time() - start_time:.2f}秒")
return R, Theta, beam
# ====== 二维截面图可视化 ======
def plot_2d_section(R, Theta, beam, offset_type):
plt.figure(figsize=(12, 8))
# 归一化波束图
beam_norm = beam / np.max(beam)
# 创建颜色映射
norm = colors.Normalize(vmin=0, vmax=1)
# 绘制热力图
plt.pcolormesh(R, Theta, beam_norm, shading='auto', cmap='viridis', norm=norm)
cbar = plt.colorbar(label='归一化幅度')
# 添加等高线
levels = np.linspace(0, 1, 11)
CS = plt.contour(R, Theta, beam_norm, levels=levels, colors='w', linewidths=0.8)
plt.clabel(CS, inline=True, fontsize=8, fmt='%1.1f')
# 添加-3dB分界线
db3_level = 0.707 # -3dB对应的幅度值
db3_contour = plt.contour(
R, Theta, beam_norm,
levels=[db3_level],
colors='red',
linewidths=2.0
)
plt.clabel(db3_contour, inline=True, fontsize=10, fmt=f'-3dB = {db3_level:.3f}')
# 添加目标位置标记
#plt.plot(r0, theta0, 'ro', markersize=8, label='目标位置')
plt.title(f'频控阵波束图 - {offset_type.capitalize()}频偏', fontsize=14)
plt.xlabel('距离 (m)', fontsize=12)
plt.ylabel('角度 (度)', fontsize=12)
plt.legend()
plt.grid(True, linestyle='--', alpha=0.6)
plt.tight_layout()
plt.savefig(f'FDA_2D_Section_{offset_type}.png', dpi=300)
plt.show()
# ====== 3D波束图可视化 ======
def plot_3d_beam(R, Theta, beam, offset_type):
print(f"生成 {offset_type} 频偏的3D波束图...")
start_time = time.time()
# 转换为直角坐标
X = R * np.sin(np.deg2rad(Theta))
Y = R * np.cos(np.deg2rad(Theta))
Z = beam / np.max(beam) # 归一化
fig = plt.figure(figsize=(14, 10))
ax = fig.add_subplot(111, projection='3d')
# 创建3D表面图 (优化性能)
surf = ax.plot_surface(
X, Y, Z,
rstride=2,
cstride=2,
cmap='viridis',
edgecolor='none',
alpha=0.8
)
# 添加颜色条
cbar = fig.colorbar(surf, shrink=0.5, aspect=10)
cbar.set_label('归一化幅度')
# 标记目标位置
x0 = r0 * np.sin(np.deg2rad(theta0))
y0 = r0 * np.cos(np.deg2rad(theta0))
#ax.scatter(x0, y0, 1.05, color='red', s=100, label='目标位置')
ax.set_title(f'3D频控阵波束图 - {offset_type.capitalize()}频偏', fontsize=16)
ax.set_xlabel('X (m)', fontsize=12)
ax.set_ylabel('Y (m)', fontsize=12)
ax.set_zlabel('归一化幅度', fontsize=12)
ax.legend()
# 优化视角
ax.view_init(elev=30, azim=45)
plt.tight_layout()
plt.savefig(f'FDA_3D_Beam_{offset_type}.png', dpi=300)
print(f"3D图生成完成，耗时: {time.time() - start_time:.2f}秒")
plt.show()
# ====== 主程序 ======
if __name__ == "__main__":
offset_types = ['sin']
for offset_type in offset_types:
print(f"\n{'=' * 40}")
print(f"处理 {offset_type} 类型的频偏")
R, Theta, beam = calculate_beam_pattern(offset_type, r0, theta0)
plot_2d_section(R, Theta, beam, offset_type)
plot_3d_beam(R, Theta, beam, offset_type)去归一化
最新发布
本文详细介绍了通信系统中时偏的概念、来源及其对性能的影响，重点讨论了时偏估计和补偿的各种方法，包括高精度时钟同步、信号处理优化、自适应技术等，并强调了在无线和光纤通信中的应用及未来研究挑战。
扫一扫
1万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
