phylmd/Orography/start_init_orog.f

MODULE start_init_orog_m

  ! From startvar.F, version 1.4
  ! 2006/01/27 15:14:22 Fairhead

  IMPLICIT NONE

  REAL, ALLOCATABLE, SAVE:: mask(:, :) ! fraction of land (iim + 1, jjm + 1)
  REAL, ALLOCATABLE, SAVE:: phis(:, :) ! surface geopotential, in m2 s-2

CONTAINS

  SUBROUTINE start_init_orog(relief, zstd_2d, zsig_2d, zgam_2d, zthe_2d, &
       zpic_2d, zval_2d)

    use conf_dat2d_m, only: conf_dat2d
    use comgeom, only: rlatu, rlonv
    use dimens_m, only: iim, jjm
    use grid_noro_m, only: grid_noro
    use indicesol, only: epsfra
    use netcdf, only: nf90_nowrite
    use netcdf95, only: nf95_open, nf95_gw_var, nf95_inq_varid, nf95_close
    use nr_util, only: pi

    REAL, intent(out):: relief(:, :) ! (iim + 1, jjm + 1) orographie moyenne

    REAL, intent(out):: zstd_2d(:, :) ! (iim + 1, jjm + 1)
    ! (deviation standard de l'orographie sous-maille)

    REAL, intent(out):: zsig_2d(:, :) ! (iim + 1, jjm + 1)
    ! (pente de l'orographie sous-maille)
    
    REAL, intent(out):: zgam_2d(:, :) ! (iim + 1, jjm + 1)
    ! (anisotropie de l'orographie sous maille)

    REAL, intent(out):: zthe_2d(:, :) ! (iim + 1, jjm + 1)
    ! (orientation de l'axe oriente dans la direction de plus grande
    ! pente de l'orographie sous maille)

    REAL, intent(out):: zpic_2d(:, :) ! (iim + 1, jjm + 1)
    ! hauteur pics de la SSO

    REAL, intent(out):: zval_2d(:, :) ! (iim + 1, jjm + 1)
    ! hauteur vallees de la SSO

    ! Local:

    INTEGER iml_rel
    INTEGER jml_rel
    INTEGER ncid, varid
    REAL, pointer:: relief_hi(:, :)
    REAL, ALLOCATABLE:: lon_rad(:), lat_rad(:)
    REAL, pointer:: lon_ini(:), lat_ini(:)

    !-----------------------------------

    print *, "Call sequence information: start_init_orog"

    if (any((/size(relief, 1), size(zstd_2d, 1), size(zsig_2d, 1), &
         size(zgam_2d, 1), size(zthe_2d, 1), size(zpic_2d, 1), &
         size(zval_2d, 1)/) /= iim + 1)) stop "start_init_orog size 1"
    if (any((/size(relief, 2), size(zstd_2d, 2), size(zsig_2d, 2), &
         size(zgam_2d, 2), size(zthe_2d, 2), size(zpic_2d, 2), &
         size(zval_2d, 2)/) /= jjm + 1)) stop "start_init_orog size 2"

    print *, 'Reading the high resolution orography...'

    call nf95_open('Relief.nc', nf90_nowrite, ncid)

    call nf95_inq_varid(ncid, "longitude", varid)
    call nf95_gw_var(ncid, varid, lon_ini)
    lon_ini = lon_ini * pi / 180. ! convert to rad
    iml_rel = size(lon_ini)

    call nf95_inq_varid(ncid, "latitude", varid)
    call nf95_gw_var(ncid, varid, lat_ini)
    lat_ini = lat_ini * pi / 180. ! convert to rad
    jml_rel = size(lat_ini)

    call nf95_inq_varid(ncid, "RELIEF", varid)
    call nf95_gw_var(ncid, varid, relief_hi)

    call nf95_close(ncid)

    ALLOCATE(lon_rad(iml_rel))
    ALLOCATE(lat_rad(jml_rel))

    CALL conf_dat2d(lon_ini, lat_ini, lon_rad, lat_rad, relief_hi , &
         interbar=.FALSE.)
    deallocate(lon_ini, lat_ini) ! pointers

    print *, 'Compute all the parameters needed for the gravity wave drag code'

    ! Interpolated fields:
    ALLOCATE(phis(iim + 1, jjm + 1))
    ALLOCATE(mask(iim + 1, jjm + 1))

    CALL grid_noro(lon_rad, lat_rad, relief_hi, rlonv, rlatu, phis, relief, &
         zstd_2d, zsig_2d, zgam_2d, zthe_2d, zpic_2d, zval_2d, mask)
    deallocate(relief_hi) ! pointer

    phis(iim + 1, :) = phis(1, :)
    phis(:, :) = phis(:, :) * 9.81

    mask(2:, 1) = mask(1, 1) ! north pole
    mask(2:, jjm + 1) = mask(1, jjm + 1) ! south pole
    mask(iim + 1, 2:jjm) = mask(1, 2:jjm) ! Greenwich
    WHERE (mask < EPSFRA)
       mask = 0.
    elsewhere (1. - mask < EPSFRA)
       mask = 1.
    endwhere

  END SUBROUTINE start_init_orog

END MODULE start_init_orog_m
1	guez	3	MODULE start_init_orog_m
2
3			! From startvar.F, version 1.4
4			! 2006/01/27 15:14:22 Fairhead
5
6			IMPLICIT NONE
7
8	guez	15	REAL, ALLOCATABLE, SAVE:: mask(:, :) ! fraction of land (iim + 1, jjm + 1)
9	guez	3	REAL, ALLOCATABLE, SAVE:: phis(:, :) ! surface geopotential, in m2 s-2
10
11			CONTAINS
12
13			SUBROUTINE start_init_orog(relief, zstd_2d, zsig_2d, zgam_2d, zthe_2d, &
14			zpic_2d, zval_2d)
15
16			use conf_dat2d_m, only: conf_dat2d
17			use comgeom, only: rlatu, rlonv
18			use dimens_m, only: iim, jjm
19	guez	39	use grid_noro_m, only: grid_noro
20	guez	3	use indicesol, only: epsfra
21	guez	43	use netcdf, only: nf90_nowrite
22			use netcdf95, only: nf95_open, nf95_gw_var, nf95_inq_varid, nf95_close
23	guez	39	use nr_util, only: pi
24	guez	3
25	guez	42	REAL, intent(out):: relief(:, :) ! (iim + 1, jjm + 1) orographie moyenne
26	guez	3
27	guez	42	REAL, intent(out):: zstd_2d(:, :) ! (iim + 1, jjm + 1)
28	guez	3	! (deviation standard de l'orographie sous-maille)
29
30	guez	42	REAL, intent(out):: zsig_2d(:, :) ! (iim + 1, jjm + 1)
31	guez	3	! (pente de l'orographie sous-maille)
32
33	guez	42	REAL, intent(out):: zgam_2d(:, :) ! (iim + 1, jjm + 1)
34	guez	3	! (anisotropie de l'orographie sous maille)
35
36	guez	42	REAL, intent(out):: zthe_2d(:, :) ! (iim + 1, jjm + 1)
37	guez	3	! (orientation de l'axe oriente dans la direction de plus grande
38			! pente de l'orographie sous maille)
39
40	guez	42	REAL, intent(out):: zpic_2d(:, :) ! (iim + 1, jjm + 1)
41			! hauteur pics de la SSO
42	guez	3
43	guez	42	REAL, intent(out):: zval_2d(:, :) ! (iim + 1, jjm + 1)
44			! hauteur vallees de la SSO
45
46	guez	3	! Local:
47
48	guez	42	INTEGER iml_rel
49			INTEGER jml_rel
50	guez	43	INTEGER ncid, varid
51			REAL, pointer:: relief_hi(:, :)
52	guez	3	REAL, ALLOCATABLE:: lon_rad(:), lat_rad(:)
53	guez	43	REAL, pointer:: lon_ini(:), lat_ini(:)
54	guez	3
55			!-----------------------------------
56
57			print *, "Call sequence information: start_init_orog"
58
59			if (any((/size(relief, 1), size(zstd_2d, 1), size(zsig_2d, 1), &
60			size(zgam_2d, 1), size(zthe_2d, 1), size(zpic_2d, 1), &
61			size(zval_2d, 1)/) /= iim + 1)) stop "start_init_orog size 1"
62			if (any((/size(relief, 2), size(zstd_2d, 2), size(zsig_2d, 2), &
63			size(zgam_2d, 2), size(zthe_2d, 2), size(zpic_2d, 2), &
64			size(zval_2d, 2)/) /= jjm + 1)) stop "start_init_orog size 2"
65
66	guez	43	print *, 'Reading the high resolution orography...'
67	guez	3
68	guez	43	call nf95_open('Relief.nc', nf90_nowrite, ncid)
69	guez	3
70	guez	43	call nf95_inq_varid(ncid, "longitude", varid)
71			call nf95_gw_var(ncid, varid, lon_ini)
72			lon_ini = lon_ini * pi / 180. ! convert to rad
73			iml_rel = size(lon_ini)
74	guez	3
75	guez	43	call nf95_inq_varid(ncid, "latitude", varid)
76			call nf95_gw_var(ncid, varid, lat_ini)
77			lat_ini = lat_ini * pi / 180. ! convert to rad
78			jml_rel = size(lat_ini)
79	guez	3
80	guez	43	call nf95_inq_varid(ncid, "RELIEF", varid)
81			call nf95_gw_var(ncid, varid, relief_hi)
82	guez	3
83	guez	43	call nf95_close(ncid)
84	guez	3
85	guez	43	ALLOCATE(lon_rad(iml_rel))
86	guez	3	ALLOCATE(lat_rad(jml_rel))
87
88			CALL conf_dat2d(lon_ini, lat_ini, lon_rad, lat_rad, relief_hi , &
89			interbar=.FALSE.)
90	guez	43	deallocate(lon_ini, lat_ini) ! pointers
91	guez	3
92			print *, 'Compute all the parameters needed for the gravity wave drag code'
93
94	guez	68	! Interpolated fields:
95	guez	3	ALLOCATE(phis(iim + 1, jjm + 1))
96	guez	15	ALLOCATE(mask(iim + 1, jjm + 1))
97	guez	3
98			CALL grid_noro(lon_rad, lat_rad, relief_hi, rlonv, rlatu, phis, relief, &
99	guez	15	zstd_2d, zsig_2d, zgam_2d, zthe_2d, zpic_2d, zval_2d, mask)
100	guez	43	deallocate(relief_hi) ! pointer
101	guez	3
102			phis(iim + 1, :) = phis(1, :)
103			phis(:, :) = phis(:, :) * 9.81
104
105	guez	15	mask(2:, 1) = mask(1, 1) ! north pole
106			mask(2:, jjm + 1) = mask(1, jjm + 1) ! south pole
107			mask(iim + 1, 2:jjm) = mask(1, 2:jjm) ! Greenwich
108			WHERE (mask < EPSFRA)
109			mask = 0.
110			elsewhere (1. - mask < EPSFRA)
111			mask = 1.
112	guez	3	endwhere
113
114			END SUBROUTINE start_init_orog
115
116			END MODULE start_init_orog_m