gridGliderData

GRIDGLIDERDATA  Grid glider trajectory data over instantaneous homogeneous regular profiles.

  Syntax:
    [DATA_GRID, META_GRID] = GRIDGLIDERDATA(DATA_PROC, META_PROC)
    [DATA_GRID, META_GRID] = GRIDGLIDERDATA(DATA_PROC, META_PROC, OPTIONS)
    [DATA_GRID, META_GRID] = GRIDGLIDERDATA(DATA_PROC, META_PROC, OPT1, VAL1, ...)

  Description:
    DATA_GRID = GRIDGLIDERDATA(DATA_PROC, META_PROC) converts glider trajectory 
    data in struct DATA_PROC to vertical instantaneous profiles defined at 
    regular intervals of depth in DATA_GRID using default option values.
    See options description below.

    DATA_PROC should be a struct in the format returned by PROCESSGLIDERDATA,
    where each field is a vector of readings of the variable with the same 
    name along the glider trajectory. At least it should have a sequence for
    each reference coordinate: time, latitude and longitude, and depth.
    It also should have a sequence of profile indices that flags each reading
    with the number of the cast it belongs to.

    META_PROC is also a struct as returned by PROCESSGLIDERDATA, and gridding
    information is added to any existing metadata of each reference coordinate
    variable or data variable in returned struct META_GRID.

    DATA_GRID is a struct with two kind of fields: bidimensional arrays with 
    profiles of gridded variables as rows, and one dimensional reference 
    coordinate sequences: LATITUDE, LONTGITUDE, DEPTH, TIME and PROFILE_INDEX.
    Coordinate sequences are selected according to preferred choices in 
    options (see below). Only variables selected in options and also present 
    in DATA_PROC are gridded. Selected variables not present in DATA_PROC are
    silently omited.

    Each cast identified in DATA_PROC is converted to an instantaneous 
    vertical profile. The position and time coordinates of the new profile are
    the mean values of the respective coordinates in the cast. All profiles 
    are defined at the same depth coordinates, computed as the depth range of
    the whole trajectory divided into regular intervals of given resolution. 
    The cast data is interpolated over the new depth grid binning the readings 
    that lay in the corresponding depth intervals.

    Options may be given in key-value pairs OPT1, VAL1... or in a struct 
    OPTIONS with field names as option keys and field values as option values.
    Recognized options are:
      PROFILE_LIST: profile index sequence choices.
        String cell array with the names of the sequence to be used as profile 
        index, in order of preference.
        Default value: {'profile_index'}
      TIME_LIST: timestamp sequence choices.
        String cell array with the names of the sequence to be used as time 
        coordinates, in order of preference.
        Default value: {'time'}
      POSITION_LIST: latitude and longitude sequence choices.
        Struct array with the names of the sequence the to be used as 
        latitude and longitude coordinates, in order of preference.
        It should have the following fields:
          LATITUDE: latitude sequence name.
          LONGITUDE: longitude sequence name.
        Default value: struct('latitude',  {'latitude'}, 
                              'longitude', {'longitude'})
      DEPTH_LIST: depth sequence choices.
        String cell array with the names of the sequence to be use as depth
        coordinate, in order of preference.
        Default value: {'depth'}
      DEPTH_STEP: depth resolution.
        Positive number setting the depth resolution for output profiles.
        Default value: 1
      VARIABLE_LIST: list of variables to be included in output profiles.
        String cell array with the names of the variables to be interpolated
        over the output profiles.
        Default value: {} (do nothing except compute profile coordinates)

  Notes:
    This function is an improved version of a previous function by Tomeu Garau
    with the same name. He is the true glider man. Main changes are:
      - Support for reference coordinate sequence selection, variables to
        interpolate, gridding options.
      - Use the mean value of readings lying in the depth interval centered at
        each new depth level with the diameter of the depth resolution
        (instead of interpolation).

  Examples:
    [data_grid, meta_grid] = gridGliderData(data_proc, meta_proc, options)

  See also:
    PROCESSGLIDERDATA

  Authors:
    Joan Pau Beltran  <joanpau.beltran@socib.cat>

