mpp/include/mpp_io_util.inc

   1 ! -*-f90-*-
   2
   3
   4 !***********************************************************************
   5 !*                   GNU Lesser General Public License
   6 !*
   7 !* This file is part of the GFDL Flexible Modeling System (FMS).
   8 !*
   9 !* FMS is free software: you can redistribute it and/or modify it under
  10 !* the terms of the GNU Lesser General Public License as published by
  11 !* the Free Software Foundation, either version 3 of the License, or (at
  12 !* your option) any later version.
  13 !*
  14 !* FMS is distributed in the hope that it will be useful, but WITHOUT
  15 !* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  16 !* FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  17 !* for more details.
  18 !*
  19 !* You should have received a copy of the GNU Lesser General Public
  20 !* License along with FMS.  If not, see <http://www.gnu.org/licenses/>.
  21 !***********************************************************************
  22
  23 !> @file
  24 !> @ingroup mpp_io_mod
  25 !> @brief Routines to retrieve data used in @ref mpp_io_mod
  26
  27   !#####################################################################
  28 ! <SUBROUTINE NAME="mpp_get_info">
  29 !   <OVERVIEW>
  30 !     Get some general information about a file.
  31 !   </OVERVIEW>
  32 !   <DESCRIPTION>
  33 !     Get some general information about a file.
  34 !   </DESCRIPTION>
  35 !   <TEMPLATE>
  36 !     call mpp_get_info( unit, ndim, nvar, natt, ntime )
  37 !   </TEMPLATE>
  38 !   <IN NAME="unit" TYPE="integer"> </IN>
  39 !   <OUT NAME="ndim" TYPE="integer"> </OUT>
  40 !   <OUT NAME="nvar" TYPE="integer"> </OUT>
  41 !   <OUT NAME="natt" TYPE="integer"> </OUT>
  42 !   <OUT NAME="ntime" TYPE="integer"> </OUT>
  43 ! </SUBROUTINE>
  44
  45     !> @brief Get some general information about a file.
  46     !!
  47     !> <br>Example usage:
  48     !! @code{.F90}
  49     !! call mpp_get_info( unit, ndim, nvar, natt, ntime )
  50     !! @endcode
  51     subroutine mpp_get_info( unit, ndim, nvar, natt, ntime )
  52
  53       integer, intent(in) :: unit
  54       integer, intent(out) :: ndim, nvar, natt, ntime
  55
  56
  57       if( .NOT.module_is_initialized )call mpp_error( FATAL, 'MPP_GET_INFO: must first call mpp_io_init.' )
  58       if( .NOT.mpp_file(unit)%opened )&
  59            call mpp_error(FATAL, 'MPP_GET_INFO: invalid unit number, file '//trim(mpp_file(unit)%name))
  60
  61       ndim = mpp_file(unit)%ndim
  62       nvar = mpp_file(unit)%nvar
  63       natt = mpp_file(unit)%natt
  64       ntime = mpp_file(unit)%time_level
  65
  66       return
  67
  68     end subroutine mpp_get_info
  69
  70   !#####################################################################
  71 ! <SUBROUTINE NAME="mpp_get_global_atts" INTERFACE="mpp_get_atts">
  72 !  <IN NAME="unit" TYPE="integer"></IN>
  73 !  <IN NAME="global_atts" TYPE="atttype" DIM="(:)"></IN>
  74 ! </SUBROUTINE>
  75 !
  76     !> @brief Copy global file attributes for use by user
  77     subroutine mpp_get_global_atts( unit, global_atts )
  78
  79       integer,       intent(in)    :: unit
  80       type(atttype), intent(inout) :: global_atts(:) !< an attribute type which is allocated from the
  81                                                      !!  calling routine
  82
  83       integer :: natt,i
  84
  85       if( .NOT.module_is_initialized )call mpp_error( FATAL, 'MPP_GET_INFO: must first call mpp_io_init.' )
  86       if( .NOT.mpp_file(unit)%opened )&
  87            call mpp_error( FATAL, 'MPP_GET_INFO: invalid unit number,file '//trim(mpp_file(unit)%name))
  88
  89       if (size(global_atts(:)).lt.mpp_file(unit)%natt) &
  90            call mpp_error(FATAL, 'MPP_GET_ATTS: atttype not dimensioned properly in calling routine, file '// &
  91            trim(mpp_file(unit)%name))
  92
  93       natt = mpp_file(unit)%natt
  94       global_atts = default_att
  95
  96       do i=1,natt
  97          global_atts(i) = mpp_file(unit)%Att(i)
  98       enddo
  99
 100       return
 101    end subroutine mpp_get_global_atts
 102
 103   !#####################################################################
 104    subroutine mpp_get_field_atts(field, name, units, longname, min, max, missing, ndim, siz, axes, atts, &
 105                                  valid, scale, add, checksum)
 106
 107      type(fieldtype),    intent(in)                            :: field
 108      character(len=*),   intent(out),                 optional :: name, units
 109      character(len=*),   intent(out),                 optional :: longname
 110      real,               intent(out),                 optional :: min,max,missing
 111      integer,            intent(out),                 optional :: ndim
 112      integer,            intent(out),   dimension(:), optional :: siz
 113      type(validtype),    intent(out),                 optional :: valid
 114      real,               intent(out),                 optional :: scale
 115      real,               intent(out),                 optional :: add
 116      integer(i8_kind), intent(out),   dimension(:), optional :: checksum
 117
 118      type(atttype),      intent(inout), dimension(:), optional :: atts
 119      type(axistype),     intent(inout), dimension(:), optional :: axes
 120
 121      integer :: n,m, check_exist
 122
 123      if (PRESENT(name)) name = field%name
 124      if (PRESENT(units)) units = field%units
 125      if (PRESENT(longname)) longname = field%longname
 126      if (PRESENT(min)) min = field%min
 127      if (PRESENT(max)) max = field%max
 128      if (PRESENT(missing)) missing = field%missing
 129      if (PRESENT(ndim)) ndim = field%ndim
 130      if (PRESENT(atts)) then
 131         atts = default_att
 132         n = size(atts(:));m=size(field%Att(:))
 133         if (n.LT.m)&
 134              call mpp_error(FATAL,'attribute array not large enough in mpp_get_field_atts, field '//&
 135              trim(field%name))
 136         do n=1,m
 137           atts(n) = field%Att(n)
 138         end do
 139      end if
 140      if (PRESENT(axes)) then
 141         axes = default_axis
 142         n = size(axes(:));m=field%ndim
 143         if (n.LT.m) &
 144              call mpp_error(FATAL,'axis array not large enough in mpp_get_field_atts, field '//&
 145              trim(field%name))
 146         do n=1,m
 147           axes(n) = field%axes(n)
 148         end do
 149      end if
 150      if (PRESENT(siz)) then
 151         siz = -1
 152         n = size(siz(:));m=field%ndim
 153         if (n.LT.m) &
 154              call mpp_error(FATAL,'size array not large enough in mpp_get_field_atts, field '//&
 155              trim(field%name))
 156         do n=1,m
 157           siz(n) = field%size(n)
 158         end do
 159      end if
 160
 161      if(PRESENT(valid)) then
 162        call mpp_get_valid(field,valid)
 163      endif
 164
 165      if(PRESENT(scale))    scale    = field%scale
 166      if(present(add))      add      = field%add
 167      if(present(checksum)) then
 168        checksum = 0
 169        check_exist = mpp_find_att(field%Att(:),"checksum")
 170        if ( check_exist >= 0 ) then
 171          if(size(checksum(:)) >size(field%checksum(:))) call mpp_error(FATAL, &
 172             & "size(checksum(:)) >size(field%checksum(:))")
 173          checksum = field%checksum(1:size(checksum(:)))
 174        endif
 175      endif
 176
 177      return
 178    end subroutine mpp_get_field_atts
 179
 180   !#####################################################################
 181    subroutine mpp_get_axis_atts( axis, name, units, longname, cartesian, &
 182                                  calendar, sense, len, natts, atts, compressed )
 183
 184      type(axistype), intent(in) :: axis
 185      character(len=*), intent(out) , optional :: name, units
 186      character(len=*), intent(out), optional :: longname, cartesian
 187      character(len=*), intent(out), optional :: compressed, calendar
 188      integer,intent(out), optional :: sense, len , natts
 189      type(atttype), intent(inout), optional, dimension(:) :: atts
 190
 191      integer :: n,m
 192
 193      if (PRESENT(name)) name = axis%name
 194      if (PRESENT(units)) units = axis%units
 195      if (PRESENT(longname)) longname = axis%longname
 196      if (PRESENT(cartesian)) cartesian = axis%cartesian
 197      if (PRESENT(compressed)) compressed = axis%compressed
 198      if (PRESENT(calendar)) calendar = axis%calendar
 199      if (PRESENT(sense)) sense = axis%sense
 200      if (PRESENT(len)) len = axis%len
 201      if (PRESENT(atts)) then
 202         atts = default_att
 203         n = size(atts(:));m=size(axis%Att(:))
 204         if (n.LT.m) &
 205              call mpp_error(FATAL,'attribute array not large enough in mpp_get_field_atts, axis '//&
 206              trim(axis%name))
 207         do n=1,m
 208           atts(n) = axis%Att(n)
 209         end do
 210      end if
 211      if (PRESENT(natts)) natts = size(axis%Att(:))
 212
 213      return
 214    end subroutine mpp_get_axis_atts
 215
 216
 217   !#####################################################################
 218     !> @brief Copy variable information from file (excluding data)
 219     subroutine mpp_get_fields( unit, variables )
 220       integer,         intent(in)    :: unit
 221       type(fieldtype), intent(inout) :: variables(:)
 222
 223       integer :: nvar,i
 224
 225       if( .NOT.module_is_initialized )call mpp_error( FATAL, 'MPP_GET_FIELDS: must first call mpp_io_init.' )
 226       if( .NOT.mpp_file(unit)%opened )call mpp_error( FATAL, 'MPP_GET_FIELDS: invalid unit number.' )
 227
 228       if (size(variables(:)).ne.mpp_file(unit)%nvar) &
 229           call mpp_error(FATAL,'MPP_GET_FIELDS: fieldtype not dimensioned properly in calling routine, file '//&
 230           trim(mpp_file(unit)%name))
 231
 232       nvar = mpp_file(unit)%nvar
 233
 234       do i=1,nvar
 235          variables(i) = mpp_file(unit)%Var(i)
 236       enddo
 237
 238       return
 239    end subroutine mpp_get_fields
 240
 241
 242
 243   !#####################################################################
 244     !> @brief Copy variable information from file (excluding data)
 245     subroutine mpp_get_axes( unit, axes, time_axis )
 246       integer, intent(in) :: unit
 247       type(axistype), intent(inout) :: axes(:)
 248       type(axistype), intent(inout), optional :: time_axis
 249       integer :: ndim,i
 250
 251       if( .NOT.module_is_initialized )call mpp_error( FATAL, 'MPP_GET_AXES: must first call mpp_io_init.' )
 252       if( .NOT.mpp_file(unit)%opened )&
 253            call mpp_error( FATAL, 'MPP_GET_AXES: invalid unit number, file '//trim(mpp_file(unit)%name))
 254
 255       if (size(axes(:)).ne.mpp_file(unit)%ndim) &
 256            call mpp_error(FATAL, 'MPP_GET_AXES: axistype not dimensioned properly in calling routine, file '//&
 257            trim(mpp_file(unit)%name))
 258
 259
 260       if (PRESENT(time_axis)) time_axis = default_axis
 261       ndim = mpp_file(unit)%ndim
 262
 263       do i=1,ndim
 264          axes(i)=mpp_file(unit)%Axis(i)
 265
 266          if (PRESENT(time_axis) &
 267              .AND. .NOT. ASSOCIATED(mpp_file(unit)%Axis(i)%data) &
 268              .AND. mpp_file(unit)%Axis(i)%type /= -1) then
 269             time_axis = mpp_file(unit)%Axis(i)
 270          endif
 271       enddo
 272
 273       return
 274    end subroutine mpp_get_axes
 275
 276   !#####################################################################
 277    !> @brief Copy variable information from file (excluding data)
 278    function mpp_get_dimension_length(unit, dimname, found)
 279      integer,           intent(in)  :: unit
 280      character(len=*),  intent(in)  :: dimname
 281      logical, optional, intent(out) :: found
 282      integer                        :: mpp_get_dimension_length
 283      logical                        :: found_dim
 284      integer                        :: i
 285
 286
 287      if( .NOT.module_is_initialized ) &
 288        call mpp_error( FATAL, 'mpp_get_dimension_length: must first call mpp_io_init.' )
 289      if( .NOT.mpp_file(unit)%opened )&
 290        call mpp_error( FATAL, 'mpp_get_dimension_length: invalid unit number, file '//trim(mpp_file(unit)%name))
 291      found_dim = .false.
 292      mpp_get_dimension_length = -1
 293      do i = 1, mpp_file(unit)%ndim
 294         if(trim(dimname) == trim(mpp_file(unit)%Axis(i)%name)) then
 295           mpp_get_dimension_length = mpp_file(unit)%Axis(i)%len
 296           found_dim = .true.
 297           exit
 298         endif
 299      enddo
 300
 301      if(present(found)) found = found_dim
 302
 303    end function mpp_get_dimension_length
 304
 305    !#####################################################################
 306     !> @brief Copy variable information from file (excluding data)
 307     subroutine mpp_get_time_axis( unit, time_axis )
 308       integer, intent(in) :: unit
 309       type(axistype), intent(inout) :: time_axis
 310
 311       if( .NOT.module_is_initialized )call mpp_error( FATAL, 'MPP_GET_AXES: must first call mpp_io_init.' )
 312       if( .NOT.mpp_file(unit)%opened )&
 313            call mpp_error( FATAL, 'MPP_GET_AXES: invalid unit number, file '//trim(mpp_file(unit)%name))
 314
 315       time_axis = mpp_file(unit)%Axis(mpp_file(unit)%recdimid)
 316
 317       return
 318    end subroutine mpp_get_time_axis
 319
 320   !####################################################################
 321   !> @brief Copy variable information from file (excluding data)
 322   function mpp_get_default_calendar( )
 323      character(len=len(default_axis%calendar)) :: mpp_get_default_calendar
 324
 325      mpp_get_default_calendar = default_axis%calendar
 326
 327   end function mpp_get_default_calendar
 328
 329    !> @brief Get file time data.
 330    !!
 331    !> Copy time information from file and convert to time_type
 332    !! <br>Example usage:
 333    !! @code{.F90}
 334    !! call mpp_get_times( unit, time_values)
 335    !! @endcode
 336    subroutine mpp_get_times( unit, time_values )
 337       integer, intent(in) :: unit
 338       real, intent(inout) :: time_values(:)
 339
 340       integer :: ntime,i
 341
 342       if( .NOT.module_is_initialized )call mpp_error( FATAL, 'MPP_GET_TIMES: must first call mpp_io_init.' )
 343       if( .NOT.mpp_file(unit)%opened )&
 344            call mpp_error(FATAL, 'MPP_GET_TIMES: invalid unit number, file '//trim(mpp_file(unit)%name))
 345
 346 ! NF_INQ_DIM returns -1 for the length of a record dimension if
 347 ! it does not exist
 348
 349       if (mpp_file(unit)%time_level == -1) then
 350           time_values = 0.0
 351           return
 352       endif
 353
 354       if (size(time_values(:)).ne.mpp_file(unit)%time_level) &
 355          call mpp_error(FATAL,'MPP_GET_TIMES: time_values not dimensioned properly in calling routine, file '//&
 356          trim(mpp_file(unit)%name))
 357
 358       ntime = mpp_file(unit)%time_level
 359
 360       do i=1,ntime
 361          time_values(i) = mpp_file(unit)%time_values(i)
 362       enddo
 363
 364       return
 365     end subroutine mpp_get_times
 366
 367   !#####################################################################
 368    function mpp_get_field_index(fields,fieldname)
 369
 370      type(fieldtype), dimension(:) :: fields
 371      character(len=*) :: fieldname
 372      integer :: mpp_get_field_index
 373
 374      integer :: n
 375
 376      mpp_get_field_index = -1
 377
 378      do n=1,size(fields(:))
 379         if (lowercase(fields(n)%name) == lowercase(fieldname)) then
 380            mpp_get_field_index = n
 381            exit
 382         endif
 383      enddo
 384
 385      return
 386    end function mpp_get_field_index
 387
 388   !#####################################################################
 389    function mpp_get_axis_index(axes,axisname)
 390
 391      type(axistype), dimension(:) :: axes
 392      character(len=*) :: axisname
 393      integer :: mpp_get_axis_index
 394
 395      integer :: n
 396
 397      mpp_get_axis_index = -1
 398
 399      do n=1,size(axes(:))
 400         if (lowercase(axes(n)%name) == lowercase(axisname)) then
 401            mpp_get_axis_index = n
 402            exit
 403         endif
 404      enddo
 405
 406      return
 407    end function mpp_get_axis_index
 408
 409   !#####################################################################
 410    function mpp_get_axis_by_name(unit,axisname)
 411
 412      integer          :: unit
 413      character(len=*) :: axisname
 414      type(axistype)   :: mpp_get_axis_by_name
 415
 416      integer :: n
 417
 418      mpp_get_axis_by_name = default_axis
 419
 420      do n=1,size(mpp_file(unit)%Axis(:))
 421         if (lowercase(mpp_file(unit)%Axis(n)%name) == lowercase(axisname)) then
 422            mpp_get_axis_by_name = mpp_file(unit)%Axis(n)
 423            exit
 424         endif
 425      enddo
 426
 427      return
 428    end function mpp_get_axis_by_name
 429
 430   !#####################################################################
 431    function mpp_get_field_size(field)
 432
 433      type(fieldtype) :: field
 434      integer :: mpp_get_field_size(4)
 435
 436      mpp_get_field_size = -1
 437
 438      mpp_get_field_size(1) = field%size(1)
 439      mpp_get_field_size(2) = field%size(2)
 440      mpp_get_field_size(3) = field%size(3)
 441      mpp_get_field_size(4) = field%size(4)
 442
 443      return
 444    end function mpp_get_field_size
 445
 446
 447   !#####################################################################
 448    function mpp_get_axis_length(axis)
 449
 450      type(axistype) :: axis
 451      integer :: mpp_get_axis_length
 452
 453      mpp_get_axis_length = axis%len
 454
 455      return
 456    end function mpp_get_axis_length
 457
 458   !#####################################################################
 459   function mpp_get_axis_bounds(axis, data, name)
 460      type(axistype), intent(in) :: axis
 461      real, dimension(:), intent(out) :: data
 462      character(len=*), optional, intent(out) :: name
 463      logical                         :: mpp_get_axis_bounds
 464
 465      if (size(data(:)).lt.axis%len+1)&
 466           call mpp_error(FATAL,'MPP_GET_AXIS_BOUNDS: data array not large enough, axis '//trim(axis%name))
 467      if (.NOT.ASSOCIATED(axis%data_bounds)) then
 468         mpp_get_axis_bounds = .false.
 469      else
 470         mpp_get_axis_bounds = .true.
 471         data(1:axis%len+1) = axis%data_bounds(:)
 472      endif
 473      if(present(name)) name = trim(axis%name_bounds)
 474
 475      return
 476   end function mpp_get_axis_bounds
 477
 478   !#####################################################################
 479    subroutine mpp_get_axis_data( axis, data )
 480
 481      type(axistype), intent(in) :: axis
 482      real, dimension(:), intent(out) :: data
 483
 484
 485      if (size(data(:)).lt.axis%len)&
 486           call mpp_error(FATAL,'MPP_GET_AXIS_DATA: data array not large enough, axis '//trim(axis%name))
 487      if (.NOT.ASSOCIATED(axis%data)) then
 488         call mpp_error(NOTE,'MPP_GET_AXIS_DATA: use mpp_get_times for record dims')
 489         data = 0.
 490      else
 491         data(1:axis%len) = axis%data
 492      endif
 493
 494      return
 495    end subroutine mpp_get_axis_data
 496
 497
 498   !#####################################################################
 499    function mpp_get_recdimid(unit)
 500 !
 501       integer, intent(in) :: unit
 502       integer  :: mpp_get_recdimid
 503
 504
 505       if( .NOT.module_is_initialized )call mpp_error( FATAL, 'MPP_GET_RECDIMID: must first call mpp_io_init.' )
 506       if( .NOT.mpp_file(unit)%opened )call mpp_error( FATAL, 'MPP_GET_RECDIMID: invalid unit number.' )
 507
 508       mpp_get_recdimid = mpp_file(unit)%recdimid
 509
 510       return
 511    end function mpp_get_recdimid
 512
 513   !#####################################################################
 514 ! <FUNCTION NAME="mpp_get_ncid">
 515 !   <OVERVIEW>
 516 !     Get netCDF ID of an open file.
 517 !   </OVERVIEW>
 518 !   <DESCRIPTION>
 519 !    This returns the <TT>ncid</TT> associated with the open file on
 520 !    <TT>unit</TT>. It is used in the instance that the user desires to
 521 !    perform netCDF calls upon the file that are not provided by the
 522 !    <TT>mpp_io_mod</TT> API itself.
 523 !   </DESCRIPTION>
 524 !   <TEMPLATE>
 525 !     mpp_get_ncid(unit)
 526 !   </TEMPLATE>
 527 !   <IN NAME="unit" TYPE="integer"> </IN>
 528 ! </FUNCTION>
 529
 530     !> @brief Get netCDF ID of an open file.
 531     !!
 532     !> This returns the <TT>ncid</TT> associated with the open file on
 533     !! <TT>unit</TT>. It is used in the instance that the user desires to
 534     !! perform netCDF calls upon the file that are not provided by the
 535     !! <TT>mpp_io_mod</TT> API itself.
 536     function mpp_get_ncid(unit)
 537       integer :: mpp_get_ncid
 538       integer, intent(in) :: unit
 539
 540       mpp_get_ncid = mpp_file(unit)%ncid
 541       return
 542     end function mpp_get_ncid
 543
 544   !#####################################################################
 545     function mpp_get_axis_id(axis)
 546       integer mpp_get_axis_id
 547       type(axistype), intent(in) :: axis
 548       mpp_get_axis_id = axis%id
 549       return
 550     end function mpp_get_axis_id
 551
 552   !#####################################################################
 553     function mpp_get_field_id(field)
 554       integer mpp_get_field_id
 555       type(fieldtype), intent(in) :: field
 556       mpp_get_field_id = field%id
 557       return
 558     end function mpp_get_field_id
 559
 560   !#####################################################################
 561     !> @brief Set the mpp_io_stack variable to be at least n LONG words long
 562     subroutine mpp_io_set_stack_size(n)
 563       integer, intent(in) :: n
 564       character(len=10) :: text
 565
 566       if( n.GT.mpp_io_stack_size .AND. allocated(mpp_io_stack) )deallocate(mpp_io_stack)
 567       if( .NOT.allocated(mpp_io_stack) )then
 568           allocate( mpp_io_stack(n) )
 569           mpp_io_stack_size = n
 570           write( text,'(i10)' )n
 571           if( pe.EQ.mpp_root_pe() )call mpp_error( NOTE, 'MPP_IO_SET_STACK_SIZE: stack size set to '//text//'.' )
 572       end if
 573
 574       return
 575     end subroutine mpp_io_set_stack_size
 576
 577   !#####################################################################
 578   !> Based on presence/absence of attributes, defines valid range or missing
 579   ! value. For details, see section 8.1 of NetCDF User Guide
 580   subroutine mpp_get_valid(f,v)
 581      type(fieldtype),intent(in)  :: f ! field
 582      type(validtype),intent(out) :: v ! validator
 583
 584      integer :: irange,imin,imax,ifill,imissing,iscale
 585      integer :: valid_T, scale_T ! types of attributes
 586
 587      v%is_range = .true.
 588      v%min = -HUGE(v%min); v%max = HUGE(v%max)
 589      if (f%natt == 0) return
 590      ! find indices of relevant attributes
 591      irange   = mpp_find_att(f%att,'valid_range')
 592      imin     = mpp_find_att(f%att,'valid_min')
 593      imax     = mpp_find_att(f%att,'valid_max')
 594      ifill    = mpp_find_att(f%att,'_FillValue')
 595      imissing = mpp_find_att(f%att,'missing_value')
 596
 597      ! find the widest type of scale and offset; note that the code
 598      ! uses assumption that NetCDF types are arranged in th order of rank,
 599      ! that is NF_BYTE < NF_SHORT < NF_INT < NF_FLOAT < NF_DOUBLE
 600      scale_T = 0
 601      iscale   = mpp_find_att(f%att,'scale_factor')
 602      if(iscale>0) scale_T = f%att(iscale)%type
 603      iscale = mpp_find_att(f%att,'add_offset')
 604      if(iscale>0) scale_T = max(scale_T,f%att(iscale)%type)
 605
 606
 607      ! examine possible range attributes
 608      valid_T = 0
 609      if (irange>0) then
 610         v%min = f%att(irange)%fatt(1)
 611         v%max = f%att(irange)%fatt(2)
 612         valid_T = f%att(irange)%type
 613      else if (imax>0.or.imin>0) then
 614         if(imax>0) then
 615            v%max = f%att(imax)%fatt(1)
 616            valid_T = max(valid_T,f%att(imax)%type)
 617         endif
 618         if(imin>0) then
 619            v%min = f%att(imin)%fatt(1)
 620            valid_T = max(valid_T,f%att(imin)%type)
 621         endif
 622      else if (imissing > 0) then
 623         v%is_range = .false.
 624         ! here we always scale, since missing_value is supposed to be in
 625         ! external representation
 626         v%min = f%att(imissing)%fatt(1)*f%scale + f%add
 627      else if (ifill>0) then
 628      !z1l ifdef is added in to be able to compile without using use_netCDF.
 629 #ifdef use_netCDF
 630         ! define min and max according to _FillValue
 631         if(f%att(ifill)%fatt(1)>0) then
 632             ! if _FillValue is positive, then it defines valid maximum
 633             v%max = f%att(ifill)%fatt(1)
 634             select case(f%type)
 635             case (NF_BYTE,NF_SHORT,NF_INT)
 636                v%max = v%max-1
 637             case (NF_FLOAT)
 638                v%max = nearest(nearest(real(v%max,4),-1.0),-1.0)
 639             case (NF_DOUBLE)
 640                v%max = nearest(nearest(real(v%max,8),-1.0),-1.0)
 641             end select
 642             ! always do the scaling, as the _FillValue is in external
 643             ! representation
 644             v%max = v%max*f%scale + f%add
 645         else
 646             ! if _FillValue is negative or zero, then it defines valid minimum
 647             v%min = f%att(ifill)%fatt(1)
 648             select case(f%type)
 649             case (NF_BYTE,NF_SHORT,NF_INT)
 650                v%min = v%min+1
 651             case (NF_FLOAT)
 652                v%min = nearest(nearest(real(v%min,4),+1.0),+1.0)
 653             case (NF_DOUBLE)
 654                v%min = nearest(nearest(real(v%min,8),+1.0),+1.0)
 655             end select
 656             ! always do the scaling, as the _FillValue is in external
 657             ! representation
 658             v%min = v%min*f%scale + f%add
 659         endif
 660 #endif
 661     endif
 662    ! If valid_range is the same type as scale_factor (actually the wider of
 663    ! scale_factor and add_offset) and this is wider than the external data, then it
 664    ! will be interpreted as being in the units of the internal (unpacked) data.
 665    ! Otherwise it is in the units of the external (packed) data.
 666     ! Note that it is not relevant if we went through the missing_data of _FillValue
 667     ! brances, because in this case all irange, imin, and imax are less then 0
 668     if(.not.((valid_T == scale_T).and.(scale_T>f%type))) then
 669        if(irange>0 .or. imin>0) then
 670           v%min = v%min*f%scale + f%add
 671        endif
 672        if(irange>0 .or. imax>0) then
 673           v%max = v%max*f%scale + f%add
 674        endif
 675     endif
 676
 677    end subroutine mpp_get_valid
 678
 679    !#####################################################################
 680    logical elemental function mpp_is_valid(x, v)
 681       real           , intent(in) :: x ! real value to be eaxmined
 682       type(validtype), intent(in) :: v ! validator
 683
 684       if (v%is_range) then
 685          mpp_is_valid = (v%min<=x).and.(x<=v%max)
 686       else
 687          mpp_is_valid = x/=v%min
 688       endif
 689     end function mpp_is_valid
 690
 691    !#####################################################################
 692    ! finds an attribute by name in the array; returns -1 if it is not
 693    ! found
 694    function mpp_find_att(atts, name)
 695      integer                   :: mpp_find_att
 696      type(atttype), intent(in) :: atts(:) ! array of attributes
 697      character(len=*)          :: name ! name of the attributes
 698
 699      integer :: i
 700
 701      mpp_find_att = -1
 702      do i = 1, size(atts)
 703         if (trim(name)==trim(atts(i)%name)) then
 704            mpp_find_att=i
 705            exit
 706         endif
 707      enddo
 708    end function mpp_find_att
 709    !#####################################################################
 710
 711    !> @brief return the name of an attribute.
 712    function mpp_get_att_name(att)
 713       type(atttype),    intent(in) :: att
 714       character(len=len(att%name)) :: mpp_get_att_name
 715
 716       mpp_get_att_name = att%name
 717       return
 718
 719    end function mpp_get_att_name
 720
 721    !#####################################################################
 722
 723    !> @brief return the type of an attribute.
 724    function mpp_get_att_type(att)
 725       type(atttype), intent(in) :: att
 726       integer                   :: mpp_get_att_type
 727
 728       mpp_get_att_type = att%type
 729       return
 730
 731    end function mpp_get_att_type
 732
 733    !#####################################################################
 734
 735    !> @brief return the length of an attribute.
 736    function mpp_get_att_length(att)
 737       type(atttype), intent(in) :: att
 738       integer                   :: mpp_get_att_length
 739
 740       mpp_get_att_length = att%len
 741
 742       return
 743
 744    end function mpp_get_att_length
 745
 746    !#####################################################################
 747
 748    !> @brief return the char value of an attribute.
 749    function mpp_get_att_char(att)
 750       type(atttype), intent(in) :: att
 751       character(len=att%len)    :: mpp_get_att_char
 752
 753       mpp_get_att_char = att%catt
 754       return
 755
 756    end function mpp_get_att_char
 757
 758    !#####################################################################
 759
 760    !> @brief return the real array value of an attribute.
 761    function mpp_get_att_real(att)
 762       type(atttype), intent(in)          :: att
 763       real, dimension(size(att%fatt(:))) :: mpp_get_att_real
 764
 765       mpp_get_att_real = att%fatt
 766       return
 767
 768    end function mpp_get_att_real
 769
 770   !#####################################################################
 771
 772    !> @brief return the real array value of an attribute.
 773    function mpp_get_att_real_scalar(att)
 774       type(atttype), intent(in)          :: att
 775       real                               :: mpp_get_att_real_scalar
 776
 777       mpp_get_att_real_scalar = att%fatt(1)
 778       return
 779
 780    end function mpp_get_att_real_scalar
 781
 782    !#####################################################################
 783    !> @brief return the name of an field
 784    function mpp_get_field_name(field)
 785       type(fieldtype), intent(in) :: field
 786       character(len=len(field%name)) :: mpp_get_field_name
 787
 788       mpp_get_field_name = field%name
 789       return
 790    end function mpp_get_field_name
 791
 792    !#####################################################################
 793    !> @brief return the  file name of corresponding unit
 794    function mpp_get_file_name(unit)
 795       integer,                  intent(in) :: unit
 796       character(len=len(mpp_file(1)%name)) :: mpp_get_file_name
 797
 798       mpp_get_file_name = mpp_file(unit)%name
 799       return
 800
 801    end function mpp_get_file_name
 802
 803    !####################################################################
 804    !> @brief return if certain file with unit is opened or not
 805    function mpp_file_is_opened(unit)
 806       integer,  intent(in) :: unit
 807       logical              :: mpp_file_is_opened
 808
 809       mpp_file_is_opened = mpp_file(unit)%opened
 810       return
 811
 812    end function mpp_file_is_opened
 813
 814    !####################################################################
 815    !> @brief return the attribute value of given field name
 816    subroutine mpp_get_field_att_text(unit, fieldname, attname, attvalue)
 817      integer,           intent(in) :: unit
 818      character(len=*),  intent(in) :: fieldname, attname
 819      character(len=*), intent(out) :: attvalue
 820      logical                       :: found_field,  found_att
 821      integer                       :: i, j, length
 822
 823      found_field = .false.
 824      found_att = .false.
 825      do i=1,mpp_file(unit)%nvar
 826         if( trim(mpp_file(unit)%Var(i)%name) == trim(fieldname)) then
 827            found_field = .true.
 828            do j=1, size(mpp_file(unit)%Var(i)%Att(:))
 829               if( trim(mpp_file(unit)%Var(i)%Att(j)%name) == trim(attname) ) then
 830                  found_att = .true.
 831                  length = mpp_file(unit)%Var(i)%Att(j)%len
 832                  if(len(attvalue) .LE. length ) call mpp_error(FATAL, &
 833                       'mpp_io_util.inc: length of attvalue is less than the length of catt')
 834                  attvalue = trim(mpp_file(unit)%Var(i)%Att(j)%catt(1:length))
 835                  exit
 836               end if
 837            end do
 838            exit
 839         end if
 840      end do
 841
 842      if(.NOT. found_field) call mpp_error(FATAL,"mpp_io_util.inc: field "//trim(fieldname)// &
 843                " does not exist in the file "//trim(mpp_file(unit)%name) )
 844      if(.NOT. found_att) call mpp_error(FATAL,"mpp_io_util.inc: attribute "//trim(attname)//" of field "&
 845                //trim(fieldname)// " does not exist in the file "//trim(mpp_file(unit)%name) )
 846
 847      return
 848
 849    end subroutine mpp_get_field_att_text
 850
 851
 852    !####################################################################
 853    !> @brief return mpp_io_nml variable io_clock_on
 854    function mpp_io_clock_on()
 855       logical :: mpp_io_clock_on
 856
 857       mpp_io_clock_on = io_clocks_on
 858       return
 859
 860    end function mpp_io_clock_on
 861
 862
 863    function mpp_attribute_exist(field,name)
 864       logical                      :: mpp_attribute_exist
 865       type(fieldtype),  intent(in) :: field ! The field that you are searching for the attribute.
 866       character(len=*), intent(in) :: name ! name of the attributes
 867
 868       if(field%natt > 0) then
 869          mpp_attribute_exist = ( mpp_find_att(field%Att(:),name) > 0 )
 870       else
 871          mpp_attribute_exist = .false.
 872       endif
 873
 874    end function mpp_attribute_exist
 875
 876 !#######################################################################
 877 subroutine mpp_dist_io_pelist(ssize,pelist)
 878   integer,              intent(in)  :: ssize   ! Stripe size for dist read
 879   integer, allocatable, intent(out) :: pelist(:)
 880   integer :: i, lsize, ioroot
 881   logical :: is_ioroot=.false.
 882
 883   ! Did you make a mistake?
 884   if(ssize < 1) call mpp_error(FATAL,'mpp_dist_io_pelist: I/O stripe size < 1')
 885
 886   is_ioroot = mpp_is_dist_ioroot(ssize,ioroot=ioroot,lsize=lsize)
 887
 888   ! Did I make a mistake?
 889   if(lsize < 1) call mpp_error(FATAL,'mpp_dist_io_pelist: size of pelist < 1')
 890
 891   allocate(pelist(lsize))
 892   do i=1,lsize
 893     pelist(i) = ioroot + i - 1
 894   enddo
 895 end subroutine mpp_dist_io_pelist
 896
 897 !#######################################################################
 898 logical function mpp_is_dist_ioroot(ssize,ioroot,lsize)
 899   integer, intent(in)  :: ssize   ! Dist io set size
 900   integer, intent(out), optional :: ioroot, lsize
 901   integer :: pe, npes, mypos, maxpe, d_ioroot, d_lsize, last_ioroot
 902   integer :: rootpe
 903
 904   if(ssize < 1) call mpp_error(FATAL,'mpp_is_dist_ioroot: I/O stripe size < 1')
 905
 906   mpp_is_dist_ioroot = .false.
 907   rootpe = mpp_root_pe()
 908   d_lsize = ssize
 909   pe = mpp_pe()
 910   mypos = modulo(pe-rootpe,ssize)  ! Which PE am I in the io group?
 911   d_ioroot = pe - mypos            ! What is the io root for the group?
 912   npes = mpp_npes()
 913   maxpe = min(d_ioroot+ssize,npes+rootpe) - 1  ! Handle end case
 914   d_lsize = maxpe - d_ioroot + 1
 915   if(mod(npes,ssize) == 1)then  ! Ensure there are no sets with 1 member
 916     last_ioroot = (npes-1) - ssize
 917     if(pe >= last_ioroot) then
 918        d_ioroot = last_ioroot
 919        d_lsize = ssize + 1
 920     endif
 921   endif
 922   if(pe == d_ioroot) mpp_is_dist_ioroot = .true.
 923   if(PRESENT(ioroot)) ioroot = d_ioroot
 924   if(PRESENT(lsize)) lsize = d_lsize
 925 end function mpp_is_dist_ioroot