mirror of
https://github.com/saitohirga/WSJT-X.git
synced 2025-07-30 12:42:27 -04:00
Improved estimates of Doppler spread. Comment the code in write_ref().
This commit is contained in:
parent
f61cb1dcbd
commit
7bd797c0e9
@ -811,88 +811,99 @@ contains
|
|||||||
end subroutine get_candidates_fst240
|
end subroutine get_candidates_fst240
|
||||||
|
|
||||||
subroutine write_ref(itone,iwave,nsps,nmax,ndown,hmod,i0,fc,fmid,w50)
|
subroutine write_ref(itone,iwave,nsps,nmax,ndown,hmod,i0,fc,fmid,w50)
|
||||||
complex cwave(nmax)
|
|
||||||
complex, allocatable :: c(:)
|
|
||||||
real,allocatable :: ss(:)
|
|
||||||
integer itone(160)
|
|
||||||
integer*2 iwave(nmax)
|
|
||||||
integer hmod
|
|
||||||
data ncall/0/
|
|
||||||
save ncall
|
|
||||||
|
|
||||||
ncall=ncall+1
|
! On "plotspec" special request, compute Doppler spread for a decoded signal
|
||||||
allocate(c(0:nmax-1))
|
|
||||||
wave=0
|
|
||||||
fsample=12000.0
|
|
||||||
nsym=160
|
|
||||||
|
|
||||||
call gen_fst240wave(itone,nsym,nsps,nmax,fsample,hmod,fc, &
|
complex, allocatable :: cwave(:) !Reconstructed complex signal
|
||||||
1,cwave,wave)
|
complex, allocatable :: g(:) !Channel gain, g(t) in QEX paper
|
||||||
cwave=cshift(cwave,-i0*ndown)
|
real,allocatable :: ss(:) !Computed power spectrum of g(t)
|
||||||
! do i=1,nmax
|
integer itone(160) !Tones for this message
|
||||||
! write(51,1000) i,iwave(i),cwave(i)
|
integer*2 iwave(nmax) !Raw Rx data
|
||||||
!1000 format(2i10,f12.6)
|
integer hmod !Modulation index
|
||||||
! enddo
|
data ncall/0/
|
||||||
|
save ncall
|
||||||
|
|
||||||
fac=1.0/32768
|
ncall=ncall+1
|
||||||
c=fac*float(iwave)*conjg(cwave)
|
nfft=2*nmax
|
||||||
call four2a(c,nmax,1,-1,1) !Forward c2c FFT
|
allocate(cwave(0:nmax-1))
|
||||||
|
allocate(g(0:nfft-1))
|
||||||
|
wave=0
|
||||||
|
fsample=12000.0
|
||||||
|
nsym=160
|
||||||
|
call gen_fst240wave(itone,nsym,nsps,nmax,fsample,hmod,fc,1,cwave,wave)
|
||||||
|
cwave=cshift(cwave,-i0*ndown)
|
||||||
|
fac=1.0/32768
|
||||||
|
g(0:nmax-1)=fac*float(iwave)*conjg(cwave)
|
||||||
|
g(nmax:)=0.
|
||||||
|
call four2a(g,nfft,1,-1,1) !Forward c2c FFT
|
||||||
|
|
||||||
df=12000.0/nmax
|
df=12000.0/nfft
|
||||||
ia=1.0/df
|
ia=1.0/df
|
||||||
smax=0.
|
smax=0.
|
||||||
do i=-ia,ia
|
do i=-ia,ia !Find smax in +/- 1 Hz around 0.
|
||||||
j=i
|
j=i
|
||||||
if(j.lt.0) j=i+nmax
|
if(j.lt.0) j=i+nfft
|
||||||
s=real(c(j))**2 + aimag(c(j))**2
|
s=real(g(j))**2 + aimag(g(j))**2
|
||||||
smax=max(s,smax)
|
smax=max(s,smax)
|
||||||
enddo
|
enddo
|
||||||
ia=10.1/df
|
|
||||||
allocate(ss(-ia:ia))
|
ia=10.1/df
|
||||||
sum1=0.
|
allocate(ss(-ia:ia)) !Allocate space for +/- 10 Hz
|
||||||
sum2=0.
|
sum1=0.
|
||||||
ns=0
|
sum2=0.
|
||||||
do i=-ia,ia
|
ns=0
|
||||||
j=i
|
do i=-ia,ia
|
||||||
if(j.lt.0) j=i+nmax
|
j=i
|
||||||
ss(i)=(real(c(j))**2 + aimag(c(j))**2)/smax
|
if(j.lt.0) j=i+nfft
|
||||||
f=i*df
|
ss(i)=(real(g(j))**2 + aimag(g(j))**2)/smax
|
||||||
if(f.ge.-4.0 .and. f.le.-2.0) then
|
f=i*df
|
||||||
sum1=sum1 + ss(i)
|
if(f.ge.-4.0 .and. f.le.-2.0) then
|
||||||
ns=ns+1
|
sum1=sum1 + ss(i) !Power between -2 and -4 Hz
|
||||||
else if(f.ge.2.0 .and. f.le.4.0) then
|
ns=ns+1
|
||||||
sum2=sum2 + ss(i)
|
else if(f.ge.2.0 .and. f.le.4.0) then
|
||||||
endif
|
sum2=sum2 + ss(i) !Power between +2 and +4 Hz
|
||||||
enddo
|
endif
|
||||||
avg=min(sum1/ns,sum2/ns)
|
enddo
|
||||||
|
avg=min(sum1/ns,sum2/ns) !Compute avg from smaller sum
|
||||||
|
|
||||||
sum1=0.
|
sum1=0.
|
||||||
do i=-ia,ia
|
do i=-ia,ia
|
||||||
f=i*df
|
f=i*df
|
||||||
if(abs(f).le.1.0) sum1=sum1 + ss(i)-avg
|
if(abs(f).le.1.0) sum1=sum1 + ss(i)-avg !Power in abs(f) < 1 Hz
|
||||||
y=0.99*ss(i) + ncall-1
|
enddo
|
||||||
write(52,1010) f,y
|
|
||||||
1010 format(f12.6,f12.6)
|
|
||||||
enddo
|
|
||||||
|
|
||||||
ia=nint(1.0/df)
|
ia=nint(1.0/df) + 1
|
||||||
sum2=0.0
|
sum2=0.0
|
||||||
i1=-999
|
xi1=-999
|
||||||
i2=-999
|
xi2=-999
|
||||||
i3=-999
|
xi3=-999
|
||||||
do i=-ia,ia
|
sum2z=0.
|
||||||
sum2=sum2 + ss(i)-avg
|
do i=-ia,ia !Find freq range that has 50% of signal power
|
||||||
if(sum2.ge.0.25*sum1 .and. i1.eq.-999) i1=i
|
sum2=sum2 + ss(i)-avg
|
||||||
if(sum2.ge.0.50*sum1 .and. i2.eq.-999) i2=i
|
if(sum2.ge.0.25*sum1 .and. xi1.eq.-999.0) then
|
||||||
if(sum2.ge.0.75*sum1) then
|
xi1=i - 1 + (sum2-0.25*sum1)/(sum2-sum2z)
|
||||||
i3=i
|
endif
|
||||||
exit
|
if(sum2.ge.0.50*sum1 .and. xi2.eq.-999.0) then
|
||||||
endif
|
xi2=i - 1 + (sum2-0.50*sum1)/(sum2-sum2z)
|
||||||
enddo
|
endif
|
||||||
fmid=i2*df
|
if(sum2.ge.0.75*sum1) then
|
||||||
w50=(i3-i1+1)*df
|
xi3=i - 1 + (sum2-0.75*sum1)/(sum2-sum2z)
|
||||||
|
exit
|
||||||
return
|
endif
|
||||||
end subroutine write_ref
|
sum2z=sum2
|
||||||
|
enddo
|
||||||
|
xdiff=sqrt(1.0+(xi3-xi1)**2) !Keep small values from fluctuating too widely
|
||||||
|
w50=xdiff*df !Compute Doppler spread
|
||||||
|
fmid=xi2*df !Frequency midpoint of signal powere
|
||||||
|
|
||||||
|
do i=-ia,ia !Save the spectrum for plotting
|
||||||
|
f=i*df
|
||||||
|
y=0.99*ss(i+nint(xi2)) + ncall-1
|
||||||
|
write(52,1010) f,y
|
||||||
|
1010 format(f12.6,f12.6)
|
||||||
|
enddo
|
||||||
|
|
||||||
|
return
|
||||||
|
end subroutine write_ref
|
||||||
|
|
||||||
end module fst240_decode
|
end module fst240_decode
|
||||||
|
Loading…
x
Reference in New Issue
Block a user