From fa080de7afc95aa1c19a6e6fc0e0708ced2eadc4 Mon Sep 17 00:00:00 2001 From: Joseph Hunkeler Date: Wed, 8 Jul 2015 20:46:52 -0400 Subject: Initial commit --- noao/mtlocal/cyber/README | 4 + noao/mtlocal/cyber/cyber.h | 187 ++++++++++++ noao/mtlocal/cyber/cykeywords.x | 102 +++++++ noao/mtlocal/cyber/cyrbits.x | 371 ++++++++++++++++++++++++ noao/mtlocal/cyber/cyrheader.x | 120 ++++++++ noao/mtlocal/cyber/cyrimage.x | 166 +++++++++++ noao/mtlocal/cyber/mkpkg | 22 ++ noao/mtlocal/cyber/pow.inc | 86 ++++++ noao/mtlocal/cyber/powd.inc | 128 ++++++++ noao/mtlocal/cyber/rdumpf.x | 283 ++++++++++++++++++ noao/mtlocal/cyber/rpft.x | 211 ++++++++++++++ noao/mtlocal/cyber/rrcopy/README | 2 + noao/mtlocal/cyber/rrcopy/Revisions | 15 + noao/mtlocal/cyber/rrcopy/mkpkg | 15 + noao/mtlocal/cyber/rrcopy/rcrbits.x | 279 ++++++++++++++++++ noao/mtlocal/cyber/rrcopy/rcrheader.x | 119 ++++++++ noao/mtlocal/cyber/rrcopy/rcrimage.x | 173 +++++++++++ noao/mtlocal/cyber/rrcopy/rrcopy.h | 41 +++ noao/mtlocal/cyber/rrcopy/rrcopy.x | 212 ++++++++++++++ noao/mtlocal/cyber/rrcopy/semicode.doc | 310 ++++++++++++++++++++ noao/mtlocal/cyber/rrcopy/t_rrcopy.x | 147 ++++++++++ noao/mtlocal/cyber/t_ldumpf.x | 220 ++++++++++++++ noao/mtlocal/cyber/t_rdumpf.x | 162 +++++++++++ noao/mtlocal/cyber/t_ridsfile.x | 516 +++++++++++++++++++++++++++++++++ noao/mtlocal/cyber/t_ridsout.x | 386 ++++++++++++++++++++++++ 25 files changed, 4277 insertions(+) create mode 100644 noao/mtlocal/cyber/README create mode 100644 noao/mtlocal/cyber/cyber.h create mode 100644 noao/mtlocal/cyber/cykeywords.x create mode 100644 noao/mtlocal/cyber/cyrbits.x create mode 100644 noao/mtlocal/cyber/cyrheader.x create mode 100644 noao/mtlocal/cyber/cyrimage.x create mode 100644 noao/mtlocal/cyber/mkpkg create mode 100644 noao/mtlocal/cyber/pow.inc create mode 100644 noao/mtlocal/cyber/powd.inc create mode 100644 noao/mtlocal/cyber/rdumpf.x create mode 100644 noao/mtlocal/cyber/rpft.x create mode 100644 noao/mtlocal/cyber/rrcopy/README create mode 100644 noao/mtlocal/cyber/rrcopy/Revisions create mode 100644 noao/mtlocal/cyber/rrcopy/mkpkg create mode 100644 noao/mtlocal/cyber/rrcopy/rcrbits.x create mode 100644 noao/mtlocal/cyber/rrcopy/rcrheader.x create mode 100644 noao/mtlocal/cyber/rrcopy/rcrimage.x create mode 100644 noao/mtlocal/cyber/rrcopy/rrcopy.h create mode 100644 noao/mtlocal/cyber/rrcopy/rrcopy.x create mode 100644 noao/mtlocal/cyber/rrcopy/semicode.doc create mode 100644 noao/mtlocal/cyber/rrcopy/t_rrcopy.x create mode 100644 noao/mtlocal/cyber/t_ldumpf.x create mode 100644 noao/mtlocal/cyber/t_rdumpf.x create mode 100644 noao/mtlocal/cyber/t_ridsfile.x create mode 100644 noao/mtlocal/cyber/t_ridsout.x (limited to 'noao/mtlocal/cyber') diff --git a/noao/mtlocal/cyber/README b/noao/mtlocal/cyber/README new file mode 100644 index 00000000..a5c88b40 --- /dev/null +++ b/noao/mtlocal/cyber/README @@ -0,0 +1,4 @@ +This directory ("dataio$cyber/") contains code for four Cyber tape readers: +LDUMPF, RDUMPF, RIDSOUT, RIDSFILE. The executable file is x_cyber; the +definition file cyber.h is common to all four readers. +The subdirectory "rrcopy" contains code for the RRCOPY reader. diff --git a/noao/mtlocal/cyber/cyber.h b/noao/mtlocal/cyber/cyber.h new file mode 100644 index 00000000..4cd4aa39 --- /dev/null +++ b/noao/mtlocal/cyber/cyber.h @@ -0,0 +1,187 @@ +# Fix length of user area +define LEN_USER_AREA 2880 + +# Offsets to elements of the IPPS raster header (RDUMPF) + +define DATA_TYPE_OFFSET 33 # Offset to data_type (nbpp) +define NCOLS_OFFSET 35 # Offset to ncols (nppr) +define NWORDS_OFFSET 37 # Offet to nwords_per_row +define NROWS_OFFSET 41 # Offset to nrows +define FIRST_PRU_OFFSET 43 # Offset to 1st pru of raster +define MIN_OFFSET 63 # Offset to data min +define MAX_OFFSET 65 # Offset to data max +define EOR_OFFSET 89 # Offset to terminating pru + +# Bit-offsets to fields in the Permanent File Table (LDUMPF) + +define PFN_OFFSET (27 * 60 + 1 - 6) # Name: Left justified +define PF_ID_OFFSET (26 * 60 + 1 - 6) # ID: Right justified +define CY_OFFSET (31 * 60 + 1 - 12) # Cycle number +define CREATE_OFFSET (32 * 60 + 1 - 18) # Creation date +define ATTACH_OFFSET (33 * 60 + 1 - 18) # Date of last attach +define ALTER_OFFSET (34 * 60 + 1 - 18) # Date of last alteration +define NCHARS_OFFSET (8 * 60 + 7) # Nchars in PF name + +# The IPPS raster descriptor structure DT (RDUMPF): + +define LEN_DT 10 + SZ_IPPS_ID + 1 + +define BITS_PIXEL Memi[$1] +define PRU_EOR Memi[$1+1] +define WRDS_PER_ROW Memi[$1+2] +define NPRU_ROW Memi[$1+3] +define PRU_ROW_ONE Memi[$1+4] +define NCOLS Memi[$1+5] +define NROWS Memi[$1+6] +define DATA_MIN Memr[P2R($1+7)] +define DATA_MAX Memr[P2R($1+8)] +define IPPS_ID Memc[P2C($1+10)] + + +# The DUMPF tape descriptor structure DMP (LDUMPF): + +define LEN_DMP 15 + SZ_PFN + SZ_PF_ID + 2 + +define CY Memi[$1] +define M_CREATE Memi[$1+1] +define D_CREATE Memi[$1+2] +define Y_CREATE Memi[$1+3] +define M_ALTER Memi[$1+4] +define D_ALTER Memi[$1+5] +define Y_ALTER Memi[$1+6] +define M_ATTACH Memi[$1+7] +define D_ATTACH Memi[$1+8] +define Y_ATTACH Memi[$1+9] +define NCHARS_PFN Memi[$1+10] +define PFN Memc[P2C($1+15)] +define ID Memc[P2C($1+50)] + +# Bit-offsets to various IDSFILE header words are defined: + +define TAPE_OFFSET ((512 + (15 - 1)) * 64 + 1) +define SCAN_OFFSET ((512 + (1 - 1)) * 64 + 1) +define ITM_OFFSET ((512 + (2 - 1)) * 64 + 1) +define NP1_OFFSET ((512 + (5 - 1)) * 64 + 1) +define NP2_OFFSET ((512 + (6 - 1)) * 64 + 1) +define BEAM_OFFSET ((512 + (7 - 1)) * 64 + 1) +define COMPANION_OFFSET ((512 + (64 - 1)) * 64 + 1) +define OLD_OFFSET ((512 + (64 - 1)) * 64 + 1) +define SMODE_OFFSET ((512 + (10 - 1)) * 64 + 1) +define UT_OFFSET ((512 + (11 - 1)) * 64 + 1) +define ST_OFFSET ((512 + (12 - 1)) * 64 + 1) +define DF_OFFSET ((512 + (16 - 1)) * 64 + 1) +define SM_OFFSET ((512 + (17 - 1)) * 64 + 1) +define QF_OFFSET ((512 + (18 - 1)) * 64 + 1) +define DC_OFFSET ((512 + (19 - 1)) * 64 + 1) +define QD_OFFSET ((512 + (20 - 1)) * 64 + 1) +define EX_OFFSET ((512 + (21 - 1)) * 64 + 1) +define BS_OFFSET ((512 + (22 - 1)) * 64 + 1) +define CA_OFFSET ((512 + (23 - 1)) * 64 + 1) +define CO_OFFSET ((512 + (24 - 1)) * 64 + 1) +define HA_OFFSET ((512 + 26) * 64 + 1) +define AIR_OFFSET ((512 + 27) * 64 + 1) +define RA_OFFSET ((512 + 12) * 64 + 1) +define DEC_OFFSET ((512 + 13) * 64 + 1) +define LAM_OFFSET ((512 + 2) * 64 + 1) +define DEL_OFFSET ((512 + 3) * 64 + 1) +define OFLAG_OFFSET ((512 + (26 - 1)) * 64 + 1) +define ALPHA1_OFFSET ((512 + (25 - 1)) * 64 + 55) +define ALPHA2_OFFSET (ALPHA1_OFFSET - 6) +define ALPHA3_OFFSET (ALPHA2_OFFSET - 6) +define IDS_ID_OFFSET (512 + 28) +define COEFF_OFFSET (512 + 37) +define IDSO_REC_OFF 1 +define IDSO_TAPE_OFF 961 + +# Definition of the control parameter descriptor structure (RIDSFILE, RIDSOUT) + +define LEN_CP 5 + SZ_FNAME + SZ_LINE + 2 + +define MAKE_IMAGE Memi[$1] +define PRINT_PIXELS Memi[$1+1] +define LONG_HEADER Memi[$1+2] +define DATA_TYPE Memi[$1+3] +define IRAF_FILE Memc[P2C($1+5)] +define REC_NUMBERS Memc[P2C($1+5+SZ_FNAME)] + +# The IDSFILE descriptor structure IDS (RIDSOUT, RIDSFILE): + +define LEN_IDS 35 + 5 + SZ_IDS_ID + 1 + +define HA Memd[P2D($1)] +define AIRMASS Memd[P2D($1+2)] +define RA Memd[P2D($1+4)] +define DEC Memd[P2D($1+6)] +define LAMBDA0 Memd[P2D($1+8)] +define DELTA_LAMBDA Memd[P2D($1+10)] +define RECORD_NUMBER Memi[$1+12] +define ITM Memi[$1+13] +define NP1 Memi[$1+14] +define NP2 Memi[$1+15] +define BEAM_NUMBER Memi[$1+16] +define COMPANION_RECORD Memi[$1+17] +define SMODE Memi[$1+18] +define UT Memi[$1+19] +define ST Memi[$1+20] +define DF_FLAG Memi[$1+21] +define SM_FLAG Memi[$1+22] +define QF_FLAG Memi[$1+23] +define DC_FLAG Memi[$1+24] +define QD_FLAG Memi[$1+25] +define EX_FLAG Memi[$1+26] +define BS_FLAG Memi[$1+27] +define CA_FLAG Memi[$1+28] +define CO_FLAG Memi[$1+29] +define OFLAG Memi[$1+30] +define COEFF Memi[$1+31] +define ALPHA_ID Memc[P2C($1+35)] +define LABEL Memc[P2C($1+40)] + + +# Definitions for the Cyber DUMPF tape reading programs RDUMPF and LDUMPF + +define NBITS_CHAR (NBITS_BYTE * SZB_CHAR) # Number of bits per char +define NBITS_CYBER_WORD 60 # Number of bits per Cyber word +define LEN_PRU 64 # Number of words per Cyber pru +define NCHARS_PRU (64 * 60 / NBITS_CHAR) # Nchars per PRU +define NBITS_PRU 3840 # Number of bits per Cyber pru +define NCHARS_NOISE (48 / NBITS_CHAR) # Nchars in a Cyber noise record +define NBITS_EOR_MARK 48 # Number of bits per eor marker +define LEN_HEADER 64 # Number of words per header +define SZ_HEADER ((64 * 60) / NBITS_CHAR) # Nchars in IPPS header +define SZ_TAPE_BLK ((512 * 60) / NBITS_CHAR) # Size in chars of tape block +define SZ_TAPE_BUFFER (SZ_TAPE_BLK + 60) # Size of tape buffer for read +define LEN_PFT 48 # Size of Permanent File Table +define SZ_IPPS_ID 127 # Max number of characters in ID +define MAX_RANGES 100 +define NOT_SET 0 # Flag for data_type not set +define BLANK 0.0 # Temporary value for blanks +define NBYTES_WORD 5 # 5 12-bit bytes per Cyber word +define NINT_CYBER_WRD 2 +define LEN_CYBER_READ (4 * 65) # Number of Cyber words read at once +define SZ_PFT ((48 * 60) / NBITS_CHAR) # Chars in Perm file table +define SZ_PFN 40 # Max characters in PF Name +define SZ_PF_ID 9 # Max characters in PF ID +define NBITS_DATE 18 # Dates are written in 18 bits +define NBITS_CY 12 # Cycle # written in 12 bits +define NBITS_DC 6 # Nbits display code character +define NCHARS_WORD 10 # Number of display code + # characters per cyber word + +define LEN_INDEX (5 * LEN_PRU) +define LEN_USER_INDEX (2 * LEN_PRU) +define LEN_IDS_RECORD (9 * LEN_PRU) +define NPIX_IDS_RECORD 1024 +define SZ_IDS_ID 64 +define NCHAR_ALPHA 3 +define START_OF_IDSFILE 6 # First PRU of IDSFILE after index +define NBITS_LRN 18 +define NBITS_HRN 18 +define NBITS_NPRU 24 +define LRN_OFFSET 19 # Offset in index to lrn +define HRN_OFFSET 1 # Offset in index to hrn +define NPRU_OFFSET 37 # Offset in index to pru ordinal +define NBITS_DC 6 # Number of bits in a display code char +define MAX_COEFF 28 # Maximum n_coeff if DF is set +define NPIX_LINE 8 # Npixels per line of text (IDSOUT) +define NLINES_PIXELS 128 # Nlines of text containing pixels diff --git a/noao/mtlocal/cyber/cykeywords.x b/noao/mtlocal/cyber/cykeywords.x new file mode 100644 index 00000000..f6a2b194 --- /dev/null +++ b/noao/mtlocal/cyber/cykeywords.x @@ -0,0 +1,102 @@ +include +include +include "cyber.h" + +define LEN_KEYWORD 8 +define UNKNOWN Memc[($1+IMU-1)*SZ_STRUCT + 1] + +# CY_STORE_KEYWORDS -- store IDS specific keywords in the IRAF image header. + +procedure cy_store_keywords (ids, im) + +pointer ids # Pointer to program data structure +pointer im # Pointer to image + +char keyword[LEN_KEYWORD] +int fd, i, nrp +int stropen() +real value + +begin + # Open image user area as a file + fd = stropen (UNKNOWN(im), LEN_USER_AREA - 1, NEW_FILE) + + # FITS keyword are formatted and appended to the image user + # area with the addcard procedures. + + call addcard_i (fd, "RECORD", RECORD_NUMBER(ids), "IDS record") + + if (SMODE(ids) != 0) { + call addcard_i (fd, "COMPANIO", COMPANION_RECORD(ids), + "Companion record") + + } + + call addcard_i (fd, "EXPOSURE", ITM(ids), "Exposure time (seconds)") + + call addcard_i (fd, "OFLAG", OFLAG(ids), "Object flag") + + call addcard_i (fd, "BEAM-NUM", BEAM_NUMBER(ids), "Beam number") + + nrp = NDIGITS_RP + value = real (LAMBDA0(ids)) + call addcard_r (fd, "W0", value, "Starting wavelength", nrp) + + value = real (DELTA_LAMBDA(ids)) + call addcard_r (fd, "WPC", value, "Wavelength per channel", nrp) + + call addcard_i (fd, "NP1", NP1(ids), "Left plot limit") + + call addcard_i (fd, "NP2", NP2(ids), "Right plot limit") + + if (IS_INDEFI (UT(ids))) + value = INDEFR + else + value = real (UT(ids) / 3600.) + call addcard_time (fd, "UT", value, "Universal time") + + if (IS_INDEFI (ST(ids))) + value = INDEFR + else + value = real (ST(ids)/ 3600.) + call addcard_time (fd, "ST", value, "Sidereal time") + + value = real (RA(ids)) + call addcard_time (fd, "RA", value, "Right ascension") + + value = real (DEC(ids)) + call addcard_time (fd, "DEC", value, "Declination") + + value = real (HA(ids)) + call addcard_time (fd, "HA", value, "Hour angle") + + value = real (AIRMASS(ids)) + call addcard_r (fd, "AIRMASS", value, "Airmass", nrp) + + # The 9 reduction flags + call addcard_i (fd, "DF-FLAG", DF_FLAG(ids), "Dispersion flag") + call addcard_i (fd, "SM-FLAG", SM_FLAG(ids), "Smoothing flag") + call addcard_i (fd, "QF-FLAG", QF_FLAG(ids), "Quartz fit flag") + call addcard_i (fd, "DC-FLAG", DC_FLAG(ids), "Dispersion corrected") + call addcard_i (fd, "QD-FLAG", QD_FLAG(ids), "Quartz division flag") + call addcard_i (fd, "EX-FLAG", EX_FLAG(ids), "Extinction flag") + call addcard_i (fd, "BS-FLAG", BS_FLAG(ids), "Beam switch flag") + call addcard_i (fd, "CA-FLAG", CA_FLAG(ids), "Calibration flag") + call addcard_i (fd, "CO-FLAG", CO_FLAG(ids), "Coincidence flag") + + # The df coeffecients are written out in the case where the df + # flag is set, and the first coefficient is nonzero. (The later + # condition is a test for IDSOUT data, where the df coeffecients + # have been applied but not stored in the header.) + + if (DF_FLAG(ids) != -1 && COEFF(ids) != 0) { + call strcpy ("DF", keyword, LEN_KEYWORD) + do i = 1, DF_FLAG(ids) { + call sprintf (keyword[3], LEN_KEYWORD, "%s") + call pargi (i) + call addcard_d (fd, keyword, Memd[COEFF(ids)+i-1], "", nrp) + } + } + + call close (fd) +end diff --git a/noao/mtlocal/cyber/cyrbits.x b/noao/mtlocal/cyber/cyrbits.x new file mode 100644 index 00000000..9f4df03f --- /dev/null +++ b/noao/mtlocal/cyber/cyrbits.x @@ -0,0 +1,371 @@ +include +include +include +include "cyber.h" + +# UNPK_12 -- Unpack 12-bit unsigned integers from an array containing +# one sixty bit cyber word in each pair of array elements. +# Each output word contains successive 12-bit pixels from the +# input array. +# Pixels are unpacked starting with element "first_word" of the input array. +# Each cyber 60-bit word contains 5 packed 12-bit pixels, the first pixel +# in the highest 12 bits. The input array contains one cyber word per two +# array elements; the cyber word occupies the lower 60 bits of each pair +# of array values. + +procedure unpk_12 (input, first_word, output, npix_unpk) + +int input[ARB] # +int first_word # +int output[npix_unpk] # +int npix_unpk # + +int n, nn, i, offset[5], off +int npix_word, ncyber_words, index +data (offset[i], i=1, 5) /15, 27, 39, 51, 63/ +int bitupk() + +begin + npix_word = 5 + if (mod (npix_unpk, npix_word) == 0) + ncyber_words = (npix_unpk) / npix_word + else + call error (0, "Incorrect number of pixels to be unpacked") + index = 1 + + i = 1 + for (n = first_word; i <= npix_unpk; n = n + 2) { + do nn = 1, npix_word { + off = (n + 1) * NBITS_INT - offset[nn] + output[i] = bitupk (input, off, 12) + if (output[i] == 7777B) + output[i] = BLANK + i = i + 1 + } + } +end + + +# UNPK-20 -- Unpack 20-bit signed integers from an array containing +# one 60-bit Cyber word per pair of array elements. Each output +# word contains successive 20-bit pixels from the input array. Pixels +# are unpacked starting with array element "first_word". Conversion +# from one's complement to two's complement is performed. Each Cyber +# word contains 3 packed 20-bit pixels, the first pixel in the highest +# 20 bits. + +procedure unpk_20 (input, first_word, output, npix_unpk) + +int input[ARB] +int output[npix_unpk], npix_unpk, first_word +int n, i, nn, off +int npix_word, ncyber_words, pix_val, offset[3] +data (offset[i], i=1, 3) /23, 43, 63/ +int bitupk() + +begin + npix_word = 3 + if (mod (npix_unpk, npix_word) == 0) + ncyber_words = npix_unpk / npix_word + else + call error (0, "Incorrect number of pixels to be unpacked") + + i = 1 + for (n = first_word; i <= npix_unpk; n = n + 2) { + do nn = 1, npix_word { + off = (n + 1) * NBITS_INT - offset[nn] + pix_val = bitupk (input, off, 20) + if (pix_val == 3777777B) + pix_val = BLANK + else if (and (pix_val, 2000000B) != 0) + # negative pixel + pix_val = -and (3777777B, not(pix_val)) + output[i] = pix_val + i = i + 1 + } + } +end + + +# UNPK_60R -- Unpack Cyber 60-bit floating point numbers from an array +# containing one 60-bit word per pair of array elements. +# The 30 most significant bits from each 60-bit word are +# unpacked and then reconstructed as a floating point number with +# with REPACK_FP. The extracted bits include an 18-bit mantissa, +# 11-bit exponent and a sign bit. This routine is used for getting +# the min and max values from the header; no 60-bit IPPS pixels are +# expected. + +procedure unpk_60r (input, first_word, fp_value, nwords) + +int input[ARB] +real fp_value[nwords] +int first_word, nwords, n, i +pointer int_buf, sp +int bitupk() + +begin + # Allocate space on stack + call smark (sp) + call salloc (int_buf, nwords, TY_INT) + + i = 1 + for (n = first_word; i <= nwords; n = n + 2) { + Memi[int_buf + i - 1] = bitupk (input[n], 31, 30) + i = i + 1 + } + + call repack_fp (Memi[int_buf], fp_value, nwords) + call sfree (sp) +end + + +# UNPK_60I -- Unpack 60-bit integers from an array containing one Cyber +# word per each NINT_CYBER_WRD elements. Each word +# of output contains only the lower 32 bits of each input word, as this +# procedure is called only for getting the reduction flags from the +# IDSFILE header. + +procedure unpk_60i (input, initial_bit_offset, output, nwords) + +char input[ARB] +int output[nwords] +int initial_bit_offset, nwords, bit_offset, n +int bitupk() +errchk bitupk + +begin + bit_offset = initial_bit_offset + + do n = 1, nwords { + output[n] = bitupk (input, bit_offset, NBITS_INT) + if (and (output[n], 2000000000B) != 0) + # negative value + output[n] = -not(output[n]) + bit_offset = bit_offset + (NINT_CYBER_WRD * NBITS_INT) + } +end + + +# CONVERT_60BIT_FP -- returns a floating point number equivalent to the Cyber +# 60-bit number input. The full 48-bit mantissa and 11-bit exponent +# is used in reconstructing the floating point value. + +double procedure convert_60bit_fp (cyber_word, bit_offset) + +int cyber_word[ARB], bit_offset +int temp1, temp2 +double float_value +int exp, lower_mantissa, upper_mantissa, i +int bitupk(), and(), not() +double tbl[255] +include "powd.inc" +errchk bitupk + +begin + # Extract cyber word in question into temp array + temp1 = bitupk (cyber_word, bit_offset, 30) + temp2 = bitupk (cyber_word, bit_offset + 30, 30) + + # Check "bit59" and complement all bits if it is set + if (bitupk (temp2, 30, 1) != 0) { + temp1 = not (temp1) + temp2 = not (temp2) + lower_mantissa = -and (temp1, 7777777777B) + upper_mantissa = -and (temp2, 777777B) + } else { + lower_mantissa = temp1 + upper_mantissa = and (temp2, 777777B) + } + + # Extract and interpret exponent; remove Cyber bias of 2000B and + # convert to two's complement if negative number + exp = bitupk (temp2, 19, 11) + if (exp > 1777B) + # "bit58" is set, positive exponent + exp = exp - 2000B + else + # negative exponent + exp = exp - 1777B + + # Reconstruct the floating point number. 30 is added to the + # exponent for the upper mantissa. The 129 is to register the data + # array matrix: tbl[1] = 2 ** -128 ==> 2 ** n = tbl[n + 129] + # float_value = mantissa * 2 ** (exp + 129) + # + # float_value = (lower_mantissa) * 2 ** (exp + 129) + + # (upper_mantissa) * 2 ** (exp + 30 + 129) + + float_value = double (lower_mantissa) * tbl[exp + 129] + + double (upper_mantissa) * tbl[exp + 30 + 129] + + return (float_value) +end + + +# UNPK_30 -- unpack Cyber 30-bit floating point numbers from an array +# containing one 60-bit Cyber word in each pair of array elements. Each +# 30-bit pixel is unpacked from this array; procedure REPACK_FP is called +# to reconstruct the floating point number. Pixels are unpacked starting +# with array element "first_word". Each Cyber word contains 2 30-bit +# pixels, the first pixel in the higher 30 bits. + +procedure unpk_30 (input, first_word, fp_value, npix) + +int input[ARB] +int first_word, npix +real fp_value[npix] +pointer int_buf, sp +int n, i, off, offset[2] +int bitupk() +data (offset[i], i = 1, 2) /33, 63/ +errchk bitupk + +begin + # Allocate buffer space, allowing for maximum of 1 extraneous pixel + call smark (sp) + call salloc (int_buf, npix + 1, TY_INT) + + i = 1 + for (n = first_word; i <= npix; n = n + 2) { + off = (n + 1) * NBITS_INT - offset[1] + Memi[int_buf + i - 1] = bitupk (input, off, 30) + off = (n + 1) * NBITS_INT - offset[2] + Memi[int_buf + i] = bitupk (input, off, 30) + i = i + 2 + } + + call repack_fp (Memi[int_buf], fp_value, npix) + call sfree (sp) +end + + +# UNPK_ID -- Unpacks ID string from input array, which contains one Cyber +# word per two array elements. The word_offset equals the number of Cyber +# words to skip before beginning to unpack. If the character string +# begins in word one of "input", word_offset = 0. The IPPS ID string +# is written in 7-bit ASCII, with eight characters per Cyber word. The lowest +# 4 bits of each 60-bit word is unused. The highest 7 bits of the first Cyber +# word contain the character count. + +procedure unpk_id (input, word_offset, output) + +int input[ARB] +int word_offset +char output[SZ_IPPS_ID] +int nbits, nchar_offset, id_offset, nchars, n +int nchars_word, ncyber_words, nn, index +int bitupk() + +begin + nbits = 7 + nchar_offset = (word_offset * NBITS_INT * NINT_CYBER_WRD) + + NBITS_CYBER_WORD - 6 + nchars = bitupk (input, nchar_offset, nbits) + ncyber_words = (nchars + 6) / 7 + index = 1 + + for (n = 1; n <= ncyber_words; n = n + 1) { + if (n == 1) { + nchars_word = 7 + id_offset = nchar_offset - 7 + } else { + nchars_word = 8 + id_offset = nchar_offset + ((n-1) * NBITS_INT * NINT_CYBER_WRD) + } + do nn = 1, nchars_word { + output[index] = bitupk (input, id_offset, nbits) + index = index + 1 + id_offset = id_offset - 7 + } + } + output[nchars+1] = EOS +end + + +# REPACK_FP -- Reconstruct a floating point number. The input to REPACK_FP +# is an array of integers containing Cyber 30-bit floating point numbers +# in the least significant bits of each array element. The exponent, mantissa +# and two bits indicating the sign are extracted and used to reassemble +# the floating point value. Cyber blanks and overflows are returned as BLANK. + +procedure repack_fp (int_value, float_value, nvalues) + +int int_value[ARB], nvalues +real float_value[nvalues] + +int i, pixel +int exp, mantissa +real tbl[255] +int bitupk(), and(), not() +include "pow.inc" + +begin + do i=1, nvalues { + pixel = int_value[i] + # Check for blanks + if (pixel == 1777000000B) { + float_value[i] = BLANK + next + } + + # Check "bit59" and complement all bits if it is set + if (and (pixel, 4000000000B) != 0) { + pixel = not (pixel) + mantissa = -and (pixel, 777777B) + } else + mantissa = and (pixel, 777777B) + + # Extract and interpret exponent: remove Cyber bias of 2000B + # and convert to two's complement if negative number + exp = bitupk (pixel, 19, 11) + if (exp > 1777B) + # "bit58" is set, positive exponent + exp = exp - 2000B + else + # negative exponent + exp = exp - 1777B + + # Reconstruct the floating point value: 30 is added to the + # exponent because only the top 18 bits of the 48-bit mantissa + # were extracted; the 129 is to register the data array index. + # float_value[i] = real(mantissa) * 2 ** (exp + 30) + # (tbl[1] = 2 ** -128) ==> (2 ** n = tbl[n + 129]). + + exp = exp + 30 + 129 + if (exp <= 0) { + #call eprintf ( + #"RDUMPF_REPACK_FP: pixel with exponent underflow seen\n") + float_value[i] = 0.0 + } else if (exp > 255) { + #call eprintf ( + #"RDUMPF_REPACK_FP: pixel with exponent overflow seen\n") + float_value[i] = MAX_REAL + } else if (exp > 0 && exp <= 255) + float_value[i] = real(mantissa) * tbl[exp] + } +end + + +# DISPLAY_CODE -- returns the ascii character equivalent to the display +# code character input. The Cyber uses the 63-character display code +# set internally, although the 64-character set is used for output. + +procedure display_code (in_char, out_char) + +char in_char, out_char +char dc[64] +int i + +data (dc[i], i=1,8) /072B, 101B, 102B, 103B, 104B, 105B, 106B, 107B/ +data (dc[i], i=9,16) /110B, 111B, 112B, 113B, 114B, 115B, 116B, 117B/ +data (dc[i], i=17,24) /120B, 121B, 122B, 123B, 124B, 125B, 126B, 127B/ +data (dc[i], i=25,32) /130B, 131B, 132B, 060B, 061B, 062B, 063B, 064B/ +data (dc[i], i=33,40) /065B, 066B, 067B, 070B, 071B, 053B, 055B, 052B/ +data (dc[i], i=41,48) /057B, 050B, 051B, 044B, 075B, 040B, 054B, 056B/ +data (dc[i], i=49,56) /043B, 133B, 135B, 045B, 042B, 137B, 041B, 046B/ +data (dc[i], i=57,64) /047B, 077B, 074B, 076B, 100B, 134B, 136B, 073B/ + +begin + out_char = dc[in_char + 1] +end diff --git a/noao/mtlocal/cyber/cyrheader.x b/noao/mtlocal/cyber/cyrheader.x new file mode 100644 index 00000000..7471b118 --- /dev/null +++ b/noao/mtlocal/cyber/cyrheader.x @@ -0,0 +1,120 @@ +include +include +include +include "cyber.h" + +# CY_READ_HEADER -- reads the IPPS header (64 60-bit words) into the 128 element +# integer array "header". Any extraneous information between the header +# and data is skipped; the tape is left positioned at the first data +# record. + +int procedure cy_read_header (rd, dt) + +int rd +pointer dt +int header[NINT_CYBER_WRD * LEN_HEADER] +int npru_skip +int read_dumpf() +errchk cy_unpk_header, read_dumpf, cy_skip_pru, unpack_cyber_record +errchk order_cyber_bits + +begin + # Read the header into array header, one cyber word per two elements + if (read_dumpf (rd, header, LEN_HEADER) == EOF) + return (EOF) + + # Unpack bit stream and fill structure dt + iferr { + call cy_unpk_header (header, dt) + npru_skip = PRU_ROW_ONE(dt) - 1 + } then { + call erract (EA_WARN) + # Position to first row of raster before posting error + if (npru_skip > 0) + call cy_skip_pru (rd, npru_skip) + call error (1, "Bad header, attempting to skip raster") + } + + # Position to first row of IPPS raster + if (npru_skip > 0) + call cy_skip_pru (rd, npru_skip) + + return (OK) +end + + +# CY_LIST_HEADER -- prints the IPPS header information. + +procedure cy_list_header (dt, file_number, raster_num) + +pointer dt +int raster_num, file_number + +begin + # Print header information from IPPS raster + call printf ("[%d.%d]%7t IPPS_ID: %s\n") + call pargi (file_number) + call pargi (raster_num) + call pargstr (IPPS_ID(dt)) + call printf ("%7t NCOLS=%d, NROWS=%d, MIN=%g, MAX=%g, NBPP=%d\n") + call pargi (NCOLS(dt)) + call pargi (NROWS(dt)) + call pargr (DATA_MIN(dt)) + call pargr (DATA_MAX(dt)) + call pargi (BITS_PIXEL(dt)) +end + + +# CY_UNPK_HEADER -- unpacks header words from the char array header +# and fills the program data structure. A few values are checked to +# make sure a valid IPPS raster is being read. Offsets to various +# header words have been defined previously. + +procedure cy_unpk_header (header, dt) + +int header[NINT_CYBER_WRD * LEN_HEADER] +pointer dt + +begin + # From array header, first the ID is unpacked + call unpk_id (header, 0, IPPS_ID(dt)) + + # An EOR marker terminates each raster + PRU_EOR(dt) = header[EOR_OFFSET] + + # The PRU containing the first data row + PRU_ROW_ONE(dt) = header[FIRST_PRU_OFFSET] + + # Most significant 30 bits of the data min are used + call unpk_60r (header, MIN_OFFSET, DATA_MIN(dt), 1) + + # Most significant 30 bits of the data max are used + call unpk_60r (header, MAX_OFFSET, DATA_MAX(dt), 1) + + # Bits per pixel is unpacked and tested + BITS_PIXEL(dt) = header[DATA_TYPE_OFFSET] + + switch (BITS_PIXEL(dt)) { + case 12,20,30,60: + ; + default: + call error (2, "Incorrect IPPS BITS_PIXEL") + } + + # Number of columns is unpacked and tested + NCOLS(dt) = header[NCOLS_OFFSET] + if (NCOLS(dt) <= 0) + call error (2, "IPPS ncols <= 0") + + # Number of Cyber words per row must be integral # of PRU's + WRDS_PER_ROW(dt) = header[NWORDS_OFFSET] + NPRU_ROW(dt) = WRDS_PER_ROW(dt) / LEN_PRU + + if (mod (WRDS_PER_ROW(dt), LEN_PRU) != 0) + call error (2, "Invalid IPPS NWPR") + + # Number of rows is unpacked and tested + NROWS(dt) = header[NROWS_OFFSET] + if (NROWS(dt) <= 0) + call error (2, "IPPS nrows <= 0") +end diff --git a/noao/mtlocal/cyber/cyrimage.x b/noao/mtlocal/cyber/cyrimage.x new file mode 100644 index 00000000..5766e838 --- /dev/null +++ b/noao/mtlocal/cyber/cyrimage.x @@ -0,0 +1,166 @@ +include +include +include +include "cyber.h" + +# READ_IPPS_ROWS -- reads the ipps raster image row by row from the tape and +# writes the output image. At the completion of READ_IMAGE, the tape is +# positioned to the header record of the next ipps raster. + +procedure read_ipps_rows (rd, out_fname, data_type, dt) + +int rd, data_type +char out_fname[SZ_FNAME] +pointer dt + +pointer sp, im, spp_buf +int npru_skip, i, stat +long clktime() + +int read_dumpf() +pointer immap(), impl2r() +errchk cy_skip_pru, immap, ipps_to_iraf, read_dumpf() + +begin + # Allocate buffer for ipps raster pixels + call smark (sp) + call salloc (spp_buf, WRDS_PER_ROW(dt) * NINT_CYBER_WRD, TY_INT) + + # Map new iraf image and set up image header + im = immap (out_fname, NEW_IMAGE, 0) + IM_LEN(im, 1) = NCOLS(dt) + IM_LEN(im, 2) = NROWS(dt) + call strcpy (IPPS_ID(dt), IM_TITLE(im), SZ_IMTITLE) + IM_MIN(im) = DATA_MIN(dt) + IM_MAX(im) = DATA_MAX(dt) + + # Set optimum image pixel type + if (data_type == NOT_SET) { + switch (BITS_PIXEL(dt)) { + case 12: + IM_PIXTYPE(im) = TY_SHORT + case 20: + IM_PIXTYPE(im) = TY_REAL + case 30: + IM_PIXTYPE(im) = TY_REAL + case 60: + IM_PIXTYPE(im) = TY_REAL + default: + call error (3, "IPPS BITS_PIXEL is incorrect") + } + } else + IM_PIXTYPE(im) = data_type + IM_LIMTIME(im) = clktime (long(0)) + + # Loop over rows to read, reorder and convert pixels. + for (i=1; i <= NROWS(dt); i=i+1) { + iferr { + stat = read_dumpf (rd, Memi[spp_buf], WRDS_PER_ROW(dt)) + } then { + call imunmap (im) + call sfree (sp) + call erract (EA_WARN) + return + } + if (stat == EOF) { + call imunmap (im) + call sfree (sp) + call eprintf ("Premature EOF in image at row # %d\n") + call pargi (i) + return + } + + call ipps_to_iraf (Memi[spp_buf], Memr[impl2r(im,i)], NCOLS(dt), + BITS_PIXEL(dt)) + } + + # Skip from present position to end of rcopy raster + npru_skip = PRU_EOR(dt) - PRU_ROW_ONE(dt) - (NPRU_ROW(dt)*NROWS(dt)) + call cy_skip_pru (rd, npru_skip) + + call imunmap (im) + call sfree (sp) +end + + +# IPPS_TO_IRAF -- performs the conversion from Cyber pixels to IRAF pixels. +# Each row of the ipps image is required to occupy an integral number of Cyber +# PRU's, so the input buffer contains pixels plus filler. The number of +# pixels unpacked from this buffer will always fill an integral number of +# Cyber words; a maximum of 4 extraneous pixels will be unpacked. + +procedure ipps_to_iraf (in_buf, iraf_real, npix, nbits_pixel) + +int in_buf[ARB] +real iraf_real[npix] +pointer iraf_int, sp +int nbits_pixel, npix, offset, npix_unpk, npix_cyber_wrd +errchk unpk_12, unpk_20 + +begin + # Allocate (maximum) space needed on the stack. + call smark (sp) + call salloc (iraf_int, npix + 4, TY_INT) + offset = 1 + + switch (nbits_pixel) { + case 12: + npix_cyber_wrd = 5 + npix_unpk = ((npix + 4) / npix_cyber_wrd) * npix_cyber_wrd + call unpk_12 (in_buf, offset, Memi[iraf_int], npix_unpk) + call achtir (Memi[iraf_int], iraf_real, npix) + + case 20: + npix_cyber_wrd = 3 + npix_unpk = ((npix + 2) / npix_cyber_wrd) * npix_cyber_wrd + call unpk_20 (in_buf, offset, Memi[iraf_int], npix_unpk) + call achtir (Memi[iraf_int], iraf_real, npix) + + case 30: + call unpk_30 (in_buf, offset, iraf_real, npix) + + case 60: + call unpk_60r (in_buf, offset, iraf_real, npix) + + default: + call error (5, "Illegal IPPS #B/P") + } + + call sfree (sp) +end + + +# CY_SKIP_IMAGE -- skips over an IPPS raster once the header has been +# read. When SKIP_IMAGE returns, the tape is positioned to the first +# record of the next tape image. + +procedure cy_skip_image (rd, dt) + +int rd +pointer dt +int npru_skip +errchk cy_skip_pru + +begin + # Calculate number of PRU's in image + npru_skip = PRU_EOR(dt) - PRU_ROW_ONE(dt) + call cy_skip_pru (rd, npru_skip) +end + +# CY_SKIP_PRU -- positions the tape by skipping the requested number of PRU's. + +procedure cy_skip_pru (rd, npru_skip) + +int rd, npru_skip +pointer sp, dummy +int read_dumpf() + +begin + call smark (sp) + call salloc (dummy, NINT_CYBER_WRD * LEN_PRU * npru_skip, TY_INT) + + if (read_dumpf (rd, Memi[dummy], npru_skip * LEN_PRU) == EOF) + call eprintf ("Unexpected EOF when skipping image\n") + + call sfree (sp) +end diff --git a/noao/mtlocal/cyber/mkpkg b/noao/mtlocal/cyber/mkpkg new file mode 100644 index 00000000..b546080e --- /dev/null +++ b/noao/mtlocal/cyber/mkpkg @@ -0,0 +1,22 @@ +# The four cyber format readers LDUMPF, RDUMPF, RIDSFILE, RIDSOUT contribute +# the following sources to the dataio package library: + +$checkout libpkg.a ../ +$update libpkg.a +$checkin libpkg.a ../ +$exit + +libpkg.a: + @rrcopy + + cykeywords.x cyber.h + cyrbits.x pow.inc powd.inc cyber.h + cyrheader.x cyber.h + cyrimage.x cyber.h + rdumpf.x cyber.h + rpft.x cyber.h + t_ldumpf.x cyber.h + t_rdumpf.x cyber.h + t_ridsfile.x cyber.h + t_ridsout.x cyber.h + ; diff --git a/noao/mtlocal/cyber/pow.inc b/noao/mtlocal/cyber/pow.inc new file mode 100644 index 00000000..b87b75bb --- /dev/null +++ b/noao/mtlocal/cyber/pow.inc @@ -0,0 +1,86 @@ + data (tbl[i], i=1, 3) /0.0 , 0.0 , 0.0 / + data (tbl[i], i=4, 6) /0.0 , 0.0 , 0.0 / + data (tbl[i], i=7, 9) /1.880791E-37, 3.761582E-37, 7.523164E-37/ + data (tbl[i], i=10, 12) /1.504633E-36, 3.009266E-36, 6.018531E-36/ + data (tbl[i], i=13, 15) /1.203706E-35, 2.407412E-35, 4.814825E-35/ + data (tbl[i], i=16, 18) /9.629650E-35, 1.925930E-34, 3.851860E-34/ + data (tbl[i], i=19, 21) /7.703720E-34, 1.540744E-33, 3.081488E-33/ + data (tbl[i], i=22, 24) /6.162976E-33, 1.232595E-32, 2.465190E-32/ + data (tbl[i], i=25, 27) /4.930381E-32, 9.860761E-32, 1.972152E-31/ + data (tbl[i], i=28, 30) /3.944305E-31, 7.888609E-31, 1.577722E-30/ + data (tbl[i], i=31, 33) /3.155444E-30, 6.310887E-30, 1.262177E-29/ + data (tbl[i], i=34, 36) /2.524355E-29, 5.048710E-29, 1.009742E-28/ + data (tbl[i], i=37, 39) /2.019484E-28, 4.038968E-28, 8.077936E-28/ + data (tbl[i], i=40, 42) /1.615587E-27, 3.231174E-27, 6.462349E-27/ + data (tbl[i], i=43, 45) /1.292470E-26, 2.584939E-26, 5.169879E-26/ + data (tbl[i], i=46, 48) /1.033976E-25, 2.067952E-25, 4.135903E-25/ + data (tbl[i], i=49, 51) /8.271806E-25, 1.654361E-24, 3.308722E-24/ + data (tbl[i], i=52, 54) /6.617445E-24, 1.323489E-23, 2.646978E-23/ + data (tbl[i], i=55, 57) /5.293956E-23, 1.058791E-22, 2.117582E-22/ + data (tbl[i], i=58, 60) /4.235165E-22, 8.470329E-22, 1.694066E-21/ + data (tbl[i], i=61, 63) /3.388132E-21, 6.776264E-21, 1.355253E-20/ + data (tbl[i], i=64, 66) /2.710505E-20, 5.421011E-20, 1.084202E-19/ + data (tbl[i], i=67, 69) /2.168404E-19, 4.336809E-19, 8.673617E-19/ + data (tbl[i], i=70, 72) /1.734723E-18, 3.469447E-18, 6.938894E-18/ + data (tbl[i], i=73, 75) /1.387779E-17, 2.775558E-17, 5.551115E-17/ + data (tbl[i], i=76, 78) /1.110223E-16, 2.220446E-16, 4.440892E-16/ + data (tbl[i], i=79, 81) /8.881784E-16, 1.776357E-15, 3.552714E-15/ + data (tbl[i], i=82, 84) /7.105427E-15, 1.421085E-14, 2.842171E-14/ + data (tbl[i], i=85, 87) /5.684342E-14, 1.136868E-13, 2.273737E-13/ + data (tbl[i], i=88, 90) /4.547474E-13, 9.094947E-13, 1.818989E-12/ + data (tbl[i], i=91, 93) /3.637979E-12, 7.275958E-12, 1.455192E-11/ + data (tbl[i], i=94, 96) /2.910383E-11, 5.820766E-11, 1.164153E-10/ + data (tbl[i], i=97, 99) /2.328306E-10, 4.656613E-10, 9.313226E-10/ + data (tbl[i], i=100, 102) / 1.862645E-9, 3.725290E-9, 7.450581E-9/ + data (tbl[i], i=103, 105) / 1.490116E-8, 2.980232E-8, 5.960464E-8/ + data (tbl[i], i=106, 108) / 1.192093E-7, 2.384186E-7, 4.768372E-7/ + data (tbl[i], i=109, 111) / 9.536743E-7, 1.907349E-6, 3.814697E-6/ + data (tbl[i], i=112, 114) / 7.629395E-6, 1.525879E-5, 3.051758E-5/ + data (tbl[i], i=115, 117) / 6.103516E-5, 1.220703E-4, 2.441406E-4/ + data (tbl[i], i=118, 120) / 4.882813E-4, 9.765625E-4, 0.001953125/ + data (tbl[i], i=121, 123) / 0.00390625, 0.0078125, 0.015625/ + data (tbl[i], i=124, 126) / 0.03125, 0.0625, 0.125/ + data (tbl[i], i=127, 129) / 0.25, 0.5, 1./ + data (tbl[i], i=130, 132) / 2., 4., 8./ + data (tbl[i], i=133, 135) / 16., 32., 64./ + data (tbl[i], i=136, 138) / 128., 256., 512./ + data (tbl[i], i=139, 141) / 1024., 2048., 4096./ + data (tbl[i], i=142, 144) / 8192., 16384., 32768./ + data (tbl[i], i=145, 147) / 65536., 131072., 262144./ + data (tbl[i], i=148, 150) / 524288., 1048576., 2097152./ + data (tbl[i], i=151, 153) / 4194304., 8388608., 16777216./ + data (tbl[i], i=154, 156) / 33554432., 67108864., 134217728./ + data (tbl[i], i=157, 159) / 268435456., 536870912., 1.073742E9/ + data (tbl[i], i=160, 162) / 2.147484E9, 4.294967E9, 8.589935E9/ + data (tbl[i], i=163, 165) / 1.717987E10, 3.435974E10, 6.871948E10/ + data (tbl[i], i=166, 168) / 1.374390E11, 2.748779E11, 5.497558E11/ + data (tbl[i], i=169, 171) / 1.099512E12, 2.199023E12, 4.398047E12/ + data (tbl[i], i=172, 174) / 8.796093E12, 1.759219E13, 3.518437E13/ + data (tbl[i], i=175, 177) / 7.036874E13, 1.407375E14, 2.814750E14/ + data (tbl[i], i=178, 180) / 5.629500E14, 1.125900E15, 2.251800E15/ + data (tbl[i], i=181, 183) / 4.503600E15, 9.007199E15, 1.801440E16/ + data (tbl[i], i=184, 186) / 3.602880E16, 7.205759E16, 1.441152E17/ + data (tbl[i], i=187, 189) / 2.882304E17, 5.764608E17, 1.152922E18/ + data (tbl[i], i=190, 192) / 2.305843E18, 4.611686E18, 9.223372E18/ + data (tbl[i], i=193, 195) / 1.844674E19, 3.689349E19, 7.378698E19/ + data (tbl[i], i=196, 198) / 1.475740E20, 2.951479E20, 5.902958E20/ + data (tbl[i], i=199, 201) / 1.180592E21, 2.361183E21, 4.722366E21/ + data (tbl[i], i=202, 204) / 9.444733E21, 1.888947E22, 3.777893E22/ + data (tbl[i], i=205, 207) / 7.555786E22, 1.511157E23, 3.022315E23/ + data (tbl[i], i=208, 210) / 6.044629E23, 1.208926E24, 2.417852E24/ + data (tbl[i], i=211, 213) / 4.835703E24, 9.671407E24, 1.934281E25/ + data (tbl[i], i=214, 216) / 3.868563E25, 7.737125E25, 1.547425E26/ + data (tbl[i], i=217, 219) / 3.094850E26, 6.189700E26, 1.237940E27/ + data (tbl[i], i=220, 222) / 2.475880E27, 4.951760E27, 9.903520E27/ + data (tbl[i], i=223, 225) / 1.980704E28, 3.961408E28, 7.922816E28/ + data (tbl[i], i=226, 228) / 1.584563E29, 3.169127E29, 6.338253E29/ + data (tbl[i], i=229, 231) / 1.267651E30, 2.535301E30, 5.070602E30/ + data (tbl[i], i=232, 234) / 1.014120E31, 2.028241E31, 4.056482E31/ + data (tbl[i], i=235, 237) / 8.112964E31, 1.622593E32, 3.245186E32/ + data (tbl[i], i=238, 240) / 6.490371E32, 1.298074E33, 2.596148E33/ + data (tbl[i], i=241, 243) / 5.192297E33, 1.038459E34, 2.076919E34/ + data (tbl[i], i=244, 246) / 4.153837E34, 8.307675E34, 1.661535E35/ + data (tbl[i], i=247, 249) / 3.323070E35, 6.646140E35, 1.329228E36/ + data (tbl[i], i=250, 252) / 2.658456E36, 5.316912E36, 1.063382E37/ + data (tbl[i], i=253, 255) / 2.126765E37, 4.253530E37, 8.507059E37/ + diff --git a/noao/mtlocal/cyber/powd.inc b/noao/mtlocal/cyber/powd.inc new file mode 100644 index 00000000..2dd8e81f --- /dev/null +++ b/noao/mtlocal/cyber/powd.inc @@ -0,0 +1,128 @@ + data (tbl(i), i=1, 2) /2.93873587706D-39, 5.87747175411D-39/ + data (tbl(i), i=3, 4) /1.17549435082D-38, 2.35098870164D-38/ + data (tbl(i), i=5, 6) /4.70197740329D-38, 9.40395480658D-38/ + data (tbl(i), i=7, 8) /1.88079096132D-37, 3.76158192263D-37/ + data (tbl(i), i=9, 10) /7.52316384526D-37, 1.50463276905D-36/ + data (tbl(i), i=11, 12) /3.00926553811D-36, 6.01853107621D-36/ + data (tbl(i), i=13, 14) /1.20370621524D-35, 2.40741243048D-35/ + data (tbl(i), i=15, 16) /4.81482486097D-35, 9.62964972194D-35/ + data (tbl(i), i=17, 18) /1.92592994439D-34, 3.85185988877D-34/ + data (tbl(i), i=19, 20) /7.70371977755D-34, 1.54074395551D-33/ + data (tbl(i), i=21, 22) /3.08148791102D-33, 6.16297582204D-33/ + data (tbl(i), i=23, 24) /1.23259516441D-32, 2.46519032882D-32/ + data (tbl(i), i=25, 26) /4.93038065763D-32, 9.86076131526D-32/ + data (tbl(i), i=27, 28) /1.97215226305D-31, 3.94430452611D-31/ + data (tbl(i), i=29, 30) /7.88860905221D-31, 1.57772181044D-30/ + data (tbl(i), i=31, 32) /3.15544362088D-30, 6.31088724177D-30/ + data (tbl(i), i=33, 34) /1.26217744835D-29, 2.52435489671D-29/ + data (tbl(i), i=35, 36) /5.04870979341D-29, 1.00974195868D-28/ + data (tbl(i), i=37, 38) /2.01948391737D-28, 4.03896783473D-28/ + data (tbl(i), i=39, 40) /8.07793566946D-28, 1.61558713389D-27/ + data (tbl(i), i=41, 42) /3.23117426779D-27, 6.46234853557D-27/ + data (tbl(i), i=43, 44) /1.29246970711D-26, 2.58493941423D-26/ + data (tbl(i), i=45, 46) /5.16987882846D-26, 1.03397576569D-25/ + data (tbl(i), i=47, 48) /2.06795153138D-25, 4.13590306277D-25/ + data (tbl(i), i=49, 50) /8.27180612553D-25, 1.65436122511D-24/ + data (tbl(i), i=51, 52) /3.30872245021D-24, 6.61744490042D-24/ + data (tbl(i), i=53, 54) /1.32348898008D-23, 2.64697796017D-23/ + data (tbl(i), i=55, 56) /5.29395592034D-23, 1.05879118407D-22/ + data (tbl(i), i=57, 58) /2.11758236814D-22, 4.23516473627D-22/ + data (tbl(i), i=59, 60) /8.47032947254D-22, 1.69406589451D-21/ + data (tbl(i), i=61, 62) /3.38813178902D-21, 6.77626357803D-21/ + data (tbl(i), i=63, 64) /1.35525271561D-20, 2.71050543121D-20/ + data (tbl(i), i=65, 66) /5.42101086243D-20, 1.08420217249D-19/ + data (tbl(i), i=67, 68) /2.16840434497D-19, 4.33680868994D-19/ + data (tbl(i), i=69, 70) /8.67361737988D-19, 1.73472347598D-18/ + data (tbl(i), i=71, 72) /3.46944695195D-18, 6.93889390391D-18/ + data (tbl(i), i=73, 74) /1.38777878078D-17, 2.77555756156D-17/ + data (tbl(i), i=75, 76) /5.55111512313D-17, 1.11022302463D-16/ + data (tbl(i), i=77, 78) /2.22044604925D-16, 4.44089209850D-16/ + data (tbl(i), i=79, 80) /8.88178419700D-16, 1.77635683940D-15/ + data (tbl(i), i=81, 82) /3.55271367880D-15, 7.10542735760D-15/ + data (tbl(i), i=83, 84) /1.42108547152D-14, 2.84217094304D-14/ + data (tbl(i), i=85, 86) /5.68434188608D-14, 1.13686837722D-13/ + data (tbl(i), i=87, 88) /2.27373675443D-13, 4.54747350886D-13/ + data (tbl(i), i=89, 90) /9.09494701773D-13, 1.81898940355D-12/ + data (tbl(i), i=91, 92) /3.63797880709D-12, 7.27595761418D-12/ + data (tbl(i), i=93, 94) /1.45519152284D-11, 2.91038304567D-11/ + data (tbl(i), i=95, 96) /5.82076609135D-11, 1.16415321827D-10/ + data (tbl(i), i=97, 98) /2.32830643654D-10, 4.65661287308D-10/ + data (tbl(i), i=99, 100) /9.31322574615D-10, 1.86264514923D-9/ + data (tbl(i), i=101, 102) /3.72529029846D-9, 7.45058059692D-9/ + data (tbl(i), i=103, 104) /1.49011611938D-8, 2.98023223877D-8/ + data (tbl(i), i=105, 106) /5.96046447754D-8, 1.19209289551D-7/ + data (tbl(i), i=107, 108) /2.38418579102D-7, 4.76837158203D-7/ + data (tbl(i), i=109, 110) /9.53674316406D-7, 1.90734863281D-6/ + data (tbl(i), i=111, 112) /3.81469726563D-6, 7.62939453125D-6/ + data (tbl(i), i=113, 114) /1.52587890625D-5, 3.05175781250D-5/ + data (tbl(i), i=115, 116) /6.10351562500D-5, 1.22070312500D-4/ + data (tbl(i), i=117, 118) /2.44140625000D-4, 4.88281250000D-4/ + data (tbl(i), i=119, 120) /9.76562500000D-4, 0.001953125/ + data (tbl(i), i=121, 122) /0.00390625, 0.0078125/ + data (tbl(i), i=123, 124) /0.015625, 0.03125/ + data (tbl(i), i=125, 126) /0.0625, 0.125/ + data (tbl(i), i=127, 128) /0.25, 0.5/ + data (tbl(i), i=129, 130) /1., 2./ + data (tbl(i), i=131, 132) /4., 8./ + data (tbl(i), i=133, 134) /16., 32./ + data (tbl(i), i=135, 136) /64., 128./ + data (tbl(i), i=137, 138) /256., 512./ + data (tbl(i), i=139, 140) /1024., 2048./ + data (tbl(i), i=141, 142) /4096., 8192./ + data (tbl(i), i=143, 144) /16384., 32768./ + data (tbl(i), i=145, 146) /65536., 131072./ + data (tbl(i), i=147, 148) /262144., 524288./ + data (tbl(i), i=149, 150) /1048576., 2097152./ + data (tbl(i), i=151, 152) /4194304., 8388608./ + data (tbl(i), i=153, 154) /16777216., 33554432./ + data (tbl(i), i=155, 156) /67108864., 134217728./ + data (tbl(i), i=157, 158) /268435456., 536870912./ + data (tbl(i), i=159, 160) /1073741824., 2147483648./ + data (tbl(i), i=161, 162) /4294967296., 8589934592./ + data (tbl(i), i=163, 164) /17179869184., 34359738368./ + data (tbl(i), i=165, 166) /68719476736., 137438953472./ + data (tbl(i), i=167, 168) /274877906944., 549755813888./ + data (tbl(i), i=169, 170) /1099511627776., 2199023255552./ + data (tbl(i), i=171, 172) /4398046511104., 8796093022208./ + data (tbl(i), i=173, 174) /17592186044416., 35184372088832./ + data (tbl(i), i=175, 176) /70368744177664., 140737488355328./ + data (tbl(i), i=177, 178) /281474976710656., 562949953421312./ + data (tbl(i), i=179, 180) /1.12589990684D15, 2.25179981369D15/ + data (tbl(i), i=181, 182) /4.50359962737D15, 9.00719925474D15/ + data (tbl(i), i=183, 184) /1.80143985095D16, 3.60287970190D16/ + data (tbl(i), i=185, 186) /7.20575940379D16, 1.44115188076D17/ + data (tbl(i), i=187, 188) /2.88230376152D17, 5.76460752303D17/ + data (tbl(i), i=189, 190) /1.15292150461D18, 2.30584300921D18/ + data (tbl(i), i=191, 192) /4.61168601843D18, 9.22337203685D18/ + data (tbl(i), i=193, 194) /1.84467440737D19, 3.68934881474D19/ + data (tbl(i), i=195, 196) /7.37869762948D19, 1.47573952590D20/ + data (tbl(i), i=197, 198) /2.95147905179D20, 5.90295810359D20/ + data (tbl(i), i=199, 200) /1.18059162072D21, 2.36118324143D21/ + data (tbl(i), i=201, 202) /4.72236648287D21, 9.44473296574D21/ + data (tbl(i), i=203, 204) /1.88894659315D22, 3.77789318630D22/ + data (tbl(i), i=205, 206) /7.55578637259D22, 1.51115727452D23/ + data (tbl(i), i=207, 208) /3.02231454904D23, 6.04462909807D23/ + data (tbl(i), i=209, 210) /1.20892581961D24, 2.41785163923D24/ + data (tbl(i), i=211, 212) /4.83570327846D24, 9.67140655692D24/ + data (tbl(i), i=213, 214) /1.93428131138D25, 3.86856262277D25/ + data (tbl(i), i=215, 216) /7.73712524553D25, 1.54742504911D26/ + data (tbl(i), i=217, 218) /3.09485009821D26, 6.18970019643D26/ + data (tbl(i), i=219, 220) /1.23794003929D27, 2.47588007857D27/ + data (tbl(i), i=221, 222) /4.95176015714D27, 9.90352031428D27/ + data (tbl(i), i=223, 224) /1.98070406286D28, 3.96140812571D28/ + data (tbl(i), i=225, 226) /7.92281625143D28, 1.58456325029D29/ + data (tbl(i), i=227, 228) /3.16912650057D29, 6.33825300114D29/ + data (tbl(i), i=229, 230) /1.26765060023D30, 2.53530120046D30/ + data (tbl(i), i=231, 232) /5.07060240091D30, 1.01412048018D31/ + data (tbl(i), i=233, 234) /2.02824096037D31, 4.05648192073D31/ + data (tbl(i), i=235, 236) /8.11296384146D31, 1.62259276829D32/ + data (tbl(i), i=237, 238) /3.24518553658D32, 6.49037107317D32/ + data (tbl(i), i=239, 240) /1.29807421463D33, 2.59614842927D33/ + data (tbl(i), i=241, 242) /5.19229685853D33, 1.03845937171D34/ + data (tbl(i), i=243, 244) /2.07691874341D34, 4.15383748683D34/ + data (tbl(i), i=245, 246) /8.30767497366D34, 1.66153499473D35/ + data (tbl(i), i=247, 248) /3.32306998946D35, 6.64613997892D35/ + data (tbl(i), i=249, 250) /1.32922799578D36, 2.65845599157D36/ + data (tbl(i), i=251, 252) /5.31691198314D36, 1.06338239663D37/ + data (tbl(i), i=253, 254) /2.12676479326D37, 4.25352958651D37/ + data tbl(255) /8.50705917302D37/ diff --git a/noao/mtlocal/cyber/rdumpf.x b/noao/mtlocal/cyber/rdumpf.x new file mode 100644 index 00000000..d9371579 --- /dev/null +++ b/noao/mtlocal/cyber/rdumpf.x @@ -0,0 +1,283 @@ +include +include "cyber.h" + +# READ_DUMPF -- Read data from an array of cyber words, passing the +# requested number of cyber words to the output buffer. +# The word array contains control words and "Zero length PRU's" as well +# as data. The control words are read and used to properly +# interpret the PRU that follows. Each read by GET_CYBER_WORDS begins +# on a control word boundry. It is expected that read_dumpf returns an +# integral number of PRU's; if not, the data remaining in the PRU is discarded, +# as the buffer pointer must be positioned to the control word preceeding +# a PRU upon entry to read_dumpf. The number of cyber words read is +# returned by read_dumpf. Pointer ip is in units of cyber words; there +# are two array elements per cyber word. + +int procedure read_dumpf (rd, out, max_cyb_words) + +int rd, max_cyb_words +int out[ARB], pru_buf[NINT_CYBER_WRD * LEN_CYBER_READ] +int nwords, nwords_read, word_count, ip +int get_cyber_words(), read_dumpf_init(), bitupk() + +begin + if (mod (max_cyb_words, LEN_PRU) != 0) + call error (0, "READ_DUMPF: Non-integral PRU requested") + + for (nwords=0; nwords < max_cyb_words; nwords = nwords + word_count) { + # If necessary, read more data into the pru buffer. This + # will be necessary when all data in the pru buffer has + # been transferred to the output buffer. + + if (ip >= nwords_read) { + # Read another chunk of cyber words; reset buffer position + nwords_read = get_cyber_words (rd, pru_buf, LEN_CYBER_READ) + ip = 1 + + if (nwords_read == EOF) + break + } + + # Get the number of 12-bit bytes in the next pru from control word. + # This number stored in lowest 9 bits of the cyber word. + word_count = bitupk (pru_buf[(NINT_CYBER_WRD * ip) - 1], 1, 9) / + NBYTES_WORD + ip = ip + 1 + + if (word_count == 0) { + # Zero length PRU to be skipped over + ip = ip + LEN_PRU + next + } + + if (mod (word_count * 60, NBITS_CHAR) != 0) + call error (0, "READ_DUMPF: Impossible PRU to CHAR Conversion") + + call amovi (pru_buf[(NINT_CYBER_WRD * ip) - 1], out[(nwords * + NINT_CYBER_WRD) + 1], word_count * NINT_CYBER_WRD) + ip = ip + LEN_PRU + } + + if (nwords == 0) + return (EOF) + else + return (nwords) + +entry read_dumpf_init () + + ip = 1 + nwords_read = 0 +end + +.help read_cyber +.nf ________________________________________________________________________ +GET_CYBER_WORDS -- Read binary chars from a file, ignoring short "noise" +records and extraneous bits inserted to fill out a byte. The requested +number of cyber words is returned, one cyber word per two array elements. +Data is read from +the file buffer and placed in the output array by UNPACK_CYBER_WORDS. The +file buffer is refilled as necessary by calling READ; the output array is +supplied by the calling procedure. Cyber noise records (48 bits) and bits +used to fill out a byte are not transferred to the output array. +Variables marking the current position and top of the buffer are initialized +by GET_CYBER_WORDS_INIT; buf pos marks the buffer position in cyber words. +The number of cyber words read is returned as the function value; the number +of output array elements filled is twice the number of cyber words read. +.endhelp ___________________________________________________________________ + +int procedure get_cyber_words (rd, out, max_cyb_words) + +int rd, max_cyb_words +int out[ARB] +char cyber_buf[SZ_TAPE_BUFFER] +int buf_pos, ncyber_words +int nwords, word_chunk, nchars_read +int read(), get_cyber_words_init() + +begin + for (nwords = 0; nwords < max_cyb_words; nwords = nwords + word_chunk) { + # See if it is necessary to transfer more data from the binary file + # to the file buffer. This will be necessary when all data from the + # file buffer has been transferred to the output array. + + if (buf_pos >= ncyber_words) { + # Read the next non-noise record into cyber_buf; reset buf_pos + repeat { + nchars_read = read (rd, cyber_buf, SZ_TAPE_BUFFER) + } until (nchars_read >= NCHARS_NOISE || nchars_read == EOF) + buf_pos = 1 + ncyber_words = (nchars_read * NBITS_CHAR) / NBITS_CYBER_WORD + + if (nchars_read == EOF) + break + } + + # The number of cyber words to output is the smaller of the number + # requested or the number of words left in the buffer. + + word_chunk = min (max_cyb_words - nwords, ncyber_words- buf_pos + 1) + call unpack_cyber_words (cyber_buf, buf_pos, out, + (NINT_CYBER_WRD * nwords) + 1, word_chunk) + buf_pos = buf_pos + word_chunk + } + + if (nwords == 0) + return (EOF) + else + return (nwords) + +entry get_cyber_words_init () + + buf_pos = 1 + ncyber_words = 0 + return (OK) +end + + + +.help unpack_cyber_words +.nf __________________________________________________________________________ +UNPACK_CYBER_WORDS -- Convert a raw data array from a 60-bit Cyber computer into +an SPP array containing one Cyber word in each 64-bit increment. +The least significant Cyber bit is bit 1 in each output word and the most +significant bit is bit 60. [MACHDEP]. + +When the Cyber outputs an array of 60 bit Cyber words it moves a stream of +bits into output bytes, filling however many bytes as necessary to output a +given number of Cyber words. The most significant bits are output first, +i.e., bit 60 of the first word is moved to bit 8 of the first output byte, +bit 59 is moved to bit 7 of the first output byte, and so on. If effect the +Cyber byte flips each 60-bit word. + +To deal with Cyber words as an ordered bit stream we must reorder the bytes +and bits so that the least significant bits are first. This function is +performed by the primitives CYBOOW and CYBOEW (order odd/even Cyber word) +for individual 60-bit Cyber words. A portable (and less efficient) version +of order_cyber_bits is available which does not use these primitives. +.endhelp _____________________________________________________________________ + + +procedure unpack_cyber_words (raw_cyber, first_cyber_word, int_array, + first_output_element, ncyber_words) + +char raw_cyber[ARB] # raw Cyber array (e.g. from tape) +int first_cyber_word # first 60-bit Cyber word to be unpacked +int int_array[ARB] # output unpacked array of Cyber words +int first_output_element # first output integer to be filled +int ncyber_words # number of Cyber words to unpack + +bool odd_word +int word, inbit, outbit + +begin + odd_word = (mod (first_cyber_word, 2) == 1) + inbit = (first_cyber_word - 1) * NBITS_CYBER_WORD + 1 + outbit = (first_output_element - 1) * NBITS_INT + 1 + + do word = 1, ncyber_words { + # Call odd or even primitive to reorder bits and move 60-bit + # ordered Cyber word to the output array. + + if (odd_word) { + call cyboow (raw_cyber, inbit, int_array, outbit) + odd_word = false + } else { + call cyboew (raw_cyber, inbit, int_array, outbit) + odd_word = true + } + + inbit = inbit + NBITS_CYBER_WORD + outbit = outbit + NINT_CYBER_WRD * NBITS_INT + } +end + + +# The portable version of order_cyber_bits follows. +#.help order_cyber_bits +#.nf __________________________________________________________________________ +#ORDER_CYBER_BITS -- Convert a raw data array from a 60-bit Cyber computer into +#an SPP bit-array. The output SPP bit-array is a bit-packed array of 60-bit +#Cyber words, i.e., bits 1-60 are word 1, bits 61-120 are word 2, and so on. +#The least significant Cyber bit is bit 1 in each output word and the most +#significant bit is bit 60. [MACHDEP]. +# +#The byte stream from the Cyber contains bits 53-60 of the first word in the +#first byte, bits 45-52 in the second byte, and so on (most significant bytes +#first). In essence we swap the order of the 7 8-bit bytes and the 4-bit half +#byte in each 60 bit word. The bits in each byte are in the correct order. +# +#Each successive pair of Cyber words fits into 15 bytes. Byte 8 contains the +#last 4 bits of word 1 in the most signficant half of the byte and the first +#4 bits of word 2 in the first half of the byte. In each 60 bit word we must +#move bit segments (bytes or half bytes) as follows (for the VAX): +# +#Odd words (from N*60 bit-offset): +# [from] [to] [nbits] +# 1 53 8 +# 9 45 8 +# 17 37 8 +# 25 29 8 +# 33 21 8 +# 41 13 8 +# 49 5 8 +# 61 1 4 +# +#Even words (from N*60 bit-offset): +# [from] [to] [nbits] +# -3 57 4 +# 5 49 8 +# 13 41 8 +# 21 33 8 +# 29 25 8 +# 37 17 8 +# 45 9 8 +# 53 1 8 +#.endhelp _____________________________________________________________________ +# +#define NBITS_PER_WORD 60 +#define NSEGMENTS 8 +# +# +#procedure order_cyber_bits (raw_cyber, first_cyber_word, bit_array, +# ncyber_words) +# +#char raw_cyber[ARB] # raw Cyber array (e.g. from tape) +#int first_cyber_word # first 60-bit Cyber word to be unpacked +#char bit_array[ARB] # output bit-array +#int ncyber_words # number of Cyber words to unpack +# +#int word, inword, inbit, outbit, temp, i +#int o_from[NSEGMENTS], o_to[NSEGMENTS], o_nbits[NSEGMENTS] +#int e_from[NSEGMENTS], e_to[NSEGMENTS], e_nbits[NSEGMENTS] +#int bitupk() +# +#data o_from / 1, 9,17,25,33,41,49,61/ # odd words +#data o_to /53,45,37,29,21,13, 5, 1/ +#data o_nbits / 8, 8, 8, 8, 8, 8, 8, 4/ +#data e_from /-3, 5,13,21,29,37,45,53/ # even words +#data e_to /57,49,41,33,25,17, 9, 1/ +#data e_nbits / 4, 8, 8, 8, 8, 8, 8, 8/ +# +#begin +# do word = 1, ncyber_words { +# inword = first_cyber_word + word - 1 +# inbit = (inword - 1) * NBITS_PER_WORD +# outbit = ( word - 1) * NBITS_PER_WORD +# +# # Move bits to the output bit array. Segment list used depends +# # on whether the word is an odd or even word. This code will work +# # even if the caller only wishes to order a single word. +# +# if (mod (inword,2) == 1) { +# do i = 1, NSEGMENTS { +# temp = bitupk (raw_cyber, inbit + o_from[i], o_nbits[i]) +# call bitpak (temp, bit_array, outbit + o_to[i], o_nbits[i]) +# } +# } else { +# do i = 1, NSEGMENTS { +# temp = bitupk (raw_cyber, inbit + e_from[i], e_nbits[i]) +# call bitpak (temp, bit_array, outbit + e_to[i], e_nbits[i]) +# } +# } +# } +#end diff --git a/noao/mtlocal/cyber/rpft.x b/noao/mtlocal/cyber/rpft.x new file mode 100644 index 00000000..f21ec44c --- /dev/null +++ b/noao/mtlocal/cyber/rpft.x @@ -0,0 +1,211 @@ +include +include "cyber.h" + +.help read_cyber +.nf ________________________________________________________________________ +READ_CYBER -- Read binary chars from a file, ignoring short "noise" records. +Data is read in chunks from the file buffer into the output buffer. The +file buffer is refilled as necessary by calling READ; the output buffer is +supplied by the calling procedure. Cyber noise records (48 bits) +are not transferred to the output buffer and so are ignored. READ_CYBER_INIT +must be called to initialize variables for CYBER_READ. Variables marking the +current position and top of the buffer are initialized. +.endhelp ___________________________________________________________________ + +int procedure read_cyber (rd, out_buffer, maxch) + +int buf_pos, buf_top +int rd, maxch +char out_buffer[ARB] +char block_buf[SZ_TAPE_BUFFER] +int nchars, chunk_size, nchars_read +int read(), read_cyber_init() + +begin + for (nchars = 0; nchars < maxch; nchars = nchars + chunk_size) { + # See if it is necessary to transfer more data from the binary file + # to the file buffer. This will be necessary when all data from the + # file buffer has been moved to the output buffer. + + if (buf_pos >= buf_top) { + # Read the next non-noise record into block_buf; reset buf_pos + repeat { + nchars_read = read (rd, block_buf, SZ_TAPE_BUFFER) + } until (nchars_read > NCHARS_NOISE || nchars_read == EOF) + buf_pos = 1 + buf_top = nchars_read + } + + # The number of chars to output is the smaller of the number of + # characters requested or the number of characters left in the + # buffer + + if (nchars_read == EOF) + break + else + chunk_size = min (maxch - nchars, buf_top - buf_pos + 1) + + # Move data to output array, increment buffer offset + call amovc (block_buf[buf_pos], out_buffer[nchars+1], chunk_size) + buf_pos = buf_pos + chunk_size + } + + if (nchars == 0) + return (EOF) + else + return (nchars) + +entry read_cyber_init () + + buf_pos = 1 + buf_top = 0 + return (OK) +end + + +.help order_cyber_bits +.nf __________________________________________________________________________ +ORDER_CYBER_BITS -- Convert a raw data array from a 60-bit Cyber computer into +an SPP bit-array. The output SPP bit-array is a bit-packed array of 60-bit +Cyber words, i.e., bits 1-60 are word 1, bits 61-120 are word 2, and so on. +The least significant Cyber bit is bit 1 in each output word and the most +significant bit is bit 60. [MACHDEP]. + +When the Cyber outputs an array of 60 bit Cyber words it moves a stream of +bits into output bytes, filling however many bytes as necessary to output a +given number of Cyber words. The most significant bits are output first, +i.e., bit 60 of the first word is moved to bit 8 of the first output byte, +bit 59 is moved to bit 7 of the first output byte, and so on. If effect the +Cyber byte flips each 60-bit word. + +To deal with Cyber words as an ordered bit stream we must reorder the bytes +and bits so that the least significant bits are first. This function is +performed by the primitives CYBOOW and CYBOEW (order odd/even Cyber word) +for individual 60-bit Cyber words. A portable (and less efficient) version +of order_cyber_bits is available which does not use these primitives. +.endhelp _____________________________________________________________________ + + +procedure order_cyber_bits (raw_cyber, first_cyber_word, bit_array, + ncyber_words) + +char raw_cyber[ARB] # raw Cyber array (e.g. from tape) +int first_cyber_word # first 60-bit Cyber word to be unpacked +char bit_array[ARB] # output bit-array +int ncyber_words # number of Cyber words to unpack + +bool odd_word +int word, inbit, outbit + +begin + odd_word = (mod (first_cyber_word, 2) == 1) + inbit = (first_cyber_word - 1) * NBITS_CYBER_WORD + 1 + outbit = 1 + + do word = 1, ncyber_words { + # Call odd or even primitive to reorder bits and move 60-bit + # ordered Cyber word to the output array. + + if (odd_word) { + call cyboow (raw_cyber, inbit, bit_array, outbit) + odd_word = false + } else { + call cyboew (raw_cyber, inbit, bit_array, outbit) + odd_word = true + } + + inbit = inbit + NBITS_CYBER_WORD + outbit = outbit + NBITS_CYBER_WORD + } +end + + +# The portable version of order_cyber_bits follows. +#.help order_cyber_bits +#.nf __________________________________________________________________________ +#ORDER_CYBER_BITS -- Convert a raw data array from a 60-bit Cyber computer into +#an SPP bit-array. The output SPP bit-array is a bit-packed array of 60-bit +#Cyber words, i.e., bits 1-60 are word 1, bits 61-120 are word 2, and so on. +#The least significant Cyber bit is bit 1 in each output word and the most +#significant bit is bit 60. [MACHDEP]. +# +#The byte stream from the Cyber contains bits 53-60 of the first word in the +#first byte, bits 45-52 in the second byte, and so on (most significant bytes +#first). In essence we swap the order of the 7 8-bit bytes and the 4-bit half +#byte in each 60 bit word. The bits in each byte are in the correct order. +# +#Each successive pair of Cyber words fits into 15 bytes. Byte 8 contains the +#last 4 bits of word 1 in the most signficant half of the byte and the first +#4 bits of word 2 in the first half of the byte. In each 60 bit word we must +#move bit segments (bytes or half bytes) as follows (for the VAX): +# +#Odd words (from N*60 bit-offset): +# [from] [to] [nbits] +# 1 53 8 +# 9 45 8 +# 17 37 8 +# 25 29 8 +# 33 21 8 +# 41 13 8 +# 49 5 8 +# 61 1 4 +# +#Even words (from N*60 bit-offset): +# [from] [to] [nbits] +# -3 57 4 +# 5 49 8 +# 13 41 8 +# 21 33 8 +# 29 25 8 +# 37 17 8 +# 45 9 8 +# 53 1 8 +#.endhelp _____________________________________________________________________ +# +#define NBITS_PER_WORD 60 +#define NSEGMENTS 8 +# +# +#procedure order_cyber_bits (raw_cyber, first_cyber_word, bit_array, +# ncyber_words) +# +#char raw_cyber[ARB] # raw Cyber array (e.g. from tape) +#int first_cyber_word # first 60-bit Cyber word to be unpacked +#char bit_array[ARB] # output bit-array +#int ncyber_words # number of Cyber words to unpack +# +#int word, inword, inbit, outbit, temp, i +#int o_from[NSEGMENTS], o_to[NSEGMENTS], o_nbits[NSEGMENTS] +#int e_from[NSEGMENTS], e_to[NSEGMENTS], e_nbits[NSEGMENTS] +#int bitupk() +# +#data o_from / 1, 9,17,25,33,41,49,61/ # odd words +#data o_to /53,45,37,29,21,13, 5, 1/ +#data o_nbits / 8, 8, 8, 8, 8, 8, 8, 4/ +#data e_from /-3, 5,13,21,29,37,45,53/ # even words +#data e_to /57,49,41,33,25,17, 9, 1/ +#data e_nbits / 4, 8, 8, 8, 8, 8, 8, 8/ +# +#begin +# do word = 1, ncyber_words { +# inword = first_cyber_word + word - 1 +# inbit = (inword - 1) * NBITS_PER_WORD +# outbit = ( word - 1) * NBITS_PER_WORD +# +# # Move bits to the output bit array. Segment list used depends +# # on whether the word is an odd or even word. This code will work +# # even if the caller only wishes to order a single word. +# +# if (mod (inword,2) == 1) { +# do i = 1, NSEGMENTS { +# temp = bitupk (raw_cyber, inbit + o_from[i], o_nbits[i]) +# call bitpak (temp, bit_array, outbit + o_to[i], o_nbits[i]) +# } +# } else { +# do i = 1, NSEGMENTS { +# temp = bitupk (raw_cyber, inbit + e_from[i], e_nbits[i]) +# call bitpak (temp, bit_array, outbit + e_to[i], e_nbits[i]) +# } +# } +# } +#end diff --git a/noao/mtlocal/cyber/rrcopy/README b/noao/mtlocal/cyber/rrcopy/README new file mode 100644 index 00000000..ec7e7c70 --- /dev/null +++ b/noao/mtlocal/cyber/rrcopy/README @@ -0,0 +1,2 @@ +This directory contains source code for rrcopy, the Cyber RCOPY tape +reader. diff --git a/noao/mtlocal/cyber/rrcopy/Revisions b/noao/mtlocal/cyber/rrcopy/Revisions new file mode 100644 index 00000000..09bf9e12 --- /dev/null +++ b/noao/mtlocal/cyber/rrcopy/Revisions @@ -0,0 +1,15 @@ +.help revisions Jun88 noao.mtlocal.cyber.rrcopy +.nf +noao$mtlocal/cyber/rrcopy/t_rrcopy.x, rcrheader.x + Fixed two places in t_rrcopy where the procedure was returning + without closing the mt file. Procedure rc_read_header was + not returning the value (OK) when NOT at EOF. These two errors, + while present all along, had not been seen until rrcopy was modified + to read tapes with more then one datafile; see below. (4-AUG-88 ShJ) + +noao$mtlocal/cyber/rrcopy/t_rrcopy.x + Added a hidden parameter "datafile" to the rrcopy task. This + allows more than one file of rcopy format data per tape. With + this "extension" to the rcopy format, many rcopy files can be + archived on a single tape. (26-JULY-88 ShJ) +.endhelp diff --git a/noao/mtlocal/cyber/rrcopy/mkpkg b/noao/mtlocal/cyber/rrcopy/mkpkg new file mode 100644 index 00000000..6d2d99d3 --- /dev/null +++ b/noao/mtlocal/cyber/rrcopy/mkpkg @@ -0,0 +1,15 @@ +# The Cyber rcopy format reader RRCOPY makes the following contributions +# to the dataio package library: + +$checkout libpkg.a ../../ +$update libpkg.a +$checkin libpkg.a ../../ +$exit + +libpkg.a: + rcrbits.x ../pow.inc rrcopy.h + rcrheader.x rrcopy.h + rcrimage.x rrcopy.h + rrcopy.x rrcopy.h + t_rrcopy.x rrcopy.h + ; diff --git a/noao/mtlocal/cyber/rrcopy/rcrbits.x b/noao/mtlocal/cyber/rrcopy/rcrbits.x new file mode 100644 index 00000000..c6525787 --- /dev/null +++ b/noao/mtlocal/cyber/rrcopy/rcrbits.x @@ -0,0 +1,279 @@ +include +include +include +include "rrcopy.h" + +# RC_UP_12 -- Unpack 12-bit unsigned integers from a stream of bits. +# Each output integer word contains successive 12-bit increments +# of the input bit stream in the least significant bit positions. +# It is assummed that the initial_bit_offset is the first bit of a +# Cyber 60-bit word containing 5 packed 12-bit pixels, the first pixel +# in the highest 12 bits. + +procedure rc_up_12 (input, initial_bit_offset, output, npix_unpk) + +char input[ARB] +int output[npix_unpk], npix_unpk +int initial_bit_offset, nbits, n, nn, bit_offset +int npix_word, ncyb_words, index +int bitupk() + +begin + nbits = 12 + npix_word = 5 + if (mod (npix_unpk, npix_word) == 0) + ncyb_words = (npix_unpk) / npix_word + else + call error (0, "Incorrect number of pixels to be unpacked") + index = 1 + + do n = 1, ncyb_words { + bit_offset = initial_bit_offset + (n * 60) + do nn = 1, npix_word { + bit_offset = bit_offset - nbits + output[index] = bitupk (input, bit_offset, nbits) + if (output[index] == 7777B) + output[index] = BLANK + index = index + 1 + } + } +end + + +# RC_UP_20 -- Unpack 20-bit signed integers from a stream of bits. +# Each output integer word contains sucessive 20-bit increments of the input. +# Conversion from one's complement to two's complement is performed. +# It is assummed that initial_bit_offset is the first bit of a Cyber +# 60-bit word containing 3 packed 20-bit pixels, the first pixel in the +# highest 20 bits. + +procedure rc_up_20 (input, initial_bit_offset, output, npix_unpk) + +char input[ARB] +int output[npix_unpk], npix_unpk +int nbits, n, index, bit_offset, initial_bit_offset +int npix_word, ncyb_words, nn, pix_val +int bitupk() + +begin + nbits = 20 + npix_word = 3 + if (mod (npix_unpk, npix_word) == 0) + ncyb_words = npix_unpk / npix_word + else + call error (0, "Incorrect number of pixels to be unpacked") + index = 1 + + do n = 1, ncyb_words { + bit_offset = initial_bit_offset + (n * 60) + do nn = 1, npix_word { + bit_offset = bit_offset - nbits + pix_val = bitupk (input, bit_offset, nbits) + if (pix_val == 3777777B) + pix_val = BLANK + else if (and (pix_val, 2000000B) != 0) + # negative pixel + pix_val = -and (3777777B, not(pix_val)) + output[index] = pix_val + index = index + 1 + } + } +end + + +# RC_UP_30 -- unpack Cyber 30-bit floating point numbers from a stream of +# bits. The input bit stream is unpacked in 30-bit increments into +# an integer array. Procedure REPACK_FP is called to reconstruct the +# floating point numbers from this array. It is assumed initial_bit_offset +# is the first bit of a Cyber 60-bit word containing 2 30-bit pixels, the +# first pixel in the higher 30 bits. + +procedure rc_up_30 (input, initial_bit_offset, fp_value, npix) + +char input[ARB] +real fp_value[npix] +pointer int_buf, sp +int initial_bit_offset, npix, bit_offset +int nbits, n +int bitupk() + +begin + # Allocate buffer space, allowing for maximum of 1 extraneous pixel + call smark (sp) + call salloc (int_buf, npix + 1, TY_INT) + + nbits = 30 + bit_offset = initial_bit_offset - 60 + + do n = 1, npix, 2 { + bit_offset = bit_offset + 90 + Memi[int_buf + n - 1] = bitupk (input, bit_offset, 30) + bit_offset = bit_offset - nbits + Memi[int_buf + n] = bitupk (input, bit_offset, 30) + } + + call rc_repack_fp (Memi[int_buf], fp_value, npix) + call sfree (sp) +end + + +# RC_UP_60R -- Unpack Cyber 60-bit floating point numbers from a stream +# of bits. The 30 most significant bits from each 60-bit word are +# unpacked into an integer array. Procedure REPACK_FP is called to +# reconstruct the floating point numbers from this array. +# An 18-bit mantissa, 11-bit exponent and a sign bit are unpacked into +# the lower 30 bits of each output word. + +procedure rc_up_60r (input, initial_bit_offset, fp_value, nwords) + +char input[ARB] +real fp_value[nwords] +int initial_bit_offset, nwords, bit_offset +pointer int_buf, sp +int n, nbits_unpk, nbits +int bitupk() + +begin + # Allocate space on stack + call smark (sp) + call salloc (int_buf, nwords, TY_INT) + + nbits = 60 + nbits_unpk = 30 + bit_offset = initial_bit_offset + 30 + + do n = 1, nwords { + Memi[int_buf + n - 1] = bitupk (input, bit_offset, nbits_unpk) + bit_offset = bit_offset + 60 + } + + call rc_repack_fp (Memi[int_buf], fp_value, nwords) + call sfree (sp) +end + + +# RC_UP_60I -- Unpack 60-bit integers from a stream of bits. Each element +# of output contains only the lower 32 bits of each input word, as this +# procedure is called only for getting NROWS, NCOLS and a few other small +# positive integer values. (A 60-bit intger is not a valid IPPS pixel type.) + +procedure rc_up_60i (input, initial_bit_offset, output, nwords) + +char input[ARB] +int output[nwords] +int initial_bit_offset, nwords, bit_offset +int n, nbits_unpk, nbits +int bitupk() + +begin + nbits_unpk = NBITS_INT + nbits = 60 + bit_offset = initial_bit_offset + + do n = 1, nwords { + output[n] = bitupk (input, bit_offset, nbits_unpk) + bit_offset = bit_offset + 60 + } +end + + +# RC_UP_ID -- Unpacks ID string from input bit stream. The IPPS ID string is +# written in 7-bit ASCII, with eight characters per Cyber word. The lowest +# 4 bits of each 60-bit word is unused. The highest 7 bits of the first Cyber +# word in the bit stream contains the character count. + +procedure rc_up_id (input, output) + +char input[SZ_HEADER] +char output[SZ_HEADER] +int nbits, nchar_offset, id_offset, nchars, n +int nchars_word, ncyb_words, nn, index +int bitupk() + +begin + nbits = 7 + nchar_offset = NBITS_CYBER_WORD - 6 + nchars = bitupk (input, nchar_offset, nbits) + ncyb_words = (nchars + 7) / 8 + index = 1 + + do n = 1, ncyb_words { + if (n == 1) { + nchars_word = 7 + id_offset = nchar_offset - 7 + } else { + nchars_word = 8 + id_offset = (n * NBITS_CYBER_WORD) - 6 + } + do nn = 1, nchars_word { + output[index] = bitupk (input, id_offset, nbits) + index = index + 1 + id_offset = id_offset - 7 + } + } + output[nchars+1] = EOS +end + + +# RC_REPACK_FP -- returns a floating point number as the function value. +# The input to REPACK_FP is an integer containing a 30-bit Cyber floating +# point number in the least significant bits. The exponent, mantissa +# and two bits indicating the sign are extracted and used to reassemble +# the floating point value. Cyber blanks and overflows are returned as BLANK. + +procedure rc_repack_fp (int_value, float_value, nvalues) + +int int_value[ARB], nvalues +real float_value[nvalues] + +int i, pixel +int exp, mantissa +real tbl[255] +int bitupk(), and(), not() +include "../pow.inc" + +begin + do i=1, nvalues { + pixel = int_value[i] + # Check for blanks + if (pixel == 1777000000B) { + float_value[i] = BLANK + next + } + + # Check "bit59" and complement all bits if it is set + if (and (pixel, 4000000000B) != 0) { + pixel = not (pixel) + mantissa = -and (pixel, 777777B) + } else + mantissa = and (pixel, 777777B) + + # Extract and interpret exponent: remove Cyber bias of 2000B + # and convert to two's complement if negative number + exp = bitupk (pixel, 19, 11) + if (exp > 1777B) + # "bit58" is set, positive exponent + exp = exp - 2000B + else + # negative exponent + exp = exp - 1777B + + # Reconstruct the floating point value: 30 is added to the + # exponent because only the top 18 bits of the 48-bit mantissa + # were extracted; the 129 is to register the data array index. + # float_value[i] = real(mantissa) * 2 ** (exp + 30) + # (tbl[1] = 2 ** -128) ==> (2 ** n = tbl[n + 129]). + + exp = exp + 30 + 129 + if (exp <= 0) { + #call eprintf ( + #"RRCOPY_RPACK_FP: Exponent underflow in following record\n") + float_value[i] = 0.0 + } else if (exp > 255) { + #call eprintf ( + #"RRCOPY_REPACK_FP: Exponent overflow in following record\n") + float_value[i] = MAX_REAL + } else if (exp > 0 && exp <= 255) + float_value[i] = double (mantissa) * tbl[exp] + } +end diff --git a/noao/mtlocal/cyber/rrcopy/rcrheader.x b/noao/mtlocal/cyber/rrcopy/rcrheader.x new file mode 100644 index 00000000..2187c6a1 --- /dev/null +++ b/noao/mtlocal/cyber/rrcopy/rcrheader.x @@ -0,0 +1,119 @@ +include +include +include +include "rrcopy.h" + +# RC_READ_HEADER -- reads the IPPS header (64 60-bit words) as +# a bit stream into consecutive elements of char array header. +# Any extraneous information between the header and data is skipped; +# the tape is left positioned at the first data record. + +int procedure rc_header_read (rd, rp) + +int rd +pointer rp +char raw_header[SZ_HEADER], header[SZ_HEADER] +int nchars_to_skip, first_word +int rc_read_cyber() +errchk rc_header_unpk, rc_skip_chars, rc_read_cyber, rc_order_cyber_bits + +begin + if (rc_read_cyber (rd, raw_header, SZ_HEADER) == EOF) + return (EOF) + + first_word = 1 + call rc_order_cyber_bits (raw_header, first_word, header, LEN_PRU) + + # Unpack bit stream and fill structure rp + iferr { + call rc_header_unpk (header, rp) + nchars_to_skip = (PRU_ROW_ONE(rp) - 1) * NBITS_PRU / NBITS_CHAR + } then { + call erract (EA_WARN) + # Position to first row of raster before posting error + if (nchars_to_skip > 0) + call rc_skip_chars (rd, nchars_to_skip) + call error (1, "Bad header, attempting to skip raster") + } + + # Position to first row of IPPS raster + if (nchars_to_skip > 0) + call rc_skip_chars (rd, nchars_to_skip) + + return (OK) +end + + +# RC_LIST_HEADER -- prints the RCOPY header information. + +procedure rc_list_header (rp, raster_num) + +pointer rp +int raster_num + +begin + # Print header information from rcopy tape + call printf ("[%d]%7t IPPS_ID: %s\n") + call pargi (raster_num) + call pargstr (IPPS_ID(rp)) + call printf ("%7t NCOLS=%d, NROWS=%d, MIN=%g, MAX=%g, NBPP=%d\n") + call pargi (NCOLS(rp)) + call pargi (NROWS(rp)) + call pargr (DATA_MIN(rp)) + call pargr (DATA_MAX(rp)) + call pargi (BITS_PIXEL(rp)) +end + + +# RC_UNPACK_HEADER -- unpacks header words from the char array header +# and fills the RCOPY data structure. A few values are checked to +# make sure a valid IPPS raster is being read. Offsets to various +# header words have been defined previously. + +procedure rc_header_unpk (header, rp) + +char header[SZ_HEADER] +pointer rp + +begin + # From the reordered array, first the ID is unpacked + call rc_up_id (header, IPPS_ID(rp)) + + # An EOR marker terminates each raster + call rc_up_60i (header, EOR_OFFSET, PRU_EOR(rp), 1) + + # The PRU containing the first data row + call rc_up_60i (header, FIRST_PRU_OFFSET, PRU_ROW_ONE(rp), 1) + + # Most significant 30 bits of the data min are used + call rc_up_60r (header, MIN_OFFSET, DATA_MIN(rp), 1) + + # Most significant 30 bits of the data max are used + call rc_up_60r (header, MAX_OFFSET, DATA_MAX(rp), 1) + + # Bits per pixel is unpacked and tested + call rc_up_60i (header, DATA_TYPE_OFFSET, BITS_PIXEL(rp), 1) + + switch (BITS_PIXEL(rp)) { + case 12,20,30,60: + ; + default: + call error (2, "Incorrect IPPS BITS_PIXEL") + } + + # Number of columns is unpacked and tested + call rc_up_60i (header, NCOLS_OFFSET, NCOLS(rp), 1) + if (NCOLS(rp) <= 0) + call error (2, "IPPS ncols <= 0") + + # Number of Cyber words per row must be integral # of PRU's + call rc_up_60i (header, NWORDS_OFFSET, WRDS_PER_ROW(rp), 1) + + if (mod (WRDS_PER_ROW(rp), LEN_PRU) != 0) + call error (2, "Invalid IPPS NWPR") + + # Number of rows is unpacked and tested + call rc_up_60i (header, NROWS_OFFSET, NROWS(rp), 1) + if (NROWS(rp) <= 0) + call error (2, "IPPS nrows <= 0") +end diff --git a/noao/mtlocal/cyber/rrcopy/rcrimage.x b/noao/mtlocal/cyber/rrcopy/rcrimage.x new file mode 100644 index 00000000..dc7ebcfb --- /dev/null +++ b/noao/mtlocal/cyber/rrcopy/rcrimage.x @@ -0,0 +1,173 @@ +include +include +include +include "rrcopy.h" + +# RC_READ_IMAGE -- reads the rcopy image row by row from the tape and writes +# the output image. At the completion of READ_IMAGE, the tape is positioned +# to the header record of the next tape raster. + +procedure rc_read_image (rd, out_fname, data_type, rp) + +int rd, data_type +char out_fname[SZ_FNAME] +pointer rp + +pointer sp, im, cyber_buf, spp_buf +int nchars_per_row, nbits_skip, nchars_to_skip, i +long clktime() + +int rc_read_cyber() +pointer immap(), impl2r() +errchk rc_skip_chars, immap, rc_ipps_to_iraf + +begin + # Allocate buffer for rcopy image pixels + call smark (sp) + nchars_per_row = WRDS_PER_ROW(rp) * NBITS_CYBER_WORD / NBITS_CHAR + call salloc (cyber_buf, nchars_per_row, TY_CHAR) + call salloc (spp_buf, nchars_per_row, TY_CHAR) + + # Map new iraf image and set up image header + im = immap (out_fname, NEW_IMAGE, 0) + IM_LEN(im, 1) = NCOLS(rp) + IM_LEN(im, 2) = NROWS(rp) + call strcpy (IPPS_ID(rp), IM_TITLE(im), SZ_IMTITLE) + IM_MIN(im) = DATA_MIN(rp) + IM_MAX(im) = DATA_MAX(rp) + + # Set optimum image pixel type + if (data_type == NOT_SET) { + switch (BITS_PIXEL(rp)) { + case 12: + IM_PIXTYPE(im) = TY_SHORT + case 20: + IM_PIXTYPE(im) = TY_REAL + case 30: + IM_PIXTYPE(im) = TY_REAL + case 60: + IM_PIXTYPE(im) = TY_REAL + default: + call error (3, "IPPS BITS_PIXEL is incorrect") + } + } else + IM_PIXTYPE(im) = data_type + IM_LIMTIME(im) = clktime (long(0)) + + # Loop over rows to read, reorder and convert pixels. + for (i=1; i <= NROWS(rp); i=i+1) { + if (rc_read_cyber (rd, Memc[cyber_buf], nchars_per_row) == EOF) + call error (4, "Unexpected EOT when reading image") + call rc_order_cyber_bits (Memc[cyber_buf], 1, Memc[spp_buf], + WRDS_PER_ROW(rp)) + call rc_ipps_to_iraf (Memc[spp_buf], Memr[impl2r(im,i)], NCOLS(rp), + BITS_PIXEL(rp)) + } + + # Skip from present position to end of rcopy raster + nbits_skip = ((PRU_EOR(rp) - PRU_ROW_ONE(rp)) * LEN_PRU - + (WRDS_PER_ROW(rp) * NROWS(rp))) * NBITS_CYBER_WORD + NBITS_EOR_MARK + + nchars_to_skip = nbits_skip / NBITS_CHAR + call rc_skip_chars (rd, nchars_to_skip) + + call imunmap (im) + call sfree (sp) +end + + +# RC_IPPS_TO_IRAF -- performs the conversion from Cyber pixels to IRAF pixels. +# Each row of the rcopy image is required to occupy an integral of Cyber +# PRU's, so the input buffer contains pixels plus filler. The entire +# buffer is converted and npix pixels are written to the output image. + +procedure rc_ipps_to_iraf (in_buf, iraf_real, npix, nbits_pixel) + +char in_buf[ARB] +real iraf_real[npix] +pointer iraf_int, sp +int nbits_pixel, npix, bit_offset, npix_unpk, npix_cyber_wrd +errchk rc_up_12, rc_up_20 + +begin + # Calculate and allocate (maximum) space needed on the stack. The + # number of pixels unpacked will always fill an integral number + # of Cyber words. A maximum of 4 extraneous pixels will be unpacked. + call smark (sp) + call salloc (iraf_int, npix + 4, TY_INT) + bit_offset = 1 + + switch (nbits_pixel) { + case 12: + npix_cyber_wrd = 5 + npix_unpk = ((npix + 4) / npix_cyber_wrd) * npix_cyber_wrd + call rc_up_12 (in_buf, bit_offset, Memi[iraf_int], npix_unpk) + call achtir (Memi[iraf_int], iraf_real, npix) + + case 20: + npix_cyber_wrd = 3 + npix_unpk = ((npix + 2) / npix_cyber_wrd) * npix_cyber_wrd + call rc_up_20 (in_buf, bit_offset, Memi[iraf_int], npix_unpk) + call achtir (Memi[iraf_int], iraf_real, npix) + + case 30: + call rc_up_30 (in_buf, bit_offset, iraf_real, npix) + + case 60: + call rc_up_60r (in_buf, bit_offset, iraf_real, npix) + + default: + call error (5, "Illegal IPPS #B/P") + } + + call sfree (sp) +end + + +# RC_SKIP_IMAGE -- skips over an RCOPY raster once the header has been +# read. When SKIP_IMAGE returns, the tape is positioned to the first +# record of the next tape image. + +procedure rc_skip_image (rd, rp) + +int rd +pointer rp +int nchars_to_skip +errchk rc_skip_chars + +begin + # Calculate number of chars in image + nchars_to_skip = ((PRU_EOR(rp) - PRU_ROW_ONE(rp)) * + LEN_PRU * NBITS_CYBER_WORD + NBITS_EOR_MARK ) / NBITS_CHAR + call rc_skip_chars (rd, nchars_to_skip) +end + +# RC_SKIP_CHARS -- positions the tape by skipping the requested number of chars. + +procedure rc_skip_chars (rd, nchars_to_skip) + +int rd, nchars_to_skip +pointer sp, dummy +int nblks_read, nchars_remaining, i +int rc_read_cyber() + +begin + call smark (sp) + call salloc (dummy, SZ_TAPE_BLK, TY_CHAR) + + # Calculate the number of full blocks to skip + nblks_read = nchars_to_skip / SZ_TAPE_BLK + nchars_remaining = mod (nchars_to_skip, SZ_TAPE_BLK) + + # Read from tape, a block at a time + for (i=1; i <= nblks_read; i=i+1) { + if (rc_read_cyber (rd, Memc[dummy], SZ_TAPE_BLK) == EOF) + call error (6, "Unexpected EOT when skipping image") + } + + # Read partial block from tape + if (rc_read_cyber (rd, Memc[dummy], nchars_remaining) == EOF) + call error (7, "Unexpected EOT when skipping image") + + call sfree (sp) +end diff --git a/noao/mtlocal/cyber/rrcopy/rrcopy.h b/noao/mtlocal/cyber/rrcopy/rrcopy.h new file mode 100644 index 00000000..9579b623 --- /dev/null +++ b/noao/mtlocal/cyber/rrcopy/rrcopy.h @@ -0,0 +1,41 @@ + +# Definitions for the Cyber RCOPY tape reader + +define NBITS_CHAR (NBITS_BYTE * SZB_CHAR) # Number of bits per char +define NBITS_CYBER_WORD 60 # Number of bits per Cyber word +define LEN_PRU 64 # Number of words per Cyber pru +define NBITS_PRU 3840 # Number of bits per Cyber pru +define NCHARS_NOISE (48 / NBITS_CHAR) # Nchars in a Cyber noise record +define NBITS_EOR_MARK 48 # Number of bits per eor marker +define SZ_HEADER ((64 * 60) / NBITS_CHAR) # Size in chars of IPPS header +define SZ_TAPE_BLK ((512 * 60) / NBITS_CHAR) # Size in chars of tape block +define SZ_BUFFER (SZ_TAPE_BLK + 100) # Size of tape buffer for read +define SZ_IPPS_ID 127 # Max number of characters in ID +define MAX_RANGES 100 +define NOT_SET 0 # Flag for data_type not set +define BLANK 0.0 # Temporary value for blanks + +# Bit-offsets to IPPS header words + +define DATA_TYPE_OFFSET (16 * 60 + 1) # Offset to data_type (nbpp) +define NCOLS_OFFSET (17 * 60 + 1) # Offset to ncols (nppr) +define NWORDS_OFFSET (18 * 60 + 1) # Offet to nwords_per_row +define NROWS_OFFSET (20 * 60 + 1) # Offset to nrows +define FIRST_PRU_OFFSET (21 * 60 + 1) # Offset to 1st pru of raster +define MIN_OFFSET (31 * 60 + 1) # Offset to data min +define MAX_OFFSET (32 * 60 + 1) # Offset to data max +define EOR_OFFSET (44 * 60 + 1) # Offset to terminating pru + +# The IPPS raster descriptor structure RP: + +define LEN_RP 10 + SZ_IPPS_ID + 1 + +define BITS_PIXEL Memi[$1] +define PRU_EOR Memi[$1+1] +define WRDS_PER_ROW Memi[$1+2] +define PRU_ROW_ONE Memi[$1+3] +define NCOLS Memi[$1+4] +define NROWS Memi[$1+5] +define DATA_MIN Memr[P2R($1+6)] +define DATA_MAX Memr[P2R($1+7)] +define IPPS_ID Memc[P2C($1+10)] diff --git a/noao/mtlocal/cyber/rrcopy/rrcopy.x b/noao/mtlocal/cyber/rrcopy/rrcopy.x new file mode 100644 index 00000000..0a6b5de7 --- /dev/null +++ b/noao/mtlocal/cyber/rrcopy/rrcopy.x @@ -0,0 +1,212 @@ +include +include "rrcopy.h" + +.help rc_read_cyber +.nf ________________________________________________________________________ +RC_READ_CYBER -- Read binary chars from a file, ignoring short "noise" records. +Data is read in chunks from the file buffer passed as a bit stream into the +output buffer. (See also: read_dumpf and get_cyber_words.) The +file buffer is refilled as necessary by calling READ; the output buffer is +supplied by the calling procedure. Cyber noise records (48 bits) +are not transferred to the output buffer and so are ignored. READ_CYBER_INIT +must be called to initialize variables for CYBER_READ. Variables marking the +current position and top of the buffer are initialized. +.endhelp ___________________________________________________________________ + +int procedure rc_read_cyber (rd, out_buffer, maxch) + +int buf_pos, buf_top +int rd, maxch +char out_buffer[ARB] +char block_buf[SZ_BUFFER] +int nchars, chunk_size, nchars_read +int read(), rc_read_cyber_init() + +begin + for (nchars = 0; nchars < maxch; nchars = nchars + chunk_size) { + # See if it is necessary to transfer more data from the binary file + # to the file buffer. This will be necessary when all data from the + # file buffer has been moved to the output buffer. + + if (buf_pos >= buf_top) { + # Read the next non-noise record into block_buf; reset buf_pos + repeat { + nchars_read = read (rd, block_buf, SZ_BUFFER) + } until (nchars_read >= NCHARS_NOISE || nchars_read == EOF) + buf_pos = 1 + buf_top = nchars_read + } + + # The number of chars to output is the smaller of the number of + # characters requested or the number of characters left in the + # buffer + + if (nchars_read == EOF) + break + else + chunk_size = min (maxch - nchars, buf_top - buf_pos + 1) + + # Move data to output array, increment buffer offset + call amovc (block_buf[buf_pos], out_buffer[nchars+1], chunk_size) + buf_pos = buf_pos + chunk_size + } + + if (nchars == 0) + return (EOF) + else + return (nchars) + +entry rc_read_cyber_init () + + buf_pos = 1 + buf_top = 0 + return (OK) +end + + +.help rc_order_cyber_bits +.nf __________________________________________________________________________ +RC_ORDER_CYBER_BITS -- Convert raw data array from a 60-bit Cyber computer into +an SPP bit-array. The output SPP bit-array is a bit-packed array of 60-bit +Cyber words, i.e., bits 1-60 are word 1, bits 61-120 are word 2, and so on. +The least significant Cyber bit is bit 1 in each output word and the most +significant bit is bit 60. [MACHDEP]. + +When the Cyber outputs an array of 60 bit Cyber words it moves a stream of +bits into output bytes, filling however many bytes as necessary to output a +given number of Cyber words. The most significant bits are output first, +i.e., bit 60 of the first word is moved to bit 8 of the first output byte, +bit 59 is moved to bit 7 of the first output byte, and so on. If effect the +Cyber byte flips each 60-bit word. + +To deal with Cyber words as an ordered bit stream we must reorder the bytes +and bits so that the least significant bits are first. This function is +performed by the primitives CYBOOW and CYBOEW (order odd/even Cyber word) +for individual 60-bit Cyber words. A portable (and less efficient) version +of order_cyber_bits is available which does not use these primitives. +.endhelp _____________________________________________________________________ + + +procedure rc_order_cyber_bits (raw_cyber, first_cyber_word, bit_array, + ncyber_words) + +char raw_cyber[ARB] # raw Cyber array (e.g. from tape) +int first_cyber_word # first 60-bit Cyber word to be unpacked +char bit_array[ARB] # output bit-array +int ncyber_words # number of Cyber words to unpack + +bool odd_word +int word, inbit, outbit + +begin + odd_word = (mod (first_cyber_word, 2) == 1) + inbit = (first_cyber_word - 1) * NBITS_CYBER_WORD + 1 + outbit = 1 + + do word = 1, ncyber_words { + # Call odd or even primitive to reorder bits and move 60-bit + # ordered Cyber word to the output array. + + if (odd_word) { + call cyboow (raw_cyber, inbit, bit_array, outbit) + odd_word = false + } else { + call cyboew (raw_cyber, inbit, bit_array, outbit) + odd_word = true + } + + inbit = inbit + NBITS_CYBER_WORD + outbit = outbit + NBITS_CYBER_WORD + } +end + + +# The portable version of order_cyber_bits follows. +#.help order_cyber_bits +#.nf __________________________________________________________________________ +#ORDER_CYBER_BITS -- Convert a raw data array from a 60-bit Cyber computer into +#an SPP bit-array. The output SPP bit-array is a bit-packed array of 60-bit +#Cyber words, i.e., bits 1-60 are word 1, bits 61-120 are word 2, and so on. +#The least significant Cyber bit is bit 1 in each output word and the most +#significant bit is bit 60. [MACHDEP]. +# +#The byte stream from the Cyber contains bits 53-60 of the first word in the +#first byte, bits 45-52 in the second byte, and so on (most significant bytes +#first). In essence we swap the order of the 7 8-bit bytes and the 4-bit half +#byte in each 60 bit word. The bits in each byte are in the correct order. +# +#Each successive pair of Cyber words fits into 15 bytes. Byte 8 contains the +#last 4 bits of word 1 in the most signficant half of the byte and the first +#4 bits of word 2 in the first half of the byte. In each 60 bit word we must +#move bit segments (bytes or half bytes) as follows (for the VAX): +# +#Odd words (from N*60 bit-offset): +# [from] [to] [nbits] +# 1 53 8 +# 9 45 8 +# 17 37 8 +# 25 29 8 +# 33 21 8 +# 41 13 8 +# 49 5 8 +# 61 1 4 +# +#Even words (from N*60 bit-offset): +# [from] [to] [nbits] +# -3 57 4 +# 5 49 8 +# 13 41 8 +# 21 33 8 +# 29 25 8 +# 37 17 8 +# 45 9 8 +# 53 1 8 +#.endhelp _____________________________________________________________________ +# +#define NBITS_PER_WORD 60 +#define NSEGMENTS 8 +# +# +#procedure order_cyber_bits (raw_cyber, first_cyber_word, bit_array, +# ncyber_words) +# +#char raw_cyber[ARB] # raw Cyber array (e.g. from tape) +#int first_cyber_word # first 60-bit Cyber word to be unpacked +#char bit_array[ARB] # output bit-array +#int ncyber_words # number of Cyber words to unpack +# +#int word, inword, inbit, outbit, temp, i +#int o_from[NSEGMENTS], o_to[NSEGMENTS], o_nbits[NSEGMENTS] +#int e_from[NSEGMENTS], e_to[NSEGMENTS], e_nbits[NSEGMENTS] +#int bitupk() +# +#data o_from / 1, 9,17,25,33,41,49,61/ # odd words +#data o_to /53,45,37,29,21,13, 5, 1/ +#data o_nbits / 8, 8, 8, 8, 8, 8, 8, 4/ +#data e_from /-3, 5,13,21,29,37,45,53/ # even words +#data e_to /57,49,41,33,25,17, 9, 1/ +#data e_nbits / 4, 8, 8, 8, 8, 8, 8, 8/ +# +#begin +# do word = 1, ncyber_words { +# inword = first_cyber_word + word - 1 +# inbit = (inword - 1) * NBITS_PER_WORD +# outbit = ( word - 1) * NBITS_PER_WORD +# +# # Move bits to the output bit array. Segment list used depends +# # on whether the word is an odd or even word. This code will work +# # even if the caller only wishes to order a single word. +# +# if (mod (inword,2) == 1) { +# do i = 1, NSEGMENTS { +# temp = bitupk (raw_cyber, inbit + o_from[i], o_nbits[i]) +# call bitpak (temp, bit_array, outbit + o_to[i], o_nbits[i]) +# } +# } else { +# do i = 1, NSEGMENTS { +# temp = bitupk (raw_cyber, inbit + e_from[i], e_nbits[i]) +# call bitpak (temp, bit_array, outbit + e_to[i], e_nbits[i]) +# } +# } +# } +#end diff --git a/noao/mtlocal/cyber/rrcopy/semicode.doc b/noao/mtlocal/cyber/rrcopy/semicode.doc new file mode 100644 index 00000000..a7ad514c --- /dev/null +++ b/noao/mtlocal/cyber/rrcopy/semicode.doc @@ -0,0 +1,310 @@ +Semicode for Cyber RCOPY Reader; Frozen at "detailed semicode" stage, Jan 84. + + +# rcopy tape descriptor + +struct rcopy { + + real data_min # minimum data value + real data_max # maximum data value + int nrows # number of rows in ipps raster + int ncols # number of columns in ipps raster + int data_type # number of bits per pixel in the ipps raster + int pru_eor # pru position of level zero eor + int wrds_per_row # number of 60-bit words per row; each row + # occupies an integral number of 64word pru's. + int pru_row_one # relative pru ordinal of the first row of raster. + char ipps_id # id string of the ipps raster +} + +procedure t_rrcopy (rcopy_file, raster_list, iraf_file) + + +begin + # get input filename and open tape drive + rcopy_fd = open (rcopy_file) + + # get output filename if it will be needed + if (make_image = yes) + outfile = fetch root name of output file + + # expand list of rasters to be read from tape + if (decode_ranges (raster_list, range, MAX_RANGES, nfiles) == ERR) + call error (0, "Illegal raster number list") + + raster_number = 0 + tape_pos = 1 + while (get_next_number (range, raster_number, != EOF) { + # position tape to first record of raster_number + if (tape_pos != raster_number) { + iferr { + stat = position_rcopy (rcopy_fd, tape_pos ,raster_number) + if (stat = EOF) + return + } then + call ERRACT (EA_FATAL) + } + iraffile = generate output filename from root + iferr { + (stat = read_rcopy(rcopy_fd, iraffile, tape_pos)) + if (stat = EOF) + return + } then { + call ERRACT (EA_WARN) + skip over rcopy raster + } + } + +end + +int procedure position_rcopy (rcopy_fd, tape_pos, raster_number) + +begin + + nrasters_skip = raster_number - tape_pos + for (i=1; i<=nrasters_skip; i=i+1) { + stat = read_header (rcopy_fd, rcopy) + if (stat = EOF) + return (EOF) + call read_header (rcopy_fd, rcopy) + call skip_image (rcopy_fd, rcopy) + } +end + +int procedure read_rcopy (rcopy_fd, iraffile, tape_pos) + +begin + + # Read ipps raster header from rcopy tape + stat = read_header (rcopy_fd, rcopy) + if (stat = EOF) + return (EOF) + + # Print header information from rcopy tape + if (print_header) + call print_header (rcopy) + + # read image data if desired + if (make_image) + call read_image(rcopy_fd, iraffile, rcopy) + + # skip image if not to be read + else + call skip_image (rcopy_fd, rcopy) + + # increment tape position marker + tape_pos = tape_pos + 1 + + +end + +int procedure read_header(rcopy_fd, rcopy) + +begin + # structure rcopy contains decoded header + # Read ipps header (64 60-bit words = 240 chars) as a + # bit stream into temp. + NBITS_CHAR = SZB_CHAR * NBITS_BYTE + NBITS_CYBER_WORD = 60 + NWORDS_CYBER_PRU = 64 + SZ_HEADER = 240 + stat = read_tape (rcopy_fd, raw_header, SZ_HEADER) + if (stat = EOF)) + return (EOF) + + # Unpack bit stream and fill structure rcopy + call unpack_header (raw_header, rcopy) + + # skip to first row of raster + if (rcopy.pru_row_one not equal to 1) { + nchars_to_skip = (rcopy.pru_row_one - 1) * NWORDS_CYBER_PRU + * NBITS_CYBER_WORD / NBITS_CHAR + call skip_chars (rcopy_fd, nchars_to_skip) + } +end + +procedure read_image (rcopy_fd, iraffile, rcopy) + +begin + + # map new iraf image and set up image header + im = immap (iraffile, NEW_IMAGE) + im.im_len[1] = rcopy.ncols + im.im_len[2] = rcopy.nrows + im.im_title = rcopy.ipps_id + im.im_min = rcopy.data_min + im.im_max = rcopy.data_max + + if (user hasn's supplied data type) { + switch (rcopy.data_type) { + case 12 : + im.im_pixtype = short + case 20 : + im.im_pixtype = real + case 30 : + im.im_pixtype = real + default: + call error (0, "error in data_type") + } + } + + # Calculate number of chars per row; this will always be an + # integral number of chars because rcopy.words_per_row + # is always divisible by 64. + NBITS_CHAR = SZB_CHAR * NBITS_BYTE + NBITS_CYBER_WORD = 60 + NWORDS_CYBER_PRU = 64 + num_chars_per_row = rcopy.wrds_per_row * NBITS_CYBER_WORD / NBITS_CHAR + + # Loop to read rcopy raster line by line into buf1, then + # convert ipps pixels to iraf pixels. + for(i=1; i<=rcopy.nrows; i=i+1) { + stat = read a single row from the internal buffer + if (stat = EOF) { + close imagefile already accumulated + call error (0, "unexpected EOF at row# ") + } + call ipps_to_iraf (buf1, Memr[impl2r(im,i)], rcopy) + } + + # Read until header of next raster into dummy buffer. + # Calculate offset in chars from present position to eor + nchars_to_skip = (((rcopy.eor - rcopy.pru_row_one) * NWORDS_CYBER_PRU + + (rcopy.words_per_row * rcopy.nrows)) * NBITS_CYBER_WORD + 48) + / NBITS_CHAR + + call skip_chars (rcopy_fd, nchars_to_skip) +end + +procedure ipps_to_iraf (buf1, buf2, rcopy) + +begin + + # convert (rcopy.ncols * rcopy.data_type ) bits from buf1. + # This is the number of bits per row that actually represent + # pixels. buf1 contains pixels plus filler. + + switch (rcopy.data_type) { + case 12: + call unpack12 (buf1, bit_offset, buf2, npix) + case 20: + call unpack20 (buf1, bit_offset, buf2, npix) + case 30: + call unpack30 (buf1, bit_offset, buf2, npix) + default: + call error (0, "illegal ipps #b/p") + } +end + +procedure skip_image(rcopy_fd, rcopy) + +begin + + # Calculate number of chars in image + nchars_to_skip = ((rcopy.eor - rcopy.pru_row_one ) * NWORDS_CYBER_PRU * + NBITS_CYBER_WORD + 48 ) / NBITS_CHAR + call skip_chars (rcopy_fd, nchars_to_skip) + +end + +procedure print_header (rcopy) + +begin + # print header information from rcopy tape + print, rcopy.ipps_id, rcopy.ncols, rcopy.nrows, rcopy.data_min, + rcopy.data_max, rcopy.ipps_data_type +end + +procedure unpack_header (raw_header, rcopy) + +begin + + get from raw_header: nrows, ncols, data_type, ipps_id, pru_eor, + wrds_per_row, pru_row_one, data_min, data_max + + if (rcopy.data_type != 12, 20, or 30) + call error (0, "invalid ipps #b/p") + if (nrows or ncols or pru_row_one !> 0) + call error (0, "invalid ipps raster") + if (wrds_per_row not divisible by 64) + call error (0, "invalid ipps raster") + +end + +procedure skip_chars (rcopy_fd, nchars_to_skip) + +begin + + # calculate the number of chars in a tape block = (512 * 60) / 16 + # This is the number of chars to be read at one time. + SZ_TAPE_BLK = 1920 + + # calculate number of blocks to skip (type int) + nblks_read = nchars_to_skip / SZ_TAPE_BLK + nchars_remaining = mod (nchars_to_skip, SZ_TAPE_BLK) + + # read chars into dummy buffer + for (i=1; i<=nblks_read; i=i+1) + stat = read_tape (rcopy_fd, dummy, SZ_TAPE_BLK) + + stat = read_tape (rcopy_fd, dummy, nchars_remaining) + + # confirm reaching eor + if (searching for eor) { + reread last 3 chars and compare with level zero eor marker + if (eor not confirmed) + call error (0, "attempted skip to eor unsuccessful") + } + +end + + +.help rrcopy 2 "Program Structure" +.sh +RCOPY Structure Chart + +.nf +procedure t_rrcopy () +# Returns when file list is satisfied or if EOT is encountered + + int procedure read_header (rcopy_fd, rcopy) + #returns EOF, ERR or OK + + int procedure read (rcopy_fd, raw_header, SZ_HEADER) + #returns EOF or OK + + int procedure unpack_header (raw_header, rcopy) + #returns ERR or OK + + int procedure skip_image (rcopy_fd, rcopy) + #returns ERR, EOF or OK + + int procedure skip_chars (rcopy_fd, nchars_to_skip) + #returns EOF, ERR or OK + + procedure read_rcopy (rcopy_fd, iraffile) + #aborts if output filename exists and noclobber is set + + int procedure read_header (rcopy_fd, rcopy) + #returns EOF, ERR or OK + + int procedure read (rcopy_fd, raw_header, SZ_HEADER) + #returns EOF or OK + + int procedure unpack_header (raw_header, rcopy) + #returns ERR or OK + + procedure print_header (rcopy) + + procedure read_image (rcopy_rd, iraffile, rcopy) + #returns EOF, or OK + + int procedure read (rcopy_rd, buf1, num_chars_per_row) + #returns EOF or OK + + int procedure ipps_to_iraf + #returns ERR or OK + + int procedure skip_image (rcopy_rd, rcopy) + #returns EOF, ERR or OK +.endhelp diff --git a/noao/mtlocal/cyber/rrcopy/t_rrcopy.x b/noao/mtlocal/cyber/rrcopy/t_rrcopy.x new file mode 100644 index 00000000..9bb83885 --- /dev/null +++ b/noao/mtlocal/cyber/rrcopy/t_rrcopy.x @@ -0,0 +1,147 @@ +include +include +include +include "rrcopy.h" + +# T_RRCOPY -- The main procedure for the RCOPY reader. The RCOPY reader +# converts IPPS rasters written in RCOPY format to IRAF images. All IPPS +# rasters on an RCOPY tape are in a single file. Each raster header must +# be read before the image can be either skipped or read. T_RRCOPY gets +# parameters from the cl and decodes the string of rasters to be read. It +# then calls READ_HEADER for each raster on the tape. The header information +# is printed if print_header=true. If make_image = true, the image is +# converted to an IRAF image by READ_IMAGE. Otherwise, the image is skipped +# with SKIP_IMAGE. T_RRCOPY terminates when the raster list is depleted or +# the tape is at EOT. +# +# Modified 26-JULY-88 to allow for multiple rcopy files on a single tape. +# This allows for rcopy format data to be archived in multiple files on +# one tape. The task is still run once per input file. The user is queried +# (hidden parameter) for the data file to be read. The tape file is actually +# datafile + 1 because of the ANSI label on each rrcopy tape. (ShJ) + +procedure t_rrcopy () + +pointer sp, rp +bool make_image, print_header, bad_header +char rcopy_file[SZ_FNAME], iraf_file[SZ_FNAME] +char out_fname[SZ_FNAME], raster_list[SZ_LINE] +int rd, ras_number, current_ras, nras, stat, tapefile +int ranges[3, MAX_RANGES], data_type, init + +bool clgetb() +char clgetc() +int get_data_type(), position_rcopy(), rc_read_cyber_init(), clgeti() +int mtopen(), decode_ranges(), get_next_number(), rc_header_read(), strlen() +int mtfile() + +begin + # Allocate space on stack for program data structure + call smark (sp) + call salloc (rp, LEN_RP, TY_STRUCT) + + # Get input filename and open tape drive to second file, skipping label + call clgstr ("rcopy_file", rcopy_file, SZ_FNAME) + if (mtfile (rcopy_file) == YES) { + tapefile = clgeti ("datafile") + 1 + call mtfname (rcopy_file, tapefile, rcopy_file, SZ_FNAME) + } + rd = mtopen (rcopy_file, READ_ONLY, SZ_BUFFER) + init = rc_read_cyber_init() + + # Get output root filename if it will be needed + make_image = clgetb ("make_image") + if (make_image) { + call clgstr ("iraf_file", iraf_file, SZ_FNAME) + data_type = get_data_type (clgetc ("data_type")) + if (data_type == ERR) + data_type = NOT_SET + } + + # Set options + print_header = clgetb ("print_header") + + # Expand list of rasters to be read from tape + call clgstr ("raster_list", raster_list, SZ_LINE) + if (decode_ranges (raster_list, ranges, MAX_RANGES, nras) == ERR) + call error (0, "Illegal raster number list") + + ras_number = 0 + current_ras = 1 + while (get_next_number (ranges, ras_number) != EOF) { + # Position tape to first record of ras_number + if (current_ras != ras_number) { + iferr (stat = position_rcopy (rd, current_ras, ras_number, rp)) + call erract (EA_FATAL) + if (stat == EOF) + break + } + + # Assume header is good + bad_header = false + iferr { + stat = rc_header_read (rd, rp) + } then { + # Error reading header; will attempt to skip raster + bad_header = true + call erract (EA_WARN) + } + + if (stat == EOF) { + call printf ("\nRCOPY tape at End of Tape\n") + break + } + + if (print_header) + call rc_list_header (rp, ras_number) + call flush (STDOUT) + + if (make_image && ! bad_header) { + # Generate output filename + call strcpy (iraf_file, out_fname, SZ_FNAME) + if (nras > 1) { + call sprintf (out_fname[strlen(out_fname)+1], SZ_FNAME, + "%03d") + call pargi (ras_number) + } + iferr (call rc_read_image (rd, out_fname, data_type, rp)) + call erract (EA_FATAL) + } else + iferr (call rc_skip_image (rd, rp)) + call erract (EA_FATAL) + + # Increment tape position + current_ras = current_ras + 1 + } + + # Return space allocated for rp, close tape unit + call close (rd) + call sfree (sp) +end + + +# POSITION_RCOPY -- Position the tape to the first +# record of the next raster to be read. Each raster header must +# be read; each image can then be skipped. + +int procedure position_rcopy (rd, current_ras, ras_number, rp) + +int rd, current_ras, ras_number +pointer rp +int nras_skip, i, stat +int rc_header_read() +errchk rc_skip_image + +begin + nras_skip = ras_number - current_ras + for (i=1; i <= nras_skip; i=i+1) { + stat = rc_header_read (rd, rp) + if (stat == EOF) { + call printf ("Cannot position RCOPY tape beyond EOF\n") + return (EOF) + } + call rc_skip_image (rd, rp) + current_ras = current_ras + 1 + } + return (OK) +end diff --git a/noao/mtlocal/cyber/t_ldumpf.x b/noao/mtlocal/cyber/t_ldumpf.x new file mode 100644 index 00000000..880c0290 --- /dev/null +++ b/noao/mtlocal/cyber/t_ldumpf.x @@ -0,0 +1,220 @@ +include +include "cyber.h" + +# T_LDUMPF -- list permanent files stored on a DUMPF tape. +# The permanent file name, owner id, cycle number, creation date, +# last attach and last alteration dates are listed for specified files +# of a Cyber DUMPF format tape. + +procedure t_ldumpf () + +pointer sp, dmp +char dumpf_file[SZ_FNAME], file_list[SZ_LINE], in_fname[SZ_FNAME] +int file_number, ranges[3, MAX_RANGES], nfiles +int read_pf_table(), get_next_number(), decode_ranges() +int mtfile() + +begin + # Allocate space for program data structure. + call smark (sp) + call salloc (dmp, LEN_DMP, TY_STRUCT) + + # Get parameters; get file_list only if dump_file is a general + # tape name. + call clgstr ("dumpf_file", dumpf_file, SZ_FNAME) + if (mtfile (dumpf_file) == YES) + call clgstr ("file_list", file_list, SZ_LINE) + else + call strcpy ("1", file_list, SZ_LINE) + if (decode_ranges (file_list, ranges, MAX_RANGES, nfiles) == ERR) + call error (0, "Illegal file number list") + + # For each file in file_list call read_pf_table. + file_number = 0 + while (get_next_number (ranges, file_number) != EOF) { + if (mtfile (dumpf_file) == YES) + call mtfname (dumpf_file, file_number + 1, in_fname, SZ_FNAME) + else + call strcpy (dumpf_file, in_fname, SZ_FNAME) + if (read_pf_table (in_fname, file_number, dmp) == EOF) { + call printf ("End of DUMPF tape\n") + call sfree (sp) + return + } + } + call sfree (sp) +end + + +# READ_PF_TABLE -- reads and prints out the Cyber permanent file information. + +int procedure read_pf_table (dumpf_file, file_num, dmp) + +char dumpf_file[SZ_FNAME] +int file_num +pointer dmp +int dump_fd, init +pointer dump_buf, spp_buf, sp +int mtopen(), read_cyber(), bitupk(), read_cyber_init() + +begin + # Allocate buffer space for the Permanent File Table + call smark (sp) + call salloc (dump_buf, SZ_PFT, TY_CHAR) + call salloc (spp_buf, SZ_PFT, TY_CHAR) + + # Open input DUMPF tape + dump_fd = mtopen (dumpf_file, READ_ONLY, SZ_TAPE_BUFFER) + init = read_cyber_init() + + # Read File Header, Permanent File Catalog and Permanent File Table + if (read_cyber (dump_fd, Memc[dump_buf], SZ_PFT) == EOF) { + call sfree (sp) + call close (dump_fd) + return (EOF) + } + + # Reorder Cyber bits into packed SPP bit array + call order_cyber_bits (Memc[dump_buf], 1, Memc[spp_buf], LEN_PFT) + + # Get number of characters in permanent file name + NCHARS_PFN(dmp) = bitupk (Memc[spp_buf], NCHARS_OFFSET, 6) + + # Unpack file information from Permanent File Table + call unpk_pf_info (Memc[spp_buf], dmp) + + # Print Permanent File information + call print_pf_info (file_num, dmp) + call flush (STDOUT) + + # Return buffer, close tape file + call sfree (sp) + call close (dump_fd) + return (OK) +end + + +# DECIPHER_DC -- An ascii character string is decoded from an input +# bit stream. An offset into the bit stream and the number of characters +# to unpack are input. + +procedure decipher_dc (inbuf, bit_offset, nchars, outbuf) + +char inbuf[ARB], outbuf[nchars + 1] +int offset, nchars, ip, op, temp, bit_offset +int bitupk() + +begin + op = 1 + offset = bit_offset + do ip = 1, nchars { + temp = bitupk (inbuf, offset, NBITS_DC) + if (temp == 55B) { + # Blank + offset = offset - NBITS_DC + next + } + else { + call display_code (temp, outbuf[op]) + op = op + 1 + offset = offset - NBITS_DC + } + } + outbuf[op] = EOS +end + + +# UNPK_PF_INFO -- unpacks words from the Permanent File Information Table +# and fills program data structure dmp. + +procedure unpk_pf_info (inbuf, dmp) + +char inbuf[SZ_PFT] +pointer dmp +int creation, attach, alteration +int bitupk() + +begin + # Extract Permanent File Name + call decipher_dc (inbuf, PFN_OFFSET, NCHARS_PFN(dmp), PFN(dmp)) + + # Extract ID + call decipher_dc (inbuf, PF_ID_OFFSET, SZ_PF_ID, ID(dmp)) + + # Extract cycle number + CY(dmp) = bitupk (inbuf, CY_OFFSET, NBITS_CY) + + # Extract creation, last_attach and last_alteration dates which are + # written in 18-bits as "yyddd". + creation = bitupk (inbuf, CREATE_OFFSET, NBITS_DATE) + call ld_get_date (creation, D_CREATE(dmp), M_CREATE(dmp), Y_CREATE(dmp)) + + attach = bitupk (inbuf, ATTACH_OFFSET, NBITS_DATE) + call ld_get_date (attach, D_ATTACH(dmp), M_ATTACH(dmp), Y_ATTACH(dmp)) + + alteration = bitupk (inbuf, ALTER_OFFSET, NBITS_DATE) + call ld_get_date (alteration, D_ALTER(dmp), M_ALTER(dmp), Y_ALTER(dmp)) +end + + +# GET_DATE -- The day, month and year are decoded from the lower 18 bits +# of the input integer. The input format is "yyddd"; three integers are +# returned as arguments: day, month, year. + +procedure ld_get_date (yyddd, day, month, year) + +int yyddd, day, month, year, ddd +int days_in_month[12], i, bitupk() + +data (days_in_month[i], i=1,9) /31, 28, 31, 30, 31, 30, 31, 31, 30/ +data (days_in_month[i], i=10, 12) /31, 30, 31/ + +begin + year = bitupk (yyddd, 10, 9) + 1900 + ddd = bitupk (yyddd, 1, 9) + + # Leap year tests + if (mod (year, 4) == 0) + days_in_month[2] = 29 + if (mod (year, 100) == 0) + days_in_month[2] = 28 + if (mod (year, 400) == 0) + days_in_month[2] = 29 + + for (i=1; i<=12; i=i+1) { + if (ddd <= days_in_month[i]) + break + else + ddd = ddd - days_in_month[i] + } + + month = i + day = ddd +end + + +# PRINT_PF_INFO -- information from the permanent file table is printed. + +procedure print_pf_info (file_num, dmp) + +pointer dmp +int file_num + +begin + call printf ("\n[%d]%6tPFN= %s, ID= %s, CY= %d\n") + call pargi (file_num) + call pargstr (PFN(dmp)) + call pargstr (ID(dmp)) + call pargi (CY(dmp)) + call printf ("%6tCreation: %d/%d/%d, Last Alteration: %d/%d/%d, ") + call pargi (M_CREATE(dmp)) + call pargi (D_CREATE(dmp)) + call pargi (Y_CREATE(dmp)) + call pargi (M_ALTER(dmp)) + call pargi (D_ALTER(dmp)) + call pargi (Y_ALTER(dmp)) + call printf ("Last Attach: %d/%d/%d\n") + call pargi (M_ATTACH(dmp)) + call pargi (D_ATTACH(dmp)) + call pargi (Y_ATTACH(dmp)) +end diff --git a/noao/mtlocal/cyber/t_rdumpf.x b/noao/mtlocal/cyber/t_rdumpf.x new file mode 100644 index 00000000..9d780b9e --- /dev/null +++ b/noao/mtlocal/cyber/t_rdumpf.x @@ -0,0 +1,162 @@ +include +include +include +include "cyber.h" + +# T_RDUMPF-- The main procedure for the DUMPF reader. Permanent files +# containing IPPS rasters are read in dumpf format and optionally +# converted to IRAF images. Each permanent file is a seperate tape file; +# the IPPS rasters are sequentially stored in the permanent file, seperated +# by "zero length PRU's". The first 48 words of each permanent file +# contain the Cyber permanent file table, catalog and file header information +# for the file. This information is listed with task LDUMPF. For each +# file in file_list, the file is opened. Then for each raster in the file, +# READ_HEADER must be called, followed by either READ_IMAGE or SKIP_IMAGE. + +procedure t_rdumpf() + +pointer sp, dt, dummy +bool make_image, print_header, bad_header +char dumpf_file[SZ_FNAME], iraf_file[SZ_FNAME], file_list[SZ_LINE] +char out_fname[SZ_FNAME], raster_list[SZ_LINE], in_fname[SZ_FNAME] +int fd, file_number, ras_number, current_ras +int stat, nfile, nras +int rasters[3, MAX_RANGES], files[3, MAX_RANGES], data_type + +bool clgetb() +char clgetc() +int get_data_type(), strlen(), mtfile() +int get_cyber_words() +int get_cyber_words_init(), read_dumpf_init(), position_dumpf() +int mtopen(), decode_ranges(), get_next_number(), cy_read_header() + +begin + call fseti (STDOUT, F_FLUSHNL, YES) + + # Allocate space for program data structure and buffers + call smark (sp) + call salloc (dt, LEN_DT, TY_STRUCT) + call salloc (dummy, NINT_CYBER_WRD * LEN_PFT, TY_INT) + + # Get paramters from cl + call clgstr ("dumpf_file", dumpf_file, SZ_FNAME) + if (mtfile (dumpf_file) == YES) + call clgstr ("file_list", file_list, SZ_LINE) + else + call strcpy ("1", file_list, SZ_LINE) + if (decode_ranges (file_list, files, MAX_RANGES, nfile) == ERR) + call error (0, "Illegal file list") + + call clgstr ("raster_list", raster_list, SZ_LINE) + if (decode_ranges (raster_list, rasters, MAX_RANGES, nras) == ERR) + call error (1, "Illegal raster list") + + print_header = clgetb ("print_header") + make_image = clgetb ("make_image") + if (make_image) { + call clgstr ("iraf_file", iraf_file, SZ_FNAME) + data_type = get_data_type (clgetc ("data_type")) + if (data_type == ERR) + data_type = NOT_SET + } + + # Expand file_list and open dumpf_file + file_number = 0 + while (get_next_number (files, file_number) != EOF) { + + # Get the file name and open file. + if (mtfile (dumpf_file) == YES) + call mtfname (dumpf_file, file_number + 1, in_fname, SZ_FNAME) + else + call strcpy (dumpf_file, in_fname, SZ_FNAME) + fd = mtopen (in_fname, READ_ONLY, SZ_TAPE_BUFFER) + + # Position to first IPPS raster in file, skipping Cyber PFT etc. + stat = get_cyber_words_init() + stat = read_dumpf_init() + + if (get_cyber_words (fd, Memi[dummy], LEN_PFT) == EOF) { + call printf ("DUMPF Tape at EOF\n") + call close (fd) + call sfree (sp) + return + } + + + ras_number = 0 + current_ras = 1 + while (get_next_number (rasters, ras_number) != EOF) { + # Position to first record of ras_number + if (current_ras != ras_number) { + iferr (stat = position_dumpf (fd, current_ras, ras_number, + dt)) + call erract (EA_FATAL) + if (stat == EOF) + break + } + + bad_header = false + iferr { + stat = cy_read_header (fd, dt) + } then { + # Error reading header; will attempt to skip raster + bad_header = true + call erract (EA_WARN) + } + + if (stat == EOF) { + call printf ("DUMPF Tape at End of File%d\n\n") + call pargi (file_number) + break + } + + if (print_header) + call cy_list_header (dt, file_number, ras_number) + + if (make_image && ! bad_header) { + call strcpy (iraf_file, out_fname, SZ_FNAME) + if (nras > 1 || nfile > 1) { + call sprintf (out_fname[strlen(out_fname)+1], SZ_FNAME, + "%d.%03d") + call pargi (file_number) + call pargi (ras_number) + } + iferr (call read_ipps_rows (fd, out_fname, data_type, dt)) + call erract (EA_FATAL) + } else + iferr (call cy_skip_image (fd, dt)) + call erract (EA_FATAL) + + current_ras = current_ras + 1 + } + call close (fd) + } + call sfree (sp) + return +end + + +# POSITION_DUMPF -- Position the tape to the first +# record of the next raster to be read. Each raster header must +# be read; each image can then be skipped. + +int procedure position_dumpf (rd, current_ras, ras_number, dt) + +int rd, current_ras, ras_number +pointer dt +int nras_skip, i +int cy_read_header() +errchk cy_skip_image + +begin + nras_skip = ras_number - current_ras + for (i=1; i <= nras_skip; i=i+1) { + if (cy_read_header (rd, dt) == EOF) { + call printf ("Cannot position DUMPF tape beyond EOF\n") + return (EOF) + } + call cy_skip_image (rd, dt) + current_ras = current_ras + 1 + } + return (OK) +end diff --git a/noao/mtlocal/cyber/t_ridsfile.x b/noao/mtlocal/cyber/t_ridsfile.x new file mode 100644 index 00000000..c7179763 --- /dev/null +++ b/noao/mtlocal/cyber/t_ridsfile.x @@ -0,0 +1,516 @@ +include +include +include +include +include "cyber.h" + + +# T_RIDSFILE __ code for the DUMPF IDSFILE reader. IDS records in an IDSFILE +# are read from a Cyber DUMPF tape and optionally converted to IRAF images. +# IDS records are not written sequentially in the IDSFILE, so, each record +# must be read and then checked against the list of "record_numbers" to +# see if the user requested the record to be read. The procedure terminates +# when the requested number of records has been read or EOF is encountered. +# The IDS trailer information is printed in either a short or long form; +# the pixel values can also be printed. + +procedure t_ridsfile() + +pointer sp, cp +char in_fname[SZ_FNAME], dumpf_file[SZ_FNAME] +int file_ordinal + +int mtfile(), clgeti(), get_data_type(), btoi() +bool clgetb() +char clgetc() + +begin + # Allocate space for the control parameter descriptor structure + call smark (sp) + call salloc (cp, LEN_CP, TY_STRUCT) + + call fseti (STDOUT, F_FLUSHNL, YES) + + # Get parameters from cl and generate input file name. If the input + # file is a tape, append the file_ordinal suffix, incremented by one + # to skip over the DUMPF tape label. + + call clgstr ("dumpf_file", dumpf_file, SZ_FNAME) + if (mtfile (dumpf_file) == YES) { + file_ordinal = clgeti ("file_ordinal") + call mtfname (dumpf_file, file_ordinal + 1, in_fname, SZ_FNAME) + } else + call strcpy (dumpf_file, in_fname, SZ_FNAME) + + LONG_HEADER(cp) = btoi (clgetb ("long_header")) + PRINT_PIXELS(cp) = btoi (clgetb ("print_pixels")) + call clgstr ("record_numbers", REC_NUMBERS(cp), SZ_LINE) + + # If an output image is to be written, get root output file name and + # output data type. + MAKE_IMAGE(cp) = btoi (clgetb ("make_image")) + if (MAKE_IMAGE(cp) == YES) { + call clgstr ("iraf_file", IRAF_FILE(cp), SZ_FNAME) + DATA_TYPE(cp) = get_data_type (clgetc ("data_type")) + if (DATA_TYPE(cp) == ERR) + DATA_TYPE(cp) = TY_REAL + } + call read_idsfile (in_fname, cp) + + call sfree (sp) +end + + +# READ_IDSFILE -- read and sort the index of record ranges. Call +# idsf_read_record for each record in each index range. + +procedure read_idsfile (in_fname, cp) + +char in_fname[SZ_FNAME] # Name of input file +pointer cp # Pointer to control parameter structure + +int records[3, MAX_RANGES], nrecs, i +pointer sp, pft, pru_buf +int fd, junk, index_buf[LEN_INDEX * NINT_CYBER_WRD], nranges, nids_read +int current_pru, n_index, next_pru, nrecords_to_read, npru_skip, n_rec +long sorted_index[LEN_INDEX] + +int mtopen(), get_cyber_words_init(), read_dumpf_init(), read_dumpf() +int get_cyber_words(), idsf_read_record(), decode_ranges() +errchk mtopen, read_dumpf, get_cyber_words, idsf_read_record +errchk sort_index, decode_ranges + +begin + # Allocate space for program data structure and buffers + call smark (sp) + call salloc (pft, NINT_CYBER_WRD * LEN_PFT, TY_INT) + call salloc (pru_buf, NINT_CYBER_WRD * LEN_PRU, TY_INT) + + # Open and initialize the tape file, and read the permanent file table + fd = mtopen (in_fname, READ_ONLY, SZ_TAPE_BUFFER) + junk = get_cyber_words_init() + junk = read_dumpf_init() + if (get_cyber_words (fd, Memi[pft], LEN_PFT) == EOF) { + call printf ("DUMPF tape at EOT\n") + call sfree (sp) + call close (fd) + return + } + + # Read and sort IDSFILE user index information. The first two + # pru's of this index are relevant. Up to 3 more pru's can + # follow, depending on the format of the idsfile. The code was + # modified 13Jan86 to read an old format tape of Paul Hintzen's + # and hopefully provide a general solution to the problem of + # different formats. + + if (read_dumpf (fd, index_buf, LEN_USER_INDEX)== EOF) { + call close (fd) + call error (1, "Unexpected EOF when reading index") + } + if (decode_ranges (REC_NUMBERS(cp), records, MAX_RANGES, junk) == ERR) + call error (2, "Error in record_numbers specification") + + call sort_index (index_buf, records, sorted_index, nranges, nrecs) + + # Loop over each range of records in the index. nids_read counts + # the number of records requested by the user that have been read. + # nrecords_to_read is the number of records in the current index range. + + nids_read = 0 + current_pru = 3 + for (n_index = 1; n_index <= nranges; n_index = n_index + 1) { + next_pru = sorted_index[n_index] / 1000 + nrecords_to_read = mod (sorted_index[n_index], 1000) + npru_skip = next_pru - current_pru + do i = 1, npru_skip { + if (read_dumpf (fd, Memi[pru_buf], LEN_PRU) == EOF) { + # At end of IDSFILE + call printf ("DUMPF tape at EOF\n") + break + } + } + + current_pru = current_pru + npru_skip + + # Loop over each record within the current range of records + for (n_rec = 1; n_rec <= nrecords_to_read; n_rec = n_rec + 1) { + if (nids_read >= nrecs) { + # No need to continue + call close (fd) + call sfree (sp) + return + } + + if (idsf_read_record (fd, records, nrecs, nids_read, + cp) == EOF) { + call close (fd) + call sfree (sp) + return + } + + + current_pru = current_pru + (LEN_IDS_RECORD / LEN_PRU) + } + } + + call close (fd) + call sfree (sp) +end + + +# IDSF_READ_RECORD -- reads a single idsrecord. If the record is in the +# set of records to be read, the record is processed and the count of requested +# records read is incremented. + +int procedure idsf_read_record (fd, records, nrecs, nids_read, cp) + +int fd # File descriptor of input file +int records[3, MAX_RANGES] # Array of ranges of records specified by user +int nrecs # Number of requested records found on tape +int nids_read # Number of requested records already read +pointer cp # Pointer to control parameter structure + +char out_fname[SZ_FNAME] +pointer sp, ids +int ids_buffer[LEN_IDS_RECORD * NINT_CYBER_WRD], this_record +int tape, scan +real pixels[NPIX_IDS_RECORD] + +bool is_in_range() +int read_dumpf(), bitupk(), strlen() +errchk read_dumpf, read_header, idsf_write_image, list_values + +begin + # Allocate space for program data structure + call smark (sp) + call salloc (ids, LEN_IDS, TY_STRUCT) + + # Read the next ids record + if (read_dumpf (fd, ids_buffer, LEN_IDS_RECORD) == EOF) { + # At end of IDSFILE + call printf ("DUMPF tape at EOF\n") + call sfree (sp) + return (EOF) + } + + scan = bitupk (ids_buffer, SCAN_OFFSET, NBITS_INT) + tape = bitupk (ids_buffer, TAPE_OFFSET, NBITS_INT) + this_record = (tape * 1000) + scan + if (is_in_range (records, this_record)) { + nids_read = nids_read + 1 + RECORD_NUMBER(ids) = this_record + iferr { + call calloc (COEFF(ids), MAX_COEFF, TY_DOUBLE) + call idsf_read_header (ids_buffer, ids) + } then { + call erract (EA_WARN) + call mfree (COEFF(ids), TY_DOUBLE) + call sfree (sp) + return (ERR) + } + + call print_header (ids, LONG_HEADER(cp)) + + if (MAKE_IMAGE(cp) == YES) { + call strcpy (IRAF_FILE(cp), out_fname, SZ_FNAME) + call sprintf (out_fname[strlen(out_fname) + 1], SZ_FNAME, ".%d") + call pargi (RECORD_NUMBER(ids)) + iferr { + + call idsf_write_image (ids_buffer, DATA_TYPE(cp), + PRINT_PIXELS(cp), out_fname, ids) + + } then { + call ERRACT (EA_WARN) + call mfree (COEFF(ids), TY_DOUBLE) + call sfree (sp) + return (ERR) + } + } + + if (PRINT_PIXELS(cp) == YES && MAKE_IMAGE(cp) == NO) { + call unpk_30 (ids_buffer, 1, pixels, NPIX_IDS_RECORD) + call list_values (pixels) + } + call mfree (COEFF(ids), TY_DOUBLE) + } + + call sfree (sp) + return (OK) +end + + +# SORT_INDEX -- Sort index information that precedes each IDSFILE. This +# index occupies 5 PRU's and points to ranges of records. Each index +# entry contains a PRU number and the low and high record numbers of the +# records that begin at the stated PRU. These three pieces of information +# are stored in a single 60-bit Cyber word. The number of records requested +# by the user that are actually in the IDSFILE is also counted. This +# number is returned as a parameter to the calling procedure. + +procedure sort_index (index_buf, records, sorted_index, nranges, nrecs_on_tape) + +int index_buf[ARB] # Buffer containing IDS index information +int records[3, MAX_RANGES] # Array of ranges of records specified by user +long sorted_index[LEN_INDEX] # Returned array of sorted index information +int nranges # Number of ranges of IDS records in IDSFILE +int nrecs_on_tape # Number of requested records actually on tape + +int i, start_pru, low_record_number, high_record_number, nrecs, j +long index[LEN_INDEX] +bool is_in_range() +int bitupk() +errchk asrtl, bitupk + +begin + nrecs_on_tape = 0 + nranges = 0 + do i = 1, NINT_CYBER_WRD * LEN_USER_INDEX, NINT_CYBER_WRD { + start_pru = bitupk (index_buf[i], NPRU_OFFSET, NBITS_NPRU) + if (start_pru == 0) + next + low_record_number = bitupk (index_buf[i], LRN_OFFSET, NBITS_LRN) + high_record_number = bitupk (index_buf[i], HRN_OFFSET, NBITS_HRN) + nrecs = high_record_number - low_record_number + 1 + nranges = nranges + 1 + index[nranges] = real (start_pru * 1000) + nrecs + + for (j=low_record_number; j<=high_record_number; j=j+1) { + if (is_in_range (records, j)) + nrecs_on_tape = nrecs_on_tape + 1 + } + } + + call asrtl (index, sorted_index, nranges) +end + + +# LIST_VALUES -- Print the ids pixel values. + +procedure list_values (pixel_buf) + +real pixel_buf[NPIX_IDS_RECORD] # Buffer containing pixels to be listed +int n_pix + +begin + for (n_pix = 1; n_pix <= NPIX_IDS_RECORD; n_pix = n_pix + 4) { + call printf ("%10.4e %10.4e %10.4e %10.4e\n") + call pargr (pixel_buf[n_pix]) + call pargr (pixel_buf[n_pix + 1]) + call pargr (pixel_buf[n_pix + 2]) + call pargr (pixel_buf[n_pix + 3]) + } + call printf ("\n") +end + +# IDSF_READ_HEADER -- Decode ids header parameters from the input buffer and +# fill the program data structure. + +procedure idsf_read_header (ids_buffer, ids) + +int ids_buffer[NINT_CYBER_WRD*LEN_IDS_RECORD] # Input IDSFILE buffer +pointer ids # Pointer to program data structure + +int n_coeff, i +char alpha[3] +int bitupk() +double convert_60bit_fp() +errchk bitupk, unpk_60i, convert_60bit_fp, unpk_id, display_code + +begin + # Get unsigned integer parameters from header + ITM(ids) = bitupk (ids_buffer, ITM_OFFSET, NBITS_INT) + NP1(ids) = bitupk (ids_buffer, NP1_OFFSET, NBITS_INT) + NP2(ids) = bitupk (ids_buffer, NP2_OFFSET, NBITS_INT) + BEAM_NUMBER(ids) = bitupk (ids_buffer, BEAM_OFFSET, NBITS_INT) + SMODE(ids) = bitupk (ids_buffer, SMODE_OFFSET, NBITS_INT) + if (SMODE(ids) != 0) { + # Determine companion record number + if (BEAM_NUMBER(ids) == 1) + COMPANION_RECORD(ids) = RECORD_NUMBER(ids) - 1 + else + COMPANION_RECORD(ids) = RECORD_NUMBER(ids) + 1 + } + UT(ids) = bitupk (ids_buffer, UT_OFFSET, NBITS_INT) + ST(ids) = bitupk (ids_buffer, ST_OFFSET, NBITS_INT) + + # The following integer parameters can be negative + call unpk_60i (ids_buffer, DF_OFFSET, DF_FLAG(ids), 1) + call unpk_60i (ids_buffer, SM_OFFSET, SM_FLAG(ids), 1) + call unpk_60i (ids_buffer, QF_OFFSET, QF_FLAG(ids), 1) + call unpk_60i (ids_buffer, DC_OFFSET, DC_FLAG(ids), 1) + call unpk_60i (ids_buffer, QD_OFFSET, QD_FLAG(ids), 1) + call unpk_60i (ids_buffer, EX_OFFSET, EX_FLAG(ids), 1) + call unpk_60i (ids_buffer, BS_OFFSET, BS_FLAG(ids), 1) + call unpk_60i (ids_buffer, CA_OFFSET, CA_FLAG(ids), 1) + call unpk_60i (ids_buffer, CO_OFFSET, CO_FLAG(ids), 1) + call unpk_60i (ids_buffer, OFLAG_OFFSET, OFLAG(ids), 1) + + # If the dispersion flag (DF) is set, get the coeffecients. The pointer + # to the coeffecient array is stored in the structure ids. + if (DF_FLAG(ids) > -1) { + n_coeff = DF_FLAG(ids) + do i = 1, n_coeff { + Memd[COEFF(ids)+i-1] = convert_60bit_fp (ids_buffer, + (COEFF_OFFSET + (i - 1)) * 64 + 1) + } + } + + # These header values converted from Cyber 60-bit floating point + HA(ids) = convert_60bit_fp (ids_buffer, HA_OFFSET) + AIRMASS(ids) = convert_60bit_fp (ids_buffer, AIR_OFFSET) + RA(ids) = convert_60bit_fp (ids_buffer, RA_OFFSET) + DEC(ids) = convert_60bit_fp (ids_buffer, DEC_OFFSET) + LAMBDA0(ids) = convert_60bit_fp (ids_buffer, LAM_OFFSET) + DELTA_LAMBDA(ids) = convert_60bit_fp (ids_buffer, DEL_OFFSET) + + # The 3 character ALPHA_ID is stored in Cyber display code + call display_code (bitupk (ids_buffer, ALPHA1_OFFSET, NBITS_DC), + alpha[1]) + call display_code (bitupk (ids_buffer, ALPHA2_OFFSET, NBITS_DC), + alpha[2]) + call display_code (bitupk (ids_buffer, ALPHA3_OFFSET, NBITS_DC), + alpha[3]) + call strcpy (alpha, ALPHA_ID(ids), NCHAR_ALPHA) + + # The ids label is written in 7-bit ascii + call unpk_id (ids_buffer, IDS_ID_OFFSET, LABEL(ids)) +end + + +# PRINT_HEADER -- print the ids header in either long or short mode. + +procedure print_header (ids, long_header) + +pointer ids # Pointer to program data structure +int long_header # Print header in long format (YES/NO)? +int i + +real value1, value2 + +begin + if (long_header == YES) { + call printf ("RECORD = %d, label = \"%s\",\n") + call pargi (RECORD_NUMBER(ids)) + call pargstr (LABEL(ids)) + + if (OFLAG(ids) == 1) { + call printf ("oflag = OBJECT, beam_number = %d, ") + call pargi (BEAM_NUMBER(ids)) + } else if (OFLAG (ids) == 0) { + call printf ("oflag = SKY, beam_number = %d, ") + call pargi (BEAM_NUMBER(ids)) + } + call printf ("alpha_ID = %s") + call pargstr (ALPHA_ID(ids)) + if (SMODE(ids) != 0) { + call printf (", companion = %d,\n") + call pargi (COMPANION_RECORD(ids)) + } else + call printf (",\n") + + call printf ("airmass = %5.3f,%24tW0 = %0.3f,") + call pargd (AIRMASS(ids)) + call pargd (LAMBDA0(ids)) + call printf (" WPC = %0.3f, ITM = %d,\n") + call pargd (DELTA_LAMBDA(ids)) + call pargi (ITM(ids)) + call printf ("NP1 = %d, NP2 = %d,") + call pargi (NP1(ids)) + call pargi (NP2(ids)) + + if (IS_INDEFI (UT(ids))) + value1 = INDEFR + else + value1 = real (UT(ids) / 3600.) + + if (IS_INDEFI (ST(ids))) + value2 = INDEFR + else + value2 = real (ST(ids) / 3600.) + call printf (" UT = %h, ST = %h,\n") + call pargr (value1) + call pargr (value2) + + call printf ("HA = %h,") + call pargd (HA(ids)) + call printf (" RA = %h, DEC = %h,\n") + call pargd (RA(ids)) + call pargd (DEC(ids)) + call printf ("df = %d, sm = %d, qf = %d, dc = %d, qd = %d, ") + call pargi (DF_FLAG(ids)) + call pargi (SM_FLAG(ids)) + call pargi (QF_FLAG(ids)) + call pargi (DC_FLAG(ids)) + call pargi (QD_FLAG(ids)) + call printf ("ex = %d, bs = %d, ca = %d, co = %d") + call pargi (EX_FLAG(ids)) + call pargi (BS_FLAG(ids)) + call pargi (CA_FLAG(ids)) + call pargi (CO_FLAG(ids)) + + # The df coeffecients are printed out in the case where the df + # flag is set, and the first coefficient is nonzero. The later + # condition is a test for IDSOUT data, where the df coeffecients + # have been applied but not stored in the header. + + if (DF_FLAG(ids) != -1 && COEFF(ids) != 0) { + call printf (",\n") + do i = 1, DF_FLAG(ids) { + call printf ("df[%d] = %10.8g") + call pargi(i) + call pargd(Memd[COEFF(ids)+i-1]) + if (i != DF_FLAG(ids)) + call printf (", ") + if (mod (i, 4) == 0) + call printf ("\n") + } + call printf ("\n") + } else + call printf ("\n") + call printf ("\n") + } else { + call printf ("RECORD = %d, label = \"%s\"\n") + call pargi (RECORD_NUMBER(ids)) + call pargstr (LABEL(ids)) + } +end + + +# IDSF_WRITE_IMAGE -- pixels are unpacked from the input buffer and written to +# a one dimensional IRAF image. + +procedure idsf_write_image (ids_buffer, data_type, print_pixels, out_fname, + ids) + +int ids_buffer[NINT_CYBER_WRD * LEN_IDS_RECORD] # Input IDSFILE buffer +int data_type # Data type of pixels to be written +int print_pixels # List pixel values (YES/NO)? +char out_fname[SZ_FNAME] # Name of output image +pointer ids # Pointer to program data structure + +pointer im, pixels +pointer impl1r(), immap() +errchk immap, unpk_30, cy_store_keywords, imunmap + +begin + # Map new iraf image and set up image header + im = immap (out_fname, NEW_IMAGE, LEN_USER_AREA) + IM_NDIM(im) = 1 + IM_LEN(im, 1) = NPIX_IDS_RECORD + call strcpy (LABEL(ids), IM_TITLE(im), SZ_IMTITLE) + IM_PIXTYPE(im) = data_type + pixels = impl1r(im) + + # Convert pixels to spp reals and write image line + call unpk_30 (ids_buffer, 1, Memr[pixels], NPIX_IDS_RECORD) + + if (print_pixels == YES) + call list_values (Memr[pixels]) + + # Write ids specific header words to iraf image header + call cy_store_keywords (ids, im) + + call imunmap (im) +end diff --git a/noao/mtlocal/cyber/t_ridsout.x b/noao/mtlocal/cyber/t_ridsout.x new file mode 100644 index 00000000..473d120f --- /dev/null +++ b/noao/mtlocal/cyber/t_ridsout.x @@ -0,0 +1,386 @@ +include +include +include +include +include +include "cyber.h" + +define MAX_WIDTH 10 # Maximum format width for pixel data + +# T_RIDSOUT -- the IDSOUT format reader, which reads a text file of IDS +# records in IDSOUT format. Each IDS record contains 133 card images: 4 +# cards of header information followed by 128 cards of pixel values and +# one blank card. The IDSOUT reader will read the text file and optionally +# convert the data to a sequence of one dimensional IRAF images. The +# header information can be printed in long or short form; the pixel values +# can also be listed. + +procedure t_ridsout () + +pointer sp, cp +int fd +char in_fname[SZ_FNAME] +int get_data_type(), clpopni(), clgfil(), btoi() +bool clgetb() +char clgetc() +errchk clpopni + +begin + # Allocate space for the control parameter descriptor structure + call smark (sp) + call salloc (cp, LEN_CP, TY_STRUCT) + + call fseti (STDOUT, F_FLUSHNL, YES) + + # Get parameters from the cl and fill control parameter structure + fd = clpopni ("idsout_file") + + LONG_HEADER(cp) = btoi (clgetb ("long_header")) + PRINT_PIXELS(cp) = btoi (clgetb ("print_pixels")) + call clgstr ("record_numbers", REC_NUMBERS(cp), SZ_LINE) + + # If an output image is to be written, get root output file name + # and output data type. + + MAKE_IMAGE(cp) = btoi (clgetb ("make_image")) + if (MAKE_IMAGE(cp) == YES) { + call clgstr ("iraf_file", IRAF_FILE(cp), SZ_FNAME) + DATA_TYPE(cp) = get_data_type (clgetc ("data_type")) + if (DATA_TYPE(cp) == ERR) + DATA_TYPE(cp) = TY_REAL + } + + while (clgfil (fd, in_fname, SZ_FNAME) != EOF) { + call read_idsout (in_fname, cp) + } + + call clpcls (fd) + call sfree (sp) +end + + +# READ_IDSOUT -- open IDSOUT text file and direct processing depending on +# user's request. + +procedure read_idsout (in_fname, cp) + +char in_fname[SZ_FNAME] # Input file name +pointer cp # Pointer to control parameter structure + +pointer sp, ids +char out_fname[SZ_FNAME] +int in, records[3, MAX_RANGES], nrecs, nrecs_read, stat +bool is_in_range() +int open(), decode_ranges(), idso_read_header(), strlen() +errchk open, read, idso_read_header, decode_ranges, idso_write_image +errchk copy_pixels + +begin + # Open text file and allocate space for descriptor structure ids + in = open (in_fname, READ_ONLY, TEXT_FILE) + call smark (sp) + call salloc (ids, LEN_IDS, TY_STRUCT) + + if (decode_ranges (REC_NUMBERS(cp), records, MAX_RANGES, nrecs) == ERR) + call error (1, "Error in record_numbers specification") + + nrecs_read = 0 + repeat { + iferr { + stat = idso_read_header (in, ids) + } then { + call erract (EA_WARN) + call eprintf ("Bad header, attempting to skip pixels\n") + call copy_pixels (in, NULL) + next + } + + if (stat == EOF) { + call printf ("\nIDSOUT file \"%s\" at EOF\n") + call pargstr (in_fname) + call close (in) + call sfree (sp) + return + } + + if (is_in_range (records, RECORD_NUMBER(ids))) { + nrecs_read = nrecs_read + 1 + call print_header (ids, LONG_HEADER(cp)) + + if (MAKE_IMAGE(cp) == YES) { + call strcpy (IRAF_FILE(cp), out_fname, SZ_FNAME) + call sprintf (out_fname[strlen(out_fname)+1], SZ_FNAME, + ".%04d") + call pargi (RECORD_NUMBER(ids)) + iferr { + call idso_write_image (in, DATA_TYPE(cp), + PRINT_PIXELS(cp), out_fname, ids) + } then { + call erract (EA_WARN) + if (PRINT_PIXELS(cp) == YES) + call copy_pixels (in, STDOUT) + else + call copy_pixels (in, NULL) + next + } + + } else if (PRINT_PIXELS(cp) == YES) + # Simply copy card image data to standard output + call copy_pixels (in, STDOUT) + + if (PRINT_PIXELS(cp) == NO && MAKE_IMAGE(cp) == NO) + call copy_pixels (in, NULL) + } else + call copy_pixels (in, NULL) + + } until (nrecs_read == nrecs) + + call sfree (sp) + call close (in) +end + + +# IDSO_READ_HEADER -- decode header parameter from IDSOUT header cards. The +# IDSO prefix implies IDSOUT format; IDSF is used for IDSFILE format. The +# first four cards of each record in the IDSOUT file are the header. + +int procedure idso_read_header (in, ids) + +int in # File descriptor of input text file +pointer ids # Pointer to program data structure + +pointer sp, temp +int record, tape, junk, ip +int fscan(), nscan(), strlen(), getline() +errchk fscan + +begin + # Allocate space on stack for temporary char storage + call smark (sp) + call salloc (temp, SZ_LINE, TY_CHAR) + + # Skip any blank lines that may preceede the first header card. + # Break out of repeat when first non-blank line is encountered. + + repeat { + if (getline (in, Memc[temp]) == EOF) { + call sfree (sp) + return (EOF) + } else { + for (ip = temp; IS_WHITE (Memc[ip]); ip = ip + 1) + ; + if (Memc[ip] == '\n' || Memc[ip] == EOS) + # Blank line + ; + else + break + } + } + + call sscan (Memc[temp]) { + # Values to be read from first card of file: + call gargi (record) + call gargi (ITM(ids)) + call gargd (LAMBDA0(ids)) + call gargd (DELTA_LAMBDA(ids)) + call gargi (NP1(ids)) + call gargi (NP2(ids)) + call gargi (BEAM_NUMBER(ids)) + call gargi (junk) + call gargi (junk) + call gargi (SMODE(ids)) + call gargi (UT(ids)) + if (nscan() < 11) + call error (2, "First header card incomplete") + } + + if (fscan (in) == EOF) { + call sfree (sp) + return (EOF) + } else { + # Values to be read from second card of file + call gargi (ST(ids)) + call gargd (RA(ids)) + call gargd (DEC(ids)) + call gargi (tape) + if (tape > 0) + RECORD_NUMBER(ids) = tape * 1000 + record + else + RECORD_NUMBER(ids) = record + call gargi (DF_FLAG(ids)) + call gargi (SM_FLAG(ids)) + call gargi (QF_FLAG(ids)) + call gargi (DC_FLAG(ids)) + call gargi (QD_FLAG(ids)) + call gargi (EX_FLAG(ids)) + call gargi (BS_FLAG(ids)) + if (nscan() < 11) + call error (3, "Second header card incomplete") + } + + if (fscan (in) == EOF) { + call sfree (sp) + return (EOF) + } else { + # Values to be read from third card of file + call gargi (CA_FLAG(ids)) + call gargi (CO_FLAG(ids)) + call gargwrd (ALPHA_ID(ids), 3) + call gargi (OFLAG(ids)) + call gargd (HA(ids)) + call gargd (AIRMASS(ids)) + if (nscan() < 6) { + call error (4, "Third header card incomplete") + } + } + + if (fscan (in) == EOF) { + call sfree (sp) + return (EOF) + } else { + call gargstr (Memc[temp], SZ_LINE) + + # Count number of characters in id string by skipping trailing + # white space. 'END' are always the last 3 characters of the + # line of text containing the ID; they are also skipped. + + for (ip = strlen(Memc[temp]) - 3; IS_WHITE (Memc[temp+ip-1]) && + ip > 0 ; ip = ip - 1) + ; + Memc[temp+ip] = EOS + call strcpy (Memc[temp], LABEL(ids), SZ_LINE) + } + + # Since this is IDSOUT format data, initialize COEFF(ids) to integer 0. + COEFF(ids) = 0 + + # Determine companion record number if appropriate + if (SMODE(ids) != 0) { + if (BEAM_NUMBER(ids) == 1) + COMPANION_RECORD(ids) = RECORD_NUMBER(ids) - 1 + else + COMPANION_RECORD(IDS) = RECORD_NUMBER(IDS) + 1 + } + + call sfree (sp) + return (OK) + +end + + +# COPY_PIXELS -- copy pixel values from input to output file. A blank line +# terminates each IDS record. + +procedure copy_pixels (in, out) + +int in # File descriptor of input text file +int out # File descriptor of output file + +pointer sp, line_buffer +bool leading_blank +int ip +int getline() + +begin + # Allocate space on stack for line_buffer + call smark (sp) + call salloc (line_buffer, SZ_LINE, TY_CHAR) + + # Ignore blank lines until first non_blank line encountered + leading_blank = true + + repeat { + if (getline (in, Memc[line_buffer]) == EOF) + call error (5, "Unexpected EOF when copying pixels") + + # Find first non-whitespace character + for (ip = line_buffer; IS_WHITE (Memc[ip]); ip = ip + 1) + ; + + if (Memc[ip] == '\n' || Memc[ip] == EOS) { + # Blank line + if (leading_blank) + ; + else + break + } else { + leading_blank = false + if (out != NULL) + call putline (out, Memc[line_buffer]) + } + } + + call sfree (sp) +end + + +# IDSO_WRITE_IMAGE -- convert card image data values to reals and write +# 1-d IRAF image. + +procedure idso_write_image (in, data_type, print_pixels, out_fname, ids) + +int in # File descriptor of input text file +int data_type # Data type of image pixels to be written +int print_pixels # Are pixel values to be listed (YES/NO)? +char out_fname[SZ_FNAME] # Name of output image +pointer ids # Pointer to program data structure + +pointer im, pixels, op +char line[SZ_LINE], temp[MAX_WIDTH] +int i, j, ip, iip, nchars +pointer immap(), impl1r() +int fscan(), nscan(), strlen(), ctowrd(), ctor() +errchk immap, imunmap, fscan, impl1r, unpk_30, cy_store_keywords + +begin + im = immap (out_fname, NEW_IMAGE, LEN_USER_AREA) + IM_NDIM(im) = 1 + IM_LEN(im, 1) = NPIX_IDS_RECORD + call strcpy (LABEL(ids), IM_TITLE(im), SZ_IMTITLE) + IM_PIXTYPE(im) = data_type + pixels = impl1r(im) + op = pixels + + # Position to first non_blank line. + repeat { + if (fscan (in) == EOF) { + call imunmap (im) + return + } else + call gargstr (line, SZ_LINE) + } until (strlen (line) > 1) + call reset_scan() + + + # Scan and convert pixel values to image pixels. A blank line + # terminates the IDS record. + + do i = 1, NLINES_PIXELS { + ip = 1 + # Extract 8 pixel values from each line of text + do j = 1, NPIX_LINE { + iip = 1 + call gargstr (line, SZ_LINE) + nchars = ctowrd (line, ip, temp, MAX_WIDTH) + nchars = ctor (temp, iip, Memr[op]) + # call gargr (Memr[op]) + op = op + 1 + } + + if (nscan () <= 1) + # Premature blank line encountered + break + + # Get next line of pixels + if (fscan (in) == EOF) + break + } + + if (print_pixels == YES) + call list_values (Memr[pixels]) + + # Write ids specific header words to iraf image header + call cy_store_keywords (ids, im) + + call imunmap (im) +end -- cgit