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 integer function l2i(l)
 174     implicit none
 175     logical, intent(in)::l
 176     if(l)then
 177         l2i = 1
 178     else
 179         l2i = 0
 180     endif
 181 end function l2i
 182
 183
 184 subroutine netcdf_var_info(ncid,varname,dims,varid,prints)
 185     implicit none
 186 !*** arguments
 187     integer, intent(in)::ncid
 188     character(len=*)::varname
 189     integer,intent(out)::dims(:),varid
 190     integer,intent(in),optional::prints
 191 !*** local
 192     integer, parameter::mdims = 256
 193     integer:: xtype, ndims, natts, dimids(mdims),i,j,attnum
 194     integer :: len
 195     character(len=nf90_max_name):: name
 196     integer:: values_int(mdims)
 197     real:: values_real(mdims)
 198     character(len=mdims):: values_char
 199     character(LEN=256):: filename, msg
 200     logical::verbose=.true.
 201
 202     if(present(prints)) verbose = prints>0
 203
 204     call check(nf90_inq_varid(ncid,trim(varname),varid),"nf90_inq_varid"//trim(varname))
 205     call check(nf90_inquire_variable(ncid, varid, name, xtype, ndims, dimids, natts),"nf90_inquire_variable")
 206     if(ndims>mdims)call crash("netcdf_var_info: increase mdims")
 207     if(ndims>size(dims))call crash("netcdf_var_info: dims too short")
 208     if(verbose)then
 209         write(msg,*)"variable ",trim(name), " xtype",xtype, "ndims",ndims, "natts",natts
 210         call message(msg)
 211     endif
 212     do i=1,ndims
 213         call check(nf90_inquire_dimension(ncid, dimids(i), name, len),"nf90_inquire_dimension")
 214         dims(i)=len
 215         if(verbose)then
 216             write(msg,*)"dimension ",i,trim(name)," length",len
 217             call message(msg)
 218         endif
 219     enddo
 220     if(.not.verbose)return
 221     do i=1,natts
 222         attnum = i
 223         call check(nf90_inq_attname(ncid, varid, attnum, name),"nf90_inq_attname")
 224         call check(nf90_inquire_attribute(ncid, varid, trim(name), xtype, len, attnum),"nf90_inquire_attribute")
 225         if(len>mdims)call crash("netcdf_var_info: increase mdims")
 226         !write(msg,*)"attribute ",i,trim(name),' type',xtype
 227         !call message(msg)
 228         select case (xtype)
 229             case (nf90_char)
 230                 call check(nf90_get_att(ncid, varid, trim(name), values_char),"nf90_get_att")
 231                 write(msg,*)"attribute ",i,trim(name)," type ",xtype," values",len," : ",trim(values_char)
 232             case (nf90_int,nf90_short,nf90_ushort,nf90_uint,nf90_int64,nf90_uint64)
 233                 call check(nf90_get_att(ncid, varid, trim(name), values_int),"nf90_get_att")
 234                 write(msg,*)"attribute ",i,trim(name)," type ",xtype," values",len," : ",(values_int(j),j=1,len)
 235             case (nf90_float,nf90_double)
 236                 call check(nf90_get_att(ncid, varid, trim(name), values_real),"nf90_get_att")
 237                 write(msg,*)"attribute ",i,trim(name)," type ",xtype," values",len," : ",(values_real(j),j=1,len)
 238             case default
 239                 write(msg,*)'attribute type ',xtype,' not supported'
 240         end select
 241         call message(msg)
 242     enddo
 243 end subroutine netcdf_var_info
 244
 245 subroutine check(ierr,errmsg)
 246     implicit none
 247     integer, intent(in)::ierr
 248     character(len=*), intent(in)::errmsg
 249     character(len=256)msg
 250     if(ierr.ne.0)then
 251         write(msg,"(a,a,i6,1x,a)")trim(errmsg)," error",ierr,trim(nf90_strerror(ierr))
 252         call crash(trim(msg))
 253     endif
 254 end subroutine check
 255
 256 end module module_netcdf