Test Failed
Push — master ( f55e77...9e8bd3 )
by Daniil
02:25 queued 21s
created

Statistics._array_to_dict()   A

Complexity

Conditions 2

Size

Total Lines 5
Code Lines 5

Duplication

Lines 0
Ratio 0 %

Importance

Changes 0
Metric Value
cc 2
eloc 5
nop 2
dl 0
loc 5
rs 10
c 0
b 0
f 0
1
"""
2
.. module:: statistics
3
   :platform: Unix
4
   :synopsis: Contains and processes statistics information for each plugin.
5
6
.. moduleauthor::Jacob Williamson <[email protected]>
7
8
"""
9
10
from savu.plugins.savers.utils.hdf5_utils import Hdf5Utils
11
from savu.data.stats.stats_utils import StatsUtils
12
from savu.core.iterate_plugin_group_utils import check_if_in_iterative_loop
13
import savu.core.utils as cu
14
15
import h5py as h5
16
import numpy as np
17
import os
18
from mpi4py import MPI
19
20
21
class Statistics(object):
22
    _pattern_list = ["SINOGRAM", "PROJECTION", "TANGENTOGRAM", "VOLUME_YZ", "VOLUME_XZ", "VOLUME_XY", "VOLUME_3D", "4D_SCAN", "SINOMOVIE"]
23
    _no_stats_plugins = ["BasicOperations", "Mipmap"]
24
    _key_list = ["max", "min", "mean", "mean_std_dev", "median_std_dev", "NRMSD"]
25
    #_savers = ["Hdf5Saver", "ImageSaver", "MrcSaver", "TiffSaver", "XrfSaver"]
26
    _has_setup = False
27
28
29
    def __init__(self):
30
        self.calc_stats = True
31
        self.stats = {'max': [], 'min': [], 'mean': [], 'std_dev': [], 'RSS': [], 'data_points': []}
32
        self.stats_before_processing = {'max': [], 'min': [], 'mean': [], 'std_dev': []}
33
        self.residuals = {'max': [], 'min': [], 'mean': [], 'std_dev': []}
34
        self._repeat_count = 0
35
        self.p_num = None
36
        self.GPU = False
37
38
    def setup(self, plugin_self, pattern=None):
39
        if not Statistics._has_setup:
40
            self._setup_class(plugin_self.exp)
41
        self.plugin_name = plugin_self.name
42
        if plugin_self.name in Statistics._no_stats_plugins:
43
            self.calc_stats = False
44
        if self.calc_stats:
45
            self.plugin = plugin_self
46
            self._pad_dims = []
47
            self._already_called = False
48
            self.p_num = Statistics.count
49
            if pattern:
50
                self.pattern = pattern
51
            else:
52
                self._set_pattern_info()
53
        if self.calc_stats:
54
            Statistics._any_stats = True
55
        self._setup_iterative()
56
57
    def _setup_iterative(self):
58
        self._iterative_group = check_if_in_iterative_loop(Statistics.exp)
59
        if self._iterative_group:
60
            if self._iterative_group.start_index == Statistics.count:
61
                Statistics._loop_counter += 1
62
                Statistics.loop_stats.append({"NRMSD": np.array([])})
63
            self.l_num = Statistics._loop_counter - 1
64
65
    @classmethod
66
    def _setup_class(cls, exp):
67
        """Sets up the statistics class for the whole plugin chain (only called once)"""
68
        try:
69
            if exp.meta_data.get("stats") == "on":
70
                cls._stats_flag = True
71
            elif exp.meta_data.get("stats") == "off":
72
                cls._stats_flag = False
73
        except KeyError:
74
            cls._stats_flag = True
75
        cls._any_stats = False
76
        cls.count = 2
77
        cls.global_stats = {}
78
        cls.loop_stats = []
79
        cls.exp = exp
80
        cls.n_plugins = len(exp.meta_data.plugin_list.plugin_list)
81
        for i in range(1, cls.n_plugins + 1):
82
            cls.global_stats[i] = np.array([])
83
        cls.global_residuals = {}
84
        cls.plugin_numbers = {}
85
        cls.plugin_names = {}
86
        cls._loop_counter = 0
87
        cls._RMSD = True
88
        cls.path = exp.meta_data['out_path']
