[Dart-dev] [5904] DART/branches/development/models/mpas_ocn: The converters work to/ from with ZERO effect on the variables in the
nancy at ucar.edu
nancy at ucar.edu
Wed Oct 24 13:51:56 MDT 2012
Revision: 5904
Author: thoar
Date: 2012-10-24 13:51:55 -0600 (Wed, 24 Oct 2012)
Log Message:
-----------
The converters work to/from with ZERO effect on the variables in the
mpas netcdf file I started with - IFF you use
&assim_model_nml
write_binary_restart_files = .true.,
if FALSE, the errors are ~ 1.e-14: roundoff in the conversion to/from ASCII.
I set the default mpas_ocn variables that we are likely to use in the
namelists, and added a netCDF attribute in the mpas_ocn files on
the variables modified by DART.
* The 'h' variable is a bit of a problem right now - I suspect it is not
sea surface height but really a layer thickness.
* mpas uses potential temperature - we have yet to convert that for the
model_mod representation.
* The wet/dry cells are not fully implemented yet.
* The vertical coordinate variables and get_state_meta_data() are not
complete. Nathan indicates a new file is on its way that will allow
this task to be completed.
Modified Paths:
--------------
DART/branches/development/models/mpas_ocn/README
DART/branches/development/models/mpas_ocn/dart_to_model.f90
DART/branches/development/models/mpas_ocn/model_mod.f90
DART/branches/development/models/mpas_ocn/model_mod.nml
DART/branches/development/models/mpas_ocn/model_to_dart.f90
DART/branches/development/models/mpas_ocn/model_to_dart.nml
DART/branches/development/models/mpas_ocn/work/input.nml
Added Paths:
-----------
DART/branches/development/models/mpas_ocn/dart_to_model.nml
-------------- next part --------------
Modified: DART/branches/development/models/mpas_ocn/README
===================================================================
--- DART/branches/development/models/mpas_ocn/README 2012-10-18 22:29:48 UTC (rev 5903)
+++ DART/branches/development/models/mpas_ocn/README 2012-10-24 19:51:55 UTC (rev 5904)
@@ -130,3 +130,47 @@
TJH - aka - SoYoung ... making changes to test the whole svn commit concept.
+
+From: Mark Petersen [mpetersen at lanl.gov]
+Sent: Thursday, September 13, 2012 3:50 PM
+To: Urban, Nathan
+Subject: restart variables
+
+Nathan,
+
+For next week, you should take a look at the restart variables in the
+Registry file at
+mpas_repo/trunk/mpas/src/core_ocean
+
+Every variable in mpas-o is listed there with a line like
+var persistent real h ( nVertLevels nCells Time ) 2 iro h state - -
+
+iro means that the variable is specified in input, restart, and output
+files. For you, you need to know what is required in the restart
+files. The great majority of these variables are grid variables, like
+
+var persistent real latVertex ( nVertices ) 0 iro latVertex mesh - -
+var persistent integer cellsOnEdge ( TWO nEdges ) 0 iro cellsOnEdge mesh - -
+etc.
+
+The important lines for you are:
+% Prognostic variables: read from input, saved in restart, and written to output
+var persistent real u (nVertLevels nEdges Time) 2 ir u state - -
+var persistent real h (nVertLevels nCells Time) 2 iro h state - -
+var persistent real rho (nVertLevels nCells Time) 2 iro rho state - -
+var persistent real temperature(nVertLevels nCells Time) 2 iro temperature state tracers dynamics
+var persistent real salinity(nVertLevels nCells Time) 2 iro salinity state tracers dynamics
+var persistent real tracer1(nVertLevels nCells Time) 2 iro tracer1 state tracers testing
+
+The variables you need to modify with DART are u,h,temperature, salinity.
+Rho is recomputed by the EOS upon start-up for our runs.
+It is in the restart file when config_vert_grid_type = 'isopycnal'
+in which case rho is specified rather than computed from the EOS.
+You can ignore tracer1 as well as it is just used for a test variable to
+test conservation. Or you can think of tracer1 as nitrogen, carbon,etc
+if that is helpful.
+
+Feel free to write me with questions next week while you are at NCAR, if
+that is helpful.
+
+Mark
Modified: DART/branches/development/models/mpas_ocn/dart_to_model.f90
===================================================================
--- DART/branches/development/models/mpas_ocn/dart_to_model.f90 2012-10-18 22:29:48 UTC (rev 5903)
+++ DART/branches/development/models/mpas_ocn/dart_to_model.f90 2012-10-24 19:51:55 UTC (rev 5904)
@@ -88,8 +88,10 @@
allocate(statevector(x_size))
write(*,*)
-write(*,*) 'dart_to_model: converting DART file ', "'"//trim(dart_to_model_input_file)//"'"
-write(*,*) 'to model analysis file ', "'"//trim(model_analysis_filename)//"'"
+write(*,*) 'dart_to_model: converting DART file ', &
+ "'"//trim(dart_to_model_input_file)//"'"
+write(*,*) ' to model analysis file ', &
+ "'"//trim(model_analysis_filename)//"'"
!----------------------------------------------------------------------
! Reads the valid time, the state, and the target time.
@@ -104,14 +106,12 @@
endif
call close_restart(iunit)
-print *, 'read state vector'
!----------------------------------------------------------------------
! update the current model state vector
! Convey the amount of time to integrate the model ...
! time_manager_nml: stop_option, stop_count increments
!----------------------------------------------------------------------
-print *, 'calling analysis file to state vector'
call statevector_to_analysis_file(statevector, model_analysis_filename, model_time)
! write time into in text format (YYYY-MM-DD_hh:mm:ss) into a file.
Added: DART/branches/development/models/mpas_ocn/dart_to_model.nml
===================================================================
--- DART/branches/development/models/mpas_ocn/dart_to_model.nml (rev 0)
+++ DART/branches/development/models/mpas_ocn/dart_to_model.nml 2012-10-24 19:51:55 UTC (rev 5904)
@@ -0,0 +1,5 @@
+&dart_to_model_nml
+ dart_to_model_input_file = 'dart_restart',
+ advance_time_present = .true.,
+ /
+
Property changes on: DART/branches/development/models/mpas_ocn/dart_to_model.nml
___________________________________________________________________
Added: svn:mime-type
+ text/plain
Added: svn:eol-style
+ native
Modified: DART/branches/development/models/mpas_ocn/model_mod.f90
===================================================================
--- DART/branches/development/models/mpas_ocn/model_mod.f90 2012-10-18 22:29:48 UTC (rev 5903)
+++ DART/branches/development/models/mpas_ocn/model_mod.f90 2012-10-24 19:51:55 UTC (rev 5904)
@@ -18,7 +18,17 @@
! be used by converters and utilities and those interfaces can be anything
! that is useful to other pieces of code.
+! Units on everything are MKS:
+!
+! u (velocity): meter / second
+! h (depth) : meter
+! rho (density) : kilograms / meter^3
+! temperature : *potential temperature* degrees C
+! salinity : PSU
+!
+! Note: the 'temperature' variable is *potential* temperature.
+
! Routines in other modules that are used here.
use types_mod, only : r4, r8, digits12, SECPERDAY, MISSING_R8, &
@@ -31,7 +41,7 @@
use location_mod, only : location_type, get_dist, query_location, &
get_close_maxdist_init, get_close_type, &
- set_location, get_location, horiz_dist_only, &
+ set_location, get_location, horiz_dist_only, &
write_location, &
vert_is_undef, VERTISUNDEF, &
vert_is_surface, VERTISSURFACE, &
@@ -54,16 +64,20 @@
open_file, file_exist, find_textfile_dims, &
file_to_text, close_file, do_nml_file, do_nml_term
-use obs_kind_mod, only : paramname_length, &
- get_raw_obs_kind_index, &
- get_raw_obs_kind_name, &
- KIND_TEMPERATURE, &
- KIND_SALINITY, &
- KIND_DRY_LAND, &
- KIND_U_CURRENT_COMPONENT, &
- KIND_V_CURRENT_COMPONENT, &
- KIND_SEA_SURFACE_HEIGHT, &
- KIND_SEA_SURFACE_PRESSURE
+use obs_kind_mod, only : paramname_length, &
+ get_raw_obs_kind_index, &
+ get_raw_obs_kind_name, &
+ KIND_VERTICAL_VELOCITY, &
+ KIND_POTENTIAL_TEMPERATURE, &
+ KIND_TEMPERATURE, &
+ KIND_SALINITY, &
+ KIND_DRY_LAND, &
+ KIND_EDGE_NORMAL_SPEED, &
+ KIND_U_CURRENT_COMPONENT, &
+ KIND_V_CURRENT_COMPONENT, &
+ KIND_SEA_SURFACE_HEIGHT, &
+ KIND_SEA_SURFACE_PRESSURE, &
+ KIND_TRACER_CONCENTRATION
use mpi_utilities_mod, only: my_task_id
@@ -78,7 +92,7 @@
! at a much greater execution-time code. they were used to interpolate
! values at arbitrary locations, not just at cell centers. with too small
! a set of basis points, the values were discontinuous at cell boundaries;
-! with too many the values were too smoothed out. we went back to
+! with too many the values were too smoothed out. we went back to
! barycentric interpolation in triangles formed by the three cell centers
! that enclosed the given point.
use get_geometry_mod
@@ -139,8 +153,9 @@
real(r8), parameter :: rcv = rgas/(cp-rgas)
! earth radius; needed to convert lat/lon to x,y,z cartesian coords.
-! FIXME: one of the example ocean files had a global attr with 6371220.0
-! instead of 1229. ??
+! FIXME: the example ocean files has a global attr with 6371220.0
+! the actual mpas code may have hardwired values (it did in the atmosphere)
+! need to check what is really going on.
real(r8), parameter :: radius = 6371229.0 ! meters
! roundoff error
@@ -154,14 +169,14 @@
type(xyz_get_close_type) :: cc_gc
type(xyz_location_type), allocatable :: cell_locs(:)
-! variables which are in the module namelist
+logical :: log_vert_interp = .false. ! if true, interpolate vertical in log space
+
+! variables which are in the module namelist
integer :: vert_localization_coord = VERTISHEIGHT
integer :: assimilation_period_days = 0
integer :: assimilation_period_seconds = 60
real(r8) :: model_perturbation_amplitude = 0.0001 ! tiny amounts
-real(r8) :: cell_size_meters = 100.0_r8 ! size of largest cell in r
logical :: output_state_vector = .false. ! output prognostic variables (if .false.)
-logical :: logp = .false. ! if true, interpolate pres in log space
integer :: debug = 0 ! turn up for more and more debug messages
integer :: xyzdebug = 0
character(len=32) :: calendar = 'Gregorian'
@@ -185,7 +200,7 @@
! if .true. edge normal winds (u) are updated based on U/V reconstructed winds
logical :: update_u_from_reconstruct = .true.
-! only if update_u_from_reconstruct is true,
+! only if update_u_from_reconstruct is true,
! if .false. use the cell center u,v reconstructed values to update edge winds
! if .true., read in the original u,v winds, compute the increments after the
! assimilation, and use only the increments to update the edge winds
@@ -226,7 +241,7 @@
character(len=NF90_MAX_NAME) :: units
character(len=NF90_MAX_NAME), dimension(NF90_MAX_VAR_DIMS) :: dimname
integer, dimension(NF90_MAX_VAR_DIMS) :: dimlens
- integer :: xtype ! netCDF variable type (NF90_double, etc.)
+ integer :: xtype ! netCDF variable type (NF90_double, etc.)
integer :: numdims ! number of dims - excluding TIME
integer :: numvertical ! number of vertical levels in variable
integer :: numcells ! number of horizontal locations (cell centers)
@@ -237,13 +252,13 @@
integer :: indexN ! location in dart state vector of last occurrence
integer :: dart_kind
character(len=paramname_length) :: kind_string
- logical :: clamping ! does variable need to be range-restricted before
+ logical :: clamping ! does variable need to be range-restricted before
real(r8) :: range(2) ! being stuffed back into MPAS analysis file.
end type progvartype
type(progvartype), dimension(max_state_variables) :: progvar
-! Grid parameters - the values will be read from an mpas analysis file.
+! Grid parameters - the values will be read from an mpas analysis file.
integer :: nCells = -1 ! Total number of cells making up the grid
integer :: nVertices = -1 ! Unique points in grid that are corners of cells
@@ -252,13 +267,12 @@
integer :: nVertLevels = -1 ! Vertical levels; count of vert cell centers
integer :: nVertLevelsP1 = -1 ! Vert levels plus 1; count of vert cell faces
integer :: vertexDegree = -1 ! Max number of cells/edges that touch any vertex
-integer :: nSoilLevels = -1 ! Number of soil layers
! scalar grid positions
! FIXME: we read in a lot of metadata about the grids. if space becomes an
! issue we could consider reading in only the x,y,z arrays for all the items
-! plus the radius, and then compute the lat/lon for locations needed by
+! plus the radius, and then compute the lat/lon for locations needed by
! get_state_meta_data() on demand. most of the computations we need to do
! are actually easier in xyz coords (no issue with poles).
@@ -274,11 +288,11 @@
real(r8), allocatable :: latEdge(:) ! edge longitudes (degrees, original radians in file)
real(r8), allocatable :: lonCell(:) ! cell center longitudes (degrees, original radians in file)
real(r8), allocatable :: latCell(:) ! cell center latitudes (degrees, original radians in file)
-real(r8), allocatable :: zGridFace(:,:) ! geometric height at cell faces (nVertLevelsP1,nCells)
-real(r8), allocatable :: zGridCenter(:,:) ! geometric height at cell centers (nVertLevels, nCells)
-real(r8), allocatable :: zGridEdge(:,:) ! geometric height at edge centers (nVertLevels, nEdges)
+real(r8), allocatable :: zGridFace(:,:) ! geometric depth at cell faces (nVertLevelsP1,nCells)
+real(r8), allocatable :: zGridCenter(:,:) ! geometric depth at cell centers (nVertLevels, nCells)
+real(r8), allocatable :: zGridEdge(:,:) ! geometric depth at edge centers (nVertLevels, nEdges)
-!real(r8), allocatable :: zEdgeFace(:,:) ! geometric height at edges faces (nVertLevelsP1,nEdges)
+real(r8), allocatable :: hZLevel(:) ! layer thicknesses - maybe - FIXME
!real(r8), allocatable :: zEdgeCenter(:,:) ! geometric height at edges faces (nVertLevels ,nEdges)
integer, allocatable :: cellsOnVertex(:,:) ! list of cell centers defining a triangle
@@ -355,7 +369,7 @@
! The triangle interpolation keeps a list of how many and which triangles
! overlap each regular lon-lat box. The number is stored in
! array triangle_num. The allocatable array
-! triangle_list lists the uniquen index
+! triangle_list lists the uniquen index
! of each overlapping triangle. The entry in
! triangle_start for a given regular lon-lat box indicates
! where the list of triangles begins in the triangle_list.
@@ -376,7 +390,7 @@
subroutine static_init_model()
! Called to do one time initialization of the model.
-!
+!
! All the grid information comes from the initialization of
! the dart_model_mod module.
@@ -430,25 +444,25 @@
call get_time(model_timestep,ss,dd) ! set_time() assures the seconds [0,86400)
-write(string1,*)'assimilation period is ',dd,' days ',ss,' seconds'
+write(string1,*)'WARNING: assimilation period is ',dd,' days ',ss,' seconds'
call error_handler(E_MSG,'static_init_model',string1,source,revision,revdate)
!---------------------------------------------------------------
! 1) get grid dimensions
-! 2) allocate space for the grids
+! 2) allocate space for the grids
! 3) read them from the analysis file
! read_grid_dims() fills in the following module global variables:
-! nCells, nVertices, nEdges, maxEdges, nVertLevels, nVertLevelsP1, vertexDegree, nSoilLevels
+! nCells, nVertices, nEdges, maxEdges, nVertLevels, nVertLevelsP1, vertexDegree
call read_grid_dims()
-allocate(latCell(nCells), lonCell(nCells))
+allocate(latCell(nCells), lonCell(nCells))
+
allocate(zGridFace(nVertLevelsP1, nCells))
allocate(zGridCenter(nVertLevels, nCells))
-
-!allocate(zEdgeFace( nVertLevelsP1, nEdges))
+allocate(zGridEdge(nVertLevels, nEdges))
+allocate(hZLevel( nVertLevels))
!allocate(zEdgeCenter(nVertLevels, nEdges))
-allocate(zGridEdge(nVertLevels, nEdges))
allocate(cellsOnVertex(vertexDegree, nVertices))
allocate(nEdgesOnCell(nCells))
@@ -456,22 +470,15 @@
allocate(cellsOnEdge(2, nEdges))
allocate(verticesOnCell(maxEdges, nCells))
allocate(edgeNormalVectors(3, nEdges))
-allocate(latEdge(nEdges), lonEdge(nEdges))
+allocate(latEdge(nEdges), lonEdge(nEdges))
allocate(xVertex(nVertices), yVertex(nVertices), zVertex(nVertices))
allocate(xEdge(nEdges), yEdge(nEdges), zEdge(nEdges))
-! this reads in latCell, lonCell, zGridFace, cellsOnVertex
+! this reads in latCell, lonCell, hZLevel, cellsOnVertex
call get_grid()
-! read in vert cell face locations and then compute vertical center locations
-do kloc=1, nCells
- do iloc=1, nVertLevels
- zGridCenter(iloc,kloc) = (zGridFace(iloc,kloc) + zGridFace(iloc+1,kloc))*0.5_r8
- enddo
-enddo
-
! FIXME: This code is supposed to check whether an edge has 2 neighbours or 1 neighbour and then
-! compute the height accordingly. HOWEVER, the array cellsOnEdge does not change with
+! compute the depth accordingly. HOWEVER, the array cellsOnEdge does not change with
! depth, but it should as an edge may have 2 neighbour cells at the top but not at depth.
do kloc=1, nEdges
do iloc=1, nVertLevels
@@ -489,7 +496,7 @@
endif
enddo
enddo
-
+
!---------------------------------------------------------------
! Compile the list of model variables to use in the creation
! of the DART state vector. Required to determine model_size.
@@ -547,7 +554,7 @@
call nc_check(nf90_inquire_variable(ncid, VarID, xtype=progvar(ivar)%xtype, &
dimids=dimIDs, ndims=numdims), 'static_init_model', 'inquire '//trim(string2))
- ! If the long_name and/or units attributes are set, get them.
+ ! If the long_name and/or units attributes are set, get them.
! They are not REQUIRED to exist but are nice to use if they are present.
if( nf90_inquire_attribute( ncid, VarID, 'long_name') == NF90_NOERR ) then
@@ -590,8 +597,6 @@
progvar(ivar)%numedges = dimlen
case ('nVertL') ! nVertLevels, nVertLevelsP1, nVertLevelsP2
progvar(ivar)%numvertical = dimlen
- case ('nSoilL') ! nSoilLevels
- progvar(ivar)%numvertical = dimlen
end select
enddo DimensionLoop
@@ -613,7 +618,7 @@
progvar(ivar)%indexN = index1 + varsize - 1
index1 = index1 + varsize ! sets up for next variable
- if ( debug > 4 ) call dump_progvar(ivar)
+ if ((debug > 4) .and. do_output()) call dump_progvar(ivar)
enddo
@@ -622,7 +627,7 @@
model_size = progvar(nfields)%indexN
-if ( debug > 0 .and. do_output()) then
+if ( (debug > 0) .and. do_output()) then
write(logfileunit,*)
write( * ,*)
write(logfileunit,'(" static_init_model: nCells, nEdges, nVertices, nVertLevels =",4(1x,i6))') &
@@ -636,21 +641,21 @@
write( * , *)'static_init_model: grid is a global grid '
else
write(logfileunit, *)'static_init_model: grid has boundaries '
- write( * , *)'static_init_model: grid has boundaries '
+ write( * , *)'static_init_model: grid has boundaries '
endif
if ( all_levels_exist_everywhere ) then
write(logfileunit, *)'static_init_model: all cells have same number of vertical levels '
write( * , *)'static_init_model: all cells have same number of vertical levels '
else
- write(logfileunit, *)'static_init_model: cells have varying number of vertical levels '
+ write(logfileunit, *)'static_init_model: cells have varying number of vertical levels '
write( * , *)'static_init_model: cells have varying number of vertical levels '
endif
endif
! do some sanity checking here:
-! if you have at least one of these wind components in the state vector,
-! you have to have them both. the subroutine will error out if only one
+! if you have at least one of these wind components in the state vector,
+! you have to have them both. the subroutine will error out if only one
! is found and not both.
if (has_uvreconstruct) then
call winds_present(z1,m1,both)
@@ -711,8 +716,8 @@
! TJH text2=string2, text3=string3)
! TJH endif
-string1 = 'fix block of required variables - detritus from atmosphere'
-call error_handler(E_ERR,'static_init_model',string1,source,revision,revdate)
+string1 = 'WARNING: fix block of required variables - detritus from atmosphere'
+call error_handler(E_MSG,'static_init_model',string1,source,revision,revdate)
allocate( ens_mean(model_size) )
@@ -732,13 +737,13 @@
integer, intent(in) :: index_in
type(location_type), intent(out) :: location
integer, optional, intent(out) :: var_type
-
+
! Local variables
integer :: nxp, nzp, iloc, vloc, nf, n
integer :: myindx
integer :: istatus
-real(r8) :: height
+real(r8) :: depth
type(location_type) :: new_location
if ( .not. module_initialized ) call static_init_model
@@ -775,37 +780,37 @@
iloc = 1 + (myindx-1) / nzp ! cell index
vloc = myindx - (iloc-1)*nzp ! vertical level index
-! the zGrid array contains the location of the cell top and bottom faces, so it has one
+! the zGrid array contains the location of the cell top and bottom faces, so it has one
! more value than the number of cells in each column. for locations of cell centers
! you have to take the midpoint of the top and bottom face of the cell.
if (progvar(nf)%numedges /= MISSING_I) then
if ( progvar(nf)%ZonHalf ) then
- height = zGridEdge(vloc,iloc)
+ depth = zGridEdge(vloc,iloc)
else
- call error_handler(E_ERR, 'get_state_meta_data', 'no support for edges at face heights', &
+ call error_handler(E_ERR, 'get_state_meta_data', 'no support for edges at face depths', &
source, revision, revdate)
endif
else
if ( progvar(nf)%ZonHalf ) then
- height = zGridCenter(vloc,iloc)
+ depth = zGridCenter(vloc,iloc)
else if (nzp <= 1) then
- height = zGridFace(1,iloc)
+ depth = zGridFace(1,iloc)
else
- height = zGridFace(vloc,iloc)
+ depth = zGridFace(vloc,iloc)
endif
endif
if (progvar(nf)%numedges /= MISSING_I) then
if (nzp <= 1) then
- location = set_location(lonEdge(iloc),latEdge(iloc), height, VERTISSURFACE)
+ location = set_location(lonEdge(iloc),latEdge(iloc), depth, VERTISSURFACE)
else
- location = set_location(lonEdge(iloc),latEdge(iloc), height, VERTISHEIGHT)
+ location = set_location(lonEdge(iloc),latEdge(iloc), depth, VERTISHEIGHT)
endif
else ! must be on cell centers
if (nzp <= 1) then
- location = set_location(lonCell(iloc),latCell(iloc), height, VERTISSURFACE)
+ location = set_location(lonCell(iloc),latCell(iloc), depth, VERTISSURFACE)
else
- location = set_location(lonCell(iloc),latCell(iloc), height, VERTISHEIGHT)
+ location = set_location(lonCell(iloc),latCell(iloc), depth, VERTISHEIGHT)
endif
endif
@@ -822,11 +827,11 @@
if(istatus == 0) location = new_location
endif
-if (debug > 12) then
+if ((debug > 12) .and. do_output()) then
write(*,'("INDEX_IN / myindx / IVAR / NX, NZ: ",2(i10,2x),3(i5,2x))') index_in, myindx, nf, nxp, nzp
write(*,'(" ILOC, KLOC: ",2(i5,2x))') iloc, vloc
- write(*,'(" LON/LAT/HGT: ",3(f12.3,2x))') lonCell(iloc), latCell(iloc), height
+ write(*,'(" LON/LAT/HGT: ",3(f12.3,2x))') lonCell(iloc), latCell(iloc), depth
endif
@@ -850,13 +855,13 @@
! ISTATUS = 99: general error in case something terrible goes wrong...
! ISTATUS = 88: this kind is not in the state vector
! ISTATUS = 11: Could not find a triangle that contains this lat/lon
-! ISTATUS = 12: Height vertical coordinate out of model range.
+! ISTATUS = 12: depth vertical coordinate out of model range.
! ISTATUS = 13: Missing value in interpolation.
! ISTATUS = 16: Don't know how to do vertical velocity for now
-! ISTATUS = 17: Unable to compute pressure values
+! ISTATUS = 17: Unable to compute pressure values
! ISTATUS = 18: altitude illegal
! ISTATUS = 19: could not compute u using RBF code
-! ISTATUS = 101: Internal error; reached end of subroutine without
+! ISTATUS = 101: Internal error; reached end of subroutine without
! finding an applicable case.
!
@@ -884,7 +889,7 @@
! to this subroutine, but this really is a kind.
obs_kind = obs_type
-if (debug > 0) then
+if ((debug > 0) .and. do_output()) then
call write_location(0,location,charstring=string1)
print *, my_task_id(), 'stt x(1), kind, loc ', x(1), obs_kind, trim(string1)
endif
@@ -897,8 +902,8 @@
else
goodkind = .false.
! exceptions if the kind isn't directly
- ! a field in the state vector:
- select case (obs_kind)
+ ! a field in the state vector:
+ select case (obs_kind)
!case (KIND_PRESSURE)
! goodkind = .true.
end select
@@ -907,7 +912,7 @@
! this kind is not in the state vector and it isn't one of the exceptions
! that we know how to handle.
if (.not. goodkind) then
- if (debug > 4) print *, 'kind rejected', obs_kind
+ if ((debug > 4) .and. do_output()) print *, 'kind rejected', obs_kind
istatus = 88
goto 100
endif
@@ -918,7 +923,7 @@
!
! compute_scalar_with_barycentric() can take up to 3 kinds at one location
! and return up to 3 values which can be combined or computed on here.
-!
+!
! to use the RBF functions, here is an example call. the third arg
! controls whether to return the meridional (.false.) or zonal (.true.) component.
! call compute_u_with_rbf(x, location, .TRUE., interp_val, istatus)
@@ -931,7 +936,7 @@
tvars(1) = ivar
call compute_scalar_with_barycentric(x, location, 1, tvars, values, istatus)
interp_val = values(1)
- if (debug > 4) print *, 'called generic compute_w_bary, kind, val, istatus: ', obs_kind, interp_val, istatus
+ if ((debug > 4) .and. do_output()) print *, 'called generic compute_w_bary, kind, val, istatus: ', obs_kind, interp_val, istatus
if (istatus /= 0) goto 100
!endif
@@ -952,12 +957,11 @@
call error_handler(E_ERR,'model_interpolate',string1,source,revision,revdate)
endif
-if (debug > 0) then
+if ((debug > 0) .and. do_output()) then
call write_location(0,location,charstring=string1)
print *, my_task_id(), 'end x(1), kind, loc, val, rc ', x(1), obs_kind, trim(string1), interp_val, istatus
endif
-
end subroutine model_interpolate
@@ -975,7 +979,7 @@
!
! Typical sequence for adding new dimensions,variables,attributes:
! NF90_OPEN ! open existing netCDF dataset
-! NF90_redef ! put into define mode
+! NF90_redef ! put into define mode
! NF90_def_dim ! define additional dimensions (if any)
! NF90_def_var ! define variables: from name, type, and dims
! NF90_put_att ! assign attribute values
@@ -1007,7 +1011,6 @@
integer :: nEdgesDimID, maxEdgesDimID
integer :: nVerticesDimID
integer :: VertexDegreeDimID
-integer :: nSoilLevelsDimID
integer :: nVertLevelsDimID
integer :: nVertLevelsP1DimID
@@ -1016,7 +1019,7 @@
integer :: ivar, VarID, mpasFileID
!----------------------------------------------------------------------
-! local variables
+! local variables
!----------------------------------------------------------------------
! we are going to need these to record the creation date in the netCDF file.
@@ -1049,7 +1052,7 @@
write(filename,*) 'ncFileID', ncFileID
!-------------------------------------------------------------------------------
-! make sure ncFileID refers to an open netCDF file,
+! make sure ncFileID refers to an open netCDF file,
! and then put into define mode.
!-------------------------------------------------------------------------------
@@ -1059,7 +1062,7 @@
!-------------------------------------------------------------------------------
! We need the dimension ID for the number of copies/ensemble members, and
-! we might as well check to make sure that Time is the Unlimited dimension.
+! we might as well check to make sure that Time is the Unlimited dimension.
! Our job is create the 'model size' dimension.
!-------------------------------------------------------------------------------
@@ -1082,7 +1085,7 @@
dimid = StateVarDimID),'nc_write_model_atts', 'state def_dim '//trim(filename))
!-------------------------------------------------------------------------------
-! Write Global Attributes
+! Write Global Attributes
!-------------------------------------------------------------------------------
call DATE_AND_TIME(crdate,crtime,crzone,values)
@@ -1113,7 +1116,7 @@
! Create a variable for the state vector
!----------------------------------------------------------------------------
- ! Define the actual (3D) state vector, which gets filled as time goes on ...
+ ! Define the actual (3D) state vector, which gets filled as time goes on ...
call nc_check(nf90_def_var(ncid=ncFileID, name='state', xtype=nf90_real, &
dimids=(/StateVarDimID,MemberDimID,unlimitedDimID/),varid=StateVarID),&
'nc_write_model_atts','state def_var '//trim(filename))
@@ -1156,10 +1159,6 @@
len = nVertLevelsP1, dimid = NVertLevelsP1DimID),'nc_write_model_atts', &
'nVertLevelsP1 def_dim '//trim(filename))
- call nc_check(nf90_def_dim(ncid=ncFileID, name='nSoilLevels', &
- len = nSoilLevels, dimid = nSoilLevelsDimID),'nc_write_model_atts', &
- 'nSoilLevels def_dim '//trim(filename))
-
!----------------------------------------------------------------------------
! Create the (empty) Coordinate Variables and the Attributes
!----------------------------------------------------------------------------
@@ -1205,17 +1204,17 @@
'nc_write_model_atts', 'zCell long_name '//trim(filename))
! Grid vertical information
- call nc_check(nf90_def_var(ncFileID,name='zgrid',xtype=nf90_double, &
- dimids=(/ nVertLevelsP1DimID, nCellsDimID /) ,varid=VarID), &
- 'nc_write_model_atts', 'zgrid def_var '//trim(filename))
- call nc_check(nf90_put_att(ncFileID, VarID, 'long_name', 'grid zgrid'), &
- 'nc_write_model_atts', 'zgrid long_name '//trim(filename))
- call nc_check(nf90_put_att(ncFileID, VarID, 'units', 'meters'), &
- 'nc_write_model_atts', 'zgrid units '//trim(filename))
- call nc_check(nf90_put_att(ncFileID, VarID, 'positive', 'up'), &
- 'nc_write_model_atts', 'zgrid units '//trim(filename))
- call nc_check(nf90_put_att(ncFileID, VarID, 'cartesian_axis', 'Z'), &
- 'nc_write_model_atts', 'zgrid cartesian_axis '//trim(filename))
+! call nc_check(nf90_def_var(ncFileID,name='zgrid',xtype=nf90_double, &
+! dimids=(/ nVertLevelsP1DimID, nCellsDimID /) ,varid=VarID), &
+! 'nc_write_model_atts', 'zgrid def_var '//trim(filename))
+! call nc_check(nf90_put_att(ncFileID, VarID, 'long_name', 'grid zgrid'), &
+! 'nc_write_model_atts', 'zgrid long_name '//trim(filename))
+! call nc_check(nf90_put_att(ncFileID, VarID, 'units', 'meters'), &
+! 'nc_write_model_atts', 'zgrid units '//trim(filename))
+! call nc_check(nf90_put_att(ncFileID, VarID, 'positive', 'up'), &
+! 'nc_write_model_atts', 'zgrid units '//trim(filename))
+! call nc_check(nf90_put_att(ncFileID, VarID, 'cartesian_axis', 'Z'), &
+! 'nc_write_model_atts', 'zgrid cartesian_axis '//trim(filename))
! Vertex Longitudes
call nc_check(nf90_def_var(ncFileID,name='lonVertex', xtype=nf90_double, &
@@ -1279,7 +1278,7 @@
! match shape of the variable to the dimension IDs
- call define_var_dims(ncFileID, ivar, MemberDimID, unlimitedDimID, myndims, mydimids)
+ call define_var_dims(ncFileID, ivar, MemberDimID, unlimitedDimID, myndims, mydimids)
! define the variable and set the attributes
@@ -1304,7 +1303,7 @@
call nc_check(nf90_enddef(ncfileID), 'prognostic enddef '//trim(filename))
!----------------------------------------------------------------------------
- ! Fill the coordinate variables that DART needs and has locally
+ ! Fill the coordinate variables that DART needs and has locally
!----------------------------------------------------------------------------
call nc_check(NF90_inq_varid(ncFileID, 'lonCell', VarID), &
@@ -1327,10 +1326,10 @@
call nc_check(nf90_put_var(ncFileID, VarID, latEdge ), &
'nc_write_model_atts', 'latEdge put_var '//trim(filename))
- call nc_check(NF90_inq_varid(ncFileID, 'zgrid', VarID), &
- 'nc_write_model_atts', 'zgrid inq_varid '//trim(filename))
- call nc_check(nf90_put_var(ncFileID, VarID, zGridFace ), &
- 'nc_write_model_atts', 'zgrid put_var '//trim(filename))
+ call nc_check(NF90_inq_varid(ncFileID, 'hZLevel', VarID), &
+ 'nc_write_model_atts', 'hZLevel inq_varid '//trim(filename))
+ call nc_check(nf90_put_var(ncFileID, VarID, hZLevel ), &
+ 'nc_write_model_atts', 'hZLevel put_var '//trim(filename))
call nc_check(NF90_inq_varid(ncFileID, 'nEdgesOnCell', VarID), &
'nc_write_model_atts', 'nEdgesOnCell inq_varid '//trim(filename))
@@ -1350,7 +1349,7 @@
!----------------------------------------------------------------------------
! Fill the coordinate variables needed for plotting only.
! DART has not read these in, so we have to read them from the input file
- ! and parrot them to the DART output file.
+ ! and parrot them to the DART output file.
!----------------------------------------------------------------------------
call nc_check(nf90_open(trim(grid_definition_filename), NF90_NOWRITE, mpasFileID), &
@@ -1403,7 +1402,7 @@
'nc_write_model_atts', 'lonVertex inq_varid '//trim(filename))
call nc_check(nf90_put_var(ncFileID, VarID, data1d ), &
'nc_write_model_atts', 'lonVertex put_var '//trim(filename))
-
+
call nc_check(nf90_inq_varid(mpasFileID, 'latVertex', VarID), &
'nc_write_model_atts', 'latVertex inq_varid ')
call nc_check(nf90_get_var(mpasFileID, VarID, data1d ), &
@@ -1429,7 +1428,7 @@
!------------------------------------------------------------------
-function nc_write_model_vars( ncFileID, state_vec, copyindex, timeindex ) result (ierr)
+function nc_write_model_vars( ncFileID, state_vec, copyindex, timeindex ) result (ierr)
! TJH 29 Aug 2011 -- all errors are fatal, so the
! return code is always '0 == normal', since the fatal errors stop execution.
@@ -1455,7 +1454,7 @@
integer :: ierr ! return value of function
integer, dimension(NF90_MAX_VAR_DIMS) :: dimIDs, mystart, mycount
-character(len=NF90_MAX_NAME) :: varname
+character(len=NF90_MAX_NAME) :: varname
integer :: i, ivar, VarID, ncNdims, dimlen
integer :: TimeDimID, CopyDimID
@@ -1479,7 +1478,7 @@
write(filename,*) 'ncFileID', ncFileID
!-------------------------------------------------------------------------------
-! make sure ncFileID refers to an open netCDF file,
+! make sure ncFileID refers to an open netCDF file,
!-------------------------------------------------------------------------------
call nc_check(nf90_inq_dimid(ncFileID, 'copy', dimid=CopyDimID), &
@@ -1501,7 +1500,7 @@
! We need to process the prognostic variables.
!----------------------------------------------------------------------------
- do ivar = 1,nfields
+ do ivar = 1,nfields
varname = trim(progvar(ivar)%varname)
string2 = trim(filename)//' '//trim(varname)
@@ -1543,7 +1542,7 @@
where(dimIDs == TimeDimID) mystart = timeindex
where(dimIDs == TimeDimID) mycount = 1
- if ( debug > 9 ) then
+ if ((debug > 9) .and. do_output()) then
write(*,*)'nc_write_model_vars '//trim(varname)//' start is ',mystart(1:ncNdims)
write(*,*)'nc_write_model_vars '//trim(varname)//' count is ',mycount(1:ncNdims)
endif
@@ -1628,7 +1627,7 @@
function get_model_size()
-! Returns the size of the model as an integer.
+! Returns the size of the model as an integer.
! Required for all applications.
integer :: get_model_size
@@ -1672,10 +1671,10 @@
ens_mean = filter_ens_mean
-if (debug > 3) then
+if ((debug > 3) .and. do_output()) then
if (do_output()) print *, 'resetting ensemble mean: '
do ivar = 1, nfields
- call dump_progvar(ivar, x=ens_mean, minmaxonly=.true.)
+ call dump_progvar(ivar, x=ens_mean, fulldump=.false.)
enddo
endif
@@ -1706,7 +1705,7 @@
! A model may choose to provide a NULL INTERFACE by returning
! .false. for the interf_provided argument. This indicates to
! the filter that if it needs to generate perturbed states, it
-! may do so by adding a perturbation to each model state
+! may do so by adding a perturbation to each model state
! variable independently. The interf_provided argument
! should be returned as .true. if the model wants to do its own
! perturbing of states.
@@ -1730,8 +1729,8 @@
! the subsequent code is a pert routine which
! can be used to add minor perturbations which can be spun up.
!
-! if all values in a field are identical (i.e. 0.0) this
-! routine will not change those values since it won't
+! if all values in a field are identical (i.e. 0.0) this
+! routine will not change those values since it won't
! make a new value outside the original min/max of that
! variable in the state vector. to handle this case you can
! remove the min/max limit lines below.
@@ -1796,10 +1795,10 @@
! "base", their indices, and their distances to the "base" ...
! For vertical distance computations, general philosophy is to convert all
-! vertical coordinates to a common coordinate. This coordinate type can be
+! vertical coordinates to a common coordinate. This coordinate type can be
! defined in the namelist with the variable "vert_localization_coord".
! But we first try a single coordinate type as the model level here.
-! FIXME: We need to add more options later.
+! FIXME: We need to add more options later.
! Vertical conversion is carried out by the subroutine vert_convert.
@@ -1847,7 +1846,7 @@
istatus1 = 1
else if (base_which /= vert_localization_coord) then
call vert_convert(ens_mean, base_obs_loc, base_obs_kind, ztypeout, istatus1)
- if(debug > 5) then
+ if ((debug > 5) .and. do_output()) then
call write_location(0,base_obs_loc,charstring=string1)
call error_handler(E_MSG, 'get_close_obs: base_obs_loc',string1,source, revision, revdate)
endif
@@ -1889,7 +1888,7 @@
dist(k) = 1.0e9_r8
else
dist(k) = get_dist(base_obs_loc, local_obs_loc, base_obs_kind, obs_kind(t_ind))
-! if (debug > 4 .and. k < 100 .and. my_task_id() == 0) then
+! if ((debug > 4) .and. (k < 100) .and. do_output()) then
! print *, 'calling get_dist'
! call write_location(0,base_obs_loc,charstring=string2)
! call error_handler(E_MSG, 'get_close_obs: base_obs_loc',string2,source, revision, revdate)
@@ -1903,7 +1902,7 @@
enddo
endif
-if (debug > 2) then
+if ((debug > 2) .and. do_output()) then
call write_location(0,base_obs_loc,charstring=string2)
print *, 'get_close_obs: nclose, base_obs_loc ', num_close, trim(string2)
endif
@@ -1915,9 +1914,9 @@
subroutine init_time(time)
-! Companion interface to init_conditions. Returns a time that is somehow
+! Companion interface to init_conditions. Returns a time that is somehow
! appropriate for starting up a long integration of the model.
-! At present, this is only used if the namelist parameter
+! At present, this is only used if the namelist parameter
! start_from_restart is set to .false. in the program perfect_model_obs.
type(time_type), intent(out) :: time
@@ -1939,13 +1938,13 @@
! Returns a model state vector, x, that is some sort of appropriate
! initial condition for starting up a long integration of the model.
-! At present, this is only used if the namelist parameter
+! At present, this is only used if the namelist parameter
! start_from_restart is set to .false. in the program perfect_model_obs.
real(r8), intent(out) :: x(:)
if ( .not. module_initialized ) call static_init_model
-
+
! this shuts up the compiler warnings about unused variables
x = 0.0_r8
@@ -1962,13 +1961,13 @@
! Does a single timestep advance of the model. The input value of
! the vector x is the starting condition and x is updated to reflect
! the changed state after a timestep. The time argument is intent
-! in and is used for models that need to know the date/time to
+! in and is used for models that need to know the date/time to
! compute a timestep, for instance for radiation computations.
! This interface is only called IF the namelist parameter
-! async is set to 0 in perfect_model_obs or filter -OR- if the
+! async is set to 0 in perfect_model_obs or filter -OR- if the
! program integrate_model is to be used to advance the model
@@ Diff output truncated at 40000 characters. @@
More information about the Dart-dev
mailing list