REM Start of BASIC! Program

FN.DEF det(a1,b1,c1,a2,b2,c2, a3,b3,c3)
 FN.RTN a1*b2*c3+a2*b3*c1+a3*b1*c2- a3*b2*c1-a2*b1*c3-a1*b3*c2
FN.END


FN.DEF ellipse(mx,my,px,py,qx,qy)
 !This draws an ellipse,
 !fitted into a parallelogram  
 !calculated with
 !Rytz's construction.

 !turn p around m with 90 degrees
 !result is p2
 l1=HYPOT(px-mx,py-my)
 w1=ATAN2(py-my,px-mx)+PI()/2
 p2x=mx+COS(w1)*l1
 p2y=my+SIN(w1)*l1

 !Middle between p2 and q
 nx=(p2x+qx)/2
 ny=(p2y+qy)/2

 !radius
 r=HYPOT(nx-mx,ny-my)

 !Direction vector p2-q, 
 rvnx=p2x-qx
 rvny=p2y-qy
 !normalized
 l3=HYPOT(rvnx,rvny)
 rvnx/=l3
 rvny/=l3

 !Point a
 ax=nx+rvnx*r
 ay=ny+rvny*r

 !Point b
 bx=nx-rvnx*r
 by=ny-rvny*r

 !Direction vector a-m,
 rvax=ax-mx
 rvay=ay-my
 !normalized
 l4=HYPOT(rvax,rvay)
 rvax/=l4
 rvay/=l4

 !Direction vector b-m,
 rvbx=bx-mx
 rvby=by-my
 !normalized
 l5=HYPOT(rvbx,rvby)
 rvbx/=l5
 rvby/=l5

 !Distance a-q
 r1=HYPOT(ax-qx,ay-qy)

 !Distance b-q
 r2=HYPOT(bx-qx,by-qy)
 GR.SCREEN sx,sy
 mx+=(sx/2)
 my+=(sy/2)
 w2=TODEGREES(ATAN2(rvay,rvax))
 GR.ROTATE.START w2,mx,my
 GR.OVAL vd,mx-r2,my-r1,mx+r2,my+r1
 GR.ROTATE.END

 FN.RTN go
FN.END

!back face np
!READ.DATA -250, 150,-150
!READ.DATA -250,-150,-150
!READ.DATA -200,-150,-150
READ.DATA -300,-200,0
READ.DATA  300,-200,0
!READ.DATA  200,-150,-150
!READ.DATA  250,-150,-150
!READ.DATA  250, 150,-150
!READ.DATA  200, 150,-150
READ.DATA  300, 200,0
READ.DATA -300, 200,0
!READ.DATA -200, 150,-150
!front face np
!READ.DATA -250, 150,150
!READ.DATA -250,-150,150
!READ.DATA -200,-150,150
!READ.DATA -200,-100,100
!READ.DATA  200,-100,100
!READ.DATA  200,-150,150
!READ.DATA  250,-150,150
!READ.DATA  250, 150,150
!READ.DATA  200, 150,150
!READ.DATA  200, 100,100
!READ.DATA -200, 100,100
!READ.DATA -200, 150,150


!center data np
!READ.DATA    0,   0,   0
!READ.DATA  300,   0,   0
!READ.DATA    0, 300,   0
!READ.DATA    0,   0, 300


!knee bone connected to thigh bone, what order to make lines nlz

READ.DATA 5,1,2,3,4,1

np=45
DIM px[np],py[np],pz[np]
DIM qx[np],qy[np],qz[np]
DIM rx[np],ry[np],rz[np]

! all twenty points being added into a "px" dimensional array
FOR i=1 TO 4
 READ.NEXT px[i],py[i],pz[i]
NEXT i

FOR i=1 TO 41
 w=i/40*pi()*2
 px[4+i]=COS(w)*300
 py[4+i]=SIN(w)*200
 pz[4+i]=0
NEXT i

! casting the points into a particular order 
nlz=2
DIM lz[nlz]
FOR i=1 TO 1
 LIST.CREATE n,llz
 lz[i]=llz
 READ.NEXT cnt
 FOR j=1 TO cnt
  READ.NEXT lin
  LIST.ADD llz,lin
 NEXT j
NEXT i

