source: LMDZ6/trunk/libf/phylmd/inlandsis/sisvat_sop.F @ 3871

Last change on this file since 3871 was 3792, checked in by evignon, 4 years ago

Ajout de INLANDSIS, nouvelle interface entre LMDZ et la neige de SISVAT
Etienne, 04/01/2021
File size: 20.7 KB
Line 
1 
2 
3      subroutine SnOptP(jjtime)
4 
5C +------------------------------------------------------------------------+
6C | MAR/SISVAT   SnOptP                                    12-08-2019  MAR |
7C |   SubRoutine SnOptP computes the Snow Pack optical Properties          |
8C +------------------------------------------------------------------------+
9C |                                                                        |
10C |   PARAMETERS:  knonv: Total Number of columns =                        |
11C |   ^^^^^^^^^^        = Total Number of continental     Grid Boxes       |
12C |                     X       Number of Mosaic Cell per Grid Box         |
13C |                                                                        |
14C |   INPUT:   isnoSV   = total Nb of Ice/Snow Layers                      |
15C |   ^^^^^    ispiSV   = 0,...,nsno: Uppermost Superimposed Ice Layer     |
16C |                                                                        |
17C |                                                                        |
18C |   INPUT:   G1snSV   : Dendricity (<0) or Sphericity (>0) of Snow Layer |
19C |   ^^^^^    G2snSV   : Sphericity (>0) or Size            of Snow Layer |
20C |            agsnSV   : Snow       Age                             [day] |
21C |            ro__SV   : Snow/Soil  Volumic Mass                  [kg/m3] |
22C |            eta_SV   : Water      Content                       [m3/m3] |
23C |            rusnSV   : Surficial  Water   Thickness   [kg/m2] .OR. [mm] |
24C |            SWS_SV   : Surficial  Water   Status                        |
25C |            dzsnSV   : Snow       Layer   Thickness                 [m] |
26C |                                                                        |
27C |            albssv   : Soil       Albedo                            [-] |
28C |            zzsnsv   : Snow       Pack    Thickness                 [m] |
29C |                                                                        |
30C |   OUTPUT:  albisv   : Snow/Ice/Water/Soil Integrated Albedo        [-] |
31C |   ^^^^^^   sEX_sv   : Verticaly Integrated Extinction Coefficient      |
32C |                                                                        |
33C |   Internal Variables:                                                  |
34C |   ^^^^^^^^^^^^^^^^^^                                                   |
35C |            SnOpSV   : Snow Grain optical Size                      [m] |
36C |            EX1_sv   : Integrated Snow Extinction (0.3--0.8micr.m)      |
37C |            EX2_sv   : Integrated Snow Extinction (0.8--1.5micr.m)      |
38C |            EX3_sv   : Integrated Snow Extinction (1.5--2.8micr.m)      |
39C |                                                                        |
40C |   METHODE:    Calcul de la taille optique des grains ? partir de       |
41C |   ^^^^^^^    -leur type decrit par les deux variables descriptives     |
42C |                    continues sur la plage -99/+99 passees en appel.    |
43C |              -la taille optique (1/10mm) des etoiles,                  |
44C |                                          des grains fins et            |
45C |                                          des jeunes faces planes       |
46C |                                                                        |
47C |   METHOD:     Computation of the optical diameter of the grains        |
48C |   ^^^^^^      described with the CROCUS formalism G1snSV / G2snSV      |
49C |                                                                        |
50C |   REFERENCE: Brun et al.      1989, J. Glaciol 35 pp. 333--342         |
51C |   ^^^^^^^^^  Brun et al.      1992, J. Glaciol 38 pp.  13-- 22         |
52C |              Eric Martin Sept.1996                                     |
53C |                                                                        |
54C |                                                                        |
55C +------------------------------------------------------------------------+
56 
57 
58 
59 
60C +--Global Variables
61C +  ================
62 
63
64      use VARphy
65      use VAR_SV
66      use VARdSV
67      use VARxSV
68      use VARySV
69      use VARtSV
70      USE surface_data, only: iflag_albzenith
71
72      IMPLICIT NONE
73
74 
75C + -- INPUT
76      integer jjtime
77
78C +--Internal Variables
79C +  ==================
80 
81      real      coalb1(knonv)                      ! weighted Coalbedo, Vis.
82      real      coalb2(knonv)                      ! weighted Coalbedo, nIR 1
83      real      coalb3(knonv)                      ! weighted Coalbedo, nIR 2
84      real      coalbm                             ! weighted Coalbedo, mean
85      real      sExt_1(knonv)                      ! Extinction Coeff., Vis.
86      real      sExt_2(knonv)                      ! Extinction Coeff., nIR 1
87      real      sExt_3(knonv)                      ! Extinction Coeff., nIR 2
88      real      SnOpSV(knonv,      nsno)           ! Snow Grain optical Size
89c #AG real      agesno
90 
91      integer   isn   ,ikl   ,isn1 
92      real      sbeta1,sbeta2,sbeta3,sbeta4,sbeta5
93      real      AgeMax,AlbMin,HSnoSV,HIceSV,doptmx,SignG1,Sph_OK
94      real      dalbed,dalbeS,dalbeW
95      real      bsegal,czemax,csegal
96      real      RoFrez,DiffRo,SignRo,SnowOK,OpSqrt
97      real      albSn1,albIc1,a_SnI1,a_SII1
98      real      albSn2,albIc2,a_SnI2,a_SII2
99      real      albSn3,albIc3,a_SnI3,a_SII3
100      real      albSno,albIce,albSnI,albSII,albWIc,albmax
101      real      doptic,Snow_H,SIce_H,SnownH,SIcenH
102      real      exarg1,exarg2,exarg3,sign_0,sExt_0
103      real      albedo_old,albCor
104      real      ro_ave,dz_ave,minalb
105      real      thetazs,thetazs0,aap,bbp,ccp,alb0      ! parameters for zenoth angle effect at Dome C
106
107
108 
109C +--Local   DATA
110C +  ============
111 
112C +--For the computation of the solar irradiance extinction in snow
113C +  ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
114      data      sbeta1/0.0192/,sbeta2/0.4000/,sbeta3/0.1098/
115      data      sbeta4/1.0000/
116      data      sbeta5/2.00e1/
117 
118C +--Snow Age Maximum (Taiga, e.g. Col de Porte)
119C +  ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
120      data      AgeMax  /60.0/
121C +...          AgeMax:  Snow Age Maximum                              [day]
122 
123      data      AlbMin  /0.94/
124C +...          AlbMin:  Albedo   Minimum / visible (0.3--0.8 micrometers)
125 
126      data      HSnoSV  /0.01/
127C +...          HSnoSV:  Snow     Thickness through witch
128C +                      Albedo is interpolated to Ice  Albedo
129      data      HIceSV  /0.10/
130C +...          HIceSV:  Snow/Ice Thickness through witch
131C +                      Albedo is interpolated to Soil Albedo
132      data      doptmx  /2.3e-3/
133C +...          doptmx:  Maximum  optical Diameter    (pi * R**2)        [m]
134C +
135      data      czeMAX  /0.173648178/  ! 80.deg (Segal et al., 1991 JAS)
136
137      data      bsegal  /4.00       /  !
138      data      albmax  /0.99       /  ! Albedo max
139
140C +-- For the computation of solar zenoth angle effect from Dome C data
141C +  ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
142      data      aap/0.00016/,bbp/-0.017/,ccp/1.2/
143
144 
145 
146C +--Snow Grain optical Size
147C +  =======================
148 
149        DO ikl=1,knonv
150         DO   isn=1,max(1,isnoSV(ikl))
151 
152          G2snSV(ikl,isn) =  max(epsi,G2snSV(ikl,isn))
153C +...    Avoid non physical Values
154 
155          SignG1          = sign(unun,G1snSV(ikl,isn))
156          Sph_OK          =  max(zero,SignG1)
157 
158          SnOpSV(ikl,isn) =   1.e-4 *
159C +...    SI:           (from 1/10 mm to m)
160 
161 
162C +--Contribution of Non Dendritic Snow
163C +  ----------------------------------
164 
165     .    (    Sph_OK *(      G2snSV(ikl,isn)*G1snSV(ikl,isn)/G1_dSV
166     .              +max(demi*G2snSV(ikl,isn),DFcdSV)
167     .                 *(unun-G1snSV(ikl,isn)                /G1_dSV))
168 
169 
170C +--Contribution of     Dendritic Snow
171C +  ----------------------------------
172 
173     .    +(1.-Sph_OK)*(     -G1snSV(ikl,isn)*DDcdSV         /G1_dSV
174     .                 +(unun+G1snSV(ikl,isn)                /G1_dSV)
175     .                  *    (G2snSV(ikl,isn)*DScdSV         /G1_dSV
176     .                 +(unun-G2snSV(ikl,isn)                /G1_dSV)
177     .                                       *DFcdSV                 )))
178          SnOpSV(ikl,isn) =  max(zero,SnOpSV(ikl,isn))
179        END DO
180      END DO
181 
182 
183C +--Snow/Ice Albedo
184C +  ===============
185
186 
187 
188C +--Uppermost effective Snow Layer
189C +  ------------------------------
190 
191        DO ikl=1,knonv
192 
193           isn   =  max(iun,isnoSV(ikl))
194 
195          SignRo = sign(unun, rocdSV - ro__SV(ikl,isn))
196          SnowOK =  max(zero,SignRo) ! Ice Density Threshold
197 
198          OpSqrt = sqrt(SnOpSV(ikl,isn))
199 
200cCA    +--Correction of snow albedo for Antarctica/Greenland
201cCA       --------------------------------------------------
202          albCor = 1.
203c #GL     albCor = 1.01
204c #AC     albCor = 1.01
205 
206          albSn1 =  0.96-1.580*OpSqrt
207          albSn1 =  max(albSn1,AlbMin)
208 
209          albSn1 =  max(albSn1,zero)
210          albSn1 =  min(albSn1*albCor,unun)
211 
212          albSn2 =  0.95-15.40*OpSqrt
213          albSn2 =  max(albSn2,zero)
214          albSn2 =  min(albSn2*albCor,unun)
215 
216          doptic =  min(SnOpSV(ikl,isn),doptmx)
217          albSn3 =  346.3*doptic -32.31*OpSqrt +0.88
218          albSn3 =  max(albSn3,zero)
219          albSn3 =  min(albSn3*albCor,unun)
220 
221          !snow albedo corection if wetsnow
222c #GL     albSn1 =  albSn1*max(0.9,(1.-1.5*eta_SV(ikl,isn)))
223c #GL     albSn2 =  albSn2*max(0.9,(1.-1.5*eta_SV(ikl,isn)))
224c #GL     albSn3 =  albSn3*max(0.9,(1.-1.5*eta_SV(ikl,isn)))
225 
226          albSno =  So1dSV*albSn1
227     .           +  So2dSV*albSn2
228     .           +  So3dSV*albSn3
229 
230cXF
231          minalb =  (aI2dSV + (aI3dSV -aI2dSV)
232     .           *  (ro__SV(ikl,isn)-ro_Ice)/(roSdSV-ro_Ice))
233          minalb =  min(aI3dSV,max(aI2dSV,minalb)) ! pure/firn albedo
234 
235          SnowOK =  SnowOK*max(zero,sign(unun,     albSno-minalb))
236          albSn1 =  SnowOK*albSn1+(1.0-SnowOK)*max(albSno,minalb)
237          albSn2 =  SnowOK*albSn2+(1.0-SnowOK)*max(albSno,minalb)
238          albSn3 =  SnowOK*albSn3+(1.0-SnowOK)*max(albSno,minalb)
239 
240c +           ro < roSdSV |          min al > aI3dSV
241c +  roSdSV < ro < rocdSV | aI2dSV < min al < aI3dSV (fct of density)
242 
243 
244C +--Snow/Ice Pack Thickness
245C +  -----------------------
246 
247          isn    =    max(min(isnoSV(ikl) ,ispiSV(ikl)),0)
248          Snow_H =  zzsnsv(ikl,isnoSV(ikl))-zzsnsv(ikl,isn)
249          SIce_H =  zzsnsv(ikl,isnoSV(ikl))
250          SnownH =  Snow_H  /  HSnoSV
251          SnownH =  min(unun,  SnownH)
252          SIcenH =  SIce_H  / (HIceSV)
253          SIcenH =  min(unun,  SIcenH)
254 
255C +       The value of SnownH is set to 1 in case of ice lenses above
256C +       1m of dry snow (ro<600kg/m3) for using CROCUS albedo
257 
258c         ro_ave =  0.
259c         dz_ave =  0.
260c         SnowOK =  1.
261c      do isn    =  isnoSV(ikl),1,-1
262c         ro_ave =  ro_ave + ro__SV(ikl,isn) * dzsnSV(ikl,isn) * SnowOK
263c         dz_ave =  dz_ave +                   dzsnSV(ikl,isn) * SnowOK
264c         SnowOK =  max(zero,sign(unun,1.-dz_ave))
265c      enddo
266 
267c         ro_ave =  ro_ave / max(dz_ave,epsi)
268c         SnowOK =  max(zero,sign(unun,600.-ro_ave))
269c         SnownH =  SnowOK + SnownH * (1. - SnowOK)
270 
271C +--Integrated Snow/Ice Albedo: Case of Water on Bare Ice
272C +  -----------------------------------------------------
273 
274          isn    =  max(min(isnoSV(ikl) ,ispiSV(ikl)),0)
275 
276          albWIc =  aI1dSV-(aI1dSV-aI2dSV)
277     .           *  exp(-(rusnSV(ikl)                      !
278     .           *  (1. -SWS_SV(ikl)                       ! 0 <=> freezing
279     .           *  (1  -min(1,iabs(isn-isnoSV(ikl)))))    ! 1 <=> isn=isnoSV
280     .           /   ru_dSV)**0.50)                        !
281c         albWIc = max(aI1dSV,min(aI2dSV,albWIc+slopSV(ikl)*
282c    .             min(5.,max(1.,dx/10000.))))
283 
284          SignRo = sign(unun,ro_Ice-5.-ro__SV(ikl,isn))    ! RoSN<920kg/m3
285          SnowOK =  max(zero,SignRo)
286 
287          albWIc = (1. - SnowOK) * albWIc + SnowOK
288     .           * (aI2dSV + (aI3dSV -aI2dSV)
289     .           * (ro__SV(ikl,isn)-ro_Ice)/(roSdSV-ro_Ice))
290 
291c +  rocdSV < ro < ro_ice | aI2dSV< al <aI3dSV (fct of density)
292c +           ro > ro_ice | aI1dSV< al <aI2dSV (fct of superficial water content
293 
294 
295C +--Integrated Snow/Ice      Albedo
296C +  -------------------------------
297 
298          a_SII1      =     albWIc      +(albSn1-albWIc)     *SnownH
299          a_SII1      = min(a_SII1       ,albSn1)
300 
301          a_SII2      =     albWIc      +(albSn2-albWIc)     *SnownH
302          a_SII2      = min(a_SII2       ,albSn2)
303 
304          a_SII3      =     albWIc      +(albSn3-albWIc)     *SnownH
305          a_SII3      = min(a_SII3       ,albSn3)
306 
307cc #AG     agesno =      min(agsnSV(ikl,isn)          ,AgeMax)
308cc #AG     a_SII1      =     a_SII1      -0.175*agesno/AgeMax
309C +...                                   Impurities: Col de Porte Parameter.
310 
311 
312!    Zenith Angle Correction (Segal et al.,         1991, JAS 48, p.1025)
313!    ----------------------- (Wiscombe & Warren, dec1980, JAS   , p.2723)
314!                            (Warren,               1982,  RG   , p.  81)
315!                            (Vignon, PhD Thesis, pp 150, conditions at Dome C)
316!                            -------------------------------------------------
317 
318         
319          dalbed = 0.0
320
321          if (iflag_albzenith .GE. 0) then
322          csegal = max(czemax                   ,coszSV(ikl))
323c #cz     dalbeS =   ((bsegal+unun)/(unun+2.0*bsegal*csegal)
324c #cz.                -       unun                          )*0.32
325c #cz.              /  bsegal
326c #cz     dalbeS = max(dalbeS,zero)
327c #cz     dalbed =     dalbeS      *       min(1,isnoSV(ikl))
328 
329          dalbeW =(0.64 - csegal  )*0.0625  ! Warren 1982, RevGeo, fig.12b
330                                            ! 0.0625 = 5% * 1/0.8,   p.81
331                                            ! 0.64   = cos(50)
332          dalbed =     dalbeW      *       min(1,isnoSV(ikl))
333
334
335! Vignon PhD thesis, Dome C conditions
336
337            if (iflag_albzenith .EQ. 1) then
338            thetazs0=-bbp/(2.0*aap)
339            alb0=(bbp**2)/4./aap-(bbp**2)/2./aap+ccp
340            thetazs=max(acos(coszSV(ikl))/pi*180.0,thetazs0) ! in degrees
341
342
343            dalbeW = aap*(thetazs**2)+bbp*thetazs+ccp-alb0       
344            dalbed =     dalbeW      *       min(1,isnoSV(ikl))
345            end if
346
347
348            end if
349 
350C +--Elsewhere    Integrated Snow/Ice Albedo
351C +  ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
352c #cp     ELSE
353            albSII =     So1dSV*a_SII1
354     .                 + So2dSV*a_SII2
355     .                 + So3dSV*a_SII3
356c #cp     END IF
357 
358 
359C +--Integrated Snow/Ice/Soil Albedo
360C +  -------------------------------
361 
362            alb1sv(ikl) =     albssv(ikl) +(a_SII1-albssv(ikl))*SIcenH
363            alb1sv(ikl) = min(alb1sv(ikl)  ,a_SII1)
364 
365            alb2sv(ikl) =     albssv(ikl) +(a_SII2-albssv(ikl))*SIcenH
366            alb2sv(ikl) = min(alb2sv(ikl)  ,a_SII2)
367 
368            alb3sv(ikl) =     albssv(ikl) +(a_SII3-albssv(ikl))*SIcenH
369            alb3sv(ikl) = min(alb3sv(ikl)  ,a_SII3)
370 
371            albisv(ikl) =     albssv(ikl) +(albSII-albssv(ikl))*SIcenH
372            albisv(ikl) = min(albisv(ikl)  ,albSII)
373 
374 
375C +--Integrated Snow/Ice/Soil Albedo: Clouds Correction! Greuell & all., 1994
376C +  --------------------------------------------------! Glob.&t Planet.Change
377                                                       ! (9):91-114
378            alb1sv(ikl) = alb1sv(ikl) + 0.05 *(cld_SV(ikl)-0.5)*SIcenH
379     .                  + dalbed      *    (1.-cld_SV(ikl))
380            alb2sv(ikl) = alb2sv(ikl) + 0.05 *(cld_SV(ikl)-0.5)*SIcenH
381     .                  + dalbed      *    (1.-cld_SV(ikl))
382            alb3sv(ikl) = alb3sv(ikl) + 0.05 *(cld_SV(ikl)-0.5)*SIcenH
383     .                  + dalbed      *    (1.-cld_SV(ikl))
384            albisv(ikl) = albisv(ikl) + 0.05 *(cld_SV(ikl)-0.5)*SIcenH
385     .                  + dalbed      *    (1.-cld_SV(ikl))
386 
387C +--Integrated Snow/Ice/Soil Albedo: Minimum snow albedo = aI1dSV
388C +  -------------------------------------------------------------
389 
390            albedo_old  = albisv(ikl)
391            albisv(ikl) = max(albisv(ikl),aI1dSV   * SIcenH
392     .                  + albssv(ikl) *(1.0        - SIcenH))
393            alb1sv(ikl) = alb1sv(ikl) - 1.0/3.0             ! 33 %
394     .                  * (albedo_old-albisv(ikl)) / So1dSV
395            alb2sv(ikl) = alb2sv(ikl) - 1.0/3.0             ! 33 %
396     .                  * (albedo_old-albisv(ikl)) / So2dSV
397            alb3sv(ikl) = alb3sv(ikl) - 1.0/3.0             ! 33 %
398     .                  * (albedo_old-albisv(ikl)) / So3dSV
399 
400 
401C +--Integrated Snow/Ice/Soil Albedo: Maximum albedo = 95%
402C +  -----------------------------------------------------
403 
404            albedo_old  = albisv(ikl)
405            albisv(ikl) = min(albisv(ikl),0.95)
406            alb1sv(ikl) = alb1sv(ikl) - 1.0/3.0             ! 33 %
407     .                  * (albedo_old-albisv(ikl)) / So1dSV
408            alb2sv(ikl) = alb2sv(ikl) - 1.0/3.0             ! 33 %
409     .                  * (albedo_old-albisv(ikl)) / So2dSV
410            alb3sv(ikl) = alb3sv(ikl) - 1.0/3.0             ! 33 %
411     .                  * (albedo_old-albisv(ikl)) / So3dSV
412
413
414        !Sea Ice/snow permanent-interractive prescription from Nemo
415        !AO_CK 20/02/2020
416
417        ! No check if coupling update since MAR and NEMO albedo are too different
418        ! and since MAR albedo is computed on properties that are not in NEMO
419        ! prescription for each time step with NEMO values
420       
421c #AO      if (LSmask(ikl) .eq. 0 .and. coupling_ao .eq. .true.) then
422c #AO       if (AOmask(ikl) .eq. 0) then 
423c #AO       albisv(ikl) =  (1.-AOmask(ikl))* albAOsisv(ikl)
424c #AO.                    +(AOmask(ikl)*albisv(ikl))
425c #AO       alb1sv(ikl) =   (1.-AOmask(ikl))* albAOsisv(ikl)
426c #AO.                    +(AOmask(ikl)*alb1sv(ikl))
427c #AO       alb2sv(ikl) =   (1.-AOmask(ikl))* albAOsisv(ikl)
428c #AO.                    +(AOmask(ikl)*alb2sv(ikl))
429c #AO       alb3sv(ikl) =   (1.-AOmask(ikl))* albAOsisv(ikl)
430c #AO.                    +(AOmask(ikl)*alb3sv(ikl))
431c #AO       endif
432c #AO      endif
433
434 
435            alb1sv(ikl) = min(max(zero,alb1sv(ikl)),albmax)
436            alb2sv(ikl) = min(max(zero,alb2sv(ikl)),albmax)
437            alb3sv(ikl) = min(max(zero,alb3sv(ikl)),albmax)
438 
439        END DO
440 
441 
442C +--Extinction Coefficient: Exponential Factor
443C +  ==========================================
444 
445        DO ikl=1,knonv
446          sExt_1(ikl)        = 1.
447          sExt_2(ikl)        = 1.
448          sExt_3(ikl)        = 1.
449          sEX_sv(ikl,nsno+1) = 1.
450 
451          coalb1(ikl) = (1.          -alb1sv(ikl))*So1dSV
452          coalb2(ikl) = (1.          -alb2sv(ikl))*So2dSV
453          coalb3(ikl) = (1.          -alb3sv(ikl))*So3dSV
454          coalbm      =  coalb1(ikl) +coalb2(ikl) +coalb3(ikl)
455          coalb1(ikl) =  coalb1(ikl)              /coalbm
456          coalb2(ikl) =  coalb2(ikl)              /coalbm
457          coalb3(ikl) =  coalb3(ikl)              /coalbm
458        END DO
459 
460cXF
461 
462        DO   isn=nsno,1,-1
463          DO ikl=1,knonv
464            sEX_sv(ikl,isn) = 1.0
465           !sEX_sv(ikl,isn) = 0.95 ! if MAR is too warm in summer
466          END DO
467        END DO
468 
469        DO ikl=1,knonv
470         DO isn=max(1,isnoSV(ikl)),1,-1
471 
472          SignRo = sign(unun, rocdSV - ro__SV(ikl,isn))
473          SnowOK =  max(zero,SignRo) ! Ice Density Threshold
474 
475          RoFrez =  1.e-3      * ro__SV(ikl,isn) * (1.0-eta_SV(ikl,isn))
476 
477          OpSqrt = sqrt(max(epsi,SnOpSV(ikl,isn)))
478          exarg1 =      SnowOK  *1.e2 *max(sbeta1*RoFrez/OpSqrt,sbeta2)
479     .            +(1.0-SnowOK)           *sbeta5
480          exarg2 =      SnowOK  *1.e2 *max(sbeta3*RoFrez/OpSqrt,sbeta4)
481     .            +(1.0-SnowOK)           *sbeta5
482          exarg3 =      SnowOK  *1.e2     *sbeta5
483     .            +(1.0-SnowOK)           *sbeta5
484 
485 
486C +--Integrated Extinction of Solar Irradiance (Normalized Value)
487C +  ============================================================
488 
489          sExt_1(ikl) = sExt_1(ikl)
490     .                          * exp(min(0.0,-exarg1 *dzsnSV(ikl,isn)))
491          sign_0      =              sign(unun,eps9   -sExt_1(ikl))
492          sExt_0      =               max(zero,sign_0)*sExt_1(ikl)
493          sExt_1(ikl) = sExt_1(ikl)                   -sExt_0
494 
495          sExt_2(ikl) = sExt_2(ikl)
496     .                          * exp(min(0.0,-exarg2 *dzsnSV(ikl,isn)))
497          sign_0      =              sign(unun,eps9   -sExt_2(ikl))
498          sExt_0      =               max(zero,sign_0)*sExt_2(ikl)
499          sExt_2(ikl) = sExt_2(ikl)                   -sExt_0
500 
501          sExt_3(ikl) = sExt_3(ikl)
502     .                          * exp(min(0.0,-exarg3 *dzsnSV(ikl,isn)))
503          sign_0      =              sign(unun,eps9   -sExt_3(ikl))
504          sExt_0      =               max(zero,sign_0)*sExt_3(ikl)
505          sExt_3(ikl) = sExt_3(ikl)                   -sExt_0
506 
507          sEX_sv(ikl,isn) = coalb1(ikl) * sExt_1(ikl)
508     .                    + coalb2(ikl) * sExt_2(ikl)
509     .                    + coalb3(ikl) * sExt_3(ikl)
510        END DO
511      END DO
512 
513      DO   isn=0,-nsol,-1
514        DO ikl=1,knonv
515          sEX_sv(ikl,isn) = 0.0
516        END DO
517      END DO
518 
519 
520      return
521      end
Note: See TracBrowser for help on using the repository browser.