Extracted fold computing from troubled cell checking for Boxplot method.

Troubled_Cell_Detector.py

@@ -5,7 +5,8 @@
 TODO: Introduce Adjusted Outer Fence method in Boxplot using global_mean
     -> Done
 TODO: Introduce overlapping cells for adjacent folds in Boxplot -> Done
-TODO: Extract fold computing from TC checking
+TODO: Extract fold computing from TC checking -> Done
+TODO: Vectorize _get_cells() in Boxplot method
 TODO: Introduce lower/upper extreme outliers in Boxplot
     (each cell is also checked for neighboring domains if existing)
 TODO: Determine max_value for Theoretical only over highest degree
@@ -344,6 +345,9 @@ class Boxplot(WaveletDetector):
         Flag whether outer fences should be adjusted using global mean.
     num_overlapping_cells : int
         Number of cells overlapping with adjacent folds.
+    folds : ndarray
+        Array with indices for elements of each fold (including
+        overlaps).
 
     """
     def _reset(self, config):
@@ -362,6 +366,15 @@ class Boxplot(WaveletDetector):
         self._whisker_len = config.pop('whisker_len', 3)
         self._adjust_outer_fences = config.pop('adjust_outer_fences', True)
         self._num_overlapping_cells = config.pop('num_overlapping_cells', 1)
+        num_folds = self._mesh.num_grid_cells//self._fold_len
+        self._folds = np.zeros([num_folds, self._fold_len
+                                + 2 * self._num_overlapping_cells]).astype(int)
+        for fold in range(num_folds):
+            self._folds[fold] = np.array(
+                [i % self._mesh.num_grid_cells for i in range(
+                    fold * self._fold_len - self._num_overlapping_cells,
+                    (fold+1) * self._fold_len + self._num_overlapping_cells)])
+        # print(self._folds)
 
     def _get_cells(self, multiwavelet_coeffs, projection):
         """Calculates troubled cells using multiwavelet coefficients.
@@ -379,37 +392,34 @@ class Boxplot(WaveletDetector):
             List of indices for all detected troubled cells.
 
         """
-        indexed_coeffs = [[multiwavelet_coeffs[0, i], i]
-                          for i in range(self._mesh.num_grid_cells)]
-        # max_multiwavelets = [multiwavelet_coeffs[0, i]
+        # indexed_coeffs = [[multiwavelet_coeffs[0, i], i]
         #                   for i in range(self._mesh.num_grid_cells)]
-        # global_mean_old = np.mean(abs(np.array(max_multiwavelets)))
-        if self._adjust_outer_fences:
-            global_mean = np.mean(abs(np.array(indexed_coeffs)[:, 0]))
-        # print(global_mean == global_mean_old)
+        coeffs = multiwavelet_coeffs[0]
+        # print(coeffs.shape)
 
         if self._mesh.num_grid_cells < self._fold_len:
             self._fold_len = self._mesh.num_grid_cells
 
         num_folds = self._mesh.num_grid_cells//self._fold_len
         troubled_cells = []
+        # troubled_cells_new = []
 
         for fold in range(num_folds):
-            indices = np.array([i % len(indexed_coeffs) for i in range(
-                fold * self._fold_len - self._num_overlapping_cells,
-                (fold+1) * self._fold_len + self._num_overlapping_cells)])
-            sorted_fold = sorted(np.array(indexed_coeffs)[indices],
-                                 key=lambda x: x[0])
+            # indexed_fold = np.array(indexed_coeffs)[self._folds[fold]]
+            # sorted_fold_old = indexed_fold[indexed_fold[:, 0].argsort()]
+
+            sorted_fold = sorted(coeffs[self._folds[fold]])
+            # print(sorted_fold == sorted_fold_old[:, 0])
 
             boundary_index = self._fold_len//4
             balance_factor = self._fold_len/4.0 - boundary_index
 
             first_quartile = (1-balance_factor) \
-                * sorted_fold[boundary_index-1][0] \
-                + balance_factor * sorted_fold[boundary_index][0]
+                * sorted_fold[boundary_index-1] \
+                + balance_factor * sorted_fold[boundary_index]
             third_quartile = (1-balance_factor) \
-                * sorted_fold[3*boundary_index-1][0]\
-                + balance_factor * sorted_fold[3*boundary_index][0]
+                * sorted_fold[3*boundary_index-1]\
+                + balance_factor * sorted_fold[3*boundary_index]
 
             lower_bound = first_quartile \
                 - self._whisker_len * (third_quartile-first_quartile)
@@ -418,22 +428,49 @@ class Boxplot(WaveletDetector):
 
             # Adjust outer fences if flag is set
             if self._adjust_outer_fences:
+                global_mean = np.mean(abs(coeffs))
                 lower_bound = min(-global_mean, lower_bound)
                 upper_bound = max(global_mean, upper_bound)
 
-            # Check for lower extreme outliers and add respective cells
-            for cell in sorted_fold:
-                if cell[0] < lower_bound:
-                    troubled_cells.append(cell[1])
-                else:
-                    break
-
-            # Check for upper extreme outliers and add respective cells
-            for cell in sorted_fold[::-1][:]:
-                if cell[0] > upper_bound:
-                    troubled_cells.append(cell[1])
-                else:
-                    break
+            # # Check for lower extreme outliers and add respective cells
+            # for cell in sorted_fold:
+            #     if cell[0] < lower_bound:
+            #         troubled_cells.append(int(cell[1]))
+            #     else:
+            #         break
+            #
+            # # Check for upper extreme outliers and add respective cells
+            # for cell in sorted_fold[::-1][:]:
+            #     if cell[0] > upper_bound:
+            #         troubled_cells.append(int(cell[1]))
+            #     else:
+            #         break
+
+            # Check for extreme outlier and add respective cells
+            for cell in self._folds[fold]:
+                if (coeffs[cell] > upper_bound) \
+                        or (coeffs[cell] < lower_bound):
+                    troubled_cells.append(int(cell))
+
+            # print(upper_bound, lower_bound)
+            # print(sorted_fold_new)
+            # print(type(sorted_fold_new))
+            # print(sorted_fold_new > upper_bound)
+            # print(sorted_fold_new < lower_bound)
+            # test =
+            # print(type(test), test)
+            # print(list(test), list(test[0]))
+
+            # troubled_cells_new += list(np.flatnonzero(np.logical_or(
+            #     sorted_fold_new > upper_bound,
+            #     sorted_fold_new < lower_bound)).astype(int))
+            # print(troubled_cells_new)
+
+        # troubled_cells_new = sorted(troubled_cells_new)
+        # print(troubled_cells_new)
+        # print(troubled_cells)
+        # print(sorted(troubled_cells) == sorted(troubled_cells_new))
+        # print(type(troubled_cells_new[0]), type(troubled_cells[0]))
 
         return sorted(troubled_cells)