WSJT-X/lib/spec9.f90
Joe Taylor 7fe6ef1f61 Many changes: all-new JT9 decoder, supposedly faster and better.
Needs thorough testing, especially for drifting signals!


git-svn-id: svn+ssh://svn.code.sf.net/p/wsjt/wsjt/branches/wsjtx@2774 ab8295b8-cf94-4d9e-aec4-7959e3be5d79
2012-11-26 21:06:41 +00:00

134 lines
3.1 KiB
Fortran

subroutine spec9(c0,npts8,nsps,fpk0,fpk,xdt,snrdb,i1SoftSymbols)
parameter (MAXFFT=31500)
complex c0(0:npts8-1)
complex c1(0:2700000)
real*4 ssym(0:7,69)
real*4 sx(0:31500-1)
complex c(0:MAXFFT-1)
integer*1 i1SoftSymbolsScrambled(207)
integer*1 i1SoftSymbols(207)
integer*1 i1
equivalence (i1,i4)
include 'jt9sync.f90'
! Fix up the data in c0()
sum=0.
do i=0,npts8-1
sum=sum + real(c0(i))**2 + aimag(c0(i))**2
enddo
rms=sqrt(sum/npts8)
fac=1.0/rms
twopi=8.0*atan(1.0)
phi=0.
dphi=twopi*500.0/1500.0
do i=0,npts8-1
phi=phi+dphi
if(phi.gt.twopi) phi=phi-twopi
if(phi.lt.-twopi) phi=phi+twopi
c1(i)=fac*cmplx(aimag(c0(i)),real(c0(i)))*cmplx(cos(phi),sin(phi))
enddo
nsps8=nsps/8
foffset=fpk0
istart=1520
call timer('peakdt9 ',0)
call peakdt9(c1,npts8,nsps8,istart,foffset,idt)
call timer('peakdt9 ',1)
istart=istart + nint(0.0625*nsps8*idt)
xdt=istart/1500.0 - 1.0
call timer('peakdf9 ',0)
call peakdf9(c1,npts8,nsps8,istart,foffset,idf)
call timer('peakdf9 ',1)
fpk=fpk0 + idf*0.1*1500.0/nsps8
foffset=foffset + idf*0.1*1500.0/nsps8
twopi=8.0*atan(1.0)
dphi=twopi*foffset/1500.0
nfft=nsps8
nsym=min((npts8-istart)/nsps8,85)
k=0
do j=1,nsym
if(isync(j).eq.1) cycle
k=k+1
ia=(j-1)*nsps8 + istart
! c(0:nfft-1)=c1(ia:ib)
do i=0,nfft-1
c(i)=0.
if(ia+i.ge.0 .and. ia+i.lt.2700000-1) c(i)=c1(ia+i)
enddo
phi=0.
do i=0,nfft-1
phi=phi + dphi
c(i)=c(i) * cmplx(cos(phi),-sin(phi))
enddo
call four2a(c,nfft,1,-1,1)
do i=0,nfft-1
sx(i)=real(c(i))**2 + aimag(c(i))**2
if(i.ge.1 .and. i.le.8) ssym(i-1,k)=sx(i)
enddo
enddo
sum=0.
sig=0.
do j=1,69
smax=0.
do i=0,7
smax=max(smax,ssym(i,j))
sum=sum+ssym(i,j)
enddo
sig=sig+smax
sum=sum-smax
enddo
ave=sum/(69*7)
call pctile(sx,nsps8,50,xmed)
ssym=ssym/ave
sig=sig/69.
df8=1500.0/nsps8
t=max(1.0,sig/xmed - 1.0)
snrdb=db(t) - db(2500.0/df8) - 5.0
m0=3
k=0
do j=1,69
smax=0.
do i=0,7
if(ssym(i,j).gt.smax) then
smax=ssym(i,j)
ipk=i
endif
enddo
do m=m0-1,0,-1 !Get bit-wise soft symbols
if(m.eq.2) then
r1=max(ssym(4,j),ssym(5,j),ssym(6,j),ssym(7,j))
r0=max(ssym(0,j),ssym(1,j),ssym(2,j),ssym(3,j))
else if(m.eq.1) then
r1=max(ssym(2,j),ssym(3,j),ssym(4,j),ssym(5,j))
r0=max(ssym(0,j),ssym(1,j),ssym(6,j),ssym(7,j))
else
r1=max(ssym(1,j),ssym(2,j),ssym(4,j),ssym(7,j))
r0=max(ssym(0,j),ssym(3,j),ssym(5,j),ssym(6,j))
endif
k=k+1
i4=nint(10.0*(r1-r0))
if(i4.lt.-127) i4=-127
if(i4.gt.127) i4=127
i4=i4+128
i1SoftSymbolsScrambled(k)=i1
enddo
enddo
call interleave9(i1SoftSymbolsScrambled,-1,i1SoftSymbols)
return
end subroutine spec9