Changeset 2293 for LMDZ5/trunk


Ignore:
Timestamp:
Jun 5, 2015, 9:16:07 PM (10 years ago)
Author:
dcugnet
Message:

Initial states creation routines have been reorganized and simplified.
As far as possible, dynamics and physics related routines have been
separated.
Some routines have been converted to fortran 90 and repeated codes sections
have been "factorized".
Array/vector arguments have become implicit in some routines to avoid usage
of "dimensions.h" ; possible for routines with explicit interfaces and if
iim and jjm can be deduced from arguments sizes.

  • dynlonlat_phylonlat/ce0l.F90 calls now phylmd/etat0phys_netcdf.F90 and dyn3d/etat0dyn_netcdf.F90 that replace phylmd/etat0_netcdf.F90. start.nc and startphy.nc creations are now independant.
  • startvar.F90 has been suppressed ; corresponding operations have been simplified and embedded in etat0*_netcdf.F90 routines as internal procedures.
  • Routines converted to fortran 90 and "factorized":
    • dyn3d_common/conf_dat_m.F90 (replaces dyn3d_common/conf_dat2d.F

and dyn3d_common/conf_dat3d.F)

  • dyn3d/dynredem.F90 (replaces dyn3d/dynredem.F)
  • dyn3d/dynetat0.F90 (replaces dyn3d/dynetat0.F)
  • phylmd/grid_noro_m.F90 (replaces dyn3d_common/grid_noro.F)
  • dynlonlat_phylonlat/grid_atob_m.F90 (replaces dyn3d_common/grid_atob.F)
  • dyn3d_common/caldyn0.F90 (replaces dyn3d_common/caldyn0.F)
  • dyn3d_common/covcont.F90 (replaces dyn3d_common/covcont.F)
  • dyn3d_common/pression.F90 (replaces dyn3d_common/pression.F)
  • phylmd/phyredem.F90 and phylmd/limit_netcdf.F90 have been slightly factorized.

TO DO:

  • little fix needed in grid_noro_m.F90 ; untouched yet to ensure results are exactly the same as before. Unsmoothed orography is used to compute "zphi", but smoothed (should be unsmoothed) one is used at poles.
  • add the dyn3dmem versions of dynredem.F90 and dynetat0.F90 (dynredem_loc.F90 and dynetat0_loc.F90, untested yet).
  • test compilation in parallel mode for a single processor.
Location:
LMDZ5/trunk/libf
Files:
10 added
11 deleted
3 edited

Legend:

Unmodified
Added
Removed

  • LMDZ5/trunk/libf/dynlonlat_phylonlat/phylmd/ce0l.F90

    r2248 r2293  
     1PROGRAM ce0l
    12!
    2 ! $Id$
    3 !
    4 !-------------------------------------------------------------------------------
    5 !
    6 PROGRAM ce0l
    73!-------------------------------------------------------------------------------
    84! Purpose: Calls etat0, creates initial states and limit_netcdf
    … …  
    1511!-------------------------------------------------------------------------------
    1612  USE control_mod
    17 #ifdef CPP_EARTH
    18 ! This prog. is designed to work for Earth
     13
     14!#ifdef CPP_EARTH
     15!! This prog. is designed to work for Earth
     16  USE etat0dyn
     17  USE netcdf, ONLY: NF90_OPEN, NF90_NOWRITE, NF90_CLOSE, NF90_NOERR
     18  USE ioipsl, ONLY: ioconf_calendar, getin, flininfo, flinopen, flinget, flinclo
     19
     20#ifdef CPP_PHYS
     21  USE etat0phys
    1922  USE dimphy
    2023  USE comgeomphy
    2124  USE infotrac
    2225  USE indice_sol_mod
     26  USE test_disvert_m, ONLY: test_disvert
     27#endif
     28! of #ifdef CPP_PHYS
     29
     30!#endif
     31!! of #ifdef CPP_EARTH
     32
     33  IMPLICIT NONE
     34
     35!-------------------------------------------------------------------------------
     36! Local variables:
     37!#ifndef CPP_EARTH
     38!  include "iniprint.h"
     39!!-------------------------------------------------------------------------------
     40!  WRITE(lunout,*)'ce0l: Earth-specific routine, needs Earth physics'
     41!#else
     42  include "dimensions.h"
     43  include "paramet.h"
     44  include "comconst.h"
     45  include "comvert.h"
     46  include "iniprint.h"
     47  include "temps.h"
     48  include "logic.h"
     49  REAL               :: masque(iip1,jjp1)             !--- CONTINENTAL MASK
     50  REAL               :: phis  (iip1,jjp1)             !--- GROUND GEOPOTENTIAL
     51  CHARACTER(LEN=256) :: modname, fmt, calnd           !--- CALENDAR TYPE
     52  LOGICAL            :: use_filtre_fft
     53  LOGICAL, PARAMETER :: interbar=.TRUE., extrap=.FALSE., oldice=.FALSE.
     54
     55!--- Local variables for ocean mask reading:
     56  INTEGER            :: nid_o2a, iml_omask, jml_omask, j
     57  INTEGER            :: fid, iret, llm_tmp, ttm_tmp, itaul(1)
     58  REAL, ALLOCATABLE  :: lon_omask(:,:), dlon_omask(:), ocemask(:,:)
     59  REAL, ALLOCATABLE  :: lat_omask(:,:), dlat_omask(:), ocetmp (:,:)
     60  REAL               :: date, lev(1)
     61!-------------------------------------------------------------------------------
     62  modname="ce0l"
     63
     64!--- Constants
     65  pi     = 4. * ATAN(1.)
     66  rad    = 6371229.
     67  daysec = 86400.
     68  omeg   = 2.*pi/daysec
     69  g      = 9.8
     70  kappa  = 0.2857143
     71  cpp    = 1004.70885
     72  jmp1   = jjm + 1
     73  preff   = 101325.
     74  pa      = 50000.
     75
     76  CALL conf_gcm( 99, .TRUE. )
     77
     78  dtvr = daysec/FLOAT(day_step)
     79  WRITE(lunout,*)'dtvr',dtvr
     80
     81  CALL iniconst()
     82#ifdef CPP_PHYS
     83  IF(pressure_exner) CALL test_disvert
     84#endif
     85  CALL inigeom()
    2386
    2487#ifdef CPP_IOIPSL
    25   USE ioipsl, ONLY: ioconf_calendar, getin
    26 #else
    27   ! if not using IOIPSL, we still need to use (a local version of) getin
    28   use ioipsl_getincom, only: getin
    29 #endif
    30 
    31 #endif
    32   IMPLICIT NONE
    33 #ifndef CPP_EARTH
    34 #include "iniprint.h"
    35   WRITE(lunout,*)'limit_netcdf: Earth-specific routine, needs Earth physics'
    36 #else
    37 !-------------------------------------------------------------------------------
    38 ! Local variables:
    39   LOGICAL, PARAMETER :: interbar=.TRUE., extrap=.FALSE., oldice=.FALSE.
    40 #include "dimensions.h"
    41 #include "paramet.h"
    42 !#include "indicesol.h"
    43 #include "iniprint.h"
    44 #include "temps.h"
    45 #include "logic.h"
    46   REAL,    DIMENSION(iip1,jjp1) :: masque
    47   CHARACTER(LEN=15)             :: calnd
    48   REAL,    DIMENSION(iip1,jjp1) :: phis ! geopotentiel au sol
    49   logical use_filtre_fft
    50 !-------------------------------------------------------------------------------
    51   CALL conf_gcm( 99, .TRUE. )
     88  calnd='gregorian'
     89  SELECT CASE(calend)
     90    CASE('earth_360d');CALL ioconf_calendar('360d');   calnd='with 360 days/year'
     91    CASE('earth_365d');CALL ioconf_calendar('noleap'); calnd='with no leap year'
     92    CASE('earth_366d');CALL ioconf_calendar('366d');   calnd='with leap years only'
     93    CASE('gregorian'); CALL ioconf_calendar('gregorian')
     94    CASE('standard');  CALL ioconf_calendar('gregorian')
     95    CASE('julian');    CALL ioconf_calendar('julian'); calnd='julian'
     96    CASE('proleptic_gregorian'); CALL ioconf_calendar('gregorian')
     97  !--- DC Bof...  => IOIPSL a mettre a jour: proleptic_gregorian /= gregorian
     98    CASE DEFAULT
     99      CALL abort_gcm('ce0l','Bad choice for calendar',1)
     100  END SELECT
     101  WRITE(lunout,*)'CHOSEN CALENDAR: Earth '//TRIM(calnd)
     102#endif
    52103
    53104  use_filtre_fft=.FALSE.
    54105  CALL getin('use_filtre_fft',use_filtre_fft)
    55   IF (use_filtre_fft) THEN
    56      write(lunout, fmt = *) 'FFT filter is not available in the ' &
    57           // 'sequential version of the dynamics.'
    58      write(lunout, fmt = *) &
    59           "Your setting of variable use_filtre_fft is not used."
     106  IF(use_filtre_fft) THEN
     107     WRITE(lunout,*)"FFT filter not available for sequential dynamics."
     108     WRITE(lunout,*)"Your setting of variable use_filtre_fft is not used."
    60109  ENDIF
    61110
     111!--- LAND MASK. TWO CASES:
     112!   1) read from ocean model    file "o2a.nc"    (coupled runs)
     113!   2) computed from topography file "Relief.nc" (masque(:,:)=-99999.)
     114! Coupled simulations (case 1) use the ocean model mask to compute the
     115! weights to ensure ocean fractions are the same for atmosphere and ocean.
     116!*******************************************************************************
     117  IF(NF90_OPEN("o2a.nc", NF90_NOWRITE, nid_o2a)/=NF90_NOERR) THEN
     118    WRITE(lunout,*)'BEWARE !! No ocean mask "o2a.nc" file found'
     119    WRITE(lunout,*)'Forced run.'
     120    masque(:,:)=-99999.
     121  ELSE
     122    iret=NF90_CLOSE(nid_o2a)
     123    WRITE(lunout,*)'BEWARE !! Ocean mask "o2a.nc" file found'
     124    WRITE(lunout,*)'Coupled run.'
     125    CALL flininfo("o2a.nc", iml_omask, jml_omask, llm_tmp, ttm_tmp, nid_o2a)
     126    IF(iml_omask/=iim .OR.jml_omask/=jjp1) THEN
     127      WRITE(lunout,*)'Mismatching dimensions for ocean mask'
     128      WRITE(lunout,*)'iim  = ',iim ,' iml_omask = ',iml_omask
     129      WRITE(lunout,*)'jjp1 = ',jjp1,' jml_omask = ',jml_omask
     130      CALL abort_gcm(modname,'',1)
     131    END IF
     132    ALLOCATE(ocemask(iim,jjp1),lon_omask(iim,jjp1),dlon_omask(iim ))
     133    ALLOCATE(ocetmp (iim,jjp1),lat_omask(iim,jjp1),dlat_omask(jjp1))
     134    CALL flinopen("o2a.nc", .FALSE.,iim,jjp1,llm_tmp,lon_omask,lat_omask,lev,  &
     135     &                                       ttm_tmp,itaul,date,dt,fid)
     136    CALL flinget(fid, "OceMask",    iim,jjp1,llm_tmp,ttm_tmp,1,1,ocetmp)
     137    CALL flinclo(fid)
     138    dlon_omask(1:iim ) = lon_omask(1:iim,1)
     139    dlat_omask(1:jjp1) = lat_omask(1,1:jjp1)
     140    ocemask = ocetmp
     141    IF(dlat_omask(1)<dlat_omask(jml_omask)) THEN
     142      DO j=1,jjp1; ocemask(:,j) = ocetmp(:,jjp1-j+1); END DO
     143    END IF
     144    DEALLOCATE(ocetmp,lon_omask,lat_omask,dlon_omask,dlat_omask)
     145    IF(prt_level>=1) THEN
     146      WRITE(fmt,"(i4,'i1)')")iim ; fmt='('//ADJUSTL(fmt)
     147      WRITE(lunout,*)'OCEAN MASK :'
     148      WRITE(lunout,fmt) NINT(ocemask)
     149    END IF
     150    masque(1:iim,:)=1.-ocemask(:,:)
     151    masque(iip1 ,:)=masque(1,:)
     152    DEALLOCATE(ocemask)
     153  END IF
     154  phis(:,:)=-99999.
     155
     156#ifdef CPP_PHYS
    62157  CALL Init_Phys_lmdz(iim,jjp1,llm,1,(/(jjm-1)*iim+2/))
    63158  WRITE(lunout,*)'---> klon=',klon
    64159  CALL InitComgeomphy
    65 
    66 #ifdef CPP_IOIPSL
    67   SELECT CASE(calend)
    68     CASE('earth_360d');CALL ioconf_calendar('360d');      calnd='a 360 jours/an'
    69     CASE('earth_365d');CALL ioconf_calendar('noleap');    calnd='a 365 jours/an'
    70     CASE('earth_366d');CALL ioconf_calendar('366d');      calnd='bissextile'
    71     CASE('gregorian'); CALL ioconf_calendar('gregorian'); calnd='gregorien'
    72     CASE('standard');  CALL ioconf_calendar('gregorian'); calnd='gregorien'
    73     CASE('julian');    CALL ioconf_calendar('julian');    calnd='julien'
    74     CASE('proleptic_gregorian'); CALL ioconf_calendar('gregorian')
    75   !--- DC Bof...  => IOIPSL a mettre a jour: proleptic_gregorian /= gregorian
    76     CASE DEFAULT
    77       CALL abort_gcm('ce0l','Mauvais choix de calendrier',1)
    78   END SELECT
    79   WRITE(lunout,*)'CALENDRIER CHOISI: Terrestre '//TRIM(calnd)
    80 #endif
    81160
    82161  IF (type_trac == 'inca') THEN
    … …  
    86165    WRITE(lunout,*)'nbtr =' , nbtr
    87166#endif
    88   END IF
    89 
    90   WRITE(lunout,'(//)')
    91   WRITE(lunout,*) '  *********************  '
    92   WRITE(lunout,*) '  ***  etat0_netcdf ***  '
    93   WRITE(lunout,*) '  *********************  '
    94   WRITE(lunout,'(//)')
    95   WRITE(lunout,*) ' interbar = ',interbar
    96   CALL etat0_netcdf(interbar,masque,phis,ok_etat0)
    97 
     167! of #ifdef INCA
     168  END IF
     169  IF(ok_etat0) THEN
     170    WRITE(lunout,'(//)')
     171    WRITE(lunout,*) '  ************************  '
     172    WRITE(lunout,*) '  ***  etat0phy_netcdf ***  '
     173    WRITE(lunout,*) '  ************************  '
     174    WRITE(lunout,'(//)')
     175    WRITE(lunout,*) ' interbar = ',interbar
     176    CALL etat0phys_netcdf(interbar,masque,phis)
     177  END IF
     178#endif
     179! of #ifdef CPP_PHYS
     180
     181  IF(ok_etat0) THEN
     182    WRITE(lunout,'(//)')
     183    WRITE(lunout,*) '  ************************  '
     184    WRITE(lunout,*) '  ***  etat0dyn_netcdf ***  '
     185    WRITE(lunout,*) '  ************************  '
     186    WRITE(lunout,'(//)')
     187    WRITE(lunout,*) ' interbar = ',interbar
     188    CALL etat0dyn_netcdf(interbar,masque,phis)
     189  END IF
     190
     191#ifdef CPP_PHYS
    98192  IF(ok_limit) THEN
    99   WRITE(lunout,'(//)')
    100   WRITE(lunout,*) '  *********************  '
    101   WRITE(lunout,*) '  ***  Limit_netcdf ***  '
    102   WRITE(lunout,*) '  *********************  '
    103   WRITE(lunout,'(//)')
    104   CALL limit_netcdf(interbar,extrap,oldice,masque)
    105   END IF
    106 
    107  
     193    WRITE(lunout,'(//)')
     194    WRITE(lunout,*) '  *********************  '
     195    WRITE(lunout,*) '  ***  Limit_netcdf ***  '
     196    WRITE(lunout,*) '  *********************  '
     197    WRITE(lunout,'(//)')
     198    CALL limit_netcdf(interbar,extrap,oldice,masque)
     199  END IF
     200#endif
     201
    108202  WRITE(lunout,'(//)')
    109203  WRITE(lunout,*) '  ***************************  '
    … …  
    113207  CALL grilles_gcm_netcdf_sub(masque,phis)
    114208
    115 #endif
     209!#endif
    116210! of #ifndef CPP_EARTH #else
    117211

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

    r2239 r2293  
    1 !
    2 ! $Id: limit_netcdf.F90 1508 2011-04-15 13:05:34Z jghattas $
    3 !-------------------------------------------------------------------------------
    4 !
    51SUBROUTINE limit_netcdf(interbar, extrap, oldice, masque)
    6 #ifndef CPP_1D
    72!
    83!-------------------------------------------------------------------------------
    … …  
    2116!  *    12/2009: D. Cugnet   (f77->f90, calendars, files from coupled runs)
    2217!-------------------------------------------------------------------------------
     18#ifndef CPP_1D
    2319  USE control_mod
    2420  USE indice_sol_mod
    … …  
    3228                   NF90_CLOBBER, NF90_ENDDEF,  NF90_UNLIMITED, NF90_FLOAT
    3329  USE inter_barxy_m, only: inter_barxy
    34   use netcdf95, only: nf95_def_var, nf95_put_att, nf95_put_var
     30  USE netcdf95,    ONLY: nf95_def_var, nf95_put_att, nf95_put_var
     31  USE grid_atob_m, ONLY: grille_m, rugosite, sea_ice
    3532#endif
    3633  IMPLICIT NONE
    3734!-------------------------------------------------------------------------------
    3835! Arguments:
    39 #include "dimensions.h"
    40 #include "paramet.h"
    41 #include "iniprint.h"
     36  include "dimensions.h"
     37  include "paramet.h"
     38  include "iniprint.h"
    4239  LOGICAL,                    INTENT(IN) :: interbar ! barycentric interpolation
    4340  LOGICAL,                    INTENT(IN) :: extrap   ! SST extrapolation flag
    … …  
    4946!-------------------------------------------------------------------------------
    5047! Local variables:
    51 #include "logic.h"
    52 #include "comvert.h"
    53 #include "comgeom2.h"
    54 #include "comconst.h"
     48  include "logic.h"
     49  include "comvert.h"
     50  include "comgeom2.h"
     51  include "comconst.h"
    5552
    5653!--- INPUT NETCDF FILES NAMES --------------------------------------------------
    57   CHARACTER(LEN=25) :: icefile, sstfile, dumstr
    58   CHARACTER(LEN=25), PARAMETER :: famipsst ='amipbc_sst_1x1.nc        ',        &
    59                                   famipsic ='amipbc_sic_1x1.nc        ',        &
    60                                   fcpldsst ='cpl_atm_sst.nc           ',        &
    61                                   fcpldsic ='cpl_atm_sic.nc           ',        &
    62                                   fhistsst ='histmth_sst.nc           ',        &
    63                                   fhistsic ='histmth_sic.nc           ',        &
    64                                   frugo    ='Rugos.nc                 ',        &
    65                                   falbe    ='Albedo.nc                ',        &
    66                                   feraisstk='sstk.nc                  ',        &
    67                                   feraici  ='ci.nc                    '
     54  CHARACTER(LEN=20) :: icefile, sstfile, dumstr, fnam
     55  CHARACTER(LEN=20), PARAMETER :: &
     56  fsst(4)=['amipbc_sst_1x1.nc   ','cpl_atm_sst.nc      ','histmth_sst.nc      '&
     57          ,'sstk.nc             ']
     58  CHARACTER(LEN=20), PARAMETER :: &
     59  fsic(4)=['amipbc_sic_1x1.nc   ','cpl_atm_sic.nc      ','histmth_sic.nc      '&
     60          ,'ci.nc               ']
     61  CHARACTER(LEN=10), PARAMETER :: &
     62  vsst(4)=['tosbcs    ','SISUTESW  ','tsol_oce  ','sstk      '], &
     63  vsic(4)=['sicbcs    ','SIICECOV  ','pourc_sic ','ci        ']
     64  CHARACTER(LEN=20), PARAMETER :: frugo='Rugos.nc            ',                &
     65                                  falbe='Albedo.nc           '
     66  CHARACTER(LEN=10), PARAMETER :: vrug='RUGOS     ', valb='ALBEDO    '
    6867  CHARACTER(LEN=10) :: varname
    6968!--- OUTPUT VARIABLES FOR NETCDF FILE ------------------------------------------
    70   REAL,   DIMENSION(klon)                :: fi_ice, verif
    71   REAL,   DIMENSION(:,:),   POINTER      :: phy_rug=>NULL(), phy_ice=>NULL()
    72   REAL,   DIMENSION(:,:),   POINTER      :: phy_sst=>NULL(), phy_alb=>NULL()
    73   REAL,   DIMENSION(:,:),   ALLOCATABLE  :: phy_bil
    74   REAL,   DIMENSION(:,:,:), ALLOCATABLE  :: pctsrf_t
    75   INTEGER                                :: nbad
     69  REAL               :: fi_ice(klon), verif(klon)
     70  REAL, POINTER      :: phy_rug(:,:)=>NULL(), phy_ice(:,:)=>NULL()
     71  REAL, POINTER      :: phy_sst(:,:)=>NULL(), phy_alb(:,:)=>NULL()
     72  REAL, ALLOCATABLE  :: phy_bil(:,:), pctsrf_t(:,:,:)
     73  INTEGER            :: nbad
    7674
    7775!--- VARIABLES FOR OUTPUT FILE WRITING -----------------------------------------
    78   INTEGER :: ierr, nid, ndim, ntim, k
    79   INTEGER, DIMENSION(2) :: dims
     76  INTEGER :: ierr, nid, ndim, ntim, k, dims(2), ix_sic, ix_sst
    8077  INTEGER :: id_tim,  id_SST,  id_BILS, id_RUG, id_ALB
    8178  INTEGER :: id_FOCE, id_FSIC, id_FTER, id_FLIC, varid_longitude, varid_latitude
    … …  
    8986  NF90_FORMAT=NF90_FLOAT
    9087#endif
    91 
    92   pi    = 4.*ATAN(1.)
    93   rad   = 6371229.
    94   daysec= 86400.
    95   omeg  = 2.*pi/daysec
    96   g     = 9.8
    97   kappa = 0.2857143
    98   cpp   = 1004.70885
    99   dtvr  = daysec/REAL(day_step)
    10088  CALL inigeom
    10189
    … …  
    10492
    10593!--- RUGOSITY TREATMENT --------------------------------------------------------
    106   IF (prt_level>1) WRITE(lunout,*) 'Traitement de la rugosite'
    107   varname='RUGOS'
    108   CALL get_2Dfield(frugo,varname,'RUG',interbar,ndays,phy_rug,mask=masque(1:iim,:))
     94  CALL msg(1,'Traitement de la rugosite')
     95  CALL get_2Dfield(frugo,vrug,'RUG',interbar,ndays,phy_rug,mask=masque(1:iim,:))
    10996
    11097!--- OCEAN TREATMENT -----------------------------------------------------------
    111   IF (prt_level>1) WRITE(lunout,*) 'Traitement de la glace oceanique'
     98  CALL msg(1,'Traitement de la glace oceanique')
    11299
    113100! Input SIC file selection
    114101! Open file only to test if available
    115   IF ( NF90_OPEN(TRIM(famipsic),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
    116      icefile=TRIM(famipsic)
    117      varname='sicbcs'
    118   ELSE IF( NF90_OPEN(TRIM(fcpldsic),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
    119      icefile=TRIM(fcpldsic)
    120      varname='SIICECOV'
    121   ELSE IF ( NF90_OPEN(TRIM(fhistsic),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
    122      icefile=TRIM(fhistsic)
    123      varname='pourc_sic'
    124   ELSE IF ( NF90_OPEN(TRIM(feraici),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
    125      icefile=TRIM(feraici)
    126      varname='ci'
    127   ELSE
     102  DO ix_sic=1,SIZE(fsic)
     103     IF ( NF90_OPEN(TRIM(fsic(ix_sic)),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
     104        icefile=fsic(ix_sic); varname=vsic(ix_sic); EXIT
     105     END IF
     106  END DO
     107  IF(ix_sic==SIZE(fsic)+1) THEN
    128108     WRITE(lunout,*) 'ERROR! No sea-ice input file was found.'
    129      WRITE(lunout,*) 'One of following files must be availible : ',trim(famipsic),', ',trim(fcpldsic),', ', &
    130                       trim(fhistsic), trim(feraici)
     109     WRITE(lunout,*) 'One of following files must be available : '
     110     DO k=1,SIZE(fsic); WRITE(lunout,*) TRIM(fsic(k)); END DO
    131111     CALL abort_gcm('limit_netcdf','No sea-ice file was found',1)
    132112  END IF
    133   ierr=NF90_CLOSE(nid)
    134   IF (prt_level>=0) WRITE(lunout,*)'Pour la glace de mer a ete choisi le fichier '//TRIM(icefile)
     113  CALL ncerr(NF90_CLOSE(nid),icefile)
     114  CALL msg(-1,'Fichier choisi pour la glace de mer:'//TRIM(icefile))
    135115
    136116  CALL get_2Dfield(icefile,varname, 'SIC',interbar,ndays,phy_ice,flag=oldice)
    … …  
    143123     pctsrf_t(:,is_ter,k)=pctsrf(:,is_ter)       ! land soil
    144124     pctsrf_t(:,is_lic,k)=pctsrf(:,is_lic)       ! land ice
    145      IF (icefile==trim(fcpldsic)) THEN           ! SIC=pICE*(1-LIC-TER)
     125     SELECT CASE(ix_sic)
     126        CASE(2)                                   ! SIC=pICE*(1-LIC-TER) (CPL)
    146127        pctsrf_t(:,is_sic,k)=fi_ice(:)*(1.-pctsrf(:,is_lic)-pctsrf(:,is_ter))
    147      ELSE IF (icefile==trim(fhistsic)) THEN      ! SIC=pICE
     128        CASE(3)                                   ! SIC=pICE            (HIST)
    148129        pctsrf_t(:,is_sic,k)=fi_ice(:)
    149      ELSE ! icefile==famipsic                    ! SIC=pICE-LIC
     130        CASE DEFAULT                              ! SIC=pICE-LIC   (AMIP,ERAI)
    150131        pctsrf_t(:,is_sic,k)=fi_ice-pctsrf_t(:,is_lic,k)
    151      END IF
     132     END SELECT
    152133     WHERE(pctsrf_t(:,is_sic,k)<=0) pctsrf_t(:,is_sic,k)=0.
    153134     WHERE(1.0-zmasq<EPSFRA)
    … …  
    167148     END WHERE
    168149     nbad=COUNT(pctsrf_t(:,is_oce,k)<0.0)
    169      IF(nbad>0) WRITE(lunout,*) 'pb sous maille pour nb point = ',nbad
    170      nbad=COUNT(abs(sum(pctsrf_t(:,:,k),dim=2)-1.0)>EPSFRA)
     150     IF(nbad>0) WRITE(lunout,*) 'pb sous maille pour nb points = ',nbad
     151     nbad=COUNT(ABS(SUM(pctsrf_t(:,:,k),DIM=2)-1.0)>EPSFRA)
    171152     IF(nbad>0) WRITE(lunout,*) 'pb sous surface pour nb points = ',nbad
    172153  END DO
    … …  
    174155
    175156!--- SST TREATMENT -------------------------------------------------------------
    176   IF (prt_level>1) WRITE(lunout,*) 'Traitement de la sst'
     157  CALL msg(1,'Traitement de la sst')
    177158
    178159! Input SST file selection
    179160! Open file only to test if available
    180   IF ( NF90_OPEN(TRIM(famipsst),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
    181      sstfile=TRIM(famipsst)
    182      varname='tosbcs'
    183   ELSE IF ( NF90_OPEN(TRIM(fcpldsst),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
    184      sstfile=TRIM(fcpldsst)
    185      varname='SISUTESW'
    186   ELSE IF ( NF90_OPEN(TRIM(fhistsst),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
    187      sstfile=TRIM(fhistsst)
    188      varname='tsol_oce'
    189   ELSE IF ( NF90_OPEN(TRIM(feraisstk),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
    190      sstfile=TRIM(feraisstk)
    191      varname='sstk'
    192   ELSE
     161  DO ix_sst=1,SIZE(fsst)
     162     IF ( NF90_OPEN(TRIM(fsst(ix_sst)),NF90_NOWRITE,nid)==NF90_NOERR ) THEN
     163       sstfile=fsst(ix_sst); varname=vsst(ix_sst); EXIT
     164     END IF
     165  END DO
     166  IF(ix_sst==SIZE(fsst)+1) THEN
    193167     WRITE(lunout,*) 'ERROR! No sst input file was found.'
    194      WRITE(lunout,*) 'One of following files must be availible : ',trim(famipsst),trim(fcpldsst),trim(fhistsst),trim(feraisstk)
     168     WRITE(lunout,*) 'One of following files must be available : '
     169     DO k=1,SIZE(fsst); WRITE(lunout,*) TRIM(fsst(k)); END DO
    195170     CALL abort_gcm('limit_netcdf','No sst file was found',1)
    196171  END IF
    197   ierr=NF90_CLOSE(nid)
    198   IF (prt_level>=0) WRITE(lunout,*)'Pour la temperature de mer a ete choisi un fichier '//TRIM(sstfile)
     172  CALL ncerr(NF90_CLOSE(nid),sstfile)
     173  CALL msg(-1,'Fichier choisi pour la temperature de mer: '//TRIM(sstfile))
    199174
    200175  CALL get_2Dfield(sstfile,varname,'SST',interbar,ndays,phy_sst,flag=extrap)
    201176
    202177!--- ALBEDO TREATMENT ----------------------------------------------------------
    203   IF (prt_level>1) WRITE(lunout,*) 'Traitement de l albedo'
    204   varname='ALBEDO'
    205   CALL get_2Dfield(falbe,varname,'ALB',interbar,ndays,phy_alb)
     178  CALL msg(1,'Traitement de l albedo')
     179  CALL get_2Dfield(falbe,valb,'ALB',interbar,ndays,phy_alb)
    206180
    207181!--- REFERENCE GROUND HEAT FLUX TREATMENT --------------------------------------
    … …  
    209183
    210184!--- OUTPUT FILE WRITING -------------------------------------------------------
    211   IF (prt_level>5) WRITE(lunout,*) 'Ecriture du fichier limit : debut'
     185  CALL msg(5,'Ecriture du fichier limit : debut')
     186  fnam="limit.nc"
    212187
    213188  !--- File creation
    214   ierr=NF90_CREATE("limit.nc",NF90_CLOBBER,nid)
    215   ierr=NF90_PUT_ATT(nid,NF90_GLOBAL,"title","Fichier conditions aux limites")
     189  CALL ncerr(NF90_CREATE(fnam,NF90_CLOBBER,nid),fnam)
     190  CALL ncerr(NF90_PUT_ATT(nid,NF90_GLOBAL,"title","Fichier conditions aux limites"),fnam)
    216191
    217192  !--- Dimensions creation
    218   ierr=NF90_DEF_DIM(nid,"points_physiques",klon,ndim)
    219   ierr=NF90_DEF_DIM(nid,"time",NF90_UNLIMITED,ntim)
    220 
    221   dims=(/ndim,ntim/)
     193  CALL ncerr(NF90_DEF_DIM(nid,"points_physiques",klon,ndim),fnam)
     194  CALL ncerr(NF90_DEF_DIM(nid,"time",NF90_UNLIMITED,ntim),fnam)
     195
     196  dims=[ndim,ntim]
    222197
    223198  !--- Variables creation
    224   ierr=NF90_DEF_VAR(nid,"TEMPS",NF90_FORMAT,(/ntim/),id_tim)
    225   ierr=NF90_DEF_VAR(nid,"FOCE", NF90_FORMAT,dims,id_FOCE)
    226   ierr=NF90_DEF_VAR(nid,"FSIC", NF90_FORMAT,dims,id_FSIC)
    227   ierr=NF90_DEF_VAR(nid,"FTER", NF90_FORMAT,dims,id_FTER)
    228   ierr=NF90_DEF_VAR(nid,"FLIC", NF90_FORMAT,dims,id_FLIC)
    229   ierr=NF90_DEF_VAR(nid,"SST",  NF90_FORMAT,dims,id_SST)
    230   ierr=NF90_DEF_VAR(nid,"BILS", NF90_FORMAT,dims,id_BILS)
    231   ierr=NF90_DEF_VAR(nid,"ALB",  NF90_FORMAT,dims,id_ALB)
    232   ierr=NF90_DEF_VAR(nid,"RUG",  NF90_FORMAT,dims,id_RUG)
     199  CALL ncerr(NF90_DEF_VAR(nid,"TEMPS",NF90_FORMAT,[ntim],id_tim),fnam)
     200  CALL ncerr(NF90_DEF_VAR(nid,"FOCE", NF90_FORMAT,dims,id_FOCE),fnam)
     201  CALL ncerr(NF90_DEF_VAR(nid,"FSIC", NF90_FORMAT,dims,id_FSIC),fnam)
     202  CALL ncerr(NF90_DEF_VAR(nid,"FTER", NF90_FORMAT,dims,id_FTER),fnam)
     203  CALL ncerr(NF90_DEF_VAR(nid,"FLIC", NF90_FORMAT,dims,id_FLIC),fnam)
     204  CALL ncerr(NF90_DEF_VAR(nid,"SST",  NF90_FORMAT,dims,id_SST),fnam)
     205  CALL ncerr(NF90_DEF_VAR(nid,"BILS", NF90_FORMAT,dims,id_BILS),fnam)
     206  CALL ncerr(NF90_DEF_VAR(nid,"ALB",  NF90_FORMAT,dims,id_ALB),fnam)
     207  CALL ncerr(NF90_DEF_VAR(nid,"RUG",  NF90_FORMAT,dims,id_RUG),fnam)
    233208  call nf95_def_var(nid, "longitude", NF90_FLOAT, ndim, varid_longitude)
    234   call nf95_def_var(nid, "latitude", NF90_FLOAT, ndim, varid_latitude)
     209  call nf95_def_var(nid, "latitude",  NF90_FLOAT, ndim, varid_latitude)
    235210
    236211  !--- Attributes creation
    237   ierr=NF90_PUT_ATT(nid,id_tim, "title","Jour dans l annee")
    238   ierr=NF90_PUT_ATT(nid,id_FOCE,"title","Fraction ocean")
    239   ierr=NF90_PUT_ATT(nid,id_FSIC,"title","Fraction glace de mer")
    240   ierr=NF90_PUT_ATT(nid,id_FTER,"title","Fraction terre")
    241   ierr=NF90_PUT_ATT(nid,id_FLIC,"title","Fraction land ice")
    242   ierr=NF90_PUT_ATT(nid,id_SST ,"title","Temperature superficielle de la mer")
    243   ierr=NF90_PUT_ATT(nid,id_BILS,"title","Reference flux de chaleur au sol")
    244   ierr=NF90_PUT_ATT(nid,id_ALB, "title","Albedo a la surface")
    245   ierr=NF90_PUT_ATT(nid,id_RUG, "title","Rugosite")
     212  CALL ncerr(NF90_PUT_ATT(nid,id_tim, "title","Jour dans l annee"),fnam)
     213  CALL ncerr(NF90_PUT_ATT(nid,id_FOCE,"title","Fraction ocean"),fnam)
     214  CALL ncerr(NF90_PUT_ATT(nid,id_FSIC,"title","Fraction glace de mer"),fnam)
     215  CALL ncerr(NF90_PUT_ATT(nid,id_FTER,"title","Fraction terre"),fnam)
     216  CALL ncerr(NF90_PUT_ATT(nid,id_FLIC,"title","Fraction land ice"),fnam)
     217  CALL ncerr(NF90_PUT_ATT(nid,id_SST ,"title","Temperature superficielle de la mer"),fnam)
     218  CALL ncerr(NF90_PUT_ATT(nid,id_BILS,"title","Reference flux de chaleur au sol"),fnam)
     219  CALL ncerr(NF90_PUT_ATT(nid,id_ALB, "title","Albedo a la surface"),fnam)
     220  CALL ncerr(NF90_PUT_ATT(nid,id_RUG, "title","Rugosite"),fnam)
    246221
    247222  call nf95_put_att(nid, varid_longitude, "standard_name", "longitude")
    … …  
    251226  call nf95_put_att(nid, varid_latitude, "units", "degrees_north")
    252227
    253   ierr=NF90_ENDDEF(nid)
     228  CALL ncerr(NF90_ENDDEF(nid),fnam)
    254229
    255230  !--- Variables saving
    256   ierr=NF90_PUT_VAR(nid,id_tim,(/(REAL(k),k=1,ndays)/))
    257   ierr=NF90_PUT_VAR(nid,id_FOCE,pctsrf_t(:,is_oce,:),(/1,1/),(/klon,ndays/))
    258   ierr=NF90_PUT_VAR(nid,id_FSIC,pctsrf_t(:,is_sic,:),(/1,1/),(/klon,ndays/))
    259   ierr=NF90_PUT_VAR(nid,id_FTER,pctsrf_t(:,is_ter,:),(/1,1/),(/klon,ndays/))
    260   ierr=NF90_PUT_VAR(nid,id_FLIC,pctsrf_t(:,is_lic,:),(/1,1/),(/klon,ndays/))
    261   ierr=NF90_PUT_VAR(nid,id_SST ,phy_sst(:,:),(/1,1/),(/klon,ndays/))
    262   ierr=NF90_PUT_VAR(nid,id_BILS,phy_bil(:,:),(/1,1/),(/klon,ndays/))
    263   ierr=NF90_PUT_VAR(nid,id_ALB ,phy_alb(:,:),(/1,1/),(/klon,ndays/))
    264   ierr=NF90_PUT_VAR(nid,id_RUG ,phy_rug(:,:),(/1,1/),(/klon,ndays/))
     231  CALL ncerr(NF90_PUT_VAR(nid,id_tim,[(REAL(k),k=1,ndays)]),fnam)
     232  CALL ncerr(NF90_PUT_VAR(nid,id_FOCE,pctsrf_t(:,is_oce,:),[1,1],[klon,ndays]),fnam)
     233  CALL ncerr(NF90_PUT_VAR(nid,id_FSIC,pctsrf_t(:,is_sic,:),[1,1],[klon,ndays]),fnam)
     234  CALL ncerr(NF90_PUT_VAR(nid,id_FTER,pctsrf_t(:,is_ter,:),[1,1],[klon,ndays]),fnam)
     235  CALL ncerr(NF90_PUT_VAR(nid,id_FLIC,pctsrf_t(:,is_lic,:),[1,1],[klon,ndays]),fnam)
     236  CALL ncerr(NF90_PUT_VAR(nid,id_SST ,phy_sst(:,:),[1,1],[klon,ndays]),fnam)
     237  CALL ncerr(NF90_PUT_VAR(nid,id_BILS,phy_bil(:,:),[1,1],[klon,ndays]),fnam)
     238  CALL ncerr(NF90_PUT_VAR(nid,id_ALB ,phy_alb(:,:),[1,1],[klon,ndays]),fnam)
     239  CALL ncerr(NF90_PUT_VAR(nid,id_RUG ,phy_rug(:,:),[1,1],[klon,ndays]),fnam)
    265240  call nf95_put_var(nid, varid_longitude, rlon)
    266241  call nf95_put_var(nid, varid_latitude, rlat)
    267242
    268   ierr=NF90_CLOSE(nid)
    269 
    270   IF (prt_level>5) WRITE(lunout,*) 'Ecriture du fichier limit : fin'
     243  CALL ncerr(NF90_CLOSE(nid),fnam)
     244
     245  CALL msg(6,'Ecriture du fichier limit : fin')
    271246
    272247  DEALLOCATE(pctsrf_t,phy_sst,phy_bil,phy_alb,phy_rug)
    … …  
    296271  USE dimphy, ONLY : klon
    297272  USE phys_state_var_mod, ONLY : pctsrf
     273  USE conf_dat_m, ONLY: conf_dat2d
    298274  USE control_mod
    299275  USE pchsp_95_m, only: pchsp_95
    … …  
    303279
    304280  IMPLICIT NONE
    305 #include "dimensions.h"
    306 #include "paramet.h"
    307 #include "comgeom2.h"
    308 #include "iniprint.h"
     281  include "dimensions.h"
     282  include "paramet.h"
     283  include "comgeom2.h"
     284  include "iniprint.h"
    309285!-----------------------------------------------------------------------------
    310286! Arguments:
    … …  
    320296! Local variables:
    321297!--- NetCDF
    322   INTEGER :: ncid, varid                  ! NetCDF identifiers
    323   CHARACTER(LEN=30)               :: dnam       ! dimension name
     298  INTEGER           :: ncid, varid        ! NetCDF identifiers
     299  CHARACTER(LEN=30) :: dnam               ! dimension name
    324300!--- dimensions
    325   INTEGER,           DIMENSION(4) :: dids       ! NetCDF dimensions identifiers
    326   REAL, ALLOCATABLE, DIMENSION(:) :: dlon_ini   ! initial longitudes vector
    327   REAL, ALLOCATABLE, DIMENSION(:) :: dlat_ini   ! initial latitudes  vector
    328   REAL, POINTER,     DIMENSION(:) :: dlon, dlat ! reordered lon/lat  vectors
     301  INTEGER           :: dids(4)            ! NetCDF dimensions identifiers
     302  REAL, ALLOCATABLE :: dlon_ini(:)        ! initial longitudes vector
     303  REAL, ALLOCATABLE :: dlat_ini(:)        ! initial latitudes  vector
     304  REAL, POINTER     :: dlon(:), dlat(:)  ! reordered lon/lat  vectors
    329305!--- fields
    330   INTEGER :: imdep, jmdep, lmdep                ! dimensions of 'champ'
    331   REAL, ALLOCATABLE, DIMENSION(:, :) :: champ   ! wanted field on initial grid
    332   REAL, ALLOCATABLE, DIMENSION(:) :: yder, timeyear
    333   REAL,              DIMENSION(iim, jjp1) :: champint  ! interpolated field
    334   REAL, ALLOCATABLE, DIMENSION(:, :, :)    :: champtime
    335   REAL, ALLOCATABLE, DIMENSION(:, :, :)    :: champan
     306  INTEGER :: imdep, jmdep, lmdep          ! dimensions of 'champ'
     307  REAL, ALLOCATABLE :: champ(:,:)         ! wanted field on initial grid
     308  REAL, ALLOCATABLE :: yder(:), timeyear(:)
     309  REAL              :: champint(iim,jjp1) ! interpolated field
     310  REAL, ALLOCATABLE :: champtime(:,:,:)
     311  REAL, ALLOCATABLE :: champan(:,:,:)
    336312!--- input files
    337   CHARACTER(LEN=20)                 :: cal_in   ! calendar
    338   CHARACTER(LEN=20)                 :: unit_sic ! attribute unit in sea-ice file
    339   INTEGER                           :: ndays_in ! number of days
     313  CHARACTER(LEN=20) :: cal_in             ! calendar
     314  CHARACTER(LEN=20) :: unit_sic          ! attribute unit in sea-ice file
     315  INTEGER           :: ndays_in          ! number of days
    340316!--- misc
    341   INTEGER :: i, j, k, l                         ! loop counters
    342   REAL, ALLOCATABLE, DIMENSION(:, :) :: work     ! used for extrapolation
    343   CHARACTER(LEN=25)                 :: title    ! for messages
    344   LOGICAL                           :: extrp    ! flag for extrapolation
    345   LOGICAL                           :: oldice   ! flag for old way ice computation
    346   REAL                              :: chmin, chmax
     317  INTEGER           :: i, j, k, l         ! loop counters
     318  REAL, ALLOCATABLE :: work(:,:)          ! used for extrapolation
     319  CHARACTER(LEN=25) :: title              ! for messages
     320  LOGICAL           :: extrp              ! flag for extrapolation
     321  LOGICAL           :: oldice             ! flag for old way ice computation
     322  REAL              :: chmin, chmax
    347323  INTEGER ierr
    348324  integer n_extrap ! number of extrapolated points
    … …  
    369345
    370346!--- GETTING SOME DIMENSIONAL VARIABLES FROM FILE -----------------------------
    371   IF (prt_level>5) WRITE(lunout,*) ' Now reading file : ',fnam
    372   ierr=NF90_OPEN(fnam, NF90_NOWRITE, ncid);             CALL ncerr(ierr, fnam)
    373   ierr=NF90_INQ_VARID(ncid, trim(varname), varid);            CALL ncerr(ierr, fnam)
    374   ierr=NF90_INQUIRE_VARIABLE(ncid, varid, dimids=dids); CALL ncerr(ierr, fnam)
     347  CALL msg(5,' Now reading file : '//TRIM(fnam))
     348  CALL ncerr(NF90_OPEN(fnam, NF90_NOWRITE, ncid),fnam)
     349  CALL ncerr(NF90_INQ_VARID(ncid, trim(varname), varid),fnam)
     350  CALL ncerr(NF90_INQUIRE_VARIABLE(ncid, varid, dimids=dids),fnam)
    375351
    376352!--- Read unit for sea-ice variable only
    377353  IF (mode=='SIC') THEN
    378      ierr=NF90_GET_ATT(ncid, varid, 'units', unit_sic)
    379      IF(ierr/=NF90_NOERR) THEN
    380         IF (prt_level>5) WRITE(lunout,*) 'No unit was given in sea-ice file. Take percentage as default value'
    381         unit_sic='X'
    382      ELSE
    383         IF (prt_level>5) WRITE(lunout,*) ' Sea-ice cover has unit=',unit_sic
    384      END IF
     354    IF(NF90_GET_ATT(ncid, varid, 'units', unit_sic)/=NF90_NOERR) THEN
     355      CALL msg(5,'No unit in sea-ice file. Take percentage as default value')
     356      unit_sic='X'
     357    ELSE
     358      CALL msg(5,'Sea-ice cover has unit='//TRIM(unit_sic))
     359    END IF
    385360  END IF
    386361
    387362!--- Longitude
    388   ierr=NF90_INQUIRE_DIMENSION(ncid, dids(1), name=dnam, len=imdep)
    389   CALL ncerr(ierr, fnam); ALLOCATE(dlon_ini(imdep), dlon(imdep))
    390   ierr=NF90_INQ_VARID(ncid, dnam, varid);                CALL ncerr(ierr, fnam)
    391   ierr=NF90_GET_VAR(ncid, varid, dlon_ini);              CALL ncerr(ierr, fnam)
    392   IF (prt_level>5) WRITE(lunout, *) 'variable ', dnam, 'dimension ', imdep
     363  CALL ncerr(NF90_INQUIRE_DIMENSION(ncid, dids(1), name=dnam, len=imdep),fnam)
     364  ALLOCATE(dlon_ini(imdep), dlon(imdep))
     365  CALL ncerr(NF90_INQ_VARID(ncid, dnam, varid), fnam)
     366  CALL ncerr(NF90_GET_VAR(ncid, varid, dlon_ini), fnam)
     367  CALL msg(5,'variable '//TRIM(dnam)//' dimension ', imdep)
    393368
    394369!--- Latitude
    395   ierr=NF90_INQUIRE_DIMENSION(ncid, dids(2), name=dnam, len=jmdep)
    396   CALL ncerr(ierr, fnam); ALLOCATE(dlat_ini(jmdep), dlat(jmdep))
    397   ierr=NF90_INQ_VARID(ncid, dnam, varid);                CALL ncerr(ierr, fnam)
    398   ierr=NF90_GET_VAR(ncid, varid, dlat_ini);              CALL ncerr(ierr, fnam)
    399   IF (prt_level>5) WRITE(lunout, *) 'variable ', dnam, 'dimension ', jmdep
     370  CALL ncerr(NF90_INQUIRE_DIMENSION(ncid, dids(2), name=dnam, len=jmdep),fnam)
     371  ALLOCATE(dlat_ini(jmdep), dlat(jmdep))
     372  CALL ncerr(NF90_INQ_VARID(ncid, dnam, varid), fnam)
     373  CALL ncerr(NF90_GET_VAR(ncid, varid, dlat_ini), fnam)
     374  CALL msg(5,'variable '//TRIM(dnam)//' dimension ', jmdep)
    400375
    401376!--- Time (variable is not needed - it is rebuilt - but calendar is)
    402   ierr=NF90_INQUIRE_DIMENSION(ncid, dids(3), name=dnam, len=lmdep)
    403   CALL ncerr(ierr, fnam); ALLOCATE(timeyear(lmdep))
    404   ierr=NF90_INQ_VARID(ncid, dnam, varid);                CALL ncerr(ierr, fnam)
     377  CALL ncerr(NF90_INQUIRE_DIMENSION(ncid, dids(3), name=dnam, len=lmdep), fnam)
     378  ALLOCATE(timeyear(lmdep))
     379  CALL ncerr(NF90_INQ_VARID(ncid, dnam, varid), fnam)
    405380  cal_in=' '
    406   ierr=NF90_GET_ATT(ncid, varid, 'calendar', cal_in)
    407   IF(ierr/=NF90_NOERR) THEN
     381  IF(NF90_GET_ATT(ncid, varid, 'calendar', cal_in)/=NF90_NOERR) THEN
    408382    SELECT CASE(mode)
    409383      CASE('RUG', 'ALB'); cal_in='360d'
    410384      CASE('SIC', 'SST'); cal_in='gregorian'
    411385    END SELECT
    412     IF (prt_level>5) WRITE(lunout, *)'ATTENTION: variable "time" sans attribut "calendrier" ' &
    413          // 'dans '//TRIM(fnam)//'. On choisit la valeur par defaut.'
    414   END IF
    415   IF (prt_level>5) WRITE(lunout, *) 'variable ', dnam, 'dimension ', lmdep, 'calendrier ', &
    416        cal_in
    417 
     386  CALL msg(5,'WARNING: missing "calendar" attribute for "time" in '&
     387     &//TRIM(fnam)//'. Choosing default value.')
     388  END IF
     389  CALL msg(5,'var, calendar, dim: '//TRIM(dnam)//' '//TRIM(cal_in), lmdep)
    418390 
    419391!--- CONSTRUCTING THE INPUT TIME VECTOR FOR INTERPOLATION --------------------
    … …  
    430402  ALLOCATE(champ(imdep, jmdep), champtime(iim, jjp1, lmdep))
    431403  IF(extrp) ALLOCATE(work(imdep, jmdep))
    432 
    433   IF (prt_level>5) WRITE(lunout, *)
    434   IF (prt_level>5) WRITE(lunout,*)'LECTURE ET INTERPOLATION HORIZ. DE ', lmdep, ' CHAMPS.'
    435   ierr=NF90_INQ_VARID(ncid, varname, varid);             CALL ncerr(ierr, fnam)
     404  CALL msg(5,'')
     405  CALL msg(5,'LECTURE ET INTERPOLATION HORIZ. DE ', lmdep, ' CHAMPS.')
     406  CALL ncerr(NF90_INQ_VARID(ncid, varname, varid), fnam)
    436407  DO l=1, lmdep
    437     ierr=NF90_GET_VAR(ncid, varid, champ, (/1, 1, l/), (/imdep, jmdep, 1/))
    438     CALL ncerr(ierr, fnam)
    439     CALL conf_dat2d(title, imdep, jmdep, dlon_ini, dlat_ini, dlon, dlat, &
    440          champ, ibar)
    441 
    442     IF (extrp) CALL extrapol(champ, imdep, jmdep, 999999., .TRUE., .TRUE., 2, &
    443          work)
     408    CALL ncerr(NF90_GET_VAR(ncid,varid,champ,[1,1,l],[imdep,jmdep,1]),fnam)
     409    CALL conf_dat2d(title, dlon_ini, dlat_ini, dlon, dlat, champ, ibar)
     410    IF(extrp) CALL extrapol(champ,imdep,jmdep,999999.,.TRUE.,.TRUE.,2,work)
    444411
    445412    IF(ibar .AND. .NOT.oldice) THEN
    446       IF(l==1 .AND. prt_level>5) THEN
    447         WRITE(lunout, *) '-------------------------------------------------------------------------'
    448         WRITE(lunout, *) 'Utilisation de l''interpolation barycentrique pour ',TRIM(title),' $$$'
    449         WRITE(lunout, *) '-------------------------------------------------------------------------'
     413      IF(l==1) THEN
     414      CALL msg(5,"----------------------------------------------------------")
     415      CALL msg(5,"$$$ Interpolation barycentrique pour "//TRIM(title)//" $$$")
     416      CALL msg(5,"----------------------------------------------------------")
    450417      END IF
    451418      IF(mode=='RUG') champ=LOG(champ)
    452       CALL inter_barxy(dlon, dlat(:jmdep-1), champ, rlonu(:iim),     &
    453                          rlatv, champint)
     419      CALL inter_barxy(dlon,dlat(:jmdep-1),champ,rlonu(:iim),rlatv,champint)
    454420      IF(mode=='RUG') THEN
    455421        champint=EXP(champint)
    … …  
    458424    ELSE
    459425      SELECT CASE(mode)
    460         CASE('RUG');       CALL rugosite(imdep, jmdep, dlon, dlat, champ,    &
    461                                     iim, jjp1, rlonv, rlatu, champint, mask)
    462         CASE('SIC');       CALL sea_ice (imdep, jmdep, dlon, dlat, champ,    &
    463                                     iim, jjp1, rlonv, rlatu, champint)
    464         CASE('SST', 'ALB'); CALL grille_m(imdep, jmdep, dlon, dlat, champ,    &
    465                                     iim, jjp1, rlonv, rlatu, champint)
     426        CASE('RUG');       CALL rugosite(dlon,dlat,champ,rlonv,rlatu,champint,mask)
     427        CASE('SIC');       CALL sea_ice (dlon,dlat,champ,rlonv,rlatu,champint)
     428        CASE('SST','ALB'); CALL grille_m(dlon,dlat,champ,rlonv,rlatu,champint)
    466429      END SELECT
    467430    END IF
    468431    champtime(:, :, l)=champint
    469432  END DO
    470   ierr=NF90_CLOSE(ncid);                                 CALL ncerr(ierr, fnam)
     433  CALL ncerr(NF90_CLOSE(ncid), fnam)
    471434
    472435  DEALLOCATE(dlon_ini, dlat_ini, dlon, dlat, champ)
    … …  
    474437
    475438!--- TIME INTERPOLATION ------------------------------------------------------
    476   IF (prt_level>5) THEN
     439  IF(prt_level>5) THEN
    477440     WRITE(lunout, *)
    478441     WRITE(lunout, *)'INTERPOLATION TEMPORELLE.'
    … …  
    508471!--- SPECIAL FILTER FOR SST: SST>271.38 --------------------------------------
    509472  IF(mode=='SST') THEN
    510     IF (prt_level>5) WRITE(lunout, *) 'Filtrage de la SST: SST >= 271.38'
     473    CALL msg(5,'Filtrage de la SST: SST >= 271.38')
    511474    WHERE(champan<271.38) champan=271.38
    512475  END IF
    … …  
    514477!--- SPECIAL FILTER FOR SIC: 0.0<SIC<1.0 -------------------------------------
    515478  IF(mode=='SIC') THEN
    516     IF (prt_level>5) WRITE(lunout, *) 'Filtrage de la SIC: 0.0 < Sea-ice < 1.0'
     479    CALL msg(5,'Filtrage de la SIC: 0.0 < Sea-ice < 1.0')
    517480
    518481    IF (unit_sic=='1') THEN
    519482       ! Nothing to be done for sea-ice field is already in fraction of 1
    520483       ! This is the case for sea-ice in file cpl_atm_sic.nc
    521        IF (prt_level>5) WRITE(lunout,*) 'Sea-ice field already in fraction of 1'
     484       CALL msg(5,'Sea-ice field already in fraction of 1')
    522485    ELSE
    523486       ! Convert sea ice from percentage to fraction of 1
    524        IF (prt_level>5) WRITE(lunout,*) 'Transformt sea-ice field from percentage to fraction of 1.'
     487       CALL msg(5,'Transformt sea-ice field from percentage to fraction of 1.')
    525488       champan(:, :, :)=champan(:, :, :)/100.
    526489    END IF
    … …  
    541504!
    542505!-------------------------------------------------------------------------------
    543 
    544506
    545507
    … …  
    620582
    621583  ELSE
    622     mnth=(/(m,m=1,nm,nd/nm)/)
     584    mnth=[(m,m=1,nm,nd/nm)]
    623585  END IF
    624586
    … …  
    634596
    635597
     598
     599!-------------------------------------------------------------------------------
     600!
     601SUBROUTINE msg(lev,str1,i,str2)
     602!
     603!-------------------------------------------------------------------------------
     604! Arguments:
     605  INTEGER,                    INTENT(IN) :: lev
     606  CHARACTER(LEN=*),           INTENT(IN) :: str1
     607  INTEGER,          OPTIONAL, INTENT(IN) :: i
     608  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: str2
     609!-------------------------------------------------------------------------------
     610  IF(prt_level>lev) THEN
     611    IF(PRESENT(str2)) THEN
     612      WRITE(lunout,*) TRIM(str1), i, TRIM(str2)
     613    ELSE IF(PRESENT(i)) THEN
     614      WRITE(lunout,*) TRIM(str1), i
     615    ELSE
     616      WRITE(lunout,*) TRIM(str1)
     617    END IF
     618  END IF
     619
     620END SUBROUTINE msg
     621!
     622!-------------------------------------------------------------------------------
     623
     624
    636625!-------------------------------------------------------------------------------
    637626!

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

    r2265 r2293  
    1 ! $Id$
    2 
    31SUBROUTINE phyredem (fichnom)
    4 
     2!
     3!-------------------------------------------------------------------------------
     4! Author: Z.X. Li (LMD/CNRS), 1993/08/18
     5!-------------------------------------------------------------------------------
     6! Purpose: Write restart state for physics.
     7!-------------------------------------------------------------------------------
    58  USE dimphy
    69  USE mod_grid_phy_lmdz
    … …  
    1922
    2023  IMPLICIT none
    21   !======================================================================
    22   ! Auteur(s) Z.X. Li (LMD/CNRS) date: 19930818
    23   ! Objet: Ecriture de l'etat de redemarrage pour la physique
    24   !======================================================================
    25   include "netcdf.inc"
     24
     25  include "iniprint.h"
    2626  include "dimsoil.h"
    2727  include "clesphys.h"
    … …  
    4949  INTEGER isoil, nsrf,isw
    5050  CHARACTER (len=7) :: str7
    51   CHARACTER (len=2) :: str2
     51  CHARACTER (len=256) :: nam, lnam
    5252  INTEGER           :: it, iiq
    5353
    … …  
    126126  CALL put_field("FSIC", "fraction glace mer", pctsrf(:, is_sic))
    127127
    128   DO nsrf = 1, nbsrf
    129      IF (nsrf.LE.99) THEN
    130         WRITE(str2, '(i2.2)') nsrf
    131         CALL put_field("TS"//str2, "Temperature de surface No."//str2, &
    132              ftsol(:, nsrf))
    133      ELSE
    134         PRINT*, "Trop de sous-mailles"
    135         call abort_gcm("phyredem", "", 1)
    136      ENDIF
    137   ENDDO
     128  IF(nbsrf>99) THEN
     129    PRINT*, "Trop de sous-mailles";  CALL abort_gcm("phyredem", "", 1)
     130  END IF
     131  IF(nsoilmx>99) THEN
     132    PRINT*, "Trop de sous-surfaces"; CALL abort_gcm("phyredem", "", 1)
     133  END IF
     134
     135  CALL put_field_srf1("TS","Temperature",ftsol(:,:))
    138136
    139137! ================== Albedo =======================================
    140138  print*,'PHYREDEM NOUVEAU'
    141   DO nsrf = 1, nbsrf
    142      DO isw=1, nsw
    143         IF (isw.LE.99 .AND. nsrf.LE.99) THEN
    144            WRITE(str7, '(i2.2, "srf", i2.2)') isw, nsrf
    145   print*,'PHYREDEM ',"A_dir_SW"//str7
    146            CALL put_field("A_dir_SW"//str7, "Albedo direct du sol bande "//str7, &
    147                 falb_dir(:, isw, nsrf))
    148            CALL put_field("A_dif_SW"//str7, "Albedo difus du sol bande "//str7, &
    149                 falb_dif(:, isw, nsrf))
    150         ELSE
    151            PRINT*, "Trop de couches"
    152            call abort_gcm("phyredem", "", 1)
    153         ENDIF
    154      ENDDO
    155   ENDDO
    156 
    157 ! ================== Tsoil =======================================
    158   DO nsrf = 1, nbsrf
    159      DO isoil=1, nsoilmx
    160         IF (isoil.LE.99 .AND. nsrf.LE.99) THEN
    161            WRITE(str7, '(i2.2, "srf", i2.2)') isoil, nsrf
    162            CALL put_field("Tsoil"//str7, "Temperature du sol No."//str7, &
    163                 tsoil(:, isoil, nsrf))
    164         ELSE
    165            PRINT*, "Trop de couches"
    166            call abort_gcm("phyredem", "", 1)
    167         ENDIF
    168      ENDDO
    169   ENDDO
    170 
    171   DO nsrf = 1, nbsrf
    172      IF (nsrf.LE.99) THEN
    173         WRITE(str2, '(i2.2)') nsrf
    174         CALL put_field("QS"//str2, "Humidite de surface No."//str2, &
    175              qsurf(:, nsrf))
    176      ELSE
    177         PRINT*, "Trop de sous-mailles"
    178         call abort_gcm("phyredem", "", 1)
    179      ENDIF
    180   END DO
    181 
    182   CALL put_field("QSOL", "Eau dans le sol (mm)", qsol)
    183 
    184   DO nsrf = 1, nbsrf
    185      IF (nsrf.LE.99) THEN
    186         WRITE(str2, '(i2.2)') nsrf
    187         CALL put_field("EVAP"//str2, "Evaporation de surface No."//str2 &
    188              , fevap(:, nsrf))
    189      ELSE
    190         PRINT*, "Trop de sous-mailles"
    191         call abort_gcm("phyredem", "", 1)
    192      ENDIF
    193   ENDDO
    194 
    195   DO nsrf = 1, nbsrf
    196      IF (nsrf.LE.99) THEN
    197         WRITE(str2, '(i2.2)') nsrf
    198         CALL put_field("SNOW"//str2, "Neige de surface No."//str2, &
    199              snow(:, nsrf))
    200      ELSE
    201         PRINT*, "Trop de sous-mailles"
    202         call abort_gcm("phyredem", "", 1)
    203      ENDIF
    204   ENDDO
     139  CALL put_field_srf2("A_dir_SW","Albedo direct",falb_dir(:,:,:))
     140  CALL put_field_srf2("A_dif_SW","Albedo diffus",falb_dif(:,:,:))
     141
     142! ================== Tsoil =========================================
     143  CALL put_field_srf2("Tsoil","Temperature",tsoil(:,:,:))
     144
     145  CALL put_field_srf1("QS"  , "Humidite",qsurf(:,:))
     146
     147  CALL put_field     ("QSOL", "Eau dans le sol (mm)", qsol)
     148
     149  CALL put_field_srf1("EVAP", "Evaporation", fevap(:,:))
     150
     151  CALL put_field_srf1("SNOW", "Neige", fevap(:,:))
    205152
    206153  CALL put_field("RADS", "Rayonnement net a la surface", radsol)
    … …  
    218165  CALL put_field("snow_f", "precipitation solide", snow_fall)
    219166
    220   DO nsrf = 1, nbsrf
    221      IF (nsrf.LE.99) THEN
    222         WRITE(str2, '(i2.2)') nsrf
    223         CALL put_field("Z0m"//str2, "rugosite de surface No."//str2, &
    224              z0m(:, nsrf))
    225         CALL put_field("Z0h"//str2, "rugosite de surface No."//str2, &
    226              z0h(:, nsrf))
    227      ELSE
    228         PRINT*, "Trop de sous-mailles"
    229         call abort_gcm("phyredem", "", 1)
    230      ENDIF
    231   ENDDO
    232 
    233   DO nsrf = 1, nbsrf
    234      IF (nsrf.LE.99) THEN
    235         WRITE(str2, '(i2.2)') nsrf
    236         CALL put_field("AGESNO"//str2, &
    237              "Age de la neige surface No."//str2, &
    238              agesno(:, nsrf))
    239      ELSE
    240         PRINT*, "Trop de sous-mailles"
    241         call abort_gcm("phyredem", "", 1)
    242      ENDIF
    243   ENDDO
     167  CALL put_field_srf1("Z0m", "rugosite", z0m(:,:))
     168
     169  CALL put_field_srf1("Z0h", "rugosite", z0h(:,:))
     170
     171  CALL put_field_srf1("AGESNO", "Age de la neige", agesno(:,:))
    244172
    245173  CALL put_field("ZMEA", "ZMEA", zmea)
    … …  
    280208
    281209  IF (iflag_pbl>1) then
    282      DO nsrf = 1, nbsrf
    283         IF (nsrf.LE.99) THEN
    284            WRITE(str2, '(i2.2)') nsrf
    285            CALL put_field("TKE"//str2, "Energ. Cineti. Turb."//str2, &
    286                 pbl_tke(:, 1:klev+1, nsrf))
    287            CALL put_field("DELTATKE"//str2, "Del TKE wk/env."//str2, &
    288                 wake_delta_pbl_tke(:, 1:klev+1, nsrf))
    289         ELSE
    290            PRINT*, "Trop de sous-mailles"
    291            call abort_gcm("phyredem", "", 1)
    292         ENDIF
    293      ENDDO
    294   ENDIF
     210    CALL put_field_srf3("TKE", "Energ. Cineti. Turb.", &
     211         pbl_tke(:,:,:))
     212    CALL put_field_srf3("DELTATKE", "Del TKE wk/env.", &
     213         wake_delta_pbl_tke(:,:,:))
     214  END IF
    295215
    296216  ! FIN TKE PBL !
    … …  
    372292  !$OMP BARRIER
    373293
     294
     295  CONTAINS
     296
     297
     298SUBROUTINE put_field_srf1(nam,lnam,field)
     299
     300  IMPLICIT NONE
     301  CHARACTER(LEN=*),  INTENT(IN) :: nam, lnam
     302  REAL,              INTENT(IN) :: field(:,:)
     303  CHARACTER(LEN=256) :: nm, lm, str
     304  DO nsrf = 1, nbsrf
     305    WRITE(str, '(i2.2)') nsrf
     306    nm=TRIM(nam)//TRIM(str)
     307    lm=TRIM(lnam)//" de surface No. "//TRIM(str)
     308    CALL put_field(nm,lm,field(:,nsrf))
     309  END DO
     310
     311END SUBROUTINE put_field_srf1
     312
     313
     314SUBROUTINE put_field_srf2(nam,lnam,field)
     315
     316  IMPLICIT NONE
     317  CHARACTER(LEN=*),  INTENT(IN) :: nam, lnam
     318  REAL,              INTENT(IN) :: field(:,:,:)
     319  CHARACTER(LEN=256) :: nm, lm, str
     320  DO nsrf = 1, nbsrf
     321    DO isoil=1, nsw
     322      WRITE(str, '(i2.2,"srf",i2.2)')isoil,nsrf
     323!      WRITE(lunout,*)"PHYREDEM ",TRIM(nam)//TRIM(str)
     324      nm=TRIM(nam)//TRIM(str)
     325      lm=TRIM(lnam)//" du sol No. "//TRIM(str)
     326      CALL put_field(nm,lm,field(:,isoil,nsrf))
     327    END DO
     328  END DO
     329
     330END SUBROUTINE put_field_srf2
     331
     332
     333SUBROUTINE put_field_srf3(nam,lnam,field)
     334
     335  IMPLICIT NONE
     336  CHARACTER(LEN=*),  INTENT(IN) :: nam, lnam
     337  REAL,              INTENT(IN) :: field(:,:,:)
     338  CHARACTER(LEN=256) :: nm, lm, str
     339  DO nsrf = 1, nbsrf
     340    WRITE(str, '(i2.2)') nsrf
     341    nm=TRIM(nam)//TRIM(str)
     342    lm=TRIM(lnam)//TRIM(str)
     343    CALL put_field(nm,lm,field(:,1:klev+1,nsrf))
     344  END DO
     345
     346END SUBROUTINE put_field_srf3
     347
     348
    374349END SUBROUTINE phyredem
Note: See TracChangeset for help on using the changeset viewer.