Satellite Trails Detection

This module contains the tools needed for satellite detection within an ACS/WFC image as published in ACS ISR 2016-01.

Note

Only tested for ACS/WFC FLT and FLC images, but it should theoretically work for any instrument.

Requires skimage version 0.11.x or later.

skimage.transform.probabilistic_hough_line() gives slightly different results from run to run, but this should not matter since detsat() only provides crude approximation of the actual trail(s).

Performance is faster when plot=False, where applicable.

Currently not supported in TEAL and PyRAF.

Examples

>>> from acstools.satdet import detsat, make_mask, update_dq

Find trail segments for a single image and extension without multiprocessing, and display plots (not shown) and verbose information:

>>> results, errors = detsat(
...     'jc8m10syq_flc.fits', chips=[4], n_processes=1,
...     plot=True, verbose=True)
1 file(s) found...
Processing jc8m10syq_flc.fits[4]...
Rescale intensity percentiles: 110.161376953, 173.693756409
Length of PHT result: 42
min(x0)=   1, min(x1)= 269, min(y0)= 827, min(y1)= 780
max(x0)=3852, max(x1)=4094, max(y0)=1611, max(y1)=1545
buf=200
topx=3896, topy=1848
...
Trail Direction: Right to Left
42 trail segment(s) detected
...
End point list:
    1. (1256, 1345), (2770, 1037)
    2. (  11, 1598), ( 269, 1545)
...
Total run time: 22.4326269627 s
>>> results[('jc8m10syq_flc.fits', 4)]
array([[[1242, 1348],
        [2840, 1023]],
       [[1272, 1341],
        [2688, 1053]],
       ...
       [[2697, 1055],
        [2967, 1000]]])
>>> errors
{}

Find trail segments for multiple images and all ACS/WFC science extensions with multiprocessing:

>>> results, errors = detsat(
...     '*_flc.fits', chips=[1, 4], n_processes=12, verbose=True)
6 file(s) found...
Using 12 processes
Number of trail segment(s) found:
  abell2744-hffpar_acs-wfc_f814w_13495_11_01_jc8n11q9q_flc.fits[1]: 0
  abell2744-hffpar_acs-wfc_f814w_13495_11_01_jc8n11q9q_flc.fits[4]: 4
  abell2744_acs-wfc_f814w_13495_51_04_jc8n51hoq_flc.fits[1]: 2
  abell2744_acs-wfc_f814w_13495_51_04_jc8n51hoq_flc.fits[4]: 34
  abell2744_acs-wfc_f814w_13495_93_02_jc8n93a7q_flc.fits[1]: 20
  abell2744_acs-wfc_f814w_13495_93_02_jc8n93a7q_flc.fits[4]: 20
  j8oc01sxq_flc.fits[1]: 0
  j8oc01sxq_flc.fits[4]: 0
  jc8m10syq_flc.fits[1]: 0
  jc8m10syq_flc.fits[4]: 38
  jc8m32j5q_flc.fits[1]: 42
  jc8m32j5q_flc.fits[4]: 12
Total run time: 34.6021330357 s
>>> results[('jc8m10syq_flc.fits', 4)]
array([[[1242, 1348],
        [2840, 1023]],
       [[1272, 1341],
        [2688, 1053]],
       ...
       [[2697, 1055],
        [2967, 1000]]])
>>> errors
{}

For a given image and extension, create a DQ mask for a satellite trail using the first segment (other segments should give similar masks) based on the results from above (plots not shown):

>>> trail_coords = results[('jc8m10syq_flc.fits', 4)]
>>> trail_segment = trail_coords[0]
>>> trail_segment
array([[1199, 1357],
       [2841, 1023]])
>>> mask = make_mask('jc8m10syq_flc.fits', 4, trail_segment,
...                  plot=True, verbose=True)
Rotation: -11.4976988695
Hit image edge at counter=26
Hit rotate edge at counter=38
Run time: 19.476323843 s

Update the corresponding DQ array using the mask from above:

>>> update_dq('jc8m10syq_flc.fits', 6, mask, verbose=True)
DQ flag value is 16384
Input... flagged NPIX=156362
Existing flagged NPIX=0
Newly... flagged NPIX=156362
jc8m10syq_flc.fits[6] updated
acstools.satdet.detsat(searchpattern, chips=[1, 4], n_processes=4, sigma=2.0, low_thresh=0.1, h_thresh=0.5, small_edge=60, line_len=200, line_gap=75, percentile=(4.5, 93.0), buf=200, plot=False, verbose=True)

Find satellite trails in the given images and extensions. The trails are calculated using Probabilistic Hough Transform.

Note

The trail endpoints found here are crude approximations. Use make_mask() to create the actual DQ mask for the trail(s) of interest.

