' rctrps.bas - jack@ord.ca - revised 10 Oct 02 - Liberty Basic v3.01
'   tests 4 different capacitors by analyzing R-C transients
'     with nominal 20 ms time constants
'   eps file option for output to a PostScript printer

' Initialize Window ---------------------------------------------------------
  nomainwin
  textbox#1.t1, 250, 0, 75, 20
  textbox#1.t2, 335, 0, 75, 20
  button#1, "C1", [b1], UL, 10, 0, 30, 20
  button#1, "C2", [b2], UL, 50, 0, 30, 20
  button#1, "C3", [b3], UL, 90, 0, 30, 20
  button#1, "C4", [b4], UL, 130, 0, 30, 20
  button#1, "fit", [b5], UL, 170, 0, 30, 20
  button#1, "eps", [b6], UL, 210, 0, 30, 20
  WindowWidth=430     ' pixel scale 0-420
  WindowHeight=374    ' pixel scale 0-345
  UpperLeftX=50: UpperLeftY=10
  open "Capacitor Test" for graphics_nsb as #1
  #1 "trapclose [quit]"
  wide=400: high=320: xm=100: ym=4096: nx=10: ny=8: np=32: kk=0
  dim x(33): dim y(33): dim z(33)
  gosub [grid]

[waitHere]
  wait

[quit]
  close #1
end

[b1] ' Select and Plot Data -------------------------------------------------
  f$="RC1"
    x(1)=0:   y(1)=4086:  x(2)=1:   y(2)=3843:  x(3)=2:   y(3)=3651
    x(4)=3:   y(4)=3468:  x(5)=4:   y(5)=3294:  x(6)=5:   y(6)=3129
    x(7)=6:   y(7)=2972:  x(8)=7:   y(8)=2823:  x(9)=8:   y(9)=2682
    x(10)=9:  y(10)=2548: x(11)=10: y(11)=2420: x(12)=11: y(12)=2299
    x(13)=12: y(13)=2184: x(14)=13: y(14)=2075: x(15)=14: y(15)=1971
    x(16)=15: y(16)=1873: x(17)=16: y(17)=1779: x(18)=17: y(18)=1690
    x(19)=18: y(19)=1606: x(20)=19: y(20)=1526: x(21)=20: y(21)=1450
    x(22)=21: y(22)=1377: x(23)=23: y(23)=1243: x(24)=25: y(24)=1123
    x(25)=27: y(25)=1014: x(26)=30: y(26)=870:  x(27)=33: y(27)=747
    x(28)=37: y(28)=610:  x(29)=41: y(29)=498:  x(30)=47: y(30)=368
    x(31)=55: y(31)=247:  x(32)=69: y(32)=125
  goto [plotdata]

[b2]
  f$="RC2"
    x(1)=0:   y(1)=4086:  x(2)=1:   y(2)=3842:  x(3)=2:   y(3)=3650
    x(4)=3:   y(4)=3468:  x(5)=4:   y(5)=3295:  x(6)=5:   y(6)=3131
    x(7)=6:   y(7)=2975:  x(8)=7:   y(8)=2827:  x(9)=8:   y(9)=2686
    x(10)=9:  y(10)=2552: x(11)=10: y(11)=2425: x(12)=11: y(12)=2305
    x(13)=12: y(13)=2190: x(14)=13: y(14)=2081: x(15)=14: y(15)=1978
    x(16)=15: y(16)=1880: x(17)=16: y(17)=1787: x(18)=17: y(18)=1698
    x(19)=18: y(19)=1614: x(20)=19: y(20)=1534: x(21)=20: y(21)=1458
    x(22)=21: y(22)=1385: x(23)=23: y(23)=1252: x(24)=25: y(24)=1131
    x(25)=27: y(25)=1022: x(26)=30: y(26)=878:  x(27)=33: y(27)=755
    x(28)=37: y(28)=617:  x(29)=41: y(29)=505:  x(30)=47: y(30)=374
    x(31)=55: y(31)=252:  x(32)=70: y(32)=123
  goto [plotdata]

[b3]
  f$="RC3"
    x(1)=0:   y(1)=4086:  x(2)=1:   y(2)=3870:  x(3)=2:   y(3)=3698
    x(4)=3:   y(4)=3536:  x(5)=4:   y(5)=3381:  x(6)=5:   y(6)=3233
    x(7)=6:   y(7)=3092:  x(8)=7:   y(8)=2957:  x(9)=8:   y(9)=2829
    x(10)=9:  y(10)=2706: x(11)=10: y(11)=2588: x(12)=11: y(12)=2476
    x(13)=12: y(13)=2368: x(14)=13: y(14)=2266: x(15)=14: y(15)=2168
    x(16)=15: y(16)=2074: x(17)=16: y(17)=1985: x(18)=17: y(18)=1899
    x(19)=19: y(19)=1739: x(20)=21: y(20)=1593: x(21)=22: y(21)=1524
    x(22)=24: y(22)=1397: x(23)=27: y(23)=1226: x(24)=29: y(24)=1124
    x(25)=32: y(25)=988:  x(26)=35: y(26)=868:  x(27)=38: y(27)=764
    x(28)=43: y(28)=618:  x(29)=48: y(29)=501:  x(30)=55: y(30)=375
    x(31)=65: y(31)=250:  x(32)=83: y(32)=124
  goto [plotdata]

[b4]
  f$="RC4"
    x(1)=0:   y(1)=4086:  x(2)=1:   y(2)=3875:  x(3)=2:   y(3)=3709
    x(4)=3:   y(4)=3562:  x(5)=4:   y(5)=3428:  x(6)=5:   y(6)=3304
    x(7)=6:   y(7)=3188:  x(8)=7:   y(8)=3078:  x(9)=8:   y(9)=2974
    x(10)=9:  y(10)=2873: x(11)=10: y(11)=2776: x(12)=11: y(12)=2681
    x(13)=13: y(13)=2494: x(14)=14: y(14)=2399: x(15)=15: y(15)=2302
    x(16)=17: y(16)=2103: x(17)=18: y(17)=2003: x(18)=19: y(18)=1904
    x(19)=21: y(19)=1721: x(20)=22: y(20)=1640: x(21)=24: y(21)=1497
    x(22)=26: y(22)=1376: x(23)=28: y(23)=1271: x(24)=31: y(24)=1135
    x(25)=34: y(25)=1018: x(26)=38: y(26)=886:  x(27)=43: y(27)=749
    x(28)=48: y(28)=635:  x(29)=55: y(29)=506:  x(30)=64: y(30)=379
    x(31)=77: y(31)=251:  x(32)=100:y(32)=125
  goto [plotdata]

[plotdata]
  kk=1: gosub [grid]
  #1 "backcolor green"
  for i=1 to np
    x1=4*x(i)+10: y1=int(high*(1-y(i)/ym)+.5)+25
    #1 "place "; x1; " "; y1: #1 "circlefilled 5"
  next i
  #1 "backcolor white"
goto [waitHere]

[grid]
  #1 "cls ; down ; color black"
  for i=0 to nx
    #1 "line "; 10+40*i; " 25 "; 10+40*i; " "; 25+high
  next i
  for i=0 to ny
    #1 "line 10 "; 25+40*i; " "; 10+wide; " "; 25+40*i
  next i
  #1.t1 "time const"
  #1.t2 "rms dev"
return

[b5] ' Fit Data and Plot Results --------------------------------------------
  if kk=1 then goto [fit] else goto [waitHere]

