Skip to content
Merged
Show file tree
Hide file tree
Changes from 1 commit
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
176 changes: 100 additions & 76 deletions src/imgproc/histogram.rs
Original file line number Diff line number Diff line change
Expand Up @@ -171,6 +171,76 @@ fn map_bin(val: f32, range: &RangeSpec, hist_size: usize) -> Option<usize> {
}
}

/// Validates that `hist_size` and `ranges` describe a well-formed set of
/// histogram dimensions, shared by `calc_hist` and `calc_back_project` so
/// the two can't drift: `hist_size[d] > 0`, `Uniform(lo, hi)` satisfies
/// `lo < hi` (also rejects NaN, since `partial_cmp` returns `None` for it),
/// and `NonUniform(boundaries)` has exactly `hist_size[d] + 1` strictly
/// increasing entries.
///
/// `fn_name` is used as the error-message prefix so callers keep their own
/// identity in the message (e.g. `"calc_hist: ..."` vs
/// `"calc_back_project: ..."`).
fn validate_hist_ranges(fn_name: &str, hist_size: &[usize], ranges: &[RangeSpec]) -> Result<()> {
for (d, &sz) in hist_size.iter().enumerate() {
if sz == 0 {
cv_bail!(
tags::IMGPROC,
InvalidInput,
"{}: hist_size[{}] must be > 0 (got 0)",
fn_name,
d
);
}
match &ranges[d] {
RangeSpec::Uniform(lo, hi) => {
if lo.partial_cmp(hi) != Some(Ordering::Less) {
cv_bail!(
tags::IMGPROC,
InvalidInput,
"{}: Uniform range [{}] must satisfy lo < hi (got {} >= {})",
fn_name,
d,
lo,
hi
);
}
}
RangeSpec::NonUniform(boundaries) => {
if boundaries.len() != sz + 1 {
Comment thread
qodo-code-review[bot] marked this conversation as resolved.
Outdated
cv_bail!(
tags::IMGPROC,
InvalidInput,
"{}: NonUniform boundaries[{}] length {} must be hist_size[{}]+1 ({})",
fn_name,
d,
boundaries.len(),
d,
sz + 1
);
}
for k in 0..boundaries.len() - 1 {
if boundaries[k].partial_cmp(&boundaries[k + 1]) != Some(Ordering::Less) {
cv_bail!(
tags::IMGPROC,
InvalidInput,
"{}: NonUniform boundaries[{}][{}] ({}) must be < boundaries[{}][{}] ({})",
fn_name,
d,
k,
boundaries[k],
d,
k + 1,
boundaries[k + 1]
);
}
}
}
}
}
Ok(())
}

// ---------------------------------------------------------------------------
// calc_hist
// ---------------------------------------------------------------------------
Expand Down Expand Up @@ -269,58 +339,7 @@ pub fn calc_hist<T: ToPrimitive + Clone + Default + Send + Sync>(
}

// Validate hist_size and ranges (prevents panics in bin mapping)
for (d, &sz) in hist_size.iter().enumerate() {
if sz == 0 {
cv_bail!(
tags::IMGPROC,
InvalidInput,
"calc_hist: hist_size[{}] must be > 0 (got 0)",
d
);
}
match &ranges[d] {
RangeSpec::Uniform(lo, hi) => {
if lo.partial_cmp(hi) != Some(Ordering::Less) {
cv_bail!(
tags::IMGPROC,
InvalidInput,
"calc_hist: Uniform range [{}] must satisfy lo < hi (got {} >= {})",
d,
lo,
hi
);
}
}
RangeSpec::NonUniform(boundaries) => {
if boundaries.len() != sz + 1 {
cv_bail!(
tags::IMGPROC,
InvalidInput,
"calc_hist: NonUniform boundaries[{}] length {} must be hist_size[{}]+1 ({})",
d,
boundaries.len(),
d,
sz + 1
);
}
for k in 0..boundaries.len() - 1 {
if boundaries[k].partial_cmp(&boundaries[k + 1]) != Some(Ordering::Less) {
cv_bail!(
tags::IMGPROC,
InvalidInput,
"calc_hist: NonUniform boundaries[{}][{}] ({}) must be < boundaries[{}][{}] ({})",
d,
k,
boundaries[k],
d,
k + 1,
boundaries[k + 1]
);
}
}
}
}
}
validate_hist_ranges("calc_hist", hist_size, ranges)?;

let total_bins: usize = hist_size.iter().product();

