Statistics
| Revision:

svn-gvsig-desktop / trunk / libraries / libRaster / src / org / gvsig / raster / grid / GridInterpolated.java @ 27361

History | View | Annotate | Download (14.6 KB)

1
/*******************************************************************************
2
                GridWrapperInterpolated
3
                Copyright (C) Victor Olaya
4
                
5
                This program is free software; you can redistribute it and/or modify
6
                it under the terms of the GNU General Public License as published by
7
                the Free Software Foundation; either version 2 of the License, or
8
                (at your option) any later version.
9

10
                This program is distributed in the hope that it will be useful,
11
                but WITHOUT ANY WARRANTY; without even the implied warranty of
12
                MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
13
                GNU General Public License for more details.
14

15
                You should have received a copy of the GNU General Public License
16
                along with this program; if not, write to the Free Software
17
                Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
18
*******************************************************************************/ 
19
package org.gvsig.raster.grid;
20

    
21
import org.gvsig.raster.buffer.BufferFactory;
22
import org.gvsig.raster.buffer.BufferInterpolation;
23
import org.gvsig.raster.buffer.RasterBuffer;
24
import org.gvsig.raster.buffer.RasterBufferInvalidAccessException;
25
import org.gvsig.raster.buffer.RasterBufferInvalidException;
26
import org.gvsig.raster.dataset.InvalidSetViewException;
27
import org.gvsig.raster.dataset.io.RasterDriverException;
28

    
29
/**
30
 * A grid wrapper that performs interpolation to calculate
31
 * cell values. This should be used when the window extent
32
 * does not 'fit' into the structure (coordinates and cellsize)
33
 * of the grid.
34
 * 
35
 * @author Victor Olaya
36
 */
