From 02a4cb95a4c608b3ee757df9598f28895cb44942 Mon Sep 17 00:00:00 2001 From: Steven Franke Date: Thu, 22 Dec 2016 17:21:23 +0000 Subject: [PATCH] Add an option on advanced time to enable MSK144 Rx filter equalizer. Other changes necessary to accommodate coexistence of static (Rx Filter) and dynamic (QSO partner-dependent) phase corrections. git-svn-id: svn+ssh://svn.code.sf.net/p/wsjt/wsjt/branches/wsjtx@7421 ab8295b8-cf94-4d9e-aec4-7959e3be5d79 --- CMakeLists.txt | 1 + Configuration.cpp | 6 +++ Configuration.hpp | 1 + Configuration.ui | 11 ++++++ lib/analytic.f90 | 79 ++++++++++++++++++++++++------------- lib/hspec.f90 | 9 +++-- lib/msk144signalquality.f90 | 19 ++++----- lib/mskrtd.f90 | 31 +++++++-------- mainwindow.cpp | 9 +++-- 9 files changed, 105 insertions(+), 61 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 571d9ed91..efa0e162d 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -421,6 +421,7 @@ set (wsjt_FSRCS lib/msk40_freq_search.f90 lib/msk144_freq_search.f90 lib/mskrtd.f90 + lib/msk144signalquality.f90 lib/msk144sim.f90 lib/mskrtd.f90 lib/options.f90 diff --git a/Configuration.cpp b/Configuration.cpp index 9750e3371..5674a1b81 100644 --- a/Configuration.cpp +++ b/Configuration.cpp @@ -540,6 +540,7 @@ private: bool twoPass_; bool x2ToneSpacing_; bool contestMode_; + bool rxEqualize_; bool realTimeDecode_; QString udp_server_name_; port_type udp_server_port_; @@ -628,6 +629,7 @@ bool Configuration::single_decode () const {return m_->single_decode_;} bool Configuration::twoPass() const {return m_->twoPass_;} bool Configuration::x2ToneSpacing() const {return m_->x2ToneSpacing_;} bool Configuration::contestMode() const {return m_->contestMode_;} +bool Configuration::rxEqualize() const {return m_->rxEqualize_;} bool Configuration::realTimeDecode() const {return m_->realTimeDecode_;} bool Configuration::split_mode () const {return m_->split_mode ();} QString Configuration::udp_server_name () const {return m_->udp_server_name_;} @@ -1057,6 +1059,7 @@ void Configuration::impl::initialize_models () ui_->cbTwoPass->setChecked(twoPass_); ui_->cbx2ToneSpacing->setChecked(x2ToneSpacing_); ui_->cbContestMode->setChecked(contestMode_); + ui_->cbRxEqualize->setChecked(rxEqualize_); ui_->cbRealTime->setChecked(realTimeDecode_); ui_->cbRealTime->setVisible(false); //Tempoary -- probably will remove this control ui_->type_2_msg_gen_combo_box->setCurrentIndex (type_2_msg_gen_); @@ -1285,6 +1288,7 @@ void Configuration::impl::read_settings () twoPass_ = settings_->value("TwoPass",true).toBool (); x2ToneSpacing_ = settings_->value("x2ToneSpacing",false).toBool (); contestMode_ = settings_->value("ContestMode",false).toBool (); + rxEqualize_ = settings_->value("RxEqualize",false).toBool (); realTimeDecode_ = settings_->value("RealTimeDecode",false).toBool (); rig_params_.poll_interval = settings_->value ("Polling", 0).toInt (); rig_params_.split_mode = settings_->value ("SplitMode", QVariant::fromValue (TransceiverFactory::split_mode_none)).value (); @@ -1384,6 +1388,7 @@ void Configuration::impl::write_settings () settings_->setValue ("TwoPass", twoPass_); settings_->setValue ("x2ToneSpacing", x2ToneSpacing_); settings_->setValue ("ContestMode", contestMode_); + settings_->setValue ("RxEqualize", rxEqualize_); settings_->setValue ("RealTimeDecode", realTimeDecode_); settings_->setValue ("UDPServer", udp_server_name_); settings_->setValue ("UDPServerPort", udp_server_port_); @@ -1779,6 +1784,7 @@ void Configuration::impl::accept () twoPass_ = ui_->cbTwoPass->isChecked (); x2ToneSpacing_ = ui_->cbx2ToneSpacing->isChecked (); contestMode_ = ui_->cbContestMode->isChecked (); + rxEqualize_ = ui_->cbRxEqualize->isChecked (); realTimeDecode_ = ui_->cbRealTime->isChecked (); frequency_calibration_intercept_ = ui_->calibration_intercept_spin_box->value (); frequency_calibration_slope_ppm_ = ui_->calibration_slope_ppm_spin_box->value (); diff --git a/Configuration.hpp b/Configuration.hpp index 22e4fa1dd..7fb1961b4 100644 --- a/Configuration.hpp +++ b/Configuration.hpp @@ -123,6 +123,7 @@ public: bool twoPass() const; bool x2ToneSpacing() const; bool contestMode() const; + bool rxEqualize() const; bool realTimeDecode() const; bool MyDx() const; bool CQMyN() const; diff --git a/Configuration.ui b/Configuration.ui index 505abe0d9..0e0ceea33 100644 --- a/Configuration.ui +++ b/Configuration.ui @@ -2295,6 +2295,16 @@ Right click for insert and delete options. + + + <html><head/><body><p>Correct for Rx filter group-delay variation.</p></body></html> + + + MSK144 Rx Equalization + + + + <html><head/><body><p>Generate Tx audio with twice the normal tone spacing. Intended for special LF/MF transmitters that use a divide-by-2 before generating RF.</p></body></html> @@ -2511,6 +2521,7 @@ soundcard changes sbBandwidth sbTxDelay cbContestMode + cbRxEqualize cbx2ToneSpacing cbRealTime diff --git a/lib/analytic.f90 b/lib/analytic.f90 index a4f5ef61f..dd9f51088 100644 --- a/lib/analytic.f90 +++ b/lib/analytic.f90 @@ -1,41 +1,42 @@ -subroutine analytic(d,npts,nfft,c,a,equalize) +subroutine analytic(d,npts,nfft,c,dpc,bseq,bdeq) ! Convert real data to analytic signal parameter (NFFTMAX=1024*1024) - real d(npts) - complex h(NFFTMAX/2) - real*8 a(5),alast(5),fp - real ac(5),a0(5) - complex c(NFFTMAX) - logical*1 equalize - data nfft0/0/ - data alast/0.0,0.0,0.0,0.0,0.0/ - data ac/1.0,0.05532,0.11438,0.12918,0.09274/ ! amp coeffs for TS2000 - data a0/0.0,0.0,-0.952,0.768,-0.565/ ! baseline phase coeffs for TS2000 - save nfft0,h,alast,a0,ac -! disable baseline phase correction for commit - should look for a file with coeffs - a0(1:5)=0.0 + real d(npts) ! passband signal + real h(NFFTMAX/2) ! real BPF magnitude + real dpc(3),dpclast(3) ! dynamic phase coeffs + real spc(3),spclast(3) ! static phase coeffs + real ac(5) ! currently unused static amp coeffs + real fp + + complex corrs(NFFTMAX/2) ! allpass static phase correction + complex corrd(NFFTMAX/2) ! allpass overall phase correction + complex c(NFFTMAX) ! analytic signal + + logical*1 bseq ! boolean static equalizer flag + logical*1 bdeq ! boolean dynamic equalizer flag + logical*1 bseqlast + + data nfft0/0/ + data bseqlast/.false./ + data spclast/0.0,0.0,0.0/ + data spc/-0.952,0.768,-0.565/ ! baseline phase coeffs for TS2000 + data ac/1.0,0.05532,0.11438,0.12918,0.09274/ ! amp coeffs for TS2000 + + save nfft0,h,spc,spclast,dpclast,ac,pi,t,beta df=12000.0/nfft nh=nfft/2 - if( nfft.ne.nfft0 .or. any(alast .ne. a) ) then + if( nfft.ne.nfft0 ) then + pi=4.0*atan(1.0) t=1.0/2000.0 beta=0.1 - pi=4.0*atan(1.0) do i=1,nh+1 ff=(i-1)*df f=ff-1500.0 - fp=f/1000.0 - h(i)=cmplx(1.0,0.0) - if( equalize ) then - phase0=a0(1)+fp*(a0(2)+fp*(a0(3)+fp*(a0(4)+fp*a0(5)))) - phase=a(1)+fp*(a(2)+fp*(a(3)+fp*(a(4)+fp*a(5)))) -! amp=ac(1)+fp*(ac(2)+fp*(ac(3)+fp*(ac(4)+fp*ac(5)))) - amp=1.0 ! no amplitude correction for now - h(i)=amp*cmplx(cos(phase),sin(phase))*cmplx(cos(phase0),sin(phase0)) - endif + h(i)=1.0 if(abs(f).gt.(1-beta)/(2*t) .and. abs(f).le.(1+beta)/(2*t)) then h(i)=h(i)*0.5*(1+cos((pi*t/beta )*(abs(f)-(1-beta)/(2*t)))) endif @@ -43,15 +44,39 @@ subroutine analytic(d,npts,nfft,c,a,equalize) nfft0=nfft endif + if( any(spclast .ne. spc) .or. any(dpclast .ne. dpc) ) then + spclast=spc + dpclast=dpc + do i=1,nh+1 + ff=(i-1)*df + f=ff-1500.0 + fp=f/1000.0 + ps=fp*fp*(spc(1)+fp*(spc(2)+fp*spc(3))) + corrs(i)=cmplx(cos(ps),sin(ps)) + pd=fp*fp*(dpc(1)+fp*(dpc(2)+fp*dpc(3))) + corrd(i)=cmplx(cos(pd),sin(pd)) + enddo + endif + fac=2.0/nfft c(1:npts)=fac*d(1:npts) c(npts+1:nfft)=0. call four2a(c,nfft,1,-1,1) !Forward c2c FFT - c(1:nh+1)=h(1:nh+1)*c(1:nh+1) + if( (.not. bseq) .and. (.not. bdeq) ) then + c(1:nh+1)=h(1:nh+1)*c(1:nh+1) + else if( bseq .and. (.not. bdeq) ) then + c(1:nh+1)=h(1:nh+1)*corrs(1:nh+1)*c(1:nh+1) + else if( (.not. bseq) .and. bdeq ) then + c(1:nh+1)=h(1:nh+1)*corrd(1:nh+1)*c(1:nh+1) + else if( bseq .and. bdeq ) then + c(1:nh+1)=h(1:nh+1)*corrs(1:nh+1)*corrd(1:nh+1)*c(1:nh+1) + endif + c(1)=0.5*c(1) !Half of DC term c(nh+2:nfft)=0. !Zero the negative frequencies call four2a(c,nfft,1,1,1) !Inverse c2c FFT - alast=a + spclast=spc + dpclast=dpc return end subroutine analytic diff --git a/lib/hspec.f90 b/lib/hspec.f90 index 9070c9487..9bf2c07ea 100644 --- a/lib/hspec.f90 +++ b/lib/hspec.f90 @@ -1,5 +1,5 @@ -subroutine hspec(id2,k,nutc0,ntrpdepth,nrxfreq,ntol,bmsk144,bcontest, & - ingain,mycall,hiscall,bshmsg,green,s,jh,line1,mygrid) +subroutine hspec(id2,k,nutc0,ntrpdepth,nrxfreq,ntol,bmsk144,bcontest, & + brxequal,ingain,mycall,hiscall,bshmsg,green,s,jh,line1,mygrid) ! Input: ! k pointer to the most recent new data @@ -8,6 +8,7 @@ subroutine hspec(id2,k,nutc0,ntrpdepth,nrxfreq,ntol,bmsk144,bcontest, & ! nrxfreq Rx audio center frequency ! ntol Decoding range is +/- ntol ! bmsk144 Boolean, true if in MSK144 mode +! brxequal Boolean, turns on equalization in MSK144 mode ! ingain Relative gain for spectra ! Output: @@ -20,7 +21,7 @@ subroutine hspec(id2,k,nutc0,ntrpdepth,nrxfreq,ntol,bmsk144,bcontest, & character*12 mycall,hiscall character*6 mygrid integer*2 id2(0:120*12000-1) - logical*1 bmsk144,bcontest,bshmsg + logical*1 bmsk144,bcontest,bshmsg,brxequal real green(0:JZ-1) real s(0:63,0:JZ-1) real x(512) @@ -83,7 +84,7 @@ subroutine hspec(id2,k,nutc0,ntrpdepth,nrxfreq,ntol,bmsk144,bcontest, & tt2=sum(float(abs(id2(k0:k0+3583)))) if(tt1.ne.0.0 .and. tt2.ne.0) then call mskrtd(id2(k-7168+1:k),nutc0,tsec,ntol,nrxfreq,ndepth, & - mycall,mygrid,hiscall,bshmsg,bcontest,line1) + mycall,mygrid,hiscall,bshmsg,bcontest,brxequal,line1) endif endif endif diff --git a/lib/msk144signalquality.f90 b/lib/msk144signalquality.f90 index 0a26ee1eb..760dd3be9 100644 --- a/lib/msk144signalquality.f90 +++ b/lib/msk144signalquality.f90 @@ -1,4 +1,5 @@ -subroutine msk144signalquality(cframe,snr,freq,t0,softbits,msg,dxcall,nbiterrors,eyeopening,trained,pcoeffs) + subroutine msk144signalquality(cframe,snr,freq,t0,softbits,msg,dxcall, & + nbiterrors,eyeopening,trained,pcoeffs) character*22 msg,msgsent character*12 dxcall @@ -35,7 +36,8 @@ subroutine msk144signalquality(cframe,snr,freq,t0,softbits,msg,dxcall,nbiterrors real d(1024) real phase(864) real twopi,freq,phi,dphi0,dphi1,dphi - real*8 x(145),y(145),pp(145),sigmay(145),a(5),chisqr,pcoeffs(5) + real*8 x(145),y(145),pp(145),sigmay(145),a(5),chisqr + real pcoeffs(3) data first/.true./ save cross_avg,abscross_avg,wt_avg,first,currently_training, & @@ -51,7 +53,7 @@ subroutine msk144signalquality(cframe,snr,freq,t0,softbits,msg,dxcall,nbiterrors training_dxcall(1:12)=' ' trained=.false. currently_training=.false. - pcoeffs(1:5)=0.0 + pcoeffs(1:3)=0.0 first=.false. endif @@ -68,7 +70,7 @@ subroutine msk144signalquality(cframe,snr,freq,t0,softbits,msg,dxcall,nbiterrors currently_training=.false. training_dxcall(1:12)=' ' trained_dxcall(1:12)=' ' - pcoeffs(1:5)=0.0 + pcoeffs(1:3)=0.0 write(*,*) 'reset to untrained state ' endif @@ -79,7 +81,7 @@ write(*,*) 'reset to untrained state ' currently_training=.true. training_dxcall=trim(dxcall) trained_dxcall(1:12)=' ' - pcoeffs(1:5)=0.0 + pcoeffs(1:3)=0.0 write(*,*) 'start training on call ',training_dxcall endif @@ -87,7 +89,7 @@ write(*,*) 'start training on call ',training_dxcall trained=.false. ! just to be sure trained_dxcall(1:12)=' ' training_dxcall=dxcall - pcoeffs(1:5)=0.0 + pcoeffs(1:3)=0.0 endif ! use decoded message to figure out how many bit errors in the frame @@ -167,8 +169,7 @@ write(*,*) 'start training on call ',training_dxcall nfft=1024 d=0 d(1:864)=waveform(0:863) - a=0 - call analytic(d,npts,nfft,canalytic,a,.false.) ! don't equalize the model + call analytic(d,npts,nfft,canalytic,pcoeffs,.false.,.false.) ! don't equalize the model call tweak1(canalytic,nfft,-freq,cmodel) call four2a(cframe(1:864),864,1,-1,1) call four2a(cmodel(1:864),864,1,-1,1) @@ -207,7 +208,7 @@ write(*,*) 'training ',navg,sqrt(chisqr),rmsdiff write(19,*) i,real(cframe(i)),imag(cframe(i)) enddo close(19) - pcoeffs=a + pcoeffs=a(3:5) trained_dxcall=dxcall training_dxcall(1:12)=' ' currently_training=.false. diff --git a/lib/mskrtd.f90 b/lib/mskrtd.f90 index 05cd40d69..d153ac3a2 100644 --- a/lib/mskrtd.f90 +++ b/lib/mskrtd.f90 @@ -1,5 +1,5 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & - bshmsg,bcontest,line) + bshmsg,bcontest,brxequal,line) ! Real-time decoder for MSK144. ! Analysis block size = NZ = 7168 samples, t_block = 0.597333 s @@ -30,11 +30,11 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & real pow(8) real softbits(144) real xmc(NPATTERNS) - real*8 pcoeffs(5) + real pcoeffs(3) - logical*1 bshmsg,bcontest + logical*1 bshmsg,bcontest,brxequal logical first - logical*1 equalized + logical*1 trained data first/.true./ data iavpatterns/ & @@ -43,14 +43,14 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & 1,1,1,1,1,0,0,0, & 1,1,1,1,1,1,1,0/ data xmc/2.0,4.5,2.5,3.5/ !Used to set time at center of averaging mask - save first,tsec0,nutc00,pnoise,nsnrlast,msglast,cdat,pcoeffs,equalized + save first,tsec0,nutc00,pnoise,nsnrlast,msglast,cdat,pcoeffs,trained if(first) then tsec0=tsec nutc00=nutc0 pnoise=-1.0 - pcoeffs(1:5)=0.0 - equalized=.false. + pcoeffs(1:3)=0.0 + trained=.false. first=.false. endif @@ -74,7 +74,7 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & fac=1.0/rms d(1:NZ)=fac*d(1:NZ) d(NZ+1:NFFT1)=0. - call analytic(d,NZ,NFFT1,cdat,pcoeffs,.true.) !Convert to analytic signal and filter + call analytic(d,NZ,NFFT1,cdat,pcoeffs,brxequal,.true.) !Convert to analytic signal and filter ! Calculate average power for each frame and for the entire block. ! If decode is successful, largest power will be taken as signal+noise. @@ -101,8 +101,6 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & if( nsuccess .eq. 1 ) then tdec=tsec+tdec - decsym=' & ' - if( equalized ) decsym=' ^ ' ipk=0 is=0 goto 900 @@ -125,7 +123,7 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & if( nsyncsuccess .eq. 0 ) cycle do ipk=1,npeaks - do is=1,3 ! With equalization, this loop may not be necessary + do is=1,3 ic0=npkloc(ipk) if(is.eq.2) ic0=max(1,ic0-1) if(is.eq.3) ic0=min(NSPM,ic0+1) @@ -133,9 +131,6 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & call msk144decodeframe(ct,softbits,msgreceived,ndecodesuccess) if(ndecodesuccess .gt. 0) then tdec=tsec+xmc(iavg)*tframe - decsym=' & ' - if( equalized ) decsym=' ^ ' - if( equalized .and. is .ne. 1 ) decsym=' ! ' !help decide if is dither is needed goto 900 endif enddo !Slicer dither @@ -165,7 +160,7 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & nsnr=nint(snr0) call msk144signalquality(ct,snr0,fest,tdec,softbits,msgreceived,hiscall, & - ncorrected,eyeopening,equalized,pcoeffs) + ncorrected,eyeopening,trained,pcoeffs) ! Dupe check. Only print if new message, or higher snr. if(msgreceived.ne.msglast .or. nsnr.gt.nsnrlast .or. tsec.lt.tsec0) then @@ -177,6 +172,10 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & if(msgreceived(1:1).ne.'<') then call fix_contest_msg(mycall(1:6),mygrid,hiscall(1:6),msgreceived) endif + decsym=' & ' + if( brxequal .and. (.not. trained) ) decsym=' ^ ' + if( brxequal .and. trained ) decsym=' $ ' + if( (.not. brxequal) .and. trained ) decsym=' @ ' write(line,1020) nutc0,nsnr,tdec,nint(fest),decsym,msgreceived, & navg,ncorrected,eyeopening,char(0) 1020 format(i6.6,i4,f5.1,i5,a3,a22,i2,i3,f5.1,a1) @@ -187,5 +186,3 @@ subroutine mskrtd(id2,nutc0,tsec,ntol,nrxfreq,ndepth,mycall,mygrid,hiscall, & end subroutine mskrtd include 'fix_contest_msg.f90' - -include 'msk144signalquality.f90' diff --git a/mainwindow.cpp b/mainwindow.cpp index 01224f8b6..d5718d878 100644 --- a/mainwindow.cpp +++ b/mainwindow.cpp @@ -67,9 +67,9 @@ extern "C" { int* minw, float* px, float s[], float* df3, int* nhsym, int* npts8); void hspec_(short int d2[], int* k, int* nutc0, int* ntrperiod, int* nrxfreq, int* ntol, - bool* bmsk144, bool* bcontest, int* ingain, char mycall[], char hiscall[], - bool* bshmsg, float green[], float s[], int* jh, char line[], - char mygrid[], int len1, int len2, int len3, int len4); + bool* bmsk144, bool* bcontest, bool* brxequalize, int* ingain, char mycall[], + char hiscall[], bool* bshmsg, float green[], float s[], int* jh, + char line[], char mygrid[], int len1, int len2, int len3, int len4); void gen4_(char* msg, int* ichk, char* msgsent, int itone[], int* itext, int len1, int len2); @@ -1284,9 +1284,10 @@ void MainWindow::fastSink(qint64 frames) QString hisCall {ui->dxCallEntry->text ()}; bool bshmsg=ui->cbShMsgs->isChecked(); bool bcontest=m_config.contestMode(); + bool brxequalize=m_config.rxEqualize(); strncpy(dec_data.params.hiscall,(hisCall + " ").toLatin1 ().constData (), 12); strncpy(dec_data.params.mygrid, (m_config.my_grid()+" ").toLatin1(),6); - hspec_(dec_data.d2,&k,&nutc0,&nTRpDepth,&m_RxFreq,&m_Ftol,&bmsk144,&bcontest, + hspec_(dec_data.d2,&k,&nutc0,&nTRpDepth,&m_RxFreq,&m_Ftol,&bmsk144,&bcontest,&brxequalize, &m_inGain,&dec_data.params.mycall[0],&dec_data.params.hiscall[0],&bshmsg, fast_green,fast_s,&fast_jh,&line[0],&dec_data.params.mygrid[0],12,12,80,6); float px = fast_green[fast_jh];