89
        if cls.path[-1] == '/':
90
            cls.path = cls.path[0:-1]
91
        cls.path = f"{cls.path}/stats"
92
        if MPI.COMM_WORLD.rank == 0:
93
            if not os.path.exists(cls.path):
94
                os.mkdir(cls.path)
95
        cls._has_setup = True
96
97
    def get_stats(self, p_num=None, stat=None, instance=-1):
98
        """Returns stats associated with a certain plugin, given the plugin number (its place in the process list).
99
100
        :param p_num: Plugin  number of the plugin whose associated stats are being fetched.
101
            If p_num <= 0, it is relative to the plugin number of the current plugin being run.
102
            E.g current plugin number = 5, p_num = -2 --> will return stats of the third plugin.
103
            By default will gather stats for the current plugin.
104
        :param stat: Specify the stat parameter you want to fetch, i.e 'max', 'mean', 'median_std_dev'.
105
            If left blank will return the whole dictionary of stats:
106
            {'max': , 'min': , 'mean': , 'mean_std_dev': , 'median_std_dev': , 'NRMSD' }
107
        :param instance: In cases where there are multiple set of stats associated with a plugin
108
            due to loops or multi-parameters, specify which set you want to retrieve, i.e 3 to retrieve the
109
            stats associated with the third run of a plugin. Pass 'all' to get a list of all sets.
110
            By default will retrieve the most recent set.
111
        """
112
        if p_num is None:
113
            p_num = self.p_num
114
        if p_num <= 0:
115
            try:
116
                p_num = self.p_num + p_num
117
            except TypeError:
118
                p_num = Statistics.count + p_num
119
        if Statistics.global_stats[p_num].ndim == 1 and instance in (None, 0, 1, -1, "all"):
120
            stats_array = Statistics.global_stats[p_num]
121
        else:
122
            if instance == "all":
123
                stats_list = [self.get_stats(p_num, stat=stat, instance=1)]
124
                n = 2
125
                if Statistics.global_stats[p_num].ndim != 1:
126
                    while n <= len(Statistics.global_stats[p_num]):
127
                        stats_list.append(self.get_stats(p_num, stat=stat, instance=n))
128
                        n += 1
129
                return stats_list
130
            if instance > 0:
131
                instance -= 1
132
            stats_array = Statistics.global_stats[p_num][instance]
133
        stats_dict = self._array_to_dict(stats_array)
134
        if stat is not None:
135
            return stats_dict[stat]
136
        else:
137
            return stats_dict
138
139
    def get_stats_from_name(self, plugin_name, n=None, stat=None, instance=-1):
140
        """Returns stats associated with a certain plugin.
141
142
        :param plugin_name: name of the plugin whose associated stats are being fetched.
143
        :param n: In a case where there are multiple instances of **plugin_name** in the process list,
144
            specify the nth instance. Not specifying will select the first (or only) instance.
145
        :param stat: Specify the stat parameter you want to fetch, i.e 'max', 'mean', 'median_std_dev'.
146
            If left blank will return the whole dictionary of stats:
147
            {'max': , 'min': , 'mean': , 'mean_std_dev': , 'median_std_dev': , 'NRMSD' }
148
        :param instance: In cases where there are multiple set of stats associated with a plugin
149
            due to iterative loops or multi-parameters, specify which set you want to retrieve, i.e 3 to retrieve the
150
            stats associated with the third run of a plugin. Pass 'all' to get a list of all sets.
151
            By default will retrieve the most recent set.
152
        """
153
        name = plugin_name
154
        if n not in (None, 0, 1):
155
            name = name + str(n)
156
        p_num = Statistics.plugin_numbers[name]
157
        return self.get_stats(p_num, stat, instance)
158
159
    def get_stats_from_dataset(self, dataset, stat=None, instance=-1):
160
        """Returns stats associated with a dataset.
161
162
        :param dataset: The dataset whose associated stats are being fetched.
163
        :param stat: Specify the stat parameter you want to fetch, i.e 'max', 'mean', 'median_std_dev'.
164
            If left blank will return the whole dictionary of stats:
165
            {'max': , 'min': , 'mean': , 'mean_std_dev': , 'median_std_dev': , 'NRMSD'}
166
        :param instance: In cases where there are multiple set of stats associated with a dataset
167
            due to iterative loops or multi-parameters, specify which set you want to retrieve, i.e 3 to retrieve the
168
            stats associated with the third run of a plugin. Pass 'all' to get a list of all sets.
169
            By default will retrieve the most recent set.
170
        """
