Реализация обнаружения линий преобразования Хафа в 2D-системе координат

Я хочу реализовать обнаружение линий в простой системе координат. Я примерно следил за статьей о том, как реализовать преобразование Хафа, но результаты, которые я получаю, весьма далеки от того, что я хочу.

Учитывая такую ​​матрицу 3 x 3:

X X X
X X X
- - -

Я хочу обнаружить строку, начинающуюся с 0,0, до 2,0. Я представляю систему координат как простой массив кортежей, первый элемент в кортеже - x, второй - y, третий - тип точки (холст или линия).

Я думал, что с помощью Hough обнаружить линию будет относительно легко, потому что определение границ - это, по сути, просто двоичное решение: либо кортеж имеет тип line, либо нет.

Я реализовал на Rust следующую программу:

use std::f32;

extern crate nalgebra as na;
use na::DMatrix;

#[derive(Debug, PartialEq, Clone)]
enum Representation {
   Canvas,
   Line,
}

fn main () {
    let image_width = 3;
    let image_height = 3;

    let grid = vec![
        (0, 0, Representation::Line), (1, 0, Representation::Line), (2, 0, Representation::Line),
        (0, 1, Representation::Canvas), (1, 1, Representation::Canvas), (2, 1, Representation::Canvas),
        (0, 2, Representation::Canvas), (1, 2, Representation::Canvas), (2, 2, Representation::Canvas),
    ];

    //let tmp:f32 = (image_width as f32 * image_width as f32) + (image_height as f32 * image_height as f32);
    let max_line_length = 3;
    let mut accumulator = DMatrix::from_element(180, max_line_length as usize, 0);

    for y in 0..image_height {
        for x in 0..image_width {
            let coords_index = (y * image_width) + x;
            let coords = grid.get(coords_index as usize).unwrap();

            // check if coords is an edge
            if coords.2 == Representation::Line {
                for angle in 0..180 {
                    let r = (x as f32) * (angle as f32).cos() + (y as f32) * (angle as f32).sin();
                    let r_scaled = scale_between(r, 0.0, 2.0, -2.0, 2.0).round() as u32;

                    accumulator[(angle as usize, r_scaled as usize)] += 1;
                }
            }
        }
    }

    let threshold = 3;

    // z = angle
    for z in 0..180 {
        for r in 0..3 {
            let val = accumulator[(z as usize, r as usize)];

            if val < threshold {
                continue;
            }

            let px = (r as f32) * (z as f32).cos();
            let py = (r as f32) * (z as f32).sin();

            let p1_px = px + (max_line_length as f32) * (z as f32).cos();
            let p1_py = py + (max_line_length as f32) * (z as f32).sin();

            let p2_px = px - (max_line_length as f32) * (z as f32).cos();
            let p2_py = px - (max_line_length as f32) * (z as f32).cos();

            println!("Found lines from {}/{} to {}/{}", p1_px.ceil(), p1_py.ceil(), p2_px.ceil(), p2_py.ceil());
        }
    }
}

fn scale_between(unscaled_num: f32, min_allowed: f32, max_allowed: f32, min: f32, max: f32) -> f32 {
    (max_allowed - min_allowed) * (unscaled_num - min) / (max - min) + min_allowed
}

Результат примерно такой:

Found lines from -1/4 to 1/1
Found lines from 2/4 to 0/0
Found lines from 2/-3 to 0/0
Found lines from -1/4 to 1/1
Found lines from 1/-3 to 0/0
Found lines from 0/4 to 1/1
...

На самом деле это довольно много, учитывая, что я хочу обнаружить только одну строку. Моя реализация явно неверна, но я не знаю, где искать, мой maths-fu недостаточно высок для дальнейшей отладки.

Я думаю, что первая часть, собственно преобразование Хафа, кажется правильной, потому что в связанной статье говорится:

for each image point p 
{
  if (p is part of an edge)
  {
    for each possible angle
    {
     r = x * cos(angle) + y * sin(angle);
     houghMatrix[angle][r]++;
    }
  }
}

