1 |
c*** Reads iobsloc |
2 |
subroutine ReadObs(obsnum,mobs,iobsloc) |
3 |
|
4 |
implicit none |
5 |
|
6 |
integer, intent(in) :: obsnum,mobs |
7 |
integer, intent(out) :: iobsloc(mobs) |
8 |
character*(100) fn |
9 |
real*8 r8seg(mobs) |
10 |
write(fn,'(A,I3.3,A)') '../inits/iobsloc-',obsnum,'.txt' |
11 |
write(*,*) fn |
12 |
|
13 |
open(unit=3,file=fn,status='old') |
14 |
read(unit=3,'(1I7)') iobsloc |
15 |
close(unit=3) |
16 |
call flush() |
17 |
end subroutine ReadObs |
18 |
|
19 |
c*** Call Model Kludge |
20 |
subroutine Model(mem) |
21 |
implicit none |
22 |
integer, intent(in) :: mem |
23 |
character*(1000) fn |
24 |
write(fn,'(A,I2.2,A)') |
25 |
& 'cd ../',mem,'/assimilate;./mitgcmuv>& LOG' |
26 |
call system(fn) |
27 |
end subroutine Model |
28 |
|
29 |
c*** Reads a pickup file... |
30 |
|
31 |
subroutine ReadPickup(mem,n,nx,ny,nz,y) |
32 |
|
33 |
implicit none |
34 |
|
35 |
integer, intent(in) :: mem,n,nx,ny,nz |
36 |
real*8, intent(out) :: y(n) |
37 |
real*8 r8seg(nx) |
38 |
integer i, irec,j,k |
39 |
character*(1024) fn |
40 |
write(fn,'(A,I2.2,A)') '../',mem,'/assimilate/pickup.in' |
41 |
open(unit=3,file=fn,status='old', |
42 |
& access='direct',recl=nx*8) |
43 |
|
44 |
ccc Set counter for state |
45 |
i =1 |
46 |
|
47 |
ccc Reading uVel |
48 |
do k = 1,nz |
49 |
do j = 1,ny |
50 |
irec = (k-1)*ny + j |
51 |
read(3,rec=irec) r8seg |
52 |
call DA_BYTESWAPR8(nx,r8seg) |
53 |
y(i:i+nx-1) = r8seg(1:nx) |
54 |
i = i+nx |
55 |
enddo |
56 |
enddo |
57 |
|
58 |
|
59 |
ccc Skip Next Two gU and gUnm1 |
60 |
ccc Read vVel |
61 |
do k = 1,nz |
62 |
do j = 1,ny |
63 |
irec = (k-1)*ny + j + nz*ny*3 |
64 |
read(3,rec=irec) r8seg |
65 |
call DA_BYTESWAPR8(nx,r8seg) |
66 |
y(i:i+nx-1) = r8seg(1:nx) |
67 |
i=i+nx |
68 |
enddo |
69 |
enddo |
70 |
|
71 |
|
72 |
ccc Skip Next two gV and gVnm1 |
73 |
ccc Read wVel |
74 |
do k = 1,nz |
75 |
do j = 1,ny |
76 |
irec = (k-1)*ny + j + nz*ny*6 |
77 |
read(3,rec=irec) r8seg |
78 |
call DA_BYTESWAPR8(nx,r8seg) |
79 |
y(i:i+nx-1) = r8seg(1:nx) |
80 |
i=i+nx |
81 |
enddo |
82 |
enddo |
83 |
|
84 |
|
85 |
ccc Read theta |
86 |
do k = 1,nz |
87 |
do j = 1,ny |
88 |
irec = (k-1)*ny + j + nz*ny*7 |
89 |
read(3,rec=irec) r8seg |
90 |
call DA_BYTESWAPR8(nx,r8seg) |
91 |
y(i:i+nx-1) = r8seg(1:nx) |
92 |
i=i+nx |
93 |
enddo |
94 |
enddo |
95 |
|
96 |
ccc Now close the file |
97 |
close(unit=3) |
98 |
ccc Open Non-Hydrostatic |
99 |
c write(fn,'(A,I2.2,A)') '../',mem,'/assimilate/pickup_nh.in' |
100 |
c open(unit=3,file=fn, |
101 |
c & status='old',access='direct',recl=nx*8) |
102 |
|
103 |
ccc Read pressure Non-hydrostatic |
104 |
c do k = 1,nz |
105 |
c do j = 1,ny |
106 |
c irec = (k-1)*ny + j |
107 |
c read(3,rec=irec) r8seg |
108 |
c call DA_BYTESWAPR8(nx,r8seg) |
109 |
c y(i:i+nx-1) = r8seg(1:nx) |
110 |
c i = i+nx |
111 |
c enddo |
112 |
c enddo |
113 |
|
114 |
c close(unit=3) |
115 |
ccc Done Reading State |
116 |
call flush() |
117 |
end subroutine ReadPickup |
118 |
|
119 |
ccc Writes a pickup |
120 |
subroutine WritePickup(mem,n,nx,ny,nz,y) |
121 |
|
122 |
implicit none |
123 |
|
124 |
integer, intent(in) :: mem,n,nx,ny,nz |
125 |
real*8, intent(in) :: y(n) |
126 |
real*8 r8seg(nx) |
127 |
integer i, irec,j,k |
128 |
character*(1000) fni, fno |
129 |
write(fni,'(A,I2.2,A)') '../',mem,'/assimilate/pickup.in' |
130 |
open(unit=3,file=fni,status='old', |
131 |
& access='direct',recl=nx*8) |
132 |
write(fni,'(A,I2.2,A)') '../',mem,'/assimilate/pickup.out' |
133 |
open(unit=4,file=fni,status='unknown', |
134 |
& access='direct',recl=nx*8) |
135 |
ccc Set counter for state |
136 |
i =1 |
137 |
|
138 |
ccc Writing uVel |
139 |
do k = 1,nz |
140 |
do j = 1,ny |
141 |
irec = (k-1)*ny + j |
142 |
r8seg(1:nx) = y(i:i+nx-1) |
143 |
call DA_BYTESWAPR8(nx,r8seg) |
144 |
write(4,rec=irec) r8seg |
145 |
i = i+nx |
146 |
enddo |
147 |
enddo |
148 |
|
149 |
ccc Copy gU and gUnm1 |
150 |
do k = 1,2*nz |
151 |
do j = 1,ny |
152 |
irec = (k-1)*ny + j+nz*ny |
153 |
read(3,rec=irec) r8seg |
154 |
write(4,rec=irec) r8seg |
155 |
enddo |
156 |
enddo |
157 |
ccc Writing vVel |
158 |
do k = 1,nz |
159 |
do j = 1,ny |
160 |
irec = (k-1)*ny + j + nz*ny*3 |
161 |
r8seg(1:nx) = y(i:i+nx-1) |
162 |
call DA_BYTESWAPR8(nx,r8seg) |
163 |
write(4,rec=irec) r8seg |
164 |
i=i+nx |
165 |
enddo |
166 |
enddo |
167 |
|
168 |
ccc copy Next two gV and gVnm1 |
169 |
do k = 1,2*nz |
170 |
do j = 1,ny |
171 |
irec = (k-1)*ny + j+nz*ny*4 |
172 |
read(3,rec=irec) r8seg |
173 |
write(4,rec=irec) r8seg |
174 |
enddo |
175 |
enddo |
176 |
ccc Write wVel and Theta |
177 |
do k = 1,2*nz |
178 |
do j = 1,ny |
179 |
irec = (k-1)*ny + j + nz*ny*6 |
180 |
r8seg(1:nx) = y(i:i+nx-1) |
181 |
call DA_BYTESWAPR8(nx,r8seg) |
182 |
write(4,rec=irec) r8seg |
183 |
i=i+nx |
184 |
enddo |
185 |
enddo |
186 |
ccc Copy Remaining fields |
187 |
do k = 1,5*nz |
188 |
do j = 1,ny |
189 |
irec = (k-1)*ny + j+nz*ny*8 |
190 |
read(3,rec=irec) r8seg |
191 |
write(4,rec=irec) r8seg |
192 |
enddo |
193 |
enddo |
194 |
do j = 1,ny |
195 |
irec = j+13*nz*ny |
196 |
read(3,rec=irec) r8seg |
197 |
write(4,rec=irec) r8seg |
198 |
enddo |
199 |
|
200 |
ccc Now close the file |
201 |
close(unit=3) |
202 |
close(unit=4) |
203 |
|
204 |
ccc Open Non-Hydrostatic |
205 |
write(fni,'(A,I2.2,A)') '../',mem,'/assimilate/pickup_nh.in' |
206 |
open(unit=3,file=fni, |
207 |
& status='old',access='direct',recl=nx*8) |
208 |
write(fni,'(A,I2.2,A)') '../',mem,'/assimilate/pickup_nh.out' |
209 |
open(unit=4,file=fni,status='unknown', |
210 |
& access='direct',recl=nx*8) |
211 |
|
212 |
ccc Write pressure Non-hydrostatic from State |
213 |
c do k = 1,nz |
214 |
c do j = 1,ny |
215 |
c irec = (k-1)*ny + j |
216 |
c r8seg(1:nx) = y(i:i+nx-1) |
217 |
c call DA_BYTESWAPR8(nx,r8seg) |
218 |
c write(4,rec=irec) r8seg |
219 |
c i = i+nx |
220 |
c enddo |
221 |
c enddo |
222 |
ccc Write out other fields |
223 |
do k = 1,2*nz |
224 |
do j = 1,ny |
225 |
irec = (k-1)*ny + j |
226 |
read(3,rec=irec) r8seg |
227 |
write(4,rec=irec) r8seg |
228 |
enddo |
229 |
enddo |
230 |
|
231 |
c close(unit=3) |
232 |
c close(unit=4) |
233 |
ccc Done Reading State |
234 |
call flush() |
235 |
end subroutine WritePickup |
236 |
|
237 |
|
238 |
|
239 |
|
240 |
subroutine loadMask(mask,nx,ny) |
241 |
implicit none |
242 |
integer, intent(in) :: nx, ny |
243 |
real*4, intent(out) :: mask(ny,nx) |
244 |
real*4 r4seg(nx) |
245 |
integer i, irec,j,k |
246 |
character*(1000) fni, fno |
247 |
|
248 |
write(fni,'(A)') '../inits/bathy61.bin' |
249 |
|
250 |
open(unit=3,file=fni, |
251 |
& status='old',access='direct',recl=nx*4) |
252 |
|
253 |
do j = 1,ny |
254 |
irec = j |
255 |
read(3,rec=irec) r4seg |
256 |
call DA_BYTESWAPR4(nx,r4seg) |
257 |
mask(j,1:nx) = r4seg(1:nx) |
258 |
enddo |
259 |
close(unit=3) |
260 |
call flush() |
261 |
end subroutine loadMask |
262 |
|
263 |
subroutine DA_BYTESWAPR8( na, arr ) |
264 |
C IN: |
265 |
C n integer - Number of 8-byte words in arr |
266 |
C IN/OUT: |
267 |
C arr real*8 - Array declared as real*4(n) |
268 |
C |
269 |
C Created: 05/05/99 adcroft@mit.edu (This is an unfortunate hack!!) |
270 |
|
271 |
implicit none |
272 |
C Arguments |
273 |
integer na |
274 |
character*(*) arr |
275 |
C Local |
276 |
integer i |
277 |
character*(1) cc |
278 |
C ------------------------------------------------------------------ |
279 |
do i=1,8*na,8 |
280 |
cc=arr(i:i) |
281 |
arr(i:i)=arr(i+7:i+7) |
282 |
arr(i+7:i+7)=cc |
283 |
cc=arr(i+1:i+1) |
284 |
arr(i+1:i+1)=arr(i+6:i+6) |
285 |
arr(i+6:i+6)=cc |
286 |
cc=arr(i+2:i+2) |
287 |
arr(i+2:i+2)=arr(i+5:i+5) |
288 |
arr(i+5:i+5)=cc |
289 |
cc=arr(i+3:i+3) |
290 |
arr(i+3:i+3)=arr(i+4:i+4) |
291 |
arr(i+4:i+4)=cc |
292 |
enddo |
293 |
C ------------------------------------------------------------------ |
294 |
return |
295 |
end |
296 |
C========================== |
297 |
|
298 |
|
299 |
C======================================================================= |
300 |
subroutine DA_BYTESWAPR4( n, arr ) |
301 |
C IN: |
302 |
C n integer - Number of 4-byte words in arr |
303 |
C IN/OUT: |
304 |
C arr real*4 - Array declared as real*4(n) |
305 |
C |
306 |
C Created: 05/05/99 adcroft@mit.edu (This is an unfortunate hack!!) |
307 |
|
308 |
implicit none |
309 |
C Arguments |
310 |
integer n |
311 |
character*(*) arr |
312 |
C Local |
313 |
integer i |
314 |
character*(1) cc |
315 |
C ------------------------------------------------------------------ |
316 |
do i=1,4*n,4 |
317 |
cc=arr(i:i) |
318 |
arr(i:i)=arr(i+3:i+3) |
319 |
arr(i+3:i+3)=cc |
320 |
cc=arr(i+1:i+1) |
321 |
arr(i+1:i+1)=arr(i+2:i+2) |
322 |
arr(i+2:i+2)=cc |
323 |
enddo |
324 |
C ------------------------------------------------------------------ |
325 |
return |
326 |
end |
327 |
C======================================================================= |
328 |
|
329 |
|
330 |
|
331 |
|
332 |
|
333 |
|