1 Rem ****** MININEC(3) **********  NOSC CODE 822 (JCL CHANGE 9) 11-26-86
2 DefInt I-K, N
3 Dim K!(6, 2), Q(14)
4 Rem ----- MAXIMUM NUMBER OF SEGMENTS (PULSES + 2 * WIRES) = 150
5 MS = 150
6 Dim X(150), Y(150), Z(150)
7 Rem ----- MAXIMUM NUMBER OF WIRES = 50
8 MW = 50
9 Dim A(50), CA(50), CB(50), CG(50), J1(50), J2(50, 2), N(50, 2), S(50)
10 Rem ----- MAXIMUM NUMBER OF LOADS = 11
11 ML = 11
12 Rem ----- MAXIMUM ORDER OF S-PARAMETER LOADS = 8
13 MA = 8
14 Dim LA(2, 11, 8), LP(11), LS(11)
15 Rem ----- MAXIMUM NUMBER OF MEDIA = 6
16 MM = 6
17 Rem ----- H MUST BE DIMENSIONED AT LEAST 6
18 Dim H(6), T(6), U(6), V(6), Z1(6), Z2(6)
19 Rem ----- MAXIMUM NUMBER OF PULSES = 50
20 MP = 50
21 Dim C%(50, 2), CI(50), CR(50), P(50), W%(50)
22 Dim ZR(50, 50), ZI(50, 50)
23 Rem ---- ARRAYS E,L & M DIMENSIONED TO MW+MP=100
24 Dim E(100), L(100), M(100)
25 COLOR 2, 0
26 GoTo 1497
27 Rem ********** KERNEL EVALUATION OF INTEGRALS I2 & I3 **********
28 If K < 0 Then GoTo 33
29 X3 = X2 + T * (V1 - X2)
30 Y3 = Y2 + T * (V2 - Y2)
31 Z3 = Z2 + T * (V3 - Z2)
32 GoTo 36
33 X3 = V1 + T * (X2 - V1)
34 Y3 = V2 + T * (Y2 - V2)
35 Z3 = V3 + T * (Z2 - V3)
36 D3 = X3 * X3 + Y3 * Y3 + Z3 * Z3
37 Rem ----- MOD FOR SMALL RADIUS TO WAVELENGTH RATIO
38 If A(P4) <= SRM Then D = Sqr(D3): GoTo 49
39 D = D3 + A2
40 If D > 0 Then D = Sqr(D)
41 Rem ----- CRITERIA FOR USING REDUCED KERNEL
42 If I6! = 0 Then GoTo 49
43 Rem ----- EXACT KERNEL CALCULATION WITH ELLIPTIC INTEGRAL
44 B = D3 / (D3 + 4 * A2)
45 W0 = C0 + B * (C1 + B * (C2 + B * (C3 + B * C4)))
46 W1 = C5 + B * (C6 + B * (C7 + B * (C8 + B * C9)))
47 V0 = (W0 - W1 * Log(B)) * Sqr(1 - B)
48 T3 = T3 + (V0 + Log(D3 / (64 * A2)) / 2) / P / A(P4) - 1 / D
49 B1 = D * W
50 Rem ----- EXP(-J*K*R)/R
51 T3 = T3 + Cos(B1) / D
52 T4 = T4 - Sin(B1) / D
53 Return
54 Rem ***** PSI(P1,P2,P3) = T1 + J * T2 **********
55 Rem ----- ENTRIES REQUIRED FOR NEAR FIELD CALCULATION
56 X1 = X0 + P1 * T5 / 2
57 Y1 = Y0 + P1 * T6 / 2
58 Z1 = Z0 + P1 * T7 / 2
59 X2 = X1 - X(P2)
60 Y2 = Y1 - Y(P2)
61 Z2 = Z1 - K * Z(P2)
62 V1 = X1 - X(P3)
63 V2 = Y1 - Y(P3)
64 V3 = Z1 - K * Z(P3)
65 GoTo 135
66 I4 = Int(P2)
67 I5 = I4 + 1
68 X2 = X0 - (X(I4) + X(I5)) / 2
69 Y2 = Y0 - (Y(I4) + Y(I5)) / 2
70 Z2 = Z0 - K * (Z(I4) + Z(I5)) / 2
71 V1 = X0 - X(P3)
72 V2 = Y0 - Y(P3)
73 V3 = Z0 - K * Z(P3)
74 GoTo 135
75 X2 = X0 - X(P2)
76 Y2 = Y0 - Y(P2)
77 Z2 = Z0 - K * Z(P2)
78 I4 = Int(P3)
79 I5 = I4 + 1
80 V1 = X0 - (X(I4) + X(I5)) / 2
81 V2 = Y0 - (Y(I4) + Y(I5)) / 2
82 V3 = Z0 - K * (Z(I4) + Z(I5)) / 2
83 GoTo 135
84 Rem ----- ENTRIES REQUIRED FOR IMPEDANCE MATRIX CALCULATION
85 Rem ----- S(M) GOES IN (X1,Y1,Z1) FOR SCALAR POTENTIAL
86 Rem ----- MOD FOR SMALL RADIUS TO WAVE LENGTH RATIO
87 FVS = 1
88 If K < 1 Then GoTo 94
89 If A(P4) > SRM Then GoTo 94
90 If (P3 = P2 + 1 And P1 = (P2 + P3) / 2) Then GoTo 91 Else GoTo 94
91 T1 = 2 * Log(S(P4) / A(P4))
92 T2 = -W * S(P4)
93 Return
94 I4 = Int(P1)
95 I5 = I4 + 1
96 X1 = (X(I4) + X(I5)) / 2
97 Y1 = (Y(I4) + Y(I5)) / 2
98 Z1 = (Z(I4) + Z(I5)) / 2
99 GoTo 113
100 Rem ----- S(M) GOES IN (X1,Y1,Z1) FOR VECTOR POTENTIAL
101 Rem ----- MOD FOR SMALL RADIUS TO WAVE LENGTH RATIO
102 FVS = 0
103 If K < 1 Then GoTo 109
104 If A(P4) >= SRM Then GoTo 109
105 If (I = J And P3 = P2 + 0.5) Then GoTo 106 Else GoTo 109
106 T1 = Log(S(P4) / A(P4))
107 T2 = -W * S(P4) / 2
108 Return
109 X1 = X(P1)
110 Y1 = Y(P1)
111 Z1 = Z(P1)
112 Rem ----- S(U)-S(M) GOES IN (X2,Y2,Z2)
113 I4 = Int(P2)
114 If I4 = P2 Then GoTo 120
115 I5 = I4 + 1
116 X2 = (X(I4) + X(I5)) / 2 - X1
117 Y2 = (Y(I4) + Y(I5)) / 2 - Y1
118 Z2 = K * (Z(I4) + Z(I5)) / 2 - Z1
119 GoTo 124
120 X2 = X(P2) - X1
121 Y2 = Y(P2) - Y1
122 Z2 = K * Z(P2) - Z1
123 Rem ----- S(V)-S(M) GOES IN (V1,V2,V3)
124 I4 = Int(P3)
125 If I4 = P3 Then GoTo 131
126 I5 = I4 + 1
127 V1 = (X(I4) + X(I5)) / 2 - X1
128 V2 = (Y(I4) + Y(I5)) / 2 - Y1
129 V3 = K * (Z(I4) + Z(I5)) / 2 - Z1
130 GoTo 135
131 V1 = X(P3) - X1
132 V2 = Y(P3) - Y1
133 V3 = K * Z(P3) - Z1
134 Rem ----- MAGNITUDE OF S(U) - S(M)
135 D0 = X2 * X2 + Y2 * Y2 + Z2 * Z2
136 Rem ----- MAGNITUDE OF S(V) - S(M)
137 If D0 > 0 Then D0 = Sqr(D0)
138 D3 = V1 * V1 + V2 * V2 + V3 * V3
139 If D3 > 0 Then D3 = Sqr(D3)
140 Rem ----- SQUARE OF WIRE RADIUS
141 A2 = A(P4) * A(P4)
142 Rem ----- MAGNITUDE OF S(V) - S(U)
143 S4 = (P3 - P2) * S(P4)
144 Rem ----- ORDER OF INTEGRATION
145 Rem ----- LTH ORDER GAUSSIAN QUADRATURE
146 T1 = 0
147 T2 = 0
148 I6! = 0
149 F2 = 1
150 L = 7
151 T = (D0 + D3) / S(P4)
152 Rem ----- CRITERIA FOR EXACT KERNEL
153 If T > 1.1 Then GoTo 165
154 If C$ = "N" Then GoTo 165
155 If J2(W%(I), 1) = J2(W%(J), 1) Then GoTo 160
156 If J2(W%(I), 1) = J2(W%(J), 2) Then GoTo 160
157 If J2(W%(I), 2) = J2(W%(J), 1) Then GoTo 160
158 If J2(W%(I), 2) = J2(W%(J), 2) Then GoTo 160
159 GoTo 165
160 If A(P4) > SRM Then GoTo 162
161 If FVS = 1 Then GoTo 91 Else GoTo 106
162 F2 = 2 * (P3 - P2)
163 I6! = (1 - Log(S4 / F2 / 8 / A(P4))) / P / A(P4)
164 GoTo 167
165 If T > 6 Then L = 3
166 If T > 10 Then L = 1
167 I5 = L + L
168 T3 = 0
169 T4 = 0
170 T = (Q(L) + 0.5) / F2
171 GoSub 28
172 T = (0.5 - Q(L)) / F2
173 GoSub 28
174 L = L + 1
175 T1 = T1 + Q(L) * T3
176 T2 = T2 + Q(L) * T4
177 L = L + 1
178 If L < I5 Then GoTo 168
179 T1 = S4 * (T1 + I6!)
180 T2 = S4 * T2
181 Return
182 Rem ********** COMPLEX SQUARE ROOT **********
183 Rem ----- W6+I*W7=SQR(Z6+I*Z7)
184 T6 = Sqr((Abs(Z6) + Sqr(Z6 * Z6 + Z7 * Z7)) / 2)
185 T7 = Abs(Z7) / 2 / T6
186 If Z6 < 0 Then GoTo 191
187 W6 = T6
188 W7 = T7
189 If Z7 < 0 Then W7 = -T7
190 Return
191 W6 = T7
192 W7 = T6
193 If Z7 < 0 Then W7 = -T6
194 Return
195 Rem ********** IMPEDANCE MATRIX CALCULATION **********
196 If FLG = 1 Then GoTo 428
197 If FLG = 2 Then GoTo 477
198 Rem ----- BEGIN MATRIX FILL TIME CALCULATION
199 OT$ = Time$
200 Q$ = "MATRIX FILL  "
201 Print
202 Print "BEGIN "; Q$
203 Rem ----- ZERO IMPEDANCE MATRIX
204 For I = 1 To N
205 For J = 1 To N
206 ZR(I, J) = 0
207 ZI(I, J) = 0
208 Next J
209 Next I
210 Rem ----- COMPUTE ROW I OF MATRIX (OBSERVATION LOOP)
211 For I = 1 To N
212 I1 = Abs(C%(I, 1))
213 I2 = Abs(C%(I, 2))
214 F4 = Sgn(C%(I, 1)) * S(I1)
215 F5 = Sgn(C%(I, 2)) * S(I2)
216 Rem ----- R(M + 1/2) - R(M - 1/2) HAS COMPONENTS (T5,T6,T7)
217 T5 = F4 * CA(I1) + F5 * CA(I2)
218 T6 = F4 * CB(I1) + F5 * CB(I2)
219 T7 = F4 * CG(I1) + F5 * CG(I2)
220 If C%(I, 1) = -C%(I, 2) Then T7 = S(I1) * (CG(I1) + CG(I2))
221 Rem ----- COMPUTE COLUMN J OF ROW I (SOURCE LOOP)
222 For J = 1 To N
223 J1 = Abs(C%(J, 1))
224 J2 = Abs(C%(J, 2))
225 F4 = Sgn(C%(J, 1))
226 F5 = Sgn(C%(J, 2))
227 F6 = 1
228 F7 = 1
229 Rem ----- IMAGE LOOP
230 For K = 1 To G Step -2
231 If C%(J, 1) <> -C%(J, 2) Then GoTo 235
232 If K < 0 Then GoTo 332
233 F6 = F4
234 F7 = F5
235 F8 = 0
236 If K < 0 Then GoTo 248
237 Rem ----- SET FLAG TO AVOID REDUNANT CALCULATIONS
238 If I1 <> I2 Then GoTo 246
239 If (CA(I1) + CB(I1)) = 0 Then GoTo 241
240 If C%(I, 1) <> C%(I, 2) Then GoTo 246
241 If J1 <> J2 Then GoTo 246
242 If (CA(J1) + CB(J1)) = 0 Then GoTo 244
243 If C%(J, 1) <> C%(J, 2) Then GoTo 246
244 If I1 = J1 Then F8 = 1
245 If I = J Then F8 = 2
246 If ZR(I, J) <> 0 Then GoTo 317
247 Rem ----- COMPUTE PSI(M,N,N+1/2)
248 P1 = 2 * W%(I) + I - 1
249 P2 = 2 * W%(J) + J - 1
250 P3 = P2 + 0.5
251 P4 = J2
252 GoSub 102
253 U1 = F5 * T1
254 U2 = F5 * T2
255 Rem ----- COMPUTE PSI(M,N-1/2,N)
256 P3 = P2
257 P2 = P2 - 0.5
258 P4 = J1
259 If F8 < 2 Then GoSub 102
260 V1 = F4 * T1
261 V2 = F4 * T2
262 Rem ----- S(N+1/2)*PSI(M,N,N+1/2) + S(N-1/2)*PSI(M,N-1/2,N)
263 X3 = U1 * CA(J2) + V1 * CA(J1)
264 Y3 = U1 * CB(J2) + V1 * CB(J1)
265 Z3 = (F7 * U1 * CG(J2) + F6 * V1 * CG(J1)) * K
266 Rem ----- REAL PART OF VECTOR POTENTIAL CONTRIBUTION
267 D1 = W2 * (X3 * T5 + Y3 * T6 + Z3 * T7)
268 X3 = U2 * CA(J2) + V2 * CA(J1)
269 Y3 = U2 * CB(J2) + V2 * CB(J1)
270 Z3 = (F7 * U2 * CG(J2) + F6 * V2 * CG(J1)) * K
271 Rem ----- IMAGINARY PART OF VECTOR POTENTIAL CONTRIBUTION
272 D2 = W2 * (X3 * T5 + Y3 * T6 + Z3 * T7)
273 Rem ----- COMPUTE PSI(M+1/2,N,N+1)
274 P1 = P1 + 0.5
275 If F8 = 2 Then P1 = P1 - 1
276 P2 = P3
277 P3 = P3 + 1
278 P4 = J2
279 If F8 <> 1 Then GoTo 283
280 U5 = F5 * U1 + T1
281 U6 = F5 * U2 + T2
282 GoTo 291
283 GoSub 87
284 If F8 < 2 Then GoTo 288
285 U1 = (2 * T1 - 4 * U1 * F5) / S(J1)
286 U2 = (2 * T2 - 4 * U2 * F5) / S(J1)
287 GoTo 314
288 U5 = T1
289 U6 = T2
290 Rem ----- COMPUTE PSI(M-1/2,N,N+1)
291 P1 = P1 - 1
292 GoSub 87
293 U1 = (T1 - U5) / S(J2)
294 U2 = (T2 - U6) / S(J2)
295 Rem ----- COMPUTE PSI(M+1/2,N-1,N)
296 P1 = P1 + 1
297 P3 = P2
298 P2 = P2 - 1
299 P4 = J1
300 GoSub 87
301 U3 = T1
302 U4 = T2
303 Rem ----- COMPUTE PSI(M-1/2,N-1,N)
304 If F8 < 1 Then GoTo 308
305 T1 = U5
306 T2 = U6
307 GoTo 311
308 P1 = P1 - 1
309 GoSub 87
310 Rem ----- GRADIENT OF SCALAR POTENTIAL CONTRIBUTION
311 U1 = U1 + (U3 - T1) / S(J1)
312 U2 = U2 + (U4 - T2) / S(J1)
313 Rem ----- SUM INTO IMPEDANCE MATRIX
314 ZR(I, J) = ZR(I, J) + K * (D1 + U1)
315 ZI(I, J) = ZI(I, J) + K * (D2 + U2)
316 Rem ----- AVOID REDUNANT CALCULATIONS
317 If J < I Then GoTo 332
318 If F8 = 0 Then GoTo 332
319 ZR(J, I) = ZR(I, J)
320 ZI(J, I) = ZI(I, J)
321 Rem ----- SEGMENTS ON SAME WIRE SAME DISTANCE APART HAVE SAME Z
322 P1 = J + 1
323 If P1 > N Then GoTo 332
324 If C%(P1, 1) <> C%(P1, 2) Then GoTo 332
325 If C%(P1, 2) = C%(J, 2) Then GoTo 328
326 If C%(P1, 2) <> -C%(J, 2) Then GoTo 332
327 If (CA(J2) + CB(J2)) <> 0 Then GoTo 332
328 P2 = I + 1
329 If P2 > N Then GoTo 332
330 ZR(P2, P1) = ZR(I, J)
331 ZI(P2, P1) = ZI(I, J)
332 Next K
333 Next J
334 PCT = I / N
335 GoSub 1599
336 Next I
337 Rem ----- END MATRIX FILL TIME CALCULATION
338 T$ = Time$
339 GoSub 1589
340 Print #3, " "
341 Print #3, "FILL MATRIX  : "; T$
342 Rem ********** ADDITION OF LOADS **********
343 If NL = 0 Then GoTo 377
344 F5 = 2 * P * F
345 For I = 1 To NL
346 If L$ = "N" Then GoTo 366
347 Rem ----- S-PARAMETER LOADS
348 U1 = 0
349 U2 = 0
350 D1 = 0
351 D2 = 0
352 S = 1
353 For J = 0 To LS(I) Step 2
354 U1 = U1 + LA(1, I, J) * S * F5 ^ J
355 D1 = D1 + LA(2, I, J) * S * F5 ^ J
356 L = J + 1
357 U2 = U2 + LA(1, I, L) * S * F5 ^ L
358 D2 = D2 + LA(2, I, L) * S * F5 ^ L
359 S = -S
360 Next J
361 J = LP(I)
362 D = D1 * D1 + D2 * D2
363 LI = (U2 * D1 - D2 * U1) / D
364 LR = (U1 * D1 + U2 * D2) / D
365 GoTo 369
366 LR = LA(1, I, 1)
367 LI = LA(2, I, 1)
368 J = LP(I)
369 F2 = 1 / M
370 If C%(J, 1) <> -C%(J, 2) Then GoTo 372
371 If K < 0 Then F2 = 2 / M
372 ZR(J, J) = ZR(J, J) + F2 * LI
373 ZI(J, J) = ZI(J, J) - F2 * LR
374 Next I
375 Rem ********** IMPEDANCE MATRIX FACTORIZATION **********
376 Rem ----- BEGIN MATRIX FACTOR TIME CALCULATION
377 OT$ = Time$
378 Q$ = "FACTOR MATRIX"
379 Print
380 Print "BEGIN "; Q$;
381 X = N
382 PCTN = X * (X - 1) * (X + X - 1)
383 For K = 1 To N - 1
384 Rem ----- SEARCH FOR PIVOT
385 T = ZR(K, K) * ZR(K, K) + ZI(K, K) * ZI(K, K)
386 I1 = K
387 For I = K + 1 To N
388 T1 = ZR(I, K) * ZR(I, K) + ZI(I, K) * ZI(I, K)
389 If T1 < T Then GoTo 392
390 I1 = I
391 T = T1
392 Next I
393 Rem ----- EXCHANGE ROWS K AND I1
394 If I1 = K Then GoTo 403
395 For J = 1 To N
396 T1 = ZR(K, J)
397 T2 = ZI(K, J)
398 ZR(K, J) = ZR(I1, J)
399 ZI(K, J) = ZI(I1, J)
400 ZR(I1, J) = T1
401 ZI(I1, J) = T2
402 Next J
403 P(K) = I1
404 Rem ----- SUBTRACT ROW K FROM ROWS K+1 TO N
405 For I = K + 1 To N
406 Rem ----- COMPUTE MULTIPLIER L(I,K)
407 T1 = (ZR(I, K) * ZR(K, K) + ZI(I, K) * ZI(K, K)) / T
408 T2 = (ZI(I, K) * ZR(K, K) - ZR(I, K) * ZI(K, K)) / T
409 ZR(I, K) = T1
410 ZI(I, K) = T2
411 Rem ----- SUBTRACT ROW K FROM ROW I
412 For J = K + 1 To N
413 ZR(I, J) = ZR(I, J) - (ZR(K, J) * T1 - ZI(K, J) * T2)
414 ZI(I, J) = ZI(I, J) - (ZR(K, J) * T2 + ZI(K, J) * T1)
415 Next J
416 Next I
417 X = N - K
418 PCT = 1 - X * (X - 1) * (X + X - 1) / PCTN
419 GoSub 1599
420 Next K
421 Rem ----- END MATRIX FACTOR TIME CALCULATION
422 T$ = Time$
423 GoSub 1589
424 Print
425 Print #3, "FACTOR MATRIX: "; T$
426 Rem ********** SOLVE **********
427 Rem ----- COMPUTE RIGHT HAND SIDE
428 For I = 1 To N
429 CR(I) = 0
430 CI(I) = 0
431 Next I
432 For J = 1 To NS
433 F2 = 1 / M
434 If C%(E(J), 1) = -C%(E(J), 2) Then F2 = 2 / M
435 CR(E(J)) = F2 * M(J)
436 CI(E(J)) = -F2 * L(J)
437 Next J
438 Rem ----- PERMUTE EXCITATION
439 For K = 1 To N - 1
440 I1 = P(K)
441 If I1 = K Then GoTo 448
442 T1 = CR(K)
443 T2 = CI(K)
444 CR(K) = CR(I1)
445 CI(K) = CI(I1)
446 CR(I1) = T1
447 CI(I1) = T2
448 Next K
449 Rem ----- FORWARD ELIMINATION
450 For I = 2 To N
451 T1 = 0
452 T2 = 0
453 For J = 1 To I - 1
454 T1 = T1 + ZR(I, J) * CR(J) - ZI(I, J) * CI(J)
455 T2 = T2 + ZR(I, J) * CI(J) + ZI(I, J) * CR(J)
456 Next J
457 CR(I) = CR(I) - T1
458 CI(I) = CI(I) - T2
459 Next I
460 Rem ----- BACK SUBSTITUTION
461 For I = N To 1 Step -1
462 T1 = 0
463 T2 = 0
464 If I = N Then GoTo 469
465 For J = I + 1 To N
466 T1 = T1 + ZR(I, J) * CR(J) - ZI(I, J) * CI(J)
467 T2 = T2 + ZR(I, J) * CI(J) + ZI(I, J) * CR(J)
468 Next J
469 T = ZR(I, I) * ZR(I, I) + ZI(I, I) * ZI(I, I)
470 T1 = CR(I) - T1
471 T2 = CI(I) - T2
472 CR(I) = (T1 * ZR(I, I) + T2 * ZI(I, I)) / T
473 CI(I) = (T2 * ZR(I, I) - T1 * ZI(I, I)) / T
474 Next I
475 FLG = 2
476 Rem ********** SOURCE DATA **********
477 Print #3, " "
478 Print #3, B$; "    SOURCE DATA     "; B$
479 PWR = 0
480 For I = 1 To NS
481 CR = CR(E(I))
482 CI = CI(E(I))
483 T = CR * CR + CI * CI
484 T1 = (L(I) * CR + M(I) * CI) / T
485 T2 = (M(I) * CR - L(I) * CI) / T
486 O2 = (L(I) * CR + M(I) * CI) / 2
487 PWR = PWR + O2
488 Print #3, "PULSE "; E(I), "VOLTAGE = ("; L(I); ","; M(I); "J)"
489 Print #3, " ", "CURRENT = ("; CR; ","; CI; "J)"
490 Print #3, " ", "IMPEDANCE = ("; T1; ","; T2; "J)"
491 Print #3, " ", "POWER = "; O2; " WATTS"
492 Next I
493 If NS > 1 Then Print #3, " "
494 If NS > 1 Then Print #3, "TOTAL POWER = "; PWR; "WATTS"
495 Return
496 Rem ********** PRINT CURRENTS **********
497 GoSub 196
498 S$ = "N"
499 Print #3, " "
500 Print #3, B$; "    CURRENT DATA    "; B$
501 For K = 1 To NW
502 If S$ = "Y" Then GoTo 507
503 Print #3, " "
504 Print #3, "WIRE NO. "; K; ":"
505 Print #3, "PULSE", "REAL", "IMAGINARY", "MAGNITUDE", "PHASE"
506 Print #3, " NO.", "(AMPS)", "(AMPS)", "(AMPS)", "(DEGREES)"
507 N1 = N(K, 1)
508 N2 = N(K, 2)
509 I = N1
510 C = C%(I, 1)
511 If (N1 = 0 And N2 = 0) Then C = K
512 If G = 1 Then GoTo 515
513 If (J1(K) = -1 And N1 > N2) Then N2 = N1
514 If J1(K) = -1 Then GoTo 525
515 E% = 1
516 GoSub 572
517 I2! = I1!
518 J2! = J1!
519 GoSub 607
520 If S$ = "N" Then Print #3, I$, I1!; Tab(29); J1!; Tab(43); S1; Tab(57); S2
521 If S$ = "Y" Then Print #1, I1!; ","; J1!; ","; S1; ","; S2
522 If N1 = 0 Then GoTo 532
523 If C = K Then GoTo 525
524 If I$ = "J" Then N1 = N1 + 1
525 For I = N1 To N2 - 1
526 I2! = CR(I)
527 J2! = CI(I)
528 GoSub 607
529 If S$ = "N" Then Print #3, I, CR(I); Tab(29); CI(I); Tab(43); S1; Tab(57); S2
530 If S$ = "Y" Then Print #1, CR(I); ","; CI(I); ","; S1; ","; S2
531 Next I
532 I = N2
533 C = C%(I, 2)
534 If (N1 = 0 And N2 = 0) Then C = K
535 If G = 1 Then GoTo 537
536 If J1(K) = 1 Then GoTo 543
537 E% = 2
538 GoSub 572
539 If (N1 = 0 And N2 = 0) Then GoTo 549
540 If N1 > N2 Then GoTo 549
541 If C = K Then GoTo 543
542 If I$ = "J" Then GoTo 549
543 I2! = CR(N2)
544 J2! = CI(N2)
545 GoSub 607
546 If S$ = "N" Then Print #3, N2, CR(N2); Tab(29); CI(N2); Tab(43); S1; Tab(57); S2
547 If S$ = "Y" Then Print #1, CR(N2); ","; CI(N2); ","; S1; ","; S2
548 If J1(K) = 1 Then GoTo 554
549 I2! = I1!
550 J2! = J1!
551 GoSub 607
552 If S$ = "N" Then Print #3, I$, I1!; Tab(29); J1!; Tab(43); S1; Tab(57); S2
553 If S$ = "Y" Then Print #1, I1!; ","; J1!; ","; S1; ","; S2
554 If S$ = "Y" Then Print #1, " 1 , 1 , 1 , 1"
555 Next K
556 If S$ = "Y" Then GoTo 569
557 Print
558 INPUT "SAVE CURRENTS TO A FILE (Y/N) ";S$
559 If S$ = "N" Then GoTo 570
560 If S$ <> "Y" Then GoTo 557
561 Print #3, " "
562 INPUT "FILENAME (NAME.OUT) ";F$
563 If Left$(Right$(F$, 4), 1) = "." Then GoTo 564 Else F$ = F$ + ".OUT"
564 If O$ > "C" Then Print #3, "FILENAME (NAME.OUT): "; F$
565 Open F$ For Output As #1
566 Print #3, " "
567 Print #1, NW; ","; PWR; ",C"
568 GoTo 501
569 Close #1
570 Return
571 Rem ----- SORT JUNCTION CURRENTS
572 I$ = "E"
573 I1! = 0!
574 J1! = 0!
575 If (C = K Or C = 0) Then GoTo 580
576 I$ = "J"
577 I1! = CR(I)
578 J1! = CI(I)
579 Rem ----- CHECK FOR OTHER OVERLAPPING WIRES
580 For J = 1 To NW
581 If J = K Then GoTo 604
582 L1 = N(J, 1)
583 L2 = N(J, 2)
584 If E% = 2 Then GoTo 590
585 CO = C%(L1, 1)
586 CT = C%(L2, 2)
587 L3 = L1
588 L4 = L2
589 GoTo 594
590 CO = C%(L2, 2)
591 CT = C%(L1, 1)
592 L3 = L2
593 L4 = L1
594 If CO = -K Then GoTo 596
595 GoTo 599
596 I1! = I1! - CR(L3)
597 J1! = J1! - CI(L3)
598 I$ = "J"
599 If CT = K Then GoTo 601
600 GoTo 604
601 I1! = I1! + CR(L4)
602 J1! = J1! + CI(L4)
603 I$ = "J"
604 Next J
605 Return
606 Rem ----- CALCULATE S1 AND S2
607 I3! = I2! * I2!
608 J3! = J2! * J2!
609 If (I3! > 0 Or J3! > 0) Then GoTo 612
610 S1 = 0!
611 GoTo 613
612 S1 = Sqr(I3! + J3!)
613 If I2! <> 0 Then GoTo 616
614 S2 = 0!
615 Return
616 S2 = Atn(J2! / I2!) / P0
617 If I2! > 0 Then Return
618 S2 = S2 + Sgn(J2!) * 180
619 Return
620 Rem ********** FAR FIELD CALCULATION **********
621 If FLG < 2 Then GoSub 196
622 O2 = PWR
623 Rem ----- TABULATE IMPEDANCE
624 If NM = 0 Then GoTo 634
625 For I = 1 To NM
626 Z6 = T(I)
627 Z7 = -V(I) / (2 * P * F * 0.00000885)
628 Rem ----- FORM IMPEDANCE=1/SQR(DIELECTRIC CONSTANT)
629 GoSub 184
630 D = W6 * W6 + W7 * W7
631 Z1(I) = W6 / D
632 Z2(I) = -W7 / D
633 Next I
634 Print #3, " "
635 Print #3, B$; "     FAR FIELD      "; B$
636 Print #3, " "
637 Rem ----- INPUT VARIABLES FOR FAR FIELD CALCULATION
638 INPUT "CALCULATE PATTERN IN DBI OR VOLTS/METER (D/V)";P$
639 If P$ = "D" Then GoTo 655
640 If P$ <> "V" Then GoTo 638
641 F1 = 1
642 Print
643 Print "PRESENT POWER LEVEL =  "; PWR; " WATTS"
644 INPUT "CHANGE POWER LEVEL (Y/N) ";A$
645 If A$ = "N" Then GoTo 650
646 If A$ <> "Y" Then GoTo 644
647 INPUT "NEW POWER LEVEL (WATTS)  ";O2
648 If O$ > "C" Then Print #3, "NEW POWER LEVEL = "; O2
649 GoTo 644
650 If (O2 < 0 Or O2 = 0) Then O2 = PWR
651 F1 = Sqr(O2 / PWR)
652 Print
653 INPUT "RADIAL DISTANCE (METERS) ";RD
654 If RD < 0 Then RD = 0
655 A$ = "ZENITH ANGLE : INITIAL,INCREMENT,NUMBER"
656 Print A$;
657 INPUT ZA,ZC,NZ
658 If NZ = 0 Then NZ = 1
659 If O$ > "C" Then Print #3, A$; ": "; ZA; ","; ZC; ","; NZ
660 A$ = "AZIMUTH ANGLE: INITIAL,INCREMENT,NUMBER"
661 Print A$;
662 INPUT AA,AC,NA
663 If NA = 0 Then NA = 1
664 If O$ > "C" Then Print #3, A$; ": "; AA; ","; AC; ","; NA
665 Print #3, " "
666 Rem ********** FILE FAR FIELD DATA **********
667 INPUT "FILE PATTERN (Y/N)";S$
668 If S$ = "N" Then GoTo 676
669 If S$ <> "Y" Then GoTo 667
670 Print #3, " "
671 INPUT "FILENAME (NAME.OUT)";F$
672 If Left$(Right$(F$, 4), 1) = "." Then GoTo 673 Else F$ = F$ + ".OUT"
673 If O$ > "C" Then Print #3, "FILENAME (NAME.OUT): "; F$
674 Open F$ For Output As #1
675 Print #1, NA * NZ; ","; O2; ","; P$
676 Print #3, " "
677 K9! = 0.016678 / PWR
678 Rem ----- PATTERN HEADER
679 Print #3, B$; "    PATTERN DATA    "; B$
680 If P$ = "V" Then GoTo 685
681 Print #3, "ZENITH", "AZIMUTH", "VERTICAL", "HORIZONTAL", "TOTAL"
682 A$ = "PATTERN (DB)"
683 Print #3, " ANGLE", " ANGLE", A$, A$, A$
684 GoTo 692
685 If RD > 0 Then Print #3, Tab(15); "RADIAL DISTANCE = "; RD; " METERS"
686 Print #3, Tab(15); "POWER LEVEL = "; PWR * F1 * F1; " WATTS"
687 Print #3, "ZENITH   AZIMUTH", "     E(THETA)     ", "     E(PHI)"
688 A$ = " MAG(V/M)    PHASE(DEG)"
689 Print #3, " ANGLE    ANGLE", A$, A$
690 If S$ = "Y" Then Print #1, RD
691 Rem ----- LOOP OVER AZIMUTH ANGLE
692 Q1 = AA
693 For I1 = 1 To NA
694 U3 = Q1 * P0
695 V1 = -Sin(U3)
696 V2 = Cos(U3)
697 Rem ----- LOOP OVER ZENITH ANGLE
698 Q2 = ZA
699 For I2 = 1 To NZ
700 U4 = Q2 * P0
701 R3 = Cos(U4)
702 T3 = -Sin(U4)
703 T1 = R3 * V2
704 T2 = -R3 * V1
705 R1 = -T3 * V2
706 R2 = T3 * V1
707 X1 = 0
708 Y1 = 0
709 Z1 = 0
710 X2 = 0
711 Y2 = 0
712 Z2 = 0
713 Rem ----- IMAGE LOOP
714 For K = 1 To G Step -2
715 For I = 1 To N
716 If K > 0 Then GoTo 718
717 If C%(I, 1) = -C%(I, 2) Then GoTo 813
718 J = 2 * W%(I) - 1 + I
719 Rem ----- FOR EACH END OF PULSE COMPUTE A CONTRIBUTION TO E-FIELD
720 For F5 = 1 To 2
721 L = Abs(C%(I, F5))
722 F3 = Sgn(C%(I, F5)) * W * S(L) / 2
723 If C%(I, 1) <> -C%(I, 2) Then GoTo 725
724 If F3 < 0 Then GoTo 812
725 If K = 1 Then GoTo 728
726 If NM <> 0 Then GoTo 747
727 Rem ----- STANDARD CASE
728 S2 = W * (X(J) * R1 + Y(J) * R2 + Z(J) * K * R3)
729 S1 = Cos(S2)
730 S2 = Sin(S2)
731 B1 = F3 * (S1 * CR(I) - S2 * CI(I))
732 B2 = F3 * (S1 * CI(I) + S2 * CR(I))
733 If C%(I, 1) = -C%(I, 2) Then GoTo 742
734 X1 = X1 + K * B1 * CA(L)
735 X2 = X2 + K * B2 * CA(L)
736 Y1 = Y1 + K * B1 * CB(L)
737 Y2 = Y2 + K * B2 * CB(L)
738 Z1 = Z1 + B1 * CG(L)
739 Z2 = Z2 + B2 * CG(L)
740 GoTo 812
741 Rem ----- GROUNDED ENDS
742 Z1 = Z1 + 2 * B1 * CG(L)
743 Z2 = Z2 + 2 * B2 * CG(L)
744 GoTo 812
745 Rem ----- REAL GROUND CASE
746 Rem ----- BEGIN BY FINDING SPECULAR DISTANCE
747 T4 = 100000!
748 If R3 = 0 Then GoTo 750
749 T4 = -Z(J) * T3 / R3
750 B9 = T4 * V2 + X(J)
751 If TB = 1 Then GoTo 755
752 B9 = B9 * B9 + (Y(J) - T4 * V1) ^ 2
753 If B9 > 0 Then B9 = Sqr(B9) Else GoTo 755
754 Rem ----- SEARCH FOR THE CORRESPONDING MEDIUM
755 J2 = NM
756 For J1 = NM To 1 Step -1
757 If B9 > U(J1) Then GoTo 759
758 J2 = J1
759 Next J1
760 Rem ----- OBTAIN IMPEDANCE AT SPECULAR POINT
761 Z4 = Z1(J2)
762 Z5 = Z2(J2)
763 Rem ----- IF PRESENT INCLUDE GROUND SCREEN IMPEDANCE IN PARALLEL
764 If NR = 0 Then GoTo 776
765 If B9 > U(1) Then GoTo 776
766 R = B9 + NR * RR
767 Z8 = W * R * Log(R / (NR * RR)) / NR
768 S8 = -Z5 * Z8
769 S9 = Z4 * Z8
770 T8 = Z4
771 T9 = Z5 + Z8
772 D = T8 * T8 + T9 * T9
773 Z4 = (S8 * T8 + S9 * T9) / D
774 Z5 = (S9 * T8 - S8 * T9) / D
775 Rem ----- FORM SQR(1-Z^2*SIN^2)
776 Z6 = 1 - (Z4 * Z4 - Z5 * Z5) * T3 * T3
777 Z7 = -(2 * Z4 * Z5) * T3 * T3
778 GoSub 184
779 Rem ----- VERTICAL REFLECTION COEFFICIENT
780 S8 = R3 - (W6 * Z4 - W7 * Z5)
781 S9 = -(W6 * Z5 + W7 * Z4)
782 T8 = R3 + (W6 * Z4 - W7 * Z5)
783 T9 = W6 * Z5 + W7 * Z4
784 D = T8 * T8 + T9 * T9
785 V8 = (S8 * T8 + S9 * T9) / D
786 V9 = (S9 * T8 - S8 * T9) / D
787 Rem ----- HORIZONTAL REFLECTION COEFFICIENT
788 S8 = W6 - R3 * Z4
789 S9 = W7 - R3 * Z5
790 T8 = W6 + R3 * Z4
791 T9 = W7 + R3 * Z5
792 D = T8 * T8 + T9 * T9
793 H8 = (S8 * T8 + S9 * T9) / D - V8
794 H9 = (S9 * T8 - S8 * T9) / D - V9
795 Rem ----- COMPUTE CONTRIBUTION TO SUM
796 S2 = W * (X(J) * R1 + Y(J) * R2 - (Z(J) - 2 * H(J2)) * R3)
797 S1 = Cos(S2)
798 S2 = Sin(S2)
799 B1 = F3 * (S1 * CR(I) - S2 * CI(I))
800 B2 = F3 * (S1 * CI(I) + S2 * CR(I))
801 W6 = B1 * V8 - B2 * V9
802 W7 = B1 * V9 + B2 * V8
803 D = CA(L) * V1 + CB(L) * V2
804 Z6 = D * (B1 * H8 - B2 * H9)
805 Z7 = D * (B1 * H9 + B2 * H8)
806 X1 = X1 - (CA(L) * W6 + V1 * Z6)
807 X2 = X2 - (CA(L) * W7 + V1 * Z7)
808 Y1 = Y1 - (CB(L) * W6 + V2 * Z6)
809 Y2 = Y2 - (CB(L) * W7 + V2 * Z7)
810 Z1 = Z1 + CG(L) * W6
811 Z2 = Z2 + CG(L) * W7
812 Next F5
813 Next I
814 Next K
815 H2 = -(X1 * T1 + Y1 * T2 + Z1 * T3) * G0
816 H1 = (X2 * T1 + Y2 * T2 + Z2 * T3) * G0
817 X4 = -(X1 * V1 + Y1 * V2) * G0
818 X3 = (X2 * V1 + Y2 * V2) * G0
819 If P$ = "D" Then GoTo 827
820 If RD = 0 Then GoTo 842
821 H1 = H1 / RD
822 H2 = H2 / RD
823 X3 = X3 / RD
824 X4 = X4 / RD
825 GoTo 842
826 Rem ----- PATTERN IN DB
827 P1 = -999
828 P2 = P1
829 P3 = P1
830 T1 = K9! * (H1 * H1 + H2 * H2)
831 T2 = K9! * (X3 * X3 + X4 * X4)
832 T3 = T1 + T2
833 Rem ----- CALCULATE VALUES IN DB
834 If T1 > 1E-30 Then P1 = 4.343 * Log(T1)
835 If T2 > 1E-30 Then P2 = 4.343 * Log(T2)
836 If T3 > 1E-30 Then P3 = 4.343 * Log(T3)
837 Print #3, Q2; Tab(15); Q1; Tab(29); P1; Tab(43); P2; Tab(57); P3
838 If S$ = "Y" Then Print #1, Q2; ","; Q1; ","; P1; ","; P2; ","; P3
839 GoTo 866
840 Rem ----- PATTERN IN VOLTS/METER
841 Rem ----- MAGNITUDE AND PHASE OF E(THETA)
842 S1 = 0
843 If (H1 = 0 And H2 = 0) Then GoTo 845
844 S1 = Sqr(H1 * H1 + H2 * H2)
845 If H1 <> 0 Then GoTo 848
846 S2 = 0
847 GoTo 851
848 S2 = Atn(H2 / H1) / P0
849 If H1 < 0 Then S2 = S2 + Sgn(H2) * 180
850 Rem ----- MAGNITUDE AND PHASE OF E(PHI)
851 S3 = 0
852 If (X3 = 0 And X4 = 0) Then GoTo 854
853 S3 = Sqr(X3 * X3 + X4 * X4)
854 If X3 <> 0 Then GoTo 857
855 S4 = 0
856 GoTo 859
857 S4 = Atn(X4 / X3) / P0
858 If X3 < 0 Then S4 = S4 + Sgn(X4) * 180
859 Print #3, USING; "###.##    "; Q2, Q1;
860 Print #3, USING; "       ##.###^^^^"; S1 * F1;
861 Print #3, USING; "   ###.##   "; S2;
862 Print #3, USING; "       ##.###^^^^"; S3 * F1;
863 Print #3, USING; "   ###.##"; S4
864 If S$ = "Y" Then Print #1, Q2; ","; Q1; ","; S1 * F1; ","; S2; ","; S3 * F1; ","; S4
865 Rem ----- INCREMENT ZENITH ANGLE
866 Q2 = Q2 + ZC
867 Next I2
868 Rem ----- INCREMENT AZIMUTH ANGLE
869 Q1 = Q1 + AC
870 Next I1
871 Close #1
872 Return
873 Rem ********** NEAR FIELD CALCULATION **********
874 Rem ----- ENSURE CURRENTS HAVE BEEN CALCULATED
875 If FLG < 2 Then GoSub 196
876 O2 = PWR
877 Print #3, " "
878 Print #3, B$; "    NEAR FIELDS     "; B$
879 Print #3, " "
880 INPUT "ELECTRIC OR MAGNETIC NEAR FIELDS (E/H) ";N$
881 If (N$ = "H" Or N$ = "E") Then GoTo 883
882 GoTo 880
883 Print
884 Rem ----- INPUT VARIABLES FOR NEAR FIELD CALCULATION
885 Print "FIELD LOCATION(S):"
886 A$ = "-COORDINATE (M): INITIAL,INCREMENT,NUMBER "
887 Print "   X"; A$;
888 INPUT XI,XC,NX
889 If NX = 0 Then NX = 1
890 If O$ > "C" Then Print #3, "X"; A$; ": "; XI; ","; XC; ","; NX
891 Print "   Y"; A$;
892 INPUT YI,YC,NY
893 If NY = 0 Then NY = 1
894 If O$ > "C" Then Print #3, "Y"; A$; ": "; YI; ","; YC; ","; NY
895 Print "   Z"; A$;
896 INPUT ZI,ZC,NZ
897 If NZ = 0 Then NZ = 1
898 If O$ > "C" Then Print #3, "Z"; A$; ": "; ZI; ","; ZC; ","; NZ
899 F1 = 1
900 Print
901 Print "PRESENT POWER LEVEL IS "; PWR; " WATTS"
902 INPUT "CHANGE POWER LEVEL (Y/N) ";A$
903 If A$ = "N" Then GoTo 908
904 If A$ <> "Y" Then GoTo 902
905 INPUT "NEW POWER LEVEL (WATTS)  ";O2
906 If O$ > "C" Then Print #3, " ": Print #3, "NEW POWER LEVEL (WATTS) = "; O2
907 GoTo 902
908 If (O2 < 0 Or O2 = 0) Then O2 = PWR
909 Rem ----- RATIO OF POWER LEVELS
910 F1 = Sqr(O2 / PWR)
911 If N$ = "H" Then F1 = F1 / S0 / 4 / P
912 Print
913 Rem ----- DESIGNATION OF OUTPUT FILE FOR NEAR FIELD DATA
914 INPUT "SAVE TO A FILE (Y/N) ";S$
915 If S$ = "N" Then GoTo 923
916 If S$ <> "Y" Then GoTo 914
917 INPUT "FILENAME (NAME.OUT)  ";F$
918 If Left$(Right$(F$, 4), 1) = "." Then GoTo 919 Else F$ = F$ + ".OUT"
919 If O$ > "C" Then Print #3, " ": Print #3, "FILENAME (NAME.OUT) "; F$
920 Open F$ For Output As #2
921 Print #2, NX * NY * NZ; ","; O2; ","; N$
922 Rem ----- LOOP OVER Z DIMENSION
923 For IZ = 1 To NZ
924 ZZ = ZI + (IZ - 1) * ZC
925 Rem ----- LOOP OVER Y DIMENSION
926 For IY = 1 To NY
927 YY = YI + (IY - 1) * YC
928 Rem ----- LOOP OVER X DIMENSION
929 For IX = 1 To NX
930 XX = XI + (IX - 1) * XC
931 Rem ----- NEAR FIELD HEADER
932 Print #3, " "
933 If N$ = "E" Then Print #3, B$; "NEAR ELECTRIC FIELDS"; B$
934 If N$ = "H" Then Print #3, B$; "NEAR MAGNETIC FIELDS"; B$
935 Print #3, Tab(10); "FIELD POINT: "; "X = "; XX; " Y = "; YY; " Z = "; ZZ
936 Print #3, "  VECTOR", "REAL", "IMAGINARY", "MAGNITUDE", "PHASE"
937 If N$ = "E" Then A$ = " V/M "
938 If N$ = "H" Then A$ = " AMPS/M "
939 Print #3, " COMPONENT  ", A$, A$, A$, " DEG"
940 A1 = 0
941 A3 = 0
942 A4 = 0
943 Rem ----- LOOP OVER THREE VECTOR COMPONENTS
944 For I = 1 To 3
945 X0 = XX
946 Y0 = YY
947 Z0 = ZZ
948 If N$ = "H" Then GoTo 958
949 T5 = 0
950 T6 = 0
951 T7 = 0
952 If I = 1 Then T5 = 2 * S0
953 If I = 2 Then T6 = 2 * S0
954 If I = 3 Then T7 = 2 * S0
955 U7 = 0
956 U8 = 0
957 GoTo 968
958 For J8 = 1 To 6
959 K!(J8, 1) = 0
960 K!(J8, 2) = 0
961 Next J8
962 J9 = 1
963 J8 = -1
964 If I = 1 Then X0 = XX + J8 * S0 / 2
965 If I = 2 Then Y0 = YY + J8 * S0 / 2
966 If I = 3 Then Z0 = ZZ + J8 * S0 / 2
967 Rem ----- LOOP OVER SOURCE SEGMENTS
968 For J = 1 To N
969 J1 = Abs(C%(J, 1))
970 J2 = Abs(C%(J, 2))
971 J3 = J2
972 If J1 > J2 Then J3 = J1
973 F4 = Sgn(C%(J, 1))
974 F5 = Sgn(C%(J, 2))
975 F6 = 1
976 F7 = 1
977 U5 = 0
978 U6 = 0
979 Rem ----- IMAGE LOOP
980 For K = 1 To G Step -2
981 If C%(J, 1) <> -C%(J, 2) Then GoTo 987
982 If K < 0 Then GoTo 1048
983 Rem ----- COMPUTE VECTOR POTENTIAL A
984 F6 = F4
985 F7 = F5
986 Rem ----- COMPUTE PSI(0,J,J+.5)
987 P1 = 0
988 P2 = 2 * J3 + J - 1
989 P3 = P2 + 0.5
990 P4 = J2
991 GoSub 75
992 U1 = T1 * F5
993 U2 = T2 * F5
994 Rem ----- COMPUTE PSI(0,J-.5,J)
995 P3 = P2
996 P2 = P2 - 0.5
997 P4 = J1
998 GoSub 66
999 V1 = F4 * T1
1000 V2 = F4 * T2
1001 Rem ----- REAL PART OF VECTOR POTENTIAL CONTRIBUTION
1002 X3 = U1 * CA(J2) + V1 * CA(J1)
1003 Y3 = U1 * CB(J2) + V1 * CB(J1)
1004 Z3 = (F7 * U1 * CG(J2) + F6 * V1 * CG(J1)) * K
1005 Rem ----- IMAGINARY PART OF VECTOR POTENTIAL CONTRIBUTION
1006 X5 = U2 * CA(J2) + V2 * CA(J1)
1007 Y5 = U2 * CB(J2) + V2 * CB(J1)
1008 Z5 = (F7 * U2 * CG(J2) + F6 * V2 * CG(J1)) * K
1009 Rem ----- MAGNETIC FIELD CALCULATION COMPLETED
1010 If N$ = "H" Then GoTo 1042
1011 D1 = (X3 * T5 + Y3 * T6 + Z3 * T7) * W2
1012 D2 = (X5 * T5 + Y5 * T6 + Z5 * T7) * W2
1013 Rem ----- COMPUTE PSI(.5,J,J+1)
1014 P1 = 0.5
1015 P2 = P3
1016 P3 = P3 + 1
1017 P4 = J2
1018 GoSub 56
1019 U1 = T1
1020 U2 = T2
1021 Rem ----- COMPUTE PSI(-.5,J,J+1)
1022 P1 = -P1
1023 GoSub 56
1024 U1 = (T1 - U1) / S(J2)
1025 U2 = (T2 - U2) / S(J2)
1026 Rem ----- COMPUTE PSI(.5,J-1,J)
1027 P1 = -P1
1028 P3 = P2
1029 P2 = P2 - 1
1030 P4 = J1
1031 GoSub 56
1032 U3 = T1
1033 U4 = T2
1034 Rem ----- COMPUTE PSI(-.5,J-1,J)
1035 P1 = -P1
1036 GoSub 56
1037 Rem ----- GRADIENT OF SCALAR POTENTIAL
1038 U5 = (U1 + (U3 - T1) / S(J1) + D1) * K + U5
1039 U6 = (U2 + (U4 - T2) / S(J1) + D2) * K + U6
1040 GoTo 1048
1041 Rem ----- COMPONENTS OF VECTOR POTENTIAL A
1042 K!(1, J9) = K!(1, J9) + (X3 * CR(J) - X5 * CI(J)) * K
1043 K!(2, J9) = K!(2, J9) + (X5 * CR(J) + X3 * CI(J)) * K
1044 K!(3, J9) = K!(3, J9) + (Y3 * CR(J) - Y5 * CI(J)) * K
1045 K!(4, J9) = K!(4, J9) + (Y5 * CR(J) + Y3 * CI(J)) * K
1046 K!(5, J9) = K!(5, J9) + (Z3 * CR(J) - Z5 * CI(J)) * K
1047 K!(6, J9) = K!(6, J9) + (Z5 * CR(J) + Z3 * CI(J)) * K
1048 Next K
1049 If N$ = "H" Then GoTo 1052
1050 U7 = U5 * CR(J) - U6 * CI(J) + U7
1051 U8 = U6 * CR(J) + U5 * CI(J) + U8
1052 Next J
1053 If N$ = "E" Then GoTo 1075
1054 Rem ----- DIFFERENCES OF VECTOR POTENTIAL A
1055 J8 = 1
1056 J9 = J9 + 1
1057 If J9 = 2 Then GoTo 964
1058 On I GoTo 1059, 1064, 1069
1059 H(3) = K!(5, 1) - K!(5, 2)
1060 H(4) = K!(6, 1) - K!(6, 2)
1061 H(5) = K!(3, 2) - K!(3, 1)
1062 H(6) = K!(4, 2) - K!(4, 1)
1063 GoTo 1097
1064 H(1) = K!(5, 2) - K!(5, 1)
1065 H(2) = K!(6, 2) - K!(6, 1)
1066 H(5) = H(5) - K!(1, 2) + K!(1, 1)
1067 H(6) = H(6) - K!(2, 2) + K!(2, 1)
1068 GoTo 1097
1069 H(1) = H(1) - K!(3, 2) + K!(3, 1)
1070 H(2) = H(2) - K!(4, 2) + K!(4, 1)
1071 H(3) = H(3) + K!(1, 2) - K!(1, 1)
1072 H(4) = H(4) + K!(2, 2) - K!(2, 1)
1073 GoTo 1097
1074 Rem ----- IMAGINARY PART OF ELECTRIC FIELD
1075 U7 = -M * U7 / S0
1076 Rem ----- REAL PART OF ELECTRIC FIELD
1077 U8 = M * U8 / S0
1078 Rem ----- MAGNITUDE AND PHASE CALCULATION
1079 S1 = 0
1080 If (U7 = 0 And U8 = 0) Then GoTo 1082
1081 S1 = Sqr(U7 * U7 + U8 * U8)
1082 S2 = 0
1083 If U8 <> 0 Then S2 = Atn(U7 / U8) / P0
1084 If U8 > 0 Then GoTo 1086
1085 S2 = S2 + Sgn(U7) * 180
1086 If I = 1 Then Print #3, "   X  ",
1087 If I = 2 Then Print #3, "   Y  ",
1088 If I = 3 Then Print #3, "   Z  ",
1089 Print #3, Tab(15); F1 * U8; Tab(29); F1 * U7; Tab(43); F1 * S1; Tab(57); S2
1090 If S$ = "Y" Then Print #2, F1 * U8; ","; F1 * U7; ","; F1 * S1; ","; S2
1091 Rem ----- CALCULATION FOR PEAK ELECTRIC FIELD
1092 S1 = S1 * S1
1093 S2 = S2 * P0
1094 A1 = A1 + S1 * Cos(2 * S2)
1095 A3 = A3 + S1 * Sin(2 * S2)
1096 A4 = A4 + S1
1097 Next I
1098 If N$ = "E" Then GoTo 1121
1099 Rem ----- MAGNETIC FIELD MAGNITUDE AND PHASE CALCULATION
1100 For I = 1 To 5 Step 2
1101 S1 = 0
1102 If (H(I) = 0 And H(I + 1) = 0) Then GoTo 1104
1103 S1 = Sqr(H(I) * H(I) + H(I + 1) * H(I + 1))
1104 S2 = 0
1105 If H(I) <> 0 Then S2 = Atn(H(I + 1) / H(I)) / P0
1106 If H(I) > 0 Then GoTo 1108
1107 S2 = S2 + Sgn(H(I + 1)) * 180
1108 If I = 1 Then Print #3, "   X  ",
1109 If I = 3 Then Print #3, "   Y  ",
1110 If I = 5 Then Print #3, "   Z  ",
1111 Print #3, Tab(15); F1 * H(I); Tab(29); F1 * H(I + 1); Tab(43); F1 * S1; Tab(57); S2
1112 If S$ = "Y" Then Print #2, F1 * H(I); ","; F1 * H(I + 1); ","; F1 * S1; ","; S2
1113 Rem ----- CALCULATION FOR PEAK MAGNETIC FIELD
1114 S1 = S1 * S1
1115 S2 = S2 * P0
1116 A1 = A1 + S1 * Cos(2 * S2)
1117 A3 = A3 + S1 * Sin(2 * S2)
1118 A4 = A4 + S1
1119 Next I
1120 Rem ----- PEAK FIELD CALCULATION
1121 PK = Sqr(A4 / 2 + Sqr(A1 * A1 + A3 * A3) / 2)
1122 Print #3, "   MAXIMUM OR PEAK FIELD = "; F1 * PK; A$
1123 If (S$ = "Y" And N$ = "E") Then Print #2, F1 * PK; ","; O2
1124 If (S$ = "Y" And N$ = "H") Then Print #2, F1 * PK; ","; O2
1125 If S$ = "Y" Then Print #2, XX; ","; YY; ","; ZZ
1126 Next IX
1127 Next IY
1128 Next IZ
1129 Close #2
1130 Return
1131 Rem ********** FREQUENCY INPUT **********
1132 Rem ----- SET FLAG
1133 Print
1134 INPUT "FREQUENCY (MHZ)";F
1135 If F = 0 Then F = 299.8
1136 If O$ > "C" Then Print #3, " ": Print #3, "FREQUENCY (MHZ):"; F
1137 W = 299.8 / F
1138 Rem -----VIRTUAL DIPOLE LENGTH FOR NEAR FIELD CALCULATION
1139 S0 = 0.001 * W
1140 Rem ----- 1 / (4 * PI * OMEGA * EPSILON)
1141 M = 4.77783352 * W
1142 Rem ----- SET SMALL RADIUS MODIFICATION CONDITION
1143 SRM = 0.0001 * W
1144 Print #3, "    WAVE LENGTH = "; W; " METERS"
1145 Rem ----- 2 PI / WAVELENGTH
1146 W = 2 * P / W
1147 W2 = W * W / 2
1148 FLG = 0
1149 Return
1150 Rem ********** GEOMETRY INPUT **********
1151 Rem ----- WHEN GEOMETRY IS CHANGED, ENVIRONMENT MUST BE CHECKED
1152 GoSub 1369
1153 Print
1154 If INFILE Then GoTo 1160
1155 INPUT "NO. OF WIRES";NW
1156 If NW = 0 Then Return
1157 If NW <= MW Then GoTo 1160
1158 Print "NUMBER OF WIRES EXCEEDS DIMENSION..."
1159 GoTo 1155
1160 If O$ > "C" Then Print #3, " ": Print #3, "NO. OF WIRES:"; NW
1161 Rem ----- INITIALIZE NUMBER OF PULSES TO ZERO
1162 N = 0
1163 For I = 1 To NW
1164 If INFILE Then GoSub 1557: GoTo 1190
1165 Print
1166 Print "WIRE NO."; I
1167 INPUT "   NO. OF SEGMENTS";S1
1168 If S1 = 0 Then GoTo 1153
1169 A$ = "   END ONE COORDINATES (X,Y,Z)"
1170 Print A$;
1171 INPUT X1,Y1,Z1
1172 If G < 0 And Z1 < 0 Then Print "Z CANNOT BE NEGATIVE": GoTo 1170
1173 A$ = "   END TWO COORDINATES (X,Y,Z)"
1174 Print A$;
1175 INPUT X2,Y2,Z2
1176 If G < 0 And Z2 < 0 Then Print "Z CANNOT BE NEGATIVE": GoTo 1174
1177 If X1 = X2 And Y1 = Y2 And Z1 = Z2 Then Print "ZERO LENGTH WIRE.": GoTo 1166
1178 A$ = "   RADIUS"
1179 Print "                     "; A$;
1180 INPUT A(I)
1181 If A(I) <= 0! Then GoTo 1179
1182 Rem ----- DETERMINE CONNECTIONS
1183 If O$ > "C" Then Print #3, " ": Print #3, "WIRE NO."; I
1184 GoSub 1299
1185 Print "CHANGE WIRE NO. "; I; " (Y/N) ";
1186 INPUT A$
1187 If A$ = "Y" Then GoTo 1165
1188 If A$ <> "N" Then GoTo 1185
1189 Rem ----- COMPUTE DIRECTION COSINES
1190 X3 = X2 - X1
1191 Y3 = Y2 - Y1
1192 Z3 = Z2 - Z1
1193 D = Sqr(X3 * X3 + Y3 * Y3 + Z3 * Z3)
1194 CA(I) = X3 / D
1195 CB(I) = Y3 / D
1196 CG(I) = Z3 / D
1197 S(I) = D / S1
1198 Rem ----- COMPUTE CONNECTIVITY DATA (PULSES N1 TO N)
1199 N1 = N + 1
1200 N(I, 1) = N1
1201 If (S1 = 1 And I1 = 0) Then N(I, 1) = 0
1202 N = N1 + S1
1203 If I1 = 0 Then N = N - 1
1204 If I2 = 0 Then N = N - 1
1205 If N > MP Then Print "PULSE NUMBER EXCEEDS DIMENSION": Close: GoTo 1155
1206 N(I, 2) = N
1207 If (S1 = 1 And I2 = 0) Then N(I, 2) = 0
1208 If N < N1 Then GoTo 1247
1209 For J = N1 To N
1210 C%(J, 1) = I
1211 C%(J, 2) = I
1212 W%(J) = I
1213 Next J
1214 C%(N1, 1) = I1
1215 C%(N, 2) = I2
1216 Rem ----- COMPUTE COORDINATES OF BREAK POINTS
1217 I1 = N1 + 2 * (I - 1)
1218 I3 = I1
1219 X(I1) = X1
1220 Y(I1) = Y1
1221 Z(I1) = Z1
1222 If C%(N1, 1) = 0 Then GoTo 1230
1223 I2 = Abs(C%(N1, 1))
1224 F3 = Sgn(C%(N1, 1)) * S(I2)
1225 X(I1) = X(I1) - F3 * CA(I2)
1226 Y(I1) = Y(I1) - F3 * CB(I2)
1227 If C%(N1, 1) = -I Then F3 = -F3
1228 Z(I1) = Z(I1) - F3 * CG(I2)
1229 I3 = I3 + 1
1230 I6 = N + 2 * I
1231 For I4 = I1 + 1 To I6
1232 J = I4 - I3
1233 X(I4) = X1 + J * X3 / S1
1234 Y(I4) = Y1 + J * Y3 / S1
1235 Z(I4) = Z1 + J * Z3 / S1
1236 Next I4
1237 If C%(N, 2) = 0 Then GoTo 1245
1238 I2 = Abs(C%(N, 2))
1239 F3 = Sgn(C%(N, 2)) * S(I2)
1240 I3 = I6 - 1
1241 X(I6) = X(I3) + F3 * CA(I2)
1242 Y(I6) = Y(I3) + F3 * CB(I2)
1243 If I = -C%(N, 2) Then F3 = -F3
1244 Z(I6) = Z(I3) + F3 * CG(I2)
1245 GoTo 1255
1246 Rem ---- SINGLE SEGMEN 0 PULSE CASE
1247 I1 = N1 + 2 * (I - 1)
1248 X(I1) = X1
1249 Y(I1) = Y1
1250 Z(I1) = Z1
1251 I1 = I1 + 1
1252 X(I1) = X2
1253 Y(I1) = Y2
1254 Z(I1) = Z2
1255 Next I
1256 Rem ********** GEOMETRY OUTPUT **********
1257 Print #3, " "
1258 Print #3, "                  **** ANTENNA GEOMETRY ****"
1259 If N > 0 Then GoTo 1264
1260 Print
1261 Print "NUMBER OF PULSES IS ZERO....RE-ENTER GEOMETRY"
1262 Print
1263 GoTo 1155
1264 K = 1
1265 J = 0
1266 For I = 1 To N
1267 I1 = 2 * W%(I) - 1 + I
1268 If K > NW Then GoTo 1279
1269 If K = J Then GoTo 1279
1270 J = K
1271 Print #3, " "
1272 Print #3, "WIRE NO. "; K; " COORDINATES", , , "CONNECTION PULSE"
1273 Print #3, "X", "Y", "Z", "RADIUS", "END1 END2  NO."
1274 If (N(K, 1) <> 0 Or N(K, 2) <> 0) Then GoTo 1279
1275 Print #3, "-", "-", "-", "    -", " -    -    0"
1276 K = K + 1
1277 If K > NW Then GoTo 1286
1278 GoTo 1270
1279 Print #3, X(I1); Tab(15); Y(I1); Tab(29); Z(I1); Tab(43); A(W%(I)); Tab(57);
1280 Print #3, USING; "###  ###   ##"; C%(I, 1), C%(I, 2), I
1281 If (I = N(K, 2) Or N(K, 1) = N(K, 2) Or C%(I, 2) = 0) Then K = K + 1
1282 If C%(I, 1) = 0 Then C%(I, 1) = W%(I)
1283 If C%(I, 2) = 0 Then C%(I, 2) = W%(I)
1284 If (K = NW And N(K, 1) = 0 And N(K, 2) = 0) Then GoTo 1270
1285 If (I = N And K < NW) Then GoTo 1270
1286 Next I
1287 Print
1288 Close 1: If INFILE Then INFILE = 0: If O$ > "C" Then GoTo 1293
1289 INPUT "    CHANGE GEOMETRY (Y/N) ";A$
1290 If A$ = "Y" Then GoTo 1153
1291 If A$ <> "N" Then GoTo 1289
1292 Rem ----- EXCITATION INPUT
1293 GoSub 1430
1294 Rem ----- LOADS/NETWORKS INPUT
1295 GoSub 1455
1296 FLG = 0
1297 Return
1298 Rem ********** CONNECTIONS **********
1299 E(I) = X1
1300 L(I) = Y1
1301 M(I) = Z1
1302 E(I + NW) = X2
1303 L(I + NW) = Y2
1304 M(I + NW) = Z2
1305 G% = 0
1306 I1 = 0
1307 I2 = 0
1308 J1(I) = 0
1309 J2(I, 1) = -I
1310 J2(I, 2) = -I
1311 If G = 1 Then GoTo 1323
1312 Rem ----- CHECK FOR GROUND CONNECTION
1313 If Z1 = 0 Then GoTo 1315
1314 GoTo 1318
1315 I1 = -I
1316 J1(I) = -1
1317 GoTo 1340
1318 If Z2 = 0 Then GoTo 1320
1319 GoTo 1323
1320 I2 = -I
1321 J1(I) = 1
1322 G% = 1
1323 If I = 1 Then GoTo 1358
1324 For J = 1 To I - 1
1325 Rem ----- CHECK FOR END1 TO END1
1326 If (X1 = E(J) And Y1 = L(J) And Z1 = M(J)) Then GoTo 1328
1327 GoTo 1333
1328 I1 = -J
1329 J2(I, 1) = J
1330 If J2(J, 1) = -J Then J2(J, 1) = J
1331 GoTo 1340
1332 Rem ----- CHECK FOR END1 TO END2
1333 If (X1 = E(J + NW) And Y1 = L(J + NW) And Z1 = M(J + NW)) Then GoTo 1335
1334 GoTo 1339
1335 I1 = J
1336 J2(I, 1) = J
1337 If J2(J, 2) = -J Then J2(J, 2) = J
1338 GoTo 1340
1339 Next J
1340 If G% = 1 Then GoTo 1358
1341 If I = 1 Then GoTo 1358
1342 For J = 1 To I - 1
1343 Rem ----- CHECK END2 TO END2
1344 If (X2 = E(J + NW) And Y2 = L(J + NW) And Z2 = M(J + NW)) Then GoTo 1346
1345 GoTo 1351
1346 I2 = -J
1347 J2(I, 2) = J
1348 If J2(J, 2) = -J Then J2(J, 2) = J
1349 GoTo 1358
1350 Rem ----- CHECK FOR END2 TO END1
1351 If (X2 = E(J) And Y2 = L(J) And Z2 = M(J)) Then GoTo 1353
1352 GoTo 1357
1353 I2 = J
1354 J2(I, 2) = J
1355 If J2(J, 1) = -J Then J2(J, 1) = J
1356 GoTo 1358
1357 Next J
1358 Print #3, "            COORDINATES", "  ", "  ", "END         NO. OF"
1359 Print #3, "   X", "   Y", "   Z", "RADIUS     CONNECTION     SEGMENTS"
1360 Print #3, X1; Tab(15); Y1; Tab(29); Z1; Tab(57); I1
1361 Print #3, X2; Tab(15); Y2; Tab(29); Z2; Tab(43); A(I); Tab(57); I2; Tab(71); S1
1362 Return
1363 Rem ********** ENVIROMENT INPUT **********
1364 Print
1365 Print "                        **** WARNING ****"
1366 Print "REDO GEOMETRY TO ENSURE PROPER GROUND CONNECTION/DISCONNECTION"
1367 Print
1368 Rem ----- INITIALIZE NUMBER OF RADIAL WIRES TO ZERO
1369 NR = 0
1370 Rem ----- SET ENVIRONMENT
1371 Print #3, " "
1372 A$ = "ENVIRONMENT (+1 FOR FREE SPACE, -1 FOR GROUND PLANE)"
1373 Print A$;
1374 INPUT G
1375 If O$ > "C" Then Print #3, A$; ": "; G
1376 If G = 1 Then GoTo 1428
1377 If G <> -1 Then GoTo 1373
1378 Rem ----- NUMBER OF MEDIA
1379 A$ = " NUMBER OF MEDIA (0 FOR PERFECTLY CONDUCTING GROUND)"
1380 Print A$;
1381 INPUT NM
1382 If NM <= MM Then GoTo 1385
1383 Print "NUMBER OF MEDIA EXCEEDS DIMENSION..."
1384 GoTo 1380
1385 If O$ > "C" Then Print #3, A$; ": "; NM
1386 Rem ----- INITIALIZE BOUNDARY TYPE
1387 TB = 1
1388 If NM = 0 Then GoTo 1428
1389 If NM = 1 Then GoTo 1396
1390 Rem ----- TYPE OF BOUNDARY
1391 A$ = " TYPE OF BOUNDARY (1-LINEAR, 2-CIRCULAR)"
1392 Print "            "; A$;
1393 INPUT TB
1394 If O$ > "C" Then Print #3, A$; ": "; TB
1395 Rem ----- BOUNDARY CONDITIONS
1396 For I = 1 To NM
1397 Print "MEDIA"; I
1398 A$ = " RELATIVE DIELECTRIC CONSTANT, CONDUCTIVITY"
1399 Print "         "; A$;
1400 INPUT T(I),V(I)
1401 If O$ > "C" Then Print #3, A$; ": "; T(I); ","; V(I)
1402 If I > 1 Then GoTo 1414
1403 If TB = 1 Then GoTo 1414
1404 A$ = " NUMBER OF RADIAL WIRES IN GROUND SCREEN"
1405 Print "            "; A$;
1406 INPUT NR
1407 If O$ > "C" Then Print #3, A$; ": "; NR
1408 If NR = 0 Then GoTo 1414
1409 A$ = " RADIUS OF RADIAL WIRES"
1410 Print "                             "; A$;
1411 INPUT RR
1412 If O$ > "C" Then Print #3, A$; ": "; RR
1413 Rem ----- INITIALIZE COORDINATE OF MEDIA INTERFACE
1414 U(I) = 1000000!
1415 Rem ----- INITIALIZE HEIGHT OF MEDIA
1416 H(I) = 0
1417 If I = NM Then GoTo 1422
1418 A$ = " X OR R COORDINATE OF NEXT MEDIA INTERFACE"
1419 Print "          "; A$;
1420 INPUT U(I)
1421 If O$ > "C" Then Print #3, A$; ": "; U(I)
1422 If I = 1 Then GoTo 1427
1423 A$ = " HEIGHT OF MEDIA"
1424 Print "                                    "; A$;
1425 INPUT H(I)
1426 If O$ > "C" Then Print #3, A$; ": "; H(I)
1427 Next I
1428 Return
1429 Rem ********** EXCITATION INPUT **********
1430 Print
1431 A$ = "NO. OF SOURCES "
1432 Print A$;
1433 INPUT NS
1434 If NS < 1 Then NS = 1
1435 If NS <= MP Then GoTo 1438
1436 Print "NO. OF SOURCES EXCEEDS DIMENSION ..."
1437 GoTo 1432
1438 If O$ > "C" Then Print #3, " ": Print #3, A$; ": "; NS
1439 For I = 1 To NS
1440 Print
1441 Print "SOURCE NO. "; I; ":"
1442 A$ = "PULSE NO., VOLTAGE MAGNITUDE, PHASE (DEGREES)"
1443 Print A$;
1444 INPUT E(I),VM,VP
1445 If E(I) <= N Then GoTo 1448
1446 Print "PULSE NUMBER EXCEEDS NUMBER OF PULSES..."
1447 GoTo 1443
1448 If O$ > "C" Then Print #3, A$; ": "; E(I); ","; VM; ","; VP
1449 L(I) = VM * Cos(VP * P0)
1450 M(I) = VM * Sin(VP * P0)
1451 Next I
1452 If FLG = 2 Then FLG = 1
1453 Return
1454 Rem ********** LOADS INPUT **********
1455 Print
1456 INPUT "NUMBER OF LOADS       ";NL
1457 If NL <= ML Then GoTo 1460
1458 Print "NUMBER OF LOADS EXCEEDS DIMENSION..."
1459 GoTo 1456
1460 If O$ > "C" Then Print #3, "NUMBER OF LOADS"; NL
1461 If NL < 1 Then GoTo 1492
1462 INPUT "S-PARAMETER (S=jw) IMPEDANCE LOAD (Y/N)";L$
1463 If L$ <> "Y" And L$ <> "N" Then GoTo 1462
1464 A$ = "PULSE NO.,RESISTANCE,REACTANCE"
1465 If L$ = "Y" Then A$ = "PULSE NO., ORDER OF S-PARAMETER FUNCTION"
1466 For I = 1 To NL
1467 Print
1468 Print "LOAD NO. "; I; ":"
1469 If L$ = "Y" Then GoTo 1476
1470 Print A$;
1471 INPUT LP(I),LA(1,I,1),LA(2,I,1)
1472 If LP(I) > N Then Print "PULSE NUMBER EXCEEDS NUMBER OF PULSES...": GoTo 1470
1473 If O$ > "C" Then Print #3, A$; ": "; LP(I); ","; LA(1, I, 1); ","; LA(2, I, 1)
1474 GoTo 1491
1475 Rem ----- S-PARAMETER LOADS
1476 Print A$;
1477 INPUT LP(I),LS(I)
1478 If LP(I) > N Then Print "PULSE NUMBER EXCEEDS NUMBER OF PULSES...": GoTo 1476
1479 If LS(I) > MA Then Print "MAXIMUM DIMENSION IS 10": GoTo 1477
1480 If O$ > "C" Then Print #3, A$; ": "; LP(I); ","; LS(I)
1481 For J = 0 To LS(I)
1482 A$ = "NUMERATOR, DENOMINATOR COEFFICIENTS OF S^"
1483 Print A$; J;
1484 INPUT LA(1,I,J),LA(2,I,J)
1485 If O$ > "C" Then Print #3, A$; J; ":"; LA(1, I, J); ","; LA(2, I, J)
1486 Next J
1487 If LS(I) > 0 Then GoTo 1491
1488 LS(I) = 1
1489 LA(1, I, 1) = 0
1490 LA(2, I, 1) = 0
1491 Next I
1492 FLG = 0
1493 Return
1494 Rem ********** MAIN PROGRAM **********
1495 Rem ----- DATA INITIALIZATION
1496 Rem ----- PI
1497 P = 4 * Atn(1)
1498 Rem ----- CHANGES DEGREES TO RADIANS
1499 P0 = P / 180
1500 B$ = "********************"
1501 Rem ----- INTRINSIC IMPEDANCE OF FREE SPACE DIVIDED BY 2 PI
1502 G0 = 29.979221
1503 Rem ---------- Q-VECTOR FOR GAUSSIAN QUADRATURE
1504 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)
1505 READ Q(13), Q(14)
1506 Data 0.288675135, 0.5, 0.430568156, 0.173927423, 0.169990522, 0.326072577
1507 Data 0.480144928, 0.050614268, 0.398333239, 0.111190517
1508 Data 0.262766205, 0.156853323, 0.091717321, 0.181341892
1509 Rem ---------- E-VECTOR FOR COEFFICIENTS OF ELLIPTIC INTEGRAL
1510 READ C0, C1, C2, C3, C4, C5, C6, C7, C8, C9
1511 Data 1.38629436112, 0.09666344259, 0.03590092383, 0.03742563713, 0.01451196212
1512 Data 0.5, 0.12498593397, 0.06880248576, 0.0332835346, 0.00441787012
1513 Rem ----- IDENTIFY OUTPUT DEVICE
1514 GoSub 1580
1515 Print #3, Tab(20); B$; B$
1516 Print #3, Tab(22); "MINI-NUMERICAL ELECTROMAGNETICS CODE"
1517 Print #3, Tab(36); "MININEC"
1518 Print #3, Tab(24); Date$; Tab(48); Time$
1519 Print #3, Tab(20); B$; B$
1520 Rem ----- FREQUENCY INPUT
1521 GoSub 1133
1522 Rem ----- ENVIRONMENT INPUT
1523 GoSub 1369
1524 Rem ----- CHECK FOR NEC-TYPE GEOMETRY INPUT
1525 GoSub 1550
1526 Rem ----- GEOMETRY INPUT
1527 GoSub 1153
1528 Rem ----- MENU
1529 Print
1530 Print B$; "    MININEC MENU    "; B$
1531 Print "   G - CHANGE GEOMETRY     C - COMPUTE/DISPLAY CURRENTS"
1532 Print "   E - CHANGE ENVIRONMENT  P - COMPUTE FAR-FIELD PATTERNS"
1533 Print "   X - CHANGE EXCITATION   N - COMPUTE NEAR-FIELDS"
1534 Print "   L - CHANGE LOADS"
1535 Print "   F - CHANGE FREQUENCY    Q - QUIT"
1536 Print B$; B$; B$
1537 INPUT "   COMMAND ";C$
1538 If C$ = "F" Then GoSub 1133
1539 If C$ = "P" Then GoSub 621
1540 If C$ = "X" Then GoSub 1430
1541 If C$ = "E" Then GoSub 1364
1542 If C$ = "G" Then GoSub 1152
1543 If C$ = "C" Then GoSub 497
1544 If C$ = "L" Then GoSub 1455
1545 If C$ = "N" Then GoSub 875
1546 If C$ <> "Q" Then GoTo 1529
1547 If O$ = "P" Then Print #3, Chr$(12) Else If O$ = "C" Then Print #3, " "
1548 Close
1549 GoTo 1617
1550 Rem ********** NEC-TYPE GEOMETRY INPUT **********
1551 Open "MININEC.INP" For Random As #1 Len = 30
1552 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$
1553 GET 1
1554 NW = CVI(S$)
1555 If NW Then INFILE = 1
1556 Return
1557 Rem ---------- GET GEOMETRY DATA FROM MININEC.INP
1558 GET 1
1559 S1 = CVI(S$)
1560 X1 = CVS(X1$)
1561 Y1 = CVS(Y1$)
1562 Z1 = CVS(Z1$)
1563 X2 = CVS(X2$)
1564 Y2 = CVS(Y2$)
1565 Z2 = CVS(Z2$)
1566 A(I) = CVS(R$)
1567 If G < 0 Then If Z1 < 0 Or Z2 < 0 Then GoSub 1572
1568 Print #3, " ": Print #3, "WIRE NO."; I
1569 If X1 = X2 And Y1 = Y2 And Z1 = Z2 Then Print "WIRE LENGTH IS ZERO.": GoTo 1547
1570 GoSub 1299
1571 Return
1572 If IZNEG Then GoTo 1576
1573 Print "NEGATIVE Z VALUE ENCOUNTERED FOR GROUND PLANE."
1574 INPUT "ABORT OR CONVERT NEGATIVE Z VALUE TO ZERO (A/C)? ";A$
1575 If A$ = "A" Then GoTo 1547 Else If A$ = "C" Then IZNEG = 1 Else GoTo 1574
1576 If Z1 < 0 Then Z1 = -Z1
1577 If Z2 < 0 Then Z2 = -Z2
1578 Return
1579 Rem ********** IDENTIFY OUTPUT DEVICE **********
1580 INPUT "OUTPUT TO CONSOLE, PRINTER, OR DISK (C/P/D)";O$
1581 If O$ = "C" Then F$ = "SCRN:": GoTo 1586
1582 If O$ = "P" Then F$ = "LPT1:": GoTo 1586
1583 If O$ <> "D" Then GoTo 1580
1584 INPUT "FILENAME (NAME.OUT)";F$
1585 If Left$(Right$(F$, 4), 1) = "." Then GoTo 1586 Else F$ = F$ + ".OUT"
1586 Open F$ For Output As #3
1587 Cls
1588 Return
1589 Rem ********** CALCULATE ELAPSED TIME **********
1590 IH = Val(Mid$(T$, 1, 2)) - Val(Mid$(OT$, 1, 2))
1591 IM = Val(Mid$(T$, 4, 2)) - Val(Mid$(OT$, 4, 2))
1592 IS=VAL(MID$(T$,7,2))-VAL(MID$(OT$,7,2))
1593 IF IS<0 THEN IS=IS+60:IM=IM-1
1594 If IM < 0 Then IM = IM + 60: IH = IH - 1
1595 If IH < 0 Then IH = IH + 24
1596 T$=":"+MID$(STR$(IS+100),3)
1597 If IH Then T$ = Mid$(Str$(IH), 2) + ":" + Mid$(Str$(IM + 100), 3) + T$ Else T$ = Mid$(Str$(IM), 2) + T$
1598 Return
1599 Rem ********** CALCULATE APPROXIMATE TIME REMAINING **********
1600 IPCT = 100 * PCT
1601 T$ = Time$
1602 IH = Val(Mid$(T$, 1, 2)) - Val(Mid$(OT$, 1, 2))
1603 If IH < 0 Then IH = IH + 24
1604 IM = Val(Mid$(T$, 4, 2)) - Val(Mid$(OT$, 4, 2))
1605 IS=VAL(MID$(T$,7,2))-VAL(MID$(OT$,7,2))
1606 IS=IS+60*(IM+60*IH)
1607 IS=IS*(1/PCT-1)
1608 IM=INT(IS/60)
1609 IS=IS MOD 60
1610 IH = Int(IM / 60)
1611 IM = IM Mod 60
1612 T$=":"+MID$(STR$(IS+100),3)
1613 If IH Then T$ = Mid$(Str$(IH), 2) + ":" + Mid$(Str$(IM + 100), 3) + T$ Else T$ = Mid$(Str$(IM), 2) + T$
1614 LOCATE CSRLIN, 1
1615 Print Q$; IPCT; "% COMPLETE - APPROX TIME REMAINING "; T$; "   ";
1616 Return
1617 End