Context Navigation

readaerosolstrato1_rrtm.F90 @ 3437

Last change on this file since 3437 was 3437, checked in by fhourdin, 5 years ago
Compiler sans XIOS Fredho
File size: 10.2 KB

Line
1	!
2	! $Id: readaerosolstrato1_rrtm.F90 2526 2016-05-26 22:13:40Z oboucher $
3	!
4
5	SUBROUTINE readaerosolstrato1_rrtm(debut)
6
7	USE netcdf95, ONLY: nf95_close, nf95_gw_var, nf95_inq_dimid, &
8	nf95_inq_varid, nf95_open
9	USE netcdf, ONLY: nf90_get_var, nf90_noerr, nf90_nowrite
10
11	USE phys_cal_mod, ONLY : mth_cur
12	USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat, klon_glo, grid2dTo1d_glo, grid_type, unstructured
13	USE mod_phys_lmdz_para
14	USE phys_state_var_mod
15	USE phys_local_var_mod
16	USE aero_mod
17	USE dimphy
18	USE YOERAD, ONLY : NLW
19	USE YOMCST
20	#ifdef CPP_XIOS
21	USE xios
22	#endif
23
24	IMPLICIT NONE
25
26	! Variable input
27	LOGICAL debut
28
29	! Variables locales
30	INTEGER n_lat ! number of latitudes in the input data
31	INTEGER n_lon ! number of longitudes in the input data
32	INTEGER n_lev ! number of levels in the input data
33	INTEGER n_month ! number of months in the input data
34	REAL, POINTER:: latitude(:)
35	REAL, POINTER:: longitude(:)
36	REAL, POINTER:: time(:)
37	REAL, POINTER:: lev(:)
38	INTEGER k, band, wave, i
39	INTEGER, SAVE :: mth_pre=1
40	!$OMP THREADPRIVATE(mth_pre)
41
42	REAL, ALLOCATABLE, DIMENSION(:,:), SAVE :: tau_aer_strat
43	!$OMP THREADPRIVATE(tau_aer_strat)
44
45	! Champs reconstitues
46	REAL, ALLOCATABLE:: tauaerstrat(:, :, :, :)
47	REAL, ALLOCATABLE:: tauaerstrat_mois(:, :, :)
48	REAL, ALLOCATABLE:: tauaerstrat_mois_glo(:, :)
49	REAL, ALLOCATABLE:: tauaerstrat_mpi(:, :)
50
51	! For NetCDF:
52	INTEGER ncid_in ! IDs for input files
53	INTEGER varid, ncerr
54
55	! Stratospheric aerosols optical properties
56	! alpha_sw_strat over the 6 bands is normalised by the 550 nm extinction coefficient
57	REAL, DIMENSION(nbands_sw_rrtm) :: alpha_sw_strat, piz_sw_strat, cg_sw_strat
58	DATA alpha_sw_strat/0.8545564, 0.8451642, 0.9821724, 0.8145110, 0.3073565, 7.7966176E-02/
59	DATA cg_sw_strat /0.6997170, 0.6810035, 0.7403592, 0.7562674, 0.6676504, 0.3478689/
60	DATA piz_sw_strat /0.9999998, 0.9999998, 1.000000000, 0.9999958, 0.9977155, 0.4510679/
61	!
62	!--diagnostics AOD in the SW
63	! alpha_sw_strat_wave is not normalised by the 550 nm extinction coefficient
64	REAL, DIMENSION(nwave_sw) :: alpha_sw_strat_wave
65	DATA alpha_sw_strat_wave/3.708007,4.125824,4.136584,3.887478,3.507738/
66	!
67	!--diagnostics AOD in the LW at 10 um (not normalised by the 550 nm ext coefficient
68	REAL :: alpha_lw_strat_wave(nwave_lw)
69	DATA alpha_lw_strat_wave/0.2746812/
70	!
71	REAL, DIMENSION(nbands_lw_rrtm) :: alpha_lw_abs_rrtm
72	DATA alpha_lw_abs_rrtm/ 8.8340312E-02, 6.9856711E-02, 6.2652975E-02, 5.7188231E-02, &
73	6.3157059E-02, 5.5072524E-02, 5.0571125E-02, 0.1349073, &
74	0.1381676, 9.6506312E-02, 5.1312990E-02, 2.4256418E-02, &
75	2.7191756E-02, 3.3862915E-02, 1.6132960E-02, 1.4275438E-02/ ! calculated with Mie_SW_LW_RRTM_V2.4 (bimodal, corrected)
76	! for r_0=/0.13E-6, 0.41E-6/ m, sigma_g=/1.26, 1.30/
77	! order: increasing wavelength!
78	!--------------------------------------------------------
79
80	IF (.not.ALLOCATED(tau_aer_strat)) ALLOCATE(tau_aer_strat(klon,klev))
81
82	!--we only read monthly strat aerosol data
83	IF (debut.OR.mth_cur.NE.mth_pre) THEN
84
85	!--only root reads the data
86	IF (is_mpi_root.AND.is_omp_root) THEN
87
88	IF (nbands_sw_rrtm.NE.6) THEN
89	print *,'nbands_sw_rrtm doit etre egal a 6 dans readaerosolstrat_rrtm'
90	STOP
91	ENDIF
92
93	CALL nf95_open("taustrat.nc", nf90_nowrite, ncid_in)
94
95	CALL nf95_inq_varid(ncid_in, "LEV", varid)
96	CALL nf95_gw_var(ncid_in, varid, lev)
97	n_lev = size(lev)
98	IF (n_lev.NE.klev) THEN
99	print *,'Le nombre de niveaux n est pas egal a klev'
100	STOP
101	ENDIF
102
103	CALL nf95_inq_varid(ncid_in, "LAT", varid)
104	CALL nf95_gw_var(ncid_in, varid, latitude)
105	n_lat = size(latitude)
106	print *, 'LAT aerosol strato=', n_lat, latitude
107
108	IF (grid_type/=unstructured) THEN
109	IF (n_lat.NE.nbp_lat) THEN
110	print *,'Le nombre de lat n est pas egal a nbp_lat'
111	STOP
112	ENDIF
113	ENDIF
114
115	CALL nf95_inq_varid(ncid_in, "LON", varid)
116	CALL nf95_gw_var(ncid_in, varid, longitude)
117	n_lon = size(longitude)
118	print *, 'LON aerosol strato=', n_lon, longitude
119
120	IF (grid_type/=unstructured) THEN
121	IF (n_lon.NE.nbp_lon) THEN
122	print *,'Le nombre de lon n est pas egal a nbp_lon'
123	STOP
124	ENDIF
125	ENDIF
126
127
128	CALL nf95_inq_varid(ncid_in, "TIME", varid)
129	CALL nf95_gw_var(ncid_in, varid, time)
130	n_month = size(time)
131	print *, 'TIME aerosol strato=', n_month, time
132	IF (n_month.NE.12) THEN
133	print *,'Le nombre de month n est pas egal a 12'
134	STOP
135	ENDIF
136
137	ALLOCATE(tauaerstrat(n_lon, n_lat, n_lev, n_month))
138	ALLOCATE(tauaerstrat_mois(n_lon, n_lat, n_lev))
139	ALLOCATE(tauaerstrat_mois_glo(klon_glo, n_lev))
140
141	!--reading stratospheric AOD at 550 nm
142	CALL nf95_inq_varid(ncid_in, "TAUSTRAT", varid)
143	ncerr = nf90_get_var(ncid_in, varid, tauaerstrat)
144	print *,'code erreur readaerosolstrato=', ncerr, varid
145
146	CALL nf95_close(ncid_in)
147
148	!---select the correct month
149	IF (mth_cur.LT.1.OR.mth_cur.GT.12) THEN
150	print *,'probleme avec le mois dans readaerosolstrat =', mth_cur
151	ENDIF
152	tauaerstrat_mois(:,:,:) = tauaerstrat(:,:,:,mth_cur)
153
154	!---reduce to a klon_glo grid
155	CALL grid2dTo1d_glo(tauaerstrat_mois,tauaerstrat_mois_glo)
156
157	ELSE
158	ALLOCATE(tauaerstrat_mois(0,0,0))
159	ENDIF !--is_mpi_root and is_omp_root
160
161	!$OMP BARRIER
162
163	!--keep memory of previous month
164	mth_pre=mth_cur
165
166	!--scatter on all proc
167
168	IF (grid_type==unstructured) THEN
169	IF (is_omp_master) THEN
170	ALLOCATE(tauaerstrat_mpi(klon_mpi,klev))
171	CALL xios_send_field("taustrat_in",tauaerstrat_mois)
172	CALL xios_recv_field("taustrat_out",tauaerstrat_mpi)
173	ELSE
174	ALLOCATE(tauaerstrat_mpi(0,0))
175	ENDIF
176	CALL scatter_omp(tauaerstrat_mpi,tau_aer_strat)
177	ELSE
178	CALL scatter(tauaerstrat_mois_glo,tau_aer_strat)
179	ENDIF
180
181	IF (is_mpi_root.AND.is_omp_root) THEN
182	!
183	DEALLOCATE(tauaerstrat)
184	DEALLOCATE(tauaerstrat_mois)
185	DEALLOCATE(tauaerstrat_mois_glo)
186	!
187	ENDIF !--is_mpi_root and is_omp_root
188
189	!$OMP BARRIER
190
191	ENDIF !--debut ou nouveau mois
192
193	!--total vertical aod at the 5 SW wavelengths
194	DO wave=1, nwave_sw
195	DO k=1, klev
196	tausum_aero(:,wave,id_STRAT_phy)=tausum_aero(:,wave,id_STRAT_phy)+ &
197	tau_aer_strat(:,k)*alpha_sw_strat_wave(wave)/alpha_sw_strat_wave(2)
198	ENDDO
199	ENDDO
200
201	!--weighted average for cg, piz and tau, adding strat aerosols on top of tropospheric ones
202	DO band=1, nbands_sw_rrtm
203	!--anthropogenic aerosols bands 1 to nbands_sw_rrtm
204	cg_aero_sw_rrtm(:,:,2,band) = ( cg_aero_sw_rrtm(:,:,2,band)piz_aero_sw_rrtm(:,:,2,band)tau_aero_sw_rrtm(:,:,2,band) + &
205	cg_sw_strat(band)piz_sw_strat(band)alpha_sw_strat(band)*tau_aer_strat(:,:) ) / &
206	MAX( piz_aero_sw_rrtm(:,:,2,band)*tau_aero_sw_rrtm(:,:,2,band) + &
207	piz_sw_strat(band)alpha_sw_strat(band)tau_aer_strat(:,:), 1.e-15 )
208	piz_aero_sw_rrtm(:,:,2,band) = ( piz_aero_sw_rrtm(:,:,2,band)*tau_aero_sw_rrtm(:,:,2,band) + &
209	piz_sw_strat(band)alpha_sw_strat(band)tau_aer_strat(:,:) ) / &
210	MAX( tau_aero_sw_rrtm(:,:,2,band) + alpha_sw_strat(band)*tau_aer_strat(:,:), 1.e-15 )
211	tau_aero_sw_rrtm(:,:,2,band) = tau_aero_sw_rrtm(:,:,2,band) + alpha_sw_strat(band)*tau_aer_strat(:,:)
212	!--natural aerosols bands 1 to nbands_sw_rrtm
213	cg_aero_sw_rrtm(:,:,1,band) = ( cg_aero_sw_rrtm(:,:,1,band)piz_aero_sw_rrtm(:,:,1,band)tau_aero_sw_rrtm(:,:,1,band) + &
214	cg_sw_strat(band)piz_sw_strat(band)alpha_sw_strat(band)*tau_aer_strat(:,:) ) / &
215	MAX( piz_aero_sw_rrtm(:,:,1,band)*tau_aero_sw_rrtm(:,:,1,band) + &
216	piz_sw_strat(band)alpha_sw_strat(band)tau_aer_strat(:,:), 1.e-15 )
217	piz_aero_sw_rrtm(:,:,1,band) = ( piz_aero_sw_rrtm(:,:,1,band)*tau_aero_sw_rrtm(:,:,1,band) + &
218	piz_sw_strat(band)alpha_sw_strat(band)tau_aer_strat(:,:) ) / &
219	MAX( tau_aero_sw_rrtm(:,:,1,band) + alpha_sw_strat(band)*tau_aer_strat(:,:),1.e-15 )
220	tau_aero_sw_rrtm(:,:,1,band) = tau_aero_sw_rrtm(:,:,1,band) + alpha_sw_strat(band)*tau_aer_strat(:,:)
221	!--no stratospheric aerosol in index 1 for these tests
222	! cg_aero_sw_rrtm(:,:,1,band) = cg_aero_sw_rrtm(:,:,1,band)
223	! piz_aero_sw_rrtm(:,:,1,band) = piz_aero_sw_rrtm(:,:,1,band)
224	! tau_aero_sw_rrtm(:,:,1,band) = tau_aero_sw_rrtm(:,:,1,band)
225	ENDDO
226
227	!--stratospheric AOD in LW
228	IF (nbands_lw_rrtm .NE. NLW) then
229	print*, 'different values for NLW (=',NLW,') and nbands_lw_rrtm (=', nbands_lw_rrtm, ')'
230	STOP
231	ENDIF
232
233	!--total vertical aod at the 1 LW wavelength
234	DO wave=1, nwave_lw
235	DO k=1, klev
236	tausum_aero(:,nwave_sw+wave,id_STRAT_phy)=tausum_aero(:,nwave_sw+wave,id_STRAT_phy)+ &
237	tau_aer_strat(:,k)*alpha_lw_strat_wave(wave)/alpha_sw_strat_wave(2)
238	ENDDO
239	ENDDO
240
241	DO band=1, nbands_lw_rrtm
242	tau_aero_lw_rrtm(:,:,2,band) = tau_aero_lw_rrtm(:,:,2,band) + alpha_lw_abs_rrtm(band)*tau_aer_strat(:,:)
243	tau_aero_lw_rrtm(:,:,1,band) = tau_aero_lw_rrtm(:,:,1,band) + alpha_lw_abs_rrtm(band)*tau_aer_strat(:,:)
244	!--no stratospheric aerosols in index 1 for these tests
245	! tau_aero_lw_rrtm(:,:,1,band) = tau_aero_lw_rrtm(:,:,1,band)
246	ENDDO
247
248	!--default SSA value if there is no aerosol
249	!--to avoid 0 values that seems to cause some problem to RRTM
250	WHERE (tau_aero_sw_rrtm.LT.1.e-14)
251	piz_aero_sw_rrtm = 1.0
252	ENDWHERE
253
254	!--in principle this should not be necessary
255	!--as these variables have min values already but just in case
256	!--put 1e-15 min value to both SW and LW AOD
257	tau_aero_sw_rrtm = MAX(tau_aero_sw_rrtm,1.e-15)
258	tau_aero_lw_rrtm = MAX(tau_aero_lw_rrtm,1.e-15)
259
260	END SUBROUTINE readaerosolstrato1_rrtm

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

Download in other formats:

Original Format