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