Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

No more centroiding in Simple Aperture Photometry plugin #1841

Merged
merged 3 commits into from
Nov 18, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 5 additions & 0 deletions CHANGES.rst
Original file line number Diff line number Diff line change
Expand Up @@ -52,6 +52,11 @@ Cubeviz
Imviz
^^^^^

- Simple Aperture Photometry plugin no longer performs centroiding.
For radial profile, curve of growth, and table reporting, the aperture
center is used instead. For centroiding, use "Recenter" feature in
the Subset Tools plugin. [#1841]
kecnry marked this conversation as resolved.
Show resolved Hide resolved

Mosviz
^^^^^^

Expand Down
12 changes: 4 additions & 8 deletions docs/imviz/export_data.rst
Original file line number Diff line number Diff line change
Expand Up @@ -51,11 +51,8 @@ The columns are as follow:

* :attr:`~photutils.aperture.ApertureStats.id`: ID number assigned to the row,
starting from 1.
* :attr:`~photutils.aperture.ApertureStats.xcentroid`,
:attr:`~photutils.aperture.ApertureStats.ycentroid`: Pixel centroids
calculated using moments. This might differ from center of the aperture.
* :attr:`~photutils.aperture.ApertureStats.sky_centroid`:
`~astropy.coordinates.SkyCoord` associated with the centroid.
* ``xcenter``, ``ycenter``: Center of the aperture (0-indexed).
* ``sky_center``: `~astropy.coordinates.SkyCoord` associated with the center.
If WCS is not available, this field is `None`.
* ``background``: The value from :guilabel:`Background value`, with unit attached.
* :attr:`~photutils.aperture.ApertureStats.sum`: Sum of flux in the aperture.
Expand Down Expand Up @@ -108,9 +105,8 @@ The columns are as follow:
.. note::

Aperture sum and statistics are done on the originally drawn aperture only.
Even though centroid is calculated, it is not used to move the aperture
to the new center. However, radial profiles (including Gaussian fitting, if any)
and curve of growth do use the centroid as zero-point on the X-axis.
You can use the :ref:`imviz-subset-plugin` plugin to center it first on the
object of interest, if you wish.

Once you have the results in a table, you can further manipulated them as
documented in :ref:`astropy:astropy-table`.
6 changes: 5 additions & 1 deletion docs/imviz/plugins.rst
Original file line number Diff line number Diff line change
Expand Up @@ -135,6 +135,10 @@ an interactively selected region. A typical workflow is as follows:
2. Draw a region over the object of interest (see :ref:`imviz_defining_spatial_regions`).
3. Select the desired image using the :guilabel:`Data` dropdown menu.
4. Select the desired region using the :guilabel:`Subset` dropdown menu.
You can use the :ref:`imviz-subset-plugin` plugin to center it first on the
object of interest using its center of mass, if you wish.
Depending on the object, it may take several iterations for re-centering
to converge, or it may never converge at all.
5. If you want to subtract background before performing photometry,
you have the following 3 options. Otherwise if your image is already
background subtracted, choose "Manual" and leave the background set at 0:
Expand Down Expand Up @@ -232,7 +236,7 @@ catalog dropdown menu.
This plugin is still under active development. As a result, the search only uses the SDSS DR17 catalog
and works best when you only have a single image loaded in a viewer.

To load a catalog from a supported `JWST ECSV catalog file <https://jwst-pipeline.readthedocs.io/en/latest/jwst/source_catalog/main.html#output-products>`_, choose "From File...".
To load a catalog from a supported `JWST ECSV catalog file <https://jwst-pipeline.readthedocs.io/en/latest/jwst/source_catalog/main.html#output-products>`_, choose "From File...".
The file must be able to be parsed by `astropy.table.Table.read` and contain a column labeled 'sky_centroid'.
Clicking :guilabel:`SEARCH` will show markers for any entry within the filtered zoom window.

Expand Down
8 changes: 0 additions & 8 deletions docs/known_bugs.rst
Original file line number Diff line number Diff line change
Expand Up @@ -143,14 +143,6 @@ to show the markers, or not at all. This is a known bug reported in
https://github.com/glue-viz/glue-jupyter/issues/243 . If you encounter this,
try a different OS/browser combo.

Simple Aperture Photometry Plugin and dithered images
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

Due to a known bug reported in https://github.com/glue-viz/glue-astronomy/issues/52 ,
aperture photometry and radial profile will report inaccurate results when you
calculate them on dithered images linked by WCS *unless* you are on the reference image
(this is usually the first loaded image).

.. _known_issues_specviz:

Specviz
Expand Down
50 changes: 28 additions & 22 deletions jdaviz/configs/imviz/plugins/aper_phot_simple/aper_phot_simple.py
Original file line number Diff line number Diff line change
Expand Up @@ -245,6 +245,8 @@ def vue_do_aper_phot(self, *args, **kwargs):

data = self._selected_data
reg = self._selected_subset
xcenter = reg.center.x
ycenter = reg.center.y

# Reset last fitted model
fit_model = None
Expand All @@ -257,6 +259,10 @@ def vue_do_aper_phot(self, *args, **kwargs):
bg = float(self.background_value)
except ValueError: # Clearer error message
raise ValueError('Missing or invalid background value')
if data.coords is not None:
sky_center = data.coords.pixel_to_world(xcenter, ycenter)
else:
sky_center = None
aperture = regions2aperture(reg)
include_pixarea_fac = False
include_counts_fac = False
Expand Down Expand Up @@ -289,12 +295,10 @@ def vue_do_aper_phot(self, *args, **kwargs):
include_flux_scale = True
phot_aperstats = ApertureStats(comp_data, aperture, wcs=data.coords, local_bkg=bg)
phot_table = phot_aperstats.to_table(columns=(
'id', 'xcentroid', 'ycentroid', 'sky_centroid', 'sum', 'sum_aper_area',
'id', 'sum', 'sum_aper_area',
'min', 'max', 'mean', 'median', 'mode', 'std', 'mad_std', 'var',
'biweight_location', 'biweight_midvariance', 'fwhm', 'semimajor_sigma',
'semiminor_sigma', 'orientation', 'eccentricity')) # Some cols excluded, add back as needed. # noqa
phot_table['xcentroid'].unit = u.pix # photutils only assumes, we make it real
phot_table['ycentroid'].unit = u.pix
rawsum = phot_table['sum'][0]

if include_pixarea_fac:
Expand Down Expand Up @@ -323,12 +327,14 @@ def vue_do_aper_phot(self, *args, **kwargs):

# Extra info beyond photutils.
phot_table.add_columns(
[bg, pixarea_fac, sum_ct, sum_ct_err, ctfac, sum_mag, flux_scale, data.label,
[xcenter * u.pix, ycenter * u.pix, sky_center,
bg, pixarea_fac, sum_ct, sum_ct_err, ctfac, sum_mag, flux_scale, data.label,
reg.meta.get('label', ''), Time(datetime.utcnow())],
names=['background', 'pixarea_tot', 'aperture_sum_counts',
'aperture_sum_counts_err', 'counts_fac', 'aperture_sum_mag', 'flux_scaling',
names=['xcenter', 'ycenter', 'sky_center', 'background', 'pixarea_tot',
'aperture_sum_counts', 'aperture_sum_counts_err', 'counts_fac',
'aperture_sum_mag', 'flux_scaling',
'data_label', 'subset_label', 'timestamp'],
indexes=[4, 6, 6, 6, 6, 6, 6, 21, 21, 21])
indexes=[1, 1, 1, 1, 3, 3, 3, 3, 3, 3, 18, 18, 18])

# Attach to app for Python extraction.
if (not hasattr(self.app, '_aper_phot_results') or
Expand All @@ -352,9 +358,9 @@ def vue_do_aper_phot(self, *args, **kwargs):
line_y_sc = bqplot.LinearScale()

if self.current_plot_type == "Curve of Growth":
self._fig.title = 'Curve of growth from source centroid'
self._fig.title = 'Curve of growth from aperture center'
x_arr, sum_arr, x_label, y_label = _curve_of_growth(
comp_data, phot_aperstats.centroid, aperture, phot_table['sum'][0],
comp_data, (xcenter, ycenter), aperture, phot_table['sum'][0],
wcs=data.coords, background=bg, pixarea_fac=pixarea_fac)
self._fig.axes = [bqplot.Axis(scale=line_x_sc, label=x_label),
bqplot.Axis(scale=line_y_sc, orientation='vertical',
Expand All @@ -370,33 +376,33 @@ def vue_do_aper_phot(self, *args, **kwargs):
label=comp.units or 'Value')]

if self.current_plot_type == "Radial Profile":
self._fig.title = 'Radial profile from source centroid'
self._fig.title = 'Radial profile from aperture center'
x_data, y_data = _radial_profile(
phot_aperstats.data_cutout, phot_aperstats.bbox, phot_aperstats.centroid,
phot_aperstats.data_cutout, phot_aperstats.bbox, (xcenter, ycenter),
raw=False)
bqplot_line = bqplot.Lines(x=x_data, y=y_data, marker='circle',
scales={'x': line_x_sc, 'y': line_y_sc},
marker_size=32, colors='gray')
else: # Radial Profile (Raw)
self._fig.title = 'Raw radial profile from source centroid'
self._fig.title = 'Raw radial profile from aperture center'
x_data, y_data = _radial_profile(
phot_aperstats.data_cutout, phot_aperstats.bbox, phot_aperstats.centroid,
phot_aperstats.data_cutout, phot_aperstats.bbox, (xcenter, ycenter),
raw=True)
bqplot_line = bqplot.Scatter(x=x_data, y=y_data, marker='circle',
scales={'x': line_x_sc, 'y': line_y_sc},
default_size=1, colors='gray')

# Fit Gaussian1D to radial profile data.
# mean is fixed at 0 because we recentered to centroid.
if self.fit_radial_profile:
fitter = LevMarLSQFitter()
y_max = y_data.max()
x_mean = x_data[np.where(y_data == y_max)].mean()
std = 0.5 * (phot_table['semimajor_sigma'][0] +
phot_table['semiminor_sigma'][0])
if isinstance(std, u.Quantity):
std = std.value
gs = Gaussian1D(amplitude=y_max, mean=0, stddev=std,
fixed={'mean': True, 'amplitude': True},
gs = Gaussian1D(amplitude=y_max, mean=x_mean, stddev=std,
fixed={'amplitude': True},
Comment on lines +399 to +405
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

seems to work reasonably well for an intentionally non-centered aperture!

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It still isn't fool-proof (it will fail if y_max returns x_data values at different peaks, for instance) but hopefully work good enough for most cases that this plugin is useful for... 🤞

bounds={'amplitude': (y_max * 0.5, y_max)})
if Version(astropy.__version__) < Version('5.2'):
fitter_kw = {}
Expand Down Expand Up @@ -433,13 +439,13 @@ def vue_do_aper_phot(self, *args, **kwargs):
continue
x = phot_table[key][0]
if (isinstance(x, (int, float, u.Quantity)) and
key not in ('xcentroid', 'ycentroid', 'sky_centroid', 'sum_aper_area',
key not in ('xcenter', 'ycenter', 'sky_center', 'sum_aper_area',
'aperture_sum_counts', 'aperture_sum_mag')):
tmp.append({'function': key, 'result': f'{x:.4e}'})
elif key == 'sky_centroid' and x is not None:
tmp.append({'function': 'RA centroid', 'result': f'{x.ra.deg:.6f} deg'})
tmp.append({'function': 'Dec centroid', 'result': f'{x.dec.deg:.6f} deg'})
elif key in ('xcentroid', 'ycentroid', 'sum_aper_area'):
elif key == 'sky_center' and x is not None:
tmp.append({'function': 'RA center', 'result': f'{x.ra.deg:.6f} deg'})
tmp.append({'function': 'Dec center', 'result': f'{x.dec.deg:.6f} deg'})
elif key in ('xcenter', 'ycenter', 'sum_aper_area'):
tmp.append({'function': key, 'result': f'{x:.1f}'})
elif key == 'aperture_sum_counts' and x is not None:
tmp.append({'function': key, 'result':
Expand All @@ -452,7 +458,7 @@ def vue_do_aper_phot(self, *args, **kwargs):
# Also display fit results
fit_tmp = []
if fit_model is not None and isinstance(fit_model, Gaussian1D):
for param in ('fwhm', 'amplitude'): # mean is fixed at 0
for param in ('mean', 'fwhm', 'amplitude'):
p_val = getattr(fit_model, param)
if isinstance(p_val, Parameter):
p_val = p_val.value
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -141,7 +141,7 @@
<jupyter-widget :widget="radial_plot" style="width: 100%; height: 480px" />
</v-row>

<div v-if="plot_available && fit_radial_profile">
<div v-if="plot_available && fit_radial_profile && current_plot_type != 'Curve of Growth'">
<j-plugin-section-header>Gaussian Fit Results</j-plugin-section-header>
<v-row no-gutters>
<v-col cols=6><U>Result</U></v-col>
Expand Down
14 changes: 7 additions & 7 deletions jdaviz/configs/imviz/tests/test_parser.py
Original file line number Diff line number Diff line change
Expand Up @@ -277,11 +277,11 @@ def test_parse_jwst_nircam_level2(self, imviz_helper):
assert_allclose(phot_plugin.flux_scaling, 0.003631)
phot_plugin.vue_do_aper_phot()
tbl = imviz_helper.get_aperture_photometry_results()
assert_quantity_allclose(tbl[0]['xcentroid'], 970.935492 * u.pix)
assert_quantity_allclose(tbl[0]['ycentroid'], 1116.967619 * u.pix)
sky = tbl[0]['sky_centroid']
assert_quantity_allclose(tbl[0]['xcenter'], 970.95 * u.pix)
assert_quantity_allclose(tbl[0]['ycenter'], 1116.05 * u.pix)
sky = tbl[0]['sky_center']
assert_allclose(sky.ra.deg, 80.48419863)
assert_allclose(sky.dec.deg, -69.494592)
assert_allclose(sky.dec.deg, -69.494608)
data_unit = u.MJy / u.sr
assert_quantity_allclose(tbl[0]['background'], 0.1741226315498352 * data_unit)
assert_quantity_allclose(tbl[0]['sum'], 4.486487e-11 * u.MJy, rtol=1e-6)
Expand Down Expand Up @@ -402,9 +402,9 @@ def test_parse_hst_drz(self, imviz_helper):
assert_allclose(phot_plugin.pixel_area, 0.0025) # Not used but still auto-populated
phot_plugin.vue_do_aper_phot()
tbl = imviz_helper.get_aperture_photometry_results()
assert_quantity_allclose(tbl[0]['xcentroid'], 1487.60825422 * u.pix, atol=2 * u.pix)
assert_quantity_allclose(tbl[0]['ycentroid'], 2573.83983184 * u.pix, atol=2 * u.pix)
sky = tbl[0]['sky_centroid']
assert_quantity_allclose(tbl[0]['xcenter'], 1488.5 * u.pix, atol=2 * u.pix)
assert_quantity_allclose(tbl[0]['ycenter'], 2576 * u.pix, atol=2 * u.pix)
sky = tbl[0]['sky_center']
assert_allclose(sky.ra.deg, 3.684062989070131, rtol=1e-3)
assert_allclose(sky.dec.deg, 10.802045612042956, rtol=1e-3)
data_unit = u.electron / u.s
Expand Down
34 changes: 19 additions & 15 deletions jdaviz/configs/imviz/tests/test_simple_aper_phot.py
Original file line number Diff line number Diff line change
Expand Up @@ -65,7 +65,7 @@ def test_plugin_wcs_dithered(self):

# Check photometry results.
assert tbl.colnames == [
'id', 'xcentroid', 'ycentroid', 'sky_centroid', 'background', 'sum',
'id', 'xcenter', 'ycenter', 'sky_center', 'background', 'sum',
'sum_aper_area', 'pixarea_tot', 'aperture_sum_counts', 'aperture_sum_counts_err',
'counts_fac', 'aperture_sum_mag', 'flux_scaling', 'min', 'max', 'mean', 'median',
'mode', 'std', 'mad_std', 'var', 'biweight_location', 'biweight_midvariance',
Expand Down Expand Up @@ -99,11 +99,11 @@ def test_plugin_wcs_dithered(self):
assert_array_equal(tbl['subset_label'], ['Subset 1', 'Subset 1'])
assert tbl['timestamp'].scale == 'utc'

# Sky is the same but xcentroid different due to dithering.
# Sky is the same but xcenter different due to dithering.
# The aperture sum is different too because mask is a little off limit in second image.
assert_quantity_allclose(tbl['xcentroid'], [4.5, 5.5] * u.pix)
assert_quantity_allclose(tbl['ycentroid'], 4.5 * u.pix)
sky = tbl['sky_centroid']
assert_quantity_allclose(tbl['xcenter'], [4.5, 5.5] * u.pix)
assert_quantity_allclose(tbl['ycenter'], 4.5 * u.pix)
sky = tbl['sky_center']
assert_allclose(sky.ra.deg, 337.518943)
assert_allclose(sky.dec.deg, -20.832083)
assert_allclose(tbl['sum'], [63.617251, 62.22684693104279])
Expand All @@ -117,9 +117,9 @@ def test_plugin_wcs_dithered(self):
tbl = self.imviz.get_aperture_photometry_results()
assert len(tbl) == 3 # New result is appended
assert tbl[-1]['id'] == 3
assert_quantity_allclose(tbl[-1]['xcentroid'], 4.5 * u.pix)
assert_quantity_allclose(tbl[-1]['ycentroid'], 2 * u.pix)
sky = tbl[-1]['sky_centroid']
assert_quantity_allclose(tbl[-1]['xcenter'], 4.5 * u.pix)
assert_quantity_allclose(tbl[-1]['ycenter'], 2 * u.pix)
sky = tbl[-1]['sky_center']
assert_allclose(sky.ra.deg, 337.51894336144454)
assert_allclose(sky.dec.deg, -20.832777499255897)
assert_quantity_allclose(tbl[-1]['sum_aper_area'], 28.274334 * (u.pix * u.pix))
Expand All @@ -139,11 +139,11 @@ def test_plugin_wcs_dithered(self):
tbl = self.imviz.get_aperture_photometry_results()
assert len(tbl) == 4 # New result is appended
assert tbl[-1]['id'] == 4
assert np.isnan(tbl[-1]['xcentroid'])
assert np.isnan(tbl[-1]['ycentroid'])
sky = tbl[-1]['sky_centroid']
assert np.isnan(sky.ra.deg)
assert np.isnan(sky.dec.deg)
assert_quantity_allclose(tbl[-1]['xcenter'], 4.5 * u.pix)
assert_quantity_allclose(tbl[-1]['ycenter'], 4.5 * u.pix)
sky = tbl[-1]['sky_center']
assert_allclose(sky.ra.deg, 337.51894336144454)
assert_allclose(sky.dec.deg, -20.832083)
assert_quantity_allclose(tbl[-1]['sum_aper_area'], 81 * (u.pix * u.pix))
assert_allclose(tbl[-1]['sum'], 0)
assert_allclose(tbl[-1]['mean'], 0)
Expand All @@ -168,7 +168,7 @@ def test_plugin_wcs_dithered(self):
# Curve of growth
phot_plugin.current_plot_type = 'Curve of Growth'
phot_plugin.vue_do_aper_phot()
assert phot_plugin._fig.title == 'Curve of growth from source centroid'
assert phot_plugin._fig.title == 'Curve of growth from aperture center'


class TestSimpleAperPhot_NoWCS(BaseImviz_WCS_NoWCS):
Expand All @@ -187,7 +187,7 @@ def test_plugin_no_wcs(self):
phot_plugin.vue_do_aper_phot()
tbl = self.imviz.get_aperture_photometry_results()
assert len(tbl) == 1 # Old table discarded due to incompatible column
assert_array_equal(tbl['sky_centroid'], None)
assert_array_equal(tbl['sky_center'], None)


def test_annulus_background(imviz_helper):
Expand Down Expand Up @@ -263,6 +263,8 @@ def test_annulus_background(imviz_helper):

# NOTE: Extracting the cutout for radial profile is aperture
# shape agnostic, so we use ellipse as representative case.
# NOTE: This test only tests the radial profile algorithm and does
# not care if the actual plugin use centroid or not.
class TestRadialProfile():
def setup_class(self):
data = np.ones((51, 51)) * u.nJy
Expand All @@ -286,6 +288,8 @@ def test_profile_imexam(self):
assert_allclose(y_arr, 1)


# NOTE: This test only tests the curve of growth algorithm and does
# not care if the actual plugin use centroid or not.
@pytest.mark.parametrize('with_unit', (False, True))
def test_curve_of_growth(with_unit):
data = np.ones((51, 51))
Expand Down