SourceXtractorPlusPlus  0.13
Please provide a description of the project.
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
measurement_images.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 os
21 import re
22 import sys
23 
24 import _SourceXtractorPy as cpp
25 
26 if sys.version_info.major < 3:
27  from StringIO import StringIO
28 else:
29  from io import StringIO
30 
31 
32 measurement_images = {}
33 
34 class FitsFile(cpp.FitsFile):
35  def __init__(self, filename):
36  super(FitsFile, self).__init__(str(filename))
37  self.hdu_list = [i for i in self.image_hdus]
38 
39  def __iter__(self):
40  return iter(self.hdu_list)
41 
42  def get_headers(self, hdu):
43  d = {}
44  headers = super(FitsFile, self).get_headers(hdu)
45 
46  try:
47  it = iter(headers)
48  while it:
49  a = next(it)
50  d[a.key()] = headers[a.key()]
51  except StopIteration:
52  pass
53 
54  return d
55 
56 class MeasurementImage(cpp.MeasurementImage):
57  """
58  A MeasurementImage is the processing unit for SourceXtractor++. Measurements and model fitting can be done
59  over one, or many, of them. It models the image, plus its associated weight file, PSF, etc.
60 
61  Parameters
62  ----------
63  fits_file : str or FitsFile object
64  The path to a FITS image, or an instance of FitsFile
65  psf_file : str
66  The path to a PSF. It can be either a FITS image, or a PSFEx model.
67  weight_file : str or FitsFile
68  The path to a FITS image with the pixel weights, or an instance of FitsFile
69  gain : float
70  Image gain. If None, `gain_keyword` will be used instead.
71  gain_keyword : str
72  Keyword for the header containing the gain.
73  saturation : float
74  Saturation value. If None, `saturation_keyword` will be used instead.
75  saturation_keyword : str
76  Keyword for the header containing the saturation value.
77  flux_scale : float
78  Flux scaling. Each pixel value will be multiplied by this. If None, `flux_scale_keyword` will be used
79  instead.
80  flux_scale_keyword : str
81  Keyword for the header containing the flux scaling.
82  weight_type : str
83  The type of the weight image. It must be one of:
84 
85  - none
86  The image itself is used to compute internally a constant variance (default)
87  - background
88  The image itself is used to compute internally a variance map
89  - rms
90  The weight image must contain a weight-map in units of absolute standard deviations
91  (in ADUs per pixel).
92  - variance
93  The weight image must contain a weight-map in units of relative variance.
94  - weight
95  The weight image must contain a weight-map in units of relative weights. The data are converted
96  to variance units.
97  weight_absolute : bool
98  If False, the weight map will be scaled according to an absolute variance map built from the image itself.
99  weight_scaling : float
100  Apply an scaling to the weight map.
101  weight_threshold : float
102  Pixels with weights beyond this value are treated just like pixels discarded by the masking process.
103  constant_background : float
104  If set a constant background of that value is assumed for the image instead of using automatic detection
105  image_hdu : int
106  For multi-extension FITS file specifies the HDU number for the image. Default 1 (primary HDU)
107  psf_hdu : int
108  For multi-extension FITS file specifies the HDU number for the psf. Defaults to the same value as image_hdu
109  weight_hdu : int
110  For multi-extension FITS file specifies the HDU number for the weight. Defaults to the same value as image_hdu
111  """
112 
113  def __init__(self, fits_file, psf_file=None, weight_file=None, gain=None,
114  gain_keyword='GAIN', saturation=None, saturation_keyword='SATURATE',
115  flux_scale=None, flux_scale_keyword='FLXSCALE',
116  weight_type='none', weight_absolute=False, weight_scaling=1.,
117  weight_threshold=None, constant_background=None,
118  image_hdu=1, psf_hdu=None, weight_hdu=None
119  ):
120  """
121  Constructor.
122  """
123  if isinstance(fits_file, FitsFile):
124  hdu_list = fits_file
125  file_path = fits_file.filename
126  else:
127  hdu_list = FitsFile(fits_file)
128  file_path = fits_file
129 
130  if isinstance(weight_file, FitsFile):
131  weight_file = weight_file.filename
132 
133  super(MeasurementImage, self).__init__(os.path.abspath(file_path),
134  os.path.abspath(psf_file) if psf_file else '',
135  os.path.abspath(weight_file) if weight_file else '')
136 
137  if image_hdu <= 0 or (weight_hdu is not None and weight_hdu <= 0) or (psf_hdu is not None and psf_hdu <= 0):
138  raise ValueError('HDU indexes start at 1')
139 
140  self.meta = {
141  'IMAGE_FILENAME': self.file,
142  'PSF_FILENAME': self.psf_file,
143  'WEIGHT_FILENAME': self.weight_file
144  }
145 
146  self.meta.update(hdu_list.get_headers(image_hdu-1))
147 
148  if gain is not None:
149  self.gain = gain
150  elif gain_keyword in self.meta:
151  self.gain = float(self.meta[gain_keyword])
152  else:
153  self.gain = 0.
154 
155  if saturation is not None:
156  self.saturation = saturation
157  elif saturation_keyword in self.meta:
158  self.saturation = float(self.meta[saturation_keyword])
159  else:
160  self.saturation = 0.
161 
162  if flux_scale is not None:
163  self.flux_scale = flux_scale
164  elif flux_scale_keyword in self.meta:
165  self.flux_scale = float(self.meta[flux_scale_keyword])
166  else:
167  self.flux_scale = 1.
168 
169  self.weight_type = weight_type
170  self.weight_absolute = weight_absolute
171  self.weight_scaling = weight_scaling
172  if weight_threshold is None:
173  self.has_weight_threshold = False
174  else:
175  self.has_weight_threshold = True
176  self.weight_threshold = weight_threshold
177 
178  if constant_background is not None:
180  self.constant_background_value = constant_background
181  else:
182  self.is_background_constant = False
183  self.constant_background_value = -1
184 
185  self.image_hdu = image_hdu
186 
187  if psf_hdu is None:
188  self.psf_hdu = image_hdu
189  else:
190  self.psf_hdu = psf_hdu
191 
192  if weight_hdu is None:
193  self.weight_hdu = image_hdu
194  else:
195  self.weight_hdu = weight_hdu
196 
197  global measurement_images
198  measurement_images[self.id] = self
199 
200  def __str__(self):
201  """
202  Returns
203  -------
204  str
205  Human readable representation for the object
206  """
207  return 'Image {}: {} / {}, PSF: {} / {}, Weight: {} / {}'.format(
208  self.id, self.meta['IMAGE_FILENAME'], self.image_hdu, self.meta['PSF_FILENAME'], self.psf_hdu,
209  self.meta['WEIGHT_FILENAME'], self.weight_hdu)
210 
211 
212 def print_measurement_images(file=sys.stderr):
213  """
214  Print a human-readable representation of the configured measurement images.
215 
216  Parameters
217  ----------
218  file : file object
219  Where to print the representation. Defaults to sys.stderr
220  """
221  print('Measurement images:', file=file)
222  for i in measurement_images:
223  im = measurement_images[i]
224  print('Image {}'.format(i), file=file)
225  print(' File: {}'.format(im.file), file=file)
226  print(' PSF: {}'.format(im.psf_file), file=file)
227  print(' Weight: {}'.format(im.weight_file), file=file)
228 
229 
230 class ImageGroup(object):
231  """
232  Models the grouping of images. Measurement can *not* be made directly on instances of this type.
233  The configuration must be "frozen" before creating a MeasurementGroup
234 
235  See Also
236  --------
237  MeasurementGroup
238  """
239 
240  def __init__(self, **kwargs):
241  """
242  Constructor. It is not recommended to be used directly. Use instead load_fits_image or load_fits_images.
243  """
244  self.__images = []
245  self.__subgroups = None
246  self.__subgroup_names = set()
247  if len(kwargs) != 1 or ('images' not in kwargs and 'subgroups' not in kwargs):
248  raise ValueError('ImageGroup only takes as parameter one of "images" or "subgroups"')
249  key = list(kwargs.keys())[0]
250  if key == 'images':
251  if isinstance(kwargs[key], list):
252  self.__images = kwargs[key]
253  else:
254  self.__images = [kwargs[key]]
255  if key == 'subgroups':
256  self.__subgroups = kwargs[key]
257  for name, _ in self.__subgroups:
258  self.__subgroup_names.add(name)
259 
260  def __len__(self):
261  """
262  See Also
263  --------
264  is_leaf
265 
266  Returns
267  -------
268  int
269  How may subgroups or images are there in this group
270  """
271  if self.__subgroups:
272  return len(self.__subgroups)
273  else:
274  return len(self.__images)
275 
276  def __iter__(self):
277  """
278  Allows to iterate on the contained subgroups or images
279 
280  See Also
281  --------
282  is_leaf
283 
284  Returns
285  -------
286  iterator
287  """
288  if self.__subgroups:
289  return self.__subgroups.__iter__()
290  else:
291  return self.__images.__iter__()
292 
293  def split(self, grouping_method):
294  """
295  Splits the group in various subgroups, applying a filter on the contained images. If the group has
296  already been split, applies the split to each subgroup.
297 
298  Parameters
299  ----------
300  grouping_method : callable
301  A callable that receives as a parameter the list of contained images, and returns
302  a list of tuples, with the grouping key value, and the list of grouped images belonging to the given key.
303 
304  See Also
305  --------
306  ByKeyword
307  ByPattern
308 
309  Raises
310  -------
311  ValueError
312  If some images have not been grouped by the callable.
313  """
314  if self.__subgroups:
315  #if we are already subgrouped, apply the split to the subgroups
316  for _, sub_group in self.__subgroups:
317  sub_group.split(grouping_method)
318  else:
319  subgrouped_images = grouping_method(self.__images)
320  if sum(len(p[1]) for p in subgrouped_images) != len(self.__images):
321  self.__subgroups = None
322  raise ValueError('Some images were not grouped')
323  self.__subgroups = []
324  for k, im_list in subgrouped_images:
325  assert k not in self.__subgroup_names
326  self.__subgroup_names.add(k)
327  self.__subgroups.append((k, ImageGroup(images=im_list)))
328  self.__images = []
329 
330  def add_images(self, images):
331  """
332  Add new images to the group.
333 
334  Parameters
335  ----------
336  images : list of, or a single, MeasurementImage
337 
338  Raises
339  ------
340  ValueError
341  If the group has been split, no new images can be added.
342  """
343  if self.__subgroups is not None:
344  raise ValueError('ImageGroup is already subgrouped')
345  if isinstance(images, MeasurementImage):
346  self.__images.append(images)
347  else:
348  self.__images.extend(images)
349 
350  def add_subgroup(self, name, group):
351  """
352  Add a subgroup to a group.
353 
354  Parameters
355  ----------
356  name : str
357  The new of the new group
358 
359  group : ImageGroup
360  """
361  if self.__subgroups is None:
362  raise Exception('ImageGroup is not subgrouped yet')
363  if name in self.__subgroup_names:
364  raise Exception('Subgroup {} alread exists'.format(name))
365  self.__subgroup_names.add(name)
366  self.__subgroups.append((name, group))
367 
368  def is_leaf(self):
369  """
370  Returns
371  -------
372  bool
373  True if the group is a leaf group
374  """
375  return self.__subgroups is None
376 
377  def __getitem__(self, name):
378  """
379  Get a subgroup.
380 
381  Parameters
382  ----------
383  name : str
384  The name of the subgroup.
385 
386  Returns
387  -------
388  ImageGroup
389  The matching group.
390 
391  Raises
392  ------
393  ValueError
394  If the group has not been split.
395  KeyError
396  If the group has not been found.
397  """
398  if self.__subgroups is None:
399  raise ValueError('ImageGroup is not subgrouped yet')
400  try:
401  return next(x for x in self.__subgroups if x[0] == name)[1]
402  except StopIteration:
403  raise KeyError('Group {} not found'.format(name))
404 
405  def print(self, prefix='', show_images=False, file=sys.stderr):
406  """
407  Print a human-readable representation of the group.
408 
409  Parameters
410  ----------
411  prefix : str
412  Print each line with this prefix. Used internally for indentation.
413  show_images : bool
414  Show the images belonging to a leaf group.
415  file : file object
416  Where to print the representation. Defaults to sys.stderr
417  """
418  if self.__subgroups is None:
419  print('{}Image List ({})'.format(prefix, len(self.__images)), file=file)
420  if show_images:
421  for im in self.__images:
422  print('{} {}'.format(prefix, im), file=file)
423  else:
424  print('{}Image sub-groups: {}'.format(prefix, ','.join(str(x) for x, _ in self.__subgroups)), file=file)
425  for name, group in self.__subgroups:
426  print('{} {}:'.format(prefix, name), file=file)
427  group.print(prefix + ' ', show_images, file)
428 
429  def __str__(self):
430  """
431  Returns
432  -------
433  str
434  A human-readable representation of the group
435  """
436  string = StringIO()
437  self.print(show_images=True, file=string)
438  return string.getvalue()
439 
440 def load_fits_image(image, psf=None, weight=None, **kwargs):
441  """Creates an image group with the images of a (possibly multi-HDU) single FITS file.
442 
443  If image is multi-hdu, psf and weight can either be multi hdu or lists of individual files.
444 
445  In any case, they are matched in order and HDUs not containing images (two dimensional arrays) are ignored.
446 
447  :param image: The filename of the FITS file containing the image(s)
448  :param psf: psf file or list of psf files
449  :param weight: FITS file for the weight image or a list of such files
450 
451  :return: A ImageGroup representing the images
452  """
453 
454  image_hdu_list = FitsFile(image)
455  image_hdu_idx = image_hdu_list.hdu_list
456 
457  # handles the PSFs
458  if isinstance(psf, list):
459  if len(psf) != len(image_hdu_idx):
460  raise ValueError("The number of psf files must match the number of images!")
461  psf_list = psf
462  psf_hdu_idx = [0] * len(psf_list)
463  else:
464  psf_list = [psf] * len(image_hdu_idx)
465  psf_hdu_idx = range(len(image_hdu_idx))
466 
467  # handles the weight maps
468  if isinstance(weight, list):
469  if len(weight) != len(image_hdu_idx):
470  raise ValueError("The number of weight files must match the number of images!")
471  weight_list = weight
472  weight_hdu_idx = [0] * len(weight_list)
473  elif weight is None:
474  weight_list = [None] * len(image_hdu_idx)
475  weight_hdu_idx = [0] * len(weight_list)
476  else:
477  weight_hdu_list = FitsFile(weight)
478  weight_hdu_idx = weight_hdu_list.hdu_list
479  weight_list = [weight_hdu_list] * len(image_hdu_idx)
480 
481  image_list = []
482  for hdu, psf_file, psf_hdu, weight_file, weight_hdu in zip(
483  image_hdu_idx, psf_list, psf_hdu_idx, weight_list, weight_hdu_idx):
484  image_list.append(MeasurementImage(image_hdu_list, psf_file, weight_file,
485  image_hdu=hdu+1, psf_hdu=psf_hdu+1, weight_hdu=weight_hdu+1, **kwargs))
486 
487  return ImageGroup(images=image_list)
488 
489 def load_fits_images(images, psfs=None, weights=None, **kwargs):
490  """Creates an image group for the given images.
491 
492  Parameters
493  ----------
494  images : list of str
495  A list of relative paths to the images FITS files. Can also be single string in which case,
496  this function acts like load_fits_image
497  psfs : list of str
498  A list of relative paths to the PSF FITS files (optional). It must match the length of image_list or be None.
499  weights : list of str
500  A list of relative paths to the weight files (optional). It must match the length of image_list or be None.
501 
502  Returns
503  -------
504  ImageGroup
505  A ImageGroup representing the images
506 
507  Raises
508  ------
509  ValueError
510  In case of mismatched list of files
511  """
512 
513  if isinstance(images, list):
514  if len(images) == 0:
515  raise ValueError("An empty list passed to load_fits_images")
516 
517  psfs = psfs or [None] * len(images)
518  weights = weights or [None] * len(images)
519 
520  if not isinstance(psfs, list) or len(psfs) != len(images):
521  raise ValueError("The number of image files and psf files must match!")
522 
523  if not isinstance(weights, list) or len(weights) != len(images):
524  raise ValueError("The number of image files and weight files must match!")
525 
526  groups = []
527  for f, p, w in zip(images, psfs, weights):
528  groups.append(load_fits_image(f, p, w, **kwargs))
529 
530  image_list = []
531  for g in groups:
532  image_list += g
533 
534  return ImageGroup(images=image_list)
535  else:
536  load_fits_image(images, psfs, weights, **kwargs)
537 
538 class ByKeyword(object):
539  """
540  Callable that can be used to split an ImageGroup by a keyword value (i.e. FILTER).
541 
542  Parameters
543  ----------
544  key : str
545  FITS header keyword (i.e. FILTER)
546 
547  See Also
548  --------
549  ImageGroup.split
550  """
551 
552  def __init__(self, key):
553  """
554  Constructor.
555  """
556  self.__key = key
557 
558  def __call__(self, images):
559  """
560  Parameters
561  ----------
562  images : list of MeasurementImage
563  List of images to group
564 
565  Returns
566  -------
567  list of tuples of str and list of MeasurementImage
568  i.e. [
569  (R, [frame_r_01.fits, frame_r_02.fits]),
570  (G, [frame_g_01.fits, frame_g_02.fits])
571  ]
572  """
573  result = {}
574  for im in images:
575  if self.__key not in im.meta:
576  raise KeyError('The image {}[{}] does not contain the key {}'.format(
577  im.meta['IMAGE_FILENAME'], im.image_hdu, self.__key
578  ))
579  if im.meta[self.__key] not in result:
580  result[im.meta[self.__key]] = []
581  result[im.meta[self.__key]].append(im)
582  return [(k, result[k]) for k in result]
583 
584 
585 class ByPattern(object):
586  """
587  Callable that can be used to split an ImageGroup by a keyword value (i.e. FILTER), applying a regular
588  expression and using the first matching group as key.
589 
590  Parameters
591  ----------
592  key : str
593  FITS header keyword
594  pattern : str
595  Regular expression. The first matching group will be used as grouping key.
596 
597  See Also
598  --------
599  ImageGroup.split
600  """
601 
602  def __init__(self, key, pattern):
603  """
604  Constructor.
605  """
606  self.__key = key
607  self.__pattern = pattern
608 
609  def __call__(self, images):
610  """
611  Parameters
612  ----------
613  images : list of MeasurementImage
614  List of images to group
615 
616  Returns
617  -------
618  list of tuples of str and list of MeasurementImage
619  """
620  result = {}
621  for im in images:
622  if self.__key not in im.meta:
623  raise KeyError('The image {}[{}] does not contain the key {}'.format(
624  im.meta['IMAGE_FILENAME'], im.image_hdu, self.__key
625  ))
626  group = re.match(self.__pattern, im.meta[self.__key]).group(1)
627  if group not in result:
628  result[group] = []
629  result[group].append(im)
630  return [(k, result[k]) for k in result]
631 
632 
633 class MeasurementGroup(object):
634  """
635  Once an instance of this class is created from an ImageGroup, its configuration is "frozen". i.e.
636  no new images can be added, or no new grouping applied.
637 
638  Parameters
639  ----------
640  image_group : ImageGroup
641  """
642 
643  _all_groups = list()
644 
645  def __init__(self, image_group, is_subgroup=False):
646  """
647  Constructor.
648  """
649  self.__images = None
650  self.__subgroups = None
651  if image_group.is_leaf():
652  self.__images = [im for im in image_group]
653  else:
654  self.__subgroups = [(n, MeasurementGroup(g, is_subgroup=True)) for n,g in image_group]
655  if not is_subgroup:
656  MeasurementGroup._all_groups.append(self)
657 
658  def __iter__(self):
659  """
660  Returns
661  -------
662  iterator
663  """
664  if self.__subgroups:
665  return self.__subgroups.__iter__()
666  else:
667  return self.__images.__iter__()
668 
669  def __getitem__(self, index):
670  """
671  The subgroup with the given name or image with the given index depending on whether this is a leaf group.
672 
673  Parameters
674  ----------
675  index : str or int
676  Subgroup name or image index
677 
678  Returns
679  -------
680  MeasurementGroup or MeasurementImage
681 
682  Raises
683  ------
684  KeyError
685  If we can't find what we want
686  """
687 
688  if self.__subgroups:
689  try:
690  return next(x for x in self.__subgroups if x[0] == index)[1]
691  except StopIteration:
692  raise KeyError('Group {} not found'.format(index))
693  else:
694  try:
695  return self.__images[index]
696  except:
697  raise KeyError('Image #{} not found'.format(index))
698 
699  def __len__(self):
700  """
701  Returns
702  -------
703  int
704  Number of subgroups, or images contained within the group
705  """
706  if self.__subgroups:
707  return len(self.__subgroups)
708  else:
709  return len(self.__images)
710 
711  def is_leaf(self):
712  """
713  Returns
714  -------
715  bool
716  True if the group is a leaf group
717  """
718  return self.__subgroups is None
719 
720  def print(self, prefix='', show_images=False, file=sys.stderr):
721  """
722  Print a human-readable representation of the group.
723 
724  Parameters
725  ----------
726  prefix : str
727  Print each line with this prefix. Used internally for indentation.
728  show_images : bool
729  Show the images belonging to a leaf group.
730  file : file object
731  Where to print the representation. Defaults to sys.stderr
732  """
733  if self.__images:
734  print('{}Image List ({})'.format(prefix, len(self.__images)), file=file)
735  if show_images:
736  for im in self.__images:
737  print('{} {}'.format(prefix, im), file=file)
738  if self.__subgroups:
739  print('{}Measurement sub-groups: {}'.format(prefix, ','.join(
740  x for x, _ in self.__subgroups)), file=file)
741  for name, group in self.__subgroups:
742  print('{} {}:'.format(prefix, name), file=file)
743  group.print(prefix + ' ', show_images, file=file)
744 
745  def __str__(self):
746  """
747  Returns
748  -------
749  str
750  A human-readable representation of the group
751  """
752  string = StringIO()
753  self.print(show_images=True, file=string)
754  return string.getvalue()
ELEMENTS_API auto join(Args &&...args) -> decltype(joinPath(std::forward< Args >(args)...))