From e7510f75d0367593711e0a0fc3df6c272524f4d2 Mon Sep 17 00:00:00 2001 From: Angus Gibson Date: Thu, 21 Apr 2016 16:24:20 +1000 Subject: [PATCH] Add UNIFORM for DIAG_REMAP_Z_GRID_DEF This is similar to specifying UNIFORM for the ALE_COORDINATE_CONFIG, however in this case the z-remapped diagnostics use the same number of levels as the model, uniformly distributed through depth. --- src/framework/MOM_diag_mediator.F90 | 42 ++++++++++++++++++++--------- 1 file changed, 29 insertions(+), 13 deletions(-) diff --git a/src/framework/MOM_diag_mediator.F90 b/src/framework/MOM_diag_mediator.F90 index 4ecd496abf..0abfa93c4a 100644 --- a/src/framework/MOM_diag_mediator.F90 +++ b/src/framework/MOM_diag_mediator.F90 @@ -246,14 +246,23 @@ subroutine set_axes_info(G, GV, param_file, diag_cs, set_vertical) " FILE:, - where is a file within\n"//& " the INPUTDIR, is\n"//& " the name of the variable that\n"//& - " contains interface positions.\n",& + " contains interface positions.\n"//& + " UNIFORM - vertical grid is uniform\n"//& + " between surface and max depth.\n",& default="") if (len_trim(string) > 0) then - if (index(trim(string),'FILE:')/=1) then - call MOM_error(FATAL,"set_axes_info: "//& - "Only the 'FILE:file,var' format has been implemented so far. "//& - "Found '"//trim(string)//"'") - else + if (trim(string) == 'UNIFORM') then + ! initialise a uniform coordinate with depth + nzi(1) = G%ke + 1 + allocate(diag_cs%zi_remap(nzi(1))) + allocate(diag_cs%zl_remap(nzi(1) - 1)) + + diag_cs%zi_remap(1) = 0 + do k = 2,nzi(1) + diag_cs%zi_remap(K) = diag_cs%zi_remap(K - 1) + G%max_depth / real(nzi(1) - 1) + enddo + elseif (index(trim(string), 'FILE:') == 1) then + ! read coordinate information from a file if (string(6:6)=='.' .or. string(6:6)=='/') then inputdir = "." filename = trim(extractWord(trim(string(6:200)), 1)) @@ -277,15 +286,22 @@ subroutine set_axes_info(G, GV, param_file, diag_cs, set_vertical) ! Log the expanded result as a comment since it cannot be read back in call log_param(param_file, mod, "! Remapping z diagnostics", & trim(inputdir)//"/"//trim(filename)//","//trim(varname)) + + ! Get interface dimensions + call field_size(filename, varname, nzi) + call assert(nzi(1) /= 0, 'set_axes_info: bad z-axis dimension size') + allocate(diag_cs%zi_remap(nzi(1))) + allocate(diag_cs%zl_remap(nzi(1) - 1)) + call MOM_read_data(filename, varname, diag_cs%zi_remap) + else + ! unsupported method + call MOM_error(FATAL,"set_axes_info: "//& + "Incorrect remapping grid specification. Only 'FILE:file,var' and"//& + "'UNIFORM' are currently supported."//& + "Found '"//trim(string)//"'") endif - ! Get interface dimensions - call field_size(filename, varname, nzi) - call assert(nzi(1) /= 0, 'set_axes_info: bad z-axis dimension size') - allocate(diag_cs%zi_remap(nzi(1))) - allocate(diag_cs%zl_remap(nzi(1) - 1)) - call MOM_read_data(filename, varname, diag_cs%zi_remap) - diag_cs%zi_remap(:) = abs( diag_cs%zi_remap(:) ) ! Always convert heights into depths + diag_cs%zi_remap(:) = abs(diag_cs%zi_remap(:)) ! Always convert heights into depths ! Calculate layer positions diag_cs%zl_remap(:) = diag_cs%zi_remap(1:nzi(1)-1) + & (diag_cs%zi_remap(2:) - diag_cs%zi_remap(:nzi(1)-1)) / 2