37 SUBROUTINE beers_cal_main(iy, id, dectime, lamdaP, lamdaF, avkdn, ldown, Temp_C, avrh, &
38 Press_hPa, Tsurf, lat, lng, alt, timezone, zenith_deg, azimuth, &
39 alb_ground, alb_bldg, emis_ground, emis_wall, &
45 INTEGER,
INTENT(in) :: iy
46 INTEGER,
INTENT(in) :: id
47 REAL(KIND(1D0)),
INTENT(in) :: lamdaP
48 REAL(KIND(1D0)),
INTENT(in) :: lamdaF
53 REAL(KIND(1D0)),
INTENT(in) :: Press_hPa
54 REAL(KIND(1D0)),
INTENT(in) :: Temp_C
55 REAL(KIND(1D0)),
INTENT(in) :: avrh
56 REAL(KIND(1D0)),
INTENT(in) :: avkdn
57 REAL(KIND(1D0)),
INTENT(in) :: ldown
58 REAL(KIND(1D0)),
INTENT(in) :: Tsurf
61 REAL(KIND(1D0)),
INTENT(in) :: zenith_deg
62 REAL(KIND(1D0)),
INTENT(in) :: azimuth
63 REAL(KIND(1D0)),
INTENT(in) :: dectime
64 REAL(KIND(1D0)),
INTENT(in) :: timezone, lat, lng, alt
65 REAL(KIND(1D0)),
INTENT(in) :: alb_bldg
66 REAL(KIND(1D0)),
INTENT(in) :: alb_ground
67 REAL(KIND(1D0)),
INTENT(in) :: emis_wall
68 REAL(KIND(1D0)),
INTENT(in) :: emis_ground
70 REAL(KIND(1D0)),
PARAMETER :: absL = 0.97
71 REAL(KIND(1D0)),
PARAMETER :: absK = 0.7
72 REAL(KIND(1D0)),
PARAMETER :: Fside = 0.22
73 REAL(KIND(1D0)),
PARAMETER :: Fup = 0.06
76 REAL(KIND(1D0)),
DIMENSION(ncolumnsDataOutBEERS - 5),
INTENT(OUT) :: dataOutLineBEERS
78 REAL(KIND(1D0)) :: t, psi
79 REAL(KIND(1D0)) :: altitude, zen
80 REAL(KIND(1D0)) :: CI, c, I0, Kt, Tw, Tg
81 REAL(KIND(1D0)) :: Ta, RH, P, radG, radD, radI
82 REAL(KIND(1D0)) :: I0et, CIuncorr
83 REAL(KIND(1D0)) :: SNDN, SNUP, DEC, DAYL
84 REAL(KIND(1D0)) :: msteg, emis_sky, ea
85 REAL(KIND(1D0)) :: shadowground, shadowwalls, shadowroof
88 REAL(KIND(1D0)) :: CIlatenight
89 REAL(KIND(1D0)) :: dectime_sunrise, zen_sunrise, I0_sunrise
94 REAL(KIND(1D0)) :: svfalfa
96 REAL(KIND(1D0)) :: Tmrt, Sstr, F_sh
98 REAL(KIND(1D0)) :: tmp, altmax
99 REAL(KIND(1D0)) :: svf_bldg_veg
100 REAL(KIND(1D0)) :: svf_ground, svf_roof
101 REAL(KIND(1D0)) :: svf_veg
102 REAL(KIND(1D0)) :: svf_aveg
103 REAL(KIND(1D0)) :: Kdown, Keast, Knorth, Ksouth, Kup2d, Kwest
104 REAL(KIND(1D0)) :: Ldown2d, Least, Lnorth, Lsouth, Lup2d, Lwest
111 REAL(KIND(1D0)),
PARAMETER :: SBC = 5.67051e-8
113 INTEGER,
PARAMETER :: onlyglobal = 1
120 INTEGER,
PARAMETER :: SOLWEIG_ldown = 1
121 INTEGER,
PARAMETER :: BEERS_tsurf = 1
165 tmp = 1 - (svf_ground + svf_veg - 1)
166 IF (tmp <= 1.e-6) tmp = 1.e-6
167 svfalfa = asin(exp(log(tmp)/2))
170 svf_bldg_veg = (svf_ground - (1 - svf_veg)*(1 - psi))
173 CALL daylen(doy, lat, dayl, dec, sndn, snup)
175 altitude = 90 - zenith_deg
178 ea = 6.107*10**((7.5*ta)/(237.3 + ta))*(rh/100)
179 msteg = 46.5*(ea/(ta + 273.15))
180 emis_sky = (1 - (1 + msteg)*exp(-((1.2 + 3.0*msteg)**0.5)))
183 IF (altitude > 0.1)
THEN
188 i0, ci, kt, i0et, ciuncorr)
192 IF (onlyglobal == 1)
THEN
202 CALL kroof(radi, radd, radg, f_sh, altitude, svf_roof, svf_veg, shadowroof, psi, alb_bldg, kdown)
208 kup2d = alb_ground*( &
209 radi*shadowground*sin(altitude*
deg2rad) &
210 + radd*svf_bldg_veg &
211 + alb_bldg*(1 - svf_bldg_veg)*(radg*(1 - f_sh) + radd*f_sh))
215 svf_ground, svf_veg, shadowground, f_sh, &
216 radi, radg, radd, azimuth, altitude, psi, t, alb_ground, alb_bldg, &
217 keast, knorth, ksouth, kwest)
219 IF (beers_tsurf == 1)
THEN
220 CALL tsurfbeers(iy, ta, rh, radi, i0, dectime, snup, altitude, zen, timezone, lat, lng, alt, &
229 lup2d = sbc*emis_ground*((shadowground*tg + ta + 273.15)**4)
234 cilatenight, dectime_sunrise, zen_sunrise, i0_sunrise)
242 IF (beers_tsurf == 1)
THEN
254 IF (dectime < (doy + 0.5) .AND. dectime > doy .AND. altitude < 1.0)
THEN
288 lup2d = sbc*emis_ground*((ta + tg + 273.15)**4)
293 IF (solweig_ldown == 1)
THEN
294 ldown2d = (svf_roof + svf_veg - 1)*emis_sky*sbc*((ta + 273.15)**4) &
295 + (2 - svf_veg - svf_aveg)*emis_wall*sbc*((ta + 273.15)**4) &
296 + (svf_aveg - svf_roof)*emis_wall*sbc*((ta + 273.15 + tw)**4) &
297 + (2 - svf_roof - svf_veg)*(1 - emis_wall)*emis_sky*sbc*((ta + 273.15)**4)
301 ldown2d = ldown2d*(1 - c) + c*( &
302 (svf_roof + svf_veg - 1)*sbc*((ta + 273.15)**4) &
303 + (2 - svf_veg - svf_aveg)*emis_wall*sbc*((ta + 273.15)**4) &
304 + (svf_aveg - svf_roof)*emis_wall*sbc*((ta + 273.15 + tw)**4) &
305 + (2 - svf_roof - svf_veg)*(1 - emis_wall)*sbc*((ta + 273.15)**4))
309 ldown2d = (svf_roof + svf_veg - 1)*ldown &
310 + (2 - svf_veg - svf_aveg)*emis_wall*sbc*((ta + 273.15)**4) &
311 + (svf_aveg - svf_roof)*emis_wall*sbc*((ta + 273.15 + tw)**4) &
312 + (2 - svf_roof - svf_veg)*(1 - emis_wall)*ldown
316 CALL lwalls(svf_ground, svf_veg, svf_aveg, &
318 altitude, ta, tw, sbc, emis_wall, &
319 emis_sky, t, ci, azimuth, ldown, svfalfa, f_sh, &
320 least, lnorth, lsouth, lwest)
323 sstr = absk*(kdown*fup + kup2d*fup + knorth*fside + keast*fside + ksouth*fside + kwest*fside) &
324 + absl*(ldown2d*fup + lup2d*fup + lnorth*fside + least*fside + lsouth*fside + lwest*fside)
325 tmrt = sqrt(sqrt((sstr/(absl*sbc)))) - 273.15
327 dataoutlinebeers = [azimuth, altitude, radg, radi, radd, &
328 kdown, kup2d, ksouth, kwest, knorth, keast, &
329 ldown2d, lup2d, lsouth, lwest, lnorth, least, &
330 tmrt, i0, ci, shadowground, shadowwalls, svf_ground, svf_roof, svf_bldg_veg, &
486 SUBROUTINE tsurfbeers(iy, Ta, RH, radI, I0, dectime, SNUP, altitude, zen, timezone, lat, lng, alt, &
492 INTEGER,
INTENT(in) :: iy
493 REAL(KIND(1D0)),
INTENT(in) :: Ta, RH, radI, I0
494 REAL(KIND(1D0)),
INTENT(in) :: dectime, SNUP, altitude, zen
495 REAL(KIND(1D0)),
INTENT(in) :: timezone, lat, lng, alt
496 REAL(KIND(1D0)),
INTENT(out) :: Tg, Tgwall, altmax
497 REAL(KIND(1D0)),
PARAMETER :: TgK = 0.37
498 REAL(KIND(1D0)),
PARAMETER :: Tstart = 3.41
499 REAL(KIND(1D0)),
PARAMETER :: TgK_wall = 0.37
500 REAL(KIND(1D0)),
PARAMETER :: Tstart_wall = -3.41
501 REAL(KIND(1D0)),
PARAMETER :: TmaxLST = 15.
502 REAL(KIND(1D0)),
PARAMETER :: TmaxLST_wall = 15.
504 REAL(KIND(1D0)) :: Ktc, notU, Tgamp, Tgampwall, radI0, corr, CI_Tg
505 REAL(KIND(1D0)) :: fifteen, sunmaximum, zen_sunmax, dectimemax, azimuth
511 DO WHILE (sunmaximum <= 90.-zen_sunmax)
512 sunmaximum = 90.-zen_sunmax
513 fifteen = fifteen + 15./1440.
514 dectimemax = floor(dectime) + 10/24.+fifteen
516 REAL(iy, KIND(1D0)), &
518 timezone, lat, lng, alt, &
524 tgamp = (tgk*altmax - tstart) + tstart
525 tgampwall = (tgk_wall*altmax - (tstart_wall)) + (tstart_wall)
526 tg = tgamp*sin((((dectime - floor(dectime)) - snup/24)/(tmaxlst/24 - snup/24))*
pi/2) + tstart
527 tgwall = tgampwall*sin((((dectime - floor(dectime)) - snup/24)/(tmaxlst_wall/24 - snup/24))*
pi/2) + (tstart_wall)
531 corr = 0.1473*log(90.-(zen*
rad2deg)) + 0.3454
532 ci_tg = (radi/radi0) + (1.-corr)
535 IF (ci_tg > 1) ci_tg = 1
537 tgwall = tgwall*ci_tg
539 IF (tgwall < 0) tgwall = 0
586 I0, CI, Kt, I0et, CIuncorr)
590 INTEGER,
INTENT(in) :: DOY
591 REAL(KIND(1D0)),
INTENT(in) :: zen
592 REAL(KIND(1D0)),
INTENT(in) :: Ta_degC
593 REAL(KIND(1D0)),
INTENT(in) :: RH_frac
594 REAL(KIND(1D0)),
INTENT(in) :: P_kPa
595 REAL(KIND(1D0)),
INTENT(in) :: radG
596 REAL(KIND(1D0)),
INTENT(in) :: lat
597 REAL(KIND(1D0)),
INTENT(out) :: I0
598 REAL(KIND(1D0)),
INTENT(out) :: CI
599 REAL(KIND(1D0)),
INTENT(out) :: Kt
600 REAL(KIND(1D0)),
INTENT(out) :: I0et
601 REAL(KIND(1D0)),
INTENT(out) :: CIuncorr
602 REAL(KIND(1D0)) :: iG, Itoa, p
603 REAL(KIND(1D0)),
DIMENSION(4) :: G
619 IF (p_kpa == -999)
THEN
627 m = 35.*cos(zen)*((1224.*(cos(zen)**2) + 1.)**(-1./2.))
628 trpg = 1.021 - 0.084*(m*(0.000949*p + 0.051))**0.5
633 g = [3.37, 2.85, 2.80, 2.64]
634 ELSE IF (lat >= 10 .AND. lat < 20)
THEN
635 g = [2.99, 3.02, 2.70, 2.93]
636 ELSE IF (lat >= 20 .AND. lat < 30)
THEN
637 g = [3.60, 3.00, 2.98, 2.93]
638 ELSE IF (lat >= 30 .AND. lat < 40)
THEN
639 g = [3.04, 3.11, 2.92, 2.94]
640 ELSE IF (lat >= 40 .AND. lat < 50)
THEN
641 g = [2.70, 2.95, 2.77, 2.71]
642 ELSE IF (lat >= 50 .AND. lat < 60)
THEN
643 g = [2.52, 3.07, 2.67, 2.93]
644 ELSE IF (lat >= 60 .AND. lat < 70)
THEN
645 g = [1.76, 2.69, 2.61, 2.61]
646 ELSE IF (lat >= 70 .AND. lat < 80)
THEN
647 g = [1.60, 1.67, 2.24, 2.63]
648 ELSE IF (lat >= 80 .AND. lat < 90)
THEN
649 g = [1.11, 1.44, 1.94, 2.02]
651 IF (doy > 335 .OR. doy <= 60)
THEN
653 ELSE IF (doy > 60 .AND. doy <= 152)
THEN
655 ELSE IF (doy > 152 .AND. doy <= 244)
THEN
657 ELSE IF (doy > 244 .AND. doy <= 335)
THEN
663 td = (b2*(((
a2*ta_degc)/(b2 + ta_degc)) + log(rh_frac)))/(
a2 - (((
a2*ta_degc)/(b2 + ta_degc)) + log(rh_frac)))
665 u = exp(0.1133 - log(ig + 1) + 0.0393*td)
666 tw = 1 - 0.077*((u*m)**0.3)
669 i0 = itoa*cos(zen)*trpg*tw*d*tar
679 corr = 0.1473*log(90 - (zen/
pi*180)) + 0.3454
682 ci = ciuncorr + (1 - corr)
683 i0et = itoa*cos(zen)*d
823 svf, svfveg, shadow, F_sh, &
824 radI, radG, radD, azimuth, altitude, psi, t, alb_ground, alb_bldg, & ! input
825 Keast, Knorth, Ksouth, Kwest)
830 REAL(KIND(1D0)),
INTENT(in) :: radI
831 REAL(KIND(1D0)),
INTENT(in) :: radG
832 REAL(KIND(1D0)),
INTENT(in) :: radD
833 REAL(KIND(1D0)),
INTENT(in) :: azimuth
834 REAL(KIND(1D0)),
INTENT(in) :: altitude
835 REAL(KIND(1D0)),
INTENT(in) :: psi
836 REAL(KIND(1D0)),
INTENT(in) :: t
837 REAL(KIND(1D0)),
INTENT(in) :: alb_bldg
838 REAL(KIND(1D0)),
INTENT(in) :: alb_ground
841 REAL(KIND(1D0)),
INTENT(in) :: shadow, F_sh, svf, svfveg
842 REAL(KIND(1D0)),
INTENT(out) :: Keast, Knorth, Ksouth, Kwest
844 REAL(KIND(1D0)) :: vikttot, aziE, aziN, aziS, aziW
845 REAL(KIND(1D0)) :: viktveg, viktwall
846 REAL(KIND(1D0)) :: KeastI, KsouthI, KwestI, KnorthI, Kuptowall
847 REAL(KIND(1D0)) :: KeastDG, KsouthDG, KwestDG, KnorthDG
848 REAL(KIND(1D0)) :: svfE, svfS, svfW, svfN, svfEveg, svfSveg, svfWveg, svfNveg, svfbuveg
849 REAL(KIND(1D0)) :: gvfalb, gvfalbnosh
852 REAL(KIND(1D0)) :: svfviktbuveg
857 azis = azimuth - 90 + t
858 aziw = azimuth - 180 + t
859 azin = azimuth - 270 + t
871 IF (azimuth > (360 - t) .OR. azimuth <= (180 - t))
THEN
876 IF (azimuth > (90 - t) .AND. azimuth <= (270 - t))
THEN
881 IF (azimuth > (180 - t) .AND. azimuth <= (360 - t))
THEN
886 IF (azimuth <= (90 - t) .OR. azimuth > (270 - t))
THEN
893 svfbuveg = svfs - (1.0 - svfsveg)*(1.0 - psi)
894 gvfalb = shadow*alb_ground
895 gvfalbnosh = (1 - shadow)*alb_ground
896 kuptowall = (gvfalb*radi*sin(altitude*
deg2rad)) &
897 + (radd*svfbuveg + alb_bldg*(1 - svfbuveg)*(radg*(1 - f_sh) + radd*f_sh))*gvfalbnosh
899 CALL kvikt_veg(svfe, svfeveg, vikttot, viktveg, viktwall)
900 svfviktbuveg = (viktwall + (viktveg)*(1 - psi))
901 keastdg = (radd*(1 - svfviktbuveg) + alb_bldg*(svfviktbuveg*(radg*(1 - f_sh) + radd*f_sh)) + kuptowall)*0.5
902 keast = keasti + keastdg
904 CALL kvikt_veg(svfs, svfsveg, vikttot, viktveg, viktwall)
905 svfviktbuveg = (viktwall + (viktveg)*(1 - psi))
906 ksouthdg = (radd*(1 - svfviktbuveg) + alb_bldg*(svfviktbuveg*(radg*(1 - f_sh) + radd*f_sh)) + kuptowall)*0.5
907 ksouth = ksouthi + ksouthdg
909 CALL kvikt_veg(svfw, svfwveg, vikttot, viktveg, viktwall)
910 svfviktbuveg = (viktwall + (viktveg)*(1 - psi))
911 kwestdg = (radd*(1 - svfviktbuveg) + alb_bldg*(svfviktbuveg*(radg*(1 - f_sh) + radd*f_sh)) + kuptowall)*0.5
912 kwest = kwesti + kwestdg
914 CALL kvikt_veg(svfn, svfnveg, vikttot, viktveg, viktwall)
915 svfviktbuveg = (viktwall + (viktveg)*(1 - psi))
916 knorthdg = (radd*(1 - svfviktbuveg) + alb_bldg*(svfviktbuveg*(radg*(1 - f_sh) + radd*f_sh)) + kuptowall)*0.5
917 knorth = knorthi + knorthdg
961 SUBROUTINE lwalls(svf, svfveg, svfaveg, &
963 altitude, Ta, Tw, SBC, emis_wall, emis_sky, t, CI, azimuth, ldown, svfalfa, F_sh_in, & !input
964 Least, Lnorth, Lsouth, Lwest)
967 REAL(KIND(1D0)),
INTENT(in) :: altitude, Ta, Tw, SBC, emis_wall, emis_sky, t, CI, azimuth, ldown
969 REAL(KIND(1D0)),
INTENT(in) :: svfalfa, svf, svfveg, svfaveg, F_sh_in
970 REAL(KIND(1D0)),
INTENT(in) :: Ldown2d, Lup2d
971 REAL(KIND(1D0)),
INTENT(out) :: Least, Lnorth, Lsouth, Lwest
973 REAL(KIND(1D0)) :: vikttot, aziE, aziN, aziS, aziW, c
975 REAL(KIND(1D0)) :: svfalfaE, svfalfaS, svfalfaW, svfalfaN
976 REAL(KIND(1D0)) :: alfaB, betaB, betasun
977 REAL(KIND(1D0)) :: Lground, Lrefl, Lsky, Lsky_allsky, Lveg, Lwallsh, Lwallsun
984 INTEGER,
PARAMETER :: SOLWEIG_ldown = 1
986 REAL(KIND(1D0)) :: viktveg, viktsky, viktrefl, viktwall
987 REAL(KIND(1D0)) :: svfE, svfS, svfW, svfN, svfEveg, svfSveg, svfWveg, svfNveg
988 REAL(KIND(1D0)) :: svfEaveg, svfSaveg, svfWaveg, svfNaveg, F_sh
1041 azin = azimuth - 90 + t
1042 azie = azimuth - 180 + t
1043 azis = azimuth - 270 + t
1045 f_sh = 2.*f_sh_in - 1.
1047 IF (solweig_ldown == 1)
THEN
1049 lsky_allsky = emis_sky*sbc*((ta + 273.15)**4)*(1 - c) + c*sbc*((ta + 273.15)**4)
1055 CALL lvikt_veg(svfe, svfeveg, svfeaveg, vikttot, &
1056 viktveg, viktsky, viktrefl, viktwall)
1058 IF (altitude > 0)
THEN
1059 alfab = atan(svfalfae)
1060 betab = atan(tan((svfalfae)*f_sh))
1061 betasun = ((alfab - betab)/2) + betab
1062 IF (azimuth > (180 - t) .AND. azimuth <= (360 - t))
THEN
1063 lwallsun = sbc*emis_wall*((ta + 273.15 + tw*sin(azie*(
pi/180)))**4)* &
1064 viktwall*(1 - f_sh)*cos(betasun)*0.5
1065 lwallsh = sbc*emis_wall*((ta + 273.15)**4)*viktwall*f_sh*0.5
1068 lwallsh = sbc*emis_wall*((ta + 273.15)**4)*viktwall*0.5
1072 lwallsh = sbc*emis_wall*((ta + 273.15)**4)*viktwall*0.5
1074 lsky = ((svfe + svfeveg - 1)*lsky_allsky)*viktsky*0.5
1075 lveg = sbc*emis_wall*((ta + 273.15)**4)*viktveg*0.5
1077 lrefl = (ldown2d + lup2d)*(viktrefl)*(1 - emis_wall)*0.5
1078 least = lsky + lwallsun + lwallsh + lveg + lground + lrefl
1081 CALL lvikt_veg(svfs, svfsveg, svfsaveg, vikttot, &
1082 viktveg, viktsky, viktrefl, viktwall)
1084 IF (altitude > 0)
THEN
1085 alfab = atan(svfalfas)
1086 betab = atan(tan((svfalfas)*f_sh))
1087 betasun = ((alfab - betab)/2) + betab
1088 IF (azimuth <= (90 - t) .OR. azimuth > (270 - t))
THEN
1089 lwallsun = sbc*emis_wall*((ta + 273.15 + tw*sin(azis*(
pi/180)))**4)* &
1090 viktwall*(1 - f_sh)*cos(betasun)*0.5
1091 lwallsh = sbc*emis_wall*((ta + 273.15)**4)*viktwall*f_sh*0.5
1094 lwallsh = sbc*emis_wall*((ta + 273.15)**4)*viktwall*0.5
1098 lwallsh = sbc*emis_wall*((ta + 273.15)**4)*viktwall*0.5
1100 lsky = ((svfs + svfsveg - 1)*lsky_allsky)*viktsky*0.5
1101 lveg = sbc*emis_wall*((ta + 273.15)**4)*viktveg*0.5
1103 lrefl = (ldown2d + lup2d)*(viktrefl)*(1 - emis_wall)*0.5
1104 lsouth = lsky + lwallsun + lwallsh + lveg + lground + lrefl
1107 CALL lvikt_veg(svfw, svfwveg, svfwaveg, vikttot, &
1108 viktveg, viktsky, viktrefl, viktwall)
1110 IF (altitude > 0)
THEN
1111 alfab = atan(svfalfaw)
1112 betab = atan(tan((svfalfaw)*f_sh))
1113 betasun = ((alfab - betab)/2) + betab
1114 IF (azimuth > (360 - t) .OR. azimuth <= (180 - t))
THEN
1115 lwallsun = sbc*emis_wall*((ta + 273.15 + tw*sin(aziw*(
pi/180)))**4)* &
1116 viktwall*(1 - f_sh)*cos(betasun)*0.5
1117 lwallsh = sbc*emis_wall*((ta + 273.15)**4)*viktwall*f_sh*0.5
1120 lwallsh = sbc*emis_wall*((ta + 273.15)**4)*viktwall*0.5
1124 lwallsh = sbc*emis_wall*((ta + 273.15)**4)*viktwall*0.5
1126 lsky = ((svfw + svfwveg - 1)*lsky_allsky)*viktsky*0.5
1127 lveg = sbc*emis_wall*((ta + 273.15)**4)*viktveg*0.5
1129 lrefl = (ldown2d + lup2d)*(viktrefl)*(1 - emis_wall)*0.5
1130 lwest = lsky + lwallsun + lwallsh + lveg + lground + lrefl
1133 CALL lvikt_veg(svfn, svfnveg, svfnaveg, vikttot, &
1134 viktveg, viktsky, viktrefl, viktwall)
1136 IF (altitude > 0)
THEN
1137 alfab = atan(svfalfan)
1138 betab = atan(tan((svfalfan)*f_sh))
1139 betasun = ((alfab - betab)/2) + betab
1140 IF (azimuth > (90 - t) .AND. azimuth <= (270 - t))
THEN
1141 lwallsun = sbc*emis_wall*((ta + 273.15 + tw*sin(azin*(
pi/180)))**4)* &
1142 viktwall*(1 - f_sh)*cos(betasun)*0.5
1143 lwallsh = sbc*emis_wall*((ta + 273.15)**4)*viktwall*f_sh*0.5
1146 lwallsh = sbc*emis_wall*((ta + 273.15)**4)*viktwall*0.5
1150 lwallsh = sbc*emis_wall*((ta + 273.15)**4)*viktwall*0.5
1152 lsky = ((svfn + svfnveg - 1)*lsky_allsky)*viktsky*0.5
1153 lveg = sbc*emis_wall*((ta + 273.15)**4)*viktveg*0.5
1155 lrefl = (ldown2d + lup2d)*(viktrefl)*(1 - emis_wall)*0.5
1156 lnorth = lsky + lwallsun + lwallsh + lveg + lground + lrefl