获取带坐标的影像中每个像元的坐标值

最新推荐文章于 2024-12-23 14:29:57 发布

陨星落云

最新推荐文章于 2024-12-23 14:29:57 发布

阅读量673

点赞数

分类专栏： gdal（python，C++）

本文链接：https://blog.csdn.net/qq_28368377/article/details/112308036

版权

gdal（python，C++）专栏收录该内容

38 篇文章

订阅专栏

# -*- coding: utf-8 -*-
"""
Created on Thu Jan  7 13:42:36 2021

@author: xiaohuihui
"""

from osgeo import gdal
import numpy as np
from numba import jit

@jit
def pixel2coord(nXSize, nYSize, c, a, b, f, d, e):
    '''
    返回：行索引，列索引，yp，xp # [nysize,nxsize,yp,xp]
    '''
    coord = np.zeros((nYSize,nXSize,4))
    for col in range(nYSize):
        for row in range(nXSize):
            xp = a * col + b * row + a * 0.5 + b * 0.5 + c
            yp = d * col + e * row + d * 0.5 + e * 0.5 + f
            coord[col,row,0] = col
            coord[col,row,1] = row
            coord[col,row,2] = yp
            coord[col,row,3] = xp         
            
    return coord


if __name__ == "__main__":
    
    #读取数据
    ds = gdal.Open(r'E:\beijing_landsat\2010-2020\LC81230322015250LGN01\LC08_L1TP_123032_20150907_20170404_01_T1_B1.TIF')

    nXSize = ds.RasterXSize #列数
    nYSize = ds.RasterYSize #行数
    
    # unravel GDAL affine transform parameters
    c, a, b, f, d, e = ds.GetGeoTransform()
    
    coord = pixel2coord(nXSize, nYSize, c, a, b, f, d, e)
    np.save("coord.npy",coord)
    
    print(coord[:5,:5,:])
    print(coord[:5,:5,0])
    print(coord[:5,:5,1])
    print(coord[:5,:5,2])
    print(coord[:5,:5,3])

结果：

[[[0.00000e+00 0.00000e+00 4.58250e+06 3.60600e+05]
  [0.00000e+00 1.00000e+00 4.58247e+06 3.60600e+05]
  [0.00000e+00 2.00000e+00 4.58244e+06 3.60600e+05]
  [0.00000e+00 3.00000e+00 4.58241e+06 3.60600e+05]
  [0.00000e+00 4.00000e+00 4.58238e+06 3.60600e+05]]

 [[1.00000e+00 0.00000e+00 4.58250e+06 3.60630e+05]
  [1.00000e+00 1.00000e+00 4.58247e+06 3.60630e+05]
  [1.00000e+00 2.00000e+00 4.58244e+06 3.60630e+05]
  [1.00000e+00 3.00000e+00 4.58241e+06 3.60630e+05]
  [1.00000e+00 4.00000e+00 4.58238e+06 3.60630e+05]]

 [[2.00000e+00 0.00000e+00 4.58250e+06 3.60660e+05]
  [2.00000e+00 1.00000e+00 4.58247e+06 3.60660e+05]
  [2.00000e+00 2.00000e+00 4.58244e+06 3.60660e+05]
  [2.00000e+00 3.00000e+00 4.58241e+06 3.60660e+05]
  [2.00000e+00 4.00000e+00 4.58238e+06 3.60660e+05]]

 [[3.00000e+00 0.00000e+00 4.58250e+06 3.60690e+05]
  [3.00000e+00 1.00000e+00 4.58247e+06 3.60690e+05]
  [3.00000e+00 2.00000e+00 4.58244e+06 3.60690e+05]
  [3.00000e+00 3.00000e+00 4.58241e+06 3.60690e+05]
  [3.00000e+00 4.00000e+00 4.58238e+06 3.60690e+05]]

 [[4.00000e+00 0.00000e+00 4.58250e+06 3.60720e+05]
  [4.00000e+00 1.00000e+00 4.58247e+06 3.60720e+05]
  [4.00000e+00 2.00000e+00 4.58244e+06 3.60720e+05]
  [4.00000e+00 3.00000e+00 4.58241e+06 3.60720e+05]
  [4.00000e+00 4.00000e+00 4.58238e+06 3.60720e+05]]]
[[0. 0. 0. 0. 0.]
 [1. 1. 1. 1. 1.]
 [2. 2. 2. 2. 2.]
 [3. 3. 3. 3. 3.]
 [4. 4. 4. 4. 4.]]
[[0. 1. 2. 3. 4.]
 [0. 1. 2. 3. 4.]
 [0. 1. 2. 3. 4.]
 [0. 1. 2. 3. 4.]
 [0. 1. 2. 3. 4.]]
[[4582500. 4582470. 4582440. 4582410. 4582380.]
 [4582500. 4582470. 4582440. 4582410. 4582380.]
 [4582500. 4582470. 4582440. 4582410. 4582380.]
 [4582500. 4582470. 4582440. 4582410. 4582380.]
 [4582500. 4582470. 4582440. 4582410. 4582380.]]
[[360600. 360600. 360600. 360600. 360600.]
 [360630. 360630. 360630. 360630. 360630.]
 [360660. 360660. 360660. 360660. 360660.]
 [360690. 360690. 360690. 360690. 360690.]
 [360720. 360720. 360720. 360720. 360720.]]