C Program PODR_PREP C Needs subroutine DOY_MMDD C Run this program on a binary data file containing open loop receiver samples in RSC-11-9 format from Parkes Observatory. C The input file is assumed to have a file name between UL0304 and UL0359; its contents are from the Voyager 2 Uranus encounter. C Program reads RSC-11-9 Parkes Original Data Records (PODRs) and stores contents in ASCII header, C data, and log files. The governing RSC-11-9 SIS is dated 1982-01-11; it describes Voyager 2 open C loop radio science data acquired at Parkes Observatory in 1986. C USAGE: podr_prep infile C C where C C infile is the file name of the input PODR file from Voyager 2 Uranus C Output: C podr_prep.hdr ASCII CSV output file with PODR record header data C podr_prep.tab ASCII output file with PODR integer data samples C podr_prep.log ASCII log file C podr_prep.txt ASCII CSV file with column headers (accompanies podr_prep.hdr)) C Output records containing PODR samples have format (2i8,20i7,a1). The first number in each output C record is the input record number, starting from 1. The second number is the number of the last C sample in the previous output record. Then there are 20 (real) data samples with values 0-255. C The final byte is the ASCII line-feed character. The sample counter resets to 0 at the beginning C of each new input record. C Notes: C (1) a "word" has 16 bits C (2) All values in podr_prep.hdr are values from the original data records C (a) source station in data records is "43" whereas Parkes was designated "49" for Voyager; C (b) the '20 Counter' values in word 26 are decimal 23 whereas the SIS claims they should be "19" C (3) POCA frequency may change even when POCA rate is 0. C This program was written for (and tested on) an LSB (little endian) platform. The PODR data C format is based on MSB architecture. Because the data are read into a byte array (iarray1) C and sign conventions are recognized at that level, the program may also work on an MSB platform. C This possibility has not been checked since MSB platforms are not available in 2024. C Exits: C C 0 normal C 2 usage (incorrect command line arguments) C 3 error on call to stat C 4 number of records in file is not an integer C 5 error opening input file C 6 error opening output file for header data C 7 error opening output file for signal 1 sample data C 8 error opening file with field headers C 9 unexpected day of year (assumes 24 or 25 for VG2 at Uranus) C 10 error reading input record C 11 error writing output record C 12 error opening log file C 13 error return from call to doy_mmdd C 14 error opening file for warning messages C 15 unexpected c_flag (program is written for 8-bit samples) C 2024-05-12 RAS: Original; adapted from podr_2prt.f integer*1 iarray1(4090) !input record as a byte array integer*2 jarray2(4090) !input byte array mapped to array of 16-bit integers integer*4 ad1_select !signal being sampled by AD1 (word 16, bits 1-2) integer*4 ad2_select !signal being sampled by AD2 (word 16, bits 3-4) integer*4 ad3_select !signal being sampled by AD3 (word 16, bits 5-6) integer*4 ad4_select !signal being sampled by AD4 (word 16, bits 7-8) character*1 colon !ASCII colon character character*1 comma !ASCII comma character double precision counter_1_phase !frequency counter 1 cumulative phase (in 2^-8 cycles) (words 17-19) double precision counter_2_phase !frequency counter 2 cumulative phase (in 2^-8 cycles) (words 20-22) integer*4 counter_1_mode_reg !frequency counter 1 mode register (word 23, bits 9-12) integer*4 counter_2_mode_reg !frequency counter 2 mode register (word 23, bits 13-16) character*1 cr !ASCII carriage return character integer*4 c_flag !conversion flag (0 = 8-bit conversion; 1 = 12-bit) (word 1, bit 4) character*24 date_time !Date and time of program execution integer*4 day_of_year !Day of year starting with 001 for January 1 (word 5, bits 1-9) integer*4 dd !day of month integer*4 debug !0 = no diagnostic output; 1 = diagnostic output character*1 dl_band !downlink band (assumed to be "X") character*1 dl_poln !downlink polarizatyion (assumed to be "R") integer*4 e_flag !master tape error flag (0 = no error; 1 = error) (word 1, bit 3) integer*4 hh !hour character*1 hyphen !ASCII hyphen character integer*4 i !miscellaneous index integer*4 ier !error return flag character*256 infile !path to RSR data file (command line) integer*4 infile_len !number of characters in infile integer*4 ios !status flag integer*4 irecl !input record length (bytes) character*1 lf !ASCII line feed character integer*4 mn !minute integer*4 mon !month integer*4 n !miscellaneous index integer*4 narg !number of command line arguments integer*4 nbrec !calculated number of records in file integer*4 n_counter_value !N counter value (word 16, bits 9-16) integer*4 nrec !record number to read (default is 1) integer*4 nwrn !number of warning messages integer*4 ohrecs !header record counter integer*4 ones_1 !three bits with value = 1 (word 28, bits 2-4) integer*4 ones_2 !three bits with value = 1 (word 28, bits 10-12) integer*4 optr !pointer into ostring integer*4 osrec !sample record counter (per input record) integer*4 osrecs !sample record counter (cumulative) character*195 ostring !output string with record header data integer*4 o_flag !0 = values in words 6-9 invalid; 1 = valid (word 1, bit 1) integer*4 overflow_flag_1 !1 = overflow; 0 = no overflow (word 28, bit 1) integer*4 overflow_flag_2 !1 = overflow; 0 = no overflow (word 28, bit 9) character*1 period !ASCII period character integer*4 pgm_status !exit status double precision poca_freq !double precison representation of POCA frequency integer*8 poca_freq_hi !POCA frequency (integer Hertz) (word 9 bits 9-16 plus word 10 plus word 11 bits 1-8) integer*8 poca_freq_lo !POCA frequency (fractional Hertz) (word 11 bits 9-16 plus word 12) double precision poca_rate !POCA rate as a float integer*4 poca_rate_mantissa !POCA rate mantissa (word 13 bits 9-16 plus word 14 bits 1-12) integer*4 poca_rate_exp !POCA rate exponent (power of 10) (word 14 bits 13-15) integer*4 poca_rate_sign !POCA rate sign (word 14 bit 16) integer*4 poca_status(8) !POCA status (word 9, bits 1-8) character*4 predict_set_id !predict set identifier (words 13-16) integer*4 record_length !PODR record length (16-bit words) (word 3) integer*4 record_number !record number in PODR file (word 2) integer*4 sample_control_reg !sample control register (word 23, bits 5-8) integer*4 sample_rate !sample rate (word 15, bits 1-16) integer*4 sampling_mode_1 !how are source channels assigned to ADCs? (word 28, bits 7-8) integer*4 sampling_mode_2 !how are source channels assigned to ADCs? (word 28, bits 15-16) integer*4 sbits !bits per sample (assumed to be 8) integer*4 sc_id !spacecraft ID (word 4, bits 1-8) integer*4 seconds_of_day !seconds past 0h on day of year (word 5, bit 16; word 6) integer*4 short_conversion_1 !1 = 8-bit; 0 = 12-bit (word 28, bit 6) integer*4 short_conversion_2 !1 = 8-bit; 0 = 12-bit (word 28, bit 14) integer*4 source_station !source station (incorrectly set to 43 for VG2U) (word 4, bits 9-16) integer*4 spares !spares (word 24, bits 1-16) integer*4 srate !sample rate for sample file integer*4 ss !second character*24 start_date_time !time tag of first record integer*4 statb(13) !file status array character*24 stop_date_time !time tag of last record with valid time tag character*256 string !miscellaneous character string integer*4 s_flag !block sequence flag (0 = continuation; 1 = start) (word 1, bit 2) character*6 tape_id !tape identifier; argument for get_vgr_info integer*4 tape_number !tape number (word 1, bits 9-16) character*1 tee !ASCII "T" character integer*4 test_mode_1 !0 = normal; 1 = test (word 28,bit 5) integer*4 test_mode_2 !0 = normal; 1 = test (word 28,bit 13) integer*4 test_signal_select !test signal selection (word 23, bits 1-4) integer*4 t_flag !compression factor (valid values are 1, 2, and 10) (word 1, bits 5-8) integer*4 version(3)/2024,05,12/ !version of program integer*4 word_25 !two bytes with zeroes (word 25, bits 1-16) integer*4 word26_hi !20 counter - high byte (word 26 bits 1-8) integer*4 word26_lo !20 counter - low byte (word 26 bits 1-8) integer*4 word_27 !two bytes with zeroes (word 27, bits 1-16) integer*4 year !4-digit year C Initializations call fdate(date_time) colon = char(58) comma = char(44) cr = char(13) debug = 0 dl_band = "X" dl_poln = "R" hyphen = char(45) irecl = 4090 lf = char(10) nrec = 1 nwrn = 0 ohrecs = 0 osrecs = 0 period = char(46) pgm_status = 0 tee = "T" call system('rm -rf podr_prep.log podr_prep.hdr podr_prep.tab podr_prep.txt') C Check command line arguments narg = iargc() if (narg .eq. 1) then call getarg(1,infile) infile_len = index(infile," ") - 1 else write(*,'("USAGE: podr_prep infile")') pgm_status = 2 go to 99 endif C Get the file size and number of records ios = stat(infile(01:infile_len),statb) if (ios .ne. 0) then pgm_status = 3 write(*,'("PODR_PREP: ios = ",i5," when calling stat")') ios go to 99 endif nbrec = statb(8)/irecl if (nbrec*irecl .ne. statb(8)) then pgm_status = 4 write(*,'("PODR_PREP: file size from stat = ",i10," bytes ")') statb(8) write(*,'(" not equal to nbrec = ",i10," records")') nbrec write(*,'(" times irecl = ",i10," bytes ")') irecl write(*,'(" May be uneven records or truncated file")') go to 99 endif C Open output log file open ( unit = 30, * file = 'podr_prep.log', * iostat = ios, * status = 'UNKNOWN') if (ios .ne. 0) then write(*,'("PODR_PREP: error opening output log file ios = ",i5)') ios pgm_status = 12 go to 99 endif C Document the run write(*,'("PODR_PREP ----- converts PODR file to ASCII",/, * "Version of: ",20x,i4.4,"-",i2.2,"-",i2.2,/, * "Today is: ",(a),/, * "Input File: ",24x,(a))') * (version(i),i=1,3),date_time,infile(1:infile_len) write(30,'("PODR_PREP ----- converts PODR file to ASCII",/, * "Version of: ",20x,i4.4,"-",i2.2,"-",i2.2,/, * "Today is: ",(a),/, * "Input File: ",24x,(a))') * (version(i),i=1,3),date_time,infile(1:infile_len) C Open input PODR file with assumed record length open ( unit = 31, * file = infile, * access = 'DIRECT', * iostat = ios, * recl = irecl, * status = 'OLD') if (ios .ne. 0) then write(*,'("PODR_2PRT: error opening input PODR file ios = ",i5)') ios pgm_status = 5 go to 99 endif C Open output file for header data open ( unit = 32, * file = 'podr_prep.hdr', * access = 'DIRECT', * iostat = ios, * recl = 195, * status = 'UNKNOWN') if (ios .ne. 0) then write(*,'("PODR_2PRT: error opening output file for header data ios = ",i5)') ios pgm_status = 6 go to 99 endif C Open output file for sample data (assumes single input signal) open ( unit = 33, * file = 'podr_prep.tab', * iostat = ios, * status = 'UNKNOWN') if (ios .ne. 0) then write(*,'("PODR_2PRT: error opening output file for sample data ios = ",i5)') ios pgm_status = 7 go to 99 endif C Open file with field headers open ( unit = 34, * file = 'podr_prep.txt', * access = 'SEQUENTIAL', * iostat = ios, * status = 'UNKNOWN') if (ios .ne. 0) then write(*,'("PODR_2PRT: error opening file with field headers ios = ",i5)') ios pgm_status = 8 go to 99 endif write(34,'("NREC, O_Flag,S_Flag,E_Flag,C_Flag,T_Flag,Tape #,", * "Record #,RECL,S/C ID,DSS ID,DOY,Secs of Day,", * "Predict Set ID,POCA Control,Control Status,Synth Pwr,Synth Lock,", * "Limit Enable,Track,Acquisition,Sweep,POCA Frequency (Hz),POCA Rate (Hz/s),", * "ADC Sample Rate,J1 Sig Select,J2 Sig Select,J3 Sig Select,J4 Sig Select,N Counter,", * "F1 Phase (2^-8 cycles),F2 Phase (2^-8 cycles),Test Sig Select,Sample Control Register,", * "F1 Mode Register,F2 Mode Register,Spares-1,Zeroes-1,20-Counter,Zeroes-2,", * "Overflow Flag,Ones,Test Mode,Short Conversion,Sampling Mode")') write(34,'("i5,i1,i1,i1,i1,i2,i2,i5,i4,i2,i2,i3,i5,a4,i1,i1,i1,i1,i1,i1,i1,i1,f16.6,e12.5,i6,i1,i1,", * "i1,i1,i3,d22.15,d22.15,i2,i2,i2,i2,i1,i1,i2,i1,i1,i1,i1,i1,i1")') close(34) C Enter file for warning messages open ( unit = 35, * file = 'podr_prep.wrn', * access = 'SEQUENTIAL', * iostat = ios, * status = 'UNKNOWN') if (ios .ne. 0) then write(*,'("PODR_2PRT: error opening file with warning messages ios = ",i5)') ios pgm_status = 14 go to 99 endif write(35,'("PODR_PREP ----- converts PODR file to ASCII",/, * "Version of: ",20x,i4.4,"-",i2.2,"-",i2.2,/, * "Today is: ",(a),/, * "Input File: ",24x,(a))') * (version(i),i=1,3),date_time,infile(1:infile_len) C Enter loop to read and unpack data; write record number into output string do 50 nrec = 1,nbrec read(31,rec=nrec,err=91) (iarray1(i),i=1,irecl) optr = 1 write(ostring(optr:optr+5),'(i5,a1)') nrec,comma optr = optr+6 C Convert byte array to int*2 to preserve unsigned byte values do i = 1,irecl if (iarray1(i) .ge. 0) then jarray2(i) = iarray1(i) else jarray2(i) = 256 + iarray1(i) endif end do C Unpack header byte 1 (first half of word 1) o_flag = int(jarray2(1)/128) s_flag = int(jarray2(1)/64) - 2*o_flag e_flag = int(jarray2(1)/32) - 4*o_flag - 2*s_flag c_flag = int(jarray2(1)/16) - 8*o_flag - 4*s_flag - 2*e_flag t_flag = mod(jarray2(1),16) if (c_flag .eq. 0) then sbits = 8 else write(*,'("unexpected c_flag = ",i3)') c_flag pgm_status = 15 go to 99 endif if (debug .eq. 1) write(*,'("O flag = ",13x,i9,/, * "S flag = ",13x,i9,/, * "E flag = ",13x,i9,/, * "C flag = ",13x,i9,/, * "T flag = ",13x,i9)') * o_flag, * s_flag, * e_flag, * c_flag, * t_flag write(ostring(optr:optr+10),'(i1,a1,i1,a1,i1,a1,i1,a1,i2,a1)') * o_flag,comma, * s_flag,comma, * e_flag,comma, * c_flag,comma, * t_flag,comma optr = optr+11 C Unpack header bytes 2-16 (second half of word 1 through word 8) tape_number = jarray2(02) record_number = 256*jarray2(03) + jarray2(04) record_length = 256*jarray2(05) + jarray2(06) sc_id = jarray2(07) source_station = jarray2(08) day_of_year = 2*jarray2(09)+jarray2(10)/128 if ((day_of_year .eq. 22) .or. (day_of_year .eq. 23) .or. * (day_of_year .eq. 24) .or. (day_of_year .eq. 25)) then year = 1986 else pgm_status = 9 write(*,'("Unexpected day or year = ",i5,5x,"at record = ",i5)') day_of_year,nrec go to 99 endif call doy_mmdd(1,year,day_of_year,mon,dd,ios) if (ios .ne. 0) then pgm_status = 13 write(*,'("Error return from call to doy_mmdd = ",i3,5x,"at record = ",i5,/, * "Mode = 1",5x,"year = ",i5,5x,"day of year = ",i5)') * ios,nrec,year,day_of_year go to 99 endif seconds_of_day = 256*256*mod(jarray2(10),128)+256*jarray2(11)+jarray2(12) hh = int(seconds_of_day/3600) mn = int((seconds_of_day-3600*hh)/60) ss = seconds_of_day - 3600*hh - 60*mn write(predict_set_id,'(4a1)') (jarray2(i),i=13,16) if (debug .eq. 1) write(*,'("tape number = ",13x,i9,/, * "record_number = ",13x,i9,/, * "record_length (words) = ",13x,i9,/, * "spacecraft ID = ",13x,i9,/, * "source station ID = ",13x,i9,/, * "day of year = ",13x,i9,/, * "seconds of day = ",13x,i9,/, * "predict set ID = ",18x,a4)') * tape_number, * record_number, * record_length, * sc_id, * source_station, * day_of_year, * seconds_of_day, * predict_set_id write(ostring(optr:optr+35),'(i3,a1,i5,a1,i4,a1,i2,a1,i2,a1,i3.3,a1,i5.5,a1,a4,a1)') * tape_number,comma, * record_number, comma, * record_length,comma, * sc_id,comma, * source_station,comma, * day_of_year,comma, * seconds_of_day,comma, * predict_set_id,comma optr = optr + 36 if (nrec .eq. 1) then write(start_date_time,'(i4.4,a1,i3.3,a1,i2.2,a1,i2.2,a1,i2.2,a1,"000000")') * year,hyphen,day_of_year,tee,hh,colon,mn,colon,ss,period endif if (o_flag .eq. 1) then write(stop_date_time,'(i4.4,a1,i3.3,a1,i2.2,a1,i2.2,a1,i2.2,a1,"000000")') * year,hyphen,day_of_year,tee,hh,colon,mn,colon,ss,period endif C Unpack bytes 17-28 (words 9-14) poca_status(1) = jarray2(17)/128 poca_status(2) = (jarray2(17)-128*poca_status(1))/64 poca_status(3) = (jarray2(17)-128*poca_status(1)-64*poca_status(2))/32 poca_status(4) = (jarray2(17)-128*poca_status(1)-64*poca_status(2)-32*poca_status(3))/16 poca_status(5) = (jarray2(17)-128*poca_status(1)-64*poca_status(2)-32*poca_status(3)-16*poca_status(4))/8 poca_status(6) = (jarray2(17)-128*poca_status(1)-64*poca_status(2)-32*poca_status(3)-16*poca_status(4)-8*poca_status(5))/4 poca_status(7) = (jarray2(17)-128*poca_status(1)-64*poca_status(2)-32*poca_status(3)-16*poca_status(4)-8*poca_status(5)-4*poca_status(6))/2 poca_status(8) = jarray2(17)-128*poca_status(1)-64*poca_status(2)-32*poca_status(3)-16*poca_status(4)-8*poca_status(5)-4*poca_status(6)-2*poca_status(7) poca_freq_hi = (jarray2(18)/16)*10000000 + * mod(jarray2(18),16)*1000000 + * (jarray2(19)/16)*100000 + * mod(jarray2(19),16)*10000 + * (jarray2(20)/16)*1000 + * mod(jarray2(20),16)*100 + * (jarray2(21)/16)*10 + * mod(jarray2(21),16) poca_freq_lo = (jarray2(22)/16)*100000 + * mod(jarray2(22),16)*10000 + * (jarray2(23)/16)*1000 + * mod(jarray2(23),16)*100 + * (jarray2(24)/16)*10 + * mod(jarray2(24),16) poca_freq = 1.0d+00*poca_freq_hi + 1.0d-06*poca_freq_lo poca_rate_mantissa = (jarray2(26)/16)*10000 + * mod(jarray2(26),16)*1000 + * (jarray2(27)/16)*100 + * mod(jarray2(27),16)*10 + * (jarray2(28)/16) poca_rate_exp = mod(jarray2(28),16)/2 if (mod(jarray2(28),2) .eq. 0) then poca_rate_sign = -1.0d+00 else poca_rate_sign = 1.0d+00 endif poca_rate = poca_rate_sign*poca_rate_mantissa*10**poca_rate_exp if (debug .eq. 1) write(*,'("POCA Status: Control Mode = ",18x,i4,/, * " Control Ready = ",18x,i4,/, * " Synthesizer Power = ",18x,i4,/, * " Synthesizer Lock = ",18x,i4,/, * " Limit Enable = ",18x,i4,/, * " Track = ",18x,i4,/, * " Acquisition = ",18x,i4,/, * " Sweep = ",18x,i4,/, * "POCA Frequency (Hertz) = ",6x,f16.6,/, * "POCA Rate (Hertz per second) = ",10x,e12.5)') * (poca_status(i),i=1,8), * poca_freq, * poca_rate write(ostring(optr:optr+45),'(8(i1,a1),f16.6,a1,e12.5,a1)') * (poca_status(i),comma,i=1,8), * poca_freq, comma, * poca_rate, comma optr = optr + 46 C Unpack bytes 29-32 (words 15-16) sample_rate = 256*jarray2(29)+jarray2(30) ad1_select = jarray2(31)/64 ad2_select = (jarray2(31)-64*ad1_select)/16 ad3_select = (jarray2(31)-64*ad1_select-16*ad2_select)/4 ad4_select = jarray2(31)-64*ad1_select-16*ad2_select-4*ad3_select n_counter_value = jarray2(32) if (debug .eq. 1) write(*,'("Sample rate per ADC = ",12x,i10,/, * "AD1 signal select = ",18x,i4,/, * "AD2 signal select = ",18x,i4,/, * "AD3 signal select = ",18x,i4,/, * "AD4 signal select = ",18x,i4,/, * "N counter value = ",18x,i4)') * sample_rate, * ad1_select, * ad2_select, * ad3_select, * ad4_select, * n_counter_value write(ostring(optr:optr+18),'(i6,a1,4(i1,a1),i3,a1)') * sample_rate,comma, * ad1_select,comma, * ad2_select,comma, * ad3_select,comma, * ad4_select,comma, * n_counter_value,comma optr = optr + 19 C Unpack bytes 33-44 (words 17-22) counter_1_phase = 4.294967296d+09*(256*jarray2(33) + jarray2(34)) + 16777216*jarray2(35) + 65536*jarray2(36) + 256*jarray2(37) + jarray2(38) counter_2_phase = 4.294967296d+09*(256*jarray2(39) + jarray2(40)) + 16777216*jarray2(41) + 65536*jarray2(42) + 256*jarray2(43) + jarray2(44) if (debug .eq. 1) write(*,'("Counter 1 Phase = ",1p1d22.15,/, * "Counter 2 Phase = ",1p1d22.15)') * counter_1_phase, * counter_2_phase write(ostring(optr:optr+45),'(1p1e22.15,a1,1p1e22.15,a1)') * counter_1_phase,comma, * counter_2_phase,comma optr = optr + 46 C Unpack bytes 45-56 (words 23-28) test_signal_select = jarray2(45)/256 sample_control_reg = mod(jarray2(45),16) counter_1_mode_reg = jarray2(46)/256 counter_2_mode_reg = mod(jarray2(46),16) spares = 256*jarray2(47)+jarray2(48) word_25 = 256*jarray2(49) + jarray2(50) word26_hi = jarray2(51) word26_lo = jarray2(52) word_27 = 256*jarray2(53) + jarray2(54) overflow_flag_1 = jarray2(55)/128 ones_1 = (jarray2(55)-128*overflow_flag_1)/16 test_mode_1 = (jarray2(55)-128*overflow_flag_1-16*ones_1)/8 short_conversion_1 = (jarray2(55)-128*overflow_flag_1-16*ones_1-8*test_mode_1)/4 sampling_mode_1 = (jarray2(55)-128*overflow_flag_1-16*ones_1-8*test_mode_1-4*short_conversion_1) if (sampling_mode_1 .eq. 0) then srate = sample_rate/1000 else if (sampling_mode_1 .eq. 1) then srate = 4*sample_rate/1000 else if (sampling_mode_1 .eq. 2) then srate = 2*sample_rate/1000 else write(*,'("Unexpected sampling_mode_1 = ",i3,5x,"Unable to compute srate")') sampling_mode_1 endif overflow_flag_2 = jarray2(56)/128 ones_2 = (jarray2(56)-128*overflow_flag_2)/16 test_mode_2 = (jarray2(56)-128*overflow_flag_2-16*ones_2)/8 short_conversion_2 = (jarray2(56)-128*overflow_flag_2-16*ones_2-8*test_mode_2)/4 sampling_mode_2 = (jarray2(56)-128*overflow_flag_2-16*ones_2-8*test_mode_2-4*short_conversion_2) C Record warning messages if expected identical values are not the same if (word26_hi .ne. word26_lo ) then write(35,'("Record = ",i5,5x,"Values for 20 Counter differ: word 26 hi = ",i4,5x,"lo = ",i4)') nrec,word26_hi,word26_lo nwrn = nwrn + 1 endif if ((overflow_flag_1 .ne. overflow_flag_2) .or. * (ones_1 .ne. ones_2) .or. * (test_mode_1 .ne. test_mode_2) .or. * (short_conversion_1 .ne. short_conversion_2) .or. * (sampling_mode_1 .ne. sampling_mode_2)) then write(35,'("Record = ",i5,5x,"Word 28 flags differ: hi = ",5i2,5x,"lo = ",5i2)') * nrec,overflow_flag_1,ones_1,test_mode_1,short_conversion_1,sampling_mode_1, * overflow_flag_2,ones_2,test_mode_2,short_conversion_2,sampling_mode_2 nwrn = nwrn + 1 endif if (debug .eq. 1) write(*,'("Counter status: test signal sel = ",18x,i4,/, * " samp ctrl reg = ",18x,i4,/, * " counter 1 reg = ",18x,i4,/, * " counter 2 reg = ",18x,i4,/, * "Spares = ",8x,i14,/, * "word_25 = ",8x,i14,/, * "20 Counter hi = ",8x,i14,/, * "20 Counter lo = ",8x,i14,/, * "word_27 = ",8x,i14,/, * "overflow_flag_1 = ",8x,i14,/, * "ones_1 = ",8x,i14,/, * "test_mode_1 = ",8x,i14,/, * "short_conversion_1 = ",8x,i14,/, * "sampling_mode_1 = ",8x,i14,/, * "overflow_flag_2 = ",8x,i14,/, * "ones_2 = ",8x,i14,/, * "test_mode_2 = ",8x,i14,/, * "short_conversion_2 = ",8x,i14,/, * "sampling_mode_2 = ",8x,i14)') * test_signal_select, * sample_control_reg, * counter_1_mode_reg, * counter_2_mode_reg, * spares, * word_25, * word26_hi, * word26_lo, * word_27, * overflow_flag_1, * ones_1, * test_mode_1, * short_conversion_1, * sampling_mode_1, * overflow_flag_2, * ones_2, * test_mode_2, * short_conversion_2, * sampling_mode_2 write(ostring(optr:optr+30),'(4(i2,a1),2(i1,a1),i2,a1,5(i1,a1),i1,a1)') * test_signal_select,comma, * sample_control_reg,comma, * counter_1_mode_reg,comma, * counter_2_mode_reg,comma, * spares,comma, * word_25,comma, * word26_hi,comma, * word_27,comma, * overflow_flag_1,comma, * ones_1,comma, * test_mode_1,comma, * short_conversion_1,comma, * sampling_mode_1,lf optr = optr + 31 C Write the unpacked header data to the output header file write(32,rec=nrec,err=92) ostring ohrecs = ohrecs + 1 C Write the sample values to the output sample file. C 20 ASCII samples per output record, repeated until all samples from input record are exhausted. osrec = 0 do n = 1,irecl-90,20 osrec = osrec + 1 osrecs = osrecs + 1 write(33,'(2i8,20i7)',err=94) nrec,n-1,(jarray2(55+n+i),i=1,20) end do C If this is the first record, add to log file if (nrec .eq. 1) then write(*,'("Input record length (bytes): ",i10,/, * "Input records in PODR file: ",i10,/, * "Sample rate (real ksps): ",i10,/, * "Bits per sample: ",i10,/, * "Spacecraft ID: ",i10,/, * "Receiving station: ",i10,/, * "Receiving band: ",9x,a1,/, * "Receiving polarization: ",9x,a1)') * irecl,nbrec,srate,sbits,sc_id,source_station,dl_band,dl_poln write(30,'("Input record length (bytes): ",i10,/, * "Input records in PODR file: ",i10,/, * "Sample rate (real ksps): ",i10,/, * "Bits per sample: ",i10,/, * "Spacecraft ID: ",i10,/, * "Receiving station: ",i10,/, * "Receiving band: ",9x,a1,/, * "Receiving polarization: ",9x,a1)') * irecl,nbrec,srate,sbits,sc_id,source_station,dl_band,dl_poln endif 50 continue go to 97 C Exits 91 continue write(*,'("Error attempting to read record nrec = ",i9)') nrec pgm_status = 10 go to 99 92 continue write(*,'("Error attempting to write record nrec = ",i9)') nrec pgm_status = 11 go to 99 94 continue write(*,'("Error attempting to write output sample record osrec = ",i9,5x,"at input record = ",i5)') osrec,nrec pgm_status = 15 go to 99 C Add data sample record count and first and last sample times to standard output and log file 97 continue write(*,'("Header output to: podr_prep.hdr",/, * " output record length (bytes): "7x,"202",/, * " output record count: ",i10,/, * "ASCII sample output to: podr_prep.tab",/, * " output record length (bytes): ",7x,"157",/, * " output sample record count: ",i10)') ohrecs,osrecs write(*,'("First time tag: ",a24)') start_date_time(01:24) write(*,'("Last valid time tag: ",a24)') stop_date_time(01:24) write(30,'("Header output to: podr_prep.hdr",/, * " output record length (bytes): "7x,"202",/, * " output record count: ",i10,/, * "ASCII sample output to: podr_prep.tab",/, * " output record length (bytes): ",7x,"157",/, * " output sample record count: ",i10)') ohrecs,osrecs write(30,'("First time tag: ",a24)') start_date_time(01:24) write(30,'("Last valid time tag: ",a24)') stop_date_time(01:24) write(*,'("Number of warning messages:",13x,i10)') nwrn write(30,'("Number of warning messages:",13x,i10)') nwrn go to 99 C Exit 99 continue close(30) close(31) close(32) close(33) close(35) call exit(pgm_status) stop end C ------------------------------------------------------------------------------------------------------------------------ C Program PODR_HIST C Run this program on the output from PODR_PREP to generate a histogram of the open loop receiver samples. C Program reads output files from podr_prep and calculates histograms for sample values. C Obtains sampling parameters from podr_prep.hdr and obtains sample values from podr_prep.tab C USAGE: podr_hist C Output: C podr_hist.tab ASCII output file -- a table with two columns: C 1) sample value (f12.1) C 2) the number of samples at that value (f12.0) C Exits: C 0 normal C 2 usage (incorrect command line arguments; none are expected C 3 error returned by get_date_time C 4 error opening PODR_PREP log file (input) C 5 error opening PODR_HIST log file (output) C 6 error opening file with ASCII sample values (input) C 7 error opening histogram file (output) C Logical unit numbers: C 29 podr_prep.log C 30 podr_hist.log C 31 podr_hist.tab C 33 podr_prep.tab C 2024-05-11 RAS: Original; based on rsc11-6_hist.f character*1 comma !ASCII comma character integer*4 debug !0 = no diagnostics; 1 = verbose diagnostic character*1 dlband !identifier for downlink band (K, S, or X) integer*4 dlpoln !identifier for downlink polarization (R or L) integer*4 dss !identifier for DSN antenna (e.g., 14, 15, 43, 65, etc.) double precision hist(256) !histogram of data samples integer*4 i !miscellaneous index integer*4 ier !error flag from get_date_time character*256 infile !name of (or path to) original PODR file integer*4 infile_len !number of characters in infile integer*4 ios !I/O status flag integer*4 j !miscellaneous index integer*4 narg !number of command line arguments integer*4 narec !number of records in ASCII file of data samples integer*4 nbrec !number of records in original PODR file integer*4 pgm_status/0/ !status upon exit from program integer*4 run_date_time(6) !date and time of program execution integer*4 samp(20) !samples in one record of input file integer*4 sbits !bits per sample integer*4 scid !identifier for spacecraft integer*4 srate !sampling rate (ksps) character*24 start_date_time !time tag for first PODR record character*24 stop_date_time !time tag for last valid PODR record real*4 value !sample value within histogram integer*4 version(3)/2024,05,11/ !version of program integer*4 iargc !Fortran library function C Check command line arguments narg = iargc() if (narg .ne. 0) then write(*,'("USAGE: No command line arguments expected")') pgm_status = 2 go to 99 end if C Initializations comma = char(44) debug = 1 call get_date_time(run_date_time,ier) if (ier .ne. 0) then pgm_status = 3 go to 99 end if do i = 1,256 hist(i) = 0.0e+00 end do C Get name of (or path to) original PODR file open ( unit = 29, * file = 'podr_prep.log', * access = 'SEQUENTIAL', * iostat = ios, * status = 'OLD') if (ios .ne. 0) then write(*,'("PODR_HIST: error opening PODR_PREP log file", * " ios = ",i5)') ios pgm_status = 4 go to 99 end if read(29,'(///,44x,(a))') infile infile_len = index(infile," ") - 1 read(29,'(/,5(38x,i12,/),2(49x,a1,/),/////,38x,i12,2(/,26x,a24))') * nbrec,srate,sbits,scid,dss,dlband,dlpoln,narec,start_date_time,stop_date_time close(29) if ((start_date_time(01:08) .eq. "1986-024") .or. (start_date_time(01:08) .eq. "1986-025")) then dss = 49 !Correct the DSS ID for VG2U at Parkes end if C Open output log file open ( unit = 30, * file = 'podr_hist.log', * access = 'sequential', * iostat = ios, * status = 'unknown') if (ios .ne. 0) then write(*,'("PODR_HIST: error opening log file", * " ios = ",i5)') ios pgm_status = 5 go to 99 end if C Document the run write(*,'("PODR_HIST: Create histogram for Real Data Samples",/, * "Version of:",29x,i4.4,"-",i2.2,"-",i2.2,/, * "Today is: ",20x,i4.4,2("-",i2.2),"T",i2.2,2(":",i2.2)," local time",/, * "PODR file name: ",19x,(a),/, * "Records in PODR file: ",15x,i10,/, * "Records in ASCII file: ",15x,i10,/, * "Sample rate (real sps): ",15x,i10,/, * "Sample bits: ",15x,i10,/, * "Spacecraft ID: ",15x,i10,/, * "DSS ID: ",15x,i10,/, * "Receive band: ",24x,a1,/, * "Receive polarization: ",24x,a1,/, * "Start date and time: ",1x,a24,/, * "Stop date and time: ",1x,a24)') * version,run_date_time,infile(1:infile_len), * nbrec,narec,srate,sbits,scid,dss,dlband,dlpoln, * start_date_time,stop_date_time write(30,'("PODR_HIST: Create histogram for Real Data Samples",/, * "Version of:",29x,i4.4,"-",i2.2,"-",i2.2,/, * "Today is: ",20x,i4.4,2("-",i2.2),"T",i2.2,2(":",i2.2)," local time",/, * "PODR file name: ",19x,(a),/, * "Records in PODR file:",19x,i10,/, * "Records in ASCII file: ",15x,i10,/, * "Sample rate (real sps): ",15x,i10,/, * "Sample bits: ",15x,i10,/, * "Spacecraft ID: ",15x,i10,/, * "DSS ID: ",15x,i10,/, * "Receive band: ",24x,a1,/, * "Receive polarization: ",24x,a1,/, * "Start date and time: ",1x,a24,/, * "Stop date and time: ",1x,a24)') * version,run_date_time,infile(1:infile_len), * nbrec,narec,srate,sbits,scid,dss,dlband,dlpoln, * start_date_time,stop_date_time C Open file with sample values open ( unit = 33, * file = 'podr_prep.tab', * iostat = ios, * status = 'old') if (ios .ne. 0) then write(*,'("PODR_HIST: error opening ASCII sample file", * " ios = ",i5)') ios pgm_status = 6 go to 99 end if C Read sample values and construct histogram do i = 1,narec read(33,'(16x,20i7)') (samp(j),j=1,20) do j=1,20 hist(samp(j)+1) = hist(samp(j)+1) + 1.0e+00 end do end do C Open file for histogram output; write histogram and close file open ( unit = 31, * file = 'podr_hist.tab', * access = 'sequential', * iostat = ios, * status = 'unknown') if (ios .ne. 0) then write(*,'("PODR_HIST: error opening output histogram file", * " ios = ",i5)') ios pgm_status = 7 go to 99 end if do i = 1,256 value = i - 1. write(31,'(f12.1,a1,f12.0)') value,comma,hist(i) end do close(31) C Exit 99 continue write(*,'("Done")') call exit(pgm_status) stop end C ------------------------------------------------------------------------------------------------------------------------ C Program PODR_POWR C Run this program on the output from PODR_PREP to generate a table of average sample power versus time. C Program reads output files from podr_prep and calculates normalized average sample power versus time. C Normalization is by the maximum power in a sample and by the number of samples in the averaging interval. C Obtains sampling parameters from podr_prep.log and obtains sample values from podr_prep.tab C USAGE: podr_powr tavg0 C where tavg0 is the requested averaging time in seconds C Output: C podr_powr.tab ASCII output file -- a table of power vs time in two columns: C 1) time in seconds from midnight of the first sample (f14.6) C 2) power averaged over tavg1 (1p1e14.6) and scaled for C maximum ADC power possible with sbits sampling C NB: tavg1 is the actual averaging time, which may C be slightly different from tavg0 C Each record ends with an ASCII carriage-return line-feed pair C podr_powr.log ASCII log file from execution of podr_powr C Exits: C 0 normal C 2 usage (incorrect command line arguments; only tavg0 is expected) C 3 error in subroutine get_date_time C 4 error opening podr_prep.log (input) file C 5 error opening log file (output) C 9 error opening sample file (input) C 10 error opening output average power file C Logical unit numbers: C 28 podr_powr.temp C 29 podr_prep.log !log file from podr_prep (input here) C 30 podr_powr.log !log file for this program C 31 podr_powr.tab !output data file from this program C 33 podr_prep.tab !sample file from podr_prep (input to this program) C 2024-05-11 RAS: Original; adapted from rsc11-6_powr.f integer*4 arecl !length of data portion of PODR record (4000 bytes) integer*4 brecl !record length in PODR (binary) file character*1 cr !ASCII carriage-return character integer*4 ddd0 !first time tag in PODR file, day of year integer*4 ddd1 !last time tag in PODR file, day of year integer*4 debug !1 = verbose output for debugging; 0 = silent character*1 dlband !identifier for downlink band (K, S, or X) integer*4 dlpoln !received polarization (R or L) integer*4 dss !identifier for DSN antenna (e.g., 49 for Parkes, after correction) double precision dt !time represented by one sample (seconds) integer*4 ffffff0 !first time tag in PODR file, microseconds integer*4 ffffff1 !last time tag in PODR file, microseconds integer*4 hh0 !first time tag in PODR file, hour integer*4 hh1 !last time tag in PODR file, hour integer*4 i !miscellaneous index integer*4 ier !subroutine return error flag character*256 infile !name of (or path to) original RSR file integer*4 infile_len !number of characters in infile integer*4 ios !I/O status flag integer*4 j !miscellaneous index integer*4 mm0 !first time tag in PODR file, minute integer*4 mm1 !last time tag in PODR file, minute integer*4 n !miscellaneous index integer*4 narec !number of records in file of ASCII data samples integer*4 narg !number of command line arguments integer*4 nbrec !records in PODR (binary) file integer*4 npts !number of points in output file integer*4 nr_avg !number of ASCII data records in each average integer*4 ns_rec !number of samples per ASCII record integer*4 ns_avg !number of samples per average integer*4 n_av !number of samples in the ASCII data file integer*4 pgm_status !status upon exit from program double precision power !accumulated (then averaged) power integer*4 run_date_time(6) !date and time of program execution integer*4 samp(20) !samples in one ASCII record double precision samp_av !average sample value over the file double precision samp_sq !average squared sample value over the file integer*4 sbits !bits per sample integer*4 scid !identifier for spacecraft double precision sd !standard deviation of sample values over the file double precision sf !full range scale factor integer*4 srate !sampling rate (ksps) integer*4 ss0 !first time tag in PODR file, integer seconds integer*4 ss1 !last time tag in PODR file, integer seconds character*80 string !miscellaneous string double precision tavg0 !requested averaging time (seconds) double precision tavg1 !actual averaging time (seconds) double precision t1 !output power times (seconds) integer*4 version(3)/2024,05,11/ !version of program integer*4 yyyy0 !from PODR_PREP command line, year of first time tag integer*4 yyyy1 !from PODR_PREP command line, year of last time tag integer*4 iargc !Fortran library function C Check command line arguments; get averaging time narg = iargc() if (narg .eq. 1) then call getarg(1,string) read(string,*) tavg0 else write(*,'("USAGE: podr_powr tavg0")') pgm_status = 2 go to 99 end if C Initializations cr = char(13) debug = 0 ns_rec = 20 pgm_status = 0 call get_date_time(run_date_time,ier) if (ier .ne. 0) then pgm_status = 3 go to 99 end if C Get name of (or path to) original PODR file open ( unit = 29, * file = 'podr_prep.log', * access = 'sequential', * iostat = ios, * status = 'old') if (ios .ne. 0) then write(*,'("PODR_POWR: error opening input PODR_PREP log file", * " ios = ",i5)') ios pgm_status = 4 go to 99 end if read(29,'(///,44x,(a))') infile infile_len = index(infile," ") - 1 C Open output log file open ( unit = 30, * file = 'podr_powr.log', * access = 'sequential', * iostat = ios, * status = 'unknown') if (ios .ne. 0) then write(*,'("PODR_POWR: error opening output log file", * " ios = ",i5)') ios pgm_status = 5 go to 99 end if C Document the run write(*,'("PODR_POWR: Calculate average sample power versus time",/, * "Version of:",29x,i4.4,"-",i2.2,"-",i2.2,/, * "Today is: ",20x,i4.4,2("-",i2.2),"T",i2.2,2(":",i2.2)," local time",/, * "PODR File name:",29x,(a))') * version,run_date_time,infile(1:infile_len) write(30,'("PODR_POWR: Calculate average sample power versus time",/, * "Version of ",29x,i4.4,"-",i2.2,"-",i2.2,/, * "Today is: ",20x,i4.4,2("-",i2.2),"T",i2.2,2(":",i2.2)," local time",/, * "PODR File name:",29x,(a))') * version,run_date_time,infile(1:infile_len) C Read sampling parameters from podr_prep log file read(29,'(6(40x,i10,/),2(49x,a1,/),//////,26x,i4,1x,i3,1x,i2,1x,i2,1x,i2,1x,i6,/, * 26x,i4,1x,i3,1x,i2,1x,i2,1x,i2,1x,i6)') * brecl,nbrec,srate,sbits,scid,dss,dlband,dlpoln, * yyyy0,ddd0,hh0,mm0,ss0,ffffff0,yyyy1,ddd1,hh1,mm1,ss1,ffffff1 arecl = brecl - 90 narec = arecl*nbrec/ns_rec if (yyyy0 .eq. 1986) dss = 49 ! Correct DSS ID for Parkes data in 1986 write(*,'("PODR record length (brecl): ",7x,i6,/, * "PODR record count (nbrec): ",7x,i6,/, * "PODR sample rate (srate, ksps): ",7x,i6,/, * "PODR samples per record (arecl): ",7x,i6,/, * "PODR sample bits (sbits): ",7x,i6,/, * "Samples per ASCII record (ns_rec): ",7x,i6,/, * "Records in ASCII sample file (narec):",4x,i9,/, * "Spacecraft number (scid): ",7x,i6,/, * "DSN antenna (dss): ",7x,i6,/, * "Downlink band (dlband): ",12x,a1,/, * "Downlink polarization (dlpoln): ",12x,a1,/, * "First time tag:",11x,i4.4,"-",i3.3,"T",2(i2.2,":"),i2.2,".",i6.6,/, * "Last time tag: ",11x,i4.4,"-",i3.3,"T",2(i2.2,":"),i2.2,".",i6.6)') * brecl,nbrec,srate,arecl,sbits,ns_rec,narec,scid,dss,dlband,dlpoln, * yyyy0,ddd0,hh0,mm0,ss0,ffffff0,yyyy1,ddd1,hh1,mm1,ss1,ffffff1 write(30,'("PODR record length (brecl): ",7x,i6,/, * "PODR record count (nbrec): ",7x,i6,/, * "PODR sample rate (srate, ksps): ",7x,i6,/, * "PODR samples per record (arecl): ",7x,i6,/, * "PODR sample bits (sbits): ",7x,i6,/, * "Samples per ASCII record (ns_rec): ",7x,i6,/, * "Records in ASCII sample file (narec):",4x,i9,/, * "Spacecraft number (scid): ",7x,i6,/, * "DSN antenna (dss): ",7x,i6,/, * "Downlink band (dlband): ",12x,a1,/, * "Downlink polarization (dlpoln): ",12x,a1,/, * "First time tag:",11x,i4.4,"-",i3.3,"T",2(i2.2,":"),i2.2,".",i6.6,/, * "Last time tag: ",11x,i4.4,"-",i3.3,"T",2(i2.2,":"),i2.2,".",i6.6)') * brecl,nbrec,srate,arecl,sbits,ns_rec,narec,scid,dss,dlband,dlpoln, * yyyy0,ddd0,hh0,mm0,ss0,ffffff0,yyyy1,ddd1,hh1,mm1,ss1,ffffff1 C Calculate number of samples to be included in each average (a multiple of 20) dt = 1.0d+00/(srate*1000.) ns_avg = int(tavg0/dt) if (ns_avg .lt. 20) then ns_avg = 20 endif nr_avg = ns_avg/20 npts = narec/nr_avg tavg1 = dt*ns_avg C Calculate full-range scale factor sf = 4**(sbits-1) write(*,'("Sample time (dt): ",4x,f11.7,/, * "Samples to average (ns_avg): ",5x,i10,/, * "Time average (requested) (tavg0): ",5x,f10.6,/, * " (actual) (tavg1): ",5x,f10.6,/, * "Full scale (sf): ",5x,f10.1)') * dt,ns_avg,tavg0,tavg1,sf write(30,'("Sample time (dt): ",4x,f11.7,/, * "Samples to average (ns_avg): ",5x,i10,/, * "Time average (requested) (tavg0): ",5x,f10.6,/, * " (actual) (tavg1): ",5x,f10.6,/, * "Full scale (sf): ",5x,f10.1)') * dt,ns_avg,tavg0,tavg1,sf C Open file with data samples open ( unit = 33, * file = 'podr_prep.tab', * iostat = ios, * status = 'old') if (ios .ne. 0) then write(*,'("PODR_POWR: error opening input ASCII sample file", * " ios = ",i5)') ios pgm_status = 9 go to 99 end if C Open file for power output open ( unit = 31, * file = 'podr_powr.tab', * access = 'sequential', * iostat = ios, * status = 'unknown') if (ios .ne. 0) then write(*,'("PODR_POWR: error opening output average power file", * " ios = ",i5)') ios pgm_status = 10 go to 99 end if C Calculate average power in samples for this file; C Rewind the file so it can be read again samp_av = 0.0d+00 samp_sq = 0.0d+00 n_av = 0 do i = 1,narec read(33,'(16x,20i7)') (samp(j),j=1,20) do j = 1,20 samp_av = samp_av + samp(j) samp_sq = samp_sq + samp(j)*samp(j) n_av = n_av + 1 end do end do samp_av = samp_av/n_av samp_sq = samp_sq/n_av sd = sqrt(samp_sq - samp_av*samp_av) write(*,'("Average sample value: ",10x,f10.3,"+/-",f7.3,/, * "Number of samples averaged: ",10x,i10)') * samp_av,sd,n_av write(30,'("Average sample value: ",10x,f10.3,"+/-",f7.3,/, * "Number of samples averaged: ",10x,i10)') * samp_av,sd,n_av rewind(33) C Calculate the time averages do i = 1,narec,nr_avg t1 = 3.6d+03*hh0 + 6.0d+01*mm0 + ss0 +1.0d-06*ffffff0 + (1.0*(i-1)/nr_avg + 0.5)*tavg1 power = 0.0d+00 n_av = 0 do n = 1,nr_avg read(33,'(16x,20i7)') (samp(j),j=1,20) do j = 1,20 power = power + (samp(j)-samp_av)**2 n_av = n_av + 1 end do end do power = power/n_av/sf write(31,'(f14.6,1p1e14.6)') t1,power end do C Exit 99 continue write(*,'("Done")') call exit(pgm_status) stop end C ------------------------------------------------------------------------------------------------------------------------ C Program PODR_LOOK C Needs subroutines FFTR and GET_DATE_TIME C Run this program on the output from PODR_PREP to generate a set of average power spectra C Program reads output files from PODR_PREP and calculates average power spectra versus time. C Obtains sampling parameters from PODR_PREP.L pair C Requires: C get_date_time (subroutine) C Exits: C 0 normal C 2 usage (incorrect command line arguments; none are expected) C 3 error in subroutine get_date_time C 4 error opening PODR_TEXT output file C 5 error opening PODR_PREP log file C 10 unexpectedly large number of records C 11 unexpectedly high sampling rate C 12 unexpectedly large number of sample bits C 13 unexpected spacecraft ID C 14 unexpected spacecraft ID and DSS ID combination C 15 unexpected band (should be "X") C Logical unit numbers: C C 28 podr_look.log C 30 podr_prep.log C 31 podr_text.txt C 2024-05-11 RAS: Original (adapted from redr_text.f) integer*4 debug !debug = 0 silent; debug 1, verbose character*1 band !received band (expected to be 'X') integer*4 dssid !identifier for DSN antenna (e.g., 14, 15, 43, 65, etc.) integer*4 i !miscellaneous index character*6 infile !name of original PODR file integer*4 ios !status flag integer*4 irecl !record length for input sample file integer*4 j !miscellaneous index integer*4 k !histogram index derived from sample value character*24 lid !logical identifer for PDS4 product real*4 maxsnrx !maximum X-band SNR from podr_look.log integer*4 narg !number of command line arguments integer*4 nrec !number of records in input sample file character*1 poln !received polarization integer*4 pgm_status/0/ !status upon exit from program integer*4 run_date_time(6) !date and time of program execution integer*4 sbits !bits per sample integer*4 scid !identifier for spacecraft integer*4 srate !sampling rate (ksps) integer*4 srate_len !length of srate (number of digits) character*21 start_date_time !first valid tme tag (yyyy-dddThh:mm:ss.sss) character*21 stop_date_time !last valid tme tag (yyyy-dddThh:mm:ss.sss) integer*4 version(3)/2024,05,11/ !version of program integer*4 iargc !Fortran library function C Initializations debug = 0 C Check command line arguments narg = iargc() if ( debug .eq. 1) write(*,'("narg: ",20x,i19)') narg if (narg .ne. 0) then write(*,'("USAGE: No command line arguments expected")') pgm_status = 2 go to 99 end if C Get today's date and time call get_date_time(run_date_time,ios) if (ios .ne. 0) then pgm_status = 3 go to 99 end if C Open output file for this program open ( unit = 31, * file = 'podr_text.txt', * access = 'sequential', * iostat = ios, * status = 'unknown') if (ios .ne. 0) then write(*,'("PODR_TEXT: error opening PODR_TEXT output file", * " ios = ",i5)') ios pgm_status = 4 go to 99 end if C open log file from podr_prep open ( unit = 30, * file = 'podr_prep.log', * access = 'sequential', * iostat = ios, * status = 'old') if (ios .ne. 0) then write(*,'("PODR_TEXT: error opening PODR_PREP log file", * " ios = ",i5)') ios pgm_status = 5 go to 99 end if C Get name of original PODR file: read(30,'(///,44x,a6)') infile if (debug .eq. 1) write(*,'("PODR tape: ",20x,(a))') infile C Get PODR record length (bytes) and file records and their lengths read(30,'(45x,i5,/,40x,i10)') irecl, nrec if (nrec .gt. 99999) then write(*,'("Unexpectedly large number of records = ",i10)') nrec pgm_status = 10 go to 99 endif if (debug .eq. 1) write(*,'("Record length: ",20x,i6,/, * "Number of records: ",20x,i6)') * irecl,nrec C Get sample rate (kbps), bits-per-sample, and their lengths; convert srate to sps read(30,'(40x,i10,/,40x,i10)') srate, sbits srate = 1000*srate if (srate .gt. 99999) then write(*,'("Unexpectedly high sampling rate = ",i10)') srate pgm_status = 11 go to 99 endif if (sbits .gt. 8) then write(*,'("Unexpectedly large sample bits = ",i10)') sbits pgm_status = 12 go to 99 endif if (debug .eq. 1) write(*,'("Sample rate: ",20x,i8,/, * "Sample bits: ",20x,i8)') * srate,sbits C Get spacecraft ID, receiving station ID, and polarization; correct DSS ID to Parkes read(30,'(48x,i2,/,48x,i2,/,49x,a1,/,49x,a1)') scid,dssid,band,poln if (band .eq. "x") then band = 'X' else if (band .ne. "X") then write(*,'("Unexpected band = ",a1)') band pgm_status = 15 go to 99 endif dssid = 49 if (scid .ne. 31 .and. scid .ne. 32 .and. scid .ne. 99) then write(*,'("Unexpected spacecraft ID = ",i10)') scid pgm_status = 13 go to 99 endif if ((scid .ne. 32) .or. (dssid .ne. 49)) then write(*,'("Unexpected scid = ",i2," and dssid = ",i2," combination")') scid,dssid pgm_status = 14 go to 99 endif if (debug .eq. 1) write(*,'("Spacecraft ID: ",20x,i5,/, * "DSN antenna number: ",20x,i5,/, * "Band: ",24x,a1,/, * "Polarization: ",24x,a1)') * scid,dssid,band,poln C Get start and end date/time: read(30,'(//////,26x,a21)') start_date_time read(30,'(26x,a21)') stop_date_time if (debug .eq. 1) write(*,'("First valid time tag: ",a21,/, * "Last valid time tag: ",a21)') * start_date_time,stop_date_time close(30) C Get max SNR from podr_look.log open ( unit = 28, * file = 'podr_look.log', * iostat = ios, * status = 'old') read(28,'(///////////////////////////////////,44x,f5.2)') maxsnrx close(28) if (debug .eq. 1) write(*,'("Max SNR:",30x,f7.2)') maxsnrx C Create product ID for the raw data product write(lid(01:24),'("vg2u_",i2.2,"xr_",a4,a3,"t",3a2)') dssid, * start_date_time(01:04),start_date_time(06:08), * start_date_time(10:11),start_date_time(13:14),start_date_time(16:17) if (debug .eq. 1) write(*,'("LID:",17x,(a))') lid C Document the run write(*,'("PODR_TEXT: Create text file of PODR information:",/, * "Version of: ",26x,i4.4,"-",i2.2,"-",i2.2,/, * "Today is: ",18x,i4.4,2("-",i2.2),"T",i2.2, * 2(":",i2.2)," local time")') * version,run_date_time write(31,'("PODR_TEXT: Create text file of PODR information:",/, * "Version of: ",25x,i4.4,"-",i2.2,"-",i2.2,/, * "Today is: ",17x,i4.4,2("-",i2.2),"T",i2.2, * 2(":",i2.2)," local time")') * version,run_date_time C Write footprint information to output file write(31,'("PDS4 Product LID:",6x,(a))') lid(1:24) if (infile(1:2) .eq. "UL") then write(31,'("PODR File Name: ",25x,(a))') infile(01:06) else write(31,'("UNEXPECTED PODR FILE NAME",(a))') infile endif write(31,'("Record Length:",29x,i4)') irecl write(31,'("Number of Records:",24x,i5)') nrec write(31,'("Sample Rate:",30x,i5)') srate write(31,'("Sample Bits:",33x,i2)') sbits write(31,'("Spacecraft ID:",31x,i2)') scid if (scid .eq. 32) write(31,'("Spacecraft Name:",22x,"Voyager 2")') write(31,'("DSN Receive Antenna:",25x,i2)') dssid write(31,'("Downlink Band: ",24x,a1)') band write(31,'("Downlink Polarization:",24x,a1)') poln write(31,'("Start Time:",15x,a21)') start_date_time write(31,'("Stop Time: ",15x,a21)') stop_date_time if (infile(1:2) .eq. "UL") then write(31,'("Target Name:",29x,"Uranus")') else write(31,'("ERROR IN FORMING TARGET NAME")') endif write(31,'("Instrument Name:",8x,"Radio Science Subsystem")') if (infile(1:2) .eq. "UL") then write(31,'("Product Type:",23x,"PODR Browse")') else write(31,'("UNABLE TO FORM PRODUCT TYPE")') endif write(31,'("Producer ID:"28x,"PDS/RMS")') go to 99 C Exit 99 continue close(31) write(*,'("Done")') call exit(pgm_status) stop end C ------------------------------------------------------------------------------------------------------------------------ C Subroutine DOY_MMDD subroutine doy_mmdd (mode, year, doy, mm, dd, status) C Converts day-of year to months and days within month (mode = 1) C Converts months and days to day-of-year (mode = -1) C Valid for years 1901 - 2099 inclusive C Exit status C 0 normal exit C 1 year is outside valid range C 2 illegal value for mode integer*4 mode integer*4 year integer*4 doy integer*4 mm integer*4 dd integer*4 status integer*4 i integer*4 month(12) /31,28,31,30,31,30,31,31,30,31,30,31/ C Initializations status = 0 if ((year .le. 1900) .or. (year .ge. 2100)) then status = 1 go to 90 end if month(2) = 28 if (mod(year,4) .eq. 0) month(2) = 29 C Enter loops for calculation if (mode .eq. 1) then dd = doy i = 1 do while (dd .gt. month(i)) dd = dd - month(i) i = i + 1 end do mm = i else if (mode. eq. -1) then doy = 0 i = 1 do while (i .lt. mm) doy = doy + month(i) i = i + 1 end do doy = doy + dd else status = 2 go to 90 end if C Exit 90 continue return end C ------------------------------------------------------------------------------------------------------------------------ C Subroutine FFTR subroutine fftr(dd,nn,ii,ier) C Forward and inverse FFT for real time samples. C C dd is a real array of length 2*nn or, equivalently, a complex array of length nn. C nn is a power of 2 (nn <= 4096) C ii sets direction of transform C ii= 1 replaces dd by its discrete Fourier transform (positive frequencies only) C ii=-1 replaces dd by its inverse discrete Fourier transform C dd is a complex array of length nn or, equivalently, a real array of length 2*nn. C nn must be an integer power of 2. C Amplitude scaling (by sf) ensures that the sum of squared time samples is twice the sum C of squared frequency values (in both directions) and that a forward transform followed by C an inverse tranform recovers the original time samples. C Exit status: C ier = 0 normal C 2 unexpected value for nn C 3 unexpected value for ii C 2021-02-07 RAS: Based on FOUR1 in Numerical Recipes (Fortran), Press et al. (1989) C Added value checks, status return, and scale factor double precision wr, wi, wpr, wpi, wtemp, theta dimension dd(2*nn) ier = 0 if ((nn .ne. 2) .and. (nn .ne. 4) .and. (nn .ne. 8) .and. * (nn .ne. 16) .and. (nn .ne. 32) .and. (nn .ne. 64) .and. * (nn .ne. 128) .and. (nn .ne. 256) .and. (nn .ne. 512) .and. * (nn .ne. 1024) .and. (nn .ne. 2048) .and. (nn .ne. 4096)) then ier = 2 return endif j = 1 n = 2*nn if (ii .eq. 1) then sf = 1./sqrt(float(n)) else if (ii .eq. -1) then sf = 2./sqrt(float(n)) else ier = 3 return endif do 11 i = 1,n,2 if (j .gt. i) then tempr = dd(j) tempi = dd(j+1) dd(j) = dd(i) dd(j+1) = dd(i+1) dd(i) = tempr dd(i+1) = tempi endif m = n/2 1 if ((m .ge. 2) .and. (j .gt. m)) then j = j - m m = m/2 go to 1 endif j = j + m 11 continue mmax = 2 C Beginning of Danielson-Lanczos section of the routine. C Outer loop executed log2(nn) times. 2 if (n .gt. mmax) then istep = 2*mmax theta = 6.28318530717959d+00/(ii*mmax) wpr = -2.0d+00*dsin(0.5d+00*theta)**2 wpi = dsin(theta) wr = 1.0d+00 wi = 0.0d+00 do 13 m = 1,mmax,2 do 12 i= m,n,istep j = i + mmax tempr = sngl(wr)*dd(j) - sngl(wi)*dd(j+1) tempi = sngl(wr)*dd(j+1) + sngl(wi)*dd(j) dd(j) = dd(i) - tempr dd(j+1) = dd(i+1) - tempi dd(i) = dd(i) + tempr dd(i+1) = dd(i+1) + tempi 12 continue wtemp = wr wr = wr*wpr - wi*wpi + wr wi = wi*wpr + wtemp*wpi + wi 13 continue mmax = istep go to 2 endif do 14 i = 1,n dd(i) = sf*dd(i) 14 continue return end C ------------------------------------------------------------------------------------------------------------------------ C Subroutine GET_DATE_TIME subroutine get_date_time(run_date_time,ier) C Subroutine calls system date/time and converts to PDS order C run_date_time(1) = YYYY C run_date_time(2) = MM C run_date_time(3) = DD C run_date_time(4) = hh C run_date_time(5) = mm C run_date_time(6) = ss C 2016-07-11 RAS: original C 2024-08-25 RAS: rewritten to use system call to 'date' when fdate was no longer available C through the HPC gfortran compiler C Declarations character*24 date_time !date and time as returned by system call fdate character*3 day !day of week (ASCII) integer*4 i !miscellaneous index character*3 month !month (ASCII) integer*4 ier !error flag integer*4 ios !error flag when opening unit 11 integer*4 run_date_time(6) !date and time in PDS order C Initializations ier = 0 call system('date > date_time.txt') open ( unit = 11, * file = 'date_time.txt', * access = 'SEQUENTIAL', * iostat = ios, * status = 'OLD') C Extraction read(11,'(a3,1x,a3,4(1x,i2),5x,i4)') day,month, * run_date_time(3), * run_date_time(4), * run_date_time(5), * run_date_time(6), * run_date_time(1) C Month conversion if (month .eq. 'Jan' ) then run_date_time(2) = 01 else if (month .eq. 'Feb' ) then run_date_time(2) = 02 else if (month .eq. 'Mar' ) then run_date_time(2) = 03 else if (month .eq. 'Apr' ) then run_date_time(2) = 04 else if (month .eq. 'May' ) then run_date_time(2) = 05 else if (month .eq. 'Jun' ) then run_date_time(2) = 06 else if (month .eq. 'Jul' ) then run_date_time(2) = 07 else if (month .eq. 'Aug' ) then run_date_time(2) = 08 else if (month .eq. 'Sep' ) then run_date_time(2) = 09 else if (month .eq. 'Oct' ) then run_date_time(2) = 10 else if (month .eq. 'Nov' ) then run_date_time(2) = 11 else if (month .eq. 'Dec' ) then run_date_time(2) = 12 else write(*,'("GET_DATE_TIME: Unexpected month from fdate = ",(a))') * month ier = 1 end if C Exit close(11) call system('rm date_time.txt') return end