171
        stats_list = [dataset.meta_data.get("stats")]
172
        n = 2
173
        while ("stats" + str(n)) in list(dataset.meta_data.get_dictionary().keys()):
174
            stats_list.append(dataset.meta_data.get("stats" + str(n)))
175
            n += 1
176
        if stat:
177
            for i in range(len(stats_list)):
178
                stats_list[i] = stats_list[i][stat]
179
        if instance in (None, 0, 1):
180
            stats = stats_list[0]
181
        elif instance == "all":
182
            stats = stats_list
183
        else:
184
            if instance >= 2:
185
                instance -= 1
186
            stats = stats_list[instance]
187
        return stats
188
189
    def set_slice_stats(self, my_slice, base_slice=None, pad=True):
190
        slice_stats_after = self.calc_slice_stats(my_slice, base_slice, pad=pad)
191
        if base_slice:
192
            slice_stats_before = self.calc_slice_stats(base_slice, pad=pad)
193
            for key in list(self.stats_before_processing.keys()):
194
                self.stats_before_processing[key].append(slice_stats_before[key])
195
        for key in list(self.stats.keys()):
196
            self.stats[key].append(slice_stats_after[key])
197
198
    def calc_slice_stats(self, my_slice, base_slice=None, pad=True):
199
        """Calculates and returns slice stats for the current slice.
200
201
        :param my_slice: The slice whose stats are being calculated.
202
        :param base_slice: Provide a base slice to calculate residuals from, to calculate RMSD.
203
        """
204
        if my_slice is not None:
205
            my_slice = self._de_list(my_slice)
206
            if pad:
207
                my_slice = self._unpad_slice(my_slice)
208
            slice_stats = {'max': np.amax(my_slice).astype('float64'), 'min': np.amin(my_slice).astype('float64'),
209
                           'mean': np.mean(my_slice), 'std_dev': np.std(my_slice), 'data_points': my_slice.size}
210
            if base_slice is not None and self._RMSD:
211
                base_slice = self._de_list(base_slice)
212
                base_slice = self._unpad_slice(base_slice)
213
                rss = self.calc_rss(my_slice, base_slice)
214
            else:
215
                rss = None
216
            slice_stats['RSS'] = rss
217
            return slice_stats
218
        return None
219
220
    def calc_rss(self, array1, array2):  # residual sum of squares # very slow needs looking at
221
        if array1.shape == array2.shape:
222
            residuals = np.subtract(array1, array2)
223
            rss = 0
224
            #for value in (np.nditer(residuals)):
225
            #    rss += value**2
226
            rss = np.sum(value**2 for value in np.nditer(residuals))
227
        else:
228
            #print("Warning: cannot calculate RSS, arrays different sizes.")
229
            rss = None
230
        return rss
231
232
    def rmsd_from_rss(self, rss, n):
233
        return np.sqrt(rss/n)
234
235
    def calc_rmsd(self, array1, array2):
236
        if array1.shape == array2.shape:
237
            rss = self.calc_rss(array1, array2)
238
            rmsd = self.rmsd_from_rss(rss, array1.size)
239
        else:
240
            print("Warning: cannot calculate RMSD, arrays different sizes.")  # need to make this an actual warning
241
            rmsd = None
242
        return rmsd
243
244
    def calc_stats_residuals(self, stats_before, stats_after):
245
        residuals = {'max': None, 'min': None, 'mean': None, 'std_dev': None}
246
        for key in list(residuals.keys()):
247
            residuals[key] = stats_after[key] - stats_before[key]
248
        return residuals
249
250
    def set_stats_residuals(self, residuals):
251
        self.residuals['max'].append(residuals['max'])
252
        self.residuals['min'].append(residuals['min'])
253
        self.residuals['mean'].append(residuals['mean'])
254
        self.residuals['std_dev'].append(residuals['std_dev'])
255
256
    def calc_volume_stats(self, slice_stats):
257
        volume_stats = np.array([max(slice_stats['max']), min(slice_stats['min']), np.mean(slice_stats['mean']),
258
                                np.mean(slice_stats['std_dev']), np.median(slice_stats['std_dev'])])
259
        if None not in slice_stats['RSS']:
260
            total_rss = sum(slice_stats['RSS'])
261
            n = sum(slice_stats['data_points'])
262
            RMSD = self.rmsd_from_rss(total_rss, n)
263
            the_range = volume_stats[0] - volume_stats[1]
264
            NRMSD = RMSD / the_range  # normalised RMSD (dividing by the range)
265
            volume_stats = np.append(volume_stats, NRMSD)
266
        else:
267
            #volume_stats = np.append(volume_stats, None)
268
            pass
269
        return volume_stats
270
271
    def _set_loop_stats(self):
272
        # NEED TO CHANGE THIS - MUST USE SLICES
273
        data_obj1 = list(self._iterative_group._ip_data_dict["iterating"].keys())[0]
274
        data_obj2 = self._iterative_group._ip_data_dict["iterating"][data_obj1]
275
        RMSD = self.calc_rmsd(data_obj1.data, data_obj2.data)
276
        the_range = self.get_stats(self.p_num, stat="max", instance=self._iterative_group._ip_iteration) -\
277
                self.get_stats(self.p_num, stat="min", instance=self._iterative_group._ip_iteration)
278
        NRMSD = RMSD/the_range
279
        Statistics.loop_stats[self.l_num]["NRMSD"] = np.append(Statistics.loop_stats[self.l_num]["NRMSD"], NRMSD)
280
281
    def set_volume_stats(self):
282
        """Calculates volume-wide statistics from slice stats, and updates class-wide arrays with these values.
283
        Links volume stats with the output dataset and writes slice stats to file.
284
        """
285
        stats = self.stats
286
        if self.GPU:
287
            comm = self.plugin.new_comm
288
        comm = self.plugin.get_communicator()
289
        combined_stats = self._combine_mpi_stats(stats, comm=comm)
290
        if not self.p_num:
291
            self.p_num = Statistics.count
292
        p_num = self.p_num
293
        name = self.plugin_name
294
        i = 2
295
        if not self._iterative_group:
296
            while name in list(Statistics.plugin_numbers.keys()):
297
                name = self.plugin_name + str(i)
298
                i += 1
299
        elif self._iterative_group._ip_iteration == 0:
300
            while name in list(Statistics.plugin_numbers.keys()):
301
                name = self.plugin_name + str(i)
302
                i += 1
303
304
        if p_num not in list(Statistics.plugin_names.keys()):
305
            Statistics.plugin_names[p_num] = name
306
        Statistics.plugin_numbers[name] = p_num
307
        if len(self.stats['max']) != 0:
308
            stats_array = self.calc_volume_stats(combined_stats)
309
            Statistics.global_residuals[p_num] = {}
310
            #before_processing = self.calc_volume_stats(self.stats_before_processing)
311
            #for key in list(before_processing.keys()):
312
            #    Statistics.global_residuals[p_num][key] = Statistics.global_stats[p_num][key] - before_processing[key]
313
314
            if len(Statistics.global_stats[p_num]) == 0:
315
                Statistics.global_stats[p_num] = stats_array
316
            else:
317
                Statistics.global_stats[p_num] = np.vstack([Statistics.global_stats[p_num], stats_array])
318
319
            stats_dict = self._array_to_dict(stats_array)
320
            self._link_stats_to_datasets(stats_dict, self._iterative_group)
321
322
        if self._iterative_group:
323
            if self._iterative_group.end_index == p_num and self._iterative_group._ip_iteration != 0:
324
                #self._set_loop_stats()
325
                pass
326
        self._write_stats_to_file(p_num, comm=comm)
327
        self._already_called = True
328
        self._repeat_count += 1
329
        if self._iterative_group:
330
            self.stats = {'max': [], 'min': [], 'mean': [], 'std_dev': [], 'RSS': [], 'data_points': []}
331
332
333
    def _combine_mpi_stats(self, slice_stats, comm=MPI.COMM_WORLD):
334
        combined_stats_list = comm.allgather(slice_stats)
335
        combined_stats = {'max': [], 'min': [], 'mean': [], 'std_dev': [], 'RSS': [], 'data_points': []}
