1 | ! (c) British Crown Copyright 2008, the Met Office. |
---|
2 | ! All rights reserved. |
---|
3 | ! $Revision: 23 $, $Date: 2011-03-31 15:41:37 +0200 (jeu. 31 mars 2011) $ |
---|
4 | ! $URL: http://cfmip-obs-sim.googlecode.com/svn/stable/v1.4.0/cosp_utils.F90 $ |
---|
5 | |
---|
6 | ! Redistribution and use in source and binary forms, with or without modification, are permitted |
---|
7 | ! provided that the following conditions are met: |
---|
8 | |
---|
9 | ! * Redistributions of source code must retain the above copyright notice, this list |
---|
10 | ! of conditions and the following disclaimer. |
---|
11 | ! * Redistributions in binary form must reproduce the above copyright notice, this list |
---|
12 | ! of conditions and the following disclaimer in the documentation and/or other materials |
---|
13 | ! provided with the distribution. |
---|
14 | ! * Neither the name of the Met Office nor the names of its contributors may be used |
---|
15 | ! to endorse or promote products derived from this software without specific prior written |
---|
16 | ! permission. |
---|
17 | |
---|
18 | ! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR |
---|
19 | ! IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND |
---|
20 | ! FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR |
---|
21 | ! CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL |
---|
22 | ! DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, |
---|
23 | ! DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER |
---|
24 | ! IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT |
---|
25 | ! OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. |
---|
26 | |
---|
27 | ! History: |
---|
28 | ! Jul 2007 - A. Bodas-Salcedo - Initial version |
---|
29 | |
---|
30 | MODULE MOD_COSP_UTILS |
---|
31 | USE MOD_COSP_CONSTANTS |
---|
32 | IMPLICIT NONE |
---|
33 | |
---|
34 | INTERFACE Z_TO_DBZ |
---|
35 | MODULE PROCEDURE Z_TO_DBZ_2D,Z_TO_DBZ_3D,Z_TO_DBZ_4D |
---|
36 | END INTERFACE |
---|
37 | |
---|
38 | INTERFACE COSP_CHECK_INPUT |
---|
39 | MODULE PROCEDURE COSP_CHECK_INPUT_1D,COSP_CHECK_INPUT_2D,COSP_CHECK_INPUT_3D |
---|
40 | END INTERFACE |
---|
41 | CONTAINS |
---|
42 | |
---|
43 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
44 | !------------------- SUBROUTINE COSP_PRECIP_MXRATIO -------------- |
---|
45 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
46 | SUBROUTINE COSP_PRECIP_MXRATIO(Npoints,Nlevels,Ncolumns,p,T,prec_frac,prec_type, & |
---|
47 | n_ax,n_bx,alpha_x,c_x,d_x,g_x,a_x,b_x,gamma1,gamma2,gamma3,gamma4, & |
---|
48 | flux,mxratio,reff) |
---|
49 | |
---|
50 | ! Input arguments, (IN) |
---|
51 | integer,intent(in) :: Npoints,Nlevels,Ncolumns |
---|
52 | real,intent(in),dimension(Npoints,Nlevels) :: p,T,flux |
---|
53 | real,intent(in),dimension(Npoints,Ncolumns,Nlevels) :: prec_frac |
---|
54 | real,intent(in) :: n_ax,n_bx,alpha_x,c_x,d_x,g_x,a_x,b_x,gamma1,gamma2,gamma3,gamma4,prec_type |
---|
55 | ! Input arguments, (OUT) |
---|
56 | real,intent(out),dimension(Npoints,Ncolumns,Nlevels) :: mxratio |
---|
57 | real,intent(inout),dimension(Npoints,Ncolumns,Nlevels) :: reff |
---|
58 | ! Local variables |
---|
59 | integer :: i,j,k |
---|
60 | REAL :: sigma,one_over_xip1,xi,rho0,rho,lambda_x,gamma_4_3_2,delta |
---|
61 | REAL :: seuil |
---|
62 | |
---|
63 | if (ok_debug_cosp) then |
---|
64 | seuil=1.e-15 |
---|
65 | else |
---|
66 | seuil=0.0 |
---|
67 | endif |
---|
68 | |
---|
69 | mxratio = 0.0 |
---|
70 | |
---|
71 | if (n_ax >= 0.0) then ! N_ax is used to control which hydrometeors need to be computed |
---|
72 | xi = d_x/(alpha_x + b_x - n_bx + 1.0) |
---|
73 | rho0 = 1.29 |
---|
74 | sigma = (gamma2/(gamma1*c_x))*(n_ax*a_x*gamma2)**xi |
---|
75 | one_over_xip1 = 1.0/(xi + 1.0) |
---|
76 | gamma_4_3_2 = 0.5*gamma4/gamma3 |
---|
77 | delta = (alpha_x + b_x + d_x - n_bx + 1.0) |
---|
78 | |
---|
79 | DO k=1,Nlevels |
---|
80 | DO j=1,Ncolumns |
---|
81 | DO i=1,Npoints |
---|
82 | if ((prec_frac(i,j,k)==prec_type).or.(prec_frac(i,j,k)==3.)) then |
---|
83 | rho = p(i,k)/(287.05*T(i,k)) |
---|
84 | mxratio(i,j,k)=(flux(i,k)*((rho/rho0)**g_x)*sigma)**one_over_xip1 |
---|
85 | mxratio(i,j,k)=mxratio(i,j,k)/rho |
---|
86 | ! Compute effective radius |
---|
87 | ! if ((reff(i,j,k) <= 0.0).and.(flux(i,k) /= 0.0)) then |
---|
88 | if ((reff(i,j,k) <= 0.0).and.(flux(i,k) > seuil)) then |
---|
89 | lambda_x = (a_x*c_x*((rho0/rho)**g_x)*n_ax*gamma1/flux(i,k))**(1./delta) |
---|
90 | reff(i,j,k) = gamma_4_3_2/lambda_x |
---|
91 | endif |
---|
92 | endif |
---|
93 | enddo |
---|
94 | enddo |
---|
95 | enddo |
---|
96 | endif |
---|
97 | END SUBROUTINE COSP_PRECIP_MXRATIO |
---|
98 | |
---|
99 | |
---|
100 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
101 | !------------------- SUBROUTINE ZERO_INT ------------------------- |
---|
102 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
103 | ELEMENTAL SUBROUTINE ZERO_INT(x,y01,y02,y03,y04,y05,y06,y07,y08,y09,y10, & |
---|
104 | y11,y12,y13,y14,y15,y16,y17,y18,y19,y20, & |
---|
105 | y21,y22,y23,y24,y25,y26,y27,y28,y29,y30) |
---|
106 | |
---|
107 | integer,intent(inout) :: x |
---|
108 | integer,intent(inout),optional :: y01,y02,y03,y04,y05,y06,y07,y08,y09,y10, & |
---|
109 | y11,y12,y13,y14,y15,y16,y17,y18,y19,y20, & |
---|
110 | y21,y22,y23,y24,y25,y26,y27,y28,y29,y30 |
---|
111 | x = 0 |
---|
112 | if (present(y01)) y01 = 0 |
---|
113 | if (present(y02)) y02 = 0 |
---|
114 | if (present(y03)) y03 = 0 |
---|
115 | if (present(y04)) y04 = 0 |
---|
116 | if (present(y05)) y05 = 0 |
---|
117 | if (present(y06)) y06 = 0 |
---|
118 | if (present(y07)) y07 = 0 |
---|
119 | if (present(y08)) y08 = 0 |
---|
120 | if (present(y09)) y09 = 0 |
---|
121 | if (present(y10)) y10 = 0 |
---|
122 | if (present(y11)) y11 = 0 |
---|
123 | if (present(y12)) y12 = 0 |
---|
124 | if (present(y13)) y13 = 0 |
---|
125 | if (present(y14)) y14 = 0 |
---|
126 | if (present(y15)) y15 = 0 |
---|
127 | if (present(y16)) y16 = 0 |
---|
128 | if (present(y17)) y17 = 0 |
---|
129 | if (present(y18)) y18 = 0 |
---|
130 | if (present(y19)) y19 = 0 |
---|
131 | if (present(y20)) y20 = 0 |
---|
132 | if (present(y21)) y21 = 0 |
---|
133 | if (present(y22)) y22 = 0 |
---|
134 | if (present(y23)) y23 = 0 |
---|
135 | if (present(y24)) y24 = 0 |
---|
136 | if (present(y25)) y25 = 0 |
---|
137 | if (present(y26)) y26 = 0 |
---|
138 | if (present(y27)) y27 = 0 |
---|
139 | if (present(y28)) y28 = 0 |
---|
140 | if (present(y29)) y29 = 0 |
---|
141 | if (present(y30)) y30 = 0 |
---|
142 | END SUBROUTINE ZERO_INT |
---|
143 | |
---|
144 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
145 | !------------------- SUBROUTINE ZERO_REAL ------------------------ |
---|
146 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
147 | ELEMENTAL SUBROUTINE ZERO_REAL(x,y01,y02,y03,y04,y05,y06,y07,y08,y09,y10, & |
---|
148 | y11,y12,y13,y14,y15,y16,y17,y18,y19,y20, & |
---|
149 | y21,y22,y23,y24,y25,y26,y27,y28,y29,y30) |
---|
150 | |
---|
151 | real,intent(inout) :: x |
---|
152 | real,intent(inout),optional :: y01,y02,y03,y04,y05,y06,y07,y08,y09,y10, & |
---|
153 | y11,y12,y13,y14,y15,y16,y17,y18,y19,y20, & |
---|
154 | y21,y22,y23,y24,y25,y26,y27,y28,y29,y30 |
---|
155 | x = 0.0 |
---|
156 | if (present(y01)) y01 = 0.0 |
---|
157 | if (present(y02)) y02 = 0.0 |
---|
158 | if (present(y03)) y03 = 0.0 |
---|
159 | if (present(y04)) y04 = 0.0 |
---|
160 | if (present(y05)) y05 = 0.0 |
---|
161 | if (present(y06)) y06 = 0.0 |
---|
162 | if (present(y07)) y07 = 0.0 |
---|
163 | if (present(y08)) y08 = 0.0 |
---|
164 | if (present(y09)) y09 = 0.0 |
---|
165 | if (present(y10)) y10 = 0.0 |
---|
166 | if (present(y11)) y11 = 0.0 |
---|
167 | if (present(y12)) y12 = 0.0 |
---|
168 | if (present(y13)) y13 = 0.0 |
---|
169 | if (present(y14)) y14 = 0.0 |
---|
170 | if (present(y15)) y15 = 0.0 |
---|
171 | if (present(y16)) y16 = 0.0 |
---|
172 | if (present(y17)) y17 = 0.0 |
---|
173 | if (present(y18)) y18 = 0.0 |
---|
174 | if (present(y19)) y19 = 0.0 |
---|
175 | if (present(y20)) y20 = 0.0 |
---|
176 | if (present(y21)) y21 = 0.0 |
---|
177 | if (present(y22)) y22 = 0.0 |
---|
178 | if (present(y23)) y23 = 0.0 |
---|
179 | if (present(y24)) y24 = 0.0 |
---|
180 | if (present(y25)) y25 = 0.0 |
---|
181 | if (present(y26)) y26 = 0.0 |
---|
182 | if (present(y27)) y27 = 0.0 |
---|
183 | if (present(y28)) y28 = 0.0 |
---|
184 | if (present(y29)) y29 = 0.0 |
---|
185 | if (present(y30)) y30 = 0.0 |
---|
186 | END SUBROUTINE ZERO_REAL |
---|
187 | |
---|
188 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
189 | !--------------------- SUBROUTINE Z_TO_DBZ_2D -------------------- |
---|
190 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
191 | SUBROUTINE Z_TO_DBZ_2D(mdi,z) |
---|
192 | real,intent(in) :: mdi |
---|
193 | real,dimension(:,:),intent(inout) :: z |
---|
194 | ! Reflectivity Z: |
---|
195 | ! Input in [m3] |
---|
196 | ! Output in dBZ, with Z in [mm6 m-3] |
---|
197 | |
---|
198 | ! 1.e18 to convert from [m3] to [mm6 m-3] |
---|
199 | z = 1.e18*z |
---|
200 | where (z > 1.0e-6) ! Limit to -60 dBZ |
---|
201 | z = 10.0*log10(z) |
---|
202 | elsewhere |
---|
203 | z = mdi |
---|
204 | end where |
---|
205 | END SUBROUTINE Z_TO_DBZ_2D |
---|
206 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
207 | !--------------------- SUBROUTINE Z_TO_DBZ_3D -------------------- |
---|
208 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
209 | SUBROUTINE Z_TO_DBZ_3D(mdi,z) |
---|
210 | real,intent(in) :: mdi |
---|
211 | real,dimension(:,:,:),intent(inout) :: z |
---|
212 | ! Reflectivity Z: |
---|
213 | ! Input in [m3] |
---|
214 | ! Output in dBZ, with Z in [mm6 m-3] |
---|
215 | |
---|
216 | ! 1.e18 to convert from [m3] to [mm6 m-3] |
---|
217 | z = 1.e18*z |
---|
218 | where (z > 1.0e-6) ! Limit to -60 dBZ |
---|
219 | z = 10.0*log10(z) |
---|
220 | elsewhere |
---|
221 | z = mdi |
---|
222 | end where |
---|
223 | END SUBROUTINE Z_TO_DBZ_3D |
---|
224 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
225 | !--------------------- SUBROUTINE Z_TO_DBZ_4D -------------------- |
---|
226 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
227 | SUBROUTINE Z_TO_DBZ_4D(mdi,z) |
---|
228 | real,intent(in) :: mdi |
---|
229 | real,dimension(:,:,:,:),intent(inout) :: z |
---|
230 | ! Reflectivity Z: |
---|
231 | ! Input in [m3] |
---|
232 | ! Output in dBZ, with Z in [mm6 m-3] |
---|
233 | |
---|
234 | ! 1.e18 to convert from [m3] to [mm6 m-3] |
---|
235 | z = 1.e18*z |
---|
236 | where (z > 1.0e-6) ! Limit to -60 dBZ |
---|
237 | z = 10.0*log10(z) |
---|
238 | elsewhere |
---|
239 | z = mdi |
---|
240 | end where |
---|
241 | END SUBROUTINE Z_TO_DBZ_4D |
---|
242 | |
---|
243 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
244 | !----------------- SUBROUTINES COSP_CHECK_INPUT_1D --------------- |
---|
245 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
246 | SUBROUTINE COSP_CHECK_INPUT_1D(vname,x,min_val,max_val) |
---|
247 | character(len=*) :: vname |
---|
248 | real,intent(inout) :: x(:) |
---|
249 | real,intent(in),optional :: min_val,max_val |
---|
250 | logical :: l_min,l_max |
---|
251 | character(len=128) :: pro_name='COSP_CHECK_INPUT_1D' |
---|
252 | |
---|
253 | l_min=.false. |
---|
254 | l_max=.false. |
---|
255 | |
---|
256 | if (present(min_val)) then |
---|
257 | ! if (x < min_val) x = min_val |
---|
258 | if (any(x < min_val)) then |
---|
259 | l_min = .true. |
---|
260 | where (x < min_val) |
---|
261 | x = min_val |
---|
262 | end where |
---|
263 | endif |
---|
264 | endif |
---|
265 | if (present(max_val)) then |
---|
266 | ! if (x > max_val) x = max_val |
---|
267 | if (any(x > max_val)) then |
---|
268 | l_max = .true. |
---|
269 | where (x > max_val) |
---|
270 | x = max_val |
---|
271 | end where |
---|
272 | endif |
---|
273 | endif |
---|
274 | |
---|
275 | if (l_min) print *,'----- WARNING: '//trim(pro_name)//': minimum value of '//trim(vname)//' set to: ',min_val |
---|
276 | if (l_max) print *,'----- WARNING: '//trim(pro_name)//': maximum value of '//trim(vname)//' set to: ',max_val |
---|
277 | END SUBROUTINE COSP_CHECK_INPUT_1D |
---|
278 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
279 | !----------------- SUBROUTINES COSP_CHECK_INPUT_2D --------------- |
---|
280 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
281 | SUBROUTINE COSP_CHECK_INPUT_2D(vname,x,min_val,max_val) |
---|
282 | character(len=*) :: vname |
---|
283 | real,intent(inout) :: x(:,:) |
---|
284 | real,intent(in),optional :: min_val,max_val |
---|
285 | logical :: l_min,l_max |
---|
286 | character(len=128) :: pro_name='COSP_CHECK_INPUT_2D' |
---|
287 | |
---|
288 | l_min=.false. |
---|
289 | l_max=.false. |
---|
290 | |
---|
291 | if (present(min_val)) then |
---|
292 | ! if (x < min_val) x = min_val |
---|
293 | if (any(x < min_val)) then |
---|
294 | l_min = .true. |
---|
295 | where (x < min_val) |
---|
296 | x = min_val |
---|
297 | end where |
---|
298 | endif |
---|
299 | endif |
---|
300 | if (present(max_val)) then |
---|
301 | ! if (x > max_val) x = max_val |
---|
302 | if (any(x > max_val)) then |
---|
303 | l_max = .true. |
---|
304 | where (x > max_val) |
---|
305 | x = max_val |
---|
306 | end where |
---|
307 | endif |
---|
308 | endif |
---|
309 | |
---|
310 | if (l_min) print *,'----- WARNING: '//trim(pro_name)//': minimum value of '//trim(vname)//' set to: ',min_val |
---|
311 | if (l_max) print *,'----- WARNING: '//trim(pro_name)//': maximum value of '//trim(vname)//' set to: ',max_val |
---|
312 | END SUBROUTINE COSP_CHECK_INPUT_2D |
---|
313 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
314 | !----------------- SUBROUTINES COSP_CHECK_INPUT_3D --------------- |
---|
315 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
316 | SUBROUTINE COSP_CHECK_INPUT_3D(vname,x,min_val,max_val) |
---|
317 | character(len=*) :: vname |
---|
318 | real,intent(inout) :: x(:,:,:) |
---|
319 | real,intent(in),optional :: min_val,max_val |
---|
320 | logical :: l_min,l_max |
---|
321 | character(len=128) :: pro_name='COSP_CHECK_INPUT_3D' |
---|
322 | |
---|
323 | l_min=.false. |
---|
324 | l_max=.false. |
---|
325 | |
---|
326 | if (present(min_val)) then |
---|
327 | ! if (x < min_val) x = min_val |
---|
328 | if (any(x < min_val)) then |
---|
329 | l_min = .true. |
---|
330 | where (x < min_val) |
---|
331 | x = min_val |
---|
332 | end where |
---|
333 | endif |
---|
334 | endif |
---|
335 | if (present(max_val)) then |
---|
336 | ! if (x > max_val) x = max_val |
---|
337 | if (any(x > max_val)) then |
---|
338 | l_max = .true. |
---|
339 | where (x > max_val) |
---|
340 | x = max_val |
---|
341 | end where |
---|
342 | endif |
---|
343 | endif |
---|
344 | |
---|
345 | if (l_min) print *,'----- WARNING: '//trim(pro_name)//': minimum value of '//trim(vname)//' set to: ',min_val |
---|
346 | if (l_max) print *,'----- WARNING: '//trim(pro_name)//': maximum value of '//trim(vname)//' set to: ',max_val |
---|
347 | END SUBROUTINE COSP_CHECK_INPUT_3D |
---|
348 | |
---|
349 | |
---|
350 | END MODULE MOD_COSP_UTILS |
---|