| 1 | MODULE aeroptproperties_mod |
|---|
| 2 | |
|---|
| 3 | IMPLICIT NONE |
|---|
| 4 | |
|---|
| 5 | CONTAINS |
|---|
| 6 | |
|---|
| 7 | SUBROUTINE aeroptproperties(ngrid,nlayer,reffrad,nueffrad, |
|---|
| 8 | & QVISsQREF3d,omegaVIS3d,gVIS3d, |
|---|
| 9 | & QIRsQREF3d,omegaIR3d,gIR3d, |
|---|
| 10 | & QREFvis3d,QREFir3d, |
|---|
| 11 | & omegaREFvis3d,omegaREFir3d) |
|---|
| 12 | use dimradmars_mod, only: nir, nsun, naerkind, |
|---|
| 13 | & radiustab, nsize, QVISsQREF, omegavis, gvis, |
|---|
| 14 | & QIRsQREF, omegaIR, gIR, QREFvis, QREFir, |
|---|
| 15 | & omegaREFvis, omegaREFir |
|---|
| 16 | use write_output_mod, only: write_output |
|---|
| 17 | |
|---|
| 18 | IMPLICIT NONE |
|---|
| 19 | c ============================================================= |
|---|
| 20 | c Aerosol Optical Properties |
|---|
| 21 | c |
|---|
| 22 | c Description: |
|---|
| 23 | c Compute the scattering parameters in each grid |
|---|
| 24 | c box, depending on aerosol grain sizes. Log-normal size |
|---|
| 25 | c distribution and Gauss-Legendre integration are used. |
|---|
| 26 | |
|---|
| 27 | c Parameters: |
|---|
| 28 | c Don't forget to set the value of varyingnueff below; If |
|---|
| 29 | c the effective variance of the distribution for the given |
|---|
| 30 | c aerosol is considered homogeneous in the atmosphere, please |
|---|
| 31 | c set varyingnueff(iaer) to .false. Resulting computational |
|---|
| 32 | c time will be much better. |
|---|
| 33 | |
|---|
| 34 | c Authors: J.-B. Madeleine, F. Forget, F. Montmessin |
|---|
| 35 | c ============================================================= |
|---|
| 36 | |
|---|
| 37 | c Local variables |
|---|
| 38 | c --------------- |
|---|
| 39 | |
|---|
| 40 | c ============================================================= |
|---|
| 41 | LOGICAL,SAVE,ALLOCATABLE :: varyingnueff(:) |
|---|
| 42 | c ============================================================= |
|---|
| 43 | |
|---|
| 44 | c Min. and max radius of the interpolation grid (in METERS) |
|---|
| 45 | REAL, SAVE, ALLOCATABLE :: refftabmin(:,:) |
|---|
| 46 | REAL, SAVE, ALLOCATABLE :: refftabmax(:,:) |
|---|
| 47 | |
|---|
| 48 | !$OMP THREADPRIVATE(varyingnueff,refftabmin,refftabmax) |
|---|
| 49 | |
|---|
| 50 | c Log of the min and max variance of the interpolation grid |
|---|
| 51 | REAL, PARAMETER :: nuefftabmin = -4.6 |
|---|
| 52 | REAL, PARAMETER :: nuefftabmax = 0. |
|---|
| 53 | c Number of effective radius of the interpolation grid |
|---|
| 54 | INTEGER, PARAMETER :: refftabsize = 100 |
|---|
| 55 | c Number of effective variances of the interpolation grid |
|---|
| 56 | INTEGER, PARAMETER :: nuefftabsize = 100 |
|---|
| 57 | c Interpolation grid indices (reff,nueff) |
|---|
| 58 | INTEGER :: grid_i,grid_j |
|---|
| 59 | c Intermediate variable |
|---|
| 60 | REAL :: var_tmp,var3d_tmp(ngrid,nlayer) |
|---|
| 61 | c Bilinear interpolation factors |
|---|
| 62 | REAL :: kx,ky,k1,k2,k3,k4 |
|---|
| 63 | c Size distribution parameters |
|---|
| 64 | REAL :: sizedistk1,sizedistk2 |
|---|
| 65 | c Pi! |
|---|
| 66 | REAL,SAVE :: pi |
|---|
| 67 | |
|---|
| 68 | !$OMP THREADPRIVATE(pi) |
|---|
| 69 | |
|---|
| 70 | c Variables used by the Gauss-Legendre integration: |
|---|
| 71 | INTEGER radius_id,gausind |
|---|
| 72 | REAL kint |
|---|
| 73 | REAL drad |
|---|
| 74 | INTEGER, PARAMETER :: ngau = 10 |
|---|
| 75 | REAL weightgaus(ngau),radgaus(ngau) |
|---|
| 76 | SAVE weightgaus,radgaus |
|---|
| 77 | |
|---|
| 78 | !$OMP THREADPRIVATE(weightgaus,radgaus) |
|---|
| 79 | |
|---|
| 80 | c DATA weightgaus/.2955242247,.2692667193,.2190863625, |
|---|
| 81 | c & .1494513491,.0666713443/ |
|---|
| 82 | c DATA radgaus/.1488743389,.4333953941,.6794095682, |
|---|
| 83 | c & .8650633666,.9739065285/ |
|---|
| 84 | DATA radgaus/0.07652652113350,0.22778585114165, |
|---|
| 85 | & 0.37370608871528,0.51086700195146, |
|---|
| 86 | & 0.63605368072468,0.74633190646476, |
|---|
| 87 | & 0.83911697181213,0.91223442826796, |
|---|
| 88 | & 0.96397192726078,0.99312859919241/ |
|---|
| 89 | |
|---|
| 90 | DATA weightgaus/0.15275338723120,0.14917298659407, |
|---|
| 91 | & 0.14209610937519,0.13168863843930, |
|---|
| 92 | & 0.11819453196154,0.10193011980823, |
|---|
| 93 | & 0.08327674160932,0.06267204829828, |
|---|
| 94 | & 0.04060142982019,0.01761400714091/ |
|---|
| 95 | |
|---|
| 96 | c Indices |
|---|
| 97 | INTEGER :: i,j,k,l,m,iaer,idomain |
|---|
| 98 | INTEGER :: ig,lg,chg |
|---|
| 99 | |
|---|
| 100 | c Local saved variables |
|---|
| 101 | c --------------------- |
|---|
| 102 | |
|---|
| 103 | c Radius axis of the interpolation grid |
|---|
| 104 | REAL,SAVE,ALLOCATABLE :: refftab(:,:,:) |
|---|
| 105 | c Variance axis of the interpolation grid |
|---|
| 106 | REAL,SAVE,ALLOCATABLE :: nuefftab(:,:,:) |
|---|
| 107 | c Volume ratio of the grid |
|---|
| 108 | REAL,SAVE,ALLOCATABLE :: logvratgrid(:,:) |
|---|
| 109 | |
|---|
| 110 | !$OMP THREADPRIVATE(refftab,nuefftab,logvratgrid) |
|---|
| 111 | |
|---|
| 112 | c Grid used to remember which calculation is done |
|---|
| 113 | LOGICAL,SAVE,ALLOCATABLE :: checkgrid(:,:,:,:) |
|---|
| 114 | c Optical properties of the grid (VISIBLE) |
|---|
| 115 | REAL,SAVE,ALLOCATABLE :: qsqrefVISgrid(:,:,:,:) |
|---|
| 116 | REAL,SAVE,ALLOCATABLE :: qextVISgrid(:,:,:,:) |
|---|
| 117 | REAL,SAVE,ALLOCATABLE :: qscatVISgrid(:,:,:,:) |
|---|
| 118 | REAL,SAVE,ALLOCATABLE :: omegVISgrid(:,:,:,:) |
|---|
| 119 | REAL,SAVE,ALLOCATABLE :: gVISgrid(:,:,:,:) |
|---|
| 120 | |
|---|
| 121 | !$OMP THREADPRIVATE(checkgrid,qsqrefVISgrid,qextVISgrid, |
|---|
| 122 | !$OMP& qscatVISgrid,omegVISgrid,gVISgrid) |
|---|
| 123 | |
|---|
| 124 | |
|---|
| 125 | c Optical properties of the grid (INFRARED) |
|---|
| 126 | REAL,SAVE,ALLOCATABLE :: qsqrefIRgrid(:,:,:,:) |
|---|
| 127 | REAL,SAVE,ALLOCATABLE :: qextIRgrid(:,:,:,:) |
|---|
| 128 | REAL,SAVE,ALLOCATABLE :: qscatIRgrid(:,:,:,:) |
|---|
| 129 | REAL,SAVE,ALLOCATABLE :: omegIRgrid(:,:,:,:) |
|---|
| 130 | REAL,SAVE,ALLOCATABLE :: gIRgrid(:,:,:,:) |
|---|
| 131 | |
|---|
| 132 | !$OMP THREADPRIVATE(qsqrefIRgrid,qextIRgrid,qscatIRgrid, |
|---|
| 133 | !$OMP& omegIRgrid,gIRgrid) |
|---|
| 134 | |
|---|
| 135 | c Optical properties of the grid (REFERENCE WAVELENGTHS) |
|---|
| 136 | REAL,SAVE,ALLOCATABLE :: qrefVISgrid(:,:,:) |
|---|
| 137 | REAL,SAVE,ALLOCATABLE :: qscatrefVISgrid(:,:,:) |
|---|
| 138 | REAL,SAVE,ALLOCATABLE :: qrefIRgrid(:,:,:) |
|---|
| 139 | REAL,SAVE,ALLOCATABLE :: qscatrefIRgrid(:,:,:) |
|---|
| 140 | REAL,SAVE,ALLOCATABLE :: omegrefVISgrid(:,:,:) |
|---|
| 141 | REAL,SAVE,ALLOCATABLE :: omegrefIRgrid(:,:,:) |
|---|
| 142 | c Firstcall |
|---|
| 143 | LOGICAL,SAVE :: firstcall = .true. |
|---|
| 144 | |
|---|
| 145 | !$OMP THREADPRIVATE(qrefVISgrid,qscatrefVISgrid,qrefIRgrid , |
|---|
| 146 | !$OMP& qscatrefIRgrid,omegrefVISgrid,omegrefIRgrid, |
|---|
| 147 | !$OMP& firstcall) |
|---|
| 148 | |
|---|
| 149 | c Variables used by the Gauss-Legendre integration: |
|---|
| 150 | REAL,SAVE,ALLOCATABLE :: normd(:,:,:,:) |
|---|
| 151 | REAL,SAVE,ALLOCATABLE :: dista(:,:,:,:,:) |
|---|
| 152 | REAL,SAVE,ALLOCATABLE :: distb(:,:,:,:,:) |
|---|
| 153 | REAL,SAVE,ALLOCATABLE :: radGAUSa(:,:,:) |
|---|
| 154 | REAL,SAVE,ALLOCATABLE :: radGAUSb(:,:,:) |
|---|
| 155 | |
|---|
| 156 | !$OMP THREADPRIVATE(normd,dista,distb,radGAUSa,radGAUSb) |
|---|
| 157 | |
|---|
| 158 | REAL,SAVE,ALLOCATABLE :: qsqrefVISa(:,:,:) |
|---|
| 159 | REAL,SAVE,ALLOCATABLE :: qrefVISa(:,:) |
|---|
| 160 | REAL,SAVE,ALLOCATABLE :: qsqrefVISb(:,:,:) |
|---|
| 161 | REAL,SAVE,ALLOCATABLE :: qrefVISb(:,:) |
|---|
| 162 | REAL,SAVE,ALLOCATABLE :: omegVISa(:,:,:) |
|---|
| 163 | REAL,SAVE,ALLOCATABLE :: omegrefVISa(:,:) |
|---|
| 164 | REAL,SAVE,ALLOCATABLE :: omegVISb(:,:,:) |
|---|
| 165 | REAL,SAVE,ALLOCATABLE :: omegrefVISb(:,:) |
|---|
| 166 | REAL,SAVE,ALLOCATABLE :: gVISa(:,:,:) |
|---|
| 167 | REAL,SAVE,ALLOCATABLE :: gVISb(:,:,:) |
|---|
| 168 | |
|---|
| 169 | !$OMP THREADPRIVATE(qsqrefVISa,qrefVISa,qsqrefVISb,qrefVISb, |
|---|
| 170 | !$OMP& omegVISa,omegrefVISa,omegVISb,omegrefVISb, |
|---|
| 171 | !$OMP& gVISa,gVISb) |
|---|
| 172 | |
|---|
| 173 | REAL,SAVE,ALLOCATABLE :: qsqrefIRa(:,:,:) |
|---|
| 174 | REAL,SAVE,ALLOCATABLE :: qrefIRa(:,:) |
|---|
| 175 | REAL,SAVE,ALLOCATABLE :: qsqrefIRb(:,:,:) |
|---|
| 176 | REAL,SAVE,ALLOCATABLE :: qrefIRb(:,:) |
|---|
| 177 | REAL,SAVE,ALLOCATABLE :: omegIRa(:,:,:) |
|---|
| 178 | REAL,SAVE,ALLOCATABLE :: omegrefIRa(:,:) |
|---|
| 179 | REAL,SAVE,ALLOCATABLE :: omegIRb(:,:,:) |
|---|
| 180 | REAL,SAVE,ALLOCATABLE :: omegrefIRb(:,:) |
|---|
| 181 | REAL,SAVE,ALLOCATABLE :: gIRa(:,:,:) |
|---|
| 182 | REAL,SAVE,ALLOCATABLE :: gIRb(:,:,:) |
|---|
| 183 | |
|---|
| 184 | !$OMP THREADPRIVATE(qsqrefIRa,qrefIRa,qsqrefIRb,qrefIRb, |
|---|
| 185 | !$OMP& omegIRa,omegrefIRa,omegIRb,omegrefIRb,gIRa,gIRb) |
|---|
| 186 | |
|---|
| 187 | REAL :: radiusm |
|---|
| 188 | REAL :: radiusr |
|---|
| 189 | |
|---|
| 190 | c Inputs |
|---|
| 191 | c ------ |
|---|
| 192 | |
|---|
| 193 | INTEGER,INTENT(IN) :: ngrid,nlayer |
|---|
| 194 | c Aerosol effective radius used for radiative transfer (meter) |
|---|
| 195 | REAL,INTENT(IN) :: reffrad(ngrid,nlayer,naerkind) |
|---|
| 196 | c Aerosol effective variance used for radiative transfer (n.u.) |
|---|
| 197 | REAL,INTENT(IN) :: nueffrad(ngrid,nlayer,naerkind) |
|---|
| 198 | |
|---|
| 199 | c Outputs |
|---|
| 200 | c ------- |
|---|
| 201 | |
|---|
| 202 | REAL,INTENT(OUT) :: QVISsQREF3d(ngrid,nlayer,nsun,naerkind) |
|---|
| 203 | REAL,INTENT(OUT) :: omegaVIS3d(ngrid,nlayer,nsun,naerkind) |
|---|
| 204 | REAL,INTENT(OUT) :: gVIS3d(ngrid,nlayer,nsun,naerkind) |
|---|
| 205 | |
|---|
| 206 | REAL,INTENT(OUT) :: QIRsQREF3d(ngrid,nlayer,nir,naerkind) |
|---|
| 207 | REAL,INTENT(OUT) :: omegaIR3d(ngrid,nlayer,nir,naerkind) |
|---|
| 208 | REAL,INTENT(OUT) :: gIR3d(ngrid,nlayer,nir,naerkind) |
|---|
| 209 | |
|---|
| 210 | REAL,INTENT(OUT) :: QREFvis3d(ngrid,nlayer,naerkind) |
|---|
| 211 | REAL,INTENT(OUT) :: QREFir3d(ngrid,nlayer,naerkind) |
|---|
| 212 | |
|---|
| 213 | REAL,INTENT(OUT) :: omegaREFvis3d(ngrid,nlayer,naerkind) |
|---|
| 214 | REAL,INTENT(OUT) :: omegaREFir3d(ngrid,nlayer,naerkind) |
|---|
| 215 | |
|---|
| 216 | c Tests |
|---|
| 217 | c ----- |
|---|
| 218 | |
|---|
| 219 | LOGICAL,SAVE :: out_qwg = .false. |
|---|
| 220 | |
|---|
| 221 | !$OMP THREADPRIVATE(out_qwg) |
|---|
| 222 | |
|---|
| 223 | INTEGER, PARAMETER :: out_iaer = 2 |
|---|
| 224 | REAL :: out_qext(ngrid,nlayer) |
|---|
| 225 | REAL :: out_omeg(ngrid,nlayer) |
|---|
| 226 | REAL :: out_g(ngrid,nlayer) |
|---|
| 227 | INTEGER :: out_nchannel |
|---|
| 228 | CHARACTER*1 :: out_str |
|---|
| 229 | |
|---|
| 230 | c Creating the effective radius and variance grid |
|---|
| 231 | c ----------------------------------------------- |
|---|
| 232 | ! AS: OK firstcall absolute |
|---|
| 233 | IF (firstcall) THEN |
|---|
| 234 | |
|---|
| 235 | c 0.0 Allocate all local saved arrays: |
|---|
| 236 | allocate(refftab(refftabsize,naerkind,2)) |
|---|
| 237 | allocate(nuefftab(nuefftabsize,naerkind,2)) |
|---|
| 238 | ! Optical properties of the grid (VISIBLE) |
|---|
| 239 | allocate(qsqrefVISgrid(refftabsize,nuefftabsize,nsun,naerkind)) |
|---|
| 240 | allocate(qextVISgrid(refftabsize,nuefftabsize,nsun,naerkind)) |
|---|
| 241 | allocate(qscatVISgrid(refftabsize,nuefftabsize,nsun,naerkind)) |
|---|
| 242 | allocate(omegVISgrid(refftabsize,nuefftabsize,nsun,naerkind)) |
|---|
| 243 | allocate(gVISgrid(refftabsize,nuefftabsize,nsun,naerkind)) |
|---|
| 244 | ! Optical properties of the grid (INFRARED) |
|---|
| 245 | allocate(qsqrefIRgrid(refftabsize,nuefftabsize,nir,naerkind)) |
|---|
| 246 | allocate(qextIRgrid(refftabsize,nuefftabsize,nir,naerkind)) |
|---|
| 247 | allocate(qscatIRgrid(refftabsize,nuefftabsize,nir,naerkind)) |
|---|
| 248 | allocate(omegIRgrid(refftabsize,nuefftabsize,nir,naerkind)) |
|---|
| 249 | allocate(gIRgrid(refftabsize,nuefftabsize,nir,naerkind)) |
|---|
| 250 | |
|---|
| 251 | allocate(qsqrefVISa(nsun,ngau,naerkind)) |
|---|
| 252 | allocate(qrefVISa(ngau,naerkind)) |
|---|
| 253 | allocate(qsqrefVISb(nsun,ngau,naerkind)) |
|---|
| 254 | allocate(qrefVISb(ngau,naerkind)) |
|---|
| 255 | allocate(omegVISa(nsun,ngau,naerkind)) |
|---|
| 256 | allocate(omegrefVISa(ngau,naerkind)) |
|---|
| 257 | allocate(omegVISb(nsun,ngau,naerkind)) |
|---|
| 258 | allocate(omegrefVISb(ngau,naerkind)) |
|---|
| 259 | allocate(gVISa(nsun,ngau,naerkind)) |
|---|
| 260 | allocate(gVISb(nsun,ngau,naerkind)) |
|---|
| 261 | |
|---|
| 262 | allocate(qsqrefIRa(nir,ngau,naerkind)) |
|---|
| 263 | allocate(qrefIRa(ngau,naerkind)) |
|---|
| 264 | allocate(qsqrefIRb(nir,ngau,naerkind)) |
|---|
| 265 | allocate(qrefIRb(ngau,naerkind)) |
|---|
| 266 | allocate(omegIRa(nir,ngau,naerkind)) |
|---|
| 267 | allocate(omegrefIRa(ngau,naerkind)) |
|---|
| 268 | allocate(omegIRb(nir,ngau,naerkind)) |
|---|
| 269 | allocate(omegrefIRb(ngau,naerkind)) |
|---|
| 270 | allocate(gIRa(nir,ngau,naerkind)) |
|---|
| 271 | allocate(gIRb(nir,ngau,naerkind)) |
|---|
| 272 | |
|---|
| 273 | allocate(qrefVISgrid(refftabsize,nuefftabsize,naerkind)) |
|---|
| 274 | allocate(qscatrefVISgrid(refftabsize,nuefftabsize,naerkind)) |
|---|
| 275 | allocate(qrefIRgrid(refftabsize,nuefftabsize,naerkind)) |
|---|
| 276 | allocate(qscatrefIRgrid(refftabsize,nuefftabsize,naerkind)) |
|---|
| 277 | allocate(omegrefVISgrid(refftabsize,nuefftabsize,naerkind)) |
|---|
| 278 | allocate(omegrefIRgrid(refftabsize,nuefftabsize,naerkind)) |
|---|
| 279 | |
|---|
| 280 | allocate(normd(refftabsize,nuefftabsize,naerkind,2)) |
|---|
| 281 | allocate(dista(refftabsize,nuefftabsize,naerkind,2,ngau)) |
|---|
| 282 | allocate(distb(refftabsize,nuefftabsize,naerkind,2,ngau)) |
|---|
| 283 | allocate(radGAUSa(ngau,naerkind,2)) |
|---|
| 284 | allocate(radGAUSb(ngau,naerkind,2)) |
|---|
| 285 | |
|---|
| 286 | allocate(checkgrid(refftabsize,nuefftabsize,naerkind,2)) |
|---|
| 287 | |
|---|
| 288 | allocate(logvratgrid(naerkind,2)) |
|---|
| 289 | |
|---|
| 290 | allocate(refftabmin(naerkind,2)) |
|---|
| 291 | allocate(refftabmax(naerkind,2)) |
|---|
| 292 | |
|---|
| 293 | allocate(varyingnueff(naerkind)) |
|---|
| 294 | |
|---|
| 295 | checkgrid(1:refftabsize,1:nuefftabsize,1:naerkind,1:2) = .false. |
|---|
| 296 | varyingnueff(1:naerkind) = .false. |
|---|
| 297 | |
|---|
| 298 | c 0.1 Pi! |
|---|
| 299 | pi = 2. * asin(1.e0) |
|---|
| 300 | |
|---|
| 301 | WRITE(*,*) "aeroptproperties: interpolation grid" |
|---|
| 302 | DO iaer = 1, naerkind ! Loop on aerosol kind |
|---|
| 303 | DO idomain = 1, 2 ! Loop on visible or infrared channel |
|---|
| 304 | |
|---|
| 305 | c 0.2 Effective radius |
|---|
| 306 | radiusm= |
|---|
| 307 | & 0.5*(radiustab(iaer,idomain,nsize(iaer,idomain))+ |
|---|
| 308 | & radiustab(iaer,idomain,1)) |
|---|
| 309 | radiusr= |
|---|
| 310 | & 0.5*(radiustab(iaer,idomain,nsize(iaer,idomain))- |
|---|
| 311 | & radiustab(iaer,idomain,1)) |
|---|
| 312 | refftabmin(iaer,idomain) = |
|---|
| 313 | & radiusm-radiusr*radgaus(ngau) |
|---|
| 314 | refftabmax(iaer,idomain) = |
|---|
| 315 | & radiustab(iaer,idomain,nsize(iaer,idomain)) |
|---|
| 316 | |
|---|
| 317 | WRITE(*,*) "Scatterer: ",iaer |
|---|
| 318 | WRITE(*,*) "Domain: ",idomain |
|---|
| 319 | WRITE(*,*) "Min radius (m): ", refftabmin(iaer,idomain) |
|---|
| 320 | WRITE(*,*) "Max radius (m): ", refftabmax(iaer,idomain) |
|---|
| 321 | |
|---|
| 322 | refftab(1,iaer,idomain) = |
|---|
| 323 | & refftabmin(iaer,idomain) |
|---|
| 324 | refftab(refftabsize,iaer,idomain) = |
|---|
| 325 | & refftabmax(iaer,idomain) |
|---|
| 326 | |
|---|
| 327 | logvratgrid(iaer,idomain) = |
|---|
| 328 | & log(refftabmax(iaer,idomain)/ |
|---|
| 329 | & refftabmin(iaer,idomain)) / |
|---|
| 330 | & float(refftabsize-1)*3. |
|---|
| 331 | do i = 2, refftabsize-1 |
|---|
| 332 | refftab(i,iaer,idomain) = refftab(i-1,iaer,idomain)* |
|---|
| 333 | & exp(1./3.*logvratgrid(iaer,idomain)) |
|---|
| 334 | enddo |
|---|
| 335 | |
|---|
| 336 | c 0.3 Effective variance |
|---|
| 337 | do i = 0, nuefftabsize-1 |
|---|
| 338 | nuefftab(i+1,iaer,idomain) = exp( nuefftabmin + |
|---|
| 339 | & i*(nuefftabmax-nuefftabmin)/(nuefftabsize-1) ) |
|---|
| 340 | enddo |
|---|
| 341 | |
|---|
| 342 | ENDDO |
|---|
| 343 | ENDDO |
|---|
| 344 | firstcall = .false. |
|---|
| 345 | ENDIF |
|---|
| 346 | |
|---|
| 347 | DO iaer = 1, naerkind ! Loop on aerosol kind |
|---|
| 348 | |
|---|
| 349 | IF ( (nsize(iaer,1).EQ.1).AND.(nsize(iaer,2).EQ.1) ) THEN |
|---|
| 350 | c================================================================== |
|---|
| 351 | c If there is one single particle size, optical |
|---|
| 352 | c properties of the considered aerosol are homogeneous |
|---|
| 353 | DO lg = 1, nlayer |
|---|
| 354 | DO ig = 1, ngrid |
|---|
| 355 | DO chg = 1, nsun |
|---|
| 356 | QVISsQREF3d(ig,lg,chg,iaer)=QVISsQREF(chg,iaer,1) |
|---|
| 357 | omegaVIS3d(ig,lg,chg,iaer)=omegaVIS(chg,iaer,1) |
|---|
| 358 | gVIS3d(ig,lg,chg,iaer)=gVIS(chg,iaer,1) |
|---|
| 359 | ENDDO |
|---|
| 360 | DO chg = 1, nir |
|---|
| 361 | QIRsQREF3d(ig,lg,chg,iaer)=QIRsQREF(chg,iaer,1) |
|---|
| 362 | omegaIR3d(ig,lg,chg,iaer)=omegaIR(chg,iaer,1) |
|---|
| 363 | gIR3d(ig,lg,chg,iaer)=gIR(chg,iaer,1) |
|---|
| 364 | ENDDO |
|---|
| 365 | QREFvis3d(ig,lg,iaer)=QREFvis(iaer,1) |
|---|
| 366 | QREFir3d(ig,lg,iaer)=QREFir(iaer,1) |
|---|
| 367 | omegaREFvis3d(ig,lg,iaer)=omegaREFvis(iaer,1) |
|---|
| 368 | omegaREFir3d(ig,lg,iaer)=omegaREFir(iaer,1) |
|---|
| 369 | ENDDO |
|---|
| 370 | ENDDO |
|---|
| 371 | ELSE ! Varying effective radius and variance |
|---|
| 372 | DO idomain = 1, 2 ! Loop on visible or infrared channel |
|---|
| 373 | c================================================================== |
|---|
| 374 | |
|---|
| 375 | c 1.1 Radius middle point and range for Gauss integration |
|---|
| 376 | radiusm= |
|---|
| 377 | & 0.5*(radiustab(iaer,idomain,nsize(iaer,idomain)) + |
|---|
| 378 | & radiustab(iaer,idomain,1)) |
|---|
| 379 | radiusr= |
|---|
| 380 | & 0.5*(radiustab(iaer,idomain,nsize(iaer,idomain)) - |
|---|
| 381 | & radiustab(iaer,idomain,1)) |
|---|
| 382 | |
|---|
| 383 | c 1.2 Interpolating data at the Gauss quadrature points: |
|---|
| 384 | DO gausind=1,ngau |
|---|
| 385 | drad=radiusr*radgaus(gausind) |
|---|
| 386 | radGAUSa(gausind,iaer,idomain)=radiusm-drad |
|---|
| 387 | |
|---|
| 388 | radius_id=minloc(abs(radiustab(iaer,idomain,:) - |
|---|
| 389 | & (radiusm-drad)),DIM=1) |
|---|
| 390 | IF ((radiustab(iaer,idomain,radius_id) - |
|---|
| 391 | & (radiusm-drad)).GT.0) THEN |
|---|
| 392 | radius_id=radius_id-1 |
|---|
| 393 | ENDIF |
|---|
| 394 | IF (radius_id.GE.nsize(iaer,idomain)) THEN |
|---|
| 395 | radius_id=nsize(iaer,idomain)-1 |
|---|
| 396 | kint = 1. |
|---|
| 397 | ELSEIF (radius_id.LT.1) THEN |
|---|
| 398 | radius_id=1 |
|---|
| 399 | kint = 0. |
|---|
| 400 | ELSE |
|---|
| 401 | kint = ( (radiusm-drad) - |
|---|
| 402 | & radiustab(iaer,idomain,radius_id) ) / |
|---|
| 403 | & ( radiustab(iaer,idomain,radius_id+1) - |
|---|
| 404 | & radiustab(iaer,idomain,radius_id) ) |
|---|
| 405 | ENDIF |
|---|
| 406 | IF (idomain.EQ.1) THEN ! VISIBLE DOMAIN ----------------- |
|---|
| 407 | DO m=1,nsun |
|---|
| 408 | qsqrefVISa(m,gausind,iaer)= |
|---|
| 409 | & (1-kint)*QVISsQREF(m,iaer,radius_id) + |
|---|
| 410 | & kint*QVISsQREF(m,iaer,radius_id+1) |
|---|
| 411 | omegVISa(m,gausind,iaer)= |
|---|
| 412 | & (1-kint)*omegaVIS(m,iaer,radius_id) + |
|---|
| 413 | & kint*omegaVIS(m,iaer,radius_id+1) |
|---|
| 414 | gVISa(m,gausind,iaer)= |
|---|
| 415 | & (1-kint)*gVIS(m,iaer,radius_id) + |
|---|
| 416 | & kint*gVIS(m,iaer,radius_id+1) |
|---|
| 417 | ENDDO |
|---|
| 418 | qrefVISa(gausind,iaer)= |
|---|
| 419 | & (1-kint)*QREFvis(iaer,radius_id) + |
|---|
| 420 | & kint*QREFvis(iaer,radius_id+1) |
|---|
| 421 | omegrefVISa(gausind,iaer)= |
|---|
| 422 | & (1-kint)*omegaREFvis(iaer,radius_id) + |
|---|
| 423 | & kint*omegaREFvis(iaer,radius_id+1) |
|---|
| 424 | ELSE ! INFRARED DOMAIN ---------------------------------- |
|---|
| 425 | DO m=1,nir |
|---|
| 426 | qsqrefIRa(m,gausind,iaer)= |
|---|
| 427 | & (1-kint)*QIRsQREF(m,iaer,radius_id) + |
|---|
| 428 | & kint*QIRsQREF(m,iaer,radius_id+1) |
|---|
| 429 | omegIRa(m,gausind,iaer)= |
|---|
| 430 | & (1-kint)*omegaIR(m,iaer,radius_id) + |
|---|
| 431 | & kint*omegaIR(m,iaer,radius_id+1) |
|---|
| 432 | gIRa(m,gausind,iaer)= |
|---|
| 433 | & (1-kint)*gIR(m,iaer,radius_id) + |
|---|
| 434 | & kint*gIR(m,iaer,radius_id+1) |
|---|
| 435 | ENDDO |
|---|
| 436 | qrefIRa(gausind,iaer)= |
|---|
| 437 | & (1-kint)*QREFir(iaer,radius_id) + |
|---|
| 438 | & kint*QREFir(iaer,radius_id+1) |
|---|
| 439 | omegrefIRa(gausind,iaer)= |
|---|
| 440 | & (1-kint)*omegaREFir(iaer,radius_id) + |
|---|
| 441 | & kint*omegaREFir(iaer,radius_id+1) |
|---|
| 442 | ENDIF |
|---|
| 443 | ENDDO |
|---|
| 444 | |
|---|
| 445 | DO gausind=1,ngau |
|---|
| 446 | drad=radiusr*radgaus(gausind) |
|---|
| 447 | radGAUSb(gausind,iaer,idomain)=radiusm+drad |
|---|
| 448 | |
|---|
| 449 | radius_id=minloc(abs(radiustab(iaer,idomain,:) - |
|---|
| 450 | & (radiusm+drad)),DIM=1) |
|---|
| 451 | IF ((radiustab(iaer,idomain,radius_id) - |
|---|
| 452 | & (radiusm+drad)).GT.0) THEN |
|---|
| 453 | radius_id=radius_id-1 |
|---|
| 454 | ENDIF |
|---|
| 455 | IF (radius_id.GE.nsize(iaer,idomain)) THEN |
|---|
| 456 | radius_id=nsize(iaer,idomain)-1 |
|---|
| 457 | kint = 1. |
|---|
| 458 | ELSEIF (radius_id.LT.1) THEN |
|---|
| 459 | radius_id=1 |
|---|
| 460 | kint = 0. |
|---|
| 461 | ELSE |
|---|
| 462 | kint = ( (radiusm+drad) - |
|---|
| 463 | & radiustab(iaer,idomain,radius_id) ) / |
|---|
| 464 | & ( radiustab(iaer,idomain,radius_id+1) - |
|---|
| 465 | & radiustab(iaer,idomain,radius_id) ) |
|---|
| 466 | ENDIF |
|---|
| 467 | IF (idomain.EQ.1) THEN ! VISIBLE DOMAIN ----------------- |
|---|
| 468 | DO m=1,nsun |
|---|
| 469 | qsqrefVISb(m,gausind,iaer)= |
|---|
| 470 | & (1-kint)*QVISsQREF(m,iaer,radius_id) + |
|---|
| 471 | & kint*QVISsQREF(m,iaer,radius_id+1) |
|---|
| 472 | omegVISb(m,gausind,iaer)= |
|---|
| 473 | & (1-kint)*omegaVIS(m,iaer,radius_id) + |
|---|
| 474 | & kint*omegaVIS(m,iaer,radius_id+1) |
|---|
| 475 | gVISb(m,gausind,iaer)= |
|---|
| 476 | & (1-kint)*gVIS(m,iaer,radius_id) + |
|---|
| 477 | & kint*gVIS(m,iaer,radius_id+1) |
|---|
| 478 | ENDDO |
|---|
| 479 | qrefVISb(gausind,iaer)= |
|---|
| 480 | & (1-kint)*QREFvis(iaer,radius_id) + |
|---|
| 481 | & kint*QREFvis(iaer,radius_id+1) |
|---|
| 482 | omegrefVISb(gausind,iaer)= |
|---|
| 483 | & (1-kint)*omegaREFvis(iaer,radius_id) + |
|---|
| 484 | & kint*omegaREFvis(iaer,radius_id+1) |
|---|
| 485 | ELSE ! INFRARED DOMAIN ---------------------------------- |
|---|
| 486 | DO m=1,nir |
|---|
| 487 | qsqrefIRb(m,gausind,iaer)= |
|---|
| 488 | & (1-kint)*QIRsQREF(m,iaer,radius_id) + |
|---|
| 489 | & kint*QIRsQREF(m,iaer,radius_id+1) |
|---|
| 490 | omegIRb(m,gausind,iaer)= |
|---|
| 491 | & (1-kint)*omegaIR(m,iaer,radius_id) + |
|---|
| 492 | & kint*omegaIR(m,iaer,radius_id+1) |
|---|
| 493 | gIRb(m,gausind,iaer)= |
|---|
| 494 | & (1-kint)*gIR(m,iaer,radius_id) + |
|---|
| 495 | & kint*gIR(m,iaer,radius_id+1) |
|---|
| 496 | ENDDO |
|---|
| 497 | qrefIRb(gausind,iaer)= |
|---|
| 498 | & (1-kint)*QREFir(iaer,radius_id) + |
|---|
| 499 | & kint*QREFir(iaer,radius_id+1) |
|---|
| 500 | omegrefIRb(gausind,iaer)= |
|---|
| 501 | & (1-kint)*omegaREFir(iaer,radius_id) + |
|---|
| 502 | & kint*omegaREFir(iaer,radius_id+1) |
|---|
| 503 | ENDIF |
|---|
| 504 | ENDDO |
|---|
| 505 | |
|---|
| 506 | c================================================================== |
|---|
| 507 | IF ( .NOT.varyingnueff(iaer) ) THEN ! CONSTANT NUEFF |
|---|
| 508 | c================================================================== |
|---|
| 509 | c 2. Compute the scattering parameters using linear |
|---|
| 510 | c interpolation over grain sizes and constant nueff |
|---|
| 511 | c --------------------------------------------------- |
|---|
| 512 | |
|---|
| 513 | DO lg = 1,nlayer |
|---|
| 514 | DO ig = 1,ngrid |
|---|
| 515 | c 2.1 Effective radius index and kx calculation |
|---|
| 516 | var_tmp=reffrad(ig,lg,iaer)/refftabmin(iaer,idomain) |
|---|
| 517 | var_tmp=log(var_tmp)*3. |
|---|
| 518 | var_tmp=var_tmp/logvratgrid(iaer,idomain)+1. |
|---|
| 519 | grid_i=floor(var_tmp) |
|---|
| 520 | IF (grid_i.GE.refftabsize) THEN |
|---|
| 521 | c WRITE(*,*) 'Warning: particle size in grid box #' |
|---|
| 522 | c WRITE(*,*) ig,' is too large to be used by the ' |
|---|
| 523 | c WRITE(*,*) 'radiative transfer; please extend the ' |
|---|
| 524 | c WRITE(*,*) 'interpolation grid to larger grain sizes.' |
|---|
| 525 | grid_i=refftabsize-1 |
|---|
| 526 | kx = 1. |
|---|
| 527 | ELSEIF (grid_i.LT.1) THEN |
|---|
| 528 | c WRITE(*,*) 'Warning: particle size in grid box #' |
|---|
| 529 | c WRITE(*,*) ig,' is too small to be used by the ' |
|---|
| 530 | c WRITE(*,*) 'radiative transfer; please extend the ' |
|---|
| 531 | c WRITE(*,*) 'interpolation grid to smaller grain sizes.' |
|---|
| 532 | grid_i=1 |
|---|
| 533 | kx = 0. |
|---|
| 534 | ELSE |
|---|
| 535 | kx = ( reffrad(ig,lg,iaer)- |
|---|
| 536 | & refftab(grid_i,iaer,idomain) ) / |
|---|
| 537 | & ( refftab(grid_i+1,iaer,idomain)- |
|---|
| 538 | & refftab(grid_i,iaer,idomain) ) |
|---|
| 539 | ENDIF |
|---|
| 540 | c 2.3 Integration |
|---|
| 541 | DO j=grid_i,grid_i+1 |
|---|
| 542 | c 2.3.1 Check if the calculation has been done |
|---|
| 543 | IF (.NOT.checkgrid(j,1,iaer,idomain)) THEN |
|---|
| 544 | c 2.3.2 Log-normal dist., r_g and sigma_g are defined |
|---|
| 545 | c in [hansen_1974], "Light scattering in planetary |
|---|
| 546 | c atmospheres", Space Science Reviews 16 527-610. |
|---|
| 547 | c Here, sizedistk1=r_g and sizedistk2=sigma_g^2 |
|---|
| 548 | sizedistk2 = log(1.+nueffrad(1,1,iaer)) |
|---|
| 549 | sizedistk1 = exp(2.5*sizedistk2) |
|---|
| 550 | sizedistk1 = refftab(j,iaer,idomain) / sizedistk1 |
|---|
| 551 | |
|---|
| 552 | normd(j,1,iaer,idomain) = 1e-30 |
|---|
| 553 | DO gausind=1,ngau |
|---|
| 554 | drad=radiusr*radgaus(gausind) |
|---|
| 555 | dista(j,1,iaer,idomain,gausind) = |
|---|
| 556 | & LOG((radiusm-drad)/sizedistk1) |
|---|
| 557 | dista(j,1,iaer,idomain,gausind) = |
|---|
| 558 | & EXP(-dista(j,1,iaer,idomain,gausind) * |
|---|
| 559 | & dista(j,1,iaer,idomain,gausind) * |
|---|
| 560 | & 0.5e0/sizedistk2)/(radiusm-drad) |
|---|
| 561 | dista(j,1,iaer,idomain,gausind) = |
|---|
| 562 | & dista(j,1,iaer,idomain,gausind) / |
|---|
| 563 | & (sqrt(2e0*pi*sizedistk2)) |
|---|
| 564 | |
|---|
| 565 | distb(j,1,iaer,idomain,gausind) = |
|---|
| 566 | & LOG((radiusm+drad)/sizedistk1) |
|---|
| 567 | distb(j,1,iaer,idomain,gausind) = |
|---|
| 568 | & EXP(-distb(j,1,iaer,idomain,gausind) * |
|---|
| 569 | & distb(j,1,iaer,idomain,gausind) * |
|---|
| 570 | & 0.5e0/sizedistk2)/(radiusm+drad) |
|---|
| 571 | distb(j,1,iaer,idomain,gausind) = |
|---|
| 572 | & distb(j,1,iaer,idomain,gausind) / |
|---|
| 573 | & (sqrt(2e0*pi*sizedistk2)) |
|---|
| 574 | |
|---|
| 575 | normd(j,1,iaer,idomain)=normd(j,1,iaer,idomain) + |
|---|
| 576 | & weightgaus(gausind) * |
|---|
| 577 | & ( |
|---|
| 578 | & distb(j,1,iaer,idomain,gausind) * pi * |
|---|
| 579 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 580 | & radGAUSb(gausind,iaer,idomain) + |
|---|
| 581 | & dista(j,1,iaer,idomain,gausind) * pi * |
|---|
| 582 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 583 | & radGAUSa(gausind,iaer,idomain) |
|---|
| 584 | & ) |
|---|
| 585 | ENDDO |
|---|
| 586 | |
|---|
| 587 | IF (normd(j,1,iaer,idomain).EQ.1e-30) THEN |
|---|
| 588 | WRITE(*,*)"normd:", normd(j,1,iaer,idomain) |
|---|
| 589 | WRITE(*,*)"Risk of division by 0 (aeroptproperties.F)" |
|---|
| 590 | WRITE(*,*)"Check the size of the interpolation grid." |
|---|
| 591 | CALL abort_physic("aeroptproperties", |
|---|
| 592 | & "Risk of division by 0",1) |
|---|
| 593 | ENDIF |
|---|
| 594 | IF (idomain.EQ.1) THEN ! VISIBLE DOMAIN ----------- |
|---|
| 595 | c 2.3.3.vis Initialization |
|---|
| 596 | qsqrefVISgrid(j,1,:,iaer)=0. |
|---|
| 597 | qextVISgrid(j,1,:,iaer)=0. |
|---|
| 598 | qscatVISgrid(j,1,:,iaer)=0. |
|---|
| 599 | omegVISgrid(j,1,:,iaer)=0. |
|---|
| 600 | gVISgrid(j,1,:,iaer)=0. |
|---|
| 601 | qrefVISgrid(j,1,iaer)=0. |
|---|
| 602 | qscatrefVISgrid(j,1,iaer)=0. |
|---|
| 603 | omegrefVISgrid(j,1,iaer)=0. |
|---|
| 604 | |
|---|
| 605 | DO gausind=1,ngau |
|---|
| 606 | DO m=1,nsun |
|---|
| 607 | c Convolution: |
|---|
| 608 | qextVISgrid(j,1,m,iaer) = |
|---|
| 609 | & qextVISgrid(j,1,m,iaer) + |
|---|
| 610 | & weightgaus(gausind) * |
|---|
| 611 | & ( |
|---|
| 612 | & qsqrefVISb(m,gausind,iaer) * |
|---|
| 613 | & qrefVISb(gausind,iaer) * |
|---|
| 614 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 615 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 616 | & distb(j,1,iaer,idomain,gausind) + |
|---|
| 617 | & qsqrefVISa(m,gausind,iaer) * |
|---|
| 618 | & qrefVISa(gausind,iaer) * |
|---|
| 619 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 620 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 621 | & dista(j,1,iaer,idomain,gausind) |
|---|
| 622 | & ) |
|---|
| 623 | qscatVISgrid(j,1,m,iaer) = |
|---|
| 624 | & qscatVISgrid(j,1,m,iaer) + |
|---|
| 625 | & weightgaus(gausind) * |
|---|
| 626 | & ( |
|---|
| 627 | & omegVISb(m,gausind,iaer) * |
|---|
| 628 | & qsqrefVISb(m,gausind,iaer) * |
|---|
| 629 | & qrefVISb(gausind,iaer) * |
|---|
| 630 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 631 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 632 | & distb(j,1,iaer,idomain,gausind) + |
|---|
| 633 | & omegVISa(m,gausind,iaer) * |
|---|
| 634 | & qsqrefVISa(m,gausind,iaer) * |
|---|
| 635 | & qrefVISa(gausind,iaer) * |
|---|
| 636 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 637 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 638 | & dista(j,1,iaer,idomain,gausind) |
|---|
| 639 | & ) |
|---|
| 640 | gVISgrid(j,1,m,iaer) = |
|---|
| 641 | & gVISgrid(j,1,m,iaer) + |
|---|
| 642 | & weightgaus(gausind) * |
|---|
| 643 | & ( |
|---|
| 644 | & omegVISb(m,gausind,iaer) * |
|---|
| 645 | & qsqrefVISb(m,gausind,iaer) * |
|---|
| 646 | & qrefVISb(gausind,iaer) * |
|---|
| 647 | & gVISb(m,gausind,iaer) * |
|---|
| 648 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 649 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 650 | & distb(j,1,iaer,idomain,gausind) + |
|---|
| 651 | & omegVISa(m,gausind,iaer) * |
|---|
| 652 | & qsqrefVISa(m,gausind,iaer) * |
|---|
| 653 | & qrefVISa(gausind,iaer) * |
|---|
| 654 | & gVISa(m,gausind,iaer) * |
|---|
| 655 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 656 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 657 | & dista(j,1,iaer,idomain,gausind) |
|---|
| 658 | & ) |
|---|
| 659 | ENDDO |
|---|
| 660 | qrefVISgrid(j,1,iaer) = |
|---|
| 661 | & qrefVISgrid(j,1,iaer) + |
|---|
| 662 | & weightgaus(gausind) * |
|---|
| 663 | & ( |
|---|
| 664 | & qrefVISb(gausind,iaer) * |
|---|
| 665 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 666 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 667 | & distb(j,1,iaer,idomain,gausind) + |
|---|
| 668 | & qrefVISa(gausind,iaer) * |
|---|
| 669 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 670 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 671 | & dista(j,1,iaer,idomain,gausind) |
|---|
| 672 | & ) |
|---|
| 673 | qscatrefVISgrid(j,1,iaer) = |
|---|
| 674 | & qscatrefVISgrid(j,1,iaer) + |
|---|
| 675 | & weightgaus(gausind) * |
|---|
| 676 | & ( |
|---|
| 677 | & omegrefVISb(gausind,iaer) * |
|---|
| 678 | & qrefVISb(gausind,iaer) * |
|---|
| 679 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 680 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 681 | & distb(j,1,iaer,idomain,gausind) + |
|---|
| 682 | & omegrefVISa(gausind,iaer) * |
|---|
| 683 | & qrefVISa(gausind,iaer) * |
|---|
| 684 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 685 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 686 | & dista(j,1,iaer,idomain,gausind) |
|---|
| 687 | & ) |
|---|
| 688 | ENDDO |
|---|
| 689 | |
|---|
| 690 | qrefVISgrid(j,1,iaer)=qrefVISgrid(j,1,iaer) / |
|---|
| 691 | & normd(j,1,iaer,idomain) |
|---|
| 692 | qscatrefVISgrid(j,1,iaer)=qscatrefVISgrid(j,1,iaer) / |
|---|
| 693 | & normd(j,1,iaer,idomain) |
|---|
| 694 | omegrefVISgrid(j,1,iaer)=qscatrefVISgrid(j,1,iaer) / |
|---|
| 695 | & qrefVISgrid(j,1,iaer) |
|---|
| 696 | DO m=1,nsun |
|---|
| 697 | qextVISgrid(j,1,m,iaer)=qextVISgrid(j,1,m,iaer) / |
|---|
| 698 | & normd(j,1,iaer,idomain) |
|---|
| 699 | qscatVISgrid(j,1,m,iaer)=qscatVISgrid(j,1,m,iaer) / |
|---|
| 700 | & normd(j,1,iaer,idomain) |
|---|
| 701 | gVISgrid(j,1,m,iaer)=gVISgrid(j,1,m,iaer) / |
|---|
| 702 | & qscatVISgrid(j,1,m,iaer) / |
|---|
| 703 | & normd(j,1,iaer,idomain) |
|---|
| 704 | |
|---|
| 705 | qsqrefVISgrid(j,1,m,iaer)=qextVISgrid(j,1,m,iaer) / |
|---|
| 706 | & qrefVISgrid(j,1,iaer) |
|---|
| 707 | omegVISgrid(j,1,m,iaer)=qscatVISgrid(j,1,m,iaer) / |
|---|
| 708 | & qextVISgrid(j,1,m,iaer) |
|---|
| 709 | ENDDO |
|---|
| 710 | ELSE ! INFRARED DOMAIN ---------- |
|---|
| 711 | c 2.3.3.ir Initialization |
|---|
| 712 | qsqrefIRgrid(j,1,:,iaer)=0. |
|---|
| 713 | qextIRgrid(j,1,:,iaer)=0. |
|---|
| 714 | qscatIRgrid(j,1,:,iaer)=0. |
|---|
| 715 | omegIRgrid(j,1,:,iaer)=0. |
|---|
| 716 | gIRgrid(j,1,:,iaer)=0. |
|---|
| 717 | qrefIRgrid(j,1,iaer)=0. |
|---|
| 718 | qscatrefIRgrid(j,1,iaer)=0. |
|---|
| 719 | omegrefIRgrid(j,1,iaer)=0. |
|---|
| 720 | |
|---|
| 721 | DO gausind=1,ngau |
|---|
| 722 | DO m=1,nir |
|---|
| 723 | c Convolution: |
|---|
| 724 | qextIRgrid(j,1,m,iaer) = |
|---|
| 725 | & qextIRgrid(j,1,m,iaer) + |
|---|
| 726 | & weightgaus(gausind) * |
|---|
| 727 | & ( |
|---|
| 728 | & qsqrefIRb(m,gausind,iaer) * |
|---|
| 729 | & qrefVISb(gausind,iaer) * |
|---|
| 730 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 731 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 732 | & distb(j,1,iaer,idomain,gausind) + |
|---|
| 733 | & qsqrefIRa(m,gausind,iaer) * |
|---|
| 734 | & qrefVISa(gausind,iaer) * |
|---|
| 735 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 736 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 737 | & dista(j,1,iaer,idomain,gausind) |
|---|
| 738 | & ) |
|---|
| 739 | qscatIRgrid(j,1,m,iaer) = |
|---|
| 740 | & qscatIRgrid(j,1,m,iaer) + |
|---|
| 741 | & weightgaus(gausind) * |
|---|
| 742 | & ( |
|---|
| 743 | & omegIRb(m,gausind,iaer) * |
|---|
| 744 | & qsqrefIRb(m,gausind,iaer) * |
|---|
| 745 | & qrefVISb(gausind,iaer) * |
|---|
| 746 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 747 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 748 | & distb(j,1,iaer,idomain,gausind) + |
|---|
| 749 | & omegIRa(m,gausind,iaer) * |
|---|
| 750 | & qsqrefIRa(m,gausind,iaer) * |
|---|
| 751 | & qrefVISa(gausind,iaer) * |
|---|
| 752 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 753 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 754 | & dista(j,1,iaer,idomain,gausind) |
|---|
| 755 | & ) |
|---|
| 756 | gIRgrid(j,1,m,iaer) = |
|---|
| 757 | & gIRgrid(j,1,m,iaer) + |
|---|
| 758 | & weightgaus(gausind) * |
|---|
| 759 | & ( |
|---|
| 760 | & omegIRb(m,gausind,iaer) * |
|---|
| 761 | & qsqrefIRb(m,gausind,iaer) * |
|---|
| 762 | & qrefVISb(gausind,iaer) * |
|---|
| 763 | & gIRb(m,gausind,iaer) * |
|---|
| 764 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 765 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 766 | & distb(j,1,iaer,idomain,gausind) + |
|---|
| 767 | & omegIRa(m,gausind,iaer) * |
|---|
| 768 | & qsqrefIRa(m,gausind,iaer) * |
|---|
| 769 | & qrefVISa(gausind,iaer) * |
|---|
| 770 | & gIRa(m,gausind,iaer) * |
|---|
| 771 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 772 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 773 | & dista(j,1,iaer,idomain,gausind) |
|---|
| 774 | & ) |
|---|
| 775 | ENDDO |
|---|
| 776 | qrefIRgrid(j,1,iaer) = |
|---|
| 777 | & qrefIRgrid(j,1,iaer) + |
|---|
| 778 | & weightgaus(gausind) * |
|---|
| 779 | & ( |
|---|
| 780 | & qrefIRb(gausind,iaer) * |
|---|
| 781 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 782 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 783 | & distb(j,1,iaer,idomain,gausind) + |
|---|
| 784 | & qrefIRa(gausind,iaer) * |
|---|
| 785 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 786 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 787 | & dista(j,1,iaer,idomain,gausind) |
|---|
| 788 | & ) |
|---|
| 789 | qscatrefIRgrid(j,1,iaer) = |
|---|
| 790 | & qscatrefIRgrid(j,1,iaer) + |
|---|
| 791 | & weightgaus(gausind) * |
|---|
| 792 | & ( |
|---|
| 793 | & omegrefIRb(gausind,iaer) * |
|---|
| 794 | & qrefIRb(gausind,iaer) * |
|---|
| 795 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 796 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 797 | & distb(j,1,iaer,idomain,gausind) + |
|---|
| 798 | & omegrefIRa(gausind,iaer) * |
|---|
| 799 | & qrefIRa(gausind,iaer) * |
|---|
| 800 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 801 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 802 | & dista(j,1,iaer,idomain,gausind) |
|---|
| 803 | & ) |
|---|
| 804 | ENDDO |
|---|
| 805 | |
|---|
| 806 | qrefIRgrid(j,1,iaer)=qrefIRgrid(j,1,iaer) / |
|---|
| 807 | & normd(j,1,iaer,idomain) |
|---|
| 808 | qscatrefIRgrid(j,1,iaer)=qscatrefIRgrid(j,1,iaer) / |
|---|
| 809 | & normd(j,1,iaer,idomain) |
|---|
| 810 | omegrefIRgrid(j,1,iaer)=qscatrefIRgrid(j,1,iaer) / |
|---|
| 811 | & qrefIRgrid(j,1,iaer) |
|---|
| 812 | DO m=1,nir |
|---|
| 813 | qextIRgrid(j,1,m,iaer)=qextIRgrid(j,1,m,iaer) / |
|---|
| 814 | & normd(j,1,iaer,idomain) |
|---|
| 815 | qscatIRgrid(j,1,m,iaer)=qscatIRgrid(j,1,m,iaer) / |
|---|
| 816 | & normd(j,1,iaer,idomain) |
|---|
| 817 | gIRgrid(j,1,m,iaer)=gIRgrid(j,1,m,iaer) / |
|---|
| 818 | & qscatIRgrid(j,1,m,iaer) / |
|---|
| 819 | & normd(j,1,iaer,idomain) |
|---|
| 820 | |
|---|
| 821 | qsqrefIRgrid(j,1,m,iaer)=qextIRgrid(j,1,m,iaer) / |
|---|
| 822 | & qrefVISgrid(j,1,iaer) |
|---|
| 823 | omegIRgrid(j,1,m,iaer)=qscatIRgrid(j,1,m,iaer) / |
|---|
| 824 | & qextIRgrid(j,1,m,iaer) |
|---|
| 825 | ENDDO |
|---|
| 826 | ENDIF ! -------------------------- |
|---|
| 827 | checkgrid(j,1,iaer,idomain) = .true. |
|---|
| 828 | ENDIF !checkgrid |
|---|
| 829 | ENDDO !grid_i |
|---|
| 830 | c 2.4 Linear interpolation |
|---|
| 831 | k1 = (1-kx) |
|---|
| 832 | k2 = kx |
|---|
| 833 | IF (idomain.EQ.1) THEN ! VISIBLE ------------------------ |
|---|
| 834 | DO m=1,nsun |
|---|
| 835 | QVISsQREF3d(ig,lg,m,iaer) = |
|---|
| 836 | & k1*qsqrefVISgrid(grid_i,1,m,iaer) + |
|---|
| 837 | & k2*qsqrefVISgrid(grid_i+1,1,m,iaer) |
|---|
| 838 | omegaVIS3d(ig,lg,m,iaer) = |
|---|
| 839 | & k1*omegVISgrid(grid_i,1,m,iaer) + |
|---|
| 840 | & k2*omegVISgrid(grid_i+1,1,m,iaer) |
|---|
| 841 | gVIS3d(ig,lg,m,iaer) = |
|---|
| 842 | & k1*gVISgrid(grid_i,1,m,iaer) + |
|---|
| 843 | & k2*gVISgrid(grid_i+1,1,m,iaer) |
|---|
| 844 | ENDDO !nsun |
|---|
| 845 | QREFvis3d(ig,lg,iaer) = |
|---|
| 846 | & k1*qrefVISgrid(grid_i,1,iaer) + |
|---|
| 847 | & k2*qrefVISgrid(grid_i+1,1,iaer) |
|---|
| 848 | omegaREFvis3d(ig,lg,iaer) = |
|---|
| 849 | & k1*omegrefVISgrid(grid_i,1,iaer) + |
|---|
| 850 | & k2*omegrefVISgrid(grid_i+1,1,iaer) |
|---|
| 851 | ELSE ! INFRARED ----------------------- |
|---|
| 852 | DO m=1,nir |
|---|
| 853 | QIRsQREF3d(ig,lg,m,iaer) = |
|---|
| 854 | & k1*qsqrefIRgrid(grid_i,1,m,iaer) + |
|---|
| 855 | & k2*qsqrefIRgrid(grid_i+1,1,m,iaer) |
|---|
| 856 | omegaIR3d(ig,lg,m,iaer) = |
|---|
| 857 | & k1*omegIRgrid(grid_i,1,m,iaer) + |
|---|
| 858 | & k2*omegIRgrid(grid_i+1,1,m,iaer) |
|---|
| 859 | gIR3d(ig,lg,m,iaer) = |
|---|
| 860 | & k1*gIRgrid(grid_i,1,m,iaer) + |
|---|
| 861 | & k2*gIRgrid(grid_i+1,1,m,iaer) |
|---|
| 862 | ENDDO !nir |
|---|
| 863 | QREFir3d(ig,lg,iaer) = |
|---|
| 864 | & k1*qrefIRgrid(grid_i,1,iaer) + |
|---|
| 865 | & k2*qrefIRgrid(grid_i+1,1,iaer) |
|---|
| 866 | omegaREFir3d(ig,lg,iaer) = |
|---|
| 867 | & k1*omegrefIRgrid(grid_i,1,iaer) + |
|---|
| 868 | & k2*omegrefIRgrid(grid_i+1,1,iaer) |
|---|
| 869 | ENDIF ! -------------------------------- |
|---|
| 870 | ENDDO !nlayer |
|---|
| 871 | ENDDO !ngrid |
|---|
| 872 | c================================================================== |
|---|
| 873 | ELSE ! VARYING NUEFF |
|---|
| 874 | c================================================================== |
|---|
| 875 | c 3. Compute the scattering parameters in each grid box |
|---|
| 876 | c using bilinear interpolation over the grain sizes |
|---|
| 877 | c and the effective variances; |
|---|
| 878 | c ----------------------------------------------------- |
|---|
| 879 | |
|---|
| 880 | DO lg = 1,nlayer |
|---|
| 881 | DO ig = 1,ngrid |
|---|
| 882 | c 3.1 Effective variance index and ky calculation |
|---|
| 883 | var_tmp=log(nueffrad(ig,lg,iaer)) |
|---|
| 884 | grid_j=floor( (nuefftabsize-1)/(nuefftabmax-nuefftabmin)* |
|---|
| 885 | & (var_tmp-nuefftabmin)+1. ) |
|---|
| 886 | IF (grid_j.GE.nuefftabsize) THEN |
|---|
| 887 | c WRITE(*,*) 'Warning: effective variance ' |
|---|
| 888 | c WRITE(*,*) 'is too large to be used by the ' |
|---|
| 889 | c WRITE(*,*) 'radiative transfer; please extend the ' |
|---|
| 890 | c WRITE(*,*) 'interpolation grid to larger values.' |
|---|
| 891 | grid_j=nuefftabsize-1 |
|---|
| 892 | ky = 1. |
|---|
| 893 | ELSEIF (grid_j.LT.1) THEN |
|---|
| 894 | c WRITE(*,*) 'Warning: effective variance ' |
|---|
| 895 | c WRITE(*,*) 'is too small to be used by the ' |
|---|
| 896 | c WRITE(*,*) 'radiative transfer; please extend the ' |
|---|
| 897 | c WRITE(*,*) 'interpolation grid to smaller values.' |
|---|
| 898 | grid_j=1 |
|---|
| 899 | ky = 0. |
|---|
| 900 | ELSE |
|---|
| 901 | ky = ( nueffrad(ig,lg,iaer)- |
|---|
| 902 | & nuefftab(grid_j,iaer,idomain) ) / |
|---|
| 903 | & ( nuefftab(grid_j+1,iaer,idomain)- |
|---|
| 904 | & nuefftab(grid_j,iaer,idomain) ) |
|---|
| 905 | ENDIF |
|---|
| 906 | c 3.2 Effective radius index and kx calculation |
|---|
| 907 | var_tmp=reffrad(ig,lg,iaer)/refftabmin(iaer,idomain) |
|---|
| 908 | var_tmp=log(var_tmp)*3. |
|---|
| 909 | var_tmp=var_tmp/logvratgrid(iaer,idomain)+1. |
|---|
| 910 | grid_i=floor(var_tmp) |
|---|
| 911 | IF (grid_i.GE.refftabsize) THEN |
|---|
| 912 | c WRITE(*,*) 'Warning: particle size in grid box #' |
|---|
| 913 | c WRITE(*,*) ig,' is too large to be used by the ' |
|---|
| 914 | c WRITE(*,*) 'radiative transfer; please extend the ' |
|---|
| 915 | c WRITE(*,*) 'interpolation grid to larger grain sizes.' |
|---|
| 916 | grid_i=refftabsize-1 |
|---|
| 917 | kx = 1. |
|---|
| 918 | ELSEIF (grid_i.LT.1) THEN |
|---|
| 919 | c WRITE(*,*) 'Warning: particle size in grid box #' |
|---|
| 920 | c WRITE(*,*) ig,' is too small to be used by the ' |
|---|
| 921 | c WRITE(*,*) 'radiative transfer; please extend the ' |
|---|
| 922 | c WRITE(*,*) 'interpolation grid to smaller grain sizes.' |
|---|
| 923 | grid_i=1 |
|---|
| 924 | kx = 0. |
|---|
| 925 | ELSE |
|---|
| 926 | kx = ( reffrad(ig,lg,iaer)- |
|---|
| 927 | & refftab(grid_i,iaer,idomain) ) / |
|---|
| 928 | & ( refftab(grid_i+1,iaer,idomain)- |
|---|
| 929 | & refftab(grid_i,iaer,idomain) ) |
|---|
| 930 | ENDIF |
|---|
| 931 | c 3.3 Integration |
|---|
| 932 | DO j=grid_i,grid_i+1 |
|---|
| 933 | DO k=grid_j,grid_j+1 |
|---|
| 934 | c 3.3.1 Check if the calculation has been done |
|---|
| 935 | IF (.NOT.checkgrid(j,k,iaer,idomain)) THEN |
|---|
| 936 | |
|---|
| 937 | c 3.3.2 Log-normal dist., r_g and sigma_g are defined |
|---|
| 938 | c in [hansen_1974], "Light scattering in planetary |
|---|
| 939 | c atmospheres", Space Science Reviews 16 527-610. |
|---|
| 940 | c Here, sizedistk1=r_g and sizedistk2=sigma_g^2 |
|---|
| 941 | sizedistk2 = log(1.+nuefftab(k,iaer,idomain)) |
|---|
| 942 | sizedistk1 = exp(2.5*sizedistk2) |
|---|
| 943 | sizedistk1 = refftab(j,iaer,idomain) / sizedistk1 |
|---|
| 944 | |
|---|
| 945 | normd(j,k,iaer,idomain) = 1e-30 |
|---|
| 946 | DO gausind=1,ngau |
|---|
| 947 | drad=radiusr*radgaus(gausind) |
|---|
| 948 | |
|---|
| 949 | dista(j,k,iaer,idomain,gausind) = |
|---|
| 950 | & LOG((radiusm-drad)/sizedistk1) |
|---|
| 951 | dista(j,k,iaer,idomain,gausind) = |
|---|
| 952 | & EXP(-dista(j,k,iaer,idomain,gausind) * |
|---|
| 953 | & dista(j,k,iaer,idomain,gausind) * |
|---|
| 954 | & 0.5e0/sizedistk2)/(radiusm-drad) |
|---|
| 955 | dista(j,k,iaer,idomain,gausind) = |
|---|
| 956 | & dista(j,k,iaer,idomain,gausind) / |
|---|
| 957 | & (sqrt(2e0*pi*sizedistk2)) |
|---|
| 958 | |
|---|
| 959 | distb(j,k,iaer,idomain,gausind) = |
|---|
| 960 | & LOG((radiusm+drad)/sizedistk1) |
|---|
| 961 | distb(j,k,iaer,idomain,gausind) = |
|---|
| 962 | & EXP(-distb(j,k,iaer,idomain,gausind) * |
|---|
| 963 | & distb(j,k,iaer,idomain,gausind) * |
|---|
| 964 | & 0.5e0/sizedistk2)/(radiusm+drad) |
|---|
| 965 | distb(j,k,iaer,idomain,gausind) = |
|---|
| 966 | & distb(j,k,iaer,idomain,gausind) / |
|---|
| 967 | & (sqrt(2e0*pi*sizedistk2)) |
|---|
| 968 | |
|---|
| 969 | normd(j,k,iaer,idomain)=normd(j,k,iaer,idomain) + |
|---|
| 970 | & weightgaus(gausind) * |
|---|
| 971 | & ( |
|---|
| 972 | & distb(j,k,iaer,idomain,gausind) * pi * |
|---|
| 973 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 974 | & radGAUSb(gausind,iaer,idomain) + |
|---|
| 975 | & dista(j,k,iaer,idomain,gausind) * pi * |
|---|
| 976 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 977 | & radGAUSa(gausind,iaer,idomain) |
|---|
| 978 | & ) |
|---|
| 979 | ENDDO |
|---|
| 980 | IF (normd(j,k,iaer,idomain).EQ.1e-30) THEN |
|---|
| 981 | WRITE(*,*)"normd:", normd(j,k,iaer,idomain) |
|---|
| 982 | WRITE(*,*)"Risk of division by 0 (aeroptproperties.F)" |
|---|
| 983 | WRITE(*,*)"Check the size of the interpolation grid." |
|---|
| 984 | CALL abort_physic("aeroptproperties", |
|---|
| 985 | & "Risk of division by 0",1) |
|---|
| 986 | ENDIF |
|---|
| 987 | IF (idomain.EQ.1) THEN ! VISIBLE DOMAIN ----------- |
|---|
| 988 | c 2.3.3.vis Initialization |
|---|
| 989 | qsqrefVISgrid(j,k,:,iaer)=0. |
|---|
| 990 | qextVISgrid(j,k,:,iaer)=0. |
|---|
| 991 | qscatVISgrid(j,k,:,iaer)=0. |
|---|
| 992 | omegVISgrid(j,k,:,iaer)=0. |
|---|
| 993 | gVISgrid(j,k,:,iaer)=0. |
|---|
| 994 | qrefVISgrid(j,k,iaer)=0. |
|---|
| 995 | qscatrefVISgrid(j,k,iaer)=0. |
|---|
| 996 | omegrefVISgrid(j,k,iaer)=0. |
|---|
| 997 | |
|---|
| 998 | DO gausind=1,ngau |
|---|
| 999 | DO m=1,nsun |
|---|
| 1000 | c Convolution: |
|---|
| 1001 | qextVISgrid(j,k,m,iaer) = |
|---|
| 1002 | & qextVISgrid(j,k,m,iaer) + |
|---|
| 1003 | & weightgaus(gausind) * |
|---|
| 1004 | & ( |
|---|
| 1005 | & qsqrefVISb(m,gausind,iaer) * |
|---|
| 1006 | & qrefVISb(gausind,iaer) * |
|---|
| 1007 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 1008 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 1009 | & distb(j,k,iaer,idomain,gausind) + |
|---|
| 1010 | & qsqrefVISa(m,gausind,iaer) * |
|---|
| 1011 | & qrefVISa(gausind,iaer) * |
|---|
| 1012 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 1013 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 1014 | & dista(j,k,iaer,idomain,gausind) |
|---|
| 1015 | & ) |
|---|
| 1016 | qscatVISgrid(j,k,m,iaer) = |
|---|
| 1017 | & qscatVISgrid(j,k,m,iaer) + |
|---|
| 1018 | & weightgaus(gausind) * |
|---|
| 1019 | & ( |
|---|
| 1020 | & omegVISb(m,gausind,iaer) * |
|---|
| 1021 | & qsqrefVISb(m,gausind,iaer) * |
|---|
| 1022 | & qrefVISb(gausind,iaer) * |
|---|
| 1023 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 1024 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 1025 | & distb(j,k,iaer,idomain,gausind) + |
|---|
| 1026 | & omegVISa(m,gausind,iaer) * |
|---|
| 1027 | & qsqrefVISa(m,gausind,iaer) * |
|---|
| 1028 | & qrefVISa(gausind,iaer) * |
|---|
| 1029 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 1030 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 1031 | & dista(j,k,iaer,idomain,gausind) |
|---|
| 1032 | & ) |
|---|
| 1033 | gVISgrid(j,k,m,iaer) = |
|---|
| 1034 | & gVISgrid(j,k,m,iaer) + |
|---|
| 1035 | & weightgaus(gausind) * |
|---|
| 1036 | & ( |
|---|
| 1037 | & omegVISb(m,gausind,iaer) * |
|---|
| 1038 | & qsqrefVISb(m,gausind,iaer) * |
|---|
| 1039 | & qrefVISb(gausind,iaer) * |
|---|
| 1040 | & gVISb(m,gausind,iaer) * |
|---|
| 1041 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 1042 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 1043 | & distb(j,k,iaer,idomain,gausind) + |
|---|
| 1044 | & omegVISa(m,gausind,iaer) * |
|---|
| 1045 | & qsqrefVISa(m,gausind,iaer) * |
|---|
| 1046 | & qrefVISa(gausind,iaer) * |
|---|
| 1047 | & gVISa(m,gausind,iaer) * |
|---|
| 1048 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 1049 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 1050 | & dista(j,k,iaer,idomain,gausind) |
|---|
| 1051 | & ) |
|---|
| 1052 | ENDDO |
|---|
| 1053 | qrefVISgrid(j,k,iaer) = |
|---|
| 1054 | & qrefVISgrid(j,k,iaer) + |
|---|
| 1055 | & weightgaus(gausind) * |
|---|
| 1056 | & ( |
|---|
| 1057 | & qrefVISb(gausind,iaer) * |
|---|
| 1058 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 1059 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 1060 | & distb(j,k,iaer,idomain,gausind) + |
|---|
| 1061 | & qrefVISa(gausind,iaer) * |
|---|
| 1062 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 1063 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 1064 | & dista(j,k,iaer,idomain,gausind) |
|---|
| 1065 | & ) |
|---|
| 1066 | qscatrefVISgrid(j,k,iaer) = |
|---|
| 1067 | & qscatrefVISgrid(j,k,iaer) + |
|---|
| 1068 | & weightgaus(gausind) * |
|---|
| 1069 | & ( |
|---|
| 1070 | & omegrefVISb(gausind,iaer) * |
|---|
| 1071 | & qrefVISb(gausind,iaer) * |
|---|
| 1072 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 1073 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 1074 | & distb(j,k,iaer,idomain,gausind) + |
|---|
| 1075 | & omegrefVISa(gausind,iaer) * |
|---|
| 1076 | & qrefVISa(gausind,iaer) * |
|---|
| 1077 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 1078 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 1079 | & dista(j,k,iaer,idomain,gausind) |
|---|
| 1080 | & ) |
|---|
| 1081 | ENDDO |
|---|
| 1082 | qrefVISgrid(j,k,iaer)=qrefVISgrid(j,k,iaer) / |
|---|
| 1083 | & normd(j,k,iaer,idomain) |
|---|
| 1084 | qscatrefVISgrid(j,k,iaer)=qscatrefVISgrid(j,k,iaer) / |
|---|
| 1085 | & normd(j,k,iaer,idomain) |
|---|
| 1086 | omegrefVISgrid(j,k,iaer)=qscatrefVISgrid(j,k,iaer) / |
|---|
| 1087 | & qrefVISgrid(j,k,iaer) |
|---|
| 1088 | DO m=1,nsun |
|---|
| 1089 | qextVISgrid(j,k,m,iaer)=qextVISgrid(j,k,m,iaer) / |
|---|
| 1090 | & normd(j,k,iaer,idomain) |
|---|
| 1091 | qscatVISgrid(j,k,m,iaer)=qscatVISgrid(j,k,m,iaer) / |
|---|
| 1092 | & normd(j,k,iaer,idomain) |
|---|
| 1093 | gVISgrid(j,k,m,iaer)=gVISgrid(j,k,m,iaer) / |
|---|
| 1094 | & qscatVISgrid(j,k,m,iaer) / |
|---|
| 1095 | & normd(j,k,iaer,idomain) |
|---|
| 1096 | |
|---|
| 1097 | qsqrefVISgrid(j,k,m,iaer)=qextVISgrid(j,k,m,iaer) / |
|---|
| 1098 | & qrefVISgrid(j,k,iaer) |
|---|
| 1099 | omegVISgrid(j,k,m,iaer)=qscatVISgrid(j,k,m,iaer) / |
|---|
| 1100 | & qextVISgrid(j,k,m,iaer) |
|---|
| 1101 | ENDDO |
|---|
| 1102 | ELSE ! INFRARED DOMAIN ---------- |
|---|
| 1103 | c 2.3.3.ir Initialization |
|---|
| 1104 | qsqrefIRgrid(j,k,:,iaer)=0. |
|---|
| 1105 | qextIRgrid(j,k,:,iaer)=0. |
|---|
| 1106 | qscatIRgrid(j,k,:,iaer)=0. |
|---|
| 1107 | omegIRgrid(j,k,:,iaer)=0. |
|---|
| 1108 | gIRgrid(j,k,:,iaer)=0. |
|---|
| 1109 | qrefIRgrid(j,k,iaer)=0. |
|---|
| 1110 | qscatrefIRgrid(j,k,iaer)=0. |
|---|
| 1111 | omegrefIRgrid(j,k,iaer)=0. |
|---|
| 1112 | |
|---|
| 1113 | DO gausind=1,ngau |
|---|
| 1114 | DO m=1,nir |
|---|
| 1115 | c Convolution: |
|---|
| 1116 | qextIRgrid(j,k,m,iaer) = |
|---|
| 1117 | & qextIRgrid(j,k,m,iaer) + |
|---|
| 1118 | & weightgaus(gausind) * |
|---|
| 1119 | & ( |
|---|
| 1120 | & qsqrefIRb(m,gausind,iaer) * |
|---|
| 1121 | & qrefVISb(gausind,iaer) * |
|---|
| 1122 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 1123 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 1124 | & distb(j,k,iaer,idomain,gausind) + |
|---|
| 1125 | & qsqrefIRa(m,gausind,iaer) * |
|---|
| 1126 | & qrefVISa(gausind,iaer) * |
|---|
| 1127 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 1128 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 1129 | & dista(j,k,iaer,idomain,gausind) |
|---|
| 1130 | & ) |
|---|
| 1131 | qscatIRgrid(j,k,m,iaer) = |
|---|
| 1132 | & qscatIRgrid(j,k,m,iaer) + |
|---|
| 1133 | & weightgaus(gausind) * |
|---|
| 1134 | & ( |
|---|
| 1135 | & omegIRb(m,gausind,iaer) * |
|---|
| 1136 | & qsqrefIRb(m,gausind,iaer) * |
|---|
| 1137 | & qrefVISb(gausind,iaer) * |
|---|
| 1138 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 1139 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 1140 | & distb(j,k,iaer,idomain,gausind) + |
|---|
| 1141 | & omegIRa(m,gausind,iaer) * |
|---|
| 1142 | & qsqrefIRa(m,gausind,iaer) * |
|---|
| 1143 | & qrefVISa(gausind,iaer) * |
|---|
| 1144 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 1145 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 1146 | & dista(j,k,iaer,idomain,gausind) |
|---|
| 1147 | & ) |
|---|
| 1148 | gIRgrid(j,k,m,iaer) = |
|---|
| 1149 | & gIRgrid(j,k,m,iaer) + |
|---|
| 1150 | & weightgaus(gausind) * |
|---|
| 1151 | & ( |
|---|
| 1152 | & omegIRb(m,gausind,iaer) * |
|---|
| 1153 | & qsqrefIRb(m,gausind,iaer) * |
|---|
| 1154 | & qrefVISb(gausind,iaer) * |
|---|
| 1155 | & gIRb(m,gausind,iaer) * |
|---|
| 1156 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 1157 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 1158 | & distb(j,k,iaer,idomain,gausind) + |
|---|
| 1159 | & omegIRa(m,gausind,iaer) * |
|---|
| 1160 | & qsqrefIRa(m,gausind,iaer) * |
|---|
| 1161 | & qrefVISa(gausind,iaer) * |
|---|
| 1162 | & gIRa(m,gausind,iaer) * |
|---|
| 1163 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 1164 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 1165 | & dista(j,k,iaer,idomain,gausind) |
|---|
| 1166 | & ) |
|---|
| 1167 | ENDDO |
|---|
| 1168 | qrefIRgrid(j,k,iaer) = |
|---|
| 1169 | & qrefIRgrid(j,k,iaer) + |
|---|
| 1170 | & weightgaus(gausind) * |
|---|
| 1171 | & ( |
|---|
| 1172 | & qrefIRb(gausind,iaer) * |
|---|
| 1173 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 1174 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 1175 | & distb(j,k,iaer,idomain,gausind) + |
|---|
| 1176 | & qrefIRa(gausind,iaer) * |
|---|
| 1177 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 1178 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 1179 | & dista(j,k,iaer,idomain,gausind) |
|---|
| 1180 | & ) |
|---|
| 1181 | qscatrefIRgrid(j,k,iaer) = |
|---|
| 1182 | & qscatrefIRgrid(j,k,iaer) + |
|---|
| 1183 | & weightgaus(gausind) * |
|---|
| 1184 | & ( |
|---|
| 1185 | & omegrefIRb(gausind,iaer) * |
|---|
| 1186 | & qrefIRb(gausind,iaer) * |
|---|
| 1187 | & pi*radGAUSb(gausind,iaer,idomain) * |
|---|
| 1188 | & radGAUSb(gausind,iaer,idomain) * |
|---|
| 1189 | & distb(j,k,iaer,idomain,gausind) + |
|---|
| 1190 | & omegrefIRa(gausind,iaer) * |
|---|
| 1191 | & qrefIRa(gausind,iaer) * |
|---|
| 1192 | & pi*radGAUSa(gausind,iaer,idomain) * |
|---|
| 1193 | & radGAUSa(gausind,iaer,idomain) * |
|---|
| 1194 | & dista(j,k,iaer,idomain,gausind) |
|---|
| 1195 | & ) |
|---|
| 1196 | ENDDO |
|---|
| 1197 | qrefIRgrid(j,k,iaer)=qrefIRgrid(j,k,iaer) / |
|---|
| 1198 | & normd(j,k,iaer,idomain) |
|---|
| 1199 | qscatrefIRgrid(j,k,iaer)=qscatrefIRgrid(j,k,iaer) / |
|---|
| 1200 | & normd(j,k,iaer,idomain) |
|---|
| 1201 | omegrefIRgrid(j,k,iaer)=qscatrefIRgrid(j,k,iaer) / |
|---|
| 1202 | & qrefIRgrid(j,k,iaer) |
|---|
| 1203 | DO m=1,nir |
|---|
| 1204 | qextIRgrid(j,k,m,iaer)=qextIRgrid(j,k,m,iaer) / |
|---|
| 1205 | & normd(j,k,iaer,idomain) |
|---|
| 1206 | qscatIRgrid(j,k,m,iaer)=qscatIRgrid(j,k,m,iaer) / |
|---|
| 1207 | & normd(j,k,iaer,idomain) |
|---|
| 1208 | gIRgrid(j,k,m,iaer)=gIRgrid(j,k,m,iaer) / |
|---|
| 1209 | & qscatIRgrid(j,k,m,iaer) / |
|---|
| 1210 | & normd(j,k,iaer,idomain) |
|---|
| 1211 | |
|---|
| 1212 | qsqrefIRgrid(j,k,m,iaer)=qextIRgrid(j,k,m,iaer) / |
|---|
| 1213 | & qrefVISgrid(j,k,iaer) |
|---|
| 1214 | omegIRgrid(j,k,m,iaer)=qscatIRgrid(j,k,m,iaer) / |
|---|
| 1215 | & qextIRgrid(j,k,m,iaer) |
|---|
| 1216 | ENDDO |
|---|
| 1217 | ENDIF ! -------------------------- |
|---|
| 1218 | checkgrid(j,k,iaer,idomain) = .true. |
|---|
| 1219 | ENDIF !checkgrid |
|---|
| 1220 | ENDDO !grid_j |
|---|
| 1221 | ENDDO !grid_i |
|---|
| 1222 | c 3.4 Bilinear interpolation |
|---|
| 1223 | k1 = (1-kx)*(1-ky) |
|---|
| 1224 | k2 = kx*(1-ky) |
|---|
| 1225 | k3 = kx*ky |
|---|
| 1226 | k4 = (1-kx)*ky |
|---|
| 1227 | IF (idomain.EQ.1) THEN ! VISIBLE ------------------------ |
|---|
| 1228 | DO m=1,nsun |
|---|
| 1229 | QVISsQREF3d(ig,lg,m,iaer) = |
|---|
| 1230 | & k1*qsqrefVISgrid(grid_i,grid_j,m,iaer) + |
|---|
| 1231 | & k2*qsqrefVISgrid(grid_i+1,grid_j,m,iaer) + |
|---|
| 1232 | & k3*qsqrefVISgrid(grid_i+1,grid_j+1,m,iaer) + |
|---|
| 1233 | & k4*qsqrefVISgrid(grid_i,grid_j+1,m,iaer) |
|---|
| 1234 | omegaVIS3d(ig,lg,m,iaer) = |
|---|
| 1235 | & k1*omegVISgrid(grid_i,grid_j,m,iaer) + |
|---|
| 1236 | & k2*omegVISgrid(grid_i+1,grid_j,m,iaer) + |
|---|
| 1237 | & k3*omegVISgrid(grid_i+1,grid_j+1,m,iaer) + |
|---|
| 1238 | & k4*omegVISgrid(grid_i,grid_j+1,m,iaer) |
|---|
| 1239 | gVIS3d(ig,lg,m,iaer) = |
|---|
| 1240 | & k1*gVISgrid(grid_i,grid_j,m,iaer) + |
|---|
| 1241 | & k2*gVISgrid(grid_i+1,grid_j,m,iaer) + |
|---|
| 1242 | & k3*gVISgrid(grid_i+1,grid_j+1,m,iaer) + |
|---|
| 1243 | & k4*gVISgrid(grid_i,grid_j+1,m,iaer) |
|---|
| 1244 | ENDDO !nsun |
|---|
| 1245 | QREFvis3d(ig,lg,iaer) = |
|---|
| 1246 | & k1*qrefVISgrid(grid_i,grid_j,iaer) + |
|---|
| 1247 | & k2*qrefVISgrid(grid_i+1,grid_j,iaer) + |
|---|
| 1248 | & k3*qrefVISgrid(grid_i+1,grid_j+1,iaer) + |
|---|
| 1249 | & k4*qrefVISgrid(grid_i,grid_j+1,iaer) |
|---|
| 1250 | omegaREFvis3d(ig,lg,iaer) = |
|---|
| 1251 | & k1*omegrefVISgrid(grid_i,grid_j,iaer) + |
|---|
| 1252 | & k2*omegrefVISgrid(grid_i+1,grid_j,iaer) + |
|---|
| 1253 | & k3*omegrefVISgrid(grid_i+1,grid_j+1,iaer) + |
|---|
| 1254 | & k4*omegrefVISgrid(grid_i,grid_j+1,iaer) |
|---|
| 1255 | ELSE ! INFRARED ----------------------- |
|---|
| 1256 | DO m=1,nir |
|---|
| 1257 | QIRsQREF3d(ig,lg,m,iaer) = |
|---|
| 1258 | & k1*qsqrefIRgrid(grid_i,grid_j,m,iaer) + |
|---|
| 1259 | & k2*qsqrefIRgrid(grid_i+1,grid_j,m,iaer) + |
|---|
| 1260 | & k3*qsqrefIRgrid(grid_i+1,grid_j+1,m,iaer) + |
|---|
| 1261 | & k4*qsqrefIRgrid(grid_i,grid_j+1,m,iaer) |
|---|
| 1262 | omegaIR3d(ig,lg,m,iaer) = |
|---|
| 1263 | & k1*omegIRgrid(grid_i,grid_j,m,iaer) + |
|---|
| 1264 | & k2*omegIRgrid(grid_i+1,grid_j,m,iaer) + |
|---|
| 1265 | & k3*omegIRgrid(grid_i+1,grid_j+1,m,iaer) + |
|---|
| 1266 | & k4*omegIRgrid(grid_i,grid_j+1,m,iaer) |
|---|
| 1267 | gIR3d(ig,lg,m,iaer) = |
|---|
| 1268 | & k1*gIRgrid(grid_i,grid_j,m,iaer) + |
|---|
| 1269 | & k2*gIRgrid(grid_i+1,grid_j,m,iaer) + |
|---|
| 1270 | & k3*gIRgrid(grid_i+1,grid_j+1,m,iaer) + |
|---|
| 1271 | & k4*gIRgrid(grid_i,grid_j+1,m,iaer) |
|---|
| 1272 | ENDDO !nir |
|---|
| 1273 | QREFir3d(ig,lg,iaer) = |
|---|
| 1274 | & k1*qrefIRgrid(grid_i,grid_j,iaer) + |
|---|
| 1275 | & k2*qrefIRgrid(grid_i+1,grid_j,iaer) + |
|---|
| 1276 | & k3*qrefIRgrid(grid_i+1,grid_j+1,iaer) + |
|---|
| 1277 | & k4*qrefIRgrid(grid_i,grid_j+1,iaer) |
|---|
| 1278 | omegaREFir3d(ig,lg,iaer) = |
|---|
| 1279 | & k1*omegrefIRgrid(grid_i,grid_j,iaer) + |
|---|
| 1280 | & k2*omegrefIRgrid(grid_i+1,grid_j,iaer) + |
|---|
| 1281 | & k3*omegrefIRgrid(grid_i+1,grid_j+1,iaer) + |
|---|
| 1282 | & k4*omegrefIRgrid(grid_i,grid_j+1,iaer) |
|---|
| 1283 | ENDIF ! -------------------------------- |
|---|
| 1284 | ENDDO !nlayer |
|---|
| 1285 | ENDDO !ngrid |
|---|
| 1286 | |
|---|
| 1287 | ENDIF ! varyingnueff |
|---|
| 1288 | c================================================================== |
|---|
| 1289 | ENDDO ! idomain |
|---|
| 1290 | |
|---|
| 1291 | ENDIF ! nsize = 1 |
|---|
| 1292 | |
|---|
| 1293 | ENDDO ! iaer (loop on aerosol kind) |
|---|
| 1294 | |
|---|
| 1295 | c=====Radiative properties - TESTS================================= |
|---|
| 1296 | IF (out_qwg) THEN |
|---|
| 1297 | c ------------------------------------------------------------- |
|---|
| 1298 | DO out_nchannel = 1, 2 |
|---|
| 1299 | c ------------------------------------------------------------- |
|---|
| 1300 | DO lg = 1, nlayer |
|---|
| 1301 | DO ig = 1, ngrid |
|---|
| 1302 | out_qext(ig,lg) = |
|---|
| 1303 | & QVISsQREF3d(ig,lg,out_nchannel,out_iaer)* |
|---|
| 1304 | & QREFvis3d(ig,lg,out_iaer) |
|---|
| 1305 | out_omeg(ig,lg) = |
|---|
| 1306 | & omegaVIS3d(ig,lg,out_nchannel,out_iaer) |
|---|
| 1307 | out_g(ig,lg) = gVIS3d(ig,lg,out_nchannel,out_iaer) |
|---|
| 1308 | ENDDO ! ig |
|---|
| 1309 | ENDDO ! lg |
|---|
| 1310 | #ifndef MESOSCALE |
|---|
| 1311 | write(out_str(1:1),'(i1.1)') out_nchannel |
|---|
| 1312 | call write_output('qextvis'//out_str,"Ext.efficiency","" |
|---|
| 1313 | & ,out_qext(:,:)) |
|---|
| 1314 | call write_output('omegvis'//out_str,"Sing.Scat.Alb.","", |
|---|
| 1315 | & out_omeg(:,:)) |
|---|
| 1316 | call write_output('gvis'//out_str,"Asym.Factor","", |
|---|
| 1317 | & out_g(:,:)) |
|---|
| 1318 | #endif |
|---|
| 1319 | c ------------------------------------------------------------- |
|---|
| 1320 | ENDDO ! out_nchannel |
|---|
| 1321 | DO out_nchannel = 2, 4 |
|---|
| 1322 | c ------------------------------------------------------------- |
|---|
| 1323 | DO lg = 1, nlayer |
|---|
| 1324 | DO ig = 1, ngrid |
|---|
| 1325 | out_qext(ig,lg) = |
|---|
| 1326 | & QIRsQREF3d(ig,lg,out_nchannel,out_iaer)* |
|---|
| 1327 | & QREFir3d(ig,lg,out_iaer) |
|---|
| 1328 | out_omeg(ig,lg) = |
|---|
| 1329 | & omegaIR3d(ig,lg,out_nchannel,out_iaer) |
|---|
| 1330 | out_g(ig,lg) = gIR3d(ig,lg,out_nchannel,out_iaer) |
|---|
| 1331 | ENDDO ! ig |
|---|
| 1332 | ENDDO ! lg |
|---|
| 1333 | #ifndef MESOSCALE |
|---|
| 1334 | write(out_str(1:1),'(i1.1)') out_nchannel |
|---|
| 1335 | call write_output('qextir'//out_str,"Ext.efficiency","", |
|---|
| 1336 | & out_qext(:,:)) |
|---|
| 1337 | call write_output('omegir'//out_str,"Sing.Scat.Alb.","", |
|---|
| 1338 | & out_omeg(:,:)) |
|---|
| 1339 | call write_output('gir'//out_str,"Asym.Factor","", |
|---|
| 1340 | & out_g(:,:)) |
|---|
| 1341 | #endif |
|---|
| 1342 | c ------------------------------------------------------------- |
|---|
| 1343 | ENDDO ! out_nchannel |
|---|
| 1344 | #ifndef MESOSCALE |
|---|
| 1345 | call write_output("omegvisref","Sing.Scat.Alb.","", |
|---|
| 1346 | & omegaREFvis3d(:,:,out_iaer)) |
|---|
| 1347 | call write_output("omegirref","Sing.Scat.Alb.","", |
|---|
| 1348 | & omegaREFir3d(:,:,out_iaer)) |
|---|
| 1349 | #endif |
|---|
| 1350 | ENDIF ! out_qwg |
|---|
| 1351 | c================================================================== |
|---|
| 1352 | |
|---|
| 1353 | END SUBROUTINE aeroptproperties |
|---|
| 1354 | |
|---|
| 1355 | END MODULE aeroptproperties_mod |
|---|
