ft8d/baseline.f90

50 lines
1.2 KiB
Fortran
Raw Normal View History

2018-03-24 07:20:34 -04:00
subroutine baseline(s,nfa,nfb,sbase)
! Fit baseline to spectrum (for FT8)
! Input: s(npts) Linear scale in power
! Output: sbase(npts) Baseline
include 'ft8_params.f90'
2018-03-24 07:20:34 -04:00
implicit real*8 (a-h,o-z)
2018-04-01 11:36:02 -04:00
real*4 s(NFFT1)
real*4 sbase(NFFT1)
2018-03-24 07:20:34 -04:00
real*4 base
real*8 x(1000),y(1000),a(5)
data nseg/10/,npct/10/
2018-04-01 16:08:21 -04:00
df=4000.0/NFFT1 !3.125 Hz
2018-03-24 07:20:34 -04:00
ia=max(1,nint(nfa/df))
ib=nint(nfb/df)
do i=ia,ib
s(i)=10.0*log10(s(i)) !Convert to dB scale
enddo
nterms=5
nlen=(ib-ia+1)/nseg !Length of test segment
i0=(ib-ia+1)/2 !Midpoint
k=0
do n=1,nseg !Loop over all segments
ja=ia + (n-1)*nlen
jb=ja+nlen-1
call pctile(s(ja),nlen,npct,base) !Find lowest npct of points
do i=ja,jb
if(s(i).le.base) then
if (k.lt.1000) k=k+1 !Save all "lower envelope" points
x(k)=i-i0
y(k)=s(i)
endif
enddo
enddo
kz=k
a=0.
call polyfit(x,y,y,kz,nterms,0,a,chisqr) !Fit a low-order polynomial
do i=ia,ib
t=i-i0
sbase(i)=a(1)+t*(a(2)+t*(a(3)+t*(a(4)+t*(a(5))))) + 0.65
! write(51,3051) i*df,s(i),sbase(i)
!3051 format(3f12.3)
enddo
return
end subroutine baseline