geotrellis使用（十四）导出定制的GeoTiff
Geotrellis系列文章链接地址http://www.cnblogs.com/shoufengwei/p/5619419.html
目录
一、前言
最近一段时间比较忙，没能继续推进Geotrellis项目开发，周末和这两天抽空又实现了一个功能——导出自定义的Tiff文件。又恰巧碰上今天这么重要的日子，当然要写点东西来纪念一下，所以就有了这篇文章，闲话莫说，进入正题。
二、需求说明
很多时候我们需要从一块（或者很多块）大的Tiff中根据需要截取一部分数据，并且需要采用某种采样方式转成特定的投影，并转成需要的数据类型。当然有人会说这个很容易，用GDAL的gdaltransform等可以很容易的实现此功能，GDAL是很强大，但是前提是你的数据不能太大并且只能处理单块栅格数据。正因为有这些问题，所以我实现了使用Geotrellis来实现该功能，下面我就为大家分析实现方法。
之前讲了很多数据处理方式，其中。
三、实现方案
1.前台界面
前台就是一个简单的地图控件，外加几个选择器。地图控件主要为了浏览区域以及手工选取想要导出的Tiff的范围，选择器主要选择目标投影方式、数据类型以及采样方式。根据用户的选择将请求采用ajax的方式发送到后台进行处理。这块不是本文的重点，不在这里具体介绍。
2.数据导入
要想处理大数据或者处理多块栅格数据就不能直接处理栅格数据，可以先将栅格数据导入到Accumulo中，当然导入之后是一块块的瓦片，这部分在geotrellis使用（三）geotrellis数据处理过程分析一文中已经进行了详细介绍，这里还是要说明的是参数layoutScheme一定要选择floating，这样在Accumulo中保存的就是原始只是切割而未经过其他处理的数据。
3.读取数据
Accumulo中已经存储了需要的数据，并且后台接收到了前台用户选择的区域范围以及投影方式、数据类型、采样方式，这样我们就可以开始实现读取需要的数据，简单的说就是从Accumulo中取出与用户选定区域有交集的数据（瓦片），实现代码如下：
val raster = reader.read[SpatialKey, Tile, TileLayerMetadata[SpatialKey]](layerId)
val masked = raster.mask(polygon)
是不是很容易，就这两句话，其中reader是AccumuloLayerReader实例，layerId表示读取的层，polygon是用户选取的范围。首先从Accumulo中读出该层数据，然后与polygon做一个mask，得到的结果就是用户想要导出的数据。
4.数据转换
4.1 得到结果瓦片数据
首先将上述得到的瓦片集做一个拼接，这样就会得到一个大的瓦片，代码也很简单，如下：
val stitch = masked.stitch
val tile = stitch.tile
4.2 得到最终结果
之后要做的就是根据采样类型、投影方式以及数据类型将上述tile进行转换，代码如下：
val rep = tile.reproject(extent, srcCRS, dstCRS, method).tile
val res = rep.convert(cellType)
其中extent是瓦片的范围，可以用上述拼接对象stitch的stitch.extent得到，srcCRS是原始数据的投影类型，可以通过上述raster对象的raster.metadata.crs得到，dstCRS表示目标投影方式，即用户想要的投影方式，同样method表示采样方法，cellType表示数据类型，这三个的获取方式我会重点介绍。
4.3 获取dstCRS、resampleMethod、cellType
首先可以肯定的是前端传来的这三个参数都是字符串型的，这就要求后台将字符串转成相应的类型。
投影方式我这里投了个懒，传递的是EPSG CODE,EPSG是投影方式的一种数字编码，具体请见：http://www.epsg.org/。总之每一个EPSG编码对应了一种投影方式，像常见的经纬度投影的编码是4326，WebMercator的编码是3857。有了这个编码之后就可以很容易的得到投影方式，代码如下：
try{
CRS.fromEpsgCode(epsg)
}
catch {
case _ => CRS.fromEpsgCode(3857)
}
采样方式和数据类型就稍显麻烦，我这里采用了反射的方式获取对应实例，我将反射的代码进行了封装，代码如下L:
def getClassName(name: String, hasEndWith$: Boolean) = {
var className = name
if (hasEndWith$) {
if (!className.endsWith("$"))
className += "$"
}
else {
if (className.endsWith("$"))
className = className.substring(0, className.length - 1)
}
className
}
def getClassFor[T: ClassTag](name: String, hasEndWith$: Boolean) = {
val className = getClassName(name, hasEndWith$)
Try[Class[_ <: T]]({
val c = Class.forName(className, false, getClass.getClassLoader).asInstanceOf[Class[_ <: T]]
val t = implicitly[ClassTag[T]].runtimeClass
if (t.isAssignableFrom(c)) c else throw new ClassCastException(t + " is not assignable from " + c)
})
}
def getInstanceByReflact[T: ClassTag](name: String) = {
val classTry = getClassFor[T](name, true) recoverWith { case _ => getClassFor[T](name, false) }
classTry flatMap{
c =>
Try{
val module = c.getDeclaredField("MODULE$") //通过获得变量"MODULE$"来初始化
module.setAccessible(true)
val t = implicitly[ClassTag[T]].runtimeClass
module.get(null) match {
case null => throw new NullPointerException
case x if !t.isInstance(x) => throw new ClassCastException(name + " is not a subtype of " + t)
case x: T => x
}
} recover{ case i: InvocationTargetException if i.getTargetException ne null => throw i.getTargetException }
}
这段代码看上去很复杂，其实是做了很多判断加强程序稳定性用的的，主要功能就是一个类的全名（包含包名）创建其实例。
有了上述代码之后就可以将前台传来的字符串直接转换为相应的实例，获取采样方式的代码如下：
val resample = s"geotrellis.raster.resample.${resampleStr}"
val method = getInstanceByReflact[ResampleMethod](resample)
if(method.isSuccess) {
println(method.get.toString)
method.get
}
else {
Bilinear
}
前台只需传入采样方式名称即可，如Bilinear，这样就会自动获取其实例。
获取数据类型同理，代码如下：
val className = getCellTypeClassName(dataTypeStr)
val cellType = getInstanceByReflact[CellType](className)
if(cellType.isSuccess) {
cellType.get
}
else {
ByteConstantNoDataCellType
}
这样就可实现获取用户想要的投影方式、采样方法以及数据类型。
四、总结
以上就是使用Geotrellis实现导出定制的GeoTiff的方法，由于时间紧，可能还有很多没有注意的细节，会在后续中进一步研究，并更新该文或者另设新篇。在这里我也呼吁大家都来撰写技术博客，一来提高自己，进行技术总结；二来真的是可以帮助别人。最近在工作中碰到很多莫名其妙的BUG，但是基本都在一番搜索之后在别人的博客中找到了解决方案，所以在此也感谢那些愿意写博客的人，希望大家都能加入进来。
