From 0b68d2c7bb1f4be4f057d9ac8cfcf9c3d9c1f624 Mon Sep 17 00:00:00 2001 From: Joe Taylor Date: Tue, 11 Sep 2012 15:31:57 +0000 Subject: [PATCH] Fix a flaw in decode1a. git-svn-id: svn+ssh://svn.code.sf.net/p/wsjt/wsjt/branches/map65@2571 ab8295b8-cf94-4d9e-aec4-7959e3be5d79 --- libm65/decode1a.f | 165 -------------------------------------------- libm65/decode1a.f90 | 156 +++++++++++++++++++++++++++++++++++++++++ libm65/ftninit.f90 | 5 +- libm65/s3avg.f90 | 4 -- mainwindow.cpp | 2 +- map65.iss | 2 +- 6 files changed, 161 insertions(+), 173 deletions(-) delete mode 100644 libm65/decode1a.f create mode 100644 libm65/decode1a.f90 diff --git a/libm65/decode1a.f b/libm65/decode1a.f deleted file mode 100644 index 112de5243..000000000 --- a/libm65/decode1a.f +++ /dev/null @@ -1,165 +0,0 @@ - subroutine decode1a(dd,newdat,f0,nflip,mode65,nfast,nfsample,xpol, - + mycall,hiscall,hisgrid,neme,ndepth,nqd,dphi,ndphi,iloop, - + nutc,nkhz,ndf,ipol,ntol,sync2,a,dt,pol,nkv,nhist,nsum,nsave, - + qual,decoded) - -! Apply AFC corrections to a candidate JT65 signal, then decode it. - - parameter (NMAX=60*96000) !Samples per 60 s - real*4 dd(4,NMAX) !92 MB: raw data from Linrad timf2 - complex cx(NMAX/64), cy(NMAX/64) !Data at 1378.125 samples/s - complex c5x(NMAX/256),c5y(NMAX/256),c5tmp(NMAX/256) !Data at 344.53125 Hz - complex c5a(512) - complex z - real s2(66,126) - real s3(64,63),sy(63) - real a(5) - logical first,xpol - character decoded*22 - character mycall*12,hiscall*12,hisgrid*6 - data first/.true./,jjjmin/1000/,jjjmax/-1000/ - data nutc0/-999/,nhz0/-9999999/ - save - -! Mix sync tone to baseband, low-pass filter, downsample to 1378.125 Hz - dt00=dt - call timer('filbig ',0) - call filbig(dd,NMAX,nfast,f0,newdat,nfsample,xpol,cx,cy,n5) -! NB: cx, cy have sample rate 96000*77125/5376000 = 1378.125 Hz - call timer('filbig ',1) - sqa=0. - sqb=0. - do i=1,n5 - sqa=sqa + real(cx(i))**2 + aimag(cx(i))**2 - if(xpol) sqb=sqb + real(cy(i))**2 + aimag(cy(i))**2 - enddo - sqa=sqa/n5 - sqb=sqb/n5 - -! Find best DF, f1, f2, DT, and pol. Start by downsampling to 344.53125 Hz - if(xpol) then - z=cmplx(cos(dphi),sin(dphi)) - cy(:n5)=z*cy(:n5) !Adjust for cable length difference - endif - call timer('fil6521 ',0) -! Add some zeros at start of c5 arrays -- empirical fix for negative DT's - nadd=1089 - c5x(:nadd)=0. - call fil6521(cx,n5,c5x(nadd+1),n6) - if(xpol) then - c5y(:nadd)=0. - call fil6521(cy,n5,c5y(nadd+1),n6) - endif - n6=n6+nadd - call timer('fil6521 ',1) - - fsample=1378.125/4. - a(5)=dt00 - i0=nint((a(5)+0.5)*fsample) - 2 + nadd - i00=i0 - if(i0.lt.1) then - write(13,*) 'i0 too small in decode1a:',i0,f0 - flush(13) - i0=1 - endif - nz=n6+1-i0 - -! We're looking only at sync tone here... so why not downsample by another -! factor of 1/8, say? Should be a significant execution speed-up. - call timer('afc65b ',0) -! Best fit for DF, f1, f2, pol - call afc65b(c5x(i0),c5y(i0),nz,nfast,fsample,nflip,ipol,xpol, - + ndphi,iloop,a,ccfbest,dtbest) - call timer('afc65b ',1) - - pol=a(4)/57.2957795 - aa=cos(pol) - bb=sin(pol) - sq0=aa*aa*sqa + bb*bb*sqb - sync2=3.7*ccfbest/sq0 - -! Apply AFC corrections to the time-domain signal -! Now we are back to using the 1378.125 Hz sample rate, enough to -! accommodate the full JT65C bandwidth. - - call timer('twkfreq ',0) - call twkfreq(cx,cy,n5,a) - call timer('twkfreq ',1) - -! Compute spectrum at best polarization for each half symbol. -! Adding or subtracting a small number (e.g., 5) to j may make it decode.\ -! NB: might want to try computing full-symbol spectra (nfft=512, even for -! submodes B and C). - - nsym=126 -! nfft=512/(nfast*mode65) - nfft=512/nfast - j=(dt00+dtbest+2.685)*1378.125 - j00=j - if(nfast.eq.2) j=j-1506 -! print*,'B',dt00,dtbest,j - if(j.lt.0) j=0 - - call timer('sh_ffts ',0) - -! Perhaps should try full-symbol-length FFTs even in B, C sub-modes? -! (Tried this, found no significant difference in decodes.) - - do k=1,nsym -! do n=1,mode65 - do n=1,1 - do i=1,nfft - j=j+1 - c5a(i)=aa*cx(j) + bb*cy(j) - enddo - call four2a(c5a,nfft,1,1,1) - if(n.eq.1) then - do i=1,66 -! s2(i,k)=real(c5a(i))**2 + aimag(c5a(i))**2 - jj=i - if(nfast.eq.1 .and. mode65.eq.2) jj=2*i-1 - if(nfast.eq.2 .and. mode65.eq.4) jj=2*i-1 - if(nfast.eq.1 .and. mode65.eq.4) jj=4*i-3 - s2(i,k)=real(c5a(jj))**2 + aimag(c5a(jj))**2 - enddo - else - do i=1,66 - s2(i,k)=s2(i,k) + real(c5a(i))**2 + aimag(c5a(i))**2 - enddo - endif - enddo - enddo - - call timer('sh_ffts ',1) - - flip=nflip - call timer('dec65b ',0) - call decode65b(s2,flip,mycall,hiscall,hisgrid,mode65,neme,ndepth, - + nqd,nkv,nhist,qual,decoded,s3,sy) - dt=dt00 + dtbest + 1.7 - if(nfast.eq.2) dt=dt00 + dtbest + 0.6 - call timer('dec65b ',1) - - if(nqd.eq.1 .and. nkv.eq.0) then - nhz=1000*nkhz + ndf - ihzdiff=min(500,ntol) - if(nutc.ne.nutc0 .or. abs(nhz-nhz0).ge.ihzdiff) syncbest=0. - if(sync2.gt.0.99999*syncbest) then - nsave=nsave+1 - nsave=mod(nsave-1,64)+1 - npol=nint(57.296*pol) - - call s3avg(nsave,mode65,nutc,nhz,xdt,npol,ntol,s3,nsum, - + nkv,decoded) - syncbest=sync2 - nhz0=nhz - endif - nutc0=nutc - endif - - write(71,3000) dt00,i00,j00,dtbest,sync2,qual,decoded - 3000 format(f7.2,2i6,3f7.2,2x,a22) - flush(71) - - return - end diff --git a/libm65/decode1a.f90 b/libm65/decode1a.f90 new file mode 100644 index 000000000..8a616e1b2 --- /dev/null +++ b/libm65/decode1a.f90 @@ -0,0 +1,156 @@ +subroutine decode1a(dd,newdat,f0,nflip,mode65,nfast,nfsample,xpol, & + mycall,hiscall,hisgrid,neme,ndepth,nqd,dphi,ndphi,iloop, & + nutc,nkhz,ndf,ipol,ntol,sync2,a,dt,pol,nkv,nhist,nsum,nsave, & + qual,decoded) + +! Apply AFC corrections to a candidate JT65 signal, then decode it. + + parameter (NMAX=60*96000) !Samples per 60 s + real*4 dd(4,NMAX) !92 MB: raw data from Linrad timf2 + complex cx(NMAX/64), cy(NMAX/64) !Data at 1378.125 samples/s + complex c5x(NMAX/256),c5y(NMAX/256) !Data at 344.53125 Hz + complex c5a(512) + complex z + real s2(66,126) + real s3(64,63),sy(63) + real a(5) + logical first,xpol + character decoded*22 + character mycall*12,hiscall*12,hisgrid*6 + data first/.true./,jjjmin/1000/,jjjmax/-1000/ + data nutc0/-999/,nhz0/-9999999/ + save + +! Mix sync tone to baseband, low-pass filter, downsample to 1378.125 Hz + dt00=dt + call timer('filbig ',0) + call filbig(dd,NMAX,nfast,f0,newdat,nfsample,xpol,cx,cy,n5) +! NB: cx, cy have sample rate 96000*77125/5376000 = 1378.125 Hz + call timer('filbig ',1) + sqa=0. + sqb=0. + do i=1,n5 + sqa=sqa + real(cx(i))**2 + aimag(cx(i))**2 + if(xpol) sqb=sqb + real(cy(i))**2 + aimag(cy(i))**2 + enddo + sqa=sqa/n5 + sqb=sqb/n5 + +! Find best DF, f1, f2, DT, and pol. Start by downsampling to 344.53125 Hz + if(xpol) then + z=cmplx(cos(dphi),sin(dphi)) + cy(:n5)=z*cy(:n5) !Adjust for cable length difference + endif + call timer('fil6521 ',0) +! Add some zeros at start of c5 arrays -- empirical fix for negative DT's + nadd=1089 + c5x(:nadd)=0. + call fil6521(cx,n5,c5x(nadd+1),n6) + if(xpol) then + c5y(:nadd)=0. + call fil6521(cy,n5,c5y(nadd+1),n6) + endif + n6=n6+nadd + call timer('fil6521 ',1) + + fsample=1378.125/4. + a(5)=dt00 + i0=nint((a(5)+0.5)*fsample) - 2 + nadd + if(i0.lt.1) then + write(13,*) 'i0 too small in decode1a:',i0,f0 + flush(13) + i0=1 + endif + nz=n6+1-i0 + +! We're looking only at sync tone here... so why not downsample by another +! factor of 1/8, say? Should be a significant execution speed-up. + call timer('afc65b ',0) +! Best fit for DF, f1, f2, pol + call afc65b(c5x(i0),c5y(i0),nz,nfast,fsample,nflip,ipol,xpol, & + ndphi,iloop,a,ccfbest,dtbest) + call timer('afc65b ',1) + + pol=a(4)/57.2957795 + aa=cos(pol) + bb=sin(pol) + sq0=aa*aa*sqa + bb*bb*sqb + sync2=3.7*ccfbest/sq0 + +! Apply AFC corrections to the time-domain signal +! Now we are back to using the 1378.125 Hz sample rate, enough to +! accommodate the full JT65C bandwidth. + + call timer('twkfreq ',0) + call twkfreq(cx,cy,n5,a) + call timer('twkfreq ',1) + +! Compute spectrum at best polarization for each half symbol. +! Adding or subtracting a small number (e.g., 5) to j may make it decode.\ +! NB: might want to try computing full-symbol spectra (nfft=512, even for +! submodes B and C). + + nsym=126 + nfft=512/nfast + j=(dt00+dtbest+2.685)*1378.125 + if(nfast.eq.2) j=j-1506 + if(j.lt.0) j=0 + + call timer('sh_ffts ',0) + +! Perhaps should try full-symbol-length FFTs even in B, C sub-modes? +! (Tried this, found no significant difference in decodes.) + + do k=1,nsym +! do n=1,mode65 + do n=1,1 + do i=1,nfft + j=j+1 + c5a(i)=aa*cx(j) + bb*cy(j) + enddo + call four2a(c5a,nfft,1,1,1) + if(n.eq.1) then + do i=1,66 +! s2(i,k)=real(c5a(i))**2 + aimag(c5a(i))**2 + jj=i + if(nfast.eq.1 .and. mode65.eq.2) jj=2*i-1 + if(nfast.eq.2 .and. mode65.eq.4) jj=2*i-1 + if(nfast.eq.1 .and. mode65.eq.4) jj=4*i-3 + s2(i,k)=real(c5a(jj))**2 + aimag(c5a(jj))**2 + enddo + else + do i=1,66 + s2(i,k)=s2(i,k) + real(c5a(i))**2 + aimag(c5a(i))**2 + enddo + endif + enddo + enddo + + call timer('sh_ffts ',1) + + flip=nflip + call timer('dec65b ',0) + call decode65b(s2,flip,mycall,hiscall,hisgrid,mode65,neme,ndepth, & + nqd,nkv,nhist,qual,decoded,s3,sy) + dt=dt00 + dtbest + 1.7 + if(nfast.eq.2) dt=dt00 + dtbest + 0.6 + call timer('dec65b ',1) + + if(nqd.eq.1 .and. decoded.eq.' ') then + nhz=1000*nkhz + ndf + ihzdiff=min(500,ntol) + if(nutc.ne.nutc0 .or. abs(nhz-nhz0).ge.ihzdiff) syncbest=0. + if(sync2.gt.0.99999*syncbest) then + nsave=nsave+1 + nsave=mod(nsave-1,64)+1 + npol=nint(57.296*pol) + + call s3avg(nsave,mode65,nutc,nhz,xdt,npol,ntol,s3,nsum,nkv,decoded) + syncbest=sync2 + nhz0=nhz + endif + nutc0=nutc + endif + + return +end subroutine decode1a diff --git a/libm65/ftninit.f90 b/libm65/ftninit.f90 index 9a9a2f950..056f4a465 100644 --- a/libm65/ftninit.f90 +++ b/libm65/ftninit.f90 @@ -23,7 +23,7 @@ subroutine ftninit(appd) character*(*) appd - character cjunk*1,firstline*30 + character firstline*30 character addpfx*8 integer junk(256) common/pfxcom/addpfx @@ -53,7 +53,8 @@ subroutine ftninit(appd) if(isuccess.ne.0) write(13,1010) firstline 1010 format('Imported FFTW wisdom: ',a30) -30 return +30 flush(13) + return 920 write(0,*) '!Error opening timer.out' stop diff --git a/libm65/s3avg.f90 b/libm65/s3avg.f90 index f9e712f8d..9d76637c4 100644 --- a/libm65/s3avg.f90 +++ b/libm65/s3avg.f90 @@ -44,10 +44,6 @@ subroutine s3avg(nsave,mode65,nutc,nhz,xdt,npol,ntol,s3,nsum,nkv,decoded) s3b=s3b + s3a(1:64,1:63,i) nsum=nsum+1 enddo - -! write(71,3001) nutc,nsum,nsave,nhz,npol,xdt -!3001 format(5i8,1f8.1) -! flush(71) decoded=' ' if(nsum.ge.2) then !Try decoding the sverage diff --git a/mainwindow.cpp b/mainwindow.cpp index d87c8d1c8..cdcfd7901 100644 --- a/mainwindow.cpp +++ b/mainwindow.cpp @@ -1,4 +1,4 @@ -//--------------------------------------------------------------- MainWindow +//-------------------------------------------------------------- MainWindow #include "mainwindow.h" #include "ui_mainwindow.h" #include "devsetup.h" diff --git a/map65.iss b/map65.iss index c07dc2afb..5bb86ffe0 100644 --- a/map65.iss +++ b/map65.iss @@ -1,6 +1,6 @@ [Setup] AppName=MAP65 -AppVerName=MAP65 Version 2.3.0 r2471 +AppVerName=MAP65 Version 2.4.0 r2570 AppCopyright=Copyright (C) 2001-2012 by Joe Taylor, K1JT DefaultDirName=c:\MAP65 DefaultGroupName=MAP65