File: C:\NOAA\NEMS_11731\src\atmos\gfs\dyn\sigio_r_module.f

1     !-------------------------------------------------------------------------------
2     module sigio_r_module
3     !$$$  Module Documentation Block
4     !
5     ! Module:    sigio_r_module  API for global spectral sigma file random I/O
6     !   Prgmmr: Iredell          Org: W/NX23     Date: 1999-01-18
7     !
8     ! Abstract: This module provides an Application Program Interface extension
9     !   for performing I/O on the sigma restart file of the global spectral model.
10     !   Functions include opening, reading, writing, and closing as well as
11     !   allocating and deallocating data buffers used in the transfers.
12     !   The I/O performed here is random.
13     !   The transfers are limited to header records, data records,
14     !   surface data records, or specific levels of upper air data records.
15     !   See the documentation for sigio_module for sequential I/O.
16     !   
17     ! Program History Log:
18     !   1999-01-18  Mark Iredell
19     !
20     ! Modules Used:
21     !   sigio_module     API for global spectral sigma file I/O
22     !
23     ! Public Variables:
24     !
25     ! Public Defined Types:
26     !   sigio_dats       Sigma file surface data fields
27     !     hs                Real(sigio_realkind)(:) pointer to spectral
28     !                       coefficients of surface height in m
29     !     ps                Real(sigio_realkind)(:) pointer to spectral
30     !                       coefficients of log of surface pressure over 1 kPa
31     !
32     !   sigio_datm       Sigma file multilevel data fields
33     !     k1                Integer(sigio_intkind) first level number
34     !     k2                Integer(sigio_intkind) last level number
35     !     t                 Real(sigio_realkind)(:,:) pointer to spectral
36     !                       coefficients of virtual temperature by level in K
37     !     d                 Real(sigio_realkind)(:,:) pointer to spectral
38     !                       coefficients of divergence by level in 1/second
39     !     z                 Real(sigio_realkind)(:,:) pointer to spectral
40     !                       coefficients of vorticity by level in 1/second
41     !     q                 Real(sigio_realkind)(:,:,:) pointer to spectral
42     !                       coefficients of tracers by tracer number and level
43     !                       in specific densities
44     !
45     !   sigio_dati       Sigma file single data field
46     !     i                 Integer(sigio_intkind) record index
47     !     f                 Real(sigio_realkind)(:) pointer to field
48     !                       
49     !   sigio_dbts       Sigma file longreal surface data fields
50     !     hs                Real(sigio_dblekind)(:) pointer to spectral
51     !                       coefficients of surface height in m
52     !     ps                Real(sigio_dblekind)(:) pointer to spectral
53     !                       coefficients of log of surface pressure over 1 kPa
54     !
55     !   sigio_dbtm       Sigma file longreal multilevel data fields
56     !     k1                Integer(sigio_intkind) first level number
57     !     k2                Integer(sigio_intkind) last level number
58     !     t                 Real(sigio_dblekind)(:,:) pointer to spectral
59     !                       coefficients of virtual temperature by level in K
60     !     d                 Real(sigio_dblekind)(:,:) pointer to spectral
61     !                       coefficients of divergence by level in 1/second
62     !     z                 Real(sigio_dblekind)(:,:) pointer to spectral
63     !                       coefficients of vorticity by level in 1/second
64     !     q                 Real(sigio_dblekind)(:,:,:) pointer to spectral
65     !                       coefficients of tracers by tracer number and level
66     !                       in specific densities
67     !
68     !   sigio_dbti       Sigma file longreal single data field
69     !     i                 Integer(sigio_intkind) record index
70     !     f                 Real(sigio_dblekind)(:) pointer to field
71     !                       
72     ! Public Subprograms:
73     !   sigio_rropen     Open sigma file for random reading
74     !     lu                Integer(sigio_intkind) input logical unit
75     !     cfname            Character(*) input filename
76     !     iret              Integer(sigio_intkind) output return code
77     !
78     !   sigio_rwopen     Open sigma file for random writing
79     !     lu                Integer(sigio_intkind) input logical unit
80     !     cfname            Character(*) input filename
81     !     iret              Integer(sigio_intkind) output return code
82     !
83     !   sigio_rxopen     Open sigma file for random reading and writing
84     !     lu                Integer(sigio_intkind) input logical unit
85     !     cfname            Character(*) input filename
86     !     iret              Integer(sigio_intkind) output return code
87     !
88     !   sigio_rclose     Close sigma file for random I/O
89     !     lu                Integer(sigio_intkind) input logical unit
90     !     iret              Integer(sigio_intkind) output return code
91     !
92     !   sigio_rrhead     Read header information with random I/O
93     !     lu                Integer(sigio_intkind) input logical unit
94     !     head              Type(sigio_head) output header information
95     !     iret              Integer(sigio_intkind) output return code
96     !
97     !   sigio_rwhead     Write header information with random I/O
98     !     lu                Integer(sigio_intkind) input logical unit
99     !     head              Type(sigio_head) input header information
100     !     iret              Integer(sigio_intkind) output return code
101     !
102     !   sigio_aldats     Allocate surface data fields
103     !     head              Type(sigio_head) input header information
104     !     dats              Type(sigio_dats) output surface data fields
105     !     iret              Integer(sigio_intkind) output return code
106     !
107     !   sigio_axdats     Deallocate surface data fields
108     !     dats              Type(sigio_dats) output surface data fields
109     !     iret              Integer(sigio_intkind) output return code
110     !
111     !   sigio_aldatm     Allocate multilevel data fields
112     !     head              Type(sigio_head) input header information
113     !     k1                Integer(sigio_intkind) input first level number
114     !     k2                Integer(sigio_intkind) input last level number
115     !     datm              Type(sigio_datm) output multilevel data fields
116     !     iret              Integer(sigio_intkind) output return code
117     !
118     !   sigio_axdatm     Deallocate multilevel data fields
119     !     datm              Type(sigio_datm) output multilevel data fields
120     !     iret              Integer(sigio_intkind) output return code
121     !
122     !   sigio_aldati     Allocate single data fields
123     !     head              Type(sigio_head) input header information
124     !     i                 Integer(sigio_intkind) input record index
125     !     dati              Type(sigio_dati) output single data field
126     !     iret              Integer(sigio_intkind) output return code
127     !
128     !   sigio_axdati     Deallocate single data fields
129     !     dati              Type(sigio_dati) output single data field
130     !     iret              Integer(sigio_intkind) output return code
131     !
132     !   sigio_rrdata     Read data fields with random I/O
133     !     lu                Integer(sigio_intkind) input logical unit
134     !     head              Type(sigio_head) input header information
135     !     data              Type(sigio_data) output data fields
136     !     iret              Integer(sigio_intkind) output return code
137     !
138     !   sigio_rwdata     Write data fields with random I/O
139     !     lu                Integer(sigio_intkind) input logical unit
140     !     head              Type(sigio_head) input header information
141     !     data              Type(sigio_data) input data fields
142     !     iret              Integer(sigio_intkind) output return code
143     !
144     !   sigio_rrohdc     Open, read header & data and close with random I/O
145     !     lu                Integer(sigio_intkind) input logical unit
146     !     cfname            Character(*) input filename
147     !     head              Type(sigio_head) output header information
148     !     data              Type(sigio_data) or type(sigio_dbta) output data fields
149     !     iret              Integer(sigio_intkind) output return code
150     !
151     !   sigio_rwohdc     Open, write header & data and close with random I/O
152     !     lu                Integer(sigio_intkind) input logical unit
153     !     cfname            Character(*) input filename
154     !     head              Type(sigio_head) input header information
155     !     data              Type(sigio_data) or type(sigio_dbta) input data fields
156     !     iret              Integer(sigio_intkind) output return code
157     !
158     !   sigio_rrdats     Read surface data fields with random I/O
159     !     lu                Integer(sigio_intkind) input logical unit
160     !     head              Type(sigio_head) input header information
161     !     dats              Type(sigio_dats) output surface data fields
162     !     iret              Integer(sigio_intkind) output return code
163     !
164     !   sigio_rwdats     Write surface data fields with random I/O
165     !     lu                Integer(sigio_intkind) input logical unit
166     !     head              Type(sigio_head) input header information
167     !     dats              Type(sigio_dats) input surface data fields
168     !     iret              Integer(sigio_intkind) output return code
169     !
170     !   sigio_rrdatm     Read multilevel data fields with random I/O
171     !     lu                Integer(sigio_intkind) input logical unit
172     !     head              Type(sigio_head) input header information
173     !     datm              Type(sigio_datm) output multilevel data fields
174     !     iret              Integer(sigio_intkind) output return code
175     !
176     !   sigio_rwdatm     Write multilevel data fields with random I/O
177     !     lu                Integer(sigio_intkind) input logical unit
178     !     head              Type(sigio_head) input header information
179     !     datm              Type(sigio_datm) input multilevel data fields
180     !     iret              Integer(sigio_intkind) output return code
181     !
182     !   sigio_rrdati     Read single data field with random I/O
183     !     lu                Integer(sigio_intkind) input logical unit
184     !     head              Type(sigio_head) input header information
185     !     dati              Type(sigio_dati) output single data field
186     !     iret              Integer(sigio_intkind) output return code
187     !
188     !   sigio_rwdati     Write single data field with random I/O
189     !     lu                Integer(sigio_intkind) input logical unit
190     !     head              Type(sigio_head) input header information
191     !     dati              Type(sigio_dati) input single data field
192     !     iret              Integer(sigio_intkind) output return code
193     !
194     !   sigio_aldbts     Allocate longreal surface data fields
195     !     head              Type(sigio_head) input header information
196     !     dbts              Type(sigio_dbts) output longreal surface data fields
197     !     iret              Integer(sigio_intkind) output return code
198     !
199     !   sigio_axdbts     Deallocate longreal surface data fields
200     !     dbts              Type(sigio_dbts) output longreal surface data fields
201     !     iret              Integer(sigio_intkind) output return code
202     !
203     !   sigio_aldbtm     Allocate longreal multilevel data fields
204     !     head              Type(sigio_head) input header information
205     !     k                 Integer(sigio_intkind) input level number
206     !     dbtm              Type(sigio_dbtm) output longreal multilevel data fields
207     !     iret              Integer(sigio_intkind) output return code
208     !
209     !   sigio_axdbtm     Deallocate longreal multilevel data fields
210     !     dbtm              Type(sigio_dbtm) output longreal multilevel data fields
211     !     iret              Integer(sigio_intkind) output return code
212     !
213     !   sigio_aldbti     Allocate longreal single data fields
214     !     head              Type(sigio_head) input header information
215     !     i                 Integer(sigio_intkind) input record index
216     !     dbti              Type(sigio_dbti) output longreal single data field
217     !     iret              Integer(sigio_intkind) output return code
218     !
219     !   sigio_axdbti     Deallocate longreal single data fields
220     !     dbti              Type(sigio_dbti) output longreal single data field
221     !     iret              Integer(sigio_intkind) output return code
222     !
223     !   sigio_rrdbta     Read longreal data fields with random I/O
224     !     lu                Integer(sigio_intkind) input logical unit
225     !     head              Type(sigio_head) input header information
226     !     dbta              Type(sigio_dbta) output longreal data fields
227     !     iret              Integer(sigio_intkind) output return code
228     !
229     !   sigio_rwdbta     Write longreal data fields with random I/O
230     !     lu                Integer(sigio_intkind) input logical unit
231     !     head              Type(sigio_head) input header information
232     !     dbta              Type(sigio_dbta) input longreal data fields
233     !     iret              Integer(sigio_intkind) output return code
234     !
235     !   sigio_rrdbts     Read longreal surface data fields with random I/O
236     !     lu                Integer(sigio_intkind) input logical unit
237     !     head              Type(sigio_head) input header information
238     !     dbts              Type(sigio_dbts) output longreal surface data fields
239     !     iret              Integer(sigio_intkind) output return code
240     !
241     !   sigio_rwdbts     Write longreal surface data fields with random I/O
242     !     lu                Integer(sigio_intkind) input logical unit
243     !     head              Type(sigio_head) input header information
244     !     dbts              Type(sigio_dbts) input longreal surface data fields
245     !     iret              Integer(sigio_intkind) output return code
246     !
247     !   sigio_rrdbtm     Read longreal multilevel data fields with random I/O
248     !     lu                Integer(sigio_intkind) input logical unit
249     !     head              Type(sigio_head) input header information
250     !     dbtm              Type(sigio_dbtm) output longreal multilevel data fields
251     !     iret              Integer(sigio_intkind) output return code
252     !
253     !   sigio_rwdbtm     Write longreal multilevel data fields with random I/O
254     !     lu                Integer(sigio_intkind) input logical unit
255     !     head              Type(sigio_head) input header information
256     !     dbtm              Type(sigio_dbtm) input longreal multilevel data fields
257     !     iret              Integer(sigio_intkind) output return code
258     !
259     !   sigio_rrdbti     Read longreal single data field with random I/O
260     !     lu                Integer(sigio_intkind) input logical unit
261     !     head              Type(sigio_head) input header information
262     !     dbti              Type(sigio_dbti) output longreal single data field
263     !     iret              Integer(sigio_intkind) output return code
264     !
265     !   sigio_rwdbti     Write longreal single data field with random I/O
266     !     lu                Integer(sigio_intkind) input logical unit
267     !     head              Type(sigio_head) input header information
268     !     dbti              Type(sigio_dbti) input longreal single data field
269     !     iret              Integer(sigio_intkind) output return code
270     !
271     ! Subprograms called:
272     !   baopenr           Byte-addressable open for reading
273     !   baopenw           Byte-addressable open for writing
274     !   baclose           Byte-addressable close
275     !   bafrindex         Byte-addressable Fortran record index
276     !   bafrread          Byte-addressable Fortran record read
277     !   bafrwrite         Byte-addressable Fortran record write
278     !
279     ! Remarks:
280     !   (1) The sigma file format follows:
281     !         ON85 label (32 bytes)
282     !         Header information record containing
283     !           real forecast hour, initial date, sigma interfaces, sigma levels,
284     !           padding to allow for 100 levels, and finally 44 identifier words
285     !           containing JCAP, LEVS, NTRAC, etc. (250 4-byte words)
286     !         Orography (NC 4-byte words, where NC=(JCAP+1)*(JCAP+2))
287     !         Log surface pressure (NC 4-byte words)
288     !         Temperature (LEVS records of NC 4-byte words)
289     !         Divergence & Vorticity interleaved (2*LEVS records of NC 4-byte words)
290     !         Tracers (LEVS*NTRAC records of NC 4-byte words)
291     !
292     !   (2) Possible return codes:
293     !          0   Successful call
294     !         -1   Open or close I/O error
295     !         -2   Header record I/O error (possible EOF)
296     !         -3   Allocation or deallocation error
297     !         -4   Data record I/O error
298     !         -5   Insufficient data dimensions allocated
299     !
300     ! Examples:
301     !   (1) Write out orography and surface pressure only from processor 0:
302     !
303     !     subroutine write_surface_fields(me,head,len,orog,lnps)
304     !     use sigio_r_module
305     !     integer,intent(in):: me
306     !     type(sigio_head),intent(in):: head
307     !     integer,intent(in):: len
308     !     real(sigio_dblekind),intent(in),target:: orog(len),lnps(len)
309     !     type(sigio_dbts) dbts
310     !     integer iret
311     !     if(me.eq.0) then
312     !       dbts%hs=>orog
313     !       dbts%ps=>lnps
314     !       call sigio_rwdbts(51,head,dbts,iret)
315     !     endif
316     !     end subroutine
317     ! 
318     ! Attributes:
319     !   Language: Fortran 90
320     !
321     !$$$
322       use sigio_module
323       implicit none
324       private
325     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
326     ! Public Variables
327     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
328     ! Public Types
329       type,public:: sigio_dats
330         real(sigio_realkind),pointer:: hs(:),ps(:)
331       end type
332       type,public:: sigio_datm
333         integer(sigio_intkind):: k1,k2
334         real(sigio_realkind),pointer:: t(:,:),d(:,:),z(:,:)
335         real(sigio_realkind),pointer:: q(:,:,:)
336       end type
337       type,public:: sigio_dati
338         integer(sigio_intkind):: i
339         real(sigio_realkind),pointer:: f(:)
340       end type
341       type,public:: sigio_dbts
342         real(sigio_dblekind),pointer:: hs(:),ps(:)
343       end type
344       type,public:: sigio_dbtm
345         integer(sigio_intkind):: k1,k2
346         real(sigio_dblekind),pointer:: t(:,:),d(:,:),z(:,:)
347         real(sigio_dblekind),pointer:: q(:,:,:)
348       end type
349       type,public:: sigio_dbti
350         integer(sigio_intkind):: i
351         real(sigio_dblekind),pointer:: f(:)
352       end type
353     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
354     ! Private Variables
355     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
356     ! Private Types
357       type sigio_head2
358         sequence
359         real(sigio_realkind):: fhour
360         integer(sigio_intkind):: idate(4)
361         real(sigio_realkind):: sisl(2*100+1)
362         real(sigio_realkind):: ext(44)
363       end type
364       type sigio_head1a
365         sequence
366         character(8):: clab8
367         integer(sigio_intkind):: ivs,nhead,ndata,reserved(3)
368       end type
369       type sigio_head3a
370         sequence
371         real(sigio_realkind) fhour
372         integer(sigio_intkind):: idate(4)
373         integer(sigio_intkind):: jcap,levs,&
374           itrun,iorder,irealf,igen,latf,lonf,&
375           latb,lonb,latr,lonr,ntrac,nvcoord,&
376           icen2,iens(2),idpp,idsl,idvc,idvm,&
377           idvt,idrun,idusr
378         real(sigio_realkind) pdryini
379         integer(sigio_intkind):: ncldt,ixgr
380         integer(sigio_intkind):: reserved(18)
381       end type
382     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
383     ! Public Subprograms
384       public sigio_rropen,sigio_rwopen,sigio_rxopen,sigio_rclose
385       public sigio_rrhead,sigio_rwhead
386       public sigio_aldats,sigio_axdats
387       public sigio_aldatm,sigio_axdatm
388       public sigio_aldati,sigio_axdati
389       public sigio_rrdata,sigio_rwdata
390       public sigio_rrohdc,sigio_rwohdc
391       public sigio_rrdats,sigio_rwdats
392       public sigio_rrdatm,sigio_rwdatm
393       public sigio_rrdati,sigio_rwdati
394       public sigio_aldbts,sigio_axdbts
395       public sigio_aldbtm,sigio_axdbtm
396       public sigio_aldbti,sigio_axdbti
397       public sigio_rrdbta,sigio_rwdbta
398       public sigio_rrdbts,sigio_rwdbts
399       public sigio_rrdbtm,sigio_rwdbtm
400       public sigio_rrdbti,sigio_rwdbti
401       interface sigio_rrohdc
402         module procedure sigio_rrohdca,sigio_rrohdcb
403       end interface
404       interface sigio_rwohdc
405         module procedure sigio_rwohdca,sigio_rwohdcb
406       end interface
407     contains
408     !-------------------------------------------------------------------------------
409       subroutine sigio_rropen(lu,cfname,iret)
410         implicit none
411         integer(sigio_intkind),intent(in):: lu
412         character*(*),intent(in):: cfname
413         integer(sigio_intkind),intent(out):: iret
414     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
415         call baopenr(lu,cfname,iret)
416         if(iret.ne.0) iret=-1
417     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
418       end subroutine
419     !-------------------------------------------------------------------------------
420       subroutine sigio_rwopen(lu,cfname,iret)
421         implicit none
422         integer(sigio_intkind),intent(in):: lu
423         character*(*),intent(in):: cfname
424         integer(sigio_intkind),intent(out):: iret
425     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
426         call baopenw(lu,cfname,iret)
427         if(iret.ne.0) iret=-1
428     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
429       end subroutine
430     !-------------------------------------------------------------------------------
431       subroutine sigio_rxopen(lu,cfname,iret)
432         implicit none
433         integer(sigio_intkind),intent(in):: lu
434         character*(*),intent(in):: cfname
435         integer(sigio_intkind),intent(out):: iret
436     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
437         call baopen(lu,cfname,iret)
438         if(iret.ne.0) iret=-1
439     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
440       end subroutine
441     !-------------------------------------------------------------------------------
442       subroutine sigio_rclose(lu,iret)
443         implicit none
444         integer(sigio_intkind),intent(in):: lu
445         integer(sigio_intkind),intent(out):: iret
446     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
447         call baclose(lu,iret)
448         if(iret.ne.0) iret=-1
449     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
450       end subroutine
451     !-------------------------------------------------------------------------------
452       subroutine sigio_rrhead(lu,head,iret)
453         implicit none
454         integer(sigio_intkind),intent(in):: lu
455         type(sigio_head),intent(inout):: head
456         integer(sigio_intkind),intent(out):: iret
457         type(sigio_head2):: head2
458         type(sigio_head1a):: head1a
459         type(sigio_head3a):: head3a
460         integer:: iskip,iread,nread
461     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
462         iret=-2
463         iskip=0
464         iread=sigio_lhead1
465         call bafrread(lu,iskip,iread,nread,head1a)
466         if(nread.lt.iread) return
467     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
468         if(head1a%clab8.eq.'GFS SIG ') then  ! modern sigma file
469           head%ivs=head1a%ivs
470           call bafrindex(lu,iskip+nread,nread,iskip)
471           iread=200
472           call bafrread(lu,iskip,iread,nread,head3a)
473           if(nread.lt.iread) return
474           head%fhour=head3a%fhour
475           head%idate=head3a%idate
476           head%jcap=head3a%jcap
477           head%levs=head3a%levs
478           head%itrun=head3a%itrun
479           head%iorder=head3a%iorder
480           head%irealf=head3a%irealf
481           head%igen=head3a%igen
482           head%latf=head3a%latf
483           head%lonf=head3a%lonf
484           head%latb=head3a%latb
485           head%lonb=head3a%lonb
486           head%latr=head3a%latr
487           head%lonr=head3a%lonr
488           head%ntrac=head3a%ntrac
489           head%nvcoord=head3a%nvcoord
490           head%icen2=head3a%icen2
491           head%iens=head3a%iens
492           head%idpp=head3a%idpp
493           head%idsl=head3a%idsl
494           head%idvc=head3a%idvc
495           head%idvm=head3a%idvm
496           head%idvt=head3a%idvt
497           head%idrun=head3a%idrun
498           head%idusr=head3a%idusr
499           head%pdryini=head3a%pdryini
500           head%ncldt=head3a%ncldt
501           head%ixgr=head3a%ixgr
502           call sigio_alhead(head,iret)
503           iskip=iskip+nread
504           iread=4*size(head%vcoord)
505           call bafrread(lu,iskip,iread,nread,head%vcoord)
506           if(nread.lt.iread) return
507           iskip=iskip+nread
508           iread=size(head%cfvars)
509           call bafrread(lu,iskip,iread,nread,head%cfvars)
510           if(nread.lt.iread) return
511     !
512           if (mod(head%idvm/10,10) == 3) then
513             iskip=iskip+nread
514             iread=4*size(head%cpi)
515             call bafrread(lu,iskip,iread,nread,head%cpi)
516             if(nread.lt.iread) return
517             iskip=iskip+nread
518             iread=4*size(head%ri)
519             call bafrread(lu,iskip,iread,nread,head%ri)
520             if(nread.lt.iread) return
521           endif
522           head%clabsig=' '
523           head%si=sigio_realfill
524           head%sl=sigio_realfill
525           head%ak=sigio_realfill
526           head%bk=sigio_realfill
527           head%pdryini=sigio_realfill
528     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
529         else
530           iskip=0
531           iread=sigio_lhead1
532           call bafrread(lu,iskip,iread,nread,head%clabsig)
533           if(nread.lt.iread) return
534           iskip=iskip+nread
535           iread=1000
536           call bafrread(lu,iskip,iread,nread,head2)
537           if(nread.lt.iread) return
538           iret=0
539           head%fhour=head2%fhour
540           head%idate=head2%idate
541           head%jcap=head2%ext(1)
542           head%levs=head2%ext(2)
543           head%itrun=head2%ext(3)
544           head%iorder=head2%ext(4)
545           head%irealf=head2%ext(5)
546           head%igen=head2%ext(6)
547           head%lonf=head2%ext(7)
548           head%latf=head2%ext(8)
549           head%lonb=head2%ext(9)
550           head%latb=head2%ext(10)
551           head%lonr=head2%ext(11)
552           head%latr=head2%ext(12)
553           head%ntrac=max(head2%ext(13),1.)
554           head%icen2=head2%ext(14)
555           head%iens=head2%ext(15:16)
556           head%idpp=head2%ext(17)
557           head%idsl=head2%ext(18)
558           head%idvc=head2%ext(19)
559           head%idvm=head2%ext(20)
560           head%idvt=head2%ext(21)
561           head%idrun=head2%ext(22)
562           head%idusr=head2%ext(23)
563           head%pdryini=head2%ext(24)
564           head%ncldt=head2%ext(25)
565           head%si=sigio_realfill
566           head%sl=sigio_realfill
567           head%ak=sigio_realfill
568           head%bk=sigio_realfill
569           if(head%idvc.eq.2) then
570             head%ak(1:head%levs+1)=head2%sisl(1:head%levs+1)
571             head%bk(1:head%levs+1)=head2%sisl(head%levs+2:2*head%levs+2)
572           else
573             head%si(1:head%levs+1)=head2%sisl(1:head%levs+1)
574             head%sl(1:head%levs)=head2%sisl(head%levs+2:2*head%levs+1)
575           endif
576           head%ivs=198410
577           if(head%idvc.eq.2) then
578             head%nvcoord=2
579             call sigio_alhead(head,iret)
580             head%vcoord(1:head%levs+1,1)=head%ak(1:head%levs+1)
581             head%vcoord(1:head%levs+1,2)=head%bk(1:head%levs+1)
582           elseif(head%idvc.eq.0.or.head%idvc.eq.1) then
583             head%nvcoord=1
584             call sigio_alhead(head,iret)
585             head%vcoord(1:head%levs+1,1)=head%si(1:head%levs+1)
586           endif
587         endif
588     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
589         call sigio_adhead(head)
590         iret=0
591       end subroutine
592     !-------------------------------------------------------------------------------
593       subroutine sigio_rwhead(lu,head,iret)
594         implicit none
595         integer(sigio_intkind),intent(in):: lu
596         type(sigio_head),intent(inout):: head
597         integer(sigio_intkind),intent(out):: iret
598         type(sigio_head2):: head2
599         type(sigio_head1a):: head1a
600         integer,allocatable:: head2a(:)
601         type(sigio_head3a):: head3a
602         integer:: iskip,iwrite,nwrite
603     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
604         iret=-2
605         call sigio_adhead(head)
606         if(head%ivs.ge.200509) then
607           head1a%clab8='GFS SIG '
608           head1a%ivs=head%ivs
609           head1a%nhead=head%nhead
610           head1a%ndata=head%ndata
611           head1a%reserved=0
612           iskip=0
613           iwrite=head%lhead(1)
614           call bafrwrite(lu,iskip,iwrite,nwrite,head1a)
615           if(nwrite.lt.iwrite) return
616           allocate(head2a(head%nhead+head%ndata))
617           head2a(:head%nhead)=head%lhead
618           head2a(head%nhead+1:)=head%ldata
619           iskip=iskip+nwrite
620           iwrite=head%lhead(2)
621           call bafrwrite(lu,iskip,iwrite,nwrite,head2a)
622           deallocate(head2a)
623           if(nwrite.lt.iwrite) return
624           head3a%fhour=head%fhour
625           head3a%idate=head%idate
626           head3a%jcap=head%jcap
627           head3a%levs=head%levs
628           head3a%itrun=head%itrun
629           head3a%iorder=head%iorder
630           head3a%irealf=head%irealf
631           head3a%igen=head%igen
632           head3a%latf=head%latf
633           head3a%lonf=head%lonf
634           head3a%latb=head%latb
635           head3a%lonb=head%lonb
636           head3a%latr=head%latr
637           head3a%lonr=head%lonr
638           head3a%ntrac=head%ntrac
639           head3a%nvcoord=head%nvcoord
640           head3a%icen2=head%icen2
641           head3a%iens=head%iens
642           head3a%idpp=head%idpp
643           head3a%idsl=head%idsl
644           head3a%idvc=head%idvc
645           head3a%idvm=head%idvm
646           head3a%idvt=head%idvt
647           head3a%idrun=head%idrun
648           head3a%idusr=head%idusr
649           head3a%pdryini=head%pdryini
650           head3a%ncldt=head%ncldt
651           head3a%ixgr=head%ixgr
652           head3a%reserved=0
653           iskip=iskip+nwrite
654           iwrite=head%lhead(3)
655           call bafrwrite(lu,iskip,iwrite,nwrite,head3a)
656           if(nwrite.lt.iwrite) return
657           iskip=iskip+nwrite
658           iwrite=head%lhead(4)
659           call bafrwrite(lu,iskip,iwrite,nwrite,head%vcoord)
660           if(nwrite.lt.iwrite) return
661           iskip=iskip+nwrite
662           iwrite=head%lhead(5)
663           call bafrwrite(lu,iskip,iwrite,nwrite,head%cfvars)
664           if(nwrite.lt.iwrite) return
665     !
666           if (mod(head%idvm/10,10) == 3) then
667             iskip=iskip+nwrite
668             iwrite=head%lhead(7)
669             call bafrwrite(lu,iskip,iwrite,nwrite,head%cpi)
670             if(nwrite.lt.iwrite) return
671             iskip=iskip+nwrite
672             iwrite=head%lhead(7)
673             call bafrwrite(lu,iskip,iwrite,nwrite,head%ri)
674             if(nwrite.lt.iwrite) return
675           endif
676     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
677         else
678           iskip=0
679           iwrite=sigio_lhead1
680           call bafrwrite(lu,iskip,iwrite,nwrite,head%clabsig)
681           if(nwrite.lt.iwrite) return
682           head2%fhour=head%fhour 
683           head2%idate=head%idate 
684           head2%sisl=0
685           if(head%idvc.eq.2) then
686             if(head%nvcoord.eq.2.and.head%vcoord(1,2).eq.1.) then
687               head2%sisl(1:head%levs+1)=head%vcoord(1:head%levs+1,1)
688               head2%sisl(head%levs+2:2*head%levs+2)=head%vcoord(1:head%levs+1,2)
689             else
690               head2%sisl(1:head%levs+1)=head%ak(1:head%levs+1)
691               head2%sisl(head%levs+2:2*head%levs+2)=head%bk(1:head%levs+1)
692             endif
693           elseif(head%idvc.eq.0.or.head%idvc.eq.1) then
694             if(head%nvcoord.eq.1.and.head%vcoord(1,1).eq.1.) then
695               head2%sisl(1:head%levs+1)=head%vcoord(1:head%levs+1,1)
696               call sigio_modpr(1,1,head%levs,head%nvcoord,head%idvc,head%idsl,&
697                                head%vcoord,iret,ps=(/1./),&
698                                pm=head2%sisl(head%levs+2:2*head%levs+1))
699             else
700               head2%sisl(1:head%levs+1)=head%si(1:head%levs+1)
701               head2%sisl(head%levs+2:2*head%levs+1)=head%sl(1:head%levs)
702             endif
703           endif
704           head2%ext(1)=head%jcap
705           head2%ext(2)=head%levs
706           head2%ext(3)=head%itrun
707           head2%ext(4)=head%iorder
708           head2%ext(5)=head%irealf
709           head2%ext(6)=head%igen
710           head2%ext(7)=head%lonf
711           head2%ext(8)=head%latf
712           head2%ext(9)=head%lonb
713           head2%ext(10)=head%latb
714           head2%ext(11)=head%lonr
715           head2%ext(12)=head%latr
716           head2%ext(13)=head%ntrac
717           head2%ext(14)=head%icen2
718           head2%ext(15:16)=head%iens
719           head2%ext(17)=head%idpp
720           head2%ext(18)=head%idsl
721           head2%ext(19)=head%idvc
722           head2%ext(20)=head%idvm
723           head2%ext(21)=head%idvt
724           head2%ext(22)=head%idrun
725           head2%ext(23)=head%idusr
726           head2%ext(24)=head%pdryini
727           head2%ext(25)=head%ncldt
728           head2%ext(26:44)=0
729           iskip=iskip+nwrite
730           iwrite=1000
731           call bafrwrite(lu,iskip,iwrite,nwrite,head2)
732           if(nwrite.lt.iwrite) return
733         endif
734     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
735         iret=0
736       end subroutine
737     !-------------------------------------------------------------------------------
738       subroutine sigio_aldats(head,dats,iret)
739         implicit none
740         type(sigio_head),intent(in):: head
741         type(sigio_dats),intent(inout):: dats
742         integer(sigio_intkind),intent(out):: iret
743         integer nc,dim1
744     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
745         call sigio_axdats(dats,iret)
746         nc=(head%jcap+1)*(head%jcap+2)
747         dim1=nc
748         allocate(dats%hs(dim1),dats%ps(dim1),stat=iret)
749         if(iret.ne.0) iret=-3
750     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
751       end subroutine
752     !-------------------------------------------------------------------------------
753       subroutine sigio_axdats(dats,iret)
754         implicit none
755         type(sigio_dats),intent(inout):: dats
756         integer(sigio_intkind),intent(out):: iret
757     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
758         deallocate(dats%hs,dats%ps,stat=iret)
759         nullify(dats%hs,dats%ps)
760         if(iret.ne.0) iret=-3
761     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
762       end subroutine
763     !-------------------------------------------------------------------------------
764       subroutine sigio_aldatm(head,k1,k2,datm,iret)
765         implicit none
766         type(sigio_head),intent(in):: head
767         integer(sigio_intkind),intent(in):: k1,k2
768         type(sigio_datm),intent(inout):: datm
769         integer(sigio_intkind),intent(out):: iret
770         integer nc,dim1,dim3q
771     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
772         call sigio_axdatm(datm,iret)
773         iret=-3
774         if(k1.lt.1.or.k2.gt.head%levs) return
775         nc=(head%jcap+1)*(head%jcap+2)
776         dim1=nc
777         dim3q=head%ntrac
778         datm%k1=k1
779         datm%k2=k2
780         allocate(datm%t(dim1,k1:k2),datm%d(dim1,k1:k2),datm%z(dim1,k1:k2),&
781                  datm%q(dim1,k1:k2,dim3q),stat=iret)
782         if(iret.ne.0) iret=-3
783     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
784       end subroutine
785     !-------------------------------------------------------------------------------
786       subroutine sigio_axdatm(datm,iret)
787         implicit none
788         type(sigio_datm),intent(inout):: datm
789         integer(sigio_intkind),intent(out):: iret
790     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
791         datm%k1=0
792         datm%k2=0
793         deallocate(datm%t,datm%d,datm%z,datm%q,stat=iret)
794         nullify(datm%t,datm%d,datm%z,datm%q)
795         if(iret.ne.0) iret=-3
796     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
797       end subroutine
798     !-------------------------------------------------------------------------------
799       subroutine sigio_aldati(head,i,dati,iret)
800         implicit none
801         type(sigio_head),intent(in):: head
802         integer(sigio_intkind),intent(in):: i
803         type(sigio_dati),intent(inout):: dati
804         integer(sigio_intkind),intent(out):: iret
805         integer dim1
806     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
807         call sigio_axdati(dati,iret)
808         iret=-3
809         if(i.lt.1.or.i.gt.head%ndata) return
810         dim1=head%ldata(i)/(4*head%irealf)
811         dati%i=i
812         allocate(dati%f(dim1),stat=iret)
813         if(iret.ne.0) iret=-3
814     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
815       end subroutine
816     !-------------------------------------------------------------------------------
817       subroutine sigio_axdati(dati,iret)
818         implicit none
819         type(sigio_dati),intent(inout):: dati
820         integer(sigio_intkind),intent(out):: iret
821     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
822         dati%i=0
823         deallocate(dati%f,stat=iret)
824         nullify(dati%f)
825         if(iret.ne.0) iret=-3
826     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
827       end subroutine
828     !-------------------------------------------------------------------------------
829       subroutine sigio_rrdata(lu,head,data,iret)
830         implicit none
831         integer(sigio_intkind),intent(in):: lu
832         type(sigio_head),intent(in):: head
833         type(sigio_data),intent(inout):: data
834         integer(sigio_intkind),intent(out):: iret
835         integer:: i,k,n
836         integer:: nc,mdim1,mdim2,mdim3q
837         integer:: iskip,iread,nread
838         type(sigio_dbta):: dbta
839     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
840         mdim1=min(size(data%hs,1),size(data%ps,1),&
841                   size(data%t,1),size(data%d,1),size(data%z,1),&
842                   size(data%q,1))
843         mdim2=min(size(data%t,2),size(data%d,2),size(data%z,2),&
844                   size(data%q,2))
845         mdim3q=size(data%q,3)
846         nc=(head%jcap+1)*(head%jcap+2)
847         iret=-5
848         if(mdim1.lt.nc.or.&
849            mdim2.lt.head%levs.or.&
850            mdim3q.lt.head%ntrac) return
851         iret=-4
852     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
853         if(head%irealf.ne.2) then
854           iskip=0
855           do i=1,head%nhead
856             call bafrindex(0,iskip,head%lhead(i),iskip)
857           enddo
858           i=1
859           iread=head%ldata(i)
860           call bafrread(lu,iskip,iread,nread,data%hs)
861           if(nread.lt.iread) return
862           i=i+1
863           iskip=iskip+nread
864           iread=head%ldata(i)
865           call bafrread(lu,iskip,iread,nread,data%ps)
866           if(nread.lt.iread) return
867           do k=1,head%levs
868             i=i+1
869             iskip=iskip+nread
870             iread=head%ldata(i)
871             call bafrread(lu,iskip,iread,nread,data%t(1,k))
872             if(nread.lt.iread) return
873           enddo
874           do k=1,head%levs
875             i=i+1
876             iskip=iskip+nread
877             iread=head%ldata(i)
878             call bafrread(lu,iskip,iread,nread,data%d(1,k))
879             if(nread.lt.iread) return
880             i=i+1
881             iskip=iskip+nread
882             iread=head%ldata(i)
883             call bafrread(lu,iskip,iread,nread,data%z(1,k))
884             if(nread.lt.iread) return
885           enddo
886           do n=1,head%ntrac
887             do k=1,head%levs
888               i=i+1
889               iskip=iskip+nread
890               iread=head%ldata(i)
891               call bafrread(lu,iskip,iread,nread,data%q(1,k,n))
892               if(nread.lt.iread) return
893             enddo
894           enddo
895           do n=1,head%nxgr
896             i=i+1
897             iskip=iskip+nread
898             iread=head%ldata(i)
899             call bafrread(lu,iskip,iread,nread,data%xgr(1,1,n))
900             if(nread.lt.iread) return
901           enddo
902           if(head%nxss.gt.0) then
903             i=i+1
904             iskip=iskip+nread
905             iread=head%ldata(i)
906             call bafrread(lu,iskip,iread,nread,data%xss)
907             if(nread.lt.iread) return
908           endif
909         else
910           call sigio_aldbta(head,dbta,iret)
911           if(iret.ne.0) return
912           call sigio_rrdbta(lu,head,dbta,iret)
913           if(iret.ne.0) return
914           data%hs(:nc)=dbta%hs(:nc)
915           data%ps(:nc)=dbta%ps(:nc)
916           data%t(:nc,:head%levs)=dbta%t(:nc,:head%levs)
917           data%d(:nc,:head%levs)=dbta%d(:nc,:head%levs)
918           data%z(:nc,:head%levs)=dbta%z(:nc,:head%levs)
919           data%q(:nc,:head%levs,:head%ntrac)=dbta%q(:nc,:head%levs,:head%ntrac)
920           data%xgr(:head%lonb,:head%latb,:head%nxgr)=&
921            dbta%xgr(:head%lonb,:head%latb,:head%nxgr)
922           data%xss(:head%nxss)=dbta%xss(:head%nxss)
923           call sigio_axdbta(dbta,iret)
924         endif
925         iret=0
926     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
927       end subroutine
928     !-------------------------------------------------------------------------------
929       subroutine sigio_rwdata(lu,head,data,iret)
930         implicit none
931         integer(sigio_intkind),intent(in):: lu
932         type(sigio_head),intent(in):: head
933         type(sigio_data),intent(in):: data
934         integer(sigio_intkind),intent(out):: iret
935         integer:: i,k,n
936         integer:: nc,mdim1,mdim2,mdim3q
937         integer:: iskip,iwrite,nwrite
938         type(sigio_dbta):: dbta
939     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
940         mdim1=min(size(data%hs,1),size(data%ps,1),&
941                   size(data%t,1),size(data%d,1),size(data%z,1),&
942                   size(data%q,1))
943         mdim2=min(size(data%t,2),size(data%d,2),size(data%z,2),&
944                   size(data%q,2))
945         mdim3q=size(data%q,3)
946         nc=(head%jcap+1)*(head%jcap+2)
947         iret=-5
948         if(mdim1.lt.nc.or.&
949            mdim2.lt.head%levs.or.&
950            mdim3q.lt.head%ntrac) return
951         iret=-4
952     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
953         if(head%irealf.ne.2) then
954           iskip=0
955           do i=1,head%nhead
956             call bafrindex(0,iskip,head%lhead(i),iskip)
957           enddo
958           i=1
959           iwrite=head%ldata(i)
960           call bafrwrite(lu,iskip,iwrite,nwrite,data%hs)
961           if(nwrite.lt.iwrite) return
962           i=i+1
963           iskip=iskip+nwrite
964           iwrite=head%ldata(i)
965           call bafrwrite(lu,iskip,iwrite,nwrite,data%ps)
966           if(nwrite.lt.iwrite) return
967           do k=1,head%levs
968             i=i+1
969             iskip=iskip+nwrite
970             iwrite=head%ldata(i)
971             call bafrwrite(lu,iskip,iwrite,nwrite,data%t(1,k))
972             if(nwrite.lt.iwrite) return
973           enddo
974           do k=1,head%levs
975             i=i+1
976             iskip=iskip+nwrite
977             iwrite=head%ldata(i)
978             call bafrwrite(lu,iskip,iwrite,nwrite,data%d(1,k))
979             if(nwrite.lt.iwrite) return
980             i=i+1
981             iskip=iskip+nwrite
982             iwrite=head%ldata(i)
983             call bafrwrite(lu,iskip,iwrite,nwrite,data%z(1,k))
984             if(nwrite.lt.iwrite) return
985           enddo
986           do n=1,head%ntrac
987             do k=1,head%levs
988               i=i+1
989               iskip=iskip+nwrite
990               iwrite=head%ldata(i)
991               call bafrwrite(lu,iskip,iwrite,nwrite,data%q(1,k,n))
992               if(nwrite.lt.iwrite) return
993             enddo
994           enddo
995           do n=1,head%nxgr
996             i=i+1
997             iskip=iskip+nwrite
998             iwrite=head%ldata(i)
999             call bafrwrite(lu,iskip,iwrite,nwrite,data%xgr(1,1,n))
1000             if(nwrite.lt.iwrite) return
1001           enddo
1002           if(head%nxss.gt.0) then
1003             i=i+1
1004             iskip=iskip+nwrite
1005             iwrite=head%ldata(i)
1006             call bafrwrite(lu,iskip,iwrite,nwrite,data%xss)
1007             if(nwrite.lt.iwrite) return
1008           endif
1009         else
1010           call sigio_aldbta(head,dbta,iret)
1011           if(iret.ne.0) return
1012           dbta%hs(:nc)=data%hs(:nc)
1013           dbta%ps(:nc)=data%ps(:nc)
1014           dbta%t(:nc,:head%levs)=data%t(:nc,:head%levs)
1015           dbta%d(:nc,:head%levs)=data%d(:nc,:head%levs)
1016           dbta%z(:nc,:head%levs)=data%z(:nc,:head%levs)
1017           dbta%q(:nc,:head%levs,:head%ntrac)=data%q(:nc,:head%levs,:head%ntrac)
1018           dbta%xgr(:head%lonb,:head%latb,:head%nxgr)=&
1019            data%xgr(:head%lonb,:head%latb,:head%nxgr)
1020           dbta%xss(:head%nxss)=data%xss(:head%nxss)
1021           call sigio_rwdbta(lu,head,dbta,iret)
1022           if(iret.ne.0) return
1023           call sigio_axdbta(dbta,iret)
1024         endif
1025         iret=0
1026     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1027       end subroutine
1028     !-------------------------------------------------------------------------------
1029       subroutine sigio_rrohdca(lu,cfname,head,data,iret)
1030         implicit none
1031         integer(sigio_intkind),intent(in):: lu
1032         character*(*),intent(in):: cfname
1033         type(sigio_head),intent(inout):: head
1034         type(sigio_data),intent(inout):: data
1035         integer(sigio_intkind),intent(out):: iret
1036     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1037         call sigio_rropen(lu,cfname,iret)
1038         if(iret.ne.0) return
1039     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1040         call sigio_rrhead(lu,head,iret)
1041         if(iret.ne.0) return
1042     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1043         call sigio_aldata(head,data,iret)
1044         if(iret.ne.0) return
1045     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1046         call sigio_rrdata(lu,head,data,iret)
1047         if(iret.ne.0) return
1048     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1049         call sigio_rclose(lu,iret)
1050         if(iret.ne.0) return
1051     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1052       end subroutine
1053     !-------------------------------------------------------------------------------
1054       subroutine sigio_rwohdca(lu,cfname,head,data,iret)
1055         implicit none
1056         integer(sigio_intkind),intent(in):: lu
1057         character*(*),intent(in):: cfname
1058         type(sigio_head),intent(inout):: head
1059         type(sigio_data),intent(in):: data
1060         integer(sigio_intkind),intent(out):: iret
1061     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1062         call sigio_rwopen(lu,cfname,iret)
1063         if(iret.ne.0) return
1064     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1065         call sigio_rwhead(lu,head,iret)
1066         if(iret.ne.0) return
1067     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1068         call sigio_rwdata(lu,head,data,iret)
1069         if(iret.ne.0) return
1070     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1071         call sigio_rclose(lu,iret)
1072         if(iret.ne.0) return
1073     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1074       end subroutine
1075     !-------------------------------------------------------------------------------
1076       subroutine sigio_rrdats(lu,head,dats,iret)
1077         implicit none
1078         integer(sigio_intkind),intent(in):: lu
1079         type(sigio_head),intent(in):: head
1080         type(sigio_dats),intent(inout):: dats
1081         integer(sigio_intkind),intent(out):: iret
1082         integer:: i
1083         integer:: nc,mdim1
1084         integer:: iskip,iread,nread
1085         type(sigio_dbts):: dbts
1086     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1087         mdim1=min(size(dats%hs,1),size(dats%ps,1))
1088         nc=(head%jcap+1)*(head%jcap+2)
1089         iret=-5
1090         if(mdim1.lt.nc) return
1091         iret=-4
1092     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1093         if(head%irealf.ne.2) then
1094           iskip=0
1095           do i=1,head%nhead
1096             call bafrindex(0,iskip,head%lhead(i),iskip)
1097           enddo
1098           i=1
1099           iread=head%ldata(i)
1100           call bafrread(lu,iskip,iread,nread,dats%hs)
1101           if(nread.lt.iread) return
1102           i=i+1
1103           iskip=iskip+nread
1104           iread=head%ldata(i)
1105           call bafrread(lu,iskip,iread,nread,dats%ps)
1106           if(nread.lt.iread) return
1107         else
1108           call sigio_aldbts(head,dbts,iret)
1109           if(iret.ne.0) return
1110           call sigio_rrdbts(lu,head,dbts,iret)
1111           if(iret.ne.0) return
1112           dats%hs(:nc)=dbts%hs(:nc)
1113           dats%ps(:nc)=dbts%ps(:nc)
1114           call sigio_axdbts(dbts,iret)
1115         endif
1116         iret=0
1117     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1118       end subroutine
1119     !-------------------------------------------------------------------------------
1120       subroutine sigio_rwdats(lu,head,dats,iret)
1121         implicit none
1122         integer(sigio_intkind),intent(in):: lu
1123         type(sigio_head),intent(in):: head
1124         type(sigio_dats),intent(in):: dats
1125         integer(sigio_intkind),intent(out):: iret
1126         integer:: i
1127         integer:: nc,mdim1
1128         integer:: iskip,iwrite,nwrite
1129         type(sigio_dbts):: dbts
1130     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1131         mdim1=min(size(dats%hs,1),size(dats%ps,1))
1132         nc=(head%jcap+1)*(head%jcap+2)
1133         iret=-5
1134         if(mdim1.lt.nc) return
1135         iret=-4
1136     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1137         if(head%irealf.ne.2) then
1138           iskip=0
1139           do i=1,head%nhead
1140             call bafrindex(0,iskip,head%lhead(i),iskip)
1141           enddo
1142           i=1
1143           iwrite=head%ldata(i)
1144           call bafrwrite(lu,iskip,iwrite,nwrite,dats%hs)
1145           if(nwrite.lt.iwrite) return
1146           i=i+1
1147           iskip=iskip+nwrite
1148           iwrite=head%ldata(i)
1149           call bafrwrite(lu,iskip,iwrite,nwrite,dats%ps)
1150           if(nwrite.lt.iwrite) return
1151         else
1152           call sigio_aldbts(head,dbts,iret)
1153           if(iret.ne.0) return
1154           dbts%hs(:nc)=dats%hs(:nc)
1155           dbts%ps(:nc)=dats%ps(:nc)
1156           call sigio_rwdbts(lu,head,dbts,iret)
1157           if(iret.ne.0) return
1158           call sigio_axdbts(dbts,iret)
1159         endif
1160         iret=0
1161     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1162       end subroutine
1163     !-------------------------------------------------------------------------------
1164       subroutine sigio_rrdatm(lu,head,datm,iret)
1165         implicit none
1166         integer(sigio_intkind),intent(in):: lu
1167         type(sigio_head),intent(in):: head
1168         type(sigio_datm),intent(inout):: datm
1169         integer(sigio_intkind),intent(out):: iret
1170         integer:: i,k,n
1171         integer:: nc,k1,k2,mdim1,ldim2,udim2,mdim3q
1172         integer:: iskip,iread,nread
1173         type(sigio_dbtm):: dbtm
1174     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1175         k1=datm%k1
1176         k2=datm%k2
1177         mdim1=min(size(datm%t,1),size(datm%d,1),size(datm%z,1),&
1178                   size(datm%q,1))
1179         ldim2=max(lbound(datm%t,2),lbound(datm%d,2),lbound(datm%z,2),&
1180                   lbound(datm%q,2))
1181         udim2=min(ubound(datm%t,2),ubound(datm%d,2),ubound(datm%z,2),&
1182                   ubound(datm%q,2))
1183         mdim3q=size(datm%q,3)
1184         nc=(head%jcap+1)*(head%jcap+2)
1185         iret=-5
1186         if(k1.lt.1.or.k2.gt.head%levs.or.&
1187            mdim1.lt.nc.or.&
1188            ldim2.gt.k1.or.udim2.lt.k2.or.&
1189            mdim3q.lt.head%ntrac) return
1190         iret=-4
1191     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1192         if(head%irealf.ne.2) then
1193           iskip=0
1194           do i=1,head%nhead
1195             call bafrindex(0,iskip,head%lhead(i),iskip)
1196           enddo
1197           i=1
1198           call bafrindex(0,iskip,head%ldata(i),iskip)
1199           i=i+1
1200           call bafrindex(0,iskip,head%ldata(i),iskip)
1201           do k=1,head%levs
1202             if(k.lt.k1.or.k.gt.k2) then
1203               i=i+1
1204               call bafrindex(0,iskip,head%ldata(i),iskip)
1205             else
1206               i=i+1
1207               iread=head%ldata(i)
1208               call bafrread(lu,iskip,iread,nread,datm%t(1,k))
1209               if(nread.lt.iread) return
1210               iskip=iskip+nread
1211             endif
1212           enddo
1213           do k=1,head%levs
1214             if(k.lt.k1.or.k.gt.k2) then
1215               i=i+1
1216               call bafrindex(0,iskip,head%ldata(i),iskip)
1217               i=i+1
1218               call bafrindex(0,iskip,head%ldata(i),iskip)
1219             else
1220               iread=head%ldata(i)
1221               call bafrread(lu,iskip,iread,nread,datm%d(1,k))
1222               if(nread.lt.iread) return
1223               iskip=iskip+nread
1224               iread=head%ldata(i)
1225               call bafrread(lu,iskip,iread,nread,datm%z(1,k))
1226               if(nread.lt.iread) return
1227               iskip=iskip+nread
1228             endif
1229           enddo
1230           do n=1,head%ntrac
1231             do k=1,head%levs
1232               if(k.lt.k1.or.k.gt.k2) then
1233                 i=i+1
1234                 call bafrindex(0,iskip,head%ldata(i),iskip)
1235               else
1236                 i=i+1
1237                 iread=head%ldata(i)
1238                 call bafrread(lu,iskip,iread,nread,datm%q(1,k,n))
1239                 if(nread.lt.iread) return
1240                 iskip=iskip+nread
1241               endif
1242             enddo
1243           enddo
1244         else
1245           call sigio_aldbtm(head,k1,k2,dbtm,iret)
1246           if(iret.ne.0) return
1247           call sigio_rrdbtm(lu,head,dbtm,iret)
1248           if(iret.ne.0) return
1249           datm%t(:nc,k1:k2)=dbtm%t(:nc,k1:k2)
1250           datm%d(:nc,k1:k2)=dbtm%d(:nc,k1:k2)
1251           datm%z(:nc,k1:k2)=dbtm%z(:nc,k1:k2)
1252           datm%q(:nc,k1:k2,:head%ntrac)=dbtm%q(:nc,k1:k2,:head%ntrac)
1253           call sigio_axdbtm(dbtm,iret)
1254         endif
1255         iret=0
1256     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1257       end subroutine
1258     !-------------------------------------------------------------------------------
1259       subroutine sigio_rwdatm(lu,head,datm,iret)
1260         implicit none
1261         integer(sigio_intkind),intent(in):: lu
1262         type(sigio_head),intent(in):: head
1263         type(sigio_datm),intent(in):: datm
1264         integer(sigio_intkind),intent(out):: iret
1265         integer:: i,k,n
1266         integer:: nc,k1,k2,mdim1,ldim2,udim2,mdim3q
1267         integer:: iskip,iwrite,nwrite
1268         type(sigio_dbtm):: dbtm
1269     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1270         k1=datm%k1
1271         k2=datm%k2
1272         mdim1=min(size(datm%t,1),size(datm%d,1),size(datm%z,1),&
1273                   size(datm%q,1))
1274         ldim2=max(lbound(datm%t,2),lbound(datm%d,2),lbound(datm%z,2),&
1275                   lbound(datm%q,2))
1276         udim2=min(ubound(datm%t,2),ubound(datm%d,2),ubound(datm%z,2),&
1277                   ubound(datm%q,2))
1278         mdim3q=size(datm%q,3)
1279         nc=(head%jcap+1)*(head%jcap+2)
1280         iret=-5
1281         if(k1.lt.1.or.k2.gt.head%levs.or.&
1282            mdim1.lt.nc.or.&
1283            ldim2.gt.k1.or.udim2.lt.k2.or.&
1284            mdim3q.lt.head%ntrac) return
1285         iret=-4
1286     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1287         if(head%irealf.ne.2) then
1288           iskip=0
1289           do i=1,head%nhead
1290             call bafrindex(0,iskip,head%lhead(i),iskip)
1291           enddo
1292           i=1
1293           call bafrindex(0,iskip,head%ldata(i),iskip)
1294           i=i+1
1295           call bafrindex(0,iskip,head%ldata(i),iskip)
1296           do k=1,head%levs
1297             if(k.lt.k1.or.k.gt.k2) then
1298               i=i+1
1299               call bafrindex(0,iskip,head%ldata(i),iskip)
1300             else
1301               i=i+1
1302               iwrite=head%ldata(i)
1303               call bafrwrite(lu,iskip,iwrite,nwrite,datm%t(1,k))
1304               if(nwrite.lt.iwrite) return
1305               iskip=iskip+nwrite
1306             endif
1307           enddo
1308           do k=1,head%levs
1309             if(k.lt.k1.or.k.gt.k2) then
1310               i=i+1
1311               call bafrindex(0,iskip,head%ldata(i),iskip)
1312               i=i+1
1313               call bafrindex(0,iskip,head%ldata(i),iskip)
1314             else
1315               iwrite=head%ldata(i)
1316               call bafrwrite(lu,iskip,iwrite,nwrite,datm%d(1,k))
1317               if(nwrite.lt.iwrite) return
1318               iskip=iskip+nwrite
1319               iwrite=head%ldata(i)
1320               call bafrwrite(lu,iskip,iwrite,nwrite,datm%z(1,k))
1321               if(nwrite.lt.iwrite) return
1322               iskip=iskip+nwrite
1323             endif
1324           enddo
1325           do n=1,head%ntrac
1326             do k=1,head%levs
1327               if(k.lt.k1.or.k.gt.k2) then
1328                 i=i+1
1329                 call bafrindex(0,iskip,head%ldata(i),iskip)
1330               else
1331                 i=i+1
1332                 iwrite=head%ldata(i)
1333                 call bafrwrite(lu,iskip,iwrite,nwrite,datm%q(1,k,n))
1334                 if(nwrite.lt.iwrite) return
1335                 iskip=iskip+nwrite
1336               endif
1337             enddo
1338           enddo
1339         else
1340           call sigio_aldbtm(head,k1,k2,dbtm,iret)
1341           if(iret.ne.0) return
1342           dbtm%t(:nc,k1:k2)=datm%t(:nc,k1:k2)
1343           dbtm%d(:nc,k1:k2)=datm%d(:nc,k1:k2)
1344           dbtm%z(:nc,k1:k2)=datm%z(:nc,k1:k2)
1345           dbtm%q(:nc,k1:k2,:head%ntrac)=datm%q(:nc,k1:k2,:head%ntrac)
1346           call sigio_rwdbtm(lu,head,dbtm,iret)
1347           if(iret.ne.0) return
1348           call sigio_axdbtm(dbtm,iret)
1349         endif
1350         iret=0
1351     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1352       end subroutine
1353     !-------------------------------------------------------------------------------
1354       subroutine sigio_rrdati(lu,head,dati,iret)
1355         implicit none
1356         integer(sigio_intkind),intent(in):: lu
1357         type(sigio_head),intent(in):: head
1358         type(sigio_dati),intent(inout):: dati
1359         integer(sigio_intkind),intent(out):: iret
1360         integer:: i,k,n
1361         integer:: mdim1
1362         integer:: mlen
1363         integer:: iskip,iread,nread
1364         type(sigio_dbti):: dbti
1365     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1366         i=dati%i
1367         mdim1=size(dati%f,1)
1368         iret=-5
1369         if(i.lt.1.or.i.gt.head%ndata) return
1370         mlen=head%ldata(i)/(4*head%irealf)
1371         if(mdim1.lt.mlen) return
1372         iret=-4
1373     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1374         if(head%irealf.ne.2) then
1375           iskip=0
1376           do i=1,head%nhead
1377             call bafrindex(0,iskip,head%lhead(i),iskip)
1378           enddo
1379           do i=1,dati%i-1
1380             call bafrindex(0,iskip,head%ldata(i),iskip)
1381           enddo
1382           i=dati%i
1383           iread=head%ldata(i)
1384           call bafrread(lu,iskip,iread,nread,dati%f)
1385           if(nread.lt.iread) return
1386         else
1387           i=dati%i
1388           call sigio_aldbti(head,i,dbti,iret)
1389           if(iret.ne.0) return
1390           call sigio_rrdbti(lu,head,dbti,iret)
1391           if(iret.ne.0) return
1392           dati%f(:mlen)=dbti%f(:mlen)
1393           call sigio_axdbti(dbti,iret)
1394         endif
1395         iret=0
1396     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1397       end subroutine
1398     !-------------------------------------------------------------------------------
1399       subroutine sigio_rwdati(lu,head,dati,iret)
1400         implicit none
1401         integer(sigio_intkind),intent(in):: lu
1402         type(sigio_head),intent(in):: head
1403         type(sigio_dati),intent(in):: dati
1404         integer(sigio_intkind),intent(out):: iret
1405         integer:: i,k,n
1406         integer:: mdim1
1407         integer:: mlen
1408         integer:: iskip,iwrite,nwrite
1409         type(sigio_dbti):: dbti
1410     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1411         i=dati%i
1412         mdim1=size(dati%f,1)
1413         iret=-5
1414         if(i.lt.1.or.i.gt.head%ndata) return
1415         mlen=head%ldata(i)/(4*head%irealf)
1416         if(mdim1.lt.mlen) return
1417         iret=-4
1418     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1419         if(head%irealf.ne.2) then
1420           iskip=0
1421           do i=1,head%nhead
1422             call bafrindex(0,iskip,head%lhead(i),iskip)
1423           enddo
1424           do i=1,dati%i-1
1425             call bafrindex(0,iskip,head%ldata(i),iskip)
1426           enddo
1427           i=dati%i
1428           iwrite=head%ldata(i)
1429           call bafrwrite(lu,iskip,iwrite,nwrite,dati%f)
1430           if(nwrite.lt.iwrite) return
1431           iret=0
1432         else
1433           i=dati%i
1434           call sigio_aldbti(head,i,dbti,iret)
1435           if(iret.ne.0) return
1436           dbti%f(:mlen)=dati%f(:mlen)
1437           call sigio_rwdbti(lu,head,dbti,iret)
1438           if(iret.ne.0) return
1439           call sigio_axdbti(dbti,iret)
1440         endif
1441     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1442       end subroutine
1443     !-------------------------------------------------------------------------------
1444       subroutine sigio_aldbts(head,dbts,iret)
1445         implicit none
1446         type(sigio_head),intent(in):: head
1447         type(sigio_dbts),intent(inout):: dbts
1448         integer(sigio_intkind),intent(out):: iret
1449         integer nc,dim1
1450     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1451         call sigio_axdbts(dbts,iret)
1452         nc=(head%jcap+1)*(head%jcap+2)
1453         dim1=nc
1454         allocate(dbts%hs(dim1),dbts%ps(dim1),stat=iret)
1455         if(iret.ne.0) iret=-3
1456     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1457       end subroutine
1458     !-------------------------------------------------------------------------------
1459       subroutine sigio_axdbts(dbts,iret)
1460         implicit none
1461         type(sigio_dbts),intent(inout):: dbts
1462         integer(sigio_intkind),intent(out):: iret
1463     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1464         deallocate(dbts%hs,dbts%ps,stat=iret)
1465         nullify(dbts%hs,dbts%ps)
1466         if(iret.ne.0) iret=-3
1467     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1468       end subroutine
1469     !-------------------------------------------------------------------------------
1470       subroutine sigio_aldbtm(head,k1,k2,dbtm,iret)
1471         implicit none
1472         type(sigio_head),intent(in):: head
1473         integer(sigio_intkind),intent(in):: k1,k2
1474         type(sigio_dbtm),intent(inout):: dbtm
1475         integer(sigio_intkind),intent(out):: iret
1476         integer nc,dim1,dim3q
1477     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1478         call sigio_axdbtm(dbtm,iret)
1479         iret=-3
1480         if(k1.lt.1.or.k2.gt.head%levs) return
1481         nc=(head%jcap+1)*(head%jcap+2)
1482         dim1=nc
1483         dim3q=head%ntrac
1484         dbtm%k1=k1
1485         dbtm%k2=k2
1486         allocate(dbtm%t(dim1,k1:k2),dbtm%d(dim1,k1:k2),dbtm%z(dim1,k1:k2),&
1487                  dbtm%q(dim1,k1:k2,dim3q),stat=iret)
1488         if(iret.ne.0) iret=-3
1489     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1490       end subroutine
1491     !-------------------------------------------------------------------------------
1492       subroutine sigio_axdbtm(dbtm,iret)
1493         implicit none
1494         type(sigio_dbtm),intent(inout):: dbtm
1495         integer(sigio_intkind),intent(out):: iret
1496     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1497         dbtm%k1=0
1498         dbtm%k2=0
1499         deallocate(dbtm%t,dbtm%d,dbtm%z,dbtm%q,stat=iret)
1500         nullify(dbtm%t,dbtm%d,dbtm%z,dbtm%q)
1501         if(iret.ne.0) iret=-3
1502     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1503       end subroutine
1504     !-------------------------------------------------------------------------------
1505       subroutine sigio_aldbti(head,i,dbti,iret)
1506         implicit none
1507         type(sigio_head),intent(in):: head
1508         integer(sigio_intkind),intent(in):: i
1509         type(sigio_dbti),intent(inout):: dbti
1510         integer(sigio_intkind),intent(out):: iret
1511         integer dim1
1512     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1513         call sigio_axdbti(dbti,iret)
1514         iret=-3
1515         if(i.lt.1.or.i.gt.head%ndata) return
1516         dim1=head%ldata(i)/(4*head%irealf)
1517         dbti%i=i
1518         allocate(dbti%f(dim1),stat=iret)
1519         if(iret.ne.0) iret=-3
1520     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1521       end subroutine
1522     !-------------------------------------------------------------------------------
1523       subroutine sigio_axdbti(dbti,iret)
1524         implicit none
1525         type(sigio_dbti),intent(inout):: dbti
1526         integer(sigio_intkind),intent(out):: iret
1527     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1528         dbti%i=0
1529         deallocate(dbti%f,stat=iret)
1530         nullify(dbti%f)
1531         if(iret.ne.0) iret=-3
1532     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1533       end subroutine
1534     !-------------------------------------------------------------------------------
1535       subroutine sigio_rrdbta(lu,head,dbta,iret)
1536         implicit none
1537         integer(sigio_intkind),intent(in):: lu
1538         type(sigio_head),intent(in):: head
1539         type(sigio_dbta),intent(inout):: dbta
1540         integer(sigio_intkind),intent(out):: iret
1541         integer:: i,k,n
1542         integer:: nc,mdim1,mdim2,mdim3q
1543         integer:: iskip,iread,nread
1544         type(sigio_data):: data
1545     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1546         mdim1=min(size(dbta%hs,1),size(dbta%ps,1),&
1547                   size(dbta%t,1),size(dbta%d,1),size(dbta%z,1),&
1548                   size(dbta%q,1))
1549         mdim2=min(size(dbta%t,2),size(dbta%d,2),size(dbta%z,2),&
1550                   size(dbta%q,2))
1551         mdim3q=size(dbta%q,3)
1552         nc=(head%jcap+1)*(head%jcap+2)
1553         iret=-5
1554         if(mdim1.lt.nc.or.&
1555            mdim2.lt.head%levs.or.&
1556            mdim3q.lt.head%ntrac) return
1557         iret=-4
1558     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1559         if(head%irealf.eq.2) then
1560           iskip=0
1561           do i=1,head%nhead
1562             call bafrindex(0,iskip,head%lhead(i),iskip)
1563           enddo
1564           i=1
1565           iread=head%ldata(i)
1566           call bafrread(lu,iskip,iread,nread,dbta%hs)
1567           if(nread.lt.iread) return
1568           i=i+1
1569           iskip=iskip+nread
1570           iread=head%ldata(i)
1571           call bafrread(lu,iskip,iread,nread,dbta%ps)
1572           if(nread.lt.iread) return
1573           do k=1,head%levs
1574             i=i+1
1575             iskip=iskip+nread
1576             iread=head%ldata(i)
1577             call bafrread(lu,iskip,iread,nread,dbta%t(1,k))
1578             if(nread.lt.iread) return
1579           enddo
1580           do k=1,head%levs
1581             i=i+1
1582             iskip=iskip+nread
1583             iread=head%ldata(i)
1584             call bafrread(lu,iskip,iread,nread,dbta%d(1,k))
1585             if(nread.lt.iread) return
1586             i=i+1
1587             iskip=iskip+nread
1588             iread=head%ldata(i)
1589             call bafrread(lu,iskip,iread,nread,dbta%z(1,k))
1590             if(nread.lt.iread) return
1591           enddo
1592           do n=1,head%ntrac
1593             do k=1,head%levs
1594               i=i+1
1595               iskip=iskip+nread
1596               iread=head%ldata(i)
1597               call bafrread(lu,iskip,iread,nread,dbta%q(1,k,n))
1598               if(nread.lt.iread) return
1599             enddo
1600           enddo
1601           do n=1,head%nxgr
1602             i=i+1
1603             iskip=iskip+nread
1604             iread=head%ldata(i)
1605             call bafrread(lu,iskip,iread,nread,dbta%xgr(1,1,n))
1606             if(nread.lt.iread) return
1607           enddo
1608           if(head%nxss.gt.0) then
1609             i=i+1
1610             iskip=iskip+nread
1611             iread=head%ldata(i)
1612             call bafrread(lu,iskip,iread,nread,dbta%xss)
1613             if(nread.lt.iread) return
1614           endif
1615         else
1616           call sigio_aldata(head,data,iret)
1617           if(iret.ne.0) return
1618           call sigio_rrdata(lu,head,data,iret)
1619           if(iret.ne.0) return
1620           dbta%hs(:nc)=data%hs(:nc)
1621           dbta%ps(:nc)=data%ps(:nc)
1622           dbta%t(:nc,:head%levs)=data%t(:nc,:head%levs)
1623           dbta%d(:nc,:head%levs)=data%d(:nc,:head%levs)
1624           dbta%z(:nc,:head%levs)=data%z(:nc,:head%levs)
1625           dbta%q(:nc,:head%levs,:head%ntrac)=data%q(:nc,:head%levs,:head%ntrac)
1626           dbta%xgr(:head%lonb,:head%latb,:head%nxgr)=&
1627            data%xgr(:head%lonb,:head%latb,:head%nxgr)
1628           dbta%xss(:head%nxss)=data%xss(:head%nxss)
1629           call sigio_axdata(data,iret)
1630         endif
1631         iret=0
1632     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1633       end subroutine
1634     !-------------------------------------------------------------------------------
1635       subroutine sigio_rwdbta(lu,head,dbta,iret)
1636         implicit none
1637         integer(sigio_intkind),intent(in):: lu
1638         type(sigio_head),intent(in):: head
1639         type(sigio_dbta),intent(in):: dbta
1640         integer(sigio_intkind),intent(out):: iret
1641         integer:: i,k,n
1642         integer:: nc,mdim1,mdim2,mdim3q
1643         integer:: iskip,iwrite,nwrite
1644         type(sigio_data):: data
1645     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1646         mdim1=min(size(dbta%hs,1),size(dbta%ps,1),&
1647                   size(dbta%t,1),size(dbta%d,1),size(dbta%z,1),&
1648                   size(dbta%q,1))
1649         mdim2=min(size(dbta%t,2),size(dbta%d,2),size(dbta%z,2),&
1650                   size(dbta%q,2))
1651         mdim3q=size(dbta%q,3)
1652         nc=(head%jcap+1)*(head%jcap+2)
1653         iret=-5
1654         if(mdim1.lt.nc.or.&
1655            mdim2.lt.head%levs.or.&
1656            mdim3q.lt.head%ntrac) return
1657         iret=-4
1658     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1659         if(head%irealf.eq.2) then
1660           iskip=0
1661           do i=1,head%nhead
1662             call bafrindex(0,iskip,head%lhead(i),iskip)
1663           enddo
1664           i=1
1665           iwrite=head%ldata(i)
1666           call bafrwrite(lu,iskip,iwrite,nwrite,dbta%hs)
1667           if(nwrite.lt.iwrite) return
1668           i=i+1
1669           iskip=iskip+nwrite
1670           iwrite=head%ldata(i)
1671           call bafrwrite(lu,iskip,iwrite,nwrite,dbta%ps)
1672           if(nwrite.lt.iwrite) return
1673           do k=1,head%levs
1674             i=i+1
1675             iskip=iskip+nwrite
1676             iwrite=head%ldata(i)
1677             call bafrwrite(lu,iskip,iwrite,nwrite,dbta%t(1,k))
1678             if(nwrite.lt.iwrite) return
1679           enddo
1680           do k=1,head%levs
1681             i=i+1
1682             iskip=iskip+nwrite
1683             iwrite=head%ldata(i)
1684             call bafrwrite(lu,iskip,iwrite,nwrite,dbta%d(1,k))
1685             if(nwrite.lt.iwrite) return
1686             i=i+1
1687             iskip=iskip+nwrite
1688             iwrite=head%ldata(i)
1689             call bafrwrite(lu,iskip,iwrite,nwrite,dbta%z(1,k))
1690             if(nwrite.lt.iwrite) return
1691           enddo
1692           do n=1,head%ntrac
1693             do k=1,head%levs
1694               i=i+1
1695               iskip=iskip+nwrite
1696               iwrite=head%ldata(i)
1697               call bafrwrite(lu,iskip,iwrite,nwrite,dbta%q(1,k,n))
1698               if(nwrite.lt.iwrite) return
1699             enddo
1700           enddo
1701           do n=1,head%nxgr
1702             i=i+1
1703             iskip=iskip+nwrite
1704             iwrite=head%ldata(i)
1705             call bafrwrite(lu,iskip,iwrite,nwrite,dbta%xgr(1,1,n))
1706             if(nwrite.lt.iwrite) return
1707           enddo
1708           if(head%nxss.gt.0) then
1709             i=i+1
1710             iskip=iskip+nwrite
1711             iwrite=head%ldata(i)
1712             call bafrwrite(lu,iskip,iwrite,nwrite,dbta%xss)
1713             if(nwrite.lt.iwrite) return
1714           endif
1715         else
1716           call sigio_aldata(head,data,iret)
1717           if(iret.ne.0) return
1718           data%hs(:nc)=dbta%hs(:nc)
1719           data%ps(:nc)=dbta%ps(:nc)
1720           data%t(:nc,:head%levs)=dbta%t(:nc,:head%levs)
1721           data%d(:nc,:head%levs)=dbta%d(:nc,:head%levs)
1722           data%z(:nc,:head%levs)=dbta%z(:nc,:head%levs)
1723           data%q(:nc,:head%levs,:head%ntrac)=dbta%q(:nc,:head%levs,:head%ntrac)
1724           data%xgr(:head%lonb,:head%latb,:head%nxgr)=&
1725            dbta%xgr(:head%lonb,:head%latb,:head%nxgr)
1726           data%xss(:head%nxss)=dbta%xss(:head%nxss)
1727           call sigio_rwdata(lu,head,data,iret)
1728           if(iret.ne.0) return
1729           call sigio_axdata(data,iret)
1730         endif
1731         iret=0
1732     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1733       end subroutine
1734     !-------------------------------------------------------------------------------
1735       subroutine sigio_rrohdcb(lu,cfname,head,dbta,iret)
1736         implicit none
1737         integer(sigio_intkind),intent(in):: lu
1738         character*(*),intent(in):: cfname
1739         type(sigio_head),intent(inout):: head
1740         type(sigio_dbta),intent(inout):: dbta
1741         integer(sigio_intkind),intent(out):: iret
1742     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1743         call sigio_rropen(lu,cfname,iret)
1744         if(iret.ne.0) return
1745     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1746         call sigio_rrhead(lu,head,iret)
1747         if(iret.ne.0) return
1748     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1749         call sigio_aldbta(head,dbta,iret)
1750         if(iret.ne.0) return
1751     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1752         call sigio_rrdbta(lu,head,dbta,iret)
1753         if(iret.ne.0) return
1754     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1755         call sigio_rclose(lu,iret)
1756         if(iret.ne.0) return
1757     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1758       end subroutine
1759     !-------------------------------------------------------------------------------
1760       subroutine sigio_rwohdcb(lu,cfname,head,dbta,iret)
1761         implicit none
1762         integer(sigio_intkind),intent(in):: lu
1763         character*(*),intent(in):: cfname
1764         type(sigio_head),intent(inout):: head
1765         type(sigio_dbta),intent(in):: dbta
1766         integer(sigio_intkind),intent(out):: iret
1767     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1768         call sigio_rwopen(lu,cfname,iret)
1769         if(iret.ne.0) return
1770     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1771         call sigio_rwhead(lu,head,iret)
1772         if(iret.ne.0) return
1773     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1774         call sigio_rwdbta(lu,head,dbta,iret)
1775         if(iret.ne.0) return
1776     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1777         call sigio_rclose(lu,iret)
1778         if(iret.ne.0) return
1779     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1780       end subroutine
1781     !-------------------------------------------------------------------------------
1782       subroutine sigio_rrdbts(lu,head,dbts,iret)
1783         implicit none
1784         integer(sigio_intkind),intent(in):: lu
1785         type(sigio_head),intent(in):: head
1786         type(sigio_dbts),intent(inout):: dbts
1787         integer(sigio_intkind),intent(out):: iret
1788         integer:: i
1789         integer:: nc,mdim1
1790         integer:: iskip,iread,nread
1791         type(sigio_dats):: dats
1792     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1793         mdim1=min(size(dbts%hs,1),size(dbts%ps,1))
1794         nc=(head%jcap+1)*(head%jcap+2)
1795         iret=-5
1796         if(mdim1.lt.nc) return
1797         iret=-4
1798     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1799         if(head%irealf.eq.2) then
1800           iskip=0
1801           do i=1,head%nhead
1802             call bafrindex(0,iskip,head%lhead(i),iskip)
1803           enddo
1804           i=1
1805           iread=head%ldata(i)
1806           call bafrread(lu,iskip,iread,nread,dbts%hs)
1807           if(nread.lt.iread) return
1808           i=i+1
1809           iskip=iskip+nread
1810           iread=head%ldata(i)
1811           call bafrread(lu,iskip,iread,nread,dbts%ps)
1812           if(nread.lt.iread) return
1813         else
1814           call sigio_aldats(head,dats,iret)
1815           if(iret.ne.0) return
1816           call sigio_rrdats(lu,head,dats,iret)
1817           if(iret.ne.0) return
1818           dbts%hs(:nc)=dats%hs(:nc)
1819           dbts%ps(:nc)=dats%ps(:nc)
1820           call sigio_axdats(dats,iret)
1821         endif
1822         iret=0
1823     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1824       end subroutine
1825     !-------------------------------------------------------------------------------
1826       subroutine sigio_rwdbts(lu,head,dbts,iret)
1827         implicit none
1828         integer(sigio_intkind),intent(in):: lu
1829         type(sigio_head),intent(in):: head
1830         type(sigio_dbts),intent(in):: dbts
1831         integer(sigio_intkind),intent(out):: iret
1832         integer:: i
1833         integer:: nc,mdim1
1834         integer:: iskip,iwrite,nwrite
1835         type(sigio_dats):: dats
1836     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1837         mdim1=min(size(dbts%hs,1),size(dbts%ps,1))
1838         nc=(head%jcap+1)*(head%jcap+2)
1839         iret=-5
1840         if(mdim1.lt.nc) return
1841         iret=-4
1842     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1843         if(head%irealf.eq.2) then
1844           iskip=0
1845           do i=1,head%nhead
1846             call bafrindex(0,iskip,head%lhead(i),iskip)
1847           enddo
1848           i=1
1849           iwrite=head%ldata(i)
1850           call bafrwrite(lu,iskip,iwrite,nwrite,dbts%hs)
1851           if(nwrite.lt.iwrite) return
1852           i=i+1
1853           iskip=iskip+nwrite
1854           iwrite=head%ldata(i)
1855           call bafrwrite(lu,iskip,iwrite,nwrite,dbts%ps)
1856           if(nwrite.lt.iwrite) return
1857         else
1858           call sigio_aldats(head,dats,iret)
1859           if(iret.ne.0) return
1860           dats%hs(:nc)=dbts%hs(:nc)
1861           dats%ps(:nc)=dbts%ps(:nc)
1862           call sigio_rwdats(lu,head,dats,iret)
1863           if(iret.ne.0) return
1864           call sigio_axdats(dats,iret)
1865         endif
1866         iret=0
1867     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1868       end subroutine
1869     !-------------------------------------------------------------------------------
1870       subroutine sigio_rrdbtm(lu,head,dbtm,iret)
1871         implicit none
1872         integer(sigio_intkind),intent(in):: lu
1873         type(sigio_head),intent(in):: head
1874         type(sigio_dbtm),intent(inout):: dbtm
1875         integer(sigio_intkind),intent(out):: iret
1876         integer:: i,k,n
1877         integer:: nc,k1,k2,mdim1,ldim2,udim2,mdim3q
1878         integer:: iskip,iread,nread
1879         type(sigio_datm):: datm
1880     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1881         k1=dbtm%k1
1882         k2=dbtm%k2
1883         mdim1=min(size(dbtm%t,1),size(dbtm%d,1),size(dbtm%z,1),&
1884                   size(dbtm%q,1))
1885         ldim2=max(lbound(dbtm%t,2),lbound(dbtm%d,2),lbound(dbtm%z,2),&
1886                   lbound(dbtm%q,2))
1887         udim2=min(ubound(dbtm%t,2),ubound(dbtm%d,2),ubound(dbtm%z,2),&
1888                   ubound(dbtm%q,2))
1889         mdim3q=size(dbtm%q,3)
1890         nc=(head%jcap+1)*(head%jcap+2)
1891         iret=-5
1892         if(k1.lt.1.or.k2.gt.head%levs.or.&
1893            mdim1.lt.nc.or.&
1894            ldim2.gt.k1.or.udim2.lt.k2.or.&
1895            mdim3q.lt.head%ntrac) return
1896         iret=-4
1897     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1898         if(head%irealf.eq.2) then
1899           iskip=0
1900           do i=1,head%nhead
1901             call bafrindex(0,iskip,head%lhead(i),iskip)
1902           enddo
1903           i=1
1904           call bafrindex(0,iskip,head%ldata(i),iskip)
1905           i=i+1
1906           call bafrindex(0,iskip,head%ldata(i),iskip)
1907           do k=1,head%levs
1908             if(k.lt.k1.or.k.gt.k2) then
1909               i=i+1
1910               call bafrindex(0,iskip,head%ldata(i),iskip)
1911             else
1912               i=i+1
1913               iread=head%ldata(i)
1914               call bafrread(lu,iskip,iread,nread,dbtm%t(1,k))
1915               if(nread.lt.iread) return
1916               iskip=iskip+nread
1917             endif
1918           enddo
1919           do k=1,head%levs
1920             if(k.lt.k1.or.k.gt.k2) then
1921               i=i+1
1922               call bafrindex(0,iskip,head%ldata(i),iskip)
1923               i=i+1
1924               call bafrindex(0,iskip,head%ldata(i),iskip)
1925             else
1926               iread=head%ldata(i)
1927               call bafrread(lu,iskip,iread,nread,dbtm%d(1,k))
1928               if(nread.lt.iread) return
1929               iskip=iskip+nread
1930               iread=head%ldata(i)
1931               call bafrread(lu,iskip,iread,nread,dbtm%z(1,k))
1932               if(nread.lt.iread) return
1933               iskip=iskip+nread
1934             endif
1935           enddo
1936           do n=1,head%ntrac
1937             do k=1,head%levs
1938               if(k.lt.k1.or.k.gt.k2) then
1939                 i=i+1
1940                 call bafrindex(0,iskip,head%ldata(i),iskip)
1941               else
1942                 i=i+1
1943                 iread=head%ldata(i)
1944                 call bafrread(lu,iskip,iread,nread,dbtm%q(1,k,n))
1945                 if(nread.lt.iread) return
1946                 iskip=iskip+nread
1947               endif
1948             enddo
1949           enddo
1950         else
1951           call sigio_aldatm(head,k1,k2,datm,iret)
1952           if(iret.ne.0) return
1953           call sigio_rrdatm(lu,head,datm,iret)
1954           if(iret.ne.0) return
1955           dbtm%t(:nc,k1:k2)=datm%t(:nc,k1:k2)
1956           dbtm%d(:nc,k1:k2)=datm%d(:nc,k1:k2)
1957           dbtm%z(:nc,k1:k2)=datm%z(:nc,k1:k2)
1958           dbtm%q(:nc,k1:k2,:head%ntrac)=datm%q(:nc,k1:k2,:head%ntrac)
1959           call sigio_axdatm(datm,iret)
1960         endif
1961         iret=0
1962     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1963       end subroutine
1964     !-------------------------------------------------------------------------------
1965       subroutine sigio_rwdbtm(lu,head,dbtm,iret)
1966         implicit none
1967         integer(sigio_intkind),intent(in):: lu
1968         type(sigio_head),intent(in):: head
1969         type(sigio_dbtm),intent(in):: dbtm
1970         integer(sigio_intkind),intent(out):: iret
1971         integer:: i,k,n
1972         integer:: nc,k1,k2,mdim1,ldim2,udim2,mdim3q
1973         integer:: iskip,iwrite,nwrite
1974         type(sigio_datm):: datm
1975     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1976         k1=dbtm%k1
1977         k2=dbtm%k2
1978         mdim1=min(size(dbtm%t,1),size(dbtm%d,1),size(dbtm%z,1),&
1979                   size(dbtm%q,1))
1980         ldim2=max(lbound(dbtm%t,2),lbound(dbtm%d,2),lbound(dbtm%z,2),&
1981                   lbound(dbtm%q,2))
1982         udim2=min(ubound(dbtm%t,2),ubound(dbtm%d,2),ubound(dbtm%z,2),&
1983                   ubound(dbtm%q,2))
1984         mdim3q=size(dbtm%q,3)
1985         nc=(head%jcap+1)*(head%jcap+2)
1986         iret=-5
1987         if(k1.lt.1.or.k2.gt.head%levs.or.&
1988            mdim1.lt.nc.or.&
1989            ldim2.gt.k1.or.udim2.lt.k2.or.&
1990            mdim3q.lt.head%ntrac) return
1991         iret=-4
1992     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
1993         if(head%irealf.eq.2) then
1994           iskip=0
1995           do i=1,head%nhead
1996             call bafrindex(0,iskip,head%lhead(i),iskip)
1997           enddo
1998           i=1
1999           call bafrindex(0,iskip,head%ldata(i),iskip)
2000           i=i+1
2001           call bafrindex(0,iskip,head%ldata(i),iskip)
2002           do k=1,head%levs
2003             if(k.lt.k1.or.k.gt.k2) then
2004               i=i+1
2005               call bafrindex(0,iskip,head%ldata(i),iskip)
2006             else
2007               i=i+1
2008               iwrite=head%ldata(i)
2009               call bafrwrite(lu,iskip,iwrite,nwrite,dbtm%t(1,k))
2010               if(nwrite.lt.iwrite) return
2011               iskip=iskip+nwrite
2012             endif
2013           enddo
2014           do k=1,head%levs
2015             if(k.lt.k1.or.k.gt.k2) then
2016               i=i+1
2017               call bafrindex(0,iskip,head%ldata(i),iskip)
2018               i=i+1
2019               call bafrindex(0,iskip,head%ldata(i),iskip)
2020             else
2021               iwrite=head%ldata(i)
2022               call bafrwrite(lu,iskip,iwrite,nwrite,dbtm%d(1,k))
2023               if(nwrite.lt.iwrite) return
2024               iskip=iskip+nwrite
2025               iwrite=head%ldata(i)
2026               call bafrwrite(lu,iskip,iwrite,nwrite,dbtm%z(1,k))
2027               if(nwrite.lt.iwrite) return
2028               iskip=iskip+nwrite
2029             endif
2030           enddo
2031           do n=1,head%ntrac
2032             do k=1,head%levs
2033               if(k.lt.k1.or.k.gt.k2) then
2034                 i=i+1
2035                 call bafrindex(0,iskip,head%ldata(i),iskip)
2036               else
2037                 i=i+1
2038                 iwrite=head%ldata(i)
2039                 call bafrwrite(lu,iskip,iwrite,nwrite,dbtm%q(1,k,n))
2040                 if(nwrite.lt.iwrite) return
2041                 iskip=iskip+nwrite
2042               endif
2043             enddo
2044           enddo
2045         else
2046           call sigio_aldatm(head,k1,k2,datm,iret)
2047           if(iret.ne.0) return
2048           datm%t(:nc,k1:k2)=dbtm%t(:nc,k1:k2)
2049           datm%d(:nc,k1:k2)=dbtm%d(:nc,k1:k2)
2050           datm%z(:nc,k1:k2)=dbtm%z(:nc,k1:k2)
2051           datm%q(:nc,k1:k2,:head%ntrac)=dbtm%q(:nc,k1:k2,:head%ntrac)
2052           call sigio_rwdatm(lu,head,datm,iret)
2053           if(iret.ne.0) return
2054           call sigio_axdatm(datm,iret)
2055         endif
2056         iret=0
2057     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2058       end subroutine
2059     !-------------------------------------------------------------------------------
2060       subroutine sigio_rrdbti(lu,head,dbti,iret)
2061         implicit none
2062         integer(sigio_intkind),intent(in):: lu
2063         type(sigio_head),intent(in):: head
2064         type(sigio_dbti),intent(inout):: dbti
2065         integer(sigio_intkind),intent(out):: iret
2066         integer:: i,k,n
2067         integer:: mdim1
2068         integer:: mlen
2069         integer:: iskip,iread,nread
2070         type(sigio_dati):: dati
2071     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2072         i=dbti%i
2073         mdim1=size(dbti%f,1)
2074         iret=-5
2075         if(i.lt.1.or.i.gt.head%ndata) return
2076         mlen=head%ldata(i)/(4*head%irealf)
2077         if(mdim1.lt.mlen) return
2078         iret=-4
2079     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2080         if(head%irealf.eq.2) then
2081           iskip=0
2082           do i=1,head%nhead
2083             call bafrindex(0,iskip,head%lhead(i),iskip)
2084           enddo
2085           do i=1,dbti%i-1
2086             call bafrindex(0,iskip,head%ldata(i),iskip)
2087           enddo
2088           i=dbti%i
2089           iread=head%ldata(i)
2090           call bafrread(lu,iskip,iread,nread,dbti%f)
2091           if(nread.lt.iread) return
2092         else
2093           i=dbti%i
2094           call sigio_aldati(head,i,dati,iret)
2095           if(iret.ne.0) return
2096           call sigio_rrdati(lu,head,dati,iret)
2097           if(iret.ne.0) return
2098           dbti%f(:mlen)=dati%f(:mlen)
2099           call sigio_axdati(dati,iret)
2100         endif
2101         iret=0
2102     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2103       end subroutine
2104     !-------------------------------------------------------------------------------
2105       subroutine sigio_rwdbti(lu,head,dbti,iret)
2106         implicit none
2107         integer(sigio_intkind),intent(in):: lu
2108         type(sigio_head),intent(in):: head
2109         type(sigio_dbti),intent(in):: dbti
2110         integer(sigio_intkind),intent(out):: iret
2111         integer:: i,k,n
2112         integer:: mdim1
2113         integer:: mlen
2114         integer:: iskip,iwrite,nwrite
2115         type(sigio_dati):: dati
2116     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2117         i=dbti%i
2118         mdim1=size(dbti%f,1)
2119         iret=-5
2120         if(i.lt.1.or.i.gt.head%ndata) return
2121         mlen=head%ldata(i)/(4*head%irealf)
2122         if(mdim1.lt.mlen) return
2123         iret=-4
2124     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2125         if(head%irealf.eq.2) then
2126           iskip=0
2127           do i=1,head%nhead
2128             call bafrindex(0,iskip,head%lhead(i),iskip)
2129           enddo
2130           do i=1,dbti%i-1
2131             call bafrindex(0,iskip,head%ldata(i),iskip)
2132           enddo
2133           i=dbti%i
2134           iwrite=head%ldata(i)
2135           call bafrwrite(lu,iskip,iwrite,nwrite,dbti%f)
2136           if(nwrite.lt.iwrite) return
2137         else
2138           i=dbti%i
2139           call sigio_aldati(head,i,dati,iret)
2140           if(iret.ne.0) return
2141           dati%f(:mlen)=dbti%f(:mlen)
2142           call sigio_rwdati(lu,head,dati,iret)
2143           if(iret.ne.0) return
2144           call sigio_axdati(dati,iret)
2145         endif
2146         iret=0
2147     ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2148       end subroutine
2149     !-------------------------------------------------------------------------------
2150     end module
2151