Recursion

本文介绍了递归算法的基本概念,通过具体实例分析了递归算法的设计思路,并探讨了递归算法的应用场景及其优缺点。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

http://www.w17x.com/AritcleDisplay.aspx?id=639

The underlying concept is to decompose a big problem into indivisible sub-problems, solve them and then combine each of these partial solutions to get the holistic solution. - Amit Saha

递归算法的定义:如果一个对象的描述中包含它本身,我们就称这个对象是递归的,这种用递归来描述的算法称为递归算法。
我们先来看看大家熟知的一个的故事:
从前有座山,山上有座庙,庙里有个老和尚在给小和尚讲故事,他说从前有座山,山上有座庙,庙里有个老和尚在给小和尚讲故事,他说……
上面的故事本身是递归的,用递归算法描述:
procedure bonze-tell-story;
begin
if 讲话被打断 then 故事结束
else begin
从前有座山,山上有座庙,庙里有个老和尚在给小和尚讲故事;
bonze-tell-story;
end
end;
从上面的递归事例不难看出,递归算法存在的两个必要条件:
(1)    必须有递归的终止条件;
(2)    过程的描述中包含它本身;
在设计递归算法中,如何将一个问题转化为递归的问题,是初学者面临的难题,下面我们通过分析汉诺塔问题,看看如何用递归算法来求解问题;
例1:汉诺塔问题,如下图,有A、B、C三根柱子。A柱子上按从小到大的顺序堆放了N个盘子,现在要把全部盘子从A柱移动到C柱,移动过程中可以借助B柱。移动时有如下要求:
(1)    一次只能移动一个盘子;
(2)    不允许把大盘放在小盘上边;
(3)    盘子只能放在三根柱子上;

算法分析:当盘子比较多的时,问题比较复杂,所以我们先分析简单的情况:
如果只有一个盘子,只需一步,直接把它从A柱移动到C柱;
如果是二个盘子,共需要移动3步:
(1)    把A柱上的小盘子移动到B柱;
(2)    把A柱上的大盘子移动到C柱;
(3)    把B柱上的大盘子移动到C柱;
如果N比较大时,需要很多步才能完成,我们先考虑是否能把复杂的移动过程转化为简单的移动过程,如果要把A柱上最大的盘子移动到C柱上去,必须先把上面的N-1个盘子从A柱移动到B柱上暂存,按这种思路,就可以把N个盘子的移动过程分作3大步:
(1)    把A柱上面的N-1个盘子移动到B柱;
(2)    把A柱上剩下的一个盘子移动到C柱;
(3)    把B柱上面的N-1个盘子移动到C柱;
其中N-1个盘子的移动过程又可按同样的方法分为三大步,这样就把移动过程转化为一个递归的过程,直到最后只剩下一个盘子,按照移动一个盘子的方法移动,递归结束。
递归过程:
procedure Hanoi(N,A,B,C:integer;);{以B柱为中转柱将N个盘子从A柱移动到C柱}
begin
if N=1 then write(A,’->’,C){把盘子直接从A移动到C}
else begin
Hanoi(N-1,A,C,B);{ 以C柱为中转柱将N-1个盘子从A柱移动到B柱}
write(A,’->’,C);{把剩下的一个盘子从A移动到C}
Hanoi(N-1,B,A,C); { 以A柱为中转柱将N-1个盘子从B柱移动到C柱}
end;
end;
从上面的例子我们可以看出,在使用递归算法时,首先弄清楚简单情况下的解法,然后弄清楚如何把复杂情况归纳为更简单的情况。
在信息学奥赛中有的问题的结构或所处理的数据本身是递归定义的,这样的问题非常适合用递归算法来求解,对于这类问题,我们把它分解为具有相同性质的若干个子问题,如果子问题解决了,原问题也就解决了。
例2求先序排列 (NOIP2001pj)
[问题描述]给出一棵二叉树的中序与后序排列。求出它的先序排列。(约定树结点用不同的大写字母表示,长度≤8)。
[样例] 输入:BADC BDCA   输出:ABCD
算法分析:我们先看看三种遍历的定义:
先序遍历是先访问根结点,再遍历左子树,最后遍历右子树;
中序遍历是先遍历左子树,再访问根结点,最后遍历右子树;
后序遍历是先遍历左子树,再遍历右子树,最后访问根结点;
从遍历的定义可知,后序排列的最后一个字符即为这棵树的根节点;在中序排列中,根结点前面的为其左子树,根结点后面的为其右子树;我们可以由后序排列求得根结点,再由根结点在中序排列的位置确定左子树和右子树,把左子树和右子树各看作一个单独的树。这样,就把一棵树分解为具有相同性质的二棵子树,一直递归下去,当分解的子树为空时,递归结束,在递归过程中,按先序遍历的规则输出求得的各个根结点,输出的结果即为原问题的解。
源程序
program noip2001_3;
var z,h : string;
procedure make(z,h:string); {z为中序排列,h为后序排列}
var s,m : integer;
begin
m:=length(h);{m为树的长度}
 write(h[m]); {输出根节点}
 s:=pos(h[m],z); {求根节点在中序排列中的位置}
 if s>1 then make(copy(z,1,s-1),copy(h,1,s-1)); {处理左子树}
 if m>s then make(copy(z,s+1,m-s),copy(h,s,m-s)); {处理右子树}
