Changeset 11993 for NEMO/trunk/src/OCE/DIA/diahth.F90
- Timestamp:
- 2019-11-28T11:20:53+01:00 (4 years ago)
- File:
-
- 1 edited
-
NEMO/trunk/src/OCE/DIA/diahth.F90 (modified) (4 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/OCE/DIA/diahth.F90
r11989 r11993 11 11 !! 3.2 ! 2009-07 (S. Masson) hc300 bugfix + cleaning + add new diag 12 12 !!---------------------------------------------------------------------- 13 #if defined key_diahth 14 !!---------------------------------------------------------------------- 15 !! 'key_diahth' : thermocline depth diag. 16 !!---------------------------------------------------------------------- 13 17 !! dia_hth : Compute varius diagnostics associated with the mixed layer 14 18 !!---------------------------------------------------------------------- … … 28 32 PUBLIC dia_hth_alloc ! routine called by nemogcm.F90 29 33 30 LOGICAL , SAVE :: l_hth!: thermocline-20d depths flag34 LOGICAL , PUBLIC, PARAMETER :: lk_diahth = .TRUE. !: thermocline-20d depths flag 31 35 32 36 ! note: following variables should move to local variables once iom_put is always used 33 37 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hth !: depth of the max vertical temperature gradient [m] 34 38 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hd20 !: depth of 20 C isotherm [m] 35 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hd26 !: depth of 26 C isotherm [m]36 39 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hd28 !: depth of 28 C isotherm [m] 37 40 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: htc3 !: heat content of first 300 m [W] 38 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: htc7 !: heat content of first 700 m [W]39 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: htc20 !: heat content of first 2000 m [W]40 41 41 42 42 !!---------------------------------------------------------------------- … … 52 52 !!--------------------------------------------------------------------- 53 53 ! 54 ALLOCATE( hth(jpi,jpj), hd20(jpi,jpj), hd26(jpi,jpj), hd28(jpi,jpj), & 55 & htc3(jpi,jpj), htc7(jpi,jpj), htc20(jpi,jpj), STAT=dia_hth_alloc ) 54 ALLOCATE( hth(jpi,jpj), hd20(jpi,jpj), hd28(jpi,jpj), htc3(jpi,jpj), STAT=dia_hth_alloc ) 56 55 ! 57 56 CALL mpp_sum ( 'diahth', dia_hth_alloc ) … … 83 82 INTEGER, INTENT( in ) :: kt ! ocean time-step index 84 83 !! 85 INTEGER :: ji, jj, jk ! dummy loop arguments 86 REAL(wp) :: zrho3 = 0.03_wp ! density criterion for mixed layer depth 87 REAL(wp) :: zrho1 = 0.01_wp ! density criterion for mixed layer depth 88 REAL(wp) :: ztem2 = 0.2_wp ! temperature criterion for mixed layer depth 89 REAL(wp) :: zztmp, zzdep ! temporary scalars inside do loop 90 REAL(wp) :: zu, zv, zw, zut, zvt ! temporary workspace 91 REAL(wp), DIMENSION(jpi,jpj) :: zabs2 ! MLD: abs( tn - tn(10m) ) = ztem2 92 REAL(wp), DIMENSION(jpi,jpj) :: ztm2 ! Top of thermocline: tn = tn(10m) - ztem2 93 REAL(wp), DIMENSION(jpi,jpj) :: zrho10_3 ! MLD: rho = rho10m + zrho3 94 REAL(wp), DIMENSION(jpi,jpj) :: zpycn ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC) 95 REAL(wp), DIMENSION(jpi,jpj) :: ztinv ! max of temperature inversion 96 REAL(wp), DIMENSION(jpi,jpj) :: zdepinv ! depth of temperature inversion 97 REAL(wp), DIMENSION(jpi,jpj) :: zrho0_3 ! MLD rho = rho(surf) = 0.03 98 REAL(wp), DIMENSION(jpi,jpj) :: zrho0_1 ! MLD rho = rho(surf) = 0.01 99 REAL(wp), DIMENSION(jpi,jpj) :: zmaxdzT ! max of dT/dz 100 REAL(wp), DIMENSION(jpi,jpj) :: zdelr ! delta rho equivalent to deltaT = 0.2 84 INTEGER :: ji, jj, jk ! dummy loop arguments 85 INTEGER :: iid, ilevel ! temporary integers 86 INTEGER, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ik20, ik28 ! levels 87 REAL(wp) :: zavt5 = 5.e-4_wp ! Kz criterion for the turbocline depth 88 REAL(wp) :: zrho3 = 0.03_wp ! density criterion for mixed layer depth 89 REAL(wp) :: zrho1 = 0.01_wp ! density criterion for mixed layer depth 90 REAL(wp) :: ztem2 = 0.2_wp ! temperature criterion for mixed layer depth 91 REAL(wp) :: zthick_0, zcoef ! temporary scalars 92 REAL(wp) :: zztmp, zzdep ! temporary scalars inside do loop 93 REAL(wp) :: zu, zv, zw, zut, zvt ! temporary workspace 94 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zabs2 ! MLD: abs( tn - tn(10m) ) = ztem2 95 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ztm2 ! Top of thermocline: tn = tn(10m) - ztem2 96 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zrho10_3 ! MLD: rho = rho10m + zrho3 97 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zpycn ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC) 98 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ztinv ! max of temperature inversion 99 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zdepinv ! depth of temperature inversion 100 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zrho0_3 ! MLD rho = rho(surf) = 0.03 101 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zrho0_1 ! MLD rho = rho(surf) = 0.01 102 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zmaxdzT ! max of dT/dz 103 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zthick ! vertical integration thickness 104 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zdelr ! delta rho equivalent to deltaT = 0.2 101 105 !!---------------------------------------------------------------------- 102 106 IF( ln_timing ) CALL timing_start('dia_hth') 103 107 104 108 IF( kt == nit000 ) THEN 105 l_hth = .FALSE.106 IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) .OR. &107 & iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) .OR. iom_use( 'mldr10_3' ) .OR. &108 & iom_use( '20d' ) .OR. iom_use( '26d' ) .OR. iom_use( '28d' ) .OR. &109 & iom_use( 'hc300' ) .OR. iom_use( 'hc700' ) .OR. iom_use( 'hc2000' ) .OR. &110 & iom_use( 'pycndep' ) .OR. iom_use( 'tinv' ) .OR. iom_use( 'depti' ) ) l_hth = .TRUE.111 109 ! ! allocate dia_hth array 112 IF( l_hth ) THEN 113 IF( dia_hth_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'dia_hth : unable to allocate standard arrays' ) 114 IF(lwp) WRITE(numout,*) 115 IF(lwp) WRITE(numout,*) 'dia_hth : diagnostics of the thermocline depth' 116 IF(lwp) WRITE(numout,*) '~~~~~~~ ' 117 IF(lwp) WRITE(numout,*) 118 ENDIF 110 IF( dia_hth_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'dia_hth : unable to allocate standard arrays' ) 111 112 IF(.NOT. ALLOCATED(ik20) ) THEN 113 ALLOCATE(ik20(jpi,jpj), ik28(jpi,jpj), & 114 & zabs2(jpi,jpj), & 115 & ztm2(jpi,jpj), & 116 & zrho10_3(jpi,jpj),& 117 & zpycn(jpi,jpj), & 118 & ztinv(jpi,jpj), & 119 & zdepinv(jpi,jpj), & 120 & zrho0_3(jpi,jpj), & 121 & zrho0_1(jpi,jpj), & 122 & zmaxdzT(jpi,jpj), & 123 & zthick(jpi,jpj), & 124 & zdelr(jpi,jpj), STAT=ji) 125 CALL mpp_sum('diahth', ji) 126 IF( ji /= 0 ) CALL ctl_stop( 'STOP', 'dia_hth : unable to allocate standard ocean arrays' ) 127 END IF 128 129 IF(lwp) WRITE(numout,*) 130 IF(lwp) WRITE(numout,*) 'dia_hth : diagnostics of the thermocline depth' 131 IF(lwp) WRITE(numout,*) '~~~~~~~ ' 132 IF(lwp) WRITE(numout,*) 119 133 ENDIF 120 134 121 IF( l_hth ) THEN 122 ! 123 IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) ) THEN 124 ! initialization 125 ztinv (:,:) = 0._wp 126 zdepinv(:,:) = 0._wp 127 zmaxdzT(:,:) = 0._wp 128 DO jj = 1, jpj 129 DO ji = 1, jpi 130 zztmp = gdepw_n(ji,jj,mbkt(ji,jj)+1) 131 hth (ji,jj) = zztmp 132 zabs2 (ji,jj) = zztmp 133 ztm2 (ji,jj) = zztmp 134 zrho10_3(ji,jj) = zztmp 135 zpycn (ji,jj) = zztmp 136 END DO 135 ! initialization 136 ztinv (:,:) = 0._wp 137 zdepinv(:,:) = 0._wp 138 zmaxdzT(:,:) = 0._wp 139 DO jj = 1, jpj 140 DO ji = 1, jpi 141 zztmp = gdepw_n(ji,jj,mbkt(ji,jj)+1) 142 hth (ji,jj) = zztmp 143 zabs2 (ji,jj) = zztmp 144 ztm2 (ji,jj) = zztmp 145 zrho10_3(ji,jj) = zztmp 146 zpycn (ji,jj) = zztmp 147 END DO 148 END DO 149 IF( nla10 > 1 ) THEN 150 DO jj = 1, jpj 151 DO ji = 1, jpi 152 zztmp = gdepw_n(ji,jj,mbkt(ji,jj)+1) 153 zrho0_3(ji,jj) = zztmp 154 zrho0_1(ji,jj) = zztmp 137 155 END DO 138 IF( nla10 > 1 ) THEN 139 DO jj = 1, jpj 140 DO ji = 1, jpi 141 zztmp = gdepw_n(ji,jj,mbkt(ji,jj)+1) 142 zrho0_3(ji,jj) = zztmp 143 zrho0_1(ji,jj) = zztmp 144 END DO 145 END DO 156 END DO 157 ENDIF 158 159 ! Preliminary computation 160 ! computation of zdelr = (dr/dT)(T,S,10m)*(-0.2 degC) 161 DO jj = 1, jpj 162 DO ji = 1, jpi 163 IF( tmask(ji,jj,nla10) == 1. ) THEN 164 zu = 1779.50 + 11.250 * tsn(ji,jj,nla10,jp_tem) - 3.80 * tsn(ji,jj,nla10,jp_sal) & 165 & - 0.0745 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_tem) & 166 & - 0.0100 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_sal) 167 zv = 5891.00 + 38.000 * tsn(ji,jj,nla10,jp_tem) + 3.00 * tsn(ji,jj,nla10,jp_sal) & 168 & - 0.3750 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_tem) 169 zut = 11.25 - 0.149 * tsn(ji,jj,nla10,jp_tem) - 0.01 * tsn(ji,jj,nla10,jp_sal) 170 zvt = 38.00 - 0.750 * tsn(ji,jj,nla10,jp_tem) 171 zw = (zu + 0.698*zv) * (zu + 0.698*zv) 172 zdelr(ji,jj) = ztem2 * (1000.*(zut*zv - zvt*zu)/zw) 173 ELSE 174 zdelr(ji,jj) = 0._wp 146 175 ENDIF 176 END DO 177 END DO 178 179 ! ------------------------------------------------------------- ! 180 ! thermocline depth: strongest vertical gradient of temperature ! 181 ! turbocline depth (mixing layer depth): avt = zavt5 ! 182 ! MLD: rho = rho(1) + zrho3 ! 183 ! MLD: rho = rho(1) + zrho1 ! 184 ! ------------------------------------------------------------- ! 185 DO jk = jpkm1, 2, -1 ! loop from bottom to 2 186 DO jj = 1, jpj 187 DO ji = 1, jpi 188 ! 189 zzdep = gdepw_n(ji,jj,jk) 190 zztmp = ( tsn(ji,jj,jk-1,jp_tem) - tsn(ji,jj,jk,jp_tem) ) / zzdep * tmask(ji,jj,jk) ! vertical gradient of temperature (dT/dz) 191 zzdep = zzdep * tmask(ji,jj,1) 192 193 IF( zztmp > zmaxdzT(ji,jj) ) THEN 194 zmaxdzT(ji,jj) = zztmp ; hth (ji,jj) = zzdep ! max and depth of dT/dz 195 ENDIF 196 197 IF( nla10 > 1 ) THEN 198 zztmp = rhop(ji,jj,jk) - rhop(ji,jj,1) ! delta rho(1) 199 IF( zztmp > zrho3 ) zrho0_3(ji,jj) = zzdep ! > 0.03 200 IF( zztmp > zrho1 ) zrho0_1(ji,jj) = zzdep ! > 0.01 201 ENDIF 202 203 END DO 204 END DO 205 END DO 147 206 148 ! Preliminary computation 149 ! computation of zdelr = (dr/dT)(T,S,10m)*(-0.2 degC) 150 DO jj = 1, jpj 151 DO ji = 1, jpi 152 IF( tmask(ji,jj,nla10) == 1. ) THEN 153 zu = 1779.50 + 11.250 * tsn(ji,jj,nla10,jp_tem) - 3.80 * tsn(ji,jj,nla10,jp_sal) & 154 & - 0.0745 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_tem) & 155 & - 0.0100 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_sal) 156 zv = 5891.00 + 38.000 * tsn(ji,jj,nla10,jp_tem) + 3.00 * tsn(ji,jj,nla10,jp_sal) & 157 & - 0.3750 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_tem) 158 zut = 11.25 - 0.149 * tsn(ji,jj,nla10,jp_tem) - 0.01 * tsn(ji,jj,nla10,jp_sal) 159 zvt = 38.00 - 0.750 * tsn(ji,jj,nla10,jp_tem) 160 zw = (zu + 0.698*zv) * (zu + 0.698*zv) 161 zdelr(ji,jj) = ztem2 * (1000.*(zut*zv - zvt*zu)/zw) 162 ELSE 163 zdelr(ji,jj) = 0._wp 164 ENDIF 165 END DO 207 CALL iom_put( "mlddzt", hth ) ! depth of the thermocline 208 IF( nla10 > 1 ) THEN 209 CALL iom_put( "mldr0_3", zrho0_3 ) ! MLD delta rho(surf) = 0.03 210 CALL iom_put( "mldr0_1", zrho0_1 ) ! MLD delta rho(surf) = 0.01 211 ENDIF 212 213 ! ------------------------------------------------------------- ! 214 ! MLD: abs( tn - tn(10m) ) = ztem2 ! 215 ! Top of thermocline: tn = tn(10m) - ztem2 ! 216 ! MLD: rho = rho10m + zrho3 ! 217 ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC) ! 218 ! temperature inversion: max( 0, max of tn - tn(10m) ) ! 219 ! depth of temperature inversion ! 220 ! ------------------------------------------------------------- ! 221 DO jk = jpkm1, nlb10, -1 ! loop from bottom to nlb10 222 DO jj = 1, jpj 223 DO ji = 1, jpi 224 ! 225 zzdep = gdepw_n(ji,jj,jk) * tmask(ji,jj,1) 226 ! 227 zztmp = tsn(ji,jj,nla10,jp_tem) - tsn(ji,jj,jk,jp_tem) ! - delta T(10m) 228 IF( ABS(zztmp) > ztem2 ) zabs2 (ji,jj) = zzdep ! abs > 0.2 229 IF( zztmp > ztem2 ) ztm2 (ji,jj) = zzdep ! > 0.2 230 zztmp = -zztmp ! delta T(10m) 231 IF( zztmp > ztinv(ji,jj) ) THEN ! temperature inversion 232 ztinv(ji,jj) = zztmp ; zdepinv (ji,jj) = zzdep ! max value and depth 233 ENDIF 234 235 zztmp = rhop(ji,jj,jk) - rhop(ji,jj,nla10) ! delta rho(10m) 236 IF( zztmp > zrho3 ) zrho10_3(ji,jj) = zzdep ! > 0.03 237 IF( zztmp > zdelr(ji,jj) ) zpycn (ji,jj) = zzdep ! > equi. delta T(10m) - 0.2 238 ! 166 239 END DO 167 168 ! ------------------------------------------------------------- ! 169 ! thermocline depth: strongest vertical gradient of temperature ! 170 ! turbocline depth (mixing layer depth): avt = zavt5 ! 171 ! MLD: rho = rho(1) + zrho3 ! 172 ! MLD: rho = rho(1) + zrho1 ! 173 ! ------------------------------------------------------------- ! 174 DO jk = jpkm1, 2, -1 ! loop from bottom to 2 175 DO jj = 1, jpj 176 DO ji = 1, jpi 177 ! 178 zzdep = gdepw_n(ji,jj,jk) 179 zztmp = ( tsn(ji,jj,jk-1,jp_tem) - tsn(ji,jj,jk,jp_tem) ) & 180 & / zzdep * tmask(ji,jj,jk) ! vertical gradient of temperature (dT/dz) 181 zzdep = zzdep * tmask(ji,jj,1) 182 183 IF( zztmp > zmaxdzT(ji,jj) ) THEN 184 zmaxdzT(ji,jj) = zztmp 185 hth (ji,jj) = zzdep ! max and depth of dT/dz 186 ENDIF 187 188 IF( nla10 > 1 ) THEN 189 zztmp = rhop(ji,jj,jk) - rhop(ji,jj,1) ! delta rho(1) 190 IF( zztmp > zrho3 ) zrho0_3(ji,jj) = zzdep ! > 0.03 191 IF( zztmp > zrho1 ) zrho0_1(ji,jj) = zzdep ! > 0.01 192 ENDIF 193 END DO 194 END DO 195 END DO 196 197 CALL iom_put( 'mlddzt', hth ) ! depth of the thermocline 198 IF( nla10 > 1 ) THEN 199 CALL iom_put( 'mldr0_3', zrho0_3 ) ! MLD delta rho(surf) = 0.03 200 CALL iom_put( 'mldr0_1', zrho0_1 ) ! MLD delta rho(surf) = 0.01 201 ENDIF 202 ! 203 ENDIF 204 ! 205 IF( iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) .OR. iom_use( 'mldr10_3' ) .OR. & 206 & iom_use( 'pycndep' ) .OR. iom_use( 'tinv' ) .OR. iom_use( 'depti' ) ) THEN 207 ! ------------------------------------------------------------- ! 208 ! MLD: abs( tn - tn(10m) ) = ztem2 ! 209 ! Top of thermocline: tn = tn(10m) - ztem2 ! 210 ! MLD: rho = rho10m + zrho3 ! 211 ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC) ! 212 ! temperature inversion: max( 0, max of tn - tn(10m) ) ! 213 ! depth of temperature inversion ! 214 ! ------------------------------------------------------------- ! 215 DO jk = jpkm1, nlb10, -1 ! loop from bottom to nlb10 216 DO jj = 1, jpj 217 DO ji = 1, jpi 218 ! 219 zzdep = gdepw_n(ji,jj,jk) * tmask(ji,jj,1) 220 ! 221 zztmp = tsn(ji,jj,nla10,jp_tem) - tsn(ji,jj,jk,jp_tem) ! - delta T(10m) 222 IF( ABS(zztmp) > ztem2 ) zabs2 (ji,jj) = zzdep ! abs > 0.2 223 IF( zztmp > ztem2 ) ztm2 (ji,jj) = zzdep ! > 0.2 224 zztmp = -zztmp ! delta T(10m) 225 IF( zztmp > ztinv(ji,jj) ) THEN ! temperature inversion 226 ztinv(ji,jj) = zztmp 227 zdepinv (ji,jj) = zzdep ! max value and depth 228 ENDIF 229 230 zztmp = rhop(ji,jj,jk) - rhop(ji,jj,nla10) ! delta rho(10m) 231 IF( zztmp > zrho3 ) zrho10_3(ji,jj) = zzdep ! > 0.03 232 IF( zztmp > zdelr(ji,jj) ) zpycn (ji,jj) = zzdep ! > equi. delta T(10m) - 0.2 233 ! 234 END DO 235 END DO 236 END DO 237 238 CALL iom_put( 'mld_dt02', zabs2 ) ! MLD abs(delta t) - 0.2 239 CALL iom_put( 'topthdep', ztm2 ) ! T(10) - 0.2 240 CALL iom_put( 'mldr10_3', zrho10_3 ) ! MLD delta rho(10m) = 0.03 241 CALL iom_put( 'pycndep' , zpycn ) ! MLD delta rho equi. delta T(10m) = 0.2 242 CALL iom_put( 'tinv' , ztinv ) ! max. temp. inv. (t10 ref) 243 CALL iom_put( 'depti' , zdepinv ) ! depth of max. temp. inv. (t10 ref) 244 ! 245 ENDIF 246 247 ! ------------------------------- ! 248 ! Depth of 20C/26C/28C isotherm ! 249 ! ------------------------------- ! 250 IF( iom_use ('20d') ) THEN ! depth of the 20 isotherm 251 ztem2 = 20. 252 CALL dia_hth_dep( ztem2, hd20 ) 253 CALL iom_put( '20d', hd20 ) 254 ENDIF 255 ! 256 IF( iom_use ('26d') ) THEN ! depth of the 26 isotherm 257 ztem2 = 26. 258 CALL dia_hth_dep( ztem2, hd26 ) 259 CALL iom_put( '26d', hd26 ) 260 ENDIF 261 ! 262 IF( iom_use ('28d') ) THEN ! depth of the 28 isotherm 263 ztem2 = 28. 264 CALL dia_hth_dep( ztem2, hd28 ) 265 CALL iom_put( '28d', hd28 ) 266 ENDIF 267 268 ! ----------------------------- ! 269 ! Heat content of first 300 m ! 270 ! ----------------------------- ! 271 IF( iom_use ('hc300') ) THEN 272 zzdep = 300. 273 CALL dia_hth_htc( zzdep, tsn(:,:,:,jp_tem), htc3 ) 274 CALL iom_put( 'hc300', rau0_rcp * htc3 ) ! vertically integrated heat content (J/m2) 275 ENDIF 276 ! 277 ! ----------------------------- ! 278 ! Heat content of first 700 m ! 279 ! ----------------------------- ! 280 IF( iom_use ('hc700') ) THEN 281 zzdep = 700. 282 CALL dia_hth_htc( zzdep, tsn(:,:,:,jp_tem), htc7 ) 283 CALL iom_put( 'hc700', rau0_rcp * htc7 ) ! vertically integrated heat content (J/m2) 284 285 ENDIF 286 ! 287 ! ----------------------------- ! 288 ! Heat content of first 2000 m ! 289 ! ----------------------------- ! 290 IF( iom_use ('hc2000') ) THEN 291 zzdep = 2000. 292 CALL dia_hth_htc( zzdep, tsn(:,:,:,jp_tem), htc20 ) 293 CALL iom_put( 'hc2000', rau0_rcp * htc20 ) ! vertically integrated heat content (J/m2) 294 ENDIF 295 ! 296 ENDIF 297 298 ! 299 IF( ln_timing ) CALL timing_stop('dia_hth') 300 ! 301 END SUBROUTINE dia_hth 302 303 SUBROUTINE dia_hth_dep( ptem, pdept ) 304 ! 305 REAL(wp), INTENT(in) :: ptem 306 REAL(wp), DIMENSION(jpi,jpj), INTENT(out) :: pdept 307 ! 308 INTEGER :: ji, jj, jk, iid 309 REAL(wp) :: zztmp, zzdep 310 INTEGER, DIMENSION(jpi,jpj) :: iktem 311 312 ! --------------------------------------- ! 313 ! search deepest level above ptem ! 314 ! --------------------------------------- ! 315 iktem(:,:) = 1 240 END DO 241 END DO 242 243 CALL iom_put( "mld_dt02", zabs2 ) ! MLD abs(delta t) - 0.2 244 CALL iom_put( "topthdep", ztm2 ) ! T(10) - 0.2 245 CALL iom_put( "mldr10_3", zrho10_3 ) ! MLD delta rho(10m) = 0.03 246 CALL iom_put( "pycndep" , zpycn ) ! MLD delta rho equi. delta T(10m) = 0.2 247 CALL iom_put( "tinv" , ztinv ) ! max. temp. inv. (t10 ref) 248 CALL iom_put( "depti" , zdepinv ) ! depth of max. temp. inv. (t10 ref) 249 250 251 ! ----------------------------------- ! 252 ! search deepest level above 20C/28C ! 253 ! ----------------------------------- ! 254 ik20(:,:) = 1 255 ik28(:,:) = 1 316 256 DO jk = 1, jpkm1 ! beware temperature is not always decreasing with depth => loop from top to bottom 317 257 DO jj = 1, jpj 318 258 DO ji = 1, jpi 319 259 zztmp = tsn(ji,jj,jk,jp_tem) 320 IF( zztmp >= ptem ) iktem(ji,jj) = jk 260 IF( zztmp >= 20. ) ik20(ji,jj) = jk 261 IF( zztmp >= 28. ) ik28(ji,jj) = jk 321 262 END DO 322 263 END DO 323 264 END DO 324 265 325 ! --------------------------- ----!326 ! Depth of ptem isotherm!327 ! --------------------------- ----!266 ! --------------------------- ! 267 ! Depth of 20C/28C isotherm ! 268 ! --------------------------- ! 328 269 DO jj = 1, jpj 329 270 DO ji = 1, jpi 330 271 ! 331 zzdep = gdepw_n(ji,jj,mbkt(ji,jj)+1) ! depth of the o cean bottom272 zzdep = gdepw_n(ji,jj,mbkt(ji,jj)+1) ! depth of the oean bottom 332 273 ! 333 iid = ik tem(ji,jj)274 iid = ik20(ji,jj) 334 275 IF( iid /= 1 ) THEN 335 zztmp =gdept_n(ji,jj,iid ) & ! linear interpolation276 zztmp = gdept_n(ji,jj,iid ) & ! linear interpolation 336 277 & + ( gdept_n(ji,jj,iid+1) - gdept_n(ji,jj,iid) ) & 337 278 & * ( 20.*tmask(ji,jj,iid+1) - tsn(ji,jj,iid,jp_tem) ) & 338 279 & / ( tsn(ji,jj,iid+1,jp_tem) - tsn(ji,jj,iid,jp_tem) + (1.-tmask(ji,jj,1)) ) 339 pdept(ji,jj) = MIN( zztmp , zzdep) * tmask(ji,jj,1) ! bound by the ocean depth280 hd20(ji,jj) = MIN( zztmp , zzdep) * tmask(ji,jj,1) ! bound by the ocean depth 340 281 ELSE 341 pdept(ji,jj) = 0._wp282 hd20(ji,jj) = 0._wp 342 283 ENDIF 343 END DO 344 END DO 345 ! 346 END SUBROUTINE dia_hth_dep 347 348 349 SUBROUTINE dia_hth_htc( pdep, ptn, phtc ) 350 ! 351 REAL(wp), INTENT(in) :: pdep ! depth over the heat content 352 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: ptn 353 REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: phtc 354 ! 355 INTEGER :: ji, jj, jk, ik 356 REAL(wp), DIMENSION(jpi,jpj) :: zthick 357 INTEGER , DIMENSION(jpi,jpj) :: ilevel 358 359 284 ! 285 iid = ik28(ji,jj) 286 IF( iid /= 1 ) THEN 287 zztmp = gdept_n(ji,jj,iid ) & ! linear interpolation 288 & + ( gdept_n(ji,jj,iid+1) - gdept_n(ji,jj,iid) ) & 289 & * ( 28.*tmask(ji,jj,iid+1) - tsn(ji,jj,iid,jp_tem) ) & 290 & / ( tsn(ji,jj,iid+1,jp_tem) - tsn(ji,jj,iid,jp_tem) + (1.-tmask(ji,jj,1)) ) 291 hd28(ji,jj) = MIN( zztmp , zzdep ) * tmask(ji,jj,1) ! bound by the ocean depth 292 ELSE 293 hd28(ji,jj) = 0._wp 294 ENDIF 295 296 END DO 297 END DO 298 CALL iom_put( "20d", hd20 ) ! depth of the 20 isotherm 299 CALL iom_put( "28d", hd28 ) ! depth of the 28 isotherm 300 301 ! ----------------------------- ! 302 ! Heat content of first 300 m ! 303 ! ----------------------------- ! 304 305 ! find ilevel with (ilevel+1) the deepest W-level above 300m (we assume we can use e3t_1d to do this search...) 306 ilevel = 0 307 zthick_0 = 0._wp 308 DO jk = 1, jpkm1 309 zthick_0 = zthick_0 + e3t_1d(jk) 310 IF( zthick_0 < 300. ) ilevel = jk 311 END DO 360 312 ! surface boundary condition 361 362 IF( .NOT. ln_linssh ) THEN ; zthick(:,:) = 0._wp ; phtc(:,:) = 0._wp 363 ELSE ; zthick(:,:) = sshn(:,:) ; phtc(:,:) = ptn(:,:,1) * sshn(:,:) * tmask(:,:,1) 313 IF( ln_linssh ) THEN ; zthick(:,:) = sshn(:,:) ; htc3(:,:) = tsn(:,:,1,jp_tem) * sshn(:,:) * tmask(:,:,1) 314 ELSE ; zthick(:,:) = 0._wp ; htc3(:,:) = 0._wp 364 315 ENDIF 365 ! 366 ilevel(:,:) = 1 367 DO jk = 2, jpkm1 368 DO jj = 1, jpj 369 DO ji = 1, jpi 370 IF( ( gdept_n(ji,jj,jk) < pdep ) .AND. ( tmask(ji,jj,jk) == 1 ) ) THEN 371 ilevel(ji,jj) = jk 372 zthick(ji,jj) = zthick(ji,jj) + e3t_n(ji,jj,jk) 373 phtc (ji,jj) = phtc (ji,jj) + e3t_n(ji,jj,jk) * ptn(ji,jj,jk) 374 ENDIF 375 ENDDO 376 ENDDO 377 ENDDO 378 ! 316 ! integration down to ilevel 317 DO jk = 1, ilevel 318 zthick(:,:) = zthick(:,:) + e3t_n(:,:,jk) 319 htc3 (:,:) = htc3 (:,:) + e3t_n(:,:,jk) * tsn(:,:,jk,jp_tem) * tmask(:,:,jk) 320 END DO 321 ! deepest layer 322 zthick(:,:) = 300. - zthick(:,:) ! remaining thickness to reach 300m 379 323 DO jj = 1, jpj 380 324 DO ji = 1, jpi 381 ik = ilevel(ji,jj) 382 zthick(ji,jj) = pdep - zthick(ji,jj) ! remaining thickness to reach depht pdep 383 phtc(ji,jj) = phtc(ji,jj) + ptn(ji,jj,ik+1) * MIN( e3t_n(ji,jj,ik+1), zthick(ji,jj) ) & 384 * tmask(ji,jj,ik+1) 385 END DO 386 ENDDO 387 ! 388 ! 389 END SUBROUTINE dia_hth_htc 325 htc3(ji,jj) = htc3(ji,jj) + tsn(ji,jj,ilevel+1,jp_tem) & 326 & * MIN( e3t_n(ji,jj,ilevel+1), zthick(ji,jj) ) * tmask(ji,jj,ilevel+1) 327 END DO 328 END DO 329 ! from temperature to heat contain 330 zcoef = rau0 * rcp 331 htc3(:,:) = zcoef * htc3(:,:) 332 CALL iom_put( "hc300", htc3 ) ! first 300m heat content 333 ! 334 IF( ln_timing ) CALL timing_stop('dia_hth') 335 ! 336 END SUBROUTINE dia_hth 337 338 #else 339 !!---------------------------------------------------------------------- 340 !! Default option : Empty module 341 !!---------------------------------------------------------------------- 342 LOGICAL , PUBLIC, PARAMETER :: lk_diahth = .FALSE. !: thermocline-20d depths flag 343 CONTAINS 344 SUBROUTINE dia_hth( kt ) ! Empty routine 345 IMPLICIT NONE 346 INTEGER, INTENT( in ) :: kt 347 WRITE(*,*) 'dia_hth: You should not have seen this print! error?', kt 348 END SUBROUTINE dia_hth 349 #endif 390 350 391 351 !!======================================================================
Note: See TracChangeset
for help on using the changeset viewer.
