实现多种光栅的反射谱与传输谱MATLAB仿真

一、光纤布拉格光栅（FBG）：耦合模理论（CMT）+ 传输矩阵法（TMM）

FBG的反射/传输谱由纤芯模与反向纤芯模的耦合决定，核心是布拉格条件（ $\lambda_B = 2n_{\text{eff}}\Lambda$ ）。MATLAB中常用传输矩阵法（TMM）或耦合模理论（CMT）仿真。

1. 理论模型

传输矩阵法（TMM）：将FBG分割为多个均匀小段，每段用2×2传输矩阵表示，通过矩阵连乘得到整体反射/传输系数。
均匀FBG的传输矩阵为：

其中， $\gamma = \sqrt{\kappa^2 - \delta^2}$ （ $\kappa$ 为耦合系数， $\delta = \beta - \pi/\Lambda$ 为失配因子）， $L$ 为光栅长度。
耦合模理论（CMT）：直接求解前向（ $A(z)$ ）与后向（ $B(z)$ ）光场的演化方程，得到反射率公式：

2. MATLAB实现案例

以下是均匀FBG的反射谱与传输谱仿真代码（基于TMM）：

% 参数设置
n_eff = 1.45;          % 有效折射率
Lambda = 530e-9;       % 光栅周期（m）
L = 10e-3;             % 光栅长度（m）
delta_n = 1e-4;        % 折射率调制深度
lambda = linspace(1500e-9, 1600e-9, 1000);  % 波长范围（m）

% 计算耦合系数与失配因子
kappa = pi * delta_n ./ lambda;  % 耦合系数
delta = 2*pi*n_eff .* (1./lambda - 1/(2*Lambda));  % 失配因子

% 传输矩阵法计算反射率与透射率
gamma = sqrt(kappa.^2 - delta.^2);
cosh_gammaL = cosh(gamma*L);
sinh_gammaL = sinh(gamma*L);

R = (sinh_gammaL.^2) ./ (cosh_gammaL.^2 - (delta./kappa).^2);  % 反射率
T = 1 ./ (cosh_gammaL.^2 - (delta./kappa).^2);  % 透射率（忽略损耗）

% 可视化
figure;
subplot(2,1,1);
plot(lambda*1e9, R, 'r', 'LineWidth', 1.5);
xlabel('波长（nm）');
ylabel('反射率');
title('均匀FBG反射谱');
grid on;

subplot(2,1,2);
plot(lambda*1e9, T, 'b', 'LineWidth', 1.5);
xlabel('波长（nm）');
ylabel('透射率');
title('均匀FBG传输谱');
grid on;

3. 扩展：啁啾FBG与非均匀FBG

对于啁啾FBG（周期随位置变化），可将光栅分割为多个小段，每段用不同的 $\Lambda_i$ 计算传输矩阵，再连乘得到整体响应。例如：

% 啁啾FBG参数
Lambda_0 = 530e-9;     % 初始周期（m）
chirp_rate = 0.1e-9/m;  % 啁啾率（周期随长度变化率）
L = 10e-3;             % 光栅长度（m）
z = linspace(0, L, 100);% 位置向量
Lambda = Lambda_0 + chirp_rate * z;  % 周期分布

% 分割为100个小段，每段长度dz = L/100
dz = L/100;
M_total = eye(2);  % 初始化总传输矩阵
for i = 1:100
    % 第i段的周期与参数
    Lambda_i = Lambda(i);
    kappa_i = pi * delta_n / lambda;  % 耦合系数（假设波长固定）
    delta_i = 2*pi*n_eff * (1/lambda - 1/(2*Lambda_i));  % 失配因子
    gamma_i = sqrt(kappa_i^2 - delta_i^2);
    
    % 第i段的传输矩阵
    M_i = [cosh(gamma_i*dz) - 1i*delta_i/gamma_i*sinh(gamma_i*dz), -1i*kappa_i/gamma_i*sinh(gamma_i*dz);
          1i*kappa_i/gamma_i*sinh(gamma_i*dz), cosh(gamma_i*dz) + 1i*delta_i/gamma_i*sinh(gamma_i*dz)];
    
    % 总传输矩阵连乘
    M_total = M_i * M_total;