336
        for single_stats in combined_stats_list:
337
            for key in list(single_stats.keys()):
338
                combined_stats[key] += single_stats[key]
339
        return combined_stats
340
341
    def _array_to_dict(self, stats_array):
342
        stats_dict = {}
343
        for i, value in enumerate(stats_array):
344
            stats_dict[Statistics._key_list[i]] = value
345
        return stats_dict
346
347
    def _set_pattern_info(self):
348
        """Gathers information about the pattern of the data in the current plugin."""
349
        out_datasets = self.plugin.get_out_datasets()
350
        try:
351
            self.pattern = self.plugin.parameters['pattern']
352
            if self.pattern == None:
353
                raise KeyError
354
        except KeyError:
355
            if not out_datasets:
356
                self.pattern = None
357
            else:
358
                patterns = out_datasets[0].get_data_patterns()
359
                for pattern in patterns:
360
                    if 1 in patterns.get(pattern)["slice_dims"]:
361
                        self.pattern = pattern
362
                        break
363
        self.calc_stats = False
364
        for dataset in out_datasets:
365
            if bool(set(Statistics._pattern_list) & set(dataset.data_info.get("data_patterns"))):
366
                self.calc_stats = True
367
368
    def _link_stats_to_datasets(self, stats_dict, iterative=False):
369
        """Links the volume wide statistics to the output dataset(s)"""
370
        out_dataset = self.plugin.get_out_datasets()[0]
371
        my_dataset = out_dataset
372
        if iterative:
373
            if "itr_clone" in out_dataset.group_name:
374
                my_dataset = list(iterative._ip_data_dict["iterating"].keys())[0]
375
        n_datasets = self.plugin.nOutput_datasets()
376
377
        i = 2
378
        group_name = "stats"
379
        #out_dataset.data_info.set([group_name], stats)
380
        while group_name in list(my_dataset.meta_data.get_dictionary().keys()):
381
            group_name = f"stats{i}"
382
            i += 1
383
        for key in list(stats_dict.keys()):
384
            my_dataset.meta_data.set([group_name, key], stats_dict[key])
385
386
    def _delete_stats_metadata(self, plugin):
387
        out_dataset = plugin.get_out_datasets()[0]
388
        out_dataset.meta_data.delete("stats")
389
390
    def _write_stats_to_file(self, p_num=None, plugin_name=None, comm=MPI.COMM_WORLD):
391
        if p_num is None:
392
            p_num = self.p_num
393
        if plugin_name is None:
394
            plugin_name = self.plugin_names[p_num]
395
        path = Statistics.path
396
        filename = f"{path}/stats.h5"
397
        stats = self.global_stats[p_num]
398
        self.hdf5 = Hdf5Utils(self.exp)
399
        self.exp._barrier(communicator=comm)
400
        if comm.rank == 0:
401
            with h5.File(filename, "a") as h5file:
402
                group = h5file.require_group("stats")
403
                if stats.shape != (0,):
404
                    if str(p_num) in list(group.keys()):
405
                        del group[str(p_num)]
406
                    dataset = group.create_dataset(str(p_num), shape=stats.shape, dtype=stats.dtype)
407
                    dataset[::] = stats[::]
408
                    dataset.attrs.create("plugin_name", plugin_name)
409
                    dataset.attrs.create("pattern", self.pattern)
410
                if self._iterative_group:
411
                    l_stats = Statistics.loop_stats[self.l_num]
412
                    group1 = h5file.require_group("iterative")
413
                    if self._iterative_group._ip_iteration == self._iterative_group._ip_fixed_iterations - 1\
414
                            and self.p_num == self._iterative_group.end_index:
415
                        dataset1 = group1.create_dataset(str(self.l_num), shape=l_stats["NRMSD"].shape, dtype=l_stats["NRMSD"].dtype)
416
                        dataset1[::] = l_stats["NRMSD"][::]
417
                        loop_plugins = []
418
                        for i in range(self._iterative_group.start_index, self._iterative_group.end_index + 1):
419
                            if i in list(self.plugin_names.keys()):
420
                                loop_plugins.append(self.plugin_names[i])
421
                        dataset1.attrs.create("loop_plugins", loop_plugins)
422
                        dataset.attrs.create("n_loop_plugins", len(loop_plugins))
