在 Matlab 中为每个纬度/经度提取 netCDF 时间序列答案

【问题标题】：Extracting netCDF time series for each lat/long in Matlab在 Matlab 中为每个纬度/经度提取 netCDF 时间序列
【发布时间】：2015-05-07 05:12:20
【问题描述】：

我目前正在处理来自气候模型的 netCDF 输出，并希望获取 netCDF 中每个纬度/经度组合的时间序列文本文件。例如，如果 netCDF 有 10 个纬度和 10 个经度，我将获得 100 个文本文件，每个文件都有一个列格式的时间序列。我对 Matlab/netCDF 语言相当熟悉，但我似乎无法理解这一点。命名文本文件并不重要；我将它们重命名为“latitude_longitude_PCP.txt”，其中PCP是经纬度位置的降水量。

任何帮助将不胜感激。谢谢。

--达伦

【问题讨论】：

很想知道为什么要将它们放在单独的文件中。您是否正在处理非常大的数据（与您机器上的内存相比），您的分析一次只能处理一个时间序列？
我正在提取单独的时间序列，以便将它们输入到水文模型中，该模型使用纬度/经度的天气时间序列。

标签： matlab netcdf

【解决方案1】：

有几种方法可以解决这个问题。

方法 1. 如果您能够将 netcdf 文件放在 THREDDS 数据服务器上，您可以使用 NetCDF 子集服务网格作为点来指定经度/纬度点并以 CSV 或 XML 格式取回数据。这是 Unidata 的 THREDDS 数据服务器的一个示例：http://thredds.ucar.edu/thredds/ncss/grid/grib/NCEP/GFS/Global_0p5deg/best/pointDataset.html

方法 2。如果您想使用 Matlab 提取特定经度/纬度位置的时间序列，您可以使用 NCTOOLBOX 中的“nj_tseries”函数，该函数位于：http://nctoolbox.github.io/nctoolbox/

方法 3。如果您真的想使用 Matlab 在 [time,lon,lat] 网格中的每个 i,j 位置编写一个 ASCII 时间序列，您可以执行以下操作（使用 NCTOOLBOX）：

url='http://thredds.ucar.edu/thredds/dodsC/grib/NCEP/GFS/Global_2p5deg/best';
nc = ncgeodataset(url);
nc.variables
var='Downward_Short-Wave_Radiation_Flux_surface_12_Hour_Average';
lon = nc.data('lon');
lat = nc.data('lat');
jd = nj_time(nc,var);
ncvar = nc.variable(var);
for j=1:length(lat)
    for i=1:length(lon)
        v=ncvar.data(:,j,i);
        outfile=sprintf('%6.2flon%6.2flat.csv',lon(i),lat(j))
        fid=fopen(outfile,'wt')
        data= [datevec(jd) v]
        fprintf(fid,'%2.2d %2.2d %2.2d %2.2d %2.2d %2.2d %7.2f\n',data')
        fclose(fid)
        disp([outfile ' created.'])
    end
end

如果你有足够的内存将所有数据读入matlab，你可以在双循环之外读取，这样会快很多。但是无论如何写 ASCII 很慢，所以可能没那么重要。

【讨论】：

【解决方案2】：

%% Create demo data 
data = reshape(1:20*30*40,[20 30 40]);
nccreate('t.nc','data','Dimensions',{'lat', 20, 'lon',30, 'time', inf});
ncwrite('t.nc', 'data',data);
ncdisp('t.nc');

%% Write timeseries to ASCII files
% Giving an idea of the size of your data can help people 
% recommend different approaches tailored to the data size. 
% For smaller data, it might be faster to read in the full 
% 3D data into memory
varInfo = ncinfo('t.nc','data');
disp(varInfo);
for latInd =1:varInfo.Size(1)
    for lonInd =1:varInfo.Size(2)
        fileName = ['t_ascii_lat',num2str(latInd),'_lon',num2str(lonInd),'.txt'];
        tSeries  = ncread('t.nc','data',[latInd, lonInd, 1],[1,1,varInfo.Size(3)]);
        dlmwrite(fileName,squeeze(tSeries));
    end
end

%% spot check
act = dlmread('t_ascii_lat10_lon29.txt');
exp = squeeze(data(10,29,:));
assert(isequal(act,exp));

【讨论】：