Fix SNR calculation for B,C,D submodes.

This commit is contained in:
Steven Franke 2020-08-04 09:15:44 -05:00
parent c91baddb2c
commit 55d83f068b
3 changed files with 20 additions and 21 deletions

View File

@ -52,7 +52,7 @@ subroutine get_fst4_bitmetrics(cd,nss,hmod,nmax,nhicoh,bitmetrics,s4,badsync)
do itone=0,3 do itone=0,3
cs(itone,k)=sum(csymb*conjg(c1(:,itone))) cs(itone,k)=sum(csymb*conjg(c1(:,itone)))
enddo enddo
s4(0:3,k)=abs(cs(0:3,k)) s4(0:3,k)=abs(cs(0:3,k))**2
enddo enddo
! Sync quality check ! Sync quality check

View File

@ -1,4 +1,4 @@
subroutine get_fst4_bitmetrics2(cd,nss,hmod,nsizes,bitmetrics,s4hmod,badsync) subroutine get_fst4_bitmetrics2(cd,nss,hmod,nsizes,bitmetrics,s4snr,badsync)
include 'fst4_params.f90' include 'fst4_params.f90'
complex cd(0:NN*nss-1) complex cd(0:NN*nss-1)
@ -15,7 +15,7 @@ subroutine get_fst4_bitmetrics2(cd,nss,hmod,nsizes,bitmetrics,s4hmod,badsync)
logical badsync logical badsync
real bitmetrics(2*NN,4) real bitmetrics(2*NN,4)
real s2(0:65535) real s2(0:65535)
real s4(0:3,NN,4),s4hmod(0:3,NN) real s4(0:3,NN,4),s4snr(0:3,NN)
data isyncword1/0,1,3,2,1,0,2,3/ data isyncword1/0,1,3,2,1,0,2,3/
data isyncword2/2,3,1,0,3,2,0,1/ data isyncword2/2,3,1,0,3,2,0,1/
data graymap/0,1,3,2/ data graymap/0,1,3,2/
@ -121,11 +121,8 @@ subroutine get_fst4_bitmetrics2(cd,nss,hmod,nsizes,bitmetrics,s4hmod,badsync)
call normalizebmet(bitmetrics(:,3),2*NN) call normalizebmet(bitmetrics(:,3),2*NN)
call normalizebmet(bitmetrics(:,4),2*NN) call normalizebmet(bitmetrics(:,4),2*NN)
! Return the s4 array corresponding to N=1/hmod. Will be used for SNR calculation ! Return the s4 array corresponding to N=1. Will be used for SNR calculation
if(hmod.eq.1) s4hmod(:,:)=s4(:,:,1) s4snr(:,:)=s4(:,:,1)
if(hmod.eq.2) s4hmod(:,:)=s4(:,:,2)
if(hmod.eq.4) s4hmod(:,:)=s4(:,:,3)
if(hmod.eq.8) s4hmod(:,:)=s4(:,:,4)
return return
end subroutine get_fst4_bitmetrics2 end subroutine get_fst4_bitmetrics2

View File

