001/**
002 * Copyright (c) 2011, The University of Southampton and the individual contributors.
003 * All rights reserved.
004 *
005 * Redistribution and use in source and binary forms, with or without modification,
006 * are permitted provided that the following conditions are met:
007 *
008 *   *  Redistributions of source code must retain the above copyright notice,
009 *      this list of conditions and the following disclaimer.
010 *
011 *   *  Redistributions in binary form must reproduce the above copyright notice,
012 *      this list of conditions and the following disclaimer in the documentation
013 *      and/or other materials provided with the distribution.
014 *
015 *   *  Neither the name of the University of Southampton nor the names of its
016 *      contributors may be used to endorse or promote products derived from this
017 *      software without specific prior written permission.
018 *
019 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
020 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
021 * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
022 * DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
023 * ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
024 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
025 * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
026 * ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
027 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
028 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
029 */
030package org.openimaj.image.processing.convolution;
031
032import org.openimaj.image.FImage;
033import org.openimaj.image.processing.algorithm.FourierTransform;
034import org.openimaj.image.processor.SinglebandImageProcessor;
035
036import edu.emory.mathcs.jtransforms.fft.FloatFFT_2D;
037
038/**
039 * {@link FImage} convolution performed in the fourier domain.
040 *
041 * @author Jonathon Hare (jsh2@ecs.soton.ac.uk)
042 *
043 */
044public class FourierConvolve implements SinglebandImageProcessor<Float, FImage> {
045        private float[][] kernel;
046
047        /**
048         * Construct the convolution operator with the given kernel
049         *
050         * @param kernel
051         *            the kernel
052         */
053        public FourierConvolve(float[][] kernel) {
054                this.kernel = kernel;
055        }
056
057        /**
058         * Construct the convolution operator with the given kernel
059         *
060         * @param kernel
061         *            the kernel
062         */
063        public FourierConvolve(FImage kernel) {
064                this.kernel = kernel.pixels;
065        }
066
067        @Override
068        public void processImage(FImage image) {
069                convolve(image, kernel, true);
070        }
071
072        /**
073         * Convolve an image with a kernel using an FFT.
074         *
075         * @param image
076         *            The image to convolve
077         * @param kernel
078         *            The kernel
079         * @param inplace
080         *            if true, then output overwrites the input, otherwise a new image
081         *            is created.
082         * @return convolved image
083         */
084        public static FImage convolve(FImage image, float[][] kernel, boolean inplace) {
085                final int cols = image.getCols();
086                final int rows = image.getRows();
087
088                final FloatFFT_2D fft = new FloatFFT_2D(rows, cols);
089
090                final float[][] preparedImage = FourierTransform.prepareData(image.pixels, rows, cols, false);
091                fft.complexForward(preparedImage);
092
093                final float[][] preparedKernel = FourierTransform.prepareData(kernel, rows, cols, false);
094                fft.complexForward(preparedKernel);
095
096                for (int y = 0; y < rows; y++) {
097                        for (int x = 0; x < cols; x++) {
098                                final float reImage = preparedImage[y][x * 2];
099                                final float imImage = preparedImage[y][1 + x * 2];
100
101                                final float reKernel = preparedKernel[y][x * 2];
102                                final float imKernel = preparedKernel[y][1 + x * 2];
103
104                                final float re = reImage * reKernel - imImage * imKernel;
105                                final float im = reImage * imKernel + imImage * reKernel;
106
107                                preparedImage[y][x * 2] = re;
108                                preparedImage[y][1 + x * 2] = im;
109                        }
110                }
111
112                fft.complexInverse(preparedImage, true);
113
114                FImage out = image;
115                if (!inplace)
116                        out = new FImage(cols, rows);
117
118                FourierTransform.unprepareData(preparedImage, out, false);
119
120                return out;
121        }
122
123        /**
124         * Convolve an image with a pre-prepared frequency domain filter. The filter
125         * must have the same height as the image and twice the width (to account for
126         * the imaginary components). Real and imaginary components should be interlaced
127         * across the rows.
128         *
129         * @param image
130         *            The image to convolve
131         * @param filter
132         *            the prepared frequency domain filter
133         * @param centered
134         *            true if the prepared filter has the highest frequency in the
135         *            centre.
136         * @return convolved image
137         */
138        public static FImage convolvePrepared(FImage image, FImage filter, boolean centered) {
139                final int cols = image.getCols();
140                final int rows = image.getRows();
141
142                final FloatFFT_2D fft = new FloatFFT_2D(rows, cols);
143
144                final float[][] preparedImage = FourierTransform.prepareData(image.pixels, rows, cols, centered);
145                fft.complexForward(preparedImage);
146
147                final float[][] preparedKernel = filter.pixels;
148
149                for (int y = 0; y < rows; y++) {
150                        for (int x = 0; x < cols; x++) {
151                                final float reImage = preparedImage[y][x * 2];
152                                final float imImage = preparedImage[y][1 + x * 2];
153
154                                final float reKernel = preparedKernel[y][x * 2];
155                                final float imKernel = preparedKernel[y][1 + x * 2];
156
157                                final float re = reImage * reKernel - imImage * imKernel;
158                                final float im = reImage * imKernel + imImage * reKernel;
159
160                                preparedImage[y][x * 2] = re;
161                                preparedImage[y][1 + x * 2] = im;
162                        }
163                }
164
165                fft.complexInverse(preparedImage, true);
166
167                final FImage out = new FImage(cols, rows);
168                FourierTransform.unprepareData(preparedImage, out, centered);
169                return out;
170        }
171}