Home

Demo

from k1lib.imports import *

2024-02-21 15:28:47,997	INFO worker.py:1458 -- Connecting to existing Ray cluster at address: 192.168.1.19:6379...
2024-02-21 15:28:48,004	INFO worker.py:1633 -- Connected to Ray cluster. View the dashboard at 127.0.0.1:8265 

This notebook serves as an introduction and main documentation for the k1lib.kop module, which mainly handles 2d optics simulation.

Intro

kop is a module within the k1lib library that's able to simulate optical phenomenons. I made it on a weekend out of frustration that most optics simulation software is really complicated, proprietary and expensive. I was still fascinated with how light behaves but without something to play with, my understanding is kinda shallow. Anyway, here're some quick examples of what you can do:

l = kop.Lens(R1=30, gp2=kop.gps["BK7"])
r = kop.Rays.parallel()
rp = l.cast(r); kop.Drawables([l, rp])

<Drawables #n=2>

l1 = kop.Lens(x=50, R1=30, gp2=kop.gps["BK7"])
l2 = kop.Lens(x=80, R1=-30, gp2=kop.gps["sapphire"])
r = kop.Rays.parallel(color="rainbow")
opsys = kop.OpticSystem([l1, l2])
rp = opsys.cast(r); kop.Drawables([l1, l2, rp])

<Drawables #n=3>

# thumbnail
r = kop.Rays.parallel(y=-4, height=13, theta=0.04, color="rainbow")
l1 = kop.Lens(x=50, R1=20, R2=300, thickness=5, gp2=kop.gps["BK7"], capture=True, surface="parabola")
l2 = kop.LineSurface(80, -20, 80, 20, gp2=kop.gps["air"], capture=True)
l3 = kop.Lens(x=100, R1=20, R2=300, thickness=5, gp2=kop.gps["BK7"])
l4 = kop.Polygon([[120, -20], [130, 20], [140, -20]], gp2=kop.gps["BK7"], capture=True, mode="mirror")
opsys = kop.OpticSystem([l1, l2, l3, l4]); rp = opsys.cast(r); kop.Drawables([l1, l2, l3, l4, rp]).img() # the .img() is not necessary btw, it just tells the system to return a png instead of the typical svg

settings.kop.display.rays.infLength = 30
l1 = kop.Lens(R1=-20, R2=20, thickness=5, gp2=kop.gps["BK7"], capture=True);
r = kop.Rays.pointToBounds(58, 50, l1.bounds(), color="rainbow", N=10, coverage=0.3)
opsys = kop.OpticSystem().add(l1); rp = opsys.cast(r, verbose=True); kop.Drawables([opsys, rp])

<Lens R1=300 R2=300 thickness=3 height=20>
0) <Line(Surface) (47.5, 10.0, 52.5, 10.0) gp1='air' gp2='BK7' capture=True>
1) <Arc(Surface) x=30.18 y=-0.0 r=20 startAngle(deg)=330.0 endAngle(deg)=390.0 gp1='BK7' gp2='air' capture=True>
2) <Line(Surface) (52.5, -10.0, 47.5, -10.0) gp1='air' gp2='BK7' capture=True>
3) <Arc(Surface) x=35.18 y=-0.0 r=20 startAngle(deg)=330.0 endAngle(deg)=390.0 gp1='air' gp2='BK7' capture=True>


Pass: 0; chosen surfaces: 0   0   0   0   0   0   0   3   3   3   
Pass: 1; chosen surfaces: 1   1   1   2   2   2   2   1   1   1   
Pass: 2; chosen surfaces: 3   3   2   0   0   0   0   3   3   3   
Pass: 3; chosen surfaces: 1   2   0   0   0   0   0   0   0   0   
Pass: 4; chosen surfaces: 2   1   0   0   0   0   0   0   0   0   
Pass: 5; chosen surfaces: 0   0   0   0   0   0   0   0   0   0   

<Drawables #n=2>

l1 = kop.Lens(R1=20, R2=30, thickness=5, gp2=kop.gps["BK7"])
l2 = kop.Lens(R1=-20, R2=30, gp2=kop.gps["BK7"], x=120, y=-3, angle=3*pi/180, mode="mirror")
l3 = kop.LineSurface(30, -20, 100, -20, gp2=kop.gps["air"], mode="mirror")
r = kop.Rays.pointToBounds(0, 0, l1.bounds(), color="rainbow", coverage=0.35, N=30)
opsys = kop.OpticSystem([l1, l2, l3]); rp = opsys.cast(r)
kop.Drawables([opsys, rp])

<Drawables #n=2>

settings.kop.display.rays.infLength = 100

So yeah, it's pretty dang powerful.

There are only a few fundamental classes and concepts you need to know before fully utilizing this:

• Rays: a collection of rays ready to interact with Surfaces
• Surface (LineSurface, ArcSurface, PolySurface): a surface that has the ability to interact and bend Rays
• RaysPath: a collection of Rays at each step after it got refracted or reflected
• OpticElement (Lens, Polygon): a collection of Surfaces encapsulated in a convenient package
• OpticSystem: a collection of OpticElements and raw Surfaces

Rays

The fundamental data structure behind this is a numpy array of shape (N, 8), where N is the number of rays total:

r = kop.Rays.parallel(); r

<Rays N=10 avgTheta=0.0deg intersectedSurface?=False ogSurface=None>

r.data | shape() | aS(print); r.data

(10, 8)

array([[ 3.06161700e-16,  5.00000000e+00,  0.00000000e+00,
                    inf,  7.00000000e+02,  0.00000000e+00,
         1.00000000e+00,  0.00000000e+00],
       [ 3.06161700e-16,  3.88888889e+00,  0.00000000e+00,
                    inf,  7.00000000e+02,  0.00000000e+00,
         1.00000000e+00,  0.00000000e+00],
       [ 3.06161700e-16,  2.77777778e+00,  0.00000000e+00,
                    inf,  7.00000000e+02,  0.00000000e+00,
         1.00000000e+00,  0.00000000e+00],
       [ 3.06161700e-16,  1.66666667e+00,  0.00000000e+00,
                    inf,  7.00000000e+02,  0.00000000e+00,
         1.00000000e+00,  0.00000000e+00],
       [ 3.06161700e-16,  5.55555556e-01,  0.00000000e+00,
                    inf,  7.00000000e+02,  0.00000000e+00,
         1.00000000e+00,  0.00000000e+00],
       [ 3.06161700e-16, -5.55555556e-01,  0.00000000e+00,
                    inf,  7.00000000e+02,  0.00000000e+00,
         1.00000000e+00,  0.00000000e+00],
       [ 3.06161700e-16, -1.66666667e+00,  0.00000000e+00,
                    inf,  7.00000000e+02,  0.00000000e+00,
         1.00000000e+00,  0.00000000e+00],
       [ 3.06161700e-16, -2.77777778e+00,  0.00000000e+00,
                    inf,  7.00000000e+02,  0.00000000e+00,
         1.00000000e+00,  0.00000000e+00],
       [ 3.06161700e-16, -3.88888889e+00,  0.00000000e+00,
                    inf,  7.00000000e+02,  0.00000000e+00,
         1.00000000e+00,  0.00000000e+00],
       [ 3.06161700e-16, -5.00000000e+00,  0.00000000e+00,
                    inf,  7.00000000e+02,  0.00000000e+00,
         1.00000000e+00,  0.00000000e+00]])

This encapsulates all the information about Rays, and has these columns: x, y, theta, length, wavelength (nm), #transforms, power (Watts), intersected?

However, you only need to pay attention to x, y, theta and wavelength. x and y specifies where the Rays are in space, and theta is the angle of the ray going counterclockwise from positive x axis. It's possible to construct any type of Rays that you want, just pass in your custom numpy array, but I've made several convenience methods to construct them:

kop.Rays.parallel(x=10, y=20, theta=pi/16, N=10, height=10, color="rainbow")

<Rays N=10 avgTheta=11.25deg intersectedSurface?=False ogSurface=None>

kop.Rays.parallelToBounds(x=-2, y=-1, bounds=[10, 0, 10, 1], N=10, color="rainbow")

<Rays N=10 avgTheta=7.13deg intersectedSurface?=False ogSurface=None>

kop.Rays.pointToBounds(x=-2, y=-1, bounds=[10, 0, 10, 1], N=10, color="rainbow")

<Rays N=10 avgTheta=7.11deg intersectedSurface?=False ogSurface=None>

Btw, before we move on, these sketches of the Rays and whatnot are designed to be informative, meaning it will try to cram a lot of info in. This might not be appropriate for presenting to another person, so you can tweak some of the display settings available:

settings.kop

