帮我写一个模拟不同密度的靶的卢瑟福散射并统计数据的matlab代码

好的，以下是一个简单的模拟代码，可以模拟不同密度的靶的卢瑟福散射，并统计数据：

% 模拟卢瑟福散射并统计数据
% 设置模拟参数
num_particles = 1000; % 粒子数量
num_targets = 3; % 靶的数量
target_density = [0.1, 0.5, 1]; % 靶的密度
target_radius = 1e-15; % 靶的半径
incident_angle = pi/2; % 入射角度
min_angle = 0; % 统计最小角度
max_angle = pi; % 统计最大角度
num_bins = 100; % 统计数据的 bin 数量

% 初始化粒子位置和速度
particle_pos = zeros(num_particles, 2);
particle_vel = zeros(num_particles, 2);
particle_pos(:, 1) = -10*target_radius; % 所有粒子都从同一位置出发
particle_vel(:, 1) = randn(num_particles, 1); % 随机速度

% 初始化统计数据
angle_bins = linspace(min_angle, max_angle, num_bins+1);
angle_counts = zeros(num_targets, num_bins);

% 开始模拟
for i = 1:num_particles
    for j = 1:num_targets
        % 计算粒子与靶的距离
        dist = particle_pos(i, 1) - target_radius*(2*j - 1 - num_targets);
        if dist > target_radius
            continue; % 粒子还没到达靶前就越过了
        end
        
        % 计算粒子与靶的散射角度
        theta = asin(dist/sqrt(dist^2 + particle_pos(i, 2)^2));
        
        % 统计数据
        bin_idx = ceil((theta - min_angle)/(max_angle - min_angle)*num_bins);
        angle_counts(j, bin_idx) = angle_counts(j, bin_idx) + 1;
        
        % 更新粒子速度和位置
        particle_vel(i, :) = [cos(theta), sin(theta); -sin(theta), cos(theta)]*particle_vel(i, :)';
        particle_pos(i, :) = particle_pos(i, :) + particle_vel(i, :)*target_radius;
    end
end

% 绘制统计图
figure;
for j = 1:num_targets
    subplot(num_targets, 1, j);
    plot(angle_bins, angle_counts(j, :), '-');
    title(sprintf('Target %d density = %.1f', j, target_density(j)));
    xlabel('Scattering angle (rad)');
    ylabel('Counts');
end

这段代码模拟了 `num_particles` 个粒子从左侧进入 `num_targets` 个密度不同的靶中，并记录了在 `min_angle` 到 `max_angle` 范围内每个 bin 的粒子数。你可以根据需要修改代码中的参数来进行模拟，并根据统计数据绘制散射图。