subroutine hazecloud(ngrid,nlayer,nq,ptimestep, & pplay,pplev,pq,pdq,pdist_sol,mu0,pfluxuv,zdqhaze, & zdqphot_prec,zdqphot_ch4,zdqconv_prec,declin) ! zdqphot_ch4,zdqconv_prec,declin,zdqhaze_col) use radinc_h, only : naerkind use comgeomfi_h implicit none !================================================================== ! Purpose ! ------- ! Production of haze in the atmosphere ! ! Inputs ! ------ ! ngrid Number of vertical columns ! nlayer Number of layers ! pplay(ngrid,nlayer) Pressure layers ! pplev(ngrid,nlayer+1) Pressure levels ! ! Outputs ! ------- ! ! Both ! ---- ! ! Authors ! ------- ! Tanguy Bertrand and Francois Forget (2014) ! !================================================================== #include "dimensions.h" #include "dimphys.h" #include "comcstfi.h" #include "surfdat.h" #include "comvert.h" #include "callkeys.h" #include "tracer.h" !----------------------------------------------------------------------- ! Arguments INTEGER ngrid, nlayer, nq ! REAL lati(ngrid) ! REAL long(ngrid) ! REAL declin LOGICAL firstcall SAVE firstcall DATA firstcall/.true./ REAL ptimestep REAL pplay(ngrid,nlayer),pplev(ngrid,nlayer+1) REAL,INTENT(IN) :: pq(ngrid,nlayer,nq) REAL,INTENT(IN) :: pdq(ngrid,nlayer,nq) ! REAL,INTENT(IN) :: mmol(nq) REAL,INTENT(IN) :: pdist_sol ! distance SUN-pluto in AU REAL,INTENT(IN) :: pfluxuv ! Lyman alpha flux at specific Ls (ph/cm/s) REAL,INTENT(IN) :: mu0(ngrid) ! cosinus of solar incident flux REAL,INTENT(IN) :: declin ! distance SUN-pluto in AU REAL,INTENT(OUT) :: zdqhaze(ngrid,nlayer,nq) ! Final tendancy ! REAL,INTENT(OUT) :: zdqhaze_col(ngrid) ! Final tendancy haze prod !----------------------------------------------------------------------- ! Local variables INTEGER l,ig,iq REAL zq_ch4(ngrid,nlayer) REAL zq_prec(ngrid,nlayer) REAL tauch4(nlayer) REAL sigch4 ! aborption cross section ch4 cm-2 mol-1 REAL flym_sol_earth ! initial flux Earth ph.m-2.s-1 REAL flym_sol_pluto ! initial Lyman alpha SOLAR flux Pluto ph.m-2.s-1 REAL flym_ipm(ngrid) ! top of atm. Incident IPM (Interplanetary) flux Lyman alpha ph.m-2.s-1 REAL fluxlym_sol(nlayer+1) ! Local Solar flux Lyman alpha ph.m-2.s-1 REAL fluxlym_ipm(nlayer+1) ! Local IPM (Interplanetary) flux Lyman alpha ph.m-2.s-1 REAL mu_ipm(ngrid) ! local Mean incident flux for IPM Lyman alpha photons REAL mu_sol(ngrid) ! local Mean incident flux for solar lyman alpha photons REAL gradflux(nlayer) ! gradient flux Lyman alpha ph.m-2.s-1 REAL zdqphot_ch4(ngrid,nlayer) ! tendancy due to photolysis of ch4 REAL zdqphot_prec(ngrid,nlayer) ! tendancy due to photolysis of ch4 REAL zdqconv_prec(ngrid,nlayer) ! tendancy due to conversion of ! prec_haze into haze REAL kch4 ! fraction of Carbon from ch4 directly dissociated ! by prec_haze CHARACTER(len=10) :: tname REAL ncratio ! ration Nitrogen/Carbon observed in tholins REAL tcon ! Time constant: conversion in aerosol REAL avogadro ! avogadro constant REAL longit REAL latit REAL valmin REAL valmax REAL valmin_dl REAL puis REAL dist REAL long2 !----------------------------------------------------------------------- !---------------- INPUT ------------------------------------------------ avogadro = 6.022141e23 sigch4 = 1.85e-17 ! aborption cross section ch4 cm-2 mol-1 !! Initial Solar flux Lyman alpha on Earth (1AU) flym_sol_earth=pfluxuv*1.e15 ! ph.m-2.s-1 -> 5e+11 ph.cm-2.s-1 !! Initial Solar flux Lyman alpha on Pluto flym_sol_pluto=flym_sol_earth/(pdist_sol*pdist_sol) ! ph.m-2.s-1 ! option decrease by 12% the initial IPM flux to account for Interplanetary H absorption: ! Fig 3 Gladstone et al 2014 : Lyalpha at Pluto flym_sol_pluto=flym_sol_pluto*0.878 !!!! Top of atm. Incident IPM (Interplanetary) flux Lyman alpha ph.m-2.s-1 ! fit Fig. 4 in Randall Gladstone et al. (2015) ! -- New version : Integration over semi sphere of Gladstone data. ! Fit of results IF (ngrid.eq.1) THEN flym_ipm(1)=75.e10 mu_ipm(1) = 0.5 !max(mu0(1), 0.5) mu_sol(1)=0.25 ELSE ! 1) get longitude/latitude (rad) of anti-subsolar point (max de mu0 - 180) longit=long((MAXLOC(mu0,DIM=1,MASK=mu0.GT.0.9))) ! rad IF (longit.GE.0) THEN longit=longit-pi ELSE longit=longit+pi ENDIF latit=-declin ! rad ! 2) Define input parameter for the fit valmin=48.74e10 ! minimum value of ipm flux in ph/m2/s valmax=115.014e10 ! maximal value of ipm flux in ph/m2/s valmin_dl=74.5e10 ! daylight minimum value puis=3.5 ! 3) Loop for each location and fit DO ig=1,ngrid ! calculation of cosinus of incident angle for IPM flux mu_sol(ig) = mu0(ig) mu_ipm(ig) = max(mu_sol(ig), 0.5) IF (mu0(ig).LT.1.e-4) THEN dist=acos(sin(lati(ig))*sin(latit)+cos(lati(ig))* & cos(latit)*cos(longit-long(ig)))*180/pi IF (dist.gt.90) dist=90 flym_ipm(ig)=(valmin_dl-valmin)/(90.**puis)*(dist)**puis & +valmin ELSE flym_ipm(ig)= mu0(ig)*(valmax-valmin_dl)+valmin_dl ENDIF ! proportional to lyman alpha solar flux (reference 2015 : 5e11 ph/cm2/s) flym_ipm(ig)=flym_ipm(ig)*pfluxuv/5. ENDDO ! nlayer ENDIF ! ngrid !--- !! Time constant of conversion in aerosol to be explore tcon= 1.e7 ! 153 ! 1.E7 !second ! tcon= 1. ! 153 ! 1.E7 !second !! Parameter of conversion precurseur to haze kch4=1. ncratio=0.5 ! boost for haze considering nitrogen contribution ! ratio n/c : =0.25 if there is 1N for 3C IF (firstcall) then write(*,*) 'hazecloud: haze parameters:' write(*,*) 'tcon, kch4, ncratio = ' , tcon, kch4, ncratio firstcall=.false. ENDIF ! note: mu0 = cos(solar zenith angle) ! max(mu0) = declin !!---------------------------------------------------------- !!---------------------------------------------------------- DO ig=1,ngrid zq_ch4(ig,:)=0. zq_prec(ig,:)=0. zdqconv_prec(ig,:)=0. zdqhaze(ig,:,:)=0. zdqphot_prec(ig,:)=0. zdqphot_ch4(ig,:)=0. tauch4(:)=0. gradflux(:)=0. fluxlym_sol(1:nlayer)=0. fluxlym_sol(nlayer+1)=flym_sol_pluto*mu_sol(ig) ! fluxlym_ipm(nlayer+1)= flym_ipm(ig) DO l=nlayer,1,-1 !! Actualisation tracer ch4 and prec_haze !IF (ngrid.eq.1) THEN !! option zq_ch4 = cte ! zq_ch4(ig,l)=0.01*0.5*16./28. ! Temporaire !ELSE zq_ch4(ig,l)=pq(ig,l,igcm_ch4_gas)+ & pdq(ig,l,igcm_ch4_gas)*ptimestep !ENDIF if (zq_ch4(ig,l).lt.0.) then zq_ch4(ig,l)=0. endif !! option zq_ch4 = cte ! zq_ch4(ig,l)=0.1*16./28. ! Temporaire zq_prec(ig,l)=pq(ig,l,igcm_prec_haze)+ & pdq(ig,l,igcm_prec_haze)*ptimestep !! Calculation optical depth ch4 in Lyman alpha for each layer l tauch4(l)=sigch4*1.e-4*avogadro* & (zq_ch4(ig,l)/(mmol(igcm_ch4_gas)*1.e-3))* & (pplev(ig,l)-pplev(ig,l+1))/g !! Calculation of Flux in each layer l if (mu_sol(ig).gt.1.e-6) then fluxlym_sol(l)=fluxlym_sol(l+1)*exp(-tauch4(l)/mu_sol(ig)) endif fluxlym_ipm(l)=fluxlym_ipm(l+1)*exp(-tauch4(l)/mu_ipm(ig)) gradflux(l)=fluxlym_sol(l+1)-fluxlym_sol(l) & + fluxlym_ipm(l+1)-fluxlym_ipm(l) !! tendancy on ch4 !! considering 1 photon destroys 1 ch4 by photolysis zdqphot_ch4(ig,l)=-mmol(igcm_ch4_gas)*1.e-3/avogadro* & gradflux(l)*g/(pplev(ig,l)-pplev(ig,l+1)) !! tendency of prec created by photolysis of ch4 zdqphot_prec(ig,l)=-zdqphot_ch4(ig,l) !! update precurseur zq_prec zq_prec(ig,l)=zq_prec(ig,l)+ & zdqphot_prec(ig,l)*ptimestep !! Conversion of prec_haze into haze !! controlled by the time constant tcon !! New tendancy for prec_haze zdqconv_prec(ig,l) = -zq_prec(ig,l)* & (1-exp(-ptimestep/tcon))/ptimestep ENDDO ! nlayer !! Final tendancy for prec_haze and haze DO iq=1,nq tname=noms(iq) !print*, 'TB17 tname=',tname,tname(1:4) if (tname(1:4).eq."haze") then zdqhaze(ig,:,iq) = -zdqconv_prec(ig,:)* & kch4*(1.+ncratio*14./12.) else if (noms(iq).eq."prec_haze") then zdqhaze(ig,:,igcm_prec_haze)= zdqphot_prec(ig,:) + & zdqconv_prec(ig,:) else if (noms(iq).eq."ch4_gas".and.(.not.hazeconservch4)) then zdqhaze(ig,:,igcm_ch4_gas)= zdqphot_ch4(ig,:) endif ENDDO ENDDO ! ngrid !! tendency kg/m2/s for haze column: ! zdqhaze_col(:)=0. ! DO ig=1,ngrid ! DO l=1,nlayer ! zdqhaze_col(ig)=zdqhaze_col(ig)+zdqhaze(ig,l,igcm_haze)* & ! (pplev(ig,l)-pplev(ig,l+1))/g ! ENDDO ! ENDDO end subroutine hazecloud