A fast fourier transform subroutine for the Honeywell X16-series computers

(1)

1975

.

v.

LT,

NAVAL SHIP RESEARCH AND DEVELOPMENT CENTER

Bethesda, Maryland 20034

A FAST FOURIER TRANSFORM SUBROUTINE FOR THE HONEYWELL X16-SERIES COMPUTERS

by

Lorraine R. Minor

APPROVED FOR PUBLIC RELEASE: DISTRIBUTION UNLIMITED

COMPUTATION AND MATHEMATICS DEPARTMENT RESEARCH AND DEVELOPMENT REPORT

February 1974 Report 4352

AUG.

(2)

The Naval Ship Research and Development Center is a U. S. Navy center for laboratory effort directed at achieving improved sea and air vehicles. It was formed in March 1967 by merging the David Taylor Model Basin at Carderock, Maryland with the Marine Engineering Laboratory at Annapolis, Maryland.

Naval Ship Research and Development Center

Bethesda, Md. 20034

*REPORT ORIGINATOR

MAJOR NSRDC ORGANIZATIONAL COMPONENTS

OFFICER-IN-CHARGE CARDE ROCK 05 SHIP PERFORMANCE DEPARTMENT 15 STRUCTURES DEPARTMENT 17 SHIP ACOUSTICS DEPARTMENT 19 MATERIALS DEPARTMENT 28 NSR DC COMMANDER 00 TECHNICAL DIRECTOR 01 OFFICER-IN-CHARGE ANNAPOLIS 04 AVIATION AND SURFACE EFFECTS DEPARTMENT 16 COMPUTATION AND MATHEMATICS DEPARTMENT 18 PROPULSION AND AUXILIARY SYSTEMS DEPARTMENT 27 CENTRAL INSTRUMENTATION DEPARTMENT 29 SYSTEMS DEVELOPMENT DEPARTMENT 11

j

(3)

UNCLASSIFIED

SECURITY CLASSIFICATION OF THIS PAGE (When Data Entered)

D, FORM 1,73

REPORT DOCUMENTATION PAGE _{BEFORE COMPLETING FORM}READ INSTRUCTIONS

1 REPORT NUMBER Report 4352

2. GOVT ACCESSION NO. 3 RECIPIENT'S CATALOG NUMBER

4 TITLE (and Subtitle)

A FAST FOURIER TRANSFORM SUBROUTINE FOR THE HONEYWELL X16-SERIES COMPUTERS

5 TYPE OF REPORT& _{PERIOD COVERED}

Final

6 PERFORMING ORG. REPORT NUMBER

7. AUTHOR(a)

Lorraine R. Minor

5 CONTRACT OR GRANT NUMBER(.)

9 PERFORMING ORGANIZATION NAME AND ADDRESS

Naval Ship Research and Development Center

Bethesda, Md. 20084

10 _{PROGRAM ELEMENT, PROJECT, TASK}

AREA a WORK UNIT NUMBERS

P.E.:

62718N; P.T.:SF11121106 W.U.: 1-1931-001

IICONTROLLING OFFICE NAME AND ADDRESS 12 REPORT DATE

Fphruary 1974 13 NUMBER OF PAGES

1

43

14 MONITORING AGENCY NAME &ADDRESS(/1 different from Controlling Office) _{13. SECURITY CLASS. (of (hi. roport)}

UNCLASSIFIED

15a. DECLASSIFICATION/DOWNGRADING SCHEDULE

16 DISTRIBUTION STATEMENT (of Milo Report)

APPROVED FOR PUBLIC RELEASE: DISTRIBUTION UNLIMITED

17 DISTRIBUTION STATEMENT (of the abstract entered in Block TO, If different from Report)

16 SUPPLEMENTARY NOTES

19 KEY WORDS (Continua on revra side If ncmar) end identify by block number)

Fast Fourier Transform, Computer program, Signal Analysis, Fortran _Program, Honeywell X16-Series Computer

20 ABSTRACT (Continua on reverse side If necersear) and Identify by block number)

This report describes a Fast Fourier Transform subroutine _{as programmed for} use with a Honeywell X16-Series computer, the DDP-516 mini-computer. This

Fast Fourier Transform, coded originally in Fortran, provides a computerized method for transforming digitized data from the time domain to the frequency

domain. When used with the DDP-516, which as a cycle time of .96 microseconds, the Fast Fourier Transform subroutine allows the computation of a transformed array of 512 complex values in .3 second.

(4)

(5)

TABLE OF CONTENTS

Page

ABSTRACT 1

ADMINISTRATIVE INFORMATION INTRODUCTION

THE FAST FOURIER TRANSFORM ALGORITHM

THEORY 1

EXAMPLE APPLICATION 5

111. THE FFT ALGORITHM AS PROGRAMMED IN FORTRAN

STEPS IN THE CODING 7

GENERATION OF SINE AND COSINE VALUES

FORTRAN PROGRAM LISTING 15

IV. THE FFT ALGORITHM AS PROGRAMMED FOR THE HONEYWELL X16-SERIES COMPUTER 18

LIMITATIONS IMPOSED BY THE HONEYWELL SYSTEM I

PROCEDURE FOR CALLING THE FFT SUBROUTINE 18

STEPS IN THE CODING 19

GENERATION OF SINE VALUES 21

BIT INVERSION ROUTINE FOR THE CASE N = 512 23

FFT AND INVT PROGRAM LISTINGS, HONEYWELL X16-SERIES COMPUTERS 24

FLOWCHARTS 35

REFERENCES 41

LIST OF FIGURES

Figure 1 The A and B Registers after Multiplication 20

Figure 2 The A and B Registers after Shift 21

Figure 3 - Generation of Sine Values . . R. . . . . . . .

...

18 . . . . -,.. . .. 14

(6)

ABSTRACT

This report describes a Fast Fourier Transform subroutine as programmed for use with a Honeywell X16-series computer, the DDP-516 mini-computer. This Fast Fourier Transform, coded originally in Fortran, provides a computerized method for transforming digitized data from the time domain to the frequency domain. When used with the DDP-516, which has a cycle time of .96 microseconds, the Fast Fourier Transform subroutine allows the computation of a transformed array of 512 complex values in .3 second.

ADMINISTRATIVE INFORMATION

The work reported was authorized and funded by the Naval Ship System Command under Program Task SF 11121106, Program Element 62718N, and performed under Work Unit 1-1931-001.

I. INTRODUCTION

In the past several years, a series of computer programs has been developed and used by the Navy in analyzing echoes from sonar transmissions. Following the perfection of these techniques, an investigation was made into the feasibility of performing such an analysis in a real-time manner aboard ship. A Honeywell DDP-516 mini-computer was procured for the study and computer programs for real-time analysis using that

machine were developed.

An important part of the analysis is the transformation of digitized data from the time domain to the frequency domain. This transformation has been accomplished using a method known as the Fourier Series, or Fourier Transform. Depending upon the analysis desired, the Fourier Transform process sometimes must be run more than once for each sonar transmission. It is therefore important that the algorithm representing the Fourier Transform and the coding of the algorithm be extremely efficient. The theory of the particular algorithm implemented on the Honeywell computer was developed and implemented as a Fortran program by

James Beard.1 The programming of the Fast Fourier Transform for use on the Honeywell DDP -516 mini-computer in connection with the shipboard real-time analysis was carried out at the Naval Ship Research and Development Center and is described within the pages of this report.

Section II of this report describes the theory underlying the algorithm; Section III discusses the Fortran coded version of the algorithm; and Section IV describes the assembly language version developed for the Honeywell X16-series computers.

(7)

II. THE FAST FOURIER TRANSFORM ALGORITHM

THEORY

The Fourier Series can be defined

N-1 1

A(n) = -

xoonk

.\71.71 _k=0 where N = 27r = [0, N-1]

Using a straight forward method of expanding this expression, N2 multiplications are required to compute the transform. Many algorithms have been designed to reduce the number of multiplications required to the order of N log, N, since such a reduction results in a reduction in computer time needed for the calculation. An algorithm which results in such a reduction is a Fast Fourier Transform, or FFT. To calculate the Fourier Series, Equation (1), the values for k and n can be represented by the following binary notation:

m-1

k = E kP 2"

p=0 M-1

n = E n /q

-q=0

where k is the binary value (0 or 1) of k for the p th binary place.

For example: M = 3 k = 5

k=k +2k +4k2

₀ 1 510 = (101)2 k = 1 + 2(0) + 4(1) k = 5 M-1 M-1

nk = E E n

k P211+P p=0 q=0

Instead of summing from p = 0. M-1 over p, we will sum from p = 1, M over M-p.

When P 1 M-p = M-1 P M M-p = 0 (1) WN .(2) (4) =

(8)

M m-I nk

= EE

nq km_p 2M+q-p p=1 q=0 M p-1 M M-1

=EE

_{nq kM-p}

2M+q-P + EE

_{nq kM-p} 2M+q-P p=1 q=0 p=1 q=p

The first sum in Equation (5) is over values of p and q such that M +

q-p < M or q < p.

The second sum contains all terms for which M + q-p M or q p.

M p-1 M M-1

EE

nq kM-p

2M+q-p + EE

nq kM-p 2M+q-p

wk = w p=1 q=0

p=1 q=p Let q-p = r If q p, then r 0 w12M+q-p 2M+r = w =

w2Mr

"N = wN2r . 2rrN2r -1

=e

=1 Equation (6) becomes M p=1

EEn k

_{q M-p} ')M+cl-p \\ink = w p=1 q=0

If X(k) and A(n) are co-ordered arrays of dimension M with subscripts k _{and nq where k}

P , P'

