subroutine syncmsk144(cdat,npts,metric,msgreceived,fest) use iso_c_binding, only: c_loc,c_size_t use packjt use hashing use timer_module, only: timer parameter (NSPM=864) character*22 msgreceived character*85 pchk_file,gen_file complex cdat(npts) !Analytic signal complex cdat2(npts) !Analytic signal complex cav(NSPM) !Analytic signal complex c(NSPM) complex ctmp(6000) complex cb(42) !Complex waveform for sync word complex csum,cfac,cca,ccb complex cc(npts) complex cc1(npts) complex cc2(npts) complex bb(6) integer s8(8),hardbits(144),hardword(128),unscrambledhardbits(128) integer*1, target:: i1Dec8BitBytes(10) integer, dimension(1) :: iloc integer*4 i4Msg6BitWords(12) !72-bit message as 6-bit words integer*4 i4Dec6BitWords(12) integer*1 decoded(80) integer*1, allocatable :: message(:) integer*1 i1hashdec integer ipeaks(10) logical ismask(6000) real cbi(42),cbq(42) real tonespec(6000) real rcw(12) real dd(npts) real pp(12) !Half-sine pulse shape real*8 dt, df, fs, twopi real softbits(144) real*8 unscrambledsoftbits(128) real lratio(128) logical first data first/.true./ data s8/0,1,1,1,0,0,1,0/ ! save first,cb,cd,pi,twopi,dt,f0,f1 save if(first) then ! These files can be found in /lib/ldpc/jtmode_codes directory pchk_file="peg-128-80-reg3.pchk" gen_file="peg-128-80-reg3.gen" call init_ldpc(trim(pchk_file)//char(0),trim(gen_file)//char(0)) ! define half-sine pulse and raised-cosine edge window pi=4d0*datan(1d0) twopi=8d0*datan(1d0) dt=1.0/12000.0 do i=1,12 angle=(i-1)*pi/12.0 pp(i)=sin(angle) rcw(i)=(1-cos(angle))/2 enddo ! define the sync word waveform s8=2*s8-1 cbq(1:6)=pp(7:12)*s8(1) cbq(7:18)=pp*s8(3) cbq(19:30)=pp*s8(5) cbq(31:42)=pp*s8(7) cbi(1:12)=pp*s8(2) cbi(13:24)=pp*s8(4) cbi(25:36)=pp*s8(6) cbi(37:42)=pp(1:6)*s8(8) cb=cmplx(cbi,cbq) first=.false. endif ! Coarse carrier frequency sync ! look for tones near 2k and 4k in the (analytic signal)**2 spectrum ! search range for coarse frequency error is +/- 100 Hz fs=12000.0 nfft=6000 !using a zero-padded fft to get 2 Hz bins df=fs/nfft ctmp=cmplx(0.0,0.0) ctmp(1:npts)=cdat**2 ctmp(1:12)=ctmp(1:12)*rcw ctmp(npts-11:npts)=ctmp(npts-11:npts)*rcw(12:1:-1) call four2a(ctmp,nfft,1,-1,1) tonespec=abs(ctmp)**2 ismask=.false. ismask(1901:2101)=.true. ! high tone search window iloc=maxloc(tonespec,ismask) ihpk=iloc(1) ah=tonespec(ihpk) ismask=.false. ismask(901:1101)=.true. ! window for low tone iloc=maxloc(tonespec,ismask) ilpk=iloc(1) al=tonespec(ilpk) if( ah .ge. al ) then ferr=(ihpk-2001)*df/2.0 tot=sum(tonespec(1901:2101)) q1=200*ah/(tot-tonespec(ihpk)) else ferr=(ilpk-1001)*df/2.0 tot=sum(tonespec(901:1101)) q1=200*al/(tot-tonespec(ilpk)) endif fdiff=(ihpk-ilpk)*df ! remove coarse freq error - should now be within a few Hz call tweak1(cdat,npts,-(1500+ferr),cdat) ! attempt frame synchronization ! correlate with sync word waveforms - the resulting complex ! correlations provide all synch information. cc=0 cc1=0 cc2=0 do i=1,npts-(56*6+41) cc1(i)=sum(cdat(i:i+41)*conjg(cb)) cc2(i)=sum(cdat(i+56*6:i+56*6+41)*conjg(cb)) enddo cc=cc1+cc2 dd=abs(cc1)*abs(cc2) ! Find 5 largest peaks do ipk=1,5 iloc=maxloc(abs(cc)) ic1=iloc(1) iloc=maxloc(dd) ic2=iloc(1) ipeaks(ipk)=ic2 dd(max(1,ic2-7):min(npts-56*6-41,ic2+7))=0.0 enddo ! See if we can find "closed brackets" - a pair of peaks that differ by 864, plus or minus ! This information is not yet used for anything do ii=1,5 do jj=ii+1,5 if( (ii .ne. jj) .and. (abs( abs(ipeaks(ii)-ipeaks(jj))-864) .le. 5) ) then write(78,*) "closed brackets: ",ii,jj,ipeaks(ii),ipeaks(jj),abs(ipeaks(ii)-ipeaks(jj)) endif enddo enddo do iav=0,1 do ipk=1,5 ! we want ic to be the index of the first sample of the message ic=ipeaks(ipk) ! This needs to be improved - it's used to protect the edges of the array from ! overruns. if( ic .lt. 12 .or. ic .gt. 2*864-12 ) then ! write(*,*) "Peak not in central section: ",ipk,is,ic cycle endif ! bb is used to place the sampling index at the center of the eye do i=1,6 io=i-3 bb(i) = sum( ( cdat(ic+io:ic+io+864:6) * conjg( cdat(ic+io+6:ic+io+6+864:6) ) )*2 ) enddo iloc=maxloc(abs(bb)) ibb=iloc(1) ! write(*,*) 'ic0: ',ic,'bb peak is at : ',ibb ! Adjust frame index to place peak of bb at desired lag ic=ic + ibb-2 ! Sanity check - recompute bb and verify that peak is now at designated lag. ! do i=1,6 ! io=i-3 ! bb(i) = sum( ( cdat(ic+io:ic+io+864:6) * conjg( cdat(ic+io+6:ic+io+6+864:6) ) )*2 ) ! enddo ! iloc=maxloc(abs(bb)) ! ibb=iloc(1) ! write(*,*) 'ic1: ',ic,'bb peak is at : ',ibb ! Average two frames on the second pass only if its incredibly easy to do ! better than nothing. Should be improved. c=cdat(ic:ic+864-1) if( iav .eq. 1 ) then id0=ic+864 id1=ic+864+863 if( id1 .le. npts ) then c=c+cdat(id0:id1) else cycle endif endif ! Estimate fine frequency error. cca=sum(c(1:1+41)*conjg(cb)) ccb=sum(c(1+56*6:1+56*6+41)*conjg(cb)) cfac=ccb*conjg(cca) ferr2=atan2(imag(cfac),real(cfac))/(twopi*56*6*dt) ! Final estimate of the carrier frequency - returned to the calling program fest=1500+ferr+ferr2 ! Remove fine frequency error and put the results in cdat2 call tweak1(c,npts,-ferr2,c) ! Estimate final frequency error and carrier phase. cca=sum(c(1:1+41)*conjg(cb)) ccb=sum(c(1+56*6:1+56*6+41)*conjg(cb)) cfac=ccb*conjg(cca) ffin=atan2(imag(cfac),real(cfac))/(twopi*56*6*dt) phase0=atan2(imag(cca+ccb),real(cca+ccb)) ! Remove the static phase error from the data cfac=cmplx(cos(phase0),sin(phase0)) c=c*conjg(cfac) ! sample to get softsamples do i=1,72 softbits(2*i-1)=imag(c(1+(i-1)*12)) softbits(2*i)=real(c(7+(i-1)*12)) enddo hardbits=0 do i=1,144 if( softbits(i) .ge. 0.0 ) then hardbits(i)=1 endif enddo ! calculated the number of sync-word bits that are incorrect ! this might come in handy some day nbadsync=sum(s8*(2*hardbits(1:8)-1)) nbadsync=nbadsync+sum(s8*(2*hardbits(57:57+7)-1)) nbadsync=16-nbadsync ! this could be used to count the number of hard errors that were corrected hardword(1:48)=hardbits(9:9+47) hardword(49:128)=hardbits(65:65+80-1) unscrambledhardbits(1:127:2)=hardword(1:64) unscrambledhardbits(2:128:2)=hardword(65:128) ! normalize the softsymbols before submitting to decoder sav=sum(softbits)/144 s2av=sum(softbits*softbits)/144 ssig=sqrt(s2av-sav*sav) softbits=softbits/ssig sigma=0.65 lratio(1:48)=softbits(9:9+47) lratio(49:128)=softbits(65:65+80-1) lratio=exp(2.0*lratio/(sigma*sigma)) unscrambledsoftbits(1:127:2)=lratio(1:64) unscrambledsoftbits(2:128:2)=lratio(65:128) max_iterations=20 max_dither=50 call ldpc_decode(unscrambledsoftbits, decoded, max_iterations, niterations, max_dither, ndither) ! if( niterations .lt. 0 ) then ! msgreceived=' ' ! return ! endif if( niterations .ge. 0.0 ) then goto 778 endif enddo enddo msgreceived=' ' goto 999 778 continue ! The decoder found a codeword - compare decoded hash with calculated ! Collapse 80 decoded bits to 10 bytes. Bytes 1-9 are the message, byte 10 is the hash do ibyte=1,10 itmp=0 do ibit=1,8 itmp=ishft(itmp,1)+iand(1,decoded((ibyte-1)*8+ibit)) enddo i1Dec8BitBytes(ibyte)=itmp enddo ! Calculate the hash using the first 9 bytes. ihashdec=nhash(c_loc(i1Dec8BitBytes),int(9,c_size_t),146) ihashdec=2*iand(ihashdec,255) ! Compare calculated hash with received byte 10 - if they agree, keep the message. i1hashdec=ihashdec if( i1hashdec .eq. i1Dec8BitBytes(10) ) then ! Good hash --- unpack 72-bit message do ibyte=1,12 itmp=0 do ibit=1,6 itmp=ishft(itmp,1)+iand(1,decoded((ibyte-1)*6+ibit)) enddo i4Dec6BitWords(ibyte)=itmp enddo call unpackmsg(i4Dec6BitWords,msgreceived) endif write(78,1001) iav,ipk,is,fest,nbadsync,phase0,niterations,ndither,i1hashdec,i1Dec8BitBytes(10),msgreceived 1001 format(i6,i6,i6,f10.1,i6,f10.2,i6,i6,i6,i6,4x,a22) 999 return end subroutine syncmsk144