Settings:                                                        
- colorD  = {'rainbow': [400, 650], 'red': [620, 750], 'ora...   ​color wavelength ranges to be used in constructing Rays                                                                                                                            
- gps     = {'BK7': (1.03961212, 0.231792344, 1.01046945, 0...   ​All builtin glass parameters of the system. All have 6 floats, for B1,B2,B3,C1,C2,C3 parameters used in the sellmeier equation: https://en.wikipedia.org/wiki/Sellmeier_equation   
- display = <Settings>                                           ​display settings                                                                                                                                                                   
  - drawable = <Settings>                                        ​generic draw settings                                                                                                                                                              
    - axes      = True                                           ​whether to add x and y axes to the sketch                                                                                                                                          
    - maxWh     = 800                                            ​when drawing a sketch, it will be rescaled so that the maximum of width and height of the final image is this number. Increase to make the sketch bigger                           
    - grid      = True                                           ​whether to add grid lines to the sketch                                                                                                                                            
    - gridColor = (255, 255, 255)                                ​                                                                                                                                                                                   
  - rays     = <Settings>                                        ​display settings of kop.Rays                                                                                                                                                       
    - showOrigin = True                                          ​whether to add a small red dot to the beginning of the mean (x,y) of a Rays or not                                                                                                 
    - infLength  = 100                                           ​length in mm to display Rays if their length is infinite                                                                                                                           
  - surface  = <Settings>                                        ​display settings for kop.Surface class                                                                                                                                             
    - showIndex = True                                           ​whether to show the index of the Surface in an OpticSystem or not                                                                                                                  
- consts  = <Settings>                                           ​magic constants throughout the sim. By default works pretty well, but you can tweak these if you need unrealistic setups, like super big focal length, etc                         
  - inchForward = 1e-06                                          ​after new Rays have been built, inch forward the origin of the new Rays by a this tiny amount so that it 'clears' the last Surface                                                 
                                                                 

To set one of them, just do something like settings.kop.display.drawable.grid = False in order to not display the grid.

Surface

A Surface is some sort of line, arc, parabola or any other continuous curves. Available types include:

• LineSurface: a Surface that's just a simple line between 2 points
• ArcSurface: a Surface that's made of a circular arc
• PolySurface: a Surface defined by any random polynomial

Let's see some examples:

kop.LineSurface(1, 1, 10, 3)

<Line(Surface) (1, 1, 10, 3) gp1=None gp2=None>

Those "gp1" and "gp2" you see in the repr() of a Surface above means "glass param". There're several default materials you can choose from:

kop.gps

{'BK7': (1.03961212,
  0.231792344,
  1.01046945,
  0.00600069867,
  0.0200179144,
  103.560653),
 'sapphire': (1.4313493,
  0.65054713,
  5.3414021,
  0.0052799261,
  0.0142382647,
  325.017834),
 'fusedSilica': (0.6961663,
  0.4079426,
  0.8974794,
  0.00467914826,
  0.0135120631,
  97.9340025),
 'MgF': (0.48755108,
  0.39875031,
  2.3120353,
  0.001882178,
  0.008951888,
  566.13559),
 'air': (0, 0, 0, 0, 0, 0)}

It's also available under settings.kop.gps These 6 numbers are the B1,B2,B3,C1,C2,C3 parameters used in the Sellmeier equation: https://en.wikipedia.org/wiki/Sellmeier_equation. You can add your custom materials here. In this LineSurface, the first point is on the left (1, 1) and second point is on the right (10, 3). Then the space on top is considered to have properties of gp1, and below is gp2. gp1 is usually that of the environment, and gp2 is usually that of the object.

If you don't have Sellmeier parameters for your material and just want a flat index of refraction for all incoming wavelengths, you can pass in that number in gp1 or gp2 directly.

Each surface has a function .castS() ("s" for single), to try to calculate the resulting Rays:

l = kop.LineSurface(30, -20, 35, 20, gp1=kop.gps["air"], gp2=kop.gps["BK7"]); l
r = kop.Rays.parallel(theta=pi/16)
r1 = l.castS(r); r1

<Rays N=10 avgTheta=364.9deg intersectedSurface?=True ogSurface=<Line(Surface) (30, -20, 35, 20) gp1='air' gp2='BK7'>>

r1.data[:,:3]

array([[33.98217994, 11.85743154,  6.36871009],
       [33.83695928, 10.69566629,  6.36871009],
       [33.69173863,  9.53390104,  6.36871009],
       [33.54651797,  8.37213579,  6.36871009],
       [33.40129731,  7.21037054,  6.36871009],
       [33.25607666,  6.04860529,  6.36871009],
       [33.110856  ,  4.88684004,  6.36871009],
       [32.96563534,  3.72507479,  6.36871009],
       [32.82041469,  2.56330954,  6.36871009],
       [32.67519403,  1.4015443 ,  6.36871009]])

The .castS() function is pretty fundamental and low level. It just calculates the resulting Rays from interacting with this singular Surface. Then, you can draw everything out (original ray, outgoing ray, surface) like this:

kop.Drawables([r, r1, l])

<Drawables #n=3>

Here, gp1 is on the left (air), and gp2 is on the right (borosilicate glass). The rays going through will be slowed down, and it looks like it bent by the right amount. What happens if I were to switch gp1 and gp2?

l = kop.LineSurface(30, -20, 35, 20, gp1=kop.gps["BK7"], gp2=kop.gps["air"]); l
r = kop.Rays.parallel(theta=pi/16)
r1 = l.castS(r); kop.Drawables([r, r1, l])

<Drawables #n=3>

Now it's travelling from BK7 to air, which widens the angle of refraction. What happens if I were to switch the first and second point of the LineSurface?

l = kop.LineSurface(35, 20, 30, -20, gp1=kop.gps["BK7"], gp2=kop.gps["air"]); l
r = kop.Rays.parallel(theta=pi/16)
r1 = l.castS(r); kop.Drawables([r, r1, l])

<Drawables #n=3>

It returns back to the original setup. Now gp1 is on the right (BK7), and gp2 is on the left (air). Just for funsies, how about let's make the 2 gps the same?

l = kop.LineSurface(35, 20, 30, -20, gp1=kop.gps["air"], gp2=kop.gps["air"]); l
r = kop.Rays.parallel(theta=pi/16)
r1 = l.castS(r); kop.Drawables([r, r1, l])

<Drawables #n=3>

It passes through without any refraction at all, which makes sense. What happens if the surface don't cover all of the Rays?

l = kop.LineSurface(30, -20, 35, 8, gp1=kop.gps["air"], gp2=kop.gps["BK7"]); l
r = kop.Rays.parallel(theta=pi/16)
r1 = l.castS(r); kop.Drawables([r, r1, l])

<Drawables #n=3>

Well then, only the rays that interact with the Surface will change anything. Let's compare the underlying numpy array for both rays:

r.data[:,:3]

array([[-0.97545161,  4.9039264 ,  0.19634954],
       [-0.75868459,  3.81416498,  0.19634954],
       [-0.54191756,  2.72440356,  0.19634954],
       [-0.32515054,  1.63464213,  0.19634954],
       [-0.10838351,  0.54488071,  0.19634954],
       [ 0.10838351, -0.54488071,  0.19634954],
       [ 0.32515054, -1.63464213,  0.19634954],
       [ 0.54191756, -2.72440356,  0.19634954],
       [ 0.75868459, -3.81416498,  0.19634954],
       [ 0.97545161, -4.9039264 ,  0.19634954]])

r1.data[:,:3]

array([[-0.97545161,  4.9039264 ,  0.19634954],
       [-0.75868459,  3.81416498,  0.19634954],
       [-0.54191756,  2.72440356,  0.19634954],
       [-0.32515054,  1.63464213,  0.19634954],
       [34.91268015,  7.51100326,  6.34974775],
       [34.70292998,  6.33640229,  6.34974775],
       [34.4931798 ,  5.16180132,  6.34974775],
       [34.28342963,  3.98720036,  6.34974775],
       [34.07367946,  2.81259939,  6.34974775],
       [33.86392928,  1.63799842,  6.34974775]])

So the last 6 rays does intersect with the Surface, so it changes dramatically, but the first 4 rays just copies that of the previous Rays. A higher level API for casting rays is .cast(), instead of .castS(). This will internally construct an OpticSystem, and will call .castS() a lot of times, until the rays don't change any more.

l = kop.LineSurface(30, -20, 35, 8, gp1=kop.gps["air"], gp2=kop.gps["BK7"]); l
r = kop.Rays.parallel(theta=pi/16)
rp = l.cast(r); rp

<RaysPath #rays=2>

.cast() will return RaysPath, which is just a collection of Rays packaged together so that it's convenient to plot. For a single Surface, you might be like "what's the point", but when dealing with multiple different Surfaces, .cast() will handle interactions with all Surfaces at the same time, instead of a single Surface that you choose. But I'm getting ahead of myself. You can also extract out individual Rays from a RaysPath:

rp[0]

<Rays N=10 avgTheta=11.25deg intersectedSurface?=False ogSurface=None>

rp[1]

<Rays N=10 avgTheta=222.79deg intersectedSurface?=True ogSurface=<Line(Surface) (30, -20, 35, 8) gp1='air' gp2='BK7'>>

By default, all Surfaces are transparent mediums. But sometimes you want a mirror instead. Then you can just switch the mode of it like this:

l = kop.LineSurface(30, -20, 35, 8, gp1=kop.gps["air"], gp2=kop.gps["BK7"], mode="mirror"); l
r = kop.Rays.parallel(theta=pi/16)
rp = l.cast(r); kop.Drawables([rp, l])

<Drawables #n=2>

Pretty neat. ArcSurface is another Surface type that's made from a segment of a circle:

kop.ArcSurface(x=1, y=1, r=30, startAngle=pi/4, endAngle=pi/2+pi/4)

<Arc(Surface) x=1 y=1 r=30 startAngle(deg)=45.0 endAngle(deg)=135.0 gp1=None gp2=None>

This is the standard definition of an ArcSurface. But the parameters required might not be convenient, so there're other methods to construct an ArcSurface:

kop.ArcSurface.from2Points(-20, 22, 20, 22, r=30)

<Arc(Surface) x=0.0 y=-0.36 r=30 startAngle(deg)=48.19 endAngle(deg)=131.81 gp1=None gp2=None>

Like LineSurface, if the first point is on the left, second point is on the right, then the top will be gp1, and bottom will be gp2. Using this, you can construct a basic lens:

a = kop.ArcSurface.from2Points(30, -20, 30, 20, r=30, gp1=kop.gps["air"], gp2=kop.gps["BK7"])
r = kop.Rays.parallel(); r1 = a.castS(r); r1

<Rays N=10 avgTheta=360.0deg intersectedSurface?=True ogSurface=<Arc(Surface) x=52.36 y=0.0 r=30 startAngle(deg)=138.19 endAngle(deg)=221.81 gp1='air' gp2='BK7'>>

kop.Drawables([r, r1, a])

<Drawables #n=3>

Here, gp1 is on the left (air), gp2 is on the right (BK7), so it focuses the light down to a single point. If I were to swap the materials around, the rays diverges instead, which makes sense:

a = kop.ArcSurface.from2Points(30, -20, 30, 20, r=30, gp1=kop.gps["BK7"], gp2=kop.gps["air"])
r = kop.Rays.parallel(); r1 = a.castS(r); kop.Drawables([r, r1, a])

<Drawables #n=3>

Finally, there's the PolySurface, which you can specify any polynomial:

kop.PolySurface(x=20, y=1, angle=-pi/4, scale=10, xmin=-1, xmax=2, coeff=(0, 0, 1))

<Poly(Surface) coeff=(0, 0, 1) x=20 y=1 angle(deg)=-45.0 scale=10 xmin=-1 xmax=2 gp1=None gp2=None>

Honestly it just looks like a Nike logo lmao. The above command creates a Surface with the equation 0 + 0*x + 1*x^2, with the polynomial's domain from -1 to 2. Then everything got rotated, scaled and translated so that it looks like the picture. Also notice how the vertex of the parabola is at (20, 1), right where we want it to be.

Just like the 2 surfaces before, you can cast Rays with this:

p = kop.PolySurface(x=20, y=-6, angle=-pi/4, scale=10, xmin=-1, xmax=2, coeff=(0, 0, 1), gp1=kop.gps["air"], gp2=kop.gps["BK7"])
r = kop.Rays.parallel(); r1 = p.castS(r); kop.Drawables([r, r1, p])

<Drawables #n=3>

Note that the above graph is not quite correct. If you use .cast(), you will get a different answer:

kop.Drawables([p.cast(r, verbose=True), p])

Raw surfaces:
0) <Poly(Surface) coeff=(0, 0, 1) x=20 y=-6 angle(deg)=-45.0 scale=10 xmin=-1 xmax=2 gp1='air' gp2='BK7'>


Pass: 0; chosen surfaces: 0   0   0   0   0   0   0   0   0   0   
Pass: 1; chosen surfaces: 0   0   0   0   0   0   0   0   0   0   
Pass: 2; chosen surfaces: 0   0   0   0   0   0   0   0   0   0   

<Drawables #n=2>

Since .castS() only tries to ray trace it once, it can't quite see the curve on the right. If you were to do .castS() twice, it would have figure it out though:

r1 = p.castS(r); r2 = p.castS(r1)
kop.Drawables([r, r1, r2, p])

<Drawables #n=4>

OpticElement

This is a collection of Surfaces that's encapsulated in a convenient package. Only 2 elements are available right now: Lens and Polygon.

kop.Lens(R1=30, angle=-pi/2)

<Lens R1=30 R2=300 thickness=3 height=20>

This doesn't have the simple single step cast function .castS(). It only has .cast():

l = kop.Lens(R1=30, gp2=kop.gps["BK7"]); r = kop.Rays.parallel()
kop.Drawables([l.cast(r), l])

<Drawables #n=2>

You can make it a little more absurd:

l = kop.Lens(R1=20, R2=-17, thickness=10, gp2=kop.gps["BK7"]); r = kop.Rays.parallel(height=15)
kop.Drawables([l.cast(r), l])

<Drawables #n=2>

The other is a Polygon, which you can make a prism out of:

p = kop.Polygon([[10, -10], [50, -10], [30, 10]]); p

<Polygon(OpticElement) #points=3 position=(30.0, -3.33)

p = kop.Polygon([[10, -10], [50, -10], [30, 10]], gp2=kop.gps["BK7"])
r = kop.Rays.parallel(); kop.Drawables([p.cast(r), p])

<Drawables #n=2>

Now isn't that beautiful? This looks a little malicious to me though, as it looks like how governments would tap into submarine fiber optic cables and whatnot. But its interesting how the output rays are perfectly parallel to the input rays.

OpticSystem

This is a collection of OpticElements and Surfaces all in a single package. This allows you to define multiple objects and have them all interact with each other:

l1 = kop.Lens(x=30, R1=30, gp2=kop.gps["BK7"])
l2 = kop.Lens(x=60, R1=-40, gp2=kop.gps["BK7"])
s = kop.LineSurface(50, -10, 48, 0, gp2=kop.gps["BK7"], mode="mirror")
r = kop.Rays.parallel()
opsys = kop.OpticSystem([l1, l2, s]); opsys

<OpticSystem #opticElements=2 #surfaces=9><Lens R1=300 R2=300 thickness=3 height=20>
0) <Line(Surface) (30.22, 10.0, 31.33, 10.0) gp1=None gp2='BK7'>
1) <Arc(Surface) x=58.5 y=0.0 r=30 startAngle(deg)=160.53 endAngle(deg)=199.47 gp1=None gp2='BK7'>
2) <Line(Surface) (31.33, -10.0, 30.22, -10.0) gp1=None gp2='BK7'>
3) <Arc(Surface) x=-268.5 y=-0.0 r=300 startAngle(deg)=358.09 endAngle(deg)=361.91 gp1=None gp2='BK7'>

