Primitives/API

class kcwidrp.primitives.ArcOffsets.ArcOffsets(action, context)[source]

Derive offset of each bar relative to reference bar.

Using cross correlation techniques, this routine calculates the relative offsets between the reference bar (loaded as a parameter as config.instrument.REFBAR). The other important config parameters is TUKEYALPHA, which controls the roll-off of the edges of the spectra that are cross-correlated. If a bright line near the edge is throwing off the cross-correlations, you can increase this parameter to lessen its impact.

See kcwidrp/configs/kcwi.cfg

Arcs must be available as context.arcs.

class kcwidrp.primitives.CalcPrelimDisp.CalcPrelimDisp(action, context)[source]

Calculate dispersion based on configuration parameters.

The parameters of the grating equation are calculated as:

  • BLUE

    • alpha = grating_angle - 13 - adjustment_angle

    • adjustment_angle = (180 for BH and 0 for all other gratings)

  • RED

    • alpha = 155.892 - grating_angle

  • BLUE & RED

    • beta = camera_angle - alpha

dispersion = cos(beta)/rho/focal_length x (pixel_scale x binning) * 1.e4

class kcwidrp.primitives.CorrectDar.CorrectDar(action, context)[source]

Correct for Differential Atmospheric Refraction

Accounts for rotator orientation, zenith angle, and parallactic angle to correct input data cube into a padded, DAR corrected output cube. Calculates the DAR correction for each wavelength slice and adjusts position in cube using scipy.nddata.shift to implement correction.

Sets flags in padded region of output cube to a value of 128.

Also corrects delta wavelength cube if it is a standard star observation.

kcwidrp.primitives.CorrectDar.atm_disper(w0, w1, airmass, temperature=10.0, pressure_pa=61100.0, humidity=50.0, co2=400.0)[source]

Calculate atmospheric dispersion at w1 relative to w0

Parameters:
  • w0 (float) – reference wavelength (Angstroms)

  • w1 (float) – offset wavelength (Angstroms)

  • airmass (float) – unitless airmass

  • temperature (float) – atmospheric temperature (C)

  • pressure_pa (float) – atmospheric pressure (Pa)

  • humidity (float) – relative humidity (%)

  • co2 (float) – Carbon-Dioxide (mu-mole/mole)

class kcwidrp.primitives.CorrectDefects.CorrectDefects(action, context)[source]

Remove known bad columns.

Looks for a defect list file in the data directory of kcwidrp based on the CCD ampmode and x and y binning. Records the defect correction in the FITS header with the following keywords:

  • BPFILE: the bad pixel file used to correct defects

  • NBPCLEAN: the number of bad pixels cleaned

Uses the following configuration parameter:

  • saveintims: if set to True write out a *_def.fits file with defects corrected. Default is False.

Updates image in returned arguments.

class kcwidrp.primitives.CorrectGain.CorrectGain(action, context)[source]

Convert raw data numbers to electrons.

Uses the ATSECn FITS header keywords to divide image into amp regions and then corrects each region with the corresponding GAINn keyword. Updates the following FITS header keywords:

  • GAINCOR: sets to True if operation performed.

  • BUNIT: sets to electron.

  • HISTORY: records the operation.

Uses the following configuration parameter:

  • saveintims: if set to True write out a *_gain.fits file with gain corrected. Default is False.

Updates image in returned arguments.

class kcwidrp.primitives.CorrectIllumination.CorrectIllumination(action, context)[source]

Correct for illumination and response variations using master flat.

Currently, gives precedence for internal flats (cflat), as they are the most universally applicable. To use other flats move the cflats aside, in which case any twilight flat (tflat) master would then have precedence followed by any dome flat (dflat) master.

class kcwidrp.primitives.CreateUncertaintyImage.CreateUncertaintyImage(action, context)[source]

Generate a standard deviation uncertainty image.

Starts with pure Poisson noise and uses astropy.nddata.StdDevUncertainty to generate the uncertainty frame. If BIASRNn header keywords present, uses these with the ATSECn keywords to apply the appropriate readnoise addition to each amplifier region.

class kcwidrp.primitives.CubeImage.CubeImage(action, context)[source]

Transform 3D data cube into 2D image.

Generates a wavelength-aligned 2d image of a 3d data cube. This is done automatically for arc lamps used for the wavelength solution as a diagnostic for the wavelength fit.

Writes out image with “_icube_2d.fits” suffix.

class kcwidrp.primitives.ExtractArcs.ExtractArcs(action, context)[source]

Use derived traces to extract arc spectra along continuum bars.

Reads in traces from continuum bars and then uses them to extract spectra along each trace. Also performs a background subtraction for each extracted spectrum. Adds a HISTORY record to the FITS header.

Uses the following configuration parameter:

  • saveintims: if set to True write out a warped version of arc image in *_warped.fits. Defaults to False.

Adds extracted arcs to returned arguments.

class kcwidrp.primitives.FindBars.FindBars(action, context)[source]

Find bars in middle row of cont bars image

Uses the middle row with surrounding rows to find the individual continuum bar locations in the image. Checks to be sure the correct number of bars has been found (see NBARS in kcwidrp/configs/kcwi.cfg).

class kcwidrp.primitives.FitCenter.FitCenter(action, context)[source]

Perform wavelength fitting on central region.

Perform a preliminary determination of the dispersion in the central region of the wavelength range. Starts with known offsets between bars and rough offset between reference bar and atlas spectrum, and the calculated dispersion.

Uses config parameter TAPERFRAC to control cross-correlation roll-off.

kcwidrp.primitives.FitCenter.pascal_shift(coefficients=None, x0=None)[source]

Shift coefficients to a new reference value (X0)

Use Pascal’s Triangle to shift the polynomial coefficients to a new reference point. Used to transfer central fits to full-wavelength fits.

Parameters:
  • coefficients (list of floats) – fit coef’s, up to 7 elements.

  • x0 (float) – New reference x value.

Returns:

Shifted coefficients

class kcwidrp.primitives.FlagSaturation.FlagSaturation(action, context)[source]

Flag saturated pixels.

Currently flags pixels with values > 60,000 counts with a value of 8 in the flags FITS extension and updates the following FITS header keywords:

  • SATFLAG: set to True if operation is performed.

  • NSATFLAG: set to the count of saturated pixels.

Uses the following configuration parameter:

  • saveintims: if set to True write out a flagged version of image in *_fsat.fits. Defaults to False.

Updates the flag extension of the image in the returned arguments.

class kcwidrp.primitives.FluxCalibrate.FluxCalibrate(action, context)[source]

Perform flux calibration.

Uses inverse sensitivity curve derived from MakeInvsens to flux calibrate input observation.

class kcwidrp.primitives.GenerateMaps.GenerateMaps(action, context)[source]

Generate geometry map images.

Uses input geometry solution and outputs the following geometry maps where the pixel value in each image indicates:

  • slice - what slice the pixel is in

  • pos - what the pixel’s location is relative to the edge of the slice

  • wave - what wavelength the pixel has

  • del - what the delta wavelength is at the pixel

Pixels not on any slice are given as a negative value (usually -1.)

class kcwidrp.primitives.GetAtlasLines.GetAtlasLines(action, context)[source]

Get relevant atlas line positions and wavelengths

Generates a list of well observed atlas lines for generating the final wavelength solution.

Uses the following configuration parameters:

  • FRACMAX: fraction of max to use for window for fitting (defaults to 0.5)

  • LINELIST: an input line list to use instead of determining one on the fly

kcwidrp.primitives.GetAtlasLines.findpeaks(x, y, wid, sth, ath, pkg=None, verbose=False)[source]

Find peaks in spectrum

Uses a gradient to determine where peaks are and then performs some cleaning of blended lines using sigma rejection.

Parameters:
  • x (list of floats) – x values of vector

  • y (list of floats) – y values of vector

  • wid (int) – boxcar smoothing width in pixels

  • sth (float) – slope threshhold

  • ath (float) – amplitude threshhold

  • pkg (float) – limit on peak wandering in pixels

  • verbose (bool) – control output messages

Returns:

  • list of peaks

  • average sigma of cleaned peaks

  • list of y values at peaks

kcwidrp.primitives.GetAtlasLines.gaus(x, a, mu, sigma)[source]

Gaussian fitting function.

Parameters:
  • x (float) – x value

  • a (float) – amplitude

  • mu (float) – x position of center of Gaussian

  • sigma (float) – Gaussian width

Returns:

Gaussian value at x position

kcwidrp.primitives.GetAtlasLines.get_line_window(y, c, thresh=0.0, logger=None, strict=False, maxwin=100, frac_max=0.5)[source]

Find a window that includes the specified amount of the line.

Get a useful window on a given line to allow finding an accurate peak for each line. Start at input center and iterate above and below until the window encompasses the specified fraction of the line maximum flux. Rejects lines that are problematic.

Parameters:
  • y (list of floats) – flux values

  • c (float) – the nominal center of the line in pixels

  • thresh (float) – threshhold for rejection (default: 0.)

  • logger (log object) – logging instance for output

  • strict (bool) – should we apply strict criterion? (default: False)

  • maxwin (int) – largest possible window in pixels (default: 100)

  • frac_max (float) – fraction of maximum flux to encompass (default: 0.5)

Returns:

  • x0 (int): lower limit of window in pixels

  • x1 (int): upper limit of window in pixels

  • count (int): number of pixels in window

class kcwidrp.primitives.MakeCube.MakeCube(action, context)[source]

Transform 2D images to 3D data cubes.

Uses thread pool to warp each slice in parallel and then assemble a 3D image cube from the slices. Rotates RED channel data by 180 degress to align with BLUE channel data.

Creates WCS keywords that reflect the new geometry of the cube based on position keywords from the telescope.

Outputs FITS 3d cube images:

  • icube - main cube with primary, uncertainty, mask, flag, and noskysub extensions

  • ocube - for nod-and-shuffle, the object frame as a cube

  • scube - for nod-and-shuffle, the sky frame as a cube

  • dcube - the delta-wavelength cube

kcwidrp.primitives.MakeCube.make_cube_helper(argument)[source]

Warp each slice.

Helper program for threaded warping of all slices.

Parameters:

argument (dict) – Dictionary of params for warping an individual slice.

Returns:

  • int: Slice number being warped

  • ndimage: warped intensity image

  • ndimage: warped uncertainty image

  • ndimage: warped mask image

  • ndimage: warped flag image

  • ndimage: warped image without sky subtraction

  • ndimage: warped nod-and-shuffle object image

  • ndimage: warped nod-and-shuffle sky image

  • ndimage: warped delta-wavelength image

class kcwidrp.primitives.MakeInvsens.MakeInvsens(action, context)[source]

Generate inverse sensitivity curve from a standard star observation.

Uses object name to determine if the observation is a standard star. Then checks that the image has been processed through DAR correction. Reads in the reference spectrum and compares it with the observed spectrum to generate the inverse sensitivity curve that can be used to flux calibrate science observations.

Outputs FITS inverse sensitivity spectrum along with diagnostic plots that include residuals between reference spectrum and calibrated observed spectrum and effective area and efficiency plots as a function of wavelength.

class kcwidrp.primitives.MakeMasterArc.MakeMasterArc(action, context)[source]

