Merge pull request #502 from joleroi/subdatastore

DataStore subset

gammapy/data/data_store.py

@@ -4,6 +4,8 @@
 import logging
 import numpy as np
 from collections import OrderedDict
+
+import subprocess
 from astropy.table import Table
 from astropy.utils import lazyproperty
 from astropy.units import Quantity
@@ -241,15 +243,21 @@ def load_many(self, obs_ids, hdu_type=None, hdu_class=None):
 
         return things
 
-    def check_integrity(self, logger):
+    def check_integrity(self, logger=None):
         """Check integrity, i.e. whether index and observation table match.
         """
+        # Todo: This is broken - remove or fix?
+        sane = True
+        if logger is None:
+            logger = logging.getLogger('default')
+
         logger.info('Checking event list files')
         available = self.check_available_event_lists(logger)
         if np.any(~available):
             logger.warning('Number of missing event list files: {}'.format(np.invert(available).sum()))
 
-            # TODO: implement better, more complete integrity checks.
+        # TODO: implement better, more complete integrity checks.
+        return sane
 
     def make_table_of_files(self, observation_table=None, filetypes=['events']):
         """Make list of files in the datastore directory.
@@ -267,7 +275,7 @@ def make_table_of_files(self, observation_table=None, filetypes=['events']):
             Table summarising info about files.
         """
         if observation_table is None:
-            observation_table = ObservationTable(self.index_table)
+            observation_table = ObservationTable(self.obs_table)
 
         data = []
         for observation in observation_table:
@@ -291,7 +299,9 @@ def check_available_event_lists(self, logger=None):
         file_available : `~numpy.ndarray`
             Boolean mask which files are available.
         """
-        observation_table = self.index_table
+        # Todo: This is broken. Remove (covered by HDUlocation class)?
+
+        observation_table = self.obs_table
         file_available = np.ones(len(observation_table), dtype='bool')
         for ii in range(len(observation_table)):
             obs_id = observation_table['OBS_ID'][ii]
@@ -304,6 +314,35 @@ def check_available_event_lists(self, logger=None):
 
         return file_available
 
+    def copy_obs(self, obs_table, outdir):
+        """Create a new `~gammapy.data.DataStore` containing a subset of observations
+
+        Parameters
+        ----------
+        obs_table : `~gammapy.data.ObservationTable`
+            Table of observation to create the subset
+        outdir : str, Path
+            Directory for the new store
+        """
+        outdir = make_path(outdir)
+        obs_ids = obs_table['OBS_ID'].data
+
+        hdutable = self.hdu_table
+        hdutable.add_index('OBS_ID')
+        subhdutable = hdutable.loc[obs_ids]
+        subobstable = self.obs_table.select_obs_id(obs_ids)
+
+        for ii in range(len(subhdutable)):
+            # Changes to the file structure could be made here
+            loc = subhdutable._location_info(ii)
+            targetdir = outdir / loc.file_dir
+            targetdir.mkdir(exist_ok=True, parents=True)
+            cmd = ['cp', str(loc.path()), str(targetdir)]
+            subprocess.call(cmd)
+
+        subhdutable.write(str(outdir/self.DEFAULT_HDU_TABLE), format='fits')
+        subobstable.write(str(outdir/self.DEFAULT_OBS_TABLE), format='fits')
+
 
 class DataStoreObservation(object):
     """IACT data store observation.

gammapy/data/tests/test_data_store.py

@@ -68,3 +68,26 @@ def test_datastore_load_all(data_manager):
     event_lists = data_store.load_all(hdu_class='events')
     assert_allclose(event_lists[0]['ENERGY'][0], 1.1156039)
     assert_allclose(event_lists[-1]['ENERGY'][0], 1.0204216)
+
+
+@requires_data('gammapy-extra')
+@requires_dependency('yaml')
+def test_datastore_subset(tmpdir, data_manager):
+    """Test creating a datastore as subset of another datastore"""
+    data_store = data_manager['hess-crab4-hd-hap-prod2']
+    selected_obs = data_store.obs_table.select_obs_id([23523, 23592])
+    storedir = tmpdir / 'substore'
+    data_store.copy_obs(selected_obs, storedir)
+
+    substore = DataStore.from_dir(storedir)
+
+    # Todo : This is broken, remove or fix?
+    # substore.check_integrity()
+
+    assert str(substore.hdu_table.base_dir) == str(storedir)
+    assert len(substore.obs_table) == 2
+
+    desired = data_store.obs(23523)
+    actual = substore.obs(23523)
+
+    assert str(actual.events) == str(desired.events)

gammapy/image/maps.py

@@ -342,13 +342,10 @@ def reproject(self, refheader, mode='interp', *args, **kwargs):
             Skymap reprojected onto ``refheader``.
         """
         from reproject import reproject_interp
-        out = reproject_interp(
-            input_data=self.to_image_hdu(),
-            output_projection=refheader,
-            *args,
-            **kwargs,
-        )
-        map = SkyMap(name=self.name, data=out[0], wcs=self.wcs, unit=self.unit, meta=self.meta)
+        out = reproject_interp(input_data=self.to_image_hdu(),
+                               output_projection=refheader, *args, **kwargs)
+        map = SkyMap(name=self.name, data=out[0], wcs=self.wcs, unit=self.unit,
+                     meta=self.meta)
         return map
 
     def show(self, viewer='mpl', **kwargs):