<Lens R1=300 R2=300 thickness=3 height=20>
4) <Line(Surface) (58.5, 10.0, 61.5, 10.0) gp1=None gp2='BK7'>
5) <Arc(Surface) x=19.77 y=-0.0 r=40 startAngle(deg)=345.52 endAngle(deg)=374.48 gp1='BK7' gp2=None>
6) <Line(Surface) (61.5, -10.0, 58.5, -10.0) gp1=None gp2='BK7'>
7) <Arc(Surface) x=-238.33 y=-0.0 r=300 startAngle(deg)=358.09 endAngle(deg)=361.91 gp1=None gp2='BK7'>

Raw surfaces:
8) <Line(Surface) (50, -10, 48, 0) gp1=None gp2='BK7'>

Once you have created an OpticSystem, you can grab any surfaces you want:

opsys[1]

<Arc(Surface) x=58.5 y=0.0 r=30 startAngle(deg)=160.53 endAngle(deg)=199.47 gp1=None gp2='BK7'>

Then you can cast any Rays you want to get a RaysPath:

settings.kop.display.rays.infLength = 30
kop.Drawables([opsys, opsys.cast(r)])

<Drawables #n=2>

There's also a verbose mode that's pretty helpful with debugging things and is a pretty nice metrology tool:

kop.Drawables([opsys, opsys.cast(r, verbose=True)])

