/[MITgcm]/MITgcm/pkg/mdsio/mdsio_read_field.F
ViewVC logotype

Contents of /MITgcm/pkg/mdsio/mdsio_read_field.F

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


Revision 1.3 - (show annotations) (download)
Tue Nov 13 19:37:44 2007 UTC (16 years, 6 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint60, checkpoint61, checkpoint59q, checkpoint59p, checkpoint59r, checkpoint59m, checkpoint59l, checkpoint59o, checkpoint59n, checkpoint59k, checkpoint61f, checkpoint61g, checkpoint61d, checkpoint61e, checkpoint61b, checkpoint61c, checkpoint61a
Changes since 1.2: +49 -33 lines
add arguments to S/R MDS_READ_FIELD and MDS_WRITE_FIELD.

1 C $Header: /u/gcmpack/MITgcm/pkg/mdsio/mdsio_read_field.F,v 1.2 2007/03/19 02:30:49 jmc Exp $
2 C $Name: $
3
4 #include "MDSIO_OPTIONS.h"
5
6 CBOP
7 C !ROUTINE: MDS_READ_FIELD
8 C !INTERFACE:
9 SUBROUTINE MDS_READ_FIELD(
10 I fName,
11 I filePrec,
12 I useCurrentDir,
13 I arrType,
14 I kSize,kLo,kHi,
15 O arr,
16 I irecord,
17 I myThid )
18
19 C !DESCRIPTION:
20 C Arguments:
21 C
22 C fName (string) :: base name for file to read
23 C filePrec (integer) :: number of bits per word in file (32 or 64)
24 C useCurrentDir(logic):: always read from the current directory (even if
25 C "mdsioLocalDir" is set)
26 C arrType (char(2)) :: declaration of "arr": either "RS" or "RL"
27 C kSize (integer) :: size of third dimension: normally either 1 or Nr
28 C kLo (integer) :: 1rst vertical level (of array "arr") to read-in
29 C kHi (integer) :: last vertical level (of array "arr") to read-in
30 C arr ( RS/RL ) :: array to read into, arr(:,:,kSize,:,:)
31 C irecord (integer) :: record number to read
32 C myIter (integer) :: time step number
33 C myThid (integer) :: thread identifier
34 C
35 C MDS_READ_FIELD first checks to see IF the file "fName" exists, then
36 C IF the file "fName.data" exists and finally the tiled files of the
37 C form "fName.xxx.yyy.data" exist. Currently, the meta-files are not
38 C read because it is difficult to parse files in fortran.
39 C The precision of the file is decsribed by filePrec, set either
40 C to floatPrec32 or floatPrec64. The precision or declaration of
41 C the array argument must be consistently described by the char*(2)
42 C string arrType, either "RS" or "RL".
43 C (kSize,kLo,kHi) allows for both 2-D and 3-D arrays to be handled, with
44 C the option to only read and fill-in a sub-set of consecutive vertical
45 C levels (from kLo to kHi) ; (kSize,kLo,kHi)=(1,1,1) implies a 2-D model
46 C field and (kSize,kLo,kHi)=(Nr,1,Nr) implies a 3-D model field.
47 C irecord is the record number to be read and must be >= 1.
48 C The file data is stored in arr *but* the overlaps are *not* updated,
49 C i.e., an exchange must be called. This is because the routine is
50 C sometimes called from within a MASTER_THID region.
51 C
52 C Created: 03/16/99 adcroft@mit.edu
53 CEOP
54
55 C !USES:
56 IMPLICIT NONE
57 C Global variables / common blocks
58 #include "SIZE.h"
59 #include "EEPARAMS.h"
60 #include "PARAMS.h"
61 #include "EESUPPORT.h"
62 #ifdef ALLOW_EXCH2
63 #include "W2_EXCH2_TOPOLOGY.h"
64 #include "W2_EXCH2_PARAMS.h"
65 #endif /* ALLOW_EXCH2 */
66 #include "MDSIO_SCPU.h"
67
68 C !INPUT PARAMETERS:
69 CHARACTER*(*) fName
70 INTEGER filePrec
71 LOGICAL useCurrentDir
72 CHARACTER*(2) arrType
73 INTEGER kSize, kLo, kHi
74 INTEGER irecord
75 INTEGER myThid
76 C !OUTPUT PARAMETERS:
77 Real arr(*)
78
79 C !FUNCTIONS
80 INTEGER ILNBLNK
81 INTEGER MDS_RECLEN
82 LOGICAL MASTER_CPU_IO
83 EXTERNAL ILNBLNK
84 EXTERNAL MDS_RECLEN
85 EXTERNAL MASTER_CPU_IO
86
87 C !LOCAL VARIABLES:
88 CHARACTER*(MAX_LEN_FNAM) dataFName,pfName
89 CHARACTER*(MAX_LEN_MBUF) msgBuf
90 LOGICAL exst
91 LOGICAL globalFile, fileIsOpen
92 LOGICAL iAmDoingIO
93 INTEGER iG,jG,bi,bj,i,j,k,nNz
94 INTEGER irec,dUnit,IL,pIL
95 INTEGER x_size,y_size,length_of_rec
96 #if defined(ALLOW_EXCH2) && !defined(MISSING_TILE_IO)
97 INTEGER iG_IO,jG_IO,npe, loc_xGlobalLo, loc_yGlobalLo
98 PARAMETER ( x_size = exch2_domain_nxt * sNx )
99 PARAMETER ( y_size = exch2_domain_nyt * sNy )
100 #else
101 PARAMETER ( x_size = Nx )
102 PARAMETER ( y_size = Ny )
103 #endif
104 Real*4 xy_buffer_r4(x_size,y_size)
105 Real*4 r4seg(sNx)
106 Real*8 r8seg(sNx)
107 Real*8 xy_buffer_r8(x_size,y_size)
108 Real*8 globalBuf(Nx,Ny)
109 #ifdef ALLOW_EXCH2
110 INTEGER iGjLoc, jGjLoc
111 c INTEGER tGy,tGx,tNy,tNx,tN
112 INTEGER tGy,tGx, tNx,tN
113 #endif /* ALLOW_EXCH2 */
114 INTEGER tNy
115
116
117 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
118
119 C Assume nothing
120 globalFile = .FALSE.
121 fileIsOpen = .FALSE.
122 IL = ILNBLNK( fName )
123 pIL = ILNBLNK( mdsioLocalDir )
124 nNz = 1 + kHi - kLo
125
126 C Only do I/O if I am the master thread (and mpi process 0 IF useSingleCpuIO):
127 iAmDoingIO = MASTER_CPU_IO(myThid)
128
129 C Only do I/O if I am the master thread
130 IF ( iAmDoingIO ) THEN
131
132 C Record number must be >= 1
133 IF (irecord .LT. 1) THEN
134 WRITE(msgBuf,'(A,I9.8)')
135 & ' MDS_READ_FIELD: argument irecord = ',irecord
136 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
137 & SQUEEZE_RIGHT , myThid)
138 WRITE(msgBuf,'(A)')
139 & ' MDS_READ_FIELD: Invalid value for irecord'
140 CALL PRINT_ERROR( msgBuf, myThid )
141 STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
142 ENDIF
143 C check for valid sub-set of levels:
144 IF ( kLo.LT.1 .OR. kHi.GT.kSize ) THEN
145 WRITE(msgBuf,'(3(A,I6))')
146 & ' MDS_READ_FIELD: arguments kSize=', kSize,
147 & ' , kLo=', kLo, ' , kHi=', kHi
148 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
149 & SQUEEZE_RIGHT , myThid)
150 WRITE(msgBuf,'(A)')
151 & ' MDS_READ_FIELD: invalid sub-set of levels'
152 CALL PRINT_ERROR( msgBuf, myThid )
153 STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
154 ENDIF
155
156 C Assign special directory
157 IF ( useCurrentDir .OR. pIL.EQ.0 ) THEN
158 pfName= fName
159 ELSE
160 WRITE(pfName,'(2a)') mdsioLocalDir(1:pIL), fName(1:IL)
161 ENDIF
162 pIL=ILNBLNK( pfName )
163
164 C Assign a free unit number as the I/O channel for this routine
165 CALL MDSFINDUNIT( dUnit, myThid )
166
167 C Check first for global file with simple name (ie. fName)
168 dataFName = fName
169 INQUIRE( file=dataFName, exist=exst )
170 IF (exst) THEN
171 IF ( debugLevel .GE. debLevA ) THEN
172 WRITE(msgBuf,'(A,A)')
173 & ' MDS_READ_FIELD: opening global file: ',dataFName(1:IL)
174 #ifndef ALLOW_ECCO
175 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
176 & SQUEEZE_RIGHT , myThid)
177 #endif
178 ENDIF
179 globalFile = .TRUE.
180 ENDIF
181
182 C If negative check for global file with MDS name (ie. fName.data)
183 IF (.NOT. globalFile) THEN
184 WRITE(dataFName,'(2a)') fName(1:IL),'.data'
185 INQUIRE( file=dataFName, exist=exst )
186 IF (exst) THEN
187 IF ( debugLevel .GE. debLevA ) THEN
188 WRITE(msgBuf,'(A,A)')
189 & ' MDS_READ_FIELD: opening global file: ',dataFName(1:IL+5)
190 #ifndef ALLOW_ECCO
191 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
192 & SQUEEZE_RIGHT , myThid)
193 #endif
194 ENDIF
195 globalFile = .TRUE.
196 ENDIF
197 ENDIF
198
199 C- endif iAmDoingIO
200 ENDIF
201
202 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
203
204 IF ( useSingleCPUIO ) THEN
205
206 C master thread of process 0, only, opens a global file
207 IF ( iAmDoingIO ) THEN
208 C If global file is visible to process 0, then open it here.
209 C Otherwise stop program.
210 IF ( globalFile) THEN
211 length_of_rec=MDS_RECLEN( filePrec, x_size*y_size, myThid )
212 OPEN( dUnit, file=dataFName, status='old',
213 & access='direct', recl=length_of_rec )
214 ELSE
215 WRITE(msgBuf,'(2A)')
216 & ' MDS_READ_FIELD: filename: ', dataFName(1:IL+5)
217 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
218 & SQUEEZE_RIGHT , myThid)
219 CALL PRINT_ERROR( msgBuf, myThid )
220 WRITE(msgBuf,'(A)')
221 & ' MDS_READ_FIELD: File does not exist'
222 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
223 & SQUEEZE_RIGHT , myThid)
224 CALL PRINT_ERROR( msgBuf, myThid )
225 STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
226 ENDIF
227 C- endif iAmDoingIO
228 ENDIF
229
230 DO k=kLo,kHi
231
232 C master thread of process 0, only, read from file
233 IF ( iAmDoingIO ) THEN
234 irec = k+1-kLo+nNz*(irecord-1)
235 IF (filePrec .EQ. precFloat32) THEN
236 READ(dUnit,rec=irec) xy_buffer_r4
237 #ifdef _BYTESWAPIO
238 CALL MDS_BYTESWAPR4( x_size*y_size, xy_buffer_r4 )
239 #endif
240 ELSEIF (filePrec .EQ. precFloat64) THEN
241 READ(dUnit,rec=irec) xy_buffer_r8
242 #ifdef _BYTESWAPIO
243 CALL MDS_BYTESWAPR8( x_size*y_size, xy_buffer_r8 )
244 #endif
245 ELSE
246 WRITE(msgBuf,'(A)')
247 & ' MDS_READ_FIELD: illegal value for filePrec'
248 CALL PRINT_ERROR( msgBuf, myThid )
249 STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
250 ENDIF
251
252 #if defined(ALLOW_EXCH2) && !defined(MISSING_TILE_IO)
253 bj=1
254 DO npe=1,nPx*nPy
255 DO bi=1,nSx
256 #ifdef ALLOW_USE_MPI
257 loc_xGlobalLo = mpi_myXGlobalLo(npe)
258 loc_yGlobalLo = mpi_myYGlobalLo(npe)
259 #else /* ALLOW_USE_MPI */
260 loc_xGlobalLo = myXGlobalLo
261 loc_yGlobalLo = myYGlobalLo
262 #endif /* ALLOW_USE_MPI */
263 tN = W2_mpi_myTileList(npe,bi)
264 IF ( exch2_mydNx(tN) .GT. x_size ) THEN
265 C- face x-size larger than glob-size : fold it
266 iGjLoc = 0
267 jGjLoc = exch2_mydNx(tN) / x_size
268 ELSEIF ( exch2_tNy(tN) .GT. y_size ) THEN
269 C- tile y-size larger than glob-size : make a long line
270 iGjLoc = exch2_mydNx(tN)
271 jGjLoc = 0
272 ELSE
273 C- default (face fit into global-IO-array)
274 iGjLoc = 0
275 jGjLoc = 1
276 ENDIF
277
278 IF (filePrec .EQ. precFloat32) THEN
279 DO J=1,sNy
280 DO I=1,sNx
281 iG = loc_xGlobalLo-1+(bi-1)*sNx+i
282 jG = loc_yGlobalLo-1+(bj-1)*sNy+j
283 iG_IO=exch2_txGlobalo(tN)+iGjLoc*(j-1)+i-1
284 jG_IO=exch2_tyGlobalo(tN)+jGjLoc*(j-1)
285 globalBuf(iG,jG) = xy_buffer_r4(iG_IO,jG_IO)
286 ENDDO
287 ENDDO
288 ELSEIF (filePrec .EQ. precFloat64) THEN
289 DO J=1,sNy
290 DO I=1,sNx
291 iG = loc_xGlobalLo-1+(bi-1)*sNx+i
292 jG = loc_yGlobalLo-1+(bj-1)*sNy+j
293 iG_IO=exch2_txGlobalo(tN)+iGjLoc*(j-1)+i-1
294 jG_IO=exch2_tyGlobalo(tN)+jGjLoc*(j-1)
295 globalBuf(iG,jG) = xy_buffer_r8(iG_IO,jG_IO)
296 ENDDO
297 ENDDO
298 ENDIF
299
300 C-- end of npe & bi loops
301 ENDDO
302 ENDDO
303 #else /* defined(ALLOW_EXCH2) && !defined(MISSING_TILE_IO) */
304 IF (filePrec .EQ. precFloat32) THEN
305 DO J=1,Ny
306 DO I=1,Nx
307 globalBuf(I,J) = xy_buffer_r4(I,J)
308 ENDDO
309 ENDDO
310 ELSEIF (filePrec .EQ. precFloat64) THEN
311 DO J=1,Ny
312 DO I=1,Nx
313 globalBuf(I,J) = xy_buffer_r8(I,J)
314 ENDDO
315 ENDDO
316 ENDIF
317 #endif /* defined(ALLOW_EXCH2) && !defined(MISSING_TILE_IO) */
318 C- endif iAmDoingIO
319 ENDIF
320 CALL SCATTER_2D(globalBuf,sharedLocalBuf,myThid)
321 IF (arrType .EQ. 'RS') THEN
322 CALL MDS_PASStoRS( sharedLocalBuf,arr,k,kSize,.TRUE.,myThid )
323 ELSEIF (arrType .EQ. 'RL') THEN
324 CALL MDS_PASStoRL( sharedLocalBuf,arr,k,kSize,.TRUE.,myThid )
325 ELSE
326 WRITE(msgBuf,'(A)')
327 & ' MDS_READ_FIELD: illegal value for arrType'
328 CALL PRINT_ERROR( msgBuf, myThid )
329 STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
330 ENDIF
331
332 ENDDO
333 c ENDDO k=kLo,kHi
334
335 IF ( iAmDoingIO ) THEN
336 CLOSE( dUnit )
337 ENDIF
338
339 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
340 C--- else .NOT.useSingleCpuIO
341 ELSE
342
343 C Only do I/O if I am the master thread
344 IF ( iAmDoingIO ) THEN
345
346 C If we are reading from a global file then we open it here
347 IF (globalFile) THEN
348 length_of_rec=MDS_RECLEN( filePrec, sNx, myThid )
349 OPEN( dUnit, file=dataFName, status='old',
350 & access='direct', recl=length_of_rec )
351 fileIsOpen=.TRUE.
352 ENDIF
353
354 C Loop over all tiles
355 DO bj=1,nSy
356 DO bi=1,nSx
357 C If we are reading from a tiled MDS file then we open each one here
358 IF (.NOT. globalFile) THEN
359 iG=bi+(myXGlobalLo-1)/sNx ! Kludge until unstructered tiles
360 jG=bj+(myYGlobalLo-1)/sNy ! Kludge until unstructered tiles
361 WRITE(dataFName,'(2A,I3.3,A,I3.3,A)')
362 & pfName(1:pIL),'.',iG,'.',jG,'.data'
363 INQUIRE( file=dataFName, exist=exst )
364 C Of course, we only open the file if the tile is "active"
365 C (This is a place-holder for the active/passive mechanism
366 IF (exst) THEN
367 IF ( debugLevel .GE. debLevA ) THEN
368 WRITE(msgBuf,'(A,A)')
369 & ' MDS_READ_FIELD: opening file: ',dataFName(1:pIL+13)
370 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
371 & SQUEEZE_RIGHT , myThid)
372 ENDIF
373 length_of_rec=MDS_RECLEN( filePrec, sNx, myThid )
374 OPEN( dUnit, file=dataFName, status='old',
375 & access='direct', recl=length_of_rec )
376 fileIsOpen=.TRUE.
377 ELSE
378 fileIsOpen=.FALSE.
379 WRITE(msgBuf,'(4A)') ' MDS_READ_FIELD: filename: ',
380 & fName(1:IL),' , ', dataFName(1:pIL+13)
381 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
382 & SQUEEZE_RIGHT , myThid)
383 CALL PRINT_ERROR( msgBuf, myThid )
384 WRITE(msgBuf,'(A)')
385 & ' MDS_READ_FIELD: Files DO not exist'
386 CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
387 & SQUEEZE_RIGHT , myThid)
388 CALL PRINT_ERROR( msgBuf, myThid )
389 STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
390 ENDIF
391 ENDIF
392
393 IF (fileIsOpen) THEN
394 tNy = sNy
395 #ifdef ALLOW_EXCH2
396 tN = W2_myTileList(bi)
397 tGy = exch2_tyGlobalo(tN)
398 tGx = exch2_txGlobalo(tN)
399 tNy = exch2_tNy(tN)
400 tNx = exch2_tNx(tN)
401 IF ( exch2_mydNx(tN) .GT. x_size ) THEN
402 C- face x-size larger than glob-size : fold it
403 iGjLoc = 0
404 jGjLoc = exch2_mydNx(tN) / x_size
405 ELSEIF ( exch2_tNy(tN) .GT. y_size ) THEN
406 C- tile y-size larger than glob-size : make a long line
407 iGjLoc = exch2_mydNx(tN)
408 jGjLoc = 0
409 ELSE
410 C- default (face fit into global-IO-array)
411 iGjLoc = 0
412 jGjLoc = 1
413 ENDIF
414 #endif /* ALLOW_EXCH2 */
415 DO k=kLo,kHi
416 DO j=1,tNy
417 IF (globalFile) THEN
418 #ifdef ALLOW_EXCH2
419 irec = 1 + ( tGx-1 + (j-1)*iGjLoc )/tNx
420 & + ( tGy-1 + (j-1)*jGjLoc )*exch2_domain_nxt
421 & + ( k-kLo + (irecord-1)*nNz
422 & )*y_size*exch2_domain_nxt
423 #else /* ALLOW_EXCH2 */
424 iG = myXGlobalLo-1 + (bi-1)*sNx
425 jG = myYGlobalLo-1 + (bj-1)*sNy
426 irec= 1 + INT(iG/sNx) + nSx*nPx*(jG+j-1)
427 & + nSx*nPx*Ny*(k-kLo)
428 & + nSx*nPx*Ny*nNz*(irecord-1)
429 #endif /* ALLOW_EXCH2 */
430 ELSE
431 irec=j + sNy*(k-kLo) + sNy*nNz*(irecord-1)
432 ENDIF
433 IF (filePrec .EQ. precFloat32) THEN
434 READ(dUnit,rec=irec) r4seg
435 #ifdef _BYTESWAPIO
436 CALL MDS_BYTESWAPR4( sNx, r4seg )
437 #endif
438 IF (arrType .EQ. 'RS') THEN
439 CALL MDS_SEG4toRS( j,bi,bj,k,kSize, r4seg, .TRUE., arr )
440 ELSEIF (arrType .EQ. 'RL') THEN
441 CALL MDS_SEG4toRL( j,bi,bj,k,kSize, r4seg, .TRUE., arr )
442 ELSE
443 WRITE(msgBuf,'(A)')
444 & ' MDS_READ_FIELD: illegal value for arrType'
445 CALL PRINT_ERROR( msgBuf, myThid )
446 STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
447 ENDIF
448 ELSEIF (filePrec .EQ. precFloat64) THEN
449 READ(dUnit,rec=irec) r8seg
450 #ifdef _BYTESWAPIO
451 CALL MDS_BYTESWAPR8( sNx, r8seg )
452 #endif
453 IF (arrType .EQ. 'RS') THEN
454 CALL MDS_SEG8toRS( j,bi,bj,k,kSize, r8seg, .TRUE., arr )
455 ELSEIF (arrType .EQ. 'RL') THEN
456 CALL MDS_SEG8toRL( j,bi,bj,k,kSize, r8seg, .TRUE., arr )
457 ELSE
458 WRITE(msgBuf,'(A)')
459 & ' MDS_READ_FIELD: illegal value for arrType'
460 CALL PRINT_ERROR( msgBuf, myThid )
461 STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
462 ENDIF
463 ELSE
464 WRITE(msgBuf,'(A)')
465 & ' MDS_READ_FIELD: illegal value for filePrec'
466 CALL PRINT_ERROR( msgBuf, myThid )
467 STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
468 ENDIF
469 C End of j loop
470 ENDDO
471 C End of k loop
472 ENDDO
473 C end if fileIsOpen
474 ENDIF
475 IF (fileIsOpen .AND. (.NOT. globalFile)) THEN
476 CLOSE( dUnit )
477 fileIsOpen = .FALSE.
478 ENDIF
479 C End of bi,bj loops
480 ENDDO
481 ENDDO
482
483 C If global file was opened then close it
484 IF (fileIsOpen .AND. globalFile) THEN
485 CLOSE( dUnit )
486 fileIsOpen = .FALSE.
487 ENDIF
488
489 C- endif iAmDoingIO
490 ENDIF
491
492 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
493 C if useSingleCpuIO / else / end
494 ENDIF
495
496 RETURN
497 END

  ViewVC Help
Powered by ViewVC 1.1.22