LIST.CREATE n,llz
lz[nlz]=llz
FOR i=1 TO 41
 LIST.ADD llz,4+i
NEXT i

!scales the size o the elements 500 increases makes smaller both eoz = goz all objects disappear
erx=0
ery=0
erz=1
eox=0
eoy=0
eoz=300

!scales the size of the elements 1000 decreases makes smaller
gox=0
goy=0
goz=-1000

! invokes the graphics
GR.OPEN 255,0,0,0,0,0
GR.SCREEN sx,sy
!GR.COLOR 255,255,255,255,1
GR.COLOR 255,0,255,0,1
GR.SET.STROKE 1

GOSUB rotate
GOSUB draw
GR.COLOR 128,0,0,255,0
GOSUB draw2
GR.RENDER

! where your finger is on the screen 
DO
 GR.CLS
 GR.TOUCH tt,tx,ty

 IF tt&!ott THEN
  otx=tx
  oty=ty


 ELSEIF tt&ott THEN
  wx+=(ty-oty)
  wy+=(tx-otx)

  GOSUB rotate
  GR.COLOR 255,0,255,0,0
  GOSUB draw
  GR.COLOR 255,0,0,255,0
  GOSUB draw2
  GR.RENDER
 ELSE
  PAUSE 100
 ENDIF

 ott=tt
 otx=tx
 oty=ty
UNTIL 0

END

!Region cut Straight plane
wts:
l=SQR(grx*grx+gry*gry+grz*grz)
grx/=l
gry/=l
grz/=l

l=SQR(erx*erx+ery*ery+erz*erz)
erx/=l
ery/=l
erz/=l

!Lot in the normal direction of the plane 
! set a value to 0, swap the other two and change the sign in one


IF erx=0 THEN
 er1x=0
 er1y=erz
 er1z=-ery
ELSEIF ery=0 THEN
 er1x=erz
 er1y=0
 er1z=-erx
ELSE
 er1x=ery
 er1y=-erx
 er1z=0
ENDIF

!1. Normalize solder
l=SQR(er1x*er1x+er1y*er1y+er1z*er1z)
er1x/=l
er1y/=l
er1z/=l

!2. Lot on normal direction level and 1st lot (cross product)
er2x=ery*er1z-erz*er1y
er2y=erx*er1z-erz*er1x
er2z=erx*er1y-ery*er1x
!Is automatically normalized
l=SQR(er2x*er2x+er2y*er2y+er2z*er2z)
er2x/=l
er2y/=l
er2z/=l

!For straight = Ov-a * Rv den 
! Calculate parameter a

d1=det(gox-eox,goy-eoy,goz-eoz, er1x,er1y,er1z, er2x,er2y,er2z)
d2=det(grx,gry,grz,er1x,er1y,er1z, er2x,er2y,er2z)
a=d1/d2

!Intersection on the line:
x=gox-grx*a
y=goy-gry*a
z=goz-grz*a
RETURN

!EndRegion

!Region Array mit Matrix
rotate:

a=wx*PI()/180
ca=COS(a)
sa=SIN(a)

!Rotations matrix x

x1=1
y1=0
z1=0

x2=0
y2=ca
z2=-sa

x3=0
y3=sa
z3=ca

FOR i=1 TO np
 !Multiply vector with matrix
 qx[i]=px[i]*x1+py[i]*y1+pz[i]*z1
 qy[i]=px[i]*x2+py[i]*y2+pz[i]*z2
 qz[i]=px[i]*x3+py[i]*y3+pz[i]*z3
NEXT i

a=wy*PI()/180
ca=COS(a)
sa=SIN(a)

!Rotations matrix y

x1=ca
y1=0
z1=-sa

x2=0
y2=1
z2=0

x3=sa
y3=0
z3=ca

FOR i=1 TO np
 !Multiply vector with matrix
 rx[i]=qx[i]*x1+qy[i]*y1+qz[i]*z1
 ry[i]=qx[i]*x2+qy[i]*y2+qz[i]*z2
 rz[i]=qx[i]*x3+qy[i]*y3+qz[i]*z3
NEXT i

RETURN
!EndRegion

!----------------------------------
draw:


