source: LMDZ5/branches/IPSLCM5A2.1_ISO/libf/phylmd/ocean_albedo.F90 @ 5396

Last change on this file since 5396 was 2227, checked in by Laurent Fairhead, 10 years ago

New ocean albedo.

To activate the new scheme, put iflag_albedo=1 in physiq.def

To activate chlorophyll concentration effect on albedo,
put ok_chlorophyll=y in def file

and download file named chlorophyll.nc
chlorophyll.nc has the same dimension as the model grid with 12 months data,
(i=lon, j=lat, L=1:12) and can be degraded from the original file of dimension
i=1:4320 , j=1:2160 , L=1:12
ada:/workgpfs/rech/psl/rpsl949/clima/chlor_seasonal_clim_seawifs.nc

For 96X96 resolution, chlorophyll.nc file is in
ada:/workgpfs/rech/psl/rpsl949/clima/chlorophyll.nc

  1. Baek
File size: 10.5 KB
Line 
1!
2! $Id$
3!
4
5!     #########
6
7subroutine ocean_albedo(knon,zrmu0,knindex,pwind,SFRWL,alb_dir_new,alb_dif_new)
8
9
10!     ##################################################################
11!
12!!****  *ALBEDO_RS14* 
13!!
14!!    PURPOSE
15!!    -------
16!       computes the direct & diffuse albedo over open water
17!
18!     
19!!**  METHOD
20!!    ------
21!
22!!    EXTERNAL
23!!    --------
24!!
25!!    IMPLICIT ARGUMENTS
26!!    ------------------
27!!
28!!     
29!!    REFERENCE
30!!    ---------
31!!
32!!     
33!!    AUTHOR
34!!    ------
35!!      R. Séférian           * Meteo-France *
36!!
37!!    MODIFICATIONS
38!!    -------------
39!!      Original    03/2014
40!                   05/2014 R. Séférian & B. Decharme :: Adaptation to spectral
41!                   computation for diffuse and direct albedo
42!                   08/2014 S. Baek :: for wider wavelength range 200-4000nm and
43!                   adaptation to LMDZ + whitecap effect by Koepke + chrolophyll
44!                   map from climatology file
45!       
46!-------------------------------------------------------------------------------
47!
48!*           DECLARATIONS
49!            ------------
50!
51USE ocean_albedo_para
52use dimphy
53!LF USE PARKIND1  ,ONLY : JPIM     ,JPRB
54use phys_state_var_mod, only : chl_con
55!
56!
57IMPLICIT NONE
58!
59!*      0.1    declarations of arguments
60!              -------------------------
61!
62
63include "clesphys.h"
64
65integer, intent(in) :: knon
66integer, dimension(klon), intent(in) :: knindex
67real, dimension(klon), intent(in) :: zrmu0,pwind
68real, dimension(klon,nsw), intent(out) ::  alb_dir_new,alb_dif_new
69real, dimension(6),intent(in) :: SFRWL
70
71
72!=== LOCAL VARIABLES
73
74REAL, parameter :: XPI=4.*atan(1.)
75
76!
77!*      0.2    declarations of local variables
78!              -------------------------
79!
80REAL, DIMENSION(klon)               :: ZCHL        ! surface chlorophyll
81REAL, DIMENSION(klon,NNWL)               :: ZDIR_ALB    ! direct  ocean surface albedo (spectral)
82REAL, DIMENSION(klon,NNWL)               :: ZSCA_ALB    ! diffuse ocean surface albedo (spectral)
83!
84INTEGER                         :: JI, JWL                  ! indexes
85REAL                            :: ZWL                      ! input parameter: wavelength and diffuse/direct fraction of light
86REAL:: ZSIG, ZREFM, ZXX2, ZR00, ZRR0, ZRRR                        ! computation variables
87REAL:: ZR22, ZUE, ZUE2, ZR11DF, ZALBT, ZFWC                       ! computation variables
88REAL:: ZCHLABS, ZAW, ZBW, ZAP, ZYLMD, ZBP550                      ! computation variables
89REAL:: ZBBP, ZNU, ZHB                                       ! computation variables
90REAL:: ZCOSZEN                                              ! Cosine of the zenith solar angle
91REAL:: ZR11, ZRW, ZRWDF, ZRDF                                   ! 4 components of the OSA
92! new damping coefficient
93REAL:: ZDAMP     
94
95!
96REAL            :: ZWORK                           ! dummy variable
97!
98!LF REAL(KIND=JPRB) :: ZHOOK_HANDLE
99!
100!-------------------------------------------------------------------------------
101!
102!
103
104
105
106!
107alb_dir_new(:,:) = 0.
108alb_dif_new(:,:) = 0.
109!
110ZDIR_ALB(:,:) = 0.
111ZSCA_ALB(:,:) = 0.
112!
113!
114
115!ZCHL(:) = CHL_CON!0.05 ! averaged global values for surface chlorophyll
116if(ok_chlorophyll)then
117  do ji=1,knon
118  ZCHL(ji)=CHL_CON(knindex(ji))
119  enddo
120else
121  ZCHL(:) = 0.05
122endif
123
124
125!
126DO JWL=1,NNWL           ! loop over the wavelength
127!
128  DO JI=1,knon  ! loop over the grid points
129
130
131    !---------------------------------------------------------------------------------
132    ! 0- Compute baseline values
133    !---------------------------------------------------------------------------------
134   
135    ! Get refractive index for the correspoding wavelength
136    ZWL=XAKWL(JWL) !!!----------  wavelength value
137    ZREFM= XAKREFM(JWL) !!!--------- refraction index value
138 
139
140    ! compute the cosine of the solar zenith angle
141!    ZCOSZEN = COS(XPI/2 - PZENITH(JI))
142     ZCOSZEN = zrmu0(knindex(JI))
143    ! Compute sigma derived from wind speed (Cox & Munk reflectance model)
144    ZSIG=SQRT(0.003+0.00512*PWIND(JI))
145
146   
147    !---------------------------------------------------------------------------------
148    ! 1- Compute direct surface albedo (ZR11)
149    !---------------------------------------------------------------------------------
150    !
151    ZXX2=SQRT(1.0-(1.0-ZCOSZEN**2)/ZREFM**2)
152    ZRR0=0.50*(((ZXX2-ZREFM*ZCOSZEN)/(ZXX2+ZREFM*ZCOSZEN))**2 +((ZCOSZEN-ZREFM*ZXX2)/(ZCOSZEN+ZREFM*ZXX2))**2)
153    ZRRR=0.50*(((ZXX2-1.34*ZCOSZEN)/(ZXX2+1.34*ZCOSZEN))**2 +((ZCOSZEN-1.34*ZXX2)/(ZCOSZEN+1.34*ZXX2))**2)
154    ZR11=ZRR0-(0.0152-1.7873*ZCOSZEN+6.8972*ZCOSZEN**2-8.5778*ZCOSZEN**3+4.071*ZSIG-7.6446*ZCOSZEN*ZSIG) &
155        & * EXP(0.1643-7.8409*ZCOSZEN-3.5639*ZCOSZEN**2-2.3588*ZSIG+10.0538*ZCOSZEN*ZSIG)*ZRR0/ZRRR
156    !
157    !---------------------------------------------------------------------------------
158    ! 2- Compute surface diffuse albedo (ZRDF)
159    !---------------------------------------------------------------------------------
160    ! Diffuse albedo from Jin et al., 2006 + estimation from diffuse fraction of
161    ! light (relying later on AOD)
162    ZRDF=-0.1482-0.012*ZSIG+0.1609*ZREFM-0.0244*ZSIG*ZREFM ! surface diffuse (Eq 5a-5b)
163   
164    !---------------------------------------------------------------------------------
165    ! *- Determine absorption and backscattering
166    ! coefficients to determine reflectance below the surface (Ro) once for all
167    !
168    ! *.1- Absorption by chlorophyll
169    ZCHLABS= XAKACHL(JWL)
170    ! *.2- Absorption by seawater
171    ZAW= XAKAW3(JWL)
172    ! *.3- Backscattering by seawater
173    ZBW= XAKBW(JWL)
174    ! *.4- Backscattering by chlorophyll
175    ZYLMD = EXP(0.014*(440.0-ZWL))
176    ZWORK= EXP(LOG(ZCHL(JI))*0.65)
177    ZAP = 0.06*ZCHLABS*ZWORK +0.2*(XAW440+0.06*ZWORK)*ZYLMD
178    ZBP550 = 0.416 * EXP(LOG(ZCHL(JI))*0.766)
179   
180    IF ( ZCHL(JI) > 2. ) THEN
181      ZNU=0.
182    ELSE
183      IF ( ZCHL(JI) > 0.02 ) THEN
184        ZWORK=LOG10(ZCHL(JI))
185        ZNU=0.5*(ZWORK-0.3)
186        ZBBP=(0.002+0.01*(0.5-0.25*ZWORK)*(ZWL/550.)**ZNU)*ZBP550
187      ELSE
188        ZBBP=0.019*(550./ZWL)*ZBP550       !ZBBPf=0.0113 at chl<=0.02
189      ENDIF
190    ENDIF
191   
192    ! Morel-Gentili(1991), Eq (12)
193    ! ZHB=h/(h+2*ZBBPf*(1.-h))       
194    ZHB=0.5*ZBW/(0.5*ZBW+ZBBP)
195   
196    !---------------------------------------------------------------------------------
197    ! 3- Compute direct water-leaving albedo (ZRW)
198    !---------------------------------------------------------------------------------
199    ! Based on Morel & Gentilli 1991 parametrization
200    ZR22=0.48168549-0.014894708*ZSIG-0.20703885*ZSIG**2
201    ! Use Morel 91 formula to compute the direct reflectance
202    ! below the surface
203    ZR00=(0.5*ZBW+ZBBP)/(ZAW+ZAP) *(0.6279-0.2227*ZHB-0.0513*ZHB**2 + (-0.3119+0.2465*ZHB)*ZCOSZEN)
204    ZRW=ZR00*(1.-ZR22)*(1.-ZR11)/(1.-ZR00*ZR22)
205
206    ZRW=ZR00*(1.-ZR22)/(1.-ZR00*ZR22)
207    !---------------------------------------------------------------------------------
208    ! 4- Compute diffuse water-leaving albedo (ZRWDF)
209    !---------------------------------------------------------------------------------
210    ! as previous water-leaving computation but assumes a uniform incidence of
211    ! shortwave at surface (ue)
212    ZUE=0.676               ! equivalent u_unif for diffuse incidence
213    ZUE2=SQRT(1.0-(1.0-ZUE**2)/ZREFM**2)
214    ZRR0=0.50*(((ZUE2-ZREFM*ZUE)/(ZUE2+ZREFM*ZUE))**2 +((ZUE-ZREFM*ZUE2)/(ZUE+ZREFM*ZUE2))**2)
215    ZRRR=0.50*(((ZUE2-1.34*ZUE)/(ZUE2+1.34*ZUE))**2 +((ZUE-1.34*ZUE2)/(ZUE+1.34*ZUE2))**2)
216    ZR11DF=ZRR0-(0.0152-1.7873*ZUE+6.8972*ZUE**2-8.5778*ZUE**3+4.071*ZSIG-7.6446*ZUE*ZSIG) &
217          & * EXP(0.1643-7.8409*ZUE-3.5639*ZUE**2-2.3588*ZSIG+10.0538*ZUE*ZSIG)*ZRR0/ZRRR
218    ! Use Morel 91 formula to compute the diffuse
219    ! reflectance below the surface
220    ZR00=(0.5*ZBW+ZBBP)/(ZAW+ZAP) *(0.6279-0.2227*ZHB-0.0513*ZHB**2 + (-0.3119+0.2465*ZHB)*ZUE)
221    ZRWDF=ZR00*(1.-ZR22)*(1.-ZR11DF)/(1.-ZR00*ZR22)
222   
223    ! original : correction for foam (Eq 16-17)
224    ZFWC=3.97e-4*PWIND(JI)**(1.59) ! Salisbury 2014 eq(2) at 37GHz, value in fraction
225    ! has to be update once we have information from wave model (discussion with G. Madec)
226     
227    ! --------------------------------------------------------------------
228    !  *- OSA estimation
229    ! --------------------------------------------------------------------
230    ! partitionning direct and diffuse albedo
231    !
232
233    ! excluding diffuse albedo ZRW on ZDIR_ALB
234      ZDIR_ALB(JI,JWL) =  XFRWL(JWL) *((1.-ZFWC) * (ZR11+ZRW) +ZFWC*XRWC(JWL))
235      ZSCA_ALB(JI,JWL) =  XFRWL(JWL) *((1.-ZFWC) * (ZRDF+ZRWDF) + ZFWC*XRWC(JWL))
236
237    !  print*,ji,ZFWC,ZDIR_ALB(JI,JWL),ZSCA_ALB(JI,JWL),pwind(ji)
238    ENDDO ! end of the loop over grid points
239
240ENDDO ! ending loop over wavelengths
241
242
243! integral for each nsw band
244
245select case(nsw)
246case(2)
247  do ji=1,knon
248   alb_dir_new(ji,1)=sum(zdir_alb(ji,1:49))/SFRWL(1) ! from 200nm to 680nm
249   alb_dir_new(ji,2)=sum(zdir_alb(ji,50:381))/SFRWL(2) ! from 690nm to 4000 nm
250
251   alb_dif_new(ji,1)=sum(zsca_alb(ji,1:49))/SFRWL(1) ! from 200nm to 680nm
252   alb_dif_new(ji,2)=sum(zsca_alb(ji,50:381))/SFRWL(2) ! from 690nm to 4000 nm
253  enddo
254case(4)
255  do ji=1,knon
256   alb_dir_new(ji,1)=sum(zdir_alb(ji,1:49))/SFRWL(1) ! from 200nm to 680nm
257   alb_dir_new(ji,2)=sum(zdir_alb(ji,50:99))/SFRWL(2) ! from 690nm to 1180 nm
258   alb_dir_new(ji,3)=sum(zdir_alb(ji,100:218))/SFRWL(3) ! from 1190nm to 2370 nm
259   alb_dir_new(ji,4)=sum(zdir_alb(ji,219:381))/SFRWL(4) ! from 2380nm to 4000 nm
260
261   alb_dif_new(ji,1)=sum(zsca_alb(ji,1:49))/SFRWL(1) ! from 200nm to 680nm
262   alb_dif_new(ji,2)=sum(zsca_alb(ji,50:99))/SFRWL(2) ! from 690nm to 1180 nm
263   alb_dif_new(ji,3)=sum(zsca_alb(ji,100:218))/SFRWL(3) ! from 1190nm to 2370 nm
264   alb_dif_new(ji,4)=sum(zsca_alb(ji,219:381))/SFRWL(4) ! from 2380nm to 4000 nm
265  enddo
266case(6)
267 do ji=1,knon
268   alb_dir_new(ji,1)=sum(zdir_alb(ji,1:5))/SFRWL(1) ! from 200nm to 240nm
269   alb_dir_new(ji,2)=sum(zdir_alb(ji,6:24))/SFRWL(2) ! from 250nm to 430 nm
270   alb_dir_new(ji,3)=sum(zdir_alb(ji,25:49))/SFRWL(3) ! from 440nm to 680 nm
271   alb_dir_new(ji,4)=sum(zdir_alb(ji,50:99))/SFRWL(4) ! from 690nm to 1180 nm
272   alb_dir_new(ji,5)=sum(zdir_alb(ji,100:218))/SFRWL(5) ! from 1190nm to 2370 nm
273   alb_dir_new(ji,6)=sum(zdir_alb(ji,219:381))/SFRWL(6) ! from 2380nm to 4000 nm
274
275   alb_dif_new(ji,1)=sum(zsca_alb(ji,1:5))/SFRWL(1) ! from 200nm to 240nm
276   alb_dif_new(ji,2)=sum(zsca_alb(ji,6:24))/SFRWL(2) ! from 250nm to 430 nm
277   alb_dif_new(ji,3)=sum(zsca_alb(ji,25:49))/SFRWL(3) ! from 440nm to 680 nm
278   alb_dif_new(ji,4)=sum(zsca_alb(ji,50:99))/SFRWL(4) ! from 690nm to 1180 nm
279   alb_dif_new(ji,5)=sum(zsca_alb(ji,100:218))/SFRWL(5) ! from 1190nm to 2370 nm
280   alb_dif_new(ji,6)=sum(zsca_alb(ji,219:381))/SFRWL(6) ! from 2380nm to 4000 nm
281 enddo
282end select
283
284
285
286END subroutine ocean_albedo
287
Note: See TracBrowser for help on using the repository browser.