python 实现光线追迹（中）：空间关系
⽂章⽬录
空间关系
变化始于相遇，所以交点是⼀切的核⼼。
相交判定⾸先考察⼀束光线能否打在某个平⾯镜上。光线被抽象成了⼀个列表[a,b,c]，平⾯镜则被抽象成为由两个点构成的线段[(x1,y1),(x2,y2)]。两条直线的交点问题属于初等数学范畴，需要先将线段转换成直线的形式，然后再求交点。但是两条直线的交点可能落在线段的外⾯，从⽽不具有判定的意义。
如果我们的光学系统中有⼤量的光学元件，那么如果有⼀种⽅法可以快速判断光线是否与光学元件有交点，将会显得更加快捷。那么，如果⼀个直线穿过某个线段，就意味着这条线段的两个端点必然在直线的两侧。
对于直线，如果为变量，则构成⼀组直线族，那么随着的不断变化，直线在坐标系中处于⼀个不断平移的过程，这个过程是连续的。假设我们的光线⽅程为，那么分布在这个光线两侧的同族直线，，其也必然分布在两侧。
假设，，那么必有，这样就得到了⼀个直线与线段是否相交的⼀个判定式。
import  numpy as  np
def  isCross (abc ,dots ):
abc = np .array (abc )                #将abc 的格式转成数组
flag1 = abc .dot (list (dots [0])+[1])    #数组的点积
flag2 = abc .dot (list (dots [1])+[1])
return  flag1*flag2我们⾮常熟悉运算符"+"，不过⽬前来说，只有数值和数组是⽀持加法运算的。所以，对于list(dots[0])+[1]这种表⽰着实让⼈有些摸不到头脑。
这个含义其实是符合⼈类直觉的。列表内的元素个数是可变的，两个列表相加可以理解为两个列表衔接在⼀起。当然，元组并不⽀持这种运算。
例如
>>> [1,2,3]+[4]
[1, 2, 3, 4]
>>> (1,2,3)+(4)
Traceback (most recent call last ):
File "<stdin>", line 1, in  <module >
TypeError : can only concatenate tuple  (not  "int") to tuple
通过加法运算符来衔接两个列表，实际上相当于新建了⼀个变量，需要开辟新的内存空间。好在对于初学者来说这样不容易出错。在numpy 中，+、-、*、/这⼏个运算符表⽰对应位置元素的运算。如果想使⽤点乘等其他运算，需要调⽤numpy 中的其他函数。>>> np .array ([1,2,3])*np .array ([4,5,6])
array ([ 4, 10, 18])
>>> np .array ([1,2,3])+np .array ([4,5,6])
array ([5, 7, 9])
>>> np .array ([1,2,3]).dot ([4,5,6])          #.dot 表点积
32
所以，。当然，我们也可以写成flag1 = abc[0]*dots[0][0]+abc[1]*dots[0][1]+c ，只是看上去不太优雅。
ax +by +c =0c c ax +by +c =00ax +by +c =11ax +by +c =20c ,c 12c 0ax +1by +1c =10ax +2by +2c =20(ax +1by +1c )(ax +02by +2c )<00flag 1=[a ,b ,c ]⋅[x ,y ,1]=00ax +0by +0c
然后，得到了flag1和flag2的值，如果⼆者异号，那么就可以断定，直线在两个点的中间。也就是说，只要flag1*flag2<0，即可说明直线与线段有焦点。
当然，从实际的⾓度出发，我们可以不去考虑光线通过镜⽚的端点或者平⾏镜⽚并掠过这种情况，从⽽只要函数返回值⼩于零，就认定⼆者相交，否则不相交。然⽽，从数学的⾓度出发，直线和线段之间可能存在三种关系：不相交、有⼀个交点、线段在直线上。
虽然这个理解没什么实际价值，但对于python学习来说，却是⾮常有意义的⼀个例⼦，代码如下，看懂了这个代码，那么就差不多可以算是⼀个初级的python程序员了。
def isCross(abc,dots):
abc = np.array(abc)
flags =[abc.dot(list(p)+[1])for p in dots]#for表达式
poses,negs,zeros =[0,0,0]
for flag in flags:#循环语句
if flag >0:#判断语句
poses +=1
elif flag <0:
negs +=1
else:
zeros +=1
return poses*negs+zeros
这个短短的程序涵盖了循环语句、判断语句以及for表达式的内容，前两者是最基础的编程知识，后者是python中⾮常亮眼的⼀种功能。
⾸先来认识⼀下运算符+=，poses += 1即为poses = poses + 1，即相当于将poses+1赋值给poses。赋值前后flag1在内存中的位置发⽣了变化，也就是说flag1已经不是原来的flag1了。在这⾥，等号也并不能读成等于，⽽是读成被赋值为。即poses被赋值为poses+1。前⾯略有提过，双等号==才表⽰真正的相等。
然后来看判断语句，对于表达式if a A elif b B else C，我们按照⼈类的语法去读即可：如果a成⽴，则执⾏A，如果b成⽴，则执⾏B，否则的话执⾏C。在上述代码中，也可以很⽅便地读成：如果flag>0，那么poses被赋值为poses+1；如果flag<0，那么negs变成negs+1；否则的话zeros变成zeros+1。
这⼏个变量顾名思义，poses表⽰正数的个数，negs表⽰负数的个数，zeros表⽰0的个数。
然后来看[abc.dot(list(p)+[1]) for p in dots]，我们⾸先使⽤⼀种没有陌⽣字符的⽅式书写：
#这是伪代码，假设dots中有n个变量，表⽰创建flag1、flag2⼀直到flagn⼀共n个变量。
flag1 = abc.dot(list(dots[0])+[1])
flag2 = abc.dot(list(dots[1])+[1])
...
flagn = abc.dot(list(dots[n])+[1])
这个表达式⾮常规律，这n个变量相当于是在dots中循环⼀遍，然后逐个赋值。for p in dots表⽰的就是将dots中的元素取出，赋值给p，然后再对p进⾏操作abc.dot(list(p)+[1])，最后将所有操作得到的值包裹再⼀个list中。最后再记⼀下这个表达形式[abc.dot(list(p)+[1]) for p in dots]，以后会经常⽤到。
ax+by+c
最后来看for flag in flags，即拿出flags中的所有元素，循环操作其下⽅的代码块。flags中的元素即为两个点分别带⼊之后的值。那么对于这两个点来说，如果⼀正⼀负，则poses*negs=1*1=1，此时代表直线和线段有⼀交点，否则这个值便为0。当
poses*negs==0时，则zeros的个数表⽰端点与直线相交的个数，zeros为0，表⽰⽆交点，为1，表⽰有⼀个端点在直线上，为2表⽰两个端点都在直线上。
射线排序
现在，我们可以判断某⼀个线段与⼀条直线是否有交点了，那么如果空间中有多个平⾯镜，光线所在的直线⼜与许多平⾯镜有交点，那么应该如何到最近的那个呢？最简单的⽅法是分别求取这些点到光源的距离，距离越近相交越早。但这样会产⽣⼀个问题，即难以判定这个最近点是否在光的传播路
径上，如果这个点在光源的后⾯，那就⽐较尴尬了。
所以，⽐较稳妥的⽅法是，按照射线的⽅向对所有点进⾏排序，那么光源后⾯的那个点，就是光线传播过程中的第⼀个交点。
刚刚我们在判定直线与线段的交点时，提到了直线族的概念。发现对于a、b取值相同的⼀组直线来说，其c值的⼤⼩与直线族的顺序是密切相关的。如Fig2-2所⽰。其到
依次递减。这启发我们需要构建出⼀组和光线想垂直的直线族[a,b,~]，则对于空间中任意⼀点，其所对应的的值即能够对射线上的点进⾏排序。考虑到a、b的值可能为0，所以不适合求倒数，故采⽤[b,-a]作为特征直线族，⽤以评价点在射线上的位置，最终代码如下。
python新手能做啥兼职def  sortByRay (abc ,points ):
ab = np .array ([abc [1],-abc [0]])        #特征直线族
pDict = {ab .dot (point ):point for  point in  points }
keyList = list (pDict .keys ())            #将pDict 的兼职转化成列表
keyList .sort (reverse =True )              #对键列表进⾏排序
return  [pDict [key ] for  key in  keyList ]
这⾥⼜涉及到了⼀个新的数据类型，即字典。在理解字典之前，我们可以先回顾⼀下列表，我们可以把列表想象成⼀组值和⾃然序列的⼀⼀对应。对于列表test = [a,b,c,d]来说，有如下的对应关系{0：a,1:b,2:c,3:d}，所以我们可以通过test[0]来索引a ，test[1]来索引b ，以此类推。那么，现在我不想⽤⾃然数来索引了，我想通过⼀个标记来索引，所以希望能够创建⼀个伪列表dic = {3:5,4:15,12:22}，于
是我们可以对此列表进⾏索引dic[3]==5,dic[4]==15,dic[12]==22。这个伪列表就可以由字典来实现。这种索引关系就叫做键值对，我们通过⼀个键来索引⼀个值。对于表达式pDict = {ab.dot(point):point for point in points}，表⽰通过point 对points 进⾏遍历，即对于每个points 中的point 都进
⾏ab.dot(point)这样的点乘操作。于是得到了由特征直线族得到的特征值与点之间的⼀⼀对应关系。pDict.keys()即可提取出字典中所有的键，pDict.values()可以提取出字典中所有的值。在此我们将所有的键提取出来之后，再将其转化为列表。
然后即可调⽤列表的排序函数，将所有的值进⾏排序。即keyList.sort()，其中reverse 参数默认为False ，此时为降序。我们选择True ，此时表⽰升序。
线弧关系
⽬前，我们已经能够精确地衡量射线与线段之间的关系了，接下来需要思考如何确定射线与透镜的位置关系。这⼀点当然也要从交点说起。⾸先，弧是圆的⼀部分，所以如果⼀条直线与弧有交点，那么必然与这段弧所在的圆有交点。⽽直线与圆的交点判定相对来说还是⾮常容易的，只要圆⼼到直线的距离⼩于半径即可。
c 1c 4(x ,y )ax +by
⽽直线和圆的交点问题则可以归结为求解⽅程组：
做代换
易得其解为：
这两个公式中隐含着对线圆交点的判定信息，即根号下的表达式。当时，即表⽰直线和圆没有交点。如果我们把交点表⽰成[(x0,y0),(x1,y1)]这样的形式，那么⼀旦没有交点，则可以返回⼀个空的列表[]。程序代码如下：
# abc 为光线参数；cir 为圆参数
# point 为光源位置，当其为[]时表⽰不考虑
def  getCrossCircle (abc =[1,1,1],cir =[0,0,2],point =()):
c =np .array (cir [0:2]+[1]).dot (abc )
a2b2 = abc [0]**2+abc [1]**2
delt = a2b2*cir [2]**2-c **2
if  delt <0: return  []          #如果⽆交点，返回空列表[]
else : delt =np .sqrt (delt )      #否则求解判别式
plusMinus = lambda  a ,b : np .array (set ([a -b ,a +b ]))  #定义函数plusMinus
yCross = plusMinus (-abc [1]*c ,abc [0]*delt )/a2b2*[1,1]+cir [1]
xCross = plusMinus (-abc [0]*c ,-abc [1]*delt )/a2b2*[1,1]+cir [0]
if  point ==[]:
return  [(xCross [i ],yCross [i ]) for  i in  [0,1]]
yFlag = (yCross -point [1])*abc [0] >= 0
xFlag = (point [0]-xCross )*abc [1] >= 0
zFlag = np .abs (xFlag )+np .abs (yFlag ) > 0
flag = yFlag *xFlag *zFlag
return  [(xCross [i ],yCross [i ])
for  i in  range (len (yFlag )) if  flag [i ]]这段程序虽然短，但信息量还是很⼤的,⽽且使⽤了⼀个lambda表达式。
plusMinus = lambda a,b : np.array([a-b,a+b])定义了⼀个名为plusMius 的函数，这个函数写成常规形式即为：
def  plusMinus (a ,b )
return  np .array ([a -b ,a +b ])
需要注意的是，lambda表达式的后⾯只能有⼀个表达式，即只能定义⼀⾏的函数。在这段代码中，我们还看到了⼀个陌⽣的运算符set ，这也是python的⼀种数据类型，集合。和我们数学上认识的集合⼀样，在集合中，不允许出现相同的值。所以，如果b==0的话，那么set(a,a)=set(a)，即起到了去重的作⽤。然后再通过np.array 将集合转换成可计算的数组数据。
ax +by +c (x −x 0)+(x −y 0)22=
=0r 2
C =c +ax +0by 0
⎩⎪⎪⎪⎨⎪⎪⎪⎧y =A +B 22−BC ±A ⋅A r +B r −C 22222
x =A +B 22−AC ±B ⋅A r +B r −C 22222δ=A r +22B r −22C 2δ<0
此外，这⾥引⼊了⽐较运算符。我们⽬前所提到的运算都是数值型的，但实际上我们所接触到的运算还有其他的类型。例如，当我们进⾏判断的时候，if delt<0：中，<；也是⼀种运算符，代表⽐较，如果delt的确⼩于0，那么将返回⼀个值True ，否则⾃然返回⼀个False 。其中，True 表⽰真，False 表⽰假，这个是符合上过英语课的⼩学⽣的直觉的。类似的运算符有>,<,>=,<=,==,!=，都可以望⽂⽣义地理解，其中!=表⽰不等于。这些都是⽐较运算符，其返回值
为True 和False 。True 和False 是⼀种不同于数值的数据类型，叫布尔型。
关系运算符不仅可以作⽤于数值类型，还可以作⽤于其他数据类型，⼀般情况下的使⽤⽅法都是符合直觉的。
>>> 1==2
False
>>> [3,3]==[3,3]
True
>>> 3>3
False
>>> 3>=3
True
然后我们再来看这个算法的逻辑，由于我们求解的是直线和圆的交点，⽽真实的光线却是射线，那么必然要考虑交点和光源的位置关系。对于直线来说，令，则可以化为参数⽅程：
则对于射线上的点，必有故代码
yFlag = (yCross -point [1])*abc [0] >= 0
xFlag = (point [0]-xCross )*abc [1] >= 0分别定义了这两个判据xFlag 、yFlag 。但是当⼆者同时为0时，说明此时，此时交点即光源，故不能算作光线与圆有交点。所以⼜有判据zFlag 。
只有当这三个判据都满⾜时，我们所得到的值才是有效的，故总判据与这三个分判据是’与’的关系，所以写为flag = yFlag*xFlag*zFlag 。此外，我们并不知道交点的个数，当判别式为0的时候，lambda表达式将只有⼀个值传出，这时的xCross和yCross中分别只有⼀个元素；如果判别式⼤于0，则分别有两个元素。这⾥⼜涉及到array 的⼀个优良特性，当维度不想等的两个变量进⾏计算时，其会⾃动对低维数据进⾏合理的扩张，例如
>>> np .array ([1,2,3])+4
array ([5, 6, 7])
>>> np .array ([[1],[2],[3]])+4
array ([[5],
[6],
[7]])
最后，⼜出现了⼀个似曾相识的表达式return [(xCross[i],yCross[i]) for i in range(len(yFlag)) if flag[i]]。这个表达式可以很容易地读出来：遍历flag，如果flag的值为真，则将对应的交点坐标放⼊列表，并返回有效的交点坐标。
这是对我们之前所使⽤的[... for ... in ...]的⼀种扩张，这种写法简洁⽽强⼤，⾮常推荐使⽤。
点弧关系ax +by +c =0c =−(ax +0b 0)y {
x =x +b ⋅t
0y =y −a ⋅t
0{
(x −x )⋅b ≥0
0(y −y )⋅a ≥0
0x =x ,y =0y 0

版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。