comparison gnv-artifacts/src/main/java/de/intevation/gnv/raster/Raster.java @ 424:21fbd254db71

Added support for converting 2D rasters into polygons. gnv-artifacts/trunk@472 c6561f87-3c4e-4783-a992-168aeb5c3f6f
author Sascha L. Teichmann <sascha.teichmann@intevation.de>
date Mon, 21 Dec 2009 18:00:54 +0000
parents
children 6e8364e766fa
comparison
equal deleted inserted replaced
423:2402173a1490 424:21fbd254db71
1 package de.intevation.gnv.raster;
2
3 /**
4 * @author Sascha L. Teichmann (sascha.teichmann@intevation.de)
5 */
6 public class Raster
7 {
8 public static class Kernel {
9
10 public static final class Coeff {
11 public int di;
12 public int dj;
13 public double c;
14
15 public Coeff() {
16 }
17
18 public Coeff(int di, int dj, double c) {
19 this.di = di;
20 this.dj = dj;
21 this.c = c;
22 }
23
24 public String toString() {
25 return String.valueOf(c);
26 }
27
28 } // class Coeff
29
30 protected Coeff [] coeffs;
31 protected int width;
32
33 public Kernel() {
34 }
35
36 public Kernel(Coeff [] coeffs, int width) {
37 this.coeffs = coeffs;
38 this.width = width;
39 }
40
41 public static double gauss(double x, double y, double sigma) {
42 double s2 = sigma*sigma;
43 return 1.0d/(2.0d*Math.PI*s2)*Math.exp(-(x*x + y*y)/(2.0d*s2));
44 }
45
46 public static Kernel createGauss() {
47 return createGauss(1.0d, 5);
48 }
49
50 public static Kernel createGauss(double sigma, int width) {
51 Coeff [] coeffs = new Coeff[width*width];
52 double sum = 0d;
53 for (int j = 0; j < width; ++j) {
54 int y = -width/2 + j;
55 for (int i = 0; i < width; ++i) {
56 int x = -width/2 + i;
57 double c = gauss(x, y, sigma);
58 coeffs[j*width + i] = new Coeff(x, y, c);
59 sum += c;
60 }
61 }
62 sum = 1.0/sum;
63 for (int i = 0; i < coeffs.length; ++i) {
64 coeffs[i].c *= sum;
65 }
66 return new Kernel(coeffs, width);
67 }
68
69 public double fold(Access access, int i, int j) {
70
71 double s = 0.0d;
72
73 double unused = 0d;
74
75 for (int k = 0; k < coeffs.length; ++k) {
76 Coeff coeff = coeffs[k];
77 double v = access.get(i + coeff.di, j + coeff.dj);
78 if (Double.isNaN(v)) {
79 unused += coeff.c;
80 }
81 else {
82 s += v * coeff.c;
83 }
84 }
85
86 return s*(1.0d - unused);
87 }
88
89 public String toString() {
90 StringBuilder sb = new StringBuilder();
91
92 int height = coeffs.length/width;
93
94 for (int j = 0; j < height; ++j) {
95 if (j > 0) {
96 sb.append(System.getProperty("line.separator"));
97 }
98 for (int i = 0; i < width; ++i) {
99 if (i > 0) {
100 sb.append("\t");
101 }
102 sb.append(coeffs[j*width + i]);
103 }
104 }
105
106 return sb.toString();
107 }
108 } // class Kernel
109
110 public interface Access {
111
112 double get(int i, int j);
113
114 } // interface Access
115
116 public interface ValueToIndex {
117
118 int toIndex(double value);
119
120 } // interface ValueToIndex;
121
122 public static class IsoClasses
123 implements ValueToIndex
124 {
125 protected double m;
126 protected double b;
127 protected int numClasses;
128
129 public IsoClasses(Raster raster, int numClasses) {
130
131 this.numClasses = numClasses;
132
133 double min = Double.MAX_VALUE;
134 double max = -Double.MAX_VALUE;
135
136 double [] src = raster.raster;
137
138 for (int i = 0; i < src.length; ++i) {
139 double x = src[i];
140 if (!Double.isNaN(x)) {
141 if (x < min) min = x;
142 if (x > max) max = x;
143 }
144 }
145
146 /* f(min) = 0, f(max) = numClasses - 1
147
148 I: 0 = m*min + b
149 II: numClasses - 1 = m*max + b
150
151 II - I:
152 numClasses - 1 = m*(max - min)
153
154 m = (numClasses - 1)/(max - min)
155 b = m*min
156 */
157
158 if (max == min) {
159 m = 0;
160 b = 0;
161 }
162 else {
163 m = (numClasses - 1)/(max - min);
164 b = m*min;
165 }
166 }
167
168 public int toIndex(double value) {
169 return Double.isNaN(value)
170 ? -1
171 : Math.min(Math.max(0, (int)Math.round(m*value + b)), numClasses-1);
172 }
173 } // class IsoClasses
174
175 protected double [] raster;
176 protected int width;
177
178 public Raster() {
179 this(new double[0], 0);
180 }
181
182 public Raster(double [] raster, int width) {
183 this.raster = raster;
184 this.width = width;
185 }
186
187 public int getWidth() {
188 return width;
189 }
190
191 public int getHeight() {
192 return raster.length / width;
193 }
194
195 public int [] toIndexed(ValueToIndex valueToIndex) {
196 int [] dst = new int[raster.length];
197 for (int i = 0; i < raster.length; ++i) {
198 dst[i] = valueToIndex.toIndex(raster[i]);
199 }
200 return dst;
201 }
202
203 public Raster create(Kernel kernel) {
204 double [] dst = new double[raster.length];
205 Raster r = new Raster(dst, width);
206 r.apply(kernel, continueBorder());
207 return r;
208 }
209
210 public void apply(Kernel kernel, Access access) {
211 int height = getHeight();
212 for (int j = 0; j < height; ++j) {
213 int row = j*width;
214 for (int i = 0; i < width; ++i) {
215 raster[row + i] = kernel.fold(access, i, j);
216 }
217 }
218 }
219
220 public Access continueBorder() {
221 return new Access() {
222 int height = getHeight()-1;
223 public double get(int i, int j) {
224 if (i < 0) i = 0;
225 else if (i >= width) i = width-1;
226 if (j < 0) j = 0;
227 else if (j > height) j = height;
228 return raster[j*width + i];
229 }
230 };
231 }
232 }
233 // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8:

http://dive4elements.wald.intevation.org