Ticket #182: def_floats.F

File def_floats.F, 21.4 KB (added by m.hadfield, 16 years ago)
Line 
1#include "cppdefs.h"
2#ifdef FLOATS
3 SUBROUTINE def_floats (ng,ldef)
4!
5!svn $Id: def_floats.F 139 2008-01-10 00:17:29Z arango $
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2008 The ROMS/TOMS Group !
8! Licensed under a MIT/X style license !
9! See License_ROMS.txt !
10!=======================================================================
11! !
12! This routine creates FLOATS NetCDF file, it defines dimensions, !
13! attributes, and variables. !
14! !
15!=======================================================================
16!
17 USE mod_param
18 USE mod_parallel
19 USE mod_floats
20# ifdef FOUR_DVAR
21 USE mod_fourdvar
22# endif
23 USE mod_grid
24 USE mod_iounits
25 USE mod_ncparam
26 USE mod_netcdf
27 USE mod_scalars
28# ifdef SEDIMENT
29 USE mod_sediment
30# endif
31# ifdef DISTRIBUTE
32!
33 USE distribute_mod, ONLY : mp_collect, mp_bcasti
34# endif
35!
36 implicit none
37!
38! Imported variable declarations.
39!
40 integer, intent(in) :: ng
41
42 logical, intent(in) :: ldef
43!
44! Local variable declarations.
45!
46 integer, parameter :: Natt = 24
47
48 logical :: got_var(-6:NV)
49
50 integer :: fltdim, i, itrc, j, l, nrec, nvd
51 integer :: recdim, status
52
53 integer :: DimIDs(29), fgrd(2), start(2), total(2)
54 integer :: Vsize(4)
55
56 integer :: def_info, def_var, wrt_info
57
58# ifdef DISTRIBUTE
59 integer :: Npts
60
61 real(r8), parameter :: Fspv = 0.0_r8
62 real(r8), dimension(Nfloats(ng)*NFV(ng)) :: Fwrk
63# endif
64
65 real(r8) :: Aval(6), Tinp(Nfloats(ng))
66
67 character (len=80) :: Vinfo(Natt)
68 character (len=80) :: dimnam, fname, ncname
69!
70!-----------------------------------------------------------------------
71! Set and report file name.
72!-----------------------------------------------------------------------
73!
74 IF (exit_flag.ne.NoError) RETURN
75 ncname=FLTname(ng)
76!
77 IF (Master) THEN
78 IF (ldef) THEN
79 WRITE (stdout,10) TRIM(ncname)
80 ELSE
81 WRITE (stdout,20) TRIM(ncname)
82 END IF
83 END IF
84!
85!=======================================================================
86! Create a new floats data file.
87!=======================================================================
88!
89! Create floats NetCDF file.
90!
91 IF (exit_flag.ne.NoError) RETURN
92 IF (ldef.and.OutThread) THEN
93 status=nf90_create(TRIM(ncname), nf90_clobber, ncFLTid(ng))
94 IF (status.ne.nf90_noerr) THEN
95 WRITE (stdout,30) TRIM(ncname)
96 exit_flag=3
97 ioerror=status
98 RETURN
99 END IF
100 END IF
101#ifdef DISTRIBUTE
102 IF (ldef) THEN
103 CALL mp_bcasti (ng, iNLM, ncFLTid(ng), 1)
104 END IF
105#endif
106!
107!-----------------------------------------------------------------------
108! Define the dimensions of staggered fields.
109!-----------------------------------------------------------------------
110!
111 IF (ldef.and.OutThread) THEN
112# ifdef SOLVE3D
113 status=nf90_def_dim(ncFLTid(ng),'s_rho' ,N(ng), &
114 & DimIDs( 9))
115 status=nf90_def_dim(ncFLTid(ng),'s_w' ,N(ng)+1, &
116 & DimIDs(10))
117 status=nf90_def_dim(ncFLTid(ng),'tracer' ,NT(ng), &
118 & DimIDs(11))
119# ifdef SEDIMENT
120 status=nf90_def_dim(ncFLTid(ng),'Nbed' ,Nbed, &
121 & DimIDs(16))
122# endif
123# ifdef ECOSIM
124 status=nf90_def_dim(ncFLTid(ng),'Nphy' ,Nphy, &
125 & DimIDs(25))
126 status=nf90_def_dim(ncFLTid(ng),'Nbac' ,Nbac, &
127 & DimIDs(26))
128 status=nf90_def_dim(ncFLTid(ng),'Ndom' ,Ndom, &
129 & DimIDs(27))
130 status=nf90_def_dim(ncFLTid(ng),'Nfec' ,Nfec, &
131 & DimIDs(28))
132# endif
133# endif
134 status=nf90_def_dim(ncFLTid(ng),'drifter' ,Nfloats(ng), &
135 & DimIDs(15))
136 status=nf90_def_dim(ncFLTid(ng),'boundary',4, &
137 & DimIDs(14))
138# ifdef FOUR_DVAR
139 status=nf90_def_dim(ncFLTid(ng),'Nstate ',NstateVar(ng), &
140 & DimIDs(29))
141# endif
142 status=nf90_def_dim(ncFLTid(ng),TRIM(ADJUSTL(Vname(5,idtime))), &
143 & nf90_unlimited,DimIDs(12))
144 recdim=DimIDs(12)
145 fltdim=DimIDs(15)
146!
147! Define dimension vectors for point variables.
148!
149 fgrd(1)=DimIDs(15)
150 fgrd(2)=DimIDs(12)
151!
152! Initialize unlimited time record dimension.
153!
154 tFLTindx(ng)=0
155!
156! Initialize local information variable arrays.
157!
158 DO i=1,Natt
159 DO j=1,80
160 Vinfo(i)(j:j)=' '
161 END DO
162 END DO
163 DO i=1,6
164 Aval(i)=0.0_r8
165 END DO
166!
167!-----------------------------------------------------------------------
168! Define time-recordless information variables.
169!-----------------------------------------------------------------------
170!
171 status=def_info(ng,ncFLTid(ng),ncname,DimIDs)
172 if (exit_flag.ne.NoError) RETURN
173!
174!-----------------------------------------------------------------------
175! Define variables and their attributes.
176!-----------------------------------------------------------------------
177!
178! Define model time.
179!
180 Vinfo( 1)=Vname(1,idtime)
181 Vinfo( 2)=Vname(2,idtime)
182 IF (INT(time_ref).eq.-2) THEN
183 Vinfo( 3)='seconds since 1968-05-23 00:00:00 GMT'
184 Vinfo( 4)='gregorian'
185 ELSE IF (INT(time_ref).eq.-1) THEN
186 Vinfo( 3)='seconds since 0001-01-01 00:00:00'
187 Vinfo( 4)='360_day'
188 ELSE IF (INT(time_ref).eq.0) THEN
189 Vinfo( 3)='seconds since 0001-01-01 00:00:00'
190 Vinfo( 4)='julian'
191 ELSE IF (time_ref.gt.0.0_r8) THEN
192 WRITE (Vinfo( 3),'(a,1x,a)') 'seconds since', TRIM(r_text)
193 Vinfo( 4)='gregorian'
194 END IF
195 Vinfo(14)=Vname(4,idtime)
196 status=def_var(ncFLTid(ng),fltVid(idtime,ng),NF_TYPE, &
197 & 1,recdim,Aval,Vinfo,ncname)
198!
199! Define floats X-grid locations.
200!
201 Vinfo( 1)='Xgrid'
202 Vinfo( 2)='x-grid floats locations'
203 Vinfo( 5)='valid_min'
204 Vinfo( 6)='valid_max'
205 Aval(2)=0.0_r8
206 Aval(3)=REAL(Lm(ng)+1,r8)
207 Vinfo(14)='Xgrid, scalar, series'
208 Vinfo(16)=Vname(1,idtime)
209 Vinfo(17)='missing_value'
210 Aval(4)=spval
211 status=def_var(ncFLTid(ng),fltVid(idXgrd,ng),NF_FOUT, &
212 & 2,fgrd,Aval,Vinfo,ncname)
213!
214! Define floats Y-grid locations.
215!
216 Vinfo( 1)='Ygrid'
217 Vinfo( 2)='Y-grid floats locations'
218 Vinfo( 5)='valid_min'
219 Vinfo( 6)='valid_max'
220 Aval(2)=0.0_r8
221 Aval(3)=REAL(Mm(ng)+1,r8)
222 Vinfo(14)='Ygrid, scalar, series'
223 Vinfo(16)=Vname(1,idtime)
224 Vinfo(17)='missing_value'
225 Aval(4)=spval
226 status=def_var(ncFLTid(ng),fltVid(idYgrd,ng),NF_FOUT, &
227 & 2,fgrd,Aval,Vinfo,ncname)
228# ifdef SOLVE3D
229!
230! Define floats Z-grid locations.
231!
232 Vinfo( 1)='Zgrid'
233 Vinfo( 2)='Z-grid floats locations'
234 Vinfo( 5)='valid_min'
235 Vinfo( 6)='valid_max'
236 Aval(2)=0.0_r8
237 Aval(3)=REAL(N(ng),r8)
238 Vinfo(14)='Zgrid, scalar, series'
239 Vinfo(16)=Vname(1,idtime)
240 Vinfo(17)='missing_value'
241 Aval(4)=spval
242 status=def_var(ncFLTid(ng),fltVid(idZgrd,ng),NF_FOUT, &
243 & 2,fgrd,Aval,Vinfo,ncname)
244# endif
245!
246! Define floats (lon,lat) or (x,y) locations.
247!
248 IF (spherical) THEN
249 Vinfo( 1)='lon'
250 Vinfo( 2)='longitude of floats trajectories'
251 Vinfo( 3)='degree_east'
252 Vinfo( 5)='valid_min'
253 Vinfo( 6)='valid_max'
254 Vinfo(14)='lon, scalar, series'
255 Vinfo(16)=Vname(1,idtime)
256 Vinfo(17)='missing_value'
257 Aval(2)=-180.0_r8
258 Aval(3)=180.0_r8
259 Aval(4)=spval
260 status=def_var(ncFLTid(ng),fltVid(idglon,ng),NF_FOUT, &
261 & 2,fgrd,Aval,Vinfo,ncname)
262 Vinfo( 1)='lat'
263 Vinfo( 2)='latitude of floats trajectories'
264 Vinfo( 3)='degree_north'
265 Vinfo( 5)='valid_min'
266 Vinfo( 6)='valid_max'
267 Vinfo(14)='lat, scalar, series'
268 Vinfo(16)=Vname(1,idtime)
269 Vinfo(17)='missing_value'
270 Aval(2)=-90.0_r8
271 Aval(3)=90.0_r8
272 Aval(4)=spval
273 status=def_var(ncFLTid(ng),fltVid(idglat,ng),NF_FOUT, &
274 & 2,fgrd,Aval,Vinfo,ncname)
275 ELSE
276 Vinfo( 1)='x'
277 Vinfo( 2)='x-location of floats trajectories'
278 Vinfo( 3)='meter'
279 Vinfo(14)='x, scalar, series'
280 Vinfo(16)=Vname(1,idtime)
281 Vinfo(17)='missing_value'
282 Aval(4)=spval
283 status=def_var(ncFLTid(ng),fltVid(idglon,ng),NF_FOUT, &
284 & 2,fgrd,Aval,Vinfo,ncname)
285 Vinfo( 1)='y'
286 Vinfo( 2)='y-location of floats trajectories'
287 Vinfo( 3)='meter'
288 Vinfo(14)='y, scalar, series'
289 Vinfo(16)=Vname(1,idtime)
290 Vinfo(17)='missing_value'
291 Aval(4)=spval
292 status=def_var(ncFLTid(ng),fltVid(idglat,ng),NF_FOUT, &
293 & 2,fgrd,Aval,Vinfo,ncname)
294 END IF
295# ifdef SOLVE3D
296!
297! Define floats depths.
298!
299 Vinfo( 1)='depth'
300 Vinfo( 2)='depth of floats trajectories'
301 Vinfo( 3)='meter'
302 Vinfo(14)='depth, scalar, series'
303 Vinfo(16)=Vname(1,idtime)
304 Vinfo(17)='missing_value'
305 Aval(4)=spval
306 status=def_var(ncFLTid(ng),fltVid(iddpth,ng),NF_FOUT, &
307 & 2,fgrd,Aval,Vinfo,ncname)
308!
309! Define density anomaly.
310!
311 Vinfo( 1)=Vname(1,idDano)
312 Vinfo( 2)=Vname(2,idDano)
313 Vinfo( 3)=Vname(3,idDano)
314 Vinfo(14)=Vname(4,idDano)
315 Vinfo(16)=Vname(1,idtime)
316 Vinfo(17)='missing_value'
317 Aval(4)=spval
318 status=def_var(ncFLTid(ng),fltVid(idDano,ng),NF_FOUT, &
319 & 2,fgrd,Aval,Vinfo,ncname)
320!
321! Define tracer type variables.
322!
323 DO itrc=1,NT(ng)
324 Vinfo( 1)=Vname(1,idTvar(itrc))
325 Vinfo( 2)=Vname(2,idTvar(itrc))
326 Vinfo( 3)=Vname(3,idTvar(itrc))
327 Vinfo(14)=Vname(4,idTvar(itrc))
328 Vinfo(16)=Vname(1,idtime)
329 Vinfo(17)='missing_value'
330 Aval(4)=spval
331# ifdef SEDIMENT
332 DO i=1,NST
333 IF (itrc.eq.idsed(i)) THEN
334 WRITE (Vinfo(19),40) 1000.0_r8*Sd50(i,ng)
335 END IF
336 END DO
337# endif
338 status=def_var(ncFLTid(ng),fltTid(itrc,ng),NF_FOUT, &
339 & 2,fgrd,Aval,Vinfo,ncname)
340 END DO
341# endif
342!
343! Initialize unlimited time record dimension.
344!
345 tFLTindx(ng)=0
346!
347!-----------------------------------------------------------------------
348! Leave definition mode.
349!-----------------------------------------------------------------------
350!
351 status=nf90_enddef(ncFLTid(ng))
352 END IF
353!
354!-----------------------------------------------------------------------
355! Write out time-recordless, information variables.
356!-----------------------------------------------------------------------
357!
358 IF (ldef) THEN
359 status=wrt_info(ng, iNLM, ncFLTid(ng), Master, ncname)
360 if (exit_flag.ne.NoError) RETURN
361 END IF
362!
363!=======================================================================
364! Open an existing floats file, check its contents, and prepare for
365! appending data.
366!=======================================================================
367!
368 IF (.not.ldef.and.OutThread) THEN
369!
370! Inquire about the contents of floats NetCDF file: Inquire about
371! the dimensions and variables. Check for consistency.
372!
373 ncname=FLTname(ng)
374 CALL opencdf (ng, 1, ncname, fname, N(ng), 0, nrec, nvd, Vsize)
375 if (exit_flag.ne.NoError) RETURN
376!
377! Open floats file for read/write.
378!
379 status=nf90_open(TRIM(ncname),nf90_write,ncFLTid(ng))
380 IF (status.ne.nf90_noerr) THEN
381 WRITE (stdout,50) TRIM(ncname)
382 exit_flag=3
383 ioerror=status
384 RETURN
385 END IF
386!
387! Initialize logical switches.
388!
389 DO i=1,NV
390 got_var(i)=.FALSE.
391 END DO
392!
393! Inquire size of the "drifter" dimension.
394!
395 status=nf90_inq_dimid(ncFLTid(ng),'drifter',fltdim)
396 IF (status.ne.nf90_noerr) THEN
397 WRITE (stdout,60) 'drifter', TRIM(ncname)
398 exit_flag=3
399 ioerror=status
400 RETURN
401 END IF
402 status=nf90_inquire_dimension(ncFLTid(ng),fltdim,dimnam, &
403 & Nfloats(ng))
404 IF (status.ne.nf90_noerr) THEN
405 WRITE (stdout,70) 'drifter', TRIM(ncname)
406 exit_flag=3
407 ioerror=status
408 RETURN
409 END IF
410!
411! Scan variable list from input NetCDF and activate switches for
412! float variables. Get variable IDs.
413!
414 DO i=1,nvars
415 IF (TRIM(varnam(i)).eq.TRIM(Vname(1,idtime))) THEN
416 got_var(idtime)=.TRUE.
417 status=nf90_inq_varid(ncFLTid(ng),TRIM(Vname(1,idtime)), &
418 & fltVid(idtime,ng))
419 END IF
420 IF (TRIM(varnam(i)).eq.'Xgrid') THEN
421 got_var(idXgrd)=.TRUE.
422 status=nf90_inq_varid(ncFLTid(ng),'Xgrid',fltVid(idXgrd,ng))
423 END IF
424 IF (TRIM(varnam(i)).eq.'Ygrid') THEN
425 got_var(idYgrd)=.TRUE.
426 status=nf90_inq_varid(ncFLTid(ng),'Ygrid',fltVid(idYgrd,ng))
427 END IF
428# ifdef SOLVE3D
429 IF (TRIM(varnam(i)).eq.'Zgrid') THEN
430 got_var(idZgrd)=.TRUE.
431 status=nf90_inq_varid(ncFLTid(ng),'Zgrid',fltVid(idZgrd,ng))
432 END IF
433# endif
434 IF (spherical) THEN
435 IF (TRIM(varnam(i)).eq.'lon') THEN
436 got_var(idglon)=.TRUE.
437 status=nf90_inq_varid(ncFLTid(ng),'lon',fltVid(idglon,ng))
438 END IF
439 IF (TRIM(varnam(i)).eq.'lat') THEN
440 got_var(idglat)=.TRUE.
441 status=nf90_inq_varid(ncFLTid(ng),'lat',fltVid(idglat,ng))
442 END IF
443 ELSE
444 IF (TRIM(varnam(i)).eq.'x') THEN
445 got_var(idglon)=.TRUE.
446 status=nf90_inq_varid(ncFLTid(ng),'x',fltVid(idglon,ng))
447 END IF
448 IF (TRIM(varnam(i)).eq.'y') THEN
449 got_var(idglat)=.TRUE.
450 status=nf90_inq_varid(ncFLTid(ng),'y',fltVid(idglat,ng))
451 END IF
452 END IF
453# ifdef SOLVE3D
454 IF (TRIM(varnam(i)).eq.'depth') THEN
455 got_var(iddpth)=.TRUE.
456 status=nf90_inq_varid(ncFLTid(ng),'depth',fltVid(iddpth,ng))
457 END IF
458 IF (TRIM(varnam(i)).eq.TRIM(Vname(1,idDano))) THEN
459 got_var(idDano)=.TRUE.
460 status=nf90_inq_varid(ncFLTid(ng),TRIM(Vname(1,idDano)), &
461 & fltVid(idDano,ng))
462 END IF
463 DO itrc=1,NT(ng)
464 IF (TRIM(varnam(i)).eq.TRIM(Vname(1,idTvar(itrc)))) THEN
465 got_var(idTvar(itrc))=.TRUE.
466 status=nf90_inq_varid(ncFLTid(ng), &
467 & TRIM(Vname(1,idTvar(itrc))), &
468 & fltTid(itrc,ng))
469 END IF
470 END DO
471# endif
472 END DO
473!
474! Check if floats variables are available in input NetCDF file.
475!
476 IF (.not.got_var(idtime)) THEN
477 WRITE (stdout,80) TRIM(Vname(1,idtime)), TRIM(ncname)
478 exit_flag=3
479 RETURN
480 END IF
481 IF (.not.got_var(idXgrd)) THEN
482 WRITE (stdout,80) 'Xgrid', TRIM(ncname)
483 exit_flag=3
484 RETURN
485 END IF
486 IF (.not.got_var(idYgrd)) THEN
487 WRITE (stdout,80) 'Ygrid', TRIM(ncname)
488 exit_flag=3
489 RETURN
490 END IF
491# ifdef SOLVE3D
492 IF (.not.got_var(idZgrd)) THEN
493 WRITE (stdout,80) 'Zgrid', TRIM(ncname)
494 exit_flag=3
495 RETURN
496 END IF
497# endif
498 IF (.not.got_var(idglon)) THEN
499 IF (spherical) THEN
500 WRITE (stdout,80) 'lon', TRIM(ncname)
501 ELSE
502 WRITE (stdout,80) 'x', TRIM(ncname)
503 END IF
504 exit_flag=3
505 RETURN
506 END IF
507 IF (.not.got_var(idglat)) THEN
508 IF (spherical) THEN
509 WRITE (stdout,80) 'lat', TRIM(ncname)
510 ELSE
511 WRITE (stdout,80) 'y', TRIM(ncname)
512 END IF
513 exit_flag=3
514 RETURN
515 END IF
516# ifdef SOLVE3D
517 IF (.not.got_var(iddpth)) THEN
518 WRITE (stdout,80) 'depth', TRIM(ncname)
519 exit_flag=3
520 RETURN
521 END IF
522 IF (.not.got_var(idDano)) THEN
523 WRITE (stdout,80) TRIM(Vname(1,idDano)), TRIM(ncname)
524 exit_flag=3
525 RETURN
526 END IF
527 DO itrc=1,NT(ng)
528 IF (.not.got_var(idTvar(itrc))) THEN
529 WRITE (stdout,80) TRIM(Vname(1,idTvar(itrc))), TRIM(ncname)
530 exit_flag=3
531 RETURN
532 END IF
533 END DO
534# endif
535!
536!-----------------------------------------------------------------------
537! Initialize floats positions to the appropriate values.
538!-----------------------------------------------------------------------
539!
540! Set-up floats time record.
541!
542 IF (frrec(ng).lt.0) THEN
543 tFLTindx(ng)=tsize
544 ELSE
545 tFLTindx(ng)=ABS(frrec(ng))
546!! tFLTindx(ng)=1+(ntstart(ng)-1)/nFLT(ng)
547 END IF
548 NrecFLT=tFLTindx(ng)
549!
550! Read in floats nondimentional horizontal positions.
551!
552 start(1)=1
553 total(1)=Nfloats(ng)
554 start(2)=tFLTindx(ng)
555 total(2)=1
556 status=nf90_get_var(ncFLTid(ng),fltVid(idXgrd,ng), &
557 & Tinp,start,total)
558 IF (status.ne.nf90_noerr) THEN
559 WRITE (stdout,90) 'Xgrid', tFLTindx(ng), TRIM(ncname)
560 exit_flag=3
561 ioerror=status
562 RETURN
563 END IF
564 DO l=1,Nfloats(ng)
565 IF ((Tinp(l).gt.REAL(Lm(ng)+1,r8)-0.5_r8).or. &
566 & (Tinp(l).lt.0.5_r8)) THEN
567 FLT(ng)%bounded(l)=.FALSE.
568 ELSE
569 FLT(ng)%bounded(l)=.TRUE.
570 DO i=0,NFT
571 FLT(ng)%track(ixgrd,i,l)=Tinp(l)
572 FLT(ng)%track(ixrhs,i,l)=0.0_r8
573 END DO
574 END IF
575 END DO
576 status=nf90_get_var(ncFLTid(ng),fltVid(idYgrd,ng), &
577 & Tinp,start,total)
578 IF (status.ne.nf90_noerr) THEN
579 WRITE (stdout,90) 'Ygrid', tFLTindx(ng), TRIM(ncname)
580 exit_flag=3
581 ioerror=status
582 RETURN
583 END IF
584 DO l=1,Nfloats(ng)
585 IF ((Tinp(l).gt.REAL(Mm(ng)+1,r8)-0.5_r8).or. &
586 & (Tinp(l).lt.0.5_r8)) THEN
587 FLT(ng)%bounded(l)=.FALSE.
588 ELSE
589 FLT(ng)%bounded(l)=.TRUE.
590 DO i=0,NFT
591 FLT(ng)%track(iygrd,i,l)=Tinp(l)
592 FLT(ng)%track(iyrhs,i,l)=0.0_r8
593 END DO
594 END IF
595 END DO
596# ifdef SOLVE3D
597 status=nf90_get_var(ncFLTid(ng),fltVid(idZgrd,ng), &
598 & Tinp,start,total)
599 IF (status.ne.nf90_noerr) THEN
600 WRITE (stdout,90) 'Zgrid', tFLTindx(ng), TRIM(ncname)
601 exit_flag=3
602 ioerror=status
603 RETURN
604 END IF
605 DO l=1,Nfloats(ng)
606 IF ((Tinp(l).gt.REAL(N(ng),r8)).or.(Tinp(l).lt.0.0_r8)) THEN
607 FLT(ng)%bounded(l)=.FALSE.
608 ELSE
609 FLT(ng)%bounded(l)=.TRUE.
610 DO i=0,NFT
611 FLT(ng)%track(izgrd,i,l)=Tinp(l)
612 FLT(ng)%track(izrhs,i,l)=0.0_r8
613 END DO
614 END IF
615 END DO
616# endif
617 END IF
618
619# ifdef DISTRIBUTE
620!
621! Broadcast restart float information to all nodes.
622!
623 IF (.not.ldef) THEN
624
625 Npts=NFV(ng)*(NFT+1)*Nfloats(ng)
626
627 CALL mp_collect (ng, iNLM, Npts, Fspv, FLT(ng)%track)
628!
629! Collect the bounded status switch.
630!
631 Fwrk=Fspv
632 DO l=1,Nfloats(ng)
633 IF (FLT(ng)%bounded(l)) THEN
634 Fwrk(l)=1.0_r8
635 END IF
636 END DO
637 CALL mp_collect (ng, iNLM, Nfloats(ng), Fspv, Fwrk)
638 DO l=1,Nfloats(ng)
639 IF (Fwrk(l).ne.Fspv) THEN
640 FLT(ng)%bounded(l)=.TRUE.
641 ELSE
642 FLT(ng)%bounded(l)=.FALSE.
643 END IF
644 END DO
645 END IF
646# endif
647!
648 10 FORMAT (6x,'DEF_FLOATS - creating floats file: ',a)
649 20 FORMAT (6x,'DEF_FLOATS - inquiring history file: ',a)
650 30 FORMAT (/,' DEF_FLOATS - unable to create floats NetCDF', &
651 & ' file: ',a)
652 40 FORMAT (1pe11.4,1x,'millimeter')
653 50 FORMAT (/,' DEF_FLOATS - unable to open floats NetCDF file: ',a)
654 60 FORMAT (/,' DEF_FLOATS - error while inquiring dimension ID ', &
655 & ' for: ',a,2x,/,13x,'in floats NetCDF file: ',a)
656 70 FORMAT (/,' DEF_FLOATS - error while inquiring size of', &
657 & ' dimension: ',a,2x,/,13x,'in floats NetCDF file: ',a)
658 80 FORMAT (/,' DEF_FLOATS - unable to find variable: ',a,2x, &
659 & ' in floats NetCDF file: ',a)
660 90 FORMAT (/,' DEF_FLOATS - error while reading variable: ',a,2x, &
661 & ' at time record = ', i6.6,/,13x, &
662 & ' in floats NetCDF file: ',a)
663 RETURN
664 END SUBROUTINE def_floats
665#else
666 SUBROUTINE def_floats
667 RETURN
668 END SUBROUTINE def_floats
669#endif