SUEWS API Site
Documentation of SUEWS source code
suews_ctrl_driver.f95
Go to the documentation of this file.
1 !========================================================================================
2 ! SUEWS driver subroutines
3 ! TS 31 Aug 2017: initial version
4 ! TS 02 Oct 2017: added as the generic wrapper
5 ! TS 03 Oct 2017: added
7  ! only the following immutable objects are imported:
8  ! 1. functions/subroutines
9  ! 2. constant variables
10  USE meteo, ONLY: qsatf, rh2qa, qa2rh
13  USE anohm_module, ONLY: anohm
16  USE estm_module, ONLY: estm
23  USE ctrl_output, ONLY: varlistall
25  use lumps_module, only: lumps_cal_qhqe
26  use evap_module, only: cal_evap
27  use rsl_module, only: rslprofile
29  use co2_module, only: co2_biogen
30  use evap_module, only: cal_evap
31  USE allocatearray, ONLY: &
32  nsurf, nvegsurf, &
38  use moist, only: avcp, avdens, lv_j_kg
40 
41  IMPLICIT NONE
42 
43 CONTAINS
44  ! ===================MAIN CALCULATION WRAPPER FOR ENERGY AND WATER FLUX===========
45  SUBROUTINE suews_cal_main( &
46  AerodynamicResistanceMethod, AH_MIN, AHProf_24hr, AH_SLOPE_Cooling, & ! input&inout in alphabetical order
47  AH_SLOPE_Heating, &
48  alb, AlbMax_DecTr, AlbMax_EveTr, AlbMax_Grass, &
49  AlbMin_DecTr, AlbMin_EveTr, AlbMin_Grass, &
50  alpha_bioCO2, alpha_enh_bioCO2, alt, avkdn, avRh, avU1, BaseT, BaseTe, &
51  BaseTMethod, &
52  BaseT_HC, beta_bioCO2, beta_enh_bioCO2, bldgH, CapMax_dec, CapMin_dec, &
53  chAnOHM, CO2PointSource, cpAnOHM, CRWmax, CRWmin, DayWat, DayWatPer, &
54  DecTreeH, Diagnose, DiagQN, DiagQS, DRAINRT, &
55  dt_since_start, dqndt, qn1_av, dqnsdt, qn1_s_av, &
56  EF_umolCO2perJ, emis, EmissionsMethod, EnEF_v_Jkm, endDLS, EveTreeH, FAIBldg, &
57  FAIDecTree, FAIEveTree, Faut, FcEF_v_kgkm, fcld_obs, FlowChange, &
58  FrFossilFuel_Heat, FrFossilFuel_NonHeat, G1, G2, G3, G4, G5, G6, GDD_id, &
59  GDDFull, Gridiv, gsModel, H_maintain, HDD_id, HumActivity_24hr, &
60  IceFrac, id, Ie_a, Ie_end, Ie_m, Ie_start, imin, &
61  InternalWaterUse_h, &
62  IrrFracPaved, IrrFracBldgs, &
63  IrrFracEveTr, IrrFracDecTr, IrrFracGrass, &
64  IrrFracBSoil, IrrFracWater, &
65  isec, it, EvapMethod, &
66  iy, kkAnOHM, Kmax, LAI_id, LAICalcYes, LAIMax, LAIMin, LAI_obs, &
67  LAIPower, LAIType, lat, lenDay_id, ldown_obs, lng, MaxConductance, MaxFCMetab, MaxQFMetab, &
68  SnowWater, MetForcingData_grid, MinFCMetab, MinQFMetab, min_res_bioCO2, &
69  NARP_EMIS_SNOW, NARP_TRANS_SITE, NetRadiationMethod, &
70  OHM_coef, OHMIncQF, OHM_threshSW, &
71  OHM_threshWD, PipeCapacity, PopDensDaytime, &
72  PopDensNighttime, PopProf_24hr, PorMax_dec, PorMin_dec, &
73  Precip, PrecipLimit, PrecipLimitAlb, Press_hPa, &
74  QF0_BEU, Qf_A, Qf_B, Qf_C, &
75  qn1_obs, qs_obs, qf_obs, &
76  RadMeltFact, RAINCOVER, RainMaxRes, resp_a, resp_b, &
77  RoughLenHeatMethod, RoughLenMomMethod, RunoffToWater, S1, S2, &
78  SatHydraulicConduct, SDDFull, SDD_id, sfr, SMDMethod, SnowAlb, SnowAlbMax, &
79  SnowAlbMin, SnowPackLimit, SnowDens, SnowDensMax, SnowDensMin, SnowfallCum, SnowFrac, &
80  SnowLimBldg, SnowLimPaved, snowFrac_obs, SnowPack, SnowProf_24hr, snowUse, SoilDepth, &
81  soilstore_id, SoilStoreCap, StabilityMethod, startDLS, state_id, StateLimit, &
82  StorageHeatMethod, StoreDrainPrm, SurfaceArea, Tair_av, tau_a, tau_f, tau_r, &
83  Tmax_id, Tmin_id, &
84  BaseT_Cooling, BaseT_Heating, Temp_C, TempMeltFact, TH, &
85  theta_bioCO2, timezone, TL, TrafficRate, TrafficUnits, &
86  TraffProf_24hr, Ts5mindata_ir, tstep, tstep_prev, veg_type, &
87  WaterDist, WaterUseMethod, WetThresh, wu_m3, &
88  WUDay_id, DecidCap_id, albDecTr_id, albEveTr_id, albGrass_id, porosity_id, &
89  WUProfA_24hr, WUProfM_24hr, xsmd, Z, z0m_in, zdm_in, &
90  datetimeLine, dataOutLineSUEWS, dataOutLineSnow, dataOutLineESTM, dataoutLineRSL, dataOutLineSOLWEIG, &!output
91  DailyStateLine)!output
92 
93  IMPLICIT NONE
94 
95  ! ########################################################################################
96  ! input variables
97  INTEGER, INTENT(IN)::AerodynamicResistanceMethod
98  INTEGER, INTENT(IN)::BaseTMethod
99  INTEGER, INTENT(IN)::Diagnose
100  INTEGER, INTENT(IN)::DiagQN
101  INTEGER, INTENT(IN)::DiagQS
102  INTEGER, INTENT(IN)::startDLS
103  INTEGER, INTENT(IN)::endDLS
104  INTEGER, INTENT(IN)::EmissionsMethod
105  INTEGER, INTENT(IN)::Gridiv
106  INTEGER, INTENT(IN)::gsModel
107  INTEGER, INTENT(IN)::id
108  INTEGER, INTENT(IN)::Ie_end
109  INTEGER, INTENT(IN)::Ie_start
110  INTEGER, INTENT(IN)::isec
111  INTEGER, INTENT(IN)::imin
112  INTEGER, INTENT(IN)::it
113  INTEGER, INTENT(IN)::EvapMethod
114  INTEGER, INTENT(IN)::iy
115  INTEGER, INTENT(IN)::LAICalcYes
116  INTEGER, INTENT(IN)::NetRadiationMethod
117  INTEGER, INTENT(IN)::OHMIncQF
118  INTEGER, INTENT(IN)::RoughLenHeatMethod
119  INTEGER, INTENT(IN)::RoughLenMomMethod
120  INTEGER, INTENT(IN)::SMDMethod
121  INTEGER, INTENT(IN)::snowUse
122  INTEGER, INTENT(IN)::StabilityMethod
123  INTEGER, INTENT(IN)::StorageHeatMethod
124  INTEGER, INTENT(IN)::tstep
125  INTEGER, INTENT(IN)::tstep_prev ! tstep size of the previous step
126  INTEGER, INTENT(in)::dt_since_start ! time since simulation starts [s]
127  INTEGER, INTENT(IN)::veg_type
128  INTEGER, INTENT(IN)::WaterUseMethod
129 
130  REAL(KIND(1D0)), INTENT(IN)::AlbMax_DecTr
131  REAL(KIND(1D0)), INTENT(IN)::AlbMax_EveTr
132  REAL(KIND(1D0)), INTENT(IN)::AlbMax_Grass
133  REAL(KIND(1D0)), INTENT(IN)::AlbMin_DecTr
134  REAL(KIND(1D0)), INTENT(IN)::AlbMin_EveTr
135  REAL(KIND(1D0)), INTENT(IN)::AlbMin_Grass
136  REAL(KIND(1D0)), INTENT(IN)::alt
137  REAL(KIND(1D0)), INTENT(IN)::avkdn
138  REAL(KIND(1D0)), INTENT(IN)::avRh
139  REAL(KIND(1D0)), INTENT(IN)::avU1
140  REAL(KIND(1D0)), INTENT(IN)::BaseT_HC
141  REAL(KIND(1D0)), INTENT(IN)::bldgH
142  REAL(KIND(1D0)), INTENT(IN)::CapMax_dec
143  REAL(KIND(1D0)), INTENT(IN)::CapMin_dec
144  REAL(KIND(1D0)), INTENT(IN)::CO2PointSource
145  REAL(KIND(1D0)), INTENT(IN)::CRWmax
146  REAL(KIND(1D0)), INTENT(IN)::CRWmin
147  REAL(KIND(1D0)), INTENT(IN)::DecTreeH
148  REAL(KIND(1D0)), INTENT(IN)::DRAINRT
149  REAL(KIND(1D0)), INTENT(IN)::EF_umolCO2perJ
150  REAL(KIND(1D0)), INTENT(IN)::EnEF_v_Jkm
151  REAL(KIND(1D0)), INTENT(IN)::EveTreeH
152  REAL(KIND(1D0)), INTENT(IN)::FAIBldg
153  REAL(KIND(1D0)), INTENT(IN)::FAIDecTree
154  REAL(KIND(1D0)), INTENT(IN)::FAIEveTree
155  REAL(KIND(1D0)), INTENT(IN)::Faut
156  REAL(KIND(1D0)), INTENT(IN)::fcld_obs
157  REAL(KIND(1D0)), INTENT(IN)::FlowChange
158  REAL(KIND(1D0)), INTENT(IN)::FrFossilFuel_Heat
159  REAL(KIND(1D0)), INTENT(IN)::FrFossilFuel_NonHeat
160  REAL(KIND(1D0)), INTENT(IN)::G1
161  REAL(KIND(1D0)), INTENT(IN)::G2
162  REAL(KIND(1D0)), INTENT(IN)::G3
163  REAL(KIND(1D0)), INTENT(IN)::G4
164  REAL(KIND(1D0)), INTENT(IN)::G5
165  REAL(KIND(1D0)), INTENT(IN)::G6
166  REAL(KIND(1D0)), INTENT(IN)::H_maintain
167  REAL(KIND(1D0)), INTENT(IN)::InternalWaterUse_h
168  REAL(KIND(1D0)), INTENT(IN)::IrrFracPaved
169  REAL(KIND(1D0)), INTENT(IN)::IrrFracBldgs
170  REAL(KIND(1D0)), INTENT(IN)::IrrFracEveTr
171  REAL(KIND(1D0)), INTENT(IN)::IrrFracDecTr
172  REAL(KIND(1D0)), INTENT(IN)::IrrFracGrass
173  REAL(KIND(1D0)), INTENT(IN)::IrrFracBSoil
174  REAL(KIND(1D0)), INTENT(IN)::IrrFracWater
175  REAL(KIND(1D0)), INTENT(IN)::Kmax
176  REAL(KIND(1D0)), INTENT(IN)::LAI_obs
177  REAL(KIND(1D0)), INTENT(IN)::lat
178  REAL(KIND(1D0)), INTENT(IN)::ldown_obs
179  REAL(KIND(1D0)), INTENT(IN)::lng
180  REAL(KIND(1D0)), INTENT(IN)::MaxFCMetab
181  REAL(KIND(1D0)), INTENT(IN)::MaxQFMetab
182  REAL(KIND(1D0)), INTENT(IN)::MinFCMetab
183  REAL(KIND(1D0)), INTENT(IN)::MinQFMetab
184  REAL(KIND(1D0)), INTENT(IN)::NARP_EMIS_SNOW
185  REAL(KIND(1D0)), INTENT(IN)::NARP_TRANS_SITE
186  REAL(KIND(1D0)), INTENT(IN)::PipeCapacity
187  REAL(KIND(1D0)), INTENT(IN)::PopDensNighttime
188  REAL(KIND(1D0)), INTENT(IN)::PorMax_dec
189  REAL(KIND(1D0)), INTENT(IN)::PorMin_dec
190  REAL(KIND(1D0)), INTENT(IN)::Precip
191  REAL(KIND(1D0)), INTENT(IN)::PrecipLimit
192  REAL(KIND(1D0)), INTENT(IN)::PrecipLimitAlb
193  REAL(KIND(1D0)), INTENT(IN)::Press_hPa
194  REAL(KIND(1D0)), INTENT(IN)::qn1_obs
195  REAL(KIND(1D0)), INTENT(IN)::qs_obs
196  REAL(KIND(1D0)), INTENT(IN)::qf_obs
197  REAL(KIND(1D0)), INTENT(IN)::RadMeltFact
198  REAL(KIND(1D0)), INTENT(IN)::RAINCOVER
199  REAL(KIND(1D0)), INTENT(IN)::RainMaxRes
200  REAL(KIND(1D0)), INTENT(IN)::RunoffToWater
201  REAL(KIND(1D0)), INTENT(IN)::S1
202  REAL(KIND(1D0)), INTENT(IN)::S2
203  REAL(KIND(1D0)), INTENT(IN)::SnowAlbMax
204  REAL(KIND(1D0)), INTENT(IN)::SnowAlbMin
205  REAL(KIND(1D0)), INTENT(IN)::SnowDensMax
206  REAL(KIND(1D0)), INTENT(IN)::SnowDensMin
207  REAL(KIND(1D0)), INTENT(IN)::SnowLimBldg
208  REAL(KIND(1D0)), INTENT(IN)::SnowLimPaved
209  REAL(KIND(1D0)), INTENT(IN)::snowFrac_obs
210  REAL(KIND(1D0)), INTENT(IN)::SurfaceArea
211  REAL(KIND(1D0)), INTENT(IN)::tau_a
212  REAL(KIND(1D0)), INTENT(IN)::tau_f
213  REAL(KIND(1D0)), INTENT(IN)::tau_r
214  REAL(KIND(1D0)), INTENT(IN)::Temp_C
215  REAL(KIND(1D0)), INTENT(IN)::TempMeltFact
216  REAL(KIND(1D0)), INTENT(IN)::TH
217  REAL(KIND(1D0)), INTENT(IN)::timezone
218  REAL(KIND(1D0)), INTENT(IN)::TL
219  REAL(KIND(1D0)), INTENT(IN)::TrafficUnits
220  REAL(KIND(1D0)), INTENT(IN)::wu_m3
221  REAL(KIND(1D0)), INTENT(IN)::xsmd
222  REAL(KIND(1D0)), INTENT(IN)::Z
223  REAL(KIND(1D0)), INTENT(IN)::z0m_in
224  REAL(KIND(1D0)), INTENT(IN)::zdm_in
225 
226  INTEGER, DIMENSION(NVEGSURF), INTENT(IN)::LAIType
227 
228  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::AH_MIN
229  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::AH_SLOPE_Cooling
230  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::AH_SLOPE_Heating
231  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::FcEF_v_kgkm
232  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::QF0_BEU
233  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::Qf_A
234  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::Qf_B
235  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::Qf_C
236  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::PopDensDaytime
237  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::BaseT_Cooling
238  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::BaseT_Heating
239  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::TrafficRate
240  REAL(KIND(1D0)), DIMENSION(3), INTENT(IN) ::Ie_a
241  REAL(KIND(1D0)), DIMENSION(3), INTENT(IN) ::Ie_m
242  REAL(KIND(1D0)), DIMENSION(3), INTENT(IN) ::MaxConductance
243  REAL(KIND(1D0)), DIMENSION(7), INTENT(IN) ::DayWat
244  REAL(KIND(1D0)), DIMENSION(7), INTENT(IN) ::DayWatPer
245  REAL(KIND(1D0)), DIMENSION(nsurf + 1), INTENT(IN) ::OHM_threshSW
246  REAL(KIND(1D0)), DIMENSION(nsurf + 1), INTENT(IN) ::OHM_threshWD
247  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::chAnOHM
248  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::cpAnOHM
249  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::emis
250  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::kkAnOHM
251  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::SatHydraulicConduct
252  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::sfr
253  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::SnowPackLimit
254  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::SoilDepth
255  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::SoilStoreCap
256  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::StateLimit
257  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::WetThresh
258  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::alpha_bioCO2
259  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::alpha_enh_bioCO2
260  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::BaseT
261  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::BaseTe
262  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::beta_bioCO2
263  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::beta_enh_bioCO2
264  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::GDDFull
265  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::LAIMax
266  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::LAIMin
267  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::min_res_bioCO2
268  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::resp_a
269  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::resp_b
270  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::SDDFull
271  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::SnowProf_24hr
272  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::theta_bioCO2
273  REAL(KIND(1D0)), DIMENSION(4, NVEGSURF), INTENT(IN) ::LAIPower
274  REAL(KIND(1D0)), DIMENSION(nsurf + 1, 4, 3), INTENT(IN) ::OHM_coef
275  REAL(KIND(1D0)), DIMENSION(NSURF + 1, NSURF - 1), INTENT(IN)::WaterDist
276  REAL(KIND(1d0)), DIMENSION(:), INTENT(IN) ::Ts5mindata_ir
277  REAL(KIND(1D0)), DIMENSION(:, :), INTENT(IN) ::MetForcingData_grid
278 
279  ! diurnal profile values for 24hr
280  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::AHProf_24hr
281  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::HumActivity_24hr
282  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::PopProf_24hr
283  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::TraffProf_24hr
284  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::WUProfA_24hr
285  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::WUProfM_24hr
286 
287  ! ########################################################################################
288 
289  ! ########################################################################################
290  ! inout variables
291  ! OHM related:
292  REAL(KIND(1d0)), INTENT(INOUT) ::qn1_av
293  REAL(KIND(1d0)), INTENT(INOUT) ::dqndt
294  REAL(KIND(1d0)), INTENT(INOUT) ::qn1_s_av
295  REAL(KIND(1d0)), INTENT(INOUT) ::dqnsdt
296 
297  ! snow related:
298  REAL(KIND(1D0)), INTENT(INOUT) ::SnowfallCum
299  REAL(KIND(1D0)), INTENT(INOUT) ::SnowAlb
300  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT)::IceFrac
301  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT)::SnowWater
302  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT)::SnowDens
303  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT)::SnowFrac
304  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT)::SnowPack
305 
306  ! water balance related:
307  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT) ::soilstore_id
308  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT) ::state_id
309  REAL(KIND(1D0)), DIMENSION(6, NSURF), INTENT(INOUT)::StoreDrainPrm
310 
311  ! phenology related:
312  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT) ::alb
313  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(INOUT)::GDD_id !Growing Degree Days (see SUEWS_DailyState.f95)
314  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(INout)::SDD_id !Senescence Degree Days (see SUEWS_DailyState.f95)
315  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(INOUT)::LAI_id !LAI for each veg surface [m2 m-2]
316  REAL(KIND(1d0)), INTENT(INout)::Tmin_id
317  REAL(KIND(1d0)), INTENT(INout)::Tmax_id
318  REAL(KIND(1d0)), INTENT(INout)::lenDay_id
319  REAL(KIND(1d0)), INTENT(INOUT)::DecidCap_id
320  REAL(KIND(1d0)), INTENT(INOUT)::albDecTr_id
321  REAL(KIND(1d0)), INTENT(INOUT)::albEveTr_id
322  REAL(KIND(1d0)), INTENT(INOUT)::albGrass_id
323  REAL(KIND(1d0)), INTENT(INOUT)::porosity_id
324 
325  ! anthropogenic heat related:
326  REAL(KIND(1d0)), DIMENSION(12), INTENT(INOUT) ::HDD_id !Heating Degree Days (see SUEWS_DailyState.f95)
327 
328  ! water use related:
329  REAL(KIND(1d0)), DIMENSION(9), INTENT(INOUT) ::WUDay_id
330 
331  ! ESTM related:
332  REAL(KIND(1d0)), INTENT(INOUT) ::Tair_av
333  ! ########################################################################################
334 
335  ! ########################################################################################
336  ! output variables
337  REAL(KIND(1D0)), DIMENSION(5), INTENT(OUT) ::datetimeLine
338  REAL(KIND(1D0)), DIMENSION(ncolumnsDataOutSUEWS - 5), INTENT(OUT) ::dataOutLineSUEWS
339  REAL(KIND(1D0)), DIMENSION(ncolumnsDataOutSnow - 5), INTENT(OUT) ::dataOutLineSnow
340  REAL(KIND(1d0)), DIMENSION(ncolumnsDataOutESTM - 5), INTENT(OUT) ::dataOutLineESTM
341  REAL(KIND(1d0)), DIMENSION(ncolumnsDataOutRSL - 5), INTENT(OUT) ::dataoutLineRSL ! RSL variable array
342  REAL(KIND(1D0)), DIMENSION(ncolumnsDataOutSol - 5), INTENT(OUT) ::dataOutLineSOLWEIG
343  REAL(KIND(1d0)), DIMENSION(ncolumnsDataOutDailyState - 5), INTENT(OUT)::DailyStateLine
344  ! ########################################################################################
345 
346  ! ########################################################################################
347  ! local variables
348  REAL(KIND(1D0))::a1
349  REAL(KIND(1D0))::a2
350  REAL(KIND(1D0))::a3
351  REAL(KIND(1D0))::AdditionalWater
352  REAL(KIND(1D0))::U10_ms
353  REAL(KIND(1D0))::azimuth
354  REAL(KIND(1D0))::chSnow_per_interval
355 
356  REAL(KIND(1D0))::dens_dry
357  REAL(KIND(1d0))::deltaLAI
358  REAL(KIND(1D0))::drain_per_tstep
359  REAL(KIND(1D0))::Ea_hPa
360  REAL(KIND(1D0))::QE_LUMPS
361  REAL(KIND(1D0))::es_hPa
362  REAL(KIND(1D0))::ev_per_tstep
363  REAL(KIND(1D0))::wu_ext
364  REAL(KIND(1D0))::Fc
365  REAL(KIND(1D0))::Fc_anthro
366  REAL(KIND(1D0))::Fc_biogen
367  REAL(KIND(1D0))::Fc_build
368  REAL(KIND(1D0))::fcld
369  REAL(KIND(1D0))::Fc_metab
370  REAL(KIND(1D0))::Fc_photo
371  REAL(KIND(1D0))::Fc_point
372  REAL(KIND(1D0))::Fc_respi
373  REAL(KIND(1D0))::Fc_traff
374  REAL(KIND(1D0))::gfunc
375  REAL(KIND(1D0))::gsc
376  REAL(KIND(1D0))::QH_LUMPS
377  REAL(KIND(1D0))::wu_int
378  REAL(KIND(1D0))::kclear
379  REAL(KIND(1D0))::kup
380  REAL(KIND(1D0))::ldown
381  REAL(KIND(1D0))::lup
382  REAL(KIND(1D0))::L_mod
383  REAL(KIND(1D0))::mwh
384  REAL(KIND(1D0))::mwstore
385  REAL(KIND(1D0))::NWstate_per_tstep
386  REAL(KIND(1D0))::planF
387  REAL(KIND(1D0))::zL
388  REAL(KIND(1D0))::q2_gkg
389  REAL(KIND(1D0))::qe
390  REAL(KIND(1D0))::qf
391  REAL(KIND(1D0))::QF_SAHP
392  REAL(KIND(1D0))::qh
393  REAL(KIND(1D0))::qh_residual
394  REAL(KIND(1D0))::qh_resist
395  REAL(KIND(1D0))::Qm
396  REAL(KIND(1D0))::QmFreez
397  REAL(KIND(1D0))::QmRain
398  REAL(KIND(1D0))::qn1
399  REAL(KIND(1D0))::qn1_S
400  REAL(KIND(1D0))::qn1_snowfree
401  REAL(KIND(1D0))::qs
402  REAL(KIND(1D0))::RA
403  REAL(KIND(1D0))::ResistSurf
404  REAL(KIND(1D0))::RH2
405  REAL(KIND(1d0))::runoffAGveg
406  REAL(KIND(1d0))::runoffAGimpervious
407  REAL(KIND(1D0))::runoff_per_tstep
408  REAL(KIND(1D0))::runoffPipes
409  REAL(KIND(1D0))::runoffSoil_per_tstep
410  REAL(KIND(1D0))::runoffwaterbody
411  REAL(KIND(1D0))::smd
412  REAL(KIND(1D0))::SoilState
413  REAL(KIND(1D0))::state_per_tstep
414  REAL(KIND(1D0))::surf_chang_per_tstep
415  REAL(KIND(1D0))::swe
416  REAL(KIND(1D0))::t2_C
417  REAL(KIND(1D0))::TSfc_C
418  REAL(KIND(1D0))::TempVeg
419  REAL(KIND(1D0))::tot_chang_per_tstep
420  REAL(KIND(1D0))::TStar
421  REAL(KIND(1D0))::tsurf
422  REAL(KIND(1D0))::UStar
423  REAL(KIND(1D0))::VPD_Pa
424  ! REAL(KIND(1D0))::wu_DecTr
425  ! REAL(KIND(1D0))::wu_EveTr
426  ! REAL(KIND(1D0))::wu_Grass
427  REAL(KIND(1D0))::z0m
428  REAL(KIND(1D0))::zdm
429  REAL(KIND(1D0))::ZENITH_deg
430  REAL(KIND(1D0))::Zh
431 
432  REAL(KIND(1D0)), DIMENSION(2)::SnowRemoval
433  REAL(KIND(1D0)), DIMENSION(NSURF)::wu_nsurf
434  REAL(KIND(1D0)), DIMENSION(NSURF)::FreezMelt
435  REAL(KIND(1d0)), DIMENSION(nsurf)::kup_ind_snow
436  REAL(KIND(1D0)), DIMENSION(NSURF)::mw_ind
437  REAL(KIND(1D0)), DIMENSION(NSURF)::Qm_freezState
438  REAL(KIND(1D0)), DIMENSION(NSURF)::Qm_melt
439  REAL(KIND(1D0)), DIMENSION(NSURF)::Qm_rain
440  REAL(KIND(1D0)), DIMENSION(NSURF)::qn1_ind_snow
441  REAL(KIND(1D0)), DIMENSION(NSURF)::rainOnSnow
442  REAL(KIND(1D0)), DIMENSION(NSURF)::runoffSoil
443  REAL(KIND(1D0)), DIMENSION(NSURF)::smd_nsurf
444  REAL(KIND(1D0)), DIMENSION(NSURF)::snowDepth
445 
446  REAL(KIND(1d0)), DIMENSION(nsurf)::Tsurf_ind_snow
447 
448  INTEGER, DIMENSION(NSURF)::snowCalcSwitch
449  INTEGER, DIMENSION(3) ::dayofWeek_id
450  INTEGER::DLS
451 
452  ! REAL(KIND(1D0))::avcp
453  ! REAL(KIND(1D0))::avdens
454  ! REAL(KIND(1D0))::lv_J_kg
455  REAL(KIND(1D0))::dq
456  REAL(KIND(1D0))::lvS_J_kg
457  REAL(KIND(1D0))::psyc_hPa
458  REAL(KIND(1D0))::RAsnow
459  REAL(KIND(1D0))::rb
460  REAL(KIND(1D0))::runoff_per_interval
461  REAL(KIND(1D0))::s_hPa
462  REAL(KIND(1D0))::sIce_hpa
463  REAL(KIND(1D0))::SoilMoistCap
464  REAL(KIND(1D0))::veg_fr
465  REAL(KIND(1D0))::VegPhenLumps
466  REAL(KIND(1D0))::VPd_hpa
467  REAL(KIND(1D0))::vsmd
468  REAL(KIND(1D0))::ZZD
469 
470  REAL(KIND(1D0)), DIMENSION(NSURF)::deltaQi
471  REAL(KIND(1D0)), DIMENSION(NSURF)::drain
472  REAL(KIND(1D0)), DIMENSION(NSURF)::FreezState
473  REAL(KIND(1D0)), DIMENSION(NSURF)::FreezStateVol
474  REAL(KIND(1D0)), DIMENSION(NSURF)::soilstoreOld
475  REAL(KIND(1D0)), DIMENSION(NSURF)::stateOld
476  REAL(KIND(1D0)), DIMENSION(NSURF)::tsurf_ind
477 
478  ! TODO: TS 25 Oct 2017
479  ! the variables are not used currently as grid-to-grid connection is NOT set up.
480  ! set these variables as zero.
481  REAL(KIND(1D0))::addImpervious = 0
482  REAL(KIND(1D0))::addPipes = 0
483  REAL(KIND(1D0))::addVeg = 0
484  REAL(KIND(1D0))::addWaterBody = 0
485  REAL(KIND(1D0)), DIMENSION(NSURF)::AddWater = 0
486  REAL(KIND(1D0)), DIMENSION(NSURF)::frac_water2runoff = 0
487 
488  ! values that are derived from tstep
489  INTEGER::nsh ! number of timesteps per hour
490  REAL(KIND(1D0))::nsh_real ! nsh in type real
491  REAL(KIND(1D0))::tstep_real ! tstep in type real
492  REAL(KIND(1D0))::dectime
493 
494  ! values that are derived from sfr (surface fractions)
495  REAL(KIND(1D0))::VegFraction
496  REAL(KIND(1D0))::ImpervFraction
497  REAL(KIND(1D0))::PervFraction
498  REAL(KIND(1D0))::NonWaterFraction
499 
500  ! snow related temporary values
501  ! REAL(KIND(1D0))::alb0
502  REAL(KIND(1D0))::alb1
503 
504  ! ########################################################################################
505  ! TS 19 Sep 2019
506  ! temporary variables to save values for inout varialbes
507  ! suffixes and denote values from last and to next tsteps, respectively
508  ! these variables are introduced to allow safe and robust iterations inccurred in this subroutine
509  ! so that these values won't updated in unexpectedly many times
510 
511  ! OHM related:
512  REAL(KIND(1D0))::qn1_av_prev, qn1_av_next
513  REAL(KIND(1D0))::dqndt_prev, dqndt_next
514  REAL(KIND(1D0))::qn1_s_av_prev, qn1_s_av_next
515  REAL(KIND(1D0))::dqnsdt_prev, dqnsdt_next
516 
517  ! snow related:
518  REAL(KIND(1D0))::SnowfallCum_prev, SnowfallCum_next
519  REAL(KIND(1D0))::SnowAlb_prev, SnowAlb_next
520 
521  REAL(KIND(1D0)), DIMENSION(NSURF)::IceFrac_prev, IceFrac_next
522  REAL(KIND(1D0)), DIMENSION(NSURF)::SnowWater_prev, SnowWater_next
523  REAL(KIND(1D0)), DIMENSION(NSURF)::SnowDens_prev, SnowDens_next
524  REAL(KIND(1D0)), DIMENSION(NSURF)::SnowFrac_prev, SnowFrac_next
525  REAL(KIND(1D0)), DIMENSION(NSURF)::SnowPack_prev, SnowPack_next
526 
527  ! water balance related:
528  REAL(KIND(1D0)), DIMENSION(NSURF) ::soilstore_id_prev, soilstore_id_next
529  REAL(KIND(1D0)), DIMENSION(NSURF) ::state_id_prev, state_id_next
530  REAL(KIND(1D0)), DIMENSION(6, NSURF)::StoreDrainPrm_prev, StoreDrainPrm_next
531 
532  ! phenology related:
533  REAL(KIND(1D0)), DIMENSION(NSURF) ::alb_prev, alb_next
534  REAL(KIND(1d0)), DIMENSION(nvegsurf)::GDD_id_prev, GDD_id_next
535  REAL(KIND(1d0)), DIMENSION(nvegsurf)::LAI_id_prev, LAI_id_next
536  REAL(KIND(1d0)), DIMENSION(nvegsurf)::SDD_id_prev, SDD_id_next
537 
538  REAL(KIND(1D0))::DecidCap_id_prev, DecidCap_id_next
539  REAL(KIND(1D0))::albDecTr_id_prev, albDecTr_id_next
540  REAL(KIND(1D0))::albEveTr_id_prev, albEveTr_id_next
541  REAL(KIND(1D0))::albGrass_id_prev, albGrass_id_next
542  REAL(KIND(1D0))::porosity_id_prev, porosity_id_next
543 
544  REAL(KIND(1d0))::Tmin_id_prev, Tmin_id_next
545  REAL(KIND(1d0))::Tmax_id_prev, Tmax_id_next
546  REAL(KIND(1d0))::lenDay_id_prev, lenDay_id_next
547 
548  ! anthropogenic heat related:
549  REAL(KIND(1d0)), DIMENSION(12)::HDD_id_prev, HDD_id_next
550 
551  ! water use related:
552  REAL(KIND(1d0)), DIMENSION(9)::WUDay_id_prev, WUDay_id_next
553 
554  REAL(KIND(1D0))::Tair_av_prev, Tair_av_next
555  ! ########################################################################################
556 
557  ! Related to RSL wind profiles
558  INTEGER, PARAMETER :: nz = 90 ! number of levels 10 levels in canopy plus 20 (3 x Zh) above the canopy
559 
560  ! flag for Tsurf convergence
561  logical:: flag_converge
562  REAL(KIND(1D0)):: Ts_iter
563  REAL(KIND(1D0)):: L_mod_iter
564  REAL(KIND(1D0)):: QH_Init
565  INTEGER::i_iter
566 
567  ! REAL(KIND(1d0)), DIMENSION(30):: psihatm_z
568  ! REAL(KIND(1d0)), DIMENSION(30):: psihath_z
569 
570  ! ########################################################################################
571  ! save initial values of inout variables
572  qn1_av_prev = qn1_av
573  dqndt_prev = dqndt
574  qn1_s_av_prev = qn1_s_av
575  dqnsdt_prev = dqnsdt
576  snowfallcum_prev = snowfallcum
577  snowalb_prev = snowalb
578  icefrac_prev = icefrac
579  snowwater_prev = snowwater
580  snowdens_prev = snowdens
581  snowfrac_prev = snowfrac
582  snowpack_prev = snowpack
583  soilstore_id_prev = soilstore_id
584  state_id_prev = state_id
585  tair_av_prev = tair_av
586  lai_id_prev = lai_id
587  gdd_id_prev = gdd_id
588  sdd_id_prev = sdd_id
589  tmin_id_prev = tmin_id
590  tmax_id_prev = tmax_id
591  lenday_id_prev = lenday_id
592  storedrainprm_prev = storedrainprm
593  decidcap_id_prev = decidcap_id
594  porosity_id_prev = porosity_id
595  alb_prev = alb
596  albdectr_id_prev = albdectr_id
597  albevetr_id_prev = albevetr_id
598  albgrass_id_prev = albgrass_id
599  hdd_id_prev = hdd_id
600  wuday_id_prev = wuday_id
601 
602  ! initialise variables
603  qn1_av_next = qn1_av
604  dqndt_next = dqndt
605  qn1_s_av_next = qn1_s_av
606  dqnsdt_next = dqnsdt
607  snowfallcum_next = snowfallcum
608  snowalb_next = snowalb
609  icefrac_next = icefrac
610  snowwater_next = snowwater
611  snowdens_next = snowdens
612  snowfrac_next = snowfrac
613  snowpack_next = snowpack
614  soilstore_id_next = soilstore_id
615  state_id_next = state_id
616  tair_av_next = tair_av
617  lai_id_next = lai_id
618  gdd_id_next = gdd_id
619  sdd_id_next = sdd_id
620  tmin_id_next = tmin_id
621  tmax_id_next = tmax_id
622  lenday_id_next = lenday_id
623  storedrainprm_next = storedrainprm
624  decidcap_id_next = decidcap_id
625  porosity_id_next = porosity_id
626  alb_next = alb
627  albdectr_id_next = albdectr_id
628  albevetr_id_next = albevetr_id
629  albgrass_id_next = albgrass_id
630  hdd_id_next = hdd_id
631  wuday_id_next = wuday_id
632 
633  !########################################################################################
634  ! main calculation starts here
635  !########################################################################################
636 
637  ! iteration is used below to get results converge
638  flag_converge = .false.
639  ts_iter = temp_c
640  l_mod_iter = 10
641  i_iter = 1
642  do while (.not. flag_converge)
643 
644  ! calculate dectime
645  CALL suews_cal_dectime( &
646  id, it, imin, isec, & ! input
647  dectime) ! output
648 
649  ! calculate tstep related VARIABLES
650  CALL suews_cal_tstep( &
651  tstep, & ! input
652  nsh, nsh_real, tstep_real) ! output
653 
654  ! calculate surface fraction related VARIABLES
655  CALL suews_cal_surf( &
656  sfr, & !input
657  vegfraction, impervfraction, pervfraction, nonwaterfraction) ! output
658 
659  ! calculate dayofweek information
660  CALL suews_cal_weekday( &
661  iy, id, lat, & !input
662  dayofweek_id) !output
663 
664  ! calculate dayofweek information
665  CALL suews_cal_dls( &
666  id, startdls, enddls, & !input
667  dls) !output
668 
669  ! calculate mean air temperature of past 24 hours
670  tair_av_next = cal_tair_av(tair_av_prev, dt_since_start, tstep, temp_c)
671 
672  !==============main calculation start=======================
673 
674  !==============surface roughness calculation=======================
675  IF (diagnose == 1) WRITE (*, *) 'Calling SUEWS_cal_RoughnessParameters...'
676  IF (diagnose == 1) print *, 'z0m_in =', z0m_in
678  roughlenmommethod, sfr, &!input
679  bldgh, evetreeh, dectreeh, &
680  porosity_id_prev, faibldg, faievetree, faidectree, &
681  z0m_in, zdm_in, z, &
682  planf, &!output
683  zh, z0m, zdm, zzd)
684 
685  !=================Calculate sun position=================
686  IF (diagnose == 1) WRITE (*, *) 'Calling NARP_cal_SunPosition...'
687  CALL narp_cal_sunposition( &
688  REAL(iy, KIND(1d0)), &!input:
689  dectime - tstep/2, &! sun position at middle of timestep before
690  timezone, lat, lng, alt, &
691  azimuth, zenith_deg)!output:
692 
693  !=================Call the SUEWS_cal_DailyState routine to get surface characteristics ready=================
694  IF (diagnose == 1) WRITE (*, *) 'Calling SUEWS_cal_DailyState...'
695  CALL suews_cal_dailystate( &
696  iy, id, it, imin, isec, tstep, tstep_prev, dt_since_start, dayofweek_id, &!input
697  tmin_id_prev, tmax_id_prev, lenday_id_prev, &
698  basetmethod, &
699  waterusemethod, ie_start, ie_end, &
700  laicalcyes, laitype, &
701  nsh_real, avkdn, temp_c, precip, baset_hc, &
702  baset_heating, baset_cooling, &
703  lat, faut, lai_obs, &
704  albmax_dectr, albmax_evetr, albmax_grass, &
705  albmin_dectr, albmin_evetr, albmin_grass, &
706  capmax_dec, capmin_dec, pormax_dec, pormin_dec, &
707  ie_a, ie_m, daywatper, daywat, &
708  baset, basete, gddfull, sddfull, laimin, laimax, laipower, &
709  decidcap_id_prev, storedrainprm_prev, lai_id_prev, gdd_id_prev, sdd_id_prev, &
710  albdectr_id_prev, albevetr_id_prev, albgrass_id_prev, porosity_id_prev, &!input
711  hdd_id_prev, &!input
712  state_id_prev, soilstore_id_prev, soilstorecap, h_maintain, &!input
713  hdd_id_next, &!output
714  tmin_id_next, tmax_id_next, lenday_id_next, &
715  albdectr_id_next, albevetr_id_next, albgrass_id_next, porosity_id_next, &!output
716  decidcap_id_next, storedrainprm_next, lai_id_next, gdd_id_next, sdd_id_next, deltalai, wuday_id_next)!output
717 
718  !=================Calculation of density and other water related parameters=================
719  IF (diagnose == 1) WRITE (*, *) 'Calling LUMPS_cal_AtmMoist...'
720  CALL cal_atmmoist( &
721  temp_c, press_hpa, avrh, dectime, &! input:
722  lv_j_kg, lvs_j_kg, &! output:
723  es_hpa, ea_hpa, vpd_hpa, vpd_pa, dq, dens_dry, avcp, avdens)
724 
725  !======== Calculate soil moisture =========
726  IF (diagnose == 1) WRITE (*, *) 'Calling SUEWS_update_SoilMoist...'
727  CALL suews_update_soilmoist( &
728  nonwaterfraction, &!input
729  soilstorecap, sfr, soilstore_id_prev, &
730  soilmoistcap, soilstate, &!output
731  vsmd, smd)
732 
733  IF (diagnose == 1) WRITE (*, *) 'Calling SUEWS_cal_WaterUse...'
734  !=================Gives the external and internal water uses per timestep=================
735  CALL suews_cal_wateruse( &
736  nsh_real, & ! input:
737  wu_m3, surfacearea, sfr, &
738  irrfracpaved, irrfracbldgs, &
739  irrfracevetr, irrfracdectr, irrfracgrass, &
740  irrfracbsoil, irrfracwater, &
741  dayofweek_id, wuprofa_24hr, wuprofm_24hr, &
742  internalwateruse_h, hdd_id_next, wuday_id_next, &
743  waterusemethod, nsh, it, imin, dls, &
744  wu_nsurf, wu_int, wu_ext)! output:
745 
746  ! ===================ANTHROPOGENIC HEAT AND CO2 FLUX======================
748  ah_min, ahprof_24hr, ah_slope_cooling, ah_slope_heating, co2pointsource, &! input:
749  dayofweek_id, dls, ef_umolco2perj, emissionsmethod, enef_v_jkm, &
750  fcef_v_kgkm, frfossilfuel_heat, frfossilfuel_nonheat, hdd_id_next, humactivity_24hr, &
751  id, imin, it, maxfcmetab, maxqfmetab, minfcmetab, minqfmetab, nsh, &
752  popdensdaytime, popdensnighttime, popprof_24hr, qf, qf0_beu, qf_a, qf_b, qf_c, &
753  qf_obs, qf_sahp, surfacearea, baset_cooling, baset_heating, &
754  temp_c, trafficrate, trafficunits, traffprof_24hr, &
755  fc_anthro, fc_build, fc_metab, fc_point, fc_traff)! output:
756 
757  ! ========================================================================
758  ! N.B.: the following parts involves snow-related calculations.
759 
760  ! ===================NET ALLWAVE RADIATION================================
761  CALL suews_cal_qn( &
762  netradiationmethod, snowuse, &!input
763  tstep, snowpack_prev, tau_a, tau_f, snowalbmax, snowalbmin, &
764  diagnose, snowfrac_obs, ldown_obs, fcld_obs, &
765  dectime, zenith_deg, ts_iter, avkdn, temp_c, avrh, ea_hpa, qn1_obs, &
766  snowalb_prev, snowfrac_prev, diagqn, &
767  narp_trans_site, narp_emis_snow, icefrac_prev, sfr, emis, &
768  alb_prev, albdectr_id_next, albevetr_id_next, albgrass_id_next, &!input
769  alb_next, &!output
770  ldown, fcld, &!output
771  qn1, qn1_snowfree, qn1_s, kclear, kup, lup, tsurf, &
772  qn1_ind_snow, kup_ind_snow, tsurf_ind_snow, tsurf_ind, &
773  alb1, snowfrac_next, snowalb_next)
774 
775  ! =================STORAGE HEAT FLUX=======================================
776  CALL suews_cal_qs( &
777  storageheatmethod, qs_obs, ohmincqf, gridiv, &!input
778  id, tstep, dt_since_start, diagnose, sfr, &
779  ohm_coef, ohm_threshsw, ohm_threshwd, &
780  soilstore_id_prev, soilstorecap, state_id_prev, snowuse, snowfrac_next, diagqs, &
781  hdd_id_next, metforcingdata_grid, ts5mindata_ir, qf, qn1, &
782  avkdn, avu1, temp_c, zenith_deg, avrh, press_hpa, ldown, &
783  bldgh, alb_next, emis, cpanohm, kkanohm, chanohm, emissionsmethod, &
784  tair_av_next, qn1_av_prev, dqndt_prev, qn1_s_av_prev, dqnsdt_prev, &
785  storedrainprm_next, &
786  qn1_s, dataoutlineestm, qs, &!output
787  qn1_av_next, dqndt_next, qn1_s_av_next, dqnsdt_next, &
788  deltaqi, a1, a2, a3)
789 
790  !==================Energy related to snow melting/freezing processes=======
791  IF (diagnose == 1) WRITE (*, *) 'Calling MeltHeat'
792 
793  CALL snow_cal_meltheat( &
794  snowuse, &!input
795  tstep, tau_r, snowdensmax, &
796  lvs_j_kg, lv_j_kg, tstep_real, radmeltfact, tempmeltfact, snowalbmax, &
797  snowdensmin, temp_c, precip, preciplimit, preciplimitalb, &
798  nsh_real, sfr, tsurf_ind, tsurf_ind_snow, state_id_prev, qn1_ind_snow, &
799  kup_ind_snow, snowwater_prev, deltaqi, alb1, &
800  snowpack_prev, snowfrac_next, snowalb_next, snowdens_prev, snowfallcum_prev, &!input
801  snowpack_next, snowfrac_next, snowalb_next, snowdens_next, snowfallcum_next, &!output
802  mwh, qm, qmfreez, qmrain, &! output
803  veg_fr, snowcalcswitch, qm_melt, qm_freezstate, qm_rain, freezmelt, &
804  freezstate, freezstatevol, rainonsnow, snowdepth, mw_ind, &
805  dataoutlinesnow)!output
806 
807  !==========================Turbulent Fluxes================================
808  IF (diagnose == 1) WRITE (*, *) 'Calling LUMPS_cal_QHQE...'
809  if (i_iter == 1) then
810  !Calculate QH and QE from LUMPS in the first iteration of each time step
811  CALL lumps_cal_qhqe( &
812  veg_type, & !input
813  snowuse, qn1, qf, qs, qm, temp_c, veg_fr, avcp, press_hpa, lv_j_kg, &
814  tstep_real, drainrt, nsh_real, &
815  precip, rainmaxres, raincover, sfr, lai_id_next, laimax, laimin, &
816  qh_lumps, & !output
817  qe_lumps, psyc_hpa, s_hpa, sice_hpa, tempveg, vegphenlumps)
818 
819  ! use LUMPS QH to do stability correction
820  qh_init = qh_lumps
821  else
822  ! use SUEWS QH to do stability correction
823  qh_init = qh
824  end if
825 
826  !============= calculate water balance =============
827  CALL suews_cal_water( &
828  diagnose, &!input
829  snowuse, nonwaterfraction, addpipes, addimpervious, addveg, addwaterbody, &
830  state_id_prev, soilstore_id_prev, sfr, storedrainprm_next, waterdist, nsh_real, &
831  drain_per_tstep, & !output
832  drain, frac_water2runoff, &
833  additionalwater, runoffpipes, runoff_per_interval, &
834  addwater, stateold, soilstoreold)
835  !============= calculate water balance end =============
836 
837  !===============Resistance Calculations=======================
838  CALL suews_cal_resistance( &
839  stabilitymethod, &!input:
840  diagnose, aerodynamicresistancemethod, roughlenheatmethod, snowuse, &
841  id, it, gsmodel, smdmethod, &
842  avdens, avcp, qh_init, zzd, z0m, zdm, &
843  avu1, temp_c, vegfraction, avkdn, &
844  kmax, &
845  g1, g2, g3, g4, &
846  g5, g6, s1, s2, &
847  th, tl, &
848  dq, xsmd, vsmd, maxconductance, laimax, lai_id_next, snowfrac_next, sfr, &
849  ustar, tstar, l_mod, &!output
850  zl, gsc, resistsurf, ra, rasnow, rb)
851 
852  !======== Evaporation and surface state_id ========
853  CALL suews_cal_qe( &
854  diagnose, snowuse, &!input
855  tstep, imin, it, evapmethod, snowcalcswitch, dayofweek_id, crwmin, crwmax, &
856  dectime, avdens, avcp, lv_j_kg, lvs_j_kg, avrh, press_hpa, temp_c, &
857  rasnow, psyc_hpa, sice_hpa, &
858  pervfraction, vegfraction, addimpervious, qn1_snowfree, qf, qs, vpd_hpa, s_hpa, &
859  resistsurf, ra, rb, snowdensmin, precip, pipecapacity, runofftowater, &
860  nonwaterfraction, wu_nsurf, addveg, addwaterbody, snowlimpaved, snowlimbldg, &
861  surfacearea, flowchange, drain, wetthresh, stateold, mw_ind, soilstorecap, rainonsnow, &
862  freezmelt, freezstate, freezstatevol, qm_melt, qm_rain, tsurf_ind, sfr, &
863  statelimit, addwater, frac_water2runoff, storedrainprm_next, snowpacklimit, snowprof_24hr, &
864  snowpack_next, snowfrac_next, snowwater_prev, icefrac_prev, snowdens_next, &! input:
865  runoff_per_interval, state_id_prev, soilstore_id_prev, &! input:
866  state_id_next, soilstore_id_next, &! output:
867  snowpack_next, snowfrac_next, snowwater_next, icefrac_next, snowdens_next, &! output
868  runoffsoil, &! output:
869  snowremoval, &
870  state_per_tstep, nwstate_per_tstep, qe, &
871  swe, chsnow_per_interval, ev_per_tstep, runoff_per_tstep, &
872  surf_chang_per_tstep, runoffpipes, mwstore, runoffwaterbody, &
873  runoffagveg, runoffagimpervious)
874  !======== Evaporation and surface state_id end========
875 
876  !============ Sensible heat flux ===============
877  IF (diagnose == 1) WRITE (*, *) 'Calling SUEWS_cal_QH...'
878  CALL suews_cal_qh( &
879  1, &
880  qn1, qf, qmrain, qe, qs, qmfreez, qm, avdens, avcp, tsurf, temp_c, ra, &
881  qh, qh_residual, qh_resist)!output
882  !============ Sensible heat flux end===============
883 
884  ! N.B.: snow-related calculations end here.
885  !===================================================
886 
887  !=== Horizontal movement between soil stores ===
888  ! Now water is allowed to move horizontally between the soil stores
889  IF (diagnose == 1) WRITE (*, *) 'Calling SUEWS_cal_HorizontalSoilWater...'
891  sfr, &! input: ! surface fractions
892  soilstorecap, &!Capacity of soil store for each surface [mm]
893  soildepth, &!Depth of sub-surface soil store for each surface [mm]
894  sathydraulicconduct, &!Saturated hydraulic conductivity for each soil subsurface [mm s-1]
895  surfacearea, &!Surface area of the study area [m2]
896  nonwaterfraction, &! sum of surface cover fractions for all except water surfaces
897  tstep_real, & !tstep cast as a real for use in calculations
898  soilstore_id_next, &! inout:!Soil moisture of each surface type [mm]
899  runoffsoil, &!Soil runoff from each soil sub-surface [mm]
900  runoffsoil_per_tstep &! output:!Runoff to deep soil per timestep [mm] (for whole surface, excluding water body)
901  )
902 
903  !========== Calculate soil moisture ============
904  IF (diagnose == 1) WRITE (*, *) 'Calling SUEWS_cal_SoilState...'
905  CALL suews_cal_soilstate( &
906  smdmethod, xsmd, nonwaterfraction, soilmoistcap, &!input
907  soilstorecap, surf_chang_per_tstep, &
908  soilstore_id_next, soilstoreold, sfr, &
909  smd, smd_nsurf, tot_chang_per_tstep, soilstate)!output
910 
911  !============ surface-level diagonostics ===============
912  ! IF (Diagnose == 1) WRITE (*, *) 'Calling SUEWS_cal_Diagnostics...'
913  ! CALL SUEWS_cal_Diagnostics( &
914  ! dectime, &!input
915  ! avU1, Temp_C, avRH, Press_hPa, &
916  ! qh, qe, &
917  ! VegFraction, z, z0m, zdm, RA, avdens, avcp, lv_J_kg, tstep_real, &
918  ! RoughLenHeatMethod, StabilityMethod, &
919  ! U10_ms, t2_C, q2_gkg, tsfc_C, RH2)!output
920 
921  !============ roughness sub-layer diagonostics ===============
922  ! IF (Diagnose == 1) WRITE (*, *) 'Calling RSLProfile...'
923  ! CALL RSLProfile( &
924  ! Zh, z0m, zdm, &
925  ! UStar, L_MOD, sfr, planF, StabilityMethod, &
926  ! avcp, lv_J_kg, &
927  ! Temp_C, avRH, Press_hPa, z, qh, qe, & ! input
928  ! T2_C, q2_gkg, U10_ms, RH2, & !output
929  ! dataoutLineRSL) ! output
930 
931  !============ calculate surface temperature ===============
932  tsfc_c = cal_tsfc(qh, avdens, avcp, ra, temp_c)
933 
934  !============ surface-level diagonostics end ===============
935 
936  ! ============ BIOGENIC CO2 FLUX =======================
937  ! CALL SUEWS_cal_BiogenCO2( &
938  ! alpha_bioCO2, alpha_enh_bioCO2, avkdn, avRh, beta_bioCO2, beta_enh_bioCO2, BSoilSurf, &! input:
939  ! ConifSurf, DecidSurf, dectime, Diagnose, EmissionsMethod, Fc_anthro, G1, G2, G3, G4, &
940  ! G5, G6, gfunc, GrassSurf, gsmodel, id, it, ivConif, ivDecid, ivGrass, Kmax, LAI_id_next, LAIMin, &
941  ! LAIMax, MaxConductance, min_res_bioCO2, nsurf, NVegSurf, Press_hPa, resp_a, &
942  ! resp_b, S1, S2, sfr, SMDMethod, SnowFrac, t2_C, Temp_C, theta_bioCO2, TH, TL, vsmd, xsmd, &
943  ! Fc, Fc_biogen, Fc_photo, Fc_respi)! output:
944 
945  ! force quit do-while, i.e., skip iteration and use NARP for Tsurf calculation
946  ! if (NetRadiationMethod < 10 .or. NetRadiationMethod > 100) exit
947 
948  ! Test if sensible heat fluxes converge in iterations
949  i_iter = i_iter + 1
950  if (abs(qh - qh_init) > 0.1) then
951  flag_converge = .false.
952  else
953  flag_converge = .true.
954  end if
955 
956  ! force quit do-while loop if not convergent after 100 iterations
957  if (i_iter > 100) exit
958 
959  ts_iter = tsfc_c
960  l_mod_iter = l_mod
961 
962  !==============main calculation end=======================
963  ENDDO ! end iteration for tsurf calculations
964 
965  !==============================================================
966  ! Calculate diagnostics: these variables are decoupled from the main SUEWS calculation
967 
968  !============ roughness sub-layer diagonostics ===============
969  IF (diagnose == 1) WRITE (*, *) 'Calling RSLProfile...'
970  CALL rslprofile( &
971  zh, z0m, zdm, &
972  l_mod, sfr, planf, stabilitymethod, &
973  avcp, lv_j_kg, avdens, &
974  avu1, temp_c, avrh, press_hpa, z, qh, qe, & ! input
975  t2_c, q2_gkg, u10_ms, rh2, & !output
976  dataoutlinersl) ! output
977 
978  ! ============ BIOGENIC CO2 FLUX =======================
979  CALL suews_cal_biogenco2( &
980  alpha_bioco2, alpha_enh_bioco2, avkdn, avrh, beta_bioco2, beta_enh_bioco2, bsoilsurf, &! input:
981  conifsurf, decidsurf, dectime, diagnose, emissionsmethod, fc_anthro, g1, g2, g3, g4, &
982  g5, g6, gfunc, grasssurf, gsmodel, id, it, ivconif, ivdecid, ivgrass, kmax, lai_id_next, laimin, &
983  laimax, maxconductance, min_res_bioco2, nsurf, nvegsurf, press_hpa, resp_a, &
984  resp_b, s1, s2, sfr, smdmethod, snowfrac, t2_c, temp_c, theta_bioco2, th, tl, vsmd, xsmd, &
985  fc, fc_biogen, fc_photo, fc_respi)! output:
986 
987  ! calculations of diagnostics end
988  !==============================================================
989 
990  !==============================================================
991  ! update inout variables with new values
992  qn1_av = qn1_av_next
993  dqndt = dqndt_next
994  qn1_s_av = qn1_s_av_next
995  dqnsdt = dqnsdt_next
996  snowfallcum = snowfallcum_next
997  snowalb = snowalb_next
998  icefrac = icefrac_next
999  snowwater = snowwater_next
1000  snowdens = snowdens_next
1001  snowfrac = snowfrac_next
1002  snowpack = snowpack_next
1003 
1004  soilstore_id = soilstore_id_next
1005  state_id = state_id_next
1006  alb = alb_next
1007  gdd_id = gdd_id_next
1008  sdd_id = sdd_id_next
1009  lai_id = lai_id_next
1010  decidcap_id = decidcap_id_next
1011  albdectr_id = albdectr_id_next
1012  albevetr_id = albevetr_id_next
1013  albgrass_id = albgrass_id_next
1014  porosity_id = porosity_id_next
1015  storedrainprm = storedrainprm_next
1016  tair_av = tair_av_next
1017  tmin_id = tmin_id_next
1018  tmax_id = tmax_id_next
1019  lenday_id = lenday_id_next
1020  hdd_id = hdd_id_next
1021  wuday_id = wuday_id_next
1022 
1023  !==============use SOLWEIG to get localised radiation flux==================
1024  if (sfr(bldgsurf) > 0) then
1025  CALL solweig_cal_main(id, it, dectime, 0.8d0, planf, avkdn, ldown, temp_c, avrh, press_hpa, tsfc_c, &
1026  lat, zenith_deg, azimuth, 1.d0, alb(1), alb(2), emis(1), emis(2), bldgh, dataoutlinesolweig)
1027  else
1028  dataoutlinesolweig = set_nan(dataoutlinesolweig)
1029  endif
1030  !==============translation of output variables into output array===========
1031  CALL suews_update_outputline( &
1032  additionalwater, alb, avkdn, u10_ms, azimuth, &!input
1033  chsnow_per_interval, dectime, &
1034  drain_per_tstep, qe_lumps, ev_per_tstep, wu_ext, fc, fc_build, fcld, &
1035  fc_metab, fc_photo, fc_respi, fc_point, fc_traff, flowchange, &
1036  qh_lumps, id, imin, wu_int, it, iy, &
1037  kup, lai_id, ldown, l_mod, lup, mwh, &
1038  mwstore, &
1039  nsh_real, nwstate_per_tstep, precip, q2_gkg, &
1040  qe, qf, qh, qh_resist, qm, qmfreez, &
1041  qmrain, qn1, qn1_s, qn1_snowfree, qs, ra, &
1042  resistsurf, rh2, runoffagimpervious, runoffagveg, &
1043  runoff_per_tstep, runoffpipes, runoffsoil_per_tstep, &
1044  runoffwaterbody, sfr, smd, smd_nsurf, snowalb, snowremoval, &
1045  state_id_next, state_per_tstep, surf_chang_per_tstep, swe, t2_c, tsfc_c, &
1046  tot_chang_per_tstep, tsurf, ustar, &
1047  wu_nsurf, &
1048  z0m, zdm, zenith_deg, &
1049  datetimeline, dataoutlinesuews)!output
1050 
1051  ! model state_id:
1052 
1053  ! daily state_id:
1054  CALL update_dailystateline( &
1055  it, imin, nsh_real, &!input
1056  gdd_id, hdd_id, lai_id, &
1057  sdd_id, &
1058  tmin_id, tmax_id, lenday_id, &
1059  decidcap_id, &
1060  albdectr_id, &
1061  albevetr_id, &
1062  albgrass_id, &
1063  porosity_id, &
1064  wuday_id, &
1065  deltalai, vegphenlumps, &
1066  snowalb, snowdens, &
1067  a1, a2, a3, &
1068  dailystateline)!out
1069 
1070  !==============translation end ================
1071 
1072  END SUBROUTINE suews_cal_main
1073  ! ================================================================================
1074 
1075  ! ===================ANTHROPOGENIC HEAT + CO2 FLUX================================
1076  SUBROUTINE suews_cal_anthropogenicemission( &
1077  AH_MIN, AHProf_24hr, AH_SLOPE_Cooling, AH_SLOPE_Heating, CO2PointSource, &! input:
1078  dayofWeek_id, DLS, EF_umolCO2perJ, EmissionsMethod, EnEF_v_Jkm, &
1079  FcEF_v_kgkm, FrFossilFuel_Heat, FrFossilFuel_NonHeat, HDD_id, HumActivity_24hr, &
1080  id, imin, it, MaxFCMetab, MaxQFMetab, MinFCMetab, MinQFMetab, nsh, &
1081  PopDensDaytime, PopDensNighttime, PopProf_24hr, QF, QF0_BEU, Qf_A, Qf_B, Qf_C, &
1082  QF_obs, QF_SAHP, SurfaceArea, BaseT_Cooling, BaseT_Heating, &
1083  Temp_C, TrafficRate, TrafficUnits, TraffProf_24hr, &
1084  Fc_anthro, Fc_build, Fc_metab, Fc_point, Fc_traff)! output:
1086  IMPLICIT NONE
1087 
1088  ! INTEGER, INTENT(in)::Diagnose
1089  INTEGER, INTENT(in)::DLS
1090  INTEGER, INTENT(in)::EmissionsMethod
1091  INTEGER, INTENT(in)::id
1092  INTEGER, INTENT(in)::it
1093  INTEGER, INTENT(in)::imin
1094  INTEGER, INTENT(in)::nsh
1095  INTEGER, DIMENSION(3), INTENT(in)::dayofWeek_id
1096 
1097  REAL(KIND(1d0)), DIMENSION(6, 2), INTENT(in)::HDD_id
1098 
1099  REAL(KIND(1d0)), DIMENSION(2), INTENT(in)::AH_MIN
1100  REAL(KIND(1d0)), DIMENSION(2), INTENT(in)::AH_SLOPE_Heating
1101  REAL(KIND(1d0)), DIMENSION(2), INTENT(in)::AH_SLOPE_Cooling
1102  REAL(KIND(1D0)), DIMENSION(2), INTENT(in)::FcEF_v_kgkm
1103  ! REAL(KIND(1d0)), DIMENSION(2), INTENT(in)::NumCapita
1104  REAL(KIND(1d0)), DIMENSION(2), INTENT(in)::PopDensDaytime
1105  REAL(KIND(1d0)), DIMENSION(2), INTENT(in)::QF0_BEU
1106  REAL(KIND(1d0)), DIMENSION(2), INTENT(in)::Qf_A
1107  REAL(KIND(1d0)), DIMENSION(2), INTENT(in)::Qf_B
1108  REAL(KIND(1d0)), DIMENSION(2), INTENT(in)::Qf_C
1109  REAL(KIND(1d0)), DIMENSION(2), INTENT(in)::BaseT_Heating
1110  REAL(KIND(1d0)), DIMENSION(2), INTENT(in)::BaseT_Cooling
1111  REAL(KIND(1d0)), DIMENSION(2), INTENT(in)::TrafficRate
1112 
1113  REAL(KIND(1d0)), DIMENSION(0:23, 2), INTENT(in)::AHProf_24hr
1114  REAL(KIND(1d0)), DIMENSION(0:23, 2), INTENT(in)::HumActivity_24hr
1115  REAL(KIND(1d0)), DIMENSION(0:23, 2), INTENT(in)::TraffProf_24hr
1116  REAL(KIND(1d0)), DIMENSION(0:23, 2), INTENT(in)::PopProf_24hr
1117 
1118  REAL(KIND(1D0)), INTENT(in)::CO2PointSource
1119  REAL(KIND(1D0)), INTENT(in)::EF_umolCO2perJ
1120  REAL(KIND(1D0)), INTENT(in)::EnEF_v_Jkm
1121  REAL(KIND(1D0)), INTENT(in)::FrFossilFuel_Heat
1122  REAL(KIND(1D0)), INTENT(in)::FrFossilFuel_NonHeat
1123  REAL(KIND(1D0)), INTENT(in)::MaxFCMetab
1124  REAL(KIND(1D0)), INTENT(in)::MaxQFMetab
1125  REAL(KIND(1D0)), INTENT(in)::MinFCMetab
1126  REAL(KIND(1D0)), INTENT(in)::MinQFMetab
1127  REAL(KIND(1D0)), INTENT(in)::PopDensNighttime
1128  REAL(KIND(1D0)), INTENT(in)::QF_obs
1129  REAL(KIND(1D0)), INTENT(in)::Temp_C
1130  REAL(KIND(1D0)), INTENT(in)::TrafficUnits
1131 
1132  ! REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::sfr
1133  ! REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::SnowFrac
1134  REAL(KIND(1D0)), INTENT(IN)::SurfaceArea
1135 
1136  REAL(KIND(1D0)), INTENT(out)::Fc_anthro
1137  REAL(KIND(1D0)), INTENT(out)::Fc_build
1138  REAL(KIND(1D0)), INTENT(out)::Fc_metab
1139  REAL(KIND(1D0)), INTENT(out)::Fc_point
1140  REAL(KIND(1D0)), INTENT(out)::Fc_traff
1141  REAL(KIND(1D0)), INTENT(out)::QF
1142  REAL(KIND(1D0)), INTENT(out)::QF_SAHP
1143 
1144  INTEGER, PARAMETER :: notUsedI = -999
1145  REAL(KIND(1D0)), PARAMETER::notUsed = -999
1146 
1147  IF (emissionsmethod == 0) THEN ! use observed qf
1148  qf = qf_obs
1149  ELSEIF ((emissionsmethod > 0 .AND. emissionsmethod <= 6) .OR. emissionsmethod >= 11) THEN
1150  CALL anthropogenicemissions( &
1151  co2pointsource, emissionsmethod, &
1152  it, imin, dls, dayofweek_id, &
1153  ef_umolco2perj, fcef_v_kgkm, enef_v_jkm, trafficunits, &
1154  frfossilfuel_heat, frfossilfuel_nonheat, &
1155  minfcmetab, maxfcmetab, minqfmetab, maxqfmetab, &
1156  popdensdaytime, popdensnighttime, &
1157  temp_c, hdd_id, qf_a, qf_b, qf_c, &
1158  ah_min, ah_slope_heating, ah_slope_cooling, &
1159  baset_heating, baset_cooling, &
1160  trafficrate, &
1161  qf0_beu, qf_sahp, &
1162  fc_anthro, fc_metab, fc_traff, fc_build, fc_point, &
1163  ahprof_24hr, humactivity_24hr, traffprof_24hr, popprof_24hr, surfacearea)
1164 
1165  ELSE
1166  CALL errorhint(73, 'RunControl.nml:EmissionsMethod unusable', notused, notused, emissionsmethod)
1167  ENDIF
1168 
1169  IF (emissionsmethod >= 1) qf = qf_sahp
1170 
1171  IF (emissionsmethod >= 0 .AND. emissionsmethod <= 6) THEN
1172  fc_anthro = 0
1173  fc_metab = 0
1174  fc_traff = 0
1175  fc_build = 0
1176  fc_point = 0
1177  ENDIF
1178 
1179  END SUBROUTINE suews_cal_anthropogenicemission
1180  ! ================================================================================
1181 
1182  !==============BIOGENIC CO2 flux==================================================
1183  SUBROUTINE suews_cal_biogenco2( &
1184  alpha_bioCO2, alpha_enh_bioCO2, avkdn, avRh, beta_bioCO2, beta_enh_bioCO2, BSoilSurf, &! input:
1185  ConifSurf, DecidSurf, dectime, Diagnose, EmissionsMethod, Fc_anthro, G1, G2, G3, G4, &
1186  G5, G6, gfunc, GrassSurf, gsmodel, id, it, ivConif, ivDecid, ivGrass, Kmax, LAI_id, LAIMin, &
1187  LAIMax, MaxConductance, min_res_bioCO2, nsurf, NVegSurf, Press_hPa, resp_a, &
1188  resp_b, S1, S2, sfr, SMDMethod, SnowFrac, t2_C, Temp_C, theta_bioCO2, TH, TL, vsmd, xsmd, &
1189  Fc, Fc_biogen, Fc_photo, Fc_respi)! output:
1191  IMPLICIT NONE
1192 
1193  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(in)::alpha_bioCO2
1194  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(in)::alpha_enh_bioCO2
1195  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(in)::beta_bioCO2
1196  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(in)::beta_enh_bioCO2
1197  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(in)::LAI_id
1198  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(in)::LAIMin
1199  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(in)::LAIMax
1200  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(in)::min_res_bioCO2
1201  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(in)::resp_a
1202  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(in)::resp_b
1203  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(in)::theta_bioCO2
1204 
1205  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::sfr
1206  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::SnowFrac
1207 
1208  REAL(KIND(1d0)), DIMENSION(3), INTENT(in)::MaxConductance
1209 
1210  INTEGER, INTENT(in)::BSoilSurf
1211  INTEGER, INTENT(in)::ConifSurf
1212  INTEGER, INTENT(in)::DecidSurf
1213  INTEGER, INTENT(in)::Diagnose
1214  INTEGER, INTENT(in)::EmissionsMethod
1215  INTEGER, INTENT(in)::GrassSurf
1216  INTEGER, INTENT(in)::gsmodel
1217  INTEGER, INTENT(in)::id
1218  INTEGER, INTENT(in)::it
1219  INTEGER, INTENT(in)::ivConif
1220  INTEGER, INTENT(in)::ivDecid
1221  INTEGER, INTENT(in)::ivGrass
1222  INTEGER, INTENT(in)::nsurf
1223  INTEGER, INTENT(in)::NVegSurf
1224  INTEGER, INTENT(in)::SMDMethod
1225 
1226  REAL(KIND(1d0)), INTENT(in)::avkdn
1227  REAL(KIND(1d0)), INTENT(in)::avRh
1228  REAL(KIND(1d0)), INTENT(in)::dectime
1229  REAL(KIND(1d0)), INTENT(in)::Fc_anthro
1230  REAL(KIND(1d0)), INTENT(in)::G1
1231  REAL(KIND(1d0)), INTENT(in)::G2
1232  REAL(KIND(1d0)), INTENT(in)::G3
1233  REAL(KIND(1d0)), INTENT(in)::G4
1234  REAL(KIND(1d0)), INTENT(in)::G5
1235  REAL(KIND(1d0)), INTENT(in)::G6
1236  REAL(KIND(1d0)), INTENT(in)::gfunc
1237  REAL(KIND(1d0)), INTENT(in)::Kmax
1238  REAL(KIND(1d0)), INTENT(in)::Press_hPa
1239  REAL(KIND(1d0)), INTENT(in)::S1
1240  REAL(KIND(1d0)), INTENT(in)::S2
1241  REAL(KIND(1d0)), INTENT(in)::t2_C
1242  REAL(KIND(1d0)), INTENT(in)::Temp_C
1243  REAL(KIND(1d0)), INTENT(in)::TH
1244  REAL(KIND(1d0)), INTENT(in)::TL
1245  REAL(KIND(1d0)), INTENT(in)::vsmd
1246  REAL(KIND(1d0)), INTENT(in)::xsmd
1247 
1248  REAL(KIND(1d0)), INTENT(out)::Fc_biogen
1249  REAL(KIND(1d0)), INTENT(out)::Fc_photo
1250  REAL(KIND(1d0)), INTENT(out)::Fc_respi
1251  REAL(KIND(1d0)), INTENT(out)::Fc
1252 
1253  REAL(KIND(1d0))::gfunc2
1254  REAL(KIND(1d0))::dq
1255  REAL(KIND(1d0))::t2
1256  REAL(KIND(1d0))::dummy1
1257  REAL(KIND(1d0))::dummy2
1258  REAL(KIND(1d0))::dummy3
1259  REAL(KIND(1d0))::dummy4
1260  REAL(KIND(1d0))::dummy5
1261  REAL(KIND(1d0))::dummy6
1262  REAL(KIND(1d0))::dummy7
1263  REAL(KIND(1d0))::dummy8
1264  REAL(KIND(1d0))::dummy9
1265  REAL(KIND(1d0))::dummy10
1266  REAL(KIND(1d0))::dummy11
1267 
1268  IF (emissionsmethod >= 11) THEN
1269 
1270  IF (gsmodel == 3 .OR. gsmodel == 4) THEN ! With modelled 2 meter temperature
1271  ! Call LUMPS_cal_AtmMoist for dq and SurfaceResistance for gfunc with 2 meter temperature
1272  ! If modelled 2 meter temperature is too different from measured air temperature then
1273  ! use temp_c
1274  IF (abs(temp_c - t2_c) > 5) THEN
1275  t2 = temp_c
1276  ELSE
1277  t2 = t2_c
1278  ENDIF
1279 
1280  CALL cal_atmmoist( &
1281  t2, press_hpa, avrh, dectime, &! input:
1282  dummy1, dummy2, &! output:
1283  dummy3, dummy4, dummy5, dummy6, dq, dummy7, dummy8, dummy9)
1284 
1285  CALL surfaceresistance( &
1286  id, it, &! input:
1287  smdmethod, snowfrac, sfr, avkdn, t2, dq, xsmd, vsmd, maxconductance, &
1288  laimax, lai_id, gsmodel, kmax, &
1289  g1, g2, g3, g4, g5, g6, th, tl, s1, s2, &
1290  gfunc2, dummy10, dummy11)! output:
1291  ENDIF
1292 
1293  ! Calculate CO2 fluxes from biogenic components
1294  IF (diagnose == 1) WRITE (*, *) 'Calling CO2_biogen...'
1295  CALL co2_biogen( &
1296  alpha_bioco2, alpha_enh_bioco2, avkdn, beta_bioco2, beta_enh_bioco2, bsoilsurf, &! input:
1297  conifsurf, decidsurf, dectime, emissionsmethod, gfunc, gfunc2, grasssurf, gsmodel, &
1298  id, it, ivconif, ivdecid, ivgrass, lai_id, laimin, laimax, min_res_bioco2, nsurf, &
1299  nvegsurf, resp_a, resp_b, sfr, snowfrac, t2, temp_c, theta_bioco2, &
1300  fc_biogen, fc_photo, fc_respi)! output:
1301  ENDIF
1302 
1303  IF (emissionsmethod >= 0 .AND. emissionsmethod <= 6) THEN
1304  fc_biogen = 0
1305  fc_photo = 0
1306  fc_respi = 0
1307  ENDIF
1308 
1309  fc = fc_anthro + fc_biogen
1310 
1311  END SUBROUTINE suews_cal_biogenco2
1312  !========================================================================
1313 
1314  !=============net all-wave radiation=====================================
1315  SUBROUTINE suews_cal_qn( &
1316  NetRadiationMethod, snowUse, &!input
1317  tstep, SnowPack_prev, tau_a, tau_f, SnowAlbMax, SnowAlbMin, &
1318  Diagnose, snowFrac_obs, ldown_obs, fcld_obs, &
1319  dectime, ZENITH_deg, Tsurf_0, avKdn, Temp_C, avRH, ea_hPa, qn1_obs, &
1320  SnowAlb_prev, snowFrac_prev, DiagQN, &
1321  NARP_TRANS_SITE, NARP_EMIS_SNOW, IceFrac, sfr, emis, &
1322  alb_prev, albDecTr_id, albEveTr_id, albGrass_id, &!input
1323  alb_next, &!output
1324  ldown, fcld, &!output
1325  qn1, qn1_snowfree, qn1_S, kclear, kup, lup, tsurf, &
1326  qn1_ind_snow, kup_ind_snow, Tsurf_ind_snow, Tsurf_ind, &
1327  alb1, snowFrac_next, SnowAlb_next)
1329 
1330  IMPLICIT NONE
1331  ! INTEGER,PARAMETER ::nsurf = 7 ! number of surface types
1332  ! INTEGER,PARAMETER ::ConifSurf = 3 !New surface classes: Grass = 5th/7 surfaces
1333  ! INTEGER,PARAMETER ::DecidSurf = 4 !New surface classes: Grass = 5th/7 surfaces
1334  ! INTEGER,PARAMETER ::GrassSurf = 5
1335 
1336  INTEGER, INTENT(in)::NetRadiationMethod
1337  INTEGER, INTENT(in)::snowUse
1338  INTEGER, INTENT(in)::Diagnose
1339  INTEGER, INTENT(in)::DiagQN
1340  INTEGER, INTENT(in)::tstep
1341 
1342  REAL(KIND(1d0)), INTENT(in)::snowFrac_obs
1343  REAL(KIND(1d0)), INTENT(in)::ldown_obs
1344  REAL(KIND(1d0)), INTENT(in)::fcld_obs
1345  REAL(KIND(1d0)), INTENT(in)::dectime
1346  REAL(KIND(1d0)), INTENT(in)::ZENITH_deg
1347  REAL(KIND(1d0)), INTENT(in)::Tsurf_0
1348  REAL(KIND(1d0)), INTENT(in)::avKdn
1349  REAL(KIND(1d0)), INTENT(in)::Temp_C
1350  REAL(KIND(1d0)), INTENT(in)::avRH
1351  REAL(KIND(1d0)), INTENT(in)::ea_hPa
1352  REAL(KIND(1d0)), INTENT(in)::qn1_obs
1353  REAL(KIND(1d0)), INTENT(in)::SnowAlb_prev
1354  REAL(KIND(1d0)), INTENT(in)::NARP_EMIS_SNOW
1355  REAL(KIND(1d0)), INTENT(in)::NARP_TRANS_SITE
1356  REAL(KIND(1d0)), INTENT(in)::tau_a, tau_f, SnowAlbMax, SnowAlbMin
1357 
1358  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in):: IceFrac
1359  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in):: sfr
1360  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in):: emis
1361 
1362  REAL(KIND(1d0)), DIMENSION(nsurf) ::alb
1363  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in) ::alb_prev
1364  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out) ::alb_next
1365  REAL(KIND(1d0)), INTENT(in) ::albDecTr_id
1366  ! REAL(KIND(1d0)), INTENT(in) ::DecidCap_id
1367  REAL(KIND(1d0)), INTENT(in) ::albEveTr_id
1368  REAL(KIND(1d0)), INTENT(in) ::albGrass_id
1369 
1370  ! REAL(KIND(1d0)), DIMENSION(6, nsurf), INTENT(inout)::StoreDrainPrm
1371 
1372  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::SnowPack_prev
1373  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::snowFrac_prev
1374  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out)::snowFrac_next
1375  REAL(KIND(1d0)), DIMENSION(nsurf)::SnowFrac
1376 
1377  REAL(KIND(1d0)), INTENT(out)::ldown
1378  REAL(KIND(1d0)), INTENT(out)::fcld
1379  REAL(KIND(1d0)), INTENT(out)::qn1
1380  REAL(KIND(1d0)), INTENT(out)::qn1_snowfree
1381  REAL(KIND(1d0)), INTENT(out)::qn1_S
1382  REAL(KIND(1d0)), INTENT(out)::kclear
1383  REAL(KIND(1d0)), INTENT(out)::kup
1384  REAL(KIND(1d0)), INTENT(out)::lup
1385  REAL(KIND(1d0)), INTENT(out)::tsurf
1386  REAL(KIND(1d0)), INTENT(out)::alb1
1387  REAL(KIND(1d0)), INTENT(out)::SnowAlb_next
1388  REAL(KIND(1d0))::alb0
1389  REAL(KIND(1d0))::SnowAlb
1390 
1391  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out) ::qn1_ind_snow
1392  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out) ::kup_ind_snow
1393  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out) ::Tsurf_ind_snow
1394  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out):: tsurf_ind
1395 
1396  REAL(KIND(1d0)), DIMENSION(nsurf):: lup_ind
1397  REAL(KIND(1d0)), DIMENSION(nsurf):: kup_ind
1398  REAL(KIND(1d0)), DIMENSION(nsurf):: qn1_ind
1399 
1400  REAL(KIND(1d0)), PARAMETER::NAN = -999
1401  INTEGER :: NetRadiationMethod_use
1402  INTEGER::AlbedoChoice, ldown_option
1403 
1404  ! translate values
1405  alb = alb_prev
1406 
1407  ! update snow albedo
1408  snowalb = update_snow_albedo( &
1409  tstep, snowpack_prev, snowalb_prev, temp_c, &
1410  tau_a, tau_f, snowalbmax, snowalbmin)
1411 
1412  CALL radmethod( &
1413  netradiationmethod, &!input
1414  snowuse, &!input
1415  netradiationmethod_use, albedochoice, ldown_option)!output
1416 
1417  snowfrac = snowfrac_prev
1418  IF (netradiationmethod_use > 0) THEN
1419 
1420  ! IF (snowUse==0) SnowFrac=snowFrac_obs
1421  IF (snowuse == 0) snowfrac = 0
1422 
1423  IF (ldown_option == 1) THEN !Observed ldown provided as forcing
1424  ldown = ldown_obs
1425  ELSE
1426  ldown = -9 !to be filled in NARP
1427  ENDIF
1428 
1429  IF (ldown_option == 2) THEN !observed cloud fraction provided as forcing
1430  fcld = fcld_obs
1431  ENDIF
1432 
1433  !write(*,*) DecidCap(id), id, it, imin, 'Calc - near start'
1434 
1435  ! Update variables that change daily and represent seasonal variability
1436  alb(decidsurf) = albdectr_id !Change deciduous albedo
1437  ! StoreDrainPrm(6, DecidSurf) = DecidCap_id !Change current storage capacity of deciduous trees
1438  ! Change EveTr and Grass albedo too
1439  alb(conifsurf) = albevetr_id
1440  alb(grasssurf) = albgrass_id
1441 
1442  IF (diagnose == 1) WRITE (*, *) 'Calling NARP...'
1443  IF (diagqn == 1) WRITE (*, *) 'NetRadiationMethodX:', netradiationmethod_use
1444  IF (diagqn == 1) WRITE (*, *) 'AlbedoChoice:', albedochoice
1445 
1446  CALL narp( &
1447  nsurf, sfr, snowfrac, alb, emis, icefrac, &! input:
1448  narp_trans_site, narp_emis_snow, &
1449  dectime, zenith_deg, tsurf_0, avkdn, temp_c, avrh, ea_hpa, qn1_obs, &
1450  snowalb, &
1451  albedochoice, ldown_option, netradiationmethod_use, diagqn, &
1452  qn1, qn1_snowfree, qn1_s, kclear, kup, ldown, lup, fcld, tsurf, &! output:
1453  qn1_ind_snow, kup_ind_snow, tsurf_ind_snow, tsurf_ind, alb0, alb1)
1454 
1455  ELSE ! NetRadiationMethod==0
1456  snowfrac = snowfrac_obs
1457  qn1 = qn1_obs
1458  qn1_snowfree = qn1_obs
1459  qn1_s = qn1_obs
1460  ldown = nan
1461  lup = nan
1462  kup = nan
1463  tsurf = nan
1464  lup_ind = nan
1465  kup_ind = nan
1466  tsurf_ind = nan
1467  qn1_ind = nan
1468  fcld = nan
1469  ENDIF
1470  snowfrac_next = snowfrac
1471 
1472  IF (ldown_option == 1) THEN
1473  fcld = nan
1474  ENDIF
1475 
1476  ! translate values
1477  alb_next = alb
1478  snowalb_next = snowalb
1479 
1480  END SUBROUTINE suews_cal_qn
1481  !========================================================================
1482 
1483  !=============storage heat flux=========================================
1484  SUBROUTINE suews_cal_qs( &
1485  StorageHeatMethod, qs_obs, OHMIncQF, Gridiv, &!input
1486  id, tstep, dt_since_start, Diagnose, sfr, &
1487  OHM_coef, OHM_threshSW, OHM_threshWD, &
1488  soilstore_id, SoilStoreCap, state_id, SnowUse, SnowFrac, DiagQS, &
1489  HDD_id, MetForcingData_grid, Ts5mindata_ir, qf, qn1, &
1490  avkdn, avu1, temp_c, zenith_deg, avrh, press_hpa, ldown, &
1491  bldgh, alb, emis, cpAnOHM, kkAnOHM, chAnOHM, EmissionsMethod, &
1492  Tair_av, qn1_av_prev, dqndt_prev, qn1_s_av_prev, dqnsdt_prev, &
1493  StoreDrainPrm, &
1494  qn1_S, dataOutLineESTM, qs, &!output
1495  qn1_av_next, dqndt_next, qn1_s_av_next, dqnsdt_next, &
1496  deltaQi, a1, a2, a3)
1498  IMPLICIT NONE
1499 
1500  INTEGER, INTENT(in) ::StorageHeatMethod
1501  INTEGER, INTENT(in) ::OHMIncQF
1502  INTEGER, INTENT(in) ::Gridiv
1503  INTEGER, INTENT(in) ::id
1504  INTEGER, INTENT(in) ::tstep ! time step [s]
1505  INTEGER, INTENT(in) ::dt_since_start ! time since simulation starts [s]
1506  INTEGER, INTENT(in) ::Diagnose
1507  ! INTEGER, INTENT(in) ::nsh ! number of timesteps in one hour
1508  INTEGER, INTENT(in) ::SnowUse ! option for snow related calculations
1509  INTEGER, INTENT(in) ::DiagQS ! diagnostic option
1510  INTEGER, INTENT(in) :: EmissionsMethod
1511 
1512  REAL(KIND(1d0)), INTENT(in)::OHM_coef(nsurf + 1, 4, 3) ! OHM coefficients
1513  REAL(KIND(1d0)), INTENT(in)::OHM_threshSW(nsurf + 1) ! OHM thresholds
1514  REAL(KIND(1d0)), INTENT(in)::OHM_threshWD(nsurf + 1) ! OHM thresholds
1515  REAL(KIND(1d0)), INTENT(in)::soilstore_id(nsurf) ! soil moisture
1516  REAL(KIND(1d0)), INTENT(in)::SoilStoreCap(nsurf) ! capacity of soil store
1517  REAL(KIND(1d0)), INTENT(in)::state_id(nsurf) ! wetness status
1518 
1519  REAL(KIND(1d0)), DIMENSION(12), INTENT(in)::HDD_id
1520  REAL(KIND(1d0)), INTENT(in)::qf
1521  REAL(KIND(1d0)), INTENT(in)::qn1
1522  REAL(KIND(1d0)), INTENT(in)::qs_obs
1523  REAL(KIND(1d0)), INTENT(in)::avkdn, avu1, temp_c, zenith_deg, avrh, press_hpa, ldown
1524  REAL(KIND(1d0)), INTENT(in)::bldgh
1525 
1526  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::sfr
1527  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::alb
1528  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::emis
1529  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::cpAnOHM
1530  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::kkAnOHM
1531  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::chAnOHM
1532  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::SnowFrac
1533 
1534  REAL(KIND(1d0)), DIMENSION(:, :), INTENT(in)::MetForcingData_grid
1535 
1536  REAL(KIND(1d0)), DIMENSION(:), INTENT(in)::Ts5mindata_ir
1537 
1538  REAL(KIND(1d0)), INTENT(in)::Tair_av
1539  REAL(KIND(1d0)), INTENT(in)::qn1_av_prev
1540  REAL(KIND(1d0)), INTENT(out)::qn1_av_next
1541  REAL(KIND(1d0)), INTENT(in)::dqndt_prev ! Rate of change of net radiation [W m-2 h-1] at t-1
1542  REAL(KIND(1d0)), INTENT(out)::dqndt_next ! Rate of change of net radiation [W m-2 h-1] at t-1
1543  REAL(KIND(1d0)), INTENT(in)::qn1_s_av_prev
1544  REAL(KIND(1d0)), INTENT(out)::qn1_s_av_next
1545  REAL(KIND(1d0)), INTENT(in)::dqnsdt_prev ! Rate of change of net radiation [W m-2 h-1] at t-1
1546  REAL(KIND(1d0)), INTENT(out)::dqnsdt_next ! Rate of change of net radiation [W m-2 h-1] at t-1
1547  ! REAL(KIND(1d0)),DIMENSION(nsh),INTENT(inout) ::qn1_store_grid
1548  ! REAL(KIND(1d0)),DIMENSION(nsh),INTENT(inout) ::qn1_S_store_grid !< stored qn1 [W m-2]
1549 
1550  ! REAL(KIND(1d0)),DIMENSION(2*nsh+1),INTENT(inout)::qn1_av_store_grid
1551  ! REAL(KIND(1d0)),DIMENSION(2*nsh+1),INTENT(inout)::qn1_S_av_store_grid !< average net radiation over previous hour [W m-2]
1552  REAL(KIND(1d0)), DIMENSION(6, nsurf), INTENT(in)::StoreDrainPrm
1553 
1554  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out)::deltaQi ! storage heat flux of snow surfaces
1555 
1556  REAL(KIND(1d0)), DIMENSION(27), INTENT(out):: dataOutLineESTM
1557  REAL(KIND(1d0)), INTENT(out)::qn1_S
1558  REAL(KIND(1d0)), INTENT(out):: qs ! storage heat flux
1559  REAL(KIND(1d0)), INTENT(out):: a1
1560  REAL(KIND(1d0)), INTENT(out):: a2
1561  REAL(KIND(1d0)), INTENT(out):: a3
1562 
1563  REAL(KIND(1d0))::Tair_mav_5d ! Tair_mav_5d=HDD(id-1,4) HDD at the begining of today (id-1)
1564  REAL(KIND(1d0))::qn1_use ! qn used in OHM calculations
1565 
1566  REAL(KIND(1d0)):: moist_surf(nsurf)
1567 
1568  ! initialise output variables
1569  !deltaQi = 0
1570  !SnowFrac = 0
1571  !qn1_S = 0
1572  dataoutlineestm = -999
1573  qs = -999
1574  a1 = -999
1575  a2 = -999
1576  a3 = -999
1577 
1578  ! calculate qn if qf should be included
1579  IF (ohmincqf == 1) THEN
1580  qn1_use = qf + qn1
1581  ELSEIF (ohmincqf == 0) THEN
1582  qn1_use = qn1
1583  ENDIF
1584 
1585  IF (storageheatmethod == 0) THEN !Use observed QS
1586  qs = qs_obs
1587 
1588  ELSEIF (storageheatmethod == 1) THEN !Use OHM to calculate QS
1589  tair_mav_5d = hdd_id(10)
1590  IF (diagnose == 1) WRITE (*, *) 'Calling OHM...'
1591  CALL ohm(qn1_use, qn1_av_prev, dqndt_prev, qn1_av_next, dqndt_next, &
1592  qn1_s, qn1_s_av_prev, dqnsdt_prev, qn1_s_av_next, dqnsdt_next, &
1593  tstep, dt_since_start, &
1594  sfr, nsurf, &
1595  tair_mav_5d, &
1596  ohm_coef, &
1597  ohm_threshsw, ohm_threshwd, &
1598  soilstore_id, soilstorecap, state_id, &
1599  bldgsurf, watersurf, &
1600  snowuse, snowfrac, &
1601  diagqs, &
1602  a1, a2, a3, qs, deltaqi)
1603 
1604  ! use AnOHM to calculate QS, TS 14 Mar 2016
1605  ELSEIF (storageheatmethod == 3) THEN
1606  IF (diagnose == 1) WRITE (*, *) 'Calling AnOHM...'
1607  ! CALL AnOHM(qn1_use,qn1_store_grid,qn1_av_store_grid,qf,&
1608  ! MetForcingData_grid,state_id/StoreDrainPrm(6,:),&
1609  ! alb, emis, cpAnOHM, kkAnOHM, chAnOHM,&
1610  ! sfr,nsurf,nsh,EmissionsMethod,id,Gridiv,&
1611  ! a1,a2,a3,qs,deltaQi)
1612  moist_surf = state_id/storedrainprm(6, :)
1613  CALL anohm( &
1614  tstep, dt_since_start, &
1615  qn1_use, qn1_av_prev, dqndt_prev, qf, &
1616  metforcingdata_grid, moist_surf, &
1617  alb, emis, cpanohm, kkanohm, chanohm, &! input
1618  sfr, nsurf, emissionsmethod, id, gridiv, &
1619  qn1_av_next, dqndt_next, &
1620  a1, a2, a3, qs, deltaqi)! output
1621 
1622  ! !Calculate QS using ESTM
1623  ELSEIF (storageheatmethod == 4 .OR. storageheatmethod == 14) THEN
1624  ! !CALL ESTM(QSestm,iMB)
1625  IF (diagnose == 1) WRITE (*, *) 'Calling ESTM...'
1626  CALL estm( &
1627  gridiv, &!input
1628  tstep, &
1629  avkdn, avu1, temp_c, zenith_deg, avrh, press_hpa, ldown, &
1630  bldgh, ts5mindata_ir, &
1631  tair_av, &
1632  dataoutlineestm, qs)!output
1633  ! CALL ESTM(QSestm,Gridiv,ir) ! iMB corrected to Gridiv, TS 09 Jun 2016
1634  ! QS=QSestm ! Use ESTM qs
1635  ENDIF
1636 
1637  END SUBROUTINE suews_cal_qs
1638  !=======================================================================
1639 
1640  !==========================drainage and runoff================================
1641  SUBROUTINE suews_cal_water( &
1642  Diagnose, &!input
1643  snowUse, NonWaterFraction, addPipes, addImpervious, addVeg, addWaterBody, &
1644  state_id, soilstore_id, sfr, StoreDrainPrm, WaterDist, nsh_real, &
1645  drain_per_tstep, & !output
1646  drain, frac_water2runoff, &
1647  AdditionalWater, runoffPipes, runoff_per_interval, &
1648  AddWater, stateOld, soilstoreOld)
1650  IMPLICIT NONE
1651  ! INTEGER,PARAMETER :: nsurf=7! number of surface types
1652  ! INTEGER,PARAMETER ::WaterSurf = 7
1653  INTEGER, INTENT(in) ::Diagnose
1654  INTEGER, INTENT(in) ::snowUse
1655 
1656  REAL(KIND(1d0)), INTENT(in)::NonWaterFraction
1657  REAL(KIND(1d0)), INTENT(in)::addPipes
1658  REAL(KIND(1d0)), INTENT(in)::addImpervious
1659  REAL(KIND(1d0)), INTENT(in)::addVeg
1660  REAL(KIND(1d0)), INTENT(in)::addWaterBody
1661  REAL(KIND(1d0)), INTENT(in)::nsh_real !nsh cast as a real for use in calculations
1662 
1663  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in) ::state_id
1664  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in) ::soilstore_id
1665  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in) ::sfr
1666  REAL(KIND(1d0)), DIMENSION(6, nsurf), INTENT(in) ::StoreDrainPrm
1667  REAL(KIND(1d0)), DIMENSION(nsurf + 1, nsurf - 1), INTENT(in)::WaterDist
1668 
1669  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out):: drain !Drainage of surface type "is" [mm]
1670  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out):: frac_water2runoff!Fraction of water going to runoff/sub-surface soil (WGWaterDist) [-]
1671  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out):: AddWater !water from other surfaces (WGWaterDist in SUEWS_ReDistributeWater.f95) [mm]
1672  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out):: stateOld
1673  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out):: soilstoreOld
1674 
1675  REAL(KIND(1d0)), INTENT(out)::drain_per_tstep
1676  REAL(KIND(1d0)), INTENT(out)::AdditionalWater
1677  REAL(KIND(1d0)), INTENT(out)::runoffPipes
1678  REAL(KIND(1d0)), INTENT(out)::runoff_per_interval
1679  INTEGER:: is
1680 
1681  ! Retain previous surface state_id and soil moisture state_id
1682  stateold = state_id !state_id of each surface [mm] for the previous timestep
1683  soilstoreold = soilstore_id !Soil moisture of each surface [mm] for the previous timestep
1684 
1685  !============= Grid-to-grid runoff =============
1686  ! Calculate additional water coming from other grids
1687  ! i.e. the variables addImpervious, addVeg, addWaterBody, addPipes
1688  !call RunoffFromGrid(GridFromFrac) !!Need to code between-grid water transfer
1689 
1690  ! Sum water coming from other grids (these are expressed as depths over the whole surface)
1691  additionalwater = addpipes + addimpervious + addveg + addwaterbody ![mm]
1692 
1693  ! Initialise runoff in pipes
1694  runoffpipes = addpipes !Water flowing in pipes from other grids. QUESTION: No need for scaling?
1695  !! CHECK p_i
1696  runoff_per_interval = addpipes !pipe plor added to total runoff.
1697 
1698  !================== Drainage ===================
1699  ! Calculate drainage for each soil subsurface (excluding water body)
1700  IF (diagnose == 1) WRITE (*, *) 'Calling Drainage...'
1701 
1702  IF (nonwaterfraction /= 0) THEN !Soil states only calculated if soil exists. LJ June 2017
1703  DO is = 1, nsurf - 1
1704 
1705  CALL drainage( &
1706  is, &! input:
1707  state_id(is), &
1708  storedrainprm(6, is), &
1709  storedrainprm(2, is), &
1710  storedrainprm(3, is), &
1711  storedrainprm(4, is), &
1712  nsh_real, &
1713  drain(is))! output
1714 
1715  ! !HCW added and changed to StoreDrainPrm(6,is) here 20 Feb 2015
1716  ! drain_per_tstep=drain_per_tstep+(drain(is)*sfr(is)/NonWaterFraction) !No water body included
1717  ENDDO
1718  drain_per_tstep = dot_product(drain(1:nsurf - 1), sfr(1:nsurf - 1))/nonwaterfraction !No water body included
1719  ELSE
1720  drain(1:nsurf - 1) = 0
1721  drain_per_tstep = 0
1722  ENDIF
1723 
1724  drain(watersurf) = 0 ! Set drainage from water body to zero
1725 
1726  ! Distribute water within grid, according to WithinGridWaterDist matrix (Cols 1-7)
1727  IF (diagnose == 1) WRITE (*, *) 'Calling ReDistributeWater...'
1728  ! CALL ReDistributeWater
1729  !Calculates AddWater(is)
1730  CALL redistributewater( &
1731  snowuse, waterdist, sfr, drain, &! input:
1732  frac_water2runoff, addwater)! output
1733 
1734  END SUBROUTINE suews_cal_water
1735  !=======================================================================
1736 
1737  !===============initialize sensible heat flux============================
1738  SUBROUTINE suews_init_qh( &
1739  avdens, avcp, h_mod, qn1, dectime, &!input
1740  H_init)!output
1742  IMPLICIT NONE
1743  ! REAL(KIND(1d0)), INTENT(in)::qh_obs
1744  REAL(KIND(1d0)), INTENT(in)::avdens
1745  REAL(KIND(1d0)), INTENT(in)::avcp
1746  REAL(KIND(1d0)), INTENT(in)::h_mod
1747  REAL(KIND(1d0)), INTENT(in)::qn1
1748  REAL(KIND(1d0)), INTENT(in)::dectime
1749  REAL(KIND(1d0)), INTENT(out)::H_init
1750 
1751  REAL(KIND(1d0)), PARAMETER::NAN = -999
1752  INTEGER, PARAMETER::notUsedI = -999
1753 
1754  ! Calculate kinematic heat flux (w'T') from sensible heat flux [W m-2] from observed data (if available) or LUMPS
1755  ! IF (qh_obs /= NAN) THEN !if(qh_obs/=NAN) qh=qh_obs !Commented out by HCW 04 Mar 2015
1756  ! H_init = qh_obs/(avdens*avcp) !Use observed value
1757  ! ELSE
1758  IF (h_mod /= nan) THEN
1759  h_init = h_mod/(avdens*avcp) !Use LUMPS value
1760  ELSE
1761  h_init = (qn1*0.2)/(avdens*avcp) !If LUMPS has had a problem, we still need a value
1762  CALL errorhint(38, 'LUMPS unable to calculate realistic value for H_mod.', h_mod, dectime, notusedi)
1763  ENDIF
1764  ! ENDIF
1765 
1766  END SUBROUTINE suews_init_qh
1767  !========================================================================
1768 
1769  !================latent heat flux and surface wetness===================
1770  ! TODO: optimise the structure of this function
1771  SUBROUTINE suews_cal_qe( &
1772  Diagnose, snowuse, &!input
1773  tstep, imin, it, EvapMethod, snowCalcSwitch, dayofWeek_id, CRWmin, CRWmax, &
1774  dectime, avdens, avcp, lv_J_kg, lvS_J_kg, avRh, Press_hPa, Temp_C, &
1775  RAsnow, psyc_hPa, sIce_hPa, &
1776  PervFraction, vegfraction, addimpervious, qn1_snowfree, qf, qs, vpd_hPa, s_hPa, &
1777  ResistSurf, RA, rb, snowdensmin, precip, PipeCapacity, RunoffToWater, &
1778  NonWaterFraction, WU_nsurf, addVeg, addWaterBody, SnowLimPaved, SnowLimBldg, &
1779  SurfaceArea, FlowChange, drain, WetThresh, stateOld, mw_ind, SoilStoreCap, rainonsnow, &
1780  freezmelt, freezstate, freezstatevol, Qm_Melt, Qm_rain, Tsurf_ind, sfr, &
1781  StateLimit, AddWater, addwaterrunoff, StoreDrainPrm, SnowPackLimit, SnowProf_24hr, &
1782  SnowPack_in, SnowFrac_in, SnowWater_in, iceFrac_in, SnowDens_in, &! input:
1783  runoff_per_interval_in, state_id_in, soilstore_id_in, &! input:
1784  state_id_out, soilstore_id_out, &! output:
1785  SnowPack_out, SnowFrac_out, SnowWater_out, iceFrac_out, SnowDens_out, &! output
1786  runoffSoil, &! output:
1787  SnowRemoval, &
1788  state_per_tstep, NWstate_per_tstep, qe, &
1789  swe, chSnow_per_interval, ev_per_tstep, runoff_per_tstep, &
1790  surf_chang_per_tstep, runoffPipes, mwstore, runoffwaterbody, &
1791  runoffAGveg, runoffAGimpervious)
1793  IMPLICIT NONE
1794 
1795  INTEGER, INTENT(in) ::Diagnose
1796  INTEGER, INTENT(in) ::snowuse
1797  INTEGER, INTENT(in) ::tstep
1798  INTEGER, INTENT(in) ::imin
1799  INTEGER, INTENT(in) ::it
1800  INTEGER, INTENT(in) ::EvapMethod !Evaporation calculated according to Rutter (1) or Shuttleworth (2)
1801 
1802  INTEGER, DIMENSION(nsurf), INTENT(in)::snowCalcSwitch
1803  INTEGER, DIMENSION(3), INTENT(in)::dayofWeek_id
1804 
1805  REAL(KIND(1d0)), INTENT(in)::CRWmin
1806  REAL(KIND(1d0)), INTENT(in)::CRWmax
1807  REAL(KIND(1d0)), INTENT(in)::dectime
1808  REAL(KIND(1d0)), INTENT(in)::lvS_J_kg
1809  REAL(KIND(1d0)), INTENT(in)::lv_j_kg
1810  REAL(KIND(1d0)), INTENT(in)::avdens
1811  REAL(KIND(1d0)), INTENT(in)::avRh
1812  REAL(KIND(1d0)), INTENT(in)::Press_hPa
1813  REAL(KIND(1d0)), INTENT(in)::Temp_C
1814  REAL(KIND(1d0)), INTENT(in)::RAsnow
1815  REAL(KIND(1d0)), INTENT(in)::psyc_hPa
1816  REAL(KIND(1d0)), INTENT(in)::avcp
1817  REAL(KIND(1d0)), INTENT(in)::sIce_hPa
1818  REAL(KIND(1d0)), INTENT(in)::PervFraction
1819  REAL(KIND(1d0)), INTENT(in)::vegfraction
1820  REAL(KIND(1d0)), INTENT(in)::addimpervious
1821  REAL(KIND(1d0)), INTENT(in)::qn1_snowfree
1822  REAL(KIND(1d0)), INTENT(in)::qf
1823  REAL(KIND(1d0)), INTENT(in)::qs
1824  REAL(KIND(1d0)), INTENT(in)::vpd_hPa
1825  REAL(KIND(1d0)), INTENT(in)::s_hPa
1826  REAL(KIND(1d0)), INTENT(in)::ResistSurf
1827  REAL(KIND(1d0)), INTENT(in)::RA
1828  REAL(KIND(1d0)), INTENT(in)::rb
1829  REAL(KIND(1d0)), INTENT(in)::snowdensmin
1830  REAL(KIND(1d0)), INTENT(in)::precip
1831  REAL(KIND(1d0)), INTENT(in)::PipeCapacity
1832  REAL(KIND(1d0)), INTENT(in)::RunoffToWater
1833  REAL(KIND(1d0)), INTENT(in)::NonWaterFraction
1834  ! REAL(KIND(1d0)), INTENT(in)::wu_EveTr!Water use for evergreen trees/shrubs [mm]
1835  ! REAL(KIND(1d0)), INTENT(in)::wu_DecTr!Water use for deciduous trees/shrubs [mm]
1836  ! REAL(KIND(1d0)), INTENT(in)::wu_Grass!Water use for grass [mm]
1837  REAL(KIND(1d0)), INTENT(in)::addVeg!Water from vegetated surfaces of other grids [mm] for whole surface area
1838  REAL(KIND(1d0)), INTENT(in)::addWaterBody!Water from water surface of other grids [mm] for whole surface area
1839  REAL(KIND(1d0)), INTENT(in)::SnowLimPaved
1840  REAL(KIND(1d0)), INTENT(in)::SnowLimBldg
1841  REAL(KIND(1d0)), INTENT(in)::SurfaceArea
1842  REAL(KIND(1d0)), INTENT(in)::FlowChange
1843 
1844  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::WU_nsurf
1845  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::drain
1846  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::WetThresh
1847  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::stateOld
1848  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::mw_ind
1849  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::SoilStoreCap
1850  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::rainonsnow
1851  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::freezmelt
1852  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::freezstate
1853  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::freezstatevol
1854  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::Qm_Melt
1855  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::Qm_rain
1856  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::Tsurf_ind
1857  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::sfr
1858  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::SnowPackLimit
1859  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::StateLimit !Limit for state_id of each surface type [mm] (specified in input files)
1860  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::AddWater
1861  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::addwaterrunoff
1862  REAL(KIND(1d0)), DIMENSION(6, nsurf), INTENT(in)::StoreDrainPrm
1863  REAL(KIND(1d0)), DIMENSION(0:23, 2), INTENT(in):: SnowProf_24hr
1864 
1865  ! Total water transported to each grid for grid-to-grid connectivity
1866  REAL(KIND(1d0)), INTENT(in)::runoff_per_interval_in
1867  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::state_id_in
1868  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::soilstore_id_in
1869  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::SnowPack_in
1870  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::SnowFrac_in
1871  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::SnowWater_in
1872  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::iceFrac_in
1873  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::SnowDens_in
1874 
1875  ! output:
1876  REAL(KIND(1d0))::runoff_per_interval_out
1877  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out)::state_id_out
1878  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out)::soilstore_id_out
1879  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out)::SnowPack_out
1880  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out)::SnowFrac_out
1881  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out)::SnowWater_out
1882  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out)::iceFrac_out
1883  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out)::SnowDens_out
1884 
1885  REAL(KIND(1d0)), DIMENSION(nsurf)::runoffSnow !Initialize for runoff caused by snowmelting
1886  REAL(KIND(1d0)), DIMENSION(nsurf)::runoff
1887  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(out)::runoffSoil
1888  REAL(KIND(1d0)), DIMENSION(nsurf)::chang
1889  REAL(KIND(1d0)), DIMENSION(nsurf)::changSnow
1890  REAL(KIND(1d0)), DIMENSION(nsurf)::snowDepth
1891  REAL(KIND(1d0)), DIMENSION(nsurf)::SnowToSurf
1892  REAL(KIND(1d0)), DIMENSION(nsurf)::ev_snow
1893  REAL(KIND(1d0)), DIMENSION(2), INTENT(out)::SnowRemoval
1894  REAL(KIND(1d0)), DIMENSION(nsurf)::evap
1895  REAL(KIND(1d0)), DIMENSION(nsurf)::rss_nsurf
1896 
1897  REAL(KIND(1d0))::p_mm!Inputs to surface water balance
1898  ! REAL(KIND(1d0)),INTENT(out)::rss
1899  REAL(KIND(1d0))::qe_surf ! latent heat flux [W m-2]
1900  REAL(KIND(1d0)), INTENT(out)::state_per_tstep
1901  REAL(KIND(1d0)), INTENT(out)::NWstate_per_tstep
1902  REAL(KIND(1d0)), INTENT(out)::qe
1903  REAL(KIND(1d0)), INTENT(out)::swe
1904  REAL(KIND(1d0))::ev
1905  REAL(KIND(1d0)), INTENT(out)::chSnow_per_interval
1906  REAL(KIND(1d0)), INTENT(out)::ev_per_tstep
1907  REAL(KIND(1d0))::qe_per_tstep
1908  REAL(KIND(1d0)), INTENT(out)::runoff_per_tstep
1909  REAL(KIND(1d0)), INTENT(out)::surf_chang_per_tstep
1910  REAL(KIND(1d0)), INTENT(out)::runoffPipes
1911  REAL(KIND(1d0)), INTENT(out)::mwstore
1912  REAL(KIND(1d0)), INTENT(out)::runoffwaterbody
1913  REAL(KIND(1d0))::runoffWaterBody_m3
1914  REAL(KIND(1d0))::runoffPipes_m3
1915  REAL(KIND(1d0)), INTENT(out)::runoffAGveg
1916  REAL(KIND(1d0)), INTENT(out)::runoffAGimpervious
1917 
1918  ! local:
1919  INTEGER:: is
1920 
1921  REAL(KIND(1d0))::runoff_per_interval
1922  REAL(KIND(1d0)), DIMENSION(nsurf)::state_id
1923  REAL(KIND(1d0)), DIMENSION(nsurf)::soilstore_id
1924  REAL(KIND(1d0)), DIMENSION(nsurf)::SnowPack
1925  REAL(KIND(1d0)), DIMENSION(nsurf)::SnowFrac
1926  REAL(KIND(1d0)), DIMENSION(nsurf)::SnowWater
1927  REAL(KIND(1d0)), DIMENSION(nsurf)::iceFrac
1928  REAL(KIND(1d0)), DIMENSION(nsurf)::SnowDens
1929 
1930  REAL(KIND(1d0)), DIMENSION(2) ::SurplusEvap !Surplus for evaporation in 5 min timestep
1931  REAL(KIND(1d0))::surplusWaterBody
1932  REAL(KIND(1d0))::pin!Rain per time interval
1933  ! REAL(KIND(1d0))::sae
1934  ! REAL(KIND(1d0))::vdrc
1935  ! REAL(KIND(1d0))::sp
1936  ! REAL(KIND(1d0))::numPM
1937  REAL(KIND(1d0))::qn_e
1938  REAL(KIND(1d0))::tlv
1939  REAL(KIND(1d0))::runoffAGimpervious_m3
1940  REAL(KIND(1d0))::runoffAGveg_m3
1941  REAL(KIND(1d0))::nsh_real
1942  REAL(KIND(1d0))::tstep_real
1943  REAL(KIND(1d0))::ev_tot
1944  REAL(KIND(1d0))::qe_tot
1945  REAL(KIND(1d0))::surf_chang_tot
1946  REAL(KIND(1d0))::runoff_tot
1947  REAL(KIND(1d0))::chSnow_tot
1948 
1949  REAL(KIND(1d0)), DIMENSION(7)::capStore ! current storage capacity [mm]
1950 
1951  runoff_per_interval = runoff_per_interval_in
1952  state_id = state_id_in
1953  soilstore_id = soilstore_id_in
1954  snowpack = snowpack_in
1955  snowfrac = snowfrac_in
1956  snowwater = snowwater_in
1957  icefrac = icefrac_in
1958  snowdens = snowdens_in
1959 
1960  tstep_real = tstep*1.d0
1961  nsh_real = 3600/tstep_real
1962 
1963  capstore = 0 !initialise capStore
1964 
1965  tlv = lv_j_kg/tstep_real !Latent heat of vapourisation per timestep
1966 
1967  pin = max(0., precip)!Initiate rain data [mm]
1968 
1969  ! Initialize the output variables
1970  qe_surf = 0
1971  ev = 0
1972  qe_tot = 0
1973  ev_tot = 0
1974  swe = 0
1975  ev_snow = 0
1976  ev_per_tstep = 0
1977  qe_per_tstep = 0
1978  surf_chang_per_tstep = 0
1979  runoff_per_tstep = 0
1980  state_per_tstep = 0
1981  nwstate_per_tstep = 0
1982  qe = 0
1983  runoffwaterbody = 0
1984  chsnow_per_interval = 0
1985  mwstore = 0
1986  runoffagveg = 0
1987  runoffagimpervious = 0
1988  surpluswaterbody = 0
1989  runoffsoil = 0
1990  runoff = 0
1991  chang = 0
1992  surplusevap = 0
1993  snowremoval = 0
1994 
1995  ! net available energy for evaporation
1996  qn_e = qn1_snowfree + qf - qs ! qn1 changed to qn1_snowfree, lj in May 2013
1997 
1998  IF (diagnose == 1) WRITE (*, *) 'Calling evap_SUEWS and SoilStore...'
1999  DO is = 1, nsurf !For each surface in turn
2000  IF (snowuse == 1 .and. snowcalcswitch(is) == 1) THEN ! snow calculation
2001  IF (sfr(is) /= 0) THEN
2002  ! IF (Diagnose == 1) WRITE (*, *) 'Calling SnowCalc...'
2003  CALL snowcalc( &
2004  tstep, imin, it, dectime, is, &!input
2005  evapmethod, crwmin, crwmax, nsh_real, lvs_j_kg, avdens, &
2006  avrh, press_hpa, temp_c, rasnow, psyc_hpa, avcp, sice_hpa, &
2007  pervfraction, vegfraction, addimpervious, &
2008  vpd_hpa, qn_e, s_hpa, resistsurf, ra, rb, tlv, snowdensmin, snowprof_24hr, precip, &
2009  pipecapacity, runofftowater, &
2010  addveg, snowlimpaved, snowlimbldg, flowchange, drain, &
2011  wetthresh, stateold, mw_ind, soilstorecap, rainonsnow, &
2012  freezmelt, freezstate, freezstatevol, &
2013  qm_melt, qm_rain, tsurf_ind, sfr, dayofweek_id, storedrainprm, snowpacklimit, &
2014  addwater, addwaterrunoff, &
2015  soilstore_id, snowpack, surplusevap, &!inout
2016  snowfrac, snowwater, icefrac, snowdens, &
2017  runoffagimpervious, runoffagveg, surpluswaterbody, &
2018  rss_nsurf, runoffsnow, & ! output
2019  runoff, runoffsoil, chang, changsnow, snowtosurf, state_id, ev_snow, &
2020  snowdepth, snowremoval, swe, ev, chsnow_tot, &
2021  ev_tot, qe_tot, runoff_tot, surf_chang_tot, &
2022  runoffpipes, mwstore, runoffwaterbody)
2023 
2024  !Actual updates here as xx_tstep variables not taken as input to snowcalc
2025  ev_per_tstep = ev_per_tstep + ev_tot
2026  qe_per_tstep = qe_per_tstep + qe_tot
2027  runoff_per_tstep = runoff_per_tstep + runoff_tot
2028  surf_chang_per_tstep = surf_chang_per_tstep + surf_chang_tot
2029  chsnow_per_interval = chsnow_per_interval + chsnow_tot
2030  ELSE
2031  snowfrac(is) = 0
2032  snowdens(is) = 0
2033  snowpack(is) = 0
2034  ENDIF
2035 
2036  !Store ev_tot for each surface
2037  evap(is) = ev_tot
2038  ELSE ! snow-free calculation
2039 
2040  capstore(is) = storedrainprm(6, is)
2041  !Calculates ev [mm]
2042  CALL cal_evap( &
2043  evapmethod, state_id(is), wetthresh(is), capstore(is), &!input
2044  vpd_hpa, avdens, avcp, qn_e, s_hpa, psyc_hpa, resistsurf, ra, rb, tlv, &
2045  rss_nsurf(is), ev, qe_surf) !output
2046 
2047  !Surface water balance and soil store updates (can modify ev, updates state_id)
2048  CALL cal_water_storage( &
2049  is, sfr, pipecapacity, runofftowater, pin, & ! input:
2050  wu_nsurf, &
2051  drain, addwater, addimpervious, nsh_real, stateold, addwaterrunoff, &
2052  pervfraction, addveg, soilstorecap, addwaterbody, flowchange, statelimit, runoffagimpervious, surpluswaterbody, &
2053  runoffagveg, runoffpipes, ev, soilstore_id, surplusevap, runoffwaterbody, &
2054  p_mm, chang, runoff, state_id)!output:
2055 
2056  !Store ev for each surface
2057  evap(is) = ev
2058 
2059  ! Sum evaporation from different surfaces to find total evaporation [mm]
2060  ev_per_tstep = ev_per_tstep + evap(is)*sfr(is)
2061 
2062  ! Sum latent heat flux from different surfaces to find total latent heat flux
2063  qe_per_tstep = qe_per_tstep + qe_surf*sfr(is)
2064 
2065  ! Sum change from different surfaces to find total change to surface state_id
2066  surf_chang_per_tstep = surf_chang_per_tstep + (state_id(is) - stateold(is))*sfr(is)
2067 
2068  ! Sum runoff from different surfaces to find total runoff
2069  runoff_per_tstep = runoff_per_tstep + runoff(is)*sfr(is)
2070 
2071  ! Calculate total state_id (including water body)
2072  state_per_tstep = state_per_tstep + state_id(is)*sfr(is)
2073 
2074  ! sum The total runoff from the area !!Check (HCW)
2075  runoff_per_interval = runoff_per_interval + (runoff(is)*sfr(is))
2076 
2077  IF (nonwaterfraction /= 0 .AND. is /= watersurf) THEN
2078  nwstate_per_tstep = nwstate_per_tstep + (state_id(is)*sfr(is)/nonwaterfraction)
2079  ENDIF
2080 
2081  changsnow(is) = 0
2082  runoffsnow(is) = 0
2083 
2084  ENDIF
2085  ENDDO !end loop over surfaces
2086 
2087  qe = qe_per_tstep
2088 
2089  ! Calculate volume of water that will move between grids
2090  ! Volume [m3] = Depth relative to whole area [mm] / 1000 [mm m-1] * SurfaceArea [m2]
2091  ! Need to use these volumes when converting back to addImpervious, AddVeg and AddWater
2092  runoffagimpervious_m3 = runoffagimpervious/1000*surfacearea
2093  runoffagveg_m3 = runoffagveg/1000*surfacearea
2094  runoffwaterbody_m3 = runoffwaterbody/1000*surfacearea
2095  runoffpipes_m3 = runoffpipes/1000*surfacearea
2096 
2097  runoff_per_interval_out = runoff_per_interval
2098  state_id_out = state_id
2099  soilstore_id_out = soilstore_id
2100  snowpack_out = snowpack
2101  snowfrac_out = snowfrac
2102  snowwater_out = snowwater
2103  icefrac_out = icefrac
2104  snowdens_out = snowdens
2105 
2106  END SUBROUTINE suews_cal_qe
2107  !========================================================================
2108 
2109  !===============sensible heat flux======================================
2110  SUBROUTINE suews_cal_qh( &
2111  QHMethod, &!input
2112  qn1, qf, QmRain, qeOut, qs, QmFreez, qm, avdens, avcp, tsurf, Temp_C, RA, &
2113  qh, qh_residual, qh_resist)!output
2114  IMPLICIT NONE
2115 
2116  INTEGER, INTENT(in) :: QHMethod ! option for QH calculation: 1, residual; 2, resistance-based
2117 
2118  REAL(KIND(1d0)), INTENT(in)::qn1
2119  REAL(KIND(1d0)), INTENT(in)::qf
2120  REAL(KIND(1d0)), INTENT(in)::QmRain
2121  REAL(KIND(1d0)), INTENT(in)::qeOut
2122  REAL(KIND(1d0)), INTENT(in)::qs
2123  REAL(KIND(1d0)), INTENT(in)::QmFreez
2124  REAL(KIND(1d0)), INTENT(in)::qm
2125  REAL(KIND(1d0)), INTENT(in)::avdens
2126  REAL(KIND(1d0)), INTENT(in)::avcp
2127  REAL(KIND(1d0)), INTENT(in)::tsurf
2128  REAL(KIND(1d0)), INTENT(in)::Temp_C
2129  REAL(KIND(1d0)), INTENT(in)::RA
2130 
2131  REAL(KIND(1d0)), INTENT(out)::qh
2132  REAL(KIND(1d0)), INTENT(out)::qh_resist
2133  REAL(KIND(1d0)), INTENT(out)::qh_residual
2134 
2135  REAL(KIND(1d0)), PARAMETER::NAN = -999
2136 
2137  ! Calculate sensible heat flux as a residual (Modified by LJ in Nov 2012)
2138  qh_residual = (qn1 + qf + qmrain) - (qeout + qs + qm + qmfreez) !qh=(qn1+qf+QmRain+QmFreez)-(qeOut+qs+Qm)
2139 
2140  ! ! Calculate QH using resistance method (for testing HCW 06 Jul 2016)
2141  ! Aerodynamic-Resistance-based method
2142  IF (ra /= 0) THEN
2143  qh_resist = avdens*avcp*(tsurf - temp_c)/ra
2144  ELSE
2145  qh_resist = nan
2146  ENDIF
2147 
2148  ! choose output QH
2149  SELECT CASE (qhmethod)
2150  CASE (1)
2151  qh = qh_residual
2152  CASE (2)
2153  qh = qh_resist
2154  END SELECT
2155 
2156  END SUBROUTINE suews_cal_qh
2157  !========================================================================
2158 
2159  !===============Resistance Calculations=======================
2160  SUBROUTINE suews_cal_resistance( &
2161  StabilityMethod, &!input:
2162  Diagnose, AerodynamicResistanceMethod, RoughLenHeatMethod, snowUse, &
2163  id, it, gsModel, SMDMethod, &
2164  avdens, avcp, QH_init, zzd, z0m, zdm, &
2165  avU1, Temp_C, VegFraction, &
2166  avkdn, Kmax, G1, G2, G3, G4, G5, G6, S1, S2, TH, TL, dq, &
2167  xsmd, vsmd, MaxConductance, LAIMax, LAI_id, SnowFrac, sfr, &
2168  UStar, TStar, L_mod, &!output
2169  zL, gsc, ResistSurf, RA, RAsnow, rb)
2171  IMPLICIT NONE
2172 
2173  INTEGER, INTENT(in)::StabilityMethod
2174  INTEGER, INTENT(in)::Diagnose
2175  INTEGER, INTENT(in)::AerodynamicResistanceMethod
2176  INTEGER, INTENT(in)::RoughLenHeatMethod
2177  INTEGER, INTENT(in)::snowUse
2178  INTEGER, INTENT(in)::id
2179  INTEGER, INTENT(in)::it !time: day of year and hour
2180  INTEGER, INTENT(in)::gsModel !Choice of gs parameterisation (1 = Ja11, 2 = Wa16)
2181  INTEGER, INTENT(in)::SMDMethod!Method of measured soil moisture
2182 
2183  ! REAL(KIND(1d0)), INTENT(in)::qh_obs
2184  REAL(KIND(1d0)), INTENT(in)::avdens
2185  REAL(KIND(1d0)), INTENT(in)::avcp
2186  REAL(KIND(1d0)), INTENT(in)::QH_init
2187  REAL(KIND(1d0)), INTENT(in)::zzd !Active measurement height (meas. height-displac. height)
2188  REAL(KIND(1d0)), INTENT(in)::z0m !Aerodynamic roughness length
2189  REAL(KIND(1d0)), INTENT(in)::zdm !Displacement height
2190  REAL(KIND(1d0)), INTENT(in)::avU1 !Average wind speed
2191  REAL(KIND(1d0)), INTENT(in)::Temp_C !Air temperature
2192  REAL(KIND(1d0)), INTENT(in)::VegFraction!Fraction of vegetation
2193  REAL(KIND(1d0)), INTENT(in)::avkdn !Average downwelling shortwave radiation
2194  REAL(KIND(1d0)), INTENT(in)::Kmax !Annual maximum hourly solar radiation
2195  REAL(KIND(1d0)), INTENT(in)::G1 !Fitted parameters related to surface res. calculations
2196  REAL(KIND(1d0)), INTENT(in)::G2 !Fitted parameters related to surface res. calculations
2197  REAL(KIND(1d0)), INTENT(in)::G3 !Fitted parameters related to surface res. calculations
2198  REAL(KIND(1d0)), INTENT(in)::G4 !Fitted parameters related to surface res. calculations
2199  REAL(KIND(1d0)), INTENT(in)::G5 !Fitted parameters related to surface res. calculations
2200  REAL(KIND(1d0)), INTENT(in)::G6 !Fitted parameters related to surface res. calculations
2201  REAL(KIND(1d0)), INTENT(in)::S1 !Fitted parameters related to surface res. calculations
2202  REAL(KIND(1d0)), INTENT(in)::S2 !Fitted parameters related to surface res. calculations
2203  REAL(KIND(1d0)), INTENT(in)::TH !Maximum temperature limit
2204  REAL(KIND(1d0)), INTENT(in)::TL !Minimum temperature limit
2205  REAL(KIND(1d0)), INTENT(in)::dq !Specific humidity deficit
2206  REAL(KIND(1d0)), INTENT(in)::xsmd !Measured soil moisture deficit
2207  REAL(KIND(1d0)), INTENT(in)::vsmd !Soil moisture deficit for vegetated surfaces only (QUESTION: what about BSoil?)
2208 
2209  REAL(KIND(1d0)), DIMENSION(3), INTENT(in) ::MaxConductance!Max conductance [mm s-1]
2210  REAL(KIND(1d0)), DIMENSION(3), INTENT(in) ::LAIMax !Max LAI [m2 m-2]
2211  REAL(KIND(1d0)), DIMENSION(3), INTENT(in) ::LAI_id !=LAI_id(id-1,:), LAI for each veg surface [m2 m-2]
2212 
2213  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::SnowFrac !Surface fraction of snow cover
2214  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in)::sfr !Surface fractions [-]
2215 
2216  REAL(KIND(1d0)), INTENT(out)::TStar !T*
2217  REAL(KIND(1d0)), INTENT(out) ::UStar !Friction velocity
2218  REAL(KIND(1d0)), INTENT(out) ::zL !
2219  REAL(KIND(1d0)), INTENT(out) ::gsc !Surface Layer Conductance
2220  REAL(KIND(1d0)), INTENT(out) ::ResistSurf!Surface resistance
2221  REAL(KIND(1d0)), INTENT(out) ::RA !Aerodynamic resistance [s m^-1]
2222  REAL(KIND(1d0)), INTENT(out) ::RAsnow !Aerodynamic resistance for snow [s m^-1]
2223  REAL(KIND(1d0)), INTENT(out) ::rb !boundary layer resistance shuttleworth
2224  REAL(KIND(1d0)), INTENT(out) ::L_mod !Obukhov length
2225  REAL(KIND(1d0)) ::gfunc !gdq*gtemp*gs*gq for photosynthesis calculations
2226  ! REAL(KIND(1d0)) ::H_init !Kinematic sensible heat flux [K m s-1] used to calculate friction velocity
2227 
2228  ! Get first estimate of sensible heat flux. Modified by HCW 26 Feb 2015
2229  ! CALL SUEWS_init_QH( &
2230  ! avdens, avcp, QH_init, qn1, dectime, &
2231  ! H_init)
2232 
2233  IF (diagnose == 1) WRITE (*, *) 'Calling STAB_lumps...'
2234  !u* and Obukhov length out
2235  CALL cal_stab( &
2236  stabilitymethod, & ! input
2237  zzd, & !Active measurement height (meas. height-displac. height)
2238  z0m, & !Aerodynamic roughness length
2239  zdm, & !zero-plane displacement
2240  avu1, & !Average wind speed
2241  temp_c, & !Air temperature
2242  qh_init, & !sensible heat flux
2243  avdens, & ! air density
2244  avcp, & ! heat capacity of air
2245  l_mod, &! output: !Obukhov length
2246  tstar, & !T*, temperature scale
2247  ustar, & !Friction velocity
2248  zl)!Stability scale
2249 
2250  IF (diagnose == 1) WRITE (*, *) 'Calling AerodynamicResistance...'
2251  CALL aerodynamicresistance( &
2252  zzd, &! input:
2253  z0m, &
2254  avu1, &
2255  l_mod, &
2256  ustar, &
2257  vegfraction, &
2258  aerodynamicresistancemethod, &
2259  stabilitymethod, &
2260  roughlenheatmethod, &
2261  ra) ! output:
2262 
2263  IF (snowuse == 1) THEN
2264  IF (diagnose == 1) WRITE (*, *) 'Calling AerodynamicResistance for snow...'
2265  CALL aerodynamicresistance( &
2266  zzd, &! input:
2267  z0m, &
2268  avu1, &
2269  l_mod, &
2270  ustar, &
2271  vegfraction, &
2272  aerodynamicresistancemethod, &
2273  stabilitymethod, &
2274  3, &
2275  rasnow) ! output:
2276  ENDIF
2277 
2278  IF (diagnose == 1) WRITE (*, *) 'Calling SurfaceResistance...'
2279  ! CALL SurfaceResistance(id,it) !qsc and surface resistance out
2280  CALL surfaceresistance( &
2281  id, it, &! input:
2282  smdmethod, snowfrac, sfr, avkdn, temp_c, dq, xsmd, vsmd, maxconductance, &
2283  laimax, lai_id, gsmodel, kmax, &
2284  g1, g2, g3, g4, g5, g6, th, tl, s1, s2, &
2285  gfunc, gsc, resistsurf)! output:
2286 
2287  IF (diagnose == 1) WRITE (*, *) 'Calling BoundaryLayerResistance...'
2288  CALL boundarylayerresistance( &
2289  zzd, &! input: !Active measurement height (meas. height- zero-plane displacement)
2290  z0m, & !Aerodynamic roughness length
2291  avu1, & !Average wind speed
2292  ustar, & ! input/output:
2293  rb) ! output:
2294 
2295  END SUBROUTINE suews_cal_resistance
2296  !========================================================================
2297 
2298  !==============Update output arrays=========================
2299  SUBROUTINE suews_update_outputline( &
2300  AdditionalWater, alb, avkdn, avU10_ms, azimuth, &!input
2301  chSnow_per_interval, dectime, &
2302  drain_per_tstep, E_mod, ev_per_tstep, ext_wu, Fc, Fc_build, fcld, &
2303  Fc_metab, Fc_photo, Fc_respi, Fc_point, Fc_traff, FlowChange, &
2304  h_mod, id, imin, int_wu, it, iy, &
2305  kup, LAI_id, ldown, l_mod, lup, mwh, &
2306  MwStore, &
2307  nsh_real, NWstate_per_tstep, Precip, q2_gkg, &
2308  qeOut, qf, qh, qh_resist, Qm, QmFreez, &
2309  QmRain, qn1, qn1_S, qn1_snowfree, qs, RA, &
2310  resistsurf, RH2, runoffAGimpervious, runoffAGveg, &
2311  runoff_per_tstep, runoffPipes, runoffSoil_per_tstep, &
2312  runoffWaterBody, sfr, smd, smd_nsurf, SnowAlb, SnowRemoval, &
2313  state_id, state_per_tstep, surf_chang_per_tstep, swe, t2_C, tskin_C, &
2314  tot_chang_per_tstep, tsurf, UStar, &
2315  wu_nsurf, &
2316  z0m, zdm, zenith_deg, &
2317  datetimeLine, dataOutLineSUEWS)!output
2318  IMPLICIT NONE
2319 
2320  REAL(KIND(1d0)), PARAMETER :: NAN = -999
2321  INTEGER, INTENT(in) :: iy
2322  ! INTEGER,INTENT(in) :: iy_prev_t
2323  INTEGER, INTENT(in) :: id
2324  ! INTEGER,INTENT(in) :: id_prev_t
2325  INTEGER, INTENT(in) :: it
2326  INTEGER, INTENT(in) :: imin
2327  ! INTEGER, INTENT(in) :: Gridiv
2328  REAL(KIND(1d0)), INTENT(in) :: AdditionalWater
2329  REAL(KIND(1d0)), INTENT(in) :: alb(nsurf)
2330  REAL(KIND(1d0)), INTENT(in) :: avkdn
2331  REAL(KIND(1d0)), INTENT(in) :: avU10_ms
2332  REAL(KIND(1d0)), INTENT(in) :: azimuth
2333  REAL(KIND(1d0)), INTENT(in) :: chSnow_per_interval
2334  REAL(KIND(1d0)), INTENT(in) :: dectime
2335  REAL(KIND(1d0)), INTENT(in) :: drain_per_tstep
2336  REAL(KIND(1d0)), INTENT(in) :: E_mod
2337  REAL(KIND(1d0)), INTENT(in) :: ev_per_tstep
2338  REAL(KIND(1d0)), INTENT(in) :: ext_wu
2339  REAL(KIND(1d0)), INTENT(in) :: Fc
2340  REAL(KIND(1d0)), INTENT(in) :: Fc_build
2341  REAL(KIND(1d0)), INTENT(in) :: Fc_metab
2342  REAL(KIND(1d0)), INTENT(in) :: Fc_photo
2343  REAL(KIND(1d0)), INTENT(in) :: Fc_respi
2344  REAL(KIND(1d0)), INTENT(in) :: Fc_point
2345  REAL(KIND(1d0)), INTENT(in) :: Fc_traff
2346  REAL(KIND(1d0)), INTENT(in) :: fcld
2347  REAL(KIND(1d0)), INTENT(in) :: FlowChange
2348  REAL(KIND(1d0)), INTENT(in) :: h_mod
2349  REAL(KIND(1d0)), INTENT(in) :: int_wu
2350  REAL(KIND(1d0)), INTENT(in) :: kup
2351  REAL(KIND(1d0)), INTENT(in) :: l_mod
2352  REAL(KIND(1d0)), INTENT(in) :: LAI_id(nvegsurf)
2353  REAL(KIND(1d0)), INTENT(in) :: ldown
2354  REAL(KIND(1d0)), INTENT(in) :: lup
2355  REAL(KIND(1d0)), INTENT(in) :: mwh
2356  REAL(KIND(1d0)), INTENT(in) :: MwStore
2357  REAL(KIND(1d0)), INTENT(in) :: nsh_real
2358  REAL(KIND(1d0)), INTENT(in) :: NWstate_per_tstep
2359  REAL(KIND(1d0)), INTENT(in) :: Precip
2360  REAL(KIND(1d0)), INTENT(in) :: q2_gkg
2361  REAL(KIND(1d0)), INTENT(in) :: qeOut
2362  REAL(KIND(1d0)), INTENT(in) :: qf
2363  REAL(KIND(1d0)), INTENT(in) :: qh
2364  REAL(KIND(1d0)), INTENT(in) :: qh_resist
2365  REAL(KIND(1d0)), INTENT(in) :: Qm
2366  REAL(KIND(1d0)), INTENT(in) :: QmFreez
2367  REAL(KIND(1d0)), INTENT(in) :: QmRain
2368  REAL(KIND(1d0)), INTENT(in) :: qn1
2369  REAL(KIND(1d0)), INTENT(in) :: qn1_S
2370  REAL(KIND(1d0)), INTENT(in) :: qn1_snowfree
2371  REAL(KIND(1d0)), INTENT(in) :: qs
2372  REAL(KIND(1d0)), INTENT(in) :: RA
2373  REAL(KIND(1d0)), INTENT(in) :: resistsurf
2374  REAL(KIND(1d0)), INTENT(in) :: RH2
2375  REAL(KIND(1d0)), INTENT(in) :: runoff_per_tstep
2376  REAL(KIND(1d0)), INTENT(in) :: runoffAGimpervious
2377  REAL(KIND(1d0)), INTENT(in) :: runoffAGveg
2378  REAL(KIND(1d0)), INTENT(in) :: runoffPipes
2379  REAL(KIND(1d0)), INTENT(in) :: runoffSoil_per_tstep
2380  REAL(KIND(1d0)), INTENT(in) :: runoffWaterBody
2381  REAL(KIND(1d0)), INTENT(in) :: sfr(nsurf)
2382  REAL(KIND(1d0)), INTENT(in) :: smd
2383  REAL(KIND(1d0)), INTENT(in) :: smd_nsurf(nsurf)
2384  REAL(KIND(1d0)), INTENT(in) :: SnowAlb
2385  REAL(KIND(1d0)), INTENT(in) :: SnowRemoval(2)
2386  REAL(KIND(1d0)), INTENT(in) :: state_id(nsurf)
2387  REAL(KIND(1d0)), INTENT(in) :: state_per_tstep
2388  REAL(KIND(1d0)), INTENT(in) :: surf_chang_per_tstep
2389  REAL(KIND(1d0)), INTENT(in) :: swe
2390  REAL(KIND(1d0)), INTENT(in) :: t2_C
2391  REAL(KIND(1d0)), INTENT(in) :: tskin_C
2392  REAL(KIND(1d0)), INTENT(in) :: tot_chang_per_tstep
2393  REAL(KIND(1d0)), INTENT(in) :: tsurf
2394  REAL(KIND(1d0)), INTENT(in) :: UStar
2395  REAL(KIND(1d0)), DIMENSION(nsurf), INTENT(in) :: wu_nsurf
2396 
2397  REAL(KIND(1d0)), INTENT(in) :: z0m
2398  REAL(KIND(1d0)), INTENT(in) :: zdm
2399  REAL(KIND(1d0)), INTENT(in) :: zenith_deg
2400 
2401  REAL(KIND(1D0)), DIMENSION(5), INTENT(OUT)::datetimeLine
2402  REAL(KIND(1d0)), DIMENSION(ncolumnsDataOutSUEWS - 5), INTENT(out) :: dataOutLineSUEWS
2403  ! REAL(KIND(1d0)),DIMENSION(ncolumnsDataOutSnow-5),INTENT(out) :: dataOutLineSnow
2404  ! REAL(KIND(1d0)),DIMENSION(ncolumnsDataOutESTM-5),INTENT(out) :: dataOutLineESTM
2405  ! INTEGER:: is
2406  REAL(KIND(1d0)):: LAI_wt
2407  REAL(KIND(1d0)):: RH2_pct ! RH2 in percentage
2408 
2409  ! the variables below with '_x' endings stand for 'exported' values
2410  REAL(KIND(1d0))::ResistSurf_x
2411  REAL(KIND(1d0))::l_mod_x
2412  REAL(KIND(1d0))::bulkalbedo
2413  REAL(KIND(1d0))::smd_nsurf_x(nsurf)
2414  REAL(KIND(1d0))::state_x(nsurf)
2415  REAL(KIND(1d0)):: wu_DecTr
2416  REAL(KIND(1d0)):: wu_EveTr
2417  REAL(KIND(1d0)):: wu_Grass
2418 
2419  !=====================================================================
2420  !====================== Prepare data for output ======================
2421  ! values outside of reasonable range are set as NAN-like numbers. TS 10 Jun 2018
2422 
2423  ! Remove non-existing surface type from surface and soil outputs ! Added back in with NANs by HCW 24 Aug 2016
2424  state_x = unpack(spread(nan, dim=1, ncopies=SIZE(sfr)), mask=(sfr < 0.00001), field=state_id)
2425  smd_nsurf_x = unpack(spread(nan, dim=1, ncopies=SIZE(sfr)), mask=(sfr < 0.00001), field=smd_nsurf)
2426 
2427  resistsurf_x = min(9999., resistsurf)
2428 
2429  l_mod_x = max(min(9999., l_mod), -9999.)
2430 
2431  ! Calculate areally-weighted LAI
2432  ! IF(iy == (iy_prev_t +1) .AND. (id-1) == 0) THEN !Check for start of next year and avoid using LAI(id-1) as this is at the start of the year
2433  ! LAI_wt=DOT_PRODUCT(LAI(id_prev_t,:),sfr(1+2:nvegsurf+2))
2434  ! ELSE
2435  ! LAI_wt=DOT_PRODUCT(LAI(id-1,:),sfr(1+2:nvegsurf+2))
2436  ! ENDIF
2437 
2438  lai_wt = dot_product(lai_id(:), sfr(1 + 2:nvegsurf + 2))
2439 
2440  ! Calculate areally-weighted albedo
2441  bulkalbedo = dot_product(alb, sfr)
2442 
2443  ! convert RH2 to a percentage form
2444  rh2_pct = rh2*100.0
2445 
2446  ! translate water use to vegetated surfaces
2447  wu_dectr = wu_nsurf(3)
2448  wu_evetr = wu_nsurf(4)
2449  wu_grass = wu_nsurf(5)
2450 
2451  !====================== update output line ==============================
2452  ! date & time:
2453  datetimeline = [ &
2454  REAL(iy, KIND(1D0)), REAL(id, KIND(1D0)), &
2455  REAL(it, KIND(1D0)), REAL(imin, KIND(1D0)), dectime]
2456  !Define the overall output matrix to be printed out step by step
2457  dataoutlinesuews = [ &
2458  avkdn, kup, ldown, lup, tsurf, &
2459  qn1, qf, qs, qh, qeout, &
2460  h_mod, e_mod, qh_resist, &
2461  precip, ext_wu, ev_per_tstep, runoff_per_tstep, tot_chang_per_tstep, &
2462  surf_chang_per_tstep, state_per_tstep, nwstate_per_tstep, drain_per_tstep, smd, &
2463  flowchange/nsh_real, additionalwater, &
2464  runoffsoil_per_tstep, runoffpipes, runoffagimpervious, runoffagveg, runoffwaterbody, &
2465  int_wu, wu_evetr, wu_dectr, wu_grass, &
2466  smd_nsurf_x(1:nsurf - 1), &
2467  state_x(1:nsurf), &
2468  zenith_deg, azimuth, bulkalbedo, fcld, &
2469  lai_wt, z0m, zdm, &
2470  ustar, l_mod, ra, resistsurf, &
2471  fc, &
2472  fc_photo, fc_respi, fc_metab, fc_traff, fc_build, fc_point, &
2473  qn1_snowfree, qn1_s, snowalb, &
2474  qm, qmfreez, qmrain, swe, mwh, mwstore, chsnow_per_interval, &
2475  snowremoval(1:2), &
2476  tskin_c, t2_c, q2_gkg, avu10_ms, rh2_pct & ! surface-level diagonostics
2477  ]
2478  ! set invalid values to NAN
2479  ! dataOutLineSUEWS = set_nan(dataOutLineSUEWS)
2480 
2481  !====================update output line end==============================
2482 
2483  END SUBROUTINE suews_update_outputline
2484  !========================================================================
2485 
2486  !==============Update output arrays=========================
2487  SUBROUTINE suews_update_output( &
2488  SnowUse, storageheatmethod, &!input
2489  ReadLinesMetdata, NumberOfGrids, &
2490  ir, gridiv, &
2491  datetimeLine, dataOutLineSUEWS, dataOutLineSnow, dataOutLineESTM, dataoutLineRSL, dataOutLineSOLWEIG, &!input
2492  dataOutSUEWS, dataOutSnow, dataOutESTM, dataOutRSL, dataOutSOLWEIG)!inout
2493  IMPLICIT NONE
2494 
2495  INTEGER, INTENT(in) ::ReadLinesMetdata
2496  INTEGER, INTENT(in) ::NumberOfGrids
2497  INTEGER, INTENT(in) ::Gridiv
2498  INTEGER, INTENT(in) ::SnowUse
2499  INTEGER, INTENT(in) ::storageheatmethod
2500  INTEGER, INTENT(in) ::ir
2501 
2502  REAL(KIND(1d0)), DIMENSION(5), INTENT(in) :: datetimeLine
2503  REAL(KIND(1d0)), DIMENSION(ncolumnsDataOutSUEWS - 5), INTENT(in) :: dataOutLineSUEWS
2504  REAL(KIND(1d0)), DIMENSION(ncolumnsDataOutESTM - 5), INTENT(in) :: dataOutLineESTM
2505  REAL(KIND(1d0)), DIMENSION(ncolumnsDataOutSnow - 5), INTENT(in) :: dataOutLineSnow
2506  REAL(KIND(1d0)), DIMENSION(ncolumnsDataOutRSL - 5), INTENT(in) :: dataoutLineRSL
2507  REAL(KIND(1d0)), DIMENSION(ncolumnsdataOutSOL - 5), INTENT(in) :: dataOutLineSOLWEIG
2508 
2509  REAL(KIND(1d0)), INTENT(inout) :: dataOutSUEWS(readlinesmetdata, ncolumnsdataoutsuews, numberofgrids)
2510  REAL(KIND(1d0)), INTENT(inout) :: dataOutSnow(readlinesmetdata, ncolumnsdataoutsnow, numberofgrids)
2511  REAL(KIND(1d0)), INTENT(inout) :: dataOutESTM(readlinesmetdata, ncolumnsdataoutestm, numberofgrids)
2512  REAL(KIND(1d0)), INTENT(inout) :: dataOutRSL(readlinesmetdata, ncolumnsdataoutrsl, numberofgrids)
2513  REAL(KIND(1d0)), INTENT(inout) :: dataOutSOLWEIG(readlinesmetdata, ncolumnsdataoutrsl, numberofgrids)
2514 
2515  !====================== update output arrays ==============================
2516  !Define the overall output matrix to be printed out step by step
2517  dataoutsuews(ir, 1:ncolumnsdataoutsuews, gridiv) = [datetimeline, set_nan(dataoutlinesuews)]
2518  dataoutrsl(ir, 1:ncolumnsdataoutrsl, gridiv) = [datetimeline, set_nan(dataoutlinersl)]
2519  dataoutsolweig(ir, 1:ncolumnsdataoutsol, gridiv) = [datetimeline, set_nan(dataoutlinesolweig)]
2520  ! ! set invalid values to NAN
2521  ! dataOutSUEWS(ir,6:ncolumnsDataOutSUEWS,Gridiv)=set_nan(dataOutSUEWS(ir,6:ncolumnsDataOutSUEWS,Gridiv))
2522 
2523  IF (snowuse == 1) THEN
2524  dataoutsnow(ir, 1:ncolumnsdataoutsnow, gridiv) = [datetimeline, set_nan(dataoutlinesnow)]
2525  END IF
2526 
2527  IF (storageheatmethod == 4) THEN
2528  dataoutestm(ir, 1:ncolumnsdataoutestm, gridiv) = [datetimeline, set_nan(dataoutlineestm)]
2529  END IF
2530 
2531  !====================update output arrays end==============================
2532 
2533  END SUBROUTINE suews_update_output
2534 
2535  ! !========================================================================
2536  ! SUBROUTINE SUEWS_cal_Diagnostics( &
2537  ! dectime, &!input
2538  ! avU1, Temp_C, avRH, Press_hPa, &
2539  ! qh, qe, &
2540  ! VegFraction, zMeas, z0m, zdm, RA, avdens, avcp, lv_J_kg, tstep_real, &
2541  ! RoughLenHeatMethod, StabilityMethod, &
2542  ! avU10_ms, t2_C, q2_gkg, tskin_C, RH2)!output
2543  ! ! TS 03 Aug 2018: added limit on q2 by restricting RH2_max to 100%
2544  ! ! TS 31 Jul 2018: removed dependence on surface variables (Tsurf, qsat)
2545  ! ! TS 26 Jul 2018: improved the calculation logic
2546  ! ! TS 05 Sep 2017: improved interface
2547  ! ! TS 20 May 2017: calculate surface-level diagonostics
2548  ! IMPLICIT NONE
2549  ! REAL(KIND(1d0)), INTENT(in) ::dectime
2550  ! REAL(KIND(1d0)), INTENT(in) ::avU1, Temp_C, avRH
2551  ! REAL(KIND(1d0)), INTENT(in) ::qh
2552  ! REAL(KIND(1d0)), INTENT(in) ::Press_hPa, qe
2553  ! REAL(KIND(1d0)), INTENT(in) :: VegFraction, z0m, RA, avdens, avcp, lv_J_kg, tstep_real
2554  ! REAL(KIND(1d0)), INTENT(in) :: zMeas! height for measurement
2555  ! REAL(KIND(1d0)), INTENT(in) :: zdm ! displacement height
2556 
2557  ! ! INTEGER,INTENT(in) :: opt ! 0 for momentum, 1 for temperature, 2 for humidity
2558  ! INTEGER, INTENT(in) :: RoughLenHeatMethod, StabilityMethod
2559 
2560  ! REAL(KIND(1d0)), INTENT(out):: avU10_ms, t2_C, q2_gkg, tskin_C, RH2
2561  ! REAL(KIND(1d0))::qa_gkg
2562  ! REAL(KIND(1d0)), PARAMETER::k = 0.4
2563 
2564  ! ! wind speed:
2565  ! CALL diagSfc( &
2566  ! 0, &
2567  ! zMeas, avU1, 0d0, 10d0, avU10_ms, &
2568  ! VegFraction, &
2569  ! z0m, zdm, avdens, avcp, lv_J_kg, &
2570  ! avU1, Temp_C, qh, &
2571  ! RoughLenHeatMethod, StabilityMethod, tstep_real, dectime)
2572 
2573  ! ! temperature at 2 m agl:
2574  ! CALL diagSfc( &
2575  ! 1, &
2576  ! zMeas, Temp_C, qh, 2d0, t2_C, &
2577  ! VegFraction, &
2578  ! z0m, zdm, avdens, avcp, lv_J_kg, &
2579  ! avU1, Temp_C, qh, &
2580  ! RoughLenHeatMethod, StabilityMethod, tstep_real, dectime)
2581 
2582  ! ! skin temperature:
2583  ! tskin_C = qh/(avdens*avcp)*RA + temp_C
2584 
2585  ! ! humidity:
2586  ! qa_gkg = RH2qa(avRH/100, Press_hPa, Temp_c)
2587  ! CALL diagSfc( &
2588  ! 2, &
2589  ! zMeas, qa_gkg, qe, 2d0, q2_gkg, &
2590  ! VegFraction, &
2591  ! z0m, zdm, avdens, avcp, lv_J_kg, &
2592  ! avU1, Temp_C, qh, &
2593  ! RoughLenHeatMethod, StabilityMethod, tstep_real, dectime)
2594  ! ! re-examine if the diagnostic RH2 > 100% ?
2595  ! RH2 = qa2RH(q2_gkg, Press_hPa, Temp_c)
2596  ! IF (RH2 > 1) THEN
2597  ! ! if so, limit RH2 to 100%
2598  ! RH2 = 1d0
2599  ! ! and adjust the diagnostic q2_gkg
2600  ! q2_gkg = RH2qa(RH2, Press_hPa, Temp_c)
2601  ! END IF
2602 
2603  ! END SUBROUTINE SUEWS_cal_Diagnostics
2604 
2605  ! calculate several surface fraction related parameters
2606  SUBROUTINE suews_cal_surf( &
2607  sfr, & !input
2608  vegfraction, ImpervFraction, PervFraction, NonWaterFraction) ! output
2609  IMPLICIT NONE
2610 
2611  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN)::sfr
2612  REAL(KIND(1D0)), INTENT(OUT)::VegFraction
2613  REAL(KIND(1D0)), INTENT(OUT)::ImpervFraction
2614  REAL(KIND(1D0)), INTENT(OUT)::PervFraction
2615  REAL(KIND(1D0)), INTENT(OUT)::NonWaterFraction
2616 
2617  vegfraction = sfr(conifsurf) + sfr(decidsurf) + sfr(grasssurf)
2618  impervfraction = sfr(pavsurf) + sfr(bldgsurf)
2619  pervfraction = 1 - impervfraction
2620  nonwaterfraction = 1 - sfr(watersurf)
2621 
2622  END SUBROUTINE suews_cal_surf
2623 
2624  ! SUBROUTINE diagSfc( &
2625  ! opt, &
2626  ! zMeas, xMeas, xFlux, zDiag, xDiag, &
2627  ! VegFraction, &
2628  ! z0m, zd, avdens, avcp, lv_J_kg, &
2629  ! avU1, Temp_C, qh, &
2630  ! RoughLenHeatMethod, StabilityMethod, tstep_real, dectime)
2631  ! ! TS 31 Jul 2018: removed dependence on surface variables (Tsurf, qsat)
2632  ! ! TS 26 Jul 2018: improved the calculation logic
2633  ! ! TS 05 Sep 2017: improved interface
2634  ! ! TS 20 May 2017: calculate surface-level diagonostics
2635 
2636  ! IMPLICIT NONE
2637  ! REAL(KIND(1d0)), INTENT(in) :: dectime
2638  ! REAL(KIND(1d0)), INTENT(in) :: qh ! sensible heat flux
2639  ! REAL(KIND(1d0)), INTENT(in) :: z0m, avdens, avcp, lv_J_kg, tstep_real
2640  ! REAL(KIND(1d0)), INTENT(in) :: avU1, Temp_C ! atmospheric level variables
2641  ! REAL(KIND(1d0)), INTENT(in) :: zDiag ! height for diagonostics
2642  ! REAL(KIND(1d0)), INTENT(in) :: zMeas! height for measurement
2643  ! REAL(KIND(1d0)), INTENT(in) :: zd ! displacement height
2644  ! REAL(KIND(1d0)), INTENT(in) :: xMeas ! measurement at height
2645  ! REAL(KIND(1d0)), INTENT(in) :: xFlux!
2646  ! REAL(KIND(1d0)), INTENT(in) :: VegFraction ! vegetation fraction
2647 
2648  ! INTEGER, INTENT(in) :: opt ! 0 for momentum, 1 for temperature, 2 for humidity
2649  ! INTEGER, INTENT(in) :: RoughLenHeatMethod, StabilityMethod
2650 
2651  ! REAL(KIND(1d0)), INTENT(out):: xDiag
2652 
2653  ! REAL(KIND(1d0)) :: L_mod
2654  ! REAL(KIND(1d0)) :: psimz0, psihzDiag, psihzMeas, psihz0, psimzDiag ! stability correction functions
2655  ! REAL(KIND(1d0)) :: z0h ! Roughness length for heat
2656  ! REAL(KIND(1d0)) :: zDiagzd! height for diagnositcs
2657  ! REAL(KIND(1d0)) :: zMeaszd
2658  ! REAL(KIND(1d0)) :: tlv, H_kms, TStar, zL, UStar
2659  ! REAL(KIND(1d0)), PARAMETER :: muu = 1.46e-5 !molecular viscosity
2660  ! REAL(KIND(1d0)), PARAMETER :: nan = -999
2661  ! REAL(KIND(1d0)), PARAMETER :: zdm = 0 ! assuming Displacement height is ZERO
2662  ! REAL(KIND(1d0)), PARAMETER::k = 0.4
2663 
2664  ! tlv = lv_J_kg/tstep_real !Latent heat of vapourisation per timestep
2665  ! zDiagzd = zDiag + z0m ! height at hgtX assuming Displacement height is ZERO; set lower limit as z0 to prevent arithmetic error, zd=0
2666 
2667  ! ! get !Kinematic sensible heat flux [K m s-1] used to calculate friction velocity
2668  ! CALL SUEWS_init_QH( &
2669  ! avdens, avcp, qh, 0d0, dectime, & ! use qh as qh_obs to initialise H_init
2670  ! H_kms)!output
2671 
2672  ! ! redo the calculation for stability correction
2673  ! CALL cal_Stab( &
2674  ! ! input
2675  ! StabilityMethod, &
2676  ! dectime, & !Decimal time
2677  ! zDiagzd, & !Active measurement height (meas. height-displac. height)
2678  ! z0m, & !Aerodynamic roughness length
2679  ! zdm, & !Displacement height
2680  ! avU1, & !Average wind speed
2681  ! Temp_C, & !Air temperature
2682  ! H_kms, & !Kinematic sensible heat flux [K m s-1] used to calculate friction velocity
2683  ! ! output:
2684  ! L_MOD, & !Obukhov length
2685  ! TStar, & !T*
2686  ! UStar, & !Friction velocity
2687  ! zL)!Stability scale
2688 
2689  ! !***************************************************************
2690  ! ! log-law based stability corrections:
2691  ! ! Roughness length for heat
2692  ! z0h = cal_z0V(RoughLenHeatMethod, z0m, VegFraction, UStar)
2693 
2694  ! ! stability correction functions
2695  ! ! momentum:
2696  ! psimzDiag = stab_psi_mom(StabilityMethod, zDiagzd/L_mod)
2697  ! ! psimz2=stab_fn_mom(StabilityMethod,z2zd/L_mod,z2zd/L_mod)
2698  ! psimz0 = stab_psi_mom(StabilityMethod, z0m/L_mod)
2699 
2700  ! ! heat and vapor: assuming both are the same
2701  ! ! psihz2=stab_fn_heat(StabilityMethod,z2zd/L_mod,z2zd/L_mod)
2702  ! psihz0 = stab_psi_heat(StabilityMethod, z0h/L_mod)
2703 
2704  ! !***************************************************************
2705  ! SELECT CASE (opt)
2706  ! CASE (0) ! wind (momentum) at hgtX=10 m
2707  ! zDiagzd = zDiag + z0m! set lower limit as z0h to prevent arithmetic error, zd=0
2708 
2709  ! ! stability correction functions
2710  ! ! momentum:
2711  ! psimzDiag = stab_psi_mom(StabilityMethod, zDiagzd/L_mod)
2712  ! psimz0 = stab_psi_mom(StabilityMethod, z0m/L_mod)
2713  ! xDiag = UStar/k*(LOG(zDiagzd/z0m) - psimzDiag + psimz0) ! Brutsaert (2005), p51, eq.2.54
2714 
2715  ! CASE (1) ! temperature at hgtX=2 m
2716  ! zMeaszd = zMeas - zd
2717  ! zDiagzd = zDiag + z0h! set lower limit as z0h to prevent arithmetic error, zd=0
2718 
2719  ! ! heat and vapor: assuming both are the same
2720  ! psihzMeas = stab_psi_heat(StabilityMethod, zMeaszd/L_mod)
2721  ! psihzDiag = stab_psi_heat(StabilityMethod, zDiagzd/L_mod)
2722  ! ! psihz0=stab_fn_heat(StabilityMethod,z0h/L_mod,z0h/L_mod)
2723  ! xDiag = xMeas + xFlux/(k*UStar*avdens*avcp)*(LOG(zMeaszd/zDiagzd) - (psihzMeas - psihzDiag)) ! Brutsaert (2005), p51, eq.2.55
2724  ! ! IF ( ABS((LOG(z2zd/z0h)-psihz2+psihz0))>10 ) THEN
2725  ! ! PRINT*, '#####################################'
2726  ! ! PRINT*, 'xSurf',xSurf
2727  ! ! PRINT*, 'xFlux',xFlux
2728  ! ! PRINT*, 'k*us*avdens*avcp',k*us*avdens*avcp
2729  ! ! PRINT*, 'k',k
2730  ! ! PRINT*, 'us',us
2731  ! ! PRINT*, 'avdens',avdens
2732  ! ! PRINT*, 'avcp',avcp
2733  ! ! PRINT*, 'xFlux/X',xFlux/(k*us*avdens*avcp)
2734  ! ! PRINT*, 'stab',(LOG(z2zd/z0h)-psihz2+psihz0)
2735  ! ! PRINT*, 'LOG(z2zd/z0h)',LOG(z2zd/z0h)
2736  ! ! PRINT*, 'z2zd',z2zd,'L_mod',L_mod,'z0h',z0h
2737  ! ! PRINT*, 'z2zd/L_mod',z2zd/L_mod
2738  ! ! PRINT*, 'psihz2',psihz2
2739  ! ! PRINT*, 'psihz0',psihz0
2740  ! ! PRINT*, 'psihz2-psihz0',psihz2-psihz0
2741  ! ! PRINT*, 'xDiag',xDiag
2742  ! ! PRINT*, '*************************************'
2743  ! ! END IF
2744 
2745  ! CASE (2) ! humidity at hgtX=2 m
2746  ! zMeaszd = zMeas - zd
2747  ! zDiagzd = zDiag + z0h! set lower limit as z0h to prevent arithmetic error, zd=0
2748 
2749  ! ! heat and vapor: assuming both are the same
2750  ! psihzMeas = stab_psi_heat(StabilityMethod, zMeaszd/L_mod)
2751  ! psihzDiag = stab_psi_heat(StabilityMethod, zDiagzd/L_mod)
2752  ! ! psihz0=stab_fn_heat(StabilityMethod,z0h/L_mod,z0h/L_mod)
2753 
2754  ! xDiag = xMeas + xFlux/(k*UStar*avdens*tlv)*(LOG(zMeaszd/zDiagzd) - (psihzMeas - psihzDiag)) ! Brutsaert (2005), p51, eq.2.56
2755 
2756  ! END SELECT
2757 
2758  ! END SUBROUTINE diagSfc
2759 
2760  !===============set variable of invalid value to NAN=====================
2761  ELEMENTAL FUNCTION set_nan(x) RESULT(xx)
2762  IMPLICIT NONE
2763  REAL(KIND(1d0)), PARAMETER::pNAN = 30000 ! 30000 to prevent water_state being filtered out as it can be large
2764  REAL(KIND(1d0)), PARAMETER::NAN = -999
2765  REAL(KIND(1d0)), INTENT(in)::x
2766  REAL(KIND(1d0))::xx
2767 
2768  IF (abs(x) > pnan) THEN
2769  xx = nan
2770  ELSE
2771  xx = x
2772  ENDIF
2773 
2774  END FUNCTION set_nan
2775  !========================================================================
2776 
2777  !===============the functions below are only for test in f2py conversion===
2778  FUNCTION square(x) RESULT(xx)
2779  IMPLICIT NONE
2780  REAL(KIND(1d0)), PARAMETER::pNAN = 9999
2781  REAL(KIND(1d0)), PARAMETER::NAN = -999
2782  REAL(KIND(1d0)), INTENT(in)::x
2783  REAL(KIND(1d0))::xx
2784 
2785  xx = x**2 + nan/pnan
2786  xx = x**2
2787 
2788  END FUNCTION square
2789 
2790  FUNCTION square_real(x) RESULT(xx)
2791  IMPLICIT NONE
2792  REAL, PARAMETER::pNAN = 9999
2793  REAL, PARAMETER::NAN = -999
2794  REAL, INTENT(in)::x
2795  REAL::xx
2796 
2797  xx = x**2 + nan/pnan
2798  xx = x**2
2799 
2800  END FUNCTION square_real
2801 
2802  SUBROUTINE output_name_n(i, name, group, aggreg, outlevel)
2803  ! used by f2py module to handle output names
2804  IMPLICIT NONE
2805  ! the dimension is potentially incorrect,
2806  ! which should be consistent with that in output module
2807  INTEGER, INTENT(in) :: i
2808  CHARACTER(len=15), INTENT(out) :: name, group, aggreg
2809  INTEGER, INTENT(out) :: outlevel
2810 
2811  INTEGER :: nVar
2812  nvar = SIZE(varlistall, dim=1)
2813  IF (i < nvar .AND. i > 0) THEN
2814  name = trim(varlistall(i)%header)
2815  group = trim(varlistall(i)%group)
2816  aggreg = trim(varlistall(i)%aggreg)
2817  outlevel = varlistall(i)%level
2818  ELSE
2819  name = ''
2820  group = ''
2821  aggreg = ''
2822  outlevel = 0
2823  END IF
2824 
2825  END SUBROUTINE output_name_n
2826 
2827  SUBROUTINE output_size(nVar)
2828  ! used by f2py module to get size of the output list
2829  IMPLICIT NONE
2830  ! the dimension is potentially incorrect,
2831  ! which should be consistent with that in output module
2832  INTEGER, INTENT(out) :: nVar
2833 
2834  nvar = SIZE(varlistall, dim=1)
2835 
2836  END SUBROUTINE output_size
2837 
2838  SUBROUTINE suews_cal_multitsteps( &
2839  MetForcingBlock, len_sim, &
2840  AerodynamicResistanceMethod, AH_MIN, AHProf_24hr, AH_SLOPE_Cooling, & ! input&inout in alphabetical order
2841  AH_SLOPE_Heating, &
2842  alb, AlbMax_DecTr, AlbMax_EveTr, AlbMax_Grass, &
2843  AlbMin_DecTr, AlbMin_EveTr, AlbMin_Grass, &
2844  alpha_bioCO2, alpha_enh_bioCO2, alt, BaseT, BaseTe, &
2845  BaseTMethod, &
2846  BaseT_HC, beta_bioCO2, beta_enh_bioCO2, bldgH, CapMax_dec, CapMin_dec, &
2847  chAnOHM, CO2PointSource, cpAnOHM, CRWmax, CRWmin, DayWat, DayWatPer, &
2848  DecTreeH, Diagnose, DiagQN, DiagQS, DRAINRT, &
2849  dt_since_start, dqndt, qn1_av, dqnsdt, qn1_s_av, &
2850  EF_umolCO2perJ, emis, EmissionsMethod, EnEF_v_Jkm, endDLS, EveTreeH, FAIBldg, &
2851  FAIDecTree, FAIEveTree, Faut, FcEF_v_kgkm, FlowChange, &
2852  FrFossilFuel_Heat, FrFossilFuel_NonHeat, G1, G2, G3, G4, G5, G6, GDD_id, &
2853  GDDFull, Gridiv, gsModel, H_maintain, HDD_id, HumActivity_24hr, &
2854  IceFrac, Ie_a, Ie_end, Ie_m, Ie_start, &
2855  InternalWaterUse_h, &
2856  IrrFracPaved, IrrFracBldgs, &
2857  IrrFracEveTr, IrrFracDecTr, IrrFracGrass, &
2858  IrrFracBSoil, IrrFracWater, &
2859  EvapMethod, &
2860  kkAnOHM, Kmax, LAI_id, LAICalcYes, LAIMax, LAIMin, &
2861  LAIPower, LAIType, lat, lng, MaxConductance, MaxFCMetab, MaxQFMetab, &
2862  SnowWater, MinFCMetab, MinQFMetab, min_res_bioCO2, &
2863  NARP_EMIS_SNOW, NARP_TRANS_SITE, NetRadiationMethod, &
2864  OHM_coef, OHMIncQF, OHM_threshSW, &
2865  OHM_threshWD, PipeCapacity, PopDensDaytime, &
2866  PopDensNighttime, PopProf_24hr, PorMax_dec, PorMin_dec, &
2867  PrecipLimit, PrecipLimitAlb, &
2868  QF0_BEU, Qf_A, Qf_B, Qf_C, &
2869  RadMeltFact, RAINCOVER, RainMaxRes, resp_a, resp_b, &
2870  RoughLenHeatMethod, RoughLenMomMethod, RunoffToWater, S1, S2, &
2871  SatHydraulicConduct, SDDFull, SDD_id, sfr, SMDMethod, SnowAlb, SnowAlbMax, &
2872  SnowAlbMin, SnowPackLimit, SnowDens, SnowDensMax, SnowDensMin, SnowfallCum, SnowFrac, &
2873  SnowLimBldg, SnowLimPaved, SnowPack, SnowProf_24hr, snowUse, SoilDepth, &
2874  soilstore_id, SoilStoreCap, StabilityMethod, startDLS, state_id, StateLimit, &
2875  StorageHeatMethod, StoreDrainPrm, SurfaceArea, Tair_av, tau_a, tau_f, tau_r, &
2876  BaseT_Cooling, BaseT_Heating, TempMeltFact, TH, &
2877  theta_bioCO2, timezone, TL, TrafficRate, TrafficUnits, &
2878  Tmin_id, Tmax_id, lenday_id, &
2879  TraffProf_24hr, Ts5mindata_ir, tstep, tstep_prev, veg_type, &
2880  WaterDist, WaterUseMethod, WetThresh, &
2881  WUDay_id, DecidCap_id, albDecTr_id, albEveTr_id, albGrass_id, porosity_id, &
2882  WUProfA_24hr, WUProfM_24hr, Z, z0m_in, zdm_in, &
2883  dataOutBlockSUEWS, dataOutBlockSnow, dataOutBlockESTM, dataOutBlockRSL, dataOutBlockSOL, &!output
2884  DailyStateBlock)
2886  IMPLICIT NONE
2887  ! input:
2888  ! met forcing block
2889  REAL(KIND(1D0)), DIMENSION(len_sim, 24), INTENT(IN) ::MetForcingBlock
2890  INTEGER, INTENT(IN) :: len_sim
2891  ! input variables
2892  INTEGER, INTENT(IN)::AerodynamicResistanceMethod
2893  INTEGER, INTENT(IN)::BaseTMethod
2894  INTEGER, INTENT(IN)::Diagnose
2895  INTEGER, INTENT(IN)::DiagQN
2896  INTEGER, INTENT(IN)::DiagQS
2897  INTEGER, INTENT(IN)::startDLS
2898  INTEGER, INTENT(IN)::endDLS
2899  INTEGER, INTENT(IN)::EmissionsMethod
2900  INTEGER, INTENT(IN)::Gridiv
2901  INTEGER, INTENT(IN)::gsModel
2902  INTEGER, INTENT(IN)::Ie_end
2903  INTEGER, INTENT(IN)::Ie_start
2904  INTEGER, INTENT(IN)::EvapMethod
2905  INTEGER, INTENT(IN)::LAICalcYes
2906  INTEGER, INTENT(IN)::NetRadiationMethod
2907  INTEGER, INTENT(IN)::OHMIncQF
2908  INTEGER, INTENT(IN)::RoughLenHeatMethod
2909  INTEGER, INTENT(IN)::RoughLenMomMethod
2910  INTEGER, INTENT(IN)::SMDMethod
2911  INTEGER, INTENT(IN)::snowUse
2912  INTEGER, INTENT(IN)::StabilityMethod
2913  INTEGER, INTENT(IN)::StorageHeatMethod
2914  INTEGER, INTENT(IN)::tstep
2915  INTEGER, INTENT(IN)::tstep_prev ! tstep size of the previous step
2916  ! dt_since_start is intentionally made as inout to keep naming consistency with the embedded subroutine
2917  INTEGER, INTENT(inout)::dt_since_start ! time since simulation starts [s]
2918  INTEGER, INTENT(IN)::veg_type
2919  INTEGER, INTENT(IN)::WaterUseMethod
2920 
2921  INTEGER, DIMENSION(NVEGSURF), INTENT(IN)::LAIType
2922 
2923  REAL(KIND(1D0)), INTENT(IN)::AlbMax_DecTr
2924  REAL(KIND(1D0)), INTENT(IN)::AlbMax_EveTr
2925  REAL(KIND(1D0)), INTENT(IN)::AlbMax_Grass
2926  REAL(KIND(1D0)), INTENT(IN)::AlbMin_DecTr
2927  REAL(KIND(1D0)), INTENT(IN)::AlbMin_EveTr
2928  REAL(KIND(1D0)), INTENT(IN)::AlbMin_Grass
2929  REAL(KIND(1D0)), INTENT(IN)::alt
2930  ! REAL(KIND(1D0)),INTENT(IN)::avkdn
2931  ! REAL(KIND(1D0)),INTENT(IN)::avRh
2932  ! REAL(KIND(1D0)),INTENT(IN)::avU1
2933  REAL(KIND(1D0)), INTENT(IN)::BaseT_HC
2934  REAL(KIND(1D0)), INTENT(IN)::bldgH
2935  REAL(KIND(1D0)), INTENT(IN)::CapMax_dec
2936  REAL(KIND(1D0)), INTENT(IN)::CapMin_dec
2937  REAL(KIND(1D0)), INTENT(IN)::CO2PointSource
2938  REAL(KIND(1D0)), INTENT(IN)::CRWmax
2939  REAL(KIND(1D0)), INTENT(IN)::CRWmin
2940  REAL(KIND(1D0)), INTENT(IN)::DecTreeH
2941  REAL(KIND(1D0)), INTENT(IN)::DRAINRT
2942  REAL(KIND(1D0)), INTENT(IN)::EF_umolCO2perJ
2943  REAL(KIND(1D0)), INTENT(IN)::EnEF_v_Jkm
2944  REAL(KIND(1D0)), INTENT(IN)::EveTreeH
2945  REAL(KIND(1D0)), INTENT(IN)::FAIBldg
2946  REAL(KIND(1D0)), INTENT(IN)::FAIDecTree
2947  REAL(KIND(1D0)), INTENT(IN)::FAIEveTree
2948  REAL(KIND(1D0)), INTENT(IN)::Faut
2949  ! REAL(KIND(1D0)),INTENT(IN)::fcld_obs
2950  REAL(KIND(1D0)), INTENT(IN)::FlowChange
2951  REAL(KIND(1D0)), INTENT(IN)::FrFossilFuel_Heat
2952  REAL(KIND(1D0)), INTENT(IN)::FrFossilFuel_NonHeat
2953  REAL(KIND(1D0)), INTENT(IN)::G1
2954  REAL(KIND(1D0)), INTENT(IN)::G2
2955  REAL(KIND(1D0)), INTENT(IN)::G3
2956  REAL(KIND(1D0)), INTENT(IN)::G4
2957  REAL(KIND(1D0)), INTENT(IN)::G5
2958  REAL(KIND(1D0)), INTENT(IN)::G6
2959  REAL(KIND(1D0)), INTENT(IN)::H_maintain
2960  REAL(KIND(1D0)), INTENT(IN)::InternalWaterUse_h
2961  REAL(KIND(1D0)), INTENT(IN)::IrrFracPaved
2962  REAL(KIND(1D0)), INTENT(IN)::IrrFracBldgs
2963  REAL(KIND(1D0)), INTENT(IN)::IrrFracEveTr
2964  REAL(KIND(1D0)), INTENT(IN)::IrrFracDecTr
2965  REAL(KIND(1D0)), INTENT(IN)::IrrFracGrass
2966  REAL(KIND(1D0)), INTENT(IN)::IrrFracBSoil
2967  REAL(KIND(1D0)), INTENT(IN)::IrrFracWater
2968  REAL(KIND(1D0)), INTENT(IN)::Kmax
2969  ! REAL(KIND(1D0)),INTENT(IN)::LAI_obs
2970  REAL(KIND(1D0)), INTENT(IN)::lat
2971  ! REAL(KIND(1D0)),INTENT(IN)::ldown_obs
2972  REAL(KIND(1D0)), INTENT(IN)::lng
2973  REAL(KIND(1D0)), INTENT(IN)::MaxFCMetab
2974  REAL(KIND(1D0)), INTENT(IN)::MaxQFMetab
2975  REAL(KIND(1D0)), INTENT(IN)::MinFCMetab
2976  REAL(KIND(1D0)), INTENT(IN)::MinQFMetab
2977  REAL(KIND(1D0)), INTENT(IN)::NARP_EMIS_SNOW
2978  REAL(KIND(1D0)), INTENT(IN)::NARP_TRANS_SITE
2979  REAL(KIND(1D0)), INTENT(IN)::PipeCapacity
2980  REAL(KIND(1D0)), INTENT(IN)::PopDensNighttime
2981  REAL(KIND(1D0)), INTENT(IN)::PorMax_dec
2982  REAL(KIND(1D0)), INTENT(IN)::PorMin_dec
2983  ! REAL(KIND(1D0)),INTENT(IN)::Precip
2984  REAL(KIND(1D0)), INTENT(IN)::PrecipLimit
2985  REAL(KIND(1D0)), INTENT(IN)::PrecipLimitAlb
2986  ! REAL(KIND(1D0)),INTENT(IN)::Press_hPa
2987  ! REAL(KIND(1D0)),INTENT(IN)::qh_obs
2988  ! REAL(KIND(1D0)),INTENT(IN)::qn1_obs
2989  ! REAL(KIND(1D0)),INTENT(IN)::qs_obs
2990  ! REAL(KIND(1D0)),INTENT(IN)::qf_obs
2991  REAL(KIND(1D0)), INTENT(IN)::RadMeltFact
2992  REAL(KIND(1D0)), INTENT(IN)::RAINCOVER
2993  REAL(KIND(1D0)), INTENT(IN)::RainMaxRes
2994  REAL(KIND(1D0)), INTENT(IN)::RunoffToWater
2995  REAL(KIND(1D0)), INTENT(IN)::S1
2996  REAL(KIND(1D0)), INTENT(IN)::S2
2997  REAL(KIND(1D0)), INTENT(IN)::SnowAlbMax
2998  REAL(KIND(1D0)), INTENT(IN)::SnowAlbMin
2999  REAL(KIND(1D0)), INTENT(IN)::SnowDensMax
3000  REAL(KIND(1D0)), INTENT(IN)::SnowDensMin
3001  REAL(KIND(1D0)), INTENT(IN)::SnowLimBldg
3002  REAL(KIND(1D0)), INTENT(IN)::SnowLimPaved
3003  ! REAL(KIND(1D0)),INTENT(IN)::snowFrac_obs
3004  REAL(KIND(1D0)), INTENT(IN)::SurfaceArea
3005  REAL(KIND(1D0)), INTENT(IN)::tau_a
3006  REAL(KIND(1D0)), INTENT(IN)::tau_f
3007  REAL(KIND(1D0)), INTENT(IN)::tau_r
3008  ! REAL(KIND(1D0)),INTENT(IN)::Temp_C
3009  REAL(KIND(1D0)), INTENT(IN)::TempMeltFact
3010  REAL(KIND(1D0)), INTENT(IN)::TH
3011  REAL(KIND(1D0)), INTENT(IN)::timezone
3012  REAL(KIND(1D0)), INTENT(IN)::TL
3013  REAL(KIND(1D0)), INTENT(IN)::TrafficUnits
3014  ! REAL(KIND(1D0)),INTENT(IN)::xsmd
3015  REAL(KIND(1D0)), INTENT(IN)::Z
3016  REAL(KIND(1D0)), INTENT(IN)::z0m_in
3017  REAL(KIND(1D0)), INTENT(IN)::zdm_in
3018 
3019  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::AH_MIN
3020  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::AH_SLOPE_Cooling
3021  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::AH_SLOPE_Heating
3022  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::FcEF_v_kgkm
3023  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::QF0_BEU
3024  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::Qf_A
3025  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::Qf_B
3026  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::Qf_C
3027  ! REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::Numcapita
3028  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::PopDensDaytime
3029  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::BaseT_Cooling
3030  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::BaseT_Heating
3031  REAL(KIND(1D0)), DIMENSION(2), INTENT(IN) ::TrafficRate
3032  REAL(KIND(1D0)), DIMENSION(3), INTENT(IN) ::Ie_a
3033  REAL(KIND(1D0)), DIMENSION(3), INTENT(IN) ::Ie_m
3034  REAL(KIND(1D0)), DIMENSION(3), INTENT(IN) ::MaxConductance
3035  REAL(KIND(1D0)), DIMENSION(7), INTENT(IN) ::DayWat
3036  REAL(KIND(1D0)), DIMENSION(7), INTENT(IN) ::DayWatPer
3037  REAL(KIND(1D0)), DIMENSION(nsurf + 1), INTENT(IN)::OHM_threshSW
3038  REAL(KIND(1D0)), DIMENSION(nsurf + 1), INTENT(IN)::OHM_threshWD
3039  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::chAnOHM
3040  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::cpAnOHM
3041  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::emis
3042  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::kkAnOHM
3043  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::SatHydraulicConduct
3044  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::sfr
3045  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::SnowPackLimit
3046  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::SoilDepth
3047  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::SoilStoreCap
3048  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::StateLimit
3049  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(IN) ::WetThresh
3050  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::alpha_bioCO2
3051  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::alpha_enh_bioCO2
3052  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::BaseT
3053  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::BaseTe
3054  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::beta_bioCO2
3055  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::beta_enh_bioCO2
3056  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::GDDFull
3057  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::LAIMax
3058  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::LAIMin
3059  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::min_res_bioCO2
3060  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::resp_a
3061  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::resp_b
3062  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::SDDFull
3063  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::SnowProf_24hr
3064  REAL(KIND(1D0)), DIMENSION(NVEGSURF), INTENT(IN) ::theta_bioCO2
3065  REAL(KIND(1D0)), DIMENSION(4, NVEGSURF), INTENT(IN) ::LAIPower
3066  REAL(KIND(1D0)), DIMENSION(nsurf + 1, 4, 3), INTENT(IN) ::OHM_coef
3067  REAL(KIND(1D0)), DIMENSION(NSURF + 1, NSURF - 1), INTENT(IN) ::WaterDist
3068  REAL(KIND(1d0)), DIMENSION(:), INTENT(IN) ::Ts5mindata_ir
3069 
3070  ! diurnal profile values for 24hr
3071  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::AHProf_24hr
3072  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::HumActivity_24hr
3073  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::PopProf_24hr
3074  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::TraffProf_24hr
3075  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::WUProfA_24hr
3076  REAL(KIND(1D0)), DIMENSION(0:23, 2), INTENT(IN) ::WUProfM_24hr
3077  ! ########################################################################################
3078 
3079  ! ########################################################################################
3080  ! inout variables
3081  ! OHM related:
3082  REAL(KIND(1d0)), INTENT(INOUT) ::qn1_av
3083  REAL(KIND(1d0)), INTENT(INOUT) ::dqndt
3084  REAL(KIND(1d0)), INTENT(INOUT) ::qn1_s_av
3085  REAL(KIND(1d0)), INTENT(INOUT) ::dqnsdt
3086 
3087  ! snow related:
3088  REAL(KIND(1D0)), INTENT(INOUT) ::SnowfallCum
3089  REAL(KIND(1D0)), INTENT(INOUT) ::SnowAlb
3090  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT) ::IceFrac
3091  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT) ::SnowWater
3092  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT) ::SnowDens
3093  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT) ::SnowFrac
3094  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT) ::SnowPack
3095 
3096  ! water balance related:
3097  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT) ::soilstore_id
3098  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT) ::state_id
3099  REAL(KIND(1D0)), DIMENSION(6, NSURF), INTENT(INOUT) ::StoreDrainPrm
3100 
3101  ! phenology related:
3102  REAL(KIND(1D0)), DIMENSION(NSURF), INTENT(INOUT) ::alb
3103  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(INOUT) ::GDD_id !Growing Degree Days (see SUEWS_DailyState.f95)
3104  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(INOUT) ::SDD_id !Senescence Degree Days (see SUEWS_DailyState.f95)
3105  REAL(KIND(1d0)), DIMENSION(nvegsurf), INTENT(INOUT)::LAI_id !LAI for each veg surface [m2 m-2]
3106  REAL(KIND(1d0)), INTENT(INOUT) :: DecidCap_id
3107  REAL(KIND(1d0)), INTENT(INOUT) :: albDecTr_id
3108  REAL(KIND(1d0)), INTENT(INOUT) :: albEveTr_id
3109  REAL(KIND(1d0)), INTENT(INOUT) :: albGrass_id
3110  REAL(KIND(1d0)), INTENT(INOUT) :: porosity_id
3111  REAL(KIND(1d0)), INTENT(INOUT) :: Tmin_id
3112  REAL(KIND(1d0)), INTENT(INOUT) :: Tmax_id
3113  REAL(KIND(1d0)), INTENT(INOUT) :: lenday_id
3114 
3115  ! anthropogenic heat related:
3116  REAL(KIND(1d0)), DIMENSION(12), INTENT(INOUT) ::HDD_id !Heating Degree Days (see SUEWS_DailyState.f95)
3117 
3118  ! water use related:
3119  REAL(KIND(1d0)), DIMENSION(9), INTENT(INOUT) ::WUDay_id
3120 
3121  ! ESTM related:
3122  REAL(KIND(1d0)), INTENT(INOUT) ::Tair_av
3123  ! ########################################################################################
3124 
3125  ! ########################################################################################
3126  ! output variables
3127  ! REAL(KIND(1D0)),DIMENSION(:,:,:),ALLOCATABLE,INTENT(OUT) ::datetimeBlock
3128  REAL(KIND(1D0)), DIMENSION(len_sim, ncolumnsDataOutSUEWS), INTENT(OUT) ::dataOutBlockSUEWS
3129  REAL(KIND(1D0)), DIMENSION(len_sim, ncolumnsDataOutSnow), INTENT(OUT) ::dataOutBlockSnow
3130  REAL(KIND(1d0)), DIMENSION(len_sim, ncolumnsDataOutESTM), INTENT(OUT) ::dataOutBlockESTM
3131  REAL(KIND(1d0)), DIMENSION(len_sim, ncolumnsDataOutRSL), INTENT(OUT) ::dataOutBlockRSL
3132  REAL(KIND(1d0)), DIMENSION(len_sim, ncolumnsdataOutSOL), INTENT(OUT) ::dataOutBlockSOL
3133  REAL(KIND(1d0)), DIMENSION(len_sim, ncolumnsDataOutDailyState), INTENT(OUT) ::DailyStateBlock
3134  ! ########################################################################################
3135 
3136  ! internal temporal iteration related variables
3137  ! INTEGER::dt_since_start ! time since simulation starts [s]
3138 
3139  ! model output blocks of the same size as met forcing block
3140 
3141  ! local variables
3142  ! length of met forcing block
3143  INTEGER :: ir
3144  ! met forcing variables
3145  INTEGER :: iy
3146  INTEGER :: id
3147  INTEGER :: it
3148  INTEGER :: imin
3149  INTEGER :: isec
3150  INTEGER, PARAMETER :: gridiv_x = 1 ! a dummy gridiv as this routine is only one grid
3151  REAL(KIND(1D0))::qn1_obs
3152  REAL(KIND(1D0))::qh_obs
3153  REAL(KIND(1D0))::qe_obs
3154  REAL(KIND(1D0))::qs_obs
3155  REAL(KIND(1D0))::qf_obs
3156  REAL(KIND(1D0))::avu1
3157  REAL(KIND(1D0))::avrh
3158  REAL(KIND(1D0))::Temp_C
3159  REAL(KIND(1D0))::Press_hPa
3160  REAL(KIND(1D0))::Precip
3161  REAL(KIND(1D0))::avkdn
3162  REAL(KIND(1D0))::snowFrac_obs
3163  REAL(KIND(1D0))::ldown_obs
3164  REAL(KIND(1D0))::fcld_obs
3165  REAL(KIND(1D0))::wu_m3
3166  REAL(KIND(1D0))::xsmd
3167  REAL(KIND(1D0))::LAI_obs
3168  REAL(KIND(1D0))::kdiff
3169  REAL(KIND(1D0))::kdir
3170  REAL(KIND(1D0))::wdir
3171 
3172  REAL(KIND(1D0)), DIMENSION(5)::datetimeLine
3173  REAL(KIND(1D0)), DIMENSION(ncolumnsDataOutSUEWS - 5)::dataOutLineSUEWS
3174  REAL(KIND(1D0)), DIMENSION(ncolumnsDataOutSnow - 5)::dataOutLineSnow
3175  REAL(KIND(1d0)), DIMENSION(ncolumnsDataOutESTM - 5)::dataOutLineESTM
3176  REAL(KIND(1d0)), DIMENSION(ncolumnsDataOutRSL - 5)::dataOutLineRSL
3177  REAL(KIND(1D0)), DIMENSION(ncolumnsDataOutSol - 5) ::dataOutLineSOLWEIG
3178  REAL(KIND(1d0)), DIMENSION(ncolumnsDataOutDailyState - 5)::DailyStateLine
3179 
3180  REAL(KIND(1D0)), DIMENSION(len_sim, ncolumnsDataOutSUEWS, 1) ::dataOutBlockSUEWS_X
3181  REAL(KIND(1D0)), DIMENSION(len_sim, ncolumnsDataOutSnow, 1) ::dataOutBlockSnow_X
3182  REAL(KIND(1d0)), DIMENSION(len_sim, ncolumnsDataOutESTM, 1) ::dataOutBlockESTM_X
3183  REAL(KIND(1d0)), DIMENSION(len_sim, ncolumnsDataOutRSL, 1) ::dataOutBlockRSL_X
3184  REAL(KIND(1d0)), DIMENSION(len_sim, ncolumnsDataOutSOL, 1) ::dataOutBlockSOL_X
3185  ! REAL(KIND(1d0)),DIMENSION(len_sim,ncolumnsDataOutDailyState,1) ::DailyStateBlock_X
3186 
3187  REAL(KIND(1D0)), DIMENSION(10, 10) ::MetForcingData_grid ! fake array as a placeholder
3188 
3189  ! CHARACTER(len=150):: FileStateInit
3190  ! CHARACTER(len=4):: year_txt
3191  ! CHARACTER(len=3):: id_text
3192  ! CHARACTER(len=2):: it_text, imin_text
3193 
3194  ! get initial dt_since_start_x from dt_since_start, dt_since_start_x is used for Qn averaging. TS 28 Nov 2018
3195  ! dt_since_start = dt_since_start
3196 
3197  DO ir = 1, len_sim, 1
3198  ! =============================================================================
3199  ! === Translate met data from MetForcingBlock to variable names used in model ==
3200  ! =============================================================================
3201  iy = int(metforcingblock(ir, 1)) !Integer variables
3202  id = int(metforcingblock(ir, 2))
3203  it = int(metforcingblock(ir, 3))
3204  imin = int(metforcingblock(ir, 4))
3205  isec = 0 ! NOT used by SUEWS but by WRF-SUEWS via the cal_main interface
3206  qn1_obs = metforcingblock(ir, 5) !Real values (kind(1d0))
3207  qh_obs = metforcingblock(ir, 6)
3208  qe_obs = metforcingblock(ir, 7)
3209  qs_obs = metforcingblock(ir, 8)
3210  qf_obs = metforcingblock(ir, 9)
3211  avu1 = metforcingblock(ir, 10)
3212  avrh = metforcingblock(ir, 11)
3213  temp_c = metforcingblock(ir, 12)
3214  press_hpa = metforcingblock(ir, 13)
3215  precip = metforcingblock(ir, 14)
3216  avkdn = metforcingblock(ir, 15)
3217  snowfrac_obs = metforcingblock(ir, 16)
3218  ldown_obs = metforcingblock(ir, 17)
3219  fcld_obs = metforcingblock(ir, 18)
3220  wu_m3 = metforcingblock(ir, 19)
3221  xsmd = metforcingblock(ir, 20)
3222  lai_obs = metforcingblock(ir, 21)
3223  kdiff = metforcingblock(ir, 22)
3224  kdir = metforcingblock(ir, 23)
3225  wdir = metforcingblock(ir, 24)
3226 
3227  ! !================================================
3228  ! ! below is for debugging
3229  ! WRITE (year_txt, '(I4)') INT(iy)
3230  ! WRITE (id_text, '(I3)') INT(id)
3231  ! WRITE (it_text, '(I4)') INT(it)
3232  ! WRITE (imin_text, '(I4)') INT(imin)
3233 
3234  ! FileStateInit = './'//TRIM(ADJUSTL(year_txt))//'_'&
3235  ! //TRIM(ADJUSTL(id_text))//'_'&
3236  ! //TRIM(ADJUSTL(it_text))//'_'&
3237  ! //TRIM(ADJUSTL(imin_text))//'_'&
3238  ! //'state_init.nml'
3239 
3240  ! OPEN (12, file=FileStateInit, position='rewind')
3241 
3242  ! write (12, *) '&state_init'
3243  ! write (12, *) 'aerodynamicresistancemethod=', aerodynamicresistancemethod
3244  ! write (12, *) 'ah_min=', ah_min
3245  ! write (12, *) 'ahprof_24hr=', ahprof_24hr
3246  ! write (12, *) 'ah_slope_cooling=', ah_slope_cooling
3247  ! write (12, *) 'ah_slope_heating=', ah_slope_heating
3248  ! write (12, *) 'alb=', alb
3249  ! write (12, *) 'albmax_dectr=', albmax_dectr
3250  ! write (12, *) 'albmax_evetr=', albmax_evetr
3251  ! write (12, *) 'albmax_grass=', albmax_grass
3252  ! write (12, *) 'albmin_dectr=', albmin_dectr
3253  ! write (12, *) 'albmin_evetr=', albmin_evetr
3254  ! write (12, *) 'albmin_grass=', albmin_grass
3255  ! write (12, *) 'alpha_bioco2=', alpha_bioco2
3256  ! write (12, *) 'alpha_enh_bioco2=', alpha_enh_bioco2
3257  ! write (12, *) 'alt=', alt
3258  ! write (12, *) 'avkdn=', avkdn
3259  ! write (12, *) 'avrh=', avrh
3260  ! write (12, *) 'avu1=', avu1
3261  ! write (12, *) 'baset=', baset
3262  ! write (12, *) 'basete=', basete
3263  ! write (12, *) 'BaseT_HC=', BaseT_HC
3264  ! write (12, *) 'beta_bioco2=', beta_bioco2
3265  ! write (12, *) 'beta_enh_bioco2=', beta_enh_bioco2
3266  ! write (12, *) 'bldgh=', bldgh
3267  ! write (12, *) 'capmax_dec=', capmax_dec
3268  ! write (12, *) 'capmin_dec=', capmin_dec
3269  ! write (12, *) 'chanohm=', chanohm
3270  ! write (12, *) 'co2pointsource=', co2pointsource
3271  ! write (12, *) 'cpanohm=', cpanohm
3272  ! write (12, *) 'crwmax=', crwmax
3273  ! write (12, *) 'crwmin=', crwmin
3274  ! write (12, *) 'daywat=', daywat
3275  ! write (12, *) 'daywatper=', daywatper
3276  ! write (12, *) 'dectreeh=', dectreeh
3277  ! write (12, *) 'diagnose=', diagnose
3278  ! write (12, *) 'diagqn=', diagqn
3279  ! write (12, *) 'diagqs=', diagqs
3280  ! write (12, *) 'drainrt=', drainrt
3281  ! write (12, *) 'dt_since_start=', dt_since_start
3282  ! write (12, *) 'dqndt=', dqndt
3283  ! write (12, *) 'qn1_av=', qn1_av
3284  ! write (12, *) 'dqnsdt=', dqnsdt
3285  ! write (12, *) 'qn1_s_av=', qn1_s_av
3286  ! write (12, *) 'ef_umolco2perj=', ef_umolco2perj
3287  ! write (12, *) 'emis=', emis
3288  ! write (12, *) 'emissionsmethod=', emissionsmethod
3289  ! write (12, *) 'enef_v_jkm=', enef_v_jkm
3290  ! write (12, *) 'enddls=', enddls
3291  ! write (12, *) 'evetreeh=', evetreeh
3292  ! write (12, *) 'faibldg=', faibldg
3293  ! write (12, *) 'faidectree=', faidectree
3294  ! write (12, *) 'faievetree=', faievetree
3295  ! write (12, *) 'faut=', faut
3296  ! write (12, *) 'fcef_v_kgkm=', fcef_v_kgkm
3297  ! write (12, *) 'fcld_obs=', fcld_obs
3298  ! write (12, *) 'flowchange=', flowchange
3299  ! write (12, *) 'frfossilfuel_heat=', frfossilfuel_heat
3300  ! write (12, *) 'frfossilfuel_nonheat=', frfossilfuel_nonheat
3301  ! write (12, *) 'g1=', g1
3302  ! write (12, *) 'g2=', g2
3303  ! write (12, *) 'g3=', g3
3304  ! write (12, *) 'g4=', g4
3305  ! write (12, *) 'g5=', g5
3306  ! write (12, *) 'g6=', g6
3307  ! write (12, *) 'gdd_id=', gdd_id
3308  ! write (12, *) 'gddfull=', gddfull
3309  ! write (12, *) 'gridiv=', gridiv
3310  ! write (12, *) 'gsmodel=', gsmodel
3311  ! write (12, *) 'hdd_id=', hdd_id
3312  ! write (12, *) 'humactivity_24hr=', humactivity_24hr
3313  ! write (12, *) 'icefrac=', icefrac
3314  ! write (12, *) 'id=', id
3315  ! write (12, *) 'ie_a=', ie_a
3316  ! write (12, *) 'ie_end=', ie_end
3317  ! write (12, *) 'ie_m=', ie_m
3318  ! write (12, *) 'ie_start=', ie_start
3319  ! write (12, *) 'imin=', imin
3320  ! write (12, *) 'internalwateruse_h=', internalwateruse_h
3321  ! write (12, *) 'IrrFracEveTr=', IrrFracEveTr
3322  ! write (12, *) 'IrrFracDecTr=', IrrFracDecTr
3323  ! write (12, *) 'irrfracgrass=', irrfracgrass
3324  ! write (12, *) 'isec=', isec
3325  ! write (12, *) 'it=', it
3326  ! write (12, *) 'evapmethod=', evapmethod
3327  ! write (12, *) 'iy=', iy
3328  ! write (12, *) 'kkanohm=', kkanohm
3329  ! write (12, *) 'kmax=', kmax
3330  ! write (12, *) 'lai_id=', lai_id
3331  ! write (12, *) 'laicalcyes=', laicalcyes
3332  ! write (12, *) 'laimax=', laimax
3333  ! write (12, *) 'laimin=', laimin
3334  ! write (12, *) 'lai_obs=', lai_obs
3335  ! write (12, *) 'laipower=', laipower
3336  ! write (12, *) 'laitype=', laitype
3337  ! write (12, *) 'lat=', lat
3338  ! write (12, *) 'lenday_id=', lenday_id
3339  ! write (12, *) 'ldown_obs=', ldown_obs
3340  ! write (12, *) 'lng=', lng
3341  ! write (12, *) 'maxconductance=', maxconductance
3342  ! write (12, *) 'maxfcmetab=', maxfcmetab
3343  ! write (12, *) 'maxqfmetab=', maxqfmetab
3344  ! write (12, *) 'snowwater=', snowwater
3345  ! ! write (12, *) 'metforcingdata_grid=', metforcingdata_grid
3346  ! write (12, *) 'minfcmetab=', minfcmetab
3347  ! write (12, *) 'minqfmetab=', minqfmetab
3348  ! write (12, *) 'min_res_bioco2=', min_res_bioco2
3349  ! write (12, *) 'narp_emis_snow=', narp_emis_snow
3350  ! write (12, *) 'narp_trans_site=', narp_trans_site
3351  ! write (12, *) 'netradiationmethod=', netradiationmethod
3352  ! write (12, *) 'ohm_coef=', ohm_coef
3353  ! write (12, *) 'ohmincqf=', ohmincqf
3354  ! write (12, *) 'ohm_threshsw=', ohm_threshsw
3355  ! write (12, *) 'ohm_threshwd=', ohm_threshwd
3356  ! write (12, *) 'pipecapacity=', pipecapacity
3357  ! write (12, *) 'popdensdaytime=', popdensdaytime
3358  ! write (12, *) 'popdensnighttime=', popdensnighttime
3359  ! write (12, *) 'popprof_24hr=', popprof_24hr
3360  ! write (12, *) 'pormax_dec=', pormax_dec
3361  ! write (12, *) 'pormin_dec=', pormin_dec
3362  ! write (12, *) 'precip=', precip
3363  ! write (12, *) 'preciplimit=', preciplimit
3364  ! write (12, *) 'preciplimitalb=', preciplimitalb
3365  ! write (12, *) 'press_hpa=', press_hpa
3366  ! write (12, *) 'qf0_beu=', qf0_beu
3367  ! write (12, *) 'qf_a=', qf_a
3368  ! write (12, *) 'qf_b=', qf_b
3369  ! write (12, *) 'qf_c=', qf_c
3370  ! write (12, *) 'qn1_obs=', qn1_obs
3371  ! write (12, *) 'qh_obs=', qh_obs
3372  ! write (12, *) 'qs_obs=', qs_obs
3373  ! write (12, *) 'qf_obs=', qf_obs
3374  ! write (12, *) 'radmeltfact=', radmeltfact
3375  ! write (12, *) 'raincover=', raincover
3376  ! write (12, *) 'rainmaxres=', rainmaxres
3377  ! write (12, *) 'resp_a=', resp_a
3378  ! write (12, *) 'resp_b=', resp_b
3379  ! write (12, *) 'roughlenheatmethod=', roughlenheatmethod
3380  ! write (12, *) 'roughlenmommethod=', roughlenmommethod
3381  ! write (12, *) 'runofftowater=', runofftowater
3382  ! write (12, *) 's1=', s1
3383  ! write (12, *) 's2=', s2
3384  ! write (12, *) 'sathydraulicconduct=', sathydraulicconduct
3385  ! write (12, *) 'sddfull=', sddfull
3386  ! write (12, *) 'sdd_id=', sdd_id
3387  ! write (12, *) 'sfr=', sfr
3388  ! write (12, *) 'smdmethod=', smdmethod
3389  ! write (12, *) 'snowalb=', snowalb
3390  ! write (12, *) 'snowalbmax=', snowalbmax
3391  ! write (12, *) 'snowalbmin=', snowalbmin
3392  ! write (12, *) 'snowpacklimit=', snowpacklimit
3393  ! write (12, *) 'snowdens=', snowdens
3394  ! write (12, *) 'snowdensmax=', snowdensmax
3395  ! write (12, *) 'snowdensmin=', snowdensmin
3396  ! write (12, *) 'snowfallcum=', snowfallcum
3397  ! write (12, *) 'snowfrac=', snowfrac
3398  ! write (12, *) 'snowlimbldg=', snowlimbldg
3399  ! write (12, *) 'snowlimpaved=', snowlimpaved
3400  ! write (12, *) 'snowfrac_obs=', snowfrac_obs
3401  ! write (12, *) 'snowpack=', snowpack
3402  ! write (12, *) 'snowprof_24hr=', snowprof_24hr
3403  ! write (12, *) 'snowuse=', snowuse
3404  ! write (12, *) 'soildepth=', soildepth
3405  ! write (12, *) 'soilstore_id=', soilstore_id
3406  ! write (12, *) 'soilstorecap=', soilstorecap
3407  ! write (12, *) 'stabilitymethod=', stabilitymethod
3408  ! write (12, *) 'startdls=', startdls
3409  ! write (12, *) 'state_id=', state_id
3410  ! write (12, *) 'statelimit=', statelimit
3411  ! write (12, *) 'storageheatmethod=', storageheatmethod
3412  ! write (12, *) 'storedrainprm=', storedrainprm
3413  ! write (12, *) 'surfacearea=', surfacearea
3414  ! write (12, *) 'tair_av=', tair_av
3415  ! write (12, *) 'tau_a=', tau_a
3416  ! write (12, *) 'tau_f=', tau_f
3417  ! write (12, *) 'tau_r=', tau_r
3418  ! write (12, *) 'tmax_id=', tmax_id
3419  ! write (12, *) 'tmin_id=', tmin_id
3420  ! write (12, *) 'BaseT_Cooling=', BaseT_Cooling
3421  ! write (12, *) 'BaseT_Heating=', BaseT_Heating
3422  ! write (12, *) 'temp_c=', temp_c
3423  ! write (12, *) 'tempmeltfact=', tempmeltfact
3424  ! write (12, *) 'th=', th
3425  ! write (12, *) 'theta_bioco2=', theta_bioco2
3426  ! write (12, *) 'timezone=', timezone
3427  ! write (12, *) 'tl=', tl
3428  ! write (12, *) 'trafficrate=', trafficrate
3429  ! write (12, *) 'trafficunits=', trafficunits
3430  ! write (12, *) 'traffprof_24hr=', traffprof_24hr
3431  ! ! write (12, *) 'ts5mindata_ir=', ts5mindata_ir
3432  ! write (12, *) 'tstep=', tstep
3433  ! write (12, *) 'tstep_prev=', tstep_prev
3434  ! write (12, *) 'veg_type=', veg_type
3435  ! write (12, *) 'waterdist=', waterdist
3436  ! write (12, *) 'waterusemethod=', waterusemethod
3437  ! write (12, *) 'wetthresh=', wetthresh
3438  ! write (12, *) 'wu_m3=', wu_m3
3439  ! write (12, *) 'wuday_id=', wuday_id
3440  ! write (12, *) 'decidcap_id=', decidcap_id
3441  ! write (12, *) 'albdectr_id=', albdectr_id
3442  ! write (12, *) 'albevetr_id=', albevetr_id
3443  ! write (12, *) 'albgrass_id=', albgrass_id
3444  ! write (12, *) 'porosity_id=', porosity_id
3445  ! write (12, *) 'wuprofa_24hr=', wuprofa_24hr
3446  ! write (12, *) 'wuprofm_24hr=', wuprofm_24hr
3447  ! write (12, *) 'xsmd=', xsmd
3448  ! write (12, *) 'z=', z
3449  ! write (12, *) 'z0m_in=', z0m_in
3450  ! write (12, *) 'zdm_in=', zdm_in
3451  ! write (12, *) '/'
3452 
3453  ! WRITE (12, *) ''
3454 
3455  ! CLOSE (12)
3456  ! !================================================
3457 
3458  CALL suews_cal_main( &
3459  aerodynamicresistancemethod, ah_min, ahprof_24hr, ah_slope_cooling, & ! input&inout in alphabetical order
3460  ah_slope_heating, &
3461  alb, albmax_dectr, albmax_evetr, albmax_grass, &
3462  albmin_dectr, albmin_evetr, albmin_grass, &
3463  alpha_bioco2, alpha_enh_bioco2, alt, avkdn, avrh, avu1, baset, basete, &
3464  basetmethod, &
3465  baset_hc, beta_bioco2, beta_enh_bioco2, bldgh, capmax_dec, capmin_dec, &
3466  chanohm, co2pointsource, cpanohm, crwmax, crwmin, daywat, daywatper, &
3467  dectreeh, diagnose, diagqn, diagqs, drainrt, &
3468  dt_since_start, dqndt, qn1_av, dqnsdt, qn1_s_av, &
3469  ef_umolco2perj, emis, emissionsmethod, enef_v_jkm, enddls, evetreeh, faibldg, &
3470  faidectree, faievetree, faut, fcef_v_kgkm, fcld_obs, flowchange, &
3471  frfossilfuel_heat, frfossilfuel_nonheat, g1, g2, g3, g4, g5, g6, gdd_id, &
3472  gddfull, gridiv, gsmodel, h_maintain, hdd_id, humactivity_24hr, &
3473  icefrac, id, ie_a, ie_end, ie_m, ie_start, imin, &
3474  internalwateruse_h, &
3475  irrfracpaved, irrfracbldgs, &
3476  irrfracevetr, irrfracdectr, irrfracgrass, &
3477  irrfracbsoil, irrfracwater, &
3478  isec, it, evapmethod, &
3479  iy, kkanohm, kmax, lai_id, laicalcyes, laimax, laimin, lai_obs, &
3480  laipower, laitype, lat, lenday_id, ldown_obs, lng, maxconductance, maxfcmetab, maxqfmetab, &
3481  snowwater, metforcingdata_grid, minfcmetab, minqfmetab, min_res_bioco2, &
3482  narp_emis_snow, narp_trans_site, netradiationmethod, &
3483  ohm_coef, ohmincqf, ohm_threshsw, &
3484  ohm_threshwd, pipecapacity, popdensdaytime, &
3485  popdensnighttime, popprof_24hr, pormax_dec, pormin_dec, &
3486  precip, preciplimit, preciplimitalb, press_hpa, &
3487  qf0_beu, qf_a, qf_b, qf_c, &
3488  qn1_obs, qs_obs, qf_obs, &
3489  radmeltfact, raincover, rainmaxres, resp_a, resp_b, &
3490  roughlenheatmethod, roughlenmommethod, runofftowater, s1, s2, &
3491  sathydraulicconduct, sddfull, sdd_id, sfr, smdmethod, snowalb, snowalbmax, &
3492  snowalbmin, snowpacklimit, snowdens, snowdensmax, snowdensmin, snowfallcum, snowfrac, &
3493  snowlimbldg, snowlimpaved, snowfrac_obs, snowpack, snowprof_24hr, snowuse, soildepth, &
3494  soilstore_id, soilstorecap, stabilitymethod, startdls, state_id, statelimit, &
3495  storageheatmethod, storedrainprm, surfacearea, tair_av, tau_a, tau_f, tau_r, &
3496  tmax_id, tmin_id, &
3497  baset_cooling, baset_heating, temp_c, tempmeltfact, th, &
3498  theta_bioco2, timezone, tl, trafficrate, trafficunits, &
3499  traffprof_24hr, ts5mindata_ir, tstep, tstep_prev, veg_type, &
3500  waterdist, waterusemethod, wetthresh, wu_m3, &
3501  wuday_id, decidcap_id, albdectr_id, albevetr_id, albgrass_id, porosity_id, &
3502  wuprofa_24hr, wuprofm_24hr, xsmd, z, z0m_in, zdm_in, &
3503  datetimeline, dataoutlinesuews, dataoutlinesnow, dataoutlineestm, dataoutlinersl, dataoutlinesolweig, &!output
3504  dailystateline)!output
3505 
3506  ! update dt_since_start_x for next iteration, dt_since_start_x is used for Qn averaging. TS 28 Nov 2018
3507  dt_since_start = dt_since_start + tstep
3508 
3509  !============ update DailyStateBlock ===============
3510  dailystateblock(ir, :) = [datetimeline, dailystateline]
3511 
3512  !============ write out results ===============
3513  ! works at each timestep
3514  CALL suews_update_output( &
3515  snowuse, storageheatmethod, &!input
3516  len_sim, 1, &
3517  ir, gridiv_x, datetimeline, dataoutlinesuews, dataoutlinesnow, dataoutlineestm, dataoutlinersl, dataoutlinesolweig, &!input
3518  dataoutblocksuews_x, dataoutblocksnow_x, dataoutblockestm_x, dataoutblockrsl_x, dataoutblocksol_x)!inout
3519 
3520  END DO
3521 
3522  dataoutblocksuews = dataoutblocksuews_x(:, :, 1)
3523  dataoutblocksnow = dataoutblocksnow_x(:, :, 1)
3524  dataoutblockestm = dataoutblockestm_x(:, :, 1)
3525  dataoutblockrsl = dataoutblockrsl_x(:, :, 1)
3526  dataoutblocksol = dataoutblocksol_x(:, :, 1)
3527  ! DailyStateBlock=DailyStateBlock_X(:,:,1)
3528 
3529  END SUBROUTINE suews_cal_multitsteps
3530 
3531  ! a wrapper of NARP_cal_SunPosition used by supy
3532  SUBROUTINE suews_cal_sunposition( &
3533  year, idectime, UTC, locationlatitude, locationlongitude, locationaltitude, & !input
3534  sunazimuth, sunzenith) !output
3535  IMPLICIT NONE
3536 
3537  REAL(KIND(1D0)), INTENT(in) :: year, idectime, UTC, &
3538  locationlatitude, locationlongitude, locationaltitude
3539  REAL(KIND(1D0)), INTENT(out) ::sunazimuth, sunzenith
3540 
3541  CALL narp_cal_sunposition( &
3542  year, idectime, utc, locationlatitude, locationlongitude, locationaltitude, &
3543  sunazimuth, sunzenith)
3544 
3545  END SUBROUTINE suews_cal_sunposition
3546 
3547  ! function func(arg) result(retval)
3548  ! implicit none
3549  ! type :: arg
3550  ! type :: retval
3551 
3552  ! end function func
3553 
3554  function cal_tair_av(tair_av_prev, dt_since_start, tstep, temp_c) result(tair_av_next)
3555  ! calculate mean air temperature of past 24 hours
3556  ! TS, 17 Sep 2019
3557  implicit none
3558  real(KIND(1D0)), intent(in) :: tair_av_prev
3559  real(KIND(1D0)), intent(in) :: temp_c
3560  integer, intent(in) :: dt_since_start
3561  integer, intent(in) :: tstep
3562 
3563  real(KIND(1D0)) ::tair_av_next
3564 
3565  real(KIND(1D0)), parameter:: len_day_s = 24*3600 ! day length in seconds
3566  real(KIND(1D0)):: len_cal_s ! length of average period in seconds
3567  real(KIND(1D0)):: temp_k ! temp in K
3568 
3569  ! determine the average period
3570  if (dt_since_start > len_day_s) then
3571  ! if simulation has been running over one day
3572  len_cal_s = len_day_s
3573  else
3574  ! if simulation has been running less than one day
3575  len_cal_s = dt_since_start + tstep
3576  end if
3577  temp_k = temp_c + 273.15
3578  tair_av_next = tair_av_prev*(len_cal_s - tstep*1.)/len_cal_s + temp_k*tstep/len_cal_s
3579 
3580  end function cal_tair_av
3581 
3582  function cal_tsfc(qh, avdens, avcp, RA, temp_c) result(tsfc_C)
3583  ! calculate surface/skin temperature
3584  ! TS, 23 Oct 2019
3585  implicit none
3586  real(KIND(1D0)), intent(in) :: qh ! sensible heat flux [W m-2]
3587  real(KIND(1D0)), intent(in) :: avdens ! air density [kg m-3]
3588  real(KIND(1D0)), intent(in) :: avcp !air heat capacity [J m-3 K-1]
3589  real(KIND(1D0)), intent(in) :: RA !Aerodynamic resistance [s m^-1]
3590  real(KIND(1D0)), intent(in) :: temp_C ! air temperature [C]
3591 
3592  real(KIND(1D0)) ::tsfc_C ! surface temperature [C]
3593 
3594  tsfc_c = qh/(avdens*avcp)*ra + temp_c
3595  end function cal_tsfc
3596 
3597 END MODULE suews_driver
subroutine suews_update_soilmoist(NonWaterFraction, SoilStoreCap, sfr, soilstore_id, SoilMoistCap, SoilState, vsmd, smd)
subroutine anohm(tstep, dt_since_start, qn1, qn1_av_prev, dqndt_prev, qf, MetForcingData_grid, moist_surf, alb, emis, cpAnOHM, kkAnOHM, chAnOHM, sfr, nsurf, EmissionsMethod, id, Gridiv, qn1_av_next, dqndt_next, a1, a2, a3, qs, deltaQi)
High level wrapper for AnOHM calculation.
subroutine suews_cal_soilstate(SMDMethod, xsmd, NonWaterFraction, SoilMoistCap, SoilStoreCap, surf_chang_per_tstep, soilstore_id, soilstoreOld, sfr, smd, smd_nsurf, tot_chang_per_tstep, SoilState)
subroutine co2_biogen(alpha_bioCO2, alpha_enh_bioCO2, avkdn, beta_bioCO2, beta_enh_bioCO2, BSoilSurf, ConifSurf, DecidSurf, dectime, EmissionsMethod, gfunc, gfunc2, GrassSurf, gsmodel, id, it, ivConif, ivDecid, ivGrass, LAI_id, LAIMin, LAIMax, min_res_bioCO2, nsurf, NVegSurf, resp_a, resp_b, sfr, SnowFrac, t2, Temp_C, theta_bioCO2, Fc_biogen, Fc_photo, Fc_respi)
subroutine suews_cal_wateruse(nsh_real, wu_m3, SurfaceArea, sfr, IrrFracPaved, IrrFracBldgs, IrrFracEveTr, IrrFracDecTr, IrrFracGrass, IrrFracBSoil, IrrFracWater, DayofWeek_id, WUProfA_24hr, WUProfM_24hr, InternalWaterUse_h, HDD_id, WUDay_id, WaterUseMethod, NSH, it, imin, DLS, wu_nsurf, wu_int, wu_ext)
subroutine redistributewater(snowUse, WaterDist, sfr, Drain, AddWaterRunoff, AddWater)
subroutine suews_cal_biogenco2(alpha_bioCO2, alpha_enh_bioCO2, avkdn, avRh, beta_bioCO2, beta_enh_bioCO2, BSoilSurf, ConifSurf, DecidSurf, dectime, Diagnose, EmissionsMethod, Fc_anthro, G1, G2, G3, G4, G5, G6, gfunc, GrassSurf, gsmodel, id, it, ivConif, ivDecid, ivGrass, Kmax, LAI_id, LAIMin, LAIMax, MaxConductance, min_res_bioCO2, nsurf, NVegSurf, Press_hPa, resp_a, resp_b, S1, S2, sfr, SMDMethod, SnowFrac, t2_C, Temp_C, theta_bioCO2, TH, TL, vsmd, xsmd, Fc, Fc_biogen, Fc_photo, Fc_respi)
integer, parameter ncolumnsdataoutestm
subroutine suews_cal_dailystate(iy, id, it, imin, isec, tstep, tstep_prev, dt_since_start, DayofWeek_id, Tmin_id_prev, Tmax_id_prev, lenDay_id_prev, BaseTMethod, WaterUseMethod, Ie_start, Ie_end, LAICalcYes, LAIType, nsh_real, avkdn, Temp_C, Precip, BaseT_HC, BaseT_Heating, BaseT_Cooling, lat, Faut, LAI_obs, AlbMax_DecTr, AlbMax_EveTr, AlbMax_Grass, AlbMin_DecTr, AlbMin_EveTr, AlbMin_Grass, CapMax_dec, CapMin_dec, PorMax_dec, PorMin_dec, Ie_a, Ie_m, DayWatPer, DayWat, BaseT, BaseTe, GDDFull, SDDFull, LAIMin, LAIMax, LAIPower, DecidCap_id_prev, StoreDrainPrm_prev, LAI_id_prev, GDD_id_prev, SDD_id_prev, albDecTr_id_prev, albEveTr_id_prev, albGrass_id_prev, porosity_id_prev, HDD_id_prev, state_id, soilstore_id, SoilStoreCap, H_maintain, HDD_id_next, Tmin_id_next, Tmax_id_next, lenDay_id_next, albDecTr_id_next, albEveTr_id_next, albGrass_id_next, porosity_id_next, DecidCap_id_next, StoreDrainPrm_next, LAI_id_next, GDD_id_next, SDD_id_next, deltaLAI, WUDay_id)
subroutine suews_init_qh(avdens, avcp, h_mod, qn1, dectime, H_init)
integer, parameter ivgrass
subroutine narp_cal_sunposition(year, idectime, UTC, locationlatitude, locationlongitude, locationaltitude, sunazimuth, sunzenith)
integer, parameter ncolumnsdataoutdailystate
subroutine suews_cal_qe(Diagnose, snowuse, tstep, imin, it, EvapMethod, snowCalcSwitch, dayofWeek_id, CRWmin, CRWmax, dectime, avdens, avcp, lv_J_kg, lvS_J_kg, avRh, Press_hPa, Temp_C, RAsnow, psyc_hPa, sIce_hPa, PervFraction, vegfraction, addimpervious, qn1_snowfree, qf, qs, vpd_hPa, s_hPa, ResistSurf, RA, rb, snowdensmin, precip, PipeCapacity, RunoffToWater, NonWaterFraction, WU_nsurf, addVeg, addWaterBody, SnowLimPaved, SnowLimBldg, SurfaceArea, FlowChange, drain, WetThresh, stateOld, mw_ind, SoilStoreCap, rainonsnow, freezmelt, freezstate, freezstatevol, Qm_Melt, Qm_rain, Tsurf_ind, sfr, StateLimit, AddWater, addwaterrunoff, StoreDrainPrm, SnowPackLimit, SnowProf_24hr, SnowPack_in, SnowFrac_in, SnowWater_in, iceFrac_in, SnowDens_in, runoff_per_interval_in, state_id_in, soilstore_id_in, state_id_out, soilstore_id_out, SnowPack_out, SnowFrac_out, SnowWater_out, iceFrac_out, SnowDens_out, runoffSoil, SnowRemoval, state_per_tstep, NWstate_per_tstep, qe, swe, chSnow_per_interval, ev_per_tstep, runoff_per_tstep, surf_chang_per_tstep, runoffPipes, mwstore, runoffwaterbody, runoffAGveg, runoffAGimpervious)
integer, parameter bsoilsurf
subroutine suews_cal_water(Diagnose, snowUse, NonWaterFraction, addPipes, addImpervious, addVeg, addWaterBody, state_id, soilstore_id, sfr, StoreDrainPrm, WaterDist, nsh_real, drain_per_tstep, drain, frac_water2runoff, AdditionalWater, runoffPipes, runoff_per_interval, AddWater, stateOld, soilstoreOld)
subroutine estm(Gridiv, tstep, avkdn, avu1, temp_c, zenith_deg, avrh, press_hpa, ldown, bldgh, Ts5mindata_ir, Tair_av, dataOutLineESTM, QS)
real(kind(1d0)) function qa2rh(qa_gkg, pres_hPa, Ta_degC)
subroutine snow_cal_meltheat(snowUse, tstep, tau_r, SnowDensMax, lvS_J_kg, lv_J_kg, tstep_real, RadMeltFact, TempMeltFact, SnowAlbMax, SnowDensMin, Temp_C, Precip, PrecipLimit, PrecipLimitAlb, nsh_real, sfr, Tsurf_ind, Tsurf_ind_snow, state_id, qn1_ind_snow, kup_ind_snow, SnowWater, deltaQi, alb1, SnowPack_in, SnowFrac_in, SnowAlb_in, SnowDens_in, SnowfallCum_in, SnowPack_out, SnowFrac_out, SnowAlb_out, SnowDens_out, SnowfallCum_out, mwh, Qm, QmFreez, QmRain, veg_fr, snowCalcSwitch, Qm_melt, Qm_freezState, Qm_rain, FreezMelt, FreezState, FreezStateVol, rainOnSnow, SnowDepth, mw_ind, dataOutLineSnow)
subroutine snowcalc(tstep, imin, it, dectime, is, EvapMethod, CRWmin, CRWmax, nsh_real, lvS_J_kg, avdens, avRh, Press_hPa, Temp_C, RAsnow, psyc_hPa, avcp, sIce_hPa, PervFraction, vegfraction, addimpervious, vpd_hPa, qn_e, s_hPa, ResistSurf, RA, rb, tlv, snowdensmin, SnowProf_24hr, precip, PipeCapacity, RunoffToWater, addVeg, SnowLimPaved, SnowLimBldg, FlowChange, drain, WetThresh, stateOld, mw_ind, SoilStoreCap, rainonsnow, freezmelt, freezstate, freezstatevol, Qm_Melt, Qm_rain, Tsurf_ind, sfr, dayofWeek_id, StoreDrainPrm, SnowPackLimit, AddWater, addwaterrunoff, soilstore_id, SnowPack, SurplusEvap, SnowFrac, SnowWater, iceFrac, SnowDens, runoffAGimpervious, runoffAGveg, surplusWaterBody, rss_nsurf, runoffSnow, runoff, runoffSoil, chang, changSnow, SnowToSurf, state_id, ev_snow, SnowDepth, SnowRemoval, swe, ev, chSnow_tot, ev_tot, qe_tot, runoff_tot, surf_chang_tot, runoffPipes, mwstore, runoffwaterbody)
subroutine radmethod(NetRadiationMethod, snowUse, NetRadiationMethod_use, AlbedoChoice, ldown_option)
subroutine suews_cal_sunposition(year, idectime, UTC, locationlatitude, locationlongitude, locationaltitude, sunazimuth, sunzenith)
elemental real(kind(1d0)) function set_nan(x)
subroutine suews_cal_anthropogenicemission(AH_MIN, AHProf_24hr, AH_SLOPE_Cooling, AH_SLOPE_Heating, CO2PointSource, dayofWeek_id, DLS, EF_umolCO2perJ, EmissionsMethod, EnEF_v_Jkm, FcEF_v_kgkm, FrFossilFuel_Heat, FrFossilFuel_NonHeat, HDD_id, HumActivity_24hr, id, imin, it, MaxFCMetab, MaxQFMetab, MinFCMetab, MinQFMetab, nsh, PopDensDaytime, PopDensNighttime, PopProf_24hr, QF, QF0_BEU, Qf_A, Qf_B, Qf_C, QF_obs, QF_SAHP, SurfaceArea, BaseT_Cooling, BaseT_Heating, Temp_C, TrafficRate, TrafficUnits, TraffProf_24hr, Fc_anthro, Fc_build, Fc_metab, Fc_point, Fc_traff)
subroutine boundarylayerresistance(zzd, z0m, avU1, UStar, rb)
subroutine suews_update_dailystate(id, datetimeline, Gridiv, NumberOfGrids, DailyStateLine, dataOutDailyState)
subroutine cal_atmmoist(Temp_C, Press_hPa, avRh, dectime, lv_J_kg, lvS_J_kg, es_hPa, Ea_hPa, VPd_hpa, VPD_Pa, dq, dens_dry, avcp, air_dens)
integer, parameter nsurf
subroutine lumps_cal_qhqe(veg_type, snowUse, qn1, qf, qs, Qm, Temp_C, Veg_Fr, avcp, Press_hPa, lv_J_kg, tstep_real, DRAINRT, nsh_real, Precip, RainMaxRes, RAINCOVER, sfr, LAI_id_prev, LAImax, LAImin, QH_LUMPS, QE_LUMPS, psyc_hPa, s_hPa, sIce_hpa, Veg_Fr_temp, VegPhenLumps)
real(kind(1d0)) lv_j_kg
subroutine snowupdate(tstep, Temp_C, tau_a, tau_f, tau_r, SnowDensMax, SnowDensMin, SnowAlbMax, SnowAlbMin, SnowPack_prev, SnowAlb_prev, SnowDens_prev, SnowAlb_next, SnowDens_next)
subroutine suews_cal_qn(NetRadiationMethod, snowUse, tstep, SnowPack_prev, tau_a, tau_f, SnowAlbMax, SnowAlbMin, Diagnose, snowFrac_obs, ldown_obs, fcld_obs, dectime, ZENITH_deg, Tsurf_0, avKdn, Temp_C, avRH, ea_hPa, qn1_obs, SnowAlb_prev, snowFrac_prev, DiagQN, NARP_TRANS_SITE, NARP_EMIS_SNOW, IceFrac, sfr, emis, alb_prev, albDecTr_id, albEveTr_id, albGrass_id, alb_next, ldown, fcld, qn1, qn1_snowfree, qn1_S, kclear, kup, lup, tsurf, qn1_ind_snow, kup_ind_snow, Tsurf_ind_snow, Tsurf_ind, alb1, snowFrac_next, SnowAlb_next)
subroutine ohm(qn1, qn1_av_prev, dqndt_prev, qn1_av_next, dqndt_next, qn1_S, qn1_s_av_prev, dqnsdt_prev, qn1_s_av_next, dqnsdt_next, tstep, dt_since_start, sfr, nsurf, Tair_mav_5d, OHM_coef, OHM_threshSW, OHM_threshWD, soilstore_id, SoilStoreCap, state_id, BldgSurf, WaterSurf, SnowUse, SnowFrac, DiagQS, a1, a2, a3, qs, deltaQi)
integer, parameter ncolumnsdataoutsol
subroutine cal_stab(StabilityMethod, zzd, z0m, zdm, avU1, Temp_C, QH_init, avdens, avcp, L_MOD, TStar, UStar, zL)
subroutine anthropogenicemissions(CO2PointSource, EmissionsMethod, it, imin, DLS, DayofWeek_id, EF_umolCO2perJ, FcEF_v_kgkm, EnEF_v_Jkm, TrafficUnits, FrFossilFuel_Heat, FrFossilFuel_NonHeat, MinFCMetab, MaxFCMetab, MinQFMetab, MaxQFMetab, PopDensDaytime, PopDensNighttime, Temp_C, HDD_id, Qf_A, Qf_B, Qf_C, AH_MIN, AH_SLOPE_Heating, AH_SLOPE_Cooling, BaseT_Heating, BaseT_Cooling, TrafficRate, QF0_BEU, QF_SAHP, Fc_anthro, Fc_metab, Fc_traff, Fc_build, Fc_point, AHProf_24hr, HumActivity_24hr, TraffProf_24hr, PopProf_24hr, SurfaceArea)
integer, parameter conifsurf
subroutine suews_cal_tstep(tstep, nsh, nsh_real, tstep_real)
subroutine suews_cal_qh(QHMethod, qn1, qf, QmRain, qeOut, qs, QmFreez, qm, avdens, avcp, tsurf, Temp_C, RA, qh, qh_residual, qh_resist)
subroutine suews_cal_resistance(StabilityMethod, Diagnose, AerodynamicResistanceMethod, RoughLenHeatMethod, snowUse, id, it, gsModel, SMDMethod, avdens, avcp, QH_init, zzd, z0m, zdm, avU1, Temp_C, VegFraction, avkdn, Kmax, G1, G2, G3, G4, G5, G6, S1, S2, TH, TL, dq, xsmd, vsmd, MaxConductance, LAIMax, LAI_id, SnowFrac, sfr, UStar, TStar, L_mod, zL, gsc, ResistSurf, RA, RAsnow, rb)
subroutine surfaceresistance(id, it, SMDMethod, SnowFrac, sfr, avkdn, Temp_C, dq, xsmd, vsmd, MaxConductance, LAIMax, LAI_id, gsModel, Kmax, G1, G2, G3, G4, G5, G6, TH, TL, S1, S2, gfunc, gsc, ResistSurf)
real(kind(1d0)) avcp
real(kind(1d0)) function update_snow_albedo(tstep, SnowPack_prev, SnowAlb_prev, Temp_C, tau_a, tau_f, SnowAlbMax, SnowAlbMin)
AnOHM: Analytical Objective Hysteresis Model.
real(kind(1d0)) function square(x)
subroutine update_dailystateline(it, imin, nsh_real, GDD_id, HDD_id, LAI_id, SDD_id, Tmin_id, Tmax_id, lenday_id, DecidCap_id, albDecTr_id, albEveTr_id, albGrass_id, porosity_id, WUDay_id, deltaLAI, VegPhenLumps, SnowAlb, SnowDens, a1, a2, a3, DailyStateLine)
real(kind(1d0)) function cal_z0v(RoughLenHeatMethod, z0m, VegFraction, UStar)
real(kind(1d0)) function stab_psi_heat(StabilityMethod, ZL)
integer, parameter ncolumnsdataoutrsl
integer, parameter ncolumnsdataoutsnow
integer, parameter grasssurf
subroutine solweig_cal_main(id, it, dectime, lamdaP, lamdaF, avkdn, ldown, Temp_C, avrh, Press_hPa, Tg, lat, zenith_deg, azimuth, scale, alb_ground, alb_bldg, emis_ground, emis_wall, heightgravity, dataOutLineSOLWEIG)
real(kind(1d0)) function, dimension(nsurf) update_snow_dens(tstep, SnowPack_prev, SnowDens_prev, tau_r, SnowDensMax, SnowDensMin)
integer, parameter ivconif
subroutine suews_cal_main(AerodynamicResistanceMethod, AH_MIN, AHProf_24hr, AH_SLOPE_Cooling, AH_SLOPE_Heating, alb, AlbMax_DecTr, AlbMax_EveTr, AlbMax_Grass, AlbMin_DecTr, AlbMin_EveTr, AlbMin_Grass, alpha_bioCO2, alpha_enh_bioCO2, alt, avkdn, avRh, avU1, BaseT, BaseTe, BaseTMethod, BaseT_HC, beta_bioCO2, beta_enh_bioCO2, bldgH, CapMax_dec, CapMin_dec, chAnOHM, CO2PointSource, cpAnOHM, CRWmax, CRWmin, DayWat, DayWatPer, DecTreeH, Diagnose, DiagQN, DiagQS, DRAINRT, dt_since_start, dqndt, qn1_av, dqnsdt, qn1_s_av, EF_umolCO2perJ, emis, EmissionsMethod, EnEF_v_Jkm, endDLS, EveTreeH, FAIBldg, FAIDecTree, FAIEveTree, Faut, FcEF_v_kgkm, fcld_obs, FlowChange, FrFossilFuel_Heat, FrFossilFuel_NonHeat, G1, G2, G3, G4, G5, G6, GDD_id, GDDFull, Gridiv, gsModel, H_maintain, HDD_id, HumActivity_24hr, IceFrac, id, Ie_a, Ie_end, Ie_m, Ie_start, imin, InternalWaterUse_h, IrrFracPaved, IrrFracBldgs, IrrFracEveTr, IrrFracDecTr, IrrFracGrass, IrrFracBSoil, IrrFracWater, isec, it, EvapMethod, iy, kkAnOHM, Kmax, LAI_id, LAICalcYes, LAIMax, LAIMin, LAI_obs, LAIPower, LAIType, lat, lenDay_id, ldown_obs, lng, MaxConductance, MaxFCMetab, MaxQFMetab, SnowWater, MetForcingData_grid, MinFCMetab, MinQFMetab, min_res_bioCO2, NARP_EMIS_SNOW, NARP_TRANS_SITE, NetRadiationMethod, OHM_coef, OHMIncQF, OHM_threshSW, OHM_threshWD, PipeCapacity, PopDensDaytime, PopDensNighttime, PopProf_24hr, PorMax_dec, PorMin_dec, Precip, PrecipLimit, PrecipLimitAlb, Press_hPa, QF0_BEU, Qf_A, Qf_B, Qf_C, qn1_obs, qs_obs, qf_obs, RadMeltFact, RAINCOVER, RainMaxRes, resp_a, resp_b, RoughLenHeatMethod, RoughLenMomMethod, RunoffToWater, S1, S2, SatHydraulicConduct, SDDFull, SDD_id, sfr, SMDMethod, SnowAlb, SnowAlbMax, SnowAlbMin, SnowPackLimit, SnowDens, SnowDensMax, SnowDensMin, SnowfallCum, SnowFrac, SnowLimBldg, SnowLimPaved, snowFrac_obs, SnowPack, SnowProf_24hr, snowUse, SoilDepth, soilstore_id, SoilStoreCap, StabilityMethod, startDLS, state_id, StateLimit, StorageHeatMethod, StoreDrainPrm, SurfaceArea, Tair_av, tau_a, tau_f, tau_r, Tmax_id, Tmin_id, BaseT_Cooling, BaseT_Heating, Temp_C, TempMeltFact, TH, theta_bioCO2, timezone, TL, TrafficRate, TrafficUnits, TraffProf_24hr, Ts5mindata_ir, tstep, tstep_prev, veg_type, WaterDist, WaterUseMethod, WetThresh, wu_m3, WUDay_id, DecidCap_id, albDecTr_id, albEveTr_id, albGrass_id, porosity_id, WUProfA_24hr, WUProfM_24hr, xsmd, Z, z0m_in, zdm_in, datetimeLine, dataOutLineSUEWS, dataOutLineSnow, dataOutLineESTM, dataoutLineRSL, dataOutLineSOLWEIG, DailyStateLine)
type(varattr), dimension(500) varlistall
real(kind(1d0)) function rh2qa(RH_dec, pres_hPa, Ta_degC)
subroutine suews_cal_horizontalsoilwater(sfr, SoilStoreCap, SoilDepth, SatHydraulicConduct, SurfaceArea, NonWaterFraction, tstep_real, soilstore_id, runoffSoil, runoffSoil_per_tstep)
subroutine suews_cal_dectime(id, it, imin, isec, dectime)
integer, parameter nvegsurf
subroutine output_size(nVar)
subroutine suews_cal_weekday(iy, id, lat, dayofWeek_id)
subroutine suews_update_outputline(AdditionalWater, alb, avkdn, avU10_ms, azimuth, chSnow_per_interval, dectime, drain_per_tstep, E_mod, ev_per_tstep, ext_wu, Fc, Fc_build, fcld, Fc_metab, Fc_photo, Fc_respi, Fc_point, Fc_traff, FlowChange, h_mod, id, imin, int_wu, it, iy, kup, LAI_id, ldown, l_mod, lup, mwh, MwStore, nsh_real, NWstate_per_tstep, Precip, q2_gkg, qeOut, qf, qh, qh_resist, Qm, QmFreez, QmRain, qn1, qn1_S, qn1_snowfree, qs, RA, resistsurf, RH2, runoffAGimpervious, runoffAGveg, runoff_per_tstep, runoffPipes, runoffSoil_per_tstep, runoffWaterBody, sfr, smd, smd_nsurf, SnowAlb, SnowRemoval, state_id, state_per_tstep, surf_chang_per_tstep, swe, t2_C, tskin_C, tot_chang_per_tstep, tsurf, UStar, wu_nsurf, z0m, zdm, zenith_deg, datetimeLine, dataOutLineSUEWS)
subroutine cal_water_storage(is, sfr, PipeCapacity, RunoffToWater, pin, WU_nsurf, drain, AddWater, addImpervious, nsh_real, stateOld, AddWaterRunoff, PervFraction, addVeg, SoilStoreCap, addWaterBody, FlowChange, StateLimit, runoffAGimpervious, surplusWaterBody, runoffAGveg, runoffPipes, ev, soilstore_id, SurplusEvap, runoffWaterBody, p_mm, chang, runoff, state_id)
subroutine cal_evap(EvapMethod, state_is, WetThresh_is, capStore_is, vpd_hPa, avdens, avcp, qn_e, s_hPa, psyc_hPa, ResistSurf, RA, rb, tlv, rss, ev, qe)
subroutine suews_update_output(SnowUse, storageheatmethod, ReadLinesMetdata, NumberOfGrids, ir, gridiv, datetimeLine, dataOutLineSUEWS, dataOutLineSnow, dataOutLineESTM, dataoutLineRSL, dataOutLineSOLWEIG, dataOutSUEWS, dataOutSnow, dataOutESTM, dataOutRSL, dataOutSOLWEIG)
subroutine narp(nsurf, sfr, SnowFrac, alb, emis, IceFrac, NARP_TRANS_SITE, NARP_EMIS_SNOW, DTIME, ZENITH_deg, tsurf_0, kdown, Temp_C, RH, Press_hPa, qn1_obs, SnowAlb, AlbedoChoice, ldown_option, NetRadiationMethod_use, DiagQN, QSTARall, QSTAR_SF, QSTAR_S, kclear, KUPall, LDOWN, LUPall, fcld, TSURFall, qn1_ind_snow, kup_ind_snow, Tsurf_ind_snow, Tsurf_ind, alb0, alb1)
subroutine output_name_n(i, name, group, aggreg, outlevel)
integer, parameter ncolumnsdataoutsuews
subroutine rslprofile(Zh, z0m, zdm, L_MOD, sfr, FAI, StabilityMethod, avcp, lv_J_kg, avdens, avU1, Temp_C, avRH, Press_hPa, zMeas, qh, qe, T2_C, q2_gkg, U10_ms, RH2, dataoutLineRSL)
subroutine suews_cal_multitsteps(MetForcingBlock, len_sim, AerodynamicResistanceMethod, AH_MIN, AHProf_24hr, AH_SLOPE_Cooling, AH_SLOPE_Heating, alb, AlbMax_DecTr, AlbMax_EveTr, AlbMax_Grass, AlbMin_DecTr, AlbMin_EveTr, AlbMin_Grass, alpha_bioCO2, alpha_enh_bioCO2, alt, BaseT, BaseTe, BaseTMethod, BaseT_HC, beta_bioCO2, beta_enh_bioCO2, bldgH, CapMax_dec, CapMin_dec, chAnOHM, CO2PointSource, cpAnOHM, CRWmax, CRWmin, DayWat, DayWatPer, DecTreeH, Diagnose, DiagQN, DiagQS, DRAINRT, dt_since_start, dqndt, qn1_av, dqnsdt, qn1_s_av, EF_umolCO2perJ, emis, EmissionsMethod, EnEF_v_Jkm, endDLS, EveTreeH, FAIBldg, FAIDecTree, FAIEveTree, Faut, FcEF_v_kgkm, FlowChange, FrFossilFuel_Heat, FrFossilFuel_NonHeat, G1, G2, G3, G4, G5, G6, GDD_id, GDDFull, Gridiv, gsModel, H_maintain, HDD_id, HumActivity_24hr, IceFrac, Ie_a, Ie_end, Ie_m, Ie_start, InternalWaterUse_h, IrrFracPaved, IrrFracBldgs, IrrFracEveTr, IrrFracDecTr, IrrFracGrass, IrrFracBSoil, IrrFracWater, EvapMethod, kkAnOHM, Kmax, LAI_id, LAICalcYes, LAIMax, LAIMin, LAIPower, LAIType, lat, lng, MaxConductance, MaxFCMetab, MaxQFMetab, SnowWater, MinFCMetab, MinQFMetab, min_res_bioCO2, NARP_EMIS_SNOW, NARP_TRANS_SITE, NetRadiationMethod, OHM_coef, OHMIncQF, OHM_threshSW, OHM_threshWD, PipeCapacity, PopDensDaytime, PopDensNighttime, PopProf_24hr, PorMax_dec, PorMin_dec, PrecipLimit, PrecipLimitAlb, QF0_BEU, Qf_A, Qf_B, Qf_C, RadMeltFact, RAINCOVER, RainMaxRes, resp_a, resp_b, RoughLenHeatMethod, RoughLenMomMethod, RunoffToWater, S1, S2, SatHydraulicConduct, SDDFull, SDD_id, sfr, SMDMethod, SnowAlb, SnowAlbMax, SnowAlbMin, SnowPackLimit, SnowDens, SnowDensMax, SnowDensMin, SnowfallCum, SnowFrac, SnowLimBldg, SnowLimPaved, SnowPack, SnowProf_24hr, snowUse, SoilDepth, soilstore_id, SoilStoreCap, StabilityMethod, startDLS, state_id, StateLimit, StorageHeatMethod, StoreDrainPrm, SurfaceArea, Tair_av, tau_a, tau_f, tau_r, BaseT_Cooling, BaseT_Heating, TempMeltFact, TH, theta_bioCO2, timezone, TL, TrafficRate, TrafficUnits, Tmin_id, Tmax_id, lenday_id, TraffProf_24hr, Ts5mindata_ir, tstep, tstep_prev, veg_type, WaterDist, WaterUseMethod, WetThresh, WUDay_id, DecidCap_id, albDecTr_id, albEveTr_id, albGrass_id, porosity_id, WUProfA_24hr, WUProfM_24hr, Z, z0m_in, zdm_in, dataOutBlockSUEWS, dataOutBlockSnow, dataOutBlockESTM, dataOutBlockRSL, dataOutBlockSOL, DailyStateBlock)
subroutine suews_cal_qs(StorageHeatMethod, qs_obs, OHMIncQF, Gridiv, id, tstep, dt_since_start, Diagnose, sfr, OHM_coef, OHM_threshSW, OHM_threshWD, soilstore_id, SoilStoreCap, state_id, SnowUse, SnowFrac, DiagQS, HDD_id, MetForcingData_grid, Ts5mindata_ir, qf, qn1, avkdn, avu1, temp_c, zenith_deg, avrh, press_hpa, ldown, bldgh, alb, emis, cpAnOHM, kkAnOHM, chAnOHM, EmissionsMethod, Tair_av, qn1_av_prev, dqndt_prev, qn1_s_av_prev, dqnsdt_prev, StoreDrainPrm, qn1_S, dataOutLineESTM, qs, qn1_av_next, dqndt_next, qn1_s_av_next, dqnsdt_next, deltaQi, a1, a2, a3)
real(kind(1d0)) function cal_tsfc(qh, avdens, avcp, RA, temp_c)
subroutine suews_cal_roughnessparameters(RoughLenMomMethod, sfr, bldgH, EveTreeH, DecTreeH, porosity_id, FAIBldg, FAIEveTree, FAIDecTree, z0m_in, zdm_in, Z, planF, Zh, z0m, zdm, ZZD)
real(kind(1d0)) function stab_psi_mom(StabilityMethod, ZL)
subroutine aerodynamicresistance(ZZD, z0m, AVU1, L_mod, UStar, VegFraction, AerodynamicResistanceMethod, StabilityMethod, RoughLenHeatMethod, RA)
integer, parameter decidsurf
integer, parameter pavsurf
real(kind(1d0)) avdens
real(kind(1d0)) function cal_tair_av(tair_av_prev, dt_since_start, tstep, temp_c)
subroutine errorhint(errh, ProblemFile, VALUE, value2, valueI)
subroutine suews_cal_surf(sfr, vegfraction, ImpervFraction, PervFraction, NonWaterFraction)
subroutine suews_cal_dls(id, startDLS, endDLS, DLS)
integer, parameter bldgsurf
integer, parameter watersurf
real(kind(1d0)) function qsatf(T, PMB)
subroutine drainage(is, state_is, StorCap, DrainEq, DrainCoef1, DrainCoef2, nsh_real, drain_is)
real function square_real(x)
integer, parameter ivdecid