Stack arc frames into master arc

Generate a master arc frame based on the instrument config parameter arc_min_nframes, which defaults to 1 for the BLUE channel and 3 for the RED channel. It is assumed that each frame is well-exposed and the combine method ‘median’ will be used to mitigate cosmic rays (especially for the RED channel). A high sigma clipping of 2.0 is used to help with the CRs.

Uses the ccdproc.combine routine to peform the stacking.

Writes out a *_marc.fits file and records a master arc frame in the proc table, no matter how many frames are combined.

class kcwidrp.primitives.MakeMasterBias.MakeMasterBias(action, context)[source]

Stack bias frames into a master bias frame.

Generate a master bias image from overscan-subtracted and trimmed bias frames (*_intb.fits) based on the instrument config parameter bias_min_nframes, which defaults to 7. The combine method for biases is ‘average’ and so cosmic rays may be present, especially in RED channel data. A high sigma clipping of 2.0 is used to help with the CRs.

Uses the ccdproc.combine routine to perform the stacking.

Writes out a *_mbias.fits file and records a master bias frame in the proc table.

class kcwidrp.primitives.MakeMasterContbars.MakeMasterContbars(action, context)[source]

Stack continuum bars frames into master continuum bars

Generate a master cont bars frame based on the instrument config parameter contbars_min_nframes, which defaults to 1 for the BLUE channel and 3 for the RED channel. It is assumed that each frame is well-exposed and the combine method ‘median’ will be used to mitigate cosmic rays (especially for the RED channel). A high sigma clipping of 2.0 is used to help with the CRs.

Uses the ccdproc.combine routine to peform the stacking.

Writes out a *_mcbars.fits file and records a master cont bars frame in the proc table, no matter how many frames are combined.

class kcwidrp.primitives.MakeMasterDark.MakeMasterDark(action, context)[source]

Stack dark frames into master dark

Generate a master dark frame based on the instrument config parameter dark_min_nframes, which defaults to 3. The combine method ‘average’ will be used. A high sigma clipping of 2.0 is used to help with the CRs.

Uses the ccdproc.combine routine to peform the stacking.

Writes out a *_mdark.fits file and records a master dark frame in the proc table, no matter how many frames are combined.

class kcwidrp.primitives.MakeMasterFlat.MakeMasterFlat(action, context)[source]

Generate illumination correction from a stacked flat image.

Uses b-spline fits along with geometry maps to generate a master image for illumination correction. If the flat is internal, accounts for vignetting along one edge. Also accounts for ledge seen in BM grating.

Depending on the type of input stack, the following files are written out and entries are made in the proc file:

  • SFLAT - a *_mflat.fits file and an MFLAT entry

  • SDOME - a *_mdome.fits file and an MDOME entry

  • STWIF - a *_mtwif.fits file and an MTWIF entry

class kcwidrp.primitives.MakeMasterObject.MakeMasterObject(action, context)[source]

Stack object frames into master object

Generate a master object frame based on the instrument config parameter object_min_nframes, which defaults to 1 for the BLUE channel and 1 for the RED channel. The combine method ‘median’ will be used to mitigate cosmic rays (especially for the RED channel) if object_min_nframes is 3 or less. If larger than 3, then the method ‘average’ will be used. A high sigma clipping of 2.0 is used to help with the CRs.

Uses the ccdproc.combine routine to peform the stacking.

Writes out a *_mobj.fits file and records a master object frame in the proc table, no matter how many frames are combined.

class kcwidrp.primitives.MakeMasterSky.MakeMasterSky(action, context)[source]

Make master sky image.

Uses b-spline fits along with geometry maps to generate a master sky image for sky subtraction.

This routine also handles the file kcwi.sky, which controls the master sky generation. This file consists of one line per image, with the first column indicating the raw object image to be sky-subtracted. The following columns can either indicate a separate image to use for sky subtraction, the filename of a mask fits image for masking object flux, or indicate that the object is a continuum source and either automatically find the object, or specify the location and width of the continuum source. Below are example one-line entries and what they mean:

  1. Skip sky subtraction for this particular object image:

    • kr230925_00075.fits skip

2. Point to a different image for the sky (this assumes the *_sky.fits image has already been generated previously:

  • kr230925_00075.fits kr230925_00076.fits

3. Indicate that a mask file should be used to mask object flux when deriving the sky model (see kcwi_masksky_ds9.py):

  • kr230925_00075.fits kr230925_00075_smsk.fits

4. Indicate that this is a bright continuum source and automatically mask the continuum source from the sky model:

  • kr230925_00075.fits cont

5. Indicate that this is a faint continuum source and specify the location and the width of the continuum source (in pixels):

  • kr230925_00075.fits cont 45.0 7.6

If no kcwi.sky file exists, or there is no entry for the input object frame, then the entire image is used to generate the sky model.

It is good practice to run all the data through first, then inspect the sky subtraction and see which frames will benefit from masking or from a dedicated sky observation.

If a sky model is generated, the routine will write out a *_sky.fits image and add a sky entry in the proc table.

class kcwidrp.primitives.NandshuffSubtractSky.NandshuffSubtractSky(action, context)[source]

Locate object and sky panels and perform a nod-and-shuffle sky subtraction.

With a nod-and-shuffle observation, the resulting image will have a sky and an object panel with interleaved exposures of the same total duration. This routine identifies those panels and does a subtraction that should eliminate the sky accurately, as long as the interleaving occurs at a rate that tracks the changes in sky throughout the exposure.

Writes out a *_obj.fits and *_sky.fits image with the panels aligned, but not subtracted. These are subsequently processed through the pipeline to allow the confirmation that resulting features are not the result of sky subtraction.

Also writes out the sky-subtracted image in a *_intk.fits image and adds an entry in the proc table.

class kcwidrp.primitives.ProcessArc.ProcessArc(action, context)[source]

Preliminary processing of arc lamp images.

Checks if the default lamp is being read in. If not, terminates processing.

class kcwidrp.primitives.ProcessBias.ProcessBias(action, context)[source]

Preliminary processing of bias images.

class kcwidrp.primitives.ProcessContbars.ProcessContbars(action, context)[source]

Preliminary processing of continuum bars images.

class kcwidrp.primitives.ProcessDark.ProcessDark(action, context)[source]

Preliminary processing of dark images.

class kcwidrp.primitives.ProcessFlat.ProcessFlat(action, context)[source]

Preliminary processing of flat images.

class kcwidrp.primitives.ProcessObject.ProcessObject(action, context)[source]

Preliminary processing of object images.

class kcwidrp.primitives.ReadAtlas.ReadAtlas(action, context)[source]

Read in atlas spectrum and derive alignment offset.

Reads in spectral atlas spectrum of observed arc lamp. Uses instrument configuration to determine convoution kernel that will match the observed spectrum resolution and applies it. Then picks a central part of the observed spectrum and tapers the ends and does a cross-correlation to determine the pixel offset that will align the observed and atlas spectra.

Uses the following configuration parameters:

  • MIDFRAC: what central fraction to use for finding the atlas offset - defaults to 0.3.

  • TAPERFRAC: how much of the ends to taper off (in case of bright lines at the ends) - defaults to 0.2.

  • ATOFF: force a specific pixel offset for the atlas - defaults to 0.

Places convolved atlas spectrum with wavelength scale in the return argument so subsequent processing can use it to fit the observed wavelength scale accurately.

class kcwidrp.primitives.RectifyImage.RectifyImage(action, context)[source]

Ensure output image has a consistent orientation.

For the BLUE channel, identifies the CCD amplifier configuration and applies the appropriate geometric transformation to produce a consistent orientation independant of amp configuration. The RED channel is already rectified, so no geometric transformation is required.

Writes out a *_int.fits image regardless of which channel is being processed and adds an entry in the proc table.

class kcwidrp.primitives.RemoveCosmicRays.RemoveCosmicRays(action, context)[source]

Remove cosmic rays and generate a flag image recording their location.

Uses astroscrappy to detect and flag cosmic rays. Updates the following FITS header keywords:

  • CRCLEAN: set to True if operation performed.

  • NCRCLEAN: set to the number of cosmic ray pixels cleaned.

Uses the following configuration parameters:

  • saveintims: if set to True write out a CR cleaned version of image in *_crr.fits. Defaults to False.

  • CRR_MINEXPTIME - exposure time below which no CR cleaning is done.

  • CRR_* - see kcwi.cfg file for parameters controlling CR cleaning.

Updates image in returned arguments with CR cleaned image and adds flags indicating where changes were made.

class kcwidrp.primitives.SolveArcs.SolveArcs(action, context)[source]

Solve individual bar arc spectra for wavelength.

For each bar, identifies the atlas lines found in GetAtlasLines.py in the observed spectrum and use the pixel positions along with the atlas wavelengths to determine the wavelength solution for the bar.

Uses the following configuration parameters:

  • FRACMAX: fraction of line maximum to use for fitting peak. Defaults to 0.5 for BLUE and 0.25 for RED.

  • LINETHRESH: the threshhold intensity for finding observed lines. Defaults to 100 for BLUE and 10 for RED (can also be set on the command line).

Outputs diagnostic plots of the fitting and stores the coefficients for later use.

class kcwidrp.primitives.SolveGeom.SolveGeom(action, context)[source]

Solve the overall geometry of the IFU.

Takes individual bar arc spectra and generates a fit for each slice, which contains five bar arc spectra each. Given that there are only five points in the ‘x’ or spatial direction and many more in the ‘y’ or wavelength direction, an asymmetric polynomial is fit with order 2 in the x and order 4 in the y directions. The wavelength coverage of the observation is recorded in parameters with the range that includes all data being in the waveall0 and waveall1 parameters, and the range that includes only good data in wavegood0 and wavegood1 parameters. The middle of the wavelength range is recorded in the wavemid parameter.

Forward and inverse transforms, along with all the parameters for the geometric fit are written out as a python pickled dictionary in a *_geom.pkl file.

class kcwidrp.primitives.StackFlats.StackFlats(action, context)[source]

Stack flat images

Generate a stacked flat frame based on one of the following instrument config parameters:

  • flat_min_nframes: defaults to 6

  • dome_min_nframes: defaults to 3

  • twiflat_min_nframes: defaults to 1

It is assumed that each frame is well-exposed and the combine method ‘average’ will be used. A high sigma clipping of 2.0 is used to help with the CRs.

Uses the ccdproc.combine routine to peform the stacking.

Writes out the following files and adds entries in the proc table depending on the input IMTYPE:

  • FLATLAMP - writes out *_sflat.fits and puts SFLAT in proc table

  • DOMEFLAT - writes out *_sdome.fits and puts SDOME in proc table

  • TWIFLAT - writes out *_stwif.fits and puts STWIF in proc table

class kcwidrp.primitives.StartBokeh.StartBokeh(action, context)[source]

Start the bokeh server for the KCWI DRP.

This is enabled through the enable_bokeh instrument configuration parameter in the kcwi.cfg file. Set to True to enable plots.

class kcwidrp.primitives.SubtractBias.SubtractBias(action, context)[source]

Subtract the master bias frame.

Reads in the master bias created by MakeMasterBias.py and performs the subtraction (after verifying amplifier configuration agreement). Records the processing in the header.

class kcwidrp.primitives.SubtractDark.SubtractDark(action, context)[source]

Subtract the master dark frame.

Checks for existence of master dark frame and, if it exists, performs subtraction and records processing in the header.

class kcwidrp.primitives.SubtractOverscan.SubtractOverscan(action, context)[source]

Determines overscan offset and subtracts it from the image.

Uses the BIASSEC header keyword to determine where to calculate the overscan offset. Subtracts the overscan offset and records the value in the header. In addition, performs a polynomial fit and uses the residuals to determine the read noise in the overscan. Records the overscan readnoise in the header as well.

class kcwidrp.primitives.SubtractScatteredLight.SubtractScatteredLight(action, context)[source]

Subtract scattered light between slices.

Uses the centeral, un-illuminated part of the image to generate a model of the scattered light. Subtract this model from the entire image.

Uses the following configuration parameter:

  • skipscat: set to True to skip scattered light subtraction. Defaults to False.

Writes out a *_intd.fits file regardless if scattered light is subtracted or not.

class kcwidrp.primitives.SubtractSky.SubtractSky(action, context)[source]

Subtract sky model from observation.

Reads in model of sky created in MakeMasterSky.py and subtracts it from observation. Which model is determined by the kcwi.sky file, if present. If not, then the model generated from the observation itself will be used.

Which master sky was used will be recorded in the header keyword SKYMAST. If a mask file was used to generate the sky model, it will be recorded in SKYMSKF. If scaling was required prior to subtraction due to differing exposure times of observation and model, the scale value is recorded in SKYSCL.

Writes out a *_intk.fits file and adds an entry to the proc file.

class kcwidrp.primitives.Template.Template(action, context)[source]

Generic template for primitive routines.

class kcwidrp.primitives.TraceBars.TraceBars(action, context)[source]

Derive bar traces

Using the output from FindBars.py, trace each bar to the top and bottom of the image.

Uses the following configuration parameter:

  • saveintims: if set to True write out a warped version of bars image in *_warped.fits. Defaults to False.

Writes out resulting spatial control points as a FITS table in *_trace.fits.

class kcwidrp.primitives.TrimOverscan.TrimOverscan(action, context)[source]

Trim overscan region from image.

Uses the data section (DSECn header keyword) to determine how to trim the image to exclude the overscan region.

Removes raw section keywords, ASECn, BSECn, CSECn, and DSECn after trimming and replaces them with the ATSECn keyword giving the image section in the trimmed image for each amplifier.

Uses the following configuration parameter:

  • saveintims: if set to True write out the trimmed image as *_trim.fits. Default is False.

If the input image is a bias frames, writes out a *_intb.fits file, otherwise, just updates the image in the returned arguments.

class kcwidrp.primitives.WavelengthCorrections.WavelengthCorrections(action, context)[source]

Perform wavelength corrections

Based on the type specified in the kcwidrp/configs/kcwi.cfg file in the radial_velocity_correction parameter, perform the requested correction. The current options are “heliocentric”, “barycentric”, and “none”. The default is “heliocentric”.

Also, if the air_to_vacuum parameter is set to True (default), apply air-to-vacuum correction.

a2v_conversion(wave)[source]

Convert air-based wavelengths to vacuum

Adapted from wave.py in: https://github.com/pypeit/PypeIt/ Formula from https://ui.adsabs.harvard.edu/abs/1996ApOpt..35.1566C/

Parameters:

wave (ndarray) – Wavelengths

Returns:

Wavelength array corrected to vacuum wavelengths

air2vac(obj, mask=False)[source]

Convert wavelengths in a cube from standard air to vacuum.

Parameters:
  • obj (astropy HDU / HDUList) – Input HDU/HDUList with 3D data.

  • mask (bool) – Set if the cube is a mask cube.

Returns:

HDU / HDUList: Trimmed FITS object with updated header. Return type matches type of obj argument.

get_wav_axis(header)[source]

Returns a NumPy array representing the wavelength axis of a cube.

Adapted from https://github.com/dbosul/cwitools.git

Parameters:

header (astropy.io.fits.Header) – header that contains wavelength or velocity axis that is specified in ‘CTYPE’ keywords in any dimension.

Returns:

Wavelength axis for this data.

Return type:

numpy.ndarray

heliocentric(obj, correction_mode, mask=False, resample=True, vcorr=None)[source]

Apply heliocentric correction to the cubes.

Note that this only works for KCWI data because the location of Keck Observatory is hard-coded in the function.

Adapted from https://github.com/dbosul/cwitools.git

Parameters:
  • obj (astropy HDU / HDUList) – Input HDU/HDUList with 3D data.

  • correction_mode (str) – “none”, “barycentric”, or “heliocentric”

  • mask (bool) – Set if the cube is a mask cube. This only works for resampled cubes.

  • resample (bool) – Resample the cube to the original wavelength grid?

  • vcorr (float) – Use a different correction velocity.

Returns:

Trimmed FITS object with updated header. vcorr (float): (if vcorr is True) Correction velocity in km/s. Return type matches type of fits_in argument.

Return type:

HDU / HDUList

Examples

To apply heliocentric correction,

>>> hdu_new = heliocentric(hdu_old)

However, this resamples the wavelengths back to the original grid. To use the new grid without resampling the data,

>>> hdu_new = heliocentric(hdu_old, resample=False)
locate_object_file(suffix)[source]

Return FITS HDU list if current file with requested suffix can be found.

Will return None if file cannot be found.

Parameters:

suffix (str) – The file suffix that you want to operate on.

Returns:

FITS HDU List from kcwi_fits_reader routine, or None if file with suffix not found.

class kcwidrp.primitives.kcwi_file_primitives.KCCDData(*args, **kwd)[source]

A container for KCWI images based on the CCDData object that adds the noskysub frame to allow parallel processing of sky-subtracted and un-sky-subtracted KCWI images.

noskysub

Optional un-sky-subtracted frame that is created at the sky subtraction stage and processed in parallel with the primary frame. Default is None.

Type:

numpy.ndarray or None

kcwidrp.primitives.kcwi_file_primitives.fix_header(ccddata)[source]

Fix header keywords for DRP use.

Update GAINn keywords for Blue channel.

Add FITS header keywords to Red channel data to make compatible with DRP.

Adds the following keywords that are not present in raw images:

  • MJD - Modified Julian Day (only for AIT data)

  • NVIDINP - from TAPLINES keyword

  • GAINMUL - set to 1

  • CCDMODE - from CDSSPEED keyword

  • AMPNUM - based on AMPMODE

  • GAINn - from red_amp_gain dictionary

kcwidrp.primitives.kcwi_file_primitives.get_master_name(tab, target_type, loc=0)[source]

Add a specific tag to an output fits filename read from a proc table.

Parameters:
  • tab (proc table) – proc table source of filename.

  • target_type (str) – suffix to add after underscore.

  • loc (int) – row within table in tab, defaults to 0.

Returns:

constructed filename from input tab entry and target_type.

Return type:

(str)

class kcwidrp.primitives.kcwi_file_primitives.ingest_file(action, context)[source]

File ingestion class for KCWI images.

Parameters:
  • action (str) – Pipeline action

  • context – Pipeline context which includes arguments for the action.

logger

KCWI pipeline logger

adjang()[source]

Return the adjustment angle in degrees for the current grating.

Returns:

180.0 for high-resolution gratings, otherwise 0.0

Return type:

(float)

ampmode()[source]

Returns value of AMPMODE FITS header keyword.

apply()[source]

Apply method for class.

Checks _pre_condition(). If True, then call _perform() and collect output. Then check _post_condition().

Returns:

output from _perform, or None if there was an exception.

atsig()[source]

Return the Gaussian sigma to use for convolving the Atlas spectrum.

Uses grating() and ifunum() to select a Gaussian sigma that when convolved with the Atlas spectrum will match the resolution of the observed spectrum.

Raises:

ValueError – if grating or ifunum canoot be determined.

Returns:

sigma in pixels ranging from 0.75 to 14.0, or 0.0 if imtype() is BIAS.

Return type:

(float)

calibration_lamp()[source]

Determine which calibration source was used for a given frame.

Examines LMPnSTAT and LMPnSHST keywords to determine which lamp was on and which shutter was open and thus providing illuminate for the frame. Returns None if not an ARCLAMP image type.

Returns:

Which calibration lamp was active for frame.

Return type:

(str, or None)

camang()[source]

Return the articulation stage camera angle.

Determines which channel to use based on camera() call, then uses keyword BARTANG for Blue channel and RARTANG for Red.

Raises:

ValueError – if the camera ID is Unknown.

Returns:

Camera angle in degrees for the given channel.

Return type:

(float)

camera()[source]

Get camera ID number.

Uses CAMERA FITS header keyword.

Returns:

0 for Blue channel, 1 for Red, and -1 for Unknown.

check_if_file_can_be_processed(imtype)[source]

For a given image type, ensure that processing can proceed.

Based on IMTYPE keyword, makes a call to proctab to see if pre-requisite images are present.

Returns:

True if processing can proceed, False if not.

Return type:

(bool)

cwave()[source]

Return the central wavelength in Angstroms.

Determines which channel to use based on camera() call, then uses keyword BCWAVE for Blue channel and RCWAVE for Red.

Raises:

ValueError – if the camera ID is Unknown.

Returns:

Central wavelength in Angstroms for the given channel.

Return type:

(float)

delta_wave_out()[source]

Return output delta lambda in Angstroms for the given grating

Calls grating(), and ybinsize() to calculate output delta lambda.

Raises:

ValueError – if the grating is unknown.

Returns:

delta lambda in Angstroms, or 0. if imtype() is BIAS.

Return type:

(float)

dich()[source]

Query if dichroic present for image.

Returns True for all Red channel data, but checks for the presence of Red keywords in the Blue channel header to determine status of dichroic.

Returns:

True if dichroic present, False if not.

Return type:

(bool)

filter()[source]

Return the filter name.

Determines which channel to use based on camera() call, then uses keyword BFILTNAM for Blue channel and returns None for Red.

Raises:

ValueError – if the camera ID is Unknown.

Returns:

Filter name for the given channel.

Return type:

(str)

get_keyword(keyword)[source]

Get a keyword from ingested headers.

Parameters:

keyword (str) – Keyword to fetch from header

Returns:

Keyword value if present, otherwise None.

grangle()[source]

Return the grating angle.

Determines which channel to use based on camera() call, then uses keyword BGRANGLE for Blue channel and RGRANGLE for Red.

Raises:

ValueError – if the camera ID is Unknown.

Returns:

Grating angle in degrees for the given channel.

Return type:

(float)

grating()[source]

Return the grating name.

Determines which channel to use based on camera() call, then uses keyword BGRATNAM for Blue channel and RGRATNAM for Red.

Raises:

ValueError – if the camera ID is Unknown.

Returns:

Grating name for the given channel.

Return type:

(str)

ifuname()[source]

Return the value of the IFUNAM FITS header keyword.

ifunum()[source]

Return the value of the IFUNUM FITS header keyword.

illum()[source]

Generate a string that characterizes the illumination for the frame.

Uses various FITS header keywords to determine the kind of illumination that was used for the frame. If a consistent picture of the illumination cannot be determined, set the return value to Test.

Returns:

Characterization of illumination for given frame.

Return type:

(str)

imtype()[source]

Return the value of the IMTYPE FITS header keyword.

map_ccd(xbin, ybin)[source]

Return CCD section variables useful for processing

Parameters:
  • xbin (int) – binning in x

  • ybin (int) – binning in y

Uses FITS keyword NVIDINP to determine how many amplifiers were used to read out the CCD. Then reads the corresponding BSECn, and DSECn keywords, where n is the amplifier number. The indices are converted to Python (0-biased, y axis first) indices and an array is constructed for each of the two useful sections of the CCD as follows:

  • Bsec[0][0] - First amp, y lower limit

  • Bsec[0][1] - First amp, y upper limit

  • Bsec[0][2] - First amp, x lower limit

  • Bsec[0][3] - First amp, x upper limit

  • Bsec[1][0] - Second amp, y lower limit

  • etc.

Bsec is the full overscan region for the given amplifier and is used to calculate and perform the overscan subtraction.

Dsec is the full CCD region for the given amplifier and is used to trim the image after overscan subtraction has been performed.

Tsec accounts for trimming the image according to Dsec.

Amps are assumed to be organized as follows:

          BLUE                          RED
(0,ny)  --------- (nx,ny)    (0,ny)  --------- (nx,ny)
        | 2 | 3 |                    | 0 | 2 |
        ---------                    ---------
        | 0 | 1 |                    | 1 | 3 |
(0,0)   --------- (nx, 0)    (0,0)   --------- (nx, 0)
Returns:

  • list: (int) y0, y1, x0, x1 for bias section

  • list: (int) y0, y1, x0, x1 for data section

  • list: (int) y0, y1, x0, x1 for trimmed section

  • list: (bool) y-direction, x-direction, True if forward, else False

namps()[source]

Return value of NVIDINP FITS header keyword.

nasmask()[source]

Query if mask was inserted for image.

Calls camera() to determine channel, and then uses BNASNAM for Blue channel or RNASNAM for Red channel.

Raises:

ValueError – if channel is Unknown.

Returns:

True if ‘Mask’ in, False if not.

Return type:

(bool)

numopen()[source]

Returns value of NUMOPEN FITS header keyword.

plotlabel()[source]

Return label string for plot titles as (str).

resolution()[source]

Return FWHM resolution in Angstroms for the given grating.

Calls cwave(), grating(), and ifunum() to calculate resolution.

Raises:

ValueError – if the grating is unknown.

Returns:

FWHM in Angstroms, or 0. if imtype() is BIAS.

Return type:

(float)

rho()[source]

Return the rho value (lines/mm divided by 1000) for the given grating.

Uses grating() to determine ingested grating.

Returns:

rho value or None if grating unknown.

Return type:

(float, or None)

shufrows()[source]

Returns value of SHUFROWS FITS header keyword.

stdlabel()[source]

Return label string for standard star plot titles as (str).

xbinsize()[source]

Return X part of BINNING keyword value as (int).

ybinsize()[source]

Return Y part of BINNING keyword value as (int).

kcwidrp.primitives.kcwi_file_primitives.kcwi_fits_reader(file)[source]

A reader for KCCDData objects.

Currently, this is a separate function, but should probably be registered as a reader similar to fits_ccddata_reader.

Parameters:

file (str) – The filename (or pathlib.Path) of the FITS file to open.

Raises:
  • FileNotFoundError – if file not found or

  • OSError – if file not accessible

Returns:

All relevant frames in a single KCCDData object and a FITS table of exposure events, if present otherwise None.

Return type:

(KCCDData, FITS table)

kcwidrp.primitives.kcwi_file_primitives.kcwi_fits_writer(ccddata, table=None, output_file=None, output_dir=None, suffix=None)[source]

A writer for KCCDData or CCDData objects.

Updates history in FITS header with pipeline version and git repo version and date.

Converts float64 data to float32.

Uses object to_hdu() method to generate hdu list and then checks if various extra frames are present (flags, noskysub) and adds them to the hdu list prior to writing out with hdu list writeto() method.

Note

Currently fits tables are not written out.

Parameters:
  • ccddata (KCCDData or CCDData) – object to write out.

  • table (FITS Table) – currently not used.

  • output_file (str) – base filename to write to.

  • output_dir (str) – directory into which to write.

  • suffix (str) – a suffix to append to output_file string.

kcwidrp.primitives.kcwi_file_primitives.parse_imsec(section=None)[source]

Parse image section FITS header keyword into useful tuples.

Take into account one-biased IRAF-style image section keywords and the possibility that a third element (strid) may be present and generate the x and y limits as well as the stride for each axis that are useful for python image slices.

Parameters:
  • section (str) – square-bracket enclosed string with colon range

  • delimiters. (specifiers and comma) –

Returns:

  • list: (int) y0, y1, x0, x1 - zero-biased (python) slice limits.

  • list: (int) y-stride, x-stride - strides for each axis.

kcwidrp.primitives.kcwi_file_primitives.plotlabel(args)[source]

Return label string for plot titles as (str).

kcwidrp.primitives.kcwi_file_primitives.read_table(input_dir=None, file_name=None)[source]

Read FITS table

Uses astropy.Table module to read in FITS table.

Raises:

FileNotFoundError – if file not found.

Returns:

table read in or None if unsuccessful.

Return type:

(FITS Table)

kcwidrp.primitives.kcwi_file_primitives.strip_fname(filename)[source]

Return pathlib.Path.stem attribute for given filename.

Parameters:

filename (str) – filename to strip.

Returns:

Path(filename).stem, or filename if error occurs.

Return type:

(str)

kcwidrp.primitives.kcwi_file_primitives.write_table(output_dir=None, table=None, names=None, comment=None, keywords=None, output_name=None, clobber=False)[source]

Write out FITS table.

Parameters:
  • output_dir (str) – output directory for table.

  • table (list of arrays) – each array in list should have the same size.

  • names (list of str) – one for each column in table.

  • comment (string) – text for FITS COMMENT header record.

  • keywords (FITS keyword dict) – for FITS header records.

  • output_name (str) – name for output table.

  • clobber (bool) – set to to True to overwrite existing table.

Raises:
  • FileExistsError – if cannot overwrite file or

  • OSError – if some other access exception occurred.

kcwidrp.scripts.kcwi_masksky_ds9.main()[source]

Creates mask image from ds9 region file.

To use this routine, process your data with default sky subtraction. Then display the target *_intf.fits file in ds9. Use region shapes to indicate non-sky pixels in image (box, circle, etc.). Write out ds9 region file (*.reg). Then run this routine:

  • python ~/kderp/devel/kcwi_masksky_ds9.py kb180101_00111_intf.fits ds9.reg

(replace paths/filenames with your local paths/filenames)

This should create kb180101_00111_smsk.fits, which will be used when you re-run the pipeline.

Parameters:
  • imagename (string) – The name of a *_intf.fits image

  • regionname (string) – The name of a ds9 region file

Returns:

None

kcwidrp.scripts.wb.wb_main()[source]

Generate summary log and group list files for BLUE channel images.

Call get_log_string to create summary entry for each image and group them according to processing group. Write out unique processing group lists in *.txt files. These files can be input to the pipeline with the -l command line parameter to allow processing of groups one at a time. For example, 2x2 Blue biases taken with the TUP amp configuration in slow readout with gainmul 10 will end up in the file bias2x2TUP010_0.txt. A master bias can be created by issuing the following command:

>>> reduce_kcwi -b -l bias2x2TUP010_0.txt

These group files are generated for biases, darks, continuum bars, arcs, flats, and all objects. The filenames are all appended with the last four characters of the STATEID header keyword, so identical configurations from different states can be distinguished.

Always good to type out the list file before processing it.

Examples

>>> wb kb*.fits > whatb.list

This will generate a summary log file along with associated group list files that can be used as inputs to the reduce_kcwi command with the -l parameter. An example of the resulting *.txt files is below:

SN2023ixf2x2MedKBlueBL4500_75fe.txt      bias2x2TUP010_0.txt
allb.txt                                 cbars2x2MedKBlueBL_4500_0.7_75fe.txt
arcs2x2MedKBlueBLFeAr4500_10.0_75fe.txt  cflat2x2MedKBlueBL_4500_0.7_75fe.txt
arcs2x2MedKBlueBLThAr4500_20.0_75fe.txt  dflat2x2MedKBlueBL_4500_14.0_75fe.txt
bd26d26062x2MedKBlueBL4500_75fe.txt

One can proceed through processing steps like this:

>>> reduce_kcwi -b -l bias2x2TUP010_0.txt
>>> reduce_kcwi -b -l cbars2x2MedKBlueBL_4500_0.7_75fe.txt
>>> reduce_kcwi -b -l arcs2x2MedKBlueBLThAr4500_20.0_75fe.txt
>>> reduce_kcwi -b -l cflat2x2MedKBlueBL_4500_0.7_75fe.txt
>>> reduce_kcwi -b -l bd26d26062x2MedKBlueBL4500_75fe.txt
>>> reduce_kcwi -b -l SN2023ixf2x2MedKBlueBL4500_75fe.txt
kcwidrp.scripts.wr.wr_main()[source]

Generate summary log and group list files for RED channel images.

Call get_log_string to create summary entry for each image and group them according to processing group. Write out unique processing group lists in *.txt files. These files can be input to the pipeline with the -l command line parameter to allow processing of groups one at a time. For example, 2x2 Red biases taken with the L2U2 amp configuration with slow readout and high gain will end up in the file bias2x2TUP01_0.txt. A master bias can be created by issuing the following command:

>>> reduce_kcwi -r -l bias2x2L2U201_0.txt

These group files are generated for biases, darks, continuum bars, arcs, flats, and all objects. The filenames are all appended with the last four characters of the STATEID header keyword, so identical configurations from different states can be distinguished.

Always good to type out the list file before processing it.

Examples

>>> wr kr*.fits > whatr.list

This will generate a summary log file along with associated group list files that can be used as inputs to the reduce_kcwi command with the -l parameter. An example of the resulting *.txt files is below:

SN2023ixf2x2MedRL8000_75fe.txt     bias2x2L2U201_0.txt
allr.txt                           cbars2x2MedRL_8000_5.0_75fe.txt
arcs2x2MedRLFeAr8000_2.5_75fe.txt  cflat2x2MedRL_8000_5.0_75fe.txt
arcs2x2MedRLThAr8000_2.5_75fe.txt  dflat2x2MedRL_8000_20.0_75fe.txt
bd26d26062x2MedRL8000_75fe.txt

One can proceed through processing steps like this:

>>> reduce_kcwi -r -l bias2x2L2U201_0.txt
>>> reduce_kcwi -r -l cbars2x2MedRL_8000_5.0_75fe.txt
>>> reduce_kcwi -r -l arcs2x2MedRlThAr8000_2.5_75fe.txt
>>> reduce_kcwi -r -l cflat2x2MedRL_8000_5.0_75fe.txt
>>> reduce_kcwi -r -l bd26d26062x2MedRL8000_75fe.txt
>>> reduce_kcwi -r -l SN2023ixf2x2MedRL8000_75fe.txt