Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

meso_dustopacity.F @ 3094

Last change on this file since 3094 was 11, checked in by aslmd, 14 years ago
spiga@svn-planeto:ajoute le modele meso-echelle martien
File size: 9.6 KB

Line
1	SUBROUTINE meso_dustopacity(ngrid,nlayer,nq,
2	$ zday,pplay,pplev,ls,pq,
3	$ tauref,tau,aerosol)
4
5	IMPLICIT NONE
6
7	c=======================================================================
8	c
9	c CAREFUL: THIS IS A VERSION TO BE USED WITH WRF !!!
10	c
11	c ... CHECK THE ****WRF lines
12	c
13	c=======================================================================
14	c subject:
15	c --------
16	c Computing aerosol optical depth (dust opacity)
17	c In each layers
18	c
19	c author: F.Forget
20	c ------
21	c update F. Montmessin (water ice scheme)
22	c and S. Lebonnois (12/06/2003) compatibility dust/ice/chemistry
23	c
24	c input:
25	c -----
26	c ngrid Number of gridpoint of horizontal grid
27	c nlayer Number of layer
28	c nq Number of tracer
29	c ls Solar longitude (Ls) , radian
30	c pplay,pplev pressure (Pa) in the middle and boundary of each layer
31	c pq Dust mixing ratio (used if tracer =T and active=T).
32	c
33	c output:
34	c -------
35	c tauref Prescribed mean column optical depth at 700 Pa
36	c tau Column total visible dust optical depth at each point
37	c aerosol aerosol(ig,l,1) is the dust optical
38	c depth in layer l, grid point ig
39	c
40	c=======================================================================
41	#include "dimensions.h"
42	#include "dimphys.h"
43	#include "callkeys.h"
44	#include "comcstfi.h"
45	#include "comgeomfi.h"
46	#include "dimradmars.h"
47	#include "yomaer.h"
48	#include "tracer.h"
49	#include "planete.h"
50
51	c-----------------------------------------------------------------------
52	c
53	c Declarations :
54	c --------------
55	c
56	c Input/Output
57	c ------------
58	INTEGER ngrid,nlayer,nq
59
60	REAL ls,zday,expfactor
61	REAL pplev(ngrid,nlayer+1),pplay(ngrid,nlayer)
62	REAL pq(ngrid,nlayer,nq)
63	REAL tauref(ngrid), tau(ngrid,naerkind)
64	REAL aerosol(ngrid,nlayer,naerkind)
65
66	c
67	c Local variables :
68	c -----------------
69	INTEGER l,ig,iq
70	!****WRF
71	INTEGER ierr
72	!****WRF
73	real topdust(ngridmx)
74	real zlsconst, zp
75	real taueq,tauS,tauN
76	real r0,reff,coefsize
77	c
78	c local saved variables
79	c ---------------------
80
81	REAL topdust0(ngridmx)
82	SAVE topdust0
83
84	LOGICAL firstcall
85	DATA firstcall/.true./
86	SAVE firstcall
87
88
89	c----------------------------------------------------------------------
90
91	c Initialisation
92	c --------------
93
94	IF (firstcall) THEN
95
96	c altitude of the top of the aerosol layer (km) at Ls=2.76rad:
97	c in the Viking year scenario
98	DO ig=1,ngrid
99	topdust0(ig)=60. -22.SIN(lati(ig))*2
100	END DO
101	firstcall=.false.
102	END IF
103
104
105
106	c -------------------------------------------------------------
107	c 1) Prescribed dust (if tracer=F or active=F)
108	c -------------------------------------------------------------
109	IF ((.not.tracer) .or. (.not.active)) THEN
110
111	c Vertical column optical depth at 700.Pa
112	c ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
113	IF(iaervar.eq.1) THEN
114	!****WRF
115	!****WRF: additional ASCII file to define dust opacity
116	OPEN(99,file='dustopacity.def',status='old',form='formatted'
117	. ,iostat=ierr)
118	IF(ierr.NE.0) THEN
119	write(,) 'No file dustopacity.def'
120	stop
121	ELSE
122	READ(99,*) tauvis
123	write(,) 'constant tau', tauvis
124	!****WRF
125	!****WRF
126	do ig=1, ngridmx
127	tauref(ig)=max(tauvis,1.e-9) ! tauvis=cste as read in file
128	end do
129	CLOSE(99)
130	ENDIF
131	ELSE IF (iaervar.eq.2) THEN ! << "Viking" Scenario>>
132
133	tauref(1) = 0.7+.3cos(ls+80.pi/180.) ! like seen by VL1
134	do ig=2,ngrid
135	tauref(ig) = tauref(1)
136	end do
137
138	ELSE IF (iaervar.eq.3) THEN ! << "MGS" scenario >>
139
140	taueq= 0.2 +(0.5-0.2) (cos(0.5(ls-4.363)))**14
141	tauS= 0.1 +(0.5-0.1) (cos(0.5(ls-4.363)))**14
142	tauN = 0.1
143	c if (peri_day.eq.150) then
144	c tauS=0.1
145	c tauN=0.1 +(0.5-0.1) (cos(0.5(ls+pi-4.363)))**14
146	c taueq= 0.2 +(0.5-0.2) (cos(0.5(ls+pi-4.363)))**14
147	c endif
148
149	do ig=1,ngrid/2 ! Northern hemisphere
150	tauref(ig)= tauN +
151	& (taueq-tauN)0.5(1+tanh((45-lati(ig)180./pi)6/60))
152	end do
153	do ig=ngrid/2+1, ngridmx ! Southern hemisphere
154	tauref(ig)= tauS +
155	& (taueq-tauS)0.5(1+tanh((45+lati(ig)180./pi)6/60))
156	end do
157
158	ELSE IF (iaervar.eq.4) THEN ! << "TES" scenario >>
159	c*************WRF
160	c
161	c 1. MY24 TES
162	c
163	call readtesassim(ngrid,nlayer,zday,pplev,tauref)
164
165	ELSE IF (iaervar.gt.99) THEN ! << input netcdf file >>
166	c*************WRF
167	c
168	c 2. customized dust opacity field
169	c ex: from assimilation
170	call meso_readtesassim(ngrid,nlayer,zday,pplev,tauref,
171	. iaervar)
172	c
173	c*************WRF
174
175
176	ELSE
177	stop 'problem with iaervar in dustopacity.F'
178	ENDIF
179
180
181	c Altitude of the top of the dust layer
182	c ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
183	zlsconst=SIN(ls-2.76)
184	if (iddist.eq.1) then
185	do ig=1,ngrid
186	topdust(ig)=topdustref ! constant dust layer top
187	end do
188
189	else if (iddist.eq.2) then ! "Viking" scenario
190	do ig=1,ngrid
191	topdust(ig)=topdust0(ig)+18.*zlsconst
192	end do
193
194	else if(iddist.eq.3) then !"MGS" scenario
195	do ig=1,ngrid
196	topdust(ig)=60.+18.*zlsconst
197	& -(32+18zlsconst)sin(lati(ig))**4
198	& - 8zlsconst(sin(lati(ig)))**5
199	end do
200	endif
201
202
203	c Optical depth in each layer :
204	c ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
205	if(iddist.ge.1) then
206	expfactor=0.
207	DO l=1,nlayer
208	DO ig=1,ngrid
209	if(pplay(ig,l).gt.700.
210	$ /(988.**(topdust(ig)/70.))) then
211	zp=(700./pplay(ig,l))**(70./topdust(ig))
212	expfactor=max(exp(0.007*(1.-max(zp,1.))),1.e-3)
213	else
214	expfactor=1.e-3
215	endif
216	aerosol(ig,l,1)= tauref(ig)/700. *
217	s (pplev(ig,l)-pplev(ig,l+1))
218	& *expfactor
219	c s max( exp(.007(1.-max(zp,1.))) , 1.E-3 )
220	ENDDO
221	ENDDO
222	c changement dans le calcul de la distribution verticale
223	c dans le cas des scenarios de poussieres assimiles
224	c if (iaervar.eq.4) THEN ! TES
225	c call zerophys(ngrid*naerkind,tau)
226	c
227	c do l=1,nlayer
228	c do ig=1,ngrid
229	c tau(ig,1)=tau(ig,1)+ aerosol(ig,l,1)
230	c end do
231	c end do
232	c do l=1,nlayer
233	c do ig=1,ngrid
234	c aerosol(ig,l,1)=aerosol(ig,l,1)*tauref(ig)/tau(ig,1)
235	c $ *(pplev(ig,1)/700)
236	c end do
237	c end do
238	c endif
239	cccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
240	else if(iddist.eq.0) then
241	c old dust vertical distribution function (pollack90)
242	DO l=1,nlayer
243	DO ig=1,ngrid
244	zp=700./pplay(ig,l)
245	aerosol(ig,l,1)= tauref(ig)/700. *
246	s (pplev(ig,l)-pplev(ig,l+1))
247	s max( exp(.03(1.-max(zp,1.))) , 1.E-3 )
248	ENDDO
249	ENDDO
250	end if
251
252	c ---------------------------------------------------------------------
253	c 2) Transported radiatively active dust (if tracer=T and active=T)
254	c ----------------------------------------------------------------------
255	ELSE IF ((tracer) .and. (active)) THEN
256	c The dust opacity is computed from q
257
258	c a) "doubleq" technique (transport of mass and number mixing ratio)
259	c ~~~~~~~~~~~~~~~~~~~
260	if(doubleq) then
261
262	call zerophys(ngridnlayernaerkind,aerosol)
263
264	c Computing effective radius :
265	do l=1,nlayer
266	do ig=1, ngrid
267	r0=
268	& (r3n_qpq(ig,l,1)/max(pq(ig,l,2),0.01))*(1./3.)
269	r0=min(max(r0,1.e-10),500.e-6)
270	reff= ref_r0 * r0
271	cc If reff is small, the transported dust mean Qext
272	c is reduced from the reference dust Qext by a factor "coefsize"
273
274	coefsize=min(max(2.52e6*reff-0.043 ,0.) ,1.)
275
276	cc It is added 1.e-8 to pq to avoid low
277
278	aerosol(ig,l,1)=aerosol(ig,l,1)+ 1.E-8 +
279	& ( 0.75Qext(1)coefsize/(rho_dust*reff))
280	& * (pq(ig,l,1))*
281	& (pplev(ig,l)-pplev(ig,l+1))/g
282	end do
283	end do
284	call zerophys(ngrid,tauref)
285
286	c b) Size bin technique (each aerosol can contribute to opacity))
287	c ~~~~~~~~~~~~~~~~~~
288	else
289	c The dust opacity is computed from q
290	call zerophys(ngridnlayernaerkind,aerosol)
291	do iq=1,dustbin
292	do l=1,nlayer
293	do ig=1,ngrid
294	cc qextrhor(iq) is (3/4)Qext/(rhoreff)
295	cc It is added 1.e-8 to pq to avoid low
296	aerosol(ig,l,1)=aerosol(ig,l,1)+
297	& qextrhor(iq)* (pq(ig,l,iq) + 1.e-8)*
298	& (pplev(ig,l)-pplev(ig,l+1))/g
299	end do
300	end do
301	end do
302	call zerophys(ngrid,tauref)
303	end if ! (doubleq)
304	END IF ! (dust scenario)
305
306	c --------------------------------------------------------------------------
307	c Column integrated visible optical depth in each point (used for diagnostic)
308	c --------------------------------------------------------------------------
309	call zerophys(ngrid*naerkind,tau)
310	do l=1,nlayer
311	do ig=1,ngrid
312	tau(ig,1)=tau(ig,1)+ aerosol(ig,l,1)
313	end do
314	end do
315
316	return
317	end
318

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

Download in other formats:

Original Format