c RDR2CSV is template for RDR output to CSV-style tables
c v0.99 added double precision function unsigned (int)
c v0.999 removed extraneous variables, updated comments, option to swap
c Fortran version may be replaced by a function programmed in the C language:
c double unsigned_ (unsigned int *s) {return (double) *s;}
c gfortran -ffixed-line-length-none -O -o rdr2csv rdr2csv.f unsigned.o kzext2.o iswap4.o spicelib.a
c ifort -132 -assume byterecl -O  -o rdr2csv rdr2csv_small.f unsigned.o kzext2.o iswap4.o spicelib_ifort.a
c smallendian changes are bracketed by @@@@
c v1.0 added changes for spare columns in RDR V2.4 1/19/10 and use of deltet for periodic term of TDB
c v1.1 reduced column width for thrsh gain, increase for flags, Earth centroid range 1/26/10
c v1.2 fixes bug in offset of 66.184 seconds between tdt and utc applied
c v1.3 fixes swap, but the function kzext2 must perform a short int swap.
	IMPLICIT NONE
	integer*4 rdr(64)
	integer*2 wdr(128)
	equivalence(rdr,wdr)
	integer*8 LMET
c the spacecraft 8-byte Long Mission Elapsed Time, or LMET,
c is given first in whole seconds and 32-bit fraction
c MET.fract are handled separately to maintain full precision
	equivalence(LMET,rdr)
c the Terrestrial Dynamic Time can be scaled to a real value  by division by 2^32
	integer*8 TDT
	equivalence(TDT,rdr(3))
        integer*4 MET !whole seconds 
	integer*4 missing /-2147483648/
	integer*4 iargc,itmp,numarg,iswap4
	integer*4 kzext2
	real*8 pi,r2d,d2r,r2de
	parameter(pi=3.141592653589794d0)
	parameter (r2d=180.d0/pi,d2r=pi/180.d0, r2de =9.d-03/pi)
	real*8 unsigned
	external unsigned
	real*8 solinc,solphs,offnadir,emission,earpw,earen
	integer*4 i,j,irec,orec,itdc,lseq,ifrm
        integer*4 subs,hz2fire
	real*8 d28, del2, utcsec
	parameter(d28=28.d0 + 1./8192.d0)! or 5/65536
c minor frame counts as fractions of 5000000.d0
	integer*4 minor(0:27)
	real*8 offset,fract,time,geoid
	real*8 c32 /4294967296.d0/
	real*8 et
	real*8 TX_EN_MJ,TX_PW,EAR_CE,EAR_PW
	real*8 sclon,sclat,scrad,radius
	real*8 gain(5),thrs(5),bkgd(5),plse(5),rnge(5),enrg(5)
	real*8 dlon(5),dlat(5),drad(5)
	integer nois(5),gflg(5)
	real*8 dscale/1.d-7/,gscale/1.d-6/, escale/1.d-6/
	real*8 pscale/1.d-3/,rscale/1.d-6/, tscale/1.d-6/
	character*28 utc, arg
        character*64 filename /'LOLARDR_00111N.DAT'/
	logical swap /.false./,head/.false./
c bitfield for status
c	integer*1 bTX,bRX1,bERX,bRX2,bRX3,bRX4,bRX5
c	data bTX/1/,bRX1/2/,bERX/4/,bRX2/8/,bRX3/16/,bRX4/32/,bRX5/64/
c initialize minor frame offsets - use to check range window?
c	minor(0)=0
c	do i=1,16
c	 minor(i)=minor(i-1)+178571
c	enddo
c	do i=17,27
c	 minor(i)=minor(i-1)+178572
c	enddo
	numarg=iargc() !flag for headers, swapping
	if(numarg.eq.0) then
	 write(*,*)'rdr2csv 1.3 - Converting reduced data record to a comma-separated-values file'
	 write(*,*)'usage: rdr2csv <rdrfile>  [[H|h](eader)] [[S|s](wap_bytes)]'
	 goto 99999
	endif
	do i=1,numarg
	 call getarg(i,arg)
	 if(arg(1:1).eq.'s' .or. arg(1:1).eq.'S') swap=.true.
	 if(arg(1:1).eq.'h' .or. arg(1:1).eq.'H') head=.true.
	enddo
c SPICE routines need this. No error trapping done.
	call ldpool('naif0009.tls')
        if(numarg.gt.0) call getarg(1,filename)
	open(unit=1,file=filename,recl=256,access='direct',status='old',err=9999)
	if(head) then
	write(*,'(a,$)')'Coordinated_Universal_Time,Frame,Mission_ET,   Subseconds,Terrestrial_Dyn_Time,'
	write(*,'(a,$)')' TX_Energy_mJ,TX_PulseW,SC_Longitude, SC_Latitude,   SC_radius,     Geoid,'
	do i=1,5
	write(*,5)i,i,i,i
	write(*,6)i,i,i,i,i,i
5	format('D',i1,'_Longitude, D',i1,'_Latitude,   D',i1,'_Radius,    D',i1,'_Range,',$)
6	format(' D',i1,'_PulseW, D',i1,'_Energy,D',i1,'_noi,D',i1,'_Thr, D',i1,'_Gn,Flg',i1',',$)
	enddo
c optional for spares
	write(*,7)
7	format('Offnadir,Emission, Sol_INC, Sol_Phs, Earth_Centr.,Earth_PW,Earth_E.')
	endif
        irec=0
1       continue
	j=0
        irec=irec+1
        read(1,rec=irec,err=999) rdr
