'apollo6.bas -  Rendezvous
SCREEN 12:  WINDOW (-5700, -3700)-(2300, 2300): MRad = 1738
CIRCLE (0, 0), MRad, 2: xt = MRad * COS(.1): yt =MRad * SIN(.1)
LINE (-xt, -yt)-(xt, yt), 2: PAINT (0, -500), 2, 2 'sunlit moon
PRINT " Note - Altitude exaggerated five times."
pi = 3.141593: LOCATE 15, 42: PRINT "(a)": LOCATE 2
PRINT "(a)Ascent stage at 18 km in 16 x 83 km orbit"
CMAlt = 112: CMR = MRad + CMAlt   'Command Module alt. & radius
CMRplot = MRad + 5 * (CMR - MRad)   'CM exaggerated plot radius
mgc = 4897: cw = SQR(mgc / CMR ^ 3) 'gravity const., CM ang vel
angc = pi + .12    'CM initial angle anticlockwise from X+ axis
ma = 2474: t = 0  'AM mass at time zero, end of take-off thrust
tl = 1364: thr = tl * .0098      'AM thrust, kg.wt and kg force
SIl = 315: ff = tl / SIl          'Specific Impulse & fuel flow
peri = 1754: apo = 1821: rv2 = mgc *2 *apo / (peri *(peri+apo))
Er = rv2 / 2 - mgc / peri   'reqd vel.sqrd & energy at perilune
h = 18: ra = 1756: an = pi + .4     'AM height, radius  & angle
w = SQR(2 * (Er + mgc / ra)) 'vel. at 18 km alt. (not perilune)
c = COS(an): s = SIN(an): x = ra * c: y = ra * s   'AM position
anc = an + .005: u = w * SIN(anc): v = -w * COS(anc)  'and vel.
gm = mgc / (ra * ra): aax = gm * c: aay = gm * s  'grav. accel.
zerot = -5500: zeroh = -3500: lookup22 =SIN(22.127 /57.3) *2000
LINE (zerot, 0)-(zerot, zeroh):LINE -(-2000, zeroh) 'lower axes
LOCATE 25, 4: PRINT "10 km": LOCATE 20, 4: PRINT "20 km"
LOCATE 15, 4: PRINT "30 km"            'height scale,lower left
FOR k = 1 TO 3: LINE (zerot, 1000 *k +zeroh)-STEP(200, 0): NEXT
sector = 0: col = 2              'ready to coast out to apolune
DO: t = t + 1 'looping every second.  First for the common code
 u = u + gmx: v = v + gmy 'velocity incremented by 1 sec. accel
 x = x + u: y = y + v     'position incremented by 1 sec. veloc
 ra2 = x *x + y *y: ra = SQR(ra2): gm = mgc / ra2 'radius & acc
 w2 = u * u + v * v: w = SQR(w2): E = w2 / 2 - mgc / ra 'veloc.
 gmx = -gm * x /ra: gmy = -gm * y / ra 'acceleration components
 angc = angc - cw      'decrement CM angle by 1 sec. of ang vel
 IF angc < 0 THEN angc = angc + pi * 2  'keep in range 0 to 2pi
 IF t MOD 5 = 0 THEN                   'plot things each 5 secs
  sc = SIN(angc): cc = COS(angc): xc = CMR * cc: yc = CMR * sc
  rn = SQR((xc - x) ^ 2 + (yc - y) ^ 2)  'AM separation from CM
  xcp = CMRplot * cc: ycp = CMRplot * sc: PSET (xcp, ycp), 15
  h = ra - MRad: angle = ATN(y / x)        'AM height and angle
  IF x < 0 THEN angle = angle + pi     'put in correct quadrant
  IF x > 0 AND y < 0 THEN angle = angle + 2 * pi         'ditto
  c = COS(angle): s = SIN(angle): rap = MRad + 5 * h 'plot rad.
  xp = rap * c: yp = rap * s: PSET (xp, yp), col 'plot position
  tp = t / 3 + zerot: lookup = (CMR - ra) * 2000 / rn
  PSET (tp, lookup + zeroh), 4                'plot lookup, red
  PSET (tp, (CMR - ra) * 100 + zeroh), 15   'height diff, white
  PSET (tp, rn * 100 + zeroh), 2       'plot slant range, green
  PSET (tp, ra *100 *(angle -angc) +zeroh),14 'hor diff, yellow
 END IF                 'end IF t MOD 5 ...  and of common code
 SELECT CASE sector          'Now, code specific to each sector
  CASE 0:              'AFTER TAKE-OFF, COASTING OUT TO APOLUNE
   IF t = 3000 THEN   'after 50 min., is at apolune, 83 km alt.
    LOCATE 5, 69: COLOR 4: PRINT "(b)": LOCATE 3
    PRINT "(b) Co-elliptic Seq. Init. CM"; CINT(rn); "km ahead"
    LINE (tp, zeroh)-(tp, -2500)'vert. line at (b), lower graph
    'at apolune, must burn to get circular orbit
    neww = SQR(mgc / ra)  'velocity required for circular orbit
    acc = thr / ma     'acceleration from engine at full thrust
    burn = (neww - w) / acc: fuel = burn * ff 'burn time & fuel
    ang = ATN(y / x): u = neww * SIN(ang): v = -neww * COS(ang)
    LINE (0, 0)-(xp, yp), 4: LINE -(xcp, ycp), 4 'radial at (b)
    PRINT"(b)Burnt"; CINT(fuel);"kg fuel in"; CINT(burn); "sec"
    PRINT " Circular orbit at"; CINT(h); "km and"; neww; "km/s"
    ma = ma-fuel: sector = 1: col = 1: DO:LOOP WHILE INKEY$ =""
   END IF
  CASE 1:       'NOTE CDH AT POINT (c). No burn in this program
   IF t = 6680 THEN
    LOCATE 18, 41: COLOR 4: PRINT "(c)": LOCATE 6
    PRINT "(c) Correction Differential Height"
    PRINT "  CM"; CINT(rn); "km ahead. Minor burn"
    LINE (tp, zeroh)-(tp, -2500) 'vert.line at (c), lower graph
    LINE (0, 0)-(xp, yp), 4: LINE -(xcp, ycp), 4 'radial at (c)
    sector = 2: DO: LOOP WHILE INKEY$ = ""
   END IF
  CASE 2:       'IN CIRCULAR ORBIT WATCHING FOR CRITICAL LOOKUP
   IF lookup > lookup22 THEN
    LINE (-3500, lookup22 + zeroh)-(tp, lookup22 + zeroh), 7
    LOCATE 5, 42: COLOR 4: PRINT "(d)": LOCATE 8'mark point (d)
    PRINT "(d) Terminal Phase Initiation"
    PRINT "  CM"; CINT(rn); "km ahead. For 83 X 112 km"
    LINE (tp, zeroh)-(tp, -2500)     'vertical line lower graph
    LINE (0, 0)-(xp, yp), 4: LINE -(xcp, ycp), 4 'radial at (d)
    peri = 1821: apo = 1849  'transfer orbit, 83 to 112 km alt.
    rw2 = mgc * 2 * apo / (peri * (peri + apo))  'reqd vel sqrd
    neww = SQR(rw2): acc = thr / ma 'new velocity & thrust acc.
    burn = (neww - w) / acc: fuel = burn * ff 'burn time & fuel
    ang = ATN(-y / x): u = neww * SIN(ang): v = neww * COS(ang)
    PRINT "  orbit, burnt"; CINT(fuel); "kg fuel": col = 2
    ma = ma - fuel: sector = 3: DO: LOOP WHILE INKEY$ = ""
   END IF
  CASE 3:             'COASTING IN TRANSFER ORBIT TO RENDEZVOUS
   IF rn < .2 THEN                    'until within .2 km of CM
    LINE (0, 0)-(xp, yp), 15: LINE -(xcp, ycp),15'radial at (e)
    LOCATE 20, 76: COLOR 15: PRINT "(e)": LOCATE 11
    PRINT "(e) Now within .2 km ie."; INT(rn * 1000) / 1000
    LINE (tp, zeroh)-(tp, -2500)          'vertical line at (e)
    sector = 4: DO: LOOP WHILE INKEY$ = ""
   END IF   'No CASE 4, so only common code runs in "overshoot"
 END SELECT
LOOP UNTIL t > 16000                 'terminate after overshoot
COLOR 4: LOCATE 26, 7: PRINT "lookup decision": LOCATE 25, 60
PRINT "Lookup Angle": LOCATE 23, 12:PRINT "(b)         (c) (d)"
COLOR 15: LOCATE 26, 60: PRINT "Height Diff.": LOCATE 23, 3
PRINT "(a)": LOCATE 23, 40: PRINT "(e)": COLOR 14: LOCATE 27,60
PRINT "Horizontal Diff.": COLOR 2: LOCATE 28, 60
PRINT "Slant Distance": DO: LOOP WHILE INKEY$ = ""