18 from __future__
import division, print_function
23 import _SourceXtractorPy
as cpp
24 from .measurement_images
import MeasurementGroup
26 from astropy
import units
as u
27 from astropy.coordinates
import SkyCoord
28 from astropy.coordinates
import Angle
40 Limit, and normalize, the range of values for a model fitting parameter. 45 limits : a tuple (min, max), or a callable that receives a source, and returns a tuple (min, max) 51 Normalized to engine space using a sigmoid function 55 engine = \ln \frac{world - min}{max-world} \\ 56 world = min + \frac{max - min}{1 + e^{engine}} 59 Normalized to engine space using an exponential sigmoid function 63 engine = \ln \left( \frac{\ln(world/min)}{\ln(max /world)} \right) \\ 64 world = min * e^\frac{ \ln(max / min) }{ (1 + e^{-engine}) } 80 Receives a source, and returns a tuple (min, max) 97 Human readable representation for the object 100 if hasattr(self.
__limits,
'__call__'):
105 RangeType.LINEAR :
'LIN',
106 RangeType.EXPONENTIAL :
'EXP' 108 res +=
',{}]'.format(type_str[self.
__type])
114 Unbounded, but normalize, value of a model fitting parameter 118 normalization_factor: float, or a callable that receives the initial value parameter value and a source, 120 The world value which will be normalized to 1 in engine coordinates 134 Receives the initial parameter value and a source, and returns the world value which will be 135 normalized to 1 in engine coordinates 144 Human readable representation for the object 155 constant_parameter_dict = {}
156 free_parameter_dict = {}
157 dependent_parameter_dict = {}
162 Print a human-readable representation of the configured model fitting parameters. 167 Where to print the representation. Defaults to sys.stderr 169 if constant_parameter_dict:
170 print(
'Constant parameters:', file=file)
171 for n
in constant_parameter_dict:
172 print(
' {}: {}'.format(n, constant_parameter_dict[n]), file=file)
173 if free_parameter_dict:
174 print(
'Free parameters:', file=file)
175 for n
in free_parameter_dict:
176 print(
' {}: {}'.format(n, free_parameter_dict[n]), file=file)
177 if dependent_parameter_dict:
178 print(
'Dependent parameters:', file=file)
179 for n
in dependent_parameter_dict:
180 print(
' {}: {}'.format(n, dependent_parameter_dict[n]), file=file)
185 Base class for all model fitting parameter types. 186 Can not be used directly. 194 Human readable representation for the object 196 return '(ID:{})'.format(self.id)
201 A parameter with a single value that remains constant. It will not be fitted. 205 value : float, or callable that receives a source and returns a float 206 Value for this parameter 213 ParameterBase.__init__(self)
215 constant_parameter_dict[self.id] = self
222 Receives a source and returns a value for the parameter 231 Human readable representation for the object 233 res = ParameterBase.__str__(self)[:-1] +
', value:' 234 if hasattr(self.
__value,
'__call__'):
243 A parameter that will be fitted by the model fitting engine. 247 init_value : float or callable that receives a source, and returns a float 248 Initial value for the parameter. 249 range : instance of Range or Unbounded 250 Defines if this parameter is unbounded or bounded, and how. 259 >>> sersic = FreeParameter(2.0, Range((1.0, 7.0), RangeType.LINEAR)) 266 ParameterBase.__init__(self)
269 free_parameter_dict[self.id] = self
276 Receives a source, and returns an initial value for the parameter. 293 Human readable representation for the object 295 res = ParameterBase.__str__(self)[:-1] +
', init:' 301 res +=
', range:' + str(self.
__range)
307 A DependentParameter is not fitted by itself, but its value is derived from another Parameters, whatever their type: 308 FreeParameter, ConstantParameter, or other DependentParameter 313 A callable that will be called with all the parameters specified in this constructor each time a new 314 evaluation is needed. 315 params : list of ParameterBase 316 List of parameters on which this DependentParameter depends. 320 >>> flux = get_flux_parameter() 321 >>> mag = DependentParameter(lambda f: -2.5 * np.log10(f) + args.mag_zeropoint, flux) 322 >>> add_output_column('mf_mag_' + band, mag) 329 ParameterBase.__init__(self)
332 dependent_parameter_dict[self.id] = self
337 Convenience function for the position parameter X and Y. 342 X coordinate, starting at the X coordinate of the centroid and linearly limited to X +/- the object radius. 344 Y coordinate, starting at the Y coordinate of the centroid and linearly limited to Y +/- the object radius. 347 X and Y are fitted on the detection image X and Y coordinates. Internally, these are translated to measurement 348 images using the WCS headers. 351 FreeParameter(
lambda o: o.get_centroid_x(),
Range(
lambda v,o: (v-o.get_radius(), v+o.get_radius()), RangeType.LINEAR)),
352 FreeParameter(
lambda o: o.get_centroid_y(),
Range(
lambda v,o: (v-o.get_radius(), v+o.get_radius()), RangeType.LINEAR))
358 Possible flux types to use as initial value for the flux parameter. 359 Right now, only isophotal is supported. 366 Convenience function for the flux parameter. 371 One of the values defined in FluxParameterType 373 Scaling of the initial flux. Defaults to 1. 378 Flux parameter, starting at the flux defined by `type`, and limited to +/- 1e3 times the initial value. 381 FluxParameterType.ISO :
'get_iso_flux' 383 return FreeParameter(
lambda o: getattr(o, func_map[type])() * scale,
Range(
lambda v,o: (v * 1E-3, v * 1E3), RangeType.EXPONENTIAL))
391 Model a Gaussian prior on a given parameter. 395 param : ParameterBase 396 Model fitting parameter 397 value : float or callable that receives a source and returns a float 399 sigma : float or callable that receives a source and returns a float 400 Standard deviation of the Gaussian 407 cpp.Id.__init__(self)
409 self.
value = value
if hasattr(value,
'__call__')
else lambda o: value
410 self.
sigma = sigma
if hasattr(sigma,
'__call__')
else lambda o: sigma
415 Add a prior to the given parameter. 419 param : ParameterBase 420 value : float or callable that receives a source and returns a float 422 sigma : float or callable that receives a source and returns a float 423 Standard deviation of the Gaussian 425 prior =
Prior(param, value, sigma)
426 prior_dict[prior.id] = prior
429 frame_models_dict = {}
434 if isinstance(x, tuple):
437 if not x.id
in frame_models_dict:
438 frame_models_dict[x.id] = []
439 frame_models_dict[x.id].append(model.id)
444 Add a model to be fitted to the given group. 448 group : MeasurementGroup 451 if not isinstance(group, MeasurementGroup):
452 raise TypeError(
'Models can only be added on MeasurementGroup, got {}'.format(type(group)))
453 if not hasattr(group,
'models'):
455 group.models.append(model)
461 Print a human-readable representation of the configured models. 465 group : MeasurementGroup 466 Print the models for this group. 468 If True, print also the parameters that belong to the model 470 Prefix each line with this string. Used internally for indentation. 472 Where to print the representation. Defaults to sys.stderr 474 if hasattr(group,
'models')
and group.models:
475 print(
'{}Models:'.format(prefix), file=file)
476 for m
in group.models:
477 print(
'{} {}'.format(prefix, m.to_string(show_params)), file=file)
479 if isinstance(x, tuple):
480 print(
'{}{}:'.format(prefix, x[0]), file=file)
484 constant_model_dict = {}
485 point_source_model_dict = {}
486 sersic_model_dict = {}
487 exponential_model_dict = {}
488 de_vaucouleurs_model_dict = {}
489 params_dict = {
"max_iterations" : 100,
"modified_chi_squared_scale" : 10,
"engine" :
"levmar" }
497 Max number of iterations for the model fitting. 499 params_dict[
"max_iterations"] = iterations
507 Sets u0, as used by the modified chi squared residual comparator, a function that reduces the effect of large 509 Refer to the SourceXtractor++ documentation for a better explanation of how residuals are computed and how 510 this value affects the model fitting. 512 params_dict[
"modified_chi_squared_scale"] = scale
520 Minimization engine for the model fitting : levmar or gsl 522 params_dict[
"engine"] = engine
527 Base class for all models. 532 class CoordinateModelBase(ModelBase):
534 Base class for positioned models with a flux. It can not be used directly. 538 x_coord : ParameterBase or float 539 X coordinate (in the detection image) 540 y_coord : ParameterBase or float 541 Y coordinate (in the detection image) 542 flux : ParameterBase or float 550 ModelBase.__init__(self)
558 Models a source as a point, spread by the PSF. 562 x_coord : ParameterBase or float 563 X coordinate (in the detection image) 564 y_coord : ParameterBase or float 565 Y coordinate (in the detection image) 566 flux : ParameterBase or float 574 CoordinateModelBase.__init__(self, x_coord, y_coord, flux)
575 global point_source_model_dict
576 point_source_model_dict[self.id] = self
580 Return a human readable representation of the model. 585 If True, include information about the parameters. 592 return 'PointSource[x_coord={}, y_coord={}, flux={}]'.format(
595 return 'PointSource[x_coord={}, y_coord={}, flux={}]'.format(
601 A model that is constant through all the image. 605 value: ParameterBase or float 606 Value to add to the value of all pixels from the model. 613 ModelBase.__init__(self)
615 global constant_model_dict
616 constant_model_dict[self.id] = self
620 Return a human readable representation of the model. 625 If True, include information about the parameters. 632 return 'ConstantModel[value={}]'.format(self.
value)
634 return 'ConstantModel[value={}]'.format(self.
value.id)
639 Base class for the Sersic, Exponential and de Vaucouleurs models. It can not be used directly. 643 x_coord : ParameterBase or float 644 X coordinate (in the detection image) 645 y_coord : ParameterBase or float 646 Y coordinate (in the detection image) 647 flux : ParameterBase or float 649 effective_radius : ParameterBase or float 650 Ellipse semi-major axis, in pixels on the detection image. 651 aspect_ratio : ParameterBase or float 653 angle : ParameterBase or float 654 Ellipse rotation, in radians 657 def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle):
661 CoordinateModelBase.__init__(self, x_coord, y_coord, flux)
669 Model a source with a Sersic profile. 673 x_coord : ParameterBase or float 674 X coordinate (in the detection image) 675 y_coord : ParameterBase or float 676 Y coordinate (in the detection image) 677 flux : ParameterBase or float 679 effective_radius : ParameterBase or float 680 Ellipse semi-major axis, in pixels on the detection image. 681 aspect_ratio : ParameterBase or float 683 angle : ParameterBase or float 684 Ellipse rotation, in radians 685 n : ParameterBase or float 689 def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle, n):
693 SersicModelBase.__init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle)
695 global sersic_model_dict
696 sersic_model_dict[self.id] = self
700 Return a human readable representation of the model. 705 If True, include information about the parameters. 712 return 'Sersic[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}, n={}]'.format(
715 return 'Sersic[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}, n={}]'.format(
721 Model a source with an exponential profile (Sersic model with an index of 1) 725 x_coord : ParameterBase or float 726 X coordinate (in the detection image) 727 y_coord : ParameterBase or float 728 Y coordinate (in the detection image) 729 flux : ParameterBase or float 731 effective_radius : ParameterBase or float 732 Ellipse semi-major axis, in pixels on the detection image. 733 aspect_ratio : ParameterBase or float 735 angle : ParameterBase or float 736 Ellipse rotation, in radians 739 def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle):
743 SersicModelBase.__init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle)
744 global exponential_model_dict
745 exponential_model_dict[self.id] = self
749 Return a human readable representation of the model. 754 If True, include information about the parameters. 761 return 'Exponential[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
764 return 'Exponential[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
770 Model a source with a De Vaucouleurs profile (Sersic model with an index of 4) 774 x_coord : ParameterBase or float 775 X coordinate (in the detection image) 776 y_coord : ParameterBase or float 777 Y coordinate (in the detection image) 778 flux : ParameterBase or float 780 effective_radius : ParameterBase or float 781 Ellipse semi-major axis, in pixels on the detection image. 782 aspect_ratio : ParameterBase or float 784 angle : ParameterBase or float 785 Ellipse rotation, in radians 788 def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle):
792 SersicModelBase.__init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle)
793 global de_vaucouleurs_model_dict
794 de_vaucouleurs_model_dict[self.id] = self
798 Return a human readable representation of the model. 803 If True, include information about the parameters. 810 return 'DeVaucouleurs[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
813 return 'DeVaucouleurs[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
819 Coordinates in right ascension and declination 839 Transform an (X, Y) in pixel coordinates on the detection image to (RA, DEC) in world coordinates. 849 global coordinate_system
850 wc = coordinate_system.image_to_world(cpp.ImageCoordinate(x-1, y-1))
856 Transform an (X, Y) in pixel coordinates on the detection image to astropy SkyCoord. 868 sky_coord = SkyCoord(ra=coord.ra*u.degree, dec=coord.dec*u.degree)
874 Transform a radius in pixels on the detection image to a radius in sky coordinates. 891 Get the separation angle in sky coordinates for two points defined in pixels on the detection image. 902 Separation in degrees 906 return c1.separation(c2).degree
911 Get the position angle in sky coordinates for two points defined in pixels on the detection image. 922 Position angle in degrees, normalized to -/+ 90 926 angle = c1.position_angle(c2).degree
929 return (angle + 90.0) % 180.0 - 90.0
934 Set the global coordinate system. This function is used internally by SourceXtractor++. 936 global coordinate_system
937 coordinate_system = cs
942 Convenience function for generating two dependent parameter with world (alpha, delta) coordinates 943 from image (X, Y) coordinates. 952 ra : DependentParameter 953 dec : DependentParameter 961 >>> x, y = get_pos_parameters() 962 >>> ra, dec = get_world_position_parameters(x, y) 963 >>> add_output_column('mf_ra', ra) 964 >>> add_output_column('mf_dec', dec) 973 Convenience function for generating five dependent parameters, in world coordinates, for the position 974 and shape of a model. 980 radius : ParameterBase 981 angle : ParameterBase 982 ratio : ParameterBase 986 ra : DependentParameter 988 dec : DependentParameter 990 rad : DependentParameter 992 angle : DependentParameter 994 ratio : DependentParameter 995 Aspect ratio. It has to be recomputed as the axis of the ellipse may have different ratios 996 in image coordinates than in world coordinates 1000 >>> flux = get_flux_parameter() 1001 >>> x, y = get_pos_parameters() 1002 >>> radius = FreeParameter(lambda o: o.get_radius(), Range(lambda v, o: (.01 * v, 100 * v), RangeType.EXPONENTIAL)) 1003 >>> angle = FreeParameter(lambda o: o.get_angle(), Range((-np.pi, np.pi), RangeType.LINEAR)) 1004 >>> ratio = FreeParameter(1, Range((0, 10), RangeType.LINEAR)) 1005 >>> add_model(group, ExponentialModel(x, y, flux, radius, ratio, angle)) 1006 >>> ra, dec, wc_rad, wc_angle, wc_ratio = get_world_parameters(x, y, radius, angle, ratio) 1007 >>> add_output_column('mf_world_angle', wc_angle) 1012 def get_major_axis(x, y, radius, angle, ratio):
1014 x2 = x + math.cos(angle) * radius
1015 y2 = y + math.sin(angle) * radius
1017 x2 = x + math.sin(angle) * radius * ratio
1018 y2 = y + math.cos(angle) * radius * ratio
1022 def get_minor_axis(x, y, radius, angle, ratio):
1024 x2 = x + math.sin(angle) * radius * ratio
1025 y2 = y + math.cos(angle) * radius * ratio
1027 x2 = x + math.cos(angle) * radius
1028 y2 = y + math.sin(angle) * radius
1032 def wc_rad_func(x, y, radius, angle, ratio):
1033 x2, y2 = get_major_axis(x, y, radius, angle, ratio)
1036 def wc_angle_func(x, y, radius, angle, ratio):
1037 x2, y2 = get_major_axis(x, y, radius, angle, ratio)
1040 def wc_ratio_func(x, y, radius, angle, ratio):
1041 minor_x, minor_y = get_minor_axis(x, y, radius, angle, ratio)
1044 major_x, major_y = get_major_axis(x, y, radius, angle, ratio)
1047 return minor_angle / major_angle
1053 return (ra, dec, wc_rad, wc_angle, wc_ratio)