pyreduce package¶

Subpackages¶

Submodules¶

pyreduce.pipeline module¶

Fluent Pipeline API for PyReduce.

Provides a cleaner interface for building and running reduction pipelines. Wraps the existing Step classes internally for backward compatibility.

Example usage:

from pyreduce.pipeline import Pipeline

# Simple: auto-discover files for an instrument result = Pipeline.from_instrument(

instrument=”UVES”, target=”HD132205”, night=”2010-04-01”, channel=”middle”, base_dir=”/data”,

).run()

# Or build manually with explicit files: result = (

Pipeline(“UVES”, output_dir, config=settings) .bias(bias_files) .flat(flat_files) .trace() .extract(science_files) .run()

)

# For multi-fiber instruments with separate illumination files: t1, cr1 = pipe.trace_raw([even_flat]) t2, cr2 = pipe.trace_raw([odd_flat]) pipe.organize(t1, cr1, t2, cr2) pipe.extract([science_file]).run()

class pyreduce.pipeline.Pipeline(instrument: Instrument | str, output_dir: str, target: str = '', channel: str = '', night: str = '', config: dict | None = None, trace_range: tuple[int, int] | None = None, plot: int = 0, plot_dir: str | None = None)[source]¶

Bases: object

Fluent API for building reduction pipelines.

STEP_CLASSES = {'bias': <class 'pyreduce.reduce.Bias'>, 'continuum': <class 'pyreduce.reduce.ContinuumNormalization'>, 'curvature': <class 'pyreduce.reduce.SlitCurvatureDetermination'>, 'finalize': <class 'pyreduce.reduce.Finalize'>, 'flat': <class 'pyreduce.reduce.Flat'>, 'freq_comb': <class 'pyreduce.reduce.LaserFrequencyCombFinalize'>, 'freq_comb_master': <class 'pyreduce.reduce.LaserFrequencyCombMaster'>, 'mask': <class 'pyreduce.reduce.Mask'>, 'norm_flat': <class 'pyreduce.reduce.NormalizeFlatField'>, 'rectify': <class 'pyreduce.reduce.RectifyImage'>, 'scatter': <class 'pyreduce.reduce.BackgroundScatter'>, 'science': <class 'pyreduce.reduce.ScienceExtraction'>, 'trace': <class 'pyreduce.reduce.Trace'>, 'wavecal': <class 'pyreduce.reduce.WavelengthCalibrationFinalize'>, 'wavecal_init': <class 'pyreduce.reduce.WavelengthCalibrationInitialize'>, 'wavecal_master': <class 'pyreduce.reduce.WavelengthCalibrationMaster'>}¶

STEP_ORDER = {'bias': 10, 'continuum': 90, 'curvature': 40, 'finalize': 100, 'flat': 20, 'freq_comb': 72, 'freq_comb_master': 70, 'mask': 5, 'norm_flat': 50, 'rectify': 75, 'scatter': 45, 'science': 80, 'trace': 30, 'wavecal': 67, 'wavecal_init': 64, 'wavecal_master': 60}¶

bias(files: list[str]) → Pipeline[source]¶: Combine bias frames into master bias.

continuum() → Pipeline[source]¶: Normalize continuum.

curvature(files: list[str] | None = None) → Pipeline[source]¶: Determine slit curvature (p1/p2).

extract(files: list[str]) → Pipeline[source]¶: Extract science spectra.

finalize() → Pipeline[source]¶: Write final output files.

flat(files: list[str]) → Pipeline[source]¶: Combine flat frames into master flat.

freq_comb() → Pipeline[source]¶: Finalize frequency comb calibration.

freq_comb_master(files: list[str]) → Pipeline[source]¶: Extract laser frequency comb spectrum.

classmethod from_files(files: dict, output_dir: str, target: str, instrument, channel: str, night: str, config: dict, trace_range=None, steps='all', plot: int = 0, plot_dir: str | None = None) → Pipeline[source]¶

Create pipeline from a files dict and run specified steps.

This provides a simpler interface similar to the legacy Reducer class.

Parameters:

files (dict) – Files for each step (bias, flat, orders, wavecal, science, etc.)
output_dir (str) – Output directory
target (str) – Target name
instrument (Instrument or str) – Instrument instance or name
channel (str) – Instrument channel
night (str) – Observation night
config (dict) – Configuration dict
trace_range (tuple, optional) – Order range to process
steps (list or "all") – Steps to run
plot (int, optional) – Plot level (0=off, 1=basic, 2=detailed). Default 0.
plot_dir (str, optional) – Directory to save plots as PNG files. If None, plots are shown interactively.

Returns:

Configured pipeline ready to run

Return type:

Pipeline

Create pipeline from instrument name with automatic file discovery.

This is the recommended entry point for running reductions. It handles loading the instrument, finding and sorting files, and setting up the pipeline with the correct configuration.

Parameters:

instrument (str) – Instrument name (e.g., “UVES”, “HARPS”, “XSHOOTER”)
target (str) – Target name or regex pattern to match in headers
night (str, optional) – Observation night (YYYY-MM-DD format or regex)
channel (str, optional) – Instrument channel (e.g., “RED”, “BLUE”, “middle”). If None, uses all available channels for the instrument.
steps (tuple, list, or "all") – Steps to run. Default “all” runs all applicable steps.
base_dir (str, optional) – Base directory for data. Default: $REDUCE_DATA or ~/REDUCE_DATA
input_dir (str, optional) – Input directory relative to base_dir. Default: from config
output_dir (str, optional) – Output directory relative to base_dir. Default: from config
configuration (dict, optional) – Configuration overrides. Default: instrument defaults
trace_range (tuple, optional) – (first, last+1) orders to process
plot (int) – Plot level (0=off, 1=basic, 2=detailed)
plot_dir (str, optional) – Directory to save plots. If None, shows interactively.

Returns:

Configured pipeline ready to call .run()

Return type:

Pipeline

Example

>>> result = Pipeline.from_instrument(
...     instrument="UVES",
...     target="HD132205",
...     night="2010-04-01",
...     channel="middle",
...     steps=("bias", "flat", "trace", "science"),
... ).run()

load(step: str, data=None) → Pipeline[source]¶

Load intermediate result instead of computing.

Parameters:

step (str) – Name of step whose output to load
data (any, optional) – Data to use directly instead of loading from disk

mask() → Pipeline[source]¶: Load or create bad pixel mask.

normalize_flat() → Pipeline[source]¶: Normalize flat field, extract blaze function.

organize(traces: list, *more) → Pipeline[source]¶

Organize traces into fiber groups based on instrument config.

Use this after trace_raw() to apply fiber grouping configuration. Can accept multiple trace lists which will be concatenated.

Parameters:

traces (list[Trace]) – Trace objects from trace_raw()
*more (additional list[Trace]) – Optional additional trace lists to concatenate

Returns:

Self for method chaining

Return type:

Pipeline

Example

>>> t1 = pipe.trace_raw([even_flat])
>>> t2 = pipe.trace_raw([odd_flat])
>>> pipe.organize(t1, t2)
>>> pipe.extract([science_file]).run()

rectify() → Pipeline[source]¶: Rectify 2D image.

property results: dict¶: Access results after run().

run(skip_existing: bool = False) → dict[source]¶

Execute all queued steps.

Parameters:: skip_existing (bool) – If True, skip steps whose output files already exist
Returns:: Results keyed by step name
Return type:: dict

scatter(files: list[str] | None = None) → Pipeline[source]¶: Fit background scatter model.

trace(files: list[str] | None = None) → Pipeline[source]¶

Trace fibers/orders on flat field.

Parameters:: files (list[str], optional) – Files to use for tracing. If not provided, uses flat from previous step.
Returns:: Self for method chaining
Return type:: Pipeline

trace_raw(files: list[str], order_centers: dict = None) → list[source]¶

Trace fibers/orders and return Trace objects without storing.

Use this for multi-file tracing workflows where you need to combine traces from multiple files before grouping.

Parameters:

files (list[str]) – Files to use for tracing.
order_centers (dict[int, float], optional) – Order number -> y-position mapping for m assignment.

Returns:

Trace objects with fiber_idx set (individual fibers, not grouped).

Return type:

list[Trace]

wavecal() → Pipeline[source]¶: Finalize wavelength calibration.

wavecal_init() → Pipeline[source]¶: Initialize wavelength solution from line atlas.

wavecal_master(files: list[str]) → Pipeline[source]¶: Extract wavelength calibration spectrum.

wavelength_calibration(files: list[str]) → Pipeline[source]¶: Full wavelength calibration (master + init + finalize).

pyreduce.reduce module¶

REDUCE script for spectrograph data

Authors¶

Ansgar Wehrhahn (ansgar.wehrhahn@physics.uu.se) Thomas Marquart (thomas.marquart@physics.uu.se) Alexis Lavail (alexis.lavail@physics.uu.se) Nikolai Piskunov (nikolai.piskunov@physics.uu.se)

Version¶

1.0 - Initial PyReduce

License¶

…

class pyreduce.reduce.BackgroundScatter(*args, **config)[source]¶

Bases: CalibrationStep

Determine the background scatter

load()[source]¶

Load scatter results from disk

Returns:: scatter – scatter coefficients
Return type:: array

run(files, trace: list[Trace], mask=None, bias=None)[source]¶

Execute the current step

This should fail if files are missing or anything else goes wrong. If the user does not want to run this step, they should not specify it in steps.

Parameters:: files (list(str)) – data files required for this step
Raises:: NotImplementedError – needs to be implemented for each step

save(scatter)[source]¶

Save scatter results to disk

Parameters:: scatter (array) – scatter coefficients

property savefile¶

Name of the scatter file

Type:: str

scatter_degree¶

Polynomial degrees for the background scatter fit, in row, column direction

Type:: tuple(int, int)

class pyreduce.reduce.Bias(*args, **config)[source]¶

Bases: Step

Calculates the master bias

degree¶

polynomial degree of the fit between exposure time and pixel values

Type:: int

load(mask)[source]¶

Load the master bias from a previous run

Parameters:

mask (array of shape (nrow, ncol)) – Bad pixel mask

Returns:

bias (masked array of shape (nrow, ncol)) – master bias data, with the bad pixel mask applied
bhead (FITS header) – header of the master bias

run(files, mask=None)[source]¶

Calculate the master bias

Parameters:

files (list(str)) – bias files
mask (array of shape (nrow, ncol), optional) – bad pixel map

Returns:

bias (masked array of shape (nrow, ncol)) – master bias data, with the bad pixel mask applied
bhead (FITS header) – header of the master bias

save(bias, bhead)[source]¶

Save the master bias to a FITS file

Parameters:

bias (array of shape (nrow, ncol)) – bias data
bhead (FITS header) – bias header

property savefile¶

Name of master bias fits file

Type:: str

class pyreduce.reduce.CalibrationStep(*args, **config)[source]¶

Bases: Step

bias_scaling¶

how to adjust for diferences between the bias and flat field exposure times

Type:: {‘number_of_files’, ‘exposure_time’, ‘mean’, ‘median’, ‘none’}

calibrate(files, mask, bias=None, norm_flat=None, traces=None, extraction_height=None)[source]¶

norm_scaling¶

how to apply the normalized flat field

Type:: {‘divide’, ‘none’}

class pyreduce.reduce.ContinuumNormalization(*args, **config)[source]¶

Bases: Step

Determine the continuum to each observation

load(norm_flat, science)[source]¶

Load the results from the continuum normalization

Returns:

heads (list(FITS header)) – FITS headers of each observation
specs (list(array of shape (ntrace, ncol))) – extracted spectra
sigmas (list(array of shape (ntrace, ncol))) – uncertainties of the extracted spectra
conts (list(array of shape (ntrace, ncol))) – continuum for each spectrum
columns (list(array of shape (ntrace, 2))) – column ranges for each spectra

run(science, norm_flat, trace: list)[source]¶

Determine the continuum to each observation Also splices the orders together

Parameters:

science (tuple) – results from science step: (heads, list[list[Spectrum]])
norm_flat (tuple) – results from the normalized flatfield step
trace (list[TraceData]) – Trace objects with wavelength polynomials

Returns:

heads (list(FITS header)) – FITS headers of each observation
specs (list(array of shape (ntrace, ncol))) – extracted spectra
sigmas (list(array of shape (ntrace, ncol))) – uncertainties of the extracted spectra
conts (list(array of shape (ntrace, ncol))) – continuum for each spectrum
columns (list(array of shape (ntrace, 2))) – column ranges for each spectra

save(heads, specs, sigmas, conts, columns)[source]¶

Save the results from the continuum normalization

Parameters:

heads (list(FITS header)) – FITS headers of each observation
specs (list(array of shape (ntrace, ncol))) – extracted spectra
sigmas (list(array of shape (ntrace, ncol))) – uncertainties of the extracted spectra
conts (list(array of shape (ntrace, ncol))) – continuum for each spectrum
columns (list(array of shape (ntrace, 2))) – column ranges for each spectra

property savefile¶

savefile name

Type:: str

class pyreduce.reduce.ExtractionStep(*args, **config)[source]¶

Bases: Step

extract_to_arrays(img, head, trace_list: list[Trace], scatter=None)[source]¶: Extract spectra and return as arrays (for wavecal compatibility).

extraction_kwargs¶

arguments for the extraction

Type:: dict

extraction_method¶

Extraction method to use

Type:: {‘simple’, ‘optimal’}

class pyreduce.reduce.Finalize(*args, **config)[source]¶

Bases: Step

Create the final output files

output_file(number, name)[source]¶: str: output file name

run(continuum, trace: list, config)[source]¶

Create the final output files

this is includes:

heliocentric corrections
creating one echelle file

Parameters:

continuum (tuple) – results from the continuum normalization
trace (list[TraceData]) – Trace objects with wavelength polynomials
config (dict) – Pipeline configuration

save(i, head, spec, sigma, cont, wave, columns)[source]¶

Save one output spectrum to disk

Parameters:

i (int) – individual number of each file
head (FITS header) – FITS header
spec (array of shape (ntrace, ncol)) – final spectrum
sigma (array of shape (ntrace, ncol)) – final uncertainties
cont (array of shape (ntrace, ncol)) – final continuum scales
wave (array of shape (ntrace, ncol)) – wavelength solution
columns (array of shape (ntrace, 2)) – columns that carry signal

Returns:

out_file – name of the output file

Return type:

str

save_config_to_header(head, config, prefix='PR')[source]¶

class pyreduce.reduce.FitsIOStep(*args, **kwargs)[source]¶

Bases: Step

load(mask)[source]¶

Load the master bias from a previous run

Parameters:

mask (array of shape (nrow, ncol)) – Bad pixel mask

Returns:

data (masked array of shape (nrow, ncol)) – master bias data, with the bad pixel mask applied
head (FITS header) – header of the master bias

save(data, head, dtype=None)[source]¶

Save the data to a FITS file

Parameters:

data (array of shape (nrow, ncol)) – bias data
head (FITS header) – bias header

class pyreduce.reduce.Flat(*args, **config)[source]¶

Bases: CalibrationStep

Calculates the master flat

load(mask)[source]¶

Load master flat from disk

Parameters:

mask (array of shape (nrow, ncol)) – Bad pixel mask

Returns:

flat (masked array of shape (nrow, ncol)) – Master flat with bad pixel map applied
fhead (FITS header) – Master flat FITS header

run(files, bias=None, mask=None)[source]¶

Calculate the master flat, with the bias already subtracted

Parameters:

files (list(str)) – flat files
bias (tuple(array of shape (nrow, ncol), FITS header), optional) – master bias and header
mask (array of shape (nrow, ncol), optional) – Bad pixel mask

Returns:

flat (masked array of shape (nrow, ncol)) – Master flat with bad pixel map applied
fhead (FITS header) – Master flat FITS header

save(flat, fhead)[source]¶

Save the master flat to a FITS file

Parameters:

flat (array of shape (nrow, ncol)) – master flat data
fhead (FITS header) – master flat header

property savefile¶

Name of master bias fits file

Type:: str

class pyreduce.reduce.LaserFrequencyCombFinalize(*args, **config)[source]¶

Bases: Step

Improve the precision of the wavelength calibration with a laser frequency comb

degree¶

polynomial degree of the wavelength fit

Type:: tuple(int, int)

dimensionality¶

Whether to use 1D or 2D polynomials

Type:: {‘1D’, ‘2D’}

lfc_peak_width¶

Width of the peaks for finding them in the spectrum

Type:: int

load()[source]¶: Load is a no-op - wavelengths are in traces.fits.

run(freq_comb_master, trace: list, wavecal: dict)[source]¶

Improve the wavelength calibration with a laser frequency comb.

Updates trace objects in-place with improved wavelength polynomial.

Parameters:

freq_comb_master (tuple) – extracted frequency comb spectrum and header
trace (list[TraceData]) – Trace objects with wavelength polynomials from wavecal
wavecal (dict[str, LineList]) – {group: linelist} from wavecal step (for diagnostics)

save(trace: list)[source]¶

Save updated traces to disk.

Parameters:: trace (list[TraceData]) – Already-updated trace objects

threshold¶

residual threshold in m/s above which to remove lines

Type:: float

class pyreduce.reduce.LaserFrequencyCombMaster(*args, **config)[source]¶

Bases: CalibrationStep, ExtractionStep

Create a laser frequency comb (or similar) master image

load()[source]¶

Load master comb from disk

Returns:

comb (masked array of shape (nrow, ncol)) – Master comb with bad pixel map applied
chead (FITS header) – Master comb FITS header

run(files, trace: list[Trace], mask=None, bias=None, norm_flat=None)[source]¶

Improve the wavelength calibration with a laser frequency comb (or similar)

Parameters:

files (list(str)) – observation files
trace (list[TraceData]) – Trace objects from trace step
mask (array of shape (nrow, ncol), optional) – Bad pixel mask
bias (tuple, optional) – results from the bias step
norm_flat (tuple, optional) – results from the norm_flat step

Returns:

comb (array of shape (ntrace, ncol)) – extracted frequency comb image
chead (Header) – FITS header of the combined image

save(comb, chead)[source]¶

Save the master comb to a FITS file

Parameters:

comb (array of shape (nrow, ncol)) – master comb data
chead (FITS header) – master comb header

property savefile¶

Name of the wavelength echelle file

Type:: str

class pyreduce.reduce.Mask(*args, **config)[source]¶

Bases: Step

Load the bad pixel mask for the given instrument/channel

load()[source]¶

Load the mask file from disk

Returns:: mask – Bad pixel mask for this setting
Return type:: array of shape (nrow, ncol)

run()[source]¶

Load the mask file from disk

Returns:: mask – Bad pixel mask for this setting
Return type:: array of shape (nrow, ncol)

class pyreduce.reduce.NormalizeFlatField(*args, **config)[source]¶

Bases: Step

Calculate the ‘normalized’ flat field image

extraction_kwargs¶

arguments for the extraction

Type:: dict

extraction_method¶

Extraction method to use

Type:: {‘normalize’}

load()[source]¶

Load normalized flat field results from disk

Returns:

norm (array of shape (nrow, ncol)) – normalized flat field
blaze (array of shape (ntrace, ncol)) – Continuum level as determined from the flat field for each order
slitfunc (list of arrays, or None) – Slit function for each order (None if not available)
slitfunc_meta (dict or None) – Metadata for slitfunc (extraction_height, osample, trace_range)

run(flat, trace: list[Trace], scatter=None)[source]¶

Calculate the ‘normalized’ flat field

Parameters:

flat (tuple(array, header)) – Master flat, and its FITS header
trace (list[TraceData]) – Trace objects from trace step
scatter (array, optional) – Background scatter model

Returns:

norm (array of shape (nrow, ncol)) – normalized flat field
blaze (array of shape (ntrace, ncol)) – Continuum level as determined from the flat field for each order
slitfunc (list of arrays) – Slit function for each order
slitfunc_meta (dict) – Metadata for slitfunc (extraction_height, osample, trace_range)

save(norm, blaze, slitfunc, slitfunc_meta)[source]¶

Save normalized flat field results to disk

Parameters:

norm (array of shape (nrow, ncol)) – normalized flat field
blaze (array of shape (ntrace, ncol)) – Continuum level as determined from the flat field for each order
slitfunc (list of arrays) – Slit function for each order
slitfunc_meta (dict) – Metadata for slitfunc (extraction_height, osample, trace_range)

property savefile¶

Name of the blaze file

Type:: str

threshold¶

Threshold of the normalized flat field (values below this are just 1)

Type:: int

class pyreduce.reduce.RectifyImage(*args, **config)[source]¶

Bases: Step

Create a 2D image of the rectified orders

filename(name)[source]¶

load(files)[source]¶

Load results from a previous execution

If this raises a FileNotFoundError, run() will be used instead For calibration steps it is preferred however to print a warning and return None. Other modules can then use a default value instead.

Raises:: NotImplementedError – Needs to be implemented for each step

run(files, trace: list[Trace], mask=None)[source]¶

Execute the current step

This should fail if files are missing or anything else goes wrong. If the user does not want to run this step, they should not specify it in steps.

Parameters:: files (list(str)) – data files required for this step
Raises:: NotImplementedError – needs to be implemented for each step

save(fname, image, wavelength, header=None)[source]¶

Save the results of this step

Parameters:: *args (obj) – things to save
Raises:: NotImplementedError – Needs to be implemented for each step

class pyreduce.reduce.ScienceExtraction(*args, **config)[source]¶

Bases: CalibrationStep, ExtractionStep

Extract the science spectra

load()[source]¶

Load all science spectra from disk.

Supports both new Spectra format (E_FMTVER >= 2) and legacy format.

Returns:

heads (list(FITS header)) – FITS headers of each observation
specs (list(array of shape (ntrace, ncol))) – extracted spectra
sigmas (list(array of shape (ntrace, ncol))) – uncertainties of the extracted spectra
slitfus (list or None) – slit functions (if available)
columns (list(array of shape (ntrace, 2))) – column ranges for each spectra

run(files, trace: list[Trace], bias=None, norm_flat=None, scatter=None, mask=None)[source]¶

Extract Science spectra from observation

Parameters:

files (list(str)) – list of observations
trace (list[TraceData]) – Trace objects from trace step
bias (tuple, optional) – results from master bias step
norm_flat (tuple, optional) – results from flat normalization
scatter (array, optional) – background scatter model
mask (array of shape (nrow, ncol), optional) – bad pixel map

Returns:

heads (list(FITS header)) – FITS headers of each observation
spectra_list (list(list[Spectrum])) – extracted spectra (one list per file)

save(fname, head, spectra: list[Spectrum])[source]¶

Save extracted spectra using Spectra format.

Parameters:

fname (str) – Original filename (used to derive output name)
head (FITS header) – FITS header
spectra (list[Spectrum]) – Extracted spectra from extract()

science_file(name)[source]¶

Name of the science file in disk, based on the input file

Parameters:: name (str) – name of the observation file
Returns:: name – science file name
Return type:: str

class pyreduce.reduce.SlitCurvatureDetermination(*args, **config)[source]¶

Bases: CalibrationStep, ExtractionStep

Determine the curvature of the slit

curvature_mode¶

Whether to use 1d or 2d polynomials

Type:: {‘1D’, ‘2D’}

curve_degree¶

Orders of the curvature to fit, currently supports only 1 and 2

Type:: int

curve_height¶

height of the 2D cutout for curvature fitting

Type:: float

extraction_height¶

extraction height for peak finding spectrum

Type:: float

fit_degree¶

Polynomial degree of the overall fit

Type:: int

load()[source]¶: Curvature is now stored in traces, not separate files.

peak_function¶

Function shape that is fit to individual peaks

Type:: str

peak_threshold¶

peak finding noise threshold

Type:: float

peak_width¶

peak width

Type:: int

run(files, trace: list[Trace], mask=None, bias=None)[source]¶

Determine the curvature of the slit

Parameters:

files (list(str)) – files to use for this
trace (list[TraceData]) – Trace objects from trace step
mask (array of shape (nrow, ncol), optional) – Bad pixel mask
bias (tuple, optional) – Master bias

Returns:

curvature – Slit curvature data including polynomial coefficients

Return type:

SlitCurvature

save(traces)[source]¶

Save curvature results by updating traces.fits.

Parameters:: traces (list[Trace]) – Traces with updated slit/slitdelta data

sigma_cutoff¶

how many sigma of bad lines to cut away

Type:: float

window_width¶

window width to search for peak in each row

Type:: float

class pyreduce.reduce.Step(instrument, channel, target, night, output_dir, trace_range, **config)[source]¶

Bases: object

Parent class for all steps

channel¶

Name of the instrument channel

Type:: str

property dependsOn¶

Steps that are required before running this step

Type:: list(str)

files¶

Input files dict, set by pipeline before load()

Type:: dict

instrument¶

Name of the instrument

Type:: str

load()[source]¶

Load results from a previous execution

If this raises a FileNotFoundError, run() will be used instead For calibration steps it is preferred however to print a warning and return None. Other modules can then use a default value instead.

Raises:: NotImplementedError – Needs to be implemented for each step

property loadDependsOn¶

Steps that are required before loading data from this step

Type:: list(str)

night¶

Date of the observation (as a string)

Type:: str

property output_dir¶

output directory, may contain tags {instrument}, {night}, {target}, {channel}

Type:: str

plot¶

Whether to plot the results or the progress of this step

Type:: bool

plot_title¶

Title used in the plots, if any

Type:: str

property prefix¶

temporary file prefix

Type:: str

run(files, *args)[source]¶

Execute the current step

This should fail if files are missing or anything else goes wrong. If the user does not want to run this step, they should not specify it in steps.

Parameters:: files (list(str)) – data files required for this step
Raises:: NotImplementedError – needs to be implemented for each step

save(*args)[source]¶

Save the results of this step

Parameters:: *args (obj) – things to save
Raises:: NotImplementedError – Needs to be implemented for each step

target¶

Name of the observation target

Type:: str

trace_range¶

First and Last(+1) trace to process

Type:: tuple(int, int)

class pyreduce.reduce.Trace(*args, **config)[source]¶

Bases: CalibrationStep

Determine the polynomial fits describing the pixel locations of each trace

border_width¶

Number of pixels at the edge of the detector to ignore

Type:: int

filter_type¶

Type of smoothing filter (boxcar, gaussian, whittaker)

Type:: str

filter_x¶

Smoothing width along x-axis (dispersion direction)

Type:: int

filter_y¶

Smoothing width along y-axis (cross-dispersion direction)

Type:: int

fit_degree¶

Polynomial degree of the fit to each order

Type:: int

get_traces_for_step(step_name: str) → dict[str, list[Trace]][source]¶

Get traces appropriate for a specific reduction step.

Uses the instrument’s fibers.use config to select traces.

Parameters:: step_name (str) – Name of the reduction step (e.g., “science”, “curvature”)
Returns:: {group_name: [traces]} for each selected group
Return type:: dict[str, list[TraceData]]

load()[source]¶

Load tracing results from FITS format.

Returns:: Trace objects with position, column_range, height, and identity.
Return type:: list[TraceData]

manual¶

Whether to use manual alignment

Type:: bool

min_cluster¶

Minimum size of each cluster to be included in further processing

Type:: int

min_width¶

Minimum width of each cluster after mergin

Type:: int, float

noise¶

Absolute background noise threshold

Type:: int

noise_relative¶

Relative background noise threshold (fraction of background)

Type:: float

run(files, mask=None, bias=None)[source]¶

Determine polynomial coefficients describing order locations

Parameters:

files (list(str)) – Observation used for order tracing (should only have one element)
mask (array of shape (nrow, ncol), optional) – Bad pixel mask
bias (tuple, optional) – Bias correction

Returns:

Trace objects with position, column_range, height, and identity.

Return type:

list[TraceData]

save()[source]¶: Save tracing results to disk in FITS format.

property savefile¶

Name of the tracing file (FITS format)

Type:: str

trace_objects: list[Trace]¶

class pyreduce.reduce.WavelengthCalibrationFinalize(*args, **config)[source]¶

Bases: Step

Perform wavelength calibration

atlas_name¶

name of the line atlas

Type:: str

correlate_cols¶

How many columns to use in the 2D cross correlation alignment. 0 means all pixels (slow).

Type:: int

degree¶

Polynomial degree of the wavelength calibration in order, column direction

Type:: tuple(int, int)

dimensionality¶

Whether to use 1d or 2d polynomials

Type:: {‘1D’, ‘2D’}

iterations¶

Number of iterations in the remove lines, auto id cycle

Type:: int

load()[source]¶

Load wavelength calibration linelists.

Wavelength data is stored in traces.fits, not returned here.

Returns:: results – {group: linelist} for each fiber group
Return type:: dict[str, LineList]

manual¶

Whether to use manual alignment instead of cross correlation

Type:: bool

medium¶

medium of the detector, vac or air

Type:: str

nstep¶

Number of detector offset steps, due to detector design

Type:: int

run(wavecal_master: dict, wavecal_init: dict, trace: list)[source]¶

Perform wavelength calibration for each fiber group.

Fits wavelength polynomials and updates trace objects in-place. Returns linelists for diagnostics.

Parameters:

wavecal_master (dict[str, tuple]) – {group: (wavecal_spec, thead)} from wavecal_master step
wavecal_init (dict[str, LineList]) – {group: linelist} from wavecal_init step
trace (list[TraceData]) – Trace objects to update with wavelength polynomials

Returns:

results – {group: linelist} for each fiber group (wavelengths are in traces)

Return type:

dict[str, LineList]

save(results: dict, trace: list)[source]¶

Save linelists and updated traces to disk.

Parameters:

results (dict[str, tuple]) – {group: (wave, linelist)} - wave polynomials and linelists
trace (list[TraceData]) – Already-updated trace objects

property savefile¶

Name of the linelist file (single-group compat)

Type:: str

savefile_for_group(group: str) → str[source]¶: Get savefile path for a specific group.

shift_window¶

fraction of columns, to allow individual orders to shift

Type:: float

threshold¶

residual threshold in m/s

Type:: float

class pyreduce.reduce.WavelengthCalibrationInitialize(*args, **config)[source]¶

Bases: Step

Create the initial wavelength solution file

load(config, wavecal_master: dict)[source]¶

Load linelists for each fiber group.

Falls back to instrument-provided wavecal file if custom not found.

run(wavecal_master: dict)[source]¶

Run iterative line matching for each fiber group.

Parameters:: wavecal_master (dict[str, tuple]) – {group: (wavecal_spec, thead)} from wavecal_master step
Returns:: results – {group: linelist} for each fiber group
Return type:: dict[str, LineList]

save(results: dict)[source]¶: Save linelists for each fiber group.

property savefile¶

Name of the linelist file (single-group compat)

Type:: str

savefile_for_group(group: str) → str[source]¶: Get savefile path for a specific group.

class pyreduce.reduce.WavelengthCalibrationMaster(*args, **config)[source]¶

Bases: CalibrationStep, ExtractionStep

Create wavelength calibration master image

load()[source]¶

Load master wavelength calibration from disk.

Returns:: results – {group: (wavecal_spec, thead)} for each fiber group
Return type:: dict[str, tuple]

run(files, trace: list[Trace], mask=None, bias=None, norm_flat=None)[source]¶

Extract wavelength calibration spectra, per fiber group.

Parameters:

files (list(str)) – wavelength calibration files
trace (list[TraceData]) – Trace objects from trace step
mask (array of shape (nrow, ncol), optional) – Bad pixel mask
bias (tuple, optional) – Master bias
norm_flat (tuple, optional) – Normalized flat field

Returns:

results – {group: (wavecal_spec, thead)} for each fiber group

Return type:

dict[str, tuple]

save(results: dict)[source]¶

Save the master wavelength calibration to FITS files.

Parameters:: results (dict[str, tuple]) – {group: (wavecal_spec, thead)} for each fiber group

property savefile¶

Name of the wavelength echelle file (single-group compat)

Type:: str

savefile_for_group(group: str) → str[source]¶: Get savefile path for a specific group.

pyreduce.reduce.main(instrument, target, night=None, channels=None, steps='all', base_dir=None, input_dir=None, output_dir=None, configuration=None, trace_range=None, skip_existing=False, plot=0, plot_dir=None, use_groups=None)[source]¶

Main entry point for REDUCE scripts.

Default values can be changed as required if reduce is used as a script. Finds input directories, and loops over observation nights and instrument channels.

Deprecated since version Use: Pipeline.from_instrument() instead.

Parameters:

instrument (str, list[str]) – instrument used for the observation (e.g. UVES, HARPS)
target (str, list[str]) – the observed star, as named in the folder structure/fits headers
night (str, list[str]) – the observation nights to reduce, as named in the folder structure. Accepts bash wildcards (i.e. *, ?), but then relies on the folder structure for restricting the nights
channels (str, list[str], dict[{instrument}:list], None, optional) – the instrument channels to use, if None will use all known channels for the current instrument. See instruments for possible options
steps (tuple(str), "all", optional) – which steps of the reduction process to perform the possible steps are: “bias”, “flat”, “trace”, “norm_flat”, “wavecal”, “science” alternatively set steps to “all”, which is equivalent to setting all steps Note that the later steps require the previous intermediary products to exist and raise an exception otherwise
base_dir (str, optional) – base data directory that Reduce should work in, is prefixxed on input_dir and output_dir (default: use settings_pyreduce.json)
input_dir (str, optional) – input directory containing raw files. Can contain placeholders {instrument}, {target}, {night}, {channel} as well as wildcards. If relative will use base_dir as root (default: use settings_pyreduce.json)
output_dir (str, optional) – output directory for intermediary and final results. Can contain placeholders {instrument}, {target}, {night}, {channel}, but no wildcards. If relative will use base_dir as root (default: use settings_pyreduce.json)
configuration (dict[str:obj], str, list[str], dict[{instrument}:dict,str], optional) – configuration file for the current run, contains parameters for different parts of reduce. Can be a path to a json file, or a dict with configurations for the different instruments. When a list, the order must be the same as instruments (default: settings_{instrument.upper()}.json)

pyreduce.reduce.wavelengths_from_traces(traces: list, ncol: int = None) → ndarray[source]¶

Compute wavelength array from trace objects.

Parameters:

traces (list[TraceData]) – Trace objects with .wave polynomial coefficients set
ncol (int, optional) – Number of columns. If not provided, uses max column_range.

Returns:

wlen – Wavelength for each pixel, or None if no wavelength data

Return type:

ndarray of shape (ntrace, ncol)

pyreduce.configuration module¶

Loads configuration files

This module loads json configuration files from disk, and combines them with the default settings, to create one dict that contains all parameters. It also checks that all parameters exists, and that no new parameters have been added by accident.

pyreduce.configuration.get_configuration_for_instrument(instrument, channel=None, channel_fallbacks=None, **kwargs)[source]¶

pyreduce.configuration.load_config(configuration, instrument, j=0, channel=None, channel_fallbacks=None)[source]¶

pyreduce.configuration.load_settings_override(config, settings_file)[source]¶

Apply settings overrides from a JSON file.

Parameters:

config (dict) – Base configuration to override
settings_file (str) – Path to JSON file with override settings

Returns:

config – Updated configuration

Return type:

dict

pyreduce.configuration.read_config(fname=None)[source]¶

Read the configuration file from disk

If no filename is given it will load the default configuration. The configuration file must be a json file.

Parameters:: fname (str, optional) – Filename of the configuration. By default the default settings.
Returns:: config – The read configuration file
Return type:: dict

pyreduce.configuration.update(dict1, dict2, check=True, name='dict1')[source]¶

Update entries in dict1 with entries of dict2 recursively, i.e. if the dict contains a dict value, values inside the dict will also be updated

Parameters:

dict1 (dict) – dict that will be updated
dict2 (dict) – dict that contains the values to update
check (bool) – If True, will check that the keys from dict2 exist in dict1 already. Except for those contained in field “instrument”

Returns:

dict1 – the updated dict

Return type:

dict

Raises:

KeyError – If dict2 contains a key that is not in dict1

pyreduce.configuration.validate_config(config)[source]¶

Test that the input configuration complies with the expected schema

Since it requires features from jsonschema 3+, it will only run if that is installed. Otherwise show a warning but continue. This is incase some other module needs an earlier, jsonschema (looking at you jwst).

If the function runs through without raising an exception, the check was succesful or skipped.

Parameters:: config (dict) – Configurations to check
Raises:: ValueError – If there is a problem with the configuration. Usually that means a setting has an unallowed value.

pyreduce.extract module¶

Module for extracting data from observations

Authors¶

Version¶

License¶

class pyreduce.extract.ProgressPlot(nrow, ncol, nslitf, title=None)[source]¶

Bases: object

close()[source]¶

get_slitf(img, spec, slitf, ycen, slitcurve=None, slitdeltas=None)[source]¶: get the slit function

get_spec(img, spec, slitf, ycen, slitcurve=None, slitdeltas=None)[source]¶: get the spectrum corrected by the slit function

plot(img, spec, slitf, model, ycen, input_mask, output_mask, trace_idx=0, left=0, right=0, unc=None, info=None, swath_idx=0, save=True, slitcurve=None, slitdeltas=None)[source]¶

wait_if_paused()[source]¶

class pyreduce.extract.Swath(nswath)[source]¶: Bases: object

pyreduce.extract.calc_scatter_correction(scatter, index)[source]¶

Calculate scatter correction by interpolating between values?

Parameters:

scatter (array of shape (degree_x, degree_y)) – 2D polynomial coefficients of the background scatter
index (tuple (array, array)) – indices of the swath within the overall image

Returns:

scatter_correction – correction for scattered light

Return type:

array of shape (swath_width, swath_height)

pyreduce.extract.calc_telluric_correction(telluric, img)[source]¶

Calculate telluric correction

If set to specific integer larger than 1 is used as the offset from the order center line. The sky is then estimated by computing median signal between this offset and the upper/lower limit of the extraction window.

Parameters:

telluric (int) – telluric correction parameter
img (array) – image of the swath

Returns:

tell – telluric correction

Return type:

array

pyreduce.extract.correct_for_curvature(img_order, curvature, xwd, inverse=False)[source]¶

Correct image for slit curvature by interpolation.

Parameters:

img_order (array[nrow, ncol]) – Image swath to correct
curvature (array[ncol, n_coeffs]) – Curvature coefficients [c0, c1, c2, …] where dx = c1*y + c2*y^2 + …
xwd (int) – Extraction full height in pixels
inverse (bool) – If True, apply inverse correction (for model reapplication)

Returns:

img_order – Corrected image

Return type:

array

pyreduce.extract.extend_traces(traces, nrow)[source]¶

Extrapolate extra traces above and below the existing ones

Parameters:

traces (array[ntrace, degree]) – trace polynomial coefficients
nrow (int) – number of rows in the image

Returns:

traces – extended traces

Return type:

array[ntrace + 2, degree]

pyreduce.extract.extract(img, traces: list[Trace], extraction_height: float = 0.5, extraction_type: str = 'optimal', **kwargs) → list[Spectrum][source]¶

Extract spectra from an image.

Parameters:

img (array[nrow, ncol](float)) – Observation to extract.
traces (list[Trace]) – Trace objects with position, column_range, and optional slit curvature.
extraction_height (float, optional) – Extraction height. Values below 2 are fractions of trace spacing, values above are pixels. If None, falls back to trace.height. (default: 0.5)
extraction_type ({"optimal", "simple"}, optional) – Extraction algorithm. (default: “optimal”)
**kwargs – Additional parameters for extraction functions (osample, lambda_sf, etc.)

Returns:

Extracted spectrum objects, one per trace.

Return type:

list[Spectrum]

pyreduce.extract.extract_normalize(img, traces: list[Trace], extraction_height: float = 0.5, **kwargs)[source]¶

Extract and normalize flat field image.

This is a specialized extraction mode for flat field processing that returns normalized images rather than Spectrum objects.

Parameters:

img (array[nrow, ncol](float)) – Flat field image to normalize.
traces (list[Trace]) – Trace objects with position, column_range, height, slit curvature.
extraction_height (float, optional) – Extraction height. If None, falls back to trace.height.
**kwargs – Additional parameters for extraction.

Returns:

im_norm (array[nrow, ncol](float)) – Normalized flat field image.
im_ordr (array[nrow, ncol](float)) – Image with just the trace regions.
blaze (array[ntrace, ncol](float)) – Extracted blaze function.
slitfunction (list) – Recovered slit functions.
column_range (array[ntrace, 2](int)) – Column ranges used.

pyreduce.extract.extract_spectrum(img, ycen, yrange, xrange, gain=1, readnoise=0, lambda_sf=0.1, lambda_sp=0, osample=1, swath_width=None, maxiter=20, reject_threshold=6, telluric=None, scatter=None, normalize=False, threshold=0, curvature=None, slitdeltas=None, plot=False, plot_title=None, im_norm=None, im_ordr=None, out_spec=None, out_sunc=None, out_slitf=None, out_mask=None, progress=None, ord_num=0, preset_slitfunc=None, **kwargs)[source]¶

Extract the spectrum of a single order from an image The order is split into several swathes of roughly swath_width length, which overlap half-half For each swath a spectrum and slitfunction are extracted overlapping sections are combined using linear weights (centrum is strongest, falling off to the edges) Here is the layout for the bins:

   1st swath    3rd swath    5th swath      ...
/============|============|============|============|============|

          2nd swath    4th swath    6th swath
       |------------|------------|------------|------------|
       |.....|
       overlap

       +     ******* 1
        +   *
         + *
          *            weights (+) previous swath, (*) current swath
         * +
        *   +
       *     +++++++ 0

Parameters:

img (array[nrow, ncol]) – observation (or similar)
ycen (array[ncol]) – order trace of the current order
yrange (tuple(int, int)) – extraction width in pixles, below and above
xrange (tuple(int, int)) – columns range to extract (low, high)
gain (float, optional) – adu to electron, amplifier gain (default: 1)
readnoise (float, optional) – read out noise factor (default: 0)
lambda_sf (float, optional) – slit function smoothing parameter, usually very small (default: 0.1)
lambda_sp (int, optional) – spectrum smoothing parameter, usually very small (default: 0)
osample (int, optional) – oversampling factor, i.e. how many subpixels to create per pixel (default: 1, i.e. no oversampling)
swath_width (int, optional) – swath width suggestion, actual width depends also on ncol, see make_bins (default: None, which will determine the width based on the order tracing)
telluric ({float, None}, optional) – telluric correction factor (default: None, i.e. no telluric correction)
scatter ({array, None}, optional) – background scatter as 2d polynomial coefficients (default: None, no correction)
normalize (bool, optional) – whether to create a normalized image. If true, im_norm and im_ordr are used as output (default: False)
threshold (int, optional) – threshold for normalization (default: 0)
curvature (array[ncol, n_coeffs], optional) – Slit curvature polynomial coefficients for this trace (default: None, i.e. vertical extraction)
slitdeltas (array[nrows_stored], optional) – Per-row residual offsets from polynomial curvature fit (default: None). Will be interpolated to match swath nrows if lengths differ.
plot (bool, optional) – wether to plot the progress, plotting will slow down the procedure significantly (default: False)
ord_num (int, optional) – current order number, just for plotting (default: 0)
im_norm (array[nrow, ncol], optional) – normalized image, only output if normalize is True (default: None)
im_ordr (array[nrow, ncol], optional) – image of the order blaze, only output if normalize is True (default: None)

Returns:

spec (array[ncol]) – extracted spectrum
slitf (array[nslitf]) – extracted slitfunction
mask (array[ncol]) – mask of the column range to use in the spectrum
unc (array[ncol]) – uncertainty on the spectrum

pyreduce.extract.fix_column_range(column_range, traces, extraction_height, nrow, ncol)[source]¶

Fix the column range, so that no pixels outside the image will be accessed (Thus avoiding errors)

Parameters:

img (array[nrow, ncol]) – image
traces (array[ntrace, degree]) – trace polynomial coefficients
extraction_height (array[ntrace]) – extraction full height in pixels
column_range (array[ntrace, 2]) – current column range
no_clip (bool, optional) – if False, new column range will be smaller or equal to current column range, otherwise it can also be larger (default: False)

Returns:

column_range (array[ntrace, 2]) – updated column range
traces (array[ntrace, degree]) – trace polynomial coefficients (may have rows removed if no valid pixels)

pyreduce.extract.fix_extraction_height(xwd, traces, cr, ncol)[source]¶