c swap all words
	if(swap) then
	 do i=1,60
	  rdr(i)=iswap4(rdr(i))
	 enddo
	 else
c@@@@ smallendian
	  itmp=rdr(4)
	  rdr(4)=rdr(3)
	  rdr(3)=itmp
c@@@@
	endif
c
c offset from 2000T12:00:00 to 2001T00:00:00 adopted by LRO is already added to MET
c       utcsec =(tdt-tai2tdt-leap)/c32 in 8-byte arithmetic
	utcsec = tdt/c32-66.184d0 ! difference between TDT and UTC
	call DELTET (utcsec, 'UTC', DEL2)! Return the value of Delta ET (ET-UTC)
c	write(0,*)'del2',del2
	et = utcsec + DEL2
        call et2utc (et,'isoc',8,utc)
c	ifrm = mod(irec,28)
c better to take 28*fract as ifrm
	fract=unsigned(rdr(2))/c32 !cast fraction of met seconds to real*8
	ifrm = int(d28*fract)
c The first two columns are DERIVED from the rdr
	write(*,'(a28,",",i3,",",$)') utc,ifrm
	MET=rdr(1)
	write(*,'(i10,",",$)') MET
	write(*,'(f13.10,",",$)') fract
c or lmet/c32
	write(*,'(f20.9,",",$)') tdt/c32
c xmt energy
	TX_EN_MJ= 1.d-6*rdr(5)
c xmt width
	TX_PW = 1.D-3*rdr(6)
	write(*,'(f12.6,",",f10.3,",",$)')TX_EN_MJ,TX_PW
	if(rdr(7).ne.missing) then
	 sclon = dscale*rdr(7)
	 if(sclon.lt.0.) sclon=sclon+360.d0
	else
	 sclon=999.d0
	endif
	if(rdr(8).ne.missing) then
	 sclat = dscale*rdr(8)
	else
	 sclat = 99.d0
	endif
	if(rdr(9).ne.-1) then
	 scrad = rscale*unsigned(rdr(9))
	else
	 scrad =  -99.d0
	endif
	geoid = rscale*rdr(10)
	write(*,1000) sclon,sclat,scrad,geoid
1000	format(f12.7,",",f12.7,",",f12.6,",",f10.5,",",$)
	do i=1,5
	 if(rdr(i*10+1).ne.missing) then
	  dlon(i) = dscale*rdr(i*10+1)
	  if(dlon(i).lt.0.) dlon(i)=dlon(i)+360.d0
	 else
	  dlon(i)=999.d0
	 endif
	 if(rdr(i*10+2).ne.missing) then
	  dlat(i) = dscale*rdr(i*10+2)
	 else
	  dlat(i)= 99.d0
	 endif
	 if(rdr(i*10+3).ne. -1) then
	  drad(i) = rscale*unsigned(rdr(i*10+3))
	 else
	  drad(i) = -99.d0
	 endif
	 if(rdr(i*10+4).ne. -1) then
	  rnge(i) = rscale*unsigned(rdr(i*10+4))
	 else
	  rnge(i) = -99.d0
	 endif
	 if(rdr(i*10+5).ge.0) then! missing or invalid value is negative
	  plse(i) = pscale*rdr(i*10+5)
	 else
	  plse(i) = -99.d0
	 endif
	 enrg(i) = escale*rdr(i*10+6)
	 nois(i) =        rdr(i*10+7)
	 thrs(i) = tscale*rdr(i*10+8)
	 gain(i) = gscale*rdr(i*10+9)
	 gflg(i) = iand(255,rdr(i*10+10))
c gflg(i) = iand(255+1024,rdr(i*10+10))
c while threshold and gain are scaled by 10**6, they do not need such precision.
	 write(*,1004)dlon(i),dlat(i),drad(i),rnge(i),plse(i),enrg(i),nois(i),thrs(i),gain(i),gflg(i)
	enddo
1004	format(f12.7,",",f12.7,",",f12.6,",",f12.6,",",f10.3,",",f10.4,",",i6,",",f6.2,",",f6.2,",",i4,",",$)
c spares, now used for angles and earth ranges
	offnadir=r2de*kzext2(wdr(121))!always found
	if(wdr(122).ne. -1) then
	 emission=r2de*kzext2(wdr(122))
	else
	 emission=-1.d0! or zero?
	endif
	if(wdr(123).ne. -1) then
	 solinc=r2de*kzext2(wdr(123))
	else
	 solinc=-1.d0
	endif
	if(wdr(124).ne. -1) then
	 solphs=r2de*kzext2(wdr(124))
	else
	 solphs=-1.d0
	endif
c time from T0 to centroid of earth pulse, generally 1-9 ms
	EAR_CE=rdr(63)/c32
c last two words
c	 wdr(127)=1.d3*rxeren(i)!earth_energy using ch. 1 calibrations
c	 wdr(128)=1.d3*(EAR_PW)! ps
c
	if(wdr(128).ne. -1) then
	 earen=pscale*kzext2(wdr(127))
	 earpw=pscale*kzext2(wdr(128))
	else
	 earen=0.0
	 earpw=99.d0
	endif
	write(*,1005)offnadir,emission,solinc,solphs,EAR_CE,earen,earpw
1005	format
     & (f8.3,",",f8.3,",",f8.3,",",f8.3,",",f13.10,",",f8.3,",",f8.3)
c
        goto 1
c normal exit
999     continue
        close(1)
	goto 99999
c error exit
9999	continue
	write(0,*) 'file not opened successfully ', filename
99999   end
	
