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

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

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


Revision 1.24 - (show annotations) (download)
Thu Dec 23 02:43:43 2010 UTC (13 years, 5 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint62v, checkpoint62u, checkpoint62t, checkpoint62s, checkpoint62r, checkpoint62q, checkpoint62w
Changes since 1.23: +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_xy.F,v 1.23 2010/09/26 02:51:38 gforget Exp $
2 C $Name: $
3
4 #include "CTRL_CPPOPTIONS.h"
5
6 subroutine ctrl_set_unpack_xy(
7 & lxxadxx, cunit, ivartype, fname, masktype, weighttype,
8 & nwetglobal, mythid)
9
10 c ==================================================================
11 c SUBROUTINE ctrl_set_unpack_xy
12 c ==================================================================
13 c
14 c o Unpack the control vector such that the land points are filled
15 c in.
16 c
17 c changed: heimbach@mit.edu 17-Jun-2003
18 c merged changes from Armin to replace write of
19 c nr * globfld2d by 1 * globfld3d
20 c (ad hoc fix to speed up global I/O)
21 c
22 c ==================================================================
23
24 implicit none
25
26 c == global variables ==
27
28 #include "EEPARAMS.h"
29 #include "SIZE.h"
30 #include "PARAMS.h"
31 #include "GRID.h"
32
33 #include "ctrl.h"
34 #include "optim.h"
35
36 c == routine arguments ==
37
38 logical lxxadxx
39 integer cunit
40 integer ivartype
41 character*( 80) fname, fnameGlobal
42 character*( 9) masktype
43 character*( 80) weighttype
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,nrec_nl
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 _RL globfld2d( snx,nsx,npx,sny,nsy,npy )
67 real*4 cbuff ( snx*nsx*npx*sny*nsy*npy )
68
69 character*(128) cfile
70 character*( 80) weightname
71
72 integer reclen,irectrue
73 integer cunit2, cunit3
74 character*(80) cfile2, cfile3
75 real*4 globfldtmp2( snx,nsx,npx,sny,nsy,npy )
76 real*4 globfldtmp3( snx,nsx,npx,sny,nsy,npy )
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 nbuffGlobal = nbuffGlobal + 1
95
96 c Initialise temporary file
97 do k = 1,nr
98 do jp = 1,nPy
99 do bj = jtlo,jthi
100 do j = jmin,jmax
101 do ip = 1,nPx
102 do bi = itlo,ithi
103 do i = imin,imax
104 globfld3d (i,bi,ip,j,bj,jp,k) = 0. _d 0
105 globmsk (i,bi,ip,j,bj,jp,k) = 0. _d 0
106 globfldtmp2(i,bi,ip,j,bj,jp) = 0.
107 globfldtmp3(i,bi,ip,j,bj,jp) = 0.
108 enddo
109 enddo
110 enddo
111 enddo
112 enddo
113 enddo
114 enddo
115
116 c-- Only the master thread will do I/O.
117 _BEGIN_MASTER( mythid )
118
119 if ( doPackDiag ) then
120 write(cfile2(1:80),'(80a)') ' '
121 write(cfile3(1:80),'(80a)') ' '
122 if ( lxxadxx ) then
123 write(cfile2(1:80),'(a,I2.2,a,I4.4,a)')
124 & 'diag_unpack_nondim_ctrl_',
125 & ivartype, '_', optimcycle, '.bin'
126 write(cfile3(1:80),'(a,I2.2,a,I4.4,a)')
127 & 'diag_unpack_dimens_ctrl_',
128 & ivartype, '_', optimcycle, '.bin'
129 else
130 write(cfile2(1:80),'(a,I2.2,a,I4.4,a)')
131 & 'diag_unpack_nondim_grad_',
132 & ivartype, '_', optimcycle, '.bin'
133 write(cfile3(1:80),'(a,I2.2,a,I4.4,a)')
134 & 'diag_unpack_dimens_grad_',
135 & ivartype, '_', optimcycle, '.bin'
136 endif
137
138 reclen = FLOAT(snx*nsx*npx*sny*nsy*npy*4)
139 call mdsfindunit( cunit2, mythid )
140 open( cunit2, file=cfile2, status='unknown',
141 & access='direct', recl=reclen )
142 call mdsfindunit( cunit3, mythid )
143 open( cunit3, file=cfile3, status='unknown',
144 & access='direct', recl=reclen )
145 endif
146
147 #ifdef ALLOW_NONDIMENSIONAL_CONTROL_IO
148 il=ilnblnk( weighttype)
149 write(weightname(1:80),'(80a)') ' '
150 write(weightname(1:80),'(a)') weighttype(1:il)
151 call MDSREADFIELD_2D_GL(
152 & weightname, ctrlprec, 'RL',
153 & 1, globfld2d, 1, mythid)
154 #endif
155
156 call MDSREADFIELD_3D_GL(
157 & masktype, ctrlprec, 'RL',
158 & Nr, globmsk, 1, mythid)
159
160 nrec_nl=int(ncvarrecs(ivartype)/Nr)
161 do irec = 1, nrec_nl
162 do k = 1,Nr
163 irectrue = (irec-1)*nr + k
164 #ifndef ALLOW_ADMTLM
165 read(cunit) filencvarindex(ivartype)
166 if (filencvarindex(ivartype) .NE. ncvarindex(ivartype))
167 & then
168 print *, 'ctrl_set_unpack_xy:WARNING: wrong ncvarindex ',
169 & filencvarindex(ivartype), ncvarindex(ivartype)
170 STOP 'in S/R ctrl_unpack'
171 endif
172 read(cunit) filej
173 read(cunit) filei
174 #endif /* ndef ALLOW_ADMTLM */
175 cbuffindex = nwetglobal(1)
176 if ( cbuffindex .gt. 0 ) then
177 #ifndef ALLOW_ADMTLM
178 read(cunit) filencbuffindex
179 if (filencbuffindex .NE. cbuffindex) then
180 print *, 'WARNING: wrong cbuffindex ',
181 & filencbuffindex, cbuffindex
182 STOP 'in S/R ctrl_unpack'
183 endif
184 read(cunit) filek
185 if (filek .NE. 1) then
186 print *, 'WARNING: wrong k ',
187 & filek, 1
188 STOP 'in S/R ctrl_unpack'
189 endif
190 cph#endif /* ndef ALLOW_ADMTLM */
191 read(cunit) (cbuff(ii), ii=1,cbuffindex)
192 #endif /* ndef ALLOW_ADMTLM */
193 endif
194 c
195 cbuffindex = 0
196 do jp = 1,nPy
197 do bj = jtlo,jthi
198 do j = jmin,jmax
199 do ip = 1,nPx
200 do bi = itlo,ithi
201 do i = imin,imax
202 if ( globmsk(i,bi,ip,j,bj,jp,1) .ne. 0. ) then
203 cbuffindex = cbuffindex + 1
204 globfld3d(i,bi,ip,j,bj,jp,k) = cbuff(cbuffindex)
205 cph(
206 globfldtmp2(i,bi,ip,j,bj,jp) = cbuff(cbuffindex)
207 cph)
208 #ifdef ALLOW_ADMTLM
209 nveccount = nveccount + 1
210 globfld3d(i,bi,ip,j,bj,jp,k) =
211 & phtmpadmtlm(nveccount)
212 cph(
213 globfldtmp2(i,bi,ip,j,bj,jp) =
214 & phtmpadmtlm(nveccount)
215 cph)
216 #endif
217 #ifndef ALLOW_SMOOTH_CORREL2D
218 #ifdef ALLOW_NONDIMENSIONAL_CONTROL_IO
219 if ( lxxadxx ) then
220 globfld3d(i,bi,ip,j,bj,jp,k) =
221 & globfld3d(i,bi,ip,j,bj,jp,k)/
222 & sqrt(globfld2d(i,bi,ip,j,bj,jp))
223 & * forcingPrecond
224 else
225 globfld3d(i,bi,ip,j,bj,jp,k) =
226 & globfld3d(i,bi,ip,j,bj,jp,k)*
227 & sqrt(globfld2d(i,bi,ip,j,bj,jp))
228 & / forcingPrecond
229 endif
230 #endif
231 #else /* ALLOW_SMOOTH_CORREL2D */
232 if ( lxxadxx ) then
233 globfld3d(i,bi,ip,j,bj,jp,k) =
234 & globfld3d(i,bi,ip,j,bj,jp,k)
235 & * forcingPrecond
236 else
237 globfld3d(i,bi,ip,j,bj,jp,k) =
238 & globfld3d(i,bi,ip,j,bj,jp,k)
239 & / forcingPrecond
240 endif
241 #endif /* ALLOW_SMOOTH_CORREL2D */
242 else
243 globfld3d(i,bi,ip,j,bj,jp,k) = 0. _d 0
244 endif
245 cph(
246 globfldtmp3(i,bi,ip,j,bj,jp) =
247 & globfld3d(i,bi,ip,j,bj,jp,k)
248 cph)
249 enddo
250 enddo
251 enddo
252 enddo
253 enddo
254 enddo
255 cph(
256 if ( doPackDiag ) then
257 write(cunit2,rec=irectrue) globfldtmp2
258 write(cunit3,rec=irectrue) globfldtmp3
259 endif
260 cph)
261 enddo
262
263 call MDSWRITEFIELD_3D_GL( fname, ctrlprec, 'RL',
264 & NR, globfld3d,
265 & irec, optimcycle, mythid)
266
267 enddo
268
269 do irec = nrec_nl*Nr+1,ncvarrecs(ivartype)
270 #ifndef ALLOW_ADMTLM
271 read(cunit) filencvarindex(ivartype)
272 if (filencvarindex(ivartype) .NE. ncvarindex(ivartype))
273 & then
274 print *, 'ctrl_set_unpack_xy:WARNING: wrong ncvarindex ',
275 & filencvarindex(ivartype), ncvarindex(ivartype)
276 STOP 'in S/R ctrl_unpack'
277 endif
278 read(cunit) filej
279 read(cunit) filei
280 #endif /* ALLOW_ADMTLM */
281 do k = 1,1
282 irectrue = irec
283 cbuffindex = nwetglobal(k)
284 if ( cbuffindex .gt. 0 ) then
285 #ifndef ALLOW_ADMTLM
286 read(cunit) filencbuffindex
287 if (filencbuffindex .NE. cbuffindex) then
288 print *, 'WARNING: wrong cbuffindex ',
289 & filencbuffindex, cbuffindex
290 STOP 'in S/R ctrl_unpack'
291 endif
292 read(cunit) filek
293 if (filek .NE. k) then
294 print *, 'WARNING: wrong k ',
295 & filek, k
296 STOP 'in S/R ctrl_unpack'
297 endif
298 cph#endif /* ALLOW_ADMTLM */
299 read(cunit) (cbuff(ii), ii=1,cbuffindex)
300 #endif /* ALLOW_ADMTLM */
301 endif
302 c
303 cbuffindex = 0
304 do jp = 1,nPy
305 do bj = jtlo,jthi
306 do j = jmin,jmax
307 do ip = 1,nPx
308 do bi = itlo,ithi
309 do i = imin,imax
310 if ( globmsk(i,bi,ip,j,bj,jp,k) .ne. 0. ) then
311 cbuffindex = cbuffindex + 1
312 globfld3d(i,bi,ip,j,bj,jp,k) = cbuff(cbuffindex)
313 cph(
314 globfldtmp2(i,bi,ip,j,bj,jp) = cbuff(cbuffindex)
315 cph)
316 #ifdef ALLOW_ADMTLM
317 nveccount = nveccount + 1
318 globfld3d(i,bi,ip,j,bj,jp,k) =
319 & phtmpadmtlm(nveccount)
320 cph(
321 globfldtmp2(i,bi,ip,j,bj,jp) =
322 & phtmpadmtlm(nveccount)
323 cph)
324 #endif
325 #ifndef ALLOW_SMOOTH_CORREL2D
326 #ifdef ALLOW_NONDIMENSIONAL_CONTROL_IO
327 if ( lxxadxx ) then
328 globfld3d(i,bi,ip,j,bj,jp,k) =
329 & globfld3d(i,bi,ip,j,bj,jp,k)/
330 & sqrt(globfld2d(i,bi,ip,j,bj,jp))
331 & * forcingPrecond
332 else
333 globfld3d(i,bi,ip,j,bj,jp,k) =
334 & globfld3d(i,bi,ip,j,bj,jp,k)*
335 & sqrt(globfld2d(i,bi,ip,j,bj,jp))
336 & / forcingPrecond
337 endif
338 #endif
339 #else /* ALLOW_SMOOTH_CORREL2D */
340 if ( lxxadxx ) then
341 globfld3d(i,bi,ip,j,bj,jp,k) =
342 & globfld3d(i,bi,ip,j,bj,jp,k)
343 & * forcingPrecond
344 else
345 globfld3d(i,bi,ip,j,bj,jp,k) =
346 & globfld3d(i,bi,ip,j,bj,jp,k)
347 & / forcingPrecond
348 endif
349 #endif /* ALLOW_SMOOTH_CORREL2D */
350 else
351 globfld3d(i,bi,ip,j,bj,jp,k) = 0. _d 0
352 endif
353 cph(
354 globfldtmp3(i,bi,ip,j,bj,jp) =
355 & globfld3d(i,bi,ip,j,bj,jp,k)
356 cph)
357 enddo
358 enddo
359 enddo
360 enddo
361 enddo
362 enddo
363 cph(
364 if ( doPackDiag ) then
365 write(cunit2,rec=irectrue) globfldtmp2
366 write(cunit3,rec=irectrue) globfldtmp3
367 endif
368 cph)
369 enddo
370
371 call MDSWRITEFIELD_2D_GL( fname, ctrlprec, 'RL',
372 & 1, globfld3d(1,1,1,1,1,1,1),
373 & irec, optimcycle, mythid)
374
375 enddo
376
377 if ( doPackDiag ) then
378 close ( cunit2 )
379 close ( cunit3 )
380 endif
381
382 _END_MASTER( mythid )
383
384 # else
385 c == local variables ==
386
387 integer bi,bj
388 integer ip,jp
389 integer i,j,k
390 integer ii
391 integer il
392 integer irec
393 integer itlo,ithi
394 integer jtlo,jthi
395
396 integer cbuffindex
397
398 _RL msk3d(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr,nSx,nSy)
399 real*8 msk2d_buf(sNx,sNy,nSx,nSy)
400 real*8 msk2d_buf_glo(Nx,Ny)
401
402 _RL fld2d(1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
403 real*8 fld2d_buf(sNx,sNy,nSx,nSy)
404 real*8 fld2d_buf_glo(Nx,Ny)
405
406 _RL fld2dDim(1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
407 _RL fld2dNodim(1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
408
409 _RL wei2d(1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
410
411 real*4 cbuff ( snx*nsx*npx*sny*nsy*npy )
412
413 character*(80) weightname
414 _RL delZnorm
415 character*(80) cfile2, cfile3
416
417 c == external ==
418
419 integer ilnblnk
420 external ilnblnk
421
422 c == end of interface ==
423
424 c-- part 1: preliminary reads and definitions
425
426 call active_read_xy(weighttype, wei2d, 1,
427 & .FALSE., .FALSE., 0 , mythid, 1)
428
429 call active_read_xyz(masktype, msk3d, 1,
430 & .FALSE., .FALSE., 0 , mythid, 1)
431
432 if ( doPackDiag ) then
433 write(cfile2(1:80),'(80a)') ' '
434 write(cfile3(1:80),'(80a)') ' '
435 il = ilnblnk( fname )
436 if ( lxxadxx ) then
437 write(cfile2(1:80),'(2a)') fname(1:il),'.unpack_ctrl_adim'
438 write(cfile3(1:80),'(2a)') fname(1:il),'.unpack_ctrl_dim'
439 else
440 write(cfile2(1:80),'(2a)') fname(1:il),'.unpack_grad_adim'
441 write(cfile3(1:80),'(2a)') fname(1:il),'.unpack_grad_dim'
442 endif
443 endif
444
445 c-- part 2: loop over records
446
447 do irec = 1, ncvarrecs(ivartype)
448
449 c-- 2.1: array <- buffer <- global buffer <- global file
450
451 #ifndef ALLOW_ADMTLM
452 _BEGIN_MASTER( mythid )
453 IF ( myProcId .eq. 0 ) THEN
454 read(cunit) filencvarindex(ivartype)
455 if (filencvarindex(ivartype) .NE. ncvarindex(ivartype))
456 & then
457 print *, 'ctrl_set_unpack_xy:WARNING: wrong ncvarindex ',
458 & filencvarindex(ivartype), ncvarindex(ivartype)
459 STOP 'in S/R ctrl_unpack'
460 endif
461 read(cunit) filej
462 read(cunit) filei
463 ENDIF
464 _END_MASTER( mythid )
465 _BARRIER
466 #endif /* ALLOW_ADMTLM */
467
468 do k = 1, 1
469
470 CALL MDS_PASS_R8toRL( msk2d_buf, msk3d,
471 & 0, 0, 1, k, Nr, 0, 0, .FALSE., myThid )
472 CALL BAR2( myThid )
473 CALL GATHER_2D_R8( msk2d_buf_glo, msk2d_buf,
474 & Nx,Ny,.FALSE.,.TRUE.,myThid)
475 CALL BAR2( myThid )
476
477 _BEGIN_MASTER( mythid )
478 cbuffindex = nwetglobal(k)
479 IF ( myProcId .eq. 0 ) THEN
480
481 #ifndef ALLOW_ADMTLM
482 if ( cbuffindex .gt. 0) then
483 read(cunit) filencbuffindex
484 read(cunit) filek
485 if (filencbuffindex .NE. cbuffindex) then
486 print *, 'WARNING: wrong cbuffindex ',
487 & filencbuffindex, cbuffindex
488 STOP 'in S/R ctrl_unpack'
489 endif
490 if (filek .NE. 1) then
491 print *, 'WARNING: wrong k ',
492 & filek, 1
493 STOP 'in S/R ctrl_unpack'
494 endif
495 read(cunit) (cbuff(ii), ii=1,cbuffindex)
496 endif
497 #endif
498
499 cbuffindex = 0
500 DO j=1,Ny
501 DO i=1,Nx
502 if (msk2d_buf_glo(i,j) .ne. 0. ) then
503 cbuffindex = cbuffindex + 1
504 fld2d_buf_glo(i,j) = cbuff(cbuffindex)
505 endif
506 ENDDO
507 ENDDO
508
509 ENDIF
510 _END_MASTER( mythid )
511 _BARRIER
512
513 CALL BAR2( myThid )
514 CALL SCATTER_2D_R8( fld2d_buf_glo, fld2d_buf,
515 & Nx,Ny,.FALSE.,.TRUE.,myThid)
516 CALL BAR2( myThid )
517 CALL MDS_PASS_R8toRL( fld2d_buf, fld2dNodim,
518 & 0, 0, 1, k, 1, 0, 0, .TRUE., myThid )
519
520 enddo !do k = 1, 1
521
522
523 c-- 2.2: normalize field if needed
524 DO bj = myByLo(myThid), myByHi(myThid)
525 DO bi = myBxLo(myThid), myBxHi(myThid)
526 DO j=1,sNy
527 DO i=1,sNx
528 if (msk3d(i,j,1,bi,bj).EQ.0. _d 0) then
529 fld2dDim(i,j,bi,bj)=0. _d 0
530 fld2dNodim(i,j,bi,bj)=0. _d 0
531 else
532 #ifdef ALLOW_ADMTLM
533 nveccount = nveccount + 1
534 fld2dNodim(i,j,bi,bj)=phtmpadmtlm(nveccount)
535 #endif
536 #ifdef ALLOW_SMOOTH_CORREL2D
537 if (lxxadxx) then
538 fld2dDim(i,j,bi,bj) =
539 & fld2dNodim(i,j,bi,bj) * forcingPrecond
540 else
541 fld2dDim(i,j,bi,bj) =
542 & fld2dNodim(i,j,bi,bj) / forcingPrecond
543 endif
544 #else
545 # ifndef ALLOW_NONDIMENSIONAL_CONTROL_IO
546 fld2dDim(i,j,bi,bj) = fld2dNodim(i,j,bi,bj)
547 # else
548 if (lxxadxx) then
549 fld2dDim(i,j,bi,bj) =
550 & fld2dNodim(i,j,bi,bj) / sqrt(wei2d(i,j,bi,bj))
551 & * forcingPrecond
552 else
553 fld2dDim(i,j,bi,bj) =
554 & fld2dNodim(i,j,bi,bj) * sqrt(wei2d(i,j,bi,bj))
555 & / forcingPrecond
556 endif
557 # endif /* ALLOW_NONDIMENSIONAL_CONTROL_IO */
558 #endif /* ALLOW_SMOOTH_CORREL2D */
559 endif
560 ENDDO
561 ENDDO
562 ENDDO
563 ENDDO
564
565 c-- 2.3:
566 if ( doPackDiag ) then
567 call WRITE_REC_3D_RL( cfile2, ctrlprec,
568 & 1, fld2dNodim, irec, 0, mythid)
569 call WRITE_REC_3D_RL( cfile3, ctrlprec,
570 & 1, fld2dDim, irec, 0, mythid)
571 endif
572
573 c-- 2.4:
574 call WRITE_REC_3D_RL( fname, ctrlprec,
575 & 1, fld2dDim, irec, 0, mythid)
576
577 enddo !do irec = 1, ncvarrecs(ivartype)
578
579 # endif /* ALLOW_PACKUNPACK_METHOD2 */
580 # endif /* EXCLUDE_CTRL_PACK */
581
582 return
583 end
584

  ViewVC Help
Powered by ViewVC 1.1.22