/[MITgcm]/MITgcm_contrib/osse/EnKF/rwObsState3.F
ViewVC logotype

Annotation of /MITgcm_contrib/osse/EnKF/rwObsState3.F

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


Revision 1.1 - (hide annotations) (download)
Tue May 4 18:19:35 2004 UTC (21 years, 2 months ago) by afe
Branch: MAIN
o EnKF stuff

1 afe 1.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    

  ViewVC Help
Powered by ViewVC 1.1.22