1use num_traits::clamp;
2
3use crate::{ImageBuffer, Pixel, Primitive};
4
5/// Approximation of Gaussian blur after
6/// Kovesi, P.: Fast Almost-Gaussian Filtering The Australian Pattern
7/// Recognition Society Conference: DICTA 2010. December 2010. Sydney.
8/// This method assumes alpha pre-multiplication for images that contain non-constant alpha.
9#[must_use]
10pub fn fast_blur<P: Pixel>(
11 image_buffer: &ImageBuffer<P, Vec<P::Subpixel>>,
12 sigma: f32,
13) -> ImageBuffer<P, Vec<P::Subpixel>> {
14 let (width, height) = image_buffer.dimensions();
15
16 if width == 0 || height == 0 {
17 return image_buffer.clone();
18 }
19 let mut samples = image_buffer.as_flat_samples().samples.to_vec();
20 let num_passes = 3;
21
22 let boxes = boxes_for_gauss(sigma, num_passes);
23
24 for radius in boxes.iter().take(num_passes) {
25 let horizontally_blurred_transposed = horizontal_fast_blur_half::<P::Subpixel>(
26 &samples,
27 width as usize,
28 height as usize,
29 (*radius - 1) / 2,
30 P::CHANNEL_COUNT as usize,
31 );
32 samples = horizontal_fast_blur_half::<P::Subpixel>(
33 &horizontally_blurred_transposed,
34 height as usize,
35 width as usize,
36 (*radius - 1) / 2,
37 P::CHANNEL_COUNT as usize,
38 );
39 }
40 ImageBuffer::from_raw(width, height, samples).unwrap()
41}
42
43fn boxes_for_gauss(sigma: f32, n: usize) -> Vec<usize> {
44 let w_ideal: f32 = f32::sqrt((12.0 * sigma.powi(2) / (n as f32)) + 1.0);
45 let mut w_l: f32 = w_ideal.floor();
46 if w_l % 2.0 == 0.0 {
47 w_l -= 1.0;
48 };
49 let w_u: f32 = w_l + 2.0;
50
51 let m_ideal: f32 = 0.25 * (n as f32) * (w_l + 3.0) - 3.0 * sigma.powi(2) * (w_l + 1.0).recip();
52
53 let m: usize = f32::round(self:m_ideal) as usize;
54
55 (0..n)
56 .map(|i: usize| if i < m { w_l as usize } else { w_u as usize })
57 .collect::<Vec<_>>()
58}
59
60fn channel_idx(channel: usize, idx: usize, channel_num: usize) -> usize {
61 channel_num * idx + channel
62}
63
64fn horizontal_fast_blur_half<P: Primitive>(
65 samples: &[P],
66 width: usize,
67 height: usize,
68 r: usize,
69 channel_num: usize,
70) -> Vec<P> {
71 let channel_size = width * height;
72
73 let mut out_samples = vec![P::from(0).unwrap(); channel_size * channel_num];
74 let mut vals = vec![0.0; channel_num];
75
76 let min_value = P::DEFAULT_MIN_VALUE.to_f32().unwrap();
77 let max_value = P::DEFAULT_MAX_VALUE.to_f32().unwrap();
78
79 for row in 0..height {
80 for (channel, value) in vals.iter_mut().enumerate().take(channel_num) {
81 *value = ((-(r as isize))..(r + 1) as isize)
82 .map(|x| {
83 extended_f(
84 samples,
85 width,
86 height,
87 x,
88 row as isize,
89 channel,
90 channel_num,
91 )
92 .to_f32()
93 .unwrap_or(0.0)
94 })
95 .sum();
96 }
97
98 for column in 0..width {
99 for (channel, channel_val) in vals.iter_mut().enumerate() {
100 let val = *channel_val / (2.0 * r as f32 + 1.0);
101 let val = clamp(val, min_value, max_value);
102 let val = P::from(val).unwrap();
103
104 let destination_row = column;
105 let destination_column = row;
106 let destination_sample_index = channel_idx(
107 channel,
108 destination_column + destination_row * height,
109 channel_num,
110 );
111 out_samples[destination_sample_index] = val;
112 *channel_val = *channel_val
113 - extended_f(
114 samples,
115 width,
116 height,
117 column as isize - r as isize,
118 row as isize,
119 channel,
120 channel_num,
121 )
122 .to_f32()
123 .unwrap_or(0.0)
124 + extended_f(
125 samples,
126 width,
127 height,
128 { column + r + 1 } as isize,
129 row as isize,
130 channel,
131 channel_num,
132 )
133 .to_f32()
134 .unwrap_or(0.0);
135 }
136 }
137 }
138
139 out_samples
140}
141
142fn extended_f<P: Primitive>(
143 samples: &[P],
144 width: usize,
145 height: usize,
146 x: isize,
147 y: isize,
148 channel: usize,
149 channel_num: usize,
150) -> P {
151 let x: usize = clamp(input:x, min:0, max:width as isize - 1) as usize;
152 let y: usize = clamp(input:y, min:0, max:height as isize - 1) as usize;
153 samples[channel_idx(channel, idx:y * width + x, channel_num)]
154}
155