sourcefinder.measure#
Source measuring routines.
Attributes#
Classes#
str(object='') -> str |
Functions#
|
This function represents equation 2.67 from Spreeuw's thesis. |
|
Calculate source positional values by fitting a 2D Gaussian. |
|
Calculate the goodness-of-fit values. |
|
Calculate source positional values using moments. |
|
Calculate source properties using moments. |
Module Contents#
- class sourcefinder.measure.MomentParam[source]#
Bases:
str,enum.Enumstr(object=’’) -> str str(bytes_or_buffer[, encoding[, errors]]) -> str
Create a new string object from the given object. If encoding or errors is specified, then the object must expose a data buffer that will be decoded using the given encoding and error handler. Otherwise, returns the result of object.__str__() (if defined) or repr(object). encoding defaults to ‘utf-8’. errors defaults to ‘strict’.
- sourcefinder.measure.find_true_peak(peak, T, epsilon, msq, maxpix)[source]#
This function represents equation 2.67 from Spreeuw’s thesis.
- sourcefinder.measure.fitgaussian(pixels, params, fixed=None, max_nfev=None, bounds={})[source]#
Calculate source positional values by fitting a 2D Gaussian.
- Parameters:
- pixelsnp.ma.MaskedArray or np.ndarray
Pixel values (with bad pixels masked)
- paramsdict
Initial fit parameters (possibly estimated using the moments() function)
- fixeddict, default: None
Parameters & their values to be kept frozen (ie, not fitted)
- max_nfevint, default: None
Maximum number of calls to the error function
- boundsdict, default: {}
Can be a dict such as the extract.ParamSet().bounds attribute generated by the extract.ParamSet().compute_bounds method, but any dict with keys from FIT_PARAMS and (lower_bound, upper_bound, bool) tuples as values will do. The boolean argument accommodates for loosening a bound when a fit becomes unfeasible because of the bound, see the scipy.optimize.Bounds documentation and source code for background.
- Returns:
- dict
peak, total, x barycenter, y barycenter, semimajor, semiminor, theta (radians)
- Raises:
- exceptions.ValueError
In case of a bad fit.
Notes
Perform a least squares fit to an elliptical Gaussian.
If a dict called fixed is passed in, then parameters specified within the dict with the same keys as in FIT_PARAMS will be “locked” in the fitting process.
- sourcefinder.measure.goodness_of_fit(masked_residuals, noise, correlation_lengths)[source]#
Calculate the goodness-of-fit values.
- Parameters:
- masked_residualsnp.ma.MaskedArray
The pixel-residuals from the fit.
- noisefloat or np.ma.MaskedArray
An estimate of the noise level. Can also be a masked numpy array matching the data for per-pixel noise estimates.
- correlation_lengthstuple of two floats
Distance along the semimajor and semiminor axes of the clean beam beyond which noise is assumed uncorrelated. Some background: Aperture synthesis imaging yields noise that is partially correlated over the entire image. This has a considerable effect on error estimates. We approximate this by considering all noise within the correlation length completely correlated and beyond that completely uncorrelated.
- Returns:
- chisq: float
- reduced_chisqfloat
Notes
We do not use the standard chi-squared formula for calculating these goodness-of-fit values, see <https://en.wikipedia.org/wiki/Goodness_of_fit#Regression_analysis>. These values are related to, but not quite the same as reduced chi-squared. The reduced chi-squared is statistically invalid for a Gaussian model from the outset (see <http://arxiv.org/abs/1012.3754>).
We attempt to provide a resolution-independent estimate of goodness-of-fit (‘reduced chi-squared’) by estimating the number of independent pixels in the data, that we have used for fitting the Gaussian model, to normalize the chi-squared value.
However, this will sometimes imply that we are fitting a fractional number of datapoints less than 1! As a result, it doesn’t really make sense to try and apply the ‘degrees-of-freedom’ correction, as this would likely result in a negative
reduced_chisqvalue.
- sourcefinder.measure.moments(data, fudge_max_pix_factor, beam, beamsize, threshold=0)[source]#
Calculate source positional values using moments.
The first moment of the distribution is the barycenter of an ellipse. The second moments are used to estimate the rotation angle and the length of the axes.
- Parameters:
- datanp.ma.MaskedArray or np.ndarray
A rectangular region of 2D observational image data, typically with the mean background subtracted. Units: spectral brightness, usually Jy/beam.
- fudge_max_pix_factorfloat
Correct for the underestimation of the peak by taking the maximum pixel value. This factor is slightly larger than 1.
- beam3-tuple
tuple of 3 floats: (semimajor axis, semiminor axis, theta). The axes are in units of pixels and theta, the position angle of the major axis wrt the positive y-axis, is in radians.
- beamsizefloat
The FWHM size of the clean beam in square pixels.
- thresholdfloat, default: 0
Source parameters like the semimajor and semiminor axes derived from moments can be underestimated if one does not take account of the threshold that was used to segment the source islands. This threshold is typically the value of the analysisthresholdmap at the position of the island pixel with the highest spectral brightness.
- Returns:
- dict
Dictionary with seven items, i.e. peak spectral brightness, flux density, barycenter position in x and y, i.e. in pixel coordinates, semimajor and semiminor axes in pixel units and the position angle of the semi-major axis measured east from local north in radians.
- Raises:
- exceptions.ValueError
If input contains NaN values.
- sourcefinder.measure.moments_enhanced(source_island, noise_island, chunkpos, posx, posy, min_width, no_pixels, threshold, noise, maxpos, maxi, fudge_max_pix_factor, beam, beamsize, force_beam, correlation_lengths, clean_bias_error, frac_flux_cal_error, Gaussian_islands_map, Gaussian_residuals_map, dummy, computed_moments, significance, chisq, reduced_chisq)[source]#
Calculate source properties using moments.
Vectorized using the guvectorize decorator. Also calculates the signal-to-noise ratio of detections, chi-squared and reduced chi-squared statistics, and fills in maps with Gaussian islands and Gaussian residuals. Uses the first moments of the distribution to determine the barycenter of an ellipse, while the second moments estimate rotation angle and axis lengths.
- Parameters:
- source_islandnp.ndarray
Selected from the 2D image data by taking pixels above the analysis threshold, with its peak above the detection threshold. Flattened to a 1D ndarray. Units: spectral brightness, typically Jy/beam.
- noise_islandnp.ndarray
Pixel values selected from the 2D RMS noise map at the positions of the island. Flattened to a 1D ndarray. Units: spectral brightness, typically Jy/beam.
- chunkposnp.ndarray
Index array of length 2 denoting the position of the top-left corner of the rectangular slice encompassing the island, relative to the top-left corner of the image.
- posxnp.ndarray
Row indices of the pixels in source_island relative to the top-left corner of the rectangular slice encompassing the island. The top-left corner corresponds to posx = 0. Derived from the 2D image data.
- posynp.ndarray
Column indices of the pixels in source_island relative to the top-left corner of the rectangular slice encompassing the island. The top-left corner corresponds to posy = 0. Derived from the 2D image data.
- min_widthint
Minimum width (in pixels) of the island, derived as the lesser of its maximum width along the x and y axes.
- no_pixelsint
Number of pixels that constitute the island.
- thresholdfloat
Threshold used for segmenting source islands, which can affect parameters like semimajor and semiminor axes. A higher threshold may lead to a larger underestimate of the Gaussian axes. If the analysis threshold is known, this underestimate can be corrected. Units: spectral brightness, typically Jy/beam.
- noisefloat
Local noise, i.e., the standard deviation of the background pixel values at the position of the island’s peak pixel value. Units: spectral brightness, typically Jy/beam.
- maxposnp.ndarray with int32 as dtype and length 2
The position of the maximum pixel value within the island, relative to the top-left corner of the rectangular slice encompassing the island. Units: pixels.
- maxifloat
Peak pixel value within the island. Pixel values represent spectral brightness (typically in Jy/beam). For clarity: source_island[maxpos] == maxi.
- fudge_max_pix_factorfloat
Correction factor for underestimation of the peak by considering the maximum pixel value.
- beamnp.ndarray
Array of three floats: (semimajor axis, semiminor axis, theta). The axes are in units of pixels and theta, the position angle of the major axis wrt the positive y-axis, is in radians.
- beamsizefloat
FWHM size of the clean beam. Units: pixels.
- force_beamint
If 1, the restoring beam is used to derive the peak spectral brightnesses of the sources, see equation 2.66 of Spreeuw’s (2010) thesis. The sources are then assumed to be unresolved and their shapes are set equal to the restoring beam. If 0, the peak spectral brightnesses are derived using equation 2.67 of that thesis and enhanced moments are used to derive the elliptical Gaussian axes.
- correlation_lengthsnp.ndarray
Array of two floats describing distances along the semi-major and semi-minor axes of the clean beam beyond which noise is assumed uncorrelated. Units: pixels. Some background: Aperture synthesis imaging yields noise that is partially correlated over the entire image. This has a considerable effect on error estimates. All noise within the correlation length is approximated as completely correlated, while noise beyond is considered uncorrelated.
- clean_bias_errorfloat
Extra source of error based on Condon (PASP 109, 166, 1997) formulae.
- frac_flux_cal_errorfloat
Extra source of error based on Condon (PASP 109, 166, 1997) formulae.
- Gaussian_islands_mapnp.ndarray
Initially a 2D np.float32 array filled with zeros, same shape as the astronomical image being processed. Computed Gaussian islands are added to this array at pixel positions above the analysis threshold.
- Gaussian_residuals_mapnp.ndarray
Initially a 2D np.float32 array filled with zeros, same shape as Gaussian_islands_map. Residuals are computed by subtracting Gaussian_islands_map from the input image data.
- dummynp.ndarray
Empty array matching the shape of computed_moments, required due to limitations in guvectorize.
- computed_momentsnp.ndarray
Array of shape (10, 2) containing moments such as peak spectral brightness (Jy/beam), flux density (Jy), barycenter (pixels), semi-major and semi-minor axes (pixels), and position angle (radians), along with corresponding errors.
- significancefloat
The significance of a detection is defined as the maximum signal-to-noise ratio across the island. Often this will be the ratio of the maximum pixel value within the source island divided by the noise at that position. But for extended sources, the noise can perhaps decrease away from the position of the peak spectral brightness more steeply than the source spectral brightness, and the maximum signal-to-noise ratio can be found at a different position.
- chisqfloat
Chi-squared statistic indicating goodness-of-fit, derived in the same way as in the measuring.goodness_of_fit method.
- reduced_chisqfloat
Reduced chi-squared statistic indicating goodness-of-fit, derived in the same way as in the measuring.goodness_of_fit method.
- Returns:
- None
Outputs are written to Gaussian_islands_map, Gaussian_residuals_map, computed_moments, significance, chisq, and reduced_chisq.
- Raises:
- ValueError
If input contains NaN values.