用python写一个平行光射向旋转椭球体后只接受一定角度反射光线的程序

好的，以下是一个基于Python的程序，用于模拟平行光射向旋转椭球体并只接受一定角度反射光线的情况：

import numpy as np

# 定义旋转矩阵
def rot_matrix(theta):
    return np.array([[np.cos(theta), 0, np.sin(theta)],
                     [0, 1, 0],
                     [-np.sin(theta), 0, np.cos(theta)]])

# 定义椭球体方程
def ellipsoid(x, y, z, a, b, c):
    return (x/a)**2 + (y/b)**2 + (z/c)**2

# 定义反射函数
def reflect(in_vec, norm_vec):
    return in_vec - 2 * np.dot(in_vec, norm_vec) * norm_vec

# 定义入射光线方向向量
inc_vec = np.array([0, 1, 0])

# 定义旋转角度
theta = np.pi/4

# 定义椭球体半轴长度
a = 5
b = 3
c = 4

# 定义椭球体中心点
center = np.array([0, 0, 0])

# 通过旋转矩阵对椭球体进行旋转
rot_mat = rot_matrix(theta)
a_vec = np.array([a, 0, 0])
b_vec = np.array([0, b, 0])
c_vec = np.array([0, 0, c])
a_vec = np.dot(rot_mat, a_vec)
b_vec = np.dot(rot_mat, b_vec)
c_vec = np.dot(rot_mat, c_vec)

# 对椭球体进行采样，计算反射光线
sample_num = 100
for i in range(sample_num):
    for j in range(sample_num):
        for k in range(sample_num):
            x = -a + 2*a/sample_num*i
            y = -b + 2*b/sample_num*j
            z = -c + 2*c/sample_num*k
            if ellipsoid(x, y, z, a, b, c) <= 1:
                norm_vec = np.array([x/a**2, y/b**2, z/c**2])
                out_vec = reflect(inc_vec, norm_vec)
                if out_vec[1] > np.tan(np.pi/6): # 只接受一定角度的反射光线
                    print('反射光线：', out_vec)

上述程序中，我们首先定义了一个旋转矩阵函数 `rot_matrix`，用于将椭球体进行旋转。然后定义了一个椭球体方程函数 `ellipsoid`，用于判断一个点是否在椭球体内部。接着定义了一个反射函数 `reflect`，用于计算入射光线经过一个法向量的反射光线方向。然后定义了一个入射光线方向向量 `inc_vec`，这里我们默认入射光线方向为 y 轴正方向。然后定义了一个旋转角度 `theta`，这里我们默认旋转角度为 pi/4。然后定义了椭球体半轴长度和中心点。接着通过旋转矩阵对椭球体进行旋转。最后通过采样的方式对椭球体进行采样，计算反射光线，只接受一定角度反射光线的光线输出。