| 1 | 
jmc | 
1.1 | 
function [nIt,rList,tim,vvA,listV]=read_StD(namF,sufx,listV); | 
| 2 | 
  | 
  | 
% [nIt,rList,tim,vvA,listV]=read_StD(namF,sufx,listV); | 
| 3 | 
  | 
  | 
% | 
| 4 | 
  | 
  | 
% read ASCII stat-Diags output files (after splitted by script extract_StD) | 
| 5 | 
  | 
  | 
% | 
| 6 | 
  | 
  | 
% input: listV = list of fileds to read in ; if ='all_flds' => read all fields | 
| 7 | 
  | 
  | 
%        namF  = prefix of all file names to read (after extract_StD) | 
| 8 | 
  | 
  | 
%        sufx  = suffix of all file names to read (after extract_StD) | 
| 9 | 
  | 
  | 
% output:  | 
| 10 | 
  | 
  | 
%  nIt      = number of time reccords | 
| 11 | 
  | 
  | 
%  rList    = list of region number | 
| 12 | 
  | 
  | 
%  tim(:,1) = iterations number ; tim(:,2) = time in simulation | 
| 13 | 
  | 
  | 
%  listV    = list of fields | 
| 14 | 
  | 
  | 
%  vvA      = 5 dims output array: | 
| 15 | 
  | 
  | 
%           ( kLev, time_rec, region_rec, [ave,std,min,max,vol], var_rec ) | 
| 16 | 
  | 
  | 
 | 
| 17 | 
  | 
  | 
% $Header:  $ | 
| 18 | 
  | 
  | 
% $Name:  $ | 
| 19 | 
  | 
  | 
 | 
| 20 | 
  | 
  | 
namfhd=[namF,'_head','.',sufx]; | 
| 21 | 
  | 
  | 
namfil=[namF,'_Iter','.',sufx]; | 
| 22 | 
  | 
  | 
 fprintf(['read ',sufx,' :']); | 
| 23 | 
  | 
  | 
rf=-1; if strcmp(char(listV),'all_flds'), rf=0; end | 
| 24 | 
  | 
  | 
 | 
| 25 | 
  | 
  | 
%- read from Header file:  | 
| 26 | 
  | 
  | 
%   frequency ; list of regions ; list of levels ; (+ list of Var, if rf=0) | 
| 27 | 
  | 
  | 
fid=fopen(namfhd,'r'); | 
| 28 | 
  | 
  | 
D=dir(namfhd); | 
| 29 | 
  | 
  | 
if size(D,1) == 1, | 
| 30 | 
  | 
  | 
 fprintf(' head'); | 
| 31 | 
  | 
  | 
else error(['file: ',namfhd,' not found']); return; end | 
| 32 | 
  | 
  | 
flag=1; l=0; k=3;  | 
| 33 | 
  | 
  | 
while flag > 0 | 
| 34 | 
  | 
  | 
 tline=fgetl(fid); | 
| 35 | 
  | 
  | 
 if ischar(tline), l=l+1;  | 
| 36 | 
  | 
  | 
   %disp(tline); | 
| 37 | 
  | 
  | 
    if tline(2:12) == ' frequency ',  | 
| 38 | 
  | 
  | 
     frq=sscanf(tline(17:end),'%f'); | 
| 39 | 
  | 
  | 
     k=k-1; | 
| 40 | 
  | 
  | 
    end | 
| 41 | 
  | 
  | 
    if tline(2:10) == ' Regions ', | 
| 42 | 
  | 
  | 
     rList=sscanf(tline(17:end),'%i'); | 
| 43 | 
  | 
  | 
     k=k-1; | 
| 44 | 
  | 
  | 
    end | 
| 45 | 
  | 
  | 
    if tline(2:15) == ' Nb of levels ',  | 
| 46 | 
  | 
  | 
     kList=sscanf(tline(17:end),'%i'); | 
| 47 | 
  | 
  | 
     k=k-1; | 
| 48 | 
  | 
  | 
    end | 
| 49 | 
  | 
  | 
    if rf >=0 & tline(2:10) == ' Fields  ',  | 
| 50 | 
  | 
  | 
      tmp=[tline(17:end),' ']; i1=0; | 
| 51 | 
  | 
  | 
      for i=1:length(tmp), | 
| 52 | 
  | 
  | 
        if isspace(tmp(i)),  | 
| 53 | 
  | 
  | 
          if i1 > 0 & i > i1, | 
| 54 | 
  | 
  | 
            if rf == 0, listV=cellstr(tmp(i1:i-1)); | 
| 55 | 
  | 
  | 
            else listV(rf+1)=cellstr(tmp(i1:i-1)); end | 
| 56 | 
  | 
  | 
            i1=i+1; rf=rf+1; | 
| 57 | 
  | 
  | 
          else i1=i+1;  | 
| 58 | 
  | 
  | 
          end | 
| 59 | 
  | 
  | 
        end | 
| 60 | 
  | 
  | 
      end | 
| 61 | 
  | 
  | 
    end | 
| 62 | 
  | 
  | 
    if tline(2:15) == ' end of header',  | 
