Selecting and extracting data from TIFF files

Question

I have some tiff files (2 for example) and a shapefile. The real format of the tiff files is geoTiff.

I want to get the points inside the shape boundary. Each point should have 4 parameters: geom(coordinates), red, green, infrared (we have tiff files with such info).

At the end of the operation I want to get a vectorized shape with points.

Could you please give me some advice on how can I do that programmatically. Are there any libs/tools? It is all in the context of geoserver/postgresql/geomajas framework.

Upd: now I have more information. My tiffs will be one of the Landsat combination: http://gif.berkeley.edu/documents/Landsat%20Band%20Information.pdf

The main question now is: to get RGB channels info in exact point, should I parse image points on RGB,or I can get this info from geoTiffMetadata (there is http://download.deegree.org/deegree3/nightly/core/javadoc/org/deegree/coverage/raster/io/imageio/GeoTiffIIOMetadataAdapter.html ). I didn`t find really useful info while debugging the code and calling the methods of GeoTiffIIOMetadataAdapter.

Answer 1

I tried different approaches(also with postgis db, 3rd party libraries), but stopped at geotools java library. So, I`ll describe how to get points from tiff file that are inside shape polygon boundary. Totally it gets about 15 seconds to do this operation for tiff file of 66 million pixels(200 mb).

public class Run {
    public static void main(String[] args) throws IOException, TransformException, FactoryException, SchemaException {
        long time1 = System.currentTimeMillis();
        Run run = new Run();

        Feature shapeFeature = run.getPolygonFromShape();
        List<MyRGB> myRGBs = run.getOnlyPointsFromShapeWithFragmentation(shapeFeature, run.getTiffMetadata());
        long time2 = System.currentTimeMillis();
        System.out.println("Wasted time = " + (time2 - time1) + " ms");
    }

    private Feature getPolygonFromShape() {
        String filename = "d:\\work\\myShapeFile.shp";
        File file = new File(filename);
        try {
            Map<String, URL> connect = new HashMap<String, URL>();
            connect.put("url", file.toURL());
            DataStore dataStore = DataStoreFinder.getDataStore(connect);
            String[] typeNames = dataStore.getTypeNames();
            String typeName = typeNames[0];
            System.out.println("Reading content " + typeName);
            FeatureSource featureSource = dataStore.getFeatureSource(typeName);
            FeatureCollection collection = featureSource.getFeatures();
            FeatureIterator iterator = collection.features();
            Feature feature = null;
            try {
                while (iterator.hasNext()) {
                    feature = iterator.next();
                }
                return feature;
            } finally {
                iterator.close();
            }

        } catch (Exception e) {
            e.printStackTrace();
        }
        return null;
    }
    public GeoTiffIIOMetadataAdapter getTiffMetadata() throws IOException {
        String filename = "d:\\work\\myTiffFile.tif";
        FileImageInputStream f = new FileImageInputStream(
                new RandomAccessFile(filename, "r"));
        // Look through ImageIO readers
        Iterator iter = ImageIO.getImageReaders(f);
        IIOMetadata imdata = null;
        GeoTiffIIOMetadataAdapter geo_data;
        GeoTiffIIOMetadataAdapter metadata = null;
        while (iter.hasNext() && imdata == null) {
            ImageReader reader = (ImageReader) iter.next();
            reader.setInput(f, true);
            String reader_name = reader.getFormatName().toLowerCase();
            if (reader_name.equalsIgnoreCase("tif")) {
                // Get Image metadata
                imdata = reader.getImageMetadata(0);
                geo_data = new GeoTiffIIOMetadataAdapter(imdata);
                if (geo_data.getGeoKeyDirectoryVersion() == 1) {
                    metadata = geo_data;
                }
            }
        }
        f.close();
        return metadata;
    }

    private List<MyRGB> getOnlyPointsFromShape(Feature shapeFeature, GeoTiffIIOMetadataAdapter tiffMetadata) throws IOException, FactoryException, TransformException {
        BufferedImage bi = ImageIO.read(new File("d:\\work\\myTiffFile.tif"));       
        CoordinateReferenceSystem sourceCRS = CRS.decode("EPSG:4326", true);
        CoordinateReferenceSystem targetCRS = CRS.decode("EPSG:32612");
        MathTransform transformToEpsg4326 = CRS.findMathTransform(targetCRS, sourceCRS);

        double lon = tiffMetadata.getModelTiePoints()[3];// 175784.99999999997 - X
        double lat = tiffMetadata.getModelTiePoints()[4];// 5842215.0   - Y
        GeometryFactory geometryFactory = new GeometryFactory(new PrecisionModel(), 32612);
        Geometry pointUpperLeft = geometryFactory.createPoint(new Coordinate(lon, lat));

        int[] pixel;
        LinkedList<MyRGB> myRGBList = new LinkedList<MyRGB>();
        MyRGB myRGB = null;
        double startLeftUpperTiffX = JTS.transform(pointUpperLeft, transformToEpsg4326).getCoordinate().x;// -115.79200437315563 - X
        double startTiffX = startLeftUpperTiffX;
        double startLeftUpperTiffY = JTS.transform(pointUpperLeft, transformToEpsg4326).getCoordinate().y;// 52.632842905984255   - Y
        double stepX = tiffMetadata.getModelPixelScales()[0];//30.0
        double stepY = tiffMetadata.getModelPixelScales()[1];//30.0
        double radius = 6378137;//earth radius

        MultiPolygon mp = (MultiPolygon) ((SimpleFeatureImpl) shapeFeature).getDefaultGeometry();
        Geometry tmpPoint = null;
        for (int y = 0; y < bi.getHeight(); y++) {
            startLeftUpperTiffX = startTiffX;
            for (int x = 0; x < bi.getWidth(); x++) {
                //create point
                tmpPoint = geometryFactory.createPoint(new Coordinate(startLeftUpperTiffX, startLeftUpperTiffY));
                if (mp.contains(tmpPoint)) {
                    pixel = bi.getRaster().getPixel(x, y, new int[3]);
                    myRGB = new MyRGB(pixel[0], pixel[1], pixel[2]);
                    myRGB.setLatitude(startLeftUpperTiffY);
                    myRGB.setLongitude(startLeftUpperTiffX);
                    myRGBList.add(myRGB);
                }
                double dLon = stepX / (radius * Math.cos(Math.PI * startLeftUpperTiffY / 180));//coordinate, and do plus meters, get new coordinate
                startLeftUpperTiffX += dLon * 180 / Math.PI;
            }
            startLeftUpperTiffY -= (180 / Math.PI) * (stepY / radius);   //have coordinate, and do plus meters, get new coordinate
        }
        System.out.println("myRGBList.size() = " + myRGBList.size());
        System.out.println("startLeftUpperTiffX = " + startLeftUpperTiffX);
        System.out.println("startLeftUpperTiffY = " + startLeftUpperTiffY);
        return myRGBList;
    }

 public class MyRGB {
        int red;
        int green;
        int infrared;
        double latitude;
        double longitude;
//.... constructors, getters and setters
}