home *** CD-ROM | disk | FTP | other *** search
- PROGRAM d12r1(input,output);
- (* driver for routine FOUR1 *)
- CONST
- nn=32;
- nn2=64; (* 2*nn *)
- TYPE
- gldarray = ARRAY [1..nn2] OF real;
- VAR
- ii,i,isign : integer;
- data,dcmp : gldarray;
-
- PROCEDURE prntft(data : gldarray; nn : integer);
- VAR
- ii,mm,n : integer;
- BEGIN
- writeln('n':4,'real(n)':13,'imag.(n)':13,'real(N-n)':12,'imag.(N-n)':13);
- writeln (0:4,data[1]:14:6,data[2]:12:6,data[1]:12:6,data[2]:12:6);
- mm := nn DIV 2;
- FOR ii := 1 to mm DO BEGIN
- n := 2*ii+1;
- writeln (((n-1) DIV 2):4,data[n]:14:6,data[n+1]:12:6,
- data[2*nn+2-n]:12:6,data[2*nn+3-n]:12:6)
- END;
- writeln (' press return to continue ...');
- readln
- END;
-
- (*$I MODFILE.PAS *)
- (*$I FOUR1.PAS *)
-
- BEGIN
- writeln('h(t) := real-valued even-function');
- writeln('h(n) := h(N-n) and real?');
- FOR ii := 1 to nn DO BEGIN
- i := 2*ii-1;
- data[i] := 1.0/(sqr((i-nn-1.0)/nn)+1.0);
- data[i+1] := 0.0
- END;
- isign := 1;
- four1(data,nn,isign);
- prntft(data,nn);
- writeln('h(t) := imaginary-valued even-function');
- writeln('h(n) := h(N-n) and imaginary?');
- FOR ii := 1 to nn DO BEGIN
- i := 2*ii-1;
- data[i+1] := 1.0/(sqr((i-nn-1.0)/nn)+1.0);
- data[i] := 0.0
- END;
- isign := 1;
- four1(data,nn,isign);
- prntft(data,nn);
- writeln('h(t) := real-valued odd-function');
- writeln('h(n) := -h(N-n) and imaginary?');
- FOR ii := 1 to nn DO BEGIN
- i := 2*ii-1;
- data[i] := ((i-nn-1.0)/nn)/(sqr((i-nn-1.0)/nn)+1.0);
- data[i+1] := 0.0
- END;
- data[1] := 0.0;
- isign := 1;
- four1(data,nn,isign);
- prntft(data,nn);
- writeln('h(t) := imaginary-valued odd-function');
- writeln('h(n) := -h(N-n) and real?');
- FOR ii := 1 to nn DO BEGIN
- i := 2*ii-1;
- data[i+1] := ((i-nn-1.0)/nn)/(sqr((i-nn-1.0)/nn)+1.0);
- data[i] := 0.0
- END;
- data[2] := 0.0;
- isign := 1;
- four1(data,nn,isign);
- prntft(data,nn);
- (* transform, inverse-transform test *)
- FOR ii := 1 to nn DO BEGIN
- i := 2*ii-1;
- data[i] := 1.0/(sqr(0.5*(i-nn-1.0)/nn)+1.0);
- dcmp[i] := data[i];
- data[i+1] := (0.25*(i-nn-1.0)/nn)
- *exp(-sqr(0.5*(i-nn-1.0)/nn));
- dcmp[i+1] := data[i+1]
- END;
- isign := 1;
- four1(data,nn,isign);
- isign := -1;
- four1(data,nn,isign);
- writeln;
- writeln('double fourier transform:':33,'original data:':23);
- writeln;
- writeln('k':4,'real h(k)':14,'imag h(k)':13,
- 'real h(k)':17,'imag h(k)':12);
- FOR ii := 1 to (nn DIV 2) DO BEGIN
- i := 2*ii-1;
- writeln(((i+1) DIV 2):4,dcmp[i]:14:6,dcmp[i+1]:12:6,
- data[i]/nn:17:6,data[i+1]/nn:12:6)
- END
- END.
-