| 63 | 
  | 
  | 
     flag=0; | 
| 64 | 
  | 
  | 
    end | 
| 65 | 
  | 
  | 
 else flag=-1; end | 
| 66 | 
  | 
  | 
end | 
| 67 | 
  | 
  | 
fclose(fid); | 
| 68 | 
  | 
  | 
if rf > 0, | 
| 69 | 
  | 
  | 
%- rename fields (consistent with script "extract_StD"): | 
| 70 | 
  | 
  | 
   listV=strrep(listV,'ETAN', 'Eta'); listV=strrep(listV,'ETANSQ','Et2'); | 
| 71 | 
  | 
  | 
   listV=strrep(listV,'THETA' ,'T');  listV=strrep(listV,'SALT'  ,'S'); | 
| 72 | 
  | 
  | 
   listV=strrep(listV,'UVEL'  ,'U');  listV=strrep(listV,'VVEL'  ,'V'); | 
| 73 | 
  | 
  | 
   listV=strrep(listV,'WVEL'  ,'W');  listV=strrep(listV,'PHIHYD','Phi'); | 
| 74 | 
  | 
  | 
   listV=strrep(listV,'UVELSQ','U2'); listV=strrep(listV,'VVELSQ','V2'); | 
| 75 | 
  | 
  | 
end | 
| 76 | 
  | 
  | 
if flag ~= 0 | k > 0, frq,rList,kList, end | 
| 77 | 
  | 
  | 
if flag ~= 0, error(['not normal end after reading ',int2str(l),' lines']); end | 
| 78 | 
  | 
  | 
if k > 0, error(['missing ',int2str(k),' flds (freq,regions,levels)']); end | 
| 79 | 
  | 
  | 
if flag ~= 0 | k > 0, frq=0; rList=0; kList=0; return; end | 
| 80 | 
  | 
  | 
nReg=size(rList,1); | 
| 81 | 
  | 
  | 
nbV=size(listV,2); | 
| 82 | 
  | 
  | 
%fprintf('\n kList=');fprintf(' %i',kList);fprintf('\n'); | 
| 83 | 
  | 
  | 
%if rf > 0, for l=1:nbV, fprintf([' ',char(listV(l))]); end; fprintf('\n'); end | 
| 84 | 
  | 
  | 
if ( rf > 0 & nbV < length(kList) ) | nbV > length(kList), | 
| 85 | 
  | 
  | 
 fprintf('\n Warning: nbV= %i but (only) %i in kList\n',nbV,length(kList)) | 
| 86 | 
  | 
  | 
end | 
| 87 | 
  | 
  | 
%- set Max Nb of levels: n3d | 
| 88 | 
  | 
  | 
n3d=max(kList); if n3d > 1, n3d=1+n3d ; end | 
| 89 | 
  | 
  | 
 | 
| 90 | 
  | 
  | 
%- load Iter: | 
| 91 | 
  | 
  | 
D=dir(namfil); | 
| 92 | 
  | 
  | 
if size(D,1) == 1, | 
| 93 | 
  | 
  | 
 fprintf(' & Iter'); var=load(namfil); | 
| 94 | 
  | 
  | 
else fprintf(['file: ',namfil,' not found => EXIT \n']); return; end | 
| 95 | 
  | 
  | 
nIt=size(var,1); | 
| 96 | 
  | 
  | 
 | 
| 97 | 
  | 
  | 
%- initialize output & save Iter: | 
| 98 | 
  | 
  | 
tim=zeros(nIt,2); vvA=zeros(n3d,nIt,nReg,5,nbV); | 
| 99 | 
  | 
  | 
msgA=' '; msgB=' '; | 
| 100 | 
  | 
  | 
tim(:,1)=var;  | 
| 101 | 
  | 
  | 
  | 
| 102 | 
  | 
  | 
dIt=ones(1,nIt); | 
| 103 | 
  | 
  | 
if nIt > 1, dIt(2:nIt)=tim(2:nIt,1)-tim(1:nIt-1,1); dIt(1)=dIt(2); end | 
| 104 | 
  | 
  | 
% fprintf(' (dIt: %i %i)',min(dIt(1:nIt)),max(dIt(1:nIt))); | 
| 105 | 
  | 
  | 
%- take the inverse : will be used to correct the volume | 
| 106 | 
  | 
  | 
%            (frq<0 : no volume correction if snap-shot) | 
| 107 | 
  | 
  | 
if frq > 0, dIt=1./max(1,dIt); else dIt=ones(1,nIt); end | 
| 108 | 
  | 
  | 
 | 
| 109 | 
  | 
  | 
%- build time:  | 
| 110 | 
  | 
  | 
delT=0; | 
| 111 | 
  | 
  | 
if nIt > 1, delT=tim(2,1)-tim(1,1); delT=max(0,delT); end | 
| 112 | 
  | 
  | 
if delT > 0, | 
| 113 | 
  | 
  | 
  delT=abs(frq)/delT; | 
| 114 | 
  | 
  | 
  delta=delT-round(delT); | 
| 115 | 
  | 
  | 
  delta=abs(delta*(tim(2,1)-tim(1,1))); | 
