/[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.1 - (show annotations) (download)
Fri Dec 29 05:25:15 2006 UTC (17 years, 6 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint58u_post, checkpoint58w_post, checkpoint58t_post, checkpoint58v_post
new S/R that replaces MDSREADFIELD (& _LOC), with 1 more argument,
 and cleaned-up (remove 1/3 of calls);

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

  ViewVC Help
Powered by ViewVC 1.1.22