1
|
|
|
# Copyright 2014 Diamond Light Source Ltd. |
2
|
|
|
# |
3
|
|
|
# Licensed under the Apache License, Version 2.0 (the "License"); |
4
|
|
|
# you may not use this file except in compliance with the License. |
5
|
|
|
# You may obtain a copy of the License at |
6
|
|
|
# |
7
|
|
|
# http://www.apache.org/licenses/LICENSE-2.0 |
8
|
|
|
# |
9
|
|
|
# Unless required by applicable law or agreed to in writing, software |
10
|
|
|
# distributed under the License is distributed on an "AS IS" BASIS, |
11
|
|
|
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. |
12
|
|
|
# See the License for the specific language governing permissions and |
13
|
|
|
# limitations under the License. |
14
|
|
|
|
15
|
|
|
""" |
16
|
|
|
.. module:: nxtomo_loader |
17
|
|
|
:platform: Unix |
18
|
|
|
:synopsis: A class for loading standard tomography data in Nexus format. |
19
|
|
|
|
20
|
|
|
.. moduleauthor:: Nicola Wadeson <[email protected]> |
21
|
|
|
|
22
|
|
|
""" |
23
|
|
|
|
24
|
|
|
import h5py |
25
|
|
|
import logging |
26
|
|
|
import numpy as np |
27
|
|
|
|
28
|
|
|
from savu.plugins.loaders.base_loader import BaseLoader |
29
|
|
|
from savu.plugins.utils import register_plugin |
30
|
|
|
|
31
|
|
|
from savu.data.data_structures.data_types.data_plus_darks_and_flats \ |
32
|
|
|
import ImageKey, NoImageKey |
33
|
|
|
|
34
|
|
|
|
35
|
|
|
@register_plugin |
36
|
|
|
class NxtomoLoader(BaseLoader): |
37
|
|
|
""" |
38
|
|
|
""" |
39
|
|
|
def __init__(self, name='NxtomoLoader'): |
40
|
|
|
super(NxtomoLoader, self).__init__(name) |
41
|
|
|
self.warnings = [] |
42
|
|
|
|
43
|
|
|
def log_warning(self, msg): |
44
|
|
|
logging.warning(msg) |
45
|
|
|
self.warnings.append(msg) |
46
|
|
|
|
47
|
|
|
def setup(self): |
48
|
|
|
exp = self.get_experiment() |
49
|
|
|
data_obj = exp.create_data_object('in_data', self.parameters['name']) |
50
|
|
|
|
51
|
|
|
data_obj.backing_file = self._get_data_file() |
52
|
|
|
|
53
|
|
|
data_obj.data = self._get_h5_entry( |
54
|
|
|
data_obj.backing_file, self.parameters['data_path']) |
55
|
|
|
|
56
|
|
|
synthetic_path = f"{'/'.join(self.parameters['data_path'].split('/')[:-1])}/synthetic/synthetic" |
57
|
|
|
if synthetic_path in data_obj.backing_file: |
58
|
|
|
if data_obj.backing_file[synthetic_path][()] == True: |
59
|
|
|
self.exp.meta_data.set("synthetic", True) |
60
|
|
|
|
61
|
|
|
self._set_dark_and_flat(data_obj) |
62
|
|
|
|
63
|
|
|
self.nFrames = self.__get_nFrames(data_obj) |
64
|
|
|
if self.nFrames > 1: |
65
|
|
|
self.__setup_4d(data_obj) |
66
|
|
|
self.__setup_3d_to_4d(data_obj, self.nFrames) |
67
|
|
|
else: |
68
|
|
|
if len(data_obj.data.shape) == 3: |
69
|
|
|
self._setup_3d(data_obj) |
70
|
|
|
else: |
71
|
|
|
self.__setup_4d(data_obj) |
72
|
|
|
data_obj.set_original_shape(data_obj.data.shape) |
73
|
|
|
self._set_rotation_angles(data_obj) |
74
|
|
|
self._set_projection_shifts(data_obj) |
75
|
|
|
|
76
|
|
|
try: |
77
|
|
|
control = self._get_h5_path( |
78
|
|
|
data_obj.backing_file, 'entry1/tomo_entry/control/data') |
79
|
|
|
data_obj.meta_data.set("control", control[...]) |
80
|
|
|
except Exception: |
81
|
|
|
self.log_warning("No Control information available") |
82
|
|
|
|
83
|
|
|
nAngles = len(data_obj.meta_data.get('rotation_angle')) |
84
|
|
|
self.__check_angles(data_obj, nAngles) |
85
|
|
|
|
86
|
|
|
self.set_data_reduction_params(data_obj) |
87
|
|
|
data_obj.data._set_dark_and_flat() |
88
|
|
|
|
89
|
|
|
def _get_h5_entry(self, filename, path): |
90
|
|
|
if path in filename: |
91
|
|
|
return filename[path] |
92
|
|
|
else: |
93
|
|
|
raise Exception("Path %s not found in %s" % (path, filename)) |
94
|
|
|
|
95
|
|
|
def __get_nFrames(self, dObj): |
96
|
|
|
if self.parameters['3d_to_4d'] is False: |
97
|
|
|
return 0 |
98
|
|
|
if self.parameters['3d_to_4d'] is True: |
99
|
|
|
try: |
100
|
|
|
# for backwards compatibility |
101
|
|
|
n_frames = eval(self.parameters["angles"], {"builtins": None, "np": np}) |
102
|
|
|
return np.array(n_frames).shape[0] |
103
|
|
|
except Exception: |
104
|
|
|
raise Exception("Please specify the angles, or the number of " |
105
|
|
|
"frames per scan (via 3d_to_4d param) in the loader.") |
106
|
|
|
if isinstance(self.parameters['3d_to_4d'], int): |
107
|
|
|
return self.parameters['3d_to_4d'] |
108
|
|
|
else: |
109
|
|
|
raise Exception("Unknown value for loader parameter '3d_to_4d', " |
110
|
|
|
"please specify an integer value") |
111
|
|
|
|
112
|
|
|
def _setup_3d(self, data_obj): |
113
|
|
|
logging.debug("Setting up 3d tomography data.") |
114
|
|
|
rot = 0 |
115
|
|
|
detY = 1 |
116
|
|
|
detX = 2 |
117
|
|
|
data_obj.set_axis_labels('rotation_angle.degrees', |
118
|
|
|
'detector_y.pixel', |
119
|
|
|
'detector_x.pixel') |
120
|
|
|
|
121
|
|
|
data_obj.add_pattern('PROJECTION', core_dims=(detX, detY), |
122
|
|
|
slice_dims=(rot,)) |
123
|
|
|
data_obj.add_pattern('SINOGRAM', core_dims=(detX, rot), |
124
|
|
|
slice_dims=(detY,)) |
125
|
|
|
data_obj.add_pattern('TANGENTOGRAM', core_dims=(rot, detY), |
126
|
|
|
slice_dims=(detX,)) |
127
|
|
|
|
128
|
|
|
def __setup_3d_to_4d(self, data_obj, n_frames): |
129
|
|
|
logging.debug("setting up 4d tomography data from 3d input.") |
130
|
|
|
from savu.data.data_structures.data_types.map_3dto4d_h5 \ |
131
|
|
|
import Map3dto4dh5 |
132
|
|
|
|
133
|
|
|
all_angles = data_obj.data.shape[0] |
134
|
|
|
if all_angles % n_frames != 0: |
135
|
|
|
self.log_warning("There are not a complete set of scans in this file, " |
136
|
|
|
"loading complete ones only") |
137
|
|
|
data_obj.data = Map3dto4dh5(data_obj, n_frames) |
138
|
|
|
data_obj.set_original_shape(data_obj.data.get_shape()) |
139
|
|
|
|
140
|
|
|
def __setup_4d(self, data_obj): |
141
|
|
|
logging.debug("setting up 4d tomography data.") |
142
|
|
|
rot = 0 |
143
|
|
|
detY = 1 |
144
|
|
|
detX = 2 |
145
|
|
|
scan = 3 |
146
|
|
|
|
147
|
|
|
data_obj.set_axis_labels('rotation_angle.degrees', 'detector_y.pixel', |
148
|
|
|
'detector_x.pixel', 'scan.number') |
149
|
|
|
|
150
|
|
|
data_obj.add_pattern('PROJECTION', core_dims=(detX, detY), |
151
|
|
|
slice_dims=(rot, scan)) |
152
|
|
|
data_obj.add_pattern('SINOGRAM', core_dims=(detX, rot), |
153
|
|
|
slice_dims=(detY, scan)) |
154
|
|
|
data_obj.add_pattern('TANGENTOGRAM', core_dims=(rot, detY), |
155
|
|
|
slice_dims=(detX, scan)) |
156
|
|
|
data_obj.add_pattern('SINOMOVIE', core_dims=(detX, rot, scan), |
157
|
|
|
slice_dims=(detY,)) |
158
|
|
|
|
159
|
|
|
def _set_dark_and_flat(self, data_obj): |
160
|
|
|
flat = self.parameters['flat'][0] |
161
|
|
|
dark = self.parameters['dark'][0] |
162
|
|
|
|
163
|
|
|
if not flat or not dark: |
164
|
|
|
self.__find_dark_and_flat(data_obj, flat=flat, dark=dark) |
165
|
|
|
else: |
166
|
|
|
self.__set_separate_dark_and_flat(data_obj) |
167
|
|
|
|
168
|
|
|
def __find_dark_and_flat(self, data_obj, flat=None, dark=None): |
169
|
|
|
ignore = self.parameters['ignore_flats'] if \ |
170
|
|
|
self.parameters['ignore_flats'] else None |
171
|
|
|
if self.parameters['image_key_path'] is None: |
172
|
|
|
image_key_path = 'dummypath/' |
173
|
|
|
else: |
174
|
|
|
image_key_path = self.parameters['image_key_path'] |
175
|
|
|
try: |
176
|
|
|
image_key = data_obj.backing_file[image_key_path][...] |
177
|
|
|
data_obj.data = \ |
178
|
|
|
ImageKey(data_obj, image_key, 0, ignore=ignore) |
179
|
|
|
except KeyError as Argument: |
180
|
|
|
self.log_warning("An image key was not found due to following error:"+str(Argument)) |
181
|
|
|
try: |
182
|
|
|
data_obj.data = NoImageKey(data_obj, None, 0) |
183
|
|
|
entry = 'entry1/tomo_entry/instrument/detector/' |
184
|
|
|
data_obj.data._set_flat_path(entry + 'flatfield') |
185
|
|
|
data_obj.data._set_dark_path(entry + 'darkfield') |
186
|
|
|
except KeyError: |
187
|
|
|
self.log_warning("Dark/flat data was not found in input file.") |
188
|
|
|
data_obj.data._set_dark_and_flat() |
189
|
|
|
if dark: |
190
|
|
|
data_obj.data.update_dark(dark) |
191
|
|
|
if flat: |
192
|
|
|
data_obj.data.update_flat(flat) |
193
|
|
|
|
194
|
|
|
def __set_separate_dark_and_flat(self, data_obj): |
195
|
|
|
try: |
196
|
|
|
image_key = data_obj.backing_file[ |
197
|
|
|
'entry1/tomo_entry/instrument/detector/image_key'][...] |
198
|
|
|
except: |
199
|
|
|
image_key = None |
200
|
|
|
data_obj.data = NoImageKey(data_obj, image_key, 0) |
201
|
|
|
self.__set_data(data_obj, 'flat', data_obj.data._set_flat_path) |
202
|
|
|
self.__set_data(data_obj, 'dark', data_obj.data._set_dark_path) |
203
|
|
|
|
204
|
|
|
def __set_data(self, data_obj, name, func): |
205
|
|
|
path, entry, scale = self.parameters[name] |
206
|
|
|
|
207
|
|
|
if path.split('/')[0] == 'Savu': |
208
|
|
|
import os |
209
|
|
|
savu_base_path = os.path.join(os.path.dirname( |
210
|
|
|
os.path.realpath(__file__)), '..', '..', '..', '..') |
211
|
|
|
path = os.path.join(savu_base_path, path.split('Savu')[1][1:]) |
212
|
|
|
|
213
|
|
|
ffile = h5py.File(path, 'r') |
214
|
|
|
try: |
215
|
|
|
image_key = \ |
216
|
|
|
ffile['entry1/tomo_entry/instrument/detector/image_key'][...] |
217
|
|
|
func(ffile[entry], imagekey=image_key) |
218
|
|
|
except: |
219
|
|
|
func(ffile[entry]) |
220
|
|
|
|
221
|
|
|
data_obj.data._set_scale(name, scale) |
222
|
|
|
|
223
|
|
|
def _set_rotation_angles(self, data_obj): |
224
|
|
|
angles = self.parameters['angles'] |
225
|
|
|
warn_ms = "No angles found so evenly distributing them between 0 and" \ |
226
|
|
|
" 180 degrees" |
227
|
|
|
if angles is None: |
228
|
|
|
angle_key = 'entry1/tomo_entry/data/rotation_angle' |
229
|
|
|
nxs_angles = self.__get_angles_from_nxs_file(data_obj, angle_key) |
230
|
|
|
if nxs_angles is None: |
231
|
|
|
self.log_warning(warn_ms) |
232
|
|
|
angles = np.linspace(0, 180, data_obj.get_shape()[0]) |
233
|
|
|
else: |
234
|
|
|
angles = nxs_angles |
235
|
|
|
else: |
236
|
|
|
try: |
237
|
|
|
angles = eval(angles) |
238
|
|
|
except Exception as e: |
239
|
|
|
logging.warning(e) |
240
|
|
|
try: |
241
|
|
|
angles = np.loadtxt(angles) |
242
|
|
|
except Exception as e: |
243
|
|
|
logging.debug(e) |
244
|
|
|
self.log_warning(warn_ms) |
245
|
|
|
angles = np.linspace(0, 180, data_obj.get_shape()[0]) |
246
|
|
|
data_obj.meta_data.set("rotation_angle", angles) |
247
|
|
|
return len(angles) |
248
|
|
|
|
249
|
|
|
def _set_projection_shifts(self, data_obj): |
250
|
|
|
proj_shifts = np.zeros((data_obj.get_shape()[0], 2)) # initialise a 2d array of projection shifts |
251
|
|
|
self.exp.meta_data.set('projection_shifts', proj_shifts) |
252
|
|
|
data_obj.meta_data.set("projection_shifts", proj_shifts) |
253
|
|
|
return len(proj_shifts) |
254
|
|
|
|
255
|
|
|
def __get_angles_from_nxs_file(self, data_obj, path): |
256
|
|
|
if path in data_obj.backing_file: |
257
|
|
|
idx = data_obj.data.get_image_key() == 0 if \ |
258
|
|
|
isinstance(data_obj.data, ImageKey) else slice(None) |
259
|
|
|
return data_obj.backing_file[path][idx] |
260
|
|
|
else: |
261
|
|
|
self.log_warning("No rotation angle entry found in input file.") |
262
|
|
|
return None |
263
|
|
|
|
264
|
|
|
def _get_data_file(self): |
265
|
|
|
data = self.exp.meta_data.get("data_file") |
266
|
|
|
return h5py.File(data, 'r') |
267
|
|
|
|
268
|
|
|
def __check_angles(self, data_obj, n_angles): |
269
|
|
|
rot_dim = data_obj.get_data_dimension_by_axis_label("rotation_angle") |
270
|
|
|
data_angles = data_obj.data.get_shape()[rot_dim] |
271
|
|
|
if data_angles != n_angles: |
272
|
|
|
if self.nFrames > 1: |
273
|
|
|
rot_angles = data_obj.meta_data.get("rotation_angle") |
274
|
|
|
try: |
275
|
|
|
full_rotations = n_angles // data_angles |
276
|
|
|
cleaned_size = full_rotations * data_angles |
277
|
|
|
if cleaned_size != n_angles: |
278
|
|
|
rot_angles = rot_angles[0:cleaned_size] |
279
|
|
|
self.log_warning( |
280
|
|
|
"the angle list has more values than expected in it") |
281
|
|
|
rot_angles = np.reshape( |
282
|
|
|
rot_angles, [full_rotations, data_angles]) |
283
|
|
|
data_obj.meta_data.set("rotation_angle", |
284
|
|
|
np.transpose(rot_angles)) |
285
|
|
|
return |
286
|
|
|
except: |
287
|
|
|
pass |
288
|
|
|
raise Exception("The number of angles %s does not match the data " |
289
|
|
|
"dimension length %s" % (n_angles, data_angles)) |
290
|
|
|
|
291
|
|
|
def executive_summary(self): |
292
|
|
|
""" Provide a summary to the user for the result of the plugin. |
293
|
|
|
|
294
|
|
|
e.g. |
295
|
|
|
- Warning, the sample may have shifted during data collection |
296
|
|
|
- Filter operated normally |
297
|
|
|
|
298
|
|
|
:returns: A list of string summaries |
299
|
|
|
""" |
300
|
|
|
if len(self.warnings) == 0: |
301
|
|
|
return ["Nothing to Report"] |
302
|
|
|
else: |
303
|
|
|
return self.warnings |
304
|
|
|
|