10 REM THIS PROGRAM MAKES BINARY HOLOGRAMS
15 REM BY:  PATRICK HAWLEY, PORT ELGIN, ON
20 DIM X(32),Y(32),S(13),U(13),AR(32,32),AI(32,32),CH$(11,15),E(16)
30 TI$="000000":OPEN4,4
40 IR=1: II=1
50 READ N,IT
60 PRINT"[147]"
70 INPUT"DO YOU WANT A HOLOGRAM PLOT (Y/N)";HP$
80 INPUT"RANDOM OR CONSTANT PHASE (R/C)";IP$
90 PQ=13
100 FOR D=1 TO N
110 FOR E=1 TO N
120 READ AR(D,E)
130 AI(D,E)=AR(D,E)*II
140 AR(D,E)=AR(D,E)*IR
150 IFIP$="C"THEN200
160 AR(D,E)=(-1+2*RND(1))*SQR(IR^2+II^2)*AR(D,E)
170 PM=INT(-1+3*RND(1))
180 IF PM=0 THEN 170
190 AI(D,E)=SQR(IR^2+II^2-AR(D,E)^2)*AI(D,E)*PM
200 NEXT E,D
220 DATA 32,1
230 DATA 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0
240 DATA 0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0
250 DATA 0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0
260 DATA 0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0
270 DATA 0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0
280 DATA 0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0
290 DATA 0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0
300 DATA 0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0
310 DATA 0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0
320 DATA 0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0
330 DATA 0,0,1,1,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,1,1,0,0
340 DATA 0,0,1,1,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,1,1,0,0
350 DATA 0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0
360 DATA 0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0
370 DATA 0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0
380 DATA 0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0
390 DATA 0,1,1,0,0,0,1,1,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,1,1,0,0,0,1,1,0
400 DATA 0,1,1,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,1,1,0
410 DATA 0,1,1,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,1,1,0
420 DATA 0,0,1,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,1,0,0
430 DATA 0,0,1,1,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,1,1,0,0
440 DATA 0,0,1,1,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,1,1,0,0
450 DATA 0,0,0,1,1,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,1,1,0,0,0
460 DATA 0,0,0,1,1,0,0,0,0,0,1,1,1,0,0,0,0,0,0,1,1,1,0,0,0,0,0,1,1,0,0,0
470 DATA 0,0,0,0,1,1,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,1,1,0,0,0,0
480 DATA 0,0,0,0,0,1,1,0,0,0,0,0,0,1,1,1,1,1,1,0,0,0,0,0,0,1,1,0,0,0,0,0
490 DATA 0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0
500 DATA 0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0
510 DATA 0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0,0,0,1,1,1,0,0,0,0,0,0,0,0
520 DATA 0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0
530 DATA 0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0
540 DATA 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0
550 :
560 GOSUB 660:REM FAST FOURIER TRANSFORM
570 GOSUB 2130:REM AMPLITUDE AND ANGLE
580 GOSUB 2370:REM NORMALIZE/TRUNCATE
590 GOSUB 2900:REM HOLOGRAM PLOT
600 GOSUB 4340:REM FIND COMPONENTS FOR INVERSE FFT FROM PLOT DATA
610 IT=-1
620 GOSUB 660:REM INVERSE FFT
630 GOSUB 4430:REM GRAY-SCALE PLOT
640 PRINT"FINISHED. TIME IS  ";TI$
650 END
660 REM 2D FFT
670 FOR RW=1 TO N
680 PRINT"  WORKING ON FFT;ROW:";RW
690 FOR RE=1 TO N
700 X(RE)=AR(RW,RE):Y(RE)=AI(RW,RE)
710 NEXT RE
720 GOSUB 890
730 FOR RE=1 TO N
740 AR(RW,RE)=X(RE):AI(RW,RE)=Y(RE)
750 NEXT RE,RW
770 FOR CL=1 TO N
780 PRINT"  WORKING ON FFT;COLUMN:";CL
790 FOR CE=1 TO N
800 X(CE)=AR(CE,CL):Y(CE)=AI(CE,CL)
810 NEXT CE
820 GOSUB 890
830 FOR CE=1 TO N
840 AR(CE,CL)=X(CE):AI(CE,CL)=Y(CE)
850 NEXT CE,CL
870 RETURN
880 :
890 REM  COMPLEX FAST FOURIER TRANSFORM
900 :
910 IF IT>0 THEN 940
920 FOR I=1 TO N
930 Y(I)=-Y(I): NEXT I
940 LG=LOG(N)/LOG(2)
950 IF LG<=1 THEN 1390
960 FOR K=2 TO LG STEP 2
970 M=2^(LG-K)
980 M4=4*M
990 FOR J=1 TO M
1000 AR=(J-1)/M4*2*(null)
1010 C1=COS(AR)
1020 S1=SIN(AR)
1030 C2=C1*C1-S1*S1
1040 S2=C1*S1+C1*S1
1050 C3=C2*C1-S2*S1
1060 S3=C2*S1+S2*C1
1070 FOR I=M4 TO N STEP M4
1080 J0=I+J-M4
1090 J1=J0+M
1100 J2=J1+M
1110 J3=J2+M
1120 R0=X(J0)+X(J2)
1130 R1=X(J0)-X(J2)
1140 I0=Y(J0)+Y(J2)
1150 I1=Y(J0)-Y(J2)
1160 R2=X(J1)+X(J3)
1170 R3=X(J1)-X(J3)
1180 I2=Y(J1)+Y(J3)
1190 I3=Y(J1)-Y(J3)
1200 X(J0)=R0+R2
1210 Y(J0)=I0+I2
1220 IF AR=0 THEN 1300
1230 X(J2)=(R1+I3)*C1+(I1-R3)*S1
1240 Y(J2)=(I1-R3)*C1-(R1+I3)*S1
1250 X(J1)=(R0-R2)*C2+(I0-I2)*S2
1260 Y(J1)=(I0-I2)*C2-(R0-R2)*S2
1270 X(J3)=(R1-I3)*C3+(I1+R3)*S3
1280 Y(J3)=(I1+R3)*C3-(R1-I3)*S3
1290 GOTO 1360
1300 X(J2)=R1+I3
1310 Y(J2)=I1-R3
1320 X(J1)=R0-R2
1330 Y(J1)=I0-I2
1340 X(J3)=R1-I3
1350 Y(J3)=I1+R3
1360 NEXT I,J,K
1390 IF LG=INT(LG/2)*2 THEN 1500
1400 FOR I=1 TO N STEP 2
1410 R0=X(I)+X(I+1)
1420 R1=X(I)-X(I+1)
1430 I0=Y(I)+Y(I+1)
1440 I1=Y(I)-Y(I+1)
1450 X(I)=R0
1460 Y(I)=I0
1470 X(I+1)=R1
1480 Y(I+1)=I1
1490 NEXT I
1500 S(13)=N/2
1510 U(13)=N
1520 FOR K=2 TO 12
1530 J=14-K
1540 S(J)=1
1550 U(J)=S(J+1)
1560 IF S(J+1)>1 THEN S(J)=INT(S(J+1)/2)
1570 NEXT K
1580 AL=S(2)
1590 JJ=0
1600 A=1: B=A: C=B: D=C: E=D: F=E
1660 FOR G=F TO U(7) STEP S(7)
1670 FOR H=G TO U(8) STEP S(8)
1680 FOR I=H TO U(9) STEP S(9)
1690 FOR J=I TO U(10) STEP S(10)
1700 FOR K=J TO U(11) STEP S(11)
1710 FOR L=K TO U(12) STEP S(12)
1720 FOR M=L TO U(13) STEP S(13)
1730 JJ=JJ+1
1740 IF JJ<=M THEN 1810
1750 T=X(JJ)
1760 X(JJ)=X(M)
1770 X(M)=T
1780 T=Y(JJ)
1790 Y(JJ)=Y(M)
1800 Y(M)=T
1810 :
1820 NEXT M,L,K,J,I,H,G
1890 F=F+S(6)
1900 IF F<=U(6) THEN 1660
1910 E=E+S(5)
1920 IF E<=U(5) THEN 1650
1930 D=D+S(4)
1940 IF D<=U(4) THEN 1640
1950 C=C+S(3)
1960 IF C<=U(S) THEN 1630
1970 B=B+S(2)
1980 IF B<=U(2) THEN 1620
1990 A=A+1
2000 IF A<=AL THEN 1610
2010 FOR SL=1 TO (N-2)/2
2020 BC=N+1-SL:FC=SL+1
2030 TX=X(FC):TY=Y(FC)
2040 X(FC)=X(BC):Y(FC)=Y(BC)
2050 X(BC)=TX:Y(BC)=TY
2060 NEXT SL
2070 IF IT>0 THEN RETURN
2080 FOR I=1 TO N
2090 X(I)=X(I)/N
2100 Y(I)=-Y(I)/N
2110 NEXTI
2120 RETURN
2130 REM CONVERT INTO AMPLITUDE & ANGLE
2140 PRINT"...FINDING AMPLITUDE & ANGLE"
2150 FOR D=1 TO N
2160 FOR E=1 TO N
2170 AP=SQR(AR(D,E)*AR(D,E)+AI(D,E)*AI(D,E))
2180 IF AP=0 THEN 2330
2190 IF AR(D,E)<>0 THEN 2230
2200 IF AI(D,E)<0 THEN PA=-(null)/2
2210 IF AI(D,E)>0 THEN PA=(null)/2
2220 GOTO 2320
2230 IF AI(D,E)<>0 THEN 2270
2240 IF AR(D,E)<0 THE NPA=(null)
2250 IF AR(D,E)>0 THENPA=0
2260 GOTO 2320
2270 PA=ATN(AI(D,E)/AR(D,E))
2280 IF AR(D,E)>0 AND AI(D,E)>0 THENPA=PA
2290 IF AR(D,E)<0 AND AI(D,E)>0 THENPA=(null)+PA
2300 IF AR(D,E)<0 AND AI(D,E)<0 THENPA=PA-(null)
2310 IF AR(D,E)>0 AND AI(D,E)<0 THENPA=PA
2320 AR(D,E)=AP: AI(D,E)=PA
2330 NEXT E,D
2350 RETURN
2360 :
2370 REM AMPLITUDE NORMALIZATION ROUTINE
2380 REM FIND LARGEST AMPLITUDE
2390 LA=AR(1,1)
2400 FOR D=1 TO N
2410 FOR E=1 TO N
2420 IFAR(D,E)>LATHENLA=AR(D,E)
2430 NEXT E,D
2450 IF IP$="R" THEN 2820
2460 PRINT"THE LARGEST AMPLITUDE IS";LA
2470 REM CHECK AMPLITUDE DISTRIBUTION
2480 P8=0: P6=0: P4=0: P2=0: P1=0
2490 N8=.8*LA:N6=.6*LA:N4=.4*LA:N2=.2*LA:N1=.1*LA
2500 PRINT" "
2510 FOR D=1 TO N
2520 FOR E=1 TO N
2530 IF AR(D,E)<=N8 THEN P8=P8+1
2540 IF AR(D,E)<=N6 THEN P6=P6+1
2550 IF AR(D,E)<=N4 THEN P4=P4+1
2560 IF AR(D,E)<=N2 THEN P2=P2+1
2570 IF AR(D,E)<=N1 THEN P1=P1+1
2580 NEXT E,D
2600 P8=INT(P8/1024*100):P6=INT(P6/1024*100):P4=INT(P4/1024*100)
2610 P2=INT(P2/1024*100):P1=INT(P1/1024*100)
2620 PRINT"80% OF THE MAX. AMP IS";N8
2630 PRINT P8;"% OF AMP'S ARE EQUAL OR SMALLER"
2640 PRINT" "
2650 PRINT"60% OF THE MAX. AMP IS";N6
2660 PRINTP6;"% OF AMP'S ARE EQUAL OR SMALLER"
2670 PRINT" "
2680 PRINT"40% OF THE MAX. AMP IS";N4
2690 PRINTP4;"% OF AMP'S ARE EQUAL OR SMALLER"
2700 PRINT" "
2710 PRINT"20% OF THE MAX. AMP IS";N2
2720 PRINTP2;"% OF AMP'S ARE EQUAL OR SMALLER"
2730 PRINT" "
2740 PRINT"10% OF THE MAX. AMP IS";N1
2750 PRINTP1;"% OF AMP'S ARE EQUAL OR SMALLER"
2760 PRINT" "
2770 PRINT"TIME IS  ";TI$
2780 PRINT" "
2790 INPUT"WHAT NO. DO YOU WANT TO NORMALIZE TO";NL
2800 PRINT" "
2810 INPUT"WHICH MEANS WHAT % OF DATA TRUNCATED";PT
2820 IFIP$="R"THENNL=LA
2830 FOR D=1 TO N
2840 FOR E=1 TO N
2850 AR(D,E)=AR(D,E)/NL: IF AR(D,E)>1 THEN AR(D,E)=1
2860 NEXT E,D
2880 RETURN
2890 :
2900 REM PLOT ROUTINE
2910 REM READ IN PLOT CHARACTERS
2920 IF HP$="N" THEN 3720
2930 PRINT"...DEFINING PLOT CHARACTERS"
2940 FOR I=1 TO 8:E(I)=2^I-1:E(I+8)=E(I)*2^(8-I):NEXT
2950 FOR I=1 TO 11: FOR J=1 TO 15: FOR K=1 TO 8
2960 CH=0: IF K<(13-I) AND K>(8-I) THEN CH=E(J)
2970 CH$(I,J)=CH$(I,J)+CHR$(CH)
2980 NEXT K,J,I
3710 :
3720 OPEN 2,4,5
3730 OPEN 6,4,6:PRINT#6,CHR$(21)
3740 REM MAIN PLOTTING LOOP
3750 TB=0
3760 NF=2: IF HP$="N" THEN NF=1
3770 FOR D=1 TO N
3780 FOR F=1 TO NF
3790 FOR E=1 TO N
3800 REM DETERMINE WHICH CHARACTER
3810 REM FIND CHARACTER HEIGHT
3820 HT=8*AR(D,E)
3830 DH=HT-INT(HT)
3840 IF DH<.125/2 THEN HT=INT(HT):GOTO3860
3850 HT=INT(HT)+1
3860 IF F=NF THEN AR(D,E)=HT/8*NL
3870 REM FIND PHASE
3880 CP=PQ*AI(D,E)/(2*(null))
3890 QP=INT(PQ/2)
3900 IFABS(CP)-INT(ABS(CP))>.5THEN3950
3910 IF CP<0 THEN CP=INT(CP):IF F=NF THEN AI(D,E)=(CP+1)*(null)/QP
3920 IF CP>=0 THEN CP=INT(CP)+1:IF F=NF THEN AI(D,E)=(CP-1)*(null)/QP
3930 IF HP$="N"THEN 4290
3940 GOTO 3980
3950 IF CP<0 THEN CP=INT(CP)-1:IF F=NF THEN AI(D,E)=(CP+1)*(null)/QP
3960 IF CP>0THEN CP=INT(CP)+2:IF F=NFTHEN AI(D,E)=(CP-1)*(null)/QP
3970 IF HP$="N" THEN 4290
3980 IF F<>1 THEN 4230
3990 IF HT=0 THEN PRINT#4,"  ";CHR$(141);:GOTO4260
4000 IF ABS(CP)>2THEN 4130
4010 IF CP>-2 THEN 4050
4020 PRINT#2,CH$(3,HT):PRINT#4,CHR$(254);CHR$(141);
4030 TB=TB+1:PRINT#2,CH$(11,HT):PRINT#4,TAB(TB);CHR$(254);CHR$(141);
4040 GOTO 4260
4050 IF CP>1THEN 4090
4060 PRINT#2,CH$(2,HT):PRINT#4,CHR$(254);CHR$(141);
4070 TB=TB+1:PRINT#2,CH$(10,HT):PRINT#4,TAB(TB);CHR$(254);CHR$(141);
4080 GOTO 4260
4090 PRINT#2,CH$(1,HT):PRINT#4,CHR$(254);CHR$(141);
4100 TB=TB+1:PRINT#2,CH$(9,HT):PRINT#4,TAB(TB);CHR$(254);CHR$(141);
4110 TB=TB+1:PRINT#2,CH$(9,HT):PRINT#4,TAB(TB);CHR$(254);CHR$(141);
4120 GOTO 4260
4130 IF CP<0 THEN 4180
4140 CP=CP+11-CP*2
4150 PRINT#4," ";CHR$(141);
4160 TB=TB+1:PRINT#2,CH$(CP,HT):PRINT#4,TAB(TB);CHR$(254);CHR$(141);
4170 GOTO4260
4180 IF CP>0 THEN 4230
4190 CP=ABS(CP)+1
4200 PRINT#2,CH$(CP,HT):PRINT#4,CHR$(254);CHR$(141);
4210 TB=TB+1:PRINT#4,TAB(TB);" ";CHR$(141);
4220 GOTO 4260
4230 IF HT=0 THEN 3990
4240 IF HT<8 THEN HT=HT+8
4250 GOTO 3990
4260 IF E=N THEN PRINT#4,CHR$(141);CHR$(13);
4270 TB=2*E:IF E=N THEN TB=0
4280 IF E<>N THEN PRINT#4,TAB(TB);
4290 NEXT E,F,D
4320 RETURN
4330 :
4340 REM FIND REAL AND IMAGINARY COMPONENTS FROM TRUNCATED & QUANTIZED DATA
4350 FOR D=1 TO N: FOR E=1 TO N
4370 AP=AR(D,E)
4380 AR(D,E)=AP*COS(AI(D,E))
4390 AI(D,E)=AP*SIN(AI(D,E))
4400 NEXT E,D
4420 RETURN
4430 :
4440 REM GRAY SCALE PICTURE ROUTINE
4450 PRINT#4," "
4460 FOR D=1 TO N
4470 PRINT#4,TAB(23);
4480 FOR E=1 TO N
4490 IN=SQR(AR(D,E)*AR(D,E)+AI(D,E)*AI(D,E))
4500 AF=1/SQR(IR*IR+II*II):IN=AF*IN
4510 IF IN>1 THEN IN=1
4520 IF IN<.1 THENPRINT#4," ";:GOTO4590
4530 IF IN<.2 THENPRINT#4,"'";:GOTO4590
4540 IF IN<.4 THENPRINT#4,":";:GOTO4590
4550 IF IN<.6 THENPRINT#4,"*";:GOTO4590
4560 IF IN<.8 THENPRINT#4,"[166]";:GOTO4590
4570 IF IN<.9 THENPRINT#4,"*[146]";:GOTO4590
4580 IF IN<=1 THENPRINT#4," [146]";
4590 IF E=N THEN PRINT#4,CHR$(13);
4600 NEXT E,D
4620 PRINT#4," "
4630 PRINT#4,"INPUT REAL AMPLITUDE:";IR
4640 PRINT#4,"INPUT IMAGINARY AMPLITUDE:";II
4650 PRINT#4,"PHASE ANGLE:";IP$
4660 PRINT#4,"NO. OF PHASE QUANTIZATION LEVELS:";PQ
4670 IF IP$="R" THEN 4700
4680 PRINT#4,"TRANSFORM NORMALIZED WITH...";NL
4690 PRINT#4,"PERCENT OF TRANSFORMED AMP'S TRUNCATED...";PT
4700 RETURN