Быстрый поиск ближайшего нечерного пикселя на изображении

StackOverflow https://stackoverflow.com/questions/307445

Вопрос

У меня есть 2D-изображение, беспорядочно и скудно разбросанное пикселями.
учитывая точку на изображении, мне нужно найти расстояние до ближайшего пикселя, которого нет в цвете фона (черном).
Какой самый быстрый способ сделать это?

Единственный метод, который я смог придумать, - это построение kd-дерева для пикселей.но я бы действительно хотел избежать такой дорогостоящей предварительной обработки.кроме того, кажется, что kd-дерево дает мне больше, чем мне нужно.Мне нужно только расстояние до чего-то, и меня не волнует, что это за что-то такое.

Это было полезно?

Решение

Как говорит Pyro, ищите по периметру квадрата, который вы продолжаете перемещать на один пиксель за раз от вашей исходной точки (то есть увеличивая ширину и высоту на два пикселя за раз). Когда вы попадаете на не черный пиксель, вы рассчитываете расстояние (это ваш первый дорогой расчет), а затем продолжаете поиск наружу, пока ширина вашего прямоугольника не станет в два раза больше расстояния до первой найденной точки (любые точки за пределами этого не могут быть ближе чем ваш оригинальный найденный пиксель). Сохраните все точки, не являющиеся черными, которые вы найдете во время этой части, а затем рассчитайте каждое из их расстояний, чтобы увидеть, ближе ли они к вашей исходной точке.

В идеальном случае вам нужно сделать только один дорогой расчет расстояния.

Обновление : поскольку здесь вы вычисляете межпиксельное расстояние (вместо произвольных значений с плавающей запятой точности), вы можете существенно ускорить этот алгоритм, используя предварительно рассчитанную таблицу поиска ( просто массив высотой на ширину), чтобы получить расстояние как функцию от x и y . Массив 100x100 стоит вам, по существу, 40 КБ памяти и занимает площадь 200x200 вокруг исходной точки, а также экономит затраты на дорогостоящий расчет расстояния (будь то пифагорейская или матричная алгебра) для каждого цветного пикселя, который вы найдете. Этот массив может даже быть предварительно рассчитан и встроен в ваше приложение как ресурс, чтобы сэкономить начальное время вычислений (это, вероятно, серьезное перерасход).

Обновление 2 . Кроме того, существуют способы оптимизации поиска по квадратному периметру. Ваш поиск должен начинаться в четырех точках, которые пересекают оси и перемещаются по одному пикселю за один раз к углам (у вас есть 8 движущихся точек поиска, которые могут легко создать больше проблем, чем оно того стоит, в зависимости от требований вашего приложения). Как только вы находите цветной пиксель, нет необходимости продолжать движение к углам, поскольку все остальные точки находятся дальше от начала координат.

После первого найденного пикселя вы можете дополнительно ограничить требуемую дополнительную область поиска до минимума, используя таблицу поиска, чтобы убедиться, что каждая искомая точка находится ближе, чем найденная точка (снова начиная с осей и останавливаясь, когда расстояние предел достигнут). Эта вторая оптимизация, вероятно, будет слишком дорогой, если вы будете рассчитывать каждое расстояние на лету.

Если ближайший пиксель находится в поле 200x200 (или любой другой размер подходит для ваших данных), вы будете искать только в пределах круга, ограниченного пикселем, делая только поиск и < > сравнения.

Другие советы

Лично я бы проигнорировал предложение MusiGenesis о таблице подстановки.

Вычисление расстояния между пикселями выполняется следующим образом нет дорого, особенно потому, что для этого начального теста вам не нужно фактическое расстояние, поэтому нет необходимости извлекать квадратный корень.Вы можете работать с расстоянием ^ 2, т.е.:

r^2 = dx^2 + dy^2

Кроме того, если вы выходите наружу по одному пикселю за раз, помните, что:

(n + 1)^2 = n^2 + 2n + 1

или если nx является текущим значением и бык является предыдущим значением:

    nx^2  = ox^2 + 2ox + 1
          = ox^2 + 2(nx - 1) + 1
          = ox^2 + 2nx - 1
=>  nx^2 += 2nx - 1 

Легко понять, как это работает:

1^2 =  0 + 2*1 - 1 =  1
2^2 =  1 + 2*2 - 1 =  4
3^2 =  4 + 2*3 - 1 =  9
4^2 =  9 + 2*4 - 1 = 16
5^2 = 16 + 2*5 - 1 = 25
etc...

Таким образом, на каждой итерации вам необходимо сохранять только некоторые промежуточные переменные, таким образом:

int dx2 = 0, dy2, r2;
for (dx = 1; dx < w; ++dx) {  // ignoring bounds checks
   dx2 += (dx << 1) - 1;
   dy2 = 0;
   for (dy = 1; dy < h; ++dy) {
       dy2 += (dy << 1) - 1;
       r2 = dx2 + dy2;
       // do tests here
   }
}

Tada!вычисление r ^ 2 только с битовыми сдвигами, сложениями и вычитаниями :)

Конечно, на любом приличном современном процессоре вычисление r ^ 2 = dx * dx + dy * dy может быть таким же быстрым, как это...

Вы не указали, как вы хотите измерить расстояние. Я приму L1 (прямолинейный), потому что это проще; возможно, эти идеи могут быть изменены для L2 (евклидова).

Если вы делаете это только для относительно небольшого количества пикселей, просто выполняйте поиск по направлению от исходного пикселя по спирали до тех пор, пока не достигнете черного цвета.

Если вы делаете это для многих / всех из них, как насчет этого: Создайте двумерный массив размером с изображение, где каждая ячейка хранит расстояние до ближайшего черного пикселя (и, если необходимо, координаты). этого пикселя). Выполните четыре развертки: слева направо, справа налево, снизу вверх и сверху вниз. Рассмотрим размах слева направо; при развертывании сохраняйте 1-D столбец, содержащий последний черный пиксель, видимый в каждой строке, и пометьте каждую ячейку в 2-D массиве расстоянием и / или координатами этого пикселя. O (N ^ 2).

В качестве альтернативы, k-d дерево является избыточным; Вы могли бы использовать quadtree. Только немного сложнее для написания кода, чем мой разметчик строк, немного больше памяти (но менее чем в два раза больше) и, возможно, быстрее.

Поиск " Поиск ближайшего соседа " первые две ссылки в Google должны вам помочь.

Если вы делаете это только для 1 пикселя на изображение, я думаю, что вам лучше всего выбрать линейный поиск с шириной в 1 пиксель за раз. Вы не можете взять первую найденную точку, если окно поиска квадратное. Вы должны быть осторожны

Да, поиск ближайшего соседа хорош, но не гарантирует, что вы найдете «ближайший». Перемещение на один пиксель каждый раз приведет к квадратному поиску - диагонали будут дальше, чем горизонталь / вертикаль. Если это важно, вам нужно проверить - продолжайте расширяться до тех пор, пока абсолютная горизонталь не будет иметь расстояние, превышающее «найденный» пиксель, а затем рассчитайте расстояния для всех не черных пикселей, которые были обнаружены.

Хорошо, это звучит интересно. Я сделал версию сулута на С ++, не знаю, поможет ли это вам. Я думаю, что это работает достаточно быстро, поскольку это почти мгновенно на матрице 800 * 600. Если у вас есть какие-либо вопросы, просто задавайте.

Извините за любые ошибки, которые я сделал, это код 10 минут ... Это итеративная версия (я тоже планировал сделать рекурсивную, но передумал). Алгоритм можно улучшить, не добавляя ни одной точки в массив точек, который находится на большем расстоянии от начальной точки, чем min_dist, но это включает в себя вычисление для каждого пикселя (несмотря на его цвет) расстояния от начальной точки.

Надеюсь, это поможет

//(c++ version)
#include<iostream>
#include<cmath>
#include<ctime>
using namespace std;
//ITERATIVE VERSION

//picture witdh&height
#define width 800
#define height 600
//indexex
int i,j;

//initial point coordinates
int x,y;
//variables to work with the array
int p,u;
//minimum dist
double min_dist=2000000000;
//array for memorising the points added
struct point{
  int x;
  int y;
} points[width*height];
double dist;
bool viz[width][height];
// direction vectors, used for adding adjacent points in the "points" array.
int dx[8]={1,1,0,-1,-1,-1,0,1};
int dy[8]={0,1,1,1,0,-1,-1,-1};
int k,nX,nY;
//we will generate an image with white&black pixels (0&1)
bool image[width-1][height-1];
int main(){
    srand(time(0));
    //generate the random pic
    for(i=1;i<=width-1;i++)
        for(j=1;j<=height-1;j++)
            if(rand()%10001<=9999) //9999/10000 chances of generating a black pixel
            image[i][j]=0;
            else image[i][j]=1;
    //random coordinates for starting x&y
    x=rand()%width;
    y=rand()%height;
    p=1;u=1;
    points[1].x=x;
    points[1].y=y;
    while(p<=u){
        for(k=0;k<=7;k++){
          nX=points[p].x+dx[k];
          nY=points[p].y+dy[k];
          //nX&nY are the coordinates for the next point
          //if we haven't added the point yet
          //also check if the point is valid
          if(nX>0&&nY>0&&nX<width&&nY<height)
          if(viz[nX][nY] == 0 ){
              //mark it as added
              viz[nX][nY]=1;
              //add it in the array
              u++;
              points[u].x=nX;
              points[u].y=nY;
              //if it's not black
              if(image[nX][nY]!=0){
              //calculate the distance
              dist=(x-nX)*(x-nX) + (y-nY)*(y-nY);
              dist=sqrt(dist);
              //if the dist is shorter than the minimum, we save it
              if(dist<min_dist)
                  min_dist=dist;
                  //you could save the coordinates of the point that has
                  //the minimum distance too, like sX=nX;, sY=nY;
              }
            }
        }
        p++;
}
    cout<<"Minimum dist:"<<min_dist<<"\n";
return 0;
}