end

% 计算反射率与透射率
r = M_total(2,1)/M_total(1,1);
R = abs(r)^2;
t = 1/M_total(1,1);
T = abs(t)^2;

二、平面透射/反射光栅：傅里叶模态法（FMM）+ 夫琅禾费衍射

平面光栅的反射/传输谱由二维周期结构的衍射决定，核心是夫琅禾费衍射（远场衍射）。MATLAB中常用傅里叶变换或傅里叶模态法（FMM）仿真。

1. 理论模型

夫琅禾费衍射：对于振幅型平面光栅（透过率调制），远场衍射光强分布为：

其中， $N$ 为缝数， $d$ 为光栅周期， $a$ 为缝宽， $\theta$ 为衍射角。
傅里叶模态法（FMM）：将光栅的折射率分布与电磁场用傅里叶级数展开，求解麦克斯韦方程组，适用于亚波长光栅（周期小于波长）。

2. MATLAB实现案例

以下是振幅型平面光栅的远场衍射光强分布仿真代码（基于夫琅禾费衍射）：

% 参数设置
lambda = 632.8e-9;     % 波长（He-Ne激光，m）
d = 1e-6;              % 光栅周期（m）
a = 0.5e-6;            % 缝宽（m）
N = 100;               % 缝数
theta = linspace(-pi/6, pi/6, 1000);  % 衍射角范围（rad）

% 计算光强分布
I0 = 1;  % 入射光强
term1 = (sin(N*pi*d*sin(theta)/lambda) ./ (N*sin(pi*d*sin(theta)/lambda))).^2;  % 多缝干涉项
term2 = (sin(pi*a*sin(theta)/lambda) ./ sin(pi*d*sin(theta)/lambda)).^2;          % 单缝衍射项
I = I0 * term1 .* term2;

% 转换为角度（deg）
theta_deg = theta * 180/pi;

% 可视化
figure;
plot(theta_deg, I, 'b', 'LineWidth', 1.5);
xlabel('衍射角（deg）');
ylabel('相对光强');
title('振幅型平面光栅远场衍射光强分布');
grid on;

3. 扩展：相位型平面光栅

对于相位型平面光栅（相位调制），需将透过率函数替换为相位函数（如 $t(x) = e^{i\phi(x)}$ ），再通过傅里叶变换计算远场衍射光强。例如，矩形相位光栅（相位差为 $\phi_0$ ）的透过率函数为：

其远场衍射光强为：

三、体全息光栅（VHG）：严格耦合波理论（RCWA）+ Kogelnik理论

体全息光栅的衍射特性由三维体积内的折射率调制决定，核心是布拉格衍射的相位匹配。MATLAB中常用严格耦合波理论（RCWA）或Kogelnik理论仿真。

1. 理论模型

Kogelnik理论：描述反射式VHG的衍射效率，考虑折射率调制度（ $\Delta n$ ）、光栅厚度（ $d$ ）、周期（ $\Lambda$ ）的影响：

其中， $\psi = \frac{\pi \Delta n d}{\lambda \cos\theta}$ ， $\xi = \frac{d\tau}{2\cos\theta}$ （ $\tau = K\cos(\phi-\theta) - 2\pi n/\lambda$ 为失配因子， $K=2\pi/\Lambda$ 为光栅矢量， $\phi$ 为光栅倾角）。
严格耦合波理论（RCWA）：将光栅区域内的电磁场用正交Legendre多项式展开，求解耦合波方程组，适用于复杂折射率分布的VHG（如倾斜VHG）。

2. MATLAB实现案例

以下是反射式体全息光栅的衍射效率仿真代码（基于Kogelnik理论）：

% 参数设置
lambda = 532e-9;       % 波长（m）
n = 1.5;               % 介质折射率
Delta_n = 0.01;        % 折射率调制度
d = 10e-6;             % 光栅厚度（m）
Lambda = 1e-6;         % 光栅周期（m）
theta_i = 30;          % 入射角（deg）
phi = 0;               % 光栅倾角（deg）

% 转换为弧度
theta_i_rad = theta_i * pi/180;
phi_rad = phi * pi/180;

