From 0c788ea2256056a797d322b8903d9c858a02ce06 Mon Sep 17 00:00:00 2001 From: Alex Tudorica Date: Tue, 19 May 2026 23:18:12 +0200 Subject: [PATCH 1/3] Improve fireball detector recall with stdpixel decontamination A bright fireball inflates stdpixel along its trail within the 256-frame FF block, pushing the threshold (avepixel + k1*stdpixel + j1) above 255. The uint8 clip kills near-saturated trail pixels (250-254). Fix: before thresholding, replace stdpixel at contaminated pixels (maxpixel >= P90 AND stdpixel > 3x background median) with the background median. Pure numpy, no Cython changes needed. Also: clip threshold to 254 (not 255) so near-saturated pixels can pass, and replace the hard-coded 2 px/frame fireball velocity cap with a configurable angular velocity (fireball_max_ang_vel, default 60 deg/s) converted to a per-sensor px/frame limit at runtime using fps and the platepar-derived deg/pixel scale. Tested on 129-event reference dataset (1755 FF files, 743 stations): FR-positive recall 99.0% (679/686), missed station rescue 74.0% (77/104), overall station recall 86.0% -> 96.4% (716/743). Zero threshold failures on real data. Line finder has zero failures in all 686 tests -- the problem was entirely upstream in the threshold stage. --- .config | 7 ++++ RMS/ConfigReader.py | 10 ++++- RMS/Routines/Grouping3D.py | 69 ++++++++++++++++++++++------------- RMS/Routines/Grouping3Dcy.pyx | 6 +-- RMS/VideoExtraction.py | 28 +++++++++++++- 5 files changed, 89 insertions(+), 31 deletions(-) diff --git a/.config b/.config index 35c5d4677..9891a1940 100644 --- a/.config +++ b/.config @@ -516,6 +516,13 @@ point_ratio_threshold: 0.7 ; Maximum number of lines which are allowed to be found on the image max_lines: 5 +; Maximum angular velocity for fireball candidates in deg/s. Converted to a per-sensor +; px/frame cap at runtime using fps and the platepar-derived deg/pixel scale, so the same +; value transfers between cameras with different FOV/resolution. Default is wider than the +; meteor detector's ang_vel_max (35 deg/s) because the fireball detector is recall-first +; and the cap exists mostly to reject birds/planes, not to constrain real fireballs. +fireball_max_ang_vel: 60.0 + [MeteorDetection] diff --git a/RMS/ConfigReader.py b/RMS/ConfigReader.py index 0d7de609f..676177383 100644 --- a/RMS/ConfigReader.py +++ b/RMS/ConfigReader.py @@ -560,6 +560,12 @@ def __init__(self): self.line_distance_const = 4 # constant that determines the influence of average point distance on the line quality self.point_ratio_threshold = 0.7# ratio of how many points must be close to the line before considering searching for another line self.max_lines = 5 # maximum number of lines + # Maximum angular velocity for fireball candidates in deg/s. Converted to a per-sensor + # px/frame cap at runtime using fps and the platepar-derived deg/pixel scale, so the same + # value transfers between cameras with different FOV/resolution. Default is wider than + # the meteor detector's ang_vel_max (35 deg/s) because the fireball detector is recall- + # first and the cap exists mostly to reject birds/planes, not to constrain real fireballs. + self.fireball_max_ang_vel = 60.0 ##### MeteorDetection @@ -1512,7 +1518,9 @@ def parseFireballDetection(config, parser): if parser.has_option(section, "max_lines"): config.max_lines = parser.getint(section, "max_lines") - + + if parser.has_option(section, "fireball_max_ang_vel"): + config.fireball_max_ang_vel = parser.getfloat(section, "fireball_max_ang_vel") diff --git a/RMS/Routines/Grouping3D.py b/RMS/Routines/Grouping3D.py index ec13f34d9..52e3ff886 100644 --- a/RMS/Routines/Grouping3D.py +++ b/RMS/Routines/Grouping3D.py @@ -138,21 +138,49 @@ def __init__(self, config): -def findCoefficients(line_list): +def findCoefficients(line_list, config=None): """ Extract coefficients from list of lines that can be consumed by RMS.VideoExtraction. - + Arguments: line_list: [list] list of detected lines - + config: [config object] configuration parameters (optional). Used to convert the + fireball_max_ang_vel limit (deg/s) into a per-sensor px/frame threshold via + fps and the average deg/pixel scale (fov_h/height + fov_w/width)/2. + Return: - coeff: [list] coefficients for each detected line in format: [first point, slope of XZ, slope of YZ, + coeff: [list] coefficients for each detected line in format: [first point, slope of XZ, slope of YZ, first frame, last frame] """ - + + # Fallback px/frame cap matches the legacy hard-coded value when no config is available + # (e.g. unit tests). Stations always pass a config in production. + max_velocity_px = 2.0 + + if config is not None: + max_ang_vel = getattr(config, 'fireball_max_ang_vel', None) + fps = getattr(config, 'fps', None) + fov_h = getattr(config, 'fov_h', None) + fov_w = getattr(config, 'fov_w', None) + height = getattr(config, 'height', None) + width = getattr(config, 'width', None) + f = getattr(config, 'f', None) + + if (max_ang_vel is not None and fps and fov_h and fov_w and height and width and f): + # deg/s -> original px/frame via the average deg/pixel scale and fps, + # then -> subsampled px/frame via the decimation factor f (thresholdAndSubsample + # divides x,y by f, so slopes are in subsampled coordinates). + scale = (fov_h/float(height) + fov_w/float(width))/2.0 + if scale > 0: + max_velocity_px = (max_ang_vel/scale)/float(fps)/float(f) + else: + log.warning("findCoefficients: fireball_max_ang_vel conversion skipped " + "(missing fps/fov/resolution/f); using fallback %.2f px/frame", + max_velocity_px) + coeff = [] - + for detected_line in line_list: - + if detected_line[0][2] < detected_line[1][2]: point1 = np.array(detected_line[0], dtype=np.float64) point2 = np.array(detected_line[1], dtype=np.float64) @@ -160,30 +188,21 @@ def findCoefficients(line_list): point1 = np.array(detected_line[1], dtype=np.float64) point2 = np.array(detected_line[0], dtype=np.float64) else: - # skip if points are on the same frame (that shouldn't happen, though) log.debug("Points on the same frame!") continue - - # difference between last point and first point that represent a line + point3 = point2 - point1 - - # slope - slopeXZ = point3[1]/point3[2] # speed on X axis - slopeYZ = point3[0]/point3[2] # speed on Y axis - - # length of velocity vector + + slopeXZ = point3[1]/point3[2] + slopeYZ = point3[0]/point3[2] + total = sqrt(slopeXZ**2 + slopeYZ**2) - - #print('Fireball slope:', total) - # ignore line if too fast - # TODO: this limit should be read from config file and calculated for FOV - # 1.6 is better estimate on upper speed limit, set to 2 for safety - if total > 2: + if total > max_velocity_px: continue - - coeff.append([point1, slopeXZ, slopeYZ, detected_line[4], detected_line[5]]) #first point, slope of XZ, slope of YZ, first frame, last frame - + + coeff.append([point1, slopeXZ, slopeYZ, detected_line[4], detected_line[5]]) + return coeff diff --git a/RMS/Routines/Grouping3Dcy.pyx b/RMS/Routines/Grouping3Dcy.pyx index f84203ee8..dd9ecce16 100644 --- a/RMS/Routines/Grouping3Dcy.pyx +++ b/RMS/Routines/Grouping3Dcy.pyx @@ -487,9 +487,9 @@ def thresholdAndSubsample(np.ndarray[UINT8_TYPE_t, ndim=3] frames, \ # Compute the threshold limit avg_std = int(float(compressed[2, y, x]) + k1*float(compressed[3, y, x])) + j1 - # Make sure the threshold limit is not above the maximum possible value - if avg_std > 255: - avg_std = 255 + # Clip threshold to 254 so that near-saturated pixels (250-254) can still pass + if avg_std > 254: + avg_std = 254 if ((max_val > min_level) and (max_val >= avg_std)): diff --git a/RMS/VideoExtraction.py b/RMS/VideoExtraction.py index 02362cebc..6bbfc3fb5 100644 --- a/RMS/VideoExtraction.py +++ b/RMS/VideoExtraction.py @@ -68,8 +68,32 @@ def findPoints(self): (y, x, z): [tuple] Coordinates of points that form the fireball, where Z is the frame number. """ + # Decontaminate stdpixel before thresholding: a bright fireball inflates stdpixel + # along its trail, pushing the threshold (avg + k1*std + j1) above 255. Replace + # stdpixel at contaminated pixels with the background median. Only copy the + # compressed array when contamination is present, so the no-op path is zero-cost. + maxpix = self.compressed[0] + stdpix = self.compressed[3] + bright_threshold = np.percentile(maxpix, 90) + bright_mask = maxpix >= bright_threshold + bg_mask = ~bright_mask + # Exclude masked-out pixels (camera borders/obstructions) from the background sample + if self.mask is not None: + bg_mask = bg_mask & (self.mask.img > 0) + if np.count_nonzero(bg_mask) > 100: + std_ref = float(np.median(stdpix[bg_mask])) + else: + std_ref = float(np.median(stdpix)) + contaminated = bright_mask & (stdpix > 3*max(std_ref, 1)) + if np.any(contaminated): + compressed = self.compressed.copy() + replacement = max(1, int(round(min(std_ref, 255)))) + compressed[3][contaminated] = np.uint8(replacement) + else: + compressed = self.compressed + # Threshold and subsample frames - length, x, y, z = Grouping3D.thresholdAndSubsample(self.frames, self.compressed, \ + length, x, y, z = Grouping3D.thresholdAndSubsample(self.frames, compressed, \ self.config.min_level, self.config.min_pixels, self.config.k1, self.config.j1, self.config.f) @@ -302,7 +326,7 @@ def executeAll(self): t = time.time() # Find the parameters of the line in the 3D point cloud - coeff = Grouping3D.findCoefficients(line_list) + coeff = Grouping3D.findCoefficients(line_list, self.config) log.debug("[" + self.filename + "] Time for finding coefficients: " + str(time.time() - t) + "s") if len(coeff) == 0: From 63b743ff442c37edde51e581ce5fd171c4070423 Mon Sep 17 00:00:00 2001 From: Alexandru Tudorica Date: Tue, 2 Jun 2026 08:30:31 +0200 Subject: [PATCH 2/3] fix: review fixes for fireball decontamination PR - Exclude masked pixels from contamination set, not just background sample - Use O(n) np.partition instead of O(n log n) np.percentile for P90 - Copy only stdpixel layer instead of full compressed array (~2.8 MB saved) - Report which specific config attributes are missing/zero in warning - Fix Cython comment: only max_val=254 is newly enabled, not 250-254 - Document fallback velocity cap (~39 deg/s) vs production default (60 deg/s) - Add from __future__ import division for Python 2 safety --- RMS/Routines/Grouping3D.py | 15 +++++++++++---- RMS/Routines/Grouping3Dcy.pyx | 2 +- RMS/VideoExtraction.py | 20 +++++++++++++++----- 3 files changed, 27 insertions(+), 10 deletions(-) diff --git a/RMS/Routines/Grouping3D.py b/RMS/Routines/Grouping3D.py index 52e3ff886..7d30f8f40 100644 --- a/RMS/Routines/Grouping3D.py +++ b/RMS/Routines/Grouping3D.py @@ -153,7 +153,8 @@ def findCoefficients(line_list, config=None): """ # Fallback px/frame cap matches the legacy hard-coded value when no config is available - # (e.g. unit tests). Stations always pass a config in production. + # (e.g. unit tests, external callers). At default 720p/25fps/f=16 this is ~39 deg/s -- + # tighter than the production default of 60 deg/s, so tests get conservative behavior. max_velocity_px = 2.0 if config is not None: @@ -165,7 +166,13 @@ def findCoefficients(line_list, config=None): width = getattr(config, 'width', None) f = getattr(config, 'f', None) - if (max_ang_vel is not None and fps and fov_h and fov_w and height and width and f): + # Check each required attribute individually for clearer diagnostics + missing = [name for name, val in [('fireball_max_ang_vel', max_ang_vel), + ('fps', fps), ('fov_h', fov_h), ('fov_w', fov_w), + ('height', height), ('width', width), ('f', f)] + if not val] + + if not missing: # deg/s -> original px/frame via the average deg/pixel scale and fps, # then -> subsampled px/frame via the decimation factor f (thresholdAndSubsample # divides x,y by f, so slopes are in subsampled coordinates). @@ -174,8 +181,8 @@ def findCoefficients(line_list, config=None): max_velocity_px = (max_ang_vel/scale)/float(fps)/float(f) else: log.warning("findCoefficients: fireball_max_ang_vel conversion skipped " - "(missing fps/fov/resolution/f); using fallback %.2f px/frame", - max_velocity_px) + "(missing/zero: %s); using fallback %.2f px/frame", + ', '.join(missing), max_velocity_px) coeff = [] diff --git a/RMS/Routines/Grouping3Dcy.pyx b/RMS/Routines/Grouping3Dcy.pyx index dd9ecce16..8ce148bf2 100644 --- a/RMS/Routines/Grouping3Dcy.pyx +++ b/RMS/Routines/Grouping3Dcy.pyx @@ -487,7 +487,7 @@ def thresholdAndSubsample(np.ndarray[UINT8_TYPE_t, ndim=3] frames, \ # Compute the threshold limit avg_std = int(float(compressed[2, y, x]) + k1*float(compressed[3, y, x])) + j1 - # Clip threshold to 254 so that near-saturated pixels (250-254) can still pass + # Clip threshold to 254 so that pixels at exactly max_val=254 can still pass if avg_std > 254: avg_std = 254 diff --git a/RMS/VideoExtraction.py b/RMS/VideoExtraction.py index 6bbfc3fb5..2b21f4da2 100644 --- a/RMS/VideoExtraction.py +++ b/RMS/VideoExtraction.py @@ -16,6 +16,8 @@ +from __future__ import division + import math import time from multiprocessing import Process, Event @@ -70,23 +72,31 @@ def findPoints(self): # Decontaminate stdpixel before thresholding: a bright fireball inflates stdpixel # along its trail, pushing the threshold (avg + k1*std + j1) above 255. Replace - # stdpixel at contaminated pixels with the background median. Only copy the - # compressed array when contamination is present, so the no-op path is zero-cost. + # stdpixel at contaminated pixels with the background median. Only the stdpixel + # layer is copied when contamination is present, so the no-op path is zero-cost. maxpix = self.compressed[0] stdpix = self.compressed[3] - bright_threshold = np.percentile(maxpix, 90) + # O(n) partition-based percentile instead of O(n log n) full sort + bright_threshold = np.partition(maxpix.ravel(), -maxpix.size // 10)[-maxpix.size // 10] bright_mask = maxpix >= bright_threshold bg_mask = ~bright_mask # Exclude masked-out pixels (camera borders/obstructions) from the background sample if self.mask is not None: - bg_mask = bg_mask & (self.mask.img > 0) + mask_valid = self.mask.img > 0 + bg_mask = bg_mask & mask_valid + else: + mask_valid = None if np.count_nonzero(bg_mask) > 100: std_ref = float(np.median(stdpix[bg_mask])) else: std_ref = float(np.median(stdpix)) contaminated = bright_mask & (stdpix > 3*max(std_ref, 1)) + # Also exclude masked-out pixels from the contamination set + if mask_valid is not None: + contaminated = contaminated & mask_valid if np.any(contaminated): - compressed = self.compressed.copy() + compressed = np.array(self.compressed) + compressed[3] = self.compressed[3].copy() replacement = max(1, int(round(min(std_ref, 255)))) compressed[3][contaminated] = np.uint8(replacement) else: From 2d8994f19e663ac5e6bdab0fccedcfc92371be20 Mon Sep 17 00:00:00 2001 From: Alex Tudorica Date: Thu, 9 Jul 2026 22:48:15 +0200 Subject: [PATCH 3/3] fix: address review on fireball stdpixel decontamination - Fix extractor crash on stations with detection_binning_factor > 1: the binned mask no longer matches the capture-resolution frames, so guard the mask by shape (mirrors maskImage()) instead of raising a broadcast error in findPoints(). - Extract the decontamination into a testable decontaminateStdpixel() helper; compute the brightness percentile over mask-valid pixels only (executeAll fills masked-out maxpixel with the image mean, which otherwise skews it on heavily-masked cameras). - Drop the redundant second copy of the stdpixel plane and correct the comment (np.array() already deep-copies all 4 planes). - Name the tuning constants (BRIGHT_PERCENTILE, CONTAMINATION_STD_FACTOR, MIN_BACKGROUND_PIXELS) at module level. - Correct the clip-254 comment and the .config/ConfigReader wording: the recall loss was unsaturated trail pixels clipped to 255 (fixed by decontamination), not saturated pixels; the velocity scale comes from the nominal [Capture] FOV, not the platepar. - Add Tests/test_fireball_decontamination.py: decontamination (clean no-op, contaminated trail, binned-mask shape guard, mask exclusion, replacement clamp) and the deg/s -> px/frame velocity conversion. --- .config | 10 +- RMS/ConfigReader.py | 9 +- RMS/Routines/Grouping3Dcy.pyx | 7 +- RMS/VideoExtraction.py | 122 +++++++++---- Tests/test_fireball_decontamination.py | 234 +++++++++++++++++++++++++ 5 files changed, 342 insertions(+), 40 deletions(-) create mode 100644 Tests/test_fireball_decontamination.py diff --git a/.config b/.config index 9891a1940..efb9969eb 100644 --- a/.config +++ b/.config @@ -517,10 +517,12 @@ point_ratio_threshold: 0.7 max_lines: 5 ; Maximum angular velocity for fireball candidates in deg/s. Converted to a per-sensor -; px/frame cap at runtime using fps and the platepar-derived deg/pixel scale, so the same -; value transfers between cameras with different FOV/resolution. Default is wider than the -; meteor detector's ang_vel_max (35 deg/s) because the fireball detector is recall-first -; and the cap exists mostly to reject birds/planes, not to constrain real fireballs. +; px/frame cap at runtime using fps and the deg/pixel scale derived from the nominal FOV +; (fov_w/fov_h in the [Capture] section above), so the same value transfers between cameras +; with different FOV/resolution. Keep those nominal FOV values correct for the cap to be +; meaningful. Default is wider than the meteor detector's ang_vel_max (35 deg/s) because the +; fireball detector is recall-first and the cap exists mostly to reject birds/planes, not to +; constrain real fireballs. fireball_max_ang_vel: 60.0 diff --git a/RMS/ConfigReader.py b/RMS/ConfigReader.py index 676177383..33f5a7c90 100644 --- a/RMS/ConfigReader.py +++ b/RMS/ConfigReader.py @@ -561,10 +561,11 @@ def __init__(self): self.point_ratio_threshold = 0.7# ratio of how many points must be close to the line before considering searching for another line self.max_lines = 5 # maximum number of lines # Maximum angular velocity for fireball candidates in deg/s. Converted to a per-sensor - # px/frame cap at runtime using fps and the platepar-derived deg/pixel scale, so the same - # value transfers between cameras with different FOV/resolution. Default is wider than - # the meteor detector's ang_vel_max (35 deg/s) because the fireball detector is recall- - # first and the cap exists mostly to reject birds/planes, not to constrain real fireballs. + # px/frame cap at runtime using fps and the deg/pixel scale derived from the nominal FOV + # (fov_w/fov_h in the [Capture] section), so the same value transfers between cameras with + # different FOV/resolution. Default is wider than the meteor detector's ang_vel_max + # (35 deg/s) because the fireball detector is recall-first and the cap exists mostly to + # reject birds/planes, not to constrain real fireballs. self.fireball_max_ang_vel = 60.0 ##### MeteorDetection diff --git a/RMS/Routines/Grouping3Dcy.pyx b/RMS/Routines/Grouping3Dcy.pyx index 8ce148bf2..849bbf59a 100644 --- a/RMS/Routines/Grouping3Dcy.pyx +++ b/RMS/Routines/Grouping3Dcy.pyx @@ -487,7 +487,12 @@ def thresholdAndSubsample(np.ndarray[UINT8_TYPE_t, ndim=3] frames, \ # Compute the threshold limit avg_std = int(float(compressed[2, y, x]) + k1*float(compressed[3, y, x])) + j1 - # Clip threshold to 254 so that pixels at exactly max_val=254 can still pass + # The threshold can exceed the uint8 range and must be clipped. With the old + # clip at 255 a saturated pixel (max_val == 255) still passed (255 >= 255) -- + # the real recall loss was unsaturated trail pixels whose threshold clipped to + # 255, and that is fixed upstream by stdpixel decontamination in + # VideoExtraction.findPoints. Clipping to 254 additionally lets max_val == 254 + # pixels through: a minor, harmless extra rescue. if avg_std > 254: avg_std = 254 diff --git a/RMS/VideoExtraction.py b/RMS/VideoExtraction.py index 2b21f4da2..169b8cdd0 100644 --- a/RMS/VideoExtraction.py +++ b/RMS/VideoExtraction.py @@ -35,6 +35,93 @@ log = getLogger("rmslogger") +# --- Fireball stdpixel decontamination tuning constants --- +# A bright fireball inflates stdpixel along its trail within the 256-frame FF block. +# The fireball threshold (avepixel + k1*stdpixel + j1) then exceeds the uint8 range and +# gets clipped, so thresholdAndSubsample keeps only the near-saturated trail core and +# drops the fainter trail edges -- exactly where recall is lost. These constants control +# how that self-contamination is detected and neutralised before thresholding. +BRIGHT_PERCENTILE = 90 # maxpixel percentile above which a pixel is a bright/trail candidate +CONTAMINATION_STD_FACTOR = 3.0 # stdpixel above this * background median counts as contaminated +MIN_BACKGROUND_PIXELS = 100 # need at least this many background pixels to trust their median + + +def decontaminateStdpixel(compressed, mask_img=None): + """ Replace stdpixel at fireball-contaminated pixels with the background median. + + A bright fireball inflates stdpixel along its trail, pushing the detection threshold + (avepixel + k1*stdpixel + j1) past the uint8 ceiling where it is clipped. + thresholdAndSubsample then keeps only the near-saturated trail core and drops the + fainter edges. Replacing the inflated stdpixel with the frame's background median lets + the whole trail threshold normally, which is the change that actually recovers recall. + + Pure numpy, no Cython, Python 2/Pi compatible. The 4-plane compressed array is copied + only when contamination is actually present, so the common (no-fireball) path is + zero-copy. + + Arguments: + compressed: [ndarray] (4, H, W) FTP-compressed frames + (maxpixel, maxframe, avepixel, stdpixel). + + Keyword arguments: + mask_img: [ndarray or None] mask image (H, W); pixels where it is 0 are excluded + from the background and contamination sets. Ignored when its shape does not + match the compressed frames (e.g. a binned mask on a full-resolution + extractor), mirroring maskImage()'s silent shape guard. + + Return: + [ndarray] compressed array with stdpixel decontaminated. The input array is + returned unchanged (same object) when no contamination is detected. + """ + + maxpix = compressed[0] + stdpix = compressed[3] + + # Only trust the mask when it matches the frame resolution. Extractor.__init__ bins + # the mask via binImageCalibration(), but the extractor runs on capture-resolution + # frames, so on stations with detection_binning_factor > 1 the shapes disagree. Skip + # the mask then instead of crashing on the boolean-and below (mirrors maskImage()). + if (mask_img is not None) and (mask_img.shape == maxpix.shape): + mask_valid = mask_img > 0 + else: + mask_valid = None + + # Brightness cut for the background sample. executeAll() fills masked-out pixels of + # maxpixel with the image mean before this runs, so restrict the percentile to + # mask-valid pixels to avoid skewing it on heavily-masked cameras. + maxpix_valid = maxpix[mask_valid] if mask_valid is not None else maxpix.ravel() + if maxpix_valid.size == 0: + return compressed + + # O(n) partition instead of an O(n log n) full sort for the percentile + n_bright = max(1, (maxpix_valid.size*(100 - BRIGHT_PERCENTILE))//100) + bright_threshold = np.partition(maxpix_valid, -n_bright)[-n_bright] + + bright_mask = maxpix >= bright_threshold + bg_mask = ~bright_mask + if mask_valid is not None: + bg_mask &= mask_valid + + if np.count_nonzero(bg_mask) > MIN_BACKGROUND_PIXELS: + std_ref = float(np.median(stdpix[bg_mask])) + else: + std_ref = float(np.median(stdpix)) + + contaminated = bright_mask & (stdpix > CONTAMINATION_STD_FACTOR*max(std_ref, 1)) + if mask_valid is not None: + contaminated &= mask_valid + + if not np.any(contaminated): + return compressed + + # np.array() deep-copies all 4 planes; only plane 3 (stdpixel) is then modified. + compressed = np.array(compressed) + replacement = max(1, int(round(min(std_ref, 255)))) + compressed[3][contaminated] = np.uint8(replacement) + + return compressed + + class Extractor(Process): """ Detects fireballs and brighter meteors on the FF files, and extracts raw frames while they are still in memory. @@ -70,37 +157,10 @@ def findPoints(self): (y, x, z): [tuple] Coordinates of points that form the fireball, where Z is the frame number. """ - # Decontaminate stdpixel before thresholding: a bright fireball inflates stdpixel - # along its trail, pushing the threshold (avg + k1*std + j1) above 255. Replace - # stdpixel at contaminated pixels with the background median. Only the stdpixel - # layer is copied when contamination is present, so the no-op path is zero-cost. - maxpix = self.compressed[0] - stdpix = self.compressed[3] - # O(n) partition-based percentile instead of O(n log n) full sort - bright_threshold = np.partition(maxpix.ravel(), -maxpix.size // 10)[-maxpix.size // 10] - bright_mask = maxpix >= bright_threshold - bg_mask = ~bright_mask - # Exclude masked-out pixels (camera borders/obstructions) from the background sample - if self.mask is not None: - mask_valid = self.mask.img > 0 - bg_mask = bg_mask & mask_valid - else: - mask_valid = None - if np.count_nonzero(bg_mask) > 100: - std_ref = float(np.median(stdpix[bg_mask])) - else: - std_ref = float(np.median(stdpix)) - contaminated = bright_mask & (stdpix > 3*max(std_ref, 1)) - # Also exclude masked-out pixels from the contamination set - if mask_valid is not None: - contaminated = contaminated & mask_valid - if np.any(contaminated): - compressed = np.array(self.compressed) - compressed[3] = self.compressed[3].copy() - replacement = max(1, int(round(min(std_ref, 255)))) - compressed[3][contaminated] = np.uint8(replacement) - else: - compressed = self.compressed + # Decontaminate stdpixel before thresholding so bright fireball trails are not + # suppressed by their own inflated stdpixel (see decontaminateStdpixel). + mask_img = self.mask.img if self.mask is not None else None + compressed = decontaminateStdpixel(self.compressed, mask_img) # Threshold and subsample frames length, x, y, z = Grouping3D.thresholdAndSubsample(self.frames, compressed, \ diff --git a/Tests/test_fireball_decontamination.py b/Tests/test_fireball_decontamination.py new file mode 100644 index 000000000..ab3968e30 --- /dev/null +++ b/Tests/test_fireball_decontamination.py @@ -0,0 +1,234 @@ +""" Unit tests for the fireball-detector recall improvements: + + - RMS.VideoExtraction.decontaminateStdpixel: the stdpixel decontamination that + stops bright fireball trails from suppressing themselves through the threshold. + - RMS.Routines.Grouping3D.findCoefficients: the FOV-aware deg/s -> px/frame + velocity cap that replaced the hard-coded "total > 2" filter. + + These are pure-numpy / pure-python paths, so no capture pipeline, Cython threshold + or Extractor process is needed. +""" + +import pytest + +np = pytest.importorskip("numpy") + +from RMS.VideoExtraction import ( + decontaminateStdpixel, + BRIGHT_PERCENTILE, + CONTAMINATION_STD_FACTOR, + MIN_BACKGROUND_PIXELS, +) +from RMS.Routines.Grouping3D import findCoefficients + + +# --------------------------------------------------------------------------- +# Helpers +# --------------------------------------------------------------------------- + +def _make_compressed(maxpix, stdpix): + """ Build a (4, H, W) uint8 FTP-compressed array from maxpixel/stdpixel planes. + maxframe (plane 1) and avepixel (plane 2) are filled with zeros; the + decontamination only reads planes 0 (maxpixel) and 3 (stdpixel). """ + + h, w = maxpix.shape + compressed = np.zeros((4, h, w), dtype=np.uint8) + compressed[0] = maxpix + compressed[3] = stdpix + return compressed + + +def _line(p1, p2, first_frame, last_frame): + """ Build a line_list entry as consumed by findCoefficients: endpoints at + indices 0/1 (each a (y, x, z) triple) and first/last frame at indices 4/5. """ + + return [p1, p2, 0, 0, first_frame, last_frame] + + +class _Cfg(object): + """ Minimal stand-in for the RMS Config object (only the attributes + findCoefficients reads). """ + + def __init__(self, **kwargs): + self.fireball_max_ang_vel = 60.0 + self.fps = 25.0 + self.fov_h = 35.0 + self.fov_w = 64.0 + self.height = 720 + self.width = 1280 + self.f = 16 + for key, val in kwargs.items(): + setattr(self, key, val) + + +# --------------------------------------------------------------------------- +# decontaminateStdpixel +# --------------------------------------------------------------------------- + +class TestDecontaminateStdpixel: + + def test_clean_image_is_noop(self): + """ With uniform stdpixel there is no contamination, so the same array object + is returned untouched (zero-copy fast path). """ + + maxpix = np.full((100, 100), 30, dtype=np.uint8) + maxpix[:20, :] = 60 # a bright but clean region (low, uniform std) + stdpix = np.full((100, 100), 10, dtype=np.uint8) + compressed = _make_compressed(maxpix, stdpix) + + result = decontaminateStdpixel(compressed) + + assert result is compressed + assert np.array_equal(result[3], stdpix) + + def test_contaminated_trail_is_replaced(self): + """ A bright trail with inflated stdpixel gets its stdpixel replaced by the + background median; every other plane and pixel is left untouched, and the + input array is not mutated. """ + + maxpix = np.full((100, 100), 20, dtype=np.uint8) + maxpix[:20, :] = 40 # bright-but-clean band (not contaminated) + maxpix[50, 25:75] = 255 # the fireball trail + stdpix = np.full((100, 100), 12, dtype=np.uint8) + stdpix[50, 25:75] = 200 # self-contaminated stdpixel on the trail + compressed = _make_compressed(maxpix, stdpix) + + result = decontaminateStdpixel(compressed) + + # A copy is returned when contamination is present + assert result is not compressed + + # Trail stdpixel is pulled down to the background median (12) + assert np.all(result[3][50, 25:75] == 12) + + # The clean bright band and background are unchanged + assert np.all(result[3][:20, :] == 12) + assert np.all(result[3][60:, :] == 12) + + # maxpixel plane is untouched + assert np.array_equal(result[0], compressed[0]) + + # The original array was not mutated in place + assert np.all(compressed[3][50, 25:75] == 200) + + def test_binned_mask_shape_mismatch_does_not_crash(self): + """ Regression for the blocking bug: on stations with detection_binning_factor > 1 + the mask is binned to a smaller shape than the capture-resolution frames. + The decontamination must skip the mask instead of raising a broadcast error. """ + + maxpix = np.full((100, 100), 20, dtype=np.uint8) + maxpix[50, 25:75] = 255 + stdpix = np.full((100, 100), 12, dtype=np.uint8) + stdpix[50, 25:75] = 200 + compressed = _make_compressed(maxpix, stdpix) + + # Half-resolution mask, as produced by binImageCalibration() at binning factor 2 + binned_mask = np.full((50, 50), 255, dtype=np.uint8) + + # Must not raise, and must still decontaminate as if no mask were present + result = decontaminateStdpixel(compressed, mask_img=binned_mask) + + assert result is not compressed + assert np.all(result[3][50, 25:75] == 12) + + def test_matching_mask_excludes_masked_region(self): + """ Pixels where the mask is 0 (camera borders/obstructions) must be excluded + from the contamination set, so a contaminated pixel under the mask keeps its + original stdpixel. """ + + maxpix = np.full((100, 100), 20, dtype=np.uint8) + maxpix[50, 25:75] = 255 # trail in the visible region + maxpix[70, 25:75] = 255 # trail in the masked-out region + stdpix = np.full((100, 100), 12, dtype=np.uint8) + stdpix[50, 25:75] = 200 + stdpix[70, 25:75] = 200 + compressed = _make_compressed(maxpix, stdpix) + + mask_img = np.full((100, 100), 255, dtype=np.uint8) + mask_img[60:, :] = 0 # mask out the bottom 40 rows + + result = decontaminateStdpixel(compressed, mask_img=mask_img) + + # Visible trail is decontaminated + assert np.all(result[3][50, 25:75] == 12) + + # Masked-out trail is left untouched + assert np.all(result[3][70, 25:75] == 200) + + def test_replacement_is_clamped_to_at_least_one(self): + """ The replacement value is always a valid uint8 >= 1 even when the background + median is 0, so a decontaminated pixel never gets a zero threshold offset. """ + + maxpix = np.full((100, 100), 20, dtype=np.uint8) + maxpix[50, 25:75] = 255 + stdpix = np.zeros((100, 100), dtype=np.uint8) # background median std == 0 + stdpix[50, 25:75] = 200 + compressed = _make_compressed(maxpix, stdpix) + + result = decontaminateStdpixel(compressed) + + assert np.all(result[3][50, 25:75] == 1) + + def test_tuning_constants_have_expected_defaults(self): + """ Pin the documented tuning constants so a future edit is a conscious choice. """ + + assert BRIGHT_PERCENTILE == 90 + assert CONTAMINATION_STD_FACTOR == 3.0 + assert MIN_BACKGROUND_PIXELS == 100 + + +# --------------------------------------------------------------------------- +# findCoefficients velocity cap +# --------------------------------------------------------------------------- + +class TestFindCoefficientsVelocityCap: + + def _slope_line(self, velocity_px): + """ A one-frame line whose total subsampled speed equals velocity_px: + dz = 1, dx = velocity_px, dy = 0 -> total = sqrt(v^2 + 0) = v. """ + + p1 = (0.0, 0.0, 0.0) + p2 = (0.0, float(velocity_px), 1.0) + return _line(p1, p2, first_frame=0, last_frame=1) + + def test_fallback_cap_without_config(self): + """ With no config the legacy ~2.0 px/frame cap applies. """ + + kept = findCoefficients([self._slope_line(1.5)], config=None) + dropped = findCoefficients([self._slope_line(2.5)], config=None) + + assert len(kept) == 1 + assert len(dropped) == 0 + + def test_default_config_widens_cap(self): + """ At the default 720p/25 fps/f=16 config, 60 deg/s converts to ~3.04 + subsampled px/frame -- wider than the legacy 2.0, so a 3.0 px/frame line is + kept while a 3.1 px/frame line is dropped. """ + + cfg = _Cfg() + + assert len(findCoefficients([self._slope_line(3.0)], config=cfg)) == 1 + assert len(findCoefficients([self._slope_line(3.1)], config=cfg)) == 0 + + def test_line_faster_than_cap_would_be_kept_under_default(self): + """ A 2.5 px/frame line is dropped by the fallback but kept by the default + config -- the exact behavior that rescues slower fireballs. """ + + assert len(findCoefficients([self._slope_line(2.5)], config=None)) == 0 + assert len(findCoefficients([self._slope_line(2.5)], config=_Cfg())) == 1 + + def test_missing_config_attribute_falls_back(self): + """ If a required attribute is missing/zero, the conversion is skipped and the + conservative 2.0 px/frame fallback is used. """ + + cfg = _Cfg(fov_h=None) + + assert len(findCoefficients([self._slope_line(2.5)], config=cfg)) == 0 + + def test_points_on_same_frame_are_skipped(self): + """ Degenerate lines with both endpoints on the same frame are ignored rather + than dividing by zero. """ + + same_frame = _line((0.0, 0.0, 5.0), (0.0, 1.0, 5.0), first_frame=5, last_frame=5) + + assert findCoefficients([same_frame], config=_Cfg()) == []