The cascade.spectral_extraction module

Module defining the spectral extraction functionality used in cascade.

directional_filters(return_valid_angle_range=False)[source]

Directional filters for smoothing and filtering.

These filters can be used in a Nagao&Matsuyama like edge preserving smoothing approach and are apropriate for dispersed spectra with a vertical dispersion direction. If the angle from vertical of the spectral trace of the dispersed light exceeds +- max(angle) radians, additional larger values need to be added to the angles list.

Parameters:

return_valid_angle_range ('bool') – optional, if True it returns the maximum angle range of the directional filters

Returns:

  • nm_mask (numpy.ndarray of ‘bool’ type) – Array containing all oriented masks used for edge preserving smooting.

  • maximum_angle_range (‘tuple’ of ‘float’) – If return_valid_angle_range is True, the maximum range of angles from vertical is returned

Notes

When adding kernels, make sure the maximum is in the central pixel

create_edge_mask(kernel, roi_mask)[source]

Create an edge mask.

Helper function for the optimal extraction task. This function creates an edge mask to mask all pixels for which the convolution kernel extends beyond the region of interest.

Parameters:
  • kernel (array_like) – Convolution kernel specific for a given instrument and observing mode, used in tasks such as replacing bad pixels and spectral extraction.

  • roi_mask (ndarray) – Mask defining the region of interest from which the speectra will be extracted.

Returns:

edge_mask (‘array_like’) – The edge mask based on the input kernel and roi_mask

determine_optimal_filter(ImageCube, Filters, ROIcube, selector)[source]

Determine optimal fileter.

Determine the optimal Filter for the image cube using a procedure similar to Nagao & Matsuyama edge preserving filtering.

Parameters:
  • ImageCube ('ndarray') – Cube of Spectral images.

  • Filters ('ndarray') – Cube of directional filters

  • ROIcube ('ndarray' of 'bool') – Region of Interests for the input ImageCube

  • selector ('list') – list containing all relevant information for each pixel within the ROI to select the sub region in the image cube and filter cube on which the filtering will be applied.

Returns:

  • selectorNumber (‘int’) – id number of the selector (pixel)

  • optimum_filter_index (‘ndarray’ of ‘int’) – array containing the optimal filter number for each pixel within the ROI in the image cube.

  • SubImageOptimalMean (‘ndarray’) – Optimal mean for each sub image in the sub image cube defined by the selector.

  • SubImageOptimaVariance (‘ndarray’) – Variance for each sub image in the sub image cube defined by the selector using the optimal Filter.

define_image_regions_to_be_filtered(ROI, filterShape)[source]

Define image regions to be filtered.

This function defines all pixels and corresponding sub-regions in the data cube to be be filtered.

Parameters:
  • ROI ('ndarray' of 'float') – Region of interest on the detector images

  • fiterShape ('tuple' of 'int') – y (dispersion direction) and x (spatial direstion) size of the directional filters.

Returns:

enumerated_sub_regions (‘list’) – enumerated definition of all regions in the spectral image cube and corresponding regions of the direction filter.

filter_image_cube(data_in, Filters, ROIcube, enumeratedSubRegions, useMultiProcesses=True, ReturnCleanedData=True)[source]

Filter image cube.

This routine filters in input data clube using an optimal choise of a directional filter and returns the filtered (smoothed) data. Optionally it also returns a cleaned dataset, where the masked data is replaced by the fitered data.

Parameters:
  • data_in ('np.ndarray' of 'float') – Input spectral data cube.

  • Filters ('ndarray' of 'float') – Directional filter kernels.

  • ROIcube ('ndarray' of 'bool') – Region of interest for each spectral image.

  • enumeratedSubRegions ('list') – List containing all information to define the regions around the pixels within the ROI used in the filtering. This list is created with the define_image_regions_to_be_filtered function.

  • useMultiProcesses ('bool' (optional|True)) – Flag to choose to use the multiprocessing module or not.

  • ReturnCleanedData ('bool' (optional|True)) – Flag to choose to return a cleaned data set or not.

