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

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

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


Revision 1.8 - (show annotations) (download)
Mon Jan 26 21:08:07 2004 UTC (20 years, 5 months ago) by afe
Branch: MAIN
Changes since 1.7: +5 -5 lines
o modified pkg/mdsio/mdsio_readfield.F -- commented out diskspace-eating
  debug output

1 C $Header: /u/u3/gcmpack/MITgcm/pkg/mdsio/mdsio_readfield.F,v 1.7 2004/01/23 21:26:55 afe Exp $
2 C $Name: $
3
4 #include "MDSIO_OPTIONS.h"
5
6 SUBROUTINE MDSREADFIELD(
7 I fName,
8 I filePrec,
9 I arrType,
10 I nNz,
11 O arr,
12 I irecord,
13 I myThid )
14 C
15 C Arguments:
16 C
17 C fName string base name for file to read
18 C filePrec integer number of bits per word in file (32 or 64)
19 C arrType char(2) declaration of "arr": either "RS" or "RL"
20 C nNz integer size of third dimension: normally either 1 or Nr
21 C arr RS/RL array to read into, arr(:,:,nNz,:,:)
22 C irecord integer record number to read
23 C myThid integer thread identifier
24 C
25 C MDSREADFIELD first checks to see if the file "fName" exists, then
26 C if the file "fName.data" exists and finally the tiled files of the
27 C form "fName.xxx.yyy.data" exist. Currently, the meta-files are not
28 C read because it is difficult to parse files in fortran.
29 C The precision of the file is decsribed by filePrec, set either
30 C to floatPrec32 or floatPrec64. The precision or declaration of
31 C the array argument must be consistently described by the char*(2)
32 C string arrType, either "RS" or "RL". nNz allows for both 2-D and
33 C 3-D arrays to be handled. nNz=1 implies a 2-D model field and
34 C nNz=Nr implies a 3-D model field. irecord is the record number
35 C to be read and must be >= 1. The file data is stored in
36 C arr *but* the overlaps are *not* updated. ie. An exchange must
37 C be called. This is because the routine is sometimes called from
38 C within a MASTER_THID region.
39 C
40 C Created: 03/16/99 adcroft@mit.edu
41
42 implicit none
43 C Global variables / common blocks
44 #include "SIZE.h"
45 #include "EEPARAMS.h"
46 #include "PARAMS.h"
47 #include "EESUPPORT.h"
48 #ifdef USE_W2
49 #include "W2_EXCH2_TOPOLOGY.h"
50 #include "W2_EXCH2_PARAMS.h"
51 #endif /* USE_W2 */
52
53 C Routine arguments
54 character*(*) fName
55 integer filePrec
56 character*(2) arrType
57 integer nNz
58 Real arr(*)
59 integer irecord
60 integer myThid
61 C Functions
62 integer ILNBLNK
63 integer MDS_RECLEN
64 C Local variables
65 character*(80) dataFName,pfName
66 integer iG,jG,irec,bi,bj,j,k,dUnit,IL,pIL
67 logical exst
68 Real*4 r4seg(sNx)
69 Real*8 r8seg(sNx)
70 logical globalFile,fileIsOpen
71 integer length_of_rec
72 character*(max_len_mbuf) msgbuf
73 integer i
74 Real*4 global_r4(Nx,Ny)
75 Real*8 global (Nx,Ny)
76 _RL local(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
77 #ifdef USE_W2
78 integer domainHeight,domainLength,tby,tgx,tny,tnx,tn
79 #endif /* USE_W2 */
80 C ------------------------------------------------------------------
81
82 C Only do I/O if I am the master thread
83 _BEGIN_MASTER( myThid )
84
85 C Record number must be >= 1
86 if (irecord .LT. 1) then
87 write(msgbuf,'(a,i9.8)')
88 & ' MDSREADFIELD: argument irecord = ',irecord
89 call print_message( msgbuf, standardmessageunit,
90 & SQUEEZE_RIGHT , mythid)
91 write(msgbuf,'(a)')
92 & ' MDSREADFIELD: Invalid value for irecord'
93 call print_error( msgbuf, mythid )
94 stop 'ABNORMAL END: S/R MDSREADFIELD'
95 endif
96
97 C Assume nothing
98 globalFile = .FALSE.
99 fileIsOpen = .FALSE.
100 IL = ILNBLNK( fName )
101 pIL = ILNBLNK( mdsioLocalDir )
102
103 C Assign special directory
104 if ( mdsioLocalDir .NE. ' ' ) then
105 write(pFname(1:80),'(2a)') mdsioLocalDir(1:pIL), fName(1:IL)
106 else
107 pFname= fName
108 endif
109 pIL=ILNBLNK( pfName )
110
111 C Assign a free unit number as the I/O channel for this routine
112 call MDSFINDUNIT( dUnit, mythid )
113
114 C Check first for global file with simple name (ie. fName)
115 dataFName = fName
116 inquire( file=dataFname, exist=exst )
117 if (exst) then
118 if ( debugLevel .GE. debLevA ) then
119 write(msgbuf,'(a,a)')
120 & ' MDSREADFIELD: opening global file: ',dataFName
121 call print_message( msgbuf, standardmessageunit,
122 & SQUEEZE_RIGHT , mythid)
123 endif
124 globalFile = .TRUE.
125 endif
126
127 C If negative check for global file with MDS name (ie. fName.data)
128 if (.NOT. globalFile) then
129 write(dataFname(1:80),'(2a)') fName(1:IL),'.data'
130 inquire( file=dataFname, exist=exst )
131 if (exst) then
132 if ( debugLevel .GE. debLevA ) then
133 write(msgbuf,'(a,a)')
134 & ' MDSREADFIELD: opening global file: ',dataFName
135 call print_message( msgbuf, standardmessageunit,
136 & SQUEEZE_RIGHT , mythid)
137 endif
138 globalFile = .TRUE.
139 endif
140 endif
141
142 if ( .not. ( globalFile .and. useSingleCPUIO ) ) then
143
144 C If we are reading from a global file then we open it here
145 if (globalFile) then
146 length_of_rec=MDS_RECLEN( filePrec, sNx, mythid )
147 open( dUnit, file=dataFName, status='old',
148 & access='direct', recl=length_of_rec )
149 fileIsOpen=.TRUE.
150 endif
151
152 #ifdef USE_W2
153
154
155 domainLength = exch2_domain_nxt
156 domainHeight = exch2_domain_nyt
157
158 write(*,fmt='(1X,A,I3,A,I3)') 'L=', domainlength,
159 & 'H=', domainheight
160 C Loop over all tiles
161 do bj=1,nSy
162 do bi=1,nSx
163 C If we are reading from a tiled MDS file then we open each one here
164 if (.NOT. globalFile) then
165 write(msgbuf,'(a)')
166 & ' MDSREADFIELD: non-global input files not
167 & implemented with exch2'
168 call print_error( msgbuf, mythid )
169 stop 'ABNORMAL END: S/R MDSREADFIELD'
170 endif
171 tn = W2_myTileList(bi)
172 tby = exch2_tbasey(tn)
173 tgx = exch2_txglobalo(tn)
174 tny = exch2_tny(tn)
175 tnx = exch2_tnx(tn)
176 if (fileIsOpen) then
177 do k=1,nNz
178 do j=1,tNy
179 c write(*,fmt='(1X,A,I3,A,I3,A,I3,A,I3,A,I3,A,I3)') 'tby=', tby,
180 c & ', tgx=', tgx,
181 c & ', tnx=',tnx, ', tny=', tny, ', j=',j,', tn=',tn
182
183 irec = domainLength*tby + (tgx-1)/tnx + 1 +
184 & domainLength*(j-1) +
185 & domainLength*domainHeight*tny*(k-1) +
186 & domainLength*domainHeight*tny*nNz*(irecord-1)
187 c write(*,fmt='(1X,A,I6,A,I3)') 'record = ',irec,',thingy=',
188 c & (tgx-1)/tnx
189 c write(*,fmt='(1X,A,I6)') 'record = ',irec
190
191 if (filePrec .eq. precFloat32) then
192 read(dUnit,rec=irec) r4seg
193 #ifdef _BYTESWAPIO
194 call MDS_BYTESWAPR4( sNx, r4seg )
195 #endif
196 if (arrType .eq. 'RS') then
197 call MDS_SEG4toRS( j,bi,bj,k,nNz, r4seg, .TRUE., arr )
198 elseif (arrType .eq. 'RL') then
199 call MDS_SEG4toRL( j,bi,bj,k,nNz, r4seg, .TRUE., arr )
200 else
201 write(msgbuf,'(a)')
202 & ' MDSREADFIELD: illegal value for arrType'
203 call print_error( msgbuf, mythid )
204 stop 'ABNORMAL END: S/R MDSREADFIELD'
205 endif
206 elseif (filePrec .eq. precFloat64) then
207 read(dUnit,rec=irec) r8seg
208 #ifdef _BYTESWAPIO
209 call MDS_BYTESWAPR8( sNx, r8seg )
210 #endif
211 if (arrType .eq. 'RS') then
212 call MDS_SEG8toRS( j,bi,bj,k,nNz, r8seg, .TRUE., arr )
213 elseif (arrType .eq. 'RL') then
214 call MDS_SEG8toRL( j,bi,bj,k,nNz, r8seg, .TRUE., arr )
215 else
216 write(msgbuf,'(a)')
217 & ' MDSREADFIELD: illegal value for arrType'
218 call print_error( msgbuf, mythid )
219 stop 'ABNORMAL END: S/R MDSREADFIELD'
220 endif
221 else
222 write(msgbuf,'(a)')
223 & ' MDSREADFIELD: illegal value for filePrec'
224 call print_error( msgbuf, mythid )
225 stop 'ABNORMAL END: S/R MDSREADFIELD'
226 endif
227 C End of j loop
228 enddo
229 C End of k loop
230 enddo
231 if (.NOT. globalFile) then
232 close( dUnit )
233 fileIsOpen = .FALSE.
234 endif
235 endif
236 C End of bi,bj loops
237 enddo
238 enddo
239
240 #else /* .not. USE_W2 */
241
242
243
244 C Loop over all tiles
245 do bj=1,nSy
246 do bi=1,nSx
247 C If we are reading from a tiled MDS file then we open each one here
248 if (.NOT. globalFile) then
249 iG=bi+(myXGlobalLo-1)/sNx ! Kludge until unstructered tiles
250 jG=bj+(myYGlobalLo-1)/sNy ! Kludge until unstructered tiles
251 write(dataFname(1:80),'(2a,i3.3,a,i3.3,a)')
252 & pfName(1:pIL),'.',iG,'.',jG,'.data'
253 inquire( file=dataFname, exist=exst )
254 C Of course, we only open the file if the tile is "active"
255 C (This is a place-holder for the active/passive mechanism
256 if (exst) then
257 if ( debugLevel .GE. debLevA ) then
258 write(msgbuf,'(a,a)')
259 & ' MDSREADFIELD: opening file: ',dataFName
260 call print_message( msgbuf, standardmessageunit,
261 & SQUEEZE_RIGHT , mythid)
262 endif
263 length_of_rec=MDS_RECLEN( filePrec, sNx, mythid )
264 open( dUnit, file=dataFName, status='old',
265 & access='direct', recl=length_of_rec )
266 fileIsOpen=.TRUE.
267 else
268 fileIsOpen=.FALSE.
269 write(msgbuf,'(3a)')
270 & ' MDSREADFIELD: filename: ',dataFName, pfName
271 call print_message( msgbuf, standardmessageunit,
272 & SQUEEZE_RIGHT , mythid)
273 write(msgbuf,'(a)')
274 & ' MDSREADFIELD: File does not exist'
275 call print_error( msgbuf, mythid )
276 stop 'ABNORMAL END: S/R MDSREADFIELD'
277 endif
278 endif
279
280 if (fileIsOpen) then
281 do k=1,nNz
282 do j=1,sNy
283 if (globalFile) then
284 iG = myXGlobalLo-1 + (bi-1)*sNx
285 jG = myYGlobalLo-1 + (bj-1)*sNy
286 irec=1 + INT(iG/sNx) + nSx*nPx*(jG+j-1) + nSx*nPx*Ny*(k-1)
287 & + nSx*nPx*Ny*nNz*(irecord-1)
288 else
289 iG = 0
290 jG = 0
291 irec=j + sNy*(k-1) + sNy*nNz*(irecord-1)
292 endif
293 if (filePrec .eq. precFloat32) then
294 read(dUnit,rec=irec) r4seg
295 #ifdef _BYTESWAPIO
296 call MDS_BYTESWAPR4( sNx, r4seg )
297 #endif
298 if (arrType .eq. 'RS') then
299 call MDS_SEG4toRS( j,bi,bj,k,nNz, r4seg, .TRUE., arr )
300 elseif (arrType .eq. 'RL') then
301 call MDS_SEG4toRL( j,bi,bj,k,nNz, r4seg, .TRUE., arr )
302 else
303 write(msgbuf,'(a)')
304 & ' MDSREADFIELD: illegal value for arrType'
305 call print_error( msgbuf, mythid )
306 stop 'ABNORMAL END: S/R MDSREADFIELD'
307 endif
308 elseif (filePrec .eq. precFloat64) then
309 read(dUnit,rec=irec) r8seg
310 #ifdef _BYTESWAPIO
311 call MDS_BYTESWAPR8( sNx, r8seg )
312 #endif
313 if (arrType .eq. 'RS') then
314 call MDS_SEG8toRS( j,bi,bj,k,nNz, r8seg, .TRUE., arr )
315 elseif (arrType .eq. 'RL') then
316 call MDS_SEG8toRL( j,bi,bj,k,nNz, r8seg, .TRUE., arr )
317 else
318 write(msgbuf,'(a)')
319 & ' MDSREADFIELD: illegal value for arrType'
320 call print_error( msgbuf, mythid )
321 stop 'ABNORMAL END: S/R MDSREADFIELD'
322 endif
323 else
324 write(msgbuf,'(a)')
325 & ' MDSREADFIELD: illegal value for filePrec'
326 call print_error( msgbuf, mythid )
327 stop 'ABNORMAL END: S/R MDSREADFIELD'
328 endif
329 C End of j loop
330 enddo
331 C End of k loop
332 enddo
333 if (.NOT. globalFile) then
334 close( dUnit )
335 fileIsOpen = .FALSE.
336 endif
337 endif
338 C End of bi,bj loops
339 enddo
340 enddo
341
342 #endif /* USE_W2 */
343
344 C If global file was opened then close it
345 if (fileIsOpen .AND. globalFile) then
346 close( dUnit )
347 fileIsOpen = .FALSE.
348 endif
349
350 endif
351 c endif ( .not. ( globalFile .and. useSingleCPUIO ) )
352
353 _END_MASTER( myThid )
354
355 if ( globalFile .and. useSingleCPUIO ) then
356
357 #ifdef USE_W2
358 write(msgbuf,'(a)')
359 & ' MDSREADFIELD: useSingleCPUIO not
360 & implemented with exch2'
361 call print_error( msgbuf, mythid )
362 stop 'ABNORMAL END: S/R MDSREADFIELD'
363 #endif /* USE_W2 */
364
365 C master thread of process 0, only, opens a global file
366 _BEGIN_MASTER( myThid )
367 #ifdef ALLOW_USE_MPI
368 IF( mpiMyId .EQ. 0 ) THEN
369 #else
370 IF ( .TRUE. ) THEN
371 #endif /* ALLOW_USE_MPI */
372 length_of_rec=MDS_RECLEN( filePrec, Nx*Ny, mythid )
373 open( dUnit, file=dataFName, status='old',
374 & access='direct', recl=length_of_rec )
375 ENDIF
376 _END_MASTER( myThid )
377
378 DO k=1,nNz
379
380 _BEGIN_MASTER( myThid )
381 #ifdef ALLOW_USE_MPI
382 IF( mpiMyId .EQ. 0 ) THEN
383 #else
384 IF ( .TRUE. ) THEN
385 #endif /* ALLOW_USE_MPI */
386 irec = k+nNz*(irecord-1)
387 if (filePrec .eq. precFloat32) then
388 read(dUnit,rec=irec) global_r4
389 #ifdef _BYTESWAPIO
390 call MDS_BYTESWAPR4( Nx*Ny, global_r4 )
391 #endif
392 DO J=1,Ny
393 DO I=1,Nx
394 global(I,J) = global_r4(I,J)
395 ENDDO
396 ENDDO
397 elseif (filePrec .eq. precFloat64) then
398 read(dUnit,rec=irec) global
399 #ifdef _BYTESWAPIO
400 call MDS_BYTESWAPR8( Nx*Ny, global )
401 #endif
402 else
403 write(msgbuf,'(a)')
404 & ' MDSREADFIELD: illegal value for filePrec'
405 call print_error( msgbuf, mythid )
406 stop 'ABNORMAL END: S/R MDSREADFIELD'
407 endif
408 ENDIF
409 _END_MASTER( myThid )
410
411 CALL SCATTER_2D(global,local,mythid)
412 if (arrType .eq. 'RS') then
413 call PASStoRS( local,arr,k,nNz,mythid )
414 elseif (arrType .eq. 'RL') then
415 call PASStoRL( local,arr,k,nNz,mythid )
416 else
417 write(msgbuf,'(a)')
418 & ' MDSREADFIELD: illegal value for arrType'
419 call print_error( msgbuf, mythid )
420 stop 'ABNORMAL END: S/R MDSREADFIELD'
421 endif
422
423 ENDDO
424 c ENDDO k=1,nNz
425
426 _BEGIN_MASTER( myThid )
427 close( dUnit )
428 _END_MASTER( myThid )
429
430 endif
431 c endif ( globalFile .and. useSingleCPUIO )
432
433 return
434 end
435
436
437 C ==================================================================
438
439 subroutine passToRS(local,arr,k,nNz,mythid)
440 implicit none
441 #include "EEPARAMS.h"
442 #include "SIZE.h"
443 _RL local(1-oLx:sNx+oLx,1-oLy:sNy+oLy,nSx,nSy)
444 integer i,j,k,bi,bj,nNz
445 _RS arr(1-oLx:sNx+oLx,1-oLy:sNy+oLy,nNz,nSx,nSy)
446 integer mythid
447 DO bj = myByLo(myThid), myByHi(myThid)
448 DO bi = myBxLo(myThid), myBxHi(myThid)
449 DO J=1-Oly,sNy+Oly
450 DO I=1-Olx,sNx+Olx
451 arr(I,J,k,bi,bj) = local(I,J,bi,bj)
452 ENDDO
453 ENDDO
454 ENDDO
455 ENDDO
456 return
457 end
458
459 subroutine passToRL(local,arr,k,nNz,mythid)
460 implicit none
461 #include "EEPARAMS.h"
462 #include "SIZE.h"
463 _RL local(1-oLx:sNx+oLx,1-oLy:sNy+oLy,nSx,nSy)
464 integer i,j,k,bi,bj,nNz
465 _RL arr(1-oLx:sNx+oLx,1-oLy:sNy+oLy,nNz,nSx,nSy)
466 integer mythid
467 DO bj = myByLo(myThid), myByHi(myThid)
468 DO bi = myBxLo(myThid), myBxHi(myThid)
469 DO J=1-Oly,sNy+Oly
470 DO I=1-Olx,sNx+Olx
471 arr(I,J,k,bi,bj) = local(I,J,bi,bj)
472 ENDDO
473 ENDDO
474 ENDDO
475 ENDDO
476 return
477 end

  ViewVC Help
Powered by ViewVC 1.1.22