1"""The ten statistical / heuristic anomaly checks (e1..e10).
2
3Each public ``check_*`` function consumes a column-oriented numpy
4:data:`Window` (plus the device ``config`` where needed) and returns a
5``dict`` whose keys are the error-coded field names the detector will
6emit. Every function handles the empty / ``None`` window case by
7returning ``{}``.
8
9The error codes are:
10
11* ``e1`` constant data (per parameter)
12* ``e2`` threshold breach, lower/upper (per parameter, ``_l``/``_u`` suffix)
13* ``e3`` data completeness (global)
14* ``e4`` dust-order violation (global)
15* ``e5`` low humidity / low t2 (single flag)
16* ``e6`` missing parameter keys (per parameter) — delegated
17* ``e7`` temperature below dew point (condensing environment)
18* ``e8`` RH-affected dust (global flag)
19* ``e9`` VOC sensor fault (v21)
20* ``e10`` RH-affected VOC (global flag)
21"""
22
23from __future__ import annotations
24
25import logging
26import math
27from typing import Any, Dict, List, Optional, Tuple
28
29import numpy as np
30
31from OzWrapper.OzAnomaly.data_provider import Window, is_empty
32from OzWrapper.OzAnomaly.preprocessing import get_ogs_list, para_check, select_columns
33from OzWrapper.OzAnomaly.thresholds import (
34 CONSTANT_CHECK_SKIP_LIST,
35 THRESHOLD_CHECK_SKIP_LIST,
36 THRESHOLD_DICT,
37)
38
39logger = logging.getLogger(__name__)
40
41# Columns the constant-data check (e1) must never flag. v25 is the VOC raw
42# channel that legitimately sits flat on many devices.
43CONSTANT_CHECK_COL_IGNORE_LIST = ["v25"]
44
45
46# ---------------------------------------------------------------------------
47# Private skip helpers — decide whether a specific sensor code should be
48# excluded from a given check, based on the OGS cartridge label(s) in config.
49# ---------------------------------------------------------------------------
50
51
[docs]
52def _skip_evaluation_constant(ogs_list: list, sensor_code: str) -> bool:
53 """Return True when the constant-data check (e1) should skip ``sensor_code``."""
54 if sensor_code not in CONSTANT_CHECK_SKIP_LIST:
55 return False
56 lb_values = {ogs_device.get("lb") for ogs_device in ogs_list}
57 return any(lb in CONSTANT_CHECK_SKIP_LIST[sensor_code] for lb in lb_values)
58
59
[docs]
60def _skip_evaluation_threshold(ogs_label: Optional[str], sensor_code: str) -> bool:
61 """Return True when the threshold-breach check (e2) should skip ``sensor_code``."""
62 if sensor_code not in THRESHOLD_CHECK_SKIP_LIST:
63 return False
64 return ogs_label in THRESHOLD_CHECK_SKIP_LIST[sensor_code]
65
66
[docs]
67def _row_count(window: Window) -> int:
68 """Total rows in the window (uses the ``t`` column which is always present)."""
69 return len(window["t"])
70
71
72# ---------------------------------------------------------------------------
73# Check functions
74# ---------------------------------------------------------------------------
75
76
[docs]
77def check_constant_instances(
78 config: Dict[str, Any], device_type: Optional[str], window: Optional[Window]
79) -> Dict[str, int]:
80 """Detect long runs of identical values per sensor column (e1)."""
81 if is_empty(window):
82 return {}
83
84 sensor_columns = select_columns(window)
85 if not sensor_columns:
86 return {}
87
88 ogs_list = get_ogs_list(config, device_type)
89 n_rows = _row_count(window)
90 constant_threshold = n_rows * 0.2
91
92 constant_instances: Dict[str, Dict[str, Any]] = {}
93
94 for column in sensor_columns:
95 if column in CONSTANT_CHECK_COL_IGNORE_LIST:
96 continue
97
98 if _skip_evaluation_constant(ogs_list, column):
99 continue
100
101 if column not in window:
102 continue
103
104 arr = window[column]
105 consecutive_count = 1
106 prev_value = arr[0]
107 total_consecutive_count = 0
108 constant_values: list = []
109
110 # NaN comparisons return False, so a run of NaN is treated as
111 # "non-constant" — same behaviour as pandas Series equality.
112 for i in range(1, n_rows):
113 v = arr[i]
114 if v == prev_value:
115 consecutive_count += 1
116 else:
117 if consecutive_count >= constant_threshold:
118 total_consecutive_count += 1
119 constant_values.append((prev_value, consecutive_count))
120 consecutive_count = 1
121 prev_value = v
122
123 if consecutive_count >= constant_threshold:
124 total_consecutive_count += 1
125 constant_values.append((prev_value, consecutive_count))
126
127 if total_consecutive_count > 0:
128 constant_instances[column] = {
129 "count": total_consecutive_count,
130 "instances": constant_values,
131 }
132
133 if not constant_instances:
134 return {}
135
136 final_dictionary: Dict[str, int] = {}
137 for key, value in constant_instances.items():
138 if key in CONSTANT_CHECK_COL_IGNORE_LIST:
139 continue
140 sum_second_element = sum(item[1] for item in value["instances"])
141 final_dictionary[f"{key}_e1"] = int((sum_second_element * 100) / n_rows)
142 return final_dictionary
143
144
[docs]
145def check_threshold_breach(
146 config: Dict[str, Any], device_type: Optional[str], window: Optional[Window]
147) -> Dict[str, int]:
148 """Flag columns where >30% of samples sit outside their bounds (e2, e5)."""
149 if is_empty(window):
150 return {}
151
152 ogs_list = get_ogs_list(config, device_type)
153 tdict = config.get("soft", {}).get("config", {}).get("tdict", {})
154 breached_records: Dict[str, int] = {}
155
156 n_rows = _row_count(window)
157
158 for col in window.keys():
159 if col == "t":
160 continue
161 if any(_skip_evaluation_threshold(ogs.get("lb"), col) for ogs in ogs_list):
162 continue
163
164 arr = window[col]
165 lower_key = f"{col}_l"
166 upper_key = f"{col}_u"
167
168 # Soft thresholds (tdict) win over the hardcoded defaults.
169 lower_threshold = tdict.get(lower_key, THRESHOLD_DICT.get(lower_key))
170 upper_threshold = tdict.get(upper_key, THRESHOLD_DICT.get(upper_key))
171
172 # NaN comparisons return False -> NaN samples are not counted as
173 # breaches, matching the pandas behaviour.
174 if lower_threshold is not None:
175 n_breached = int(np.sum(arr < lower_threshold))
176 if n_breached > n_rows * 0.3:
177 breached_records[f"{col}_e2_l"] = int((n_breached * 100) / n_rows)
178
179 if upper_threshold is not None:
180 n_breached = int(np.sum(arr > upper_threshold))
181 if n_breached > n_rows * 0.3:
182 breached_records[f"{col}_e2_u"] = int((n_breached * 100) / n_rows)
183
184 # Low humidity marker — t2 is preferred, falls back to hum. Any non-zero
185 # number of breaches counts (not the 30% rule).
186 if "t2" in window:
187 t2_threshold = tdict.get("t2", THRESHOLD_DICT.get("t2"))
188 if t2_threshold is not None:
189 n_low = int(np.sum(window["t2"] < t2_threshold))
190 if n_low > 0:
191 breached_records["t2_e5"] = int((n_low * 100) / n_rows)
192 elif "hum" in window:
193 hum_threshold = tdict.get("hum", THRESHOLD_DICT.get("hum"))
194 if hum_threshold is not None:
195 n_low = int(np.sum(window["hum"] < hum_threshold))
196 if n_low > 0:
197 breached_records["hum_e5"] = int((n_low * 100) / n_rows)
198
199 return breached_records
200
201
[docs]
202def check_data_completeness(
203 config: Dict[str, Any], window: Optional[Window]
204) -> Dict[str, int]:
205 """Report the shortfall from expected sample count (e3)."""
206 if is_empty(window):
207 return {}
208
209 try:
210 interval = config["cloud"]["interval"]
211 except (KeyError, TypeError):
212 logger.warning("cloud.interval missing from config; skipping e3 check")
213 return {}
214
215 n_rows = _row_count(window)
216 data_target = 86400 / (interval / 1000)
217 percentage_data = round((n_rows * 100 / data_target), 2)
218 if percentage_data < 70:
219 return {"e3": 100 - int(percentage_data)}
220 return {}
221
222
[docs]
223def check_dust_order(window: Optional[Window]) -> Dict[str, float]:
224 """Flag violations of the dust-size ordering invariant (e4)."""
225 if is_empty(window):
226 return {}
227
228 required_cols = ["z1", "z2", "z3", "z4"]
229 available_cols = [col for col in required_cols if col in window]
230 if not available_cols:
231 return {}
232
233 n_rows = _row_count(window)
234
235 # "keep" = rows where ALL available dust cols are non-NaN. Mirrors
236 # pandas' df.dropna(subset=available_cols).
237 keep = np.ones(n_rows, dtype=bool)
238 for c in available_cols:
239 keep &= ~np.isnan(window[c])
240
241 # Build the violation mask on the same column set the original code
242 # branches over. NaN comparisons yield False, so NaN rows can't appear
243 # in the violation mask anyway — but we still AND with `keep` so
244 # the count we divide matches the dropna-then-count semantics.
245 has_z1 = "z1" in window
246 has_z2 = "z2" in window
247 has_z3 = "z3" in window
248 has_z4 = "z4" in window
249
250 if has_z1 and has_z2 and has_z3 and has_z4:
251 z1, z2, z3, z4 = window["z1"], window["z2"], window["z3"], window["z4"]
252 violation = (z4 < z2) | (z2 < z1) | (z1 < z3)
253 elif has_z1 and has_z2 and has_z3 and not has_z4:
254 z1, z2, z3 = window["z1"], window["z2"], window["z3"]
255 violation = (z2 < z1) | (z1 < z3)
256 elif has_z1 and has_z2 and not has_z3 and not has_z4:
257 z1, z2 = window["z1"], window["z2"]
258 violation = z2 < z1
259 else:
260 violation = np.zeros(n_rows, dtype=bool)
261
262 n_violations = int(np.sum(keep & violation))
263 if n_violations > 0:
264 e4 = (n_violations * 100) / n_rows
265 if e4 > 30:
266 return {"e4": e4}
267 return {}
268
269
[docs]
270def check_missing_parameters(
271 window: Optional[Window], config: Dict[str, Any]
272) -> Dict[str, int]:
273 """Flag parameters expected by config but absent from the window (e6)."""
274 return para_check(window, config)
275
276
[docs]
277def check_temp_dewpoint(window: Optional[Window]) -> Dict[str, int]:
278 """Flag a likely condensing environment (e7)."""
279 if is_empty(window):
280 return {}
281
282 temp_col = "t1" if "t1" in window else ("temp" if "temp" in window else None)
283 if not temp_col:
284 return {}
285
286 dp_col = "dp1" if "dp1" in window else ("dp" if "dp" in window else None)
287 if not dp_col:
288 return {}
289
290 temp_arr = window[temp_col]
291 dp_arr = window[dp_col]
292 n_rows = _row_count(window)
293
294 keep = ~np.isnan(temp_arr) & ~np.isnan(dp_arr)
295 err = temp_arr < (dp_arr + 2)
296 n_error = int(np.sum(keep & err))
297
298 if n_error > 0:
299 return {"e7": int((n_error * 100) / n_rows)}
300 return {}
301
302
[docs]
303def detect_pm_voc(
304 window: Optional[Window],
305 mode: str,
306 rh_col: str = "hum",
307 cols: Optional[list] = None,
308 degree: int = 2,
309 rh_threshold: Optional[Tuple[float, float]] = None,
310 r2_threshold: float = 0.3,
311 spike_rh_min: float = 80,
312 spike_z_thresh: float = 1.0,
313 cluster_dist: float = 1.0,
314 min_cluster_points: int = 3,
315) -> Dict[str, Any]:
316 """Test whether humidity drives PM (mode="pm") or VOC (mode="voc") signal.
317
318 Pure numpy implementation — no scikit-learn. Polynomial fit via
319 ``np.polyfit``; R² computed manually; spike detection via z-score
320 threshold; cluster decision via O(n²) distance check.
321 """
322 if is_empty(window):
323 return {}
324 if rh_col not in window:
325 return {}
326
327 if mode == "pm":
328 err_code = "e8"
329 cols = cols or ["z1", "z2", "z3", "z4"]
330 rh_threshold = rh_threshold or (75, 100)
331 elif mode == "voc":
332 err_code = "e10"
333 cols = cols or ["v21"]
334 rh_threshold = rh_threshold or (50, 100)
335 else:
336 return {}
337
338 rh_arr = window[rh_col]
339 sub_mask = (rh_arr >= rh_threshold[0]) & (rh_arr <= rh_threshold[1])
340 if not bool(sub_mask.any()):
341 return {}
342
343 result: Dict[str, Any] = {err_code: {}}
344 overall_flag = False
345
346 for col in cols:
347 if col not in window:
348 continue
349
350 col_arr = window[col]
351 # Rows in RH range AND non-NaN in both rh and target column.
352 keep = sub_mask & ~np.isnan(rh_arr) & ~np.isnan(col_arr)
353 if not bool(keep.any()):
354 continue
355
356 X = rh_arr[keep]
357 y = col_arr[keep]
358
359 if X.size < degree + 1:
360 # Not enough points to fit a polynomial of this degree.
361 continue
362
363 # Polynomial fit via numpy — no sklearn.
364 try:
365 coeffs = np.polyfit(X, y, degree)
366 y_pred = np.polyval(coeffs, X)
367 except (np.linalg.LinAlgError, ValueError) as exc:
368 logger.warning("polyfit failed for %s: %s", col, exc)
369 continue
370
371 # Manual R² (guard against zero-variance targets).
372 ss_res = float(np.sum((y - y_pred) ** 2))
373 ss_tot = float(np.sum((y - np.mean(y)) ** 2))
374 r2 = 0.0 if ss_tot == 0 else 1 - (ss_res / ss_tot)
375
376 # Z-score spikes — fall back to 1e-6 std to avoid /0 when all y equal.
377 mean_y = float(np.mean(y))
378 std_y = float(np.std(y)) or 1e-6
379 z_scores = np.abs((y - mean_y) / std_y)
380
381 spike_mask = (z_scores > spike_z_thresh) & (X >= spike_rh_min)
382 rh_anom = X[spike_mask]
383 gas_anom = y[spike_mask]
384
385 # Distance-based clustering — a point is clustered if it has
386 # at least `min_cluster_points` neighbours within `cluster_dist`
387 # on both axes (including itself).
388 spikes_flag = False
389 if rh_anom.size >= min_cluster_points:
390 clusters = 0
391 for i in range(rh_anom.size):
392 neighbors = int(np.sum(
393 (np.abs(rh_anom - rh_anom[i]) <= cluster_dist)
394 & (np.abs(gas_anom - gas_anom[i]) <= cluster_dist)
395 ))
396 if neighbors >= min_cluster_points:
397 clusters += 1
398 spikes_flag = clusters > 0
399
400 flag = (r2 >= r2_threshold) or spikes_flag
401 if flag:
402 result[err_code][col] = {
403 "flag": flag,
404 "r2": round(r2, 2),
405 "spikes": spikes_flag,
406 }
407 overall_flag = True
408
409 if overall_flag:
410 result[err_code]["overall_flag"] = True
411 return result
412 return {}
413
414
[docs]
415def check_rh_affected_dust(window: Optional[Window]) -> Dict[str, Any]:
416 """Humidity vs dust (e8). Wraps :func:`detect_pm_voc` with ``mode='pm'``."""
417 return detect_pm_voc(window, mode="pm", degree=3)
418
419
[docs]
420def check_voc_sensor_fault(window: Optional[Window]) -> Dict[str, int]:
421 """Flag a stuck VOC sensor (e9).
422
423 If the max-min range of ``v21`` across the window is below 3 mV, the
424 sensor is treated as faulty.
425 """
426 if is_empty(window):
427 return {}
428
429 col = "v21"
430 threshold = 3
431 if col not in window:
432 return {}
433
434 arr = window[col]
435 arr = arr[~np.isnan(arr)] # numpy.min/max propagate NaN; pandas skips by default
436 if arr.size == 0:
437 return {}
438
439 diff = float(arr.max() - arr.min())
440 if diff < threshold:
441 return {f"{col}_e9": 1}
442 return {}
443
444
[docs]
445def check_rh_affected_voc(window: Optional[Window]) -> Dict[str, Any]:
446 """Humidity vs VOC (e10). Wraps :func:`detect_pm_voc` with ``mode='voc'``."""
447 return detect_pm_voc(window, mode="voc", degree=2)