Convert fractional extraction height to pixel range.

Fractions (< 2) are multiplied by the minimum distance to neighboring traces.

Parameters:

xwd (array[ntrace]) – extraction full height per trace
traces (array[ntrace, degree]) – trace polynomial coefficients
cr (array[ntrace, 2]) – column range to use
ncol (int) – number of columns in image

Returns:

xwd – updated extraction full height in pixels

Return type:

array[ntrace]

pyreduce.extract.fix_parameters(xwd, cr, traces, nrow, ncol, ntrace, ignore_column_range=False)[source]¶

Fix extraction width and column range, so that all pixels used are within the image. I.e. the column range is cut so that the everything is within the image

Parameters:

xwd (float) – Total extraction height. Split evenly above/below trace. Values below 3 are fractions of trace spacing.
cr (2-tuple(int), array) – Column range, either one value for all traces, or the whole array
traces (array) – polynomial coefficients that describe each trace
nrow (int) – Number of rows in the image
ncol (int) – Number of columns in the image
ntrace (int) – Number of traces in the image
ignore_column_range (bool, optional) – if true does not change the column range, however this may lead to problems with the extraction, by default False

Returns:

xwd (array) – fixed extraction width
cr (array) – fixed column range
traces (array) – the same traces as before

pyreduce.extract.get_y_scale(ycen, xrange, extraction_height, nrow)[source]¶

Calculate the y limits of the order for C extraction code.

Parameters:

ycen (array[ncol]) – order trace
xrange (tuple(int, int)) – column range
extraction_height (int) – extraction full height in pixels
nrow (int) – number of rows in the image, defines upper edge

Returns:

ylow, yhigh – lower and upper y bound for extraction (pixels below/above trace) These satisfy: ylow + yhigh + 1 = extraction_height

Return type:

int, int

pyreduce.extract.make_bins(swath_width, xlow, xhigh, ycen)[source]¶

Create bins for the swathes Bins are roughly equally sized, have roughly length swath width (if given) and overlap roughly half-half with each other

Parameters:

swath_width ({int, None}) – initial value for the swath_width, bins will have roughly that size, but exact value may change if swath_width is None, determine a good value, from the data
xlow (int) – lower bound for x values
xhigh (int) – upper bound for x values
ycen (array[ncol]) – center of the order trace

Returns:

nbin (int) – number of bins
bins_start (array[nbin]) – left(beginning) side of the bins
bins_end (array[nbin]) – right(ending) side of the bins

pyreduce.extract.model(spec, slitf)[source]¶

pyreduce.extract.optimal_extraction(img, traces, extraction_height, column_range, curvature=None, slitdeltas=None, plot=False, plot_title=None, **kwargs)[source]¶

Use optimal extraction to get spectra

This functions just loops over the traces, the actual work is done in extract_spectrum

Parameters:

img (array[nrow, ncol]) – image to extract
traces (array[ntrace, degree]) – trace polynomial coefficients
extraction_height (array[ntrace]) – extraction full height in pixels
column_range (array[ntrace, 2]) – column range to use
curvature (array[ntrace, ncol, n_coeffs] or None) – Slit curvature polynomial coefficients (default: None for vertical extraction)
slitdeltas (array[ntrace, nrows] or None) – Per-row residual offsets from curvature fit (default: None)
**kwargs – other parameters for the extraction (see extract_spectrum)

Returns:

spectrum (array[ntrace, ncol]) – extracted spectrum
slitfunction (array[ntrace, nslitf]) – recovered slitfunction
uncertainties (array[ntrace, ncol]) – uncertainties on the spectrum

pyreduce.extract.plot_comparison(original, traces, spectrum, slitf, extraction_height, column_range, title=None)[source]¶

pyreduce.extract.simple_extraction(img, traces, extraction_height, column_range, gain=1, readnoise=0, dark=0, plot=False, plot_title=None, curvature=None, collapse_function='median', **kwargs)[source]¶

Use simple extraction to get a spectrum Simple extraction takes the sum/mean/median orthogonal to the trace for extraction_height pixels

This extraction makes a few rough assumptions and does not provide the most accurate results, but rather a good approximation

Parameters:

img (array[nrow, ncol]) – image to extract
traces (array[ntrace, degree]) – trace polynomial coefficients
extraction_height (array[ntrace]) – extraction full height in pixels
column_range (array[ntrace, 2]) – column range to use
gain (float, optional) – adu to electron, amplifier gain (default: 1)
readnoise (float, optional) – read out noise (default: 0)
dark (float, optional) – dark current noise (default: 0)
plot (bool, optional) – wether to plot the results (default: False)

Returns:

spectrum (array[ntrace, ncol]) – extracted spectrum
uncertainties (array[ntrace, ncol]) – uncertainties on extracted spectrum

pyreduce.extract.validate_traces_for_extraction(traces: list[Trace], extraction_height: float | ndarray, nrow: int, ncol: int) → None[source]¶

Validate traces and mark invalid ones.

Checks if each trace’s extraction aperture fits within the image bounds. Invalid traces are marked with trace.invalid = “reason”.

Parameters:

traces (list[Trace]) – Trace objects to validate. Modified in-place.
extraction_height (float or array) – Extraction height(s) in pixels.
nrow (int) – Number of rows in image.
ncol (int) – Number of columns in image.

pyreduce.trace module¶

Find clusters of pixels with signal and fit polynomial traces.

Note on terminology: - “trace”: A single polynomial fit to a cluster of pixels (e.g., one fiber) - “spectral order”: A group of traces at similar wavelengths (e.g., all fibers in one echelle order)

The main function trace() detects and fits individual traces, returning Trace objects. Use group_fibers() to merge traces into fiber groups according to instrument config.

pyreduce.trace.best_fit(x, y)[source]¶

pyreduce.trace.calculate_mean_cluster_thickness(x, y)[source]¶

pyreduce.trace.combine(i, j, x, y, merge, mct, nrow, ncol, deg, threshold)[source]¶

pyreduce.trace.create_merge_array(x, y, mean_cluster_thickness, nrow, ncol, deg, threshold)[source]¶

pyreduce.trace.delete(i, x, y, merge)[source]¶

pyreduce.trace.determine_overlap_rating(xi, yi, xj, yj, mean_cluster_thickness, nrow, ncol, deg=2)[source]¶

pyreduce.trace.fit(x, y, deg, regularization=0)[source]¶

pyreduce.trace.fit_polynomials_to_clusters(x, y, clusters, degree, regularization=0)[source]¶

Fits a polynomial of degree opower to points x, y in cluster clusters

Parameters:

x (dict(int: array)) – x coordinates seperated by cluster
y (dict(int: array)) – y coordinates seperated by cluster
clusters (list(int)) – cluster labels, equivalent to x.keys() or y.keys()
degree (int) – degree of polynomial fit

Returns:

traces – coefficients of polynomial fit for each cluster

Return type:

dict(int, array[degree+1])

pyreduce.trace.group_fibers(traces: list[Trace], fibers_config, degree: int = 4) → list[Trace][source]¶

Merge individual fiber traces into groups according to config.

Takes traces with fiber_idx set (individual fibers) and returns NEW traces with group set (merged/grouped result). The input traces are not modified.

Parameters:

traces (list[Trace]) – Individual fiber traces with fiber_idx set and group=None. Can have m set (from order_centers) or None (to be assigned later).
fibers_config (FibersConfig) – Configuration specifying groups or bundles.
degree (int, optional) – Polynomial degree for refitted traces (used with “average” merge).

Returns:

New Trace objects with: - group: set to group name (e.g., “A”, “B”, “bundle_1”) - fiber_idx: None (merged, no individual identity) - m: preserved from input traces - pos: new polynomial from merging - column_range: intersection of member traces - height: from config or computed from member traces

Returns empty list if no grouping config is provided.

Return type:

list[Trace]

pyreduce.trace.merge_clusters(img, x, y, n_clusters, manual=True, deg=2, auto_merge_threshold=0.9, merge_min_threshold=0.1, plot_title=None)[source]¶

Merge clusters that belong together

Parameters:

img (array[nrow, ncol]) – the image the order trace is based on
x (dict(int, array(int))) – x coordinates of cluster points
y (dict(int, array(int))) – y coordinates of cluster points
n_clusters (array(int)) – cluster numbers
manual (bool, optional) – if True ask before merging clusters (default: True)
deg (int, optional) – polynomial degree for fitting (default: 2)
auto_merge_threshold (float, optional) – overlap threshold for automatic merging (default: 0.9)
merge_min_threshold (float, optional) – minimum overlap to consider merging (default: 0.1)
plot_title (str, optional) – title for plots

Returns:

x (dict(int: array)) – x coordinates of clusters, key=cluster id
y (dict(int: array)) – y coordinates of clusters, key=cluster id
n_clusters (array(int)) – cluster labels

pyreduce.trace.plot_trace_pair(i, j, x, y, img, deg, title='')[source]¶: Plot two trace candidates for merge decision

pyreduce.trace.plot_traces(im, traces, title=None)[source]¶

Plot traces and image.

Parameters:

im (ndarray) – Input image
traces (list[Trace]) – Trace objects to plot
title (str, optional) – Plot title

pyreduce.trace.select_traces_for_step(traces: list[Trace], fibers_config, step_name: str) → dict[str, list[Trace]][source]¶

Select which traces to use for a given reduction step.

Looks up fibers_config.use[step_name] to determine selection mode.

Parameters:

traces (list[Trace]) – All trace objects
fibers_config (FibersConfig or None) – Fiber configuration (may be None for single-fiber instruments)
step_name (str) – Name of the reduction step (e.g., “science”, “curvature”)

Returns:

selected – {group_name: [traces]} for each selected group

Return type:

dict[str, list[Trace]]

pyreduce.trace.trace(im, min_cluster=None, min_width=None, filter_x=0, filter_y=None, filter_type='boxcar', noise=0, noise_relative=0, degree=4, border_width=None, degree_before_merge=2, regularization=0, closing_shape=(5, 5), opening_shape=(2, 2), plot=False, plot_title=None, manual=True, auto_merge_threshold=0.9, merge_min_threshold=0.1, sigma=0, debug_dir=None, order_centers: dict[int, float] | None = None, fibers_per_order: int | None = None)[source]¶

Identify and trace orders, returning Trace objects.

Parameters:

im (array[nrow, ncol]) – order definition image
min_cluster (int, optional) – minimum cluster size in pixels (default: 500)
filter_x (int, optional) – Smoothing width along x-axis/dispersion direction (default: 0, no smoothing). Useful for noisy data or thin fiber traces.
filter_y (int, optional) – Smoothing width along y-axis/cross-dispersion direction (default: auto). Used to estimate local background. For thin closely-spaced traces, use small values.
filter_type (str, optional) – Type of smoothing filter: “boxcar” (default), “gaussian”, or “whittaker”. Boxcar is a uniform moving average. Whittaker preserves edges better.
noise (float, optional) – Absolute noise threshold added to background (default: 0).
noise_relative (float, optional) – Relative noise threshold as fraction of background (default: 0). If both noise and noise_relative are 0, defaults to 0.001 (0.1%).
degree (int, optional) – polynomial degree of the order fit (default: 4)
border_width (int or list of 4 int, optional) – Pixels to ignore at image edges for order tracing. If int, same value applied to all edges. If list: [top, bottom, left, right] for per-side control.
plot (bool, optional) – wether to plot the final order fits (default: False)
manual (bool, optional) – wether to manually select clusters to merge (strongly recommended) (default: True)
debug_dir (str, optional) – if set, write intermediate images (filtered, background, mask) to this directory
order_centers (dict[int, float], optional) – Mapping of order number (m) -> y-position at detector center. If provided, traces are assigned m values by matching to these centers. Otherwise, m remains None (to be assigned later from wavecal obase).
fibers_per_order (int, optional) – Number of fiber traces per spectral order. When set and order_centers is None, consecutive traces (sorted by y) are grouped into orders of this size. Used for instruments like HARPSpol where a Wollaston prism splits each order into multiple beams.

Returns:

Trace objects with: - m: assigned from order_centers if provided, else None - fiber_idx: 1 for single-fiber, or sequential within each order for multi-fiber - group: None (not yet grouped) - pos: polynomial coefficients - column_range: valid column range - height: computed from neighbor distances (None for single trace)

Return type:

list[Trace]

pyreduce.trace.update_merge_array(merge, x, y, j, mean_cluster_thickness, nrow, ncol, deg, threshold)[source]¶

pyreduce.trace.whittaker_smooth(y, lam, axis=0)[source]¶

Whittaker smoother (optimal filter).

Solves: min sum((y - z)^2) + lam * sum((z[i] - z[i-1])^2)

Parameters:

y (array) – Input data (1D or 2D)
lam (float) – Smoothing parameter (higher = smoother)
axis (int) – Axis along which to smooth (for 2D arrays)

Returns:

z – Smoothed data

Return type:

array

pyreduce.trace_model module¶

Trace data model for PyReduce.

This module defines the Trace dataclass and I/O functions for storing trace positions, curvature, and wavelength calibration in FITS format.

The Trace dataclass consolidates what was previously scattered across separate files (traces.npz, curve.npz, wavecal.npz) into a single structure.

Bases: object

Container for a single trace’s geometry and calibration data.

A trace represents a single spectral order (or fiber within an order) on the detector.

m¶

Spectral order number (diffraction order). This is the physical order number from the grating equation, not a sequential index. In echelle spectrographs, higher order numbers correspond to shorter wavelengths.

The order number is assigned in one of three ways:

From order_centers.yaml (preferred): If the instrument provides an order_centers_{channel}.yaml file with known order positions, traces are matched to these centers during detection and assigned the corresponding order numbers immediately.
From wavelength calibration: If no order_centers file exists, m is initially None. During wavelength calibration, the linelist file provides obase (the base order number). Each trace is then assigned m = obase + trace_index.
Sequential fallback: For legacy files or MOSAIC mode where order identity cannot be determined, m may remain None or be assigned sequentially from 0.

The order number is critical for 2D wavelength calibration, which fits a polynomial in both pixel position (x) and order number (m). When evaluating wavelengths via Trace.wlen(), the trace’s m value is used as the second coordinate in the 2D polynomial.

Type:: int | None

group¶

Group identifier, or None if trace is ungrouped. When set, indicates this trace is the result of grouping/merging fibers for this order. There should be exactly one trace per (m, group). String for named groups (‘A’, ‘B’, ‘cal’), int for bundle indices.

Mutually exclusive with fiber_idx. A trace has either: - group set (merged/grouped result) and fiber_idx=None, or - fiber_idx set (individual fiber) and group=None, or - both None (ungrouped single-fiber instrument)

Type:: str | int | None

fiber_idx¶

Physical fiber index (1-indexed). For multi-fiber instruments where fibers are tracked individually (not merged). Used for per-fiber wavelength calibration. There should be exactly one trace per (m, fiber_idx).

Mutually exclusive with group. See group docstring for details.

