1
|
|
|
import collections |
|
|
|
|
2
|
|
|
|
3
|
|
|
from hawc_hal.serialize import Serialization |
4
|
|
|
from hawc_hal.region_of_interest import get_roi_from_dict |
5
|
|
|
|
6
|
|
|
from threeML.exceptions.custom_exceptions import custom_warnings |
|
|
|
|
7
|
|
|
|
8
|
|
|
from ..healpix_handling import SparseHealpix, DenseHealpix |
9
|
|
|
from data_analysis_bin import DataAnalysisBin |
|
|
|
|
10
|
|
|
|
11
|
|
|
|
12
|
|
|
def from_hdf5_file(map_tree_file, roi): |
|
|
|
|
13
|
|
|
""" |
14
|
|
|
Create a MapTree object from a HDF5 file and a ROI. Do not use this directly, use map_tree_factory instead. |
15
|
|
|
|
16
|
|
|
:param map_tree_file: |
17
|
|
|
:param roi: |
18
|
|
|
:return: |
19
|
|
|
""" |
20
|
|
|
|
21
|
|
|
# Read the data frames contained in the file |
22
|
|
|
with Serialization(map_tree_file) as serializer: |
23
|
|
|
|
24
|
|
|
analysis_bins_df, _ = serializer.retrieve_pandas_object('/analysis_bins') |
25
|
|
|
meta_df, _ = serializer.retrieve_pandas_object('/analysis_bins_meta') |
26
|
|
|
_, roi_meta = serializer.retrieve_pandas_object('/ROI') |
27
|
|
|
|
28
|
|
|
# Let's see if the file contains the definition of an ROI |
29
|
|
|
if len(roi_meta) > 0: |
|
|
|
|
30
|
|
|
|
31
|
|
|
# Yes. Let's build it |
32
|
|
|
file_roi = get_roi_from_dict(roi_meta) |
33
|
|
|
|
34
|
|
|
# Now let's check that the ROI the user has provided (if any) is compatible with the one contained |
35
|
|
|
# in the file (i.e., either they are the same, or the user-provided one is smaller) |
36
|
|
|
if roi is not None: |
37
|
|
|
|
38
|
|
|
# Let's test with a nside=1024 (the highest we will use in practice) |
39
|
|
|
active_pixels_file = file_roi.active_pixels(1024) |
40
|
|
|
active_pixels_user = roi.active_pixels(1024) |
41
|
|
|
|
42
|
|
|
# This verifies that active_pixels_file is a superset (or equal) to the user-provided set |
43
|
|
|
assert set(active_pixels_file) >= set(active_pixels_user), \ |
44
|
|
|
"The ROI you provided (%s) is not a subset " \ |
45
|
|
|
"of the one contained in the file (%s)" % (roi, file_roi) |
46
|
|
|
|
47
|
|
|
else: |
48
|
|
|
|
49
|
|
|
# The user has provided no ROI, but the file contains one. Let's issue a warning |
50
|
|
|
custom_warnings.warn("You did not provide any ROI but the map tree %s contains " |
51
|
|
|
"only data within the ROI %s. " |
52
|
|
|
"Only those will be used." % (map_tree_file, file_roi)) |
53
|
|
|
|
54
|
|
|
# Make a copy of the file ROI and use it as if the user provided that one |
55
|
|
|
roi = get_roi_from_dict(file_roi.to_dict()) |
56
|
|
|
|
57
|
|
|
# Get the name of the analysis bins |
58
|
|
|
|
59
|
|
|
bin_names = analysis_bins_df.index.levels[0] |
60
|
|
|
|
61
|
|
|
# Loop over them and build a DataAnalysisBin instance for each one |
62
|
|
|
|
63
|
|
|
data_analysis_bins = collections.OrderedDict() |
64
|
|
|
|
65
|
|
|
for bin_name in bin_names: |
66
|
|
|
|
67
|
|
|
this_df = analysis_bins_df.loc[bin_name] |
68
|
|
|
this_meta = meta_df.loc[bin_name] |
69
|
|
|
|
70
|
|
|
if roi is not None: |
71
|
|
|
|
72
|
|
|
# Get the active pixels for this plane |
73
|
|
|
active_pixels_user = roi.active_pixels(this_meta['nside']) |
74
|
|
|
|
75
|
|
|
# Read only the pixels that the user wants |
76
|
|
|
observation_hpx_map = SparseHealpix(this_df.loc[active_pixels_user, 'observation'].values, |
77
|
|
|
active_pixels_user, this_meta['nside']) |
78
|
|
|
background_hpx_map = SparseHealpix(this_df.loc[active_pixels_user, 'background'].values, |
79
|
|
|
active_pixels_user, this_meta['nside']) |
80
|
|
|
|
81
|
|
|
else: |
82
|
|
|
|
83
|
|
|
# Full sky |
84
|
|
|
observation_hpx_map = DenseHealpix(this_df.loc[:, 'observation'].values) |
85
|
|
|
background_hpx_map = DenseHealpix(this_df.loc[:, 'background'].values) |
86
|
|
|
|
87
|
|
|
# This signals the DataAnalysisBin that we are dealing with a full sky map |
88
|
|
|
active_pixels_user = None |
89
|
|
|
|
90
|
|
|
# Let's now build the instance |
91
|
|
|
this_bin = DataAnalysisBin(bin_name, |
92
|
|
|
observation_hpx_map=observation_hpx_map, |
93
|
|
|
background_hpx_map=background_hpx_map, |
94
|
|
|
active_pixels_ids=active_pixels_user, |
95
|
|
|
n_transits=this_meta['n_transits'], |
96
|
|
|
scheme='RING' if this_meta['scheme'] == 0 else 'NEST') |
97
|
|
|
|
98
|
|
|
data_analysis_bins[bin_name] = this_bin |
99
|
|
|
|
100
|
|
|
return data_analysis_bins |
101
|
|
|
|
The coding style of this project requires that you add a docstring to this code element. Below, you find an example for methods:
If you would like to know more about docstrings, we recommend to read PEP-257: Docstring Conventions.