Parameters:
  • searchpattern (str) – Search pattern for input FITS images, as accepted by glob.glob().
  • chips (list) – List of extensions for science data, as accepted by astropy.io.fits. The default values of [1, 4] are tailored for ACS/WFC.
  • n_processes (int) – Number of processes for multiprocessing, which is only useful if you are processing a lot of images or extensions. If 1 is given, no multiprocessing is done.
  • sigma (float, optional) – The size of a Gaussian filter to use before edge detection. The default is 2, which is good for almost all images.
  • low_thresh (float, optional) – The lower threshold for hysteresis linking of edge pieces. This should be between 0 and 1, and less than h_thresh.
  • h_thresh (float, optional) – The upper threshold for hysteresis linking of edge pieces. This should be between 0 and 1, and greater than low_thresh.
  • small_edge (int, optional) – Size of perimeter of small objects to remove in edge image. This significantly reduces noise before doing Hough Transform. If it is set too high, you will remove the edge of the satellite you are trying to find.
  • line_len (int, optional) – Minimum line length for Probabilistic Hough Transform to fit.
  • line_gap (int, optional) – The largest gap in points allowed for the Probabilistic Hough Transform.
  • percentile (tuple of float, optional) – The percent boundaries to scale the image to before creating edge image.
  • buf (int, optional) – How close to the edge of the image the satellite trail has to be to be considered a trail.
  • plot (bool, optional) – Make plots of edge image, Hough space transformation, and rescaled image. This is only applicable if n_processes=1.
  • verbose (bool, optional) – Print extra information to the terminal, mostly for debugging. In multiprocessing mode, info from individual process is not printed.
Returns:

  • results (dict) – Dictionary mapping (filename, ext) to an array of endpoints of line segments in the format of [[x0, y0], [x1, y1]] (if found) or an empty array (if not). These are the segments that have been identified as making up part of a satellite trail.
  • errors (dict) – Dictionary mapping (filename, ext) to the error message explaining why processing failed.

Raises:

ImportError – Missing scipy or skimage>=0.11 packages.

acstools.satdet.make_mask(filename, ext, trail_coords, sublen=75, subwidth=200, order=3, sigma=4, pad=10, plot=False, verbose=False)

Create DQ mask for an image for a given satellite trail. This mask can be added to existing DQ data using update_dq().

Note

Unlike detsat(), multiprocessing is not available for this function.

Parameters:
  • filename (str) – FITS image filename.
  • ext (int, str, or tuple) – Extension for science data, as accepted by astropy.io.fits.
  • trail_coords (ndarray) – One of the trails returned by detsat(). This must be in the format of [[x0, y0], [x1, y1]].
  • sublen (int, optional) – Length of strip to use as the fitting window for the trail.
  • subwidth (int, optional) – Width of box to fit trail on.
  • order (int, optional) – The order of the spline interpolation for image rotation. See skimage.transform.rotate().
  • sigma (float, optional) – Sigma of the satellite trail for detection. If points are a given sigma above the background in the subregion then it is marked as a satellite. This may need to be lowered for resolved trails.
  • pad (int, optional) – Amount of extra padding in pixels to give the satellite mask.
  • plot (bool, optional) – Plot the result.
  • verbose (bool, optional) – Print extra information to the terminal, mostly for debugging.
Returns:

mask – Boolean array marking the satellite trail with True.

Return type:

ndarray

Raises:
  • ImportError – Missing scipy or skimage>=0.11 packages.
  • IndexError – Invalid subarray indices.
  • ValueError – Image has no positive values, trail subarray too small, or trail profile not found.
acstools.satdet.update_dq(filename, ext, mask, dqval=16384, verbose=True)

Update the given image and DQ extension with the given satellite trails mask and flag.

Parameters:
  • filename (str) – FITS image filename to update.
  • ext (int, str, or tuple) – DQ extension, as accepted by astropy.io.fits, to update.
  • mask (ndarray) – Boolean mask, with True marking the satellite trail(s). This can be the result(s) from make_mask().
  • dqval (int, optional) – DQ value to use for the trail. Default value of 16384 is tailored for ACS/WFC.
  • verbose (bool, optional) – Print extra information to the terminal.

Version

acstools.satdet.__version__ = u'0.3.3'

unicode(object=’‘) -> unicode object unicode(string[, encoding[, errors]]) -> unicode object

Create a new Unicode object from the given encoded string. encoding defaults to the current default string encoding. errors can be ‘strict’, ‘replace’ or ‘ignore’ and defaults to ‘strict’.

acstools.satdet.__vdate__ = u'11-Jul-2017'

unicode(object=’‘) -> unicode object unicode(string[, encoding[, errors]]) -> unicode object

Create a new Unicode object from the given encoded string. encoding defaults to the current default string encoding. errors can be ‘strict’, ‘replace’ or ‘ignore’ and defaults to ‘strict’.

Author

acstools.satdet.__author__ = u'David Borncamp, Pey Lian Lim'

unicode(object=’‘) -> unicode object unicode(string[, encoding[, errors]]) -> unicode object

Create a new Unicode object from the given encoded string. encoding defaults to the current default string encoding. errors can be ‘strict’, ‘replace’ or ‘ignore’ and defaults to ‘strict’.