Ignore:
Timestamp:
Aug 6, 2013, 3:33:18 PM (11 years ago)
Author:
lguez
Message:

Changed names of variables ema_work1 and ema_work2 to more meaningful
sig1 and w01. Same change in (re)startphy.nc. phyetat0 tries to find
old names ema_work1 and ema_work2 if new names sig1 and w01 are not
found, so the program can run with an old restartphy.nc. restartphy.nc
is modified compared to the previous SVN revision because of the change of
names but the data content is not modified (this can be checked with
max_diff_nc.sh -i).

File:
1 moved

Legend:

Unmodified
Added
Removed

  • LMDZ5/trunk/libf/phylmd/phyetat0.F90

    r1826 r1827  
    1 !
    21! $Id$
    3 !
    4 c
    5 c
    6       SUBROUTINE phyetat0 (fichnom,
    7      .           clesphy0,
    8      .           tabcntr0)
    9 
    10       USE dimphy
    11       USE mod_grid_phy_lmdz
    12       USE mod_phys_lmdz_para
    13       USE iophy
    14       USE ocean_cpl_mod,    ONLY : ocean_cpl_init
    15       USE fonte_neige_mod,  ONLY : fonte_neige_init
    16       USE pbl_surface_mod,  ONLY : pbl_surface_init
    17       USE surface_data,     ONLY : type_ocean
    18       USE phys_state_var_mod
    19       USE iostart
    20       USE write_field_phy
    21       USE infotrac
    22       USE traclmdz_mod,    ONLY : traclmdz_from_restart
    23       USE carbon_cycle_mod,ONLY :
    24      &     carbon_cycle_tr, carbon_cycle_cpl, co2_send
    25       USE indice_sol_mod
    26 
    27       IMPLICIT none
    28 c======================================================================
    29 c Auteur(s) Z.X. Li (LMD/CNRS) date: 19930818
    30 c Objet: Lecture de l'etat initial pour la physique
    31 c======================================================================
    32 #include "dimensions.h"
    33 #include "netcdf.inc"
    34 #include "dimsoil.h"
    35 #include "clesphys.h"
    36 #include "temps.h"
    37 #include "thermcell.h"
    38 #include "compbl.h"
    39 c======================================================================
    40       CHARACTER*(*) fichnom
    41 
    42 c les variables globales lues dans le fichier restart
    43 
    44       REAL tsoil(klon,nsoilmx,nbsrf)
    45       REAL tslab(klon), seaice(klon)
    46       REAL qsurf(klon,nbsrf)
    47       REAL qsol(klon)
    48       REAL snow(klon,nbsrf)
    49       REAL evap(klon,nbsrf)
    50       real fder(klon)
    51       REAL frugs(klon,nbsrf)
    52       REAL agesno(klon,nbsrf)
    53       REAL run_off_lic_0(klon)
    54       REAL fractint(klon)
    55       REAL trs(klon,nbtr)
    56 
    57       CHARACTER*6 ocean_in
    58       LOGICAL ok_veget_in
    59 
    60       INTEGER        longcles
    61       PARAMETER    ( longcles = 20 )
    62       REAL clesphy0( longcles )
    63 c
    64       REAL xmin, xmax
    65 c
    66       INTEGER nid, nvarid
    67       INTEGER ierr, i, nsrf, isoil ,k
    68       INTEGER length
    69       PARAMETER (length=100)
    70       INTEGER it, iiq
    71       REAL tab_cntrl(length), tabcntr0(length)
    72       CHARACTER*7 str7
    73       CHARACTER*2 str2
    74       LOGICAL :: found
    75 
    76 c FH1D
    77 c     real iolat(jjm+1)
    78       real iolat(jjm+1-1/(iim*jjm))
    79 c
    80 c Ouvrir le fichier contenant l'etat initial:
    81 c
    82 
    83      
    84       CALL open_startphy(fichnom)
    85      
    86 
    87 c
    88 c Lecture des parametres de controle:
    89 c
    90       CALL get_var("controle",tab_cntrl)
    91        
    92 c
     2
     3SUBROUTINE phyetat0 (fichnom, clesphy0, tabcntr0)
     4
     5  USE dimphy
     6  USE mod_grid_phy_lmdz
     7  USE mod_phys_lmdz_para
     8  USE iophy
     9  USE ocean_cpl_mod,    ONLY : ocean_cpl_init
     10  USE fonte_neige_mod,  ONLY : fonte_neige_init
     11  USE pbl_surface_mod,  ONLY : pbl_surface_init
     12  USE surface_data,     ONLY : type_ocean
     13  USE phys_state_var_mod
     14  USE iostart
     15  USE write_field_phy
     16  USE infotrac
     17  USE traclmdz_mod,    ONLY : traclmdz_from_restart
     18  USE carbon_cycle_mod, ONLY : carbon_cycle_tr, carbon_cycle_cpl, co2_send
     19  USE indice_sol_mod
     20
     21  IMPLICIT none
     22  !======================================================================
     23  ! Auteur(s) Z.X. Li (LMD/CNRS) date: 19930818
     24  ! Objet: Lecture de l'etat initial pour la physique
     25  !======================================================================
     26  include "dimensions.h"
     27  include "netcdf.inc"
     28  include "dimsoil.h"
     29  include "clesphys.h"
     30  include "temps.h"
     31  include "thermcell.h"
     32  include "compbl.h"
     33  !======================================================================
     34  CHARACTER*(*) fichnom
     35
     36  ! les variables globales lues dans le fichier restart
     37
     38  REAL tsoil(klon, nsoilmx, nbsrf)
     39  REAL tslab(klon), seaice(klon)
     40  REAL qsurf(klon, nbsrf)
     41  REAL qsol(klon)
     42  REAL snow(klon, nbsrf)
     43  REAL evap(klon, nbsrf)
     44  real fder(klon)
     45  REAL frugs(klon, nbsrf)
     46  REAL agesno(klon, nbsrf)
     47  REAL run_off_lic_0(klon)
     48  REAL fractint(klon)
     49  REAL trs(klon, nbtr)
     50
     51  CHARACTER*6 ocean_in
     52  LOGICAL ok_veget_in
     53
     54  INTEGER        longcles
     55  PARAMETER    ( longcles = 20 )
     56  REAL clesphy0( longcles )
     57
     58  REAL xmin, xmax
     59
     60  INTEGER nid, nvarid
     61  INTEGER ierr, i, nsrf, isoil , k
     62  INTEGER length
     63  PARAMETER (length=100)
     64  INTEGER it, iiq
     65  REAL tab_cntrl(length), tabcntr0(length)
     66  CHARACTER*7 str7
     67  CHARACTER*2 str2
     68  LOGICAL :: found
     69
     70  ! FH1D
     71  !     real iolat(jjm+1)
     72  real iolat(jjm+1-1/(iim*jjm))
     73
     74  ! Ouvrir le fichier contenant l'etat initial:
     75
     76  CALL open_startphy(fichnom)
     77
     78  ! Lecture des parametres de controle:
     79
     80  CALL get_var("controle", tab_cntrl)
     81
    9382!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    94 ! FH 2008/05/09 On elimine toutes les clefs physiques dans la dynamique
    95 ! Les constantes de la physiques sont lues dans la physique seulement.
    96 ! Les egalites du type
    97 !             tab_cntrl( 5 )=clesphy0(1)
    98 ! sont remplacees par
    99 !             clesphy0(1)=tab_cntrl( 5 )
    100 ! On inverse aussi la logique.
    101 ! On remplit les tab_cntrl avec les parametres lus dans les .def
     83  ! FH 2008/05/09 On elimine toutes les clefs physiques dans la dynamique
     84  ! Les constantes de la physiques sont lues dans la physique seulement.
     85  ! Les egalites du type
     86  !             tab_cntrl( 5 )=clesphy0(1)
     87  ! sont remplacees par
     88  !             clesphy0(1)=tab_cntrl( 5 )
     89  ! On inverse aussi la logique.
     90  ! On remplit les tab_cntrl avec les parametres lus dans les .def
    10291!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    10392
    104          DO i = 1, length
    105            tabcntr0( i ) = tab_cntrl( i )
    106          ENDDO
    107 c
    108          tab_cntrl(1)=dtime
    109          tab_cntrl(2)=radpas
    110 
    111 c co2_ppm : value from the previous time step
    112          IF (carbon_cycle_tr .OR. carbon_cycle_cpl) THEN
    113             co2_ppm = tab_cntrl(3)
    114             RCO2    = co2_ppm * 1.0e-06  * 44.011/28.97
    115 c ELSE : keep value from .def
    116          END IF
    117 
    118 c co2_ppm0 : initial value of atmospheric CO2 (from create_etat0_limit.e .def)
    119          co2_ppm0   = tab_cntrl(16)
    120 
    121          solaire_etat0      = tab_cntrl(4)
    122          tab_cntrl(5)=iflag_con
    123          tab_cntrl(6)=nbapp_rad
    124 
    125          if (cycle_diurne) tab_cntrl( 7) =1.
    126          if (soil_model) tab_cntrl( 8) =1.
    127          if (new_oliq) tab_cntrl( 9) =1.
    128          if (ok_orodr) tab_cntrl(10) =1.
    129          if (ok_orolf) tab_cntrl(11) =1.
    130          if (ok_limitvrai) tab_cntrl(12) =1.
    131 
    132 
    133       itau_phy = tab_cntrl(15)
    134 
    135        
    136       clesphy0(1)=tab_cntrl( 5 )
    137       clesphy0(2)=tab_cntrl( 6 )
    138       clesphy0(3)=tab_cntrl( 7 )
    139       clesphy0(4)=tab_cntrl( 8 )
    140       clesphy0(5)=tab_cntrl( 9 )
    141       clesphy0(6)=tab_cntrl( 10 )
    142       clesphy0(7)=tab_cntrl( 11 )
    143       clesphy0(8)=tab_cntrl( 12 )
    144 
    145 c
    146 c Lecture des latitudes (coordonnees):
    147 c
    148       CALL get_field("latitude",rlat)
    149 
    150 c
    151 c Lecture des longitudes (coordonnees):
    152 c
    153       CALL get_field("longitude",rlon)
    154 
    155 C
    156 C
    157 C Lecture du masque terre mer
    158 C
    159       CALL get_field("masque",zmasq,found)
    160       IF (.NOT. found) THEN
    161         PRINT*, 'phyetat0: Le champ <masque> est absent'
    162         PRINT *, 'fichier startphy non compatible avec phyetat0'
    163       ENDIF
    164 
    165        
    166 C Lecture des fractions pour chaque sous-surface
    167 C
    168 C initialisation des sous-surfaces
    169 C
    170       pctsrf = 0.
    171 C
    172 C fraction de terre
    173 C
    174 
    175       CALL get_field("FTER",pctsrf(:,is_ter),found)
    176       IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FTER> est absent'
    177 
    178 C
    179 C fraction de glace de terre
    180 C
    181       CALL get_field("FLIC",pctsrf(:,is_lic),found)
    182       IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FLIC> est absent'
    183 
    184 C
    185 C fraction d'ocean
    186 C
    187       CALL get_field("FOCE",pctsrf(:,is_oce),found)
    188       IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FOCE> est absent'
    189 
    190 C
    191 C fraction glace de mer
    192 C
    193       CALL get_field("FSIC",pctsrf(:,is_sic),found)
    194       IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FSIC> est absent'
    195 
    196 C
    197 C  Verification de l'adequation entre le masque et les sous-surfaces
    198 C
    199       fractint( 1 : klon) = pctsrf(1 : klon, is_ter)
    200      $    + pctsrf(1 : klon, is_lic)
    201       DO i = 1 , klon
    202         IF ( abs(fractint(i) - zmasq(i) ) .GT. EPSFRA ) THEN
    203             WRITE(*,*) 'phyetat0: attention fraction terre pas ',
    204      $          'coherente ', i, zmasq(i), pctsrf(i, is_ter)
    205      $          ,pctsrf(i, is_lic)
    206             WRITE(*,*) 'Je force la coherence zmasq=fractint'
    207             zmasq(i) = fractint(i)
    208         ENDIF
    209       END DO
    210       fractint (1 : klon) =  pctsrf(1 : klon, is_oce)
    211      $    + pctsrf(1 : klon, is_sic)
    212       DO i = 1 , klon
    213         IF ( abs( fractint(i) - (1. - zmasq(i))) .GT. EPSFRA ) THEN
    214             WRITE(*,*) 'phyetat0 attention fraction ocean pas ',
    215      $          'coherente ', i, zmasq(i) , pctsrf(i, is_oce)
    216      $          ,pctsrf(i, is_sic)
    217             WRITE(*,*) 'Je force la coherence zmasq=fractint'
    218             zmasq(i) = fractint(i)
    219         ENDIF
    220       END DO
    221 
    222 C
    223 c Lecture des temperatures du sol:
    224 c
    225 
    226        CALL get_field("TS",ftsol(:,1),found)
    227        IF (.NOT. found) THEN
    228          PRINT*, 'phyetat0: Le champ <TS> est absent'
    229          PRINT*, '          Mais je vais essayer de lire TS**'
    230          DO nsrf = 1, nbsrf
    231            IF (nsrf.GT.99) THEN
    232              PRINT*, "Trop de sous-mailles"
    233              CALL abort
    234            ENDIF
    235            WRITE(str2,'(i2.2)') nsrf
    236            CALL get_field("TS"//str2,ftsol(:,nsrf))
    237 
    238            xmin = 1.0E+20
    239            xmax = -1.0E+20
    240            DO i = 1, klon
    241               xmin = MIN(ftsol(i,nsrf),xmin)
    242               xmax = MAX(ftsol(i,nsrf),xmax)
    243            ENDDO
    244            PRINT*,'Temperature du sol TS**:', nsrf, xmin, xmax
    245          ENDDO
    246       ELSE
    247          PRINT*, 'phyetat0: Le champ <TS> est present'
    248          PRINT*, '          J ignore donc les autres temperatures TS**'
    249          xmin = 1.0E+20
    250          xmax = -1.0E+20
    251          DO i = 1, klon
    252             xmin = MIN(ftsol(i,1),xmin)
    253             xmax = MAX(ftsol(i,1),xmax)
    254          ENDDO
    255          PRINT*,'Temperature du sol <TS>', xmin, xmax
    256          DO nsrf = 2, nbsrf
    257          DO i = 1, klon
    258             ftsol(i,nsrf) = ftsol(i,1)
    259          ENDDO
    260          ENDDO
    261       ENDIF
    262 
    263 c
    264 c Lecture des temperatures du sol profond:
    265 c
    266       DO nsrf = 1, nbsrf
    267         DO isoil=1, nsoilmx
    268           IF (isoil.GT.99 .AND. nsrf.GT.99) THEN
    269             PRINT*, "Trop de couches ou sous-mailles"
    270             CALL abort
    271           ENDIF
    272           WRITE(str7,'(i2.2,"srf",i2.2)') isoil, nsrf
    273          
    274           CALL get_field('Tsoil'//str7,tsoil(:,isoil,nsrf),found)
    275           IF (.NOT. found) THEN
    276             PRINT*, "phyetat0: Le champ <Tsoil"//str7//"> est absent"
    277             PRINT*, "          Il prend donc la valeur de surface"
    278             DO i=1, klon
    279                tsoil(i,isoil,nsrf)=ftsol(i,nsrf)
    280             ENDDO
    281           ENDIF
    282         ENDDO
    283       ENDDO
    284 c
    285 c Lecture de l'humidite de l'air juste au dessus du sol:
    286 c
    287 
    288       CALL get_field("QS",qsurf(:,1),found)
    289       IF (.NOT. found) THEN
    290          PRINT*, 'phyetat0: Le champ <QS> est absent'
    291          PRINT*, '          Mais je vais essayer de lire QS**'
    292          DO nsrf = 1, nbsrf
    293            IF (nsrf.GT.99) THEN
    294              PRINT*, "Trop de sous-mailles"
    295              CALL abort
    296            ENDIF
    297            WRITE(str2,'(i2.2)') nsrf
    298            CALL get_field("QS"//str2,qsurf(:,nsrf))
    299            xmin = 1.0E+20
    300            xmax = -1.0E+20
    301            DO i = 1, klon
    302               xmin = MIN(qsurf(i,nsrf),xmin)
    303               xmax = MAX(qsurf(i,nsrf),xmax)
    304            ENDDO
    305            PRINT*,'Humidite pres du sol QS**:', nsrf, xmin, xmax
    306          ENDDO
    307       ELSE
    308          PRINT*, 'phyetat0: Le champ <QS> est present'
    309          PRINT*, '          J ignore donc les autres humidites QS**'
    310          xmin = 1.0E+20
    311          xmax = -1.0E+20
    312          DO i = 1, klon
    313             xmin = MIN(qsurf(i,1),xmin)
    314             xmax = MAX(qsurf(i,1),xmax)
    315          ENDDO
    316          PRINT*,'Humidite pres du sol <QS>', xmin, xmax
    317          DO nsrf = 2, nbsrf
    318            DO i = 1, klon
    319              qsurf(i,nsrf) = qsurf(i,1)
    320            ENDDO
    321          ENDDO
    322       ENDIF
    323 
    324 C
    325 C Eau dans le sol (pour le modele de sol "bucket")
    326 C
    327       CALL get_field("QSOL",qsol,found)
    328       IF (.NOT. found) THEN
    329         PRINT*, 'phyetat0: Le champ <QSOL> est absent'
    330         PRINT*, '          Valeur par defaut nulle'
    331           qsol(:)=0.
    332       ENDIF
    333 
    334       xmin = 1.0E+20
    335       xmax = -1.0E+20
    336       DO i = 1, klon
    337         xmin = MIN(qsol(i),xmin)
    338         xmax = MAX(qsol(i),xmax)
    339       ENDDO
    340       PRINT*,'Eau dans le sol (mm) <QSOL>', xmin, xmax
    341 
    342 c
    343 c Lecture de neige au sol:
    344 c
    345 
    346       CALL get_field("SNOW",snow(:,1),found)
    347       IF (.NOT. found) THEN
    348         PRINT*, 'phyetat0: Le champ <SNOW> est absent'
    349         PRINT*, '          Mais je vais essayer de lire SNOW**'
    350         DO nsrf = 1, nbsrf
    351           IF (nsrf.GT.99) THEN
    352             PRINT*, "Trop de sous-mailles"
    353             CALL abort
    354           ENDIF
    355           WRITE(str2,'(i2.2)') nsrf
    356           CALL get_field( "SNOW"//str2,snow(:,nsrf))
    357           xmin = 1.0E+20
    358           xmax = -1.0E+20
    359           DO i = 1, klon
    360             xmin = MIN(snow(i,nsrf),xmin)
    361             xmax = MAX(snow(i,nsrf),xmax)
    362           ENDDO
    363           PRINT*,'Neige du sol SNOW**:', nsrf, xmin, xmax
    364         ENDDO
    365       ELSE
    366          PRINT*, 'phyetat0: Le champ <SNOW> est present'
    367          PRINT*, '          J ignore donc les autres neiges SNOW**'
    368          xmin = 1.0E+20
    369          xmax = -1.0E+20
    370          DO i = 1, klon
    371             xmin = MIN(snow(i,1),xmin)
    372             xmax = MAX(snow(i,1),xmax)
    373          ENDDO
    374          PRINT*,'Neige du sol <SNOW>', xmin, xmax
    375          DO nsrf = 2, nbsrf
    376          DO i = 1, klon
    377             snow(i,nsrf) = snow(i,1)
    378          ENDDO
    379          ENDDO
    380       ENDIF
    381 c
    382 c Lecture de albedo de l'interval visible au sol:
    383 c
    384       CALL get_field("ALBE",falb1(:,1),found)
    385       IF (.NOT. found) THEN
    386          PRINT*, 'phyetat0: Le champ <ALBE> est absent'
    387          PRINT*, '          Mais je vais essayer de lire ALBE**'
    388          DO nsrf = 1, nbsrf
    389            IF (nsrf.GT.99) THEN
    390              PRINT*, "Trop de sous-mailles"
    391              CALL abort
    392            ENDIF
    393            WRITE(str2,'(i2.2)') nsrf
    394            CALL get_field("ALBE"//str2,falb1(:,nsrf))
    395            xmin = 1.0E+20
    396            xmax = -1.0E+20
    397            DO i = 1, klon
    398               xmin = MIN(falb1(i,nsrf),xmin)
    399               xmax = MAX(falb1(i,nsrf),xmax)
    400            ENDDO
    401            PRINT*,'Albedo du sol ALBE**:', nsrf, xmin, xmax
    402          ENDDO
    403       ELSE
    404          PRINT*, 'phyetat0: Le champ <ALBE> est present'
    405          PRINT*, '          J ignore donc les autres ALBE**'
    406          xmin = 1.0E+20
    407          xmax = -1.0E+20
    408          DO i = 1, klon
    409             xmin = MIN(falb1(i,1),xmin)
    410             xmax = MAX(falb1(i,1),xmax)
    411          ENDDO
    412          PRINT*,'Neige du sol <ALBE>', xmin, xmax
    413          DO nsrf = 2, nbsrf
    414            DO i = 1, klon
    415             falb1(i,nsrf) = falb1(i,1)
    416            ENDDO
    417          ENDDO
    418       ENDIF
    419 
    420 c
    421 c Lecture de albedo au sol dans l'interval proche infra-rouge:
    422 c
    423       CALL get_field("ALBLW",falb2(:,1),found)
    424       IF (.NOT. found) THEN
    425          PRINT*, 'phyetat0: Le champ <ALBLW> est absent'
    426          PRINT*, '          Mais je vais prendre ALBE**'
    427          DO nsrf = 1, nbsrf
    428            DO i = 1, klon
    429              falb2(i,nsrf) = falb1(i,nsrf)
    430            ENDDO
    431          ENDDO
    432       ELSE
    433          PRINT*, 'phyetat0: Le champ <ALBLW> est present'
    434          PRINT*, '          J ignore donc les autres ALBLW**'
    435          xmin = 1.0E+20
    436          xmax = -1.0E+20
    437          DO i = 1, klon
    438             xmin = MIN(falb2(i,1),xmin)
    439             xmax = MAX(falb2(i,1),xmax)
    440          ENDDO
    441          PRINT*,'Neige du sol <ALBLW>', xmin, xmax
    442          DO nsrf = 2, nbsrf
    443            DO i = 1, klon
    444              falb2(i,nsrf) = falb2(i,1)
    445            ENDDO
    446          ENDDO
    447       ENDIF
    448 c
    449 c Lecture de evaporation: 
    450 c
    451       CALL get_field("EVAP",evap(:,1),found)
    452       IF (.NOT. found) THEN
    453          PRINT*, 'phyetat0: Le champ <EVAP> est absent'
    454          PRINT*, '          Mais je vais essayer de lire EVAP**'
    455          DO nsrf = 1, nbsrf
    456            IF (nsrf.GT.99) THEN
    457              PRINT*, "Trop de sous-mailles"
    458              CALL abort
    459            ENDIF
    460            WRITE(str2,'(i2.2)') nsrf
    461            CALL get_field("EVAP"//str2, evap(:,nsrf))
    462            xmin = 1.0E+20
    463            xmax = -1.0E+20
    464            DO i = 1, klon
    465               xmin = MIN(evap(i,nsrf),xmin)
    466               xmax = MAX(evap(i,nsrf),xmax)
    467            ENDDO
    468            PRINT*,'evap du sol EVAP**:', nsrf, xmin, xmax
    469          ENDDO
    470       ELSE
    471          PRINT*, 'phyetat0: Le champ <EVAP> est present'
    472          PRINT*, '          J ignore donc les autres EVAP**'
    473          xmin = 1.0E+20
    474          xmax = -1.0E+20
    475          DO i = 1, klon
    476             xmin = MIN(evap(i,1),xmin)
    477             xmax = MAX(evap(i,1),xmax)
    478          ENDDO
    479          PRINT*,'Evap du sol <EVAP>', xmin, xmax
    480          DO nsrf = 2, nbsrf
    481          DO i = 1, klon
    482             evap(i,nsrf) = evap(i,1)
    483          ENDDO
    484          ENDDO
    485       ENDIF
    486 c
    487 c Lecture precipitation liquide:
    488 c
    489       CALL get_field("rain_f",rain_fall)
    490       xmin = 1.0E+20
    491       xmax = -1.0E+20
    492       DO i = 1, klon
    493          xmin = MIN(rain_fall(i),xmin)
    494          xmax = MAX(rain_fall(i),xmax)
    495       ENDDO
    496       PRINT*,'Precipitation liquide rain_f:', xmin, xmax
    497 c
    498 c Lecture precipitation solide:
    499 c
    500       CALL get_field("snow_f",snow_fall)
    501       xmin = 1.0E+20
    502       xmax = -1.0E+20
    503       DO i = 1, klon
    504          xmin = MIN(snow_fall(i),xmin)
    505          xmax = MAX(snow_fall(i),xmax)
    506       ENDDO
    507       PRINT*,'Precipitation solide snow_f:', xmin, xmax
    508 c
    509 c Lecture rayonnement solaire au sol:
    510 c
    511       CALL get_field("solsw",solsw,found)
    512       IF (.NOT. found) THEN
    513          PRINT*, 'phyetat0: Le champ <solsw> est absent'
    514          PRINT*, 'mis a zero'
    515          solsw(:) = 0.
    516       ENDIF
    517       xmin = 1.0E+20
    518       xmax = -1.0E+20
    519       DO i = 1, klon
    520          xmin = MIN(solsw(i),xmin)
    521          xmax = MAX(solsw(i),xmax)
    522       ENDDO
    523       PRINT*,'Rayonnement solaire au sol solsw:', xmin, xmax
    524 c
    525 c Lecture rayonnement IF au sol:
    526 c
    527       CALL get_field("sollw",sollw,found)
    528       IF (.NOT. found) THEN
    529          PRINT*, 'phyetat0: Le champ <sollw> est absent'
    530          PRINT*, 'mis a zero'
    531          sollw = 0.
    532       ENDIF
    533       xmin = 1.0E+20
    534       xmax = -1.0E+20
    535       DO i = 1, klon
    536          xmin = MIN(sollw(i),xmin)
    537          xmax = MAX(sollw(i),xmax)
    538       ENDDO
    539       PRINT*,'Rayonnement IF au sol sollw:', xmin, xmax
    540      
    541 c
    542 c Lecture derive des flux:
    543 c
    544       CALL get_field("fder",fder,found)
    545       IF (.NOT. found) THEN
    546          PRINT*, 'phyetat0: Le champ <fder> est absent'
    547          PRINT*, 'mis a zero'
    548          fder = 0.
    549       ENDIF
    550       xmin = 1.0E+20
    551       xmax = -1.0E+20
    552       DO i = 1, klon
    553          xmin = MIN(fder(i),xmin)
    554          xmax = MAX(fder(i),xmax)
    555       ENDDO
    556       PRINT*,'Derive des flux fder:', xmin, xmax
    557 
    558 c
    559 c Lecture du rayonnement net au sol:
    560 c
    561       CALL get_field("RADS",radsol)
    562       xmin = 1.0E+20
    563       xmax = -1.0E+20
    564       DO i = 1, klon
    565          xmin = MIN(radsol(i),xmin)
    566          xmax = MAX(radsol(i),xmax)
    567       ENDDO
    568       PRINT*,'Rayonnement net au sol radsol:', xmin, xmax
    569 c
    570 c Lecture de la longueur de rugosite
    571 c
    572 c
    573       CALL get_field("RUG",frugs(:,1),found)
    574       IF (.NOT. found) THEN
    575          PRINT*, 'phyetat0: Le champ <RUG> est absent'
    576          PRINT*, '          Mais je vais essayer de lire RUG**'
    577          DO nsrf = 1, nbsrf
    578            IF (nsrf.GT.99) THEN
    579              PRINT*, "Trop de sous-mailles"
    580              CALL abort
    581            ENDIF
    582            WRITE(str2,'(i2.2)') nsrf
    583            CALL get_field("RUG"//str2,frugs(:,nsrf))
    584            xmin = 1.0E+20
    585            xmax = -1.0E+20
    586            DO i = 1, klon
    587               xmin = MIN(frugs(i,nsrf),xmin)
    588               xmax = MAX(frugs(i,nsrf),xmax)
    589            ENDDO
    590            PRINT*,'rugosite du sol RUG**:', nsrf, xmin, xmax
    591          ENDDO
    592       ELSE
    593          PRINT*, 'phyetat0: Le champ <RUG> est present'
    594          PRINT*, '          J ignore donc les autres RUG**'
    595          xmin = 1.0E+20
    596          xmax = -1.0E+20
    597          DO i = 1, klon
    598             xmin = MIN(frugs(i,1),xmin)
    599             xmax = MAX(frugs(i,1),xmax)
    600          ENDDO
    601          PRINT*,'rugosite <RUG>', xmin, xmax
    602          DO nsrf = 2, nbsrf
    603          DO i = 1, klon
    604             frugs(i,nsrf) = frugs(i,1)
    605          ENDDO
    606          ENDDO
    607       ENDIF
    608 
    609 c
    610 c Lecture de l'age de la neige:
    611 c
    612       CALL get_field("AGESNO",agesno(:,1),found)
    613       IF (.NOT. found) THEN
    614          PRINT*, 'phyetat0: Le champ <AGESNO> est absent'
    615          PRINT*, '          Mais je vais essayer de lire AGESNO**'
    616          DO nsrf = 1, nbsrf
    617            IF (nsrf.GT.99) THEN
    618              PRINT*, "Trop de sous-mailles"
    619              CALL abort
    620            ENDIF
    621            WRITE(str2,'(i2.2)') nsrf
    622            CALL get_field("AGESNO"//str2,agesno(:,nsrf),found)
    623            IF (.NOT. found) THEN
    624               PRINT*, "phyetat0: Le champ <AGESNO"//str2//"> est absent"
    625               agesno = 50.0
    626            ENDIF
    627            xmin = 1.0E+20
    628            xmax = -1.0E+20
    629            DO i = 1, klon
    630               xmin = MIN(agesno(i,nsrf),xmin)
    631               xmax = MAX(agesno(i,nsrf),xmax)
    632            ENDDO
    633            PRINT*,'Age de la neige AGESNO**:', nsrf, xmin, xmax
    634          ENDDO
    635       ELSE
    636          PRINT*, 'phyetat0: Le champ <AGESNO> est present'
    637          PRINT*, '          J ignore donc les autres AGESNO**'
    638          xmin = 1.0E+20
    639          xmax = -1.0E+20
    640          DO i = 1, klon
    641             xmin = MIN(agesno(i,1),xmin)
    642             xmax = MAX(agesno(i,1),xmax)
    643          ENDDO
    644          PRINT*,'Age de la neige <AGESNO>', xmin, xmax
    645          DO nsrf = 2, nbsrf
    646          DO i = 1, klon
    647             agesno(i,nsrf) = agesno(i,1)
    648          ENDDO
    649          ENDDO
    650       ENDIF
    651 
    652 c
    653       CALL get_field("ZMEA", zmea)
    654       xmin = 1.0E+20
    655       xmax = -1.0E+20
    656       DO i = 1, klon
    657          xmin = MIN(zmea(i),xmin)
    658          xmax = MAX(zmea(i),xmax)
    659       ENDDO
    660       PRINT*,'OROGRAPHIE SOUS-MAILLE zmea:', xmin, xmax
    661 c
    662 c
    663       CALL get_field("ZSTD",zstd)
    664       xmin = 1.0E+20
    665       xmax = -1.0E+20
    666       DO i = 1, klon
    667          xmin = MIN(zstd(i),xmin)
    668          xmax = MAX(zstd(i),xmax)
    669       ENDDO
    670       PRINT*,'OROGRAPHIE SOUS-MAILLE zstd:', xmin, xmax
    671 c
    672 c
    673       CALL get_field("ZSIG",zsig)
    674       xmin = 1.0E+20
    675       xmax = -1.0E+20
    676       DO i = 1, klon
    677          xmin = MIN(zsig(i),xmin)
    678          xmax = MAX(zsig(i),xmax)
    679       ENDDO
    680       PRINT*,'OROGRAPHIE SOUS-MAILLE zsig:', xmin, xmax
    681 c
    682 c
    683       CALL get_field("ZGAM",zgam)
    684       xmin = 1.0E+20
    685       xmax = -1.0E+20
    686       DO i = 1, klon
    687          xmin = MIN(zgam(i),xmin)
    688          xmax = MAX(zgam(i),xmax)
    689       ENDDO
    690       PRINT*,'OROGRAPHIE SOUS-MAILLE zgam:', xmin, xmax
    691 c
    692 c
    693       CALL get_field("ZTHE",zthe)
    694       xmin = 1.0E+20
    695       xmax = -1.0E+20
    696       DO i = 1, klon
    697          xmin = MIN(zthe(i),xmin)
    698          xmax = MAX(zthe(i),xmax)
    699       ENDDO
    700       PRINT*,'OROGRAPHIE SOUS-MAILLE zthe:', xmin, xmax
    701 c
    702 c
    703       CALL get_field("ZPIC",zpic)
    704       xmin = 1.0E+20
    705       xmax = -1.0E+20
    706       DO i = 1, klon
    707          xmin = MIN(zpic(i),xmin)
    708          xmax = MAX(zpic(i),xmax)
    709       ENDDO
    710       PRINT*,'OROGRAPHIE SOUS-MAILLE zpic:', xmin, xmax
    711 c
    712       CALL get_field("ZVAL",zval)
    713       xmin = 1.0E+20
    714       xmax = -1.0E+20
    715       DO i = 1, klon
    716          xmin = MIN(zval(i),xmin)
    717          xmax = MAX(zval(i),xmax)
    718       ENDDO
    719       PRINT*,'OROGRAPHIE SOUS-MAILLE zval:', xmin, xmax
    720 c
    721 c
    722       CALL get_field("RUGSREL",rugoro)
    723       xmin = 1.0E+20
    724       xmax = -1.0E+20
    725       DO i = 1, klon
    726          xmin = MIN(rugoro(i),xmin)
    727          xmax = MAX(rugoro(i),xmax)
    728       ENDDO
    729       PRINT*,'Rugosite relief (ecart-type) rugsrel:', xmin, xmax
    730 c
    731 c
    732      
    733 c
    734       ancien_ok = .TRUE.
    735 
    736       CALL get_field("TANCIEN",t_ancien,found)
    737       IF (.NOT. found) THEN
    738          PRINT*, "phyetat0: Le champ <TANCIEN> est absent"
    739          PRINT*, "Depart legerement fausse. Mais je continue"
    740          ancien_ok = .FALSE.
    741       ENDIF
    742 
    743 
    744       CALL get_field("QANCIEN",q_ancien,found)
    745       IF (.NOT. found) THEN
    746          PRINT*, "phyetat0: Le champ <QANCIEN> est absent"
    747          PRINT*, "Depart legerement fausse. Mais je continue"
    748          ancien_ok = .FALSE.
    749       ENDIF
    750 
    751       CALL get_field("UANCIEN",u_ancien,found)
    752       IF (.NOT. found) THEN
    753          PRINT*, "phyetat0: Le champ <UANCIEN> est absent"
    754          PRINT*, "Depart legerement fausse. Mais je continue"
    755          ancien_ok = .FALSE.
    756       ENDIF
    757 
    758       CALL get_field("VANCIEN",v_ancien,found)
    759       IF (.NOT. found) THEN
    760          PRINT*, "phyetat0: Le champ <VANCIEN> est absent"
    761          PRINT*, "Depart legerement fausse. Mais je continue"
    762          ancien_ok = .FALSE.
    763       ENDIF
    764 
    765       clwcon=0.
    766       CALL get_field("CLWCON",clwcon,found)
    767       IF (.NOT. found) THEN
    768          PRINT*, "phyetat0: Le champ CLWCON est absent"
    769          PRINT*, "Depart legerement fausse. Mais je continue"
    770       ENDIF
    771       xmin = 1.0E+20
    772       xmax = -1.0E+20
    773       xmin = MINval(clwcon)
    774       xmax = MAXval(clwcon)
    775       PRINT*,'Eau liquide convective (ecart-type) clwcon:', xmin, xmax
    776 c
    777       rnebcon = 0.
    778       CALL get_field("RNEBCON",rnebcon,found)
    779       IF (.NOT. found) THEN
    780          PRINT*, "phyetat0: Le champ RNEBCON est absent"
    781          PRINT*, "Depart legerement fausse. Mais je continue"
    782       ENDIF
    783       xmin = 1.0E+20
    784       xmax = -1.0E+20
    785       xmin = MINval(rnebcon)
    786       xmax = MAXval(rnebcon)
    787       PRINT*,'Nebulosite convective (ecart-type) rnebcon:', xmin, xmax
    788 
    789 c
    790 c Lecture ratqs
    791 c
    792       ratqs=0.
    793       CALL get_field("RATQS",ratqs,found)
    794       IF (.NOT. found) THEN
    795          PRINT*, "phyetat0: Le champ <RATQS> est absent"
    796          PRINT*, "Depart legerement fausse. Mais je continue"
    797       ENDIF
    798       xmin = 1.0E+20
    799       xmax = -1.0E+20
    800       xmin = MINval(ratqs)
    801       xmax = MAXval(ratqs)
    802       PRINT*,'(ecart-type) ratqs:', xmin, xmax
    803 c
    804 c Lecture run_off_lic_0
    805 c
    806       CALL get_field("RUNOFFLIC0",run_off_lic_0,found)
    807       IF (.NOT. found) THEN
    808          PRINT*, "phyetat0: Le champ <RUNOFFLIC0> est absent"
    809          PRINT*, "Depart legerement fausse. Mais je continue"
    810          run_off_lic_0 = 0.
    811       ENDIF
    812       xmin = 1.0E+20
    813       xmax = -1.0E+20
    814       xmin = MINval(run_off_lic_0)
    815       xmax = MAXval(run_off_lic_0)
    816       PRINT*,'(ecart-type) run_off_lic_0:', xmin, xmax
    817 
    818 
    819 c Lecture de l'energie cinetique turbulente
    820 c
    821 
    822       IF (iflag_pbl>1) then
    823         DO nsrf = 1, nbsrf
    824           IF (nsrf.GT.99) THEN
    825             PRINT*, "Trop de sous-mailles"
    826             CALL abort
    827           ENDIF
    828           WRITE(str2,'(i2.2)') nsrf
    829           CALL get_field("TKE"//str2,pbl_tke(:,1:klev+1,nsrf),found)
    830           IF (.NOT. found) THEN
    831             PRINT*, "phyetat0: <TKE"//str2//"> est absent"
    832             pbl_tke(:,:,nsrf)=1.e-8
    833           ENDIF
    834           xmin = 1.0E+20
    835           xmax = -1.0E+20
    836           DO k = 1, klev+1
    837             DO i = 1, klon
    838               xmin = MIN(pbl_tke(i,k,nsrf),xmin)
    839               xmax = MAX(pbl_tke(i,k,nsrf),xmax)
    840             ENDDO
    841           ENDDO
    842           PRINT*,'Temperature du sol TKE**:', nsrf, xmin, xmax
    843         ENDDO
    844       ENDIF
    845 c
    846 c zmax0
    847       CALL get_field("ZMAX0",zmax0,found)
    848       IF (.NOT. found) THEN
    849         PRINT*, "phyetat0: Le champ <ZMAX0> est absent"
    850         PRINT*, "Depart legerement fausse. Mais je continue"
    851         zmax0=40.
    852       ENDIF
    853       xmin = 1.0E+20
    854       xmax = -1.0E+20
    855       xmin = MINval(zmax0)
    856       xmax = MAXval(zmax0)
    857       PRINT*,'(ecart-type) zmax0:', xmin, xmax
    858 c
    859 c           f0(ig)=1.e-5
    860 c f0
    861       CALL get_field("F0",f0,found)
    862       IF (.NOT. found) THEN
    863          PRINT*, "phyetat0: Le champ <f0> est absent"
    864          PRINT*, "Depart legerement fausse. Mais je continue"
    865          f0=1.e-5
    866       ENDIF
    867       xmin = 1.0E+20
    868       xmax = -1.0E+20
    869       xmin = MINval(f0)
    870       xmax = MAXval(f0)
    871       PRINT*,'(ecart-type) f0:', xmin, xmax
    872 c
    873 c ema_work1
    874 c
    875       CALL get_field("EMA_WORK1",ema_work1,found)
    876       IF (.NOT. found) THEN
    877         PRINT*, "phyetat0: Le champ <EMA_WORK1> est absent"
    878         PRINT*, "Depart legerement fausse. Mais je continue"
    879         ema_work1=0.
    880       ELSE
     93  DO i = 1, length
     94     tabcntr0( i ) = tab_cntrl( i )
     95  ENDDO
     96
     97  tab_cntrl(1)=dtime
     98  tab_cntrl(2)=radpas
     99
     100  ! co2_ppm : value from the previous time step
     101  IF (carbon_cycle_tr .OR. carbon_cycle_cpl) THEN
     102     co2_ppm = tab_cntrl(3)
     103     RCO2    = co2_ppm * 1.0e-06  * 44.011/28.97
     104     ! ELSE : keep value from .def
     105  END IF
     106
     107  ! co2_ppm0 : initial value of atmospheric CO2 (from create_etat0_limit.e .def)
     108  co2_ppm0   = tab_cntrl(16)
     109
     110  solaire_etat0      = tab_cntrl(4)
     111  tab_cntrl(5)=iflag_con
     112  tab_cntrl(6)=nbapp_rad
     113
     114  if (cycle_diurne) tab_cntrl( 7) =1.
     115  if (soil_model) tab_cntrl( 8) =1.
     116  if (new_oliq) tab_cntrl( 9) =1.
     117  if (ok_orodr) tab_cntrl(10) =1.
     118  if (ok_orolf) tab_cntrl(11) =1.
     119  if (ok_limitvrai) tab_cntrl(12) =1.
     120
     121  itau_phy = tab_cntrl(15)
     122
     123  clesphy0(1)=tab_cntrl( 5 )
     124  clesphy0(2)=tab_cntrl( 6 )
     125  clesphy0(3)=tab_cntrl( 7 )
     126  clesphy0(4)=tab_cntrl( 8 )
     127  clesphy0(5)=tab_cntrl( 9 )
     128  clesphy0(6)=tab_cntrl( 10 )
     129  clesphy0(7)=tab_cntrl( 11 )
     130  clesphy0(8)=tab_cntrl( 12 )
     131
     132  ! Lecture des latitudes (coordonnees):
     133
     134  CALL get_field("latitude", rlat)
     135
     136  ! Lecture des longitudes (coordonnees):
     137
     138  CALL get_field("longitude", rlon)
     139
     140  ! Lecture du masque terre mer
     141
     142  CALL get_field("masque", zmasq, found)
     143  IF (.NOT. found) THEN
     144     PRINT*, 'phyetat0: Le champ <masque> est absent'
     145     PRINT *, 'fichier startphy non compatible avec phyetat0'
     146  ENDIF
     147
     148  ! Lecture des fractions pour chaque sous-surface
     149
     150  ! initialisation des sous-surfaces
     151
     152  pctsrf = 0.
     153
     154  ! fraction de terre
     155
     156  CALL get_field("FTER", pctsrf(:, is_ter), found)
     157  IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FTER> est absent'
     158
     159  ! fraction de glace de terre
     160
     161  CALL get_field("FLIC", pctsrf(:, is_lic), found)
     162  IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FLIC> est absent'
     163
     164  ! fraction d'ocean
     165
     166  CALL get_field("FOCE", pctsrf(:, is_oce), found)
     167  IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FOCE> est absent'
     168
     169  ! fraction glace de mer
     170
     171  CALL get_field("FSIC", pctsrf(:, is_sic), found)
     172  IF (.NOT. found) PRINT*, 'phyetat0: Le champ <FSIC> est absent'
     173
     174  !  Verification de l'adequation entre le masque et les sous-surfaces
     175
     176  fractint( 1 : klon) = pctsrf(1 : klon, is_ter)  &
     177       + pctsrf(1 : klon, is_lic)
     178  DO i = 1 , klon
     179     IF ( abs(fractint(i) - zmasq(i) ) .GT. EPSFRA ) THEN
     180        WRITE(*, *) 'phyetat0: attention fraction terre pas ',  &
     181             'coherente ', i, zmasq(i), pctsrf(i, is_ter) &
     182             , pctsrf(i, is_lic)
     183        WRITE(*, *) 'Je force la coherence zmasq=fractint'
     184        zmasq(i) = fractint(i)
     185     ENDIF
     186  END DO
     187  fractint (1 : klon) =  pctsrf(1 : klon, is_oce)  &
     188       + pctsrf(1 : klon, is_sic)
     189  DO i = 1 , klon
     190     IF ( abs( fractint(i) - (1. - zmasq(i))) .GT. EPSFRA ) THEN
     191        WRITE(*, *) 'phyetat0 attention fraction ocean pas ',  &
     192             'coherente ', i, zmasq(i) , pctsrf(i, is_oce) &
     193             , pctsrf(i, is_sic)
     194        WRITE(*, *) 'Je force la coherence zmasq=fractint'
     195        zmasq(i) = fractint(i)
     196     ENDIF
     197  END DO
     198
     199  ! Lecture des temperatures du sol:
     200
     201  CALL get_field("TS", ftsol(:, 1), found)
     202  IF (.NOT. found) THEN
     203     PRINT*, 'phyetat0: Le champ <TS> est absent'
     204     PRINT*, '          Mais je vais essayer de lire TS**'
     205     DO nsrf = 1, nbsrf
     206        IF (nsrf.GT.99) THEN
     207           PRINT*, "Trop de sous-mailles"
     208           CALL abort
     209        ENDIF
     210        WRITE(str2, '(i2.2)') nsrf
     211        CALL get_field("TS"//str2, ftsol(:, nsrf))
     212
    881213        xmin = 1.0E+20
    882214        xmax = -1.0E+20
    883         DO k = 1, klev
    884           DO i = 1, klon
    885             xmin = MIN(ema_work1(i,k),xmin)
    886             xmax = MAX(ema_work1(i,k),xmax)
    887           ENDDO
    888         ENDDO
    889         PRINT*,'ema_work1:', xmin, xmax
    890       ENDIF
    891 c
    892 c ema_work2
    893 c
    894       CALL get_field("EMA_WORK2",ema_work2,found)
    895       IF (.NOT. found) THEN
    896         PRINT*, "phyetat0: Le champ <EMA_WORK2> est absent"
    897         PRINT*, "Depart legerement fausse. Mais je continue"
    898         ema_work2=0.
    899       ELSE
     215        DO i = 1, klon
     216           xmin = MIN(ftsol(i, nsrf), xmin)
     217           xmax = MAX(ftsol(i, nsrf), xmax)
     218        ENDDO
     219        PRINT*, 'Temperature du sol TS**:', nsrf, xmin, xmax
     220     ENDDO
     221  ELSE
     222     PRINT*, 'phyetat0: Le champ <TS> est present'
     223     PRINT*, '          J ignore donc les autres temperatures TS**'
     224     xmin = 1.0E+20
     225     xmax = -1.0E+20
     226     DO i = 1, klon
     227        xmin = MIN(ftsol(i, 1), xmin)
     228        xmax = MAX(ftsol(i, 1), xmax)
     229     ENDDO
     230     PRINT*, 'Temperature du sol <TS>', xmin, xmax
     231     DO nsrf = 2, nbsrf
     232        DO i = 1, klon
     233           ftsol(i, nsrf) = ftsol(i, 1)
     234        ENDDO
     235     ENDDO
     236  ENDIF
     237
     238  ! Lecture des temperatures du sol profond:
     239
     240  DO nsrf = 1, nbsrf
     241     DO isoil=1, nsoilmx
     242        IF (isoil.GT.99 .AND. nsrf.GT.99) THEN
     243           PRINT*, "Trop de couches ou sous-mailles"
     244           CALL abort
     245        ENDIF
     246        WRITE(str7, '(i2.2, "srf", i2.2)') isoil, nsrf
     247
     248        CALL get_field('Tsoil'//str7, tsoil(:, isoil, nsrf), found)
     249        IF (.NOT. found) THEN
     250           PRINT*, "phyetat0: Le champ <Tsoil"//str7//"> est absent"
     251           PRINT*, "          Il prend donc la valeur de surface"
     252           DO i=1, klon
     253              tsoil(i, isoil, nsrf)=ftsol(i, nsrf)
     254           ENDDO
     255        ENDIF
     256     ENDDO
     257  ENDDO
     258
     259  ! Lecture de l'humidite de l'air juste au dessus du sol:
     260
     261  CALL get_field("QS", qsurf(:, 1), found)
     262  IF (.NOT. found) THEN
     263     PRINT*, 'phyetat0: Le champ <QS> est absent'
     264     PRINT*, '          Mais je vais essayer de lire QS**'
     265     DO nsrf = 1, nbsrf
     266        IF (nsrf.GT.99) THEN
     267           PRINT*, "Trop de sous-mailles"
     268           CALL abort
     269        ENDIF
     270        WRITE(str2, '(i2.2)') nsrf
     271        CALL get_field("QS"//str2, qsurf(:, nsrf))
    900272        xmin = 1.0E+20
    901273        xmax = -1.0E+20
    902         DO k = 1, klev
    903           DO i = 1, klon
    904             xmin = MIN(ema_work2(i,k),xmin)
    905             xmax = MAX(ema_work2(i,k),xmax)
    906           ENDDO
    907         ENDDO
    908         PRINT*,'ema_work2:', xmin, xmax
    909       ENDIF
    910 c
    911 c wake_deltat
    912 c
    913       CALL get_field("WAKE_DELTAT",wake_deltat,found)
    914       IF (.NOT. found) THEN
    915         PRINT*, "phyetat0: Le champ <WAKE_DELTAT> est absent"
    916         PRINT*, "Depart legerement fausse. Mais je continue"
    917         wake_deltat=0.
    918       ELSE
     274        DO i = 1, klon
     275           xmin = MIN(qsurf(i, nsrf), xmin)
     276           xmax = MAX(qsurf(i, nsrf), xmax)
     277        ENDDO
     278        PRINT*, 'Humidite pres du sol QS**:', nsrf, xmin, xmax
     279     ENDDO
     280  ELSE
     281     PRINT*, 'phyetat0: Le champ <QS> est present'
     282     PRINT*, '          J ignore donc les autres humidites QS**'
     283     xmin = 1.0E+20
     284     xmax = -1.0E+20
     285     DO i = 1, klon
     286        xmin = MIN(qsurf(i, 1), xmin)
     287        xmax = MAX(qsurf(i, 1), xmax)
     288     ENDDO
     289     PRINT*, 'Humidite pres du sol <QS>', xmin, xmax
     290     DO nsrf = 2, nbsrf
     291        DO i = 1, klon
     292           qsurf(i, nsrf) = qsurf(i, 1)
     293        ENDDO
     294     ENDDO
     295  ENDIF
     296
     297  ! Eau dans le sol (pour le modele de sol "bucket")
     298
     299  CALL get_field("QSOL", qsol, found)
     300  IF (.NOT. found) THEN
     301     PRINT*, 'phyetat0: Le champ <QSOL> est absent'
     302     PRINT*, '          Valeur par defaut nulle'
     303     qsol(:)=0.
     304  ENDIF
     305
     306  xmin = 1.0E+20
     307  xmax = -1.0E+20
     308  DO i = 1, klon
     309     xmin = MIN(qsol(i), xmin)
     310     xmax = MAX(qsol(i), xmax)
     311  ENDDO
     312  PRINT*, 'Eau dans le sol (mm) <QSOL>', xmin, xmax
     313
     314  ! Lecture de neige au sol:
     315
     316  CALL get_field("SNOW", snow(:, 1), found)
     317  IF (.NOT. found) THEN
     318     PRINT*, 'phyetat0: Le champ <SNOW> est absent'
     319     PRINT*, '          Mais je vais essayer de lire SNOW**'
     320     DO nsrf = 1, nbsrf
     321        IF (nsrf.GT.99) THEN
     322           PRINT*, "Trop de sous-mailles"
     323           CALL abort
     324        ENDIF
     325        WRITE(str2, '(i2.2)') nsrf
     326        CALL get_field( "SNOW"//str2, snow(:, nsrf))
    919327        xmin = 1.0E+20
    920328        xmax = -1.0E+20
    921         DO k = 1, klev
    922           DO i = 1, klon
    923             xmin = MIN(wake_deltat(i,k),xmin)
    924             xmax = MAX(wake_deltat(i,k),xmax)
    925           ENDDO
    926         ENDDO
    927         PRINT*,'wake_deltat:', xmin, xmax
    928       ENDIF
    929 c
    930 c wake_deltaq
    931 c   
    932       CALL get_field("WAKE_DELTAQ",wake_deltaq,found)
    933       IF (.NOT. found) THEN
    934         PRINT*, "phyetat0: Le champ <WAKE_DELTAQ> est absent"
    935         PRINT*, "Depart legerement fausse. Mais je continue"
    936         wake_deltaq=0.
    937       ELSE
     329        DO i = 1, klon
     330           xmin = MIN(snow(i, nsrf), xmin)
     331           xmax = MAX(snow(i, nsrf), xmax)
     332        ENDDO
     333        PRINT*, 'Neige du sol SNOW**:', nsrf, xmin, xmax
     334     ENDDO
     335  ELSE
     336     PRINT*, 'phyetat0: Le champ <SNOW> est present'
     337     PRINT*, '          J ignore donc les autres neiges SNOW**'
     338     xmin = 1.0E+20
     339     xmax = -1.0E+20
     340     DO i = 1, klon
     341        xmin = MIN(snow(i, 1), xmin)
     342        xmax = MAX(snow(i, 1), xmax)
     343     ENDDO
     344     PRINT*, 'Neige du sol <SNOW>', xmin, xmax
     345     DO nsrf = 2, nbsrf
     346        DO i = 1, klon
     347           snow(i, nsrf) = snow(i, 1)
     348        ENDDO
     349     ENDDO
     350  ENDIF
     351
     352  ! Lecture de albedo de l'interval visible au sol:
     353
     354  CALL get_field("ALBE", falb1(:, 1), found)
     355  IF (.NOT. found) THEN
     356     PRINT*, 'phyetat0: Le champ <ALBE> est absent'
     357     PRINT*, '          Mais je vais essayer de lire ALBE**'
     358     DO nsrf = 1, nbsrf
     359        IF (nsrf.GT.99) THEN
     360           PRINT*, "Trop de sous-mailles"
     361           CALL abort
     362        ENDIF
     363        WRITE(str2, '(i2.2)') nsrf
     364        CALL get_field("ALBE"//str2, falb1(:, nsrf))
    938365        xmin = 1.0E+20
    939366        xmax = -1.0E+20
    940         DO k = 1, klev
    941           DO i = 1, klon
    942             xmin = MIN(wake_deltaq(i,k),xmin)
    943             xmax = MAX(wake_deltaq(i,k),xmax)
    944           ENDDO
    945         ENDDO
    946         PRINT*,'wake_deltaq:', xmin, xmax
    947       ENDIF
    948 c
    949 c wake_s
    950 c
    951       CALL get_field("WAKE_S",wake_s,found)
    952       IF (.NOT. found) THEN
    953         PRINT*, "phyetat0: Le champ <WAKE_S> est absent"
    954         PRINT*, "Depart legerement fausse. Mais je continue"
    955         wake_s=0.
    956       ENDIF
    957       xmin = 1.0E+20
    958       xmax = -1.0E+20
    959       xmin = MINval(wake_s)
    960       xmax = MAXval(wake_s)
    961       PRINT*,'(ecart-type) wake_s:', xmin, xmax
    962 c
    963 c wake_cstar
    964 c
    965       CALL get_field("WAKE_CSTAR",wake_cstar,found)
    966       IF (.NOT. found) THEN
    967          PRINT*, "phyetat0: Le champ <WAKE_CSTAR> est absent"
    968          PRINT*, "Depart legerement fausse. Mais je continue"
    969          wake_cstar=0.
    970       ENDIF
    971       xmin = 1.0E+20
    972       xmax = -1.0E+20
    973       xmin = MINval(wake_cstar)
    974       xmax = MAXval(wake_cstar)
    975       PRINT*,'(ecart-type) wake_cstar:', xmin, xmax
    976 c
    977 c wake_pe
    978 c
    979       CALL get_field("WAKE_PE",wake_pe,found)
    980       IF (.NOT. found) THEN
    981          PRINT*, "phyetat0: Le champ <WAKE_PE> est absent"
    982          PRINT*, "Depart legerement fausse. Mais je continue"
    983          wake_pe=0.
    984       ENDIF
    985       xmin = 1.0E+20
    986       xmax = -1.0E+20
    987       xmin = MINval(wake_pe)
    988       xmax = MAXval(wake_pe)
    989       PRINT*,'(ecart-type) wake_pe:', xmin, xmax
    990 c
    991 c wake_fip
    992 c
    993       CALL get_field("WAKE_FIP",wake_fip,found)
    994       IF (.NOT. found) THEN
    995          PRINT*, "phyetat0: Le champ <WAKE_FIP> est absent"
    996          PRINT*, "Depart legerement fausse. Mais je continue"
    997          wake_fip=0.
    998       ENDIF
    999       xmin = 1.0E+20
    1000       xmax = -1.0E+20
    1001       xmin = MINval(wake_fip)
    1002       xmax = MAXval(wake_fip)
    1003       PRINT*,'(ecart-type) wake_fip:', xmin, xmax
    1004 c
    1005 c  thermiques
    1006 c
    1007 
    1008       CALL get_field("FM_THERM",fm_therm,found)
    1009       IF (.NOT. found) THEN
    1010          PRINT*, "phyetat0: Le champ <fm_therm> est absent"
    1011          PRINT*, "Depart legerement fausse. Mais je continue"
    1012          fm_therm=0.
    1013       ENDIF
    1014       xmin = 1.0E+20
    1015       xmax = -1.0E+20
    1016       xmin = MINval(fm_therm)
    1017       xmax = MAXval(fm_therm)
    1018       PRINT*,'(ecart-type) fm_therm:', xmin, xmax
    1019 
    1020       CALL get_field("ENTR_THERM",entr_therm,found)
    1021       IF (.NOT. found) THEN
    1022          PRINT*, "phyetat0: Le champ <entr_therm> est absent"
    1023          PRINT*, "Depart legerement fausse. Mais je continue"
    1024          entr_therm=0.
    1025       ENDIF
    1026       xmin = 1.0E+20
    1027       xmax = -1.0E+20
    1028       xmin = MINval(entr_therm)
    1029       xmax = MAXval(entr_therm)
    1030       PRINT*,'(ecart-type) entr_therm:', xmin, xmax
    1031 
    1032       CALL get_field("DETR_THERM",detr_therm,found)
    1033       IF (.NOT. found) THEN
    1034          PRINT*, "phyetat0: Le champ <detr_therm> est absent"
    1035          PRINT*, "Depart legerement fausse. Mais je continue"
    1036          detr_therm=0.
    1037       ENDIF
    1038       xmin = 1.0E+20
    1039       xmax = -1.0E+20
    1040       xmin = MINval(detr_therm)
    1041       xmax = MAXval(detr_therm)
    1042       PRINT*,'(ecart-type) detr_therm:', xmin, xmax
    1043 
    1044 
    1045 
    1046 c
    1047 c Read and send field trs to traclmdz
    1048 c
    1049       IF (type_trac == 'lmdz') THEN
    1050          DO it=1,nbtr
    1051             iiq=niadv(it+2)
    1052             CALL get_field("trs_"//tname(iiq),trs(:,it),found)
    1053             IF (.NOT. found) THEN
    1054                PRINT*,
    1055      $           "phyetat0: Le champ <trs_"//tname(iiq)//"> est absent"
    1056                PRINT*, "Depart legerement fausse. Mais je continue"
    1057                trs(:,it) = 0.
    1058             ENDIF
    1059             xmin = 1.0E+20
    1060             xmax = -1.0E+20
    1061             xmin = MINval(trs(:,it))
    1062             xmax = MAXval(trs(:,it))
    1063             PRINT*,"(ecart-type) trs_"//tname(iiq)//" :", xmin, xmax
    1064 
    1065          END DO
    1066          CALL traclmdz_from_restart(trs)
    1067 
    1068          IF (carbon_cycle_cpl) THEN
    1069             ALLOCATE(co2_send(klon), stat=ierr)
    1070             IF (ierr /= 0) CALL abort_gcm
    1071      &           ('phyetat0','pb allocation co2_send',1)
    1072             CALL get_field("co2_send",co2_send,found)
    1073             IF (.NOT. found) THEN
    1074                PRINT*,"phyetat0: Le champ <co2_send> est absent"
    1075                PRINT*,"Initialisation uniforme a co2_ppm=",co2_ppm
    1076                co2_send(:) = co2_ppm
    1077             END IF
    1078          END IF
    1079       END IF
    1080 
    1081 
    1082 c on ferme le fichier
    1083       CALL close_startphy
    1084 
    1085       CALL init_iophy_new(rlat,rlon)
    1086        
    1087 
    1088 c
    1089 c Initialize module pbl_surface_mod
    1090 c
    1091       CALL pbl_surface_init(qsol, fder, snow, qsurf,
    1092      $     evap, frugs, agesno, tsoil)
    1093 
    1094 c Initialize module ocean_cpl_mod for the case of coupled ocean
    1095       IF ( type_ocean == 'couple' ) THEN
    1096          CALL ocean_cpl_init(dtime, rlon, rlat)
    1097       ENDIF
    1098 c
    1099 c Initilialize module fonte_neige_mod     
    1100 c
    1101       CALL fonte_neige_init(run_off_lic_0)
    1102 
    1103 
    1104       RETURN
    1105       END
     367        DO i = 1, klon
     368           xmin = MIN(falb1(i, nsrf), xmin)
     369           xmax = MAX(falb1(i, nsrf), xmax)
     370        ENDDO
     371        PRINT*, 'Albedo du sol ALBE**:', nsrf, xmin, xmax
     372     ENDDO
     373  ELSE
     374     PRINT*, 'phyetat0: Le champ <ALBE> est present'
     375     PRINT*, '          J ignore donc les autres ALBE**'
     376     xmin = 1.0E+20
     377     xmax = -1.0E+20
     378     DO i = 1, klon
     379        xmin = MIN(falb1(i, 1), xmin)
     380        xmax = MAX(falb1(i, 1), xmax)
     381     ENDDO
     382     PRINT*, 'Neige du sol <ALBE>', xmin, xmax
     383     DO nsrf = 2, nbsrf
     384        DO i = 1, klon
     385           falb1(i, nsrf) = falb1(i, 1)
     386        ENDDO
     387     ENDDO
     388  ENDIF
     389
     390  ! Lecture de albedo au sol dans l'interval proche infra-rouge:
     391
     392  CALL get_field("ALBLW", falb2(:, 1), found)
     393  IF (.NOT. found) THEN
     394     PRINT*, 'phyetat0: Le champ <ALBLW> est absent'
     395     PRINT*, '          Mais je vais prendre ALBE**'
     396     DO nsrf = 1, nbsrf
     397        DO i = 1, klon
     398           falb2(i, nsrf) = falb1(i, nsrf)
     399        ENDDO
     400     ENDDO
     401  ELSE
     402     PRINT*, 'phyetat0: Le champ <ALBLW> est present'
     403     PRINT*, '          J ignore donc les autres ALBLW**'
     404     xmin = 1.0E+20
     405     xmax = -1.0E+20
     406     DO i = 1, klon
     407        xmin = MIN(falb2(i, 1), xmin)
     408        xmax = MAX(falb2(i, 1), xmax)
     409     ENDDO
     410     PRINT*, 'Neige du sol <ALBLW>', xmin, xmax
     411     DO nsrf = 2, nbsrf
     412        DO i = 1, klon
     413           falb2(i, nsrf) = falb2(i, 1)
     414        ENDDO
     415     ENDDO
     416  ENDIF
     417
     418  ! Lecture de evaporation: 
     419
     420  CALL get_field("EVAP", evap(:, 1), found)
     421  IF (.NOT. found) THEN
     422     PRINT*, 'phyetat0: Le champ <EVAP> est absent'
     423     PRINT*, '          Mais je vais essayer de lire EVAP**'
     424     DO nsrf = 1, nbsrf
     425        IF (nsrf.GT.99) THEN
     426           PRINT*, "Trop de sous-mailles"
     427           CALL abort
     428        ENDIF
     429        WRITE(str2, '(i2.2)') nsrf
     430        CALL get_field("EVAP"//str2, evap(:, nsrf))
     431        xmin = 1.0E+20
     432        xmax = -1.0E+20
     433        DO i = 1, klon
     434           xmin = MIN(evap(i, nsrf), xmin)
     435           xmax = MAX(evap(i, nsrf), xmax)
     436        ENDDO
     437        PRINT*, 'evap du sol EVAP**:', nsrf, xmin, xmax
     438     ENDDO
     439  ELSE
     440     PRINT*, 'phyetat0: Le champ <EVAP> est present'
     441     PRINT*, '          J ignore donc les autres EVAP**'
     442     xmin = 1.0E+20
     443     xmax = -1.0E+20
     444     DO i = 1, klon
     445        xmin = MIN(evap(i, 1), xmin)
     446        xmax = MAX(evap(i, 1), xmax)
     447     ENDDO
     448     PRINT*, 'Evap du sol <EVAP>', xmin, xmax
     449     DO nsrf = 2, nbsrf
     450        DO i = 1, klon
     451           evap(i, nsrf) = evap(i, 1)
     452        ENDDO
     453     ENDDO
     454  ENDIF
     455
     456  ! Lecture precipitation liquide:
     457
     458  CALL get_field("rain_f", rain_fall)
     459  xmin = 1.0E+20
     460  xmax = -1.0E+20
     461  DO i = 1, klon
     462     xmin = MIN(rain_fall(i), xmin)
     463     xmax = MAX(rain_fall(i), xmax)
     464  ENDDO
     465  PRINT*, 'Precipitation liquide rain_f:', xmin, xmax
     466
     467  ! Lecture precipitation solide:
     468
     469  CALL get_field("snow_f", snow_fall)
     470  xmin = 1.0E+20
     471  xmax = -1.0E+20
     472  DO i = 1, klon
     473     xmin = MIN(snow_fall(i), xmin)
     474     xmax = MAX(snow_fall(i), xmax)
     475  ENDDO
     476  PRINT*, 'Precipitation solide snow_f:', xmin, xmax
     477
     478  ! Lecture rayonnement solaire au sol:
     479
     480  CALL get_field("solsw", solsw, found)
     481  IF (.NOT. found) THEN
     482     PRINT*, 'phyetat0: Le champ <solsw> est absent'
     483     PRINT*, 'mis a zero'
     484     solsw(:) = 0.
     485  ENDIF
     486  xmin = 1.0E+20
     487  xmax = -1.0E+20
     488  DO i = 1, klon
     489     xmin = MIN(solsw(i), xmin)
     490     xmax = MAX(solsw(i), xmax)
     491  ENDDO
     492  PRINT*, 'Rayonnement solaire au sol solsw:', xmin, xmax
     493
     494  ! Lecture rayonnement IF au sol:
     495
     496  CALL get_field("sollw", sollw, found)
     497  IF (.NOT. found) THEN
     498     PRINT*, 'phyetat0: Le champ <sollw> est absent'
     499     PRINT*, 'mis a zero'
     500     sollw = 0.
     501  ENDIF
     502  xmin = 1.0E+20
     503  xmax = -1.0E+20
     504  DO i = 1, klon
     505     xmin = MIN(sollw(i), xmin)
     506     xmax = MAX(sollw(i), xmax)
     507  ENDDO
     508  PRINT*, 'Rayonnement IF au sol sollw:', xmin, xmax
     509
     510  ! Lecture derive des flux:
     511
     512  CALL get_field("fder", fder, found)
     513  IF (.NOT. found) THEN
     514     PRINT*, 'phyetat0: Le champ <fder> est absent'
     515     PRINT*, 'mis a zero'
     516     fder = 0.
     517  ENDIF
     518  xmin = 1.0E+20
     519  xmax = -1.0E+20
     520  DO i = 1, klon
     521     xmin = MIN(fder(i), xmin)
     522     xmax = MAX(fder(i), xmax)
     523  ENDDO
     524  PRINT*, 'Derive des flux fder:', xmin, xmax
     525
     526  ! Lecture du rayonnement net au sol:
     527
     528  CALL get_field("RADS", radsol)
     529  xmin = 1.0E+20
     530  xmax = -1.0E+20
     531  DO i = 1, klon
     532     xmin = MIN(radsol(i), xmin)
     533     xmax = MAX(radsol(i), xmax)
     534  ENDDO
     535  PRINT*, 'Rayonnement net au sol radsol:', xmin, xmax
     536
     537  ! Lecture de la longueur de rugosite
     538
     539  CALL get_field("RUG", frugs(:, 1), found)
     540  IF (.NOT. found) THEN
     541     PRINT*, 'phyetat0: Le champ <RUG> est absent'
     542     PRINT*, '          Mais je vais essayer de lire RUG**'
     543     DO nsrf = 1, nbsrf
     544        IF (nsrf.GT.99) THEN
     545           PRINT*, "Trop de sous-mailles"
     546           CALL abort
     547        ENDIF
     548        WRITE(str2, '(i2.2)') nsrf
     549        CALL get_field("RUG"//str2, frugs(:, nsrf))
     550        xmin = 1.0E+20
     551        xmax = -1.0E+20
     552        DO i = 1, klon
     553           xmin = MIN(frugs(i, nsrf), xmin)
     554           xmax = MAX(frugs(i, nsrf), xmax)
     555        ENDDO
     556        PRINT*, 'rugosite du sol RUG**:', nsrf, xmin, xmax
     557     ENDDO
     558  ELSE
     559     PRINT*, 'phyetat0: Le champ <RUG> est present'
     560     PRINT*, '          J ignore donc les autres RUG**'
     561     xmin = 1.0E+20
     562     xmax = -1.0E+20
     563     DO i = 1, klon
     564        xmin = MIN(frugs(i, 1), xmin)
     565        xmax = MAX(frugs(i, 1), xmax)
     566     ENDDO
     567     PRINT*, 'rugosite <RUG>', xmin, xmax
     568     DO nsrf = 2, nbsrf
     569        DO i = 1, klon
     570           frugs(i, nsrf) = frugs(i, 1)
     571        ENDDO
     572     ENDDO
     573  ENDIF
     574
     575  ! Lecture de l'age de la neige:
     576
     577  CALL get_field("AGESNO", agesno(:, 1), found)
     578  IF (.NOT. found) THEN
     579     PRINT*, 'phyetat0: Le champ <AGESNO> est absent'
     580     PRINT*, '          Mais je vais essayer de lire AGESNO**'
     581     DO nsrf = 1, nbsrf
     582        IF (nsrf.GT.99) THEN
     583           PRINT*, "Trop de sous-mailles"
     584           CALL abort
     585        ENDIF
     586        WRITE(str2, '(i2.2)') nsrf
     587        CALL get_field("AGESNO"//str2, agesno(:, nsrf), found)
     588        IF (.NOT. found) THEN
     589           PRINT*, "phyetat0: Le champ <AGESNO"//str2//"> est absent"
     590           agesno = 50.0
     591        ENDIF
     592        xmin = 1.0E+20
     593        xmax = -1.0E+20
     594        DO i = 1, klon
     595           xmin = MIN(agesno(i, nsrf), xmin)
     596           xmax = MAX(agesno(i, nsrf), xmax)
     597        ENDDO
     598        PRINT*, 'Age de la neige AGESNO**:', nsrf, xmin, xmax
     599     ENDDO
     600  ELSE
     601     PRINT*, 'phyetat0: Le champ <AGESNO> est present'
     602     PRINT*, '          J ignore donc les autres AGESNO**'
     603     xmin = 1.0E+20
     604     xmax = -1.0E+20
     605     DO i = 1, klon
     606        xmin = MIN(agesno(i, 1), xmin)
     607        xmax = MAX(agesno(i, 1), xmax)
     608     ENDDO
     609     PRINT*, 'Age de la neige <AGESNO>', xmin, xmax
     610     DO nsrf = 2, nbsrf
     611        DO i = 1, klon
     612           agesno(i, nsrf) = agesno(i, 1)
     613        ENDDO
     614     ENDDO
     615  ENDIF
     616
     617  CALL get_field("ZMEA", zmea)
     618  xmin = 1.0E+20
     619  xmax = -1.0E+20
     620  DO i = 1, klon
     621     xmin = MIN(zmea(i), xmin)
     622     xmax = MAX(zmea(i), xmax)
     623  ENDDO
     624  PRINT*, 'OROGRAPHIE SOUS-MAILLE zmea:', xmin, xmax
     625
     626  CALL get_field("ZSTD", zstd)
     627  xmin = 1.0E+20
     628  xmax = -1.0E+20
     629  DO i = 1, klon
     630     xmin = MIN(zstd(i), xmin)
     631     xmax = MAX(zstd(i), xmax)
     632  ENDDO
     633  PRINT*, 'OROGRAPHIE SOUS-MAILLE zstd:', xmin, xmax
     634
     635  CALL get_field("ZSIG", zsig)
     636  xmin = 1.0E+20
     637  xmax = -1.0E+20
     638  DO i = 1, klon
     639     xmin = MIN(zsig(i), xmin)
     640     xmax = MAX(zsig(i), xmax)
     641  ENDDO
     642  PRINT*, 'OROGRAPHIE SOUS-MAILLE zsig:', xmin, xmax
     643
     644  CALL get_field("ZGAM", zgam)
     645  xmin = 1.0E+20
     646  xmax = -1.0E+20
     647  DO i = 1, klon
     648     xmin = MIN(zgam(i), xmin)
     649     xmax = MAX(zgam(i), xmax)
     650  ENDDO
     651  PRINT*, 'OROGRAPHIE SOUS-MAILLE zgam:', xmin, xmax
     652
     653  CALL get_field("ZTHE", zthe)
     654  xmin = 1.0E+20
     655  xmax = -1.0E+20
     656  DO i = 1, klon
     657     xmin = MIN(zthe(i), xmin)
     658     xmax = MAX(zthe(i), xmax)
     659  ENDDO
     660  PRINT*, 'OROGRAPHIE SOUS-MAILLE zthe:', xmin, xmax
     661
     662  CALL get_field("ZPIC", zpic)
     663  xmin = 1.0E+20
     664  xmax = -1.0E+20
     665  DO i = 1, klon
     666     xmin = MIN(zpic(i), xmin)
     667     xmax = MAX(zpic(i), xmax)
     668  ENDDO
     669  PRINT*, 'OROGRAPHIE SOUS-MAILLE zpic:', xmin, xmax
     670
     671  CALL get_field("ZVAL", zval)
     672  xmin = 1.0E+20
     673  xmax = -1.0E+20
     674  DO i = 1, klon
     675     xmin = MIN(zval(i), xmin)
     676     xmax = MAX(zval(i), xmax)
     677  ENDDO
     678  PRINT*, 'OROGRAPHIE SOUS-MAILLE zval:', xmin, xmax
     679
     680  CALL get_field("RUGSREL", rugoro)
     681  xmin = 1.0E+20
     682  xmax = -1.0E+20
     683  DO i = 1, klon
     684     xmin = MIN(rugoro(i), xmin)
     685     xmax = MAX(rugoro(i), xmax)
     686  ENDDO
     687  PRINT*, 'Rugosite relief (ecart-type) rugsrel:', xmin, xmax
     688
     689  ancien_ok = .TRUE.
     690
     691  CALL get_field("TANCIEN", t_ancien, found)
     692  IF (.NOT. found) THEN
     693     PRINT*, "phyetat0: Le champ <TANCIEN> est absent"
     694     PRINT*, "Depart legerement fausse. Mais je continue"
     695     ancien_ok = .FALSE.
     696  ENDIF
     697
     698  CALL get_field("QANCIEN", q_ancien, found)
     699  IF (.NOT. found) THEN
     700     PRINT*, "phyetat0: Le champ <QANCIEN> est absent"
     701     PRINT*, "Depart legerement fausse. Mais je continue"
     702     ancien_ok = .FALSE.
     703  ENDIF
     704
     705  CALL get_field("UANCIEN", u_ancien, found)
     706  IF (.NOT. found) THEN
     707     PRINT*, "phyetat0: Le champ <UANCIEN> est absent"
     708     PRINT*, "Depart legerement fausse. Mais je continue"
     709     ancien_ok = .FALSE.
     710  ENDIF
     711
     712  CALL get_field("VANCIEN", v_ancien, found)
     713  IF (.NOT. found) THEN
     714     PRINT*, "phyetat0: Le champ <VANCIEN> est absent"
     715     PRINT*, "Depart legerement fausse. Mais je continue"
     716     ancien_ok = .FALSE.
     717  ENDIF
     718
     719  clwcon=0.
     720  CALL get_field("CLWCON", clwcon, found)
     721  IF (.NOT. found) THEN
     722     PRINT*, "phyetat0: Le champ CLWCON est absent"
     723     PRINT*, "Depart legerement fausse. Mais je continue"
     724  ENDIF
     725  xmin = 1.0E+20
     726  xmax = -1.0E+20
     727  xmin = MINval(clwcon)
     728  xmax = MAXval(clwcon)
     729  PRINT*, 'Eau liquide convective (ecart-type) clwcon:', xmin, xmax
     730
     731  rnebcon = 0.
     732  CALL get_field("RNEBCON", rnebcon, found)
     733  IF (.NOT. found) THEN
     734     PRINT*, "phyetat0: Le champ RNEBCON est absent"
     735     PRINT*, "Depart legerement fausse. Mais je continue"
     736  ENDIF
     737  xmin = 1.0E+20
     738  xmax = -1.0E+20
     739  xmin = MINval(rnebcon)
     740  xmax = MAXval(rnebcon)
     741  PRINT*, 'Nebulosite convective (ecart-type) rnebcon:', xmin, xmax
     742
     743  ! Lecture ratqs
     744
     745  ratqs=0.
     746  CALL get_field("RATQS", ratqs, found)
     747  IF (.NOT. found) THEN
     748     PRINT*, "phyetat0: Le champ <RATQS> est absent"
     749     PRINT*, "Depart legerement fausse. Mais je continue"
     750  ENDIF
     751  xmin = 1.0E+20
     752  xmax = -1.0E+20
     753  xmin = MINval(ratqs)
     754  xmax = MAXval(ratqs)
     755  PRINT*, '(ecart-type) ratqs:', xmin, xmax
     756
     757  ! Lecture run_off_lic_0
     758
     759  CALL get_field("RUNOFFLIC0", run_off_lic_0, found)
     760  IF (.NOT. found) THEN
     761     PRINT*, "phyetat0: Le champ <RUNOFFLIC0> est absent"
     762     PRINT*, "Depart legerement fausse. Mais je continue"
     763     run_off_lic_0 = 0.
     764  ENDIF
     765  xmin = 1.0E+20
     766  xmax = -1.0E+20
     767  xmin = MINval(run_off_lic_0)
     768  xmax = MAXval(run_off_lic_0)
     769  PRINT*, '(ecart-type) run_off_lic_0:', xmin, xmax
     770
     771  ! Lecture de l'energie cinetique turbulente
     772
     773  IF (iflag_pbl>1) then
     774     DO nsrf = 1, nbsrf
     775        IF (nsrf.GT.99) THEN
     776           PRINT*, "Trop de sous-mailles"
     777           CALL abort
     778        ENDIF
     779        WRITE(str2, '(i2.2)') nsrf
     780        CALL get_field("TKE"//str2, pbl_tke(:, 1:klev+1, nsrf), found)
     781        IF (.NOT. found) THEN
     782           PRINT*, "phyetat0: <TKE"//str2//"> est absent"
     783           pbl_tke(:, :, nsrf)=1.e-8
     784        ENDIF
     785        xmin = 1.0E+20
     786        xmax = -1.0E+20
     787        DO k = 1, klev+1
     788           DO i = 1, klon
     789              xmin = MIN(pbl_tke(i, k, nsrf), xmin)
     790              xmax = MAX(pbl_tke(i, k, nsrf), xmax)
     791           ENDDO
     792        ENDDO
     793        PRINT*, 'Temperature du sol TKE**:', nsrf, xmin, xmax
     794     ENDDO
     795  ENDIF
     796
     797  ! zmax0
     798  CALL get_field("ZMAX0", zmax0, found)
     799  IF (.NOT. found) THEN
     800     PRINT*, "phyetat0: Le champ <ZMAX0> est absent"
     801     PRINT*, "Depart legerement fausse. Mais je continue"
     802     zmax0=40.
     803  ENDIF
     804  xmin = 1.0E+20
     805  xmax = -1.0E+20
     806  xmin = MINval(zmax0)
     807  xmax = MAXval(zmax0)
     808  PRINT*, '(ecart-type) zmax0:', xmin, xmax
     809
     810  !           f0(ig)=1.e-5
     811  ! f0
     812  CALL get_field("F0", f0, found)
     813  IF (.NOT. found) THEN
     814     PRINT*, "phyetat0: Le champ <f0> est absent"
     815     PRINT*, "Depart legerement fausse. Mais je continue"
     816     f0=1.e-5
     817  ENDIF
     818  xmin = 1.0E+20
     819  xmax = -1.0E+20
     820  xmin = MINval(f0)
     821  xmax = MAXval(f0)
     822  PRINT*, '(ecart-type) f0:', xmin, xmax
     823
     824  ! sig1 or ema_work1
     825
     826  CALL get_field("sig1", sig1, found)
     827  IF (.NOT. found) CALL get_field("EMA_WORK1", sig1, found)
     828  IF (.NOT. found) THEN
     829     PRINT*, "phyetat0: Le champ sig1 est absent"
     830     PRINT*, "Depart legerement fausse. Mais je continue"
     831     sig1=0.
     832  ELSE
     833     xmin = 1.0E+20
     834     xmax = -1.0E+20
     835     DO k = 1, klev
     836        DO i = 1, klon
     837           xmin = MIN(sig1(i, k), xmin)
     838           xmax = MAX(sig1(i, k), xmax)
     839        ENDDO
     840     ENDDO
     841     PRINT*, 'sig1:', xmin, xmax
     842  ENDIF
     843
     844  ! w01 or ema_work2
     845
     846  CALL get_field("w01", w01, found)
     847  IF (.NOT. found) CALL get_field("EMA_WORK2", w01, found)
     848  IF (.NOT. found) THEN
     849     PRINT*, "phyetat0: Le champ w01 est absent"
     850     PRINT*, "Depart legerement fausse. Mais je continue"
     851     w01=0.
     852  ELSE
     853     xmin = 1.0E+20
     854     xmax = -1.0E+20
     855     DO k = 1, klev
     856        DO i = 1, klon
     857           xmin = MIN(w01(i, k), xmin)
     858           xmax = MAX(w01(i, k), xmax)
     859        ENDDO
     860     ENDDO
     861     PRINT*, 'w01:', xmin, xmax
     862  ENDIF
     863
     864  ! wake_deltat
     865
     866  CALL get_field("WAKE_DELTAT", wake_deltat, found)
     867  IF (.NOT. found) THEN
     868     PRINT*, "phyetat0: Le champ <WAKE_DELTAT> est absent"
     869     PRINT*, "Depart legerement fausse. Mais je continue"
     870     wake_deltat=0.
     871  ELSE
     872     xmin = 1.0E+20
     873     xmax = -1.0E+20
     874     DO k = 1, klev
     875        DO i = 1, klon
     876           xmin = MIN(wake_deltat(i, k), xmin)
     877           xmax = MAX(wake_deltat(i, k), xmax)
     878        ENDDO
     879     ENDDO
     880     PRINT*, 'wake_deltat:', xmin, xmax
     881  ENDIF
     882
     883  ! wake_deltaq
     884
     885  CALL get_field("WAKE_DELTAQ", wake_deltaq, found)
     886  IF (.NOT. found) THEN
     887     PRINT*, "phyetat0: Le champ <WAKE_DELTAQ> est absent"
     888     PRINT*, "Depart legerement fausse. Mais je continue"
     889     wake_deltaq=0.
     890  ELSE
     891     xmin = 1.0E+20
     892     xmax = -1.0E+20
     893     DO k = 1, klev
     894        DO i = 1, klon
     895           xmin = MIN(wake_deltaq(i, k), xmin)
     896           xmax = MAX(wake_deltaq(i, k), xmax)
     897        ENDDO
     898     ENDDO
     899     PRINT*, 'wake_deltaq:', xmin, xmax
     900  ENDIF
     901
     902  ! wake_s
     903
     904  CALL get_field("WAKE_S", wake_s, found)
     905  IF (.NOT. found) THEN
     906     PRINT*, "phyetat0: Le champ <WAKE_S> est absent"
     907     PRINT*, "Depart legerement fausse. Mais je continue"
     908     wake_s=0.
     909  ENDIF
     910  xmin = 1.0E+20
     911  xmax = -1.0E+20
     912  xmin = MINval(wake_s)
     913  xmax = MAXval(wake_s)
     914  PRINT*, '(ecart-type) wake_s:', xmin, xmax
     915
     916  ! wake_cstar
     917
     918  CALL get_field("WAKE_CSTAR", wake_cstar, found)
     919  IF (.NOT. found) THEN
     920     PRINT*, "phyetat0: Le champ <WAKE_CSTAR> est absent"
     921     PRINT*, "Depart legerement fausse. Mais je continue"
     922     wake_cstar=0.
     923  ENDIF
     924  xmin = 1.0E+20
     925  xmax = -1.0E+20
     926  xmin = MINval(wake_cstar)
     927  xmax = MAXval(wake_cstar)
     928  PRINT*, '(ecart-type) wake_cstar:', xmin, xmax
     929
     930  ! wake_pe
     931
     932  CALL get_field("WAKE_PE", wake_pe, found)
     933  IF (.NOT. found) THEN
     934     PRINT*, "phyetat0: Le champ <WAKE_PE> est absent"
     935     PRINT*, "Depart legerement fausse. Mais je continue"
     936     wake_pe=0.
     937  ENDIF
     938  xmin = 1.0E+20
     939  xmax = -1.0E+20
     940  xmin = MINval(wake_pe)
     941  xmax = MAXval(wake_pe)
     942  PRINT*, '(ecart-type) wake_pe:', xmin, xmax
     943
     944  ! wake_fip
     945
     946  CALL get_field("WAKE_FIP", wake_fip, found)
     947  IF (.NOT. found) THEN
     948     PRINT*, "phyetat0: Le champ <WAKE_FIP> est absent"
     949     PRINT*, "Depart legerement fausse. Mais je continue"
     950     wake_fip=0.
     951  ENDIF
     952  xmin = 1.0E+20
     953  xmax = -1.0E+20
     954  xmin = MINval(wake_fip)
     955  xmax = MAXval(wake_fip)
     956  PRINT*, '(ecart-type) wake_fip:', xmin, xmax
     957
     958  !  thermiques
     959
     960  CALL get_field("FM_THERM", fm_therm, found)
     961  IF (.NOT. found) THEN
     962     PRINT*, "phyetat0: Le champ <fm_therm> est absent"
     963     PRINT*, "Depart legerement fausse. Mais je continue"
     964     fm_therm=0.
     965  ENDIF
     966  xmin = 1.0E+20
     967  xmax = -1.0E+20
     968  xmin = MINval(fm_therm)
     969  xmax = MAXval(fm_therm)
     970  PRINT*, '(ecart-type) fm_therm:', xmin, xmax
     971
     972  CALL get_field("ENTR_THERM", entr_therm, found)
     973  IF (.NOT. found) THEN
     974     PRINT*, "phyetat0: Le champ <entr_therm> est absent"
     975     PRINT*, "Depart legerement fausse. Mais je continue"
     976     entr_therm=0.
     977  ENDIF
     978  xmin = 1.0E+20
     979  xmax = -1.0E+20
     980  xmin = MINval(entr_therm)
     981  xmax = MAXval(entr_therm)
     982  PRINT*, '(ecart-type) entr_therm:', xmin, xmax
     983
     984  CALL get_field("DETR_THERM", detr_therm, found)
     985  IF (.NOT. found) THEN
     986     PRINT*, "phyetat0: Le champ <detr_therm> est absent"
     987     PRINT*, "Depart legerement fausse. Mais je continue"
     988     detr_therm=0.
     989  ENDIF
     990  xmin = 1.0E+20
     991  xmax = -1.0E+20
     992  xmin = MINval(detr_therm)
     993  xmax = MAXval(detr_therm)
     994  PRINT*, '(ecart-type) detr_therm:', xmin, xmax
     995
     996  ! Read and send field trs to traclmdz
     997
     998  IF (type_trac == 'lmdz') THEN
     999     DO it=1, nbtr
     1000        iiq=niadv(it+2)
     1001        CALL get_field("trs_"//tname(iiq), trs(:, it), found)
     1002        IF (.NOT. found) THEN
     1003           PRINT*,  &
     1004                "phyetat0: Le champ <trs_"//tname(iiq)//"> est absent"
     1005           PRINT*, "Depart legerement fausse. Mais je continue"
     1006           trs(:, it) = 0.
     1007        ENDIF
     1008        xmin = 1.0E+20
     1009        xmax = -1.0E+20
     1010        xmin = MINval(trs(:, it))
     1011        xmax = MAXval(trs(:, it))
     1012        PRINT*, "(ecart-type) trs_"//tname(iiq)//" :", xmin, xmax
     1013
     1014     END DO
     1015     CALL traclmdz_from_restart(trs)
     1016
     1017     IF (carbon_cycle_cpl) THEN
     1018        ALLOCATE(co2_send(klon), stat=ierr)
     1019        IF (ierr /= 0) CALL abort_gcm &
     1020             ('phyetat0', 'pb allocation co2_send', 1)
     1021        CALL get_field("co2_send", co2_send, found)
     1022        IF (.NOT. found) THEN
     1023           PRINT*, "phyetat0: Le champ <co2_send> est absent"
     1024           PRINT*, "Initialisation uniforme a co2_ppm=", co2_ppm
     1025           co2_send(:) = co2_ppm
     1026        END IF
     1027     END IF
     1028  END IF
     1029
     1030  ! on ferme le fichier
     1031  CALL close_startphy
     1032
     1033  CALL init_iophy_new(rlat, rlon)
     1034
     1035  ! Initialize module pbl_surface_mod
     1036
     1037  CALL pbl_surface_init(qsol, fder, snow, qsurf, &
     1038       evap, frugs, agesno, tsoil)
     1039
     1040  ! Initialize module ocean_cpl_mod for the case of coupled ocean
     1041  IF ( type_ocean == 'couple' ) THEN
     1042     CALL ocean_cpl_init(dtime, rlon, rlat)
     1043  ENDIF
     1044
     1045  ! Initilialize module fonte_neige_mod     
     1046
     1047  CALL fonte_neige_init(run_off_lic_0)
     1048
     1049END SUBROUTINE phyetat0
Note: See TracChangeset for help on using the changeset viewer.