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

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

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


Revision 1.23 - (show annotations) (download)
Thu Dec 23 02:43:43 2010 UTC (13 years, 4 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint62v, checkpoint62u, checkpoint62t, checkpoint62s, checkpoint62r, checkpoint62q, checkpoint62w
Changes since 1.22: +10 -10 lines
change arg. list of S/R MDSIO_PASS_R4/8toRL/S.

1 C $Header: /u/gcmpack/MITgcm/pkg/ctrl/ctrl_set_unpack_xyz.F,v 1.22 2010/09/26 02:51:38 gforget Exp $
2 C $Name: $
3
4 #include "CTRL_CPPOPTIONS.h"
5
6 subroutine ctrl_set_unpack_xyz(
7 & lxxadxx, cunit, ivartype, fname, masktype, weighttype,
8 & weightfld, nwetglobal, mythid)
9
10 c ==================================================================
11 c SUBROUTINE ctrl_set_unpack_xyz
12 c ==================================================================
13 c
14 c o Unpack the control vector such that land points are filled in.
15 c
16 c o Use a more precise nondimensionalization that depends on (x,y)
17 c Added weighttype to the argument list so that I can geographically
18 c vary the nondimensionalization.
19 c gebbie@mit.edu, 18-Mar-2003
20 c
21 c ==================================================================
22
23 implicit none
24
25 c == global variables ==
26
27 #include "EEPARAMS.h"
28 #include "SIZE.h"
29 #include "PARAMS.h"
30 #include "GRID.h"
31
32 #include "ctrl.h"
33 #include "optim.h"
34
35 c == routine arguments ==
36
37 logical lxxadxx
38 integer cunit
39 integer ivartype
40 character*( 80) fname
41 character*( 9) masktype
42 character*( 80) weighttype
43 _RL weightfld( nr,nsx,nsy )
44 integer nwetglobal(nr)
45 integer mythid
46
47 #ifndef EXCLUDE_CTRL_PACK
48 # ifndef ALLOW_PACKUNPACK_METHOD2
49 c == local variables ==
50
51 integer bi,bj
52 integer ip,jp
53 integer i,j,k
54 integer ii
55 integer il
56 integer irec
57 integer itlo,ithi
58 integer jtlo,jthi
59 integer jmin,jmax
60 integer imin,imax
61
62 integer cbuffindex
63
64 _RL globmsk ( snx,nsx,npx,sny,nsy,npy,nr )
65 _RL globfld3d( snx,nsx,npx,sny,nsy,npy,nr )
66 #ifdef CTRL_UNPACK_PRECISE
67 _RL weightfld3d( snx,nsx,npx,sny,nsy,npy,nr )
68 #endif
69 real*4 cbuff ( snx*nsx*npx*sny*nsy*npy )
70 real*4 globfldtmp2( snx,nsx,npx,sny,nsy,npy )
71 real*4 globfldtmp3( snx,nsx,npx,sny,nsy,npy )
72
73 character*(128) cfile
74 character*(80) weightname
75
76 _RL delZnorm
77 integer reclen, irectrue
78 integer cunit2, cunit3
79 character*(80) cfile2, cfile3
80
81 c == external ==
82
83 integer ilnblnk
84 external ilnblnk
85
86 cc == end of interface ==
87
88 jtlo = 1
89 jthi = nsy
90 itlo = 1
91 ithi = nsx
92 jmin = 1
93 jmax = sny
94 imin = 1
95 imax = snx
96
97 #ifdef CTRL_DELZNORM
98 delZnorm = 0.
99 do k = 1, Nr
100 delZnorm = delZnorm + delR(k)/FLOAT(Nr)
101 enddo
102 #endif
103
104 c Initialise temporary file
105 do k = 1,nr
106 do jp = 1,nPy
107 do bj = jtlo,jthi
108 do j = jmin,jmax
109 do ip = 1,nPx
110 do bi = itlo,ithi
111 do i = imin,imax
112 globfld3d (i,bi,ip,j,bj,jp,k) = 0. _d 0
113 globmsk (i,bi,ip,j,bj,jp,k) = 0. _d 0
114 globfldtmp2(i,bi,ip,j,bj,jp) = 0.
115 globfldtmp3(i,bi,ip,j,bj,jp) = 0.
116 enddo
117 enddo
118 enddo
119 enddo
120 enddo
121 enddo
122 enddo
123
124 c-- Only the master thread will do I/O.
125 _BEGIN_MASTER( mythid )
126
127 #ifdef CTRL_DELZNORM
128 do k = 1, nr
129 print *, 'ph-delznorm ', k, delZnorm, delR(k)
130 print *, 'ph-weight ', weightfld(k,1,1)
131 enddo
132 #endif
133
134 if ( doPackDiag ) then
135 write(cfile2(1:80),'(80a)') ' '
136 write(cfile3(1:80),'(80a)') ' '
137 if ( lxxadxx ) then
138 write(cfile2(1:80),'(a,I2.2,a,I4.4,a)')
139 & 'diag_unpack_nondim_ctrl_',
140 & ivartype, '_', optimcycle, '.bin'
141 write(cfile3(1:80),'(a,I2.2,a,I4.4,a)')
142 & 'diag_unpack_dimens_ctrl_',
143 & ivartype, '_', optimcycle, '.bin'
144 else
145 write(cfile2(1:80),'(a,I2.2,a,I4.4,a)')
146 & 'diag_unpack_nondim_grad_',
147 & ivartype, '_', optimcycle, '.bin'
148 write(cfile3(1:80),'(a,I2.2,a,I4.4,a)')
149 & 'diag_unpack_dimens_grad_',
150 & ivartype, '_', optimcycle, '.bin'
151 endif
152
153 reclen = FLOAT(snx*nsx*npx*sny*nsy*npy*4)
154 call mdsfindunit( cunit2, mythid )
155 open( cunit2, file=cfile2, status='unknown',
156 & access='direct', recl=reclen )
157 call mdsfindunit( cunit3, mythid )
158 open( cunit3, file=cfile3, status='unknown',
159 & access='direct', recl=reclen )
160 endif
161
162 #ifdef CTRL_UNPACK_PRECISE
163 il=ilnblnk( weighttype)
164 write(weightname(1:80),'(80a)') ' '
165 write(weightname(1:80),'(a)') weighttype(1:il)
166
167 call MDSREADFIELD_3D_GL(
168 & weightname, ctrlprec, 'RL',
169 & Nr, weightfld3d, 1, mythid)
170 #endif
171
172 call MDSREADFIELD_3D_GL(
173 & masktype, ctrlprec, 'RL',
174 & Nr, globmsk, 1, mythid)
175
176 do irec = 1, ncvarrecs(ivartype)
177 #ifndef ALLOW_ADMTLM
178 read(cunit) filencvarindex(ivartype)
179 if (filencvarindex(ivartype) .NE. ncvarindex(ivartype))
180 & then
181 print *, 'ctrl_set_unpack_xyz:WARNING: wrong ncvarindex ',
182 & filencvarindex(ivartype), ncvarindex(ivartype)
183 STOP 'in S/R ctrl_unpack'
184 endif
185 read(cunit) filej
186 read(cunit) filei
187 #endif /* ALLOW_ADMTLM */
188 do k = 1, Nr
189 irectrue = (irec-1)*nr + k
190 if ( doZscaleUnpack ) then
191 delZnorm = (delR(1)/delR(k))**delZexp
192 else
193 delZnorm = 1. _d 0
194 endif
195 cbuffindex = nwetglobal(k)
196 if ( cbuffindex .gt. 0 ) then
197 #ifndef ALLOW_ADMTLM
198 read(cunit) filencbuffindex
199 if (filencbuffindex .NE. cbuffindex) then
200 print *, 'WARNING: wrong cbuffindex ',
201 & filencbuffindex, cbuffindex
202 STOP 'in S/R ctrl_unpack'
203 endif
204 read(cunit) filek
205 if (filek .NE. k) then
206 print *, 'WARNING: wrong k ',
207 & filek, k
208 STOP 'in S/R ctrl_unpack'
209 endif
210 cph#endif /* ALLOW_ADMTLM */
211 read(cunit) (cbuff(ii), ii=1,cbuffindex)
212 #endif /* ALLOW_ADMTLM */
213 endif
214 c
215 cbuffindex = 0
216 do jp = 1,nPy
217 do bj = jtlo,jthi
218 do j = jmin,jmax
219 do ip = 1,nPx
220 do bi = itlo,ithi
221 do i = imin,imax
222 if ( globmsk(i,bi,ip,j,bj,jp,k) .ne. 0. ) then
223 cbuffindex = cbuffindex + 1
224 globfld3d(i,bi,ip,j,bj,jp,k) = cbuff(cbuffindex)
225 cph(
226 globfldtmp2(i,bi,ip,j,bj,jp) = cbuff(cbuffindex)
227 cph)
228 #ifdef ALLOW_ADMTLM
229 nveccount = nveccount + 1
230 globfld3d(i,bi,ip,j,bj,jp,k) =
231 & phtmpadmtlm(nveccount)
232 cph(
233 globfldtmp2(i,bi,ip,j,bj,jp) =
234 & phtmpadmtlm(nveccount)
235 cph)
236 #endif
237 #ifndef ALLOW_SMOOTH_CORREL3D
238 #ifdef ALLOW_NONDIMENSIONAL_CONTROL_IO
239 if ( lxxadxx ) then
240 globfld3d(i,bi,ip,j,bj,jp,k) = delZnorm
241 & * globfld3d(i,bi,ip,j,bj,jp,k)
242 # ifdef CTRL_UNPACK_PRECISE
243 & / sqrt(weightfld3d(i,bi,ip,j,bj,jp,k))
244 # else
245 & / sqrt(weightfld(k,bi,bj))
246 # endif
247 else
248 globfld3d(i,bi,ip,j,bj,jp,k) = 1/delZnorm
249 & * globfld3d(i,bi,ip,j,bj,jp,k)
250 # ifdef CTRL_UNPACK_PRECISE
251 & * sqrt(weightfld3d(i,bi,ip,j,bj,jp,k))
252 # else
253 & * sqrt(weightfld(k,bi,bj))
254 # endif
255 endif
256 #endif /* ALLOW_NONDIMENSIONAL_CONTROL_IO */
257 #endif /* ALLOW_SMOOTH_CORREL3D */
258 else
259 globfld3d(i,bi,ip,j,bj,jp,k) = 0. _d 0
260 endif
261 cph(
262 globfldtmp3(i,bi,ip,j,bj,jp) =
263 & globfld3d(i,bi,ip,j,bj,jp,k)
264 cph)
265 enddo
266 enddo
267 enddo
268 enddo
269 enddo
270 enddo
271 c
272 if ( doPackDiag ) then
273 write(cunit2,rec=irectrue) globfldtmp2
274 write(cunit3,rec=irectrue) globfldtmp3
275 endif
276 c
277 enddo
278
279 call MDSWRITEFIELD_3D_GL( fname, ctrlprec, 'RL',
280 & Nr, globfld3d,
281 & irec, optimcycle, mythid)
282
283 enddo
284
285 if ( doPackDiag ) then
286 close ( cunit2 )
287 close ( cunit3 )
288 endif
289
290 _END_MASTER( mythid )
291
292 # else
293 c == local variables ==
294
295 integer bi,bj
296 integer ip,jp
297 integer i,j,k
298 integer ii
299 integer il
300 integer irec
301 integer itlo,ithi
302 integer jtlo,jthi
303
304 integer cbuffindex
305
306 _RL msk3d(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
307 real*8 msk2d_buf(sNx,sNy,nSx,nSy)
308 real*8 msk2d_buf_glo(Nx,Ny)
309
310 _RL fld3d(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
311 real*8 fld2d_buf(sNx,sNy,nSx,nSy)
312 real*8 fld2d_buf_glo(Nx,Ny)
313
314 _RL fld3dDim(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
315 _RL fld3dNodim(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
316
317 #ifdef CTRL_PACK_PRECISE
318 _RL wei3d(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
319 #endif
320
321 real*4 cbuff ( snx*nsx*npx*sny*nsy*npy )
322
323 character*(80) weightname
324 _RL delZnorm
325 character*(80) cfile2, cfile3
326
327 c == external ==
328
329 integer ilnblnk
330 external ilnblnk
331
332 c == end of interface ==
333
334 c-- part 1: preliminary reads and definitions
335
336 #ifdef CTRL_PACK_PRECISE
337 call active_read_xyz(weighttype, wei3d, 1,
338 & .FALSE., .FALSE., 0 , mythid, 1)
339 #endif
340
341 call active_read_xyz(masktype, msk3d, 1,
342 & .FALSE., .FALSE., 0 , mythid, 1)
343
344 if ( doPackDiag ) then
345 write(cfile2(1:80),'(80a)') ' '
346 write(cfile3(1:80),'(80a)') ' '
347 il = ilnblnk( fname )
348 if ( lxxadxx ) then
349 write(cfile2(1:80),'(2a)') fname(1:il),'.unpack_ctrl_adim'
350 write(cfile3(1:80),'(2a)') fname(1:il),'.unpack_ctrl_dim'
351 else
352 write(cfile2(1:80),'(2a)') fname(1:il),'.unpack_grad_adim'
353 write(cfile3(1:80),'(2a)') fname(1:il),'.unpack_grad_dim'
354 endif
355 endif
356
357 c-- part 2: loop over records
358
359 do irec = 1, ncvarrecs(ivartype)
360
361 c-- 2.1: array <- buffer <- global buffer <- global file
362
363 #ifndef ALLOW_ADMTLM
364 _BEGIN_MASTER( mythid )
365 IF ( myProcId .eq. 0 ) THEN
366 read(cunit) filencvarindex(ivartype)
367 if (filencvarindex(ivartype) .NE. ncvarindex(ivartype))
368 & then
369 print *, 'ctrl_set_unpack_xyz:WARNING: wrong ncvarindex ',
370 & filencvarindex(ivartype), ncvarindex(ivartype)
371 STOP 'in S/R ctrl_unpack'
372 endif
373 read(cunit) filej
374 read(cunit) filei
375 ENDIF
376 _END_MASTER( mythid )
377 _BARRIER
378 #endif /* ALLOW_ADMTLM */
379
380 do k = 1, nr
381
382 CALL MDS_PASS_R8toRL( msk2d_buf, msk3d,
383 & 0, 0, 1, k, Nr, 0, 0, .FALSE., myThid )
384 CALL BAR2( myThid )
385 CALL GATHER_2D_R8( msk2d_buf_glo, msk2d_buf,
386 & Nx,Ny,.FALSE.,.TRUE.,myThid)
387 CALL BAR2( myThid )
388
389 _BEGIN_MASTER( mythid )
390 cbuffindex = nwetglobal(k)
391 IF ( myProcId .eq. 0 ) THEN
392
393 #ifndef ALLOW_ADMTLM
394 if ( cbuffindex .gt. 0) then
395 read(cunit) filencbuffindex
396 read(cunit) filek
397 if (filencbuffindex .NE. cbuffindex) then
398 print *, 'WARNING: wrong cbuffindex ',
399 & filencbuffindex, cbuffindex
400 STOP 'in S/R ctrl_unpack'
401 endif
402 if (filek .NE. k) then
403 print *, 'WARNING: wrong k ',
404 & filek, k
405 STOP 'in S/R ctrl_unpack'
406 endif
407 read(cunit) (cbuff(ii), ii=1,cbuffindex)
408 endif
409 #endif
410
411 cbuffindex = 0
412 DO j=1,Ny
413 DO i=1,Nx
414 if (msk2d_buf_glo(i,j) .ne. 0. ) then
415 cbuffindex = cbuffindex + 1
416 fld2d_buf_glo(i,j) = cbuff(cbuffindex)
417 endif
418 ENDDO
419 ENDDO
420
421 ENDIF
422 _END_MASTER( mythid )
423 _BARRIER
424
425 CALL BAR2( myThid )
426 CALL SCATTER_2D_R8( fld2d_buf_glo, fld2d_buf,
427 & Nx,Ny,.FALSE.,.TRUE.,myThid)
428 CALL BAR2( myThid )
429 CALL MDS_PASS_R8toRL( fld2d_buf, fld3dNodim,
430 & 0, 0, 1, k, Nr, 0, 0, .TRUE., myThid )
431
432 enddo !do k = 1, nr
433
434 c-- 2.2: normalize field if needed
435 DO bj = myByLo(myThid), myByHi(myThid)
436 DO bi = myBxLo(myThid), myBxHi(myThid)
437 DO k=1,Nr
438 if ( doZscalePack ) then
439 delZnorm = (delR(1)/delR(k))**delZexp
440 else
441 delZnorm = 1. _d 0
442 endif
443 DO j=1,sNy
444 DO i=1,sNx
445 if (msk3d(i,j,k,bi,bj).EQ.0. _d 0) then
446 fld3dDim(i,j,k,bi,bj)=0. _d 0
447 fld3dNodim(i,j,k,bi,bj)=0. _d 0
448 else
449 #ifdef ALLOW_ADMTLM
450 nveccount = nveccount + 1
451 fld3dNodim(i,j,k,bi,bj)=phtmpadmtlm(nveccount)
452 #endif
453 #ifdef ALLOW_SMOOTH_CORREL3D
454 fld3dDim(i,j,k,bi,bj)=fld3dNodim(i,j,k,bi,bj)
455 #else
456 # ifndef ALLOW_NONDIMENSIONAL_CONTROL_IO
457 fld3dDim(i,j,k,bi,bj) = fld3dNodim(i,j,k,bi,bj)
458 # else
459 if (lxxadxx) then
460 fld3dDim(i,j,k,bi,bj) =
461 & fld3dNodim(i,j,k,bi,bj) * delZnorm
462 # ifdef CTRL_PACK_PRECISE
463 & / sqrt(wei3d(i,j,k,bi,bj))
464 # else
465 & / sqrt(weightfld(k,bi,bj))
466 # endif
467 else
468 fld3dDim(i,j,k,bi,bj) =
469 & fld3dNodim(i,j,k,bi,bj) / delZnorm
470 # ifdef CTRL_PACK_PRECISE
471 & * sqrt(wei3d(i,j,k,bi,bj))
472 # else
473 & * sqrt(weightfld(k,bi,bj))
474 # endif
475 endif
476 # endif /* ALLOW_NONDIMENSIONAL_CONTROL_IO */
477 #endif /* ALLOW_SMOOTH_CORREL3D */
478 endif
479 ENDDO
480 ENDDO
481 ENDDO
482 ENDDO
483 ENDDO
484
485 c-- 2.3:
486 if ( doPackDiag ) then
487 c error: twice the same one
488 call WRITE_REC_3D_RL( cfile2, ctrlprec,
489 & Nr, fld3dNodim, irec, 0, mythid)
490 call WRITE_REC_3D_RL( cfile3, ctrlprec,
491 & Nr, fld3dDim, irec, 0, mythid)
492 endif
493
494 c-- 2.4:
495 call WRITE_REC_3D_RL( fname, ctrlprec,
496 & Nr, fld3dDim, irec, 0, mythid)
497
498
499 enddo !do irec = 1, ncvarrecs(ivartype)
500
501 # endif /* ALLOW_PACKUNPACK_METHOD2 */
502 # endif /* EXCLUDE_CTRL_PACK */
503
504 return
505 end
506

  ViewVC Help
Powered by ViewVC 1.1.22