source: TOOLS/CMIP6_FORCING/AER_OPTICS/refr_ind_sulfate_MACC_v2.f @ 3385

Last change on this file since 3385 was 3385, checked in by oboucher, 4 years ago

Uploading package to prepare tropospheric aerosol optics LUT in LMDz fo rCMIP6

File size: 3.3 KB
Line 
1      PROGRAM ammonium_sulfate
2!--based on Tang and Munkelwitz JGR 1994 and Tang JGR 1997
3      IMPLICIT NONE
4      REAL xms  !--solute mass fraction
5      REAL aw !-water activity (RH here)
6      REAL molality !-molality
7      REAL rhop !-dry droplet density
8      REAL nr !--real refractive index
9      REAL RR, Rwater, Ramm !--molar refraction
10      PARAMETER (Rwater=3.717, Ramm=23.50) 
11      REAL VV !--molal volume
12      REAL Mwater, Mamm 
13      PARAMETER (Mwater=18.02, Mamm=132.14)
14      REAL y !--mole fraction
15      PARAMETER (rhop=1.769) !--g cm-3
16      REAL B0, B1, B2, B3, B4 
17      PARAMETER (B0=110.65495, B1=-367.59197,B2=504.62934) 
18      PARAMETER (B3=-315.43839,B4=67.70824)
19      REAL C1, C2, C3, C4 
20      PARAMETER (C1=-2.715E-3,C2=3.113E-5,C3=-2.336E-6,C4=1.412E-8)
21      REAL A0, A1, A2, A3
22      PARAMETER (A0=0.9971, A1=5.92E-3, A2=-5.036E-6, A3=1.024E-8)
23c-------------------------
24      INTEGER nbre_rh,IRH,wve
25      INTEGER nbre_dry,ii
26      PARAMETER(nbre_rh=12,nbre_dry=4)
27      REAL DELTA
28      REAL RH_tab(nbre_rh)
29      DATA RH_tab/0.,10.,20.,30.,40.,50.,60.,70.,80.,85.,90.,95./
30      INTEGER Nwvmax
31      PARAMETER(Nwvmax=100)
32c
33c--wavelength in m
34      REAL lambda_min, lambda_max, lambda(Nwvmax)
35      PARAMETER (lambda_min=0.2E-6, lambda_max=5.0E-6)
36
37      REAL n_rr(nbre_rh,Nwvmax)
38      REAL n_ii(nbre_rh,Nwvmax)
39      REAL ratio(nbre_rh)      !--diameter growth factor
40      REAL rhod(nbre_rh)       !--density
41      REAL n_r_exact(nbre_rh)  !--refractive index
42c
43      DO wve=1,Nwvmax
44        lambda(wve)=
45     .  lambda_min+FLOAT(wve-1)*(lambda_max-lambda_min)/FLOAT(Nwvmax-1)
46      ENDDO 
47c
48c--xms is the varying variable as per Tang and Munkelwitz 1994
49      OPEN (unit=4,file='ri_sul_v2')
50      OPEN (unit=10,file='growth_sul_v2')
51c
52      DO IRH=1,nbre_rh
53c
54      IF (IRH.LE.nbre_dry) THEN
55c
56      ratio(IRH)=1.0
57      rhod(IRH)=rhop
58      n_r_exact(IRH)=1.521
59c
60      ELSE
61c
62      aw=rh_tab(IRH)/100.   !--water activity
63      molality=B0+aw*(B1+aw*(B2+aw*(B3+aw*B4)))   !--molality
64      xms=(molality*Mamm/1000.)/(molality*Mamm/1000.+1.)*100. !--solute mass fraction
65c
66      rhod(IRH)=A0+xms*(A1+xms*(A2+xms*A3))       !--density
67      ratio(IRH)=(rhop/rhod(IRH))**(1./3.) * (xms/100.)**(-1./3.)  !--growth factor
68      y=Mwater*xms/100./(Mamm*(1.-xms/100.)+xms/100.*Mwater)
69      VV=1./rhod(IRH)*(Mwater+(Mamm-Mwater)*y)
70      RR=Rwater+(Ramm-Rwater)*y
71      nr=(2.*RR/VV+1.)/(1.-RR/VV)
72      n_r_exact(IRH)=sqrt(nr)
73c
74      ENDIF
75c
76      ENDDO
77c
78      write(10,20)'RH', (rh_tab(IRH),IRH=1,nbre_rh)
79      write(10,10)"growth factor",(ratio(IRH),IRH=1,nbre_rh)
80      write(10,10)"density factor",(rhod(IRH)/rhod(1),IRH=1,nbre_rh)
81      write(10,10)'rho @ RH', (rhod(IRH),IRH=1,nbre_rh)
82      write(10,10)'nr @ RH', (n_r_exact(IRH),IRH=1,nbre_rh)
8310    FORMAT(A14, F6.3,11(',',F6.3))
8420    FORMAT(A14, F6.1,11(',',F6.1))
85
86c
87c-----Computation of n_r at different wavelengths using
88c     n_r(lambda)=n_r(0.589)-0.03(lambda-0.589) from 
89c     Palmer and Williams (1975) and  Kent et al. (1983).
90c     Above computed n_r_exact is assumed to be at 0.589 nm
91c
92      DO IRH=1,nbre_rh
93      DO wve=1,Nwvmax
94c
95      n_rr(IRH,wve)=n_r_exact(irh)-0.03*(lambda(wve)*1.E6-0.589)
96      n_ii(IRH,wve)=0.0
97c
98      WRITE(4,*) RH_tab(IRH), lambda(wve), n_rr(IRH,wve), n_ii(IRH,wve)
99c
100      ENDDO
101      ENDDO
102c
10363    FORMAT(F10.6,E13.6)
104c
105      END
Note: See TracBrowser for help on using the repository browser.