Returns:

  • optimalFilterIndex (‘MaskedArray’ of ‘int’) – Index number of the optimal filter for each pixel in the input data.

  • filteredImage (‘MaskedArray’ of ‘float’) – The filtered data set.

  • filteredImageVariance (‘MaskedArray’ of ‘float’) – The variance of the filtered data set.

  • cleanedData (‘MaskedArray’ of ‘float’) – The cleaned data set.

iterative_bad_pixel_flagging(dataset, ROIcube, Filters, enumeratedSubRegions, sigmaLimit=4.0, maxNumberOfIterations=12, fractionalAcceptanceLimit=0.005, useMultiProcesses=True, maxNumberOfCPUs=2)[source]

Flag bad pixels.

This routine flags the bad pixels found in the input dataset and creates a cleaned dataset.

Parameters:
  • dataset ('SpectralDataTimeSeries')

  • ROIcube ('ndarray' of 'bool') – Region of interest for each spectral image.

  • Filters ('ndarray' of 'float') – Directional filter kernels.

  • enumeratedSubRegions ('list') – List containing all information to define the regions around the pixels within the ROI used in the filtering. This list is created with the define_image_regions_to_be_filtered function.

  • sigmaLimit ('float' (optional|4.0)) – Standard diviation limit used to identify bad pixels

  • maxNumberOfIterations ('int' (optional|12)) – The maximum number of iterations used in bad pixel masking

  • fractionalAcceptanceLimit ('float' (optional|0.005)) – Fractional number of bad pixels still found in the dataset below which the iteration can be terminated.

  • useMultiProcesses ('bool' (optional|True)) – Use multiprocessing or not.

  • max_number_of_cpus ('int' (optional|12)) – Maximum number of CPU’s which can be used.

Returns:

  • dataset (‘SpectralDataTimeSeries’) – Input data set containing the observed spectral time series. After completinion of this function, the dataset is returned with an updated mask with all diviating pixels flaged. Tho indicate this the flag isSigmaCliped=True is set.

  • filteredDataset (‘SpectralDataTimeSeries’) – The optimal filterd dataset. Includied in this dataset is the optimal Filter Index, i.e. the index indicating the used optimal filter for each pixel. To indicate that this is a filtered dataset the flag isFilteredData=True is set.

  • cleanedDataset – The cleaned dataset. To indicte that this is a cleaned dataset, the isCleanedData=True is set.

Notes

In the current version no double progress bar is used as in some IDE’s double progress bars do not work properly.

extract_spectrum(dataset, ROICube, extractionProfile=None, optimal=False, verbose=False, verboseSaveFile=None)[source]

Extract 1D spectra.

This routine extracts the spectrum both optimally as well as using an aperture.

Parameters:
  • dataset ('SpectralDataTimeSeries') – Input spectral dataset

  • ROIcube ('ndarray' of 'bool') – Region of interest for each spectral image.

  • extractionProfile ('MaskedArray', optinonal) – Normalized extraction profile. Has to be set if optimal=True

  • verbose ('bool', optional) – If true diagnostic plots will be generated.

  • verboseSaveFile ('str', optional) – If not None, verbose output will be saved to the specified file.

Returns:

extracted1dDataset (SpectralDataTimeSeries’) – Extracted 1D spectral timeseries dataset

Raises:

ValueError

create_extraction_profile(fiteredSpectralDataset, ROI=None)[source]

Create an extraction profile.

Parameters:
  • fiteredSpectralDataset ('SpectralDataTimeSeries') – Input filtered dataset on which the extraction profile is based

  • ROI ('ndarray' of 'bool', optional)

Returns:

spectralExtractionProfile (‘MaskedArray’)

determine_relative_source_position(spectralImageCube, ROICube, refIntegration, upsampleFactor=111, AngleOversampling=2)[source]

Determine the shift of the spectra relative to the first integration.

