| [57] | 1 | subroutine lwu (kdlon,kflev |
|---|
| 2 | & ,dp,plev,tlay,aerosol |
|---|
| 3 | & ,aer_t,co2_u,co2_up |
|---|
| 4 | & ,tautotal,omegtotal,gtotal) |
|---|
| 5 | |
|---|
| 6 | c---------------------------------------------------------------------- |
|---|
| 7 | c LWU computes - co2: longwave effective absorber amounts including |
|---|
| 8 | c pressure and temperature effects |
|---|
| 9 | c - aerosols: amounts for every band |
|---|
| 10 | c transmission for bandes 1 and 2 of co2 |
|---|
| 11 | c---------------------------------------------------------------------- |
|---|
| 12 | |
|---|
| 13 | c----------------------------------------------------------------------- |
|---|
| 14 | c ATTENTION AUX UNITES: |
|---|
| 15 | c le facteur 10*g fait passer des kg m-2 aux g cm-2 |
|---|
| 16 | c----------------------------------------------------------------------- |
|---|
| 17 | c! modif diffusion |
|---|
| 18 | c! on ne change rien a la bande CO2 : les quantites d'absorbant CO2 |
|---|
| 19 | c! sont multipliees par 1.66 |
|---|
| 20 | c! pview= 1/cos(teta0)=1.66 |
|---|
| 21 | c----------------------------------------------------------------------- |
|---|
| 22 | |
|---|
| 23 | implicit none |
|---|
| 24 | |
|---|
| 25 | #include "dimensions.h" |
|---|
| 26 | #include "dimphys.h" |
|---|
| 27 | #include "dimradmars.h" |
|---|
| 28 | #include "comcstfi.h" |
|---|
| 29 | |
|---|
| 30 | #include "yomaer.h" |
|---|
| 31 | #include "yomlw.h" |
|---|
| 32 | |
|---|
| 33 | #include "fisice.h" |
|---|
| 34 | #include "callkeys.h" |
|---|
| 35 | |
|---|
| 36 | #include "aerice.h" |
|---|
| 37 | |
|---|
| 38 | c---------------------------------------------------------------------- |
|---|
| 39 | c 0.1 arguments |
|---|
| 40 | c --------- |
|---|
| 41 | c inputs: |
|---|
| 42 | c ------- |
|---|
| 43 | integer kdlon ! part of ngrid |
|---|
| 44 | integer kflev ! part of nalyer |
|---|
| 45 | |
|---|
| 46 | real dp (ndlo2,kflev) ! layer pressure thickness (Pa) |
|---|
| 47 | real plev (ndlo2,kflev+1) ! level pressure (Pa) |
|---|
| 48 | real tlay (ndlo2,kflev) ! layer temperature (K) |
|---|
| 49 | real aerosol (ndlo2,kflev,naerkind) ! aerosol extinction optical depth |
|---|
| 50 | c at reference wavelength "longrefvis" set |
|---|
| 51 | c in dimradmars.h , in each layer, for one of |
|---|
| 52 | c the "naerkind" kind of aerosol optical properties. |
|---|
| 53 | |
|---|
| 54 | |
|---|
| 55 | c outputs: |
|---|
| 56 | c -------- |
|---|
| 57 | real aer_t (ndlo2,nuco2,kflev+1) ! transmission (aer) |
|---|
| 58 | real co2_u (ndlo2,nuco2,kflev+1) ! absorber amounts (co2) |
|---|
| 59 | real co2_up (ndlo2,nuco2,kflev+1) ! idem scaled by the pressure (co2) |
|---|
| 60 | |
|---|
| 61 | real tautotal(ndlo2,kflev,nir) ! \ Total single scattering |
|---|
| 62 | real omegtotal(ndlo2,kflev,nir) ! > properties (Addition of the |
|---|
| 63 | real gtotal(ndlo2,kflev,nir) ! / NAERKIND aerosols properties) |
|---|
| 64 | |
|---|
| 65 | c---------------------------------------------------------------------- |
|---|
| 66 | c 0.2 local arrays |
|---|
| 67 | c ------------ |
|---|
| 68 | |
|---|
| 69 | integer jl,jk,jkl,ja,n |
|---|
| 70 | |
|---|
| 71 | real aer_u (ndlon,nir,nflev+1) ! absorber amounts (aer) extinction |
|---|
| 72 | real co2c ! co2 concentration (pa/pa) |
|---|
| 73 | real pview ! cosecant of viewing angle |
|---|
| 74 | real pref ! reference pressure (1013 mb = 101325 Pa) |
|---|
| 75 | real tx,tx2 |
|---|
| 76 | real phi (ndlon,nuco2) |
|---|
| 77 | real psi (ndlon,nuco2) |
|---|
| 78 | real plev2 (ndlon,nflev+1) |
|---|
| 79 | real zzz |
|---|
| 80 | |
|---|
| 81 | real ray,coefsize,coefsizew,coefsizeg |
|---|
| 82 | |
|---|
| 83 | c************************************************************************ |
|---|
| 84 | c---------------------------------------------------------------------- |
|---|
| 85 | c 0.3 Initialisation |
|---|
| 86 | c ------------- |
|---|
| 87 | |
|---|
| 88 | pview = 1.66 |
|---|
| 89 | co2c = 0.95 |
|---|
| 90 | pref = 101325. |
|---|
| 91 | |
|---|
| 92 | do jk=1,nlaylte+1 |
|---|
| 93 | do jl=1,kdlon |
|---|
| 94 | plev2(jl,jk)=plev(jl,jk)*plev(jl,jk) |
|---|
| 95 | enddo |
|---|
| 96 | enddo |
|---|
| 97 | |
|---|
| 98 | c---------------------------------------------------------------------- |
|---|
| 99 | c Computing TOTAL single scattering parameters by adding properties of |
|---|
| 100 | c all the NAERKIND kind of aerosols in each IR band |
|---|
| 101 | |
|---|
| 102 | call zerophys(ndlon*kflev*nir,tautotal) |
|---|
| 103 | call zerophys(ndlon*kflev*nir,omegtotal) |
|---|
| 104 | call zerophys(ndlon*kflev*nir,gtotal) |
|---|
| 105 | |
|---|
| 106 | do n=1,naerkind |
|---|
| 107 | do ja=1,nir |
|---|
| 108 | do jk=1,nlaylte |
|---|
| 109 | do jl = 1,kdlon |
|---|
| 110 | c TEST : to account for the varying sol/ir optical depth of ice with varying crystal size |
|---|
| 111 | c : and for varying w,g with varying crystal size |
|---|
| 112 | if (activice.and.n.eq.naerkind) then |
|---|
| 113 | ray=min( max(rice(jl,jk)*1.e+6, 1.),10.) |
|---|
| 114 | if (ja.eq.1.or.ja.eq.2) then |
|---|
| 115 | coefsize=(-0.00382*ray**3.+0.0503*ray**2.+0.03531) |
|---|
| 116 | & /QIRsQREF(ja,n) |
|---|
| 117 | coefsizew=(-0.011*ray**2.+0.1824*ray-0.1283) |
|---|
| 118 | & /omegaIR(ja,n) |
|---|
| 119 | coefsizeg=(-0.00122*ray**3.+0.009161*ray**2.+ |
|---|
| 120 | & 0.1182*ray-0.096877)/gIR(ja,n) |
|---|
| 121 | elseif (ja.eq.3) then |
|---|
| 122 | coefsize=(-0.00324*ray**3.+0.0419*ray**2.+0.0459) |
|---|
| 123 | & /QIRsQREF(ja,n) |
|---|
| 124 | coefsizew=(-0.01292*ray**2.+0.1963*ray-0.06566) |
|---|
| 125 | & /omegaIR(ja,n) |
|---|
| 126 | coefsizeg=(0.00271*ray**3.-0.05959*ray**2.+ |
|---|
| 127 | & 0.4411*ray-0.2724)/gIR(ja,n) |
|---|
| 128 | elseif (ja.eq.4) then |
|---|
| 129 | coefsize=(-0.0003823*ray**3.+0.0104*ray**2.+0.005966) |
|---|
| 130 | & /QIRsQREF(ja,n) |
|---|
| 131 | coefsizew=(-0.002*ray**3.+0.02623*ray**2.-0.014465) |
|---|
| 132 | & /omegaIR(ja,n) |
|---|
| 133 | coefsizeg=(-0.0017192*ray**3.+0.0259*ray**2.- |
|---|
| 134 | & 0.027692*ray+0.016099)/gIR(ja,n) |
|---|
| 135 | endif |
|---|
| 136 | if (coefsize.le.0.or.coefsizew.le.0.or.coefsizeg.le.0) |
|---|
| 137 | & stop'pb dans lwu avec prp opt glace' |
|---|
| 138 | |
|---|
| 139 | tautotal(jl,jk,ja)=tautotal(jl,jk,ja) + |
|---|
| 140 | & QIRsQREF(ja,n)*aerosol(jl,jk,n) |
|---|
| 141 | & *coefsize |
|---|
| 142 | omegtotal(jl,jk,ja) = omegtotal(jl,jk,ja) + |
|---|
| 143 | & QIRsQREF(ja,n)*aerosol(jl,jk,n)*omegaIR(ja,n) |
|---|
| 144 | & *coefsize*coefsizew |
|---|
| 145 | gtotal(jl,jk,ja) = gtotal(jl,jk,ja) + |
|---|
| 146 | & QIRsQREF(ja,n)*aerosol(jl,jk,n)*omegaIR(ja,n)*gIR(ja,n) |
|---|
| 147 | & *coefsize*coefsizew*coefsizeg |
|---|
| 148 | else |
|---|
| 149 | tautotal(jl,jk,ja)=tautotal(jl,jk,ja) + |
|---|
| 150 | & QIRsQREF(ja,n)*aerosol(jl,jk,n) |
|---|
| 151 | omegtotal(jl,jk,ja) = omegtotal(jl,jk,ja) + |
|---|
| 152 | & QIRsQREF(ja,n)*aerosol(jl,jk,n)*omegaIR(ja,n) |
|---|
| 153 | gtotal(jl,jk,ja) = gtotal(jl,jk,ja) + |
|---|
| 154 | & QIRsQREF(ja,n)*aerosol(jl,jk,n)*omegaIR(ja,n)*gIR(ja,n) |
|---|
| 155 | endif |
|---|
| 156 | enddo |
|---|
| 157 | enddo |
|---|
| 158 | enddo |
|---|
| 159 | enddo |
|---|
| 160 | do ja=1,nir |
|---|
| 161 | do jk=1,nlaylte |
|---|
| 162 | do jl = 1,kdlon |
|---|
| 163 | gtotal(jl,jk,ja)=gtotal(jl,jk,ja)/omegtotal(jl,jk,ja) |
|---|
| 164 | omegtotal(jl,jk,ja)=omegtotal(jl,jk,ja)/tautotal(jl,jk,ja) |
|---|
| 165 | enddo |
|---|
| 166 | enddo |
|---|
| 167 | enddo |
|---|
| 168 | |
|---|
| 169 | c---------------------------------------------------------------------- |
|---|
| 170 | c 1.0 cumulative (aerosol) amounts (for every band) |
|---|
| 171 | c ---------------------------- |
|---|
| 172 | |
|---|
| 173 | jk=nlaylte+1 |
|---|
| 174 | do ja=1,nir |
|---|
| 175 | do jl = 1 , kdlon |
|---|
| 176 | aer_u(jl,ja,jk)=0. |
|---|
| 177 | enddo |
|---|
| 178 | enddo |
|---|
| 179 | |
|---|
| 180 | do jk=1,nlaylte |
|---|
| 181 | jkl=nlaylte+1-jk |
|---|
| 182 | do ja=1,nir |
|---|
| 183 | do jl=1,kdlon |
|---|
| 184 | aer_u(jl,ja,jkl)=aer_u(jl,ja,jkl+1)+ tautotal(jl,jkl,ja) |
|---|
| 185 | enddo |
|---|
| 186 | enddo |
|---|
| 187 | enddo |
|---|
| 188 | |
|---|
| 189 | c---------------------------------------------------------------------- |
|---|
| 190 | c 1.0 bands 1 and 2 of co2 |
|---|
| 191 | c -------------------- |
|---|
| 192 | |
|---|
| 193 | jk=nlaylte+1 |
|---|
| 194 | do ja=1,nuco2 |
|---|
| 195 | do jl = 1 , kdlon |
|---|
| 196 | co2_u(jl,ja,jk)=0. |
|---|
| 197 | co2_up(jl,ja,jk)=0. |
|---|
| 198 | aer_t(jl,ja,jk)=1. |
|---|
| 199 | enddo |
|---|
| 200 | enddo |
|---|
| 201 | |
|---|
| 202 | do jk=1,nlaylte |
|---|
| 203 | jkl=nlaylte+1-jk |
|---|
| 204 | do ja=1,nuco2 |
|---|
| 205 | do jl=1,kdlon |
|---|
| 206 | |
|---|
| 207 | c introduces temperature effects on absorber(co2) amounts |
|---|
| 208 | c ------------------------------------------------------- |
|---|
| 209 | tx = sign(min(abs(tlay(jl,jkl)-tref),70.) |
|---|
| 210 | . ,tlay(jl,jkl)-tref) |
|---|
| 211 | tx2=tx*tx |
|---|
| 212 | phi(jl,ja)=at(1,ja)*tx+bt(1,ja)*tx2 |
|---|
| 213 | psi(jl,ja)=at(2,ja)*tx+bt(2,ja)*tx2 |
|---|
| 214 | phi(jl,ja)=exp(phi(jl,ja)/cst_voigt(2,ja)) |
|---|
| 215 | psi(jl,ja)=exp(2.*psi(jl,ja)) |
|---|
| 216 | |
|---|
| 217 | c cumulative absorber(co2) amounts |
|---|
| 218 | c -------------------------------- |
|---|
| 219 | co2_u(jl,ja,jkl)=co2_u(jl,ja,jkl+1) |
|---|
| 220 | . + pview/(10*g)*phi(jl,ja)*dp(jl,jkl)*co2c |
|---|
| 221 | |
|---|
| 222 | co2_up(jl,ja,jkl)=co2_up(jl,ja,jkl+1) |
|---|
| 223 | . + pview/(10*g*2*pref)*psi(jl,ja) |
|---|
| 224 | . * (plev2(jl,jkl)-plev2(jl,jkl+1))*co2c |
|---|
| 225 | |
|---|
| 226 | |
|---|
| 227 | c (aerosol) transmission |
|---|
| 228 | c ---------------------- |
|---|
| 229 | c on calcule directement les transmissions pour les aerosols. |
|---|
| 230 | c on multiplie le Qext par 1-omega dans la bande du CO2. |
|---|
| 231 | c et pourquoi pas d'abord? hourdin@lmd.ens.fr |
|---|
| 232 | |
|---|
| 233 | zzz=pview*(1.-omegtotal(jl,jkl,ja)) |
|---|
| 234 | aer_t(jl,ja,jkl)=exp(-zzz*aer_u(jl,ja,jkl)) |
|---|
| 235 | |
|---|
| 236 | enddo |
|---|
| 237 | enddo |
|---|
| 238 | enddo |
|---|
| 239 | |
|---|
| 240 | c---------------------------------------------------------------------- |
|---|
| 241 | return |
|---|
| 242 | end |
|---|