0001 function [data_grid, meta_grid] = gridGliderData(data_proc, meta_proc, varargin)
0002 %GRIDGLIDERDATA  Grid glider trajectory data over instantaneous homogeneous regular profiles.
0003 %
0004 %  Syntax:
0005 %    [DATA_GRID, META_GRID] = GRIDGLIDERDATA(DATA_PROC, META_PROC)
0006 %    [DATA_GRID, META_GRID] = GRIDGLIDERDATA(DATA_PROC, META_PROC, OPTIONS)
0007 %    [DATA_GRID, META_GRID] = GRIDGLIDERDATA(DATA_PROC, META_PROC, OPT1, VAL1, ...)
0008 %
0009 %  Description:
0010 %    DATA_GRID = GRIDGLIDERDATA(DATA_PROC, META_PROC) converts glider trajectory
0011 %    data in struct DATA_PROC to vertical instantaneous profiles defined at
0012 %    regular intervals of depth in DATA_GRID using default option values.
0013 %    See options description below.
0014 %
0015 %    DATA_PROC should be a struct in the format returned by PROCESSGLIDERDATA,
0016 %    where each field is a vector of readings of the variable with the same
0017 %    name along the glider trajectory. At least it should have a sequence for
0018 %    each reference coordinate: time, latitude and longitude, and depth.
0019 %    It also should have a sequence of profile indices that flags each reading
0020 %    with the number of the cast it belongs to.
0021 %
0022 %    META_PROC is also a struct as returned by PROCESSGLIDERDATA, and gridding
0023 %    information is added to any existing metadata of each reference coordinate
0024 %    variable or data variable in returned struct META_GRID.
0025 %
0026 %    DATA_GRID is a struct with two kind of fields: bidimensional arrays with
0027 %    profiles of gridded variables as rows, and one dimensional reference
0028 %    coordinate sequences: LATITUDE, LONTGITUDE, DEPTH, TIME and PROFILE_INDEX.
0029 %    Coordinate sequences are selected according to preferred choices in
0030 %    options (see below). Only variables selected in options and also present
0031 %    in DATA_PROC are gridded. Selected variables not present in DATA_PROC are
0032 %    silently omited.
0033 %
0034 %    Each cast identified in DATA_PROC is converted to an instantaneous
0035 %    vertical profile. The position and time coordinates of the new profile are
0036 %    the mean values of the respective coordinates in the cast. All profiles
0037 %    are defined at the same depth coordinates, computed as the depth range of
0038 %    the whole trajectory divided into regular intervals of given resolution.
0039 %    The cast data is interpolated over the new depth grid binning the readings
0040 %    that lay in the corresponding depth intervals.
0041 %
0042 %    Options may be given in key-value pairs OPT1, VAL1... or in a struct
0043 %    OPTIONS with field names as option keys and field values as option values.
0044 %    Recognized options are:
0045 %      PROFILE_LIST: profile index sequence choices.
0046 %        String cell array with the names of the sequence to be used as profile
0047 %        index, in order of preference.
0048 %        Default value: {'profile_index'}
0049 %      TIME_LIST: timestamp sequence choices.
0050 %        String cell array with the names of the sequence to be used as time
0051 %        coordinates, in order of preference.
0052 %        Default value: {'time'}
0053 %      POSITION_LIST: latitude and longitude sequence choices.
0054 %        Struct array with the names of the sequence the to be used as
0055 %        latitude and longitude coordinates, in order of preference.
0056 %        It should have the following fields:
0057 %          LATITUDE: latitude sequence name.
0058 %          LONGITUDE: longitude sequence name.
0059 %        Default value: struct('latitude',  {'latitude'},
0060 %                              'longitude', {'longitude'})
0061 %      DEPTH_LIST: depth sequence choices.
0062 %        String cell array with the names of the sequence to be use as depth
0063 %        coordinate, in order of preference.
0064 %        Default value: {'depth'}
0065 %      DEPTH_STEP: depth resolution.
0066 %        Positive number setting the depth resolution for output profiles.
0067 %        Default value: 1
0068 %      VARIABLE_LIST: list of variables to be included in output profiles.
0069 %        String cell array with the names of the variables to be interpolated
0070 %        over the output profiles.
0071 %        Default value: {} (do nothing except compute profile coordinates)
0072 %
0073 %  Notes:
0074 %    This function is an improved version of a previous function by Tomeu Garau
0075 %    with the same name. He is the true glider man. Main changes are:
0076 %      - Support for reference coordinate sequence selection, variables to
0077 %        interpolate, gridding options.
0078 %      - Use the mean value of readings lying in the depth interval centered at
0079 %        each new depth level with the diameter of the depth resolution
0080 %        (instead of interpolation).
0081 %
0082 %  Examples:
0083 %    [data_grid, meta_grid] = gridGliderData(data_proc, meta_proc, options)
0084 %
0085 %  See also:
0086 %    PROCESSGLIDERDATA
0087 %
0088 %  Authors:
0089 %    Joan Pau Beltran  <joanpau.beltran@socib.cat>
0090 
0091 %  Copyright (C) 2013-2016
0092 %  ICTS SOCIB - Servei d'observacio i prediccio costaner de les Illes Balears
0093 %  <http://www.socib.es>
0094 %
0095 %  This program is free software: you can redistribute it and/or modify
0096 %  it under the terms of the GNU General Public License as published by
0097 %  the Free Software Foundation, either version 3 of the License, or
0098 %  (at your option) any later version.
0099 %
0100 %  This program is distributed in the hope that it will be useful,
0101 %  but WITHOUT ANY WARRANTY; without even the implied warranty of
0102 %  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
0103 %  GNU General Public License for more details.
0104 %
0105 %  You should have received a copy of the GNU General Public License
0106 %  along with this program.  If not, see <http://www.gnu.org/licenses/>.
0107 
0108   error(nargchk(2, 16, nargin, 'struct'));  
0109   
0110   
0111   %% Set gridding options from default values and extra arguments.
0112   % Set default option values.
0113   options = struct();
0114   options.profile_list = {'profile_index'};
0115   options.time_list = {'time'};
0116   options.position_list.latitude = 'latitude';
0117   options.position_list.longitude = 'longitude';
0118   options.depth_list = {'depth'};
0119   options.depth_step = 1;
0120   options.variable_list = {};
0121   % Parse option key-value pairs in any accepted call signature.
0122   if isscalar(varargin) && isstruct(varargin{1})
0123     % Options passed as a single option struct argument:
0124     % field names are option keys and field values are option values.
0125     option_key_list = fieldnames(varargin{1});
0126     option_val_list = struct2cell(varargin{1});
0127   elseif mod(numel(varargin), 2) == 0
0128     % Options passed as key-value argument pairs.
0129     option_key_list = varargin(1:2:end);
0130     option_val_list = varargin(2:2:end);
0131   else
0132     error('glider_toolbox:processGliderData:InvalidOptions', ...
0133           'Invalid optional arguments (neither key-value pairs nor struct).');
0134   end
0135   % Overwrite default options with values given in extra arguments.
0136   for opt_idx = 1:numel(option_key_list)
0137     opt = lower(option_key_list{opt_idx});
0138     val = option_val_list{opt_idx};
0139     if isfield(options, opt)
0140       options.(opt) = val;
0141     else
0142       error('glider_toolbox:gridGliderData:InvalidOption', ...
0143             'Invalid option: %s.', opt);
0144     end
0145   end
0146 
0147   
0148   %% Get list of sequences in trajectory data.
0149   sequence_list = fieldnames(data_proc);
0150   
0151   
0152   %% Select profile index sequence and coordinate sequences.
0153   % Profile index and coordinate variables (latitude, longitude, depth and time)
0154   % are mandatory.
0155   profile_available = false;
0156   time_available = false;
0157   position_available = false;
0158   depth_available = false;
0159   % Select profile index sequence.
0160   profile_choice_list = cellstr(options.profile_list);
0161   for profile_choice_idx = 1:numel(profile_choice_list)
0162     profile_sequence = profile_choice_list{profile_choice_idx};
0163     if ismember(profile_sequence, sequence_list) ...
0164         && ~all(isnan(data_proc.(profile_sequence)))
0165       fprintf('Selected profile index coordinate sequence:\n');
0166       fprintf('  profile: %s\n', profile_sequence);
0167       profile_available = true;
0168       break
0169     end
0170   end
0171   % Select time coordinate sequence.
0172   time_choice_list = cellstr(options.time_list);
0173   for time_choice_idx = 1:numel(time_choice_list)
0174     time_sequence = time_choice_list{time_choice_idx};
0175     if ismember(time_sequence, sequence_list) ...
0176         && any(data_proc.(time_sequence) > 0)
0177       fprintf('Selected time coordinate sequence:\n');
0178       fprintf('  time: %s\n', time_sequence);
0179       time_available = true;
0180       break
0181     end
0182   end
0183   % Select position coordinate sequences.
0184   position_choice_list = options.position_list;
0185   for position_choice_idx = 1:numel(position_choice_list)
0186     latitude_sequence = position_choice_list(position_choice_idx).latitude;
0187     longitude_sequence = position_choice_list(position_choice_idx).longitude;
0188     if all(ismember({latitude_sequence longitude_sequence}, sequence_list)) ...
0189         && ~all(isnan(data_proc.(latitude_sequence))) ...
0190         && ~all(isnan(data_proc.(longitude_sequence))) 
0191       fprintf('Selected position coordinate sequences:\n');
0192       fprintf('  longitude: %s\n', latitude_sequence);
0193       fprintf('  latitude : %s\n', longitude_sequence);
0194       position_available = true;
0195       break
0196     end
0197   end
0198   % Select depth sequence.
0199   depth_choice_list = cellstr(options.depth_list);
0200   for depth_choice_idx = 1:numel(depth_choice_list)
0201     depth_sequence = depth_choice_list{depth_choice_idx};
0202     if ismember(depth_sequence, sequence_list) ...
0203         && ~all(isnan(data_proc.(depth_sequence)))
0204       fprintf('Selected depth coordinate sequence:\n');
0205       fprintf('  depth: %s\n', depth_sequence);
0206       depth_available = true;
0207       break
0208     end
0209   end
0210   % Check all required inputs are present.
0211   coordinate_sequence_available = ...
0212     [profile_available time_available position_available depth_available];
0213   if ~all(coordinate_sequence_available)
0214     coordinate_sequences = {'profile' 'time' 'position' 'depth'};
0215     miss_coords = coordinate_sequences(~coordinate_sequence_available);
0216     miss_coords_str = [sprintf('%s, ', miss_coords{1:end-1}) miss_coords{end}];
0217     error('glider_toolbox:processGliderData:MissingCoordinateSequence', ...
0218           'Missing coordinate sequences in data set: %s.', miss_coords_str);
0219   end
0220 
0221 
0222   %% Select variables to grid.
0223   variable_list_choice = cellstr(options.variable_list);
0224   variable_name_list = intersect(sequence_list, variable_list_choice);
0225   fprintf('Selected variables to interpolate:\n');
0226   fprintf('  %s\n', variable_name_list{:});
0227 
0228   
0229   %% Store variables as columns in a single array to accelerate binning below.
0230   profile = data_proc.(profile_sequence);
0231   time = data_proc.(time_sequence);
0232   latitude = data_proc.(latitude_sequence);
0233   longitude = data_proc.(longitude_sequence);
0234   depth = data_proc.(depth_sequence);
0235   num_variables = numel(variable_name_list);
0236   num_instants = numel(time);
0237   variables = nan(num_instants, num_variables);
0238   for variable_name_idx = 1:num_variables
0239     variable_name = variable_name_list{variable_name_idx};
0240     variables(:, variable_name_idx) = data_proc.(variable_name)(:);
0241   end
0242 
0243   
0244   %% Compute number of casts.
0245   num_casts = fix(max(profile));
0246   profile_range = (1:num_casts);
0247 
0248 
0249   %% Compute depth intervals.
0250   depth_resolution = options.depth_step;
0251   depth_min = round(min(depth) / depth_resolution) * depth_resolution;
0252   depth_max = round(max(depth) / depth_resolution) * depth_resolution;
0253   depth_range = depth_min : depth_resolution : depth_max;
0254   num_levels = numel(depth_range);
0255   
0256   
0257   %% Initialize output.
0258   data_grid.depth = depth_range(:);
0259   data_grid.profile_index = profile_range(:);
0260   data_grid.time = nan(num_casts, 1);
0261   data_grid.longitude = nan(num_casts, 1);
0262   data_grid.latitude = nan(num_casts, 1);
0263   for variable_name_idx = 1:numel(variable_name_list)
0264     variable_name = variable_name_list{variable_name_idx};
0265     data_grid.(variable_name) = nan(num_casts, num_levels);
0266   end
0267   
0268 
0269   %% Compute profile coordinates and profile data.
0270   % Spatial and temporal coordinates are the mean values among cast readings.
0271   % Selected variable data is interpolated at selected depth levels.
0272   %{
0273   for cast_idx = 1:num_casts
0274     cast_select = (profile == cast_idx);
0275     cast_lat = latitude(cast_select);
0276     cast_lon = longitude(cast_select);
0277     cast_depth = depth(cast_select);
0278     cast_time = time(cast_select);
0279     data_grid.time(cast_idx) = nanmean(cast_time);
0280     data_grid.latitude(cast_idx) = nanmean(cast_lat);
0281     data_grid.longitude(cast_idx) = nanmean(cast_lon);
0282     for variable_name_idx = 1:numel(variable_name_list)
0283       variable_name = variable_name_list{variable_name_idx};
0284       cast_variable = data_proc.(variable_name)(cast_select);
0285       cast_valid = ~(isnan(cast_depth(:)) | isnan(cast_variable(:)));
0286       if sum(cast_valid) > 2
0287         data_grid.(variable_name)(cast_idx, :) = ...
0288           interp1(cast_depth(cast_valid), cast_variable(cast_valid), ...
0289           depth_range(:));
0290       end
0291     end
0292   end
0293   %}
0294   %%{
0295   % Spatial and temporal coordinates are the mean values among cast readings.
0296   % Selected variable data is binned taking the mean values of readings in depth
0297   % intervals centered at selected depth levels.
0298   % For better performance, compute variable data in single array and move it to
0299   % output struct at the end.
0300   fprintf('Gridding variables with settings:\n');
0301   fprintf('  depth level min : %d\n', depth_min);
0302   fprintf('  depth level max : %d\n', depth_max);
0303   fprintf('  depth level step: %d\n', depth_resolution);
0304   fprintf('  number of depth levels: %d\n', num_levels);
0305   fprintf('  number of profiles    : %d\n', num_casts);
0306   fprintf('  number of variables   : %d\n', num_variables);
0307   data_grid_variables = nan(num_casts, num_levels, num_variables);
0308   for cast_idx = 1:num_casts
0309     cast_select = (profile == cast_idx);
0310     cast_lat = latitude(cast_select);
0311     cast_lon = longitude(cast_select);
0312     cast_depth = depth(cast_select);
0313     cast_time = time(cast_select);
0314     cast_variables = variables(cast_select, :);
0315     data_grid.time(cast_idx) = nanmean(cast_time);
0316     data_grid.latitude(cast_idx) = nanmean(cast_lat);
0317     data_grid.longitude(cast_idx) = nanmean(cast_lon);
0318     if ~isempty(cast_variables) % Speed up when there are no variables.
0319       data_grid_variables(cast_idx, :, :) = ...
0320         cell2mat(arrayfun(@(d) nanmean(cast_variables(abs(cast_depth-d)<=0.5*depth_resolution, :), 1), ...
0321                           depth_range(:), 'UniformOutput', false));
0322     end
0323   end
0324   % Move binned variable data to output struct.
0325   for variable_name_idx = 1:num_variables
0326     variable_name = variable_name_list{variable_name_idx};
0327     data_grid.(variable_name) = data_grid_variables(:, :, variable_name_idx);
0328   end
0329   %%}
0330 
0331   
0332   %% Add gridding metadata:
0333   meta_grid.profile_index = meta_proc.(profile_sequence);
0334   meta_grid.profile_index.grid_sources = profile_sequence;
0335   meta_grid.profile_index.grid_resolution = 1;
0336   meta_grid.profile_index.grid_min = min(profile_range);
0337   meta_grid.profile_index.grid_max = max(profile_range);
0338   meta_grid.depth = meta_proc.(depth_sequence);
0339   meta_grid.depth.grid_sources = depth_sequence;
0340   meta_grid.depth.grid_resolution = depth_resolution;
0341   meta_grid.depth.grid_min = depth_min;
0342   meta_grid.depth.grid_max = depth_max;
0343   meta_grid.time = meta_proc.(time_sequence);
0344   meta_grid.time.grid_sources = time_sequence;
0345   meta_grid.time.grid_coordinates = {'profile_index'};
0346   meta_grid.time.grid_method = {'mean'};
0347   meta_grid.longitude = meta_proc.(longitude_sequence);
0348   meta_grid.longitude.grid_sources = longitude_sequence;
0349   meta_grid.longitude.grid_coordinates = {'profile_index'};
0350   meta_grid.longitude.grid_method = {'mean'};
0351   meta_grid.latitude = meta_proc.(latitude_sequence);
0352   meta_grid.latitude.grid_sources = latitude_sequence;
0353   meta_grid.latitude.grid_coordinates = {'profile_index'};
0354   meta_grid.latitude.grid_method = {'mean'};
0355   for variable_name_idx = 1:numel(variable_name_list)
0356     variable_name = variable_name_list{variable_name_idx};
0357     meta_grid.(variable_name) = meta_proc.(variable_name);
0358     meta_grid.(variable_name).grid_sources = variable_name;
0359     meta_grid.(variable_name).grid_coordinates = {'profile_index' 'depth'};
0360     meta_grid.(variable_name).grid_method = {'index' 'mean'};
0361   end
0362   
0363 end