Parameters:
  • spectralImageCube ('ndarray') – Input spectral image data cube. Fist dimintion is dispersion direction, second dimintion is cross dispersion direction and the last dimension is time.

  • ROICube ('ndarray' of 'bool') – Cube containing the Region of interest for each integration. If not given, it is assumed that the mask of the spectralImageCube contains the region of interest.

  • refIntegration ('int') – Index number of of integration to be taken as reference.

  • upsampleFactor ('int', optional) – integer factor by which to upsample image for FFT analysis to get sub-pixel accuracy. Default value is 111

  • AngleOversampling ('int', optional) – Oversampling factor for angle determination, Default value 2

Returns:

relativeSourcePosition (‘collections.OrderedDict’) – relative rotation angle, scaling and x and y position as a function of time.

Raises:

ValueError – In case refIntegration exceeds number of integrations

Notes

Note that for sub-pixel registration to work correctly it should be performed on cleaned data i.e. bad pixels have been identified and corrected using the tools in this module.

warp_polar(image, center=None, *, radius=None, AngleOversampling=1, output_shape=None, scaling='linear', multichannel=False, **kwargs)[source]

Remap image to polor or log-polar coordinates space.

Parameters:
  • image (ndarray) – Input image. Only 2-D arrays are accepted by default. If multichannel=True, 3-D arrays are accepted and the last axis is interpreted as multiple channels.

  • center (tuple (row, col), optional) – Point in image that represents the center of the transformation (i.e., the origin in cartesian space). Values can be of type float. If no value is given, the center is assumed to be the center point of the image.

  • radius (float, optional) – Radius of the circle that bounds the area to be transformed.

  • AngleOversampling (int) – Oversample factor for number of angles

  • output_shape (tuple (row, col), optional)

  • scaling ({'linear', 'log'}, optional) – Specify whether the image warp is polar or log-polar. Defaults to ‘linear’.

  • multichannel (bool, optional) – Whether the image is a 3-D array in which the third axis is to be interpreted as multiple channels. If set to False (default), only 2-D arrays are accepted.

  • **kwargs (keyword arguments) – Passed to transform.warp.

Returns:

warped (ndarray) – The polar or log-polar warped image.

Examples

Perform a basic polar warp on a grayscale image: >>> from skimage import data >>> from skimage.transform import warp_polar >>> image = data.checkerboard() >>> warped = warp_polar(image) Perform a log-polar warp on a grayscale image: >>> warped = warp_polar(image, scaling=’log’) Perform a log-polar warp on a grayscale image while specifying center, radius, and output shape: >>> warped = warp_polar(image, (100,100), radius=100, … output_shape=image.shape, scaling=’log’) Perform a log-polar warp on a color image: >>> image = data.astronaut() >>> warped = warp_polar(image, scaling=’log’, multichannel=True)

highpass(shape)[source]

Return highpass filter to be multiplied with fourier transform.

Parameters:

shape ('ndarray' of 'int') – Input shape of 2d filter

Returns:

filter – high pass filter

_log_polar_mapping(output_coords, k_angle, k_radius, center)[source]

Inverse mapping function to convert from cartesion to polar coordinates.

Parameters:
  • output_coords (ndarray) – (M, 2) array of (col, row) coordinates in the output image

  • k_angle (float) – Scaling factor that relates the intended number of rows in the output image to angle: k_angle = nrows / (2 * np.pi)

  • k_radius (float) – Scaling factor that relates the radius of the circle bounding the area to be transformed to the intended number of columns in the output image: k_radius = width / np.log(radius)

  • center (tuple (row, col)) – Coordinates that represent the center of the circle that bounds the area to be transformed in an input image.

Returns:

coords (ndarray) – (M, 2) array of (col, row) coordinates in the input image that correspond to the output_coords given as input.

_determine_relative_source_shift(reference_image, image, referenceROI=None, ROI=None, upsampleFactor=111, space='real')[source]

Determine the relative shift of the spectral images.

This routine determine the relative shift between a reference spectral image and another spectral image.

