Прямое БПФ изображения и обратное БПФ изображения, чтобы получить тот же результат

#c #c #image-processing #fft #fftw

#c #c #обработка изображений #БПФ #fftw #c

Вопрос:

Я пытаюсь выполнить БПФ изображения, используя библиотеку изhttp://www.fftw.org / чтобы я мог выполнить свертку в частотной области. Но я не могу понять, как заставить это работать. Чтобы понять, как это сделать, я пытаюсь перенаправить БПФ изображения в виде массива pixelcolors, а затем обратное БПФ, чтобы получить тот же массив pixelcolors. Вот что я делаю:

 fftw_plan planR, planG, planB;
fftw_complex *inR, *inG, *inB, *outR, *outG, *outB, *resultR, *resultG, *resultB;

//Allocate arrays.
inR = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * width * width);
inG = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * width * width);
inB = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * width * width);

outR = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * width * width);
outG = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * width * width);
outB = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * width * width);

resultR = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * width * width);
resultG = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * width * width);
resultB = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * width * width);

//Fill in arrays with the pixelcolors.
for (int y = 0; y < height; y  ) {
    for (int x = 0; x < width; x  ) {
        int currentIndex = ((y * width)   (x)) * 3;
        inR[y * width   x][0] = pixelColors[currentIndex];
        inG[y * width   x][0] = pixelColors[currentIndex   1];
        inB[y * width   x][0] = pixelColors[currentIndex   2];
    }
}

//Forward plans.
planR = fftw_plan_dft_2d(width, width, inR, outR, FFTW_FORWARD, FFTW_MEASURE);
planG = fftw_plan_dft_2d(width, width, inG, outG, FFTW_FORWARD, FFTW_MEASURE);
planB = fftw_plan_dft_2d(width, width, inB, outB, FFTW_FORWARD, FFTW_MEASURE);

//Forward FFT.
fftw_execute(planR);
fftw_execute(planG);
fftw_execute(planB);

//Backward plans.
planR = fftw_plan_dft_2d(width, width, outR, resultR, FFTW_BACKWARD, FFTW_MEASURE);
planG = fftw_plan_dft_2d(width, width, outG, resultG, FFTW_BACKWARD, FFTW_MEASURE);
planB = fftw_plan_dft_2d(width, width, outB, resultB, FFTW_BACKWARD, FFTW_MEASURE);

//Backward fft
fftw_execute(planR);
fftw_execute(planG);
fftw_execute(planB);

//Overwrite the pixelcolors with the result.
for (int y = 0; y < height; y  ) {
    for (int x = 0; x < width; x  ) {
        int currentIndex = ((y * width)   (x)) * 3;
        pixelColors[currentIndex] = resultR[y * width   x][0];
        pixelColors[currentIndex   1] = resultG[y * width   x][0];
        pixelColors[currentIndex   2] = resultB[y * width   x][0];
    }
}
  

Может кто-нибудь, пожалуйста, показать мне пример того, как переслать БПФ изображение, а затем выполнить обратное БПФ изображение, используя FFTW, чтобы получить тот же результат? Я просмотрел множество примеров, показывающих, как использовать FFTW для БПФ, но я не могу понять, как это применимо к моей ситуации, когда у меня есть массив цветов пикселей, представляющих изображение.

Ответ №1:

Одна важная вещь, на которую следует обратить внимание, когда вы выполняете прямое БПФ, за которым следует обратное БПФ, заключается в том, что это обычно приводит к тому, что к конечному результату применяется коэффициент масштабирования N, т. Е. Результирующие значения пикселей изображения необходимо разделить на N, чтобы соответствовать исходным значениям пикселей. (N — размер БПФ.) Таким образом, ваш цикл вывода, вероятно, должен выглядеть примерно так:

 //Overwrite the pixelcolors with the result.
for (int y = 0; y < height; y  ) {
    for (int x = 0; x < width; x  ) {
        int currentIndex = ((y * width)   (x)) * 3;
        pixelColors[currentIndex] = resultR[y * width   x][0] / (width * height);
        pixelColors[currentIndex   1] = resultG[y * width   x][0] / (width * height);
        pixelColors[currentIndex   2] = resultB[y * width   x][0] / (width * height);
    }
}
  

Также обратите внимание, что вы, вероятно, захотите выполнить БПФ от реального к сложному, за которым последует преобразование от сложного к реальному IFFT (несколько более эффективное с точки зрения как памяти, так и производительности). На данный момент, похоже, что вы выполняете сложное преобразование в обоих направлениях, что нормально, но вы неправильно заполняете входные массивы. Если вы собираетесь придерживаться complex-to-complex, то вы, вероятно, захотите изменить свой цикл ввода на что-то вроде этого:

 //Fill in arrays with the pixelcolors.
for (int y = 0; y < height; y  ) {
    for (int x = 0; x < width; x  ) {
        int currentIndex = ((y * width)   (x)) * 3;
        inR[y * width   x][0] = (double)pixelColors[currentIndex];
        inR[y * width   x][1] = 0.0;
        inG[y * width   x][0] = (double)pixelColors[currentIndex   1];
        inG[y * width   x][1] = 0.0;
        inB[y * width   x][0] = (double)pixelColors[currentIndex   2];
        inB[y * width   x][1] = 0.0;
    }
}
  

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

Еще одна вещь, на которую следует обратить внимание: когда вы в конечном итоге получите эту работу, вы обнаружите, что производительность ужасна — создание плана занимает много времени относительно времени, затраченного на фактическое БПФ. Идея в том, что вы создаете план только один раз, но используете его для выполнения многих БПФ. Итак, вы захотите отделить создание плана от фактического кода БПФ и поместить его в процедуру инициализации или конструктор или что-то еще.

Ответ №2:

Но если вы используете realToComplex или complextorealфункцию, обратите внимание на тот факт, что изображение будет сохранено в матрице размеров [высота x (ширина / 2 1)], и если вы хотите выполнить некоторые промежуточные вычисления в частотной области, они станут немного сложнее…

Ответ №3:

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

Комментарии:

1. Действительно. Другой вариант — использовать FFTW_ESTIMATE вместо FFTW_MEASURE , тогда массивы не будут перезаписаны.