home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
Antennas
/
Antennas_CD-ROM_Walnut_Creek_September_1996.iso
/
mininec
/
amiga
/
smn9.asc
< prev
next >
Wrap
Text File
|
1996-06-30
|
48KB
|
1,620 lines
5 REM GEOMETRY MODIFIED 17 OCT 86 R.P.HAVILAND
10 REM ****** MININEC(3) ********** NOSC CODE 822 (JCL) 4-86 WITH REVS 1-9
20 DEFINT I-N,MA,MI,ML,MM,MP,FLG,FVS
30 DIM K!(6,2),Q(14)
40 REM ----- MAXIMUM NUMBER OF SEGMENTS (PULSES + 2 * WIRES) = 105
50 MS=105
60 DIM X(105),Y(105),Z(105)
70 REM ----- MAXIMUM NUMBER OF WIRES = 35
80 MW=35
90 DIM A(35),CA(35),CB(35),CG(35),J1(35),J2(35,2),N(35,2),S(35)
100 REM ----- MAXIMUM NUMBER OF LOADS = 8
110 ML=8
120 REM ----- MAXIMUM ORDER OF S-PARAMETER LOADS = 8
130 MA=8
140 DIM LA(2,8,8),LP(8),LS(8)
150 REM ----- MAXIMUM NUMBER OF MEDIA = 6
160 MM=6
170 REM ----- H MUST BE DIMENSIONED AT LEAST 6
180 DIM H(6),T(6),U(6),V(6),Z1(6),Z2(6)
190 REM ----- MAXIMUM NUMBER OF PULSES = 35
200 MP=35
210 DIM C%(35,2),CI(35),CR(35),P(35),W%(35)
220 DIM ZR(35,35),ZI(35,35)
230 REM ---- ARRAYS E,L & M DIMENSIONED TO MW+MP=70
240 DIM E(70),L(70),M(70)
250 COLOR 2,0
260 GOTO 14870
270 REM ********** KERNEL EVALUATION OF INTEGRALS I2 & I3 **********
280 IF K<0 THEN 330
290 X3=X2+T*(V1-X2)
300 Y3=Y2+T*(V2-Y2)
310 Z3=Z2+T*(V3-Z2)
320 GOTO 360
330 X3=V1+T*(X2-V1)
340 Y3=V2+T*(Y2-V2)
350 Z3=V3+T*(Z2-V3)
360 D3=X3*X3+Y3*Y3+Z3*Z3
370 REM ----- MOD FOR SMALL RADIUS TO WAVELENGTH RATIO
380 IF A(P4)<=SRM THEN D=SQR(D3):GOTO 490
390 D=D3+A2
400 IF D>0 THEN D=SQR(D)
410 REM ----- CRITERIA FOR USING REDUCED KERNEL
420 IF I6!=0 THEN 490
430 REM ----- EXACT KERNEL CALCULATION WITH ELLIPTIC INTEGRAL
440 B=D3/(D3+4*A2)
450 W0=C0+B*(C1+B*(C2+B*(C3+B*C4)))
460 W1=C5+B*(C6+B*(C7+B*(C8+B*C9)))
470 V0=(W0-W1*LOG(B))*SQR(1-B)
480 T3=T3+(V0+LOG(D3/(64*A2))/2)/P/A(P4)-1/D
490 B1=D*W
500 REM ----- EXP(-J*K*R)/R
510 T3=T3+COS(B1)/D
520 T4=T4-SIN(B1)/D
530 RETURN
540 REM ***** PSI(P1,P2,P3) = T1 + J * T2 **********
550 REM ----- ENTRIES REQUIRED FOR NEAR FIELD CALCULATION
560 X1=X0+P1*T5/2
570 Y1=Y0+P1*T6/2
580 Z1=Z0+P1*T7/2
590 X2=X1-X(P2)
600 Y2=Y1-Y(P2)
610 Z2=Z1-K*Z(P2)
620 V1=X1-X(P3)
630 V2=Y1-Y(P3)
640 V3=Z1-K*Z(P3)
650 GOTO 1350
660 I4=INT(P2)
670 I5=I4+1
680 X2=X0-(X(I4)+X(I5))/2
690 Y2=Y0-(Y(I4)+Y(I5))/2
700 Z2=Z0-K*(Z(I4)+Z(I5))/2
710 V1=X0-X(P3)
720 V2=Y0-Y(P3)
730 V3=Z0-K*Z(P3)
740 GOTO 1350
750 X2=X0-X(P2)
760 Y2=Y0-Y(P2)
770 Z2=Z0-K*Z(P2)
780 I4=INT(P3)
790 I5=I4+1
800 V1=X0-(X(I4)+X(I5))/2
810 V2=Y0-(Y(I4)+Y(I5))/2
820 V3=Z0-K*(Z(I4)+Z(I5))/2
830 GOTO 1350
840 REM ----- ENTRIES REQUIRED FOR IMPEDANCE MATRIX CALCULATION
850 REM ----- S(M) GOES IN (X1,Y1,Z1) FOR SCALAR POTENTIAL
860 REM ----- MOD FOR SMALL RADIUS TO WAVE LENGTH RATIO
870 FVS=1
880 IF K<1 THEN 940
890 IF A(P4)>SRM THEN 940
900 IF (P3=P2+1 AND P1=(P2+P3)/2) THEN 910 ELSE 940
910 T1=2*LOG(S(P4)/A(P4))
920 T2=-W*S(P4)
930 RETURN
940 I4=INT(P1)
950 I5=I4+1
960 X1=(X(I4)+X(I5))/2
970 Y1=(Y(I4)+Y(I5))/2
980 Z1=(Z(I4)+Z(I5))/2
990 GOTO 1130
1000 REM ----- S(M) GOES IN (X1,Y1,Z1) FOR VECTOR POTENTIAL
1010 REM ----- MOD FOR SMALL RADIUS TO WAVE LENGTH RATIO
1020 FVS=0
1030 IF K<1 THEN 1090
1040 IF A(P4)>=SRM THEN 1090
1050 IF (I=J AND P3=P2+.5) THEN 1060 ELSE 1090
1060 T1=LOG(S(P4)/A(P4))
1070 T2=-W*S(P4)/2
1080 RETURN
1090 X1=X(P1)
1100 Y1=Y(P1)
1110 Z1=Z(P1)
1120 REM ----- S(U)-S(M) GOES IN (X2,Y2,Z2)
1130 I4=INT(P2)
1140 IF I4=P2 THEN 1200
1150 I5=I4+1
1160 X2=(X(I4)+X(I5))/2-X1
1170 Y2=(Y(I4)+Y(I5))/2-Y1
1180 Z2=K*(Z(I4)+Z(I5))/2-Z1
1190 GOTO 1240
1200 X2=X(P2)-X1
1210 Y2=Y(P2)-Y1
1220 Z2=K*Z(P2)-Z1
1230 REM ----- S(V)-S(M) GOES IN (V1,V2,V3)
1240 I4=INT(P3)
1250 IF I4=P3 THEN 1310
1260 I5=I4+1
1270 V1=(X(I4)+X(I5))/2-X1
1280 V2=(Y(I4)+Y(I5))/2-Y1
1290 V3=K*(Z(I4)+Z(I5))/2-Z1
1300 GOTO 1350
1310 V1=X(P3)-X1
1320 V2=Y(P3)-Y1
1330 V3=K*Z(P3)-Z1
1340 REM ----- MAGNITUDE OF S(U) - S(M)
1350 D0=X2*X2+Y2*Y2+Z2*Z2
1360 REM ----- MAGNITUDE OF S(V) - S(M)
1370 IF D0>0 THEN D0=SQR(D0)
1380 D3=V1*V1+V2*V2+V3*V3
1390 IF D3>0 THEN D3=SQR(D3)
1400 REM ----- SQUARE OF WIRE RADIUS
1410 A2=A(P4)*A(P4)
1420 REM ----- MAGNITUDE OF S(V) - S(U)
1430 S4=(P3-P2)*S(P4)
1440 REM ----- ORDER OF INTEGRATION
1450 REM ----- LTH ORDER GAUSSIAN QUADRATURE
1460 T1=0
1470 T2=0
1480 I6!=0
1490 F2=1
1500 L=7
1510 T=(D0+D3)/S(P4)
1520 REM ----- CRITERIA FOR EXACT KERNEL
1530 IF T>1.1 THEN 1650
1540 IF C$="N" THEN 1650
1550 IF J2(W%(I),1)=J2(W%(J),1) THEN 1600
1560 IF J2(W%(I),1)=J2(W%(J),2) THEN 1600
1570 IF J2(W%(I),2)=J2(W%(J),1) THEN 1600
1580 IF J2(W%(I),2)=J2(W%(J),2) THEN 1600
1590 GOTO 1650
1600 IF A(P4)>SRM THEN 1620
1610 IF FVS=1 THEN 910 ELSE 1060
1620 F2=2*(P3-P2)
1630 I6!=(1-LOG(S4/F2/8/A(P4)))/P/A(P4)
1640 GOTO 1670
1650 IF T>6 THEN L=3
1660 IF T>10 THEN L=1
1670 I5=L+L
1680 T3=0
1690 T4=0
1700 T=(Q(L)+.5)/F2
1710 GOSUB 280
1720 T=(.5-Q(L))/F2
1730 GOSUB 280
1740 L=L+1
1750 T1=T1+Q(L)*T3
1760 T2=T2+Q(L)*T4
1770 L=L+1
1780 IF L<I5 THEN 1680
1790 T1=S4*(T1+I6!)
1800 T2=S4*T2
1810 RETURN
1820 REM ********** COMPLEX SQUARE ROOT **********
1830 REM ----- W6+I*W7=SQR(Z6+I*Z7)
1840 T6=SQR((ABS(Z6)+SQR(Z6*Z6+Z7*Z7))/2)
1850 T7=ABS(Z7)/2/T6
1860 IF Z6<0 THEN 1910
1870 W6=T6
1880 W7=T7
1890 IF Z7<0 THEN W7=-T7
1900 RETURN
1910 W6=T7
1920 W7=T6
1930 IF Z7<0 THEN W7=-T6
1940 RETURN
1950 REM ********** IMPEDANCE MATRIX CALCULATION **********
1960 IF FLG=1 THEN 4270
1970 IF FLG=2 THEN 4760
1980 REM ----- BEGIN MATRIX FILL TIME CALCULATION
1990 OT$=TIME$
2000 Q$="MATRIX FILL "
2010 CLS
2020 PRINT "BEGIN ";Q$
2030 REM ----- ZERO IMPEDANCE MATRIX
2040 FOR I=1 TO N
2050 FOR J=1 TO N
2060 ZR(I,J)=0
2070 ZI(I,J)=0
2080 NEXT J
2090 NEXT I
2100 REM ----- COMPUTE ROW I OF MATRIX (OBSERVATION LOOP)
2110 FOR I=1 TO N
2120 I1=ABS(C%(I,1))
2130 I2=ABS(C%(I,2))
2140 F4=SGN(C%(I,1))*S(I1)
2150 F5=SGN(C%(I,2))*S(I2)
2160 REM ----- R(M + 1/2) - R(M - 1/2) HAS COMPONENTS (T5,T6,T7)
2170 T5=F4*CA(I1)+F5*CA(I2)
2180 T6=F4*CB(I1)+F5*CB(I2)
2190 T7=F4*CG(I1)+F5*CG(I2)
2200 IF C%(I,1)=-C%(I,2) THEN T7=S(I1)*(CG(I1)+CG(I2))
2210 REM ----- COMPUTE COLUMN J OF ROW I (SOURCE LOOP)
2220 FOR J=1 TO N
2230 J1=ABS(C%(J,1))
2240 J2=ABS(C%(J,2))
2250 F4=SGN(C%(J,1))
2260 F5=SGN(C%(J,2))
2270 F6=1
2280 F7=1
2290 REM ----- IMAGE LOOP
2300 FOR K=1 TO G STEP -2
2310 IF C%(J,1)<>-C%(J,2) THEN 2350
2320 IF K<0 THEN 3320
2330 F6=F4
2340 F7=F5
2350 F8=0
2360 IF K<0 THEN 2480
2370 REM ----- SET FLAG TO AVOID REDUNANT CALCULATIONS
2380 IF I1<>I2 THEN 2460
2390 IF (CA(I1)+CB(I1))=0 THEN 2410
2400 IF C%(I,1)<>C%(I,2) THEN 2460
2410 IF J1<>J2 THEN 2460
2420 IF (CA(J1)+CB(J1))=0 THEN 2440
2430 IF C%(J,1)<>C%(J,2) THEN 2460
2440 IF I1=J1 THEN F8=1
2450 IF I=J THEN F8=2
2460 IF ZR(I,J)<>0 THEN 3170
2470 REM ----- COMPUTE PSI(M,N,N+1/2)
2480 P1=2*W%(I)+I-1
2490 P2=2*W%(J)+J-1
2500 P3=P2+.5
2510 P4=J2
2520 GOSUB 1020
2530 U1=F5*T1
2540 U2=F5*T2
2550 REM ----- COMPUTE PSI(M,N-1/2,N)
2560 P3=P2
2570 P2=P2-.5
2580 P4=J1
2590 IF F8<2 THEN GOSUB 1020
2600 V1=F4*T1
2610 V2=F4*T2
2620 REM ----- S(N+1/2)*PSI(M,N,N+1/2) + S(N-1/2)*PSI(M,N-1/2,N)
2630 X3=U1*CA(J2)+V1*CA(J1)
2640 Y3=U1*CB(J2)+V1*CB(J1)
2650 Z3=(F7*U1*CG(J2)+F6*V1*CG(J1))*K
2660 REM ----- REAL PART OF VECTOR POTENTIAL CONTRIBUTION
2670 D1=W2*(X3*T5+Y3*T6+Z3*T7)
2680 X3=U2*CA(J2)+V2*CA(J1)
2690 Y3=U2*CB(J2)+V2*CB(J1)
2700 Z3=(F7*U2*CG(J2)+F6*V2*CG(J1))*K
2710 REM ----- IMAGINARY PART OF VECTOR POTENTIAL CONTRIBUTION
2720 D2=W2*(X3*T5+Y3*T6+Z3*T7)
2730 REM ----- COMPUTE PSI(M+1/2,N,N+1)
2740 P1=P1+.5
2750 IF F8=2 THEN P1=P1-1
2760 P2=P3
2770 P3=P3+1
2780 P4=J2
2790 IF F8<>1 THEN 2830
2800 U5=F5*U1+T1
2810 U6=F5*U2+T2
2820 GOTO 2910
2830 GOSUB 870
2840 IF F8<2 THEN 2880
2850 U1=(2*T1-4*U1*F5)/S(J1)
2860 U2=(2*T2-4*U2*F5)/S(J1)
2870 GOTO 3140
2880 U5=T1
2890 U6=T2
2900 REM ----- COMPUTE PSI(M-1/2,N,N+1)
2910 P1=P1-1
2920 GOSUB 870
2930 U1=(T1-U5)/S(J2)
2940 U2=(T2-U6)/S(J2)
2950 REM ----- COMPUTE PSI(M+1/2,N-1,N)
2960 P1=P1+1
2970 P3=P2
2980 P2=P2-1
2990 P4=J1
3000 GOSUB 870
3010 U3=T1
3020 U4=T2
3030 REM ----- COMPUTE PSI(M-1/2,N-1,N)
3040 IF F8<1 THEN 3080
3050 T1=U5
3060 T2=U6
3070 GOTO 3110
3080 P1=P1-1
3090 GOSUB 870
3100 REM ----- GRADIENT OF SCALAR POTENTIAL CONTRIBUTION
3110 U1=U1+(U3-T1)/S(J1)
3120 U2=U2+(U4-T2)/S(J1)
3130 REM ----- SUM INTO IMPEDANCE MATRIX
3140 ZR(I,J)=ZR(I,J)+K*(D1+U1)
3150 ZI(I,J)=ZI(I,J)+K*(D2+U2)
3160 REM ----- AVOID REDUNANT CALCULATIONS
3170 IF J<I THEN 3320
3180 IF F8=0 THEN 3320
3190 ZR(J,I)=ZR(I,J)
3200 ZI(J,I)=ZI(I,J)
3210 REM ----- SEGMENTS ON SAME WIRE SAME DISTANCE APART HAVE SAME Z
3220 P1=J+1
3230 IF P1>N THEN 3320
3240 IF C%(P1,1)<>C%(P1,2) THEN 3320
3250 IF C%(P1,2)=C%(J,2) THEN 3280
3260 IF C%(P1,2)<>-C%(J,2) THEN 3320
3270 IF (CA(J2)+CB(J2))<>0 THEN 3320
3280 P2=I+13290 IF P2>N THEN 3320
3300 ZR(P2,P1)=ZR(I,J)
3310 ZI(P2,P1)=ZI(I,J)
3320 NEXT K
3330 NEXT J
3340 PCT=I/N
3350 GOSUB 15890
3360 NEXT I
3370 REM ----- END MATRIX FILL TIME CALCULATION
3380 T$=TIME$
3390 GOSUB 15790
3400 PRINT #3," "
3410 PRINT #3,"FILL MATRIX : ";T$
3420 REM ********** ADDITION OF LOADS **********
3430 IF NL=0 THEN 3760
3440 F5=2*P*F
3450 FOR I=1 TO NL
3460 IF L$="N" THEN 3650
3470 REM ----- S-PARAMETER LOADS
3480 U1=0
3490 U2=0
3500 D1=0
3510 D2=0
3520 S=1
3530 FOR J=0 TO LS(I) STEP 2
3540 U1=U1+LA(1,I,J)*S*F5^J
3550 D1=D1+LA(2,I,J)*S*F5^J
3560 L=J+1
3570 U2=U2+LA(1,I,L)*S*F5^L
3580 D2=D2+LA(2,I,L)*S*F5^L
3585 S=-S
3590 NEXT J
3600 J=LP(I)
3610 D=D1*D1+D2*D2
3620 LI=(U2*D1-D2*U1)/D
3630 LR=(U1*D1+U2*D2)/D
3640 GOTO 3680
3650 LR=LA(1,I,1)
3660 LI=LA(2,I,1)
3670 J=LP(I)
3680 F2=1/M
3690 IF C%(J,1)<>-C%(J,2) THEN 3710
3700 IF K<0 THEN F2=2/M
3710 ZR(J,J)=ZR(J,J)+F2*LI
3720 ZI(J,J)=ZI(J,J)-F2*LR
3730 NEXT I
3740 REM ********** IMPEDANCE MATRIX FACTORIZATION **********
3750 REM ----- BEGIN MATRIX FACTOR TIME CALCULATION
3760 OT$=TIME$
3770 Q$="FACTOR MATRIX"
3780 CLS
3790 PRINT "BEGIN ";Q$;
3800 X=N
3810 PCTN=X*(X-1)*(X+X-1)
3820 FOR K=1 TO N-1
3830 REM ----- SEARCH FOR PIVOT
3840 T=ZR(K,K)*ZR(K,K)+ZI(K,K)*ZI(K,K)
3850 I1=K
3860 FOR I=K+1 TO N
3870 T1=ZR(I,K)*ZR(I,K)+ZI(I,K)*ZI(I,K)
3880 IF T1<T THEN 3910
3890 I1=I
3900 T=T1
3910 NEXT I
3920 REM ----- EXCHANGE ROWS K AND I1
3930 IF I1=K THEN 4020
3940 FOR J=1 TO N
3950 T1=ZR(K,J)
3960 T2=ZI(K,J)
3970 ZR(K,J)=ZR(I1,J)
3980 ZI(K,J)=ZI(I1,J)
3990 ZR(I1,J)=T1
4000 ZI(I1,J)=T2
4010 NEXT J
4020 P(K)=I1
4030 REM ----- SUBTRACT ROW K FROM ROWS K+1 TO N
4040 FOR I=K+1 TO N
4050 REM ----- COMPUTE MULTIPLIER L(I,K)
4060 T1=(ZR(I,K)*ZR(K,K)+ZI(I,K)*ZI(K,K))/T
4070 T2=(ZI(I,K)*ZR(K,K)-ZR(I,K)*ZI(K,K))/T
4080 ZR(I,K)=T1
4090 ZI(I,K)=T2
4100 REM ----- SUBTRACT ROW K FROM ROW I
4110 FOR J=K+1 TO N
4120 ZR(I,J)=ZR(I,J)-(ZR(K,J)*T1-ZI(K,J)*T2)
4130 ZI(I,J)=ZI(I,J)-(ZR(K,J)*T2+ZI(K,J)*T1)
4140 NEXT J
4150 NEXT I
4160 X=N-K
4170 PCT=1-X*(X-1)*(X+X-1)/PCTN
4180 GOSUB 15890
4190 NEXT K
4200 REM ----- END MATRIX FACTOR TIME CALCULATION
4210 T$=TIME$
4220 GOSUB 15790
4230 PRINT
4240 PRINT #3, "FACTOR MATRIX: ";T$
4250 REM ********** SOLVE **********
4260 REM ----- COMPUTE RIGHT HAND SIDE
4270 FOR I=1 TO N
4280 CR(I)=0
4290 CI(I)=0
4300 NEXT I
4310 FOR J=1 TO NS
4320 F2=1/M
4330 IF C%(E(J),1)=-C%(E(J),2) THEN F2=2/M
4340 CR(E(J))=F2*M(J)
4350 CI(E(J))=-F2*L(J)
4360 NEXT J
4370 REM ----- PERMUTE EXCITATION
4380 FOR K=1 TO N-1
4390 I1=P(K)
4400 IF I1=K THEN 4470
4410 T1=CR(K)
4420 T2=CI(K)
4430 CR(K)=CR(I1)
4440 CI(K)=CI(I1)
4450 CR(I1)=T1
4460 CI(I1)=T2
4470 NEXT K
4480 REM ----- FORWARD ELIMINATION
4490 FOR I=2 TO N
4500 T1=0
4510 T2=0
4520 FOR J=1 TO I-1
4530 T1=T1+ZR(I,J)*CR(J)-ZI(I,J)*CI(J)
4540 T2=T2+ZR(I,J)*CI(J)+ZI(I,J)*CR(J)
4550 NEXT J
4560 CR(I)=CR(I)-T1
4570 CI(I)=CI(I)-T2
4580 NEXT I
4590 REM ----- BACK SUBSTITUTION
4600 FOR I=N TO 1 STEP -1
4610 T1=0
4620 T2=0
4630 IF I=N THEN 4680
4640 FOR J=I+1 TO N
4650 T1=T1+ZR(I,J)*CR(J)-ZI(I,J)*CI(J)
4660 T2=T2+ZR(I,J)*CI(J)+ZI(I,J)*CR(J)
4670 NEXT J
4680 T=ZR(I,I)*ZR(I,I)+ZI(I,I)*ZI(I,I)
4690 T1=CR(I)-T1
4700 T2=CI(I)-T2
4710 CR(I)=(T1*ZR(I,I)+T2*ZI(I,I))/T
4720 CI(I)=(T2*ZR(I,I)-T1*ZI(I,I))/T
4730 NEXT I
4740 FLG=2
4750 REM ********** SOURCE DATA **********
4760 PRINT #3," "
4770 PRINT #3,B$;" SOURCE DATA ";B$
4780 PWR=0
4790 FOR I=1 TO NS
4800 CR=CR(E(I))
4810 CI=CI(E(I))
4820 T=CR*CR+CI*CI
4830 T1=(L(I)*CR+M(I)*CI)/T
4840 T2=(M(I)*CR-L(I)*CI)/T
4850 O2=(L(I)*CR+M(I)*CI)/2
4860 PWR=PWR+O2
4870 PRINT #3,"PULSE ";E(I),"VOLTAGE = (";L(I);",";M(I);"J)"
4880 PRINT #3," ","CURRENT = (";CR;",";CI;"J)"
4890 PRINT #3," ","IMPEDANCE = (";T1;",";T2;"J)"
4900 PRINT #3," ","POWER = ";O2;" WATTS"
4910 NEXT I
4920 IF NS>1 THEN PRINT #3," "
4930 IF NS>1 THEN PRINT #3,"TOTAL POWER = ";PWR;"WATTS"
4940 RETURN
4950 REM ********** PRINT CURRENTS **********
4960 GOSUB 1960
4970 S$="N"
4980 PRINT #3, " "
4990 PRINT #3,B$;" CURRENT DATA ";B$
5000 FOR K=1 TO NW
5010 IF S$="Y" THEN 5060
5020 PRINT #3, " "
5030 PRINT #3, "WIRE NO. ";K;":"
5040 PRINT #3, "PULSE","REAL","IMAGINARY","MAGNITUDE","PHASE"
5050 PRINT #3, " NO.","(AMPS)","(AMPS)","(AMPS)","(DEGREES)"
5060 N1=N(K,1)
5070 N2=N(K,2)
5080 I=N1
5090 C=C%(I,1)
5100 IF (N1=0 AND N2=0) THEN C=K
5110 IF G=1 THEN 5140
5120 IF (J1(K)=-1 AND N1>N2) THEN N2=N1
5130 IF J1(K)=-1 THEN 5240
5140 E%=1
5150 GOSUB 5710
5160 I2!=I1!
5170 J2!=J1!
5180 GOSUB 6060
5190 IF S$="N" THEN PRINT #3, I$,I1!;TAB(29);J1!;TAB(43);S1;TAB(57);S2
5200 IF S$="Y" THEN PRINT #1,I1!;",";J1!;",";S1;",";S2
5210 IF N1=0 THEN 5310
5220 IF C=K THEN 5240
5230 IF I$="J" THEN N1=N1+1
5240 FOR I=N1 TO N2-1
5250 I2!=CR(I)
5260 J2!=CI(I)
5270 GOSUB 6060
5280 IF S$="N" THEN PRINT #3, I,CR(I);TAB(29);CI(I);TAB(43);S1;TAB(57);S2
5290 IF S$="Y" THEN PRINT #1,CR(I);",";CI(I);",";S1;",";S2
5300 NEXT I
5310 I=N2
5320 C=C%(I,2)
5330 IF (N1=0 AND N2=0) THEN C=K
5340 IF G=1 THEN 5360
5350 IF J1(K)=1 THEN 5420
5360 E%=2
5370 GOSUB 5710
5380 IF (N1=0 AND N2=0) THEN 5480
5390 IF N1>N2 THEN 5480
5400 IF C=K THEN 5420
5410 IF I$="J" THEN 5480
5420 I2!=CR(N2)
5430 J2!=CI(N2)
5440 GOSUB 6060
5450 IF S$="N" THEN PRINT #3, N2,CR(N2);TAB(29);CI(N2);TAB(43);S1;TAB(57);S2
5460 IF S$="Y" THEN PRINT #1,CR(N2);",";CI(N2);",";S1;",";S2
5470 IF J1(K)=1 THEN 5530
5480 I2!=I1!
5490 J2!=J1!
5500 GOSUB 6060
5510 IF S$="N" THEN PRINT #3,I$,I1!;TAB(29);J1!;TAB(43);S1;TAB(57);S2
5520 IF S$="Y" THEN PRINT #1,I1!;",";J1!;",";S1;",";S2
5530 IF S$="Y" THEN PRINT #1," 1 , 1 , 1 , 1"
5540 NEXT K
5550 IF S$="Y" THEN 5680
5560 PRINT
5570 INPUT "SAVE CURRENTS TO A FILE (Y/N) ";S$
5580 IF S$="N" THEN 5690
5590 IF S$<>"Y" THEN 5560
5600 PRINT #3," "
5610 INPUT "FILE PATH +FILENAME (.OUT IS ADDED) ";F$
5620 IF LEFT$(RIGHT$(F$,4),1)="." THEN 5630 ELSE F$=F$+".OUT"
5630 IF O$>"C" THEN PRINT #3,"FILENAME (NAME.OUT): ";F$
5640 OPEN F$ FOR OUTPUT AS #1
5650 PRINT #3," "
5660 PRINT #1,NW;",";PWR;",C"
5670 GOTO 5000
5680 CLOSE #1
5690 RETURN
5700 ----- SORT JUNCTION CURRENTS
5710 I$="E"
5720 I1!=0!
5730 J1!=0!
5740 IF (C=K OR C=0) THEN 5790
5750 I$="J"
5760 I1!=CR(I)
5770 J1!=CI(I)
5780 REM ----- CHECK FOR OTHER OVERLAPPING WIRES
5790 FOR J=1 TO NW
5800 IF J=K GOTO 6030
5810 L1=N(J,1)
5820 L2=N(J,2)
5830 IF E%=2 THEN 5890
5840 CO=C%(L1,1)
5850 CT=C%(L2,2)
5860 L3=L1
5870 L4=L2
5880 GOTO 5930
5890 CO=C%(L2,2)
5900 CT=C%(L1,1)
5910 L3=L2
5920 L4=L1
5930 IF CO=-K THEN 5950
5940 GOTO 5980
5950 I1!=I1!-CR(L3)
5960 J1!=J1!-CI(L3)
5970 I$="J"
5980 IF CT=K THEN 6000
5990 GOTO 6030
6000 I1!=I1!+CR(L4)
6010 J1!=J1!+CI(L4)
6020 I$="J"
6030 NEXT J
6040 RETURN
6050 REM ----- CALCULATE S1 AND S2
6060 I3!=I2!*I2!
6070 J3!=J2!*J2!
6080 IF (I3!>0 OR J3!>0) THEN 6110
6090 S1=0!
6100 GOTO 6120
6110 S1=SQR(I3!+J3!)
6120 IF I2!><0 THEN 6150
6130 S2=0!
6140 RETURN
6150 S2=ATN(J2!/I2!)/P0
6160 IF I2!>0 THEN RETURN
6170 S2=S2+SGN(J2!)*180
6180 RETURN
6190 REM ********** FAR FIELD CALCULATION **********
6200 IF FLG<2 THEN GOSUB 1960
6210 O2=PWR
6220 REM ----- TABULATE IMPEDANCE
6230 IF NM=0 THEN 6330
6240 FOR I=1 TO NM
6250 Z6=T(I)
6260 Z7=-V(I)/(2*P*F*8.85E-06)
6270 REM ----- FORM IMPEDANCE=1/SQR(DIELECTRIC CONSTANT)
6280 GOSUB 1840
6290 D=W6*W6+W7*W7
6300 Z1(I)=W6/D
6310 Z2(I)=-W7/D
6320 NEXT I
6330 PRINT #3," "
6340 PRINT #3,B$;" FAR FIELD ";B$
6350 PRINT #3," "
6360 REM ----- INPUT VARIABLES FOR FAR FIELD CALCULATION
6370 INPUT "CALCULATE PATTERN IN DBI OR VOLTS/METER (D/V)";P$
6380 IF P$="D" THEN 6540
6390 IF P$<>"V" THEN 6370
6400 F1=1
6410 PRINT
6420 PRINT "PRESENT POWER LEVEL = ";PWR;" WATTS"
6430 INPUT "CHANGE POWER LEVEL (Y/N) ";A$
6440 IF A$="N" THEN 6490
6450 IF A$<>"Y" THEN 6430
6460 INPUT "NEW POWER LEVEL (WATTS) ";O2
6470 IF O$>"C" THEN PRINT #3,"NEW POWER LEVEL = ";O2
6480 GOTO 6430
6490 IF (O2<0 OR O2=0) THEN O2=PWR
6500 F1=SQR(O2/PWR)
6510 PRINT
6520 INPUT "RADIAL DISTANCE (METERS) ";RD
6530 IF RD<0 THEN RD=0
6540 A$="ZENITH ANGLE : INITIAL,INCREMENT,NUMBER"
6550 PRINT A$;
6560 INPUT ZA,ZC,NZ
6570 IF NZ=0 THEN NZ=1
6580 IF O$>"C" THEN PRINT #3,A$;": ";ZA;",";ZC;",";NZ
6590 A$="AZIMUTH ANGLE: INITIAL,INCREMENT,NUMBER"
6600 PRINT A$;
6610 INPUT AA,AC,NA
6620 IF NA=0 THEN NA=1
6630 IF O$>"C" THEN PRINT #3,A$;": ";AA;",";AC;",";NA
6640 PRINT #3," "
6650 REM ********** FILE FAR FIELD DATA **********
6660 INPUT "FILE PATTERN (Y/N)";S$
6670 IF S$="N" THEN 6750
6680 IF S$<>"Y" THEN 6660
6690 PRINT #3," "
6700 INPUT "FILENAME (NAME.OUT)";F$
6710 IF LEFT$(RIGHT$(F$,4),1)="." THEN 6720 ELSE F$=F$+".OUT"
6720 IF O$>"C" THEN PRINT #3,"FILENAME (NAME.OUT): ";F$
6730 OPEN F$ FOR OUTPUT AS #1
6740 PRINT #1,NA*NZ;",";O2;",";P$
6750 PRINT #3, " "
6760 K9!=.016678/PWR
6770 REM ----- PATTERN HEADER
6780 PRINT #3,B$;" PATTERN DATA ";B$
6790 IF P$="V" GOTO 684
6800 PRINT #3,"ZENITH","AZIMUTH","VERTICAL","HORIZONTAL","TOTAL"
6810 A$="PATTERN (DB)"
6820 PRINT #3," ANGLE"," ANGLE",A$,A$,A$
6830 GOTO 6910
6840 IF RD>0 THEN PRINT #3,TAB(15);"RADIAL DISTANCE = ";RD;" METERS"
6850 PRINT #3,TAB(15);"POWER LEVEL = ";PWR*F1*F1;" WATTS"
6860 PRINT #3,"ZENITH AZIMUTH"," E(THETA) "," E(PHI)"
6870 A$=" MAG(V/M) PHASE(DEG)"
6880 PRINT #3," ANGLE ANGLE",A$,A$
6890 IF S$="Y" THEN PRINT #1,RD
6900 REM ----- LOOP OVER AZIMUTH ANGLE
6910 Q1=AA
6920 FOR I1=1 TO NA
6930 U3=Q1*P0
6940 V1=-SIN(U3)
6950 V2=COS(U3)
6960 REM ----- LOOP OVER ZENITH ANGLE
6970 Q2=ZA
6980 FOR I2=1 TO NZ
6990 U4=Q2*P0
7000 R3=COS(U4)
7010 T3=-SIN(U4)
7020 T1=R3*V2
7030 T2=-R3*V1
7040 R1=-T3*V2
7050 R2=T3*V1
7060 X1=0
7070 Y1=0
7080 Z1=0
7090 X2=0
7100 Y2=0
7110 Z2=0
7120 REM ----- IMAGE LOOP
7130 FOR K=1 TO G STEP -2
7140 FOR I=1 TO N
7150 IF K>0 THEN 7170
7160 IF C%(I,1)=-C%(I,2) THEN 8110
7170 J=2*W%(I)-1+I
7180 REM ----- FOR EACH END OF PULSE COMPUTE A CONTRIBUTION TO E-FIELD
7190 FOR F5=1 TO 2
7200 L=ABS(C%(I,F5))
7210 F3=SGN(C%(I,F5))*W*S(L)/2
7220 IF C%(I,1)<>-C%(I,2) THEN 7240
7230 IF F3<0 THEN 8100
7240 IF K=1 THEN 7270
7250 IF NM<>0 THEN 7460
7260 REM ----- STANDARD CASE
7270 S2=W*(X(J)*R1+Y(J)*R2+Z(J)*K*R3)
7280 S1=COS(S2)
7290 S2=SIN(S2)
7300 B1=F3*(S1*CR(I)-S2*CI(I))
7310 B2=F3*(S1*CI(I)+S2*CR(I))
7320 IF C%(I,1)=-C%(I,2) THEN 7410
7330 X1=X1+K*B1*CA(L)
7340 X2=X2+K*B2*CA(L)
7350 Y1=Y1+K*B1*CB(L)
7360 Y2=Y2+K*B2*CB(L)
7370 Z1=Z1+B1*CG(L)
7380 Z2=Z2+B2*CG(L)
7390 GOTO 8100
7400 REM ----- GROUNDED ENDS
7410 Z1=Z1+2*B1*CG(L)
7420 Z2=Z2+2*B2*CG(L)
7430 GOTO 8100
7440 REM ----- REAL GROUND CASE
7450 REM ----- BEGIN BY FINDING SPECULAR DISTANCE
7460 T4=100000!
7470 IF R3=0 THEN 7490
7480 T4=-Z(J)*T3/R3
7490 B9=T4*V2+X(J)
7500 IF TB=1 THEN 7530
7510 B9=B9*B9+(Y(J)-T4*V1)^2
7515 IF B9>0 THEN B9=SQR(B9) ELSE 7530
7520 REM ----- SEARCH FOR THE CORRESPONDING MEDIUM
7530 J2=NM
7540 FOR J1=NM TO 1 STEP -1
7550 IF B9>U(J1) THEN 757
7560 J2=J1
7570 NEXT J1
7580 REM ----- OBTAIN IMPEDANCE AT SPECULAR POINT
7590 Z4=Z1(J2)
7600 Z5=Z2(J2)
7610 REM ----- IF PRESENT INCLUDE GROUND SCREEN IMPEDANCE IN PARALLEL
7620 IF NR=0 THEN 7740
7630 IF B9>U(1) THEN 7740
7640 R=B9+NR*RR
7650 Z8=W*R*LOG(R/(NR*RR))/NR
7660 S8=-Z5*Z8
7670 S9=Z4*Z8
7680 T8=Z4
7690 T9=Z5+Z8
7700 D=T8*T8+T9*T9
7710 Z4=(S8*T8+S9*T9)/D
7720 Z5=(S9*T8-S8*T9)/D
7730 REM ----- FORM SQR(1-Z^2*SIN^2)
7740 Z6=1-(Z4*Z4-Z5*Z5)*T3*T3
7750 Z7=-(2*Z4*Z5)*T3*T3
7760 GOSUB 1840
7770 REM ----- VERTICAL REFLECTION COEFFICIENT
7780 S8=R3-(W6*Z4-W7*Z5)
7790 S9=-(W6*Z5+W7*Z4)
7800 T8=R3+(W6*Z4-W7*Z5)
7810 T9=W6*Z5+W7*Z4
7820 D=T8*T8+T9*T9
7830 V8=(S8*T8+S9*T9)/D
7840 V9=(S9*T8-S8*T9)/D
7850 REM ----- HORIZONTAL REFLECTION COEFFICIENT
7860 S8=W6-R3*Z4
7870 S9=W7-R3*Z5
7880 T8=W6+R3*Z4
7890 T9=W7+R3*Z5
7900 D=T8*T8+T9*T9
7910 H8=(S8*T8+S9*T9)/D-V8
7920 H9=(S9*T8-S8*T9)/D-V9
7930 REM ----- COMPUTE CONTRIBUTION TO SUM
7940 S2=W*(X(J)*R1+Y(J)*R2-(Z(J)-2*H(J2))*R3)
7950 S1=COS(S2)
7960 S2=SIN(S2)
7970 B1=F3*(S1*CR(I)-S2*CI(I))
7980 B2=F3*(S1*CI(I)+S2*CR(I))
7990 W6=B1*V8-B2*V9
8000 W7=B1*V9+B2*V8
8010 D=CA(L)*V1+CB(L)*V2
8020 Z6=D*(B1*H8-B2*H9)
8030 Z7=D*(B1*H9+B2*H8)
8040 X1=X1-(CA(L)*W6+V1*Z6)
8050 X2=X2-(CA(L)*W7+V1*Z7)
8060 Y1=Y1-(CB(L)*W6+V2*Z6)
8070 Y2=Y2-(CB(L)*W7+V2*Z7)
8080 Z1=Z1+CG(L)*W6
8090 Z2=Z2+CG(L)*W7
8100 NEXT F5
8110 NEXT I
8120 NEXT K
8130 H2=-(X1*T1+Y1*T2+Z1*T3)*G0
8140 H1=(X2*T1+Y2*T2+Z2*T3)*G0
8150 X4=-(X1*V1+Y1*V2)*G0
8160 X3=(X2*V1+Y2*V2)*G0
8170 IF P$="D" THEN 8240
8180 IF RD=0 THEN 8390
8190 H1=H1/RD
8191 H2=H2/RD
8200 X3=X3/RD
8210 X4=X4/RD
8220 GOTO 8390
8230 REM ----- PATTERN IN DB
8240 P1=-999
8250 P2=P1
8260 P3=P1
8270 T1=K9!*(H1*H1+H2*H2)
8280 T2=K9!*(X3*X3+X4*X4)
8290 T3=T1+T2
8300 REM ----- CALCULATE VALUES IN DB
8310 IF T1>1E-30 THEN P1=4.343*LOG(T1)
8320 IF T2>1E-30 THEN P2=4.343*LOG(T2)
8330 IF T3>1E-30 THEN P3=4.343*LOG(T3)
8340 PRINT #3,Q2;TAB(15);Q1;TAB(29);P1;TAB(43);P2;TAB(57);P3
8350 IF S$="Y" THEN PRINT #1,Q2;",";Q1;",";P1;",";P2;",";P3
8360 GOTO 8630
8370 REM ----- PATTERN IN VOLTS/METER
8380 REM ----- MAGNITUDE AND PHASE OF E(THETA)
8390 S1=0
8400 IF (H1=0 AND H2=0) THEN 8420
8410 S1=SQR(H1*H1+H2*H2)
8420 IF H1><0 THEN 8450
8430 S2=0
8440 GOTO 8480
8450 S2=ATN(H2/H1)/P0
8460 IF H1<0 THEN S2=S2+SGN(H2)*180
8470 REM ----- MAGNITUDE AND PHASE OF E(PHI)
8480 S3=0
8490 IF (X3=0 AND X4=0) THEN 8510
8500 S3=SQR(X3*X3+X4*X4)
8510 IF X3><0 THEN 8540
8520 S4=0
8530 GOTO 8560
8540 S4=ATN(X4/X3)/P0
8550 IF X3<0 THEN S4=S4+SGN(X4)*180
8560 PRINT #3,USING "###.## ";Q2,Q1;
8570 PRINT #3,USING " ##.###^^^^";S1*F1;
8580 PRINT #3,USING " ###.## ";S2;
8590 PRINT #3,USING " ##.###^^^^";S3*F1;
8600 PRINT #3,USING " ###.##";S4
8610 IF S$="Y" THEN PRINT #1,Q2;",";Q1;",";S1*F1;",";S2;",";S3*F1;","S4
8620 REM ----- INCREMENT ZENITH ANGLE
8630 Q2=Q2+ZC
8640 NEXT I2
8650 REM ----- INCREMENT AZIMUTH ANGLE
8660 Q1=Q1+AC
8670 NEXT I1
8680 CLOSE #1
8690 RETURN
8700 REM ********** NEAR FIELD CALCULATION **********
8710 REM ----- ENSURE CURRENTS HAVE BEEN CALCULATED
8720 IF FLG<2 THEN GOSUB 1960
8730 O2=PWR
8740 PRINT #3," "
8750 PRINT #3,B$;" NEAR FIELDS ";B$
8760 PRINT #3," "
8770 INPUT "ELECTRIC OR MAGNETIC NEAR FIELDS (E/H) ";N$
8780 IF(N$="H" OR N$="E") GOTO 8800
8790 GOTO 8770
8800 PRINT
8810 REM ----- INPUT VARIABLES FOR NEAR FIELD CALCULATION
8820 PRINT "FIELD LOCATION(S):"
8830 A$="-COORDINATE (M): INITIAL,INCREMENT,NUMBER "
8840 PRINT " X";A$;
8850 INPUT XI,XC,NX
8860 IF NX=0 THEN NX=1
8870 IF O$>"C" THEN PRINT #3,"X";A$;": ";XI;",";XC;",";NX
8880 PRINT " Y";A$;
8890 INPUT YI,YC,NY
8900 IF NY=0 THEN NY=1
8910 IF O$>"C" THEN PRINT #3,"Y";A$;": ";YI;",";YC;",";NY
8920 PRINT " Z";A$;
8930 INPUT ZI,ZC,NZ
8940 IF NZ=0 THEN NZ=1
8950 IF O$>"C" THEN PRINT #3,"Z";A$;": ";ZI;",";ZC;",";NZ
8960 F1=1
8970 PRINT
8980 PRINT "PRESENT POWER LEVEL IS ";PWR;" WATTS"
8990 INPUT "CHANGE POWER LEVEL (Y/N) ";A$
9000 IF A$="N" THEN 9050
9010 IF A$<>"Y" THEN 8990
9020 INPUT "NEW POWER LEVEL (WATTS) ";O2
9030 IF O$>"C" THEN PRINT #3," ":PRINT #3,"NEW POWER LEVEL (WATTS) = ";O2
9040 GOTO 8990
9050 IF (O2<0 OR O2=0) THEN O2=PWR
9060 REM ----- RATIO OF POWER LEVELS
9070 F1=SQR(O2/PWR)
9080 IF N$="H" THEN F1=F1/S0/4/P
9090 PRINT
9100 REM ----- DESIGNATION OF OUTPUT FILE FOR NEAR FIELD DATA
9110 INPUT "SAVE TO A FILE (Y/N) ";S$
9120 IF S$="N" THEN 9200
9130 IF S$<>"Y" THEN 9110
9140 INPUT "FILENAME (NAME.OUT) ";F$
9150 IF LEFT$(RIGHT$(F$,4),1)="." THEN 916 ELSE F$=F$+".OUT"
9160 IF O$>"C" THEN PRINT #3," ":PRINT #3,"FILENAME (NAME.OUT) ";F$
9170 OPEN F$ FOR OUTPUT AS #2
9180 PRINT #2,NX*NY*NZ;",";O2;",";N$
9190 REM ----- LOOP OVER Z DIMENSION
9200 FOR IZ=1 TO NZ
9205 ZZ=ZI+(IZ-1)*ZC
9210 REM ----- LOOP OVER Y DIMENSION
9220 FOR IY=1 TO NY
9235 YY=YI+(IY-1)*YC
9230 REM ----- LOOP OVER X DIMENSION
9240 FOR IX=1 TO NX
9255 XX=XI+(IX-1)*XC
9250 REM ----- NEAR FIELD HEADER
9260 PRINT #3," "
9270 IF N$="E" THEN PRINT #3,B$;"NEAR ELECTRIC FIELDS";B$
9280 IF N$="H" THEN PRINT #3,B$;"NEAR MAGNETIC FIELDS";B$
9290 PRINT #3,TAB(10);"FIELD POINT: ";"X = ";XX;" Y = ";YY;" Z = ";ZZ
9300 PRINT #3," VECTOR","REAL","IMAGINARY","MAGNITUDE","PHASE"
9310 IF N$="E" THEN A$=" V/M "
9320 IF N$="H" THEN A$=" AMPS/M "
9330 PRINT #3," COMPONENT ",A$,A$,A$," DEG"
9340 A1=0
9350 A3=0
9360 A4=0
9370 REM ----- LOOP OVER THREE VECTOR COMPONENTS
9380 FOR I=1 TO 3
9390 X0=XX
9400 Y0=YY
9410 Z0=ZZ
9420 IF N$="H" THEN 952
9430 T5=0
9440 T6=0
9450 T7=0
9460 IF I=1 THEN T5=2*S0
9470 IF I=2 THEN T6=2*S0
9480 IF I=3 THEN T7=2*S0
9490 U7=0
9500 U8=0
9510 GOTO 9620
9520 FOR J8=1 TO 6
9530 K!(J8,1)=0
9540 K!(J8,2)=0
9550 NEXT J8
9560 J9=1
9570 J8=-1
9580 IF I=1 THEN X0=XX+J8*S0/2
9590 IF I=2 THEN Y0=YY+J8*S0/2
9600 IF I=3 THEN Z0=ZZ+J8*S0/2
9610 REM ----- LOOP OVER SOURCE SEGMENTS
9620 FOR J=1 TO N
9630 J1=ABS(C%(J,1))
9640 J2=ABS(C%(J,2))
9650 J3=J2
9660 IF J1>J2 THEN J3=J1
9670 F4=SGN(C%(J,1))
9680 F5=SGN(C%(J,2))
9690 F6=1
9700 F7=1
9710 U5=0
9720 U6=0
9730 REM ----- IMAGE LOOP
9740 FOR K=1 TO G STEP -2
9750 IF C%(J,1)<>-C%(J,2) THEN 9810
9760 IF K<0 THEN 10420
9770 REM ----- COMPUTE VECTOR POTENTIAL A
9780 F6=F4
9790 F7=F5
9800 REM ----- COMPUTE PSI(0,J,J+.5)
9810 P1=0
9820 P2=2*J3+J-1
9830 P3=P2+.5
9840 P4=J2
9850 GOSUB 750
9860 U1=T1*F5
9870 U2=T2*F5
9880 REM ----- COMPUTE PSI(0,J-.5,J)
9890 P3=P2
9900 P2=P2-.5
9910 P4=J1
9920 GOSUB 660
9930 V1=F4*T1
9940 V2=F4*T2
9950 REM ----- REAL PART OF VECTOR POTENTIAL CONTRIBUTION
9960 X3=U1*CA(J2)+V1*CA(J1)
9970 Y3=U1*CB(J2)+V1*CB(J1)
9980 Z3=(F7*U1*CG(J2)+F6*V1*CG(J1))*K
9990 REM ----- IMAGINARY PART OF VECTOR POTENTIAL CONTRIBUTION
10000 X5=U2*CA(J2)+V2*CA(J1)
10010 Y5=U2*CB(J2)+V2*CB(J1)
10020 Z5=(F7*U2*CG(J2)+F6*V2*CG(J1))*K
10030 REM ----- MAGNETIC FIELD CALCULATION COMPLETED
10040 IF N$="H" THEN 10360
10050 D1=(X3*T5+Y3*T6+Z3*T7)*W2
10060 D2=(X5*T5+Y5*T6+Z5*T7)*W2
10070 REM ----- COMPUTE PSI(.5,J,J+1)
10080 P1=.5
10090 P2=P3
10100 P3=P3+1
10110 P4=J2
10120 GOSUB 560
10130 U1=T1
10140 U2=T2
10150 REM ----- COMPUTE PSI(-.5,J,J+1)
10160 P1=-P1
10170 GOSUB 560
10180 U1=(T1-U1)/S(J2)
10190 U2=(T2-U2)/S(J2)
10200 REM ----- COMPUTE PSI(.5,J-1,J)
10210 P1=-P1
10220 P3=P2
10230 P2=P2-1
10240 P4=J1
10250 GOSUB 560
10260 U3=T1
10270 U4=T2
10280 REM ----- COMPUTE PSI(-.5,J-1,J)
10290 P1=-P1
10300 GOSUB 560
10310 REM ----- GRADIENT OF SCALAR POTENTIAL
10320 U5=(U1+(U3-T1)/S(J1)+D1)*K+U5
10330 U6=(U2+(U4-T2)/S(J1)+D2)*K+U6
10340 GOTO 10420
10350 REM ----- COMPONENTS OF VECTOR POTENTIAL A
10360 K!(1,J9)=K!(1,J9)+(X3*CR(J)-X5*CI(J))*K
10370 K!(2,J9)=K!(2,J9)+(X5*CR(J)+X3*CI(J))*K
10380 K!(3,J9)=K!(3,J9)+(Y3*CR(J)-Y5*CI(J))*K
10390 K!(4,J9)=K!(4,J9)+(Y5*CR(J)+Y3*CI(J))*K
10400 K!(5,J9)=K!(5,J9)+(Z3*CR(J)-Z5*CI(J))*K
10410 K!(6,J9)=K!(6,J9)+(Z5*CR(J)+Z3*CI(J))*K
10420 NEXT K
10430 IF N$="H" THEN 1046
10440 U7=U5*CR(J)-U6*CI(J)+U7
10450 U8=U6*CR(J)+U5*CI(J)+U8
10460 NEXT J
10470 IF N$="E" THEN 10690
10480 REM ----- DIFFERENCES OF VECTOR POTENTIAL A
10490 J8=1
10500 J9=J9+1
10510 IF J9=2 THEN 9580
10520 ON I GOTO 10530,10580,10630
10530 H(3)=K!(5,1)-K!(5,2)
10540 H(4)=K!(6,1)-K!(6,2)
10550 H(5)=K!(3,2)-K!(3,1)
10560 H(6)=K!(4,2)-K!(4,1)
10570 GOTO 10910
10580 H(1)=K!(5,2)-K!(5,1)
10590 H(2)=K!(6,2)-K!(6,1)
10600 H(5)=H(5)-K!(1,2)+K!(1,1)
10610 H(6)=H(6)-K!(2,2)+K!(2,1)
10620 GOTO 10910
10630 H(1)=H(1)-K!(3,2)+K!(3,1)
10640 H(2)=H(2)-K!(4,2)+K!(4,1)
10650 H(3)=H(3)+K!(1,2)-K!(1,1)
10660 H(4)=H(4)+K!(2,2)-K!(2,1)
10670 GOTO 10910
10680 REM ----- IMAGINARY PART OF ELECTRIC FIELD
10690 U7=-M*U7/S0
10700 REM ----- REAL PART OF ELECTRIC FIELD
10710 U8=M*U8/S0
10720 REM ----- MAGNITUDE AND PHASE CALCULATION
10730 S1=0
10740 IF (U7=0 AND U8=0) THEN 10760
10750 S1=SQR(U7*U7+U8*U8)
10760 S2=0
10770 IF U8<>0 THEN S2=ATN(U7/U8)/P0
10780 IF U8>0 THEN 10800
10790 S2=S2+SGN(U7)*180
10800 IF I=1 THEN PRINT #3," X ",
10810 IF I=2 THEN PRINT #3," Y ",
10820 IF I=3 THEN PRINT #3," Z ",
10830 PRINT #3,TAB(15);F1*U8;TAB(29);F1*U7;TAB(43);F1*S1;TAB(57);S2
10840 IF S$="Y" THEN PRINT #2,F1*U8;",";F1*U7;",";F1*S1;",";S2
10850 REM ----- CALCULATION FOR PEAK ELECTRIC FIELD
10860 S1=S1*S1
10870 S2=S2*P0
10880 A1=A1+S1*COS(2*S2)
10890 A3=A3+S1*SIN(2*S2)
10900 A4=A4+S1
10910 NEXT I
10920 IF N$="E" THEN 11150
10930 REM ----- MAGNETIC FIELD MAGNITUDE AND PHASE CALCULATION
10940 FOR I=1 TO 5 STEP 2
10950 S1=0
10960 IF (H(I)=0 AND H(I+1)=0) THEN 10980
10970 S1=SQR(H(I)*H(I)+H(I+1)*H(I+1))
10980 S2=0
10990 IF H(I)<>0 THEN S2=ATN(H(I+1)/H(I))/P0
11000 IF H(I)>0 THEN 11020
11010 S2=S2+SGN(H(I+1))*180
11020 IF I=1 THEN PRINT #3," X ",
11030 IF I=3 THEN PRINT #3," Y ",
11040 IF I=5 THEN PRINT #3," Z ",
11050 PRINT #3,TAB(15);F1*H(I);TAB(29);F1*H(I+1);TAB(43);F1*S1;TAB(57);S2
11060 IF S$="Y" THEN PRINT #2,F1*H(I);",";F1*H(I+1);",";F1*S1;",";S2
11070 REM ----- CALCULATION FOR PEAK MAGNETIC FIELD
11080 S1=S1*S1
11090 S2=S2*P0
11100 A1=A1+S1*COS(2*S2)
11110 A3=A3+S1*SIN(2*S2)
11120 A4=A4+S1
11130 NEXT I
11140 REM ----- PEAK FIELD CALCULATION
11150 PK=SQR(A4/2+SQR(A1*A1+A3*A3)/2)
11160 PRINT #3," MAXIMUM OR PEAK FIELD = ";F1*PK;A$
11170 IF (S$="Y" AND N$="E") THEN PRINT #2,F1*PK;",";O2
11180 IF (S$="Y" AND N$="H") THEN PRINT #2,F1*PK;",";O2
11190 IF S$="Y" THEN PRINT #2,XX;",";YY;",";ZZ
11220 NEXT IX
11250 NEXT IY
11280 NEXT IZ
11290 CLOSE #2
11300 RETURN
11310 REM ********** FREQUENCY INPUT **********
11320 REM ----- SET FLAG
11330 PRINT
11340 INPUT "FREQUENCY (MHZ)";F
11350 IF F=0 THEN F=299.8
11360 IF O$>"C" THEN PRINT #3, " ":PRINT #3, "FREQUENCY (MHZ):";F
11370 W=299.8/F
11380 REM -----VIRTUAL DIPOLE LENGTH FOR NEAR FIELD CALCULATION
11390 S0=.001*W
11400 REM ----- 1 / (4 * PI * OMEGA * EPSILON)
11410 M=4.77783352#*W
11420 REM ----- SET SMALL RADIUS MODIFICATION CONDITION
11430 SRM=.0001*W
11440 PRINT #3, " WAVE LENGTH = ";W;" METERS"
11450 REM ----- 2 PI / WAVELENGTH
11460 W=2*P/W
11470 W2=W*W/2
11480 FLG=0
11490 RETURN
11500 REM ********** GEOMETRY INPUT **********
11510 REM ----- WHEN GEOMETRY IS CHANGED, ENVIRONMENT MUST BE CHECKED
11520 GOSUB 13590
11530 PRINT
11540 IF INFILE THEN 11600
11550 INPUT "NO. OF WIRES";NW
11560 IF NW=0 THEN RETURN
11570 IF NW<=MW THEN 11600
11580 PRINT "NUMBER OF WIRES EXCEEDS DIMENSION..."
11590 GOTO 11550
11600 IF O$>"C" THEN PRINT #3," ":PRINT #3,"NO. OF WIRES:";NW
11610 REM ----- INITIALIZE NUMBER OF PULSES TO ZERO
11620 N=0
11630 FOR I=1 TO NW
11640 IF INFILE THEN GOSUB 15470:GOTO 11900
11650 PRINT
11660 PRINT "WIRE NO.";I
11670 INPUT " NO. OF SEGMENTS";S1
11680 IF S1=0 THEN 11530
11690 A$=" END ONE COORDINATES (X,Y,Z)"
11700 PRINT A$;
11710 INPUT X1,Y1,Z1
11720 IF G<0 AND Z1<0 THEN PRINT "Z CANNOT BE NEGATIVE":GOTO 11700
11730 A$=" END TWO COORDINATES (X,Y,Z)"
11740 PRINT A$;
11750 INPUT X2,Y2,Z2
11760 IF G<0 AND Z2<0 THEN PRINT "Z CANNOT BE NEGATIVE":GOTO 11740
11770 IF X1=X2 AND Y1=Y2 AND Z1=Z2 THEN PRINT"ZERO LENGTH WIRE.":GOTO 11660
11780 A$=" RADIUS"
11790 PRINT " "A$;
11800 INPUT A(I)
11810 IF A(I)<=0! THEN 11790
11820 REM ----- DETERMINE CONNECTIONS11830 IF O$>"C" THEN PRINT #3," ":PRINT #3,"WIRE NO.";I
11840 GOSUB 12890
11850 PRINT "CHANGE WIRE NO. ";I;" (Y/N) ";
11860 INPUT A$
11870 IF A$="Y" THEN 11650
11880 IF A$<>"N" THEN 11850
11890 REM ----- COMPUTE DIRECTION COSINES
11900 X3=X2-X1
11910 Y3=Y2-Y1
11920 Z3=Z2-Z1
11930 D=SQR(X3*X3+Y3*Y3+Z3*Z3)
11940 CA(I)=X3/D
11950 CB(I)=Y3/D
11960 CG(I)=Z3/D
11970 S(I)=D/S1
11980 REM ----- COMPUTE CONNECTIVITY DATA (PULSES N1 TO N)
11990 N1=N+1
12000 N(I,1)=N1
12010 IF (S1=1 AND I1=0) THEN N(I,1)=0
12020 N=N1+S1
12030 IF I1=0 THEN N=N-1
12040 IF I2=0 THEN N=N-1
12050 IF N>MP THEN PRINT "PULSE NUMBER EXCEEDS DIMENSION":CLOSE:GOTO 11550
12060 N(I,2)=N
12070 IF (S1=1 AND I2=0) THEN N(I,2)=0
12080 IF N<N1 THEN 12442
12090 FOR J=N1 TO N
12100 C%(J,1)=I
12110 C%(J,2)=I
12120 W%(J)=I
12130 NEXT J
12140 C%(N1,1)=I1
12150 C%(N,2)=I2
12160 REM ----- COMPUTE COORDINATES OF BREAK POINTS
12170 I1=N1+2*(I-1)
12180 I3=I1
12190 X(I1)=X1
12200 Y(I1)=Y1
12210 Z(I1)=Z1
12220 IF C%(N1,1)=0 THEN 12300
12230 I2=ABS(C%(N1,1))
12240 F3=SGN(C%(N1,1))*S(I2)
12250 X(I1)=X(I1)-F3*CA(I2)
12260 Y(I1)=Y(I1)-F3*CB(I2)
12270 IF C%(N1,1)=-I THEN F3=-F3
12280 Z(I1)=Z(I1)-F3*CG(I2)
12290 I3=I3+1
12300 I6=N+2*I
12310 FOR I4=I1+1 TO I6
12320 J=I4-I3
12330 X(I4)=X1+J*X3/S1
12340 Y(I4)=Y1+J*Y3/S1
12350 Z(I4)=Z1+J*Z3/S1
12360 NEXT I4
12370 IF C%(N,2)=0 THEN 12450
12380 I2=ABS(C%(N,2))
12390 F3=SGN(C%(N,2))*S(I2)
12400 I3=I6-1
12410 X(I6)=X(I3)+F3*CA(I2)
12420 Y(I6)=Y(I3)+F3*CB(I2)
12430 IF I=-C%(N,2) THEN F3=-F3
12440 Z(I6)=Z(I3)+F3*CG(I2)
12441 GOTO 12450
12442 I1=N1-2*(I-1): REM SINGLE SEGMENT/PULSE CASE
12443 X(I1)=X1
12444 Y(I1)=Y1
12445 Z(I1)=Z1
12446 I1=I1+1
12447 X(I1)=X2
12448 Y(I1)=Y2
12449 Z(I1)=Z2
12450 NEXT I
12460 REM ********** GEOMETRY OUTPUT **********
12470 PRINT #3, " "
12480 PRINT #3, " **** ANTENNA GEOMETRY ****"
12490 IF N>0 THEN 12540
12500 PRINT
12510 PRINT "NUMBER OF PULSES IS ZERO....RE-ENTER GEOMETRY"
12520 PRINT
12530 GOTO 11550
12540 K=1
12550 J=0
12560 FOR I=1 TO N
12570 I1=2*W%(I)-1+I
12580 IF K>NW THEN 12690
12590 IF K=J THEN 12690
12600 J=K
12610 PRINT #3," "
12620 PRINT #3,"WIRE NO. ";K;" COORDINATES",,,"CONNECTION PULSE"
12630 PRINT #3,"X","Y","Z","RADIUS","END1 END2 NO."
12640 IF (N(K,1)><0 OR N(K,2)><0) THEN 12690
12650 PRINT #3,"-","-","-"," -"," - - 0"
12660 K=K+1
12670 IF K>NW THEN 12760
12680 GOTO 12600
12690 PRINT #3,X(I1);TAB(15);Y(I1);TAB(29);Z(I1);TAB(43);A(W%(I));TAB(57);
12700 PRINT #3, USING "### ### ##";C%(I,1),C%(I,2),I
12710 IF (I=N(K,2) OR N(K,1)=N(K,2) OR C%(I,2)=0) THEN K=K+1
12720 IF C%(I,1)=0 THEN C%(I,1)=W%(I)
12730 IF C%(I,2)=0 THEN C%(I,2)=W%(I)
12740 IF (K=NW AND N(K,1)=0 AND N(K,2)=0) THEN 12600
12750 IF (I=N AND K<NW) THEN 12600
12760 NEXT I
12770 PRINT
12780 CLOSE 1:IF INFILE THEN INFILE=0:IF O$>"C" THEN 12830
12790 INPUT " CHANGE GEOMETRY (Y/N) ";A$
12800 IF A$="Y" THEN 11530
12810 IF A$<>"N" THEN 12790
12820 REM ----- EXCITATION INPUT
12830 GOSUB 14200
12840 REM ----- LOADS/NETWORKS INPUT
12850 GOSUB 14450
12860 FLG=0
12870 RETURN
12880 REM ********** CONNECTIONS **********
12890 E(I)=X1
12900 L(I)=Y1
12910 M(I)=Z1
12920 E(I+NW)=X2
12930 L(I+NW)=Y2
12940 M(I+NW)=Z2
12950 G%=0
12960 I1=0
12970 I2=0
12980 J1(I)=0
12990 J2(I,1)=-I
13000 J2(I,2)=-I
13010 IF G=1 THEN 13130
13020 REM ----- CHECK FOR GROUND CONNECTION
13030 IF Z1=0 THEN 13050
13040 GOTO 13080
13050 I1=-I
13060 J1(I)=-1
13070 GOTO 13300
13080 IF Z2=0 THEN 13100
13090 GOTO 13130
13100 I2=-I
13110 J1(I)=1
13120 G%=1
13130 IF I=1 THEN 13480
13140 FOR J=1 TO I-1
13150 REM ----- CHECK FOR END1 TO END1
13160 IF (X1=E(J) AND Y1=L(J) AND Z1=M(J)) THEN 13180
13170 GOTO 13230
13180 I1=-J
13190 J2(I,1)=J
13200 IF J2(J,1)=-J THEN J2(J,1)=J
13210 GOTO 13300
13220 REM ----- CHECK FOR END1 TO END2
13230 IF (X1=E(J+NW) AND Y1=L(J+NW) AND Z1=M(J+NW)) THEN 13250
13240 GOTO 13290
13250 I1=J
13260 J2(I,1)=J
13270 IF J2(J,2)=-J THEN J2(J,2)=J
13280 GOTO 13300
13290 NEXT J
13300 IF G%=1 THEN 13480
13310 IF I=1 THEN 13480
13320 FOR J=1 TO I-1
13330 REM ----- CHECK END2 TO END2
13340 IF (X2=E(J+NW) AND Y2=L(J+NW) AND Z2=M(J+NW)) THEN 13360
13350 GOTO 13410
13360 I2=-J
13370 J2(I,2)=J
13380 IF J2(J,2)=-J THEN J2(J,2)=J
13390 GOTO 13480
13400 REM ----- CHECK FOR END2 TO END1
13410 IF (X2=E(J) AND Y2=L(J) AND Z2=M(J)) THEN 13430
13420 GOTO 13470
13430 I2=J
13440 J2(I,2)=J
13450 IF J2(J,1)=-J THEN J2(J,1)=J
13460 GOTO 13480
13470 NEXT J
13480 PRINT #3," COORDINATES"," "," ","END NO. OF"
13490 PRINT #3," X"," Y"," Z","RADIUS CONNECTION SEGMENTS"
13500 PRINT #3,X1;TAB(15);Y1;TAB(29);Z1;TAB(57);I1
13510 PRINT #3,X2;TAB(15);Y2;TAB(29);Z2;TAB(43);A(I);TAB(57);I2;TAB(71);S1
13520 RETURN
13530 REM ********** ENVIROMENT INPUT **********
13540 PRINT
13550 PRINT " **** WARNING ****"
13560 PRINT "REDO GEOMETRY TO ENSURE PROPER GROUND CONNECTION/DISCONNECTION"
13570 PRINT
13580 REM ----- INITIALIZE NUMBER OF RADIAL WIRES TO ZERO
13590 NR=0
13600 REM ----- SET ENVIRONMENT
13610 PRINT #3," "
13620 A$="ENVIRONMENT (+1 FOR FREE SPACE, -1 FOR GROUND PLANE)"
13630 PRINT A$;
13640 INPUT G
13650 IF O$>"C" THEN PRINT #3,A$;": ";G
13660 IF G=1 THEN 14180
13670 IF G<>-1 THEN 13630
13680 REM ----- NUMBER OF MEDIA
13690 A$=" NUMBER OF MEDIA (0 FOR PERFECTLY CONDUCTING GROUND)"
13700 PRINT A$;
13710 INPUT NM
13720 IF NM<=MM THEN 13750
13730 PRINT "NUMBER OF MEDIA EXCEEDS DIMENSION..."
13740 GOTO 13700
13750 IF O$>"C" THEN PRINT #3,A$;": ";NM
13760 REM ----- INITIALIZE BOUNDARY TYPE
13770 TB=1
13780 IF NM=0 THEN 14180
13790 IF NM=1 THEN 13860
13800 REM ----- TYPE OF BOUNDARY
13810 A$=" TYPE OF BOUNDARY (1-LINEAR, 2-CIRCULAR)"
13820 PRINT " ";A$;
13830 INPUT TB
13840 IF O$>"C" THEN PRINT #3,A$;": ";TB
13850 REM ----- BOUNDARY CONDITIONS
13860 FOR I=1 TO NM
13870 PRINT "MEDIA";I
13880 A$=" RELATIVE DIELECTRIC CONSTANT, CONDUCTIVITY"
13890 PRINT " ";A$;
13900 INPUT T(I),V(I)
13910 IF O$>"C" THEN PRINT #3,A$;": ";T(I)","V(I)
13920 IF I>1 THEN 14040
13930 IF TB=1 THEN 14040
13940 A$=" NUMBER OF RADIAL WIRES IN GROUND SCREEN"
13950 PRINT " ";A$;
13960 INPUT NR
13970 IF O$>"C" THEN PRINT #3,A$;": ";NR
13980 IF NR=0 THEN 14040
13990 A$=" RADIUS OF RADIAL WIRES"
14000 PRINT " ";A$;
14010 INPUT RR
14020 IF O$>"C" THEN PRINT #3,A$;": ";RR
14030 REM ----- INITIALIZE COORDINATE OF MEDIA INTERFACE
14040 U(I)=1000000!
14050 REM ----- INITIALIZE HEIGHT OF MEDIA
14060 H(I)=0
14070 IF I=NM THEN 14120
14080 A$=" X OR R COORDINATE OF NEXT MEDIA INTERFACE"
14090 PRINT " ";A$;
14100 INPUT U(I)
14110 IF O$>"C" THEN PRINT #3,A$;": ";U(I)
14120 IF I=1 THEN 14170
14130 A$=" HEIGHT OF MEDIA"
14140 PRINT " ";A$;
14150 INPUT H(I)
14160 IF O$>"C" THEN PRINT #3,A$;": ";H(I)
14170 NEXT I
14180 RETURN
14190 REM ********** EXCITATION INPUT **********
14200 PRINT
14210 A$="NO. OF SOURCES "
14220 PRINT A$;
14230 INPUT NS
14240 IF NS<1 THEN NS=1
14250 IF NS<=MP THEN 14280
14260 PRINT "NO. OF SOURCES EXCEEDS DIMENSION ..."
14270 GOTO 14220
14280 IF O$>"C" THEN PRINT #3," ":PRINT #3, A$;": ";NS
14290 FOR I=1 TO NS
14300 PRINT
14310 PRINT "SOURCE NO. ";I;":"
14320 A$="PULSE NO., VOLTAGE MAGNITUDE, PHASE (DEGREES)"
14330 PRINT A$;
14340 INPUT E(I),VM,VP
14350 IF E(I)<=N THEN 14380
14360 PRINT "PULSE NUMBER EXCEEDS NUMBER OF PULSES..."
14370 GOTO 14330
14380 IF O$>"C" THEN PRINT #3,A$;": ";E(I)","VM","VP
14390 L(I)=VM*COS(VP*P0)
14400 M(I)=VM*SIN(VP*P0)
14410 NEXT I
14420 IF FLG=2 THEN FLG=1
14430 RETURN
14440 REM ********** LOADS INPUT **********
14450 PRINT
14460 INPUT "NUMBER OF LOADS ";NL
14470 IF NL<=ML THEN 14500
14480 PRINT "NUMBER OF LOADS EXCEEDS DIMENSION..."
14490 GOTO 14460
14500 IF O$>"C" THEN PRINT #3,"NUMBER OF LOADS";NL
14510 IF NL<1 THEN 14820
14520 INPUT "S-PARAMETER (S=jW) IMPEDANCE LOAD (Y/N)";L$
14530 IF L$<>"Y" AND L$<>"N" THEN 14520
14540 A$="PULSE NO.,RESISTANCE,REACTANCE"
14550 IF L$="Y" THEN A$= "PULSE NO., ORDER OF S-PARAMETER FUNCTION"
14560 FOR I=1 TO NL
14570 PRINT
14580 PRINT "LOAD NO. ";I;":"
14590 IF L$="Y" THEN 14660
14600 PRINT A$;
14610 INPUT LP(I),LA(1,I,1),LA(2,I,1)
14620 IF LP(I)>N THEN PRINT "PULSE NUMBER EXCEEDS NUMBER OF PULSES...": GOTO 14600
14630 IF O$>"C" THEN PRINT #3,A$;": ";LP(I);",";LA(1,I,1);",";LA(2,I,1)
14640 GOTO 14810
14650 REM ----- S-PARAMETER LOADS
14660 PRINT A$;
14670 INPUT LP(I),LS(I)
14680 IF LP(I)>N THEN PRINT "PULSE NUMBER EXCEEDS NUMBER OF PULSES...": GOTO 14660
14690 IF LS(I)>MA THEN PRINT "MAXIMUM DIMENSION IS 10":GOTO 14670
14700 IF O$>"C" THEN PRINT #3,A$;": ";LP(I);",";LS(I)
14710 FOR J=0 TO LS(I)
14720 A$="NUMERATOR, DENOMINATOR COEFFICIENTS OF S^"
14730 PRINT A$;J;
14740 INPUT LA(1,I,J),LA(2,I,J)
14750 IF O$>"C" THEN PRINT #3,A$;J;":";LA(1,I,J);",";LA(2,I,J)
14760 NEXT J
14770 IF LS(I)>0 THEN 14810
14780 LS(I)=1
14790 LA(1,I,1)=0
14800 LA(2,I,1)=0
14810 NEXT I
14820 FLG=0
14830 RETURN
14840 REM ********** MAIN PROGRAM **********
14850 REM ----- DATA INITIALIZATION
14860 REM ----- PI
14870 P=4*ATN(1)
14880 REM ----- CHANGES DEGREES TO RADIANS
14890 P0=P/180
14900 B$="********************"
14910 REM ----- INTRINSIC IMPEDANCE OF FREE SPACE DIVIDED BY 2 PI
14920 G0=29.979221#
14930 REM ---------- Q-VECTOR FOR GAUSSIAN QUADRATURE
14940 READ Q(1),Q(2),Q(3),Q(4),Q(5),Q(6),Q(7),Q(8),Q(9),Q(10),Q(11),Q(12)
14950 READ Q(13),Q(14)
14960 DATA .288675135,.5,.430568156,.173927423,.169990522,.326072577
14970 DATA .480144928,.050614268,.398333239,.111190517
14980 DATA .262766205,.156853323,.091717321,.181341892
14990 REM ---------- E-VECTOR FOR COEFFICIENTS OF ELLIPTIC INTEGRAL
15000 READ C0,C1,C2,C3,C4,C5,C6,C7,C8,C9
15010 DATA 1.38629436112,.09666344259,.03590092383,.03742563713,.01451196212
15020 DATA .5,.1249859397,.06880248576,.0332355346,.00441787012
15030 REM ----- IDENTIFY OUTPUT DEVICE
15040 GOSUB 15700
15050 PRINT #3,TAB(20);B$;B$
15060 PRINT #3,TAB(22);"MINI-NUMERICAL ELECTROMAGNETICS CODE"
15070 PRINT #3,TAB(36);"MININEC"
15080 PRINT #3,TAB(24);DATE$;TAB(48);TIME$
15090 PRINT #3,TAB(20);B$;B$
15100 REM ----- FREQUENCY INPUT
15110 GOSUB 11330
15120 REM ----- ENVIRONMENT INPUT
15130 GOSUB 13590
15140 REM ----- CHECK GEOMETRY INPUT
15141 INPUT "GEOMETRY FROM FILE, Y/N "; NA$
15142 IF NA$ <> "Y" THEN NA$="": GOTO 15170
15143 INPUT " ENTER FILEPATH + NAME OF FILE (.GEO IS ADDED)"; NA$: NA$=NA$+".GEO"
15144 OPEN NA$ AS #1 LEN=30
15150 GOSUB 15420
15160 REM ----- GEOMETRY INPUT
15170 GOSUB 11530
15180 REM ----- MENU
15190 PRINT
15200 PRINT B$;" MININEC MENU ";B$
15210 PRINT " G - CHANGE GEOMETRY C - COMPUTE/DISPLAY CURRENTS"
15220 PRINT " E - CHANGE ENVIRONMENT P - COMPUTE FAR-FIELD PATTERNS"
15230 PRINT " X - CHANGE EXCITATION N - COMPUTE NEAR-FIELDS"
15240 PRINT " L - CHANGE LOADS"
15250 PRINT " F - CHANGE FREQUENCY Q - QUIT"
15260 PRINT B$;B$;B$
15270 INPUT " COMMAND ";C$
15280 IF C$="F" THEN GOSUB 11330
15290 IF C$="P" THEN GOSUB 6200
15300 IF C$="X" THEN GOSUB 14200
15310 IF C$="E" THEN GOSUB 13540
15320 IF C$="G" THEN GOSUB 11520
15330 IF C$="C" THEN GOSUB 4960
15340 IF C$="L" THEN GOSUB 14450
15350 IF C$="N" THEN GOSUB 8720
15360 IF C$<>"Q" THEN 15190
15370 IF O$="P" THEN PRINT #3, CHR$(12) ELSE IF O$="C" THEN PRINT #3, " "
15380 CLOSE
15390 GOTO 16070
15400 REM ********** NEC-TYPE GEOMETRY INPUT **********
15410 OPEN "MININEC.INP" AS #1 LEN=30
15420 FIELD #1,2 AS S$,4 AS X1$,4 AS Y1$,4 AS Z1$,4 AS X2$,4 AS Y2$,4 AS Z2$,4 AS R$
15430 GET 1
15440 NW=CVI(S$)
15450 IF NW THEN INFILE=1
15460 RETURN
15470 REM ---------- GET GEOMETRY DATA FROM MININEC.INP ETC
15480 GET 1
15490 S1=CVI(S$)
15500 X1=CVS(X1$)
15510 Y1=CVS(Y1$)
15520 Z1=CVS(Z1$)
15530 X2=CVS(X2$)
15540 Y2=CVS(Y2$)
15550 Z2=CVS(Z2$)
15560 A(I)=CVS(R$)
15570 IF G<0 THEN IF Z1<0 OR Z2<0 THEN GOSUB 15620
15580 PRINT #3," ":PRINT #3,"WIRE NO.";I
15590 IF X1=X2 AND Y1=Y2 AND Z1=Z2 THEN PRINT"WIRE LENGTH IS ZERO.":GOTO 15370
15600 GOSUB 12890
15610 RETURN
15620 IF IZNEG THEN 15660
15630 PRINT"NEGATIVE Z VALUE ENCOUNTERED FOR GROUND PLANE."
15640 INPUT "ABORT OR CONVERT NEGATIVE Z VALUE TO ZERO (A/C)? ";A$
15650 IF A$="A" THEN 15370 ELSE IF A$="C" THEN IZNEG=1 ELSE 15640
15660 IF Z1<0 THEN Z1=-Z1
15670 IF Z2<0 THEN Z2=-Z2
15680 RETURN
15690 REM ********** IDENTIFY OUTPUT DEVICE **********
15700 INPUT "OUTPUT TO CONSOLE, PRINTER, OR DISK (C/P/D)";O$
15710 IF O$="C" THEN F$="SCRN:":GOTO 15760
15720 IF O$="P" THEN F$="LPT1:":GOTO 15760
15730 IF O$<>"D" THEN 15700
15740 INPUT "ENTER FILEPATH + FILENAME (.OUT IS ADDED)";F$
15750 IF LEFT$(RIGHT$(F$,4),1)="." THEN 15760 ELSE F$=F$+".OUT"
15760 OPEN F$ FOR OUTPUT AS #3
15770 CLS
15780 RETURN
15790 REM ********** CALCULATE ELAPSED TIME **********
15800 IH=VAL(MID$(T$,1,2))-VAL(MID$(OT$,1,2))
15810 IM=VAL(MID$(T$,4,2))-VAL(MID$(OT$,4,2))
15820 IS=VAL(MID$(T$,7,2))-VAL(MID$(OT$,7,2))
15830 IF IS<0 THEN IS=IS+60:IM=IM-1
15840 IF IM<0 THEN IM=IM+60:IH=IH-1
15850 IF IH<0 THEN IH=IH+24
15860 T$=":"+MID$(STR$(IS+100),3)
15870 IF IH THEN T$=MID$(STR$(IH),2)+":"+MID$(STR$(IM+100),3)+T$ ELSE T$=MID$(STR$(IM),2)+T$
15880 RETURN
15890 REM ********** CALCULATE APPROXIMATE TIME REMAINING **********
15900 IPCT=100*PCT
15910 T$=TIME$
15920 IH=VAL(MID$(T$,1,2))-VAL(MID$(OT$,1,2))
15930 IF IH<0 THEN IH=IH+24
15940 IM=VAL(MID$(T$,4,2))-VAL(MID$(OT$,4,2))
15950 IS=VAL(MID$(T$,7,2))-VAL(MID$(OT$,7,2))
15960 IS=IS+60*(IM+60*IH)
15970 IS=IS*(1/PCT-1)
15980 IM=INT(IS/60)
15990 IS=IS MOD 60
16000 IH=INT(IM/60)
16010 IM=IM MOD 60
16020 T$=":"+MID$(STR$(IS+100),3)
16030 IF IH THEN T$=MID$(STR$(IH),2)+":"+MID$(STR$(IM+100),3)+T$ ELSE T$=MID$(STR$(IM),2)+T$
16040 LOCATE CSRLIN,1
16050 PRINT Q$;IPCT;"% COMPLETE - APPROX TIME REMAINING "T$" ";
16060 RETURN
16070 END
**