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.processor.SinglebandImageProcessor;
034import org.openimaj.math.matrix.MatrixUtils;
035
036import Jama.SingularValueDecomposition;
037
038/**
039 * Base class for implementation of classes that perform convolution operations
040 * on @link{FImage}s as a @link{SinglebandImageProcessor}, with the kernel
041 * itself formed from and @link{FImage}.
042 * 
043 * @author Jonathon Hare (jsh2@ecs.soton.ac.uk)
044 */
045public class FConvolution implements SinglebandImageProcessor<Float, FImage> {
046        /** The kernel */
047        public FImage kernel;
048
049        private ConvolveMode mode;
050
051        interface ConvolveMode {
052                public void convolve(FImage f);
053
054                class OneD implements ConvolveMode {
055                        private float[] kernel;
056                        private boolean rowMode;
057
058                        OneD(FImage image) {
059                                if (image.height == 1) {
060                                        this.rowMode = true;
061                                        this.kernel = image.pixels[0];
062
063                                }
064                                else {
065                                        this.rowMode = false;
066                                        this.kernel = new float[image.height];
067                                        for (int i = 0; i < image.height; i++)
068                                                this.kernel[i] = image.pixels[i][0];
069                                }
070                        }
071
072                        @Override
073                        public void convolve(FImage f) {
074                                if (this.rowMode)
075                                        FImageConvolveSeparable.convolveHorizontal(f, kernel);
076                                else
077                                        FImageConvolveSeparable.convolveVertical(f, kernel);
078                        }
079
080                }
081
082                class Separable implements ConvolveMode {
083                        private float[] row;
084                        private float[] col;
085
086                        Separable(SingularValueDecomposition svd) {
087
088                                final int nrows = svd.getU().getRowDimension();
089
090                                this.row = new float[nrows];
091                                this.col = new float[nrows];
092
093                                final float factor = (float) Math.sqrt(svd.getS().get(0, 0));
094                                for (int i = 0; i < nrows; i++) {
095                                        this.row[i] = (float) svd.getU().get(i, 0) * factor;
096                                        this.col[i] = (float) svd.getV().get(i, 0) * factor;
097                                }
098                        }
099
100                        @Override
101                        public void convolve(FImage f) {
102                                FImageConvolveSeparable.convolveHorizontal(f, row);
103                                FImageConvolveSeparable.convolveVertical(f, col);
104                        }
105                }
106
107                class BruteForce implements ConvolveMode {
108                        protected FImage kernel;
109
110                        BruteForce(FImage kernel) {
111                                this.kernel = kernel;
112                        }
113
114                        @Override
115                        public void convolve(FImage image) {
116                                final int kh = kernel.height;
117                                final int kw = kernel.width;
118                                final int hh = kh / 2;
119                                final int hw = kw / 2;
120                                final FImage clone = image.newInstance(image.width, image.height);
121                                for (int y = hh; y < image.height - (kh - hh); y++) {
122                                        for (int x = hw; x < image.width - (kw - hw); x++) {
123                                                float sum = 0;
124                                                for (int j = 0, jj = kh - 1; j < kh; j++, jj--) {
125                                                        for (int i = 0, ii = kw - 1; i < kw; i++, ii--) {
126                                                                final int rx = x + i - hw;
127                                                                final int ry = y + j - hh;
128
129                                                                sum += image.pixels[ry][rx] * kernel.pixels[jj][ii];
130                                                        }
131                                                }
132                                                clone.pixels[y][x] = sum;
133                                        }
134                                }
135                                image.internalAssign(clone);
136                        }
137                }
138        }
139
140        /**
141         * Construct the convolution operator with the given kernel
142         * 
143         * @param kernel
144         *            the kernel
145         */
146        public FConvolution(FImage kernel) {
147                this.kernel = kernel;
148                setup(false);
149        }
150
151        /**
152         * Construct the convolution operator with the given kernel
153         * 
154         * @param kernel
155         *            the kernel
156         */
157        public FConvolution(float[][] kernel) {
158                this.kernel = new FImage(kernel);
159                setup(false);
160        }
161
162        /**
163         * Set brute-force convolution; disables kernel separation and other
164         * optimisations.
165         * 
166         * @param brute
167         */
168        public void setBruteForce(boolean brute) {
169                setup(brute);
170        }
171
172        private void setup(boolean brute) {
173                if (brute) {
174                        this.mode = new ConvolveMode.BruteForce(this.kernel);
175                        return;
176                }
177                if (this.kernel.width == 1 || this.kernel.height == 1) {
178                        this.mode = new ConvolveMode.OneD(kernel);
179                }
180                else {
181                        MatrixUtils.matrixFromFloat(this.kernel.pixels);
182                        final SingularValueDecomposition svd = new SingularValueDecomposition(
183                                        MatrixUtils.matrixFromFloat(this.kernel.pixels));
184                        if (svd.rank() == 1)
185                                this.mode = new ConvolveMode.Separable(svd);
186                        else
187                                this.mode = new ConvolveMode.BruteForce(this.kernel);
188                }
189        }
190
191        /*
192         * (non-Javadoc)
193         * 
194         * @see
195         * org.openimaj.image.processor.ImageProcessor#processImage(org.openimaj
196         * .image.Image)
197         */
198        @Override
199        public void processImage(FImage image) {
200                mode.convolve(image);
201        }
202
203        /**
204         * Return the kernel response at the x,y in the given image.
205         * 
206         * This method will throw an array index out of bounds if x,y requests
207         * pixels outside the image bounds
208         * 
209         * @param x
210         * @param y
211         * @param image
212         * @return the kernel response at the given coordinates
213         */
214        public float responseAt(int x, int y, FImage image) {
215                float sum = 0;
216                final int kh = kernel.height;
217                final int kw = kernel.width;
218                final int hh = kh / 2;
219                final int hw = kw / 2;
220
221                for (int j = 0, jj = kh - 1; j < kh; j++, jj--) {
222                        for (int i = 0, ii = kw - 1; i < kw; i++, ii--) {
223                                final int rx = x + i - hw;
224                                final int ry = y + j - hh;
225
226                                sum += image.pixels[ry][rx] * kernel.pixels[jj][ii];
227                        }
228                }
229                return sum;
230        }
231}