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

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

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


Revision 1.14 - (hide annotations) (download)
Thu Apr 27 12:50:39 2006 UTC (18 years ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint58e_post
Changes since 1.13: +2 -1 lines
o supressing admtlm-related vector output for now
  (such ad admtlm_vector, admtlm_eigen)

1 heimbach 1.2
2     #include "CTRL_CPPOPTIONS.h"
3    
4     subroutine ctrl_set_pack_xyz(
5 heimbach 1.4 & cunit, ivartype, fname, masktype, weighttype,
6 heimbach 1.2 & 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 heimbach 1.4 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 heimbach 1.2 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 heimbach 1.11 character*( 9) masktype
40 heimbach 1.4 character*( 80) weighttype
41 heimbach 1.2 _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 heimbach 1.4 #ifdef CTRL_PACK_PRECISE
63     _RL weightfld3d( snx,nsx,npx,sny,nsy,npy,nr )
64     #endif
65 heimbach 1.9 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 heimbach 1.4
69     character*(80) weightname
70 heimbach 1.2
71 heimbach 1.9 _RL delZnorm
72     integer reclen, irectrue
73     integer cunit2, cunit3
74     character*(80) cfile2, cfile3
75    
76 heimbach 1.2 c == external ==
77    
78     integer ilnblnk
79     external ilnblnk
80    
81     c == end of interface ==
82    
83     jtlo = 1
84     jthi = nsy
85     itlo = 1
86     ithi = nsx
87     jmin = 1
88     jmax = sny
89     imin = 1
90     imax = snx
91    
92 heimbach 1.9 #ifdef CTRL_DELZNORM
93     delZnorm = 0.
94     do k = 1, Nr
95     delZnorm = delZnorm + delR(k)/FLOAT(Nr)
96     enddo
97     #endif
98    
99 heimbach 1.2 c Initialise temporary file
100     do k = 1,nr
101     do jp = 1,nPy
102     do bj = jtlo,jthi
103     do j = jmin,jmax
104     do ip = 1,nPx
105     do bi = itlo,ithi
106     do i = imin,imax
107 heimbach 1.9 globfld3d (i,bi,ip,j,bj,jp,k) = 0. _d 0
108     globmsk (i,bi,ip,j,bj,jp,k) = 0. _d 0
109     globfldtmp2(i,bi,ip,j,bj,jp) = 0.
110     globfldtmp3(i,bi,ip,j,bj,jp) = 0.
111 heimbach 1.2 enddo
112     enddo
113     enddo
114     enddo
115     enddo
116     enddo
117     enddo
118    
119     c-- Only the master thread will do I/O.
120     _BEGIN_MASTER( mythid )
121    
122 heimbach 1.9 if ( doPackDiag ) then
123     write(cfile2(1:80),'(80a)') ' '
124     write(cfile3(1:80),'(80a)') ' '
125     if ( lxxadxx ) then
126     write(cfile2(1:80),'(a,I2.2,a,I4.4,a)')
127     & 'diag_pack_nonout_ctrl_',
128     & ivartype, '_', optimcycle, '.bin'
129     write(cfile3(1:80),'(a,I2.2,a,I4.4,a)')
130     & 'diag_pack_dimout_ctrl_',
131     & ivartype, '_', optimcycle, '.bin'
132     else
133     write(cfile2(1:80),'(a,I2.2,a,I4.4,a)')
134     & 'diag_pack_nonout_grad_',
135     & ivartype, '_', optimcycle, '.bin'
136     write(cfile3(1:80),'(a,I2.2,a,I4.4,a)')
137     & 'diag_pack_dimout_grad_',
138     & ivartype, '_', optimcycle, '.bin'
139     endif
140    
141     reclen = FLOAT(snx*nsx*npx*sny*nsy*npy*4)
142     call mdsfindunit( cunit2, mythid )
143     open( cunit2, file=cfile2, status='unknown',
144     & access='direct', recl=reclen )
145     call mdsfindunit( cunit3, mythid )
146     open( cunit3, file=cfile3, status='unknown',
147     & access='direct', recl=reclen )
148     endif
149    
150 heimbach 1.4 #ifdef CTRL_PACK_PRECISE
151     il=ilnblnk( weighttype)
152     write(weightname(1:80),'(80a)') ' '
153     write(weightname(1:80),'(a)') weighttype(1:il)
154    
155     call MDSREADFIELD_3D_GL(
156     & weightname, ctrlprec, 'RL',
157     & Nr, weightfld3d, 1, mythid)
158     #endif
159    
160 heimbach 1.2 call MDSREADFIELD_3D_GL(
161     & masktype, ctrlprec, 'RL',
162     & Nr, globmsk, 1, mythid)
163    
164     do irec = 1, ncvarrecs(ivartype)
165    
166     call MDSREADFIELD_3D_GL( fname, ctrlprec, 'RL',
167     & Nr, globfld3d, irec, mythid)
168    
169 heimbach 1.13 #ifndef ALLOW_ADMTLM
170 heimbach 1.2 write(cunit) ncvarindex(ivartype)
171     write(cunit) 1
172     write(cunit) 1
173 heimbach 1.13 #endif
174 heimbach 1.2 do k = 1, nr
175 heimbach 1.9 irectrue = (irec-1)*nr + k
176 heimbach 1.10 if ( doZscalePack ) then
177 heimbach 1.12 cph delZnorm = SQRT(delR(1)/delR(k))
178     delZnorm = delR(1)/delR(k)
179 heimbach 1.10 else
180     delZnorm = 1. _d 0
181     endif
182 heimbach 1.2 cbuffindex = 0
183     do jp = 1,nPy
184     do bj = jtlo,jthi
185     do j = jmin,jmax
186     do ip = 1,nPx
187     do bi = itlo,ithi
188     do i = imin,imax
189 heimbach 1.9 if (globmsk(i,bi,ip,j,bj,jp,k) .ne. 0. ) then
190 heimbach 1.2 cbuffindex = cbuffindex + 1
191 heimbach 1.9 cph(
192     globfldtmp3(i,bi,ip,j,bj,jp) =
193     & globfld3d(i,bi,ip,j,bj,jp,k)
194     cph)
195 heimbach 1.2 #ifdef ALLOW_NONDIMENSIONAL_CONTROL_IO
196     if (lxxadxx) then
197 heimbach 1.10 cbuff(cbuffindex) = delZnorm
198     & * globfld3d(i,bi,ip,j,bj,jp,k)
199 heimbach 1.4 # ifdef CTRL_PACK_PRECISE
200 heimbach 1.9 & * sqrt(weightfld3d(i,bi,ip,j,bj,jp,k))
201 heimbach 1.4 # else
202 heimbach 1.9 & * sqrt(weightfld(k,bi,bj))
203 heimbach 1.4 # endif
204 heimbach 1.2 else
205 heimbach 1.10 cbuff(cbuffindex) = delZnorm
206     & * globfld3d(i,bi,ip,j,bj,jp,k)
207 heimbach 1.4 # ifdef CTRL_PACK_PRECISE
208 heimbach 1.9 & / sqrt(weightfld3d(i,bi,ip,j,bj,jp,k))
209 heimbach 1.4 # else
210 heimbach 1.9 & / sqrt(weightfld(k,bi,bj))
211 heimbach 1.4 # endif
212 heimbach 1.2 endif
213 heimbach 1.9 cph(
214     globfldtmp2(i,bi,ip,j,bj,jp) = cbuff(cbuffindex)
215     cph)
216 heimbach 1.4 #else /* ALLOW_NONDIMENSIONAL_CONTROL_IO undef */
217 heimbach 1.2 cbuff(cbuffindex) = globfld3d(i,bi,ip,j,bj,jp,k)
218 heimbach 1.4 #endif /* ALLOW_NONDIMENSIONAL_CONTROL_IO */
219 heimbach 1.13 #ifdef ALLOW_ADMTLM
220     nveccount = nveccount + 1
221     phtmpadmtlm(nveccount) = cbuff(cbuffindex)
222     #endif
223 heimbach 1.2 endif
224     enddo
225     enddo
226     enddo
227     enddo
228     enddo
229     enddo
230     c --> check cbuffindex.
231     if ( cbuffindex .gt. 0) then
232 heimbach 1.13 #ifndef ALLOW_ADMTLM
233 heimbach 1.2 write(cunit) cbuffindex
234     write(cunit) k
235 heimbach 1.14 cph#endif
236     write(cunit) (cbuff(ii), ii=1,cbuffindex)
237 heimbach 1.13 #endif
238 heimbach 1.2 endif
239 heimbach 1.9 c
240     if ( doPackDiag ) then
241     write(cunit2,rec=irectrue) globfldtmp2
242     write(cunit3,rec=irectrue) globfldtmp3
243     endif
244     c
245 heimbach 1.2 enddo
246     c
247     c -- end of irec loop --
248     enddo
249    
250 heimbach 1.9 if ( doPackDiag ) then
251     close ( cunit2 )
252     close ( cunit3 )
253     endif
254    
255 heimbach 1.2 _END_MASTER( mythid )
256    
257     return
258     end
259    

  ViewVC Help
Powered by ViewVC 1.1.22