Changeset 5100 for LMDZ6/branches/Amaury_dev/libf/dynphy_lonlat
- Timestamp:
- Jul 23, 2024, 7:00:20 AM (13 months ago)
- Location:
- LMDZ6/branches/Amaury_dev/libf/dynphy_lonlat/phylmd
- Files:
-
- 1 edited
- 1 moved
-
ce0l.F90 (modified) (5 diffs)
-
limit_netcdf.f90 (moved) (moved from LMDZ6/branches/Amaury_dev/libf/dynphy_lonlat/phylmd/limit_netcdf.F90) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/branches/Amaury_dev/libf/dynphy_lonlat/phylmd/ce0l.F90
r5099 r5100 21 21 USE etat0phys, ONLY: etat0phys_netcdf 22 22 USE limit, ONLY: limit_netcdf 23 USE netcdf, ONLY: NF90_OPEN, NF90_NOWRITE, NF90_CLOSE, nf90_noerr, &23 USE netcdf, ONLY: nf90_open, nf90_nowrite, nf90_close, nf90_noerr, & 24 24 nf90_inquire_dimension, nf90_inq_dimid, NF90_INQ_VARID, nf90_get_var 25 25 USE infotrac, ONLY: init_infotrac … … 154 154 ! weights to ensure ocean fractions are the same for atmosphere and ocean. 155 155 !******************************************************************************* 156 IF( NF90_OPEN("o2a.nc", NF90_NOWRITE, nid_o2a)==nf90_noerr) THEN157 iret= NF90_CLOSE(nid_o2a)156 IF(nf90_open("o2a.nc", nf90_nowrite, nid_o2a)==nf90_noerr) THEN 157 iret=nf90_close(nid_o2a) 158 158 WRITE(lunout,*)'BEWARE !! Ocean mask "o2a.nc" file found' 159 159 WRITE(lunout,*)'Coupled run.' … … 188 188 masque(iip1 ,:)=masque(1,:) 189 189 DEALLOCATE(ocemask) 190 ELSE IF( NF90_OPEN("startphy0.nc", NF90_NOWRITE, nid_sta)==nf90_noerr) THEN190 ELSE IF(nf90_open("startphy0.nc", nf90_nowrite, nid_sta)==nf90_noerr) THEN 191 191 WRITE(lunout,*)'BEWARE !! File "startphy0.nc" found.' 192 192 WRITE(lunout,*)'Getting the land mask from a previous run.' … … 196 196 WRITE(lunout,*)'Mismatching dimensions for land mask' 197 197 WRITE(lunout,*)'nphys = ',nphys ,' klon = ',klon 198 iret= NF90_CLOSE(nid_sta)198 iret=nf90_close(nid_sta) 199 199 CALL abort_gcm(modname,'',1) 200 200 END IF … … 202 202 iret=NF90_INQ_VARID(nid_sta,'masque',nid_msk) 203 203 iret=nf90_get_var(nid_sta,nid_msk,masktmp) 204 iret= NF90_CLOSE(nid_sta)204 iret=nf90_close(nid_sta) 205 205 CALL gr_fi_dyn(1,klon,iip1,jjp1,masktmp,masque) 206 206 IF(prt_level>=1) THEN -
LMDZ6/branches/Amaury_dev/libf/dynphy_lonlat/phylmd/limit_netcdf.f90
r5099 r5100 1 1 MODULE limit 2 2 3 !******************************************************************************* 4 ! Author : L. Fairhead, 27/01/94 5 !------------------------------------------------------------------------------- 6 ! Purpose: Boundary conditions files building for new model using climatologies. 7 ! Both grids have to be regular. 8 !------------------------------------------------------------------------------- 9 ! Note: This routine is designed to work for Earth 10 !------------------------------------------------------------------------------- 11 ! Modification history: 12 ! * 23/03/1994: Z. X. Li 13 ! * 09/1999: L. Fairhead (netcdf reading in LMDZ.3.3) 14 ! * 07/2001: P. Le Van 15 ! * 11/2009: L. Guez (ozone day & night climatos, see etat0_netcdf.F90) 16 ! * 12/2009: D. Cugnet (f77->f90, calendars, files from coupled runs) 17 !------------------------------------------------------------------------------- 18 19 USE ioipsl, ONLY : flininfo, flinopen, flinget, flinclo 20 USE assert_eq_m, ONLY : assert_eq 21 USE cal_tools_m, ONLY : year_len, mid_month 22 USE conf_dat_m, ONLY : conf_dat2d, conf_dat3d 23 USE dimphy, ONLY : klon, zmasq 24 USE geometry_mod, ONLY : longitude_deg, latitude_deg 25 USE phys_state_var_mod, ONLY : pctsrf 26 USE control_mod, ONLY : anneeref 27 USE init_ssrf_m, ONLY : start_init_subsurf 28 29 INTEGER, PARAMETER :: ns = 256 30 CHARACTER(LEN = ns), PARAMETER :: & 31 fsst(5) = ['amipbc_sst_1x1.nc ', 'amip_sst_1x1.nc ', 'cpl_atm_sst.nc '& 32 , 'histmth_sst.nc ', 'sstk.nc '], & 33 fsic(5) = ['amipbc_sic_1x1.nc ', 'amip_sic_1x1.nc ', 'cpl_atm_sic.nc '& 34 , 'histmth_sic.nc ', 'ci.nc '], & 35 vsst(5) = ['tosbcs ', 'tos ', 'SISUTESW ', 'tsol_oce ', 'sstk '], & 36 vsic(5) = ['sicbcs ', 'sic ', 'SIICECOV ', 'pourc_sic ', 'ci '], & 37 frugo = 'Rugos.nc ', falbe = 'Albedo.nc ', frelf = 'Relief.nc ', & 38 vrug = 'RUGOS ', valb = 'ALBEDO ', vrel = 'RELIEF ', & 39 DegK(11) = ['degK ', 'degree_K ', 'degreeK ', 'deg_K '& 40 , 'degsK ', 'degrees_K ', 'degreesK ', 'degs_K '& 41 , 'degree_kelvin ', 'degrees_kelvin', 'K '], & 42 DegC(10) = ['degC ', 'degree_C ', 'degreeC ', 'deg_C '& 43 , 'degsC ', 'degrees_C ', 'degreesC ', 'degs_C '& 44 , 'degree_Celsius', 'celsius '], & 45 Perc(2) = ['% ', 'percent '], & 46 Frac(2) = ['1.0 ', '1 '] 47 48 CONTAINS 49 50 !------------------------------------------------------------------------------- 51 52 SUBROUTINE limit_netcdf(masque, phis, extrap) 53 54 !------------------------------------------------------------------------------- 55 ! Author : L. Fairhead, 27/01/94 56 !------------------------------------------------------------------------------- 57 ! Purpose: Boundary conditions files building for new model using climatologies. 58 ! Both grids have to be regular. 59 !------------------------------------------------------------------------------- 60 ! Note: This routine is designed to work for Earth 61 !------------------------------------------------------------------------------- 62 ! Modification history: 63 ! * 23/03/1994: Z. X. Li 64 ! * 09/1999: L. Fairhead (netcdf reading in LMDZ.3.3) 65 ! * 07/2001: P. Le Van 66 ! * 11/2009: L. Guez (ozone day & night climatos, see etat0_netcdf.F90) 67 ! * 12/2009: D. Cugnet (f77->f90, calendars, files from coupled runs) 68 ! * 04/2016: D. Cugnet (12/14 recs SST/SIC files: cyclic/interannual runs) 69 ! * 05/2017: D. Cugnet (linear time interpolation for BCS files) 70 !------------------------------------------------------------------------------- 71 USE indice_sol_mod 72 USE netcdf, ONLY : nf90_open, nf90_create, nf90_close, & 73 nf90_def_dim, nf90_def_var, nf90_put_var, nf90_put_att, & 74 nf90_noerr, nf90_nowrite, nf90_global, & 75 nf90_clobber, nf90_enddef, nf90_unlimited, nf90_float, & 76 nf90_64bit_offset 77 USE lmdz_cppkeys_wrapper, ONLY : nf90_format 78 USE inter_barxy_m, ONLY : inter_barxy 79 USE netcdf95, ONLY : nf95_def_var, nf95_put_att, nf95_put_var 80 USE comconst_mod, ONLY : pi 81 USE phys_cal_mod, ONLY : calend 82 IMPLICIT NONE 83 !------------------------------------------------------------------------------- 84 ! Arguments: 85 include "iniprint.h" 86 include "dimensions.h" 87 include "paramet.h" 88 REAL, DIMENSION(iip1, jjp1), INTENT(INOUT) :: masque ! land mask 89 REAL, DIMENSION(iip1, jjp1), INTENT(INOUT) :: phis ! ground geopotential 90 LOGICAL, INTENT(IN) :: extrap ! SST extrapolation flag 91 !------------------------------------------------------------------------------- 92 ! Local variables: 93 include "comgeom2.h" 94 95 !--- INPUT NETCDF FILES AND VARIABLES NAMES ------------------------------------ 96 CHARACTER(LEN = ns) :: icefile, sstfile, fnam, varname 97 98 !--- OUTPUT VARIABLES FOR NETCDF FILE ------------------------------------------ 99 REAL :: fi_ice(klon) 100 REAL, POINTER :: phy_rug(:, :) => NULL(), phy_ice(:, :) => NULL() 101 REAL, POINTER :: phy_sst(:, :) => NULL(), phy_alb(:, :) => NULL() 102 REAL, ALLOCATABLE :: phy_bil(:, :), pctsrf_t(:, :, :) 103 INTEGER :: nbad 104 105 !--- VARIABLES FOR OUTPUT FILE WRITING ----------------------------------------- 106 INTEGER :: nid, ndim, ntim, k, dims(2), ix_sic, ix_sst 107 INTEGER :: id_tim, id_SST, id_BILS, id_RUG, id_ALB 108 INTEGER :: id_FOCE, id_FSIC, id_FTER, id_FLIC, varid_longitude, varid_latitude 109 INTEGER :: ndays !--- Depending on the output calendar 110 CHARACTER(LEN = ns) :: str 111 112 !--- INITIALIZATIONS ----------------------------------------------------------- 113 CALL inigeom 114 115 !--- MASK, GROUND GEOPOT. & SUBSURFACES COMPUTATION (IN CASE ok_etat0==.FALSE.) 116 IF(ALL(masque==-99999.)) THEN 117 CALL start_init_orog0(rlonv, rlatu, phis, masque) 118 CALL gr_dyn_fi(1, iip1, jjp1, klon, masque, zmasq) !--- To physical grid 119 ALLOCATE(pctsrf(klon, nbsrf)) 120 CALL start_init_subsurf(.FALSE.) 121 !--- TO MATCH EXACTLY WHAT WOULD BE DONE IN etat0phys_netcdf 122 WHERE(masque(:, :)<EPSFRA) masque(:, :) = 0. 123 WHERE(1. - masque(:, :)<EPSFRA) masque(:, :) = 1. 124 END IF 125 126 !--- Beware: anneeref (from gcm.def) is used to determine output time sampling 127 ndays = year_len(anneeref) 128 129 !--- RUGOSITY TREATMENT -------------------------------------------------------- 130 CALL msg(0, ""); CALL msg(0, " *** TRAITEMENT DE LA RUGOSITE ***") 131 CALL get_2Dfield(frugo, vrug, 'RUG', ndays, phy_rug, mask = masque(1:iim, :)) 132 133 !--- OCEAN TREATMENT ----------------------------------------------------------- 134 CALL msg(0, ""); CALL msg(0, " *** TRAITEMENT DE LA GLACE OCEANIQUE ***") 135 136 ! Input SIC file selection 137 ! Open file only to test if available 138 DO ix_sic = 1, SIZE(fsic) 139 IF (nf90_open(TRIM(fsic(ix_sic)), nf90_nowrite, nid)==nf90_noerr) THEN 140 icefile = fsic(ix_sic); varname = vsic(ix_sic); EXIT 141 END IF 142 END DO 143 IF(ix_sic==SIZE(fsic) + 1) THEN 144 WRITE(lunout, *) 'ERROR! No sea-ice input file was found.' 145 WRITE(lunout, *) 'One of following files must be available : ' 146 DO k = 1, SIZE(fsic); WRITE(lunout, *) TRIM(fsic(k)); 147 END DO 148 CALL abort_physic('limit_netcdf', 'No sea-ice file was found', 1) 149 END IF 150 CALL ncerr(nf90_close(nid), icefile) 151 CALL msg(0, 'Fichier choisi pour la glace de mer:' // TRIM(icefile)) 152 153 CALL get_2Dfield(icefile, varname, 'SIC', ndays, phy_ice) 154 155 ALLOCATE(pctsrf_t(klon, nbsrf, ndays)) 156 DO k = 1, ndays 157 fi_ice = phy_ice(:, k) 158 WHERE(fi_ice>=1.0) fi_ice = 1.0 159 WHERE(fi_ice<EPSFRA) fi_ice = 0.0 160 pctsrf_t(:, is_ter, k) = pctsrf(:, is_ter) ! land soil 161 pctsrf_t(:, is_lic, k) = pctsrf(:, is_lic) ! land ice 162 SELECT CASE(ix_sic) 163 CASE(3) ! SIC=pICE*(1-LIC-TER) (CPL) 164 pctsrf_t(:, is_sic, k) = fi_ice(:) * (1. - pctsrf(:, is_lic) - pctsrf(:, is_ter)) 165 CASE(4) ! SIC=pICE (HIST) 166 pctsrf_t(:, is_sic, k) = fi_ice(:) 167 CASE DEFAULT ! SIC=pICE-LIC (AMIP,ERAI) 168 pctsrf_t(:, is_sic, k) = fi_ice - pctsrf_t(:, is_lic, k) 169 END SELECT 170 WHERE(pctsrf_t(:, is_sic, k)<=0) pctsrf_t(:, is_sic, k) = 0. 171 WHERE(1.0 - zmasq<EPSFRA) 172 pctsrf_t(:, is_sic, k) = 0.0 173 pctsrf_t(:, is_oce, k) = 0.0 174 ELSEWHERE 175 WHERE(pctsrf_t(:, is_sic, k)>=1.0 - zmasq) 176 pctsrf_t(:, is_sic, k) = 1.0 - zmasq 177 pctsrf_t(:, is_oce, k) = 0.0 178 ELSEWHERE 179 pctsrf_t(:, is_oce, k) = 1.0 - zmasq - pctsrf_t(:, is_sic, k) 180 WHERE(pctsrf_t(:, is_oce, k)<EPSFRA) 181 pctsrf_t(:, is_oce, k) = 0.0 182 pctsrf_t(:, is_sic, k) = 1.0 - zmasq 183 END WHERE 184 END WHERE 185 END WHERE 186 nbad = COUNT(pctsrf_t(:, is_oce, k)<0.0) 187 IF(nbad>0) WRITE(lunout, *) 'pb sous maille pour nb points = ', nbad 188 nbad = COUNT(ABS(SUM(pctsrf_t(:, :, k), DIM = 2) - 1.0)>EPSFRA) 189 IF(nbad>0) WRITE(lunout, *) 'pb sous surface pour nb points = ', nbad 190 END DO 191 DEALLOCATE(phy_ice) 192 193 !--- SST TREATMENT ------------------------------------------------------------- 194 CALL msg(0, ""); CALL msg(0, " *** TRAITEMENT DE LA SST ***") 195 196 ! Input SST file selection 197 ! Open file only to test if available 198 DO ix_sst = 1, SIZE(fsst) 199 IF (nf90_open(TRIM(fsst(ix_sst)), nf90_nowrite, nid)==nf90_noerr) THEN 200 sstfile = fsst(ix_sst); varname = vsst(ix_sst); EXIT 201 END IF 202 END DO 203 IF(ix_sst==SIZE(fsst) + 1) THEN 204 WRITE(lunout, *) 'ERROR! No sst input file was found.' 205 WRITE(lunout, *) 'One of following files must be available : ' 206 DO k = 1, SIZE(fsst); WRITE(lunout, *) TRIM(fsst(k)); 207 END DO 208 CALL abort_physic('limit_netcdf', 'No sst file was found', 1) 209 END IF 210 CALL ncerr(nf90_close(nid), sstfile) 211 CALL msg(0, 'Fichier choisi pour la temperature de mer: ' // TRIM(sstfile)) 212 213 CALL get_2Dfield(sstfile, varname, 'SST', ndays, phy_sst, flag = extrap) 214 215 !--- ALBEDO TREATMENT ---------------------------------------------------------- 216 CALL msg(0, ""); CALL msg(0, " *** TRAITEMENT DE L'ALBEDO ***") 217 CALL get_2Dfield(falbe, valb, 'ALB', ndays, phy_alb) 218 219 !--- REFERENCE GROUND HEAT FLUX TREATMENT -------------------------------------- 220 ALLOCATE(phy_bil(klon, ndays)); phy_bil = 0.0 221 222 !--- OUTPUT FILE WRITING ------------------------------------------------------- 223 CALL msg(0, ""); CALL msg(0, ' *** Ecriture du fichier limit : debut ***') 224 fnam = "limit.nc" 225 226 !--- File creation 227 CALL ncerr(nf90_create(fnam, IOR(nf90_clobber, nf90_64bit_offset), nid), fnam) 228 CALL ncerr(nf90_put_att(nid, nf90_global, "title", "Fichier conditions aux limites"), fnam) 229 str = 'File produced using ce0l executable.' 230 str = TRIM(str) // NEW_LINE(' ') // 'Sea Ice Concentration built from' 231 SELECT CASE(ix_sic) 232 CASE(1); str = TRIM(str) // ' Amip mid-month boundary condition (BCS).' 233 CASE(2); str = TRIM(str) // ' Amip monthly mean observations.' 234 CASE(3); str = TRIM(str) // ' IPSL coupled model outputs.' 235 CASE(4); str = TRIM(str) // ' LMDZ model outputs.' 236 CASE(5); str = TRIM(str) // ' ci.nc file.' 237 END SELECT 238 str = TRIM(str) // NEW_LINE(' ') // 'Sea Surface Temperature built from' 239 SELECT CASE(ix_sst) 240 CASE(1); str = TRIM(str) // ' Amip mid-month boundary condition (BCS).' 241 CASE(2); str = TRIM(str) // ' Amip monthly mean observations.' 242 CASE(3); str = TRIM(str) // ' IPSL coupled model outputs.' 243 CASE(4); str = TRIM(str) // ' LMDZ model outputs.' 244 CASE(5); str = TRIM(str) // ' sstk.nc file.' 245 END SELECT 246 CALL ncerr(nf90_put_att(nid, nf90_global, "history", TRIM(str)), fnam) 247 248 !--- Dimensions creation 249 CALL ncerr(nf90_def_dim(nid, "points_physiques", klon, ndim), fnam) 250 CALL ncerr(nf90_def_dim(nid, "time", nf90_unlimited, ntim), fnam) 251 252 dims = [ndim, ntim] 253 254 !--- Variables creation 255 CALL ncerr(nf90_def_var(nid, "TEMPS", nf90_format, [ntim], id_tim), fnam) 256 CALL ncerr(nf90_def_var(nid, "FOCE", nf90_format, dims, id_FOCE), fnam) 257 CALL ncerr(nf90_def_var(nid, "FSIC", nf90_format, dims, id_FSIC), fnam) 258 CALL ncerr(nf90_def_var(nid, "FTER", nf90_format, dims, id_FTER), fnam) 259 CALL ncerr(nf90_def_var(nid, "FLIC", nf90_format, dims, id_FLIC), fnam) 260 CALL ncerr(nf90_def_var(nid, "SST", nf90_format, dims, id_SST), fnam) 261 CALL ncerr(nf90_def_var(nid, "BILS", nf90_format, dims, id_BILS), fnam) 262 CALL ncerr(nf90_def_var(nid, "ALB", nf90_format, dims, id_ALB), fnam) 263 CALL ncerr(nf90_def_var(nid, "RUG", nf90_format, dims, id_RUG), fnam) 264 call nf95_def_var(nid, "longitude", nf90_float, ndim, varid_longitude) 265 call nf95_def_var(nid, "latitude", nf90_float, ndim, varid_latitude) 266 267 !--- Attributes creation 268 CALL ncerr(nf90_put_att(nid, id_tim, "title", "Jour dans l annee"), fnam) 269 CALL ncerr(nf90_put_att(nid, id_tim, "calendar", calend), fnam) 270 CALL ncerr(nf90_put_att(nid, id_FOCE, "title", "Fraction ocean"), fnam) 271 CALL ncerr(nf90_put_att(nid, id_FSIC, "title", "Fraction glace de mer"), fnam) 272 CALL ncerr(nf90_put_att(nid, id_FTER, "title", "Fraction terre"), fnam) 273 CALL ncerr(nf90_put_att(nid, id_FLIC, "title", "Fraction land ice"), fnam) 274 CALL ncerr(nf90_put_att(nid, id_SST, "title", "Temperature superficielle de la mer"), fnam) 275 CALL ncerr(nf90_put_att(nid, id_BILS, "title", "Reference flux de chaleur au sol"), fnam) 276 CALL ncerr(nf90_put_att(nid, id_ALB, "title", "Albedo a la surface"), fnam) 277 CALL ncerr(nf90_put_att(nid, id_RUG, "title", "Rugosite"), fnam) 278 279 call nf95_put_att(nid, varid_longitude, "standard_name", "longitude") 280 call nf95_put_att(nid, varid_longitude, "units", "degrees_east") 281 282 call nf95_put_att(nid, varid_latitude, "standard_name", "latitude") 283 call nf95_put_att(nid, varid_latitude, "units", "degrees_north") 284 285 CALL ncerr(nf90_enddef(nid), fnam) 286 287 !--- Variables saving 288 CALL ncerr(nf90_put_var(nid, id_tim, [(REAL(k), k = 1, ndays)]), fnam) 289 CALL ncerr(nf90_put_var(nid, id_FOCE, pctsrf_t(:, is_oce, :), [1, 1], [klon, ndays]), fnam) 290 CALL ncerr(nf90_put_var(nid, id_FSIC, pctsrf_t(:, is_sic, :), [1, 1], [klon, ndays]), fnam) 291 CALL ncerr(nf90_put_var(nid, id_FTER, pctsrf_t(:, is_ter, :), [1, 1], [klon, ndays]), fnam) 292 CALL ncerr(nf90_put_var(nid, id_FLIC, pctsrf_t(:, is_lic, :), [1, 1], [klon, ndays]), fnam) 293 CALL ncerr(nf90_put_var(nid, id_SST, phy_sst(:, :), [1, 1], [klon, ndays]), fnam) 294 CALL ncerr(nf90_put_var(nid, id_BILS, phy_bil(:, :), [1, 1], [klon, ndays]), fnam) 295 CALL ncerr(nf90_put_var(nid, id_ALB, phy_alb(:, :), [1, 1], [klon, ndays]), fnam) 296 CALL ncerr(nf90_put_var(nid, id_RUG, phy_rug(:, :), [1, 1], [klon, ndays]), fnam) 297 call nf95_put_var(nid, varid_longitude, longitude_deg) 298 call nf95_put_var(nid, varid_latitude, latitude_deg) 299 300 CALL ncerr(nf90_close(nid), fnam) 301 302 CALL msg(0, ""); CALL msg(0, ' *** Ecriture du fichier limit : fin ***') 303 304 DEALLOCATE(pctsrf_t, phy_sst, phy_bil, phy_alb, phy_rug) 305 306 307 !=============================================================================== 308 309 CONTAINS 310 311 !=============================================================================== 312 313 314 !------------------------------------------------------------------------------- 315 316 SUBROUTINE get_2Dfield(fnam, varname, mode, ndays, champo, flag, mask) 317 318 !----------------------------------------------------------------------------- 319 ! Comments: 320 ! There are two assumptions concerning the NetCDF files, that are satisfied 321 ! with files that are conforming NC convention: 322 ! 1) The last dimension of the variables used is the time record. 323 ! 2) Dimensional variables have the same names as corresponding dimensions. 324 !----------------------------------------------------------------------------- 325 USE netcdf, ONLY : nf90_open, NF90_INQ_VARID, NF90_INQUIRE_VARIABLE, & 326 nf90_close, nf90_inq_dimid, nf90_inquire_dimension, nf90_get_var, & 327 nf90_get_att 328 USE pchsp_95_m, only : pchsp_95 329 USE pchfe_95_m, only : pchfe_95 330 USE arth_m, only : arth 331 USE indice_sol_mod 332 333 IMPLICIT NONE 334 include "dimensions.h" 335 include "paramet.h" 336 include "comgeom2.h" 337 !----------------------------------------------------------------------------- 338 ! Arguments: 339 CHARACTER(LEN = *), INTENT(IN) :: fnam ! NetCDF file name 340 CHARACTER(LEN = *), INTENT(IN) :: varname ! NetCDF variable name 341 CHARACTER(LEN = *), INTENT(IN) :: mode ! RUG, SIC, SST or ALB 342 INTEGER, INTENT(IN) :: ndays ! current year number of days 343 REAL, POINTER, DIMENSION(:, :) :: champo ! output field = f(t) 344 LOGICAL, OPTIONAL, INTENT(IN) :: flag ! extrapol. (SST) old ice (SIC) 345 REAL, OPTIONAL, DIMENSION(iim, jjp1), INTENT(IN) :: mask 346 !------------------------------------------------------------------------------ 347 ! Local variables: 348 !--- NetCDF 349 INTEGER :: ncid, varid ! NetCDF identifiers 350 CHARACTER(LEN = ns) :: dnam ! dimension name 351 !--- dimensions 352 INTEGER :: dids(4) ! NetCDF dimensions identifiers 353 REAL, ALLOCATABLE :: dlon_ini(:) ! initial longitudes vector 354 REAL, ALLOCATABLE :: dlat_ini(:) ! initial latitudes vector 355 REAL, POINTER :: dlon(:), dlat(:) ! reordered lon/lat vectors 356 !--- fields 357 INTEGER :: imdep, jmdep, lmdep ! dimensions of 'champ' 358 REAL, ALLOCATABLE :: champ(:, :) ! wanted field on initial grid 359 REAL, ALLOCATABLE :: yder(:), timeyear(:) 360 REAL :: champint(iim, jjp1) ! interpolated field 361 REAL, ALLOCATABLE :: champtime(:, :, :) 362 REAL, ALLOCATABLE :: champan(:, :, :) 363 !--- input files 364 CHARACTER(LEN = ns) :: fnam_m, fnam_p ! previous/next files names 365 CHARACTER(LEN = ns) :: cal_in ! calendar 366 CHARACTER(LEN = ns) :: units ! attribute "units" in sic/sst file 367 INTEGER :: ndays_in ! number of days 368 REAL :: value ! mean/max value near equator 369 !--- misc 370 INTEGER :: i, j, k, l ! loop counters 371 REAL, ALLOCATABLE :: work(:, :) ! used for extrapolation 372 CHARACTER(LEN = ns) :: title, mess ! for messages 373 LOGICAL :: is_bcs ! flag for BCS data 374 LOGICAL :: extrp ! flag for extrapolation 375 LOGICAL :: ll 376 REAL :: chmin, chmax, timeday, al 377 INTEGER ierr, idx 378 integer n_extrap ! number of extrapolated points 379 logical skip 380 381 !------------------------------------------------------------------------------ 382 !---Variables depending on keyword 'mode' ------------------------------------- 383 NULLIFY(champo) 384 385 SELECT CASE(mode) 386 CASE('RUG'); title = 'Rugosite' 387 CASE('SIC'); title = 'Sea-ice' 388 CASE('SST'); title = 'SST' 389 CASE('ALB'); title = 'Albedo' 390 END SELECT 391 extrp = .FALSE.; IF(PRESENT(flag).AND.mode=='SST') extrp = flag 392 is_bcs = (mode=='SIC'.AND.ix_sic==1).OR.(mode=='SST'.AND.ix_sst==1) 393 idx = INDEX(fnam, '.nc') - 1 394 395 !--- GETTING SOME DIMENSIONAL VARIABLES FROM FILE ----------------------------- 396 CALL msg(5, ' Now reading file : ' // TRIM(fnam)) 397 CALL ncerr(nf90_open(fnam, nf90_nowrite, ncid), fnam) 398 CALL ncerr(NF90_INQ_VARID(ncid, trim(varname), varid), fnam) 399 CALL ncerr(NF90_INQUIRE_VARIABLE(ncid, varid, dimids = dids), fnam) 400 401 !--- Longitude 402 CALL ncerr(nf90_inquire_dimension(ncid, dids(1), name = dnam, len = imdep), fnam) 403 ALLOCATE(dlon_ini(imdep), dlon(imdep)) 404 CALL ncerr(NF90_INQ_VARID(ncid, dnam, varid), fnam) 405 CALL ncerr(nf90_get_var(ncid, varid, dlon_ini), fnam) 406 CALL msg(5, 'variable ' // TRIM(dnam) // ' dimension ', imdep) 407 408 !--- Latitude 409 CALL ncerr(nf90_inquire_dimension(ncid, dids(2), name = dnam, len = jmdep), fnam) 410 ALLOCATE(dlat_ini(jmdep), dlat(jmdep)) 411 CALL ncerr(NF90_INQ_VARID(ncid, dnam, varid), fnam) 412 CALL ncerr(nf90_get_var(ncid, varid, dlat_ini), fnam) 413 CALL msg(5, 'variable ' // TRIM(dnam) // ' dimension ', jmdep) 414 415 !--- Time (variable is not needed - it is rebuilt - but calendar is) 416 CALL ncerr(nf90_inquire_dimension(ncid, dids(3), name = dnam, len = lmdep), fnam) 417 ALLOCATE(timeyear(lmdep + 2)) 418 CALL ncerr(NF90_INQ_VARID(ncid, dnam, varid), fnam) 419 cal_in = ' ' 420 IF(nf90_get_att(ncid, varid, 'calendar', cal_in)/=nf90_noerr) THEN 421 SELECT CASE(mode) 422 CASE('RUG', 'ALB'); cal_in = '360_day' 423 CASE('SIC', 'SST'); cal_in = 'gregorian' 424 END SELECT 425 CALL msg(0, 'WARNING: missing "calendar" attribute for "time" in '& 426 // TRIM(fnam) // '. Choosing default value.') 427 END IF 428 CALL strclean(cal_in) !--- REMOVE (WEIRD) NULL CHARACTERS 429 CALL msg(0, 'var, calendar, dim: ' // TRIM(dnam) // ' ' // TRIM(cal_in), lmdep) 430 431 !--- Determining input file number of days, depending on calendar 432 ndays_in = year_len(anneeref, cal_in) 433 434 !--- Rebuilding input time vector (field from input file might be unreliable) 435 IF(lmdep==12) THEN 436 timeyear = mid_month(anneeref, cal_in) 437 CALL msg(0, 'Monthly input file(s) for ' // TRIM(title) // '.') 438 ELSE IF(lmdep==ndays_in) THEN 439 timeyear = [(REAL(k) - 0.5, k = 0, ndays_in + 1)] 440 CALL msg(0, 'Daily input file (no time interpolation).') 441 ELSE 442 WRITE(mess, '(a,i3,a,i3,a)')'Mismatching input file: found', lmdep, & 443 ' records, 12/', ndays_in, ' (monthly/daily needed).' 444 CALL abort_physic('mid_month', TRIM(mess), 1) 445 END IF 446 447 !--- GETTING THE FIELD AND INTERPOLATING IT ---------------------------------- 448 ALLOCATE(champ(imdep, jmdep), champtime(iim, jjp1, lmdep + 2)) 449 IF(extrp) ALLOCATE(work(imdep, jmdep)) 450 CALL msg(5, '') 451 CALL msg(5, 'READ AND INTERPOLATE HORIZONTALLY ', lmdep, ' FIELDS.') 452 CALL ncerr(NF90_INQ_VARID(ncid, varname, varid), fnam) 453 DO l = 1, lmdep 454 CALL ncerr(nf90_get_var(ncid, varid, champ, [1, 1, l], [imdep, jmdep, 1]), fnam) 455 CALL conf_dat2d(title, dlon_ini, dlat_ini, dlon, dlat, champ, .TRUE.) 456 457 !--- FOR SIC/SST FIELDS ONLY 458 IF(l==1.AND.is_in(mode, ['SIC', 'SST'])) THEN 459 460 !--- DETERMINE THE UNIT: READ FROM FILE OR ASSUMED USING FIELD VALUES 461 ierr = nf90_get_att(ncid, varid, 'units', units) 462 IF(ierr==nf90_noerr) THEN !--- ATTRIBUTE "units" FOUND IN THE FILE 463 CALL strclean(units) 464 IF(mode=='SIC'.AND.is_in(units, Perc)) units = "%" 465 IF(mode=='SIC'.AND.is_in(units, Frac)) units = "1" 466 IF(mode=='SST'.AND.is_in(units, DegC)) units = "C" 467 IF(mode=='SST'.AND.is_in(units, DegK)) units = "K" 468 ELSE !--- CHECK THE FIELD VALUES 469 IF(mode=='SIC') value = MAXVAL(champ(:, :)) 470 IF(mode=='SST') value = SUM(champ(:, jmdep / 2), DIM = 1) / REAL(imdep) 471 IF(mode=='SIC') THEN; units = "1"; IF(value>= 10.) units = "%"; 472 END IF 473 IF(mode=='SST') THEN; units = "C"; IF(value>=100.) units = "K"; 474 END IF 475 END IF 476 CALL msg(0, 'INPUT FILE ' // TRIM(title) // ' UNIT IS: "' // TRIM(units) // '".') 477 IF(ierr/=nf90_noerr) CALL msg(0, 'WARNING ! UNIT TO BE CHECKED ! ' & 478 // 'No "units" attribute, so only based on the fields values.') 479 480 !--- CHECK VALUES ARE IN THE EXPECTED RANGE 481 SELECT CASE(units) 482 CASE('%'); ll = ANY(champ>100.0 + EPSFRA); str = 'percentages > 100.' 483 CASE('1'); ll = ANY(champ> 1.0 + EPSFRA); str = 'fractions > 1.' 484 CASE('C'); ll = ANY(champ<-100.).OR.ANY(champ> 60.); str = '<-100 or >60 DegC' 485 CASE('K'); ll = ANY(champ< 180.).OR.ANY(champ>330.); str = '<180 or >330 DegK' 486 CASE DEFAULT; CALL abort_physic(mode, 'Unrecognized ' // TRIM(title) & 487 // ' unit: ' // TRIM(units), 1) 488 END SELECT 489 490 !--- DROPPED FOR BCS DATA (FRACTIONS CAN BE HIGHER THAN 1) 491 IF(ll.AND.ix_sic/=1.AND.mode=='SIC') & 492 CALL abort_physic(mode, 'unrealistic ' // TRIM(mode) // ' found: ' // TRIM(str), 1) 493 494 END IF 495 496 IF(extrp) CALL extrapol(champ, imdep, jmdep, 999999., .TRUE., .TRUE., 2, work) 497 IF(l==1) THEN 498 CALL msg(5, "--------------------------------------------------------") 499 CALL msg(5, "$$$ Barycentric interpolation for " // TRIM(title) // " $$$") 500 CALL msg(5, "--------------------------------------------------------") 501 END IF 502 IF(mode=='RUG') champ = LOG(champ) 503 CALL inter_barxy(dlon, dlat(:jmdep - 1), champ, rlonu(:iim), rlatv, champint) 504 IF(mode=='RUG') THEN 505 champint = EXP(champint) 506 WHERE(NINT(mask)/=1) champint = 0.001 507 END IF 508 champtime(:, :, l + 1) = champint 509 END DO 510 CALL ncerr(nf90_close(ncid), fnam) 511 512 !--- FIRST RECORD: LAST ONE OF PREVIOUS YEAR (CURRENT YEAR IF UNAVAILABLE) 513 fnam_m = fnam(1:idx) // '_m.nc' 514 IF(nf90_open(fnam_m, nf90_nowrite, ncid)==nf90_noerr) THEN 515 CALL msg(0, 'Reading previous year file ("' // TRIM(fnam_m) // '") last record for ' // TRIM(title)) 516 CALL ncerr(NF90_INQ_VARID(ncid, varname, varid), fnam_m) 517 CALL ncerr(NF90_INQUIRE_VARIABLE(ncid, varid, dimids = dids), fnam_m) 518 CALL ncerr(nf90_inquire_dimension(ncid, dids(3), len = l), fnam_m) 519 CALL ncerr(nf90_get_var(ncid, varid, champ, [1, 1, l], [imdep, jmdep, 1]), fnam_m) 520 CALL ncerr(nf90_close(ncid), fnam_m) 521 CALL conf_dat2d(title, dlon_ini, dlat_ini, dlon, dlat, champ, .TRUE.) 522 IF(extrp) CALL extrapol(champ, imdep, jmdep, 999999., .TRUE., .TRUE., 2, work) 523 IF(mode=='RUG') champ = LOG(champ) 524 CALL inter_barxy(dlon, dlat(:jmdep - 1), champ, rlonu(:iim), rlatv, champint) 525 IF(mode=='RUG') THEN 526 champint = EXP(champint) 527 WHERE(NINT(mask)/=1) champint = 0.001 528 END IF 529 champtime(:, :, 1) = champint 530 ELSE 531 CALL msg(0, 'Using current year file ("' // TRIM(fnam) // '") last record for ' // TRIM(title)) 532 champtime(:, :, 1) = champtime(:, :, lmdep + 1) 533 END IF 534 535 !--- LAST RECORD: FIRST ONE OF NEXT YEAR (CURRENT YEAR IF UNAVAILABLE) 536 fnam_p = fnam(1:idx) // '_p.nc' 537 IF(nf90_open(fnam_p, nf90_nowrite, ncid)==nf90_noerr) THEN 538 CALL msg(0, 'Reading next year file ("' // TRIM(fnam_p) // '") first record for ' // TRIM(title)) 539 CALL ncerr(NF90_INQ_VARID(ncid, varname, varid), fnam_p) 540 CALL ncerr(nf90_get_var(ncid, varid, champ, [1, 1, 1], [imdep, jmdep, 1]), fnam_p) 541 CALL ncerr(nf90_close(ncid), fnam_p) 542 CALL conf_dat2d(title, dlon_ini, dlat_ini, dlon, dlat, champ, .TRUE.) 543 IF(extrp) CALL extrapol(champ, imdep, jmdep, 999999., .TRUE., .TRUE., 2, work) 544 IF(mode=='RUG') champ = LOG(champ) 545 CALL inter_barxy(dlon, dlat(:jmdep - 1), champ, rlonu(:iim), rlatv, champint) 546 IF(mode=='RUG') THEN 547 champint = EXP(champint) 548 WHERE(NINT(mask)/=1) champint = 0.001 549 END IF 550 champtime(:, :, lmdep + 2) = champint 551 ELSE 552 CALL msg(0, 'Using current year file ("' // TRIM(fnam) // '") first record for ' // TRIM(title)) 553 champtime(:, :, lmdep + 2) = champtime(:, :, 2) 554 END IF 555 DEALLOCATE(dlon_ini, dlat_ini, dlon, dlat, champ) 556 IF(extrp) DEALLOCATE(work) 557 558 !--- TIME INTERPOLATION ------------------------------------------------------ 559 IF(prt_level>0) THEN 560 IF(ndays/=ndays_in) THEN 561 WRITE(lunout, *)'DIFFERENT YEAR LENGTHS:' 562 WRITE(lunout, *)' In the input file: ', ndays_in 563 WRITE(lunout, *)' In the output file: ', ndays 564 END IF 565 IF(lmdep==ndays_in) THEN 566 WRITE(lunout, *)'NO TIME INTERPOLATION.' 567 WRITE(lunout, *)' Daily input file.' 568 ELSE 569 IF(is_bcs) WRITE(lunout, *)'LINEAR TIME INTERPOLATION.' 570 IF(.NOT.is_bcs) WRITE(lunout, *)'SPLINES TIME INTERPOLATION.' 571 WRITE(lunout, *)' Input time vector: ', timeyear 572 WRITE(lunout, *)' Output time vector: from 0.5 to ', ndays - 0.5 573 END IF 574 END IF 575 ALLOCATE(champan(iip1, jjp1, ndays)) 576 577 IF(lmdep==ndays_in) THEN !--- DAILY DATA: NO TIME INTERPOLATION 578 DO l = 1, lmdep 579 champan(1:iim, :, l) = champtime(:, :, l + 1) 580 END DO 581 ELSE IF(is_bcs) THEN !--- BCS DATA: LINEAR TIME INTERPOLATION 582 l = 1 583 DO k = 1, ndays 584 timeday = (REAL(k) - 0.5) * REAL(ndays_in) / ndays 585 IF(timeyear(l + 1)<timeday) l = l + 1 586 al = (timeday - timeyear(l)) / (timeyear(l + 1) - timeyear(l)) 587 champan(1:iim, :, k) = champtime(1:iim, :, l) + al * (champtime(1:iim, :, l + 1) - champtime(1:iim, :, l)) 588 END DO 589 ELSE !--- AVE DATA: SPLINE TIME INTERPOLATION 590 skip = .false. 591 n_extrap = 0 592 ALLOCATE(yder(lmdep + 2)) 593 DO j = 1, jjp1 594 DO i = 1, iim 595 yder = pchsp_95(timeyear, champtime(i, j, :), ibeg = 2, iend = 2, & 596 vc_beg = 0., vc_end = 0.) 597 CALL pchfe_95(timeyear, champtime(i, j, :), yder, skip, & 598 arth(0.5, real(ndays_in) / ndays, ndays), champan(i, j, :), ierr) 599 if (ierr < 0) call abort_physic("get_2Dfield", "", 1) 600 n_extrap = n_extrap + ierr 601 END DO 602 END DO 603 IF(n_extrap /= 0) WRITE(lunout, *) "get_2Dfield pchfe_95: n_extrap = ", n_extrap 604 DEALLOCATE(yder) 605 END IF 606 champan(iip1, :, :) = champan(1, :, :) 607 DEALLOCATE(champtime, timeyear) 608 609 !--- Checking the result 610 DO j = 1, jjp1 611 CALL minmax(iip1, champan(1, j, 10), chmin, chmax) 612 IF (prt_level>5) WRITE(lunout, *)' ', TRIM(title), ' at time 10 ', chmin, chmax, j 613 END DO 614 615 !--- SPECIAL FILTER FOR SST: SST>271.38 -------------------------------------- 616 IF(mode=='SST') THEN 617 SELECT CASE(units) 618 CASE("K"); CALL msg(0, 'SST field is already in kelvins.') 619 CASE("C"); CALL msg(0, 'SST field converted from celcius degrees to kelvins.') 620 champan(:, :, :) = champan(:, :, :) + 273.15 621 END SELECT 622 CALL msg(0, 'Filtering SST: Sea Surface Temperature >= 271.38') 623 WHERE(champan<271.38) champan = 271.38 624 END IF 625 626 !--- SPECIAL FILTER FOR SIC: 0.0<SIC<1.0 ------------------------------------- 627 IF(mode=='SIC') THEN 628 SELECT CASE(units) 629 CASE("1"); CALL msg(0, 'SIC field already in fraction of 1') 630 CASE("%"); CALL msg(0, 'SIC field converted from percentage to fraction of 1.') 631 champan(:, :, :) = champan(:, :, :) / 100. 632 END SELECT 633 CALL msg(0, 'Filtering SIC: 0.0 <= Sea-ice <=1.0') 634 WHERE(champan>1.0) champan = 1.0 635 WHERE(champan<0.0) champan = 0.0 636 END IF 637 638 !--- DYNAMICAL TO PHYSICAL GRID ---------------------------------------------- 639 ALLOCATE(champo(klon, ndays)) 640 DO k = 1, ndays 641 CALL gr_dyn_fi(1, iip1, jjp1, klon, champan(1, 1, k), champo(:, k)) 642 END DO 643 DEALLOCATE(champan) 644 645 END SUBROUTINE get_2Dfield 646 647 !------------------------------------------------------------------------------- 648 649 650 !------------------------------------------------------------------------------- 651 652 SUBROUTINE start_init_orog0(lon_in, lat_in, phis, masque) 653 654 !------------------------------------------------------------------------------- 655 USE grid_noro_m, ONLY : grid_noro0 656 IMPLICIT NONE 657 !=============================================================================== 658 ! Purpose: Compute "phis" just like it would be in start_init_orog. 659 !=============================================================================== 660 ! Arguments: 661 REAL, INTENT(IN) :: lon_in(:), lat_in(:) ! dim (iml) (jml) 662 REAL, INTENT(INOUT) :: phis(:, :), masque(:, :) ! dim (iml,jml) 663 !------------------------------------------------------------------------------- 664 ! Local variables: 665 CHARACTER(LEN = ns) :: modname = "start_init_orog0" 666 INTEGER :: fid, llm_tmp, ttm_tmp, iml, jml, iml_rel, jml_rel, itau(1) 667 REAL :: lev(1), date, dt, deg2rad 668 REAL, ALLOCATABLE :: lon_rad(:), lon_ini(:), lon_rel(:, :), relief_hi(:, :) 669 REAL, ALLOCATABLE :: lat_rad(:), lat_ini(:), lat_rel(:, :) 670 !------------------------------------------------------------------------------- 671 iml = assert_eq(SIZE(lon_in), SIZE(phis, 1), SIZE(masque, 1), TRIM(modname) // " iml") 672 jml = assert_eq(SIZE(lat_in), SIZE(phis, 2), SIZE(masque, 2), TRIM(modname) // " jml") 673 IF(iml/=iip1) CALL abort_gcm(TRIM(modname), 'iml/=iip1', 1) 674 IF(jml/=jjp1) CALL abort_gcm(TRIM(modname), 'jml/=jjp1', 1) 675 pi = 2.0 * ASIN(1.0); deg2rad = pi / 180.0 676 IF(ANY(phis/=-99999.)) RETURN !--- phis ALREADY KNOWN 677 678 !--- HIGH RESOLUTION OROGRAPHY 679 CALL flininfo(frelf, iml_rel, jml_rel, llm_tmp, ttm_tmp, fid) 680 681 ALLOCATE(lat_rel(iml_rel, jml_rel), lon_rel(iml_rel, jml_rel)) 682 CALL flinopen(frelf, .FALSE., iml_rel, jml_rel, llm_tmp, lon_rel, lat_rel, & 683 lev, ttm_tmp, itau, date, dt, fid) 684 ALLOCATE(relief_hi(iml_rel, jml_rel)) 685 CALL flinget(fid, vrel, iml_rel, jml_rel, llm_tmp, ttm_tmp, 1, 1, relief_hi) 686 CALL flinclo(fid) 687 688 !--- IF ANGLES ARE IN DEGREES, THEY ARE CONVERTED INTO RADIANS 689 ALLOCATE(lon_ini(iml_rel), lat_ini(jml_rel)) 690 lon_ini(:) = lon_rel(:, 1); IF(MAXVAL(lon_rel)>pi) lon_ini = lon_ini * deg2rad 691 lat_ini(:) = lat_rel(1, :); IF(MAXVAL(lat_rel)>pi) lat_ini = lat_ini * deg2rad 692 693 !--- FIELDS ARE PROCESSED TO BE ON STANDARD ANGULAR DOMAINS 694 ALLOCATE(lon_rad(iml_rel), lat_rad(jml_rel)) 695 CALL conf_dat2d(vrel, lon_ini, lat_ini, lon_rad, lat_rad, relief_hi, .FALSE.) 696 DEALLOCATE(lon_ini, lat_ini) 697 698 !--- COMPUTING SURFACE GEOPOTENTIAL USING ROUTINE grid_noro0 699 WRITE(lunout, *) 700 WRITE(lunout, *)'*** Compute surface geopotential ***' 701 702 !--- CALL OROGRAPHY MODULE (REDUCED VERSION) TO COMPUTE FIELDS 703 CALL grid_noro0(lon_rad, lat_rad, relief_hi, lon_in, lat_in, phis, masque) 704 phis = phis * 9.81 705 phis(iml, :) = phis(1, :) 706 DEALLOCATE(relief_hi, lon_rad, lat_rad) 707 708 END SUBROUTINE start_init_orog0 709 710 !------------------------------------------------------------------------------- 711 712 713 !------------------------------------------------------------------------------- 714 715 SUBROUTINE msg(lev, str1, i, str2) 716 717 !------------------------------------------------------------------------------- 718 ! Arguments: 719 INTEGER, INTENT(IN) :: lev 720 CHARACTER(LEN = *), INTENT(IN) :: str1 721 INTEGER, OPTIONAL, INTENT(IN) :: i 722 CHARACTER(LEN = *), OPTIONAL, INTENT(IN) :: str2 723 !------------------------------------------------------------------------------- 724 IF(prt_level>=lev) THEN 725 IF(PRESENT(str2)) THEN 726 WRITE(lunout, *) TRIM(str1), i, TRIM(str2) 727 ELSE IF(PRESENT(i)) THEN 728 WRITE(lunout, *) TRIM(str1), i 729 ELSE 730 WRITE(lunout, *) TRIM(str1) 731 END IF 732 END IF 733 734 END SUBROUTINE msg 735 736 !------------------------------------------------------------------------------- 737 738 739 !------------------------------------------------------------------------------- 740 741 SUBROUTINE ncerr(ncres, fnam) 742 743 !------------------------------------------------------------------------------- 744 ! Purpose: NetCDF errors handling. 745 !------------------------------------------------------------------------------- 746 USE netcdf, ONLY : nf90_noerr, NF90_STRERROR 747 IMPLICIT NONE 748 !------------------------------------------------------------------------------- 749 ! Arguments: 750 INTEGER, INTENT(IN) :: ncres 751 CHARACTER(LEN = *), INTENT(IN) :: fnam 752 !------------------------------------------------------------------------------- 753 IF(ncres/=nf90_noerr) THEN 754 WRITE(lunout, *)'Problem with file ' // TRIM(fnam) // ' in routine limit_netcdf.' 755 CALL abort_physic('limit_netcdf', NF90_STRERROR(ncres), 1) 756 END IF 757 758 END SUBROUTINE ncerr 759 760 !------------------------------------------------------------------------------- 761 762 763 !------------------------------------------------------------------------------- 764 765 SUBROUTINE strclean(s) 766 767 !------------------------------------------------------------------------------- 768 IMPLICIT NONE 769 !------------------------------------------------------------------------------- 770 ! Purpose: Remove tail null characters from the input string. 771 !------------------------------------------------------------------------------- 772 ! Parameters: 773 CHARACTER(LEN = *), INTENT(INOUT) :: s 774 !------------------------------------------------------------------------------- 775 ! Local variable: 776 INTEGER :: k 777 !------------------------------------------------------------------------------- 778 k = LEN_TRIM(s); DO WHILE(ICHAR(s(k:k))==0); s(k:k) = ' '; k = LEN_TRIM(s); 779 END DO 780 781 END SUBROUTINE strclean 782 783 !------------------------------------------------------------------------------- 784 785 786 !------------------------------------------------------------------------------- 787 788 FUNCTION is_in(s1, s2) RESULT(res) 789 790 !------------------------------------------------------------------------------- 791 IMPLICIT NONE 792 !------------------------------------------------------------------------------- 793 ! Purpose: Check wether s1 is present in the s2(:) list (case insensitive). 794 !------------------------------------------------------------------------------- 795 ! Arguments: 796 CHARACTER(LEN = *), INTENT(IN) :: s1, s2(:) 797 LOGICAL :: res 798 !------------------------------------------------------------------------------- 799 res = .FALSE.; DO k = 1, SIZE(s2); res = res.OR.strLow(s1)==strLow(s2(k)); 800 END DO 801 802 END FUNCTION is_in 803 804 !------------------------------------------------------------------------------- 805 806 807 !------------------------------------------------------------------------------- 808 809 ELEMENTAL FUNCTION strLow(s) RESULT(res) 810 811 !------------------------------------------------------------------------------- 812 IMPLICIT NONE 813 !------------------------------------------------------------------------------- 814 ! Purpose: Lower case conversion. 815 !------------------------------------------------------------------------------- 816 ! Arguments: 817 CHARACTER(LEN = *), INTENT(IN) :: s 818 CHARACTER(LEN = ns) :: res 819 !------------------------------------------------------------------------------- 820 ! Local variable: 821 INTEGER :: k, ix 822 !------------------------------------------------------------------------------- 823 res = s 824 DO k = 1, LEN(s); ix = IACHAR(s(k:k)) 825 IF(64<ix.AND.ix<91) res(k:k) = ACHAR(ix + 32) 826 END DO 827 828 END FUNCTION strLow 829 830 !------------------------------------------------------------------------------- 831 832 END SUBROUTINE limit_netcdf 833 834 END MODULE limit 835 3 836 !******************************************************************************* 4 ! Author : L. Fairhead, 27/01/94 5 !------------------------------------------------------------------------------- 6 ! Purpose: Boundary conditions files building for new model using climatologies. 7 ! Both grids have to be regular. 8 !------------------------------------------------------------------------------- 9 ! Note: This routine is designed to work for Earth 10 !------------------------------------------------------------------------------- 11 ! Modification history: 12 ! * 23/03/1994: Z. X. Li 13 ! * 09/1999: L. Fairhead (netcdf reading in LMDZ.3.3) 14 ! * 07/2001: P. Le Van 15 ! * 11/2009: L. Guez (ozone day & night climatos, see etat0_netcdf.F90) 16 ! * 12/2009: D. Cugnet (f77->f90, calendars, files from coupled runs) 17 !------------------------------------------------------------------------------- 18 19 USE ioipsl, ONLY: flininfo, flinopen, flinget, flinclo 20 USE assert_eq_m, ONLY: assert_eq 21 USE cal_tools_m, ONLY: year_len, mid_month 22 USE conf_dat_m, ONLY: conf_dat2d, conf_dat3d 23 USE dimphy, ONLY: klon, zmasq 24 USE geometry_mod, ONLY: longitude_deg, latitude_deg 25 USE phys_state_var_mod, ONLY: pctsrf 26 USE control_mod, ONLY: anneeref 27 USE init_ssrf_m, ONLY: start_init_subsurf 28 29 INTEGER, PARAMETER :: ns=256 30 CHARACTER(LEN=ns), PARAMETER :: & 31 fsst(5)=['amipbc_sst_1x1.nc ','amip_sst_1x1.nc ','cpl_atm_sst.nc '& 32 ,'histmth_sst.nc ','sstk.nc '], & 33 fsic(5)=['amipbc_sic_1x1.nc ','amip_sic_1x1.nc ','cpl_atm_sic.nc '& 34 ,'histmth_sic.nc ','ci.nc '], & 35 vsst(5)=['tosbcs ','tos ','SISUTESW ','tsol_oce ','sstk '], & 36 vsic(5)=['sicbcs ','sic ','SIICECOV ','pourc_sic ','ci '], & 37 frugo='Rugos.nc ', falbe='Albedo.nc ', frelf='Relief.nc ', & 38 vrug='RUGOS ', valb='ALBEDO ', vrel='RELIEF ', & 39 DegK(11)=['degK ','degree_K ','degreeK ','deg_K '& 40 ,'degsK ','degrees_K ','degreesK ','degs_K '& 41 ,'degree_kelvin ','degrees_kelvin','K '], & 42 DegC(10)=['degC ','degree_C ','degreeC ','deg_C '& 43 ,'degsC ','degrees_C ','degreesC ','degs_C '& 44 ,'degree_Celsius','celsius '], & 45 Perc(2) =['% ','percent '], & 46 Frac(2) =['1.0 ','1 '] 47 48 CONTAINS 49 50 !------------------------------------------------------------------------------- 51 52 SUBROUTINE limit_netcdf(masque, phis, extrap) 53 54 !------------------------------------------------------------------------------- 55 ! Author : L. Fairhead, 27/01/94 56 !------------------------------------------------------------------------------- 57 ! Purpose: Boundary conditions files building for new model using climatologies. 58 ! Both grids have to be regular. 59 !------------------------------------------------------------------------------- 60 ! Note: This routine is designed to work for Earth 61 !------------------------------------------------------------------------------- 62 ! Modification history: 63 ! * 23/03/1994: Z. X. Li 64 ! * 09/1999: L. Fairhead (netcdf reading in LMDZ.3.3) 65 ! * 07/2001: P. Le Van 66 ! * 11/2009: L. Guez (ozone day & night climatos, see etat0_netcdf.F90) 67 ! * 12/2009: D. Cugnet (f77->f90, calendars, files from coupled runs) 68 ! * 04/2016: D. Cugnet (12/14 recs SST/SIC files: cyclic/interannual runs) 69 ! * 05/2017: D. Cugnet (linear time interpolation for BCS files) 70 !------------------------------------------------------------------------------- 71 #ifndef CPP_1D 72 USE indice_sol_mod 73 USE netcdf, ONLY: NF90_OPEN, NF90_CREATE, NF90_CLOSE, & 74 NF90_DEF_DIM, NF90_DEF_VAR, NF90_PUT_VAR, NF90_PUT_ATT, & 75 nf90_noerr, NF90_NOWRITE, NF90_GLOBAL, & 76 NF90_CLOBBER, NF90_ENDDEF, NF90_UNLIMITED, NF90_FLOAT, & 77 NF90_64BIT_OFFSET 78 USE lmdz_cppkeys_wrapper, ONLY: nf90_format 79 USE inter_barxy_m, ONLY: inter_barxy 80 USE netcdf95, ONLY: nf95_def_var, nf95_put_att, nf95_put_var 81 USE comconst_mod, ONLY: pi 82 USE phys_cal_mod, ONLY: calend 83 IMPLICIT NONE 84 !------------------------------------------------------------------------------- 85 ! Arguments: 86 include "iniprint.h" 87 include "dimensions.h" 88 include "paramet.h" 89 REAL, DIMENSION(iip1,jjp1), INTENT(INOUT) :: masque ! land mask 90 REAL, DIMENSION(iip1,jjp1), INTENT(INOUT) :: phis ! ground geopotential 91 LOGICAL, INTENT(IN) :: extrap ! SST extrapolation flag 92 !------------------------------------------------------------------------------- 93 ! Local variables: 94 include "comgeom2.h" 95 96 !--- INPUT NETCDF FILES AND VARIABLES NAMES ------------------------------------ 97 CHARACTER(LEN=ns) :: icefile, sstfile, fnam, varname 98 99 !--- OUTPUT VARIABLES FOR NETCDF FILE ------------------------------------------ 100 REAL :: fi_ice(klon) 101 REAL, POINTER :: phy_rug(:,:)=>NULL(), phy_ice(:,:)=>NULL() 102 REAL, POINTER :: phy_sst(:,:)=>NULL(), phy_alb(:,:)=>NULL() 103 REAL, ALLOCATABLE :: phy_bil(:,:), pctsrf_t(:,:,:) 104 INTEGER :: nbad 105 106 !--- VARIABLES FOR OUTPUT FILE WRITING ----------------------------------------- 107 INTEGER :: nid, ndim, ntim, k, dims(2), ix_sic, ix_sst 108 INTEGER :: id_tim, id_SST, id_BILS, id_RUG, id_ALB 109 INTEGER :: id_FOCE, id_FSIC, id_FTER, id_FLIC, varid_longitude, varid_latitude 110 INTEGER :: ndays !--- Depending on the output calendar 111 CHARACTER(LEN=ns) :: str 112 113 !--- INITIALIZATIONS ----------------------------------------------------------- 114 CALL inigeom 115 116 !--- MASK, GROUND GEOPOT. & SUBSURFACES COMPUTATION (IN CASE ok_etat0==.FALSE.) 117 IF(ALL(masque==-99999.)) THEN 118 CALL start_init_orog0(rlonv,rlatu,phis,masque) 119 CALL gr_dyn_fi(1,iip1,jjp1,klon,masque,zmasq) !--- To physical grid 120 ALLOCATE(pctsrf(klon,nbsrf)) 121 CALL start_init_subsurf(.FALSE.) 122 !--- TO MATCH EXACTLY WHAT WOULD BE DONE IN etat0phys_netcdf 123 WHERE( masque(:,:)<EPSFRA) masque(:,:)=0. 124 WHERE(1.-masque(:,:)<EPSFRA) masque(:,:)=1. 125 END IF 126 127 !--- Beware: anneeref (from gcm.def) is used to determine output time sampling 128 ndays=year_len(anneeref) 129 130 !--- RUGOSITY TREATMENT -------------------------------------------------------- 131 CALL msg(0,""); CALL msg(0," *** TRAITEMENT DE LA RUGOSITE ***") 132 CALL get_2Dfield(frugo,vrug,'RUG',ndays,phy_rug,mask=masque(1:iim,:)) 133 134 !--- OCEAN TREATMENT ----------------------------------------------------------- 135 CALL msg(0,""); CALL msg(0," *** TRAITEMENT DE LA GLACE OCEANIQUE ***") 136 137 ! Input SIC file selection 138 ! Open file only to test if available 139 DO ix_sic=1,SIZE(fsic) 140 IF ( NF90_OPEN(TRIM(fsic(ix_sic)),NF90_NOWRITE,nid)==nf90_noerr ) THEN 141 icefile=fsic(ix_sic); varname=vsic(ix_sic); EXIT 142 END IF 143 END DO 144 IF(ix_sic==SIZE(fsic)+1) THEN 145 WRITE(lunout,*) 'ERROR! No sea-ice input file was found.' 146 WRITE(lunout,*) 'One of following files must be available : ' 147 DO k=1,SIZE(fsic); WRITE(lunout,*) TRIM(fsic(k)); END DO 148 CALL abort_physic('limit_netcdf','No sea-ice file was found',1) 149 END IF 150 CALL ncerr(NF90_CLOSE(nid),icefile) 151 CALL msg(0,'Fichier choisi pour la glace de mer:'//TRIM(icefile)) 152 153 CALL get_2Dfield(icefile,varname, 'SIC',ndays,phy_ice) 154 155 ALLOCATE(pctsrf_t(klon,nbsrf,ndays)) 156 DO k=1,ndays 157 fi_ice=phy_ice(:,k) 158 WHERE(fi_ice>=1.0 ) fi_ice=1.0 159 WHERE(fi_ice<EPSFRA) fi_ice=0.0 160 pctsrf_t(:,is_ter,k)=pctsrf(:,is_ter) ! land soil 161 pctsrf_t(:,is_lic,k)=pctsrf(:,is_lic) ! land ice 162 SELECT CASE(ix_sic) 163 CASE(3) ! SIC=pICE*(1-LIC-TER) (CPL) 164 pctsrf_t(:,is_sic,k)=fi_ice(:)*(1.-pctsrf(:,is_lic)-pctsrf(:,is_ter)) 165 CASE(4) ! SIC=pICE (HIST) 166 pctsrf_t(:,is_sic,k)=fi_ice(:) 167 CASE DEFAULT ! SIC=pICE-LIC (AMIP,ERAI) 168 pctsrf_t(:,is_sic,k)=fi_ice-pctsrf_t(:,is_lic,k) 169 END SELECT 170 WHERE(pctsrf_t(:,is_sic,k)<=0) pctsrf_t(:,is_sic,k)=0. 171 WHERE(1.0-zmasq<EPSFRA) 172 pctsrf_t(:,is_sic,k)=0.0 173 pctsrf_t(:,is_oce,k)=0.0 174 ELSEWHERE 175 WHERE(pctsrf_t(:,is_sic,k)>=1.0-zmasq) 176 pctsrf_t(:,is_sic,k)=1.0-zmasq 177 pctsrf_t(:,is_oce,k)=0.0 178 ELSEWHERE 179 pctsrf_t(:,is_oce,k)=1.0-zmasq-pctsrf_t(:,is_sic,k) 180 WHERE(pctsrf_t(:,is_oce,k)<EPSFRA) 181 pctsrf_t(:,is_oce,k)=0.0 182 pctsrf_t(:,is_sic,k)=1.0-zmasq 183 END WHERE 184 END WHERE 185 END WHERE 186 nbad=COUNT(pctsrf_t(:,is_oce,k)<0.0) 187 IF(nbad>0) WRITE(lunout,*) 'pb sous maille pour nb points = ',nbad 188 nbad=COUNT(ABS(SUM(pctsrf_t(:,:,k),DIM=2)-1.0)>EPSFRA) 189 IF(nbad>0) WRITE(lunout,*) 'pb sous surface pour nb points = ',nbad 190 END DO 191 DEALLOCATE(phy_ice) 192 193 !--- SST TREATMENT ------------------------------------------------------------- 194 CALL msg(0,""); CALL msg(0," *** TRAITEMENT DE LA SST ***") 195 196 ! Input SST file selection 197 ! Open file only to test if available 198 DO ix_sst=1,SIZE(fsst) 199 IF ( NF90_OPEN(TRIM(fsst(ix_sst)),NF90_NOWRITE,nid)==nf90_noerr ) THEN 200 sstfile=fsst(ix_sst); varname=vsst(ix_sst); EXIT 201 END IF 202 END DO 203 IF(ix_sst==SIZE(fsst)+1) THEN 204 WRITE(lunout,*) 'ERROR! No sst input file was found.' 205 WRITE(lunout,*) 'One of following files must be available : ' 206 DO k=1,SIZE(fsst); WRITE(lunout,*) TRIM(fsst(k)); END DO 207 CALL abort_physic('limit_netcdf','No sst file was found',1) 208 END IF 209 CALL ncerr(NF90_CLOSE(nid),sstfile) 210 CALL msg(0,'Fichier choisi pour la temperature de mer: '//TRIM(sstfile)) 211 212 CALL get_2Dfield(sstfile,varname,'SST',ndays,phy_sst,flag=extrap) 213 214 !--- ALBEDO TREATMENT ---------------------------------------------------------- 215 CALL msg(0,""); CALL msg(0," *** TRAITEMENT DE L'ALBEDO ***") 216 CALL get_2Dfield(falbe,valb,'ALB',ndays,phy_alb) 217 218 !--- REFERENCE GROUND HEAT FLUX TREATMENT -------------------------------------- 219 ALLOCATE(phy_bil(klon,ndays)); phy_bil=0.0 220 221 !--- OUTPUT FILE WRITING ------------------------------------------------------- 222 CALL msg(0,""); CALL msg(0,' *** Ecriture du fichier limit : debut ***') 223 fnam="limit.nc" 224 225 !--- File creation 226 CALL ncerr(NF90_CREATE(fnam,IOR(NF90_CLOBBER,NF90_64BIT_OFFSET),nid),fnam) 227 CALL ncerr(NF90_PUT_ATT(nid,NF90_GLOBAL,"title","Fichier conditions aux limites"),fnam) 228 str='File produced using ce0l executable.' 229 str=TRIM(str)//NEW_LINE(' ')//'Sea Ice Concentration built from' 230 SELECT CASE(ix_sic) 231 CASE(1); str=TRIM(str)//' Amip mid-month boundary condition (BCS).' 232 CASE(2); str=TRIM(str)//' Amip monthly mean observations.' 233 CASE(3); str=TRIM(str)//' IPSL coupled model outputs.' 234 CASE(4); str=TRIM(str)//' LMDZ model outputs.' 235 CASE(5); str=TRIM(str)//' ci.nc file.' 236 END SELECT 237 str=TRIM(str)//NEW_LINE(' ')//'Sea Surface Temperature built from' 238 SELECT CASE(ix_sst) 239 CASE(1); str=TRIM(str)//' Amip mid-month boundary condition (BCS).' 240 CASE(2); str=TRIM(str)//' Amip monthly mean observations.' 241 CASE(3); str=TRIM(str)//' IPSL coupled model outputs.' 242 CASE(4); str=TRIM(str)//' LMDZ model outputs.' 243 CASE(5); str=TRIM(str)//' sstk.nc file.' 244 END SELECT 245 CALL ncerr(NF90_PUT_ATT(nid,NF90_GLOBAL,"history",TRIM(str)),fnam) 246 247 !--- Dimensions creation 248 CALL ncerr(NF90_DEF_DIM(nid,"points_physiques",klon,ndim),fnam) 249 CALL ncerr(NF90_DEF_DIM(nid,"time",NF90_UNLIMITED,ntim),fnam) 250 251 dims=[ndim,ntim] 252 253 !--- Variables creation 254 CALL ncerr(NF90_DEF_VAR(nid,"TEMPS",nf90_format,[ntim],id_tim),fnam) 255 CALL ncerr(NF90_DEF_VAR(nid,"FOCE", nf90_format,dims,id_FOCE),fnam) 256 CALL ncerr(NF90_DEF_VAR(nid,"FSIC", nf90_format,dims,id_FSIC),fnam) 257 CALL ncerr(NF90_DEF_VAR(nid,"FTER", nf90_format,dims,id_FTER),fnam) 258 CALL ncerr(NF90_DEF_VAR(nid,"FLIC", nf90_format,dims,id_FLIC),fnam) 259 CALL ncerr(NF90_DEF_VAR(nid,"SST", nf90_format,dims,id_SST),fnam) 260 CALL ncerr(NF90_DEF_VAR(nid,"BILS", nf90_format,dims,id_BILS),fnam) 261 CALL ncerr(NF90_DEF_VAR(nid,"ALB", nf90_format,dims,id_ALB),fnam) 262 CALL ncerr(NF90_DEF_VAR(nid,"RUG", nf90_format,dims,id_RUG),fnam) 263 call nf95_def_var(nid, "longitude", NF90_FLOAT, ndim, varid_longitude) 264 call nf95_def_var(nid, "latitude", NF90_FLOAT, ndim, varid_latitude) 265 266 !--- Attributes creation 267 CALL ncerr(NF90_PUT_ATT(nid,id_tim, "title","Jour dans l annee"),fnam) 268 CALL ncerr(NF90_PUT_ATT(nid,id_tim, "calendar",calend),fnam) 269 CALL ncerr(NF90_PUT_ATT(nid,id_FOCE,"title","Fraction ocean"),fnam) 270 CALL ncerr(NF90_PUT_ATT(nid,id_FSIC,"title","Fraction glace de mer"),fnam) 271 CALL ncerr(NF90_PUT_ATT(nid,id_FTER,"title","Fraction terre"),fnam) 272 CALL ncerr(NF90_PUT_ATT(nid,id_FLIC,"title","Fraction land ice"),fnam) 273 CALL ncerr(NF90_PUT_ATT(nid,id_SST ,"title","Temperature superficielle de la mer"),fnam) 274 CALL ncerr(NF90_PUT_ATT(nid,id_BILS,"title","Reference flux de chaleur au sol"),fnam) 275 CALL ncerr(NF90_PUT_ATT(nid,id_ALB, "title","Albedo a la surface"),fnam) 276 CALL ncerr(NF90_PUT_ATT(nid,id_RUG, "title","Rugosite"),fnam) 277 278 call nf95_put_att(nid, varid_longitude, "standard_name", "longitude") 279 call nf95_put_att(nid, varid_longitude, "units", "degrees_east") 280 281 call nf95_put_att(nid, varid_latitude, "standard_name", "latitude") 282 call nf95_put_att(nid, varid_latitude, "units", "degrees_north") 283 284 CALL ncerr(NF90_ENDDEF(nid),fnam) 285 286 !--- Variables saving 287 CALL ncerr(NF90_PUT_VAR(nid,id_tim,[(REAL(k),k=1,ndays)]),fnam) 288 CALL ncerr(NF90_PUT_VAR(nid,id_FOCE,pctsrf_t(:,is_oce,:),[1,1],[klon,ndays]),fnam) 289 CALL ncerr(NF90_PUT_VAR(nid,id_FSIC,pctsrf_t(:,is_sic,:),[1,1],[klon,ndays]),fnam) 290 CALL ncerr(NF90_PUT_VAR(nid,id_FTER,pctsrf_t(:,is_ter,:),[1,1],[klon,ndays]),fnam) 291 CALL ncerr(NF90_PUT_VAR(nid,id_FLIC,pctsrf_t(:,is_lic,:),[1,1],[klon,ndays]),fnam) 292 CALL ncerr(NF90_PUT_VAR(nid,id_SST ,phy_sst(:,:),[1,1],[klon,ndays]),fnam) 293 CALL ncerr(NF90_PUT_VAR(nid,id_BILS,phy_bil(:,:),[1,1],[klon,ndays]),fnam) 294 CALL ncerr(NF90_PUT_VAR(nid,id_ALB ,phy_alb(:,:),[1,1],[klon,ndays]),fnam) 295 CALL ncerr(NF90_PUT_VAR(nid,id_RUG ,phy_rug(:,:),[1,1],[klon,ndays]),fnam) 296 call nf95_put_var(nid, varid_longitude, longitude_deg) 297 call nf95_put_var(nid, varid_latitude, latitude_deg) 298 299 CALL ncerr(NF90_CLOSE(nid),fnam) 300 301 CALL msg(0,""); CALL msg(0,' *** Ecriture du fichier limit : fin ***') 302 303 DEALLOCATE(pctsrf_t,phy_sst,phy_bil,phy_alb,phy_rug) 304 305 306 !=============================================================================== 307 308 CONTAINS 309 310 !=============================================================================== 311 312 313 !------------------------------------------------------------------------------- 314 315 SUBROUTINE get_2Dfield(fnam, varname, mode, ndays, champo, flag, mask) 316 317 !----------------------------------------------------------------------------- 318 ! Comments: 319 ! There are two assumptions concerning the NetCDF files, that are satisfied 320 ! with files that are conforming NC convention: 321 ! 1) The last dimension of the variables used is the time record. 322 ! 2) Dimensional variables have the same names as corresponding dimensions. 323 !----------------------------------------------------------------------------- 324 USE netcdf, ONLY: NF90_OPEN, NF90_INQ_VARID, NF90_INQUIRE_VARIABLE, & 325 NF90_CLOSE, nf90_inq_dimid, nf90_inquire_dimension, nf90_get_var, & 326 NF90_GET_ATT 327 USE pchsp_95_m, only: pchsp_95 328 USE pchfe_95_m, only: pchfe_95 329 USE arth_m, only: arth 330 USE indice_sol_mod 331 332 IMPLICIT NONE 333 include "dimensions.h" 334 include "paramet.h" 335 include "comgeom2.h" 336 !----------------------------------------------------------------------------- 337 ! Arguments: 338 CHARACTER(LEN=*), INTENT(IN) :: fnam ! NetCDF file name 339 CHARACTER(LEN=*), INTENT(IN) :: varname ! NetCDF variable name 340 CHARACTER(LEN=*), INTENT(IN) :: mode ! RUG, SIC, SST or ALB 341 INTEGER, INTENT(IN) :: ndays ! current year number of days 342 REAL, POINTER, DIMENSION(:, :) :: champo ! output field = f(t) 343 LOGICAL, OPTIONAL, INTENT(IN) :: flag ! extrapol. (SST) old ice (SIC) 344 REAL, OPTIONAL, DIMENSION(iim, jjp1), INTENT(IN) :: mask 345 !------------------------------------------------------------------------------ 346 ! Local variables: 347 !--- NetCDF 348 INTEGER :: ncid, varid ! NetCDF identifiers 349 CHARACTER(LEN=ns) :: dnam ! dimension name 350 !--- dimensions 351 INTEGER :: dids(4) ! NetCDF dimensions identifiers 352 REAL, ALLOCATABLE :: dlon_ini(:) ! initial longitudes vector 353 REAL, ALLOCATABLE :: dlat_ini(:) ! initial latitudes vector 354 REAL, POINTER :: dlon(:), dlat(:) ! reordered lon/lat vectors 355 !--- fields 356 INTEGER :: imdep, jmdep, lmdep ! dimensions of 'champ' 357 REAL, ALLOCATABLE :: champ(:,:) ! wanted field on initial grid 358 REAL, ALLOCATABLE :: yder(:), timeyear(:) 359 REAL :: champint(iim,jjp1) ! interpolated field 360 REAL, ALLOCATABLE :: champtime(:,:,:) 361 REAL, ALLOCATABLE :: champan(:,:,:) 362 !--- input files 363 CHARACTER(LEN=ns) :: fnam_m, fnam_p ! previous/next files names 364 CHARACTER(LEN=ns) :: cal_in ! calendar 365 CHARACTER(LEN=ns) :: units ! attribute "units" in sic/sst file 366 INTEGER :: ndays_in ! number of days 367 REAL :: value ! mean/max value near equator 368 !--- misc 369 INTEGER :: i, j, k, l ! loop counters 370 REAL, ALLOCATABLE :: work(:,:) ! used for extrapolation 371 CHARACTER(LEN=ns) :: title, mess ! for messages 372 LOGICAL :: is_bcs ! flag for BCS data 373 LOGICAL :: extrp ! flag for extrapolation 374 LOGICAL :: ll 375 REAL :: chmin, chmax, timeday, al 376 INTEGER ierr, idx 377 integer n_extrap ! number of extrapolated points 378 logical skip 379 380 !------------------------------------------------------------------------------ 381 !---Variables depending on keyword 'mode' ------------------------------------- 382 NULLIFY(champo) 383 384 SELECT CASE(mode) 385 CASE('RUG'); title='Rugosite' 386 CASE('SIC'); title='Sea-ice' 387 CASE('SST'); title='SST' 388 CASE('ALB'); title='Albedo' 389 END SELECT 390 extrp=.FALSE.; IF(PRESENT(flag).AND.mode=='SST') extrp=flag 391 is_bcs=(mode=='SIC'.AND.ix_sic==1).OR.(mode=='SST'.AND.ix_sst==1) 392 idx=INDEX(fnam,'.nc')-1 393 394 !--- GETTING SOME DIMENSIONAL VARIABLES FROM FILE ----------------------------- 395 CALL msg(5,' Now reading file : '//TRIM(fnam)) 396 CALL ncerr(NF90_OPEN(fnam, NF90_NOWRITE, ncid),fnam) 397 CALL ncerr(NF90_INQ_VARID(ncid, trim(varname), varid),fnam) 398 CALL ncerr(NF90_INQUIRE_VARIABLE(ncid, varid, dimids=dids),fnam) 399 400 !--- Longitude 401 CALL ncerr(nf90_inquire_dimension(ncid, dids(1), name=dnam, len=imdep),fnam) 402 ALLOCATE(dlon_ini(imdep), dlon(imdep)) 403 CALL ncerr(NF90_INQ_VARID(ncid, dnam, varid), fnam) 404 CALL ncerr(nf90_get_var(ncid, varid, dlon_ini), fnam) 405 CALL msg(5,'variable '//TRIM(dnam)//' dimension ', imdep) 406 407 !--- Latitude 408 CALL ncerr(nf90_inquire_dimension(ncid, dids(2), name=dnam, len=jmdep),fnam) 409 ALLOCATE(dlat_ini(jmdep), dlat(jmdep)) 410 CALL ncerr(NF90_INQ_VARID(ncid, dnam, varid), fnam) 411 CALL ncerr(nf90_get_var(ncid, varid, dlat_ini), fnam) 412 CALL msg(5,'variable '//TRIM(dnam)//' dimension ', jmdep) 413 414 !--- Time (variable is not needed - it is rebuilt - but calendar is) 415 CALL ncerr(nf90_inquire_dimension(ncid, dids(3), name=dnam, len=lmdep), fnam) 416 ALLOCATE(timeyear(lmdep+2)) 417 CALL ncerr(NF90_INQ_VARID(ncid, dnam, varid), fnam) 418 cal_in=' ' 419 IF(NF90_GET_ATT(ncid, varid, 'calendar', cal_in)/=nf90_noerr) THEN 420 SELECT CASE(mode) 421 CASE('RUG', 'ALB'); cal_in='360_day' 422 CASE('SIC', 'SST'); cal_in='gregorian' 423 END SELECT 424 CALL msg(0,'WARNING: missing "calendar" attribute for "time" in '& 425 //TRIM(fnam)//'. Choosing default value.') 426 END IF 427 CALL strclean(cal_in) !--- REMOVE (WEIRD) NULL CHARACTERS 428 CALL msg(0,'var, calendar, dim: '//TRIM(dnam)//' '//TRIM(cal_in), lmdep) 429 430 !--- Determining input file number of days, depending on calendar 431 ndays_in=year_len(anneeref, cal_in) 432 433 !--- Rebuilding input time vector (field from input file might be unreliable) 434 IF(lmdep==12) THEN 435 timeyear=mid_month(anneeref, cal_in) 436 CALL msg(0,'Monthly input file(s) for '//TRIM(title)//'.') 437 ELSE IF(lmdep==ndays_in) THEN 438 timeyear=[(REAL(k)-0.5,k=0,ndays_in+1)] 439 CALL msg(0,'Daily input file (no time interpolation).') 440 ELSE 441 WRITE(mess,'(a,i3,a,i3,a)')'Mismatching input file: found',lmdep, & 442 ' records, 12/',ndays_in,' (monthly/daily needed).' 443 CALL abort_physic('mid_month',TRIM(mess),1) 444 END IF 445 446 !--- GETTING THE FIELD AND INTERPOLATING IT ---------------------------------- 447 ALLOCATE(champ(imdep, jmdep), champtime(iim, jjp1, lmdep+2)) 448 IF(extrp) ALLOCATE(work(imdep, jmdep)) 449 CALL msg(5,'') 450 CALL msg(5,'READ AND INTERPOLATE HORIZONTALLY ', lmdep, ' FIELDS.') 451 CALL ncerr(NF90_INQ_VARID(ncid, varname, varid), fnam) 452 DO l=1, lmdep 453 CALL ncerr(nf90_get_var(ncid,varid,champ,[1,1,l],[imdep,jmdep,1]),fnam) 454 CALL conf_dat2d(title, dlon_ini, dlat_ini, dlon, dlat, champ, .TRUE.) 455 456 !--- FOR SIC/SST FIELDS ONLY 457 IF(l==1.AND.is_in(mode,['SIC','SST'])) THEN 458 459 !--- DETERMINE THE UNIT: READ FROM FILE OR ASSUMED USING FIELD VALUES 460 ierr=NF90_GET_ATT(ncid, varid, 'units', units) 461 IF(ierr==nf90_noerr) THEN !--- ATTRIBUTE "units" FOUND IN THE FILE 462 CALL strclean(units) 463 IF(mode=='SIC'.AND.is_in(units,Perc)) units="%" 464 IF(mode=='SIC'.AND.is_in(units,Frac)) units="1" 465 IF(mode=='SST'.AND.is_in(units,DegC)) units="C" 466 IF(mode=='SST'.AND.is_in(units,DegK)) units="K" 467 ELSE !--- CHECK THE FIELD VALUES 468 IF(mode=='SIC') value=MAXVAL(champ(:,:)) 469 IF(mode=='SST') value= SUM(champ(:,jmdep/2),DIM=1)/REAL(imdep) 470 IF(mode=='SIC') THEN; units="1"; IF(value>= 10.) units="%"; END IF 471 IF(mode=='SST') THEN; units="C"; IF(value>=100.) units="K"; END IF 472 END IF 473 CALL msg(0,'INPUT FILE '//TRIM(title)//' UNIT IS: "'//TRIM(units)//'".') 474 IF(ierr/=nf90_noerr) CALL msg(0,'WARNING ! UNIT TO BE CHECKED ! ' & 475 //'No "units" attribute, so only based on the fields values.') 476 477 !--- CHECK VALUES ARE IN THE EXPECTED RANGE 478 SELECT CASE(units) 479 CASE('%'); ll=ANY(champ>100.0+EPSFRA); str='percentages > 100.' 480 CASE('1'); ll=ANY(champ> 1.0+EPSFRA); str='fractions > 1.' 481 CASE('C'); ll=ANY(champ<-100.).OR.ANY(champ> 60.); str='<-100 or >60 DegC' 482 CASE('K'); ll=ANY(champ< 180.).OR.ANY(champ>330.); str='<180 or >330 DegK' 483 CASE DEFAULT; CALL abort_physic(mode, 'Unrecognized '//TRIM(title) & 484 //' unit: '//TRIM(units),1) 485 END SELECT 486 487 !--- DROPPED FOR BCS DATA (FRACTIONS CAN BE HIGHER THAN 1) 488 IF(ll.AND.ix_sic/=1.AND.mode=='SIC') & 489 CALL abort_physic(mode,'unrealistic '//TRIM(mode)//' found: '//TRIM(str), 1) 490 491 END IF 492 493 IF(extrp) CALL extrapol(champ,imdep,jmdep,999999.,.TRUE.,.TRUE.,2,work) 494 IF(l==1) THEN 495 CALL msg(5,"--------------------------------------------------------") 496 CALL msg(5,"$$$ Barycentric interpolation for "//TRIM(title)//" $$$") 497 CALL msg(5,"--------------------------------------------------------") 498 END IF 499 IF(mode=='RUG') champ=LOG(champ) 500 CALL inter_barxy(dlon,dlat(:jmdep-1),champ,rlonu(:iim),rlatv,champint) 501 IF(mode=='RUG') THEN 502 champint=EXP(champint) 503 WHERE(NINT(mask)/=1) champint=0.001 504 END IF 505 champtime(:, :, l+1)=champint 506 END DO 507 CALL ncerr(NF90_CLOSE(ncid), fnam) 508 509 !--- FIRST RECORD: LAST ONE OF PREVIOUS YEAR (CURRENT YEAR IF UNAVAILABLE) 510 fnam_m=fnam(1:idx)//'_m.nc' 511 IF(NF90_OPEN(fnam_m,NF90_NOWRITE,ncid)==nf90_noerr) THEN 512 CALL msg(0,'Reading previous year file ("'//TRIM(fnam_m)//'") last record for '//TRIM(title)) 513 CALL ncerr(NF90_INQ_VARID(ncid, varname, varid),fnam_m) 514 CALL ncerr(NF90_INQUIRE_VARIABLE(ncid, varid, dimids=dids),fnam_m) 515 CALL ncerr(nf90_inquire_dimension(ncid, dids(3), len=l), fnam_m) 516 CALL ncerr(nf90_get_var(ncid,varid,champ,[1,1,l],[imdep,jmdep,1]),fnam_m) 517 CALL ncerr(NF90_CLOSE(ncid), fnam_m) 518 CALL conf_dat2d(title, dlon_ini, dlat_ini, dlon, dlat, champ, .TRUE.) 519 IF(extrp) CALL extrapol(champ,imdep,jmdep,999999.,.TRUE.,.TRUE.,2,work) 520 IF(mode=='RUG') champ=LOG(champ) 521 CALL inter_barxy(dlon,dlat(:jmdep-1),champ,rlonu(:iim),rlatv,champint) 522 IF(mode=='RUG') THEN 523 champint=EXP(champint) 524 WHERE(NINT(mask)/=1) champint=0.001 525 END IF 526 champtime(:, :, 1)=champint 527 ELSE 528 CALL msg(0,'Using current year file ("'//TRIM(fnam)//'") last record for '//TRIM(title)) 529 champtime(:, :, 1)=champtime(:, :, lmdep+1) 530 END IF 531 532 !--- LAST RECORD: FIRST ONE OF NEXT YEAR (CURRENT YEAR IF UNAVAILABLE) 533 fnam_p=fnam(1:idx)//'_p.nc' 534 IF(NF90_OPEN(fnam_p,NF90_NOWRITE,ncid)==nf90_noerr) THEN 535 CALL msg(0,'Reading next year file ("'//TRIM(fnam_p)//'") first record for '//TRIM(title)) 536 CALL ncerr(NF90_INQ_VARID(ncid, varname, varid),fnam_p) 537 CALL ncerr(nf90_get_var(ncid,varid,champ,[1,1,1],[imdep,jmdep,1]),fnam_p) 538 CALL ncerr(NF90_CLOSE(ncid), fnam_p) 539 CALL conf_dat2d(title, dlon_ini, dlat_ini, dlon, dlat, champ, .TRUE.) 540 IF(extrp) CALL extrapol(champ,imdep,jmdep,999999.,.TRUE.,.TRUE.,2,work) 541 IF(mode=='RUG') champ=LOG(champ) 542 CALL inter_barxy(dlon,dlat(:jmdep-1),champ,rlonu(:iim),rlatv,champint) 543 IF(mode=='RUG') THEN 544 champint=EXP(champint) 545 WHERE(NINT(mask)/=1) champint=0.001 546 END IF 547 champtime(:, :, lmdep+2)=champint 548 ELSE 549 CALL msg(0,'Using current year file ("'//TRIM(fnam)//'") first record for '//TRIM(title)) 550 champtime(:, :, lmdep+2)=champtime(:, :, 2) 551 END IF 552 DEALLOCATE(dlon_ini, dlat_ini, dlon, dlat, champ) 553 IF(extrp) DEALLOCATE(work) 554 555 !--- TIME INTERPOLATION ------------------------------------------------------ 556 IF(prt_level>0) THEN 557 IF(ndays/=ndays_in) THEN 558 WRITE(lunout,*)'DIFFERENT YEAR LENGTHS:' 559 WRITE(lunout,*)' In the input file: ',ndays_in 560 WRITE(lunout,*)' In the output file: ',ndays 561 END IF 562 IF(lmdep==ndays_in) THEN 563 WRITE(lunout, *)'NO TIME INTERPOLATION.' 564 WRITE(lunout, *)' Daily input file.' 565 ELSE 566 IF( is_bcs) WRITE(lunout, *)'LINEAR TIME INTERPOLATION.' 567 IF(.NOT.is_bcs) WRITE(lunout, *)'SPLINES TIME INTERPOLATION.' 568 WRITE(lunout, *)' Input time vector: ', timeyear 569 WRITE(lunout, *)' Output time vector: from 0.5 to ', ndays-0.5 570 END IF 571 END IF 572 ALLOCATE(champan(iip1, jjp1, ndays)) 573 574 IF(lmdep==ndays_in) THEN !--- DAILY DATA: NO TIME INTERPOLATION 575 DO l=1,lmdep 576 champan(1:iim,:,l)=champtime(:,:,l+1) 577 END DO 578 ELSE IF(is_bcs) THEN !--- BCS DATA: LINEAR TIME INTERPOLATION 579 l=1 580 DO k=1, ndays 581 timeday = (REAL(k)-0.5)*REAL(ndays_in)/ndays 582 IF(timeyear(l+1)<timeday) l=l+1 583 al=(timeday-timeyear(l))/(timeyear(l+1)-timeyear(l)) 584 champan(1:iim,:,k) = champtime(1:iim,:,l)+al*(champtime(1:iim,:,l+1)-champtime(1:iim,:,l)) 585 END DO 586 ELSE !--- AVE DATA: SPLINE TIME INTERPOLATION 587 skip = .false. 588 n_extrap = 0 589 ALLOCATE(yder(lmdep+2)) 590 DO j=1, jjp1 591 DO i=1, iim 592 yder = pchsp_95(timeyear, champtime(i, j, :), ibeg=2, iend=2, & 593 vc_beg=0., vc_end=0.) 594 CALL pchfe_95(timeyear, champtime(i, j, :), yder, skip, & 595 arth(0.5, real(ndays_in) / ndays, ndays), champan(i, j, :), ierr) 596 if (ierr < 0) call abort_physic("get_2Dfield", "", 1) 597 n_extrap = n_extrap + ierr 598 END DO 599 END DO 600 IF(n_extrap /= 0) WRITE(lunout,*) "get_2Dfield pchfe_95: n_extrap = ", n_extrap 601 DEALLOCATE(yder) 602 END IF 603 champan(iip1, :, :)=champan(1, :, :) 604 DEALLOCATE(champtime, timeyear) 605 606 !--- Checking the result 607 DO j=1, jjp1 608 CALL minmax(iip1, champan(1, j, 10), chmin, chmax) 609 IF (prt_level>5) WRITE(lunout, *)' ',TRIM(title),' at time 10 ', chmin, chmax, j 610 END DO 611 612 !--- SPECIAL FILTER FOR SST: SST>271.38 -------------------------------------- 613 IF(mode=='SST') THEN 614 SELECT CASE(units) 615 CASE("K"); CALL msg(0,'SST field is already in kelvins.') 616 CASE("C"); CALL msg(0,'SST field converted from celcius degrees to kelvins.') 617 champan(:, :, :)=champan(:, :, :)+273.15 618 END SELECT 619 CALL msg(0,'Filtering SST: Sea Surface Temperature >= 271.38') 620 WHERE(champan<271.38) champan=271.38 621 END IF 622 623 !--- SPECIAL FILTER FOR SIC: 0.0<SIC<1.0 ------------------------------------- 624 IF(mode=='SIC') THEN 625 SELECT CASE(units) 626 CASE("1"); CALL msg(0,'SIC field already in fraction of 1') 627 CASE("%"); CALL msg(0,'SIC field converted from percentage to fraction of 1.') 628 champan(:, :, :)=champan(:, :, :)/100. 629 END SELECT 630 CALL msg(0,'Filtering SIC: 0.0 <= Sea-ice <=1.0') 631 WHERE(champan>1.0) champan=1.0 632 WHERE(champan<0.0) champan=0.0 633 END IF 634 635 !--- DYNAMICAL TO PHYSICAL GRID ---------------------------------------------- 636 ALLOCATE(champo(klon, ndays)) 637 DO k=1, ndays 638 CALL gr_dyn_fi(1, iip1, jjp1, klon, champan(1, 1, k), champo(:, k)) 639 END DO 640 DEALLOCATE(champan) 641 642 END SUBROUTINE get_2Dfield 643 644 !------------------------------------------------------------------------------- 645 646 647 !------------------------------------------------------------------------------- 648 649 SUBROUTINE start_init_orog0(lon_in,lat_in,phis,masque) 650 651 !------------------------------------------------------------------------------- 652 USE grid_noro_m, ONLY: grid_noro0 653 IMPLICIT NONE 654 !=============================================================================== 655 ! Purpose: Compute "phis" just like it would be in start_init_orog. 656 !=============================================================================== 657 ! Arguments: 658 REAL, INTENT(IN) :: lon_in(:), lat_in(:) ! dim (iml) (jml) 659 REAL, INTENT(INOUT) :: phis(:,:), masque(:,:) ! dim (iml,jml) 660 !------------------------------------------------------------------------------- 661 ! Local variables: 662 CHARACTER(LEN=ns) :: modname="start_init_orog0" 663 INTEGER :: fid, llm_tmp,ttm_tmp, iml,jml, iml_rel,jml_rel, itau(1) 664 REAL :: lev(1), date, dt, deg2rad 665 REAL, ALLOCATABLE :: lon_rad(:), lon_ini(:), lon_rel(:,:), relief_hi(:,:) 666 REAL, ALLOCATABLE :: lat_rad(:), lat_ini(:), lat_rel(:,:) 667 !------------------------------------------------------------------------------- 668 iml=assert_eq(SIZE(lon_in),SIZE(phis,1),SIZE(masque,1),TRIM(modname)//" iml") 669 jml=assert_eq(SIZE(lat_in),SIZE(phis,2),SIZE(masque,2),TRIM(modname)//" jml") 670 IF(iml/=iip1) CALL abort_gcm(TRIM(modname),'iml/=iip1',1) 671 IF(jml/=jjp1) CALL abort_gcm(TRIM(modname),'jml/=jjp1',1) 672 pi=2.0*ASIN(1.0); deg2rad=pi/180.0 673 IF(ANY(phis/=-99999.)) RETURN !--- phis ALREADY KNOWN 674 675 !--- HIGH RESOLUTION OROGRAPHY 676 CALL flininfo(frelf, iml_rel, jml_rel, llm_tmp, ttm_tmp, fid) 677 678 ALLOCATE(lat_rel(iml_rel,jml_rel),lon_rel(iml_rel,jml_rel)) 679 CALL flinopen(frelf, .FALSE., iml_rel, jml_rel, llm_tmp, lon_rel, lat_rel, & 680 lev, ttm_tmp, itau, date, dt, fid) 681 ALLOCATE(relief_hi(iml_rel,jml_rel)) 682 CALL flinget(fid, vrel, iml_rel, jml_rel, llm_tmp, ttm_tmp, 1, 1, relief_hi) 683 CALL flinclo(fid) 684 685 !--- IF ANGLES ARE IN DEGREES, THEY ARE CONVERTED INTO RADIANS 686 ALLOCATE(lon_ini(iml_rel),lat_ini(jml_rel)) 687 lon_ini(:)=lon_rel(:,1); IF(MAXVAL(lon_rel)>pi) lon_ini=lon_ini*deg2rad 688 lat_ini(:)=lat_rel(1,:); IF(MAXVAL(lat_rel)>pi) lat_ini=lat_ini*deg2rad 689 690 !--- FIELDS ARE PROCESSED TO BE ON STANDARD ANGULAR DOMAINS 691 ALLOCATE(lon_rad(iml_rel),lat_rad(jml_rel)) 692 CALL conf_dat2d(vrel, lon_ini, lat_ini, lon_rad, lat_rad, relief_hi, .FALSE.) 693 DEALLOCATE(lon_ini,lat_ini) 694 695 !--- COMPUTING SURFACE GEOPOTENTIAL USING ROUTINE grid_noro0 696 WRITE(lunout,*) 697 WRITE(lunout,*)'*** Compute surface geopotential ***' 698 699 !--- CALL OROGRAPHY MODULE (REDUCED VERSION) TO COMPUTE FIELDS 700 CALL grid_noro0(lon_rad, lat_rad, relief_hi, lon_in, lat_in, phis, masque) 701 phis = phis * 9.81 702 phis(iml,:) = phis(1,:) 703 DEALLOCATE(relief_hi,lon_rad,lat_rad) 704 705 END SUBROUTINE start_init_orog0 706 707 !------------------------------------------------------------------------------- 708 709 710 !------------------------------------------------------------------------------- 711 712 SUBROUTINE msg(lev,str1,i,str2) 713 714 !------------------------------------------------------------------------------- 715 ! Arguments: 716 INTEGER, INTENT(IN) :: lev 717 CHARACTER(LEN=*), INTENT(IN) :: str1 718 INTEGER, OPTIONAL, INTENT(IN) :: i 719 CHARACTER(LEN=*), OPTIONAL, INTENT(IN) :: str2 720 !------------------------------------------------------------------------------- 721 IF(prt_level>=lev) THEN 722 IF(PRESENT(str2)) THEN 723 WRITE(lunout,*) TRIM(str1), i, TRIM(str2) 724 ELSE IF(PRESENT(i)) THEN 725 WRITE(lunout,*) TRIM(str1), i 726 ELSE 727 WRITE(lunout,*) TRIM(str1) 728 END IF 729 END IF 730 731 END SUBROUTINE msg 732 733 !------------------------------------------------------------------------------- 734 735 736 !------------------------------------------------------------------------------- 737 738 SUBROUTINE ncerr(ncres,fnam) 739 740 !------------------------------------------------------------------------------- 741 ! Purpose: NetCDF errors handling. 742 !------------------------------------------------------------------------------- 743 USE netcdf, ONLY : nf90_noerr, NF90_STRERROR 744 IMPLICIT NONE 745 !------------------------------------------------------------------------------- 746 ! Arguments: 747 INTEGER, INTENT(IN) :: ncres 748 CHARACTER(LEN=*), INTENT(IN) :: fnam 749 !------------------------------------------------------------------------------- 750 IF(ncres/=nf90_noerr) THEN 751 WRITE(lunout,*)'Problem with file '//TRIM(fnam)//' in routine limit_netcdf.' 752 CALL abort_physic('limit_netcdf',NF90_STRERROR(ncres),1) 753 END IF 754 755 END SUBROUTINE ncerr 756 757 !------------------------------------------------------------------------------- 758 759 760 !------------------------------------------------------------------------------- 761 762 SUBROUTINE strclean(s) 763 764 !------------------------------------------------------------------------------- 765 IMPLICIT NONE 766 !------------------------------------------------------------------------------- 767 ! Purpose: Remove tail null characters from the input string. 768 !------------------------------------------------------------------------------- 769 ! Parameters: 770 CHARACTER(LEN=*), INTENT(INOUT) :: s 771 !------------------------------------------------------------------------------- 772 ! Local variable: 773 INTEGER :: k 774 !------------------------------------------------------------------------------- 775 k=LEN_TRIM(s); DO WHILE(ICHAR(s(k:k))==0); s(k:k)=' '; k=LEN_TRIM(s); END DO 776 777 END SUBROUTINE strclean 778 779 !------------------------------------------------------------------------------- 780 781 782 !------------------------------------------------------------------------------- 783 784 FUNCTION is_in(s1,s2) RESULT(res) 785 786 !------------------------------------------------------------------------------- 787 IMPLICIT NONE 788 !------------------------------------------------------------------------------- 789 ! Purpose: Check wether s1 is present in the s2(:) list (case insensitive). 790 !------------------------------------------------------------------------------- 791 ! Arguments: 792 CHARACTER(LEN=*), INTENT(IN) :: s1, s2(:) 793 LOGICAL :: res 794 !------------------------------------------------------------------------------- 795 res=.FALSE.; DO k=1,SIZE(s2); res=res.OR.strLow(s1)==strLow(s2(k)); END DO 796 797 END FUNCTION is_in 798 799 !------------------------------------------------------------------------------- 800 801 802 !------------------------------------------------------------------------------- 803 804 ELEMENTAL FUNCTION strLow(s) RESULT(res) 805 806 !------------------------------------------------------------------------------- 807 IMPLICIT NONE 808 !------------------------------------------------------------------------------- 809 ! Purpose: Lower case conversion. 810 !------------------------------------------------------------------------------- 811 ! Arguments: 812 CHARACTER(LEN=*), INTENT(IN) :: s 813 CHARACTER(LEN=ns) :: res 814 !------------------------------------------------------------------------------- 815 ! Local variable: 816 INTEGER :: k, ix 817 !------------------------------------------------------------------------------- 818 res=s 819 DO k=1,LEN(s); ix=IACHAR(s(k:k)) 820 IF(64<ix.AND.ix<91) res(k:k)=ACHAR(ix+32) 821 END DO 822 823 END FUNCTION strLow 824 825 !------------------------------------------------------------------------------- 826 827 #endif 828 ! of #ifndef CPP_1D 829 END SUBROUTINE limit_netcdf 830 831 END MODULE limit 832 833 !******************************************************************************* 834 837
Note: See TracChangeset
for help on using the changeset viewer.
