10 implicit none ;
private
19 real,
allocatable,
dimension(:) :: coordinateresolution
22 real,
allocatable,
dimension(:) :: target_density
25 real,
allocatable,
dimension(:) :: max_interface_depths
28 real,
allocatable,
dimension(:) :: max_layer_thickness
34 public init_coord_hycom, set_hycom_params, build_hycom1_column, end_coord_hycom
39 subroutine init_coord_hycom(CS, nk, coordinateResolution, target_density, interp_CS)
41 integer,
intent(in) :: nk
42 real,
dimension(nk),
intent(in) :: coordinateresolution
43 real,
dimension(nk+1),
intent(in) :: target_density
46 if (
associated(cs))
call mom_error(fatal,
"init_coord_hycom: CS already associated!")
48 allocate(cs%coordinateResolution(nk))
49 allocate(cs%target_density(nk+1))
52 cs%coordinateResolution(:) = coordinateresolution(:)
53 cs%target_density(:) = target_density(:)
54 cs%interp_CS = interp_cs
55 end subroutine init_coord_hycom
58 subroutine end_coord_hycom(CS)
62 if (.not.
associated(cs))
return
63 deallocate(cs%coordinateResolution)
64 deallocate(cs%target_density)
65 if (
allocated(cs%max_interface_depths))
deallocate(cs%max_interface_depths)
66 if (
allocated(cs%max_layer_thickness))
deallocate(cs%max_layer_thickness)
68 end subroutine end_coord_hycom
71 subroutine set_hycom_params(CS, max_interface_depths, max_layer_thickness, interp_CS)
73 real,
dimension(:),
optional,
intent(in) :: max_interface_depths
74 real,
dimension(:),
optional,
intent(in) :: max_layer_thickness
77 if (.not.
associated(cs))
call mom_error(fatal,
"set_hycom_params: CS not associated")
79 if (
present(max_interface_depths))
then
80 if (
size(max_interface_depths) /= cs%nk+1) &
81 call mom_error(fatal,
"set_hycom_params: max_interface_depths inconsistent size")
82 allocate(cs%max_interface_depths(cs%nk+1))
83 cs%max_interface_depths(:) = max_interface_depths(:)
86 if (
present(max_layer_thickness))
then
87 if (
size(max_layer_thickness) /= cs%nk) &
88 call mom_error(fatal,
"set_hycom_params: max_layer_thickness inconsistent size")
89 allocate(cs%max_layer_thickness(cs%nk))
90 cs%max_layer_thickness(:) = max_layer_thickness(:)
93 if (
present(interp_cs)) cs%interp_CS = interp_cs
94 end subroutine set_hycom_params
97 subroutine build_hycom1_column(CS, eqn_of_state, nz, depth, h, T, S, p_col, &
98 z_col, z_col_new, zScale, h_neglect, h_neglect_edge)
100 type(
eos_type),
pointer :: eqn_of_state
101 integer,
intent(in) :: nz
102 real,
intent(in) :: depth
103 real,
dimension(nz),
intent(in) :: t
104 real,
dimension(nz),
intent(in) :: s
105 real,
dimension(nz),
intent(in) :: h
106 real,
dimension(nz),
intent(in) :: p_col
107 real,
dimension(nz+1),
intent(in) :: z_col
108 real,
dimension(CS%nk+1),
intent(inout) :: z_col_new
109 real,
optional,
intent(in) :: zscale
111 real,
optional,
intent(in) :: h_neglect
114 real,
optional,
intent(in) :: h_neglect_edge
120 real,
dimension(nz) :: rho_col
121 real,
dimension(CS%nk) :: h_col_new
126 logical :: maximum_depths_set
127 logical :: maximum_h_set
129 maximum_depths_set =
allocated(cs%max_interface_depths)
130 maximum_h_set =
allocated(cs%max_layer_thickness)
132 z_scale = 1.0 ;
if (
present(zscale)) z_scale = zscale
139 rho_col(k) = min( rho_col(k), rho_col(k+1) )
144 call build_and_interpolate_grid(cs%interp_CS, rho_col, nz, h(:), z_col, &
145 cs%target_density, cs%nk, h_col_new, z_col_new, h_neglect, h_neglect_edge)
150 stretching = z_col(nz+1) / depth
152 nominal_z = nominal_z + (z_scale * cs%coordinateResolution(k-1)) * stretching
153 z_col_new(k) = max( z_col_new(k), nominal_z )
154 z_col_new(k) = min( z_col_new(k), z_col(nz+1) )
157 if (maximum_depths_set .and. maximum_h_set)
then ;
do k=2,cs%nk
160 z_col_new(k) = min(z_col_new(k), cs%max_interface_depths(k), &
161 z_col_new(k-1) + cs%max_layer_thickness(k-1))
162 enddo ;
elseif (maximum_depths_set)
then ;
do k=2,cs%nk
163 z_col_new(k) = min(z_col_new(k), cs%max_interface_depths(k))
164 enddo ;
elseif (maximum_h_set)
then ;
do k=2,cs%nk
165 z_col_new(k) = min(z_col_new(k), z_col_new(k-1) + cs%max_layer_thickness(k-1))
167 end subroutine build_hycom1_column