<p><b>duda</b> 2012-03-27 15:47:50 -0600 (Tue, 27 Mar 2012)</p><p>BRANCH COMMIT<br>
<br>
Add code for reading and writing variable attributes at the stream level.<br>
<br>
<br>
M src/framework/mpas_io.F<br>
M src/framework/mpas_io_streams.F<br>
M src/framework/mpas_attlist.F<br>
</p><hr noshade><pre><font color="gray">Modified: branches/omp_blocks/io/src/framework/mpas_attlist.F
===================================================================
--- branches/omp_blocks/io/src/framework/mpas_attlist.F        2012-03-27 17:15:51 UTC (rev 1726)
+++ branches/omp_blocks/io/src/framework/mpas_attlist.F        2012-03-27 21:47:50 UTC (rev 1727)
@@ -31,13 +31,13 @@
end interface mpas_get_att
- !!!!! PRIVATE !!!!!
+ !!!!! PRIVATE? !!!!!
- integer, parameter, private :: ATT_INT = 1
- integer, parameter, private :: ATT_INTA = 2
- integer, parameter, private :: ATT_REAL = 3
- integer, parameter, private :: ATT_REALA = 4
- integer, parameter, private :: ATT_TEXT = 5
+ integer, parameter :: ATT_INT = 1
+ integer, parameter :: ATT_INTA = 2
+ integer, parameter :: ATT_REAL = 3
+ integer, parameter :: ATT_REALA = 4
+ integer, parameter :: ATT_TEXT = 5
contains
Modified: branches/omp_blocks/io/src/framework/mpas_io.F
===================================================================
--- branches/omp_blocks/io/src/framework/mpas_io.F        2012-03-27 17:15:51 UTC (rev 1726)
+++ branches/omp_blocks/io/src/framework/mpas_io.F        2012-03-27 21:47:50 UTC (rev 1727)
@@ -110,11 +110,11 @@
!!!!!!!! PRIVATE !!!!!!!!
- integer, parameter :: ATT_INT = 1
- integer, parameter :: ATT_INTA = 2
- integer, parameter :: ATT_REAL = 3
- integer, parameter :: ATT_REALA = 4
- integer, parameter :: ATT_TEXT = 5
+! integer, parameter :: ATT_INT = 1
+! integer, parameter :: ATT_INTA = 2
+! integer, parameter :: ATT_REAL = 3
+! integer, parameter :: ATT_REALA = 4
+! integer, parameter :: ATT_TEXT = 5
type decomphandle_type
integer :: field_type
@@ -3152,7 +3152,7 @@
end if
end if
- pio_ierr = PIO_put_att(handle % pio_file, varid, attName, attValue)
+ pio_ierr = PIO_put_att(handle % pio_file, varid, attName, trim(attValue))
if (pio_ierr /= PIO_noerr) then
if (present(ierr)) ierr = MPAS_IO_ERR_PIO
return
Modified: branches/omp_blocks/io/src/framework/mpas_io_streams.F
===================================================================
--- branches/omp_blocks/io/src/framework/mpas_io_streams.F        2012-03-27 17:15:51 UTC (rev 1726)
+++ branches/omp_blocks/io/src/framework/mpas_io_streams.F        2012-03-27 21:47:50 UTC (rev 1727)
@@ -173,7 +173,9 @@
return
end if
+ call put_get_field_atts(stream % fileHandle, stream % ioDirection, trim(field % fieldname), field % attList)
+
!
! Set field pointer and type in fieldList
!
@@ -289,7 +291,9 @@
return
end if
+ call put_get_field_atts(stream % fileHandle, stream % ioDirection, trim(field % fieldname), field % attList)
+
!
! Set field pointer and type in fieldList
!
@@ -405,7 +409,9 @@
return
end if
+ call put_get_field_atts(stream % fileHandle, stream % ioDirection, trim(field % fieldname), field % attList)
+
!
! Set field pointer and type in fieldList
!
@@ -523,7 +529,9 @@
return
end if
+ call put_get_field_atts(stream % fileHandle, stream % ioDirection, trim(field % fieldname), field % attList)
+
!
! Set field pointer and type in fieldList
!
@@ -590,14 +598,16 @@
call MPAS_streamAddField_generic(stream, trim(field % fieldName), MPAS_IO_DOUBLE, dimNames, dimSizes, &
field % hasTimeDimension, isDecomposed, totalDimSize, globalDimSize, indices, io_err)
+
+ deallocate(indices)
+ deallocate(dimSizes)
+ deallocate(dimNames)
if (io_err /= MPAS_STREAM_NOERR) then
if (present(ierr)) ierr = MPAS_IO_ERR
return
end if
- deallocate(indices)
- deallocate(dimSizes)
- deallocate(dimNames)
+ call put_get_field_atts(stream % fileHandle, stream % ioDirection, trim(field % fieldname), field % attList)
!
@@ -715,7 +725,9 @@
return
end if
+ call put_get_field_atts(stream % fileHandle, stream % ioDirection, trim(field % fieldname), field % attList)
+
!
! Set field pointer and type in fieldList
!
@@ -831,7 +843,9 @@
return
end if
+ call put_get_field_atts(stream % fileHandle, stream % ioDirection, trim(field % fieldname), field % attList)
+
!
! Set field pointer and type in fieldList
!
@@ -949,7 +963,9 @@
return
end if
+ call put_get_field_atts(stream % fileHandle, stream % ioDirection, trim(field % fieldname), field % attList)
+
!
! Set field pointer and type in fieldList
!
@@ -2252,4 +2268,55 @@
end subroutine mergeArrays
+
+ subroutine put_get_field_atts(handle, ioDirection, fieldname, attList)
+
+ implicit none
+
+ type (MPAS_IO_Handle_type), intent(inout) :: handle
+ integer, intent(in) :: ioDirection
+ character (len=*), intent(in) :: fieldname
+ type (att_list_type), pointer :: attList
+
+ type (att_list_type), pointer :: att_cursor
+
+ if (.not. associated(attList)) return
+
+ att_cursor => attList
+ if (ioDirection == MPAS_IO_WRITE) then
+ do while (associated(att_cursor))
+ select case (att_cursor % attType)
+ case (ATT_INT)
+ call MPAS_io_put_att(handle, trim(att_cursor % attName), att_cursor % attValueInt, fieldname)
+ case (ATT_INTA)
+ call MPAS_io_put_att(handle, trim(att_cursor % attName), att_cursor % attValueIntA, fieldname)
+ case (ATT_REAL)
+ call MPAS_io_put_att(handle, trim(att_cursor % attName), att_cursor % attValueReal, fieldname)
+ case (ATT_REALA)
+ call MPAS_io_put_att(handle, trim(att_cursor % attName), att_cursor % attValueRealA, fieldname)
+ case (ATT_TEXT)
+ call MPAS_io_put_att(handle, trim(att_cursor % attName), att_cursor % attValueText, fieldname)
+ end select
+ att_cursor => att_cursor % next
+ end do
+ else
+ do while (associated(att_cursor))
+ select case (att_cursor % attType)
+ case (ATT_INT)
+ call MPAS_io_get_att(handle, trim(att_cursor % attName), att_cursor % attValueInt, fieldname)
+ case (ATT_INTA)
+ call MPAS_io_get_att(handle, trim(att_cursor % attName), att_cursor % attValueIntA, fieldname)
+ case (ATT_REAL)
+ call MPAS_io_get_att(handle, trim(att_cursor % attName), att_cursor % attValueReal, fieldname)
+ case (ATT_REALA)
+ call MPAS_io_get_att(handle, trim(att_cursor % attName), att_cursor % attValueRealA, fieldname)
+ case (ATT_TEXT)
+ call MPAS_io_get_att(handle, trim(att_cursor % attName), att_cursor % attValueText, fieldname)
+ end select
+ att_cursor => att_cursor % next
+ end do
+ end if
+
+ end subroutine put_get_field_atts
+
end module mpas_io_streams
</font>
</pre>