[fit]
  test=0: kk=2: t=x(np)/log(np): dt=t/4: ymn=-4: ymx=4
  gosub [lsq]: h=lsy
  while test=0                                         ' vary time constant
    t=t+dt: gosub [lsq]
    if (lsy>h) then
      dt=0-dt/2
      if abs(dt)<t/100000 then test=1: lsy=h: t=t+2*dt ' best time constant
    end if
    h=lsy
  wend
  gosub [lsq]                                          ' best fit
  #1 "backcolor blue"                                  ' plot y(i) vs z(i)
  for i=1 to np
    x1=int(z(i)*wide+.5)+10: y1=int(high*(1-y(i)/ym)+.5)+25
    #1 "place "; x1; " "; y1: #1 "circlefilled 5"
  next i
  y1=int(high*(1-lb/ym)+.5)+25: y2=int(high*(1-(la+lb)/ym)+.5)+25
  #1 "color red ; line 10 "; y1; " 410 "; y2           ' fitted line
  #1 "place 10 "; y2
  for i=1 to 100                                       ' fitted curve
    y1=int(high*(1-(la*exp(0-i/t)+lb)/ym)+.5)+25
    #1 "goto "; 10+4*i; " "; y1
  next i
  #1 "color black ; backcolor red"
  for i=1 to np                                        ' scatter plot
    yy=y(i)-la*z(i)-lb
    if lsy>3 then yy=yy/20: ymn=-80: ymx=80
    yy=25+high/2-40*yy
    x1=int(wide*z(i)+.5)+10: y1=int(yy+.5)
    #1 "place "; x1; " "; y1: #1 "circlefilled 5"
  next i
  #1 "backcolor white"
  #1.t1 "t = "; left$(str$(t),6)
  #1.t2 "d = "; left$(str$(lsy),6)
goto [waitHere]

[lsq]                                                  ' fit y=la*x+lb
  l0=0: l1=0: l2=0: l3=0: l4=0: ls=0
  for i=1 to np                                        ' convert to exp
    z(i)=exp(0-x(i)/t)
  next i
  for i=2 to np                                        ' start analysis ..
    l1=l1+z(i): l4=l4+y(i)                             ' .. from 2nd point
  next i
  l1=l1/(np-1): l4=l4/(np-1)
  for i=2 to np
    l0=z(i)-l1: l2=l2+l0*(y(i)-l4): l3=l3+l0*l0
  next i
  la=l2/l3: lb=l4-la*l1
  for i=2 to np
    l0=y(i)-la*z(i)-lb: ls=ls+l0*l0
  next i
  lsy=(ls/(np-3))^.5                                   ' deviation
return

[b6] ' Write PostScript .eps File -------------------------------------------
  if kk=2 then goto [eps] else goto [waitHere]

[eps]
'Set axis labels, title
  xt1$="TIME  milliseconds": yt1$="A/D READING  lsb"
  xt2$="EXP(-t/T)": yt2$="DEVIATION  lsb"
'Set origin of hardcopy graph in points (size in points = size in pixels)
  xorg=120: yorg=340
  'note bounding box: lower left (xorg-60, yorg-140)
  '                   upper right (xorg+wide+60, yorg+high+60)
