Changeset 2947 for LMDZ5/trunk


Ignore:
Timestamp:
Jul 13, 2017, 6:43:58 PM (7 years ago)
Author:
dcugnet
Message:
  • The "units" attribute of SST variable is now checked to allow celcius degrees to kelvins conversion if required.
  • Few cosmetic changes for output text file to make it more readable (and help to detect if ce0l has really done what was expected).
File:
1 edited

Legend:

Unmodified
Added
Removed
  • LMDZ5/trunk/libf/dynphy_lonlat/phylmd/limit_netcdf.F90

    r2942 r2947  
    128128
    129129!--- RUGOSITY TREATMENT --------------------------------------------------------
    130   CALL msg(1,'Traitement de la rugosite')
     130  CALL msg(0,""); CALL msg(0," *** TRAITEMENT DE LA RUGOSITE ***")
    131131  CALL get_2Dfield(frugo,vrug,'RUG',ndays,phy_rug,mask=masque(1:iim,:))
    132132
    133133!--- OCEAN TREATMENT -----------------------------------------------------------
    134   CALL msg(1,'Traitement de la glace oceanique')
     134  CALL msg(0,""); CALL msg(0," *** TRAITEMENT DE LA GLACE OCEANIQUE ***")
    135135
    136136! Input SIC file selection
     
    148148  END IF
    149149  CALL ncerr(NF90_CLOSE(nid),icefile)
    150   CALL msg(-1,'Fichier choisi pour la glace de mer:'//TRIM(icefile))
     150  CALL msg(0,'Fichier choisi pour la glace de mer:'//TRIM(icefile))
    151151
    152152  CALL get_2Dfield(icefile,varname, 'SIC',ndays,phy_ice)
     
    191191
    192192!--- SST TREATMENT -------------------------------------------------------------
    193   CALL msg(1,'Traitement de la sst')
     193  CALL msg(0,""); CALL msg(0," *** TRAITEMENT DE LA SST ***")
    194194
    195195! Input SST file selection
     
    207207  END IF
    208208  CALL ncerr(NF90_CLOSE(nid),sstfile)
    209   CALL msg(-1,'Fichier choisi pour la temperature de mer: '//TRIM(sstfile))
     209  CALL msg(0,'Fichier choisi pour la temperature de mer: '//TRIM(sstfile))
    210210
    211211  CALL get_2Dfield(sstfile,varname,'SST',ndays,phy_sst,flag=extrap)
    212212
    213213!--- ALBEDO TREATMENT ----------------------------------------------------------
    214   CALL msg(1,'Traitement de l albedo')
     214  CALL msg(0,""); CALL msg(0," *** TRAITEMENT DE L'ALBEDO ***")
    215215  CALL get_2Dfield(falbe,valb,'ALB',ndays,phy_alb)
    216216
     
    219219
    220220!--- OUTPUT FILE WRITING -------------------------------------------------------
    221   CALL msg(5,'Ecriture du fichier limit : debut')
     221  CALL msg(0,""); CALL msg(0,' *** Ecriture du fichier limit : debut ***')
    222222  fnam="limit.nc"
    223223
     
    298298  CALL ncerr(NF90_CLOSE(nid),fnam)
    299299
    300   CALL msg(6,'Ecriture du fichier limit : fin')
     300  CALL msg(0,""); CALL msg(0,' *** Ecriture du fichier limit : fin ***')
    301301
    302302  DEALLOCATE(pctsrf_t,phy_sst,phy_bil,phy_alb,phy_rug)
     
    362362  CHARACTER(LEN=20) :: fnam_m, fnam_p     ! previous/next files names
    363363  CHARACTER(LEN=20) :: cal_in             ! calendar
    364   CHARACTER(LEN=20) :: unit_sic           ! attribute unit in sea-ice file
     364  CHARACTER(LEN=20) :: unit_sic           ! attribute "units" in sea-ice file
     365  CHARACTER(LEN=20) :: unit_sst           ! attribute "units" in sst     file
    365366  INTEGER           :: ndays_in           ! number of days
    366367!--- misc
     
    395396  CALL ncerr(NF90_INQUIRE_VARIABLE(ncid, varid, dimids=dids),fnam)
    396397
    397 !--- Read unit for sea-ice variable only
     398!--- Read unit for sea-ice and sst only
    398399  IF (mode=='SIC') THEN
    399400    ierr=NF90_GET_ATT(ncid, varid, 'units', unit_sic); CALL strclean(unit_sic)
     401    IF(ierr/=NF90_NOERR) THEN; unit_sic='%'
     402      CALL msg(0,'No unit in sea-ice file. Take percentage as default value')
     403    ELSE IF(TRIM(unit_sic)=="%") THEN
     404      CALL msg(0,'Sea-ice cover is a PERCENTAGE.')
     405    ELSE IF(ANY(unit_sic==["1.0","1  "])) THEN; unit_sic="1"
     406      CALL msg(0,'Sea-ice cover is a FRACTION.')
     407    ELSE
     408      CALL abort_physic('SIC','Unrecognized sea-ice unit: '//TRIM(unit_sic),1)
     409    END IF
     410  END IF
     411  IF (mode=='SST') THEN
     412    ierr=NF90_GET_ATT(ncid, varid, 'units', unit_sst); CALL strclean(unit_sst)
    400413    IF(ierr/=NF90_NOERR) THEN
    401       CALL msg(5,'No unit in sea-ice file. Take percentage as default value')
    402       unit_sic='X'
    403     ELSE IF(ANY(["%  ","1.0","1  "]==unit_sic)) THEN
    404       CALL msg(5,'Sea-ice cover has unit='//TRIM(unit_sic))
     414      CALL msg(0,'No unit in sst file. Take default: kelvins.')
     415      unit_sst='X'
     416    ELSE IF(ANY(unit_sst==["degC  ","DegC  "])) THEN; unit_sst='C'
     417      CALL msg(0,'Sea-surface temperature is in CELCIUS DEGREES.')
     418    ELSE IF(ANY(unit_sst==["K     ","Kelvin"])) THEN; unit_sst='K'
     419      CALL msg(0,'Sea-surface temperature is in KELVINS.')
    405420    ELSE
    406       CALL abort_physic('SIC','Unrecognized unit for sea-ice file: '&
    407         //TRIM(unit_sic),1)
     421      CALL abort_physic('SST','Unrecognized sst unit: '//TRIM(unit_sst),1)
    408422    END IF
    409423  END IF
     
    433447      CASE('SIC', 'SST'); cal_in='gregorian'
    434448    END SELECT
    435     CALL msg(5,'WARNING: missing "calendar" attribute for "time" in '&
     449    CALL msg(0,'WARNING: missing "calendar" attribute for "time" in '&
    436450     &//TRIM(fnam)//'. Choosing default value.')
    437451  END IF
    438452  CALL strclean(cal_in)                     !--- REMOVE (WEIRD) NULL CHARACTERS
    439   CALL msg(5,'var, calendar, dim: '//TRIM(dnam)//' '//TRIM(cal_in), lmdep)
     453  CALL msg(0,'var, calendar, dim: '//TRIM(dnam)//' '//TRIM(cal_in), lmdep)
    440454 
    441455!--- Determining input file number of days, depending on calendar
     
    469483        IF(ANY(champ>1.0+EPSFRA)) &
    470484          CALL abort_physic('SIC','Found sea-ice fractions greater than 1.')
    471       ELSE IF(TRIM(unit_sic)=="X".OR.TRIM(unit_sic)=="%") THEN
     485      ELSE IF(TRIM(unit_sic)=="%") THEN
    472486        IF(ANY(champ>100.0+EPSFRA)) &
    473487          CALL abort_physic('SIC','Found sea-ice percentages greater than 100.')
     
    542556  IF(prt_level>0) THEN
    543557     IF(ndays/=ndays_in) THEN
    544         WRITE(lunout, *)
    545558        WRITE(lunout,*)'DIFFERENT YEAR LENGTHS:'
    546         WRITE(lunout,*)' In the input file: ',ndays_in
     559        WRITE(lunout,*)' In the  input file: ',ndays_in
    547560        WRITE(lunout,*)' In the output file: ',ndays
    548561     END IF
    549562     IF(lmdep==ndays_in) THEN
    550         WRITE(lunout, *)
    551563        WRITE(lunout, *)'NO TIME INTERPOLATION.'
    552564        WRITE(lunout, *)' Daily input file.'
     
    602614!--- SPECIAL FILTER FOR SST: SST>271.38 --------------------------------------
    603615  IF(mode=='SST') THEN
     616    IF(TRIM(unit_sst)=="K") THEN
     617       ! Nothing to be done if the sst field is already in kelvins
     618      CALL msg(0,'SST field is already in kelvins.')
     619    ELSE
     620       ! Convert sst field from celcius degrees to kelvins
     621      CALL msg(0,'SST field converted from celcius degrees to kelvins.')
     622      champan=champan+273.15
     623    END IF
    604624    CALL msg(0,'Filtering SST: SST >= 271.38')
    605625    WHERE(champan<271.38) champan=271.38
     
    609629  IF(mode=='SIC') THEN
    610630    CALL msg(0,'Filtering SIC: 0.0 < Sea-ice < 1.0')
    611     IF(TRIM(unit_sic)=="1".OR.TRIM(unit_sic)=="1.0") THEN
    612        ! Nothing to be done for sea-ice field is already in fraction of 1
     631    IF(TRIM(unit_sic)=="1") THEN
     632       ! Nothing to be done if the sea-ice field is already in fraction of 1
    613633       ! This is the case for sea-ice in file cpl_atm_sic.nc
    614634       CALL msg(0,'Sea-ice field already in fraction of 1')
     
    709729  CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: str2
    710730!-------------------------------------------------------------------------------
    711   IF(prt_level>lev) THEN
     731  IF(prt_level>=lev) THEN
    712732    IF(PRESENT(str2)) THEN
    713733      WRITE(lunout,*) TRIM(str1), i, TRIM(str2)
Note: See TracChangeset for help on using the changeset viewer.