/[MITgcm]/MITgcm_contrib/afe/osse_MkII/piv/datamaker.m
ViewVC logotype

Contents of /MITgcm_contrib/afe/osse_MkII/piv/datamaker.m

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


Revision 1.1 - (show annotations) (download)
Fri Jul 22 17:53:22 2005 UTC (20 years ago) by afe
Branch: MAIN
CVS Tags: HEAD
added datafile maker

1 nx=120;ny=31;nz=29;n=nx*ny*nz;
2 %xdata=nx;ydata=ny;zdata=nz;ndata=xdata*ydata*zdata
3 xdata=90;ydata=21;zdata=29;ndata=xdata*ydata*zdata;
4 xdata=90;ydata=21;zdata=29;ndata=xdata*ydata*zdata;
5 %xdata=120;ydata=31;zdata=29;ndata=xdata*ydata*zdata;
6 %xout=90;yout=21;zout=29;nout=xout*yout*zout;
7 %iter=68;
8 z=15;
9 vars=['U' 'V' 'W' 'P' 'T'];
10 vars=['U' 'V' 'W' 'T'];
11 vars=['U' 'V' 'T'];
12 U=find(vars=='U');
13 V=find(vars=='V');
14 %W=find(vars=='W');
15 %P=find(vars=='P');
16 T=find(vars=='T');
17 varnum= size(vars,2);
18
19 if(0)
20 filename=sprintf('%s','../da/00/assimilate/pickup.in')
21 %filename=sprintf('%s%03i%s','../da/00/assimilate/pickup.',iter,'.out');
22 %filename=sprintf('%s','../da/foo/00/assimilate/pickup.0000016820.001.001.data')
23 fid=fopen(filename,'r','ieee-be');
24 foo=fread(fid,13*ndata+ydata*xdata,'float64');
25 i=0;tru(:,:,:,U)=reshape(foo(i*ndata+1:(i+1)*ndata),[xdata ydata zdata]);
26 i=3;tru(:,:,:,V)=reshape(foo(i*ndata+1:(i+1)*ndata),[xdata ydata zdata]);
27 i=6;tru(:,:,:,W)=reshape(foo(i*ndata+1:(i+1)*ndata),[xdata ydata zdata]);
28 i=7;tru(:,:,:,T)=reshape(foo(i*ndata+1:(i+1)*ndata),[xdata ydata zdata]);
29 fclose(fid);
30 end
31
32 interval=100;
33 %last=16800+200;
34 last=17800;
35 %last=16800+20*100;
36 first=16800;
37 %interval=10;
38 %last=1000;
39 %first=200;
40
41 %for k=first:interval:last
42 for k=1:200
43 iter=16800+(k-1)*100;
44
45 if(0)
46 filename=sprintf('%s%010i%s','pickup.',k,'.001.001.data')
47 fid=fopen(filename,'r','ieee-be');
48 foo=fread(fid,13*ndata+ydata*xdata,'float64');
49 i=0;tru(:,:,:,U)=reshape(foo(i*ndata+1:(i+1)*ndata),[xdata ydata zdata]);
50 i=3;tru(:,:,:,V)=reshape(foo(i*ndata+1:(i+1)*ndata),[xdata ydata zdata]);
51 %i=6;tru(:,:,:,W)=reshape(foo(i*ndata+1:(i+1)*ndata),[xdata ydata zdata]);
52 i=7;tru(:,:,:,T)=reshape(foo(i*ndata+1:(i+1)*ndata),[xdata ydata zdata]);
53 fclose(fid);
54
55 Up=tru(:,:,[5:5:25],U);
56 Vp=tru(:,:,[5:5:25],V);
57 Up(find(Up==0))=NaN;
58 Vp(find(Vp==0))=NaN;
59 Tp=tru(:,:,:,T);
60 Tp(:,[1:5 7:19 21],:)=NaN;
61
62 edge=20;
63 step=1;
64 method='nearest';
65 [u v]=cyl2cartuv(Up,Vp, -edge:step:edge, -edge:step:edge,method);
66 t=cyl2cart(Tp, -edge:step:edge, -edge:step:edge,method);
67 end
68
69 %clear u,v
70 u=0.0015*GU1(:,:,(5*(k-1)+1):(5*(k-1)+5));
71 v=0.0015*GV1(:,:,(5*(k-1)+1):(5*(k-1)+5));
72
73 %u=ones(27,37,5)./100;
74 %v=ones(27,37,5)./100;
75
76
77 if(1)
78 iobsloc=[];
79 datavect=[];
80
81 loccount=0;
82 [xout yout zout]=size(u)
83 for j=1:zout
84 u(:,:,j)=flipud(u(:,:,j));
85 for i=1:yout
86 for h=1:xout
87 % loccount=loccount+1;
88 % if u(h,i,j) ~= nan
89 if ~isnan(u(h,i,j))
90 datavect(end+1)=u(h,i,j);
91 iobsloc(end+1)=((j*5)-1)*xout*yout+(i-1)*xout+h;
92 end % if u(h,i,j) ~= nan
93 end % h=1:xout
94 end % i=1:yout
95 end % j=1:zout
96
97 [xout yout zout]=size(v);
98 %for j=1:1
99 for j=1:zout
100 v(:,:,j)=flipud(v(:,:,j));
101 for i=1:yout
102 for h=1:xout
103 % loccount=loccount+1;
104 % if v(h,i,j) ~= nan
105 if ~isnan(v(h,i,j))
106 datavect(end+1)=v(h,i,j);
107 iobsloc(end+1)=xout*yout*zdata+((j*5)-1)*xout*yout+(i-1)*xout+h;
108 end % if u(h,i,j) ~= nan
109 end % h=1:xout
110 end % i=1:yout
111 end % j=1:zout
112 end
113
114 if(0)
115 [xout yout zout]=size(t);
116 for j=1:zout
117 for i=1:yout
118 for h=1:xout
119 % loccount=loccount+1;
120 % if t(h,i,j) ~= nan
121 if ~isnan(t(h,i,j))
122 datavect(end+1)=t(h,i,j);
123 iobsloc(end+1)=xout*yout*zdata*3+((j)-1)*xout*yout+(i-1)*xout+h;
124 end % if u(h,i,j) ~= nan
125 end % h=1:xout
126 end % i=1:yout
127 end % j=1:zout
128 end
129
130 [meep locsize]=size(iobsloc)
131 [meep vectsize]=size(datavect)
132
133 if (locsize ~= vectsize)
134 locsize
135 vectsize
136 error('sizes wrong!')
137 end
138
139 mobsize=2*xout*yout*zdata+2*zdata*xout;
140
141 iobsloc(end+1:mobsize)=999999999;
142 datavect(end+1:mobsize)=999999999;
143 if(1)
144 %filename=sprintf('%s%010i%s','iobsloc.',iter,'.txt')
145 filename=sprintf('%s%03i%s','iobsloc.',k,'.txt')
146 fid = fopen(filename,'w');
147 fprintf(fid,'%i\n',iobsloc);
148 fclose(fid);
149
150 dtype='float64';
151 %filename=sprintf('%s%010i%s','obs.',iter,'.dat')
152 filename=sprintf('%s%03i%s','obs.',k,'.dat')
153 fid = fopen(filename,'w','ieee-be');
154 fwrite(fid,datavect,dtype);
155 fclose(fid);
156 end
157
158
159
160
161 %for j=1:size(iobsloc)
162 % sq=iobsloc(j)
163 % zig=sq
164
165
166 end

  ViewVC Help
Powered by ViewVC 1.1.22