18 from __future__
import division, print_function
25 import _SourceXtractorPy
as cpp
29 warnings.warn(
'Could not import pyston: running outside sourcextractor?', ImportWarning)
30 from astropy
import units
as u
31 from astropy.coordinates
import SkyCoord
33 from .measurement_images
import MeasurementGroup
43 Limit, and normalize, the range of values for a model fitting parameter.
48 limits : a tuple (min, max), or a callable that receives a source, and returns a tuple (min, max)
54 Normalized to engine space using a sigmoid function
58 engine = \ln \frac{world - min}{max-world} \\
59 world = min + \frac{max - min}{1 + e^{engine}}
62 Normalized to engine space using an exponential sigmoid function
66 engine = \ln \left( \frac{\ln(world/min)}{\ln(max /world)} \right) \\
67 world = min * e^\frac{ \ln(max / min) }{ (1 + e^{-engine}) }
76 self.
__callable__callable = limits
if hasattr(limits,
'__call__')
else lambda v, o: limits
84 o : object being fitted
88 The minimum acceptable value for the range
97 o : object being fitted
101 The maximum acceptable value for the range
118 Human readable representation for the object
121 if hasattr(self.
__limits__limits,
'__call__'):
126 RangeType.LINEAR :
'LIN',
127 RangeType.EXPONENTIAL :
'EXP'
129 res +=
',{}]'.format(type_str[self.
__type__type])
135 Unbounded, but normalize, value of a model fitting parameter
139 normalization_factor: float, or a callable that receives the initial value parameter value and a source,
141 The world value which will be normalized to 1 in engine coordinates
149 if hasattr(normalization_factor,
'__call__'):
152 self.
__callable__callable =
lambda v, o: normalization_factor
159 o : object being fitted
164 The world value which will be normalized to 1 in engine coordinates
173 Human readable representation for the object
184 constant_parameter_dict = {}
185 free_parameter_dict = {}
186 dependent_parameter_dict = {}
191 Print a human-readable representation of the configured model fitting parameters.
196 Where to print the representation. Defaults to sys.stderr
198 if constant_parameter_dict:
199 print(
'Constant parameters:', file=file)
200 for n
in constant_parameter_dict:
201 print(
' {}: {}'.format(n, constant_parameter_dict[n]), file=file)
202 if free_parameter_dict:
203 print(
'Free parameters:', file=file)
204 for n
in free_parameter_dict:
205 print(
' {}: {}'.format(n, free_parameter_dict[n]), file=file)
206 if dependent_parameter_dict:
207 print(
'Dependent parameters:', file=file)
208 for n
in dependent_parameter_dict:
209 print(
' {}: {}'.format(n, dependent_parameter_dict[n]), file=file)
214 Base class for all model fitting parameter types.
215 Can not be used directly.
223 Human readable representation for the object
225 return '(ID:{})'.format(self.id)
230 A parameter with a single value that remains constant. It will not be fitted.
234 value : float, or callable that receives a source and returns a float
235 Value for this parameter
242 ParameterBase.__init__(self)
244 self.
__callable__callable = value
if hasattr(value,
'__call__')
else lambda o: value
245 constant_parameter_dict[self.id] = self
251 o : object being fitted
256 Value of the constant parameter
265 Human readable representation for the object
267 res = ParameterBase.__str__(self)[:-1] +
', value:'
268 if hasattr(self.
__value__value,
'__call__'):
271 res += str(self.
__value__value)
277 A parameter that will be fitted by the model fitting engine.
281 init_value : float or callable that receives a source, and returns a float
282 Initial value for the parameter.
283 range : instance of Range or Unbounded
284 Defines if this parameter is unbounded or bounded, and how.
293 >>> sersic = FreeParameter(2.0, Range((1.0, 7.0), RangeType.LINEAR))
300 ParameterBase.__init__(self)
302 self.
__init_call__init_call = init_value
if hasattr(init_value,
'__call__')
else lambda o: init_value
304 free_parameter_dict[self.id] = self
310 o : object being fitted
315 Initial value for the free parameter
332 Human readable representation for the object
334 res = ParameterBase.__str__(self)[:-1] +
', init:'
340 res +=
', range:' + str(self.
__range__range)
346 A DependentParameter is not fitted by itself, but its value is derived from another Parameters, whatever their type:
347 FreeParameter, ConstantParameter, or other DependentParameter
352 A callable that will be called with all the parameters specified in this constructor each time a new
353 evaluation is needed.
354 params : list of ParameterBase
355 List of parameters on which this DependentParameter depends.
359 >>> flux = get_flux_parameter()
360 >>> mag = DependentParameter(lambda f: -2.5 * np.log10(f) + args.mag_zeropoint, flux)
361 >>> add_output_column('mf_mag_' + band, mag)
368 ParameterBase.__init__(self)
370 self.
paramsparams = [p.id
for p
in params]
371 dependent_parameter_dict[self.id] = self
376 Convenience function for the position parameter X and Y.
381 X coordinate, starting at the X coordinate of the centroid and linearly limited to X +/- the object radius.
383 Y coordinate, starting at the Y coordinate of the centroid and linearly limited to Y +/- the object radius.
386 X and Y are fitted on the detection image X and Y coordinates. Internally, these are translated to measurement
387 images using the WCS headers.
389 param_range =
Range(
lambda v,o: (v-o.radius, v+o.radius), RangeType.LINEAR)
398 Possible flux types to use as initial value for the flux parameter.
399 Right now, only isophotal is supported.
406 Convenience function for the flux parameter.
411 One of the values defined in FluxParameterType
413 Scaling of the initial flux. Defaults to 1.
418 Flux parameter, starting at the flux defined by `type`, and limited to +/- 1e3 times the initial value.
421 FluxParameterType.ISO :
'isophotal_flux'
423 return FreeParameter(
lambda o: getattr(o, attr_map[type]) * scale,
Range(
lambda v,o: (v * 1E-3, v * 1E3), RangeType.EXPONENTIAL))
431 Model a Gaussian prior on a given parameter.
435 param : ParameterBase
436 Model fitting parameter
437 value : float or callable that receives a source and returns a float
439 sigma : float or callable that receives a source and returns a float
440 Standard deviation of the Gaussian
447 cpp.Id.__init__(self)
449 self.
valuevalue = value
if hasattr(value,
'__call__')
else lambda o: value
450 self.
sigmasigma = sigma
if hasattr(sigma,
'__call__')
else lambda o: sigma
455 Add a prior to the given parameter.
459 param : ParameterBase
460 value : float or callable that receives a source and returns a float
462 sigma : float or callable that receives a source and returns a float
463 Standard deviation of the Gaussian
465 prior =
Prior(param, value, sigma)
466 prior_dict[prior.id] = prior
469 frame_models_dict = {}
474 if isinstance(x, tuple):
477 if not x.id
in frame_models_dict:
478 frame_models_dict[x.id] = []
479 frame_models_dict[x.id].append(model.id)
484 Add a model to be fitted to the given group.
488 group : MeasurementGroup
491 if not isinstance(group, MeasurementGroup):
492 raise TypeError(
'Models can only be added on MeasurementGroup, got {}'.format(type(group)))
493 if not hasattr(group,
'models'):
495 group.models.append(model)
501 Print a human-readable representation of the configured models.
505 group : MeasurementGroup
506 Print the models for this group.
508 If True, print also the parameters that belong to the model
510 Prefix each line with this string. Used internally for indentation.
512 Where to print the representation. Defaults to sys.stderr
514 if hasattr(group,
'models')
and group.models:
515 print(
'{}Models:'.format(prefix), file=file)
516 for m
in group.models:
517 print(
'{} {}'.format(prefix, m.to_string(show_params)), file=file)
519 if isinstance(x, tuple):
520 print(
'{}{}:'.format(prefix, x[0]), file=file)
524 constant_model_dict = {}
525 point_source_model_dict = {}
526 sersic_model_dict = {}
527 exponential_model_dict = {}
528 de_vaucouleurs_model_dict = {}
530 params_dict = {
"max_iterations": 200,
"modified_chi_squared_scale": 10,
"engine":
"",
"use_iterative_fitting":
True,
"meta_iterations": 5,
531 "deblend_factor": 0.95,
"meta_iteration_stop": 0.0001}
539 Max number of iterations for the model fitting.
541 params_dict[
"max_iterations"] = iterations
549 Sets u0, as used by the modified chi squared residual comparator, a function that reduces the effect of large
551 Refer to the SourceXtractor++ documentation for a better explanation of how residuals are computed and how
552 this value affects the model fitting.
554 params_dict[
"modified_chi_squared_scale"] = scale
562 Minimization engine for the model fitting : levmar or gsl
564 params_dict[
"engine"] = engine
570 use_iterative_fitting : boolean
571 use iterative model fitting or legacy
573 params_dict[
"use_iterative_fitting"] = use_iterative_fitting
579 meta_iterations : int
580 number of meta iterations on the whole group (when using iterative model fitting)
582 params_dict[
"meta_iterations"] = meta_iterations
590 params_dict[
"deblend_factor"] = deblend_factor
598 params_dict[
"meta_iteration_stop"] = meta_iteration_stop
605 Base class for all models.
610 class CoordinateModelBase(ModelBase):
612 Base class for positioned models with a flux. It can not be used directly.
616 x_coord : ParameterBase or float
617 X coordinate (in the detection image)
618 y_coord : ParameterBase or float
619 Y coordinate (in the detection image)
620 flux : ParameterBase or float
628 ModelBase.__init__(self)
636 Models a source as a point, spread by the PSF.
640 x_coord : ParameterBase or float
641 X coordinate (in the detection image)
642 y_coord : ParameterBase or float
643 Y coordinate (in the detection image)
644 flux : ParameterBase or float
652 CoordinateModelBase.__init__(self, x_coord, y_coord, flux)
653 global point_source_model_dict
654 point_source_model_dict[self.id] = self
658 Return a human readable representation of the model.
663 If True, include information about the parameters.
670 return 'PointSource[x_coord={}, y_coord={}, flux={}]'.format(
673 return 'PointSource[x_coord={}, y_coord={}, flux={}]'.format(
679 A model that is constant through all the image.
683 value: ParameterBase or float
684 Value to add to the value of all pixels from the model.
691 ModelBase.__init__(self)
693 global constant_model_dict
694 constant_model_dict[self.id] = self
698 Return a human readable representation of the model.
703 If True, include information about the parameters.
710 return 'ConstantModel[value={}]'.format(self.
valuevalue)
712 return 'ConstantModel[value={}]'.format(self.
valuevalue.id)
717 Base class for the Sersic, Exponential and de Vaucouleurs models. It can not be used directly.
721 x_coord : ParameterBase or float
722 X coordinate (in the detection image)
723 y_coord : ParameterBase or float
724 Y coordinate (in the detection image)
725 flux : ParameterBase or float
727 effective_radius : ParameterBase or float
728 Ellipse semi-major axis, in pixels on the detection image.
729 aspect_ratio : ParameterBase or float
731 angle : ParameterBase or float
732 Ellipse rotation, in radians
735 def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle):
739 CoordinateModelBase.__init__(self, x_coord, y_coord, flux)
747 Model a source with a Sersic profile.
751 x_coord : ParameterBase or float
752 X coordinate (in the detection image)
753 y_coord : ParameterBase or float
754 Y coordinate (in the detection image)
755 flux : ParameterBase or float
757 effective_radius : ParameterBase or float
758 Ellipse semi-major axis, in pixels on the detection image.
759 aspect_ratio : ParameterBase or float
761 angle : ParameterBase or float
762 Ellipse rotation, in radians
763 n : ParameterBase or float
767 def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle, n):
771 SersicModelBase.__init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle)
773 global sersic_model_dict
774 sersic_model_dict[self.id] = self
778 Return a human readable representation of the model.
783 If True, include information about the parameters.
790 return 'Sersic[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}, n={}]'.format(
793 return 'Sersic[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}, n={}]'.format(
799 Model a source with an exponential profile (Sersic model with an index of 1)
803 x_coord : ParameterBase or float
804 X coordinate (in the detection image)
805 y_coord : ParameterBase or float
806 Y coordinate (in the detection image)
807 flux : ParameterBase or float
809 effective_radius : ParameterBase or float
810 Ellipse semi-major axis, in pixels on the detection image.
811 aspect_ratio : ParameterBase or float
813 angle : ParameterBase or float
814 Ellipse rotation, in radians
817 def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle):
821 SersicModelBase.__init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle)
822 global exponential_model_dict
823 exponential_model_dict[self.id] = self
827 Return a human readable representation of the model.
832 If True, include information about the parameters.
839 return 'Exponential[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
842 return 'Exponential[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
848 Model a source with a De Vaucouleurs profile (Sersic model with an index of 4)
852 x_coord : ParameterBase or float
853 X coordinate (in the detection image)
854 y_coord : ParameterBase or float
855 Y coordinate (in the detection image)
856 flux : ParameterBase or float
858 effective_radius : ParameterBase or float
859 Ellipse semi-major axis, in pixels on the detection image.
860 aspect_ratio : ParameterBase or float
862 angle : ParameterBase or float
863 Ellipse rotation, in radians
866 def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle):
870 SersicModelBase.__init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle)
871 global de_vaucouleurs_model_dict
872 de_vaucouleurs_model_dict[self.id] = self
876 Return a human readable representation of the model.
881 If True, include information about the parameters.
888 return 'DeVaucouleurs[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
891 return 'DeVaucouleurs[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
896 ComputeGraphModel model
900 models: string or list of strings, corresponding to path to Onnx format models,
901 (specifying more than one allows using the most efficient model based on render size.)
902 x_coord : ParameterBase or float
903 X coordinate (in the detection image)
904 y_coord : ParameterBase or float
905 Y coordinate (in the detection image)
906 flux : ParameterBase or float
908 params : Dictionary of String and ParameterBase or float
909 Dictionary of custom parameters for the ONNX model
912 def __init__(self, models, x_coord, y_coord, flux, params={}):
916 CoordinateModelBase.__init__(self, x_coord, y_coord, flux)
918 ratio_name =
"_aspect_ratio"
919 angle_name =
"_angle"
920 scale_name =
"_scale"
922 for k
in params.keys():
923 if not isinstance(params[k], ParameterBase):
926 aspect_ratio = params[ratio_name]
if ratio_name
in params.keys()
else 1.0
927 angle = params[angle_name]
if angle_name
in params.keys()
else 0.0
928 scale = params[scale_name]
if scale_name
in params.keys()
else 1.0
934 params.pop(ratio_name,
None)
935 params.pop(angle_name,
None)
936 params.pop(scale_name,
None)
939 self.
modelsmodels = models
if isinstance(models, list)
else [models]
941 global onnx_model_dict
942 onnx_model_dict[self.id] = self
946 Return a human readable representation of the model.
951 If True, include information about the parameters.
958 return 'ComputeGraph[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
961 return 'ComputeGraph[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
967 Coordinates in right ascension and declination
987 Transform an (X, Y) in pixel coordinates on the detection image to (RA, DEC) in world coordinates.
998 wc_alpha = pyston.image_to_world_alpha(x - 1, y - 1)
999 wc_delta = pyston.image_to_world_delta(x - 1, y - 1)
1005 Transform an (X, Y) in pixel coordinates on the detection image to astropy SkyCoord.
1017 sky_coord = SkyCoord(ra=coord.ra*u.degree, dec=coord.dec*u.degree)
1023 Transform a radius in pixels on the detection image to a radius in sky coordinates.
1040 Get the separation angle in sky coordinates for two points defined in pixels on the detection image.
1051 Separation in degrees
1055 return c1.separation(c2).degree
1060 Get the position angle in sky coordinates for two points defined in pixels on the detection image.
1071 Position angle in degrees, normalized to -/+ 90
1075 angle = c1.position_angle(c2).degree
1078 return (angle + 90.0) % 180.0 - 90.0
1083 Convenience function for generating two dependent parameter with world (alpha, delta) coordinates
1084 from image (X, Y) coordinates.
1093 ra : DependentParameter
1094 dec : DependentParameter
1102 >>> x, y = get_pos_parameters()
1103 >>> ra, dec = get_world_position_parameters(x, y)
1104 >>> add_output_column('mf_ra', ra)
1105 >>> add_output_column('mf_dec', dec)
1114 Convenience function for generating five dependent parameters, in world coordinates, for the position
1115 and shape of a model.
1121 radius : ParameterBase
1122 angle : ParameterBase
1123 ratio : ParameterBase
1127 ra : DependentParameter
1129 dec : DependentParameter
1131 rad : DependentParameter
1133 angle : DependentParameter
1135 ratio : DependentParameter
1136 Aspect ratio. It has to be recomputed as the axis of the ellipse may have different ratios
1137 in image coordinates than in world coordinates
1141 >>> flux = get_flux_parameter()
1142 >>> x, y = get_pos_parameters()
1143 >>> radius = FreeParameter(lambda o: o.radius, Range(lambda v, o: (.01 * v, 100 * v), RangeType.EXPONENTIAL))
1144 >>> angle = FreeParameter(lambda o: o.angle, Range((-np.pi, np.pi), RangeType.LINEAR))
1145 >>> ratio = FreeParameter(1, Range((0, 10), RangeType.LINEAR))
1146 >>> add_model(group, ExponentialModel(x, y, flux, radius, ratio, angle))
1147 >>> ra, dec, wc_rad, wc_angle, wc_ratio = get_world_parameters(x, y, radius, angle, ratio)
1148 >>> add_output_column('mf_world_angle', wc_angle)
1153 def get_major_axis(x, y, radius, angle, ratio):
1155 x2 = x + math.cos(angle) * radius
1156 y2 = y + math.sin(angle) * radius
1158 x2 = x + math.sin(angle) * radius * ratio
1159 y2 = y + math.cos(angle) * radius * ratio
1163 def get_minor_axis(x, y, radius, angle, ratio):
1165 x2 = x + math.sin(angle) * radius * ratio
1166 y2 = y + math.cos(angle) * radius * ratio
1168 x2 = x + math.cos(angle) * radius
1169 y2 = y + math.sin(angle) * radius
1173 def wc_rad_func(x, y, radius, angle, ratio):
1174 x2, y2 = get_major_axis(x, y, radius, angle, ratio)
1177 def wc_angle_func(x, y, radius, angle, ratio):
1178 x2, y2 = get_major_axis(x, y, radius, angle, ratio)
1181 def wc_ratio_func(x, y, radius, angle, ratio):
1182 minor_x, minor_y = get_minor_axis(x, y, radius, angle, ratio)
1185 major_x, major_y = get_major_axis(x, y, radius, angle, ratio)
1188 return minor_angle / major_angle
1194 return ra, dec, wc_rad, wc_angle, wc_ratio