[GeoTools]GridCoverage Statistics

GeoTools에서 GridCoverage(Raster)의 단순통계를 계산하는 Java 샘플 코드입니다.
 - JAI의 RectIter 인터페이스를 이용하여 계산
 - AWT의 RenderedImage인터페이스를 이용하여 계산
 - 단일 밴드에 대한 샘플입니다.

▣ Code Snippet

import java.awt.image.Raster;
import java.awt.image.RenderedImage;
import java.text.DecimalFormat;

import javax.media.jai.PlanarImage;
import javax.media.jai.iterator.RectIter;
import javax.media.jai.iterator.RectIterFactory;

import org.geotools.coverage.grid.GridCoverage2D;
import org.geotools.coverage.grid.GridEnvelope2D;

■ JAI(Java Advanced Imaging): RectIter Interface

static void calculateStatistics1(GridCoverage2D inputGc) {
    int    count = 0;
    double minVal = Double.MAX_VALUE;
    double maxVal = Double.MIN_VALUE;
    double sumOfVals = 0.0;
    double sumOfSqrs = 0.0;
    
    double noDataValue = (Double) inputGc.getProperty("GC_NODATA");
    
    PlanarImage renderedImage = (PlanarImage) inputGc.getRenderedImage();
    java.awt.Rectangle bounds = renderedImage.getBounds();
    
    RectIter rectIter = RectIterFactory.create(renderedImage, bounds);
    rectIter.startLines();
    while (!rectIter.finishedLines()) {
        rectIter.startPixels();
        while (!rectIter.finishedPixels()) {
            double sampleVal = rectIter.getSampleDouble(0);

            if (Double.compare(sampleVal, noDataValue) != 0) {
                count++;
                sumOfVals += sampleVal;
                minVal = Math.min(minVal, sampleVal);
                maxVal = Math.max(maxVal, sampleVal);
                sumOfSqrs += Math.pow(sampleVal, 2);
            }

            rectIter.nextPixel();
        }
        rectIter.nextLine();
    }

    // Print output
    final DecimalFormat DF = new DecimalFormat("#.###############");
    System.out.println("Count: " + count);
    System.out.println("NoData: " + DF.format(noDataValue));
    System.out.println("Minimum: " + DF.format(minVal));
    System.out.println("Maximum: " + DF.format(maxVal));
    System.out.println("Sum: " + DF.format(sumOfVals));
    System.out.println("Mean: " + DF.format(sumOfVals / count));

    // Population Standard Deviation
    double variance = (sumOfSqrs - Math.pow(sumOfVals, 2) / count) / count;
    System.out.println("Variance: " + DF.format(variance));
    System.out.println("Standard Deviation: " + DF.format(Math.sqrt(variance)));
}

■ Java AWT: RenderedImage Interface

static void calculateStatistics2(GridCoverage2D inputGc) {
    int    count = 0;
    double minVal = Double.MAX_VALUE;
    double maxVal = Double.MIN_VALUE;
    double sumOfVals = 0.0;
    double sumOfSqrs = 0.0;

    double noDataValue = (Double) inputGc.getProperty("GC_NODATA");

    GridEnvelope2D ge2D = inputGc.getGridGeometry().getGridRange2D();
    RenderedImage renderedImage = inputGc.getRenderedImage();

    for (int tileX = 0; tileX < renderedImage.getNumXTiles(); tileX++) {
        for (int tileY = 0; tileY < renderedImage.getNumYTiles(); tileY++) {
            Raster tileRs = renderedImage.getTile(tileX, tileY);

            java.awt.Rectangle bounds = tileRs.getBounds();
            for (int dy = bounds.y, drow = 0; drow < bounds.height; dy++, drow++) {
                for (int dx = bounds.x, dcol = 0; dcol < bounds.width; dx++, dcol++) {
                    if (ge2D.contains(dx, dy)) {
                        double sampleVal = tileRs.getSampleDouble(dx, dy, 0);

                        if (Double.compare(sampleVal, noDataValue) != 0) {
                            count++;
                            sumOfVals += sampleVal;
                            minVal = Math.min(minVal, sampleVal);
                            maxVal = Math.max(maxVal, sampleVal);
                            sumOfSqrs += Math.pow(sampleVal, 2);
                        }                         
                    }
                }
            }
        }
    }

    // Print output
    final DecimalFormat DF = new DecimalFormat("#.###############");
    System.out.println("Count: " + count);
    System.out.println("NoData: " + DF.format(noDataValue));
    System.out.println("Minimum: " + DF.format(minVal));
    System.out.println("Maximum: " + DF.format(maxVal));
    System.out.println("Sum: " + DF.format(sumOfVals));
    System.out.println("Mean: " + DF.format(sumOfVals / count));

    // Population Standard Deviation
    double variance = (sumOfSqrs - Math.pow(sumOfVals, 2) / count) / count;
    System.out.println("Variance: " + DF.format(variance));
    System.out.println("Standard Deviation: " + DF.format(Math.sqrt(variance)));
}

■ Usage

protected static final Logger LOGGER = Logging.getLogger(GridCoverageStatistics.class);

public static void main(String[] args) {
    String gridFile = "C:/dev/GeoServer 2.1.1/data_dir/data/sf/sfdem.tif";
    
    try {
        AbstractGridFormat format = GridFormatFinder.findFormat(gridFile);
        AbstractGridCoverage2DReader reader = format.getReader(gridFile);
        
        GridCoverage2D inputGc = reader.read(null);

        calculateStatistics1(inputGc);
        calculateStatistics2(inputGc);
    } catch (IllegalArgumentException e) {
        LOGGER.log(Level.FINER, e.getMessage(), e);
    } catch (IOException e) {
        LOGGER.log(Level.FINER, e.getMessage(), e);
    }
}

■ Output

* Statistics1
Count: 292317
NoData: -9999999933815813000000000000000000000
Minimum: 1066
Maximum: 1840
Sum: 395700554
Mean: 1353.6693178980354
Variance: 31343.402316475993
Standard Deviation: 177.04067983510456

* Statistics2
Count: 292317
NoData: -9999999933815813000000000000000000000
Minimum: 1066
Maximum: 1840
Sum: 395700554
Mean: 1353.6693178980354
Variance: 31343.402316475993
Standard Deviation: 177.04067983510456

▣ 참고
 - GeoTools Coverage
 - JAI RectIter Interface
 - Java AWT RenderedImage