DECLARE SUB font (s$, x!, y!)
DECLARE FUNCTION lax! (x AS LONG, y AS LONG, z AS LONG)
DECLARE FUNCTION lay! (x AS LONG, y AS LONG, z AS LONG)
DECLARE SUB prod (a AS ANY, b AS ANY, c AS ANY)
DECLARE SUB campo (x AS LONG, y AS LONG, z AS LONG, a1 AS ANY)
' il centro della Terra Š in O - asse parallelo a y
TYPE vet
  x AS DOUBLE
  y AS DOUBLE
  z AS DOUBLE
END TYPE
TYPE partic
  vel AS vet
  mas AS DOUBLE
  car AS DOUBLE
  ene AS DOUBLE
  x AS LONG
  y AS LONG
  z AS LONG
END TYPE

DIM SHARED proiezione AS INTEGER, scala AS SINGLE 'proiezione: 0 lato, 1 alto
CONST mu0 = .000001257#     'permeabilit… magnetica
CONST mom = 8D+22           'momento magnetico terrestre - porlo = 0 per
CONST rt = 6000000          'raggio terrestre             ["spegnere" il campo
CONST enne = 30            'numero di particelle
CONST elec = 1.6D-19        'carica elementare
CONST elem = 9.1D-31        'massa dell'elettrone
CONST g = .000000000066#    'costante di gravitazione universale - porla = 0
CONST mt = 6D+24            'massa terrestre        [per "spegnere" la gravit…
CONST passo = .0003             'secondi per operazione

proiezione = 0: traccia = 0
DIM par(enne) AS partic, opar(enne) AS partic
DIM bb AS vet, aus AS vet, fza AS vet, ipl AS vet


SCREEN 12
font "Schema del campo magnetico terrestre in sezione verticale lungo un meridiano", 10, 20
PALETTE 9, 256 ^ 2 * 60 + 256 * 20
RANDOMIZE TIMER
scala = 1 / 400000' un pixel= 400 km
FOR i = 0 TO enne
  IF RND > .5 THEN par(i).car = -elec ELSE par(i).car = elec
  IF par(i).car < 0 THEN par(i).mas = elem ELSE par(i).mas = 1835 * elem
  par(i).vel.z = RND * 1000000000 + 100000
  IF par(i).car < o THEN par(i).vel.z = 5 * par(i).vel.z
  par(i).vel.x = RND * 8000 - 4000
  par(i).vel.y = RND * 8000 - 4000
  par(i).z = (-RND * 50000000 - 30000000) / 200000 / scala
  par(i).y = -RND * 20000000 + 10000000
  par(i).x = -RND * 20000000 + 10000000
  par(i).ene = (par(i).vel.x ^ 2 + par(i).vel.y ^ 2 + par(i).vel.z ^ 2) * par(i).mas / 2
  par(i).ene = par(i).ene - par(i).mas * mt * g / SQR(par(i).x ^ 2 + par(i).y ^ 2 + par(i).z ^ 2)
NEXT i
scala = 1 / 200000' un pixel= 200 km
CIRCLE (320, 240), INT(rt * scala), 9: PAINT (320, 240), 9
DIM z AS LONG, y AS LONG, x AS LONG
FOR z = -30000000 TO 30000000 STEP 1500000
  FOR y = -30000000 TO 30000000 STEP 1500000
    a$ = INKEY$
    IF SQR(z ^ 2 + y ^ 2) > 6000000 THEN
      campo 0, y, z, bb
      IF bb.y < 0 THEN dy = -2 ELSE dy = 2
      IF bb.z < 0 THEN dz = -2 ELSE dz = 2
      LINE (lax(x, y, z), lay(x, y, z))-(lax(x, y, z) + dz * LOG(ABS(bb.z) * 10 ^ 7 + 1), lay(x, y, z) + dy * LOG(ABS(bb.y) * 10 ^ 7 + 1)), 15
    END IF
    IF a$ = CHR$(27) THEN EXIT FOR
  NEXT y
  IF a$ = CHR$(27) THEN EXIT FOR
NEXT z
scala = 1 / 400000' un pixel=400 km
WHILE INKEY$ = "": WEND: CLS
font "Interazione delle particelle del vento solare con il campo magnetico terrestre", 10, 10
DRAW "c15"
font ". Elettroni", 550, 420
DRAW "c12"
font ". Protoni", 550, 440
DRAW "c14"
font "  Sole", 10, 470
DRAW "bm 10,465 ne4nf4r17"



proiezione = 2
CIRCLE (320, 240), INT(rt * scala), 9: PAINT (320, 240), 9
DO
  FOR i = 0 TO enne
              'calcolo forze e impulsi (su tempo doppio al passo)
  campo par(i).x, par(i).y, par(i).z, bb
 
  aus.x = par(i).vel.x * par(i).car
  aus.y = par(i).vel.y * par(i).car
  aus.z = par(i).vel.z * par(i).car
  prod aus, bb, fza
  ipl.x = fza.x * 2 * passo
  ipl.y = fza.y * 2 * passo
  ipl.z = fza.z * 2 * passo
             'gravit…, facoltativa
  f1 = g * mt * par(i).mas / (par(i).x ^ 2 + par(i).y ^ 2 + par(i).z ^ 2)
  fza.x = -f1 * par(i).x / SQR(par(i).x ^ 2 + par(i).y ^ 2 + par(i).z ^ 2)
  fza.y = -f1 * par(i).y / SQR(par(i).x ^ 2 + par(i).y ^ 2 + par(i).z ^ 2)
  fza.z = -f1 * par(i).z / SQR(par(i).x ^ 2 + par(i).y ^ 2 + par(i).z ^ 2)
  ipl.x = ipl.x + fza.x * 2 * passo
  ipl.y = ipl.y + fza.y * 2 * passo
  ipl.z = ipl.z + fza.z * 2 * passo
             
              'movimento
  opar(i) = par(i)
  par(i).x = par(i).x + par(i).vel.x * passo
  par(i).y = par(i).y + par(i).vel.y * passo
  par(i).z = par(i).z + par(i).vel.z * passo

              'applicazione impulso delle forze (eventuale controllo energetico)
  par(i).vel.x = par(i).vel.x + ipl.x / par(i).mas
  par(i).vel.y = par(i).vel.y + ipl.y / par(i).mas
  par(i).vel.z = par(i).vel.z + ipl.z / par(i).mas
  'energia
  energia = par(i).ene + par(i).mas * mt * g / SQR(par(i).x ^ 2 + par(i).y ^ 2 + par(i).z ^ 2)
  velocita2 = (2 * energia / par(i).mas)
  rap = SQR(velocita2 / (par(i).vel.x ^ 2 + par(i).vel.y ^ 2 + par(i).vel.z ^ 2))
  par(i).vel.x = par(i).vel.x * rap
  par(i).vel.y = par(i).vel.y * rap
  par(i).vel.z = par(i).vel.z * rap
              'movimento
  par(i).x = par(i).x + par(i).vel.x * passo
  par(i).y = par(i).y + par(i).vel.y * passo
  par(i).z = par(i).z + par(i).vel.z * passo
             
              'collisione!
  IF SQR(par(i).x ^ 2 + par(i).y ^ 2 + par(i).z ^ 2) < 6000000 THEN
    PALETTE 9, 60: PALETTE 0, 30
    par(i).vel.z = RND * 30000 + 10000
    IF par(i).car < o THEN par(i).vel.z = 10 * par(i).vel.z
    par(i).vel.x = RND * 800 - 400
    par(i).vel.y = RND * 800 - 400
    par(i).z = (-RND * 20000000 - 50000000) / 200000 / scala
    par(i).y = -RND * 30000000 + 15000000
    par(i).x = -RND * 30000000 + 15000000
    par(i).ene = (par(i).vel.x ^ 2 + par(i).vel.y ^ 2 + par(i).vel.z ^ 2) * par(i).mas / 2
    par(i).ene = par(i).ene - par(i).mas * mt * g / SQR(par(i).x ^ 2 + par(i).y ^ 2 + par(i).z ^ 2)
    PALETTE 9, 256 ^ 2 * 60 + 256 * 20: PALETTE 0, 0
  END IF
  IF SQR(par(i).x ^ 2 + par(i).y ^ 2 + par(i).z ^ 2) > 500 / scala THEN
    par(i).vel.z = RND * 1000000000 + 100000
    'IF par(i).car < o THEN par(i).vel.z = 10 * par(i).vel.z
    par(i).vel.x = RND * 800 - 400
    par(i).vel.y = RND * 800 - 400
    par(i).z = (-RND * 20000000 - 50000000) / 200000 / scala
    par(i).y = -RND * 30000000 + 15000000
    par(i).x = -RND * 30000000 + 15000000
    par(i).ene = (par(i).vel.x ^ 2 + par(i).vel.y ^ 2 + par(i).vel.z ^ 2) * par(i).mas / 2
    par(i).ene = par(i).ene - par(i).mas * mt * g / SQR(par(i).x ^ 2 + par(i).y ^ 2 + par(i).z ^ 2)
  END IF
              'disegno
  IF par(i).car > 0 THEN col = 12 ELSE col = 15
  PSET (lax(par(i).x, par(i).y, par(i).z), lay(par(i).x, par(i).y, par(i).z)), col
  IF (lax(opar(i).x, opar(i).y, opar(i).z) <> lax(par(i).x, par(i).y, par(i).z) OR lay(opar(i).x, opar(i).y, opar(i).z) <> lay(par(i).x, par(i).y, par(i).z)) AND traccia = 0 THEN
    IF proiezione = 0 THEN
      IF SQR(opar(i).z ^ 2 + opar(i).y ^ 2) < 6000000 THEN PSET (lax(opar(i).x, opar(i).y, opar(i).z), lay(opar(i).x, opar(i).y, opar(i).z)), 9 ELSE PSET (lax(opar(i).x, opar(i).y, opar(i).z), lay(opar(i).x, opar(i).y, opar(i).z)), 0
    END IF
    IF proiezione = 1 THEN
      IF SQR(opar(i).z ^ 2 + opar(i).x ^ 2) < 6000000 THEN PSET (lax(opar(i).x, opar(i).y, opar(i).z), lay(opar(i).x, opar(i).y, opar(i).z)), 9 ELSE PSET (lax(opar(i).x, opar(i).y, opar(i).z), lay(opar(i).x, opar(i).y, opar(i).z)), 0
    END IF
    IF proiezione = 2 THEN
      IF SQR((lax(opar(i).x, opar(i).y, opar(i).z) - 320) ^ 2 + (lay(opar(i).x, opar(i).y, opar(i).z) - 240) ^ 2) < 6000000 * scala THEN PSET (lax(opar(i).x, opar(i).y, opar(i).z), lay(opar(i).x, opar(i).y, opar(i).z)), 9 ELSE PSET (lax(opar(i).x,  _
opar(i).y, opar(i).z), lay(opar(i).x, opar(i).y, opar(i).z)), 0
    END IF
  END IF
  NEXT i
