一、应用场景
最近在求非平坦地区的太阳高度角,其中有一个参数就是纬度,因此就想到了把纬度作为值来创建二维数组的想法。
二、参考文献
《读取遥感影像的信息》 《使用python读取tiff文件中的经纬度,并将数据以excel表的形式输出(详细步骤)》
三、代码展示
import os
import numpy as np
from osgeo import gdal
import glob
import time
import math
import copy
DEM_Path = "D:\\DATA\\SEBAL\\0. InitialData\\DEM1984.tif"
dataset = gdal.Open(DEM_Path)
width = dataset.RasterXSize
height = dataset.RasterYSize
geo_information = dataset.GetGeoTransform()
proj = dataset.GetProjection()
filevalue= np.array(dataset.ReadAsArray(0, 0, width, height), dtype=float)
LongitudeValue=np.zeros((height,width))
LatitudeValue=np.zeros((height,width))
for x in range(height):
for y in range(width):
LongitudeValue[x][y]=geo_information[0] + y * geo_information[1] + x * geo_information[2]
LatitudeValue[x][y]=geo_information[3] + y * geo_information[4] + x * geo_information[5]
最后得到的数组和图像分别是这样的: 经度:
纬度:
四、反思与注意事项
1、遥感影像的长度和宽度以及二维数组的行列要对应正确,比如for循环的时候,如果写反了会出现长度不够的错误 2、 创建0矩阵的时候也是要注意参数中的行列对应,如果不行就交换一下 3、这个算法的核心是掌握地理信息的六个参数: geo_information(0):左上像素左上角的x坐标。 geo_information(1):w - e像素分辨率 / 像素宽度。 geo_information(2):行旋转(通常为零)。 geo_information(3):左上像素左上角的y坐标。 geo_information(4):列旋转(通常为零)。 geo_information(5):n - s像素分辨率 / 像素高度(北半球上图像为负值)
|