Я застрял в отображении и фильтрации, что согласно статье:

1. Каждая точка в пространстве Хафа задается углом a и расстоянием r. Используя эти значения, можно вычислить одну точку p (x, y) линии как px = r * cos (угол) py = r * sin (угол).

2. Максимальная длина строки ограничена sqrt (imagewidth2 + imageheight2).

3. Точка p, угол a линии и максимальная длина линии maxLength могут использоваться для вычисления двух других точек линии. Максимальная длина здесь гарантирует, что обе вычисляемые точки лежат за пределами фактического изображения, в результате чего, если линия проведена между этими двумя точками, линия в любом случае идет от границы изображения к границе изображения и никогда не обрезается. где-то внутри изображения.

4. Эти две точки p1 и p2 вычисляются по формуле: p1_x = px + maxLength * cos (угол); p1_y = py + maxLength * sin (угол); p2_x = px - maxLength * cos (угол); p2_y = py - maxLength * sin (угол);

5. ...

ИЗМЕНИТЬ

Обновленная версия, которая учитывает размер изображения, как предложено @RaymoAisla

use std::f32;

extern crate nalgebra as na;
use na::DMatrix;

fn main () {
    let image_width = 3;
    let image_height = 3;

    let mut grid = DMatrix::from_element(image_width as usize, image_height as usize, 0);
    grid[(0, 0)] = 1;
    grid[(1, 0)] = 1;
    grid[(2, 0)] = 1;

    let accu_width = 7;
    let accu_height = 3;
    let max_line_length = 3;

    let mut accumulator = DMatrix::from_element(accu_width as usize, accu_height as usize, 0);


    for y in 0..image_height {
        for x in 0..image_width {
            let coords = (x, y);
            let is_edge = grid[coords] == 1;

            if !is_edge {
                continue;
            }

            for i in 0..7 {
                let angle = i * 30;

                let r = (x as f32) * (angle as f32).cos() + (y as f32) * (angle as f32).sin();
                let r_scaled = scale_between(r, 0.0, 2.0, -2.0, 2.0).round() as u32;

                accumulator[(i as usize, r_scaled as usize)] += 1;

                println!("angle: {}, r: {}, r_scaled: {}", angle, r, r_scaled);
            }
        }
    }

    let threshold = 3;

    // z = angle index
    for z in 0..7 {
        for r in 0..3 {
            let val = accumulator[(z as usize, r as usize)];

            if val < threshold {
                continue;
            }

            let px = (r as f32) * (z as f32).cos();
            let py = (r as f32) * (z as f32).sin();

            let p1_px = px + (max_line_length as f32) * (z as f32).cos();
            let p1_py = py + (max_line_length as f32) * (z as f32).sin();

            let p2_px = px - (max_line_length as f32) * (z as f32).cos();
            let p2_py = px - (max_line_length as f32) * (z as f32).cos();

            println!("Found lines from {}/{} to {}/{} - val: {}", p1_px.ceil(), p1_py.ceil(), p2_px.ceil(), p2_py.ceil(), val);
        }
    }
}

fn scale_between(unscaled_num: f32, min_allowed: f32, max_allowed: f32, min: f32, max: f32) -> f32 {
    (max_allowed - min_allowed) * (unscaled_num - min) / (max - min) + min_allowed
}

Сообщаемый результат теперь:

angle: 0, r: 0, r_scaled: 1
angle: 30, r: 0, r_scaled: 1
angle: 60, r: 0, r_scaled: 1
angle: 90, r: 0, r_scaled: 1
angle: 120, r: 0, r_scaled: 1
angle: 150, r: 0, r_scaled: 1
angle: 180, r: 0, r_scaled: 1
...
Found lines from 3/4 to -1/-1
Found lines from -3/1 to 2/2 

Я нарисовал линии в системе координат, они очень далеко от линии, которую я ожидал. Интересно, конвертация обратно в баллы все еще отключена.