1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
|
Aperture photometry
===================
There are four aperture functions available:
================== =========================
Function Sums data within...
================== =========================
`sep.sum_circle` circle(s)
`sep.sum_circann` circular annulus/annuli
`sep.sum_ellipse` ellipse(s)
`sep.sum_ellipann` elliptical annulus/annuli
================== =========================
The follow examples demonstrate options for circular aperture
photometry. The other functions behave similarly.
.. code-block:: python
# sum flux in circles of radius=3.0
flux, fluxerr, flag = sep.sum_circle(data, objs['x'], objs['y'], 3.0)
# x, y and r can be arrays and obey numpy broadcasting rules.
# Here, r is an array instead of a single number:
flux, fluxerr, flag = sep.sum_circle(data, objs['x'], objs['y'],
3.0 * np.ones(len(objs)))
# use a different subpixel sampling (default is 5; 0 means "exact")
flux, fluxerr, flag = sep.sum_circle(data, objs['x'], objs['y'], 3.0,
subpix=0)
**Error calculation**
In the default modes illustrated above, the uncertainty ``fluxerr`` is
not calculated and values of 0 are simply returned. The uncertainty can be
flexibly and efficiently calculated, depending on the characteristics
of your data. The presence of an ``err`` or ``var`` keyword indicates
a per-pixel noise, while the presense of a ``gain`` keyword indicates
that the Poisson uncertainty on the total sum should be included. Some
illustrative examples:
.. code-block:: python
# Specify a per-pixel "background" error and a gain. This is suitable
# when the data have been background subtracted.
flux, fluxerr, flag = sep.sum_circle(data, objs['x'], objs['y'], 3.0,
err=bkg.globalrms, gain=1.0)
# Variance can be passed instead of error, with identical results.
flux, fluxerr, flag = sep.sum_circle(data, objs['x'], objs['y'], 3.0,
var=bkg.globalrms**2, gain=1.0)
# Error or variance can be arrays, indicating that the background error
# is variable.
bkgrms = bkg.rms() # array, same shape as data
flux, fluxerr, flag = sep.sum_circle(data, objs['x'], objs['y'], 3.0,
err=bkgrms, gain=1.0)
# If your uncertainty array already includes Poisson noise from the object,
# leave gain as None (default):
flux, fluxerr, flag = sep.sum_circle(data, objs['x'], objs['y'], 3.0,
err=error_array)
# If your data represent raw counts (it is not background-subtracted),
# set only gain to get the poisson error:
flux, fluxerr, flag = sep.sum_circle(data, objs['x'], objs['y'], 3.0,
gain=1.0)
The error is calculated as
.. math::
\sigma_F^2 = \sum_i \sigma_i^2 + F/g
where the sum is over pixels in the aperture, :math:`\sigma_i` is the
noise in each pixel, :math:`F` is the sum in the aperture and
:math:`g` is the gain. The last term is not added if ``gain`` is
`None`.
**Masking**
Apply a mask (same shape as data). Pixels where the mask is True are
"corrected" to the average value within the aperture.
.. code-block:: python
flux, fluxerr, flag = sep.sum_circle(data, objs['x'], objs['y'], 3.0,
mask=mask)
**Local background subtraction**
The `~sep.sum_circle` and `~sep.sum_ellipse` functions have options
for performing local background subtraction. For example, to subtract the
background calculated in an annulus between 6 and 8 pixel radius:
.. code-block:: python
flux, fluxerr, flag = sep.sum_circle(data, objs['x'], objs['y'], 3.0,
mask=mask, bkgann=(6., 8.))
Pixels in the background annulus are not subsampled and any masked
pixels in the annulus are completely igored rather than corrected.
The inner and outer radii can also be arrays. The error in the background
is included in the reported error.
Equivalent of FLUX_AUTO (e.g., MAG_AUTO) in Source Extractor
------------------------------------------------------------
This is a two-step process. First we calculate the Kron radius for each
object, then we perform elliptical aperture photometry within that radius:
.. code-block:: python
kronrad, krflag = sep.kron_radius(data, x, y, a, b, theta, 6.0)
flux, fluxerr, flag = sep.sum_ellipse(data, x, y, a, b, theta, 2.5*kronrad,
subpix=1)
flag |= krflag # combine flags into 'flag'
This specific example is the equilvalent of setting ``PHOT_AUTOPARAMS
2.5, 0.0`` in Source Extractor (note the 2.5 in the argument to
``sep.sum_ellipse``). The second Source Extractor parameter is a
minimum diameter. To replicate Source Extractor behavior for non-zero
values of the minimum diameter, one would put in logic to use circular
aperture photometry if the Kron radius is too small. For example:
.. code-block:: python
r_min = 1.75 # minimum diameter = 3.5
use_circle = kronrad * np.sqrt(a * b) < r_min
cflux, cfluxerr, cflag = sep.sum_circle(data, x[use_circle], y[use_circle],
r_min, subpix=1)
flux[use_circle] = cflux
fluxerr[use_circle] = cfluxerr
flag[use_circle] = cflag
.. warning::
Caution should be used when calculating Kron radii in crowded fields.
In almost all cases, one would want to pass in a segmentation map to
mask out nearby objects, as described below in
:ref:`segmentation masking`.
Equivalent of FLUX_RADIUS in Source Extractor
---------------------------------------------
In Source Extractor, the FLUX_RADIUS parameter gives the radius of a
circle enclosing a desired fraction of the total flux. For example,
with the setting ``PHOT_FLUXFRAC 0.5``, FLUX_RADIUS will give the
radius of a circle containing half the "total flux" of the object. For
the definition of "total flux", Source Extractor uses its measurement
of FLUX_AUTO, which is taken through an elliptical aperture (see
above). Thus, with the setting ``PHOT_FLUXFRAC 1.0``, you would find
the circle containing the same flux as whatever ellipse Source
Extractor used for ``FLUX_AUTO``.
Given a previous calculation of ``flux`` as above, calculate the
radius for a flux fraction of 0.5:
.. code-block:: python
r, flag = sep.flux_radius(data, objs['x'], objs['y'], 6.*objs['a'], 0.5,
normflux=flux, subpix=5)
And for multiple flux fractions:
.. code-block:: python
r, flag = sep.flux_radius(data, objs['x'], objs['y'], 6.*objs['a'],
[0.5, 0.6], normflux=flux, subpix=5)
Equivalent of XWIN_IMAGE, YWIN_IMAGE in Source Extractor
--------------------------------------------------------
Source Extractor's XWIN_IMAGE, YWIN_IMAGE parameters can be used for
more accurate object centroids than the default X_IMAGE, Y_IMAGE.
Here, the ``winpos`` function provides this behavior. To match Source
Extractor exactly, the right ``sig`` parameter (giving a description
of the effective width) must be used for each object. Source
Extractor uses ``2. / 2.35 * (half-light radius)`` where the
half-light radius is calculated using ``flux_radius`` with a fraction
of 0.5 and a normalizing flux of ``FLUX_AUTO``. The equivalent here is:
.. code-block:: python
sig = 2. / 2.35 * r # r from sep.flux_radius() above, with fluxfrac = 0.5
xwin, ywin, flag = sep.winpos(data, objs['x'], objs['y'], sig)
.. _segmentation masking:
Segmentation-masked image measurements
--------------------------------------
SourceExtractor provides a mechanism for measuring the "AUTO" and
"FLUX_RADIUS" parameters for a given object including a mask for
neighboring sources. In addition to the mask, setting the SourceExtractor
parameter ``MASK_TYPE=CORRECT`` further fills the masked pixels of a given
source with "good" pixel values reflected opposite of the masked pixels.
The ``SEP`` photometry and measurement functions provide an option for
simple masking without reflection or subtracting neighbor flux.
For example, using a segmentation array provided by ``sep.extract``, we
can compute the masked ``flux_radius`` that could otherwise be
artificially large due to flux from nearby sources:
.. code-block:: python
# list of object id numbers that correspond to the segments
seg_id = np.arange(1, len(objs)+1, dtype=np.int32)
r, flag = sep.flux_radius(data, objs['x'], objs['y'], 6.*objs['a'],
0.5, seg_id=seg_id, seg=seg,
normflux=flux, subpix=5)
To enforce that a given measurement **only** includes pixels within a
segment, provide negative values in the ``seg_id`` list. Otherwise the
mask for a given object will be pixels with
``(seg == 0) | (seg_id == id_i)``.
The following functions include the segmentation masking:
``sum_circle``, ``sum_circann``, ``sum_ellipse``, ``sum_ellipann``,
``flux_radius`` , and ``kron_radius`` (``winpos`` **currently does not**).
Masking image regions
---------------------
Create a boolean array with elliptical regions set to True:
.. code-block:: python
mask = np.zeros(data.shape, dtype=np.bool)
sep.mask_ellipse(mask, objs['x'], objs['y'], obs['a'], objs['b'],
objs['theta'], r=3.)
|
