1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30 package org.openimaj.image.pixel.statistics;
31
32 import java.util.Arrays;
33
34 import org.apache.commons.math.DimensionMismatchException;
35 import org.apache.commons.math.stat.descriptive.MultivariateSummaryStatistics;
36 import org.openimaj.image.MBFImage;
37 import org.openimaj.image.pixel.sampling.FLineSampler;
38 import org.openimaj.math.geometry.line.Line2d;
39 import org.openimaj.math.geometry.point.Point2d;
40 import org.openimaj.math.geometry.point.Point2dImpl;
41 import org.openimaj.math.util.FloatArrayStatsUtils;
42 import org.openimaj.util.array.ArrayUtils;
43
44 import Jama.Matrix;
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61 @SuppressWarnings("deprecation")
62 public class MBFStatisticalPixelProfileModel implements PixelProfileModel<MBFImage> {
63 private MultivariateSummaryStatistics statistics;
64 private int nsamples;
65 private FLineSampler sampler;
66
67 private double[] mean;
68 private Matrix invCovar;
69
70
71
72
73
74
75
76
77
78
79 public MBFStatisticalPixelProfileModel(int nsamples, FLineSampler sampler) {
80 this.nsamples = nsamples;
81 this.statistics = new MultivariateSummaryStatistics(nsamples, true);
82 this.sampler = sampler;
83 }
84
85 private float[] normaliseSamples(float[] samples) {
86 final float sum = FloatArrayStatsUtils.sum(samples);
87
88 for (int i = 0; i < samples.length; i++) {
89 samples[i] /= sum;
90 }
91
92 return samples;
93 }
94
95 @Override
96 public void updateModel(MBFImage image, Line2d line) {
97 final float[] samples = extractNormalisedStacked(image, line);
98
99 try {
100 statistics.addValue(ArrayUtils.convertToDouble(samples));
101 } catch (final DimensionMismatchException e) {
102 throw new RuntimeException(e);
103 }
104
105 invCovar = null;
106 mean = null;
107 }
108
109 private float[][] extract(MBFImage image, Line2d line, int numSamples) {
110 final float[][] samples = new float[image.numBands()][];
111
112 for (int i = 0; i < image.numBands(); i++) {
113 samples[i] = sampler.extractSamples(line, image.getBand(i), numSamples);
114 }
115
116 return samples;
117 }
118
119
120
121
122
123
124
125
126 private float[] extractNormalisedStacked(MBFImage image, Line2d line) {
127 final float[] samples = new float[nsamples * image.numBands()];
128
129 for (int i = 0; i < image.numBands(); i++) {
130 final float[] s = sampler.extractSamples(line, image.getBand(i), nsamples);
131 System.arraycopy(s, 0, samples, i * nsamples, nsamples);
132 }
133
134 return normaliseSamples(samples);
135 }
136
137
138
139
140 public double[] getMean() {
141 if (mean == null) {
142 mean = statistics.getMean();
143 invCovar = new Matrix(statistics.getCovariance().getData()).inverse();
144 }
145
146 return mean;
147 }
148
149
150
151
152 public Matrix getCovariance() {
153 if (mean == null) {
154 mean = statistics.getMean();
155 invCovar = new Matrix(statistics.getCovariance().getData()).inverse();
156 }
157 return new Matrix(statistics.getCovariance().getData());
158 }
159
160
161
162
163 public Matrix getInverseCovariance() {
164 if (mean == null) {
165 mean = statistics.getMean();
166 invCovar = new Matrix(statistics.getCovariance().getData()).inverse();
167 }
168 return invCovar;
169 }
170
171
172
173
174
175
176
177
178
179
180 public float computeMahalanobis(float[] vector) {
181 if (mean == null) {
182 mean = statistics.getMean();
183 try {
184 invCovar = new Matrix(statistics.getCovariance().getData()).inverse();
185 } catch (final RuntimeException e) {
186 invCovar = Matrix.identity(nsamples, nsamples);
187 }
188 }
189
190 final double[] meanCentered = new double[mean.length];
191 for (int i = 0; i < mean.length; i++) {
192 meanCentered[i] = vector[i] - mean[i];
193 }
194
195 final Matrix mct = new Matrix(new double[][] { meanCentered });
196 final Matrix mc = mct.transpose();
197
198 final Matrix dist = mct.times(invCovar).times(mc);
199
200 return (float) dist.get(0, 0);
201 }
202
203
204
205
206
207
208
209
210
211
212
213 public float computeMahalanobis(MBFImage image, Line2d line) {
214 final float[] samples = extractNormalisedStacked(image, line);
215 return computeMahalanobis(samples);
216 }
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235 public float[] computeMahalanobisWindowed(MBFImage image, Line2d line, int numSamples) {
236 final float[][] samples = extract(image, line, numSamples);
237 return computeMahalanobisWindowed(samples);
238 }
239
240 @Override
241 public Point2d computeNewBest(MBFImage image, Line2d line, int numSamples) {
242 final float[] resp = computeMahalanobisWindowed(image, line, numSamples);
243
244 final int minIdx = ArrayUtils.minIndex(resp);
245 final int offset = (numSamples - nsamples) / 2;
246
247 if (resp[offset] == resp[minIdx])
248
249 return line.calculateCentroid();
250
251
252
253 line = this.sampler.getSampleLine(line, image.getBand(0), numSamples);
254
255 final float x = line.begin.getX();
256 final float y = line.begin.getY();
257 final float dxStep = (line.end.getX() - x) / (numSamples - 1);
258 final float dyStep = (line.end.getY() - y) / (numSamples - 1);
259
260 return new Point2dImpl(x + (minIdx + offset) * dxStep, y + (minIdx + offset) * dyStep);
261 }
262
263 @Override
264 public float computeMovementDistance(MBFImage image, Line2d line, int numSamples, Point2d pt) {
265 final Line2d sampleLine = sampler.getSampleLine(line, image.getBand(0), numSamples);
266
267 return (float) (2 * Line2d.distance(sampleLine.calculateCentroid(), pt) / sampleLine.calculateLength());
268 }
269
270
271
272
273
274
275
276
277
278
279
280
281
282 public float[] computeMahalanobisWindowed(float[][] vector) {
283 final int maxShift = vector.length - nsamples + 1;
284
285 final float[] responses = new float[maxShift];
286 float[] samples = new float[nsamples * vector.length];
287 for (int i = 0; i < maxShift; i++) {
288 for (int j = 0; j < vector.length; j++) {
289 System.arraycopy(vector[j], i, samples, nsamples * j, nsamples);
290 }
291 samples = normaliseSamples(samples);
292 responses[i] = computeMahalanobis(samples);
293 }
294
295 return responses;
296 }
297
298 @Override
299 public String toString() {
300 return "\nPixelProfileModel[\n" +
301 "\tcount = " + statistics.getN() + "\n" +
302 "\tmean = " + Arrays.toString(statistics.getMean()) + "\n" +
303 "\tcovar = " + statistics.getCovariance() + "\n" +
304 "]";
305 }
306
307
308
309
310 public int getNumberSamples() {
311 return nsamples;
312 }
313
314
315
316
317 public FLineSampler getSampler() {
318 return sampler;
319 }
320
321 @Override
322 public float computeCost(MBFImage image, Line2d line) {
323 return computeMahalanobis(image, line);
324 }
325 }