| 116 | 
  | 
  | 
  if delta <= 0.5, delT=round(delT); end | 
| 117 | 
  | 
  | 
end | 
| 118 | 
  | 
  | 
fprintf(': nReg= %i ; delta.T= %6.1f ; freq= %8.0f',nReg,delT,frq); | 
| 119 | 
  | 
  | 
if delT > 0, fprintf(' (= %8.3f it)',frq/delT); end | 
| 120 | 
  | 
  | 
fprintf('\n') | 
| 121 | 
  | 
  | 
 | 
| 122 | 
  | 
  | 
if delT > 0, | 
| 123 | 
  | 
  | 
 tim(:,2)=delT*tim(:,1); | 
| 124 | 
  | 
  | 
else | 
| 125 | 
  | 
  | 
 tim(:,2)=tim(:,1); | 
| 126 | 
  | 
  | 
end | 
| 127 | 
  | 
  | 
 fprintf(' + var:'); var=load(namfil); | 
| 128 | 
  | 
  | 
 | 
| 129 | 
  | 
  | 
%--------- | 
| 130 | 
  | 
  | 
for nv=1:nbV, | 
| 131 | 
  | 
  | 
  namV=char(listV(nv)); | 
| 132 | 
  | 
  | 
 [vv1,nk,msg]=read_StD_1v(namF,namV,sufx,nReg,nIt,dIt); | 
| 133 | 
  | 
  | 
  if nk > 0, fprintf('_%i',nk); elseif nk < 0, msgA=sprintf([msgA,msg]); | 
| 134 | 
  | 
  | 
  else msgB=sprintf([msgB,msg]) ; end | 
| 135 | 
  | 
  | 
  if nk > n3d, | 
| 136 | 
  | 
  | 
    error([' too many levels(=',int2str(nk),' > ',int2str(n3d), ... | 
| 137 | 
  | 
  | 
           ' reading field: ',namV]); | 
| 138 | 
  | 
  | 
  end | 
| 139 | 
  | 
  | 
  if nk > 0, vvA(1:nk,:,:,:,nv)=vv1; end | 
| 140 | 
  | 
  | 
end | 
| 141 | 
  | 
  | 
fprintf(' <= end \n'); | 
| 142 | 
  | 
  | 
if length(msgA) ~= 1, fprintf(['  not found:',msgA,'\n']) ; end | 
| 143 | 
  | 
  | 
if length(msgB) ~= 1, fprintf(['  dim Pb:',msgB,'\n']) ; end | 
| 144 | 
  | 
  | 
if n3d == 0, n3d=1; vvA=zeros(n3d,nIt,nReg,5,nbV); end | 
| 145 | 
  | 
  | 
 | 
| 146 | 
  | 
  | 
return | 
| 147 | 
  | 
  | 
 | 
| 148 | 
  | 
  | 
function [vvm,nk,msg] = read_StD_1v(namF,nmV,sufx,nReg,nit,dit) | 
| 149 | 
  | 
  | 
 msg=' '; nk=0; fprintf([' ',nmV]); | 
| 150 | 
  | 
  | 
 namfil=[namF,'_',nmV,'.',sufx]; D=dir(namfil); | 
| 151 | 
  | 
  | 
 if size(D,1) == 1, | 
| 152 | 
  | 
  | 
   var=load(namfil); | 
| 153 | 
  | 
  | 
   nk=1+max(var(:,1)); | 
| 154 | 
  | 
  | 
   if size(var,1) == nk*nReg*nit, | 
| 155 | 
  | 
  | 
     vv1=reshape(var(:,2:6),[nk*nReg nit 5]); | 
| 156 | 
  | 
  | 
%--  put back Integrated Volume (was cumulated over the time period) in 5: | 
| 157 | 
  | 
  | 
     vv1(:,:,5)=vv1(:,:,5).*(ones(nk*nReg,1)*dit); | 
| 158 | 
  | 
  | 
     vv1=reshape(vv1,[nk nReg nit 5]); | 
| 159 | 
  | 
  | 
     vvm=permute(vv1,[1 3 2 4]); | 
| 160 | 
  | 
  | 
   else | 
| 161 | 
  | 
  | 
    %fprintf(['\n ERROR in Dim, var=',nmV,' :']); fprintf(' %i',size(var));  | 
| 162 | 
  | 
  | 
    %fprintf(' ; nk,nReg,nit= %i %i %i\n',nk,nReg,nit); | 
| 163 | 
  | 
  | 
     msg=sprintf([' in ',nmV,' : %i <> %ix%ix%i'],size(var,1),nk,nReg,nit); nk=0; | 
| 164 | 
  | 
  | 
     vvm=zeros(1,nit,nReg,5); | 
| 165 | 
  | 
  | 
     return | 
| 166 | 
  | 
  | 
   end | 
| 167 | 
  | 
  | 
 else msg=sprintf([msg,namfil,' ']); nk=-1; vvm=zeros(1,nit,nReg,5); | 
| 168 | 
  | 
  | 
 end | 
| 169 | 
  | 
  | 
return |