02142-503-M 2142-503-M = RL-2043 THE UNIVERSITY OF MICHIGAN COLLEGE OF HtN G DEPARTMENT OF ILICTRICAL ENGINEERN4G Radiation LabmiAtr y ANALYSES AND DIGITAL PROGRAMS FOR COMPUTING THE MONOPULSE POINTING ERRORS ASSOCIATED WITH THREE-DIMENSIONAL CONICAL AND OGIVAL RADOMES WITH OR WITHOUT A SURROUNDING (WEAK) PLASMA By T. B. A. Senior and V. V. Lepa February 1970 Prepared for Harry mond Laboratories Microwave Branch, 250 Attn: Dr. Irvin Pollin, AMXDO-RDB Washington, DC 20438 Ann Arbor, Michigan, 48108 201 East Catherine Street

THE UNIVERSITY OF MICHIGAN TABLE OF CONTENTS I INTRODUCTION 1 II RAY TRACING PRESCRIPTION 3 m INCIDENT FIELD SAMPLING 18 IV MONOPULSE PLATE AND POINTING ERROR 25 V PLASMA MODIFICATION 36 VI COMPUTER PROGRAM AND NUMERICAL DISCUSSION 43 VH CONCLUSIONS 2 APPENDIX: COMPUTER PROGRAM PRINT-OUTS

THE UNIVERSITY OF MICHIGAN I. INTRODUCTION In a preceding memorandum (No. 02142-502-M) entitled "The Monopulse Pointing Error Associated With a Two-Dimensional Conical or Ogival Radome With or Without a Surrounding (Weak) Plasma", hereinafter referred to as S-L, we described the motivation and bckground of the present study, and discussed the procedures which are appropriate to the determination of the pointing error in the two-dimensional case. In addition the analytical steps were described, and a summary and print-out of the computer program were prsented. It was felt desirable to study the two-dimensional case in some detail before attempting the more general and complex problem posed by a threedimensional geometry. This was undoubtedly a wise decision, and the knowledge gained from a consideration of a problem whose geometry is sufficiently simple to permit a direct interpretation of the results has nabled us to arriw at a threedimensional problem which is more efficient than might otherwise have been the case. There is, moreover a certain degree of simplicity between the analyses in the two cases, and a considerable similarity in the concepts. The present memorandum is devoted entirely to the three-dimensional problem in which a plane wave of arbitrary polarization is incident at an arbitrary direction on a rotationally-symmetric radome whose surfaces are described by either linear or quadratic equations of stated form, corresponding to conical and ogival radomes respectively. Wherever appropriate, reference is made to the previous memorandum, S-L, for details which are common to the analyses for the two- and three-dimensional problems. There are, however, numerous features which are unique to the latter case, and these are described in Sections II through V. A discussion of the numerical procedures, including the choice of tolerances, ray sampling areas, etc, is given in Section VI, and the Appendix contains a listing of the computer program. It should be noted that there are actually two 1

