1! This is part of the netCDF package. 2! Copyright 2006 University Corporation for Atmospheric Research/Unidata. 3! See COPYRIGHT file for conditions of use. 4 5! This program tests netCDF-4 variable functions from fortran. 6 7! Ed Hartnett, 2009 8 9program f90tst_vars3 10 use typeSizes 11 use netcdf 12 implicit none 13 14 ! This is the name of the data file we will create. 15 character (len = *), parameter :: FILE_NAME = "f90tst_vars3.nc" 16 17 ! We are writing 2D data, a 6 x 12 grid. 18 integer, parameter :: MAX_DIMS = 2 19 integer, parameter :: NX = 6, NY = 12 20 integer :: data_out(NY, NX), data_in(NY, NX) 21 integer :: data_out_1d(NX), data_in_1d(NX) 22 23 ! We need these ids and other gunk for netcdf. 24 integer :: ncid, varid1, varid2, varid3, varid4, varid5, dimids(MAX_DIMS) 25 integer :: chunksizes(MAX_DIMS), chunksizes_in(MAX_DIMS) 26 integer :: x_dimid, y_dimid 27 integer :: nvars, ngatts, ndims, unlimdimid, file_format 28 integer :: x, y 29 integer, parameter :: DEFAULT_CACHE_NELEMS = 10000, DEFAULT_CACHE_SIZE = 1000000 30 integer, parameter :: DEFAULT_CACHE_PREEMPTION = 22 31 integer, parameter :: DEFLATE_LEVEL = 4 32 integer (kind = EightByteInt), parameter :: TOE_SAN_VALUE = 2147483648_EightByteInt 33 character (len = *), parameter :: VAR1_NAME = "Chon-Ji" 34 character (len = *), parameter :: VAR2_NAME = "Tan-Gun" 35 character (len = *), parameter :: VAR3_NAME = "Toe-San" 36 character (len = *), parameter :: VAR4_NAME = "Won-Hyo" 37 character (len = *), parameter :: VAR5_NAME = "Yul-Guk" 38 integer, parameter :: CACHE_SIZE = 8, CACHE_NELEMS = 571, CACHE_PREEMPTION = 66 39 40 ! Information read back from the file to check correctness. 41 integer :: varid1_in, varid2_in, varid3_in, varid4_in, varid5_in 42 integer :: xtype_in, ndims_in, natts_in, dimids_in(MAX_DIMS) 43 character (len = nf90_max_name) :: name_in 44 integer :: endianness_in, deflate_level_in 45 logical :: shuffle_in, fletcher32_in, contiguous_in 46 integer (kind = EightByteInt) :: toe_san_in 47 integer :: cache_size_in, cache_nelems_in, cache_preemption_in 48 49 print *, '' 50 print *,'*** Testing definition of netCDF-4 vars from Fortran 90.' 51 52 ! Create some pretend data. 53 do x = 1, NX 54 do y = 1, NY 55 data_out(y, x) = (x - 1) * NY + (y - 1) 56 end do 57 end do 58 do x = 1, NX 59 data_out_1d(x) = x 60 end do 61 62 ! Create the netCDF file. 63 call check(nf90_create(FILE_NAME, nf90_netcdf4, ncid, cache_nelems = DEFAULT_CACHE_NELEMS, & 64 cache_size = DEFAULT_CACHE_SIZE, cache_preemption = DEFAULT_CACHE_PREEMPTION)) 65 66 ! Define the dimensions. 67 call check(nf90_def_dim(ncid, "x", NX, x_dimid)) 68 call check(nf90_def_dim(ncid, "y", NY, y_dimid)) 69 dimids = (/ y_dimid, x_dimid /) 70 71 ! Define some variables. 72 chunksizes = (/ NY, NX /) 73 call check(nf90_def_var(ncid, VAR1_NAME, NF90_INT, dimids, varid1, chunksizes = chunksizes, & 74 shuffle = .TRUE., fletcher32 = .TRUE., endianness = nf90_endian_big, deflate_level = DEFLATE_LEVEL)) 75 call check(nf90_def_var(ncid, VAR2_NAME, NF90_INT, dimids, varid2, contiguous = .TRUE.)) 76 call check(nf90_def_var(ncid, VAR3_NAME, NF90_INT64, varid3)) 77 call check(nf90_def_var(ncid, VAR4_NAME, NF90_INT, x_dimid, varid4, contiguous = .TRUE.)) 78 call check(nf90_def_var(ncid, VAR5_NAME, NF90_INT, dimids, varid5, chunksizes = chunksizes, & 79 deflate_level = DEFLATE_LEVEL, cache_size = CACHE_SIZE, cache_nelems = CACHE_NELEMS, & 80 cache_preemption = CACHE_PREEMPTION)) 81 82 call check(nf90_inquire_variable(ncid, varid1, cache_size = cache_size_in, & 83 cache_nelems = cache_nelems_in, cache_preemption = cache_preemption_in)) 84! if (cache_size_in .ne. DEFAULT_CACHE_SIZE .or. cache_nelems_in .ne. & 85! DEFAULT_CACHE_NELEMS .or. cache_preemption_in .ne. DEFAULT_CACHE_PREEMPTION) stop 101 86 call check(nf90_inquire_variable(ncid, varid5, cache_size = cache_size_in, & 87 cache_nelems = cache_nelems_in, cache_preemption = cache_preemption_in)) 88 if (cache_size_in .ne. CACHE_SIZE .or. cache_nelems_in .ne. & 89 CACHE_NELEMS .or. cache_preemption_in .ne. CACHE_PREEMPTION) stop 1 90 91 ! Write the pretend data to the file. 92 call check(nf90_put_var(ncid, varid1, data_out)) 93 call check(nf90_put_var(ncid, varid2, data_out)) 94 call check(nf90_put_var(ncid, varid3, TOE_SAN_VALUE)) 95 call check(nf90_put_var(ncid, varid4, data_out_1d)) 96 call check(nf90_put_var(ncid, varid5, data_out)) 97 98 ! Close the file. 99 call check(nf90_close(ncid)) 100 101 ! Reopen the file. 102 call check(nf90_open(FILE_NAME, nf90_nowrite, ncid)) 103 104 ! Check some stuff out. 105 call check(nf90_inquire(ncid, ndims, nvars, ngatts, unlimdimid, file_format)) 106 if (ndims /= 2 .or. nvars /= 5 .or. ngatts /= 0 .or. unlimdimid /= -1 .or. & 107 file_format /= nf90_format_netcdf4) stop 2 108 109 ! Get varids. 110 call check(nf90_inq_varid(ncid, VAR1_NAME, varid1_in)) 111 call check(nf90_inq_varid(ncid, VAR2_NAME, varid2_in)) 112 call check(nf90_inq_varid(ncid, VAR3_NAME, varid3_in)) 113 call check(nf90_inq_varid(ncid, VAR4_NAME, varid4_in)) 114 call check(nf90_inq_varid(ncid, VAR5_NAME, varid5_in)) 115 116 ! Check variable 1. 117 call check(nf90_inquire_variable(ncid, varid1_in, name_in, xtype_in, ndims_in, dimids_in, & 118 natts_in, chunksizes = chunksizes_in, endianness = endianness_in, fletcher32 = fletcher32_in, & 119 deflate_level = deflate_level_in, shuffle = shuffle_in)) 120 if (name_in .ne. VAR1_NAME .or. xtype_in .ne. NF90_INT .or. ndims_in .ne. MAX_DIMS .or. & 121 natts_in .ne. 0 .or. dimids_in(1) .ne. dimids(1) .or. dimids_in(2) .ne. dimids(2)) stop 3 122 if (chunksizes_in(1) /= chunksizes(1) .or. chunksizes_in(2) /= chunksizes(2)) & 123 stop 4 124 if (endianness_in .ne. nf90_endian_big) stop 5 125 126 ! Check variable 2. 127 call check(nf90_inquire_variable(ncid, varid2_in, name_in, xtype_in, ndims_in, dimids_in, & 128 natts_in, contiguous = contiguous_in, endianness = endianness_in, fletcher32 = fletcher32_in, & 129 deflate_level = deflate_level_in, shuffle = shuffle_in)) 130 if (name_in .ne. VAR2_NAME .or. xtype_in .ne. NF90_INT .or. ndims_in .ne. MAX_DIMS .or. & 131 natts_in .ne. 0 .or. dimids_in(1) .ne. dimids(1) .or. dimids_in(2) .ne. dimids(2)) stop 6 132 if (deflate_level_in .ne. 0 .or. .not. contiguous_in .or. fletcher32_in .or. shuffle_in) stop 7 133 134 ! Check variable 3. 135 call check(nf90_inquire_variable(ncid, varid3_in, name_in, xtype_in, ndims_in, dimids_in, & 136 natts_in, contiguous = contiguous_in, endianness = endianness_in, fletcher32 = fletcher32_in, & 137 deflate_level = deflate_level_in, shuffle = shuffle_in)) 138 if (name_in .ne. VAR3_NAME .or. xtype_in .ne. NF90_INT64 .or. ndims_in .ne. 0 .or. & 139 natts_in .ne. 0) stop 8 140 if (deflate_level_in .ne. 0 .or. .not. contiguous_in .or. fletcher32_in .or. shuffle_in) stop 9 141 142 ! Check variable 4. 143 call check(nf90_inquire_variable(ncid, varid4_in, name_in, xtype_in, ndims_in, dimids_in, & 144 natts_in, contiguous = contiguous_in, endianness = endianness_in, fletcher32 = fletcher32_in, & 145 deflate_level = deflate_level_in, shuffle = shuffle_in)) 146 if (name_in .ne. VAR4_NAME .or. xtype_in .ne. NF90_INT .or. ndims_in .ne. 1 .or. & 147 natts_in .ne. 0 .or. dimids_in(1) .ne. x_dimid) stop 10 148 if (deflate_level_in .ne. 0 .or. .not. contiguous_in .or. fletcher32_in .or. shuffle_in) stop 11 149 150 ! Check the data. 151 call check(nf90_get_var(ncid, varid1_in, data_in)) 152 do x = 1, NX 153 do y = 1, NY 154 if (data_out(y, x) .ne. data_in(y, x)) stop 12 155 end do 156 end do 157 call check(nf90_get_var(ncid, varid2_in, data_in)) 158 do x = 1, NX 159 do y = 1, NY 160 if (data_out(y, x) .ne. data_in(y, x)) stop 13 161 end do 162 end do 163 call check(nf90_get_var(ncid, varid3_in, toe_san_in)) 164 if (toe_san_in .ne. TOE_SAN_VALUE) stop 14 165 call check(nf90_get_var(ncid, varid4_in, data_in_1d)) 166 do x = 1, NX 167 if (data_out_1d(x) .ne. data_in_1d(x)) stop 15 168 end do 169 call check(nf90_get_var(ncid, varid5_in, data_in)) 170 do x = 1, NX 171 do y = 1, NY 172 if (data_out(y, x) .ne. data_in(y, x)) stop 12 173 end do 174 end do 175 176 ! Close the file. 177 call check(nf90_close(ncid)) 178 179 print *,'*** SUCCESS!' 180 181contains 182! This subroutine handles errors by printing an error message and 183! exiting with a non-zero status. 184 subroutine check(errcode) 185 use netcdf 186 implicit none 187 integer, intent(in) :: errcode 188 189 if(errcode /= nf90_noerr) then 190 print *, 'Error: ', trim(nf90_strerror(errcode)) 191 stop 2 192 endif 193 end subroutine check 194end program f90tst_vars3 195 196