Changeset 173 for LMDZ.3.3/branches/rel-LF
- Timestamp:
- Jan 18, 2001, 4:43:58 PM (24 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ.3.3/branches/rel-LF/libf/dyn3d/create_limit.F
r136 r173 1 c $Header$ 1 2 PROGRAM create_limit 3 USE startvar 4 USE ioipsl 2 5 IMPLICIT none 3 6 c … … 24 27 #include "comgeom2.h" 25 28 #include "comconst.h" 26 c 29 #include "dimphy.h" 30 #include "indicesol.h" 27 31 c----------------------------------------------------------------------- 28 INTEGER KIDIA, KFDIA, KLON, KLEV 29 PARAMETER (KIDIA=1,KFDIA=iim*(jjm-1)+2, 30 . KLON=KFDIA-KIDIA+1,KLEV=llm) 31 c----------------------------------------------------------------------- 32 REAL phy_nat(klon,360), phy_nat0(klon) 32 REAL phy_nat(klon,360) 33 real phy_nat0(klon) 33 34 REAL phy_alb(klon,360) 34 35 REAL phy_sst(klon,360) 35 36 REAL phy_bil(klon,360) 36 37 REAL phy_rug(klon,360) 37 REAL phy_ice(klon ,360)38 REAL phy_ice(klon) 38 39 CPB 39 REAL phy_icet(klon,360) 40 REAL phy_oce(klon,360) 40 c REAL phy_icet(klon,360) 41 c REAL phy_oce(klon,360) 42 real pctsrf_t(klon,nbsrf,360) 43 real pctsrf(klon,nbsrf) 41 44 REAL verif 42 45 c … … 89 92 INTEGER id_FOCE, id_FSIC, id_FTER, id_FLIC 90 93 91 INTEGER i, j, k, l 94 INTEGER i, j, k, l, ji 95 c declarations pour lecture glace de mer 96 INTEGER :: iml_lic, jml_lic, llm_tmp, ttm_tmp, iret 97 INTEGER :: itaul(1), fid 98 REAL :: lev(1), date, dt 99 REAL, ALLOCATABLE, DIMENSION(:,:) :: lon_lic, lat_lic 100 REAL, ALLOCATABLE, DIMENSION(:) :: dlon_lic, dlat_lic 101 REAL, ALLOCATABLE, DIMENSION (:,:) :: fraclic 102 REAL :: flic_tmp(iip1, jjp1) 103 REAL :: champint(iim, jjp1) 92 104 c Diverses variables locales 93 105 REAL time 106 ! pour la lecture du fichier masque ocean 107 integer :: nid_o2a 108 logical :: couple = .false. 109 INTEGER :: iml_omask, jml_omask 110 REAL, ALLOCATABLE, DIMENSION(:,:) :: lon_omask, lat_omask 111 REAL, ALLOCATABLE, DIMENSION(:) :: dlon_omask, dlat_omask 112 REAL, ALLOCATABLE, DIMENSION (:,:) :: ocemask, ocetmp 113 real, dimension(klon) :: ocemask_fi 114 94 115 95 116 INTEGER longcles … … 99 120 INTEGER ncid,varid,ndimid(4),dimid 100 121 character*30 namedim 122 CHARACTER*80 :: varname 101 123 102 124 c initialisations: … … 122 144 C Traitement du relief au sol 123 145 c 124 write(*,*) 'Traitement du relief au sol pour fabriquer masque' 125 ierr = NF_OPEN('Relief.nc', NF_NOWRITE, ncid) 126 127 if (ierr.ne.0) then 128 print *, NF_STRERROR(ierr) 129 STOP 130 ENDIF 131 132 ierr = NF_INQ_VARID(ncid,'RELIEF',varid) 133 if (ierr.ne.0) then 134 print *, NF_STRERROR(ierr) 135 STOP 136 ENDIF 137 ierr = NF_INQ_VARDIMID (ncid,varid,ndimid) 138 if (ierr.ne.0) then 139 print *, NF_STRERROR(ierr) 140 STOP 141 ENDIF 142 ierr = NF_INQ_DIM(ncid,ndimid(1), namedim, imdep) 143 if (ierr.ne.0) then 144 print *, NF_STRERROR(ierr) 145 STOP 146 ENDIF 147 print*,'variable ', namedim, 'dimension ', imdep 148 ierr = NF_INQ_VARID(ncid,namedim,dimid) 149 if (ierr.ne.0) then 150 print *, NF_STRERROR(ierr) 151 STOP 152 ENDIF 153 ierr = NF_GET_VAR_REAL(ncid,dimid,dlon_msk) 154 if (ierr.ne.0) then 155 print *, NF_STRERROR(ierr) 156 STOP 157 ENDIF 158 ierr = NF_INQ_DIM(ncid,ndimid(2), namedim, jmdep) 159 if (ierr.ne.0) then 160 print *, NF_STRERROR(ierr) 161 STOP 162 ENDIF 163 print*,'variable ', namedim, 'dimension ', jmdep 164 ierr = NF_INQ_VARID(ncid,namedim,dimid) 165 if (ierr.ne.0) then 166 print *, NF_STRERROR(ierr) 167 STOP 168 ENDIF 169 ierr = NF_GET_VAR_REAL(ncid,dimid,dlat_msk) 170 if (ierr.ne.0) then 171 print *, NF_STRERROR(ierr) 172 STOP 173 ENDIF 174 ierr = NF_GET_VAR_REAL(ncid,varid,champ_msk) 175 if (ierr.ne.0) then 176 print *, NF_STRERROR(ierr) 177 STOP 178 ENDIF 179 180 c 181 CALL mask_c_o(imdep, jmdep, dlon_msk, dlat_msk,champ_msk, 182 . iim, jjp1, rlonv, rlatu, champint) 183 CALL gr_int_dyn(champint, masque, iim, jjp1) 146 write(*,*) 'Fabrication masque' 147 148 varname = 'masque' 149 masque(:,:) = 0.0 150 CALL startget(varname, iip1, jjp1, rlonv, rlatu, masque, 0.0) 151 pctsrf=0. 152 varname = 'zmasq' 153 zmasq(:) = 0. 154 CALL startget(varname,iip1,jjp1,rlonv,rlatu,klon,zmasq,0.0) 155 WHERE (zmasq(1 : klon) .LE. EPSFRA) 156 zmasq(1 : klon) = 0. 157 END WHERE 158 ! WRITE(*,*)zmasq 159 184 160 IF ( fracterre ) THEN 185 161 DO i = 1, iim … … 195 171 DO i = 1, iim 196 172 DO j = 1, jjp1 197 mask(i,j) = champint(i,j)173 mask(i,j) = masque(i,j) 198 174 ENDDO 199 175 ENDDO 200 176 CALL gr_dyn_fi(1, iip1, jjp1, klon, masque, phy_nat0) 201 ierr = NF_CLOSE(ncid) 177 C 178 C En cas de simulation couplee, lecture du masque ocean issu du modele ocean 179 C utilise pour calculer les poids et pour assurer l'adequation entre les 180 C fractions d'ocean vu par l'atmosphere et l'ocean 181 C 182 183 write(*,*)'Essai de lecture masque ocean' 184 iret = nf_open("o2a.nc", NF_NOWRITE, nid_o2a) 185 if (iret .ne. 0) then 186 write(*,*)'ATTENTION!! pas de fichier o2a.nc trouve' 187 write(*,*)'Run force' 188 else 189 couple = .true. 190 iret = nf_close(nid_o2a) 191 call flininfo("o2a.nc", iml_omask, jml_omask, llm_tmp, ttm_tmp 192 $ , nid_o2a) 193 if (iml_omask /= iim .or. jml_omask /= jjp1) then 194 write(*,*)'Dimensions non compatibles pour masque ocean' 195 write(*,*)'iim = ',iim,' iml_omask = ',iml_omask 196 write(*,*)'jjp1 = ',jjp1,' jml_omask = ',jml_omask 197 stop 198 endif 199 ALLOCATE(lat_omask(iml_omask, jml_omask), stat=iret) 200 ALLOCATE(lon_omask(iml_omask, jml_omask), stat=iret) 201 ALLOCATE(dlon_omask(iml_omask), stat=iret) 202 ALLOCATE(dlat_omask(jml_omask), stat=iret) 203 ALLOCATE(ocemask(iml_omask, jml_omask), stat=iret) 204 ALLOCATE(ocetmp(iml_omask, jml_omask), stat=iret) 205 CALL flinopen("o2a.nc", .FALSE., iml_omask, jml_omask, llm_tmp 206 $ , lon_omask, lat_omask, lev, ttm_tmp, itaul, date, dt, fid) 207 CALL flinget(fid, 'OceMask', iml_omask, jml_omask, llm_tmp, 208 $ ttm_tmp, 1, 1, ocetmp) 209 CALL flinclo(fid) 210 dlon_omask(1 : iml_omask) = lon_omask(1 : iml_omask, 1) 211 dlat_omask(1 : jml_omask) = lat_omask(1 , 1 : jml_omask) 212 ocemask = ocetmp 213 if (dlat_omask(1) < dlat_omask(jml_omask)) then 214 do j = 1, jml_omask 215 ocemask(:,j) = ocetmp(:,jml_omask-j+1) 216 enddo 217 endif 218 C 219 C passage masque ocean a la grille physique 220 C 221 ocemask_fi(1) = ocemask(1,1) 222 do j = 2, jjm 223 do i = 1, iim 224 ocemask_fi((j-2)*iim + i + 1) = ocemask(i,j) 225 enddo 226 enddo 227 ocemask_fi(klon) = ocemask(1,jjp1) 228 zmasq = 1. - ocemask_fi 229 endif 230 231 232 C 233 C lecture du fichier glace de terre pour fixer la fraction de terre 234 C et de glace de terre 235 C 236 CALL flininfo("landiceref.nc", iml_lic, jml_lic,llm_tmp, ttm_tmp 237 $ , fid) 238 ALLOCATE(lat_lic(iml_lic, jml_lic), stat=iret) 239 ALLOCATE(lon_lic(iml_lic, jml_lic), stat=iret) 240 ALLOCATE(dlon_lic(iml_lic), stat=iret) 241 ALLOCATE(dlat_lic(jml_lic), stat=iret) 242 ALLOCATE(fraclic(iml_lic, jml_lic), stat=iret) 243 CALL flinopen("landiceref.nc", .FALSE., iml_lic, jml_lic, llm_tmp 244 $ , lon_lic, lat_lic, lev, ttm_tmp, itaul, date, dt, fid) 245 CALL flinget(fid, 'landice', iml_lic, jml_lic, llm_tmp, ttm_tmp 246 $ , 1, 1, fraclic) 247 CALL flinclo(fid) 248 C 249 C interpolation sur la grille T du modele 250 C 251 WRITE(*,*) 'dimensions de landice iml_lic, jml_lic : ', 252 $ iml_lic, jml_lic 253 c 254 C sil les coordonnees sont en degres, on les transforme 255 C 256 IF( MAXVAL( lon_lic(:,:) ) .GT. 2.0 * asin(1.0) ) THEN 257 lon_lic(:,:) = lon_lic(:,:) * 2.0* ASIN(1.0) / 180. 258 ENDIF 259 IF( maxval( lat_lic(:,:) ) .GT. 2.0 * asin(1.0)) THEN 260 lat_lic(:,:) = lat_lic(:,:) * 2.0 * asin(1.0) / 180. 261 ENDIF 262 263 dlon_lic(1 : iml_lic) = lon_lic(1 : iml_lic, 1) 264 dlat_lic(1 : jml_lic) = lat_lic(1 , 1 : jml_lic) 265 C 266 CALL grille_m(iml_lic, jml_lic, dlon_lic, dlat_lic, fraclic 267 $ ,iim, jjp1, 268 $ rlonv, rlatu, flic_tmp(1 : iim, 1 : jjp1)) 269 c$$$ flic_tmp(1 : iim, 1 : jjp1) = champint(1: iim, 1 : jjp1) 270 flic_tmp(iip1, 1 : jjp1) = flic_tmp(1 , 1 : jjp1) 271 C 272 C passage sur la grille physique 273 C 274 CALL gr_dyn_fi(1, iip1, jjp1, klon, flic_tmp, 275 $ pctsrf(1:klon, is_lic)) 276 C adequation avec le maque terre/mer 277 WHERE (pctsrf(1 : klon, is_lic) .LE. EPSFRA ) 278 pctsrf(1 : klon, is_lic) = 0. 279 END WHERE 280 WHERE (zmasq( 1 : klon) .LE. EPSFRA) 281 pctsrf(1 : klon, is_lic) = 0. 282 END WHERE 283 pctsrf(1 : klon, is_ter) = zmasq(1 : klon) 284 DO ji = 1, klon 285 IF (zmasq(ji) .GT. EPSFRA) THEN 286 IF ( pctsrf(ji, is_lic) .GE. zmasq(ji)) THEN 287 pctsrf(ji, is_lic) = zmasq(ji) 288 pctsrf(ji, is_ter) = 0. 289 ELSE 290 pctsrf(ji,is_ter) = zmasq(ji) - pctsrf(ji, is_lic) 291 ENDIF 292 ENDIF 293 END DO 294 202 295 c 203 296 c … … 295 388 ENDDO 296 389 ENDDO 390 c write(70,*)champtime 297 391 c 298 392 DO l = 1, lmdep … … 453 547 ENDDO 454 548 c 455 WRITE(*,*) 'phy_nat'549 c WRITE(*,*) 'phy_nat' 456 550 c WRITE(*,'(72f4.1)') phy_nat0(1:klon) 457 551 c 458 552 DO k = 1, 360 459 553 CALL gr_dyn_fi(1, iip1, jjp1, klon, 460 . champan(1,1,k), phy_ice (1,k))554 . champan(1,1,k), phy_ice) 461 555 IF ( newlmt) THEN 462 556 463 557 CPB en attendant de mettre fraction de terre 464 558 c 465 WHERE(phy_ice(1:klon, k) .GT. 1.) phy_ice(1 : klon, k) = 1. 466 WHERE(phy_ice(1:klon, k) .LT. 0.) phy_ice(1 : klon, k) = 0. 467 WRITE(*,*) 'phy_ice : ', k 468 c WRITE(*,'(72f4.1)') phy_ice(1 : klon, k) 559 WHERE(phy_ice(1:klon) .GT. 1.) phy_ice(1 : klon) = 1. 560 WHERE(phy_ice(1:klon) .LT. EPSFRA) phy_ice(1 : klon) = 0. 469 561 c 470 471 472 DO i = 1, klon473 phy_nat(i,k) = phy_nat0(i)474 IF (phy_nat0(i) .GE. 0.5 ) THEN475 IF(phy_ice(i,k) .GE. 1.e-5) THEN476 IF ( phy_ice(i,k) .LE. phy_nat(i,k)) THEN477 phy_oce(i,k) = 1. - phy_nat(i,k)478 phy_icet(i,k) = phy_ice(i,k)479 phy_ice(i,k) = 0.480 phy_nat(i,k)= phy_nat(i,k)- phy_icet(i,k)481 ELSE482 phy_oce(i,k) = 1. - phy_ice(i,k)483 phy_icet(i,k) = phy_nat(i,k)484 phy_nat(i,k) = 0.485 phy_ice(i,k) = phy_ice(i,k)486 $ - phy_icet(i,k)487 ENDIF488 ELSE489 phy_icet(i,k) = 0.490 phy_ice(i,k) = 0.491 phy_oce(i,k) = 1. - phy_nat(i,k)492 ENDIF493 ELSE494 phy_oce(i,k) = 1. - phy_nat(i,k)495 IF(phy_ice(i,k) .GE. 1.e-5) THEN496 IF( phy_ice(i,k) .LE. phy_oce(i,k) ) THEN497 phy_icet(i,k) = 0.498 phy_oce(i,k) = phy_oce(i,k) - phy_ice(i,k)499 ELSE500 phy_icet(i,k)=phy_ice(i,k) - phy_oce(i,k)501 phy_ice(i,k) = phy_oce(i,k)502 phy_oce(i,k) = 0.503 phy_nat(i,k) = phy_nat(i,k)-phy_icet(i,k)504 ENDIF505 ELSE506 phy_icet(i,k) = 0.507 phy_ice(i,k) = 0.508 ENDIF509 ENDIF510 verif = phy_nat(i,k) + phy_icet(i,k)+ phy_ice(i,k)511 $ + phy_oce(i,k)512 513 $ 514 562 IF (fracterre ) THEN 563 c WRITE(*,*) 'passe dans cas fracterre' 564 pctsrf_t(:,is_ter,k) = pctsrf(:,is_ter) 565 pctsrf_t(:,is_lic,k) = pctsrf(:,is_lic) 566 DO i = 1, klon 567 pctsrf_t(i,is_sic,k) = (1. - pctsrf_t(i,is_lic,k) - 568 . pctsrf_t(i,is_ter,k)) * phy_ice(i) 569 pctsrf_t(i,is_oce,k) = 1. - pctsrf_t(i,is_lic,k) - 570 . pctsrf_t(i,is_ter,k) - pctsrf_t(i,is_sic,k) 571 if (pctsrf_t(i,is_oce,k) .lt. 0.) then 572 WRITE(*,*) 'pb sous maille au point : i,k ' 573 $ , i,k,pctsrf_t(:,is_oce,k) 574 ENDIF 575 END DO 576 ELSE 577 DO i = 1, klon 578 pctsrf_t(i,is_ter,k) = pctsrf(i,is_ter) 579 IF (NINT(pctsrf(i,is_ter)).EQ.1 ) THEN 580 pctsrf_t(i,is_sic,k) = 0. 581 pctsrf_t(i,is_oce,k) = 0. 582 IF(phy_ice(i) .GE. 1.e-5) THEN 583 pctsrf_t(i,is_lic,k) = phy_ice(i) 584 pctsrf_t(i,is_ter,k) = pctsrf_t(i,is_ter,k) 585 . - pctsrf_t(i,is_lic,k) 586 ELSE 587 pctsrf_t(i,is_lic,k) = 0. 588 ENDIF 589 ELSE 590 pctsrf_t(i,is_lic,k) = 0. 591 IF(phy_ice(i) .GE. 1.e-5) THEN 592 pctsrf_t(i,is_ter,k) = 0. 593 pctsrf_t(i,is_sic,k) = phy_ice(i) 594 pctsrf_t(i,is_oce,k) = 1. - pctsrf_t(i,is_sic,k) 595 ELSE 596 pctsrf_t(i,is_sic,k) = 0. 597 pctsrf_t(i,is_oce,k) = 1. 598 ENDIF 599 ENDIF 600 verif = pctsrf_t(i,is_ter,k) + 601 . pctsrf_t(i,is_oce,k) + 602 . pctsrf_t(i,is_sic,k) + 603 . pctsrf_t(i,is_lic,k) 604 IF ( verif .LT. 1. - 1.e-5 .OR. 605 $ verif .GT. 1 + 1.e-5) THEN 606 WRITE(*,*) 'pb sous maille au point : i,k,verif ' 515 607 $ , i,k,verif 516 ENDIF 517 END DO 518 ELSE 519 DO i = 1, klon 520 phy_nat(i,k) = phy_nat0(i) 521 IF (NINT(phy_nat0(i)).EQ.1 ) THEN 522 IF(phy_ice(i,k) .GE. 1.e-5) THEN 523 phy_icet(i,k) = phy_ice(i,k) 524 phy_ice(i,k) = 0. 525 phy_nat(i,k) = phy_nat(i,k) - phy_icet(i,k) 526 phy_oce(i,k) = 0. 527 ELSE 528 phy_icet(i,k) = 0. 529 phy_ice(i,k) = 0. 530 phy_oce(i,k) = 0. 531 ENDIF 532 ELSE 533 IF(phy_ice(i,k) .GE. 1.e-5) THEN 534 phy_icet(i,k) = 0. 535 phy_nat(i,k) = 0. 536 phy_oce(i,k) = 1. - phy_ice(i,k) 537 ELSE 538 phy_icet(i,k) = 0. 539 phy_ice(i,k) = 0. 540 phy_oce(i,k) = 1. 541 ENDIF 542 ENDIF 543 verif = phy_nat(i,k) + phy_icet(i,k)+ phy_ice(i,k) 544 $ + phy_oce(i,k) 545 IF ( verif .LT. 1. - 1.e-5 .OR. 546 $ verif .GT. 1 + 1.e-5) THEN 547 WRITE(*,*) 'pb sous maille au point : i,k,verif ' 548 $ , i,k,verif 549 ENDIF 550 END DO 551 ENDIF 608 ENDIF 609 END DO 610 ENDIF 552 611 ELSE 553 DO i = 1, klon 554 phy_nat(i,k) = phy_nat0(i) 555 IF ( (phy_ice(i,k) - 0.5).GE.1.e-5 ) THEN 556 IF (NINT(phy_nat0(i)).EQ.0) THEN 557 phy_nat(i,k) = 3.0 558 ELSE 559 phy_nat(i,k) = 2.0 560 ENDIF 612 DO i = 1, klon 613 phy_nat(i,k) = phy_nat0(i) 614 IF ( (phy_ice(i) - 0.5).GE.1.e-5 ) THEN 615 IF (NINT(phy_nat0(i)).EQ.0) THEN 616 phy_nat(i,k) = 3.0 617 ELSE 618 phy_nat(i,k) = 2.0 561 619 ENDIF 562 END DO 620 ENDIF 621 END DO 563 622 ENDIF 564 623 ENDDO … … 892 951 IF (newlmt ) THEN 893 952 ierr = NF_PUT_VARA_DOUBLE (nid,id_FOCE,debut,epais 894 $ ,p hy_oce(1,k))953 $ ,pctsrf_t(1,is_oce,k)) 895 954 ierr = NF_PUT_VARA_DOUBLE (nid,id_FSIC,debut,epais 896 $ ,p hy_ice(1,k))955 $ ,pctsrf_t(1,is_sic,k)) 897 956 ierr = NF_PUT_VARA_DOUBLE (nid,id_FTER,debut,epais 898 $ ,p hy_nat(1,k))957 $ ,pctsrf_t(1,is_ter,k)) 899 958 ierr = NF_PUT_VARA_DOUBLE (nid,id_FLIC,debut,epais 900 $ ,p hy_icet(1,k))959 $ ,pctsrf_t(1,is_lic,k)) 901 960 ELSE 902 961 ierr = NF_PUT_VARA_DOUBLE (nid,id_NAT,debut,epais … … 912 971 IF (newlmt ) THEN 913 972 ierr = NF_PUT_VARA_REAL (nid,id_FOCE,debut,epais 914 $ ,p hy_oce(1,k))973 $ ,pctsrf_t(1,is_oce,k)) 915 974 ierr = NF_PUT_VARA_REAL (nid,id_FSIC,debut,epais 916 $ ,p hy_ice(1,k))975 $ ,pctsrf_t(1,is_sic,k)) 917 976 ierr = NF_PUT_VARA_REAL (nid,id_FTER,debut,epais 918 $ ,p hy_nat(1,k))977 $ ,pctsrf_t(1,is_ter,k)) 919 978 ierr = NF_PUT_VARA_REAL (nid,id_FLIC,debut,epais 920 $ ,p hy_icet(1,k))979 $ ,pctsrf_t(1,is_lic,k)) 921 980 ELSE 922 981 ierr = NF_PUT_VARA_REAL (nid,id_NAT,debut,epais
Note: See TracChangeset
for help on using the changeset viewer.
