diff --git a/Makefile.win b/Makefile.win index 24e29f09c..f5e39bbc4 100644 --- a/Makefile.win +++ b/Makefile.win @@ -11,29 +11,23 @@ all: MAP65.EXE OBJS2C = init_rs.o encode_rs.o decode_rs.o jtaudio.o -F2PYONLY = ftn_init ftn_quit audio_init spec getfile azdist0 astro0 map65a0 +F2PYONLY = ftn_init ftn_quit audio_init getfile azdist0 astro0 \ + spec map65a0 SRCS2F90 = a2d.f90 astro0.f90 audio_init.f90 azdist0.f90 \ - decode1.f90 decode2.f90 decode3.f90 ftn_init.f90 \ - ftn_quit.f90 get_fname.f90 getfile.f90 \ - i1tor4.f90 pix2d65.f90 rfile.f90 savedata.f90 spec.f90 \ - wsjtgen.f90 runqqq.f90 fivehz.f90 \ - flushqqq.f90 sysqqq.f90 map65a0.f90 rfile3a.f90 + decode1.f90 ftn_init.f90 ftn_quit.f90 wsjtgen.f90 \ + runqqq.f90 fivehz.f90 flushqqq.f90 sysqqq.f90 map65a0.f90 \ + rfile.f90 rfile3a.f90 spec.f90 map65a.f90 -SRCS2F77 = map65a.f wsjt1.f bzap.f avesp2.f flatten.f \ - indexx.f flat2.f gen65.f chkmsg.f \ - gentone.f gencwid.f \ - set.f db.f pctile.f sort.f ssort.f ps.f smooth.f \ - peakup.f avemsg65.f decode65.f demod64a.f \ - encode65.f extract.f flat1.f four2.f getpfx1.f \ - getpfx2.f getsnr.f graycode.f grid2k.f interleave63.f k2grid.f \ - limit.f lpf1.f deep65.f morse.f nchar.f packcall.f packgrid.f \ - packmsg.f packtext.f setup65.f short65.f slope.f spec2d65.f \ - sync65.f unpackcall.f unpackgrid.f unpackmsg.f unpacktext.f \ - xcor.f xfft.f xfft2.f wsjt65.f astro.f azdist.f coord.f dcoord.f \ +SRCS2F77 = indexx.f gen65.f chkmsg.f \ + gentone.f gencwid.f set.f db.f pctile.f sort.f ssort.f \ + avemsg65.f demod64a.f encode65.f extract.f four2a.f getpfx1.f \ + getpfx2.f graycode.f grid2k.f interleave63.f k2grid.f \ + deep65.f morse.f nchar.f packcall.f packgrid.f \ + packmsg.f packtext.f setup65.f unpackcall.f unpackgrid.f \ + unpackmsg.f unpacktext.f astro.f azdist.f coord.f dcoord.f \ deg2grid.f dot.f ftsky.f geocentric.f GeoDist.f grid2deg.f \ moon2.f MoonDop.f sun.f toxyz.f pfxdump.f \ - ftpeak65.f fil651.f fil652.f fil653.f symsync65.f \ symspec.f ccf65.f trimlist.f display.f chkhist.f decode1a.f \ filbig.f fil659.f fil658.f fil6521.f twkfreq.f decode65b.f \ afc65b.f fchisq.f ccf2.f diff --git a/a2d.f90 b/a2d.f90 index 5cbef8d5f..18665eb6e 100644 --- a/a2d.f90 +++ b/a2d.f90 @@ -10,21 +10,17 @@ subroutine a2d(iarg) ! JTaudio goes into a test-and-sleep loop. write(*,1000) -1000 format('Using PortAudio.') - idevin=ndevin +1000 format('Using Linrad for input, PortAudio for output.') idevout=ndevout call padevsub(numdevs,ndefin,ndefout,nchin,nchout) - write(*,1002) ndefin,ndefout -1002 format(/'Default Input:',i3,' Output:',i3) - write(*,1004) idevin,idevout -1004 format('Requested Input:',i3,' Output:',i3) - if(idevin.lt.0 .or. idevin.ge.numdevs) idevin=ndefin + write(*,1002) ndefout +1002 format(/'Default Output:',i3) + write(*,1004) idevout +1004 format('Requested Output:',i3) if(idevout.lt.0 .or. idevout.ge.numdevs) idevout=ndefout - if(idevin.eq.0 .and. idevout.eq.0) then - idevin=ndefin - idevout=ndefout - endif + if(idevout.eq.0) idevout=ndefout + idevin=0 ierr=jtaudio(idevin,idevout,y1,y2,NMAX,iwrite,iwave,nwave, & 11025,NSPB,TRPeriod,TxOK,ndebug,Transmitting, & Tsec,ngo,nmode,tbuf,ibuf,ndsec) @@ -32,7 +28,7 @@ subroutine a2d(iarg) print*,'Error ',ierr,' in JTaudio, cannot continue.' else write(*,1006) -1006 format('Audio streams terminated normally.') +1006 format('Audio output stream terminated normally.') endif return end subroutine a2d diff --git a/audio_init.F90 b/audio_init.F90 index d93a4b244..71bb696ea 100644 --- a/audio_init.F90 +++ b/audio_init.F90 @@ -10,16 +10,10 @@ subroutine audio_init(ndin,ndout) include 'gcom1.f90' include 'gcom2.f90' - nmode=1 - if(mode(1:4).eq.'JT65') then - nmode=2 - if(mode(5:5).eq.'A') mode65=1 - if(mode(5:5).eq.'B') mode65=2 - if(mode(5:5).eq.'C') mode65=4 - endif - if(mode.eq.'Echo') nmode=3 - if(mode.eq.'JT6M') nmode=4 - ndevin=ndin + nmode=2 + if(mode(5:5).eq.'A') mode65=1 + if(mode(5:5).eq.'B') mode65=2 + if(mode(5:5).eq.'C') mode65=4 ndevout=ndout TxOK=0 Transmitting=0 diff --git a/decode1.F90 b/decode1.F90 index fac24fdb2..4c30ed2ba 100644 --- a/decode1.F90 +++ b/decode1.F90 @@ -1,4 +1,3 @@ - !---------------------------------------------------- decode1 subroutine decode1(iarg) @@ -24,22 +23,11 @@ subroutine decode1(iarg) ns0=999999 10 continue - if(mode(1:4).eq.'JT65') then - if(rxdone) then - call savedata - rxdone=.false. - endif - else - if(ntr.ne.ntr0 .and. monitoring.gt.0) then - if(ntr.ne.TxFirst .or. (lauto.eq.0)) call savedata - ntr0=ntr - endif - endif - if(ndecoding.gt.0) then ndecdone=0 - call decode2 - ndecdone=1 +! call decode2 + call map65a + ndecdone=2 if(mousebutton.eq.0) ndecoding0=ndecoding ndecoding=0 endif diff --git a/decode3.F90 b/decode3.F90 index e28f464ff..fdb2245cc 100644 --- a/decode3.F90 +++ b/decode3.F90 @@ -55,12 +55,13 @@ subroutine decode3(d2,jz,istart,filename) enddo jz=min(60*11025,jz+nzero) endif - call wsjt1(d2d,jz,istart,samfacin,FileID,ndepth,MinSigdB, & - NQRN,DFTolerance,MouseButton,NClearAve, & - nMode,NFreeze,NAFC,NZap,mode65,idf, & - MyCall,HisCall,HisGrid,neme,nsked,ntx2,s2, & - ps0,npkept,lumsg,basevb,rmspower,nslim2,psavg,ccf,Nseg, & - MouseDF,NAgain,LDecoded,nspecial,ndf,ss1,ss2) +! call wsjt1(d2d,jz,istart,samfacin,FileID,ndepth,MinSigdB, & +! NQRN,DFTolerance,MouseButton,NClearAve, & +! nMode,NFreeze,NAFC,NZap,mode65,idf, & +! MyCall,HisCall,HisGrid,neme,nsked,ntx2,s2, & +! ps0,npkept,lumsg,basevb,rmspower,nslim2,psavg,ccf,Nseg, & +! MouseDF,NAgain,LDecoded,nspecial,ndf,ss1,ss2) + basevb=-999. close(23) if(basevb.le.-98.0) go to 999 diff --git a/filbig.f b/filbig.f index 529817aa9..9c881cc59 100644 --- a/filbig.f +++ b/filbig.f @@ -13,7 +13,7 @@ C f_stop = 750 Hz. real halfpulse(8) !Impulse response of filter (one side) complex cfilt(NFFT2) !Filter (complex; imag = 0) real rfilt(NFFT2) !Filter (real) - integer*8 plan1,plan2,plan3,plan4,plan5 + integer plan1,plan2,plan3,plan4,plan5 logical first include 'fftw3.f' equivalence (rfilt,cfilt) diff --git a/four2a.f b/four2a.f index 10b2d5a3b..0c9208983 100644 --- a/four2a.f +++ b/four2a.f @@ -21,7 +21,7 @@ C The transform will be real and returned to the input array. complex a(nfft) complex aa(32768) integer nn(NPMAX),ns(NPMAX),nf(NPMAX),nl(NPMAX) - integer*8 plan(NPMAX) + real*8 plan(NPMAX) !Should be i*8 data nplan/0/ include 'fftw3.f' save diff --git a/ftn_quit.f90 b/ftn_quit.f90 index 551149eef..ea38d7b74 100644 --- a/ftn_quit.f90 +++ b/ftn_quit.f90 @@ -1,6 +1,8 @@ - !------------------------------------------------ ftn_quit subroutine ftn_quit + include 'gcom1.f90' + ngo=0 call four2a(a,-1,1,1,1) +! @@@ Should also terminate the plans created in filbig! return end subroutine ftn_quit diff --git a/gcom2.f90 b/gcom2.f90 index 7c10b0ce9..0addd90a7 100644 --- a/gcom2.f90 +++ b/gcom2.f90 @@ -31,6 +31,7 @@ integer minsigdb !Decoder threshold setting GUI integer nclearave !Set to 1 to clear JT65 avg GUI,Decoder integer nfreeze !Is Freeze checked? GUI integer nafc !Is AFC checked? GUI +integer newspec !New spectra in ss(4,322,NSMAX) GUI,Decoder integer nmode !Which WSJT mode? GUI,Decoder integer mode65 !JT65 sub-mode (A/B/C ==> 1/2/4) GUI,SoundIn,Decoder integer nclip !Clipping level GUI @@ -89,7 +90,7 @@ common/gcom2/ps0(431),psavg(450),s2(64,3100),ccf(-5:540), & green(500),ngreen,dgain,iter,ndecoding,ndecoding0,mousebutton, & ndecdone,npingtime,ierr,lauto,mantx,nrestart,ntr,nmsg,nsave,nadd5, & dftolerance,LDecoded,rxdone,monitoring,nzap,nsavecum,minsigdb, & - nclearave,nfreeze,nafc,nmode,mode65,nclip,ndebug,nblank,nport, & + nclearave,nfreeze,nafc,newspec,nmode,mode65,nclip,ndebug,nblank,nport, & mousedf,neme,nsked,naggressive,ntx2,nslim2,nagain,nsavelast, & shok,sendingsh,d2a(661500),d2b(661500),b(60000),jza,jzb,ntime, & idinterval,msmax,lenappdir,idf,ndiskdat,nlines,nflat,ntxreq,ntxnow, & diff --git a/jtaudio.c b/jtaudio.c index b6f76475c..fc1b2241d 100644 --- a/jtaudio.c +++ b/jtaudio.c @@ -238,19 +238,18 @@ int padevsub_(int *numdev, int *ndefin, int *ndefout, } - printf("\nAudio Input Output Device Name\n"); - printf("Device Channels Channels\n"); - printf("------------------------------------------------------------------\n"); + printf("\nAudio Output Device Name\n"); + printf("Device Channels\n"); + printf("----------------------------------------------------------\n"); for( i=0; imaxInputChannels; nchout[i]=pdi->maxOutputChannels; - printf(" %2d %2d %2d %s\n",i,nchin[i],nchout[i],pdi->name); + if(nchin[i]>0) + printf(" %2d %2d %s\n",i,nchout[i],pdi->name); } Pa_Terminate(); diff --git a/limit.f b/limit.f deleted file mode 100644 index aaa2927c2..000000000 --- a/limit.f +++ /dev/null @@ -1,31 +0,0 @@ - subroutine limit(x,jz) - - real x(jz) - logical noping - common/limcom/ nslim2 - - noping=.false. - xlim=1.e30 - if(nslim2.eq.1) xlim=3.0 - if(nslim2.ge.2) xlim=1.0 - if(nslim2.ge.3) noping=.true. - - sq=0. - do i=1,jz - sq=sq+x(i)*x(i) - enddo - rms=sqrt(sq/jz) - rms0=14.5 - x1=xlim*rms0 - fac=1.0/xlim - if(fac.lt.1.0) fac=1.0 - if(noping .and. rms.gt.20.0) fac=0.01 !Crude attempt at ping excision - - do i=1,jz - if(x(i).lt.-x1) x(i)=-x1 - if(x(i).gt.x1) x(i)=x1 - x(i)=fac*x(i) - enddo - - return - end diff --git a/map65.py b/map65.py index 54367d054..ed8296ba8 100644 --- a/map65.py +++ b/map65.py @@ -251,11 +251,12 @@ def avetextkey(event=NONE): def decode(event=NONE): if Audio.gcom2.ndecoding==0: #If already busy, ignore request Audio.gcom2.nagain=1 - Audio.gcom2.npingtime=0 #Decode whole record n=1 Audio.gcom2.mousebutton=0 if Audio.gcom2.ndecoding0==4: n=4 Audio.gcom2.ndecoding=n #Standard decode, full file (d2a) +# if Audio.gcom2.ndecoding: +# Audio.map65a0() # @@@ Temporary @@@ #------------------------------------------------------ decode_include def decode_include(event=NONE): @@ -490,7 +491,7 @@ def cleartext(): def ModeJT65(): global slabel,isync,textheight,itol cleartext() - lab2.configure(text='FileID Sync dB DT DF W') + lab2.configure(text='FileID Sync dB DT DF W') lab4.configure(fg='gray85') lab5.configure(fg='gray85') Audio.gcom1.trperiod=60 @@ -1179,9 +1180,6 @@ def update(): bdecode.configure(bg='gray85',activebackground='gray95') if Audio.gcom2.ndecoding: #Set button bg=light_blue while decoding bdecode.configure(bg='#66FFFF',activebackground='#66FFFF') -# print 'A' - Audio.map65a0() # @@@ Temporary @@@ -# print 'B' tx1.configure(bg='white') tx2.configure(bg='white') @@ -1225,7 +1223,7 @@ def update(): t='Receiving' msg7.configure(text=t,bg=bgcolor) - if Audio.gcom2.ndecdone==1 or g.cmap != cmap0: + if Audio.gcom2.ndecdone>0 or g.cmap != cmap0: if Audio.gcom2.ndecdone==1: if isync==-99 or isync==99: isync=isync_save @@ -1255,7 +1253,10 @@ def update(): avetext.insert(END,lines[0]) avetext.insert(END,lines[1]) # avetext.configure(state=DISABLED) + Audio.gcom2.ndecdone=0 + + if Audio.gcom2.ndecdone==2: try: f=open(appdir+'/bandmap.txt',mode='r') lines=f.readlines() @@ -1265,12 +1266,11 @@ def update(): bmtext.configure(state=NORMAL) bmtext.insert(END,' Freq DF Pol UTC\n') bmtext.insert(END,'--------------------------------------------\n') - for i in range(len(lines)): bmtext.insert(END,lines[i]) bmtext.see(END) - Audio.gcom2.ndecdone=2 + Audio.gcom2.ndecdone=3 if g.cmap != cmap0: im.putpalette(g.palette) diff --git a/map65a.f b/map65a.f deleted file mode 100644 index a8447b783..000000000 --- a/map65a.f +++ /dev/null @@ -1,288 +0,0 @@ - subroutine map65a - -C Processes timf2 data from Linrad to find and decode JT65 signals. - - parameter (NMAX=60*96000) !Samples per 60 s file - parameter (MAXMSG=1000) !Size of decoded message list - integer*2 id(4,NMAX) !46 MB: raw data from Linrad timf2 - parameter (NFFT=32768) !Half symbol = 17833 samples; - real ss(4,322,NFFT) !169 MB: half-symbol spectra - real savg(4,NFFT) - real tavg(-50:50) !Temp for finding local base level - real base(4) !Local basel level at 4 pol'ns - real tmp (200) !Temp storage for pctile sorting - real short(3,NFFT) !SNR dt ipol for potential shorthands - real sig(MAXMSG,30) !Parameters of detected signals - real a(5) - character*22 msg(MAXMSG) - character*3 shmsg0(4),shmsg - character arg*12,infile*11,outfile*11 - integer indx(MAXMSG),nsiz(MAXMSG) - logical done(MAXMSG) - integer rfile3 - character decoded*22,blank*22,cbad*1 - data blank/' '/ - data shmsg0/'ATT','RO ','RRR','73 '/ - - tskip=0. - fselect=0. -! fselect=103.0 - nmin=1 - infile='061111.0745' - -C Initialize some constants - -! open(22,file='kvasd.dat',access='direct',recl=1024, -! + status='unknown') - open(23,file='CALL3.TXT',status='unknown') - -! nbytes=8*(4*96000+9000) !Empirical, for 061111_0744.dat.48 -! nskip=8*nint(96000*(tskip+4.09375)) -! n=rfile3(infile,id,nskip) !Skip to start of minute -! if(n.ne.nskip) go to 9999 - - df=96000.0/NFFT !df = 96000/NFFT = 2.930 Hz - fa=0.0 - fb=60000.0 - ia=nint((fa+23000.0)/df + 1.0) ! 23000 = 48000 - 25000 - ib=nint((fb+23000.0)/df + 1.0) - ftol=0.020 !Frequency tolerance (kHz) - kk=0 - nkk=1 - - do nfile=1,nmin -! n=rfile3(infile,id,8*NMAX) !Read 60 s of data (approx 46 MB) - n=8*NMAX - call rfile3a(infile,id,n,ierr) - newdat=1 - nz=n/8 - read(infile(8:11),*) nutc - if(fselect.gt.0.0) then - nfilt=1 !nfilt=2 is faster for selected freq - freq=fselect+1.600 - nflip=-1 !May need to try both +/- 1 - ipol=4 !Try all four? - dt=2.314240 !Not needed? - call decode1a(id,newdat,nfilt,freq,nflip,ipol,sync2, - + a,dt,pol,nkv,nhist,qual,decoded) - write(11,1010) nutc,nsync1,nsync2,dt,ndf,decoded, - + nkv,nqual - 1010 format(i4.4,i5,i4,f5.1,i5,2x,a22,2i3) - if(nfile.eq.1) go to 999 - endif - - nfilt=1 - do i=1,NFFT - short(1,i)=0. - short(2,i)=0. - short(3,i)=0. - enddo - - call symspec(id,nz,ss,savg) - - freq0=-999. - sync10=-999. - fshort0=-999. - sync20=-999. - ntry=0 - do i=ia,ib !Search over freq range - freq=0.001*((i-1)*df - 23000) + 100.0 - -C Find the local base level for each polarization; update every 10 bins. - if(mod(i-ia,10).eq.0) then - do jp=1,4 - do ii=-50,50 - tavg(ii)=savg(jp,i+ii) - enddo - call pctile(tavg,tmp,101,50,base(jp)) - enddo - endif - -C Find max signal at this frequency - smax=0. - do jp=1,4 - if(savg(jp,i)/base(jp).gt.smax) smax=savg(jp,i)/base(jp) - enddo - - if(smax.gt.1.1) then - ntry=ntry+1 -C Look for JT65 sync patterns and shorthand square-wave patterns. - call ccf65(ss(1,1,i),sync1,ipol,dt,flipk, - + syncshort,snr2,ipol2,dt2) - - shmsg=' ' -C Is there a shorthand tone above threshold? - if(syncshort.gt.1.0) then - -C### Do shorthand AFC here (or maybe after finding a pair?) ### - - short(1,i)=syncshort - short(2,i)=dt2 - short(3,i)=ipol2 -C Check to see if lower tone of shorthand pair was found. - do j=2,4 - i0=i-nint(j*53.8330078/df) -C Should this be i0 +/- 1, or just i0? -C Should we also insist that difference in DT be either 1.5 or -1.5 s? - if(short(1,i0).gt.1.0) then - fshort=0.001*((i0-1)*df - 23000) + 100.0 - -C Keep only the best candidate within ftol. - if(fshort-fshort0.le.ftol .and. sync2.gt.sync20 - + .and. nkk.eq.2) kk=kk-1 - if(fshort-fshort0.gt.ftol .or. - + sync2.gt.sync20) then - kk=kk+1 - sig(kk,1)=nfile - sig(kk,2)=nutc - sig(kk,3)=fshort - sig(kk,4)=syncshort - sig(kk,5)=dt2 - sig(kk,6)=45*(ipol2-1)/57.2957795 - sig(kk,7)=0 - sig(kk,8)=snr2 - sig(kk,9)=0 - sig(kk,10)=0 -! sig(kk,11)=rms0 - sig(kk,12)=savg(ipol2,i) - sig(kk,13)=0 - sig(kk,14)=0 - sig(kk,15)=0 - sig(kk,16)=0 -! sig(kk,17)=0 - sig(kk,18)=0 - msg(kk)=shmsg0(j) - fshort0=fshort - sync20=sync2 - nkk=2 - endif - endif - enddo - endif - -C Is sync1 above threshold? - if(sync1.gt.1.0) then - -C Keep only the best candidate within ftol. -C (Am I deleting any good decodes by doing this? Any harm in omitting -C these statements??) - if(freq-freq0.le.ftol .and. sync1.gt.sync10 .and. - + nkk.eq.1) kk=kk-1 - - if(freq-freq0.gt.ftol .or. sync1.gt.sync10) then - nflip=nint(flipk) - - call decode1a(id,newdat,nfilt,freq,nflip,ipol, - + sync2,a,dt,pol,nkv,nhist,qual,decoded) - - kk=kk+1 - sig(kk,1)=nfile - sig(kk,2)=nutc - sig(kk,3)=freq - sig(kk,4)=sync1 - sig(kk,5)=dt - sig(kk,6)=pol - sig(kk,7)=flipk - sig(kk,8)=sync2 - sig(kk,9)=nkv - sig(kk,10)=qual -! sig(kk,11)=rms0 - sig(kk,12)=savg(ipol,i) - sig(kk,13)=a(1) - sig(kk,14)=a(2) - sig(kk,15)=a(3) - sig(kk,16)=a(4) -! sig(kk,17)=a(5) - sig(kk,18)=nhist - msg(kk)=decoded - freq0=freq - sync10=sync1 - nkk=1 - endif - endif - endif - enddo - -! write(*,1010) - -C Trim the list and produce a sorted index and sizes of groups. -C (Should trimlist remove all but best SNR for given UTC and message content?) - call trimlist(sig,kk,indx,nsiz,nz) - - do i=1,kk - done(i)=.false. - enddo - j=0 - ilatest=-1 - do n=1,nz - ifile0=0 - do m=1,nsiz(n) - i=indx(j+m) - ifile=sig(i,1) - if(ifile.gt.ifile0 .and.msg(i).ne.blank) then - ilatest=i - ifile0=ifile - endif - enddo - i=ilatest - if(i.ge.1) then - if(.not.done(i)) then - done(i)=.true. - nutc=sig(i,2) - freq=sig(i,3) - sync1=sig(i,4) - dt=sig(i,5) - npol=nint(57.2957795*sig(i,6)) - flip=sig(i,7) - sync2=sig(i,8) - nkv=sig(i,9) - nqual=min(sig(i,10),10.0) -! rms0=sig(i,11) - nsavg=sig(i,12) !Was used for diagnostic ... - do k=1,5 - a(k)=sig(i,12+k) - enddo - nhist=sig(i,18) - decoded=msg(i) - - if(flip.lt.0.0) then - do i=22,1,-1 - if(decoded(i:i).ne.' ') go to 10 - enddo - stop 'Error in message format' - 10 if(i.le.18) decoded(i+2:i+4)='OOO' - endif - nkHz=nint(freq-1.600) - f0=144.0+0.001*nkHz - ndf=nint(1000.0*(freq-1.600-nkHz)) - ndf0=nint(a(1)) - ndf1=nint(a(2)) - ndf2=nint(a(3)) - nsync1=sync1 - nsync2=nint(10.0*log10(sync2)) - 40 !### empirical ### - cbad=' ' - - if(abs(f0-144.103).lt.0.001) then - write(11,1010) nutc,nsync1,nsync2,dt,ndf,decoded, - + nkv,nqual - endif - - write(19,1012) f0,ndf,npol,nutc,decoded - 1012 format(f7.3,i5,i4,i5.4,2x,a22) - - write(26,1014) f0,ndf,ndf0,ndf1,ndf2,dt,npol,nsync1, - + nsync2,nutc,decoded,nkv,nqual,nhist - 1014 format(f7.3,i5,3i3,f5.1,i5,i3,i4,i5.4,2x,a22,3i3) - - endif - endif - j=j+nsiz(n) - enddo - call display(nutc) -! if(nfile.ge.1) go to 999 - 100 continue - enddo - - 999 call four2a(cx,-1,1,1,1) !Destroy the FFTW plans - - 9999 end diff --git a/map65a.f90 b/map65a.f90 new file mode 100644 index 000000000..22001ee78 --- /dev/null +++ b/map65a.f90 @@ -0,0 +1,296 @@ +subroutine map65a + +! Processes timf2 data from Linrad to find and decode JT65 signals. + + parameter (NSMAX=60*96000) !Samples per 60 s file + parameter (MAXMSG=1000) !Size of decoded message list + integer*2 id(4,NSMAX) !46 MB: raw data from Linrad timf2 + parameter (NFFT=32768) !Half symbol = 17833 samples; + real savg(4,NFFT) + real tavg(-50:50) !Temp for finding local base level + real base(4) !Local basel level at 4 pol'ns + real tmp (200) !Temp storage for pctile sorting + real short(3,NFFT) !SNR dt ipol for potential shorthands + real sig(MAXMSG,30) !Parameters of detected signals + real a(5) + character*22 msg(MAXMSG) + character*3 shmsg0(4),shmsg + character arg*12,infile*11,outfile*11 + integer indx(MAXMSG),nsiz(MAXMSG) + logical done(MAXMSG) + integer rfile3 + character decoded*22,blank*22,cbad*1 + common/spcom/ss(4,322,NFFT) !169 MB: half-symbol spectra + data blank/' '/ + data shmsg0/'ATT','RO ','RRR','73 '/ + + include 'gcom2.f90' + + tskip=0. + fselect=0. + fselect=104.5303 + nmin=1 + infile='061111.0745' + +! Initialize some constants + +! open(22,file='kvasd.dat',access='direct',recl=1024, +! + status='unknown') + open(23,file='CALL3.TXT',status='unknown') + +! nbytes=8*(4*96000+9000) !Empirical, for 061111_0744.dat.48 +! nskip=8*nint(96000*(tskip+4.09375)) +! n=rfile3(infile,id,nskip) !Skip to start of minute +! if(n.ne.nskip) go to 9999 + + df=96000.0/NFFT !df = 96000/NFFT = 2.930 Hz + fa=0.0 + fb=60000.0 + ia=nint((fa+23000.0)/df + 1.0) ! 23000 = 48000 - 25000 + ib=nint((fb+23000.0)/df + 1.0) + ftol=0.020 !Frequency tolerance (kHz) + kk=0 + nkk=1 + + do nfile=1,nmin +! n=rfile3(infile,id,8*NSMAX) !Read 60 s of data (approx 46 MB) + n=8*NSMAX + call rfile3a(infile,id,n,ierr) + newdat=1 + nz=n/8 + read(infile(8:11),*) nutc + if(fselect.gt.0.0) then + +! nfilt=2 should be faster (but doesn't work right?) + nfilt=1 !nfilt=2 is faster for selected freq +! freq=fselect+1.600 + freq=fselect + nflip=-1 !May need to try both +/- 1 + ipol=4 !Try all four? + dt=2.314240 !Not needed? + + call decode1a(id,newdat,nfilt,freq,nflip,ipol,sync2, & + a,dt,pol,nkv,nhist,qual,decoded) + nsync2=nint(10.0*log10(sync2)) - 40 !### empirical ### + nw=0 +! Insert 'OOO' if flip<0. + write(11,1010) nutc,nsync1,nsync2,dt,ndf,nw,decoded, & + nkv,nqual +1010 format(i4.4,i3,i4,f5.1,i5,i3,2x,a22,2i3) + ndecdone=1 + endif + + nfilt=1 + do i=1,NFFT + short(1,i)=0. + short(2,i)=0. + short(3,i)=0. + enddo + + call symspec(id,nz,ss,savg) + newspec=1 + + freq0=-999. + sync10=-999. + fshort0=-999. + sync20=-999. + ntry=0 + do i=ia,ib !Search over freq range + freq=0.001*((i-1)*df - 23000) + 100.0 + +! Find the local base level for each polarization; update every 10 bins. + if(mod(i-ia,10).eq.0) then + do jp=1,4 + do ii=-50,50 + tavg(ii)=savg(jp,i+ii) + enddo + call pctile(tavg,tmp,101,50,base(jp)) + enddo + endif + +! Find max signal at this frequency + smax=0. + do jp=1,4 + if(savg(jp,i)/base(jp).gt.smax) smax=savg(jp,i)/base(jp) + enddo + + if(smax.gt.1.1) then + ntry=ntry+1 +! Look for JT65 sync patterns and shorthand square-wave patterns. + call ccf65(ss(1,1,i),sync1,ipol,dt,flipk, & + syncshort,snr2,ipol2,dt2) + + shmsg=' ' +! Is there a shorthand tone above threshold? + if(syncshort.gt.1.0) then + +! ### Do shorthand AFC here (or maybe after finding a pair?) ### + + short(1,i)=syncshort + short(2,i)=dt2 + short(3,i)=ipol2 +! Check to see if lower tone of shorthand pair was found. + do j=2,4 + i0=i-nint(j*53.8330078/df) +! Should this be i0 +/- 1, or just i0? +! Should we also insist that difference in DT be either 1.5 or -1.5 s? + if(short(1,i0).gt.1.0) then + fshort=0.001*((i0-1)*df - 23000) + 100.0 + +! Keep only the best candidate within ftol. + if(fshort-fshort0.le.ftol .and. sync2.gt.sync20 & + .and. nkk.eq.2) kk=kk-1 + if(fshort-fshort0.gt.ftol .or. & + sync2.gt.sync20) then + kk=kk+1 + sig(kk,1)=nfile + sig(kk,2)=nutc + sig(kk,3)=fshort + sig(kk,4)=syncshort + sig(kk,5)=dt2 + sig(kk,6)=45*(ipol2-1)/57.2957795 + sig(kk,7)=0 + sig(kk,8)=snr2 + sig(kk,9)=0 + sig(kk,10)=0 +! sig(kk,11)=rms0 + sig(kk,12)=savg(ipol2,i) + sig(kk,13)=0 + sig(kk,14)=0 + sig(kk,15)=0 + sig(kk,16)=0 +! sig(kk,17)=0 + sig(kk,18)=0 + msg(kk)=shmsg0(j) + fshort0=fshort + sync20=sync2 + nkk=2 + endif + endif + enddo + endif + +! Is sync1 above threshold? + if(sync1.gt.1.0) then + +! Keep only the best candidate within ftol. +! (Am I deleting any good decodes by doing this? Any harm in omitting +! these statements??) + if(freq-freq0.le.ftol .and. sync1.gt.sync10 .and. & + nkk.eq.1) kk=kk-1 + + if(freq-freq0.gt.ftol .or. sync1.gt.sync10) then + nflip=nint(flipk) + + call decode1a(id,newdat,nfilt,freq,nflip,ipol, & + sync2,a,dt,pol,nkv,nhist,qual,decoded) + kk=kk+1 + sig(kk,1)=nfile + sig(kk,2)=nutc + sig(kk,3)=freq + sig(kk,4)=sync1 + sig(kk,5)=dt + sig(kk,6)=pol + sig(kk,7)=flipk + sig(kk,8)=sync2 + sig(kk,9)=nkv + sig(kk,10)=qual +! sig(kk,11)=rms0 + sig(kk,12)=savg(ipol,i) + sig(kk,13)=a(1) + sig(kk,14)=a(2) + sig(kk,15)=a(3) + sig(kk,16)=a(4) +! sig(kk,17)=a(5) + sig(kk,18)=nhist + msg(kk)=decoded + freq0=freq + sync10=sync1 + nkk=1 + endif + endif + endif + enddo + +! write(*,1010) + +! Trim the list and produce a sorted index and sizes of groups. +! (Should trimlist remove all but best SNR for given UTC and message content?) + call trimlist(sig,kk,indx,nsiz,nz) + + do i=1,kk + done(i)=.false. + enddo + j=0 + ilatest=-1 + do n=1,nz + ifile0=0 + do m=1,nsiz(n) + i=indx(j+m) + ifile=sig(i,1) + if(ifile.gt.ifile0 .and.msg(i).ne.blank) then + ilatest=i + ifile0=ifile + endif + enddo + i=ilatest + if(i.ge.1) then + if(.not.done(i)) then + done(i)=.true. + nutc=sig(i,2) + freq=sig(i,3) + sync1=sig(i,4) + dt=sig(i,5) + npol=nint(57.2957795*sig(i,6)) + flip=sig(i,7) + sync2=sig(i,8) + nkv=sig(i,9) + nqual=min(sig(i,10),10.0) +! rms0=sig(i,11) + nsavg=sig(i,12) !Was used for diagnostic ... + do k=1,5 + a(k)=sig(i,12+k) + enddo + nhist=sig(i,18) + decoded=msg(i) + + if(flip.lt.0.0) then + do i=22,1,-1 + if(decoded(i:i).ne.' ') go to 10 + enddo + stop 'Error in message format' +10 if(i.le.18) decoded(i+2:i+4)='OOO' + endif + nkHz=nint(freq-1.600) + f0=144.0+0.001*nkHz + ndf=nint(1000.0*(freq-1.600-nkHz)) + ndf0=nint(a(1)) + ndf1=nint(a(2)) + ndf2=nint(a(3)) + nsync1=sync1 + nsync2=nint(10.0*log10(sync2)) - 40 !### empirical ### + cbad=' ' + +! if(abs(f0-144.103).lt.0.001) then +! write(11,1010) nutc,nsync1,nsync2,dt,ndf,decoded, +! + nkv,nqual +! endif + + write(19,1012) f0,ndf,npol,nutc,decoded +1012 format(f7.3,i5,i4,i5.4,2x,a22) + + write(26,1014) f0,ndf,ndf0,ndf1,ndf2,dt,npol,nsync1, & + nsync2,nutc,decoded,nkv,nqual,nhist +1014 format(f7.3,i5,3i3,f5.1,i5,i3,i4,i5.4,2x,a22,3i3) + + endif + endif + j=j+nsiz(n) + enddo + call display(nutc) +! if(nfile.ge.1) go to 999 +100 continue + enddo + +999 return +end subroutine map65a diff --git a/map65a0.f90 b/map65a0.f90 index 05b275d18..8eee53ab8 100644 --- a/map65a0.f90 +++ b/map65a0.f90 @@ -1,9 +1,5 @@ subroutine map65a0 - include 'gcom2.f90' - call map65a - ndecdone=1 - return end subroutine map65a0 diff --git a/padevsub.c b/padevsub.c deleted file mode 100644 index dded911f6..000000000 --- a/padevsub.c +++ /dev/null @@ -1,36 +0,0 @@ -// #include -// #include -#include "portaudio.h" - -int __stdcall PADEVSUB(int *numdev, int *ndefin, int *ndefout, - int nchin[], int nchout[]) -{ - int i,j; - int numDevices; - const PaDeviceInfo *pdi; - PaError err; - - Pa_Initialize(); - numDevices = Pa_CountDevices(); - *numdev=numDevices; - if( numDevices < 0 ) { - err = numDevices; - goto error; - } - - for( i=0; imaxInputChannels; - nchout[i]=pdi->maxOutputChannels; - printf("Audio device %d: In=%d Out=%d %s\n",i,nchin[i],nchout[i],pdi->name); - } - - Pa_Terminate(); - return 0; - - error: - Pa_Terminate(); - return err; -} diff --git a/spec.f90 b/spec.f90 index d56cd5312..195cf22db 100644 --- a/spec.f90 +++ b/spec.f90 @@ -1,180 +1,55 @@ -!---------------------------------------------------- spec subroutine spec(brightness,contrast,logmap,ngain,nspeed,a) -! Called by SpecJT. -! Probably should use the "!f2py intent(...)" structure here. + parameter (NX=750,NY=130,NTOT=NX*NY,NFFT=32768) ! Input: - parameter (NX=750,NY=130,NTOT=NX*NY) integer brightness,contrast !Display parameters integer ngain !Digital gain for input audio integer nspeed !Scrolling speed index + ! Output: integer*2 a(NTOT) !Pixel values for NX x NY array - real a0(NTOT) !Save the last NY spectra +! real a0(NTOT) !Save the last NY spectra integer nstep(5) integer b0,c0 - real x(4096) !Data for FFT - complex c(0:2048) !Complex spectrum - real ss(2048) !Power spectrum - logical first + real s(NFFT) + common/spcom/ss(4,322,NFFT) !169 MB: half-symbol spectra include 'gcom1.f90' include 'gcom2.f90' include 'gcom3.f90' include 'gcom4.f90' data jz/0/ !Number of spectral lines available data nstep/15,10,5,2,1/ !Integration limits - data first/.true./ - - equivalence (x,c) save - if(first) then - istep=2205 - nfft=4096 - nh=nfft/2 - do i=1,nh - ss(i)=0. - enddo - df=11025.0/nfft - fac=2.0/10000. - nsum=0 - iread=0 - first=.false. - b0=-999 - c0=-999 - logmap0=-999 - nspeed0=-999 - ncall=0 - jza=0 - rms=0. - endif + df=96000.0/nfft + b0=-999 + c0=-999 + logmap0=-999 + nspeed0=-999 + nmode=2 !JT65 mode - nmode=1 - if(mode(1:4).eq.'JT65') nmode=2 - if(mode.eq.'Echo') nmode=3 - if(mode.eq.'JT6M') nmode=4 - - nlines=0 - newdat=0 - npts=iwrite-iread - if(ndiskdat.eq.1) then - npts=jzc/2048 - npts=2048*npts - kread=0 - endif - if(npts.lt.0) npts=npts+nmax - if(npts.lt.nfft) go to 900 !Not enough data available - -10 continue - if(ndiskdat.eq.1) then -! Data read from disk - k=kread + nadd=nstep(nspeed) + nlines=322/nadd + j=0 + print*,'A',nspeed,nadd,nlines + do k=1,nlines do i=1,nfft - k=k+1 - x(i)=0.4*d2c(k) - enddo - kread=kread+istep !Update pointer - else -! Real-time data - dgain=2.0*10.0**(0.005*ngain) - k=iread - do i=1,nfft - k=k+1 - if(k.gt.nmax) k=k-nmax - x(i)=0.5*dgain*y1(k) - enddo - iread=iread+istep !Update pointer - if(iread.gt.nmax) iread=iread-nmax - endif - - sum=0. !Get ave, rms of data - do i=1,nfft - sum=sum+x(i) - enddo - ave=sum/nfft - sq=0. - do i=1,nfft - d=x(i)-ave - sq=sq+d*d - x(i)=fac*d - enddo - rms1=sqrt(sq/nfft) - if(rms.eq.0) rms=rms1 - rms=0.25*rms1 + 0.75*rms - - if(ndiskdat.eq.0) then - level=0 !Compute S-meter level - if(rms.gt.0.0) then !Scale 0-100, steps = 0.4 dB - dB=20.0*log10(rms/800.0) - level=50 + 2.5*dB - if(level.lt.0) level=0 - if(level.gt.100) level=100 - endif - endif - - call xfft2(x,nfft) - - do i=1,nh !Accumulate power spectrum - ss(i)=ss(i) + real(c(i))**2 + aimag(c(i))**2 - enddo - nsum=nsum+1 - - if(nsum.ge.nstep(nspeed)) then !Integrate for specified time - nlines=nlines+1 - do i=NTOT,NX+1,-1 !Move spectra up one row - a0(i)=a0(i-NX) ! (will be "down" on display) - enddo - if(ndiskdat.eq.1 .and. nlines.eq.1) then - do i=1,NX - a0(i)=255 + s(i)=0. + do n=1,nadd + j=j+1 + s(i)=s(i) + ss(1,j,i) + ss(3,j,i) enddo - do i=NTOT,NX+1,-1 - a0(i)=a0(i-NX) - enddo - endif - - if(nflat.gt.0) call flat2(ss,nh,nsum) - - ia=1 - if(nfrange.eq.2000) then - i0=182 + nint((nfmid-1500)/df) - if(i0.lt.0) ia=1-i0 - else if(nfrange.eq.4000) then - i0=nint(nfmid/df - (NX+2.0)) - if(i0.lt.0) ia=1-i0/2 - endif - do i=ia,NX !Insert new data in top row - if(nfrange.eq.2000) then - a0(i)=5*ss(i+i0)/nsum - else if(nfrange.eq.4000) then - smax=max(ss(2*i+i0),ss(2*i+i0-1)) - a0(i)=5*smax/nsum - endif + s(i)=s(i)/nadd enddo - nsum=0 - newdat=1 !Flag for new spectrum available - do i=1,nh !Zero the accumulating array - ss(i)=0. - enddo - if(jz.lt.NY) jz=jz+1 - endif - - if(ndiskdat.eq.1) then - npts=jzc-kread - else - npts=iwrite-iread - if(npts.lt.0) npts=npts+nmax - endif - - if(npts.ge.4096) go to 10 + print*,'B',k ! Compute pixel values - iz=NX - logmap=0 - if(brightness.ne.b0 .or. contrast.ne.c0 .or. logmap.ne.logmap0 .or. & - nspeed.ne.nspeed0 .or. nlines.gt.1) then +! iz=NX +! logmap=0 +! if(brightness.ne.b0 .or. contrast.ne.c0 .or. logmap.ne.logmap0 .or. & +! nspeed.ne.nspeed0 .or. nlines.gt.1) then iz=NTOT gain=40*sqrt(nstep(nspeed)/5.0) * 5.0**(0.01*contrast) gamma=1.3 + 0.01*contrast @@ -183,17 +58,20 @@ subroutine spec(brightness,contrast,logmap,ngain,nspeed,a) c0=contrast logmap0=logmap nspeed0=nspeed - endif + ia=10001 + ib=10750 - do i=1,iz - n=0 - if(a0(i).gt.0.0 .and. logmap.eq.1) n=gain*log10(0.001*a0(i)) + offset + 20 - if(a0(i).gt.0.0 .and. logmap.eq.0) n=(0.01*a0(i))**gamma + offset - n=min(252,max(0,n)) - a(i)=n + do i=ia,ib + n=0 + if(s(i).gt.0.0 .and. logmap.eq.1) n=gain*log10(0.001*s(i)) & + + offset + 20 + if(s(i).gt.0.0 .and. logmap.eq.0) n=(0.01*s(i))**gamma + offset + n=min(252,max(0,n)) + a(i)=n + enddo + print*,'C' enddo + print*,'D' -900 continue - if(ndiskdat.eq.1) ndecoding=4 return end subroutine spec diff --git a/specjt.py b/specjt.py index 7bc3c0814..f5147d7fd 100644 --- a/specjt.py +++ b/specjt.py @@ -247,16 +247,12 @@ def update(): contrast=sc2.get() logm=logmap.get() g0=sc3.get() - -# Don't calculate spectra for waterfall while decoding - if Audio.gcom2.ndecoding==0 and \ - (Audio.gcom2.monitoring or Audio.gcom2.ndiskdat): - Audio.spec(brightness,contrast,logm,g0,nspeed,a) #Call Fortran routine spec - newdat=Audio.gcom1.newdat #True if new data available - else: - newdat=0 - if newdat or brightness!=b0 or contrast!=c0 or logm!=logm0: + newspec=Audio.gcom2.newspec #True if new data available + if newspec: + Audio.spec(brightness,contrast,logm,g0,nspeed,a) #Call Fortran routine spec + + if newspec or brightness!=b0 or contrast!=c0 or logm!=logm0: if brightness==b0 and contrast==c0 and logm==logm0: n=Audio.gcom2.nlines box=(0,0,NX,300-n) #Define region @@ -275,7 +271,7 @@ def update(): c0=contrast logm0=logm - if newdat: + if newspec: if Audio.gcom2.monitoring: if minsep.get() and newMinute: draw.line((0,0,749,0),fill=128) #Draw the minute separator @@ -290,6 +286,7 @@ def update(): #For some reason, top two lines are invisible, so we move down 2 graph1.create_image(0,0+2,anchor='nw',image=pim) newMinute=0 + Audio.gcom2.newspec=0 if (Audio.gcom2.mousedf != mousedf0 or Audio.gcom2.dftolerance != tol0): df_mark() @@ -312,7 +309,7 @@ def update(): df_mark() nmark0=nmark.get() - if newdat: Audio.gcom2.ndiskdat=0 + if newspec: Audio.gcom2.ndiskdat=0 Audio.gcom2.nlines=0 Audio.gcom2.nflat=nflat.get() frange=nfr.get()*2000 diff --git a/wsjt1.F b/wsjt1.F deleted file mode 100644 index c320c0ba6..000000000 --- a/wsjt1.F +++ /dev/null @@ -1,194 +0,0 @@ - subroutine wsjt1(d,jz0,istart,samfacin,FileID,ndepth,MinSigdB, - + NQRN,DFTolerance,MouseButton,NClearAve, - + Mode,NFreeze,NAFC,NZap,mode65,idf, - + MyCall,HisCall,HisGrid,neme,nsked,ntx2,s2, - + ps0,npkept,lumsg,basevb,rmspower,nslim2,psavg,ccf,Nseg, - + MouseDF,NAgain,LDecoded,nspecial,ndf,ss1,ss2) - - parameter (NP2=1024*1024) - - integer*2 d(jz0) !Buffer for raw one-byte data - integer istart !Starting location in original d() array - character FileID*40 !Name of file being processed - integer MinSigdB !Minimum ping strength, dB - integer NQRN !QRN rejection parameter - integer DFTolerance !Defines DF search range - integer NSyncOK !Set to 1 if JT65 file synchronized OK - character*12 mycall - character*12 hiscall - character*6 hisgrid - real ps0(431) !Spectrum of best ping - integer npkept !Number of pings kept and decoded - integer lumsg !Logical unit for decoded.txt - real basevb !Baseline signal level, dB - integer nslim2 !Minimum strength for single-tone pings, dB - real psavg(450) !Average spectrum of the whole file - integer Nseg !First or second Tx sequence? - integer MouseDF !Freeze position for DF - logical pick !True if this is a mouse-picked ping - logical stbest !True if the best decode was Single-Tone - logical STfound !True if at least one ST decode - logical LDecoded !True if anything was decoded - real s2(64,3100) !2D spectral array - real ccf(-5:540) !X-cor function in JT65 mode (blue line) - real red(512) - real ss1(-224:224) !Magenta curve (for JT65 shorthands) - real ss2(-224:224) !Orange curve (for JT65 shorthands) - real yellow(216) - real yellow0(216) - real fzap(200) - integer resample - real*8 samfacin,samratio - real dat2(NP2) - character msg3*3 - character cfile6*6 - logical lcum - integer indx(100) - character*90 line - common/avecom/dat(NP2),labdat,jza,modea - common/ccom/nline,tping(100),line(100) - common/limcom/ nslim2a - common/clipcom/ nclip - save - - lcum=.true. - jz=jz0 - modea=Mode - nclip=NQRN-5 - nslim2a=nclip - MinWidth=40 !Minimum width of pings, ms - call zero(psavg,450) - rewind 11 - rewind 12 - - do i=1,40 - if(FileID(i:i).eq.'.') go to 3 - enddo - i=4 - 3 ia=max(1,i-6) - cfile6=FileID(ia:i-1) - - nline=0 - ndiag=0 -! If file "/wsjt.reg" exists, set ndiag=1 - open(16,file='/wsjt.reg',status='old',err=4) - ndiag=1 - close(16) - - 4 if(jz.gt.655360) jz=655360 - - sum=0. - do j=1,jz !Convert raw data from i*2 to real, remove DC - dat(j)=0.1*d(j) - sum=sum + dat(j) - enddo - ave=sum/jz - samratio=1.d0/samfacin - if(samratio.eq.1.d0) then - do j=1,jz - dat(j)=dat(j)-ave - enddo - else - do j=1,jz - dat2(j)=dat(j)-ave - enddo - -#if (USE_PORTAUDIO==1) || defined(Win32) - ierr=resample(dat2,dat,samratio,jz) - if(ierr.ne.0) print*,'Resample error.',samratio -#endif - - endif - - if(ndiag.ne.0 .and. nclip.lt.0) then -C Intentionally degrade SNR by -nclip dB. - sq=0. - do i=1,jz - sq=sq + dat(i)**2 - enddo - p0=sq/jz - p1=p0*10.0**(-0.1*nclip) - dnoise=sqrt(4*(p1-p0)) - idum=-1 - do i=1,jz - dat(i)=dat(i) + dnoise*gran(idum) - enddo - endif - - sq=0. - do j=1,jz !Compute power level for whole array - sq=sq + dat(j)**2 - enddo - avesq=sq/jz - basevb=dB(avesq) - 44 !Base power level to send back to GUI - if(avesq.eq.0) go to 900 - - nz=600 - nstep=jz/nz - sq=0. - k=0 - do j=1,nz - sum=0. - do n=1,nstep - k=k+1 - sum=sum+dat(k)**2 - enddo - sum=sum/nstep - sq=sq + (sum-avesq)**2 - enddo - rmspower=sqrt(sq/nz) - - pick=.false. - if(istart.gt.1) pick=.true. !This is a mouse-picked decoding - if(.not.pick .and. (basevb.lt.-15.0 .or. basevb.gt.20.0)) goto 900 - nchan=64 !Save 64 spectral channels - nstep=221 !Set step size to ~20 ms - nz=jz/nstep - 1 !# of spectra to compute - if(.not.pick) then - MouseButton=0 - jza=jz - labdat=labdat+1 - endif - tbest=0. - NsyncOK=0 - -! Only JT65 mode is supported. -! Check for a JT65 shorthand message - nstest=0 - if(ntx2.ne.1) call short65(dat,jz,NFreeze,MouseDF, - + DFTolerance,mode65,nspecial,nstest,dfsh,iderrsh, - + idriftsh,snrsh,ss1,ss2,nwsh,idfsh) -! Lowpass filter and decimate by 2 - call lpf1(dat,jz,jz2,MouseDF,MouseDF2) - idf=mousedf-mousedf2 - jz=jz2 - nadd=1 - fzap(1)=0. - if(nzap.eq.1) call avesp2(dat,jz,nadd,mode,NFreeze,MouseDF2, - + DFTolerance,fzap) - if(nzap.eq.1.and.nstest.eq.0) call bzap(dat,jz,nadd,mode,fzap) - - i=index(MyCall,char(0)) - if(i.le.0) i=index(MyCall,' ') - mycall=MyCall(1:i-1)//' ' - i=index(HisCall,char(0)) - if(i.le.0) i=index(HisCall,' ') - hiscall=HisCall(1:i-1)//' ' - -! Offset data by about 1 s. - if(jz.ge.126*2048) call wsjt65(dat(4097),jz-4096,cfile6, - + NClearAve,MinSigdB,DFTolerance,NFreeze,NAFC,mode65,Nseg, - + MouseDF2,NAgain,ndepth,neme,nsked,idf,idfsh, - + mycall,hiscall,hisgrid,lumsg,lcum,nspecial,ndf, - + nstest,dfsh,snrsh, - + NSyncOK,ccf,psavg,ndiag,nwsh) - - 900 LDecoded = ((NSyncOK.gt.0) .or. npkept.gt.0) - end file 11 - call flushqqq(11) - call flushqqq(12) - call flushqqq(21) - - return - end - diff --git a/xcor.f b/xcor.f deleted file mode 100644 index 44c7c10c2..000000000 --- a/xcor.f +++ /dev/null @@ -1,84 +0,0 @@ - subroutine xcor(s2,ipk,nsteps,nsym,lag1,lag2, - + ccf,ccf0,lagpk,flip,fdot) - -C Computes ccf of a row of s2 and the pseudo-random array pr. Returns -C peak of the CCF and the lag at which peak occurs. For JT65, the -C CCF peak may be either positive or negative, with negative implying -C the "OOO" message. - - parameter (NHMAX=1024) !Max length of power spectra - parameter (NSMAX=320) !Max number of half-symbol steps - real s2(NHMAX,NSMAX) !2d spectrum, stepped by half-symbols - real a(NSMAX),a2(NSMAX) - real ccf(-5:540) - include 'prcom.h' - common/clipcom/ nclip - data lagmin/0/ !Silence g77 warning - save - - df=11025.0/4096. - dtstep=0.5/df - fac=dtstep/(60.0*df) - - do j=1,nsteps - ii=nint((j-nsteps/2)*fdot*fac)+ipk - a(j)=s2(ii,j) - enddo - -C If requested, clip the spectrum that will be cross correlated. - nclip=0 !Turn it off - if(nclip.gt.0) then - call pctile(a,a2,nsteps,50,base) - alow=a2(nint(nsteps*0.16)) - ahigh=a2(nint(nsteps*0.84)) - rms=min(base-alow,ahigh-base) - clip=4.0-nclip - atop=base+clip*rms - abot=base-clip*rms - do i=1,nsteps - if(nclip.lt.4) then - a(i)=min(a(i),atop) - a(i)=max(a(i),abot) - else - if(a(i).ge.base) then - a(i)=1.0 - else - a(i)=-1.0 - endif - endif - enddo - endif - - ccfmax=0. - ccfmin=0. - do lag=lag1,lag2 - x=0. - do i=1,nsym - j=2*i-1+lag - if(j.ge.1 .and. j.le.nsteps) x=x+a(j)*pr(i) - enddo - ccf(lag)=2*x !The 2 is for plotting scale - if(ccf(lag).gt.ccfmax) then - ccfmax=ccf(lag) - lagpk=lag - endif - - if(ccf(lag).lt.ccfmin) then - ccfmin=ccf(lag) - lagmin=lag - endif - enddo - - ccf0=ccfmax - flip=1.0 - if(-ccfmin.gt.ccfmax) then - do lag=lag1,lag2 - ccf(lag)=-ccf(lag) - enddo - lagpk=lagmin - ccf0=-ccfmin - flip=-1.0 - endif - - return - end