Index: LMDZ5/trunk/libf/dynphy_lonlat/phylmd/limit_netcdf.F90 =================================================================== --- LMDZ5/trunk/libf/dynphy_lonlat/phylmd/limit_netcdf.F90 (revision 2738) +++ LMDZ5/trunk/libf/dynphy_lonlat/phylmd/limit_netcdf.F90 (revision 2766) @@ -17,12 +17,12 @@ !------------------------------------------------------------------------------- - USE ioipsl, ONLY: flininfo, flinopen, flinget, flinclo, & - ioconf_calendar, ioget_calendar, lock_calendar, ioget_mon_len, ioget_year_len + USE ioipsl, ONLY: flininfo, flinopen, flinget, flinclo USE assert_eq_m, ONLY: assert_eq + USE cal_tools_m, ONLY: year_len, mid_month USE conf_dat_m, ONLY: conf_dat2d, conf_dat3d USE dimphy, ONLY: klon, zmasq - USE geometry_mod, ONLY: longitude_deg, latitude_deg + USE geometry_mod, ONLY: longitude_deg, latitude_deg USE phys_state_var_mod, ONLY: pctsrf - USE control_mod, ONLY: anneeref + USE control_mod, ONLY: anneeref USE init_ssrf_m, ONLY: start_init_subsurf @@ -123,5 +123,5 @@ !--- Beware: anneeref (from gcm.def) is used to determine output time sampling - ndays=ioget_year_len(anneeref) + ndays=year_len(anneeref) !--- RUGOSITY TREATMENT -------------------------------------------------------- @@ -375,9 +375,13 @@ !--- Read unit for sea-ice variable only IF (mode=='SIC') THEN - IF(NF90_GET_ATT(ncid, varid, 'units', unit_sic)/=NF90_NOERR) THEN + ierr=NF90_GET_ATT(ncid, varid, 'units', unit_sic); CALL strclean(unit_sic) + IF(ierr/=NF90_NOERR) THEN CALL msg(5,'No unit in sea-ice file. Take percentage as default value') unit_sic='X' + ELSE IF(ANY(["% ","1.0","1 "]==unit_sic)) THEN + CALL msg(5,'Sea-ice cover has unit='//TRIM(unit_sic)) ELSE - CALL msg(5,'Sea-ice cover has unit='//TRIM(unit_sic)) + CALL abort_physic('SIC','Unrecognized unit for sea-ice file: '& + //TRIM(unit_sic),1) END IF END IF @@ -418,5 +422,5 @@ !--- Rebuilding input time vector (field from input file might be unreliable) IF(lmdep==12) THEN - timeyear=mid_month(anneeref, cal_in, ndays_in) + timeyear=mid_month(anneeref, cal_in) CALL msg(0,'Monthly input file(s) for '//TRIM(title)//'.') ELSE IF(lmdep==ndays_in) THEN @@ -426,5 +430,5 @@ WRITE(mess,'(a,i3,a,i3,a)')'Mismatching input file: found',lmdep, & ' records, 12/',ndays_in,' (monthly/daily needed).' - CALL abort_physic('mid_months',TRIM(mess),1) + CALL abort_physic('mid_month',TRIM(mess),1) END IF @@ -433,13 +437,23 @@ IF(extrp) ALLOCATE(work(imdep, jmdep)) CALL msg(5,'') - CALL msg(5,'LECTURE ET INTERPOLATION HORIZ. DE ', lmdep, ' CHAMPS.') + CALL msg(5,'READ AND INTERPOLATE HORIZONTALLY ', lmdep, ' FIELDS.') CALL ncerr(NF90_INQ_VARID(ncid, varname, varid), fnam) DO l=1, lmdep CALL ncerr(NF90_GET_VAR(ncid,varid,champ,[1,1,l],[imdep,jmdep,1]),fnam) + !--- Check whether values are acceptable for SIC, depending on unit. + IF(mode=='SIC') THEN + IF(TRIM(unit_sic)=="1".OR.TRIM(unit_sic)=="1.0") THEN + IF(ANY(champ>1.0+EPSFRA)) & + CALL abort_physic('SIC','Found sea-ice fractions greater than 1.') + ELSE IF(TRIM(unit_sic)=="X".OR.TRIM(unit_sic)=="%") THEN + IF(ANY(champ>100.0+EPSFRA)) & + CALL abort_physic('SIC','Found sea-ice percentages greater than 100.') + END IF + END IF CALL conf_dat2d(title, dlon_ini, dlat_ini, dlon, dlat, champ, .TRUE.) IF(extrp) CALL extrapol(champ,imdep,jmdep,999999.,.TRUE.,.TRUE.,2,work) IF(l==1) THEN CALL msg(5,"----------------------------------------------------------") - CALL msg(5,"$$$ Interpolation barycentrique pour "//TRIM(title)//" $$$") + CALL msg(5,"$$$ Barycentrique interpolation for "//TRIM(title)//" $$$") CALL msg(5,"----------------------------------------------------------") END IF @@ -504,16 +518,16 @@ IF(ndays/=ndays_in) THEN WRITE(lunout, *) - WRITE(lunout,*)'DIFFERENTES LONGEURS D ANNEES:' - WRITE(lunout,*)' Dans le fichier d entree: ',ndays_in - WRITE(lunout,*)' Dans le fichier de sortie: ',ndays + WRITE(lunout,*)'DIFFERENT YEAR LENGTHS:' + WRITE(lunout,*)' In the input file: ',ndays_in + WRITE(lunout,*)' In the output file: ',ndays END IF IF(lmdep==ndays_in) THEN WRITE(lunout, *) - WRITE(lunout, *)'PAS D INTERPOLATION TEMPORELLE.' - WRITE(lunout, *)' Fichier journalier en entree.' + WRITE(lunout, *)'NO TIME INTERPOLATION.' + WRITE(lunout, *)' Daily input file.' ELSE - WRITE(lunout, *)'INTERPOLATION TEMPORELLE.' - WRITE(lunout, *)' Vecteur temps en entree: ', timeyear - WRITE(lunout, *)' Vecteur temps en sortie de 0 a ', ndays-1 + WRITE(lunout, *)'TIME INTERPOLATION.' + WRITE(lunout, *)' Input time vector: ', timeyear + WRITE(lunout, *)' Output time vector from 0 to ', ndays-1 END IF END IF @@ -542,10 +556,10 @@ DO j=1, jjp1 CALL minmax(iip1, champan(1, j, 10), chmin, chmax) - IF (prt_level>5) WRITE(lunout, *)' ',TRIM(title),' au temps 10 ', chmin, chmax, j + IF (prt_level>5) WRITE(lunout, *)' ',TRIM(title),' at time 10 ', chmin, chmax, j END DO !--- SPECIAL FILTER FOR SST: SST>271.38 -------------------------------------- IF(mode=='SST') THEN - CALL msg(5,'Filtrage de la SST: SST >= 271.38') + CALL msg(0,'Filtering SST: SST >= 271.38') WHERE(champan<271.38) champan=271.38 END IF @@ -553,16 +567,14 @@ !--- SPECIAL FILTER FOR SIC: 0.01.0) champan=1.0 @@ -647,80 +659,4 @@ !------------------------------------------------------------------------------- ! -FUNCTION year_len(y,cal_in) -! -!------------------------------------------------------------------------------- - IMPLICIT NONE -!------------------------------------------------------------------------------- -! Arguments: - INTEGER :: year_len - INTEGER, INTENT(IN) :: y - CHARACTER(LEN=*), INTENT(IN) :: cal_in -!------------------------------------------------------------------------------- -! Local variables: - CHARACTER(LEN=20) :: cal_out ! calendar (for outputs) -!------------------------------------------------------------------------------- -!--- Getting the input calendar to reset at the end of the function - CALL ioget_calendar(cal_out) - -!--- Unlocking calendar and setting it to wanted one - CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_in)) - -!--- Getting the number of days in this year - year_len=ioget_year_len(y) - -!--- Back to original calendar - CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_out)) - -END FUNCTION year_len -! -!------------------------------------------------------------------------------- - - -!------------------------------------------------------------------------------- -! -FUNCTION mid_month(y,cal_in,ndays_out) -! -!------------------------------------------------------------------------------- - IMPLICIT NONE -!------------------------------------------------------------------------------- -! Arguments: - INTEGER, INTENT(IN) :: y ! year - CHARACTER(LEN=*), INTENT(IN) :: cal_in ! calendar - INTEGER, INTENT(IN) :: ndays_out ! days number - REAL, DIMENSION(14) :: mid_month ! mid-bins times -!------------------------------------------------------------------------------- -! Local variables: - CHARACTER(LEN=99) :: mess ! error message - CHARACTER(LEN=20) :: cal_out ! output calendar - INTEGER, DIMENSION(14) :: tlen ! months lengths (days) - INTEGER :: m ! months counter - INTEGER :: nd ! number of days -!------------------------------------------------------------------------------- - !--- Save the input calendar to reset it at the end of the function - CALL ioget_calendar(cal_out) - - !--- Unlock calendar and set it to "cal_in" - CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_in)) - - !--- Get the length of each month - tlen(1 )=ioget_mon_len(y-1,12) - DO m=1,12; tlen(m+1)=ioget_mon_len(y,m); END DO - tlen(14)=ioget_mon_len(y+1, 1) - - !--- Mid-bins times - mid_month(1)=-0.5*REAL(tlen(1)) - DO m=2,14; mid_month(m)=mid_month(m-1)+0.5*REAL(tlen(m-1)+tlen(m)); END DO - - !--- Back to original calendar - CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_out)) - -END FUNCTION mid_month -! -!------------------------------------------------------------------------------- - - - -!------------------------------------------------------------------------------- -! SUBROUTINE msg(lev,str1,i,str2) !