pyklip.instruments.utils package
Submodules
pyklip.instruments.utils.global_centroid module
- class pyklip.instruments.utils.global_centroid.Consumer(task_queue, result_queue)[source]
Bases:
Process
- pyklip.instruments.utils.global_centroid.fitallrelcen(cubes, ivars, r1=15, r2=50, smooth=True, maxcpus=1)[source]
Function fitallrelcen. Fit for the relative centroids between all pairs of frames at the central wavelength using the PSF in an annulus around the center. Return the best-fit relative offets.
- Parameters:
cubes – all image cubes in the data set, shape (ncube, nwv, ny, nx)
ivars – inverse variance frames corresponding to cubes
r1 – int, minimum separation in lenslets from the image center for annular reference region. Default 15.
r2 – int, maximum separation in lenslets from the image center for annular reference region. Default 50.
smooth – whether to smooth the reference image before fitting
maxcpus – int, maximum number of CPUs to allocate for parallelization. Default 1/2 of the available CPUs.
- Returns:
1D ndarray of the relative centers in x ysol : 1D ndarray of the relative centers in y
- Return type:
xsol
- pyklip.instruments.utils.global_centroid.fitcen(cube, ivar, lam, spotsep=None, guess_center_loc=None, i1=1, i2=-1, r1=15, r2=35, spot_dx=4, astrogrid='XYdiag', smooth=True)[source]
Function fitcen. Fit for the center of a CHARIS data cube using the satellite spots by maximizing the agreement between scaled cube slices around the spot locations. If no spot locations are provided, use only the diffraction pattern itself in an annulus around the image center.
- Parameters:
cube – 3D ndarray, CHARIS data cube
ivar – 3D ndarray, inverse variance of the CHARIS data cube
lam – 1D ndarray, wavelengths in microns
spotsep – float or None. If float, separation of the satellite spots in units of lambda/D. If None, only use the diffraction pattern in an annulus between r1 and r2.
guess_center_loc – manually specify initial location of image center if necessary, in [x, y] format
i1 – int, first slice of the data cube to use. Default 1 (skip slice 0)
i2 – int, high limit of slices to use. Default -1 (skip last slice)
r1 – float, minimum separation from approximate center for the annulus of the diffraction pattern to use in centroiding. Default 15
r2 – float, maximum separation from approximate center for the annulus of the diffraction pattern to use in centroiding. Default 35
spot_dx – float, radius around spot location to cut out in order to match the spot location as a function of wavelength. Default 4
astrogrid – astrogrid status read from the header, determines the pattern of the diffraction spots
smooth – whether to smooth image before fitting
- Returns:
- list of floats
p[0] is the angle of spots in radians p[1] is the separation in lambda/D p[2] - p[-1] are the coefficients of the polynomial fit to the centroid
- Return type:
p
- pyklip.instruments.utils.global_centroid.fitcen_parallel(infiles, cubes, ivars, prihdrs, astrogrid_status=None, astrogrid_sep=None, smooth_coef=True, guess_center_loc=None, smooth_cubes=True, maxcpus=1)[source]
Function fitcen_parallel. Centroid a series of CHARIS data cubes in parallel using fitcen. By default, get the wavelengths and astrogrid parameters from the headers. This might fail on early CHARIS data before the standardization of headers.
- Parameters:
infiles – input files
cubes – all image cubes in the data set corresponding to filenames, shape (ncube, nwv, ny, nx)
ivars – inverse variance frames corresponding to cubes
prihdrs – primary headers for the cubes
astrogrid_status – None or list of astrogrid configurations for SCExAO. If None, try to read the astrogrid configuration from the header. If this fails, assume there is no astrogrid and centroid using the general diffraction pattern. Default None.
astrogrid_sep – None or list of astrogrid spot separations in units of lambda/D. If None, try to read from the header. If that fails, centroid using the general diffraction pattern. Default None.
smooth_coef – boolean. smooth the nonlinear coefficients of the centroid fit (the terms proportional to lambda and lambda^2) over the sequence of cubes? Default True.
guess_center_loc – manually specify initial location of image center if necessary, in [x, y] format
smooth_cubes – whether to smooth the data before fitting
maxcpus – int, maximum number of CPUs to use in parallelization
- Returns:
[centroid_params, x, y, mask]
- centroid_params2D array of centroid parameters, first dimension is the number of files.
Second dimension is the length of the wavelength-dependent model.
x : x-coordinates of the centroid at the middle wavelength y : y-coordinates of the centroid at the middle wavelength mask : 1D boolean array, True if astrogrid was on, or None if the astrogrid was never on.
- pyklip.instruments.utils.global_centroid.fitrelcen(image1, image2, x, y, method='Powell')[source]
Function fitrelcen fits for the offset between two images without any wavelength dependence by calling _resid_xy, minimizing the sum of the squared differences between the images (after optimizing over the wavelength-dependent relative normalization).
- Parameters:
image1 – 2D ndarray, first image
image2 – 2D ndarray, second image
x – ndarray, x coordinates of pixels/lenslets to use
y – ndarray, y coordinates of pixels/lenslets to use
method – method passed to scipy.optimize.minimize. Default ‘Powell’.
- Returns:
- two floating point numbers giving the best-fit offset between
image1 and image2
- Return type:
xc, yc
- pyklip.instruments.utils.global_centroid.get_sats_satf(p, cube, lam, astrogrid='XYdiag')[source]
retrieves the pixel locations of all four satellite spots at each wavelength, and the negative sum of the four spot intensities at each wavelength
- Parameters:
p – list of floats, the coefficients of the polynomial fit to the centroid
cube – ndarray, data cube for which the centroid is fitted
lam – ndarray, wavelengths for the datacube
astrogrid – astrogrid status read from the header, determines the pattern of the diffraction spots
- Returns:
pixel locations (in [x,y] format) of all four satellite spots at each wavelength, shape (wvs, 4, 2) satf: float, peak fluxes (interpolated pixel value) at the fitted spot locations
- Return type:
sats
- pyklip.instruments.utils.global_centroid.polyfit(x, y, order=2, clip=2.5, niter=5, mask=None, return_y=True)[source]
Smooth a series of points with a polynomial, iteratively clipping outliers.
- Parameters:
x – 1D ndarray of x coordinates for the polynomial fit
y – 1D ndarray of y coordinates for the polynomial fit
order – int, order of the polynomial to fit. Default 2.
clip – float, number of sigma outliers to clip. Default 2.5.
niter – int, number of iterations of sigma clipping. Default 5.
mask – boolean ndarray or None: mask each y value? Default None
return_y – boolean, return smoothed y values (as opposed to coefficients of the polynomial fit)? Default True.
- Returns:
1D array, if return_y=True coef : array of the polynomial coefficients if return_y=False
- Return type:
y_smoothed
- pyklip.instruments.utils.global_centroid.recen(p, cube, lam, sats, satf, n=None, scale=False, head=None, outfile=None, mask=None, data_HDU1=True)[source]
Function recen recenters the input data according to the offset parameters given in the argument p. Optionally scale the data by wavelength to undo the scaling of diffraction. Return the recentered cube.
- Parameters:
p – list of floats, the coefficients of the polynomial fit to the centroid
cube – 3D ndarray, data cube
lam – 1D array of wavelengths
sats – fitted satellite spot indices before recentering
satf – fitted satellite spot fluxes
n – integer, spatial dimension of recentered cube. Default original cube size
scale – boolean, rescale by wavelength? Default False
head – fits header for output file. Default None
outfile – string or None, name of output file. Default None
mask – boolean lenslet mask
data_HDU1 – boolean, write data to HDU1 and leave HDU0 with no data? Default True.
- Returns:
3D ndarray, nlam x n x n, recentered (optionally scaled by wavelength) data cube.
- pyklip.instruments.utils.global_centroid.specphotcal(infiles, cubes, prihdrs, cencoef, aperture=1.0)[source]
Function specphotcal. Computes approximate photometry from the satellite spots (using aperture photometry) and scale each wavelength to this photometric value. This should crudely put the cubes in units of contrast, though it omits the scaling of satellite spot intensity with 1/lambda^2.
- Parameters:
infiles – list of file names with CHARIS data cubes
cubes – all image cubes in the data set corresponding to filenames, shape (ncube, nwv, ny, nx)
prihdrs – primary headers for the cubes
cencoef – 2D ndarray with coefficeitns
aperture – float, radius of aperture for photometry in units of lambda/D
- Returns:
photocalibration coefficients
- Return type:
all_phot
pyklip.instruments.utils.nair module
- pyklip.instruments.utils.nair.get_coeff_mathar(i, P, T, H, wvrange=1)[source]
Calculate the coefficients for the polynomial series expansion of index of refraction (Mathar (2008)) ***Only valid for between 1.3 and 2.5 microns! and 2.8 through 4.2 microns!
- Inputs:
i: degree of expansion in wavenumber P: Pressure in Pa T: Temperature in Kelvin H: relative humiditiy in % (i.e. between 0 and 100) wvrange (int): 1 = (1.3-2.5 um), 2 = (2.8 - 4.2 um)
- Returns:
Coefficient [cm^-i]
- Return type:
coeff
- pyklip.instruments.utils.nair.nMathar(wv, P, T, H=10)[source]
Calculate the index of refraction as given by Mathar (2008): http://arxiv.org/pdf/physics/0610256v2.pdf ***Only valid for between 1.3 and 2.5 microns!
- Inputs:
wv: wavelength in microns P: Pressure in Pa T: Temperature in Kelvin H: relative humiditiy in % (i.e. between 0 and 100)
- Returns:
index of refratoin
- Return type:
n
- pyklip.instruments.utils.nair.nRoe(wv, P, T, H=10)[source]
Compute n for air from the formula in Henry Roe’s PASP paper: http://arxiv.org/pdf/astro-ph/0201273v1.pdf which in turn is referenced from Allen’s Astrophysical Quantities.
- Inputs:
wv: wavelength in microns P: pressure in Pascal T: temperature in Kelvin H: relative humidity in % (0-100)
- Returns:
index of refraction of air
- Return type:
n
pyklip.instruments.utils.radonCenter module
- pyklip.instruments.utils.radonCenter.samplingRegion(size_window, theta=[45, 135], m=0.2, M=0.8, step=1, decimals=2, ray=False)[source]
This function returns all the coordinates of the sampling region, the center of the region is (0,0) When applying to matrices, don’t forget to SHIFT THE CENTER! Input:
size_window: the radius of the sampling region. The whole region should thus have a length of 2*size_window+1. theta: the angle range of the sampling region, default: [45, 135] for the anti-diagonal and diagonal directions. m: the minimum fraction of size_window, default: 0.2 (i.e., 20%). In this way, the saturated region can be excluded. M: the maximum fraction of size_window, default: 0.8 (i.e., 80%). Just in case if there’s some star along the diagonals. step: the seperation between sampling dots (units: pixel), default value is 1pix. decimals: the precisoin of the sampling dots (units: pixel), default value is 0.01pix. ray: only half of the line?
- Output: (xs, ys)
xs: x indecies, flattend. ys: y indecies, flattend.
Example
1. If you call “xs, ys = samplingRegion(5)”, you will get: xs: array([-2.83, -2.12, -1.41, -0.71, 0.71, 1.41, 2.12, 2.83, 2.83, 2.12, 1.41, 0.71, -0.71, -1.41, -2.12, -2.83] ys: array([-2.83, -2.12, -1.41, -0.71, 0.71, 1.41, 2.12, 2.83, -2.83, -2.12, -1.41, -0.71, 0.71, 1.41, 2.12, 2.83])) 2. For “radonCenter.samplingRegion(5, ray=True)”, you will get: xs: array([ 0.71, 1.41, 2.12, 2.83, -0.71, -1.41, -2.12, -2.83]) ys: array([ 0.71, 1.41, 2.12, 2.83, 0.71, 1.41, 2.12, 2.83])
- pyklip.instruments.utils.radonCenter.searchCenter(image, x_ctr_assign, y_ctr_assign, size_window, m=0.2, M=0.8, size_cost=5, theta=[45, 135], ray=False, smooth=2, decimals=2, output_cost=False)[source]
This function searches the center in a grid, calculate the cost function of Radon Transform (Pueyo et al., 2015), then interpolate the cost function, get the center which corresponds to the maximum value in the cost function.
- Input:
image: 2d array. x_ctr_assign: the assigned x-center, or starting x-position; for STIS, the “CRPIX1” header is suggested. x_ctr_assign: the assigned y-center, or starting y-position; for STIS, the “CRPIX2” header is suggested. size_window: half width of the sampling region; size_window = image.shape[0]/2 is suggested.
m & M: The sampling region will be (-M*size_window, -m*size_window)U(m*size_window, M*size_window).
size_cost: search the center within +/- size_cost pixels, i.e., a square region. theta: the angle range of the sampling region; default: [45, 135] for the anti-diagonal and diagonal directions. ray: is the theta a line or a ray? Default: line. smooth: smooth the cost function, for one pixel, replace it by the average of its +/- smooth neighbours; defualt = 2. decimals: the precision of the centers; default = 2 for a precision of 0.01.
- Output:
x_cen, y_cen
- pyklip.instruments.utils.radonCenter.smoothCostFunction(costFunction, halfWidth=0)[source]
smoothCostFunction will smooth the function within +/- halfWidth, i.e., to replace the value with the average within +/- halfWidth pixel. This function can be genrally used to smooth any 2D matrix. Input:
costFunction: original cost function, a matrix. halfWdith: the half width of the smoothing region, default = 0 pix.
- Output:
newFunction: smoothed cost function.
pyklip.instruments.utils.xcorr module
- pyklip.instruments.utils.xcorr.find_best_shift(frame, ref_frame, guess_center=None, inner_mask=5, outer_mask=None)[source]
Finds the best shift based on xcorr
- Parameters:
frame – frame to find offset of. Shape of (Y, X)
ref_frame – reference frame to align frame to. Shape of (Y, X)
- Returns:
best shift to shift frame to match reference frame
- Return type:
(dx, dy)