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
|
|
|
from savu.data.stats.statistics import Statistics |
31
|
|
|
|
32
|
|
|
from savu.data.data_structures.data_types.data_plus_darks_and_flats \ |
33
|
|
|
import ImageKey, NoImageKey |
34
|
|
|
|
35
|
|
|
|
36
|
|
|
@register_plugin |
37
|
|
|
class NxtomoLoader(BaseLoader): |
38
|
|
|
""" |
39
|
|
|
""" |
40
|
|
|
def __init__(self, name='NxtomoLoader'): |
41
|
|
|
super(NxtomoLoader, self).__init__(name) |
42
|
|
|
self.warnings = [] |
43
|
|
|
|
44
|
|
|
def log_warning(self, msg): |
45
|
|
|
logging.warning(msg) |
46
|
|
|
self.warnings.append(msg) |
47
|
|
|
|
48
|
|
|
def setup(self): |
49
|
|
|
exp = self.get_experiment() |
50
|
|
|
data_obj = exp.create_data_object('in_data', self.parameters['name']) |
51
|
|
|
|
52
|
|
|
data_obj.backing_file = self._get_data_file() |
53
|
|
|
|
54
|
|
|
data_obj.data = self._get_h5_entry( |
55
|
|
|
data_obj.backing_file, self.parameters['data_path']) |
56
|
|
|
exp.meta_data.set('data_path', self.parameters['data_path']) |
57
|
|
|
exp.meta_data.set('image_key_path', self.parameters['image_key_path']) |
58
|
|
|
|
59
|
|
|
self._setup_after_pre_run(data_obj) |
60
|
|
|
|
61
|
|
|
synthetic_path = f"{'/'.join(self.parameters['data_path'].split('/')[:-1])}/synthetic/synthetic" |
62
|
|
|
if synthetic_path in data_obj.backing_file: |
63
|
|
|
if data_obj.backing_file[synthetic_path][()] == True: |
64
|
|
|
self.exp.meta_data.set("synthetic", True) |
65
|
|
|
|
66
|
|
|
self._set_dark_and_flat(data_obj) |
67
|
|
|
|
68
|
|
|
self.nFrames = self.__get_nFrames(data_obj) |
69
|
|
|
if self.nFrames > 1: |
70
|
|
|
self.__setup_4d(data_obj) |
71
|
|
|
self.__setup_3d_to_4d(data_obj, self.nFrames) |
72
|
|
|
else: |
73
|
|
|
if len(data_obj.data.shape) == 3: |
74
|
|
|
self._setup_3d(data_obj) |
75
|
|
|
else: |
76
|
|
|
self.__setup_4d(data_obj) |
77
|
|
|
data_obj.set_original_shape(data_obj.data.shape) |
78
|
|
|
self._set_rotation_angles(data_obj) |
79
|
|
|
self._set_projection_shifts(data_obj) |
80
|
|
|
|
81
|
|
|
try: |
82
|
|
|
control = self._get_h5_path( |
83
|
|
|
data_obj.backing_file, 'entry1/tomo_entry/control/data') |
84
|
|
|
data_obj.meta_data.set("control", control[...]) |
85
|
|
|
except Exception: |
86
|
|
|
self.log_warning("No Control information available") |
87
|
|
|
|
88
|
|
|
nAngles = len(data_obj.meta_data.get('rotation_angle')) |
89
|
|
|
self.__check_angles(data_obj, nAngles) |
90
|
|
|
|
91
|
|
|
self.set_data_reduction_params(data_obj) |
92
|
|
|
data_obj.data._set_dark_and_flat() |
93
|
|
|
|
94
|
|
|
def _setup_after_pre_run(self, data_obj): |
95
|
|
|
|
96
|
|
|
fsplit = self.exp.meta_data["data_path"].split("/") |
97
|
|
|
fsplit[-1] = "stats" |
98
|
|
|
stats_path = "/".join(fsplit) |
99
|
|
|
if stats_path in data_obj.backing_file: |
100
|
|
|
print("STATS FOUND IN INPUT DATA") |
101
|
|
|
stats_obj = Statistics() |
102
|
|
|
stats_obj.p_num = 1 |
103
|
|
|
stats_obj.pattern = "PROJECTION" |
104
|
|
|
stats = data_obj.backing_file[stats_path] |
105
|
|
|
stats_obj.stats_key = list(stats.attrs.get("stats_key")) |
106
|
|
|
stats_dict = stats_obj._array_to_dict(stats) |
107
|
|
|
Statistics.global_stats[1] = [stats_dict] |
108
|
|
|
Statistics.plugin_names[1] = "raw_data" |
109
|
|
|
for key in list(stats_dict.keys()): |
110
|
|
|
data_obj.meta_data.set(["stats", key], stats_dict[key]) |
111
|
|
|
stats_obj._write_stats_to_file(p_num=1, plugin_name="raw_data") |
112
|
|
|
|
113
|
|
|
fsplit = self.exp.meta_data["data_path"].split("/") |
114
|
|
|
fsplit[-1] = "preview" |
115
|
|
|
preview_path = "/".join(fsplit) |
116
|
|
|
if preview_path in data_obj.backing_file: |
117
|
|
|
print("PREVIEW FOUND IN INPUT DATA") |
118
|
|
|
preview_str = data_obj.backing_file[preview_path][()] |
119
|
|
|
preview = preview_str.split(",") |
120
|
|
|
for i, crop in enumerate(self.parameters["preview"]): |
121
|
|
|
if crop != ":": # Take preview dimensions from user parameter where they exist |
122
|
|
|
preview[i] = crop |
123
|
|
|
self.parameters["preview"] = preview |
124
|
|
|
|
125
|
|
|
def _get_h5_entry(self, filename, path): |
126
|
|
|
if path in filename: |
127
|
|
|
return filename[path] |
128
|
|
|
else: |
129
|
|
|
raise Exception("Path %s not found in %s" % (path, filename)) |
130
|
|
|
|
131
|
|
|
def __get_nFrames(self, dObj): |
132
|
|
|
if self.parameters['3d_to_4d'] is False: |
133
|
|
|
return 0 |
134
|
|
|
if self.parameters['3d_to_4d'] is True: |
135
|
|
|
try: |
136
|
|
|
# for backwards compatibility |
137
|
|
|
n_frames = eval(self.parameters["angles"], {"builtins": None, "np": np}) |
138
|
|
|
return np.array(n_frames).shape[0] |
139
|
|
|
except Exception: |
140
|
|
|
raise Exception("Please specify the angles, or the number of " |
141
|
|
|
"frames per scan (via 3d_to_4d param) in the loader.") |
142
|
|
|
if isinstance(self.parameters['3d_to_4d'], int): |
143
|
|
|
return self.parameters['3d_to_4d'] |
144
|
|
|
else: |
145
|
|
|
raise Exception("Unknown value for loader parameter '3d_to_4d', " |
146
|
|
|
"please specify an integer value") |
147
|
|
|
|
148
|
|
|
def _setup_3d(self, data_obj): |
149
|
|
|
logging.debug("Setting up 3d tomography data.") |
150
|
|
|
rot = 0 |
151
|
|
|
detY = 1 |
152
|
|
|
detX = 2 |
153
|
|
|
data_obj.set_axis_labels('rotation_angle.degrees', |
154
|
|
|
'detector_y.pixel', |
155
|
|
|
'detector_x.pixel') |
156
|
|
|
|
157
|
|
|
data_obj.add_pattern('PROJECTION', core_dims=(detX, detY), |
158
|
|
|
slice_dims=(rot,)) |
159
|
|
|
data_obj.add_pattern('SINOGRAM', core_dims=(detX, rot), |
160
|
|
|
slice_dims=(detY,)) |
161
|
|
|
data_obj.add_pattern('TANGENTOGRAM', core_dims=(rot, detY), |
162
|
|
|
slice_dims=(detX,)) |
163
|
|
|
|
164
|
|
|
def __setup_3d_to_4d(self, data_obj, n_frames): |
165
|
|
|
logging.debug("setting up 4d tomography data from 3d input.") |
166
|
|
|
from savu.data.data_structures.data_types.map_3dto4d_h5 \ |
167
|
|
|
import Map3dto4dh5 |
168
|
|
|
|
169
|
|
|
all_angles = data_obj.data.shape[0] |
170
|
|
|
if all_angles % n_frames != 0: |
171
|
|
|
self.log_warning("There are not a complete set of scans in this file, " |
172
|
|
|
"loading complete ones only") |
173
|
|
|
data_obj.data = Map3dto4dh5(data_obj, n_frames) |
174
|
|
|
data_obj.set_original_shape(data_obj.data.get_shape()) |
175
|
|
|
|
176
|
|
|
def __setup_4d(self, data_obj): |
177
|
|
|
logging.debug("setting up 4d tomography data.") |
178
|
|
|
rot = 0 |
179
|
|
|
detY = 1 |
180
|
|
|
detX = 2 |
181
|
|
|
scan = 3 |
182
|
|
|
|
183
|
|
|
data_obj.set_axis_labels('rotation_angle.degrees', 'detector_y.pixel', |
184
|
|
|
'detector_x.pixel', 'scan.number') |
185
|
|
|
|
186
|
|
|
data_obj.add_pattern('PROJECTION', core_dims=(detX, detY), |
187
|
|
|
slice_dims=(rot, scan)) |
188
|
|
|
data_obj.add_pattern('SINOGRAM', core_dims=(detX, rot), |
189
|
|
|
slice_dims=(detY, scan)) |
190
|
|
|
data_obj.add_pattern('TANGENTOGRAM', core_dims=(rot, detY), |
191
|
|
|
slice_dims=(detX, scan)) |
192
|
|
|
data_obj.add_pattern('SINOMOVIE', core_dims=(detX, rot, scan), |
193
|
|
|
slice_dims=(detY,)) |
194
|
|
|
|
195
|
|
|
def _set_dark_and_flat(self, data_obj): |
196
|
|
|
flat = self.parameters['flat'][0] |
197
|
|
|
dark = self.parameters['dark'][0] |
198
|
|
|
|
199
|
|
|
if not flat or not dark: |
200
|
|
|
self.__find_dark_and_flat(data_obj, flat=flat, dark=dark) |
201
|
|
|
else: |
202
|
|
|
self.__set_separate_dark_and_flat(data_obj) |
203
|
|
|
|
204
|
|
|
def __find_dark_and_flat(self, data_obj, flat=None, dark=None): |
205
|
|
|
ignore = self.parameters['ignore_flats'] if \ |
206
|
|
|
self.parameters['ignore_flats'] else None |
207
|
|
|
if self.parameters['image_key_path'] is None: |
208
|
|
|
image_key_path = 'dummypath/' |
209
|
|
|
else: |
210
|
|
|
image_key_path = self.parameters['image_key_path'] |
211
|
|
|
try: |
212
|
|
|
image_key = data_obj.backing_file[image_key_path][...] |
213
|
|
|
data_obj.data = \ |
214
|
|
|
ImageKey(data_obj, image_key, 0, ignore=ignore) |
215
|
|
|
except KeyError as Argument: |
216
|
|
|
self.log_warning("An image key was not found due to following error:"+str(Argument)) |
217
|
|
|
try: |
218
|
|
|
data_obj.data = NoImageKey(data_obj, None, 0) |
219
|
|
|
entry = 'entry1/tomo_entry/instrument/detector/' |
220
|
|
|
data_obj.data._set_flat_path(entry + 'flatfield') |
221
|
|
|
data_obj.data._set_dark_path(entry + 'darkfield') |
222
|
|
|
except KeyError: |
223
|
|
|
self.log_warning("Dark/flat data was not found in input file.") |
224
|
|
|
data_obj.data._set_dark_and_flat() |
225
|
|
|
if dark: |
226
|
|
|
data_obj.data.update_dark(dark) |
227
|
|
|
if flat: |
228
|
|
|
data_obj.data.update_flat(flat) |
229
|
|
|
|
230
|
|
|
def __set_separate_dark_and_flat(self, data_obj): |
231
|
|
|
try: |
232
|
|
|
image_key = data_obj.backing_file[ |
233
|
|
|
'entry1/tomo_entry/instrument/detector/image_key'][...] |
234
|
|
|
except: |
235
|
|
|
image_key = None |
236
|
|
|
data_obj.data = NoImageKey(data_obj, image_key, 0) |
237
|
|
|
self.__set_data(data_obj, 'flat', data_obj.data._set_flat_path) |
238
|
|
|
self.__set_data(data_obj, 'dark', data_obj.data._set_dark_path) |
239
|
|
|
|
240
|
|
|
def __set_data(self, data_obj, name, func): |
241
|
|
|
path, entry, scale = self.parameters[name] |
242
|
|
|
|
243
|
|
|
if path.split('/')[0] == 'Savu': |
244
|
|
|
import os |
245
|
|
|
savu_base_path = os.path.join(os.path.dirname( |
246
|
|
|
os.path.realpath(__file__)), '..', '..', '..', '..') |
247
|
|
|
path = os.path.join(savu_base_path, path.split('Savu')[1][1:]) |
248
|
|
|
|
249
|
|
|
ffile = h5py.File(path, 'r') |
250
|
|
|
try: |
251
|
|
|
image_key = \ |
252
|
|
|
ffile['entry1/tomo_entry/instrument/detector/image_key'][...] |
253
|
|
|
func(ffile[entry], imagekey=image_key) |
254
|
|
|
except: |
255
|
|
|
func(ffile[entry]) |
256
|
|
|
|
257
|
|
|
data_obj.data._set_scale(name, scale) |
258
|
|
|
|
259
|
|
|
def _set_rotation_angles(self, data_obj): |
260
|
|
|
angles = self.parameters['angles'] |
261
|
|
|
warn_ms = "No angles found so evenly distributing them between 0 and" \ |
262
|
|
|
" 180 degrees" |
263
|
|
|
if angles is None: |
264
|
|
|
angle_key = 'entry1/tomo_entry/data/rotation_angle' |
265
|
|
|
nxs_angles = self.__get_angles_from_nxs_file(data_obj, angle_key) |
266
|
|
|
if nxs_angles is None: |
267
|
|
|
self.log_warning(warn_ms) |
268
|
|
|
angles = np.linspace(0, 180, data_obj.get_shape()[0]) |
269
|
|
|
else: |
270
|
|
|
angles = nxs_angles |
271
|
|
|
else: |
272
|
|
|
try: |
273
|
|
|
angles = eval(angles) |
274
|
|
|
except Exception as e: |
275
|
|
|
logging.warning(e) |
276
|
|
|
try: |
277
|
|
|
angles = np.loadtxt(angles) |
278
|
|
|
except Exception as e: |
279
|
|
|
logging.debug(e) |
280
|
|
|
self.log_warning(warn_ms) |
281
|
|
|
angles = np.linspace(0, 180, data_obj.get_shape()[0]) |
282
|
|
|
data_obj.meta_data.set("rotation_angle", angles) |
283
|
|
|
return len(angles) |
284
|
|
|
|
285
|
|
|
def _set_projection_shifts(self, data_obj): |
286
|
|
|
proj_shifts = np.zeros((data_obj.get_shape()[0], 2)) # initialise a 2d array of projection shifts |
287
|
|
|
self.exp.meta_data.set('projection_shifts', proj_shifts) |
288
|
|
|
data_obj.meta_data.set("projection_shifts", proj_shifts) |
289
|
|
|
return len(proj_shifts) |
290
|
|
|
|
291
|
|
|
def __get_angles_from_nxs_file(self, data_obj, path): |
292
|
|
|
if path in data_obj.backing_file: |
293
|
|
|
idx = data_obj.data.get_image_key() == 0 if \ |
294
|
|
|
isinstance(data_obj.data, ImageKey) else slice(None) |
295
|
|
|
return data_obj.backing_file[path][idx] |
296
|
|
|
else: |
297
|
|
|
self.log_warning("No rotation angle entry found in input file.") |
298
|
|
|
return None |
299
|
|
|
|
300
|
|
|
def _get_data_file(self): |
301
|
|
|
data = self.exp.meta_data.get("data_file") |
302
|
|
|
return h5py.File(data, 'r') |
303
|
|
|
|
304
|
|
|
def __check_angles(self, data_obj, n_angles): |
305
|
|
|
rot_dim = data_obj.get_data_dimension_by_axis_label("rotation_angle") |
306
|
|
|
data_angles = data_obj.data.get_shape()[rot_dim] |
307
|
|
|
if data_angles != n_angles: |
308
|
|
|
if self.nFrames > 1: |
309
|
|
|
rot_angles = data_obj.meta_data.get("rotation_angle") |
310
|
|
|
try: |
311
|
|
|
full_rotations = n_angles // data_angles |
312
|
|
|
cleaned_size = full_rotations * data_angles |
313
|
|
|
if cleaned_size != n_angles: |
314
|
|
|
rot_angles = rot_angles[0:cleaned_size] |
315
|
|
|
self.log_warning( |
316
|
|
|
"the angle list has more values than expected in it") |
317
|
|
|
rot_angles = np.reshape( |
318
|
|
|
rot_angles, [full_rotations, data_angles]) |
319
|
|
|
data_obj.meta_data.set("rotation_angle", |
320
|
|
|
np.transpose(rot_angles)) |
321
|
|
|
return |
322
|
|
|
except: |
323
|
|
|
pass |
324
|
|
|
raise Exception("The number of angles %s does not match the data " |
325
|
|
|
"dimension length %s" % (n_angles, data_angles)) |
326
|
|
|
|
327
|
|
|
def executive_summary(self): |
328
|
|
|
""" Provide a summary to the user for the result of the plugin. |
329
|
|
|
|
330
|
|
|
e.g. |
331
|
|
|
- Warning, the sample may have shifted during data collection |
332
|
|
|
- Filter operated normally |
333
|
|
|
|
334
|
|
|
:returns: A list of string summaries |
335
|
|
|
""" |
336
|
|
|
if len(self.warnings) == 0: |
337
|
|
|
return ["Nothing to Report"] |
338
|
|
|
else: |
339
|
|
|
return self.warnings |
340
|
|
|
|