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

Contents of /MITgcm_contrib/gael/matlab_class/gcmfaces_IO/read_nctiles.m

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


Revision 1.7 - (show annotations) (download)
Fri Mar 13 16:16:09 2015 UTC (10 years, 4 months ago) by gforget
Branch: MAIN
CVS Tags: checkpoint65x, checkpoint65r, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, HEAD
Changes since 1.6: +42 -9 lines
contribution from M. Buckley: allow for
-specification of kk without specification of tt
-a range of contiguous indicies for kk  rather than just a single level

1 function [fld]=read_nctiles(fileName,fldName,varargin);
2 %usage: fld=read_nctiles(fileName); reads full field using fileName as field name
3 %usage: fld=read_nctiles(fileName,fldName); reads full field (all depths, all times)
4 %usage: fld=read_nctiles(fileName,fldName,tt); reads 3D or 2D field, at time index tt (all depths)
5 %usage: fld=read_nctiles(fileName,fldName,tt,kk); reads 3D field, at depth index(es) kk, at time index tt
6
7 gcmfaces_global;
8 nz=length(mygrid.RC);
9
10 if nargin==1;
11 tmp1=['/' fileName];
12 tmp2=strfind(tmp1,filesep);
13 fldName=tmp1(tmp2(end)+1:end);
14 end;
15 if nargin>2;
16 tt=varargin{1};
17 else;
18 tt=[];
19 end;
20 if nargin>3;
21 kk=varargin{2};
22 if min(kk)<1 | max(kk)>nz
23 error(['kk must be an integer between 1 and ' num2str(nz)])
24 return
25 end
26 if length(kk)>1
27 if ~isempty(find(diff(kk)~=1))
28 error(['kk must be a contigous range of indices'])
29 return
30 end
31 end
32 else;
33 kk=[];
34 end;
35
36 %A) determine map of tile indices (if not already done)
37
38 global nctiles;
39
40 test1=isempty(nctiles);
41 test2=0;
42 if ~test1;
43 tmp1=dir([fileName '*']);
44 test2=(length(tmp1)~=length(nctiles.no));
45 end;
46
47 if test1|test2;
48 %build map of tile indices
49 if (length(dir([fileName '*']))==mygrid.nFaces);
50 nctiles.map=NaN*mygrid.XC;
51 for ff=1:mygrid.nFaces; nctiles.map{ff}(:)=ff; end;
52 else;
53 fileIn=sprintf('%s.%04d.nc',fileName,1);
54 nc=netcdf.open(fileIn,0);
55 vv = netcdf.inqVarID(nc,fldName);
56 [varname,xtype,dimids,natts]=netcdf.inqVar(nc,vv);
57 tileSize=zeros(1,2);
58 [tmp,tileSize(1)] = netcdf.inqDim(nc,dimids(1));
59 [tmp,tileSize(2)] = netcdf.inqDim(nc,dimids(2));
60 netcdf.close(nc);
61 nctiles.map=gcmfaces_loc_tile(tileSize(1),tileSize(2));
62 end;
63 %determine tile list
64 nctiles.no=unique(convert2vector(nctiles.map));
65 nctiles.no=nctiles.no(~isnan(nctiles.no));
66 %determine location of each tile
67 for ff=1:length(nctiles.no);
68 for gg=1:mygrid.nFaces;
69 [tmpi,tmpj]=find(nctiles.map{gg}==ff);
70 if ~isempty(tmpi);
71 nctiles.f{ff}=gg;
72 nctiles.i{ff}=[min(tmpi(:)):max(tmpi(:))];
73 nctiles.j{ff}=[min(tmpj(:)):max(tmpj(:))];
74 end;
75 end;
76 end;
77 end;
78
79 %B) the file read operation itself
80
81 for ff=1:length(nctiles.no);
82
83 %read one tile
84 fileIn=sprintf('%s.%04d.nc',fileName,ff);
85 nc=netcdf.open(fileIn,0);
86
87 vv = netcdf.inqVarID(nc,fldName);
88 [varname,xtype,dimids,natts]=netcdf.inqVar(nc,vv);
89
90 [dimname,siz(1)] = netcdf.inqDim(nc,dimids(1));
91 [dimname,siz(2)] = netcdf.inqDim(nc,dimids(2));
92 siz=[siz length(mygrid.RC)];
93
94 if ~isempty(tt);
95 t0=tt(1)-1;
96 nt=tt(end)-tt(1)+1;
97 if length(dimids)==3;
98 start=[0 0 t0];
99 count=[siz(1) siz(2) nt];
100 elseif isempty(kk);
101 start=[0 0 0 t0];
102 count=[siz(1) siz(2) siz(3) nt];
103 else;
104 start=[0 0 kk(1)-1 t0];
105 count=[siz(1) siz(2) length(kk) nt];
106 end;
107 else
108 [dimname,nt] = netcdf.inqDim(nc,dimids(end));
109 if length(dimids)==2;
110 start=[0 0];
111 count=[siz(1) siz(2)];
112 elseif length(dimids)==3;
113 start=[0 0 0];
114 count=[siz(1) siz(2) nt];
115 elseif isempty(kk);
116 start=[0 0 0 0];
117 count=[siz(1) siz(2) siz(3) nt];
118 else;
119 start=[0 0 kk(1)-1 0];
120 count=[siz(1) siz(2) length(kk) nt];
121 end;
122 end;
123 fldTile=netcdf.getVar(nc,vv,start,count);
124 fldTile=squeeze(fldTile);
125 netcdf.close(nc);
126
127 %initialize fld (full gcmfaces object)
128 if ff==1;
129 siz=[1 1 size(fldTile,3) size(fldTile,4)];
130 fld=NaN*repmat(mygrid.XC,siz);
131 end;
132
133 %place tile in fld
134 fld{nctiles.f{ff}}(nctiles.i{ff},nctiles.j{ff},:,:)=fldTile;
135
136 end;%for ff=1:mygrid.nFaces;
137
138
139

  ViewVC Help
Powered by ViewVC 1.1.22