理论地震学程序通常涉及地震波在介质中传播的理论计算和模拟。以下是一个基于MATLAB编写的简单理论地震学程序的示例，该程序用于计算和绘制地震波在两种不同介质中的传播情况。
示例程序：基本地震学理论（斯奈尔定律）
该程序使用斯奈尔定律来计算固-固介质界面上的P波和S波反射及透射角。斯奈尔定律描述了入射波和反射波之间的角度关系。
```matlab
function [a1, a2, b1, b2] = snel（i, a, vp1, vp2, vs1, vs2）
% i： 输入波形类型，1表示P波，2表示S波，3表示SH波
% a： 输入入射角
% vp1, vp2： 第一种和第二种介质的P波速度
% vs1, vs2： 第一种和第二种介质的S波速度
% 计算反射和透射角
% 这里只是一个简单的示例，实际应用中可能需要更复杂的计算
a1 = a； % P波入射角等于反射角
a2 = a； % S波入射角等于反射角
b1 = （vp1^2 - vs1^2） / （vp1^2 + vs1^2） * a； % P波透射角
b2 = （vp2^2 - vs2^2） / （vp2^2 + vs2^2） * a； % S波透射角
end
```
使用示例
输入参数
`i`：输入波形类型，1表示P波，2表示S波，3表示SH波。
`a`：输入入射角（以度为单位）。
`vp1`：第一种介质的P波速度（单位：m/s）。
`vp2`：第二种介质的P波速度（单位：m/s）。
`vs1`：第一种介质的S波速度（单位：m/s）。
`vs2`：第二种介质的S波速度（单位：m/s）。
调用函数
```matlab
[a1, a2, b1, b2] = snel（1, 30, 3000, 2500, 1500, 1000）；
```
输出结果
`a1` 和 `a2`：分别为P波和S波的反射角。
`b1` 和 `b2`：分别为P波和S波的透射角。
建议
参数选择：在实际应用中，需要根据具体的地质条件和介质特性选择合适的参数值。
复杂计算：对于更复杂的情况，如非均匀介质或各向异性介质，可能需要更复杂的数学模型和数值方法。
验证和测试：在编写和使用程序时，应进行充分的验证和测试，确保计算结果的准确性和可靠性。
这个示例程序提供了一个基本的框架，可以根据具体需求进行扩展和优化。