end;
begin
 readln(z);
 readln(h);
 make(z,h);
end.
递归算法不仅仅是用于求解递归描述的问题,在其它很多问题中也可以用到递归思想,如回溯法、分治法、动态规划法等算法中都可以使用递归思想来实现,从而使编写的程序更加简洁。
比如上期回溯法所讲的例2《数的划分问题》,若用递归来求解,程序非常短小且效率很高,源程序如下
var
n,k:integer;
tol:longint;
procedure make(sum,t,d:integer);
var i:integer;
begin
if d=k then inc(tol)
else for i:=t to sum div 2 do make(sum-i,i,d+1);
end;
begin
readln(n,k);
tol:=0;
make(n,1,1);
writeln(tol);
end.

有些问题本身是递归定义的,但它并不适合用递归算法来求解,如斐波那契(Fibonacci)数列,它的递归定义为:
F(n)=1          (n=1,2)
F(n)=F(n-2)+F(n-1) (n>2)
用递归过程描述为:
Funtion fb(n:integer):integer;
Begin
if n<3 then fb:=1
else fb:=fb(n-1)+fb(n-2);
End;
上面的递归过程,调用一次产生二个新的调用,递归次数呈指数增长,时间复杂度为O(2n),把它改为非递归:
x:=1;y:=1;
for i:=3 to n do
begin
z:=y;y:=x+y;x:=z;
end;
修改后的程序,它的时间复杂度为O(n)。
我们在编写程序时是否使用递归算法,关键是看问题是否适合用递归算法来求解。由于递归算法编写的程序逻辑性强,结构清晰,正确性易于证明,程序调试也十分方便,在NOIP中,数据的规模一般也不大,只要问题适合用递归算法求解,我们还是可以大胆地使用递归算法。