@ -49,7 +49,7 @@ contains
complex, allocatable :: cframe(:) complex, allocatable :: cframe(:)
complex, allocatable :: c_bigfft(:) !Complex waveform complex, allocatable :: c_bigfft(:) !Complex waveform
real llr(240),llra(240),llrb(240),llrc(240),llrd(240) real llr(240),llra(240),llrb(240),llrc(240),llrd(240)
real candidates(100,4) real candidates(200,4)
real bitmetrics(320,4) real bitmetrics(320,4)
real s4(0:3,NN) real s4(0:3,NN)
real minsync real minsync
@ -253,6 +253,7 @@ contains
call four2a(c_bigfft,nfft1,1,-1,0) !r2c call four2a(c_bigfft,nfft1,1,-1,0) !r2c
! call blank2(nfa,nfb,nfft1,c_bigfft,iwave) ! call blank2(nfa,nfb,nfft1,c_bigfft,iwave)
nhicoh=0
if(hmod.eq.1) then if(hmod.eq.1) then
if(fMHz.lt.2.0) then if(fMHz.lt.2.0) then
nsyncoh=8 ! Use N=8 for sync nsyncoh=8 ! Use N=8 for sync
@ -277,7 +278,7 @@ contains
if(hmod.eq.1) then if(hmod.eq.1) then
if(ntrperiod.eq.15) minsync=1.15 if(ntrperiod.eq.15) minsync=1.15
if(ntrperiod.gt.15) minsync=1.20 if(ntrperiod.gt.15) minsync=1.25
elseif(hmod.gt.1) then elseif(hmod.gt.1) then
minsync=1.2 minsync=1.2
endif endif
@ -410,7 +411,7 @@ contains
ns4=count(hbits(229:244).eq.(/1,1,1,0,0,1,0,0,1,0,1,1,0,0,0,1/)) ns4=count(hbits(229:244).eq.(/1,1,1,0,0,1,0,0,1,0,1,1,0,0,0,1/))
ns5=count(hbits(305:320).eq.(/0,0,0,1,1,0,1,1,0,1,0,0,1,1,1,0/)) ns5=count(hbits(305:320).eq.(/0,0,0,1,1,0,1,1,0,1,0,0,1,1,1,0/))
nsync_qual=ns1+ns2+ns3+ns4+ns5 nsync_qual=ns1+ns2+ns3+ns4+ns5
! if(nsync_qual.lt. 46) cycle !### Value ?? ### if(nsync_qual.lt. 46) cycle !### Value ?? ###
scalefac=2.83 scalefac=2.83
llra( 1: 60)=bitmetrics( 17: 76, 1) llra( 1: 60)=bitmetrics( 17: 76, 1)
llra( 61:120)=bitmetrics( 93:152, 1) llra( 61:120)=bitmetrics( 93:152, 1)
@ -457,7 +458,7 @@ contains
if(itry.gt.nblock) then if(itry.gt.nblock) then
llr=llra llr=llra
if(nblock.gt.1) then if(nblock.gt.1) then
if(hmod.eq.1) llr=llrd if(hmod.eq.1) llr=llrc
if(hmod.eq.2) llr=llrb if(hmod.eq.2) llr=llrb
if(hmod.eq.4) llr=llrc if(hmod.eq.4) llr=llrc
if(hmod.eq.8) llr=llrd if(hmod.eq.8) llr=llrd
@ -550,7 +551,7 @@ contains
endif endif
xsig=0 xsig=0
do i=1,NN do i=1,NN
xsig=xsig+s4(itone(i),i)**2 xsig=xsig+s4(itone(i),i)
enddo enddo
arg=600.0*(xsig/base)-1.0 arg=600.0*(xsig/base)-1.0
if(arg.gt.0.0) then if(arg.gt.0.0) then
@ -737,7 +738,7 @@ contains
complex c_bigfft(0:nfft1/2) !Full length FFT of raw data complex c_bigfft(0:nfft1/2) !Full length FFT of raw data
integer hmod !Modulation index (submode) integer hmod !Modulation index (submode)
integer im(1) !For maxloc integer im(1) !For maxloc
real candidates(100,4) !Candidate list real candidates(200,4) !Candidate list
real, allocatable :: s(:) !Low resolution power spectrum real, allocatable :: s(:) !Low resolution power spectrum
real, allocatable :: s2(:) !CCF of s() with 4 tones real, allocatable :: s2(:) !CCF of s() with 4 tones
real xdb(-3:3) !Model 4-tone CCF peaks real xdb(-3:3) !Model 4-tone CCF peaks
@ -794,17 +795,18 @@ contains
! Find candidates, using the CLEAN algorithm to remove a model of each one ! Find candidates, using the CLEAN algorithm to remove a model of each one
! from s2() after it has been found. ! from s2() after it has been found.
pval=99.99 pval=99.99
do while(ncand.lt.100) do while(ncand.lt.200)
im=maxloc(s2(ia:ib)) im=maxloc(s2(ia:ib))
iploc=ia+im(1)-1 !Index of CCF peak iploc=ia+im(1)-1 !Index of CCF peak
pval=s2(iploc) !Peak value pval=s2(iploc) !Peak value
if(pval.lt.minsync) exit if(pval.lt.minsync) exit
do i=-3,+3 !Remove 0.9 of a model CCF at ! do i=-3,+3 !Remove 0.9 of a model CCF at
k=iploc+2*hmod*i !this frequency from s2() ! k=iploc+2*hmod*i !this frequency from s2()
if(k.ge.ia .and. k.le.ib) then ! if(k.ge.ia .and. k.le.ib) then
s2(k)=max(0.,s2(k)-0.9*pval*xdb(i)) ! s2(k)=max(0.,s2(k)-0.9*pval*xdb(i))
endif ! endif
enddo ! enddo
s2(max(1,iploc-2*hmod*3):min(nnw,iploc+2*hmod*3))=0.0
ncand=ncand+1 ncand=ncand+1
candidates(ncand,1)=df2*iploc !Candidate frequency candidates(ncand,1)=df2*iploc !Candidate frequency
candidates(ncand,2)=pval !Rough estimate of SNR candidates(ncand,2)=pval !Rough estimate of SNR