Index: /issm/trunk-jpl/src/m/contrib/defleurian/netCDF/export_netCDF.m =================================================================== --- /issm/trunk-jpl/src/m/contrib/defleurian/netCDF/export_netCDF.m (revision 26760) +++ /issm/trunk-jpl/src/m/contrib/defleurian/netCDF/export_netCDF.m (revision 26761) @@ -1,19 +1,20 @@ function export_netCDF(md,filename) - -%Now going on Real treatment +%verbosity of the code, 0 is no messages, 5 is chatty + verbose = 5; if exist(filename), - disp(sprintf('File %s allready exist', filename)); - prompt = 'Give a new name or "delete" to replace: '; - newname = input(prompt,'s'); - if strcmp(newname,'delete') - delete(filename) - else - disp(sprintf('New file name is %s ', newname)); - filename=newname - end - end + delete(filename) + % disp(sprintf('File %s allready exist', filename)); + % prompt = 'Give a new name or "delete" to replace: '; + % newname = input(prompt,'s'); + % if strcmp(newname,'delete') + % delete(filename) + % else + % disp(sprintf('New file name is %s ', newname)); + % filename=newname + % end + end %open file and write description mode = netcdf.getConstant('NC_NETCDF4'); - mode = bitor(mode,netcdf.getConstant('NC_NOCLOBBER'));%NOCLOBBER to avoid overwrite + mode = bitor(mode,netcdf.getConstant('NC_NOCLOBBER')); %NOCLOBBER to avoid overwrite ncid = netcdf.create(filename,mode); netcdf.putAtt(ncid,netcdf.getConstant('NC_GLOBAL'),'description',['Results for run ' md.miscellaneous.name]); @@ -21,138 +22,286 @@ %gather geometry and timestepping as dimensions - resfields=fieldnames(md.results); - Duration=size(eval(['md.results. ' resfields{1} ]),2); + if isempty(fieldnames(md.results)), + %results as no field so no time is present + Duration = 0; + else + resfields = fieldnames(md.results); + Duration = size(eval(['md.results. ' resfields{1} ]),2); + end if Duration>0, - StepNum=Duration; + StepNum = Duration; else StepNum=1; - end - - dimlist=[2,md.mesh.numberofelements,md.mesh.numberofvertices,size(md.mesh.elements,2)]; - - %define netcdf dimensions - DimSize(1).index=netcdf.defDim(ncid,'Time',StepNum); + end + DimSize(1).index=netcdf.defDim(ncid,'Time',StepNum); %time is the first dimension [DimSize(1).name,DimSize(1).value]=netcdf.inqDim(ncid,DimSize(1).index); DimValue(1)=DimSize(1).value; - DimSize(2).index=netcdf.defDim(ncid,'UnLim',netcdf.getConstant('NC_UNLIMITED')); + DimSize(2).index=netcdf.defDim(ncid,'UnLim',netcdf.getConstant('NC_UNLIMITED')); % we add an unlimited dimension if needed [DimSize(2).name,DimSize(2).value]=netcdf.inqDim(ncid,DimSize(2).index); DimValue(2)=DimSize(2).value; + % adding mesh related dimensions + dimlist=[2,40,md.mesh.numberofelements,md.mesh.numberofvertices,size(md.mesh.elements,2)]; + dimnames=["DictDummy" "StringLength" "EltNum" "VertNum" "VertPerElt"]; + if verbose > 0, + disp('===Creating dimensions ==='); + end + %define netcdf dimensions for i=1:5 + % do not add the dimension if it exists already if sum(dimlist(i) == DimValue) == 0 - DimSize(i+2).index=netcdf.defDim(ncid,['Dimension' num2str(i+2)],dimlist(i)); + DimSize(i+2).index=netcdf.defDim(ncid,dimnames(i),dimlist(i)); [DimSize(i+2).name,DimSize(i+2).value]=netcdf.inqDim(ncid,DimSize(i+2).index); DimValue(i+2)=DimSize(i+2).value; end end - - typelist=[{'numeric'} {'logical'} {'string'} {'char'} {'cell'}]; - + issmclasses = fieldnames(md)'; + typelist={'half', 'single','double','int8','int16'... + ,'int32','int64','uint8','uint16','uint32'... + ,'uint64','logical','char','string'}; %all malab types that are 0D + + for cl=1:length(issmclasses), + subclasses=fieldnames(md.(issmclasses{cl}))'; + for sc=1:length(subclasses), + if sum(strcmp(class(md.(issmclasses{cl}).(subclasses{sc})), typelist)) == 0, + issmclasses = [issmclasses class(md.(issmclasses{cl}).(subclasses{sc}))]; + end + end + end %get all model classes and create respective groups groups=fieldnames(md); + if verbose > 0, + disp('===Creating and populating groups==='); + end for i=1:length(groups), - disp(sprintf('group name in tree %s ',groups{i})); + if verbose >1, + disp(sprintf('===Now treating %s===',groups{i})); + end + if strcmp(groups{i}, 'qmu'), + disp('qmu is skipped until it is more stable'); + continue + end groupID=netcdf.defGrp(ncid,groups{i}); %In each group gather the fields of the class - groupfields=fieldnames(md.(groups{i})); - for j=1:length(groupfields), - Var=md.(groups{i}).(groupfields{j}); + fields=fieldnames(md.(groups{i})); + if isempty(fields), + disp(sprintf("WARNING: md.%s as no fields, we skip it.",groups{i})) + continue + end + %looping on fields in each group + for j=1:length(fields), + Var=md.(groups{i}).(fields{j}); + %treatment for lists if isa(Var,'cell') - Stdlist=false; + Stdlist=false; %first assume it is not a standard list if length(Var) == 0 - Stdlist=true; + Stdlist=true; %It is empty and so standard (for us) else for k=1:length(typelist) if isa(Var{1},typelist{k}) - Stdlist=true; + Stdlist=true; %if the list is of a known type (to matlab) if not it is probably some exotic ISSM stuff end end end - - netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),'classtype',class(md.(groups{i}))); - if(Stdlist) - disp(sprintf('=====Field name in tree %s ',groupfields{j})); - [DimSize,DimValue]=DefCreateVar(ncid,Var,groupID,groupfields{j},DimSize,DimValue); - else - listsize=length(Var); - subgroupID=netcdf.defGrp(groupID,groupfields{j}); - netcdf.putAtt(subgroupID,netcdf.getConstant('NC_GLOBAL'),'classtype',class(Var)); - for l=1:listsize - if isprop(Var{l},'name') - lname=Var{l}.name; - elseif isprop(Var{l},'step') - lname=Var{l}.step - else - lname=[class(Var{l}) int2str(l)]; - end - listgroupID=netcdf.defGrp(subgroupID,lname); - netcdf.putAtt(listgroupID,netcdf.getConstant('NC_GLOBAL'),'classtype',class(Var{l})); - subfields=fieldnames(Var{l}); - for m=1:length(subfields) - if ~strcmp(subfields{m},'outlog') - [DimSize,DimValue]=DefCreateVar(ncid,Var{l}.(subfields{m}),listgroupID,subfields{m},DimSize,DimValue); + %print the issm class as a classtype attribute + klass = class(md.(groups{i})); + klasstring = strcat(klass, '.',klass); + netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),'classtype',klasstring); + if(Stdlist) % this is a standard or empty list just proceed + if verbose > 4, + disp(sprintf("=££=creating var for %s.%s with classtype : %s",groups{i}, fields{j}, klasstring)) + end + if ~isempty(Var) && isa(Var{1}, 'char'), % we have a char array, pad it to a given length + Var=char(Var)'; + end + [DimSize,DimValue,varid]=CreateVar(ncid,Var,groupID,fields{j},DimSize,DimValue); + if ~isempty(varid), + FillVar(Var,groupID,varid); + end + + else % this is a list of fields, specific treatment needed (perhaps) + if verbose > 4, + disp(sprintf("=??=we have a list of fields for %s.%s with classtype : %s",groups{i}, fields{j}, klasstring)); + end + if strcmp(groups{i}, 'outputdefinition'), + listsize=length(Var); + for k=1:listsize, + subgroupname=md.(groups{i}).(fields{j}){k}.definitionstring; + subgroupID=netcdf.defGrp(groupID,subgroupname); + klass=class(md.(groups{i}).(fields{j}){k}); + klasstring = strcat(klass, '.',klass); + netcdf.putAtt(subgroupID,netcdf.getConstant('NC_GLOBAL'),'classtype',klasstring); + subfields=fieldnames(md.(groups{i}).(fields{j}){k}); + for l=1:length(subfields) + if verbose > 4, + disp(sprintf("=--=creating var for %s.%s[%i].%s",groups{i}, fields{j}, k, subfields{l})); + end + Var = md.(groups{i}).(fields{j}){k}.(subfields{l}); + [DimSize,DimValue,varid]=CreateVar(ncid,Var,subgroupID,subfields{l},DimSize,DimValue); + if ~isempty(varid), + FillVar(Var,subgroupID,varid); + end end end - end - end - elseif isa(Var,'struct') && ~strcmp(groupfields{j},'bamg') - classtype=class(md.(groups{i})); - if strcmp(classtype,'struct') - classtype=groups{i}; - end - netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),'classtype',classtype); - if length(Var)>1 - listsize=length(Var); - subgroupID=netcdf.defGrp(groupID,groupfields{j}); - classtype=class(Var); - if strcmp(classtype,'struct') - classtype=groups{i}; - end - netcdf.putAtt(subgroupID,netcdf.getConstant('NC_GLOBAL'),'classtype',classtype); - for l=1:listsize - if isfield(Var(l),'step') - lname=[int2str(Var(l).step)]; - else - lname=[class(Var(l)) int2str(l)]; + else + disp(sprintf("WARNING: unknown treatment for md.%s",groups{i})); + end + end + elseif sum(strcmp(class(Var), typelist))==1, %this is a standard matlab class with no subgrouping + if verbose > 4, + disp(sprintf("====creating var for %s.%s", groups{i}, fields{j})) + end + klass=class(md.(groups{i})); + klasstring = strcat(klass, '.',klass); + netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),'classtype',klasstring); + if sum(numel(Var) == size(Var)) == 0, %this is a 2D array or more (and not a vector with dimension 2 = 1) + Var = Var'; + end + + [DimSize,DimValue,varid]=CreateVar(ncid,Var,groupID,fields{j},DimSize,DimValue); + if ~isempty(varid), + FillVar(Var,groupID,varid); + end + + + elseif isa(Var,'struct') % structures need special treatment + if strcmp(groups{i}, 'results'), + klasstring='results.results'; + netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),'classtype',klasstring); + subgroupname=fields{j}; + subgroupID=netcdf.defGrp(groupID,subgroupname); + klasstring = 'results.solutionstep'; + netcdf.putAtt(subgroupID,netcdf.getConstant('NC_GLOBAL'),'classtype',klasstring); + subfields=fieldnames(md.(groups{i}).(fields{j})); + if isempty(subfields), + disp(sprintf("WARNING: md.%s.%s as no subfields, we skip it.",groups{i}, fields{j})); + continue + end + for k=1:length(subfields), + if verbose > 4, + disp(sprintf("=@@=creating var for %s.%s.%s",groups{i}, fields{j}, subfields{k})); end - classtype=class(Var(l)); - if strcmp(classtype,'struct') - classtype=groups{i}; + Var = md.(groups{i}).(fields{j}).(subfields{k}); + [DimSize,DimValue,varid]=CreateVar(ncid,Var,subgroupID,subfields{k},DimSize,DimValue); + if ~isempty(varid), + FillVar(Var,subgroupID,varid); end - listgroupID=netcdf.defGrp(subgroupID,lname); - netcdf.putAtt(listgroupID,netcdf.getConstant('NC_GLOBAL'),'classtype',classtype); - subfields=fieldnames(Var(l)); - for m=1:length(subfields) - if ~strcmp(subfields{m},'outlog') - [DimSize,DimValue]=DefCreateVar(ncid,Var(l).(subfields{m}),listgroupID,subfields{m},DimSize,DimValue); + end + elseif isempty(fieldnames(md.(groups{i}).(fields{j}))) % this is an empty struct, jus treat it as normal + klass=class(md.(groups{i})); + klasstring = strcat(klass, '.',klass); + netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),'classtype',klasstring); + if verbose > 4, + disp(sprintf("=[]=creating var for %s.%s",groups{i}, fields{j})); + end + + [DimSize,DimValue,varid]=CreateVar(ncid,Var,groupID,fields{j},DimSize,DimValue); + if ~isempty(varid), + FillVar(Var,groupID,varid); + end + + else + disp(sprintf("WARNING, md.%s.%s is not treated as it does not fall in one of the existing cases with class '%s'.",groups{i}, fields{j}, class(md.(groups{i}).(fields{j})))) + end + % elseif strcmp(fields{j},'bamg') + % klass=class(md.(groups{i})); + % % if strcmp(klass,'struct') + % % klass=groups{i}; + % % end + % netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),'classtype',klass); + % if length(Var)>1 + % listsize=length(Var); + % subgroupID=netcdf.defGrp(groupID,fields{j}); + % klass=class(md.(groups{i}).(fields{j})); + % netcdf.putAtt(subgroupID,netcdf.getConstant('NC_GLOBAL'),'classtype',klass); + % for k=1:listsize + % if isfield(Var(k),'step') + % lname=[int2str(Var(k).step)]; + % else + % lname=[class(Var(k)) int2str(k)]; + % end + + % listgroupID=netcdf.defGrp(subgroupID,lname); + % klass=class(md.(groups{i}).(fields{j}){k}); + % netcdf.putAtt(listgroupID,netcdf.getConstant('NC_GLOBAL'),'classtype',klass); + % subfields=fieldnames(Var(k)); + % for l=1:length(subfields) + % if ~strcmp(subfields{l},'outlog') + % if verbose > 4, + % disp(sprintf("=[]=creating var for %s.%s[%i].%s", groups{i}, fields{j}, k, subfields{l})) + % end + % Var=md.(groups{i}).(fields{j}){k}.(subfields{l}); + % [DimSize,DimValue,varid]=CreateVar(ncid,Var,listgroupID,subfields{l},DimSize,DimValue); + % if ~isempty(varid), + % FillVar(Var,listgroupID,varid); + % end + % end + % end + % end + % else + % subgroupID=netcdf.defGrp(groupID,fields{j}); + % klass=class(md.(groups{i})); + % netcdf.putAtt(subgroupID,netcdf.getConstant('NC_GLOBAL'),'classtype',klass); + % subfields=fieldnames(Var); + % for k=1:length(subfields) + % if ~strcmp(subfields{k},'outlog') + % if verbose > 4, + % disp(sprintf("====creating var for %s.%s.%s", groups{i}, fields{j}, subfields{k})) + % end + % Var=md.(groups{i}).(fields{j}).(subfields{k}); + % [DimSize,DimValue,varid]=CreateVar(ncid,Var,subgroupID,subfields{k},DimSize,DimValue); + % if ~isempty(varid), + % FillVar(Var,subgroupID,varid); + % end + % end + % end + % end + elseif sum(strcmp(class(Var), issmclasses)) == 1, % that is an issm class + if strcmp(class(Var), 'solution'), + if verbose > 4, + disp(sprintf("=$$=creating var for %s.%s",groups{i}, fields{j})) + disp("NEED treatment") + end + elseif strcmp(class(Var), 'dict'), %we have potential for a dict in py not to sure what it translates to here. + if verbose > 4, + disp(sprintf("=WW=creating var for %s.%s",groups{i}, fields{j})) + disp("NEED Treatment") + end + + else + klass=class(md.(groups{i})); + klasstring = strcat(klass, '.',klass); + netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),'classtype',klasstring); + subgroupID=netcdf.defGrp(groupID,fields{j}); + klass=class(md.(groups{i}).(fields{j})); + klasstring = strcat(klass, '.',klass); + netcdf.putAtt(subgroupID,netcdf.getConstant('NC_GLOBAL'),'classtype',klasstring); + subfields=fieldnames(Var); + for k=1:length(subfields), + if sum(strcmp(subfields{k},["outlog" "errlog"])) == 0, + if verbose > 4, + disp(sprintf("+==+creating var for %s.%s.%s",groups{i}, fields{j}, subfields{k})) + end + Var=md.(groups{i}).(fields{j}).(subfields{k}); + [DimSize,DimValue,varid]=CreateVar(ncid,Var,subgroupID,subfields{k},DimSize,DimValue); + if ~isempty(varid), + FillVar(Var,subgroupID,varid); end end end - else - subgroupID=netcdf.defGrp(groupID,groupfields{j}); - classtype=class(Var); - if strcmp(classtype,'struct') - classtype=groups{i}; - end - netcdf.putAtt(subgroupID,netcdf.getConstant('NC_GLOBAL'),'classtype',classtype); - subfields=fieldnames(Var); - for m=1:length(subfields) - if ~strcmp(subfields{m},'outlog') - [DimSize,DimValue]=DefCreateVar(ncid,Var.(subfields{m}),subgroupID,subfields{m},DimSize,DimValue); - end - end + end else - netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),'classtype',class(md.(groups{i}))); - [DimSize,DimValue]=DefCreateVar(ncid,Var,groupID,groupfields{j},DimSize,DimValue); + disp(sprintf("WARNING, md.%s.%s is not treated as it does not fall in one of the existing cases with class '%s'.",groups{i}, fields{j}, class(Var))) end end - end - netcdf.close(ncid); + end + netcdf.close(ncid); end -function [DimSize,DimValue]=DefCreateVar(ncid,Var,groupID,field,DimSize,DimValue,last,md,midfield) - varclass=class(Var); +function [DimSize,DimValue,varid]=CreateVar(ncid,Var,groupID,field,DimSize,DimValue) +% Grab dimensions varsize=size(Var); varlength=length(Var); + % treating scalar string or bool as atribute if isa(Var,'logical'), if Var, @@ -160,11 +309,58 @@ else, LogicString='False'; - end + end netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),field,LogicString); + varid=[]; + elseif isa(Var,'char'), - netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),field,Var); + if strcmp(field,'name'), % it looks like netCDF does not like attributes that are called "name" + field = 'varname'; + end + if size(Var,1) <= 1 %that is a single string or empty + netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),field,Var); + varid=[]; + else % that is a character array + [dims,DimSize,DimValue]=GetDims(ncid,Var,DimSize,DimValue); + varid = netcdf.defVar(groupID,field,'NC_CHAR',dims); + end + elseif isa(Var,'double'), %dealing with arrays + if all(mod(Var, 1) == 0, 'all') %those are actually integers, + [dims,DimSize,DimValue]=GetDims(ncid,Var,DimSize,DimValue); + varid = netcdf.defVar(groupID,field,'NC_INT64',dims); + else + [dims,DimSize,DimValue]=GetDims(ncid,Var,DimSize,DimValue); + varid = netcdf.defVar(groupID,field,'NC_DOUBLE',dims); + end + elseif isa(Var,'cell'), [dims,DimSize,DimValue]=GetDims(ncid,Var,DimSize,DimValue); - varid = netcdf.defVar(groupID,field,'NC_DOUBLE',dims); + varid = netcdf.defVar(groupID,field,'NC_CHAR',dims); + + elseif isa(Var,'struct'), + if isempty(fieldnames(Var)), + netcdf.putAtt(groupID,netcdf.getConstant('NC_GLOBAL'),field,'emptystruct'); + varid=[]; + else + %Start by getting the structure fields and size + [dims,DimSize,DimValue]=GetDims(ncid,Var,DimSize,DimValue); + varid = netcdf.defVar(groupID,field,'NC_CHAR',dims); + end + else + disp(sprintf('no support for class %s of field %s',class(Var),field)); + varid=[]; + end + return +end + + +function FillVar(Var,groupID,varid) +% Grab dimensions + varsize=size(Var); + varlength=length(Var); + % treating scalar string or bool as atribute + if isa(Var,'double'), %dealing with arrays + if all(mod(Var, 1) == 0, 'all') %those are actually integers, + Var = int64(Var); + end if length(Var)==0, netcdf.putVar(groupID,varid,NaN); @@ -172,24 +368,26 @@ netcdf.putVar(groupID,varid,Var); end + elseif isa(Var,'char'), % at this point this should be a character array + netcdf.putVar(groupID,varid,Var); elseif isa(Var,'cell'), - [dims,DimSize,DimValue]=GetDims(ncid,Var,DimSize,DimValue); - %dirty hack to be able to pass strings - varid = netcdf.defVar(groupID,field,'NC_CHAR',dims); - if length(Var)==0, - netcdf.putVar(groupID,varid,0,9,'emptycell') - else - for i=1:length(Var), - if length(Var)>1, - endpoint=[1,min(length(Var{i}),40)]; - startpoint=[1 0]; - else - endpoint=min(length(Var{i}),40); - startpoint=0; - end - if length(Var{i})>40, - netcdf.putVar(groupID,varid,startpoint,extent,Var{i}(1:40)) - disp(sprintf('some variable have been truncated')); - else - netcdf.putVar(groupID,varid,startpoint,endpoint,Var{i}) + if ~isempty(Var), + if length(Var)==0, + netcdf.putVar(groupID,varid,0,9,'emptycell'); + else + for i=1:length(Var), + if length(Var)>1, + count=[min(length(Var{i}),40), 1]; + startpoint=[0 i-1]; + else + count=min(length(Var{i}),40); + startpoint=0; + end + + if length(Var{i})>40, + netcdf.putVar(groupID,varid,startpoint,count,Var{i}(1:40)); + disp(sprintf('some variable have been truncated')); + else + netcdf.putVar(groupID,varid,startpoint,count,Var{i}); + end end end @@ -197,43 +395,37 @@ elseif isa(Var,'struct'), %Start by getting the structure fields and size - locfields=fieldnames(Var); - [dims,DimSize,DimValue]=GetDims(ncid,Var,DimSize,DimValue); - varid = netcdf.defVar(groupID,field,'NC_CHAR',dims); - if length(locfields)==0, - netcdf.putVar(groupID,varid,[0,0],[1,11],'emptystruct') - else - for i=1:length(locfields), - for j=1:2, - if j==1, - CharVar=locfields{i}; - if length(CharVar)==0 - CharVar='emptystruct'; - end - startpoint=[i-1,0,0]; + locfields=fieldnames(Var) + for i=1:length(locfields), + for j=1:2, + if j==1, + CharVar=locfields{i}; + if length(CharVar)==0 + CharVar='emptystruct'; + end + startpoint=[i-1,0,0]; + else + if isa(Var.(locfields{i}),'char'), + CharVar=Var.(locfields{i}); else - if isa(Var.(locfields{i}),'char'), - CharVar=Var.(locfields{i}); - else - CharVar=num2str(Var.(locfields{i})); - end - if length(CharVar)==0 - CharVar='emptystruct'; - end - startpoint=[i-1,1,0]; - end - - extent=[1,1,min(length(CharVar),40)]; - if length(CharVar)>40, - netcdf.putVar(groupID,varid,startpoint,extent,CharVar(1:40)) - disp(sprintf('some variable have been truncated')); - else - netcdf.putVar(groupID,varid,startpoint,extent,CharVar) - end + CharVar=num2str(Var.(locfields{i})); + end + if length(CharVar)==0 + CharVar='emptystruct'; + end + startpoint=[i-1,1,0]; + end + + extent=[1,1,min(length(CharVar),40)]; + if length(CharVar)>40, + netcdf.putVar(groupID,varid,startpoint,extent,CharVar(1:40)); + disp(sprintf('some variable have been truncated')); + else + netcdf.putVar(groupID,varid,startpoint,extent,CharVar); end end end else - disp(sprintf('no support for class %s of field %s',varclass,field)); - end + disp(sprintf('no support for class %s',class(Var))); + end return end @@ -241,44 +433,42 @@ function [dims,DimSize,DimValue]=GetDims(ncid,Var,DimSize,DimValue) dims=[]; - if isa(Var,'cell'), - varsize=size(Var'); - elseif isa(Var,'struct') + % if isa(Var,'cell'), + % varsize=size(Var'); + % else + if isa(Var,'struct') varsize=length(fieldnames(Var)); else varsize=size(Var); end - dim=sum(varsize>1); - if dim>0 - for i=1:dim - indsize=find(varsize(i)==DimValue); - if length(indsize)>0 - dims=[dims DimSize(indsize).index]; - else - indsize=length(DimSize)+1; - DimSize(indsize).index=netcdf.defDim(ncid,['Dimension' num2str(indsize)],varsize(i)); - [DimSize(indsize).name,DimSize(indsize).value]=netcdf.inqDim(ncid,DimSize(indsize).index); - DimValue(indsize)=DimSize(indsize).value; - dims=[dims DimSize(indsize).index]; + + % dim=sum(varsize>1); + dim=ndims(Var); + if dim>0, + for i=1:dim, + if size(Var, i) >1 + indsize=find(varsize(i)==DimValue); + if length(indsize)>0 + dims=[dims DimSize(indsize).index]; + else + indsize=length(DimSize)+1; + DimSize(indsize).index=netcdf.defDim(ncid,['DimNum' num2str(indsize)],varsize(i)); + [DimSize(indsize).name,DimSize(indsize).value]=netcdf.inqDim(ncid,DimSize(indsize).index); + DimValue(indsize)=DimSize(indsize).value; + dims=[dims DimSize(indsize).index]; + end end end end - %if we have a cell variable we need to add a stringlength dimension + %if we have an empty cell variable we need to add a stringlength for the no data + % if isa(Var,'cell'), + % if isempty(dims) + % dims=[DimSize(4).index]; + % else + % dims=[DimSize(4).index dims]; + % end + % end + % struct also need an extra dimension 2, but only if non empty if isa(Var,'struct'), - if DimValue(3)~=2 - if DimValue(2)~=2 - dims=[dims DimSize(1).index]; - else - dims=[dims DimSize(2).index]; - end - else - dims=[dims DimSize(3).index]; - end - end - if isa(Var,'cell') || isa(Var,'struct'), - if DimValue(2)~=40 - dims=[dims DimSize(1).index]; - else - dims=[dims DimSize(2).index]; - end + dims=[dims DimSize(3).index DimSize(4).index]; end end Index: /issm/trunk-jpl/src/m/contrib/defleurian/netCDF/export_netCDF.py =================================================================== --- /issm/trunk-jpl/src/m/contrib/defleurian/netCDF/export_netCDF.py (revision 26760) +++ /issm/trunk-jpl/src/m/contrib/defleurian/netCDF/export_netCDF.py (revision 26761) @@ -60,5 +60,5 @@ def export_netCDF(md, filename): # {{{ #verbosity of the code, 0 is no messages, 5 is chatty - verbose = 0 + verbose = 5 if path.exists(filename): print('File {} allready exist'.format(filename)) @@ -88,6 +88,6 @@ dimindex = 2 #add mesh related dimension that we know are needed - dimlist = [2, md.mesh.numberofelements, md.mesh.numberofvertices, np.shape(md.mesh.elements)[1]] - dimnames = ['DictDummy', 'EltNum', 'VertNum', 'VertPerElt'] + dimlist = [2, 40, md.mesh.numberofelements, md.mesh.numberofvertices, np.shape(md.mesh.elements)[1]] + dimnames = ['DictDummy', 'StringLength', 'EltNum', 'VertNum', 'VertPerElt'] if verbose > 0: print('===Creating dimensions ===') @@ -121,12 +121,13 @@ # looping on fields in each group for field in fields: + Var = md.__dict__[group].__dict__[field] # Special treatment for list fields - if type(md.__dict__[group].__dict__[field]) == list: + if type(Var) == list: StdList = False - if len(md.__dict__[group].__dict__[field]) == 0: + if len(Var) == 0: StdList = True #this is an empty list else: #returns False for exotic types (typicaly results) - StdList = type(md.__dict__[group].__dict__[field][0]) in typelist + StdList = type(Var[0]) in typelist klass = type(md.__dict__[group]).__module__ + '.' + type(md.__dict__[group]).__name__ NCgroup.__setattr__('classtype', klass) @@ -134,5 +135,4 @@ if verbose > 4: print("=££=creating var for {}.{} with classtype : {}".format(group, field, klass)) - Var = md.__dict__[group].__dict__[field] Var = SqueezeVar(Var) DimDict, ncvar = CreateVar(NCData, Var, field, NCgroup, DimDict) @@ -141,6 +141,6 @@ else: # this is a list of fields, specific treatment needed (usually results or outputdefinitions) if verbose > 4: - print("list of fields happens for {}.{} with classtype : {}".format(group, field, klass)) - Listsize = len(md.__dict__[group].__dict__[field]) + print("=??=we have a list of fields for {}.{} with classtype : {}".format(group, field, klass)) + Listsize = len(Var) if group == 'results': #for results we reshape the datas following time rather than subgrouping Subgroup = NCgroup.createGroup(str(field)) @@ -172,5 +172,5 @@ StackedVar.update(Var) if verbose > 4: - print("=$$=creating var for {}.{}.{}".format(group, field, subfield)) + print("=@@=creating var for {}.{}.{}".format(group, field, subfield)) print("last index of the list is {}".format(lastindex)) StackedVar = SqueezeVar(StackedVar.finalize(int(lastindex))) @@ -289,5 +289,5 @@ FillVar(ncvar, Var) else: - print("WARNING, md.{}.{} is not treated as it does not fall in one of the existing cases.".format(str(group), str(field))) + print("WARNING, md.{}.{} is not treated as it does not fall in one of the existing cases.".format(group, field)) NCData.close() @@ -329,5 +329,6 @@ if field == 'name': # it looks like netCDF does not like attributes that are called "name" field = 'varname' - Group.__setattr__(str(field).swapcase(), str(var)) + #Group.__setattr__(str(field).swapcase(), str(var)) + Group.__setattr__(str(field), str(var)) ncvar = None # numpy array of strings @@ -337,5 +338,6 @@ if field == 'name': field = 'varname' - Group.__setattr__(str(field).swapcase(), str(var[0])) + #Group.__setattr__(str(field).swapcase(), str(var[0])) + Group.__setattr__(str(field), str(var[0])) ncvar = None else: Index: /issm/trunk-jpl/src/m/io/loadvars.py =================================================================== --- /issm/trunk-jpl/src/m/io/loadvars.py (revision 26760) +++ /issm/trunk-jpl/src/m/io/loadvars.py (revision 26761) @@ -7,5 +7,5 @@ from re import findall, split import shelve -from netCDF4 import Dataset +from netCDF4 import Dataset, chartostring import numpy as np from model import * @@ -100,5 +100,6 @@ #that is the current treatment #here we have a more NC approach with time being a dimension - listtype = split(r'\.', classtype[mod][0])[0] + listtype = split(r'\.', classtype[mod][0])[1] + print(listtype) if len(NCFile.dimensions['Time']) == 1: nvdict['md'].__dict__[classtree[mod][0]].__dict__[classtree[mod][1]] = getattr(classtype[mod][1], listtype)() @@ -137,8 +138,8 @@ Tree = nvdict['md'].__dict__[classtree[mod][0]].__dict__[classtree[mod][1]][:] #}}} - elif classtype[mod][0] == 'massfluxatgate': #this is for output definitions {{{ + elif classtype[mod][0] == 'massfluxatgate.massfluxatgate': #this is for output definitions {{{ defname = split('Output|[0-9]+', classtree[mod][1])[1] + 's' defindex = int(findall('[0-9]+', classtree[mod][1])[0]) - nvdict['md'].__dict__[classtree[mod][0]].__dict__[defname].append(getattr(classtype[mod][1], classtype[mod][0])()) + nvdict['md'].__dict__[classtree[mod][0]].__dict__[defname].append(getattr(classtype[mod][1], 'massfluxatgate')()) Tree = nvdict['md'].__dict__[classtree[mod][0]].__dict__[defname][defindex - 1] #}}} @@ -223,5 +224,5 @@ else: if vardim == 0: #that is a scalar - if str(varval[0]) == '': #no value + if str(varval[0]) == '' or str(varval[0]) == '--': #no value Tree.__dict__[str(var)] = [] elif varval[0] == 'True': #treatin bool @@ -233,4 +234,6 @@ elif vardim == 1: #that is a vector + if debug: + print(" for variable {} type is {}".format(str(var), varval.dtype)) if varval.dtype == str: if varval.shape[0] == 1: @@ -240,4 +243,7 @@ else: Tree.__dict__[str(var)] = [str(vallue) for vallue in varval[:]] + elif varval.dtype == "|S1": #that is for matlab chararcter arrays + stringlist = chartostring(varval[:]) + Tree.__dict__[str(var)] = [stringlist.tolist(), ] else: try: @@ -253,4 +259,6 @@ elif vardim == 2: #dealling with dict + if debug: + print(" for variable {} type is {}".format(str(var), varval.dtype)) if varval.dtype == str: #that is for dictionaries if any(varval[:, 0] == 'toolkit'): #toolkit definition have to be first @@ -263,4 +271,8 @@ strings2 = [str(arg[1]) for arg in varval] Tree.__dict__[str(var)] = OrderedDict(list(zip(strings1, strings2))) + elif varval.dtype == "|S1": #that is for matlab chararcter arrays + stringlist = chartostring(varval[:, :]) + stringlist = [string.strip() for string in stringlist] + Tree.__dict__[str(var)] = stringlist else: if type(Tree) == list: @@ -282,7 +294,11 @@ print(" ==> treating attribute {}".format(attr)) if attr != 'classtype': #classtype is for treatment, don't get it back - attribute = str(attr).swapcase() #there is a reason for swapcase, no sure what it isanymore - if attr == 'VARNAME': + # attribute = str(attr).swapcase() #there is a reason for swapcase, no sure what it isanymore + # if attr == 'VARNAME': + # attribute = 'name' + if attr == 'varname': attribute = 'name' + else: + attribute = attr if type(Tree) == list: if debug: