source: LMDZ.3.3/trunk/libf/dyn3d/modif_etat0.F @ 3651

Last change on this file since 3651 was 96, checked in by lmdzadmin, 24 years ago

Modification d'un etat0 physique pour rajouter un masque. P. Braconnot
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 8.1 KB
Line 
1      PROGRAM modif_etat0
2     
3      USE startvar
4      USE ioipsl
5      !
6      IMPLICIT NONE
7      !
8#include "dimensions.h"
9#include "paramet.h"
10      !
11      !
12#include "dimphy.h"
13      !
14#include "comgeom2.h"
15#include "comvert.h"
16#include "comconst.h"
17#include "indicesol.h"
18#include "dimsoil.h"
19      !
20      REAL :: latfi(klon), lonfi(klon)
21      REAL :: orog(iip1,jjp1), rugo(iip1,jjp1), masque(iip1,jjp1),
22     . psol(iip1, jjp1), phis(iip1, jjp1)
23      REAL :: tsol(klon,nbsrf), qsol(klon,nbsrf), sn(klon,nbsrf)
24      REAL :: evap(klon,nbsrf), albe(klon,nbsrf)
25      real :: rain_fall(klon), snow_fall(klon)
26      real :: sollw(klon), solsw(klon)
27      REAL :: radsol(klon),deltat(klon), rugmer(klon), agesno(klon)
28      REAL :: zmea(iip1*jjp1), zstd(iip1*jjp1)
29      REAL :: zsig(iip1*jjp1), zgam(iip1*jjp1), zthe(iip1*jjp1)
30      REAL :: zpic(iip1*jjp1), zval(iip1*jjp1), rugsrel(iip1*jjp1)
31      REAL :: pctsrf(klon, nbsrf)
32      REAL :: tsoil(klon, nsoilmx, nbsrf)
33      REAL :: t_ancien(klon, klev), q_ancien(klon, klev)
34      REAL :: tabcntrl0(100)
35      LOGICAL :: ancien_ok
36      ! declarations pour lecture glace de mer
37      INTEGER :: iml_lic, jml_lic, llm_tmp, ttm_tmp, iret
38      INTEGER :: itau(1), fid
39      REAL :: lev(1), date, dt
40      REAL, ALLOCATABLE, DIMENSION(:,:) :: lon_lic, lat_lic
41      REAL, ALLOCATABLE, DIMENSION(:)  :: dlon_lic, dlat_lic
42      REAL, ALLOCATABLE, DIMENSION (:,:) :: fraclic
43      REAL :: flic_tmp(iip1, jjp1)
44      REAL :: champint(iim, jjp1)
45      !
46      CHARACTER*80 :: varname
47      !
48#include "comdissnew.h"
49#include "control.h"
50#include "serre.h"
51#include "clesph0.h"
52
53      INTEGER  ::        longcles
54      PARAMETER      ( longcles  = 20 )
55      REAL :: clesphy0 ( longcles       )
56      REAL ::phystep,co2_ppm,solaire
57      REAL :: unskap
58      INTEGER :: radpas
59
60      CHARACTER*80 :: visu_file
61      INTEGER :: visuid
62      !
63      ! indices de boucle
64      INTEGER :: ji
65      !   Constantes
66      !
67      pi     = 4. * ATAN(1.)
68      rad    = 6371229.
69      omeg   = 4.* ASIN(1.)/(24.*3600.)
70      g      = 9.8
71      daysec = 86400.
72      kappa  = 0.2857143
73      cpp    = 1004.70885
74      !
75      preff     = 101325.
76      unskap = 1./kappa
77      !
78      jmp1    = jjm + 1
79      !
80      !    Construct a grid
81      !
82
83      CALL defrun_new(99,.TRUE.,clesphy0)
84
85      dtvr   = daysec/FLOAT(day_step)
86      print*,'dtvr',dtvr
87
88      CALL inicons0()
89      CALL inigeom()
90      !
91      CALL inifilr()
92      !
93      !
94      !! 1. READ orography
95      !
96      varname = 'relief'
97      ! This line needs to be replaced by a call to restget to get the values in the restart file
98      orog(:,:) = 0.0
99       CALL startget(varname, iip1, jjp1, rlonv, rlatu, orog, 0.0)
100      !
101      WRITE(*,*) 'OUT OF GET VARIABLE : Relief'
102      WRITE(*,'(49I1)') INT(orog(:,:))
103      !
104      varname = 'rugosite'
105      ! This line needs to be replaced by a call to restget to get the values in the restart file
106      rugo(:,:) = 0.0
107       CALL startget(varname, iip1, jjp1, rlonv, rlatu, rugo, 0.0)
108      !
109      WRITE(*,*) 'OUT OF GET VARIABLE : Rugosite'
110      WRITE(*,'(49I1)') INT(rugo(:,:)*10)
111      !
112      varname = 'masque'
113      ! This line needs to be replaced by a call to restget to get the values in the restart file
114      masque(:,:) = 0.0
115       CALL startget(varname, iip1, jjp1, rlonv, rlatu, masque, 0.0)
116      !
117      WRITE(*,*) 'MASQUE construit : Masque'
118      WRITE(*,*) masque(:,:)
119      WRITE(*,*) 'orographie'
120      WRITE(*,*) orog(:,:)
121      !
122      !! 2. create masque and READ land iice fraction
123      !!
124      ! lit bande startphy a modifier
125      !
126     
127      CALL phyetat0("startphy.nc", phystep, co2_ppm, solaire, latfi, 
128     $    lonfi, pctsrf, tsol, tsoil, deltat, qsol, sn,
129     $    albe, evap, rain_fall, snow_fall, solsw, sollw,
130     $    radsol, rugmer, agesno, clesphy0,
131     $    zmea, zstd, zsig, zgam, zthe, zpic, zval,
132     $    rugsrel, tabcntrl0, t_ancien, q_ancien, ancien_ok)
133C
134C on initialise les sous surfaces
135C
136      pctsrf=0.
137      !cree le masque a partir du fichier relief
138      varname = 'zmasq'
139      zmasq(:) = 0.
140      CALL startget(varname,iip1,jjp1,rlonv,rlatu,klon,zmasq,0.0)
141      WHERE (zmasq(1 : klon) .LE. EPSFRA)
142          zmasq(1 : klon) = 0.
143      END WHERE
144      WRITE(*,*)zmasq
145      varname = 'zmea'
146      zmea(:) = 0.0
147      CALL startget(varname,iip1,jjp1,rlonv,rlatu,klon,zmea,0.0)
148      varname = 'zstd'
149      zstd(:) = 0.0
150      CALL startget(varname,iip1,jjp1,rlonv,rlatu,klon,zstd,0.0)
151      varname = 'zsig'
152      zsig(:) = 0.0
153      CALL startget(varname,iip1,jjp1,rlonv,rlatu,klon,zsig,0.0)
154      varname = 'zgam'
155      zgam(:) = 0.0
156      CALL startget(varname,iip1,jjp1,rlonv,rlatu,klon,zgam,0.0)
157      varname = 'zthe'
158      zthe(:) = 0.0
159      CALL startget(varname,iip1,jjp1,rlonv,rlatu,klon,zthe,0.0)
160      varname = 'zpic'
161      zpic(:) = 0.0
162      CALL startget(varname,iip1,jjp1,rlonv,rlatu,klon,zpic,0.0)
163      varname = 'zval'
164      zval(:) = 0.0
165      CALL startget(varname,iip1,jjp1,rlonv,rlatu,klon,zval,0.0)
166      rugsrel(:) = 0.0
167C
168C lecture du fichier glace de terre pour fixer la fraction de terre
169C et de glace de terre
170C
171      CALL flininfo("landiceref.nc", iml_lic, jml_lic,llm_tmp, ttm_tmp
172     $    , fid)
173      ALLOCATE(lat_lic(iml_lic, jml_lic), stat=iret)
174      ALLOCATE(lon_lic(iml_lic, jml_lic), stat=iret)
175      ALLOCATE(dlon_lic(iml_lic), stat=iret)
176      ALLOCATE(dlat_lic(jml_lic), stat=iret)
177      ALLOCATE(fraclic(iml_lic, jml_lic), stat=iret)
178      CALL flinopen("landiceref.nc", .FALSE., iml_lic, jml_lic, llm_tmp
179     $    , lon_lic, lat_lic, lev, ttm_tmp, itau, date, dt, fid)
180      CALL flinget(fid, 'landice', iml_lic, jml_lic, llm_tmp, ttm_tmp
181     $    , 1, 1, fraclic)
182      CALL flinclo(fid)
183C
184C interpolation sur la grille T du modele
185C
186      WRITE(*,*) 'dimensions de landice iml_lic, jml_lic : ',
187     $    iml_lic, jml_lic
188c
189C sil les coordonnees sont en degres, on les transforme
190C
191      IF( MAXVAL( lon_lic(:,:) ) .GT. 2.0 * asin(1.0) )  THEN
192          lon_lic(:,:) = lon_lic(:,:) * 2.0* ASIN(1.0) / 180.
193      ENDIF
194      IF( maxval( lat_lic(:,:) ) .GT. 2.0 * asin(1.0)) THEN
195          lat_lic(:,:) = lat_lic(:,:) * 2.0 * asin(1.0) / 180.
196      ENDIF
197
198      dlon_lic(1 : iml_lic) = lon_lic(1 : iml_lic, 1)
199      dlat_lic(1 : jml_lic) = lat_lic(1 , 1 : jml_lic)
200C
201      CALL grille_m(iml_lic, jml_lic, dlon_lic, dlat_lic, fraclic
202     $    ,iim, jjp1,
203     $    rlonv, rlatu, flic_tmp(1 : iim, 1 : jjp1))
204c$$$      flic_tmp(1 : iim, 1 : jjp1) = champint(1: iim, 1 : jjp1)
205      flic_tmp(iip1, 1 : jjp1) = flic_tmp(1 , 1 : jjp1)
206C
207C passage sur la grille physique
208C
209      CALL gr_dyn_fi(1, iip1, jjp1, klon, flic_tmp,
210     $    pctsrf(1:klon, is_lic))
211C adequation avec le maque terre/mer
212      WHERE (pctsrf(1 : klon, is_lic) .LE. EPSFRA )
213          pctsrf(1 : klon, is_lic) = 0.
214      END WHERE
215      WHERE (zmasq( 1 : klon) .LE. EPSFRA)
216          pctsrf(1 : klon, is_lic) = 0.
217      END WHERE
218      pctsrf(1 : klon, is_ter) = zmasq(1 : klon)
219      DO ji = 1, klon
220        IF (zmasq(ji) .GT. EPSFRA) THEN
221            IF ( pctsrf(ji, is_lic) .GE. zmasq(ji)) THEN
222                pctsrf(ji, is_lic) = zmasq(ji)
223                pctsrf(ji, is_ter) = 0.
224            ELSE
225                pctsrf(ji,is_ter) = zmasq(ji) - pctsrf(ji, is_lic)
226            ENDIF
227        ENDIF
228      END DO
229C
230C sous surface ocean et glace de mer (pour demarrer on met glace de mer a 0)
231C
232      pctsrf(1 : klon, is_oce) = (1. - zmasq(1 : klon))
233      WHERE (pctsrf(1 : klon, is_oce) .LT. EPSFRA)
234          pctsrf(1 : klon, is_oce) = 0.
235      END WHERE
236C
237C verif que somme des sous surface = 1
238C
239      ji=count( (abs( sum(pctsrf(1 : klon, 1 : nbsrf), dim = 2)) - 1.0 )
240     $    .GT. EPSFRA)
241      IF (ji .NE. 0) THEN
242          WRITE(*,*) 'pb repartition sous maille pour ',ji,' points'
243      ENDIF
244C
245C Ecriture etat initial physique
246C
247      WRITE(*,*) 'zmasq avant phyredem'
248      WRITE(*,*) zmasq
249
250      call phyredem("restartphy.nc",phystep,radpas, co2_ppm, solaire,
251     $    latfi, lonfi, pctsrf, tsol, tsoil, deltat, qsol, sn,
252     $    albe, evap, rain_fall, snow_fall, solsw, sollw,
253     $    radsol, rugmer,  agesno,
254     $    zmea, zstd, zsig, zgam, zthe, zpic, zval, rugsrel,
255     $    t_ancien, q_ancien)
256
257      CALL histclo
258      !
259      END 
Note: See TracBrowser for help on using the repository browser.