歡迎來到Linux教程網
Linux教程網
Linux教程網
Linux教程網
您现在的位置: Linux教程網 >> UnixLinux >  >> Linux編程 >> Linux編程

在Python中使用GDAL為衛片制作縮略圖

衛片和航片文件一般都很大,

有時候需要制作一個縮略圖,放在導航窗體中做底圖

這個文件不用太清晰,只要反映整體面貌即可

在ArcMap,PCI和Erdas中,可以直接截屏,或者使用重投影等功能很方便地完成,

最近在用python和GDAL,於是乎,查了查資料,並不復雜,

首先要確定python安裝了numpy庫,代碼

  1. import gdal  
  2. import numpy  
  3. from gdalconst import *  
  4. dataset = gdal.Open("F://數據//Raster//earth.img")  
  5. width = dataset.RasterXSize  
  6. height = dataset.RasterYSize  
  7. bw = 0.1  
  8. bh = 0.1  
  9. x = int(width*bw)  
  10. y = int(height*bh)  
  11. datas = []  
  12. for i in range(3):  
  13.     band = dataset.GetRasterBand(i+1)  
  14.     data = band.ReadAsArray(0,0,width,height,x,y)  
  15.     datas.append(numpy.reshape(data,(1,-1)))  
  16. datas = numpy.concatenate(datas)  
  17. driver = gdal.GetDriverByName("GTiff")  
  18. tods = driver.Create("F://數據//Raster//tpix1.tif",x,y,3,options=["INTERLEAVE=PIXEL"])  
  19. tods.WriteRaster(0,0,x,y,datas.tostring(),x,y,band_list=[1,2,3])  

以上的這個方法,用到了numpy的矩陣運算,不過使用band的ReadAsArray

方法取得的數據,采用的是最鄰近采樣,雖然縮略圖的要求不高,但是這個

采樣方式總讓人不爽,查了資料,使用gdal.ReprojectImage

可以設置采樣方式,代碼如下:

  1. # -*- coding: cp936 -*-   
  2. import gdal  
  3. import numpy  
  4. from gdalconst import *  
  5. dataset = gdal.Open("F://數據//Raster//earth.img")  
  6. width = dataset.RasterXSize  
  7. height = dataset.RasterYSize  
  8. bw = 0.1  
  9. bh = 0.1  
  10. x = int(width*bw)  
  11. y = int(height*bh)  
  12. driver = gdal.GetDriverByName("GTiff")  
  13. tods = driver.Create("F://數據//Raster//tpix1.tif",x,y,3,options=["INTERLEAVE=PIXEL"])  
  14. tods.SetGeoTransform(dataset.GetGeoTransform())  
  15. tods.SetProjection(dataset.GetProjection())  
  16. gdal.ReprojectImage(dataset,tods,dataset.GetProjection(),tods.GetProjection(),GRA_Cubic)  

但這個生成的不是縮略圖,而是從大圖中裁剪出來的一部分

Copyright © Linux教程網 All Rights Reserved