get_POLCOMS_netCDF.m 4.63 KB
Newer Older
1
function ncdata = get_POLCOMS_netCDF(files, varlist)
2 3 4
% Read temperature and salinity from NetCDF model output files and
% interpolate onto the open boundaries in Mobj.
%
5
% function struct = get_POLCOMS_netCDF(Mobj, files, varlist)
6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21
%
% DESCRIPTION:
%    Extract variables in varlist to a struct from the files given in
%    files.
%
% INPUT:
%   files   = Cell array of NetCDF file(s).
%   varlist = Cell array of variables names.
%
% OUTPUT:
%    Struct in which the field names are the variable names supplied in var
%    list. Each field has a data and units array which contain the data and
%    units variables from the NetCDF. If no units are found, it is left
%    blank for that variable.
%
% EXAMPLE USAGE
22
%    S = get_POLCOMS_netCDF({'/tmp/2000.nc', '/tmp/2001.nc', {'temp', 'salt'})
23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45
%
% NOTES:
%
%   - If you supply multiple files, there are a few assumptions:
%
%       - Variables are only appended if there are 3 or 4 dimensions; fewer
%       than that, and the values are assumed to be static across all the
%       given files (e.g. longitude, latitude etc.). The last dimension
%       is assumed to be time.
%       - The order of the files given should be chronological.
% 
%   - This has been tested on NetCDF files generated from the PML
%   POLCOMS-ERSEM daily mean model output
%
% Author(s):
%    Pierre Cazenave (Plymouth Marine Laboratory)
%
% Revision history
%    2013-02-08 First version based on restart_FVCOM_AMM.m and
%    get_AMM_tsobc.m.
%
%==========================================================================

46
subname = 'get_POLCOMS_netCDF';
47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64

global ftbverbose;
if ftbverbose
    fprintf('\n')
    fprintf(['begin : ' subname '\n'])
end

% Get the results. Check we have a cell array, and if we don't, assume it's
% a file name.

if iscell(files)
    todo = length(files);
else
    todo = 1;
end

for ii = 1:todo

65 66 67 68 69
    if iscell(files)
        ftn = files{ii};
    else
        ftn = files;
    end
70

71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115
    if ftbverbose
        % Strip path from filename for the verbose output.
        [~, basename, ext] = fileparts(ftn);
        tmp_fn = [basename, ext];

        if todo == 1
            fprintf('%s: extracting file %s... ', subname, tmp_fn)
        else
            fprintf('%s: extracting file %s (%i of %i)... ', subname, tmp_fn, ii, todo)
        end
    end

    nc = netcdf.open(ftn, 'NOWRITE');

    for var = 1:numel(varlist)

        getVar = varlist{var};
        varid = netcdf.inqVarID(nc, getVar);

        data = double(netcdf.getVar(nc, varid, 'single'));
        if ii == 1
            ncdata.(getVar).data = data;
        else
            if ndims(data) < 3
                if strcmpi(varlist{var}, 'time')
                    % If the dimension is time, we need to be a bit more
                    % clever since we'll need a concatenated time series
                    % (in which values are continuous and from which we
                    % can extract a sensible time). As such, we need to add
                    % the maximum of the existing time. On the first
                    % iteration, we should save ourselves the base time
                    % (from the units of the time variable).
                    ncdata.(getVar).data = [ncdata.(getVar).data; data + max(ncdata.(getVar).data)];
                else
                    % This should be a fixed set of values (probably lon or
                    % lat) in which case we don't need to append them, so
                    % just replace the existing values with those in the
                    % current NetCDF file.
                    ncdata.(getVar).data = data;
                end
            elseif ndims(data) == 3 || ndims(data) == 4
                % Concatenate along the last dimension and hope/assume it's
                % a time dimension.
                ncdata.(getVar).data = cat(ndims(data), ncdata.(getVar).data, data);
            else
116
                error('Unsupported number of dimensions in PML POLCOMS-ERSEM data')
117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146
            end
        end
        % Try to get some units (important for the calculation of MJD).
        try
            if ii == 1
                units = netcdf.getAtt(nc, varid, 'units');
            else
                % Leave the units values alone so we always use the values
                % from the first file. This is particularly important for
                % the time calculation later on which is dependent on
                % knowing the time origin of the first file.
                continue
            end
        catch
            units = [];
        end
        ncdata.(getVar).units = units;
    end

    netcdf.close(nc)

    if ftbverbose
        fprintf('done.\n')
    end

end

if ftbverbose
    fprintf(['end   : ' subname '\n'])
end