From 68de5dd4bdab42648352f2affa361c0c71254069 Mon Sep 17 00:00:00 2001
From: Joe Taylor <k1jt@arrl.org>
Date: Thu, 4 Oct 2012 19:03:39 +0000
Subject: [PATCH] Better symspec now.

git-svn-id: svn+ssh://svn.code.sf.net/p/wsjt/wsjt/branches/wsjtx@2646 ab8295b8-cf94-4d9e-aec4-7959e3be5d79
---
 getfile.cpp     |  1 -
 lib/jt9sim.f90  |  3 ++
 lib/symspec.f90 | 88 +++++++++++++++++++++----------------------------
 lib/timf2.f90   |  1 -
 mainwindow.cpp  |  4 +--
 5 files changed, 41 insertions(+), 56 deletions(-)

diff --git a/getfile.cpp b/getfile.cpp
index ae2e528b6..e286e00f8 100644
--- a/getfile.cpp
+++ b/getfile.cpp
@@ -19,7 +19,6 @@ void getfile(QString fname, int ntrperiod)
     fread(jt9com_.d2,1,44,fp);
     int nrd=fread(jt9com_.d2,2,npts,fp);
     fclose(fp);
-//    for(int i=0; i<npts; i++) jt9com_.d2[i]/=100;
   }
 }
 
diff --git a/lib/jt9sim.f90 b/lib/jt9sim.f90
index 0b49cb6ad..99aabc815 100644
--- a/lib/jt9sim.f90
+++ b/lib/jt9sim.f90
@@ -63,6 +63,8 @@ program jt9sim
   if(minutes.eq.30) nsps=252000
   if(nsps.eq.0) stop 'Bad value for minutes.'
   ihdr=0                             !Temporary ###
+
+  if(msg0(1:3).eq.'sin') read(msg0(4:),*) sinfreq
   
   open(12,file='msgs.txt',status='old')
 
@@ -111,6 +113,7 @@ program jt9sim
         k=12000                             !Start at t = 1 s
         do isym=1,85
            freq=f + itone(isym)*baud
+           if(msg0(1:3).eq.'sin') freq=sinfreq
            dphi=twopi*freq*dt
            do i=1,nsps
               phi=phi + dphi