nq = [0 or 11, and the summation in Equation (1) is made from 0 to 1 for each kp, then Equation (1)

becomes

1 1 i 1

A(nm_i, nsi-2, . . ., n1, no)

=E E... E

E

X(km_, , km_2,

k0=0 k1=0 kM-2 =0 kM-1 =0 M p-1 (8)

EE

nq kM-p 2M+q-p , ko) WN p=1 q=0 (5) N q q

(9)

This can be rewritten M-1 s-1 1 k0

E

nq2q

1k 2M-s

E n 2q

E w m_s

q

,

A(nm , nm_2, . . .

, ni,

no) = E w

_-N 0,0

q0

k0=0 kM-s =0

1k

2M-I_no

E w

_"N X_{(kM-1' kM-2'} _{1 'k 0'}

k =0

M-1

Now suppose we let

X0 = X0 (kM-1, km_2, . . ., k2, 1(1, 1(0) 1 k 2M-1n =X1 (no, km_2, ,k2, ki , ko)

= E

wNm--1 0 kM-1=° s-1 M-s

okNm, 2

n 2q

= (no, n1, . .

-1'

n

s1'

_kM-s-

"

_q'₌₀ q _Xs-1 k_M-s=0 M-1

ko E nq 2q

= Xm= Xm (no, n1, . . . , WN

q0

E

xm_i (no, n1, . .

'M-2'

1(0) 1(0=0

But A(nm, , nm-2, , no) from Equation (9) is equal to

Xm(no, n1, . . . , nm) from Equation (10).

This X is the Fourier Series desired, except that the bit pattern for the address of the elements is reversed.

To avoid this undesirable ordering, the element positions in the original representation (Equation (10)) are reversed so that the correct bit pattern for the Xm addresses of the elements is achieved (Equation (11)). (Note that both the symbol X in Equation (10) and the symbol Y in Equation (11) refer to the same array.)

.

(10)

Since 1 _2k1

E

nci2g Y

E

w

q(ko,

, no) k1 =0,

Ya =-Y0 (ko, k1,..

km_.1) k 2m-1 n

= Y1 (ko, k1,

, 2, no) '=

E

WNM-1 '

Y (

k 0.k 0, 1' "" k M-2"km-i) k =0 0. M-I 2M-s

E

n 2g

1 km-s q=0

Ys (ko, k..

, no)

= E

_{Yll tk}

_s-

_0'

kM-s=0 km:Ls, nsL2 . no),

M-1 ko

E

n 2g

.

= _"I'M no)

=. E.

g=° _(ko,

WN .

k =0₀

and Ym , n1,, no) = A(n) for the Fourier transform.

EXAMPLE APPLICATION

As an example of the application of the Fourier Transform algorithm, let us suppose that M = 3 'so that N = 2M = S. Using the conventional .method of computing the Fourier Series and the expansion of WN pro-vided on page 2, the following relationship is achieved

A(3) = X(0)w0 + X(1)w3 +.X(2)w6 + X(3)w9 + X(4)w12 X(5)w15 + X(6)w18 + X(7)w21

= X(0) t X(1)w3 + X(2)j + X(3)w1 X(4) - X(5)w3 - X(6)j - X(7)w1 Using the FFT algorithm and working backwards, the following, notation Is derivved 2

A(3) = Y3(3) -='Y3(01 1)

2 k

E

n 2g

1 0 q Y3(011)

= E w

Y2(k011) =0, k0 = Y2(011) + Y2(111)w(1+2+0)

= Y(010+ Y2(l11)w3

= . (11) . . n1, (nm_i, . + - - - -(ko, nl, no)

(11)

Y3(011) -= Y1(004) + Y1(01 1)w2(1+2) + w3 1Y1(10,1) + Y1(1 T )w20 +211

,=& Y1(001) + Y1(01 1)w6,+ .Y1(101)w3 + Y1(1 1.1)W9

Yi(001) + jY1(011) w3Y1(101)+ wlY1(111)

1 4k, no Since ,

E w

Yo (ko, Icy k2) Y1 (kW kl nO) k =0₂

Y3(011) = Y0(000) + Y0(001)W4 + [VOW) + Y0(01 1)w4 11 + [Y0(100) + (101)81.141wa

(Y0(1 10) + Y0(1141)w4 liwl

= Y0(000) + Y0(001)w4 + fY0(010), + jY0(011)w4 + Y0(100)w + Y00,01)W7 + Y0(110)w1 + Y (111)w5

Making the final bit inversion, Yo, (ko, k1,, k2) = X (k2, k1 ko)

Y3(01 1) = X(000) - X(1 oo) + j:X(01 .Q) - jX(1 10) + w3 Woo + w7X(1,01) wl X(aa 1) + w5 X(111)

= X(0) + X(1)w3 + X(2)j + X(3)w - X(4) + X(5)w7 - X(6)j + X(7)w5 = X(0) + ,X(1)w3 + X(2)j + X(3)w1 - X(4) , X(5)w3' X(6)j, = N(7)wl = + +

-

-- -

(12)

-Expansion of W, for N = 8: . 2rr rr WN = = e 0

_{e0 =l}

. rr -1

=e

= cos - j sin 4.7r = WN1 271. 2 _{= e-11 4} = cos 77 - j sin

-= cos j sin

-

-j

37r

T

37 37r = e _{= cos} _{-j sin} _-j 2 2

8arrbrr

brr

- T

T

w 8a+b = w N8a wNb = c

ej

= e-j 8a+b = w b

III. THE FFT ALGORITHM AS PROGRAMMED IN FORTRAN

STEPS IN THE CODING

The computation of the Fourier Series of a complex array can be divided into three distinct steps. No complex variables are used in this program. All real values are stored in contiguous locations with all the imaginary values following, again in contiguous locations. Complex computations are handled directly by the coding. The first step performs the bit pattern inversion and the computation of the values for the Y1 array. The second step computes the values for the Y2 array. The third step computes the values for the arrays Y3 through Yrvi, and divides the elements of the Ym array by the square root of N. Each step is explained separately in the pages that follow.

w 4 -N

. 477

3-] 4

57r

= cos Tr -j sin rr = -1 = -WNo

W N5 = 4 = cos 51r -j sin s rr= = _n,, 1 4 2 . 67r WN6 -J

=e

_{= cos} 3rr _{-j sin T = j = -wN}3rr 2 wN7 J- 4

=e

77r 7 7 = cos -4- -j sin 4

y j

_w 3 w 7-.1 e N

(13)

A complete listing of the Fortran coding is provided on pages 15-17. In the following discussion, the first location of an array has an index of zero. Since Fortran does not allow a zero index, the indexes used in the coding are incremented by I.

Step One

This portion of the coding performs the bit inversion on the input array X0 = X0 (km_i , . . , ko) and then computes the Y1 array using the following equation:

(ko, , , km no) =

E

w_2, kM-1 `10

kl,

. .

, km,)

2m-i

(12)

k =0

M-1

This coding is complex, since the computation of the value and the bit reversal of the address of the value are performed in one operation and the results (the Y/ array) occupy the same locations as the input array Xo.

The value in the variable 11 controls the number of times this step is repeated. At each iteration, II is incremented by 2. The bit inverse oft! is determined and stored in 12. If the values of 11 and 12 are equal, control passes to Statement 13 in the coding. Since M = 9

and 11 = index (even number) (Example: 010 000 010)

12 = inverse of 11 (Example: 010 000 010) let J1 = 11 + N/2

13 = inverse of .11 11 + 1

J3 = .11 + I =13 + N/2 J3 is its own inverse.

If X is the input array, and Y is the array after reordering for bit inverse, then

Y(11) = X(11)

Y(J3)=X(J3)

Y(J1)= X(13) Y(13) = X(JI)

Note that for II and JI, no = kM-1 = 0

and for 13 and 13, no = km-1 = I

Now Y is the Yo in Equation (12). The values for Y1 may now be determined as follows: Y (11) = Y(I1)w0 + Y(13)0

= X(I1) + X(J1)

Y1(13) . Y(II)w00 Y(I3)w2566

= X(11) - X(J1) (Example: 110 000 010) (Example: 010 000 011) (Example: 110 000 Or k) (ko, . , 1 = +

(14)

-Y1 (J1) = X(13) + X(J3)

Y (J3) = X(I3) - X(J3)

When 12 is greater than 11, control passes to Statement 14 in the coding. In this case

and the following values can be computed: Y (J2) = X(14) + X(J4) Y (J3) = X(14) - X(J4) Y1 (11) = X(12) + X(J2) Y1(14) = X(12) - X(J2) Y1 (J1) = X(13) + X(J3) Y (J4) = X(13) - X(J3) Y (12) = X(I1) + X(J1) Y1(13) = X(II) - X(J1)

When 11 is greater than the inverse 12, the computation will already have been performed and no further coding is required.

Step One is repeated until 11 is equal to N/2. Since at this point all values greater than N/2 will have been computed, control will pass to Step Two.

Step Two

The second step of the coding is by far the simplest. The operations compute the Y2 array:

11 = index (even number) (Example: 000 000 010)

12 = inverse of Ii (Example: 010 000 000) 13 = 12 + 1 (Example: 010 000 001) 14 = II + 1 (Example: 000 000 011) JI = 11 + N/2 = inverse of 13 (Example: 100 000 010) J2 = 12 + N/2 = inverse of14 (Example: 110 000 000) J3 = 13 + N/2 (Example: 110 000 001) and J4 = 14 + N/2 = inverse ofJ3 Y(II) = X(12) (Example: 100 000 011) Y(12) = X(I1) Y(13) = X(J1) Y(14) = X(J2) Y(J1) = X(13) Y(J2) = X(14) Y(J3) = X(J4) Y(J4) = X(J3)

(15)

-1,1 = index

= 11 + .11 = i + 2

J2. = l + 3

ho = 0, n1 = 0 Y2(11) = Y1(111) + Y1(J1)

= L,. n= G Y2(12) = Y1(12) + Y (J2).w`,

no = 0, ni 4'11 Y2(J1) ='Y1(1.1) + Yi(J1)w2M-2(2)

n₀ n Y2(J2)='Yt(12).+Y102)w2M-2(3k 2kir2M-2 wk/M-2' =

eM

2 . = 2)k 72M-2-M . krr

=e

= cos (ir k) -j sin k)

-2 k = 0 wik2m = k = 1 k2M-2

_{=_j}

k2M-2 k k. = 3 wk 2M-2 = Y2(11) = Y101) + Y1(J1),

Y2(1.2) = .Y/(12) -j Y1(J2) Y2(J1) Y1(11)

Y1(n)

Y2(J2) = Y1(12) + j Y1(.12)

The real or imaginary values for Y2(1t) and Y2(J1) are the sum or the difference of the real or imaginary values in the Yt array. For the Y2(12) and Y2(.12), computations, however, this 'is not the case. For example,

when no =1, nt 7 of, (Example:, 000 GOO 0001), (Example: Q00000 001) (Example: 000

we

,1340) (Example:- 000 goo 011) 1 Y(k(1, k r, . k , 1.10) =

E

M-2 M-2

E

_n92q q =0 _k no) k_M-2=0 1 1 Ifwe let For for = = 1 =-1 = 2 =

-k 1 12 I

(16)

let Y(I2) + j YM(12) be the 12 complex element of the array

and Y(J2) + j YM(J2) be the J2 complex element of the array.

Then (Y2(12) +j YM2(12)) 1Y1(12) +j YM1(12)] - j [Y1(J2) +j YM1(12)1 = [Y1(12) + YMI (J2)1 + j [-Y, (J2) + YM (12)1 By the same process

(Y2(12) + j Y2(j2)) = [Y1(12) - YM, (J2)1 + j [Y(J2) + YM, (12)1

The Fortran coding handles this operation in a straightforward manner. By incrementing the Index 11 by 4 at each iteration until N-1 has been reached, all elements of the array will be processed.

Step Three

The third step handles the computation of the values for the arrays Y3 through Ym where

s-1

1

km_s2MsE

n 42g

Ys (ko, ki , ...km_s_i, ns_1, ... no) =

E wN

q=0

Ys- 1 (k0 '

km_s,ns_2 ... no)

(13) km-, =0

KJH2 = KJH+ 2

Step Three is divided into two separate loops, the first to accommodate the following special case. Suppose rig = 0 for q = 0, , S 3.

and s takes on the values 3, 4, .

.. N.

In the coding, 13 is the variable representing s in the Equation (13). The other variables are defined as follows:

KI1 = 2M-13 KJJ = 213-1 NN

=2's

KJH = 213-2 JM = KJH+ 1 . . _s, .

(17)

s-1 2M-s

E

n 2g

Then for n5-2 = 0 and ns-1 = 0 _WN q=0 w 2 (0)

M-s s-1 2M-s

E

n 2g 2M-s2s-2 1 and ns_i. = 0 _WN q=0 WN = WN2M-2 _

-s-1 2M-s

E

n 2g 2M-s2s-1 ns-2 = 0 and ns-1 = 1 WN q=0 WN ,M-1 s-1 2M-s

E

n 2q 2M-s(2s-2 + 2s-1) ,M-2 2M-1 = WN- WN- = ns-2 = 1 and ns-1 = 1 WN q=0 WN

For this case, the coding parallels that of Step Two (remembering that S=13) in which 11 is incremented by 2s at each iteration of the loop

12 = 11 + J1 = 11

J2 = 11 + 's-2 + 2s-1

The second loop performs the remainder of the computations necessary to compute values for the Ys

array. Note that in the first loop, bits s-2 and s-1 are used, while in this loop, the bits 0, s-3, s are

used. Consequently Loop 2 is executed s-2 times.

We define

11 = 1

J1 = 11 + 2s-1

12 = + 2s-2

J2 = 12 + 2s-1

and continue to use 2s as the increment for the index II within the loop. The angle for eachindex is

s-1

,M-s V

,q

Z--d

WN

q=0 _. _{The bits q} _{0, q = 1,} . .

. q = s-

1 do not change when the index is incremented by2s, hence in a particular loop, the angle for a given index remains constant. Since w.e start with 11 = I,

s-1 2M-s n 2g

2M-s WN

q=0

= WN and the index into the sine array for all values of 11 in this loop is /M-s.

7s,2

+

(18)

The value in J1 is 11 + 2s-1 and it follows that s-1 2M's

E

2g Ms WN q=0 _w 2

If JI = 11 + 2s-1 then nq =.1 for q = Q and, q =.s-1, and we have

s-1 2Ms

E

n 2g

q0

₌

2Ms:

2s-1;) Ms M-1 = w 2 .+ 2

Ms N

w 2 M's, , N

,Ms

N But Nt,/ 2 7- W W N 2" = 4 N

The angles for the values 12 and J2 can be found in the same manner. The angle for 12 is jWN2M-s'

and the angle for J2 is -jWN2

This second loop is repeated s -2 times. For each iteration the yalue in 11. is incremented by I and the angle ,aSsociated With each index is incremented by

As an example suppose s 5. Then

11 =000 000 001 J1 = 000 010 001 12 = 000 001 001 J2 = 000 011 001 Then

Y5(11)+j Y154(I1) = Y4(11)+j Y

Y(1l)

[Y4(J1)+ YM_I 4 (J1)1w24

= Y4(11)+S Y700+ [Y4(.11)+,JYT(J1)11(cs(16)-cs(112)11

= y400+

111)cs (TO+ YT(31) 0042)+i LyT(l1)+NT(J1)cs,06)--y401)cs(112),]

The values for J1,12, and J2 are computed in the same manner.

After the second loop has; been repeated s - 2 times, s is incremented and the entire Step Three process is repeated. This process continues until s M, at which time the entire Ys4 array is divided by \TIC] and the transformation is complete.,

=

W2Ms

2M-s.

=

(19)

GENERATION OF SINE AND COSINE VALUES

The FFT algorithm uses the factorWx

. 27rx

where N for x = 0, 1, . . . ,N -1

Since e-Y = cos y - j sin y wx= cos (21,7

-j sin (

A first' glance would seem to indicate a need for the computation of N-values of both sine and cosine,

27r

at increments of However, further investigation reveals that such is not the case. By using the elementary trigonometric relationship

cos x = sin x)

2

and the fact that

for 27rk and

e-3 N

= wk Therefore for x = + k 4 27ry 7T N 2

where x is the location in the array of the sine values. Thus a separate calculation of cosine values is

unnecessary.

Moreover, not all of the sine values need be calculated, as the followingrelationships show:

27(.11+4 k)

N,

x = + K W X = e 2Nirk) = e \ Tr -

7r7r

But . e = cos 7.; -j sin T.= -j wx _j wk

The following set of equivalences is valid:

wx = wx for x = 0, 1, . N -1 4 N j w \x 4/ for x = , + 1 , . . . , -1 N N N 3N

w =

-w(x. 7)

for x = -1,

+ 1, .... ,--4 -1

wx = w _

2p)

ior x = -473N 3N

+, ,

. . . . , N - I N 27r

Therefore only sine values in increments of need be calculated. These values may be generated by coding, or they may be inserted within the program as a constant block.

for y = [27r IN

x)]

N 1 4 27r k we have cos 27rx = sin 2irx + N J.. . ir = e

-.

(20)

FORTRAN PROGRAM LISTING SUBROUTINE FFT(XR,XI,CS,M) DIMENSION XR(1),XI(1),CS(11 IF(M .EQ. MS) 1,2' 2 MS=M N=2"M N2=N/2 NA=N/4

C "" OPTIMAL PROCESSING FINITE FOURIER TRANSFORM "'""

C "" THE PROGRAM PERFORMS THE COMPUTATION

A(J+1) = (SUM(K=0,11(X(K41)*WN(KYJ)))/SORT(N)v J=0,N-.1

WHERE WN(L)=EXP(2."I*LIN)

IC _{AND N IS THE DIMENSION OF THE ARRAYS}

XR AND XI

C "" N MUST BE A POWER OF 2

M IS THE INTEGER SUCH THAT N = 2"M

CS IS A SCRATCH ARRAY OF MINIMUM DIMENSION' Nf4

XR AND XI ARE THE REAL AND IMAGINARY COMPONENT _ARRAYS. _FOR BOTH INPUT AND OUTPUT-- THE INPUT ARRAY IS OVERSTORED

C YYYY _{'GENERATE SIN-COS ARRAY CS ""}

CS(I) = SIN(DTHET) WHERE OTMET = 2.YPIYI/N TINC=6.2631853/FLOAT(N) CS(1) = SIN(TINC) DO 3 11=2,N4 TI=I1 CS(I1) = SIN(TIYTINC) 3 CONTINUE

REVERSE BIT AND FORM /ST ARRAY IN SEQUENCE

C

4""

_{ADDRESSES AND BIT-REVERSED MATES ARE}

(11912), (I3,J1), (149J2), (J304) C "" STORAGE SCHEME IS AS FOLLOWS

I19 j1 INTO 12,13 J2 INTO I1, 14 J3 INTO J19J4 I4, J14 INTO J203 .

C **1.* _{THE IF STATEMENTS PREVENT DUPLICATION' Of COMPUTATION WHEN} 12 = Ii OR WHEN 12 ,LT. I1 INDICATES THAT THE COMPUTATION HAS ALREADY BEEN PERFORMED

***1- WHEN I2= Ii, LT ALSO FOLLOWS THAT 14 = 139 J2 = J19

AND J4 =J3 SO THAT Ii, I3v Jit AND J3 FORM A cuNDAmtNTAt SET WHICH IS DUPLICATED BY 12, J29 I41 AND j4

***411, 1 DO 12 I1=1,N2,2 12 = IRVB(I11)+1 IF(I2 I1) 12,13,14 13 13. = 11+1 J1 I1 4. N2 J3' = 13 + N2 XRI3 = XR(I3) = XI(I3) XR(I3 )= XP(I1) XR1J11

Xr(I3), = (,XI(I1) - XI1J1)) C C -=

(21)

-XR(I1) XR(I1) + XR(J1)

XI(I1) XI(I1) + XI(J1)

XR(J1) = XRI3 + XR(J3)

XI(J1) = XII3 + XI(J3)

XR(J3) = h(RI3 - XR(J3)

XI(J3) = XII3 - XI(J3)

GO TO 12 14 13=12+1 14 = II + 1 J1 = I1 + N2 J2 = 12+ N2 J3 = 13 + N2 J4 = 14 + N2 XRI2=XR(I2) XII2 XI(I2) XR(I2) = XR(I1) + XR(J1) XI(I2) = XI(I1) + XI(J1)

XRI3= XR(I3)

XII3 XI(I3)

XR(I3) = XR(I1) - XR(J1) XI(I3) = XI(I1) - XI(J1)

XR(I1) = XRI2 + XR(J2) XI(I1) = XII2 + XI(J2)

XRI4 = XR(I4) XII4 = XI(I4)

XR(I4) = XRI2 - XR(J2) XI(I4) = XII2 - XI(J2) XR(J1) = XRI3 + XR(J3)

XI(J1) - = XII3 + XI(J3) XRJ4 = XR(J4)

XIJ4 XI(J4)

XR(J4) = XPI3 - XR(J3) XI(J4) = XII3 XI(J3)

XR(J2)= XRI4 + XRJ4

XI(J2) XII4 + XIJ4

XR(J3) = XRI4 - XRJ4 XI(J3) = XII4 - XIJ4

12 CONTINUE

C FORM REMAINING ARRAYS IN SEQUENCE C K IS INDEX OF ARRAY SEQUENCE

C "" KJJ = 2"(K-1), KII = 2"(M-K), NN = 2"K

C ALL SPECIAL CASES(C OR S = 0. AND +/- 1.0) ARE TREATED SEPARATELY C 13 IS USED FOR K BELOW TO SIMPLIFY THE COMPILED OBJECT PROGRAM C

.."

LOOP 2 DO 16 I1=19N.4 12 = I1 + 1 J1 = 12 + 1 J2 = Ji 1 X = XR(J1) Y = XI(J1) XR(J1) = XR(I1) -X XI(J1) = XI(I1) -Y XR(I1) = XR(I1) + X XI(I1) = XI(I1) + Y X = XI(J2) Y= -XR(J2) XR(J2) = XR(I2) - X XI(J2) = XI(I2) - Y XR(I2) = XR(I2) + X XI(I2) = XI(I2) + = = = -= -+

(22)

-16 CONTINUE C" **** _{END LOOP Z} XII = N4 KjJ = 2 NN = 4 DO 15 13 = XJH = KJJ XJJ = NN NN = NN + NN DO 17 I1 = 19NON J1 = It + KJJ 12 = I1 + KJH J2 = 12 + KJ.) -X = -XR(J1) Y = XI(J1) XR(J1) = XR1I14 X XI(J1) = XI(I1) -Y XR(I1) = XR(I1) + X XI(I1) = XI(I1) + Y X = XI(J2) = -XR(J2) llt(J2) = XR(I2) - X XI(J2) = XI(I2) XR(I2) = XR(I2) + X XI(I2) = XI(I2)

r

17 CONTINUE KII = KII,2 NN2 = KIT M2 = N. - XII DO 19 J = 2,KJHI C = CS(M2) CS(NN2) DO 1

II=

J4 = II + KJJ X = XR(J1)."C + xr(Ji).*s Y = XI(J1),C - XR(J1)*S XR(Ji) = XR(I1) - X XI(J1) = XI(I1) - Y XR(I1) = XR(I1) + X VI(I1) = XI(I1) + 12 = I1 + KJH J2 = 12 + KJJ X = XR(J2)*(-S) + V1(J2)*C Y = XI(J2),(-S) - YR(J2)*C XR(J2) = XR(I2) XI(J2) = XI(I2) -XR(I2) = -XR(I2) +. V XI(I2) = XI(I2) + 18 CONTINUE NN2 = NN2 + M2 = N2 - XII 19 CONTINUE 15 CONTINUE' X = SORT(FLOAT(N)) 00 20 = 1,NN XR(J) = XR(J),X XI(J) = 20 CONTINUE RETURN'. ND, 3,M -Y -+ S = J,N,NN -- X -J XI(J)f)(

(23)

IV. THE FFT ALGORITHM AS PROGRAMMED FOR THE HONEYWELL X16SERIES COMPUTERS

LIMITATIONS IMPOSED BY THE HONEYWELL SYSTEM

The Honeywell X16-Series computer for which this program was written uses paged, or sectored

memory for addressing, each sector consisting of 512 words. Any instruction in a given sector may access any other location within the same sector or any location in the zero sector without more cycle time or memory being used. However, if a sector other than the current one or the zero sector is to be used, the loader must store the address being accessed in the zero sector and change the instruction to use indirect addressing through the zero-sector address. For each instruction so changed, another word in the zero sector must be used. This is true even when the location being referenced is one that has already been referenced. Thus, even if only one location is referenced indirectly, a number of words in the zero sector will be needed. Indirect addressing adds one cycle to the instruction execution time. To minimize the amount of cycle time and memory required, the user should construct his routines to fit within one sector of memory (512 words). If he must use more memory, it is best for him to arrange for the indirect addressing within the program itself.

A Honeywell X16-series computer has one index register which can be accessed either in the usual way or as

memory location 000000. All addition, subtraction, loading, and storing operations are handled through an A register. Multiplication and division are handled through both the A and B registers. The B register can be accessed with shift instructions or by interchanging contents with the A register. No floating-point operations are used.

For specific information on the Honeywell computer and its assemblylanguage (DAP), consult the

Honeywell Series. 16 Programmers Reference Manual.2

PROCEDURE FOR CALLING THE FFT SUBROUTINE

1 N-1 2irnk

The FFT subroutine will compute the transform A(n) =

E X (k) e

N

for the complex integer array X of length 512. The complex array to be transformed must be stored in 1024 contiguous core locations with the first element of the array stored in the first address of a sector other than Sector Zero. All 512 real elements are stored in one sector with the 512 imaginary elements in the sector

immediately following. All elements should be of a magnitude less than or equal to 4096 to prevent overflow

during the calculations.

In the calling program, a location should be defined by the following statement:

ADD DAC ARR + 512, 1

where ARR is the address of the first element in the real array. Location ADD will contain the address of the first element in the complex array with bit 15 set to allow indexing. The following two instructions are needed to call the FFT subroutine:

LDA ADD

CALL FFT

The resulting array is stored in the same memory locations as those used for the input array.

2Programmers Reference Manual for the H316 and DDP-516 General Purpose Computer, Honeywell Computer Control Division, Framingham, Mass (April 1969).

(24)

STEPS IN THE CODING

The steps of the Honeywell version of the FFT subroutine parallel those of the Fortran version

described in Section hIll. A listing of the Honeywell version is provided on pages 24 -34.

The program goes through nine loops. All values computed in the even-numbered loops are divided by 2. The values computed in Loop 9 are divided by \ [T.. Thus after the computations are complete, the results have been divided by N/Tff.

As soon as the FFT starts executing, it sets up the following locations which will be referenced through indirect addressing:

XRN Contains the address of the real array + 512, with index bit set. XIN Contains the address of the imaginary array + 512, with index bit set. AXR Contains the address of the real array.

AX1 Contains the address of the imaginary array.

Step One

Step One consists of the group of instructions Li -1 to L2. The variable LOPC is used to control the number of iterations through the step.

At each iteration of LOPC, TEST is set to -1 and the index register is loaded with the contents of AXR. Values of the real array are then processed in the following manner. The index II contains the contents of LOPC with the index bit set. The instruction LDA* 11 loads the contents of the address AXR+LOPC into the A register. When the computations of the real variables are completed. TEST is incremented to zero and the procedure repeated for the imaginary values, this time with AXI in the index register.

When TEST is incremented the next time (to +1), control passes to the section of the program which increments LOPC by 2. When the value in LOPC exceeds 256, control passes to Step Two.

Step Two

Step Two consists of the instructions in the locations L2 to L3. The following variables are used:

AR Contains the address AXR, with the index bit set.

AR2 Contains the address AXR+2, with the index bit set.

Al Contains the address AXI, with the index bit set.

Al2 Contains the address AXI+2, with the index bit set.

For the first iteration, the index register has a value I such that mod(I, 4) = 0, and the real and imaginary values for (I) and X(I +2) are computed. The index register is then incremented by 1, giving a

value I such that mod(I, 4) = I. Since this value of I is odd, the values for X(1+1) and X(I +3) are computed.

repeating the instructions beginning at Location REP3. The index register is again incremented by 1, resulting in an index register value that is even, so control passes to Location FL2, which increments the index register by two more. Once again the value I of the index register is a value such that mod(I, 4) = 0. If the register value I is less than 512, control returns to the section which repeats the iterations. The results of all Step.

(25)

Step 'Three

Step Three consists of two loops. The computation of the first loop parallels, thatof Step Two,. The index register is incremented by KJH instead of by 1 before the second half of the loop is begun. LOPC, which is used to control the /lumber of iterations, is incremented by NN until it equals 512 at Which time control passes to the second loop in the step.

The procedure used in the second loop corresponds in most respects to that of the 'Fortran-coded version described in Section The values for II and J1 are computed and the index is incremented by KR-I. The temporary values (TEM1 and TEM2) used in the computation of 12 and J2 are computed using the same coding

as for Ii and J1. Then 'TEM1 is set equal to TEM2 and TEM2 is set to minus TEM1 before the actual

calculations of 12 and J2 are begun:

The following paragraphs explain the techniques used in the Honeywell assembly language version which 'differ from those used in the Fortran. version,

Suppose ACS is the address of the first location of the sine array and Kit is the index, for the angle of interest with KI,=Ku - 1. Then the values AS equal to sin (K11) is located at address ACS+ KI (the ,first element of the sine array has- an index (KII) equal to 1). Thevalue ,AC equal to cOs (KM) is located at address ACS +

-KI -2.

4

-[

cos x sin

(

-4

cos (KM = sin

(

-

KII)

which is 'found at location ACS + Kit

-4 4

N

Since KI = KII - we have ACS + - Kit-1 = ACS +

- Kt

.2.]

The .sine array (CS) contains values derived from multiplying the sines by 2'12, or 4096. When these CS values are used by the program in multiplication, the quotients (in the combined A and B registers) reflect this multiplication factor, as is shown in Figure L.

A REGISTER B REGISTER

L'O 00000XXXXXXXXX, XXXXXXXXXXXXXXX

1, 2 16 2..

bits of interest 212 offset

Figure 1 The A and B Registers afterMurtiplicatinn

111.

=

-

- 1.

(26)

-The condition of the A and B Registers subsequent to the performance of the instruction LLS 3 is shown in

Figure 2.

A REGISTER B REGISTER

H000XXXXXXXXXXXXXXX

XXXXXXXXXXXX000

2 16 2 16

Figure 2 The A and B Registers after Shift

The A register now contains the result desired (in which the 212 offset has been eliminated).

The magnitude of the values used as input to the FFT subroutine should be such that overflow does not occur. A check is made after the performance of each calculation that might result in an overflow. If the overflow indicator is set (indicating that an overflow has occurred), the sign of the value contained in the A register will be the opposite of what it should be: If the result should have been negative, it is changed to 1000008' the largest possible negative value; if the result should have been positive, it is changed to 777778' the largest possible positive value. This procedure tends to minimize any error caused by overflow.

After 13 is incremented, a check is made to determine whether the value in 13 is odd. If it is, the real and imaginary arrays are divided by 2 before computations for the new value of 13 are begun. When 13 has a value of 10, the real and imaginary arrays are divided by the value SQ2, and the calculation is complete.

The value SQ2 is equal to N/7.2 4096. When a value is divided by \727, the B register is cleared and the value is placed in the A register. An arithmetic long shift of 3 is executed. The value in the combined A and B register is now the value to be divided by VT multiplied by 4096. Division by SQ2 yields the

desired result.

GENERATION OF SINE VALUES

The 127 sine values defined for the CS array (the array referred to by the assembly language coded listing on page 31) were obtained using the Fortran statements listed in Figure 3, page 22. Note that the values have been multiplied by 212 (which is 4096) and then converted to integer format. These integer values were printed in octal.

i000XXXXXXXXXXXXXXX

(27)

PROGRAM GENR(OUTPUT,TAPE6=OUTPUT) DIMENSION CC(13),SS(13),CS(120) DIMENSION IS(126) M=9 N=512 SCALE=2**12 SSS=SORT(.5) CCS=SSS SORT2=CCS SS(3)=SORT2 CC(3)=SS(3) DO 3 I1=40 CCS=SOPT(.54(1.+CCS)) CC(I1)=CCS SSS=.5.SSS/CCS SS(I1)=SSS 3 CONTINUE N4=N/4 CS(1)=SS(M) CS(N4-l)=CC(M) I1=M-1 I2=I1-1 CS(2)=SS(I1) CS(N4-2)=CC(I1) CS(3)=SS(I2)+CC(M)-CC(I2),SS(M) CS(N4-3)=CC(I2).CC(M)+SS(I2)*SS(M) IH=M-4 00 L. I=2,IH I1=2*.I ILL=I1+1 IHH=2.I1-1 J1=N4-I1 J2=M-I CS(I1)=SS(J2) SSS=CS(I1) CCC=CC(J2) OS(J1)=OCC DO 4 I3=ILL,IHH J3=N4-I3 14=13-11 J4=N4-I4 CS(I3)=CCC,CS(I4)*SSS*CS(J4) 4 CS(J3)=CCC+CS(J4)-SSS+CS(I4) N8=N/8 CS(N8)=SORT2 DO 1110 J=1,127 IS(J)=CS(J)*SCALE 1110 CONTINUE WRITE(6,1003)(IS(J),J=1,127) 1003 FORMAT(11H1SINE TAeLE/(1608/)) STOP END

(28)

BIT-INVERSION ROUTINE FOR THE CASE N = 512

The Fast Fourier Transform routine uses a bit-inversion process for an index I

_{u -1,}

i _jm-i)

where jp is the binary value of j for the p th binary place. The result after inversion is J = (j,_, , jm-2, . .

For the case where M = 9 and N = 512, the binary values can be treated in groups of three. Each group of three is inverted separately by using an 8-place table lookup;

1 = 0 = 0002 inverts to J = 0 = 0002

I = 1 = 0012 inverts to J = 4 = 1002

etc.

After being converted, each group of three bits is shifted, according to its position in I, and the results are summed. The right-most group of three binary values in I should become the left-most group in J by multiplying by 64 (shifting six positions to the left). The middle group remains the same and is multiplied by 8 (shifted three positions to the left). The left-most group of binary values is not shifted prior to being added

to the sum of the other two bit patterns.

For the Honeywell X16 computer series, the FFT program requests bit inversion with the statement CALL INVT where the A register contains the value to be inverted. When control returns to FFT, the A register will contain the inverted bit pattern. The listing for INVT is included within the FFT listing which is given on the following pages.

(29)

FFT AND INVT PROGRAM LISTINGS, HONEYWELL X16-SERIES COMPUTERS

0001 0002 0003 SWAP FFT SUPP INvT PEL 0004 00000 TI ASS 1 0005 oaool 12 ASS 1 0006 00002 13 ASS 1 0007 00003 I4 ASS I 000P 00004 J1 ASS 1 0009 00005 J2 ASS 1 0010 00006 J3 ASS I 0011 00007 J4 ASS 1 0012 00010 IFS! ASS 1 0013 00011 TEMI ASS 1 0014 00012 TFm7 ASS 1 0015 00013 TEm3 ASS 1 0016 00014 040000 4AS1 OCT 040000 0017 00015 000011 M DEC 9 001A 00016 000400 N2 DEC 256 0019 00017 001000 N nEC 51? 0020 00020 000200 N4 nEc 12P 0021 00021 KIT ASS 1 0022 00022 KJJ ASS I 0023 00023 NN ASS 1 0024 00024 KJm ASS 1 0025 00025 ApJ1 ASS 1 0026 00026 AIJI ASS 1 0027 00027 0 000622 ACS DAC CS 002P 00030 013240 c02 OCT 13240 0029 00031 C 000000 AxR PAC 00 0030 00032 0 000000 Ax! DAC 0031 00033 10Ar ASS 1 0032 00034 xi ASS 1 0033 00035 0 000000 AP OAc 41.0 0034 00036 0 000000 Ap2 nair 0035 0036 0037 00037 00040 00041 0 000000 0 000000 AT DAr

.0

Al2 nAc 44 AC ASS 1 003S 00042 AS ASS 1 0039 00041 xRN ASS I 0040 00044 wTN ASS I 0041 00045 100000 NFGm OCT 100000 0047 00046 077777 POSm OCT 077777 0043 0044

AxP CONTAINS ADDRESS AxT CONTAINS APORFSS

OF REAL APRA),

OF IMAGINARY ARRAY

0045 CONTAINS AXP.512 WITH INoEX SET

0046 xIN CONTAINS AxI.512 wITH INnEx SET

0047 00047 0 000000 FFT nAc 0048 00050 0 04 00043 STA xAN 0049 00051 0 06 00017 Ann N 0050 00052 0 04 00044 STA XI-N 0051 00053 0 07 00014 cU mAS1 0057 00054 0 07 00017 SUP N

(30)

0053

nonss

a

04 00032 STA AxT 0054 0055 00056, 00057 o 07 0 04 00017 00031 SUP STA N' Axp 0056 00060 0 02 01067 LOA =0

0057 _:0 _{11 PFPFORmS PEORDPRING. &ND COMPUTES 'FIRST}

ARRAY.

0058 _{tORC COUNTER FOR Loop}

0059 4 _{TNCREmENT BY 2}

0060,

0061 00061 0 04 00033 11

UNTIL GREATER THAN 25A

STA LoPC 0062 00062 0 10 00000 CALL INVT 0063 00063 0 II 00033 CAS LOPC 0064 00064 0 01 P0132 JMR GI 0065 00065 0 01 00074 jmp 0066 00066 0 02 00033 FL1 LOA. LOFT 0067 00067 10 06 01066 40,0 0068 0069 00070 00071 0 11 0 01 00016, 00222 CAS JMP N2 L2 0070 00072 0 01 00222 jMP L2 0071 00073 0

ol

00061 JMR LI

0072 _* _{INVERSE AND' ORIGINAL APE' FOML}

0073 00074 0' 02 00033 FO IDA tOPC 0074 00075 0

as

00014 ERA mASI 0075 00076 0 04 00000 STA II 0076 00077 0 06 00016 ADD N2 0077 00100 0.04 00004 STA JI 0078 00101 0 06 01065 ADD =1 007(7 00102 0 04 00006 STA J3 0080 00103 0 017 00016 SUB AI? 0081 00104 0 P4 00002 STA 13 0082 00105 0 02 01064 LOA =-1 0081 00106 0 OA 00010 _STA _TEST 0084 00107 0 35 00031 LOX AXR 0085 rOmPUTE VALUES 0086 010110

-o

02' 0,0002 pppl LOA* 13 0087 00111 b 04 00011 STA TEml 0088 0011? -0 0? 00000 LOA* 11 0089 00111 -0 07 00004 1SUR* J1 0090 00114 -0' 04 00002 STA* 13 0091 00115 -0 02 00000 L0A* Ti 0092 00116 -0 06 00004 ADO* J1 .0093 00117 -0 04 00000 sTA* 0094 00120 0 02 00011 LOA TPA1 0095 00121 -0 06 00006 Ann* J3 0096 00122 04 00004 STA* J1 0097 00123 0 02 00011 tDA TEml 0098 00124 -0 07 00006 SUP* J3 0099 00125 -0 04 00006 STA* J3 0100 00126 o 12 00010 IRS TEST 0101 00127 0 01 00066 jmP ELI 0102 00130 0 35 00032 LOX AX! 0103 00131 "0 01 00110 jmP PERI

0104 INVERSE GREATER THAN ORIGINAL

0105 0013? 0 Os 00014 GT EPA mAS1 0106 00133 n 04 00001 STA 12 010.7 00134 0 06 n0016 /ion N? 0108 0109 00135 00136 0 04 0 06 00005 01069 STA, ADD J2 =1 011,0 00137 0 04 00006 _s14 _J3 0111 0112' 0113 00140 00141 0014? 0 07 0 04 0 02 00016 00002 00031 SUR, STA, LOA AT 13 LoPr, FO =2 * II -0

(31)

0114 00143 0 05 00014 ERA mAS1 0115 00144 0 04 00000 STA 11 0116 00145 0 06 00016 ADO N2 0117 00146 0 04 00004 STA Ji 0118 00147 0 06 01065 ADD *1 0119 00150 0 04 00007 STA J4 0120 00151 0 07 00016 SUP N2 0121 0015? 0 04 00001 STA 14 0122 00153 0 02 01064 LOA *-1 0123 00154 0 04 00010 STA TEST 0124 00155 0 35 00031 10* Axp 0125 romPuTE VALUES 0126 00156 -0 02 00001 PEP? I0A* 12 0127 00157 0 04 c0011 STA TEml 0128 00160 -0 02 00000 L0A* 0129 00161 -0 06 00004 ADn* JI 0130 00162 -0 04 00001 514* I? 0131 00163 -0 02 00002

LOA

13 0132 00164 0 04 00012 STA TEMP 0133 00165 -0 02 00000 104* 11 0134 00166 -0 07 00004

SUP

J1 0135 00167 -0 04 00002 STA* 13 0136 00170 0 02 00011 LOA TFm1 0137 00171 -0 06 00005 ADD* J2 0138 00172 -0 04 00000 STA. Ti 0139 00173 -0 02 00001 LOA* 14 0140 00174 0 04 00013 STA TFm3 0141 00175 0 0? 00011 iDA TEmi 014? 00176 -0 07 00005 SUP*

j?

0143 00177 -0 04 00003 STA* T4 0144 00200 0 02 00012 LOA TEm2 0145 00201 -0 06 00006 400* J3 0146 00202 -0 04 00004 STA* JI 0147 00201 -0 07 00007 L04* J4 0148 00204 0 04 00011 sTe TFml 0149 00205 0 02 00012 Loa 1Fm2 0150 00206 -o 07 00006 SUP* J3 0151 00207 -0 04 00007 STA. j4 0152 00210 0 0? 00013 LOA TEm3 0153 00211 0 06 00011 ADO TFmi 0154 00712 -0 04 00005

STA

J2

0155 00213 n 02 nOn13 IDA TEm3

0156 00214 0 07 00011 SUP TEmi

0157 00715 -0 04 00006 ST** J3

0158 TEST -1 INT1Fx PEG CONTAINS AxP okonPFSS

0159 TEST n INDEX PEG coNTAINS AX! *DnPFSS

0160 00216 0 12 00010 TRS TEST 0161 00717 0 01 00066 JMP ELI 0162 00220 0 35 00032 10x Ax! 0163 00721 0 01 00156 jmp PEP? 0164 LOOP ? 0165 AP CONTAINS AxP.1

0166 APP CONTAINS AxP.2.1

0167 AT CONTAINS Ag1.I

0168 AT? CoN1AINS 4x1.2.1

0169 ALL RESULTS APE

urvrnrn HY 2

0170 00222 0 15 01067 L2 lox .0

0171 00223 0 02 00031 IDA AXP

0172 00224 0 06 00014 ADO mAS1

0173 00225 0 04 00035 STA AR

(32)

0174 00226

06 01066

Ann, =2

0175 00727

04 00036( STA AR?

0176 00230

0 02 00032'

1.DA AXT

0177 00231

06 00014

Ann mAsi

0178 0073?

0 04 0003/

STA A/

0179 00231

0 06 01066

ADD =2 0100 ₀₀₂₃₄

04 00040

_STA _Al2

0181 _{sFT HP TEml} ANn TE1M2

0182 00235 o 62 00016

L2A LOA* AR?

4183 00736

M4 ₆₀₀₁₁ STA TFml

0184 00237

02 00040

1DA. AT 0185 00740

04 00012

STA Tfm? 0186 * COmPUTf vALUIFS

0187 00241

02 00035

PEP1 LOA* AP

0188 00242

0 07 00011

SuR '1EM

0189 00243

0405 77

Ams 1

0190 00244

-0 04 00036

ST.A AQ2

0191 00745

,0 02 00035

tOA* AQ

0192 00246

0 06 00011

Ann TEM1

0191 00247

0405 77

ARs 1

0194 00254

-0 04 00035

STA* AP 0195 00251

-0 02 00017

LOA* AT

0196 0025?

0 07 00012

SUR TEm2 0197 002511

0405 71

1

.0198 00254

-0

04 00046

STA Al2

0199 00255

-0 02 00037'

ICA* AT

Woo 00256

0 06 0001?

Ann Tfm2

0201 00757

0405 77

APS 1

0202 00260

-0 04 00037

STA*, Ay

0201 TNCREmFNT ANn CHECK FOP SFCONn HALF

0204 ₀₀₂₆₁

0 07 00006

LDA

0205

0026?

0 06 01065

ADO

0206 00263

101100 SLN

0207 00264

01 00274

JMP FL?'

0208 00265

0 OA 00000

STA 10

0209

sFT up Fop FcoND.

HALF

0210 00266, ,0

op 40040

(OA,* AT?

0211 00267

0 04 00011

STA TEml

021? 00270

-0 02 0003*

IDA* AR?

0211 00271

140407 TCA

0214 0077?

0 04 nonl?

STA TEN?

0215 00273

0 64 00241 JMP PEP3

0216

* TNC X PEG RN 7 0217 coop IINTTL GF 512.

0218 00274

0 06 010,66 FL?

Ann ,m2

0219 00775

o 04 noono STA 0

0220 00276

q 11 00017

CAS N

0221 00777

0 01 00102

JMP _L3

0222 onloo,

o 01 00302

jmp L3

0223 00301

4 01 00735

JP

LPA 0224 a 4:n0P5 3 To 19

0225 00102

n 0? 60017 L.3

!.DA

0226 00303

0405 76

ARS

0227 00304

,0

04 00021

STA 4(TT,

0228 00305

4 02 01066

LOA =2 0229 0030.6

04 00022

STA

0230 00307

0 02 01063

LOA =4 .0231 00310

0 04 00023

STA NN

0232 00111

6 02 111062 LOA =3

4213 P0117

0' 04, 001002 STA 13 0 0 0 0 -0 0 -0 ARS 0 =1 0

*

0

(33)

0234 0235 0236 0237 00313 00114 00315 00316 0 02 00022 04 00024 0 02 00021 0405 77 L31 1DA STA LOA

AS

KJJ KJH KIT 1 0238 00117 0 04 00021 sTA KIT 0239 00120 0 OP 00021 iDa NN 0240 00321 0 04 00022 STA KJJ 0241 00322 0 06 00023 Ann NN 024? 00121 0 04 00023 STA NN 0243 00324 0 OP 00022 LOA Kjj 0244 00325 0 06 00014 ADO mAS1 0245 00326 0 06 00031 ADO AXR 0246 00327 0 04 00025 STA APJ1 0247 00130 0 07 00031 SUP AxR 0248 00131 0 ON 00032 Arm *XT 0249 0013? Al 04 00026 STA ATJ1 0250 00333 0 35 01067 LOX so 0281 00134 o 19 00033 sTx LOPC

u25? * FIRST PART

0253 00335 -0 02 00625 P17 LuA* ARJ1 0254 0255 00336 00337 0 04 00011 -0 0? 00026 STA LOA' Tpm, ATJ1 0256 00140 0 04 0001? STA TFm2 0257 00141 -0 02 00035 PEP4 LOA* AP 0258 00142 0 07 00011 cu8 TFml 0259 00143 100001 SRC 0260 00144 0 10 00613 JST ovF? 0261 00145 -0 04 00n2c cTA* APJ1 026? 00346 n 06 00011 Am) TFH1 0263 00347 0 06 00011 ADO TFml 0264 00150 100001 SRC 0265 00351 0 10 00613 JST OvF2 0266 0035? -0 04 00035 STA. AR 0267 00353 -0 OP 00037 LIJA* AT 0268 00354 0 07 0001? SUR TFm? 0269 00155 100001 SP4c 0270 00356 10 00613 JST nvF2 0271 00357 -0 04 00026 STA* ATJ1 027? 00160 0 06 00012 ADn TEmP 0273 00161 0 06 00012 ADD TEm2 0274 0036? 100001 SPC 0275 00363 0 10 00613 JST OVF2 0276 00364 -0 04 00037 cTA* AT

0277 rHFCK FOR SFcOND TTMF THROUGH LOOP

0278 00365 0 02 00000 LOA 0

0279 00166 0 03 00024 ANA KJM

0280 00367 100040 SZF

0281 00370 0 01 00402 jMR FL17

0282 SET uP FOR SECOND HALF

0281 00371 0 0? 00000 10A 0 0284 00172 0 06 00024 Ann KJM 028C 00373 0 04 00000 STA 0 0286 00374 -0 0? 00026 LOA* ATJ1 0287 00375 n 04 00011 STA TFml 0288 00376 .0 02 00025 LOA* Apjl 0289 00177 140407 TCA 0290 00400 0 04 0001? STA TFm? 0291 00401 0 01 00341 jMR PFP4 0292 CHFCK TO SEE IF FINISHED 0293 0040? 0 OP 00011 FL17 LOA LoPC *

(34)

-0294 00403 06 00023 *OD NN 0295 00404 0 11 00017 CAS N 0296 00405 0 01 00412 Jm0 L16 0297 00406 0 01 0041? JmP LL6 029R 00407 0 04 00033 'Sin it0pc 0299 00410 15 00033 tOx LOPC 0300 0301 00411 0 0] 00335 Jmi, sECON0 R17 PART 0302 0041? 4 35 01065 06 LOX 2c1 0301 00413 0 15 00033 Six LoPC 0304 00414 0 0? 41467 LOA =0 0305 00415 006201 TAR 0306 00416 0 0? 00021 tOA KIT 0307 00417 0 07 01465 SUR :1 030R 00420 4 04 00034 STA KT 0309 00421 0 02 00434 R16 104 KT 0310 004?? 0 06 00027 ADP ACS 0311 00423 0 n* 00042 STA AS 031? 00424 4 017 00034 SUA Kt 0313 00425 0 07 00034 SUR xf 0314 00426 0 06 00020 ADn N4 0315 04427 0 07 01066 SUR ,7.2 0316 00410 0 04 00041 5TA AC

0317 rpmPHTE TEml ANO TFm2,

031R 00431 00n25 p16,A 1DA* ARM

0319 00432 -4 16 00041 ImPy* AC 0320 00433 0411 75 ILLS 3 0321 004314 n 04 00011 STA TEml 0322 00435 -4 0? 00026 LOA* AIJI 0323 00436 -0 16 0004? mRy* AS 0324 00437 0411 75 LLS 3 0325 00440 0 06 04411 400 TEml 0326 00441 100001 SPC 0327 0044? 0 10

owa

Jsr

OvF2 032R 00443 n 04 00011 STA TFml 03291 00444 _0 0? 00025 LIDA*

ARA

033n 00445 =0 16 60042 KAPy* AS 0331 00446 0411 75 ILS 3 0332 00447 0 04 0001? sTA TFm2 0333 00460 =0 132 60026 LOA* AfJ1 0334 00451 -0 16 00041 mPm* AC 0335 00452 4411 75 ILS 0336 00453 07 00012 SUR TFm2 b337 04454 100001 SPe 033R 00455 0 10 00613 JS1 0vF? 033R 00456 0 04 00012 TEm2

0340 CHECK FOR SFCON'O, HALF

0341 00457 0 0? 00000 ILDA 0 0342 00460 0 03 00024 .ANA KJH 0343 00461 100040 SZE 0344 00462 0, 01 00516

j"

CHG 0345 4 rnmPUTE vALUF5 0346, 0347 00463 00464 '0 02 00035 0 07 00011 CONT LOA* SUP AR TEml 034R 0046s 141001 SRC 0349 00466 _O 10 00613 JST OvE2 0350 00467 -0 04 00025 STA* ARJI 0351 00470 0 06 00011 ADO TEml 035? 00471 oc 00011 AOn TEmi 0353 00472 .10000T SRC 0 0 * -0 0

(35)

JST ovF2 STAG AR LOA* AT SUR TEMP SRC J5T OvE2 ATJ1 ADO TEMP ADD TFm2 5RC J5T OvF2 STA* AT

cHECK FOR SECOND HALF

LOA 0 A00 mJH sTA 0 ANA XJH 52E. jmp 1816A JMP EL16

tmANGE TEm1 ANn tEm2

Cr,H LDA TEm1 TCA sTA TEml LDA TEN2 STA TEN1 tOA TFm3 STA TEMP JmP CnNT

4 noNE ALL OF lonp" FL16 LOA Ann KJJ 616 0 CAS IN JMP JmP

0.2

jMP Pl6A LOA LnPc ADO =1 STA

Lon

cAs KJR YAP OUT jmp _OUT LOA 1JOPC LDA .K1 ADO KTT STA KT YAP P16: TNCRfmFN4 ANO TEST 11

toti

/1

AOn STA 13 SLN _imp FIN LOX -512

nivInF By Two IF 13 NOW 000, pHF7 IDA,* XPN, ARS 1 STA* xpN LOA* xl" 0354 00473 0 10 006131 0355 00474 -0 04 00035 0356 00475 -0 02 00037 0357 00476 0 07 00012' 0358 00477 100001 0359 00500 0 10 00613 0360 00501 -0 04 00026 0361 00502 0 06 00012 0362 00503 0 06 00012 0363 00504 100001 0364 00505 0 10 00613 0365 00506 =0 104 0.0037 10366 0367 00507 0 12 40000 0368 00510 0 06 00024 0369 00511 0 04 00000 0370 00512 03 00024 0371 00513 100040 0372. 00514 0 01 00431 0373 0374 00515, 0 01 00526 0175 00516 0 02 00011 0376 00517 740407 0377 005?0 0 04 00011 037A 00521 0 02 00012 0379 00522 0 04 00011 03A0 00523 t 02 00011 0381 00524 04 00012 0382 00525 0 11 00463 0383 0384 00526 02 00000 10385 00527 OA 000?? 0386 00530 04 00000 03A7 00531 11 00017 0388 00532 01 00535 0389 00533 ,0 01 00535 0390 00534 01 00431 0391 00635 0,2 00033 039? 00536 06 01065 0393 00537 04 00033 0394 00540 11 00024 0395 00541 01 00550 0396 00542

0 ni

00550 0397 00543 35 00031 0398 00544 0 0? 00034 0399 00545 06 00021 0400 00546 04 00034 0401 00547 0 01 00421 0402 0403 00550 0 02 00002 0404 00551 0 06 01065 0405 00552 n 04 00002 0406 0407 00551, 00554 101100 0 01 00566 0406 00555 0 35 01061 0409 0410 00556,

.0 02 00041

0411 00557 0405 77 0412 00560

-0 04

00041 0413 nn561 sn O2 00044 0 0 0 0 0 0 0 *.3 0 0 0 0

(36)

EOLE*E2924E*Se.t9,2610*ESSEEeSIEZ 100 t410,0126E2000,[162991.109160E5t 13U otE100 1E900

etetoo

':9.E90u liII00 5E900 Lt0I00 ,E900 130 4Z*155EUttkeI.EI2161([01,014591 59/000 EE900 05t0 (01000 2E900 . t9000 1E900 .4E5000 0E900 S5,000 LE900 ELE000 92900 01(000 52900 92E000 t2900 **I000 (4900 E00190LES6WOELE,01(69226,11429 lOu Sj 290000 2E900 6,470 24A0

am'

(1900 10 0- JE900 847,0 NO3N[ *OA e.)34 519000 e0 o-

40000

LttO *awl' E4A0

(1900 IO Of 41900 9,t0 wSUd VO1 9,000 20 a to ozsoo dm' 293N o

Is

00,101 *4' Jvu idAj 000000 0 SNU113bdU0 aU1.183A,U 11A dw1'. ' 4,000 10.6- cit)s dwr VLSOO 10 0 U Sul 00000 21 0 Nix *viS vir000 '0

Of 0 LI 0E000 AIO ZUS

beg! Si. 1047'0 NIX 4001 14000 20 0- tAVI loe000' U: L90[0. eu

a

0, E1000

*vis

Nax

a-

0E000 AIU

dos

41 u E Se1 Si. 4401 147000 20 WVI 10200,0 Oa Val Aiu9 1.9010 20 0

xal

190LO St 0

Jultau

AR AdOS

t2)

EL1 dWr t1COU to

a to EtEoo

ale

(El

awr ELSUu IU u

It

SIU00

-

SIO 0

vol

LI NIA dut)OU 0 0

'IV

El * dWr 14HS 9SS00,10 0 SeI 00000 21 0,000

'0

0-

50,0 Lt. .StiV I.

seuao

stezoo v912U 4 0'1200 ceoeoo E/LIO0 299100 109100 u2SI00 kEtIO0 55E100 ES900 I5900 25,0 05900 /t900 9,900 91900 tt900 0,900 Et900 15,0 1,900 04,900 91900 -Stt0 51900

tt*D

*1900 0,470 (1900 EttO .11,0

1900

U0-

6EtO 11900 01900 trE40 1.0900 lEt0 90900 9E4,0 50900 SEt0 ,0900 47E00 C0900 EE*0 e0900 eEtO 10400 lEt0 00900 00,0 LL900 bEt0

9L00

dEt0 5L900 12470 tis00 92*0 LLyuu tEt0

elsoa

EE*0 1/900 E2t0 01900 12,0 /9t00 02,0 99500 bl,0 V1,0 59S00 LIVE) a9SOU 91,0

E900

SI,0 e9s00 ,14,0 VO1 NbX -= 0 U

(37)

00653 002404 00654 002464 00655 002543 00656 002623 00657 002702 0453 00660 007760 ncT 2760,3037.3115.3173.3251.3327.3404 00661 003017 00662 003115 00663 003173 00664 003251 0066C 003127 00666 003404 0454 00667 001461 OCT 346193536,3612,3667.3742,4016.4071 00670 003536 00671 003612 00672 003667 00673 003742 00674 004016 00675 004071 0455 00676 004144 OCT 4144,4217.4271.434104415.4466.4537 00677 004?17 00700 004271 00701 004341 00702 004415 00703 004466 00704 004537 0456 00705 004607 OCT 4607,4660,47?7,477705046,5115,5163 00706 004660 00707 004727 00710 004777 00711 005046 00712 005115 00713 005163 0457 00714 009231 nC- 9231.5276,5343,541099454,5520.5563 00715 005276 00716 005143 00717 005410 00720 005454 00721 005920 0072? 009561 045P 00723 005626 OCT 5626,5670.5712.5774.6035.6076.6136 00724 005670 00725 00571? 00726 005774 00727 006035 00730 006076 00731 006136 0459 0073? 006175 OCT 6175.6235.6273.6331,6367,64?4,6461 00733 006235 00734 006273 00735 006331 00736 006367 00737 006424 00740 0460 00741 006461 006515 OCT Ac15,6551.6604.6637.6671,67?2,6753 00742 006551 00743 006604 00744 006637 00749 006671 00746 006722

(38)

00747 006753 0461 00750 007004 OCT 7004,7034,7063.7112.7140,7166,7213. 00751 007034 00757 007063 00753 ₀₀₇₁₁₂ 00754 007140 06755 007166 00756 007713 0462 00757 007240 OCT _{7240,7764,7310.733197355.7377,7420} 00760 007264 00761 ₀₀₇₃₁₀ 00762 007133 00763 ₀₀₇₃₅₅ 00764 007377 00765 007420 0463 00766 007441 OCT

_{7.41,7461,75n0os17.7515.7c51.7570}

00767 007461 00770 007500 00771 007517 00772 007535 00773 007553 00774 007570 0464 00775 007605 OCT _{7605.7671.7614.7647.7661,7677.7703} 00776 007671 00777 007634 01000 007647 01001 007661 01002 007672 01003 007703 0465 01004 007713 OCT 7713,7773,7732.7741.7747,7754.7760 01005 007773 01006 007732 01007 007741 01010 007747 01011 007754 01017 007760 0466 01013 007764 oCT _{7764,7770.7773,7775,7776,7777} 01014 007770 01015 007773 01016 007775 01017 007776 01070 ₀₀₇₇₇₇ 0467 01021 000000 TABI OCT

0,4,7.6,1,5,1.7

01022 000004 01023 000002 01024 000006 01025 000001 01026 ₀₀₀₀₀₅ 01077 000003 01030 000007

04614 01031 _{n 07 01021 ATAR LOA} _TAAL

046P 01032 _{TEMP BSS}

0470 01033 0 000000 TNyT nAc

0471 4 _{(-ONO-PT HIGH ORDER BITS Awl}

STORE 0472 01034 0400 72 LBL 6 0473 01035 0 06 01031 *on ATAB 0474 01036 0 04 01037 STA 41.1 0475 01037 0 02 00000 LOA AA 0476 01040 _{0 04 01032} _STA _TEMP

(39)

OW WARNING OR ERROR FLAGS

DAP-16 MOD 2 REV'. n 06-28-71

AC 000041 ACS 000027 AT 0.00037 ATP 001140

AIJII 000076 AR 000035 ARP 4100036 APJ1 000025

AS 000042 ATA8- 001031 AxT 00.0032 AXP 000031

CMG 000516 CONT 000463 CS 000622 E0 000074

FET 000047 FIN 0100566 ELI 000066 1116 00032.6

FL17 000402 EL? 000274 GT 000132 Ooonon 12 000001 13 000002 I4 000003 INvf 001 031 J1 00000.4 J2 000005 J3 000006 J4 0010007 KT 000034 KIT nonnal KJm. 000024 KJJ 000022 LI 000061 L16 000412 LP 000222' LPA 400235 L3 000302 L33 000113 LOPC 000033 000015 m4S1 000014 ni 001017 Np . 000016 "Nk

loon2o

NEC,? 000620 NFOm 000045 NN 1001123 poncsn

nVE2 000613 ROSPA 000046 R16 001421 R16A 000431

R17 000335 'REP1I 000110 PEP? 110156 44E423 000P41

RER4 000341 SHET 000556 SO2 000030 000574

TARL 001021 TEmi 000011' TEm2 000012 '1 EM 000013

TEMO' 001012 TEST nonnlin xTN 000044 )(RN 000043 A047.8 0479 0480 0481 0482 0483 0484 0485 0486 0487 0488 0489 0490 0491 0492 0491 0494 049s 0446 0497 0498 11041 0104? 01043 01044 01045 01046 01047 01050 01051

nios2

01053 0105'. 01055. 01036 010S7 01060 01061 01067 01063 01064 01065 01066 01067 0 02 01067 0410 75 0 06 01031 0 04 01045 0 02 00000 v415 75 06 01032, 0 04 0103? 0 0? 01067' 0410 75 0 06 01031 04 01055 0 02 0.000V 0415 77 0 4)6 01037 -V 01 01013 177000 000001 -0000(14 177777 000001 000002 0,00000 LOA tl tLL 3 ADO . LIAR STA 0.1 LOA **

cmIFT AND Ann To Icml,

ALS 3

A00 TEMP

SIA TEMP

.CONVERT LOW 'ORDER

Rrrs

LOA an

LLL

3-ADOI ATAR

STA 4,1

LOA

SHIFT ANn ADD TEMP' RESULTS TN. A REG AL5 6 ADO TEMP jm0A TAIVT NO 11 0

(40)

FLOWCHARTS FL1 ADD 2 TO LOOP COUNTER ENTER STORE ADDRESSES OF ARRAYS LOOP COUNTER =0 L1 STORE LOOP COUNTER IN A FIND INVERSE + OR GT EG SUBROUTINE FFT STORE ADDRESSES 12, J2, 13, 13 11, J2, 14, J4 STORE ADDRESSES IN 11, 13, J1, J3 TEST = -1 INDEX REGISTER = ADDRESS OF REAL REP1 COMPUTE 11, 13, J1, J3 LOOP 1 TEST = -1 INDEX REGISTER = ADDRESS OF IMAGINARY INDEX REGISTER = ADDRESS OF REAL REP2 COMPUTE VALUES INDEX REGISTER = ADDRESS OF IMAGINARY INVR -LOOP COUNTER YES = YES

(41)

.1.2 L2A REP3, COMPUTE XRUI, ,St DIVIDE BY 2 FL2, INDEX = SET UP XR, XR +2, XI, XI +2 SET 'UP TEM1 AND TEM2 ,ADD 1 TO REGISTER YES ADD TWO TO INDEX NO LOOP 2 SET UP NEW TEM1 AND TEM2

Y

0 XR(I), XI(I) NO

(42)

L3 COMPUTE Kit, KJJ, NN 13 = 3 L33 COMPUTE KJH, KU, KJJ, NN LOPC = 0 SET UP XR KJJ, 1 XI + KJJ, 1 R17 COMPUTE TEM1 & TEM2 REP4 COMPUTE X(J1I, X111) ADD NN TO LOPC ADD KJH TO INDEX LOOP 3 COMPUTE NEW TEM1 & TEM2 NO _{INDEX =} LOPC

(43)

L16y

INDEX = 1 LOPC = INDEX COMPUTE KI R16 COMPUTE AS & AC R16A COMPUTE TEM1 AND TEM2 NO CHG CONT COMPUTE X(J1) & X(I1) ADD KJH TO INDEX LOOP 3 (Continu FL16 ADD KJJ TO INDEX ed) NO YES TEM1 = TEM2

(44)

4

NZ

7 ADD 1 TO LOPC INDEX = LOPC OUT ADD 1 TO 13 LOOP 3 (Continued) SHFT SOD DIVIDE XR & XI BY 2 DIVIDE XR & XI BY RETURN = YES

(45)

INVT C ENTER 4 PICK UP THREE HIGH ORDER BITS FORM ADDRESS OF INVERSE

i

FIND INVERSE SAVE AS THREE LOW ORDER BITS OF INVERSE 4 PICK UP THREE MIDDLE BITS

i

FORM ADDRESS OF INVERSE

i

FIND INVERSE

i

SHIFT LEFT 3

i

)

ADD TO LOW ORDER BITS AND SAVE

i

SUBROUTINE INVT

Y

PICK UP THREE LOW ORDER BITS

i

FORM ADDRESS OF INVERSE

i

FIND INVERSE

i

SHIFT LEFT 6

i

ADD TO SAVED VALUE AS HIGH ORDER BITS

i

C

RETURN

)

(46)

REFERENCES

I. Beard, "Internal Reordering and Multivariate Indexing with the Fast Fourier Transform," Applied Research Laboratories Report ARL-TR-71-24 (1971).

2. Programmers Reference Manual for the H316 and DDP-5I6 General Purpose Computer, Honeywell Computer Control Division, Framingham, Mass (April 1969).

(47)

Copies

INITIAL DISTRIBUTION

CENTER DISTRIBUTION

3 CNO Copies Code

1 0P090

1 0P96

1 00 CAPT P. W. Nelson

1 OP 98T/R. Burns ₁ ₁₈₀₀ _{Gene H. Gleissner}

1 NRL 2627 Lib ₄ ₁₈₀₂

1 USNA Lib 1 Daniel Shanks

1 NAVPGSCOL Lib

1 Francois N. Frenkiel

1 Hans J. Lugt

1 NROTC & NAVADMINU, MIT 1 Feodor Theilheimer

1 NAVWARCOL 1 1805 Elizabeth H. Cuthill

1 NAVELEXSYSCOM 1 1830 Herbert M. Ernst

9 NAVSEA 1 1840 John W. Wrench, Jr

2 SEA 037 2 SEA 09G32

1 1850 Thomas Corin

2 SEA 06H 1 1860 Robert A. Sulit

1 SEA 06H1-1 1 1880 Abel W. Camara

1 SEA 06H1-3 1 SEA 06H2

1 1890 Gilbert R. Gray

2 189 Central Depository

1 NAVAIRDEVCEN ADT Lib 4 1892 Lorraine R. Minor

2 NAVUSEACEN 2 1926 David J. Vendittis

1 1311 Lib 1 13111 Lib

3 1931 Theo Kooij

1 NOL 730 Lib

1 1932 Quentin E. Dolechek 1 NPTLAB NUSC Lib

2 1933 James P. Lee 1 NLONLAB NUSC 1 1935 Gerald R. Clark 1 1935 William G. Niner 12 DDC 30 5614 Reports Distribution 3 U TEXAS ARLJSIGNAL PHYS DIV

1 WHO! OCEAN ENGR DEPT

1 5641 Main Library