Changeset 27875 for issm/trunk/src/m/contrib/musselman/read_netCDF.m

Timestamp:

08/10/23 18:34:53 (20 months ago)

Author:

musselman

Message:

Added support in matlab read/write_netCDF files for multidimensional structure arrays.

File:

• issm/trunk/src/m/contrib/musselman/read_netCDF.m (modified) (7 diffs)

issm/trunk/src/m/contrib/musselman/read_netCDF.m

@@ -17 +17 @@
 
 function model_copy = read_netCDF(filename)
-    fprintf('NetCDF42C v1.1.13\n');
+    fprintf('NetCDF42C v1.1.14\n');
 
     % make a model framework to fill that is in the scope of this file
@@ -109 +109 @@
         for variable = variables
             [varname, xtype, dimids, numatts] = netcdf.inqVar(group_location_in_file, variable);
-            copy_variable_data_to_new_model(group_location_in_file,varname, xtype);
+
+            % keep an eye out for nested structs:
+            if strcmp(varname, 'this_is_a_nested')
+                is_nested = true;
+                copy_nested_struct(group_location_in_file)
+            else
+                is_nested = false;
+                copy_variable_data_to_new_model(group_location_in_file,varname, xtype);
+            end
         end
     catch ME
@@ -117 +125 @@
     % try to find groups in current level, if it doesn't work it's because there is nothing there
     try
