"""
Overview
--------
This module provides the implementation of the :class:`Table`
Classes and Inheritance Structure
----------------------------------------------
.. inheritance-diagram::
DoDENCLUEDeblendingTask
mary::
Table
Summary
---------
.. autosummary::
:nosignatures:
DoDENCLUEDeblendingTask
do_denclue_deblending
do_denclue_of_cluster
gauss_laplace_down_sampling
Module API
-----------
"""
from __future__ import division, absolute_import, print_function
import numpy as np
from .deblending import DeblendedProducts
from ....core.clustering.density_based.denclue import DENCLUE
from ....pipeline_manager.analysis_tasks import AnalysisTask
from ....core.image_processing.features import BlobDetectionGaussLaplace
from ....core.geometry.distance import dist_eval
__author__ = 'andrea tramacere'
[docs]def denclue_plot(denclue_debl,parent_cluster,target_coordinates=None,target_coordinates_w_array=None):
from asterism.plotting.plot_tools import show_image,plot_contour
import pylab as plt
if target_coordinates is not None:
x_t=target_coordinates[:,0]
y_t=target_coordinates[:,1]
image,off_set_x,off_set_y,masked_pixels=parent_cluster.get_cluster_Image(border=1,bkg=0)
im,ax=show_image(image.array)
if denclue_debl.attractors_position_array is not None:
x=denclue_debl.attractors_position_array[:,0]
y=denclue_debl.attractors_position_array[:,1]
dx=off_set_x
dy=off_set_y
ax.scatter(x-dx,y-dy,marker='o', edgecolor='black', linewidth='0',facecolor='white',s=100)
plot_contour(ax,parent_cluster.contour_shape[:,0]-dx,parent_cluster.contour_shape[:,1]-dy,'black')
if target_coordinates_w_array is not None:
ax.scatter(x_t-dx,y_t-dy,marker='+',color='black',s=10)
for cl in denclue_debl.clusters_list:
cl.add_contour_points()
if cl.contour_shape is not None:
plot_contour(ax,cl.contour_shape[:,0]-dx,cl.contour_shape[:,1]-dy,'w')
plt.plot(cl.x_c-dx,cl.y_c-dy,'+',ms=10)
print("|dbs cl=%d, x_c=%f, y_c=%f"%(cl.ID,cl.x_c-dx,cl.y_c-dy))
plt.show()
[docs]def plot_final_clusters(parent_cluster,deblended_clusters_list):
from asterism.plotting.plot_tools import plot_contour,show_image
import pylab as plt
image,off_set_x,off_set_y,masked_pixels=parent_cluster.get_cluster_Image(border=1,bkg=0)
dx=off_set_x
dy=off_set_y
im,ax=show_image(image.array)
for cl in deblended_clusters_list:
if cl.contour_shape is not None:
plot_contour(ax,cl.contour_shape[:,0]-dx,cl.contour_shape[:,1]-dy,'w')
plt.show()
[docs]def get_pixel_matching_gauss_laplace_lm(parent_cluster,h,local_maxima_xy,off_x,off_y):
lm_pix_coords=np.copy(local_maxima_xy)
for ID,lm in enumerate(local_maxima_xy):
#print ('x,y ',lm[0],lm[1],lm[0]+off_x,lm[1]+off_y)
d=dist_eval(parent_cluster.cartesian.coords,x_c=lm[0]+off_x,y_c=lm[1]+off_y,metric='euclidean')
msk1=d<h
if parent_cluster.flux[msk1].size>0:
lm_id=np.argmax(parent_cluster.flux[msk1])
#move back to cluster image coords, shift to cluster coords is done later
#print ('f1',parent_cluster.flux[msk1].max())
lm_pix_coords[ID]=parent_cluster.cartesian.coords[msk1][lm_id]-[off_x,off_y]
return lm_pix_coords
[docs]def gauss_laplace_down_sampling(parent_cluster,h,verbose=False,gl_th_rel=0.1,h_grid=3):
"""
Performs the down sampling of the target coordinates to use in the Denclue.
* the cluster image is extracted
* local maxima are extracted using the :class:`BlobDetectionGaussLaplace`
* if number of local maxima is leq 1 then
* target_coordinates are set to None
* target_coordinates_w_array are set to None
* if number of local maxima is geq 2 then
* target_coordinates are obtained from a grid with size = max(1,h*0.5)
* target_coordinates_w_array are obtained from the image flux at grid coordinates
* coordinates and w_array corresponding to local maxima are added to the target arrays
Parameters
----------
parent_cluster :
h : float
sets the side of the grid and the max_sigma=h, and min_sigma=h*0.5 used
in the :class:`BlobDetectionGaussLaplace`
Returns
-------
target_coordinates :
target_coordinates_w_array :
"""
bkg=parent_cluster.flux.min()
image,off_x,off_y,masked=parent_cluster.get_cluster_image_array(bkg=bkg)
h2=2.0*h
h1=max(1.0,h2*0.5)
#h_grid=max(1,np.int(h2*0.5))
if verbose == True:
print("|downsampling start")
print('|h=',h,)
print('|h1,h2=',h1, h2)
print('|h_grid=', h_grid)
blob_log_det=BlobDetectionGaussLaplace(min_sigma=h1,
max_sigma=h2,
num_sigma=2,
threshold_rel=gl_th_rel,
overlap=0.5)
local_maxima_xy,local_maxima_flux=blob_log_det.apply(image,coords_in_xy=True,verbose=verbose)
#print ("local maxima",local_maxima_xy,local_maxima_flux)
#print ("pixel matched local maxima",local_maxima_xy,local_maxima_flux)
# Local maxima added to target point final list
if local_maxima_flux is not None:
local_maxima_xy = get_pixel_matching_gauss_laplace_lm(parent_cluster, h, local_maxima_xy, off_x, off_y)
if verbose == True:
print('|local maxima', local_maxima_flux.size)
r=local_maxima_xy[:,1].astype(np.int16)
c=local_maxima_xy[:,0].astype(np.int16)
m=masked[r,c]==False
r=r[m]
c=c[m]
x=c
y=r
else:
x=None
y=None
#grid to extract target points with step given by kernel width
step=max(1.0,np.int(h_grid))
row_sel=np.arange(0,image.shape[0],step,dtype=np.int16)
col_sel=np.arange(0,image.shape[1],step,dtype=np.int16)
mgrid=np.zeros(image.shape, dtype=np.bool)
for row in row_sel:
mgrid[row][col_sel]=True
#target points are extracted from grid vertices
#and appended to targets from local maxima
y_grid,x_grid=np.where(np.logical_not(masked)*mgrid)
if local_maxima_flux is not None:
x_t=np.append(x,x_grid).astype(np.int16)
y_t=np.append(y,y_grid).astype(np.int16)
else:
x_t=x_grid
y_t=y_grid
#fluxes are of target points are extracted
target_coordinates_w_array=image[np.int16(y_t),np.int16(x_t)]
#cluster image offset is added to target coordinates
x_t=x_t+off_x
y_t=y_t+off_y
target_coordinates=np.column_stack((x_t.astype(np.int16),y_t.astype(np.int16)))
if verbose == True:
print("|downsampling stop")
return target_coordinates,target_coordinates_w_array,local_maxima_xy,local_maxima_flux
[docs]def do_denclue_of_cluster(parent_cluster,
min_size,
kernel,
eps,
h_frac,
h_min,
attr_dbs_K,
attr_dbs_eps,
mask_unchanged,
k_table_size,
R_max_frac,
R_max_kern_th,
gl_th_rel=0.1,
h_max=None,
plot=False,
digitize_attractors=False,
gl_downsampling=False,
h_grid=3,
verbose=False):
"""
This function sets-up the parameters for the :class:`.DENCLUE` partitioning of the parent clusters.
The following features are
Parameters
----------
parent_cluster
min_size
kernel
eps
h_frac
attr_dbs_K
attr_dbs_eps
mask_unchanged
k_table_size
R_max_frac
R_max_kern_th
plot
digitize_attractors
gl_downsampling
Returns
-------
"""
print ("|debelnding for cluster",parent_cluster.ID, "r_cluter",parent_cluster.r_cluster,"r_max",parent_cluster.r_max,'size',parent_cluster.cl_len,min_size)
if parent_cluster.cl_len>min_size:
h=h_frac*np.sqrt(parent_cluster.r_max**2+parent_cluster.r_cluster**2)
if h_min is not None:
h=max(h,h_min)
if h_max is not None:
h=min(h,h_max)
if R_max_frac is not None:
R_max_method='fixed'
R_max=R_max_frac*parent_cluster.r_cluster
elif R_max_kern_th is not None:
R_max_method='kern_rel_th'
R_max=None
else:
R_max_method='auto'
R_max=None
if verbose==True:
print ("|eps",eps)
print ("|h",h)
denclue_debl=DENCLUE(eps=eps,h=h,kernel=kernel,k_table_size=k_table_size,R_max_method=R_max_method,R_max_th_rel=R_max_kern_th,R_max=R_max)
if gl_downsampling == True:
target_coordinates,target_coordinates_w_array,local_maxima_xy,local_maxima_flux=gauss_laplace_down_sampling(parent_cluster,h,gl_th_rel,h_grid=h_grid)
else:
target_coordinates=None
target_coordinates_w_array=None
denclue_debl.run(parent_cluster.cartesian.coords,
dbscan_eps=attr_dbs_eps,
dbscan_K=attr_dbs_K,
weight_array=parent_cluster.pts_weight,
mask_unchanged=mask_unchanged,
R_max=R_max,
make_convolved_image=False,
digitize_attractors=digitize_attractors,
target_coordinates=target_coordinates,
target_coordinates_w_array=target_coordinates_w_array,
verbose=verbose)
print("|denclue done")
if plot is True:
denclue_plot(denclue_debl,parent_cluster,target_coordinates=target_coordinates,target_coordinates_w_array=target_coordinates_w_array)
children_clusters_list=denclue_debl.clusters_list
if len(children_clusters_list)==1:
children_clusters_list=[]
else:
children_clusters_list=[]
return children_clusters_list
[docs]def do_denclue_deblending(clusters_list,
eps=0.1,
h_frac=0.20,
h_min=1.0,
h_max=None,
kernel='gauss',
gl_downsampling=False,
h_grid=3,
gl_th_rel=0.1,
attr_dbs_eps=1.0,
attr_dbs_K=4.0,
min_size=9,
mask_unchanged=0.5,
k_table_size=None,
R_max_frac=1.0,
R_max_kern_th=None,
plot=False,
verbose=False,
digitize_attractors=False):
"""
This function implements the top-level algorithm for the Denclue-based deblending:
* each parent cluster in the `cluster_list` is partitioned by the :func:`do_denclue_of_cluster` function.
* the parent cluster with his children clusters are used to build :class:`DeblendedProducts` object
* a list of class:`DeblendedProducts` object is returned
Parameters
----------
clusters_list
eps
h_frac
kernel
gl_downsampling
attr_dbs_eps
attr_dbs_K
children_min_frac_integ_flux
children_min_frac_peak_flux
children_min_frac_size
children_bright_frac_peak_flux
min_size
mask_unchanged
k_table_size
R_max_frac
R_max_kern_th
plot
denclue_catalog
verbose
do_denclue_deblending
out_catalog
digitize_attractors
Returns
-------
"""
deblended_prod_list=[]
for parent_cluster in clusters_list:
print("|-")
denclue_clusters=do_denclue_of_cluster(parent_cluster,
min_size,
kernel,
eps,
h_frac,
h_min,
attr_dbs_K,
attr_dbs_eps,
mask_unchanged,
k_table_size,
gl_th_rel=gl_th_rel,
h_max=h_max,
R_max_frac=R_max_frac,
R_max_kern_th=R_max_kern_th,
plot=plot,
digitize_attractors=digitize_attractors,
gl_downsampling=gl_downsampling,
h_grid=h_grid,
verbose=verbose)
print("|-")
#if plot==True:
# plot_final_clusters(parent_cluster,denclue_clusters)
deblended_prod_list.append(DeblendedProducts(parent=parent_cluster,
children_list=denclue_clusters,))
return deblended_prod_list
[docs]class DoDENCLUEDeblendingTask(AnalysisTask):
"""
"""
def __init__(self, name='denclue_deblending', func=do_denclue_deblending, parser=None, process=None):
super(DoDENCLUEDeblendingTask,self).__init__(name=name,func=func,parser=parser,process=process)
self.parameters_list.add_par('-gl_downsampling',help='if True, a GuassLaplace downsampling of the target coordiantes is performed ',action='store_true')
self.parameters_list.add_par('-h_grid', type=np.int, help=' grid mesh in pixels for the downsampling', default=3)
self.parameters_list.add_par('-gl_th_rel', type=np.float, help='relative threshold for the GaussLaplace local maxima detection', default=0.1)
self.parameters_list.add_par('-eps', type=np.float ,help='sets the stop threshold for the recursive attractor update rule',default=0.01)
self.parameters_list.add_par('-digitize_attractors', help='if True, will approximate the final attractor position to the coordinate of the closest image pixel ',action='store_true')
self.parameters_list.add_par('-kernel',type=str, help='smoothing kernel',default='gauss')
self.parameters_list.add_par('-k_table_size',type=np.int, help='if not None, the kernel is precomputed in a look-up table over a grid with size equal to k_table_size', default=None)
self.parameters_list.add_par('-h_frac', type=np.float ,help='sets the width of the kernel as h=h_frac*sqrt(r_max^2+r_cluster^2)',default=0.1 )
self.parameters_list.add_par('-h_min', type=np.float ,help='kernel width obtained by h_frac can not be lower than h_min',default=1.0 )
self.parameters_list.add_par('-h_max', type=np.float ,help='kernel width obtained by h_frac can not be lerger than h_max',default=None )
self.parameters_list.add_par('-mask_unchanged', type=np.float ,help='th value to flag and attractor as unchanged',default=None )
self.parameters_list.add_par('-min_size',type=np.int, help=' sets the minimum size in pixels to perform a deblending', default=9)
self.parameters_list.add_par('-attr_dbs_eps', type=np.float ,help='eps for dbscan of attractors',default=1.0 )
self.parameters_list.add_par('-attr_dbs_K', type=np.float ,help='K for dbscan of attractors',default=4.0 )
self.parameters_list.add_par('-verbose',help='set verbosity on',action='store_true')
self.parameters_list.add_par('-plot',help='plot ',action='store_true')
group_method=self.parameters_list.add_mex_group('R_max_method')
group_method.add_par('-R_max_frac', type=np.float ,help='max kernel influence radius as fraction of parent cluster r_max',default=None )
group_method.add_par('-R_max_kern_th', type=np.float, help='max kernel influence radius when kernel<=th with th relative to kenerl max value',default=None)
self.func=func
