gvsig-raster / org.gvsig.raster.principalcomponents / trunk / org.gvsig.raster.principalcomponents / org.gvsig.raster.principalcomponents.algorithm / src / main / java / org / gvsig / raster / principalcomponents / algorithm / PCAStatisticsProcess.java @ 2125
History | View | Annotate | Download (10.2 KB)
1 |
/* gvSIG. Geographic Information System of the Valencian Government
|
---|---|
2 |
*
|
3 |
* Copyright (C) 2007-2008 Infrastructures and Transports Department
|
4 |
* of the Valencian Government (CIT)
|
5 |
*
|
6 |
* This program is free software; you can redistribute it and/or
|
7 |
* modify it under the terms of the GNU General Public License
|
8 |
* as published by the Free Software Foundation; either version 2
|
9 |
* of the License, or (at your option) any later version.
|
10 |
*
|
11 |
* This program is distributed in the hope that it will be useful,
|
12 |
* but WITHOUT ANY WARRANTY; without even the implied warranty of
|
13 |
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
14 |
* GNU General Public License for more details.
|
15 |
*
|
16 |
* You should have received a copy of the GNU General Public License
|
17 |
* along with this program; if not, write to the Free Software
|
18 |
* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
|
19 |
* MA 02110-1301, USA.
|
20 |
*
|
21 |
*/
|
22 |
package org.gvsig.raster.principalcomponents.algorithm; |
23 |
|
24 |
|
25 |
import java.awt.geom.Rectangle2D; |
26 |
import java.util.List; |
27 |
|
28 |
import org.cresques.Messages; |
29 |
import org.gvsig.fmap.dal.coverage.RasterLocator; |
30 |
import org.gvsig.fmap.dal.coverage.dataset.Buffer; |
31 |
import org.gvsig.fmap.dal.coverage.datastruct.Extent; |
32 |
import org.gvsig.fmap.dal.coverage.datastruct.NoData; |
33 |
import org.gvsig.fmap.dal.coverage.exception.CloneException; |
34 |
import org.gvsig.fmap.dal.coverage.exception.ProcessInterruptedException; |
35 |
import org.gvsig.fmap.dal.coverage.store.RasterDataStore; |
36 |
import org.gvsig.fmap.dal.exception.CloseException; |
37 |
import org.gvsig.raster.algorithm.process.DataProcess; |
38 |
import org.gvsig.raster.algorithm.process.ProcessException; |
39 |
import org.gvsig.raster.roi.ROI; |
40 |
import org.slf4j.Logger; |
41 |
import org.slf4j.LoggerFactory; |
42 |
|
43 |
import Jama.EigenvalueDecomposition; |
44 |
import Jama.Matrix; |
45 |
|
46 |
/**
|
47 |
* This Process computes a <code>PCStatsDataStructure</code> Object from bands of two images in a study zone.
|
48 |
* <code>PCStatsDataStructure</code> is a object that contains the needed statistics to Principal Components decomposition
|
49 |
* of images.
|
50 |
*
|
51 |
* @author Nacho Brodin (nachobrodin@gmail.com)
|
52 |
*
|
53 |
*/
|
54 |
public class PCAStatisticsProcess extends DataProcess { |
55 |
public static final String LAYER_ZONA_ESTUDIO = "LAYER_ZONA_ESTUDIO"; |
56 |
public static String BANDS = "BANDS"; |
57 |
public static String PATH = "PATH"; |
58 |
public static String ROI = "ROI"; |
59 |
|
60 |
public static String RASTER_STORE = "RASTER_STORE"; |
61 |
public static String STATS_RESULT = "STATS_RESULT"; |
62 |
|
63 |
private boolean[] bands = null; |
64 |
private RasterDataStore store = null; |
65 |
private Buffer buffer = null; |
66 |
private double autovalors[] = null; |
67 |
private Matrix coVarMatrix = null; |
68 |
private Matrix autoVectorMatrix = null; |
69 |
private double coVar[][] = null; |
70 |
private static Logger logger = LoggerFactory.getLogger(PCAStatisticsProcess.class.getName()); |
71 |
private double[][] maxminmean = null; |
72 |
private Extent extentResult = null; |
73 |
private String path = null; |
74 |
private List<ROI> rois = null; |
75 |
|
76 |
private boolean[] statsByBandCalculated = null; |
77 |
|
78 |
public static void registerParameters() { |
79 |
registerInputParameter(RASTER_STORE, RasterDataStore.class, PrincipalComponentsAlgorithmLibrary.PC_STATS_PROCESS_LABEL); |
80 |
registerInputParameter(BANDS, Boolean[].class, PrincipalComponentsAlgorithmLibrary.PC_STATS_PROCESS_LABEL); |
81 |
registerInputParameter(PATH, String.class, PrincipalComponentsAlgorithmLibrary.PC_STATS_PROCESS_LABEL);
|
82 |
registerInputParameter(ROI, List.class, PrincipalComponentsAlgorithmLibrary.PC_STATS_PROCESS_LABEL);
|
83 |
|
84 |
registerOutputParameter(STATS_RESULT, PCStatsDataStructure.class, PrincipalComponentsAlgorithmLibrary.PC_STATS_PROCESS_LABEL); |
85 |
registerOutputParameter(PATH, String.class, PrincipalComponentsAlgorithmLibrary.PC_STATS_PROCESS_LABEL);
|
86 |
registerOutputParameter(BANDS, Boolean[].class, PrincipalComponentsAlgorithmLibrary.PC_STATS_PROCESS_LABEL); |
87 |
registerOutputParameter(ROI, List.class, PrincipalComponentsAlgorithmLibrary.PC_STATS_PROCESS_LABEL);
|
88 |
} |
89 |
|
90 |
@SuppressWarnings("unchecked") |
91 |
@Override
|
92 |
public void init() { |
93 |
store = (RasterDataStore)getParam(RASTER_STORE); |
94 |
bands = (boolean[])getParam(BANDS); |
95 |
path = getStringParam(PATH); |
96 |
rois = getParam(ROI) != null ? (List<ROI>)getParam(ROI) : null; |
97 |
} |
98 |
|
99 |
public double[] getAutoValors() { |
100 |
return autovalors;
|
101 |
} |
102 |
|
103 |
/**
|
104 |
* @return Autovectors Matrix
|
105 |
*/
|
106 |
public Matrix getAutoVectorMatrix() {
|
107 |
return autoVectorMatrix;
|
108 |
} |
109 |
|
110 |
/**
|
111 |
* @return Varianze-covarianze Matrix
|
112 |
*/
|
113 |
public Matrix getcoVarMatrix() {
|
114 |
return coVarMatrix;
|
115 |
} |
116 |
|
117 |
/**
|
118 |
* Calculate the Varianze-covarianze Matrix from bands of a Grid
|
119 |
* @throws InterruptedException
|
120 |
* @throws ProcessInterruptedException
|
121 |
*/
|
122 |
private double[][] covarianceOptimize() throws InterruptedException, ProcessInterruptedException { |
123 |
//NoData
|
124 |
NoData nd = RasterLocator.getManager().getDataStructFactory().createDefaultNoData(1, Buffer.TYPE_DOUBLE); |
125 |
|
126 |
double nodata = nd.getValue().doubleValue();
|
127 |
double dSum = 0; |
128 |
int iValues = 0; |
129 |
double[][] coV = new double[buffer.getBandCount()][buffer.getBandCount()]; |
130 |
|
131 |
double valorBandai = 0, valorBandaj = 0; |
132 |
maxminmean = new double[buffer.getBandCount()][4]; |
133 |
statsByBandCalculated = new boolean[buffer.getBandCount()]; |
134 |
for (int i = 0; i < statsByBandCalculated.length; i++) { |
135 |
statsByBandCalculated[i] = false;
|
136 |
} |
137 |
|
138 |
loadDataStoreStatistics(); |
139 |
|
140 |
int percent = 0; |
141 |
int maxCounter = 0; |
142 |
for (int i = buffer.getBandCount(); i >= 0; i--) { |
143 |
maxCounter += i; |
144 |
} |
145 |
|
146 |
Buffer roiBuf = createOutputBuffer(buffer.getWidth(), buffer.getHeight(), 1, Buffer.TYPE_BYTE); |
147 |
|
148 |
for (int i = 0; i < buffer.getBandCount(); i++) { |
149 |
double meanBufferToOperatei = getMeanValue(i, buffer);
|
150 |
for (int j = i; j < buffer.getBandCount(); j++) { |
151 |
iValues = 0;
|
152 |
dSum = 0;
|
153 |
double meanBufferToOperatej = getMeanValue(j, buffer);
|
154 |
|
155 |
for (int row = 0; row < buffer.getHeight(); row++) { |
156 |
for (int col = 0; col < buffer.getWidth(); col++) { |
157 |
if (i == 0 && j == 0) { |
158 |
roiBuf.setElem(row, col, 0, isInsideOfROI(col, row, rois, extentResult) ? (byte)1: (byte)0); |
159 |
} |
160 |
if(roiBuf.getElemByte(row, col, 0) == 1) { |
161 |
valorBandai = getData(buffer, row, col, i) - meanBufferToOperatei; |
162 |
valorBandaj = getData(buffer, row, col, j) - meanBufferToOperatej; |
163 |
|
164 |
dSum += valorBandai * valorBandaj; |
165 |
iValues++; |
166 |
} |
167 |
updatePercent((percent * 100) / maxCounter, 100); |
168 |
} |
169 |
} |
170 |
// Asigno el valor a la matriz
|
171 |
if (iValues > 1) |
172 |
coV[i][j] = dSum / (double)(iValues);
|
173 |
else
|
174 |
coV[i][j] = nodata; |
175 |
|
176 |
percent ++; |
177 |
updatePercent((percent * 100) / maxCounter, 100); |
178 |
} |
179 |
|
180 |
} |
181 |
|
182 |
roiBuf.dispose(); |
183 |
|
184 |
for (int i = 0; i < buffer.getBandCount(); i++) { |
185 |
for (int j = 0; j < buffer.getBandCount(); j++) { |
186 |
if(j < i)
|
187 |
coV[i][j] = coV[j][i]; |
188 |
} |
189 |
} |
190 |
|
191 |
return coV;
|
192 |
} |
193 |
|
194 |
// * Devuelve la media de los valores de los pixeles de un buffer incluidos en la zona de estudio
|
195 |
/**
|
196 |
* Get the mean of values of pixels from a buffer in a study zone
|
197 |
* @param bufferToOperate
|
198 |
* @return
|
199 |
*/
|
200 |
private double getMeanValue(int band, Buffer bufferToOperate) { |
201 |
if(statsByBandCalculated[band])
|
202 |
return maxminmean[band][2]; |
203 |
|
204 |
double mean = 0; |
205 |
double variance = 0; |
206 |
double min = Double.MAX_VALUE; |
207 |
double max = Double.MIN_VALUE; |
208 |
int valuesCounter = 0; |
209 |
for (int j = 0; j < bufferToOperate.getHeight(); j++) { |
210 |
for (int i = 0; i < bufferToOperate.getWidth(); i++) { |
211 |
if (isInsideOfROI(i, j, rois, extentResult)) {
|
212 |
double value = getData(bufferToOperate, j, i, band);
|
213 |
|
214 |
if(value < min) {
|
215 |
min = value; |
216 |
} |
217 |
|
218 |
if(value > max) {
|
219 |
max = value; |
220 |
} |
221 |
|
222 |
mean += value; |
223 |
variance += value * value; |
224 |
valuesCounter++; |
225 |
} |
226 |
} |
227 |
|
228 |
} |
229 |
mean = (mean / valuesCounter); |
230 |
variance = variance / (double) valuesCounter - mean * mean;
|
231 |
maxminmean[band][0] = max;
|
232 |
maxminmean[band][1] = min;
|
233 |
maxminmean[band][2] = mean;
|
234 |
maxminmean[band][3] = variance;
|
235 |
statsByBandCalculated[band] = true;
|
236 |
return mean;
|
237 |
} |
238 |
|
239 |
private void loadDataStoreStatistics() { |
240 |
if(store.getStatistics().isCalculated()) {
|
241 |
double[] max = store.getStatistics().getMax(); |
242 |
double[] min = store.getStatistics().getMin(); |
243 |
double[] mean = store.getStatistics().getMean(); |
244 |
double[] variance = store.getStatistics().getVariance(); |
245 |
int iBand = 0; |
246 |
for (int i = 0; i < mean.length; i++) { |
247 |
if(bands[i]) {
|
248 |
maxminmean[iBand][0] = max[i];
|
249 |
maxminmean[iBand][1] = min[i];
|
250 |
maxminmean[iBand][2] = mean[i];
|
251 |
maxminmean[iBand][3] = variance[i];
|
252 |
statsByBandCalculated[iBand] = true;
|
253 |
iBand ++; |
254 |
} |
255 |
} |
256 |
} |
257 |
} |
258 |
|
259 |
@Override
|
260 |
public void process() throws ProcessInterruptedException, ProcessException { |
261 |
insertLineLog(Messages.getText("preparing_buffers"));
|
262 |
|
263 |
try {
|
264 |
store = store.cloneDataStore(); |
265 |
} catch (CloneException e1) {
|
266 |
throw new ProcessException("Error cloning the sources", e1); |
267 |
} |
268 |
|
269 |
//Calculo de datos de entrada
|
270 |
extentResult = getExtentResult(rois, store); |
271 |
Rectangle2D sourcePxBBox = getSourcePxBox(extentResult, store);
|
272 |
|
273 |
buffer = createSourceBuffer(store, sourcePxBBox, bands); |
274 |
|
275 |
insertLineLog(Messages.getText("calc_stats"));
|
276 |
|
277 |
try {
|
278 |
coVar = covarianceOptimize(); |
279 |
} catch (InterruptedException e) { |
280 |
logger.error(Messages.getText("error_covarianza_optimize"), e);
|
281 |
} |
282 |
// Calculo de autovectores:
|
283 |
coVarMatrix = new Matrix(coVar);
|
284 |
EigenvalueDecomposition eigenvalueDecomp = new EigenvalueDecomposition(coVarMatrix);
|
285 |
//Autovectores
|
286 |
autoVectorMatrix = eigenvalueDecomp.getV(); |
287 |
|
288 |
// Autovalores
|
289 |
autovalors = eigenvalueDecomp.getRealEigenvalues(); |
290 |
|
291 |
PCStatsDataStructure resultStatistics = new PCStatsDataStructure(autoVectorMatrix, autovalors, coVarMatrix);
|
292 |
resultStatistics.setMaxMinMean(maxminmean); |
293 |
|
294 |
//Object[] dataROI = new Object[]{vectROI, (int)shift[0].getX(), (int)shift[0].getY()};
|
295 |
addOutputValue(STATS_RESULT, resultStatistics); |
296 |
addOutputValue(PATH, path); |
297 |
addOutputValue(BANDS, bands); |
298 |
addOutputValue(ROI, rois); |
299 |
try {
|
300 |
store.close(); |
301 |
} catch (CloseException e) {
|
302 |
throw new ProcessException("Error closing the sources", e); |
303 |
} |
304 |
} |
305 |
|
306 |
|
307 |
public String getTitle() { |
308 |
return null; |
309 |
} |
310 |
} |