Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
for i in range(number_apertures-1):
if sens == 'clock':
xx[i+1] = cosangle*xx[i] + sinangle*yy[i]
yy[i+1] = cosangle*yy[i] - sinangle*xx[i]
elif sens == 'counterclock':
xx[i+1] = cosangle*xx[i] - sinangle*yy[i]
yy[i+1] = cosangle*yy[i] + sinangle*xx[i]
xx += centerx
yy += centery
rad = fwhm/2.
if exclude_negative_lobes:
xx = np.concatenate(([xx[0]], xx[2:-1]))
yy = np.concatenate(([yy[0]], yy[2:-1]))
apertures = photutils.CircularAperture((xx, yy), r=rad) # Coordinates (X,Y)
fluxes = photutils.aperture_photometry(array, apertures, method='exact')
fluxes = np.array(fluxes['aperture_sum'])
if array2 is not None:
fluxes2 = photutils.aperture_photometry(array2, apertures,
method='exact')
fluxes2 = np.array(fluxes2['aperture_sum'])
if use2alone:
fluxes = np.concatenate(([fluxes[0]], fluxes2[:]))
else:
fluxes = np.concatenate((fluxes, fluxes2))
f_source = fluxes[0].copy()
fluxes = fluxes[1:]
n2 = fluxes.shape[0]
backgr_apertures_std = fluxes.std(ddof=1)
def aperture(r, pos):
""" Creates circular & rectangular apertures of given size """
circ = CircularAperture(pos, r)
rect = RectangularAperture(pos, r, r, 0.0)
return circ, rect
'false positive fraction.')
# Initialize a numpy array in which we will store the integrated flux of all reference apertures
ap_phot = np.zeros(num_ap)
# Loop over all reference apertures and measure the integrated flux
for i, theta in enumerate(ap_theta):
# Compute the position of the current aperture in polar coordinates and convert to Cartesian
x_tmp = center[1] + (x_pos - center[1]) * math.cos(theta) - \
(y_pos - center[0]) * math.sin(theta)
y_tmp = center[0] + (x_pos - center[1]) * math.sin(theta) + \
(y_pos - center[0]) * math.cos(theta)
# Place a circular aperture at this position and sum up the flux inside the aperture
aperture = CircularAperture((x_tmp, y_tmp), size)
phot_table = aperture_photometry(image, aperture, method='exact')
ap_phot[i] = phot_table['aperture_sum']
# Define shortcuts to the signal and the noise aperture sums
signal_aperture = ap_phot[0]
noise_apertures = ap_phot[1:]
# Compute the "signal", that is, the numerator of the signal-to-noise ratio: According to
# eq. (8) in Mawet et al. (2014), this is given by the difference between the integrated flux
# in the signal aperture and the mean of the integrated flux in the noise apertures
signal = signal_aperture - np.mean(noise_apertures)
# Compute the "noise", that is, the denominator of the signal-to-noise-ratio: According to
# eq. (8) in Mawet et al. (2014), this is given by the standard deviation of the integrated flux
# in the noise apertures times a correction factor to account for the small sample statistics.
# NOTE: `ddof=1` is a necessary argument for np.std() in order to compute the *unbiased*
def aperture(r, pos):
circ = CircularAperture(pos, r)
rect = RectangularAperture(pos, r, r, 0.0)
return circ, rect
depending on whether its center is in or out of the aperture.
'subpixel': A pixel is divided into subpixels and the center of each
subpixel is tested (as above).
'exact': (default) The exact overlap between the aperture and each pixel is
calculated.
"""
n_obj = len(yc)
flux = np.zeros((n_obj))
for i, (y, x) in enumerate(zip(yc, xc)):
if mean:
ind = circle(y, x, (ap_factor*fwhm)/2)
values = array[ind]
obj_flux = np.mean(values)
else:
aper = photutils.CircularAperture((x, y), (ap_factor*fwhm)/2)
obj_flux = photutils.aperture_photometry(array, aper,
method='exact')
obj_flux = np.array(obj_flux['aperture_sum'])
flux[i] = obj_flux
if verbose:
print('Coordinates of object {} : ({},{})'.format(i, y, x))
print('Object Flux = {:.2f}'.format(flux[i]))
return flux
if init_rad is None:
init_rad = fwhm
if debug:
_, ax = plt.subplots(figsize=(6, 6))
ax.imshow(array, origin='lower', interpolation='nearest',
alpha=0.5, cmap='gray')
for i in range(n_annuli):
y = centery + init_rad + separation * i
rad = dist(centery, centerx, y, centerx)
yy, xx = find_coords(rad, fwhm, init_angle, fin_angle)
yy += centery
xx += centerx
apertures = photutils.CircularAperture((xx, yy), fwhm/2)
fluxes = photutils.aperture_photometry(array, apertures)
fluxes = np.array(fluxes['aperture_sum'])
noise_ann = np.std(fluxes)
noise.append(noise_ann)
vector_radd.append(rad)
if debug:
for j in range(xx.shape[0]):
# Circle takes coordinates as (X,Y)
aper = plt.Circle((xx[j], yy[j]), radius=fwhm/2, color='r',
fill=False, alpha=0.8)
ax.add_patch(aper)
cent = plt.Circle((xx[j], yy[j]), radius=0.8, color='r',
fill=True, alpha=0.5)
ax.add_patch(cent)
depending on whether its center is in or out of the aperture.
'subpixel': A pixel is divided into subpixels and the center of each
subpixel is tested (as above).
'exact': (default) The exact overlap between the aperture and each pixel is
calculated.
"""
n_obj = len(yc)
flux = np.zeros((n_obj))
for i, (y, x) in enumerate(zip(yc, xc)):
if mean:
ind = circle(y, x, (ap_factor*fwhm)/2)
values = array[ind]
obj_flux = np.mean(values)
else:
aper = photutils.CircularAperture((x, y), (ap_factor*fwhm)/2)
obj_flux = photutils.aperture_photometry(array, aper,
method='exact')
obj_flux = np.array(obj_flux['aperture_sum'])
flux[i] = obj_flux
if verbose:
print('Coordinates of object {} : ({},{})'.format(i, y, x))
print('Object Flux = {:.2f}'.format(flux[i]))
return flux
psf_e, _ = psf.generate_aureole(psf_range=2 * psf_size)
### Hybrid radius (in pixel)
try:
hybrid_r = config.starhalo.hybrid_r
except:
hybrid_r = 12
### Inner PSF: from stacking stars
inner_psf = copy.deepcopy(median_psf)
inner_psf /= np.sum(inner_psf) # Normalize
inner_size = inner_psf.shape
inner_cen = [int(x / 2) for x in inner_size]
##### flux_inn is the flux inside an annulus, we use this to scale inner and outer parts
flux_inn = compute_Rnorm(inner_psf, None, inner_cen, R=hybrid_r, display=False, mask_cross=False)[1]
##### We only remain the stacked PSF inside hybrid radius.
aper = CircularAperture(inner_cen, hybrid_r).to_mask()
mask = aper.to_image(inner_size) == 0
inner_psf[mask] = np.nan
### Make new empty PSF
outer_cen = (int(psf_size / 2), int(psf_size / 2))
new_psf = np.zeros((int(psf_size), int(psf_size)))
new_psf[outer_cen[0] - inner_cen[0]:outer_cen[0] + inner_cen[0] + 1,
outer_cen[1] - inner_cen[1]:outer_cen[1] + inner_cen[1] + 1] = inner_psf
### Outer PSF: from model
outer_psf = psf_e.drawImage(nx=psf_size, ny=psf_size, scale=config.lowres.pixel_scale, method="no_pixel").array
outer_psf /= np.sum(outer_psf) # Normalize
##### flux_out is the flux inside an annulus, we use this to scale inner and outer parts
flux_out = compute_Rnorm(outer_psf, None, outer_cen,
R=hybrid_r, display=False, mask_cross=False)[1]
# Determine the maximum value in the box and the minium value for plotting
if vmin is None: vmin = max(np.nanmin(box), 0.)
if vmax is None: vmax = 0.5 * (np.nanmax(box) + vmin)
# Set the normalization
norm = ImageNormalize(stretch=SqrtStretch())
# Make the plot
plt.figure(figsize=(8,2.5))
plt.imshow(box, origin='lower', norm=norm, interpolation='nearest', vmin=vmin, vmax=vmax, cmap="viridis")
if radius is None: plt.plot(x_peaks, y_peaks, ls='none', color='white', marker='+', ms=40, lw=10, mew=4)
else:
positions = (x_peaks, y_peaks)
apertures = CircularAperture(positions, r=radius)
apertures.plot(color='green', lw=1.5, alpha=0.5)
plt.xlim(0, box.shape[1]-1)
plt.ylim(0, box.shape[0]-1)
if title is not None: plt.title(title)
plt.show()