THE UNIVERSITY OF MICHIGAN separate programs, one pertaining to the conical radome and the other to the ogival one. This separation was found convenient because of the need for numerical iteration to locate the ray intercepts when the radome is ogival. In the two-dimensional case it was possible to run the program for a wide variety of incidence angles and or a number of radome/plasma/polarization ombinations, and to examin the nau of th rsulting pointig errors. Because of muh greatr nth (dhe, ai ti of the threedimeonal program, such detailed checkn has not ben possible in the preet case. Indeed, the funds that were available have allowed only those tests that were required to de-bug the program and to give confidence that it will work in any given case, together with a few runs of the complete program for single Incidence and polarization angles. The conditions under which these particular runs were carried out, and the results that were obtained, are described in Section VI. It is a pleasure to acknowledge the assistance of Dr. K-H Hsu with the programming. 2

THE UNIVERSITY OF MICHIGAN II RAY TRACING PRESCRIPTIONS The radomes considered here are assumed to be of single layer construction and formed from a material which can be treated as a homogeneous, isotropic and lossless dielectric whose permeability is the same as that of free space. The electromagnetic properties can therefore be represented by a real refractive index n, a typical value of which is 2.5 appropriate to fiberglas in the GHz range of frequencies. Two particular radome configurations are treated, namely, ogival and conical. In either instane it is assumed that outer and inner surfaces are both ogival or both conical, and formed by the revolution of a straight line about some axis, or by the revolution of an arc of a circle about a chord. Although the choice does permit the radome thickness to be non-uniform, it will be appreciated that the type of thickness variation that is encompassed is distinctly limited. The radome surface themselves are described in terms of a rectangular Cartesian coordinate system (x, y, z) whose z axis is the axis of the radome. For computational purposes it is convenient to choose the origin of coordinates at a small but non-zero distance I to the left of the radome tip (see Fig. 1). In the ogival case, the outer surface can then be defined as Pouter /A2+B-( 2 -A (1) o ter valid for z > C-B, where p -Yx gy; A, B and C are positive real numerical constants in terms of which the maximum diameter of the radome is 2 (' A2+B2 -A) occurring Mt x = C, the (longtudnal) raus of curv re is 3

THE UNIVERSITY OF MICHIGAN R "A+B2 ad I = C-B. The overall length of the radome is 2B, but in practice interest is confined to that portion of the interior extending to at most the position of the maximum diameter, i.e. to z < C. The inner radome surface is defined in a similar manner. For the conical radome the definitions of the radome surfaces are more straightforward, and for the outer surface we have p a(z-t ) (2) valid for z> I, with a similar definition for the inner surface. The electromagnetic field incident on the radome is assumed to be a uniform plane wave of arbitrary polarization propagating in an arbitrary direction. It is convenient to write the electric vector as E= ( m )ik(e x y +nz)(3) where (1, m, n) and (l,mln) are two sets of direction cosines measured with respect to (x, y, z), and such that (l, m, n) (11' Iin ni), IlI +ml+nnl = 0 (4) to make E perpendicular to the direction of propagation, as required. We note that the magnetic vector corresponding to E is H = Y(2' m2 n2)ek(lx + my nz) (5) where Y is the intrinsic admittance of free space. Since (12, m 2,n2) = (1, m, n)^ (11, ml n) (5) Of course, t2+m2n = 1 1= n 2 2 2 4

THE UNIVERSITY OF MICHIGAN aknowledge of 5 of the 6 drection cosines l,m,n, I,ml, nserves to specify the incident field in its entirety. In practice, however, we shall regard all 6 as input parameters chosen in accordance with Eq. (4). We now turn to a consideration of ray tracing. A general discussion of the philosophy and limitations of ray tracing as applied to this type of problem was given in S-L, and will not be enlarged upon here. Instead, we shall confine ourselves merely to a description of the 'mechanics' of carrying through this procedure for a three-dimensional geometry of the sort shown in Fig. 1. The incident field (3) was taken to have zero phase wavefront passing through the origin of the coordinate system. This is the wavefront (or phase plane) which will be sampled. Choose x, y (0) such that the point x x(0) 1 (0) (0) y(), z=z() lies in this plane. Then z) - - (Ix0)~+ my) ) (7) The ray through this point is perpendicular to the phase front and has equations x(~)+ (z ) yy()+ m (z.(0) ) (8) This will intercept (say) the lower surface of the outer radome. The z coordinate, (1) z=z of the point of interception is given by the equation p [{x(~+ (z-z(~) + {y )+ m (z() ) p(z) (9) where p=p(z) is the equation of the radome surface. For a conical radome (see Eq. (2) ), the determination of z(1) from (9) requires only the solution of a quadratic equation, but for an ogival radome (see Eq. (1) ) it is necessary to resort to iteration. Having found z(l, the remaning coordinates x(1) and y(l) of the intercept point follow from (8). The distance 5

THE UNIVERSITY OF MICHIGAN d -= {(1).x))2+(1).L 0) )2+(z(l)-() )2 }1/2 (10) is then computed and stored. Compute the direction cosines m ) 1 (X. e. ) (11) n ul+(8p/8 p 1+ (a/Oz)2 of the outward normal to the radome at the point of ntercept. This vector and the ray direction establish the plane of incidence at the point, and because the manner of transmission and reflection at the Interface depends on the inclination of the electric vector to the plane, we now perform the first of an almost endless series of rotations of coordinates. Choose (x1, Y1' z) such that nm x (XYlZ) 1 ml n, (12) 12 m2 The incident field then takes the form ikx fkXt E yle, H= zle, (13) and the ray that strikes the radome is traveling in the positive x direction. Snce A A- A I,A ^A xI = x+my +nm and the normal direction I n +m nz spelify the plane of incidence, a rotation of these new coordinates about the x1 axis can produce a coordinate system x2, y, 2 z with x2x and z2 normal to this plane. But this normal has direction cosines +(In m,,n,(, m,n)t +(m n-n m, n l- n, 1 m-m I), n n n n n n n n n 6

THE UNIVERSITY OF MICHIGAN and hence we can take (2' Y 2) = 1 0 0 A 0 -B 0 x B Y1 A Z1 (14) with A ( n 1 n1) {(n 1 mnml+nnnl) +(1n 2+m m2 nn)2 B =y(l 2+m mn n mn 1+ mn n1 2+m nm2m+n n2 -1/2 n 2 n 2 n 2 n 1 n i n n2 n 2 n2 (15) and y7 = 1. Although either sign for 7 would be acoeptable, it is necessary to choose one and to know the implications of that choice in order that we may later deduce angles correctly from a knowledge of their sines or cosines. To remove this arbitrariness, we impose the requirement that Y2 ( 'm n'n)> 0 (16) 2 n n ni.e. that Y2 makes an acute angle with the outward normal. t is a trivial matter to show that (16) implies y= 1. (17) These new coordinates are related to the original coordinates x, y, z via the matrix equation x2 =M 2 1 x (18) where. 4 A I 1 m n ' L D1i (I1 m nl) ' 2 = 0 12 m2 n2 o and A and B are given by Eq. (15) and (17), and the takes the form 0 A -B B A/ incident electric vector now 7

THE UNIVERSITY OF MICHIGAN E = (AY2-B2)e. (19( Recalling that Y2 lies in the plane of incidence, A and -B therefore represent the electric field components in and perpendioular to the plane of incidence, and we can invoke the Fresnel reflection and transmission formulae to determine the strengths of the reflected and transmitted flekb. Consider the elementary problem illustrated in Fig. 2. Medium 1 is free space and medium 2 is a pure dielectric with refractive index 2 relative to free space. A plane wave is incident at an angle a to the normal to the interfacoe, and from S&ell's law, the direction which the transmitted wave makes with the normal is then 3, where sin a (90) ain t= a-,. (20) n If the incident field is polarized with its electric vector perpendicular to the plane of incidence, the reflection and transmission coefficients are respectively 1 1-r L 2r1 R s= -am, T1 -— =- (21) 12 1+r 1 11+r (see Section LI of S-L) with r...; (22) 3cosI3 whereas for the electric vector in the plate of Incidenoe 12 i +r i1' 2 3 with noso =a n 2r (24) The above coefficients are, of course, amplitude factors only, relating to the electric field, and in using these formulae, particularly those for the parallel 8

THE UNIVERSITY OF MICHIGAN case, it is essential that the orientation of the electric vectors in the reflected and transmitted waves be taken as shown in Fig. 2. We also note in passing that if the plane wave were to have originated in medium 2 and to have struck the interface with medium 1 (fre space) at an angle P to the normal, the reflection and transmission coefficients corresponding to the above would be l-r' r 2 R' T (25) 21 1F'' 2= 1 +' (25) R 1-.1 2T (26) 21 '- r ' T21 1i+r' where r' and r are as given in (22) and (24) respectively. Reverting now to the ray tracing procedure through the radome, we observe that the angle a appearing in the above formulae is that which the ray direction makes with the inward normal at the point of intercept, and hence cos a = - x2-(nx M y +n ) - (n +mm ann ) (27) n n n with in a Y2' (1 +m y 4n A). n n n 1 (28) Equations (27) and (28) specify a uniquely, and using (20), the angle B oan then be found. The reflected and transmitted rays both lie n the plane of incidence, and though a ray which is reflected at its initial impact with the radome is of no real ooncern to us, it is convenient to describe the tracing procedure for such a ray as well uince the analysis is almost directly applicable to reflection and transmission at any surface of the radome. 9

THE UNIVERSITY OF MICHIGAN To describe the reflected ray, we introduce new coordinates (x3, Y3 z3) obtained from x2, Y22 z) by rotation about the z2 axis and such that 23 is in the direction of the reflected ray (see Fig. 3a). Clearly with x 3 = 3.a -oos2a 3 = -sin2a 3 \ 0 (29) (30) osi 0 0 ) implying '3 M3 M2 M1 (31) nd hence from Eq. (19), (21) and (23), the electric vector associated with the reflected ray is kx3) E (ARi -BR )e -12 3 12 3 (32) The final step is to rotate the coordinates yet again so as to cast this result into a form as compact as possessed by the electric vector of the incident ray in the (x. Y1, zl) coordinates. This can be achieved by a rotation about the x3 axis to produce coordinates (x4,y4, z4): where with x =S0 40: (33) 0 ARI 74 BRI2 12 -. w4 0 -L AR12 74 )21/ (q4>0). (34) (35) Y4 {(AR" )2 (B 74 R 12 + (BR12 10

THE UNIVERSITY OF MICHIGAN implying M 3 M2 x - <36) in which cae the electric vector becomes ikx4 E Y4 y4 e. (37) The procedure for the transmitted ray is very similar. By rotating the ooordinate system about the z2 axis, a new set of coordinates (x3 y3, z8 ) can be established having x3 in the direction of the transmitted ray (see Fig. 3b). We have X3 2 5 (38) with /oos(c-) -sin(-3) 0\ ( sn(,i-) cos(AB) 0 (39) 0 0 1/ implying Al M2 MM (40) '4 ~3=2=1 -The electric vector associated with the transmitted ray is then E=(AT12tA-BT I )e 3 (41) -- s 12 Y3 BT~2 t3 )e (43) and if we now iftts the coordinates again to produce (x, y,, z' where M4 14X (42) with 11

THE UNIVERSITY OF MICHIGAN 1 0 0 AT BTL 0 AT12 12 and 4 {(ATI)2 (BT2)2 (74> ), (44) Implying M M4 t3 Mx (45) the electric vector in the transmitted ray becomes E A e i k (48) This restores the ray to the form it had prior to its Intercept with the radome. Although the above discussion has been phrased in terms of the initial impact of a ray with the outer surface of the radome, the procedure is that which must be applied at any interface. To illustrate this fact, consider the ray (45) transmitted through the outer surface. Under almost all circumstances, this ray will sooner or later strike the inner surface, and knowing the ray direction, XI, in terms of the coordinates (x,y, z), together with a point on the ray (namely, inner surface can be found in the same way as before. The distance = { 2)-x( ) )2+(y2)-y(l) )2+(z(2).z(l) )2 ) 1/2 (47) is the co ted n computed, and n d ed to d to give the accumulated optical istance from the referene pase plane. From now on te a c e the aal through in the 12

THE UNIVERSITY OF MICHIGAN manner already described apart from such trivial changes as are caused by the fact that the ray is incident in the denser medium. Specifically, in Eqs. (30) and (39) a must be replaced by v, and vice versa, and the reflection and transmisbion coefficients given by (21) and (23) must be replaced by those in (25) ad (26). These changes apart, the same analysis is applicable. The above procedure enables us to traoe a ray from the outside of the radome through the layer and Into the interior space, keeping track of the direction of the ray and of the electric vector, as well as of the amplitude and phase (or, in effect, the optical distanoe) associated with the ray. At the inner surface of the radome a reflected ray is, of course, produced, and this is directed back to the outer surface where transmission and reflection again occur. Only the reflected wave is of Interest, and we remark that the analysis appropriate to this is identical to that corresponding to reflection at the inner radome surface providing that the normal direction now used is the inward-directed one. hb the analysis, therefore, the direction cosines In, m, n of the outward normal must be reversed in sign. This reflected ray can now be traced back to the inner radome surface; and so on. And if, at some stage, a ray crosses the interior of the radome to strike the upper surface, vstitis the upper surface initially rather than the lower, the nature of the analysis is in no way affected. By repeated applications of the above procedures we can trace the progress of any one ray from its starting point on the (reference) phase plane, and can follow all of those rays which are spawned by reflection at the radome surfaces. There are an nfinity of such rays, but inasmuch as the amplitude of any one is reduced at each reflection, it is in practice sufficient to establish an amplitude level (or tolerance) and to ignore any ray once its amplitude has fallen below this level. We can likewise ignore any ray that passes through the monopulse plane in the region outside the monopulse plate, but even so we may still be faced with a relatively large number of significant ray contributions all stemming from the single ray incident on the radome. These rays provide our first (partial) sampling of the field distribution over the monopulse plate. To complete the 13

THE UNIVERSITY OF MICHIGAN determination of this field we must now sample the incident phase front at other points, and continue this process until we achieve an adequate specification of the field to which the monopulse is subject. lUihiwuires us to sample the incident phase front at a sufficiently dense number of points over a region which at least inciudes the origin of all rays capable of yielding any contribution to the monopulse excitation exceeding the tolerance level. It is convenient to carry out this sampling at a rect angular lattice of points whose separation (or stepping distance) A will be taken equal to that which was found desirable in the twodimensional problem (see S-L, pp. 42-43 ), and the sampling procedure is described in the next section. 14

A K A z - --- — v lrlIL - - - —.41 y &'iiase Plane Fig. 1:- Coordinate System. 15

E I Medium E E T2 I MediMum 1 I TL /12 Medium 1. I \21 E E 2 12 Medium 2 121 1 Medium 1 12 I R2 21 E E E E Fig. 2: Geometry for Reflection and Transmission at an Interface. Medium 16

\ <Y3 xa * Y2 ' Normal (a) (b) Fig. 3: Coordinate Transformations for Refleotion and Transmission. 17

THE UNIVERSITY OF MICHIGAN Im INCIDENT FIELD SAMPLING The monopulse plate pivots about a point on the (z) axis of symmetry of the radome (see Figs. 1 and 4). Let zm be the coordinate of the pivot point, and let zn be the corresponding coordinate of the radome tip, i.e. the intercept of the outer radome surface with with z axis. The monpulse plate itself is often somewhat irregular in shape, but for the sake of the following analysis it is convenient to speak in terms of a disc of radius a centered at the point zUZm. It is natural to choose 2a to equal the maximum diameter (or dimension) of the monpulse plate, but under most circumstances it would appear adequate to choose a somewhat smaller value such that the disc just encompasses the monopulse slots. A reduction in the value of a is, of course, reflected in a decrease in the number of sampling points and, hence, in the running time of the program. The equation of the incident phase front is Ix+my + nz0. (48) We are required to compute the field distribution over a disc of radius a lying in the center of the monopulse plane when the latter is aligned parallel to the above phase front. The projection of this disc on to the wavefront is therefore itself a disc of radius a. Since any line perpendicular to the wavefront and passing through the point (x',y', z') has equations (49) I m n the projection of the center of the disc has coordinates (xo,, zo ) where x -zm, yo-mnzm, z= (l-2)zm, (50) and the projection of the cone tip is likewise (xt Yt, zt ) where t=-nz, yt -mnz, zt (1-n2z.(51) t n t n t 18

THE UNIVERSITY OF MICHIGAN The equation of the line passing through (xo. yo. z0) and (xt y, z t) is therefore x-x y-y z-z o o 0 X-Xt *'Yt z-zt and inserting the preceding expressions for x,xt, etc, (51) reduces te In m n n z y = n —n z (52) 1-n 1-n2 This is the reference line for our sampling procedure. It is clearly just the projection of the z axis on to the incident wavefront, and therefore passes through the origin of coordinates as well as through the points (x o yo zo) and ( Yt, zt) It is now a trivial matter to locate the projection of the monopulse 'disc' on to the phase front. Since the projection is also a disc of radius a, moving down the reference line a distance a from the point (xo, yY z) shows the point A (see Fig. 5) to have coordinates x = -in (z + - ), a y -m n(m+ + --- ), (53) z -(l-n2)(zm+). l1l-n2' Likewise, B has coordinates x =-n (z - a m jl/.n2') y =-mn(zm- ( Lm ), (54) z (1-n2)(z a 19

THE UNIVERSITY OF MICHIGAN But any line parallel to AB has equations (see (52)) In mn x = - z+ C y = - 2 z+C2 (55) 1-n 1-n and lies in the phase front if I C+mC2 =0. (56) By considering a sphere of radius a centered on the point (x, y, z) and requiring that the above line have only one intersection with the sphere, i.e. be tangent to it, we can specify C1 appropriate to the points C and D and, hence, determine the z coordinates of these points. I then follows that the coordinates of C are x -lnz + ma y -mnz - (57) m T.2 z (1-n2) (m+ a- ), whilst those of D are ma x = -tnz - a am y = -mnz + a-, (58) m ~/1_n2' Z (1-n2)(Zm+ a1. ) ) The reference or starting point for the sampling procedure that has been developed is the point A in Fig. 5, the coordinates of which are 20

THE UNIVERSITY OF MICHIGAN x =-in +^ &- + mT,.a-N 4 \ M,Q y = -mn (m + ) - (59) V 1-(n2 ( 2 with M = N = 0. Sampling therefore starts at A and proceeds by stepping up the line AB a distance A each time. The oordinates of the sucoessive points are given by Eqs. (59) with N = 1, 2,3... and M = 0, and amplig continues beyond the point B (for which N = [ ) until the oorrespondn rays no longer ntercept the radome. For angles of incidence outside the "bakward oone", the last contributing sampling point will be just beyond T, corresponding to N = [ a + (m -z) /.1-n Having completed this traverse, the next sampling point is on the line AF a distance A towards E. This is the starting point for a second traverse parallel to the first, and the coordinates of the successive points are now given by (59) with N 1, 2, 3,... and M = 1. The third traverse is to the tghtof AB by a distance A ( so that M =-1), and so on, running alternatively to the iet; (M = 1,2, 3,.. Jndttotharijft (M = -1, -2, -3,...) of the central line AT. In each case sampling along a given traverse oeases when no further intercepts with the radome are obtained. The final two traverses are cartled out with M = ~ [|] and proceed beyond the points 0 or 0 respectively (for which N = [ ] ) until terminated in the same manner as before. The region of the wavefront over which the sampling actually extends in any practical situation depends, amongst other things, on the incidence angles. It is, of course, obvious that if a ray emanating from a given point does not intercept [x] denotes the integer which is Just larger than (or equal to) x. 21

THE UNIVERSITY OF MICHIGAN the radome surface, it cannot possibly produce any contribution to the field distribution over the monopulse plate, but if a ray does intercept, it does not follow automatically that it will excite the monopulse. Only by following through the laborious and time consuming process described in Section II can we discover whether it was necessary to have considered that sampling point in the first place. Neverthelsss, any ray which intercepts the radome is a pottial (and, in general, actual) soure of monopulse e citation, and since the time spent in searchi for, and oomputing, nitial intercepts with the outer radome surface is a very small fraction of the total computation time, it is prudent (if not necessary) to carry out the wavefront sampling over a region which is large enough to encompass all rays which could oonevably intercept the radome. Although the above procedure, in which the computer itself decides when sampling of the wavefront along a particular traverse shall oease, is aesthetically a satisfying one, the actual numerical program has been written on the basis of a pre-specified sampling area. The area selected is rectangular in shape and extends from EF up to a parallel line through B or through the tip projection T, whichever is the farthest. In consequence, M runs from 0 to + [a], and N from 0 to max. ([]. [^ ta+(zm n) 7] 22

THE UNIVERSITY OF MICHIGAN Z=m * m z / / / / / / / / / / / B / / / / / f1 A Fig. 4: Projection of the Equivalent Monopulse Disc. 23

THE UNIVERSITY OF MICHIGAN I 0 T ~B C /I" '.N — Fig. 5:- Geometry for Samplin the Wavefront. 24

THE UNIVERSITY OF MICHIGAN IV MONOPULSE PLATE AND POINTING ERROR The rays which emanate from the incident wave front are all parallel and of constant amplitude (unity), and if the wave front is sampled in the manner just described, these rays provide for us a uniform sampling of the uniform field that would exist in the absenoe of the radome. When each ray strikes the radome, however, reflection and refraction occur, generating additional rays, and this takes place whenever the original ray or any of its subsidiaries Intercepts a radome surface. It is therefore not surprising to find that the field within the radome is markedly non-uniform, and our knowledge of it is obtained from the information carried by those rays which penetrate into the interior. This information consists of the amplitude, optical distance of travel, direction of propagation and polarization of each ray, i. e. As, d (, m,n), (11, mil ni) respectively, at some position along the path of the s'th ray. If the radome were not present, each ray within th radome would be only a continuation of one emanating from the phase plane, and its direction cosines (1, m, n) would be the same as those of the normal to the wave front (1, m, n), but because of reflection and refraction at the radome walls, the actual values of (1, m, n)a can be quite different from those of the incident ray. The amplitude As must be real, positive and less than unity. In contrast to the situation in the two-dimensional case where the amplitude could be positive or negative, any 180~ change of phase (equivalent to a reversal in direction of the electric vector) is now incorporated into the polarization direction cosines (11, m, n1) The optical distance d. is measured in inches, and can be converted to an electrical phase by multiplying by k = 2/X, where X is the free space wavelength (also measured in inches). The ensemble of all these rays constitutes a non-uniform sampling over any given plane of the non-uniform field existing withi e radome, and it is this field which excites the monopulse. The monopulse plate configuration 25

THE UNIVERSITY OF MICHIGAN used in this program is that described * by Dr. I. Pollin and shown in Fig. 6. The overall dimensions of the pate are 9 x 9 inches, but since the plate has its corners chopped off, the maximum cross dimension is about 10. 5 inches. The plate contains a total of 32 horizontal slots symmetrically located about x and y axes through its center. Each slot is 1.060 inches long and 0.140 inches wide, corresponding to 0. 442 X and 0. 0583 X, respectively, at 5 GHz. Apparently the slots have been chosen to be of resonant size at this frequency of operation. The electric vector in the incident field is specified by the direction cosines (11, mn,,), and since these can be chosen at will providing E remains in the plane perpendicular to the propagation vector, there is no loas of generality in assuming that any motion of the monopulse plate is such as to keep the slots with their long dimension horizontal (appropriate to a field vertically polarized with respect to a horizontal plane, x = 0). Any displacement of the monopulse plate can then be treated as a rotation about a line through its center parallel to the x axis, followed by a rotation (or tilting) about a horizontal line. Such displacement is, of course, brought about by the voltages induced in the slots by the field incident upon them, and for the reasons discussed in S-L, these voltages are computed for a monopulse plate aligned parallel to the wavefront of the plane wave incident on the radome. The above mentioned angles of rotation are now specified by the direction cosines (, m, n), and with the notation shown in Fig. 7, * have sin a m -co a sin a, n os a cos a. (60) x x y x y Private communications, dated 13 November 1968 and 26 May1969 26

THE UNIVERSITY OF MICHIGAN Similarly, if the electric vector lies in the -"'n y" plane and makes an angle 2 with the x" axis, the direction cosines of this vector with respect to the original coordinate system (x, y, z) are (cos a coo Q2, sina sina oos +oosa sin 1,-ina cosoa coso +sina sin Q). x x y y x y y (61) In the computer program the description of the monopulse plate is fed in as a statement of the center point locations of each slot, and of their length and width. The symmetrical layout of the slots is employed to reduce the number of parametesatobe specified, and only the fiformation concerning the slots in the first quadrant, Q = 1, is explicitly introduced. For the other three quadrants (mirror images of the first), the slot locations are deduced by appropriate changes in the signs of the coordinates. The convention used in the program for designating the quadrants Q = 1,... 4, and the slots 8(1, 1),..., 8(4, 8) is shown in Fig. 6. The first number in the slot designation specifies the quadrant the slot is in, whilst the second specifies the slot number based on the first quadrant as reference. Note that the slot numbers in the second through fourth quadrants are those of the mirror-image slots in the first quadrant. The voltages induced in the slots are produced by the field incident on the monopulse plate, and it is to be expected that this field will be of a highly non-uniform nature. The response function for a receiving antenna, e. g. a slot, is generally stated for the case of an incident plane wave, and it might now appear that to use such information in the present case would require us to match the sampled data provided by the ensemble of ray contributions to a set of plane waves impinging on the slots at the appropriate * Although the present computer program has the monopulse plate data inserted as part of the main program, it would not be a difficult task to modify it to have the plate specifications read in as input data. 27

THE UNIVERSITY OF MICHIGAN angles of incidence. Fortunately, this process is not necessary, and we can avoid the matching by assuming each ray to excite a slot directly. This is in line with the procedure developed in the two-dimensional case. Sctions IV and VI of 8-L contain a detailed discussion of the effects of amplitude taper, polar diagram, slot weighting, etc., and lead us to characterize the response of a slot by the following parameters. (W Slot size: The dimensions chosen for the slot are not the physical dimensions shown in Fig. 7. The reasons for this are two-fold. Firstly, examinations of the surface currents excited on a conducting sheet by a narrow resonant slot shows that a strong field exists beyond the immediate confines of the slot, and the actual area 'excited' may be much larger than the area of the slot. Without a detailed knowledge of the actual slot characteristics, it is not possible to estimate the effective dimensions of a slot, but it was found in the twodimensional study that the number and locations 6f slots in the monopulse plate affect only the smaller 'wiggles' in the pointing error curve, whereas the dominant features of the curve are almost exclusively determined by the radome and plasma (If present). The second reason for choosing increased slot dimensions is to keep to a minimum the number of sampling rays necessary for an accurate pointing error computation. Experience has shown that of order ten (or more) primary rays are required to 'hit' each slot. By 'primary' rays we here mean those rays which reach the monopulse plate after a single passage through the radome layer, and without reflection at any surface. In the two-dimensional case it was found that for a spacing A = 0. 2 inches of the successive points at which the incident wave front was sampled, the number of rays striking a X/2 slot was adequate for the pointing error determination, but it is most unlikely that this same sampling D. K. Reynolds, "Surface Current and Charge Measurements on Flat Metal Sheets, Harvard University Cruft Laboratory Tech. Report No. TR-53 (p. 68), 1 August 1948. 28

THE UNIVERSITY OF MICHIGAN distance would produce more than 2 or 3 primary rays which strike a 1.060 by 0.140 inch slot. To decrease the sampling distance would be most undesirable fromi the #itCatof view of the computation time, and thus the only alternative is to increase the area over which a ray striking the monopulse plate can contribute. It is fortunate that such an increase in the slot dimensions is a physically-justifiable procedure. (b) Amplitude taper: The concept of an amplitude taper has been taken over directly from the two-dimensional study where it was found that its use removes the discontinuity in induoed voltage that would otherwise occur when any one ray hits just inside, rather than just outside, a slot. Here again the concept is not only a physically reasonable one, but also enables us to decrease the sampling Oate. For a slot whose mid point is located at the origin of the x"', y coordinates (see right hand lower corner of Fig. 6), the following taper factor is assumed: ' ytI,'<w Ix91< *t fg = coo ( v ) COs ( I Y < H < (62) = 0 otherwise where w and I are the effective width and effective length, respectively, of the slot. (c) Pattern factor: This takes into account the direction of arrival of a ray with respect to the normal to the slot. In general, the pattern factor is a rather complicated function involving the polarization as well as the direction of arrival, but it is convenient to treat these two effects separately. In the two-dimensional analysis, a sin x/x pattern factor was employed appropriate to a uniformly excited aperture, but if the aperture dimensions are no more than X/2 and the incidence is not too far from normal to the slot, a sin x/x factor is indistinguishable from a cos x one. We shall therefore take the pattern factor to be 29

THE UNIVERSITY OF MICHIGAN fP = cos i (63) where 6i is the angle between the direction of the incoming ray and the normal to the monopulse plate. We note that this can be written as fP (I,m,n) ~ (1, m',n'), (64) where ( ', m', n' ) are the direction cosines of the ray in the (, y, z) coordinates, and also note that this choioe of factor is not at all unreasonable. 1n the horizontal plane (see Fig. 6) a slot (bein the dual of a dipole) will have a cosine pattern, and though the pattern would be a contant in a vertical plane were the conducting plate infinite in size, the fact that the actual plate dimensions are finite will produce some reduction in the pattern in this direction also. (d) Polarization Sensititty: It is obvious that the slots will be polarization sensitive. The induced signal will be a maximum when the incoming ray has its electric vector in the vertical direction, and for other polarizations the variation in signal strength will be proportional to fI = cos 0 (65) where 0 is the angle which the electric vector makes with the plaite containing the direction of propagation and the y" axis (see Fig. 7). Hence, cos f= (l.,m'l n) ( {(l.mt'n')y "t where W, mt, n' ) and (I', m!,,n') are the direction cosixes of the electric and propagation vectors, respectively, and since, from Fig. 7, = = { (mn' n') + '(my +n )), we have fpol. ((mm +nn')+l'(mmi+nn. (66) m I1 l l 30

THE UNIVERSITY OF MICHIGAN The signal (in amplitude and phase) which a single ray induces in, say, the slot S(1, 3) is now obtained by multiplying the amplitude of that ray by all of the abMWe factors, and by converting the optical distance of travel into an electrical distance. We thus have V =,3 AliSi e (67) V (1, 3)A fPIPO' e ikds (67) and the total signal induced in this slot is V(1,3)= V(1, 3), (68) where the summation extends over all rays that strike the slot. The signals induced in the remaining 31 slots are computed in a similar manner, and are designated V(l, 1),..., V(4,8). Were the slots excited by a uniform plane wave, i. e. were the radome and plasma absent, the signals induced in all 32 slots would be the same apart from such small discrepancies as are associated with the finite sampling rate. With the radome and, perhaps, plasma present, however, the field that strikes the monopulse plate is no longer uniform or planar, and the slot signals will now be different. The pointing error ascribed to the monopulse system will then depend on the way in which the various slot signals are combined and utilized. It is obvious that with a total of 32 slots there are numerous ways in which the slots can be interconnected, with each ray equivalent to some form of addition (or subtraction) with varying (amplitude) weightings and phase shifts associated with the outputs from the individual slots. Although some information about the amplitude weightings has been supplied to us, an uncertainty about the electrical connections and the implied phase shifts has made it desirable to print out the computed values for the signals (ee, for example, (68) ) induced in each slot. For any given design of the electrical system, it is then possible to find the appropriate pointing error with only a trivial computation. 31

THE UNIVERSITY OF MICHIGAN For illustrative purposes, however, a very elementary form of monopulse system has been considered in which the monopulse plate is treated as four antennas, each occupying a quadrant. The signals obtained from the quadrants are then 8 V- n V(1,n), na1 8 V= E4 V(4, n), n=1 and these can be combined to give til12 i123 V23= V+V3 = IV23 e.if. (69) V34= V3+V4 341 e i41 41 = 41 41 e from which a pointing error can be deduced in the same manner as in the two-dimensional program (see Section IV of S-L). The 'vertical' pointing error is then = sin'1 4(-1 (t 34- 12) 0) and the 'horizontal' error is H r sin' { 4 '^ '41l 34>} (71) 32

THE UNIVERSITY OF MICHIGAN where EV and EH are measured in radians. The quantities xV and iH are the mean, normalized (with respect to wavelength) slot displacements from the y" and x" axes respectively, and for the geometry shown in Fig. 6 with X = 2.4, we have V =. 8698, H = O 8810. (72) Concerning the signs of EV and CH we observe that if V> 0, then /34> q12 ' This in turn implies that the dominant rays have arrived in directions 'above' the plane xSO, and must travel further to reach the lower slots compared with the upper, leading to an upward tilting of the monopulse plate in its effort to find a constant-phase plane. Hence, in the notation of Fig. 7, a positive value-of V represents a deerease in az, whereas a positive value of EH represents an increase in a. 1n other words, as viewed from the incident wave front, the monopulse points too high or too low according as V> 0 or < 0, respectively, and points to the right or left according as ~H > 0 or < 0, respectively. 33

* 14( (Ty' 1.500 4.65 (Tyl (Typ.) 3x450(Typ.) II i - M _ i C.(1,72) Cm D PoPt 6'I -C- -v -7-. '.o.* 4g.W.dm (3.A) *,. ( 4,1) I..I *Fg:McanclSeiicto4o h Monpule4lat. A'C Q=1 t t I..AA 4.34 3 y a.,. ** I ,AA A' 2'

THE UNIVERSITY OF MICHIGAN Xt xmx' Z" -.0- was vp z I yr y JI f y'=y" Fig. 7: Monopulse Plate Displacements. 35

THE UNIVERSITY OF MICHIGAN V PLASMA MODIFICATION Under some conditions of operation, the radome may be surrounded by a weak plasma, and even though the plasma can be treated as a lossless dielectric for purposes of analysis, it will still produce some effect on the field that penetrates into the interior of the radome. This will, in turn, lead to a ohange in the pointing error of the monopulse system, and the determination of this change is the main objective of the present study. A detailed examination of the effect of the plasma on the rays passing through it was given in Section V of S-L. It was there noted that because the maximum departure of the equivalent refractive index from its free space value, unity, is less than 0.02 for the particular electron density profiles supplied to us, it is possible to avoid the time-consuming process of numerical integration that would be required to follow precisely the progress of each ray through the layer. The alternative procedure is to simulate the plasma by a form of jump discontinuity in which each ray is displaced a certain distance A 1 along the tangent to the radome At:te plane of incidence, and the optical distance along the ray is increased by 61. In terms of the angle a between the incident ray direction and the (inward) normal ^ 11.61 x 10 12 t I tan a sec2 a (73) 61 SL - A 1cos 2acoseca, (74) where t is the thickness of the layer and I is a normalized integrated electron density. For the plasma characteristics that have been specified (see Section V of S-L ), both t and I are functions of the coordinate z(1) of the radome intercept, but are independent of x(l) and y(l) In view of the detailed discussion given in S-L, we shall here confine ourselves to a description of the way in which we can implement the above 36

THE UNIVERSITY OF MICHIGAN procedure in the three-dimensional case. Since the essence of the procedure is to rely as much as possible on the computational techniques developed for the bare radome, the analysis follows rather closely that developed in Section n of the present Memorandum. Statig with a point x((~), y() z(O) on the incident phase front, we first locate the point at which the ray intercepts the radome. Let this point have coordinates x, y() z). We now perform two notations of coordinate to produoe the system X21y2, Z2 (see Fig. 3), and if, for convenience, we take the origin of the new coordinates at the intercept point, we have (see Eq. 18) x2 L(x-x(1)) (y-y() 4n(z-z1) ), Y2= (A11+B2)(x-xO) )+(Aml+Bm2)(y-y() ))Anl+Bn2)(z-z(l)), (75 z2(-Bl+A12)X-x() )+(-Bmi+Am2) (-y(l))+(-B+An2)(z-z()), where A and B are given by Eqs. (15) and (17). Conversely, x-x(1)= Ix2+.(AIi+BlJY2+(-Bl +A 2)z2, y-y()= mx2+(Aml+Bm2)y2+(-BmI+Am2)z2, 7 z-z() nx2+(Anl+Bn2)y2+(.Bnl+An2 2 In the x2y2 plane, the geometry is as shown in Fig. 8. By virtue of the definitions of x2,y2, z2 in terms of x,y, z, the normal lies in this plane, and the angle. is auch that (see Eq. 16) cos0 Y2' fx+mn+n z) -A(Jnl+mnml nnnl) +B(ln/2+mnm2 n 2 sr (t7 where r {(Inll+nml+nnl)2 +(l2+mnm24nnn2) 2, 7(J 37

THE UNIVERSITY OF MICHIGAN Hence /l r2' tan.. and the line representing the Intersection of the tangent plane with the plane z2'0 now has equation ~l-r2 Y2= - X2 f7 If we step a distance A1 along this linethe displacement being in the direction of increasing x2 (so as to have a positive compont in the direction of the incident ray), we arrive at the point X2= A i', Y2 = A 1 1'r (80) and here we now draw a line parallel to the original normal to the radome. The equation of this line is clearly Y2- Al l- rz = -cot (X2- Al r) i.e. A1 Y2=" (x2-), (81) and the point at which it intercept* the radome is our new, displaced intercept where penetration into the radome is presumeds t take place. The corresponding ray has the direction that it originally had (of course), but its optical distance of travel is increased by an amount 6 (see Eq. 74) over and above its value dol (see Eq. 10) appropriate to the undisplaced intercept. It will be appreciated that because of the curvature of the radome surface, the normal to the radome at the new intercept point will not, in general, lie in the x2y2 plane. It is therefore necessary to establish a further set of coordinates 38

THE UNIVERSITY OF MICHIGAN tailored to the plane of intercept at the new point, and to proceed as we did in transforming from the c, yj. zb) system to the (x2 y2 z2) system with the bare radome. With the above formulation of the method, the determination of the point at which the line (81) intercepts the radome (always, of course, in the plane Z2=O) requires us to express the equation p = p(z) of the outer surface of the radome in terms of the (x2,y2, z2) coordinates. This is not necessarily a convenient thing to do, primarily because all of the coordinate transformations demanded by the bare radome analysis are merely rotations of base vectors, whereas we arhere concerned with an actual coordinate system which is displaced in origin as well as being rotated relative to the original (x, y, z) system. The simplest way to avoid this difficulty is to express Eq. (81) in terms of (x, y, z). From Eq. (75) and bearing in mind that the line represented by (81) lies in the plane z2=0, its equations now take the form (All+B2+.2r ) (x-x )) + (Aml+bm 2+ m ) (-y(l ) 1 2 1+bm2+_ )+ +B+ (z-z (1) &1 (82) (-B I +Al (-x() )+(-Bml+Am2)0y-y(l))+(-Bn )(z-z 0. These can be reduced to { nl-I r 2(A+n2 B)4 x^ -X rp2-r(lA+I 2B4 (z-z(l))= A(Bm3 ) + ni~r-r (X^*A.2B)y-Y) )4wr r(m+M 2B) ) (B I(+A8 2) (83) 39

THE UNIVERSITY OF MICHIGAN but the form is still not a convenient one: it is by no means apparent that the line does have direction cosines 1n, mn, nn, which it must have since it is parallel to the radome normal at the point x(1), y) z) The final step in the reduction process is one of those on which hours can (and were) 'wasted'. We first note that because (, m, n), Qll, ml, nl) and(l2, m2, n2) are the direction cosines of three mutually orthogonal vectors, 12 12 12 12+t2+ rt-T2^+ 2^1 5 1 (84) and lm+l ml+ + m2 ln+l n +2n2= mn-n +m ll 2n2 0. Hence, from -Eq. (15), n7-r An2B)=n nn(+n2)+n(1 lnl+ 2n2)+^n(mlnl+m2n2)} =-2 =n ni r (lI+m mn n n) n n n n ' and by a further application of the identities (84), it can be shown that (1 nm m-+n n) +(ll+m nm l nn2)+ + nm m 2nn2)2= 1 i.e. i l+m m+n n=- -P. n n n Similarly for the other factors in the Eqs. (83), leading to the following equations specifying the line: x l —x+ 1 (z-z ))+.A (Bm-Am (85) y=y(l + m (z-z_ ) + A1 (B(l-A I+ n 40

THE UNIVERSITY OF MICHIGAN In this form, the equations are quite similar to those describing the original ray emanating from the point x(O), y(O) z() on the phase front. The point of interception of that ray with the radome specifies zd) and hence, from (8), x1) and yl). The z coordinate, z = z(2), of the displasd intercept is given by (1) iX(1)+ (Z{I (l))q.Al(Bml-Am~J 2 +y {mn(ZZ (BI 2} {P fi), (86) where A1 is given its value appropriate to the original location z=z, and where p = p(z) is the equation of the outer radome surface; and having determined z(, the x) and y(2) follow from the Eqs. (85). The ray is now presumed to enter the radome at the displaced point x(2) y(), z(2) rather than at the point x(l) y) z(l) that it would have done in the abseoe of the plasma, and the phase (or optical distance) that is associated with it is d + 61 (see ol Eqs. 1(10) and (74) ) instead of d. Ray tracing then proceeds as it did in the case of the bare radome. 41

THE UNIVERSITY OF M ICHIGAN Original intercept X2 Tangent "o- Di\ 'line'.\. Dinter Y2 \ I Oute: surfi Fig. 8: Geometry for Displacement of the Intercept. laced mep p ace 42

THE UNIVERSITY OF MICHIGAN VI COMPUTER PROGRAM AND NUMERICAL DISCUSSION An examination of the basic analyses presented in Sections II through V suggest that the computer program for ray tracing through a three-dimensional radome will be a fairly long and complicated one, and this is indeed the ease. The program consists of a main program, four subprograms and 17 subroutines. The FORTRAN listing of the complete program* is presented in the Appendix to this Memorandum. In essence, the MAIN program controls the input and output data, initiates a sampling ray at the &ilent phase front, takes into account the plasma (if present), and guides a ray into the interior of the radome. If the ray hits the monopulse plate, the results are recorded or neglected depending on whether the ray intercepts a slot or not. Another possibility for a ray inside the radome is that it hit the upper surface. In such a case, the SECOND program (a subprogram) takes over and follows the rays through reflection and refraction at the upper surface. The THIRD program picks up As noted in the introduction, the programs for the conical and ogival radomes are separate. Although the general format is the same for each, the following descriptions in this section pertain to the conical radome program. We assume here that the incident illumination is from below the axis of the radome. 43

THE UNIVERSITY OF MICHIGAN the downward rays, some of which will either hit the monopulse plate or the lower surface. The rays that hit the monopulse plate are 'home' and are recorded, whilst the others are picked up by the FOURTH program and, if the situation calls for it, the FIFTH (and final) program that is similar to the THIRD. Throughout this ray tracing process, the rays resulting from multiple bounces within the dielectric layer are picked up by a subprogram appropriate to rays in that particular region. It is to be expected that in most cases the computations will be performed by the MAIN, SECOND and THIRD programs and only for angles of incidence greater than (about) 50~ will the FOURTH and FIFTH programs come into play. Subroutines are used to perform the mathematical, vector, and other operations required in the ray tracing procedure. In total there are 17 subroutines and these carry out the following operations: 1. SUBROUTINE REPEAT - traces a ray inside the dielectric layer; the transmitted portion of the ray is picked up by one of the subprograms. 2. SUBROUTINE PINT - determines the point of intersection of a ray with the outer (and inner) radome surfaces. 3. SUBROUTINE PNORM - computes the direction cosines for the normal to the surface at the point where a ray intersects. 4. SUBROUTINE ANGLES - computes the angles a and, as described in Fig. 3b. 5. SUBROUTINE NEWCOR - sets up a new coordinate system at the point of intersection of a ray with a surface. 6. SUBROUTINE TR12 - computes the transmission and reflection coefficients for a ray gotg, from air into the dielectric. 7. SUBROUTINE TR 21 - computes the transmission and reflection coefficients for the ray going from the dielectric into air. 8. SUBROUTINE CONAB - computes the quantities A and B as defed as defined in Eqs. (15) and (17). 44

THE UNIVERSITY OF MICHIGAN 9. FUNCTION DIST - computes the distance between two points. 10. SUBROUTINE TRANSM - determines the direction cosines for the transmitted ray. 11. SUBROUTINE REFLCT - determines the direction cosines for the reflected ray. 12. SUBROUTINE PLASMA - adjusts the path uad optical distance traveled by a ray when it goes through a plasma layer. 13. SUBROUTINE RECORD - computes the intersection point of a ray with the onopilse plane. 14. SUBROUTINE SLOTS - determines which slot, if any, is intercepted by a ray. 15. SUBROUTINE APSLOT - computes the signal indued in the slot by a ray. 16. SUBROUTINE VECTOP - computes a vector cross product. 17. SUBROUTINE VECTDT - computes a vector dot product. In the theoretiol analysis as well as the computer program it is convenient to specify the direction of propagation and the polarization of the incident wave in terms of the direction cosines (I, m, n) and(l, ml, nl) respectively. Since it is not always easy to visualize the particular wave which is implied by given values of these quantities, we have also employed the same angles a and ac x y that were used (see Section IV) to describe the rotation of the monopulse plate to describe the direction of propagation. Such a choice is natural inasmuch as we have specified that the initial alignment of the monopulse plate shall be perpendicular to the direction of hiaidence outside the radome and, hence, parallel to the external wavefront. In addition, the polarization of the incoming wave can be described in terms of the angle Q referred to in Section IV. 45

THE UNIVERSITY OF MICHIGAN To illustrate the application of these angles, consider the following situation. For a plane wave with propagation vector in the positive z direction (axial incidence), the phase front lies in the xy plane. The pbllaruiion angle PSis then the angle between the E vector and the x axis, and ax and ay are zero. We now imagine any arbitrary alignent of the wave front to be taied by an initial rotation about the x axis through an angle a, followed by a rotation or tilting of the new plane through an angle ax away from the vertical (see Fig. 7), and if we do this, the pertinent polarization angle Qis the angle between the E vector and the x" axis (see p. 27). Thus, for example, a vertically-polarized wave incident at an angle of 20~ below the z axis has a =200, a =0, Q= 0, x y whereas for a horizontally-polarized wave incident in the horizontal plane but from a direction 15~ to the left as viewed from the monopulse plate, a =0, a = 15~, Q= 90~. x y The expressions for the direction cosines of the propagation and electric vectors in terms of a ay and fPare given in Eqs. (60) and (61) respectively, a (ANGLE X), ay (ANGLE Y) and f2 (OMEGA) to be employed as input data. The other input parameters for the program a the stepping or sampling distance A (STEP OF DATA), radome refractive index n (ANEX), preence (or otherwise) of the plasma (KPLASM = 1 or 0), equivalent radius a of the monopulse plate (CONSTANT A), and the number of increments (M(MAX), NMAX)) over which the phase front sampling extends. It may be noted that the tolerance level i..-. the level at which any ray amplitude is regarded s so small as to be negligible, is not an input parameter, and has been set at 0.01 in accordance with out findings in the two-dimensional study. 46

THE UNIVERSITY OF MICHIGAN Concerning the choice of values for CONSTANT A, M(MAX) and N(MAX), we observe that together with (1, m,n) these parameters define a set of coordinate points on the incident phase front from which the sampling rays originate. The procedure for specifying the sampling area was discussed in Section IM, and in essence the area that must be covered is the projected area of the monopulse plate extended upwards to include the projection of the frontal portion of the radome. For computational purposes, however, it has been found convenient to specify a rectangular area which includes the above region. The rays that emanate from sampling points that do not produce intercepts with the radome are soon abandoned, and their initial inclusion in the computational sequence does not add significantly to the overall computation time. Indeed, the execution time required to find whether a ray does intercept the radome is only a small fraction (less thea one percent) of that entailed in tracking an intercepting ray through all of its reflections and refractions to the monopulse plate. The coordinates of a general point within the rectangular sampling area of the incident phase front are given by the Eqs. (59), and sampling is carried out at all such points having M=M.., M max' 'max N 0,... Nma max with x = N =max( 2a] +(z6-z ) 1 ]) (87) In order to test the program and thereby obtain some assurance that it is operating correctly, a very small number of test runs of the conical radome program have been made for a specific angle of incidence of a plane wave on the conical radome described* by Dr. I. Pollin. The outer and inner radome Private communication, 13 November 1968. 47

THE UNIVERSITY OF MICHIGAN surfaces were defined by Eq. (2) with a = 6, = 1 (outer), 1 a =, X=4. 191 (inner). The quantity zn ( I ) above is therefore 1, and zm = 43. All dimensions are,of course, ininches. The frequency employed was 4. 918 GHz for which X = 2.4, though we may note that other runs have been made for X = i (f = 3.757 GHz). For a fiberglass radome, the refractive index n has been assumed to be 2.5. To test the computer program we chose the case a 15~, a = 5~, Q= 10~ x y a = 5.25, no plasma and took A = 0.4. This value of A is somewhat larger than is desirable for an accurate estimate of the pointing error. In the two-dimensional program we found that a stepping distance A = 0.2 was necessary to obtain a pointing error that was not unduly sensitive to the particular choice of initial sampling point. For lack of any other information it would seem prudent to choose A = 0.2 in the three-dimensional program as well, but since the objectives of the present test were to check the computer program, and to provide spaififc outputs against which to compare any subsequent runs of this program elsewhere, the larger value of A was entirely adequate. We may note that the choice A= 0.4 reduces the number of sampling rays by a factor of 4, and decreases the running time (and cost) by the same amount. To a person looking in the direction of the negative z axis, the above incident field is approaching the radome from a direction 150 below the horizontal and 5~ to the left. The selection of a value for the polarization angle f' not equal 48

THE UNIVERSITY OF MICHIGAN to zero was made only to produce a case of complete generality and to avoid any suggestion of symmetry. With the above input data we have, from Eq. (87), M = 13, N = 37, max max implying M max(1+2Nm) = 999 sampling rays of which about 500 will intercept the radome. The computer prinltd out a) the induced signal in each of the 32 slots, b) the combined signal for each of the four quadrants, c) the signals obtained by adding thoX in each of the two upper quadrants, and in each of the two lower quadrants, d) the signals obtained by adding those in each of the two right-hand quadrants and in each of the two left-hand quadrants, e) the pointing errors e -6.5 andH = 9.1 milliradians for the vertical and horizontal planes respectively. These pointing errors are about twice as large as were obtained with a similar case for the two-dimensional radome, but this is not unreasonable since the rays now undergo many more reflections and refractions, leading to perturbations of the wave front in a second plane. The signs of the pointing errors are also reasonable. According to definition, the pointing error is the angle between the direction of arrival of the signal outside the radome and the direction in which the monopulse plate points, i. e. the normal to the monopulse plate. If eV > 0 (eV< 0 ) the monopulse points too high (too low), and for eH > 0 (e < 0) the plt. is pointing too far to the right (or left). With ax= 15~, H H the rays which strike the monopulse have come through the lower surfaces of the radome arthe dominant contributions are provided by those rays which then proceed directly to the monopulse plate; in addition, however, there are contributions from rays reflected off the upper portion of the inner radome 49

THE UNIVERSITY OF MICHIGAN surface, and because of these it might be thought that the radome would tend to point above the true direction (i. e. give eV > 0 ). This is not the case. The action of the monopulse plate is to rotate in such a way as to take up a position in a (mean) constant-phase plane. The rays which are reflected off the upper radome surface must travel a greater distance than the direct rays, and for the chosen angle of incidence, these rays are incident primarily on the upper slots. Accordingly, the phase df the signal induced in the upper slots should exceed that of the lower slot signals, and this was confirmed by an examination of the actual slot signals. To achieve a constant-phase plane, the monopulse must now point below the direction of the main (direct) rays, corresponding to a value of EV less than zero. As a second test, the previous case was now run with the fiberglas radome replaced with a dimensionally-equivalent styrofoam one for which the appropriate refractive index is "s 1.015. Since '3 is so little different from unity, the radome should have a correspondingly small effect on the field passing through it, and the test was run with the expectation that any small pointing errors computed would be due primarily to the low rate at which the wavefront is sampled, It was found that V = -0.72 and ~H = -0.077 mlliradlan. Compared with 2 milliradians, EH is certainly so small as to be entirely negligible, but ~V is almost 40 percent of the maximum pointing error that is tolerable, and is therefore too large to be ignored. Part of this is undoubtedly a real effect indicative of the actual field perturbation produced by the radome, but our belief is that a significant fraction of it is due to an insufficient sampling rate. This could be verified in the same manner as we did in the two-dimensional program by systematically displacing the initial sampling point a distance A/ 10 in each of two directions, and thereby generating a pointing error curve which is an oscillatory function of the 50

THE UNIVERSITY OF MICHIGAN displacement. A must then be chosen so that the total variation of this curve is within an acceptable range - say, less than 0. 2 milliradians, and it was on this basis that we were led to choose A = 0. 2 in the two-dimensional program. The third test run that was made was for the same situation as in the first run, i.e. for a fiberglas radome with ax = 159, y = 5~ and nf 10~, but with the plasma present. The pointing errors that were computed are V = -7. 1 and fH 8.5 milliradians. The change in pointing error produced by the plasma is quite small, amounting to only 0.6 milliradians in each direction.

THE UNIVERSITY OF MICHIGAN VII CONCLUSIONS In this Memorandum we have been concerned with the determination of the pointing errors associated with a general type of monopulse system mounted within a (three-dimensional) oonical or oglval radome that may have a weak plasma surrounding it. n oonoept and formulation the method that we have developed is analogous to that previously discussed (Memorandum 02142-502-M, and referred to as 8-L' ) in relation to a 'two-dimnsional' radome, and it is suggested that the reader acquaint himself with this more simple case before attempting to understand the more complex analysis and procedures required by the three-dimensional geometry. We have omitted from the present Memorandum such details as are common to the two oase, and have ineluded here only thosefatures ae f rs aare ique te tonal tomety. Nevertheless, it has been our Intention to include full details of the complicated analytical techniques which the geometry entails, and print-outs of the resulting computer programs mai fbefmsidi'Lthe AWpIidix. We believe that these programs constitute a unique and versatile tool for assessing the Influence of the radome (and, If present, plasma) on the behavior of a monopulse system, and because of the considerable experience gained from the study of the two-dimensional problem, the programs themselves are rather efficient. In spite of this, the running times of the program are not inconsiderable, and the contract funds available have permitted only such test runsm as were essential to the checking out of the program. As an indication of these running times, we may note that each of the test runs 1 and 3 discussed in Section VI took between 3 and 4 minutes on an IBM 360/B7eoospitvr, and were the sampling distance reduced to no more than the 0. 2 value required for an accurate pointing error determination, these times would increase by at least a factor 4. The times will be greater still at wide angles, 52

THE UNIVERSITY OF MICHIGAN and yet greater for the oglval radome where the computation of the ray intercepts is more troublesome. To compute the pointing errors for a variety of incidence angles (in two directions) and polarizations would therefore require no small amount of computer time. In conclusion, it should be emphasized that because the entire approach has been based on ray theory, the results obtained are only approximate. It is our belief, based primarily on the detailed studies of the twodimensional geometry,,tat the values for the pointing errors obtained with these programs will be aoourate to within 1 millirdian, but prudence would dictate that before complete reliance is placed on these data, some attempt be made to verify these conclusions experimentally. We await suoh confirmation with considerable optimism in the outcome. 53

APPENDIX COMPUTER PROGRAM PRINT-OUTS 54

MAIN PROGRAM (Conical) D____U I M F N I W 1. ( 3),r (3),Nk ( 3 j).,lP~Lty3~ I XL I( 3) XIl I( ' ) XNI (73)I? X2()XrA2 3)XN?(3) 4 XL4(3),XM4(3),XN4( 3),~...-... 5 XL'i 3) XM15 (3) XN5 (3) __ __ __ X L 6 (3 ) XM 6_(3),X N 6 3) _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ I) I MINS ION AilPVAL (L 8 CII:.,'IL H-X V/ALUI-( 4,8)SVALUL.... C ( Vi 1)pL FX V., V1,,3, \'4, V1.? V?3, \/34, V4 I ~IAL KflU- M ___________ C(''ON \P LI IF/ /I( C, AAXNI- ANF )A IA P/ 3.,1 4 1 57)2 P7/ lA 'I A I )F-LIT'A /l / I I) F LT/'06 59/ PF-AL SLlTY(4, ) /3.78,-0. 78, v- 'A 8,Li.78, 5.30'6-?.05,? *0562.Li~ 3 3"2 *3,- * 305 3?. 3,2 I 3'i, 3.8 2-3 1.3 3 3,38 C )6I~04, )= *5 659AX -3, 65 )X=37, 65PLSI'A=0,k' AN9XZ-S.F H t 'P I. F I 6, ) 3I0 I ) TAI\,FxAG YI1H; )TP,1AM A IA COPEP.9ALItI II ANOL- Xk=A' C,LI.. F YLhA0-iGFLLY ',L.2LLMFA,5 IF A1I PL'~ *) IO )) O T.4 11F(,3f0/ ) _\_;___ 9Ai. F vM A9O.)9O~9MvlAX 1 30C F~U'A ( 1 TH PL-ASMA' I__I_____ ( 31 1 V I 11F(5,30 - ( C C, C, ) (1) C.L):(;I (k I( A )-C IS(CGL)+SIN.GA) -. I f: I1I:AA 1 ( 3/ 2 IL,1I 0.3 55 _____ _~~

aA-' li (^A t.\ I'lli)__ __ - - - - - -. ~1.9'/- 1!,'t,1- '::.\ A A LI+ I-i,/ J ( A.<\ 7 1 (. 1.. );() T n'jV li- I M>-^\ t"1 ] 1 1 1 -I 1,,..1.1..1......\. LU ( I, ~J ) = (, _.... ~._____________________..___ '! =l.( 1 ):-I.(1 /I.('m /I.(3 ('. — ( Fh.'-Cl':aA ) r. - = 1;;. A\ A;:::i)j- 1 | +4^ 3 * ~4 fl, (:,.;=((-?:;:( (...... _,, ~..._.:: (. 3 /= -" ' p? - 7. 7 7 P,,'D1?:771 T'r Tb 72.T 2r~,n, I '1i-:(,. T h ~(n).bllc, ~-rii."A' mi,(,) \I1vl'-RSI-CI']t:lI FO'R ZF. I)')o _________ *fS(^1 (,4() _ __________ ______________ ___ ________ _ 22 (1 ' Ii c-..A?? Cr, i It.:lt 779 (-.^ C^ /?. X sr; 2= U Ak' ( 1. -. (3)':. (3)) Hi' KK T T"~4=1 t!.,I,,. ~... ~.. -... 7 =i i.....,.,/7 X-!.-j!=77. ':;. T _ _. i,~ ~ ~ ~ ~ ~ ~ i:,H2I ''l, FA >' TC:) 4 A-, N: () — J)',- T,,,:1SB P V) tS(r i V"" —.. ( 2 )::.( 'q ):::XC2;'+XU:-1 ( 1) 'L =- (:: L (2..Y___________________________________,-(,''~ (I. ( 1.)) * IF 0:, 001 ( );ul% TI)? 7 T,,i:=7 (7'-:1.(3 )/t ( ). 7. T 'i 2- A T 8,i:,(I.(2) ) *I. Fo?f01 )r, Tb '...........-...... / T = 7; ( - y ( 3 ):!/I. ( ) ') o -( ii'!;.A *(VM) IN^T-LSTA n )(;R TFN0 8 I-I (. ___. _Z__________________________'.._.!_L__ )____TO_.... 7.7,L~..j Jj(~jj Tb' 99 ~ li~ I Y9 s- = ( *(.(-. (3I) Tll 87)./,.I ',, > 1 'i ' ____?.n, ^) '[ 1 I *.,|t.A ',; -1 (, 1),, ( ^ C - X l L (2 l...i.-Q.,...-_..../',j'l-........._ ()) { t T. g It- (,,,", (. ( 1T) f'. 1^. 0 )t':f 'TtO ~'; ' 1; "T 1 =71 y ^.(3 ~/.(1 )....,, I r' i i,;*! '.,-t, - 7 '.(.1..............~..)......T.. i 'Jri ( 1,.,, / - T, — - ft(. 1,1 r. j, ( iT,5 -56 _ _____.

SUB3-PROGRAMS - SECOND, THIRD, FOURTH, FIFTH (Conical) CI 00I I' T T ~,IIFii ')I1I (>i. IF (I).AAX(,Y(~/.HLXIY1I I,1KTT(4 I ( I C. T '. J p (/ 1 * KH /,r~) (~G I' I ( /.1 * * T I *I C. I t4 I 1 bj" C - I HI \1 A T' F pI C 1~) H\1I 11)r' AL 1) J* - H('I'II I 1\ T II 'TR1T. t N CA LI. F. 1\H.' ( \ 1, s E C,/ 1, DLI)!, DN CI- I 1) Cfl ( I T\5 A 1\I T.A f $ 1 ANI) (;APIA C. A. ('I 11\,A I (, I qt S~ 9C - C.A l. A i I C&1.`S(l_ u1,G AMAALI-PHAMf-~TAOKI ~CK) I-(K I C. K (; 1)(r 1 TiF~ (KPI.SO'A *I)* 1 )(; *I'H ('2 2 C-AL-L I1 AS — A (Xi, Y1 1.,L,,N I),(AcA A AL~L T_',JILLA ~ 'Z ~.J~, A T TH,2 1 0 T F ( )) * H1. () ';f 6 T (Th19 in WI C AY U III PASS F~ - P 1 fiN " ) I IN I~1 TI) 2 Co~ (',A IP I? A(1 AI IA F kTA R 1 2P T 1'.P, Pt1. 2 VTI?V) A I' I- A, IGA I X1 Xi)1AXN C C P 'I I HI CT IH CSI 1\1F S (I1F TRkA NSMI\ TTT FI) RA Y 1 XI 1,r' XH\'I., (AiI A k XL I X!A NJ, XNJI2 -IK ~JLL2 —TA-ss lHiu~ j1 i EILs-R C CF4.)) 4 J~i\)-)~JHIT'FP.SFCT PIIINT (>A11I PI-4'T( FLTAA,1Y,,~LQY.Z2KCO I( T IC CT. 1 )(';I TI) 9q )AI=)AHlI ST (X 1 Y 1 vZ 1 A~E X,9X2,9Y.2Z.9 2) _....,... r L L [.(Jeo"'"X(-? X',Y?-/2,II)NtDi\I)GA,,1A LA LI f ~ F x I 2 X"? (,,;A,A Tl A G VA -) -i -.. —I-.'.-...C ('T i CV * T * c)() T I QL A__ _ _ _ _ L I? (.A A l P I_________I____ U A 'T ( C( (A, —,T? l ) *c + ( IA8):T?21 V ):c P /V; Ar=S 'P T ( ( 2 - C ~I11pr.((AW, IAFTA,, ALIAT A 1PIT TVXLXA, 1 xP I /., X 4,I 4.vX I1 -~ * I. 'I * Zr- J( )T. I 'II3 7 (. I~j'57

J -1$= \/A I- II i-, f _____________ V A I P1 I. ( AVPP IX _ A x 1 '"'1-I C dlii- \1F 1CAII 1.r.~~ ( XLI X"/ I X, ~I T YA X? AA, N3 PANi A AA iv 7 ATL T 7 F I) P XO TI. 44A AAXY _2 1X1 J;,\d, 7IK-ILl 7 NI) 4 I) =I 14)~'\A - i ( P A C.I H lf:~H) T~ 7~' t I F) I HII F) T 44 Z, A D'1 LTA.,-:rLlH.KI I.\ T1 T;i. 7- =Cdf I X Ab o PPI, Al-)P_ L 2TI- (AkS AII) T.AL)T (1.01 )(h TO 994 I F di* L * 0)6 To 40 0 3 =,%l 1 I lT I i- ' IL - - GH TO 9' d I..4 1=G194O4.4A1 t 1. l1X d F AL-6XOH 1IAI- H __ _ 1.-.T T_ _ __ _ _ __ _ __ _ _ __ _ __ _ _ __ _ __ _ _ __ _ __ _ _. At G "), F ALt (.VIJ.

IH -I T 'I A~h ( A T I P~( ( "3), 1kF Aj ( _____________________________________ V V" +\ V*" _______C.1'.\A.L i Vit V0~ 1 A=( VS( /3 j1-1 2= I A? 1\ ( A IA( (i A I G V2)4 k F A L (V ) ))... ~ P~Il T '~4 A 1I A H? (Al r-)A( ( Wl2 )t, PFA I ( \I?4 TH 14 I A1 Ai~2 A 1M A(;A( V4 1. ),P F A L (V\4 1) (Q flfi4 U A I( X / 'X, ('1. t V/, I3 v / 1. 1 i F ( C' A) 0 V, V Is.A, PHI 4o 9 n I IF, ( 'I, I4 fl (v's V? V4 _________ ' L..K(~i2!)XLV12\)11?APHI 12, 1 4'3 /2A,PHI2 \/I'01 ~ 8/~ PHI.3t&PHI 1]2) (42f -P-i)PI I 4 I4-1P5.5) _________ (I ('ViP 2 284 L1I- PIT - 59

x I xl ix y~I-~fH)X XFY.,2,Y272AATHX, _ I I- I )lY SI 1 K( L T 4 Al t1-\S T I' 'H I' ( )T,\ ( 3 ) _, H(3,X. 3 ~L t3 X 4 (3) XMI (3)..... *i X J ( 7 X O.I7 ( X, N 7(3 C'lr'JLI-X V/AI I) (4, ), S\/ALI IF I-AL VC. 1\1i 7FI f= 7F I) CA LL I) I1\1T ( 1Di- TKr A AX?,Y?7.?, XI-X,Y 5,.5,KTIC, 1) __ K (KC. iT. 1 G fII T H 1. 7......... I ) 7 = I ' f: X ~ A A, _ _ _ _ ___ _ _ _ _ A>\I. P P'k ( x5, Y' 7 '5,I- IH.- 1) j 1f\ MI 1.0 T=1 1-3 C7ALL1 FC, HlK\ (XL,\t M1 X1i, XN CA, C, XLA,\ XA 1 T ( K I C. le (;GT. 1 )(;I I Tol 1. 7 A.: A kI=S( T (, (A T 12? P ): C,2 (;B~ 1?V ):c 2 -1=KA K C A I;AI CALL - L CT G A (f AL-P HA, R IPR l2 VXLI,X MI, X\1lI 1 q\ A kIk,XL 2,XI.I2,9X\12)- -...... C-A 1.1 hCHkD(LI,9Xi.? tX5,qY5, 75,9X69,Y6 t76,P)1)K) n CAl I '11H I R I) ( XL?, XI,?, XI'?, X 5 Y'5, Z 5, AA, T NDEX, DFLTA, I9 LF- I 1 P1, ',L, LSIZ FN1) 1, LE D2,9S L-T X-,S LLLTY)T......... >1 R )'; (ITYIJI4 JJSP L I -(41, I n-l () C) T) 5 1 C, =C, IR A K(C 1T'') ~ ' i._ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ L'F1 (IIIJ =VALUIFI.3 (IIIJJ )+,.';AL U E................. 1, /o Ak Y 1. 3i, 1, XI' 'A -. 14..T A 1. 'l-AJAX+Y5.X.L3.

v 7 ' 'T / V N Gf T ' 3 t '9e77,7 1~7.y i.AI-1 ( ' i * II) )X I Y i-,7 X Y X1 I) I) A h, -J 7 I 17 Ai (17 A Pi TF IN )( Ii 311GTH 7 1 P1H I,AI TI I, I SI JJJ, ANH PSH\'),L(i S_)Y ( -11\1 KC.1 1\) 1 3 IAK\ )1 I.IA1A A (. AkS ( APT) S *C.S LI. (.(1 )(t 0 _________ Al IIIiC. i~yl r I A PT R XYAYvX1 P,X7 \1 IA IX I-JT,I~ VAI U IF AI TI JSLAUI " C'. H"= l( (. III)QI AR = H ( ) TKC.Hi\0 ).O T - 1.7.... Si K12 SI. AK H T) A. I= IA \ FP T I. C S V1 _ _ __ _ A P I~ A 1.P.- 1. TI 1 ~~; N ZF I l) HA1. I I (I_ IIII- Y J V LF(IIJ~)+VLJ 1 k ( e f ) - * IL. LOT1 Y'd 4 178 17Xf-L 73X F d1 A01 'v A N R RJA.NIF (3) XN (3 W<L3 3) TX'(3),X_ A N3Y(73)I 5y5Z IC * L I( ), X&7 ( 'T ) XN (3 72 I X 1. Ii I T( TCI E. YrjC',fl T C, 4Q. I..~ C,, 1T GI I 1 TI)O 1 71 A C,~(l. Z ' )1 '):c;:c2+ ( (m:c k I2V )::2 r, -F CT (1A,GI', Al1. PHAA, R IP 2 N)RI?/X L]1,X MI X X\1 I 61

A V \'?. 2, X i 2, X iM1)\__ _ _ _ _ __ _ _ _ _ (i O i 1) 2)4 A f I.A1 \ 1S1.AF~IT(A IIT, X 6,Y 6,XI-2X M?, 1~ I Tv\ SI1. T Y I I I, I3, IJ AMPiP SL )....... C. S ( 1 H ) fif] T C,1 A K_ __ _ __ __ _ ___I_ VK A II I=' C,,-'IPI" ( AK ~ li'.) - AMPI I '\/AlA I-=CPX (S\W/ ALA1P JI t T. TJ... JF 12 C l.. A - - A P. (G t.kt l.PA. i.PT-.T...Pq. 1 ~A$~, XL3 0, X~1, XI3 1-(7 L X9 7h~(' 7I 23V _I2 GI Tr.) 1 I1 CA I P- 1Tj1TH ( X1.7 X XM7,t X N 7, X'i, Y,, 7 5, AA, I N) F.X,F )L T AD OL Tf~ P ATH, ____IL_ 1LE I11, L,.L a! IuJ~LI.ZLzfI-, LN ) 's I)2 f &I)IX4SU. IYT T f C A kS CAkPT ) - LT. 0.0r)1. F) 6 TI( 1 7 C.ALL A) PSL OT (APTvX7,9Y7,9XL 7,9XM7, — 9 1 IL TY',. STY, I II,J33, AMPSL_ _______1 Ii~ 1~ (V<(11\ T 5 K C, T 15 r ' C.N K; T ( K I ed Q o!- I- I I;-.. J X-LA PLA(-IP L___________________________ lJ' I \1AC..F II JJJ + 1A II 1', (~ Ik P. (APT )L LF. 0.01 ) GO TI 1 7 V(JJ C- L 2160 TI) 1 7.... 1 ~I A,11I i, PATH, bAMPI., I SI(N,~l 7FNID1,7NI2 51Th SIC. Ii /1 4,V),SL(ITY (4, A I: Y-JK!1),X L(f'..fX LL.L.N..3-.L., _______________ ____62

(}f~n P -X VAI,, h (4,t ),- VAI-^'i 1 '- 'I K(;I< n' r '.;, VAL.iL://l/K/CuN,'!X\t1XN/;,[/AIEX..X.,1=1 ____1 1 I S 1.... A.) G LLLL6AL_..X.-_.__._____5_K 'n / 7- T,'t-?7 i-nlt)? ] I - ( K I. T 1) (I 1.7 A'IHi= A I ll-+i)T SI' ( X?, Y, /.2? 1., X5, Y5 7.5),1l 1)7= 1 M!I' 8X;XA^.......... - A ' A. l i'if ) i ( X, 5 Y', /'5, S )l1 L Z t, N) i j 1( I = 1,:..'-,,... -..'.... n, ( T )=-,l ( I( ) _______ ___ t__A I. C,(l'At ( XIi, X'i,I)fI, A, (;:,r(AMAA ) CA L. TH IC.Hm ( X. L, X Mt XN,;A, A (;, 9X I 1 X l, X l ].) CAt.L Ai, SI.-F F (XL,INX,- AIt A, ALPHA, IFTAtKJCK)..... -.......... IF(K ICK.(;1. 1. ),I ( TH 1 7 CA l. IL.? ( A,PH A, ETAr21PT1.?PRV 2 ) T................. 1. (.A^, XL? X_, xiX2 )... C AI I C C I, ) (, XL, X S,,- 5, XZ, YX,,,6, )1), K) ______-1 _ (1..... I(: J..lJLLL______________________._____________________________ T(7;.IT. Z - I G )(;'r Ti ) T) 30 - 1i )'?.................. 3.n C.A. -I FTH(XI.?, X 2,X\?, XS,Y,Z 5,AA, IN D X, )FLTA^, D)FLTR, 1 PA H,A ii R,L, ISIG'Zh -,, ZEND2,.SL TX,SLJ TY. -...............-. i.,i 1 n 1 2 C.\LL APS1.f 1 lT(A(.PTX^,Y6,XL?,tXM?, -1, '(I,.SL!..i' 1. X1..7, Xi,, X, ', X... ).... 1 t1....:....:,SL fLTY,Ii,JJJ, MP.SL)............................... ' ];1 1 vSL.CSKO,", rl r'i,.: Ar Cl I( t, T T I. I AU)= ALIE( I I I, JJ).S VA.. J.~.SV.E,,l..i... 1,,i, ) ' — 2.:," LL 't S-:. (.GA,_G b_, A LPJAA.B. ET.A_ _T.2P. T 12V_:XL1. X M1, X I',l_ X2r*,i,tI-ji +1. - T (iI PA I'X, t)F L TA, 1)FL'TI, AAt XS YS Z5,X: X M3 X L3 'X X3" PAT H, i t PTX p /,X (,Xi-7,X1<7,X7,Y7Z7,lt ), t K KLJJ)..... /' ' I i T,:, -j, 'i i-( '. 1r.! )(;11 T i 1 3 M-:' (:. * - r ))(J -r i... 1 (............... __............. ___. IT- (/'7.t T. /TF i")(;l' 't'll 3) V (' 1 1 1 ) T(/ < C'. I. 7- i -rl ( y.'Xf( xl, x1\i!' X' Y 7/- tAA, TNI)FX, FIL:TA, Fl I. t, 1 '1 I, IT L, I SIG. Z7 i EN t 1, 7Z iI')2, SL')TX, SLC OTY). I. F i 't 1 /. 63_ __

1)3 =1-1n4 PAT F(A 1P S (~ '\ T ) T, L, * 0*)1 )(( 1 1 7 C. AI I APSLUT A 01)PT,7,Yi, X LY, X v17, 1 SL.iY I I)TY,? I II W, I IjA M P SI.) A s=K C.' I \i \* I: 1I C.;. S f '.sJ' A K( C I. A ~.A PSL1. K. 1-.. —.-. --- "'.A Lt:CP L X ( AFIPR L,All IA0 \1\ II T!, dI V VALF (I II J dId + SV ALIFI 4 1F (Ak ( Fe'P 2 GOL. 0.0 )7 O. T 0 I e'- I I'lI"\ SI L YX (4,8 ) L I IT Y ( 4, L I '-'SI TI! A N (3),XL (3), M.(. 3),XN 3),- 9...._. 1 L 1( ) X'' (3 q Xi\1(3 3i Y I.I '3),X'33, xi33 ~~X4'/ 1 ),X 4 ( 3) X1\14 ( 3 XI _ 'Fi' X f-I\ ( 3 ) X...-.-7.. ---.- -. il P L FXY /A L IF (4 V VA LI I F I)-A l. K (. I Ij'i CL L P T(IE L T~ A-AX2,Y 2,Z2,X L,9X5-,Y.5-,25.,K IC t1) T T(.C. *G T. 1 )(;I- T f) 17 r Al. I i (~,5L5DfZ )... LL ArL(Xi. 11\, GANIA, AIPHA, kFTA,0,K ICK) CI 2( Al.IPHA, 4FTA, klP, T12p, R 12\/,T12\) A - _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ A1~.'-\AtT'IYLPXsA?,XN2). -.-. —. V........ K "II dFC.IPrGI)A 1,XA P,XY57ixS6,,,)) ( * 3. ).'; 1 T! G O_ __T H__1__7 12; l H,R1 P M V, L M r1 Cf11 tP i;ITv X 6 Y, XI, X vI? ~IH IUT Y T1J I I,1.I, AH'SL I r-(.I.I C I ));I f)4 '. ',- hi t64'-[I

SV 1 V =C -III.u __________ ~' LK'(II I1., J JJ) =AV A L! 11(1_ v-JJ,UJ dI: ~SVALUJL....-__ _________ IC) 1 ( 1 I' II 1 A 2A V1 x,, N T, A 1. * -f A T FX I..1 T - A N j PAJ (1iNI Y71 17 DI11 1i I KK ( 0S P )(*L. l 1 ) TI) 17 T, A P~S.T ( A, PI' T, XYJX, l IL,_____ _______-______ I S L I ~X, S 1)T Y T T I, I J, AMPSL T:J. I I*i (JcIlr T ( 53 1 C.KI= I\ C. I K HII ) I A,.K I I A KiP ~C e, i I ) \I \L U IE-= C!iPL X ( A MP RL, A ivP IF I"I...~.. \V ALIli- ( (II I IJJJ =V/Al(F- (II, JJJ )+,SVALIJIE 14 TiF(APS(AS',PT).1LF. 0.01)(1 Tn 17 (Ti Tr)? n 65

5uBiiUwn'Iluonicaii ________ ~'1i WPHIVI 'T1\1-I., I III'AT-( I 10 IO rHPB A 1 2,,? __ 1, I',Xf~, xm-,A TH,AMPT, AM~k IXI4 NX4 2 ~Y ~/,HI L JII ______ \I A " \"4( 3 )X04( 3 1 4k-At. XL o1 iL PA iI 1)A Il- +)lTT (X29?,Y,72,ANFX, X1,YL, 71) 0t P1) = I 1\'I'HX A A C.AL L P'Pk'I X I Y 1, 71,D L P7,DN) PP 4 9 11,3 -.. rA I I_ )jI A k(A d,X\ 1~IN~ A...., AA______________________________________ I F (K G I (K C1T. L5 1) ). T A 1 A,(). CA~. 121 ALHA -TA A.P,TA2L PHA, tL9KJVT21V) (A,.,A HR=S I RT GA;*.R 1. P ) 44*c?+( GR* R? I V xc 2 T____ P_ T(j *LF.?)(11 Tn I11 (P AF 12.)Ai 11 P PA P T G A b1AkPB C.AlLI P-l)FLC(ACl~tIFTA k?1. Pq k2iV,/ XL6? 1,X6, Xr\16 G A L I. CA XI A (- -, 2, XI?,2-,(,~k,(~I CA LL P '( LES~Y.LT.,i) A A I X AI. ]LA 1AL2IIA,- 1. Z29K I CK IF K IA (;T e 1 ) C,) T1) IO CLL XP21(ALPHA ZETA9kDLP,2Z 2V.T1V. DN.~..~ T A(KAI=S~K T ( (()(T71) )1-2+ If) T1\ )* CP!~- 1 =. A P H4A A ETA R2 P T2 P R21V-921 A, A,XL'.,Xi~, XH4? 1 '- C A, XI.' XIA ~i X XN ' I- I I.\? g CII I( T i Ilt- PITPLAXOYOZO,0P1LXitYi,71,KIC9,jj),,J ___ f _ I'iL (3) _________ ___________________________________________________.1 - 1. 1 ) 1).C ~fLL (L 3 03 =MMADE LIDEL-X5+2. 4X4lO-X2'~ZO"'ZO

7 7 + 9C' 2 ( I F //-I.L I L II. L' 2 * L. D E L ) G~. C A....,. T - /1 T >. 7 7. ) (I T 1( 4 IF _(,ILI * 1 l 1. ), (O T I I 6 3 7 7 7=7 71 3 1x )+(/7)c)L)/(3 Y1 =Y'+ ( 7/-O)~L(? JL. __ — (GH TH '5 __________ I 51 = C I -Xe):O 2 Y 1 -YO )*? 1-(:* TI S =S()R T ( 1)PS T) 'NX IR S1 Iiz 1t.O I Il'Il i F- P r1rM ( P X 1.,Y.I, 71,)D. ) DN I P 1( XI.':XI +Y1cY] II(1 =1/P 1 1k?CI =V I- 17) = -I~ 1)7~ / P (,2=\!FC.A CH (' 12 _ ________.ik=AI.AL Gi~G~: AM\=AL L*KIC- E.-A` L SL — ~ -C H - I ( )I-' L-.....-.-... —.... -... '-2 2 1 I. 1' 1CY 11 P- (I * I )(I T ( 2 _ _ _ _ 67 __ _ _ _ _ _

I (- S ( R A\(3) * (;1.:.I a)(;f T(I) 4 it- ( \(1.L0 *ANDs. AR(;2-.I.I.) UTf2. K 1 CK = 21Al PHA=/'kS TN'A16?) I AK A =fl ARW!f 1k i lJ 1 '1 J NF TR I~ 2 (A LPH A.?H TA- ~ h _ _ __.. ~.. ~. AlX~ =Al PH<)SA-M:T A )/(Sl F 1 ~=? / ( +(;x) (YAIXCl(.kFTA)IC( I S(AAp 0-~~ V 1 *-(y )/( +GY) SiIlK f~lhI I T \)i — TP;121(AL-PHARFATARPTPR\/,TV) Al1 FAI=AI.P HA-k1-.ITA ________(A Fd YL C fLSALPMAJJRCOS (-BTA(.Y=Artl-X*-C,1S (lFTA ) /CAS(AL-PHA) 1~2.4;/C, i.+;X T hl~I 1' PH1-:WrF1R ML 9fv\1.A.,F- XL MXIAXN.)___ ___ ) I~"SI f U' 3 ) M( 3 ),N (3 ),XL (3 ),XMC3 ),XN C3) till,~ '' I 1-' =Al ________________"1. I.1.~I t J..3. YCI )CHISC1 A.k T ):: T S I H(A L1A T )xcv( 1) Y 1.) A!3T I L1 / I )=YI C T - Yi4' )-l*T/ Y7C ("A.I-L-V —TC (X t X- l.- -- _ _ __ 8 __ _ _ __ _ _ __. - -- - — l_ _

SIAkOIH II 'F P ~"I.CT (A, i, A L )I HAPP) NV, 1-._,MN,(;A MlAq XL Xi N A3 X1 L IM(3) Nii4 hI S? I =Y L 3 Y1 - -YL( 1)=-C)OS(tIL.PH?)'*L( I I)+SIN(1A.LPHF2) *IM(:\' (1) =-'S I \(AL- PH? )' (T ) CD)S (AL PH? )*:IA( 01 I, 33 lI, V4~ x(i- = ( ~A'k P *Y () I KH)cR\/~cYN (I) /(AMA 1) 7' dN J-f~ 1_ ( 3) X (=C? -I C 3 I ),C? 7 = C /' / C1). z2) Y=Y2+XL ( X 3)*(7-? )S)TI)IST(X?,Y? Z?,9 1.OXYZ) (IH Ilii3...... ~ Al I (3 ) M 3 9l( 3 =1 1 ~( -L ( 2 );j2 T i ______3 C,d T GI I I T~ ~ Ri~~.~'r T, IL(1S. XXYFGIIJJSGA -—.I T T '. (T 0,) I *~Tr? 2 Y ~ Yf IG.(.)C TO 4 ___________69_______________

T T T~ C ( r.7;iC _______________'~>~ 4* ~ cc Y=AkS ( YY) Tc X 1. IX3 *AND, Y.LT*.Y4l~ Tl' (IIo ~ 10T~(X.L1. X?.ANII._ Y_.LT. Y?)G TO 't _ TF1:(X *c F. X? e ANr'. Y.c;F.e Y2)GO TO 3(0 ________ —Y 1,~yAN) Y_ F - Y 2(hO T 1 T f(~X CG[. X l )(;II Tf) 13 II( ) f 20 (Te Y (FIT n 1 1I l(X *I-F. X l *A ND&l Y. LF a Y3 )(() T(1 15 TIX * C, F XI)G O TO0 37 (;F1 Pi 20 (0 1fJJ= 0 6. -J 'i J N S1 kL~ (I I I I rF A PS IOT (A MP, X,9 Y, X L XM, S L TX t S L OT Y ITI J J IJ, AM1P SL vI hN~Ir(I\' S L (ITX (.4 _______Q — (.4 C IH-.F\1 /kRk/kG;LC(rgL., L f.1. A II _____ ' i =. r I 1,~ _2 4( = HH / 2. I ( I. ITS CXY Y,k I - T T, II 33 C(L ) i T~-( 33 ~ ) ThTO 3-.- - -—.- ---- - 7=~1 VI T IT,Ja1 ), 'fn I fSI; QY=~Hiy ( IT It k,3j) *C()S(CGL) — I X. f. [W(i A A\Mo. PY * L F a~l WA I TO P. 70

I *- H AC.="~ 1 V -C ) (XI. I. _____________ is________ (;.j I ~ ~ ')'~ H JH I Y IPI I M1; 1I ) p1. A StIAA ( Xl.,Y 1,7-1,,.M N1,Otp1.,tf.G. GG JI~. __ I AI1-PHAIAH TA q MJ 07,I)FTAAAtPATH,1O) P 1-AL I.,,"iNI )\ L. 7 7 1 =7. 1. C,~ A= (A F(7 1.I.T. 0 OR~ (.,G Tf) 99 ______ h-/..I _T. in. )GO TIt iF( 7 1 I-To 48. )(;G) To 3). 1 01=~1 2? 7 * 7-(fl.5)+O.OO 1? CoI F=Z o3CE o-5 I J-09 — 2. T Fno2:c n-, FooR/ - - 'ONI=1. 7* 1 47 0. FORh_____ G(l TI) 5 31 rCfli'=I=1.47*10,.Fns 5 C-f,, Wi \IF__________F).Z T f P ~F F:::CAN P:cT *rAN (Aly L~(L.PHA )*)-, 'i-I=Io-* 2 )::/\I P-HA )/S I( ALPHA) ~.=~i11: (rG:h 2 ) -(; A*f\ 2) /)N(3 '= ( A GIlTA):c;..( IkMY~~ D.'2 L:[ - Tf%:cA A AA+ 2 *X' (BXI3c —D-(. 1IJ.~Y1B-DN.L-2iiiDNA3..... C IA AI:A-(O (3:f~ (1 +N (2 ) 1\)N 2) /IN (3 ) 0T) 3 F (C4 *IT. T. C),()T) 99 (7 = (c24',)/ C r 7(C71 *r.r/?; T ( r1 C71 I. (.~,1.. C Y +I iP ( 2,c 1 /)1 (o 3) G P (;AL. Gi'8 N, ',IJ,(;, (~ ~AA ) 71

C_____AI AI(r. A|S Af(I F((I, l)N, AI,,A,AILPHA, IFTA r) K] CK ) (KIICK, GT. I)(;r TH 99 I o n..,..,..,.,.. I.f.,..... _ _. _...... _... ^......... ^...,._...... ^...... 'i lif '499q ~0 1 ( = _. ~ 9 9 -.I(J- - -.. -. ~. _...............- _. ^._...-^ _ _ -...-..-... -....-.. I Q= I i..... h L'd 379 LI. I1S PRN IT................ mmim..... ' -..,i,, m *..-... r l 72

OGIVA L PROGRAM._-Main-Program~ Sub.-Programs- rSubroutines._____ DIMENSION _L(3),M(3),N(3),DN(3),v 1 XL1(3),XM1(3),XN1(3)9 2 XL2(3),tXM2(3),tXN2(3), 3 X L3( 3) X M3(3),XN 3(3) __4XL(3,X43,X43) 5 X L5(3),XM5 (3)XN5(3), 6XL6(3),9XM6(3).,XN6( 3)____ DIMENSION AMPVAL(498) COMPLEX VALUE (4,98),SVALUE_______ COMPLEX VlV2,V3qV4,V12,V23qV34-,V41 REAL LMN REAL KCON. _____ COMMON VALUE/R/KCONMAXN/RR/ANEX____________________ COMMON 7RRR/BGL,CGL, L COMMON /RRRR/AAtCCtBOUTBIN_ ____ DATA P1/3.1415927/ DATA DELTA/1.0/,DELT/0.689/ _ ____ REAL SLOTX(4, 8)13.85,v3.859-3. 859-3.8593.6593.659-3.65,-3.659 1 _______2.3592.359-2.359-2.3592.1592.159-2. 159-2.15, 2 2.3592. 359-2. 359-2.3590.85,0.85,-O.85,-0.85, 3 _ 0.65,0.65,-0.o65,-'0.65,tO.85,0.8_5_,_-0.o85,P-0.85/ _ ----. ----. REAL SLOTY(4,8) /0.789-0.*7 8,-0.7 8,90.*7 8,2.e3 0 52.*3 05 -2.*3 05,2.*3 0 5 1 0.78,-0.78,-0.78,0.78,2.305,-2.3059-2.305,2.305, 2 3. 839,-3. 83,v-3.a83,p3. 83909.78p,r-0. 78, —O78,v0.78,v 3 2.3059-2. 3059-2.305,2.305,3.83,-3o83,-3.83,3o83/ C _C SAMPLE INPUT DATA C &INPUT ANGLEX=15.,ANGLEY=5.,OMEGA=10., C DSTEP=O.4, DA=5. 25,MMAX=13,NMAX=37,KPLSMA=0, ANEX=2. 5,9 C AA~=140.,PCC=49.,vBOUT=48.,vBIN=45.832 &END _ _ C _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ NAMELIS T /1 NPUT/ANGLEX, ANGLEYOMEGA9 1 DSTEPDAMMAXNMAXKPLSMAANEXt 2. AA,CC,9BOUT9B I N 1 READ(1,INPUT) 3003 FORMAT(lHl,'ANGLEX~',F6.2/'ANGLEY=',F6.2/'OMEGA=',F6.2/ ___ 'STEP OF DATA='ItF6. 2/ 'CONSTANTA='1,0F6.*2/ 2'M(MAX)='15~/ 'N(MAX)=',15/) WRITE (6,93003 )ANGLEX,ANGLEYOMEGA,DSTEP, DA,tMMAXNMAX __ WR I TE (6,v3333) AA,9CC,9BOUT9B I N 3333 FORMAT( 'A = 'F1O.2/'C = 'FlO.2/ -..1 'B(OUTER SURFACE)=',F1O.2/ ____ 2 'R(INNER SURFACE)=',F1O.2) IF(KPLSMA.EQ. 0)GO TO 15 WRI TE(6,r3004) 3004 FORMAT('IWITH PLASMA') GO TO 16 15 WRITE(6,93005) 3005_FORMAT('$NOPLASMA') 16 CONTINUE V l=(0. 90.) ---—.. ----_ _ V2=(0. 90.) V 3 = (0.,0. )....- - -. -- - - -. - - -_ __ _ _ _ -- --. - - -.. V4= (0.,90.) BGL =ANGL EX* PI /180. CGL=ANGLEY*PI /180. ROMEGA=OMEGA*PI /180.. o ~ L(1 )=SIN(BGL) L( 2)=-COS (BGL)*SIN( CGL)' L( 3) =COS (BGL) *C OS(CGL) _ _ _ _ _ _ _ _ _ _ _ _ _ _ 7 3 _ _ _ _ _ _ _ _ _ __ _ _ _ _ _ _ _ _ _ _ _

M( 1)=COS(BGL ).*COS( ROMEGA) M(2)=SIN(BGL)*COS(ROMEGA)*SIN(CGL)+SIN(ROMEGA)*COS(CGL) __M(3)=-SIN( BGL)*COS(ROMEGA)*COS(CGL)+SIN(ROMEGA)*SIN(CGL) W RI TE(6,p3001)L( 1) L( 2),L(3) __3001 FORMAT(3X,2HL=Fl0.3,3X,2HM=FlO.3,3X,2HN=FJO.3 WRI TE( 6,v 3002)M( 1),M( 2),M(3) 3002 FORMAT(2X_,3HL1=Fl10.3,?X,3HM1=FlO.3,2X,3HN1=F1O.3// CALL VECTCP(LMN) KCON=2.*PI /2.146 _ INDEX=l DO 111 1=1,4 _ _ _ DO 111 J=198 111 VALUE (IJ)=0.,90.) DEL TA=CC~-B OUT DELTB=CC-BIN _ ______ _ R A DI US=S QR T ( A A*A A +BOU T*BOUT )AA Z EN DI= C C -RAD IUS*S I N( BG L) _ ____ ZEND2=CC+RADIUS*SIN( BGL) _________SON2SO0RTC(I.L (3)L(3)) MM=2*-MMAX+1l ~_ _ _ _ _ _ _ _ _ _ _ _ _ __ _ _ _ _ _ _ _ _ _ _ _ _ _ _ DO 555 IIM=1,MM______ I M= I I M-MMAX-1 X C 1=1I M *D ST E P / SQN2 __ __ ____ DO 555 IIN=1,NMAX _ _ _ _ _ _ _ I N = I I N - i1_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ X C 2 =4 3. 0 +( D A-!I N *D S TE PS )SN 2 XO=-L( 1)*L(3)*XC2-XC1*L(2)___ YO=-.L(2)*L(3)*XC2+XC1*L(l1) ZO=-(L(1)*XO+L(2)*YO)/L(3)__ I F(ABS (L (1)) LE* 0. 001) GO TO 3 Z INT=ZO-XO*L (3) /L ( 1) GO TO 2 3 IF(ABS(L(2)) *LE* 0.001)GO TO 82___ ZINT=ZO-YO*L (3) /L( 2) 2 IF(ZINT *GT. DELTA *AND* ZINT.L.DLBG O9 IF(ZITNT *GT. DELTA)GO TO 80 ISIGN= 1 _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ ~~ G ~ O TO 81- __ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 82 IFCL(3i oEQ. 1.)GO TO 83 _ GO TO 99 83 IF(XO.EQ. 0. *AND* YO.EQ._0.)GO TO 99 I F(XO. GT. 0.)GO TO 87 ___ ~GO TO 80__ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 87 IISIGN,=1l' — GO TO 81 80 ISIGN=0 81 CONTINUE I F (I SI GN.EQ. 0) GO TO 6O Z____ IFIN =Z E ND02 GO TO 61 60 ZFIN=ZEND1 ___ 61 CONTINUE C FIND INTERSECTION POINT CALL PI NT ( XO,vYO, ZO, LXl1,V1,9Z1,9KIC I SI GN, 1,DLDZ) — I_ iF ( KTICL-.9GLT.AJ)-GO TO 99 IF(Z1.GE. ZEND2)GO TO 99 I F (ZI *G E. 4 8 o)GO TO0 9 9 C. COMPUTE DISTANCE PATH=DI ST (XOYOZ~qO,.Xjy1,Z1) C FIND NORMAL DIRECTION AT INTERSECT POINT I -------.74 ~

CALL PNORM(Xl,Yl,ZltDLDZ,DN) C FIND CONSTANTS A $ B AND GAMA CALL CONAB(M,N,DN,GA,GB,GAMA) __ __ C FIND ANGLES ALPHA AND BETA CALL ANGLES(L,DN,GAMA', ALPHA,BETA,0,KICK) _ - IF(KICK.GT. 1)GO TO 99 IF(KPLSMA.EQ. 1)GO TO 62 GO TO 63 62 CALL PLASMA(Xl,Yl,ZlL,MM,N,DNtGAGBGAMA, - 1 ALPHA,BETA,PATH,DLDZ,IQ) IF (IQ.EQ. O)GO TO 63 GO TO 99 C RAY WILL _PASS FROM MEDIUM 1 TO 2 ____ - C COMPUTE T12 AND R12 63 CALL TR12(ALPHA,BETA,R12P,T12P,Rl2V,T12V) __ GAMAB=SQRT ((GA*T12P)**2+ (GB*Tl2V)**2) AMP =GAMAB -. CALL NEWCOR(LM,N,GAGBXL1,XM1,XN1) C COMPUTE DIRECTION COSINES OF TRANSMITTED RAY CALL TRANSM(GA, GB, ALPHA, BETA, T12P T12V, 1 XL1,XM1,XNl,GAMAB,XL2, XM2,XN2) C RAY WILL PASS THE INNER SURFACE C COMPUTE INNER INTERSECT POINT CALL PINT(X1,Y1,ZItXL2,X2,Y2,Z2,KIC,0,O,DLDZ) IF(KIC.GT. 1)GO TO 99 ____ --- PATH=PATH+DIST(X,YlZ 1, ANEX,X2, Y2,Z2) PAT=PATH DD=PATH CALL PNORM(X2,Y2,Z2,DLDZ,DN) __ CALL CONAB(XM2,XN2,DN,GAtGB,GAMA) ________ CALL ANGLES(XL2,DNtGAMA,BETAALPHA l1KICK) IF(KICK.GT. 1)GO TO 99 C COMPUTE T21 AND R21 CALL TR21(ALPHA,BETA,R21P,T21P,R21VT21V) GAMAB=SORT((GA*T21P)**2+(GB*T21V)**2) GAMABB=SQRT((GA*R21P)**2+(GB*R21V)**2) ____ AMPT=AMP*GAMAB AMPR=AMP*GAMABB TAMP=AMPT. _ CALL NEWCOR(XL2,XM2,XN2,GA,GB,XL3,XM3,XN3) C COMPUTE DIRECTION COSINES OF TRANSMITTED RAY CALL 1RANSM(GA,GB,BETA,ALPHA,T21P,T21V,XL3,XM3,XN3, 1 GAMAB,XL4,XM4,XN4)____________________ JJ=1 CALL RECORD(L,XL4,X2,Y2tZ2,X3,Y3,Z3,DDT,K) IF(K.EO. 1)GO TO 30 IF(Z3.LT. ZFIN)GO TO 31 GO TO 32 ____ 30 DDT=DDT+DD _____________________________________________________ * CALL APSLOT(AMPT,X3,Y3,XL4,XM4, 1 SLOTX,SLOTY, I I I,JJJ,AMPSL) IF(JJJ.EQ. 9)GO TO 51 AKCON=KCON*DDT CSKD=COS (AKCON) SNKD=SIN(AKCON) AMPRL=AMPSL*CSKD AMP I M=AMPSL*SNKD SVALUE=CMPLX ( AMPRLAMPIM) VALUE(III,JJJ)=VALUE( I I,JJJ)+SVALUE 51 CONTINUE GO TO 32 31 CALL SECOND(XL4,XM4.XN4,X29Y2,Z2, INDEXt 1 DFLTA.DELTR.PAT.TAMP. I. T T N.7FNnl 7 F:Nn7

32 CALL REFLCT(GAGBBETAR21PR2lVXL3,XM3, XN3,GAMABB, 40 JJ=JJ+l CALL REPEAT (INDEXDELTA, EL TB, X'2 V2 2XL 5,X5_,N 5,TH__ 1 AMPTAMPRXL4,XM4,pXN4,pX3,vY3, Z3DDKL,3) R AM P=AM PR TAMP =AMPT _ PAT=PATH 'IF(K.EQ. 1)GO TO 41 ________F(K EQ. 0)GO TO 99 IF(Z3.LT. ZFIN)GO TO 44 -GO TO 42 - ~ - —. — - - I.. 4-4 CALL SECOND( XL4,XM4,XN4,X2,v2,Z2, INDEX, DELTADELTBPATTAMPL, 1 ISIGNZEND1,ZEND2)__ GO TO 42 _ _ _ 4 1 DD T=DD+ PATH _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ CtA L-L — A-P-S —L-OTT —AM PT,9X 3,Y3,9X L4,9X M4, 1 SLOTX,9SLOTY,9I IIJJJ,9AMPSL) IF(JJJ.EQ. 9)GO TO 91 A KGON =KC ON* DDT SNKD=SIN( AKCON) CSKD=COS( AKCON) AMPR L=AMPSL*CSKD AMPIM=AMPSL*SNKD _ _ S VAL UE= C MP L X( AMP R LA MP IM) 91VALUE(IIIJJJ)=VALUE(III,_JJJ)+SVALUE _ 91CONTINUE 42_ IF(ABS(TAMP) *LE. 0.01 )GO TO 99_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ IF(JJ.LT. 20)GO TO 40 99 CONTINUE 555 CONTINUE WR ITE (6 6,001) 9001 FORMAT(1H1,2X,4HSLOT,3X, llX,4HREAL,6X,9HIMAGINARV, ___ 1 6X99HAMPLITUDE/) DO 432 1=1,4 DO 432 J=198 432 AMPVAL(IJ)=CABS(VALUE(1vJ)) DO 431 1=1,4 DO 431 J=198 WRI TE(6,t9000)I,tJVALUE(1IJ),AMPVAL (I J) 9000 FORMAT(215,3E15.5) 431 CONTINUE DO 433 1=1,8 V1=Vl+VALUE( 1,1) V2=V24-VALUE( 2, 1) ___ V3=V3+VALUE ( 3,j IJ_____________________________ V4=V4+VALUE( 4, I 433 CONTINUE V1 A=CABS (Vi) V2A=CABS (V2) V3A=CABS (V3) __ _ _V4A=CABS (V4)__ _ _ __ _ _ _ _ _ _ _ _ _ _ _ _ _ _ ___ _ _ _ _ _ PH11=ATAN2(AIMAG(V1),REAL(V1)) PHI2=ATAN2(AIMAG(V2),REAL(V2))_ _ PH I3=A TAN2 (AlIMAG (V3),?REAL (V3)) PHI4 L=AT AN21( A TM AGr-IvIL4)-,RE AII %I. V12=VI+V2 _ _ _ V? 3 = V 2 + V 3__ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ V 34=V3+V4 V41 =V4+VL V12A=CABS(V12) V?3A=CABS (V23) V34A=CABS (V34) - _ _ _ ~7 6 _ _ _

___ ~V41A=CABS (V41)________________________ PH112=ATAN2(AIMAG(Vl2),REAL(Vl2)) _PHI 23=ATAN2(AIMAG(V23),vREAL (V23)) PH134=ATAN2(AIMAG(V34),REAL(V34)) PHI141 =ATAN2 (A I MAG (V4 1),R EA L (V4 1) WRI TE(6,P9004) 9004 FORMAT(1X//lOX,'(VlV2,V3,V4)t) WRI TE(6,9006) 9006 FORMAT(l11X,4HREAL,6X,9HIMAGINARY,5X,1OHABSv VALUE,1OX,5HPHASE) WRITE(6,,9002) V1,VIAPHI1, 1 V 2,9V 2A,9P H!2, — 2 V3,?V3APHI3, 3 V4,V4APHI4 _ _____________ WR ITE (6,t9 005) 9005 FORMAT(lX//lOX,'(Vl2,V23,V34,V41)')_ WRI TE(6,r9002)V12,Vl2APHI12, I V23,vV23APHI23,. 2 V 349V 3 4A,9P H 134, _____ 3 _V41,V41APHI41__________________ VVV=1./4./PI VV1=VVV/0.8698*( PHI34-PHI12) VV2=VVV/0. 8810* (PHI41-PHI23) ERROR V=ARSTN (VVI) ERROR.H=ARS IN( VV2) -_______WR I TE (6v9 003) ERRORV,ERRORH 9003 FORMAT.('***','ERRORV=',El2.5,5X,'ERRORH=',El2.5) 9002 FORMAT(4E15.5) GO TO 1 END 77

SUBROUTINE SECOND( XLXMXNvX2,Y2,Z2, INDEX, 1 DELTADEL-TBPATHBAMPLIlSIGNZEND1,ZEND2, 2SLn TX,'SLOTY) D IMENS ION SL OTX (49,8),9SLOTY (4 v8) DIMENSION DN(3),XL(3)tXM(3),XN(3), lXL1(3)tXMl(3),XNl(3)., 2XL2( 3) vXM2 (3) XN2( 3), __-_ 3XL3( 3) vXM3( 3),XN3( 3) 4XL4( 3) IXM4(3),XN4( 3) 5,XL7( 3),XM7( 3),XN7( 3) DIMENSION L(3) COMPL EX VAL UE (4,v8),SVALUE REAL L__ __ ___ _ _ REAL KCON COMMON VALUE/R/KCONMAXN/RR/ANEX COMMON /RRRR/AAC,CCBOUTBIN JJ=1 IF(ISIGN *E~o 0)GO TO 60 ZFIN=ZEND1J_ GO TO 61 60 ZFIN=ZEND2 61 CONTINUE CALL PINT( X2,Y2,Z2,XLX5,Y5,Z5,KIC, 1,0,DLDZ)_ IF(KIC.GT. 1)GO TO 17 PATH=PATH+DIST(X2,Y2,Pz2,1.,9X5Y5,Z5) CALL PNORM( X5,Y5,Z5,DLDZDN) DO 10 1=1,3 10 D)N(T)=-DN(I) CALL CONAB(XMXND-NGAGBGAMA) CALL NEWCOR(XLXMXNGAGBXLI,XM1,XN1) CALL ANGLES( XLDNGAMAALPHA,BETA,0,KICK) IF(KICK.GT. 1)GO TO 17 CALL 1R12(ALPHABETAR12PT12PR12VT12V) GAMABB=SQRT((GA*Rl2P)**2+(GB*Rl2V)**2) GAMAB=SQRT(C(GA*T12P) **2+ (GB*T12V )**2) AMP T=BAMP*GAMAB AM PR=BAMP*GAMAaBB ___ _ CALL REFLCT(GAGBALPHARl2PR12VXL1,XM1,XN1, 1 GAMABB,XL2,XM2,XN2) CAL L R ECORD ( L,XL 2,9X5,IY5,pZ 5,pX6,rY6,pZ 6,9DDtK) IF(K.EQ. 1)GO TO 11 IF(Z6 *LT. ZFIN)GO TO 30 30GO TO 12 _ _ __ 30CALL THIRD(XL2,XM2tXN2,X5,Y5,Z5, INDEXiDELTADELTB, 1 PA TH, AMPRiL,I SI GNZ END1,tZ END2, SLOTXvSLOTY) - GO TO 12 11 DD=OD+PATH CALL APSLOT(CAMPTX6,Y6,XL2,XM2, 1 SLOTX,9SLOTY I I I,9JJJtAMPSL) IF(JJJ.EQ. 9)GOTO5 AKCON=KC ON*DD CSKD=COS(CAKC ON) SNKD=S IN(CAKCON) AMPRL=AMPSL*CSKD AMPI M=AMPSL*SNKD SVALUJE=CMPLX(CAMPRL,v'AMPITM) VALUE ( I II 1,JJJ ) =VALUE ( II I,pJJJ )+SVALJE. 51 CONTINUE 12 CALL TRANSM(GAGBALPHABETATl2PT12VXL1,XM11,XN1,p I GAMABXL3,XM3,XN3) 20 JJ=JJ+1l__ CALL REPEAT(CINDEX,DELTADELTB, X5vY5,Z5,XL3,XM3,XN3,qPATHt 1 AMPTAMPRXL7,tXM7,tXN7,oX7,pY7,tZ7,PDDKKL,-JJ)_ 78 I

IFIKK.EQ.v 1)GO TO 13 IF(KK.EQO. O)GO TO 17. IF(Z7 *LT* ZFIN)GO TO 31 GO TO 14___ 31 CALL THIRD(XL7,XM7,XN7,X5,Y5,Z5, INDEX,DELTADELTB, 1 PATHAMPTL, ISIGNZEND1,ZEND2,SLOTXSLOTY) GO TO 14 13 DD=DD+PATH IF(ABS(AMPT).LT. 0.01 )GO TO 17 CALL APSL0T(AMPTX7,Y7,XL7,XM7, _____ __ 1 SLOTX,9SL OTY, I II,9JJJtAMPSL) IF(JJJ.EQ. 9)GO TO 91 A KCON =KGCON*DD --- —-- CSKD=COS (AKC ON) SNKD=SIN(AKCON) AMPR L=AMPSL*CSKD__ _ __ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ AMP IM=AMPSL*SNKD SVA L UE=C MP L X( AM PR LAMP IM) VALUE( III,JJJ) =VALUJE( III,JJJ)+SVALUE 91 CONTINUE 14 IF(ABS(AMPT).LE. 0.01 )GO TO 17 IF(JJ *GE. 20)GO TO 17 _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ GO TO 20 17 RETURN END SUBROUTINE THIRD( XLXM, XNX2,Y2,Z2, INDEX, 1 DELTA,DELTB,PATHBAMPL, ISIGNZENDlZEND2, 2SLO TX, SLOTY) _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ DIMENSION SLOTX(4,8),9SLOTY(4,8) DIMENSION DN(3),XL(3),XM(3),XN( 3),. LXL-1(3),XM1(3),XNl( 3), 2X L2( 3)v,XM2 (3), XN2 (3),9 3XL3( 3),XM3( 3),XN3( 3)t 4X L4 ( 3),qX M4 (3),XN4(3) ___ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 5,9X L7 (3),9X M7(3),XN7 (3) DIMENSION L(3) COMPL EX VAL UE (4,8), SVALUE REAL L REAL KCON COMMON VALUE/R/KCONMAXN/RR/ANEX COMMON /RRRR/AACCBOUTBIN JJ=1 IF( ISIGN.EO. 0)GO TO 40 ZFIN=Z END? GO TO 41 40 ZFIN=ZEND1 _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 41 CONTINUE CALLPNTX,22XLXY5ZI,,,Lz IF(KIC.GTo 1)GO TO 17 PATH=PATH+DIST(X2,PY2,pZ2,1l.X5,PY5,Z5)CALL PNORM(X5,Y5,Z5,DLDZDN) 0 0 1 0 1= 1, 3_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ 10 ON (I)=DN (I) CALL CONAB (XM, XN, ONGA, GBGAMA)__ CALL NEWCOR(XLXMXNGAGBXL1,XMltXN1) CALL ANGLES-,(xlnNr,GAMAAL PHABErTA,0,KICK) _ _ _ IF(KICK.GT. 1)GO TO 17 CALL iRi?2(ALPHA, BETA, R12P, T12PR12VTI2V) G'AMABB=SORT((GA*Rl2P)**2+(GB*Rl2V)**2) GAMAB=SQRT( (GA*T12P)**2+(GB*T12V)**2) AMP T=BAMP*GAMAB AMPR=BAMP*GAMABB CALL REFLCT( GAGBALPHARl2P,Rl2VXL1,XM1,XNl,

1 GAMABBXL2,XM2,XN2) CALL RECORD( UXL2,X5,V5,Z5,X6,Y6,Z6,DDK) TF (K * EQ. 1 )GO TO 1 1 TF(Z6.LT. ZFIN)GO TO 30 GO 10 12 30 CALL FOURTH( XL2,XM2,XN2,X5,Y5,Z5, INDEXqDELTAvDELTBPATH, 1 AMPRL, ISIGNZEND1,ZEND2,SLOTXSLOTY)___ _____ __ GO TO 12 11 DD=DD+PATH CALL APSL0T(AMPTX6qY6,XL2tXM2, 1 SLOTXSL0TYIIIJJJAMPSL) IF(JJJ E~o. 9)GO TO 91 AKCON=KC ON*DD CSKD=Cns C AKC ON) SNKD=SI NCAKCON) AMPRL=AMPSL*CSKD AMP IM=AMPSL*SNKD SVALUE=CMPLX( AMPRLAMPIM) - -- VALUE(IiIJJJ)=VALUE(II1,JJJ)+SVALUJE_________________ 12CALL TRANSMCGAGBALPH-ABETAT12PT12VXL1 X.M1.,XNI1, 1 GAMABgXL3,XM3,XN3) 20 JJ=JJ+1 CALL REPEAT( INDEXDELTADELTB, X5,Y5vZ5,XL3vXM3,XN3,PATH, 1 AMPTAMPRXL7,vXM7,tXN7,X7,Y7,pZ7,DDKKLJJ) _____ IF(KK.EQ. 1)GO TO 13 IFCKK.E~o O)GO TO 17 IF(Z7.LTo ZFIN)GO TO 31 GO TO 14 31 CALL FOURTH(XL7,XM7,XN7,X5,Y5,Z5, INDEXPDELTADELTBPATH, 1 AMPTLISIGNZEND1,ZEND2,tSLOTXSLOTY) _____________ GO TO 14 13 DD=DD+PATH IF(AB.S(AMPT).LT* 0.01 )GO TO 17 CALL APSLOT(AMPTX7,Y7,XL7,XM7, 1 SLOTXSLOTY, I II,JJJAMPSL) I F(JJJ. EQ* 9 )GO- TO 51 _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ A KCON =KC ON* DD CSKD=COS (AKC ON) SNKD=SINC AKOON) AMPR L=AMPSL*CSKD AMP IM=AMPSL*SNKD SVALU)E=CMPLX (AMPRLiAMPIM)____ _____________________ VALUE ( I II, JJJ) =VALUEC I II, JJJ )+SVALUE 51 CONTINUE 14 IF(ABSCAMPT).LE* 0.01 )GO TO 17 IF (JJ.GE. 20) GO TO 17 GO TO 20 17 RE TURN _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ END a _ _ / 80 __ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _

SUBROUTINE FOURTH( XLXMXNX2,tY2,PZ2, INDEX, 1 DELTADELTF3,PATHBAMPL, ISIGNqZEND1,ZEND2, 2SLO TX, SLOTY) DIMENSION SLOTX(4,8),SLOTY(4,8) DIMENSION DN(3),XL(3),XM(3),XN(3), 1XL1( 3),9XM1C3) 9XNI( 3), 2XL2( 3),XM2(3),XN2(3), I_ 3 X L3 3) X M 3 3),XN3 3) 4X L4( 3) X M4(3),XN 4 3) 5tXL7( 3),XM7( 3),XN7( 3) DIMENSION L(3) COMPLEX VALUEC(498),9SVALUE REAL L REAL KCON COMMON VALUE/R/KCONMAXN/RR/ANEX COMMON /RRkR/AACCBOUTqBIN JJ=1 IF( ISIGN.E~o 0)GO TO 60 ZFIN=ZEND1__ ___ _ _ _ _ _ _ _ _ _ _ _ _ _ GO TO 61 60 ZFIN=ZEND? 61 CONTINUE CALL PINT( X2,Y2,Z2,XLX5,V5,Z5,KIC,1,0,DLDZ) IF(KIC oGT. 1)GO TO 17 PA TH=PATH+DIST( X2,rY2,PZ2,.t aX5vY5,Z5)_________ CALL PNORM( X5,Y5,Z5,DLDZDN) DO 10 1=1,3 10 DN (I )=-DN (I) CALL CONAB ( XMtXN, DN, GAGB, GAMA) CALL NEWCOR( XLXMXNGAGBXL1,XMlXN1) CALL ANGLES( XLDNGAMAALPHABETA,0,PKICK)_______ _ IF(KICK.G7T. 1)GO TO 17 CALL 1R12(ALPHABETARl2PT12PR12VT12V) GAMABB=SQRT(CGA*Rl2P)**2+(GB*Rl2V)**2) GAMAB=SQRT((GA*Tl2P)**2+(GB*Tl2V)**2) AMP T=BAMP*GAMAB AMPR=BAMP*GAMABB- _ ____ CALL REFLCT(CGAGBALPHAR12PRl2V;XL1,XM1,XNl, 1 GAMABBXL2,XM2,XN2) CALL RECORD( LXL2,X5,Y5,15,X6,Y6,Z6,DDK) IF(K.E~o 1)GO TO 11 IF(Z6.LT. ZFPIN)GO, TO 30 GO TO 1?2____ 30 CALL FIFTH( XL2,XM2,XN2,X5,Y5,Z5, INDEX9DELTAvDELTB, 1 PATHAMPRL, ISIGNZENDlZEND2,SLOTXSLOTY). GO TO 12 11 DD=DD+PATH CALL APSLOT( AMPTX6,Y6,XL2,XM2, 1 SLOTXSLOTY,1IIIJJJAMPSL) ___ _______________ _ IF(JJJ.EQ. 9)GO TO 51 A KCON=KC ON*DD C SK D= COS C AKO ON) SNKD=SIN(CAKCON) AMPRL=AMPSL*CSKD AMPI M=AMPSL*SNKD__ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ __ _ _ SVALUEj=C-MPLX (AMPRL,AMP TM) VALIJEC 111, JJJ )=VALLJE (III, JJJ )+SVALUJE 51 CONTINUE 12 CALL lRANSM(GA, GB, ALPHA, BETA qT12PT12V, XLl,)XM1,iXN1,v 1 GAMABXL3,XM3,XN3) 20 JJ=jj+1l_ _ _ _ _ _ _ _ _ _ _ _ _ _ CALL REPEAT(INDEXDELTADELTB, X5,Y5,159XL3,XM3,XN3,PATH, 1 AMPTAMPRXL7,tXM7,XN7,,X7,V7,Z71DDFKKLJJ)

IF(KK.E~o J)GO TO 13 IF(KK.EO. O)G.() TO 17 lF(Z7.LT., 7FIN)GO TO 31 GO0 TO 14 31 CALL FlFTH(XL7,XM7,XN7,X5?Y5vZ5, INDEXDELTADELTB, 1 PATHAMPTL, ISIGNZEND1,ZEND2,SLOTXSLOTY) GO TO 14 13 DD=DD+PATH IF(ABS(AMPT).LT. 0.01 )GO TO 17 CALL APSLOT(CAMPT, X7,vY7, XL7,XM7, __ __ ___ 1 SLOTX,9SLOTYqI IItJJJtAMPSL) IF(JJJ.EQ. 9)GO TO 9.1 A KGON =KC ON* DD CSKD=COS (AKCON) SNKD=SINC AKCON) AMPR L=AMPSL*C SKD _ AMP IM=AMPSL*SNKD S VA LUJE = CM PL X( AMP R L,AMPIM) VALLJE(IlIIJJJ)=VALUE(II1,JJJ)+SVALUE 91 CONTINUE 14 IF(ABS(AMPT).LE. 0.01 )GO TO 17 IF(JJ.GE. 20)GO TO 17 __ ______ __ _____ GO TO 2 0 17 REIURN END SUBROUTI NE FT FTH( XLXM XN, X2,Y2,Z2, INDEX, 1 DELTADELTBPATHBAMPL, ISIGNZEND1,ZEND2, 2SLOTXqSLOTY) ___ _ _ _ _ _ _ _ DIMENSION SLOTX(4,8),SLOTY(4,8) DIMENSION DN(3),XL(3),XM(3),XN(3;), 1XL1(3),XM1(3),XNI(3), 2 XL2(3),XM2(3) X N2(3) 3 XL3 3),XM( 3) X N3 3), 4 XL4 ( 3),XM4 (3 ),XN4 ( 3) _ _ _ _ _ _ _ _ _ __ _ _ _ 5,X L7 ( 3),XM 7C 3) XN 7( 3) DIMENSION L(3) COMPLEX VALUEC(4,98),SVALUJE REAL L REAL KCON COMMON VALUE/R/KCONMAXN/RR/ANEX __ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ COMMON /RRRR/AACCBOUTBIN JJ=1 CALL PINT (X2,Y2,12, XLX5,Y5,Z5,KIC,1,0,DLDZ) IF(KIC *GT. 1)GO TO 17 PA TH=PA TH+D IST (X2,9Y2,tZ 2,91., X5, Y5,tZ 5) CALL PNORM(X5lY5,Z5,DLDZDN)__ _ _ _ _ _ _ _ _ _ _ _ _ _ _ DO 10 1=1,3 10 DN(I)=-DNCI) CALL CONAB( XMXNDNGAGBGAMA) CA LL NEWCOR (XL,9XMvXNIGA, GBiXLL19xM1,9xN1) CALL ANGLES(CXLDNGAMAALPHABETA,0,KICK) IF(KICK.GT. 1)GO TO 17 _ _ _ _ _ _ _ _ _ _ _ _ _ _ CALL lRi?2(ALPHABETAR12PTl2PR).2VT12V) GAMABB=SQRT( (GA*R12P)**2+(GB*Rl2V)**2)GAMAB=SQRT( (GA*Tl2P)**2+(GB*T12V)**2) AMP T=BAMP*GAMAB AMPR=BAMP*GAMARB TF(AMPR.LE* 0.01 )GO TO 17 _ _ _ _ _ _ _ _ __ _ _ _ _ _ _ _ _ _ CALL REFLCT CGAGBALPHAR12PR12VXL1,tXM1,XN1,l 1 GAMABBgXL2,XM2,XN2) CALL RECORD( LXL2,X5,Y5,15,X6,Y6,16,DDK) IFCK.E~o 1)GO TO 11 GO TO 12 -- - /82

11 DD=DD+PATH CALL APSLOT (AMPT,9X6,9Y6,9XL2,vXM2,1 SLOTXSLOTY, IglJJJAMPSL) IF(JJJ.EO. 9)GO TO 91 AKCON=KC ON*DD CSKD=COS (AKC ON) SNKD=SI N (AKCON) AMPRL =AMPSL*CSKD AMPI M=AMPSL"SNKD SVALUJE=CMPLX( AMPRLAMPIM) VALUJE(IIIvJJJ)=VALUJE(IIItJJJ)+SVALUJE 9)1 CONTI NUE 12 CALL IRANSM( GA, GBALPHA, BETA, Tl2PT12V, XL1, XM1,XN1, _ 1 GAMABXL3,XM39XN3') 20 JJ=JJ+1 CALL REPEAT( INDEXDELTADELTB, X5,pY5,Z5,XL39XM3,XN3,PATH,9 1 AMP TAMPRXL7,XM7,PXN7,tX7,PY7,PZ7,PDDKKLJJ) IF(KK.EC). 1)GO TO 13 TF(KK.EO. O)GO TO 17__ _ _ _ _ _ __ _ _ _ _ _ GOnITO 14 13 DD=DD+PATH IF(ABS(AMPT).LT. 0.01 )GO) TO 17 CALL APSLOT (AMPT,9X7, Y7,XL7,9XM7,9 1 SLOTX,9SLOTY,9 I I IJJJ,9AMPSL) TF(JJJ oE~o 9)GO TO 51 ____ _ A KCON =KC ON* DD CSKD=COS ( AKCON) SNKD=SIN( AKCON) AMPRL=AMPSL*CSKD AMP IM=AMPSL*SNKD SVALI-JE=CMPLX( AMPRLAMPIM) ____ _ VALU)E(IIIPJJJ)=VALUJE(III,tJJJ)+SVALUE 51 CON TINUE 14 IF(ABS(AMPT) *LE* 0.01 )GO, TO 17 IF(JJ *GE* 10)GO TO 17 GO TO 20 17 RETURN _ _ _ _ _ _ _ _ __ _ END 83 -

SUBROUTINE REPEAT(INDEXDELTADELTB, X29Y29Z2 1,XL5,XM5,XN5, PATHAMPTAMPRXL4,XM4,XN4, 2 X3 Y3, Z3qDD, K, LJJ) DIMENSION DN(3),L(3), IXL2( 3),XM2(3),XN2(3), 2XL3(3)gXM3(3)tXN3(3), 3XL4( 3),XM4(3),XN4( 3), 4XL5 (3),XM5 (3), XN5( 3),0 5XL6(3),XM6C3),XN6(3) REAL L COMMON /RR/ANEX COMMON /RRRR/AACCBOUJTBIN CALL P INT (X2 Y2,9Z 2, XL5,9X1,9YlZ 1,KIC,11,DLDZ) TF(KIC.GT. 1)GO TO 10 PATH=PATH+DIST(X2,PY2,?Z2,ANEXXlYlZl) CALL PNORM ( X 1,vY 1,Z'1,DLDZ,9DN) DO 48 1=1,3 48 DN(I)=-DN(I) __CALL CONAB (XM5, XN5,DNGAGBGAMA)_________________ CALL NEWCOR( XL5,XM5,XN5,GAGBXL6,XM6,XN6) CALL ANGLES(XL5,DNGAMABETA,ALPHA,1,KICK) IF(KICK.GT. 1)GO TO 10 CALL TR21 (ALPHABETAR2lPT21PR21VT21V,) GAMABB=SQRT((GA*R21P)**2+(GB*R21V)**2) IF(JJ.LE. 2)GYO TO 11 _ _ _ _ _ __ _ _ _ _ _ AMP=AMPR*GAMABB GO TO 12 11 AMP=AMPT*GAMABB 12 CONTINUE CALL REFLCT (GAGB, BETAR2lP,R2lVXL6,XM6,XN6, 1 GAMABBXL?,gXM2,,.XN2) ___ CALL PINT(XlY1,ZlXL2,X2,Y2,Z2,KIC,0,ODLDZ) IF(KIC.GT. I)GO TO 10 PA1TH= PA TH+D I ST ( X 1,pY1,PZ 1 i,ANEXp,X2,pY2r,Z2) CALL PNORM( X2, Y2, Z2,DLDZ,DN) CALL CONAB( XM2,XN2,DNGAGBGAMA) CALL ANGLES( XL2,DNGAMA, BETAALPHA, 1,KICK)_____________ IF(KICK.GT. 1)GO TO 10 CALL TR21 (ALPHA, BETAR?1PT2lPR2IVT21V) GAMA=SORT ((GA*T21P) **2+( GB*T21V) **2) GAMAB=SQRT((GA*R21P)**2+(GB,*R21V)**2) AMP T=AMP*GAMA AMPR=AMP*GAMAB ______ __ __ CALL NEWCOR(XL2,XM2,XN2,tGAGBXL3,pXM3,tXN3) CALL IRANSM(GAGBBETAALPHAT21PT21V,XL3,PXM3,XN3-,1 GAMAXL4,XM4,XN4) CALL REFLCT(GAGBBETAR21PR2lVXL3,XM3, XN3, 1 GAMABXL5tXM5,XN5) CALL RECORD(LXL4,oX2,tY2, Z2,PX3,Y3,PZ3,DDK)____ R ETURN 10 K=0 R E TUR N END SUBROUTINE PNORM(X1,V1,vZltDLDZ,DN) DIMENSION DN(3)__ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ PL=SORT( 1+DLDZ*DLDZ) P2=SQR T( X1*X1+Y1*Yl) DN( 1) =X1/P1/P2 DN( 2)=Y1/Pl/P2 DN(3)=-DLDZ/Pl RE TURN _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ END SUBROUTINE CONAB( MN,DN A, BPGAMA) 84

DIMENSION M(3),N(3),DN(3) REAL MN Gl=VECTDT (DNM) G2=VEC TDT(DONN) G3=SORT( Gl*Gl+G2*G2) G7AMA=-l */G3 A.=-G7AMA*Gl B=-GAMA*'G2 'RETURN END _ SUBROUTINE ANGLES(LDNGAMAALPHABETAlI IKICK) DI MENSInON L (3),ON (3) COMMON /RR/ANEX REAL L P1=3.1415927 ARG1 =-VECTDT (LDN)' ARG2=-1 */GAMA IF(ABS(ARG2).GT. 1)GO TO 10 GD TO 11 10 KICK=2 GO TO 23 1 1 I E lh.EO. 1)GO TO 2__ _ _ _ _ _ _ _ _ _ _ __ _ _ _ _ ARG3=ARG2 /ANEX GO TO 3 2 ARG3=ARG2*ANEX I F(ABS (ARG3). GE. le)GO TO 4 GO TO 3 4 KICK=2 ___ GO TO 23 3 CONTINUE TECARG1 *GE* 0. *AND* ARG2.GE.-D-O.-) GO TO 21 KI CK=2 GO TO 23 21 ALPHA=ARSIN(ARG2) _____ _ BETA=ARSIN( ARG3) KI CK=0 23 RETURN END SUBROUTINE TRL2( ALPHABETARPTPRVTV) COMMON /RR/ANEX _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ ALBT=ALPHA-BETA GX=ANEX*COS (ALPHA) /COS(CBETA) RP=(l1.-GX)/( 1.+GX) TP=2. / Cl+GX) GY=ANEX*COS C BETA) /COS( ALPHA) RV=( l1.-GY)/( l1.+GY) _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ IV=2./(l1.+GY) RE TURN FND SUBROUTINE TR21(ALPHABETARPTPRVTV) COMMON /RR/ANEX ALB T=ALPHA-BETA__ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ G X=AN E X*C OS ( AL PHA) /COS(BE TA) GY=ANEX*COS (BETA) /COS( ALPHA)RP = -( 1.,- G X )/( l.o+ G X T P =2,*G X/ 1+ G X RV=-(l1.-GY)/(l1.+GY) TV=2o*GY/(1lo+GY) _ _ _ __ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ RE TURN END SUBROUTINE NEWCOR( LMNABXLXMXN) DIMENSION L(3),M(3),PN(3),XL(3),PXM(3),PXN(,3)-.. REAL LMtN,I85

DO 30 1=1,3 XL (1 =L (I) XM( I )=A*M( I ) +B*N ( I) 30 XN(I)=-B*M(I)+A*N(I) RE TURN END SUBROUTINE TR ANSM(A,BALPHABETATPTVLMNGAMAXLXMI,.XN) DIMENSION L(3),M(3),N(3),XL(3),XM(3),XN(3) DIMENSION YL(3) p,YM( 3),YN(3) REAL LMN ALB T=ALPHA-BETA DO 32 1=1,3 YL(I)=COS(ALBT)*L(I)-SIN(ALBT)*M(I)_ YM( I)=SIN(ALBT)*L( I)+COS(ALBT)*M( I) 32 YN(I)=N(I) DO 33 1=1,3 XL (I)=YL (I) 33 XM(I)=(A*TP*YM(I)-B*TV*YN(I) )/GAMA CALL VECTCP(XLXMXN)___ RETURN END SUBROUTINE REFLCT( ABALPHARPRVPLMNGAMAXLXMXN) DIMENSION L (3),M( 3),NC 3)XL(3),XM(3),XN(3)4 DIMENSION YL (3),YM( 3),YN( 3) REAL LMN ALPH2=ALPHA*2. 1)0 3? 1=1,P3 YL( I) =-COS( ALPH2)#L( I) +SIN( ALPH2 )*M( I) YM ( I ) =-S IN( ALPH2)~L( I )-COS (ALPH2)*M(I) 3? YN(I)=N(I) DO 33 I=1,3 — _ - - - - -- - - _ XL (I) =YL (I) 33 XM(I)=(-A*RP*YM(I)-B*RV*YN(l))/GAMA CALL VECTCP(XLXMXN) RE TURN END SUBROUTINE RECORD( L,XLX2,Y2,Z2,XYZDDK) ________ D I ME N SION L (3),L (3) REAL L C,1=L(1)*XL(1)/XL(3)+L(2)*XL(2)/XL(3) C2=C1*72 C3=-L Cl)*,X2-L( 2)*Y2+L( 3)*43* C 4 =C 2+ C3__ _ _ __ _ _ _ __ _ _ _ _ _ _ _ __ _ _ _ C5=C1+L( 3) Z=C4/C5 X=X?+XL (1)/XL (3) *(CZ-Z2) Y=Y2+XL (2) /XL (3 )*( Z-Z2) C.C=X*X+Y*Y+ (Z-43. ) **2 DD=DIST(X2,Y?,Z2,.pl 0,XYZ) __ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ IF(CC.LE* 30.25)GO TO 2 K=3 GO TO 3 2 K=1 3 RElURN E N D_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ SIJBROtJIINE VECTCP(LvMqN) DIMENSION L (3), M( 3),N(3) REAL LMN N(1) =L (2 )*M( 3)-L (3) *M( 2) -- ---------- - 86~

N(3)=L( I)*M(2)-L(2)*M( 1 RETUJRN ENn FUNCTION VECTOT(LM) DIMENSION L (3) M( 3),N( 3) REAL LMN DO 1 1=1,3 1 (I )=L( I ) M I VFCTDT=N(1)+N(2)+N(3) RETURN END SUBROUTINE SLOTS(XXYYB IF(XX.GT, 0. ) GO TO 2 IF(YVY GT. 0.)GO TO 3 111I=3 GO TO 5 3 111=4 G,1O TO 5 2 F (YY.GT. 0.)GO TO 4 111=2 GO TO 5 4 III=l 5 CONTINUE CCC=COS ( BGL) DDD=COS(CGL) X3=4. 5*CCC X2=3.0*CCC Xl=l. *CCC Y1=1.5*DDD Y2=3. O*DDD Y3=4o5*DDD X=ABS(XX) Y=ABS(YY) IF(X.LT* X3.AND. V.LT GO TO 30 10 IF(X.LT* X2 *AND* V.LT IF( X GE* X2.AND* V_ *GE IF(X *GE* X2.AND* Y.LE IF(X. GT. Xl)GO TO 13 JIJJ=8 GO TO 20 13 JJJ=5 GO TO 20 12 TF(Y.GT. Y1)GO T0 14 GL, IIItJJJCGL) - — -- --- ----- --— c` --- —-L -— * Y3)GO TO 10 * V2)GO TO 11 *_Y2)GO TO 30 * Y2)GO TO 12 J JJ=1 GO TO 20 14 JJJ=2 GO TO 20 11 TF(X * LE. IF(X.GE. IF(X *GE. j1JJ=7 GO0 TO 20 17 JJJ==3 GO TO 20 16 JJJ=4 GO TO 20 Xl.AND* V_ LE*_Y1)GO TO 15 Xl *AND* V *GE* Y1)GO TO 16 X1)GO TO 17 - - 15 20 30,I J J = 6 RE TURN JJJ=9 - -- e87

R ETURN END SUB~ROUTINE APSLOT(AMPXYXL,XMSLOTXSLOTYIIIJJJAMPSL) DIMENSION SLOTX(4,98),SLOTY(4,8),XL(3),vXM(3) DIMENSION L(3).CnmmnN /RRR/BGL,CGLL RFAL L P1=3. 1415927 H =I. 0 6 W=1.060 HH=H*'COS (BGL) WW=W*COS(CGL) Hn2=HH/2. wn?=WW/2. C A LL S LOT S ( X,9 Y BGLII I, J JJ CGL) IF(JJJ.EO. 9)GO TO 3 CX=SLnTX (III, JJJ ) *COS(BfGL) CY=SLO1Y(11II JJJ)*COS(CGL) PX=AB~S ( X-CX)__ __ PY=ABS ( Y-CY) IF(PX.LE* H02 *AND, PY.LE. W02-)GO TO 2 JJJ=9 (;O TO 3 2 AMPSL =AMP*COS (PY*P I/WW) *COS (PX*P I/HH) AC1=VECTDT(XLPL-) _ _ _ B1=L(2)*L(2)+L(3)*L(3) B2=SOR T( BI) IA3=L (2 )*XL (2) +L (3) *XL (3)R4=L(2)*XM(2)+L(3)*XM( 3) AC?=(XM(1)*B3+XL( 1)*B4)/B2 AMPS L=AMPSL*ACI1*AC2 __ __ _ 3 REIURN END -~88

I I.- I - I r - f I'? I.F -, - - i I - v If 7 7' fr 'J-i. 7 r -r If'I' l S~I MKP 11 TJ i\i- PI '\T (XO, YO, 7 0,91)1,, XI J Y 1, 7 1It K I C- JJJ T II-?,.) LH7 1)1 (- \I S I I N' I) L )-) I1171=/fl- C C 4 '7 4 I I T I H) I (C TII 170 1)) 11 10 K=0 I 1 A A HRA A AA+K C,1I=0.KA A/1)1.(3)/L(3 C.?= (X.+CC,+:C.C4)/AA 0 (- AA ) 1=B*K + CC*CCCA1 3) =2.? 1. *'C3+ C A ll. 7 P (1IN T ( 7-0,D), CCZ,9K JKIC) T(K IC,.(;T. 0)17W 1T7 99 IH( KJ *FO 0)17I TO 99 7 7 =CC7 (I)__ ___ __ IF (O).LT. 0. )I I TH 6 0 (I.=S(.1 'I O) -AA XX=XO+(I Z-Z0)~:clnL(I)/ML( YY=YO+ (77-70 n:'i;2 /L (3) (OF!. =0A-XY IF(-~ABSOHJF.LT, 0.0)l G01T50 '0 K K=K K + TIF( KK *FO, 2)171 TO 70 (;() TH -) 70 71 =AMAX I (T7(1),T7(2)) 7 S=Afr HivIi (17(Tl), 17( 2) 7 1=7 108 PO 71=7 L K 11C7= 0 YY=- 71-CC) f I )7 YY/XX -.. i~ F 1 1kR 1" 09 KI C= 2 F-K) RI~(11 1\ E PA SPi A( X IY -P,71,Z vL r MNO! DNC, AC,,C, A mAA ALP HA RFT AIPATH, OLOZI1,) 01'FiSI(NL (3),M (3nN 3,Oi3) D( 5) I~i'-115101"CC7 (4), T7 (4) -'HAL L,F89

C, i.llh /R!-'/AA,CC,.OOT, kTiN VA'A 1 l/ C l. F-il. f 1 H-1.? / I t- (71.LT. 0 1. 48.. 7 4R. )&() To 99 I F( 1.I T. 1~. )C'iG T I 1 I- (7 1. LT. 1 5. )(1 1TO'? I (/71. LT. 48. )(;l T 3 ' 1' 1 - 7.. 1,: (0.5)+-. 001 F (7.1. -.,. b. ). ) T( 11. C( N I =72. T ' ) 1 0. F O) I') l ) )bn 1 Nix 1=2.?3:o. HFl 8n (; ) )..........27 FF=0. 01.27 1-O.nO ClONuI =1. 1 47:1.n. F n C 10n n 3 TFF=O. 02? 4' Z 7 1 ': '(O.R)-0. 0(958 IF(71.LT..);. ) T't 31 Ci,*I =.-);; 1 n. FOR CIfA.(1 31 Cf I =1. 1.47 0. FOR S C. 'llI1 F I I F 0F LI =Cl f i T F -C)i'C, I" T A ' ( A l PHA ) /C() S ( ALPH A ) L ':"2 I I- =-) FL.: C, (2.) S. ALPHA) /SIN( ALPHA) PA 1H=PA TH+I) l.. ___)___ P 1 =: iF 1 ': ( (;G:: M ( 2 )-GA^'N ( 2 )) /I)N ( 3 ) p,2=i)L1 ( C, i ( ] )-GA,:N( 1 ) ) /1)\ ( 3) Xk =X +PI1 VY=Y1 -P? A [An:IA -I 11 CcC=cC-z1 ________ A A F = A A,:: A A + IAk ' [kt A C, =-A A R + CCC C CCC Cn =(). S/O)N,(3 ) /ON( 3) /AA C2=(r),[(1 ):Xf+n)(2(2)':YR)/)N(3)-CCC C 2= C 2 / A A C, ' = X P X H + Y ':YB -A A A A + A K C C.,=C3/?./AA f( h ) =C 1:C 1. I)( 4 ) =. 2 C 1 =: C? )(3) =?. C1:-C,+C:C+1. i(? )=.::C? C3-?. CCC i( 1 ) =C':C3+A __ _ C A I. 7 P NT (7-, )D,CZ,KJ,KIC) IF(KIC, CT, r))(; T.I 99 TF (K J.F0. n)GC Tn 99 < K =0 / 7 =CC7 (1 ______ "'= AA -^ ( 77-CC ) 2 T (0. T. nT. )C) TO 60 X=X' +), ( 1 /N ( 3)((7O)-A71 Y v = Y P + () ' ( 1 ) / I) J ( 3 ) Z ( Z 7 7 z ), Yv=YF.i,+I\( 2) /[)N(3 )'*( ZZ-Z71) X Y = RSR 1 ( X X XX -YY:yYY) FI:Il=.-X Y IFr(ARS(I)FL).LT. 0.01)G0 T(- 50 KK KK+= I ) I/ ( KK ) =C.C (I ) - - - - -... - - 1 90

