针对模式结果经常需要提取站点或者指定区域的数据，而wrfout作为曲面网格，不能像再分析数据那样进行规则的提取
slp(:,{50:70},{110:120}) ;直线网格的变量在有经纬度信息的情况下可以直接读取需要的范围
wrfout需要用到函数wrf_user_ll_to_ij，可参考官网wrf_user_ll_to_ij (ucar.edu)
有两种使用方法，一种是输入一组经纬度，输出的是离这个点最近的格点的index
loc = wrf_user_ll_to_ij(a,lon_station, lat_station,True)
output_pm25 = pm25(:,0,loc(1),loc(0)) ;输出最近格点上的PM2.5浓度
注意loc的0维是longitude（x axis）的index，loc的1维是latitude（y axis）的index，根据自己数据的维度（例如我的就是time*level*latitude*longitude）进行使用，不要弄反了
另一种是输入两组经纬度数据，输出的index可以框出一个区域来
loc = wrf_user_ll_to_ij(a,(/lon_min,lon_max/), (/lat_min,lat_max/),True)output_pm25 = wgt_areaave(pm25(:,0,loc(1,0):loc(1,1),loc(0,0):loc(0,1)),1.0,1.0,0)
同上，loc的0维是longitude的最小最大index，1维是latitude的最小最大index
最近还遇到一个特殊情况，数据里没有经纬度信息，只有投影信息wrf_ll_to_ij函数可以通过给定投影信息找到给定一组经纬度的对应格点，用法见官网wrf_ll_to_ij (ucar.edu)。但是在使用这个函数时，返回的格点数并不是整数，提取出来的变量值也很奇怪，不知道是不是有什么隐藏的用法，我还没有仔细研究，毕竟是个少数情况。
最后是有大佬根据投影信息给弄出了经纬度信息，我直接选取了最原始的办法，算站点到格点的最小距离选临近格点的值，再根据local_min函数找出最小值的坐标
dis = (lon_cmaq-station_lon(istation))^2 + (lat_cmaq-station_lat(istation))^2 ;计算距离
nmin = local_min(dis, True, 0.0) ;找出最小值及其坐标
output = pm_cmaq(:,nmin@yi,nmin@xi)
