Line data Source code
1 : !> \file mo_meteo_spatial_tools.f90
2 : !> \brief \copybrief mo_meteo_spatial_tools
3 : !> \details \copydetails mo_meteo_spatial_tools
4 :
5 : !> \brief Spatial aggegation or disaggregation of meteorological input data.
6 : !> \details This module contains two subroutines to upscale and downscale, respectively,
7 : !! the level-2 meterological inputs to a required Level-1 hydrological spatial resolution.
8 : !> \authors Rohini Kumar
9 : !> \date Jan 2013
10 : !> \copyright Copyright 2005-\today, the mHM Developers, Luis Samaniego, Sabine Attinger: All rights reserved.
11 : !! mHM is released under the LGPLv3+ license \license_note
12 : !> \ingroup f_meteo
13 : MODULE mo_meteo_spatial_tools
14 :
15 : ! This module provides routines for spatial aggegation or disaggregation of meteorological input data.
16 :
17 : USE mo_kind, ONLY : i4, dp
18 :
19 : IMPLICIT NONE
20 :
21 : PUBLIC :: spatial_aggregation ! Spatial aggregation or upscaling
22 : PUBLIC :: spatial_disaggregation ! Spatial disaggregation or downscaling
23 :
24 :
25 : ! ------------------------------------------------------------------
26 :
27 : ! NAME
28 : ! spatial_aggregation
29 :
30 : ! PURPOSE
31 : !> \brief Spatial aggregation of meterological variables
32 :
33 : !> \details Aggregate (or upscale) the given level-2 meteorological data to the
34 : !> required level-1 spatial resolution for the mHM run.
35 :
36 : ! HISTORY
37 : !> \authors Rohini Kumar
38 :
39 : !> \date Jan 2013
40 :
41 : ! Modifications:
42 : ! Rohini Kumar Nov 2013 - data1 changed from intent(inout) to intent(out)
43 : ! RK, MZ, DS May 2014 - added mask2
44 : ! Robert Schweppe Jun 2018 - refactoring and reformatting
45 :
46 :
47 : INTERFACE spatial_aggregation
48 : MODULE PROCEDURE spatial_aggregation_3d, spatial_aggregation_4d
49 : END INTERFACE spatial_aggregation
50 :
51 : ! ------------------------------------------------------------------
52 :
53 : ! NAME
54 : ! spatial_disaggregation
55 :
56 : ! PURPOSE
57 : !> \brief Spatial disaggregation of meterological variables
58 :
59 : !> \details Disaggregate (or downscale) the given level-2 meteorological data to the
60 : !> required level-1 spatial resolution for the mHM run.
61 :
62 : ! INTENT(IN)
63 : !> \param[in] "real(dp), dimension(:, :, :) :: data2" Level-2 data
64 : !> \param[in] "real(dp) :: cellsize2" Level-2 resolution
65 : !> \param[in] "real(dp) :: cellsize1" Level-1 resolution
66 : !> \param[in] "logical, dimension(:, :) :: mask1" Level-1 mask
67 : !> \param[in] "logical, dimension(:, :) :: mask2" Level-2 mask
68 :
69 : ! INTENT(OUT)
70 : !> \param[out] "real(dp), dimension(:, :, :) :: data1" Level-1 data
71 :
72 : ! HISTORY
73 : !> \authors Rohini Kumar
74 :
75 : !> \date Jan 2013
76 :
77 : ! Modifications:
78 : ! Rohini Kumar Nov 2013 - data1 changed from intent(inout) to intent(out)
79 : ! RK, MZ, DS May 2014 - added mask2
80 :
81 :
82 : INTERFACE spatial_disaggregation
83 : MODULE PROCEDURE spatial_disaggregation_3d, spatial_disaggregation_4d
84 : end INTERFACE spatial_disaggregation
85 :
86 : ! ------------------------------------------------------------------
87 :
88 : PRIVATE
89 :
90 : ! ------------------------------------------------------------------
91 :
92 : CONTAINS
93 :
94 0 : subroutine spatial_aggregation_3d(data2, cellsize2, cellsize1, mask1, mask2, data1)
95 :
96 : use mo_common_constants, only : nodata_dp
97 :
98 : implicit none
99 :
100 : ! Level-2 data
101 : real(dp), dimension(:, :, :), intent(in) :: data2
102 :
103 : ! Level-2 resolution
104 : real(dp), intent(in) :: cellsize2
105 :
106 : ! Level-1 resolution
107 : real(dp), intent(in) :: cellsize1
108 :
109 : ! Level-1 mask
110 : logical, dimension(:, :), intent(in) :: mask1
111 :
112 : ! Level-2 mask
113 : logical, dimension(:, :), intent(in) :: mask2
114 :
115 : ! Level-1 data
116 : real(dp), dimension(:, :, :), allocatable, intent(out) :: data1
117 :
118 : ! No. of rows and cols at Level-1
119 : integer(i4) :: nr2, nc2
120 :
121 : ! No. of rows and cols at Level-1
122 : integer(i4) :: nr1, nc1
123 :
124 0 : real(dp) :: cellFactor
125 :
126 0 : integer(i4), dimension(:, :), allocatable :: nTCells
127 :
128 : integer(i4) :: nTimeSteps
129 :
130 : integer(i4) :: i, j, ic, jc, t
131 :
132 :
133 : ! get number of rows and cols at level-2 from mask2
134 : ! and the total time steps
135 0 : nr2 = size(data2, 1)
136 0 : nc2 = size(data2, 2)
137 0 : nTimeSteps = size(data2, 3)
138 :
139 : ! get number of rows and cols at level-1 from mask1
140 0 : nr1 = size(mask1, 1)
141 0 : nc1 = size(mask1, 2)
142 :
143 : !-----------------------------------------------------------------------
144 : ! Allocate and initalize nTCells which comprises
145 : ! of number of L2 cells that belongs to a given L1 cell
146 : ! NOTE:: 1) cell size of L1 > L2 (see CellFactor)
147 : ! 2) nTCells is estimated over the valid masked domain only
148 : !-----------------------------------------------------------------------
149 :
150 : ! cellFactor = level-1 resolution (hydro) / level-2 resolution (meteo)
151 0 : cellFactor = cellsize1 / cellsize2
152 :
153 : ! nTCells calculations
154 0 : allocate(nTCells(nr1, nc1))
155 0 : nTCells(:, :) = 0
156 :
157 0 : do j = 1, nc2
158 0 : jc = ceiling(real(j, dp) / cellFactor)
159 0 : do i = 1, nr2
160 0 : ic = ceiling(real(i, dp) / cellFactor)
161 0 : if(.not. mask2(i, j)) cycle
162 0 : nTCells(ic, jc) = nTcells(ic, jc) + 1
163 : end do
164 : end do
165 :
166 :
167 : ! allocate and initalize L1_data
168 0 : allocate(data1(nr1, nc1, nTimeSteps))
169 0 : data1(:, :, :) = 0.0_dp
170 :
171 : ! time loop
172 0 : do t = 1, nTimeSteps
173 :
174 : ! perform spatial aggregation
175 0 : do j = 1, nc2
176 0 : jc = ceiling(real(j, dp) / cellFactor)
177 0 : do i = 1, nr2
178 0 : ic = ceiling(real(i, dp) / cellFactor)
179 :
180 : ! only in valid masked area
181 0 : if(.not. mask2(i, j)) cycle
182 0 : data1(ic, jc, t) = data1(ic, jc, t) + data2(i, j, t)
183 :
184 : end do
185 : end do
186 :
187 : ! perform spatial average only over valid masked domain
188 : ! out of the masked domain nTCells(:,:) = 0
189 0 : where(mask1)
190 0 : data1(:, :, t) = data1(:, :, t) / real(nTcells(:, :), dp)
191 : elsewhere
192 0 : data1(:, :, t) = nodata_dp
193 : endwhere
194 :
195 : end do
196 :
197 : ! free space
198 0 : deallocate(nTCells)
199 :
200 0 : end subroutine spatial_aggregation_3d
201 :
202 0 : subroutine spatial_aggregation_4d(data2, cellsize2, cellsize1, mask1, mask2, data1)
203 :
204 0 : use mo_common_constants, only : nodata_dp
205 :
206 : implicit none
207 :
208 : ! Level-2 data
209 : real(dp), dimension(:, :, :, :), intent(in) :: data2
210 :
211 : ! Level-2 resolution
212 : real(dp), intent(in) :: cellsize2
213 :
214 : ! Level-1 resolution
215 : real(dp), intent(in) :: cellsize1
216 :
217 : ! Level-1 mask
218 : logical, dimension(:, :), intent(in) :: mask1
219 :
220 : ! Level-2 mask
221 : logical, dimension(:, :), intent(in) :: mask2
222 :
223 : ! Level-1 data
224 : real(dp), dimension(:, :, :, :), allocatable, intent(out) :: data1
225 :
226 : ! No. of rows and cols at Level-1
227 : integer(i4) :: nr2, nc2
228 :
229 : ! No. of rows and cols at Level-1
230 : integer(i4) :: nr1, nc1
231 :
232 0 : real(dp) :: cellFactor
233 :
234 0 : integer(i4), dimension(:, :), allocatable :: nTCells
235 :
236 : integer(i4) :: nMonths, nHours
237 :
238 : integer(i4) :: i, j, ic, jc, t, h
239 :
240 :
241 : ! get number of rows and cols at level-2 from mask2
242 : ! and the total time steps
243 0 : nr2 = size(data2, 1)
244 0 : nc2 = size(data2, 2)
245 0 : nMonths = size(data2, 3)
246 0 : nHours = size(data2, 4)
247 :
248 : ! get number of rows and cols at level-1 from mask1
249 0 : nr1 = size(mask1, 1)
250 0 : nc1 = size(mask1, 2)
251 :
252 : !-----------------------------------------------------------------------
253 : ! Allocate and initalize nTCells which comprises
254 : ! of number of L2 cells that belongs to a given L1 cell
255 : ! NOTE:: 1) cell size of L1 > L2 (see CellFactor)
256 : ! 2) nTCells is estimated over the valid masked domain only
257 : !-----------------------------------------------------------------------
258 :
259 : ! cellFactor = level-1 resolution (hydro) / level-2 resolution (meteo)
260 0 : cellFactor = cellsize1 / cellsize2
261 :
262 : ! nTCells calculations
263 0 : allocate(nTCells(nr1, nc1))
264 0 : nTCells(:, :) = 0
265 :
266 0 : do j = 1, nc2
267 0 : jc = ceiling(real(j, dp) / cellFactor)
268 0 : do i = 1, nr2
269 0 : ic = ceiling(real(i, dp) / cellFactor)
270 0 : if(.not. mask2(i, j)) cycle
271 0 : nTCells(ic, jc) = nTcells(ic, jc) + 1
272 : end do
273 : end do
274 :
275 :
276 : ! allocate and initalize L1_data
277 0 : allocate(data1(nr1, nc1, nMonths, nHours))
278 0 : data1(:, :, :, :) = 0.0_dp
279 :
280 : ! time loop
281 0 : do t = 1, nMonths
282 0 : do h = 1, nHours
283 :
284 : ! perform spatial aggregation
285 0 : do j = 1, nc2
286 0 : jc = ceiling(real(j, dp) / cellFactor)
287 0 : do i = 1, nr2
288 0 : ic = ceiling(real(i, dp) / cellFactor)
289 :
290 : ! only in valid masked area
291 0 : if(.not. mask2(i, j)) cycle
292 0 : data1(ic, jc, t, h) = data1(ic, jc, t, h) + data2(i, j, t, h)
293 :
294 : end do
295 : end do
296 :
297 : ! perform spatial average only over valid masked domain
298 : ! out of the masked domain nTCells(:,:) = 0
299 0 : where(mask1)
300 0 : data1(:, :, t, h) = data1(:, :, t, h) / real(nTcells(:, :), dp)
301 : elsewhere
302 0 : data1(:, :, t, h) = nodata_dp
303 : endwhere
304 :
305 : end do
306 : end do
307 :
308 : ! free space
309 0 : deallocate(nTCells)
310 :
311 0 : end subroutine spatial_aggregation_4d
312 :
313 18 : subroutine spatial_disaggregation_3d(data2, cellsize2, cellsize1, mask1, mask2, data1)
314 :
315 0 : use mo_common_constants, only : nodata_dp
316 :
317 : implicit none
318 :
319 : ! Level-2 data
320 : real(dp), dimension(:, :, :), intent(in) :: data2
321 :
322 : ! Level-2 resolution
323 : real(dp), intent(in) :: cellsize2
324 :
325 : ! Level-1 resolution
326 : real(dp), intent(in) :: cellsize1
327 :
328 : ! Level-1 mask
329 : logical, dimension(:, :), intent(in) :: mask1
330 :
331 : ! Level-2 mask
332 : logical, dimension(:, :), intent(in) :: mask2
333 :
334 : ! Level-1 data
335 : real(dp), dimension(:, :, :), allocatable, intent(out) :: data1
336 :
337 : ! No. of rows and cols at Level-1
338 : integer(i4) :: nr1, nc1
339 :
340 9 : real(dp) :: cellFactor
341 :
342 : integer(i4) :: nTimeSteps
343 :
344 : integer(i4) :: i, j, t, ic, jc
345 :
346 :
347 : ! get number of rows and cols at level-2 from mask2
348 9 : nr1 = size(mask1, 1)
349 9 : nc1 = size(mask1, 2)
350 :
351 : ! cellFactor = level-2 resolution (meteo) / level-1 resolution (hydro)
352 9 : cellFactor = cellsize2 / cellsize1
353 :
354 : ! total time steps
355 9 : nTimeSteps = size(data2, 3)
356 :
357 : ! allocate and initalize L1_data
358 45 : allocate(data1(nr1, nc1, nTimeSteps))
359 1320474 : data1(:, :, :) = nodata_dp
360 :
361 : ! over the time loop
362 5628 : do t = 1, nTimeSteps
363 :
364 : ! spatial disaggregation
365 106770 : do j = 1, nc1
366 101142 : jc = ceiling(real(j, dp) / cellFactor)
367 1320465 : do i = 1, nr1
368 1213704 : ic = ceiling(real(i, dp) / cellFactor)
369 : ! only over the valid masked area
370 1213704 : if(.not. mask2(ic, jc)) cycle
371 1314846 : data1(i, j, t) = data2(ic, jc, t)
372 : end do
373 : end do
374 :
375 : end do
376 :
377 9 : end subroutine spatial_disaggregation_3d
378 :
379 0 : subroutine spatial_disaggregation_4d(data2, cellsize2, cellsize1, mask1, mask2, data1)
380 :
381 9 : use mo_common_constants, only : nodata_dp
382 :
383 : implicit none
384 :
385 : ! Level-2 data
386 : real(dp), dimension(:, :, :, :), intent(in) :: data2
387 :
388 : ! Level-2 resolution
389 : real(dp), intent(in) :: cellsize2
390 :
391 : ! Level-1 resolution
392 : real(dp), intent(in) :: cellsize1
393 :
394 : ! Level-1 mask
395 : logical, dimension(:, :), intent(in) :: mask1
396 :
397 : ! Level-2 mask
398 : logical, dimension(:, :), intent(in) :: mask2
399 :
400 : ! Level-1 data
401 : real(dp), dimension(:, :, :, :), allocatable, intent(out) :: data1
402 :
403 : ! No. of rows and cols at Level-1
404 : integer(i4) :: nr1, nc1
405 :
406 0 : real(dp) :: cellFactor
407 :
408 : integer(i4) :: nMonths, nHours
409 :
410 : integer(i4) :: i, j, t, ic, jc, h
411 :
412 :
413 : ! get number of rows and cols at level-2 from mask2
414 0 : nr1 = size(mask1, 1)
415 0 : nc1 = size(mask1, 2)
416 :
417 : ! cellFactor = level-2 resolution (meteo) / level-1 resolution (hydro)
418 0 : cellFactor = cellsize2 / cellsize1
419 :
420 : ! time axis
421 0 : nMonths = size(data2, 3)
422 0 : nHours = size(data2, 4)
423 :
424 : ! allocate and initalize L1_data
425 0 : allocate(data1(nr1, nc1, nMonths, nHours))
426 0 : data1(:, :, :, :) = nodata_dp
427 :
428 : ! over the time loop
429 0 : do t = 1, nMonths
430 0 : do h = 1, nHours
431 :
432 : ! spatial disaggregation
433 0 : do j = 1, nc1
434 0 : jc = ceiling(real(j, dp) / cellFactor)
435 0 : do i = 1, nr1
436 0 : ic = ceiling(real(i, dp) / cellFactor)
437 : ! only over the valid masked area
438 0 : if(.not. mask2(ic, jc)) cycle
439 0 : data1(i, j, t, h) = data2(ic, jc, t, h)
440 : end do
441 : end do
442 : end do
443 : end do
444 :
445 0 : end subroutine spatial_disaggregation_4d
446 :
447 : END MODULE mo_meteo_spatial_tools
|