LOOP UNTIL INKEY$ = CHR$(27)

SUB campo (x AS LONG, y AS LONG, z AS LONG, a1 AS vet)
  IF ABS(x) < 100 AND ABS(z) < 100 THEN
    IF y > 0 THEN alfa = 0 ELSE alfa = 3.141592
    ELSE alfa = ATN(y / SQR(z ^ 2 + x ^ 2))
  END IF
  bt = mu0 * mom * COS(alfa) / 4 / 3.141592 / SQR(x ^ 2 + y ^ 2 + z ^ 2) ^ 3
  br = mu0 * mom * SIN(alfa) / 2 / 3.141592 / SQR(x ^ 2 + y ^ 2 + z ^ 2) ^ 3
  a1.y = SIN(alfa) * br - COS(alfa) * bt
  bxz = SIN(alfa) * bt - COS(alfa) * br
  IF ABS(x) < 100 AND ABS(z) < 100 THEN
    a1.z = 0: a1.x = 0
    ELSE a1.x = bxz * x / SQR(z ^ 2 + x ^ 2): a1.z = bxz * z / SQR(z ^ 2 + x ^ 2)
  END IF
END SUB

SUB font (s$, x, y)
ee = 0
FOR e = 1 TO LEN(s$)
a$ = MID$(s$, e, 1)
ee = ee + 1
'si inizia a scrivere la lettera dal pixel in basso a sinistra
'Lettere maiuscole e numeri: h=10
'Lettere minuscole: h=6
'Spessore=6
xx = x
yy = y
e$ = "bm+" + STR$(xx + ee * 8) + "," + STR$(yy)
DRAW "bm0,0"
DRAW e$
IF a$ = "0" THEN a$ = "O"
IF a$ = "A" THEN DRAW "m+3,-9 m+3,+9 bm-2,-3 l3"
IF a$ = "R" OR a$ = "P" OR a$ = "B" THEN DRAW "u9 r4 f1 d2 g1 nl4  "
IF a$ = "B" THEN DRAW "bd1br1d3bd1bl1l3"
IF a$ = "Q" OR a$ = "O" OR a$ = "C" OR a$ = "G" THEN DRAW "bu1nf1u7e1r3f1bd7g1l3"
IF a$ = "D" THEN DRAW "u9r4f1d7g1l4"
IF a$ = "E" THEN DRAW "nr5u5nr4u4r5"
IF a$ = "F" THEN DRAW "u5nr4u4r5"
IF a$ = "G" THEN DRAW "r3e1u2l1"
IF a$ = "W" THEN DRAW "br1 nm-2,-9 m+2,-9  m+2,+9 m+2,-9": ee = ee + .25
IF a$ = ")" THEN DRAW "br1nd1e1u1e1u4h1u1h1u1": ee = ee - .375
IF a$ = "(" THEN DRAW "br1nd1h1u1h1u4e1u1e1u1": ee = ee - .375
IF a$ = "H" THEN DRAW "u9d4r5u4d9"
IF a$ = "I" THEN DRAW "r2l1u9r1l2 ": ee = ee - .375
IF a$ = "L" THEN DRAW "nu9r5"
IF a$ = "M" THEN DRAW "u9f3e3d9": ee = ee + .125
IF a$ = "N" THEN DRAW "u9 m+5,+9 u9"
IF a$ = "O" OR a$ = "Q" THEN DRAW "nh1r3e1nu7"
IF a$ = "." THEN DRAW "u1": ee = ee - .625
IF a$ = "Q" THEN DRAW "bd1bl1d1u2h1"
IF a$ = "J" THEN DRAW "bu1f1r3e1u8"
IF a$ = "K" THEN DRAW "u9d5nm+5,-5u1nm+5,+5"
IF a$ = "R" THEN DRAW "l2 m+3,+5"
IF a$ = "S" THEN DRAW "bu1f1r3e1u3h1l3h1u2e1r3f1"
IF a$ = "T" THEN DRAW "br3u9 nr3nl3"
IF a$ = "U" THEN DRAW "bu9 d8 m+1,+1r3m+1,-1 u8"
IF a$ = "V" THEN DRAW "bu9 d1 m+3,+8m+3,-9"
IF a$ = "Y" THEN DRAW "br3u4 nm+3,-5m-3,-5 "
IF a$ = "Z" THEN DRAW "nr5m+5,-9l5"
IF a$ = "X" THEN DRAW "m+5,-9bd9m-5,-9"
IF a$ = "-" THEN DRAW "bu5 r5"
IF a$ = "?" THEN DRAW "br3 u1 bu2 u1e2u1h2l2g1"
IF a$ = "!" THEN DRAW "u1 bu2 u6": ee = ee - .625
IF a$ = "," THEN DRAW "br2g1": ee = ee - .5
IF a$ = ";" THEN DRAW "br2ng1bu3u1": ee = ee - .5
IF a$ = ":" THEN DRAW "u1 bu3 u1": ee = ee - .625
IF a$ = "q" THEN DRAW "bu2u3e1r3f1d9u5g1l3h1"
IF a$ = "w" THEN DRAW "bu1nu5f1r1e1nu5f1r1e1nu5": ee = ee + .125
IF a$ = "e" THEN DRAW "bu1u4e1r3f1d2l5d2f1r3e1"
IF a$ = "r" THEN DRAW "u6d1r1e1r2f1"
IF a$ = "t" THEN DRAW "br5bu1g1l3h1u6nr3u2"
IF a$ = "y" THEN DRAW "bu1nu5f1r4nu6d3g1l3"
IF a$ = "u" THEN DRAW "bu1nu5f1r3e1u5"
IF a$ = "i" THEN DRAW "u6bu2u1": ee = ee - .625
IF a$ = "'" THEN DRAW "bm+0,-9d1": ee = ee - .625
IF a$ = "o" THEN DRAW "bu1u4e1r3f1d4g1l3h1"
IF a$ = "p" THEN DRAW "nd4u5e1r3f1d4g1l4"
IF a$ = "a" THEN DRAW "bu1u1e1r4d2nf1g1l3nh1 r3e1u4h1l3g1"
IF a$ = "s" THEN DRAW "bu1f1r3e1u1h1l3h1u1e1r3f1"
IF a$ = "d" THEN DRAW "bu1u4e1r4nu3d5g1l3h1"
IF a$ = "f" THEN DRAW "u5nr3u3e1r3f1"
IF a$ = "g" THEN DRAW "bu1u4e1r3f1d5l4nh1r4d3g1l3"
IF a$ = "h" THEN DRAW "nu9u5r1e1r2f1d5"
IF a$ = "j" THEN DRAW "bd3f1r3e1u9bu2u1"
IF a$ = "k" THEN DRAW "nu9u3r1ne2nf3": ee = ee - .125
IF a$ = "l" THEN DRAW "nu9": ee = ee - .625
IF a$ = "m" THEN DRAW "u6d1r1e1r2f1nd5r1e1r2f1d5": ee = ee + .625
IF a$ = "n" THEN DRAW "u6d1r1e1r2f1d5"
IF a$ = "b" THEN DRAW "bu1nu8f1r3e1u4h1l3"
IF a$ = "v" THEN DRAW "bu6d1f1d1f1d1f1e1u1e1u1e1u1": ee = ee + .125
IF a$ = "c" THEN DRAW "br5bu1g1l3h1u4e1r3f1"
IF a$ = "x" THEN DRAW "e3nh3l1ne3f3"
IF a$ = "z" THEN DRAW "nr5u1e5l5"
IF a$ = "•" THEN DRAW "bu1u4e1r3f1d4g1l3h1u4e1r2bu2h1"
IF a$ = "…" THEN DRAW "bu1u1e1r4d2nf1g1l3nh1r3e1u4h1l3ng1r2bu2h1"
IF a$ = "—" THEN DRAW "bu1nu5f1r3e1u5bl2bu2h1"
IF a$ = "Š" THEN DRAW "bu1u4e1r3f1d2l5d2f1r3e1bu7bl3h1"
IF a$ = "Ť" THEN DRAW "u6bu2h1": ee = ee - .625
IF a$ = "1" THEN DRAW "r2l1u9g1": ee = ee - .375
IF a$ = "2" THEN DRAW "nr5u1e5u2h1l3g1"
IF a$ = "3" THEN DRAW "bu1f1r3e1u3h1n l2 e1u2h1l3g1"
IF a$ = "4" THEN DRAW "bu2nr5u1e1u1e1u1e1e1d9"
IF a$ = "5" THEN DRAW "bu1f1r3e1u4h1l4nd1u3r5"
IF a$ = "6" THEN DRAW "bu1nu5f1r3e1u3h1l2g1l1u4e1r3f1"
IF a$ = "7" THEN DRAW "br2u3e1u2e1u1e1l5"
IF a$ = "8" THEN DRAW "bu1nf1u3e1nr3h1u2e1r3f1d2g1f1d3g1l3"
IF a$ = "9" THEN DRAW "bu1f1r3e1u7h1l3g1d3f1r2e1r1"
IF a$ = "" THEN ee = ee - .25
NEXT e


END SUB

FUNCTION lax (x AS LONG, y AS LONG, z AS LONG)
  IF proiezione = 0 THEN lax = INT(scala * z + 320)
  IF proiezione = 1 THEN lax = INT(scala * z + 320)
  IF proiezione = 2 THEN lax = INT(scala * (z + x * (SQR(2) / 2))) + 320
   'xx(a) = (z(a) + x(a) / (2 * SQR(2)))  'CAVALIERA
   'yy(a) = (x(a) / (2 * SQR(2)) + y(a))

END FUNCTION

FUNCTION lay (x AS LONG, y AS LONG, z AS LONG)
  IF proiezione = 0 THEN lay = -INT(scala * y) + 240
  IF proiezione = 1 THEN lay = -INT(scala * x) + 240
  IF proiezione = 2 THEN lay = -INT(scala * (y + x * (SQR(2) / 2))) + 240
END FUNCTION

SUB prod (a AS vet, b AS vet, c AS vet)
  c.x = a.y * b.z - b.y * a.z
  c.y = a.z * b.x - b.z * a.x
  c.z = a.x * b.y - b.x * a.y
END SUB