37

    
38
public class GridInterpolated extends GridReader {
39

    
40
        public static final int INTERPOLATION_NearestNeighbour = BufferInterpolation.INTERPOLATION_NearestNeighbour;
41
        public static final int INTERPOLATION_Bilinear = BufferInterpolation.INTERPOLATION_Bilinear;
42
        public static final int INTERPOLATION_InverseDistance = BufferInterpolation.INTERPOLATION_InverseDistance;
43
        public static final int INTERPOLATION_BicubicSpline = BufferInterpolation.INTERPOLATION_BicubicSpline;
44
        public static final int INTERPOLATION_BSpline = BufferInterpolation.INTERPOLATION_BSpline;
45
        
46
        double m_dXMin , m_dYMax; //coordinates of the layer, not the window
47
        double m_dCellSizeX; //cellsize of the layer, not the window
48
        double m_dCellSizeY; //cellsize of the layer, not the window
49
        int m_iInterpolationMethod = INTERPOLATION_BSpline;
50

    
51
        /**
52
         * Crea un objeto lector a partir de un buffer de datos y el extent de la extensi?n
53
         * completa y de la ventana accedida.
54
         * @param rb Buffer de datos
55
         * @param layerExtent extent de la capa completa
56
         * @param windowExtent Extent
57
         * @param bands N?mero de bandas del origen
58
         */
59
        public GridInterpolated(        RasterBuffer rb, 
60
                                                                GridExtent layerExtent,
61
                                                                GridExtent windowExtent,
62
                                                                int[] bands){
63
                super(rb, layerExtent, windowExtent, bands);
64
                m_dXMin = layerExtent.minX();
65
                m_dYMax = layerExtent.maxY();
66
                
67
                m_dCellSizeX = layerExtent.getCellSizeX();
68
                m_dCellSizeY = layerExtent.getCellSizeY();
69
        }
70
        
71
        /**
72
         * Crea un objeto lector a partir de una fuente de datos y el extent de la extensi?n
73
         * completa y de la ventana accedida.
74
         * @param ds Fuente de datos
75
         * @param layerExtent extent de la capa completa
76
         * @param windowExtent Extent
77
         * @param bands N?mero de bandas del origen
78
         */
79
        public GridInterpolated(        BufferFactory ds, 
80
                                                                GridExtent layerExtent,
81
                                                                GridExtent windowExtent,
82
                                                                int[] bands){
83
                super(ds, layerExtent, windowExtent, bands);
84
                init();
85
        }
86
        
87
        /**
88
         * Inicializaci?n de la fuente de datos y carga del buffer.
89
         */
90
        private void init(){
91
                m_dXMin = bufferFactory.getDataSource().getExtent().minX();
92
                m_dYMax = bufferFactory.getDataSource().getExtent().maxY();
93
                
94
                m_dCellSizeX = bufferFactory.getDataSource().getAffineTransform(0).getScaleX();
95
                m_dCellSizeY = bufferFactory.getDataSource().getAffineTransform(0).getScaleY();
96
                
97
                bufferFactory.setDrawableBands(bands);
98
        }
99
                
100
        /**
101
         * Asigna el m?todo de interpolaci?n 
102
         * @param iMethod
103
         */
104
        public void setInterpolationMethod(int iMethod){
105
                m_iInterpolationMethod = iMethod;
106
        }
107

    
108
        /*
109
         *  (non-Javadoc)
110
         * @see org.gvsig.fmap.grid.GridReader#getCellValueAsByte(int, int)
111
         */
112
        public byte getCellValueAsByte(int x, int y) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException, InterruptedException {
113
                return  (byte)_getValueAt(x, y);
114
        }
115
        
116
        /*
117
         *  (non-Javadoc)
118
         * @see org.gvsig.fmap.grid.GridReader#getCellValueAsShort(int, int)
119
         */
120
        public short getCellValueAsShort(int x, int y) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException {
121
                return (short) _getValueAt(x, y);
122
        }
123
        
124
        /*
125
         *  (non-Javadoc)
126
         * @see org.gvsig.fmap.grid.GridReader#getCellValueAsInt(int, int)
127
         */
128
        public int getCellValueAsInt(int x, int y) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException  {                
129
                return (int) _getValueAt(x, y);
130
        }
131

    
132
        /*
133
         *  (non-Javadoc)
134
         * @see org.gvsig.fmap.grid.GridReader#getCellValueAsFloat(int, int)
135
         */
136
        public float getCellValueAsFloat(int x, int y)  throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException {
137
                return (float) _getValueAt(x, y);
138
        }
139
        
140
        /*
141
         *  (non-Javadoc)
142
         * @see org.gvsig.fmap.grid.GridReader#getCellValueAsDouble(int, int)
143
         */
144
        public double getCellValueAsDouble(int x, int y) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException {
145
                return _getValueAt(x, y);
146
        }
147
        
148
        /*
149
         *  (non-Javadoc)
150
         * @see org.gvsig.fmap.grid.GridReader#getCellValueAsDouble(int, int)
151
         */
152
        public double getCellValue(int x, int y) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException {
153
                return _getValueAt(x, y);
154
        }
155
        
156
        private double _getValueAt(int x, int y) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException {
157
                double dX = windowExtent.minX() + windowExtent.getCellSize() * (x + 0.5);
158
                double dY = windowExtent.maxY() - windowExtent.getCellSize() * (y + 0.5);
159
                return getValueAt(dX, dY);
160
        }
161

    
162
        
163
        /** Devuelve el valor interpolado para unas coordenadas de la imagen no enteras
164
         *  la distancia se toma al vertice superior izquierdo del pixel
165
         * @throws InterruptedException 
166
         * */
167
        public double _getValueAt(double x, double y) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException {
168
                int x_int, y_int;
169
                x_int= (int)Math.floor(x);
170
                y_int= (int)Math.floor(y);
171
                double dX = windowExtent.minX() + windowExtent.getCellSize() * (x_int + (x-x_int));
172
                double dY = windowExtent.maxY() - windowExtent.getCellSize() * (y_int + (y-y_int));
173
                return getValueAt(dX, dY);
174
        }
175
        
176
        
177
        private double getValueAt(double xPosition, double yPosition) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException {
178
                int        x, y;
179
                double        dx, dy;
180
                double dValue;
181

    
182
                x        = (int) Math.floor(xPosition = (xPosition - m_dXMin) / m_dCellSizeX);
183
                y        = (int) Math.floor(yPosition = (yPosition - m_dYMax) / m_dCellSizeY);
184

    
185
                dValue = _getCellValueInLayerCoords(x, y);
186
                
187
                if(        !isNoDataValue(dValue) ) {
188

    
189
                        dx        = xPosition - x;
190
                        dy        = yPosition - y;
191
                        
192
                        switch( m_iInterpolationMethod ) {
193
                        case INTERPOLATION_NearestNeighbour:
194
                                dValue = _getValueNearestNeighbour (x, y, dx, dy);
195
                                break;
196

    
197
                        case INTERPOLATION_Bilinear:
198
                                dValue        = _getValueBiLinear (x, y, dx, dy);
199
                                break;
200

    
201
                        case INTERPOLATION_InverseDistance:
202
                                dValue        = _getValueInverseDistance(x, y, dx, dy);
203
                                break;
204

    
205
                        case INTERPOLATION_BicubicSpline:
206
                                dValue        = _getValueBiCubicSpline (x, y, dx, dy);
207
                                break;
208

    
209
                        case INTERPOLATION_BSpline:
210
                                dValue        = _getValueBSpline (x, y, dx, dy);
211
                                break;
212
                        }
213
                } else {
214
                        dValue = getNoDataValue();
215
                }
216

    
217
                return dValue;
218
        }
219

    
220
        private double _getValueNearestNeighbour(int x, int y, double dx, double dy) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException  {
221
                x += (int)(0.5 + dx);
222
                y += (int)(0.5 + dy);
223

    
224
                return _getCellValueInLayerCoords(x, y);
225
        }
226

    
227
        private  double _getValueBiLinear(int x, int y, double dx, double dy) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException  {
228
                double        z = 0.0, n = 0.0, d;
229
                double dValue;
230
                
231
                dValue = _getCellValueInLayerCoords(x, y);
232
                if (!isNoDataValue(dValue)){
233
                         d = (1.0 - dx) * (1.0 - dy);
234
                         z += d * dValue; 
235
                         n += d;
236
                }
237
                
238
                dValue = _getCellValueInLayerCoords(x + 1, y);
239
                if (!isNoDataValue(dValue)){
240
                         d = (dx) * (1.0 - dy);
241
                         z += d * dValue; 
242
                         n += d;
243
                }
244
                
245
                dValue = _getCellValueInLayerCoords(x, y + 1);
246
                if (!isNoDataValue(dValue)){
247
                         d = (1.0 - dx) * (dy);
248
                         z += d * dValue; 
249
                         n += d;
250
                }
251
                
252
                dValue = _getCellValueInLayerCoords(x + 1, y + 1);
253
                if (!isNoDataValue(dValue)){
254
                         d = (dx) * (dy);
255
                         z += d * dValue; 
256
                         n += d;
257
                }
258

    
259
                if( n > 0.0 ){
260
                        return( z / n );
261
                }
262

    
263
                return( getNoDataValue() );
264
        }
265

    
266
        private double _getValueInverseDistance(int x, int y, double dx, double dy) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException {
267
                double        z = 0.0, n = 0.0, d;
268
                double dValue;
269
                
270
                if( dx > 0.0 || dy > 0.0 ){
271

    
272
                        dValue = _getCellValueInLayerCoords(x, y);
273
                        if (!isNoDataValue(dValue)){
274
                                d = 1.0 / Math.sqrt(dx*dx + dy*dy); 
275
                                z += d * dValue; 
276
                                n += d;
277
                        }
278
                        
279
                        dValue = _getCellValueInLayerCoords(x + 1, y);
280
                        if (!isNoDataValue(dValue)){
281
                                d = 1.0 / Math.sqrt((1.0-dx)*(1.0-dx) + dy*dy); 
282
                                z += d * dValue; 
283
                                n += d;
284
                        }
285
                        
286
                        dValue = _getCellValueInLayerCoords(x, y + 1);
287
                        if (!isNoDataValue(dValue)){
288
                                d = 1.0 / Math.sqrt(dx*dx + (1.0-dy)*(1.0-dy)); 
289
                                z += d * dValue; 
290
                                n += d;
291
                        }
292
                        
293
                        dValue = _getCellValueInLayerCoords(x + 1, y + 1);
294
                        if (!isNoDataValue(dValue)){
295
                                d = 1.0 / Math.sqrt((1.0-dx)*(1.0-dx) + (1.0-dy)*(1.0-dy)); 
296
                                z += d * dValue; 
297
                                n += d;
298
                        }
299

    
300
                        if( n > 0.0 ){
301
                                return( z / n );
302
                        }
303
                }else
304
                        return( _getCellValueInLayerCoords(x, y) );
305
                
306
                return( getNoDataValue());
307
        }
308

    
309
        private double _getValueBiCubicSpline(int x, int y, double dx, double dy) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException  {
310
                int                i;
311
                double        a0, a2, a3, b1, b2, b3, c[], z_xy[][];
312

    
313
                c = new double[4];
314
                z_xy = new double[4][4];
315
                
316
                if( _get4x4Submatrix(x, y, z_xy) ){
317
                        
318
                        for(i=0; i<4; i++){
319
                                a0                = z_xy[0][i] - z_xy[1][i];
320
                                a2                = z_xy[2][i] - z_xy[1][i];
321
                                a3                = z_xy[3][i] - z_xy[1][i];
322

    
323
                                b1                = -a0 / 3.0 + a2       - a3 / 6.0;
324
                                b2                =  a0 / 2.0 + a2 / 2.0;
325
                                b3                = -a0 / 6.0 - a2 / 2.0 + a3 / 6.0;
326

    
327
                                c[i]        = z_xy[1][i] + b1 * dx + b2 * dx*dx + b3 * dx*dx*dx;
328
                        }
329

    
330
                        a0                = c[0] - c[1];
331
                        a2                = c[2] - c[1];
332
                        a3                = c[3] - c[1];
333

    
334
                        b1                = -a0 / 3.0 + a2       - a3 / 6.0;
335
                        b2                =  a0 / 2.0 + a2 / 2.0;
336
                        b3                = -a0 / 6.0 - a2 / 2.0 + a3 / 6.0;
337

    
338
                        return( c[1] + b1 * dy + b2 * dy*dy + b3 * dy*dy*dy );
339
                }
340

    
341
                return( _getValueBiLinear(x, y, dx, dy) );
342
        }
343

    
344
        private double _getValueBSpline(int x, int y, double dx, double dy) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException {
345
                int                i, ix, iy;
346
                double        z, px, py, Rx[], Ry[], z_xy[][];
347
                
348
                Rx = new double[4];
349
                Ry = new double[4];
350
                z_xy = new double [4][4];
351

    
352
                if( _get4x4Submatrix(x, y, z_xy) ) {
353
                        for(i=0, px=-1.0-dx, py=-1.0-dy; i<4; i++, px++, py++) {
354
                                Rx[i]        = 0.0;
355
                                Ry[i]        = 0.0;
356

    
357
                                if( (z = px + 2.0) > 0.0 )
358
                                        Rx[i]        +=        z*z*z;
359
                                if( (z = px + 1.0) > 0.0 )
360
                                        Rx[i]        += -4.0 * z*z*z;
361
                                if( (z = px + 0.0) > 0.0 )
362
                                        Rx[i]        +=  6.0 * z*z*z;
363
                                if( (z = px - 1.0) > 0.0 )
364
                                        Rx[i]        += -4.0 * z*z*z;
365
                                if( (z = py + 2.0) > 0.0 )
366
                                        Ry[i]        +=        z*z*z;
367
                                if( (z = py + 1.0) > 0.0 )
368
                                        Ry[i]        += -4.0 * z*z*z;
369
                                if( (z = py + 0.0) > 0.0 )
370
                                        Ry[i]        +=  6.0 * z*z*z;
371
                                if( (z = py - 1.0) > 0.0 )
372
                                        Ry[i]        += -4.0 * z*z*z;
373

    
374
                                Rx[i]        /= 6.0;
375
                                Ry[i]        /= 6.0;
376
                        }
377

    
378
                        for(iy=0, z=0.0; iy<4; iy++) {
379
                                for(ix=0; ix<4; ix++){
380
                                        z        += z_xy[ix][iy] * Rx[ix] * Ry[iy];
381
                                }
382
                        }
383
                        return( z );
384
                }
385
                return( _getValueBiLinear(x, y, dx, dy) );
386
        }
387

    
388
        private boolean _get4x4Submatrix(int x, int y, double z_xy[][]) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException {
389
                int        ix, iy, px, py;
390
                double dValue;
391

    
392
                for(iy=0, py=y-1; iy<4; iy++, py++) {
393
                        for(ix=0, px=x-1; ix<4; ix++, px++) {
394
                                dValue = _getCellValueInLayerCoords(px, py);
395
                                if (isNoDataValue(dValue)) {
396
                                        return false;
397
                                } else {
398
                                        z_xy[ix][iy] = dValue;
399
                                }
400
                        }
401
                }
402
                return( true );
403
        }
404

    
405
        private double _getCellValueInLayerCoords(int x, int y) throws RasterBufferInvalidAccessException, RasterBufferInvalidException, InterruptedException {
406
                int w, h;
407
                int pX = 0, pY = 0;
408
                if (bufferFactory != null) {
409
                        w = bufferFactory.getSourceWidth();
410
                        h = bufferFactory.getSourceHeight();
411
                        try {
412
                                if (y >= h || x >= w || x < 0 || y < 0)
413
                                        return getNoDataValue();
414
                                bufferFactory.setAreaOfInterest(x, y, 1, 1);
415
                        } catch (InvalidSetViewException e) {
416
                                throw new RasterBufferInvalidAccessException("");
417
                        } catch (RasterDriverException e) {
418
                                throw new RasterBufferInvalidException("");
419
                        } catch (InterruptedException e) {
420
                                // La cancelaci?n de la lectura de un pixel no es significativa
421
                        }
422
                        rasterBuf = (RasterBuffer) bufferFactory.getRasterBuf();
423
                } else {
424
                        w = rasterBuf.getWidth();
425
                        h = rasterBuf.getHeight();
426
                        pX = x;
427
                        pY = y;
428
                }
429

    
430
                if (x >= 0 && y >= 0 && x < w && y < h) {
431
                        switch (dataType) {
432
                                case RasterBuffer.TYPE_DOUBLE:
433
                                        return rasterBuf.getElemDouble(pY, pX, bandToOperate);
434
                                case RasterBuffer.TYPE_INT:
435
                                        return (double) rasterBuf.getElemInt(pY, pX, bandToOperate);
436
                                case RasterBuffer.TYPE_FLOAT:
437
                                        return (double) rasterBuf.getElemFloat(pY, pX, bandToOperate);
438
                                case RasterBuffer.TYPE_BYTE:
439
                                        return ((byte) (rasterBuf.getElemByte(pY, pX, bandToOperate)) & 0xff);
440
                                case RasterBuffer.TYPE_SHORT:
441
                                case RasterBuffer.TYPE_USHORT:
442
                                        return (double) rasterBuf.getElemShort(pY, pX, bandToOperate);
443
                        }
444
                }
445
                return getNoDataValue();
446
        }
447

    
448
        public byte[] getBandsValuesAsByte(int x, int y) {
449
                //TODO: FUNCIONALIDAD: getBandsValuesAsByte:Obtener los valores interpolados en todas las bandas
450
                return null;
451
        }
452

    
453
        public short[] getBandsValuesAsShort(int x, int y) {
454
                //TODO: FUNCIONALIDAD: getBandsValuesAsShort:Obtener los valores interpolados en todas las bandas
455
                return null;
456
        }
457

    
458
        public int[] getBandsValuesAsInt(int x, int y) {
459
                //TODO: FUNCIONALIDAD: getBandsValuesAsInt:Obtener los valores interpolados en todas las bandas
460
                return null;
461
        }
462

    
463
        public float[] getBandsValuesAsFloat(int x, int y) {
464
                //TODO: FUNCIONALIDAD: getBandsValuesAsFloat:Obtener los valores interpolados en todas las bandas
465
                return null;
466
        }
467

    
468
        public double[] getBandsValuesAsDouble(int x, int y) {
469
                //TODO: FUNCIONALIDAD: getBandsValuesAsDouble:Obtener los valores interpolados en todas las bandas
470
                return null;
471
        }
472
}