%% 光谱法测量大深径微孔深度计算程序
% 该程序通过光谱反射法测量微孔深度
% 使用傅里叶变换光谱分析技术计算微孔深度
% 数据来源于外部文件输入
%% 清除环境并设置参数
clear all;
close all;
clc;
%% 参数设置
refractive_index = 1; % 微孔内介质折射率
filename = 'C:\Users\m\Desktop\新建文件夹\USB2H090911_11-28-50-828.txt'; % 输入文件名
%% 从文件读取光谱数据
try
data = readtable(filename);
% 检查数据格式
if height(data) < 100
error('数据点过少，请检查输入文件');
end
% 假设第一列是波长，第二列是光强
if width(data) >= 2
lambda = data{:,1} * 1e-9; % 转换为米
I_measured = data{:,2};
else
error('输入文件至少需要包含波长和光强两列数据');
end
catch ME
error(['读取文件时出错: ', ME.message]);
end
%% 光谱预处理
% 数据检查和处理
if any(diff(lambda) <= 0)
% 确保波长是严格递增的
[lambda, idx] = sort(lambda);
I_measured = I_measured(idx);
end
% 基线校正 - 使用二次多项式拟合
x = lambda - mean(lambda); % 中心化
x_scaled = x / std(x); % 缩放
p = polyfit(x_scaled, I_measured, 1); % 使用二次多项式
baseline = polyval(p, x_scaled);
I_corrected = I_measured - baseline + 0.5; % 调整到合理范围
% 归一化
I_normalized = (I_corrected - min(I_corrected)) / (max(I_corrected) - min(I_corrected));
%% 傅里叶变换光谱分析
% 转换到波数空间 (单位: m^-1)
wave_number = 1 ./ lambda;
wave_number = wave_number - min(wave_number); % 从0开始
% 插值到均匀波数网格
wave_number_uniform = linspace(min(wave_number), max(wave_number), length(lambda));
I_uniform = interp1(wave_number, I_normalized, wave_number_uniform, 'spline');
% 傅里叶变换
FT_data = fftshift(fft(I_uniform));
FT_magnitude = abs(FT_data);
% 创建距离轴 (单位: m)
sample_spacing = (max(wave_number_uniform) - min(wave_number_uniform)) / (length(wave_number_uniform) - 1);
distance_axis = (-length(FT_magnitude)/2 : length(FT_magnitude)/2-1) * (1 / (sample_spacing * length(FT_magnitude)));
%% 深度计算
% 寻找主峰位置 (排除零频附近区域)
n = length(FT_magnitude);
center_idx = floor(n/2) + 1; % 零频位置
search_range = max(10, floor(n/20)); % 搜索范围，至少10个点
% 创建排除零频附近的搜索向量
left_part = FT_magnitude(1:center_idx-search_range);
right_part = FT_magnitude(center_idx+search_range:end);
search_vector = (left_part:right_part); % 修正拼接符号
% 找到最大值及其索引
[~, peak_idx] = max(search_vector);
% 转换为全局索引
if peak_idx <= length(left_part)
peak_idx = peak_idx;
else
peak_idx = peak_idx - length(left_part) + center_idx + search_range;
end
% 计算深度 (单位: m)
%% 深度计算
% 寻找主峰位置 (排除零频附近区域)
n = length(FT_magnitude);
center_idx = floor(n/2) + 1; % 零频位置
search_range = max(10, floor(n/20)); % 搜索范围，至少10个点
% 只考虑正频率部分
positive_side = center_idx+search_range:n;
[~, max_idx] = max(FT_magnitude(positive_side));
peak_idx = positive_side(max_idx);
% 计算深度 (单位: m)
% 深度 = 光程差/2 = (距离轴值 * 折射率)/2
calculated_depth = abs(distance_axis(peak_idx)) * refractive_index / 2;
% 确保深度为正数
calculated_depth = abs(calculated_depth);
%% 结果显示
figure('Position', [100, 100, 1200, 800]);
% 绘制原始光谱
subplot(2, 2, 1);
plot(lambda*1e9, I_measured, 'r-', 'LineWidth', 1.5);
hold on;
plot(lambda*1e9, baseline, 'b--', 'LineWidth', 1);
xlabel('波长 (nm)');
ylabel('反射强度');
title('测量光谱与基线');
legend('测量数据', '拟合基线');
grid on;
% 绘制预处理后的光谱
subplot(2, 2, 2);
plot(lambda*1e9, I_normalized, 'g-', 'LineWidth', 1.5);
xlabel('波长 (nm)');
ylabel('归一化强度');
title('预处理后的光谱');
grid on;
% 绘制傅里叶变换结果
subplot(2, 2, 3);
plot(distance_axis*1e6, FT_magnitude, 'k-', 'LineWidth', 1.5);
xlabel('光程差 (\mum)');
ylabel('幅度');
title('傅里叶变换光谱');
grid on;
xlim([-100 100]); % 聚焦感兴趣区域
% 在FT图上标记峰值位置
hold on;
plot(distance_axis(peak_idx)*1e6, FT_magnitude(peak_idx), 'ro', 'MarkerSize', 8, 'MarkerFaceColor', 'r');
text(distance_axis(peak_idx)*1e6, FT_magnitude(peak_idx), [' 峰值位置: ' num2str(distance_axis(peak_idx)*1e6, '%.2f') ' \mum'], 'FontSize', 8);
hold off;
% 绘制结果汇总
subplot(2, 2, 4);
% 绘制微孔示意图
rectangle('Position', [0.3, 0.1, 0.1, 0.8], 'FaceColor', [0.8 0.8 0.8]);
line([0.3, 0.4], [0.9, 0.9], 'Color', 'k', 'LineWidth', 2);
line([0.3, 0.3], [0.1, 0.9], 'Color', 'k', 'LineWidth', 2);
line([0.4, 0.4], [0.1, 0.9], 'Color', 'k', 'LineWidth', 2);
text(0.25, 0.45, ['计算深度: ' num2str(calculated_depth*1e6, '%.2f') ' \mum'], 'FontSize', 10);
axis off;
title('测量结果汇总');
%% 输出结果
fprintf('计算微孔深度: %.2f \mu m \n', calculated_depth*1e6);修改代码使计算结果正确
最新发布
在使用MATLAB绘制以日期为x轴的数据图时，遇到了xlim设置范围的错误。文章详细介绍了错误信息及原因，提供了正确的xlim设置方法，即使用datetime类型变量的索引。
扫一扫
5614
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
