Simulating a mirror

[Rod]

Ok, my current project is light and lenses. As usual my maths is letting me down. I want a mirror placed on screen that will properly reflect light rays shone on it. I have a sort of solution for a ray coming right to left but not a foolproof solution. Any help or pointers welcome.

I will post the project as it is now (broken) so you can get a feel for where I am headed.

Here is the link to the project code and .zip




rX=current x position of ray
rY=current y position of ray
'rXD,rYD define the normalised vector of the ray between -1 to 1

'"mirr" is a rotatable mirror on screen which has one face.
'oA defines the orientation of the mirror.
'at 0o it faces left at 45o left and up at 90o up.
'the mirror is restricted to eight orientations.

'I move the ray thus
rX=rX+rXD
rY=rY+rYD

'If I hit a mirror I need to reflect correctly
'or stop if it is the back of the mirror
case "mirr"
   if oA=0 then rXD=rXD*-1
   if oA=45 thenrXD=0 : rYD=-1
   if oA=90 then nul
   if oA=135 then stopped=1
   if oA=180 then stopped=1
   if oA=225 then stopped=1
   if oA=270 then rXD=0 : rYD=-1
   if oA=315 then rXD=0 : rYD=1

'The above sort of works for an incoming ray from left to right


[Rod]

I think my answer is here, no atan2 just vectors. But how to code?

stackoverflow.com/questions/30970103/2d-line-reflection-on-a-mirror

[tsh73]

I used vector library I have. It reflects at any angle
(and uses dot product to check side of the mirror, had to google - it was really long time ago)

'Mirror: (tsh73 Jun 2023)
'make ray reflect off mirror
'1. make that mirror
'2. rotate
'3. reflect center
'4. simplify direction check via dot product
' + show mirror angle
'*****************************************************
'nomainwin

global pi
pi = acs(-1)   'global for atan2, dePi
pi2=pi/2

desiredWidth = 400
desiredHeight = 400
mirrR=50
mirrA=0 '0 is left, 90 up
mirrAA=mirrA*pi/180-pi2 'same in radians, rotated as needed

gosub [ajustWindow]
'now, center window
   UpperLeftX = 1'(DisplayWidth - WindowWidth)/2
   UpperLeftY = 1'(DisplayHeight - WindowHeight)/2

open "Mirror" for graphics_nsb_nf as #gr
#gr "trapclose [quit]"
#gr "down; fill white; flush"
#gr "home; posxy cx cy"
#gr "setfocus"
#gr "when characterInput [key]"

gosub [msg]
gosub [drawMirror]

#gr "when leftButtonDown [start]"
#gr "when leftButtonMove [moving]"
#gr "when leftButtonUp [end]"
#gr "when rightButtonUp [ray]"

wait

[start]
mx=MouseX:my=MouseY
'print mx, my
if sqr((mx-cx)^2+(my-cy)^2)>mirrR then wait
''or else
started=1
mirrAA0=mirrAA
a0=vectAngle(vect$(mx-cx, my-cy))
'print "started", a0
wait

[moving]
if not(started) then wait
mx=MouseX:my=MouseY
a1=vectAngle(vect$(mx-cx, my-cy))
'print "moving", a1
mirrAA=mirrAA0+a1-a0
gosub [drawMirror]
wait

[end]
if not(started) then wait
started=0
'mirrAA=mirrA*pi/180-pi2
mirrA=(mirrAA+pi2)*180/pi
#gr "place ";2*cx-70;" 20"
#gr "\";using("###.#",mirrA);"   "
wait

[ray]
mx=MouseX:my=MouseY
'ray from mx, my to cx, cy
#gr "place ";mx;" ";my
#gr "goto "; base$
'now check if ray falls on right side of the mirror
ray$=vectSub$(base$, vect$(mx,my))
'if right side of the mirror - reflect
if vectDotProduct(ray$, mirrNormal$) <0 then
rayT$=vectTangent$(ray$,mirrVect$)
rayN$=vectNorm$(ray$,mirrVect$)
' #gr "color blue"
' #gr "line ";base$;" "; vectAdd$(base$, vectScale$(0.1,rayT$))
' #gr "color red"
' #gr "line ";base$;" "; vectAdd$(base$, vectScale$(0.1,rayN$))
rayRefl$=vectSub$(rayT$, rayN$)
#gr "color green"
#gr "line ";base$;" "; vectAdd$(base$, rayRefl$)
end if
#gr "color black"
wait