'Set axis tic option (0 for grid, 1 for tics)
  tic=0

  open f$+".eps" for output as #2
  #2 "%!PS-Adobe-3.0 EPSF-2.0"
  #2 "%%Creator: Jack Ord"
  #2 "%%Title: "; f$+".eps"
  #2 "%%CreationDate: 12 Nov 99"
  #2 "%%Pages: 0"
  #2 "%%BoundingBox: 60 200 580 720"
  #2 "/scaleHelv {/Helvetica findfont exch scalefont setfont} def"
  #2 "/scaleSymb {/Symbol findfont exch scalefont setfont} def"
  #2 "/centerx"
  #2 "{dup stringwidth pop 2 div neg 0 rmoveto show} def"
  #2 "/justify"
  #2 "{dup stringwidth pop 6 add neg 0 rmoveto show} def"
  #2 "/centery"
  #2 "{dup stringwidth pop 2 div neg 0 exch rmoveto"
  #2 "90 rotate show -90 rotate} def"
  #2 "/tau {14 scaleSymb (\164) show 14 scaleHelv} def"
  #2 "14 scaleHelv .5 setlinewidth"
  if tic=0 then zt=high else zt=9
  for i=0 to nx: x=wide*i/nx+xorg                      ' x axis
    if i>0 and i<nx then
      #2 "newpath "; x; " "; yorg; " moveto 0 "; zt; " rlineto stroke"
    end if
    #2 x; " "; yorg; " 16 sub moveto (";
    #2 xm*i/nx; : #2 ") centerx"
    #2 x; " "; yorg+high+4; " moveto (";
    #2 i/nx; : #2 ") centerx"
  next i
  if tic=1 then                                        ' tic option
    for i=1 to nx-1: x=wide*i/nx+xorg
      #2 "newpath "; x; " "; yorg + high; " moveto 0 "; 0-zt; " rlineto stroke"
    next i
    for i=1 to ny-1: y=high*i/ny+yorg
      #2 "newpath "; xorg + wide; " "; y; " moveto "; 0-zt; " 0 rlineto stroke"
    next i
  end if
  for i=0 to ny: y=high*i/ny+yorg                      ' y axis
    if tic=0 then zt=wide
    if i>0 and i<ny then
      #2 "newpath "; xorg; " "; y; " moveto "; zt; " 0 rlineto stroke"
    end if
    #2 xorg; " "; y-5; " moveto (";
    #2 ym*i/ny; : #2 ") justify"
    #2 xorg+wide+4; " "; y-5; " moveto (";
    #2 ymn+(ymx-ymn)*i/ny; : #2 ") show"
  next i
  x=xorg+wide/2: y=yorg-36                             ' x axis caption
  #2 x; " "; y; " moveto ("; xt1$; ") centerx"
  x=xorg+wide/2-24: y=yorg+high+24
  #2 x; " "; y; " moveto (exp\(-t/) show tau (\)) show"
  x=xorg: y=yorg+high/2                                ' y axis caption
  #2 x; " "; y; " moveto ("; str$(ym); ") stringwidth pop neg 14 sub 0 rmoveto (";
  #2 yt1$; ") centery"
  x=xorg+wide: y=yorg+high/2
  #2 x; " "; y; " moveto ("; str$(ymn); ") stringwidth pop 20 add 0 rmoveto (";
  #2 yt2$; ") centery"
  #2 "newpath "; xorg; " "; yorg; " moveto "; wide; " 0 rlineto 0 "
  #2 high; " rlineto "; 0-wide; " 0 rlineto closepath 1.2 setlinewidth"
  #2 "stroke 12 scaleHelv .9 setlinewidth"

  gray=1
  for i=1 to np                                        ' plot data vs t
    x1=x(i): y1=y(i): gosub [pplot]
  next i
  gray=.75
  for i=1 to np                                        ' plot data vs exp(-t/T)
    x1=100*z(i): y1=y(i): gosub [pplot]
  next i
  gray=.25
  for i=1 to np                                        ' plot deviation
    y1=y(i)-la*z(i)-lb
    if lsy>3 then y1=y1/20
    x1=100*z(i): y1=int(ym/2*(1+y1/4)+.5)
    gosub [pplot]
  next i
  x1=0: y1=lb: gosub [lstrt]                           ' plot fitted line
  x1=100: y1=la+lb: gosub [lcont]
  x1=0: gosub [lstrt]
  for i=1 to 100                                       ' plot fitted curve
    x1=i: y1=la*exp(0-i/t)+lb: gosub [lcont]
  next i
  z$="120 280": lbl$="Transient "+f$: gosub [label]    ' title and comments
  z$="120 260": lbl$="Time Constant = "+left$(str$(t),6)+" milliseconds": gosub [label]
  z$="120 240": lbl$="Deviation (rms) = "+left$(str$(lsy),6)+" lsb": gosub [label]
  #2 "stroke showpage": #2 "%%EOF": CLOSE #2
  kk=0
goto [waitHere]

[label]
  #2 "stroke newpath "; z$; " moveto ("; lbl$; ") show"
return

[lstrt]
  gosub [sc]
  #2 "stroke newpath "; z$; " moveto"
return

[lcont]
  gosub [sc]
  #2 z$; " lineto"
return

[pplot]
  gosub [sc]
  #2 "stroke newpath "; z$; " 4 0 360 arc"
  #2 "gsave "; gray; " setgray fill grestore 0 setgray stroke"
return

[sc]
  zx=xorg+int(10*wide*x1/xm)/10
  zy=yorg+int(10*high*y1/ym)/10
  z$=str$(zx)+" "+str$(zy)
return