Expand Down Expand Up @@ -477,6 +496,7 @@ pub fn calc_back_project<T: ToPrimitive + Clone + Default + Send + Sync>(
dims
);
}
validate_hist_ranges("calc_back_project", hist_size, ranges)?;

let rows = images[0].rows;
let cols = images[0].cols;
Expand Down Expand Up @@ -1199,6 +1219,30 @@ fn clip_and_redistribute(tile_hist: &mut [i32], clip_limit: i32, hist_size: usiz
}
}

/// Computes the pair of tile indices and blend weight for CLAHE's bilinear
/// tile interpolation along one axis.
///
/// `coord_f` is the pixel's tile-space coordinate, already offset by -0.5 so
/// tile centers land on integers. Returns `(idx1, idx2, weight1, weight2)`
/// where `weight1 + weight2 == 1.0` and both weights are in `[0, 1]`.
///
/// The fraction is derived from the *unclamped* floor of `coord_f` so it
/// always lands in `[0, 1)` - a genuine convex blend. Only the indices used
/// for LUT lookup are clamped to the valid tile range; clamping the index
/// *before* computing the fraction (as an earlier version of this code did)
/// produces fractions outside `[0, 1]` at the first/last tile, extrapolating
/// beyond that tile's LUT instead of blending within it.
#[inline(always)]
pub(crate) fn tile_interp_weights(coord_f: f64, num_tiles: usize) -> (usize, usize, f64, f64) {
let idx1_raw = coord_f.floor() as i32;
let weight2 = coord_f - idx1_raw as f64;
let weight1 = 1.0 - weight2;
let max_idx = num_tiles as i32 - 1;
let idx1 = idx1_raw.clamp(0, max_idx) as usize;
let idx2 = (idx1_raw + 1).clamp(0, max_idx) as usize;
(idx1, idx2, weight1, weight2)
}

#[allow(clippy::too_many_arguments)]
fn interpolate_tiles_u8(
src: &Matrix<u8>,
Expand All @@ -1219,21 +1263,11 @@ fn interpolate_tiles_u8(

let process_row = |y: usize, dst_row: &mut [u8]| {
let tyf = y as f64 * inv_th - 0.5;
let ty1 = (tyf.floor() as i32).max(0);
let ty2 = (ty1 + 1).min(tiles_y as i32 - 1);
let ya = tyf - ty1 as f64;
let ya1 = 1.0 - ya;
let ty1 = ty1 as usize;
let ty2 = ty2 as usize;
let (ty1, ty2, ya1, ya) = tile_interp_weights(tyf, tiles_y);

for (x, out_pixel) in dst_row.iter_mut().enumerate() {
let txf = x as f64 * inv_tw - 0.5;
let tx1 = (txf.floor() as i32).max(0);
let tx2 = (tx1 + 1).min(tiles_x as i32 - 1);
let xa = txf - tx1 as f64;
let xa1 = 1.0 - xa;
let tx1 = tx1 as usize;
let tx2 = tx2 as usize;
let (tx1, tx2, xa1, xa) = tile_interp_weights(txf, tiles_x);

let src_val = *src.get(y, x, 0).unwrap_or(&0) as usize;
let bin = (src_val >> bit_shift).min(hist_size - 1);
Expand Down Expand Up @@ -1286,21 +1320,11 @@ fn interpolate_tiles_u16(

let process_row = |y: usize, dst_row: &mut [u16]| {
let tyf = y as f64 * inv_th - 0.5;
let ty1 = (tyf.floor() as i32).max(0);
let ty2 = (ty1 + 1).min(tiles_y as i32 - 1);
let ya = tyf - ty1 as f64;
let ya1 = 1.0 - ya;
let ty1 = ty1 as usize;
let ty2 = ty2 as usize;
let (ty1, ty2, ya1, ya) = tile_interp_weights(tyf, tiles_y);

for (x, out_pixel) in dst_row.iter_mut().enumerate() {
let txf = x as f64 * inv_tw - 0.5;
let tx1 = (txf.floor() as i32).max(0);
let tx2 = (tx1 + 1).min(tiles_x as i32 - 1);
let xa = txf - tx1 as f64;
let xa1 = 1.0 - xa;
let tx1 = tx1 as usize;
let tx2 = tx2 as usize;
let (tx1, tx2, xa1, xa) = tile_interp_weights(txf, tiles_x);

let src_val = *src.get(y, x, 0).unwrap_or(&0) as usize;
let bin = (src_val >> bit_shift).min(hist_size - 1);
Expand Down
98 changes: 98 additions & 0 deletions src/imgproc/tests.rs
Original file line number Diff line number Diff line change
Expand Up @@ -1556,6 +1556,66 @@ mod imgproc_tests {
.is_err());
}

#[test]
fn test_calc_back_project_invalid_uniform_range_error() {
// Regression test: calc_back_project used to skip the range
// validation calc_hist already had (lo < hi, NaN rejection,
// non-uniform boundary shape/ordering), silently producing a
// plausible-but-wrong projection instead of an error.
let img = Matrix::from_vec(1, 4, 1, vec![0u8, 1, 2, 3]);
let hist = Matrix::from_vec(4, 1, 1, vec![10.0, 20.0, 30.0, 40.0]);

// NaN bounds: lo.partial_cmp(hi) is None, must be rejected.
assert!(calc_back_project(
&[&img],
&[0],
&[4],
&hist,
&[RangeSpec::Uniform(f32::NAN, f32::NAN)],
1.0,
)
.is_err());

// lo >= hi.
assert!(calc_back_project(
&[&img],
&[0],
&[4],
&hist,
&[RangeSpec::Uniform(4.0, 0.0)],
1.0,
)
.is_err());
}

#[test]
fn test_calc_back_project_invalid_nonuniform_boundaries_error() {
let img = Matrix::from_vec(1, 4, 1, vec![0u8, 1, 2, 3]);
let hist = Matrix::from_vec(4, 1, 1, vec![10.0, 20.0, 30.0, 40.0]);

// Wrong length: needs hist_size[0] + 1 = 5 boundaries, only 3 given.
assert!(calc_back_project(
&[&img],
&[0],
&[4],
&hist,
&[RangeSpec::NonUniform(vec![0.0, 2.0, 4.0])],
1.0,
)
.is_err());

// Not strictly increasing.
assert!(calc_back_project(
&[&img],
&[0],
&[4],
&hist,
&[RangeSpec::NonUniform(vec![0.0, 2.0, 1.0, 3.0, 4.0])],
1.0,
)
.is_err());
}

#[test]
fn test_compare_hist_correl_identical() {
let h1 = Matrix::from_vec(4, 1, 1, vec![1.0, 2.0, 3.0, 4.0]);
Expand Down Expand Up @@ -1968,4 +2028,42 @@ mod imgproc_tests {
let dst = clahe.apply_u8(&img).unwrap();
assert_eq!(dst.data.len(), 256);
}

#[test]
fn test_clahe_tile_interp_weights_are_bounded() {
// Regression test: the interpolation fraction used to be derived
// *after* clamping the tile index, which produced weights outside
// [0, 1] (extrapolation) at the top/left edges - e.g. at y=0 with
// tile_rows=4, tyf=-0.5 used to yield ya=-0.5, ya1=1.5 instead of a
// convex blend. Check the invariant holds across the coordinate
// range that actually occurs (coord_f = pixel*inv_tile_extent - 0.5,
// for pixel in 0..src_extent), for several tile counts.
use crate::imgproc::histogram::tile_interp_weights;

for num_tiles in [1usize, 2, 3, 5] {
let tile_extent = 4usize; // pixels per tile
let inv = 1.0 / tile_extent as f64;
let src_extent = num_tiles * tile_extent;

for pixel in 0..src_extent {
let coord_f = pixel as f64 * inv - 0.5;
let (idx1, idx2, w1, w2) = tile_interp_weights(coord_f, num_tiles);

assert!(
(0.0..=1.0).contains(&w1) && (0.0..=1.0).contains(&w2),
"weights out of [0,1] at pixel={pixel}, num_tiles={num_tiles}: w1={w1}, w2={w2}"
);
assert!(
(w1 + w2 - 1.0).abs() < 1e-12,
"weights don't sum to 1 at pixel={pixel}: w1={w1}, w2={w2}"
);
assert!(idx1 < num_tiles && idx2 < num_tiles);
}
}

// The exact case from the bug report: y=0, tile_rows=4 -> tyf=-0.5.
let (idx1, idx2, w1, w2) = tile_interp_weights(-0.5, 2);
assert_eq!((idx1, idx2), (0, 0)); // both clamp to the first tile
assert!((w1 - 0.5).abs() < 1e-12 && (w2 - 0.5).abs() < 1e-12);
}
}