利用griddata进行插值

利用griddata进行插值因为最近在做算法优化,所以对数据统一性有一定要求,在最近的研究中主要用一个简单的最近邻插值对数据集进行降尺度处理。主要运用到的函数时scipy里面的griddata第一步:导入相关库importxarrayasxrfromscipy.interpolateimportgriddata#插值函数importnumpyasnp第二步:给出插值到的经纬度信息(目标经纬度)mask_tmp=xr.open_dataset(‘G:/China/temperatu

大家好,又见面了,我是你们的朋友全栈君。

因为最近在做算法优化,所以对数据统一性有一定要求,在最近的研究中主要用一个简单的最近邻插值对数据集进行降尺度处理。

主要运用到的函数时scipy里面的
griddata

griddata函数讲解

`scipy.interpolate.griddata(points, values, xi, method='linear', fill_value=nan, rescale=False)`
  • points:2-D ndarray of floats with shape (n, D), 或 length D tuple of 1-D ndarrays with shape (n,).
  • values:ndarray of float 或 complex, shape (n,)
  • xi:2-D ndarray of floats with shape (m, D), 或 length D tuple of ndarrays broadcastable to the same shape.
  • method:{‘linear’, ‘nearest’, ‘cubic’}, 可选参数

nearest:返回最接近插值点的数据点的值
linear:线性插值
cubic:三次样条插值

第一步:导入相关库

import xarray as xr    
from scipy.interpolate import griddata    # 插值函数
import numpy as np

第二步:给出插值到的经纬度信息(目标经纬度)

mask_tmp = xr.open_dataset('G:/China/temperature_max/nc/2000/tmx_2000_1.nc')
# 待插值的标准经纬度
mask_tmp_lon = mask_tmp.lon.values  # (7680,) 一维
mask_tmp_lat = mask_tmp.lat.values  # (4717,) 一维
mask_LON,mask_LAT = np.meshgrid(mask_tmp_lon,mask_tmp_lat)  # (4717, 7680) 生成经纬度网格

第三步:待插值数据

rad = xr.open_dataset('D:/Users/62692/Desktop/rad.nc')
# 辐射数据经纬度
rad_lon = rad.lon.values    # 辐射数据经度 (641,)
rad_lat = rad.lat.values    # 辐射数据纬度 (394,)
rad_LON, rad_LAT = np.meshgrid(rad_lon,rad_lat)   # (394, 641)
rad_LON = rad_LON.ravel().reshape(-1,1)  # 展平 (252554, 1)
rad_LAT = rad_LAT.ravel().reshape(-1,1)  # 展平 (252554, 1)
rad_values = rad['rad'].values     # 需要插值的辐射数据 (394, 641)
points = np.concatenate([rad_LON,rad_LAT],axis = 1)   # (252554, 2)

第四步:插值

# 插值
data = griddata(points, rad_values.ravel(),(mask_LON,mask_LAT),method='nearest')  # 用最近邻插值即可
# rad_values.ravel() (252554,)
# 这里用最邻近主要考虑到辐射数据的连续性变化,对于线性插值或者三次插值并没有多大影响

汇总成函数

''' Created on 1 23, 2022 @author: GongHaixing 将一个文件夹里面所有的nc文件进行插值 '''
def interp2D(maskpath,mask_lon='lon',mask_lat='lat',inputpath='', outputpath='',data_lon='lon',data_lat='lat',variable='',interp_method='nearest',save=True):
    """输入插值目标的相关信息以及需要插值的数据 :maskpath: 需要插值到对应数据的数据路径 :mask_lon: 标准数据的经度名称,比如:x,lon :mask_lat: 标准数据的纬度名称,比如:y,lat :inputpath: 需要做插值处理的nc文件所在的目录 :outputpath: 插值完nc文件保存的路径,注意要是'/' :data_lon: 需要做插值数据经度名称,比如:'x','lon' :data_lat: 需要做插值数据经度名称,比如:'y','lat' :variable:需要做插值数据变量的名称,比如:'tmp','ndvi' :interp_method: griddata的插值方法,比如:'nearest','linear','cubic' :save:是否对文件进行存储 """
    #导入相关库
    import xarray as xr
    import os
    from scipy.interpolate import griddata    # 插值函数
    import numpy as np

    ### 目标插值
    mask_data = xr.open_dataset(maskpath)
    mask_data_lon = mask_data[mask_lon].values
    mask_data_lat = mask_data[mask_lat].values
    mask_LON,mask_LAT = np.meshgrid(mask_data_lon,mask_data_lat)
    mask_LON1 = np.array(mask_LON)
    mask_LAT1 = np.array(mask_LAT)

    ### 插值对象
    os.chdir(inputpath)                            # 给出nc文件所在的目录(路径)
    files = os.listdir()
    savepath = outputpath
    for file in files:
        inputfile = xr.open_dataset(file)
        inputfile_lon = inputfile[data_lon].values                                 # 数据的经纬度
        inputfile_lat = inputfile[data_lat].values
        inputfile_LON, inputfile_LAT = np.meshgrid(inputfile_lon,inputfile_lat)
        inputfile_LON = inputfile_LON.ravel().reshape(-1,1)
        inputfile_LAT = inputfile_LAT.ravel().reshape(-1,1)
        inputfile_values = np.array(inputfile[variable].values,dtype=np.float32)     # 需要插值的数据
        points = np.concatenate([inputfile_LON,inputfile_LAT],axis = 1)
        # 插值
        print('开始对'+file+'进行插值')
        inputfile_interp = griddata(points, inputfile_values.ravel(),(mask_LON1,mask_LAT1),method=interp_method).astype(np.float32)  # 用最近邻插值即可,不然数据会nan
        outfile= xr.Dataset(
                      data_vars = { 
   variable:(('lat','lon'),inputfile_interp)},
                      coords = { 
   
                                'lon':('lon',np.array(mask_data_lon)),
                                'lat':('lat',np.array(mask_data_lat))}
                                )
        if save == True:
            outfile.to_netcdf(outputpath+'/'+file)
            print(file+'已经插值成功,且已经保存到'+outputpath+'路径下')
        else:
            print(file+'已经插值成功,但是我没有保存文件')
