home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
HAM Radio 1
/
HamRadio.cdr
/
tech
/
design2
/
sabin.asc
< prev
next >
Wrap
Text File
|
1986-10-29
|
12KB
|
289 lines
10 PRINT:PRINT "WAVEFORM ANALYSIS:"
20 PRINT "W.E. SABIN - EDN JUNE 1983 - PAGE 243"
30 INPUT "ENTER EXPONENT M ";M
40 N=2^M: PI=3.1415928159#:E1=2
50 DIM X(N+1,5)
60 PRINT: INPUT "1 TO READ A DATA FILE, 2 TO GO ON";A
70 ON A GOTO 80,110
80 INPUT "FILE TO READ";F$
81 OPEN "I",1,F$
82 FOR I=1 TO N
83 INPUT #1,X(I,0)
84 NEXT I
85 CLOSE
100 PRINT "DATA LOADED":PRINT:GOTO 280
110 REM ****************************************************************
120 REM * X(I,0) REAL, X(I,1) IMAGINARY *
130 REM * EVALUATE REAL, IMAGINARY IN LINES 190 - 269 *
140 REM * FOR AUTOCORRELATION, AUTOSPEC. USE X(I,0), X(I,1) *
150 REM * FOR CROSS SPECTRUM, CROSS CORRELATION *
160 REM * AND CONVOLUTION, *
170 REM * USE X(I,0),X(I,1) AND X(I,2),X(I,3) *
180 REM ****************************************************************
190 REM
200 FOR I=0 TO N
210 IF I<N/2 THEN X(I,0)=1 ELSE X(I,0)=0
220 X(I,1)=0
230 NEXT I
240 REM
250 REM
260 REM
270 PRINT
280 PRINT: PRINT "ITEMS 1-7 FOR X(I,0), X(I,1) ONLY":PRINT
290 PRINT "TYPE 1 FOR FORWARD TRANSFORM"
300 PRINT "TYPE 2 FOR INVERSE TRANSFORM"
310 PRINT "TYPE 3 TO LIST X REAL, X IMAG"
320 PRINT "TYPE 4 FOR SINE/COSINE"
330 PRINT "TYPE 5 FOR MAG AND PHASE"
340 PRINT "TYPE 6 FOR SMOOTHING"
350 PRINT "TYPE 7 FOR WINDOWING"
360 PRINT "TYPE 8 FOR POWER SPECTRUM"
370 PRINT "TYPE 9 FOR CORRELATION"
380 PRINT "TYPE 10 FOR CONVOLUTION"
390 PRINT "TYPE 11 FOR TO SAVE DATA IN FILE"
400 PRINT "TYPE 12 TO END"
410 PRINT "TYPE 13 TO EXCHANGE SEQ 1 & 2"
420 PRINT "TYPE 14 FOR DEQSEQ(1)=SEQ(1)*SEQ(2)"
430 INPUT X
440 ON X GOTO 470,1410,500,600,600,1450,1690,1790,2000,2310,1230,2770,460,450
450 FOR I=1 TO N:X(I,4)=X(I,0)*X(I,2)-X(I,1)*X(I,3)451 X(I,5)=X(I,0)*X(I,3)+X(I,2)*X(I,1):X(I,0)=X(I,4)452 X(I,1)=X(I,5):NEXT:GOTO 270
460 FOR I=1 TO N:X(I,4)=X(I,0):X(I,5)=X(I,1):X(I,0)=X(I,2):X(I,1)=X(I,3)461 X(I,2)=X(I,4):X(I,3)=X(I,5):NEXT:GOTO 270
470 REM ******** COMPUTE FORWARD TRANSFORM **********
480 D=0:GOSUB 2480
490 GOTO 270
500 REM ******** DISPLAY EXPONENTIALS **************
510 PRINT:PRINT "N";TAB(6);"XR";TAB(25);"XI":PRINT
520 FOR J=O0 TO N STEP 10:PRINT
530 FOR K=1 TO 10
540 IF J+K>N GOTO 570
550 G=(2-E1)*N/2
560 PRINT J+K-G-1;TAB(6)X(J+K,0);TAB(25)X(J+K,1)
570 NEXT K
580 PRINT: PRINT "PRESS 'Q' TO QUIT, ANY OTHER KEY TO CONTINUE"
581 D$=INPUT$(1): IF D$="Q" THEN GOTO 270
582 PRINT:NEXT J
590 GOTO 270
600 REM ********** CALCULATE SINE, COSINE ****************
610 D=0:GOSUB 2480
620 X(1,2)=X(1,0):X(1,4)=X(1,1)
630 FOR I=2 TO N/2
640 X(I,2)=X(I,0)+X(N+2-I,0)
650 X(I,3)=X(N+2-I,1)-X(I,1)
660 X(I,4)=X(N+2-I,1)+X(I,1)
670 X(I,5)=X(I,0)-X(N+2-I,0)
680 NEXT I
690 REM ******* PRINT SIN/COS OF FIND MAG ***********
700 IF X=5 GOTO 880
710 PRINT:PRINT "REAL PART SINE, COSINE":PRINT
720 FOR J=0 TO N/2 STEP 10
730 PRINT "N" TAB(6)"COSINE" TAB(25)"SINE":PRINT
740 FOR K=1 TO 10
750 IF J+K>N/2 GOTO 270
760 PRINT J+K-1 TAB(6)X(J+K,2) TAB(25)X(J+K,3)
770 NEXT K
780 GET E$:PRINT:NEXT J
790 PRINT:PRINT "IMAG PART, SINE, COSINE:P":PRINT
800 FOR J=0 TO N/2 STEP 10
810 PRINT "N" TAB(6)"COSINE"TAB(325)"SINE":PRINT
820 FOR K=1 TO 10
830 IF J+K>N/2 GOTO 850
840 PRINT J+K-1 TAB(6)X(J+K),4)TAB(25,)X(J+K,5)
850 NEXT K
860 PRINT:INPUT "PRESS <RET> TO CONTINUE:";D$:PRINT:NEXT J:GOTO 270
870 REM ******** EMD SINE COSINE ROUTINE **********
880 REM ****** ** START AMPLITUDE, PHASE **********
890 FOR I=1 TO N/2
900 V=SQR(X(I,2)^2+X(I,3)^2)
910 EW=DQSQR(X(I,4)^2+X(I,5)^2)
920 IF ABS(X(I,2))<1E-12 THEN X(I,2)=1E-12
930 IF ABS(X(I,4))<1E-12 THEN X(I,4)=1E-12
940 Y=ATN(X(I,3)/X(I,2))*57.2958
950 IF X(I,2)<0 AND X(I,3)>0 THEN Y=Y+180
960 IF X(I,2)<0 AND X(I,3)<0 THEN Y=Y-180
970 X(I,2)=V:X(I,3)=Y
980 Y=ATN(X(I,5)/X(I,4))*57.2958
990 IF X(I,4)<0 AND X(I,5)>0 THEN Y=Y+180
1000 IF X(I,4)<0 AND X(I,5)<0 THEN Y=Y-180
1010 X(I,4)=W:X(I,5)=Y
1020 NEXT I
1030 PRINT CHR$(26):PRINT:PRINT "MAGNITUDE AND PHASE":PRINT
1040 PRINT "REAL PART":PRINT
1050 FOR J=0 TO N/2 STEP 10
1060 PRINT "N" TAB(6)"MAG" TAB(25) "PHASE":PRINT
1070 FOR K=1 TO 10
1080 IF J+K>N/2 GOTO 1100
1090 PRINT J+K-1 TAB(6)X(J+K,2) TAB(25)X(J+K,3)
1100 NEXT K
1110 PRINT:INPUT "PRESS <RET> TO CONTINUE";D$:PRINT:NEXT J
1120 PRINT: PRINT "IMAGINARY PART::":PRINT
1130 FOR J=0 TO N/2 STEP 10
1140 PRINT "N" TAB(6) "MAG" TAB(25) "PHASE":PRINT
1150 FOR K=1 TO 10
1160 IF J+K>N/2 GOTO 1180
1170 PRINT J+K-1 TAB(6) X(J+K,4) TAB(25) X(J+K,5)
1180 NEXT K
1190 PRINT:INPUT "PRESS <RET> TO CONTINUE";D$:PRINT:NEXT J
1200 GOTO 270
1210 REM ******* END MAG,PHASE ******
1220 REM ****** OUTPUT DATA FILE ****
1230 PRINT CHR$(26):PRINT:PRINT "INSTRUCTIONS TO SAVE DATA IN DATA FILE"
1240 PRINT "SAVE X OR MAGNITUDE OR PHASE":PRINT
1250 PRINT "X IS X(N),X(K) , SPEC.,CONV.,CORR.":PRINT
1260 PRINT "0=X(I,0) REAL"
1270 PRINT "1=X(I,1) IMAGINARY"
1280 PRINT "2=MAGNITUDE, REAL PART"
1290 PRINT "3=PHASE, REAL PART"
1300 PRINT "4=MAGNITUDE, IMAGINARY PART"
1310 PRINT "5=PHASE, IMAGINARY PART"
1320 INPUT R
1330 Q=1
1340 IF R>1 THEN Q=2
1350 DIM A(N/Q)
1360 FOR I=1 TO N/Q:A(I)=X(I,R):NEXT
1370 INPUT "NAME OF FILE TO SAVE DATA";F$
1371 OPEN "O",1,F$
1372 FOR I=1 TO N/Q
1373 PRINT #1,A(I)
1374 NEXT I
1375 CLOSE #1
1380 PRINT CHR$(13):PRINT "DATA FILED"
1390 GOTO 270
1400 REM ********* END DATA OUTPUT ************
1410 REM ********* INVERSE DFS ************
1420 D=1:GOSUB 2480
1430 GOTO 270
1440 REM ******** SMOOTHING ************
1450 PRINT CHR$(26):PRINT:PRINT "SEQUENCE SMOOTHING":PRINT
1460 PRINT "TYPE 1 FOR LOW PASS"
1470 PRINT "TYPE 2 FOR HIGH PASS"
1480 INPUT Z
1490 ON Z GOTO 1500, 1590
1500 X(1,5)=.25*X(N,0)+.5*X(I1,0)+.25*X(2,0)
1510 X(N,5)=.25*X(N-1,0)+.5*X(1,0)+.25*X(2,0)
1520 FOR J=2 TO N-1:X(J,5)=.25*X(J-1,0)+.5*X(J,0)+.25*X(J+1,0):NEXT
1530 FOR J=1 TO N:X(J,0)=X(J,5):NEXT
1540 X(1,T5)=.25*X(N,1)+.5*JX(1,1)+.25*X(2,1)
1550 X(N,5)=.25*X(N-1,1)+.5*X(N,1)+.25*X(1,1)
1560 FOR J=2 TO N-1:X(J,5)=.25*X(J-1,0)+.5*X(J,1)+.25*X(J+1,1):NEXT J
1570 FOR J=1 TO N:X(J,1)=X(J,5):NEXT J
1580 GOTO 270
1590 X(1,5)=-.25*X(N,0)+.5*X(1,0)-.25*X(2,0)
1600 X(N,5)=-.25*X(N-1,0)+.5*X(N,0)-.25*X(1,0)
1610 FOR J=2 TO N-1:X(J,5)=-.25*X(J-1,1)+.5*X(J,0)-.25*X(J+1,0):NEXT
1620 FOR J=1 TO N:X(J,0)=X(J,5):NEXT J
1630 X(1,5)=-.25*X(N,1)+.5*X(1,1)-.25*X(2,1)
1640 X(N,5)=-.25*X(N-1,1)+.5*X(N,1)-.25*X(2,1)
1650 FOR J=2 TO N-1:X(J,5)=-.25*X(J-1,1)+.5*X(J,1)-.25*X(J+1,1):NEXT J
1660 FOR J=1 TO N:X(J,1)=X(J,5):NEXT J
1670 GOTO 270
1680 REM ******** WINDOWING ************
1690 PRINT CHR$(26):PRINT:PRINT: PRINT "SEQUENCE WINDOW":PRINT
1700 PRINT "TYPE 1 FOR HANNING"
1710 PRINT "TYPE 2 FOR HAMMING"
1720 INPUT Q
1730 IF Q=2 THEN Q=.8519
1740 FOR I=1 TO N
1750 X(I,0)=X(I,0)*(1-Q*COS(2*PI*(I-1)/N))
1760 X(I,1)=X(I,1)*(1-Q*COS (2*PI*(I-1)/N))
1770 NEXT I
1780 GOTO 270
1790 REM ******** POWER SPECTRUM ***********
1800 REM ******* AUTO SPEC. USE X(I,0),X(I,1) ***********
1810 REM ****** CROSS SPEC. USE X(I,0),X(I,1) AND X(I,2),X(I,3) ********
1820 PRINT CHR$(26):PRINT:PRINT "TWO-SIDED POWER SPECTRUM":PRINT
1830 PRINT "TYPE 1 FOR SPECTRUM OF SEQUENCE ONE"
1840 PRINT "TYPE 2 FOR SPECTRUM OF SEQUENCE TWO"
1850 PRINT "TYPE 3 FOR CROSS SPECTRUM"
1860 INPUT F:PRINT
1870 ON F GOTO 1890,1880,1920
1880 FOR I=1 TO N: X(I,0)=X(I,2):X(I,1)=X(I,3):NEXT
1890 D=0:GOSUB 2480
1900 FOR I=1 TO N:X(I,0)=X(I,0)*X(I,0)+X(I,1)*X(I,1):X(I,1)=0:NEXT
1910 PRINT:PRINT "AUTOSPECTRUM":PRINT:GOTO 1980
1920 FOR I=1 TO N:X(I,4)=X(I,0):X(I,5)=X(I,1):X(I,0)=X(I,2):X(I,1)=X(I,3):NEXT
1930 D=0:GOSUB 2480
1940 FOR I=1 TO N:X(I,2)=X(I,0):X(I,3)=X(I,1):X(I,0)=X(I,4):X(I,1)=X(I,5):NEXT
1950 D=0:GOSUB 2480
1960 FOR I=1 TO N:X(I,4)=X(I,0)*X(I,2)+X(I,1)*X(I,3)1961 X(I,5)=X(I,1)*X(I,2)-X(I,0)*X(I,3):X(I,0)=X(I,4):X(I,1)=X(I,5):NEXT
1970 PRINT:PRINT "CROSS SPECTRUM":PRINT
1980 PRINT "TYPE 3 TO LIST POWER SPECTRUM":PRINT:GOTO 290
1990 REM ******* CORRELATION *************
2000 PRINT CHR$(26):PRINT:PRINT "CORRELATION"
2010 PRINT: PRINT "FOR LINEAR CORRELATION, DOUBLE THE VALUE OF N AND FILL IN ";"ZEROS FROM N/2 TO N-1 IN X(N) BEFORE PROCEEDING"
2020 PRINT:PRINT "TYPE 1 FOR AUTOCORRELATION OF SEQUENCE X(N) IN X(I,0),X(I,1)"
2030 PRINT "TYPE 2 FOR AUTOCORRELATION OF SEQUENCE X(N) IN X(I,2),X(I,3)."
2040 PRINT "TYPE 3 FOR CROSS-CORRELATION"
2050 INPUT C
2060 PRINT:PRINT "TYPE 1 FOR CORRELATION:"
2070 PRINT "TYPE 2 FOR COVARIANCE."
2080 INPUT Q
2090 PRINT:PRINT "TYPE 1 IF LINEAR"
2100 PRINT "TYPE 2 IF CIRCULAR"
2110 INPUT E1
2120 IF Q=2 THEN GOSUB 2280
2130 ON C GOTO 2150,2140,2160
2140 FOR I=1 TO N:X(I,0)=X(I,2):X(I,1)=X(I,3):NEXT
2150 D=0:GOSUB 2480:GOSUB 2200:D=1:GOSUB 2480:GOSUB 2210:GOTO 2170
2160 D=0:GOSUB 2480:GOSUB 2260:D=0:GOSUB 2480:GOSUB 2270:D=1:GOSUB 2480: GOSUB 2210
2170 PRINT: ON Q GOTO 2180,2190
2180 PRINT "TYPE 3 TO LIST CORRELATION":PRINT:GOTO 290
2190 PRINT "TYPE 3 TO LIST COVARIANCE":PRINT:GOTO 290
2200 FOR I=1 TO N:X(I,0)=X(I,0)^2+X(I,1)^2:X(I,1)=0:NEXT
2210 IF E1=2 GOTO 2250
2220 FOR I=1 TO N:X(I,0)=X(I,0)*2:X(I,1)=X(I,1)*2:NEXT
2230 FOR I=1 TO N/2:X(I+N/2,4)=X(I,0):X(I,4)=X(I+N/2,0):X(I+N/2,5)=X(I,1):X(I,5)=X(I+N/2,1):NEXT
2240 FOR I=1 TO N:X(I,0)=X(I,4):X(I,1)=X(I,5):NEXT
2250 RETURN
2260 FOR I=1 TO N:X(I,4)=X(I,0):X(I,5)=X(I,1):X(I,0)=X(I,2):X(I,1)=X(I,3):NEXT:RETURN
2270 FOR I=1 TO N:X(I,2)=X(I,0)*X(I,4)+X(I,1)*X(I,5):X(I,3)=X(I,0)*X(I,5)-X(I,1)*X(I,4):X(I,0)=X(I,2):X(I,1)=X(I,3):NEXT:RETURN
2280 U=N/(3-E1):AA=0:BB=0:CC=0:DD=0
2290 FOR I=1 TO U:AA=AA+X(I,0):BB=BB+X(I,1):CC=CC+X(I,2):DD=DD+X(I,3):NEXT
2300 FOR I=1 TO U:X(I,0)=X(I,0)-AA/U:X(I,1)=X(I,1)-BB/U:X(I,2)=X(I,2)-CC/U:X(I,3)=X(I,3)=-DD/U:NEXT:RETURN
2310 REM ********* CONVOLUTION ************
2320 PRINT CHR$(26):PRINT:PRINT "CONVOLUTION":PRINT
2330 PRINT "SEQUENCE 1 IN X(I,0),X(I,1)"
2340 PRINT "SEQUENCE 2 IN X(I,3),X(I,4)":PRINT
2350 PRINT "FOR LINEAR CONVOLUTION, DOUBLE THE VALUE OF N AND ";"ARGUMENT WITH ZEROS IN BOTH SEQUENCES"
2360 PRINT:PRINT "TYPE 1 IF LINEAR"
2370 PRINT "TYPE 2 IF CIRCULAR"
2380 INPUT QQ
2390 D=0:GOSUB 2480:GOSUB 2440:GOSUB 2480:GOSUB 2450:D=1:GOSUB 2480:GOSUB 2460
2400 PRINT:PRINT TYPE 1 TO "TYPE 1 TO MULTIPLY FOR N"
2410 GET A$:PRINT A$
2420 IF A$="1" THEN: FOR I=1 TO N:X(I,0)=X(I,0)*N:X(I,1)=X(I,1)*N:NEXT
2430 PRINT:PRINT "TYPE 3 TO LIST CONVOLUTION OF X(1,N) AND X2(N).":PRINT: GOTO 290
2440 FOR I=1 TO N:X(I,4)=X(I,0):X(I,5)=X(I,1):X(I,0)=X(I,2):X(I,1)=X(I,3):NEXT:RETURN
2450 FOR I=1 TO N:X(I,2)=X(I,0)*X(I,4)-X(I,1)*X(I,5):X(I,3)=X(I,0)*X(I,5)+X(I,4)*X(I,1):X(I,0)=X(I,2):X(I,1)=X(I,3):NEXT:RETURN
2460 IF QQ=1 THEN: FOR I=1 TO N:X(I,0)=2*X(I,0):X(I,1)=2*X(I,1):NEXT:RETURN
2470 RETURN
2480 REM ***** FFT ROUTINE, COMPLEX DATA ARRAY ******
2490 REM ****** X(I,0) REAL, X(I,1) IMAGINARY ******
2500 REM ****** D=0, FORWARD. D-=1, REVERSE **********
2510 N2=N/2:N1=N-1:J=1
2520 FOR I=1 TO N1
2530 IF I>=J THEN 2550
2540 T1=X(J,0):T2=X(J,1):X(J,0)=X(I,0):X(J,1)=X(I,1):X(I,0)=T1:X(I,1)=T2
2550 K=N2
2560 IF K>=J THEN 2590
2570 J=J-K:K=K/2
2580 GOTO 2560
2590 J=J+K
2600 NEXT I
2610 S1=-1
2620 IF D=0 THEN 2640
2630 S1=1
2640 FOR L=1 TO M
2650 L1=2^L:L2=L1/2:U1=1:U2=0:W1=COS(PI/L2):W2=S1*SIN(PI/L2)
2660 FOR J=1 TO L2
2670 FOR I=J TO N STEP L1
2680 I1=I+L2
2690 V1=X(I1,0)*U1-X(I1,1)*U2:V2=X(I1,1)*U1+X(I1,0)*U2
2700 X(I1,0)=X(I,0)-V1:X(I1,1)=X(I,1)-V2:X(I,0)=X(I,0)+V1:X(I,1)=X(I,1)+V2
2710 NEXT I
2720 U3=U1:U4=U2:U1=U3*W1-U4*W2:U2=U4*W1+U3*W2
2730 NEXT J,L
2740 IF D=1 THEN 2760
2750 FOR I=1 TO N:X(I,0)=X(I,0)/N:X(I,1)=X(I,1)/N:NEXT
2760 RETURN
2770 PRINT:PRINT "END":END