Type:: int | None

pos¶

y(x) trace position polynomial coefficients, shape (deg+1,). Coefficients in numpy.polyval order (highest power first).

Type:: np.ndarray

column_range¶

Valid x range [start, end) for this trace.

Type:: tuple[int, int]

height¶

Extraction aperture height in pixels. None to use settings default.

Type:: float | None

slit¶

Slit curvature coefficients, shape (deg_y+1, deg_x+1). Evaluates to x_offset = P(y) where P’s coefficients vary with x. slit[i, :] are coefficients for the y^i term as a function of x.

Type:: np.ndarray | None

slitdelta¶

Per-row slit correction, shape (height_pixels,). Residual offsets beyond polynomial fit.

Type:: np.ndarray | None

wave¶

Wavelength polynomial coefficients. Can be: - 1D array, shape (deg+1,): per-trace polynomial, wavelength = polyval(x) - 2D array, shape (deg_x+1, deg_m+1): global 2D polynomial shared across

all traces. Wavelength = polyval2d(x, m) where m is this trace’s order.

Type:: np.ndarray | None

column_range: tuple[int, int]¶

fiber_idx: int | None = None¶

group: str | int | None = None¶

height: float | None = None¶

invalid: str | None = None¶

m: int | None¶

pos: ndarray¶

slit: ndarray | None = None¶

slit_at_x(x: float | ndarray) → ndarray | None[source]¶

Evaluate slit polynomial coefficients at position x.

Parameters:: x (float or np.ndarray) – Column position(s) to evaluate at.
Returns:: Polynomial coefficients for y_offset = c0 + c1*y + c2*y^2 + … Shape (deg_y+1,) for scalar x, or (len(x), deg_y+1) for array x. Returns None if no slit curvature is set.
Return type:: np.ndarray or None

slitdelta: ndarray | None = None¶

wave: ndarray | None = None¶

wlen(x: ndarray) → ndarray | None[source]¶

Evaluate wavelength polynomial at column positions.

Parameters:: x (np.ndarray) – Column positions to evaluate at.
Returns:: Wavelength values at each x position. Returns None if no wavelength calibration is set.
Return type:: np.ndarray or None

y_at_x(x: ndarray) → ndarray[source]¶

Evaluate trace y-position at column positions.

Parameters:: x (np.ndarray) – Column positions to evaluate at.
Returns:: Y positions of the trace center at each x.
Return type:: np.ndarray

pyreduce.trace_model.load_traces(path: str | Path) → tuple[list[Trace], Header][source]¶

Load traces from a FITS file.

Also supports loading legacy NPZ format for backwards compatibility.

Parameters:

path (str | Path) – Input file path (.fits or .npz).

Returns:

traces (list[Trace]) – Loaded traces.
header (fits.Header) – FITS header (empty for NPZ files).

pyreduce.trace_model.save_traces(path: str | Path, traces: list[Trace], header: Header = None, steps: list[str] = None) → None[source]¶

Save traces to a FITS binary table.

Parameters:

path (str | Path) – Output file path.
traces (list[Trace]) – Traces to save.
header (fits.Header, optional) – FITS header to include. If None, a minimal header is created.
steps (list[str], optional) – Pipeline steps that have been run (stored in E_STEPS header).

Raises:

ValueError – If traces have duplicate (group, m) keys.

pyreduce.spectra module¶

Spectrum data model for PyReduce.

This module defines the Spectrum, ExtractionParams, and Spectra classes for storing extracted spectral data.

Replaces the legacy Echelle class with cleaner semantics: - NaN masking instead of COLUMNS + MASK redundancy - Per-trace metadata (m, fiber, extraction_height) - Un-normalized spectra with separate continuum

class pyreduce.spectra.ExtractionParams(osample: int, lambda_sf: float, lambda_sp: float, swath_width: int | None = None)[source]¶

Bases: object

Global extraction parameters (same for all traces in a file).

osample¶

Oversampling factor for slit function.

Type:: int

lambda_sf¶

Slitfunction smoothing parameter.

Type:: float

lambda_sp¶

Spectrum smoothing parameter.

Type:: float

swath_width¶

Swath width for extraction, or None for automatic.

Type:: int | None

classmethod from_header(header: Header) → ExtractionParams | None[source]¶

Read extraction parameters from FITS header.

Parameters:: header (fits.Header) – Header to read from.
Returns:: Extraction parameters, or None if not present.
Return type:: ExtractionParams or None

lambda_sf: float¶

lambda_sp: float¶

osample: int¶

swath_width: int | None = None¶

to_header(header: Header) → None[source]¶

Write extraction parameters to FITS header.

Parameters:: header (fits.Header) – Header to write to.

class pyreduce.spectra.Spectra(header: Header, data: list[Spectrum], params: ExtractionParams | None = None)[source]¶

Bases: object

Container for multiple spectra from one observation.

Replaces the legacy Echelle class.

header¶

FITS header with observation metadata.

Type:: fits.Header

data¶

Individual spectra, one per trace.

Type:: list[Spectrum]

params¶

Global extraction parameters (stored in header).

Type:: ExtractionParams | None

data: list[Spectrum]¶

get_arrays() → dict[str, ndarray][source]¶

Get stacked arrays for all spectra.

Returns:: Dictionary with keys ‘spec’, ‘sig’, ‘wave’, ‘cont’, ‘m’, ‘group’, ‘fiber_idx’. Arrays are stacked along axis 0 (one row per trace).
Return type:: dict

header: Header¶

property ncol: int¶: Number of columns (pixels) in spectra.

property ntrace: int¶: Number of traces/spectra.

params: ExtractionParams | None = None¶

static read(fname: str | Path, raw: bool = False, continuum_normalization: bool = False, barycentric_correction: bool = True, radial_velocity_correction: bool = True) → Spectra[source]¶

Read spectra from a FITS file.

Supports both new format (E_FMTVER >= 2) and legacy format.

Parameters:

fname (str or Path) – Input file path.
raw (bool) – If True, skip all corrections.
continuum_normalization (bool) – Apply continuum normalization (default False for new format).
barycentric_correction (bool) – Apply barycentric correction to wavelength.
radial_velocity_correction (bool) – Apply radial velocity correction to wavelength.

Returns:

Loaded spectra.

Return type:

Spectra

save(fname: str | Path, steps: list[str] = None) → None[source]¶

Save spectra to a FITS file.

Parameters:

fname (str or Path) – Output file path.
steps (list[str], optional) – Pipeline steps that have been run.

select(m: int | None = None, group: str | int | None = None, fiber_idx: int | None = None) → list[Spectrum][source]¶

Filter spectra by order, group, and/or fiber index.

Parameters:

m (int, optional) – Select spectra with this spectral order number.
group (str or int, optional) – Select spectra with this group identifier.
fiber_idx (int, optional) – Select spectra with this fiber index.

Returns:

Matching spectra.

Return type:

list[Spectrum]

Bases: object

Output of extraction for one trace.

m¶

Spectral order number. None if unknown.

Type:: int | None

group¶

Group identifier (e.g., ‘A’, ‘B’, ‘cal’ or bundle index). None if trace has not been through fiber grouping.

Type:: str | int | None

fiber_idx¶

Fiber index within group (1-indexed). None if unknown.

Type:: int | None

spec¶

Flux values, un-normalized. NaN for masked pixels.

Type:: np.ndarray

sig¶

Uncertainty values. NaN for masked pixels.

Type:: np.ndarray

wave¶

Wavelength values (evaluated from polynomial). Same length as spec.

Type:: np.ndarray | None

cont¶

Continuum values (full array, not polynomial). Same length as spec.

Type:: np.ndarray | None

slitfu¶

Slit function (shape depends on osample: height * osample + 1).

Type:: np.ndarray | None

extraction_height¶

Extraction aperture used for this trace.

Type:: float | None

cont: ndarray | None = None¶

extraction_height: float | None = None¶

fiber_idx: int | None = None¶

classmethod from_trace(trace: Trace, spec: np.ndarray, sig: np.ndarray, **kwargs) → Spectrum[source]¶

Factory method that copies identity from Trace.

Parameters:

trace (Trace) – Source trace for identity (m, group, fiber_idx).
spec (np.ndarray) – Extracted spectrum.
sig (np.ndarray) – Spectrum uncertainty.
**kwargs – Additional fields (wave, cont, slitfu, extraction_height).

Returns:

New spectrum with identity copied from trace.

Return type:

Spectrum

group: str | int | None = None¶

m: int | None¶

property mask: ndarray¶: Boolean mask where True indicates invalid (masked) pixels.

normalized() → tuple[ndarray, ndarray][source]¶

Return continuum-normalized spectrum and uncertainty.

Returns:

spec_norm (np.ndarray) – Spectrum divided by continuum.
sig_norm (np.ndarray) – Uncertainty divided by continuum.

Raises:

ValueError – If no continuum is available.

sig: ndarray¶

slitfu: ndarray | None = None¶

spec: ndarray¶

wave: ndarray | None = None¶

pyreduce.spectra.read(fname: str | Path, **kwargs) → Spectra[source]¶

Read spectra from a file.

Alias for Spectra.read().

pyreduce.wavelength_calibration module¶

Wavelength Calibration by comparison to a reference spectrum Loosely bases on the IDL wavecal function

class pyreduce.wavelength_calibration.AlignmentPlot(ax, obs, lines, offset=(0, 0), plot_title=None)[source]¶

Bases: object

Makes a plot which can be clicked to align the two spectra, reference and observed

connect()[source]¶: connect the click event with the appropiate function

make_ref_image()[source]¶: create and show the reference plot, with the two spectra

on_click(event)[source]¶: On click offset the reference by the distance between click positions

class pyreduce.wavelength_calibration.LineAtlas(element, medium='vac', search_dirs=None)[source]¶

Bases: object

load_fits(fname)[source]¶

class pyreduce.wavelength_calibration.LineList(lines=None, obase=None)[source]¶

Bases: object

add_line(wave, order, pos, width, height, flag)[source]¶

append(linelist)[source]¶

dtype = dtype((numpy.record, [(('wlc', 'WLC'), '>f8'), (('wll', 'WLL'), '>f8'), (('posc', 'POSC'), '>f8'), (('posm', 'POSM'), '>f8'), (('xfirst', 'XFIRST'), '>i2'), (('xlast', 'XLAST'), '>i2'), (('approx', 'APPROX'), 'O'), (('width', 'WIDTH'), '>f8'), (('height', 'HEIGHT'), '>f8'), (('order', 'ORDER'), '>i2'), ('flag', '?')]))¶

classmethod from_list(wave, order, pos, width, height, flag)[source]¶

classmethod load(filename)[source]¶

save(filename)[source]¶

class pyreduce.wavelength_calibration.WavelengthCalibration(threshold=100, degree=(6, 6), iterations=3, dimensionality='2D', nstep=0, correlate_cols=0, shift_window=0.01, manual=False, polarim=False, lfc_peak_width=3, closing=5, atlas_name=None, atlas_search_dirs=None, medium='vac', plot=True, plot_title=None, element=None)[source]¶

Bases: object

Wavelength Calibration Module

Takes an observed wavelength image and the reference linelist and returns the wavelength at each pixel

align(obs, lines)[source]¶

Align the observation with the reference spectrum Either automatically using cross correlation or manually (visually)

Parameters:

obs (array[nrow, ncol]) – observed wavelength calibration spectrum (e.g. obs=ThoriumArgon)
lines (struct_array) – reference line data
manual (bool, optional) – wether to manually align the spectra (default: False)
plot (bool, optional) – wether to plot the alignment (default: False)

Returns:

offset – offset in order and column

Return type:

tuple(int, int)

align_manual(obs, lines)[source]¶

Open an AlignmentPlot window for manual selection of the alignment

Parameters:

obs (array of shape (ntrace, ncol)) – observed image
lines (recarray of shape (nlines,)) – reference linelist

Returns:

offset – offset in order and column to be applied to each line in the linelist

Return type:

tuple(int, int)

apply_alignment_offset(lines, offset, select=None)[source]¶

Apply an offset to the linelist

Parameters:

lines (recarray of shape (nlines,)) – reference linelist
offset (tuple(int, int)) – offset in (order, column)
select (array of shape(nlines,), optional) – Mask that defines which lines the offset applies to

Returns:

lines – linelist with offset applied

Return type:

recarray of shape (nlines,)

atlas_name¶

Name of the line atlas used for calibration

Type:: str

atlas_search_dirs¶

Directories to search for atlas files (before the default)

Type:: list

auto_id(obs, wave_img, lines)[source]¶

Automatically identify peaks that are close to known lines

Parameters:

obs (array of shape (ntrace, ncol)) – observed spectrum
wave_img (array of shape (ntrace, ncol)) – wavelength solution image
lines (struc_array) – line data
threshold (int, optional) – difference threshold between line positions in m/s, until which a line is considered identified (default: 1)
plot (bool, optional) – wether to plot the new lines

Returns:

lines – line data with new flags

Return type:

struct_array

build_2d_solution(lines, plot=False)[source]¶

Create a 2D polynomial fit to flagged lines.

Parameters:

lines (struc_array) – line data
plot (bool, optional) – whether to plot the solution (default: False)

Returns:

coef – 2d polynomial coefficients

Return type:

array[degree_x, degree_y]

build_step_solution(lines, plot=False)[source]¶

Fit the least squares fit to the wavelength points, with additional free parameters for detector gaps, e.g. due to stitching.

The exact method of the fit depends on the dimensionality. Either way we are using the usual polynomial fit for the wavelength, but the x points are modified beforehand by shifting them some amount, at specific indices. We assume that the stitching effects are distributed evenly and we know how many steps we expect (this is set as “nstep”).

Parameters:

lines (np.recarray) – linedata
plot (bool, optional) – whether to plot results or not, by default False

Returns:

coefficients of the best fit

Return type:

coef

calculate_AIC(lines, wave_solution)[source]¶

calculate_residual(wave_solution, lines)[source]¶

Calculate all residuals of all given lines

Residual = (Wavelength Solution - Expected Wavelength) / Expected Wavelength * speed of light

Parameters:

wave_solution (array of shape (degree_x, degree_y)) – polynomial coefficients of the wavelength solution (in numpy format)
lines (recarray of shape (nlines,)) – contains the position of the line on the detector (posm), the order (order), and the expected wavelength (wll)

Returns:

residual – Residual of each line in m/s

Return type:

array of shape (nlines,)

closing¶

grey closing range for the input image

Type:: int

correlate_cols¶

How many columns to use in the 2D cross correlation alignment. 0 means all pixels (slow).

Type:: int

create_image_from_lines(lines)[source]¶

Create a reference image based on a line list Each line will be approximated by a Gaussian Space inbetween lines is 0

Parameters:: lines (recarray of shape (nlines,)) – line data
Returns:: img – New reference image
Return type:: array of shape (ntrace, ncol)

degree¶

polynomial degree of the wavelength fit in (pixel, order) direction

Type:: tuple(int, int)

property dimensionality¶

Whether to use 1d or 2d fit

Type:: {“1D”, “2D”}

evaluate_solution(pos, order, solution)[source]¶

Evaluate the 1d or 2d wavelength solution at the given pixel positions and orders

Parameters:

pos (array) – pixel position on the detector (i.e. x axis)
order (array) – order of each point
solution (array of shape (ntrace, ndegree) or (degree_x, degree_y)) – polynomial coefficients. For mode=1D, one set of coefficients per order. For mode=2D, the first dimension is for the positions and the second for the orders
mode (str, optional) – Wether to interpret the solution as 1D or 2D polynomials, by default “1D”

Returns:

result – Evaluated polynomial

Return type:

array

Raises:

ValueError – If pos and order have different shapes, or mode is of the wrong value

evaluate_step_solution(pos, order, solution)[source]¶

execute(obs, lines)[source]¶

Perform the whole wavelength calibration procedure with the current settings

Parameters:

obs (array of shape (ntrace, ncol)) – observed image
lines (recarray of shape (nlines,)) – reference linelist

Returns:

wave_img – Wavelength solution for each pixel

Return type:

array of shape (ntrace, ncol)

Raises:

NotImplementedError – If polarimitry flag is set

f(x, poly_coef, step_coef_pos, step_coef_diff)[source]¶

fit_lines(obs, lines)[source]¶

Determine exact position of each line on the detector based on initial guess

This fits a Gaussian to each line, and uses the peak position as a new solution

Parameters:

obs (array of shape (ntrace, ncol)) – observed wavelength calibration image
lines (recarray of shape (nlines,)) – reference line data

Returns:

lines – Updated line information (posm is changed)

Return type:

recarray of shape (nlines,)

g(x, step_coef_pos, step_coef_diff)[source]¶

iterations¶

Number of iterations in the remove residuals, auto id, loop

Type:: int

lfc_peak_width¶

Laser Frequency Peak width (for scipy.signal.find_peaks)

Type:: int

make_wave(wave_solution, plot=False)[source]¶

Expand polynomial wavelength solution into full image

Parameters:

wave_solution (array of shape(degree,)) – polynomial coefficients of wavelength solution
plot (bool, optional) – wether to plot the solution (default: False)

Returns:

wave_img – wavelength solution for each point in the spectrum

Return type:

array of shape (ntrace, ncol)

manual¶

Whether to manually align the reference instead of using cross correlation

Type:: bool

medium¶

Medium of the detector, vac or air

Type:: str

ncol¶

Number of columns in the observation

Type:: int

normalize(obs, lines)[source]¶

Normalize the observation and reference list in each order individually Copies the data if the image, but not of the linelist

Parameters:

obs (array of shape (ntrace, ncol)) – observed image
lines (recarray of shape (nlines,)) – reference linelist

Returns:

obs (array of shape (ntrace, ncol)) – normalized image
lines (recarray of shape (nlines,)) – normalized reference linelist

nstep¶

Whether to fit for pixel steps (offsets) in the detector

Type:: bool

ntrace¶

Number of orders in the observation

Type:: int

plot¶

Whether to plot the results. Set to 2 to plot during all steps.

Type:: int

plot_residuals(lines, coef, title='Residuals')[source]¶

plot_results(wave_img, obs)[source]¶

polarim¶

Whether to use polarimetric orders instead of the usual ones. I.e. Each pair of two orders represents the same data. Not Supported yet

Type:: bool

reject_lines(lines, plot=False)[source]¶

Reject the largest outlier one by one until all residuals are lower than the threshold

Parameters:

lines (recarray of shape (nlines,)) – Line data with pixel position, and expected wavelength
threshold (float, optional) – upper limit for the residual, by default 100
degree (tuple, optional) – polynomial degree of the wavelength solution (pixel, column) (default: (6, 6))
plot (bool, optional) – Wether to plot the results (default: False)

Returns:

lines – Line data with updated flags

Return type:

recarray of shape (nlines,)

reject_outlier(residual, lines)[source]¶

Reject the strongest outlier

Parameters:

residual (array of shape (nlines,)) – residuals of all lines
lines (recarray of shape (nlines,)) – line data

Returns:

lines (struct_array) – line data with one more flagged line
residual (array of shape (nlines,)) – residuals of each line, with outliers masked (including the new one)

shift_window¶

Fraction if the number of columns to use in the alignment of individual orders. Set to 0 to disable

Type:: float

property step_mode¶

threshold¶

Residual threshold in m/s above which to remove lines

Type:: float

class pyreduce.wavelength_calibration.WavelengthCalibrationComb(threshold=100, degree=(6, 6), iterations=3, dimensionality='2D', nstep=0, correlate_cols=0, shift_window=0.01, manual=False, polarim=False, lfc_peak_width=3, closing=5, atlas_name=None, atlas_search_dirs=None, medium='vac', plot=True, plot_title=None, element=None)[source]¶

Bases: WavelengthCalibration

execute(comb, wave, lines=None)[source]¶

Perform the whole wavelength calibration procedure with the current settings

Parameters:

obs (array of shape (ntrace, ncol)) – observed image
lines (recarray of shape (nlines,)) – reference linelist

Returns:

wave_img – Wavelength solution for each pixel

Return type:

array of shape (ntrace, ncol)

Raises:

NotImplementedError – If polarimitry flag is set

class pyreduce.wavelength_calibration.WavelengthCalibrationInitialize(degree=2, plot=False, plot_title='Wavecal Initial', resid_delta=2000, match_tolerance=1.0, iterations=3, edge_margin=10, width_min=1.0, width_max=8.0, cutoff=0.01, smoothing=0, atlas_name='thar', atlas_search_dirs=None, medium='vac', element=None)[source]¶

Bases: WavelengthCalibration

execute(spectrum, wave_range) → LineList[source]¶

Perform the whole wavelength calibration procedure with the current settings

Parameters:

obs (array of shape (ntrace, ncol)) – observed image
lines (recarray of shape (nlines,)) – reference linelist

Returns:

wave_img – Wavelength solution for each pixel

Return type:

array of shape (ntrace, ncol)

Raises:

NotImplementedError – If polarimitry flag is set

get_cutoff(spectrum)[source]¶

identify_lines_for_order(spectrum, atlas, wave_range, order) → LineList[source]¶

Identify spectral lines via iterative peak matching (IDL algorithm).

Detects peaks in the observed spectrum, matches them to atlas lines using a polynomial wavelength solution, and iteratively refines the fit with outlier rejection.

Parameters:

spectrum (array) – observed spectrum for one order
atlas (LineAtlas) – reference line atlas
wave_range (2-tuple) – initial wavelength guess (begin, end) in Angstrom
order (int) – order index

Returns:

matched lines for this order

Return type:

LineList

normalize(spectrum)[source]¶

Normalize the observation and reference list in each order individually Copies the data if the image, but not of the linelist

Parameters:

obs (array of shape (ntrace, ncol)) – observed image
lines (recarray of shape (nlines,)) – reference linelist

Returns:

obs (array of shape (ntrace, ncol)) – normalized image
lines (recarray of shape (nlines,)) – normalized reference linelist

pyreduce.wavelength_calibration.polyfit(x, y, deg)[source]¶

pyreduce.slit_curve module¶

Calculate slit curvature based on a reference spectrum with high SNR, e.g. Wavelength calibration image

Authors¶

Nikolai Piskunov Ansgar Wehrhahn

Version¶

0.9 - NP - IDL Version 1.0 - AW - Python Version

License¶

class pyreduce.slit_curve.Curvature(traces, curve_height=0.5, extraction_height=0.2, trace_range=None, window_width=9, peak_threshold=10, peak_width=1, fit_degree=2, sigma_cutoff=3, mode='1D', plot=False, plot_title=None, peak_function='gaussian', curve_degree=2)[source]¶

Bases: object

eval(peaks, order, fitted_coeffs)[source]¶

Evaluate fitted curvature coefficients at given positions.

Parameters:

peaks (array) – Column positions to evaluate at
order (array) – Order indices (same shape as peaks)
fitted_coeffs (array) – Fitted coefficients from fit() method

Returns:

coeffs – Evaluated curvature coefficients at each position

Return type:

array of shape (len(peaks), curve_degree)

eval_legacy(peaks, order, fitted_coeffs)[source]¶: Evaluate and return legacy p1, p2 format for backward compatibility.

execute(original, compute_slitdeltas=True)[source]¶

Execute curvature determination using row-tracking method.

Parameters:

original (array of shape (nrow, ncol)) – Original image
compute_slitdeltas (bool) – Whether to compute slitdeltas from fit residuals (default: True)

Returns:

Curvature data with keys: coeffs, slitdeltas, degree, fitted_coeffs, fit_degree.

Return type:

dict

fit(peaks, all_coeffs)[source]¶

Fit smooth polynomial to curvature coefficients.

Parameters:

peaks (list of arrays) – Peak columns for each order
all_coeffs (list of arrays) – Curvature coefficients for each order. Each entry has shape (n_peaks, curve_degree).

Returns:

fitted_coeffs – For 1D mode: shape (n_orders, curve_degree, fit_degree + 1) For 2D mode: shape (curve_degree, …) with polyfit2d coefficients

Return type:

array

property mode¶

property n¶

property ntrace¶

plot_comparison(original, coeffs_array, peaks, slitdeltas=None)[source]¶

Plot comparison of curvature model vs data.

Parameters:

original (array) – Original image
coeffs_array (array of shape (ntrace, ncol, curve_degree + 1)) – Curvature coefficients at each point
peaks (list of arrays) – Peak columns for each order
slitdeltas (array of shape (ntrace, nrows), optional) – Per-row residual offsets. If provided, plotted as white lines offset from the polynomial (red lines).

plot_results(ncol, plot_peaks, plot_vec, plot_coeffs, fitted_coeffs)[source]¶

Plot curvature fitting results.

Parameters:

ncol (int) – Number of columns in image
plot_peaks (list of arrays) – Peak columns for each order
plot_vec (list of arrays) – Middle spectrum for each order
plot_coeffs (list of arrays) – Raw curvature coefficients for each order
fitted_coeffs (array) – Fitted coefficients from fit() method

pyreduce.slit_curve.gaussian(x, A, mu, sig)[source]¶: A: height mu: offset from central line sig: standard deviation

pyreduce.slit_curve.lorentzian(x, A, x0, mu)[source]¶: A: height x0: offset from central line mu: width of lorentzian

pyreduce.combine_frames module¶

Combine several fits files into one master frame

Used to create master bias and master flat

pyreduce.combine_frames.calculate_probability(buffer, window, method='sum')[source]¶

Construct a probability function based on buffer data.

Parameters:

buffer (array of shape (nx, ny)) – buffer
window (int) – size of the running window
method ({"sum", "median"}, optional) – which method to use to average the probabilities (default: “sum”) “sum” is much faster, but “median” is more resistant to outliers

Returns:

weights – probabilities

Return type:

array of shape (nx, ny - 2 * window)

pyreduce.combine_frames.combine_bias(files, instrument, channel, extension=None, plot=False, plot_title=None, science_observation_time=None, **kwargs)[source]¶

Combine bias frames, determine read noise, reject bad pixels. Read noise calculation only valid if both lists yield similar noise.

Parameters:

files (list(str)) – bias files to combine
instrument (str) – instrument channel for arminfo
extension ({int, str}, optional) – fits extension to use (default: 1)
xr (2-tuple(int), optional) – x range to use (default: None, i.e. whole image)
yr (2-tuple(int), optional) – y range to use (default: None, i.e. whole image)
dtype (np.dtype, optional) – datatype of the combined bias frame (default: float32)

Returns:

bias image and header

Return type:

bias, bhead

pyreduce.combine_frames.combine_calibrate(files, instrument, channel, mask=None, bias=None, bhead=None, norm=None, bias_scaling='exposure_time', norm_scaling='divide', plot=False, plot_title=None, traces=None, extraction_height=None, **kwargs)[source]¶

Combine the input files and then calibrate the image with the bias and normalized flat field if provided

Parameters:

files (list) – list of file names to load
instrument (Instrument) – PyReduce instrument object with load_fits method
channel (str) – descriptor of the instrument channel
mask (array) – 2D Bad Pixel Mask to apply to the master image
bias (tuple(bias, bhead), optional) – bias correction to apply to the combiend image, if bias has 3 dimensions it is used as polynomial coefficients scaling with the exposure time, by default None
norm_flat (tuple(norm, blaze), optional) – normalized flat to divide the combined image with after the bias subtraction, by default None
bias_scaling (str, optional) – defines how the bias is subtracted, by default “exposure_time”
plot (bool, optional) – whether to plot the results, by default False
plot_title (str, optional) – Name to put on the plot, by default None

Returns:

orig (array) – combined image with calibrations applied
thead (Header) – header of the combined image

Raises:

ValueError – Unrecognised bias_scaling option

pyreduce.combine_frames.combine_frames(files, instrument, channel, extension=None, threshold=3.5, window=50, dtype=<class 'numpy.float32'>, **kwargs)[source]¶

Subroutine to correct cosmic rays blemishes, while adding otherwise similar images.

combine_frames co-adds a group of FITS files with 2D images of identical dimensions. In the process it rejects cosmic ray, detector defects etc. It is capable of handling images that have strip pattern (e.g. echelle spectra) using the REDUCE arminfo conventions to figure out image orientation and useful pixel ranges. It can handle many frames. Special cases: 1 file in the list (the input is returned as output) and 2 files (straight sum is returned).

If the image orientation is not predominantly vertical, the image is rotated 90 degrees (and rotated back afterwards).

Open all FITS files in the list. Loop through the rows. Read next row from each file into a row buffer mBuff[nCol, nFil]. Optionally correct the data for non-linearity.

calc_probability:

Go through the row creating "probability" vector. That is for column iCol take the median of
the part of the row mBuff[iCol-win:iCol+win,iFil] for each file and divide these medians by the
mean of them computed across the stack of files. In other words:
>>> filwt[iFil] = median(mBuff[iCol-win:iCol+win,iFil])
>>> norm_filwt = mean(filwt)
>>> prob[iCol,iFil] = (norm_filtwt>0)?filwt[iCol]/norm_filwt:filwt[iCol]

This is done for all iCol in the range of [win:nCol-win-1]. It is then linearly extrapolated to
the win zones of both ends. E.g. for iCol in [0:win-1] range:
>>> prob[iCol,iFil]=2*prob[win,iFil]-prob[2*win-iCol,iFil]

For the other end ([nCol-win:nCol-1]) it is similar:
>>> prob[iCol,iFil]=2*prob[nCol-win-1,iFil]-prob[2*(nCol-win-1)-iCol,iFil]

fix_bad_pixels:

Once the probailities are constructed we can do the fitting, measure scatter and detect outliers.
We ignore negative or zero probabilities as it should not happen. For each iCol with (some)
positive probabilities we compute tha ratios of the original data to the probabilities and get
the mean amplitude of these ratios after rejecting extreme values:
>>> ratio = mBuff[iCol,iFil]/prob[iCol,iFil]
>>> amp = (total(ratio)-min(ratio)-max(ratio))/(nFil-2)
>>> mFit[iCol,iFil] = amp*prob[iCol,iFil]

Note that for iFil whereprob[iCol,iFil] is zero we simply set mFit to zero. The scatter (noise)
consists readout noise and shot noise of the model (fit) co-added in quadratures:
>>> sig=sqrt(rdnoise*rdnoise + abs(mFit[iCol,iFil]/gain))

and the outliers are defined as:
>>> iBad=where(mBuff-mFit gt thresh*sig)

>>> Bad values are replaced from the fit:
>>> mBuff[iBad]=mFit[iBad]

and mBuff is summed across the file dimension to create an output row.

Parameters:

files (list(str)) – list of fits files to combine
instrument (str) – instrument id for arminfo
channel (str) – instrument channel
extension (int, optional) – fits extension to load (default: 1)
threshold (float, optional) – threshold for bad pixels (default: 3.5)
window (int, optional) – horizontal window size (default: 50)
mask (array(bool), optional) – mask for the fits image (default: None)
xr (int, optional) – xrange (default: None)
yr (int, optional) – yrange (default: None)
debug (bool, optional) – show debug plot of noise distribution (default: False)
dtype (np.dtype, optional) – datatype of the combined image (default float32)

Returns:

combined image data, header

Return type:

combined_data, header

pyreduce.combine_frames.combine_frames_simple(files, instrument, channel, extension=None, dtype=<class 'numpy.float32'>, **kwargs)[source]¶

Simple addition of similar images.

Parameters:

files (list(str)) – list of fits files to combine
instrument (str) – instrument id for modinfo
channel (str) – instrument channel
extension (int, optional) – fits extension to load (default: 1)
dtype (np.dtype, optional) – datatype of the combined image (default float32)

Returns:

combined image data, header

Return type:

combined_data, header

pyreduce.combine_frames.combine_polynomial(files, instrument, channel, mask, degree=1, plot=False, plot_title=None)[source]¶

Combine the input files by fitting a polynomial of the pixel value versus the exposure time of each pixel

Parameters:

files (list) – list of file names
instrument (Instrument) – PyReduce instrument object with load_fits method
channel (str) – channel identifier for this instrument
mask (array) – bad pixel mask to apply to the coefficients
degree (int, optional) – polynomial degree of the fit, by default 1
plot (bool, optional) – whether to plot the results, by default False
plot_title (str, optional) – Title of the plot, by default None

Returns:

bias (array) – 3d array with the coefficients for each pixel
bhead (Header) – combined FITS header of the coefficients

pyreduce.combine_frames.fix_bad_pixels(probability, buffer, readnoise, gain, threshold)[source]¶

find and fix bad pixels

Parameters:

probability (array(float)) – probabilities
buffer (array(int)) – image buffer
readnoise (float) – readnoise of current amplifier
gain (float) – gain of current amplifier
threshold (float) – sigma threshold between observation and fit for bad pixels

Returns:

input buffer, with bad pixels fixed

Return type:

array(int)

pyreduce.combine_frames.running_median(arr, size)[source]¶

Calculate the running median of a 2D sequence

Parameters:

seq (2d array [n, l]) – n datasets of length l
size (int) – number of elements to consider for each median

Returns:

running median

Return type:

2d array [n, l-size]

pyreduce.combine_frames.running_sum(arr, size)[source]¶

Calculate the running sum over the 2D sequence

Parameters:

arr (array[n, l]) – sequence to calculate running sum over, n datasets of length l
size (int) – number of elements to sum

Returns:

running sum

Return type:

2D array

pyreduce.continuum_normalization module¶

Find the continuum level

Currently only splices orders together First guess of the continuum is provided by the flat field

class pyreduce.continuum_normalization.Plot_Normalization(wsort, sB, new_wave, contB, iteration=0, title=None)[source]¶

Bases: object

close()[source]¶

plot(wsort, sB, new_wave, contB, iteration)[source]¶

pyreduce.continuum_normalization.continuum_normalize(spec, wave, cont, sigm, iterations=10, smooth_initial=100000.0, smooth_final=5000000.0, scale_vert=1, plot=True, plot_title=None)[source]¶

Fit a continuum to a spectrum by slowly approaching it from the top. We exploit here that the continuum varies only on large wavelength scales, while individual lines act on much smaller scales

TODO automatically find good parameters for smooth_initial and smooth_final TODO give variables better names

Parameters:

spec (masked array of shape (ntrace, ncol)) – Observed input spectrum, masked values describe column ranges
wave (masked array of shape (ntrace, ncol)) – Wavelength solution of the spectrum
cont (masked array of shape (ntrace, ncol)) – Initial continuum guess, for example based on the blaze
sigm (masked array of shape (ntrace, ncol)) – Uncertainties of the spectrum
iterations (int, optional) – Number of iterations of the algorithm, note that runtime roughly scales with the number of iterations squared (default: 10)
smooth_initial (float, optional) – Smoothing parameter in the initial runs, usually smaller than smooth_final (default: 1e5)
smooth_final (float, optional) – Smoothing parameter of the final run (default: 5e6)
scale_vert (float, optional) – Vertical scale of the spectrum. Usually 1 if a previous normalization exists (default: 1)
plot (bool, optional) – Wether to plot the current status and results or not (default: True)

Returns:

cont – New continuum

Return type:

masked array of shape (ntrace, ncol)

pyreduce.continuum_normalization.splice_orders(spec, wave, cont, sigm, scaling=True, plot=False, plot_title=None)[source]¶

Splice orders together so that they form a continous spectrum This is achieved by linearly combining the overlaping regions

Parameters:

spec (array[ntrace, ncol]) – Spectrum to splice, with seperate orders
wave (array[ntrace, ncol]) – Wavelength solution for each point
cont (array[ntrace, ncol]) – Continuum, blaze function will do fine as well
sigm (array[ntrace, ncol]) – Errors on the spectrum
scaling (bool, optional) – If true, the spectrum/continuum will be scaled to 1 (default: False)
plot (bool, optional) – If true, will plot the spliced spectrum (default: False)

Raises:

NotImplementedError – If neighbouring orders dont overlap

Returns:

spec, wave, cont, sigm – spliced spectrum

Return type:

array[ntrace, ncol]

pyreduce.estimate_background_scatter module¶

Module that estimates the background scatter

pyreduce.estimate_background_scatter.estimate_background_scatter(img, traces, extraction_height=0.1, scatter_degree=4, sigma_cutoff=2, border_width=10, plot=False, plot_title=None)[source]¶

Estimate the background by fitting a 2d polynomial to inter-trace data

Inter-trace data is all pixels minus the traces +- the extraction width

Parameters:

img (array[nrow, ncol]) – (flat) image data
traces (list[Trace]) – Trace objects with pos, column_range attributes
extraction_height (float, optional) – extraction full height, values below 2 are considered fractional (default: 0.1)
scatter_degree (int, optional) – polynomial degree of the 2d fit for the background scatter (default: 4)
plot (bool, optional) – wether to plot the fitted polynomial and the data or not (default: False)

Returns:

coeff – 2D polynomial coefficients for background scatter

Return type:

array

pyreduce.rectify module¶

pyreduce.rectify.merge_images(images, wave, column_range, extraction_height)[source]¶

pyreduce.rectify.rectify_image(img, traces, extraction_height, trace_range)[source]¶

Rectify image by extracting and straightening each trace.

Parameters:

img (array) – Input image
traces (list[Trace]) – Trace objects with pos, column_range, and optional slit curvature
extraction_height (float) – Extraction height (fraction if < 3, else pixels)
trace_range (tuple) – (start, end) indices of traces to process

Returns:

images (dict) – Rectified images keyed by trace index
column_range (array) – Column ranges for each trace
extraction_height (array) – Extraction heights for each trace

pyreduce.clipnflip module¶

Module that

Clips remove pre- and overscan regions
Flips orients the image so that orders are roughly horizontal

pyreduce.clipnflip.clipnflip(image, header, xrange=None, yrange=None, orientation=None, transpose=None)[source]¶

Process an image and associated FITS header already in memory as follows: 1. Trim image to desired subregion: newimage = image(xlo:xhi,ylo:yhi) 2. Transform to standard orientation (red at top, orders run left to right)

Parameters:

image (array[nrow, ncol]) – raw image to be processed.
header (fits.header, dict) – FITS header of the image
xrange (tuple(int, int), optional) – column - range to keep in the image (default: data from header/instrument)
yrange (tuple(int, int), optional) – row - range to keep in the image (default: data from header/instrument)
orientation (int, optional) – number of counterclockwise 90 degrees rotation to apply to the image (default: data from header/instrument)
transpose (bool, optional) – if True the image will be transposed before rotation (default: data from header/instrument)

Returns:

image – clipped and flipped image

Return type:

array[yrange, xrange]

Raises:

NotImplementedError – nonlinear images are not supported yet

pyreduce.cwrappers module¶

Wrapper for REDUCE C functions

This module provides access to the extraction algorithms in the C libraries and sanitizes the input parameters.

Mask convention:

FITS files and Python: numpy convention (1/True = bad/masked pixel)
C code uses REDUCE convention: 1 = good pixel, 0 = bad pixel

All conversions happen at the boundary in this module.

pyreduce.cwrappers.create_spectral_model(ncols: int, nrows: int, osample: int, xi: [('x', <class 'ctypes.c_int'>), ('y', <class 'ctypes.c_int'>), ('w', <class 'ctypes.c_double'>)], spec: ~numpy.ndarray, slitfunc: ~numpy.ndarray)[source]¶

pyreduce.cwrappers.extract_with_slitfunc(img: ndarray, ycen: ndarray, slitfunc: ndarray, slitfunc_meta: dict, yrange: tuple[int, int], osample: int, p1: ndarray | float = 0, p2: ndarray | float = 0, lambda_sp: float = 0, gain: float = 1, maxiter: int = 1, reject_threshold: float = 0) → tuple[ndarray, ndarray, ndarray, ndarray, ndarray, ndarray][source]¶

Extract spectrum using a preset slit function (single-pass extraction).

This function validates and adapts the slit function to match the current extraction parameters, then calls the C extraction code with the preset slit function (skipping the sL solve step).

Parameters:

img (array[nrows, ncols]) – Image to extract from
ycen (array[ncols]) – Order center positions (fractional pixel)
slitfunc (array[ny_src]) – Preset slit function from a previous extraction (e.g., from norm_flat)
slitfunc_meta (dict) – Metadata about the slit function source, must contain: - “osample”: oversampling factor used when slitfunc was computed - “extraction_height”: extraction full height used
yrange (tuple(int, int)) – Target extraction bounds: pixels (below, above) the trace center
osample (int) – Target oversampling factor for slit function
p1 (array[ncols] or float) – Linear curvature coefficient
p2 (array[ncols] or float) – Quadratic curvature coefficient
lambda_sp (float) – Spectrum smoothing parameter (0 = no smoothing)
gain (float) – Detector gain for uncertainty estimation
maxiter (int) – Maximum iterations (default 1 for single-pass)
reject_threshold (float) – Outlier rejection threshold (default 0 = disabled)

Returns:

Same as slitfunc_curved()

Return type:

sp, sl, model, unc, mask, info

Raises:

ValueError – If extraction height is larger than source (can’t extrapolate)

pyreduce.cwrappers.slitfunc(img, ycen, lambda_sp=0, lambda_sf=0.1, osample=1)[source]¶

Decompose image into spectrum and slitfunction

This is for horizontal straight orders only, for curved orders use slitfunc_curved instead

Parameters:

img (array[n, m]) – image to decompose, should just contain a small part of the overall image
ycen (array[n]) – traces the center of the order along the image, relative to the center of the image?
lambda_sp (float, optional) – smoothing parameter of the spectrum (the default is 0, which no smoothing)
lambda_sf (float, optional) – smoothing parameter of the slitfunction (the default is 0.1, which )
osample (int, optional) – Subpixel ovsersampling factor (the default is 1, which no oversampling)

Returns:

spectrum, slitfunction, model, spectrum uncertainties

Return type:

sp, sl, model, unc

pyreduce.cwrappers.slitfunc_curved(img, ycen, p1, p2, lambda_sp, lambda_sf, osample, yrange, maxiter=20, gain=1, reject_threshold=6, preset_slitfunc=None)[source]¶

Decompose an image into a spectrum and a slitfunction, image may be curved

Parameters:

img (array[n, m]) – input image
ycen (array[n]) – traces the center of the order
p1 (array[n]) – 1st order curvature of the order along the image, set to 0 if order straight
p2 (array[n]) – 2nd order curvature of the order along the image, set to 0 if order straight
osample (int) – Subpixel ovsersampling factor (the default is 1, no oversampling)
lambda_sp (float) – smoothing factor spectrum (the default is 0, no smoothing)
lambda_sl (float) – smoothing factor slitfunction (the default is 0.1, small smoothing)
yrange (array[2]) – number of pixels below and above the central line that have been cut out
maxiter (int, optional) – maximumim number of iterations, by default 20
gain (float, optional) – gain of the image, by default 1
reject_threshold (float, optional) – outlier rejection threshold in sigma, by default 6. Set to 0 to disable.
preset_slitfunc (array[ny], optional) – If provided, use this slit function instead of solving for it. Size must be osample * (nrows + 1) + 1.

Returns:

spectrum, slitfunction, model, spectrum uncertainties

Return type:

sp, sl, model, unc

pyreduce.cwrappers.xi_zeta_tensors(ncols: int, nrows: int, ycen: ndarray, yrange, osample: int, p1: ndarray, p2: ndarray)[source]¶

pyreduce.datasets module¶

Download example datasets for PyReduce.

Downloads tarballs from the PyReduce server and extracts them to $REDUCE_DATA (or ~/REDUCE_DATA by default).

pyreduce.datasets.HARPS(local_dir=None)[source]¶: Download HARPS example dataset (target: HD109200).

pyreduce.datasets.JWST_MIRI(local_dir=None)[source]¶: Download JWST/MIRI example dataset (simulated with MIRIsim).

pyreduce.datasets.JWST_NIRISS(local_dir=None)[source]¶: Download JWST/NIRISS example dataset (simulated with awesimsoss).

pyreduce.datasets.KECK_NIRSPEC(local_dir=None)[source]¶: Download Keck/NIRSPEC example dataset (target: GJ1214).

pyreduce.datasets.LICK_APF(local_dir=None)[source]¶: Download Lick APF example dataset (target: KIC05005618).

pyreduce.datasets.MCDONALD(local_dir=None)[source]¶: Download McDonald Observatory example dataset.

pyreduce.datasets.UVES(local_dir=None)[source]¶: Download UVES example dataset (target: HD132205).