import numpy as np &#39;&#39;&#39; ---L、U计算模块--- 包含: Calc(...), 给定孔径系数、视场系数、物面、镜面、入瞳信息,光线类型,利用迭代递推模块进行中间计算,返回各面L、U、I、PA信息。 &#39;&#39;&#39; def Calc(K_eta, K_W, object, lens, pupil, image, raytype): """ 作用:计算某一孔径、视场光线在逐个面的(L,U)坐标数据 输入:K_eta为孔径取点系数,K_W为视场取点系数 """ #初始化(L,U)坐标列表 L0 = [] L1 = [] U = [] I0 = [] I1 = [] PA = [] A = pupil.semi_dia #入瞳半径 eta = K_eta * A #像高 #物在无穷远且光线与光轴平行时的情形 if(object.thickness == None and K_W == 0): h0 = K_eta * pupil.semi_dia L0, L1, U, I0, I1, PA = lens_recursion(0, 0, h0, lens, raytype, 1) #其余情形 -> 转化为轴上点在有限远处 else: if((object.thickness == None and K_W != 0)): U = -K_W * object.semi_dia * np.pi/180 L = pupil.thickness + eta/np.tan(U) elif(object.semi_dia == 0): L = -object.thickness U = np.arcsin(K_eta * A/np.sqrt(A**2 + Lp**2)) elif(object.semi_dia != 0): Lp = pupil.thickness #第一面到入瞳的距离 y = K_W * object.semi_dia #物方线视场 U = np.arctan((y-eta)/(Lp+object.thickness)) L = Lp + eta/((y-eta)/(Lp+object.thickness)) L0, L1, U, I0, I1, PA = lens_recursion(L, U, 0, lens, raytype, 0) return L0, L1, U, I0, I1, PA &#39;&#39;&#39; ---迭代递推模块--- 包含: lens_recursion(...), 给定入射光线L、U,对各镜面进行光线追迹并返回各面L、U数据 astigmatism_recursion(...), 给定各面L、U、I、PA参数,求子午场曲,弧矢场曲,像散 &#39;&#39;&#39; EPS = 1e-5 #PA校验精度 def lens_recursion(L, U, H, lens, raytype, isParallel): """ 函数作用: 给定轴上光线入射初值, 对每个镜面进行递推 输入: 初始物方截距(L)、物方倾斜角(U)、镜面数据(lens)、光线类型(raytype = &#39;d&#39;, &#39;F&#39;, &#39;C&#39;) 是否平行于主光轴(isParallel, 1平行, 0倾斜) 输出: 光线经过每面镜面后的L, L&#39;, U&#39;, I, I&#39;, PA数组 """ snumber = len(surface) recursion = 0 sinI = 0; sinII = 0 I = 0; II = 0 UU = 0 LL = 0 L_data = [L] LL_data = [] UU_data = [U] I_data = [] II_data = [] PA_data = [] if(raytype == &#39;d&#39;): #d光,折射率取nd for recursion in range(lensnumber): if(recursion == 0 and isParallel): sinI = H / lens[recursion].radius else: sinI = (L - lens[recursion].radius) / lens[recursion].radius * np.sin(U) sinII = lens[recursion].nd0 / lens[recursion].nd1 * sinI if(np.abs(sinI) < 1 and np.abs(sinII) < 1): I = np.arcsin(sinI); II = np.arcsin(sinII) I_data.append(I); II_data.append(II) UU = U + I - II LL = lens[recursion].radius + lens[recursion].radius * sinII / np.sin(UU) LL_data.append(LL); UU_data.append(UU) #PA校对法 PA1 = L * np.sin(U) / np.cos((I-U)/2) PA2 = LL * np.sin(UU) / np.cos((II-UU)/2) if(np.abs(PA1-PA2) > EPS and isParallel == 0): print("Recursion not right.") print(f&#39;PA1 = {PA1}, PA2 = {PA2}&#39;) PA_data.append(PA2) #过渡公式 U = UU L = LL - lens[recursion].thickness #舍去最后一个过渡的物方截距 if(recursion != lensnumber - 1): L_data.append(L) if(raytype == &#39;F&#39;): #F光,折射率取nF for recursion in range(lensnumber): if(recursion == 0 and isParallel): sinI = H / lens[recursion].radius else: sinI = (L - lens[recursion].radius) / lens[recursion].radius * np.sin(U) sinII = lens[recursion].nF0 / lens[recursion].nF1 * sinI if(np.abs(sinI) < 1 and np.abs(sinII) < 1): I = np.arcsin(sinI); II = np.arcsin(sinII) I_data.append(I); II_data.append(II) UU = U + I - II LL = lens[recursion].radius + lens[recursion].radius * sinII / np.sin(UU) LL_data.append(LL); UU_data.append(UU) #PA校对法 PA1 = L * np.sin(U) / np.cos((I-U)/2) PA2 = LL * np.sin(UU) / np.cos((II-UU)/2) if(np.abs(PA1-PA2) > EPS and isParallel == 0): print("Recursion not right.") print(f&#39;PA1 = {PA1}, PA2 = {PA2}&#39;) PA_data.append(PA2) #过渡公式 U = UU L = LL - lens[recursion].thickness #舍去最后一个过渡的物方截距 if(recursion != lensnumber - 1): L_data.append(L) if(raytype == &#39;C&#39;): #C光,折射率取nC for recursion in range(lensnumber): if(recursion == 0 and isParallel): sinI = H / lens[recursion].radius else: sinI = (L - lens[recursion].radius) / lens[recursion].radius * np.sin(U) sinII = lens[recursion].nC0 / lens[recursion].nC1 * sinI if(np.abs(sinI) < 1 and np.abs(sinII) < 1): I = np.arcsin(sinI); II = np.arcsin(sinII) I_data.append(I); II_data.append(II) UU = U + I - II LL = lens[recursion].radius + lens[recursion].radius * sinII / np.sin(UU) LL_data.append(LL); UU_data.append(UU) #PA校对法 PA1 = L * np.sin(U) / np.cos((I-U)/2) PA2 = LL * np.sin(UU) / np.cos((II-UU)/2) if(np.abs(PA1-PA2) > EPS and isParallel == 0): print("Recursion not right.") print(f&#39;PA1 = {PA1}, PA2 = {PA2}&#39;) PA_data.append(PA2) #过渡公式 U = UU L = LL - lens[recursion].thickness #舍去最后一个过渡的物方截距 if(recursion != lensnumber - 1): L_data.append(L) #print(f&#39;LEN_DATA\n{len(L_data)}, {len(LL_data)}, {len(UU_data)}, {len(I_data)}, {len(II_data)}, {len(PA_data)}&#39;) return L_data, LL_data, UU_data, I_data, II_data, PA_data def astigmatism_recursion(L_data, LL_data, UU_data, I_data, II_data, PA_data, lens, raytype, isParallel): """ 函数作用: 求像散, 场曲 输入: 由lens_recursion算得的各参数 输出: 子午场曲(xt), 弧矢场曲(xs), 像散(delta_x) """ lensnumber = len(lens) for recursion in range(lensnumber): if(raytype == &#39;d&#39;): x = PA_data[recursion] ** 2 / (2 * lens[recursion].radius) #求t、s初值 if(recursion == 0 and isParallel == 0): t = s = (L_data[0] - x) / np.cos(UU_data[0]) if(recursion == 1 and isParallel == 1): t = s = (L_data[1] - x) / np.cos(UU_data[1]) temp_tt = (lens[recursion].nd1 * np.cos(II_data[recursion]) - lens[recursion].nd0 * np.cos(I_data[recursion]))\ / lens[recursion].radius + lens[recursion].nd0 * np.cos(I_data[recursion]) ** 2 / t tt = 1 / temp_tt * lens[recursion].nd1 * np.cos(II_data[recursion]) ** 2 #t&#39; temp_ss = (lens[recursion].nd1 * np.cos(II_data[recursion]) - lens[recursion].nd0 * np.cos(I_data[recursion]))\ / lens[recursion].radius + lens[recursion].nd0 / s ss = 1 / temp_ss * lens[recursion].nd1 #s&#39; if(recursion < lensnumber-1): #过渡公式 xx = PA_data[recursion + 1] ** 2 / (2 * lens[recursion+1].radius) #x(i+1) D = (lens[recursion].thickness - x + xx) / np.cos(UU_data[recursion+1]) t = tt - D s = ss - D if(recursion == lensnumber-1): #最后一面 lt = tt * np.cos(UU_data[recursion + 1]) + x ls = ss * np.cos(UU_data[recursion + 1]) + x xt = lt - LL_data[-1] xs = ls - LL_data[-1] delta_x = xt - xs if(raytype == &#39;C&#39;): x = PA_data[recursion] ** 2 / (2 * lens[recursion].radius) #求t、s初值 if(recursion == 0 and isParallel == 0): t = s = (L_data[0] - x) / np.cos(UU_data[0]) if(recursion == 1 and isParallel == 1): t = s = (L_data[1] - x) / np.cos(UU_data[1]) temp_tt = (lens[recursion].nC1 * np.cos(II_data[recursion]) - lens[recursion].nC0 * np.cos(I_data[recursion]))\ / lens[recursion].radius + lens[recursion].nC0 * np.cos(I_data[recursion]) ** 2 / t tt = 1 / temp_tt * lens[recursion].nC1 * np.cos(II_data[recursion]) ** 2 #t&#39; temp_ss = (lens[recursion].nC1 * np.cos(II_data[recursion]) - lens[recursion].nC0 * np.cos(I_data[recursion]))\ / lens[recursion].radius + lens[recursion].nC0 / s ss = 1 / temp_ss * lens[recursion].nC1 #s&#39; if(recursion < lensnumber-1): #过渡公式 xx = PA_data[recursion + 1] ** 2 / (2 * lens[recursion+1].radius) #x(i+1) D = (lens[recursion].thickness - x + xx) / np.cos(UU_data[recursion+1]) t = tt - D s = ss - D if(recursion == lensnumber-1): #最后一面 lt = tt * np.cos(UU_data[recursion + 1]) + x ls = ss * np.cos(UU_data[recursion + 1]) + x xt = lt - LL_data[recursion] xs = ls - LL_data[recursion] delta_x = xt - xs if(raytype == &#39;F&#39;): x = PA_data[recursion] ** 2 / (2 * lens[recursion].radius) #求t、s初值 if(recursion == 0 and isParallel == 0): t = s = (L_data[0] - x) / np.cos(UU_data[0]) if(recursion == 1 and isParallel == 1): t = s = (L_data[1] - x) / np.cos(UU_data[1]) temp_tt = (lens[recursion].nF1 * np.cos(II_data[recursion]) - lens[recursion].nF0 * np.cos(I_data[recursion]))\ / lens[recursion].radius + lens[recursion].nF0 * np.cos(I_data[recursion]) ** 2 / t tt = 1 / temp_tt * lens[recursion].nF1 * np.cos(II_data[recursion]) ** 2 #t&#39; temp_ss = (lens[recursion].nF1 * np.cos(II_data[recursion]) - lens[recursion].nF0 * np.cos(I_data[recursion]))\ / lens[recursion].radius + lens[recursion].nF0 / s ss = 1 / temp_ss * lens[recursion].nF1 #s&#39; if(recursion < lensnumber-1): #过渡公式 xx = PA_data[recursion + 1] ** 2 / (2 * lens[recursion+1].radius) #x(i+1) D = (lens[recursion].thickness - x + xx) / np.cos(UU_data[recursion+1]) t = tt - D s = ss - D if(recursion == lensnumber-1): #最后一面 lt = tt * np.cos(UU_data[recursion + 1]) + x ls = ss * np.cos(UU_data[recursion + 1]) + x xt = lt - LL_data[recursion] xs = ls - LL_data[recursion] delta_x = xt - xs return xt, xs, delta_x &#39;&#39;&#39; ---数据处理模块--- 包含: data_processing(...), 总输出函数,给定物面、镜面、入瞳、像面信息,向主函数返回所有待求量 Calc_ff(...), 计算计算像方焦距、焦点位置、主点位置 Calc_ll(...), 计算像距、像高、像散、场曲 Calc_lp(...), 计算出瞳距 Calc_sp(...), 计算全孔径和0.7孔径时的实际像距、球差、位置色差 Calc_coma(...), 计算全视场、0.7视场、正负全孔径、正负0.7孔径的彗差、两个视场dFC光实际像高,绝对畸变、相对畸变、倍率色差 &#39;&#39;&#39; def data_processing(object, lens, pupil, image): &#39;&#39;&#39; 共输出38个数据 data1:f&#39;, l&#39;, lH&#39;, lp&#39;, y0&#39;, y0, 0.707&#39;, xt&#39;, xs&#39;, delta_x&#39;, F光近轴像位置, C光近轴像位置;共12个数据 data2:轴上点全孔径的 d光像点, F光像点, C光像点, 球差, 位置色差;共5个数据 data3:轴上点0.707孔径的 d光像点, F光像点, C光像点, 球差, 位置色差;共5个数据 data4:轴上点0孔径位置色差;共1个数据 data5:全视场的 全孔径子午彗差, 0.7孔径子午彗差, d光实际像高, F光实际像高, C光实际像高, 绝对畸变, 相对畸变, 倍率色差;共8个数据 data6:0.707视场的 全孔径子午彗差, 0.7孔径子午彗差, d光实际像高, F光实际像高, C光实际像高, 绝对畸变, 相对畸变, 倍率色差;共8个数据 &#39;&#39;&#39; oo = 1e-20 #定义近轴小量 ff, l_H = Calc_ff(object, lens, pupil, image, oo) ll, ll_F, ll_C, y0, y0_707, xt, xs, delta_x = Calc_ll(object, lens, pupil, image, oo, l_H, ff) l_p = Calc_lp(lens, pupil, image, oo) data1 = [ff, ll, l_H, l_p, y0, y0_707, xt, xs, delta_x, ll_F, ll_C] data2, data3, data4 = Calc_sp(object, lens, pupil, image, ll, oo) data5, data6 = Calc_coma(object, lens, pupil, image, oo, y0, y0_707) return data1, data2, data3, data4, data5, data6 #1、计算像方焦距、焦点位置、主点位置 def Calc_ff(object, lens, pupil, image, oo): #定义一个近轴无穷远物面 object_ = Surface(0, None, object.nd1, object.nF1, object.nC1, None) #计算平行入射的近轴光线 L0 = Calc(oo, 0, object_, lens, pupil, image, &#39;d&#39;)[0] L1 = Calc(oo, 0, object_, lens, pupil, image, &#39;d&#39;)[1] #计算L0,L1的乘积 L0_product = L1_product = L0F_product = L1F_product = L0C_product = L1C_product = 1 for i in range(len(L0)-1): L0_product = L0_product * L0[i+1] for i in range(len(L1)): L1_product = L1_product * L1[i] #利用正切计算法计算焦距 ff = L1_product/L0_product #计算像方焦点和主点位置 l_F = L1[-1] l_H = ff - l_F return ff, l_H 这个是正确计算的代码,但是写的比较混乱,下方代码计算结果不正确import tkinter as tk from tkinter import ttk, messagebox, filedialog import json import os import math import numpy as np class OpticalSystem: def __init__(self): self.surfaces = [] # 存储光学表面对象 self.entrance_pupil_diameter = None # 入瞳直径 (mm) self.entrance_pupil_position = None # 入瞳位置 (相对第一个面顶点) (mm) self.object_infinite = True # 物在无穷远 self.field_angle = None # 半视场角 (度) self.object_distance = None # 物距 (有限远) (mm) self.object_height = None # 物高 (有限远) (mm) self.aperture_angle = None # 孔径角 (有限远) (度) self.light_type = "d" # 色光类型,默认d光 self.ray_paths = [] # 存储光线路径数据 # 计算结果存储 self.results = { "focal_length": None, # 焦距 f&#39; "ideal_image_distance": None, # 理想像距 l&#39; "actual_image_position": None, # 实际像位置 "image_principal_plane": None, # 像方主面位置 lH&#39; "exit_pupil_distance": None, # 出瞳距 lp&#39; "ideal_image_height": None, # 理想像高 y0&#39; "spherical_aberration": None, # 球差 "longitudinal_chromatic": None, # 位置色差 "tangential_field_curvature": None, # 子午场曲 xt&#39; "sagittal_field_curvature": None, # 弧矢场曲 xs&#39; "astigmatism": None, # 像散 &Delta;xts&#39; "actual_image_height": None, # 实际像高 "relative_distortion": None, # 相对畸变 "absolute_distortion": None, # 绝对畸变 "lateral_chromatic": None, # 倍率色差 "tangential_coma": None # 子午慧差 } class Surface: def __init__(self, r=float(&#39;inf&#39;), d=0.0, nd=1.0, vd=0.0): self.r = r # 曲率半径 (mm) self.d = d # 厚度/间隔 (mm) self.nd = nd # d光折射率 self.vd = vd # 阿贝数 def to_dict(self): """将表面数据转换为字典""" return { "r": self.r, "d": self.d, "nd": self.nd, "vd": self.vd } @classmethod def from_dict(cls, data): """从字典创建表面对象""" return cls( r=data.get(&#39;r&#39;, float(&#39;inf&#39;)), d=data.get(&#39;d&#39;, 0.0), nd=data.get(&#39;nd&#39;, 1.0), vd=data.get(&#39;vd&#39;, 0.0) ) class Calculate: EPS = 1e-5 # PA校验精度 @staticmethod def recursion(L, U, H, surfaces, P, light_type=&#39;d&#39;): """ 光线追迹递归计算 :param L: 物方截距 :param U: 物方孔径角 :param H: 光线高度 :param surfaces: 光学表面列表 :param P: PA校验标志 :param light_type: 色光类型 :return: 光线路径数据 """ # 数据初始化 snumber = len(surfaces) sinI = 0; sinII = 0 I = 0; II = 0 UU = 0 LL = 0 L_data = [L] LL_data = [] UU_data = [U] I_data = [] II_data = [] PA_data = [] # 前一个面的折射率(初始为空气) n_prev = 1.0 for i, surf in enumerate(surfaces): # 根据色光类型选择折射率 if light_type == &#39;d&#39;: n_current = surf.nd elif light_type == &#39;f&#39;: # F光折射率计算 (nd + (nf - nc)/(vd * (nf - nc)/ (nd - 1))) n_current = surf.nd + 1/(surf.vd * (1.0/(surf.nd - 1.0))) elif light_type == &#39;c&#39;: # C光折射率计算 n_current = surf.nd - 1/(surf.vd * (1.0/(surf.nd - 1.0))) else: n_current = surf.nd # 默认为d光 # 计算入射角 if i == 0 and P: # 第一面且平行光入射 sinI = H / surf.r if surf.r != 0 else 0 sinII = (n_prev / n_current) * sinI else: # 普通情况 sinI = (L_data[-1] - surf.r) / surf.r * np.sin(U) if surf.r != 0 else 0 sinII = (n_prev / n_current) * sinI # 确保sin值在合理范围内 sinI = np.clip(sinI, -1.0, 1.0) sinII = np.clip(sinII, -1.0, 1.0) I = np.arcsin(sinI) II = np.arcsin(sinII) I_data.append(I) II_data.append(II) # 计算折射后角度 UU = U + I - II UU_data.append(UU) # 计算折射后截距 if abs(sinII) > 1e-10: # 避免除以0 LL = surf.r + surf.r * sinII / np.sin(UU) else: LL = surf.r # 垂直入射 LL_data.append(LL) # PA校对法 if abs(U) > 1e-5 and abs(UU) > 1e-5: # 避免除以0 PA1 = L_data[-1] * np.sin(U) / np.cos((I - U)/2) PA2 = LL * np.sin(UU) / np.cos((II - UU)/2) if abs(PA1 - PA2) > Calculate.EPS and P == 0: print(f"Warning: PA校验不通过! 表面 {i+1}: PA1={PA1:.6f}, PA2={PA2:.6f}") PA_data.append(PA2) else: PA_data.append(0.0) # 无法计算PA值 # 过渡到下一个面 U = UU if i < snumber - 1: # 最后一个面不需要过渡 L = LL - surf.d L_data.append(L) # 更新折射率为当前面的折射率 n_prev = n_current return L_data, LL_data, UU_data, I_data, II_data, PA_data @staticmethod def calc(system, K_eta=0.0, K_W=0.0): """ 计算某一孔径、视场光线在逐个面的(L,U)坐标数据 :param system: OpticalSystem对象 :param K_eta: 孔径取点系数 (0-1) :param K_W: 视场取点系数 (0-1) :return: 光线路径数据 """ # 初始化(L,U)坐标列表 L0 = [] U_data = [] I0 = [] I1 = [] PA = [] A = system.entrance_pupil_diameter / 2 # 入瞳半径 eta = K_eta * A # 实际像高系数 # 物在无穷远且光线与光轴平行时的情形 if system.object_infinite and K_W == 0: h0 = K_eta * A L0, LL, U_data, I0, I1, PA = Calculate.recursion( 0, 0, h0, system.surfaces, True, system.light_type ) else: # 物在无穷远但视场不为0 if system.object_infinite: W = K_W * system.field_angle * np.pi / 180 # 转换为弧度 U_angle = -W L_start = system.entrance_pupil_position + eta / np.tan(U_angle) else: # 有限远物距 Lp = system.entrance_pupil_position # 第一面到入瞳的距离 y = K_W * system.object_height # 物方线视场 if K_W == 0: # 轴上点 dist = np.sqrt(Lp**2 + A**2) U_angle = np.arcsin(K_eta * A / dist) L_start = -system.object_distance else: # 轴外点 U_angle = np.arctan((y - eta) / (Lp + system.object_distance)) L_start = Lp + eta / np.tan(U_angle) L0, LL, U_data, I0, I1, PA = Calculate.recursion( L_start, U_angle, 0, system.surfaces, False, system.light_type ) return L0, LL, U_data, I0, I1, PA @staticmethod def calculate_focal_length(system): """计算焦距""" # 计算近轴光线 _, LL, U_data, _, _, _ = Calculate.calc(system, K_eta=0.0001, K_W=0.0) last_U = U_data[-1] last_L = LL[-1] # 焦距 f&#39; = h / tan(U&#39;) if abs(last_U) > 1e-5: focal_length = system.entrance_pupil_diameter / 2 / np.tan(last_U) return focal_length return float(&#39;inf&#39;) # 平行光出射 @staticmethod def calculate_ideal_image_distance(system): """计算理想像距""" # 计算主光线 L0, LL, U_data, _, _, _ = Calculate.calc(system, K_eta=0.0, K_W=0.0) last_L = LL[-1] return last_L @staticmethod def calculate_spherical_aberration(system): """计算球差""" # 计算边缘光线 _, LL_edge, U_edge, _, _, _ = Calculate.calc(system, K_eta=1.0, K_W=0.0) # 计算近轴光线 _, LL_paraxial, U_paraxial, _, _, _ = Calculate.calc(system, K_eta=0.0001, K_W=0.0) # 球差 &delta;L&#39; = L&#39;边缘 - L&#39;近轴 return LL_edge[-1] - LL_paraxial[-1]请说说第二段代码运算逻辑哪里有问题
最新发布
06-27
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值