-        % search for nested groups in the current level to feed back to this function
-        groups = netcdf.inqGrps(group_location_in_file);
-        if not(isempty(groups))
-            for group = groups
-                %disp('found nested group!!')
-                group_id = netcdf.inqNcid(group_location_in_file, netcdf.inqGrpName(group));
-                %disp(netcdf.inqGrpNameFull(group_id))
-                walk_nested_groups(group);
-            end
-        end
-    catch
-    end
-end
+        % if it's not a nested struct, keep searching for subgroups
+        if isnested
+            % do nothing
+        else
+            % search for nested groups in the current level to feed back to this function
+            groups = netcdf.inqGrps(group_location_in_file);
+            if not(isempty(groups))
+                for group = groups
+                    %disp('found nested group!!')
+                    group_id = netcdf.inqNcid(group_location_in_file, netcdf.inqGrpName(group));
+                    %disp(netcdf.inqGrpNameFull(group_id))
+                    walk_nested_groups(group);
+                end
+            end
+        end
+    catch % no nested groups here
+    end
+end
+
+
+
+function copy_nested_struct(group_location_in_file)
+    global model_copy;
+    global NCData;
+    %{
+        A common multidimensional struct array is the 1xn md.results.TransientSolution struct. 
+        The process to recreate is as follows:
+            1. Get the name of the struct from group variable name_of_struct
+            2. Get the fieldnames from the subgroups 
+            3. Recreate the struct with fieldnames 
+            4. Populate the fields with their respective values
+    %}
+
+    % step 1
+    varid = netcdf.inqVarID(group_location_in_file, 'name_of_struct');
+    name_of_struct = netcdf.getVar(group_location_in_file, varid)';
+
+    % step 2
+    subgroups = netcdf.inqGrps(group_location_in_file); % numerical cell array with ID's of subgroups
+    % get single subgroup's data
+    single_subgroup_ID = subgroups(1);
+    subgroup_varids = netcdf.inqVarIDs(single_subgroup_ID);
+    fieldnames = {};
+    for variable = subgroup_varids
+        [varname, xtype, dimids, numatts] = netcdf.inqVar(single_subgroup_ID, variable);
+        fieldnames{end+1} = varname;
+    end
+
+    % step 3
+    address_in_model = strrep(netcdf.inqGrpNameFull(group_location_in_file), '/', '');
+    % we cannot assign a variable to represent this object as MATLAB treats all variables as copies
+    % and not pointers to the same memory address
+    % this means that if address_in_model is more than 1 item, we need to modify the code. For now, 
+    % we just hope this will do
+    
+    model_copy.(address_in_model).(name_of_struct) = struct();
+    
+    % for every fieldname in the subgroup, create an empty field
+    for fieldname = string(fieldnames)
+        model_copy.(address_in_model).(name_of_struct).(fieldname) = {};
+    end
+
+    % use repmat to make the struct array multidimensional along the fields axis
+    number_of_dimensions = numel(subgroups);
+    model_copy.(address_in_model).(name_of_struct) = repmat(model_copy.(address_in_model).(name_of_struct), 1, number_of_dimensions);
+    
+    % step 4
+    % for every layer of the multidimensional struct array, populate the fields
+    for current_layer = 1:number_of_dimensions
+        % choose subgroup
+        current_layer_subgroup_ID = subgroups(current_layer);
+        % get all vars
+        current_layer_subgroup_varids = netcdf.inqVarIDs(current_layer_subgroup_ID);
+        % get individual vars and set fields at layer current_layer
+        for varid = current_layer_subgroup_varids
+            [varname, xtype, dimids, numatts] = netcdf.inqVar(current_layer_subgroup_ID, varid);
+            data = netcdf.getVar(current_layer_subgroup_ID, varid);
+
+            % netcdf uses Row Major Order but MATLAB uses Column Major Order so we need to transpose all arrays w/ more than 1 dim
+            if all(size(data)~=1) || xtype == 2
+                data = data.';
+            end
+            
+            % set the field
+            model_copy.(address_in_model).(name_of_struct)(current_layer).(varname) = data;
+            %address_to_struct_in_model = setfield(address_to_struct_in_model(current_layer), varname, data)
+        end
+        model_copy.(address_in_model).(name_of_struct)(current_layer);
+        fprintf("Successfully saved layer %s to multidimension struct array\n", num2str(current_layer))
+    end
+    fprintf('Successfully recreated multidimensional structure array %s in md.%s\n', name_of_struct, address_in_model)
+end
+
+
+
 
 %{
@@ -142 +232 @@
     %disp(varname)
     % this is an inversion band-aid
-    if strcmp(varname, 'inversion_class_name')
+    if strcmp(varname, 'inversion_class_name') || strcmp(varname, 'name_of_struct') || strcmp(varname, 'solution')
         % we don't need this
     else
@@ -149 +239 @@
             %disp(netcdf.inqGrpNameFull(group_location_in_file))
             %disp(class(netcdf.inqGrpNameFull(group_location_in_file)))
-            adress_to_attr = strrep(netcdf.inqGrpNameFull(group_location_in_file), '/', '.');
-            
-            data = netcdf.getVar(group_location_in_file, netcdf.inqVarID(group_location_in_file, varname));
-            
-    
-            % matlab needs to know that ' ' = char()
+            address_to_attr = strrep(netcdf.inqGrpNameFull(group_location_in_file), '/', '.');
+            varid = netcdf.inqVarID(group_location_in_file, varname);
+            data = netcdf.getVar(group_location_in_file, varid);
+            
+    
+            % if we have an empty string
             if xtype == 2 && isempty(all(data))
                 data = cell(char());
+            % if we have an empty cell-char array
             elseif numel(data) == 1 && xtype == 3 && data == -32767
                 data = cell(char());
             end
-            % band-aid for cell-char-arrays:
+            % band-aid for some cell-char-arrays:
             if xtype == 2 && strcmp(data, 'default')
                 data = {'default'};
@@ -168 +259 @@
             if all(size(data)~=1) || xtype == 2
                 data = data.';
+            end
+
+            % if we have a list of strings
+            if xtype == 2
+                try
+                    if strcmp(netcdf.getAtt(group, varid, "type_is"), 'cell_array_of_strings')
+                        data = cellstr(data)
+                    end
+                catch
+                    % no attr found so we pass
+                end
             end
             
@@ -174 +276 @@
             %xtype
             if xtype == 10
-                arg_to_eval = ['model_copy', adress_to_attr, '.', varname, ' = ' , 'double(data);'];
+                arg_to_eval = ['model_copy', address_to_attr, '.', varname, ' = ' , 'double(data);'];
                 eval(arg_to_eval);
                 %disp('saved int64 as int16')
             else
-                arg_to_eval = ['model_copy', adress_to_attr, '.', varname, ' = ' , 'data;'];
+                arg_to_eval = ['model_copy', address_to_attr, '.', varname, ' = ' , 'data;'];
                 eval(arg_to_eval);
             end