'ORB3.BAS: third module of ORBWEAVE.BAS
'by Peter Bird, UCLA, August 1993-March 2000.
'(c) Copyright 1993, 2000 by Peter Bird and the Regents of the University of California.
DECLARE SUB ABG2XY (Ins!, TempVec!(), Outs!, Visible%, x!, Y!)
DECLARE SUB Blanker (r1%, R2%, c1%, c2%)
DECLARE SUB Cross (Ins, V1!(), V2!(), Outs, VCross!())
DECLARE SUB Divide (X1#, Y1#, z1#, X2#, Y2#, z2#, X3#, Y3#, z3#, NSlice%)
DECLARE SUB ElementOut (X1#, Y1#, z1#, X2#, Y2#, z2#, X3#, Y3#, z3#)
DECLARE SUB GetFileName (Text$, LineN%, NewName$)
DECLARE SUB LineOfNumbers (Ins!, Rec$, Outs!, AreNumbers%)
DECLARE SUB Pixels (Ins!, x!, Y!, Outs!, col%, row%)
DECLARE SUB RaiseOne (Ins!, VCol%(), VRow%())
DECLARE SUB MouseHide ()
DECLARE SUB MouseShow ()
DECLARE SUB Unitize (Mods, x#, Y#, Z#)
DECLARE SUB WaitForKey (a$)
DECLARE SUB XY2ABG (Ins!, x!, Y!, Outs!, Outside%, TempVec!())
'
OPTION BASE 1
DIM WinOut!(3), WinRight!(3), WinUp!(3)
DIM Cart2!(3, 3)
DIM FIPoint!(4)
DIM Icon1%(30), Icon2%(30)
DIM scales!(2, 2)
DIM UnScale!(2, 2)
COMMON SHARED Cart2!(), Colored%, CursorOn%, DoIcon%, FEP%, FIPoint!(), HiColor%, HiCol%, HiRow%, Icon1%(), Icon2%(), LastPutC%, LastPutR%, R2C!, scales!(), UnScale!(), WindowHeight!, WinOut!(), WinRight!(), WinUp!()

SUB CloseOut (Ins, GridLoaded%, Remind%)
'terminates program (when no grid is loaded)
	IF GridLoaded% THEN
		BEEP
	ELSE  'avoid orphan binary files
		IF LEN(DIR$("ELEMENTS.BIN")) > 0 THEN KILL "ELEMENTS.BIN"
		IF LEN(DIR$("BADNODES.BIN")) > 0 THEN
			CLOSE #3
			KILL "BADNODES.BIN"
		END IF
		IF Remind% THEN
			CLS
			LOCATE 10, 10
			PRINT "Since you have created or modified a grid, remember that"
			LOCATE 12, 10
			PRINT "it is necessary to renumber the nodes with utility program"
			LOCATE 14, 10
			PRINT "ORBNUMBR before using the grid in finite element program SHELLS."
			LOCATE 24, 60: PRINT "Press any key ...";
			CALL WaitForKey(a$)
		END IF
		CLS
		SYSTEM
	END IF
END SUB

SUB Divide (X1#, Y1#, z1#, X2#, Y2#, z2#, X3#, Y3#, z3#, NSlice%)
'Accepts 3 unit vectors of triangle corners as input, subdivides the
'triangle NSlice% times into 4 triangles, and outputs results to a file.
'NOTE that Divide calls itself recursively, so all the arguments should
'be simple numbers passed by value, not by address.
'To avoid complications, output is dumped into an open binary file (#5).
xe1# = .5# * (X2# + X3#)
xe2# = .5# * (X3# + X1#)
xe3# = .5# * (X1# + X2#)
ye1# = .5# * (Y2# + Y3#)
ye2# = .5# * (Y3# + Y1#)
ye3# = .5# * (Y1# + Y2#)
ze1# = .5# * (z2# + z3#)
ze2# = .5# * (z3# + z1#)
ze3# = .5# * (z1# + z2#)
CALL Unitize(Mods, xe1#, ye1#, ze1#)
CALL Unitize(Mods, xe2#, ye2#, ze2#)
CALL Unitize(Mods, xe3#, ye3#, ze3#)
IF NSlice% > 1 THEN
	MSlice% = NSlice% - 1
	CALL Divide(X1#, Y1#, z1#, xe2#, ye2#, ze2#, xe3#, ye3#, ze3#, MSlice%)
	CALL Divide(X2#, Y2#, z2#, xe1#, ye1#, ze1#, xe3#, ye3#, ze3#, MSlice%)
	CALL Divide(X3#, Y3#, z3#, xe1#, ye1#, ze1#, xe2#, ye2#, ze2#, MSlice%)
	CALL Divide(xe1#, ye1#, ze1#, xe2#, ye2#, ze2#, xe3#, ye3#, ze3#, MSlice%)
ELSEIF NSlice% = 1 THEN
	CALL ElementOut(X1#, Y1#, z1#, xe2#, ye2#, ze2#, xe3#, ye3#, ze3#)
	CALL ElementOut(X2#, Y2#, z2#, xe1#, ye1#, ze1#, xe3#, ye3#, ze3#)
	CALL ElementOut(X3#, Y3#, z3#, xe1#, ye1#, ze1#, xe2#, ye2#, ze2#)
	CALL ElementOut(xe1#, ye1#, ze1#, xe2#, ye2#, ze2#, xe3#, ye3#, ze3#)
ELSE 'Nslice% =0
	CALL ElementOut(X1#, Y1#, z1#, X2#, Y2#, z2#, X3#, Y3#, z3#)
END IF
END SUB

SUB DrawElement (BotF!, ColorIn%, Contour%, DFC!, EMemo%(), EQCM!(), i%, IData%, Logs%, Mus!(), NODES%(), Redo%, UseColor%, NodeABG!())
'Draw ONE element (i%), unless it wraps over the horizon.
'Switch Contour% = T will contour the interior using data from:
'      EQCM!(IData%, NODES%(1-3, i%))  [ if IData% = 1, 2, 3, 4 ]; OR
'      Mus!(IData%-4, i%)              [ if IData% = 5, 6, 7 ].
'Or, use switch ColorIn% = T to raise the color index of each pixel in the triangle by one.
DIM CList%(6), ECol%(3), ERow%(3), RList%(6), TempVec!(3), Xt!(3), Yt!(3)
K1% = NODES%(1, i%)
TempVec!(1) = NodeABG!(1, K1%): TempVec!(2) = NodeABG!(2, K1%): TempVec!(3) = NodeABG!(3, K1%)
CALL ABG2XY(Ins, TempVec!(), Outs, Visible%, Xt!(1), Yt!(1))
IF Visible% THEN  ELSE EXIT SUB
R21! = Xt!(1) * Xt!(1) + Yt!(1) * Yt!(1)
K2% = NODES%(2, i%)
TempVec!(1) = NodeABG!(1, K2%): TempVec!(2) = NodeABG!(2, K2%): TempVec!(3) = NodeABG!(3, K2%)
CALL ABG2XY(Ins, TempVec!(), Outs, Visible%, Xt!(2), Yt!(2))
IF Visible% THEN  ELSE EXIT SUB
R22! = Xt!(2) * Xt!(2) + Yt!(2) * Yt!(2)
K3% = NODES%(3, i%)
TempVec!(1) = NodeABG!(1, K3%): TempVec!(2) = NodeABG!(2, K3%): TempVec!(3) = NodeABG!(3, K3%)
CALL ABG2XY(Ins, TempVec!(), Outs, Visible%, Xt!(3), Yt!(3))
IF Visible% THEN  ELSE EXIT SUB
R23! = Xt!(3) * Xt!(3) + Yt!(3) * Yt!(3)
IF (R21! < R2C!) OR (R22! < R2C!) OR (R23! < R2C!) THEN
	FOR K% = 1 TO 3
		CALL Pixels(Ins, Xt!(K%), Yt!(K%), Outs, ECol%(K%), ERow%(K%))
	NEXT K%
	IF Redo% THEN  'blank triangular area before redrawing
		MouseHide
		LINE (ECol%(1), ERow%(1))-(ECol%(2), ERow%(2)), 1
		LINE (ECol%(2), ERow%(2))-(ECol%(3), ERow%(3)), 1
		LINE (ECol%(3), ERow%(3))-(ECol%(1), ERow%(1)), 1
		SeedC% = (ECol%(1) + ECol%(2) + ECol%(3)) / 3
		SeedR% = (ERow%(1) + ERow%(2) + ERow%(3)) / 3
		PAINT (SeedC%, SeedR%), 1, 1
		MouseShow
	END IF
	IF Contour% THEN
		  MouseHide
		  IF (IData% > 4) THEN
			f1! = Mus!(IData% - 4, i%)
			f2! = Mus!(IData% - 4, i%)
			f3! = Mus!(IData% - 4, i%)
		  ELSE
			f1! = EQCM!(IData%, K1%)
			f2! = EQCM!(IData%, K2%)
			f3! = EQCM!(IData%, K3%)
		  END IF
		  IF Logs% THEN
			IF f1! > 0! THEN
				f1! = LOG(f1!) / LOG(10#)
			ELSE
				f1! = BotF!
			END IF
			IF f2! > 0! THEN
				f2! = LOG(f2!) / LOG(10#)
			ELSE
				f2! = BotF!
			END IF
			IF f3! > 0! THEN
				f3! = LOG(f3!) / LOG(10#)
			ELSE
				f3! = BotF!
			END IF
		  END IF
		  fLow! = f1!
		  IF f2! < fLow! THEN fLow! = f2!
		  IF f3! < fLow! THEN fLow! = f3!
		  fHigh! = f1!
		  IF f2! > fHigh! THEN fHigh! = f2!
		  IF f3! > fHigh! THEN fHigh! = f3!
		  ILOW% = INT(2.01 + (fLow! - BotF!) / DFC!)
		  IHIGH% = INT(1.99 + (fHigh! - BotF!) / DFC!)
		  IF (ILOW% >= IHIGH%) OR (ILOW% > HiColor%) OR (IHIGH% < 2) THEN
'                UNIFORM COLOR; NO CONTOURS
			  CList%(1) = ECol%(1): CList%(2) = ECol%(2): CList%(3) = ECol%(3)
			  RList%(1) = ERow%(1): RList%(2) = ERow%(2): RList%(3) = ERow%(3)
			  NPNTS% = 3
			  index% = IHIGH%
			  GOSUB Shade
		  ELSE
'                CONTOUR(S) MUST BE LOCATED
			  NUMCON% = IHIGH% - ILOW%
'                FIRST, ROTATE TO STANDARD VIEW
			  IF (f1! = fLow!) THEN
				  XLOW% = ECol%(1)
				  YLOW% = ERow%(1)
				  IF (f2! <= f3!) THEN
					  FMID! = f2!
					  XMID% = ECol%(2)
					  YMID% = ERow%(2)
					  XHIGH% = ECol%(3)
					  YHIGH% = ERow%(3)
				  ELSE
					  FMID! = f3!
					  XMID% = ECol%(3)
					  YMID% = ERow%(3)
					  XHIGH% = ECol%(2)
					  YHIGH% = ERow%(2)
				  END IF
			  ELSEIF (f2! = fLow!) THEN
				  XLOW% = ECol%(2)
				  YLOW% = ERow%(2)
				  IF (f1! <= f3!) THEN
					  FMID! = f1!
					  XMID% = ECol%(1)
					  YMID% = ERow%(1)
					  XHIGH% = ECol%(3)
					  YHIGH% = ERow%(3)
				  ELSE
					  FMID! = f3!
					  XMID% = ECol%(3)
					  YMID% = ERow%(3)
					  XHIGH% = ECol%(1)
					  YHIGH% = ERow%(1)
				  END IF
			  ELSE
'                    (F3 IS THE LOW VALUE)
				  XLOW% = ECol%(3)
				  YLOW% = ERow%(3)
				  IF (f1! <= f2!) THEN
					  FMID! = f1!
					  XMID% = ECol%(1)
					  YMID% = ERow%(1)
					  XHIGH% = ECol%(2)
					  YHIGH% = ERow%(2)
				  ELSE
					  FMID! = f2!
					  XMID% = ECol%(2)
					  YMID% = ERow%(2)
					  XHIGH% = ECol%(1)
					  YHIGH% = ERow%(1)
				  END IF
			  END IF
			  XLOW! = XLOW%
			  YLOW! = YLOW%
			  XMID! = XMID%
			  YMID! = YMID%
			  XHIGH! = XHIGH%
			  YHIGH! = YHIGH%
			  FOR J% = 1 TO NUMCON%
				  IIN% = ILOW% + J% - 1
				  IOUT% = IIN% + 1
				  FCON! = BotF! + (IIN% - 1) * DFC!
				  SBASE! = (FCON! - fLow!) / (fHigh! - fLow!)
				  XBASE% = CINT(XLOW! + SBASE! * (XHIGH! - XLOW!))
				  YBASE% = CINT(YLOW! + SBASE! * (YHIGH! - YLOW!))
				  IF (FCON! < FMID!) THEN
					inner% = (NOT 0) 'true
				  ELSE
					inner% = 0       'false
				  END IF
				  IF (inner%) THEN
					  SFAR! = (FCON! - fLow!) / (FMID! - fLow!)
					  XFAR% = CINT(XLOW! + SFAR! * (XMID! - XLOW!))
					  YFAR% = CINT(YLOW! + SFAR! * (YMID! - YLOW!))
				  ELSE
					  SFAR! = (FCON! - FMID!) / (fHigh! - FMID!)
					  XFAR% = CINT(XMID! + SFAR! * (XHIGH! - XMID!))
					  YFAR% = CINT(YMID! + SFAR! * (YHIGH! - YMID!))
				  END IF
				  IF (J% = 1) THEN
'                          FIRST POLYGON
					  CList%(1) = XLOW%
					  RList%(1) = YLOW%
					  CList%(2) = XBASE%
					  RList%(2) = YBASE%
					  CList%(3) = XFAR%
					  RList%(3) = YFAR%
					  IF (inner%) THEN
						  NPNTS% = 3
					  ELSE
						  NPNTS% = 4
						  CList%(4) = XMID%
						  RList%(4) = YMID%
					  END IF
					  index% = IIN%
					  GOSUB Shade
				  ELSE
'                          ONE MORE POLYGON IN A SERIES
					  CList%(1) = XOBASE%
					  RList%(1) = YOBASE%
					  CList%(2) = XBASE%
					  RList%(2) = YBASE%
					  CList%(3) = XFAR%
					  RList%(3) = YFAR%
					  IF (inner% = INOLD%) THEN
						  NPNTS% = 4
						  CList%(4) = XOFAR%
						  RList%(4) = YOFAR%
					  ELSE
						  NPNTS% = 5
						  CList%(4) = XMID%
						  RList%(4) = YMID%
						  CList%(5) = XOFAR%
						  RList%(5) = YOFAR%
					  END IF
					  index% = IIN%
					  GOSUB Shade
				  END IF
				  IF (J% = NUMCON%) THEN
'                          LAST POLYGON TO FILL IN
					  CList%(1) = XBASE%
					  RList%(1) = YBASE%
					  CList%(2) = XHIGH%
					  RList%(2) = YHIGH%
					  IF (inner%) THEN
						  NPNTS% = 4
						  CList%(3) = XMID%
						  RList%(3) = YMID%
						  CList%(4) = XFAR%
						  RList%(4) = YFAR%
					  ELSE
						  NPNTS% = 3
						  CList%(3) = XFAR%
						  RList%(3) = YFAR%
					  END IF
					  index% = IOUT%
					  GOSUB Shade
				  END IF
				  XOBASE% = XBASE%
				  YOBASE% = YBASE%
				  XOFAR% = XFAR%
				  YOFAR% = YFAR%
				  INOLD% = inner%
			  NEXT J%
		  END IF
		  MouseShow
	END IF 'Contour%
	IF Colored% THEN COLOR UseColor%
	IF Colored% THEN
		LINE (ECol%(1), ERow%(1))-(ECol%(2), ERow%(2))
		LINE (ECol%(2), ERow%(2))-(ECol%(3), ERow%(3))
		LINE (ECol%(3), ERow%(3))-(ECol%(1), ERow%(1))
	ELSE
		LINE (ECol%(1), ERow%(1))-(ECol%(2), ERow%(2)), 0
		LINE (ECol%(2), ERow%(2))-(ECol%(3), ERow%(3)), 0
		LINE (ECol%(3), ERow%(3))-(ECol%(1), ERow%(1)), 0
		LINE (ECol%(1), ERow%(1))-(ECol%(2), ERow%(2)), , , 3855
		LINE (ECol%(2), ERow%(2))-(ECol%(3), ERow%(3)), , , 3855
		LINE (ECol%(3), ERow%(3))-(ECol%(1), ERow%(1)), , , 3855
	END IF 'Colored% or not
	'add a schematic triangle at the center point, to click on
	ColC% = ((ECol%(1) + ECol%(2) + ECol%(3)) / 3) - 3
	RowC% = ((ERow%(1) + ERow%(2) + ERow%(3)) / 3) - 3
	IF DoIcon% THEN
		IF ColC% >= 0 AND ColC% <= LastPutC% THEN
			IF RowC% >= 0 AND RowC% <= LastPutR% THEN
				IF (EMemo%(i%) AND 1) THEN 'X for flipped element
					IF UseColor% THEN
						PUT (ColC%, RowC%), Icon2%, PSET
					ELSE
						PUT (ColC%, RowC%), Icon2%, XOR
					END IF
				ELSE 'triangle for normal element
					IF UseColor% THEN
						PUT (ColC%, RowC%), Icon1%, PSET
					ELSE
						PUT (ColC%, RowC%), Icon1%, XOR
					END IF
				END IF
			END IF
		END IF
	END IF 'DoIcon%
	IF ColorIn% THEN
		IF (EMemo%(i%) AND 1) THEN  'element flipped; reorder before call
			J% = ECol%(2): ECol%(2) = ECol%(3): ECol%(3) = J%
			J% = ERow%(2): ERow%(2) = ERow%(3): ERow%(3) = J%
		END IF
		CALL RaiseOne(Ins, ECol%(), ERow%())
	END IF 'ColorIn%
END IF
EXIT SUB
'==============================
Shade:    '(Ins, NPNTS%, CList%(), RList%(), Mods, INDEX% )
	IF index% < 2 THEN index% = 0
	IF index% > HiColor% THEN index% = 0
	'draw "safety net" around whole element
	LINE (ECol%(1), ERow%(1))-(ECol%(2), ERow%(2)), index%
	LINE (ECol%(2), ERow%(2))-(ECol%(3), ERow%(3)), index%
	LINE (ECol%(3), ERow%(3))-(ECol%(1), ERow%(1)), index%
	SeedC! = 0!: SeedR! = 0!
	R2! = 0!
	area! = 0!
	FOR jt% = 1 TO NPNTS%
		SeedC! = SeedC! + CSNG(CList%(jt%))
		SeedR! = SeedR! + CSNG(RList%(jt%))
		IF jt% < NPNTS% THEN
			jp% = jt% + 1
		ELSE
			jp% = 1
		END IF
		c1% = CList%(jt%): c2% = CList%(jp%)
		r1% = RList%(jt%): R2% = RList%(jp%)
		LINE (c1%, r1%)-(c2%, R2%), index%
		t! = CSNG(c1% - c2%) ^ 2 + CSNG(r1% - R2%) ^ 2
		IF t! > R2! THEN R2! = t!
		area! = area! + .5 * CSNG(r1% + R2%) * CSNG(c2% - c1%)
	NEXT jt%
	IF area! < 0! THEN area! = -area!
	Longest! = SQR(R2!)
	IF Longest! < 1! THEN Longest! = 1!
	Wide! = area! / Longest!
	NPNTS! = CSNG(NPNTS%)
	SeedC% = CINT(SeedC! / NPNTS!)
	SeedR% = CINT(SeedR! / NPNTS!)
	Hue% = POINT(SeedC%, SeedR%)
	IF (Hue% <> index%) THEN
		IF Wide! > 1.9 THEN
			PAINT (SeedC%, SeedR%), index%, index%
		END IF
	END IF
RETURN
'===============================================================
END SUB

SUB ElementOut (X1#, Y1#, z1#, X2#, Y2#, z2#, X3#, Y3#, z3#)
'Outputs an element defined by 3 unit radius vectors to its vertices.
'The output (binary) file must already be open, as device !5.
X1! = CSNG(X1#)
Y1! = CSNG(Y1#)
z1! = CSNG(z1#)
X2! = CSNG(X2#)
Y2! = CSNG(Y2#)
z2! = CSNG(z2#)
X3! = CSNG(X3#)
Y3! = CSNG(Y3#)
z3! = CSNG(z3#)
PUT #5, , X1!
PUT #5, , Y1!
PUT #5, , z1!
PUT #5, , X2!
PUT #5, , Y2!
PUT #5, , z2!
PUT #5, , X3!
PUT #5, , Y3!
PUT #5, , z3!
END SUB

SUB GetBase (Ins, CD$, Outs, BinaryN$, NeedMenu%, NeedToDraw%, UsingBase%)
'gets name of basemap file, and creates binary file if needed.
FALSE% = 0
TRUE% = NOT FALSE%
NeedToDraw% = 0
B1:
	FileN$ = CD$ + "\*.DIG"
	IF LEN(DIR$(FileN$)) = 0 THEN
		CALL Blanker(16, 24, 2, 79)
		LOCATE 19, 10
		PRINT "There are no .DIG files in the current directory.";
		LOCATE 20, 10
		PRINT "You may wish to change directories with command D.";
		LOCATE 21, 10
		PRINT "Press Enter to continue."
		CALL WaitForKey(a$)
	ELSE
		CALL Blanker(16, 24, 2, 79)
		LOCATE 16, 10: PRINT "The following .DIG files are available in the current directory,";
		LOCATE 17, 3: FILES FileN$
		NeedMenu% = NOT 0 '(file listing may scroll menu up!)
		CALL Blanker(23, 24, 2, 79)
		LOCATE 23, 10
		PRINT "If you do not want any of these basemap files, press ESCape now.";
		LOCATE 24, 10
		PRINT "If you do, press Enter to continue ...";
		CALL WaitForKey(a$)
		IF a$ = CHR$(27) THEN
			UsingBase% = FALSE%
		ELSE
			CALL Blanker(23, 24, 3, 78)
			a$ = "Enter an existing filename (1-8 characters, no extension): "
			CALL GetFileName(a$, 23, NewName$)
			BaseMapN$ = CD$ + "\" + NewName$ + ".DIG"
			IF LEN(DIR$(BaseMapN$)) = 0 THEN
				BEEP
				GOTO B1
			ELSE
				BinaryN$ = CD$ + "\" + NewName$ + ".BIN"
				IF (LEN(DIR$(BinaryN$)) = 0) OR (NewName$ = "BADNODES") THEN 'create a binary file
					' Nature of endmark must be inferred:
					SeeEWidth% = 3      'default, if none found
					SeeEndMark$ = "***" 'default, if none found
					CLOSE #2
					OPEN BaseMapN$ FOR INPUT AS #2
					PastTitles% = FALSE%
					DO UNTIL EOF(2)
						LINE INPUT #2, Rec$
						CALL LineOfNumbers(Ins, Rec$, Outs, AreNumbers%)
						IF AreNumbers% THEN
							PastTitles% = TRUE%
						ELSE
							IF PastTitles% THEN
								SeeEndMark$ = Rec$
								EXIT DO
							END IF
						END IF
					LOOP
					SeeEWidth% = LEN(SeeEndMark$)
					IF SeeEWidth% > 4 THEN
						SeeEWidth% = 4
						SeeEndMark$ = LEFT$(SeeEndMark$, 4)
					END IF
					CLOSE #2
					'convert to binary (Alpha,Beta,Gamma) format:
					OPEN BaseMapN$ FOR INPUT AS #2
					OPEN BinaryN$ FOR BINARY AS #3
					np& = 0
					CALL Blanker(23, 24, 2, 79)
					LOCATE 23, 19: PRINT "Making a binary-file copy to speed access...";
					LOCATE 24, 26: PRINT "Number of points read =";
					DO UNTIL EOF(2) 'loop on line segments
						'work past any titles (any number of titles):
						DO UNTIL EOF(2) 'loop on title lines
							LINE INPUT #2, Rec$
							CALL LineOfNumbers(Ins, Rec$, Outs, AreNumbers%)
							IF AreNumbers% THEN EXIT DO
						LOOP
						'Rec$ should now contain first numbers:
						DO  'loop on number pairs
							IF LEFT$(Rec$, SeeEWidth%) = SeeEndMark$ THEN
								Alpha! = 99! 'impossible value used as endmark
								Beta! = 99!
								Gamma! = 99!
								PUT #3, , Alpha!
								PUT #3, , Beta!
								PUT #3, , Gamma!
								EXIT DO
							ELSE
								np& = np& + 1&
								LOCATE 24, 50: PRINT np&;
								Lon! = VAL(MID$(Rec$, 2, 12))
								Lat! = VAL(MID$(Rec$, 15, 12))
								Lon! = Lon! * .017453293#
								Lat! = Lat! * .017453292#
								CosLat! = COS(Lat!)
								Alpha! = CosLat! * COS(Lon!)
								Beta! = CosLat! * SIN(Lon!)
								Gamma! = SIN(Lat!)
								PUT #3, , Alpha!
								PUT #3, , Beta!
								PUT #3, , Gamma!
							END IF
							'prepare to loop by reading next line:
							IF EOF(2) THEN EXIT DO
							LINE INPUT #2, Rec$
						LOOP  'until this line segment is completed
					LOOP 'until all line segments are read
					CLOSE #2
					CLOSE #3
				END IF
			END IF
			UsingBase% = NOT 0
			NeedToDraw% = NOT 0
		END IF
	END IF
END SUB

SUB GetCD (Outs, CD$)
'uses DOS shell to determine current directory
CLS
LOCATE 12, 30
PRINT "One moment, please ..."
SHELL "CHDIR > TEMP9375.284"
OPEN "TEMP9375.284" FOR INPUT AS #1
LINE INPUT #1, CD$
CLOSE #1
KILL "TEMP9375.284"
CD$ = UCASE$(CD$)
IF RIGHT$(CD$, 1) = "\" THEN CD$ = LEFT$(CD$, LEN(CD$) - 1)
END SUB

SUB GetElement (Ins, NodeABG!(), NODES%(), NUMEL%, x!, Y!, Outs, J%, R2MIN!)
'finds triangular continuum element closest to (X!,Y!)
DIM TempVec!(3), TempV2!(3), TempV3!(3)
R2MIN! = 3.3E+38
FOR i% = 1 TO NUMEL%  'search for closest element centroid (in projection)
	n1% = NODES%(1, i%): n2% = NODES%(2, i%): n3% = NODES%(3, i%)
	TempVec!(1) = NodeABG!(1, n1%): TempVec!(2) = NodeABG!(2, n1%): TempVec!(3) = NodeABG!(3, n1%)
	CALL ABG2XY(Ins, TempVec!(), Outs, Visible%, X1!, Y1!)
	IF Visible% THEN
		TempV2!(1) = NodeABG!(1, n2%): TempV2!(2) = NodeABG!(2, n2%): TempV2!(3) = NodeABG!(3, n2%)
		CALL ABG2XY(Ins, TempV2!(), Outs, Visible%, X2!, Y2!)
		IF Visible% THEN
			TempV3!(1) = NodeABG!(1, n3%): TempV3!(2) = NodeABG!(2, n3%): TempV3!(3) = NodeABG!(3, n3%)
			CALL ABG2XY(Ins, TempV3!(), Outs, Visible%, X3!, Y3!)
			IF Visible% THEN
				XC! = (X1! + X2! + X3!) / 3!
				YC! = (Y1! + Y2! + Y3!) / 3!
				DX! = x! - XC!: DY! = Y! - YC!
				R2! = DX! * DX! + DY! * DY!
				IF R2! < R2MIN! THEN
					R2MIN! = R2!: J% = i%
				END IF
			END IF
		END IF
	END IF
NEXT i%
END SUB

SUB GetNewCD (Mods, CD$, NSlice%)
	LOCATE 23, 5
	PRINT "Current directory is now set to:";
	Blanks$ = "                                        "
	Pad% = 40 - LEN(CD$)
	a$ = CD$ + RIGHT$(Blanks$, Pad%)
	LOCATE 23, 38
	PRINT a$;
	LOCATE 24, 5
	PRINT "(Return/Enter to accept, or type preferred path)";
	LOCATE 23, 38
	CALL WaitForKey(a$)
	IF a$ = CHR$(13) THEN
	ELSE
		LOCATE 23, 38
		PRINT Blanks$
		LOCATE 23, 38
		PRINT a$;
		a$ = UCASE$(a$)
		INPUT ; "", b$
		TrialCD$ = a$ + b$
		IF RIGHT$(TrialCD$, 1) = ":" THEN TrialCD$ = TrialCD$ + "\"
		CHDIR TrialCD$
		IF MID$(TrialCD$, 2, 1) = ":" THEN CHDRIVE TrialCD$
		CD$ = UCASE$(TrialCD$)
		IF NSlice% >= 0 THEN 'avoid orphan binary files
			NSlice% = -99
			IF LEN(DIR$("ELEMENTS.BIN")) > 0 THEN KILL "ELEMENTS.BIN"
		END IF
	END IF
END SUB

SUB Icosahedron (Ins, NSlice%)
'Generate an icosahedron:
'there are 12 vertices, 30 edges, and 20 triangles.
'Then, subdivide each face NSlice% times, and output triangular elements.
CONST S# = 1.107148719#
DIM Lat#(12), Lon#(12)
DIM V1#(3), V2#(3), V3#(3)
'---------------------------------------------------------------
' generate basic form with a vertex (a 5-fold axis) up; highest symmetry axis
Lat#(1) = 1.570796327#  '90d
Lon#(1) = 0#
FOR i% = 2 TO 6
	Lat#(i%) = Lat#(1) - S#
	Lon#(i%) = (i% - 2#) * 1.256637061#  '72d
NEXT i%
FOR i% = 7 TO 11
	Lat#(i%) = -Lat#(1) + S#
	Lon#(i%) = (i% - 7#) * 1.256637061# + .628318531#  '72d, 36d
NEXT i%
Lat#(12) = -Lat#(1)
Lon#(12) = 0#
DIM ABG#(3, 12)  'Cartesian (alpha, beta, gamma) coordinates of these vertices.
FOR i% = 1 TO 12
	ABG#(1, i%) = COS(Lat#(i%)) * COS(Lon#(i%))
	ABG#(2, i%) = COS(Lat#(i%)) * SIN(Lon#(i%))
	ABG#(3, i%) = SIN(Lat#(i%))
NEXT i%
'-------------------------------------------------------
'find all 20 faces and subdivide each into four spherical triangles;
IF LEN(DIR$("ELEMENTS.BIN")) > 0 THEN KILL "ELEMENTS.BIN"
OPEN "ELEMENTS.BIN" FOR BINARY AS #5 LEN = 3600
FOR i% = 1 TO 10
	FOR J% = (i% + 1) TO 11
		FOR K% = (J% + 1) TO 12
			dot1# = ABG#(1, i%) * ABG#(1, J%) + ABG#(2, i%) * ABG#(2, J%) + ABG#(3, i%) * ABG#(3, J%)
			dot2# = ABG#(1, J%) * ABG#(1, K%) + ABG#(2, J%) * ABG#(2, K%) + ABG#(3, J%) * ABG#(3, K%)
			dot3# = ABG#(1, K%) * ABG#(1, i%) + ABG#(2, K%) * ABG#(2, i%) + ABG#(3, K%) * ABG#(3, i%)
			IF (dot1# > .3#) AND (dot2# > .3#) AND (dot3# > .3#) THEN
				X1# = ABG#(1, i%)
				X2# = ABG#(1, J%)
				X3# = ABG#(1, K%)
				Y1# = ABG#(2, i%)
				Y2# = ABG#(2, J%)
				Y3# = ABG#(2, K%)
				z1# = ABG#(3, i%)
				z2# = ABG#(3, J%)
				z3# = ABG#(3, K%)
				'Note: Divide will call itself, recursively.
				'Therefore, all inputs are simple numbers (not vectors),
				'to go on stack as values, not addresses.
				CALL Divide(X1#, Y1#, z1#, X2#, Y2#, z2#, X3#, Y3#, z3#, NSlice%)
			END IF
		NEXT K%
	NEXT J%
NEXT i%
CLOSE #5
END SUB

SUB KolorMenu (Ins, Blue%, Brown%, Green%, Red%, White%, Yellow%, Mods, BaseC%, CircleC%, DoIcon%, ElementC%, FaultC1%, FaultC2%, NodeC1%, NodeC2%)
DIM Kolors%(0 TO 6)
FALSE% = 0
TRUE% = NOT FALSE%
	Kolors%(0) = 0
	Kolors%(1) = Blue%
	Kolors%(2) = Green%
	Kolors%(3) = Brown%
	Kolors%(4) = Red%
	Kolors%(5) = Yellow%
	Kolors%(6) = White%
	CALL Blanker(16, 24, 2, 79)
	LOCATE 16, 5: PRINT "CODE: 0=black; 1=blue; 2=green; 3=brown; 4=red; 5=yellow; 6=white.";
K1:
	FOR i% = 0 TO 6
		IF CircleC% = Kolors%(i%) THEN ColorI% = i%
	NEXT i%
	LOCATE 18, 57: PRINT USING "#"; ColorI%;
	LOCATE 18, 10: INPUT ; "Color of outer circle representing the sphere: ", a$
	IF a$ = "" THEN GOTO K2
	KolorI% = VAL(a$)  'avoid redo-from'start scrolling!
	IF KolorI% < 0 OR KolorI% > 6 THEN
		BEEP
		CALL Blanker(18, 19, 2, 79)
		GOTO K1
	END IF
	CircleC% = Kolors%(KolorI%)
K2:
	FOR i% = 0 TO 6
		IF BaseC% = Kolors%(i%) THEN ColorI% = i%
	NEXT i%
	LOCATE 19, 32: PRINT USING "#"; ColorI%;
	LOCATE 19, 10: INPUT ; "Color of the basemap: ", a$
	IF a$ = "" THEN GOTO K3
	KolorI% = VAL(a$)  'avoid redo-from'start scrolling!
	IF KolorI% < 0 OR KolorI% > 6 THEN
		BEEP
		CALL Blanker(19, 20, 2, 79)
		GOTO K2
	END IF
	BaseC% = Kolors%(KolorI%)
K3:
	FOR i% = 0 TO 6
		IF NodeC1% = Kolors%(i%) THEN ColorI% = i%
	NEXT i%
	LOCATE 20, 26: PRINT USING "#"; ColorI%;
	LOCATE 20, 10: INPUT ; "Color of nodes: ", a$
	IF a$ = "" THEN GOTO K4
	KolorI% = VAL(a$)  'avoid redo-from'start scrolling!
	IF KolorI% < 0 OR KolorI% > 6 THEN
		BEEP
		CALL Blanker(20, 21, 2, 79)
		GOTO K3
	END IF
	NodeC1% = Kolors%(KolorI%)
	IF NodeC1% <> Red% THEN
		NodeC2% = Red%
	ELSE
		NodeC2% = Yellow%
	END IF
K4:
	FOR i% = 0 TO 6
		IF ElementC% = Kolors%(i%) THEN ColorI% = i%
	NEXT i%
	LOCATE 21, 34: PRINT USING "#"; ColorI%;
	LOCATE 21, 10: INPUT ; "Color of element sides: ", a$
	IF a$ = "" THEN GOTO K5
	KolorI% = VAL(a$)  'avoid redo-from'start scrolling!
	IF KolorI% < 0 OR KolorI% > 6 THEN
		BEEP
		CALL Blanker(21, 22, 2, 79)
		GOTO K4
	END IF
	ElementC% = Kolors%(KolorI%)
K5:
	FOR i% = 0 TO 6
		IF FaultC1% = Kolors%(i%) THEN ColorI% = i%
	NEXT i%
	LOCATE 22, 27: PRINT USING "#"; ColorI%;
	LOCATE 22, 10: INPUT ; "Color of faults: ", a$
	IF a$ = "" THEN GOTO K6
	KolorI% = VAL(a$)  'avoid redo-from'start scrolling!
	IF KolorI% < 0 OR KolorI% > 6 THEN
		BEEP
		CALL Blanker(22, 23, 2, 79)
		GOTO K5
	END IF
	FaultC1% = Kolors%(KolorI%)
	IF FaultC1% <> Yellow% THEN
		FaultC2% = Yellow%
	ELSE
		FaultC2% = Red%
	END IF
K6:
	IF DoIcon% THEN a$ = "yes" ELSE a$ = "no"
	LOCATE 24, 71: PRINT a$;
	LOCATE 24, 6: PRINT "Should the element-center icons (";
	IF Colored% THEN COLOR Green%
	PRINT CHR$(17);
	IF Colored% THEN COLOR White%
	PRINT " or ";
	IF Colored% THEN COLOR Yellow%
	PRINT "X";
	IF Colored% THEN COLOR White%
	INPUT ; ") be displayed (Yes/No)?: ", a$
	IF a$ = "" THEN GOTO K7
	IF LEFT$(a$, 1) = "N" OR LEFT$(a$, 1) = "n" THEN
		DoIcon% = FALSE%
	ELSE
		DoIcon% = TRUE%
	END IF
K7:
	CALL Blanker(16, 24, 2, 79)
END SUB

SUB LineOfNumbers (Ins, Rec$, Outs, AreNumbers%)
'Evaluates one line of a .dig file to decide if it contains a number pair.
	FALSE% = 0
	TRUE% = NOT 0
	a$ = LEFT$(Rec$, 2)
	IF (a$ = " +") OR (a$ = " -") THEN
		AreNumbers% = TRUE%
	ELSE
		IF (a$ = "  ") THEN
			b$ = LEFT$(Rec$, 3)
			b3$ = RIGHT$(b$, 1)
			IF (LEN(b3$) < 1) THEN b3$ = "x"
			b3% = ASC(b3$) 'look for digits 0, 1, ... 9
			IF (b3% > 47) AND (b3% < 58) THEN
				AreNumbers% = TRUE%
			ELSE
				AreNumbers% = FALSE%
			END IF
		ELSE
			AreNumbers% = FALSE%
		END IF
	END IF
END SUB

SUB ModeParms (Ins, BestMode%, Outs, background%, Blue%, Brown%, Foreground%, Green%, HotColor%, Lines%, Red%, White%, Yellow%)
'Table of parameters describing the different screen modes
'and colors (if any) in a standard format.
'Colored%, HiColor%, HiCol%, and HiRow% are passed thru COMMON SHARED.
DIM ColArray&(0 TO 15)
FALSE% = 0
TRUE% = NOT FALSE%
SELECT CASE BestMode%
	   CASE 2
			 Colored% = FALSE%
			 HiColor% = 1
			 HiRow% = 199
			 HiCol% = 639
			 Lines% = 25
			 HotColor% = 1
			 Brown% = 1
			 Green% = 1
			 Red% = 1
			 Yellow% = 1
			 White% = 1
			 Blue% = 1
	   CASE 3
			 Colored% = FALSE%
			 HiColor% = 1
			 HiRow% = 347
			 HiCol% = 719
			 Lines% = 25
			 HotColor% = 1
			 Brown% = 1
			 Green% = 1
			 Red% = 1
			 Yellow% = 1
			 White% = 1
			 Blue% = 1
	   CASE 9
			 Colored% = TRUE%
			 HiColor% = 15
			 ColArray&(15) = 63&
			 ColArray&(14) = 55&
			 ColArray&(13) = 28&
			 ColArray&(12) = 36&
			 ColArray&(11) = 38&
			 ColArray&(10) = 20&
			 ColArray&(9) = 54&
			 ColArray&(8) = 18&
			 ColArray&(7) = 16&
			 ColArray&(6) = 11&
			 ColArray&(5) = 9&
			 ColArray&(4) = 8&
			 ColArray&(3) = 45&
			 ColArray&(2) = 40&
			 ColArray&(1) = 56&
			 ColArray&(0) = 0&
			 PALETTE USING ColArray&(0)
			 Brown% = 10
			 Green% = 8
			 Red% = 12
			 Yellow% = 9
			 White% = 15
			 Blue% = 5
			 HotColor% = Yellow%
			 Foreground% = White%
			 background% = 20 'for this mode, uses color (not attribute)
			 HiRow% = 349
			 HiCol% = 639
			 Lines% = 25
	   CASE 10
			 Colored% = FALSE%
			 HiColor% = 1
			 HiRow% = 349
			 HiCol% = 639
			 Lines% = 25
			 HotColor% = 1
			 Brown% = 1
			 Green% = 1
			 Red% = 1
			 Yellow% = 1
			 White% = 1
			 Blue% = 1
	   CASE 11
			 Colored% = FALSE%
			 HiColor% = 1
			 HiRow% = 479
			 HiCol% = 639
			 Lines% = 30
			 HotColor% = 1
			 Brown% = 1
			 Green% = 1
			 Red% = 1
			 Yellow% = 1
			 White% = 1
			 Blue% = 1
	   CASE 12
			 Colored% = TRUE%
			 HiColor% = 15
			 ColArray&(15) = 4144959
			 ColArray&(14) = 3421236
			 ColArray&(13) = 1976383
			 ColArray&(12) = 63&
			 ColArray&(11) = 207167
			 ColArray&(10) = 332840
			 ColArray&(9) = 16191&
			 ColArray&(8) = 16128&
			 ColArray&(7) = 6400&
			 ColArray&(6) = 3289600
			 ColArray&(5) = 4128768
			 ColArray&(4) = 2031616
			 ColArray&(3) = 4128831
			 ColArray&(2) = 1966110
			 ColArray&(1) = 986895
			 ColArray&(0) = 0&
			 PALETTE USING ColArray&(0)
			 Brown% = 10
			 Green% = 8
			 Red% = 12
			 Yellow% = 9
			 White% = 15
			 Blue% = 5
			 HotColor% = Yellow%
			 Foreground% = White%
			 background% = Brown% 'for this mode, uses attribute (not color)
			 HiRow% = 479
			 HiCol% = 639
			 Lines% = 30
END SELECT
END SUB

SUB Nearest (Ins, NodeABG!(), NUMNOD%, x!, Y!, Outs, NumNear%, NearOnes%())
'identifies up to 4 nodes (in same spot) nearest to (X!,Y!)
DIM TempVec!(3)
CALL XY2ABG(Ins, x!, Y!, Outs, Outside%, TempVec!())
IF Outside% THEN
	NumNear% = 0
	FOR i% = 1 TO 4: NearOnes%(i%) = 0: NEXT i%
ELSE
	R2toN! = 9.99E+29
	FOR i% = 1 TO NUMNOD%
		R2! = (TempVec!(1) - NodeABG!(1, i%)) ^ 2 + (TempVec!(2) - NodeABG!(2, i%)) ^ 2 + (TempVec!(3) - NodeABG!(3, i%)) ^ 2
		IF R2! < R2toN! THEN
			AN! = NodeABG!(1, i%)
			BN! = NodeABG!(2, i%)
			GN! = NodeABG!(3, i%)
			R2toN! = R2!
		END IF
	NEXT i%
	NumNear% = 0
	NearOnes%(1) = 0: NearOnes%(2) = 0: NearOnes%(3) = 0: NearOnes%(4) = 0:
	FOR i% = 1 TO NUMNOD%
		IF NodeABG!(1, i%) = AN! THEN
			IF NodeABG!(2, i%) = BN! THEN
				IF NodeABG!(3, i%) = GN! THEN
					NumNear% = NumNear% + 1
					NearOnes%(NumNear%) = i%
				END IF
			END IF
		END IF
		IF NumNear% = 4 THEN GOTO DidIt:
	NEXT i%
END IF
DidIt:
END SUB

SUB SecondScreen (Ins, BestMode%, Foreground%, HotColor%, Outs, PageNow%)
'choice between (COMMON) FEP% = 1 (Shells), 2 (Restore2/NeoKinema), 3 (Restore3)
TheTop:
SCREEN BestMode%: CLS
PageNow% = 0
IF Colored% THEN COLOR HotColor%
LOCATE 4, 33
PRINT "DECISION POINT:";
LOCATE 5, 33
PRINT "===============";
LOCATE 7, 2
PRINT "Programs Shells, NeoKinema, and Restore use different file formats,";
LOCATE 8, 2
PRINT "because in Shells/NeoKinema/Restore2 there are no per-element data, but in";
LOCATE 9, 2
PRINT "Restore3 there are three real numbers per element.  Also, .feg files for";
LOCATE 10, 2
PRINT "NeoKinema and Restore2/3 are not permitted to have faults.";
LOCATE 15, 2
PRINT "Enter 1 for Shells format (faults, but no per-element data),";
LOCATE 16, 2
PRINT "   or 2 for Restore2/NeoKinema format (no faults or per-element data).";
LOCATE 17, 2
PRINT "   or 3 for Restore3 format (per-element data, but no faults).";
LOCATE 19, 30
INPUT "SELECTION (1, 2, or 3): ", FEP%
IF FEP% < 1 OR FEP% > 3 THEN
	BEEP
	GOTO TheTop
END IF
END SUB

SUB SetBins (Ins, Ask%, Contour%, EQCM!(), IData%, NUMEL%, NUMNOD%, Mus!(), Mods, OldBotF!, OldIData%, OldTopF!, Outs, BotF!, DFC!, Logs%, NeedToDraw%, TopF!)
'Figure out best assignment of (HiColor% - 1) colors to nodal data values
'Logs% indicates that base-10 logarithmic scale should be used.
'If IData% = 1, 2, 3, 4, then data is taken from EQCM!(IData%, nodes).
'If IData% = 5, 6, 7, then data is taken from Mus!(IData%-4, elements).
FALSE% = 0
TRUE% = NOT FALSE%
Q2:  CALL Blanker(23, 24, 2, 79)
	IF (IData% > 4) THEN
		MaxF! = Mus!(IData% - 4, 1)
		MinF! = MaxF!
		FOR i% = 1 TO NUMEL%
			IF (Mus!(IData% - 4, i%) > MaxF!) THEN MaxF! = Mus!(IData% - 4, i%)
			IF (Mus!(IData% - 4, i%) < MinF!) THEN MinF! = Mus!(IData% - 4, i%)
		NEXT i%
	ELSE
		MaxF! = EQCM!(IData%, 1)
		MinF! = MaxF!
		FOR i% = 1 TO NUMNOD%
			IF (EQCM!(IData%, i%) > MaxF!) THEN MaxF! = EQCM!(IData%, i%)
			IF (EQCM!(IData%, i%) < MinF!) THEN MinF! = EQCM!(IData%, i%)
		NEXT i%
	END IF
	IF MaxF! > MinF! THEN
		IF Ask% THEN
			LOCATE 23, 15: INPUT ; "Do you want colors to span the data range? (Yes/no): ", a$
			IF (LEFT$(a$, 1) = "N") OR (LEFT$(a$, 1) = "n") THEN
				GetF1F2% = TRUE%
			ELSE
				GetF1F2% = FALSE%
			END IF
		ELSE
			GetF1F2% = FALSE%
		END IF
	ELSE
		LOCATE 23, 15: PRINT "All values are identical: "; MaxF!
		GetF1F2% = TRUE%
	END IF
	IF GetF1F2% THEN
		LOCATE 24, 15: INPUT ; "Enter lowest colored value: ", a$
		MinF! = VAL(a$)
		CALL Blanker(24, 24, 2, 79)
		LOCATE 24, 15: INPUT ; "Enter highest colored value: ", a$
		MaxF! = VAL(a$)
		IF MaxF! <= MinF! THEN
			BEEP
			GOTO Q2
		END IF
	END IF
	IF (MinF! > 0!) AND (MaxF! <= .001) THEN
		Logs% = TRUE%
		BotF! = LOG(MinF!) / LOG(10#)
		TopF! = LOG(MaxF!) / LOG(10#)
		DFC! = (TopF! - BotF!) / CSNG(HiColor% - 1)
	ELSE
		Logs% = FALSE%
		BotF! = MinF!
		TopF! = MaxF!
	END IF
	DFC! = (TopF! - BotF!) / CSNG(HiColor% - 1)
	IF (NOT Contour%) OR (IData% <> OldIData%) OR (TopF! <> OldTopF!) OR (BotF! <> OldBotF!) THEN
		NeedToDraw% = TRUE%
	END IF
	OldIData% = IData%
	OldBotF! = BotF!
	OldTopF! = TopF!
END SUB

SUB Unitize (Mods, x#, Y#, Z#)
'Converts any vector in 3-component double-precision form to a unit vector.
R# = SQR(x# * x# + Y# * Y# + Z# * Z#)
IF R# > 0 THEN
	x# = x# / R#
	Y# = Y# / R#
	Z# = Z# / R#
ELSE
	x# = 1#
	Y# = 0#
	Z# = 0#
END IF
END SUB

SUB WinFrame (Ins, WinLat!, WinLon!, Outs, WinRight!(), WinOut!(), WinUp!())
'determines unit vectors of window frame
WinOut!(1) = COS(WinLat!) * COS(WinLon!)
WinOut!(2) = COS(WinLat!) * SIN(WinLon!)
WinOut!(3) = SIN(WinLat!)
WinRight!(1) = COS(WinLon! + 1.570796)
WinRight!(2) = SIN(WinLon! + 1.570796)
WinRight!(3) = 0!
CALL Cross(Ins, WinOut!(), WinRight!(), Outs, WinUp!())
END SUB