SCALE-RM
mod_cnvlanduse.f90
Go to the documentation of this file.
1 !-------------------------------------------------------------------------------
9 !-------------------------------------------------------------------------------
10 module mod_cnvlanduse
11  !-----------------------------------------------------------------------------
12  !
13  !++ used modules
14  !
15  use scale_precision
16  use scale_stdio
17  use scale_prof
18  use scale_grid_index
19  use scale_tracer
20  !-----------------------------------------------------------------------------
21  implicit none
22  private
23  !-----------------------------------------------------------------------------
24  !
25  !++ Public procedure
26  !
27  public :: cnvlanduse_setup
28  public :: cnvlanduse
29 
30  !-----------------------------------------------------------------------------
31  !
32  !++ Public parameters & variables
33  !
34  logical, public :: cnvlanduse_donothing
35  logical, public :: cnvlanduse_useglccv2 = .false.
36  logical, public :: cnvlanduse_uselu100m = .false.
37  logical, public :: cnvlanduse_usejibis = .false.
38 
39  !-----------------------------------------------------------------------------
40  !
41  !++ Private procedure
42  !
43  private :: cnvlanduse_glccv2
44  private :: cnvlanduse_lu100m
45  private :: cnvlanduse_jibis
46 
47  !-----------------------------------------------------------------------------
48  !
49  !++ Private parameters & variables
50  !
51  real(RP), private :: cnvlanduse_unittile_ddeg = 0.0_rp ! dx for unit tile [deg]
52  real(RP), private :: cnvlanduse_oversampling_factor = 2.0_rp ! factor of min. dx against the unit tile
53 
54  !-----------------------------------------------------------------------------
55 contains
56  !-----------------------------------------------------------------------------
58  subroutine cnvlanduse_setup
59  use scale_process, only: &
60  prc_mpistop
61  use scale_const, only: &
62  d2r => const_d2r
63  use scale_comm, only: &
64  comm_horizontal_min
65  use scale_grid_real, only: &
66  real_dlat, &
67  real_dlon
68  implicit none
69 
70  character(len=H_SHORT) :: cnvlanduse_name = 'NONE' ! keep backward compatibility
71 
72  namelist / param_cnvlanduse / &
73  cnvlanduse_name, &
74  cnvlanduse_useglccv2, &
75  cnvlanduse_uselu100m, &
76  cnvlanduse_usejibis, &
77  cnvlanduse_unittile_ddeg, &
78  cnvlanduse_oversampling_factor
79 
80  real(RP) :: drad(ia,ja)
81  real(RP) :: drad_min
82 
83  integer :: ierr
84  !---------------------------------------------------------------------------
85 
86  if( io_l ) write(io_fid_log,*)
87  if( io_l ) write(io_fid_log,*) '++++++ Module[convert landuseindex] / Categ[preprocess] / Origin[SCALE-RM]'
88  !--- read namelist
89  rewind(io_fid_conf)
90  read(io_fid_conf,nml=param_cnvlanduse,iostat=ierr)
91  if( ierr < 0 ) then !--- missing
92  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
93  elseif( ierr > 0 ) then !--- fatal error
94  write(*,*) 'xxx Not appropriate names in namelist PARAM_CNVLANDUSE. Check!'
95  call prc_mpistop
96  endif
97  if( io_nml ) write(io_fid_nml,nml=param_cnvlanduse)
98 
99  select case(cnvlanduse_name)
100  case('NONE')
101  ! do nothing
102  case('GLCCv2')
103  cnvlanduse_useglccv2 = .true.
104  cnvlanduse_uselu100m = .false.
105  cnvlanduse_usejibis = .false.
106  case('LU100M')
107  cnvlanduse_useglccv2 = .false.
108  cnvlanduse_uselu100m = .true.
109  cnvlanduse_usejibis = .false.
110  case('COMBINE')
111  cnvlanduse_useglccv2 = .true.
112  cnvlanduse_uselu100m = .true.
113  cnvlanduse_usejibis = .true.
114  case default
115  write(*,*) 'xxx Unsupported TYPE:', trim(cnvlanduse_name)
116  call prc_mpistop
117  endselect
118 
119  cnvlanduse_donothing = .true.
120 
121  if ( cnvlanduse_useglccv2 ) then
122  cnvlanduse_donothing = .false.
123  if( io_l ) write(io_fid_log,*) '*** Use GLCC ver.2, global 30 arcsec. data'
124  if ( cnvlanduse_uselu100m ) then
125  if( io_l ) write(io_fid_log,*) '*** Use KSJ landuse 100m data for Japan region'
126  if( io_l ) write(io_fid_log,*) '*** Overwrite Japan region'
127  if ( cnvlanduse_usejibis ) then
128  if( io_l ) write(io_fid_log,*) '*** Use J-IBIS map 100m data for Japan region'
129  if( io_l ) write(io_fid_log,*) '*** Overwrite Japan region (PFT only)'
130  endif
131  endif
132  elseif ( cnvlanduse_uselu100m ) then
133  cnvlanduse_donothing = .false.
134  if( io_l ) write(io_fid_log,*) '*** Use KSJ landuse 100m data, Japan region only'
135  if ( cnvlanduse_usejibis ) then
136  if( io_l ) write(io_fid_log,*) '*** Use J-IBIS map 100m data for Japan region'
137  if( io_l ) write(io_fid_log,*) '*** Overwrite Japan region (PFT only)'
138  endif
139  endif
140 
141  if ( cnvlanduse_donothing ) then
142  if( io_l ) write(io_fid_log,*) '*** Do nothing for landuse index'
143  else
144  drad(:,:) = min( real_dlat(:,:), real_dlon(:,:) )
145  call comm_horizontal_min( drad_min, drad(:,:) )
146 
147  if ( cnvlanduse_unittile_ddeg > 0.0_rp ) then
148  cnvlanduse_oversampling_factor = ( drad_min / d2r ) / cnvlanduse_unittile_ddeg
149  endif
150  cnvlanduse_oversampling_factor = max( 1.0_rp, cnvlanduse_oversampling_factor )
151  cnvlanduse_unittile_ddeg = ( drad_min / d2r ) / cnvlanduse_oversampling_factor
152 
153  if( io_l ) write(io_fid_log,*) '*** The size of tile [deg] = ', cnvlanduse_unittile_ddeg
154  if( io_l ) write(io_fid_log,*) '*** oversampling factor = ', cnvlanduse_oversampling_factor
155  endif
156 
157  return
158  end subroutine cnvlanduse_setup
159 
160  !-----------------------------------------------------------------------------
162  subroutine cnvlanduse
163  use scale_process, only: &
164  prc_mpistop
165  use scale_landuse, only: &
166  landuse_calc_fact, &
167  landuse_write
168  implicit none
169  !---------------------------------------------------------------------------
170 
171  if ( cnvlanduse_donothing ) then
172  if( io_l ) write(io_fid_log,*)
173  if( io_l ) write(io_fid_log,*) '++++++ SKIP CONVERT LANDUSE DATA ++++++'
174  else
175  if( io_l ) write(io_fid_log,*)
176  if( io_l ) write(io_fid_log,*) '++++++ START CONVERT LANDUSE DATA ++++++'
177 
178  if ( cnvlanduse_useglccv2 ) then
179  call cnvlanduse_glccv2
180  endif
181 
182  if ( cnvlanduse_uselu100m ) then
183  call cnvlanduse_lu100m
184  endif
185 
186  if ( cnvlanduse_usejibis ) then
187  call cnvlanduse_jibis
188  endif
189 
190  ! calculate landuse factors
191  call landuse_calc_fact
192 
193  if( io_l ) write(io_fid_log,*) '++++++ END CONVERT LANDUSE DATA ++++++'
194 
195  ! output landuse file
196  call landuse_write
197  endif
198 
199  return
200  end subroutine cnvlanduse
201 
202  !-----------------------------------------------------------------------------
204  subroutine cnvlanduse_glccv2
205  use scale_process, only: &
206  prc_mpistop
207  use scale_const, only: &
208  radius => const_radius, &
209  pi => const_pi, &
210  eps => const_eps, &
211  d2r => const_d2r
212  use scale_landuse, only: &
213  landuse_frac_land, &
214  landuse_frac_lake, &
215  landuse_frac_urban, &
216  landuse_pft_mosaic, &
217  landuse_pft_nmax, &
218  landuse_frac_pft, &
219  landuse_index_pft
220  use scale_grid_real, only: &
221  real_laty, &
222  real_lonx
223  implicit none
224 
225  character(len=H_LONG) :: glccv2_in_catalogue = ''
226  character(len=H_LONG) :: glccv2_in_dir = ''
227  real(RP) :: limit_urban_fraction = 1.0_rp
228 
229  namelist / param_cnvlanduse_glccv2 / &
230  glccv2_in_catalogue, &
231  glccv2_in_dir, &
232  limit_urban_fraction
233 
234  ! data catalogue list
235  integer, parameter :: tile_nlim = 100
236  integer :: tile_nmax
237  real(RP) :: tile_lats (tile_nlim)
238  real(RP) :: tile_late (tile_nlim)
239  real(RP) :: tile_lons (tile_nlim)
240  real(RP) :: tile_lone (tile_nlim)
241  character(len=H_LONG) :: tile_fname(tile_nlim)
242 
243  ! GLCCv2 data
244  integer, parameter :: isize_orig = 3600
245  integer(1) :: tile_landuse_orig(isize_orig,isize_orig)
246  real(RP) :: tile_dlat_orig, tile_dlon_orig
247 
248  ! GLCCv2 data (oversampling)
249  integer :: ios
250  integer :: isize
251  integer(1), allocatable :: tile_landuse(:,:)
252  real(RP), allocatable :: tile_lath (:)
253  real(RP), allocatable :: tile_lonh (:)
254  real(RP) :: tile_dlat, tile_dlon
255  real(RP) :: area, area_fraction
256 
257  integer :: iloc
258  integer :: jloc
259  real(RP) :: ifrac_l ! fraction for iloc
260  real(RP) :: jfrac_b ! fraction for jloc
261 
262  real(RP) :: real_lonx_mod(0:ia,ja)
263  real(RP) :: domain_lats, domain_late
264  real(RP) :: domain_lons, domain_lone
265  integer :: domain_lonsloc(2), domain_loneloc(2)
266  logical :: check_idl
267 
268  ! tentative: LANDUSE_PFT_nmax is assumed to be 4
269  real(RP) :: categ_sum(ia,ja,-2:landuse_pft_nmax)
270 
271  integer :: lookuptable(1:25)
272  data lookuptable / 0, & ! 1 Urban and Built-Up Land -> 0 urban
273  7, & ! 2 Dryland Cropland and Pasture -> 7 Dryland Cropland and Pasture
274  8, & ! 3 Irrigated Cropland and Pasture -> 8 Irrigated Cropland and Pasture
275  9, & ! 4 Mixed Cropland and Pasture -> 9 Mixed Cropland and Pasture
276  5, & ! 5 Cropland/Grassland Mosaic -> 5 Cropland/Grassland Mosaic
277  6, & ! 6 Cropland/Woodland Mosaic -> 6 Cropland/Woodland Mosaic
278  2, & ! 7 Grassland -> 2 Grassland
279  3, & ! 8 Shrubland -> 3 Shrubland
280  4, & ! 9 Mixed Shrubland/Grassland -> 4 Mixed Shrubland/Grassland
281  4, & ! 10 Savanna -> 4 Mixed Shrubland/Grassland
282  11, & ! 11 Deciduous Broadleaf Forest -> 11 Deciduous Broadleaf Forest
283  12, & ! 12 Deciduous Needleleaf Forest -> 12 Deciduous Needleleaf Forest
284  13, & ! 13 Evergreen Broadleaf Forest -> 13 Deciduous Broadleaf Forest
285  14, & ! 14 Evergreen Needleleaf Forest -> 14 Deciduous Needleleaf Forest
286  15, & ! 15 Mixed Forest -> 15 Mixed Forest
287  -2, & ! 16 Water Bodies -> -2 Lake/River
288  10, & ! 17 Herbaceous Wetland -> 10 Paddy
289  10, & ! 18 Wooded Wetland -> 10 Paddy
290  1, & ! 19 Barren or Sparsely Vegetated -> 1 Bare Ground
291  1, & ! 20 Herbaceous Tundra -> 1 Bare Ground
292  1, & ! 21 Wooded Tundra -> 1 Bare Ground
293  1, & ! 22 Mixed Tundra -> 1 Bare Ground
294  1, & ! 23 Bare Ground Tundra -> 1 Bare Ground
295  1, & ! 24 Snow or Ice -> 1 Bare Ground
296  -1 / ! 25+ Sea Surface -> -1 Sea Surface
297 
298  real(RP) :: categ_pftsum, allsum
299  real(RP) :: pft (landuse_pft_nmax)
300  integer :: pft_idx(landuse_pft_nmax)
301  real(RP) :: temp
302  integer :: temp_idx
303 
304  character(len=H_LONG) :: fname
305 
306  real(RP) :: zerosw
307  logical :: hit_lat, hit_lon
308  integer :: index
309  integer :: fid, ierr
310  integer :: i, j, ii, jj, iii, jjj, t, p, pp
311  !---------------------------------------------------------------------------
312 
313  if( io_l ) write(io_fid_log,*)
314  if( io_l ) write(io_fid_log,*) '+++ Module[GLCCv2]/Categ[CNVLANDUSE]'
315 
316  !--- read namelist
317  rewind(io_fid_conf)
318  read(io_fid_conf,nml=param_cnvlanduse_glccv2,iostat=ierr)
319  if( ierr < 0 ) then !--- missing
320  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
321  elseif( ierr > 0 ) then !--- fatal error
322  write(*,*) 'xxx Not appropriate names in namelist PARAM_CNVLANDUSE_GLCCv2. Check!'
323  call prc_mpistop
324  endif
325  if( io_nml ) write(io_fid_nml,nml=param_cnvlanduse_glccv2)
326 
327  do p = -2, landuse_pft_nmax
328  do j = 1, ja
329  do i = 1, ia
330  categ_sum(i,j,p) = 0.0_rp
331  enddo
332  enddo
333  enddo
334 
335  real_lonx_mod(:,:) = mod( real_lonx(:,:)+3.0_dp*pi, 2.0_dp*pi ) - pi
336 
337  domain_lats = minval(real_laty(:,:))
338  domain_late = maxval(real_laty(:,:))
339  domain_lons = minval(real_lonx_mod(:,:))
340  domain_lone = maxval(real_lonx_mod(:,:))
341 
342  domain_lonsloc = minloc(real_lonx_mod(:,:))
343  domain_loneloc = maxloc(real_lonx_mod(:,:))
344 
345  check_idl = .false.
346  if ( domain_lons < real_lonx_mod(0 ,domain_lonsloc(2)) &
347  .OR. domain_lone > real_lonx_mod(ia,domain_loneloc(2)) ) then
348  check_idl = .true.
349  domain_lons = minval(real_lonx_mod(:,:),mask=(real_lonx_mod>0.0_rp))
350  domain_lone = maxval(real_lonx_mod(:,:),mask=(real_lonx_mod<0.0_rp))
351  endif
352 
353  ios = nint( 30.0_rp / 60.0_rp / 60.0_rp / cnvlanduse_unittile_ddeg - 0.5_rp ) + 1
354  isize = isize_orig * ios
355 
356  allocate( tile_landuse(isize,isize) )
357  allocate( tile_lath(0:isize) )
358  allocate( tile_lonh(0:isize) )
359 
360  if( io_l ) write(io_fid_log,*) '*** Oversampling orig = ', isize_orig, ', use = ', isize
361 
362  tile_dlat_orig = 30.0_rp / 60.0_rp / 60.0_rp * d2r
363  tile_dlon_orig = 30.0_rp / 60.0_rp / 60.0_rp * d2r
364  if( io_l ) write(io_fid_log,*) '*** TILE_DLAT :', tile_dlat_orig/d2r
365  if( io_l ) write(io_fid_log,*) '*** TILE_DLON :', tile_dlon_orig/d2r
366 
367  tile_dlat = tile_dlat_orig / ios
368  tile_dlon = tile_dlon_orig / ios
369  if( io_l ) write(io_fid_log,*) '*** TILE_DLAT (OS) :', tile_dlat/d2r
370  if( io_l ) write(io_fid_log,*) '*** TILE_DLON (OS) :', tile_dlon/d2r
371 
372  !---< READ from external files >---
373 
374  ! catalogue file
375  fname = trim(glccv2_in_dir)//'/'//trim(glccv2_in_catalogue)
376 
377  if( io_l ) write(io_fid_log,*)
378  if( io_l ) write(io_fid_log,*) '+++ Input catalogue file:', trim(fname)
379 
380  fid = io_get_available_fid()
381  open( fid, &
382  file = trim(fname), &
383  form = 'formatted', &
384  status = 'old', &
385  iostat = ierr )
386 
387  if ( ierr /= 0 ) then
388  write(*,*) 'xxx catalogue file not found!', trim(fname)
389  call prc_mpistop
390  endif
391 
392  do t = 1, tile_nlim
393  read(fid,*,iostat=ierr) index, tile_lats(t), tile_late(t), & ! South->North
394  tile_lons(t), tile_lone(t), & ! WEST->EAST
395  tile_fname(t)
396  if ( ierr /= 0 ) exit
397 
398  if ( tile_lons(t) >= 180.0_rp ) then
399  tile_lons(t) = tile_lons(t) - 360.0_rp
400  tile_lone(t) = tile_lone(t) - 360.0_rp
401  endif
402  if ( tile_lons(t) < -180.0_rp ) tile_lons(t) = tile_lons(t) + 360.0_rp
403  if ( tile_lone(t) < -180.0_rp ) tile_lone(t) = tile_lone(t) + 360.0_rp
404 
405  enddo
406 
407  tile_nmax = t - 1
408  close(fid)
409 
410  ! data file
411  do t = 1, tile_nmax
412  hit_lat = .false.
413  hit_lon = .false.
414 
415  if ( ( tile_lats(t)*d2r >= domain_lats .AND. tile_lats(t)*d2r < domain_late ) &
416  .OR. ( tile_late(t)*d2r >= domain_lats .AND. tile_late(t)*d2r < domain_late ) ) then
417  hit_lat = .true.
418  endif
419 
420  if ( ( domain_lats >= tile_lats(t)*d2r .AND. domain_lats < tile_late(t)*d2r ) &
421  .OR. ( domain_late >= tile_lats(t)*d2r .AND. domain_late < tile_late(t)*d2r ) ) then
422  hit_lat = .true.
423  endif
424 
425  if ( check_idl ) then
426  if ( ( tile_lons(t)*d2r >= domain_lons .AND. tile_lons(t)*d2r < pi ) &
427  .OR. ( tile_lons(t)*d2r >= -pi .AND. tile_lons(t)*d2r < domain_lone ) &
428  .OR. ( tile_lone(t)*d2r >= domain_lons .AND. tile_lone(t)*d2r < pi ) &
429  .OR. ( tile_lone(t)*d2r >= -pi .AND. tile_lone(t)*d2r < domain_lone ) ) then
430  hit_lon = .true.
431  endif
432  else
433  if ( ( tile_lons(t)*d2r >= domain_lons .AND. tile_lons(t)*d2r < domain_lone ) &
434  .OR. ( tile_lone(t)*d2r >= domain_lons .AND. tile_lone(t)*d2r < domain_lone ) ) then
435  hit_lon = .true.
436  endif
437  endif
438 
439  if ( ( domain_lons >= tile_lons(t)*d2r .AND. domain_lons < tile_lone(t)*d2r ) &
440  .OR. ( domain_lone >= tile_lons(t)*d2r .AND. domain_lone < tile_lone(t)*d2r ) ) then
441  hit_lon = .true.
442  endif
443 
444  if ( hit_lat .AND. hit_lon ) then
445  fname = trim(glccv2_in_dir)//'/'//trim(tile_fname(t))
446 
447  if( io_l ) write(io_fid_log,*)
448  if( io_l ) write(io_fid_log,*) '+++ Input data file :', trim(fname)
449  if( io_l ) write(io_fid_log,*) '*** Domain (LAT) :', domain_lats/d2r, domain_late/d2r
450  if( io_l ) write(io_fid_log,*) '*** (LON) :', domain_lons/d2r, domain_lone/d2r
451  if ( check_idl ) then
452  if( io_l ) write(io_fid_log,*) '*** (Date line exists within the domain)'
453  endif
454  if( io_l ) write(io_fid_log,*) '*** Tile (LAT) :', tile_lats(t), tile_late(t)
455  if( io_l ) write(io_fid_log,*) '*** (LON) :', tile_lons(t), tile_lone(t)
456 
457  fid = io_get_available_fid()
458  open( fid, &
459  file = trim(fname), &
460  form = 'unformatted', &
461  access = 'direct', &
462  status = 'old', &
463  recl = isize_orig*isize_orig*1, &
464  iostat = ierr )
465 
466  if ( ierr /= 0 ) then
467  write(*,*) 'xxx data file not found!'
468  call prc_mpistop
469  endif
470 
471  read(fid,rec=1) tile_landuse_orig(:,:)
472  close(fid)
473 
474  ! oversampling
475  do jj = 1, isize_orig
476  do ii = 1, isize_orig
477  do j = 1, ios
478  do i = 1, ios
479  jjj = (jj-1) * ios + j
480  iii = (ii-1) * ios + i
481 
482  tile_landuse(iii,jjj) = tile_landuse_orig(ii,jj)
483  enddo
484  enddo
485  enddo
486  enddo
487 
488  tile_lath(0) = tile_lats(t) * d2r
489  do jj = 1, isize
490  tile_lath(jj) = tile_lath(jj-1) + tile_dlat
491 ! if( IO_L ) write(IO_FID_LOG,*) jj, TILE_LATH(jj)
492  enddo
493 
494  tile_lonh(0) = tile_lons(t) * d2r
495  do ii = 1, isize
496  tile_lonh(ii) = tile_lonh(ii-1) + tile_dlon
497 ! if( IO_L ) write(IO_FID_LOG,*) ii, TILE_LONH(ii)
498  enddo
499 
500  ! match and calc fraction
501  do jj = 1, isize
502  do ii = 1, isize
503 
504  iloc = 1 ! Z_sfc(1,1) is used for dummy grid
505  ifrac_l = 1.0_rp
506 
507  jloc = 1 ! Z_sfc(1,1) is used for dummy grid
508  jfrac_b = 1.0_rp
509 
510  if ( tile_lath(jj ) < domain_lats &
511  .OR. tile_lath(jj-1) > domain_late ) then
512  cycle
513  endif
514 
515  if ( check_idl ) then
516  if ( tile_lonh(ii ) < domain_lons &
517  .AND. tile_lonh(ii-1) > domain_lone ) then
518  cycle
519  endif
520  else
521  if ( tile_lonh(ii ) < domain_lons &
522  .OR. tile_lonh(ii-1) > domain_lone ) then
523  cycle
524  endif
525  endif
526 
527  jloop: do j = js-1, je+1
528  iloop: do i = is-1, ie+1
529  if ( tile_lonh(ii-1) >= real_lonx_mod(i-1,j ) &
530  .AND. tile_lonh(ii-1) < real_lonx_mod(i ,j ) &
531  .AND. tile_lath(jj-1) >= real_laty(i ,j-1) &
532  .AND. tile_lath(jj-1) < real_laty(i ,j ) ) then
533 
534  iloc = i
535  ifrac_l = min( real_lonx_mod(i,j)-tile_lonh(ii-1), tile_dlon ) / tile_dlon
536 
537  jloc = j
538  jfrac_b = min( real_laty(i,j)-tile_lath(jj-1), tile_dlat ) / tile_dlat
539  exit jloop
540 
541  endif
542 
543  if ( real_lonx_mod(i-1,j) >= real_lonx_mod(i ,j ) &
544  .AND. tile_lath(jj-1) >= real_laty(i ,j-1) &
545  .AND. tile_lath(jj-1) < real_laty(i ,j ) ) then ! across the IDL
546 
547  if ( tile_lonh(ii-1) >= real_lonx_mod(i-1,j) &
548  .AND. tile_lonh(ii-1) < pi ) then
549 
550  iloc = i
551  ifrac_l = min( real_lonx_mod(i,j)-tile_lonh(ii-1)+2.0_rp*pi, tile_dlon ) / tile_dlon
552 
553  jloc = j
554  jfrac_b = min( real_laty(i,j)-tile_lath(jj-1), tile_dlat ) / tile_dlat
555  exit jloop
556 
557  elseif( tile_lonh(ii-1) >= -pi &
558  .AND. tile_lonh(ii-1) < real_lonx_mod(i ,j) ) then
559 
560  iloc = i
561  ifrac_l = min( real_lonx_mod(i,j)-tile_lonh(ii-1), tile_dlon ) / tile_dlon
562 
563  jloc = j
564  jfrac_b = min( real_laty(i,j)-tile_lath(jj-1), tile_dlat ) / tile_dlat
565  exit jloop
566 
567  endif
568 
569  endif
570  enddo iloop
571  enddo jloop
572 
573  if( iloc == 1 .AND. jloc == 1 ) cycle
574 
575  area = radius * radius * tile_dlon * ( sin(tile_lath(jj))-sin(tile_lath(jj-1)) )
576 
577  pp = min( max( int(tile_landuse(ii,jj),kind=4), 0 ), 25 )
578  p = lookuptable(pp)
579 
580  !if( IO_L ) write(IO_FID_LOG,*) ii, jj, iloc, jloc, p, pp
581 
582  area_fraction = ( ifrac_l) * ( jfrac_b) * area
583  categ_sum(iloc ,jloc ,p) = categ_sum(iloc ,jloc ,p) + area_fraction
584 
585  area_fraction = (1.0_rp-ifrac_l) * ( jfrac_b) * area
586  categ_sum(iloc+1,jloc ,p) = categ_sum(iloc+1,jloc ,p) + area_fraction
587 
588  area_fraction = ( ifrac_l) * (1.0_rp-jfrac_b) * area
589  categ_sum(iloc ,jloc+1,p) = categ_sum(iloc ,jloc+1,p) + area_fraction
590 
591  area_fraction = (1.0_rp-ifrac_l) * (1.0_rp-jfrac_b) * area
592  categ_sum(iloc+1,jloc+1,p) = categ_sum(iloc+1,jloc+1,p) + area_fraction
593  enddo
594  enddo
595 
596  endif
597  enddo ! tile loop
598 
599  do j = js, je
600  do i = is, ie
601 ! area = RADIUS * RADIUS * (REAL_LONX_mod(i,j)-REAL_LONX_mod(i-1,j)) * ( sin(REAL_LATY(i,j))-sin(REAL_LATY(i,j-1)) )
602 ! allsum = categ_sum(i,j,-2) + categ_sum(i,j,-1) + categ_sum(i,j,0) + categ_pftsum
603 ! if ( abs(allsum/area-1.0_RP) > EPS ) then
604 ! if( IO_L ) write(IO_FID_LOG,*) i,j,allsum/area
605 ! endif
606 
607  do p = 1, landuse_pft_nmax
608  pft(p) = categ_sum(i,j,p)
609  pft_idx(p) = p
610  enddo
611 
612  ! bubble sort
613  do p = 1, landuse_pft_nmax-1
614  do pp = p+1, landuse_pft_nmax
615  if ( pft(pp) > pft(p) ) then
616  temp = pft(p)
617  pft(p) = pft(pp)
618  pft(pp) = temp
619 
620  temp_idx = pft_idx(p)
621  pft_idx(p) = pft_idx(pp)
622  pft_idx(pp) = temp_idx
623  endif
624  enddo
625  enddo
626 
627  categ_pftsum = sum( pft(:) )
628 
629  ! land fraction : 1 - ocean / total
630  allsum = categ_sum(i,j,-2) + categ_sum(i,j,-1) + categ_sum(i,j,0) + categ_pftsum
631  zerosw = 0.5_rp - sign( 0.5_rp, allsum-eps )
632  landuse_frac_land(i,j) = ( allsum-categ_sum(i,j,-1) ) * ( 1.0_rp-zerosw ) / ( allsum-zerosw )
633 
634  ! lake fraction : lake / ( total - ocean )
635  allsum = categ_sum(i,j,-2) + categ_sum(i,j,0) + categ_pftsum
636  zerosw = 0.5_rp - sign( 0.5_rp, allsum-eps )
637  landuse_frac_lake(i,j) = categ_sum(i,j,-2) * ( 1.0_rp-zerosw ) / ( allsum-zerosw )
638 
639  ! urban fraction : urban / ( total - ocean - lake )
640  allsum = categ_sum(i,j,0) + categ_pftsum
641  zerosw = 0.5_rp - sign( 0.5_rp, allsum-eps )
642  landuse_frac_urban(i,j) = categ_sum(i,j,0) * ( 1.0_rp-zerosw ) / ( allsum-zerosw )
643 
644  ! PFT fraction : PFT / sum( PFT(1:mosaic) )
645  allsum = sum( pft(1:landuse_pft_mosaic) )
646  if ( abs(allsum) > eps ) then
647  do p = 1, landuse_pft_mosaic
648  landuse_frac_pft(i,j,p) = pft(p) / allsum
649  landuse_index_pft(i,j,p) = pft_idx(p)
650  enddo
651  ! if no second PFT, set to same as PFT1
652  if ( abs(landuse_frac_pft(i,j,1)-1.0_rp) <= eps ) then
653  landuse_frac_pft(i,j,:) = 0.0_rp
654  landuse_frac_pft(i,j,1) = 1.0_rp
655  landuse_index_pft(i,j,:) = pft_idx(1)
656  endif
657  else ! if no PFT, set to bare ground
658  landuse_frac_pft(i,j,:) = 0.0_rp
659  landuse_frac_pft(i,j,1) = 1.0_rp
660  landuse_index_pft(i,j,:) = 1
661  endif
662 
663  enddo
664  enddo
665 
666  if ( limit_urban_fraction < 1.0_rp ) then
667  do j = js, je
668  do i = is, ie
669  if ( landuse_frac_urban(i,j) == 1.0_rp ) then ! if no PFT, set to grassland
670  landuse_frac_pft(i,j,:) = 0.0_rp
671  landuse_frac_pft(i,j,1) = 1.0_rp
672  landuse_index_pft(i,j,:) = 2
673  endif
674  landuse_frac_urban(i,j) = min( landuse_frac_urban(i,j), limit_urban_fraction )
675  enddo
676  enddo
677  endif
678 
679  return
680  end subroutine cnvlanduse_glccv2
681 
682  !-----------------------------------------------------------------------------
684  subroutine cnvlanduse_lu100m
685  use scale_process, only: &
686  prc_mpistop
687  use scale_const, only: &
688  radius => const_radius, &
689  pi => const_pi, &
690  eps => const_eps, &
691  d2r => const_d2r
692  use scale_landuse, only: &
693  landuse_frac_land, &
694  landuse_frac_lake, &
695  landuse_frac_urban, &
696  landuse_pft_mosaic, &
697  landuse_pft_nmax, &
698  landuse_frac_pft, &
699  landuse_index_pft
700  use scale_grid_real, only: &
701  real_laty, &
702  real_lonx
703  implicit none
704 
705  character(len=H_LONG) :: lu100m_in_catalogue = ''
706  character(len=H_LONG) :: lu100m_in_dir = ''
707  real(RP) :: limit_urban_fraction = 1.0_rp
708 
709  namelist / param_cnvlanduse_lu100m / &
710  lu100m_in_catalogue, &
711  lu100m_in_dir, &
712  limit_urban_fraction
713 
714  ! data catalogue list
715  integer, parameter :: tile_nlim = 1000
716  integer :: tile_nmax
717  real(RP) :: tile_lats (tile_nlim)
718  real(RP) :: tile_late (tile_nlim)
719  real(RP) :: tile_lons (tile_nlim)
720  real(RP) :: tile_lone (tile_nlim)
721  character(len=H_LONG) :: tile_fname(tile_nlim)
722 
723  ! LU100M data
724  integer, parameter :: isize_orig = 800
725  real(SP) :: tile_landuse_orig(isize_orig,isize_orig)
726  real(RP) :: tile_dlat_orig, tile_dlon_orig
727 
728  ! LU100M data (oversampling)
729  integer :: ios
730  integer :: isize
731  real(RP), allocatable :: tile_landuse(:,:)
732  real(RP), allocatable :: tile_lath (:)
733  real(RP), allocatable :: tile_lonh (:)
734  real(RP) :: tile_dlat, tile_dlon
735  real(RP) :: area, area_fraction
736 
737  integer :: iloc
738  integer :: jloc
739  real(RP) :: ifrac_l ! fraction for iloc
740  real(RP) :: jfrac_b ! fraction for jloc
741 
742  real(RP) :: real_lonx_mod(0:ia,ja)
743  real(RP) :: domain_lats, domain_late
744  real(RP) :: domain_lons, domain_lone
745  integer :: domain_lonsloc(2), domain_loneloc(2)
746  logical :: check_idl
747 
748  ! tentative: LANDUSE_PFT_nmax is assumed to be 4
749  real(RP) :: categ_sum(ia,ja,-2:landuse_pft_nmax)
750 
751  integer :: lookuptable(0:16)
752  data lookuptable / -1, & ! 0 missing -> -1 Sea Surface
753  10, & ! 1 paddy -> 10 Paddy
754  9, & ! 2 cropland -> 9 Mixed Cropland and Pasture
755  1, & ! 3 UNDEF -> 1 Bare Ground
756  1, & ! 4 UNDEF -> 1 Bare Ground
757  11, & ! 5 forest -> 11 Deciduous Broadleaf Forest
758  1, & ! 6 bareground -> 1 Bare Ground
759  0, & ! 7 urban building -> 0 Urban and Built-up Land
760  1, & ! 8 UNDEF -> 1 Bare Ground
761  0, & ! 9 motorway -> 0 Urban and Built-up Land
762  0, & ! 10 urban ground -> 0 Urban and Built-up Land
763  -2, & ! 11 lake,river -> -2 Lake/River
764  1, & ! 12 UNDEF -> 1 Bare Ground
765  1, & ! 13 UNDEF -> 1 Bare Ground
766  1, & ! 14 seashore -> 1 Bare Ground
767  -1, & ! 15 ocean -> -1 Sea Surface
768  2 / ! 16 golf course -> 2 Grassland
769 
770  real(RP) :: categ_pftsum, allsum
771  real(RP) :: pft (landuse_pft_nmax)
772  integer :: pft_idx(landuse_pft_nmax)
773  real(RP) :: temp
774  integer :: temp_idx
775  real(RP) :: frac, mask
776 
777  character(len=H_LONG) :: fname
778 
779  real(RP) :: zerosw
780  logical :: hit_lat, hit_lon
781  integer :: index
782  integer :: fid, ierr
783  integer :: i, j, ii, jj, iii, jjj, t, p, pp
784  !---------------------------------------------------------------------------
785 
786  if( io_l ) write(io_fid_log,*)
787  if( io_l ) write(io_fid_log,*) '+++ Module[LU100M]/Categ[CNVLANDUSE]'
788 
789  !--- read namelist
790  rewind(io_fid_conf)
791  read(io_fid_conf,nml=param_cnvlanduse_lu100m,iostat=ierr)
792  if( ierr < 0 ) then !--- missing
793  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
794  elseif( ierr > 0 ) then !--- fatal error
795  write(*,*) 'xxx Not appropriate names in namelist PARAM_CNVLANDUSE_LU100M. Check!'
796  call prc_mpistop
797  endif
798  if( io_nml ) write(io_fid_nml,nml=param_cnvlanduse_lu100m)
799 
800  do p = -2, landuse_pft_nmax
801  do j = 1, ja
802  do i = 1, ia
803  categ_sum(i,j,p) = 0.0_rp
804  enddo
805  enddo
806  enddo
807 
808  real_lonx_mod(:,:) = mod( real_lonx(:,:)+3.0_dp*pi, 2.0_dp*pi ) - pi
809 
810  domain_lats = minval(real_laty(:,:))
811  domain_late = maxval(real_laty(:,:))
812  domain_lons = minval(real_lonx_mod(:,:))
813  domain_lone = maxval(real_lonx_mod(:,:))
814 
815  domain_lonsloc = minloc(real_lonx_mod(:,:))
816  domain_loneloc = maxloc(real_lonx_mod(:,:))
817 
818  check_idl = .false.
819  if ( domain_lons < real_lonx_mod(0 ,domain_lonsloc(2)) &
820  .OR. domain_lone > real_lonx_mod(ia,domain_loneloc(2)) ) then
821  check_idl = .true.
822  domain_lons = minval(real_lonx_mod(:,:),mask=(real_lonx_mod>0.0_rp))
823  domain_lone = maxval(real_lonx_mod(:,:),mask=(real_lonx_mod<0.0_rp))
824  endif
825 
826  ios = nint( 5.0_rp / 60.0_rp / 100.0_rp / cnvlanduse_unittile_ddeg - 0.5_rp ) + 1
827  isize = isize_orig * ios
828 
829  allocate( tile_landuse(isize,isize) )
830  allocate( tile_lath(0:isize) )
831  allocate( tile_lonh(0:isize) )
832 
833  if( io_l ) write(io_fid_log,*) '*** Oversampling orig = ', isize_orig, ', use = ', isize
834 
835  tile_dlat_orig = 5.0_rp / 60.0_rp / 100.0_rp * d2r
836  tile_dlon_orig = 7.5_rp / 60.0_rp / 100.0_rp * d2r
837  if( io_l ) write(io_fid_log,*) '*** TILE_DLAT :', tile_dlat_orig/d2r
838  if( io_l ) write(io_fid_log,*) '*** TILE_DLON :', tile_dlon_orig/d2r
839 
840  tile_dlat = tile_dlat_orig / ios
841  tile_dlon = tile_dlon_orig / ios
842  if( io_l ) write(io_fid_log,*) '*** TILE_DLAT (OS) :', tile_dlat/d2r
843  if( io_l ) write(io_fid_log,*) '*** TILE_DLON (OS) :', tile_dlon/d2r
844 
845  !---< READ from external files >---
846 
847  ! catalogue file
848  fname = trim(lu100m_in_dir)//'/'//trim(lu100m_in_catalogue)
849 
850  if( io_l ) write(io_fid_log,*)
851  if( io_l ) write(io_fid_log,*) '+++ Input catalogue file:', trim(fname)
852 
853  fid = io_get_available_fid()
854  open( fid, &
855  file = trim(fname), &
856  form = 'formatted', &
857  status = 'old', &
858  iostat = ierr )
859 
860  if ( ierr /= 0 ) then
861  write(*,*) 'xxx catalogue file not found!', trim(fname)
862  call prc_mpistop
863  endif
864 
865  do t = 1, tile_nlim
866  read(fid,*,iostat=ierr) index, tile_lats(t), tile_late(t), & ! South->North
867  tile_lons(t), tile_lone(t), & ! WEST->EAST
868  tile_fname(t)
869  if ( ierr /= 0 ) exit
870 
871  if ( tile_lons(t) >= 180.0_rp ) then
872  tile_lons(t) = tile_lons(t) - 360.0_rp
873  tile_lone(t) = tile_lone(t) - 360.0_rp
874  endif
875  if ( tile_lons(t) < -180.0_rp ) tile_lons(t) = tile_lons(t) + 360.0_rp
876  if ( tile_lone(t) < -180.0_rp ) tile_lone(t) = tile_lone(t) + 360.0_rp
877 
878  enddo
879 
880  tile_nmax = t - 1
881  close(fid)
882 
883  ! data file
884  do t = 1, tile_nmax
885  hit_lat = .false.
886  hit_lon = .false.
887 
888  if ( ( tile_lats(t)*d2r >= domain_lats .AND. tile_lats(t)*d2r < domain_late ) &
889  .OR. ( tile_late(t)*d2r >= domain_lats .AND. tile_late(t)*d2r < domain_late ) ) then
890  hit_lat = .true.
891  endif
892 
893  if ( ( domain_lats >= tile_lats(t)*d2r .AND. domain_lats < tile_late(t)*d2r ) &
894  .OR. ( domain_late >= tile_lats(t)*d2r .AND. domain_late < tile_late(t)*d2r ) ) then
895  hit_lat = .true.
896  endif
897 
898  if ( check_idl ) then
899  if ( ( tile_lons(t)*d2r >= domain_lons .AND. tile_lons(t)*d2r < pi ) &
900  .OR. ( tile_lons(t)*d2r >= -pi .AND. tile_lons(t)*d2r < domain_lone ) &
901  .OR. ( tile_lone(t)*d2r >= domain_lons .AND. tile_lone(t)*d2r < pi ) &
902  .OR. ( tile_lone(t)*d2r >= -pi .AND. tile_lone(t)*d2r < domain_lone ) ) then
903  hit_lon = .true.
904  endif
905  else
906  if ( ( tile_lons(t)*d2r >= domain_lons .AND. tile_lons(t)*d2r < domain_lone ) &
907  .OR. ( tile_lone(t)*d2r >= domain_lons .AND. tile_lone(t)*d2r < domain_lone ) ) then
908  hit_lon = .true.
909  endif
910  endif
911 
912  if ( ( domain_lons >= tile_lons(t)*d2r .AND. domain_lons < tile_lone(t)*d2r ) &
913  .OR. ( domain_lone >= tile_lons(t)*d2r .AND. domain_lone < tile_lone(t)*d2r ) ) then
914  hit_lon = .true.
915  endif
916 
917  if ( hit_lat .AND. hit_lon ) then
918  fname = trim(lu100m_in_dir)//'/'//trim(tile_fname(t))
919 
920  if( io_l ) write(io_fid_log,*)
921  if( io_l ) write(io_fid_log,*) '+++ Input data file :', trim(fname)
922  if( io_l ) write(io_fid_log,*) '*** Domain (LAT) :', domain_lats/d2r, domain_late/d2r
923  if( io_l ) write(io_fid_log,*) '*** (LON) :', domain_lons/d2r, domain_lone/d2r
924  if ( check_idl ) then
925  if( io_l ) write(io_fid_log,*) '*** (Date line exists within the domain)'
926  endif
927  if( io_l ) write(io_fid_log,*) '*** Tile (LAT) :', tile_lats(t), tile_late(t)
928  if( io_l ) write(io_fid_log,*) '*** (LON) :', tile_lons(t), tile_lone(t)
929 
930  fid = io_get_available_fid()
931  open( fid, &
932  file = trim(fname), &
933  form = 'unformatted', &
934  access = 'direct', &
935  status = 'old', &
936  recl = isize_orig*isize_orig*4, &
937  iostat = ierr )
938 
939  if ( ierr /= 0 ) then
940  write(*,*) 'xxx data file not found!'
941  call prc_mpistop
942  endif
943 
944  read(fid,rec=1) tile_landuse_orig(:,:)
945  close(fid)
946 
947  ! oversampling
948  do jj = 1, isize_orig
949  do ii = 1, isize_orig
950  do j = 1, ios
951  do i = 1, ios
952  jjj = (jj-1) * ios + j
953  iii = (ii-1) * ios + i
954 
955  tile_landuse(iii,jjj) = real( TILE_LANDUSE_orig(ii,jj), kind=rp )
956  enddo
957  enddo
958  enddo
959  enddo
960 
961  tile_lath(0) = tile_lats(t) * d2r
962  do jj = 1, isize
963  tile_lath(jj) = tile_lath(jj-1) + tile_dlat
964 ! if( IO_L ) write(IO_FID_LOG,*) jj, TILE_LATH(jj)
965  enddo
966 
967  tile_lonh(0) = tile_lons(t) * d2r
968  do ii = 1, isize
969  tile_lonh(ii) = tile_lonh(ii-1) + tile_dlon
970 ! if( IO_L ) write(IO_FID_LOG,*) ii, TILE_LONH(ii)
971  enddo
972 
973  ! match and calc fraction
974  do jj = 1, isize
975  do ii = 1, isize
976 
977  iloc = 1 ! Z_sfc(1,1) is used for dummy grid
978  ifrac_l = 1.0_rp
979 
980  jloc = 1 ! Z_sfc(1,1) is used for dummy grid
981  jfrac_b = 1.0_rp
982 
983  if ( tile_lath(jj ) < domain_lats &
984  .OR. tile_lath(jj-1) > domain_late ) then
985  cycle
986  endif
987 
988  if ( check_idl ) then
989  if ( tile_lonh(ii ) < domain_lons &
990  .AND. tile_lonh(ii-1) > domain_lone ) then
991  cycle
992  endif
993  else
994  if ( tile_lonh(ii ) < domain_lons &
995  .OR. tile_lonh(ii-1) > domain_lone ) then
996  cycle
997  endif
998  endif
999 
1000  jloop: do j = js-1, je+1
1001  iloop: do i = is-1, ie+1
1002  if ( tile_lonh(ii-1) >= real_lonx_mod(i-1,j ) &
1003  .AND. tile_lonh(ii-1) < real_lonx_mod(i ,j ) &
1004  .AND. tile_lath(jj-1) >= real_laty(i ,j-1) &
1005  .AND. tile_lath(jj-1) < real_laty(i ,j ) ) then
1006 
1007  iloc = i
1008  ifrac_l = min( real_lonx_mod(i,j)-tile_lonh(ii-1), tile_dlon ) / tile_dlon
1009 
1010  jloc = j
1011  jfrac_b = min( real_laty(i,j)-tile_lath(jj-1), tile_dlat ) / tile_dlat
1012  exit jloop
1013 
1014  endif
1015 
1016  if ( real_lonx_mod(i-1,j) >= real_lonx_mod(i ,j ) &
1017  .AND. tile_lath(jj-1) >= real_laty(i ,j-1) &
1018  .AND. tile_lath(jj-1) < real_laty(i ,j ) ) then ! across the IDL
1019 
1020  if ( tile_lonh(ii-1) >= real_lonx_mod(i-1,j) &
1021  .AND. tile_lonh(ii-1) < pi ) then
1022 
1023  iloc = i
1024  ifrac_l = min( real_lonx_mod(i,j)-tile_lonh(ii-1)+2.0_rp*pi, tile_dlon ) / tile_dlon
1025 
1026  jloc = j
1027  jfrac_b = min( real_laty(i,j)-tile_lath(jj-1), tile_dlat ) / tile_dlat
1028  exit jloop
1029 
1030  elseif( tile_lonh(ii-1) >= -pi &
1031  .AND. tile_lonh(ii-1) < real_lonx_mod(i ,j) ) then
1032 
1033  iloc = i
1034  ifrac_l = min( real_lonx_mod(i,j)-tile_lonh(ii-1), tile_dlon ) / tile_dlon
1035 
1036  jloc = j
1037  jfrac_b = min( real_laty(i,j)-tile_lath(jj-1), tile_dlat ) / tile_dlat
1038  exit jloop
1039 
1040  endif
1041 
1042  endif
1043  enddo iloop
1044  enddo jloop
1045 
1046  if( iloc == 1 .AND. jloc == 1 ) cycle
1047 
1048  area = radius * radius * tile_dlon * ( sin(tile_lath(jj))-sin(tile_lath(jj-1)) )
1049 
1050  pp = int( max( tile_landuse(ii,jj), 0.0_rp ) )
1051  p = lookuptable(pp)
1052 
1053  !if( IO_L ) write(IO_FID_LOG,*) ii, jj, iloc, jloc, p, pp
1054 
1055  area_fraction = ( ifrac_l) * ( jfrac_b) * area
1056  categ_sum(iloc ,jloc ,p) = categ_sum(iloc ,jloc ,p) + area_fraction
1057 
1058  area_fraction = (1.0_rp-ifrac_l) * ( jfrac_b) * area
1059  categ_sum(iloc+1,jloc ,p) = categ_sum(iloc+1,jloc ,p) + area_fraction
1060 
1061  area_fraction = ( ifrac_l) * (1.0_rp-jfrac_b) * area
1062  categ_sum(iloc ,jloc+1,p) = categ_sum(iloc ,jloc+1,p) + area_fraction
1063 
1064  area_fraction = (1.0_rp-ifrac_l) * (1.0_rp-jfrac_b) * area
1065  categ_sum(iloc+1,jloc+1,p) = categ_sum(iloc+1,jloc+1,p) + area_fraction
1066  enddo
1067  enddo
1068 
1069  endif
1070  enddo ! tile loop
1071 
1072  do j = js, je
1073  do i = is, ie
1074 ! area = RADIUS * RADIUS * (REAL_LONX_mod(i,j)-REAL_LONX_mod(i-1,j)) * ( sin(REAL_LATY(i,j))-sin(REAL_LATY(i,j-1)) )
1075 ! allsum = categ_sum(i,j,-2) + categ_sum(i,j,-1) + categ_sum(i,j,0) + categ_pftsum
1076 ! if ( abs(allsum/area-1.0_RP) > EPS ) then
1077 ! if( IO_L ) write(IO_FID_LOG,*) i,j,allsum/area
1078 ! endif
1079 
1080  do p = 1, landuse_pft_nmax
1081  pft(p) = categ_sum(i,j,p)
1082  pft_idx(p) = p
1083  enddo
1084 
1085  ! bubble sort
1086  do p = 1, landuse_pft_nmax-1
1087  do pp = p+1, landuse_pft_nmax
1088  if ( pft(pp) > pft(p) ) then
1089  temp = pft(p)
1090  pft(p) = pft(pp)
1091  pft(pp) = temp
1092 
1093  temp_idx = pft_idx(p)
1094  pft_idx(p) = pft_idx(pp)
1095  pft_idx(pp) = temp_idx
1096  endif
1097  enddo
1098  enddo
1099 
1100  categ_pftsum = sum( pft(:) )
1101 
1102  ! overwrite existing map (GLCCv2) by LU100
1103  allsum = categ_sum(i,j,-2) + categ_sum(i,j,0) + categ_pftsum
1104  mask = 0.5_rp + sign( 0.5_rp, allsum-eps ) ! if any land data is found, overwrite
1105 
1106  ! land fraction : 1 - ocean / total
1107  allsum = categ_sum(i,j,-2) + categ_sum(i,j,-1) + categ_sum(i,j,0) + categ_pftsum
1108  zerosw = 0.5_rp - sign( 0.5_rp, allsum-eps )
1109  frac = ( 1.0_rp-zerosw ) - categ_sum(i,j,-1) * ( 1.0_rp-zerosw ) / ( allsum-zerosw )
1110  landuse_frac_land(i,j) = ( mask ) * frac & ! overwrite
1111  + ( 1.0_rp-mask ) * landuse_frac_land(i,j) ! keep existing value
1112 
1113  ! lake fraction : lake / ( total - ocean )
1114  allsum = categ_sum(i,j,-2) + categ_sum(i,j,0) + categ_pftsum
1115  zerosw = 0.5_rp - sign( 0.5_rp, allsum-eps )
1116  frac = categ_sum(i,j,-2) * ( 1.0_rp-zerosw ) / ( allsum-zerosw )
1117  landuse_frac_lake(i,j) = ( mask ) * frac & ! overwrite
1118  + ( 1.0_rp-mask ) * landuse_frac_lake(i,j) ! keep existing value
1119 
1120  ! urban fraction : urban / ( total - ocean - lake )
1121  allsum = categ_sum(i,j,0) + categ_pftsum
1122  zerosw = 0.5_rp - sign( 0.5_rp, allsum-eps )
1123  frac = categ_sum(i,j,0) * ( 1.0_rp-zerosw ) / ( allsum-zerosw )
1124  landuse_frac_urban(i,j) = ( mask ) * frac & ! overwrite
1125  + ( 1.0_rp-mask ) * landuse_frac_urban(i,j) ! keep existing value
1126 
1127  ! PFT fraction : PFT / sum( PFT(1:mosaic) )
1128  allsum = sum( pft(1:landuse_pft_mosaic) )
1129  if ( abs(allsum) > eps ) then
1130  do p = 1, landuse_pft_mosaic
1131  landuse_frac_pft(i,j,p) = ( mask ) * pft(p) / allsum & ! overwrite
1132  + ( 1.0_rp-mask ) * landuse_frac_pft(i,j,p) ! keep existing value
1133  landuse_index_pft(i,j,p) = ( mask ) * pft_idx(p) & ! overwrite
1134  + ( 1.0_rp-mask ) * landuse_index_pft(i,j,p) ! keep existing value
1135  enddo
1136  ! if no second PFT, set to same as PFT1
1137  if ( abs(landuse_frac_pft(i,j,1)-1.0_rp) <= eps ) then
1138  landuse_frac_pft(i,j,:) = 0.0_rp
1139  landuse_frac_pft(i,j,1) = 1.0_rp
1140  landuse_index_pft(i,j,:) = pft_idx(1)
1141  endif
1142  else ! if no PFT, set to bare ground
1143  do p = 1, landuse_pft_mosaic
1144  landuse_frac_pft(i,j,p) = ( mask ) * 0.0_rp & ! overwrite
1145  + ( 1.0_rp-mask ) * landuse_frac_pft(i,j,p) ! keep existing value
1146  landuse_index_pft(i,j,p) = ( mask ) * 1 & ! overwrite
1147  + ( 1.0_rp-mask ) * landuse_index_pft(i,j,p) ! keep existing value
1148  enddo
1149  landuse_frac_pft(i,j,1) = ( mask ) * 1.0_rp & ! overwrite
1150  + ( 1.0_rp-mask ) * landuse_frac_pft(i,j,1) ! keep existing value
1151  endif
1152 
1153  enddo
1154  enddo
1155 
1156  if ( limit_urban_fraction < 1.0_rp ) then
1157  do j = js, je
1158  do i = is, ie
1159  if ( landuse_frac_urban(i,j) == 1.0_rp ) then ! if no PFT, set to grassland
1160  landuse_frac_pft(i,j,:) = 0.0_rp
1161  landuse_frac_pft(i,j,1) = 1.0_rp
1162  landuse_index_pft(i,j,:) = 2
1163  endif
1164  landuse_frac_urban(i,j) = min( landuse_frac_urban(i,j), limit_urban_fraction )
1165  enddo
1166  enddo
1167  endif
1168 
1169  return
1170  end subroutine cnvlanduse_lu100m
1171 
1172  !-----------------------------------------------------------------------------
1174  subroutine cnvlanduse_jibis
1175  implicit none
1176  !---------------------------------------------------------------------------
1177 
1178  return
1179  end subroutine cnvlanduse_jibis
1180 
1181 end module mod_cnvlanduse
scale_grid_index::is
integer, public is
start point of inner domain: x, local
Definition: scale_grid_index.F90:60
scale_grid_index::je
integer, public je
end point of inner domain: y, local
Definition: scale_grid_index.F90:63
scale_landuse::landuse_frac_pft
real(rp), dimension(:,:,:), allocatable, public landuse_frac_pft
fraction of PFT for each mosaic
Definition: scale_landuse.F90:46
scale_process::prc_mpistop
subroutine, public prc_mpistop
Abort MPI.
Definition: scale_process.F90:267
scale_stdio::io_l
logical, public io_l
output log or not? (this process)
Definition: scale_stdio.F90:61
scale_const::const_radius
real(rp), public const_radius
radius of the planet [m]
Definition: scale_const.F90:46
scale_stdio
module STDIO
Definition: scale_stdio.F90:12
scale_landuse::landuse_calc_fact
subroutine, public landuse_calc_fact
Definition: scale_landuse.F90:156
mod_cnvlanduse::cnvlanduse
subroutine, public cnvlanduse
Driver.
Definition: mod_cnvlanduse.f90:163
scale_const::const_d2r
real(rp), public const_d2r
degree to radian
Definition: scale_const.F90:35
scale_landuse::landuse_frac_urban
real(rp), dimension(:,:), allocatable, public landuse_frac_urban
urban fraction
Definition: scale_landuse.F90:41
mod_cnvlanduse::cnvlanduse_useglccv2
logical, public cnvlanduse_useglccv2
Definition: mod_cnvlanduse.f90:35
scale_grid_index
module grid index
Definition: scale_grid_index.F90:14
scale_stdio::io_nml
logical, public io_nml
output log or not? (for namelist, this process)
Definition: scale_stdio.F90:62
scale_tracer
module TRACER
Definition: scale_tracer.F90:15
scale_grid_index::ia
integer, public ia
of whole cells: x, local, with HALO
Definition: scale_grid_index.F90:55
scale_stdio::io_get_available_fid
integer function, public io_get_available_fid()
search & get available file ID
Definition: scale_stdio.F90:313
scale_grid_real
module GRID (real space)
Definition: scale_grid_real.F90:14
scale_landuse
module LANDUSE
Definition: scale_landuse.F90:11
scale_grid_real::real_dlon
real(rp), dimension(:,:), allocatable, public real_dlon
delta longitude
Definition: scale_grid_real.F90:52
scale_landuse::landuse_write
subroutine, public landuse_write
Write landuse data.
Definition: scale_landuse.F90:271
scale_comm
module COMMUNICATION
Definition: scale_comm.F90:23
scale_grid_index::js
integer, public js
start point of inner domain: y, local
Definition: scale_grid_index.F90:62
scale_process
module PROCESS
Definition: scale_process.F90:14
scale_landuse::landuse_pft_mosaic
integer, public landuse_pft_mosaic
number of PFT mosaic
Definition: scale_landuse.F90:43
scale_grid_real::real_dlat
real(rp), dimension(:,:), allocatable, public real_dlat
delta latitude
Definition: scale_grid_real.F90:53
scale_landuse::landuse_index_pft
integer, dimension(:,:,:), allocatable, public landuse_index_pft
index of PFT for each mosaic
Definition: scale_landuse.F90:47
scale_const
module CONSTANT
Definition: scale_const.F90:14
mod_cnvlanduse
module Convert LandUseIndex
Definition: mod_cnvlanduse.f90:10
scale_prof
module profiler
Definition: scale_prof.F90:10
scale_grid_index::ie
integer, public ie
end point of inner domain: x, local
Definition: scale_grid_index.F90:61
mod_cnvlanduse::cnvlanduse_setup
subroutine, public cnvlanduse_setup
Setup.
Definition: mod_cnvlanduse.f90:59
scale_const::const_eps
real(rp), public const_eps
small number
Definition: scale_const.F90:36
scale_landuse::landuse_frac_lake
real(rp), dimension(:,:), allocatable, public landuse_frac_lake
lake fraction
Definition: scale_landuse.F90:40
mod_cnvlanduse::cnvlanduse_uselu100m
logical, public cnvlanduse_uselu100m
Definition: mod_cnvlanduse.f90:36
scale_precision
module PRECISION
Definition: scale_precision.F90:13
mod_cnvlanduse::cnvlanduse_usejibis
logical, public cnvlanduse_usejibis
Definition: mod_cnvlanduse.f90:37
scale_const::const_pi
real(rp), public const_pi
pi
Definition: scale_const.F90:34
scale_landuse::landuse_pft_nmax
integer, public landuse_pft_nmax
number of plant functional type(PFT)
Definition: scale_landuse.F90:44
scale_stdio::io_fid_conf
integer, public io_fid_conf
Config file ID.
Definition: scale_stdio.F90:55
scale_stdio::io_fid_log
integer, public io_fid_log
Log file ID.
Definition: scale_stdio.F90:56
scale_precision::rp
integer, parameter, public rp
Definition: scale_precision.F90:43
mod_cnvlanduse::cnvlanduse_donothing
logical, public cnvlanduse_donothing
Definition: mod_cnvlanduse.f90:34
scale_grid_real::real_lonx
real(rp), dimension(:,:), allocatable, public real_lonx
longitude at staggered point (uy) [rad,0-2pi]
Definition: scale_grid_real.F90:46
scale_landuse::landuse_frac_land
real(rp), dimension(:,:), allocatable, public landuse_frac_land
land fraction
Definition: scale_landuse.F90:39
scale_stdio::io_fid_nml
integer, public io_fid_nml
Log file ID (only for output namelist)
Definition: scale_stdio.F90:57
scale_grid_real::real_laty
real(rp), dimension(:,:), allocatable, public real_laty
latitude at staggered point (xv) [rad,-pi,pi]
Definition: scale_grid_real.F90:50
scale_grid_index::ja
integer, public ja
of whole cells: y, local, with HALO
Definition: scale_grid_index.F90:56
Generated by   doxygen 1.8.14