如果有人想要每个经纬度网格单元中平方米数的 netcdf 。
这可能不是最干净的解决方案，但会使用 xarray 在每个网格中创建 m2 的 netcdf (earth_m2.nc)。
The gridsize()函数改编自另一个堆栈溢出问题。
然后，我们可以创建一个虚拟数组，并使用每个位置的经度距离创建一个 m2s 的地球场。
"""
This will create a global grid of the approximate size of each grid square.
"""
import numpy as np
import xarray as xr
import matplotlib.pyplot as plt
def gridsize(lat1):
#https://en.wikipedia.org/wiki/Haversine_formula
#https://stackoverflow.com/questions/639695/how-to-convert-latitude-or-longitude-to-meters/11172685#11172685
lon1=200
import math
lat2=lat1
lon2=lon1+1
R = 6378.137 # // Radius of earth in km
dLat = lat2 * np.pi / 180 - lat1 * np.pi / 180
dLon = lon2 * np.pi / 180 - lon1 * np.pi / 180
a = np.sin(dLat/2) * np.sin(dLat/2) + np.cos(lat1 * np.pi / 180) * np.cos(lat2 * np.pi / 180) * np.sin(dLon/2) * np.sin(dLon/2)
c = 2 * np.arctan2(np.sqrt(a), np.sqrt(1-a))
d = R * c
return d * 1000 #; // meters
boxlo,boxla=np.array(np.meshgrid(np.arange(-179.5,179.5,1),np.arange(-89.5,89.5,1)))
sizes=np.ones(boxlo.shape)
grid=gridsize(boxla)
grid_nc=xr.DataArray(grid,coords={'lat':boxla[:,1],'lon':boxlo[1,:]},dims=['lat','lon'])
lat_size=110567 #in m
grid_nc['m2']=grid_nc*lat_size
grid_nc=grid_nc['m2']
grid_nc.to_netcdf('earth_m2.nc')
plt.pcolormesh(boxlo[1,:],boxla[:,1],grid_nc)
plt.colorbar()
plt.show()