[drawMirror]
base$=vect$(cx, cy)
#gr "place ";base$
#gr "circlefilled ";mirrR 'clear
#gr "circle ";mirrR
mirrPt1$=vect$(cx+mirrR*cos(mirrAA),cx+mirrR*sin(mirrAA))
mirrPt2$=vect$(cx+mirrR*cos(mirrAA+pi),cx+mirrR*sin(mirrAA+pi))
#gr "line ";mirrPt1$;" ";mirrPt2$
'#gr "size 3; set ";mirrPt2$ 'mark endpoint
'#gr "size 1"
'get and draw a normal
mirrVect$=vectSub$(mirrPt2$,mirrPt1$)
mirrNormal$=vectScale$(0.1,vectRotate$(mirrVect$,pi/2))
#gr "place ";base$
#gr "goto "; vectAdd$(base$, mirrNormal$)
return

[msg]
#gr "place 20 20"
#gr "\"
#gr "\LMB drag in a circle rotates mirror"
#gr "\RMB shoots ray"
#gr "\ESC clears screen"
return

[key]
if Inkey$=chr$(27) then
#gr "cls"
gosub [msg]
gosub [drawMirror]
end if
   wait

[quit]
   close #gr
   end

return
'-------------------------------------------------
[ajustWindow]
   UpperLeftX = 20
   UpperLeftY = 20
   WindowWidth = 200   '100 seems to be too much - works different
   WindowHeight = 100
'    MENU #gr, "dummy"
   open "Ajusting..." for graphics_nsb_nf as #gr
'   graphics
'   graphics_nsb
'   graphics_nsb_nf

   #gr, "home ; down ; posxy x y"
   'x, y give us width, height
   width = 2*x : height = 2*y
   close #gr

   slackX = 200-width
   slackY = 100-height

   WindowWidth = desiredWidth + slackX
   WindowHeight = desiredHeight + slackY

return
'------------------------------
'GetPixelValue$ returns a string with the RGB values of the pixel
'in coordinates x and y in window/graphicbox names handle$ (e.g, "#main.graph")
function GetPixelValue$(x, y, handle$)

'Grab a 1*1 bitmap
   #handle$, "getbmp gpv "; x; " "; y; " "; 1; " "; 1

'Save in a bmp file
   bmpsave "gpv", "getpvaluetemp.bmp"

