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

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

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


Revision 1.4 - (hide annotations) (download)
Tue Dec 30 00:13:35 2008 UTC (15 years, 5 months ago) by jahn
Branch: MAIN
CVS Tags: checkpoint61l, checkpoint61m, checkpoint61j, checkpoint61k, checkpoint61h, checkpoint61i
Changes since 1.3: +2 -10 lines
move buffers to common block to save some memory

1 jahn 1.4 C $Header: /u/gcmpack/MITgcm/pkg/mdsio/mdsio_read_field.F,v 1.3 2007/11/13 19:37:44 jmc Exp $
2 jmc 1.1 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 jmc 1.3 I kSize,kLo,kHi,
15 jmc 1.1 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 jmc 1.3 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 jmc 1.1 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 jmc 1.3 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 jmc 1.1 C The precision of the file is decsribed by filePrec, set either
40 jmc 1.3 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 jmc 1.1 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 jmc 1.3 INTEGER kSize, kLo, kHi
74 jmc 1.1 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 jmc 1.3 INTEGER iG,jG,bi,bj,i,j,k,nNz
94     INTEGER irec,dUnit,IL,pIL
95 jahn 1.4 INTEGER length_of_rec
96 jmc 1.1 #if defined(ALLOW_EXCH2) && !defined(MISSING_TILE_IO)
97 jmc 1.2 INTEGER iG_IO,jG_IO,npe, loc_xGlobalLo, loc_yGlobalLo
98 jmc 1.1 #endif
99     Real*4 r4seg(sNx)
100     Real*8 r8seg(sNx)
101     #ifdef ALLOW_EXCH2
102 jmc 1.2 INTEGER iGjLoc, jGjLoc
103     c INTEGER tGy,tGx,tNy,tNx,tN
104     INTEGER tGy,tGx, tNx,tN
105 jmc 1.1 #endif /* ALLOW_EXCH2 */
106     INTEGER tNy
107    
108    
109     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
110    
111     C Assume nothing
112     globalFile = .FALSE.
113     fileIsOpen = .FALSE.
114     IL = ILNBLNK( fName )
115     pIL = ILNBLNK( mdsioLocalDir )
116 jmc 1.3 nNz = 1 + kHi - kLo
117 jmc 1.1
118     C Only do I/O if I am the master thread (and mpi process 0 IF useSingleCpuIO):
119     iAmDoingIO = MASTER_CPU_IO(myThid)
120    
121     C Only do I/O if I am the master thread
122     IF ( iAmDoingIO ) THEN
123    
124     C Record number must be >= 1
125     IF (irecord .LT. 1) THEN
126     WRITE(msgBuf,'(A,I9.8)')
127     & ' MDS_READ_FIELD: argument irecord = ',irecord
128     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
129     & SQUEEZE_RIGHT , myThid)
130     WRITE(msgBuf,'(A)')
131     & ' MDS_READ_FIELD: Invalid value for irecord'
132     CALL PRINT_ERROR( msgBuf, myThid )
133     STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
134     ENDIF
135 jmc 1.3 C check for valid sub-set of levels:
136     IF ( kLo.LT.1 .OR. kHi.GT.kSize ) THEN
137     WRITE(msgBuf,'(3(A,I6))')
138     & ' MDS_READ_FIELD: arguments kSize=', kSize,
139     & ' , kLo=', kLo, ' , kHi=', kHi
140     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
141     & SQUEEZE_RIGHT , myThid)
142     WRITE(msgBuf,'(A)')
143     & ' MDS_READ_FIELD: invalid sub-set of levels'
144     CALL PRINT_ERROR( msgBuf, myThid )
145     STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
146     ENDIF
147 jmc 1.1
148     C Assign special directory
149     IF ( useCurrentDir .OR. pIL.EQ.0 ) THEN
150     pfName= fName
151     ELSE
152     WRITE(pfName,'(2a)') mdsioLocalDir(1:pIL), fName(1:IL)
153     ENDIF
154     pIL=ILNBLNK( pfName )
155    
156     C Assign a free unit number as the I/O channel for this routine
157     CALL MDSFINDUNIT( dUnit, myThid )
158    
159     C Check first for global file with simple name (ie. fName)
160     dataFName = fName
161     INQUIRE( file=dataFName, exist=exst )
162     IF (exst) THEN
163     IF ( debugLevel .GE. debLevA ) THEN
164     WRITE(msgBuf,'(A,A)')
165     & ' MDS_READ_FIELD: opening global file: ',dataFName(1:IL)
166     #ifndef ALLOW_ECCO
167     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
168     & SQUEEZE_RIGHT , myThid)
169     #endif
170     ENDIF
171     globalFile = .TRUE.
172     ENDIF
173    
174     C If negative check for global file with MDS name (ie. fName.data)
175     IF (.NOT. globalFile) THEN
176     WRITE(dataFName,'(2a)') fName(1:IL),'.data'
177     INQUIRE( file=dataFName, exist=exst )
178     IF (exst) THEN
179     IF ( debugLevel .GE. debLevA ) THEN
180     WRITE(msgBuf,'(A,A)')
181     & ' MDS_READ_FIELD: opening global file: ',dataFName(1:IL+5)
182     #ifndef ALLOW_ECCO
183     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
184     & SQUEEZE_RIGHT , myThid)
185     #endif
186     ENDIF
187     globalFile = .TRUE.
188     ENDIF
189     ENDIF
190    
191     C- endif iAmDoingIO
192     ENDIF
193    
194     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
195    
196     IF ( useSingleCPUIO ) THEN
197    
198     C master thread of process 0, only, opens a global file
199     IF ( iAmDoingIO ) THEN
200     C If global file is visible to process 0, then open it here.
201     C Otherwise stop program.
202     IF ( globalFile) THEN
203     length_of_rec=MDS_RECLEN( filePrec, x_size*y_size, myThid )
204     OPEN( dUnit, file=dataFName, status='old',
205     & access='direct', recl=length_of_rec )
206     ELSE
207     WRITE(msgBuf,'(2A)')
208     & ' MDS_READ_FIELD: filename: ', dataFName(1:IL+5)
209     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
210     & SQUEEZE_RIGHT , myThid)
211     CALL PRINT_ERROR( msgBuf, myThid )
212     WRITE(msgBuf,'(A)')
213     & ' MDS_READ_FIELD: File does not exist'
214     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
215     & SQUEEZE_RIGHT , myThid)
216     CALL PRINT_ERROR( msgBuf, myThid )
217     STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
218     ENDIF
219     C- endif iAmDoingIO
220     ENDIF
221    
222 jmc 1.3 DO k=kLo,kHi
223 jmc 1.1
224     C master thread of process 0, only, read from file
225     IF ( iAmDoingIO ) THEN
226 jmc 1.3 irec = k+1-kLo+nNz*(irecord-1)
227 jmc 1.1 IF (filePrec .EQ. precFloat32) THEN
228     READ(dUnit,rec=irec) xy_buffer_r4
229     #ifdef _BYTESWAPIO
230     CALL MDS_BYTESWAPR4( x_size*y_size, xy_buffer_r4 )
231     #endif
232 jmc 1.2 ELSEIF (filePrec .EQ. precFloat64) THEN
233     READ(dUnit,rec=irec) xy_buffer_r8
234     #ifdef _BYTESWAPIO
235     CALL MDS_BYTESWAPR8( x_size*y_size, xy_buffer_r8 )
236     #endif
237     ELSE
238     WRITE(msgBuf,'(A)')
239     & ' MDS_READ_FIELD: illegal value for filePrec'
240     CALL PRINT_ERROR( msgBuf, myThid )
241     STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
242     ENDIF
243    
244 jmc 1.1 #if defined(ALLOW_EXCH2) && !defined(MISSING_TILE_IO)
245 jmc 1.2 bj=1
246     DO npe=1,nPx*nPy
247     DO bi=1,nSx
248     #ifdef ALLOW_USE_MPI
249     loc_xGlobalLo = mpi_myXGlobalLo(npe)
250     loc_yGlobalLo = mpi_myYGlobalLo(npe)
251     #else /* ALLOW_USE_MPI */
252     loc_xGlobalLo = myXGlobalLo
253     loc_yGlobalLo = myYGlobalLo
254     #endif /* ALLOW_USE_MPI */
255     tN = W2_mpi_myTileList(npe,bi)
256     IF ( exch2_mydNx(tN) .GT. x_size ) THEN
257     C- face x-size larger than glob-size : fold it
258     iGjLoc = 0
259     jGjLoc = exch2_mydNx(tN) / x_size
260     ELSEIF ( exch2_tNy(tN) .GT. y_size ) THEN
261     C- tile y-size larger than glob-size : make a long line
262     iGjLoc = exch2_mydNx(tN)
263     jGjLoc = 0
264     ELSE
265     C- default (face fit into global-IO-array)
266     iGjLoc = 0
267     jGjLoc = 1
268     ENDIF
269    
270     IF (filePrec .EQ. precFloat32) THEN
271 jmc 1.1 DO J=1,sNy
272     DO I=1,sNx
273 jmc 1.2 iG = loc_xGlobalLo-1+(bi-1)*sNx+i
274     jG = loc_yGlobalLo-1+(bj-1)*sNy+j
275     iG_IO=exch2_txGlobalo(tN)+iGjLoc*(j-1)+i-1
276     jG_IO=exch2_tyGlobalo(tN)+jGjLoc*(j-1)
277 jmc 1.1 globalBuf(iG,jG) = xy_buffer_r4(iG_IO,jG_IO)
278     ENDDO
279     ENDDO
280 jmc 1.2 ELSEIF (filePrec .EQ. precFloat64) THEN
281     DO J=1,sNy
282     DO I=1,sNx
283     iG = loc_xGlobalLo-1+(bi-1)*sNx+i
284     jG = loc_yGlobalLo-1+(bj-1)*sNy+j
285     iG_IO=exch2_txGlobalo(tN)+iGjLoc*(j-1)+i-1
286     jG_IO=exch2_tyGlobalo(tN)+jGjLoc*(j-1)
287     globalBuf(iG,jG) = xy_buffer_r8(iG_IO,jG_IO)
288     ENDDO
289     ENDDO
290     ENDIF
291    
292     C-- end of npe & bi loops
293 jmc 1.1 ENDDO
294 jmc 1.2 ENDDO
295 jmc 1.1 #else /* defined(ALLOW_EXCH2) && !defined(MISSING_TILE_IO) */
296 jmc 1.2 IF (filePrec .EQ. precFloat32) THEN
297 jmc 1.1 DO J=1,Ny
298     DO I=1,Nx
299     globalBuf(I,J) = xy_buffer_r4(I,J)
300     ENDDO
301     ENDDO
302     ELSEIF (filePrec .EQ. precFloat64) THEN
303     DO J=1,Ny
304     DO I=1,Nx
305     globalBuf(I,J) = xy_buffer_r8(I,J)
306     ENDDO
307     ENDDO
308 jmc 1.2 ENDIF
309 jmc 1.1 #endif /* defined(ALLOW_EXCH2) && !defined(MISSING_TILE_IO) */
310     C- endif iAmDoingIO
311     ENDIF
312     CALL SCATTER_2D(globalBuf,sharedLocalBuf,myThid)
313     IF (arrType .EQ. 'RS') THEN
314 jmc 1.3 CALL MDS_PASStoRS( sharedLocalBuf,arr,k,kSize,.TRUE.,myThid )
315 jmc 1.1 ELSEIF (arrType .EQ. 'RL') THEN
316 jmc 1.3 CALL MDS_PASStoRL( sharedLocalBuf,arr,k,kSize,.TRUE.,myThid )
317 jmc 1.1 ELSE
318     WRITE(msgBuf,'(A)')
319     & ' MDS_READ_FIELD: illegal value for arrType'
320     CALL PRINT_ERROR( msgBuf, myThid )
321     STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
322     ENDIF
323    
324     ENDDO
325 jmc 1.3 c ENDDO k=kLo,kHi
326 jmc 1.1
327     IF ( iAmDoingIO ) THEN
328     CLOSE( dUnit )
329     ENDIF
330    
331     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
332     C--- else .NOT.useSingleCpuIO
333     ELSE
334    
335     C Only do I/O if I am the master thread
336     IF ( iAmDoingIO ) THEN
337    
338     C If we are reading from a global file then we open it here
339     IF (globalFile) THEN
340     length_of_rec=MDS_RECLEN( filePrec, sNx, myThid )
341     OPEN( dUnit, file=dataFName, status='old',
342     & access='direct', recl=length_of_rec )
343     fileIsOpen=.TRUE.
344     ENDIF
345    
346     C Loop over all tiles
347     DO bj=1,nSy
348     DO bi=1,nSx
349     C If we are reading from a tiled MDS file then we open each one here
350     IF (.NOT. globalFile) THEN
351     iG=bi+(myXGlobalLo-1)/sNx ! Kludge until unstructered tiles
352     jG=bj+(myYGlobalLo-1)/sNy ! Kludge until unstructered tiles
353     WRITE(dataFName,'(2A,I3.3,A,I3.3,A)')
354     & pfName(1:pIL),'.',iG,'.',jG,'.data'
355     INQUIRE( file=dataFName, exist=exst )
356     C Of course, we only open the file if the tile is "active"
357     C (This is a place-holder for the active/passive mechanism
358     IF (exst) THEN
359     IF ( debugLevel .GE. debLevA ) THEN
360     WRITE(msgBuf,'(A,A)')
361     & ' MDS_READ_FIELD: opening file: ',dataFName(1:pIL+13)
362     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
363     & SQUEEZE_RIGHT , myThid)
364     ENDIF
365     length_of_rec=MDS_RECLEN( filePrec, sNx, myThid )
366     OPEN( dUnit, file=dataFName, status='old',
367     & access='direct', recl=length_of_rec )
368     fileIsOpen=.TRUE.
369     ELSE
370     fileIsOpen=.FALSE.
371     WRITE(msgBuf,'(4A)') ' MDS_READ_FIELD: filename: ',
372     & fName(1:IL),' , ', dataFName(1:pIL+13)
373     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
374     & SQUEEZE_RIGHT , myThid)
375     CALL PRINT_ERROR( msgBuf, myThid )
376     WRITE(msgBuf,'(A)')
377     & ' MDS_READ_FIELD: Files DO not exist'
378     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
379     & SQUEEZE_RIGHT , myThid)
380     CALL PRINT_ERROR( msgBuf, myThid )
381     STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
382     ENDIF
383     ENDIF
384    
385     IF (fileIsOpen) THEN
386     tNy = sNy
387     #ifdef ALLOW_EXCH2
388 jmc 1.2 tN = W2_myTileList(bi)
389     tGy = exch2_tyGlobalo(tN)
390     tGx = exch2_txGlobalo(tN)
391     tNy = exch2_tNy(tN)
392     tNx = exch2_tNx(tN)
393     IF ( exch2_mydNx(tN) .GT. x_size ) THEN
394     C- face x-size larger than glob-size : fold it
395     iGjLoc = 0
396     jGjLoc = exch2_mydNx(tN) / x_size
397     ELSEIF ( exch2_tNy(tN) .GT. y_size ) THEN
398     C- tile y-size larger than glob-size : make a long line
399     iGjLoc = exch2_mydNx(tN)
400     jGjLoc = 0
401     ELSE
402     C- default (face fit into global-IO-array)
403     iGjLoc = 0
404     jGjLoc = 1
405     ENDIF
406 jmc 1.1 #endif /* ALLOW_EXCH2 */
407 jmc 1.3 DO k=kLo,kHi
408 jmc 1.1 DO j=1,tNy
409     IF (globalFile) THEN
410     #ifdef ALLOW_EXCH2
411 jmc 1.2 irec = 1 + ( tGx-1 + (j-1)*iGjLoc )/tNx
412     & + ( tGy-1 + (j-1)*jGjLoc )*exch2_domain_nxt
413 jmc 1.3 & + ( k-kLo + (irecord-1)*nNz
414 jmc 1.2 & )*y_size*exch2_domain_nxt
415 jmc 1.1 #else /* ALLOW_EXCH2 */
416     iG = myXGlobalLo-1 + (bi-1)*sNx
417     jG = myYGlobalLo-1 + (bj-1)*sNy
418     irec= 1 + INT(iG/sNx) + nSx*nPx*(jG+j-1)
419 jmc 1.3 & + nSx*nPx*Ny*(k-kLo)
420 jmc 1.1 & + nSx*nPx*Ny*nNz*(irecord-1)
421     #endif /* ALLOW_EXCH2 */
422     ELSE
423 jmc 1.3 irec=j + sNy*(k-kLo) + sNy*nNz*(irecord-1)
424 jmc 1.1 ENDIF
425     IF (filePrec .EQ. precFloat32) THEN
426     READ(dUnit,rec=irec) r4seg
427     #ifdef _BYTESWAPIO
428     CALL MDS_BYTESWAPR4( sNx, r4seg )
429     #endif
430     IF (arrType .EQ. 'RS') THEN
431 jmc 1.3 CALL MDS_SEG4toRS( j,bi,bj,k,kSize, r4seg, .TRUE., arr )
432 jmc 1.1 ELSEIF (arrType .EQ. 'RL') THEN
433 jmc 1.3 CALL MDS_SEG4toRL( j,bi,bj,k,kSize, r4seg, .TRUE., arr )
434 jmc 1.1 ELSE
435     WRITE(msgBuf,'(A)')
436     & ' MDS_READ_FIELD: illegal value for arrType'
437     CALL PRINT_ERROR( msgBuf, myThid )
438     STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
439     ENDIF
440     ELSEIF (filePrec .EQ. precFloat64) THEN
441     READ(dUnit,rec=irec) r8seg
442     #ifdef _BYTESWAPIO
443     CALL MDS_BYTESWAPR8( sNx, r8seg )
444     #endif
445     IF (arrType .EQ. 'RS') THEN
446 jmc 1.3 CALL MDS_SEG8toRS( j,bi,bj,k,kSize, r8seg, .TRUE., arr )
447 jmc 1.1 ELSEIF (arrType .EQ. 'RL') THEN
448 jmc 1.3 CALL MDS_SEG8toRL( j,bi,bj,k,kSize, r8seg, .TRUE., arr )
449 jmc 1.1 ELSE
450     WRITE(msgBuf,'(A)')
451     & ' MDS_READ_FIELD: illegal value for arrType'
452     CALL PRINT_ERROR( msgBuf, myThid )
453     STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
454     ENDIF
455     ELSE
456     WRITE(msgBuf,'(A)')
457     & ' MDS_READ_FIELD: illegal value for filePrec'
458     CALL PRINT_ERROR( msgBuf, myThid )
459     STOP 'ABNORMAL END: S/R MDS_READ_FIELD'
460     ENDIF
461     C End of j loop
462     ENDDO
463     C End of k loop
464     ENDDO
465     C end if fileIsOpen
466     ENDIF
467     IF (fileIsOpen .AND. (.NOT. globalFile)) THEN
468     CLOSE( dUnit )
469     fileIsOpen = .FALSE.
470     ENDIF
471     C End of bi,bj loops
472     ENDDO
473     ENDDO
474    
475     C If global file was opened then close it
476     IF (fileIsOpen .AND. globalFile) THEN
477     CLOSE( dUnit )
478     fileIsOpen = .FALSE.
479     ENDIF
480    
481     C- endif iAmDoingIO
482     ENDIF
483    
484     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
485     C if useSingleCpuIO / else / end
486     ENDIF
487    
488     RETURN
489     END

  ViewVC Help
Powered by ViewVC 1.1.22