pyreduce.datasets.XSHOOTER(local_dir=None)[source]¶: Download XSHOOTER example dataset (target: Ux-Ori).

pyreduce.datasets.get_data_dir()[source]¶

Get the default data directory.

Returns $REDUCE_DATA if set, otherwise ~/REDUCE_DATA

pyreduce.datasets.get_dataset(name, local_dir=None)[source]¶

Download and extract a dataset.

Parameters:

name (str) – Name of the dataset (e.g., “UVES”, “HARPS”)
local_dir (str, optional) – Directory to save data (default: $REDUCE_DATA or ~/REDUCE_DATA)

Returns:

Directory where the data was saved

Return type:

str

pyreduce.echelle module¶

Deprecated: use pyreduce.spectra instead.

This module is provided for backward compatibility only. All functionality has been moved to pyreduce.spectra.

pyreduce.util module¶

Collection of various useful and/or reoccuring functions across PyReduce

pyreduce.util.air2vac(wl_air)[source]¶: Convert wavelengths in air to vacuum wavelength Author: Nikolai Piskunov

pyreduce.util.bezier_interp(x_old, y_old, x_new)[source]¶

Bezier interpolation, based on the scipy methods

This mostly sanitizes the input by removing masked values and duplicate entries Note that in case of duplicate entries (in x_old) the results are not well defined as only one of the entries is used and the other is discarded

Parameters:

x_old (array[n]) – old x values
y_old (array[n]) – old y values
x_new (array[m]) – new x values

Returns:

y_new – new y values

Return type:

array[m]

pyreduce.util.bottom(f, order=1, iterations=40, eps=0.001, poly=False, weight=1, **kwargs)[source]¶

bottom tries to fit a smooth curve to the lower envelope of 1D data array f. Filter size “filter” together with the total number of iterations determine the smoothness and the quality of the fit. The total number of iterations can be controlled by limiting the maximum number of iterations (iter) and/or by setting the convergence criterion for the fit (eps) 04-Nov-2000 N.Piskunov wrote. 09-Nov-2011 NP added weights and 2nd derivative constraint as LAM2

Parameters:

f (Callable) – Function to fit
filter (int) – Smoothing parameter of the optimal filter (or polynomial degree of poly is True)
iter (int) – maximum number of iterations [def: 40]
eps (float) – convergence level [def: 0.001]
mn (float) – minimum function values to be considered [def: min(f)]
mx (float) – maximum function values to be considered [def: max(f)]
lam2 (float) – constraint on 2nd derivative
weight (array(float)) – vector of weights.

pyreduce.util.cutout_image(img, ymin, ymax, xmin, xmax)[source]¶

Cut a section of an image out

Parameters:

img (array) – image
ymin (array[ncol](int)) – lower y value
ymax (array[ncol](int)) – upper y value
xmin (int) – lower x value
xmax (int) – upper x value

Returns:

cutout – selection of the image

Return type:

array[height, ncol]

pyreduce.util.find_first_index(arr, value)[source]¶: find the first element equal to value in the array arr

pyreduce.util.gaussbroad(x, y, hwhm)[source]¶

Apply gaussian broadening to x, y data with half width half maximum hwhm

Parameters:

x (array(float)) – x values
y (array(float)) – y values
hwhm (float > 0) – half width half maximum

Returns:

broadened y values

Return type:

array(float)

pyreduce.util.gaussfit(x, y)[source]¶

Fit a simple gaussian to data

gauss(x, a, mu, sigma) = a * exp(-z**2/2) with z = (x - mu) / sigma

Parameters:

x (array(float)) – x values
y (array(float)) – y values

Returns:

fitted values for x, fit paramters (a, mu, sigma)

Return type:

gauss(x), parameters

pyreduce.util.gaussfit2(x, y)[source]¶

Fit a gaussian(normal) curve to data x, y

gauss = A * exp(-(x-mu)**2/(2*sig**2)) + offset

Parameters:

x (array[n]) – x values
y (array[n]) – y values

Returns:

popt – coefficients of the gaussian: A, mu, sigma**2, offset

Return type:

array[4]

pyreduce.util.gaussfit3(x, y)[source]¶

A very simple (and relatively fast) gaussian fit gauss = A * exp(-(x-mu)**2/(2*sig**2)) + offset

Parameters:

x (array of shape (n,)) – x data
y (array of shape (n,)) – y data

Returns:

popt – Parameters A, mu, sigma**2, offset

Return type:

list of shape (4,)

pyreduce.util.gaussval2(x, a, mu, sig, const)[source]¶

pyreduce.util.helcorr(obs_long, obs_lat, obs_alt, ra2000, dec2000, jd, system='barycentric')[source]¶

calculates heliocentric Julian date, barycentric and heliocentric radial velocity corrections, using astropy functions

Parameters:

obs_long (float) – Longitude of observatory (degrees, western direction is positive)
obs_lat (float) – Latitude of observatory (degrees)
obs_alt (float) – Altitude of observatory (meters)
ra2000 (float) – Right ascension of object for epoch 2000.0 (hours)
dec2000 (float) – Declination of object for epoch 2000.0 (degrees)
jd (float) – Julian date for the middle of exposure in MJD
system ({"barycentric", "heliocentric"}, optional) – reference system of the result, barycentric: around earth-sun gravity center, heliocentric: around sun, usually barycentric is preferred (default: “barycentric)

Returns:

correction (float) – radial velocity correction due to barycentre offset
hjd (float) – Heliocentric Julian date for middle of exposure

pyreduce.util.interpolate_masked(masked)[source]¶

Interpolate masked values, from non masked values

Parameters:: masked (masked_array) – masked array to interpolate on
Returns:: interpolated – interpolated non masked array
Return type:: array

pyreduce.util.is_interactive_plot_mode()[source]¶: Return True if plots are shown interactively (block mode).

pyreduce.util.log_version()[source]¶: For Debug purposes

pyreduce.util.make_index(ymin, ymax, xmin, xmax, zero=0)[source]¶

Create an index (numpy style) that will select part of an image with changing position but fixed height

The user is responsible for making sure the height is constant, otherwise it will still work, but the subsection will not have the desired format

Parameters:

ymin (array[ncol](int)) – lower y border
ymax (array[ncol](int)) – upper y border
xmin (int) – leftmost column
xmax (int) – rightmost colum
zero (bool, optional) – if True count y array from 0 instead of xmin (default: False)

Returns:

index – numpy index for the selection of a subsection of an image

Return type:

tuple(array[height, width], array[height, width])

pyreduce.util.middle(f, param, x=None, iterations=40, eps=0.001, poly=False, weight=1, lambda2=-1, mn=None, mx=None)[source]¶

middle tries to fit a smooth curve that is located along the “middle” of 1D data array f. Filter size “filter” together with the total number of iterations determine the smoothness and the quality of the fit. The total number of iterations can be controlled by limiting the maximum number of iterations (iter) and/or by setting the convergence criterion for the fit (eps) 04-Nov-2000 N.Piskunov wrote. 09-Nov-2011 NP added weights and 2nd derivative constraint as LAM2

Parameters:

f (Callable) – Function to fit
filter (int) – Smoothing parameter of the optimal filter (or polynomial degree of poly is True)
iter (int) – maximum number of iterations [def: 40]
eps (float) – convergence level [def: 0.001]
mn (float) – minimum function values to be considered [def: min(f)]
mx (float) – maximum function values to be considered [def: max(f)]
lam2 (float) – constraint on 2nd derivative
weight (array(float)) – vector of weights.

pyreduce.util.opt_filter(y, par, par1=None, weight=None, lambda2=-1, maxiter=100)[source]¶

Optimal filtering of 1D and 2D arrays. Uses tridiag in 1D case and sprsin and linbcg in 2D case. Written by N.Piskunov 8-May-2000

Parameters:

f (array) – 1d or 2d array
xwidth (int) – filter width (for 2d array width in x direction (1st index)
ywidth (int) – (for 2d array only) filter width in y direction (2nd index) if ywidth is missing for 2d array, it set equal to xwidth
weight (array(float)) – an array of the same size(s) as f containing values between 0 and 1
lambda1 (float) – regularization parameter
maxiter (int) – maximum number of iteration for filtering of 2d array

pyreduce.util.plot2d(x, y, z, coeff, title=None)[source]¶

pyreduce.util.plot_traces(im, traces, column_range=None, ax=None, imshow_kwargs=None, **line_kwargs)[source]¶

Plot image with polynomial traces overlaid.

Parameters:

im (array[nrow, ncol]) – 2D image to display
traces (array or list) – Polynomial coefficients for traces. Either a 2D array of shape (n_traces, degree+1) or a list of 1D coefficient arrays. Coefficients are in numpy polyval order (highest degree first).
column_range (array[n_traces, 2], optional) – Column range [start, end] for each trace. If None, traces span full width.
ax (matplotlib.axes.Axes, optional) – Axes to plot on. If None, creates new figure.
imshow_kwargs (dict, optional) – Keyword arguments passed to imshow (e.g., vmin, vmax, cmap)
**line_kwargs – Additional keyword arguments passed to plot for trace lines (e.g., color, linewidth, alpha)

Returns:

ax – The axes with the plot

Return type:

matplotlib.axes.Axes

pyreduce.util.polyfit1d(x, y, degree=1, regularization=0)[source]¶

pyreduce.util.polyfit2d(x, y, z, degree=1, max_degree=None, scale=True, plot=False, plot_title=None)[source]¶

A simple 2D plynomial fit to data x, y, z The polynomial can be evaluated with numpy.polynomial.polynomial.polyval2d

Parameters:

x (array[n]) – x coordinates
y (array[n]) – y coordinates
z (array[n]) – data values
degree (int, optional) – degree of the polynomial fit (default: 1)
max_degree ({int, None}, optional) – if given the maximum combined degree of the coefficients is limited to this value
scale (bool, optional) – Wether to scale the input arrays x and y to mean 0 and variance 1, to avoid numerical overflows. Especially useful at higher degrees. (default: True)
plot (bool, optional) – wether to plot the fitted surface and data (slow) (default: False)

Returns:

coeff – the polynomial coefficients in numpy 2d format, i.e. coeff[i, j] for x**i * y**j

Return type:

array[degree+1, degree+1]

pyreduce.util.polyfit2d_2(x, y, z, degree=1, x0=None, loss='arctan', method='trf', plot=False)[source]¶

pyreduce.util.polyscale2d(coeff, scale_x, scale_y, copy=True)[source]¶

pyreduce.util.polyshift2d(coeff, offset_x, offset_y, copy=True)[source]¶

pyreduce.util.polyvander2d(x, y, degree)[source]¶

pyreduce.util.remove_bias(img, ihead, bias, bhead, nfiles=1)[source]¶

pyreduce.util.resample(array, new_size)[source]¶

pyreduce.util.safe_interpolation(x_old, y_old, x_new=None, fill_value=0)[source]¶

‘Safe’ interpolation method that should avoid the common pitfalls of spline interpolation

masked arrays are compressed, i.e. only non masked entries are used remove NaN input in x_old and y_old only unique x values are used, corresponding y values are ‘random’ if all else fails, revert to linear interpolation

Parameters:

x_old (array of size (n,)) – x values of the data
y_old (array of size (n,)) – y values of the data
x_new (array of size (m, ) or None, optional) – x values of the interpolated values if None will return the interpolator object (default: None)

Returns:

y_new – if x_new was given, return the interpolated values otherwise return the interpolator object

Return type:

array of size (m, ) or interpolator

pyreduce.util.set_plot_dir(path)[source]¶: Set directory for saving plots. If None, plots will not be saved.

pyreduce.util.set_plot_show(mode, plot_level=0)[source]¶

Set plot display mode: ‘off’, ‘block’, or ‘defer’.

Parameters:

mode (str) – One of ‘off’, ‘block’, or ‘defer’
plot_level (int) – Current plot level. If >= 2, warns that progress plots won’t work with defer/off.

pyreduce.util.show_all()[source]¶: Show all deferred plots. Call at end of pipeline when using defer mode.

pyreduce.util.show_or_save(name='plot')[source]¶

Save plot to file and/or show interactively based on settings.

Parameters:: name (str) – Base name for the saved file (without extension)

pyreduce.util.start_logging(log_file='log.log')[source]¶

Start logging to log file and command line

Parameters:: log_file (str, optional) – name of the logging file (default: “log.log”)

pyreduce.util.swap_extension(fname, ext, path=None)[source]¶: exchange the extension of the given file with a new one

pyreduce.util.top(f, order=1, iterations=40, eps=0.001, poly=False, weight=1, lambda2=-1, mn=None, mx=None)[source]¶

top tries to fit a smooth curve to the upper envelope of 1D data array f. Filter size “filter” together with the total number of iterations determine the smoothness and the quality of the fit. The total number of iterations can be controlled by limiting the maximum number of iterations (iter) and/or by setting the convergence criterion for the fit (eps) 04-Nov-2000 N.Piskunov wrote. 09-Nov-2011 NP added weights and 2nd derivative constraint as LAM2

Parameters:

f (Callable) – Function to fit
filter (int) – Smoothing parameter of the optimal filter (or polynomial degree of poly is True)
iter (int) – maximum number of iterations [def: 40]
eps (float) – convergence level [def: 0.001]
mn (float) – minimum function values to be considered [def: min(f)]
mx (float) – maximum function values to be considered [def: max(f)]
lam2 (float) – constraint on 2nd derivative
weight (array(float)) – vector of weights.

pyreduce.util.vac2air(wl_vac)[source]¶: Convert vacuum wavelengths to wavelengths in air Author: Nikolai Piskunov

pyreduce.cli module¶

Click-based CLI for PyReduce.

Usage:: uv run reduce –help uv run reduce bias INSTRUMENT –files bias/*.fits –output output/ uv run reduce run reduction.yaml

pyreduce.cli.main()[source]¶: Entry point for the CLI.

Module contents¶

class pyreduce.TqdmLoggingHandler(level=0)[source]¶

Bases: Handler

emit(record)[source]¶

Do whatever it takes to actually log the specified logging record.

This version is intended to be implemented by subclasses and so raises a NotImplementedError.

pyreduce package¶

Subpackages¶

Submodules¶

pyreduce.pipeline module¶

pyreduce.reduce module¶

Authors¶

Version¶

License¶

pyreduce.configuration module¶

pyreduce.extract module¶

Authors¶

Version¶

License¶

pyreduce.trace module¶

pyreduce.trace_model module¶

pyreduce.spectra module¶

pyreduce.wavelength_calibration module¶

pyreduce.slit_curve module¶

Authors¶

Version¶

License¶

pyreduce.combine_frames module¶

pyreduce.continuum_normalization module¶

pyreduce.estimate_background_scatter module¶

pyreduce.rectify module¶

pyreduce.clipnflip module¶

pyreduce.cwrappers module¶

pyreduce.datasets module¶

pyreduce.echelle module¶

pyreduce.util module¶

pyreduce.cli module¶

Module contents¶

PyReduce

Navigation

Related Topics