From 6a59223b707d57afe99202bc4bb70b96abb8604e Mon Sep 17 00:00:00 2001
From: Joe Taylor <joe@princeton.edu>
Date: Tue, 27 Feb 2024 11:39:48 -0500
Subject: [PATCH] Add quick-and-dirty measurement of sync width (Doppler
 spread).

---
 lib/superfox/sfox_sync.f90 | 14 +++++++++++++-
 lib/superfox/sfoxtest.f90  | 13 ++++++++-----
 2 files changed, 21 insertions(+), 6 deletions(-)

diff --git a/lib/superfox/sfox_sync.f90 b/lib/superfox/sfox_sync.f90
index 17baf1d25..09bf23869 100644
--- a/lib/superfox/sfox_sync.f90
+++ b/lib/superfox/sfox_sync.f90
@@ -1,4 +1,4 @@
-subroutine sfox_sync(iwave,fsample,isync,f,t)
+subroutine sfox_sync(iwave,fsample,isync,f,t,fwidth)
 
   use sfox_mod
   parameter (NSTEP=8)
@@ -130,11 +130,23 @@ subroutine sfox_sync(iwave,fsample,isync,f,t)
   enddo
   ipeak2=maxloc(s2)
   ipk=ipeak2(1)-1
+
   dxi=0.
   if(ipk.gt.1 .and. ipk.lt.nfft/4) then
      call peakup(s2(ipk-1),s2(ipk),s2(ipk+1),dxi)
   endif
   f=(ipk+dxi)*df2 + bw/2.0
 
+  call pctile(s2(ipk-100:ipk+100),201,48,base)
+  s2=s2-base
+!  do i=ipk-100,ipk+100
+!     write(51,3051) i*df2,s2(i)
+!3051 format(2f15.3)
+!  enddo
+  smax=maxval(s2(ipk-10:ipk+10))
+  w=count(s2(ipk-10:ipk+10).gt.0.5*smax)
+  fwidth=0.
+  if(w.gt.4.0) fwidth=sqrt(w*w - 4*4)*df2
+  
   return
 end subroutine sfox_sync
diff --git a/lib/superfox/sfoxtest.f90 b/lib/superfox/sfoxtest.f90
index e23c24b8f..4b68add1c 100644
--- a/lib/superfox/sfoxtest.f90
+++ b/lib/superfox/sfoxtest.f90
@@ -162,6 +162,7 @@ program sfoxtest
      nworst=0
      sqt=0.
      sqf=0.
+     sumw=0.
 
      do ifile=1,nfiles
         xnoise=0.
@@ -203,7 +204,7 @@ program sfoxtest
            t=xdt
         else
            call timer('sync    ',0)
-           call sfox_sync(iwave,fsample,isync,f,t) ! Find signal freq and DT
+           call sfox_sync(iwave,fsample,isync,f,t,fwidth) !Find freq, DT, width
            call timer('sync    ',1)
         endif
         ferr=f-f1
@@ -215,6 +216,7 @@ program sfoxtest
            ngoodsync=ngoodsync+1
            sqt=sqt + terr*terr
            sqf=sqf + ferr*ferr
+           sumw=sumw+fwidth
         endif
 
         a=0.
@@ -257,14 +259,15 @@ program sfoxtest
      fgoodsync=float(ngoodsync)/nfiles
      fgood=float(ngood)/nfiles
      if(isnr.eq.isnr0) write(*,1300)
-1300 format('    SNR  Eb/No  iters fsync  fgood  averr   rmsf   rmst'/  &
-            '-------------------------------------------------------')
+1300 format('    SNR  Eb/No  iters fsync  fgood  averr   rmsf   rmst  avew'/  &
+            '-------------------------------------------------------------')
      ave_harderr=float(ntot)/nfiles
      rmst=sqrt(sqt/ngoodsync)
      rmsf=sqrt(sqf/ngoodsync)
      ebno=snr-10*log10(baud/2500*mm0*KK/NN)
-     write(*,1310) snr,ebno,nfiles,fgoodsync,fgood,ave_harderr,rmsf,rmst
-1310 format(f7.2,f7.2 i6,2f7.4,f7.1,f7.2,f7.4)
+     avew=sumw/ngoodsync
+     write(*,1310) snr,ebno,nfiles,fgoodsync,fgood,ave_harderr,rmsf,rmst,avew
+1310 format(f7.2,f7.2 i6,2f7.4,f7.1,f7.2,f7.4,f6.1)
      if(fgood.le.0.5 .and. fgood0.gt.0.5) then
         threshold=isnr + 1 - (fgood0-0.50)/(fgood0-fgood+0.000001)
      endif