Java 调用gdal API(二)——栅格裁剪

gdal可以说是GIS数据处理比较好的工具之一，虽然也提供了Java API，但是官方文档确实太过简单，用起来确实太难受，每次都需要去参考对应的C++api，然后在对应使用。
因此小编决定从这篇文章开始，将自己以前以及今后用到的API，都记录下，以帮助更多从事java GIS服务开发者。
今天主要来记录栅格数据的裁剪功能，gdal中提供了多种裁剪方法，今天主要介绍以下两种：
第一中是通过gdal自身的Warp方法，其实现有以下几个接口；
public static int Warp(Dataset dstDS, Dataset[] object_list_count, WarpOptions warpAppOptions, ProgressCallback callback)；public static int Warp(Dataset dstDS, Dataset[] object_list_count, WarpOptions warpAppOptions)；public static Dataset Warp(String dest, Dataset[] object_list_count, WarpOptions warpAppOptions, ProgressCallback callback)；public static Dataset Warp(String dest, Dataset[] object_list_count, WarpOptions warpAppOptions)；其实这个几个接口基本差不多
1、dstDs/dest 其实就是目标数据，两种方式，要么传入目标路径，要么是目标Dataset
2、object_list_count 这个是源数据的Dataset 数据
3、warpAppOptions 这个参数很重要，是我们裁剪需要传递的重要参数，后面代码会专门解释
4、callback 处理进度
下面我们结合最后一个方法，来处理一个栅格裁剪
 

1. private static void cutTif() throws FactoryException, IOException, ParseException {
2.         Dataset dataSet = gdal.Open("E:\\数据\\影像数据\\tiff\\xian18.tif");
3.         //通过最大最小经纬度生成一个面矢量
4.         createCutShapeRegion(108.583740,33.928528,108.883740,34.260242);
5.         Vector vector = new Vector();
6.         vector.add("-cutline");
7.         vector.add("D:\\data\\cut2.shp");
8.         vector.add("-crop_to_cutline");
9.         Dataset ds = gdal.Warp("D:\\data\\xian18_1.tif",new Dataset[]{dataSet},new WarpOptions(vector));
10.         if (ds != null){
11.             System.out.println("成功！");
12.         }
13.     }

复制代码

1. createCutShapeRegion是用通过bbox来生成shape数据的，我这边是用矩形来裁剪的，当然你也可以用多边形，通过传入多边形的顶点来生层一个面矢量，记住一定是面矢量，<br>

复制代码

1. D:\\data\\cut2.shp是createCutShapeRegion生成的矢量文件，最主要的两个参数一个是-cutline 一个是-crop_to_cutline，其中-cutline后面跟裁剪的矢量文件，-crop_to_cutline不需要具体值，这些具体参数可参考gdal的warp的参数介绍，其实和安装版gdalwarp的参数是基本一致的。<br><br>但是这里说明下，并不是gdalwarp的所有参数都在gdal的java api中进行了适配，就比如gdalwarp这个命令中除了上面使用shp作为裁剪之外，还有一个-te参数，这个参数后面跟一个经纬度的范围<br>主要作用就是将原始数据裁剪为制定的范围，就比如我们使用安装的gdal中可以这样调用裁剪<br>gdalwarp.exe -te 108.583740 33.928528 108.883740 34.260242  E:\\数据\\影像数据\\tiff\\xian18.tif D:\\data\\xian18_2.tif<br>然而这个-te参数如果使用代码来调用，就会报没有-te这个option，就比如下面的这个例子其实是跑不通的<br>

复制代码

1. public static void main(String[] args) throws NoSuchFieldException, IllegalAccessException, FactoryException, IOException, ParseException {
2.         gdal.AllRegister();
3.         Dataset dataSet = gdal.Open("E:\\数据\\影像数据\\tiff\\xian18.tif");
4.         Vector vector = new Vector();
5.         vector.add("-te");
6.         vector.add("108.583740 33.928528 108.883740 34.260242"); //裁剪范围
7.         vector.add("-crop_to_cutline");
8.         Dataset ds = gdal.Warp("D:\\data\\xian18_1.tif",new Dataset[]{dataSet},new WarpOptions(vector));
9.     }

复制代码

[code][/code]上述代码是跑不通的，虽然安装板的gdalwarp可以使用-te来裁剪，而且效率比上面使用shp要更快，但是java代码的api却没有集成这个-te参数，当然我不知道是我调用的问题，还是确实没有适配
如果有了解的朋友，麻烦告知一声。
 
除了上述使用gdalwarp的接口外，我们还可以通过使用driver来快速裁剪出一个tif,具体代码如下
[code]private static void clipRaster(Dataset dataset, Envelope boundingBox){        //获取当前数据范围        ReferencedEnvelope originalEnvelope = (ReferencedEnvelope) GdalDataAnalysis.readBboxFromFile(dataset);        //判断是否相交，不想交则直接返回        if (!originalEnvelope.intersects(boundingBox)){            return;        }        //求请求与原始数据的范围        ReferencedEnvelope interEnvelope = originalEnvelope.intersection(boundingBox);        //原始数据上偏移位置        int offset_x = (int) ((interEnvelope.getMinX() - originalEnvelope.getMinX()) / GdalDataAnalysis.getResoulution(dataset,true));        int offset_y = (int) ((interEnvelope.getMaxY() - originalEnvelope.getMaxY()) / GdalDataAnalysis.getResoulution(dataset,false));        //需要多少个原始像素        int block_xsize = (int) ((interEnvelope.getMaxX() - interEnvelope.getMinX()) / GdalDataAnalysis.getResoulution(dataset,true));        int block_ysize = (int)Math.abs((interEnvelope.getMaxY() - interEnvelope.getMinY()) / GdalDataAnalysis.getResoulution(dataset,false));        //实际瓦片的长度，裁剪数据的分辨率和原始图像一致        int image_Xbuf = block_xsize;        int image_Ybuf =  block_ysize;        //写入瓦片的偏移位置        int imageOffSetX = 0;        int imageOffSetY = 0;        //循环读写        Driver driver = gdal.GetDriverByName("GTiff");        int dataType = dataset.GetRasterBand(1).getDataType();        int bands = dataset.getRasterCount();        byte[] pixelArray = new byte[image_Xbuf*image_Ybuf];        Dataset dsDataset = driver.Create("D:\\data\\xian18_1.tif",image_Xbuf,image_Ybuf,bands, dataType);        for (int i = 1;i