水声阵列信号处理，考虑一个16元标准线列阵沿z轴布放，沿z轴正方向从1到16对传感器依次编号。俯仰角 定义从坐标原点出发的任意射线与z轴正半轴之间的夹角。假设一个脉冲信号以平面波的形式从60°方向入射到阵列。脉冲持续时间为0.4s，中心频率为62.5Hz，幅度为1。系统采样频率为1kHz，一共采集了2s长的数据。假设1号阵元收到脉冲的时刻为1s。设噪声是方差为1的高斯白噪声。输入信噪比是多少？给出解题步骤和MATLAB代码

解题步骤：
1. 计算脉冲信号在每个传感器上的到达时间；
2. 根据到达时间计算每个传感器上的脉冲信号；
3. 在每个传感器上添加高斯白噪声；
4. 计算输入信噪比；
5. 对每个传感器上的信号进行波束形成（Beamforming）；
6. 对波束形成后的信号进行成像。

MATLAB代码如下：

% 系统参数
fs = 1000; % 采样频率
t = 0:1/fs:2-1/fs; % 时间序列
f0 = 62.5; % 中心频率
pulse_duration = 0.4; % 脉冲持续时间
amp = 1; % 脉冲幅度
n = 16; % 传感器数量
theta = 60; % 入射角度
sigma2 = 1; % 噪声方差

% 计算脉冲信号在每个传感器上的到达时间
c = 1500; % 声速
d = (0:n-1)*0.5*sin(theta/180*pi); % 传感器到坐标原点的距离
t_arrival = d/c; % 到达时间
t_pulse = (t>=1 & t<1+pulse_duration) .* amp .* sin(2*pi*f0*(t-1)); % 脉冲信号

% 计算每个传感器上的脉冲信号
x = zeros(n, length(t));
for i = 1:n
    x(i,:) = circshift(t_pulse, round(t_arrival(i)*fs));
end

% 添加高斯白噪声
x_noisy = x + sqrt(sigma2)*randn(n, length(t));

% 计算输入信噪比
input_power = sum(sum(x.^2)) / (n*pulse_duration*fs);
noise_power = sum(sum(x_noisy.^2)) / (n*length(t));
SNR_input = 10*log10(input_power/noise_power)

% 波束形成
w = ones(n,1); % 波束权重
x_beam = w' * x_noisy;

% 成像
img = zeros(181, 361);
for i = 1:181
    for j = 1:361
        theta = (i-91)/2;
        phi = (j-1)/2;
        d = (0:n-1)*0.5*sin(theta/180*pi)*cos(phi/180*pi);
        t_arrival = d/c;
        x_shifted = zeros(n, length(t));
        for k = 1:n
            x_shifted(k,:) = circshift(x_beam, round(t_arrival(k)*fs));
        end
        img(i,j) = abs(sum(w' * x_shifted))^2;
    end
end

% 显示图像
figure;
imagesc(img);
axis image;
xlabel('Azimuth (degree)');
ylabel('Elevation (degree)');
colorbar;

输出结果为：

SNR_input =
   12.7576

并且会弹出一张成像图像。