subroutine aeroptproperties(ngrid,nlayer,reffrad,nueffrad, & QVISsQREF3d,omegaVIS3d,gVIS3d, & QIRsQREF3d,omegaIR3d,gIR3d, & QREFvis3d,QREFir3d) use radinc_h, only: L_NSPECTI,L_NSPECTV,naerkind,nsizemax use radcommon_h, only: QVISsQREF,omegavis,gvis,QIRsQREF,omegair,gir use radcommon_h, only: radiustab,nsize,qrefvis,qrefir implicit none !================================================================== ! ! Purpose ! ------- ! Compute the scattering parameters in each grid ! box, depending on aerosol grain sizes. ! ! Notes ! ----- ! Don't forget to set the value of varyingnueff below; If ! the effective variance of the distribution for the given ! aerosol is considered homogeneous in the atmosphere, please ! set varyingnueff(iaer) to .false. Resulting computational ! time will be much better. ! ! Authors ! ------- ! J.-B. Madeleine, F. Montmessin ! Slightly modified and converted to F90 by R. Wordsworth (2009) ! Varying nueff section removed by R. Wordsworth for simplicity ! !================================================================== #include "dimensions.h" #include "dimphys.h" #include "callkeys.h" ! Local variables ! --------------- ! ============================================================= LOGICAL, PARAMETER :: varyingnueff(naerkind) = .false. ! ============================================================= ! Radius axis used for integration DOUBLE PRECISION :: radiusint(nsizemax+1) ! Min. and max radii of the interpolation grid (in METERS) REAL, PARAMETER :: refftabmin = 2e-8 ! REAL, PARAMETER :: refftabmax = 35e-6 REAL, PARAMETER :: refftabmax = 1e-3 ! CHANGED BY RDW ! Log of the min and max variance of the interpolation grid REAL, PARAMETER :: nuefftabmin = -4.6 REAL, PARAMETER :: nuefftabmax = 0. ! Number of effective radii of the interpolation grid INTEGER, PARAMETER :: refftabsize = 200 ! Number of effective variances of the interpolation grid ! INTEGER, PARAMETER :: nuefftabsize = 100 INTEGER, PARAMETER :: nuefftabsize = 1 ! CHANGED BY RDW ! Interpolation grid indices (reff,nueff) INTEGER :: grid_i,grid_j ! Volume ratio of the look-up table (different in VIS and IR) DOUBLE PRECISION :: vrat ! r_g and sigma_g for the log-normal distribution ! as defined by [hansen_1974] REAL :: r_g,sigma_g ! Error function used for integration DOUBLE PRECISION :: derf ! Density function f(x)dx of the log-normal distribution REAL :: dfi DOUBLE PRECISION :: dfi_tmp(nsizemax+1) ! Intermediate variable REAL :: var_tmp,var3d_tmp(ngridmx,nlayermx) ! Bilinear interpolation factors REAL :: kx,ky,k1,k2,k3,k4 ! Indices INTEGER :: i,j,k,l,m,iaer,idomain INTEGER :: ig,lg,chg ! Local saved variables ! --------------------- ! Radius axis of the interpolation grid DOUBLE PRECISION,SAVE :: refftab(refftabsize) ! Variance axis of the interpolation grid DOUBLE PRECISION,SAVE :: nuefftab(nuefftabsize) ! Volume ratio of the grid DOUBLE PRECISION,SAVE :: logvratgrid,vratgrid ! Grid used to remember which calculation is done LOGICAL,SAVE :: checkgrid(refftabsize,nuefftabsize,naerkind,2) & = .false. ! Optical properties of the grid (VISIBLE) REAL,SAVE :: epVISgrid(refftabsize,nuefftabsize,L_NSPECTV,naerkind) REAL,SAVE :: omegVISgrid(refftabsize,nuefftabsize,L_NSPECTV,naerkind) REAL,SAVE :: gVISgrid(refftabsize,nuefftabsize,L_NSPECTV,naerkind) ! Optical properties of the grid (INFRARED) REAL,SAVE :: epIRgrid(refftabsize,nuefftabsize,L_NSPECTI,naerkind) REAL,SAVE :: omegIRgrid(refftabsize,nuefftabsize,L_NSPECTI,naerkind) REAL,SAVE :: gIRgrid(refftabsize,nuefftabsize,L_NSPECTI,naerkind) ! Optical properties of the grid (REFERENCE WAVELENGTHS) REAL,SAVE :: eprefVISgrid(refftabsize,nuefftabsize,naerkind) REAL,SAVE :: eprefIRgrid(refftabsize,nuefftabsize,naerkind) ! Firstcall LOGICAL,SAVE :: firstcall = .true. ! Inputs ! ------ INTEGER :: ngrid,nlayer ! Aerosol effective radius used for radiative transfer (meter) REAL :: reffrad(ngridmx,nlayermx,naerkind) ! Aerosol effective variance used for radiative transfer (n.u.) REAL :: nueffrad(ngridmx,nlayermx,naerkind) ! Outputs ! ------- REAL :: QVISsQREF3d(ngridmx,nlayermx,L_NSPECTV,naerkind) REAL :: omegaVIS3d(ngridmx,nlayermx,L_NSPECTV,naerkind) REAL :: gVIS3d(ngridmx,nlayermx,L_NSPECTV,naerkind) REAL :: QIRsQREF3d(ngridmx,nlayermx,L_NSPECTI,naerkind) REAL :: omegaIR3d(ngridmx,nlayermx,L_NSPECTI,naerkind) REAL :: gIR3d(ngridmx,nlayermx,L_NSPECTI,naerkind) REAL :: QREFvis3d(ngridmx,nlayermx,naerkind) REAL :: QREFir3d(ngridmx,nlayermx,naerkind) DO iaer = 1, naerkind ! Loop on aerosol kind IF ( (nsize(iaer,1).EQ.1).AND.(nsize(iaer,2).EQ.1) ) THEN !================================================================== ! If there is one single particle size, optical ! properties of the considered aerosol are homogeneous DO lg = 1, nlayer DO ig = 1, ngrid DO chg = 1, L_NSPECTV QVISsQREF3d(ig,lg,chg,iaer)=QVISsQREF(chg,iaer,1) omegaVIS3d(ig,lg,chg,iaer)=omegaVIS(chg,iaer,1) gVIS3d(ig,lg,chg,iaer)=gVIS(chg,iaer,1) ENDDO DO chg = 1, L_NSPECTI QIRsQREF3d(ig,lg,chg,iaer)=QIRsQREF(chg,iaer,1) omegaIR3d(ig,lg,chg,iaer)=omegaIR(chg,iaer,1) gIR3d(ig,lg,chg,iaer)=gIR(chg,iaer,1) ENDDO QREFvis3d(ig,lg,iaer)=QREFvis(iaer,1) QREFir3d(ig,lg,iaer)=QREFir(iaer,1) ENDDO ENDDO if (firstcall) then print*,'Optical properties of the aerosol are homogenous for:' print*,'iaer = ',iaer endif ELSE ! Varying effective radius and variance DO idomain = 1, 2 ! Loop on visible or infrared channel !================================================================== ! 1. Creating the effective radius and variance grid ! -------------------------------------------------- IF (firstcall) THEN ! 1.1 Effective radius refftab(1) = refftabmin refftab(refftabsize) = refftabmax logvratgrid = log(refftabmax/refftabmin) / & float(refftabsize-1)*3. vratgrid = exp(logvratgrid) do i = 2, refftabsize-1 refftab(i) = refftab(i-1)*vratgrid**(1./3.) enddo ! 1.2 Effective variance do i = 0, nuefftabsize-1 nuefftab(i+1) = exp( nuefftabmin + & i*(nuefftabmax-nuefftabmin)/(nuefftabsize-1) ) enddo firstcall = .false. ENDIF ! 2. Compute the scattering parameters using linear ! interpolation over grain sizes and constant nueff ! --------------------------------------------------- ! 2.1 Initialization vrat = log(radiustab(iaer,idomain,nsize(iaer,idomain)) / & radiustab(iaer,idomain,1)) / & float(nsize(iaer,idomain)-1)*3. vrat = exp(vrat) radiusint(1) = 1.e-9 DO i = 2,nsize(iaer,idomain) radiusint(i) = ( (2.*vrat) / (vrat+1.) )**(1./3.) * & radiustab(iaer,idomain,i-1) ENDDO radiusint(nsize(iaer,idomain)+1) = 1.e-2 DO lg = 1,nlayer DO ig = 1,ngrid ! 2.1 Effective radius index and kx calculation var_tmp=reffrad(ig,lg,iaer)/refftabmin var_tmp=log(var_tmp)*3. var_tmp=var_tmp/logvratgrid+1. grid_i=floor(var_tmp) IF (grid_i.GE.refftabsize) THEN WRITE(*,*) 'Warning: Aerosol particle size in grid box #' WRITE(*,*) ig,' is too large to be used by the ' WRITE(*,*) 'radiative transfer; please extend the ' WRITE(*,*) 'interpolation grid to larger sizes.' grid_i=refftabsize-1 kx = 1. ELSEIF (grid_i.LT.1) THEN WRITE(*,*) 'Warning: Aerosol particle size in grid box #' WRITE(*,*) ig,' is too small to be used by the ' WRITE(*,*) 'radiative transfer; please extend the ' WRITE(*,*) 'interpolation grid to smaller sizes.' grid_i=1 kx = 0. ELSE kx = ( reffrad(ig,lg,iaer)-refftab(grid_i) ) / & ( refftab(grid_i+1)-refftab(grid_i) ) ENDIF ! 2.3 Integration DO j=grid_i,grid_i+1 ! 2.3.1 Check if the calculation has been completed IF (.NOT.checkgrid(j,1,iaer,idomain)) THEN ! 2.3.2 Compute r_g and sigma_g for the log-normal ! distribution as defined by [hansen_1974], "Light ! scattering in planetary atmospheres", Space ! Science Reviews 16 527-610, p558 sigma_g = log(1.+nueffrad(1,1,iaer)) r_g = exp(2.5*sigma_g) sigma_g = sqrt(sigma_g) r_g = refftab(j) / r_g IF (idomain.EQ.1) THEN ! VISIBLE DOMAIN ----------- ! 2.3.3.vis Initialization epVISgrid(j,1,:,iaer)=0. omegVISgrid(j,1,:,iaer)=0. gVISgrid(j,1,:,iaer)=0. eprefVISgrid(j,1,iaer)=0. ! 2.3.4.vis Log-normal distribution DO l=1,nsize(iaer,idomain)+1 dfi_tmp(l) = log(radiusint(l)/r_g) / & sqrt(2.)/sigma_g ENDDO DO l=1,nsize(iaer,idomain) dfi = 0.5*( derf(dfi_tmp(l+1)) - & derf(dfi_tmp(l)) ) DO m=1,L_NSPECTV epVISgrid(j,1,m,iaer) = & epVISgrid(j,1,m,iaer) & + QVISsQREF(m,iaer,l)*dfi omegVISgrid(j,1,m,iaer) = & omegVISgrid(j,1,m,iaer) & + omegaVIS(m,iaer,l)*dfi gVISgrid(j,1,m,iaer) = & gVISgrid(j,1,m,iaer) & + gVIS(m,iaer,l)*dfi ENDDO !L_NSPECTV eprefVISgrid(j,1,iaer) = & eprefVISgrid(j,1,iaer) & + QREFvis(iaer,l)*dfi ENDDO !nsize ELSE ! INFRARED DOMAIN ---------- ! 2.3.3.ir Initialization epIRgrid(j,1,:,iaer)=0. omegIRgrid(j,1,:,iaer)=0. gIRgrid(j,1,:,iaer)=0. eprefIRgrid(j,1,iaer)=0. ! 2.3.4.ir Log-normal distribution DO l=1,nsize(iaer,idomain)+1 dfi_tmp(l) = log(radiusint(l)/r_g) / & sqrt(2.)/sigma_g ENDDO DO l=1,nsize(iaer,idomain) dfi = 0.5*( derf(dfi_tmp(l+1)) - & derf(dfi_tmp(l)) ) DO m=1,L_NSPECTI epIRgrid(j,1,m,iaer) = & epIRgrid(j,1,m,iaer) & + QIRsQREF(m,iaer,l)*dfi omegIRgrid(j,1,m,iaer) = & omegIRgrid(j,1,m,iaer) & + omegaIR(m,iaer,l)*dfi gIRgrid(j,1,m,iaer) = & gIRgrid(j,1,m,iaer) & + gIR(m,iaer,l)*dfi ENDDO !L_NSPECTI eprefIRgrid(j,1,iaer) = & eprefIRgrid(j,1,iaer) & + QREFir(iaer,l)*dfi ENDDO !nsize ENDIF ! -------------------------- checkgrid(j,1,iaer,idomain) = .true. ENDIF !checkgrid ENDDO !grid_i ! 2.4 Linear interpolation k1 = (1-kx) k2 = kx IF (idomain.EQ.1) THEN ! VISIBLE ------------------------ DO m=1,L_NSPECTV QVISsQREF3d(ig,lg,m,iaer) = & k1*epVISgrid(grid_i,1,m,iaer) + & k2*epVISgrid(grid_i+1,1,m,iaer) omegaVIS3d(ig,lg,m,iaer) = & k1*omegVISgrid(grid_i,1,m,iaer) + & k2*omegVISgrid(grid_i+1,1,m,iaer) gVIS3d(ig,lg,m,iaer) = & k1*gVISgrid(grid_i,1,m,iaer) + & k2*gVISgrid(grid_i+1,1,m,iaer) ENDDO !L_NSPECTV QREFvis3d(ig,lg,iaer) = & k1*eprefVISgrid(grid_i,1,iaer) + & k2*eprefVISgrid(grid_i+1,1,iaer) ELSE ! INFRARED ----------------------- DO m=1,L_NSPECTI QIRsQREF3d(ig,lg,m,iaer) = & k1*epIRgrid(grid_i,1,m,iaer) + & k2*epIRgrid(grid_i+1,1,m,iaer) omegaIR3d(ig,lg,m,iaer) = & k1*omegIRgrid(grid_i,1,m,iaer) + & k2*omegIRgrid(grid_i+1,1,m,iaer) gIR3d(ig,lg,m,iaer) = & k1*gIRgrid(grid_i,1,m,iaer) + & k2*gIRgrid(grid_i+1,1,m,iaer) ENDDO !L_NSPECTI QREFir3d(ig,lg,iaer) = & k1*eprefIRgrid(grid_i,1,iaer) + & k2*eprefIRgrid(grid_i+1,1,iaer) ENDIF ! -------------------------------- ENDDO !nlayermx ENDDO !ngridmx ENDDO ! idomain ENDIF ! nsize = 1 ENDDO ! iaer (loop on aerosol kind) ! open(116,file='QIRsQREF3dO.dat') ! write(116,*) QIRsQREF3d ! close(116) ! open(117,file='omegaIR3dO.dat') ! write(117,*) omegaIR3d ! close(117) ! open(118,file='gIR3dO.dat') ! write(118,*) gIR3d ! close(118) ! stop RETURN END subroutine aeroptproperties