<Lens R1=300 R2=300 thickness=3 height=20>
0) <Line(Surface) (30.22, 10.0, 31.33, 10.0) gp1='air' gp2='BK7'>
1) <Arc(Surface) x=58.5 y=0.0 r=30 startAngle(deg)=160.53 endAngle(deg)=199.47 gp1='air' gp2='BK7'>
2) <Line(Surface) (31.33, -10.0, 30.22, -10.0) gp1='air' gp2='BK7'>
3) <Arc(Surface) x=-268.5 y=-0.0 r=300 startAngle(deg)=358.09 endAngle(deg)=361.91 gp1='air' gp2='BK7'>

<Lens R1=300 R2=300 thickness=3 height=20>
4) <Line(Surface) (58.5, 10.0, 61.5, 10.0) gp1='air' gp2='BK7'>
5) <Arc(Surface) x=19.77 y=-0.0 r=40 startAngle(deg)=345.52 endAngle(deg)=374.48 gp1='BK7' gp2='air'>
6) <Line(Surface) (61.5, -10.0, 58.5, -10.0) gp1='air' gp2='BK7'>
7) <Arc(Surface) x=-238.33 y=-0.0 r=300 startAngle(deg)=358.09 endAngle(deg)=361.91 gp1='air' gp2='BK7'>

