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/imred/quadred/src/ccdproc/generic/ccdred.h | 155 +++++ noao/imred/quadred/src/ccdproc/generic/cor.x | 695 ++++++++++++++++++++++ noao/imred/quadred/src/ccdproc/generic/corinput.x | 436 ++++++++++++++ noao/imred/quadred/src/ccdproc/generic/mkpkg | 12 + noao/imred/quadred/src/ccdproc/generic/proc.x | 678 +++++++++++++++++++++ 5 files changed, 1976 insertions(+) create mode 100644 noao/imred/quadred/src/ccdproc/generic/ccdred.h create mode 100644 noao/imred/quadred/src/ccdproc/generic/cor.x create mode 100644 noao/imred/quadred/src/ccdproc/generic/corinput.x create mode 100644 noao/imred/quadred/src/ccdproc/generic/mkpkg create mode 100644 noao/imred/quadred/src/ccdproc/generic/proc.x (limited to 'noao/imred/quadred/src/ccdproc/generic') diff --git a/noao/imred/quadred/src/ccdproc/generic/ccdred.h b/noao/imred/quadred/src/ccdproc/generic/ccdred.h new file mode 100644 index 00000000..ef41f592 --- /dev/null +++ b/noao/imred/quadred/src/ccdproc/generic/ccdred.h @@ -0,0 +1,155 @@ +# CCDRED Data Structures and Definitions + +# The CCD structure: This structure is used to communicate processing +# parameters between the package procedures. It contains pointers to +# data, calibration image IMIO pointers, scaling parameters, and the +# correction flags. The corrections flags indicate which processing +# operations are to be performed. The subsection parameters do not +# include a step size. A step size is assumed. If arbitrary subsampling +# is desired this would be the next generalization. + +define LEN_CCD 75 # Length of CCD structure + +# CCD data coordinates +define CCD_C1 Memi[$1] # CCD starting column +define CCD_C2 Memi[$1+1] # CCD ending column +define CCD_L1 Memi[$1+2] # CCD starting line +define CCD_L2 Memi[$1+3] # CCD ending line + +# Input data +define IN_IM Memi[$1+4] # Input image pointer +define IN_C1 Memi[$1+5] # Input data starting column +define IN_C2 Memi[$1+6] # Input data ending column +define IN_L1 Memi[$1+7] # Input data starting line +define IN_L2 Memi[$1+8] # Input data ending line +define IN_NSEC Memi[$1+71] # Number of input pieces +define IN_SEC Memi[$1+72] # Pointer to sections (c1,c2,l1,l2)xn + +# Output data +define OUT_IM Memi[$1+9] # Output image pointer +define OUT_C1 Memi[$1+10] # Output data starting column +define OUT_C2 Memi[$1+11] # Output data ending column +define OUT_L1 Memi[$1+12] # Output data starting line +define OUT_L2 Memi[$1+13] # Output data ending line +define OUT_SEC Memi[$1+73] # Pointer to sections (c1,c2,l1,l2)xn + +# Zero level data +define ZERO_IM Memi[$1+14] # Zero level image pointer +define ZERO_C1 Memi[$1+15] # Zero level data starting column +define ZERO_C2 Memi[$1+16] # Zero level data ending column +define ZERO_L1 Memi[$1+17] # Zero level data starting line +define ZERO_L2 Memi[$1+18] # Zero level data ending line + +# Dark count data +define DARK_IM Memi[$1+19] # Dark count image pointer +define DARK_C1 Memi[$1+20] # Dark count data starting column +define DARK_C2 Memi[$1+21] # Dark count data ending column +define DARK_L1 Memi[$1+22] # Dark count data starting line +define DARK_L2 Memi[$1+23] # Dark count data ending line + +# Flat field data +define FLAT_IM Memi[$1+24] # Flat field image pointer +define FLAT_C1 Memi[$1+25] # Flat field data starting column +define FLAT_C2 Memi[$1+26] # Flat field data ending column +define FLAT_L1 Memi[$1+27] # Flat field data starting line +define FLAT_L2 Memi[$1+28] # Flat field data ending line + +# Illumination data +define ILLUM_IM Memi[$1+29] # Illumination image pointer +define ILLUM_C1 Memi[$1+30] # Illumination data starting column +define ILLUM_C2 Memi[$1+31] # Illumination data ending column +define ILLUM_L1 Memi[$1+32] # Illumination data starting line +define ILLUM_L2 Memi[$1+33] # Illumination data ending line + +# Fringe data +define FRINGE_IM Memi[$1+34] # Fringe image pointer +define FRINGE_C1 Memi[$1+35] # Fringe data starting column +define FRINGE_C2 Memi[$1+36] # Fringe data ending column +define FRINGE_L1 Memi[$1+37] # Fringe data starting line +define FRINGE_L2 Memi[$1+38] # Fringe data ending line + +# Trim section +define TRIM_C1 Memi[$1+39] # Trim starting column +define TRIM_C2 Memi[$1+40] # Trim ending column +define TRIM_L1 Memi[$1+41] # Trim starting line +define TRIM_L2 Memi[$1+42] # Trim ending line + +# Bias section +define BIAS_C1 Memi[$1+43] # Bias starting column +define BIAS_C2 Memi[$1+44] # Bias ending column +define BIAS_L1 Memi[$1+45] # Bias starting line +define BIAS_L2 Memi[$1+46] # Bias ending line +define BIAS_SEC Memi[$1+74] # Multiple bias sections + +define READAXIS Memi[$1+47] # Read out axis (1=cols, 2=lines) +define CALCTYPE Memi[$1+48] # Calculation data type +define NBADCOLS Memi[$1+49] # Number of column interpolation regions +define BADCOLS Memi[$1+50] # Pointer to col interpolation regions +define NBADLINES Memi[$1+51] # Number of line interpolation regions +define BADLINES Memi[$1+52] # Pointer to line interpolation regions +define OVERSCAN_VEC Memi[$1+53] # Pointer to overscan vector +define DARKSCALE Memr[P2R($1+54)] # Dark count scale factor +define FRINGESCALE Memr[P2R($1+55)] # Fringe scale factor +define FLATSCALE Memr[P2R($1+56)] # Flat field scale factor +define ILLUMSCALE Memr[P2R($1+57)] # Illumination scale factor +define MINREPLACE Memr[P2R($1+58)] # Minimum replacement value +define MEAN Memr[P2R($1+59)] # Mean of output image +define COR Memi[$1+60] # Overall correction flag +define CORS Memi[$1+61+($2-1)] # Individual correction flags + +# Individual components of input, output, and bias section pieces. +define IN_SC1 Memi[IN_SEC($1)+4*$2-4] +define IN_SC2 Memi[IN_SEC($1)+4*$2-3] +define IN_SL1 Memi[IN_SEC($1)+4*$2-2] +define IN_SL2 Memi[IN_SEC($1)+4*$2-1] +define OUT_SC1 Memi[OUT_SEC($1)+4*$2-4] +define OUT_SC2 Memi[OUT_SEC($1)+4*$2-3] +define OUT_SL1 Memi[OUT_SEC($1)+4*$2-2] +define OUT_SL2 Memi[OUT_SEC($1)+4*$2-1] +define BIAS_SC1 Memi[BIAS_SEC($1)+4*$2-4] +define BIAS_SC2 Memi[BIAS_SEC($1)+4*$2-3] +define BIAS_SL1 Memi[BIAS_SEC($1)+4*$2-2] +define BIAS_SL2 Memi[BIAS_SEC($1)+4*$2-1] + +# The correction array contains the following elements with array indices +# given by the macro definitions. + +define NCORS 10 # Number of corrections + +define FIXPIX 1 # Fix bad pixels +define TRIM 2 # Trim image +define OVERSCAN 3 # Apply overscan correction +define ZEROCOR 4 # Apply zero level correction +define DARKCOR 5 # Apply dark count correction +define FLATCOR 6 # Apply flat field correction +define ILLUMCOR 7 # Apply illumination correction +define FRINGECOR 8 # Apply fringe correction +define FINDMEAN 9 # Find the mean of the output image +define MINREP 10 # Check and replace minimum value + +# The following definitions identify the correction values in the correction +# array. They are defined in terms of bit fields so that it is possible to +# add corrections to form unique combination corrections. Some of +# these combinations are implemented as compound operations for efficiency. + +define O 001B # overscan +define Z 002B # zero level +define D 004B # dark count +define F 010B # flat field +define I 020B # Illumination +define Q 040B # Fringe + +# The following correction combinations are recognized. + +define ZO 003B # zero level + overscan +define DO 005B # dark count + overscan +define DZ 006B # dark count + zero level +define DZO 007B # dark count + zero level + overscan +define FO 011B # flat field + overscan +define FZ 012B # flat field + zero level +define FZO 013B # flat field + zero level + overscan +define FD 014B # flat field + dark count +define FDO 015B # flat field + dark count + overscan +define FDZ 016B # flat field + dark count + zero level +define FDZO 017B # flat field + dark count + zero level + overscan +define QI 060B # fringe + illumination diff --git a/noao/imred/quadred/src/ccdproc/generic/cor.x b/noao/imred/quadred/src/ccdproc/generic/cor.x new file mode 100644 index 00000000..0dc21310 --- /dev/null +++ b/noao/imred/quadred/src/ccdproc/generic/cor.x @@ -0,0 +1,695 @@ +include "ccdred.h" + + +.help cor Feb87 noao.imred.ccdred +.nf ---------------------------------------------------------------------------- +cor -- Process CCD image lines + +These procedures are the heart of the CCD processing. They do the desired +set of processing operations on the image line data as efficiently as +possible. They are called by the PROC procedures. There are four procedures +one for each readout axis and one for short and real image data. +Some sets of operations are coded as single compound operations for efficiency. +To keep the number of combinations managable only the most common +combinations are coded as compound operations. The combinations +consist of any set of line overscan, column overscan, zero level, dark +count, and flat field and any set of illumination and fringe +correction. The corrections are applied in place to the output vector. + +The column readout procedure is more complicated in order to handle +zero level and flat field corrections specified as one dimensional +readout corrections instead of two dimensional calibration images. +Column readout format is probably extremely rare and the 1D readout +corrections are used only for special types of data. +.ih +SEE ALSO +proc, ccdred.h +.endhelp ----------------------------------------------------------------------- + + +# COR1 -- Correct image lines with readout axis 1 (lines). + +procedure cor1s (cors, out, overscan, zero, dark, flat, illum, + fringe, n, darkscale, flatscale, illumscale, frgscale) + +int cors[ARB] # Correction flags +short out[n] # Output data +real overscan # Overscan value +short zero[n] # Zero level correction +short dark[n] # Dark count correction +short flat[n] # Flat field correction +short illum[n] # Illumination correction +short fringe[n] # Fringe correction +int n # Number of pixels +real darkscale # Dark count scale factor +real flatscale # Flat field scale factor +real illumscale # Illumination scale factor +real frgscale # Fringe scale factor + +int i, op + +begin + op = cors[OVERSCAN] + cors[ZEROCOR] + cors[DARKCOR] + cors[FLATCOR] + switch (op) { + case O: # overscan + do i = 1, n + out[i] = out[i] - overscan + case Z: # zero level + do i = 1, n + out[i] = out[i] - zero[i] + + case ZO: # zero level + overscan + do i = 1, n + out[i] = out[i] - overscan - zero[i] + + case D: # dark count + do i = 1, n + out[i] = out[i] - darkscale * dark[i] + case DO: # dark count + overscan + do i = 1, n + out[i] = out[i] - overscan - darkscale * dark[i] + case DZ: # dark count + zero level + do i = 1, n + out[i] = out[i] - zero[i] - darkscale * dark[i] + case DZO: # dark count + zero level + overscan + do i = 1, n + out[i] = out[i] - overscan - zero[i] - darkscale * dark[i] + + case F: # flat field + do i = 1, n + out[i] = out[i] * flatscale / flat[i] + case FO: # flat field + overscan + do i = 1, n + out[i] = (out[i] - overscan) * flatscale / flat[i] + case FZ: # flat field + zero level + do i = 1, n + out[i] = (out[i] - zero[i]) * flatscale / flat[i] + case FZO: # flat field + zero level + overscan + do i = 1, n + out[i] = (out[i] - overscan - zero[i]) * flatscale / + flat[i] + case FD: # flat field + dark count + do i = 1, n + out[i] = (out[i] - darkscale * dark[i]) * flatscale / flat[i] + case FDO: # flat field + dark count + overscan + do i = 1, n + out[i] = (out[i] - overscan - darkscale * dark[i]) * + flatscale / flat[i] + case FDZ: # flat field + dark count + zero level + do i = 1, n + out[i] = (out[i] - zero[i] - darkscale * dark[i]) * + flatscale / flat[i] + case FDZO: # flat field + dark count + zero level + overscan + do i = 1, n + out[i] = (out[i] - overscan - zero[i] - + darkscale * dark[i]) * flatscale / flat[i] + } + + # Often these operations will not be performed so test for no + # correction rather than go through the switch. + + op = cors[ILLUMCOR] + cors[FRINGECOR] + if (op != 0) { + switch (op) { + case I: # illumination + do i = 1, n + out[i] = out[i] * illumscale / illum[i] + case Q: # fringe + do i = 1, n + out[i] = out[i] - frgscale * fringe[i] + case QI: # fringe + illumination + do i = 1, n + out[i] = out[i]*illumscale/illum[i] - frgscale*fringe[i] + } + } +end + + +# COR2 -- Correct lines for readout axis 2 (columns). This procedure is +# more complex than when the readout is along the image lines because the +# zero level and/or flat field corrections may be single readout column +# vectors. + +procedure cor2s (line, cors, out, overscan, zero, dark, flat, illum, + fringe, n, zeroim, flatim, darkscale, flatscale, illumscale, frgscale) + +int line # Line to be corrected +int cors[ARB] # Correction flags +short out[n] # Output data +real overscan[n] # Overscan value +short zero[n] # Zero level correction +short dark[n] # Dark count correction +short flat[n] # Flat field correction +short illum[n] # Illumination correction +short fringe[n] # Fringe correction +int n # Number of pixels +pointer zeroim # Zero level IMIO pointer (NULL if 1D vector) +pointer flatim # Flat field IMIO pointer (NULL if 1D vector) +real darkscale # Dark count scale factor +real flatscale # Flat field scale factor +real illumscale # Illumination scale factor +real frgscale # Fringe scale factor + +short zeroval +real flatval +int i, op + +begin + op = cors[OVERSCAN] + cors[ZEROCOR] + cors[DARKCOR] + cors[FLATCOR] + switch (op) { + case O: # overscan + do i = 1, n + out[i] = out[i] - overscan[i] + case Z: # zero level + if (zeroim != NULL) + do i = 1, n + out[i] = out[i] - zero[i] + else { + zeroval = zero[line] + do i = 1, n + out[i] = out[i] - zeroval + } + + case ZO: # zero level + overscan + if (zeroim != NULL) + do i = 1, n + out[i] = out[i] - overscan[i] - zero[i] + else { + zeroval = zero[line] + do i = 1, n + out[i] = out[i] - overscan[i] - zeroval + } + + case D: # dark count + do i = 1, n + out[i] = out[i] - darkscale * dark[i] + case DO: # dark count + overscan + do i = 1, n + out[i] = out[i] - overscan[i] - darkscale * dark[i] + case DZ: # dark count + zero level + if (zeroim != NULL) + do i = 1, n + out[i] = out[i] - zero[i] - darkscale * dark[i] + else { + zeroval = zero[line] + do i = 1, n + out[i] = out[i] - zeroval - darkscale * dark[i] + } + case DZO: # dark count + zero level + overscan + if (zeroim != NULL) + do i = 1, n + out[i] = out[i] - overscan[i] - zero[i] - + darkscale * dark[i] + else { + zeroval = zero[line] + do i = 1, n + out[i] = out[i] - overscan[i] - zeroval - + darkscale * dark[i] + } + + case F: # flat field + if (flatim != NULL) { + do i = 1, n + out[i] = out[i] * flatscale / flat[i] + } else { + flatval = flatscale / flat[line] + do i = 1, n + out[i] = out[i] * flatval + } + case FO: # flat field + overscan + if (flatim != NULL) { + do i = 1, n + out[i] = (out[i] - overscan[i]) * flatscale / flat[i] + } else { + flatval = flatscale / flat[line] + do i = 1, n + out[i] = (out[i] - overscan[i]) * flatval + } + case FZ: # flat field + zero level + if (flatim != NULL) { + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - zero[i]) * flatscale / flat[i] + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - zeroval) * flatscale / flat[i] + } + } else { + flatval = flatscale / flat[line] + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - zero[i]) * flatval + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - zeroval) * flatval + } + } + case FZO: # flat field + zero level + overscan + if (flatim != NULL) { + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - overscan[i] - zero[i]) * + flatscale / flat[i] + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - overscan[i] - zeroval) * + flatscale / flat[i] + } + } else { + flatval = flatscale / flat[line] + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - overscan[i] - zero[i]) * flatval + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - overscan[i] - zeroval) * flatval + } + } + case FD: # flat field + dark count + if (flatim != NULL) { + do i = 1, n + out[i] = (out[i] - darkscale * dark[i]) * flatscale/flat[i] + } else { + flatval = flatscale / flat[line] + do i = 1, n + out[i] = (out[i] - darkscale * dark[i]) * flatval + } + case FDO: # flat field + dark count + overscan + if (flatim != NULL) { + do i = 1, n + out[i] = (out[i] - overscan[i] - darkscale * dark[i]) * + flatscale / flat[i] + } else { + flatval = flatscale / flat[line] + do i = 1, n + out[i] = (out[i] - overscan[i] - darkscale * dark[i]) * + flatval + } + case FDZ: # flat field + dark count + zero level + if (flatim != NULL) { + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - zero[i] - darkscale * dark[i]) * + flatscale / flat[i] + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - zeroval - darkscale * dark[i]) * + flatscale / flat[i] + } + } else { + flatval = flatscale / flat[line] + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - zero[i] - darkscale * dark[i]) * + flatval + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - zeroval - darkscale * dark[i]) * + flatval + } + } + case FDZO: # flat field + dark count + zero level + overscan + if (flatim != NULL) { + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - overscan[i] - zero[i] - + darkscale * dark[i]) * flatscale / flat[i] + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - overscan[i] - zeroval - + darkscale * dark[i]) * flatscale / flat[i] + } + } else { + flatval = flatscale / flat[line] + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - overscan[i] - zero[i] - + darkscale * dark[i]) * flatval + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - overscan[i] - zeroval - + darkscale * dark[i]) * flatval + } + } + } + + # Often these operations will not be performed so test for no + # correction rather than go through the switch. + + op = cors[ILLUMCOR] + cors[FRINGECOR] + if (op != 0) { + switch (op) { + case I: # illumination + do i = 1, n + out[i] = out[i] * illumscale / illum[i] + case Q: # fringe + do i = 1, n + out[i] = out[i] - frgscale * fringe[i] + case QI: # fringe + illumination + do i = 1, n + out[i] = out[i]*illumscale/illum[i] - frgscale*fringe[i] + } + } +end + +# COR1 -- Correct image lines with readout axis 1 (lines). + +procedure cor1r (cors, out, overscan, zero, dark, flat, illum, + fringe, n, darkscale, flatscale, illumscale, frgscale) + +int cors[ARB] # Correction flags +real out[n] # Output data +real overscan # Overscan value +real zero[n] # Zero level correction +real dark[n] # Dark count correction +real flat[n] # Flat field correction +real illum[n] # Illumination correction +real fringe[n] # Fringe correction +int n # Number of pixels +real darkscale # Dark count scale factor +real flatscale # Flat field scale factor +real illumscale # Illumination scale factor +real frgscale # Fringe scale factor + +int i, op + +begin + op = cors[OVERSCAN] + cors[ZEROCOR] + cors[DARKCOR] + cors[FLATCOR] + switch (op) { + case O: # overscan + do i = 1, n + out[i] = out[i] - overscan + case Z: # zero level + do i = 1, n + out[i] = out[i] - zero[i] + + case ZO: # zero level + overscan + do i = 1, n + out[i] = out[i] - overscan - zero[i] + + case D: # dark count + do i = 1, n + out[i] = out[i] - darkscale * dark[i] + case DO: # dark count + overscan + do i = 1, n + out[i] = out[i] - overscan - darkscale * dark[i] + case DZ: # dark count + zero level + do i = 1, n + out[i] = out[i] - zero[i] - darkscale * dark[i] + case DZO: # dark count + zero level + overscan + do i = 1, n + out[i] = out[i] - overscan - zero[i] - darkscale * dark[i] + + case F: # flat field + do i = 1, n + out[i] = out[i] * flatscale / flat[i] + case FO: # flat field + overscan + do i = 1, n + out[i] = (out[i] - overscan) * flatscale / flat[i] + case FZ: # flat field + zero level + do i = 1, n + out[i] = (out[i] - zero[i]) * flatscale / flat[i] + case FZO: # flat field + zero level + overscan + do i = 1, n + out[i] = (out[i] - overscan - zero[i]) * flatscale / + flat[i] + case FD: # flat field + dark count + do i = 1, n + out[i] = (out[i] - darkscale * dark[i]) * flatscale / flat[i] + case FDO: # flat field + dark count + overscan + do i = 1, n + out[i] = (out[i] - overscan - darkscale * dark[i]) * + flatscale / flat[i] + case FDZ: # flat field + dark count + zero level + do i = 1, n + out[i] = (out[i] - zero[i] - darkscale * dark[i]) * + flatscale / flat[i] + case FDZO: # flat field + dark count + zero level + overscan + do i = 1, n + out[i] = (out[i] - overscan - zero[i] - + darkscale * dark[i]) * flatscale / flat[i] + } + + # Often these operations will not be performed so test for no + # correction rather than go through the switch. + + op = cors[ILLUMCOR] + cors[FRINGECOR] + if (op != 0) { + switch (op) { + case I: # illumination + do i = 1, n + out[i] = out[i] * illumscale / illum[i] + case Q: # fringe + do i = 1, n + out[i] = out[i] - frgscale * fringe[i] + case QI: # fringe + illumination + do i = 1, n + out[i] = out[i]*illumscale/illum[i] - frgscale*fringe[i] + } + } +end + + +# COR2 -- Correct lines for readout axis 2 (columns). This procedure is +# more complex than when the readout is along the image lines because the +# zero level and/or flat field corrections may be single readout column +# vectors. + +procedure cor2r (line, cors, out, overscan, zero, dark, flat, illum, + fringe, n, zeroim, flatim, darkscale, flatscale, illumscale, frgscale) + +int line # Line to be corrected +int cors[ARB] # Correction flags +real out[n] # Output data +real overscan[n] # Overscan value +real zero[n] # Zero level correction +real dark[n] # Dark count correction +real flat[n] # Flat field correction +real illum[n] # Illumination correction +real fringe[n] # Fringe correction +int n # Number of pixels +pointer zeroim # Zero level IMIO pointer (NULL if 1D vector) +pointer flatim # Flat field IMIO pointer (NULL if 1D vector) +real darkscale # Dark count scale factor +real flatscale # Flat field scale factor +real illumscale # Illumination scale factor +real frgscale # Fringe scale factor + +real zeroval +real flatval +int i, op + +begin + op = cors[OVERSCAN] + cors[ZEROCOR] + cors[DARKCOR] + cors[FLATCOR] + switch (op) { + case O: # overscan + do i = 1, n + out[i] = out[i] - overscan[i] + case Z: # zero level + if (zeroim != NULL) + do i = 1, n + out[i] = out[i] - zero[i] + else { + zeroval = zero[line] + do i = 1, n + out[i] = out[i] - zeroval + } + + case ZO: # zero level + overscan + if (zeroim != NULL) + do i = 1, n + out[i] = out[i] - overscan[i] - zero[i] + else { + zeroval = zero[line] + do i = 1, n + out[i] = out[i] - overscan[i] - zeroval + } + + case D: # dark count + do i = 1, n + out[i] = out[i] - darkscale * dark[i] + case DO: # dark count + overscan + do i = 1, n + out[i] = out[i] - overscan[i] - darkscale * dark[i] + case DZ: # dark count + zero level + if (zeroim != NULL) + do i = 1, n + out[i] = out[i] - zero[i] - darkscale * dark[i] + else { + zeroval = zero[line] + do i = 1, n + out[i] = out[i] - zeroval - darkscale * dark[i] + } + case DZO: # dark count + zero level + overscan + if (zeroim != NULL) + do i = 1, n + out[i] = out[i] - overscan[i] - zero[i] - + darkscale * dark[i] + else { + zeroval = zero[line] + do i = 1, n + out[i] = out[i] - overscan[i] - zeroval - + darkscale * dark[i] + } + + case F: # flat field + if (flatim != NULL) { + do i = 1, n + out[i] = out[i] * flatscale / flat[i] + } else { + flatval = flatscale / flat[line] + do i = 1, n + out[i] = out[i] * flatval + } + case FO: # flat field + overscan + if (flatim != NULL) { + do i = 1, n + out[i] = (out[i] - overscan[i]) * flatscale / flat[i] + } else { + flatval = flatscale / flat[line] + do i = 1, n + out[i] = (out[i] - overscan[i]) * flatval + } + case FZ: # flat field + zero level + if (flatim != NULL) { + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - zero[i]) * flatscale / flat[i] + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - zeroval) * flatscale / flat[i] + } + } else { + flatval = flatscale / flat[line] + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - zero[i]) * flatval + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - zeroval) * flatval + } + } + case FZO: # flat field + zero level + overscan + if (flatim != NULL) { + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - overscan[i] - zero[i]) * + flatscale / flat[i] + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - overscan[i] - zeroval) * + flatscale / flat[i] + } + } else { + flatval = flatscale / flat[line] + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - overscan[i] - zero[i]) * flatval + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - overscan[i] - zeroval) * flatval + } + } + case FD: # flat field + dark count + if (flatim != NULL) { + do i = 1, n + out[i] = (out[i] - darkscale * dark[i]) * flatscale/flat[i] + } else { + flatval = flatscale / flat[line] + do i = 1, n + out[i] = (out[i] - darkscale * dark[i]) * flatval + } + case FDO: # flat field + dark count + overscan + if (flatim != NULL) { + do i = 1, n + out[i] = (out[i] - overscan[i] - darkscale * dark[i]) * + flatscale / flat[i] + } else { + flatval = flatscale / flat[line] + do i = 1, n + out[i] = (out[i] - overscan[i] - darkscale * dark[i]) * + flatval + } + case FDZ: # flat field + dark count + zero level + if (flatim != NULL) { + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - zero[i] - darkscale * dark[i]) * + flatscale / flat[i] + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - zeroval - darkscale * dark[i]) * + flatscale / flat[i] + } + } else { + flatval = flatscale / flat[line] + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - zero[i] - darkscale * dark[i]) * + flatval + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - zeroval - darkscale * dark[i]) * + flatval + } + } + case FDZO: # flat field + dark count + zero level + overscan + if (flatim != NULL) { + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - overscan[i] - zero[i] - + darkscale * dark[i]) * flatscale / flat[i] + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - overscan[i] - zeroval - + darkscale * dark[i]) * flatscale / flat[i] + } + } else { + flatval = flatscale / flat[line] + if (zeroim != NULL) { + do i = 1, n + out[i] = (out[i] - overscan[i] - zero[i] - + darkscale * dark[i]) * flatval + } else { + zeroval = zero[line] + do i = 1, n + out[i] = (out[i] - overscan[i] - zeroval - + darkscale * dark[i]) * flatval + } + } + } + + # Often these operations will not be performed so test for no + # correction rather than go through the switch. + + op = cors[ILLUMCOR] + cors[FRINGECOR] + if (op != 0) { + switch (op) { + case I: # illumination + do i = 1, n + out[i] = out[i] * illumscale / illum[i] + case Q: # fringe + do i = 1, n + out[i] = out[i] - frgscale * fringe[i] + case QI: # fringe + illumination + do i = 1, n + out[i] = out[i]*illumscale/illum[i] - frgscale*fringe[i] + } + } +end + diff --git a/noao/imred/quadred/src/ccdproc/generic/corinput.x b/noao/imred/quadred/src/ccdproc/generic/corinput.x new file mode 100644 index 00000000..07afaa41 --- /dev/null +++ b/noao/imred/quadred/src/ccdproc/generic/corinput.x @@ -0,0 +1,436 @@ +include +include "ccdred.h" + + +# CORINPUT -- Get an input image line, fix the bad pixels, and trim. +# Return the corrected input line in the output array. + +procedure corinputs (in, line, ccd, output, ncols) + +pointer in # Input IMIO pointer +int line # Corrected output line +pointer ccd # CCD pointer +short output[ncols] # Output data (returned) +int ncols # Number of output columns + +int i, inline +pointer inbuf, imgl2s() + +begin + # Determine the input line in terms of the trimmed output line. + if (IN_SEC(ccd) == NULL) + inline = IN_L1(ccd) + line - 1 + else { + do i = 1, IN_NSEC(ccd) { + if (line < OUT_SL1(ccd,i) || line > OUT_SL2(ccd,i)) + next + inline = IN_SL1(ccd,i) + line - OUT_SL1(ccd,i) + break + } + } + + # If there are bad lines call a procedure to fix them. Otherwise + # read the image line directly. + + if (NBADLINES(ccd) != 0) + call lfixs (in, inline, Mems[BADLINES(ccd)], IM_LEN(in,1), + IM_LEN(in,2), NBADLINES(ccd), inbuf) + else + inbuf = imgl2s (in, inline) + + # IF there are bad columns call a procedure to fix them. + if (NBADCOLS(ccd) != 0) + call cfixs (inline, Mems[BADCOLS(ccd)], IM_LEN(in,1), + IM_LEN(in,2), NBADCOLS(ccd), Mems[inbuf]) + + # Move the pixels to the output line. + if (IN_SEC(ccd) == NULL) + call amovs (Mems[inbuf+IN_C1(ccd)-OUT_C1(ccd)], output, ncols) + else { + do i = 1, IN_NSEC(ccd) { + if (inline < IN_SL1(ccd,i) || inline > IN_SL2(ccd,i)) + next + call amovs (Mems[inbuf+IN_SC1(ccd,i)-OUT_C1(ccd)], + output[OUT_SC1(ccd,i)], OUT_SC2(ccd,i)-OUT_SC1(ccd,i)+1) + } + } +end + + +# CFIX -- Interpolate across bad columns defined in the bad column array. + +procedure cfixs (line, badcols, ncols, nlines, nbadcols, data) + +int line # Line to be fixed +short badcols[2, nlines, nbadcols] # Bad column array +int ncols # Number of columns +int nlines # Number of lines +int nbadcols # Number of bad column regions +short data[ncols] # Data to be fixed + +short val +real del +int i, j, col1, col2 + +begin + do i = 1, nbadcols { + col1 = badcols[1, line, i] + if (col1 == 0) # No bad columns + return + col2 = badcols[2, line, i] + if (col1 == 1) { # Bad first column + val = data[col2+1] + do j = col1, col2 + data[j] = val + } else if (col2 == ncols) { # Bad last column + val = data[col1-1] + do j = col1, col2 + data[j] = val + } else { # Interpolate + del = (data[col2+1] - data[col1-1]) / (col2 - col1 + 2) + val = data[col1-1] + del + do j = col1, col2 + data[j] = val + (j - col1) * del + } + } +end + + +# LFIX -- Get image line and replace bad pixels by interpolation from +# neighboring lines. Internal buffers are used to keep the last fixed +# line and the next good line. They are allocated with LFIXINIT and +# freed with LFIXFREE. + +procedure lfixs (im, line, badlines, ncols, nlines, nbadlines, data) + +pointer im # IMIO pointer +int line # Line to be obtained and fixed +short badlines[2,nlines,nbadlines] # Bad line region array +int ncols # Number of columns in image +int nlines # Number of lines in images +int nbadlines # Number of bad line regions +pointer data # Data line pointer (returned) + +real wt1, wt2 +int i, nextgood, lastgood, col1, col2 +pointer imgl2s() + +pointer lastbuf, nextbuf +common /lfixcom/ lastbuf, nextbuf, lastgood + +begin + # If this line has bad pixels replace them. Otherwise just + # read the line. + + if (badlines[1, line, 1] != 0) { + # Save the last line which has already been fixed. + if (line != 1) + call amovs (Mems[data], Mems[lastbuf], ncols) + + # Determine the next line with no bad line pixels. Note that + # this requirement is overly strict since the bad columns + # may not be the same in neighboring lines. + + nextgood = 0 + do i = line+1, nlines { + if (badlines[1, i, 1] == 0) { + nextgood = i + break + } + } + + # If the next good line is not the same as previously + # read the data line and store it in a buffer. + + if ((nextgood != lastgood) && (nextgood != 0)) { + data = imgl2s (im, nextgood) + call amovs (Mems[data], Mems[nextbuf], ncols) + lastgood = nextgood + } + + # Get the data line. + data = imgl2s (im, line) + + # Interpolate the bad columns. At the ends of the image use + # extension otherwise use linear interpolation. + + if (line == 1) { # First line is bad + do i = 1, nbadlines { + col1 = badlines[1,line,i] - 1 + if (col1 == -1) + break + col2 = badlines[2,line,i] + call amovs (Mems[nextbuf+col1], Mems[data+col1], + col2-col1) + } + } else if (nextgood == 0) { # Last line is bad + do i = 1, nbadlines { + col1 = badlines[1,line,i] - 1 + if (col1 == -1) + break + col2 = badlines[2,line,i] + call amovs (Mems[lastbuf+col1], Mems[data+col1], + col2-col1) + } + } else { # Interpolate + wt1 = 1. / (nextgood - line + 1) + wt2 = 1. - wt1 + do i = 1, nbadlines { + col1 = badlines[1,line,i] - 1 + if (col1 == -1) + break + col2 = badlines[2,line,i] - 1 + call awsus (Mems[nextbuf+col1], Mems[lastbuf+col1], + Mems[data+col1], col2-col1+1, wt1, wt2) + } + } + } else + data = imgl2s (im, line) +end + + +# LFIXINIT -- Allocate internal buffers. + +procedure lfixinits (im) + +pointer im # IMIO pointer + +int lastgood +pointer lastbuf, nextbuf +common /lfixcom/ lastbuf, nextbuf, lastgood + +begin + call malloc (lastbuf, IM_LEN(im,1), TY_SHORT) + call malloc (nextbuf, IM_LEN(im,1), TY_SHORT) + lastgood=0 +end + +# LFIXFREE -- Free memory when the last line has been obtained. + +procedure lfixfrees () + +int lastgood +pointer lastbuf, nextbuf +common /lfixcom/ lastbuf, nextbuf, lastgood + +begin + call mfree (lastbuf, TY_SHORT) + call mfree (nextbuf, TY_SHORT) +end + +# CORINPUT -- Get an input image line, fix the bad pixels, and trim. +# Return the corrected input line in the output array. + +procedure corinputr (in, line, ccd, output, ncols) + +pointer in # Input IMIO pointer +int line # Corrected output line +pointer ccd # CCD pointer +real output[ncols] # Output data (returned) +int ncols # Number of output columns + +int i, inline +pointer inbuf, imgl2r() + +begin + # Determine the input line in terms of the trimmed output line. + if (IN_SEC(ccd) == NULL) + inline = IN_L1(ccd) + line - 1 + else { + do i = 1, IN_NSEC(ccd) { + if (line < OUT_SL1(ccd,i) || line > OUT_SL2(ccd,i)) + next + inline = IN_SL1(ccd,i) + line - OUT_SL1(ccd,i) + break + } + } + + # If there are bad lines call a procedure to fix them. Otherwise + # read the image line directly. + + if (NBADLINES(ccd) != 0) + call lfixr (in, inline, Mems[BADLINES(ccd)], IM_LEN(in,1), + IM_LEN(in,2), NBADLINES(ccd), inbuf) + else + inbuf = imgl2r (in, inline) + + # IF there are bad columns call a procedure to fix them. + if (NBADCOLS(ccd) != 0) + call cfixr (inline, Mems[BADCOLS(ccd)], IM_LEN(in,1), + IM_LEN(in,2), NBADCOLS(ccd), Memr[inbuf]) + + # Move the pixels to the output line. + if (IN_SEC(ccd) == NULL) + call amovr (Memr[inbuf+IN_C1(ccd)-OUT_C1(ccd)], output, ncols) + else { + do i = 1, IN_NSEC(ccd) { + if (inline < IN_SL1(ccd,i) || inline > IN_SL2(ccd,i)) + next + call amovr (Memr[inbuf+IN_SC1(ccd,i)-OUT_C1(ccd)], + output[OUT_SC1(ccd,i)], OUT_SC2(ccd,i)-OUT_SC1(ccd,i)+1) + } + } +end + + +# CFIX -- Interpolate across bad columns defined in the bad column array. + +procedure cfixr (line, badcols, ncols, nlines, nbadcols, data) + +int line # Line to be fixed +short badcols[2, nlines, nbadcols] # Bad column array +int ncols # Number of columns +int nlines # Number of lines +int nbadcols # Number of bad column regions +real data[ncols] # Data to be fixed + +real val +real del +int i, j, col1, col2 + +begin + do i = 1, nbadcols { + col1 = badcols[1, line, i] + if (col1 == 0) # No bad columns + return + col2 = badcols[2, line, i] + if (col1 == 1) { # Bad first column + val = data[col2+1] + do j = col1, col2 + data[j] = val + } else if (col2 == ncols) { # Bad last column + val = data[col1-1] + do j = col1, col2 + data[j] = val + } else { # Interpolate + del = (data[col2+1] - data[col1-1]) / (col2 - col1 + 2) + val = data[col1-1] + del + do j = col1, col2 + data[j] = val + (j - col1) * del + } + } +end + + +# LFIX -- Get image line and replace bad pixels by interpolation from +# neighboring lines. Internal buffers are used to keep the last fixed +# line and the next good line. They are allocated with LFIXINIT and +# freed with LFIXFREE. + +procedure lfixr (im, line, badlines, ncols, nlines, nbadlines, data) + +pointer im # IMIO pointer +int line # Line to be obtained and fixed +short badlines[2,nlines,nbadlines] # Bad line region array +int ncols # Number of columns in image +int nlines # Number of lines in images +int nbadlines # Number of bad line regions +pointer data # Data line pointer (returned) + +real wt1, wt2 +int i, nextgood, lastgood, col1, col2 +pointer imgl2r() + +pointer lastbuf, nextbuf +common /lfixcom/ lastbuf, nextbuf, lastgood + +begin + # If this line has bad pixels replace them. Otherwise just + # read the line. + + if (badlines[1, line, 1] != 0) { + # Save the last line which has already been fixed. + if (line != 1) + call amovr (Memr[data], Memr[lastbuf], ncols) + + # Determine the next line with no bad line pixels. Note that + # this requirement is overly strict since the bad columns + # may not be the same in neighboring lines. + + nextgood = 0 + do i = line+1, nlines { + if (badlines[1, i, 1] == 0) { + nextgood = i + break + } + } + + # If the next good line is not the same as previously + # read the data line and store it in a buffer. + + if ((nextgood != lastgood) && (nextgood != 0)) { + data = imgl2r (im, nextgood) + call amovr (Memr[data], Memr[nextbuf], ncols) + lastgood = nextgood + } + + # Get the data line. + data = imgl2r (im, line) + + # Interpolate the bad columns. At the ends of the image use + # extension otherwise use linear interpolation. + + if (line == 1) { # First line is bad + do i = 1, nbadlines { + col1 = badlines[1,line,i] - 1 + if (col1 == -1) + break + col2 = badlines[2,line,i] + call amovr (Memr[nextbuf+col1], Memr[data+col1], + col2-col1) + } + } else if (nextgood == 0) { # Last line is bad + do i = 1, nbadlines { + col1 = badlines[1,line,i] - 1 + if (col1 == -1) + break + col2 = badlines[2,line,i] + call amovr (Memr[lastbuf+col1], Memr[data+col1], + col2-col1) + } + } else { # Interpolate + wt1 = 1. / (nextgood - line + 1) + wt2 = 1. - wt1 + do i = 1, nbadlines { + col1 = badlines[1,line,i] - 1 + if (col1 == -1) + break + col2 = badlines[2,line,i] - 1 + call awsur (Memr[nextbuf+col1], Memr[lastbuf+col1], + Memr[data+col1], col2-col1+1, wt1, wt2) + } + } + } else + data = imgl2r (im, line) +end + + +# LFIXINIT -- Allocate internal buffers. + +procedure lfixinitr (im) + +pointer im # IMIO pointer + +int lastgood +pointer lastbuf, nextbuf +common /lfixcom/ lastbuf, nextbuf, lastgood + +begin + call malloc (lastbuf, IM_LEN(im,1), TY_REAL) + call malloc (nextbuf, IM_LEN(im,1), TY_REAL) + lastgood=0 +end + +# LFIXFREE -- Free memory when the last line has been obtained. + +procedure lfixfreer () + +int lastgood +pointer lastbuf, nextbuf +common /lfixcom/ lastbuf, nextbuf, lastgood + +begin + call mfree (lastbuf, TY_REAL) + call mfree (nextbuf, TY_REAL) +end + diff --git a/noao/imred/quadred/src/ccdproc/generic/mkpkg b/noao/imred/quadred/src/ccdproc/generic/mkpkg new file mode 100644 index 00000000..0f12b368 --- /dev/null +++ b/noao/imred/quadred/src/ccdproc/generic/mkpkg @@ -0,0 +1,12 @@ +# Make CCDRED Package. + +$checkout libpkg.a ../ +$update libpkg.a +$checkin libpkg.a ../ +$exit + +libpkg.a: + cor.x ccdred.h + corinput.x ccdred.h + proc.x ccdred.h + ; diff --git a/noao/imred/quadred/src/ccdproc/generic/proc.x b/noao/imred/quadred/src/ccdproc/generic/proc.x new file mode 100644 index 00000000..0251f4f8 --- /dev/null +++ b/noao/imred/quadred/src/ccdproc/generic/proc.x @@ -0,0 +1,678 @@ +include +include "ccdred.h" + + +.help proc Feb87 noao.imred.ccdred +.nf ---------------------------------------------------------------------------- +proc -- Process CCD images + +These are the main CCD reduction procedures. There is one for each +readout axis (lines or columns) and one for short and real image data. +They apply corrections for bad pixels, overscan levels, zero levels, +dark counts, flat field response, illumination response, and fringe +effects. The image is also trimmed if it was mapped with an image +section. The mean value for the output image is computed when the flat +field or illumination image is processed to form the scale factor for +these calibrations in order to avoid reading through these image a +second time. + +The processing information and parameters are specified in the CCD +structure. The processing operations to be performed are specified by +the correction array CORS in the ccd structure. There is one array +element for each operation with indices defined symbolically by macro +definitions (see ccdred.h); i.e. FLATCOR. The value of the array +element is an integer bit field in which the bit set is the same as the +array index; i.e element 3 will have the third bit set for an operation +with array value 2**(3-1)=4. If an operation is not to be performed +the bit is not set and the array element has the numeric value zero. +Note that the addition of several correction elements gives a unique +bit field describing a combination of operations. For efficiency the +most common combinations are implemented as separate units. + +The CCD structure also contains the correction or calibration data +consisting either pointers to data, IMIO pointers for the calibration +images, and scale factors. + +The processing is performed line-by-line. The procedure CORINPUT is +called to get an input line. This procedure trims and fixes bad pixels by +interpolation. The output line and lines from the various calibration +images are read. The image vectors as well as the overscan vector and +the scale factors are passed to the procedure COR (which also +dereferences the pointer data into simple arrays and variables). That +procedure does the actual corrections apart from bad pixel +corrections. + +The final optional step is to add each corrected output line to form a +mean. This adds efficiency since the operation is done only if desired +and the output image data is already in memory so there is no I/O +penalty. + +SEE ALSO + ccdred.h, cor, fixpix, setfixpix, setoverscan, settrim, + setzero, setdark, setflat, setillum, setfringe +.endhelp ---------------------------------------------------------------------- + + + +# PROC1 -- Process CCD images with readout axis 1 (lines). + +procedure proc1s (ccd) + +pointer ccd # CCD structure + +int i, line, nlin, ncols, nlines, findmean, rep +int c1, c2, l1, l2 +real overscan, darkscale, flatscale, illumscale, frgscale, mean +short minrep +pointer in, out, zeroim, darkim, flatim, illumim, fringeim +pointer outbuf, overscan_vec, zerobuf, darkbuf, flatbuf, illumbuf, fringebuf + +real asums() +pointer imgl2s(), impl2s(), ccd_gls() + +begin + # Initialize. If the correction image is 1D then just get the + # data once. + + in = IN_IM(ccd) + out = OUT_IM(ccd) + nlin = IM_LEN(in,2) + ncols = OUT_C2(ccd) - OUT_C1(ccd) + 1 + nlines = OUT_L2(ccd) - OUT_L1(ccd) + 1 + + if (CORS(ccd, FIXPIX) == YES) + call lfixinits (in) + + findmean = CORS(ccd, FINDMEAN) + if (findmean == YES) + mean = 0. + rep = CORS(ccd, MINREP) + if (rep == YES) + minrep = MINREPLACE(ccd) + + overscan_vec = OVERSCAN_VEC(ccd) + + if (CORS(ccd, ZEROCOR) == 0) { + zeroim = NULL + zerobuf = 1 + } else if (IM_LEN(ZERO_IM(ccd),2) == 1) { + zeroim = NULL + zerobuf = ccd_gls (ZERO_IM(ccd), ZERO_C1(ccd), ZERO_C2(ccd), 1) + } else + zeroim = ZERO_IM(ccd) + + if (CORS(ccd, DARKCOR) == 0) { + darkim = NULL + darkbuf = 1 + } else if (IM_LEN(DARK_IM(ccd),2) == 1) { + darkim = NULL + darkbuf = ccd_gls (DARK_IM(ccd), DARK_C1(ccd), DARK_C2(ccd), 1) + darkscale = FLATSCALE(ccd) + } else { + darkim = DARK_IM(ccd) + darkscale = DARKSCALE(ccd) + } + + if (CORS(ccd, FLATCOR) == 0) { + flatim = NULL + flatbuf = 1 + } else if (IM_LEN(FLAT_IM(ccd),2) == 1) { + flatim = NULL + flatbuf = ccd_gls (FLAT_IM(ccd), FLAT_C1(ccd), FLAT_C2(ccd), 1) + flatscale = FLATSCALE(ccd) + } else { + flatim = FLAT_IM(ccd) + flatscale = FLATSCALE(ccd) + } + + if (CORS(ccd, ILLUMCOR) == 0) { + illumim = NULL + illumbuf = 1 + } else { + illumim = ILLUM_IM(ccd) + illumscale = ILLUMSCALE(ccd) + } + + if (CORS(ccd, FRINGECOR) == 0) { + fringeim = NULL + fringebuf = 1 + } else { + fringeim = FRINGE_IM(ccd) + frgscale = FRINGESCALE(ccd) + } + + # For each line read lines from the input. Procedure CORINPUT + # replaces bad pixels by interpolation and applies a trim to the + # input. Get lines from the output image and from the zero level, + # dark count, flat field, illumination, and fringe images. + # Call COR1 to do the actual pixel corrections. Finally, add the + # output pixels to a sum for computing the mean. + # We must copy data outside of the output data section. + + do line = 2 - OUT_L1(ccd), 0 + call amovs ( + Mems[imgl2s(in,IN_L1(ccd)+line-1)+IN_C1(ccd)-OUT_C1(ccd)], + Mems[impl2s(out,OUT_L1(ccd)+line-1)], IM_LEN(out,1)) + + do line = 1, nlines { + outbuf = impl2s (out, OUT_L1(ccd)+line-1) + call corinputs (in, line, ccd, Mems[outbuf], IM_LEN(out,1)) + + outbuf = outbuf + OUT_C1(ccd) - 1 + if (overscan_vec != NULL) + overscan = Memr[overscan_vec+line-1] + if (zeroim != NULL) + zerobuf = ccd_gls (zeroim, ZERO_C1(ccd), ZERO_C2(ccd), + ZERO_L1(ccd)+line-1) + if (darkim != NULL) + darkbuf = ccd_gls (darkim, DARK_C1(ccd), DARK_C2(ccd), + DARK_L1(ccd)+line-1) + if (flatim != NULL) + flatbuf = ccd_gls (flatim, FLAT_C1(ccd), FLAT_C2(ccd), + FLAT_L1(ccd)+line-1) + if (illumim != NULL) + illumbuf = ccd_gls (illumim, ILLUM_C1(ccd), ILLUM_C2(ccd), + ILLUM_L1(ccd)+line-1) + if (fringeim != NULL) + fringebuf = ccd_gls (fringeim, FRINGE_C1(ccd), FRINGE_C2(ccd), + FRINGE_L1(ccd)+line-1) + + if (OUT_SEC(ccd) == NULL) { + call cor1s (CORS(ccd,1), Mems[outbuf], + overscan, Mems[zerobuf], Mems[darkbuf], + Mems[flatbuf], Mems[illumbuf], Mems[fringebuf], ncols, + darkscale, flatscale, illumscale, frgscale) + } else { + do i = 1, IN_NSEC(ccd) { + l1 = OUT_SL1(ccd,i) + l2 = OUT_SL2(ccd,i) + if (line < l1 || line > l2) + next + c1 = OUT_SC1(ccd,i) - 1 + c2 = OUT_SC2(ccd,i) - 1 + ncols = c2 - c1 + 1 + if (overscan_vec != NULL) + overscan = Memr[overscan_vec+(i-1)*nlin+line-l1] + + call cor1s (CORS(ccd,1), Mems[outbuf+c1], + overscan, Mems[zerobuf+c1], Mems[darkbuf+c1], + Mems[flatbuf+c1], Mems[illumbuf+c1], + Mems[fringebuf+c1], ncols, + darkscale, flatscale, illumscale, frgscale) + } + } + + if (rep == YES) + call amaxks (Mems[outbuf], minrep, Mems[outbuf], ncols) + if (findmean == YES) + mean = mean + asums (Mems[outbuf], ncols) + } + + do line = nlines+1, IM_LEN(out,2)-OUT_L1(ccd)+1 + call amovs ( + Mems[imgl2s(in,IN_L1(ccd)+line-1)+IN_C1(ccd)-OUT_C1(ccd)], + Mems[impl2s(out,OUT_L1(ccd)+line-1)], IM_LEN(out,1)) + + # Compute the mean from the sum of the output pixels. + if (findmean == YES) + MEAN(ccd) = mean / ncols / nlines + + if (CORS(ccd, FIXPIX) == YES) + call lfixfrees () +end + + +# PROC2 -- Process CCD images with readout axis 2 (columns). + +procedure proc2s (ccd) + +pointer ccd # CCD structure + +int line, ncols, nlines, findmean, rep +real darkscale, flatscale, illumscale, frgscale, mean +short minrep +pointer in, out, zeroim, darkim, flatim, illumim, fringeim +pointer outbuf, overscan_vec, zerobuf, darkbuf, flatbuf, illumbuf, fringebuf + +real asums() +pointer imgl2s(), impl2s(), imgs2s(), ccd_gls() + +begin + # Initialize. If the correction image is 1D then just get the + # data once. + + in = IN_IM(ccd) + out = OUT_IM(ccd) + ncols = OUT_C2(ccd) - OUT_C1(ccd) + 1 + nlines = OUT_L2(ccd) - OUT_L1(ccd) + 1 + + if (CORS(ccd, FIXPIX) == YES) + call lfixinits (in) + + findmean = CORS(ccd, FINDMEAN) + if (findmean == YES) + mean = 0. + rep = CORS(ccd, MINREP) + if (rep == YES) + minrep = MINREPLACE(ccd) + + overscan_vec = OVERSCAN_VEC(ccd) + + if (CORS(ccd, ZEROCOR) == 0) { + zeroim = NULL + zerobuf = 1 + } else if (IM_LEN(ZERO_IM(ccd),1) == 1) { + zeroim = NULL + zerobuf = imgs2s (ZERO_IM(ccd), 1, 1, ZERO_L1(ccd), ZERO_L2(ccd)) + } else + zeroim = ZERO_IM(ccd) + + if (CORS(ccd, DARKCOR) == 0) { + darkim = NULL + darkbuf = 1 + } else if (IM_LEN(DARK_IM(ccd),1) == 1) { + darkim = NULL + darkbuf = imgs2s (DARK_IM(ccd), 1, 1, DARK_L1(ccd), DARK_L2(ccd)) + darkscale = DARKSCALE(ccd) + } else { + darkim = DARK_IM(ccd) + darkscale = DARKSCALE(ccd) + } + + if (CORS(ccd, FLATCOR) == 0) { + flatim = NULL + flatbuf = 1 + } else if (IM_LEN(FLAT_IM(ccd),1) == 1) { + flatim = NULL + flatbuf = imgs2s (FLAT_IM(ccd), 1, 1, FLAT_L1(ccd), FLAT_L2(ccd)) + flatscale = FLATSCALE(ccd) + } else { + flatim = FLAT_IM(ccd) + flatscale = FLATSCALE(ccd) + } + + if (CORS(ccd, ILLUMCOR) == 0) { + illumim = NULL + illumbuf = 1 + } else { + illumim = ILLUM_IM(ccd) + illumscale = ILLUMSCALE(ccd) + } + + if (CORS(ccd, FRINGECOR) == 0) { + fringeim = NULL + fringebuf = 1 + } else { + fringeim = FRINGE_IM(ccd) + frgscale = FRINGESCALE(ccd) + } + + # For each line read lines from the input. Procedure CORINPUT + # replaces bad pixels by interpolation and applies a trim to the + # input. Get lines from the output image and from the zero level, + # dark count, flat field, illumination, and fringe images. + # Call COR2 to do the actual pixel corrections. Finally, add the + # output pixels to a sum for computing the mean. + # We must copy data outside of the output data section. + + do line = 2 - OUT_L1(ccd), 0 + call amovs ( + Mems[imgl2s(in,IN_L1(ccd)+line-1)+IN_C1(ccd)-OUT_C1(ccd)], + Mems[impl2s(out,OUT_L1(ccd)+line-1)], IM_LEN(out,1)) + + do line = 1, nlines { + outbuf = impl2s (out, OUT_L1(ccd)+line-1) + call corinputs (in, line, ccd, Mems[outbuf], IM_LEN(out,1)) + + outbuf = outbuf + OUT_C1(ccd) - 1 + if (zeroim != NULL) + zerobuf = ccd_gls (zeroim, ZERO_C1(ccd), ZERO_C2(ccd), + ZERO_L1(ccd)+line-1) + if (darkim != NULL) + darkbuf = ccd_gls (darkim, DARK_C1(ccd), DARK_C2(ccd), + DARK_L1(ccd)+line-1) + if (flatim != NULL) + flatbuf = ccd_gls (flatim, FLAT_C1(ccd), FLAT_C2(ccd), + FLAT_L1(ccd)+line-1) + if (illumim != NULL) + illumbuf = ccd_gls (illumim, ILLUM_C1(ccd), ILLUM_C2(ccd), + ILLUM_L1(ccd)+line-1) + if (fringeim != NULL) + fringebuf = ccd_gls (fringeim, FRINGE_C1(ccd), FRINGE_C2(ccd), + FRINGE_L1(ccd)+line-1) + + call cor2s (line, CORS(ccd,1), Mems[outbuf], + Memr[overscan_vec], Mems[zerobuf], Mems[darkbuf], + Mems[flatbuf], Mems[illumbuf], Mems[fringebuf], ncols, + zeroim, flatim, darkscale, flatscale, illumscale, frgscale) + + if (rep == YES) + call amaxks (Mems[outbuf], minrep, Mems[outbuf], ncols) + if (findmean == YES) + mean = mean + asums (Mems[outbuf], ncols) + } + + do line = nlines+1, IM_LEN(out,2)-OUT_L1(ccd)+1 + call amovs ( + Mems[imgl2s(in,IN_L1(ccd)+line-1)+IN_C1(ccd)-OUT_C1(ccd)], + Mems[impl2s(out,OUT_L1(ccd)+line-1)], IM_LEN(out,1)) + + # Compute the mean from the sum of the output pixels. + if (findmean == YES) + MEAN(ccd) = mean / ncols / nlines + + if (CORS(ccd, FIXPIX) == YES) + call lfixfrees () +end + +# PROC1 -- Process CCD images with readout axis 1 (lines). + +procedure proc1r (ccd) + +pointer ccd # CCD structure + +int i, line, nlin, ncols, nlines, findmean, rep +int c1, c2, l1, l2 +real overscan, darkscale, flatscale, illumscale, frgscale, mean +real minrep +pointer in, out, zeroim, darkim, flatim, illumim, fringeim +pointer outbuf, overscan_vec, zerobuf, darkbuf, flatbuf, illumbuf, fringebuf + +real asumr() +pointer imgl2r(), impl2r(), ccd_glr() + +begin + # Initialize. If the correction image is 1D then just get the + # data once. + + in = IN_IM(ccd) + out = OUT_IM(ccd) + nlin = IM_LEN(in,2) + ncols = OUT_C2(ccd) - OUT_C1(ccd) + 1 + nlines = OUT_L2(ccd) - OUT_L1(ccd) + 1 + + if (CORS(ccd, FIXPIX) == YES) + call lfixinitr (in) + + findmean = CORS(ccd, FINDMEAN) + if (findmean == YES) + mean = 0. + rep = CORS(ccd, MINREP) + if (rep == YES) + minrep = MINREPLACE(ccd) + + overscan_vec = OVERSCAN_VEC(ccd) + + if (CORS(ccd, ZEROCOR) == 0) { + zeroim = NULL + zerobuf = 1 + } else if (IM_LEN(ZERO_IM(ccd),2) == 1) { + zeroim = NULL + zerobuf = ccd_glr (ZERO_IM(ccd), ZERO_C1(ccd), ZERO_C2(ccd), 1) + } else + zeroim = ZERO_IM(ccd) + + if (CORS(ccd, DARKCOR) == 0) { + darkim = NULL + darkbuf = 1 + } else if (IM_LEN(DARK_IM(ccd),2) == 1) { + darkim = NULL + darkbuf = ccd_glr (DARK_IM(ccd), DARK_C1(ccd), DARK_C2(ccd), 1) + darkscale = FLATSCALE(ccd) + } else { + darkim = DARK_IM(ccd) + darkscale = DARKSCALE(ccd) + } + + if (CORS(ccd, FLATCOR) == 0) { + flatim = NULL + flatbuf = 1 + } else if (IM_LEN(FLAT_IM(ccd),2) == 1) { + flatim = NULL + flatbuf = ccd_glr (FLAT_IM(ccd), FLAT_C1(ccd), FLAT_C2(ccd), 1) + flatscale = FLATSCALE(ccd) + } else { + flatim = FLAT_IM(ccd) + flatscale = FLATSCALE(ccd) + } + + if (CORS(ccd, ILLUMCOR) == 0) { + illumim = NULL + illumbuf = 1 + } else { + illumim = ILLUM_IM(ccd) + illumscale = ILLUMSCALE(ccd) + } + + if (CORS(ccd, FRINGECOR) == 0) { + fringeim = NULL + fringebuf = 1 + } else { + fringeim = FRINGE_IM(ccd) + frgscale = FRINGESCALE(ccd) + } + + # For each line read lines from the input. Procedure CORINPUT + # replaces bad pixels by interpolation and applies a trim to the + # input. Get lines from the output image and from the zero level, + # dark count, flat field, illumination, and fringe images. + # Call COR1 to do the actual pixel corrections. Finally, add the + # output pixels to a sum for computing the mean. + # We must copy data outside of the output data section. + + do line = 2 - OUT_L1(ccd), 0 + call amovr ( + Memr[imgl2r(in,IN_L1(ccd)+line-1)+IN_C1(ccd)-OUT_C1(ccd)], + Memr[impl2r(out,OUT_L1(ccd)+line-1)], IM_LEN(out,1)) + + do line = 1, nlines { + outbuf = impl2r (out, OUT_L1(ccd)+line-1) + call corinputr (in, line, ccd, Memr[outbuf], IM_LEN(out,1)) + + outbuf = outbuf + OUT_C1(ccd) - 1 + if (overscan_vec != NULL) + overscan = Memr[overscan_vec+line-1] + if (zeroim != NULL) + zerobuf = ccd_glr (zeroim, ZERO_C1(ccd), ZERO_C2(ccd), + ZERO_L1(ccd)+line-1) + if (darkim != NULL) + darkbuf = ccd_glr (darkim, DARK_C1(ccd), DARK_C2(ccd), + DARK_L1(ccd)+line-1) + if (flatim != NULL) + flatbuf = ccd_glr (flatim, FLAT_C1(ccd), FLAT_C2(ccd), + FLAT_L1(ccd)+line-1) + if (illumim != NULL) + illumbuf = ccd_glr (illumim, ILLUM_C1(ccd), ILLUM_C2(ccd), + ILLUM_L1(ccd)+line-1) + if (fringeim != NULL) + fringebuf = ccd_glr (fringeim, FRINGE_C1(ccd), FRINGE_C2(ccd), + FRINGE_L1(ccd)+line-1) + + if (OUT_SEC(ccd) == NULL) { + call cor1r (CORS(ccd,1), Memr[outbuf], + overscan, Memr[zerobuf], Memr[darkbuf], + Memr[flatbuf], Memr[illumbuf], Memr[fringebuf], ncols, + darkscale, flatscale, illumscale, frgscale) + } else { + do i = 1, IN_NSEC(ccd) { + l1 = OUT_SL1(ccd,i) + l2 = OUT_SL2(ccd,i) + if (line < l1 || line > l2) + next + c1 = OUT_SC1(ccd,i) - 1 + c2 = OUT_SC2(ccd,i) - 1 + ncols = c2 - c1 + 1 + if (overscan_vec != NULL) + overscan = Memr[overscan_vec+(i-1)*nlin+line-l1] + + call cor1r (CORS(ccd,1), Memr[outbuf+c1], + overscan, Memr[zerobuf+c1], Memr[darkbuf+c1], + Memr[flatbuf+c1], Memr[illumbuf+c1], + Memr[fringebuf+c1], ncols, + darkscale, flatscale, illumscale, frgscale) + } + } + + if (rep == YES) + call amaxkr (Memr[outbuf], minrep, Memr[outbuf], ncols) + if (findmean == YES) + mean = mean + asumr (Memr[outbuf], ncols) + } + + do line = nlines+1, IM_LEN(out,2)-OUT_L1(ccd)+1 + call amovr ( + Memr[imgl2r(in,IN_L1(ccd)+line-1)+IN_C1(ccd)-OUT_C1(ccd)], + Memr[impl2r(out,OUT_L1(ccd)+line-1)], IM_LEN(out,1)) + + # Compute the mean from the sum of the output pixels. + if (findmean == YES) + MEAN(ccd) = mean / ncols / nlines + + if (CORS(ccd, FIXPIX) == YES) + call lfixfreer () +end + + +# PROC2 -- Process CCD images with readout axis 2 (columns). + +procedure proc2r (ccd) + +pointer ccd # CCD structure + +int line, ncols, nlines, findmean, rep +real darkscale, flatscale, illumscale, frgscale, mean +real minrep +pointer in, out, zeroim, darkim, flatim, illumim, fringeim +pointer outbuf, overscan_vec, zerobuf, darkbuf, flatbuf, illumbuf, fringebuf + +real asumr() +pointer imgl2r(), impl2r(), imgs2r(), ccd_glr() + +begin + # Initialize. If the correction image is 1D then just get the + # data once. + + in = IN_IM(ccd) + out = OUT_IM(ccd) + ncols = OUT_C2(ccd) - OUT_C1(ccd) + 1 + nlines = OUT_L2(ccd) - OUT_L1(ccd) + 1 + + if (CORS(ccd, FIXPIX) == YES) + call lfixinitr (in) + + findmean = CORS(ccd, FINDMEAN) + if (findmean == YES) + mean = 0. + rep = CORS(ccd, MINREP) + if (rep == YES) + minrep = MINREPLACE(ccd) + + overscan_vec = OVERSCAN_VEC(ccd) + + if (CORS(ccd, ZEROCOR) == 0) { + zeroim = NULL + zerobuf = 1 + } else if (IM_LEN(ZERO_IM(ccd),1) == 1) { + zeroim = NULL + zerobuf = imgs2r (ZERO_IM(ccd), 1, 1, ZERO_L1(ccd), ZERO_L2(ccd)) + } else + zeroim = ZERO_IM(ccd) + + if (CORS(ccd, DARKCOR) == 0) { + darkim = NULL + darkbuf = 1 + } else if (IM_LEN(DARK_IM(ccd),1) == 1) { + darkim = NULL + darkbuf = imgs2r (DARK_IM(ccd), 1, 1, DARK_L1(ccd), DARK_L2(ccd)) + darkscale = DARKSCALE(ccd) + } else { + darkim = DARK_IM(ccd) + darkscale = DARKSCALE(ccd) + } + + if (CORS(ccd, FLATCOR) == 0) { + flatim = NULL + flatbuf = 1 + } else if (IM_LEN(FLAT_IM(ccd),1) == 1) { + flatim = NULL + flatbuf = imgs2r (FLAT_IM(ccd), 1, 1, FLAT_L1(ccd), FLAT_L2(ccd)) + flatscale = FLATSCALE(ccd) + } else { + flatim = FLAT_IM(ccd) + flatscale = FLATSCALE(ccd) + } + + if (CORS(ccd, ILLUMCOR) == 0) { + illumim = NULL + illumbuf = 1 + } else { + illumim = ILLUM_IM(ccd) + illumscale = ILLUMSCALE(ccd) + } + + if (CORS(ccd, FRINGECOR) == 0) { + fringeim = NULL + fringebuf = 1 + } else { + fringeim = FRINGE_IM(ccd) + frgscale = FRINGESCALE(ccd) + } + + # For each line read lines from the input. Procedure CORINPUT + # replaces bad pixels by interpolation and applies a trim to the + # input. Get lines from the output image and from the zero level, + # dark count, flat field, illumination, and fringe images. + # Call COR2 to do the actual pixel corrections. Finally, add the + # output pixels to a sum for computing the mean. + # We must copy data outside of the output data section. + + do line = 2 - OUT_L1(ccd), 0 + call amovr ( + Memr[imgl2r(in,IN_L1(ccd)+line-1)+IN_C1(ccd)-OUT_C1(ccd)], + Memr[impl2r(out,OUT_L1(ccd)+line-1)], IM_LEN(out,1)) + + do line = 1, nlines { + outbuf = impl2r (out, OUT_L1(ccd)+line-1) + call corinputr (in, line, ccd, Memr[outbuf], IM_LEN(out,1)) + + outbuf = outbuf + OUT_C1(ccd) - 1 + if (zeroim != NULL) + zerobuf = ccd_glr (zeroim, ZERO_C1(ccd), ZERO_C2(ccd), + ZERO_L1(ccd)+line-1) + if (darkim != NULL) + darkbuf = ccd_glr (darkim, DARK_C1(ccd), DARK_C2(ccd), + DARK_L1(ccd)+line-1) + if (flatim != NULL) + flatbuf = ccd_glr (flatim, FLAT_C1(ccd), FLAT_C2(ccd), + FLAT_L1(ccd)+line-1) + if (illumim != NULL) + illumbuf = ccd_glr (illumim, ILLUM_C1(ccd), ILLUM_C2(ccd), + ILLUM_L1(ccd)+line-1) + if (fringeim != NULL) + fringebuf = ccd_glr (fringeim, FRINGE_C1(ccd), FRINGE_C2(ccd), + FRINGE_L1(ccd)+line-1) + + call cor2r (line, CORS(ccd,1), Memr[outbuf], + Memr[overscan_vec], Memr[zerobuf], Memr[darkbuf], + Memr[flatbuf], Memr[illumbuf], Memr[fringebuf], ncols, + zeroim, flatim, darkscale, flatscale, illumscale, frgscale) + + if (rep == YES) + call amaxkr (Memr[outbuf], minrep, Memr[outbuf], ncols) + if (findmean == YES) + mean = mean + asumr (Memr[outbuf], ncols) + } + + do line = nlines+1, IM_LEN(out,2)-OUT_L1(ccd)+1 + call amovr ( + Memr[imgl2r(in,IN_L1(ccd)+line-1)+IN_C1(ccd)-OUT_C1(ccd)], + Memr[impl2r(out,OUT_L1(ccd)+line-1)], IM_LEN(out,1)) + + # Compute the mean from the sum of the output pixels. + if (findmean == YES) + MEAN(ccd) = mean / ncols / nlines + + if (CORS(ccd, FIXPIX) == YES) + call lfixfreer () +end + -- cgit