[Dart-dev] [3811] DART/trunk/models/wrf: First of several expected major updates to the WRF model_mod and
nancy at ucar.edu
nancy at ucar.edu
Wed Apr 15 14:10:13 MDT 2009
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/dart-dev/attachments/20090415/65c1cdc2/attachment-0001.html
-------------- next part --------------
Added: DART/trunk/models/wrf/WRF_DART_utilities/dart_to_wrf.f90
===================================================================
--- DART/trunk/models/wrf/WRF_DART_utilities/dart_to_wrf.f90 (rev 0)
+++ DART/trunk/models/wrf/WRF_DART_utilities/dart_to_wrf.f90 2009-04-15 20:10:12 UTC (rev 3811)
@@ -0,0 +1,354 @@
+! Data Assimilation Research Testbed -- DART
+! Copyright 2004-2007, Data Assimilation Research Section
+! University Corporation for Atmospheric Research
+! Licensed under the GPL -- www.gpl.org/licenses/gpl.html
+
+PROGRAM dart_to_wrf
+
+! <next few lines under version control, do not edit>
+! $URL$
+! $Id$
+! $Revision$
+! $Date$
+
+use types_mod, only : r8
+use time_manager_mod, only : time_type, write_time, read_time, get_date, set_date, operator(-), &
+ get_time, print_time, set_calendar_type, GREGORIAN, julian_day
+use utilities_mod, only : get_unit, file_exist, open_file, close_file, &
+ error_handler, E_ERR, E_MSG, initialize_utilities, &
+ register_module, logfileunit, nmlfileunit, timestamp, &
+ find_namelist_in_file, check_namelist_read, &
+ nc_check
+use assim_model_mod, only : open_restart_read, open_restart_write, aread_state_restart, &
+ awrite_state_restart
+use model_mod, only : max_state_variables, &
+ num_state_table_columns, read_wrf_dimensions, &
+ get_number_of_wrf_variables, &
+ get_variable_size_from_file, wrf_dom, &
+ fill_default_state_table, &
+ trans_1Dto3D, trans_1Dto2D,&
+ set_wrf_date
+
+use netcdf
+
+implicit none
+
+! version controlled file description for error handling, do not edit
+character(len=128), parameter :: &
+ source = "$URL$", &
+ revision = "$Revision$", &
+ revdate = "$Date$"
+
+
+!-----------------------------------------------------------------------
+! Model namelist parameters with default values.
+!-----------------------------------------------------------------------
+
+logical :: output_state_vector = .false. ! state vs. prognostic format
+logical :: default_state_variables = .true. ! use default state list?
+character(len=129) :: wrf_state_variables(num_state_table_columns,max_state_variables) = 'NULL'
+integer :: num_moist_vars = 3
+integer :: num_domains = 1
+integer :: calendar_type = GREGORIAN
+integer :: assimilation_period_seconds = 21600
+logical :: surf_obs = .true.
+logical :: soil_data = .true.
+logical :: h_diab = .false.
+logical :: allow_obs_below_vol = .false.
+character(len = 72) :: adv_mod_command = './wrf.exe'
+real (kind=r8) :: center_search_half_length = 500000.0_r8
+integer :: center_spline_grid_scale = 10
+integer :: vert_localization_coord = 3 ! 1,2,3 == level,pressure,height
+! candidates for including in the WRF netcdf files:
+logical :: polar = .false. ! wrap over the poles
+logical :: periodic_x = .false. ! wrap in longitude or x
+logical :: periodic_y = .false. ! used for single column model, wrap in y
+!JPH -- single column model flag
+logical :: scm = .false. ! using the single column model
+
+
+namelist /model_nml/ output_state_vector, num_moist_vars, &
+ num_domains, calendar_type, surf_obs, soil_data, h_diab, &
+ default_state_variables, wrf_state_variables, &
+ adv_mod_command, assimilation_period_seconds, &
+ allow_obs_below_vol, vert_localization_coord, &
+ center_search_half_length, center_spline_grid_scale, &
+ polar, periodic_x, periodic_y, scm
+
+
+!-------------------------------------------------------------
+
+type(wrf_dom) :: wrf
+
+real(r8), pointer :: dart(:)
+real(r8), pointer :: wrf_var_3d(:,:,:), wrf_var_2d(:,:)
+type(time_type) :: dart_time(2)
+integer :: number_dart_values, ndays, &
+ year, month, day, hour, minute, second
+integer :: ndims, idims(2), dimids(2)
+integer :: i, ivtype, ind, dart_ind, my_index
+character(len=80) :: varname
+character(len=19) :: timestring
+character(len=1) :: idom
+
+logical, parameter :: debug = .true.
+integer :: ncid(50), io, var_id, id, iunit, dart_unit
+integer :: var_element_list(max_state_variables)
+
+write(*,*) 'DART to WRF' !(.false./F)?'
+
+call initialize_utilities('dart_to_wrf')
+call register_module(source, revision, revdate)
+
+! Begin by reading the namelist input
+call find_namelist_in_file("input.nml", "model_nml", iunit)
+read(iunit, nml = model_nml, iostat = io)
+call check_namelist_read(iunit, io, "model_nml")
+
+! Record the namelist values used for the run ...
+write(nmlfileunit, nml=model_nml)
+write( * , nml=model_nml)
+
+call set_calendar_type(calendar_type)
+
+call error_handler(E_MSG,'dart_to_wrf', &
+ 'Converting a dart_state_vector to a WRF netcdf file', &
+ source, revision, revdate)
+
+allocate(wrf%dom(num_domains))
+
+! get default state variable table if asked
+if ( default_state_variables ) then
+ wrf_state_variables = 'NULL'
+ call fill_default_state_table(wrf_state_variables)
+endif
+
+if ( debug ) then
+ print*,'WRF state vector table'
+ print*,'default_state_variables = ',default_state_variables
+ print*,wrf_state_variables
+endif
+
+
+! open wrf data netCDF file 'wrfinput_d0x'
+! we get sizes of the WRF geometry and resolution
+
+! approach is to do everything in a loop, without any regard to the actual
+! order. just use the namelist to read and unroll one at a time.
+
+! big loop over domains, just like in static_init
+number_dart_values = 0
+WRFDomains : do id = 1,num_domains
+
+ write(idom,'(I1)') id
+
+ if(file_exist('wrfinput_d0'//idom)) then
+
+ call nc_check( nf90_open('wrfinput_d0'//idom, NF90_WRITE, ncid(id)), &
+ 'dart_to_wrf','open wrfinput_d0'//idom )
+
+ else
+
+ call error_handler(E_ERR,'dart_to_wrf', &
+ 'Please put wrfinput_d0'//idom//' in the work directory.', source, revision,revdate)
+
+ endif
+
+! read WRF dimensions
+ call read_wrf_dimensions(ncid(id),wrf%dom(id)%bt, wrf%dom(id)%bts, &
+ wrf%dom(id)%sn, wrf%dom(id)%sns, &
+ wrf%dom(id)%we, wrf%dom(id)%wes, &
+ wrf%dom(id)%sls)
+
+
+! get the number of wrf variables wanted in this domain's state
+ wrf%dom(id)%number_of_wrf_variables = get_number_of_wrf_variables(id,wrf_state_variables,var_element_list)
+
+ if ( debug ) then
+ print*,'Domain ',id,' number of wrf variables is ',wrf%dom(id)%number_of_wrf_variables
+ endif
+
+! allocate and store the table locations of the variables valid on this domain
+ allocate(wrf%dom(id)%var_index_list(wrf%dom(id)%number_of_wrf_variables))
+ wrf%dom(id)%var_index_list = var_element_list(1:wrf%dom(id)%number_of_wrf_variables)
+
+! allocate var size
+ allocate(wrf%dom(id)%var_size(3,wrf%dom(id)%number_of_wrf_variables))
+
+! allocate stagger
+ allocate(wrf%dom(id)%stagger(wrf%dom(id)%number_of_wrf_variables))
+
+! add to the dart state vector size
+ do ind = 1,wrf%dom(id)%number_of_wrf_variables
+
+ ! actual location in state variable table
+ my_index = wrf%dom(id)%var_index_list(ind)
+
+ ! get stagger and variable size
+ call get_variable_size_from_file(ncid(id),id, &
+ wrf_state_variables(1,my_index), &
+ wrf%dom(id)%bt, wrf%dom(id)%bts, &
+ wrf%dom(id)%sn, wrf%dom(id)%sns, &
+ wrf%dom(id)%we, wrf%dom(id)%wes, &
+ wrf%dom(id)%stagger(ind), &
+ wrf%dom(id)%var_size(:,ind))
+
+ if ( debug ) then
+ print*,'variable size ',trim(wrf_state_variables(1,my_index)),' ',wrf%dom(id)%var_size(:,ind)
+ endif
+
+ number_dart_values = number_dart_values &
+ + (wrf%dom(id)%var_size(1,ind) &
+ * wrf%dom(id)%var_size(2,ind) &
+ * wrf%dom(id)%var_size(3,ind))
+ enddo
+
+enddo WRFDomains
+
+if ( debug ) then
+ print*,'dart vector size ',number_dart_values
+endif
+
+! allocate dart state vector
+allocate(dart(number_dart_values))
+
+! open dart file
+dart_unit = open_restart_read("dart_wrf_vector")
+
+! read dart vector
+call aread_state_restart(dart_time(2), dart, dart_unit, dart_time(1))
+
+! record wrf.info
+iunit = get_unit()
+open(unit = iunit, file = 'wrf.info')
+call write_time(iunit, dart_time(1))
+call write_time(iunit, dart_time(2))
+call get_date(dart_time(2), year, month, day, hour, minute, second)
+write (iunit,FMT='(I4,5I3.2)') year, month, day, hour, minute, second
+
+write (iunit,*) num_domains
+write (iunit,*) adv_mod_command
+close(iunit)
+
+! set new times
+call set_wrf_date(timestring, year, month, day, hour, minute, second)
+ndays = julian_day(year, month, day)
+
+! loop through domains again and pull each variable from the state
+dart_ind = 1
+WRFDomains2 : do id = 1,num_domains
+
+ do ind = 1,wrf%dom(id)%number_of_wrf_variables
+
+ ! actual location in state variable table
+ my_index = wrf%dom(id)%var_index_list(ind)
+
+ if ( debug ) then
+ write(*,*),'Rolling up variable ',trim(wrf_state_variables(1,my_index))
+ endif
+
+ ! get stagger and variable size
+ call nc_check( nf90_inq_varid(ncid(id),wrf_state_variables(1,my_index), &
+ var_id), 'dart_to_wrf', &
+ 'inq_var_id '//wrf_state_variables(1,my_index))
+
+ if ( wrf%dom(id)%var_size(3,ind) == 1 ) then
+
+ if ( debug ) then
+ write(*,*),trim(wrf_state_variables(1,my_index)),' is 2D'
+ write(*,*),'size ',wrf%dom(id)%var_size(:,ind)
+ endif
+
+ allocate(wrf_var_2d(wrf%dom(id)%var_size(1,ind),wrf%dom(id)%var_size(2,ind)))
+
+ wrf_var_2d = 0.0_r8
+
+ call trans_1Dto2D( dart(dart_ind:),wrf_var_2d, &
+ wrf%dom(id)%var_size(1,ind),wrf%dom(id)%var_size(2,ind))
+
+ call nc_check( nf90_put_var(ncid(id), var_id, wrf_var_2d), &
+ 'dart_to_wrf','put_var '//wrf_state_variables(1,my_index) )
+
+
+ deallocate(wrf_var_2d)
+
+ else
+
+ if ( debug ) then
+ write(*,*),trim(wrf_state_variables(1,my_index)),' is 3D'
+ write(*,*),'size ',wrf%dom(id)%var_size(:,ind)
+ endif
+
+ allocate(wrf_var_3d(wrf%dom(id)%var_size(1,ind),wrf%dom(id)%var_size(2,ind),wrf%dom(id)%var_size(3,ind)))
+
+ wrf_var_3d = 0.0_r8
+
+ call trans_1Dto3D( dart(dart_ind:),wrf_var_3d, &
+ wrf%dom(id)%var_size(1,ind), &
+ wrf%dom(id)%var_size(2,ind),wrf%dom(id)%var_size(3,ind))
+
+ ! better way to set 0, or even necessary with WRF functionality?
+! trim(wrf_state_variables(1,my_index)) == 'QICE' .or. &
+! trim(wrf_state_variables(1,my_index)) == 'QSNOW' .or. &
+! trim(wrf_state_variables(1,my_index)) == 'QGRAUP' .or. &
+! trim(wrf_state_variables(1,my_index)) == 'QNDROP' .or. &
+! trim(wrf_state_variables(1,my_index)) == 'QNGRAUPEL' .or. &
+! trim(wrf_state_variables(1,my_index)) == 'QNSNOW' .or. &
+! trim(wrf_state_variables(1,my_index)) == 'QNRAIN' .or. &
+! trim(wrf_state_variables(1,my_index)) == 'QNICE' .or. &
+! trim(wrf_state_variables(1,my_index)) == 'QNICE' ) then
+ if ( trim(wrf_state_variables(1,my_index)) == 'QVAPOR' .or. &
+ trim(wrf_state_variables(1,my_index)) == 'QRAIN' .or. &
+ trim(wrf_state_variables(1,my_index)) == 'QCLOUD' ) then
+
+ if ( debug ) then
+ write(*,*),'Setting 0 lower bound on ', &
+ trim(wrf_state_variables(1,my_index))
+ endif
+
+ wrf_var_3d = max(0.0_r8,wrf_var_3d)
+
+ endif
+
+ call nc_check( nf90_put_var(ncid(id), var_id, wrf_var_3d), &
+ 'dart_to_wrf','put_var '//wrf_state_variables(1,my_index) )
+
+ deallocate(wrf_var_3d)
+
+ endif
+
+ dart_ind = dart_ind &
+ + (wrf%dom(id)%var_size(1,ind) &
+ * wrf%dom(id)%var_size(2,ind) &
+ * wrf%dom(id)%var_size(3,ind))
+
+ enddo
+
+ ! output times too
+ call nc_check( nf90_inq_varid(ncid(id), "Times", var_id), 'dart_to_wrf', &
+ 'inq_varid Times' )
+ call nc_check( nf90_put_var(ncid(id), var_id, timestring), 'dart_to_wrf', &
+ 'put_var Times' )
+ call nc_check( nf90_put_att(ncid(id), nf90_global, "START_DATE", timestring), &
+ 'dart_to_wrf','put_att START_DATE' )
+ call nc_check( nf90_put_att(ncid(id), nf90_global, "JULYR", year), &
+ 'dart_to_wrf','put_att JULYR' )
+ call nc_check( nf90_put_att(ncid(id), nf90_global, "JULDAY", ndays), &
+ 'dart_to_wrf','put_att JULDY' )
+
+ ! at this point we are done with this domain
+ call nc_check ( nf90_sync(ncid(id)),'dart_to_wrf','sync wrfinput' )
+ call nc_check ( nf90_close(ncid(id)),'dart_to_wrf','close wrfinput' )
+
+enddo WRFDomains2
+
+deallocate(wrf%dom)
+deallocate(dart)
+
+write(logfileunit,*)'FINISHED dart_to_wrf.'
+write(logfileunit,*)
+
+call timestamp(source,revision,revdate,'end') ! That closes the log file, too.
+
+end PROGRAM dart_to_wrf
+
Property changes on: DART/trunk/models/wrf/WRF_DART_utilities/dart_to_wrf.f90
___________________________________________________________________
Name: svn:keywords
+ Date Revision Author HeadURL Id
Added: DART/trunk/models/wrf/WRF_DART_utilities/wrf_to_dart.f90
===================================================================
--- DART/trunk/models/wrf/WRF_DART_utilities/wrf_to_dart.f90 (rev 0)
+++ DART/trunk/models/wrf/WRF_DART_utilities/wrf_to_dart.f90 2009-04-15 20:10:12 UTC (rev 3811)
@@ -0,0 +1,319 @@
+! Data Assimilation Research Testbed -- DART
+! Copyright 2004-2007, Data Assimilation Research Section
+! University Corporation for Atmospheric Research
+! Licensed under the GPL -- www.gpl.org/licenses/gpl.html
+
+PROGRAM wrf_to_dart
+
+! <next few lines under version control, do not edit>
+! $URL$
+! $Id$
+! $Revision$
+! $Date$
+
+use types_mod, only : r8
+use time_manager_mod, only : time_type, write_time, read_time, get_date, set_date, operator(-), &
+ get_time, print_time, set_calendar_type, GREGORIAN, julian_day
+use utilities_mod, only : get_unit, file_exist, open_file, close_file, &
+ error_handler, E_ERR, E_MSG, initialize_utilities, &
+ register_module, logfileunit, nmlfileunit, timestamp, &
+ find_namelist_in_file, check_namelist_read, &
+ nc_check
+use assim_model_mod, only : open_restart_read, open_restart_write, aread_state_restart, &
+ awrite_state_restart
+use model_mod, only : max_state_variables, &
+ num_state_table_columns, read_wrf_dimensions, &
+ get_number_of_wrf_variables, &
+ get_variable_size_from_file, wrf_dom, &
+ fill_default_state_table, &
+ trans_3Dto1D, trans_2Dto1D, &
+ get_wrf_date
+
+use netcdf
+
+implicit none
+
+! version controlled file description for error handling, do not edit
+character(len=128), parameter :: &
+ source = "$URL$", &
+ revision = "$Revision$", &
+ revdate = "$Date$"
+
+
+!-----------------------------------------------------------------------
+! Model namelist parameters with default values.
+!-----------------------------------------------------------------------
+
+logical :: output_state_vector = .false. ! state vs. prognostic format
+logical :: default_state_variables = .true. ! use default state list?
+character(len=129) :: wrf_state_variables(num_state_table_columns,max_state_variables) = 'NULL'
+integer :: num_moist_vars = 3
+integer :: num_domains = 1
+integer :: calendar_type = GREGORIAN
+integer :: assimilation_period_seconds = 21600
+logical :: surf_obs = .true.
+logical :: soil_data = .true.
+logical :: h_diab = .false.
+logical :: allow_obs_below_vol = .false.
+character(len = 72) :: adv_mod_command = './wrf.exe'
+real (kind=r8) :: center_search_half_length = 500000.0_r8
+integer :: center_spline_grid_scale = 10
+integer :: vert_localization_coord = 3 ! 1,2,3 == level,pressure,height
+! candidates for including in the WRF netcdf files:
+logical :: polar = .false. ! wrap over the poles
+logical :: periodic_x = .false. ! wrap in longitude or x
+logical :: periodic_y = .false. ! used for single column model, wrap in y
+!JPH -- single column model flag
+logical :: scm = .false. ! using the single column model
+
+
+namelist /model_nml/ output_state_vector, num_moist_vars, &
+ num_domains, calendar_type, surf_obs, soil_data, h_diab, &
+ default_state_variables, wrf_state_variables, &
+ adv_mod_command, assimilation_period_seconds, &
+ allow_obs_below_vol, vert_localization_coord, &
+ center_search_half_length, center_spline_grid_scale, &
+ polar, periodic_x, periodic_y, scm
+
+
+!-------------------------------------------------------------
+
+type(wrf_dom) :: wrf
+
+real(r8), pointer :: dart(:)
+real(r8), allocatable :: wrf_var_3d(:,:,:), wrf_var_2d(:,:)
+type(time_type) :: dart_time(2)
+integer :: number_dart_values, ndays, &
+ year, month, day, hour, minute, second
+integer :: ndims, idims(2), dimids(2)
+integer :: i, ivtype, ind, dart_ind, my_index
+character(len=80) :: varname
+character(len=19) :: timestring
+character(len=1) :: idom
+
+logical, parameter :: debug = .false.
+integer :: ncid(50), io, var_id, id, iunit, dart_unit
+integer :: var_element_list(max_state_variables)
+
+write(*,*) 'WRF to DART' !(.false./F)?'
+
+call initialize_utilities('wrf_to_dart')
+call register_module(source, revision, revdate)
+
+! Begin by reading the namelist input
+call find_namelist_in_file("input.nml", "model_nml", iunit)
+read(iunit, nml = model_nml, iostat = io)
+call check_namelist_read(iunit, io, "model_nml")
+
+! Record the namelist values used for the run ...
+write(nmlfileunit, nml=model_nml)
+write( * , nml=model_nml)
+
+call set_calendar_type(calendar_type)
+
+call error_handler(E_MSG,'wrf_to_dart', &
+ 'Converting a WRF netcdf file to a dart state vector', &
+ source, revision, revdate)
+
+allocate(wrf%dom(num_domains))
+
+! get default state variable table if asked
+if ( default_state_variables ) then
+ wrf_state_variables = 'NULL'
+ call fill_default_state_table(wrf_state_variables)
+endif
+
+if ( debug ) then
+ print*,'WRF state vector table'
+ print*,'default_state_variables = ',default_state_variables
+ print*,wrf_state_variables
+endif
+
+
+! open wrf data netCDF file 'wrfinput_d0x'
+! we get sizes of the WRF geometry and resolution
+
+! approach is to do everything in a loop, without any regard to the actual
+! order. just use the namelist to read and unroll one at a time.
+
+! big loop over domains, just like in static_init
+number_dart_values = 0
+WRFDomains : do id = 1,num_domains
+
+ write(idom,'(I1)') id
+
+ if(file_exist('wrfinput_d0'//idom)) then
+
+ call nc_check( nf90_open('wrfinput_d0'//idom, NF90_NOWRITE, ncid(id)), &
+ 'wrf_to_dart','open wrfinput_d0'//idom )
+
+ else
+
+ call error_handler(E_ERR,'wrf_to_dart', &
+ 'Please put wrfinput_d0'//idom//' in the work directory.', source, revision,revdate)
+
+ endif
+
+! read WRF dimensions
+ call read_wrf_dimensions(ncid(id),wrf%dom(id)%bt, wrf%dom(id)%bts, &
+ wrf%dom(id)%sn, wrf%dom(id)%sns, &
+ wrf%dom(id)%we, wrf%dom(id)%wes, &
+ wrf%dom(id)%sls)
+
+
+! get the number of wrf variables wanted in this domain's state
+ wrf%dom(id)%number_of_wrf_variables = get_number_of_wrf_variables(id,wrf_state_variables,var_element_list)
+
+ if ( debug ) then
+ print*,'Domain ',id,' number of wrf variables is ',wrf%dom(id)%number_of_wrf_variables
+ endif
+
+! allocate and store the table locations of the variables valid on this domain
+ allocate(wrf%dom(id)%var_index_list(wrf%dom(id)%number_of_wrf_variables))
+ wrf%dom(id)%var_index_list = var_element_list(1:wrf%dom(id)%number_of_wrf_variables)
+
+! allocate var size
+ allocate(wrf%dom(id)%var_size(3,wrf%dom(id)%number_of_wrf_variables))
+
+! allocate stagger
+ allocate(wrf%dom(id)%stagger(wrf%dom(id)%number_of_wrf_variables))
+
+! figure out the dart state vector size
+ do ind = 1,wrf%dom(id)%number_of_wrf_variables
+
+ ! actual location in state variable table
+ my_index = wrf%dom(id)%var_index_list(ind)
+
+ ! get stagger and variable size
+ call get_variable_size_from_file(ncid(id),id, &
+ wrf_state_variables(1,my_index), &
+ wrf%dom(id)%bt, wrf%dom(id)%bts, &
+ wrf%dom(id)%sn, wrf%dom(id)%sns, &
+ wrf%dom(id)%we, wrf%dom(id)%wes, &
+ wrf%dom(id)%stagger(ind), &
+ wrf%dom(id)%var_size(:,ind))
+ if ( debug ) then
+ print*,'variable size ',trim(wrf_state_variables(1,my_index)),' ',wrf%dom(id)%var_size(:,ind)
+ endif
+
+ number_dart_values = number_dart_values &
+ + (wrf%dom(id)%var_size(1,ind) &
+ * wrf%dom(id)%var_size(2,ind) &
+ * wrf%dom(id)%var_size(3,ind))
+ enddo
+
+enddo WRFDomains
+
+if ( debug ) then
+ print*,'dart vector size ',number_dart_values
+endif
+
+! allocate dart state vector
+allocate(dart(number_dart_values))
+
+
+! loop through domains again and append each variable to the state
+dart_ind = 1
+WRFDomains2 : do id = 1,num_domains
+
+ do ind = 1,wrf%dom(id)%number_of_wrf_variables
+
+ ! actual location in state variable table
+ my_index = wrf%dom(id)%var_index_list(ind)
+
+ ! get stagger and variable size
+ call nc_check( nf90_inq_varid(ncid(id),wrf_state_variables(1,my_index), &
+ var_id), 'wrf_to_dart', &
+ 'inq_var_id '//wrf_state_variables(1,my_index))
+
+ if ( wrf%dom(id)%var_size(3,ind) == 1 ) then
+
+ allocate(wrf_var_2d(wrf%dom(id)%var_size(1,ind),wrf%dom(id)%var_size(2,ind)))
+
+ call nc_check( nf90_get_var(ncid(id), var_id, wrf_var_2d), &
+ 'wrf_to_dart','get_var '//wrf_state_variables(1,my_index) )
+
+ call trans_2Dto1D( dart(dart_ind:),wrf_var_2d, &
+ wrf%dom(id)%var_size(1,ind),wrf%dom(id)%var_size(2,ind))
+
+ deallocate(wrf_var_2d)
+
+ else
+
+ allocate(wrf_var_3d(wrf%dom(id)%var_size(1,ind),wrf%dom(id)%var_size(2,ind),wrf%dom(id)%var_size(3,ind)))
+
+ call nc_check( nf90_get_var(ncid(id), var_id, wrf_var_3d), &
+ 'wrf_to_dart','get_var '//wrf_state_variables(1,my_index) )
+
+ call trans_3Dto1D( dart(dart_ind:),wrf_var_3d, &
+ wrf%dom(id)%var_size(1,ind), &
+ wrf%dom(id)%var_size(2,ind),wrf%dom(id)%var_size(3,ind))
+
+ deallocate(wrf_var_3d)
+
+ endif
+
+ dart_ind = dart_ind &
+ + (wrf%dom(id)%var_size(1,ind) &
+ * wrf%dom(id)%var_size(2,ind) &
+ * wrf%dom(id)%var_size(3,ind))
+
+ enddo
+
+enddo WRFDomains2
+
+!---
+! output
+
+if(debug) write(*,*) ' state output '
+call nc_check( nf90_inq_varid(ncid(1), "Times", var_id), 'wrf_to_dart', &
+ 'inq_varid Times' )
+call nc_check( nf90_inquire_variable(ncid(1), var_id, varname, xtype=ivtype, &
+ ndims=ndims, dimids=dimids), 'wrf_to_dart', &
+ 'inquire_variable Times' )
+do i=1,ndims
+ call nc_check( nf90_inquire_dimension(ncid(1), dimids(i), &
+ len=idims(i)),'wrf_to_dart','inquire_dimensions Times' )
+ if(debug) write(*,*) ' dimension ',i,idims(i)
+enddo
+
+call nc_check( nf90_get_var(ncid(1), var_id, timestring, &
+ start = (/ 1, idims(2) /)), 'wrf_to_dart','get_var Times' )
+
+call get_wrf_date(timestring, year, month, day, hour, minute, second)
+dart_time(1) = set_date(year, month, day, hour, minute, second)
+
+call print_time(dart_time(1),str='Time from wrfinput_d0x:')
+
+iunit = get_unit()
+if(file_exist('wrf.info')) then
+ open(unit = iunit, file = 'wrf.info')
+ dart_time(1) = read_time(iunit)
+ close(iunit)
+endif
+
+call print_time(dart_time(1),str='Time written to dart vector file:')
+
+! open DART data file
+
+dart_unit = open_restart_write("dart_wrf_vector")
+
+call awrite_state_restart(dart_time(1), dart, dart_unit)
+
+if(debug) write(*,*) ' returned from state output '
+
+do id=1,num_domains
+ call nc_check ( nf90_sync(ncid(id)),'wrf_to_dart','sync wrfinput' )
+ call nc_check ( nf90_close(ncid(id)),'wrf_to_dart','close wrfinput' )
+enddo
+
+deallocate(wrf%dom)
+deallocate(dart)
+
+write(logfileunit,*)'FINISHED wrf_to_dart.'
+write(logfileunit,*)
+
+call timestamp(source,revision,revdate,'end') ! That closes the log file, too.
+
+end PROGRAM wrf_to_dart
+
Property changes on: DART/trunk/models/wrf/WRF_DART_utilities/wrf_to_dart.f90
___________________________________________________________________
Name: svn:keywords
+ Date Revision Author HeadURL Id
Modified: DART/trunk/models/wrf/model_mod.f90
===================================================================
--- DART/trunk/models/wrf/model_mod.f90 2009-04-13 17:56:46 UTC (rev 3810)
+++ DART/trunk/models/wrf/model_mod.f90 2009-04-15 20:10:12 UTC (rev 3811)
@@ -53,7 +53,8 @@
KIND_ICE_NUMBER_CONCENTRATION, KIND_GEOPOTENTIAL_HEIGHT, &
KIND_POTENTIAL_TEMPERATURE, KIND_SOIL_MOISTURE, &
KIND_VORTEX_LAT, KIND_VORTEX_LON, &
- KIND_VORTEX_PMIN, KIND_VORTEX_WMAX
+ KIND_VORTEX_PMIN, KIND_VORTEX_WMAX,&
+ get_raw_obs_kind_index
!nc -- module_map_utils split the declarations of PROJ_* into a separate module called
@@ -99,12 +100,25 @@
! Here is the appropriate place for other users to make additional routines
! contained within model_mod available for public use:
public :: get_number_domains, &
- wrf_static_data_for_dart, &
get_wrf_static_data, &
model_pressure, &
- pres_to_zk
+ pres_to_zk, &
+ fill_default_state_table, &
+ read_wrf_dimensions, &
+ get_number_of_wrf_variables, &
+ get_variable_size_from_file, &
+ trans_3Dto1D, trans_1Dto3D, &
+ trans_2Dto1D, trans_1Dto2D, &
+ get_wrf_date, set_wrf_date
+! public parameters
+public :: max_state_variables, &
+ num_state_table_columns
+! types
+public :: wrf_dom, wrf_static_data_for_dart
+
+
!-----------------------------------------------------------------------
! version controlled file description for error handling, do not edit
character(len=128), parameter :: &
@@ -112,6 +126,10 @@
revision = "$Revision$", &
revdate = "$Date$"
+! miscellaneous
+integer, parameter :: max_state_variables = 100
+integer, parameter :: num_state_table_columns = 5
+
!-----------------------------------------------------------------------
! Model namelist parameters with default values.
!
@@ -121,6 +139,8 @@
!-----------------------------------------------------------------------
logical :: output_state_vector = .false. ! output prognostic variables
+logical :: default_state_variables = .true. ! use default state list?
+character(len=129) :: wrf_state_variables(num_state_table_columns,max_state_variables) = 'NULL'
integer :: num_moist_vars = 3
integer :: num_domains = 1
integer :: calendar_type = GREGORIAN
@@ -146,15 +166,16 @@
!JPH -- single column model flag
logical :: scm = .false. ! using the single column model
-real(r8), allocatable :: ens_mean(:)
-
namelist /model_nml/ output_state_vector, num_moist_vars, &
num_domains, calendar_type, surf_obs, soil_data, h_diab, &
+ default_state_variables, wrf_state_variables, &
adv_mod_command, assimilation_period_seconds, &
allow_obs_below_vol, vert_localization_coord, &
center_search_half_length, center_spline_grid_scale, &
polar, periodic_x, periodic_y, scm
+real(r8), allocatable :: ens_mean(:)
+
character(len = 20) :: wrf_nml_file = 'namelist.input'
!-----------------------------------------------------------------------
@@ -162,22 +183,12 @@
! Private definition of domain map projection use by WRF
!nc -- added in CASSINI and CYL according to module_map_utils convention
-!JPH -- change variable name from map_sphere to more general map_latlon
+!JPH -- change variable name from map_sphere to more specific map_latlon
integer, parameter :: map_latlon = 0, map_lambert = 1, map_polar_stereo = 2, map_mercator = 3
integer, parameter :: map_cassini = 6, map_cyl = 5
! Private definition of model variable types
-integer, parameter :: TYPE_U = 1, TYPE_V = 2, TYPE_W = 3, &
- TYPE_GZ = 4, TYPE_T = 5, TYPE_MU = 6, &
- TYPE_TSK = 7, TYPE_QV = 8, TYPE_QC = 9, &
- TYPE_QR = 10, TYPE_QI = 11, TYPE_QS = 12, &
- TYPE_QG = 13, TYPE_QNICE = 14, TYPE_U10 = 15, &
- TYPE_V10 = 16, TYPE_T2 = 17, TYPE_TH2 = 18, &
- TYPE_Q2 = 19, TYPE_PS = 20, TYPE_TSLB = 21, &
- TYPE_SMOIS = 22, TYPE_SH2O = 23, TYPE_HDIAB = 24
-integer, parameter :: num_model_var_types = 24
-
real (kind=r8), PARAMETER :: kappa = 2.0_r8/7.0_r8 ! gas_constant / cp
real (kind=r8), PARAMETER :: ts0 = 300.0_r8 ! Base potential temperature for all levels.
@@ -217,8 +228,10 @@
integer, dimension(:,:), pointer :: var_index
integer, dimension(:,:), pointer :: var_size
integer, dimension(:), pointer :: var_type
+ integer, dimension(:), pointer :: var_index_list
integer, dimension(:), pointer :: dart_kind
integer, dimension(:,:), pointer :: land
+ character(len=129),dimension(:),pointer :: description, units, stagger
integer, dimension(:,:,:,:), pointer :: dart_ind
@@ -231,6 +244,8 @@
type(wrf_dom) :: wrf
+! JPH move map stuff into common (can move back into S/R later?)
+real(r8) :: stdlon,truelat1,truelat2 !,latinc,loninc
contains
@@ -238,7 +253,7 @@
subroutine static_init_model()
-! INitializes class data for WRF???
+! Initializes class data for WRF
integer :: ncid
integer :: io, iunit
@@ -247,9 +262,9 @@
character (len=1) :: idom
logical, parameter :: debug = .false.
integer :: var_id, ind, i, j, k, id, dart_index, model_type
+integer :: my_index
+integer :: var_element_list(max_state_variables)
-integer :: proj_code
-real(r8) :: stdlon,truelat1,truelat2,latinc,loninc
character(len=129) :: errstring
!----------------------------------------------------------------------
@@ -268,8 +283,23 @@
allocate(wrf%dom(num_domains))
-wrf%dom(:)%n_moist = num_moist_vars
+! get default state variable table if asked
+if ( default_state_variables ) then
+ wrf_state_variables = 'NULL'
+ call fill_default_state_table(wrf_state_variables)
+endif
+if ( debug ) then
+ print*,'WRF state vector table'
+ print*,'default_state_variables = ',default_state_variables
+ print*,wrf_state_variables
+endif
+
+!---------------------------
+! next block to be obsolete
+!---------------------------
+wrf%dom(:)%n_moist = num_moist_vars
+
if( num_moist_vars > 7) then
write(*,'(''num_moist_vars = '',i3)')num_moist_vars
call error_handler(E_ERR,'static_init_model', &
@@ -279,6 +309,10 @@
wrf%dom(:)%surf_obs = surf_obs
wrf%dom(:)%soil_data = soil_data
+!---------------------------
+! end obsolete
+!---------------------------
+
if ( debug ) then
if ( output_state_vector ) then
write(*,*)'netcdf file in state vector format'
@@ -287,6 +321,7 @@
endif
endif
+! set calendar type
call set_calendar_type(calendar_type)
! Store vertical localization coordinate
@@ -306,15 +341,15 @@
call error_handler(E_ERR,'static_init_model', errstring, source, revision,revdate)
endif
-dart_index = 1
-
! the agreement amongst the dart/wrf users was that there was no need to
! read the wrf namelist, since the only thing it was extracting was the
! timestep, which is part of the wrf input netcdf file.
! call read_dt_from_wrf_nml()
-do id=1,num_domains
+dart_index = 1
+WRFDomains : do id=1,num_domains
+
write( idom , '(I1)') id
! only print this once, no matter how many parallel tasks are running
@@ -338,681 +373,157 @@
if(debug) write(*,*) ' ncid is ',ncid
-! get wrf grid dimensions
+!-------------------------------------------------------
+! read WRF dimensions
+!-------------------------------------------------------
+ call read_wrf_dimensions(ncid,wrf%dom(id)%bt, wrf%dom(id)%bts, &
+ wrf%dom(id)%sn, wrf%dom(id)%sns, &
+ wrf%dom(id)%we, wrf%dom(id)%wes, &
+ wrf%dom(id)%sls)
- call nc_check( nf90_inq_dimid(ncid, "bottom_top", var_id), &
- 'static_init_model','inq_dimid bottom_top')
- call nc_check( nf90_inquire_dimension(ncid, var_id, name, wrf%dom(id)%bt), &
- 'static_init_model','inquire_dimension '//trim(name))
+!-------------------------------------------------------
+! read WRF file attributes
+!-------------------------------------------------------
+ call read_wrf_file_attributes(ncid,id)
- call nc_check( nf90_inq_dimid(ncid, "bottom_top_stag", var_id), &
- 'static_init_model','inq_dimid bottom_top_stag') ! reuse var_id, no harm
- call nc_check( nf90_inquire_dimension(ncid, var_id, name, wrf%dom(id)%bts), &
- 'static_init_model','inquire_dimension '//trim(name))
+!-------------------------------------------------------
+! assign boundary condition flags
+!-------------------------------------------------------
- call nc_check( nf90_inq_dimid(ncid, "south_north", var_id), &
- 'static_init_model','inq_dimid south_north')
- call nc_check( nf90_inquire_dimension(ncid, var_id, name, wrf%dom(id)%sn), &
- 'static_init_model','inquire_dimension '//trim(name))
+ call assign_boundary_conditions(id)
- call nc_check( nf90_inq_dimid(ncid, "south_north_stag", var_id), &
- 'static_init_model','inq_dimid south_north_stag') ! reuse var_id, no harm
- call nc_check( nf90_inquire_dimension(ncid, var_id, name, wrf%dom(id)%sns), &
- 'static_init_model','inquire_dimension '//trim(name))
+!-------------------------------------------------------
+! read static data
+!-------------------------------------------------------
- call nc_check( nf90_inq_dimid(ncid, "west_east", var_id), &
- 'static_init_model','inq_dimid west_east')
- call nc_check( nf90_inquire_dimension(ncid, var_id, name, wrf%dom(id)%we), &
- 'static_init_model','inquire_dimension '//trim(name))
+ call read_wrf_static_data(ncid,id)
- call nc_check( nf90_inq_dimid(ncid, "west_east_stag", var_id), &
- 'static_init_model','inq_dimid west_east_stag') ! reuse var_id, no harm
- call nc_check( nf90_inquire_dimension(ncid, var_id, name, wrf%dom(id)%wes), &
- 'static_init_model','inquire_dimension '//trim(name))
- call nc_check( nf90_inq_dimid(ncid, "soil_layers_stag", var_id), &
- 'static_init_model','inq_dimid soil_layers_stag') ! reuse var_id, no harm
- call nc_check( nf90_inquire_dimension(ncid, var_id, name, wrf%dom(id)%sls), &
- 'static_init_model','inquire_dimension '//trim(name))
+!-------------------------------------------------------
+! next block set up map
+!-------------------------------------------------------
- if(debug) then
- write(*,*) ' dimensions bt, sn, we are ',wrf%dom(id)%bt, &
- wrf%dom(id)%sn, wrf%dom(id)%we
- write(*,*) ' staggered bt, sn, we are ',wrf%dom(id)%bts, &
- wrf%dom(id)%sns,wrf%dom(id)%wes
- endif
+ call setup_map_projection(id)
-! get meta data and static data we need
+!-------------------------------------------------------
+! end block set up map
+!-------------------------------------------------------
- call nc_check( nf90_get_att(ncid, nf90_global, 'DX', wrf%dom(id)%dx), &
- 'static_init_model', 'get_att DX')
- call nc_check( nf90_get_att(ncid, nf90_global, 'DY', wrf%dom(id)%dy), &
- 'static_init_model', 'get_att DY')
- call nc_check( nf90_get_att(ncid, nf90_global, 'DT', wrf%dom(id)%dt), &
- 'static_init_model', 'get_att DT')
- if (do_output()) print*,'dt from wrfinput file is: ', wrf%dom(id)%dt
- if(debug) write(*,*) ' dx, dy, dt are ',wrf%dom(id)%dx, &
- wrf%dom(id)%dy, wrf%dom(id)%dt
+! get the number of wrf variables wanted in this domain's state
+ wrf%dom(id)%number_of_wrf_variables = get_number_of_wrf_variables(id,wrf_state_variables,var_element_list)
- call nc_check( nf90_get_att(ncid, nf90_global, 'MAP_PROJ', wrf%dom(id)%map_proj), &
- 'static_init_model', 'get_att MAP_PROJ')
- if(debug) write(*,*) ' map_proj is ',wrf%dom(id)%map_proj
+! allocate and store the table locations of the variables valid on this domain
+ allocate(wrf%dom(id)%var_index_list(wrf%dom(id)%number_of_wrf_variables))
+ wrf%dom(id)%var_index_list = var_element_list(1:wrf%dom(id)%number_of_wrf_variables)
- call nc_check( nf90_get_att(ncid, nf90_global, 'CEN_LAT', wrf%dom(id)%cen_lat), &
- 'static_init_model', 'get_att CEN_LAT')
- if(debug) write(*,*) ' cen_lat is ',wrf%dom(id)%cen_lat
+! allocation for wrf variable types
+ allocate(wrf%dom(id)%var_type(wrf%dom(id)%number_of_wrf_variables))
- call nc_check( nf90_get_att(ncid, nf90_global, 'CEN_LON', wrf%dom(id)%cen_lon), &
- 'static_init_model', 'get_att CEN_LON')
- if(debug) write(*,*) ' cen_lon is ',wrf%dom(id)%cen_lon
+! allocation for dart kinds
+ allocate(wrf%dom(id)%dart_kind(wrf%dom(id)%number_of_wrf_variables))
- call nc_check( nf90_get_att(ncid, nf90_global, 'TRUELAT1', truelat1), &
- 'static_init_model', 'get_att TRUELAT1')
- if(debug) write(*,*) ' truelat1 is ',truelat1
+! allocation of var size
+ allocate(wrf%dom(id)%var_size(3,wrf%dom(id)%number_of_wrf_variables))
- call nc_check( nf90_get_att(ncid, nf90_global, 'TRUELAT2', truelat2), &
- 'static_init_model', 'get_att TRUELAT2')
- if(debug) write(*,*) ' truelat2 is ',truelat2
+! allocation for wrf variable metadata
+ allocate(wrf%dom(id)%stagger(wrf%dom(id)%number_of_wrf_variables))
+ allocate(wrf%dom(id)%description(wrf%dom(id)%number_of_wrf_variables))
+ allocate(wrf%dom(id)%units(wrf%dom(id)%number_of_wrf_variables))
- call nc_check( nf90_get_att(ncid, nf90_global, 'STAND_LON', stdlon), &
- 'static_init_model', 'get_att STAND_LON')
- if(debug) write(*,*) ' stdlon is ',stdlon
+! build the variable indices
+! this accounts for the fact that some variables might not be on all domains
+ do ind = 1,wrf%dom(id)%number_of_wrf_variables
-!nc -- fill in the boundary conditions (periodic_x and polar) here. This will
@@ Diff output truncated at 40000 characters. @@
More information about the Dart-dev
mailing list