'Open the file for string input and get it's full contents
   open "getpvaluetemp.bmp" for input as #gpv
   s$ = input$(#gpv, lof(#gpv))
   close #gpv

'Check if user's display is 32-bit, and read the red-green-blue values
'If display 16 bit, then colors are masked. So some last (3 for red, 2 for green, 3 for blue) bits always 0
'That means that you did not get 255 255 255 for white - (248 252 248) instead. You have to experiment
'otherwise function returns nothing (support for other display types could be added (?))
   bpp =  asc(mid$(s$, 29, 1))
   select case bpp
   case 32
       red = asc(mid$(s$, 69, 1))
       green = asc(mid$(s$, 68, 1))
       blue = asc(mid$(s$, 67, 1))
  case 24  'from Richard, for LBB
       red = asc(mid$(s$, 57, 1))
       green = asc(mid$(s$, 56, 1))
       blue = asc(mid$(s$, 55, 1))
   case 16
       bytes = asc(mid$( s$, 67, 1)) + 256*asc(mid$( s$, 68, 1))
       red =  (bytes AND 63488) /256       '0xF800
       green =  (bytes AND 2016) / 32 * 4  '0x7E0
       blue =  (bytes AND 31) * 8          '0x1F
   end select

'concatenate the return value, delete temporary file and free memory
   GetPixelValue$ = using("###",red)+using("####",green)+using("####",blue)
   kill "getpvaluetemp.bmp"
   unloadbmp "gpv"
end function

'===================================
function vect$(x,y)
   vect$=x;" ";y
end function

function vectX(v$)
   vectX=val(word$(v$,1))
end function

function vectY(v$)
   vectY=val(word$(v$,2))
end function

function vectLen(v$)
   x=val(word$(v$,1))
   y=val(word$(v$,2))
   vectLen=sqr(x*x+y*y)
end function

function vectUnit$(v$)
   x=val(word$(v$,1))
   y=val(word$(v$,2))
   vectLen=sqr(x*x+y*y)
   vectUnit$=x/vectLen;" ";y/vectLen
end function

function vectAdd$(v1$,v2$)
   x1=val(word$(v1$,1))
   y1=val(word$(v1$,2))
   x2=val(word$(v2$,1))
   y2=val(word$(v2$,2))
   vectAdd$=x1+x2;" ";y1+y2
end function

function vectSub$(v1$,v2$)
   x1=val(word$(v1$,1))
   y1=val(word$(v1$,2))
   x2=val(word$(v2$,1))
   y2=val(word$(v2$,2))
   vectSub$=x1-x2;" ";y1-y2
end function

function vectDotProduct(v1$,v2$)
   x1=val(word$(v1$,1))
   y1=val(word$(v1$,2))
   x2=val(word$(v2$,1))
   y2=val(word$(v2$,2))
   vectDotProduct=x1*x2+y1*y2
end function

function vectScale$(a,v$)   'a * vector v$
   x=val(word$(v$,1))
   y=val(word$(v$,2))
   vectScale$=a*x;" ";a*y
end function

function vectTangent$(v$,base$)
   n$=vectUnit$(base$)
   vectTangent$=vectScale$(vectDotProduct(n$,v$),n$)
end function

function vectNorm$(v$,base$)
   vectNorm$=vectSub$(v$,vectTangent$(v$,base$))
end function

function vectAngle(v$)
   x=val(word$(v$,1))
   y=val(word$(v$,2))
   vectAngle=atan2(y,x)
end function

function vectFromPolar$(rho, phi)
   vectFromPolar$=rho*cos(phi);" ";rho*sin(phi)
end function

function vectRotate$(v$,alpha)
   x=val(word$(v$,1))
   y=val(word$(v$,2))
   rho=sqr(x*x+y*y)
   phi=atan2(y,x)+alpha
   vectRotate$=rho*cos(phi);" ";rho*sin(phi)
end function

function dePi$(x)   'pure aestetics
   dePi$=x/pi;"Pi"
end function

'---------------------------
function atan2(y,x)
   'pi = asn(1) * 2    'global
   if x <> 0 then arctan = atn(y/x)

   select case
       case x > 0
       atan2 = arctan

       case y>=0 and x<0
       atan2 = pi + arctan

       case y<0 and x<0
       atan2 = arctan - pi

       case y>0 and x=0
       atan2 = pi / 2

       case y<0 and x=0
       atan2 = pi / -2
   end select
end function


[Rod]

Thank you for that code. It will take me a while to work through. I realise now that I will have to use 2d vector maths for ALL of the objects if I intend they rotate. So mirror first.

[Rod]

Well some progress. It appears to hang together. This will draw a line representing a mirror and bounce a ray if it is on the correct side of the mirror.

I need to test hitting somewhere other than the centre of the mirror but one step at a time.

   rad=57.29577951
   WindowWidth = 220
   WindowHeight = 220
   open "mirror" for graphics_nsb as #m
   #m "trapclose [exit]"
   #m "down ; fill black ; color white"
   'draw a line about a central position 100,100
   'use angle A and calculate the two opposite points
   'on the circle to draw a line representing the mirror
   'use a radius of 50 to get a 100 pixel mirror
   'then draw a short red line at 90o to the mirror line
   'this represents the mirrors reflective vector
   'we will draw at 45o increments and shine a random ray at the mirror
   for a = 0 to 315 step 45
       'point 1
       x1=100-(50*sin(a/rad))
       y1=100-(50*cos(a/rad))
       'point 2 180o opposite
       x2=100-(50*sin((a+180)/rad))
       y2=100-(50*cos((a+180)/rad))
       '90o to A, radius 10
       x3=100-(10*sin((a+90)/rad))
       y3=100-(10*cos((a+90)/rad))

       'get 90o mirror vector to central point
       vx=100-x3
       vy=100-y3
       'get vector length
       vl=sqr(vx*vx+vy*vy)
       'calculate mirror unit vector (norm)
       'unit vector is simply vector/length to give
       'a direction between 1 and -1 for both x and y
       mnx=0-vx/vl
       mny=0-vy/vl


       'pick a random starting point on screen for the ray
       rx=int(rnd(0)*220)
       ry=int(rnd(0)*220)
       'calculate ray vector to central point
       vx=100-rx
       vy=100-ry
       'get vector length
       vl=sqr(vx*vx+vy*vy)
       'calculate ray's unit vector (norm)
       rnx=vx/vl
       rny=vy/vl

       'get dot product of mirror norm and ray norm
       'if it is - the ray is correct side of mirror
       'if it is + the ray is on wrong side of mirror
       rd=mnx*rnx+mny*rny

       'now calculate the reflection vector
       'using formula v=d-2*(n.d)*n where d = ray norm and n = mirror norm and n.d is dot product
       refx=rnx-2*rd*mnx
       refy=rny-2*rd*mny
       fx=100+refx*vl
       fy=100+refy*vl
       print "Angle ";a
       print "Dot Product";using("##.##",rd)
       print "Mirror ",using("##.##",mnx),using("##.##",mny)
       print "Ray ",using("##.##",rnx),using("##.##",rny)
       print "Reflection ",using("##.##",refx),using("##.##",refy)

       #m "cls ; fill black ; color green ; line ";rx;" ";ry;" 100 100"
       if rd<0 then #m "color blue ; line 100 100 ";fx;" ";fy
       #m "color white ; line ";x1;" ";y1;" ";x2;" ";y2;" ; color red ; line 100 100 ";x3;" ";y3

       input a$
   next
   wait

   [exit]
   close #m
   end



[tsh73]

Read a Wikipedia on refraction.
Added a refraction, now cyan thing is glass (n=1.5) half-sphere
'Mirror: (tsh73 Jun 2023)
'make ray reflect off mirror
'1. make that mirror
'2. rotate
'3. reflect center
'4. simplify direction check via dot product
'    + show mirror angle
'5. Add refraction law - now cyan thing is glass half-sphere , n=1.5
'   full internal refraction possible then going from glass to air
'*****************************************************
'nomainwin

global pi
pi = acs(-1)   'global for atan2, dePi
pi2=pi/2

desiredWidth = 400
desiredHeight = 400
mirrR=50
mirrA=0    '0 is left, 90 up
mirrAA=mirrA*pi/180-pi2    'same in radians, rotated as needed

gosub [ajustWindow]
'now, center window
   UpperLeftX = 1'(DisplayWidth - WindowWidth)/2
   UpperLeftY = 1'(DisplayHeight - WindowHeight)/2

open "Mirror" for graphics_nsb_nf as #gr
#gr "trapclose [quit]"
#gr "down; fill white; flush"
#gr "home; posxy cx cy"
#gr "setfocus"
#gr "when characterInput [key]"

gosub [msg]
gosub [drawMirror]

#gr "when leftButtonDown [start]"
#gr "when leftButtonMove [moving]"
#gr "when leftButtonUp [end]"
#gr "when rightButtonUp [ray]"

wait

[start]
   mx=MouseX:my=MouseY
   'print mx, my
   if sqr((mx-cx)^2+(my-cy)^2)>mirrR then wait
   ''or else
   started=1
   mirrAA0=mirrAA
   a0=vectAngle(vect$(mx-cx, my-cy))
   'print "started", a0
   wait

[moving]
   if not(started) then wait
   mx=MouseX:my=MouseY
   a1=vectAngle(vect$(mx-cx, my-cy))
   'print "moving", a1
   mirrAA=mirrAA0+a1-a0
   gosub [drawMirror]
   wait

[end]
   if not(started) then wait
   started=0
   'mirrAA=mirrA*pi/180-pi2
   mirrA=(mirrAA+pi2)*180/pi
   gosub [clear]
   #gr "place ";2*cx-70;" 20"
   #gr "\";using("###.#",mirrA);"   "
   wait

[ray]
   mx=MouseX:my=MouseY
   'ray from mx, my to cx, cy
   #gr "place ";mx;" ";my
   #gr "goto "; base$
   'now check if ray falls on right side of the mirror
   ray$=vectSub$(base$, vect$(mx,my))
   'if right side of the mirror - reflect
   if vectDotProduct(ray$, mirrNormal$) <0 then
       rayT$=vectTangent$(ray$,mirrVect$)
       rayN$=vectNorm$(ray$,mirrVect$)
       ' #gr "color blue"
       ' #gr "line ";base$;" "; vectAdd$(base$, vectScale$(0.1,rayT$))
       ' #gr "color red"
       ' #gr "line ";base$;" "; vectAdd$(base$, vectScale$(0.1,rayN$))
       rayRefl$=vectSub$(rayT$, rayN$)
       rayRScaled$=vectScale$(sqr(2*cx*cx),vectUnit$(rayRefl$))
       REM print cx, sqr(2*cx*cx)
       REM print rayRefl$
       REM print vectUnit$(rayRefl$)
       REM print, rayRScaled$
       #gr "color red"
       #gr "line ";base$;" "; vectAdd$(base$, rayRScaled$)

       'now, Snell's law / law of refraction
       n=1.5    'glass
       'sin (angle of incidence) /sin(angle of refraction) =1/n
       '=>sin(angle of refraction) = sin (angle of incidence)/n
       'btw if n<1 (say from glass to air), could get |sin|>1 then full intenal reflection (no refracted ray at all)
       AoI = vectAngle(ray$) - vectAngle(mirrNormal$)'angle of incidence
       print AoI*180/pi-180
       AoI = AoI - pi
       sinAoR=sin(AoI)/n
       if abs(sinAoR)>1 then
           print , "Full internal reflection"
       else
           AoR=asn(sinAoR)
           print , AoR*180/pi
           #gr "color blue"
           rayRefr$=vectRotate$(mirrNormal$, pi+AoR)
           rayRefr$=vectScale$(sqr(2*cx*cx),vectUnit$(rayRefr$))
           #gr "line ";base$;" "; vectAdd$(base$, vectScale$(2.5,rayRefr$))
       end if
   else
       'from glass to air
       rayT$=vectTangent$(ray$,mirrVect$)
       rayN$=vectNorm$(ray$,mirrVect$)
       ' #gr "color blue"
       ' #gr "line ";base$;" "; vectAdd$(base$, vectScale$(0.1,rayT$))
       ' #gr "color red"
       ' #gr "line ";base$;" "; vectAdd$(base$, vectScale$(0.1,rayN$))
       rayRefl$=vectSub$(rayT$, rayN$)
       rayRScaled$=vectScale$(sqr(2*cx*cx),vectUnit$(rayRefl$))
       REM print cx, sqr(2*cx*cx)
       REM print rayRefl$
       REM print vectUnit$(rayRefl$)
       REM print, rayRScaled$
       '#gr "color red"
       '#gr "line ";base$;" "; vectAdd$(base$, rayRScaled$)

       'now, Snell's law / law of refraction
       n=1/ 1.5    'glass to air
       'sin (angle of incidence) /sin(angle of refraction) =1/n
       '=>sin(angle of refraction) = sin (angle of incidence)/n
       'btw if n<1 (say from glass to air), could get |sin|>1 then full intenal reflection (no refracted ray at all)
       AoI = vectAngle(ray$) - vectAngle(mirrNormal$)'angle of incidence
       print AoI*180/pi-180
       AoI = AoI - pi
       sinAoR=sin(AoI)/n
       if abs(sinAoR)>1 then
           print , "Full internal reflection"
           #gr "color red"
           #gr "line ";base$;" "; vectAdd$(base$, rayRScaled$)
       else
           AoR=asn(sinAoR)
           print , AoR*180/pi
           #gr "color blue"
           rayRefr$=vectRotate$(mirrNormal$, 0-AoR)
           rayRefr$=vectScale$(sqr(2*cx*cx),vectUnit$(rayRefr$))
           #gr "line ";base$;" "; vectAdd$(base$, vectScale$(2.5,rayRefr$))
       end if
   end if

   #gr "color black"
   wait

[drawMirror]
   base$=vect$(cx, cy)
   #gr "place ";base$
   #gr "backcolor white"
   #gr "circlefilled ";mirrR    'clear
   #gr "backcolor cyan"
   #gr "piefilled ";2*mirrR;" ";2*mirrR;" ";mirrA-90;" 180"
   #gr "backcolor white"
   #gr "circle ";mirrR
   mirrPt1$=vect$(cx+mirrR*cos(mirrAA),cx+mirrR*sin(mirrAA))
   mirrPt2$=vect$(cx+mirrR*cos(mirrAA+pi),cx+mirrR*sin(mirrAA+pi))
   #gr "line ";mirrPt1$;" ";mirrPt2$
   '#gr "size 3; set ";mirrPt2$ 'mark endpoint
   '#gr "size 1"
   'get and draw a normal
   mirrVect$=vectSub$(mirrPt2$,mirrPt1$)
   mirrNormal$=vectScale$(0.1,vectRotate$(mirrVect$,pi/2))
   #gr "place ";base$
   #gr "goto "; vectAdd$(base$, mirrNormal$)
return

[msg]
   #gr "place 20 20"
   #gr "\"
   #gr "\LMB drag in a circle rotates mirror"
   #gr "\RMB shoots ray"
   #gr "\ESC clears screen"
return

[key]
   if Inkey$=chr$(27) then gosub [clear]
   wait

[clear]
   #gr "cls"
   gosub [msg]
   gosub [drawMirror]
return

[quit]
   close #gr
   end

return
'-------------------------------------------------
[ajustWindow]
   UpperLeftX = 20
   UpperLeftY = 20
   WindowWidth = 200   '100 seems to be too much - works different
   WindowHeight = 100
'    MENU #gr, "dummy"
   open "Ajusting..." for graphics_nsb_nf as #gr
'   graphics
'   graphics_nsb
'   graphics_nsb_nf

   #gr, "home ; down ; posxy x y"
   'x, y give us width, height
   width = 2*x : height = 2*y
   close #gr

   slackX = 200-width
   slackY = 100-height

   WindowWidth = desiredWidth + slackX
   WindowHeight = desiredHeight + slackY

return
'------------------------------
'GetPixelValue$ returns a string with the RGB values of the pixel
'in coordinates x and y in window/graphicbox names handle$ (e.g, "#main.graph")
function GetPixelValue$(x, y, handle$)

'Grab a 1*1 bitmap
   #handle$, "getbmp gpv "; x; " "; y; " "; 1; " "; 1

'Save in a bmp file
   bmpsave "gpv", "getpvaluetemp.bmp"

'Open the file for string input and get it's full contents
   open "getpvaluetemp.bmp" for input as #gpv
   s$ = input$(#gpv, lof(#gpv))
   close #gpv

'Check if user's display is 32-bit, and read the red-green-blue values
'If display 16 bit, then colors are masked. So some last (3 for red, 2 for green, 3 for blue) bits always 0
'That means that you did not get 255 255 255 for white - (248 252 248) instead. You have to experiment
'otherwise function returns nothing (support for other display types could be added (?))
   bpp =  asc(mid$(s$, 29, 1))
   select case bpp
   case 32
       red = asc(mid$(s$, 69, 1))
       green = asc(mid$(s$, 68, 1))
       blue = asc(mid$(s$, 67, 1))
  case 24  'from Richard, for (censored)
       red = asc(mid$(s$, 57, 1))
       green = asc(mid$(s$, 56, 1))
       blue = asc(mid$(s$, 55, 1))
   case 16
       bytes = asc(mid$( s$, 67, 1)) + 256*asc(mid$( s$, 68, 1))
       red =  (bytes AND 63488) /256       '0xF800
       green =  (bytes AND 2016) / 32 * 4  '0x7E0
       blue =  (bytes AND 31) * 8          '0x1F
   end select

'concatenate the return value, delete temporary file and free memory
   GetPixelValue$ = using("###",red)+using("####",green)+using("####",blue)
   kill "getpvaluetemp.bmp"
   unloadbmp "gpv"
end function

'===================================
function vect$(x,y)
   vect$=x;" ";y
end function

function vectX(v$)
   vectX=val(word$(v$,1))
end function

function vectY(v$)
   vectY=val(word$(v$,2))
end function

function vectLen(v$)
   x=val(word$(v$,1))
   y=val(word$(v$,2))
   vectLen=sqr(x*x+y*y)
end function

function vectUnit$(v$)
   x=val(word$(v$,1))
   y=val(word$(v$,2))
   vectLen=sqr(x*x+y*y)
   vectUnit$=x/vectLen;" ";y/vectLen
end function

function vectAdd$(v1$,v2$)
   x1=val(word$(v1$,1))
   y1=val(word$(v1$,2))
   x2=val(word$(v2$,1))
   y2=val(word$(v2$,2))
   vectAdd$=x1+x2;" ";y1+y2
end function

function vectSub$(v1$,v2$)
   x1=val(word$(v1$,1))
   y1=val(word$(v1$,2))
   x2=val(word$(v2$,1))
   y2=val(word$(v2$,2))
   vectSub$=x1-x2;" ";y1-y2
end function

function vectDotProduct(v1$,v2$)
   x1=val(word$(v1$,1))
   y1=val(word$(v1$,2))
   x2=val(word$(v2$,1))
   y2=val(word$(v2$,2))
   vectDotProduct=x1*x2+y1*y2
end function

function vectScale$(a,v$)   'a * vector v$
   x=val(word$(v$,1))
   y=val(word$(v$,2))
   vectScale$=a*x;" ";a*y
end function

function vectTangent$(v$,base$)
   n$=vectUnit$(base$)
   vectTangent$=vectScale$(vectDotProduct(n$,v$),n$)
end function

function vectNorm$(v$,base$)
   vectNorm$=vectSub$(v$,vectTangent$(v$,base$))
end function

function vectAngle(v$)
   x=val(word$(v$,1))
   y=val(word$(v$,2))
   vectAngle=atan2(y,x)
end function

function vectFromPolar$(rho, phi)
   vectFromPolar$=rho*cos(phi);" ";rho*sin(phi)
end function

function vectRotate$(v$,alpha)
   x=val(word$(v$,1))
   y=val(word$(v$,2))
   rho=sqr(x*x+y*y)
   phi=atan2(y,x)+alpha
   vectRotate$=rho*cos(phi);" ";rho*sin(phi)
end function

function dePi$(x)   'pure aestetics
   dePi$=x/pi;"Pi"
end function

'---------------------------
function atan2(y,x)
   'pi = asn(1) * 2    'global
   if x <> 0 then arctan = atn(y/x)

   select case
       case x > 0
       atan2 = arctan

       case y>=0 and x<0
       atan2 = pi + arctan

       case y<0 and x<0
       atan2 = arctan - pi

       case y>0 and x=0
       atan2 = pi / 2

       case y<0 and x=0
       atan2 = pi / -2
   end select
end function


[marshawn]

who's the fairest of them all?

[Rod]

I am still persevering, it is slow going for me. There are two vector formula I think I need. Reflection and Refraction

Reflection seems the easier of the two.

cosI = - dot(normal,incident)
reflection=incident + 2 * cosI * normal

Where normal is the mirror's normalised vector, incident is the incoming ray's normalised vector and reflection is the reflected ray's normalised vector. cosI is not cos() it is the dot product of the two mirror and incidents vectors. So far so good that is working. Here is a short demo.


       WindowWidth = 220
       WindowHeight = 220
       open "mirror" for graphics_nsb as #m
       #m "trapclose [exit]"
       #m "down ; fill black ; color white"
       rad=57.29577951
       a=0
       'point 1
       x1=100-(50*cos(a/rad))
       y1=100-(50*sin(a/rad))
       'point 2 180o opposite
       x2=100-(50*cos((a+180)/rad))
       y2=100-(50*sin((a+180)/rad))
       '90o to A, radius 10
       x3=100-(10*cos((a+90)/rad))
       y3=100-(10*sin((a+90)/rad))

       'get 90o mirror vector to central point
       vx=100-x3
       vy=100-y3
       'get vector length
       vl=sqr(vx*vx+vy*vy)
       'calculate mirror unit vector (norm)
       'unit vector is simply vector/length to give
       'a direction between 1 and -1 for both x and y
       mnx=0-vx/vl
       mny=0-vy/vl


       'pick a random starting point on screen for the ray
       rx=int(rnd(0)*220)
       ry=int(rnd(0)*22)
       'calculate ray vector to central point
       vx=100-rx
       vy=100-ry
       'get vector length
       vl=sqr(vx*vx+vy*vy)
       'calculate ray's unit vector (norm)
       rnx=vx/vl
       rny=vy/vl

       'get dot product of mirror norm and ray norm
       'if it is - the ray is correct side of mirror
       'if it is + the ray is on wrong side of mirror
       rd=mnx*rnx+mny*rny

       'now calculate the reflection vector
       'using formula v=d-2*(n.d)*n where d = ray norm and n = mirror norm
       refx=rnx-2*rd*mnx
       refy=rny-2*rd*mny

       'multiply by length of vector and offset by central point
       fx=100+refx*vl
       fy=100+refy*vl
       print "Angle ";a
       print "Dot Product";using("##.##",rd)
       print "Mirror ",using("##.##",mnx),using("##.##",mny)
       print "Ray ",using("##.##",rnx),using("##.##",rny)
       print "Reflection ",using("##.##",refx),using("##.##",refy)

       #m "cls ; fill black ; color green ; line ";rx;" ";ry;" 100 100"
       if rd<0 then #m "color blue ; line 100 100 ";fx;" ";fy
       #m "color white ; line ";x1;" ";y1;" ";x2;" ";y2;" ; color red ; line 100 100 ";x3;" ";y3
       wait

       [exit]
       close #m
       end


Now refraction, first of all I am really appreciative of the support I have had so far. It has spurred me on. But atan2! Not for me I seek a vector solution where all that trig is hidden. Here is the pseudo code.

n = n1 / n2
cosI = - dot(normal,incident)
sinT2 = n * n * (1.0 - cosI * cosI)
if sinT2 > 1 then return invalid (total internal reflection)
cosT = sqrt(1.0 - sinT2)
refraction = n * incident + (n * cosI - cosT) * normal

n1 and n2 are the refractive index of the medium the ray passes through, so air 1 and glass 1.5
cosI again is the dot product of the mirror and ray normals. Total internal reflection is handled.

So my head is down and coding is in progress, I may be a while.

[plus]

Oh man I might have to try another lecture on ATan2 ;-))

