use crate::{BlurError, KernelShape};
pub fn lens_kernel(
shape: KernelShape,
n: f32,
m: f32,
k: f32,
rotation: f32,
) -> Result<Vec<f32>, BlurError> {
if shape.width == 0 || shape.height == 0 {
return Err(BlurError::ZeroBaseSize);
}
if shape.width.is_multiple_of(2) {
return Err(BlurError::OddKernel(shape.width));
}
if shape.height.is_multiple_of(2) {
return Err(BlurError::OddKernel(shape.height));
}
assert!(k.abs() <= 1., "Roundness contract is not satisfied");
let eps = 1f32 / shape.width as f32;
let radius = (shape.width as f32 - 1.) / 2. - 1.;
let mut new_buffer = vec![0f32; shape.width * shape.height];
for (j, row) in new_buffer.chunks_exact_mut(shape.width).enumerate() {
for (i, dst) in row.iter_mut().enumerate() {
let x = (i as f32 - 1.) / radius - 1.;
let y = (j as f32 - 1.) / radius - 1.;
let r = x.hypot(y);
let m = f32::cos((2f32 * f32::asin(k) + std::f32::consts::PI * m) / (2. * n))
/ f32::cos(
(2. * f32::asin(k * f32::cos(n * (f32::atan2(y, x) + rotation)))
+ std::f32::consts::PI * m)
/ (2. * n),
);
*dst = if m >= r + eps { 1. } else { 0. };
}
}
let sum = new_buffer.iter().sum::<f32>();
if sum != 0. {
let recip = 1. / sum;
new_buffer.iter_mut().for_each(|x| *x *= recip);
}
Ok(new_buffer)
}