Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

scops.f90 @ 5408

Last change on this file since 5408 was 5268, checked in by abarral, 8 weeks ago
.f90 <-> .F90 depending on cpp key use
File size: 10.9 KB

Line
1	! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2	! Copyright (c) 2009, British Crown Copyright, the Met Office
3	! All rights reserved.
4	!
5	! Redistribution and use in source and binary forms, with or without modification, are
6	! permitted provided that the following conditions are met:
7	!
8	! 1. Redistributions of source code must retain the above copyright notice, this list of
9	! conditions and the following disclaimer.
10	!
11	! 2. 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
13	! materials provided with the distribution.
14	!
15	! 3. Neither the name of the copyright holder nor the names of its contributors may be
16	! used to endorse or promote products derived from this software without specific prior
17	! written permission.
18	!
19	! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY
20	! EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF
21	! MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL
22	! THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
23	! SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT
24	! OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
25	! INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
26	! LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
27	! OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
28	!
29	! History
30	! May 2015 - D. Swales - Modified for COSPv2.0
31	! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
32	module mod_scops
33	USE COSP_KINDS, ONLY: wp
34	USE MOD_RNG!, ONLY: rng_state,get_rng
35	use mod_cosp_error, ONLY: errorMessage
36
37	implicit none
38
39	integer,parameter :: default_overlap = 3 ! Used when invalid overlap assumption is provided.
40
41	contains
42	subroutine scops(npoints,nlev,ncol,rngs,cc,conv,overlap,frac_out,ncolprint)
43	INTEGER :: npoints, & ! Number of model points in the horizontal
44	nlev, & ! Number of model levels in column
45	ncol, & ! Number of subcolumns
46	overlap ! Overlap type (1=max, 2=rand, 3=max/rand)
47	type(rng_state),dimension(npoints) :: rngs
48	INTEGER, parameter :: huge32 = 2147483647
49	INTEGER, parameter :: i2_16 = 65536
50	INTEGER :: i,j,ilev,ibox,ncolprint,ilev2
51
52	REAL(WP), dimension(npoints,nlev) :: &
53	cc, & ! Input cloud cover in each model level (fraction)
54	! NOTE: This is the HORIZONTAL area of each
55	! grid box covered by clouds
56	conv, & ! Input convective cloud cover in each model level (fraction)
57	! NOTE: This is the HORIZONTAL area of each
58	! grid box covered by convective clouds
59	tca ! Total cloud cover in each model level (fraction)
60	! with extra layer of zeroes on top
61	! in this version this just contains the values input
62	! from cc but with an extra level
63	REAL(WP),intent(inout), dimension(npoints,ncol,nlev) :: &
64	frac_out ! Boxes gridbox divided up into equivalent of BOX in
65	! original version, but indexed by column then row, rather than
66	! by row then column
67	REAL(WP), dimension(npoints,ncol) :: &
68	threshold, & ! pointer to position in gridbox
69	maxocc, & ! Flag for max overlapped conv cld
70	maxosc, & ! Flag for max overlapped strat cld
71	boxpos, & ! ordered pointer to position in gridbox
72	threshold_min ! minimum value to define range in with new threshold is chosen.
73	REAL(WP), dimension(npoints) :: &
74	ran ! vector of random numbers
75
76	! Test for valid input overlap assumption
77	if (overlap .ne. 1 .and. overlap .ne. 2 .and. overlap .ne. 3) then
78	overlap=default_overlap
79	call errorMessage('ERROR(scops): Invalid overlap assumption provided. Using default overlap assumption (max/ran)')
80	endif
81
82	boxpos = spread(([(i, i=1,ncol)]-0.5)/ncol,1,npoints)
83
84	! #######################################################################
85	! Initialize working variables
86	! #######################################################################
87
88	! Initialize frac_out to zero
89	frac_out(1:npoints,1:ncol,1:nlev)=0.0
90
91	! Assign 2d tca array using 1d input array cc
92	tca(1:npoints,1:nlev) = cc(1:npoints,1:nlev)
93
94	if (ncolprint.ne.0) then
95	write (6,'(a)') 'frac_out_pp_rev:'
96	do j=1,npoints,1000
97	write(6,'(a10)') 'j='
98	write(6,'(8I10)') j
99	write (6,'(8f5.2)') ((frac_out(j,ibox,ilev),ibox=1,ncolprint),ilev=1,nlev)
100	enddo
101	write (6,'(a)') 'ncol:'
102	write (6,'(I3)') ncol
103	endif
104	if (ncolprint.ne.0) then
105	write (6,'(a)') 'last_frac_pp:'
106	do j=1,npoints,1000
107	write(6,'(a10)') 'j='
108	write(6,'(8I10)') j
109	write (6,'(8f5.2)') (tca(j,1))
110	enddo
111	endif
112
113	! #######################################################################
114	! ALLOCATE CLOUD INTO BOXES, FOR NCOLUMNS, NLEVELS
115	! frac_out is the array that contains the information
116	! where 0 is no cloud, 1 is a stratiform cloud and 2 is a
117	! convective cloud
118	! #######################################################################
119
120	! Loop over vertical levels
121	DO ilev = 1,nlev
122
123	! Initialise threshold
124	IF (ilev.eq.1) then
125	! If max overlap
126	IF (overlap.eq.1) then
127	! Select pixels spread evenly across the gridbox
128	threshold(1:npoints,1:ncol)=boxpos(1:npoints,1:ncol)
129	ELSE
130	DO ibox=1,ncol
131	!include 'congvec.f90'
132	ran(1:npoints) = get_rng(RNGS)
133	! select random pixels from the non-convective
134	! part the gridbox ( some will be converted into
135	! convective pixels below )
136	threshold(1:npoints,ibox) = conv(1:npoints,ilev)+(1-conv(1:npoints,ilev))*ran(npoints)
137	enddo
138	ENDIF
139	IF (ncolprint.ne.0) then
140	write (6,'(a)') 'threshold_nsf2:'
141	do j=1,npoints,1000
142	write(6,'(a10)') 'j='
143	write(6,'(8I10)') j
144	write (6,'(8f5.2)') (threshold(j,ibox),ibox=1,ncolprint)
145	enddo
146	ENDIF
147	ENDIF
148
149	IF (ncolprint.ne.0) then
150	write (6,'(a)') 'ilev:'
151	write (6,'(I2)') ilev
152	ENDIF
153
154	DO ibox=1,ncol
155	! All versions
156	!maxocc(1:npoints,ibox) = merge(1,0,boxpos(1:npoints,ibox) .le. conv(1:npoints,ilev))
157	!maxocc(1:npoints,ibox) = merge(1,0, conv(1:npoints,ilev) .gt. boxpos(1:npoints,ibox))
158	do j=1,npoints
159	if (boxpos(j,ibox).le.conv(j,ilev)) then
160	maxocc(j,ibox) = 1
161	else
162	maxocc(j,ibox) = 0
163	end if
164	enddo
165
166	! Max overlap
167	if (overlap.eq.1) then
168	threshold_min(1:npoints,ibox) = conv(1:npoints,ilev)
169	maxosc(1:npoints,ibox) = 1
170	endif
171
172	! Random overlap
173	if (overlap.eq.2) then
174	threshold_min(1:npoints,ibox) = conv(1:npoints,ilev)
175	maxosc(1:npoints,ibox) = 0
176	endif
177	! Max/Random overlap
178	if (overlap.eq.3) then
179	! DS2014 START: The bounds on tca are not valid when ilev=1.
180	!threshold_min(1:npoints,ibox) = max(conv(1:npoints,ilev),min(tca(1:npoints,ilev-1),tca(1:npoints,ilev)))
181	!maxosc(1:npoints,ibox) = merge(1,0,threshold(1:npoints,ibox) .lt. &
182	! min(tca(1:npoints,ilev-1),tca(1:npoints,ilev)) .and. &
183	! (threshold(1:npoints,ibox).gt.conv(1:npoints,ilev)))
184	if (ilev .ne. 1) then
185	threshold_min(1:npoints,ibox) = max(conv(1:npoints,ilev),min(tca(1:npoints,ilev-1),tca(1:npoints,ilev)))
186	maxosc(1:npoints,ibox) = merge(1,0,threshold(1:npoints,ibox) .lt. &
187	min(tca(1:npoints,ilev-1),tca(1:npoints,ilev)) .and. &
188	(threshold(1:npoints,ibox).gt.conv(1:npoints,ilev)))
189	else
190	threshold_min(1:npoints,ibox) = max(conv(1:npoints,ilev),min(0._wp,tca(1:npoints,ilev)))
191	maxosc(1:npoints,ibox) = merge(1,0,threshold(1:npoints,ibox) .lt. &
192	min(0._wp,tca(1:npoints,ilev)) .and. &
193	(threshold(1:npoints,ibox).gt.conv(1:npoints,ilev)))
194	endif
195	endif
196
197	! Reset threshold
198	!include 'congvec.f90'
199	ran(1:npoints) = get_rng(RNGS)
200
201	threshold(1:npoints,ibox)= maxocc(1:npoints,ibox)*(boxpos(1:npoints,ibox)) + &
202	(1-maxocc(1:npoints,ibox))((maxosc(1:npoints,ibox))(threshold(1:npoints,ibox)) + &
203	(1-maxosc(1:npoints,ibox))*(threshold_min(1:npoints,ibox)+ &
204	(1-threshold_min(1:npoints,ibox))*ran(1:npoints)))
205
206	! Fill frac_out with 1's where tca is greater than the threshold
207	frac_out(1:npoints,ibox,ilev) = merge(1,0,tca(1:npoints,ilev).gt.threshold(1:npoints,ibox))
208
209	! Code to partition boxes into startiform and convective parts goes here
210	where(threshold(1:npoints,ibox).le.conv(1:npoints,ilev) .and. conv(1:npoints,ilev).gt.0.) frac_out(1:npoints,ibox,ilev)=2
211	ENDDO ! ibox
212
213
214	! Set last_frac to tca at this level, so as to be tca from last level next time round
215	if (ncolprint.ne.0) then
216	do j=1,npoints ,1000
217	write(6,'(a10)') 'j='
218	write(6,'(8I10)') j
219	write (6,'(a)') 'last_frac:'
220	write (6,'(8f5.2)') (tca(j,ilev))
221	write (6,'(a)') 'conv:'
222	write (6,'(8f5.2)') (conv(j,ilev),ibox=1,ncolprint)
223	write (6,'(a)') 'max_overlap_cc:'
224	write (6,'(8f5.2)') (maxocc(j,ibox),ibox=1,ncolprint)
225	write (6,'(a)') 'max_overlap_sc:'
226	write (6,'(8f5.2)') (maxosc(j,ibox),ibox=1,ncolprint)
227	write (6,'(a)') 'threshold_min_nsf2:'
228	write (6,'(8f5.2)') (threshold_min(j,ibox),ibox=1,ncolprint)
229	write (6,'(a)') 'threshold_nsf2:'
230	write (6,'(8f5.2)') (threshold(j,ibox),ibox=1,ncolprint)
231	write (6,'(a)') 'frac_out_pp_rev:'
232	write (6,'(8f5.2)') ((frac_out(j,ibox,ilev2),ibox=1,ncolprint),ilev2=1,nlev)
233	enddo
234	endif
235
236	enddo ! Loop over nlev
237
238	! END
239	end subroutine scops
240	end module mod_scops

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format