da_detsurtyp.inc
References to this file elsewhere.
1 subroutine da_detsurtyp (snow,ice,landsea, vegtyp, soiltyp, &
2 is, ie, js, je, &
3 i, j, dx, dy, dxm, dym, &
4 isflg,ob_vegtyp,ob_soiltyp, &
5 seap, icep, lndp, snop)
6
7 !---------------------------------------------------------------------------
8 ! Purpose: determine surface type at obs locations.
9 !
10 ! METHOD: using the background information
11 ! 1. get the background landmask/snow/sea-ice
12 ! 2. calculate percentage of sea/ice/land/snow
13 ! 3. determine surface type at obs location
14 ! 4. find nearest grid vegtype and soil type
15 !---------------------------------------------------------------------------
16
17 implicit none
18
19 integer, intent(in) :: is, ie, js, je
20 integer, intent(in) :: i, j
21 real , intent(in) :: dx, dxm, dy, dym
22 real , intent(in) :: snow(is:ie,js:je) ! Input variable
23 real , intent(in) :: ice(is:ie,js:je) ! Input variable
24 real , intent(in) :: landsea(is:ie,js:je) ! Input variable
25 integer, intent(in) :: vegtyp(is:ie,js:je) ! Input variable
26 integer, intent(in) :: soiltyp(is:ie,js:je) ! Input variable
27 integer, intent(out) :: isflg ! Output variable
28 real, intent(out) :: ob_vegtyp ! Output variable
29 real, intent(out) :: ob_soiltyp ! Output variable
30 real, intent(out) :: seap, icep, lndp, snop ! percentage of surface type
31
32 ! isflg - surface flag
33 ! 0 sea
34 ! 1 sea ice
35 ! 2 land
36 ! 3 snow
37 ! 4 mixed predominately sea
38 ! 5 mixed predominately sea ice
39 ! 6 mixed predominately land
40 ! 7 mixed predominately snow
41
42 ! local variables
43 integer :: n, xbflag(4) ! surface type at xb location
44 ! 0:sea 1:sea-ice 2:land 3:snow
45 real :: w(4),minw ! weight at 4 xb locations
46
47 ! 1.0 determine surface type of xb at 4 location around obs
48 !-------------------------------------------------------
49 if ( nint(landsea(i ,j )) == 0 ) xbflag(1) = 0 ! sea
50 if ( nint(landsea(i+1,j )) == 0 ) xbflag(2) = 0 ! sea
51 if ( nint(landsea(i ,j+1)) == 0 ) xbflag(3) = 0 ! sea
52 if ( nint(landsea(i+1,j+1)) == 0 ) xbflag(4) = 0 ! sea
53
54 if ( nint(landsea(i ,j )) == 1 ) xbflag(1) = 2 ! land/snow/sea-ice
55 if ( nint(landsea(i+1,j )) == 1 ) xbflag(2) = 2 ! land/snow/sea-ice
56 if ( nint(landsea(i ,j+1)) == 1 ) xbflag(3) = 2 ! land/snow/sea-ice
57 if ( nint(landsea(i+1,j+1)) == 1 ) xbflag(4) = 2 ! land/snow/sea-ice
58
59 if ( nint(snow(i ,j )) == 1 ) xbflag(1) = 3 ! snow
60 if ( nint(snow(i+1,j )) == 1 ) xbflag(2) = 3 ! snow
61 if ( nint(snow(i ,j+1)) == 1 ) xbflag(3) = 3 ! snow
62 if ( nint(snow(i+1,j+1)) == 1 ) xbflag(4) = 3 ! snow
63
64 if ( nint(ice(i ,j )) == 1 ) xbflag(1) = 1 ! sea-ice
65 if ( nint(ice(i+1,j )) == 1 ) xbflag(2) = 1 ! sea-ice
66 if ( nint(ice(i ,j+1)) == 1 ) xbflag(3) = 1 ! sea-ice
67 if ( nint(ice(i+1,j+1)) == 1 ) xbflag(4) = 1 ! sea-ice
68
69 ! 2.0 determine surface type percentage at obs location
70 !------------------------------------------------------
71 ! (i,j+1)
72 ! -----------------(i+1,j+1)
73 ! | w2 |
74 ! |dym w1 |
75 ! | obs |
76 ! |------- * |
77 ! |dy w4 | w3 |
78 ! | dx | dxm |
79 ! |----------------
80 ! (i,j) (i+1,j)
81 !
82 !--------------------------------------------------------
83
84 seap = 0.
85 icep = 0.
86 snop = 0.
87 lndp = 0.
88 w(1) = dym*dxm ! weight for point (i,j)
89 w(2) = dym*dx ! weight for point (i+1,j)
90 w(3) = dy *dxm ! weight for point (i,j+1)
91 w(4) = dy *dx ! weight for point (i+1,j+1)
92
93 do n = 1, 4
94 if (xbflag(n) == 0) seap = seap+w(n)
95 if (xbflag(n) == 1) icep = icep+w(n)
96 if (xbflag(n) == 2) lndp = lndp+w(n)
97 if (xbflag(n) == 3) snop = snop+w(n)
98 end do
99
100 ! fo2d = dym*(dxm*fm2d(i,j ) + dx*fm2d(i+1,j )) &
101 ! + dy *(dxm*fm2d(i,j+1) + dx*fm2d(i+1,j+1))
102
103 ! 3.0 determine final surface flag at obs location
104 !-----------------------------------------
105 if (seap >= 0.99) isflg = 0
106 if (icep >= 0.99) isflg = 1
107 if (lndp >= 0.99) isflg = 2
108 if (snop >= 0.99) isflg = 3
109 if ( .not. (seap >= 0.99) .and. &
110 .not. (icep >= 0.99) .and. &
111 .not. (lndp >= 0.99) .and. &
112 .not. (snop >= 0.99) ) then
113 if (seap > lndp) then
114 if (seap > icep) then
115 if (seap > snop) then
116 isflg = 4
117 else
118 isflg = 7
119 end if
120 else
121 if (icep > snop) then
122 isflg = 5
123 else
124 isflg = 7
125 end if
126 end if
127 else
128 if (lndp > icep) then
129 if (lndp > snop) then
130 isflg = 6
131 else
132 isflg = 7
133 end if
134 else
135 if (icep > snop) then
136 isflg = 5
137 else
138 isflg = 7
139 end if
140 end if
141 end if
142 end if
143
144 ! 4.0 find nearest grid vegtype and soil type
145 ! at obs location
146 !-----------------------------------------
147 minw=min(w(1),w(2),w(3),w(4))
148 if (minw == w(1)) then
149 ob_vegtyp = float(vegtyp (i+1,j+1))
150 ob_soiltyp = float(soiltyp(i+1,j+1))
151 else if (minw == w(2)) then
152 ob_vegtyp = float(vegtyp (i,j+1))
153 ob_soiltyp = float(soiltyp(i,j+1))
154 else if (minw == w(3)) then
155 ob_vegtyp = float(vegtyp (i+1,j))
156 ob_soiltyp = float(soiltyp(i+1,j))
157 else if (minw == w(4)) then
158 ob_vegtyp = float(vegtyp (i,j))
159 ob_soiltyp = float(soiltyp(i,j))
160 end if
161
162 return
163
164 end subroutine da_detsurtyp
165