I couldn't stand not knowing how to get the angle of one point is to another. 

[Rod]

But you don't need to know the angle, it is implicit in the vector and so vector maths can proceed without disclosing it. While there seems a lot of code in my short demo, most is reusable so the normal and the dot product are identical for reflection or refraction. It comes down to two lines for reflection. Neat, compact, fast.

refx=rnx-2*rd*mnx
refy=rny-2*rd*mny

I was thinking about my old maths teacher, he did try but rugby, fencing and girls just filled my head with anything but maths. Possible to old and slow for lecturing now.

But if I do succeed then I will do something to promote vectors. Lots of coders use them without even knowing x=x+xd y=y+yd. if x>limit then xd=xd*-1

[tsh73]

cosI is not cos() it is the dot product of the two mirror and incidents vectors

Actually
a * b = |a|*|b|*cos θ
so for norm vectors (|a|=1, |b|=1) it is indeed cos() all right.

[Rod]

Sure its true and it was obvious with the variable naming convention but for readers with as much maths savvy as I have there is a difference between cos(I) and cosI that I wanted to clarify.

[Rod]

Ok, this seems to draw reflective and refractive lines but not totally proven. I need to build a more robust test. Also I now need to think about lenses as they are not flat! So I need to find where I have hit the lens and calculate the normal for the line that runs from the focal point to the intersect. Small steps.

