matlab/gmt/rdnctiles_bytile.m

6802a9bca5 Ed H*0001 function [tlist] = rdnctiles_bytile(fall,vit, dlev)
b2dac0a6c3 Ed H*0002
6802a9bca5 Ed H*0003 % Function [tlist] = rdnctiles_bytile(fall,vit, dlev)
b2dac0a6c3 Ed H*0004 %
                0005 % INPUTS
5e40d2151d Ed H*0006 %   fall   cell array of file names
6802a9bca5 Ed H*0007 %   vit    struct containing variable and time information
                0008 %            vit.tdname : "time" dim name (DEF: 'T')
                0009 %            vit.tdname : "time" coord var name (DEF: 'iter')
                0010 %            vit.vars.(vname) : "time" values for each var
                0011 %
                0012 %   dlev   debug level
b2dac0a6c3 Ed H*0013 %
                0014 % OUTPUTS
6802a9bca5 Ed H*0015 %   tlist  struct array of tile data
b2dac0a6c3 Ed H*0016 %
                0017 %
6802a9bca5 Ed H*0018 %  This function has minimal input checking since it is meant to be
b2dac0a6c3 Ed H*0019 %  called by a wrapper function that ensures proper inputs.
                0020 %
                0021 %  Ed Hill
                0022 tlist = struct('gtn',{});
                0023 for fi = 1:length(fall)
                0024   if dlev > 10
                0025     disp(['  Opening : ' char(fall{fi}) ]);
                0026   end
                0027   nc = netcdf(fall{fi},'read');
                0028
                0029   % Get the global tile number
                0030   if isempty(nc.tile_number) || not(isscalar(nc.tile_number(:)))
                0031     error(['No scalar "tile_number" global attribute was found' ...
                0032            ' in file "' fall{fi} '"']);
                0033   end
                0034   gti = nc.tile_number(:);
6802a9bca5 Ed H*0035   itile = find([tlist(:).gtn] == gti);
b2dac0a6c3 Ed H*0036   % tlist
                0037   % [tlist(:).gtn]
                0038   % pause
6802a9bca5 Ed H*0039   if isempty(itile)
                0040     itile = length(tlist) + 1;
                0041     tlist(itile).gtn = gti;
b2dac0a6c3 Ed H*0042   end
                0043
                0044   %  Get all the variables
                0045   if dlev > 10
                0046     fprintf(1,'    reading : ');
                0047   end
6802a9bca5 Ed H*0048   vread = fields(vit.vars);
b2dac0a6c3 Ed H*0049   for iv = 1:length(vread)
                0050     if isempty(nc{char(vread{iv})})
6802a9bca5 Ed H*0051       % disp(['    warning: no var "',vread{iv},'" in "',fall{fi},'"']);
b2dac0a6c3 Ed H*0052       continue
                0053     end
6802a9bca5 Ed H*0054     if dlev > 10
                0055       fprintf(1,[' ' char(vread{iv}) ]);
                0056     end
b2dac0a6c3 Ed H*0057
6802a9bca5 Ed H*0058     if isempty(vit.vars.(vread{iv}))
                0059       % Read all of the variable at once since no time levels are used
b2dac0a6c3 Ed H*0060       tmpv =  nc{vread{iv}}(:);
                0061       sz = size(tmpv);
                0062       nd = length(sz);
6802a9bca5 Ed H*0063       tlist(itile).var.(char(vread{iv})) = permute(tmpv,[nd:-1:1]);
b2dac0a6c3 Ed H*0064     else
6802a9bca5 Ed H*0065       % Get the rank of the time dimension
                0066       trank = 0;
b2dac0a6c3 Ed H*0067       dnames = ncnames(dim(nc{vread{iv}}));
6802a9bca5 Ed H*0068       if strcmp( ncnames(recdim(nc)), vit.tdname )
                0069         m = regexp(dnames, [ '^' vit.tdname '$'] );
b2dac0a6c3 Ed H*0070         for i = 1:length(m)
                0071           if not(isempty(m{i}))
6802a9bca5 Ed H*0072             trank = i;
b2dac0a6c3 Ed H*0073             break
                0074           end
                0075         end
                0076       end
6802a9bca5 Ed H*0077       if trank == 0
                0078         error(['no time dim was found for variable ''' ...
                0079                vread{iv} '''']);
                0080       end
                0081
                0082       % get the corresponding file-local indicies and global-assembly
                0083       % indicies along the time dimension
                0084       loc_times = nc{vit.tvname}(:);
15db930790 Ed H*0085       [v,ind1,ind2] = intersect( vit.vars.(vread{iv}), loc_times );
6802a9bca5 Ed H*0086
                0087       % only read the desired time values based on:
                0088       %   the local  "kt" indicies and
                0089       %   the global "tk" indicies
b2dac0a6c3 Ed H*0090       indstr = '';
                0091       for i = 1:length(dnames)
                0092         if i > 1
                0093          indstr = [ indstr ',' ];
6802a9bca5 Ed H*0094         end
                0095         if i == trank
                0096           indstr = [ indstr 'kt' ];
                0097         else
                0098           indstr = [ indstr ':' ];
                0099         end
b2dac0a6c3 Ed H*0100       end
                0101       rindstr = fliplr(indstr);
15db930790 Ed H*0102       for jj = 1:length(ind1)
                0103         kt = ind2(jj);
b2dac0a6c3 Ed H*0104         eval([ 'tmpv =  nc{vread{iv}}(' indstr ');' ]);
                0105         sz = size(tmpv);
                0106         nd = length(sz);
15db930790 Ed H*0107         tk = ind1(jj);
6802a9bca5 Ed H*0108         comm = [ 'tlist(itile).var.(char(vread{iv}))(' ...
4e7bfad2e0 Ed H*0109                  rindstr ') = permute(tmpv,[nd:-1:1]);' ];
                0110         eval(comm);
b2dac0a6c3 Ed H*0111       end
                0112     end
                0113   end
                0114   if dlev > 10
                0115     fprintf(1,'\n');
                0116   end
                0117
                0118   nc = close(nc);
                0119 end
                0120