GCC Code Coverage Report


Directory: ./
File: rad/readaerosolstrato1_rrtm.f90
Date: 2022-01-11 19:19:34
Exec Total Coverage
Lines: 0 79 0.0%
Branches: 0 174 0.0%

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