SUEWS API Site
Documentation of SUEWS source code
Functions/Subroutines
suews_ctrl_init.f95 File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine overallruncontrol
 
subroutine readcoeff (FileName, nlines, ncolumns, HeaderFile, Coeff)
 
subroutine numberrows (FileN, SkipHeaderLines)
 
subroutine initializesurfacecharacteristics (Gridiv, rr)
 
subroutine initialstate (GridName, year_int, Gridiv, NumberOfGrids)
 
real(kind(1d0)) function normalizevegchar (VegCol, Gridiv)
 
subroutine nextinitial (GridName, year_int)
 
subroutine suews_initializemetdata (lunit, ReadlinesMetdata_read)
 
subroutine checkinitial
 

Function/Subroutine Documentation

◆ checkinitial()

subroutine checkinitial

Definition at line 2256 of file suews_ctrl_init.f95.

2257 !Check the parameters in InitialConditions file.
2258 !Modified by HCW 04 Mar 2014, changed soilstore_id(is) checks to use names given in InitialConditions
2259 !Added by LJ in 8/2/2013
2260
2261 USE allocatearray
2262 USE data_in
2263 USE defaultnotused
2264 USE initialcond
2265 USE snowmod
2266 USE time
2267
2268 IMPLICIT NONE
2269
2270 !real(kind(1d0)):: pTol !Precision tolerance for range checks
2271
2272 IF (temp_c0 < (temp_c - 10) .OR. temp_c0 > (temp_c + 10)) THEN
2273 CALL errorhint(37, 'Temp_C0 very different to Tair.', temp_c0, temp_c, notusedi)
2274 END IF
2275
2276 !Check more thoroughly if LAI values are OK. Need to treat different hemispheres as well as tropics separately.
2277 IF (lat > 40) THEN
2278 IF ((laiinitialevetr > laimin(conifsurf - 2) + 1 .AND. (id < 60 .OR. id > 330)) .OR. &
2279 (laiinitialevetr < laimax(conifsurf - 2) - 1 .AND. (id > 130 .AND. id < 244))) THEN
2280 CALL errorhint(37, 'Check LAIinitialEveTr in InitialConditions file', laiinitialevetr, laimin(conifsurf - 2), notusedi)
2281 END IF
2282 IF ((laiinitialdectr > laimin(decidsurf - 2) + 1 .AND. (id < 60 .OR. id > 330)) .OR. &
2283 (laiinitialdectr < laimax(decidsurf - 2) - 1 .AND. (id > 130 .AND. id < 244))) THEN
2284 CALL errorhint(37, 'Check LAIinitialDecTr in InitialConditions file', laiinitialdectr, laimin(decidsurf - 2), notusedi)
2285 END IF
2286 IF ((laiinitialgrass > laimin(grasssurf - 2) + 1 .AND. (id < 60 .OR. id > 330)) .OR. &
2287 (laiinitialgrass < laimax(grasssurf - 2) - 1 .AND. (id > 130 .AND. id < 244))) THEN
2288 CALL errorhint(37, 'Check LAIinitialGrass in InitialConditions file', laiinitialgrass, laimin(grasssurf - 2), notusedi)
2289 END IF
2290
2291 ELSEIF (lat < -40) THEN
2292 IF ((laiinitialevetr < laimax(conifsurf - 2) - 1 .AND. (id < 60 .OR. id > 330)) .OR. &
2293 (laiinitialevetr > laimin(conifsurf - 2) + 1 .AND. (id > 130 .AND. id < 244))) THEN
2294 CALL errorhint(37, 'Check LAIinitialEveTr in InitialConditions file', laiinitialevetr, laimax(conifsurf - 2), notusedi)
2295 END IF
2296 IF ((laiinitialdectr > laimax(decidsurf - 2) - 1 .AND. (id < 60 .OR. id > 330)) .OR. &
2297 (laiinitialdectr > laimin(decidsurf - 2) + 1 .AND. (id > 130 .AND. id < 244))) THEN
2298 CALL errorhint(37, 'Check LAIinitialDecTr in InitialConditions file', laiinitialdectr, laimax(decidsurf - 2), notusedi)
2299 END IF
2300 IF ((laiinitialgrass < laimax(grasssurf - 2) - 1 .AND. (id < 60 .OR. id > 330)) .OR. &
2301 (laiinitialgrass > laimin(grasssurf - 2) + 1 .AND. (id > 130 .AND. id < 244))) THEN
2302 CALL errorhint(37, 'Check LAIinitialGrass in InitialConditions file', laiinitialgrass, laimax(grasssurf - 2), notusedi)
2303 END IF
2304
2305 ELSEIF (lat < 10 .AND. lat > -10) THEN
2306
2307 IF (laiinitialevetr < laimax(conifsurf - 2) - 0.5) THEN
2308 CALL errorhint(37, 'Check LAIinitialEveTr in InitialConditions file', laiinitialevetr, laimax(conifsurf - 2), notusedi)
2309 END IF
2310 IF (laiinitialdectr < laimax(decidsurf - 2) - 0.5) THEN
2311 CALL errorhint(37, 'Check LAIinitialDecTr in InitialConditions file', laiinitialdectr, laimax(decidsurf - 2), notusedi)
2312 END IF
2313 IF (laiinitialgrass < laimax(grasssurf - 2) - 0.5) THEN
2314 CALL errorhint(37, 'Check LAIinitialGrass in InitialConditions file', laiinitialgrass, laimax(grasssurf - 2), notusedi)
2315 END IF
2316
2317 END IF
2318
2319 !Soilstore check
2321 CALL errorhint(37, 'InitialCond: Check initial condition of building soil store.', &
2323 END IF
2325 CALL errorhint(37, 'InitialCond: Check initial condition of paved soil store.', &
2327 END IF
2329 CALL errorhint(37, 'InitialCond: Check initial condition of conif soil store.', &
2331 END IF
2333 CALL errorhint(37, 'InitialCond: Check initial condition of deciduous soil store.', &
2335 END IF
2337 CALL errorhint(37, 'InitialCond: Check initial condition of bare soil soil store.', &
2339 END IF
2341 CALL errorhint(37, 'InitialCond: Check initial condition of grass soil store.', &
2343 END IF
2344
2345 !Snow stuff
2346 IF (snowuse == 1) THEN
2348 CALL errorhint(37, 'InitialCond: SnowWaterBldgsState', snowwaterbldgsstate, snowpackbldgs, notusedi)
2349 END IF
2351 CALL errorhint(37, 'InitialCond: SnowWaterPavedState', snowwaterpavedstate, snowpackpaved, notusedi)
2352 END IF
2354 CALL errorhint(37, 'InitialCond: SnowWaterEveTrstate', snowwaterevetrstate, snowpackevetr, notusedi)
2355 END IF
2357 CALL errorhint(37, 'InitialCond: SnowWaterDecTrState', snowwaterdectrstate, snowpackdectr, notusedi)
2358 END IF
2360 CALL errorhint(37, 'InitialCond: SnowWaterGrassState', snowwatergrassstate, snowpackgrass, notusedi)
2361 END IF
2363 CALL errorhint(37, 'InitialCond: SnowWaterGrassUnirState', snowwaterbsoilstate, snowpackbsoil, notusedi)
2364 END IF
2365 END IF
2366
integer, parameter bldgsurf
real(kind(1d0)), dimension(nsurf) soilstorecap_surf
integer, parameter conifsurf
real(kind(1d0)), dimension(nvegsurf) laimax
real(kind(1d0)), dimension(nvegsurf) laimin
integer, parameter grasssurf
integer, parameter bsoilsurf
integer, parameter pavsurf
integer, parameter decidsurf
real(kind(1d0)) lat
real(kind(1d0)) temp_c
real(kind(1d0)) snowwaterbldgsstate
real(kind(1d0)) snowpackbsoil
real(kind(1d0)) snowpackgrass
real(kind(1d0)) snowwatergrassstate
real(kind(1d0)) soilstoredectrstate
real(kind(1d0)) laiinitialgrass
real(kind(1d0)) soilstorepavedstate
real(kind(1d0)) snowpackdectr
real(kind(1d0)) laiinitialevetr
real(kind(1d0)) snowpackevetr
real(kind(1d0)) soilstorebldgsstate
real(kind(1d0)) snowwaterbsoilstate
real(kind(1d0)) soilstoregrassstate
real(kind(1d0)) snowwaterdectrstate
real(kind(1d0)) temp_c0
real(kind(1d0)) snowwaterpavedstate
real(kind(1d0)) laiinitialdectr
real(kind(1d0)) snowpackbldgs
real(kind(1d0)) soilstoreevetrstate
real(kind(1d0)) snowwaterevetrstate
real(kind(1d0)) soilstorebsoilstate
real(kind(1d0)) snowpackpaved
real(kind(1d0)) crwmax
integer id
subroutine errorhint(errh, ProblemFile, VALUE, value2, valueI)

References allocatearray::bldgsurf, allocatearray::bsoilsurf, allocatearray::conifsurf, snowmod::crwmax, allocatearray::decidsurf, errorhint(), allocatearray::grasssurf, time::id, initialcond::laiinitialdectr, initialcond::laiinitialevetr, initialcond::laiinitialgrass, allocatearray::laimax, allocatearray::laimin, data_in::lat, defaultnotused::notusedi, allocatearray::pavsurf, initialcond::snowpackbldgs, initialcond::snowpackbsoil, initialcond::snowpackdectr, initialcond::snowpackevetr, initialcond::snowpackgrass, initialcond::snowpackpaved, data_in::snowuse, initialcond::snowwaterbldgsstate, initialcond::snowwaterbsoilstate, initialcond::snowwaterdectrstate, initialcond::snowwaterevetrstate, initialcond::snowwatergrassstate, initialcond::snowwaterpavedstate, initialcond::soilstorebldgsstate, initialcond::soilstorebsoilstate, allocatearray::soilstorecap_surf, initialcond::soilstoredectrstate, initialcond::soilstoreevetrstate, initialcond::soilstoregrassstate, initialcond::soilstorepavedstate, data_in::temp_c, and initialcond::temp_c0.

Referenced by suews_translate().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ initializesurfacecharacteristics()

subroutine initializesurfacecharacteristics ( integer  Gridiv,
integer  rr 
)

Definition at line 456 of file suews_ctrl_init.f95.

457
458 USE allocatearray
460 USE data_in
462 USE initial
463 USE sues_data
464
465 IMPLICIT NONE
466
467 INTEGER :: Gridiv, & !Row of SurfaceChar where input information will be stored
468 rr !Row of SiteSelect that matches current grid and year
469 INTEGER :: iii, &
470 ii
471
472 !-------------------------------------------------------------------------------------------
473
474 ! Initialise row of SurfaceChar
475 surfacechar(gridiv, :) = -999
476
477 ! Transfer data in SiteSelect to SurfaceChar
478 surfacechar(gridiv, 1:ncolumnssiteselect) = siteselect(rr, 1:ncolumnssiteselect) !Cols in same order as in SiteSelect.txt
479
480 ! ======== Retrieve information from other input files via codes ========
481
482 ! ---- Find code for Paved surface (Impervious) ----
484 ! Transfer characteristics to SurfaceChar for Paved surface
505 surfacechar(gridiv, c_cpanohm(pavsurf)) = nonveg_coeff(iv5, ci_cpanohm) ! heat capacity, AnOHM TS
506 surfacechar(gridiv, c_kkanohm(pavsurf)) = nonveg_coeff(iv5, ci_kkanohm) ! heat conductivity, AnOHM TS
507 surfacechar(gridiv, c_chanohm(pavsurf)) = nonveg_coeff(iv5, ci_chanohm) ! bulk transfer coef., AnOHM TS
508
509 ! Use SoilCode for Paved to find code for soil characteristics
510 CALL codematchsoil(gridiv, c_soiltcode(pavsurf))
511 ! Transfer soil characteristics to SurfaceChar
520
521 ! Get OHM characteristics for Paved
522 CALL codematchohm(gridiv, pavsurf, 'SWet') !Summer wet
523 ! Transfer OHM characteristics to SurfaceChar
527 CALL codematchohm(gridiv, pavsurf, 'SDry') !Summer dry
528 ! Transfer OHM characteristics to SurfaceChar
532 CALL codematchohm(gridiv, pavsurf, 'WWet') !Winter wet
533 ! Transfer OHM characteristics to SurfaceChar
537 CALL codematchohm(gridiv, pavsurf, 'WDry') !Winter dry
538 ! Transfer OHM characteristics to SurfaceChar
542
543 ! Get water distribution (within grid) for Paved
545 ! Transfer distribution to SurfaceChar
555
556 ! ---- Find code for Bldgs surface (Impervious) ----
558 ! Transfer characteristics to SurfaceChar for Bldgs surface
579 surfacechar(gridiv, c_cpanohm(bldgsurf)) = nonveg_coeff(iv5, ci_cpanohm) ! heat capacity, AnOHM TS
580 surfacechar(gridiv, c_kkanohm(bldgsurf)) = nonveg_coeff(iv5, ci_kkanohm) ! heat conductivity, AnOHM TS
581 surfacechar(gridiv, c_chanohm(bldgsurf)) = nonveg_coeff(iv5, ci_chanohm) ! bulk transfer coef., AnOHM TS
582
583 ! Use SoilCode for Bldgs to find code for soil characteristics
584 CALL codematchsoil(gridiv, c_soiltcode(bldgsurf))
585 ! Transfer soil characteristics to SurfaceChar
594 !Get OHM characteristics for Bldgs
595 CALL codematchohm(gridiv, bldgsurf, 'SWet') !Summer wet
596 ! Transfer OHM characteristics to SurfaceChar
600 CALL codematchohm(gridiv, bldgsurf, 'SDry') !Summer dry
601 ! Transfer OHM characteristics to SurfaceChar
605 CALL codematchohm(gridiv, bldgsurf, 'WWet') !Winter wet
606 ! Transfer OHM characteristics to SurfaceChar
610 CALL codematchohm(gridiv, bldgsurf, 'WDry') !Winter dry
611 ! Transfer OHM characteristics to SurfaceChar
615
616 ! Get water distribution (within grid) for Bldgs
618 ! Transfer distribution to SurfaceChar
628
629 ! ---- Find code for EveTr surface (Pervious) ----
630 CALL codematchveg(rr, c_evetrcode)
631 ! Transfer characteristics to SurfaceChar for EveTr surface
632 ! All surfaces (1-nsurf)
645 ! Veg surfaces only (1-nvegsurf)
660 ! OHM codes
667 ! ESTM code
669 ! AnOHM TS
670 surfacechar(gridiv, c_cpanohm(conifsurf)) = veg_coeff(iv5, cp_cpanohm) ! heat capacity, AnOHM TS
671 surfacechar(gridiv, c_kkanohm(conifsurf)) = veg_coeff(iv5, cp_kkanohm) ! heat conductivity, AnOHM TS
672 surfacechar(gridiv, c_chanohm(conifsurf)) = veg_coeff(iv5, cp_chanohm) ! bulk transfer coef., AnOHM TS
673
675
676 ! ---- Find code for Biogenic CO2 Method ----
678 ! Transfer Biogenic CO2 characteristics to SurfaceChar
687
688 ! Use SoilCode for EveTr to find code for soil characteristics
689 CALL codematchsoil(gridiv, c_soiltcode(conifsurf))
690 ! Transfer soil characteristics to SurfaceChar
699 !Get OHM characteristics for Conif
700 CALL codematchohm(gridiv, conifsurf, 'SWet') !Summer wet
701 ! Transfer OHM characteristics to SurfaceChar
705 CALL codematchohm(gridiv, conifsurf, 'SDry') !Summer dry
706 ! Transfer OHM characteristics to SurfaceChar
710 CALL codematchohm(gridiv, conifsurf, 'WWet') !Winter wet
711 ! Transfer OHM characteristics to SurfaceChar
715 CALL codematchohm(gridiv, conifsurf, 'WDry') !Winter dry
716 ! Transfer OHM characteristics to SurfaceChar
720
721 ! Get water distribution (within grid) for EveTr
723 ! Transfer distribution to SurfaceChar
733
734 ! ---- Find code for DecTr surface (Pervious) ----
735 CALL codematchveg(rr, c_dectrcode)
736 ! Transfer characteristics to SurfaceChar for DecTr surface
737 ! All surfaces (1-nsurf)
750 ! Veg surfaces only (1-nvegsurf)
765 ! OHM codes
772 ! ESTM code
774 ! AnOHM TS
775 surfacechar(gridiv, c_cpanohm(decidsurf)) = veg_coeff(iv5, cp_cpanohm) ! heat capacity, AnOHM TS
776 surfacechar(gridiv, c_kkanohm(decidsurf)) = veg_coeff(iv5, cp_kkanohm) ! heat conductivity, AnOHM TS
777 surfacechar(gridiv, c_chanohm(decidsurf)) = veg_coeff(iv5, cp_chanohm) ! bulk transfer coef., AnOHM TS
778
780
781 ! ---- Find code for Biogenic CO2 Method ----
783 ! Transfer Biogenic CO2 characteristics to SurfaceChar
792
793 ! Use SoilCode for DecTr to find code for soil characteristics
794 CALL codematchsoil(gridiv, c_soiltcode(decidsurf))
795 ! Transfer soil characteristics to SurfaceChar
804 !Get OHM characteristics for Decid
805 CALL codematchohm(gridiv, decidsurf, 'SWet') !Summer wet
806 ! Transfer OHM characteristics to SurfaceChar
810 CALL codematchohm(gridiv, decidsurf, 'SDry') !Summer dry
811 ! Transfer OHM characteristics to SurfaceChar
815 CALL codematchohm(gridiv, decidsurf, 'WWet') !Winter wet
816 ! Transfer OHM characteristics to SurfaceChar
820 CALL codematchohm(gridiv, decidsurf, 'WDry') !Winter dry
821 ! Transfer OHM characteristics to SurfaceChar
825
826 ! Get water distribution (within grid) for DecTr
828 ! Transfer distribution to SurfaceChar
838
839 ! ---- Find code for Grass surface (Pervious) ----
840 CALL codematchveg(rr, c_grasscode)
841 ! Transfer characteristics to SurfaceChar for Grass surface
842 ! All surfaces (1-nsurf)
855 ! Veg surfaces only (1-nvegsurf)
870 ! OHM codes
877 ! ESTM code
879 ! AnOHM TS
880 surfacechar(gridiv, c_cpanohm(grasssurf)) = veg_coeff(iv5, cp_cpanohm) ! heat capacity, AnOHM TS
881 surfacechar(gridiv, c_kkanohm(grasssurf)) = veg_coeff(iv5, cp_kkanohm) ! heat conductivity, AnOHM TS
882 surfacechar(gridiv, c_chanohm(grasssurf)) = veg_coeff(iv5, cp_chanohm) ! bulk transfer coef., AnOHM TS
883
885
886 ! ---- Find code for Biogenic CO2 Method ----
888 ! Transfer Biogenic CO2 characteristics to SurfaceChar
897
898 ! Use SoilCode for Grass to find code for soil characteristics
899 CALL codematchsoil(gridiv, c_soiltcode(grasssurf))
900 ! Transfer soil characteristics to SurfaceChar
909 !Get OHM characteristics for Grass
910 CALL codematchohm(gridiv, grasssurf, 'SWet') !Summer wet
911 ! Transfer OHM characteristics to SurfaceChar
915 CALL codematchohm(gridiv, grasssurf, 'SDry') !Summer dry
916 ! Transfer OHM characteristics to SurfaceChar
920 CALL codematchohm(gridiv, grasssurf, 'WWet') !Winter wet
921 ! Transfer OHM characteristics to SurfaceChar
925 CALL codematchohm(gridiv, grasssurf, 'WDry') !Winter dry
926 ! Transfer OHM characteristics to SurfaceChar
930
931 ! Get water distribution (within grid) for Grass
933 ! Transfer distribution to SurfaceChar
943
944 ! ---- Find code for BSoil surface (Impervious) ----
946 ! Transfer characteristics to SurfaceChar for BSoil surface
947 ! All surfaces (1-nsurf)
968 surfacechar(gridiv, c_cpanohm(bsoilsurf)) = nonveg_coeff(iv5, ci_cpanohm) ! heat capacity, AnOHM TS
969 surfacechar(gridiv, c_kkanohm(bsoilsurf)) = nonveg_coeff(iv5, ci_kkanohm) ! heat conductivity, AnOHM TS
970 surfacechar(gridiv, c_chanohm(bsoilsurf)) = nonveg_coeff(iv5, ci_chanohm) ! bulk transfer coef., AnOHM TS
971
972 ! Use SoilCode for BSoil to find code for soil characteristics
973 CALL codematchsoil(gridiv, c_soiltcode(bsoilsurf))
974 ! Transfer soil characteristics to SurfaceChar
983
984 ! Get OHM characteristics for BSoil
985 CALL codematchohm(gridiv, bsoilsurf, 'SWet') !Summer wet
986 ! Transfer OHM characteristics to SurfaceChar
990 CALL codematchohm(gridiv, bsoilsurf, 'SDry') !Summer dry
991 ! Transfer OHM characteristics to SurfaceChar
995 CALL codematchohm(gridiv, bsoilsurf, 'WWet') !Winter wet
996 ! Transfer OHM characteristics to SurfaceChar
1000 CALL codematchohm(gridiv, bsoilsurf, 'WDry') !Winter dry
1001 ! Transfer OHM characteristics to SurfaceChar
1005
1006 ! Get water distribution (within grid) for Bare soil
1008 ! Transfer distribution to SurfaceChar
1018
1019 ! ---- Find code for Water surface (Water) ----
1020 CALL codematchwater(rr, c_watercode)
1021 ! Transfer characteristics to SurfaceChar for Water surface
1022 ! All surfaces (1-nsurf)
1033 ! Water surface only
1035 ! OHM codes
1042 ! ESTM code
1044 ! AnOHM TS
1045 surfacechar(gridiv, c_cpanohm(watersurf)) = water_coeff(iv5, cw_cpanohm) ! heat capacity, AnOHM TS
1046 surfacechar(gridiv, c_kkanohm(watersurf)) = water_coeff(iv5, cw_kkanohm) ! heat conductivity, AnOHM TS
1047 surfacechar(gridiv, c_chanohm(watersurf)) = water_coeff(iv5, cw_chanohm) ! bulk transfer coef., AnOHM TS
1048 ! Get OHM characteristics for Water
1049 CALL codematchohm(gridiv, watersurf, 'SWet') !Summer wet
1050 ! Transfer OHM characteristics to SurfaceChar
1054 CALL codematchohm(gridiv, watersurf, 'SDry') !Summer dry
1055 ! Transfer OHM characteristics to SurfaceChar
1059 CALL codematchohm(gridiv, watersurf, 'WWet') !Winter wet
1060 ! Transfer OHM characteristics to SurfaceChar
1064 CALL codematchohm(gridiv, watersurf, 'WDry') !Winter dry
1065 ! Transfer OHM characteristics to SurfaceChar
1069
1070 ! Get water distribution (within grid) for Water
1072 ! Transfer distribution to SurfaceChar
1082
1083 ! ---- Find code for Snow surface (Snow) ----
1084 CALL codematchsnow(rr, c_snowcode)
1085 ! Transfer characteristics to SurfaceChar for Snow surface
1090 !SurfaceChar(gridiv,c_SnowAlb) = Snow_Coeff(iv5,cs_SnowAlb)
1107
1108 ! ESTM code
1110 ! SurfaceChar(Gridiv,c_CpAnOHM(nsurf+1)) = Snow_Coeff(iv5,cs_CpAnOHM) ! heat capacity, AnOHM TS
1111 ! SurfaceChar(Gridiv,c_KkAnOHM(nsurf+1)) = Snow_Coeff(iv5,cs_KkAnOHM) ! heat conductivity, AnOHM TS
1112 ! SurfaceChar(Gridiv,c_ChAnOHM(nsurf+1)) = Snow_Coeff(iv5,cs_ChAnOHM) ! bulk transfer coef., AnOHM TS
1113 ! Get OHM characteristics for Snow
1114 CALL codematchohm(gridiv, (nsurf + 1), 'SWet') !Summer wet
1115 ! Transfer OHM characteristics to SurfaceChar
1119 CALL codematchohm(gridiv, (nsurf + 1), 'SDry') !Summer dry
1120 ! Transfer OHM characteristics to SurfaceChar
1124 CALL codematchohm(gridiv, (nsurf + 1), 'WWet') !Winter wet
1125 ! Transfer OHM characteristics to SurfaceChar
1129 CALL codematchohm(gridiv, (nsurf + 1), 'WDry') !Winter dry
1130 ! Transfer OHM characteristics to SurfaceChar
1134
1135 !Transfer ESTM characteristics to SurfaceChar
1136 DO iii = 1, (nsurf + 1)
1137 IF (surfacechar(gridiv, c_estmcode(iii)) /= 0) THEN !If ESTM Code not equal to zero, use code as normal
1138 CALL codematchestm(gridiv, iii)
1154 !Extra characteristics for Bldg surfaces
1155 IF (iii == bldgsurf) THEN
1192 END IF
1193 !If ESTM Code equals zero, use codes and surface fractions from SiteSelect.txt for Paved and Bldgs
1194 ELSEIF (iii == pavsurf .AND. surfacechar(gridiv, c_estmcode(iii)) == 0) THEN
1195 DO ii = 1, 3 !for the 3x Paved ESTM classes
1196 CALL codematchestm_class(gridiv, iii, ii)
1212 END DO
1213 ELSEIF (iii == bldgsurf .AND. surfacechar(gridiv, c_estmcode(iii)) == 0) THEN
1214 DO ii = 1, 5 !for the 5x Bldgs ESTM classes
1215 CALL codematchestm_class(gridiv, iii, ii)
1231 !Extra characteristics for Bldgs surface
1268 END DO
1269 END IF
1270 END DO
1271
1272 ! ---- Find code for Surface conductances ----
1274 ! Transfer conductance characteristics to SurfaceChar
1287
1288 ! ---- Find code for Anthropogenic heat ----
1290 ! Transfer Anthropogenic heat characteristics to SurfaceChar
1329
1330 ! ---- Find code for Irrigation ----
1332 ! Transfer Irrigation characteristics to SurfaceChar
1342
1343 ! ---- Find code for Hourly Profiles ----
1344 ! Energy use (weekdays)
1345 CALL codematchprof(gridiv, c_enprofwd)
1347 ! Energy use (weekends)
1348 CALL codematchprof(gridiv, c_enprofwe)
1350 ! Water use profile (manual, weekdays)
1351 CALL codematchprof(gridiv, c_wprofmanuwd)
1353 ! Water use profile (manual, weekends)
1354 CALL codematchprof(gridiv, c_wprofmanuwe)
1356 ! Water use profile (automatic, weekdays)
1357 CALL codematchprof(gridiv, c_wprofautowd)
1359 ! Water use profile (automatic, weekends)
1360 CALL codematchprof(gridiv, c_wprofautowe)
1362 ! Snow clearing profile (weekdays)
1363 CALL codematchprof(gridiv, c_snowprofwd)
1365 ! Snow clearing profile (weekends)
1366 CALL codematchprof(gridiv, c_snowprofwe)
1368 !Human activity (weekdays)
1369 CALL codematchprof(gridiv, c_co2mwd)
1371 !Human activity (weekends)
1372 CALL codematchprof(gridiv, c_co2mwe)
1374 !Traffic (weekdays)
1375 CALL codematchprof(gridiv, c_traffprofwd)
1377 !Traffic (weekends)
1378 CALL codematchprof(gridiv, c_traffprofwe)
1380 !Population (weekdays)
1381 CALL codematchprof(gridiv, c_popprofwd)
1383 !Population (weekends)
1384 CALL codematchprof(gridiv, c_popprofwe)
1386
1387 ! TS 05 Jul 2018: No longer needed as interpolation is done through specific subroutines at each required instant
1388 ! the below is commented out by TS 05 Jul 2018
1389 ! ! ---- Interpolate Hourly Profiles to model timestep and normalise
1390 ! TstepProfiles(Gridiv,:,:) = -999 !Initialise TstepProfiles
1391 ! ! Energy use
1392 ! CALL SUEWS_InterpHourlyProfiles(Gridiv,cTP_EnUseWD,c_HrProfEnUseWD)
1393 ! CALL SUEWS_InterpHourlyProfiles(Gridiv,cTP_EnUseWE,c_HrProfEnUseWE)
1394 !
1395 ! ! For energy use, normalise so the AVERAGE of the multipliers is equal to 1
1396 ! TstepProfiles(Gridiv,cTP_EnUseWD,:) = TstepProfiles(Gridiv,cTP_EnUseWD,:) / SUM(TstepProfiles(Gridiv,cTP_EnUseWD,:))*24*nsh_real
1397 ! TstepProfiles(Gridiv,cTP_EnUseWE,:) = TstepProfiles(Gridiv,cTP_EnUseWE,:) / SUM(TstepProfiles(Gridiv,cTP_EnUseWE,:))*24*nsh_real
1398 !
1399 ! ! Water use
1400 ! CALL SUEWS_InterpHourlyProfiles(Gridiv,cTP_WUManuWD,c_HrProfWUManuWD)
1401 ! CALL SUEWS_InterpHourlyProfiles(Gridiv,cTP_WUManuWE,c_HrProfWUManuWE)
1402 ! CALL SUEWS_InterpHourlyProfiles(Gridiv,cTP_WUAutoWD,c_HrProfWUAutoWD)
1403 ! CALL SUEWS_InterpHourlyProfiles(Gridiv,cTP_WUAutoWE,c_HrProfWUAutoWE)
1404 ! ! For water use, normalise so the SUM of the multipliers is equal to 1 (profile is multiplied by daily water use)
1405 ! TstepProfiles(Gridiv,cTP_WUManuWD,:) = TstepProfiles(Gridiv,cTP_WUManuWD,:) / SUM(TstepProfiles(Gridiv,cTP_WUManuWD,:))
1406 ! TstepProfiles(Gridiv,cTP_WUManuWE,:) = TstepProfiles(Gridiv,cTP_WUManuWE,:) / SUM(TstepProfiles(Gridiv,cTP_WUManuWE,:))
1407 ! TstepProfiles(Gridiv,cTP_WUAutoWD,:) = TstepProfiles(Gridiv,cTP_WUAutoWD,:) / SUM(TstepProfiles(Gridiv,cTP_WUAutoWD,:))
1408 ! TstepProfiles(Gridiv,cTP_WUAutoWE,:) = TstepProfiles(Gridiv,cTP_WUAutoWE,:) / SUM(TstepProfiles(Gridiv,cTP_WUAutoWE,:))
1409 !
1410 ! ! Human activity for CO2 calculations
1411 ! CALL SUEWS_InterpHourlyProfiles(Gridiv,cTP_HumActivityWD,c_HrProfHumActivityWD)
1412 ! CALL SUEWS_InterpHourlyProfiles(Gridiv,cTP_HumActivityWE,c_HrProfHumActivityWE)
1413 !
1414 ! ! For human activity, check values are between 1 (night) and 2 (day)
1415 ! IF(ANY(TstepProfiles(Gridiv,cTP_HumActivityWD,:) < 1 .OR. TstepProfiles(Gridiv,cTP_HumActivityWD,:) > 2)) THEN
1416 ! CALL ErrorHint(70,'Profile value for human activity (WD) exceeds allowed range 1-2.',NotUsed,NotUsed,notUsedI)
1417 ! ENDIF
1418 ! IF(ANY(TstepProfiles(Gridiv,cTP_HumActivityWE,:) < 1 .OR. TstepProfiles(Gridiv,cTP_HumActivityWE,:) > 2)) THEN
1419 ! CALL ErrorHint(70,'Profile value for human activity (WE) exceeds allowed range 1-2.',NotUsed,NotUsed,notUsedI)
1420 ! ENDIF
1421 !
1422 ! CALL SUEWS_InterpHourlyProfiles(Gridiv,cTP_TraffProfWD,c_HrProfTraffWD)
1423 ! CALL SUEWS_InterpHourlyProfiles(Gridiv,cTP_TraffProfWE,c_HrProfTraffWE)
1424 ! ! For traffic, normalise so the AVERAGE of the multipliers is equal to 1
1425 ! TstepProfiles(Gridiv,cTP_TraffProfWD,:) = TstepProfiles(Gridiv,cTP_TraffProfWD,:) &
1426 ! / SUM(TstepProfiles(Gridiv,cTP_TraffProfWD,:))*24*nsh_real
1427 ! TstepProfiles(Gridiv,cTP_TraffProfWE,:) = TstepProfiles(Gridiv,cTP_TraffProfWE,:) &
1428 ! / SUM(TstepProfiles(Gridiv,cTP_TraffProfWE,:))*24*nsh_real
1429 !
1430 ! ! Population for CO2 calculations
1431 ! CALL SUEWS_InterpHourlyProfiles(Gridiv,cTP_PopProfWD,c_HrProfPopWD)
1432 ! CALL SUEWS_InterpHourlyProfiles(Gridiv,cTP_PopProfWE,c_HrProfPopWE)
1433
integer, dimension(5) c_wall_thick4_bldgs
integer, dimension(nvegsurf) c_gddfull
real(kind(1d0)), dimension(:, :), allocatable biogen_coeff
integer, dimension(5) c_internal_rhocp4_bldgs
integer, dimension(nsurfincsnow) c_surf_thick5
real(kind(1d0)), dimension(:, :), allocatable estmcoefficients_coeff
integer, dimension(nsurf) c_kkanohm
integer, dimension(5) c_surf_k3_bldgs
integer, parameter ncolumnssiteselect
integer, dimension(5) c_wall_thick1_bldgs
real(kind(1d0)), dimension(:, :), allocatable snow_coeff
integer, dimension(nsurfincsnow) c_surf_rhocp2
integer, dimension(24) c_hrproftraffwd
real(kind(1d0)), dimension(:, :), allocatable siteselect
real(kind(1d0)), dimension(:, :), allocatable anthropogenic_coeff
integer, dimension(nsurf) c_cpanohm
integer, dimension(nsurfincsnow) c_ohmcode_wwet
integer, dimension(5) c_ch_iwall_bldgs
integer, dimension(5) c_internal_rhocp1_bldgs
real(kind(1d0)), dimension(:, :), allocatable wgwaterdist_coeff
integer, dimension(24) c_hrprofsnowcwe
integer, dimension(5) c_wall_k1_bldgs
integer, dimension(nsurfincsnow) c_ohmthresh_sw
integer, dimension(nsurfincsnow) c_surf_k1
integer, dimension(3) c_ie_a
integer c_tcriticcooling_we
integer, dimension(5) c_internal_thick3_bldgs
integer, dimension(nsurfincsnow) c_surf_thick4
integer, dimension(5) c_internal_k3_bldgs
integer, dimension(3) c_surf_rhocp2_paved
integer, dimension(nvegsurf) c_sddfull
integer, dimension(nvegsurf) c_beta_bioco2
integer, dimension(nsurf) c_wgtorunoff
integer, dimension(nsurfincsnow) c_ohmcode_sdry
integer, dimension(nsurf) c_wgtograss
integer, dimension(nsurfincsnow) c_a2_swet
integer, dimension(5) c_surf_k1_bldgs
integer, dimension(nsurf) c_wgtobldgs
integer, dimension(nsurfincsnow) c_a1_wwet
integer, dimension(nvegsurf) c_gsmax
integer, dimension(nsurf) c_soilstcap
integer, dimension(nsurf) c_soildens
integer, dimension(nsurfincsnow) c_a1_swet
integer, dimension(nsurfincsnow) c_ohmcode_swet
integer, dimension(nsurf) c_soildepth
integer, dimension(5) c_surf_thick3_bldgs
integer, dimension(nsurf) c_wgtoevetr
integer, dimension(5) c_surf_rhocp2_bldgs
integer c_ahslopecooling_we
integer, dimension(24) c_hrprofenusewe
integer, dimension(7) c_daywatper
integer, dimension(5) c_surf_rhocp5_bldgs
integer, dimension(nsurf) c_soilinfrate
integer, dimension(5) c_ch_ibld_bldgs
integer, dimension(5) c_surf_rhocp4_bldgs
integer, dimension(5) c_internal_thick1_bldgs
integer, dimension(24) c_hrprofwuautowd
integer, dimension(nvegsurf) c_laimax
integer, dimension(nsurfincsnow) c_a2_sdry
integer, dimension(nsurf) c_statelimit
integer, dimension(3) c_ie_m
integer c_tcriticheating_we
integer, dimension(nvegsurf) c_alpha_bioco2
integer, dimension(3) c_surf_k4_paved
integer, dimension(nsurf) c_drcoef2
integer, dimension(5) c_ch_iroof_bldgs
integer, dimension(nsurf) c_ksat
integer c_ahslopecooling_wd
integer, dimension(nsurf) c_obssmmax
integer, dimension(nvegsurf) c_leafgp2
integer, dimension(nsurf) c_chanohm
integer, dimension(5) c_surf_thick2_bldgs
integer, dimension(5) c_internal_rhocp2_bldgs
integer, dimension(5) c_wall_k3_bldgs
integer, dimension(5) c_wall_k5_bldgs
integer, dimension(nsurf) c_snowlimpat
integer, dimension(nsurfincsnow) c_ohmthresh_wd
integer, dimension(nsurfincsnow) c_a3_wwet
integer, dimension(5) c_surf_rhocp1_bldgs
integer, dimension(3) c_surf_k5_paved
integer, dimension(nsurf) c_stormin
integer c_tcriticheating_wd
integer, dimension(7) c_daywat
integer, dimension(5) c_wall_rhocp1_bldgs
integer, dimension(24) c_hrprofsnowcwd
integer, dimension(5) c_wall_rhocp3_bldgs
integer, dimension(nsurfincsnow) c_a3_sdry
integer c_ahslopeheating_we
integer, dimension(5) c_surf_k2_bldgs
integer c_tcriticcooling_wd
integer, dimension(nsurfincsnow) c_surf_k5
integer, dimension(nsurf) c_wgtowater
real(kind(1d0)), dimension(:, :), allocatable water_coeff
integer, dimension(24) c_hrprofpopwe
integer, dimension(nsurfincsnow) c_surf_thick3
integer, dimension(nsurfincsnow) c_ohmcode_wdry
integer, dimension(3) c_surf_k1_paved
real(kind(1d0)), dimension(:, :), allocatable nonveg_coeff
integer, dimension(5) c_surf_k5_bldgs
integer, dimension(nsurf) c_wgtosoilstore
integer, dimension(nvegsurf) c_laimin
integer, dimension(3) c_surf_rhocp4_paved
integer, dimension(nvegsurf) c_alpha_enh_bioco2
real(kind(1d0)), dimension(:, :), allocatable conductance_coeff
integer, dimension(24) c_hrprofpopwd
integer, dimension(nsurfincsnow) c_surf_thick1
integer, dimension(nsurf) c_wgtopaved
integer, dimension(5) c_internal_thick4_bldgs
integer, dimension(nsurfincsnow) c_a3_wdry
integer, dimension(24) c_hrprofwumanuwe
integer, dimension(nvegsurf) c_min_res_bioco2
integer, dimension(nvegsurf) c_resp_b
integer, dimension(5) c_internal_rhocp5_bldgs
integer, dimension(nsurf) c_albmin
integer, dimension(5) c_internal_k5_bldgs
integer c_frfossilfuel_heat
integer, dimension(nsurf) c_obssmdepth
integer, dimension(24) c_hrproftraffwe
integer, dimension(5) c_wall_rhocp5_bldgs
integer, dimension(24) c_hrprofwumanuwd
integer, dimension(5) c_surf_thick5_bldgs
integer, dimension(5) c_wall_k2_bldgs
integer, dimension(nsurfincsnow) c_a3_swet
real(kind(1d0)), dimension(:, :), allocatable surfacechar
integer c_ahslopeheating_wd
integer, parameter ivgrass
integer, dimension(3) c_surf_thick4_paved
integer, dimension(nsurf) c_stormax
integer, dimension(nsurfincsnow) c_surf_rhocp3
integer, dimension(nsurf) c_snowlimrem
integer, dimension(3) c_surf_rhocp1_paved
integer c_frfossilfuel_nonheat
integer, dimension(nvegsurf) c_baset
integer, dimension(nsurf) c_emis
integer, dimension(nvegsurf) c_resp_a
integer, dimension(nsurfincsnow) c_surf_thick2
integer, dimension(5) c_nroom_bldgs
integer, dimension(5) c_internal_k2_bldgs
integer, dimension(5) c_internal_k4_bldgs
integer, dimension(5) c_surf_k4_bldgs
integer, dimension(24) c_hrprofhumactivitywe
integer, parameter watersurf
real(kind(1d0)), dimension(:, :), allocatable profiles_coeff
integer, dimension(3) c_surf_thick1_paved
integer, dimension(nsurfincsnow) c_estmcode
real(kind(1d0)), dimension(:, :), allocatable soil_coeff
integer, dimension(nsurfincsnow) c_surf_k3
integer, dimension(nvegsurf) c_biogenco2code
integer, dimension(3) c_surf_rhocp5_paved
integer, parameter ivdecid
integer, dimension(5) c_surf_thick4_bldgs
integer, dimension(3) c_surf_thick2_paved
integer, dimension(5) c_surf_thick1_bldgs
integer, dimension(nsurfincsnow) c_surf_k2
integer, dimension(3) c_surf_rhocp3_paved
integer, dimension(nsurfincsnow) c_surf_rhocp5
integer, dimension(nsurf) c_wgtodectr
integer, dimension(5) c_internal_rhocp3_bldgs
integer, dimension(nsurf) c_soiltcode
real(kind(1d0)), dimension(:, :), allocatable ohmcoefficients_coeff
integer, dimension(nvegsurf) c_beta_enh_bioco2
integer, dimension(nvegsurf) c_porositymax
integer, parameter nsurf
integer, dimension(5) c_em_ibld_bldgs
integer, dimension(3) c_surf_k3_paved
integer, dimension(nvegsurf) c_theta_bioco2
integer, dimension(nvegsurf) c_porositymin
integer, dimension(nvegsurf) c_leafgp1
integer, dimension(3) c_surf_thick5_paved
integer, dimension(nsurfincsnow) c_surf_k4
real(kind(1d0)), dimension(:, :), allocatable irrigation_coeff
integer, dimension(nsurf) c_obssnrfrac
integer, dimension(5) c_wall_thick5_bldgs
integer, dimension(5) c_alb_ibld_bldgs
integer, dimension(24) c_hrprofwuautowe
integer, parameter ivconif
integer, dimension(3) c_surf_thick3_paved
integer, dimension(nsurfincsnow) c_surf_rhocp1
integer, dimension(24) c_hrprofhumactivitywd
integer, dimension(5) c_wall_rhocp4_bldgs
integer, dimension(nsurf) c_wgtobsoil
integer, dimension(nsurfincsnow) c_a1_wdry
integer, dimension(nsurfincsnow) c_a2_wdry
integer, dimension(nvegsurf) c_laieq
real(kind(1d0)), dimension(:, :), allocatable veg_coeff
integer, dimension(5) c_wall_k4_bldgs
integer, dimension(nsurfincsnow) c_a1_sdry
integer, dimension(nsurfincsnow) c_surf_rhocp4
integer, dimension(nvegsurf) c_leafop1
integer, dimension(5) c_surf_rhocp3_bldgs
integer, dimension(nsurf) c_wetthresh
integer, dimension(5) c_internal_thick2_bldgs
integer, dimension(5) c_internal_k1_bldgs
integer, dimension(nsurfincsnow) c_a2_wwet
integer, dimension(nvegsurf) c_leafop2
integer, dimension(5) c_internal_thick5_bldgs
integer, dimension(24) c_hrprofenusewd
integer, dimension(nsurf) c_dreq
integer, dimension(nsurf) c_albmax
integer, dimension(5) c_wall_rhocp2_bldgs
integer, dimension(3) c_surf_k2_paved
integer, dimension(5) c_wall_thick2_bldgs
integer, dimension(nsurf) c_drcoef1
integer, dimension(nvegsurf) c_basete
integer, dimension(5) c_wall_thick3_bldgs
integer, dimension(24) cpr_hours
subroutine codematchanthropogenic(rr, CodeCol)
subroutine codematchirrigation(rr, CodeCol)
subroutine codematchohm(Gridiv, is, SWWD)
subroutine codematchsoil(Gridiv, SurfaceCharCodeCol)
subroutine codematchwater(rr, CodeCol)
subroutine codematchnonveg(rr, CodeCol)
subroutine codematchdist(rr, CodeCol, codeColSameSurf)
subroutine codematchprof(Gridiv, SurfaceCharCodeCol)
subroutine codematchestm_class(Gridiv, is, ii)
subroutine codematchestm(Gridiv, is)
subroutine codematchbiogen(Gridiv, SurfaceCharCodeCol)
subroutine codematchsnow(rr, CodeCol)
subroutine codematchveg(rr, CodeCol)
subroutine codematchconductance(rr, CodeCol)

References allocatearray::anthropogenic_coeff, allocatearray::biogen_coeff, allocatearray::bldgsurf, allocatearray::bsoilsurf, allocatearray::c_a1_sdry, allocatearray::c_a1_swet, allocatearray::c_a1_wdry, allocatearray::c_a1_wwet, allocatearray::c_a2_sdry, allocatearray::c_a2_swet, allocatearray::c_a2_wdry, allocatearray::c_a2_wwet, allocatearray::c_a3_sdry, allocatearray::c_a3_swet, allocatearray::c_a3_wdry, allocatearray::c_a3_wwet, allocatearray::c_ahmin_wd, allocatearray::c_ahmin_we, allocatearray::c_ahslopecooling_wd, allocatearray::c_ahslopecooling_we, allocatearray::c_ahslopeheating_wd, allocatearray::c_ahslopeheating_we, allocatearray::c_alb_ibld, allocatearray::c_alb_ibld_bldgs, allocatearray::c_albmax, allocatearray::c_albmin, allocatearray::c_alpha_bioco2, allocatearray::c_alpha_enh_bioco2, allocatearray::c_baset, allocatearray::c_baset_hc, allocatearray::c_basete, allocatearray::c_beta_bioco2, allocatearray::c_beta_enh_bioco2, allocatearray::c_biogenco2code, colnamesinputfiles::c_bldgscode, colnamesinputfiles::c_bsoilcode, allocatearray::c_ch_ibld, allocatearray::c_ch_ibld_bldgs, allocatearray::c_ch_iroof, allocatearray::c_ch_iroof_bldgs, allocatearray::c_ch_iwall, allocatearray::c_ch_iwall_bldgs, allocatearray::c_chanohm, allocatearray::c_co2mwd, allocatearray::c_co2mwe, allocatearray::c_co2pointsource, colnamesinputfiles::c_condcode, allocatearray::c_cpanohm, allocatearray::c_daywat, allocatearray::c_daywatper, colnamesinputfiles::c_dectrcode, allocatearray::c_drcoef1, allocatearray::c_drcoef2, allocatearray::c_dreq, allocatearray::c_ef_umolco2perj, allocatearray::c_em_ibld, allocatearray::c_em_ibld_bldgs, allocatearray::c_emis, allocatearray::c_enef_v_jkm, allocatearray::c_enprofwd, allocatearray::c_enprofwe, allocatearray::c_estmcode, colnamesinputfiles::c_evetrcode, allocatearray::c_faut, allocatearray::c_fcef_v_kgkmwd, allocatearray::c_fcef_v_kgkmwe, allocatearray::c_frfossilfuel_heat, allocatearray::c_frfossilfuel_nonheat, allocatearray::c_frpddwe, allocatearray::c_gddfull, colnamesinputfiles::c_grasscode, allocatearray::c_gsg1, allocatearray::c_gsg2, allocatearray::c_gsg3, allocatearray::c_gsg4, allocatearray::c_gsg5, allocatearray::c_gsg6, allocatearray::c_gskmax, allocatearray::c_gsmax, allocatearray::c_gsmodel, allocatearray::c_gss1, allocatearray::c_gss2, allocatearray::c_gsth, allocatearray::c_gstl, allocatearray::c_h_maintain, allocatearray::c_hrprofenusewd, allocatearray::c_hrprofenusewe, allocatearray::c_hrprofhumactivitywd, allocatearray::c_hrprofhumactivitywe, allocatearray::c_hrprofpopwd, allocatearray::c_hrprofpopwe, allocatearray::c_hrprofsnowcwd, allocatearray::c_hrprofsnowcwe, allocatearray::c_hrproftraffwd, allocatearray::c_hrproftraffwe, allocatearray::c_hrprofwuautowd, allocatearray::c_hrprofwuautowe, allocatearray::c_hrprofwumanuwd, allocatearray::c_hrprofwumanuwe, allocatearray::c_ie_a, allocatearray::c_ie_m, allocatearray::c_ieend, allocatearray::c_iestart, allocatearray::c_internal_k1, allocatearray::c_internal_k1_bldgs, allocatearray::c_internal_k2, allocatearray::c_internal_k2_bldgs, allocatearray::c_internal_k3, allocatearray::c_internal_k3_bldgs, allocatearray::c_internal_k4, allocatearray::c_internal_k4_bldgs, allocatearray::c_internal_k5, allocatearray::c_internal_k5_bldgs, allocatearray::c_internal_rhocp1, allocatearray::c_internal_rhocp1_bldgs, allocatearray::c_internal_rhocp2, allocatearray::c_internal_rhocp2_bldgs, allocatearray::c_internal_rhocp3, allocatearray::c_internal_rhocp3_bldgs, allocatearray::c_internal_rhocp4, allocatearray::c_internal_rhocp4_bldgs, allocatearray::c_internal_rhocp5, allocatearray::c_internal_rhocp5_bldgs, allocatearray::c_internal_thick1, allocatearray::c_internal_thick1_bldgs, allocatearray::c_internal_thick2, allocatearray::c_internal_thick2_bldgs, allocatearray::c_internal_thick3, allocatearray::c_internal_thick3_bldgs, allocatearray::c_internal_thick4, allocatearray::c_internal_thick4_bldgs, allocatearray::c_internal_thick5, allocatearray::c_internal_thick5_bldgs, allocatearray::c_intwu, colnamesinputfiles::c_irrcode, allocatearray::c_kkanohm, allocatearray::c_ksat, allocatearray::c_laieq, allocatearray::c_laimax, allocatearray::c_laimin, allocatearray::c_leafgp1, allocatearray::c_leafgp2, allocatearray::c_leafop1, allocatearray::c_leafop2, allocatearray::c_maxfcmetab, allocatearray::c_maxqfmetab, allocatearray::c_min_res_bioco2, allocatearray::c_minfcmetab, allocatearray::c_minqfmetab, allocatearray::c_nroom, allocatearray::c_nroom_bldgs, allocatearray::c_obssmdepth, allocatearray::c_obssmmax, allocatearray::c_obssnrfrac, allocatearray::c_ohmcode_sdry, allocatearray::c_ohmcode_swet, allocatearray::c_ohmcode_wdry, allocatearray::c_ohmcode_wwet, allocatearray::c_ohmthresh_sw, allocatearray::c_ohmthresh_wd, colnamesinputfiles::c_pavedcode, allocatearray::c_popprofwd, allocatearray::c_popprofwe, allocatearray::c_porositymax, allocatearray::c_porositymin, allocatearray::c_qf_a1, allocatearray::c_qf_a2, allocatearray::c_qf_b1, allocatearray::c_qf_b2, allocatearray::c_qf_c1, allocatearray::c_qf_c2, colnamesinputfiles::c_qfcode, allocatearray::c_resp_a, allocatearray::c_resp_b, allocatearray::c_sddfull, allocatearray::c_snowalbmax, allocatearray::c_snowalbmin, colnamesinputfiles::c_snowcode, allocatearray::c_snowcrwmax, allocatearray::c_snowcrwmin, allocatearray::c_snowemis, allocatearray::c_snowlimpat, allocatearray::c_snowlimrem, allocatearray::c_snowplimalb, allocatearray::c_snowplimsnow, colnamesinputfiles::c_snowprofwd, colnamesinputfiles::c_snowprofwe, allocatearray::c_snowrmfactor, allocatearray::c_snowsdmax, allocatearray::c_snowsdmin, allocatearray::c_snowtau_a, allocatearray::c_snowtau_f, allocatearray::c_snowtau_r, allocatearray::c_snowtmfactor, allocatearray::c_soildens, allocatearray::c_soildepth, allocatearray::c_soilinfrate, allocatearray::c_soilstcap, allocatearray::c_soiltcode, allocatearray::c_statelimit, allocatearray::c_stormax, allocatearray::c_stormin, allocatearray::c_surf_k1, allocatearray::c_surf_k1_bldgs, allocatearray::c_surf_k1_paved, allocatearray::c_surf_k2, allocatearray::c_surf_k2_bldgs, allocatearray::c_surf_k2_paved, allocatearray::c_surf_k3, allocatearray::c_surf_k3_bldgs, allocatearray::c_surf_k3_paved, allocatearray::c_surf_k4, allocatearray::c_surf_k4_bldgs, allocatearray::c_surf_k4_paved, allocatearray::c_surf_k5, allocatearray::c_surf_k5_bldgs, allocatearray::c_surf_k5_paved, allocatearray::c_surf_rhocp1, allocatearray::c_surf_rhocp1_bldgs, allocatearray::c_surf_rhocp1_paved, allocatearray::c_surf_rhocp2, allocatearray::c_surf_rhocp2_bldgs, allocatearray::c_surf_rhocp2_paved, allocatearray::c_surf_rhocp3, allocatearray::c_surf_rhocp3_bldgs, allocatearray::c_surf_rhocp3_paved, allocatearray::c_surf_rhocp4, allocatearray::c_surf_rhocp4_bldgs, allocatearray::c_surf_rhocp4_paved, allocatearray::c_surf_rhocp5, allocatearray::c_surf_rhocp5_bldgs, allocatearray::c_surf_rhocp5_paved, allocatearray::c_surf_thick1, allocatearray::c_surf_thick1_bldgs, allocatearray::c_surf_thick1_paved, allocatearray::c_surf_thick2, allocatearray::c_surf_thick2_bldgs, allocatearray::c_surf_thick2_paved, allocatearray::c_surf_thick3, allocatearray::c_surf_thick3_bldgs, allocatearray::c_surf_thick3_paved, allocatearray::c_surf_thick4, allocatearray::c_surf_thick4_bldgs, allocatearray::c_surf_thick4_paved, allocatearray::c_surf_thick5, allocatearray::c_surf_thick5_bldgs, allocatearray::c_surf_thick5_paved, allocatearray::c_tcriticcooling_wd, allocatearray::c_tcriticcooling_we, allocatearray::c_tcriticheating_wd, allocatearray::c_tcriticheating_we, allocatearray::c_theta_bioco2, allocatearray::c_trafficunits, allocatearray::c_traffprofwd, allocatearray::c_traffprofwe, allocatearray::c_wall_k1, allocatearray::c_wall_k1_bldgs, allocatearray::c_wall_k2, allocatearray::c_wall_k2_bldgs, allocatearray::c_wall_k3, allocatearray::c_wall_k3_bldgs, allocatearray::c_wall_k4, allocatearray::c_wall_k4_bldgs, allocatearray::c_wall_k5, allocatearray::c_wall_k5_bldgs, allocatearray::c_wall_rhocp1, allocatearray::c_wall_rhocp1_bldgs, allocatearray::c_wall_rhocp2, allocatearray::c_wall_rhocp2_bldgs, allocatearray::c_wall_rhocp3, allocatearray::c_wall_rhocp3_bldgs, allocatearray::c_wall_rhocp4, allocatearray::c_wall_rhocp4_bldgs, allocatearray::c_wall_rhocp5, allocatearray::c_wall_rhocp5_bldgs, allocatearray::c_wall_thick1, allocatearray::c_wall_thick1_bldgs, allocatearray::c_wall_thick2, allocatearray::c_wall_thick2_bldgs, allocatearray::c_wall_thick3, allocatearray::c_wall_thick3_bldgs, allocatearray::c_wall_thick4, allocatearray::c_wall_thick4_bldgs, allocatearray::c_wall_thick5, allocatearray::c_wall_thick5_bldgs, colnamesinputfiles::c_watercode, allocatearray::c_waterdepth, allocatearray::c_wetthresh, colnamesinputfiles::c_wgbldgscode, colnamesinputfiles::c_wgbsoilcode, colnamesinputfiles::c_wgdectrcode, colnamesinputfiles::c_wgevetrcode, colnamesinputfiles::c_wggrasscode, colnamesinputfiles::c_wgpavedcode, allocatearray::c_wgtobldgs, allocatearray::c_wgtobsoil, allocatearray::c_wgtodectr, allocatearray::c_wgtoevetr, allocatearray::c_wgtograss, allocatearray::c_wgtopaved, allocatearray::c_wgtorunoff, allocatearray::c_wgtosoilstore, allocatearray::c_wgtowater, colnamesinputfiles::c_wgwatercode, colnamesinputfiles::c_wprofautowd, colnamesinputfiles::c_wprofautowe, colnamesinputfiles::c_wprofmanuwd, colnamesinputfiles::c_wprofmanuwe, colnamesinputfiles::ca_ahmin_wd, colnamesinputfiles::ca_ahmin_we, colnamesinputfiles::ca_ahslopecooling_wd, colnamesinputfiles::ca_ahslopecooling_we, colnamesinputfiles::ca_ahslopeheating_wd, colnamesinputfiles::ca_ahslopeheating_we, colnamesinputfiles::ca_baset_hc, colnamesinputfiles::ca_co2mwd, colnamesinputfiles::ca_co2mwe, colnamesinputfiles::ca_co2pointsource, colnamesinputfiles::ca_ef_umolco2perj, colnamesinputfiles::ca_enef_v_jkm, colnamesinputfiles::ca_enprofwd, colnamesinputfiles::ca_enprofwe, colnamesinputfiles::ca_fcef_v_kgkmwd, colnamesinputfiles::ca_fcef_v_kgkmwe, colnamesinputfiles::ca_frfossilfuel_heat, colnamesinputfiles::ca_frfossilfuel_nonheat, colnamesinputfiles::ca_frpddwe, colnamesinputfiles::ca_maxfcmetab, colnamesinputfiles::ca_maxqfmetab, colnamesinputfiles::ca_minfcmetab, colnamesinputfiles::ca_minqfmetab, colnamesinputfiles::ca_popprofwd, colnamesinputfiles::ca_popprofwe, colnamesinputfiles::ca_qf_a1, colnamesinputfiles::ca_qf_a2, colnamesinputfiles::ca_qf_b1, colnamesinputfiles::ca_qf_b2, colnamesinputfiles::ca_qf_c1, colnamesinputfiles::ca_qf_c2, colnamesinputfiles::ca_tcriticcooling_wd, colnamesinputfiles::ca_tcriticcooling_we, colnamesinputfiles::ca_tcriticheating_wd, colnamesinputfiles::ca_tcriticheating_we, colnamesinputfiles::ca_trafficunits, colnamesinputfiles::ca_traffprofwd, colnamesinputfiles::ca_traffprofwe, colnamesinputfiles::cb_alpha, colnamesinputfiles::cb_alpha_enh, colnamesinputfiles::cb_beta, colnamesinputfiles::cb_min_r, colnamesinputfiles::cb_resp_a, colnamesinputfiles::cb_resp_b, colnamesinputfiles::cb_theta, colnamesinputfiles::cc_gsg1, colnamesinputfiles::cc_gsg2, colnamesinputfiles::cc_gsg3, colnamesinputfiles::cc_gsg4, colnamesinputfiles::cc_gsg5, colnamesinputfiles::cc_gsg6, colnamesinputfiles::cc_gskmax, colnamesinputfiles::cc_gsmodel, colnamesinputfiles::cc_gss1, colnamesinputfiles::cc_gss2, colnamesinputfiles::cc_gsth, colnamesinputfiles::cc_gstl, colnamesinputfiles::ce_alb_ibld, colnamesinputfiles::ce_ch_ibld, colnamesinputfiles::ce_ch_iroof, colnamesinputfiles::ce_ch_iwall, colnamesinputfiles::ce_em_ibld, colnamesinputfiles::ce_internal_k1, colnamesinputfiles::ce_internal_k2, colnamesinputfiles::ce_internal_k3, colnamesinputfiles::ce_internal_k4, colnamesinputfiles::ce_internal_k5, colnamesinputfiles::ce_internal_rhocp1, colnamesinputfiles::ce_internal_rhocp2, colnamesinputfiles::ce_internal_rhocp3, colnamesinputfiles::ce_internal_rhocp4, colnamesinputfiles::ce_internal_rhocp5, colnamesinputfiles::ce_internal_thick1, colnamesinputfiles::ce_internal_thick2, colnamesinputfiles::ce_internal_thick3, colnamesinputfiles::ce_internal_thick4, colnamesinputfiles::ce_internal_thick5, colnamesinputfiles::ce_nroom, colnamesinputfiles::ce_surf_k1, colnamesinputfiles::ce_surf_k2, colnamesinputfiles::ce_surf_k3, colnamesinputfiles::ce_surf_k4, colnamesinputfiles::ce_surf_k5, colnamesinputfiles::ce_surf_rhocp1, colnamesinputfiles::ce_surf_rhocp2, colnamesinputfiles::ce_surf_rhocp3, colnamesinputfiles::ce_surf_rhocp4, colnamesinputfiles::ce_surf_rhocp5, colnamesinputfiles::ce_surf_thick1, colnamesinputfiles::ce_surf_thick2, colnamesinputfiles::ce_surf_thick3, colnamesinputfiles::ce_surf_thick4, colnamesinputfiles::ce_surf_thick5, colnamesinputfiles::ce_wall_k1, colnamesinputfiles::ce_wall_k2, colnamesinputfiles::ce_wall_k3, colnamesinputfiles::ce_wall_k4, colnamesinputfiles::ce_wall_k5, colnamesinputfiles::ce_wall_rhocp1, colnamesinputfiles::ce_wall_rhocp2, colnamesinputfiles::ce_wall_rhocp3, colnamesinputfiles::ce_wall_rhocp4, colnamesinputfiles::ce_wall_rhocp5, colnamesinputfiles::ce_wall_thick1, colnamesinputfiles::ce_wall_thick2, colnamesinputfiles::ce_wall_thick3, colnamesinputfiles::ce_wall_thick4, colnamesinputfiles::ce_wall_thick5, colnamesinputfiles::ci_albmax, colnamesinputfiles::ci_albmin, colnamesinputfiles::ci_chanohm, colnamesinputfiles::ci_cpanohm, colnamesinputfiles::ci_drcoef1, colnamesinputfiles::ci_drcoef2, colnamesinputfiles::ci_dreq, colnamesinputfiles::ci_emis, colnamesinputfiles::ci_estmcode, colnamesinputfiles::ci_kkanohm, colnamesinputfiles::ci_ohmcode_sdry, colnamesinputfiles::ci_ohmcode_swet, colnamesinputfiles::ci_ohmcode_wdry, colnamesinputfiles::ci_ohmcode_wwet, colnamesinputfiles::ci_ohmthresh_sw, colnamesinputfiles::ci_ohmthresh_wd, colnamesinputfiles::ci_snowlimpat, colnamesinputfiles::ci_snowlimrem, colnamesinputfiles::ci_soiltcode, colnamesinputfiles::ci_statelimit, colnamesinputfiles::ci_stormax, colnamesinputfiles::ci_stormin, colnamesinputfiles::ci_wetthresh, colnamesinputfiles::cir_daywat1, colnamesinputfiles::cir_daywat7, colnamesinputfiles::cir_daywatper1, colnamesinputfiles::cir_daywatper7, colnamesinputfiles::cir_faut, colnamesinputfiles::cir_h_maintain, colnamesinputfiles::cir_ie_a1, colnamesinputfiles::cir_ie_a3, colnamesinputfiles::cir_ie_m1, colnamesinputfiles::cir_ie_m3, colnamesinputfiles::cir_ieend, colnamesinputfiles::cir_iestart, colnamesinputfiles::cir_intwu, colnamesinputfiles::co_a1, colnamesinputfiles::co_a2, colnamesinputfiles::co_a3, codematchanthropogenic(), codematchbiogen(), codematchconductance(), codematchdist(), codematchestm(), codematchestm_class(), codematchirrigation(), codematchnonveg(), codematchohm(), codematchprof(), codematchsnow(), codematchsoil(), codematchveg(), codematchwater(), allocatearray::conductance_coeff, allocatearray::conifsurf, colnamesinputfiles::cp_albmax, colnamesinputfiles::cp_albmin, colnamesinputfiles::cp_baset, colnamesinputfiles::cp_basete, colnamesinputfiles::cp_biogenco2code, colnamesinputfiles::cp_chanohm, colnamesinputfiles::cp_cpanohm, colnamesinputfiles::cp_drcoef1, colnamesinputfiles::cp_drcoef2, colnamesinputfiles::cp_dreq, colnamesinputfiles::cp_emis, colnamesinputfiles::cp_estmcode, colnamesinputfiles::cp_gddfull, colnamesinputfiles::cp_gsmax, colnamesinputfiles::cp_kkanohm, colnamesinputfiles::cp_laieq, colnamesinputfiles::cp_laimax, colnamesinputfiles::cp_laimin, colnamesinputfiles::cp_leafgp1, colnamesinputfiles::cp_leafgp2, colnamesinputfiles::cp_leafop1, colnamesinputfiles::cp_leafop2, colnamesinputfiles::cp_ohmcode_sdry, colnamesinputfiles::cp_ohmcode_swet, colnamesinputfiles::cp_ohmcode_wdry, colnamesinputfiles::cp_ohmcode_wwet, colnamesinputfiles::cp_ohmthresh_sw, colnamesinputfiles::cp_ohmthresh_wd, colnamesinputfiles::cp_porositymax, colnamesinputfiles::cp_porositymin, colnamesinputfiles::cp_sddfull, colnamesinputfiles::cp_snowlimpat, colnamesinputfiles::cp_soiltcode, colnamesinputfiles::cp_statelimit, colnamesinputfiles::cp_stormax, colnamesinputfiles::cp_stormin, colnamesinputfiles::cp_wetthresh, colnamesinputfiles::cpr_hours, colnamesinputfiles::cs_estmcode, colnamesinputfiles::cs_ohmcode_sdry, colnamesinputfiles::cs_ohmcode_swet, colnamesinputfiles::cs_ohmcode_wdry, colnamesinputfiles::cs_ohmcode_wwet, colnamesinputfiles::cs_ohmthresh_sw, colnamesinputfiles::cs_ohmthresh_wd, colnamesinputfiles::cs_snowalbmax, colnamesinputfiles::cs_snowalbmin, colnamesinputfiles::cs_snowcrwmax, colnamesinputfiles::cs_snowcrwmin, colnamesinputfiles::cs_snowemis, colnamesinputfiles::cs_snowplimalb, colnamesinputfiles::cs_snowplimsnow, colnamesinputfiles::cs_snowrmfactor, colnamesinputfiles::cs_snowsdmax, colnamesinputfiles::cs_snowsdmin, colnamesinputfiles::cs_snowtau_a, colnamesinputfiles::cs_snowtau_f, colnamesinputfiles::cs_snowtau_r, colnamesinputfiles::cs_snowtmfactor, colnamesinputfiles::cso_ksat, colnamesinputfiles::cso_obssmdepth, colnamesinputfiles::cso_obssmmax, colnamesinputfiles::cso_obssnrfrac, colnamesinputfiles::cso_soildens, colnamesinputfiles::cso_soildepth, colnamesinputfiles::cso_soilinfrate, colnamesinputfiles::cso_soilstcap, colnamesinputfiles::cw_albmax, colnamesinputfiles::cw_albmin, colnamesinputfiles::cw_chanohm, colnamesinputfiles::cw_cpanohm, colnamesinputfiles::cw_drcoef1, colnamesinputfiles::cw_drcoef2, colnamesinputfiles::cw_dreq, colnamesinputfiles::cw_emis, colnamesinputfiles::cw_estmcode, colnamesinputfiles::cw_kkanohm, colnamesinputfiles::cw_ohmcode_sdry, colnamesinputfiles::cw_ohmcode_swet, colnamesinputfiles::cw_ohmcode_wdry, colnamesinputfiles::cw_ohmcode_wwet, colnamesinputfiles::cw_ohmthresh_sw, colnamesinputfiles::cw_ohmthresh_wd, colnamesinputfiles::cw_statelimit, colnamesinputfiles::cw_stormax, colnamesinputfiles::cw_stormin, colnamesinputfiles::cw_waterdepth, colnamesinputfiles::cw_wetthresh, colnamesinputfiles::cwg_tobldgs, colnamesinputfiles::cwg_tobsoil, colnamesinputfiles::cwg_todectr, colnamesinputfiles::cwg_toevetr, colnamesinputfiles::cwg_tograss, colnamesinputfiles::cwg_topaved, colnamesinputfiles::cwg_torunoff, colnamesinputfiles::cwg_tosoilstore, colnamesinputfiles::cwg_towater, allocatearray::decidsurf, allocatearray::estmcoefficients_coeff, allocatearray::grasssurf, allocatearray::irrigation_coeff, initial::iv5, allocatearray::ivconif, allocatearray::ivdecid, allocatearray::ivgrass, allocatearray::ncolumnssiteselect, allocatearray::nonveg_coeff, allocatearray::nsurf, allocatearray::ohmcoefficients_coeff, allocatearray::pavsurf, allocatearray::profiles_coeff, allocatearray::siteselect, allocatearray::snow_coeff, allocatearray::soil_coeff, allocatearray::surfacechar, allocatearray::veg_coeff, allocatearray::water_coeff, allocatearray::watersurf, and allocatearray::wgwaterdist_coeff.

Referenced by suews_program().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ initialstate()

subroutine initialstate ( character(len=20)  GridName,
integer  year_int,
integer  Gridiv,
integer  NumberOfGrids 
)

Definition at line 1445 of file suews_ctrl_init.f95.

1446 ! Last modified HCW 13 Jan 2017 - Major changes to InitialConditions file
1447 ! Last modified HCW 24 May 2016 - Removed unused argument year_txt
1448 ! Last modified HCW 03 Jul 2015 - Added initial conditions albEveTr0 and albGrass0
1449 ! Last modified HCW 03 Dec 2014
1450 !------------------------------------------------------------------------
1451
1452 USE allocatearray
1453 USE data_in
1455 USE defaultnotused
1456 USE filename
1457 USE gis_data
1458 USE initialcond
1459 USE mod_z
1460 USE resist
1461 USE snowmod
1462 USE sues_data
1463 USE time
1464 USE initialcond
1467
1468 IMPLICIT NONE
1469
1470 CHARACTER(len=20) :: GridName !Name of the evaluated grid
1471 CHARACTER(len=4) :: year_txt
1472 INTEGER :: NumberOfGrids
1473
1474 CHARACTER(len=150) :: fileInit !Initial conditions filename
1475 INTEGER :: DaysSinceRain, Gridiv, & !number of days since rain, grid number,
1476 gamma1, gamma2 !switches related to cooling and heating degree days
1477 INTEGER :: wd, seas, date, mb, & !weekday information, season, date, month
1478 year_int, switch = 0, & !year as an integer, switch related to previous day
1479 id_next !next day,counter in irrigation calculations
1480
1481 REAL(KIND(1D0)) :: PavedState, BldgsState, EveTrState, DecTrState, GrassState, BSoilState, WaterState, &
1482 SnowFracPaved, SnowFracBldgs, SnowFracEveTr, SnowFracDecTr, &
1483 SnowFracGrass, SnowFracBSoil, SnowFracWater, &
1484 SnowDensPaved, SnowDensBldgs, SnowDensEveTr, SnowDensDecTr, &
1485 SnowDensGrass, SnowDensBSoil, SnowDensWater
1486
1487 INTEGER :: LeavesOutInitially !Allows for quick setting of veg-related initial conditions for full leaf-out (1) or leaf-off (0)
1488 INTEGER :: SnowInitially !Allows for quick setting of snow-related initial conditions for no snow initially (0)
1489
1490 INTEGER :: GridsInitialised = 0 ! Number of grids initialised at start of model run
1491 INTEGER :: YearsInitialised = 0 ! Number of years initialised at start of model run
1492 INTEGER :: dayofWeek_id(3)
1493 REAL(KIND(1D0)) :: NormalizeVegChar !Function
1494
1495 ! Define InitialConditions namelist ---------------------------------------
1496 namelist /initialconditions/ dayssincerain, &
1497 temp_c0, &
1498 !ID_Prev,& !Now calculated from met forcing file
1499 leavesoutinitially, &
1500 gdd_1_0, &
1501 gdd_2_0, &
1505 albevetr0, &
1506 albdectr0, &
1507 albgrass0, &
1508 decidcap0, &
1509 porosity0, &
1510 pavedstate, &
1511 bldgsstate, &
1512 evetrstate, &
1513 dectrstate, &
1514 grassstate, &
1515 bsoilstate, &
1516 waterstate, &
1523 snowinitially, &
1531 snowpackpaved, &
1532 snowpackbldgs, &
1533 snowpackevetr, &
1534 snowpackdectr, &
1535 snowpackgrass, &
1536 snowpackbsoil, &
1537 snowpackwater, &
1538 snowfracpaved, &
1539 snowfracbldgs, &
1540 snowfracevetr, &
1541 snowfracdectr, &
1542 snowfracgrass, &
1543 snowfracbsoil, &
1544 snowfracwater, &
1545 snowdenspaved, &
1546 snowdensbldgs, &
1547 snowdensevetr, &
1548 snowdensdectr, &
1549 snowdensgrass, &
1550 snowdensbsoil, &
1551 snowdenswater, &
1552 snowalb0 !,&
1553 ! BoInit ! removed as no longer needed by AnOHM
1554
1555 ! Initialise namelist to NAN ----------------------------------------------
1556 dayssincerain = int(nan)
1557 temp_c0 = nan
1558 leavesoutinitially = int(nan)
1559 gdd_1_0 = nan
1560 gdd_2_0 = nan
1564 albevetr0 = nan
1565 albdectr0 = nan
1566 albgrass0 = nan
1567 decidcap0 = nan
1568 porosity0 = nan
1569 pavedstate = nan
1570 bldgsstate = nan
1571 evetrstate = nan
1572 dectrstate = nan
1573 grassstate = nan
1574 bsoilstate = nan
1575 waterstate = nan
1582 snowinitially = int(nan)
1597 snowfracpaved = nan
1598 snowfracbldgs = nan
1599 snowfracevetr = nan
1600 snowfracdectr = nan
1601 snowfracgrass = nan
1602 snowfracbsoil = nan
1603 snowfracwater = nan
1604 snowdenspaved = nan
1605 snowdensbldgs = nan
1606 snowdensevetr = nan
1607 snowdensdectr = nan
1608 snowdensgrass = nan
1609 snowdensbsoil = nan
1610 snowdenswater = nan
1611 snowalb0 = nan
1612 ! BoInit=NAN
1613
1614 WRITE (year_txt, '(I4)') year_int !Get year as a text string
1615
1616 ! Define InitialConditions file -------------------------------------------
1617 fileinit = trim(fileinputpath)//trim("InitialConditions")//trim(gridname)//'.nml'
1618 ! On very first InitialConditions for each grid, can use one initial conditions file specified for all grids
1619 IF (multipleinitfiles == 0 .AND. yearsinitialised == 0) THEN
1620 fileinit = trim(fileinputpath)//trim("InitialConditions")//trim(filecode)//'_'//trim(year_txt)//'.nml'
1621 gridsinitialised = gridsinitialised + 1
1622 IF (gridsinitialised == numberofgrids) THEN
1623 yearsinitialised = yearsinitialised + 1
1624 gridsinitialised = 0 !reset GridsInitialised
1625 END IF
1626 END IF
1627 !write(*,*) TRIM(FileInit)
1628
1629 ! Open, read and close InitialConditions file -----------------------------
1630 OPEN (56, file=trim(fileinit), err=600, status='old')
1631 READ (56, iostat=ios_out, nml=initialconditions, err=601)
1632 CLOSE (56)
1633
1634 ! Write InitialConditions to FileChoices ----------------------------------
1635 filechoices = trim(fileoutputpath)//trim(filecode)//'_FileChoices.txt'
1636 OPEN (12, file=filechoices, position='append')
1637 WRITE (12, *) '----- '//trim("InitialConditions")//trim(gridname)//'.nml'//' -----'
1638 WRITE (12, nml=initialconditions)
1639 CLOSE (12)
1640
1641 !--------------------------------------------------------------------------
1642 ! Check initial conditions and assign values if not provided --------------
1643
1644 ! Calculate previous day --------------------------------------------------
1645 id_prev = int(metforcingdata(1, 2, gridiv)) - 1
1646
1647 ! If no. days since rainfall unknown, set to zero -------------------------
1648 IF (dayssincerain == int(nan)) dayssincerain = 0
1649
1650 ! If average temperature for previous day unknown, use average of first day
1651 IF (temp_c0 == nan) temp_c0 = sum(metforcingdata(1:(24*nsh), 12, gridiv))/(24*nsh)
1652
1653 ! Set vegetation-related initial conditions -------------------------------
1654 ! If LeavesOutInitially is -999, don't use and check all required conditions have been provided
1655 IF (leavesoutinitially == int(nan)) THEN
1656 IF (gdd_1_0 == nan .OR. gdd_2_0 == nan) THEN
1657 CALL errorhint(36, 'Specify values for GDD_1_0 and GDD_2_0.', notused, notused, notusedi)
1658 END IF
1659 IF (laiinitialevetr == nan .OR. laiinitialdectr == nan .OR. laiinitialgrass == nan) THEN
1660 CALL errorhint(36, 'Specify initial values for LAI for all vegetated surface types.', notused, notused, notusedi)
1661 END IF
1662 IF (albevetr0 == nan .OR. albdectr0 == nan .OR. albgrass0 == nan) THEN
1663 CALL errorhint(36, 'Specify initial values for albedo for all vegetated surface types.', notused, notused, notusedi)
1664 END IF
1665 IF (decidcap0 == nan) THEN
1666 CALL errorhint(36, 'Specify DecidCap0.', notused, notused, notusedi)
1667 END IF
1668 IF (porosity0 == nan) THEN
1669 CALL errorhint(36, 'Specify Porosity0.', notused, notused, notusedi)
1670 END IF
1671 ELSEIF (leavesoutinitially == 1) THEN !If leaves out, set to summertime values using SUEWS_Veg.txt
1673 gdd_2_0 = 0
1674 laiinitialevetr = surfacechar(gridiv, c_laimax(ivconif)) !Max LAI
1677 albevetr0 = surfacechar(gridiv, c_albmax(conifsurf)) !Max albedo
1680 decidcap0 = surfacechar(gridiv, c_stormax(decidsurf)) !Max storage capacity (DecTr only)
1681 porosity0 = surfacechar(gridiv, c_porositymin(ivdecid)) !Min porosity (DecTr only)
1682 ELSEIF (leavesoutinitially == 0) THEN !If leaves off, set to wintertime values using SUEWS_Veg.txt
1683 gdd_1_0 = 0
1685 laiinitialevetr = surfacechar(gridiv, c_laimin(ivconif)) !Min LAI
1688 albevetr0 = surfacechar(gridiv, c_albmin(conifsurf)) !Min albedo
1691 decidcap0 = surfacechar(gridiv, c_stormin(decidsurf)) !Min storage capacity (DecTr only)
1692 porosity0 = surfacechar(gridiv, c_porositymax(ivdecid)) !Max porosity (DecTr only)
1693 ELSE
1694 CALL errorhint(36, 'LeavesOutInitially must be 0, 1, or -999 (or omitted from InitialConditions namelist)', &
1696 END IF
1697
1698 ! If surface wetness states unknown, set to zero --------------------------
1699 IF (pavedstate == nan) pavedstate = 0
1700 IF (bldgsstate == nan) bldgsstate = 0
1701 IF (evetrstate == nan) evetrstate = 0
1702 IF (dectrstate == nan) dectrstate = 0
1703 IF (grassstate == nan) grassstate = 0
1704 IF (bsoilstate == nan) bsoilstate = 0
1705 ! except for water surface - set using WaterDepth in SUEWS_Water.txt
1706 IF (waterstate == nan) waterstate = surfacechar(gridiv, c_waterdepth)
1707
1708 ! Check initial soil moisture states are provided -------------------------
1711 CALL errorhint(36, 'Initial soil moisture must be provided for all surface types except water.', notused, notused, notusedi)
1712 END IF
1713
1714 ! Set snow-related initial conditions -------------------------------
1715 ! If snow part not used, or no snow initially, set all snow-related initial conditions to zero
1716 IF (snowuse == 0 .OR. snowinitially == 0) THEN
1724 snowpackpaved = 0
1725 snowpackbldgs = 0
1726 snowpackevetr = 0
1727 snowpackdectr = 0
1728 snowpackgrass = 0
1729 snowpackbsoil = 0
1730 snowpackwater = 0
1731 snowfracpaved = 0
1732 snowfracbldgs = 0
1733 snowfracevetr = 0
1734 snowfracdectr = 0
1735 snowfracgrass = 0
1736 snowfracbsoil = 0
1737 snowfracwater = 0
1738 snowdenspaved = 0
1739 snowdensbldgs = 0
1740 snowdensevetr = 0
1741 snowdensdectr = 0
1742 snowdensgrass = 0
1743 snowdensbsoil = 0
1744 snowdenswater = 0
1745 snowalb0 = 0
1746 ELSEIF (snowinitially == int(nan)) THEN !Check all required snow-related conditions are provided
1749 snowwaterwaterstate == nan) THEN
1750 CALL errorhint(36, 'Specify SnowWater state for all 7 surface types.', notused, notused, notusedi)
1751 END IF
1752 IF (snowpackpaved == nan .OR. snowpackbldgs == nan .OR. snowpackevetr == nan .OR. &
1753 snowpackdectr == nan .OR. snowpackgrass == nan .OR. snowpackbsoil == nan .OR. &
1754 snowpackwater == nan) THEN
1755 CALL errorhint(36, 'Specify SnowPack for all 7 surface types.', notused, notused, notusedi)
1756 END IF
1757 IF (snowfracpaved == nan .OR. snowfracbldgs == nan .OR. snowfracevetr == nan .OR. &
1758 snowfracdectr == nan .OR. snowfracgrass == nan .OR. snowfracbsoil == nan .OR. &
1759 snowfracwater == nan) THEN
1760 CALL errorhint(36, 'Specify SnowFrac for all 7 surface types.', notused, notused, notusedi)
1761 END IF
1762 IF (snowdenspaved == nan .OR. snowdensbldgs == nan .OR. snowdensevetr == nan .OR. &
1763 snowdensdectr == nan .OR. snowdensgrass == nan .OR. snowdensbsoil == nan .OR. &
1764 snowdenswater == nan) THEN
1765 CALL errorhint(36, 'Specify SnowDens for all 7 surface types.', notused, notused, notusedi)
1766 END IF
1767 IF (snowalb0 == nan) THEN
1768 CALL errorhint(36, 'Specify SnowAlb0.', notused, notused, notusedi)
1769 END IF
1770 ELSE
1771 CALL errorhint(36, 'SnowInitially must be 0 or -999 (or omitted from InitialConditions namelist)', &
1773 END IF
1774
1775 ! removed as no longer needed, TS 30 Jan 2018
1776 ! ! If AnOHM option selected, check initial Bowen ratio is provided ---------
1777 ! IF(StorageHeatMethod==3 .AND. BoInit == NAN) THEN
1778 ! CALL ErrorHint(36,'Specify BoInit for AnOHM calculations.', notUsed,notUsed,notUsedI)
1779 ! ENDIF
1780
1781 ! -------------------------------------------------------------------------
1782 ! -------------------------------------------------------------------------
1783
1784 ! Previous day DOY number (needed in file allocations)
1785 IF (id_prev >= 364) id_prev = 0 !If previous day is larger than 364, set this to zero
1786
1787 ! Save initial conditions to ModelDailyState array ------------------------
1794 modeldailystate(gridiv, cmds_gddmin) = 90 !QUESTION: Going to check for minimum GDD
1795 modeldailystate(gridiv, cmds_gddmax) = -90 !QUESTION: Going to check for maximum GDD
1801
1802 modeldailystate(gridiv, cmds_snowfallcum) = 0 !!Check this
1803
1804 modeldailystate(gridiv, cmds_dayssincerain) = real(dayssincerain, kind(1d0))
1806 ! Assume that the temperature has been the same for the previous days
1810
1811 ! -- Anthropogenic heat flux initializations --
1812 ! Need to get BaseT_HC from SurfaceChar, as info not transferred until SUEWS_Translate called
1814
1815 IF (emissionsmethod >= 0) THEN
1816 !Calculations related to heating and cooling degree days (BaseT_HC is used always)
1817 IF ((temp_c0 - baset_hc) >= 0) THEN !Cooling
1818 gamma2 = 1
1819 ELSE
1820 gamma2 = 0
1821 END IF
1822 IF ((baset_hc - temp_c0) >= 0) THEN !Heating
1823 gamma1 = 1
1824 ELSE
1825 gamma1 = 0
1826 END IF
1827 modeldailystate(gridiv, cmds_hdd1) = gamma1*(baset_hc - temp_c0) ! Heating
1828 modeldailystate(gridiv, cmds_hdd2) = gamma2*(temp_c0 - baset_hc) ! Cooling
1829 END IF
1830
1831 ! -- Save snow density and snow albedo info in InitialConditions to ModelDailyState array --
1832 modeldailystate(gridiv, cmds_snowdens(pavsurf)) = snowdenspaved
1833 modeldailystate(gridiv, cmds_snowdens(bldgsurf)) = snowdensbldgs
1834 modeldailystate(gridiv, cmds_snowdens(conifsurf)) = snowdensevetr
1835 modeldailystate(gridiv, cmds_snowdens(decidsurf)) = snowdensdectr
1836 modeldailystate(gridiv, cmds_snowdens(grasssurf)) = snowdensgrass
1837 modeldailystate(gridiv, cmds_snowdens(bsoilsurf)) = snowdensbsoil
1838 modeldailystate(gridiv, cmds_snowdens(watersurf)) = snowdenswater
1839
1841
1842 ! -------------------------------------------------------------------------
1843
1844 ! Saving to ModelOutputData array -----------------------------------------
1845
1846 ! -- Initial wetness status of each surface (above ground) --
1847 modeloutputdata(0, cmod_state(pavsurf), gridiv) = pavedstate
1848 modeloutputdata(0, cmod_state(bldgsurf), gridiv) = bldgsstate
1849 modeloutputdata(0, cmod_state(conifsurf), gridiv) = evetrstate
1850 modeloutputdata(0, cmod_state(decidsurf), gridiv) = dectrstate
1851 modeloutputdata(0, cmod_state(grasssurf), gridiv) = grassstate
1852 modeloutputdata(0, cmod_state(bsoilsurf), gridiv) = bsoilstate
1853 modeloutputdata(0, cmod_state(watersurf), gridiv) = waterstate
1854
1855 ! -- Initial soil stores for each surface (below ground) --
1862 modeloutputdata(0, cmod_soilstate(watersurf), gridiv) = 0 ! No soil layer for water surface
1863
1864 ! -- Initial liquid (melted) water for each surface --
1872
1873 ! -- Initial snow water equivalent for each surface --
1881
1882 ! -- Initial fraction of snow on each surface --
1883 modeloutputdata(0, cmod_snowfrac(pavsurf), gridiv) = snowfracpaved
1884 modeloutputdata(0, cmod_snowfrac(bldgsurf), gridiv) = snowfracbldgs
1885 modeloutputdata(0, cmod_snowfrac(conifsurf), gridiv) = snowfracevetr
1886 modeloutputdata(0, cmod_snowfrac(decidsurf), gridiv) = snowfracdectr
1887 modeloutputdata(0, cmod_snowfrac(grasssurf), gridiv) = snowfracgrass
1888 modeloutputdata(0, cmod_snowfrac(bsoilsurf), gridiv) = snowfracbsoil
1889 modeloutputdata(0, cmod_snowfrac(watersurf), gridiv) = snowfracwater
1890
1891 icefrac = 0.2 !Estimated fraction of ice. Should be improved in the future
1892
1893 ! initialise surface temperatures for ESTM_ext
1894 IF (storageheatmethod == 5 .OR. netradiationmethod > 1000) THEN
1895
1896 tsfc_roof_grids(gridiv, :) = temp_c0
1897 tsfc_wall_grids(gridiv, :) = temp_c0
1898 END IF
1899 tsfc_surf_grids(gridiv, :) = temp_c0
1900 tin_surf_grids(gridiv, :) = temp_c0
1901
1902 ! At this point translate arrays to variables (needed for SUEWS_cal_RoughnessParameters)
1903 IF (diagnose == 1) print *, 'calling in initial state: SUEWS_Translate'
1904 CALL suews_translate(gridiv, 0, 0)
1905
1906 !Calculation of roughness parameters (N.B. uses porosity)
1907 IF (diagnose == 1) print *, 'calling in initial state: SUEWS_cal_RoughnessParameters'
1909 roughlenmommethod, sfr_surf, & !input
1912 z0m_in, zdm_in, z, &
1913 fai, pai, & !output
1914 zh, z0m, zdm, zzd)
1915
1916 !=============================================================================
1917 ! If the run start day is at previous year, then calculate the number of days
1918 ! in that year.
1919
1920 !First we need to know if the previous day given in initial conditions (id_prev) is
1921 !on previous year as this is needed in the initialization of DayofWeek matrix.
1922 !In this case switch is set to one for date calculations.
1923 IF (id_prev == 0) THEN !If id_prev = 0, means that the first modelled day is 1 Jan
1924 year_int = year_int - 1 !1) find the number of days on that year
1925 CALL leapyearcalc(year_int, id_prev) !2) set switch to 1 so that the code knows to change back to current year (switch=0)
1926 switch = 1
1927 END IF
1928
1929 CALL day2month(id_prev, mb, date, seas, year_int, lat) !Calculate date information (mb = month, date = day,...)
1930 CALL day_of_week(date, mb, year_int, wd) !Calculate weekday of the previous day (wd) (1=Sun, ..., 7=Sat)
1931
1932 !After the day in previous year switch is changed back to zero:
1933 ! ie not previous day anymore
1934 !Also the size of DayofWeek is from 0:NdaysinYear meaning
1935 !that in zero slot is the previous day information
1936 IF (switch == 1) THEN
1937 year_int = year_int + 1
1938 id_prev = 0
1939 switch = 0
1940 END IF
1941
1942 ! DayofWeek(id_prev,1)=wd ! day of week
1943 ! DayofWeek(id_prev,2)=mb ! month
1944 ! DayofWeek(id_prev,3)=seas ! season (summer=1, winter=2) needed for accumulation
1945
1946 ! in case next day goes to next year calculate again the date information for DayofWeek matrix.
1947 id_next = id_prev + 1
1948 IF (id_next > nofdaysthisyear) THEN
1949 id_next = 1
1950 year_int = year_int + 1
1951 switch = 1
1952 CALL errorhint(43, 'switch- years', notused, notused, notusedi)
1953 END IF
1954
1955 CALL day2month(id_next, mb, date, seas, year_int, lat) !Calculate real date from doy
1956 CALL day_of_week(date, mb, year_int, wd) !Calculate weekday (1=Sun, ..., 7=Sat)
1957
1958 IF (switch == 1) THEN
1959 iy = iy - 1
1960 switch = 0
1961 END IF
1962
1963 ! DayofWeek(id_next,1)=wd ! day of week
1964 ! DayofWeek(id_next,2)=mb ! month
1965 ! DayofWeek(id_next,3)=seas ! season
1966
1967 !=============================================================================
1968
1969 !id=id_prev
1970 !it= 23 !!LastTimeOfDay
1971 IF (id_prev >= startdls .AND. id_prev <= enddls) THEN !Summertime
1972 dls = 1
1973 ELSE
1974 dls = 0
1975 END IF
1976
1977 ! -----------------------------------------------------------------------
1978 ! Calculate daily water use if modelled (i.e. if WaterUseMethod = 0).
1979 ! Calculated from previous day information given in InitialConditions file
1980 CALL suews_cal_weekday( &
1981 iy, id, lat, & !input
1982 dayofweek_id) !output
1983
1984 state_surf = [pavedstate, bldgsstate, evetrstate, dectrstate, grassstate, bsoilstate, waterstate]
1987 CALL update_wateruse( &
1988 id, waterusemethod, dayofweek_id, lat, faut, hdd_id, & !input
1991 wuday_id) !output
1992
1993 ! ---- AnOHM TS ---------------------
1994 ! initialize Bowen ratio
1995 ! Bo_grids(0,:)=2.
1996 ! mAH_grids(0,:)=25.
1997
1998 ! -----------------------------------
1999 CALL suews_translateback(gridiv, 0, 0)
2000
2001 RETURN
2002
2003600 CALL errorhint(47, trim(fileinit), notused, notused, notusedi)
2004601 CALL errorhint(48, trim(fileinit), notused, notused, ios_out)
2005
real(kind(1d0)), dimension(:, :), allocatable modeldailystate
real(kind(1d0)), dimension(:, :), allocatable tsfc_wall_grids
real(kind(1d0)), dimension(:, :), allocatable tin_surf_grids
real(kind(1d0)), dimension(nsurf) state_surf
real(kind(1d0)) porosity_id
real(kind(1d0)), dimension(:, :), allocatable tsfc_roof_grids
real(kind(1d0)), dimension(nsurf) sfr_surf
real(kind(1d0)), dimension(nsurf) soilstore_surf
real(kind(1d0)), dimension(:, :), allocatable tsfc_surf_grids
integer, dimension(nsurf) cmod_snowfrac
integer, dimension(nsurf) cmds_snowdens
integer, dimension(nsurf) cmod_snowwaterstate
real(kind(1d0)), dimension(:, :, :), allocatable metforcingdata
integer, dimension(nsurf) cmod_snowpack
real(kind(1d0)), dimension(9) wuday_id
integer, dimension(nsurf) cmod_soilstate
real(kind(1d0)), dimension(:, :, :), allocatable modeloutputdata
integer, dimension(nsurf) cmod_state
real(kind(1d0)), dimension(12) hdd_id
real(kind(1d0)), dimension(nsurf) icefrac
subroutine update_wateruse(id, WaterUseMethod, DayofWeek_id, lat, FrIrriAuto, HDD_id, state_id, soilstore_id, SoilStoreCap, H_maintain, Ie_a, Ie_m, Ie_start, Ie_end, DayWatPer, DayWat, WUDay_id)
integer netradiationmethod
character(len=20) filecode
character(len=150) fileoutputpath
character(len=150) filechoices
real(kind(1d0)) baset_hc
integer emissionsmethod
integer startdls
integer diagnose
integer multipleinitfiles
integer waterusemethod
integer storageheatmethod
character(len=150) fileinputpath
integer roughlenmommethod
real(kind(1d0)) nan
real(kind(1d0)) notused
real(kind(1d0)) faibldg
real(kind(1d0)) bldgh
real(kind(1d0)) evetreeh
real(kind(1d0)) dectreeh
real(kind(1d0)) faievetree
real(kind(1d0)) faidectree
real(kind(1d0)) porosity0
real(kind(1d0)) decidcap0
real(kind(1d0)) albdectr0
real(kind(1d0)) albevetr0
real(kind(1d0)) snowwaterwaterstate
real(kind(1d0)) albgrass0
real(kind(1d0)) snowpackwater
real(kind(1d0)) snowalb0
real(kind(1d0)) gdd_2_0
real(kind(1d0)) gdd_1_0
real(kind(1d0)) zdm
real(kind(1d0)) z
real(kind(1d0)) z0m
real(kind(1d0)) zdm_in
real(kind(1d0)) zzd
real(kind(1d0)) z0m_in
subroutine suews_cal_roughnessparameters(RoughLenMomMethod, sfr_surf, bldgH, EveTreeH, DecTreeH, porosity_dectr, FAIBldg, FAIEveTree, FAIDecTree, z0m_in, zdm_in, Z, FAI, PAI, Zh, z0m, zdm, ZZD)
real(kind(1d0)) faut
real(kind(1d0)) zh
real(kind(1d0)) fai
real(kind(1d0)), dimension(7) daywat
real(kind(1d0)) pai
real(kind(1d0)), dimension(3) ie_a
real(kind(1d0)), dimension(7) daywatper
real(kind(1d0)), dimension(3) ie_m
real(kind(1d0)) h_maintain
integer dls
integer nofdaysthisyear
integer iy
real(kind(1d0)) function normalizevegchar(VegCol, Gridiv)
subroutine suews_translateback(Gridiv, ir, irMax)
subroutine suews_translate(Gridiv, ir, iMB)
subroutine suews_cal_weekday(iy, id, lat, dayofWeek_id)
subroutine leapyearcalc(year_int, nroDays)
subroutine day_of_week(DATE, MONTH, YEAR, DOW)
subroutine day2month(b, mb, md, seas, year, latitude)

References initialcond::albdectr0, initialcond::albevetr0, initialcond::albgrass0, data_in::baset_hc, gis_data::bldgh, allocatearray::bldgsurf, allocatearray::bsoilsurf, allocatearray::c_albmax, allocatearray::c_albmin, allocatearray::c_baset_hc, allocatearray::c_gddfull, allocatearray::c_laimax, allocatearray::c_laimin, allocatearray::c_porositymax, allocatearray::c_porositymin, allocatearray::c_sddfull, allocatearray::c_stormax, allocatearray::c_stormin, allocatearray::c_waterdepth, colnamesmodeldailystate::cmds_albdectr, colnamesmodeldailystate::cmds_albevetr, colnamesmodeldailystate::cmds_albgrass, colnamesmodeldailystate::cmds_dayssincerain, colnamesmodeldailystate::cmds_decidcap, colnamesmodeldailystate::cmds_gdd1_0, colnamesmodeldailystate::cmds_gdd2_0, colnamesmodeldailystate::cmds_gddmax, colnamesmodeldailystate::cmds_gddmin, colnamesmodeldailystate::cmds_hdd1, colnamesmodeldailystate::cmds_hdd2, colnamesmodeldailystate::cmds_id_prev, colnamesmodeldailystate::cmds_laiinitialdectr, colnamesmodeldailystate::cmds_laiinitialevetr, colnamesmodeldailystate::cmds_laiinitialgrass, colnamesmodeldailystate::cmds_porosity, colnamesmodeldailystate::cmds_snowalb, allocatearray::cmds_snowdens, colnamesmodeldailystate::cmds_snowfallcum, colnamesmodeldailystate::cmds_tempc, colnamesmodeldailystate::cmds_tempcold1, colnamesmodeldailystate::cmds_tempcold2, colnamesmodeldailystate::cmds_tempcold3, allocatearray::cmod_snowfrac, allocatearray::cmod_snowpack, allocatearray::cmod_snowwaterstate, allocatearray::cmod_soilstate, allocatearray::cmod_state, allocatearray::conifsurf, day2month(), day_of_week(), sues_data::daywat, sues_data::daywatper, initialcond::decidcap0, allocatearray::decidsurf, gis_data::dectreeh, data_in::diagnose, time::dls, data_in::emissionsmethod, data_in::enddls, errorhint(), gis_data::evetreeh, sues_data::fai, gis_data::faibldg, gis_data::faidectree, gis_data::faievetree, sues_data::faut, data_in::filechoices, data_in::filecode, data_in::fileinputpath, data_in::fileoutputpath, initialcond::gdd_1_0, initialcond::gdd_2_0, allocatearray::grasssurf, sues_data::h_maintain, allocatearray::hdd_id, allocatearray::icefrac, time::id, initialcond::id_prev, sues_data::ie_a, sues_data::ie_end, sues_data::ie_m, sues_data::ie_start, defaultnotused::ios_out, allocatearray::ivconif, allocatearray::ivdecid, allocatearray::ivgrass, time::iy, initialcond::laiinitialdectr, initialcond::laiinitialevetr, initialcond::laiinitialgrass, data_in::lat, leapyearcalc(), allocatearray::metforcingdata, allocatearray::modeldailystate, allocatearray::modeloutputdata, data_in::multipleinitfiles, defaultnotused::nan, data_in::netradiationmethod, time::nofdaysthisyear, defaultnotused::notused, defaultnotused::notusedi, sues_data::nsh, sues_data::pai, allocatearray::pavsurf, initialcond::porosity0, allocatearray::porosity_id, data_in::roughlenmommethod, allocatearray::sfr_surf, initialcond::snowalb0, initialcond::snowpackbldgs, initialcond::snowpackbsoil, initialcond::snowpackdectr, initialcond::snowpackevetr, initialcond::snowpackgrass, initialcond::snowpackpaved, initialcond::snowpackwater, data_in::snowuse, initialcond::snowwaterbldgsstate, initialcond::snowwaterbsoilstate, initialcond::snowwaterdectrstate, initialcond::snowwaterevetrstate, initialcond::snowwatergrassstate, initialcond::snowwaterpavedstate, initialcond::snowwaterwaterstate, allocatearray::soilstore_surf, initialcond::soilstorebldgsstate, initialcond::soilstorebsoilstate, allocatearray::soilstorecap_surf, initialcond::soilstoredectrstate, initialcond::soilstoreevetrstate, initialcond::soilstoregrassstate, initialcond::soilstorepavedstate, data_in::startdls, allocatearray::state_surf, data_in::storageheatmethod, resist_module::suews_cal_roughnessparameters(), suews_cal_weekday(), suews_translate(), suews_translateback(), allocatearray::surfacechar, initialcond::temp_c0, allocatearray::tin_surf_grids, allocatearray::tsfc_roof_grids, allocatearray::tsfc_surf_grids, allocatearray::tsfc_wall_grids, dailystate_module::update_wateruse(), allocatearray::watersurf, data_in::waterusemethod, allocatearray::wuday_id, mod_z::z, mod_z::z0m, mod_z::z0m_in, mod_z::zdm, mod_z::zdm_in, sues_data::zh, and mod_z::zzd.

Referenced by suews_program().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ nextinitial()

subroutine nextinitial ( character(len=15)  GridName,
integer  year_int 
)

Definition at line 2035 of file suews_ctrl_init.f95.

2036 ! Last modified HCW 24 May 2016
2037 ! Year of InitialConditions output file fixed. _EndofRun appended to files where the run finishes before the year end
2038 ! Last modified LJ 06 Jul 2015
2039 ! Initial conditions of SnowAlb added, densSnow changed to SnowDens
2040 ! Last modified HCW 03 Jul 2015
2041 ! Added initial conditions albEveTr0 and albGrass0
2042 ! Modified by HCW 21 Nov 2014
2043 ! Last day of year is not anymore the number of days on that year, but rather
2044 ! id == 1. Thus nofDaysThisYear was changed to 1. LJ 9/4/2015
2045 !------------------------------------------------------------------------
2046
2047 USE allocatearray
2050 USE data_in
2051 USE defaultnotused
2052 USE initial
2053 USE sues_data
2054 USE snowmod
2055 USE time
2056 USE initialcond
2057
2058 IMPLICIT NONE
2059
2060 CHARACTER(len=15) :: GridName
2061 CHARACTER(len=4) :: year_txt2
2062 INTEGER :: year_int2
2063 INTEGER :: year_int
2064 INTEGER :: ID_Prev_Out ! ID_Prev written to next Initial Conditions file
2065 INTEGER :: nofDaysThisYear_ForOutput !Added HCW 13 Jan 2017
2066
2067 year = year_int !HCW added 21 Nov 2014
2068
2069 ! Modified by HCW 24 May 2016
2070 IF (id == 1 .AND. iy == (year + 1)) THEN !if id = 1 and this is the first row of next year
2071 year_int2 = int(year + 1)
2072 WRITE (year_txt2, '(I4)') year_int2
2073 OPEN (57, &
2074 file=trim(fileinputpath)//trim("InitialConditions")//trim(gridname)//'_'//trim(adjustl(year_txt2))//'.nml', &
2075 err=200)
2076 nofdaysthisyear_foroutput = nofdaysthisyear
2077 ELSE
2078 year_int2 = int(year) !End of Run but not end of year
2079 WRITE (year_txt2, '(I4)') year_int2
2080 OPEN (57, &
2081 file=trim(fileinputpath)//trim("InitialConditions")//trim(gridname)//'_'//trim(adjustl(year_txt2))//'_EndofRun.nml', &
2082 err=201)
2083 nofdaysthisyear_foroutput = id - 1
2084 END IF
2085 id_prev_out = (id - 1)
2086
2087 !! If last time of day, then DailyState variables will have been updated so can write out arrays for id rather than id-1
2088 !if(it==23 .and. imin == (nsh_real-1)/nsh_real*60) then !!LastTimeofday
2089 ! id=id+1
2090 !endif
2091 WRITE (57, *) '&InitialConditions'
2092 WRITE (57, *) 'DaysSinceRain=', int(hdd_id(12))
2093 WRITE (57, *) 'Temp_C0=', hdd_id(9)
2094 !WRITE(57,*)'ID_Prev=',ID_Prev_Out !No longer included in initial conditions (HCW 13 Jan 2017)
2095 WRITE (57, *) 'GDD_1_0=', gdd_id(1)
2096 WRITE (57, *) 'GDD_2_0=', gdd_id(2)
2097 WRITE (57, *) 'LAIinitialEveTr=', lai_id(ivconif)
2098 WRITE (57, *) 'LAIinitialDecTr=', lai_id(ivdecid)
2099 WRITE (57, *) 'LAIinitialGrass=', lai_id(ivgrass)
2100 WRITE (57, *) 'AlbEveTr0=', albevetr_id
2101 WRITE (57, *) 'AlbDecTr0=', albdectr_id
2102 WRITE (57, *) 'AlbGrass0=', albgrass_id
2103 WRITE (57, *) 'DecidCap0=', decidcap_id
2104 WRITE (57, *) 'Porosity0=', porosity_id
2105 WRITE (57, *) 'SoilStorePavedState=', soilstore_surf(pavsurf)
2106 WRITE (57, *) 'SoilStoreBldgsState=', soilstore_surf(bldgsurf)
2107 WRITE (57, *) 'SoilStoreEveTrState=', soilstore_surf(conifsurf)
2108 WRITE (57, *) 'SoilStoreDecTrState=', soilstore_surf(decidsurf)
2109 WRITE (57, *) 'SoilStoreGrassState=', soilstore_surf(grasssurf)
2110 WRITE (57, *) 'SoilStoreBSoilState=', soilstore_surf(bsoilsurf)
2111 WRITE (57, *) 'PavedState=', state_surf(pavsurf)
2112 WRITE (57, *) 'BldgsState=', state_surf(bldgsurf)
2113 WRITE (57, *) 'EveTrState=', state_surf(conifsurf)
2114 WRITE (57, *) 'DecTrState=', state_surf(decidsurf)
2115 WRITE (57, *) 'GrassState=', state_surf(grasssurf)
2116 WRITE (57, *) 'BSoilState=', state_surf(bsoilsurf)
2117 WRITE (57, *) 'WaterState=', state_surf(watersurf)
2118 ! Only write snow variables if snow part is running
2119 IF (snowuse == 1) THEN
2120 WRITE (57, *) 'SnowWaterPavedState=', snowwater(pavsurf)
2121 WRITE (57, *) 'SnowWaterBldgsState=', snowwater(bldgsurf)
2122 WRITE (57, *) 'SnowWaterEveTrState=', snowwater(conifsurf)
2123 WRITE (57, *) 'SnowWaterDecTrState=', snowwater(decidsurf)
2124 WRITE (57, *) 'SnowWaterGrassState=', snowwater(grasssurf)
2125 WRITE (57, *) 'SnowWaterBSoilState=', snowwater(bsoilsurf)
2126 WRITE (57, *) 'SnowWaterWaterState=', snowwater(watersurf)
2127 WRITE (57, *) 'SnowPackPaved=', snowpack(pavsurf)
2128 WRITE (57, *) 'SnowPackBldgs=', snowpack(bldgsurf)
2129 WRITE (57, *) 'SnowPackEveTr=', snowpack(conifsurf)
2130 WRITE (57, *) 'SnowPackDecTr=', snowpack(decidsurf)
2131 WRITE (57, *) 'SnowPackGrass=', snowpack(grasssurf)
2132 WRITE (57, *) 'SnowPackBSoil=', snowpack(bsoilsurf)
2133 WRITE (57, *) 'SnowPackWater=', snowpack(watersurf)
2134 WRITE (57, *) 'SnowFracPaved=', snowfrac(pavsurf)
2135 WRITE (57, *) 'SnowFracBldgs=', snowfrac(bldgsurf)
2136 WRITE (57, *) 'SnowFracEveTr=', snowfrac(conifsurf)
2137 WRITE (57, *) 'SnowFracDecTr=', snowfrac(decidsurf)
2138 WRITE (57, *) 'SnowFracGrass=', snowfrac(grasssurf)
2139 WRITE (57, *) 'SnowFracBSoil=', snowfrac(bsoilsurf)
2140 WRITE (57, *) 'SnowFracWater=', snowfrac(watersurf)
2141 WRITE (57, *) 'SnowDensPaved=', snowdens(pavsurf)
2142 WRITE (57, *) 'SnowDensBldgs=', snowdens(bldgsurf)
2143 WRITE (57, *) 'SnowDensEveTr=', snowdens(conifsurf)
2144 WRITE (57, *) 'SnowDensDecTr=', snowdens(decidsurf)
2145 WRITE (57, *) 'SnowDensGrass=', snowdens(grasssurf)
2146 WRITE (57, *) 'SnowDensBSoil=', snowdens(bsoilsurf)
2147 WRITE (57, *) 'SnowDensWater=', snowdens(watersurf)
2148 WRITE (57, *) 'SnowAlb0=', snowalb
2149 END IF
2150 ! WRITE(57,*)'BoInit=',BoInit
2151 WRITE (57, *) '/'
2152 CLOSE (57)
2153
2154 IF (it == 23 .AND. imin == (nsh_real - 1)/nsh_real*60) THEN
2155 id = id - 1
2156 END IF
2157
2158 RETURN
2159
2160200 CALL errorhint(49, trim("InitialConditions")//trim(gridname)// &
2161 '_'//trim(adjustl(year_txt2))//'.nml', notused, notused, notusedi)
2162201 CALL errorhint(49, trim("InitialConditions")//trim(gridname)// &
2163 '_'//trim(adjustl(year_txt2))//'EoR.nml', notused, notused, notusedi)
2164
real(kind(1d0)), dimension(nsurf) snowpack
real(kind(1d0)), dimension(nvegsurf) lai_id
real(kind(1d0)), dimension(nvegsurf) gdd_id
real(kind(1d0)) decidcap_id
real(kind(1d0)) albevetr_id
real(kind(1d0)), dimension(nsurf) snowwater
real(kind(1d0)) albdectr_id
real(kind(1d0)), dimension(nsurf) snowdens
real(kind(1d0)) albgrass_id
real(kind(1d0)), dimension(nsurf) snowfrac
real(kind(1d0)) year
real(kind(1d0)) snowalb
real(kind(1d0)) nsh_real
integer it
integer imin

References allocatearray::albdectr_id, allocatearray::albevetr_id, allocatearray::albgrass_id, allocatearray::bldgsurf, allocatearray::bsoilsurf, allocatearray::conifsurf, allocatearray::decidcap_id, allocatearray::decidsurf, errorhint(), data_in::fileinputpath, allocatearray::gdd_id, allocatearray::grasssurf, allocatearray::hdd_id, time::id, time::imin, time::it, allocatearray::ivconif, allocatearray::ivdecid, allocatearray::ivgrass, time::iy, allocatearray::lai_id, time::nofdaysthisyear, defaultnotused::notused, defaultnotused::notusedi, sues_data::nsh_real, allocatearray::pavsurf, allocatearray::porosity_id, snowmod::snowalb, allocatearray::snowdens, allocatearray::snowfrac, allocatearray::snowpack, data_in::snowuse, allocatearray::snowwater, allocatearray::soilstore_surf, allocatearray::state_surf, allocatearray::watersurf, and data_in::year.

Referenced by suews_program().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ normalizevegchar()

real(kind(1d0)) function normalizevegchar ( integer, dimension(nvegsurf)  VegCol,
integer  Gridiv 
)

Definition at line 2009 of file suews_ctrl_init.f95.

2010
2011 USE allocatearray
2013
2014 IMPLICIT NONE
2015
2016 INTEGER, DIMENSION(nvegsurf) :: VegCol !Must be column numbers defined for veg surfaces only
2017 INTEGER :: Gridiv
2018 REAL(KIND(1D0)) :: NormVegResult
2019 IF (surfacechar(gridiv, c_frevetr) + &
2020 surfacechar(gridiv, c_frdectr) + &
2021 surfacechar(gridiv, c_frgrass) == 0.) THEN ! prevent arithmetic error under a full impervious scenario
2022 normvegresult = 0.
2023 ELSE
2024 normvegresult = (surfacechar(gridiv, vegcol(ivconif))*surfacechar(gridiv, c_frevetr) + &
2025 surfacechar(gridiv, vegcol(ivdecid))*surfacechar(gridiv, c_frdectr) + &
2026 surfacechar(gridiv, vegcol(ivgrass))*surfacechar(gridiv, c_frgrass))/ &
2027 (surfacechar(gridiv, c_frevetr) + surfacechar(gridiv, c_frdectr) + surfacechar(gridiv, c_frgrass))
2028 END IF
2029
2030 RETURN

References colnamesinputfiles::c_frdectr, colnamesinputfiles::c_frevetr, colnamesinputfiles::c_frgrass, allocatearray::ivconif, allocatearray::ivdecid, allocatearray::ivgrass, and allocatearray::surfacechar.

◆ numberrows()

subroutine numberrows ( character(len=50)  FileN,
integer  SkipHeaderLines 
)

Definition at line 405 of file suews_ctrl_init.f95.

406
407 USE data_in
409 USE initial
410
411 IMPLICIT NONE
412
413 CHARACTER(len=50) :: FileN
414 INTEGER :: SkipHeaderLines, RunNumber
415 INTEGER :: SkipCounter
416 INTEGER :: ios
417
418 WRITE (*, *) filen
419 OPEN (39, file=trim(fileinputpath)//trim(filen), err=204, status='old')
420
421 IF (skipheaderlines > 0) THEN
422 DO skipcounter = 1, skipheaderlines
423 READ (39, *, err=205)
424 !write(*,*) SkipCounter, SkipHeaderLines
425 END DO
426 END IF
427
428 nlines = 0 !Initialize nlines
429 DO
430 READ (39, *, iostat=ios) runnumber
431 IF (ios < 0 .OR. runnumber == -9) EXIT !IF(RunNumber==-9) EXIT
432 nlines = nlines + 1
433 END DO
434 !write(*,*) 'nlines read: ',nlines
435 CLOSE (39)
436
437 RETURN
438
439204 CALL errorhint(47, trim(fileinputpath)//trim(filen), notused, notused, notusedi)
440205 CALL errorhint(48, trim(fileinputpath)//trim(filen), notused, notused, notusedi)
441

References errorhint(), data_in::fileinputpath, initial::nlines, defaultnotused::notused, and defaultnotused::notusedi.

Referenced by overallruncontrol().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ overallruncontrol()

subroutine overallruncontrol

Definition at line 5 of file suews_ctrl_init.f95.

6 ! Last modified:
7 ! MH 21 Jun 2017 - Added anthropogenic CO2 parameters and changed AnthropogenicHeat to Anthropogenic
8 ! MH 16 Jun 2017 - Added biogenic CO2 parameters
9 ! HCW 21 Apr 2017 - Added new method for precip disaggregation
10 ! HCW 13 Jan 2017 - Changes to RunControl and InitialConditions
11 ! HCW 04 Nov 2016 - minor bug fix in LAImin/LAImax warnings related to 3 veg surface types out cf 7 surface types
12 ! LJ 27 Jan 2016 - Removal of tabs, cleaning of the code
13 ! HCW 06 Mar 2015 - Removed options 10,20,30 (NARPOutput) for NetRadiationMethod
14 ! HCW 06 Feb 2015 - File ID numbers changed so they are unique
15 ! HCW 19 Dec 2014
16 ! To Do:
17 ! - Holidays.txt input file needs to be read in and coded into model
18 ! - Add column header checks for SiteSelect
19 !-------------------------------------------------------------------------
20
23 USE data_in
25 USE filename
26 USE initial
27 USE gis_data
28 USE mod_z
29 USE resist
30 USE snowmod
31 USE sues_data
32 USE time
33
34 IMPLICIT NONE
35
36 INTEGER :: iv, i, SkipCounter, iFile !iv and i, ii are integers used in do loops
37 CHARACTER(len=50) :: FileN
38 INTEGER, PARAMETER :: nFile = 13
39 CHARACTER(len=50), DIMENSION(nFile) :: &
40 FileNames = [character(len=50) :: &
41 'SUEWS_NonVeg.txt', 'SUEWS_Veg.txt', 'SUEWS_Water.txt', 'SUEWS_Snow.txt', &
42 'SUEWS_Soil.txt', 'SUEWS_Conductance.txt', 'SUEWS_OHMCoefficients.txt', &
43 'SUEWS_ESTMCoefficients.txt', 'SUEWS_AnthropogenicEmission.txt', 'SUEWS_Irrigation.txt', &
44 'SUEWS_Profiles.txt', 'SUEWS_WithinGridWaterDist.txt', 'SUEWS_BiogenCO2.txt']
45
46 ! ---- Namelist for RunControl.nml ----
47 namelist /runcontrol/ filecode, &
50 tstep, &
61 cbluse, &
62 snowuse, &
63 ! SOLWEIGuse, &
64 diagmethod, &
70 smdmethod, &
73 ohmincqf, &
78 rainamongn, &
81 kdownzen, &
83 ! ncMode, &
84 ! nRow, &
85 ! nCol, &
86 diagnose, &
89 diagqn, &
90 diagqs
91
92 ! -------------------------------------
93
94 !Initialise namelist with default values
98 disaggmethod = 1 ! linear disaggregation of averages
99 disaggmethodestm = 1 ! linear disaggregation of averages
100 raindisaggmethod = 100 ! even distribution among all subintervals
101 rainamongn = -999 ! no default setting for number of rainy subintervals
102 multrainamongn = -999 ! no default setting for number of rainy subintervals
103 multrainamongnupperi = -999 ! no default setting for rain intensity upper bound
104 kdownzen = 1 ! use zenith angle by default
105
106 suppresswarnings = 0 ! write warnings file
107 resolutionfilesin = 0 ! Set to zero so that if not found, automatically set to Tstep below
108
109 ! Set Diagnose switch to off (0). If Diagnose = 1 is set in RunControl, model progress will be printed
110 diagnose = 0
113 diagqn = 0
114 diagqs = 0
115
116 filecode = 'none'
117 !smithFile='Smith1966.grd'
118
119 !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
120 !Read in the RunControl.nml file
121 OPEN (55, file='RunControl.nml', err=200, status='old') !Change with needs
122 READ (55, nml=runcontrol, err=201)
123 CLOSE (55)
124
125 ! set tstep_prev as tstep ! TS 04 Jul 2018
126 ! tstep_prev is NOT used by SUEWS but by WRF-SUEWS through the main interface
128
129 IF (diagnose == 1) WRITE (*, *) 'Diagnosis switched on (model progress will be printed to screen)...'
130
131 !Check for problems with FileCode
132 IF (filecode == 'none') CALL errorhint(26, trim("RunControl.nml FileCode is missing"), notused, notused, notusedi)
133
134 IF (resolutionfilesin == 0) resolutionfilesin = tstep !If ResolutionFilesIn not found, automatically set to Tstep
135
136 !-----------------------------------------------------------------------
137
138 !Write RunControl information to FileChoices.txt
139 filechoices = trim(fileoutputpath)//trim(filecode)//'_FileChoices.txt'
140 OPEN (12, file=filechoices, err=203)
141 WRITE (12, *) '----- RunControl -----'
142 WRITE (12, nml=runcontrol)
143 CLOSE (12)
144
145 ! !Determine what should be done with respect to radiation
146 ! ! TODO: this can be wrapped into a subroutine, TS 20 Oct 2017
147 ! AlbedoChoice=0
148 ! ldown_option=0
149 ! IF(NetRadiationMethod==0)THEN !Observed Q* from the met input file will be used
150 ! IF(SnowUse==1) THEN !If snow is modelled, NARP is needed for surface temperature
151 ! NetRadiationMethod=3000
152 ! ldown_option=3 !Ldown will be modelled
153 ! !NetRadiationMethod=NetRadiationMethod/1000
154 ! ENDIF
155 !
156 ! ELSEIF(NetRadiationMethod>0)THEN !Modelled Q* is used (NARP)
157 ! AlbedoChoice=-9
158 ! IF(NetRadiationMethod<10) THEN
159 ! AlbedoChoice=0
160 ! IF(NetRadiationMethod==1)ldown_option=1
161 ! IF(NetRadiationMethod==2)ldown_option=2
162 ! IF(NetRadiationMethod==3)ldown_option=3
163 !
164 ! ELSEIF(NetRadiationMethod>=100.AND.NetRadiationMethod<1000) THEN
165 ! AlbedoChoice=1
166 ! IF(NetRadiationMethod==100)ldown_option=1
167 ! IF(NetRadiationMethod==200)ldown_option=2
168 ! IF(NetRadiationMethod==300)ldown_option=3
169 ! NetRadiationMethod=NetRadiationMethod/100
170 ! ENDIF
171 !
172 ! !If bad NetRadiationMethod value
173 ! IF(NetRadiationMethod>3.OR. AlbedoChoice==-9)THEN
174 ! WRITE(*,*) 'NetRadiationMethod=',NetRadiationMethod
175 ! WRITE(*,*) 'Value not usable'
176 ! STOP
177 ! ENDIF
178 ! ENDIF
179
180 ! Adjust input for precip downscaling using different intensities (HCW 21 Apr 2017)
181 IF (raindisaggmethod == 102) THEN
182 DO i = 1, 5
184 END DO
185 END IF
186
187 !------------------------------------------------------------------
188 !Print run information on the screen
189 WRITE (*, *) '--------------------------------------------------------'
190 WRITE (*, *) "LUMPS/SUEWS - relevant references"
191 WRITE (*, *) "LUMPS - Grimmond and Oke (2002) JAM, 41, 79-810"
192 WRITE (*, *) "OHM - Grimmond and Oke (1999) JAM, 38, 922-940"
193 WRITE (*, *) "NARP - Offerle et al. (2003) JAM"
194 WRITE (*, *) "SUES - Evaporation Grimmond & Oke (1991) WRR"
195 WRITE (*, *) "Water Balance Model Grimmond et al. (1986) WRR"
196 WRITE (*, *) "NARP - Long wave improvements (Loridan et al. 2011 JAMC)"
197 WRITE (*, *) "SUEWS - Anthropogenic heat, etc (Jarvi et al. 2011 JH)"
198 WRITE (*, *) "SUEWS - Snow module included (Jarvi et al. 2014 GMD)"
199 WRITE (*, *) "SUEWS - v2016a release (Ward et al. 2016 UC)"
200 WRITE (*, *) '--------------------------------------------------------'
201
202 !=======================================================================
203 !======================== Read input files =============================
204 ! This part reads the input files derived from the SiteInfo spreadsheet
205
206 WRITE (*, *) 'Reading the following input files:'
207
208 !=======================SUEWS_SiteSelect.txt============================
209 filen = 'SUEWS_SiteSelect.txt'
210 CALL numberrows(filen, skipheadersiteinfo) !Find number of rows in input file
213 !Read input file
214 OPEN (21, file=trim(fileinputpath)//trim(filen), err=300, status='old')
215 DO skipcounter = 1, (skipheadersiteinfo - 1)
216 READ (21, *) !Skip lines before header
217 END DO
218 READ (21, *) (headersiteselect_file(iv), iv=1, ncolumnssiteselect) !Get header
219
220 DO i = 1, nlinessiteselect
221 READ (21, *) (siteselect(i, iv), iv=1, ncolumnssiteselect)
222 !write(*,*) (SiteSelect(i,iv),iv=1,ncolumnsSiteSelect)
223 END DO
224 CLOSE (21)
225
226 !call InputHeaderCheck(FileN) !! Need to add column checks for SiteSelect.txt
227
228 ! FileNames = (/'SUEWS_NonVeg.txt', 'SUEWS_Veg.txt', 'SUEWS_Water.txt', 'SUEWS_Snow.txt', &
229 ! 'SUEWS_Soil.txt', 'SUEWS_Conductance.txt', 'SUEWS_OHMCoefficients.txt', &
230 ! 'SUEWS_ESTMCoefficients.txt', 'SUEWS_AnthropogenicEmission.txt', 'SUEWS_Irrigation.txt', &
231 ! 'SUEWS_Profiles.txt', 'SUEWS_WithinGridWaterDist.txt', 'SUEWS_BiogenCO2.txt'/)
232
233 DO ifile = 1, nfile
234 CALL numberrows(filenames(ifile), skipheadersiteinfo) !Find number of rows in input file
235 SELECT CASE (ifile)
236 CASE (1)
240 CASE (2)
242 ALLOCATE (veg_coeff(nlinesveg, ncolumnsveg))
243 CALL readcoeff(filenames(ifile), nlinesveg, ncolumnsveg, headerveg_file, veg_coeff)
244 CASE (3)
248 CASE (4)
252 CASE (5)
256 CASE (6)
260 CASE (7)
265 CASE (8)
270 CASE (9)
274 CASE (10)
278 CASE (11)
282 CASE (12)
286 CASE (13)
290 END SELECT
291 END DO
292
293 !=======================================================================
294 !=======================================================================
295
296 !-----------------------------------------------------------------------
297 !SUEWS run information
298 inputmetformat = 10 !Input met data file in LUMPS format(1) or SUEWS format(10)
299 laicalcyes = 1 !Use observed(0) or modelled(1) LAI
300 evapmethod = 2 !Evaporation calculated according to Rutter(1) or Shuttleworth(2)
301 writedailystate = 1 !Daily state file written
302 tstepcount = 0
303
304 t_interval = 3600 !Number of seconds in an hour
305
306 !Calculate nsh (number of steps per hour) from model timestep (tstep) set in in RunControl
307 nsh_real = t_interval/real(tstep, kind(1d0))
308
309 ! Check nsh is an integer
310 IF (nsh_real == int(nsh_real)) THEN
311 nsh = int(nsh_real)
312 ELSE
313 CALL errorhint(39, &
314 'TSTEP must divide into t_INTERVAL exactly.', real(tstep, kind(1d0)), real(t_interval, kind(1d0)), notusedi)
315 END IF
316
317 ! Check nsh is reasonable
318 IF (nsh_real < 6 .OR. nsh_real > 60) THEN
319 CALL errorhint(39, 'TSTEP is too small or too large.', real(tstep, kind(1d0)), real(t_interval, kind(1d0)), notusedi)
320 END IF
321
322 ! Cast integer nsh as nsh_real for use in calculations
323 nsh_real = real(nsh, kind(1d0))
324 ! Cast integer tstep as tstep_real for use in calculations
325 tstep_real = real(tstep, kind(1d0))
326 ! get integer nsd from nsh for use in AnOHM checking, 20160708 TS
327 nsd = 24*nsh
328
329 !! Check this is still valid for v2016a
330 halftimestep = real(tstep_real)/2/(24*3600) !Used in NARP_cal_SunPosition to get sunpos in the middle of timestep
331
332 RETURN
333
334 !-------Possible problems-----------------------------------------------
335200 CALL errorhint(47, 'RunControl.nml', notused, notused, notusedi)
336201 CALL errorhint(48, 'RunControl.nml', notused, notused, notusedi)
337
338203 CALL errorhint(47, trim(filechoices), notused, notused, notusedi)
339
340300 CALL errorhint(48, trim(filen), notused, notused, notusedi)
341 !-----------------------------------------------------------------------
342
343 !pause
344
integer, parameter ncolumnssoil
character(len=20), dimension(ncolumnsnonveg) headernonveg_file
character(len=20), dimension(ncolumnswgwaterdist) headerwgwaterdist_file
integer, parameter ncolumnsconductance
integer, parameter ncolumnsohmcoefficients
character(len=20), dimension(ncolumnswater) headerwater_file
character(len=20), dimension(ncolumnssnow) headersnow_file
character(len=20), dimension(ncolumnsveg) headerveg_file
integer, parameter ncolumnsprofiles
integer, parameter ncolumnswater
integer, parameter ncolumnsanthropogenic
integer, parameter ncolumnsveg
integer, parameter ncolumnsnonveg
integer, parameter ncolumnswgwaterdist
character(len=20), dimension(ncolumnssoil) headersoil_file
character(len=20), dimension(ncolumnsanthropogenic) headeranthropogenic_file
character(len=20), dimension(ncolumnsohmcoefficients) headerohmcoefficients_file
integer, parameter ncolumnsestmcoefficients
character(len=20), dimension(ncolumnsprofiles) headerprofiles_file
character(len=20), dimension(ncolumnsirrigation) headerirrigation_file
character(len=20), dimension(ncolumnsbiogen) headerbiogen_file
character(len=20), dimension(ncolumnssiteselect) headersiteselect_file
character(len=20), dimension(ncolumnsconductance) headercond_file
integer, parameter ncolumnsbiogen
character(len=20), dimension(ncolumnsestmcoefficients) headerestmcoefficients_file
integer, parameter ncolumnsirrigation
integer, parameter ncolumnssnow
integer evapmethod
integer resolutionfilesin
integer resolutionfilesinestm
integer multiplemetfiles
integer ohmincqf
integer smdmethod
real(kind(1d0)), dimension(5) multrainamongnupperi
integer laicalcyes
integer diagnosedisagg
integer, dimension(5) multrainamongn
integer inputmetformat
integer keeptstepfilesin
integer disaggmethodestm
integer resolutionfilesout
integer suppresswarnings
integer multiplelayoutfiles
integer diagmethod
integer writeoutoption
integer skipheadersiteinfo
integer basetmethod
integer keeptstepfilesout
integer rainamongn
integer disaggmethod
integer multipleestmfiles
integer kdownzen
integer writedailystate
integer raindisaggmethod
integer diagnosedisaggestm
integer nlinesirrigation
integer nlinesohmcoefficients
integer nlinessiteselect
integer nlinesveg
integer nlinesprofiles
integer nlineswgwaterdist
integer nlineswater
integer nlinesestmcoefficients
integer nlinesconductance
integer nlinesbiogen
integer nlinessoil
integer nlinesanthropogenic
integer nlinessnow
integer nlinesnonveg
integer stabilitymethod
real(kind(1d0)) halftimestep
real(kind(1d0)) tstep_real
integer roughlenheatmethod
real(kind(1d0)) tstepcount
subroutine readcoeff(FileName, nlines, ncolumns, HeaderFile, Coeff)
subroutine numberrows(FileN, SkipHeaderLines)

References allocatearray::anthropogenic_coeff, data_in::basetmethod, allocatearray::biogen_coeff, data_in::cbluse, allocatearray::conductance_coeff, data_in::diagmethod, data_in::diagnose, data_in::diagnosedisagg, data_in::diagnosedisaggestm, data_in::diagqn, data_in::diagqs, data_in::disaggmethod, data_in::disaggmethodestm, data_in::emissionsmethod, errorhint(), allocatearray::estmcoefficients_coeff, data_in::evapmethod, data_in::filechoices, data_in::filecode, data_in::fileinputpath, data_in::fileoutputpath, sues_data::halftimestep, allocatearray::headeranthropogenic_file, allocatearray::headerbiogen_file, allocatearray::headercond_file, allocatearray::headerestmcoefficients_file, allocatearray::headerirrigation_file, allocatearray::headernonveg_file, allocatearray::headerohmcoefficients_file, allocatearray::headerprofiles_file, allocatearray::headersiteselect_file, allocatearray::headersnow_file, allocatearray::headersoil_file, allocatearray::headerveg_file, allocatearray::headerwater_file, allocatearray::headerwgwaterdist_file, data_in::inputmetformat, allocatearray::irrigation_coeff, data_in::kdownzen, data_in::keeptstepfilesin, data_in::keeptstepfilesout, data_in::laicalcyes, data_in::multipleestmfiles, data_in::multipleinitfiles, data_in::multiplelayoutfiles, data_in::multiplemetfiles, data_in::multrainamongn, data_in::multrainamongnupperi, allocatearray::ncolumnsanthropogenic, allocatearray::ncolumnsbiogen, allocatearray::ncolumnsconductance, allocatearray::ncolumnsestmcoefficients, allocatearray::ncolumnsirrigation, allocatearray::ncolumnsnonveg, allocatearray::ncolumnsohmcoefficients, allocatearray::ncolumnsprofiles, allocatearray::ncolumnssiteselect, allocatearray::ncolumnssnow, allocatearray::ncolumnssoil, allocatearray::ncolumnsveg, allocatearray::ncolumnswater, allocatearray::ncolumnswgwaterdist, data_in::netradiationmethod, initial::nlines, initial::nlinesanthropogenic, initial::nlinesbiogen, initial::nlinesconductance, initial::nlinesestmcoefficients, initial::nlinesirrigation, initial::nlinesnonveg, initial::nlinesohmcoefficients, initial::nlinesprofiles, initial::nlinessiteselect, initial::nlinessnow, initial::nlinessoil, initial::nlinesveg, initial::nlineswater, initial::nlineswgwaterdist, allocatearray::nonveg_coeff, defaultnotused::notused, defaultnotused::notusedi, sues_data::nsd, sues_data::nsh, sues_data::nsh_real, numberrows(), allocatearray::ohmcoefficients_coeff, data_in::ohmincqf, allocatearray::profiles_coeff, data_in::rainamongn, data_in::raindisaggmethod, readcoeff(), data_in::resolutionfilesin, data_in::resolutionfilesinestm, data_in::resolutionfilesout, sues_data::roughlenheatmethod, data_in::roughlenmommethod, allocatearray::siteselect, data_in::skipheadersiteinfo, data_in::smdmethod, allocatearray::snow_coeff, data_in::snowuse, allocatearray::soil_coeff, sues_data::stabilitymethod, data_in::storageheatmethod, data_in::suppresswarnings, sues_data::t_interval, sues_data::tstep, sues_data::tstep_prev, sues_data::tstep_real, time::tstepcount, allocatearray::veg_coeff, allocatearray::water_coeff, data_in::waterusemethod, allocatearray::wgwaterdist_coeff, data_in::writedailystate, and data_in::writeoutoption.

Referenced by suews_program().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ readcoeff()

subroutine readcoeff ( character(len=*), intent(in)  FileName,
integer, intent(in)  nlines,
integer, intent(in)  ncolumns,
character(len=*), dimension(ncolumns), intent(out)  HeaderFile,
real(kind(1d0)), dimension(nlines, ncolumns), intent(out)  Coeff 
)

Definition at line 348 of file suews_ctrl_init.f95.

349
350 USE data_in
352
353 IMPLICIT NONE
354 !----------------------------------------------------------------------
355 ! dummy arguments
356 !----------------------------------------------------------------------
357 CHARACTER(len=*), INTENT(in) :: FileName
358 INTEGER, INTENT(in) :: nlines, ncolumns
359 CHARACTER(len=*), INTENT(out) :: HeaderFile(ncolumns)
360 REAL(KIND(1D0)), INTENT(out) :: Coeff(nlines, ncolumns)
361
362 !----------------------------------------------------------------------
363 ! local variables
364 !----------------------------------------------------------------------
365 INTEGER :: SkipCounter, iv, i, ii
366
367 !Read input file
368 OPEN (22, file=trim(fileinputpath)//trim(filename), err=301, status='old')
369
370 DO skipcounter = 1, skipheadersiteinfo - 1
371 READ (22, *) !Skip lines before header
372 END DO
373 READ (22, *) (headerfile(iv), iv=1, ncolumns) !Get header
374
375 DO i = 1, nlines
376 READ (22, *) (coeff(i, iv), iv=1, ncolumns)
377 !write(*,*) (NonVeg_Coeff(i,iv),iv=1,ncolumnsNonVeg)
378 END DO
379 CLOSE (22)
380
382
383 ! Check codes are unique
384 DO i = 1, nlines
385 DO ii = i + 1, nlines
386 IF (coeff(i, 1) == coeff(ii, 1) .AND. i /= ii) THEN
387 WRITE (*, *) 'Code', coeff(i, 1), 'in ', trim(filename), ' not unique!'
388 CALL errorhint(60, filename, coeff(i, 1), notused, notusedi)
389 END IF
390 END DO
391 END DO
392
393 RETURN
394
395301 CALL errorhint(48, trim(filename), notused, notused, notusedi)
396
subroutine inputheadercheck(FileName)

References errorhint(), data_in::fileinputpath, inputheadercheck(), defaultnotused::notused, defaultnotused::notusedi, and data_in::skipheadersiteinfo.

Referenced by overallruncontrol().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ suews_initializemetdata()

subroutine suews_initializemetdata ( integer  lunit,
integer  ReadlinesMetdata_read 
)

Definition at line 2172 of file suews_ctrl_init.f95.

2173
2174 USE allocatearray
2175 USE data_in
2176 USE sues_data
2177 USE time
2178 USE defaultnotused
2179 USE initial
2180
2181 IMPLICIT NONE
2182
2183 INTEGER :: lunit, i, iyy !,RunNumber,NSHcounter
2184 INTEGER :: ReadlinesMetdata_read !,RunNumber,NSHcounter
2185 REAL(KIND(1D0)), DIMENSION(24) :: MetArray
2186 REAL(KIND(1D0)) :: imin_prev, ih_prev, iday_prev, tstep_met, iy_only !For checks on temporal resolution of met data
2187
2188 ! initialisation
2189 iy_only = 1
2190 ih_prev = 1
2191 imin_prev = 1
2192 iday_prev = 1
2193
2194 !---------------------------------------------------------------
2195
2196 !Open the file for reading and read the actual data
2197 WRITE (*, *) filemet
2198 OPEN (lunit, file=trim(filemet), status='old', err=314)
2199 CALL skipheader(lunit, skipheadermet)
2200
2201 ! Skip to the right place in the met file, depending on how many chunks have been read already
2202 IF (skippedlines > 0) THEN
2203 DO iyy = 1, skippedlines
2204 READ (lunit, *)
2205 END DO
2206 END IF
2207
2208 ! Read in next chunk of met data and fill MetForcingData array with data for every timestep
2209 !NSHcounter = 1
2210 !write(*,*) 'ReadlinesMetdata:',ReadlinesMetdata
2211 DO i = 1, readlinesmetdata_read
2212 CALL metread(lunit, metarray, inputmetformat, ldown_option, netradiationmethod, &
2214 !DO iv=1,NSH
2215 ! MetForcingData(NSHcounter,1:24,GridCounter) = MetArray
2216 ! NSHcounter = NSHcounter + 1
2217 !ENDDO
2218 metforcingdata(i, 1:24, gridcounter) = metarray
2219 ! Check timestamp of met data file matches TSTEP specified in RunControl
2220 IF (i == 1) THEN
2221 imin_prev = metarray(4)
2222 ih_prev = metarray(3)
2223 iday_prev = metarray(2)
2224 iy_only = metarray(1)
2225 ELSEIF (i == 2) THEN
2226 tstep_met = ((metarray(4) + 60*metarray(3)) - (imin_prev + 60*ih_prev))*60 !tstep in seconds
2227 IF (tstep_met /= tstep_real .AND. metarray(2) == iday_prev) THEN
2228 CALL errorhint(39, 'TSTEP in RunControl does not match TSTEP of met data (DOY).', real(tstep, kind(1d0)), tstep_met, &
2229 int(metarray(2)))
2230 END IF
2231 END IF
2232
2233 ! Check file only contains a single year --------------------------------------------
2234 ! Very last data point is allowed to be (should be) timestamped with following year
2235 IF (metarray(1) /= iy_only) THEN
2236 IF (metarray(1) == iy_only + 1 .AND. metarray(2) == 1 .AND. metarray(3) == 0 .AND. metarray(4) == 0) THEN
2237 !write(*,*) 'end of year - no problem'
2238 ELSE
2239 CALL errorhint(3, 'Problem in SUEWS_Initial: multiple years found in met forcing file.', &
2240 metarray(1), notused, notusedi)
2241 END IF
2242 END IF
2243
2244 END DO
2245
2246 CLOSE (lunit)
2247
2248 RETURN
2249
2250314 CALL errorhint(11, trim(filemet), notused, notused, ios_out)
2251
integer ldown_option
character(len=150) filemet
integer skipheadermet
integer gridcounter
integer skippedlines
real(kind(1d0)) soildensity
real(kind(1d0)) soilrocks
real(kind(1d0)) smcap
real(kind(1d0)) soildepthmeas
subroutine metread(lfn, MetArray, InputmetFormat, ldown_option, NetRadiationMethod, SnowUse, SMDMethod, SoilDepthMeas, SoilRocks, SoilDensity, SmCap)
subroutine skipheader(lfn, skip)

References errorhint(), data_in::filemet, initial::gridcounter, data_in::inputmetformat, defaultnotused::ios_out, data_in::ldown_option, allocatearray::metforcingdata, metread(), data_in::netradiationmethod, defaultnotused::notused, defaultnotused::notusedi, skipheader(), data_in::skipheadermet, initial::skippedlines, sues_data::smcap, data_in::smdmethod, data_in::snowuse, sues_data::soildensity, sues_data::soildepthmeas, sues_data::soilrocks, sues_data::tstep, and sues_data::tstep_real.

Referenced by suews_program().

Here is the call graph for this function:
Here is the caller graph for this function: