/[MITgcm]/MITgcm/verification/natl_box_forward/matlab/writegcm.m
ViewVC logotype

Diff of /MITgcm/verification/natl_box_forward/matlab/writegcm.m

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

revision 1.1 by heimbach, Tue Feb 5 20:34:50 2002 UTC revision 1.1.2.1 by heimbach, Tue Feb 5 20:34:50 2002 UTC
# Line 0  Line 1 
1    function writegcm(fnam,fld,tme,nPx,nPy,prec,oflag)
2    
3    % Function writegcm(fnam,fld,tme,nPx,nPy,prec,oflag)
4    % write array fld in MIT GCM UV format
5    %
6    % INPUTS
7    % fnam      output file name
8    % fld       Nx*Ny*Nz*Nt input array
9    % tme       model integration times in days has dimension Nt
10    % nPx, nPy  number of processes in x and y direction (default 1,1)
11    % prec      numeric precision (default 'real*4')
12    % oflag     overwrite flag specifies field to be overwritten
13    %           oflag=0 means create new file (the default)
14    %           oflag>0 must have dimension Nt
15    %
16    % SEE ALSO
17    % readgcm
18    
19    if nargin < 2, error('please specify array and output file name'); end
20    [Nx Ny Nz Nt]=size(fld);
21    if nargin < 7, oflag=0; end
22    if oflag > 0
23      if length(oflag) ~= Nt, error('length(oflag) must equal Nt'); end
24    end
25    if nargin < 6, prec='real*4'; end
26    if nargin < 5, nPy=1; end
27    if nargin < 4, nPx=1; end
28    if nargin < 3, tme=1:Nt; end
29    if length(tme) ~= Nt, error('length(tme) must equal Nt'); end
30    
31    sNx=Nx/nPx;                                 % number of X points in sub-grid
32    sNy=Ny/nPy;                                 % number of Y points in sub-grid
33    if (sNx-floor(sNx))~=0 | (sNy-floor(sNy))~=0
34      error('Nx/nPx and Ny/nPy must be integer')
35    end
36    
37    switch prec
38      case {'float32', 'real*4'}
39        rlength=(sNx*sNy*Nz+1)*4;
40      case {'float64', 'real*8'}
41        rlength=(sNx*sNy*Nz+1)*8;
42    end
43    
44    if oflag == 0
45    
46      fid=fopen(fnam,'w','ieee-be');
47      for t=1:Nt
48        for j=1:nPy
49          jx=((j-1)*sNy+1):(j*sNy);
50          for i=1:nPx
51            ix=((i-1)*sNx+1):(i*sNx);
52            fwrite(fid,tme(t),prec);
53            fwrite(fid,fld(ix,jx,:,t),prec);
54          end
55        end
56      end
57    
58    else
59    
60      fid=fopen(fnam,'r+','ieee-be');
61      for t=1:Nt
62        if fseek(fid,nPx*nPy*(oflag(t)-1)*rlength,'bof')
63          error('unexpected end of file reached')
64        end
65        for j=1:nPy
66          jx=((j-1)*sNy+1):(j*sNy);
67          for i=1:nPx
68            ix=((i-1)*sNx+1):(i*sNx);
69            fwrite(fid,tme(t),prec);
70            fwrite(fid,fld(ix,jx,:,t),prec);
71          end
72        end
73      end
74    
75    end
76    
77    fid=fclose(fid);

Legend:
Removed from v.1.1  
changed lines
  Added in v.1.1.2.1

  ViewVC Help
Powered by ViewVC 1.1.22