0 ignored issues
show
introduced by
The variable dataset does not seem to be defined in case stats.shape != TupleNode on line 403 is False. Are you sure this can never be the case?
Loading history...
423
        self.exp._barrier(communicator=comm)
424
425
    def write_slice_stats_to_file(self, slice_stats=None, p_num=None, comm=MPI.COMM_WORLD):
426
        """Writes slice statistics to a h5 file. Placed in the stats folder in the output directory."""
427
        if not slice_stats:
428
            slice_stats = self.stats
429
        if not p_num:
430
            p_num = self.count
431
            plugin_name = self.plugin_name
432
        else:
433
            plugin_name = self.plugin_names[p_num]
434
        combined_stats = self._combine_mpi_stats(slice_stats)
435
        slice_stats_arrays = {}
436
        datasets = {}
437
        path = Statistics.path
438
        filename = f"{path}/stats_p{p_num}_{plugin_name}.h5"
439
        self.hdf5 = Hdf5Utils(self.plugin.exp)
440
        with h5.File(filename, "a", driver="mpio", comm=comm) as h5file:
441
            i = 2
442
            group_name = "/stats"
443
            while group_name in h5file:
444
                group_name = f"/stats{i}"
445
                i += 1
446
            group = h5file.create_group(group_name, track_order=None)
447
            for key in list(combined_stats.keys()):
448
                slice_stats_arrays[key] = np.array(combined_stats[key])
449
                datasets[key] = self.hdf5.create_dataset_nofill(group, key, (len(slice_stats_arrays[key]),), slice_stats_arrays[key].dtype)
450
                datasets[key][::] = slice_stats_arrays[key]
451
452
    def _unpad_slice(self, slice1):
453
        """If data is padded in the slice dimension, removes this pad."""
454
        out_datasets = self.plugin.get_out_datasets()
455
        if len(out_datasets) == 1:
456
            out_dataset = out_datasets[0]
457
        else:
458
            for dataset in out_datasets:
459
                if self.pattern in list(dataset.data_info.get(["data_patterns"]).keys()):
460
                    out_dataset = dataset
461
                    break
462
        slice_dims = out_dataset.get_slice_dimensions()
0 ignored issues
show
introduced by
The variable out_dataset does not seem to be defined for all execution paths.
Loading history...
463
        if self.plugin.pcount == 0:
464
            self._slice_list, self._pad = self._get_unpadded_slice_list(slice1, slice_dims)
465
        if self._pad:
466
            #for slice_dim in slice_dims:
467
            slice_dim = slice_dims[0]
468
            temp_slice = np.swapaxes(slice1, 0, slice_dim)
469
            temp_slice = temp_slice[self._slice_list[slice_dim]]
470
            slice1 = np.swapaxes(temp_slice, 0, slice_dim)
471
        return slice1
472
473
    def _get_unpadded_slice_list(self, slice1, slice_dims):
474
        """Creates slice object(s) to un-pad slices in the slice dimension(s)."""
475
        slice_list = list(self.plugin.slice_list[0])
476
        pad = False
477
        if len(slice_list) == len(slice1.shape):
478
            #for i in slice_dims:
479
            i = slice_dims[0]
480
            slice_width = self.plugin.slice_list[0][i].stop - self.plugin.slice_list[0][i].start
481
            if slice_width != slice1.shape[i]:
482
                pad = True
483
                pad_width = (slice1.shape[i] - slice_width) // 2  # Assuming symmetrical padding
484
                slice_list[i] = slice(pad_width, pad_width + 1, 1)
485
            return tuple(slice_list), pad
486
        else:
487
            return self.plugin.slice_list[0], pad
488
489
    def _de_list(self, slice1):
490
        """If the slice is in a list, remove it from that list."""
491
        if type(slice1) == list:
492
            if len(slice1) != 0:
493
                slice1 = slice1[0]
494
                slice1 = self._de_list(slice1)
495
        return slice1
496
497
498
    @classmethod
499
    def _count(cls):
500
        cls.count += 1
501
502
    @classmethod
503
    def _post_chain(cls):
504
        if cls._any_stats & cls._stats_flag:
505
            stats_utils = StatsUtils()
506
            stats_utils.generate_figures(f"{cls.path}/stats.h5", cls.path)
507