/[MITgcm]/MITgcm_contrib/gael/matlab_class/gcmfaces_IO/process2nctiles.m
ViewVC logotype

Diff of /MITgcm_contrib/gael/matlab_class/gcmfaces_IO/process2nctiles.m

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.9 by gforget, Sun Aug 3 19:49:58 2014 UTC revision 1.13 by gforget, Wed Mar 16 15:02:19 2016 UTC
# Line 24  gcmfaces_global; Line 24  gcmfaces_global;
24  doClim=0;  doClim=0;
25    
26  %directory names  %directory names
27  listDirs={'STATE/','TRSP/'};  listDirs={'STATE/','TRSP/'};%BUDG?
28  filAvailDiag=[dirModel 'available_diagnostics.log'];  filAvailDiag=[dirModel 'available_diagnostics.log'];
29  filReadme=[dirModel 'README'];  filReadme=[dirModel 'README'];
30  dirOut=[dirModel 'nctiles_tmp/'];  dirOut=[dirModel 'nctiles_tmp/'];
31    %dirOut=[dirModel 'nctiles_post_tmp/'];
32  if ~isdir(dirOut); mkdir(dirOut); end;  if ~isdir(dirOut); mkdir(dirOut); end;
33    
34  %search in subdirectories  %search in subdirectories
35  subDir=[];  subDir=[];
36    diagsDir='diags/';
37    %diagsDir='diags_post/';
38    %diagsDir='diags_interp/';
39  for ff=1:length(listDirs);  for ff=1:length(listDirs);
40  tmp1=dir([dirModel 'diags/' listDirs{ff} fileModel '*']);  tmp1=dir([dirModel diagsDir listDirs{ff} fileModel '*']);
41  if ~isempty(tmp1); subDir=listDirs{ff}; end;  if ~isempty(tmp1); subDir=listDirs{ff}; end;
42  end;  end;
43    
44    if isempty(subDir);
45    tmp1=dir([dirModel diagsDir fileModel '/' fileModel '*']);
46    if ~isempty(tmp1); subDir=[fileModel '/']; end;
47    end;
48    
49  if isempty(subDir);  if isempty(subDir);
50    error(['file ' fileModel ' was not found']);    error(['file ' fileModel ' was not found']);
51  else;  else;
52    dirIn=[dirModel 'diags/' subDir];    dirIn=[dirModel diagsDir subDir];
53    nn=length(dir([dirIn fileModel '*data']));    nn=length(dir([dirIn fileModel '*data']));
54    fprintf('%s (%d files) was found in \n %s \n',fileModel,nn,dirIn);    fprintf('%s (%d files) was found in \n %s \n',fileModel,nn,dirIn);
55  end;  end;
# Line 87  tim=datenum(tim)-datenum([1992 1 0]); Line 96  tim=datenum(tim)-datenum([1992 1 0]);
96  timUnits='days since 1992-1-1 0:0:0';  timUnits='days since 1992-1-1 0:0:0';
97    
98  %get time step axis  %get time step axis
99  [listTimes]=diags_list_times({[dirModel 'diags/STATE/']},{'state_2d_set1'});  [listTimes]=diags_list_times({dirIn},{fileModel});
100    
101  %get units and long name from available_diagnostics.log  %get units and long name from available_diagnostics.log
102  [avail_diag]=read_avail_diag(filAvailDiag,nameDiag);  [avail_diag]=read_avail_diag(filAvailDiag,nameDiag);
# Line 96  timUnits='days since 1992-1-1 0:0:0'; Line 105  timUnits='days since 1992-1-1 0:0:0';
105  [rdm]=read_readme(filReadme);  [rdm]=read_readme(filReadme);
106  disp(rdm');  disp(rdm');
107    
108  %convert to MITgcm format (90x1170 array)  %set output directory/file name
109  myFile=[dirOut nameDiag];%first instance is for subdirectory name  myFile=[dirOut nameDiag];%first instance is for subdirectory name
110  if ~isdir(myFile); mkdir(myFile); end;  if ~isdir(myFile); mkdir(myFile); end;
111  myFile=[myFile '/' nameDiag];%second instance is for file name base  myFile=[myFile filesep nameDiag];%second instance is for file name base
112    
113  %get grid params  %get grid params
114  [grid_diag]=set_grid_diag(avail_diag);  [grid_diag]=set_grid_diag(avail_diag);
115    
116  %apply mask, and convert to land mask  %apply mask, and convert to land mask
117  msk=grid_diag.msk;  if ~isempty(mygrid.RAC);
118  if length(size(myDiag{1}))==3;    msk=grid_diag.msk;
119    msk=repmat(msk,[1 1 size(myDiag{1},3)]);    if length(size(myDiag{1}))==3;
120  end;      msk=repmat(msk(:,:,1),[1 1 size(myDiag{1},3)]);
121  for kk=1:size(myDiag{1},4);    else;
122    myDiag(:,:,:,kk)=myDiag(:,:,:,kk).*msk;      msk=repmat(msk,[1 1 1 size(myDiag{1},4)]);
123  end;    end;
124  clear msk;    myDiag=myDiag.*msk;
125  %    clear msk;
126  land=isnan(grid_diag.msk);    %
127      land=isnan(grid_diag.msk);
128    end;
129    
130  %set 'coord' attribute  %set 'coord' attribute
131  if avail_diag.nr~=1;  if avail_diag.nr~=1;
# Line 130  if doClim; Line 141  if doClim;
141    avail_diag.longNameDiag=[avail_diag.longNameDiag ' (climatology) '];    avail_diag.longNameDiag=[avail_diag.longNameDiag ' (climatology) '];
142  end;  end;
143    
144  %create netcdf file using write2nctiles (works only with old matlab, thus far ...)  %create netcdf file using write2nctiles
145  doCreate=1;  doCreate=1;
146  dimlist=write2nctiles(myFile,myDiag,doCreate,{'tileNo',tileNo},...  dimlist=write2nctiles(myFile,myDiag,doCreate,{'tileNo',tileNo},...
147      {'fldName',nameDiag},{'longName',avail_diag.longNameDiag},...      {'fldName',nameDiag},{'longName',avail_diag.longNameDiag},...
# Line 153  end; Line 164  end;
164  doCreate=0;  doCreate=0;
165    
166  %now add fields  %now add fields
 write2nctiles(myFile,tim,doCreate,{'tileNo',tileNo},...  
   {'fldName','tim'},{'longName','time'},...  
   {'units',timUnits},{'dimIn',dim.tim});  
 write2nctiles(myFile,listTimes,doCreate,{'tileNo',tileNo},...  
   {'fldName','step'},{'longName','final time step number'},...  
   {'units','1'},{'dimIn',dim.tim});  
167  write2nctiles(myFile,grid_diag.lon,doCreate,{'tileNo',tileNo},...  write2nctiles(myFile,grid_diag.lon,doCreate,{'tileNo',tileNo},...
168    {'fldName','lon'},{'units','degrees_east'},{'dimIn',dim.twoD});    {'fldName','lon'},{'units','degrees_east'},{'dimIn',dim.twoD});
169  write2nctiles(myFile,grid_diag.lat,doCreate,{'tileNo',tileNo},...  write2nctiles(myFile,grid_diag.lat,doCreate,{'tileNo',tileNo},...
170    {'fldName','lat'},{'units','degrees_north'},{'dimIn',dim.twoD});    {'fldName','lat'},{'units','degrees_north'},{'dimIn',dim.twoD});
 write2nctiles(myFile,grid_diag.msk,doCreate,{'tileNo',tileNo},...  
   {'fldName','land'},{'units','1'},{'longName','land mask'},{'dimIn',dim.threeD});  
 write2nctiles(myFile,grid_diag.RAC,doCreate,{'tileNo',tileNo},...  
   {'fldName','area'},{'units','m^2'},{'longName','grid cell area'},{'dimIn',dim.twoD});  
171  if isfield(grid_diag,'dep');  if isfield(grid_diag,'dep');
172      write2nctiles(myFile,grid_diag.dep,doCreate,{'tileNo',tileNo},...      write2nctiles(myFile,grid_diag.dep,doCreate,{'tileNo',tileNo},...
173        {'fldName','dep'},{'units','m'},{'dimIn',dim.dep});        {'fldName','dep'},{'units','m'},{'dimIn',dim.dep});
174    end;
175    write2nctiles(myFile,tim,doCreate,{'tileNo',tileNo},...
176      {'fldName','tim'},{'longName','time'},...
177      {'units',timUnits},{'dimIn',dim.tim});
178    if ~isempty(mygrid.RAC);
179      write2nctiles(myFile,listTimes,doCreate,{'tileNo',tileNo},...
180        {'fldName','timstep'},{'longName','final time step number'},...
181        {'units','1'},{'dimIn',dim.tim});
182      write2nctiles(myFile,grid_diag.msk,doCreate,{'tileNo',tileNo},...
183        {'fldName','land'},{'units','1'},{'longName','land mask'},{'dimIn',dim.threeD});
184      write2nctiles(myFile,grid_diag.RAC,doCreate,{'tileNo',tileNo},...
185        {'fldName','area'},{'units','m^2'},{'longName','grid cell area'},{'dimIn',dim.twoD});
186      if isfield(grid_diag,'dz');
187      write2nctiles(myFile,grid_diag.dz,doCreate,{'tileNo',tileNo},...      write2nctiles(myFile,grid_diag.dz,doCreate,{'tileNo',tileNo},...
188        {'fldName','thic'},{'units','m'},{'dimIn',dim.dep});        {'fldName','thic'},{'units','m'},{'dimIn',dim.dep});
189      end;
190  end;  end;
191    
192    clear myDiag;
193    
194  end;%for vv=1:length(listFlds);  end;%for vv=1:length(listFlds);
195    
196  function [meta]=read_meta(fileName);  function [meta]=read_meta(fileName);
Legend:
Removed from v.1.9  
changed lines
  Added in v.1.13
  ViewVC Help
Powered by ViewVC 1.1.22