[PATCH 4/5] libipa: histogram: Fix interQuantileMean() for small ranges

[Laurent Pinchart]

Hi Stefan,

Thank you for the patch.

On Mon, Mar 24, 2025 at 06:07:39PM +0100, Stefan Klug wrote:
> The interQuantileMean() is supposed to return a weighted mean value
> between two quantiles. This works for reasonably fine histograms, but
> fails for coarse histograms and small quantile ranges because the weight
> is always taken from the lower border of the bin.
> 
> Fix that by rewriting the algorithm to calculate a lower and upper bound
> for every (partial) bin that goes into the mean calculation and weight
> the bins by the middle of these bounds.
> 
> Signed-off-by: Stefan Klug <stefan.klug at ideasonboard.com>
> ---
>  src/ipa/libipa/histogram.cpp | 20 +++++++++++---------
>  1 file changed, 11 insertions(+), 9 deletions(-)
> 
> diff --git a/src/ipa/libipa/histogram.cpp b/src/ipa/libipa/histogram.cpp
> index c19a4cbbf3cd..31f017af3458 100644
> --- a/src/ipa/libipa/histogram.cpp
> +++ b/src/ipa/libipa/histogram.cpp
> @@ -153,22 +153,24 @@ double Histogram::interQuantileMean(double lowQuantile, double highQuantile) con
>  	double lowPoint = quantile(lowQuantile);
>  	/* Proportion of pixels which lies below highQuantile */
>  	double highPoint = quantile(highQuantile, static_cast<uint32_t>(lowPoint));

Those two variables can now be const. You can write

	ASSERT(highQuantile > lowQuantile);
	
	/* Proportion of pixels which lies below lowQuantile and highQuantile. */
	const double lowPoint = quantile(lowQuantile);
        const double highPoint = quantile(highQuantile, static_cast<uint32_t>(lowPoint));

> -	double sumBinFreq = 0, cumulFreq = 0;
> +	double sumBinFreq = 0;
> +	double cumulFreq = 0;
> +

Let's document the algorithm (and see if I understand it correctly :-)).

	/*
	 * Calculate the mean pixel value between the low and high points by
	 * summing all the pixels between the two points, and dividing the sum
	 * by the number of pixels. Given the discrete nature of the histogram
	 * data, the sum of the pixels is approximated by accummulating the
	 * product of the bin values (calculated as the mid point of the bin) by
	 * the number of pixels they contain, for each bin in the internal.
	 */

> +	for (int bin = std::floor(lowPoint); bin < std::ceil(highPoint); bin++) {

It looks like bin can be unsigned.

> +		double lowBound = std::max(static_cast<double>(bin), lowPoint);

I think you can also write

		double lowBound = std::max<double>(bin, lowPoint);

Same for the next line. Up to you. Oh, and you can make them const too.

> +		double highBound = std::min(static_cast<double>(bin + 1), highPoint);

If I understand the code correctly, this is only meaningful for the
first and last iterations. I can't easily find a better construct that
wouldn't need to be run for each iteration, so this seems fine.

>  
> -	for (double p_next = floor(lowPoint) + 1.0;
> -	     p_next <= ceil(highPoint);
> -	     lowPoint = p_next, p_next += 1.0) {
> -		int bin = floor(lowPoint);
>  		double freq = (cumulative_[bin + 1] - cumulative_[bin])
> -			* (std::min(p_next, highPoint) - lowPoint);
> +			* (highBound - lowBound);

	 	/*
		 * The low and high quantile may not lie at bin boundaries, so
		 * the first and last bins need to be weighted accordingly. The
		 * best available approximation is to multiply the number of
		 * pixels by the partial bin width.
		 */
		const double freq = (cumulative_[bin + 1] - cumulative_[bin])
				  * (highBound - lowBound);

>  
>  		/* Accumulate weighted bin */
> -		sumBinFreq += bin * freq;
> +		sumBinFreq += 0.5 * (highBound + lowBound) * freq;

I wondered for a moment where the 0.5 came from. I think

		sumBinFreq += (highBound + lowBound) / 2 * freq;

would better reflect the intent.

> +
>  		/* Accumulate weights */
>  		cumulFreq += freq;

I wonder if we should rename sumBinFreq to sumPixelValues and numPixels.

Reviewed-by: Laurent Pinchart <laurent.pinchart at ideasonboard.com>

>  	}
> -	/* add 0.5 to give an average for bin mid-points */
> -	return sumBinFreq / cumulFreq + 0.5;
> +
> +	return sumBinFreq / cumulFreq;
>  }
>  
>  } /* namespace ipa */

-- 
Regards,

Laurent Pinchart