WindowWidth = 220
       WindowHeight = 220
       open "mirror" for graphics_nsb as #m
       #m "trapclose [exit]"
       #m "down ; fill black ; color white"
       rad=57.29577951
       a=0
       for a= 0 to 315  step 45
       'point 1
       x1=100-(50*cos(a/rad))
       y1=100-(50*sin(a/rad))
       'point 2 180o opposite
       x2=100-(50*cos((a+180)/rad))
       y2=100-(50*sin((a+180)/rad))
       '90o to A, radius 10
       x3=100-(10*cos((a+90)/rad))
       y3=100-(10*sin((a+90)/rad))

       'get 90o mirror vector to central point
       vx=100-x3
       vy=100-y3
       'get vector length
       vl=sqr(vx*vx+vy*vy)
       'calculate mirror unit vector (norm)
       'unit vector is simply vector/length to give
       'a direction between 1 and -1 for both x and y
       mnx=0-vx/vl
       mny=0-vy/vl


       'pick a random starting point on screen for the ray
       rx=int(rnd(0)*220)
       ry=int(rnd(0)*220)
       'calculate ray vector to central point
       vx=100-rx
       vy=100-ry
       'get vector length
       vl=sqr(vx*vx+vy*vy)
       'calculate ray's unit vector (norm)
       rnx=vx/vl
       rny=vy/vl

       'get dot product of mirror norm and ray norm
       'if it is - the ray is correct side of mirror
       'if it is + the ray is on wrong side of mirror
       rd=mnx*rnx+mny*rny

       'now calculate the reflection vector
       'using formula v=d-2*(n.d)*n where d = ray norm and n = mirror norm
       refx=rnx-2*rd*mnx
       refy=rny-2*rd*mny

       'multiply by length of vector and offset by central point
       fx=100+refx*vl
       fy=100+refy*vl


       'calculate the refraction vector
       n=1/1.5 'refractive index of air and glass
       cosI = 0-rd ' - dot(normal,incident)
       sinT2 = n * n * (1.0 - cosI * cosI)
       '>1 means total internal reflection so only act on <=1
       if sinT2 <=1 then
           cosT = sqr(1.0 - sinT2)
           refrx= n * rnx + (n * cosI - cosT) * mnx
           refry= n * rny + (n * cosI - cosT) * mny
       end if

       'multiply by length of vector and offset by central point
       frx=100+refrx*vl
       fry=100+refry*vl

       print "Angle ";a
       print "Dot Product";using("##.##",rd)
       print "SinT2 ";using("##.##",sinT2)
       print "Mirror ",using("##.##",mnx),using("##.##",mny)
       print "Ray ",using("##.##",rnx),using("##.##",rny)
       print "Reflection ",using("##.##",refx),using("##.##",refy)
       print "Refraction ",using("##.##",refrx),using("##.##",refry)


       #m "cls ; fill black ; color green ; line ";rx;" ";ry;" 100 100"
       if rd<0 then #m "color blue ; line 100 100 ";fx;" ";fy
       if sinT2<=1 and rd<0 then #m "color cyan ; line 100 100 ";frx;" ";fry
       #m "color white ; line ";x1;" ";y1;" ";x2;" ";y2;" ; color red ; line 100 100 ";x3;" ";y3
       input a$
       next
       wait

       [exit]
       close #m
       end