I'm working with Pygrib trying to get surface temperatures for particular lat/lon coordinates using the NBM grib data (available here if it helps).
I've been stuck trying to get an index value to use with representative data for a particular latitude and longitude. I was able to derive an index, but the problem is the latitude and longitude appear to have 2 coordinates each. I'll use Miami, FL (25.7617° N, 80.1918° W) as an example to illustrate this. Formatted to be minimum reproducible IF a grib file is provided.
def get_grib_data(self, gribfile, shortName):
grbs = pygrib.open(gribfile)
# Temp needs level specified
if shortName == '2t':
grib_param = grbs.select(shortName=shortName, level=2)
# Convention- use short name for less than 5 chars
# Else, use name
elif len(shortName) < 5:
grib_param = grbs.select(shortName=shortName)
else:
grib_param = grbs.select(name=shortName)
data_values = grib_param[0].values
# Need varying returns depending on parameter
grbs.close()
if shortName == '2t':
return data_values, grib_param
else:
return data_values
# This function is used to find the closest lat/lon value to the entered one
def closest(self, coordinate, value):
ab_array = np.abs(coordinate - value)
smallest_difference_index = np.amin(ab_array)
ind = np.unravel_index(np.argmin(ab_array, axis=None), ab_array.shape)
return ind
def get_local_value(data, j, in_lats, in_lons, lats, lons):
lat_ind = closest(lats, in_lats[j])
lon_ind = closest(lons, in_lons[j])
print(lat_ind[0])
print(lat_ind[1])
print(lon_ind[0])
print(lon_ind[1])
if len(lat_ind) > 1 or len(lon_ind) > 1:
lat_ind = lat_ind[0]
lon_ind = lon_ind[0]
dtype = data[lat_ind][lon_ind]
else:
dtype = data[lat_ind][lon_ind]
return dtype
if __name__ == '__main__':
tfile = # Path to grib file
temps, param = grib_data.get_grib_data(tfile, '2t')
lats, lons = param[0].latlons()
j = 0
in_lats = [25.7617, 0 , 0]
in_lons = [-80.198, 0, 0]
temp = grib_data.get_local_value(temps, j, in_lats, in_lons, lats, lons)
When I do the print listed, I get the following for indices:
lat_ind[0]: 182
lat_ind[1]: 1931
lon_ind[0]: 1226
lon_ind[1]: 1756
So if my lat/lon were 1 dimensional, I would just do temp = data[lat[0]][lon[0]], but in this case that would give non-representative data. How would I go about handling the fact that lat/lon are in 2 coordinates? I have verified that lats[lat_ind[0][lat_ind1] gives the input latitude and the same for longitude.
from How to Deal with Lat/Lon Arrays with Multiple Dimensions?
No comments:
Post a Comment