femwind/fortran/module_netcdf.f90

   1 module module_netcdf
   2
   3 use netcdf
   4 use module_utils
   5
   6 integer::netcdf_msglevel = 0
   7
   8 contains
   9
  10 ! from https://github.com/openwfm/wrf-fire/blob/master/standalone/wrf_netcdf.F
  11 subroutine ncopen(filename,mode,ncid)
  12 !*** purpose: open netcdf file wrapper with an informative error message
  13 implicit none
  14 !*** arguments
  15 character(len=*), intent(in):: filename
  16 integer, intent(in)::mode
  17 integer, intent(out):: ncid
  18 !*** executable
  19 call check(nf90_open(trim(filename),mode,ncid),"Cannot open file "//trim(filename))
  20 print *,"Opened netcdf file ",trim(filename)," as ",ncid," mode ",mode
  21 end subroutine ncopen
  22
  23 subroutine ncclose(ncid)
  24 !*** purpose: open netcdf file wrapper with informative error message
  25 implicit none
  26 !*** arguments
  27 integer, intent(in):: ncid
  28 print *,"Closing netcdf file ",ncid
  29 call check(nf90_close(ncid),"Cannot close netcdf file ")
  30 end subroutine ncclose
  31
  32 real function netcdf_read_att(ncid,name)
  33 ! read real global scalar attribute
  34 !*** arguments
  35 integer, intent(in)::ncid
  36 character(len=*), intent(in)::name
  37 !*** local
  38 real::value
  39 integer::xtype,len,attnum
  40 character(len=256)::msg
  41 !*** executable
  42 call check(nf90_inquire_attribute(ncid, nf90_global, trim(name), xtype, len, attnum),"nf90_inquire_attribute")
  43 if(xtype.ne.nf90_float.or.len.ne.1)then
  44    write(msg,*)"netcdf_read_att can read only float scalar but ",trim(name)," has xtype=",xtype," len=",len
  45    call crash(trim(msg))
  46 endif
  47 call check(nf90_get_att(ncid,nf90_global, trim(name), value),"nf90_get_att")
  48 write(msg,*)"netcdf_read_att returning ",trim(name),"=",value
  49 call message(msg)
  50 netcdf_read_att = value
  51 end function netcdf_read_att
  52
  53 integer function netcdf_read_int_wrf(ncid,name,istep)
  54 implicit none
  55 !*** Read one integer
  56 !*** arguments
  57     integer, intent(in)::ncid                      ! open netcdf file
  58     character(LEN=*),intent(in)::name              ! variable name
  59     integer, intent(in)::istep                     ! index in unlimited dimension (timestep number)
  60 !*** local
  61     integer::ia(1)                     ! variable to store
  62     integer::ierr,varid
  63 !*** executable
  64         print *,"netcdf_read_int_wrf reading variable ",trim(name)," time step ",istep
  65         call check(nf90_inq_varid(ncid, trim(name), varid), &
  66             "netcdf_read_int_wrf/nf90_inq_varid:"//trim(name))
  67         call check(nf90_get_var(ncid, varid, ia, start = (/istep/), count = (/1/)), &
  68              "netcdf_read_int_wrf/nf90_get_var:"//trim(name))
  69         print *,"netcdf_read_int_wrf:", trim(name), " = ",ia
  70         netcdf_read_int_wrf = ia(1)
  71 end function netcdf_read_int_wrf
  72
  73 subroutine netcdf_write_int(ncid,ia,varname)
  74 implicit none
  75 !*** arguments
  76     integer, intent(in)::                         &
  77     ncid,                                         & ! open netcdf file
  78     ia                                              ! variable to write
  79     character(LEN=*),intent(in):: varname
  80 !*** local
  81     integer::varid,ival
  82     character(len=256)::msg
  83 !*** executable
  84     write(msg,*)'netcdf_write_int: varname=',varname,' value=',ia
  85     call message(msg)
  86         call check(nf90_inq_varid(ncid, trim(varname), varid), &
  87             "netcdf_write_int/nf90_inq_varid:"//trim(varname))
  88         ival = ia
  89         call check(nf90_put_var(ncid, varid, ival), &
  90             "netcdf_write_int/nf90_put_var:"//trim(varname))
  91 end subroutine netcdf_write_int
  92
  93 subroutine netcdf_write_array(ncid,a,name)
  94     implicit none
  95
  96 !*** arguments
  97     integer, intent(in)::ncid                    ! open netcdf file
  98     real,intent(in),dimension(:,:,:)::a
  99     character(LEN=*),intent(in):: name
 100
 101 !*** local
 102     integer,dimension(4)::star,cnts
 103     integer::i,j,k,varid,ends(4),dims(4),n(3)
 104     real,dimension(:,:,:,:),allocatable::at
 105     character(len=256) msg
 106
 107     ! get idx
 108     n=shape(a)
 109     call netcdf_var_info(ncid,name,dims,varid,netcdf_msglevel)
 110     star   = (/1,1,1,1/)
 111     ends   = (/dims(1),dims(2),dims(3),1/)
 112     ends   = min(ends,dims)
 113     cnts = ends - star + 1
 114
 115     ! transpose a -> at
 116     allocate(at(star(1):ends(1),star(2):ends(2),star(3):ends(3),1))
 117     do k=star(3),ends(3)
 118         do j=star(2),ends(2)
 119             do i=star(1),ends(1)
 120                 at(i,j,k,1)=a(i,k,j)
 121             enddo
 122         enddo
 123     enddo
 124
 125     if(netcdf_msglevel>=0) &
 126       write(msg,*)"writing ",trim(name),n(1),star(1),ends(1),n(2),star(2),ends(2),n(3),star(3),ends(3)
 127     call message(msg)
 128
 129     ! write to file
 130     call check(nf90_put_var(ncid, varid, at, start = star, count = cnts),"nf90_put_var:"//trim(name))
 131
 132     deallocate(at)
 133
 134 end subroutine netcdf_write_array
 135
 136
 137 subroutine netcdf_write_2d(ncid,a,name,iframe)
 138     implicit none
 139 !*** purpose
 140 !   write a 2d array to netcdf file
 141
 142 !*** arguments
 143     integer, intent(in)::ncid                    ! open netcdf file
 144     real,intent(in),dimension(:,:)::a
 145     character(LEN=*),intent(in):: name
 146     integer, intent(in)::iframe                  ! time frame to write in
 147
 148 !*** local
 149     integer,dimension(3)::star,cnts
 150     integer::i,j,k,varid,ends(3),dims(3),n(2)
 151     character(len=256) msg
 152
 153     ! get idx
 154     n=shape(a)
 155     call netcdf_var_info(ncid,name,dims,varid,netcdf_msglevel)
 156     write(msg,*)"array ",trim(name)," shape ",n," NetCDF dimensions ",dims
 157     call message(msg)
 158
 159     if(dims(1).lt.n(1).or.dims(2).lt.n(2))call crash("array shape too large")
 160     star   = (/1,1,iframe/)
 161     ends   = (/n(1),n(2),iframe/)
 162     if(iframe.gt.dims(3))call crash('netcdf_write_2d: frame not in file')
 163     cnts = ends - star + 1
 164
 165     write(msg,*)"writing ",trim(name)," from ",star," to ",ends
 166     call message(msg)
 167
 168     ! write to file
 169     call check(nf90_put_var(ncid, varid, a, start = star, count = cnts),"nf90_put_var:"//trim(name))
 170
 171 end subroutine netcdf_write_2d
 172
 173
 174 subroutine netcdf_write_wrf_3d(ncid,a,name,iframe)
 175     implicit none
 176 !*** purpose
 177 !   write a 3d array to netcdf file with WRF ordering
 178
 179 !*** arguments
 180     integer, intent(in)::ncid                    ! open netcdf file
 181     real,intent(in),dimension(:,:,:)::a
 182     character(LEN=*),intent(in):: name
 183     integer, intent(in)::iframe                  ! time frame to write in
 184
 185 !*** local
 186     integer:: i,j,k,n(3)
 187     real, allocatable:: aw(:,:,:)
 188 !*** executable
 189     n = shape(a)
 190     allocate(aw(n(1),n(3),n(2)))
 191     do i = 1,n(1)
 192        do k = 1,n(3)
 193           do j = 1,n(2)
 194              aw(i,k,j) = a(i,j,k)
 195           enddo
 196        enddo
 197     enddo
 198     call netcdf_write_3d(ncid,aw,name,iframe)
 199     deallocate(aw)
 200 end subroutine netcdf_write_wrf_3d
 201
 202 subroutine netcdf_write_3d(ncid,a,name,iframe)
 203     implicit none
 204 !*** purpose
 205 !   write a 3d array to netcdf file
 206
 207 !*** arguments
 208     integer, intent(in)::ncid                    ! open netcdf file
 209     real,intent(in),dimension(:,:,:)::a
 210     character(LEN=*),intent(in):: name
 211     integer, intent(in)::iframe                  ! time frame to write in
 212
 213 !*** local
 214     integer,dimension(4)::star,cnts,ends,dims
 215     integer::i,j,k,varid,n(3)
 216     character(len=256) msg
 217
 218     ! get idx
 219     n=shape(a)
 220     call netcdf_var_info(ncid,name,dims,varid,netcdf_msglevel)
 221     write(msg,*)"array ",trim(name)," shape ",n," NetCDF dimensions ",dims
 222     call message(msg)
 223
 224     if(dims(1).lt.n(1).or.dims(2).lt.n(2).or.dims(3).lt.n(3))call crash("array shape too large")
 225     star   = (/1,1,1,iframe/)
 226     ends   = (/n(1),n(2),n(3),iframe/)
 227     if(iframe.gt.dims(4))call crash('netcdf_write_2d: frame not in file')
 228     cnts = ends - star + 1
 229
 230     write(msg,*)"writing ",trim(name)," from ",star," to ",ends
 231     call message(msg)
 232
 233     ! write to file
 234     call check(nf90_put_var(ncid, varid, a, start = star, count = cnts),"nf90_put_var:"//trim(name))
 235
 236 end subroutine netcdf_write_3d
 237
 238 integer function l2i(l)
 239     implicit none
 240     logical, intent(in)::l
 241     if(l)then
 242         l2i = 1
 243     else
 244         l2i = 0
 245     endif
 246 end function l2i
 247
 248
 249 subroutine netcdf_var_info(ncid,varname,dims,varid,prints)
 250     implicit none
 251 !*** arguments
 252     integer, intent(in)::ncid
 253     character(len=*)::varname
 254     integer,intent(out)::dims(:),varid
 255     integer,intent(in),optional::prints
 256 !*** local
 257     integer, parameter::mdims = 256
 258     integer:: xtype, ndims, natts, dimids(mdims),i,j,attnum
 259     integer :: len
 260     character(len=nf90_max_name):: name
 261     integer:: values_int(mdims)
 262     real:: values_real(mdims)
 263     character(len=mdims):: values_char
 264     character(LEN=256):: filename, msg
 265     logical::verbose=.true.
 266
 267     if(present(prints)) verbose = prints>0
 268
 269     call check(nf90_inq_varid(ncid,trim(varname),varid),"nf90_inq_varid"//trim(varname))
 270     call check(nf90_inquire_variable(ncid, varid, name, xtype, ndims, dimids, natts),"nf90_inquire_variable")
 271     if(ndims>mdims)call crash("netcdf_var_info: increase mdims")
 272     if(ndims>size(dims))call crash("netcdf_var_info: dims too short")
 273     if(verbose)then
 274         write(msg,*)"variable ",trim(name), " xtype",xtype, "ndims",ndims, "natts",natts
 275         call message(msg)
 276     endif
 277     do i=1,ndims
 278         call check(nf90_inquire_dimension(ncid, dimids(i), name, len),"nf90_inquire_dimension")
 279         dims(i)=len
 280         if(verbose)then
 281             write(msg,*)"dimension ",i,trim(name)," length",len
 282             call message(msg)
 283         endif
 284     enddo
 285     if(.not.verbose)return
 286     do i=1,natts
 287         attnum = i
 288         call check(nf90_inq_attname(ncid, varid, attnum, name),"nf90_inq_attname")
 289         call check(nf90_inquire_attribute(ncid, varid, trim(name), xtype, len, attnum),"nf90_inquire_attribute")
 290         if(len>mdims)call crash("netcdf_var_info: increase mdims")
 291         !write(msg,*)"attribute ",i,trim(name),' type',xtype
 292         !call message(msg)
 293         select case (xtype)
 294             case (nf90_char)
 295                 call check(nf90_get_att(ncid, varid, trim(name), values_char),"nf90_get_att")
 296                 write(msg,*)"attribute ",i,trim(name)," type ",xtype," values",len," : ",trim(values_char)
 297             case (nf90_int,nf90_short,nf90_ushort,nf90_uint,nf90_int64,nf90_uint64)
 298                 call check(nf90_get_att(ncid, varid, trim(name), values_int),"nf90_get_att")
 299                 write(msg,*)"attribute ",i,trim(name)," type ",xtype," values",len," : ",(values_int(j),j=1,len)
 300             case (nf90_float,nf90_double)
 301                 call check(nf90_get_att(ncid, varid, trim(name), values_real),"nf90_get_att")
 302                 write(msg,*)"attribute ",i,trim(name)," type ",xtype," values",len," : ",(values_real(j),j=1,len)
 303             case default
 304                 write(msg,*)'attribute type ',xtype,' not supported'
 305         end select
 306         call message(msg)
 307     enddo
 308 end subroutine netcdf_var_info
 309
 310 subroutine check(ierr,errmsg)
 311     implicit none
 312     integer, intent(in)::ierr
 313     character(len=*), intent(in)::errmsg
 314     character(len=256)msg
 315     if(ierr.ne.0)then
 316         write(msg,"(a,a,i6,1x,a)")trim(errmsg)," error",ierr,trim(nf90_strerror(ierr))
 317         call crash(trim(msg))
 318     endif
 319 end subroutine check
 320
 321 end module module_netcdf