SourceXtractorPlusPlus  0.16
Please provide a description of the project.
model_fitting.py
Go to the documentation of this file.
1 # -*- coding: utf-8 -*-
2 
3 # Copyright © 2019 Université de Genève, LMU Munich - Faculty of Physics, IAP-CNRS/Sorbonne Université
4 #
5 # This library is free software; you can redistribute it and/or modify it under
6 # the terms of the GNU Lesser General Public License as published by the Free
7 # Software Foundation; either version 3.0 of the License, or (at your option)
8 # any later version.
9 #
10 # This library is distributed in the hope that it will be useful, but WITHOUT
11 # ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
12 # FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
13 # details.
14 #
15 # You should have received a copy of the GNU Lesser General Public License
16 # along with this library; if not, write to the Free Software Foundation, Inc.,
17 # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18 from __future__ import division, print_function
19 
20 import math
21 import sys
22 import warnings
23 from enum import Enum
24 
25 import _SourceXtractorPy as cpp
26 try:
27  import pyston
28 except ImportError:
29  warnings.warn('Could not import pyston: running outside sourcextractor?', ImportWarning)
30 from astropy import units as u
31 from astropy.coordinates import SkyCoord
32 
33 from .measurement_images import MeasurementGroup
34 
35 
36 class RangeType(Enum):
37  LINEAR = 1
38  EXPONENTIAL = 2
39 
40 
41 class Range(object):
42  r"""
43  Limit, and normalize, the range of values for a model fitting parameter.
44 
45 
46  Parameters
47  ----------
48  limits : a tuple (min, max), or a callable that receives a source, and returns a tuple (min, max)
49  type : RangeType
50 
51  Notes
52  -----
53  RangeType.LINEAR
54  Normalized to engine space using a sigmoid function
55 
56  .. math::
57 
58  engine = \ln \frac{world - min}{max-world} \\
59  world = min + \frac{max - min}{1 + e^{engine}}
60 
61  RangeType.EXPONENTIAL
62  Normalized to engine space using an exponential sigmoid function
63 
64  .. math::
65 
66  engine = \ln \left( \frac{\ln(world/min)}{\ln(max /world)} \right) \\
67  world = min * e^\frac{ \ln(max / min) }{ (1 + e^{-engine}) }
68 
69  """
70 
71  def __init__(self, limits, type):
72  """
73  Constructor.
74  """
75  self.__limits__limits = limits
76  self.__callable__callable = limits if hasattr(limits, '__call__') else lambda v, o: limits
77  self.__type__type = type
78 
79  def get_min(self, v, o):
80  """
81  Parameters
82  ----------
83  v : initial value
84  o : object being fitted
85 
86  Returns
87  -------
88  The minimum acceptable value for the range
89  """
90  return self.__callable__callable(v, o)[0]
91 
92  def get_max(self, v, o):
93  """
94  Parameters
95  ----------
96  v : initial value
97  o : object being fitted
98 
99  Returns
100  -------
101  The maximum acceptable value for the range
102  """
103  return self.__callable__callable(v, o)[1]
104 
105  def get_type(self):
106  """
107  Returns
108  -------
109  RangeType
110  """
111  return self.__type__type
112 
113  def __str__(self):
114  """
115  Returns
116  -------
117  str
118  Human readable representation for the object
119  """
120  res = '['
121  if hasattr(self.__limits__limits, '__call__'):
122  res += 'func'
123  else:
124  res += '{},{}'.format(self.__limits__limits[0], self.__limits__limits[1])
125  type_str = {
126  RangeType.LINEAR : 'LIN',
127  RangeType.EXPONENTIAL : 'EXP'
128  }
129  res += ',{}]'.format(type_str[self.__type__type])
130  return res
131 
132 
133 class Unbounded(object):
134  """
135  Unbounded, but normalize, value of a model fitting parameter
136 
137  Parameters
138  ----------
139  normalization_factor: float, or a callable that receives the initial value parameter value and a source,
140  and returns a float
141  The world value which will be normalized to 1 in engine coordinates
142  """
143 
144  def __init__(self, normalization_factor=1):
145  """
146  Constructor.
147  """
148  self.__normalization_factor__normalization_factor = normalization_factor
149  if hasattr(normalization_factor, '__call__'):
150  self.__callable__callable = normalization_factor
151  else:
152  self.__callable__callable = lambda v, o: normalization_factor
153 
154  def get_normalization_factor(self, v, o):
155  """
156  Parameters
157  ----------
158  v : initial value
159  o : object being fitted
160 
161  Returns
162  -------
163  float
164  The world value which will be normalized to 1 in engine coordinates
165  """
166  return self.__callable__callable(v, o)
167 
168  def __str__(self):
169  """
170  Returns
171  -------
172  str
173  Human readable representation for the object
174  """
175  res = '['
176  if hasattr(self.__normalization_factor__normalization_factor, '__call__'):
177  res += 'func'
178  else:
179  res += '{}'.format(self.__normalization_factor__normalization_factor)
180  res += ']'
181  return res
182 
183 
184 constant_parameter_dict = {}
185 free_parameter_dict = {}
186 dependent_parameter_dict = {}
187 
188 
189 def print_parameters(file=sys.stderr):
190  """
191  Print a human-readable representation of the configured model fitting parameters.
192 
193  Parameters
194  ----------
195  file : file object
196  Where to print the representation. Defaults to sys.stderr
197  """
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)
210 
211 
212 class ParameterBase(cpp.Id):
213  """
214  Base class for all model fitting parameter types.
215  Can not be used directly.
216  """
217 
218  def __str__(self):
219  """
220  Returns
221  -------
222  str
223  Human readable representation for the object
224  """
225  return '(ID:{})'.format(self.id)
226 
227 
229  """
230  A parameter with a single value that remains constant. It will not be fitted.
231 
232  Parameters
233  ----------
234  value : float, or callable that receives a source and returns a float
235  Value for this parameter
236  """
237 
238  def __init__(self, value):
239  """
240  Constructor.
241  """
242  ParameterBase.__init__(self)
243  self.__value__value = value
244  self.__callable__callable = value if hasattr(value, '__call__') else lambda o: value
245  constant_parameter_dict[self.id] = self
246 
247  def get_value(self, o):
248  """
249  Parameters
250  ----------
251  o : object being fitted
252 
253  Returns
254  -------
255  float
256  Value of the constant parameter
257  """
258  return self.__callable__callable(o)
259 
260  def __str__(self):
261  """
262  Returns
263  -------
264  str
265  Human readable representation for the object
266  """
267  res = ParameterBase.__str__(self)[:-1] + ', value:'
268  if hasattr(self.__value__value, '__call__'):
269  res += 'func'
270  else:
271  res += str(self.__value__value)
272  return res + ')'
273 
274 
276  """
277  A parameter that will be fitted by the model fitting engine.
278 
279  Parameters
280  ----------
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.
285 
286  See Also
287  --------
288  Unbounded
289  Range
290 
291  Examples
292  --------
293  >>> sersic = FreeParameter(2.0, Range((1.0, 7.0), RangeType.LINEAR))
294  """
295 
296  def __init__(self, init_value, range=Unbounded()):
297  """
298  Constructor.
299  """
300  ParameterBase.__init__(self)
301  self.__init_value__init_value = init_value
302  self.__init_call__init_call = init_value if hasattr(init_value, '__call__') else lambda o: init_value
303  self.__range__range = range
304  free_parameter_dict[self.id] = self
305 
306  def get_init_value(self, o):
307  """
308  Parameters
309  ----------
310  o : object being fitted
311 
312  Returns
313  -------
314  float
315  Initial value for the free parameter
316  """
317  return self.__init_call__init_call(o)
318 
319  def get_range(self):
320  """
321  Returns
322  -------
323  Unbounded or Range
324  """
325  return self.__range__range
326 
327  def __str__(self):
328  """
329  Returns
330  -------
331  str
332  Human readable representation for the object
333  """
334  res = ParameterBase.__str__(self)[:-1] + ', init:'
335  if hasattr(self.__init_value__init_value, '__call__'):
336  res += 'func'
337  else:
338  res += str(self.__init_value__init_value)
339  if self.__range__range:
340  res += ', range:' + str(self.__range__range)
341  return res + ')'
342 
343 
345  """
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
348 
349  Parameters
350  ----------
351  func : callable
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.
356 
357  Examples
358  --------
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)
362  """
363 
364  def __init__(self, func, *params):
365  """
366  Constructor.
367  """
368  ParameterBase.__init__(self)
369  self.funcfunc = func
370  self.paramsparams = [p.id for p in params]
371  dependent_parameter_dict[self.id] = self
372 
373 
375  """
376  Convenience function for the position parameter X and Y.
377 
378  Returns
379  -------
380  x : FreeParameter
381  X coordinate, starting at the X coordinate of the centroid and linearly limited to X +/- the object radius.
382  y : FreeParameter
383  Y coordinate, starting at the Y coordinate of the centroid and linearly limited to Y +/- the object radius.
384  Notes
385  -----
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.
388  """
389  param_range = Range(lambda v,o: (v-o.radius, v+o.radius), RangeType.LINEAR)
390  return (
391  FreeParameter(lambda o: o.centroid_x, param_range),
392  FreeParameter(lambda o: o.centroid_y, param_range)
393  )
394 
395 
397  """
398  Possible flux types to use as initial value for the flux parameter.
399  Right now, only isophotal is supported.
400  """
401  ISO = 1
402 
403 
404 def get_flux_parameter(type=FluxParameterType.ISO, scale=1):
405  """
406  Convenience function for the flux parameter.
407 
408  Parameters
409  ----------
410  type : int
411  One of the values defined in FluxParameterType
412  scale : float
413  Scaling of the initial flux. Defaults to 1.
414 
415  Returns
416  -------
417  flux : FreeParameter
418  Flux parameter, starting at the flux defined by `type`, and limited to +/- 1e3 times the initial value.
419  """
420  attr_map = {
421  FluxParameterType.ISO : 'isophotal_flux'
422  }
423  return FreeParameter(lambda o: getattr(o, attr_map[type]) * scale, Range(lambda v,o: (v * 1E-3, v * 1E3), RangeType.EXPONENTIAL))
424 
425 
426 prior_dict = {}
427 
428 
429 class Prior(cpp.Id):
430  """
431  Model a Gaussian prior on a given parameter.
432 
433  Parameters
434  ----------
435  param : ParameterBase
436  Model fitting parameter
437  value : float or callable that receives a source and returns a float
438  Mean of the Gaussian
439  sigma : float or callable that receives a source and returns a float
440  Standard deviation of the Gaussian
441  """
442 
443  def __init__(self, param, value, sigma):
444  """
445  Constructor.
446  """
447  cpp.Id.__init__(self)
448  self.paramparam = param.id
449  self.valuevalue = value if hasattr(value, '__call__') else lambda o: value
450  self.sigmasigma = sigma if hasattr(sigma, '__call__') else lambda o: sigma
451 
452 
453 def add_prior(param, value, sigma):
454  """
455  Add a prior to the given parameter.
456 
457  Parameters
458  ----------
459  param : ParameterBase
460  value : float or callable that receives a source and returns a float
461  Mean of the Gaussian
462  sigma : float or callable that receives a source and returns a float
463  Standard deviation of the Gaussian
464  """
465  prior = Prior(param, value, sigma)
466  prior_dict[prior.id] = prior
467 
468 
469 frame_models_dict = {}
470 
471 
472 def _set_model_to_frames(group, model):
473  for x in group:
474  if isinstance(x, tuple):
475  _set_model_to_frames(x[1], model)
476  else:
477  if not x.id in frame_models_dict:
478  frame_models_dict[x.id] = []
479  frame_models_dict[x.id].append(model.id)
480 
481 
482 def add_model(group, model):
483  """
484  Add a model to be fitted to the given group.
485 
486  Parameters
487  ----------
488  group : MeasurementGroup
489  model : ModelBase
490  """
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'):
494  group.models = []
495  group.models.append(model)
496  _set_model_to_frames(group, model)
497 
498 
499 def print_model_fitting_info(group, show_params=False, prefix='', file=sys.stderr):
500  """
501  Print a human-readable representation of the configured models.
502 
503  Parameters
504  ----------
505  group : MeasurementGroup
506  Print the models for this group.
507  show_params : bool
508  If True, print also the parameters that belong to the model
509  prefix : str
510  Prefix each line with this string. Used internally for indentation.
511  file : file object
512  Where to print the representation. Defaults to sys.stderr
513  """
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)
518  for x in group:
519  if isinstance(x, tuple):
520  print('{}{}:'.format(prefix, x[0]), file=file)
521  print_model_fitting_info(x[1], show_params, prefix + ' ', file=file)
522 
523 
524 constant_model_dict = {}
525 point_source_model_dict = {}
526 sersic_model_dict = {}
527 exponential_model_dict = {}
528 de_vaucouleurs_model_dict = {}
529 onnx_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}
532 
533 
534 def set_max_iterations(iterations):
535  """
536  Parameters
537  ----------
538  iterations : int
539  Max number of iterations for the model fitting.
540  """
541  params_dict["max_iterations"] = iterations
542 
543 
545  """
546  Parameters
547  ----------
548  scale : float
549  Sets u0, as used by the modified chi squared residual comparator, a function that reduces the effect of large
550  deviations.
551  Refer to the SourceXtractor++ documentation for a better explanation of how residuals are computed and how
552  this value affects the model fitting.
553  """
554  params_dict["modified_chi_squared_scale"] = scale
555 
556 
557 def set_engine(engine):
558  """
559  Parameters
560  ----------
561  engine : str
562  Minimization engine for the model fitting : levmar or gsl
563  """
564  params_dict["engine"] = engine
565 
566 def use_iterative_fitting(use_iterative_fitting):
567  """
568  Parameters
569  ----------
570  use_iterative_fitting : boolean
571  use iterative model fitting or legacy
572  """
573  params_dict["use_iterative_fitting"] = use_iterative_fitting
574 
575 def set_meta_iterations(meta_iterations):
576  """
577  Parameters
578  ----------
579  meta_iterations : int
580  number of meta iterations on the whole group (when using iterative model fitting)
581  """
582  params_dict["meta_iterations"] = meta_iterations
583 
584 def set_deblend_factor(deblend_factor):
585  """
586  Parameters
587  ----------
588 
589  """
590  params_dict["deblend_factor"] = deblend_factor
591 
592 def set_meta_iteration_stop(meta_iteration_stop):
593  """
594  Parameters
595  ----------
596 
597  """
598  params_dict["meta_iteration_stop"] = meta_iteration_stop
599 
600 
601 
602 
603 class ModelBase(cpp.Id):
604  """
605  Base class for all models.
606  """
607  pass
608 
609 
610 class CoordinateModelBase(ModelBase):
611  """
612  Base class for positioned models with a flux. It can not be used directly.
613 
614  Parameters
615  ----------
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
621  Total flux
622  """
623 
624  def __init__(self, x_coord, y_coord, flux):
625  """
626  Constructor.
627  """
628  ModelBase.__init__(self)
629  self.x_coordx_coord = x_coord if isinstance(x_coord, ParameterBase) else ConstantParameter(x_coord)
630  self.y_coordy_coord = y_coord if isinstance(y_coord, ParameterBase) else ConstantParameter(y_coord)
631  self.fluxflux = flux if isinstance(flux, ParameterBase) else ConstantParameter(flux)
632 
633 
635  """
636  Models a source as a point, spread by the PSF.
637 
638  Parameters
639  ----------
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
645  Total flux
646  """
647 
648  def __init__(self, x_coord, y_coord, flux):
649  """
650  Constructor.
651  """
652  CoordinateModelBase.__init__(self, x_coord, y_coord, flux)
653  global point_source_model_dict
654  point_source_model_dict[self.id] = self
655 
656  def to_string(self, show_params=False):
657  """
658  Return a human readable representation of the model.
659 
660  Parameters
661  ----------
662  show_params: bool
663  If True, include information about the parameters.
664 
665  Returns
666  -------
667  str
668  """
669  if show_params:
670  return 'PointSource[x_coord={}, y_coord={}, flux={}]'.format(
671  self.x_coordx_coord, self.y_coordy_coord, self.fluxflux)
672  else:
673  return 'PointSource[x_coord={}, y_coord={}, flux={}]'.format(
674  self.x_coordx_coord.id, self.y_coordy_coord.id, self.fluxflux.id)
675 
676 
678  """
679  A model that is constant through all the image.
680 
681  Parameters
682  ----------
683  value: ParameterBase or float
684  Value to add to the value of all pixels from the model.
685  """
686 
687  def __init__(self, value):
688  """
689  Constructor.
690  """
691  ModelBase.__init__(self)
692  self.valuevalue = value if isinstance(value, ParameterBase) else ConstantParameter(value)
693  global constant_model_dict
694  constant_model_dict[self.id] = self
695 
696  def to_string(self, show_params=False):
697  """
698  Return a human readable representation of the model.
699 
700  Parameters
701  ----------
702  show_params: bool
703  If True, include information about the parameters.
704 
705  Returns
706  -------
707  str
708  """
709  if show_params:
710  return 'ConstantModel[value={}]'.format(self.valuevalue)
711  else:
712  return 'ConstantModel[value={}]'.format(self.valuevalue.id)
713 
714 
716  """
717  Base class for the Sersic, Exponential and de Vaucouleurs models. It can not be used directly.
718 
719  Parameters
720  ----------
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
726  Total flux
727  effective_radius : ParameterBase or float
728  Ellipse semi-major axis, in pixels on the detection image.
729  aspect_ratio : ParameterBase or float
730  Ellipse ratio.
731  angle : ParameterBase or float
732  Ellipse rotation, in radians
733  """
734 
735  def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle):
736  """
737  Constructor.
738  """
739  CoordinateModelBase.__init__(self, x_coord, y_coord, flux)
740  self.effective_radiuseffective_radius = effective_radius if isinstance(effective_radius, ParameterBase) else ConstantParameter(effective_radius)
741  self.aspect_ratioaspect_ratio = aspect_ratio if isinstance(aspect_ratio, ParameterBase) else ConstantParameter(aspect_ratio)
742  self.angleangle = angle if isinstance(angle, ParameterBase) else ConstantParameter(angle)
743 
744 
746  """
747  Model a source with a Sersic profile.
748 
749  Parameters
750  ----------
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
756  Total flux
757  effective_radius : ParameterBase or float
758  Ellipse semi-major axis, in pixels on the detection image.
759  aspect_ratio : ParameterBase or float
760  Ellipse ratio.
761  angle : ParameterBase or float
762  Ellipse rotation, in radians
763  n : ParameterBase or float
764  Sersic index
765  """
766 
767  def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle, n):
768  """
769  Constructor.
770  """
771  SersicModelBase.__init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle)
772  self.nn = n if isinstance(n, ParameterBase) else ConstantParameter(n)
773  global sersic_model_dict
774  sersic_model_dict[self.id] = self
775 
776  def to_string(self, show_params=False):
777  """
778  Return a human readable representation of the model.
779 
780  Parameters
781  ----------
782  show_params: bool
783  If True, include information about the parameters.
784 
785  Returns
786  -------
787  str
788  """
789  if show_params:
790  return 'Sersic[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}, n={}]'.format(
791  self.x_coordx_coord, self.y_coordy_coord, self.fluxflux, self.effective_radiuseffective_radius, self.aspect_ratioaspect_ratio, self.angleangle, self.nn)
792  else:
793  return 'Sersic[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}, n={}]'.format(
794  self.x_coordx_coord.id, self.y_coordy_coord.id, self.fluxflux.id, self.effective_radiuseffective_radius.id, self.aspect_ratioaspect_ratio.id, self.angleangle.id, self.nn.id)
795 
796 
798  """
799  Model a source with an exponential profile (Sersic model with an index of 1)
800 
801  Parameters
802  ----------
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
808  Total flux
809  effective_radius : ParameterBase or float
810  Ellipse semi-major axis, in pixels on the detection image.
811  aspect_ratio : ParameterBase or float
812  Ellipse ratio.
813  angle : ParameterBase or float
814  Ellipse rotation, in radians
815  """
816 
817  def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle):
818  """
819  Constructor.
820  """
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
824 
825  def to_string(self, show_params=False):
826  """
827  Return a human readable representation of the model.
828 
829  Parameters
830  ----------
831  show_params: bool
832  If True, include information about the parameters.
833 
834  Returns
835  -------
836  str
837  """
838  if show_params:
839  return 'Exponential[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
840  self.x_coordx_coord, self.y_coordy_coord, self.fluxflux, self.effective_radiuseffective_radius, self.aspect_ratioaspect_ratio, self.angleangle)
841  else:
842  return 'Exponential[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
843  self.x_coordx_coord.id, self.y_coordy_coord.id, self.fluxflux.id, self.effective_radiuseffective_radius.id, self.aspect_ratioaspect_ratio.id, self.angleangle.id)
844 
845 
847  """
848  Model a source with a De Vaucouleurs profile (Sersic model with an index of 4)
849 
850  Parameters
851  ----------
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
857  Total flux
858  effective_radius : ParameterBase or float
859  Ellipse semi-major axis, in pixels on the detection image.
860  aspect_ratio : ParameterBase or float
861  Ellipse ratio.
862  angle : ParameterBase or float
863  Ellipse rotation, in radians
864  """
865 
866  def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle):
867  """
868  Constructor.
869  """
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
873 
874  def to_string(self, show_params=False):
875  """
876  Return a human readable representation of the model.
877 
878  Parameters
879  ----------
880  show_params: bool
881  If True, include information about the parameters.
882 
883  Returns
884  -------
885  str
886  """
887  if show_params:
888  return 'DeVaucouleurs[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
889  self.x_coordx_coord, self.y_coordy_coord, self.fluxflux, self.effective_radiuseffective_radius, self.aspect_ratioaspect_ratio, self.angleangle)
890  else:
891  return 'DeVaucouleurs[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
892  self.x_coordx_coord.id, self.y_coordy_coord.id, self.fluxflux.id, self.effective_radiuseffective_radius.id, self.aspect_ratioaspect_ratio.id, self.angleangle.id)
893 
895  """
896  ComputeGraphModel model
897 
898  Parameters
899  ----------
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
907  Total flux
908  params : Dictionary of String and ParameterBase or float
909  Dictionary of custom parameters for the ONNX model
910  """
911 
912  def __init__(self, models, x_coord, y_coord, flux, params={}):
913  """
914  Constructor.
915  """
916  CoordinateModelBase.__init__(self, x_coord, y_coord, flux)
917 
918  ratio_name = "_aspect_ratio"
919  angle_name = "_angle"
920  scale_name = "_scale"
921 
922  for k in params.keys():
923  if not isinstance(params[k], ParameterBase):
924  params[k] = ConstantParameter(params[k])
925 
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
929 
930  self.aspect_ratioaspect_ratio = aspect_ratio if isinstance(aspect_ratio, ParameterBase) else ConstantParameter(aspect_ratio)
931  self.angleangle = angle if isinstance(angle, ParameterBase) else ConstantParameter(angle)
932  self.scalescale = scale if isinstance(scale, ParameterBase) else ConstantParameter(scale)
933 
934  params.pop(ratio_name, None)
935  params.pop(angle_name, None)
936  params.pop(scale_name, None)
937 
938  self.paramsparams = params
939  self.modelsmodels = models if isinstance(models, list) else [models]
940 
941  global onnx_model_dict
942  onnx_model_dict[self.id] = self
943 
944  def to_string(self, show_params=False):
945  """
946  Return a human readable representation of the model.
947 
948  Parameters
949  ----------
950  show_params: bool
951  If True, include information about the parameters.
952 
953  Returns
954  -------
955  str
956  """
957  if show_params:
958  return 'ComputeGraph[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
959  self.x_coordx_coord, self.y_coordy_coord, self.fluxflux, self.effective_radius, self.aspect_ratioaspect_ratio, self.angleangle)
960  else:
961  return 'ComputeGraph[x_coord={}, y_coord={}, flux={}, effective_radius={}, aspect_ratio={}, angle={}]'.format(
962  self.x_coordx_coord.id, self.y_coordy_coord.id, self.fluxflux.id, self.effective_radius.id, self.aspect_ratioaspect_ratio.id, self.angleangle.id)
963 
964 
966  """
967  Coordinates in right ascension and declination
968 
969  Parameters
970  ----------
971  ra : float
972  Right ascension
973  dec : float
974  Declination
975  """
976 
977  def __init__(self, ra, dec):
978  """
979  Constructor.
980  """
981  self.rara = ra
982  self.decdec = dec
983 
984 
986  """
987  Transform an (X, Y) in pixel coordinates on the detection image to (RA, DEC) in world coordinates.
988  Parameters
989  ----------
990  x : float
991  y : float
992 
993  Returns
994  -------
995  WorldCoordinate
996  """
997  # -1 as we go FITS -> internal
998  wc_alpha = pyston.image_to_world_alpha(x - 1, y - 1)
999  wc_delta = pyston.image_to_world_delta(x - 1, y - 1)
1000  return WorldCoordinate(wc_alpha, wc_delta)
1001 
1002 
1003 def get_sky_coord(x, y):
1004  """
1005  Transform an (X, Y) in pixel coordinates on the detection image to astropy SkyCoord.
1006 
1007  Parameters
1008  ----------
1009  x : float
1010  y : float
1011 
1012  Returns
1013  -------
1014  SkyCoord
1015  """
1016  coord = pixel_to_world_coordinate(x, y)
1017  sky_coord = SkyCoord(ra=coord.ra*u.degree, dec=coord.dec*u.degree)
1018  return sky_coord
1019 
1020 
1021 def radius_to_wc_angle(x, y, rad):
1022  """
1023  Transform a radius in pixels on the detection image to a radius in sky coordinates.
1024 
1025  Parameters
1026  ----------
1027  x : float
1028  y : float
1029  rad : float
1030 
1031  Returns
1032  -------
1033  Radius in degrees
1034  """
1035  return (get_separation_angle(x, y, x+rad, y) + get_separation_angle(x, y, x, y+rad)) / 2.0
1036 
1037 
1038 def get_separation_angle(x1, y1, x2, y2):
1039  """
1040  Get the separation angle in sky coordinates for two points defined in pixels on the detection image.
1041 
1042  Parameters
1043  ----------
1044  x1 : float
1045  y1 : float
1046  x2 : float
1047  y2 : float
1048 
1049  Returns
1050  -------
1051  Separation in degrees
1052  """
1053  c1 = get_sky_coord(x1, y1)
1054  c2 = get_sky_coord(x2, y2)
1055  return c1.separation(c2).degree
1056 
1057 
1058 def get_position_angle(x1, y1, x2, y2):
1059  """
1060  Get the position angle in sky coordinates for two points defined in pixels on the detection image.
1061 
1062  Parameters
1063  ----------
1064  x1
1065  y1
1066  x2
1067  y2
1068 
1069  Returns
1070  -------
1071  Position angle in degrees, normalized to -/+ 90
1072  """
1073  c1 = get_sky_coord(x1, y1)
1074  c2 = get_sky_coord(x2, y2)
1075  angle = c1.position_angle(c2).degree
1076 
1077  # return angle normalized to range: -90 <= angle < 90
1078  return (angle + 90.0) % 180.0 - 90.0
1079 
1080 
1082  """
1083  Convenience function for generating two dependent parameter with world (alpha, delta) coordinates
1084  from image (X, Y) coordinates.
1085 
1086  Parameters
1087  ----------
1088  x : ParameterBase
1089  y : ParameterBase
1090 
1091  Returns
1092  -------
1093  ra : DependentParameter
1094  dec : DependentParameter
1095 
1096  See Also
1097  --------
1098  get_pos_parameters
1099 
1100  Examples
1101  --------
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)
1106  """
1107  ra = DependentParameter(lambda x,y: pixel_to_world_coordinate(x, y).ra, x, y)
1108  dec = DependentParameter(lambda x,y: pixel_to_world_coordinate(x, y).dec, x, y)
1109  return (ra, dec)
1110 
1111 
1112 def get_world_parameters(x, y, radius, angle, ratio):
1113  """
1114  Convenience function for generating five dependent parameters, in world coordinates, for the position
1115  and shape of a model.
1116 
1117  Parameters
1118  ----------
1119  x : ParameterBase
1120  y : ParameterBase
1121  radius : ParameterBase
1122  angle : ParameterBase
1123  ratio : ParameterBase
1124 
1125  Returns
1126  -------
1127  ra : DependentParameter
1128  Right ascension
1129  dec : DependentParameter
1130  Declination
1131  rad : DependentParameter
1132  Radius as degrees
1133  angle : DependentParameter
1134  Angle in degrees
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
1138 
1139  Examples
1140  --------
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)
1149  """
1150  ra = DependentParameter(lambda x,y: pixel_to_world_coordinate(x, y).ra, x, y)
1151  dec = DependentParameter(lambda x,y: pixel_to_world_coordinate(x, y).dec, x, y)
1152 
1153  def get_major_axis(x, y, radius, angle, ratio):
1154  if ratio <= 1:
1155  x2 = x + math.cos(angle) * radius
1156  y2 = y + math.sin(angle) * radius
1157  else:
1158  x2 = x + math.sin(angle) * radius * ratio
1159  y2 = y + math.cos(angle) * radius * ratio
1160 
1161  return x2, y2
1162 
1163  def get_minor_axis(x, y, radius, angle, ratio):
1164  if ratio <= 1:
1165  x2 = x + math.sin(angle) * radius * ratio
1166  y2 = y + math.cos(angle) * radius * ratio
1167  else:
1168  x2 = x + math.cos(angle) * radius
1169  y2 = y + math.sin(angle) * radius
1170 
1171  return x2, y2
1172 
1173  def wc_rad_func(x, y, radius, angle, ratio):
1174  x2, y2 = get_major_axis(x, y, radius, angle, ratio)
1175  return get_separation_angle(x, y, x2, y2)
1176 
1177  def wc_angle_func(x, y, radius, angle, ratio):
1178  x2, y2 = get_major_axis(x, y, radius, angle, ratio)
1179  return get_position_angle(x, y, x2, y2)
1180 
1181  def wc_ratio_func(x, y, radius, angle, ratio):
1182  minor_x, minor_y = get_minor_axis(x, y, radius, angle, ratio)
1183  minor_angle = get_separation_angle(x,y, minor_x, minor_y)
1184 
1185  major_x, major_y = get_major_axis(x, y, radius, angle, ratio)
1186  major_angle = get_separation_angle(x,y, major_x, major_y)
1187 
1188  return minor_angle / major_angle
1189 
1190  wc_rad = DependentParameter(wc_rad_func, x, y, radius, angle, ratio)
1191  wc_angle = DependentParameter(wc_angle_func, x, y, radius, angle, ratio)
1192  wc_ratio = DependentParameter(wc_ratio_func, x, y, radius, angle, ratio)
1193 
1194  return ra, dec, wc_rad, wc_angle, wc_ratio
1195 
def __init__(self, models, x_coord, y_coord, flux, params={})
def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle)
def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle)
def __init__(self, init_value, range=Unbounded())
def __init__(self, x_coord, y_coord, flux)
def __init__(self, param, value, sigma)
def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle)
def __init__(self, x_coord, y_coord, flux, effective_radius, aspect_ratio, angle, n)
def __init__(self, normalization_factor=1)
def set_meta_iterations(meta_iterations)
def set_deblend_factor(deblend_factor)
def get_position_angle(x1, y1, x2, y2)
def get_world_parameters(x, y, radius, angle, ratio)
def print_parameters(file=sys.stderr)
def use_iterative_fitting(use_iterative_fitting)
def get_separation_angle(x1, y1, x2, y2)
def get_flux_parameter(type=FluxParameterType.ISO, scale=1)
def add_prior(param, value, sigma)
def print_model_fitting_info(group, show_params=False, prefix='', file=sys.stderr)
def set_meta_iteration_stop(meta_iteration_stop)