Line data Source code
1 : !> \file mo_restart.f90
2 : !> \brief \copybrief mo_restart
3 : !> \details \copydetails mo_restart
4 :
5 : !> \brief reading and writing states, fluxes and configuration for restart of mHM.
6 : !> \details routines are seperated for reading and writing variables for:
7 : !! - states and fluxes, and
8 : !! - configuration.
9 : !!
10 : !! Reading of L11 configuration is also seperated from the rest,
11 : !! since it is only required when routing is activated.
12 : !> \authors Stephan Thober
13 : !> \date Jul 2013
14 : !> \copyright Copyright 2005-\today, the mHM Developers, Luis Samaniego, Sabine Attinger: All rights reserved.
15 : !! mHM is released under the LGPLv3+ license \license_note
16 : !> \ingroup f_mhm
17 : MODULE mo_restart
18 :
19 : use mo_common_constants, only : soilHorizonsVarName, landCoverPeriodsVarName, LAIVarName
20 :
21 : IMPLICIT NONE
22 :
23 : PRIVATE
24 :
25 : PUBLIC :: read_restart_states ! read restart files for state variables from a given path
26 : PUBLIC :: write_restart_files ! write restart files for configuration to a given path
27 :
28 : CONTAINS
29 :
30 : !> \brief write restart files for each domain
31 : !> \details write restart files for each domain. For each domain
32 : !! three restart files are written. These are xxx_states.nc,
33 : !! xxx_L11_config.nc, and xxx_config.nc (xxx being the three digit
34 : !! domain index). If a variable is added here, it should also be added
35 : !! in the read restart routines below.
36 : !> \changelog
37 : !! - Stephan Thober Aug 2015
38 : !! - moved write of routing states to mRM
39 : !! - David Schaefer Nov 2015
40 : !! - mo_netcdf
41 : !! - Stephan Thober Nov 2016
42 : !! - moved processMatrix to common variables
43 : !! - Zink M. Demirel C. Mar 2017
44 : !! - Added Jarvis soil water stress function at SM process(3)
45 : !! - Robert Schweppe Feb 2018
46 : !! - Removed all L0 references
47 : !! - Robert Schweppe Jun 2018
48 : !! - refactoring and reformatting
49 : !! - Stephan Thober Dec 2022
50 : !! - added grid info for level0
51 : !> \authors Stephan Thober
52 : !> \date Jun 2014
53 19 : subroutine write_restart_files(OutFile)
54 :
55 : use mo_common_constants, only : nodata_dp
56 : use mo_common_restart, only : write_grid_info
57 : use mo_common_variables, only : level0, level1, nLCoverScene, domainMeta, LC_year_start, LC_year_end
58 : use mo_global_variables, only : L1_Inter, L1_Throughfall, L1_aETCanopy, L1_aETSealed, L1_aETSoil, L1_baseflow, &
59 : L1_fastRunoff, L1_infilSoil, L1_melt, L1_percol, L1_preEffect, L1_rain, &
60 : L1_runoffSeal, L1_satSTW, L1_sealSTW, L1_slowRunoff, L1_snow, L1_snowPack, &
61 : L1_soilMoist, L1_total_runoff, L1_unsatSTW
62 :
63 : use mo_kind, only : dp, i4
64 : use mo_message, only : message
65 : use mo_mpr_global_variables, only : nLAI, nSoilHorizons_mHM, HorizonDepth_mHM, LAIBoundaries
66 : use mo_mpr_restart, only : write_eff_params
67 : use mo_netcdf, only : NcDataset, NcDimension, NcVariable
68 : use mo_string_utils, only : num2str
69 :
70 : implicit none
71 :
72 : character(256) :: Fname
73 :
74 : ! Output Path for each domain
75 : character(256), dimension(:), intent(in) :: OutFile
76 :
77 : integer(i4) :: iDomain, domainID
78 :
79 : integer(i4) :: ii
80 :
81 : ! start index at level 1
82 : integer(i4) :: s1
83 :
84 : ! end index at level 1
85 : integer(i4) :: e1
86 :
87 : ! mask at level 1
88 9 : logical, dimension(:, :), allocatable :: mask1
89 :
90 : ! dummy variable
91 9 : real(dp), dimension(:, :, :), allocatable :: dummy_3D
92 9 : real(dp), dimension(:), allocatable :: dummy_1D
93 :
94 : integer(i4) :: max_extent
95 :
96 : type(NcDataset) :: nc
97 :
98 : type(NcDimension) :: rows1, cols1, soil1, lcscenes, lais
99 :
100 : type(NcVariable) :: var
101 :
102 :
103 : ! get maximum extent of one dimension 2 or 3
104 9 : max_extent = max(nSoilHorizons_mHM, nLCoverScene, nLAI)
105 :
106 25 : domain_loop : do iDomain = 1, domainMeta%nDomains
107 16 : domainID = domainMeta%indices(iDomain)
108 :
109 : ! write restart file for iDomain
110 16 : Fname = trim(OutFile(iDomain))
111 : ! print a message
112 16 : call message(" Writing Restart-file: ", trim(adjustl(Fname)), " ...")
113 :
114 16 : nc = NcDataset(fname, "w")
115 :
116 16 : call write_grid_info(level1(iDomain), "1", nc)
117 16 : call write_grid_info(level0(domainMeta%L0DataFrom(iDomain)), "0", nc)
118 :
119 16 : rows1 = nc%getDimension("nrows1")
120 16 : cols1 = nc%getDimension("ncols1")
121 :
122 : ! write the dimension to the file and also save bounds
123 48 : allocate(dummy_1D(nSoilHorizons_mHM+1))
124 16 : dummy_1D(1) = 0.0_dp
125 48 : dummy_1D(2:nSoilHorizons_mHM+1) = HorizonDepth_mHM(:)
126 16 : soil1 = nc%setCoordinate(trim(soilHorizonsVarName), nSoilHorizons_mHM, dummy_1D, 2_i4)
127 16 : deallocate(dummy_1D)
128 48 : allocate(dummy_1D(nLCoverScene+1))
129 48 : dummy_1D(1:nLCoverScene) = LC_year_start(:)
130 : ! this is done because bounds are always stored as real so e.g.
131 : ! 1981-1990,1991-2000 is thus saved as 1981.0-1991.0,1991.0-2001.0
132 : ! it is translated back into ints correctly during reading
133 16 : dummy_1D(nLCoverScene+1) = LC_year_end(nLCoverScene) + 1
134 16 : lcscenes = nc%setCoordinate(trim(landCoverPeriodsVarName), nLCoverScene, dummy_1D, 0_i4)
135 16 : deallocate(dummy_1D)
136 : ! write the dimension to the file
137 16 : lais = nc%setCoordinate(trim(LAIVarName), nLAI, LAIBoundaries, 0_i4)
138 :
139 : ! for appending and intialization
140 80 : allocate(dummy_3D(rows1%getLength(), cols1%getLength(), max_extent))
141 64 : allocate(mask1(rows1%getLength(), cols1%getLength()))
142 16 : s1 = level1(iDomain)%iStart
143 16 : e1 = level1(iDomain)%iEnd
144 1316 : mask1 = level1(iDomain)%mask
145 :
146 48 : var = nc%setVariable("L1_Inter", "f64", (/rows1, cols1/))
147 16 : call var%setFillValue(nodata_dp)
148 16 : call var%setData(unpack(L1_inter(s1 : e1), mask1, nodata_dp))
149 16 : call var%setAttribute("long_name", "Interception storage at level 1")
150 :
151 48 : var = nc%setVariable("L1_snowPack", "f64", (/rows1, cols1/))
152 16 : call var%setFillValue(nodata_dp)
153 16 : call var%setData(unpack(L1_snowPack(s1 : e1), mask1, nodata_dp))
154 16 : call var%setAttribute("long_name", "Snowpack at level 1")
155 :
156 48 : var = nc%setVariable("L1_sealSTW", "f64", (/rows1, cols1/))
157 16 : call var%setFillValue(nodata_dp)
158 16 : call var%setData(unpack(L1_sealSTW(s1 : e1), mask1, nodata_dp))
159 16 : call var%setAttribute("long_name", "Retention storage of impervious areas at level 1")
160 :
161 48 : do ii = 1, nSoilHorizons_mHM
162 48 : dummy_3D(:, :, ii) = unpack(L1_soilMoist(s1 : e1, ii), mask1, nodata_dp)
163 : end do
164 :
165 64 : var = nc%setVariable("L1_soilMoist", "f64", (/rows1, cols1, soil1/))
166 16 : call var%setFillValue(nodata_dp)
167 16 : call var%setData(dummy_3D(:, :, 1 : nSoilHorizons_mHM))
168 16 : call var%setAttribute("long_name", "soil moisture at level 1")
169 :
170 48 : var = nc%setVariable("L1_unsatSTW", "f64", (/rows1, cols1/))
171 16 : call var%setFillValue(nodata_dp)
172 16 : call var%setData(unpack(L1_unsatSTW(s1 : e1), mask1, nodata_dp))
173 16 : call var%setAttribute("long_name", "upper soil storage at level 1")
174 :
175 48 : var = nc%setVariable("L1_satSTW", "f64", (/rows1, cols1/))
176 16 : call var%setFillValue(nodata_dp)
177 16 : call var%setData(unpack(L1_satSTW(s1 : e1), mask1, nodata_dp))
178 16 : call var%setAttribute("long_name", "groundwater storage at level 1")
179 :
180 48 : do ii = 1, nSoilHorizons_mHM
181 48 : dummy_3D(:, :, ii) = unpack(L1_aETSoil(s1 : e1, ii), mask1, nodata_dp)
182 : end do
183 :
184 64 : var = nc%setVariable("L1_aETSoil", "f64", (/rows1, cols1, soil1/))
185 16 : call var%setFillValue(nodata_dp)
186 16 : call var%setData(dummy_3D(:, :, 1 : nSoilHorizons_mHM))
187 16 : call var%setAttribute("long_name", "soil actual ET at level 1")
188 :
189 48 : var = nc%setVariable("L1_aETCanopy", "f64", (/rows1, cols1/))
190 16 : call var%setFillValue(nodata_dp)
191 16 : call var%setData(unpack(L1_aETCanopy(s1 : e1), mask1, nodata_dp))
192 16 : call var%setAttribute("long_name", "canopy actual ET at level 1")
193 :
194 48 : var = nc%setVariable("L1_aETSealed", "f64", (/rows1, cols1/))
195 16 : call var%setFillValue(nodata_dp)
196 16 : call var%setData(unpack(L1_aETSealed(s1 : e1), mask1, nodata_dp))
197 16 : call var%setAttribute("long_name", "sealed actual ET at level 1")
198 :
199 48 : var = nc%setVariable("L1_baseflow", "f64", (/rows1, cols1/))
200 16 : call var%setFillValue(nodata_dp)
201 16 : call var%setData(unpack(L1_baseflow(s1 : e1), mask1, nodata_dp))
202 16 : call var%setAttribute("long_name", "baseflow at level 1")
203 :
204 48 : do ii = 1, nSoilHorizons_mHM
205 48 : dummy_3D(:, :, ii) = unpack(L1_infilSoil(s1 : e1, ii), mask1, nodata_dp)
206 : end do
207 :
208 64 : var = nc%setVariable("L1_infilSoil", "f64", (/rows1, cols1, soil1/))
209 16 : call var%setFillValue(nodata_dp)
210 16 : call var%setData(dummy_3D(:, :, 1 : nSoilHorizons_mHM))
211 16 : call var%setAttribute("long_name", "soil in-exfiltration at level 1")
212 :
213 48 : var = nc%setVariable("L1_fastRunoff", "f64", (/rows1, cols1/))
214 16 : call var%setFillValue(nodata_dp)
215 16 : call var%setData(unpack(L1_fastRunoff(s1 : e1), mask1, nodata_dp))
216 16 : call var%setAttribute("long_name", "fast runoff")
217 :
218 48 : var = nc%setVariable("L1_percol", "f64", (/rows1, cols1/))
219 16 : call var%setFillValue(nodata_dp)
220 16 : call var%setData(unpack(L1_percol(s1 : e1), mask1, nodata_dp))
221 16 : call var%setAttribute("long_name", "percolation at level 1")
222 :
223 48 : var = nc%setVariable("L1_melt", "f64", (/rows1, cols1/))
224 16 : call var%setFillValue(nodata_dp)
225 16 : call var%setData(unpack(L1_melt(s1 : e1), mask1, nodata_dp))
226 16 : call var%setAttribute("long_name", "snow melt at level 1")
227 :
228 48 : var = nc%setVariable("L1_preEffect", "f64", (/rows1, cols1/))
229 16 : call var%setFillValue(nodata_dp)
230 16 : call var%setData(unpack(L1_preEffect(s1 : e1), mask1, nodata_dp))
231 16 : call var%setAttribute("long_name", "effective precip. depth (snow melt + rain) at level 1")
232 :
233 48 : var = nc%setVariable("L1_rain", "f64", (/rows1, cols1/))
234 16 : call var%setFillValue(nodata_dp)
235 16 : call var%setData(unpack(L1_rain(s1 : e1), mask1, nodata_dp))
236 16 : call var%setAttribute("long_name", "rain (liquid water) at level 1")
237 :
238 48 : var = nc%setVariable("L1_runoffSeal", "f64", (/rows1, cols1/))
239 16 : call var%setFillValue(nodata_dp)
240 16 : call var%setData(unpack(L1_runoffSeal(s1 : e1), mask1, nodata_dp))
241 16 : call var%setAttribute("long_name", "runoff from impervious area at level 1")
242 :
243 48 : var = nc%setVariable("L1_slowRunoff", "f64", (/rows1, cols1/))
244 16 : call var%setFillValue(nodata_dp)
245 16 : call var%setData(unpack(L1_slowRunoff(s1 : e1), mask1, nodata_dp))
246 16 : call var%setAttribute("long_name", "slow runoff at level 1")
247 :
248 48 : var = nc%setVariable("L1_snow", "f64", (/rows1, cols1/))
249 16 : call var%setFillValue(nodata_dp)
250 16 : call var%setData(unpack(L1_snow(s1 : e1), mask1, nodata_dp))
251 16 : call var%setAttribute("long_name", "snow (solid water) at level 1")
252 :
253 48 : var = nc%setVariable("L1_Throughfall", "f64", (/rows1, cols1/))
254 16 : call var%setFillValue(nodata_dp)
255 16 : call var%setData(unpack(L1_Throughfall(s1 : e1), mask1, nodata_dp))
256 16 : call var%setAttribute("long_name", "throughfall at level 1")
257 :
258 48 : var = nc%setVariable("L1_total_runoff", "f64", (/rows1, cols1/))
259 16 : call var%setFillValue(nodata_dp)
260 16 : call var%setData(unpack(L1_total_runoff(s1 : e1), mask1, nodata_dp))
261 16 : call var%setAttribute("long_name", "total runoff at level 1")
262 :
263 16 : call write_eff_params(mask1, s1, e1, rows1, cols1, soil1, lcscenes, lais, nc)
264 :
265 16 : call nc%close()
266 :
267 57 : deallocate(dummy_3D, mask1)
268 : end do domain_loop
269 :
270 9 : end subroutine write_restart_files
271 :
272 :
273 : !> \brief reads fluxes and state variables from file
274 : !> \details read fluxes and state variables from given
275 : !! restart directory and initialises all state variables
276 : !! that are initialized in the subroutine initialise,
277 : !! contained in module mo_startup.
278 : !> \changelog
279 : !! - Stephan Thober Aug 2015
280 : !! - moved read of routing states to mRM
281 : !! - David Schaefer Nov 2015
282 : !! - mo_netcdf
283 : !! - Stephan Thober Nov 2016
284 : !! - moved processMatrix to common variables
285 : !! - Robert Schweppe Jun 2018
286 : !! - refactoring and reformatting
287 : !! - Sebastian Müller Mar 2023
288 : !! - compatibility layer for 2D/3D data
289 : !! - move reading of nLAI to mo_startup (needed beforehand)
290 : !> \authors Stephan Thober
291 : !> \date Apr 2013
292 1 : subroutine read_restart_states(iDomain, domainID, InFile)
293 :
294 9 : use mo_common_variables, only : LC_year_end, LC_year_start, level1, nLCoverScene, processMatrix
295 : use mo_global_variables, only : L1_Inter, L1_Throughfall, L1_aETCanopy, &
296 : L1_aETSealed, L1_aETSoil, L1_baseflow, L1_fastRunoff, L1_infilSoil, L1_melt, &
297 : L1_percol, L1_preEffect, L1_rain, L1_runoffSeal, L1_satSTW, L1_sealSTW, &
298 : L1_slowRunoff, L1_snow, L1_snowPack, L1_soilMoist, L1_total_runoff, L1_unsatSTW
299 : use mo_kind, only : dp, i4
300 : use mo_mpr_global_variables, only : L1_HarSamCoeff, &
301 : L1_PrieTayAlpha, L1_aeroResist, L1_alpha, L1_degDay, L1_degDayInc, L1_degDayMax, &
302 : L1_degDayNoPre, L1_fAsp, L1_fRoots, L1_fSealed, L1_jarvis_thresh_c1, &
303 : L1_kBaseFlow, L1_kPerco, L1_kSlowFlow, L1_karstLoss, L1_kfastFlow, L1_maxInter, &
304 : L1_petLAIcorFactor, L1_sealedThresh, L1_soilMoistExp, L1_soilMoistFC, &
305 : L1_soilMoistSat, L1_surfResist, L1_tempThresh, L1_unsatThresh, L1_wiltingPoint, &
306 : nLAI, nSoilHorizons_mHM, &
307 : ! neutron count
308 : L1_No_Count, L1_bulkDens, L1_latticeWater, L1_COSMICL3
309 :
310 : use mo_netcdf, only : NcDataset, NcDimension, NcVariable
311 : use mo_string_utils, only : num2str
312 : use mo_message, only: message, error_message
313 :
314 : implicit none
315 :
316 : !> number of domain
317 : integer(i4), intent(in) :: iDomain
318 : !> ID of domain
319 : integer(i4), intent(in) :: domainID
320 : !> Input Path including trailing slash
321 : character(256), intent(in) :: InFile
322 :
323 : character(256) :: Fname
324 : ! variable rank
325 : integer(i4) :: var_rank
326 : ! loop index
327 : integer(i4) :: ii, jj
328 : ! start index at level 1
329 : integer(i4) :: s1
330 : ! end index at level 1
331 : integer(i4) :: e1
332 : ! mask at level 1
333 1 : logical, dimension(:, :), allocatable :: mask1
334 : ! dummy, 2 dimension
335 1 : real(dp), dimension(:, :), allocatable :: dummyD2
336 : ! dummy, 3 dimension
337 1 : real(dp), dimension(:, :, :), allocatable :: dummyD3
338 : ! dummy, 4 dimension
339 1 : real(dp), dimension(:, :, :, :), allocatable :: dummyD4
340 :
341 : type(NcDataset) :: nc
342 : type(NcVariable) :: var
343 :
344 :
345 1 : Fname = trim(InFile)
346 1 : call message(' Reading states from ', trim(adjustl(Fname)),' ...')
347 :
348 : ! get domain information at level 1
349 4 : allocate(mask1 (level1(iDomain)%nrows, level1(iDomain)%ncols))
350 65 : mask1 = level1(iDomain)%mask
351 1 : s1 = level1(iDomain)%iStart
352 1 : e1 = level1(iDomain)%iEnd
353 :
354 1 : nc = NcDataset(fname, "r")
355 :
356 1 : if (nc%hasVariable('L1_Inter')) then
357 : !-------------------------------------------
358 : ! STATE VARIABLES (optionally)
359 : !-------------------------------------------
360 :
361 : ! Interception
362 1 : var = nc%getVariable("L1_Inter")
363 1 : call var%getData(dummyD2)
364 1 : L1_inter(s1 : e1) = pack(dummyD2, mask1)
365 :
366 : ! Snowpack
367 1 : var = nc%getVariable("L1_snowPack")
368 1 : call var%getData(dummyD2)
369 1 : L1_snowPack(s1 : e1) = pack(dummyD2, mask1)
370 :
371 : ! Retention storage of impervious areas
372 1 : var = nc%getVariable("L1_sealSTW")
373 1 : call var%getData(dummyD2)
374 1 : L1_sealSTW(s1 : e1) = pack(dummyD2, mask1)
375 :
376 : ! upper soil storage
377 1 : var = nc%getVariable("L1_unsatSTW")
378 1 : call var%getData(dummyD2)
379 1 : L1_unsatSTW(s1 : e1) = pack(dummyD2, mask1)
380 :
381 : ! groundwater storage
382 1 : var = nc%getVariable("L1_satSTW")
383 1 : call var%getData(dummyD2)
384 1 : L1_satSTW(s1 : e1) = pack(dummyD2, mask1)
385 :
386 : ! Soil moisture of each horizon
387 1 : var = nc%getVariable("L1_soilMoist")
388 1 : call var%getData(dummyD3)
389 3 : do ii = 1, nSoilHorizons_mHM
390 3 : L1_soilMoist(s1 : e1, ii) = pack(dummyD3(:, :, ii), mask1)
391 : end do
392 :
393 : !-------------------------------------------
394 : ! FLUXES (optionally)
395 : !-------------------------------------------
396 :
397 : ! soil actual ET
398 1 : var = nc%getVariable("L1_aETSoil")
399 1 : call var%getData(dummyD3)
400 3 : do ii = 1, nSoilHorizons_mHM
401 3 : L1_aETSoil(s1 : e1, ii) = pack(dummyD3(:, :, ii), mask1)
402 : end do
403 :
404 : ! canopy actual ET
405 1 : var = nc%getVariable("L1_aETCanopy")
406 1 : call var%getData(dummyD2)
407 1 : L1_aETCanopy(s1 : e1) = pack(dummyD2, mask1)
408 :
409 : ! sealed area actual ET
410 1 : var = nc%getVariable("L1_aETSealed")
411 1 : call var%getData(dummyD2)
412 1 : L1_aETSealed(s1 : e1) = pack(dummyD2, mask1)
413 :
414 : ! baseflow
415 1 : var = nc%getVariable("L1_baseflow")
416 1 : call var%getData(dummyD2)
417 1 : L1_baseflow(s1 : e1) = pack(dummyD2, mask1)
418 :
419 : ! soil in-exfiltration
420 1 : var = nc%getVariable("L1_infilSoil")
421 1 : call var%getData(dummyD3)
422 3 : do ii = 1, nSoilHorizons_mHM
423 3 : L1_infilSoil(s1 : e1, ii) = pack(dummyD3(:, :, ii), mask1)
424 : end do
425 :
426 : ! fast runoff
427 1 : var = nc%getVariable("L1_fastRunoff")
428 1 : call var%getData(dummyD2)
429 1 : L1_fastRunoff(s1 : e1) = pack(dummyD2, mask1)
430 :
431 : ! snow melt
432 1 : var = nc%getVariable("L1_melt")
433 1 : call var%getData(dummyD2)
434 1 : L1_melt(s1 : e1) = pack(dummyD2, mask1)
435 :
436 : ! percolation
437 1 : var = nc%getVariable("L1_percol")
438 1 : call var%getData(dummyD2)
439 1 : L1_percol(s1 : e1) = pack(dummyD2, mask1)
440 :
441 : ! effective precip. depth (snow melt + rain)
442 1 : var = nc%getVariable("L1_preEffect")
443 1 : call var%getData(dummyD2)
444 1 : L1_preEffect(s1 : e1) = pack(dummyD2, mask1)
445 :
446 : ! rain (liquid water)
447 1 : var = nc%getVariable("L1_rain")
448 1 : call var%getData(dummyD2)
449 1 : L1_rain(s1 : e1) = pack(dummyD2, mask1)
450 :
451 : ! runoff from impervious area
452 1 : var = nc%getVariable("L1_runoffSeal")
453 1 : call var%getData(dummyD2)
454 1 : L1_runoffSeal(s1 : e1) = pack(dummyD2, mask1)
455 :
456 : ! slow runoff
457 1 : var = nc%getVariable("L1_slowRunoff")
458 1 : call var%getData(dummyD2)
459 1 : L1_slowRunoff(s1 : e1) = pack(dummyD2, mask1)
460 :
461 : ! snow (solid water)
462 1 : var = nc%getVariable("L1_snow")
463 1 : call var%getData(dummyD2)
464 1 : L1_snow(s1 : e1) = pack(dummyD2, mask1)
465 :
466 : ! throughfall
467 1 : var = nc%getVariable("L1_Throughfall")
468 1 : call var%getData(dummyD2)
469 1 : L1_Throughfall(s1 : e1) = pack(dummyD2, mask1)
470 :
471 : ! total runoff
472 1 : var = nc%getVariable("L1_total_runoff")
473 1 : call var%getData(dummyD2)
474 1 : L1_total_runoff(s1 : e1) = pack(dummyD2, mask1)
475 : end if
476 :
477 : !-------------------------------------------
478 : ! EFFECTIVE PARAMETERS
479 : !-------------------------------------------
480 1 : var = nc%getVariable("L1_fSealed")
481 1 : call var%getData(dummyD3)
482 3 : do ii = 1, nLCoverScene
483 3 : L1_fSealed(s1 : e1, 1, ii) = pack(dummyD3(:, :, ii), mask1)
484 : end do
485 :
486 : ! exponent for the upper reservoir
487 1 : var = nc%getVariable("L1_alpha")
488 1 : var_rank = var%getrank()
489 0 : select case(var_rank)
490 : case(2)
491 : ! distribute over all land cover scenes
492 0 : call var%getData(dummyD2)
493 0 : do ii = 1, nLCoverScene
494 0 : L1_alpha(s1 : e1, 1, ii) = pack(dummyD2, mask1)
495 : end do
496 : case(3)
497 1 : call var%getData(dummyD3)
498 3 : do ii = 1, nLCoverScene
499 3 : L1_alpha(s1 : e1, 1, ii) = pack(dummyD3(:, :, ii), mask1)
500 : end do
501 : case default
502 1 : call error_message("Restart: L1_alpha rank needs to be 2 or 3")
503 : end select
504 :
505 : ! increase of the Degree-day factor per mm of increase in precipitation
506 1 : var = nc%getVariable("L1_degDayInc")
507 1 : call var%getData(dummyD3)
508 3 : do ii = 1, nLCoverScene
509 3 : L1_degDayInc(s1 : e1, 1, ii) = pack(dummyD3(:, :, ii), mask1)
510 : end do
511 :
512 : ! maximum degree-day factor
513 1 : var = nc%getVariable("L1_degDayMax")
514 1 : call var%getData(dummyD3)
515 3 : do ii = 1, nLCoverScene
516 3 : L1_degDayMax(s1 : e1, 1, ii) = pack(dummyD3(:, :, ii), mask1)
517 : end do
518 :
519 : ! degree-day factor with no precipitation
520 1 : var = nc%getVariable("L1_degDayNoPre")
521 1 : call var%getData(dummyD3)
522 3 : do ii = 1, nLCoverScene
523 3 : L1_degDayNoPre(s1 : e1, 1, ii) = pack(dummyD3(:, :, ii), mask1)
524 : end do
525 :
526 : ! degree-day factor
527 1 : var = nc%getVariable("L1_degDay")
528 1 : var_rank = var%getrank()
529 0 : select case(var_rank)
530 : case(2)
531 : ! distribute over all land cover scenes
532 0 : call var%getData(dummyD2)
533 0 : do ii = 1, nLCoverScene
534 0 : L1_degDay(s1 : e1, 1, ii) = pack(dummyD2, mask1)
535 : end do
536 : case(3)
537 1 : call var%getData(dummyD3)
538 3 : do ii = 1, nLCoverScene
539 3 : L1_degDay(s1 : e1, 1, ii) = pack(dummyD3(:, :, ii), mask1)
540 : end do
541 : case default
542 1 : call error_message("Restart: L1_degDay rank needs to be 2 or 3")
543 : end select
544 :
545 : ! Karstic percolation loss
546 1 : var = nc%getVariable("L1_karstLoss")
547 1 : call var%getData(dummyD2)
548 1 : L1_karstLoss(s1 : e1, 1, 1) = pack(dummyD2, mask1)
549 :
550 : ! Fraction of roots in soil horizons
551 1 : var = nc%getVariable("L1_fRoots")
552 1 : call var%getData(dummyD4)
553 3 : do jj = 1, nLCoverScene
554 7 : do ii = 1, nSoilHorizons_mHM
555 6 : L1_fRoots(s1 : e1, ii, jj) = pack(dummyD4(:, :, ii, jj), mask1)
556 : end do
557 : end do
558 :
559 : ! Maximum interception
560 1 : var = nc%getVariable("L1_maxInter")
561 1 : call var%getData(dummyD3)
562 13 : do ii = 1, nLAI
563 13 : L1_maxInter(s1 : e1, ii, 1) = pack(dummyD3(:, :, ii), mask1)
564 : end do
565 :
566 : ! fast interflow recession coefficient ("L1_kFastFlow" or "L1_kfastFlow")
567 1 : if (nc%hasvariable("L1_kFastFlow")) then
568 0 : var = nc%getVariable("L1_kFastFlow")
569 : else
570 1 : var = nc%getVariable("L1_kfastFlow")
571 : end if
572 1 : call var%getData(dummyD3)
573 3 : do ii = 1, nLCoverScene
574 3 : L1_kfastFlow(s1 : e1, 1, ii) = pack(dummyD3(:, :, ii), mask1)
575 : end do
576 :
577 : ! slow interflow recession coefficient
578 1 : var = nc%getVariable("L1_kSlowFlow")
579 1 : var_rank = var%getrank()
580 0 : select case(var_rank)
581 : case(2)
582 : ! distribute over all land cover scenes
583 0 : call var%getData(dummyD2)
584 0 : do ii = 1, nLCoverScene
585 0 : L1_kSlowFlow(s1 : e1, 1, ii) = pack(dummyD2, mask1)
586 : end do
587 : case(3)
588 1 : call var%getData(dummyD3)
589 3 : do ii = 1, nLCoverScene
590 3 : L1_kSlowFlow(s1 : e1, 1, ii) = pack(dummyD3(:, :, ii), mask1)
591 : end do
592 : case default
593 1 : call error_message("Restart: L1_kSlowFlow rank needs to be 2 or 3")
594 : end select
595 :
596 : ! baseflow recession coefficient
597 1 : var = nc%getVariable("L1_kBaseFlow")
598 1 : var_rank = var%getrank()
599 0 : select case(var_rank)
600 : case(2)
601 : ! distribute over all land cover scenes
602 0 : call var%getData(dummyD2)
603 0 : do ii = 1, nLCoverScene
604 0 : L1_kBaseFlow(s1 : e1, 1, ii) = pack(dummyD2, mask1)
605 : end do
606 : case(3)
607 1 : call var%getData(dummyD3)
608 3 : do ii = 1, nLCoverScene
609 3 : L1_kBaseFlow(s1 : e1, 1, ii) = pack(dummyD3(:, :, ii), mask1)
610 : end do
611 : case default
612 1 : call error_message("Restart: L1_kBaseFlow rank needs to be 2 or 3")
613 : end select
614 :
615 : ! percolation coefficient
616 1 : var = nc%getVariable("L1_kPerco")
617 1 : var_rank = var%getrank()
618 0 : select case(var_rank)
619 : case(2)
620 : ! distribute over all land cover scenes
621 0 : call var%getData(dummyD2)
622 0 : do ii = 1, nLCoverScene
623 0 : L1_kPerco(s1 : e1, 1, ii) = pack(dummyD2, mask1)
624 : end do
625 : case(3)
626 1 : call var%getData(dummyD3)
627 3 : do ii = 1, nLCoverScene
628 3 : L1_kPerco(s1 : e1, 1, ii) = pack(dummyD3(:, :, ii), mask1)
629 : end do
630 : case default
631 1 : call error_message("Restart: L1_kPerco rank needs to be 2 or 3")
632 : end select
633 :
634 : ! Soil moisture below which actual ET is reduced linearly till PWP
635 : ! for processCase(3) = 1
636 1 : var = nc%getVariable("L1_soilMoistFC")
637 1 : call var%getData(dummyD4)
638 3 : do jj = 1, nLCoverScene
639 7 : do ii = 1, nSoilHorizons_mHM
640 6 : L1_soilMoistFC(s1 : e1, ii, jj) = pack(dummyD4(:, :, ii, jj), mask1)
641 : end do
642 : end do
643 :
644 : ! Saturation soil moisture for each horizon [mm]
645 1 : var = nc%getVariable("L1_soilMoistSat")
646 1 : call var%getData(dummyD4)
647 3 : do jj = 1, nLCoverScene
648 7 : do ii = 1, nSoilHorizons_mHM
649 6 : L1_soilMoistSat(s1 : e1, ii, jj) = pack(dummyD4(:, :, ii, jj), mask1)
650 : end do
651 : end do
652 :
653 : ! Exponential parameter to how non-linear is the soil water retention
654 1 : var = nc%getVariable("L1_soilMoistExp")
655 1 : call var%getData(dummyD4)
656 3 : do jj = 1, nLCoverScene
657 7 : do ii = 1, nSoilHorizons_mHM
658 6 : L1_soilMoistExp(s1 : e1, ii, jj) = pack(dummyD4(:, :, ii, jj), mask1)
659 : end do
660 : end do
661 :
662 3 : if (any(processMatrix(3, 1) == (/2, 3/))) then
663 : ! jarvis critical value for normalized soil water content
664 0 : var = nc%getVariable("L1_jarvis_thresh_c1")
665 0 : call var%getData(dummyD2)
666 0 : L1_jarvis_thresh_c1(s1 : e1, 1, 1) = pack(dummyD2, mask1)
667 : end if
668 :
669 : ! Threshold temperature for snow/rain
670 1 : var = nc%getVariable("L1_tempThresh")
671 1 : call var%getData(dummyD3)
672 3 : do ii = 1, nLCoverScene
673 3 : L1_tempThresh(s1 : e1, 1, ii) = pack(dummyD3(:, :, ii), mask1)
674 : end do
675 :
676 : ! Threshold water depth controlling fast interflow
677 1 : var = nc%getVariable("L1_unsatThresh")
678 1 : call var%getData(dummyD2)
679 1 : L1_unsatThresh(s1 : e1, 1, 1) = pack(dummyD2, mask1)
680 :
681 : ! Threshold water depth for surface runoff in sealed surfaces
682 1 : var = nc%getVariable("L1_sealedThresh")
683 1 : call var%getData(dummyD2)
684 1 : L1_sealedThresh(s1 : e1, 1, 1) = pack(dummyD2, mask1)
685 :
686 : ! Permanent wilting point
687 1 : var = nc%getVariable("L1_wiltingPoint")
688 1 : call var%getData(dummyD4)
689 3 : do jj = 1, nLCoverScene
690 7 : do ii = 1, nSoilHorizons_mHM
691 6 : L1_wiltingPoint(s1 : e1, ii, jj) = pack(dummyD4(:, :, ii, jj), mask1)
692 : end do
693 : end do
694 :
695 : ! different parameters dependent on PET formulation
696 1 : select case (processMatrix(5, 1))
697 : case(-1) ! PET is input
698 :
699 : ! PET correction factor due to LAI
700 0 : var = nc%getVariable("L1_petLAIcorFactor")
701 0 : call var%getData(dummyD4)
702 0 : do jj = 1, nLCoverScene
703 0 : do ii = 1, nLAI
704 0 : L1_petLAIcorFactor(s1 : e1, ii, jj) = pack(dummyD4(:, :, ii, jj), mask1)
705 : end do
706 : end do
707 :
708 : case(0) ! PET is input
709 :
710 : ! PET correction factor due to terrain aspect
711 0 : var = nc%getVariable("L1_fAsp")
712 0 : call var%getData(dummyD2)
713 0 : L1_fAsp(s1 : e1, 1, 1) = pack(dummyD2, mask1)
714 :
715 : case(1) ! Hargreaves-Samani
716 :
717 : ! PET correction factor due to terrain aspect
718 1 : var = nc%getVariable("L1_fAsp")
719 1 : call var%getData(dummyD2)
720 1 : L1_fAsp(s1 : e1, 1, 1) = pack(dummyD2, mask1)
721 :
722 : ! Hargreaves Samani coeffiecient
723 1 : var = nc%getVariable("L1_HarSamCoeff")
724 1 : call var%getData(dummyD2)
725 1 : L1_HarSamCoeff(s1 : e1, 1, 1) = pack(dummyD2, mask1)
726 :
727 : case(2) ! Priestely-Taylor
728 :
729 : ! Priestley Taylor coeffiecient (alpha)
730 0 : var = nc%getVariable("L1_PrieTayAlpha")
731 0 : call var%getData(dummyD3)
732 0 : do ii = 1, nLAI
733 0 : L1_PrieTayAlpha(s1 : e1, ii, 1) = pack(dummyD3(:, :, ii), mask1)
734 : end do
735 :
736 : case(3) ! Penman-Monteith
737 :
738 : ! aerodynamical resitance
739 0 : var = nc%getVariable("L1_aeroResist")
740 0 : call var%getData(dummyD4)
741 0 : do jj = 1, nLCoverScene
742 0 : do ii = 1, nLAI
743 0 : L1_aeroResist(s1 : e1, ii, jj) = pack(dummyD4(:, :, ii, jj), mask1)
744 : end do
745 : end do
746 :
747 : ! bulk surface resitance
748 0 : var = nc%getVariable("L1_surfResist")
749 0 : call var%getData(dummyD3)
750 1 : do ii = 1, nLAI
751 0 : L1_surfResist(s1 : e1, ii, 1) = pack(dummyD3(:, :, ii), mask1)
752 : end do
753 :
754 : end select
755 :
756 :
757 : ! neutron count
758 1 : select case (processMatrix(10, 1))
759 : case(1) ! deslet
760 : ! N0 count
761 0 : var = nc%getVariable("L1_No_Count")
762 0 : call var%getData(dummyD2)
763 0 : L1_No_Count(s1:e1, 1, 1) = pack(dummyD2, mask1)
764 :
765 : ! Bulk density
766 0 : var = nc%getVariable("L1_bulkDens")
767 0 : call var%getData(dummyD4)
768 0 : do jj = 1, nLCoverScene
769 0 : do ii = 1, nSoilHorizons_mHM
770 0 : L1_bulkDens(s1:e1, ii, jj) = pack(dummyD4(:, :, ii, jj), mask1)
771 : end do
772 : end do
773 :
774 : ! Lattice water
775 0 : var = nc%getVariable("L1_latticeWater")
776 0 : call var%getData(dummyD4)
777 0 : do jj = 1, nLCoverScene
778 0 : do ii = 1, nSoilHorizons_mHM
779 0 : L1_latticeWater(s1:e1, ii, jj) = pack(dummyD4(:, :, ii, jj), mask1)
780 : end do
781 : end do
782 :
783 : case(2) ! COSMIC
784 : ! N0 count
785 0 : var = nc%getVariable("L1_No_Count")
786 0 : call var%getData(dummyD2)
787 0 : L1_No_Count(s1 : e1, 1, 1) = pack(dummyD2, mask1)
788 :
789 : ! Bulk density
790 0 : var = nc%getVariable("L1_bulkDens")
791 0 : call var%getData(dummyD4)
792 0 : do jj = 1, nLCoverScene
793 0 : do ii = 1, nSoilHorizons_mHM
794 0 : L1_bulkDens(s1:e1, ii, jj) = pack(dummyD4(:, :, ii, jj), mask1)
795 : end do
796 : end do
797 :
798 : ! Lattice water
799 0 : var = nc%getVariable("L1_latticeWater")
800 0 : call var%getData(dummyD4)
801 0 : do jj = 1, nLCoverScene
802 0 : do ii = 1, nSoilHorizons_mHM
803 0 : L1_latticeWater(s1: e1, ii, jj) = pack(dummyD4(:, :, ii, jj), mask1)
804 : end do
805 : end do
806 :
807 : ! COSMIC L3 parameter
808 0 : var = nc%getVariable("L1_COSMICL3")
809 0 : call var%getData(dummyD4)
810 1 : do jj = 1, nLCoverScene
811 0 : do ii = 1, nSoilHorizons_mHM
812 0 : L1_COSMICL3(s1:e1, ii, jj) = pack(dummyD4(:, :, ii, jj), mask1)
813 : end do
814 : end do
815 :
816 : end select
817 :
818 : ! close file
819 1 : call nc%close()
820 :
821 1 : end subroutine read_restart_states
822 :
823 : END MODULE mo_restart
|