! -*- Mode: f90 -*- PROGRAM lmdgrille !!! !!! Ecris les grilles LMDZ au format OASIS !! !! Attention : !! pour OASIS 3, le nord est en haut (j=jpai) et le sud en bas (j=1) !! pour OASI3-MCT, il faut conserver le sens du modèle !! !! Unite IEEE : 32 !! USE declare USE dimensions USE mod_prih USE fliocom USE getincom USE errioipsl USE mod_inipar IMPLICIT none !! INTEGER (kind=il) :: jai, jaj, jk, jc REAL (kind=rl), DIMENSION (:, :), ALLOCATABLE :: xlont, xlatt, xlonu, xlatu, xlonv, xlatv, xlonf, xlatf REAL (kind=rl), DIMENSION (:, :), ALLOCATABLE :: xlat_o2a REAL (kind=rl), DIMENSION (:, :), ALLOCATABLE :: xmas, xare, xtmp, xfra REAL (kind=rl), DIMENSION (:, :, :), ALLOCATABLE :: xclon, xclat REAL (kind=rl), DIMENSION (:, :), ALLOCATABLE :: xaw, xae, yan, yas INTEGER (kind=il) :: id_restart, id_restartphy, id_hist, id_grid, id_coo, no = 32, narg = 0, iargc REAL (kind=rl) :: ra, zzz, xdel, ydel CHARACTER (len=180) :: cfname, clmd, c_restart, c_restartphy, c_hist CHARACTER (len=80) :: cmodel CHARACTER (len=8) :: clon, clat, cmsk, csrf, clarg REAL (kind=rl), DIMENSION (:, :), ALLOCATABLE :: flux_iceberg INTEGER :: ia_sk=1_il, ja_sk=1_il INTEGER :: nlmd REAL (kind=rl), DIMENSION (:), ALLOCATABLE :: tab1d !! INTEGER (kind=il) :: il_grid, il_mask, il_area, il_frac, id INTEGER (kind=il) :: ierr ! LOGICAL :: l_build = .FALSE. !< Construction a partir de rien LOGICAL :: l_restart = .FALSE. !< Construction a partir de restart et restartphy LOGICAL :: l_read = .FALSE. !< Construction a partir de grids.nc, areas.nc, masks.nc !LOGICAL :: lec_msk = .TRUE., lec_srf = .TRUE., lec_coo = .TRUE., lec_fra = .TRUE. LOGICAL :: la_retourn_lat !< Retournement nor/sud des champs lmdz CHARACTER (LEN=1) :: cl_rk !! CALL inipar !! IF ( jpai .GT. 100_il ) ia_sk = 2 IF ( jpai .GT. 200_il ) ia_sk = 5 IF ( jpai .GT. 500_il ) ia_sk = 10 IF ( jpaj .GT. 100_il ) ia_sk = 2 IF ( jpaj .GT. 200_il ) ja_sk = 5 IF ( jpaj .GT. 500_il ) ja_sk = 10 !! ALLOCATE (xlont (jpai , jpaj ), STAT=ierr) ; CALL chk_allo (ierr, 'xlont', lreset=.TRUE., crout='lmdz') ALLOCATE (xlatt (jpai , jpaj ), STAT=ierr) ; CALL chk_allo (ierr, 'xlatt') ALLOCATE (xlonu (jpaiu, jpaju), STAT=ierr) ; CALL chk_allo (ierr, 'xlonu') ALLOCATE (xlatu (jpaiu, jpaju), STAT=ierr) ; CALL chk_allo (ierr, 'xlatu') ALLOCATE (xlonv (jpaiv, jpajv), STAT=ierr) ; CALL chk_allo (ierr, 'xlonv') ALLOCATE (xlatv (jpaiv, jpajv), STAT=ierr) ; CALL chk_allo (ierr, 'xlatv') ALLOCATE (xmas (jpai , jpaj ), STAT=ierr) ; CALL chk_allo (ierr, 'xmas') ALLOCATE (xare (jpai , jpaj ), STAT=ierr) ; CALL chk_allo (ierr, 'xare') ALLOCATE (xfra (jpai , jpaj ), STAT=ierr) ; CALL chk_allo (ierr, 'xfra') ALLOCATE (xtmp (jpai , jpaj )) ; CALL chk_allo (ierr, 'xtmp') ALLOCATE (xclon (jpai , jpaj, 4), STAT=ierr) ; CALL chk_allo (ierr, 'xclon') ALLOCATE (xclat (jpai , jpaj, 4), STAT=ierr) ; CALL chk_allo (ierr, 'xclat') ALLOCATE (xaw (jpai , jpaj), STAT=ierr) ; CALL chk_allo (ierr, 'xaw') ALLOCATE (xae (jpai , jpaj), STAT=ierr) ; CALL chk_allo (ierr, 'xae') ALLOCATE (yan (jpai , jpaj), STAT=ierr) ; CALL chk_allo (ierr, 'yan') ALLOCATE (yas (jpai , jpaj), STAT=ierr) ; CALL chk_allo (ierr, 'yas') ALLOCATE (flux_iceberg (jpai, jpaj)) ; CALL chk_allo (ierr, 'flux_iceberg') !! ra = r_earth !! clon = 'nav_lon' ; clat = 'nav_lat' ; cmsk = 'msks' ; csrf = 'aire' ; clmd = 'grillelmd.nc' !! !CALL ipsldbg (.TRUE.) !! WRITE (UNIT=cmodel, FMT='("LMDZ ", 1I3.3, "x", 1I3.3 )') jpai, jpaj WRITE ( nout, *) cmodel !! !! Lecture de la ligne de commande narg = iargc() IF ( narg > 0 ) THEN CALL getarg ( 1, clarg) SELECT CASE (TRIM(clarg)) CASE ('-b') l_build = .TRUE. CASE ('-l') l_restart = .TRUE. c_restart = 'restart.nc' c_restartphy = 'restartphy.nc' !-$$ c_hist = 'histmth.nc' !-$$ csrf = 'AIRE' c_hist = 'gridarea_zoom_correct.nc' csrf = 'area' WRITE ( nout, *) 'Lecture LMDZ dans restart : ', TRIM (c_restart) CASE ('-r') l_read = .TRUE. clmd = 'o2a.diag.nc' cmsk = 'OceMask' ; csrf = 'Surface' WRITE ( nout, *) 'Relecture du masque calcule par mozaic' CASE Default IF ( INDEX ( TRIM(clarg), '.def' ) /= 0 ) THEN WRITE (nout,*) 'Lecture des parametres dans ', TRIM(clarg) CALL getin_name ( TRIM (clarg) ) ELSE WRITE (nout,*) 'Lecture des parametres dans run.def' ENDIF END SELECT END IF !! WRITE ( nout,*) 'Dimension atmopshere : ', jpai, jpaj, jpan WRITE ( nout,*) 'Dimension ocean : ', jpoi, jpoj, jpon !-$$ WRITE ( nout,*) 'Lecture mask : ', lec_msk !-$$ WRITE ( nout,*) 'Lecture coord : ', lec_coo !! !! WRITE ( unit = cl_rk, fmt = '(1I1)' ) rk_in cfname = 'lmd.grids' // TRIM (c_suffix) // '.r' // cl_rk OPEN ( unit = no, file = TRIM(cfname), form = 'unformatted', action = 'write', position = 'rewind', status = 'unknown' ) ! !CALL ipsldbg ( .TRUE. ) !! IF ( l_build ) THEN WRITE (nout,*) 'Construction de la grille ex nihilo' ! Lon, lat - regulier uniquement ! Contruction avec le Sud en bas de la carte DO jai = 1, jpai xlont (jai,:) = -180.0_rl + (REAL (jai-1, rl) ) / REAL (jpai , rl) * 360.0_rl xaw (jai,:) = -180.0_rl + (REAL (jai-1, rl)-0.5_rl ) / REAL (jpai , rl) * 360.0_rl xae (jai,:) = -180.0_rl + (REAL (jai-1, rl)+0.5_rl ) / REAL (jpai , rl) * 360.0_rl END DO DO jaj = 1, jpaj xlatt (:,jaj) = 90.0_rl - ( REAL (jaj-1, rl) ) / REAL (jpaj-1, rl) * 180.0_rl yas (:,jaj) = 90.0_rl - ( REAL (jaj-1, rl)+0.5_rl ) / REAL (jpaj-1, rl) * 180.0_rl yan (:,jaj) = 90.0_rl - ( REAL (jaj-1, rl)-0.5_rl ) / REAL (jpaj-1, rl) * 180.0_rl END DO ! yas = MIN ( 90.0_rl, MAX (-90.0_rl, yas)) yan = MIN ( 90.0_rl, MAX (-90.0_rl, yan)) ! xclon (:,:,1) = xae xclon (:,:,2) = xaw xclon (:,:,3) = xaw xclon (:,:,4) = xae ! xclat (:,:,1) = yas xclat (:,:,2) = yas xclat (:,:,3) = yan xclat (:,:,4) = yan DO jk = 1, 4 xclon(:,:,jk) = clo_lon ( xclon(:,:,jk), xlont(:,:) ) END DO ! Masque et fractions xmas = 0.0_rl xfra = 1.0_rl ! Surfaces CALL surface ( 'globe' ) END IF IF ( l_restart ) THEN !! Lecture restart LMDZ : le sud est en bas nlmd = jpait*(jpajv-1) + 2 WRITE (nout,*) 'Lecture restart LMD ' WRITE (nout,*) 'Dimension 1D : ', nlmd ALLOCATE ( tab1d (nlmd)) ; CALL chk_allo (ierr, 'tab1d') CALL flioopfd ( TRIM(c_restart) , id_restart) CALL flioopfd ( TRIM(c_restartphy), id_restartphy) ! WRITE (nout,*) 'Lecture longitude ', jpait*(jpajv-1), (/jpait, jpajv-1/) CALL fliogetv (id_restartphy, 'longitude', tab1d (:)) xlont (:,2:jpai-1) = RESHAPE ( tab1d (2:nlmd-1), (/jpait, jpajv-1/) ) xlont (:, 1 ) = xlont (:,jpaj/2) xlont (:,jpaj) = xlont (:,jpaj/2) WRITE (nout,*) 'Lecture latitude ', jpait*(jpajv-1), (/jpait, jpajv-1/) CALL fliogetv (id_restartphy, 'latitude' , tab1d (:)) xlatt (:,2:jpai-1) = RESHAPE ( tab1d (2:nlmd-1), (/jpait, jpajv-1/) ) xlatt (:,1) = 90.0_rl ; xlatt (:,jpaj) = -90.0_rl WRITE (nout,*) 'Lecture FOCE ', jpait*(jpajv-1), (/jpait, jpajv-1/) CALL fliogetv (id_restartphy, 'FOCE' , tab1d (:)) xfra (:,2:jpai-1) = RESHAPE ( tab1d (1:nlmd-1), (/jpait, jpajv-1/) ) WRITE (nout,*) 'Lecture FSIC ', jpait*(jpajv-1), (/jpait, jpajv-1/) CALL fliogetv (id_restartphy, 'FSIC' , tab1d (:)) xfra (:,2:jpai-1) = xfra (:,2:jpai-1) & & + RESHAPE ( tab1d (1:nlmd-1), (/jpait, jpajv-1/) ) xfra (:, 1 ) = 1.0_rl xfra (:, jpaj) = 0.0_rl xmas (:,:) = 1.0_rl WHERE ( xfra (:,:) .GT. 0.0_rl) xmas (:,:) = 0.0_rl CALL fliogetv (id_restart, 'rlonu', xlonu(:,1) ) CALL fliogetv (id_restart, 'rlatu', xlatu(1,:) ) CALL fliogetv (id_restart, 'rlonv', xlonv(:,1) ) CALL fliogetv (id_restart, 'rlatv', xlatv(1,:) ) ! Eventuellement conversion IF ( MAXVAL(ABS(xlonu)) .LT. 10.0_rl ) THEN xlonu (:,:) = SPREAD ( xlonu(:,1), DIM=2, ncopies=jpaju) * dar xlatu (:,:) = SPREAD ( xlatu(1,:), DIM=1, ncopies=jpaiu) * dar xlonv (:,:) = SPREAD ( xlonv(:,1), DIM=2, ncopies=jpajv) * dar xlatv (:,:) = SPREAD ( xlatv(1,:), DIM=1, ncopies=jpaiv) * dar END IF ! xaw (2:jpai ,:) = xlonu (1:jpai-1,1:jpaj) xaw (1 ,:) = xlonu (jpai ,1:jpaj) - 360.0_rl xae (1:jpai ,:) = xlonu (1:jpai ,1:jpaj) ! yan (1:jpai,2:jpaj ) = xlatv(1:jpai,1:jpaj-1) yan ( : , 1 ) = 90.0_rl yas (1:jpai,1:jpaj-1) = xlatv (1:jpai, 1:jpaj-1) yas ( : ,jpaj ) = -90.0_rl xclon (:,:,1) = xae xclon (:,:,2) = xaw xclon (:,:,3) = xaw xclon (:,:,4) = xae ! xclat (:,:,1) = yas xclat (:,:,2) = yas xclat (:,:,3) = yan xclat (:,:,4) = yan ! Surfaces WRITE (nout,*) 'Lecture hist LMD ', TRIM(c_hist) CALL flioopfd ( TRIM(c_hist), id_hist) ! CALL fliogetv (id_hist, TRIM(csrf), xare(:,:), start=(/1,1,1/), count=(/jpai,jpaj,1/) ) !! Bug surface dans le zoom de Yves et Rng !CALL surface ( 'pole' ) ENDIF IF ( l_read ) THEN ! Le sud n'est pas forcement au bon endroit CALL flioopfd ('grids' // TRIM(c_suffix) // '.nc', id) CALL fliogetv (id, 'tlmd.lon', xlont) CALL fliogetv (id, 'tlmd.lat', xlatt) CALL fliogetv (id, 'tlmd.clo', xclon) CALL fliogetv (id, 'tlmd.cla', xclat) CALL flioclo (id) CALL flioopfd ('areas' // TRIM(c_suffix) // '.nc', id) CALL fliogetv (id, 'tlmd.srf', xare) CALL flioclo (id) WRITE ( nout, *) 'Lecture des fractions dans ' // TRIM(clmd) // ' variable ' // TRIM(cmsk) ! Les fractions ne sont pas forcément dans le même sens que les autres fichiers ... ! CALL flioopfd ( TRIM (clmd), id) CALL fliogetv ( id, TRIM (cmsk), xfra) WRITE ( UNIT = nout, fmt = *) 'Msk lmd ', TRIM ( cmsk) CALL prihre ( xfra, 1.0_rl, nout) IF ( TRIM (cmsk) == 'phis' ) THEN xtmp = xmas xmas = 0.0_rl WHERE ( xtmp .GT. 0.1_rl ) xmas = 1.0_rl xfra = xtmp END IF IF ( TRIM (cmsk) == 'OceMask' ) THEN xmas = 0.0_rl WHERE ( xfra .EQ. 0.0_rl ) xmas = 1.0_rl ALLOCATE ( xlat_o2a(jpai,jpaj)) CALL fliogetv ( id, 'nav_lat', xlat_o2a) ! IF ( xlat_o2a(jpai/2,1) .LT. xlat_o2a(jpai/2,jpaj) .AND. xlatt(jpai/2,1) .GT. xlatt(jpai/2,jpaj) & .OR. xlat_o2a(jpai/2,1) .GT. xlat_o2a(jpai/2,jpaj) .AND. xlatt(jpai/2,1) .LT. xlatt(jpai/2,jpaj) ) THEN xfra ( :,:) = xfra (:, jpaj:1:-1) xmas ( :,:) = xmas (:, jpaj:1:-1) END IF END IF END IF !! !! Inversion eventuelle des indices J !! IF ( la_nortop ) THEN WRITE (unit = nout, fmt = *) "Cas nord en haut" IF ( xlatt(jpai/2,1) .GT. xlatt(jpai/2,jpaj) ) THEN WRITE (nout,*) "Il faut retourner les indices J" la_retourn_lat = .TRUE. ELSE WRITE (nout,*) "Indices J dans le bon ordre" la_retourn_lat = .FALSE. END IF ELSE WRITE (unit = nout, fmt = *) "Cas nord en bas" IF ( xlatt(jpai/2,1) .LT. xlatt(jpai/2,jpaj) ) THEN WRITE (nout,*) "Il faut retourner les indices J" la_retourn_lat = .TRUE. ELSE WRITE (nout,*) "Indices J dans le bon ordre" la_retourn_lat = .FALSE. END IF END IF IF ( la_retourn_lat ) THEN WRITE (nout,*) "Retournement des indices J" xlont (:,:) = xlont (1:jpai , jpaj :1:-1) xlonu (:,:) = xlonu (1:jpaiu, jpaju:1:-1) xlonv (:,:) = xlonv (1:jpaiv, jpajv:1:-1) xlatt (:,:) = xlatt (1:jpai , jpaj :1:-1) xlatu (:,:) = xlatu (1:jpaiu, jpaju:1:-1) xlatv (:,:) = xlatv (1:jpaiv, jpajv:1:-1) xmas (:,:) = xmas (1:jpai , jpaj :1:-1) xfra (:,:) = xfra (1:jpai , jpaj :1:-1) xare (:,:) = xare (1:jpai , jpaj :1:-1) xclon (:,:,:) = xclon (1:jpai , jpaj :1:-1,:) xclat (:,:,:) = xclat (1:jpai , jpaj :1:-1,:) xaw (:,:) = xaw (1:jpai , jpaj :1:-1) xae (:,:) = xae (1:jpai , jpaj :1:-1) yas (:,:) = yas (1:jpai , jpaj :1:-1) yan (:,:) = yan (1:jpai , jpaj :1:-1) xclon (:,:,1) = xaw (:,:) xclon (:,:,2) = xae (:,:) xclon (:,:,3) = xae (:,:) xclon (:,:,4) = xaw (:,:) ! xclat (:,:,1) = yan (:,:) xclat (:,:,2) = yan (:,:) xclat (:,:,3) = yas (:,:) xclat (:,:,4) = yas (:,:) END IF WRITE ( unit = nout, fmt = *) 'Longitudes T' CALL prihre ( xlont, 1._rl, nout) WRITE ( unit = no) 't'//TRIM(camod)//'.lon' WRITE ( unit = no) REAL ( xlont(1:jpai,1:jpaj), kind = rk_in ) ! WRITE ( unit = nout, fmt = *) 'Latitudes T' CALL prihre ( xlatt, 1._rl, nout) WRITE ( unit = no) 't'//TRIM(camod)//'.lat' WRITE ( unit = no) REAL (xlatt(1:jpai,1:jpaj), kind = rk_in ) ! WRITE ( unit = nout, fmt = *) 'Corners' WRITE ( unit = no) 't'//TRIM(camod)//'.clo' WRITE ( unit = no) REAL (xclon, kind = rk_in ) WRITE ( unit = no) 't'//TRIM(camod)//'.cla' WRITE ( unit = no) REAL (xclat, kind = rk_in ) DO jc = 1, 4 WRITE ( unit = nout, fmt = *) 'Coins - Longitudes ', jc CALL prihre ( xclon(:,:,jc), 1._rl, nout) WRITE ( unit = nout, fmt = *) 'Coins - Latitudes ', jc CALL prihre ( xclat(:,:,jc), 1._rl, nout) END DO WRITE ( unit = nout, fmt = *) 'Fractions (%)' CALL prihre ( xfra, 100.0_rl, nout) WRITE ( unit = nout, fmt = *) 'Surface ' CALL prihre (xare, 0.0_rl, nout) zzz = SUM ( xare) / ( ra * ra * rpi * 4.0_rl ) WRITE (UNIT = nout, FMT = *) 'Surface normalisee : ', zzz !! !! Ecriture OASIS NetCDF CALL fliocrfd ( 'lmd.grids' // TRIM(c_suffix), (/'jpia ', 'jpja ', 'corners'/), (/jpai, jpaj, 4/), il_grid, mode='REPLACE') CALL flioputa ( il_grid, '?', 'Model', cmodel) CALL flioputa ( il_grid, '?', 'Comment', TRIM(c_comment) ) CALL fliodefv ( il_grid, 'tlmd.lon', (/1, 2/), v_t=flio_r8 ) CALL flioputa ( il_grid, 'tlmd.lon', 'units' , 'degrees_east' ) CALL flioputa ( il_grid, 'tlmd.lon', 'title' , 'lmdz-t longitudes' ) CALL flioputa ( il_grid, 'tlmd.lon', 'grid_type' , 'P' ) CALL flioputa ( il_grid, 'tlmd.lon', 'overlap' , 0_i_4 ) CALL fliodefv ( il_grid, 'tlmd.lat', (/1, 2/), v_t=flio_r8 ) CALL flioputa ( il_grid, 'tlmd.lat', 'units' , 'degrees_north' ) CALL flioputa ( il_grid, 'tlmd.lat', 'title' ,'lmdz-t latitudes' ) CALL flioputa ( il_grid, 'tlmd.lat', 'grid_type' ,'P' ) CALL flioputa ( il_grid, 'tlmd.lat', 'overlap' , 0_i_4 ) CALL fliodefv ( il_grid, 'tlmd.clo', (/1, 2, 3/), v_t=flio_r8 ) CALL flioputa ( il_grid, 'tlmd.clo', 'units' , 'degrees_east' ) CALL flioputa ( il_grid, 'tlmd.clo', 'title' ,'Longitudes for t-cell corner' ) CALL flioputa ( il_grid, 'tlmd.clo', 'grid_type' ,'P' ) CALL flioputa ( il_grid, 'tlmd.clo', 'overlap' , 0_i_4 ) CALL fliodefv ( il_grid, 'tlmd.cla', (/1, 2, 3/), v_t=flio_r8 ) CALL flioputa ( il_grid, 'tlmd.cla', 'units' , 'degrees_north' ) CALL flioputa ( il_grid, 'tlmd.cla', 'title' ,'Latitudes for t-cell corner' ) CALL flioputa ( il_grid, 'tlmd.cla', 'grid_type' ,'P' ) CALL flioputa ( il_grid, 'tlmd.cla', 'overlap' , 0_i_4 ) CALL fliodefv ( il_grid, 'aone.lon', (/1, 2/), v_t=flio_r8 ) CALL flioputa ( il_grid, 'aone.lon', 'units' , 'degrees_east' ) CALL flioputa ( il_grid, 'aone.lon', 'title' , 'lmdz-t longitudes' ) CALL flioputa ( il_grid, 'aone.lon', 'grid_type' , 'P' ) CALL flioputa ( il_grid, 'aone.lon', 'overlap' , 0_i_4 ) CALL fliodefv ( il_grid, 'aone.lat', (/1, 2/), v_t=flio_r8 ) CALL flioputa ( il_grid, 'aone.lat', 'units' , 'degrees_north' ) CALL flioputa ( il_grid, 'aone.lat', 'title' ,'lmdz-t latitudes' ) CALL flioputa ( il_grid, 'aone.lat', 'grid_type' ,'P' ) CALL flioputa ( il_grid, 'aone.lat', 'overlap' , 0_i_4 ) CALL fliodefv ( il_grid, 'aone.clo', (/1, 2, 3/), v_t=flio_r8 ) CALL flioputa ( il_grid, 'aone.clo', 'units' , 'degrees_east' ) CALL flioputa ( il_grid, 'aone.clo', 'title' ,'Longitudes for t-cell corner' ) CALL flioputa ( il_grid, 'aone.clo', 'grid_type' ,'P' ) CALL flioputa ( il_grid, 'aone.clo', 'overlap' , 0_i_4 ) CALL fliodefv ( il_grid, 'aone.cla', (/1, 2, 3/), v_t=flio_r8 ) CALL flioputa ( il_grid, 'aone.cla', 'units' , 'degrees_north' ) CALL flioputa ( il_grid, 'aone.cla', 'title' ,'Latitudes for t-cell corner' ) CALL flioputa ( il_grid, 'aone.cla', 'grid_type' ,'P' ) CALL flioputa ( il_grid, 'aone.cla', 'overlap' , 0_i_4 ) CALL fliodefv ( il_grid, 'afra.lon', (/1, 2/), v_t=flio_r8 ) CALL flioputa ( il_grid, 'afra.lon', 'units' , 'degrees_east' ) CALL flioputa ( il_grid, 'afra.lon', 'title' , 'lmdz-t longitudes' ) CALL flioputa ( il_grid, 'afra.lon', 'grid_type' , 'P' ) CALL flioputa ( il_grid, 'afra.lon', 'overlap' , 0_i_4 ) CALL fliodefv ( il_grid, 'afra.lat', (/1, 2/), v_t=flio_r8 ) CALL flioputa ( il_grid, 'afra.lat', 'units' , 'degrees_north' ) CALL flioputa ( il_grid, 'afra.lat', 'title' ,'lmdz-t latitudes' ) CALL flioputa ( il_grid, 'afra.lat', 'grid_type' ,'P' ) CALL flioputa ( il_grid, 'afra.lat', 'overlap' , 0_i_4 ) CALL fliodefv ( il_grid, 'afra.clo', (/1, 2, 3/), v_t=flio_r8 ) CALL flioputa ( il_grid, 'afra.clo', 'units' , 'degrees_east' ) CALL flioputa ( il_grid, 'afra.clo', 'title' ,'Longitudes for t-cell corner' ) CALL flioputa ( il_grid, 'afra.clo', 'grid_type' ,'P' ) CALL flioputa ( il_grid, 'afra.clo', 'overlap' , 0_i_4 ) CALL fliodefv ( il_grid, 'afra.cla', (/1, 2, 3/), v_t=flio_r8 ) CALL flioputa ( il_grid, 'afra.cla', 'units' , 'degrees_north' ) CALL flioputa ( il_grid, 'afra.cla', 'title' ,'Latitudes for t-cell corner' ) CALL flioputa ( il_grid, 'afra.cla', 'grid_type' ,'P' ) CALL flioputa ( il_grid, 'afra.cla', 'overlap' , 0_i_4 ) CALL flioputv ( il_grid, 'tlmd.lon', xlont ) CALL flioputv ( il_grid, 'tlmd.lat', xlatt ) CALL flioputv ( il_grid, 'tlmd.clo', xclon ) CALL flioputv ( il_grid, 'tlmd.cla', xclat ) CALL flioputv ( il_grid, 'aone.lon', xlont ) CALL flioputv ( il_grid, 'aone.lat', xlatt ) CALL flioputv ( il_grid, 'aone.clo', xclon ) CALL flioputv ( il_grid, 'aone.cla', xclat ) CALL flioputv ( il_grid, 'afra.lon', xlont ) CALL flioputv ( il_grid, 'afra.lat', xlatt ) CALL flioputv ( il_grid, 'afra.clo', xclon ) CALL flioputv ( il_grid, 'afra.cla', xclat ) ! CALL flioclo ( il_grid) ! CLOSE ( unit = no) !! WRITE ( unit = nout, fmt = *) 'Masque atm ', i_4 WRITE ( unit = cl_rk, fmt = '(1I1)' ) i_4 cfname = 'lmd.masks' // TRIM (c_suffix) // '.i' // cl_rk OPEN ( unit = no, file = cfname, form = 'unformatted', action = 'write', position = 'rewind', status = 'unknown' ) CALL prihin ( INT (xmas,KIND=il), 1_il, nout) WRITE ( unit = no) 't'//TRIM(camod)//'.msk' WRITE ( unit = no) INT ( xmas, KIND = i_4) WRITE ( unit = no) 'u'//TRIM(camod)//'.msk' WRITE ( unit = no) INT ( xmas, KIND = i_4) WRITE ( unit = no) 'v'//TRIM(camod)//'.msk' WRITE ( unit = no) INT ( xmas, KIND = i_4) CLOSE ( unit = no) !! WRITE ( unit = nout, fmt = *) 'Masque atm ', i_8 WRITE ( unit = cl_rk, fmt = '(1I1)' ) i_8 cfname = 'lmd.masks' // TRIM (c_suffix) // '.i' // cl_rk OPEN ( unit = no, file = cfname, form = 'unformatted', action = 'write', position = 'rewind', status = 'unknown' ) CALL prihin ( INT (xmas,KIND=il), 1_il, nout) WRITE ( unit = no) 't'//TRIM(camod)//'.msk' WRITE ( unit = no) INT ( xmas, KIND = i_8) WRITE ( unit = no) 'u'//TRIM(camod)//'.msk' WRITE ( unit = no) INT ( xmas, KIND = i_8) WRITE ( unit = no) 'v'//TRIM(camod)//'.msk' WRITE ( unit = no) INT ( xmas, KIND = i_8) CLOSE ( unit = no) !! CALL fliocrfd ( 'lmd.masks' // TRIM(c_suffix), (/'jpia', 'jpja' /), (/jpai, jpaj/), il_mask, mode='REPLACE' ) CALL flioputa ( il_mask, '?', 'Model', cmodel ) CALL flioputa ( il_mask, '?', 'Comment', TRIM(c_comment) ) CALL fliodefv ( il_mask, 'tlmd.msk', (/1, 2/), v_t=flio_i4 ) CALL flioputa ( il_mask, 'tlmd.msk', 'units' , 'n' ) CALL flioputa ( il_mask, 'tlmd.msk', 'land_value', 1_i_4) CALL flioputa ( il_mask, 'tlmd.msk', 'sea_value' , 0_i_4) CALL flioputa ( il_mask, 'tlmd.msk', 'title' , 'lmdz-t land-sea mask' ) CALL fliodefv ( il_mask, 'aone.msk', (/1, 2/), v_t=flio_i4 ) CALL flioputa ( il_mask, 'aone.msk', 'units' , 'n' ) CALL flioputa ( il_mask, 'aone.msk', 'land_value', 1_i_4) CALL flioputa ( il_mask, 'aone.msk', 'sea_value' , 0_i_4) CALL flioputa ( il_mask, 'aone.msk', 'title' , 'lmdz-t land-sea mask, forced to not masked' ) CALL fliodefv ( il_mask, 'afra.msk', (/1, 2/), v_t=flio_i4 ) CALL flioputa ( il_mask, 'afra.msk', 'units' , 'n' ) CALL flioputa ( il_mask, 'afra.msk', 'land_value', 1_i_4) CALL flioputa ( il_mask, 'afra.msk', 'sea_value' , 0_i_4) CALL flioputa ( il_mask, 'afra.msk', 'title' , 'lmdz-t land-sea mask' ) !! CALL flioputv ( il_mask, 'tlmd.msk', xmas ) CALL flioputv ( il_mask, 'aone.msk', 0.0*xmas ) CALL flioputv ( il_mask, 'afra.msk', xmas ) CALL flioclo ( il_mask ) !! WRITE ( unit = cl_rk, fmt = '(1I1)' ) rk_in cfname = 'lmd.areas' // TRIM (c_suffix) // '.r' // cl_rk OPEN ( unit = no, file = cfname, form = 'unformatted', action = 'write', position = 'rewind', status = 'unknown' ) ! WRITE (unit = no) 't'//TRIM(camod)//'.srf' WRITE (unit = no) REAL ( xare, KIND = rk_in) ! CLOSE ( unit = no) !! CALL fliocrfd ( 'lmd.areas' // TRIM(c_suffix), (/'jpia', 'jpja' /), (/jpai, jpaj/), il_area, mode='REPLACE' ) CALL flioputa ( il_area, '?', 'Model', cmodel ) CALL flioputa ( il_area, '?', 'Comment', TRIM(c_comment) ) CALL fliodefv ( il_area, 'tlmd.srf', (/1, 2/), v_t=flio_r8 ) CALL flioputa ( il_area, 'tlmd.srf', 'units' , 'm^2' ) CALL flioputa ( il_area, 'tlmd.srf', 'title', 'lmdz-t mesh surface') CALL fliodefv ( il_area, 'aone.srf', (/1, 2/), v_t=flio_r8 ) CALL flioputa ( il_area, 'aone.srf', 'units' , 'm^2' ) CALL flioputa ( il_area, 'aone.srf', 'title', 'lmdz-t mesh surface, set to 1') CALL fliodefv ( il_area, 'afra.srf', (/1, 2/), v_t=flio_r8 ) CALL flioputa ( il_area, 'afra.srf', 'units' , 'm^2' ) CALL flioputa ( il_area, 'afra.srf', 'title', 'lmdz-t mesh surface, set to 1') CALL flioputv ( il_area, 'tlmd.srf', xare ) CALL flioputv ( il_area, 'aone.srf', 1.0_rl + 0.0_rl*xare ) CALL flioputv ( il_area, 'afra.srf', xare*xfra ) CALL flioclo ( il_area) ! !! CALL fliocrfd ( 'lmd.fracs' // TRIM(c_suffix), (/'jpia', 'jpja' /), (/jpai, jpaj/), il_frac, mode='REPLACE' ) CALL flioputa ( il_frac, '?', 'Model', cmodel ) CALL flioputa ( il_frac, '?', 'Comment', TRIM(c_comment) ) CALL fliodefv ( il_frac, 'tlmd.fra', (/1, 2/), v_t=flio_r8 ) CALL flioputa ( il_frac, 'tlmd.fra', 'units' , 'm^2' ) CALL flioputa ( il_frac, 'tlmd.fra', 'title', 'lmdz-t mesh fraction') CALL fliodefv ( il_frac, 'aone.fra', (/1, 2/), v_t=flio_r8 ) CALL flioputa ( il_frac, 'aone.fra', 'units' , 'm^2' ) CALL flioputa ( il_frac, 'aone.fra', 'title', 'lmdz-t mesh fraction, set to 1') CALL fliodefv ( il_frac, 'afra.fra', (/1, 2/), v_t=flio_r8 ) CALL flioputa ( il_frac, 'afra.fra', 'units' , 'm^2' ) CALL flioputa ( il_frac, 'afra.fra', 'title', 'lmdz-t mesh fraction, set to 1') CALL flioputv ( il_frac, 'tlmd.fra', xfra ) CALL flioputv ( il_frac, 'aone.fra', 1.0_rl + 0.0_rl*xfra ) CALL flioputv ( il_frac, 'afra.fra', xfra ) CALL flioclo ( il_frac) ! IF ( l_read ) THEN !! Sortie fichier bathy OPEN ( unit = 27, file = 'bathy.lmdz', FORM = 'formatted', STATUS = 'unknown', ACTION = 'write', POSITION = 'rewind' ) DO jaj = jpaj, 1, -1 WRITE ( unit = 27, FMT = '(92I1)' ) INT ( xmas(:,jaj), KIND = i_4) END DO CLOSE ( unit = 27) !! Sortie fichier bathy flux iceberg flux_iceberg = 0.0 DO jai = 1, jpai DO jaj = 1, jpaj IF ( xmas(jai,jaj) == 0 ) THEN flux_iceberg(jai,jaj ) = 2.0 flux_iceberg(jai,jaj+1) = 2.0 flux_iceberg(jai,jaj+2) = 1.0 flux_iceberg(jai,jaj+3) = 1.0 EXIT END IF END DO END DO !! OPEN ( unit = 27, file = 'flux_iceberg', FORM = 'formatted', STATUS = 'unknown', ACTION = 'write', POSITION = 'rewind' ) DO jaj = jpaj, 1, -1 WRITE ( unit = 27, FMT = '(96F3.0)' ) flux_iceberg(:,jaj) END DO END IF !! STOP CONTAINS ELEMENTAL FUNCTION clo_lon ( zlon, zlon0) ! Find closest longitude IMPLICIT NONE REAL (kind=rl) :: clo_lon REAL (kind=rl), INTENT ( in) :: zlon, zlon0 REAL (kind=rl) :: zz, zp, zm, z0 INTEGER (kind=il) :: jn z0 = zlon DO jn = 1_il, 2_il zz = ABS ( (z0 ) - zlon0 ) zp = ABS ( (z0 + 360.0_rl) - zlon0 ) zm = ABS ( (z0 - 360.0_rl) - zlon0 ) IF ( zp < MIN ( zm, zz) ) THEN z0 = z0 + 360.0_rl ELSE IF ( zm < MIN ( zp, zz) ) THEN z0 = z0 - 360.0_rl END IF END DO clo_lon = z0 END FUNCTION clo_lon SUBROUTINE surface ( c_type ) IMPLICIT NONE CHARACTER (len=*), OPTIONAL :: c_type CHARACTER (len=8) :: cl_type IF ( PRESENT (c_type)) THEN cl_type = TRIM (c_type) ELSE cl_type = 'globe' END IF IF ( TRIM(c_type) == 'globe') THEN WRITE (nout,*) 'Calcul des surfaces atm partout' xare (:,:) = (xae(:,:)-xaw(:,:))*(yan(:,:)-yas(:,:)) * COS(rad*xlatt(:,:)) * rad * rad * ra * ra END IF WRITE (nout,*) 'Calcul des surfaces atm aux poles' xare (:, 1 ) = (xae(:, 1)-xaw(:, 1 ))*(1.0_rl - ABS(SIN(rad*yas(:, 1 )))) * rad * ra * ra xare (:,jpaj) = (xae(:,jpaj)-xaw(:,jpaj))*(1.0_rl - ABS(SIN(rad*yan(:,jpaj)))) * rad * ra * ra END SUBROUTINE surface SUBROUTINE lmd1d_2d ( ztab1, ztab2 ) REAL (kind=rl), INTENT (in) , DIMENSION (:) :: ztab1 REAL (kind=rl), INTENT (out), DIMENSION (:,:) :: ztab2 INTEGER :: jn, jn_min, jn_max jn_min = 9999999 jn_max = -9999999 ztab2 ( :, 1 ) = ztab1 ( 1) ztab2 ( :, jpai) = ztab1 ( nlmd) DO jaj = 2, jpai-1 DO jai = 1, jpai jn = 1 + (jaj-2)*jpai + jai jn_min = MIN ( jn, jn_min ) jn_max = MIN ( jn, jn_max ) ztab2 (jai, jaj) = ztab1 ( jn ) END DO END DO WRITE (nout,*) 'lmd1d_2d : ', jn_min, jn_max END SUBROUTINE lmd1d_2d END PROGRAM lmdgrille