other/tign_history_api/tign_history.F90

   1 !*** tign_history ***
   2 !
   3 ! A simple api for reading/writing level function histories in wrfinput files.
   4 !
   5 ! The convention in this code is that two extra variables will be added to
   6 ! the wrfinput_d## files.
   7 !
   8 ! 1. float tign_g(nx,ny) contains the fire history as the time in seconds from
   9 !    the start of the simulation that each point ignited.
  10 !
  11 ! The api takes care of all the particulars of reading/writing... it should
  12 ! be invisible to the caller.
  13 !
  14 ! subroutine:
  15 !   write_tign_history, write (or add) level function history data to
  16 !                      the input file
  17 !   read_tign_history, read all level function history data from the
  18 !                     input file
  19 !   get_grid_info, get any relevant information (such as grid size) from
  20 !                  input file
  21
  22 module tign_history
  23 use netcdf
  24 implicit none
  25 private
  26
  27 ! Various parameters describing the layout of wrf files
  28 integer,parameter::FILELEN=16,XTYPE=nf90_float,INVALID=-9999999,write_time=1
  29 character(len=FILELEN),parameter::FILEFMT="(A10,I02.2)"
  30 character(len=FILELEN),parameter::FILEBSE='wrfinput_d'
  31 character(len=nf90_max_name),parameter::   &
  32               XNAME='west_east_subgrid',   &
  33               YNAME='south_north_subgrid', &
  34               TNAME='Time',                &
  35               XATM='west_east_stag',       &
  36               YATM='south_north_stag',     &
  37               tign_NAME='TIGN_G',          &
  38               DTNAME='DT',                 &
  39               DXNAME='DX',                 &
  40               DYNAME='DY'
  41
  42 public::write_tign_history,read_tign_history,get_grid_info
  43
  44 contains
  45
  46
  47 subroutine write_tign_history(idom,nx,ny,tign)
  48 implicit none
  49
  50 ! The subroutine writes a level function to a wrfinput file in the
  51 ! current directory.
  52 !
  53 ! idom : Input integer describing the the domain number that we will
  54 !        write the array to.  This is only to determine the file name
  55 !        as printf 'wrfinput_d%02i' idom.
  56 !
  57 ! nx,ny : Input size of level function (must be the same size as the fire
  58 !         grid in the input file.  See get_grid_info to get this from the file.
  59 !
  60 ! tign(nx,ny) : The level function history from time(1) to time(ntime).
  61
  62 integer,intent(in)::idom,nx,ny
  63 real,dimension(nx,ny),intent(in)::tign
  64 integer::lx,ly,vhist,ncid,i
  65 call get_grid_info(idom,nx=lx,ny=ly)
  66 if(lx.ne.nx.or.ly.ne.ny)then
  67   print*,'invalid input level function size'
  68   call abort()
  69 endif
  70 call create_hist(idom,vhist)
  71 ncid=open_file(idom,nf90_write)
  72 call check(nf90_put_var(ncid,vhist,tign,start=(/1,1,write_time/),count=(/nx,ny,1/)))
  73 call check(nf90_close(ncid))
  74 end subroutine write_tign_history
  75
  76 subroutine read_tign_history(idom,nx,ny,tign)
  77 implicit none
  78
  79 ! The subroutine reads a level function from a wrfinput file in the
  80 ! current directory.
  81 !
  82 ! idom : Input integer describing the the domain number that we will
  83 !        read the array from.  This is only to determine the file name
  84 !        as printf 'wrfinput_d%02i' idom.
  85 !
  86 ! nx,ny : Input size of level function (must be the same size as the fire
  87 !         grid in the input file.  See get_grid_info to get this from the file.
  88 !
  89 ! tign(nx,ny) : The fire history read from the input file.
  90
  91 integer,intent(in)::idom,nx,ny
  92 real,dimension(nx,ny),intent(out)::tign
  93 integer::lx,ly,vhist,vtime,ncid,ltime
  94 real,dimension(1)::arr
  95 call get_grid_info(idom,nx=lx,ny=ly)
  96 if(lx.ne.nx.or.ly.ne.ny)then
  97   print*,'invalid input array size'
  98   call abort()
  99 endif
 100 call create_hist(idom,vhist)
 101 ncid=open_file(idom,nf90_nowrite)
 102 call check(nf90_get_var(ncid,vhist,tign,start=(/1,1,write_time/),count=(/nx,ny,1/)))
 103 call check(nf90_close(ncid))
 104 end subroutine read_tign_history
 105
 106 subroutine get_grid_info(idom,nx,ny,dx,dy,dt,sr_x,sr_y)
 107 implicit none
 108
 109 ! This subroutine inquires a wrfinput file in the current directory about
 110 ! information relevant to the computation and manipulation of tign history
 111 ! arrays.  The return values are all optional so they can be used all at
 112 ! once or individually.  Calls to this subroutine should be made with the
 113 ! key=value syntax as more arguments may be added in the future to account
 114 ! for additional information required.
 115 !
 116 ! idom : Input integer describing the the domain number that we will
 117 !        read the array from.  This is only to determine the file name
 118 !        as printf 'wrfinput_d%02i' idom.
 119 !
 120 ! Optional outputs:
 121 !
 122 ! nx,ny : The dimensions of the fire grid in the given file.
 123 !
 124 ! dx,dy : The grid resolution of the fire grid in meters.
 125 !
 126 ! dt : The atmospheric time step in seconds.
 127 !
 128 ! sr_x,sr_y : The atmospheric/fire grid refinement factor.
 129
 130 integer,intent(in)::idom
 131 integer,optional,intent(out)::nx,ny,sr_x,sr_y
 132 real,optional,intent(out)::dx,dy,dt
 133
 134 integer::ncid,n,m,lrx,lry
 135
 136 ncid=open_file(idom,nf90_nowrite)
 137
 138 ! get subgrid refinement needed for several of the computations.
 139 call get_dim(ncid,XNAME,n)
 140 call get_dim(ncid,XATM,m)
 141 lrx=n/m
 142 call get_dim(ncid,YNAME,n)
 143 call get_dim(ncid,YATM,m)
 144 lry=n/m
 145
 146 if(present(nx))then
 147   call get_dim(ncid,XNAME,n)
 148   call get_dim(ncid,XATM,m)
 149   nx=n-lrx  ! correct for extra memory in fire grid arrays
 150 endif
 151 if(present(ny))then
 152   call get_dim(ncid,YNAME,n)
 153   call get_dim(ncid,YATM,m)
 154   ny=n-lry  ! correct for extra memory in fire grid arrays
 155 endif
 156 if(present(dx))then
 157   dx=get_attr(ncid,DXNAME)
 158 endif
 159 if(present(dy))then
 160   dy=get_attr(ncid,DYNAME)
 161 endif
 162 if(present(dt))then
 163   dt=get_attr(ncid,DTNAME)
 164 endif
 165 if(present(sr_x))sr_x=lrx
 166 if(present(sr_y))sr_y=lry
 167 call check(nf90_close(ncid))
 168 end subroutine get_grid_info
 169
 170 character(len=FILELEN) function get_file_name(idom)
 171 integer,intent(in)::idom
 172 write(get_file_name,FILEFMT) FILEBSE,idom
 173 end function get_file_name
 174
 175 integer function open_file(idom,rw)
 176 implicit none
 177 integer,intent(in)::idom,rw
 178 call check(nf90_open(trim(get_file_name(idom)),rw,open_file))
 179 end function open_file
 180
 181 subroutine get_dim(ncid,dname,dlen,dimid)
 182 implicit none
 183 integer,intent(in)::ncid
 184 character(len=*),intent(in)::dname
 185 integer,optional,intent(out)::dlen,dimid
 186 integer::did
 187 call check(nf90_inq_dimid(ncid,dname,did))
 188 if(present(dimid))dimid=did
 189 if(present(dlen))then
 190   call check(nf90_inquire_dimension(ncid,did,len=dlen))
 191 endif
 192 end subroutine get_dim
 193
 194 real function get_attr(ncid,aname) result(aval)
 195 implicit none
 196 integer,intent(in)::ncid
 197 character(len=*),intent(in)::aname
 198 call check(nf90_get_att(ncid,nf90_global,aname,aval))
 199 end function get_attr
 200
 201 subroutine create_hist(idom,ihist)
 202 implicit none
 203 integer,intent(in)::idom
 204 integer,intent(out)::ihist
 205 integer::ncid,ix,iy,ih,err
 206 ncid=open_file(idom,nf90_write)
 207 call check(nf90_redef(ncid))
 208 call get_dim(ncid,XNAME,dimid=ix)
 209 call get_dim(ncid,YNAME,dimid=iy)
 210 call get_dim(ncid,TNAME,dimid=ih)
 211 err=nf90_def_var(ncid,tign_NAME,XTYPE,(/ix,iy,ih/),ihist)
 212 call check(nf90_enddef(ncid))
 213 call check(nf90_inq_varid(ncid,tign_NAME,ihist))
 214 call check(nf90_close(ncid))
 215 end subroutine create_hist
 216
 217 subroutine check(ncout)
 218 implicit none
 219 integer,intent(in)::ncout
 220 if(ncout.ne.nf90_noerr)then
 221   print*,'ERROR in netcdf call, ierr=',ncout
 222   print*,nf90_strerror(ncout)
 223   call abort()
 224 endif
 225 end subroutine check
 226
 227 end module tign_history
 228
 229 #ifdef TESTING
 230 program tign_history_main
 231 use tign_history
 232 implicit none
 233 integer,parameter::i=1
 234 integer::nx,ny,ncid,sr_x,sr_y
 235 real::dx,dy,dt
 236 call get_grid_info(i,nx=nx,ny=ny,dx=dx,dy=dy,dt=dt,sr_x=sr_x,sr_y=sr_y)
 237 print*,nx,ny,dx,dy,dt,sr_x,sr_y
 238 end program tign_history_main
 239 #endif