Python叠加两幅栅格图像的实现方法

yipeiwu_com5年前Python基础

目的

现有两幅栅格图像,一个是某地区道路栅格图,一个是某地区土地利用类型图,需要将道路叠加到土地利用类型图中,即叠加后,重合的像元值以道路图为准,其余的像元值仍是土地利用类型图原有的像元值。

图1 道路信息图

图2 土地利用类型图

图3 结果图

具体实现

from gdalconst import *
from osgeo import gdal
import osr
import sys
import copy

#叠加两个栅格图像(一个道路栅格图,一个土地利用类型图),两幅图像重叠的像元值都是第一个图像的值,
#未重叠的像元值还是土地利用类型图上的值,最终结果便是土地利用类型图上面多了道路信息。

roadFile = 'E:\\Exercise\\test\\grasstest\\road_rastercalc.tif'
landuseFile = 'E:\\Exercise\\test\\grasstest\\landuse.tif'
roadDs = gdal.Open(roadFile, GA_ReadOnly)
landuseDs = gdal.Open(landuseFile, GA_ReadOnly)
if roadDs is None:
  print 'Can not open ', roadFile
  sys.exit(1)

geotransform = roadDs.GetGeoTransform()
projection=roadDs.GetProjection()
cols = roadDs.RasterXSize
rows = roadDs.RasterYSize
roadBand = roadDs.GetRasterBand(1)
roadData = roadBand.ReadAsArray(0,0,cols,rows)
roadNoData = roadBand.GetNoDataValue()

landuseBand = landuseDs.GetRasterBand(1)
landuseData = landuseBand.ReadAsArray(0,0,cols,rows)
landuseNoData = landuseBand.GetNoDataValue()


result = landuseData

for i in range(0,rows):
  for j in range(0,cols):
    if(abs(roadData[i,j] - 20) < 0.0001):
      result[i,j] = 20
    if((abs(landuseData[i,j] - landuseNoData)>0.0001) and (abs(roadData[i,j] - roadNoData) < 0.0001)):
      result[i,j] = landuseData[i,j]
    if((abs(landuseData[i,j] - landuseNoData)<0.0001) and (abs(roadData[i,j] - roadNoData) < 0.0001)):
      result[i,j] = landuseNoData
#write result to disk
resultPath = 'E:\\Exercise\\test\\grasstest\\result_landuse.tif'

format = "GTiff"  
driver = gdal.GetDriverByName(format)
ds = driver.Create(resultPath, cols, rows, 1, GDT_Float32)
ds.SetGeoTransform(geotransform)
ds.SetProjection(projection)
ds.GetRasterBand(1).SetNoDataValue(landuseNoData)
ds.GetRasterBand(1).WriteArray(result)  
ds = None

print 'ok---------'

以上这篇Python叠加两幅栅格图像的实现方法就是小编分享给大家的全部内容了,希望能给大家一个参考,也希望大家多多支持【听图阁-专注于Python设计】。

相关文章

python实现关键词提取的示例讲解

新人小菜鸟又来写博客啦!!!没人表示不开心~~(>_<)~~ 今天我来弄一个简单的关键词提取的代码 文章内容关键词的提取分为三大步: (1) 分词 (2) 去停用词 (...

使用 Python 清理收藏夹里已失效的网站

使用 Python 清理收藏夹里已失效的网站

失效的书签们 我们日常浏览网站的时候,时不时会遇到些新奇的东西( 你懂的.jpg ),于是我们就默默的点了个收藏或者加书签。然而当我们面对成百上千的书签和收藏夹的时候,总会头疼不已……...

浅谈python对象数据的读写权限

面向对象的编程语言在写大型程序的的时候,往往比面向过程的语言用起来更方便,安全。其中原因之一在于:类机制。 类,对众多的数据进行分类,封装,让一个数据对象成为一个完整的个体,贴近现实生活...

Flask框架学习笔记之路由和反向路由详解【图文与实例】

Flask框架学习笔记之路由和反向路由详解【图文与实例】

本文实例讲述了Flask框架学习笔记之路由和反向路由。分享给大家供大家参考,具体如下: #-*- coding:utf-8 -*- from flask import Flask,...

Python常见字典内建函数用法示例

本文实例讲述了Python常见字典内建函数用法。分享给大家供大家参考,具体如下: 1、len(mapping)      &n...