% 计算光栅矢量与入射波矢
K = 2*pi/Lambda;                          % 光栅矢量大小
k0 = 2*pi/lambda;                         % 波数
kx_i = k0*n*sin(theta_i_rad)*cos(phi_rad);  % 入射波矢x分量
ky_i = k0*n*sin(theta_i_rad)*sin(phi_rad);  % 入射波矢y分量
kz_i = k0*n*cos(theta_i_rad);              % 入射波矢z分量

% 计算失配因子
tau = K*cos(phi_rad - theta_i_rad) - 2*pi*n/lambda;  % 失配因子
xi = d*tau/(2*cos(theta_i_rad));                     % 失配项
psi = pi*Delta_n*d/(lambda*cos(theta_i_rad));       % 相位项

% 计算衍射效率
eta = (sinh(psi)^2) / (cosh(psi)^2 - (xi/psi)^2);

% 输出结果
fprintf('反射式体全息光栅衍射效率：%.2f%%\n', eta*100);

3. 扩展：透射式体全息光栅

对于透射式体全息光栅，只需将Kogelnik理论中的相位匹配条件调整为：

其中， $\vec{k}_d$ 为衍射波矢， $\vec{k}_i$ 为入射波矢， $\vec{K}$ 为光栅矢量。衍射效率公式修改为：

四、多种光栅的MATLAB仿真总结

光栅类型	核心理论模型	MATLAB实现方法	关键参数
光纤布拉格光栅（FBG）	耦合模理论（CMT）	传输矩阵法（TMM）、耦合模方程求解	周期 $\Lambda$ 、长度 $L$ 、折射率调制 $\Delta n$
平面透射/反射光栅	夫琅禾费衍射、FMM	傅里叶变换、传输矩阵法	周期 $d$ 、缝宽 $a$ 、缝数 $N$
体全息光栅（VHG）	Kogelnik理论、RCWA	严格耦合波方程求解、Kogelnik公式	周期 $\Lambda$ 、厚度 $d$ 、折射率调制度 $\Delta n$

参考代码多种光栅的反射谱与传输谱MATLAB仿真 www.youwenfan.com/contentcsr/100486.html

五、注意事项与优化技巧

参数校准：仿真前需通过实验测量（如光谱分析仪）校准光栅参数（如 $\Lambda$ 、 $\Delta n$ ），确保模型准确性。
数值稳定性：
- 传输矩阵法中，矩阵连乘顺序（从基底到表面）需正确，避免溢出。
- 对于高耦合系数（ $\kappa \gg \delta$ ）， $\gamma$ 可能为虚数，需用双曲函数（如 $\cosh$ 、 $\sinh$ ）代替三角函数。
偏振处理：
- 倾斜光栅与亚波长光栅的反射/传输谱具有偏振依赖性，需分别计算TE/TM模的结果（如平面光栅的s光与p光衍射效率）。
可视化：用MATLAB Plot或Origin绘制反射谱/传输谱，标注布拉格波长、3 dB带宽、峰值反射率等关键参数，便于分析。

六、应用案例

光纤传感：通过FBG的反射谱偏移（ $\Delta\lambda_B = 2n_{\text{eff}}\Delta\Lambda + 2\Lambda\Delta n_{\text{eff}}$ ）测量温度/应变，仿真时需考虑温度系数（ $d n_{\text{eff}}/dT$ ）与应变系数（ $d\Lambda/d\varepsilon$ ）。
光谱仪：用平面透射光栅的衍射分光特性，仿真时需设计光栅周期（如1 μm）与焦距（如100 mm），以覆盖可见光范围（400-800 nm）。
窄线宽激光器：用VHG的窄光谱带宽（~100 pm）作为激光腔的选频元件，仿真时需优化光栅厚度（ $d$ ）与折射率调制度（ $\Delta n$ ）以提高衍射效率。

七、总结

多种光栅的反射谱/传输谱MATLAB仿真需结合光栅类型选择理论模型（如CMT、RCWA、FMM），并通过代码实现（如传输矩阵法、傅里叶变换）。关键是理解光栅的工作机制（如模式耦合、布拉格衍射），并通过参数优化（如周期、长度、折射率）满足应用需求（如传感、通信、激光）。