Context Navigation

MISR_simulator.F90 @ 5186

Last change on this file since 5186 was 3358, checked in by idelkadi, 6 years ago
Implementation de la nouvelle version COSPv2 dans LMDZ. Pour compiler avec makelmdz_fcma utiliser l'option "-cosp2 true"
File size: 13.0 KB

Line
1	! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2	! Copyright (c) 2009, Roger Marchand, version 1.2
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_MISR_SIMULATOR
33	use cosp_kinds, only: wp
34	use MOD_COSP_STATS, ONLY: hist2D
35	use mod_cosp_config, ONLY: R_UNDEF,numMISRHgtBins,numMISRTauBins,misr_histHgt, &
36	misr_histTau
37	implicit none
38
39	! Parameters
40	real(wp),parameter :: &
41	misr_taumin = 0.3_wp, & ! Minimum optical depth for joint-histogram
42	tauchk = -1.*log(0.9999999) ! Lower limit on optical depth
43
44	contains
45
46	! ######################################################################################
47	! SUBROUTINE misr_subcolumn
48	! ######################################################################################
49	SUBROUTINE MISR_SUBCOLUMN(npoints,ncol,nlev,dtau,zfull,at,sunlit,tauOUT, &
50	dist_model_layertops,box_MISR_ztop)
51	! INPUTS
52	INTEGER, intent(in) :: &
53	npoints, & ! Number of horizontal gridpoints
54	ncol, & ! Number of subcolumns
55	nlev ! Number of vertical layers
56	INTEGER, intent(in),dimension(npoints) :: &
57	sunlit ! 1 for day points, 0 for night time
58	REAL(WP),intent(in),dimension(npoints,ncol,nlev) :: &
59	dtau ! Optical thickness
60	REAL(WP),intent(in),dimension(npoints,nlev) :: &
61	zfull, & ! Height of full model levels (i.e. midpoints), [nlev] is bottom
62	at ! Temperature (K)
63
64	! OUTPUTS
65	REAL(WP),intent(out),dimension(npoints,ncol) :: &
66	box_MISR_ztop, & ! Cloud-top height in each column
67	tauOUT ! Optical depth in each column
68	REAL(WP),intent(out),dimension(npoints,numMISRHgtBins) :: &
69	dist_model_layertops !
70
71	! INTERNAL VARIABLES
72	INTEGER :: ilev,j,loop,ibox,thres_crossed_MISR
73	INTEGER :: iMISR_ztop
74	REAL(WP) :: cloud_dtau,MISR_penetration_height,ztest
75
76	! ############################################################################
77	! Initialize
78	box_MISR_ztop(1:npoints,1:ncol) = 0._wp
79
80	do j=1,npoints
81
82	! Estimate distribution of Model layer tops
83	dist_model_layertops(j,:)=0
84	do ilev=1,nlev
85	! Define location of "layer top"
86	if(ilev.eq.1 .or. ilev.eq.nlev) then
87	ztest=zfull(j,ilev)
88	else
89	ztest=0.5_wp*(zfull(j,ilev)+zfull(j,ilev-1))
90	endif
91
92	! Find MISR layer that contains this level
93	! NOTE the first MISR level is "no height" level
94	iMISR_ztop=2
95	do loop=2,numMISRHgtBins
96	if ( ztest .gt. 1000*misr_histHgt(loop+1) ) then
97	iMISR_ztop=loop+1
98	endif
99	enddo
100
101	dist_model_layertops(j,iMISR_ztop) = dist_model_layertops(j,iMISR_ztop)+1
102	enddo
103
104	! For each GCM cell or horizontal model grid point
105	do ibox=1,ncol
106	! Compute optical depth as a cummulative distribution in the vertical (nlev).
107	tauOUT(j,ibox)=sum(dtau(j,ibox,1:nlev))
108
109	thres_crossed_MISR=0
110	do ilev=1,nlev
111	! If there a cloud, start the counter and store this height
112	if(thres_crossed_MISR .eq. 0 .and. dtau(j,ibox,ilev) .gt. 0.) then
113	! First encountered a "cloud"
114	thres_crossed_MISR = 1
115	cloud_dtau = 0
116	endif
117
118	if( thres_crossed_MISR .lt. 99 .and. thres_crossed_MISR .gt. 0 ) then
119	if( dtau(j,ibox,ilev) .eq. 0.) then
120	! We have come to the end of the current cloud layer without yet
121	! selecting a CTH boundary. Restart cloud tau counter
122	cloud_dtau=0
123	else
124	! Add current optical depth to count for the current cloud layer
125	cloud_dtau=cloud_dtau+dtau(j,ibox,ilev)
126	endif
127
128	! If the cloud is continuous but optically thin (< 1) from above the
129	! current layer cloud top to the current level then MISR will like
130	! see a top below the top of the current layer.
131	if( dtau(j,ibox,ilev).gt.0 .and. (cloud_dtau-dtau(j,ibox,ilev)) .lt. 1) then
132	if(dtau(j,ibox,ilev) .lt. 1 .or. ilev.eq.1 .or. ilev.eq.nlev) then
133	! MISR will likely penetrate to some point within this layer ... the middle
134	MISR_penetration_height=zfull(j,ilev)
135	else
136	! Take the OD = 1.0 level into this layer
137	MISR_penetration_height=0.5_wp*(zfull(j,ilev)+zfull(j,ilev-1)) - &
138	0.5_wp*(zfull(j,ilev-1)-zfull(j,ilev+1))/dtau(j,ibox,ilev)
139	endif
140	box_MISR_ztop(j,ibox)=MISR_penetration_height
141	endif
142
143	! Check for a distinctive water layer
144	if(dtau(j,ibox,ilev) .gt. 1 .and. at(j,ilev) .gt. 273 ) then
145	! Must be a water cloud, take this as CTH level
146	thres_crossed_MISR=99
147	endif
148
149	! If the total column optical depth is "large" than MISR can't see
150	! anything else. Set current point as CTH level
151	if(sum(dtau(j,ibox,1:ilev)) .gt. 5) then
152	thres_crossed_MISR=99
153	endif
154	endif
155	enddo
156
157	! Check to see if there was a cloud for which we didn't
158	! set a MISR cloud top boundary
159	if( thres_crossed_MISR .eq. 1) then
160	! If the cloud has a total optical depth of greater
161	! than ~ 0.5 MISR will still likely pick up this cloud
162	! with a height near the true cloud top
163	! otherwise there should be no CTH
164	if(sum(dtau(j,ibox,1:nlev)) .gt. 0.5) then
165	! keep MISR detected CTH
166	elseif(sum(dtau(j,ibox,1:nlev)) .gt. 0.2) then
167	! MISR may detect but wont likley have a good height
168	box_MISR_ztop(j,ibox)=-1
169	else
170	! MISR not likely to even detect.
171	! so set as not cloudy
172	box_MISR_ztop(j,ibox)=0
173	endif
174	endif
175	enddo ! loop of subcolumns
176
177	enddo ! loop of gridpoints
178
179	! Modify MISR CTH for satellite spatial / pattern matcher effects
180	! Code in this region added by roj 5/2006 to account
181	! for spatial effect of the MISR pattern matcher.
182	! Basically, if a column is found between two neighbors
183	! at the same CTH, and that column has no hieght or
184	! a lower CTH, THEN misr will tend to but place the
185	! odd column at the same height as it neighbors.
186
187	! This setup assumes the columns represent a about a 1 to 4 km scale
188	! it will need to be modified significantly, otherwise
189	! ! DS2015: Add loop over gridpoints and index accordingly.
190	! if(ncol.eq.1) then
191	! ! Adjust based on neightboring points.
192	! do j=2,npoints-1
193	! if(box_MISR_ztop(j-1,1) .gt. 0 .and. &
194	! box_MISR_ztop(j+1,1) .gt. 0 .and. &
195	! abs(box_MISR_ztop(j-1,1)-box_MISR_ztop(j+1,1)) .lt. 500 .and. &
196	! box_MISR_ztop(j,1) .lt. box_MISR_ztop(j+1,1)) then
197	! box_MISR_ztop(j,1) = box_MISR_ztop(j+1,1)
198	! endif
199	! enddo
200	! else
201	! ! Adjust based on neighboring subcolumns.
202	! do j=1,npoints
203	! do ibox=2,ncol-1
204	! if(box_MISR_ztop(j,ibox-1) .gt. 0 .and. &
205	! box_MISR_ztop(j,ibox+1) .gt. 0 .and. &
206	! abs(box_MISR_ztop(j,ibox-1)-box_MISR_ztop(j,ibox+1)) .lt. 500 .and. &
207	! box_MISR_ztop(j,ibox) .lt. box_MISR_ztop(j,ibox+1)) then
208	! box_MISR_ztop(j,ibox) = box_MISR_ztop(j,ibox+1)
209	! endif
210	! enddo
211	! enddo
212	! endif
213	! ! DS2015 END
214
215	! Fill dark scenes
216	do j=1,numMISRHgtBins
217	where(sunlit .ne. 1) dist_model_layertops(1:npoints,j) = R_UNDEF
218	enddo
219
220	end SUBROUTINE MISR_SUBCOLUMN
221
222	! ######################################################################################
223	! SUBROUTINE misr_column
224	! ######################################################################################
225	SUBROUTINE MISR_COLUMN(npoints,ncol,box_MISR_ztop,sunlit,tau,MISR_cldarea,MISR_mean_ztop,fq_MISR_TAU_v_CTH)
226
227	! INPUTS
228	INTEGER, intent(in) :: &
229	npoints, & ! Number of horizontal gridpoints
230	ncol ! Number of subcolumns
231	INTEGER, intent(in),dimension(npoints) :: &
232	sunlit ! 1 for day points, 0 for night time
233	REAL(WP),intent(in),dimension(npoints,ncol) :: &
234	box_MISR_ztop, & ! Cloud-top height in each column
235	tau ! Column optical thickness
236
237	! OUTPUTS
238	REAL(WP),intent(inout),dimension(npoints) :: &
239	MISR_cldarea, & ! Fraction area covered by clouds
240	MISR_mean_ztop ! Mean cloud top height MISR would observe
241	REAL(WP),intent(inout),dimension(npoints,7,numMISRHgtBins) :: &
242	fq_MISR_TAU_v_CTH ! Joint histogram of cloud-cover and tau
243
244	! INTERNAL VARIABLES
245	INTEGER :: j
246	LOGICAL,dimension(ncol) :: box_cloudy
247	real(wp),dimension(npoints,ncol) :: tauWRK,box_MISR_ztopWRK
248	! ############################################################################
249
250	! Compute column quantities and joint-histogram
251	MISR_cldarea(1:npoints) = 0._wp
252	MISR_mean_ztop(1:npoints) = 0._wp
253	fq_MISR_TAU_v_CTH(1:npoints,1:7,1:numMISRHgtBins) = 0._wp
254	tauWRK(1:npoints,1:ncol) = tau(1:npoints,1:ncol)
255	box_MISR_ztopWRK(1:npoints,1:ncol) = box_MISR_ztop(1:npoints,1:ncol)
256	do j=1,npoints
257
258	! Subcolumns that are cloudy(true) and not(false)
259	box_cloudy(1:ncol) = merge(.true.,.false.,tau(j,1:ncol) .gt. tauchk)
260
261	! Fill optically thin clouds with fill value
262	where(.not. box_cloudy(1:ncol)) tauWRK(j,1:ncol) = -999._wp
263	where(box_MISR_ztopWRK(j,1:ncol) .eq. 0) box_MISR_ztopWRK(j,1:ncol)=-999._wp
264
265	! Compute joint histogram and column quantities for points that are sunlit and cloudy
266	if (sunlit(j) .eq. 1) then
267	! Joint histogram
268	call hist2D(tauWRK(j,1:ncol),box_MISR_ztopWRK(j,1:ncol),ncol,misr_histTau,numMISRTauBins,&
269	1000*misr_histHgt,numMISRHgtBins,fq_MISR_TAU_v_CTH(j,1:numMISRTauBins,1:numMISRHgtBins))
270	fq_MISR_TAU_v_CTH(j,1:numMISRTauBins,1:numMISRHgtBins) = &
271	100._wp*fq_MISR_TAU_v_CTH(j,1:numMISRTauBins,1:numMISRHgtBins)/ncol
272
273	! Column cloud area
274	MISR_cldarea(j)=real(count(box_MISR_ztopWRK(j,1:ncol) .ne. -999.))/ncol
275
276	! Column cloud-top height
277	if ( count(box_MISR_ztopWRK(j,1:ncol) .ne. -999.) .ne. 0 ) then
278	MISR_mean_ztop(j) = sum(box_MISR_ztopWRK(j,1:ncol),box_MISR_ztopWRK(j,1:ncol) .ne. -999.)/ &
279	count(box_MISR_ztopWRK(j,1:ncol) .ne. -999.)
280	else
281	MISR_mean_ztop(j) = R_UNDEF
282	endif
283
284	else
285	MISR_cldarea(j) = R_UNDEF
286	MISR_mean_ztop(npoints) = R_UNDEF
287	endif
288	enddo
289
290	end SUBROUTINE MISR_COLUMN
291
292	end MODULE MOD_MISR_SIMULATOR

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

Download in other formats:

Original Format