001 //$HeadURL: https://svn.wald.intevation.org/svn/deegree/base/branches/2.3_testing/src/org/deegree/processing/raster/filter/Convolve.java $
002 /*----------------------------------------------------------------------------
003 This file is part of deegree, http://deegree.org/
004 Copyright (C) 2001-2009 by:
005 Department of Geography, University of Bonn
006 and
007 lat/lon GmbH
008
009 This library is free software; you can redistribute it and/or modify it under
010 the terms of the GNU Lesser General Public License as published by the Free
011 Software Foundation; either version 2.1 of the License, or (at your option)
012 any later version.
013 This library is distributed in the hope that it will be useful, but WITHOUT
014 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
015 FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
016 details.
017 You should have received a copy of the GNU Lesser General Public License
018 along with this library; if not, write to the Free Software Foundation, Inc.,
019 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
020
021 Contact information:
022
023 lat/lon GmbH
024 Aennchenstr. 19, 53177 Bonn
025 Germany
026 http://lat-lon.de/
027
028 Department of Geography, University of Bonn
029 Prof. Dr. Klaus Greve
030 Postfach 1147, 53001 Bonn
031 Germany
032 http://www.geographie.uni-bonn.de/deegree/
033
034 e-mail: info@deegree.org
035 ----------------------------------------------------------------------------*/
036 package org.deegree.processing.raster.filter;
037
038 /**
039 * This class offeres some convenience methods for convolving a data matrix. <br>
040 * Convolution filtering is often used to reduce the effects of noise in images or to sharpen the
041 * detail in rasters. Convolution filtering is a form of spatial filtering that computes each output
042 * sample by multiplying elements of a kernel with the samples surrounding a particular source
043 * sample. <br>
044 * Convolution filtering operates on a group of input values surrounding a center pixel. The
045 * adjoining values provide important information about values trends in the area of the values
046 * being processed. <br>
047 * Convolution filtering moves across the source raster, cell by cell, placing resulting values into
048 * the destination raster. The resulting value of each source cell depends on the group of values
049 * surrounding the source cell. Using the cell values of the source cell's neighbors, the
050 * convolution process calculates the spatial frequency activity in the area, making it possible to
051 * filter the values based on the spatial frequency of the area. <br>
052 * Convolution filtering uses a convolve kernel, containing an array of convolution coefficient
053 * values, called key elements. The array is not restricted to any particular size, and does not
054 * even have to be square. The kernel can be 1 x 1, 3 x 3, 5 x 5, M x N, and so on. A larger kernel
055 * size affords a more precise filtering operation by increasing the number of neighboring cells
056 * used in the calculation. However, the kernel cannot be bigger in any dimension than the source
057 * data. Also, the larger the kernel, the more computations that are required to be performed. For
058 * example, given a 640 x 480 raster and a 3 x 3 kernel, the convolve operation requires over five
059 * million total multiplications and additions. <br>
060 * The convolution filtering operation computes each output sample by multiplying the key elements
061 * of the kernel with the samples surrounding a particular source cell. For each destination cell,
062 * the kernel is rotated 180 degrees and its key element is placed over the source pixel
063 * corresponding with the destination pixel. The key elements are multiplied with the source cell
064 * value under them, and the resulting products are summed together to produce the destination
065 * sample value. <br>
066 * The selection of the weights for the key elements determines the nature of the filtering action,
067 * such as high-pass or low-pass. If the values of the key elements are the reciprocal of the number
068 * of key elements in the kernel (for example, 1/9 for a 3 x 3 kernel), the result is a conventional
069 * low-pass averaging process. If the weights are altered, certain cells in the kernel will have an
070 * increased or decreased influence in the average. <br>
071 * (modified from JAI documentation)
072 *
073 * @version $Revision: 18195 $
074 * @author <a href="mailto:poth@lat-lon.de">Andreas Poth</a>
075 * @author last edited by: $Author: mschneider $
076 *
077 * @version 1.0. $Revision: 18195 $, $Date: 2009-06-18 17:55:39 +0200 (Do, 18. Jun 2009) $
078 *
079 * @since 2.0
080 */
081 public class Convolve {
082
083 /**
084 * convolves the passed float matrix by applying the passed kernel. A kernel for lowpass
085 * filtering may looks like:
086 *
087 * <pre>
088 * float[][] kernel = new float[5][5];
089 * for ( int i = 0; i < kernel.length; i++ ) {
090 * for ( int j = 0; j < kernel[i].length; j++ ) {
091 * kernel[i][j] = 1;
092 * }
093 * }
094 * <pre>
095 * which results in a strong smoothing of the raster.
096 * <BR>
097 * A kernel for highpass filtering may looks like:
098 * <pre>
099 * float[][] kernel = new float[3][3];
100 * for ( int i = 0; i < kernel.length; i++ ) {
101 * for ( int j = 0; j < kernel[i].length; j++ ) {
102 * kernel[i][j] = -1;
103 * }
104 * }
105 * kernel[kernel.length/2][kernel[0].length/2] = 2;
106 * </pre>
107 * Notice that the sum of all kernel values does not must be
108 * == 1! Correct weighting will be done by this method.
109 *
110 * @param source to apply the kernel to.
111 * @param kernel to apply
112 * @return the raster after applying the given kernel.
113 * @throws RasterFilterException if the kernel does not have an odd number of cells.
114 *
115 */
116 public static float[][] perform( float[][] source, float[][] kernel )
117 throws RasterFilterException {
118
119 if ( kernel.length % 2 == 0 || kernel[0].length % 2 == 0 ) {
120 throw new RasterFilterException( "A Kernel must have an odd number of cells in x- and y-direction" );
121 }
122
123 kernel = rotateKernel( kernel );
124
125 float[][] dest = new float[source.length][source[0].length];
126
127 int ww = kernel.length;
128 int hh = kernel[0].length;
129 int xOrigin = ww / 2;
130 int yOrigin = hh / 2;
131 for ( int y = 0; y < source.length; y++ ) {
132 for ( int x = 0; x < source[y].length; x++ ) {
133 float g = 0;
134 float v = 0;
135 for ( int i = -xOrigin; i < -xOrigin + ww; i++ ) {
136 for ( int j = -yOrigin; j < -yOrigin + ww; j++ ) {
137 if ( y + i >= 0 && y + i < source.length && x + j >= 0 && x + j < source[y].length ) {
138 v += source[y + i][x + j] * kernel[xOrigin + i][yOrigin + j];
139 g += kernel[xOrigin + i][yOrigin + j];
140 }
141 }
142 }
143 dest[y][x] = v / g;
144 }
145 }
146
147 return dest;
148 }
149
150 /**
151 * rotates a convolution kernel by 180°
152 *
153 * @param kernel
154 * @return the rotated kernel.
155 */
156 private static float[][] rotateKernel( float[][] kernel ) {
157 for ( int i = 0; i < kernel.length / 2; i++ ) {
158 for ( int j = 0; j < kernel[i].length / 2; j++ ) {
159 float v = kernel[i][j];
160 kernel[i][j] = kernel[kernel.length - i - 1][kernel[i].length - j - 1];
161 kernel[kernel.length - i - 1][kernel[i].length - j - 1] = v;
162 }
163 }
164 return kernel;
165 }
166
167 }