diff --git a/datcom.f90 b/datcom.f90 index 06cb1da86..c278c6389 100644 --- a/datcom.f90 +++ b/datcom.f90 @@ -1,5 +1,5 @@ parameter (NSMAX=60*96000) !Samples per 60 s file -integer*2 id !46 MB: raw data from Linrad timf2 +real*4 dd !92 MB: raw data from Linrad timf2 character*80 fname80 -common/datcom/id(4,NSMAX,2),nutc,newdat2,kbuf,kxp,kk,kkdone,nlost, & +common/datcom/dd(4,NSMAX,2),nutc,newdat2,kbuf,kxp,kk,kkdone,nlost, & nlen,fname80 diff --git a/decode1.F90 b/decode1.F90 index 041346f3e..b0768ff3b 100644 --- a/decode1.F90 +++ b/decode1.F90 @@ -41,7 +41,7 @@ subroutine decode1(iarg) n=Tsec if((ndiskdat.eq.1 .or. ndecoding.eq.0) .and. ((kkk-kkdone).gt.32768)) then - call symspec(id,kbuf,kk,kkdone,nutc,newdat) + call symspec(dd,kbuf,kk,kkdone,nutc,newdat) call sleep_msec(10) endif diff --git a/decode1a.f b/decode1a.f index 82a1c8c2d..1a1908f8b 100644 --- a/decode1a.f +++ b/decode1a.f @@ -1,4 +1,4 @@ - subroutine decode1a(id,newdat,freq,nflip, + subroutine decode1a(dd,newdat,freq,nflip, + mycall,hiscall,hisgrid,neme,ndepth,nqd,dphi,ndphi, + ipol,sync2,a,dt,pol,nkv,nhist,qual,decoded) @@ -7,7 +7,7 @@ C to decode it. parameter (NFFT1=77760,NFFT2=2430) parameter (NMAX=60*96000) !Samples per 60 s - integer*2 id(4,NMAX) !46 MB: raw data from Linrad timf2 + real*4 dd(4,NMAX) !92 MB: raw data from Linrad timf2 complex c2x(NMAX/4), c2y(NMAX/4) !After 1/4 filter and downsample complex c3x(NMAX/16),c3y(NMAX/16) !After 1/16 filter and downsample complex c4x(NMAX/64),c4y(NMAX/64) !After 1/64 filter and downsample @@ -29,7 +29,7 @@ C Mix sync tone to baseband, low-pass filter, and decimate by 64 dt00=dt C If freq=125.0 kHz, f0=48000 Hz. f0=1000*(freq-77.0) !Freq of sync tone (0-96000 Hz) - call filbig(id,NMAX,f0,newdat,cx,cy,n5) + call filbig(dd,NMAX,f0,newdat,cx,cy,n5) joff=0 sqa=0. sqb=0. diff --git a/filbig.f b/filbig.f index da1fd570e..5773ccd55 100644 --- a/filbig.f +++ b/filbig.f @@ -1,12 +1,12 @@ - subroutine filbig(id,nmax,f0,newdat,c4a,c4b,n4) + subroutine filbig(dd,nmax,f0,newdat,c4a,c4b,n4) C Filter and downsample complex data for X and Y polarizations, -C stored in array id(4,nmax). Output is downsampled from 96000 Hz +C stored in array dd(4,nmax). Output is downsampled from 96000 Hz C to 1500 Hz, and the low-pass filter has f_cutoff = 375 Hz and C f_stop = 750 Hz. parameter (NFFT1=5376000,NFFT2=77175) - integer*2 id(4,nmax) !Input data + real*4 dd(4,nmax) !Input data complex c4a(NFFT2),c4b(NFFT2) !Output data complex ca(NFFT1),cb(NFFT1) !FFTs of input real*8 df @@ -71,8 +71,8 @@ C If we just have a new f0, continue with the existing ca and cb. if(newdat.ne.0) then nz=min(nmax,NFFT1) do i=1,nz - ca(i)=cmplx(float(int(id(1,i))),float(int(id(2,i)))) - cb(i)=cmplx(float(int(id(3,i))),float(int(id(4,i)))) + ca(i)=cmplx(dd(1,i),dd(2,i)) + cb(i)=cmplx(dd(3,i),dd(4,i)) enddo if(nmax.lt.NFFT1) then diff --git a/getfile2.F90 b/getfile2.F90 index ecf17c7e4..fbd84f2c7 100644 --- a/getfile2.F90 +++ b/getfile2.F90 @@ -11,6 +11,7 @@ subroutine getfile2(fname,len) include 'gcom1.f90' include 'gcom2.f90' include 'gcom4.f90' + integer*2 id(4,NSMAX) 1 if(ndecoding.eq.0) go to 2 #ifdef CVF @@ -34,7 +35,15 @@ subroutine getfile2(fname,len) kbuf=1 call cs_lock('getfile2a') +!### +! NB: not really necessary to read whole file at once. Save memory! call rfile3a(fname,id,n,ierr) + do i=1,NSMAX + dd(1,i,1)=id(1,i) + dd(2,i,1)=id(2,i) + enddo +!### + call cs_unlock if(ierr.ne.0) then print*,'Error opening or reading file: ',fname,ierr @@ -45,8 +54,7 @@ subroutine getfile2(fname,len) ka=0.1*NSMAX kb=0.8*NSMAX do k=ka,kb - sq=sq + float(int(id(1,k,1)))**2 + float(int(id(2,k,1)))**2 + & - float(int(id(3,k,1)))**2 + float(int(id(4,k,1)))**2 + sq=sq + dd(1,k,1)**2 + dd(2,k,1)**2 + dd(3,k,1)**2 + dd(4,k,1)**2 enddo sqave=174*sq/(kb-ka+1) rxnoise=10.0*log10(sqave) - 48.0 diff --git a/map65.py b/map65.py index 72b3a93ee..ce3a6df4f 100644 --- a/map65.py +++ b/map65.py @@ -1,4 +1,4 @@ -#--------------------------------------------------------------------- MAP65 +#-------------------------------------------------------------------- MAP65 # $Date$ $Revision$ # from Tkinter import * diff --git a/map65a.F90 b/map65a.F90 index 2cf1526e3..44cb2dff5 100644 --- a/map65a.F90 +++ b/map65a.F90 @@ -201,7 +201,7 @@ subroutine map65a(newdat) nkm.eq.1) km=km-1 if(freq-freq0.gt.ftol .or. sync1.gt.sync10) then nflip=nint(flipk) - call decode1a(id(1,1,kbuf),newdat,freq,nflip, & + call decode1a(dd(1,1,kbuf),newdat,freq,nflip, & mycall,hiscall,hisgrid,neme,ndepth,nqd,dphi,ndphi, & ipol,sync2,a,dt,pol,nkv,nhist,qual,decoded) @@ -395,8 +395,10 @@ subroutine map65a(newdat) call display(nkeep,ncsmin) ndecdone=2 - if(nsave.gt.0 .and. ndiskdat.eq.0) call savetf2(id(1,1,kbuf), & - fnamedate,savedir) +!### Temporarily disable the optional saving of raw data +! if(nsave.gt.0 .and. ndiskdat.eq.0) call savetf2(id(1,1,kbuf), & +! fnamedate,savedir) +!### 999 close(23) ndphi=0 diff --git a/recvpkt.F90 b/recvpkt.F90 index 64b1b5100..dcc337e00 100644 --- a/recvpkt.F90 +++ b/recvpkt.F90 @@ -1,19 +1,23 @@ subroutine recvpkt(iarg) -! Receive timf2 packets from Linrad and stuff data into array id(). +! Receive timf2 packets from Linrad and stuff data into array dd(). ! (This routine runs in a background thread and will never return.) parameter (NSZ=2*60*96000) - real*8 d8(NSZ) integer*1 userx_no,iusb integer*2 nblock,nblock0 logical first,synced - real*8 center_freq,buf8 + real*8 center_freq,d8,buf8 + complex*16 c16,buf16(87) + integer*2 jd(4) + real*4 xd(4) common/plrscom/center_freq,msec,fqso,iptr,nblock,userx_no,iusb,buf8(174) include 'datcom.f90' include 'gcom1.f90' include 'gcom2.f90' - equivalence (id,d8) + equivalence (jd,d8) + equivalence (xd,c16) + equivalence (buf8,buf16) data nblock0/0/,kb/1/,ns00/99/,first/.true./ data sqave/0.0/,u/0.001/,rxnoise/0.0/,pctblank/0.0/,kbuf/1/,lost_tot/0/ data multicast0/-99/ @@ -36,6 +40,9 @@ subroutine recvpkt(iarg) 10 if(multicast.ne.multicast0) go to 1 call recv_pkt(center_freq) + iz=174 + if(nfloat.ne.0) iz=87 + ! Should receive a new packet every 174/96000 = 0.0018125 s nsec=mod(Tsec,86400.d0) !Time according to MAP65 nseclr=msec/1000 !Time according to Linrad @@ -60,8 +67,8 @@ subroutine recvpkt(iarg) if(transmitting.eq.1) ntx=1 ! Test for buffer full - if((kb.eq.1 .and. (k+174).gt.NSMAX) .or. & - (kb.eq.2 .and. (k+174).gt.2*NSMAX)) go to 20 + if((kb.eq.1 .and. (k+iz).gt.NSMAX) .or. & + (kb.eq.2 .and. (k+iz).gt.2*NSMAX)) go to 20 if(.not.first) then ! Check for lost packets @@ -72,10 +79,12 @@ subroutine recvpkt(iarg) nb0=nblock0 if(nb0.lt.0) nb0=nb0+65536 lost_tot=lost_tot + lost ! Insert zeros for the lost data. - do i=1,174*lost - k=k+1 - d8(k)=0 - enddo +!### +! do i=1,iz*lost +! k=k+1 +! d8(k)=0 +! enddo +!### endif endif first=.false. @@ -87,23 +96,42 @@ subroutine recvpkt(iarg) ! Move data into Rx buffer and compute average signal level. sq=0. - do i=1,174 + do i=1,iz k=k+1 - d8(k)=buf8(i) k2=k n=1 if(k.gt.NSMAX) then k2=k2-NSMAX n=2 endif - x1=id(1,k2,n) - x2=id(2,k2,n) - x3=id(3,k2,n) - x4=id(4,k2,n) - sq=sq + x1*x1 + x2*x2 + x3*x3 + x4*x4 + + if(nfloat.eq.0) then + d8=buf8(i) + x1=jd(1) + x2=jd(2) + x3=jd(3) + x4=jd(4) + dd(1,k2,n)=x1 + dd(2,k2,n)=x2 + dd(3,k2,n)=x3 + dd(4,k2,n)=x4 + sq=sq + x1*x1 + x2*x2 + x3*x3 + x4*x4 + else + c16=buf16(i) + x1=xd(1) + x2=xd(2) + x3=xd(3) + x4=xd(4) + dd(1,k2,n)=x1 + dd(2,k2,n)=x2 + dd(3,k2,n)=x3 + dd(4,k2,n)=x4 + sq=sq + x1*x1 + x2*x2 + x3*x3 + x4*x4 + endif enddo + sq=sq/(2.0*iz) sqave=sqave + u*(sq-sqave) - rxnoise=10.0*log10(sqave) - 48.0 + rxnoise=10.0*log10(sqave) - 20.0 ! Was -48.0 kxp=k 20 if(nsec.ne.nsec0) then diff --git a/symspec.f90 b/symspec.f90 index bba6b8147..4cb46ee5d 100644 --- a/symspec.f90 +++ b/symspec.f90 @@ -1,9 +1,9 @@ -subroutine symspec(id,kbuf,kk,kkdone,nutc,newdat) +subroutine symspec(dd,kbuf,kk,kkdone,nutc,newdat) ! Compute spectra at four polarizations, using half-symbol steps. parameter (NSMAX=60*96000) - integer*2 id(4,NSMAX,2) + real*4 dd(4,NSMAX,2) complex z real*8 ts,hsym include 'spcom.f90' @@ -41,10 +41,10 @@ subroutine symspec(id,kbuf,kk,kkdone,nutc,newdat) sq=0. do i=1,n1 !Find power in each block k=k+1 - x1=id(1,k,kbuf) - x2=id(2,k,kbuf) - x3=id(3,k,kbuf) - x4=id(4,k,kbuf) + x1=dd(1,k,kbuf) + x2=dd(2,k,kbuf) + x3=dd(3,k,kbuf) + x4=dd(4,k,kbuf) sq=sq + x1*x1 + x2*x2 + x3*x3 + x4*x4 enddo if(sq.lt.n1*10000.) then !Find power in good blocks @@ -65,19 +65,19 @@ subroutine symspec(id,kbuf,kk,kkdone,nutc,newdat) sq=0. do i=1,n1 k=k+1 - x1=id(1,k,kbuf) - x2=id(2,k,kbuf) - x3=id(3,k,kbuf) - x4=id(4,k,kbuf) + x1=dd(1,k,kbuf) + x2=dd(2,k,kbuf) + x3=dd(3,k,kbuf) + x4=dd(4,k,kbuf) sq=sq + x1*x1 + x2*x2 + x3*x3 + x4*x4 enddo ! If power in this block is excessive, blank it. if(sq.gt.1.5*sqave) then do i=k-n1+1,k - id(1,i,kbuf)=0 - id(2,i,kbuf)=0 - id(3,i,kbuf)=0 - id(4,i,kbuf)=0 + dd(1,i,kbuf)=0 + dd(2,i,kbuf)=0 + dd(3,i,kbuf)=0 + dd(4,i,kbuf)=0 enddo nclip=nclip+1 endif @@ -94,11 +94,11 @@ subroutine symspec(id,kbuf,kk,kkdone,nutc,newdat) i1=ts+2*hsym !Next starting sample pointer ts=ts+hsym !OK, update the exact sample pointer do i=1,npts !Copy data to FFT arrays - xr=fac*id(1,i0+i,kbuf) - xi=fac*id(2,i0+i,kbuf) + xr=fac*dd(1,i0+i,kbuf) + xi=fac*dd(2,i0+i,kbuf) cx(i)=cmplx(xr,xi) - yr=fac*id(3,i0+i,kbuf) - yi=fac*id(4,i0+i,kbuf) + yr=fac*dd(3,i0+i,kbuf) + yi=fac*dd(4,i0+i,kbuf) cy(i)=cmplx(yr,yi) enddo