from interp2D import *
import xarray as xr
import os
import pandas as pd
from scipy.interpolate import griddata    # 插值函数


maskpath = 'H:/China/temperature_max/nc/2000/tmx_2000_1.nc'
mask_lon='lon'
mask_lat='lat'
inputpath='H:/China/LAI/nc'
outputpath='H:/China/LAI/nc_1km'
data_lon='lon'
data_lat='lat'
variable='LAI'
interp_method='nearest'
save=True

interp2D(maskpath=maskpath,mask_lon='lon',mask_lat='lat',inputpath=inputpath, outputpath=outputpath,data_lon='lon',data_lat='lat',variable=variable,interp_method='nearest',save=True)

结果对比

插值前(10km)

在这里插入图片描述

插值后(1km)

在这里插入图片描述

版权声明:本文内容由互联网用户自发贡献,该文观点仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请联系我们举报,一经查实,本站将立刻删除。

发布者:全栈程序员-站长,转载请注明出处:https://javaforall.net/141272.html原文链接:https://javaforall.net

(0)
全栈程序员-站长的头像全栈程序员-站长


相关推荐

  • 如何将ER图转化为关系模型(超详细、含例题)「建议收藏」

    如何将ER图转化为关系模型(超详细、含例题)「建议收藏」本文转载自:http://blog.csdn.net/HaoDaWang/article/details/78098937?locationNum=4&fps=1超级感谢博主分享本篇博文中的例题只涉及到N:M和N:M:P的转换,大家看完后可以再去看下百度文库里面的这个例题,这个例题几个转换准则都用到了。https://wenku.baidu.com/view/33849a73a…

    2022年9月22日
    1
  • 2021年7月整理–简单方法 暴力激活成功教程WIFI密码

    2021年7月整理–简单方法 暴力激活成功教程WIFI密码2021年7月整理–简单方法暴力激活成功教程WIFI密码很多人都面临过短期租房、短期出差、住院而没有WIFI可用等境遇,有的是宽带太多办不起、有的是临时一阵子不值得折腾、有的是运营商不给扯线等等原因。然后就用手机下载了WIFI智能钥匙等APP,然后发现卵用么有,根本没有人共享自家WIFI密码给你用。以下是按步骤整理的软件和详细教程笔记本电脑+软件暴力激活成功教程出的密码我亲身用这个软件解开N多个密码此软件是家用路由器安全审计工具,切勿用作非法用途!!!….

    2022年10月13日
    3
  • ES集群搭建「建议收藏」

    ES集群搭建「建议收藏」1.设置内核参数:/etc/sysctl.confvm.max_map_count=6553602.设置文件描述符限制:/etc/security/limits.conf*softnofile65536*hardnofile131072*softnproc2048*hardnproc4096&&1和2的配置,可以等到启动时候报错,再修改(可以清晰地知道1、2在处理什么问题)3.es不允许root用户启动,修改用户和分组chownuser1.

    2022年10月13日
    2
  • lldp 命令「建议收藏」

    lldp 命令「建议收藏」LLDP协议介绍:简单说来,LLDP是一种邻近发现协议。它为以太网网络设备,如交换机、路由器和无线局域网接入点定义了一种标准的方法,使其可以向网络中其他节点公告自身的存在,并保存各个邻近设备的发现信息。例如设备配置和设备识别等详细信息都可以用该协议进行公告。具体来说,LLDP定义了一个通用公告信息集、一个传输公告的协议和一种用来存储所收到的公告信息的方法。要公告自身信息的设备…

    2022年5月28日
    69
  • 数据库严格模式_mysql关系型数据库

    数据库严格模式_mysql关系型数据库设置严格模式为当前会话设置严格模式:执行SETsql_mode=’STRICT_TRANS_TABLES’或者SETsql_mode=’STRICT_ALL_TABLES’。全局设置严格模式:执行SETglobalsql_mode=’STRICT_TRANS_TABLES’或者SETglobalsql_mode=’STRICT_ALL_TABLES’。或者在MyS…

    2025年7月10日
    7
  • 深入理解 HashMap

    深入理解 HashMap什么是HashMap?​ HashMap是基于哈希表的Map接口是实现的。此实现提供所有可选操作,并允许使用null做为值(key)和键(value)。HashMap不保证映射的顺序,特别是它不保证该顺序恒久不变。此实现假定哈希函数将元素适当的分布在各个桶之间,可作为基本操作(get和put)提供稳定的性能。在jdk1.7中的HashMap是基于数组+链表实现的,在jdk1….

    2025年10月19日
    5

发表回复

您的邮箱地址不会被公开。 必填项已用 * 标注

关注全栈程序员社区公众号