Raw surfaces:
8) <Line(Surface) (50, -10, 48, 0) gp1='air' gp2='BK7'>


Pass: 0; chosen surfaces: 1   1   1   1   1   1   1   1   1   1   
Pass: 1; chosen surfaces: 3   3   3   3   3   3   3   3   3   3   
Pass: 2; chosen surfaces: 5   5   5   5   5   8   8   8   8   8   
Pass: 3; chosen surfaces: 7   7   7   7   7   3   3   3   3   3   
Pass: 4; chosen surfaces: 0   0   0   0   0   1   1   1   1   1   
Pass: 5; chosen surfaces: 0   0   0   0   0   0   0   0   0   0   

<Drawables #n=2>

You can clearly see what ray chose what surface on each pass, and you can also see the surface index numbers in the sketch as well.

Color

I just want to touch on this topic a little bit. Each ray can have a different wavelength, and there're convenience methods to construct rays with a particular wavelength:

settings.kop.display.rays.infLength = 100
kop.Rays.parallel(color="rainbow")

<Rays N=10 avgTheta=0.0deg intersectedSurface?=False ogSurface=None>

kop.Rays.parallel(color="yellow")

<Rays N=10 avgTheta=0.0deg intersectedSurface?=False ogSurface=None>

These color codes are configurable here:

settings.kop.colorD

{'rainbow': [400, 650],
 'red': [620, 750],
 'orange': [590, 620],
 'yellow': [570, 590],
 'green': [495, 570],
 'blue': [450, 495],
 'purple': [400, 450]}

These are the wavelength ranges for each color code. Speaking of wavelength:

kop.Rays.parallel(color=500) # 500nm wavelength

<Rays N=10 avgTheta=0.0deg intersectedSurface?=False ogSurface=None>

kop.Rays.parallel(color=[500, 600]) # 500nm to 600nm wavelengths, spreaded out

<Rays N=10 avgTheta=0.0deg intersectedSurface?=False ogSurface=None>

kop.Rays.parallel(color=[600, 500]) # reverse of that 

<Rays N=10 avgTheta=0.0deg intersectedSurface?=False ogSurface=None>

Pretty nice. As the index of refraction depends on the wavelength, lets see that classic prism setup that splits white light apart:

settings.kop.display.rays.infLength = 30
p = kop.Polygon([[10, -10], [30, -10], [20, 10]], gp2=kop.gps["BK7"])
r = kop.Rays.parallel(height=0, color="rainbow"); kop.Drawables([p.cast(r), p])

<Drawables #n=2>

What's the angle spread of the outward rays?

angles = p.cast(r)[-1].data[:,2]; angles

array([-0.66288376, -0.65645809, -0.65128918, -0.6470508 , -0.6435176 ,
       -0.64052918, -0.63796854, -0.63574849, -0.63380284, -0.63208056])

angles | toMin() & toMax() | ~aS(lambda x,y: y-x) | op()*180/pi

1.7648932897648684

1.76 degrees, that's not too big, smaller than what textbook diagrams usually depict. Gotta say, I'm loving the versatility of this thing.