! -*- Mode: f90 -*- PROGRAM cote USE declare USE modeles USE mod_lecdata USE formula USE math USE mod_rectifie USE mod_lbc USE mod_norma USE bords USE mod_inter USE mod_prih USE mod_wri_wei USE mod_inipar USE fliocom USE getincom USE errioipsl !!! IMPLICIT NONE !!! !!! A partir d'un fichier de poids, complete avec les rivieres et !!! les poids pour le run off !! On met a zero les poids sur les points non cotiers, et on renormalise !! INTEGER (kind=il), PARAMETER :: jpr = 52 ! Nombre de rivieres REAL (kind=rl) :: zsuma, zsumo REAL (kind=rl), DIMENSION (:), ALLOCATABLE :: za REAL (kind=rl), DIMENSION (:), ALLOCATABLE :: zo !! CHARACTER (LEN = 1) :: clriv CHARACTER (LEN = 30) :: clname INTEGER (kind=il) :: ja, jo, jo2, jo3, jn, jno, joi, joj, jai, jaj, jr, kn ! indices de boucle INTEGER (kind=il) :: jai_p, jai_m, jaj_p, jaj_m, n1 INTEGER (kind=il) :: joi_p, joi_m, joj_p, joj_m REAL (kind=rl) :: z_d_1, z_d_2, zzmin INTEGER (kind=il) :: nriv = 32, nlmd = 35 INTEGER (kind=il), DIMENSION (:, :), ALLOCATABLE :: jx ! Point trouve sous une maille atm INTEGER (kind=il) :: kom, ko LOGICAL, DIMENSION (:), ALLOCATABLE :: lacot, laland, laoce, lanoroute LOGICAL, DIMENSION (:), ALLOCATABLE :: locot, looce INTEGER (kind=il), DIMENSION (:), ALLOCATABLE :: ktemp REAL(kind=rl), DIMENSION (:), ALLOCATABLE :: wtemp, d_cote INTEGER (kind=il), DIMENSION (:,:), ALLOCATABLE :: ktemp2 REAL(kind=rl), DIMENSION (:,:), ALLOCATABLE :: wtemp2 INTEGER (kind=il), DIMENSION (:), ALLOCATABLE :: iamsk INTEGER (kind=il), DIMENSION (:), ALLOCATABLE :: iomsk INTEGER (kind=il) :: jovoid, jot, inum CHARACTER (len=80) :: cfname4, cfname8 LOGICAL :: ltrouve, ltrouve2 INTEGER (kind=il), DIMENSION (:), ALLOCATABLE :: itarget REAL (kind=rl), DIMENSION (:), ALLOCATABLE :: wtarget INTEGER (kind=il), DIMENSION (:,:), ALLOCATABLE :: i2a INTEGER (kind=il) :: no, ntrouve INTEGER (kind=il), DIMENSION (:), ALLOCATABLE :: jo_n LOGICAL :: l_fulldiag = .TRUE. CHARACTER (len=180) :: c_date, c_time, c_zone, c_tmp CHARACTER (len=80) :: clarg INTEGER :: narg REAL(kind=rl), DIMENSION (:), ALLOCATABLE :: tmask_atl, tmask_noclose, tmask_nomed, tmask_med REAL(kind=rl), DIMENSION (:, :), ALLOCATABLE :: z2o CHARACTER (len=40) :: cl_bassin INTEGER (kind=il) :: il_ncid !! INTEGER (kind=il), DIMENSION (2) :: jind INTEGER (kind=il) :: ierr !! !! Recherche voisins proches INTEGER (kind=il) :: jp_nb INTEGER (kind=il), DIMENSION (:), ALLOCATABLE :: jo_tab INTEGER (kind=il) :: jo_nb !! INTEGER :: nid_it !! CHARACTER (len = 50) :: cldiag_a2o, clw_a2o, clw_a2o_cdf, clw_a2o_mct CHARACTER (len = 25) :: mo_name, a2o_name CHARACTER (len = 1) :: c_i, c_r !! Lecture de la ligne de comande narg = iargc() IF ( narg > 0 ) THEN CALL getarg ( 1, clarg) SELECT CASE (TRIM(clarg)) 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 !! CALL inipar CALL alloc_modeles !! ALLOCATE (za (jpan) , STAT=ierr) ; CALL chk_allo (ierr, 'za', lreset=.TRUE., crout='cotes') ALLOCATE (zo (jpon) , STAT=ierr) ; CALL chk_allo (ierr, 'zo') ALLOCATE (jx (jpo2a,2) , STAT=ierr) ; CALL chk_allo (ierr, 'jx') ALLOCATE (lacot (jpan) , STAT=ierr) ; CALL chk_allo (ierr, 'lacot') ; lacot = .FALSE. ALLOCATE (laland (jpan) , STAT=ierr) ; CALL chk_allo (ierr, 'laland') ; laland = .FALSE. ALLOCATE (laoce (jpan) , STAT=ierr) ; CALL chk_allo (ierr, 'laoce') ; laoce = .FALSE. ! ALLOCATE (locot (jpon) , STAT=ierr) ; CALL chk_allo (ierr, 'locot') ; locot = .FALSE. ALLOCATE (looce (jpon) , STAT=ierr) ; CALL chk_allo (ierr, 'looce') ; looce = .FALSE. ! ALLOCATE (ktemp (jpa2o) , STAT=ierr) ; CALL chk_allo (ierr, 'ktemp') ALLOCATE (wtemp (jpa2o) , STAT=ierr) ; CALL chk_allo (ierr, 'wtemp') ALLOCATE (iamsk (jpan) , STAT=ierr) ; CALL chk_allo (ierr, 'iamsk') ALLOCATE (iomsk (jpon) , STAT=ierr) ; CALL chk_allo (ierr, 'iomsk') ALLOCATE (d_cote (jpon) , STAT=ierr) ; CALL chk_allo (ierr, 'd_cote') ! ALLOCATE (itarget (jpan) , STAT=ierr) ; CALL chk_allo (ierr, 'itarget') ALLOCATE (wtarget (jpan) , STAT=ierr) ; CALL chk_allo (ierr, 'wtarget') ALLOCATE (i2a (jpai, jpaj), STAT=ierr) ; CALL chk_allo (ierr, 'i2a') ! ALLOCATE (tmask_atl (jpon), STAT=ierr) ; CALL chk_allo (ierr, 'tmask_atl') ALLOCATE (tmask_noclose (jpon), STAT=ierr) ; CALL chk_allo (ierr, 'tmask_noclose') ALLOCATE (tmask_nomed (jpon), STAT=ierr) ; CALL chk_allo (ierr, 'tmasl_nomed') ALLOCATE (tmask_med (jpon), STAT=ierr) ; CALL chk_allo (ierr, 'tmask_med') ALLOCATE (z2o (jpoi, jpoj), STAT=ierr) ; CALL chk_allo (ierr, 'z2o') ! jp_nb = jpon ALLOCATE (jo_tab (jpon)) ; CALL chk_allo (ierr, 'jo_tab') ! IF (l_limit_iosize) l_fulldiag = .FALSE. !! !! Initialisation ngrd = 11 ; nsrf = 12 ; nmsk = 13 ; nchk = 14 ; ndeb = 9 ; nbug1 = 10 nwei4o2a = 81 ; nwei8o2a = 82 ; nwei4a2o = 83 ; nwei8a2o = 84 ! nout = 6 IF (nout /= 6 .AND. nout /= 0 ) & & OPEN (unit = nout, file = 'cote' // TRIM(c_suffix) // '.out', action = 'WRITE', status = 'unknown', position = 'REWIND' ) !! !! Diagnostics files !! IF (nout /= 6 ) & & OPEN (unit = nout, file = 'cotes' // TRIM(c_suffix) // '.out', position = 'REWIND', action = 'WRITE', status = 'replace' ) !! !! Read coordinates of all models, and weights !! CALL lecdata (lecweights = .TRUE., leco2amask = .TRUE. ) !! Lecture masques de bassin IF ( TRIM (cotyp) == 'orca4' .OR. & & TRIM (cotyp) == 'orca2' .OR. TRIM (cotyp) == 'orca2.0' .OR. TRIM (cotyp) == 'orca2.1' .OR. TRIM (cotyp) == 'orca2.3' .OR. & & TRIM (cotyp) == 'orca1' & & ) THEN cl_bassin = TRIM (cotyp) !-$$ IF (c_period /= 'none' ) cl_bassin = TRIM (cl_bassin) // TRIM (c_period) cl_bassin = TRIM (cl_bassin) // '.nc' WRITE (nout,*) 'Reading bassin masks in : ', TRIM (cl_bassin) CALL flioopfd (TRIM (cl_bassin), il_ncid ) CALL fliogetv (il_ncid, 'mask_atl' , z2o ) tmask_atl = RESHAPE (z2o, (/ jpon /) ) CALL fliogetv (il_ncid, 'mask_noclose', z2o ) tmask_noclose = RESHAPE (z2o, (/ jpon /) ) CALL fliogetv (il_ncid, 'mask_nomed' , z2o ) tmask_nomed = RESHAPE (z2o, (/ jpon /) ) CALL flioclo (il_ncid) ELSE tmask_atl = REAL (1-iomskt, rl) tmask_noclose = REAL (1-iomskt, rl) tmask_nomed = REAL (1-iomskt, rl) END IF ! WRITE (unit=nout, fmt=*) 'iomskt' CALL prihin (RESHAPE(iomskt ,(/jpoi,jpoj/)), 1_il, nout) WRITE (unit=nout, fmt=*) 'tmask_atl' CALL prihin (RESHAPE(NINT(tmask_atl) ,(/jpoi,jpoj/)), 1_il, nout) WRITE (unit=nout, fmt=*) 'tmask_med' !CALL prihin (RESHAPE(NINT(tmask_med) ,(/jpoi,jpoj/)), 1_il, nout) WRITE (unit=nout, fmt=*) 'tmask_noclose' !CALL prihin (RESHAPE(NINT(tmask_noclose),(/jpoi,jpoj/)), 1_il, nout) ! ! Set correct periodicity to grids CALL lbc (xolont, noperio, 'T', jpoi, jpoj) CALL lbc (xolonu, noperio, 'U', jpoi, jpoj) CALL lbc (xolonv, noperio, 'V', jpoi, jpoj) CALL lbc (xolonf, noperio, 'F', jpoi, jpoj) CALL lbc (xolatt, noperio, 'T', jpoi, jpoj) CALL lbc (xolatu, noperio, 'U', jpoi, jpoj) CALL lbc (xolatv, noperio, 'V', jpoi, jpoj) CALL lbc (xolatf, noperio, 'F', jpoi, jpoj) CALL lbc (xosrft, noperio, 'T', jpoi, jpoj) CALL lbc (iomskt, noperio, 'T', jpoi, jpoj) ! CALL lbc (xalont, naperio, 'T', jpai, jpaj) CALL lbc (xalonu, naperio, 'U', jpai, jpaj) CALL lbc (xalonv, naperio, 'V', jpai, jpaj) CALL lbc (xalatt, naperio, 'T', jpai, jpaj) CALL lbc (xalatu, naperio, 'U', jpai, jpaj) CALL lbc (xalatv, naperio, 'V', jpai, jpaj) CALL lbc (xasrft, naperio, 'T', jpai, jpaj) CALL lbc (iamskt, naperio, 'T', jpai, jpaj) !! CALL rectifie !! CALL fliocrfd ('itarget' // TRIM(c_suffix), (/'x', 'y'/), (/jpai, jpaj/), nid_it, mode='REPLACE') CALL flioputa (nid_it, "?", "title" , "itarget" ) CALL flioputa (nid_it, '?', 'Comment', TRIM(c_comment) ) CALL DATE_AND_TIME (c_date, c_time, c_zone ) CALL flioputa (nid_it, "?", "history" , "Created: "//c_date(1:4)//"-"//c_date(5:6)//"-"// & & c_date(7:8)//" "//c_time(1:2)//"h"//c_time(3:4)//" GMT"//TRIM(c_zone) ) CALL flioputa (nid_it, "?", "method" , "MOSAIC" ) CALL flioputa (nid_it, "?", "source_grid" , "curvilinear" ) CALL flioputa (nid_it, "?", "dest_grid" , "curvilinear" ) CALL flioputa (nid_it, "?", "Institution" , "IPSL" ) CALL flioputa (nid_it, "?", "Model" , "IPSL CM6" ) CALL GET_ENVIRONMENT_VARIABLE ( NAME="HOSTNAME", VALUE=c_tmp, TRIM_NAME=.TRUE., STATUS=ierr) IF ( ierr == 0 ) THEN CALL flioputa (nid_it, "?", "HOSTNAME" , TRIM(c_tmp) ) ELSE WRITE (nout,*) 'Environment variable not found : $HOSTNAME' END IF CALL GET_ENVIRONMENT_VARIABLE ( NAME="LOGNAME" , VALUE=c_tmp, TRIM_NAME=.TRUE., STATUS=ierr) IF ( ierr == 0 ) THEN CALL flioputa (nid_it, "?", "LOGNAME" , TRIM(c_tmp) ) ELSE WRITE (nout,*) 'Environment variable not found : $LOGNAME' END IF CALL fliodefv (nid_it, 'lon', (/1, 2/), units="degree_north", v_t=flio_r4 ) CALL fliodefv (nid_it, 'lat', (/1, 2/), units="degree_east", v_t=flio_r4 ) CALL fliodefv (nid_it, 'o2amask', (/1,2/), v_t=flio_r4) CALL fliodefv (nid_it, 'laland', (/1,2/), v_t=flio_i4) CALL fliodefv (nid_it, 'lacot' , (/1,2/), v_t=flio_i4) CALL fliodefv (nid_it, 'laoce' , (/1,2/), v_t=flio_i4) CALL fliodefv (nid_it, 'iamskp', (/1,2/), v_t=flio_i4) CALL flioputv (nid_it, 'lon' , RESHAPE ( xalont, (/jpai, jpaj/)) ) CALL flioputv (nid_it, 'lat' , RESHAPE ( xalatt, (/jpai, jpaj/)) ) CALL flioputv (nid_it, 'o2amask', RESHAPE (o2amask, (/jpai, jpaj/)) ) CALL flioputv (nid_it, 'iamskp' , RESHAPE (iamskp , (/jpai, jpaj/)) ) !! !! Compute edges of grid boxes !! CALL bord_a CALL bord_o CALL bord_oa !! mo_name = "mozaic" a2o_name = "a2o" !! WRITE (nout,*) 'o2amask (>0, <=0, >1) ', COUNT (o2amask > 0.0_rl), COUNT (o2amask <= 0.0_rl) & & , COUNT (o2amask >= 1.0_rl) i2a = 0_il WHERE (RESHAPE(o2amask,(/jpai,jpaj/)) .GT. 0.0_rl ) i2a = 1_il END WHERE WRITE (nout, *) 'i2a (o2amask > 0)' CALL prihin (i2a , 1_il, nout) !! !! Verif poids avant reduction !! WRITE (unit=nout,fmt=*) 'Verif poids globale' za = 0.0_rl WRITE (nout, *) 'Comptage xasrft ' , COUNT ( xasrft < 1.0E-10_rl ) WRITE (nout, *) 'Comptage o2amask ' , COUNT ( o2amask < 0.5_rl ) WHERE ( xasrft * o2amask /= 0.0_rl ) za = 1.0_rl / (xasrft * o2amask) END WHERE WRITE (nout, *) 'Somme za ', SUM ( za) WRITE (nout, *) 'Somme xasrft ', SUM ( xasrft) zsuma = DOT_PRODUCT (za, xasrft * o2amask) !! WRITE (nout, *) 'jma2o for coastal runoff : ', jma2o !! Atm mask interpolated to oce za = REAL (1_il-iamskt, KIND = rl ) CALL inter_a2o (za, a2omask ) !! 1 on Atm interpolated to oce za = REAL (1_il , KIND = rl ) CALL inter_a2o (za, a2ofull ) !! CALL verif ('Apres lecture') !! CALL bilan ("Poids ", c_case="complet", c_int_atm="Local", c_int_oce="Local") CALL bilan ("Poids ", c_case="ocean" , c_int_atm="Local", c_int_oce="Local") !! !! Remise des poids en integre DO jn = 1, jpa2o DO jo = 1, jpon ja = ka2o (jn, jo) IF (ja > 0) THEN IF ( o2amask (ja) .GT. 0.0_rl) & & wa2o (jn, jo) = wa2o (jn, jo) * xosrft (jo) / ( xasrft (ja) * o2amask (ja)) END IF END DO END DO ! CALL bilan ("Integre ", c_case="complet", c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Integre ", c_case="ocean" , c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Integre ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") !! E_DRYRUN: IF (l_dryrun) THEN WRITE (nout,*) "Remplissage bidon des champs pour tester les I/O" jma2o = jma2or nvo = jma2or WRITE (nout,*) 'jma2o pour test : ', jma2o CALL RANDOM_SEED CALL RANDOM_NUMBER (wa2o) ; ka2o = NINT (wa2o*REAL(jpan,rl) ) DO jo = 1, jpon IF ( iomskt (jo) == 1 ) THEN wa2o (:, jo) = 0.0_rl ka2o (:, jo) = 0_il END IF END DO DO jo = 1, jpon DO jno = 1, jpa2o ja = ka2o (jno, jo) IF (ja > 0) THEN IF (iamskt (ja) == 1) THEN wa2o (jno, jo) = 0.0_rl ; ka2o (jno, jo) = 0_il END IF END IF END DO END DO ELSE L_lcoast: IF (lcoast) THEN !! Dataset 4 : idem 2, mais en interpolant seulement des points cotes atmosphere vers !! les points cotes ocean. WRITE (nout,*) 'Demarrage L_coast' !! WRITE (nout,*) 'Calculs points terre/oce/cote sur grille atmosphere' WHERE (o2amask .LT. epsfrac ) laland = .TRUE. WHERE (o2amask .GE. epsfrac .AND. o2amask .LE. (1.0_rl-epsfrac) ) lacot = .TRUE. WHERE ( o2amask .GT. (1.0_rl-epsfrac) ) laoce = .TRUE. WRITE (UNIT=nout, FMT=*) 'Points cotiers ', COUNT(lacot) WRITE (nout,*) 'Ecriture ' i2a = 0 WHERE (RESHAPE(laland,(/jpai,jpaj/))) i2a=1 CALL flioputv (nid_it, "laland", i2a) i2a = 0 WHERE (RESHAPE(laoce ,(/jpai,jpaj/))) i2a=1 CALL flioputv (nid_it, "laoce", i2a) i2a = 0 WHERE (RESHAPE(lacot ,(/jpai,jpaj/))) i2a=1 CALL flioputv (nid_it, "lacot", i2a) WRITE (nout,*) 'Determination des points cotes ocean ' ALLOCATE (ktemp2(jpoi,jpoj)) ktemp2 = 0_il locot = .FALSE. looce = .FALSE. IF ( lcoast ) THEN DO jo = 1, jpon joi = moi (jo) ; joj = moj (jo) ! Indices 2D joi_p = joi+1 ; joi_m = joi-1 joj_p = MIN(joj+1, jpoj) ; joj_m = MAX (joj-1, 1) IF ( iomskt (jo) == 0_il .AND. iomskp (jo) == 0_il ) THEN IF ( iomskt (m1o (joi_p, joj )) == 1 & .OR. iomskt (m1o (joi_m, joj )) == 1 & .OR. iomskt (m1o (joi , joj_p)) == 1 & .OR. iomskt (m1o (joi , joj_m)) == 1 & .OR. iomskt (m1o (joi_p, joj_p)) == 1 & .OR. iomskt (m1o (joi_p, joj_m)) == 1 & .OR. iomskt (m1o (joi_m, joj_p)) == 1 & .OR. iomskt (m1o (joi_m, joj_m)) == 1 ) THEN locot (jo) = .TRUE. ktemp2 (joi, joj) = 1_il ELSE looce (jo) = .TRUE. END IF END IF END DO ELSE locot (:) = ( iomskt == 0 ) ENDIF CALL lbc (locot, noperio, 'T', jpoi, jpoj) CALL lbc (looce, noperio, 'T', jpoi, jpoj) WRITE (nout,*) "iomskt" CALL prihin (RESHAPE(iomskt,(/jpoi,jpoj/)), 1, nout) WRITE (nout,*) "Locot" CALL prihin (ktemp2, 1, nout) DEALLOCATE (ktemp2) WRITE (nout,*) 'Verif coherence (tout doit etre nul)' WRITE (nout,*) COUNT ( laland .AND. lacot) WRITE (nout,*) COUNT ( lacot .AND. laoce) WRITE (nout,*) COUNT ( laoce .AND. laland) WRITE (nout,*) COUNT ( .NOT. ( laland .OR. lacot .OR. laoce) ) !-$$ iamsk = 0 ; WHERE (lacot) iamsk = 1 !-$$ WRITE(nout,*) 'Fraction lacot' !-$$ CALL prihin(RESHAPE(iamsk,(/jpai,jpaj/)),1_il,nout) !-$$ !-$$ iamsk = 0 ; WHERE (laland) iamsk = 1 !-$$ WRITE(nout,*) 'Points avec terre laland' !-$$ CALL prihin(RESHAPE(iamsk,(/jpai,jpaj/)),1_il,nout) !-$$ !-$$ iamsk = 0 ; WHERE (laoce) iamsk = 1 !-$$ WRITE(nout,*) 'Fraction laoce' !-$$ CALL prihin(RESHAPE(iamsk,(/jpai,jpaj/)),1_il,nout) !-$$ !-$$ WRITE(nout,*) 'Periodicite atm' !-$$ CALL prihin(RESHAPE(iamskp,(/jpai,jpaj/)),1_il,nout) !! WRITE (nout,*) "iamskt" CALL prihin (RESHAPE(iamskt,(/jpai,jpaj/)), 1, nout) WRITE (nout,*) "Lacot" ALLOCATE (ktemp2 (jpai, jpaj)) ktemp2 = 0_il ; WHERE (RESHAPE(lacot, (/jpai, jpaj/))) ktemp2 = 1_il CALL prihin (ktemp2, 1, nout) !! WRITE (unit = nout, fmt = *) 'Nombre de points cotes atmosphere : ', COUNT (lacot) WRITE (unit = nout, fmt = *) 'Nombre de points cotes ocean : ', COUNT (locot) !! CALL bilan ("Avec cote ", c_case="complet", c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Avec cote ", c_case="ocean" , c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Avec cote ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") ! ! Mise a zero des poids sur les points non cotiers DO jo = 1, jpon IF (.NOT. locot (jo) ) THEN wa2o (:, jo) = 0.0e0_rl ka2o (:, jo) = 0_il ELSE DO jn = 1, jpa2o ja = ka2o (jn, jo) IF (ja > 0 ) THEN IF (.NOT. lacot (ja) ) THEN wa2o (jn, jo) = 0.0e0_rl ka2o (jn, jo) = 0_il END IF END IF END DO END IF END DO CALL verif ('Apres remise a zero') !! !! Elimination des points redondant par periodicite DO jo = 1, jpon IF (iomskp (jo) .EQ. 1_il ) THEN wa2o (:,jo) = 0.0e0_rl ka2o (:,jo) = 0_il END IF END DO !! !! Controle !! CALL verif ('Apres elimination de periodicite') !! CALL bilan ("Apres remise a zero ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") !! !! Renormalize ! Cas 1 RenormA : DO ja = 1, jpan IF (.NOT. lacot (ja) ) CYCLE RenormA z_d_1 = 0.0e0_rl kom = 0 RenormO : DO jo = 1, jpon IF (.NOT. locot (jo) ) CYCLE RenormO DO jn = 1, jpa2o IF ( ka2o (jn,jo) == ja ) THEN z_d_1 = z_d_1 + wa2o (jn,jo) kom = kom + 1 jx (kom,:) = (/ jn, jo /) !!$ WRITE (unit = nout, fmt = '(6I6,16F18.6) ') ja, jo, jn, kom, jx (kom, :), & !!$ wa2o (jn,jo), wa2o (jx (kom, 1), jx (kom, 2)), z_d_1, & !!$ xosrft (jo), xasrft (ja) END IF END DO ENDDO RenormO IF (kom /= 0 ) THEN IF (z_d_1 == 0.0_rl ) THEN WRITE (nout, *) 'z_d_1 nul ', ja, xasrft (ja), lacot (ja), iamskt (ja), & & wa2o (:,ja), ka2o (:,ja) STOP ELSE DO ko = 1, kom wa2o (jx (ko,1), jx (ko,2)) = wa2o (jx (ko,1), jx (ko,2)) / z_d_1 END DO ENDIF END IF ENDDO RenormA !! !! Calcul des distances à la cote des points oceans !-$$ WRITE (nout,*) 'Distance cotes' !-$$ d_cote = rpi !-$$ Cal_dist_1: DO jo = 1, jpon !-$$ IF ( iomskt (jo) == 1 .OR. locot (jo) ) THEN !-$$ d_cote (jo) = 0.0E0_rl !-$$ ELSE !-$$ Cal_dist_2:DO jo2 = 1, jpon !-$$ IF ( locot (jo2) ) THEN !-$$ d_cote (jo) = MIN (d_cote (jo), REAL (geodist(xolont(jo), xolatt(jo), xolont(jo2), xolatt(jo2)), KIND=rl)) !-$$ END IF !-$$ END DO Cal_dist_2 !-$$ END IF !-$$ END DO Cal_dist_1 !-$$ d_cote = d_cote * r_earth !CALL prihre (RESHAPE(d_cote,(/jpoi,jpoj/)), 1.0E-5_rl, nout) WRITE (nout,*) 'Distance maxi a la cote (ocean) : ', MAXVAL (d_cote) CALL fliocrfd ("dist_cote" // TRIM(c_suffix), (/ "x", "y"/), (/jpoi, jpoj/), n1, MODE="REPLACE") CALL flioputa (n1, '?', 'Comment', TRIM(c_comment) ) CALL fliodefv (n1, "dist_cote", (/1, 2/)) CALL flioputv (n1, "dist_cote", RESHAPE (d_cote, (/jpoi, jpoj/) )) CALL flioclo (n1) ! !! Controle !! CALL verif ('Premiere normalisation') CALL gather_wei ('Premiere etape') CALL verif ('Premier gather') !! CALL bilan ("Renorm ", c_case="complet", c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Renorm ", c_case="ocean" , c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Renorm ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") !! !! Traitement des points atmospheres non atteint !! total: IF (ltotal) THEN WRITE (nout,*) 'Va chercher tout les points terre (ltotal=TRUE)' SeekAtm_1_1: DO ja = 1, jpan !-$$ WRITE(nout,*) 'Point atm : ', ja, m2ai(ja), m2aj(ja), itarget(ja), iamskp(ja), laland(ja), lacot(ja), laoce(ja) IF (itarget (ja) .GE. 1_il) CYCLE SeekAtm_1_1 IF (iamskp (ja) .EQ. 1_il) CYCLE SeekAtm_1_1 IF (.NOT. laland (ja) ) CYCLE SeekAtm_1_1 zzmin = REAL (cartdist(0.0_rl, 90.0_rl, 0.0_rl, -90.0_rl), KIND=rl) ltrouve = .FALSE. ; jot = 0_il !-$$ WRITE(nout,*) 'Point atm : ', ja, m2ai(ja), m2aj(ja) SeekOce_1_1: DO jo = 1, jpon ! Cherche point ocean le plus proche IF (.NOT. locot (jo) ) CYCLE SeekOce_1_1 IF (iomskp(jo) .EQ. 1_il) CYCLE SeekOce_1_1 z_d_1 = REAL (cartdist(xalont(ja), xalatt(ja), xolont(jo), xolatt(jo)), KIND=rl) IF (z_d_1 .LT. zzmin) THEN zzmin = z_d_1 ltrouve = .TRUE. ; jot = jo END IF ENDDO SeekOce_1_1 Trouve_1_1: IF (ltrouve) THEN ! Calcul du poids z_d_1 = 1.0_rl ! Rangement des poids IF (nvo (jot) < jpa2o ) THEN ! Reste un poids libre nvo (jot) = nvo (jot) + 1 wa2o (nvo (jot), jot) = z_d_1 ka2o (nvo (jot), jot) = ja !-$$ WRITE (nout,*) ' Vers point oce ', jot, m2oi(jot), m2oj(jot), nvo(jot) ELSE ! Trouve12 WRITE(nout,*) ' Plus de place ', ja, jot, nvo (jot), z_d_1 STOP ENDIF ELSE ! Trouve !! Message d'erreur WRITE (unit = nout, fmt = '("Pas de voisin proche pour ", 1I5,2I4)' ) ja, m2ai(ja), m2aj(ja) END IF Trouve_1_1 END DO SeekAtm_1_1 !! CALL verif ('ltotal avant gather') CALL bilan ("Ltotal avant gather ", c_case="terre" , c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Ltotal avant gather ", c_case="terre100", c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Ltotal avant gather ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") !! CALL gather_wei ('apres total') CALL verif ( 'Cas total', l_non_a=.TRUE.) !! WRITE (unit = nout, fmt = * ) 'Nombre de voisins apres completion : ', jma2o !! CALL bilan ("Fin ltotal ", c_case="terre" , c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Fin ltotal ", c_case="terre100", c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Fin ltotal ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") !! END IF Total !! Total_dist: IF (ltotal_dist) THEN !! On limite la recherche à une certaine distance de la cote. WRITE (nout,*) 'Cas ltotal_dist' SeekAtm_2_1: DO ja = 1, jpan IF (itarget (ja) .GE. 1_il) CYCLE SeekAtm_2_1 IF (iamskp (ja) .EQ. 1_il) CYCLE SeekAtm_2_1 IF (.NOT. laland (ja) ) CYCLE SeekAtm_2_1 zzmin = dist_max ltrouve = .FALSE. ; jot = 0_il !$$$ WRITE(nout,*) 'Point atm ', ja, m2ai(ja), m2aj(ja) SeekOce_2_1: DO jo = 1, jpon ! Cherche point ocean le plus proche IF (.NOT. locot (jo) ) CYCLE SeekOce_2_1 IF (iomskp(jo) .EQ. 1_il) CYCLE SeekOce_2_1 z_d_1 = r_earth * REAL (geodist(xalont(ja), xalatt(ja), xolont(jo), xolatt(jo)), KIND=rl) IF (z_d_1 .LT. zzmin) THEN zzmin = z_d_1 ltrouve = .TRUE. ; jot = jo END IF ENDDO SeekOce_2_1 Trouve_2_1: IF (ltrouve) THEN ! Calcul du poids z_d_1 = 1.0_rl IF (nvo (jot) < jpa2o ) THEN ! Reste un poids libre nvo (jot) = nvo (jot) + 1 wa2o (nvo (jot), jot) = z_d_1 ka2o (nvo (jot), jot) = ja !$$$ WRITE (nout,*) ' Vers point oce ', jot, m2oi(jot), m2oj(jot), nvo(jot) ELSE WRITE(nout,*) ' Plus de place ', ja, jot, nvo (jot), z_d_1 STOP ENDIF END IF Trouve_2_1 END DO SeekAtm_2_1 !! CALL verif ('Ltotal_dist avant gather') CALL bilan ("Ltotal_dist avant gather ", c_case="terre" , c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Ltotal_dist avant gather ", c_case="terre100", c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Ltotal_dist avant gather ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") !! CALL gather_wei ('apres total_dist') CALL verif ( 'Cas total_dist') !! WRITE (unit = nout, fmt = * ) 'Nombre de voisins apres completion : ', jma2o !! CALL bilan ("Fin ltotal_dist ", c_case="terre" , c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Fin ltotal_dist ", c_case="terre100", c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Fin ltotal_dist ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") !! END IF Total_dist !! Total_dist_2: IF (ltotal_dist_2) THEN !! On limite la recherche à une certaine distance de la cote. !! On etale sur les point oceans cotiers proches WRITE (nout,*) '-------------------' WRITE (nout,*) 'Cas ltotal_dist_2' WRITE (nout,*) 'Distance recherche : ', dist_max/1.0E3_rl, 'km' WRITE (nout,*) 'Distance etalement ocean : ', dist_max_voisin/1.0E3_rl, 'km' !! SeekAtm_3_1: DO ja = 1, jpan IF (itarget (ja) .GE. 1_il) CYCLE SeekAtm_3_1 IF (iamskp (ja) .EQ. 1_il) CYCLE SeekAtm_3_1 IF (.NOT. laland (ja) ) CYCLE SeekAtm_3_1 zzmin = dist_max ltrouve = .FALSE. ; jot = 0_il ; ntrouve = 0 !WRITE(nout,*) 'Point atm : ', ja, m2ai(ja), m2aj(ja) SeekOce_3_1: DO jo = 1, jpon ! Cherche point ocean le plus proche IF (.NOT. locot (jo) ) CYCLE SeekOce_3_1 IF (iomskp (jo) .EQ. 1_il) CYCLE SeekOce_3_1 z_d_1 = r_earth * REAL (geodist (xalont (ja), xalatt (ja), xolont (jo), xolatt (jo)), KIND=rl) IF (z_d_1 .LT. zzmin) THEN zzmin = z_d_1 ; ltrouve = .TRUE. ; jot = jo ; ntrouve = ntrouve+1 END IF ENDDO SeekOce_3_1 WRITE(nout,*) 'Point atm ', ja, m2ai(ja), m2aj(ja), jot, ntrouve Trouve_3_1: IF (ltrouve) THEN ! On cherche les points proches du point ocean trouve z_d_1 = 0.0_rl ; jo_tab = 0 ; jo_nb = 0 SeekOce_3_2: DO jo2 = 1, jpon IF (.NOT. locot (jo2) ) CYCLE SeekOce_3_2 IF (iomskp (jo2) .EQ. 1_il) CYCLE SeekOce_3_2 IF ( r_earth * REAL (geodist (xolont (jot), xolatt (jot), xolont (jo2), xolatt (jo2)), KIND=rl) & & < dist_max_voisin) THEN IF (jo_nb >= jp_nb) THEN WRITE (nout,*) 'Trop de point proches ', ja, jot, jo_tab STOP ELSE jo_nb = jo_nb + 1 z_d_1 = z_d_1 + xosrft (jo2) jo_tab (jo_nb) = jo2 END IF END IF ! Calcul du poids END DO SeekOce_3_2 IF (jo_nb == 0) THEN WRITE (nout,*) 'Pas de voisin au voisin ! ', ja, jot, m2oi(jot), m2oj(jot), xolont(jot), xolatt(jot) WRITE (nout,*) 'Point atm ', ja, m2ai(ja), m2aj(ja) WRITE (nout,*) r_earth * REAL (geodist (xolont (jot), xolatt (jot), xolont (jot), xolatt (jot)), KIND=rl) STOP ELSE !-$$ WRITE(nout,*) 'Point atm ', ja, m2ai(ja), m2aj(ja), jot, jo_nb, jo_tab (1:jo_nb) DO jn = 1, jo_nb jo3 = jo_tab (jn) IF (nvo (jo3) < jpa2o ) THEN ! Reste un poids libre nvo (jo3) = nvo (jo3) + 1 wa2o (nvo (jo3), jo3) = xosrft (jo3) / z_d_1 ka2o (nvo (jo3), jo3) = ja !$$$ WRITE (nout,*) ' Vers point oce ', jot, m2oi(jot), m2oj(jot), nvo(jot) ELSE WRITE(nout,*) ' Plus de place ', ja, jo3, nvo (jo3), z_d_1 STOP ENDIF END DO END IF !-$$ ELSE !-$$ WRITE (nout,*) 'Pas de voisins proche trouve ' !-$$ WRITE (nout,*) 'Point atm ', ja, m2ai(ja), m2aj(ja) !-$$ WRITE (nout,*) dist_max, zzmin, z_d_1 !-$$ !STOP END IF Trouve_3_1 END DO SeekAtm_3_1 !! CALL verif ('Ltotal_dist_2 avant gather') CALL bilan ("Ltotal_dist_2 avant gather ", c_case="terre" , c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Ltotal_dist_2 avant gather ", c_case="terre100", c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Ltotal_dist_2 avant gather ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") !! CALL gather_wei ('apres total_dist_2') CALL verif ( 'Cas total_dist_2') !! WRITE (unit = nout, fmt = * ) 'Nombre de voisins apres completion : ', jma2o !! CALL bilan ("Fin ltotal_dist_2 ", c_case="terre" , c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Fin ltotal_dist_2 ", c_case="terre100", c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Fin ltotal_dist_2 ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") !! END IF Total_dist_2 !! Total_dist_3: IF (ltotal_dist_3) THEN !! On limite la recherche à une certaine distance de la cote. !! On etale sur les point oceans proches, cote ou pas WRITE (nout,*) '-------------------' WRITE (nout,*) 'Cas ltotal_dist_3' WRITE (nout,*) 'Distance recherche : ', dist_max/1.0E3_rl, 'km' WRITE (nout,*) 'Distance etalement ocean : ', dist_max_voisin/1.0E3_rl, 'km' !! SeekAtm_4_1: DO ja = 1, jpan IF (itarget (ja) .GE. 1_il) CYCLE SeekAtm_4_1 IF (iamskp (ja) .EQ. 1_il) CYCLE SeekAtm_4_1 IF (.NOT. laland (ja) ) CYCLE SeekAtm_4_1 zzmin = dist_max ltrouve = .FALSE. ; jot = 0_il !-$$ WRITE(nout,*) 'Point atm : ', ja, m2ai(ja), m2aj(ja) SeekOce_4_1: DO jo = 1, jpon ! Cherche point ocean le plus proche IF (iomskp (jo) .EQ. 1_il .OR. iomskt (jo) .EQ. 1_il ) CYCLE SeekOce_4_1 z_d_1 = r_earth * REAL (geodist (xalont (ja), xalatt (ja), xolont (jo), xolatt (jo)), KIND=rl) IF (z_d_1 .LT. zzmin) THEN zzmin = z_d_1 ; ltrouve = .TRUE. ; jot = jo END IF ENDDO SeekOce_4_1 !-$$ WRITE(nout,*) 'Point atm ', ja, m2ai(ja), m2aj(ja), jot Trouve_4_1: IF (ltrouve) THEN ! On cherche les points proches du point ocean trouve, le long de la cote et vers le large z_d_1 = 0.0_rl ; jo_tab = 0 ; jo_nb = 0 SeekOce_4_2: DO jo2 = 1, jpon IF (iomskp (jo2) .EQ. 1_il .OR. iomskt (jo2) .EQ. 1_il) CYCLE SeekOce_4_2 IF ( .NOT. l_large .AND. .NOT. locot (jo2)) CYCLE SeekOce_4_2 z_d_2 = r_earth * REAL (geodist (xolont (jot), xolatt (jot), xolont (jo2), xolatt (jo2)), KIND=rl) IF ( ( locot (jo2) .AND. z_d_2 < dist_max_cote ) & .OR. ( .NOT. locot (jo2) .AND. z_d_2 < dist_max_large ) ) THEN IF (jo_nb > jp_nb) THEN WRITE (nout,*) 'Trop de point proches ', ja, jot, jo_tab STOP ELSE jo_nb = jo_nb + 1 z_d_1 = z_d_1 + xosrft (jo2) jo_tab (jo_nb) = jo2 END IF END IF ! Calcul du poids END DO SeekOce_4_2 IF (jo_nb == 0) THEN WRITE (nout,*) 'Pas de voisin au voisin ! ', ja, jot, m2oi(jot), m2oj(jot), xolont(jot), xolatt(jot) WRITE (nout,*) 'Point atm ', ja, m2ai(ja), m2aj(ja) WRITE (nout,*) r_earth * REAL (geodist (xolont (jot), xolatt (jot), xolont (jot), xolatt (jot)), KIND=rl) STOP ELSE !-$$ WRITE(nout,*) 'Point atm ', ja, m2ai(ja), m2aj(ja), jot, jo_nb, jo_tab (1:jo_nb) DO jn = 1, jo_nb jo3 = jo_tab (jn) IF (nvo (jo3) < jpa2o ) THEN ! Reste un poids libre nvo (jo3) = nvo (jo3) + 1 wa2o (nvo (jo3), jo3) = xosrft (jo3) / z_d_1 ka2o (nvo (jo3), jo3) = ja !$$$ WRITE (nout,*) ' Vers point oce ', jot, m2oi(jot), m2oj(jot), nvo(jot) ELSE WRITE(nout,*) ' Plus de place ', ja, jo3, nvo (jo3), z_d_1 STOP ENDIF END DO END IF ELSE WRITE (nout,*) 'Pas de voisins proche trouve ' WRITE (nout,*) 'Point atm ', ja, m2ai(ja), m2aj(ja) WRITE (nout,*) dist_max, zzmin, z_d_1 STOP END IF Trouve_4_1 END DO SeekAtm_4_1 !! CALL verif ('Ltotal_dist_3 avant gather') CALL bilan ("Ltotal_dist_3 avant gather ", c_case="terre" , c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Ltotal_dist_3 avant gather ", c_case="terre100", c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Ltotal_dist_3 avant gather ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") !! CALL gather_wei ('apres total_dist_3') CALL verif ( 'Cas total_dist_3') !! WRITE (unit = nout, fmt = * ) 'Nombre de voisins apres completion : ', jma2o !! CALL bilan ("Fin ltotal_dist_3 ", c_case="terre" , c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Fin ltotal_dist_3 ", c_case="terre100", c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Fin ltotal_dist_3 ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") !! END IF Total_dist_3 !! Essai: IF (lessai) THEN WRITE (nout,*) '--- Methode essai ----' no = COUNT ( laland .AND. (itarget == 0_il)) WRITE (nout,*) 'Nombre de points atm non vises : ', no WRITE (nout,*) 'jma2o, jpa2o, jma2o+no ', jma2o, jpa2o, jma2o+no IF ( jma2o + no .GT. jpa2o ) THEN WRITE (nout, *) 'Pas assez de voisins pour la methode essai' STOP END IF !! Liste des points oceans no = COUNT ( looce .AND. iomskp == 0_il) WRITE (nout,*) 'Nombre de points oce cibles : ', no ALLOCATE ( jo_n (no), STAT=ierr) ; CALL chk_allo (ierr, 'jo_n') !! no = 0 ; z_d_1 = 0.0_rl DO jo = 1, jpon IF ( .NOT. looce (jo) ) CYCLE IF ( iomskp (jo) == 1 ) CYCLE no = no + 1 jo_n (no) = jo z_d_1 = z_d_1 + xosrft (jo) END DO no = COUNT ( looce .AND. iomskp == 0_il) WRITE (nout,*) 'Surface : ', z_d_1 SeekAtm_5_1: DO ja = 1, jpan IF (itarget (ja) .GE. 1_il) CYCLE SeekAtm_5_1 IF (iamskp (ja) .EQ. 1_il) CYCLE SeekAtm_5_1 IF (.NOT. laland (ja) ) CYCLE SeekAtm_5_1 DO kn = 1, no jo = jo_n (kn) IF (nvo (jo) < jpa2o ) THEN ! Reste un poids libre nvo (jo) = nvo (jo) + 1 wa2o (nvo (jo), jo) = xosrft (jo) / z_d_1 ka2o (nvo (jo), jo) = ja ELSE WRITE(nout,*) ' Plus de place ', ja, jot, nvo (jot), z_d_1 STOP ENDIF END DO END DO SeekAtm_5_1 !! CALL gather_wei ('apres essai') !! CALL verif ( 'Cas essai') CALL bilan ("Lessai avant gather ", c_case="terre" , c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Lessai avant gather ", c_case="terre100", c_int_atm="Integral", c_int_oce="Integral") CALL bilan ("Lessai avant gather ", c_case="cote" , c_int_atm="Integral", c_int_oce="Integral") !! WRITE (unit = nout, fmt = * ) 'Nombre de voisins apres completion : ', jma2o END IF Essai !! !! Traitement des points atmosphere en bord de cote !! !-$$ IfNear: IF (lnear) THEN !-$$ WRITE (nout,*) 'Extension de 1 point a l''interieur, vers le point ocean le plus proche (lnear=TRUE)' !-$$ SeekAtm1: DO ja = 1, jpan !-$$ IF (.NOT. laland_bord(ja) ) CYCLE SeekAtm1 !-$$ IF (iamskp(ja) .EQ. 1_il ) CYCLE SeekAtm1 !-$$ IF (.NOT. laland(ja) ) CYCLE SeekAtm1 !-$$ zzmin = REAL (geodist(0.0_rl, 90.0_rl, 0.0_rl, -90.0_rl), KIND=rl) !-$$ ltrouve = .FALSE. ; jot = 0_il !-$$ SeekOce1: DO jo = 1, jpon !-$$ IF (.NOT. locot(jo) ) CYCLE SeekOce1 !-$$ IF (iomskp(jo) .EQ. 1_il ) CYCLE SeekOce1 !-$$ z_d_1 = REAL (geodist(xalont(ja), xalatt(ja), xolont(jo), xolatt(jo)), KIND=rl) !-$$ IF (z_d_1 .LT. zzmin) THEN !-$$ zzmin = z_d_1 !-$$ ltrouve = .TRUE. ; jot = jo !-$$ END IF !-$$ ENDDO SeekOce1 !-$$ Trouve11: IF (ltrouve) THEN !-$$ WRITE (nout,*) 'Point oce proche ', jot, m2oi(jot), m2oj(jot) !-$$ !! Calcul du poids !-$$ z_d_1 = xasrft(ja)/xosrft(jot) !-$$ !! Rangement des poids !-$$ ! On cherche si ce point atm est deja vise par ce point ocean (normalement non) !-$$ ltrouve2 = .FALSE. !-$$ L11: DO jn = 1, jpa2o !-$$ IF (ka2o(jn,jot) == ja ) THEN !-$$ wa2o(jn,jot) = wa2o(jn,jot) + z_d_1 !-$$ WRITE (nout,*) 'Ajout au point oce ', jo, m2oi(jo), m2oj(jo) !-$$ ltrouve2 = .TRUE. !-$$ EXIT L11 !-$$ ENDIF !-$$ ENDDO L11 !-$$ Trouve12: IF (.NOT. ltrouve2 ) THEN !-$$ IF (nvo(jot) < jpa2o ) THEN ! Reste un poids libre !-$$ nvo(jot) = nvo(jot) + 1 !-$$ wa2o(nvo(jot),jot) = z_d_1 !-$$ ka2o(nvo(jot),jot) = ja !-$$ WRITE (nout,*) 'Vers point oce ', jot, m2oi(jot), m2oj(jot) !-$$ ELSE ! Trouve2 !-$$ WRITE(nout,*) 'Plus de place' !-$$ STOP !-$$ ENDIF !-$$ END IF Trouve12 !-$$ ELSE ! Trouve !-$$ !! Message d'erreur !-$$ WRITE (unit = nout, fmt = '("Pas de voisin proche pour ", 1I5,2I4)' ) ja, m2ai(ja), m2aj(ja) !-$$ END IF Trouve11 !-$$ END DO SeekAtm1 !-$$ !! !-$$ CALL verif ( 'Cas Near' ) !-$$ !! !-$$ WRITE (unit = nout, fmt = * ) 'Nombre de voisins apres completion : ', jma2o !-$$ !! !-$$ END IF IfNear !! !-$$ IfNei: IF (lnei) THEN !-$$ WRITE(nout,*) 'Route le run-off vers les autres voisins mouilles (lnei=TRUE)' !-$$ i2a = 0_il !-$$ WHERE (RESHAPE(o2amask,(/jpai,jpaj/)) .GT. 0.0_rl ) i2a = 1_il !-$$ SeekAtm2: DO ja = 1_il, jpan !-$$ IF (.NOT. laland_bord(ja) ) CYCLE SeekAtm2 !-$$ IF (iamskp(ja) .EQ. 1_il ) CYCLE SeekAtm2 !-$$ IF (.NOT. laland(ja) ) CYCLE SeekAtm2 !-$$ !! Nombre de voisins !-$$ jai = m2ai(ja) ; jaj = m2aj(ja) ! Indices 2D !-$$ jai_p = mai(jai+1_il) ; jai_m = mai(jai-1_il) !-$$ jaj_p = maj(jaj+1_il) ; jaj_m = maj(jaj-1_il) !-$$ inum = & !-$$ & i2a (jai_p, jaj ) & !-$$ & + i2a (jai_m, jaj ) & !-$$ & + i2a (jai , jaj_p) & !-$$ & + i2a (jai , jaj_m) & !-$$ & + i2a (jai_p, jaj_p) & !-$$ & + i2a (jai_p, jaj_m) & !-$$ & + i2a (jai_m, jaj_p) & !-$$ & + i2a (jai_m, jaj_m) !-$$ IF (inum == 0) WRITE(nout,*) "erreur ", ja !-$$ !! !-$$ IF (i2a(jai_p, jaj ) == 1_il ) CALL complet_run (m1a(jai_p, jaj ), m1a(jai,jaj), inum) !-$$ IF (i2a(jai_m, jaj ) == 1_il ) CALL complet_run (m1a(jai_m, jaj ), m1a(jai,jaj), inum) !-$$ IF (i2a(jai , jaj_p) == 1_il ) CALL complet_run (m1a(jai , jaj_p), m1a(jai,jaj), inum) !-$$ IF (i2a(jai , jaj_m) == 1_il ) CALL complet_run (m1a(jai , jaj_m), m1a(jai,jaj), inum) !-$$ IF (i2a(jai_p, jaj_p) == 1_il ) CALL complet_run (m1a(jai_p, jaj_p), m1a(jai,jaj), inum) !-$$ IF (i2a(jai_p, jaj_m) == 1_il ) CALL complet_run (m1a(jai_p, jaj_m), m1a(jai,jaj), inum) !-$$ IF (i2a(jai_m, jaj_p) == 1_il ) CALL complet_run (m1a(jai_m, jaj_p), m1a(jai,jaj), inum) !-$$ IF (i2a(jai_m, jaj_m) == 1_il ) CALL complet_run (m1a(jai_m, jaj_m), m1a(jai,jaj), inum) !-$$ END DO SeekAtm2 !-$$ !! !-$$ CALL verif ('Cas Nei') !-$$ !! !-$$ WRITE (unit = nout, fmt = * ) 'Nombre de voisins apres completion : ', jma2o !-$$ !! !-$$ END IF IfNei !! CALL int_wei !! !! Controle !! CALL verif (' Normalisation') !! CALL bilan ("Apres lcoast ", c_case="terre" , c_int_atm="Integral", c_int_oce="Local") CALL bilan ("Apres lcoast ", c_case="terre100", c_int_atm="Integral", c_int_oce="Local") CALL bilan ("Apres lcoast ", c_case="cote" , c_int_atm="Integral", c_int_oce="Local") !! WRITE (unit = nout, fmt = * ) 'Nombre de voisins : ', jma2o !! CALL gather_wei ('apres normalisation') CALL verif ('apres gather' ) !! Atm mask interpolated to oce za = REAL (1_il-iamskt, KIND = rl ) CALL inter_a2o (za, a2omask ) !! 1 on Atm interpolated to oce za = REAL (1_il , KIND = rl ) CALL inter_a2o (za, a2ofull ) !! !! Checking weights WRITE (unit = nout, fmt = '("Neighbors a2o : ", 3I9 )' ) & & MINVAL (nvo, nvo /=0), MAXVAL (nvo), jma2o WRITE (unit = nout, fmt = '("Index a2o : ", 2I9 )' ) & & MINVAL (ka2o, ka2o /= 0), MAXVAL (ka2o) jind = MINLOC(wa2o, wa2o /= 0.0_rl) ; jo = jind(2) ; jn = jind(1) joi = m2oi(jo) ; joj = m2oj(jo) WRITE (unit = nout, fmt = '("Weights a2o mini : ", 4I6, 1E13.4, 3I6)' ) & & jo, jn, joi, joj, wa2o(jn,jo), ka2o(jn,jo), m2ai(ka2o(jn,jo)), m2aj(ka2o(jn,jo)) jind = MAXLOC(wa2o, wa2o /= 0.0_rl) ; jo = jind(2) ; jn = jind(1) joi = m2oi(jo) ; joj = m2oj(jo) WRITE (unit = nout, fmt = '("Weights a2o maxi : ", 4I6, 1E13.4, 3I6)' ) & & jo, jn, joi, joj, wa2o(jn,jo), ka2o(jn,jo), m2ai(ka2o(jn,jo)), m2aj(ka2o(jn,jo)) ! jind(1:1) = INT(MAXLOC (nvo),il) ; jo = jind(1) WRITE (unit = nout, fmt = *) jo, m2oi(jo), m2oj(jo), ka2o(1:nvo(jo),jo), & & (m2ai(ka2o(jn,jo)), m2aj(ka2o(jn,jo)), jn = 1, nvo(jo)) !! clweight = "WEIGHTS4" ; cladress = "ADRESSE4" WRITE (unit = nout, fmt = *) cladress, " ", clweight, " ", jpa2o, jma2o !! !! Controle !! WRITE (nout, *) 'Verif avant ecriture' WRITE (nout, *) 'wa2o ', MINVAL (wa2o), MAXVAL (wa2o), COUNT (wa2o > 0.0_rl) WRITE (nout, *) 'ka2o ', MINVAL (ka2o), MAXVAL (ka2o), COUNT (ka2o > 0 ) WRITE (nout, *) 'xasrft ', MINVAL (xasrft), MAXVAL (xasrft) WRITE (nout, *) 'o2amask ', MINVAL (o2amask), MAXVAL (o2amask), SUM (o2amask) WRITE (nout, *) 'xosrft ', MINVAL (xosrft), MAXVAL (xosrft) WRITE (nout, *) 'iomskt ', MINVAL (iomskt), MAXVAL (iomskt), COUNT (iomskt == 0) WRITE (nout, *) 'iomskp ', MINVAL (iomskp), MAXVAL (iomskp), COUNT (iomskp == 0) WRITE (nout, *) 'iamskt ', MINVAL (iamskt), MAXVAL (iamskt), COUNT (iamskt == 0) WRITE (nout, *) 'locot ', COUNT (locot) WRITE (nout, *) 'lacot ', COUNT (lacot) !! ENDIF L_LCOAST !! END IF E_DRYRUN !! !! WRITE (nout, *) ' ' WRITE (nout, *) ' ' WRITE (nout, *) 'Eventuels traitement specifiques' WRITE (nout, *) ' ' WRITE (nout, *) ' ' !! Extension ?? !! Traitements specifiques, tres specifiques ... !! IF ( TRIM(c_suffix) == "_runoffNord09" ) THEN WRITE (nout, *) 'Traitement specifique _runoffNord09' DO jo = 1, jpon IF ( tmask_atl(jo) == 1 .AND. xolatt(jo) > 40.0_rl ) wa2o(:,jo) = wa2o(:,jo) * 0.9 END DO END IF ! IF ( TRIM(c_suffix) == "_runoffNord07" ) THEN WRITE (nout, *) 'Traitement specifique _runoffNord07' DO jo = 1, jpon IF ( tmask_atl(jo) == 1 .AND. xolatt(jo) > 40.0_rl ) wa2o(:,jo) = wa2o(:,jo) * 0.7 END DO END IF IF ( TRIM(c_suffix) == "_runoffNord00" ) THEN WRITE (nout, *) 'Traitement specifique _runoffNord00' DO jo = 1, jpon IF ( tmask_atl(jo) == 1 .AND. xolatt(jo) > 40.0_rl ) wa2o(:,jo) = wa2o(:,jo) * 0.0 END DO END IF !! !! !! Verif nombre de voisins !! IF ( jma2or /= jma2o ) THEN WRITE (nout, *) 'Erreur nombre de voisins dans run.def : ' WRITE (nout, *) ' jma2o calcule / lu : ', jma2o, jma2or !STOP ENDIF !! !! Ecriture des poids !! WRITE (unit = c_r, fmt = '(1i1)' ) rk_out ! IF ( l_wei_i4) THEN WRITE (unit = c_i, fmt = '(1i1)' ) i_4 cfname4 = TRIM(mo_name) // ".wa2o.runoff.i" // TRIM(c_i) // ".r" // TRIM(c_r) // TRIM(c_suffix) OPEN (unit = nwei4a2o, file = TRIM(cfname4), form = 'unformatted', status = 'replace',& & action = 'write') WRITE (nout, *) 'Poids i_4 a -> o: ', TRIM(cfname4) END IF ! IF (l_wei_i8) THEN WRITE (unit = c_i, fmt = '(1i1)' ) i_8 cfname8 = TRIM(mo_name) // ".wa2o.runoff.i" // TRIM(c_i) // ".r" // TRIM(c_r) // TRIM(c_suffix) OPEN (unit = nwei8a2o, file = TRIM(cfname8), form = 'unformatted', status = 'replace',& & action = 'write') WRITE (nout, *) 'Poids i_8 a -> o : ', TRIM(cfname8) END IF ! cldiag_a2o = TRIM(a2o_name) // '.runcoa.diag' clw_a2o = TRIM(mo_name) // '.wa2o.runcoa' clw_a2o_mct= TRIM(mo_name) // '.wa2o_mct.runcoa' ! CALL wri_weights_a2o (cldiag_a2o, clw_a2o, clw_a2o_mct, "4", l_fulldiag, lo_src_grid_frac=.FALSE., lo_dst_grid_frac=.FALSE., & & co_omsk=TRIM(cotes_omsk), co_amsk=TRIM(cotes_amsk) ) ! IF (lriv) THEN !! Dataset 3 : cas particulier du run off des rivieres !! ka2o = 0 ; wa2o = 0.0e0_rl ; nvo = 0 ; nva = 0 ; jma2o = 1 OPEN (unit = nriv, file = '/home/clim/ocean/CPL/DATA/riviere.dat', status = 'old', action = 'read', & & form = 'formatted', position = 'rewind' ) jovoid = 1 DO jr = 1, jpr READ (unit = nriv, fmt = '(1A1,2I3,1I5,1A30)' ) clriv, joi, joj, nlmd, clname clname = ADJUSTL (TRIM (ADJUSTL (clname))) jai = m2ai (nlmd) ; jaj = m2aj (nlmd) jaj = jpaj - jaj + 1 ja = m1a(jai,jaj) IF (joi == 99 .AND. joj == 99 ) THEN !! On ramene sur la ligne j=1 de l'ocean pour faire un bilan plus tard jovoid = jovoid + 1 WRITE (unit = nout, fmt = & & '(1A, " (", 1A30, ") ", 1I3, " : ", 4I5, " : ", 2F6.1 )' ) & & clriv, clname, jr, nlmd, ja, jai, jaj, xalont (ja), xalatt (ja) joi = jovoid ; joj = 1 jo = m1o (joi, joj) nva (ja) = nva (ja) + 1 nvo (jo) = nvo (jo) + 1 ka2o (nvo (jo), jo) = ja wa2o (nvo (jo), jo) = 1.0e0_rl ELSE jo = m1o (joi, joj) nva (ja) = nva (ja) + 1 nvo (jo) = nvo (jo) + 1 ka2o (nvo (jo), jo) = ja IF (lint_atm .AND. lint_oce ) THEN !! Atm envoie : kg.s-1, Oce reçoit : kg.s-1 wa2o(nvo(jo),jo) = 1.0e0_rl ELSE !! Atm envoie kg/m ^2/s. Oce reçoit kg/m^2/s wa2o (nvo (jo), jo) = xasrft (ja) / xosrft (jo) END IF WRITE (unit = nout, fmt = & & '(1A, " (", 1A30, ") ", 1I3, " : ", 4I5, " -> ", 2I4, 1I2, " : ", 1E12.4, " ", 2F6.1, " -> ", 2F6.1)' ) & & clriv, clname, jr, nlmd, ja, jai, jaj, joi, joj, nvo(jo), wa2o (nvo(jo),jo), xalont(ja), xalatt(ja), & & clo_lon(xolont(jo),xalont(ja)), xolatt (jo) END IF END DO CLOSE (unit = nriv) jma2o = MAXVAL (nvo) WRITE (unit = nout, fmt = '("Voisins atm : ", 2I9 )' ) & & MINVAL(nva), MAXVAL(nva) WRITE (unit = nout, fmt = '("Voisins ocean : ", 2I9 )' ) & & MINVAL(nvo), MAXVAL(nvo) WRITE (unit = nout, fmt = '("Point bidons : ", 2I9 )' ) jovoid !! END IF !! !! Verif nombre de voisins !! IF ( jma2or /= jma2o ) THEN WRITE (nout, *) 'Erreur nombre de voisins dans run.def : ' WRITE (nout, *) ' jma2o calcule / lu : ', jma2o, jma2or !STOP ENDIF !! !! Ecriture des poids !! !! Ouverture des fichiers poids !! WRITE (unit = c_r, fmt = '(1i1)' ) rk_out ! IF (l_wei_i4) THEN WRITE (unit = c_i, fmt = '(1i1)' ) i_4 cfname4 = TRIM(mo_name) // ".wa2o.rivflu" // TRIM(c_suffix) // ".i" // TRIM(c_i) // ".r" // TRIM(c_r) OPEN (unit = nwei4a2o, file = TRIM(cfname4), form = 'unformatted', status = 'replace',& & action = 'write') WRITE (nout, *) 'Poids i4 a -> o: ', TRIM(cfname4) END IF ! IF (l_wei_i4) THEN WRITE (unit = c_i, fmt = '(1i1)' ) i_8 cfname8 = TRIM (mo_name) // ".wa2o.rivflu" // TRIM(c_suffix) // ".i" // TRIM(c_i) // ".r" // TRIM(c_r) OPEN (unit = nwei8a2o, file = TRIM(cfname8), form = 'unformatted', status = 'replace',& & action = 'write') WRITE (nout, *) 'Poids i8 a -> o : ', TRIM(cfname8) END IF ! cldiag_a2o = TRIM (a2o_name) // '.rivflu.diag' clw_a2o = TRIM (mo_name) // '.wa2o.rivflu' clw_a2o_mct = "rmp_"//TRIM(camod_t)//"_to_"//TRIM(comod_t)//"_MOSAIC_rivflu" ! CALL wri_weights_a2o (cldiag_a2o, clw_a2o, clw_a2o_mct, "3", l_fulldiag, lo_src_grid_frac=.FALSE., lo_dst_grid_frac=.FALSE., & & co_omsk='noperio', co_amsk='full' ) !! IF (l_wei_i4) CLOSE (unit = nwei8a2o) IF (l_wei_i8) CLOSE (unit = nwei4a2o) !! WRITE (nout,*) 'Appel flioclo' CALL flioclo !! STOP !! CONTAINS !! SUBROUTINE diag_itarget (clmess, lwrite) CHARACTER (len=*), INTENT (in) :: clmess LOGICAL, INTENT (in) :: lwrite INTEGER, SAVE :: nappel = 0 CHARACTER (len=2) :: cc !! nappel = nappel + 1 !! itarget = 0_il ; wtarget = 0.0_rl DO jo = 1, jpon DO jn = 1, jpa2o ja = ka2o (jn, jo) IF (ja /= 0 ) THEN itarget (ja) = itarget (ja) + 1_il wtarget (ja) = wtarget (ja) + wa2o (jn, jo) ENDIF END DO END DO !! !CALL lbc (itarget, naperio, 'T', jpai, jpaj) !CALL lbc (wtarget, naperio, 'T', jpai, jpaj) !! IF (lwrite) THEN !$$$ WRITE(nout,*) 'itarget, etape '// clmess !$$$ CALL prihin (RESHAPE(itarget,(/jpai,jpaj/)),1_il,nout) WRITE (unit = cc, fmt = '(1I2.2)' ) nappel CALL fliodefv (nid_it, 'itarget' // cc, (/1, 2/), v_t=flio_i4) CALL fliodefv (nid_it, 'ntarget' // cc, (/1, 2/), v_t=flio_i4) CALL fliodefv (nid_it, 'wtarget' // cc, (/1, 2/), v_t=flio_r4) CALL fliodefv (nid_it, 'non_target' // cc, (/1, 2/), v_t=flio_i4) CALL flioputa (nid_it, 'itarget' // cc, 'long_name ', clmess // ' itarget' // cc) CALL flioputv (nid_it, 'itarget' // cc, MIN (1, RESHAPE (itarget, (/jpai, jpaj/))) ) CALL flioputa (nid_it, 'ntarget' // cc, 'long_name ', clmess // ' ntarget' // cc) CALL flioputv (nid_it, 'ntarget' // cc, RESHAPE (itarget, (/jpai, jpaj/)) ) CALL flioputa (nid_it, 'wtarget' // cc, 'long_name ', clmess // ' wtarget' // cc) CALL flioputv (nid_it, 'wtarget' // cc, RESHAPE (wtarget, (/jpai, jpaj/)) ) wtarget = 0 WHERE (itarget .EQ. 0 .AND. (laland .OR. lacot)) wtarget = 1.0_rl CALL flioputa (nid_it, 'non_target' // cc, 'long_name ', clmess // ' non_target' // cc) CALL flioputv (nid_it, 'non_target' //cc, RESHAPE (NINT (wtarget), (/jpai, jpaj/)) ) END IF !! END SUBROUTINE diag_itarget !! SUBROUTINE int_wei !! IMPLICIT NONE WRITE(nout,*) 'Debut normalisation' !! IF (.NOT. lint_atm) THEN WRITE(nout,*) 'Multiplication par les surfaces atm' DO jo = 1, jpon DO jn = 1, jma2o ja = ka2o (jn, jo) IF (ja > 0 ) wa2o (jn,jo) = wa2o (jn,jo) * xasrft (ja) END DO ENDDO ENDIF !! IF (.NOT. lint_oce) THEN WRITE (nout,*) 'Division par les surfaces oce' DO jo = 1, jpon DO jn = 1, jma2o wa2o (jn,jo) = wa2o (jn,jo) / xosrft (jo) END DO ENDDO ENDIF END SUBROUTINE int_wei !! SUBROUTINE complet_run (ka1, ka2, kn) INTEGER (kind=il), INTENT(in) :: ka1, ka2 ! Index of atm boxes INTEGER (kind=il), INTENT(in) :: kn ! Number of neighbors INTEGER (kind=il) :: jt !! DO jo = 1, jpon IF (iomskt(jo) == 1 ) CYCLE IF (iomskp(jo) .EQ. 1_il ) CYCLE DO jn = 1, jpa2o IF (ka2o(jn,jo) == ka1 ) THEN nvo(jo) = nvo(jo)+1 IF (nvo(jo) > jpa2o ) THEN WRITE(nout,*) 'Pas assez de voisins pour la completion ' WRITE(nout,fmt='(1I7,2I4,2F10.1)' ) ka1, m2ai(ka1), m2aj(ka1), & & xalont(ka1), xalatt(ka1) WRITE(nout,fmt='(1I7,2I4,2F10.1)' ) ka2, m2ai(ka2), m2aj(ka2), & & xalont(ka2), xalatt(ka2) WRITE(nout,fmt='(1I7,2I4,2F10.1,2I8)') jo, m2oi(jo), m2oj(jo), & & xolont(jo), xolatt(jo), jn, nvo(jo) DO jt = 1, jpa2o WRITE(nout,fmt=*) jt, ka2o(jt,jo), xalont(ka2o(jt,jo)), xalatt(ka2o(jt,jo)), & & lacot(ka2o(jt,jo)) ENDDO END IF ka2o (nvo (jo), jo) = ka2 wa2o (nvo (jo), jo) = wa2o(jn,jo) * xasrft(ka2) / xasrft(ka1) / REAL(kn,KIND=rl) END IF END DO ENDDO END SUBROUTINE complet_run !! SUBROUTINE gather_wei (clmess) CHARACTER (len=*), INTENT (in) :: clmess INTEGER (kind=il) :: jn1, jn2, jo2 !! Eliminate redundancy DO jo = 1, jpon DO jn1 = 1, jpa2o-1 IF (ka2o (jn1,jo) == 0 ) CYCLE DO jn2 = jn1+1, jpa2o IF (ka2o (jn1, jo) == ka2o (jn2, jo) ) THEN WRITE(nout,*) 'Redondance ', jo, jn1, jn2, ka2o(jn1,jo), ka2o (jn2,jo), wa2o(jn1,jo), wa2o(jn2,jo) wa2o (jn1, jo) = wa2o (jn1, jo) + wa2o (jn2, jo) ka2o (jn2, jo) = 0 END IF END DO END DO END DO !! !! Seek for strangeness DO jo = 1, jpon DO jn = 1, jpa2o ja = ka2o (jn, jo) IF (ja == 0) CYCLE IF (wa2o (jn, jo) == 0.0E0 ) THEN WRITE (nout, *) 'Null weight : ', jo, jn, ja, locot (jo), lacot (ja) END IF END DO END DO !! !! Periodicity !DO jn = 1, jpa2o ! CALL lbc (wa2o (jn,:), noperio, 'T', jpoi, jpoj) ! CALL lbc (ka2o (jn,:), noperio, 'T', jpoi, jpoj) !END DO CALL diag_itarget (' gather/perio - ' // clmess, lwrite=.TRUE.) !! !! Gather weights DO jo = 1, jpon ktemp = ka2o (:,jo) ; wtemp = wa2o (:,jo) ka2o (:,jo) = 0_il ; wa2o (:,jo) = 0.0_rl kn = 1_il DO jo2 = 1, jpa2o IF (ktemp (jo2) /= 0_il ) THEN ka2o (kn, jo) = ktemp (jo2) wa2o (kn, jo) = wtemp (jo2) kn = kn + 1_il END IF END DO END DO nvo = COUNT (ka2o /= 0_il, DIM = 1) jma2o = MAXVAL (nvo) ! !! !! Periodicity !DO jn = 1, jpa2o ! CALL lbc (wa2o(jn,:), noperio, 'T', jpoi, jpoj) ! CALL lbc (ka2o(jn,:), noperio, 'T', jpoi, jpoj) !END DO !! nvo = COUNT (ka2o /= 0_il, DIM = 1) jma2o = MAXVAL (nvo) !! END SUBROUTINE gather_wei !! SUBROUTINE verif (clmess, l_non_a) !! CHARACTER (len=*), INTENT (in) :: clmess INTEGER :: ja LOGICAL, INTENT (in), OPTIONAL :: l_non_a !! nvo = COUNT (ka2o /= 0_il, DIM = 1) jma2o = MAXVAL (nvo) CALL diag_itarget (clmess, lwrite=.TRUE.) !! WRITE (nout, *) ' --- Verifs -- : ', TRIM (clmess) WRITE (nout, *) 'wa2o ', MINVAL (wa2o), MAXVAL (wa2o), COUNT (wa2o > 0.0_rl) WRITE (nout, *) 'ka2o ', MINVAL (ka2o), MAXVAL (ka2o), COUNT (ka2o > 0 ) WRITE (nout, *) 'xasrft ', MINVAL (xasrft), MAXVAL (xasrft) WRITE (nout, *) 'o2amask ', MINVAL (o2amask), MAXVAL (o2amask), SUM (o2amask) WRITE (nout, *) 'xosrft ', MINVAL (xosrft), MAXVAL (xosrft) WRITE (nout, *) 'iomskt ', MINVAL (iomskt), MAXVAL (iomskt), COUNT (iomskt == 0) WRITE (nout, *) 'iomskp ', MINVAL (iomskp), MAXVAL (iomskp), COUNT (iomskp == 0) WRITE (nout, *) 'iamskt ', MINVAL (iamskt), MAXVAL (iamskt), COUNT (iamskt == 0) WRITE (nout, *) 'locot ', COUNT (locot) WRITE (nout, *) 'lacot ', COUNT (lacot) WRITE (nout, *) 'Points non attribues : ', COUNT (itarget .EQ. 0 .AND. (laland .OR. lacot)) IF (PRESENT (l_non_a)) THEN DO ja = 1, jpan IF (itarget (ja) == 0 .AND. iamskp(ja) == 0 .AND. (laland (ja) .OR. lacot(ja) )) THEN !-$$ WRITE (UNIT=nout,FMT='("Point atm non attribues : ", 3I7, 2F7.1, 1I4)' ) ja, m2ai(ja), m2aj (ja), xalont(ja), xalatt(ja), itarget(ja) WRITE (nout,*) ja, m2ai(ja), m2aj (ja), xalont(ja), xalatt(ja), itarget(ja), laland(ja), lacot(ja), laoce(ja) ENDIF END DO END IF END SUBROUTINE verif !! !! SUBROUTINE noroute !! Determine les points atmospheres non routes lanoroute = .TRUE. DO jn = 1, jpa2o DO jo = 1, jpon ja = ka2o ( jo, jn) IF ( ja /= 0 .AND. wa2o (jo,jn) > 0.0_rl ) THEN lanoroute = .FALSE. END IF END DO END DO END SUBROUTINE noroute !! SUBROUTINE bilan (clmess, c_case, c_int_atm, c_int_oce) !! CHARACTER (len=*), INTENT (in) :: clmess CHARACTER (len=*), INTENT (in) :: c_case ! Type de champ d'entree CHARACTER (len=*), INTENT (in) :: c_int_atm ! Type d'integrale sur l'atm CHARACTER (len=*), INTENT (in) :: c_int_oce ! Type d'integrale sur l'ocean !! WRITE (nout,*) "Bilan :" // TRIM (clmess) // ", cas :", TRIM (c_case), ", c_int_atm :", TRIM (c_int_atm),& & ", c_int_oce :", TRIM(c_int_oce), "." !! za = 0.0_rl !! SELECT CASE (TRIM (c_case)) CASE ('complet') za = 1.0_rl * REAL (1-iamskp, KIND=rl) SELECT CASE (TRIM (c_int_atm)) CASE ('Integral') za = za * xasrft zsuma = SUM (za) CASE ('Local') zsuma = DOT_PRODUCT (za, xasrft) CASE default WRITE (nout,*) 'Erreur type de bilan atm' ; STOP END SELECT !! CASE ('terre') WHERE (laland .OR. lacot) za = 1.0_rl * REAL (1-iamskp, KIND=rl) END WHERE SELECT CASE (TRIM (c_int_atm)) CASE ('Integral') za = za * xasrft zsuma = SUM (za) CASE ('Local') zsuma = DOT_PRODUCT (za, xasrft*(1.0-o2amask)) CASE default WRITE (nout,*) 'Erreur type de bilan atm' ; STOP END SELECT !! CASE ('terre100') WHERE (laland) za = 1.0_rl * REAL (1-iamskp, KIND=rl) END WHERE SELECT CASE (TRIM (c_int_atm)) CASE ('Integral') za = za * xasrft zsuma = SUM (za) CASE ('Local') zsuma = DOT_PRODUCT (za, xasrft) CASE default WRITE (nout,*) 'Erreur type de bilan atm' ; STOP END SELECT !! CASE ('ocean') WHERE (lacot) za = 1.0_rl * REAL (1-iamskp, KIND=rl) END WHERE SELECT CASE (TRIM (c_int_atm)) CASE ('Integral') za = za * xasrft zsuma = SUM (za) CASE ('Local') zsuma = DOT_PRODUCT (za, xasrft*o2amask) CASE default WRITE (nout,*) TRIM(clmess), ' : Erreur type de bilan atm' ; STOP END SELECT !! CASE ('cote') WHERE (lacot) za = 1.0_rl * REAL (1-iamskp, KIND=rl) END WHERE SELECT CASE (TRIM (c_int_atm)) CASE ('Integral') za = za * xasrft zsuma = SUM (za) CASE ('Local') zsuma = DOT_PRODUCT (za, xasrft*o2amask) CASE default WRITE (nout,*) 'Erreur type de bilan atm' ; STOP END SELECT !! CASE default WRITE (nout,*)'Erreur choix de cas' ; STOP STOP END SELECT !! CALL inter_a2o (za, zo) !! SELECT CASE (TRIM (c_int_oce)) CASE ('Integral') zsumo = DOT_PRODUCT (zo, REAL(1-iomskp, KIND=rl)) CASE ('Local') zsumo = DOT_PRODUCT (zo*xosrft,REAL(1-iomskp, KIND=rl)) CASE default WRITE (nout,*) TRIM (clmess), ' : Erreur type de bilan ' ; STOP END SELECT !! WRITE (nout, fmt='(" bilan atm: ", 1E15.6, ", bilan oce: ", 1E15.6)' ) zsuma, zsumo !! RETURN !! END SUBROUTINE bilan !! END PROGRAM cote