Parameters:
  • reference_image ('ndarray or np.ma.MaskedArray' of 'float') – Reference spectral image

  • image ('ndarray or np.ma.MaskedArray' of 'float') – Spectral image

  • referenceROI (ndarray' of 'bool', optional)

  • ROI ('ndarray' of 'bool', optional)

  • upsampleFactor ('int', optional) – Default value is 111

  • space ('str', optional) – Default value is ‘real’

Returns:

  • relativeImageShiftY – relative shift compared to the reference image in the dispersion direction of the light (from top to bottom, shortest wavelength should be at row 0. Note that this shift is defined such that shifting a spectral image by this amound will place the trace at the exact same position as that of the reference image

  • relativeImageShiftX – relative shift compared to the reference image in the cross-dispersion direction of the light (from top to bottom, shortest wavelength should be at row 0. Note that this shift is defined such that shifting a spectral image by this amound will place the trace at the exact same position as that of the reference image.

register_telescope_movement(cleanedDataset, ROICube=None, nreferences=6, mainReference=4, upsampleFactor=111, AngleOversampling=2, verbose=False, verboseSaveFile=None, maxNumberOfCPUs=2, useMultiProcesses=True)[source]

Register the telescope movement.

Parameters:
  • cleanedDataset ('SpectralDataTimeSeries') – Input dataset. Note that for image registration to work properly, bad pixels need ti be removed (cleaned) first. This routine checks if a cleaned dataset is used by checking for the isCleanedData flag.

  • ROICube ('ndarray' of 'bool', optional) – Cube containing the Region of interest for each integration. If not given, it is assumed that the mask of the cleanedDataset contains the region of interest.

  • nreferences ('int', optional) – Default is 6.

  • mainReference ('int', optional) – Default is 4.

  • upsampleFactor ('int, optional) – Upsample factor of FFT images to determine relative movement at sub-pixel level. Default is 111

  • AngleOversampling ('int, optional) – Upsampling factor of the angle in the to polar coordinates transformed FFT images to determine the relative rotation adn scale change. Default is 2.

  • verbose ('bool', optional) – If true diagnostic plots will be generated. Default is False

  • verboseSaveFile ('str', optional) – If not None, verbose output will be saved to the specified file.

  • max_number_of_cpus ('int', optional) – Maxiumum bumber of cpu’s to be used.

Returns:

spectralMovement (‘OrderedDict’) – Ordered dict containing the relative rotation, scaling, and movement in the dispersion and cross dispersion direction.

Raises:

ValueError, TypeError – Errors are raised if certain data is not present of from the wrong type.

_determine_relative_rotation_and_scale(reference_image, referenceROI, image, ROI, upsampleFactor=111, AngleOversampling=2)[source]

Determine rotation and scalng changes.

This routine determines the relative rotation and scale change between an reference spectral image and another spectral image.

Parameters:
  • reference_image ('ndarray or np.ma.MaskedArray' of 'float') – Reference image

  • referenceROI ('ndarray' of 'float')

  • image ('ndarray or np.ma.MaskedArray' of 'float') – Image for which the rotation and translation relative to the reference image will be determined

  • ROI ('ndarray' of 'float') – Region of interest.

  • upsampleFactor ('int', optional) – Upsampling factor of FFT image used to determine sub-pixel shift. By default set to 111.

  • AngleOversampling ('int', optional) – Upsampling factor of the FFT image in polar coordinates for the determination of sub-degree rotation. Set by default to 2.

Returns:

  • relative_rotation – Relative rotation angle in degrees. The angle is defined such that the image needs to be rotated by this angle to have the same orientation as the reference spectral image

  • relative_scaling – Relative image scaling

_derotate_image(image, angle, ROI=None, order=3)[source]

Derotate image.

Parameters:
  • image ('2-D ndarray' of 'float') – Input image to be de-rotated by ‘angle’ degrees.

  • ROI ('2-D ndarray' of 'bool') – Region of interest (default None)

  • angle ('float') – Rotaton angle in degrees

  • order ('int') – Order of the used interpolation in the rotation function of the skimage package.

Returns:

derotatedImage (‘2-D ndarray’ of ‘float’) – The zero padded and derotated image.

_pad_to_size(image, h, w)[source]

Zero pad the input image to an image of hight h and width w.

Parameters:
  • image ('2-D ndarray' of 'float') – Input image to be zero-padded to size (h, w).

  • h ('int') – Hight (number of rows) of output image.

  • w ('int') – Width (number of columns) of output image.

Returns:

padded_image (‘2-D ndarray’ of ‘float’)

_pad_region_of_interest_to_square(image, ROI=None)[source]

Pad ROI to square.

Zero pad the extracted Region Of Interest of a larger image such that the resulting image is squared.

Parameters:
  • image ('2-D ndarray' of 'float') – Input image to be de-rotated by ‘angle’ degrees.

  • ROI ('2-D ndarray' of 'bool') – Region of interest (default None)

Returns:

padded_image (‘2-D ndarray’ of ‘float’)

correct_wavelength_for_source_movent(datasetIn, spectral_movement, useScale=False, useCrossDispersion=False, verbose=False, verboseSaveFile=None)[source]

Correct wavelengths for source movement.

This routine corrects the wavelength cube attached to the spectral image data cube for source (telescope) movements

Parameters:
  • datasetIn ('SpectralDataTimeSeries') – Input dataset for which the waveength will be corrected for telescope movement

  • spectral_movement ('OrderedDict') –

    Ordered dict containing the relative rotation, scaling,

    and movement in the dispersion and cross dispersion direction.

    useScale’bool’, optional

    If set the scale parameter is used to correct the wavelength. Default is False.

  • useCrossDispersion ('bool', optional) – If set the coress dispersion movement is used to correct the wavelength. Default is False.

  • verbose ('bool', optional) – If true diagnostic plots will be generated. Default is False

  • verboseSaveFile ('str', optional) – If not None, verbose output will be saved to the specified file.

Returns:

dataset_out (‘SpectralDataTimeSeries’) – The flag isMovementCorrected=True is set to indicate that this dataset is corrected

Notes

Scaling changes are not corrected at the moment. Note that the used rotation and translation to correct the wavelengths is the relative source movement defined such that shifting the observed spectral image by these angles and shifts the position would be identical to the reference image. The correction of the wavelength using the reference spectral image is hence in the oposite direction.

rebin_to_common_wavelength_grid(dataset, referenceIndex, nrebin=None, verbose=False, verboseSaveFile=None, return_weights=False)[source]

Rebin the spectra to single wavelength per row.

Parameters:
  • dataset ('SpectralDataTimeSeries') – Input dataset

  • referenceIndex ('int') – Exposure index number which will be used as reference defining the uniform wavelength grid.

  • nrebin ('float', optional) – rebinning factor for the new wavelength grid compare to the old.

  • verbose ('bool', optional) – If True, diagnostic plots will be created

  • verboseSaveFile ('str', optional) – If not None, verbose output will be saved to the specified file.

  • return_weights ('bool', optional) – If set, returns weights used in rebinning.

Returns:

rebinnedDataset (‘SpectralDataTimeSeries’) – Output to common wavelength grid rebinned dataset

determine_center_of_light_posision(cleanData, ROI=None, verbose=False, quantileCut=0.5, orderTrace=2)[source]

Determine the center of light position.

This routine determines the center of light position (cross-dispersion) of the dispersed light. The center of light is defined in a similar way as the center of mass. This routine also fits a polynomial to the spectral trace.

Parameters:
  • cleanData ('maskedArray' or 'ndarray') – Input data

  • ROI ('ndarray' of 'bool', optional) – Region of interest

  • verbose ('bool') – Default is False

  • quantileCut ('float', optional) – Default is 0.5

  • orderTrace ('int') – Default is 2

Returns:

  • total_light (‘ndarray’) – Total summed signal on the detector as function of wavelength.

  • idx (‘int’)

  • COL_pos (‘ndarray’) – Center of light position of the dispersed light.

  • ytrace (‘ndarray’) – Spectral trace position in fraction of pixels in the dispersion direction

  • xtrace (‘ndarray’) – Spectral trace position in fraction of pixels in the cross dispersion direction

determine_absolute_cross_dispersion_position(cleanedDataset, initialTrace, ROI=None, verbose=False, verboseSaveFile=None, quantileCut=0.5, orderTrace=2)[source]

Determine the initial cross dispersion position.

This routine updates the initial spectral trace for positional shifts in the cross dispersion direction for the first exposure of the the time series observation.

Parameters:
  • cleanedDataset ('SpectralDataTimeSeries')

  • initialTrace ('OrderedDict') – input spectral trace.

  • ROI ('ndarray' of 'bool') – Region of interest.

  • verbose ('bool', optional) – If true diagnostic plots will be generated. Default is False

  • verboseSaveFile ('str', optional) – If not None, verbose output will be saved to the specified file.

  • quantileCut ('float', optional) – Default is 0.5

  • orderTrace ('int', optional) – Default is 2

Returns:

  • newShiftedTrace (‘OrderedDict’) – To the observed source poisiton shifted spectral trace

  • newFittedTrace (‘OrderedDict’) – Trace determined by fit to the center of light position.

  • initialCrossDispersionShift (‘float’) – Shift between initial guess for spectral trace position and fitted trace position of the first spectral image.

combine_scan_samples(datasetIn, scanDictionary, verbose=False)[source]

Combine all (scan) samples.

This routine creates a new SpectralDataTimeSeries of integration averaged spectra from time series data per sample-up-the-ramp.

Parameters:
  • datasetIn ('SpectralDataTimeSeries') – Input dataset

  • scanDictionary ('dict') – Dictionary containg relevant information about the scans

  • verbose ('bool', optional) – If True, diagnostic plots will be created (default False).

Returns:

combinedDataset (‘SpectralDataTimeSeries’) – Output dataset with average signals per integration

Raises:

AttributeError

sigma_clip_data(datasetIn, sigma, nfilter)[source]

Perform sigma clip on science data to flag bad data.

Parameters:
  • datasetIn ('SpectralDataTimeSeries') – Input dataset

  • sigma (float) – Sigma value of sigmaclip.

  • nfilter ('int') – Filter length for sigma clip.

Returns:

datasetOut (‘SpectralDataTimeSeries’) – Input data set containing the observed spectral time series. After completinion of this function, the dataset is returned with an updated mask with all diviating pixels flaged. Tho indicate this the flag isSigmaCliped=True is set.

sigma_clip_data_cosmic(data, sigma)[source]

Sigma clip of time series data cube allong the time axis.

Parameters:
  • data (ndarray) – Input data to be cliped, last axis of data to be assumed the time

  • sigma (float) – Sigma value of sigmaclip

Returns:

sigmaClipMask (ndarray of ‘bool’) – Updated mask for input data with bad data points flagged (=1)

create_cleaned_dataset(datasetIn, ROIcube, kernel, stdvKernelTime)[source]

Create a cleaned dataset to be used in regresion analysis.

Parameters:
  • datasetIn ('SpectralDataTimeSeries') – Input dataset

  • ROIcube ('ndarray' of 'bool') – Region of interest.

  • kernel ('ndarray') – Instrument convolution kernel

  • stdvKernelTime ('float') – Standeard devistion in time direction used in convolution.

Returns:

cleanedDataset (SpectralDataTimeSeries) – A cleaned version of the spectral timeseries data of the transiting exoplanet system

compressROI(ROI, compressMask)[source]

Remove masked wavelengths from ROI.

Parameters:
  • ROI ('ndarray') – Region of interest on detector.

  • compressMask ('ndarray') – Compression mask indicating all valid data.

Returns:

compressedROI (‘ndarray’) – Row (wavelength) compressed region of interest.

compressSpectralTrace(spectralTrace, compressMask)[source]

Remove masked wavelengths from spectral trace.

Parameters:
  • spectralTrace ('dict') – Spectral trace of the dispersed light on the detector.

  • compressMask ('ndarray') – Compression mask indicating all valid data.

Returns:

compressedsSpectralTrace (‘dict’) – Row (wavelength) compressed spectral trace.

compressDataset(datasetIn, ROI)[source]

Remove all flaged wavelengths from data set.

Parameters:
  • datasetIn ('SpectralDataset') – Spectral dataset.

  • ROI ('ndarray') – Region of interest.

Returns:

compressedDataset (SpectralDataset’) – Row (wavelength) compressed dataset.

correct_initial_wavelength_shift(referenceDataset, cascade_configuration, *otherDatasets)[source]

Determine if there is an initial wavelength shift and correct.

Parameters:
  • referenceDataset ('cascade.data_model.SpectralDataTimeSeries') – Dataset who’s wavelength is used as refernce of the wavelength correction.

  • cascade_configuration ('cascade.initialize.initialize.configurator') – Singleton containing the confifuration parameters of cascade.

  • **otherDatasets ('cascade.data_model.SpectralDataTimeSeries') – Optional. Other datasets assumed to have the same walengths as the reference dataset and which will be corrected simultaneously with the reference.

Returns:

  • referenceDataset (‘list’ of ‘cascade.data_model.SpectralDataTimeSeries’) – Dataset with corrected wavelengths.

  • otherDatasets_list (‘list’ of ‘cascade.data_model.SpectralDataTimeSeries’) – Optinal output.

  • modeled_observations (‘list’ of ‘ndarray’)

  • stellar_model (‘list’ of ‘ndarray’)

  • corrected_observations (‘list’ of ‘ndarray’)

renormalize_spatial_scans(referenceDataset, *otherDatasets)[source]

bla.

Parameters:
  • referenceDataset (TYPE) – DESCRIPTION.

  • *otherDatasets (TYPE) – DESCRIPTION.

Returns:

TYPE – DESCRIPTION.

_define_band_limits(wave)[source]

Define the limits of the rebin intervals.

Parameters:

wave ('ndarray' of 'float') – Wavelengths (1D or 2D array)

Returns:

limits (‘tuple’) – lower and upper band limits

Raises:
  • TypeError – When the input wavelength array has more than 2 dimensions

  • ValueError – When the wavelength is not defined for each data point

_define_covariance_matrix(x_stddev=None, y_stddev=None, theta=None)[source]

Define covariance matrix.

Define 2D covariance matrix based on standard deviation in x and y and rotation angle

_define_rebin_weights(lr0, ur0, lr, ur)[source]

Define the weights used in spectral rebin.

Define the summation weights (i.e. the fractions of the original intervals used in the new interval after rebinning)

Parameters:
  • lr0 ('ndarray') – lower limits wavelength band of new wavelength grid

  • ur0 ('ndarray') – upper limits wavelength band of new wavelength grid

  • lr ('ndarray') – lower limits

  • ur ('ndarray') – upper limits

Returns:

weights (‘ndarray’ of ‘float’) – Summation weights used to rebin to new wavelength grid

Raises:

TypeError – When the input wavelength array has more than 2 dimensions.

_linear_polar_mapping(output_coords, k_angle, k_radius, center)[source]

Inverse mapping function to convert from cartesion to polar coordinates.

Parameters:
  • output_coords (ndarray) – (M, 2) array of (col, row) coordinates in the output image

  • k_angle (float) – Scaling factor that relates the intended number of rows in the output image to angle: k_angle = nrows / (2 * np.pi)

  • k_radius (float) – Scaling factor that relates the radius of the circle bounding the area to be transformed to the intended number of columns in the output image: k_radius = ncols / radius

  • center (tuple (row, col)) – Coordinates that represent the center of the circle that bounds the area to be transformed in an input image.

Returns:

coords (ndarray) – (M, 2) array of (col, row) coordinates in the input image that correspond to the output_coords given as input.

_rebin_spectra(spectra, errors, weights)[source]

Rebin spectra.

Parameters:
  • spectra ('ndarray') – Input spectral data values

  • errors ('ndarray') – Input error on spectral data values

  • weights ('ndarray') – rebin weights

Returns:

  • newSpectra (‘ndarray’) – Output spectra

  • newErrors (‘ndarray’) – Output error