source: LMDZ.3.3/branches/rel-LF/libf/phylmd/albedo.F @ 1546

Last change on this file since 1546 was 580, checked in by Laurent Fairhead, 20 years ago

Valeur de fmagic pour les runs IPCC
LF

  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Author Date Id Revision
File size: 5.4 KB
Line 
1c
2c $Header$
3c
4      SUBROUTINE alboc(rjour,rlat,albedo)
5      IMPLICIT none
6c======================================================================
7c Auteur(s): Z.X. Li (LMD/CNRS) (adaptation du GCM du LMD)
8c Date: le 16 mars 1995
9c Objet: Calculer l'albedo sur l'ocean
10c Methode: Integrer numeriquement l'albedo pendant une journee
11c
12c Arguments;
13c rjour (in,R)  : jour dans l'annee (a compter du 1 janvier)
14c rlat (in,R)   : latitude en degre
15c albedo (out,R): albedo obtenu (de 0 a 1)
16c======================================================================
17#include "dimensions.h"
18#include "dimphy.h"
19#include "YOMCST.h"
20c
21      REAL fmagic ! un facteur magique pour regler l'albedo
22ccc      PARAMETER (fmagic=0.7)
23cccIM => a remplacer 
24c        PARAMETER (fmagic=1.32)
25       PARAMETER (fmagic=1.0)
26c       PARAMETER (fmagic=0.7)
27      INTEGER npts ! il controle la precision de l'integration
28      PARAMETER (npts=120) ! 120 correspond a l'interval 6 minutes
29c
30      REAL rlat(klon), rjour, albedo(klon)
31      REAL zdist, zlonsun, zpi, zdeclin
32      REAL rmu,alb, srmu, salb, fauxo, aa, bb
33      INTEGER i, k
34cccIM
35      LOGICAL ancien_albedo
36      PARAMETER(ancien_albedo=.FALSE.)
37c     SAVE albedo
38c
39      IF ( ancien_albedo ) THEN
40c
41      zpi = 4. * ATAN(1.)
42c
43c Calculer la longitude vraie de l'orbite terrestre:
44      CALL orbite(rjour,zlonsun,zdist)
45c
46c Calculer la declinaison du soleil (qui varie entre + et - R_incl):
47      zdeclin = ASIN(SIN(zlonsun*zpi/180.0)*SIN(R_incl*zpi/180.0))
48c
49      DO 999 i=1,klon
50      aa = SIN(rlat(i)*zpi/180.0) * SIN(zdeclin)
51      bb = COS(rlat(i)*zpi/180.0) * COS(zdeclin)
52c
53c Midi local (angle du temps = 0.0):
54      rmu = aa + bb * COS(0.0)
55      rmu = MAX(0.0, rmu)
56      fauxo = (1.47-ACOS(rmu))/.15
57      alb = 0.03+0.630/(1.+fauxo*fauxo)
58      srmu = rmu
59      salb = alb * rmu
60c
61c Faire l'integration numerique de midi a minuit (le facteur 2
62c prend en compte l'autre moitie de la journee):
63      DO k = 1, npts
64         rmu = aa + bb * COS(FLOAT(k)/FLOAT(npts)*zpi)
65         rmu = MAX(0.0, rmu)
66         fauxo = (1.47-ACOS(rmu))/.15
67         alb = 0.03+0.630/(1.+fauxo*fauxo)
68         srmu = srmu + rmu * 2.0
69         salb = salb + alb*rmu * 2.0
70      ENDDO
71      IF (srmu .NE. 0.0) THEN
72         albedo(i) = salb / srmu * fmagic
73      ELSE ! nuit polaire (on peut prendre une valeur quelconque)
74         albedo(i) = fmagic
75      ENDIF
76  999 CONTINUE
77c
78c nouvel albedo
79c
80      ELSE
81c
82      zpi = 4. * ATAN(1.)
83c
84c Calculer la longitude vraie de l'orbite terrestre:
85      CALL orbite(rjour,zlonsun,zdist)
86c
87c Calculer la declinaison du soleil (qui varie entre + et - R_incl):
88      zdeclin = ASIN(SIN(zlonsun*zpi/180.0)*SIN(R_incl*zpi/180.0))
89c
90      DO 1999 i=1,klon
91      aa = SIN(rlat(i)*zpi/180.0) * SIN(zdeclin)
92      bb = COS(rlat(i)*zpi/180.0) * COS(zdeclin)
93c
94c Midi local (angle du temps = 0.0):
95      rmu = aa + bb * COS(0.0)
96      rmu = MAX(0.0, rmu)
97cIM cf. PB  alb = 0.058/(rmu + 0.30)
98c     alb = 0.058/(rmu + 0.30) * 1.5
99      alb = 0.058/(rmu + 0.30) * 1.2
100c     alb = 0.058/(rmu + 0.30) * 1.3
101      srmu = rmu
102      salb = alb * rmu
103c
104c Faire l'integration numerique de midi a minuit (le facteur 2
105c prend en compte l'autre moitie de la journee):
106      DO k = 1, npts
107         rmu = aa + bb * COS(FLOAT(k)/FLOAT(npts)*zpi)
108         rmu = MAX(0.0, rmu)
109cIM cf. PB      alb = 0.058/(rmu + 0.30)
110c        alb = 0.058/(rmu + 0.30) * 1.5
111         alb = 0.058/(rmu + 0.30) * 1.2
112c        alb = 0.058/(rmu + 0.30) * 1.3
113         srmu = srmu + rmu * 2.0
114         salb = salb + alb*rmu * 2.0
115      ENDDO
116      IF (srmu .NE. 0.0) THEN
117         albedo(i) = salb / srmu * fmagic
118      ELSE ! nuit polaire (on peut prendre une valeur quelconque)
119         albedo(i) = fmagic
120      ENDIF
1211999  CONTINUE
122      ENDIF
123      RETURN
124      END
125c=====================================================================
126      SUBROUTINE alboc_cd(rmu0,albedo)
127      IMPLICIT none
128c======================================================================
129c Auteur(s): Z.X. Li (LMD/CNRS)
130c date: 19940624
131c Calculer l'albedo sur l'ocean en fonction de l'angle zenithal moyen
132c Formule due a Larson and Barkstrom (1977) Proc. of the symposium
133C on radiation in the atmosphere, 19-28 August 1976, science Press,
134C 1977 pp 451-453, ou These de 3eme cycle de Sylvie Joussaume.
135c
136c Arguments
137c rmu0    (in): cosinus de l'angle solaire zenithal
138c albedo (out): albedo de surface de l'ocean
139c======================================================================
140#include "dimensions.h"
141#include "dimphy.h"
142      REAL rmu0(klon), albedo(klon)
143c
144      REAL fmagic ! un facteur magique pour regler l'albedo
145ccc      PARAMETER (fmagic=0.7)
146cccIM => a remplacer 
147c        PARAMETER (fmagic=1.32)
148       PARAMETER (fmagic=1.0)
149c       PARAMETER (fmagic=0.7)
150c
151      REAL fauxo
152      INTEGER i
153cccIM
154      LOGICAL ancien_albedo
155      PARAMETER(ancien_albedo=.FALSE.)
156c     SAVE albedo
157c
158      IF ( ancien_albedo ) THEN
159c
160      DO i = 1, klon
161c
162         rmu0(i) = MAX(rmu0(i),0.0)
163c
164         fauxo = ( 1.47 - ACOS( rmu0(i) ) )/0.15
165         albedo(i) = fmagic*( .03 + .630/( 1. + fauxo*fauxo))
166         albedo(i) = MAX(MIN(albedo(i),0.60),0.04)
167      ENDDO
168c
169c nouvel albedo
170c
171      ELSE
172c
173      DO i = 1, klon
174         rmu0(i) = MAX(rmu0(i),0.0)
175cIM:orig albedo(i) = 0.058/(rmu0(i) + 0.30)
176         albedo(i) = fmagic * 0.058/(rmu0(i) + 0.30)
177         albedo(i) = MAX(MIN(albedo(i),0.60),0.04)
178      ENDDO
179c
180      ENDIF
181c
182      RETURN
183      END
184c========================================================================
Note: See TracBrowser for help on using the repository browser.