diff --git a/lib/symspec.f90 b/lib/symspec.f90
index 53715d455..90434b3ef 100644
--- a/lib/symspec.f90
+++ b/lib/symspec.f90
@@ -26,20 +26,14 @@ subroutine symspec(k,ntrperiod,nsps,ndiskdat,nb,nbslider,pxdb,s,f0a,df3,    &
   real*4 s(NSMAX),w(NFFT1),w3(MAXFFT3)
   real*4 stmp(NFFT2/2)
   real*4 x0(NFFT1),x1(NFFT1)
-  real*4 x2(NFFT2)
-  complex cx2(0:NFFT2/2)
-  complex cx2a(NFFT2A)
-  complex z,zfac
-  complex zsumx
-  complex cx(MAXFFT3)
-  complex cx00(NFFT1)
-  complex cx0(0:1023),cx1(0:1023)
+  real*4 x2(NFFT1+105)
+  real*4 xx(NMAX)
+  complex cx(0:MAXFFT3-1)
   logical*1 lstrong(0:1023)               !Should be (0:512)
   integer*2 id2
   complex c0
   common/jt9com/id2(NMAX),ss(184,NSMAX),savg(NSMAX),c0(NDMAX),     &
        nutc,npts8,junk(20)
-  equivalence (x2,cx2)
   data rms/999.0/,k0/99999999/,ntrperiod0/0/,nfft3z/0/
   save
 
@@ -58,22 +52,22 @@ subroutine symspec(k,ntrperiod,nsps,ndiskdat,nb,nbslider,pxdb,s,f0a,df3,    &
   if(nfft3.ne.nfft3z) then
      pi=4.0*atan(1.0)
      do i=1,nfft3
-        w3(i)=(sin(i*pi/nfft3))**2             !Window for nfft3
+        w3(i)=2.0*(sin(i*pi/nfft3))**2             !Window for nfft3
      enddo
      stmp=0.
      nfft3z=nfft3
   endif
 
   if(k.lt.k0) then
-     ja=-2*jstep
+     ja=-3*jstep
      savg=0.
      ihsym=0
      k1=0
      k8=0
-     if(ndiskdat.eq.0) id2(k+1:)=0.        !### Should not be needed ??? ###
+     x2=0.
   endif
   k0=k
-
+ 
   nzap=0
   sigmas=1.5*(10.0**(0.01*nbslider)) + 0.7
   peaklimit=sigmas*max(10.0,rms)
@@ -81,77 +75,69 @@ subroutine symspec(k,ntrperiod,nsps,ndiskdat,nb,nbslider,pxdb,s,f0a,df3,    &
   px=0.
   df2=12000.0/NFFT2
 
-!  nwindow=2
-  nwindow=0                                    !### No windowing ###
+  nwindow=2
+!  nwindow=0                                    !### No windowing ###
   kstep1=NFFT1
   if(nwindow.ne.0) kstep1=NFFT1/2
-  fac=1.0/(NFFT1*NFFT2)
+  fac=2.0/NFFT1
   nblks=(k-k1)/kstep1
   do nblk=1,nblks
+     j=k1+1
      do i=1,NFFT1
-        x0(i)=fac*id2(k1+i)
+        x0(i)=id2(k1+i)
      enddo
      call timf2(x0,k,NFFT1,nwindow,nb,peaklimit,faclim,x1,   &
           slimit,lstrong,px,nzap)
-!     x1=x0
-     x2=x1
-     call four2a(x2,NFFT2,1,-1,0)              !Second forward FFT, r2c
-
-     i0=nint(1000.0/df2)
-     f0a=i0*df2
-     cx2a(1:NFFT2A/2)=cx2(i0:NFFT2A/2+i0-1)
-     cx2a(NFFT2A/2+1:NFFT2A)=cx2(i0-1-NFFT2A/2:i0-1)
-     call four2a(cx2a,NFFT2A,1,1,1)
-
-     c0(k8+1:k8+NFFT2A)=cx2a
-     npts8=k8+NFFT2A
-
-!###                                   Test for gliches at multiples of 128
-!     if(k8.lt.1000) then
-!        do i=k8+1,k8+NFFT2A
-!           write(82,4002) i,c0(i)
-!4002       format(i8,2e12.3)
-!        enddo
-!     endif
-!###
 
+! Mix at 1500 Hz, lowpass at +/-750 Hz, and downsample to 1500 Hz complex.
+     x2(106:105+kstep1)=x1(1:kstep1)
+     call fil3(x2,kstep1+105,c0(k8+1),n2)
+     x2(1:105)=x1(kstep1-104:kstep1)   !Save 105 trailing samples
      k1=k1+kstep1
      k8=k8+kstep1/8
   enddo
 
+  npts8=k8
   ja=ja+jstep                         !Index of first sample
-  if(ja.lt.0) go to 999
-  do i=1,nfft3                          !Copy data into cx
-     cx(i)=c0(ja+i)
+  if(ja.lt.0 .or. npts8.lt.ja+nfft3) go to 999
+  do i=0,nfft3-1                      !Copy data into cx
+     cx(i)=c0(ja+i+1)
   enddo
 
 !###
 !  if(nzap/178.lt.50 .and. (ndiskdat.eq.0 .or. ihsym.lt.280)) then
-     nsum=nblks*kstep1 - nzap
-     if(nsum.le.0) nsum=1
-     rms=sqrt(0.5*px/nsum)
+  nsum=nblks*kstep1 - nzap
+  if(nsum.le.0) nsum=1
+  rms=sqrt(0.5*px/nsum)
 !  endif
   pxdb=0.
-  pxdb=20.0*log10(rms) + 166.0           !### ??? ###
+  if(rms.gt.0.0) pxdb=20.0*log10(rms)
   if(pxdb.gt.60.0) pxdb=60.0
-!  write(71,*) kstep1,nzap,nsum,px,rms,pxdb
-!  flush(71)
-
 !###
 
 
   ihsym=ihsym+1
-  call four2a(cx,nfft3,1,-1,1)           !Third forward FFT (X)
+  cx(0:nfft3-1)=w3(1:nfft3)*cx(0:nfft3-1)  !Apply window w3
+  call four2a(cx,nfft3,1,1,1)           !Third forward FFT (X)
 
   n=min(184,ihsym)
   df3=1500.0/nfft3
+  i0=nint(-500.0/df3)
   iz=min(NSMAX,nint(1000.0/df3))
+  fac=(1.0/nfft3)**2
   do i=1,iz
-     sx=real(cx(i))**2 + aimag(cx(i))**2  
+     j=i0+i-1
+     if(j.lt.0) j=j+nfft3
+     sx=fac*(real(cx(j))**2 + aimag(cx(j))**2)
      ss(n,i)=sx
      savg(i)=savg(i) + sx
      s(i)=sx
   enddo
 
-999 return
+999 continue
+!  write(71,3003) k,nsum,nzap,px,rms,pxdb
+!3003 format(3i10,3f12.2)
+!  flush(71)
+
+  return
 end subroutine symspec
diff --git a/lib/timf2.f90 b/lib/timf2.f90
index eecbb2f3d..0debfba0d 100644
--- a/lib/timf2.f90
+++ b/lib/timf2.f90
@@ -147,7 +147,6 @@ subroutine timf2(x0,k,nfft,nwindow,nb,peaklimit,faclim,x1,     &
   do i=0,kstep-1
      px=px + xw(i)*xw(i)
   enddo
-  px=px/kstep
 
   x1(0:kstep-1)=xw(0:kstep-1) + xs(0:kstep-1)     !Recombine weak + strong
 
diff --git a/mainwindow.cpp b/mainwindow.cpp
index b0b5e3268..b57af5fe4 100644
--- a/mainwindow.cpp
+++ b/mainwindow.cpp
@@ -1,4 +1,4 @@
-//----------------------------------------------------------- MainWindow
+//------------------------------------------------------------ MainWindow
 #include "mainwindow.h"
 #include "ui_mainwindow.h"
 #include "devsetup.h"
@@ -314,7 +314,6 @@ void MainWindow::dataSink(int k)
   static uchar lstrong[1024];
   static float slimit;
 
-
   if(m_diskData) {
     ndiskdat=1;
     jt9com_.ndiskdat=1;
@@ -719,7 +718,6 @@ void MainWindow::diskDat()                                   //diskDat()
 {
   int kstep=m_nsps/2;
   m_diskData=true;
-
   for(int n=1; n<=184; n++) {              // Do the half-symbol FFTs
     int k=(n+1)*kstep;
     dataSink(k);