2

背景

我有一个网格大小的 NETCDF4 文件0.125x0.125。纬度从90-90,经度从0360。因此,完整的表格大小为1441 x 2880(纬度 x 经度)。

我正在获取我的位置坐标(以度为单位的纬度)并试图找到我所在的位置cell......

为了计算cell我所处的位置,我这样做:

'''
This function takes an array and a value.
It finds the item in the array for which the value most closely matches, and returns the index.
'''
def GetNearestIndex(arrayIn, value):
    distance = math.fabs(arrayIn[0] - value);
    nearest = 0;
    index = 0
    for idx, val in enumerate(arrayIn):
        delta_distance = math.fabs(val - value)
        if delta_distance < distance:
            nearest = val
            index = idx
            distance = delta_distance
    return index

#Lats and Longs arrays from the NETCDF4 dataset
lats = dataset.variables['latitude'][:]
longs = dataset.variables['longitude'][:]

#GetNearestIndex finds the item in the array for which the value most closely matches, and returns the index.
nearestLatitudeIndex = Utilities.GetNearestIndex(lats,  myLat)
nearestLongitudeIndex = Utilities.GetNearestIndex(longs, myLon%360)

所以给定我的 NECTDF4 数据集,如果我的位置是[31.351621, -113.305864](lat, lon),我发现我与cell[31.375, 246.75] (lat, lon) 匹配。将计算出的 lat 和 lon 插入GetNearestIndex,然后我将获得我所在单元格的“地址”(x,y)。

现在我知道cell我离哪个最近,我从 NETCDF 文件中获取值,然后我可以说“你所在位置的温度是 X”。

问题是,我不知道我这样做是否正确,所以我的问题是:

问题

  1. 如何正确确定我所在的单元格并获取 x 和 y 索引?

  2. 如何验证我的计算是否正确?

  3. myLon%360 是从 myLon 转换为从 0 到 360 度的网格的正确方法吗?网格单元大小无关紧要吗?

4

1 回答 1

0

我没有时间检查/测试您的方法,但我会使用 numpy,这是 NetCDF4 所基于的库。geo_idx()下面,是我用于 lat/lon 度的常规网格,lons 介于 -180 和 180 之间。这种方法避免了循环遍历所有 lat/lon 数组。

  import numpy as np
  import netCDF4

  def geo_idx(dd, dd_array):
    """
     - dd - the decimal degree (latitude or longitude)
     - dd_array - the list of decimal degrees to search.
     search for nearest decimal degree in an array of decimal degrees and return the index.
     np.argmin returns the indices of minium value along an axis.
     so subtract dd from all values in dd_array, take absolute value and find index of minium.
   """
    geo_idx = (np.abs(dd_array - dd)).argmin()
    return geo_idx

使用和测试

  # to test
  in_lat = 31.351621
  in_lon = -113.305864

  nci = netCDF4.Dataset(infile)
  lats = nci.variables['lat'][:]
  lons = nci.variables['lon'][:]
  # since lons are 0 thru 360, convert to -180 thru 180
  converted_lons = lons - ( lons.astype(np.int32) / 180) * 360

  lat_idx = geo_idx(in_lat, lats)
  lon_idx = geo_idx(in_lon, converted_lons)
  print lats[lat_idx]
  print converted_lons[lon_idx] 
于 2017-05-05T10:59:04.247 回答