Index: /issm/trunk-jpl/src/m/archive/arch.py =================================================================== --- /issm/trunk-jpl/src/m/archive/arch.py (revision 21165) +++ /issm/trunk-jpl/src/m/archive/arch.py (revision 21166) @@ -228,5 +228,5 @@ elif data_type==3: data_size='{0}x{1}'.format(rows,cols) - data_type_str='vector' + data_type_str='vector/matrix' result=OrderedDict() @@ -249,5 +249,5 @@ 1 : string 2 : double (scalar) - 3 : vector (of type double) + 3 : vector or matrix (of type double) """ Index: /issm/trunk-jpl/src/m/archive/archdisp.m =================================================================== --- /issm/trunk-jpl/src/m/archive/archdisp.m (revision 21165) +++ /issm/trunk-jpl/src/m/archive/archdisp.m (revision 21166) @@ -43,5 +43,5 @@ archive_data{end+1}=fread(fid,1,'double','ieee-be'); elseif field_type==3 - field_types{end+1}='vector'; + field_types{end+1}='vector/matrix'; rows=fread(fid,1,'int','ieee-be'); cols=fread(fid,1,'int','ieee-be'); Index: /issm/trunk-jpl/src/m/archive/archread.m =================================================================== --- /issm/trunk-jpl/src/m/archive/archread.m (revision 21165) +++ /issm/trunk-jpl/src/m/archive/archread.m (revision 21166) @@ -76,5 +76,5 @@ elseif isscalar(format) code=2; % scalar - elseif isvector(format) + elseif isvector(format) or ismatrix(format) code=3; % vector else Index: /issm/trunk-jpl/src/m/archive/archwrite.m =================================================================== --- /issm/trunk-jpl/src/m/archive/archwrite.m (revision 21165) +++ /issm/trunk-jpl/src/m/archive/archwrite.m (revision 21166) @@ -51,8 +51,8 @@ elseif isscalar(format) code=2; % scalar - elseif isvector(format) - code=3; % vector + elseif (isvector(format) || ismatrix(format)) + code=3; % vector or matrix else - error('Error! Please ensure arguments are strings, scalars, or vectors.'); + error('Error! Please ensure arguments are strings, scalars, or vectors/matrixes.'); end end%}}} Index: /issm/trunk-jpl/test/Data/convertmattoarch.m =================================================================== --- /issm/trunk-jpl/test/Data/convertmattoarch.m (revision 21166) +++ /issm/trunk-jpl/test/Data/convertmattoarch.m (revision 21166) @@ -0,0 +1,47 @@ +function []=convertmattoarch(matfile,archfile); +% convertmattoarch -- function to convert mat-format file to arch-format file +% +% usage: +% convertmattoarch(matfile,archfile); +% where: +% matfile name of mat-format file +% archfile name of arch-format file (optional) +% + + if ~exist('matfile','var') || isempty(matfile) + help convertmattoarch + error('convertmattoarch usage error.'); + end + + [pathstr,name,ext]=fileparts(matfile); + if isempty(ext) + ext='.mat'; + end + matfile=fullfile(pathstr,[name ext]); + + if ~exist('archfile','var') || isempty(archfile) + archfile=fullfile(pathstr,[name '.arch']); + end + + if exist(archfile,'file') + delete(archfile); + end + + a=load(matfile,'-mat'); + disp(sprintf('mat-format file ''%s'' read.',matfile)); + fnames=fieldnames(a); + + for i=1:length(fnames) + if isstruct(a.(fnames{i})) || iscell(a.(fnames{i})) + warning('field ''%s'' is of class ''%s'' and will not be written.',fnames{i},class(a.(fnames{i}))); + else + % matlab writes the dimensions reversed and matrices transposed into binary files, so compensate for that + archwrite(archfile,fnames{i},transpose(a.(fnames{i}))); + disp(sprintf('field ''%s'' of class ''%s'' and size [%dx%d] written.',... + fnames{i},class(a.(fnames{i})),size(a.(fnames{i}),1),size(a.(fnames{i}),2))); + end + end + disp(sprintf('arch-format file ''%s'' written.',archfile)); + +end + Index: sm/trunk-jpl/test/Data/convertmattonc.m =================================================================== --- /issm/trunk-jpl/test/Data/convertmattonc.m (revision 21165) +++ (revision ) @@ -1,50 +1,0 @@ -function []=convertmattonc(matfile,ncfile); -% convertmattonc -- function to convert mat-format file to nc-format file -% -% usage: -% convertmattonc(matfile,ncfile); -% where: -% matfile name of mat-format file -% ncfile name of nc-format file (optional) -% - - if ~exist('matfile','var') || isempty(matfile) - help convertmattonc - error('convertmattonc usage error.'); - end - - [pathstr,name,ext]=fileparts(matfile); - if isempty(ext) - ext='.mat'; - end - matfile=fullfile(pathstr,[name ext]); - - if ~exist('ncfile','var') || isempty(ncfile) - ncfile=fullfile(pathstr,[name '.nc']); - end - - if exist(ncfile,'file') - delete(ncfile); - end - - a=load(matfile,'-mat'); - disp(sprintf('mat-format file ''%s'' read.',matfile)); - fnames=fieldnames(a); - - for i=1:length(fnames) - if isstruct(a.(fnames{i})) || iscell(a.(fnames{i})) - warning('field ''%s'' is of class ''%s'' and will not be written.',fnames{i},class(a.(fnames{i}))); - else - % matlab writes the dimensions reversed and matrices transposed into netcdf, so compensate for that - nccreate(ncfile,fnames{i},... - 'Dimensions',{[fnames{i} '_2'] size(a.(fnames{i}),2) [fnames{i} '_1'] size(a.(fnames{i}),1)},... - 'Format','classic'); - ncwrite(ncfile,fnames{i},transpose(a.(fnames{i}))); - disp(sprintf('field ''%s'' of class ''%s'' and size [%dx%d] written.',... - fnames{i},class(a.(fnames{i})),size(a.(fnames{i}),1),size(a.(fnames{i}),2))); - end - end - disp(sprintf('nc-format file ''%s'' written.',ncfile)); - -end - Index: /issm/trunk-jpl/test/Data/loadarch.m =================================================================== --- /issm/trunk-jpl/test/Data/loadarch.m (revision 21166) +++ /issm/trunk-jpl/test/Data/loadarch.m (revision 21166) @@ -0,0 +1,67 @@ +function [s]=loadarch(archfile); +%LOADARCH - Function to read the variables from an arch-formatted file +% +% Usage: +% s=loadarch(archfile); +% where: +% archfile name of an arch-format file + + if ~exist('archfile','var') || isempty(archfile) + help loadarch + error('loadarch usage error.'); + end + + [pathstr,name,ext]=fileparts(archfile); + if isempty(ext) + ext='.arch'; + end + archfile=fullfile(pathstr,[name ext]); + + % Fields for our structure 's' + variable_names={}; + variable_sizes={}; + variable_types={}; + variable_data={}; + + % Read file and gather variable names and sizes + fid=fopen(archfile,'rb'); + while ~feof(fid) + [reclen,count]=fread(fid,1,'int','ieee-be'); + if count == 0, % reached eof + break; + end + % Read variable traits + fread(fid,1,'int','ieee-be'); + variable_name_len=fread(fid,1,'int','ieee-be'); + variable_names{end+1}=char(fread(fid,variable_name_len,'char','ieee-be')'); + variable_reclen=fread(fid,1,'int','ieee-be'); + variable_type=fread(fid,1,'int','ieee-be'); + % set related data type + if variable_type==2 + variable_types{end+1}='double'; + variable_sizes{end+1}='1x1'; + variable_data{end+1}=fread(fid,1,'double','ieee-be'); % read data + elseif variable_type==3 + variable_types{end+1}='double vector / matrix'; + rows=fread(fid,1,'int','ieee-be'); + cols=fread(fid,1,'int','ieee-be'); + variable_data{end+1}=fread(fid,[rows,cols],'double','ieee-be'); + variable_size_cat=strcat(num2str(rows),'x',num2str(cols)); + variable_sizes{end+1}=variable_size_cat; + else + fclose(fid); + error('Error: Encountered invalid data type while loading arch information.'); + end + end + fclose(fid); + disp(sprintf('arch-format file ''%s'' read.',archfile)); + + % Relate the variable fields to the structure + % Access fields by doing (for example): s(1).Name -> 'VelocityXObservation' + s=struct('FileName',archfile,'Name',variable_names,'Size',variable_sizes,'Type',variable_types,'Data',variable_data); + % Let user know data has been successfully read + for i=1:numel(variable_names) + disp(sprintf('field ''%s'' of class ''%s'' and size [%s] read.',... + variable_names{i},variable_types{i},variable_sizes{i})); + end +end Index: sm/trunk-jpl/test/Data/loadnc.m =================================================================== --- /issm/trunk-jpl/test/Data/loadnc.m (revision 21165) +++ (revision ) @@ -1,32 +1,0 @@ -function [s]=loadnc(ncfile); -% loadnc -- function to read the variables from an nc-format file -% -% usage: -% s=loadnc(ncfile); -% where: -% ncfile name of nc-format file -% - - if ~exist('ncfile','var') || isempty(ncfile) - help loadnc - error('loadnc usage error.'); - end - - [pathstr,name,ext]=fileparts(ncfile); - if isempty(ext) - ext='.nc'; - end - ncfile=fullfile(pathstr,[name ext]); - - a=ncinfo(ncfile); - disp(sprintf('nc-format file ''%s'' read.',ncfile)); - - for i=1:length(a.Variables) - % matlab reads the dimensions reversed and matrices transposed from netcdf, so compensate for that - s.(a.Variables(i).Name)=transpose(ncread(ncfile,a.Variables(i).Name)); - disp(sprintf('field ''%s'' of class ''%s'' and size [%dx%d] read.',... - a.Variables(i).Name,class(s.(a.Variables(i).Name)),size(s.(a.Variables(i).Name),1),size(s.(a.Variables(i).Name),2))); - end - -end - Index: /issm/trunk-jpl/test/Par/79North.par =================================================================== --- /issm/trunk-jpl/test/Par/79North.par (revision 21165) +++ /issm/trunk-jpl/test/Par/79North.par (revision 21166) @@ -2,11 +2,11 @@ %Geometry and observation -x=transpose(ncread('../Data/79North.nc','x')); -y=transpose(ncread('../Data/79North.nc','y')); -vx=transpose(ncread('../Data/79North.nc','vx')); -vy=transpose(ncread('../Data/79North.nc','vy')); -index=transpose(ncread('../Data/79North.nc','index')); -surface=transpose(ncread('../Data/79North.nc','surface')); -thickness=transpose(ncread('../Data/79North.nc','thickness')); +x=archread('../Data/79North.nc','x'); +y=archread('../Data/79North.nc','y'); +vx=archread('../Data/79North.nc','vx'); +vy=archread('../Data/79North.nc','vy'); +index=archread('../Data/79North.nc','index'); +surface=archread('../Data/79North.nc','surface'); +thickness=archread('../Data/79North.nc','thickness'); md.initialization.vx =InterpFromMeshToMesh2d(index,x,y,vx,md.mesh.x,md.mesh.y); md.initialization.vy =InterpFromMeshToMesh2d(index,x,y,vy,md.mesh.x,md.mesh.y); Index: /issm/trunk-jpl/test/Par/GiaBenchmarksAB.par =================================================================== --- /issm/trunk-jpl/test/Par/GiaBenchmarksAB.par (revision 21165) +++ /issm/trunk-jpl/test/Par/GiaBenchmarksAB.par (revision 21166) @@ -27,9 +27,9 @@ %Initial velocity -x = transpose(ncread('../Data/SquareSheetConstrained.nc','x')); -y = transpose(ncread('../Data/SquareSheetConstrained.nc','y')); -vx = transpose(ncread('../Data/SquareSheetConstrained.nc','vx')); -vy = transpose(ncread('../Data/SquareSheetConstrained.nc','vy')); -index = transpose(ncread('../Data/SquareSheetConstrained.nc','index')); +x = archread('../Data/SquareSheetConstrained.nc','x'); +y = archread('../Data/SquareSheetConstrained.nc','y'); +vx = archread('../Data/SquareSheetConstrained.nc','vx'); +vy = archread('../Data/SquareSheetConstrained.nc','vy'); +index = archread('../Data/SquareSheetConstrained.nc','index'); md.initialization.vx=InterpFromMeshToMesh2d(index,x,y,vx,md.mesh.x,md.mesh.y); Index: /issm/trunk-jpl/test/Par/GiaBenchmarksCD.par =================================================================== --- /issm/trunk-jpl/test/Par/GiaBenchmarksCD.par (revision 21165) +++ /issm/trunk-jpl/test/Par/GiaBenchmarksCD.par (revision 21166) @@ -26,9 +26,9 @@ %Initial velocity -x = transpose(ncread('../Data/SquareSheetConstrained.nc','x')); -y = transpose(ncread('../Data/SquareSheetConstrained.nc','y')); -vx = transpose(ncread('../Data/SquareSheetConstrained.nc','vx')); -vy = transpose(ncread('../Data/SquareSheetConstrained.nc','vy')); -index = transpose(ncread('../Data/SquareSheetConstrained.nc','index')); +x = archread('../Data/SquareSheetConstrained.nc','x'); +y = archread('../Data/SquareSheetConstrained.nc','y'); +vx = archread('../Data/SquareSheetConstrained.nc','vx'); +vy = archread('../Data/SquareSheetConstrained.nc','vy'); +index = archread('../Data/SquareSheetConstrained.nc','index'); md.initialization.vx=InterpFromMeshToMesh2d(index,x,y,vx,md.mesh.x,md.mesh.y); Index: /issm/trunk-jpl/test/Par/ISMIPE.par =================================================================== --- /issm/trunk-jpl/test/Par/ISMIPE.par (revision 21165) +++ /issm/trunk-jpl/test/Par/ISMIPE.par (revision 21166) @@ -2,5 +2,5 @@ disp(' creating thickness'); -data=transpose(ncread('../Data/ISMIPE.nc','data')); +data=archread('../Data/ISMIPE.nc','data'); md.geometry.surface=zeros(md.mesh.numberofvertices,1); md.geometry.base=zeros(md.mesh.numberofvertices,1); Index: /issm/trunk-jpl/test/Par/Pig.par =================================================================== --- /issm/trunk-jpl/test/Par/Pig.par (revision 21165) +++ /issm/trunk-jpl/test/Par/Pig.par (revision 21166) @@ -2,11 +2,11 @@ %Geometry and observation -x = transpose(ncread('../Data/Pig.nc','x')); -y = transpose(ncread('../Data/Pig.nc','y')); -vx_obs = transpose(ncread('../Data/Pig.nc','vx_obs')); -vy_obs = transpose(ncread('../Data/Pig.nc','vy_obs')); -index = transpose(ncread('../Data/Pig.nc','index')); -surface = transpose(ncread('../Data/Pig.nc','surface')); -thickness = transpose(ncread('../Data/Pig.nc','thickness')); +x = archread('../Data/Pig.nc','x'); +y = archread('../Data/Pig.nc','y'); +vx_obs = archread('../Data/Pig.nc','vx_obs'); +vy_obs = archread('../Data/Pig.nc','vy_obs'); +index = archread('../Data/Pig.nc','index'); +surface = archread('../Data/Pig.nc','surface'); +thickness = archread('../Data/Pig.nc','thickness'); md.inversion.vx_obs =InterpFromMeshToMesh2d(index,x,y,vx_obs,md.mesh.x,md.mesh.y); md.inversion.vy_obs =InterpFromMeshToMesh2d(index,x,y,vy_obs,md.mesh.x,md.mesh.y); Index: /issm/trunk-jpl/test/Par/SquareSheetConstrained.par =================================================================== --- /issm/trunk-jpl/test/Par/SquareSheetConstrained.par (revision 21165) +++ /issm/trunk-jpl/test/Par/SquareSheetConstrained.par (revision 21166) @@ -13,9 +13,9 @@ %Initial velocity -x = transpose(ncread('../Data/SquareSheetConstrained.nc','x')); -y = transpose(ncread('../Data/SquareSheetConstrained.nc','y')); -vx = transpose(ncread('../Data/SquareSheetConstrained.nc','vx')); -vy = transpose(ncread('../Data/SquareSheetConstrained.nc','vy')); -index = transpose(ncread('../Data/SquareSheetConstrained.nc','index')); +x = archread('../Data/SquareSheetConstrained.nc','x'); +y = archread('../Data/SquareSheetConstrained.nc','y'); +vx = archread('../Data/SquareSheetConstrained.nc','vx'); +vy = archread('../Data/SquareSheetConstrained.nc','vy'); +index = archread('../Data/SquareSheetConstrained.nc','index'); md.initialization.vx=InterpFromMeshToMesh2d(index,x,y,vx,md.mesh.x,md.mesh.y); Index: /issm/trunk-jpl/test/Par/SquareSheetShelf.par =================================================================== --- /issm/trunk-jpl/test/Par/SquareSheetShelf.par (revision 21165) +++ /issm/trunk-jpl/test/Par/SquareSheetShelf.par (revision 21166) @@ -16,9 +16,9 @@ %Initial velocity -x = transpose(ncread('../Data/SquareSheetShelf.nc','x')); -y = transpose(ncread('../Data/SquareSheetShelf.nc','y')); -vx = transpose(ncread('../Data/SquareSheetShelf.nc','vx')); -vy = transpose(ncread('../Data/SquareSheetShelf.nc','vy')); -index = transpose(ncread('../Data/SquareSheetShelf.nc','index')); +x = archread('../Data/SquareSheetShelf.nc','x'); +y = archread('../Data/SquareSheetShelf.nc','y'); +vx = archread('../Data/SquareSheetShelf.nc','vx'); +vy = archread('../Data/SquareSheetShelf.nc','vy'); +index = archread('../Data/SquareSheetShelf.nc','index'); md.initialization.vx=InterpFromMeshToMesh2d(index,x,y,vx,md.mesh.x,md.mesh.y); md.initialization.vy=InterpFromMeshToMesh2d(index,x,y,vy,md.mesh.x,md.mesh.y); Index: /issm/trunk-jpl/test/Par/SquareShelf.par =================================================================== --- /issm/trunk-jpl/test/Par/SquareShelf.par (revision 21165) +++ /issm/trunk-jpl/test/Par/SquareShelf.par (revision 21166) @@ -13,9 +13,9 @@ %Initial velocity and pressure -x = transpose(ncread('../Data/SquareShelf.nc','x')); -y = transpose(ncread('../Data/SquareShelf.nc','y')); -vx = transpose(ncread('../Data/SquareShelf.nc','vx')); -vy = transpose(ncread('../Data/SquareShelf.nc','vy')); -index = transpose(ncread('../Data/SquareShelf.nc','index')); +x = archread('../Data/SquareShelf.nc','x'); +y = archread('../Data/SquareShelf.nc','y'); +vx = archread('../Data/SquareShelf.nc','vx'); +vy = archread('../Data/SquareShelf.nc','vy'); +index = archread('../Data/SquareShelf.nc','index'); md.initialization.vx=InterpFromMeshToMesh2d(index,x,y,vx,md.mesh.x,md.mesh.y); md.initialization.vy=InterpFromMeshToMesh2d(index,x,y,vy,md.mesh.x,md.mesh.y); Index: /issm/trunk-jpl/test/Par/SquareShelf2.par =================================================================== --- /issm/trunk-jpl/test/Par/SquareShelf2.par (revision 21165) +++ /issm/trunk-jpl/test/Par/SquareShelf2.par (revision 21166) @@ -13,9 +13,9 @@ %Initial velocity and pressure -x = transpose(ncread('../Data/SquareShelf.nc','x')); -y = transpose(ncread('../Data/SquareShelf.nc','y')); -vx = transpose(ncread('../Data/SquareShelf.nc','vx')); -vy = transpose(ncread('../Data/SquareShelf.nc','vy')); -index = transpose(ncread('../Data/SquareShelf.nc','index')); +x = archread('../Data/SquareShelf.nc','x'); +y = archread('../Data/SquareShelf.nc','y'); +vx = archread('../Data/SquareShelf.nc','vx'); +vy = archread('../Data/SquareShelf.nc','vy'); +index = archread('../Data/SquareShelf.nc','index'); md.initialization.vx=InterpFromMeshToMesh2d(index,x,y,vx,md.mesh.x,md.mesh.y); md.initialization.vy=InterpFromMeshToMesh2d(index,x,y,vy,md.mesh.x,md.mesh.y); Index: /issm/trunk-jpl/test/Par/SquareShelfConstrained.par =================================================================== --- /issm/trunk-jpl/test/Par/SquareShelfConstrained.par (revision 21165) +++ /issm/trunk-jpl/test/Par/SquareShelfConstrained.par (revision 21166) @@ -14,9 +14,9 @@ %Initial velocity -x = transpose(ncread([issmdir() '/test/Data/SquareShelfConstrained.nc'],'x')); -y = transpose(ncread([issmdir() '/test/Data/SquareShelfConstrained.nc'],'y')); -vx = transpose(ncread([issmdir() '/test/Data/SquareShelfConstrained.nc'],'vx')); -vy = transpose(ncread([issmdir() '/test/Data/SquareShelfConstrained.nc'],'vy')); -index = transpose(ncread([issmdir() '/test/Data/SquareShelfConstrained.nc'],'index')); +x = archread([issmdir() '/test/Data/SquareShelfConstrained.nc'],'x'); +y = archread([issmdir() '/test/Data/SquareShelfConstrained.nc'],'y'); +vx = archread([issmdir() '/test/Data/SquareShelfConstrained.nc'],'vx'); +vy = archread([issmdir() '/test/Data/SquareShelfConstrained.nc'],'vy'); +index = archread([issmdir() '/test/Data/SquareShelfConstrained.nc'],'index'); md.initialization.vx=InterpFromMeshToMesh2d(index,x,y,vx,md.mesh.x,md.mesh.y); md.initialization.vy=InterpFromMeshToMesh2d(index,x,y,vy,md.mesh.x,md.mesh.y);