| 1 | Subroutine aeropacity(ngrid,nlayer,pplay,pplev,pt, & |
|---|
| 2 | aerosol,reffrad,nueffrad,QREFvis3d,tau_col) |
|---|
| 3 | |
|---|
| 4 | use radinc_h, only : L_TAUMAX,naerkind |
|---|
| 5 | use aerosol_mod |
|---|
| 6 | |
|---|
| 7 | implicit none |
|---|
| 8 | #include "YOMCST.h" |
|---|
| 9 | |
|---|
| 10 | !================================================================== |
|---|
| 11 | ! |
|---|
| 12 | ! Purpose |
|---|
| 13 | ! ------- |
|---|
| 14 | ! Compute aerosol optical depth in each gridbox. |
|---|
| 15 | ! |
|---|
| 16 | ! Authors |
|---|
| 17 | ! ------- |
|---|
| 18 | ! F. Forget |
|---|
| 19 | ! F. Montmessin (water ice scheme) |
|---|
| 20 | ! update J.-B. Madeleine (2008) |
|---|
| 21 | ! dust removal, simplification by Robin Wordsworth (2009) |
|---|
| 22 | ! Generic n-layer aerosol - J. Vatant d'Ollone (2020) |
|---|
| 23 | ! |
|---|
| 24 | ! Input |
|---|
| 25 | ! ----- |
|---|
| 26 | ! ngrid Number of horizontal gridpoints |
|---|
| 27 | ! nlayer Number of layers |
|---|
| 28 | ! pplev Pressure (Pa) at each layer boundary |
|---|
| 29 | ! reffrad(ngrid,nlayer,naerkind) Aerosol effective radius |
|---|
| 30 | ! QREFvis3d(ngrid,nlayer,naerkind) \ 3d extinction coefficients |
|---|
| 31 | ! QREFir3d(ngrid,nlayer,naerkind) / at reference wavelengths |
|---|
| 32 | ! |
|---|
| 33 | ! Output |
|---|
| 34 | ! ------ |
|---|
| 35 | ! aerosol Aerosol optical depth in layer l, grid point ig |
|---|
| 36 | ! tau_col Total column optical depth at grid point ig |
|---|
| 37 | ! |
|---|
| 38 | !======================================================================= |
|---|
| 39 | |
|---|
| 40 | INTEGER,INTENT(IN) :: ngrid ! number of atmospheric columns |
|---|
| 41 | INTEGER,INTENT(IN) :: nlayer ! number of atmospheric layers |
|---|
| 42 | REAL,INTENT(IN) :: pplay(ngrid,nlayer) ! mid-layer pressure (Pa) |
|---|
| 43 | REAL,INTENT(IN) :: pplev(ngrid,nlayer+1) ! inter-layer pressure (Pa) |
|---|
| 44 | REAL,INTENT(IN) :: pt(ngrid,nlayer) ! mid-layer temperatre (K) |
|---|
| 45 | REAL,INTENT(OUT) :: aerosol(ngrid,nlayer,naerkind) ! aerosol optical depth |
|---|
| 46 | REAL,INTENT(IN) :: reffrad(ngrid,nlayer,naerkind) ! aerosol effective radius |
|---|
| 47 | REAL,INTENT(IN) :: nueffrad(ngrid,nlayer,naerkind) ! aerosol effective variance |
|---|
| 48 | REAL,INTENT(IN) :: QREFvis3d(ngrid,nlayer,naerkind) ! extinction coefficient in the visible |
|---|
| 49 | REAL,INTENT(OUT):: tau_col(ngrid) !column integrated visible optical depth |
|---|
| 50 | |
|---|
| 51 | real aerosol0, obs_tau_col_aurora, pm |
|---|
| 52 | real pcloud_deck, cloud_slope |
|---|
| 53 | |
|---|
| 54 | real dp_strato(ngrid) |
|---|
| 55 | real dp_tropo(ngrid) |
|---|
| 56 | real dp_layer(ngrid) |
|---|
| 57 | |
|---|
| 58 | INTEGER l,ig,iq,iaer,ia,ll |
|---|
| 59 | |
|---|
| 60 | LOGICAL,SAVE :: firstcall=.true. |
|---|
| 61 | !$OMP THREADPRIVATE(firstcall) |
|---|
| 62 | REAL CBRT |
|---|
| 63 | EXTERNAL CBRT |
|---|
| 64 | |
|---|
| 65 | ! for fixed dust profiles |
|---|
| 66 | real topdust, expfactor, zp |
|---|
| 67 | REAL taudusttmp(ngrid) ! Temporary dust opacity used before scaling |
|---|
| 68 | REAL tauh2so4tmp(ngrid) ! Temporary h2so4 opacity used before scaling |
|---|
| 69 | |
|---|
| 70 | real CLFtot |
|---|
| 71 | |
|---|
| 72 | ! for venus clouds |
|---|
| 73 | real :: p_bot2,p_bot,p_top,p_top2,h_bot2,h_bot,h_top,h_top2 |
|---|
| 74 | real :: mode_dens,h_lay,nmode1,nmode2,nmode2p,nmode3,nmodeuv |
|---|
| 75 | |
|---|
| 76 | ! First call only |
|---|
| 77 | IF (firstcall) THEN |
|---|
| 78 | |
|---|
| 79 | ! identify tracers |
|---|
| 80 | write(*,*) "Active aerosols found in aeropacity:" |
|---|
| 81 | |
|---|
| 82 | if (iaero_venus1.ne.0) then |
|---|
| 83 | print*,'iaero_venus1= ',iaero_venus1 |
|---|
| 84 | endif |
|---|
| 85 | if (iaero_venus2.ne.0) then |
|---|
| 86 | print*,'iaero_venus2= ',iaero_venus2 |
|---|
| 87 | endif |
|---|
| 88 | if (iaero_venus2p.ne.0) then |
|---|
| 89 | print*,'iaero_venus2p= ',iaero_venus2p |
|---|
| 90 | endif |
|---|
| 91 | if (iaero_venus3.ne.0) then |
|---|
| 92 | print*,'iaero_venus3= ',iaero_venus3 |
|---|
| 93 | endif |
|---|
| 94 | if (iaero_venusUV.ne.0) then |
|---|
| 95 | print*,'iaero_venusUV= ',iaero_venusUV |
|---|
| 96 | endif |
|---|
| 97 | |
|---|
| 98 | firstcall=.false. |
|---|
| 99 | ENDIF ! of IF (firstcall) |
|---|
| 100 | |
|---|
| 101 | |
|---|
| 102 | ! --------------------------------------------------------- |
|---|
| 103 | |
|---|
| 104 | !================================================================== |
|---|
| 105 | ! Venus clouds (4 modes) |
|---|
| 106 | ! S. Lebonnois (jan 2016) |
|---|
| 107 | !================================================================== |
|---|
| 108 | ! distributions from Haus et al, 2016 |
|---|
| 109 | ! mode 1 2 2p 3 |
|---|
| 110 | ! r (microns) 0.30 1.05 1.40 3.65 |
|---|
| 111 | ! sigma 1.56 1.29 1.23 1.28 |
|---|
| 112 | ! reff (microns) 0.49 1.23 1.56 4.25 |
|---|
| 113 | ! nueff 0.21 0.067 0.044 0.063 |
|---|
| 114 | ! (nueff=exp(ln^2 sigma)-1) |
|---|
| 115 | ! |
|---|
| 116 | ! p_bot <=> zb ; p_top <=> zb+zc ; h_bot <=> Hlo ; h_top <=> Hup |
|---|
| 117 | ! p<p_top: N(i+1)=N(i)*(p(i+1)/p(i))**(h_lay(i)/h_top) |
|---|
| 118 | ! h_lay=R(T(i)+T(i+1))/2/g (in m) |
|---|
| 119 | ! R=8.31/mu (mu in kg/mol) |
|---|
| 120 | ! p>p_bot: N(i-1)=N(i)*(p(i)/p(i-1))**(h_lay(i)/h_bot) |
|---|
| 121 | ! N is in m-3 |
|---|
| 122 | ! |
|---|
| 123 | ! dTau = Qext*[pi*reff**2*exp(-3*ln(1+nueff))]*N*h_lay*(-dp)/p |
|---|
| 124 | |
|---|
| 125 | ! Mode 1 |
|---|
| 126 | if (iaero_venus1 .ne.0) then |
|---|
| 127 | iaer=iaero_venus1 |
|---|
| 128 | |
|---|
| 129 | ! 1. Initialization |
|---|
| 130 | aerosol(1:ngrid,1:nlayer,iaer)=0.0 |
|---|
| 131 | ! two scales heights for below and above max density (cf Haus et al 16) |
|---|
| 132 | p_bot2 = 2.0e5 ! Pa 2.0e5 |
|---|
| 133 | h_bot2 = 5.0e3 ! m |
|---|
| 134 | p_bot = 1.2e5 ! 1.2e5 |
|---|
| 135 | h_bot = 1.0e3 |
|---|
| 136 | nmode1 = 1.935e8 ! m-3 |
|---|
| 137 | p_top = 1.0e4 |
|---|
| 138 | h_top = 3.5e3 |
|---|
| 139 | p_top2 = 4.5e2 |
|---|
| 140 | h_top2 = 2.0e3 |
|---|
| 141 | |
|---|
| 142 | ! 2. Opacity calculation |
|---|
| 143 | |
|---|
| 144 | DO ig=1,ngrid |
|---|
| 145 | ! determine the approximate middle of the mode layer |
|---|
| 146 | ll=1 |
|---|
| 147 | DO l=1,nlayer-1 ! high p to low p |
|---|
| 148 | IF (pplay(ig,l) .gt. (p_top+p_bot)/2) ll=l |
|---|
| 149 | ENDDO |
|---|
| 150 | ! from there go down and up for profile N |
|---|
| 151 | mode_dens = nmode1 ! m-3 |
|---|
| 152 | DO l=ll,1,-1 |
|---|
| 153 | h_lay=RD*(pt(ig,l+1)+pt(ig,l))/2./RG |
|---|
| 154 | IF (pplay(ig,l) .le. p_bot) THEN |
|---|
| 155 | mode_dens = nmode1 ! m-3 |
|---|
| 156 | ELSEIF (pplay(ig,l) .gt. p_bot .and. pplay(ig,l) .le. p_bot2) THEN |
|---|
| 157 | mode_dens = mode_dens*(pplay(ig,l+1)/pplay(ig,l))**(h_lay/h_bot) |
|---|
| 158 | ELSEIF (pplay(ig,l) .gt. p_bot2) THEN |
|---|
| 159 | mode_dens = mode_dens*(pplay(ig,l+1)/pplay(ig,l))**(h_lay/h_bot2) |
|---|
| 160 | ENDIF |
|---|
| 161 | mode_dens=max(1.e-30,mode_dens) |
|---|
| 162 | aerosol(ig,l,iaer) = QREFvis3d(ig,l,iaer)* & |
|---|
| 163 | RPI*(reffrad(ig,l,iaer))**2*exp(-3.*log(1+nueffrad(ig,l,iaer)))* & |
|---|
| 164 | mode_dens*h_lay*(pplev(ig,l)-pplev(ig,l+1))/pplay(ig,l) |
|---|
| 165 | ! aerosol(ig,l,iaer) = mode_dens |
|---|
| 166 | ENDDO |
|---|
| 167 | mode_dens = nmode1 ! m-3 |
|---|
| 168 | DO l=ll+1,nlayer-1 |
|---|
| 169 | h_lay=RD*(pt(ig,l-1)+pt(ig,l))/2./RG |
|---|
| 170 | IF (pplay(ig,l) .gt. p_top) THEN |
|---|
| 171 | mode_dens = nmode1 ! m-3 |
|---|
| 172 | ELSEIF (pplay(ig,l) .le. p_top .and. pplay(ig,l) .gt. p_top2) THEN |
|---|
| 173 | mode_dens = mode_dens*(pplay(ig,l)/pplay(ig,l-1))**(h_lay/h_top) |
|---|
| 174 | ELSEIF (pplay(ig,l) .le. p_top2) THEN |
|---|
| 175 | mode_dens = mode_dens*(pplay(ig,l)/pplay(ig,l-1))**(h_lay/h_top2) |
|---|
| 176 | ENDIF |
|---|
| 177 | mode_dens=max(1.e-30,mode_dens) |
|---|
| 178 | aerosol(ig,l,iaer) = QREFvis3d(ig,l,iaer)* & |
|---|
| 179 | RPI*(reffrad(ig,l,iaer))**2*exp(-3.*log(1+nueffrad(ig,l,iaer)))* & |
|---|
| 180 | mode_dens*h_lay*(pplev(ig,l)-pplev(ig,l+1))/pplay(ig,l) |
|---|
| 181 | ! aerosol(ig,l,iaer) = mode_dens |
|---|
| 182 | ENDDO |
|---|
| 183 | ENDDO |
|---|
| 184 | |
|---|
| 185 | end if ! mode 1 |
|---|
| 186 | |
|---|
| 187 | ! Mode 2 |
|---|
| 188 | if (iaero_venus2 .ne.0) then |
|---|
| 189 | iaer=iaero_venus2 |
|---|
| 190 | |
|---|
| 191 | ! 1. Initialization |
|---|
| 192 | aerosol(1:ngrid,1:nlayer,iaer)=0.0 |
|---|
| 193 | p_bot = 1.0e4 |
|---|
| 194 | h_bot = 3.0e3 |
|---|
| 195 | nmode2 = 1.00e8 ! m-3 |
|---|
| 196 | p_top = 8.0e3 |
|---|
| 197 | h_top = 3.5e3 |
|---|
| 198 | p_top2 = 4.5e2 |
|---|
| 199 | h_top2 = 2.0e3 |
|---|
| 200 | |
|---|
| 201 | ! 2. Opacity calculation |
|---|
| 202 | |
|---|
| 203 | DO ig=1,ngrid |
|---|
| 204 | ! determine the approximate middle of the mode layer |
|---|
| 205 | ll=1 |
|---|
| 206 | DO l=1,nlayer-1 ! high p to low p |
|---|
| 207 | IF (pplay(ig,l) .gt. (p_top+p_bot)/2) ll=l |
|---|
| 208 | ENDDO |
|---|
| 209 | ! from there go down and up for profile N |
|---|
| 210 | mode_dens = nmode2 ! m-3 |
|---|
| 211 | DO l=ll,1,-1 |
|---|
| 212 | h_lay=RD*(pt(ig,l+1)+pt(ig,l))/2./RG |
|---|
| 213 | IF (pplay(ig,l) .le. p_bot) THEN |
|---|
| 214 | mode_dens = nmode2 ! m-3 |
|---|
| 215 | ELSEIF (pplay(ig,l) .gt. p_bot) THEN |
|---|
| 216 | mode_dens = mode_dens*(pplay(ig,l+1)/pplay(ig,l))**(h_lay/h_bot) |
|---|
| 217 | ENDIF |
|---|
| 218 | mode_dens=max(1.e-30,mode_dens) |
|---|
| 219 | aerosol(ig,l,iaer) = QREFvis3d(ig,l,iaer)* & |
|---|
| 220 | RPI*(reffrad(ig,l,iaer))**2*exp(-3.*log(1+nueffrad(ig,l,iaer)))* & |
|---|
| 221 | mode_dens*h_lay*(pplev(ig,l)-pplev(ig,l+1))/pplay(ig,l) |
|---|
| 222 | ! aerosol(ig,l,iaer) = mode_dens |
|---|
| 223 | ENDDO |
|---|
| 224 | mode_dens = nmode2 ! m-3 |
|---|
| 225 | DO l=ll+1,nlayer-1 |
|---|
| 226 | h_lay=RD*(pt(ig,l-1)+pt(ig,l))/2./RG |
|---|
| 227 | IF (pplay(ig,l) .gt. p_top) THEN |
|---|
| 228 | mode_dens = nmode2 ! m-3 |
|---|
| 229 | ELSEIF (pplay(ig,l) .le. p_top .and. pplay(ig,l) .gt. p_top2) THEN |
|---|
| 230 | mode_dens = mode_dens*(pplay(ig,l)/pplay(ig,l-1))**(h_lay/h_top) |
|---|
| 231 | ELSEIF (pplay(ig,l) .le. p_top2) THEN |
|---|
| 232 | mode_dens = mode_dens*(pplay(ig,l)/pplay(ig,l-1))**(h_lay/h_top2) |
|---|
| 233 | ENDIF |
|---|
| 234 | mode_dens=max(1.e-30,mode_dens) |
|---|
| 235 | aerosol(ig,l,iaer) = QREFvis3d(ig,l,iaer)* & |
|---|
| 236 | RPI*(reffrad(ig,l,iaer))**2*exp(-3.*log(1+nueffrad(ig,l,iaer)))* & |
|---|
| 237 | mode_dens*h_lay*(pplev(ig,l)-pplev(ig,l+1))/pplay(ig,l) |
|---|
| 238 | ! aerosol(ig,l,iaer) = mode_dens |
|---|
| 239 | ENDDO |
|---|
| 240 | ENDDO |
|---|
| 241 | |
|---|
| 242 | end if ! mode 2 |
|---|
| 243 | |
|---|
| 244 | ! Mode 2p |
|---|
| 245 | if (iaero_venus2p .ne.0) then |
|---|
| 246 | iaer=iaero_venus2p |
|---|
| 247 | |
|---|
| 248 | ! 1. Initialization |
|---|
| 249 | aerosol(1:ngrid,1:nlayer,iaer)=0.0 |
|---|
| 250 | p_bot = 1.2e5 ! 1.2e5 |
|---|
| 251 | h_bot = 0.1e3 |
|---|
| 252 | nmode2p = 5.0e7 ! m-3 |
|---|
| 253 | p_top = 2.3e4 |
|---|
| 254 | h_top = 1.0e3 |
|---|
| 255 | |
|---|
| 256 | ! 2. Opacity calculation |
|---|
| 257 | |
|---|
| 258 | DO ig=1,ngrid |
|---|
| 259 | ! determine the approximate middle of the mode layer |
|---|
| 260 | ll=1 |
|---|
| 261 | DO l=1,nlayer-1 ! high p to low p |
|---|
| 262 | IF (pplay(ig,l) .gt. (p_top+p_bot)/2) ll=l |
|---|
| 263 | ENDDO |
|---|
| 264 | ! from there go down and up for profile N |
|---|
| 265 | mode_dens = nmode2p ! m-3 |
|---|
| 266 | DO l=ll,1,-1 |
|---|
| 267 | h_lay=RD*(pt(ig,l+1)+pt(ig,l))/2./RG |
|---|
| 268 | IF (pplay(ig,l) .le. p_bot) THEN |
|---|
| 269 | mode_dens = nmode2p ! m-3 |
|---|
| 270 | ELSEIF (pplay(ig,l) .gt. p_bot) THEN |
|---|
| 271 | mode_dens = mode_dens*(pplay(ig,l+1)/pplay(ig,l))**(h_lay/h_bot) |
|---|
| 272 | ENDIF |
|---|
| 273 | mode_dens=max(1.e-30,mode_dens) |
|---|
| 274 | aerosol(ig,l,iaer) = QREFvis3d(ig,l,iaer)* & |
|---|
| 275 | RPI*(reffrad(ig,l,iaer))**2*exp(-3.*log(1+nueffrad(ig,l,iaer)))* & |
|---|
| 276 | mode_dens*h_lay*(pplev(ig,l)-pplev(ig,l+1))/pplay(ig,l) |
|---|
| 277 | ! aerosol(ig,l,iaer) = mode_dens |
|---|
| 278 | ENDDO |
|---|
| 279 | mode_dens = nmode2p ! m-3 |
|---|
| 280 | DO l=ll+1,nlayer-1 |
|---|
| 281 | h_lay=RD*(pt(ig,l-1)+pt(ig,l))/2./RG |
|---|
| 282 | IF (pplay(ig,l) .gt. p_top) THEN |
|---|
| 283 | mode_dens = nmode2p ! m-3 |
|---|
| 284 | ELSEIF (pplay(ig,l) .le. p_top) THEN |
|---|
| 285 | mode_dens = mode_dens*(pplay(ig,l)/pplay(ig,l-1))**(h_lay/h_top) |
|---|
| 286 | ENDIF |
|---|
| 287 | mode_dens=max(1.e-30,mode_dens) |
|---|
| 288 | aerosol(ig,l,iaer) = QREFvis3d(ig,l,iaer)* & |
|---|
| 289 | RPI*(reffrad(ig,l,iaer))**2*exp(-3.*log(1+nueffrad(ig,l,iaer)))* & |
|---|
| 290 | mode_dens*h_lay*(pplev(ig,l)-pplev(ig,l+1))/pplay(ig,l) |
|---|
| 291 | ! aerosol(ig,l,iaer) = mode_dens |
|---|
| 292 | ENDDO |
|---|
| 293 | ENDDO |
|---|
| 294 | |
|---|
| 295 | end if ! mode 2p |
|---|
| 296 | |
|---|
| 297 | ! Mode 3 |
|---|
| 298 | if (iaero_venus3 .ne.0) then |
|---|
| 299 | iaer=iaero_venus3 |
|---|
| 300 | |
|---|
| 301 | ! 1. Initialization |
|---|
| 302 | aerosol(1:ngrid,1:nlayer,iaer)=0.0 |
|---|
| 303 | p_bot = 1.2e5 ! 1.2e5 |
|---|
| 304 | h_bot = 0.5e3 |
|---|
| 305 | nmode3 = 1.4e7 ! m-3 |
|---|
| 306 | p_top = 3.9e4 |
|---|
| 307 | h_top = 1.0e3 |
|---|
| 308 | |
|---|
| 309 | ! 2. Opacity calculation |
|---|
| 310 | |
|---|
| 311 | DO ig=1,ngrid |
|---|
| 312 | ! determine the approximate middle of the mode layer |
|---|
| 313 | ll=1 |
|---|
| 314 | DO l=1,nlayer-1 ! high p to low p |
|---|
| 315 | IF (pplay(ig,l) .gt. (p_top+p_bot)/2) ll=l |
|---|
| 316 | ENDDO |
|---|
| 317 | ! from there go down and up for profile N |
|---|
| 318 | mode_dens = nmode3 ! m-3 |
|---|
| 319 | DO l=ll,1,-1 |
|---|
| 320 | h_lay=RD*(pt(ig,l+1)+pt(ig,l))/2./RG |
|---|
| 321 | IF (pplay(ig,l) .le. p_bot) THEN |
|---|
| 322 | mode_dens = nmode3 ! m-3 |
|---|
| 323 | ELSEIF (pplay(ig,l) .gt. p_bot) THEN |
|---|
| 324 | mode_dens = mode_dens*(pplay(ig,l+1)/pplay(ig,l))**(h_lay/h_bot) |
|---|
| 325 | ENDIF |
|---|
| 326 | mode_dens=max(1.e-30,mode_dens) |
|---|
| 327 | aerosol(ig,l,iaer) = QREFvis3d(ig,l,iaer)* & |
|---|
| 328 | RPI*(reffrad(ig,l,iaer))**2*exp(-3.*log(1+nueffrad(ig,l,iaer)))* & |
|---|
| 329 | mode_dens*h_lay*(pplev(ig,l)-pplev(ig,l+1))/pplay(ig,l) |
|---|
| 330 | ! aerosol(ig,l,iaer) = mode_dens |
|---|
| 331 | ENDDO |
|---|
| 332 | mode_dens = nmode3 ! m-3 |
|---|
| 333 | DO l=ll+1,nlayer-1 |
|---|
| 334 | h_lay=RD*(pt(ig,l-1)+pt(ig,l))/2./RG |
|---|
| 335 | IF (pplay(ig,l) .gt. p_top) THEN |
|---|
| 336 | mode_dens = nmode3 ! m-3 |
|---|
| 337 | ELSEIF (pplay(ig,l) .le. p_top) THEN |
|---|
| 338 | mode_dens = mode_dens*(pplay(ig,l)/pplay(ig,l-1))**(h_lay/h_top) |
|---|
| 339 | ENDIF |
|---|
| 340 | mode_dens=max(1.e-30,mode_dens) |
|---|
| 341 | aerosol(ig,l,iaer) = QREFvis3d(ig,l,iaer)* & |
|---|
| 342 | RPI*(reffrad(ig,l,iaer))**2*exp(-3.*log(1+nueffrad(ig,l,iaer)))* & |
|---|
| 343 | mode_dens*h_lay*(pplev(ig,l)-pplev(ig,l+1))/pplay(ig,l) |
|---|
| 344 | ! aerosol(ig,l,iaer) = mode_dens |
|---|
| 345 | ENDDO |
|---|
| 346 | ENDDO |
|---|
| 347 | |
|---|
| 348 | end if ! mode 3 |
|---|
| 349 | |
|---|
| 350 | ! UV absorber |
|---|
| 351 | if (iaero_venusUV .ne.0) then |
|---|
| 352 | iaer=iaero_venusUV |
|---|
| 353 | |
|---|
| 354 | ! 1. Initialization |
|---|
| 355 | aerosol(1:ngrid,1:nlayer,iaer)=0.0 |
|---|
| 356 | p_bot = 3.3e4 ! 58 km |
|---|
| 357 | h_bot = 1.0e3 |
|---|
| 358 | nmodeuv = 1.00e7 !m-3 |
|---|
| 359 | p_top = 3.7e3 ! 70 km |
|---|
| 360 | h_top = 1.0e3 |
|---|
| 361 | |
|---|
| 362 | ! 2. Opacity calculation |
|---|
| 363 | |
|---|
| 364 | DO ig=1,ngrid |
|---|
| 365 | ! determine the approximate middle of the mode layer |
|---|
| 366 | ll=1 |
|---|
| 367 | DO l=1,nlayer-1 ! high p to low p |
|---|
| 368 | IF (pplay(ig,l) .gt. (p_top+p_bot)/2) ll=l |
|---|
| 369 | ENDDO |
|---|
| 370 | ! from there go down and up for profile N |
|---|
| 371 | mode_dens = nmodeuv ! m-3 |
|---|
| 372 | DO l=ll,1,-1 |
|---|
| 373 | h_lay=RD*(pt(ig,l+1)+pt(ig,l))/2./RG |
|---|
| 374 | IF (pplay(ig,l) .le. p_bot) THEN |
|---|
| 375 | mode_dens = nmodeuv ! m-3 |
|---|
| 376 | ELSEIF (pplay(ig,l) .gt. p_bot) THEN |
|---|
| 377 | mode_dens = mode_dens*(pplay(ig,l+1)/pplay(ig,l))**(h_lay/h_bot) |
|---|
| 378 | ENDIF |
|---|
| 379 | mode_dens=max(1.e-30,mode_dens) |
|---|
| 380 | ! normalized to 0.35 microns (peak of absorption) |
|---|
| 381 | aerosol(ig,l,iaer) = QREFvis3d(ig,l,iaer)*mode_dens |
|---|
| 382 | ENDDO |
|---|
| 383 | mode_dens = nmodeuv ! m-3 |
|---|
| 384 | DO l=ll+1,nlayer-1 |
|---|
| 385 | h_lay=RD*(pt(ig,l-1)+pt(ig,l))/2./RG |
|---|
| 386 | IF (pplay(ig,l) .gt. p_top) THEN |
|---|
| 387 | mode_dens = nmodeuv ! m-3 |
|---|
| 388 | ELSEIF (pplay(ig,l) .le. p_top) THEN |
|---|
| 389 | mode_dens = mode_dens*(pplay(ig,l)/pplay(ig,l-1))**(h_lay/h_top) |
|---|
| 390 | ENDIF |
|---|
| 391 | mode_dens=max(1.e-30,mode_dens) |
|---|
| 392 | ! normalized to 0.35 microns (peak of absorption) |
|---|
| 393 | aerosol(ig,l,iaer) = QREFvis3d(ig,l,iaer)*mode_dens |
|---|
| 394 | ENDDO |
|---|
| 395 | ENDDO |
|---|
| 396 | |
|---|
| 397 | ! 3. Re-normalize to Haus et al 2015 total column optical depth |
|---|
| 398 | tau_col(:)=0.0 |
|---|
| 399 | DO l=1,nlayer |
|---|
| 400 | DO ig=1,ngrid |
|---|
| 401 | tau_col(ig) = tau_col(ig) & |
|---|
| 402 | + aerosol(ig,l,iaer) |
|---|
| 403 | ENDDO |
|---|
| 404 | ENDDO |
|---|
| 405 | DO ig=1,ngrid |
|---|
| 406 | DO l=1,nlayer-1 |
|---|
| 407 | aerosol(ig,l,iaer)=aerosol(ig,l,iaer)*0.205/tau_col(ig) |
|---|
| 408 | ENDDO |
|---|
| 409 | ENDDO |
|---|
| 410 | |
|---|
| 411 | end if ! UV absorber |
|---|
| 412 | |
|---|
| 413 | !================================================================== |
|---|
| 414 | ! if (is_master) then |
|---|
| 415 | ! ig=10 |
|---|
| 416 | ! do l=1,nlayer |
|---|
| 417 | ! write(82,*) l,pplay(ig,l),aerosol(ig,l,1) |
|---|
| 418 | ! write(83,*) l,pplay(ig,l),aerosol(ig,l,2) |
|---|
| 419 | ! write(84,*) l,pplay(ig,l),aerosol(ig,l,3) |
|---|
| 420 | ! write(85,*) l,pplay(ig,l),aerosol(ig,l,4) |
|---|
| 421 | ! write(86,*) l,pplay(ig,l),aerosol(ig,l,5) |
|---|
| 422 | ! enddo |
|---|
| 423 | ! endif |
|---|
| 424 | ! stop |
|---|
| 425 | !================================================================== |
|---|
| 426 | |
|---|
| 427 | ! -------------------------------------------------------------------------- |
|---|
| 428 | ! Column integrated visible optical depth in each point (used for diagnostic) |
|---|
| 429 | |
|---|
| 430 | tau_col(:)=0.0 |
|---|
| 431 | do iaer = 1, naerkind |
|---|
| 432 | do l=1,nlayer |
|---|
| 433 | do ig=1,ngrid |
|---|
| 434 | tau_col(ig) = tau_col(ig) + aerosol(ig,l,iaer) |
|---|
| 435 | end do |
|---|
| 436 | end do |
|---|
| 437 | end do |
|---|
| 438 | |
|---|
| 439 | do ig=1,ngrid |
|---|
| 440 | do l=1,nlayer |
|---|
| 441 | do iaer = 1, naerkind |
|---|
| 442 | if(aerosol(ig,l,iaer).gt.1.e3)then |
|---|
| 443 | print*,'WARNING: aerosol=',aerosol(ig,l,iaer) |
|---|
| 444 | print*,'at ig=',ig,', l=',l,', iaer=',iaer |
|---|
| 445 | print*,'QREFvis3d=',QREFvis3d(ig,l,iaer) |
|---|
| 446 | print*,'reffrad=',reffrad(ig,l,iaer) |
|---|
| 447 | endif |
|---|
| 448 | end do |
|---|
| 449 | end do |
|---|
| 450 | end do |
|---|
| 451 | |
|---|
| 452 | do ig=1,ngrid |
|---|
| 453 | if(tau_col(ig).gt.1.e3)then |
|---|
| 454 | print*,'WARNING: tau_col=',tau_col(ig) |
|---|
| 455 | print*,'at ig=',ig |
|---|
| 456 | print*,'aerosol=',aerosol(ig,:,:) |
|---|
| 457 | print*,'QREFvis3d=',QREFvis3d(ig,:,:) |
|---|
| 458 | print*,'reffrad=',reffrad(ig,:,:) |
|---|
| 459 | endif |
|---|
| 460 | end do |
|---|
| 461 | return |
|---|
| 462 | end subroutine aeropacity |
|---|
| 463 | |
|---|
