| << Chapter < Page | Chapter >> Page > |
1 ; v:\ece420\54x\dspclib\c_fft_given.asm
2 ; dgs - 9/14/2001
3 .copy "v:\ece420\54x\dspclib\core.inc"
4
5 .global _bit_rev_data
6 .global _fft_data
7 .global _window
8
9 .global _bit_rev_fft
10
11 .sect ".data"
12
13 .align 4*N
14 _bit_rev_data .space 16*2*N ; Input to _bit_rev_fft
15
16 .align 4*N
17 _fft_data .space 16*2*N ; FFT output buffer
18
19
20 ; Copy in the Hamming window
21 _window ; The Hamming window
22 .copy "window.asm"
23
24 .sect ".text"
25
26 _bit_rev_fft
27 ENTER_ASM
28
29 call bit_rev ; Do the bit-reversal.
30
31 call fft ; Do the FFT
32
33 LEAVE_ASM
34 RET
35
36 bit_rev:
37 STM #_bit_rev_data,AR3 ; AR3 -> original input
38 STM #_fft_data,AR7 ; AR7 -> data processing buffer
39 MVMM AR7,AR2 ; AR2 -> bit-reversed data
40 STM #K_FFT_SIZE-1,BRC
41 RPTBD bit_rev_end-1
42 STM #K_FFT_SIZE,AR0 ; AR0 = 1/2 size of circ buffer
43 MVDD *AR3+,*AR2+
44 MVDD *AR3-,*AR2+
45 MAR *AR3+0B
46 bit_rev_end:
47 NOP
48 RET
49
50 ; Copy the actual FFT subroutine.
51 fft_data .set _fft_data ; FFT code needs this.
52 .copy "v:\ece420\54x\dsplib\fft.asm"
53
54
55 ; If you need any more assembly subroutines, make sure you name them
56 ; _name, and include a ".global _name" directive at the top. Also,
57 ; don't forget to use ENTER_ASM at the beginning, and LEAVE_ASM
58 ; and RET at the end!
1 #define N 1024 /* Number of FFT points */
2 #define logN 10
1 /*****************************************************************/
2 /* lab4fft.c */
3 /* Douglas L. Jones */
4 /* University of Illinois at Urbana-Champaign */
5 /* January 19, 1992 */
6 /* Changed for use w/ short integers and lookup table for ECE420 */
7 /* Matt Kleffner */
8 /* February 10, 2004 */
9 /* */
10 /* fft: in-place radix-2 DIT DFT of a complex input */
11 /* */
12 /* Permission to copy and use this program is granted */
13 /* as long as this header is included. */
14 /* */
15 /* WARNING: */
16 /* This file is intended for educational use only, since most */
17 /* manufacturers provide hand-tuned libraries which typically */
18 /* include the fastest fft routine for their DSP/processor */
19 /* architectures. High-quality, open-source fft routines */
20 /* written in C (and included in MATLAB) can be found at */
21 /* http://www.fftw.org */
22 /* */
23 /* #defines expected in lab4.h */
24 /* N: length of FFT: must be a power of two */
25 /* logN: N = 2**logN */
26 /* */
27 /* 16-bit-limited input/output (must be defined elsewhere) */
28 /* real: integer array of length N with real part of data */
29 /* imag: integer array of length N with imag part of data */
30 /* */
31 /* sinetables.h must */
32 /* 1) #define Nt to an equal or greater power of two than N */
33 /* 2) contain the following integer arrays with */
34 /* element magnitudes bounded by M = 2**15-1: */
35 /* costable: M*cos(-2*pi*n/Nt), n=0,1,...,Nt/2-1 */
36 /* sintable: M*sin(-2*pi*n/Nt), n=0,1,...,Nt/2-1 */
37 /* */
38 /*****************************************************************/
39
40 #include "lab4.h"
41 #include "sinetables.h"
42
43 extern int real[N];
44 extern int imag[N];
45
46 void fft(void)
47 {
48 int i,j,k,n1,n2,n3;
49 int c,s,a,t,Wr,Wi;
50
51 j = 0; /* bit-reverse */
52 n2 = N >> 1;
53 for (i=1; i < N - 1; i++)
54 {
55 n1 = n2;
56 while ( j >= n1 )
57 {
58 j = j - n1;
59 n1 = n1 >> 1;
60 }
61 j = j + n1;
62
63 if (i < j)
64 {
65 t = real[i];
66 real[i] = real[j];
67 real[j] = t;
68 t = imag[i];
69 imag[i] = imag[j];
70 imag[j] = t;
71 }
72 }
73
74 /* FFT */
75 n2 = 1; n3 = Nt;
76
77 for (i=0; i < logN; i++)
78 {
79 n1 = n2; /* n1 = 2**i */
80 n2 = n2 + n2; /* n2 = 2**(i+1) */
81 n3 = n3 >> 1; /* cos/sin arg of -6.283185307179586/n2 */
82 a = 0;
83
84 for (j=0; j < n1; j++)
85 {
86 c = costable[a];
87 s = sintable[a];
88 a = a + n3;
89
90 for (k=j; k < N; k=k+n2)
91 {
92 /* Code for standard 32-bit hardware, */
93 /* with real,imag limited to 16 bits */
94 /*
95 Wr = (c*real[k+n1] - s*imag[k+n1]) >> 15;
96 Wi = (s*real[k+n1] + c*imag[k+n1]) >> 15;
97 real[k+n1] = (real[k] - Wr) >> 1;
98 imag[k+n1] = (imag[k] - Wi) >> 1;
99 real[k] = (real[k] + Wr) >> 1;
100 imag[k] = (imag[k] + Wi) >> 1;
101 */
102 /* End standard 32-bit code */
103
104 /* Code for TI TMS320C54X series */
105
106 Wr = ((long int)(c*real[k+n1]) - (long int)(s*imag[k+n1])) >> 15;
107 Wi = ((long int)(s*real[k+n1]) + (long int)(c*imag[k+n1])) >> 15;
108 real[k+n1] = ((long int)real[k] - (long int)Wr) >> 1;
109 imag[k+n1] = ((long int)imag[k] - (long int)Wi) >> 1;
110 real[k] = ((long int)real[k] + (long int)Wr) >> 1;
111 imag[k] = ((long int)imag[k] + (long int)Wi) >> 1;
112
113 /* End code for TMS320C54X series */
114
115 /* Intrinsic code for TMS320C54X series */
116 /*
117 Wr = _ssub(_smpy(c, real[k+n1]), _smpy(s, imag[k+n1]));
118 Wi = _sadd(_smpy(s, real[k+n1]), _smpy(c, imag[k+n1]));
119 real[k+n1] = _sshl(_ssub(real[k], Wr),-1);
120 imag[k+n1] = _sshl(_ssub(imag[k], Wi),-1);
121 real[k] = _sshl(_sadd(real[k], Wr),-1);
122 imag[k] = _sshl(_sadd(imag[k], Wi),-1);
123 */
124 /* End intrinsic code for TMS320C54X series */
125 }
126 }
127 }
128 return;
129 }
Notification Switch
Would you like to follow the 'Digital signal processing laboratory (ece 420)' conversation and receive update notifications?
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|