9 INTEGER,
PARAMETER ::
nz = 30
16 L_MOD, sfr_surf, FAI, PAI, &
17 StabilityMethod, RA_h, &
18 avcp, lv_J_kg, avdens, &
19 avU1, Temp_C, avRH, Press_hPa, zMeas, qh, qe, & ! input
20 T2_C, q2_gkg, U10_ms, RH2, & !output
35 REAL(KIND(1D0)),
DIMENSION(nsurf),
INTENT(in) :: sfr_surf
36 REAL(KIND(1D0)),
INTENT(in) :: zMeas
37 REAL(KIND(1D0)),
INTENT(in) :: avU1
38 REAL(KIND(1D0)),
INTENT(in) :: Temp_C
39 REAL(KIND(1D0)),
INTENT(in) :: avRH
40 REAL(KIND(1D0)),
INTENT(in) :: Press_hPa
41 REAL(KIND(1D0)),
INTENT(in) :: L_MOD
42 REAL(KIND(1D0)),
INTENT(in) :: RA_h
43 REAL(KIND(1D0)),
INTENT(in) :: avcp
44 REAL(KIND(1D0)),
INTENT(in) :: lv_J_kg
45 REAL(KIND(1D0)),
INTENT(in) :: avdens
46 REAL(KIND(1D0)),
INTENT(in) :: qh
47 REAL(KIND(1D0)),
INTENT(in) :: qe
48 REAL(KIND(1D0)),
INTENT(in) :: Zh
49 REAL(KIND(1D0)),
INTENT(in) :: z0m
50 REAL(KIND(1D0)),
INTENT(in) :: z0v
51 REAL(KIND(1D0)),
INTENT(in) :: zdm
52 REAL(KIND(1D0)),
INTENT(in) :: FAI
53 REAL(KIND(1D0)),
INTENT(in) :: PAI
57 INTEGER,
INTENT(in) :: StabilityMethod
58 INTEGER,
INTENT(in) :: DiagMethod
60 REAL(KIND(1D0)),
INTENT(out) :: T2_C
61 REAL(KIND(1D0)),
INTENT(out) :: q2_gkg
62 REAL(KIND(1D0)),
INTENT(out) :: U10_ms
63 REAL(KIND(1D0)),
INTENT(out) :: RH2
65 INTEGER,
PARAMETER :: nz = 30
67 REAL(KIND(1D0)),
PARAMETER :: cd_tree = 1.2
68 REAL(KIND(1D0)),
PARAMETER :: a_tree = 0.05
69 REAL(KIND(1D0)),
PARAMETER :: kappa = 0.40
71 REAL(KIND(1D0)),
PARAMETER :: beta_N = 0.40
72 REAL(KIND(1D0)),
PARAMETER :: pi = 4.*atan(1.0), r = 0.1
73 REAL(KIND(1D0)),
PARAMETER :: a1 = 4.,
a2 = -0.1,
a3 = 1.5, a4 = -1.
75 REAL(KIND(1D0)),
PARAMETER :: porosity_evetr = 0.32
81 REAL(KIND(1D0)),
INTENT(out),
DIMENSION(ncolumnsDataOutRSL - 5) :: dataoutLineRSL
82 REAL(KIND(1D0)),
DIMENSION(nz) :: psihatm_z
83 REAL(KIND(1D0)),
DIMENSION(nz) :: psihath_z
86 REAL(KIND(1D0)),
DIMENSION(nz) :: zarray
87 REAL(KIND(1D0)),
DIMENSION(nz) :: dataoutLineURSL
88 REAL(KIND(1D0)),
DIMENSION(nz) :: dataoutLineTRSL
89 REAL(KIND(1D0)),
DIMENSION(nz) :: dataoutLineqRSL
91 REAL(KIND(1D0)) :: z0_RSL
92 REAL(KIND(1D0)) :: zd_RSL
98 REAL(KIND(1D0)) :: Scc
99 REAL(KIND(1D0)) :: psimz, psimz0, psimza, psihz, psihz0, psihza
101 REAL(KIND(1D0)) :: beta
102 REAL(KIND(1D0)) :: elm
104 REAL(KIND(1D0)) :: fx
105 REAL(KIND(1D0)) :: t_h, q_h
106 REAL(KIND(1D0)) :: TStar_RSL
107 REAL(KIND(1D0)) :: UStar_RSL
108 REAL(KIND(1D0)) :: UStar_heat
111 REAL(KIND(1D0)) :: L_MOD_RSL
115 REAL(KIND(1D0)) :: zH_RSL
116 REAL(KIND(1D0)) :: dz_above
117 REAL(KIND(1D0)) :: dz_can
122 REAL(KIND(1D0)) :: qa_gkg, qStar_RSL
148 IF (diagmethod == 0)
THEN
151 ELSEIF (diagmethod == 1)
THEN
154 ELSEIF (diagmethod == 2)
THEN
162 pai > 0.1 .AND. pai < 0.68 .AND. &
178 ELSE IF (zh <= 10)
THEN
190 IF (dz_can > 2) zarray(1) = 1.999
193 nz_above = nz - nz_can
194 dz_above = (zmeas - zh)/nz_above
195 DO i = nz_can + 1, nz
196 zarray(i) = zh + (i - nz_can)*dz_above
216 psihatm_z = 0.*zarray
217 psihath_z = 0.*zarray
221 nz_above, zarray(nz_can + 1:nz), &
222 zh, l_mod, sfr_surf, fai, pai, &
223 psihatm_z(nz_can + 1:nz), psihath_z(nz_can + 1:nz), &
225 lc, beta, zd_rsl, z0_rsl, elm, scc, fx)
276 psimz0 =
stab_psi_mom(stabilitymethod, z0_rsl/l_mod_rsl)
277 psimza =
stab_psi_mom(stabilitymethod, (zmeas - zd_rsl)/l_mod_rsl)
278 psihza =
stab_psi_heat(stabilitymethod, (zmeas - zd_rsl)/l_mod_rsl)
280 ustar_rsl = avu1*kappa/(log((zmeas - zd_rsl)/z0_rsl) - psimza + psimz0 + psihatm_z(nz))
284 ustar_rsl = max(0.001, ustar_rsl)
286 IF ((zmeas - zd_rsl)/l_mod_rsl < -
neut_limit) ustar_rsl = max(0.15, ustar_rsl)
291 ustar_heat = max(0.15, ustar_rsl)
295 ustar_heat = 1/(kappa*ra_h)*(log((zmeas - zd_rsl)/z0v) - psihza + psihz0)
297 tstar_rsl = -1.*(qh/(avcp*avdens))/ustar_heat
299 qstar_rsl = 10.**(-10)
301 qstar_rsl = -1.*(qe/lv_j_kg*avdens)/ustar_heat
303 qa_gkg =
rh2qa(avrh/100, press_hpa, temp_c)
306 psimz =
stab_psi_mom(stabilitymethod, (zarray(z) - zd_rsl)/l_mod_rsl)
307 psihz =
stab_psi_heat(stabilitymethod, (zarray(z) - zd_rsl)/l_mod_rsl)
308 dataoutlineursl(z) = (log((zarray(z) - zd_rsl)/z0_rsl) - psimz + psimz0 + psihatm_z(z))/kappa
309 dataoutlinetrsl(z) = (log((zarray(z) - zd_rsl)/(zmeas - zd_rsl)) - psihz + psihza + psihath_z(z) - psihath_z(nz))/kappa
310 dataoutlineqrsl(z) = dataoutlinetrsl(z)
318 t_h = scc*tstar_rsl/(beta*fx)
319 q_h = scc*qstar_rsl/(beta*fx)
321 dataoutlineursl(z) = dataoutlineursl(nz_can)*exp(beta*(zarray(z) - zh_rsl)/elm)
322 dataoutlinetrsl(z) = dataoutlinetrsl(nz_can) + (t_h*exp(beta*fx*(zarray(z) - zh_rsl)/elm) - t_h)/tstar_rsl
323 dataoutlineqrsl(z) = dataoutlineqrsl(nz_can) + (q_h*exp(beta*fx*(zarray(z) - zh_rsl)/elm) - q_h)/qstar_rsl
330 IF (zarray(z) <= zd_rsl) zarray(z) = 1.01*(zd_rsl + z0_rsl)
331 psimz =
stab_psi_mom(stabilitymethod, (zarray(z) - zd_rsl)/l_mod_rsl)
332 psihz =
stab_psi_heat(stabilitymethod, (zarray(z) - zd_rsl)/l_mod_rsl)
333 dataoutlineursl(z) = (log((zarray(z) - zd_rsl)/z0_rsl) - psimz + psimz0)/kappa
334 dataoutlinetrsl(z) = (log((zarray(z) - zd_rsl)/(zmeas - zd_rsl)) - psihz + psihza)/kappa
335 dataoutlineqrsl(z) = dataoutlinetrsl(z)
340 dataoutlineursl = dataoutlineursl*ustar_rsl
341 dataoutlinetrsl = dataoutlinetrsl*tstar_rsl + temp_c
342 dataoutlineqrsl = (dataoutlineqrsl*qstar_rsl + qa_gkg/1000.)*1000.
345 dataoutlinersl = [zarray, dataoutlineursl, dataoutlinetrsl, dataoutlineqrsl, &
352 beta, zd_rsl, z0_rsl, elm, scc, fx, &
353 ustar_rsl, ustar_heat, tstar_rsl, fai, pai, merge(1.d0, 0.d0, flag_rsl) &
362 t2_c =
interp_z(2d0, zarray, dataoutlinetrsl)
363 q2_gkg =
interp_z(2d0, zarray, dataoutlineqrsl)
364 u10_ms =
interp_z(10d0, zarray, dataoutlineursl)
367 t2_c =
interp_z(2d0 + zd_rsl + z0_rsl, zarray, dataoutlinetrsl)
368 q2_gkg =
interp_z(2d0 + zd_rsl + z0_rsl, zarray, dataoutlineqrsl)
369 u10_ms =
interp_z(10d0 + zd_rsl + z0_rsl, zarray, dataoutlineursl)
372 rh2 =
qa2rh(q2_gkg, press_hpa, t2_c)
380 sfr_paved, sfr_bldg, sfr_evetr, sfr_dectr, sfr_grass, sfr_bsoil, sfr_water, &
382 StabilityMethod, RA_h, &
383 avcp, lv_J_kg, avdens, &
384 avU1, Temp_C, avRH, Press_hPa, zMeas, qh, qe, & ! input
385 T2_C, q2_gkg, U10_ms, RH2, & !output
400 REAL(KIND(1D0)),
INTENT(IN) :: sfr_paved
401 REAL(KIND(1D0)),
INTENT(IN) :: sfr_bldg
402 REAL(KIND(1D0)),
INTENT(IN) :: sfr_evetr
403 REAL(KIND(1D0)),
INTENT(IN) :: sfr_dectr
404 REAL(KIND(1D0)),
INTENT(IN) :: sfr_grass
405 REAL(KIND(1D0)),
INTENT(IN) :: sfr_bsoil
406 REAL(KIND(1D0)),
INTENT(IN) :: sfr_water
407 REAL(KIND(1D0)),
DIMENSION(nsurf) :: sfr_surf
409 REAL(KIND(1D0)),
INTENT(in) :: zMeas
410 REAL(KIND(1D0)),
INTENT(in) :: avU1
411 REAL(KIND(1D0)),
INTENT(in) :: Temp_C
412 REAL(KIND(1D0)),
INTENT(in) :: avRH
413 REAL(KIND(1D0)),
INTENT(in) :: Press_hPa
414 REAL(KIND(1D0)),
INTENT(in) :: L_MOD
415 REAL(KIND(1D0)),
INTENT(in) :: RA_h
416 REAL(KIND(1D0)),
INTENT(in) :: avcp
417 REAL(KIND(1D0)),
INTENT(in) :: lv_J_kg
418 REAL(KIND(1D0)),
INTENT(in) :: avdens
419 REAL(KIND(1D0)),
INTENT(in) :: qh
420 REAL(KIND(1D0)),
INTENT(in) :: qe
421 REAL(KIND(1D0)),
INTENT(in) :: Zh
422 REAL(KIND(1D0)),
INTENT(in) :: z0m
423 REAL(KIND(1D0)),
INTENT(in) :: z0v
424 REAL(KIND(1D0)),
INTENT(in) :: zdm
425 REAL(KIND(1D0)),
INTENT(in) :: FAI
426 REAL(KIND(1D0)),
INTENT(in) :: PAI
430 INTEGER,
INTENT(in) :: StabilityMethod
431 INTEGER,
INTENT(in) :: DiagMethod
433 REAL(KIND(1D0)),
INTENT(out) :: T2_C
434 REAL(KIND(1D0)),
INTENT(out) :: q2_gkg
435 REAL(KIND(1D0)),
INTENT(out) :: U10_ms
436 REAL(KIND(1D0)),
INTENT(out) :: RH2
438 INTEGER,
PARAMETER :: nz = 30
440 REAL(KIND(1D0)),
PARAMETER :: cd_tree = 1.2
441 REAL(KIND(1D0)),
PARAMETER :: a_tree = 0.05
442 REAL(KIND(1D0)),
PARAMETER :: kappa = 0.40
444 REAL(KIND(1D0)),
PARAMETER :: beta_N = 0.40
445 REAL(KIND(1D0)),
PARAMETER :: pi = 4.*atan(1.0), r = 0.1
446 REAL(KIND(1D0)),
PARAMETER :: a1 = 4.,
a2 = -0.1,
a3 = 1.5, a4 = -1.
448 REAL(KIND(1D0)),
PARAMETER :: porosity_evetr = 0.32
454 REAL(KIND(1D0)),
INTENT(out),
DIMENSION(ncolumnsDataOutRSL - 5) :: dataoutLineRSL
455 REAL(KIND(1D0)),
DIMENSION(nz) :: psihatm_z
456 REAL(KIND(1D0)),
DIMENSION(nz) :: psihath_z
459 REAL(KIND(1D0)),
DIMENSION(nz) :: zarray
460 REAL(KIND(1D0)),
DIMENSION(nz) :: dataoutLineURSL
461 REAL(KIND(1D0)),
DIMENSION(nz) :: dataoutLineTRSL
462 REAL(KIND(1D0)),
DIMENSION(nz) :: dataoutLineqRSL
464 REAL(KIND(1D0)) :: z0_RSL
465 REAL(KIND(1D0)) :: zd_RSL
468 REAL(KIND(1D0)) :: Lc
471 REAL(KIND(1D0)) :: Scc
472 REAL(KIND(1D0)) :: psimz, psimz0, psimza, psihz, psihz0, psihza
474 REAL(KIND(1D0)) :: beta
475 REAL(KIND(1D0)) :: elm
477 REAL(KIND(1D0)) :: fx
478 REAL(KIND(1D0)) :: t_h, q_h
479 REAL(KIND(1D0)) :: TStar_RSL
480 REAL(KIND(1D0)) :: UStar_RSL
481 REAL(KIND(1D0)) :: UStar_heat
484 REAL(KIND(1D0)) :: L_MOD_RSL
488 REAL(KIND(1D0)) :: zH_RSL
489 REAL(KIND(1D0)) :: dz_above
490 REAL(KIND(1D0)) :: dz_can
495 REAL(KIND(1D0)) :: qa_gkg, qStar_RSL
520 sfr_surf = [sfr_paved, sfr_bldg, sfr_evetr, sfr_dectr, sfr_grass, sfr_bsoil, sfr_water]
522 IF (diagmethod == 0)
THEN
525 ELSEIF (diagmethod == 1)
THEN
528 ELSEIF (diagmethod == 2)
THEN
535 pai > 0.1 .AND. pai < 0.68 .AND. &
551 ELSE IF (zh <= 10)
THEN
563 IF (dz_can > 2) zarray(1) = 1.999
566 nz_above = nz - nz_can
567 dz_above = (zmeas - zh)/nz_above
568 DO i = nz_can + 1, nz
569 zarray(i) = zh + (i - nz_can)*dz_above
589 psihatm_z = 0.*zarray
590 psihath_z = 0.*zarray
594 nz_above, zarray(nz_can + 1:nz), &
595 zh, l_mod, sfr_surf, fai, pai, &
596 psihatm_z(nz_can + 1:nz), psihath_z(nz_can + 1:nz), &
598 lc, beta, zd_rsl, z0_rsl, elm, scc, fx)
649 psimz0 =
stab_psi_mom(stabilitymethod, z0_rsl/l_mod_rsl)
650 psimza =
stab_psi_mom(stabilitymethod, (zmeas - zd_rsl)/l_mod_rsl)
651 psihza =
stab_psi_heat(stabilitymethod, (zmeas - zd_rsl)/l_mod_rsl)
653 ustar_rsl = avu1*kappa/(log((zmeas - zd_rsl)/z0_rsl) - psimza + psimz0 + psihatm_z(nz))
657 ustar_rsl = max(0.001, ustar_rsl)
659 IF ((zmeas - zd_rsl)/l_mod_rsl < -
neut_limit) ustar_rsl = max(0.15, ustar_rsl)
664 ustar_heat = max(0.15, ustar_rsl)
668 ustar_heat = 1/(kappa*ra_h)*(log((zmeas - zd_rsl)/z0v) - psihza + psihz0)
670 tstar_rsl = -1.*(qh/(avcp*avdens))/ustar_heat
672 qstar_rsl = 10.**(-10)
674 qstar_rsl = -1.*(qe/lv_j_kg*avdens)/ustar_heat
676 qa_gkg =
rh2qa(avrh/100, press_hpa, temp_c)
679 psimz =
stab_psi_mom(stabilitymethod, (zarray(z) - zd_rsl)/l_mod_rsl)
680 psihz =
stab_psi_heat(stabilitymethod, (zarray(z) - zd_rsl)/l_mod_rsl)
681 dataoutlineursl(z) = (log((zarray(z) - zd_rsl)/z0_rsl) - psimz + psimz0 + psihatm_z(z))/kappa
682 dataoutlinetrsl(z) = (log((zarray(z) - zd_rsl)/(zmeas - zd_rsl)) - psihz + psihza + psihath_z(z) - psihath_z(nz))/kappa
683 dataoutlineqrsl(z) = dataoutlinetrsl(z)
691 t_h = scc*tstar_rsl/(beta*fx)
692 q_h = scc*qstar_rsl/(beta*fx)
694 dataoutlineursl(z) = dataoutlineursl(nz_can)*exp(beta*(zarray(z) - zh_rsl)/elm)
695 dataoutlinetrsl(z) = dataoutlinetrsl(nz_can) + (t_h*exp(beta*fx*(zarray(z) - zh_rsl)/elm) - t_h)/tstar_rsl
696 dataoutlineqrsl(z) = dataoutlineqrsl(nz_can) + (q_h*exp(beta*fx*(zarray(z) - zh_rsl)/elm) - q_h)/qstar_rsl
703 IF (zarray(z) <= zd_rsl) zarray(z) = 1.01*(zd_rsl + z0_rsl)
704 psimz =
stab_psi_mom(stabilitymethod, (zarray(z) - zd_rsl)/l_mod_rsl)
705 psihz =
stab_psi_heat(stabilitymethod, (zarray(z) - zd_rsl)/l_mod_rsl)
706 dataoutlineursl(z) = (log((zarray(z) - zd_rsl)/z0_rsl) - psimz + psimz0)/kappa
707 dataoutlinetrsl(z) = (log((zarray(z) - zd_rsl)/(zmeas - zd_rsl)) - psihz + psihza)/kappa
708 dataoutlineqrsl(z) = dataoutlinetrsl(z)
713 dataoutlineursl = dataoutlineursl*ustar_rsl
714 dataoutlinetrsl = dataoutlinetrsl*tstar_rsl + temp_c
715 dataoutlineqrsl = (dataoutlineqrsl*qstar_rsl + qa_gkg/1000.)*1000.
718 dataoutlinersl = [zarray, dataoutlineursl, dataoutlinetrsl, dataoutlineqrsl, &
725 beta, zd_rsl, z0_rsl, elm, scc, fx, &
726 ustar_rsl, ustar_heat, tstar_rsl, fai, pai, merge(1.d0, 0.d0, flag_rsl) &
735 t2_c =
interp_z(2d0, zarray, dataoutlinetrsl)
736 q2_gkg =
interp_z(2d0, zarray, dataoutlineqrsl)
737 u10_ms =
interp_z(10d0, zarray, dataoutlineursl)
740 t2_c =
interp_z(2d0 + zd_rsl + z0_rsl, zarray, dataoutlinetrsl)
741 q2_gkg =
interp_z(2d0 + zd_rsl + z0_rsl, zarray, dataoutlineqrsl)
742 u10_ms =
interp_z(10d0 + zd_rsl + z0_rsl, zarray, dataoutlineursl)
745 rh2 =
qa2rh(q2_gkg, press_hpa, t2_c)
751 REAL(kind(1d0)),
INTENT(in) :: z_x
752 REAL(kind(1d0)),
DIMENSION(nz),
INTENT(in) :: z
753 REAL(kind(1d0)),
DIMENSION(nz),
INTENT(in) :: v
756 REAL(kind(1d0)) :: v_x
759 REAL(kind(1d0)) :: slope
760 REAL(kind(1d0)) :: dz
761 REAL(kind(1d0)),
DIMENSION(nz) :: dif
766 INTEGER,
PARAMETER ::
nz = 30
772 idx_x = maxloc(dif, 1, abs(dif) < 1.d-6)
773 idx_low = maxloc(dif, 1, dif < 0.)
774 idx_high = minloc(dif, 1, dif > 0.)
781 dz = z(idx_high) - z(idx_low)
782 slope = (v(idx_high) - v(idx_low))/dz
783 v_x = v(idx_low) + (z_x - z(idx_low))*slope
790 REAL(kind(1d0)),
INTENT(in) :: lc
791 REAL(kind(1d0)),
INTENT(in) :: beta
794 REAL(kind(1d0)) :: elm
801 psihatm_top, psihatm_mid, &
802 z_top, z_mid, z_btm, &
804 zh_RSL, zd_RSL, L_MOD, beta, elm, Lc) &
809 INTEGER,
INTENT(in) :: stabilitymethod
810 REAL(kind(1d0)),
INTENT(in) :: psihatm_top
811 REAL(kind(1d0)),
INTENT(in) :: z_top
812 REAL(kind(1d0)),
INTENT(in) :: psihatm_mid
813 REAL(kind(1d0)),
INTENT(in) :: z_mid
814 REAL(kind(1d0)),
INTENT(in) :: z_btm
815 REAL(kind(1d0)),
INTENT(in) :: zh_rsl
816 REAL(kind(1d0)),
INTENT(in) :: zd_rsl
817 REAL(kind(1d0)),
INTENT(in) :: lc
818 REAL(kind(1d0)),
INTENT(in) :: beta
819 REAL(kind(1d0)),
INTENT(in) :: elm
820 REAL(kind(1d0)),
INTENT(in) :: l_mod
821 REAL(kind(1d0)),
INTENT(in) :: cm
822 REAL(kind(1d0)),
INTENT(in) :: c2
825 REAL(kind(1d0)) :: psihatm_btm
828 REAL(kind(1d0)) :: phim_btm
829 REAL(kind(1d0)) :: phim_mid
830 REAL(kind(1d0)) :: slope_top
831 REAL(kind(1d0)) :: slope_btm
832 REAL(kind(1d0)) :: z_plus
833 REAL(kind(1d0)) :: psihatm_plus
835 REAL(kind(1d0)),
PARAMETER :: tol = 0.5
836 REAL(kind(1d0)),
PARAMETER :: kappa = 0.40
837 REAL(kind(1d0)) :: dz_above
839 dz_above = z_mid - z_btm
841 phim_mid =
stab_phi_mom(stabilitymethod, (z_mid - zd_rsl)/l_mod)
842 phim_btm =
stab_phi_mom(stabilitymethod, (z_btm - zd_rsl)/l_mod)
844 psihatm_btm = psihatm_mid + dz_above/2.*phim_mid*(cm*exp(-1.*c2*beta*(z_mid - zd_rsl)/elm)) &
846 psihatm_btm = psihatm_btm + dz_above/2.*phim_btm*(cm*exp(-1.*c2*beta*(z_btm - zd_rsl)/elm)) &
848 IF (dz_above/zd_rsl < 1e-2)
THEN
852 IF (abs(psihatm_btm) > 1e-3)
THEN
854 slope_top = (psihatm_top - psihatm_mid)/(z_top - z_mid)
855 slope_btm = (psihatm_mid - psihatm_btm)/(z_mid - z_btm)
856 IF (abs(slope_btm) > (1 + tol)*abs(slope_top))
THEN
857 z_plus = (z_mid + z_btm)/2
860 psihatm_top, psihatm_mid, &
861 z_top, z_mid, z_plus, &
863 zh_rsl, zd_rsl, l_mod, beta, elm, lc)
866 psihatm_mid, psihatm_plus, &
867 z_mid, z_plus, z_btm, &
869 zh_rsl, zd_rsl, l_mod, beta, elm, lc)
876 psihath_top, psihath_mid, &
877 z_top, z_mid, z_btm, &
879 zh_RSL, zd_RSL, L_MOD, beta, elm, Lc) &
884 INTEGER,
INTENT(in) :: stabilitymethod
885 REAL(kind(1d0)),
INTENT(in) :: psihath_top
886 REAL(kind(1d0)),
INTENT(in) :: z_top
887 REAL(kind(1d0)),
INTENT(in) :: psihath_mid
888 REAL(kind(1d0)),
INTENT(in) :: z_mid
889 REAL(kind(1d0)),
INTENT(in) :: z_btm
890 REAL(kind(1d0)),
INTENT(in) :: zh_rsl
891 REAL(kind(1d0)),
INTENT(in) :: zd_rsl
892 REAL(kind(1d0)),
INTENT(in) :: lc
893 REAL(kind(1d0)),
INTENT(in) :: beta
894 REAL(kind(1d0)),
INTENT(in) :: elm
895 REAL(kind(1d0)),
INTENT(in) :: l_mod
896 REAL(kind(1d0)),
INTENT(in) :: ch
897 REAL(kind(1d0)),
INTENT(in) :: c2h
900 REAL(kind(1d0)) :: psihath_btm
903 REAL(kind(1d0)) :: phih_btm
904 REAL(kind(1d0)) :: phih_mid
906 REAL(kind(1d0)) :: slope_top
907 REAL(kind(1d0)) :: slope_btm
908 REAL(kind(1d0)) :: z_plus
909 REAL(kind(1d0)) :: psihath_plus
911 REAL(kind(1d0)),
PARAMETER :: tol = 0.1
912 REAL(kind(1d0)),
PARAMETER :: kappa = 0.40
913 REAL(kind(1d0)) :: dz_above
915 dz_above = z_mid - z_btm
917 phih_mid =
stab_phi_heat(stabilitymethod, (z_mid - zd_rsl)/l_mod)
918 phih_btm =
stab_phi_heat(stabilitymethod, (z_btm - zd_rsl)/l_mod)
920 psihath_btm = psihath_mid + dz_above/2.*phih_mid*(ch*exp(-1.*c2h*beta*(z_mid - zd_rsl)/elm)) &
922 psihath_btm = psihath_btm + dz_above/2.*phih_btm*(ch*exp(-1.*c2h*beta*(z_btm - zd_rsl)/elm)) &
924 IF (dz_above/zd_rsl < 1e-2)
THEN
928 IF (abs(psihath_btm) > 1e-3)
THEN
930 slope_top = (psihath_top - psihath_mid)/(z_top - z_mid)
931 slope_btm = (psihath_mid - psihath_btm)/(z_mid - z_btm)
932 IF (abs(slope_btm) > (1 + tol)*abs(slope_top))
THEN
933 z_plus = (z_mid + z_btm)/2
936 psihath_top, psihath_mid, &
937 z_top, z_mid, z_plus, &
939 zh_rsl, zd_rsl, l_mod, beta, elm, lc)
942 psihath_mid, psihath_plus, &
943 z_mid, z_plus, z_btm, &
945 zh_rsl, zd_rsl, l_mod, beta, elm, lc)
951 FUNCTION cal_phim_hat(StabilityMethod, z, zh_RSL, L_MOD, beta, lc)
RESULT(phim_hat)
953 INTEGER,
INTENT(in) :: stabilitymethod
954 REAL(kind(1d0)),
INTENT(in) :: z
955 REAL(kind(1d0)),
INTENT(in) :: zh_rsl
956 REAL(kind(1d0)),
INTENT(in) :: l_mod
957 REAL(kind(1d0)),
INTENT(in) :: beta
958 REAL(kind(1d0)),
INTENT(in) :: lc
959 REAL(kind(1d0)) :: phim_hat
960 REAL(kind(1d0)) :: zd_rsl
962 REAL(kind(1d0)) :: elm
963 REAL(kind(1d0)) :: c2
964 REAL(kind(1d0)) :: cm
970 CALL cal_cm(stabilitymethod, zh_rsl, zd_rsl, lc, beta, l_mod, c2, cm)
972 phim_hat = 1 - cm*exp(-1.*c2*beta*(z - zd_rsl)/elm)
976 SUBROUTINE cal_cm(StabilityMethod, zh_RSL, zd_RSL, Lc, beta, L_MOD, c2, cm)
979 INTEGER,
INTENT(in) :: StabilityMethod
981 REAL(KIND(1D0)),
INTENT(in) :: zh_RSL
982 REAL(KIND(1D0)),
INTENT(in) :: zd_RSL
983 REAL(KIND(1D0)),
INTENT(in) :: Lc
984 REAL(KIND(1D0)),
INTENT(in) :: beta
985 REAL(KIND(1D0)),
INTENT(in) :: L_MOD
988 REAL(KIND(1D0)),
INTENT(out) :: c2
989 REAL(KIND(1D0)),
INTENT(out) :: cm
993 REAL(KIND(1D0)) :: phi_hatmZh
994 REAL(KIND(1D0)) :: phim_zh
995 REAL(KIND(1D0)) :: phim_zhdz
996 REAL(KIND(1D0)) :: dphi
999 REAL(KIND(1D0)),
PARAMETER :: kappa = 0.40
1000 REAL(KIND(1D0)),
PARAMETER :: dz = 0.1
1002 phim_zh =
stab_phi_mom(stabilitymethod, (zh_rsl - zd_rsl)/l_mod)
1003 phim_zhdz =
stab_phi_mom(stabilitymethod, (zh_rsl - zd_rsl + dz)/l_mod)
1005 dphi = (phim_zhdz - phim_zh)/dz
1006 IF (phim_zh /= 0.)
THEN
1007 phi_hatmzh = kappa/(2.*beta*phim_zh)
1013 IF (phi_hatmzh >= 1.)
THEN
1019 c2 = (kappa*(3.-(2.*beta**2.*lc/phim_zh*dphi)))/(2.*beta*phim_zh - kappa)
1025 cm = (1.-phi_hatmzh)*exp(c2/2.)
1029 SUBROUTINE cal_ch(StabilityMethod, zh_RSL, zd_RSL, Lc, beta, L_MOD, Scc, f, c2h, ch)
1032 INTEGER,
INTENT(in) :: StabilityMethod
1034 REAL(KIND(1D0)),
INTENT(in) :: zh_RSL
1035 REAL(KIND(1D0)),
INTENT(in) :: zd_RSL
1036 REAL(KIND(1D0)),
INTENT(in) :: Scc
1037 REAL(KIND(1D0)),
INTENT(in) :: f
1038 REAL(KIND(1D0)),
INTENT(in) :: Lc
1039 REAL(KIND(1D0)),
INTENT(in) :: beta
1040 REAL(KIND(1D0)),
INTENT(in) :: L_MOD
1043 REAL(KIND(1D0)),
INTENT(out) :: ch
1044 REAL(KIND(1D0)),
INTENT(out) :: c2h
1047 REAL(KIND(1D0)) :: phih_zh
1048 REAL(KIND(1D0)) :: phih_zhdz
1049 REAL(KIND(1D0)) :: dphih
1050 REAL(KIND(1D0)) :: phi_hathZh
1052 REAL(KIND(1D0)),
PARAMETER :: kappa = 0.40
1053 REAL(KIND(1D0)),
PARAMETER :: dz = 0.1
1055 phih_zh =
stab_phi_heat(stabilitymethod, (zh_rsl - zd_rsl)/l_mod)
1056 phih_zhdz =
stab_phi_heat(stabilitymethod, (zh_rsl - zd_rsl + 1.)/l_mod)
1058 dphih = phih_zhdz - phih_zh
1059 IF (phih_zh /= 0.)
THEN
1060 phi_hathzh = kappa*scc/(2.*beta*phih_zh)
1065 IF (phi_hathzh >= 1.)
THEN
1071 c2h = (kappa*scc*(2.+f - (dphih*2.*beta**2.*lc/phih_zh)))/(2.*beta*phih_zh - kappa*scc)
1077 ch = (1.-phi_hathzh)*exp(c2h/2.)
1238 REAL(kind(1d0)),
INTENT(in) :: zh_rsl
1239 REAL(kind(1d0)),
INTENT(in) :: lc
1240 REAL(kind(1d0)),
INTENT(in) :: beta
1243 REAL(kind(1d0)) :: zd_rsl
1245 zd_rsl = zh_rsl - (beta**2.)*lc
1251 FUNCTION cal_z0_rsl(StabilityMethod, zH_RSL, zd_RSL, beta, L_MOD_RSL, psihatm_Zh)
RESULT(z0_RSL)
1255 INTEGER,
INTENT(in) :: stabilitymethod
1256 REAL(kind(1d0)),
INTENT(in) :: zh_rsl
1257 REAL(kind(1d0)),
INTENT(in) :: zd_rsl
1258 REAL(kind(1d0)),
INTENT(in) :: l_mod_rsl
1260 REAL(kind(1d0)),
INTENT(in) :: beta
1261 REAL(kind(1d0)),
INTENT(in) :: psihatm_zh
1264 REAL(kind(1d0)) :: z0_rsl
1267 REAL(kind(1d0)) :: psimzh, psimz0, z0_rsl_x
1268 REAL(kind(1d0)) :: err
1271 REAL(kind(1d0)),
PARAMETER :: kappa = 0.40
1275 psimzh =
stab_psi_mom(stabilitymethod, (zh_rsl - zd_rsl)/l_mod_rsl)
1283 DO WHILE ((err > 0.001) .AND. (it < 10))
1285 psimz0 =
stab_psi_mom(stabilitymethod, z0_rsl_x/l_mod_rsl)
1286 z0_rsl = (zh_rsl - zd_rsl)*exp(-1.*kappa/beta)*exp(-1.*psimzh + psimz0)*exp(psihatm_zh)
1287 err = abs(z0_rsl_x - z0_rsl)
1292 z0_rsl = merge(z0_rsl, 1d-2, z0_rsl > 1d-2)
1297 StabilityMethod, nz_above, z_array, zh, L_MOD, sfr_surf, FAI, PAI, & !input
1298 psihatm_array, psihath_array, zH_RSL, L_MOD_RSL, Lc, beta, zd_RSL, z0_RSL, elm, Scc, fx)
1301 INTEGER,
INTENT(in) :: StabilityMethod
1302 INTEGER,
INTENT(IN) :: nz_above
1303 REAL(KIND(1D0)),
DIMENSION(nz_above),
INTENT(in) :: z_array
1304 REAL(KIND(1D0)),
DIMENSION(nz_above),
INTENT(out) :: psihatm_array
1305 REAL(KIND(1D0)),
DIMENSION(nz_above),
INTENT(out) :: psihath_array
1306 REAL(KIND(1D0)),
INTENT(in) :: zh
1307 REAL(KIND(1D0)),
INTENT(in) :: FAI
1309 REAL(KIND(1D0)),
INTENT(in) :: PAI
1311 REAL(KIND(1D0)),
INTENT(in) :: L_MOD
1312 REAL(KIND(1D0)),
DIMENSION(nsurf),
INTENT(in) :: sfr_surf
1317 REAL(KIND(1D0)),
INTENT(out) :: L_MOD_RSL
1318 REAL(KIND(1D0)),
INTENT(out) :: zH_RSL
1321 REAL(KIND(1D0)),
INTENT(out) :: Lc
1322 REAL(KIND(1D0)),
INTENT(out) :: beta
1323 REAL(KIND(1D0)),
INTENT(out) :: zd_RSL
1324 REAL(KIND(1D0)),
INTENT(out) :: z0_RSL
1325 REAL(KIND(1D0)),
INTENT(out) :: elm
1326 REAL(KIND(1D0)),
INTENT(out) :: Scc
1327 REAL(KIND(1D0)),
INTENT(out) :: fx
1331 REAL(KIND(1D0)) :: sfr_tr
1332 REAL(KIND(1D0)) :: psihatm_Zh
1335 REAL(KIND(1D0)) :: lc_over_L
1336 REAL(KIND(1D0)) :: z_top, z_mid, z_btm
1337 REAL(KIND(1D0)) :: psihatm_top, psihatm_mid, psihatm_btm
1338 REAL(KIND(1D0)) :: psihath_top, psihath_mid, psihath_btm
1339 REAL(KIND(1D0)) :: cm, ch, c2m, c2h
1342 REAL(KIND(1D0)) :: Lc_min
1343 REAL(KIND(1D0)) :: dim_bluffbody
1345 REAL(KIND(1D0)),
PARAMETER :: cd_tree = 1.2
1346 REAL(KIND(1D0)),
PARAMETER :: a_tree = 0.05
1348 REAL(KIND(1D0)),
PARAMETER :: beta_N = 0.40
1349 REAL(KIND(1D0)),
PARAMETER :: pi = 4.*atan(1.0)
1351 REAL(KIND(1D0)),
PARAMETER :: planF_low = 1e-6
1352 REAL(KIND(1D0)),
PARAMETER :: kappa = 0.40
1354 REAL(KIND(1D0)),
PARAMETER :: r = 0.1
1355 REAL(KIND(1D0)),
PARAMETER :: a1 = 4.,
a2 = -0.1,
a3 = 1.5, a4 = -1.
1356 REAL(KIND(1D0)),
PARAMETER :: Zh_min = 0.4
1357 REAL(KIND(1D0)),
PARAMETER :: porosity_evetr = 0.32
1363 zh_rsl = max(zh, zh_min)
1382 lc = (1.-pai)/fai*zh_rsl
1388 dim_bluffbody = zh_rsl*pai/fai
1389 lc_min = dim_bluffbody*pai**(-0.5)/3.
1390 lc = merge(lc, lc_min, lc > lc_min)
1393 lc_over_l = lc/l_mod
1395 IF (lc_over_l > 0)
THEN
1396 lc_over_l = min(2., lc_over_l)
1398 lc_over_l = max(-2., lc_over_l)
1401 l_mod_rsl = lc/lc_over_l
1405 beta =
cal_beta_rsl(stabilitymethod, pai, sfr_tr, lc_over_l)
1408 scc = 0.5 + 0.3*tanh(2.*lc_over_l)
1409 fx = 0.5*((1.+4.*r*scc)**0.5) - 0.5
1417 CALL cal_ch(stabilitymethod, zh_rsl, zd_rsl, lc, beta, l_mod_rsl, scc, fx, c2h, ch)
1418 CALL cal_cm(stabilitymethod, zh_rsl, zd_rsl, lc, beta, l_mod_rsl, c2m, cm)
1423 psihatm_array(nz_above) = 0
1424 psihatm_array(nz_above - 1) = 0
1428 psihath_array(nz_above) = 0
1429 psihath_array(nz_above - 1) = 0
1431 DO iz = nz_above, 3, -1
1433 z_mid = z_array(iz - 1)
1434 z_btm = z_array(iz - 2)
1438 psihatm_top, psihatm_mid, &
1439 z_top, z_mid, z_btm, &
1441 zh_rsl, zd_rsl, l_mod, beta, elm, lc)
1442 psihatm_array(iz - 2) = psihatm_btm
1443 psihatm_top = psihatm_mid
1444 psihatm_mid = psihatm_btm
1448 psihath_top, psihath_mid, &
1449 z_top, z_mid, z_btm, &
1451 zh_rsl, zd_rsl, l_mod, beta, elm, lc)
1452 psihath_top = psihath_mid
1453 psihath_mid = psihath_btm
1454 psihath_array(iz - 2) = psihath_btm
1459 psihatm_zh = psihatm_array(1)
1460 z0_rsl =
cal_z0_rsl(stabilitymethod, zh_rsl, zd_rsl, beta, l_mod_rsl, psihatm_zh)
1464 FUNCTION cal_beta_rsl(StabilityMethod, PAI, sfr_tr, lc_over_L)
RESULT(beta)
1469 INTEGER,
INTENT(in) :: stabilitymethod
1470 REAL(kind(1d0)),
INTENT(in) :: pai
1471 REAL(kind(1d0)),
INTENT(in) :: sfr_tr
1472 REAL(kind(1d0)),
INTENT(in) :: lc_over_l
1475 REAL(kind(1d0)) :: beta
1478 REAL(kind(1d0)) :: betahf
1479 REAL(kind(1d0)) :: betanl
1481 REAL(kind(1d0)),
PARAMETER :: kappa = 0.4
1482 REAL(kind(1d0)),
PARAMETER ::
a1 = 4.,
a2 = -0.1,
a3 = 1.5, a4 = -1.
1486 REAL(kind(1d0)) :: betan2
1494 betan2 = 0.30*sfr_tr/pai + (pai - sfr_tr)/pai*0.4
1499 betahf =
cal_beta_lc(stabilitymethod, betan2, lc_over_l)
1501 betanl =
cal_beta_lc(stabilitymethod, kappa/2., lc_over_l)
1503 IF (lc_over_l >
a2)
THEN
1506 beta = betanl + ((betahf - betanl)/(1.+
a1*abs(lc_over_l -
a2)**
a3))
1509 IF (beta > 0.5)
THEN
1515 FUNCTION cal_beta_lc(stabilityMethod, beta0, lc_over_l)
RESULT(beta_x)
1520 INTEGER,
INTENT(in) :: stabilitymethod
1521 REAL(kind(1d0)),
INTENT(in) :: beta0
1522 REAL(kind(1d0)),
INTENT(in) :: lc_over_l
1523 REAL(kind(1d0)) :: beta_x
1525 REAL(kind(1d0)) :: phim, err, beta_x0
1535 DO WHILE ((err > 0.001) .AND. (it < 20))
1536 beta_x0 = beta0/phim
1537 phim =
stab_phi_heat(stabilitymethod, (beta_x0**2)*lc_over_l)
1541 err = abs(beta_x - beta_x0)
integer, parameter bldgsurf
integer, parameter conifsurf
integer, parameter ncolumnsdataoutrsl
integer, parameter decidsurf
real(kind(1d0)) function stab_psi_mom(stabilitymethod, zl)
subroutine cal_stab(stabilitymethod, zzd, z0m, zdm, avu1, temp_c, qh_init, avdens, avcp, l_mod, tstar, ustar, zl)
real(kind(1d0)) function stab_phi_mom(stabilitymethod, zl)
real(kind(1d0)) function stab_phi_heat(stabilitymethod, zl)
real(kind(1d0)), parameter neut_limit
real(kind(1d0)) function stab_psi_heat(stabilitymethod, zl)
real(kind(1d0)) function rh2qa(rh_dec, pres_hpa, ta_degc)
real(kind(1d0)) function qa2rh(qa_gkg, pres_hpa, ta_degc)
subroutine suews_cal_roughnessparameters(roughlenmommethod, faimethod, sfr_surf, surfacearea, bldgh, evetreeh, dectreeh, porosity_dectr, faibldg, faievetree, faidectree, z0m_in, zdm_in, z, fai, pai, zh, z0m, zdm, zzd)
real(kind(1d0)) function cal_beta_rsl(stabilitymethod, pai, sfr_tr, lc_over_l)
subroutine rslprofile_dts(diagmethod, zh, z0m, zdm, z0v, l_mod, sfr_paved, sfr_bldg, sfr_evetr, sfr_dectr, sfr_grass, sfr_bsoil, sfr_water, fai, pai, stabilitymethod, ra_h, avcp, lv_j_kg, avdens, avu1, temp_c, avrh, press_hpa, zmeas, qh, qe, t2_c, q2_gkg, u10_ms, rh2, dataoutlinersl)
real(kind(1d0)) function interp_z(z_x, z, v)
subroutine cal_cm(stabilitymethod, zh_rsl, zd_rsl, lc, beta, l_mod, c2, cm)
real(kind(1d0)) function cal_phim_hat(stabilitymethod, z, zh_rsl, l_mod, beta, lc)
subroutine cal_ch(stabilitymethod, zh_rsl, zd_rsl, lc, beta, l_mod, scc, f, c2h, ch)
subroutine rslprofile(diagmethod, zh, z0m, zdm, z0v, l_mod, sfr_surf, fai, pai, stabilitymethod, ra_h, avcp, lv_j_kg, avdens, avu1, temp_c, avrh, press_hpa, zmeas, qh, qe, t2_c, q2_gkg, u10_ms, rh2, dataoutlinersl)
real(kind(1d0)) function cal_zd_rsl(zh_rsl, beta, lc)
real(kind(1d0)) function cal_beta_lc(stabilitymethod, beta0, lc_over_l)
real(kind(1d0)) function cal_z0_rsl(stabilitymethod, zh_rsl, zd_rsl, beta, l_mod_rsl, psihatm_zh)
recursive real(kind(1d0)) function cal_psim_hat(stabilitymethod, psihatm_top, psihatm_mid, z_top, z_mid, z_btm, cm, c2, zh_rsl, zd_rsl, l_mod, beta, elm, lc)
real(kind(1d0)) function cal_elm_rsl(beta, lc)
subroutine rsl_cal_prms(stabilitymethod, nz_above, z_array, zh, l_mod, sfr_surf, fai, pai, psihatm_array, psihath_array, zh_rsl, l_mod_rsl, lc, beta, zd_rsl, z0_rsl, elm, scc, fx)
recursive real(kind(1d0)) function cal_psih_hat(stabilitymethod, psihath_top, psihath_mid, z_top, z_mid, z_btm, ch, c2h, zh_rsl, zd_rsl, l_mod, beta, elm, lc)