Line data Source code
1 : !> \file mo_mrm_read_data.f90
2 : !> \brief \copybrief mo_mrm_read_data
3 : !> \details \copydetails mo_mrm_read_data
4 :
5 : !> \brief mRM reading routines
6 : !> \details This module contains all routines to read mRM data from file.
7 : !> \authors Stephan Thober
8 : !> \date Aug 2015
9 : !> \copyright Copyright 2005-\today, the mHM Developers, Luis Samaniego, Sabine Attinger: All rights reserved.
10 : !! mHM is released under the LGPLv3+ license \license_note
11 : !> \ingroup f_mrm
12 : module mo_mrm_read_data
13 : use mo_kind, only : i4, dp
14 : use mo_message, only: message, error_message
15 :
16 : implicit none
17 :
18 : public :: mrm_read_L0_data
19 : public :: mrm_read_discharge
20 : public :: mrm_read_total_runoff
21 : public :: mrm_read_bankfull_runoff
22 : private
23 : contains
24 : ! ------------------------------------------------------------------
25 :
26 : ! NAME
27 : ! mrm_read_L0_data
28 :
29 : ! PURPOSE
30 : !> \brief read L0 data from file
31 :
32 : !> \details With the exception of L0_mask, L0_elev, and L0_LCover, all
33 : !> L0 variables are read from file. The former three are only read if they
34 : !> are not provided as variables.
35 :
36 : ! INTENT(IN)
37 : !> \param[in] "logical :: do_reinit"
38 : !> \param[in] "logical :: do_readlatlon"
39 : !> \param[in] "logical :: do_readlcover"
40 :
41 : ! HISTORY
42 : !> \authors Juliane Mai, Matthias Zink, and Stephan Thober
43 :
44 : !> \date Aug 2015
45 :
46 : ! Modifications:
47 : ! Stephan Thober Sep 2015 - added L0_mask, L0_elev, and L0_LCover
48 : ! Robert Schweppe Jun 2018 - refactoring and reformatting
49 : ! Stephan Thober Jun 2018 - including varying celerity functionality
50 :
51 12 : subroutine mrm_read_L0_data(do_reinit, do_readlatlon, do_readlcover)
52 :
53 : use mo_append, only : append
54 : use mo_common_constants, only : nodata_i4
55 : use mo_common_read_data, only : read_dem, read_lcover
56 : use mo_common_types, only: Grid
57 : use mo_common_variables, only : L0_LCover, dirMorpho, level0, domainMeta, processMatrix
58 : use mo_mpr_file, only: file_slope, uslope
59 : use mo_mrm_file, only : file_facc, file_fdir, &
60 : file_gaugeloc, ufacc, ufdir, ugaugeloc
61 : use mo_mrm_global_variables, only : L0_InflowGaugeLoc, L0_fAcc, L0_fDir, L0_gaugeLoc, domain_mrm
62 : use mo_read_latlon, only : read_latlon
63 : use mo_read_spatial_data, only : read_spatial_data_ascii
64 : use mo_string_utils, only : num2str
65 :
66 : implicit none
67 :
68 : logical, intent(in) :: do_reinit
69 :
70 : logical, intent(in) :: do_readlatlon
71 :
72 : logical, intent(in) :: do_readlcover
73 :
74 : integer(i4) ::domainID, iDomain
75 :
76 : integer(i4) :: iVar
77 :
78 : integer(i4) :: iGauge
79 :
80 : character(256) :: fname
81 :
82 : integer(i4) :: nunit
83 :
84 12 : integer(i4), dimension(:, :), allocatable :: data_i4_2d
85 :
86 12 : integer(i4), dimension(:, :), allocatable :: dataMatrix_i4
87 :
88 12 : logical, dimension(:, :), allocatable :: mask_2d
89 :
90 12 : logical, dimension(:, :), allocatable :: mask_global
91 :
92 : type(Grid), pointer :: level0_iDomain => null()
93 :
94 :
95 : ! ************************************************
96 : ! READ SPATIAL DATA FOR EACH DOMAIN
97 : ! ************************************************
98 12 : if (do_reinit) then
99 0 : call read_dem()
100 : end if
101 :
102 12 : if (do_readlcover .and. processMatrix(8, 1) .eq. 1) then
103 9 : call read_lcover()
104 3 : else if (do_readlcover .and. ((processMatrix(8, 1) .eq. 2) .or. (processMatrix(8, 1) .eq. 3))) then
105 0 : allocate(dataMatrix_i4(count(mask_global), 1))
106 0 : dataMatrix_i4 = nodata_i4
107 0 : call append(L0_LCover, dataMatrix_i4)
108 : ! free memory
109 0 : deallocate(dataMatrix_i4)
110 : end if
111 :
112 36 : do iDomain = 1, domainMeta%nDomains
113 24 : domainID = domainMeta%indices(iDomain)
114 :
115 24 : level0_iDomain => level0(domainMeta%L0DataFrom(iDomain))
116 :
117 : ! ToDo: check if change is correct
118 : ! check whether L0 data is shared
119 24 : if (domainMeta%L0DataFrom(iDomain) < iDomain) then
120 : !
121 : call message(' Using data of domain ', &
122 3 : trim(adjustl(num2str(domainMeta%indices(domainMeta%L0DataFrom(iDomain))))), ' for domain: ',&
123 3 : trim(adjustl(num2str(domainID))), '...')
124 3 : cycle
125 : !
126 : end if
127 : !
128 21 : call message(' Reading data for domain: ', trim(adjustl(num2str(domainID))), ' ...')
129 :
130 25 : if (do_readlatlon) then
131 : ! read lat lon coordinates of each domain
132 4 : call read_latlon(iDomain, "lon_l0", "lat_l0", "level0", level0_iDomain)
133 : end if
134 :
135 : ! read fAcc, fDir, gaugeLoc
136 117 : do iVar = 1, 4
137 21 : select case (iVar)
138 : case(1) ! flow accumulation
139 21 : fName = trim(adjustl(dirMorpho(iDomain))) // trim(adjustl(file_facc))
140 21 : nunit = ufacc
141 : case(2) ! flow direction
142 21 : fName = trim(adjustl(dirMorpho(iDomain))) // trim(adjustl(file_fdir))
143 21 : nunit = ufdir
144 : case(3) ! location of gauging stations
145 21 : fName = trim(adjustl(dirMorpho(iDomain))) // trim(adjustl(file_gaugeloc))
146 21 : nunit = ugaugeloc
147 : case(4)
148 21 : fName = trim(adjustl(dirMorpho(iDomain))) // trim(adjustl(file_slope))
149 105 : nunit = uslope
150 : end select
151 :
152 84 : if (iVar .le. 2) then
153 : !
154 : ! reading and transposing
155 : call read_spatial_data_ascii(trim(fName), nunit, &
156 : level0_iDomain%nrows, level0_iDomain%ncols, level0_iDomain%xllcorner, &
157 42 : level0_iDomain%yllcorner, level0_iDomain%cellsize, data_i4_2d, mask_2d)
158 :
159 : ! put global nodata value into array (probably not all grid cells have values)
160 4975668 : data_i4_2d = merge(data_i4_2d, nodata_i4, mask_2d)
161 : end if
162 :
163 84 : if (iVar .eq. 3) then
164 21 : if (domain_mrm(iDomain)%nGauges .ge. 1_i4) then
165 : !
166 : ! reading and transposing
167 : call read_spatial_data_ascii(trim(fName), nunit, &
168 : level0_iDomain%nrows, level0_iDomain%ncols, level0_iDomain%xllcorner, &
169 15 : level0_iDomain%yllcorner, level0_iDomain%cellsize, data_i4_2d, mask_2d)
170 :
171 : ! put global nodata value into array (probably not all grid cells have values)
172 1805742 : data_i4_2d = merge(data_i4_2d, nodata_i4, mask_2d)
173 :
174 : else
175 : ! set to nodata, but not ommit because data association between arrays and domains might break
176 682092 : data_i4_2d = merge(nodata_i4, nodata_i4, level0_iDomain%mask)
177 : end if
178 : end if
179 :
180 : ! put data into global L0 variable
181 : select case (iVar)
182 : case(1) ! flow accumulation
183 21 : call append(L0_fAcc, pack(data_i4_2d, level0_iDomain%mask))
184 : case(2) ! flow direction
185 : ! rotate flow direction and any other variable with directions
186 : ! NOTE: ONLY when ASCII files are read
187 21 : call rotate_fdir_variable(data_i4_2d)
188 : ! append
189 21 : call append(L0_fDir, pack(data_i4_2d, level0_iDomain%mask))
190 : case(3) ! location of evaluation and inflow gauging stations
191 : ! evaluation gauges
192 : ! Input data check
193 37 : do iGauge = 1, domain_mrm(iDomain)%nGauges
194 : ! If gaugeId is found in gauging location file?
195 209111 : if (.not. any(data_i4_2d .EQ. domain_mrm(iDomain)%gaugeIdList(iGauge))) then
196 : call error_message('***ERROR: Gauge ID "', trim(adjustl(num2str(domain_mrm(iDomain)%gaugeIdList(iGauge)))), &
197 0 : '" not found in ', raise=.false.)
198 : call error_message(' Gauge location input file: ', &
199 0 : trim(adjustl(dirMorpho(iDomain))) // trim(adjustl(file_gaugeloc)))
200 : end if
201 : end do
202 :
203 21 : call append(L0_gaugeLoc, pack(data_i4_2d, level0_iDomain%mask))
204 :
205 : ! inflow gauges
206 : ! if no inflow gauge for this subdomain exists still matirx with dim of subdomain has to be paded
207 21 : if (domain_mrm(iDomain)%nInflowGauges .GT. 0_i4) then
208 : ! Input data check
209 2 : do iGauge = 1, domain_mrm(iDomain)%nInflowGauges
210 : ! If InflowGaugeId is found in gauging location file?
211 55318 : if (.not. any(data_i4_2d .EQ. domain_mrm(iDomain)%InflowGaugeIdList(iGauge))) then
212 : call error_message('***ERROR: Inflow Gauge ID "', &
213 : trim(adjustl(num2str(domain_mrm(iDomain)%InflowGaugeIdList(iGauge)))), &
214 0 : '" not found in ', raise=.false.)
215 : call error_message(' Gauge location input file: ', &
216 0 : trim(adjustl(dirMorpho(iDomain))) // trim(adjustl(file_gaugeloc)))
217 : end if
218 : end do
219 : end if
220 :
221 105 : call append(L0_InflowGaugeLoc, pack(data_i4_2d, level0_iDomain%mask))
222 :
223 : end select
224 : !
225 : ! deallocate arrays
226 84 : if (allocated(data_i4_2d)) deallocate(data_i4_2d)
227 105 : if (allocated(mask_2d)) deallocate(mask_2d)
228 : !
229 : end do
230 : end do
231 :
232 12 : end subroutine mrm_read_L0_data
233 : ! ---------------------------------------------------------------------------
234 :
235 : ! NAME
236 : ! mrm_read_discharge
237 :
238 : ! PURPOSE
239 : !> \brief Read discharge timeseries from file
240 :
241 : !> \details Read Observed discharge at the outlet of a catchment
242 : !> and at the inflow of a catchment. Allocate global runoff
243 : !> variable that contains the simulated runoff after the simulation.
244 :
245 : ! HISTORY
246 : !> \authors Matthias Zink & Stephan Thober
247 :
248 : !> \date Aug 2015
249 :
250 : ! Modifications:
251 : ! Robert Schweppe Jun 2018 - refactoring and reformatting
252 :
253 13 : subroutine mrm_read_discharge
254 :
255 12 : use mo_append, only : paste
256 : use mo_common_constants, only : nodata_dp
257 : use mo_common_mHM_mRM_variables, only : evalPer, nTstepDay, opti_function, optimize, simPer
258 : use mo_common_variables, only : domainMeta
259 : use mo_mrm_file, only : udischarge
260 : use mo_mrm_global_variables, only : InflowGauge, gauge, mRM_runoff, nGaugesLocal, &
261 : nInflowGaugesTotal, nMeasPerDay, &
262 : riv_temp_pcs
263 : use mo_read_timeseries, only : read_timeseries
264 : use mo_string_utils, only : num2str
265 :
266 : implicit none
267 :
268 : integer(i4) :: iGauge
269 :
270 : integer(i4) :: iDomain
271 :
272 : integer(i4) :: maxTimeSteps
273 :
274 : ! file name of file to read
275 : character(256) :: fName
276 :
277 : integer(i4), dimension(3) :: start_tmp, end_tmp
278 :
279 13 : real(dp), dimension(:), allocatable :: data_dp_1d
280 :
281 13 : logical, dimension(:), allocatable :: mask_1d
282 :
283 :
284 : !----------------------------------------------------------
285 : ! INITIALIZE RUNOFF
286 : !----------------------------------------------------------
287 51 : maxTimeSteps = maxval(simPer(1 : domainMeta%nDomains)%julEnd - simPer(1 : domainMeta%nDomains)%julStart + 1) * nTstepDay
288 52 : allocate(mRM_runoff(maxTimeSteps, nGaugesLocal))
289 358584 : mRM_runoff = nodata_dp
290 :
291 : ! READ GAUGE DATA
292 35 : do iGauge = 1, nGaugesLocal
293 : ! get domain id
294 22 : iDomain = gauge%domainId(iGauge)
295 : ! get start and end dates
296 88 : start_tmp = (/evalPer(iDomain)%yStart, evalPer(iDomain)%mStart, evalPer(iDomain)%dStart/)
297 88 : end_tmp = (/evalPer(iDomain)%yEnd, evalPer(iDomain)%mEnd, evalPer(iDomain)%dEnd /)
298 : ! evaluation gauge
299 22 : fName = trim(adjustl(gauge%fname(iGauge)))
300 : call read_timeseries(trim(fName), udischarge, &
301 : start_tmp, end_tmp, optimize, opti_function, &
302 22 : data_dp_1d, mask = mask_1d, nMeasPerDay = nMeasPerDay)
303 9540 : data_dp_1d = merge(data_dp_1d, nodata_dp, mask_1d)
304 22 : call paste(gauge%Q, data_dp_1d, nodata_dp)
305 35 : deallocate (data_dp_1d)
306 : ! TODO-RIV-TEMP: read temperature at gauge
307 : end do
308 : !
309 : ! inflow gauge
310 : !
311 : ! in mhm call InflowGauge%Q has to be initialized -- dummy allocation with period of domain 1 and initialization
312 13 : if (nInflowGaugesTotal .EQ. 0) then
313 58 : allocate(data_dp_1d(maxval(simPer(:)%julEnd - simPer(:)%julStart + 1)))
314 8043 : data_dp_1d = nodata_dp
315 12 : call paste(InflowGauge%Q, data_dp_1d, nodata_dp)
316 : else
317 3 : do iGauge = 1, nInflowGaugesTotal
318 : ! get domain id
319 2 : iDomain = InflowGauge%domainId(iGauge)
320 : ! get start and end dates
321 8 : start_tmp = (/simPer(iDomain)%yStart, simPer(iDomain)%mStart, simPer(iDomain)%dStart/)
322 8 : end_tmp = (/simPer(iDomain)%yEnd, simPer(iDomain)%mEnd, simPer(iDomain)%dEnd /)
323 : ! inflow gauge
324 2 : fName = trim(adjustl(InflowGauge%fname(iGauge)))
325 : call read_timeseries(trim(fName), udischarge, &
326 : start_tmp, end_tmp, optimize, opti_function, &
327 2 : data_dp_1d, mask = mask_1d, nMeasPerDay = nMeasPerDay)
328 1094 : if (.NOT. (all(mask_1d))) then
329 0 : call error_message('***ERROR: Nodata values in inflow gauge time series. File: ', trim(fName), raise=.false.)
330 0 : call error_message(' During simulation period from ', num2str(simPer(iDomain)%yStart) &
331 0 : , ' to ', num2str(simPer(iDomain)%yEnd))
332 : end if
333 1096 : data_dp_1d = merge(data_dp_1d, nodata_dp, mask_1d)
334 2 : call paste(InflowGauge%Q, data_dp_1d, nodata_dp)
335 3 : deallocate (data_dp_1d)
336 : end do
337 : end if
338 :
339 13 : end subroutine mrm_read_discharge
340 :
341 : ! ---------------------------------------------------------------------------
342 :
343 : ! NAME
344 : ! mrm_read_total_runoff
345 :
346 : ! PURPOSE
347 : !> \brief read simulated runoff that is to be routed
348 :
349 : !> \details read spatio-temporal field of total runoff that has been
350 : !> simulated by a hydrologic model or land surface model. This
351 : !> total runoff will then be aggregated to the level 11 resolution
352 : !> and then routed through the stream network.
353 :
354 : ! INTENT(IN)
355 : !> \param[in] "integer(i4) :: iDomain" domain id
356 :
357 : ! HISTORY
358 : !> \authors Stephan Thober
359 :
360 : !> \date Sep 2015
361 :
362 : ! Modifications:
363 : ! Stephan Thober Feb 2016 - refactored deallocate statements
364 : ! Stephan Thober Sep 2016 - added ALMA convention
365 : ! Robert Schweppe Jun 2018 - refactoring and reformatting
366 :
367 0 : subroutine mrm_read_total_runoff(iDomain)
368 :
369 13 : use mo_append, only : append
370 : use mo_constants, only : HourSecs
371 : use mo_common_constants, only : nodata_dp
372 : use mo_common_mHM_mRM_variables, only : simPer, timestep
373 : use mo_common_variables, only : ALMA_convention, level1
374 : use mo_mrm_global_variables, only : L1_total_runoff_in, dirTotalRunoff, filenameTotalRunoff, &
375 : varnameTotalRunoff
376 : use mo_read_nc, only : read_nc
377 :
378 : implicit none
379 :
380 : ! domain id
381 : integer(i4), intent(in) :: iDomain
382 :
383 : integer(i4) :: tt
384 :
385 : integer(i4) :: nTimeSteps
386 :
387 : ! tell nc file to read daily or hourly values
388 : integer(i4) :: nctimestep
389 :
390 : ! read data from file
391 0 : real(dp), dimension(:, :, :), allocatable :: L1_data
392 :
393 0 : real(dp), dimension(:, :), allocatable :: L1_data_packed
394 :
395 :
396 0 : if (timestep .eq. 1) nctimestep = -4 ! hourly input
397 0 : if (timestep .eq. 24) nctimestep = -1 ! daily input
398 0 : call read_nc(trim(dirTotalRunoff(iDomain)), level1(iDomain)%nrows, level1(iDomain)%ncols, &
399 0 : varnameTotalRunoff, level1(iDomain)%mask, L1_data, target_period = simPer(iDomain), &
400 0 : nctimestep = nctimestep, filename = filenameTotalRunoff)
401 : ! pack variables
402 0 : nTimeSteps = size(L1_data, 3)
403 0 : allocate(L1_data_packed(level1(iDomain)%nCells, nTimeSteps))
404 0 : do tt = 1, nTimeSteps
405 0 : L1_data_packed(:, tt) = pack(L1_data(:, :, tt), mask = level1(iDomain)%mask)
406 : end do
407 : ! free space immediately
408 0 : deallocate(L1_data)
409 :
410 : ! convert if ALMA conventions have been given
411 0 : if (ALMA_convention) then
412 : ! convert from kg m-2 s-1 to mm TS-1
413 : ! 1 kg m-2 -> 1 mm depth
414 : ! multiply with time to obtain per timestep
415 0 : L1_data_packed = L1_data_packed * timestep * HourSecs
416 : end if
417 :
418 : ! append
419 0 : call append(L1_total_runoff_in, L1_data_packed(:, :), nodata_dp)
420 :
421 : !free space
422 0 : deallocate(L1_data_packed)
423 :
424 0 : end subroutine mrm_read_total_runoff
425 :
426 0 : subroutine mrm_read_bankfull_runoff(iDomain)
427 : ! ---------------------------------------------------------------------------
428 :
429 : ! NAME
430 : ! mrm_read_bankfull_runoff
431 :
432 : !> \brief reads the bankfull runoff for approximating the channel widths
433 :
434 : !> \details reads the bankfull runoff, which can be calculated with
435 : !> the script in mhm/post_proc/bankfull_discharge.py
436 :
437 : ! INTENT(IN)
438 : !> \param[in] "integer(i4) :: iDomain" domain id
439 :
440 : ! INTENT(INOUT)
441 : ! None
442 :
443 : ! INTENT(OUT)
444 : ! None
445 :
446 : ! INTENT(IN), OPTIONAL
447 : ! None
448 :
449 : ! INTENT(INOUT), OPTIONAL
450 : ! None
451 :
452 : ! INTENT(OUT), OPTIONAL
453 : ! None
454 :
455 : ! RETURN
456 : ! None
457 :
458 : ! RESTRICTIONS
459 : !> \note The file read in must contain a double precision float variable with the name
460 : !> "Q_bkfl".
461 :
462 : ! EXAMPLE
463 : ! None
464 :
465 : ! LITERATURE
466 : ! None
467 :
468 : ! HISTORY
469 : ! \author Lennart Schueler
470 : ! \date May 2018
471 :
472 0 : use mo_mrm_global_variables, only: level11
473 : use mo_read_nc, only: read_const_nc
474 : use mo_mrm_global_variables, only: &
475 : dirBankfullRunoff, & ! directory of bankfull_runoff file for each domain
476 : L11_bankfull_runoff_in ! bankfull runoff at L1
477 :
478 : implicit none
479 :
480 : ! input variables
481 : integer(i4), intent(in) :: iDomain
482 :
483 : ! local variables
484 0 : real(dp), dimension(:,:), allocatable :: L11_data ! read data from file
485 0 : real(dp), dimension(:), allocatable :: L11_data_packed
486 :
487 0 : call read_const_nc(trim(dirBankfullRunoff(iDomain)), &
488 0 : level11(iDomain)%nrows, &
489 : level11(iDomain)%ncols, &
490 0 : "Q_bkfl", L11_data)
491 :
492 0 : allocate(L11_data_packed(level11(iDomain)%nCells))
493 0 : L11_data_packed(:) = pack(L11_data(:,:), mask=level11(iDomain)%mask)
494 :
495 : ! append
496 0 : if (allocated(L11_bankfull_runoff_in)) then
497 0 : L11_bankfull_runoff_in = [L11_bankfull_runoff_in, L11_data_packed]
498 : else
499 0 : allocate(L11_bankfull_runoff_in(size(L11_data_packed)))
500 0 : L11_bankfull_runoff_in = L11_data_packed
501 : end if
502 :
503 0 : deallocate(L11_data)
504 0 : deallocate(L11_data_packed)
505 :
506 0 : end subroutine mrm_read_bankfull_runoff
507 :
508 : ! NAME
509 : ! rotate_fdir_variable
510 :
511 : ! PURPOSE
512 : !> \brief TODO: add description
513 :
514 : !> \details TODO: add description
515 :
516 : ! INTENT(INOUT)
517 : !> \param[inout] "integer(i4), dimension(:, :) :: x"
518 :
519 : ! HISTORY
520 : !> \authors L. Samaniego & R. Kumar
521 :
522 : !> \date Jun 2018
523 :
524 : ! Modifications:
525 : ! Robert Schweppe Jun 2018 - refactoring and reformatting
526 :
527 21 : subroutine rotate_fdir_variable(x)
528 :
529 0 : use mo_common_constants, only : nodata_i4
530 :
531 : implicit none
532 :
533 : integer(i4), dimension(:, :), intent(INOUT) :: x
534 :
535 :
536 : ! 0 -1 0 |
537 : integer(i4) :: i, j
538 :
539 :
540 : !-------------------------------------------------------------------
541 : ! NOTE:
542 : !
543 : ! Since the DEM was transposed from (lat,lon) to (lon,lat), i.e.
544 : ! new DEM = transpose(DEM_old), then
545 : !
546 : ! the flow direction X (which was read) for every i, j cell needs
547 : ! to be rotated as follows
548 : !
549 : ! X(i,j) = [R] * {uVector}
550 : !
551 : ! with
552 : ! {uVector} denoting the fDir_old (read) in vector form, and
553 : ! e.g. dir 8 => {-1, -1, 0 }
554 : ! [R] denting a full rotation matrix to transform the flow
555 : ! direction into the new coordinate system (lon,lat).
556 : !
557 : ! [R] = [rx][rz]
558 : !
559 : ! with
560 : ! | 1 0 0 |
561 : ! [rx] =| 0 cos T sin T | = elemental rotation along x axis
562 : ! | 0 -sin T cos T |
563 : !
564 : ! | cos F sin F 0 |
565 : ! [rz] =| -sin F cos F 0 | = elemental rotation along z axis
566 : ! | 0 0 1 |
567 : !
568 : ! and T = pi, F = - pi/2
569 : ! thus
570 : ! | 0 -1 0 |
571 : ! [R] =| -1 0 0 |
572 : ! | 0 0 -1 |
573 : ! making all 8 directions the following transformation were
574 : ! obtained.
575 : !-------------------------------------------------------------------
576 5973 : do i = 1, size(x, 1)
577 2485077 : do j = 1, size(x, 2)
578 2479104 : if (x(i, j) .eq. nodata_i4) cycle
579 5952 : select case (x(i, j))
580 : case(1)
581 122693 : x(i, j) = 4
582 : case(2)
583 73066 : x(i, j) = 2
584 : case(4)
585 123704 : x(i, j) = 1
586 : case(8)
587 84531 : x(i, j) = 128
588 : case(16)
589 159247 : x(i, j) = 64
590 : case(32)
591 102046 : x(i, j) = 32
592 : case(64)
593 163544 : x(i, j) = 16
594 : case(128)
595 926405 : x(i, j) = 8
596 : end select
597 : end do
598 : end do
599 :
600 21 : end subroutine rotate_fdir_variable
601 :
602 : end module mo_mrm_read_data
|