NCL處理nc就是會更方便,但是避免不了一些用不了NCL的情況
寫了一個簡單用matlab提取nc文件站點數據的腳本,通過判斷經緯度與格點經緯度距離,選取最近格點值。
后面也附上了提取wrfout中站點數據的NCL腳本,如果是再分析資料的話,wrf系列函數不知道是否可以,可能還是需要通過距離判斷。
matlab code:
clc;???%清屏
clear all;?%清空
%%
%nc批量讀取的數據
datadir1='C:\Users\sunxy\Desktop\';?????????????????%指定批量數據所在的文件夾
filelist1=dir([datadir1,'200901.nc']);???????%指定批量數據的類型
%%
%讀取數據
lat_station=32.5;lon_station=90.8;
?filename1=[datadir1,filelist1.name];
ncid1=netcdf.open(filename1,'NC_NOWRITE');??????????%打開nc文件???
lat=ncread(filename1,'latitude');?????????????%讀入變量lat
lon=ncread(filename1,'longitude');?????????????%讀入變量lon
?[lat_2D lon_2D]=meshgrid(lat,lon);
distance=sqrt((lat_2D-lat_station).^2+(lon_2D-lon_station).^2);???%算出站點與每個格點距離
?[ilon ilat]=find(distance==min(min(distance)));????????%尋找距離最近格點
?lat_back=lat(ilat);??????????????????%輸出經緯度檢查是否為希望讀取站點
?lon_back=lon(ilon,:);
?blh= ncread(filename1,'blh');
?blh_station=blh(ilon,ilat,:);
?netcdf.close(ncid1);??????????????????%關閉nc文件
%?
% % ncdisp('C:\Users\sunxy\Desktop\200901.nc','/','full')
NCL code:
??row = numAsciiRow("/datadir1/station.txt")
??col = numAsciiCol("/datadir1/station.txt")
??f_h = asciiread("/datadir1/station.txt",(/row,col/),"float")
??station=f_h(:,0)
??lon=f_h(:,1)
??lat=f_h(:,2)
?station_num=dimsizes(station)
??;read wrfout?u v data
?DATADir=(/"/datadir1/zhuyan/wrfout/2020/expt1/18/"/)
??FILE= systemfunc (" ls -1 " + DATADir + "wrfout_d01_2019*")
?filenum=dimsizes(FILE)
??a = addfiles(FILE + ".nc", "r")
??u?= wrf_user_getvar(a,"U10",-1)
?output_u=new((/station_num,filenum+5/), float)
?output_u(:,0)=station
??output_u(:,1)=lon
??output_u(:,2)=lat
??do istation=0,station_num-1
???loc?= wrf_user_ll_to_ij(a,lon(istation), lat(istation), True)-1
???output_u(istation,3)=loc(0)
???output_u(istation,4)=loc(1)
???output_u(istation,5::)=u(:,loc(1),loc(0))
??end do
??opt??=True
??M=filenum+5
??fWidth = 7??; specify the format width
??fDec= 2?; specify the number to the right of decimal point
??fmtf= M + "f" + fWidth + "." + fDec???; same as "7f7.2"
?opt@fout?="output_u_202001.txt"
??write_matrix (output_u, fmtf, opt)