全球高分辨率(3小时,10公里)地表太阳辐射数据集(1983-2018)数据单像元数值获取
全球高分辨率(3小时,10公里)地表太阳辐射数据集(1983-2018)数据单像元数值获取

数据格式是nc,本文采用Python进行数据提取(仅单个像元单年份或多年份的提取,多个像元暂时不知道如何提取,也由于暂时用不到,因此没有深入研究,欢迎各位老师指教)。
import os
import numpy as nup
from netCDF4 import Dataset
directory_name = "C:\\Users\\mhb\\Desktop\\xinjian"
#上面一行是文件目录所在位置
ff = os.listdir(directory_name)
for item in ff:
nc_file = directory_name + "\\" +item
nc_obj = Dataset(nc_file)
abc = nc_obj['daily_global_radiation'][:]
print(abc[1164][1193])
#上面一行是要提取的数据所在该nc文件的“列+行”号,此处存疑,仅个人推测是列行号而非行列号,待各位测试后反馈。
需要更改的部分仅目录位置和需要提取的像元所在数据列行号。
先说如何提取,再说推算的过程(经过与@利未的博客_CSDN博客-领域博主的交流,确定了以下推算方法,在此特别感谢!)。
+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
下述abs=取绝对值。
令a=abs(经度*10)&忽略小数位数取整,
东经的情况下:x=a,例如E 119.359°,则x=1993;
西经的情况下:x=3600-a-1,例如W 170.312°,则x=1896;
令b=纬度,
北纬的情况下:y=abs(b-89.95)*10,例如N 26.297°,则y=636;
南纬的情况下:y=abs(b+89.95)*10,例如S 26.297°,则y=1162。
在pycharm中将最后一行改为[y][x]的形式,即:
print(abc[1164][1193])
运行,输出最后结果,如下:

+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
像元所在数据列行号是根据软件Panoply确定的,具体的下载、安装及配置参考以下文章。

双击后配置经纬度,此处更改为X轴对应经度,Y轴为纬度。

可以得到下图及行列号对应的数组,我们要做的就是根据经纬度和来确定单个像元的位置。

上面可以看到,纬度共180°被分为了1800行,经度共360°被分为了3600列。官网上下载的坐标信息显示左上角纬度89.95°,经度0.05°;右下角纬度-88.95°,经度359.95°。

因此根据所在区域的经纬度可以进行以下换算
令a=abs(经度*10)&忽略小数位数取整,
东经的情况下:x=a,例如E 119.359°,则x=1993;
西经的情况下:x=3600-a-1,例如W 170.312°,则x=1896;
令b=纬度,
北纬的情况下:y=abs(b-89.95)*10,例如N 26.297°,则y=636;
南纬的情况下:y=abs(b+89.95)*10,例如S 26.297°,则y=1162。
完毕
更多推荐
https://data.tpdc.ac.cn/zh-hans/data/be562de3-6367-402f-956d-59f7c21ad294/?q=%E5%85%A8%E7%90%83%E9%AB%98%E5%88%86%E8%BE%A8%E7%8E%87%EF%BC%883%E5%B0%8F%E6%97%B6%EF%BC%8C10%E5%85%AC%E9%87%8C%EF%BC%89%E5%9C%B0%E8%A1%A8%E5%A4%AA%E9%98%B3%E8%BE%90%E5%B0%84%E6%95%B0%E6%8D%AE%E9%9B%86%EF%BC%881983-2018%EF%BC%89
https://blog.csdn.net/ch206265/article/details/103516883

所有评论(0)