module aeropacity_mod implicit none contains Subroutine aeropacity(ngrid,nlayer,nq,pplay,pplev,zzlev,pt,pq, & dtau_aer,reffrad,nueffrad, QREFvis3d,QREFir3d,tau_col) use radinc_h, only : L_TAUMAX,naerkind use aerosol_mod, only: iaero_haze, i_haze USE tracer_h, only: noms,rho_n2,rho_ice,rho_q,mmol,micro_indx use comcstfi_mod, only: g, pi, mugaz, avocado use geometry_mod, only: latitude use callkeys_mod, only: kastprof, callmufi use mp2m_diagnostics implicit none !================================================================== ! ! Purpose ! ------- ! Compute aerosol optical depth in each gridbox. ! ! Authors ! ------- ! F. Forget ! F. Montmessin (water ice scheme) ! update J.-B. Madeleine (2008) ! dust removal, simplification by Robin Wordsworth (2009) ! ! Input ! ----- ! ngrid Number of horizontal gridpoints ! nlayer Number of layers ! nq Number of tracers ! pplev Pressure (Pa) at each layer boundary ! pq Aerosol mixing ratio ! reffrad(ngrid,nlayer,naerkind) Aerosol effective radius ! QREFvis3d(ngrid,nlayer,naerkind) \ 3d extinction coefficients ! QREFir3d(ngrid,nlayer,naerkind) / at reference wavelengths ! ! Output ! ------ ! dtau_aer Aerosol optical depth in layer l, grid point ig ! tau_col Total column optical depth at grid point ig ! !======================================================================= INTEGER,INTENT(IN) :: ngrid ! number of atmospheric columns INTEGER,INTENT(IN) :: nlayer ! number of atmospheric layers INTEGER,INTENT(IN) :: nq ! number of tracers REAL,INTENT(IN) :: pplay(ngrid,nlayer) ! mid-layer pressure (Pa) REAL,INTENT(IN) :: pplev(ngrid,nlayer+1) ! inter-layer pressure (Pa) REAL,INTENT(IN) :: zzlev(ngrid,nlayer) ! Altitude at the layer boundaries. REAL,INTENT(IN) :: pq(ngrid,nlayer,nq) ! tracers (.../kg_of_air) REAL,INTENT(IN) :: pt(ngrid,nlayer) ! mid-layer temperature (K) REAL,INTENT(OUT) :: dtau_aer(ngrid,nlayer,naerkind) ! aerosol optical depth REAL,INTENT(IN) :: reffrad(ngrid,nlayer,naerkind) ! aerosol effective radius REAL,INTENT(IN) :: nueffrad(ngrid,nlayer,naerkind) ! aerosol effective variance REAL,INTENT(IN) :: QREFvis3d(ngrid,nlayer,naerkind) ! extinction coefficient in the visible REAL,INTENT(IN) :: QREFir3d(ngrid,nlayer,naerkind) REAL,INTENT(OUT):: tau_col(ngrid) !column integrated visible optical depth real aerosol0, obs_tau_col_aurora, pm real pcloud_deck, cloud_slope real dp_strato(ngrid) real dp_tropo(ngrid) real dp_layer(ngrid) ! Microphysical tracers real sig real m0as(ngrid,nlayer) real m0af(ngrid,nlayer) INTEGER l,ig,iq,iaer,ia LOGICAL,SAVE :: firstcall=.true. !$OMP THREADPRIVATE(firstcall) REAL CBRT EXTERNAL CBRT INTEGER,SAVE :: i_n2ice=0 ! n2 ice !$OMP THREADPRIVATE(i_n2ice) CHARACTER(LEN=20) :: tracername ! to temporarily store text real CLFtot integer igen_ice,igen_gas ! to store the index of generic tracer logical dummy_bool ! dummy boolean just in case we need one ! for venus clouds real :: p_bot,p_top,h_bot,h_top,mode_dens,h_lay ! identify tracers IF (firstcall) THEN ia =0 write(*,*) "Tracers found in aeropacity:" do iq=1,nq tracername=noms(iq) if (tracername.eq."n2_ice") then i_n2ice=iq write(*,*) "i_n2ice=",i_n2ice endif enddo firstcall=.false. ENDIF ! of IF (firstcall) ! --------------------------------------------------------- !================================================================== ! Haze aerosols !================================================================== if (iaero_haze.ne.0) then if (callmufi) then ! Convert intensive microphysical tracers to extensive [m-2] m0as(:,:) = pq(:,:,micro_indx(1)) * (pplev(:,1:nlayer) - pplev(:,2:nlayer+1)) / g m0af(:,:) = pq(:,:,micro_indx(3)) * (pplev(:,1:nlayer) - pplev(:,2:nlayer+1)) / g ! Spherical aerosols sig = 0.2 dtau_aer(:,:,1) = m0as(:,:) * QREFvis3d(:,:,1) * pi * mp2m_rc_sph(:,:)**2 * exp(2*sig**2) ! Fractal aerosols sig = 0.35 dtau_aer(:,:,2) = m0af(:,:) * QREFvis3d(:,:,2) * pi * mp2m_rc_fra(:,:)**2 * exp(2*sig**2) ! write(*,*) 'dtau_as :', MINVAL(dtau_aer(:,:,1)), '-', MAXVAL(dtau_aer(:,:,1)) ! write(*,*) 'dtau_af :', MINVAL(dtau_aer(:,:,2)), '-', MAXVAL(dtau_aer(:,:,2)) else do iaer = 1, naerkind ! 1. Initialization dtau_aer(1:ngrid,1:nlayer,iaer)=0.0 ! 2. Opacity calculation DO ig = 1, ngrid DO l = 1, nlayer-1 ! to stop the rad tran bug ! If fractal, radius doit etre equivalent sphere radius ! Eq. 2.37 - Madeleine's PhD (2011). aerosol0 = & ( 0.75 * QREFvis3d(ig,l,iaer) / & ( rho_q(i_haze) * reffrad(ig,l,iaer) ) ) * & ( pq(ig,l,i_haze) + 1.E-10 ) * & ( pplev(ig,l) - pplev(ig,l+1) ) / g aerosol0 = max(aerosol0,1.e-10) aerosol0 = min(aerosol0,L_TAUMAX) dtau_aer(ig,l,iaer) = aerosol0 ENDDO ENDDO !QREF est le meme dans toute la colonne par def si size uniforme !print*, 'TB17: QREFvis3d=',QREFvis3d(1,:,1) !print*, 'TB17: rho_q=',rho_q(i_haze) !print*, 'TB17: reffrad=',reffrad(1,:,1) !print*, 'TB17: pq=',pq(1,:,i_haze) !print*, 'TB17: deltap=',pplev(1,1) - pplev(1,nlayer) enddo ! end iaer = 1, naerkind endif ! end callmufi endif ! if haze aerosols ! -------------------------------------------------------------------------- ! Column integrated visible optical depth in each point (used for diagnostic) tau_col(:)=0.0 do iaer = 1, naerkind do l=1,nlayer do ig=1,ngrid tau_col(ig) = tau_col(ig) + dtau_aer(ig,l,iaer) end do end do end do do ig=1,ngrid do l=1,nlayer do iaer = 1, naerkind if(dtau_aer(ig,l,iaer).gt.1.e3)then print*,'WARNING: dtau_aer=',dtau_aer(ig,l,iaer) print*,'at ig=',ig,', l=',l,', iaer=',iaer print*,'QREFvis3d=',QREFvis3d(ig,l,iaer) print*,'reffrad=',reffrad(ig,l,iaer) endif end do end do end do do ig=1,ngrid if(tau_col(ig).gt.1.e3)then print*,'WARNING: tau_col=',tau_col(ig) print*,'at ig=',ig print*,'dtau_aer=',dtau_aer(ig,:,:) print*,'QREFvis3d=',QREFvis3d(ig,:,:) print*,'reffrad=',reffrad(ig,:,:) endif end do end subroutine aeropacity end module aeropacity_mod