FOR i=1 TO nlz
 LIST.SIZE lz[i],nlin
 LIST.GET lz[i],1,p1
 grx=rx[p1]-gox
 gry=ry[p1]-goy
 grz=rz[p1]-goz
 GOSUB wts
 x2=x
 y2=y
 FOR j=2 TO nlin
  LIST.GET lz[i],j,p2

  grx=rx[p2]-gox
  gry=ry[p2]-goy
  grz=rz[p2]-goz

  GOSUB wts

  GR.LINE vd,x+sx/2,y+sy/2,x2+sx/2,y2+sy/2

  x2=x
  y2=y

 NEXT j
NEXT i


RETURN

draw2:

grx=(rx[1]+rx[3])/2-gox
gry=(ry[1]+ry[3])/2-goy
grz=(rz[1]+rz[3])/2-goz
GOSUB wts
mx=x
my=y

grx=(rx[1]+rx[2])/2-gox
gry=(ry[1]+ry[2])/2-goy
grz=(rz[1]+rz[2])/2-goz
GOSUB wts
px1=x
py1=y

grx=(rx[2]+rx[3])/2-gox
gry=(ry[2]+ry[3])/2-goy
grz=(rz[2]+rz[3])/2-goz
GOSUB wts
qx=x
qy=y

grx=(rx[3]+rx[4])/2-gox
gry=(ry[3]+ry[4])/2-goy
grz=(rz[3]+rz[4])/2-goz
GOSUB wts
q2x=x
q2y=y

grx=(rx[4]+rx[1])/2-gox
gry=(ry[4]+ry[1])/2-goy
grz=(rz[4]+rz[1])/2-goz
GOSUB wts
q3x=x
q3y=y



grx=rx[1]-gox
gry=ry[1]-goy
grz=rz[1]-goz
GOSUB wts
s1x=x
s1y=y

grx=rx[2]-gox
gry=ry[2]-goy
grz=rz[2]-goz
GOSUB wts
s2x=x
s2y=y

grx=rx[3]-gox
gry=ry[3]-goy
grz=rz[3]-goz
GOSUB wts
s3x=x
s3y=y

grx=rx[4]-gox
gry=ry[4]-goy
grz=rz[4]-goz
GOSUB wts
s4x=x
s4y=y



a1=q2y-py1
b1=px1-q2x
c1=q2x*py1-px1*q2y

a2=q3y-qy
b2=qx-q3x
c2=q3x*qy-qx*q3y

Parallel = A1 * B2 - A2 * B1

IF Parallel <> 0 THEN
 spx=(B1*C2-B2*C1)/Parallel
 spy=(A2*C1-A1*C2)/Parallel
ELSE

ENDIF

mx=spx
my=spy

gr.color 255,255,0,0,0
vd=ellipse(mx,my,px1,py1,qx,qy)





t1x=(s1x+s2x)/2
t1y=(s1y+s2y)/2

t2x=(s2x+s3x)/2
t2y=(s2y+s3y)/2

t3x=(s3x+s4x)/2
t3y=(s3y+s4y)/2

t4x=(s4x+s1x)/2
t4y=(s4y+s1y)/2


a1=t3y-t1y
b1=t1x-t3x
c1=t3x*t1y-t1x*t3y

a2=t4y-t2y
b2=t2x-t4x
c2=t4x*t2y-t2x*t4y



Parallel = A1 * B2 - A2 * B1

IF Parallel <> 0 THEN
 spx=(B1*C2-B2*C1)/Parallel
 spy=(A2*C1-A1*C2)/Parallel
ELSE

ENDIF

mx=spx
my=spy

gr.color 255,255,0,0,1
gr.line vd,px1+sx/2,py1+sy/2,q2x+sx/2,q2y+sy/2
gr.line vd,qx+sx/2,qy+sy/2,q3x+sx/2,q3y+sy/2
gr.line vd,px1+sx/2,py1+sy/2,0,0


gr.color 255,0,0,255,0
gr.line vd,t1x+sx/2,t1y+sy/2,t3x+sx/2,t3y+sy/2
gr.line vd,t2x+sx/2,t2y+sy/2,t4x+sx/2,t4y+sy/2
gr.line vd,mx+sx/2,my+sy/2,0,0

vd=ellipse(mx,my,px1,py1,qx,qy)

RETURN