I F- ( KK * `0 2 ) G 1' TA C) / S~=AfA~ N\) (T/ (I T/ ( ) 1F(ZS *(T. 1. )GII TI) 9 0 TF(71 G *T. 7 1 )GH T 0 91__ (ii 7 =/L. 90 /8=/ S 8 1 Xl1 = XK+DN( I) ID()(7-1 Y 1 =YR+0I)N ( 2 )1/\ (3)7-i)_ _ X X=SOk 1 ( AAKK- ( Il- CC) 2 YY= ( 71. -C C) LI!)) =YY /X X C~A! I-I PNf~ HRM (X I, Y 1, 7 I 1 F)L07,OIN CA LL C, NA R (N, N, 1P\ GN (A, Gk,_9G A MA CALI-L A NCLF S(LIDr1,PA MiA,9A L PHAKRF TAO0 I I'(K ICK.(T. I)CI TI) 99 T F 1 kfT'iIST K PF I I)R\t FIv T ARf)I IH I 1)1S F lT ( XY,70, KNFX,Xl Y1JK C, 0 1 1I ( XIN -X) (:2 l -YO F ) HU T?4 (71-0) I F 17 4)CC H C1 I ) I I I -- N~iS I (I F) \ T 7 (4' 1F(S),R U R(),O TF(IFk.Fo. O)GOi TI) 1?? CIN 11 Ii IFI I t4 C 7 (T)f*__ I1F(A FSk MU))T I (1))LT, O.01) Ml TO 2 1 21? 2 4 22 K=K+l C7 () =PHOTR C 2 A F(A S (kR lU)ITJ3))L T. 0O0 1)C4) CATA 25 24 K=K+ I C 7 (K) = PkI(IITR (3 (' I1A 27 25~K + C Y ( K) =oT k ( 4) TO 2 TO 26 - --- ---- 9 1

27 K\=IK=1 II- ( K F ) l T Cl 99? =1. KK?M '17 (1T) =7 0I+C/. (IT:) II-(KK * VA. 2 )Gil TA 3() 1 F (KK *, 4)GP TO 40 30 KJI=0 GT* DO GO TO 3 1 3I KJ=K.I+1 rC.7 (K<J ) =T 7 ( 1. 3? IF 1~7(2) (IT, DO)6( TI) 33 33 KJ =Kj+)I.Cfl7 ( K T.T7.( 3 4 IF(KJ.F C 0)601 T0.99C C(T ) 1 48 40 KI=fl0 I( 17-( ).GT. DDflGF TOl 4 (4) I0 f? 4 1 Ki =J+lI ( t-C Z ( K.T,) T 7Tr 4 4? iF(17(2) oT i. ) c)n T(C) 43 C;rin I 44 4.3 KJ =K.)+ + C7 (K.I =T7(2) 44 To(T.M OT-,4 (;I) 1I I (., 41) Kc,I=K J +l C.C7 ( KJ ) =TZ ( 4AS I1- II(4).GT. [)I) )G O TO 47 G I In -its 47 Kj =KJ+I UC7(Kl) =TZ(4) K I C.=n I -li C N' '4', rTIC=2 '-I -- -- -.........-.. _.I- _...-. -....-......... - L. — --— ~- -- ~ ~ -# 4. --- SUBROUTINE -POLRT -- - - - ---- -- This subroutine computes real and complex in IBM Scientific Subroutine Package (SSP). Programmer's Manual, Form 120-0205-n. 221 LINFS PRII\TED roots of a real polynomial and is available Ref. sstem/36 0 Scientific Subroutine Package,. -- -- -- ---— ' --- —— — --- 92