Я уверен, что это можно было бы сделать лучше, но вот некоторый код, который выполняет поиск по периметру квадрата вокруг центрального пикселя, сначала исследуя центр и двигаясь к углам.Если пиксель не найден, периметр (радиус) расширяется до тех пор, пока не будет достигнуто либо ограничение по радиусу, либо найден пиксель.Первой реализацией был цикл, выполняющий простую спираль вокруг центральной точки, но, как уже отмечалось, он не находит абсолютно ближайший пиксель.Создание SomeBigObjCStruct внутри цикла было очень медленным - удаление его из цикла сделало его достаточно хорошим, и был использован спиральный подход.Но в любом случае, вот эта реализация - будьте осторожны, тестирование практически не проводилось.

Все это делается с помощью сложения и вычитания целых чисел.

- (SomeBigObjCStruct *)nearestWalkablePoint:(SomeBigObjCStruct)point {    

typedef struct _testPoint { // using the IYMapPoint object here is very slow
    int x;
    int y;
} testPoint;

// see if the point supplied is walkable
testPoint centre;
centre.x = point.x;
centre.y = point.y;

NSMutableData *map = [self getWalkingMapDataForLevelId:point.levelId];

// check point for walkable (case radius = 0)
if(testThePoint(centre.x, centre.y, map) != 0) // bullseye
    return point;

// radius is the distance from the location of point. A square is checked on each iteration, radius units from point.
// The point with y=0 or x=0 distance is checked first, i.e. the centre of the side of the square. A cursor variable
// is used to move along the side of the square looking for a walkable point. This proceeds until a walkable point
// is found or the side is exhausted. Sides are checked until radius is exhausted at which point the search fails.
int radius = 1;

BOOL leftWithinMap = YES, rightWithinMap = YES, upWithinMap = YES, downWithinMap = YES;

testPoint leftCentre, upCentre, rightCentre, downCentre;
testPoint leftUp, leftDown, rightUp, rightDown;
testPoint upLeft, upRight, downLeft, downRight;

leftCentre = rightCentre = upCentre = downCentre = centre;

int foundX = -1;
int foundY = -1;

while(radius < 1000) {

    // radius increases. move centres outward
    if(leftWithinMap == YES) {

        leftCentre.x -= 1; // move left

        if(leftCentre.x < 0) {

            leftWithinMap = NO;
        }
    }

    if(rightWithinMap == YES) {

        rightCentre.x += 1; // move right

        if(!(rightCentre.x < kIYMapWidth)) {

            rightWithinMap = NO;
        }
    }

    if(upWithinMap == YES) {

        upCentre.y -= 1; // move up

        if(upCentre.y < 0) {

            upWithinMap = NO;
        }
    }

    if(downWithinMap == YES) {

        downCentre.y += 1; // move down

        if(!(downCentre.y < kIYMapHeight)) {

            downWithinMap = NO;
        }
    }

    // set up cursor values for checking along the sides of the square
    leftUp = leftDown = leftCentre;
    leftUp.y -= 1;
    leftDown.y += 1;
    rightUp = rightDown = rightCentre;
    rightUp.y -= 1;
    rightDown.y += 1;
    upRight = upLeft = upCentre;
    upRight.x += 1;
    upLeft.x -= 1;
    downRight = downLeft = downCentre;
    downRight.x += 1;
    downLeft.x -= 1;

    // check centres
    if(testThePoint(leftCentre.x, leftCentre.y, map) != 0) {

        foundX = leftCentre.x;
        foundY = leftCentre.y;
        break;
    }
    if(testThePoint(rightCentre.x, rightCentre.y, map) != 0) {

        foundX = rightCentre.x;
        foundY = rightCentre.y;
        break;
    }
    if(testThePoint(upCentre.x, upCentre.y, map) != 0) {

        foundX = upCentre.x;
        foundY = upCentre.y;
        break;
    }
    if(testThePoint(downCentre.x, downCentre.y, map) != 0) {

        foundX = downCentre.x;
        foundY = downCentre.y;
        break;
    }

    int i;

    for(i = 0; i < radius; i++) {

        if(leftWithinMap == YES) {
            // LEFT Side - stop short of top/bottom rows because up/down horizontal cursors check that line
            // if cursor position is within map
            if(i < radius - 1) {

                if(leftUp.y > 0) {
                    // check it
                    if(testThePoint(leftUp.x, leftUp.y, map) != 0) {
                        foundX = leftUp.x;
                        foundY = leftUp.y;
                        break;
                    }
                    leftUp.y -= 1; // moving up
                }
                if(leftDown.y < kIYMapHeight) {
                    // check it
                    if(testThePoint(leftDown.x, leftDown.y, map) != 0) {
                        foundX = leftDown.x;
                        foundY = leftDown.y;
                        break;
                    }
                    leftDown.y += 1; // moving down
                }
            }
        }

        if(rightWithinMap == YES) {
            // RIGHT Side
            if(i < radius - 1) {

                if(rightUp.y > 0) {

                    if(testThePoint(rightUp.x, rightUp.y, map) != 0) {
                        foundX = rightUp.x;
                        foundY = rightUp.y;
                        break;
                    }
                    rightUp.y -= 1; // moving up
                }
                if(rightDown.y < kIYMapHeight) {

                    if(testThePoint(rightDown.x, rightDown.y, map) != 0) {
                        foundX = rightDown.x;
                        foundY = rightDown.y;
                        break;
                    }
                    rightDown.y += 1; // moving down
                }
            }
        }

        if(upWithinMap == YES) {
            // UP Side
            if(upRight.x < kIYMapWidth) {

                if(testThePoint(upRight.x, upRight.y, map) != 0) {
                    foundX = upRight.x;
                    foundY = upRight.y;
                    break;
                }
                upRight.x += 1; // moving right
            }
            if(upLeft.x > 0) {

                if(testThePoint(upLeft.x, upLeft.y, map) != 0) {
                    foundX = upLeft.x;
                    foundY = upLeft.y;
                    break;
                }
                upLeft.y -= 1; // moving left
            }
        }

        if(downWithinMap == YES) {
            // DOWN Side
            if(downRight.x < kIYMapWidth) {

                if(testThePoint(downRight.x, downRight.y, map) != 0) {
                    foundX = downRight.x;
                    foundY = downRight.y;
                    break;
                }
                downRight.x += 1; // moving right
            }
            if(downLeft.x > 0) {

                if(testThePoint(upLeft.x, upLeft.y, map) != 0) {
                    foundX = downLeft.x;
                    foundY = downLeft.y;
                    break;
                }
                downLeft.y -= 1; // moving left
            }
        }
    }

    if(foundX != -1 && foundY != -1) {
        break;
    }

    radius++;
}

// build the return object
if(foundX != -1 && foundY != -1) {

    SomeBigObjCStruct *foundPoint = [SomeBigObjCStruct mapPointWithX:foundX Y:foundY levelId:point.levelId];
    foundPoint.z = [self zWithLevelId:point.levelId];
    return foundPoint;
}
return nil;

}

