/[MITgcm]/MITgcm/pkg/ctrl/ctrl_set_pack_xyz.F
ViewVC logotype

Contents of /MITgcm/pkg/ctrl/ctrl_set_pack_xyz.F

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


Revision 1.9 - (show annotations) (download)
Tue Jan 4 22:02:31 2005 UTC (19 years, 5 months ago) by heimbach
Branch: MAIN
Changes since 1.8: +75 -10 lines
o Add ctrlvec diagnostics in pack/unpack for nondimensional I/O
o May be enabled via doPackDiag

1
2 #include "CTRL_CPPOPTIONS.h"
3
4 subroutine ctrl_set_pack_xyz(
5 & cunit, ivartype, fname, masktype, weighttype,
6 & weightfld, lxxadxx, mythid)
7
8 c ==================================================================
9 c SUBROUTINE ctrl_set_pack_xyz
10 c ==================================================================
11 c
12 c o Compress the control vector such that only ocean points are
13 c written to file.
14 c
15 c o Use a more precise nondimensionalization that depends on (x,y)
16 c Added weighttype to the argument list so that I can geographically
17 c vary the nondimensionalization.
18 c gebbie@mit.edu, 18-Mar-2003
19 c
20 c ==================================================================
21
22 implicit none
23
24 c == global variables ==
25
26 #include "EEPARAMS.h"
27 #include "SIZE.h"
28 #include "PARAMS.h"
29 #include "GRID.h"
30
31 #include "ctrl.h"
32 #include "optim.h"
33
34 c == routine arguments ==
35
36 integer cunit
37 integer ivartype
38 character*( 80) fname
39 character* (5) masktype
40 character*( 80) weighttype
41 _RL weightfld( nr,nsx,nsy )
42 logical lxxadxx
43 integer mythid
44
45 c == local variables ==
46
47 integer bi,bj
48 integer ip,jp
49 integer i,j,k
50 integer ii
51 integer il
52 integer irec
53 integer itlo,ithi
54 integer jtlo,jthi
55 integer jmin,jmax
56 integer imin,imax
57
58 integer cbuffindex
59
60 _RL globmsk ( snx,nsx,npx,sny,nsy,npy,nr )
61 _RL globfld3d( snx,nsx,npx,sny,nsy,npy,nr )
62 #ifdef CTRL_PACK_PRECISE
63 _RL weightfld3d( snx,nsx,npx,sny,nsy,npy,nr )
64 #endif
65 real*4 cbuff ( snx*nsx*npx*sny*nsy*npy )
66 real*4 globfldtmp2( snx,nsx,npx,sny,nsy,npy )
67 real*4 globfldtmp3( snx,nsx,npx,sny,nsy,npy )
68
69 character*(80) weightname
70
71 #ifdef CTRL_DELZNORM
72 _RL delZnorm
73 #endif
74 integer reclen, irectrue
75 integer cunit2, cunit3
76 character*(80) cfile2, cfile3
77
78 c == external ==
79
80 integer ilnblnk
81 external ilnblnk
82
83 c == end of interface ==
84
85 jtlo = 1
86 jthi = nsy
87 itlo = 1
88 ithi = nsx
89 jmin = 1
90 jmax = sny
91 imin = 1
92 imax = snx
93
94 #ifdef CTRL_DELZNORM
95 delZnorm = 0.
96 do k = 1, Nr
97 delZnorm = delZnorm + delR(k)/FLOAT(Nr)
98 enddo
99 #endif
100
101 c Initialise temporary file
102 do k = 1,nr
103 do jp = 1,nPy
104 do bj = jtlo,jthi
105 do j = jmin,jmax
106 do ip = 1,nPx
107 do bi = itlo,ithi
108 do i = imin,imax
109 globfld3d (i,bi,ip,j,bj,jp,k) = 0. _d 0
110 globmsk (i,bi,ip,j,bj,jp,k) = 0. _d 0
111 globfldtmp2(i,bi,ip,j,bj,jp) = 0.
112 globfldtmp3(i,bi,ip,j,bj,jp) = 0.
113 enddo
114 enddo
115 enddo
116 enddo
117 enddo
118 enddo
119 enddo
120
121 c-- Only the master thread will do I/O.
122 _BEGIN_MASTER( mythid )
123
124 if ( doPackDiag ) then
125 write(cfile2(1:80),'(80a)') ' '
126 write(cfile3(1:80),'(80a)') ' '
127 if ( lxxadxx ) then
128 write(cfile2(1:80),'(a,I2.2,a,I4.4,a)')
129 & 'diag_pack_nonout_ctrl_',
130 & ivartype, '_', optimcycle, '.bin'
131 write(cfile3(1:80),'(a,I2.2,a,I4.4,a)')
132 & 'diag_pack_dimout_ctrl_',
133 & ivartype, '_', optimcycle, '.bin'
134 else
135 write(cfile2(1:80),'(a,I2.2,a,I4.4,a)')
136 & 'diag_pack_nonout_grad_',
137 & ivartype, '_', optimcycle, '.bin'
138 write(cfile3(1:80),'(a,I2.2,a,I4.4,a)')
139 & 'diag_pack_dimout_grad_',
140 & ivartype, '_', optimcycle, '.bin'
141 endif
142
143 reclen = FLOAT(snx*nsx*npx*sny*nsy*npy*4)
144 call mdsfindunit( cunit2, mythid )
145 open( cunit2, file=cfile2, status='unknown',
146 & access='direct', recl=reclen )
147 call mdsfindunit( cunit3, mythid )
148 open( cunit3, file=cfile3, status='unknown',
149 & access='direct', recl=reclen )
150 endif
151
152 #ifdef CTRL_PACK_PRECISE
153 il=ilnblnk( weighttype)
154 write(weightname(1:80),'(80a)') ' '
155 write(weightname(1:80),'(a)') weighttype(1:il)
156
157 call MDSREADFIELD_3D_GL(
158 & weightname, ctrlprec, 'RL',
159 & Nr, weightfld3d, 1, mythid)
160 #endif
161
162 call MDSREADFIELD_3D_GL(
163 & masktype, ctrlprec, 'RL',
164 & Nr, globmsk, 1, mythid)
165
166 do irec = 1, ncvarrecs(ivartype)
167
168 call MDSREADFIELD_3D_GL( fname, ctrlprec, 'RL',
169 & Nr, globfld3d, irec, mythid)
170
171 write(cunit) ncvarindex(ivartype)
172 write(cunit) 1
173 write(cunit) 1
174 do k = 1, nr
175 irectrue = (irec-1)*nr + k
176 cbuffindex = 0
177 do jp = 1,nPy
178 do bj = jtlo,jthi
179 do j = jmin,jmax
180 do ip = 1,nPx
181 do bi = itlo,ithi
182 do i = imin,imax
183 if (globmsk(i,bi,ip,j,bj,jp,k) .ne. 0. ) then
184 cbuffindex = cbuffindex + 1
185 cph(
186 globfldtmp3(i,bi,ip,j,bj,jp) =
187 & globfld3d(i,bi,ip,j,bj,jp,k)
188 cph)
189 #ifdef ALLOW_NONDIMENSIONAL_CONTROL_IO
190 if (lxxadxx) then
191 cbuff(cbuffindex) =
192 & globfld3d(i,bi,ip,j,bj,jp,k)
193 # ifdef CTRL_PACK_PRECISE
194 & * sqrt(weightfld3d(i,bi,ip,j,bj,jp,k))
195 # else
196 & * sqrt(weightfld(k,bi,bj))
197 # endif
198 else
199 cbuff(cbuffindex) =
200 & globfld3d(i,bi,ip,j,bj,jp,k)
201 # ifdef CTRL_PACK_PRECISE
202 & / sqrt(weightfld3d(i,bi,ip,j,bj,jp,k))
203 # else
204 & / sqrt(weightfld(k,bi,bj))
205 # endif
206 endif
207 cph(
208 globfldtmp2(i,bi,ip,j,bj,jp) = cbuff(cbuffindex)
209 cph)
210 #else /* ALLOW_NONDIMENSIONAL_CONTROL_IO undef */
211 cbuff(cbuffindex) = globfld3d(i,bi,ip,j,bj,jp,k)
212 #endif /* ALLOW_NONDIMENSIONAL_CONTROL_IO */
213 endif
214 enddo
215 enddo
216 enddo
217 enddo
218 enddo
219 enddo
220 c --> check cbuffindex.
221 if ( cbuffindex .gt. 0) then
222 write(cunit) cbuffindex
223 write(cunit) k
224 write(cunit) (cbuff(ii), ii=1,cbuffindex)
225 endif
226 c
227 if ( doPackDiag ) then
228 write(cunit2,rec=irectrue) globfldtmp2
229 write(cunit3,rec=irectrue) globfldtmp3
230 endif
231 c
232 enddo
233 c
234 c -- end of irec loop --
235 enddo
236
237 if ( doPackDiag ) then
238 close ( cunit2 )
239 close ( cunit3 )
240 endif
241
242 _END_MASTER( mythid )
243
244 return
245 end
246

  ViewVC Help
Powered by ViewVC 1.1.22