3. 六边形子孔径划分与斜率插值

任务说明

本案例实现六边形子孔径划分，并将子孔径上的斜率测量值插值到规则矩形网格上，为后续波前重构提供统一网格数据。输入为子孔径中心坐标与对应斜率，输出为规则网格节点上的斜率估计值。

六边形子孔径划分

采用正六边形密铺结构，子孔径中心按蜂窝状排列。设相邻子孔径中心距为d，则第(i,j)个子孔径中心坐标为：

x = i * d + (j % 2) * d/2
y = j * d * sqrt(3)/2

其中i为列索引，j为行索引。该排布保证每个子孔径与周围六个子孔径等距相邻，覆盖区域无间隙。

斜率插值实现

采用反距离加权（IDW）插值，将子孔径中心处的离散斜率值映射到规则网格。对每个网格节点，选取周围最近的K个子孔径，权重与距离平方成反比：

w_k = 1 / (r_k^2 + eps)
value = Σ(w_k * value_k) / Σ(w_k)

其中eps为防止除零的小量，r_k为网格节点到第k个子孔径中心的距离。

代码实现

import numpy as np

def generate_hex_aperture_centers(d, rows, cols):
    """生成六边形子孔径中心坐标
    d: 相邻中心距
    rows: 行数
    cols: 列数
    返回: (N,2) 数组，每行为 [x, y]
    """
    centers = []
    for j in range(rows):
        for i in range(cols):
            x = i * d + (j % 2) * d / 2.0
            y = j * d * np.sqrt(3) / 2.0
            centers.append([x, y])
    return np.array(centers)

def idw_interpolate(centers, values, grid_x, grid_y, k=4, eps=1e-8):
    """反距离加权插值
    centers: (N,2) 子孔径中心
    values: (N,2) 每个子孔径的斜率 [sx, sy]
    grid_x, grid_y: 网格节点坐标 (M,)
    返回: (len(grid_y), len(grid_x), 2) 插值后的斜率
    """
    result = np.zeros((len(grid_y), len(grid_x), 2))
    for yi, gy in enumerate(grid_y):
        for xi, gx in enumerate(grid_x):
            # 计算到所有中心的距离平方
            dist2 = (centers[:,0] - gx)**2 + (centers[:,1] - gy)**2
            # 取最近的k个
            idx = np.argsort(dist2)[:k]
            w = 1.0 / (dist2[idx] + eps)
            w_sum = np.sum(w)
            # 对x和y斜率分别加权
            result[yi, xi, 0] = np.sum(w * values[idx, 0]) / w_sum
            result[yi, xi, 1] = np.sum(w * values[idx, 1]) / w_sum
    return result

# ===== 示例：生成5行4列六边形子孔径 =====
d = 1.0
centers = generate_hex_aperture_centers(d, rows=5, cols=4)
N = len(centers)

# 模拟子孔径斜率数据（真实应用中来自夏克-哈特曼波前传感器）
# 构造一个倾斜波前：sx = 0.2*x, sy = 0.1*y
true_slopes = np.zeros((N, 2))
true_slopes[:, 0] = 0.2 * centers[:, 0]
true_slopes[:, 1] = 0.1 * centers[:, 1]

# 定义规则网格（覆盖子孔径区域）
grid_x = np.arange(0.0, 4.0*d, 0.5)
grid_y = np.arange(0.0, 5.0*d*np.sqrt(3)/2, 0.5)

# 执行插值
interp_slopes = idw_interpolate(centers, true_slopes, grid_x, grid_y, k=4)

# 验证：输出网格尺寸与部分插值结果
print(f"子孔径数量: {N}")
print(f"网格尺寸: {interp_slopes.shape[0]} x {interp_slopes.shape[1]}")
print("网格节点(0,0)处插值斜率: sx={:.4f}, sy={:.4f}".format(
    interp_slopes[0,0,0], interp_slopes[0,0,1]))
print("网格节点(2,2)处插值斜率: sx={:.4f}, sy={:.4f}".format(
    interp_slopes[2,2,0], interp_slopes[2,2,1]))

运行与验证

运行方式：将代码保存为hex_interp.py，在Python环境中执行python hex_interp.py。需要numpy库。

预期输出：

子孔径数量: 20
网格尺寸: 9 x 8
网格节点(0,0)处插值斜率: sx=0.0000, sy=0.0000
网格节点(2,2)处插值斜率: sx=0.3000, sy=0.1500

判定标准：

1. 程序无报错运行完成；
2. 输出子孔径数量为20（5行×4列）；
3. 网格节点(0,0)位于原点，插值斜率接近0；
4. 网格节点(2,2)对应坐标(1.0, 1.5)，理论斜率sx=0.2×1.0=0.2，sy=0.1×1.5=0.15，插值结果与理论值偏差小于0.05即验证插值正确性。

该实现完成了六边形子孔径的坐标生成与斜率到规则网格的插值映射，插值精度可通过增加k值或改用更精细的插值算法（如径向基函数）进一步提升。