Вы можете комбинировать множество способов ускорить это.

  • Способ ускорить поиск пикселей - это использовать то, что я называю картой пространственного поиска.По сути, это карта с уменьшенной дискретизацией (скажем, размером 8x8 пикселей, но это компромисс) пикселей в этом блоке.Значениями могут быть "пикселы не установлены", "частичные пикселы установлены", "все пикселы установлены".Таким образом, при одном чтении можно определить, заполнен ли блок / ячейка, частично заполнена или пуста.
  • сканирование прямоугольника вокруг центра может быть не идеальным, поскольку есть много пикселей / ячеек, которые находятся далеко друг от друга.Я использую алгоритм рисования окружностей (Bresenham), чтобы уменьшить накладные расходы.
  • считывание необработанных значений пикселей может выполняться горизонтальными пакетами, например байтом (для ячейки размером 8x8 или кратным ему), dword или long.Это должно снова дать вам серьезное ускорение.
  • вы также можете использовать несколько уровней "карт пространственного поиска", это опять же компромисс.

Для вычисления расстояния можно использовать упомянутую таблицу подстановки, но это компромисс между пропускной способностью (кэша) и скоростью вычисления (я не знаю, как это работает, например, на графическом процессоре).

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

Лицензировано под: CC-BY-SA с атрибуция
Не связан с StackOverflow
scroll top