Index: /trunk/LMDZ.TITAN/libf/phytitan/comcstfi.h =================================================================== --- /trunk/LMDZ.TITAN/libf/phytitan/comcstfi.h (revision 1516) +++ /trunk/LMDZ.TITAN/libf/phytitan/comcstfi.h (revision 1517) @@ -2,8 +2,7 @@ c INCLUDE comcstfi.h - COMMON/comcstfi/ - * pi,rad,g,r,cpp,rcp,dtphys,daysec,mugaz,omeg + COMMON/comcstfi/rad,g,r,rcp,dtphys,daysec,mugaz,omeg - REAL pi,rad,g,r,cpp,rcp,dtphys,daysec,mugaz,omeg + REAL rad,g,r,rcp,dtphys,daysec,mugaz,omeg c----------------------------------------------------------------------- Index: /trunk/LMDZ.TITAN/libf/phytitan/dyn1d/rcm1d.F =================================================================== --- /trunk/LMDZ.TITAN/libf/phytitan/dyn1d/rcm1d.F (revision 1516) +++ /trunk/LMDZ.TITAN/libf/phytitan/dyn1d/rcm1d.F (revision 1517) @@ -5,4 +5,5 @@ use comgeomphy USE phys_state_var_mod + USE comconst_mod, ONLY: cpp,t0_venus,nu_venus use cpdet_mod, only: ini_cpdet use moyzon_mod, only: plevmoy @@ -76,7 +77,5 @@ character*2 str2 -c normalement dans dyn3d/comconst_mod.F90 - COMMON/cpdetvenus/cppdyn,nu_venus,t0_venus - REAL cppdyn,nu_venus,t0_venus + real pi c======================================================================= @@ -103,5 +102,4 @@ mugaz=28. ! Masse molaire de l'atm (g.mol-1) cpp=1.039e3 - cppdyn=cpp r= 8.314511E+0 *1000.E+0/mugaz rcp= r/cpp @@ -193,5 +191,7 @@ CALL init_phys_lmdz(iim,jjm,llm,1,(/1/)) - call initcomgeomphy + CALL iniphysiq(1,1,llm,daysec,day0,dtphys, + . lati,long,area,cufi,cvfi,rad,g,r,cpp,1) + call ini_cpdet @@ -205,7 +205,4 @@ c par la pression de surface ---> phisfi(1)=0.E+0 - - CALL iniphysiq(1,1,llm,daysec,day0,dtphys, - . lati,long,area,cufi,cvfi,rad,g,r,cpp,1) c Initialisation pour prendre en compte les vents en 1-D Index: /trunk/LMDZ.TITAN/libf/phytitan/iniinterp_horiz.F =================================================================== --- /trunk/LMDZ.TITAN/libf/phytitan/iniinterp_horiz.F (revision 1517) +++ /trunk/LMDZ.TITAN/libf/phytitan/iniinterp_horiz.F (revision 1517) @@ -0,0 +1,179 @@ +C +C $Header$ +C + subroutine iniinterp_horiz (imo,jmo,imn,jmn ,kllm, + & rlonuo,rlatvo,rlonun,rlatvn, + & ktotal,iik,jjk,jk,ik,intersec,airen) + + implicit none + + + +c --------------------------------------------------------- +c Prepare l' interpolation des variables d'une grille LMDZ +c dans une autre grille LMDZ en conservant la quantite +c totale pour les variables intensives (/m2) : ex : Pression au sol +c +c (Pour chaque case autour d'un point scalaire de la nouvelle +c grille, on calcule la surface (en m2)en intersection avec chaque +c case de l'ancienne grille , pour la future interpolation) +c +c on calcule aussi l' aire dans la nouvelle grille +c +c +c Auteur: F.Forget 01/1995 +c ------- +c +c --------------------------------------------------------- +c Declarations: +c ============== +c +c ARGUMENTS +c """"""""" +c INPUT + integer imo, jmo ! dimensions ancienne grille + integer imn,jmn ! dimensions nouvelle grille + integer kllm ! taille du tableau des intersections + real rlonuo(imo+1) ! Latitude et + real rlatvo(jmo) ! longitude des + real rlonun(imn+1) ! bord des + real rlatvn(jmn) ! cases "scalaires" (input) + +c OUTPUT + integer ktotal ! nombre totale d'intersections reperees + integer iik(kllm), jjk(kllm),jk(kllm),ik(kllm) + real intersec(kllm) ! surface des intersections (m2) + real airen (imn+1,jmn+1) ! aire dans la nouvelle grille + + + + +c Autres variables +c """""""""""""""" + integer i,j,ii,jj,k + real a(imo+1),b(imo+1),c(jmo+1),d(jmo+1) + real an(imn+1),bn(imn+1),cn(jmn+1),dn(jmn+1) + real aa, bb,cc,dd + real pi + + pi = 2.*ASIN( 1. ) + + + +c On repere les frontieres des cases : +c =================================== +c +c Attention, on ruse avec des latitudes = 90 deg au pole. + + +c ANcienne grile +c """""""""""""" + a(1) = - rlonuo(imo+1) + b(1) = rlonuo(1) + do i=2,imo+1 + a(i) = rlonuo(i-1) + b(i) = rlonuo(i) + end do + + d(1) = pi/2 + do j=1,jmo + c(j) = rlatvo(j) + d(j+1) = rlatvo(j) + end do + c(jmo+1) = -pi/2 + + +c Nouvelle grille +c """"""""""""""" + an(1) = - rlonun(imn+1) + bn(1) = rlonun(1) + do i=2,imn+1 + an(i) = rlonun(i-1) + bn(i) = rlonun(i) + end do + + dn(1) = pi/2 + do j=1,jmn + cn(j) = rlatvn(j) + dn(j+1) = rlatvn(j) + end do + cn(jmn+1) = -pi/2 + +c Calcul de la surface des cases scalaires de la nouvelle grille +c ============================================================== + do ii=1,imn + 1 + do jj = 1,jmn+1 + airen(ii,jj) = (bn(ii)-an(ii))*(sin(dn(jj))-sin(cn(jj))) + end do + end do + +c Calcul de la surface des intersections +c ====================================== + +c boucle sur la nouvelle grille +c """""""""""""""""""""""""""" + ktotal = 0 + do jj = 1,jmn+1 + do j=1, jmo+1 + if((cn(jj).lt.d(j)).and.(dn(jj).gt.c(j)))then + do ii=1,imn + 1 + do i=1, imo +1 + if ( ((an(ii).lt.b(i)).and.(bn(ii).gt.a(i))) + & .or. ((an(ii).lt.b(i)-2*pi).and.(bn(ii).gt.a(i)-2*pi) + & .and.(b(i)-2*pi.lt.-pi) ) + & .or. ((an(ii).lt.b(i)+2*pi).and.(bn(ii).gt.a(i)+2*pi) + & .and.(a(i)+2*pi.gt.pi) ) + & )then + ktotal = ktotal +1 + iik(ktotal) =ii + jjk(ktotal) =jj + ik(ktotal) =i + jk(ktotal) =j + + dd = min(d(j), dn(jj)) + cc = cn(jj) + if (cc.lt. c(j))cc=c(j) + if((an(ii).lt.b(i)-2*pi).and. + & (bn(ii).gt.a(i)-2*pi)) then + bb = min(b(i)-2*pi,bn(ii)) + aa = an(ii) + if (aa.lt.a(i)-2*pi) aa=a(i)-2*pi + else if((an(ii).lt.b(i)+2*pi).and. + & (bn(ii).gt.a(i)+2*pi)) then + bb = min(b(i)+2*pi,bn(ii)) + aa = an(ii) + if (aa.lt.a(i)+2*pi) aa=a(i)+2*pi + else + bb = min(b(i),bn(ii)) + aa = an(ii) + if (aa.lt.a(i)) aa=a(i) + end if + intersec(ktotal)=(bb-aa)*(sin(dd)-sin(cc)) + end if + end do + end do + end if + end do + end do + + + +c TEST INFO +c DO k=1,ktotal +c ii = iik(k) +c jj = jjk(k) +c i = ik(k) +c j = jk(k) +c if ((ii.eq.10).and.(jj.eq.10).and.(i.eq.10).and.(j.eq.10))then +c if (jj.eq.2.and.(ii.eq.1))then +c write(*,*) '**************** jj=',jj,'ii=',ii +c write(*,*) 'i,j =',i,j +c write(*,*) 'an,bn,cn,dn', an(ii), bn(ii), cn(jj),dn(jj) +c write(*,*) 'a,b,c,d', a(i), b(i), c(j),d(j) +c write(*,*) 'intersec(k)',intersec(k) +c write(*,*) 'airen(ii,jj)=',airen(ii,jj) +c end if +c END DO + + return + end Index: /trunk/LMDZ.TITAN/libf/phytitan/interp_horiz.F =================================================================== --- /trunk/LMDZ.TITAN/libf/phytitan/interp_horiz.F (revision 1517) +++ /trunk/LMDZ.TITAN/libf/phytitan/interp_horiz.F (revision 1517) @@ -0,0 +1,154 @@ +c +c $Id: interp_horiz.F 1403 2010-07-01 09:02:53Z fairhead $ +c + subroutine interp_horiz (varo,varn,imo,jmo,imn,jmn,lm, + & rlonuo,rlatvo,rlonun,rlatvn) + +c=========================================================== +c Interpolation Horizontales des variables d'une grille LMDZ +c (des points SCALAIRES au point SCALAIRES) +c dans une autre grille LMDZ en conservant la quantite +c totale pour les variables intensives (/m2) : ex : Pression au sol +c +c Francois Forget (01/1995) +c=========================================================== + + IMPLICIT NONE + +c Declarations: +c ============== +c +c ARGUMENTS +c """"""""" + + integer imo, jmo ! dimensions ancienne grille (input) + integer imn,jmn ! dimensions nouvelle grille (input) + + real rlonuo(imo+1) ! Latitude et + real rlatvo(jmo) ! longitude des + real rlonun(imn+1) ! bord des + real rlatvn(jmn) ! cases "scalaires" (input) + + integer lm ! dimension verticale (input) + real varo (imo+1, jmo+1,lm) ! var dans l'ancienne grille (input) + real varn (imn+1,jmn+1,lm) ! var dans la nouvelle grille (output) + +c Autres variables +c """""""""""""""" + real airetest(imn+1,jmn+1) + integer ii,jj,l + + real airen (imn+1,jmn+1) ! aire dans la nouvelle grille +c Info sur les ktotal intersection entre les cases new/old grille + integer kllm, k, ktotal + parameter (kllm = 400*200*10) + integer iik(kllm), jjk(kllm),jk(kllm),ik(kllm) + real intersec(kllm) + real R + real totn, tots + + logical firstcall, firsttest, aire_ok + save firsttest + data firsttest /.true./ + data aire_ok /.true./ + + + + + +c initialisation +c -------------- +c Si c'est le premier appel, on prepare l'interpolation +c en calculant pour chaque case autour d'un point scalaire de la +c nouvelle grille, la surface de intersection avec chaque +c case de l'ancienne grille. + + + call iniinterp_horiz (imo,jmo,imn,jmn ,kllm, + & rlonuo,rlatvo,rlonun,rlatvn, + & ktotal,iik,jjk,jk,ik,intersec,airen) + + do l=1,lm + do jj =1 , jmn+1 + do ii=1, imn+1 + varn(ii,jj,l) =0. + end do + end do + end do + +c Interpolation horizontale +c ------------------------- +c boucle sur toute les ktotal intersections entre les cases +c de l'ancienne et la nouvelle grille +c + PRINT *, 'ktotal 1 = ', ktotal + + do k=1,ktotal + do l=1,lm + varn(iik(k),jjk(k),l) = varn(iik(k),jjk(k),l) + & + varo(ik(k), jk(k),l)*intersec(k)/airen(iik(k),jjk(k)) + end do + end do + +c Une seule valeur au pole pour les variables ! : +c ----------------------------------------------- + do l=1, lm + totn =0. + tots =0. + do ii =1, imn+1 + totn = totn + varn(ii,1,l) + tots = tots + varn (ii,jmn+1,l) + end do + do ii =1, imn+1 + varn(ii,1,l) = totn/REAL(imn+1) + varn(ii,jmn+1,l) = tots/REAL(imn+1) + end do + end do + + +c--------------------------------------------------------------- +c TEST TEST TEST TEST TEST TEST TEST TEST TEST TEST +!! if (.not.(firsttest)) goto 99 +!! firsttest = .false. +!! ! write (*,*) 'INTERP. HORIZ. : TEST SUR LES AIRES:' +!! do jj =1 , jmn+1 +!! do ii=1, imn+1 +!! airetest(ii,jj) =0. +!! end do +!! end do +!! PRINT *, 'ktotal = ', ktotal +!! PRINT *, 'jmn+1 =', jmn+1, 'imn+1', imn+1 +!! +!! do k=1,ktotal +!! airetest(iik(k),jjk(k))= airetest(iik(k),jjk(k)) +intersec(k) +!! end DO +!! +!! +!! PRINT *, 'fin boucle' +!! do jj =1 , jmn+1 +!! do ii=1, imn+1 +!! r = airen(ii,jj)/airetest(ii,jj) +!! if ((r.gt.1.001).or.(r.lt.0.999)) then +!! ! write (*,*) '********** PROBLEME D'' AIRES !!!', +!! ! & ' DANS L''INTERPOLATION HORIZONTALE' +!! ! write(*,*)'ii,jj,airen,airetest', +!! ! & ii,jj,airen(ii,jj),airetest(ii,jj) +!! aire_ok = .false. +!! end if +!! end do +!! end do +!! ! if (aire_ok) write(*,*) 'INTERP. HORIZ. : AIRES OK' +!! 99 continue + +c FIN TEST FIN TEST FIN TEST FIN TEST FIN TEST FIN TEST FIN TEST +c--------------------------------------------------------------- + + + + + + + + + return + end Index: /trunk/LMDZ.VENUS/libf/phyvenus/comcstfi.h =================================================================== --- /trunk/LMDZ.VENUS/libf/phyvenus/comcstfi.h (revision 1516) +++ /trunk/LMDZ.VENUS/libf/phyvenus/comcstfi.h (revision 1517) @@ -2,7 +2,7 @@ ! INCLUDE comcstfi.h - COMMON/comcstfi/rad,g,r,cpp,rcp,dtphys,daysec,mugaz,omeg + COMMON/comcstfi/rad,g,r,rcp,dtphys,daysec,mugaz,omeg - REAL rad,g,r,cpp,rcp,dtphys,daysec,mugaz,omeg + REAL rad,g,r,rcp,dtphys,daysec,mugaz,omeg !----------------------------------------------------------------------- Index: unk/LMDZ.VENUS/libf/phyvenus/dyn1d/comgeom.h =================================================================== --- /trunk/LMDZ.VENUS/libf/phyvenus/dyn1d/comgeom.h (revision 1516) +++ (revision ) @@ -1,1 +1,0 @@ -link ../../dyn3d_common/comgeom.h Index: unk/LMDZ.VENUS/libf/phyvenus/dyn1d/comgeom2.h =================================================================== --- /trunk/LMDZ.VENUS/libf/phyvenus/dyn1d/comgeom2.h (revision 1516) +++ (revision ) @@ -1,1 +1,0 @@ -link ../../dyn3d_common/comgeom2.h Index: unk/LMDZ.VENUS/libf/phyvenus/dyn1d/mod_interface_dyn_phys.F90 =================================================================== --- /trunk/LMDZ.VENUS/libf/phyvenus/dyn1d/mod_interface_dyn_phys.F90 (revision 1516) +++ (revision ) @@ -1,1 +1,0 @@ -link ../../dynlonlat_phylonlat/mod_interface_dyn_phys.F90 Index: /trunk/LMDZ.VENUS/libf/phyvenus/dyn1d/rcm1d.F =================================================================== --- /trunk/LMDZ.VENUS/libf/phyvenus/dyn1d/rcm1d.F (revision 1516) +++ /trunk/LMDZ.VENUS/libf/phyvenus/dyn1d/rcm1d.F (revision 1517) @@ -6,4 +6,5 @@ USE phys_state_var_mod use chemparam_mod + USE comconst_mod, ONLY: cpp,t0_venus,nu_venus use cpdet_mod, only: ini_cpdet use moyzon_mod, only: tmoy @@ -11,4 +12,5 @@ . aps,bps,scaleheight,pseudoalt, . disvert_type,pressure_exner + use conc, only: rho IMPLICIT NONE @@ -78,7 +80,4 @@ character*2 str2 -c normalement dans dyn3d/comconst_mod.F90 - COMMON/cpdetvenus/cppdyn,nu_venus,t0_venus - REAL cppdyn,nu_venus,t0_venus real pi @@ -110,5 +109,4 @@ cpp=1.0e3 ! cpp=9.0e2 ! version constante - cppdyn=cpp r= 8.314511E+0 *1000.E+0/mugaz rcp= r/cpp @@ -201,5 +199,7 @@ CALL init_phys_lmdz(iim,jjm,llm,1,(/1/)) - call initcomgeomphy + CALL iniphysiq(1,1,llm,daysec,day0,dtphys, + . lati,long,area,cufi,cvfi,rad,g,r,cpp,1) + call ini_cpdet @@ -213,7 +213,4 @@ c par la pression de surface ---> phisfi(1)=0.E+0 - - CALL iniphysiq(1,1,llm,daysec,day0,dtphys, - . lati,long,area,cufi,cvfi,rad,g,r,cpp,1) c Initialisation pour prendre en compte les vents en 1-D @@ -399,5 +396,10 @@ s du,dv,dtemp,dq,dpsurf) -c print*,"DT APRES PHYSIQ=",day,time +c calcul de rho + rho = 0. +c print*,rho + + +c print*,"DT APRES PHYSIQ=",day,time,dtime c print*,dtemp c print*,temp Index: /trunk/LMDZ.VENUS/libf/phyvenus/iniinterp_horiz.F =================================================================== --- /trunk/LMDZ.VENUS/libf/phyvenus/iniinterp_horiz.F (revision 1517) +++ /trunk/LMDZ.VENUS/libf/phyvenus/iniinterp_horiz.F (revision 1517) @@ -0,0 +1,179 @@ +C +C $Header$ +C + subroutine iniinterp_horiz (imo,jmo,imn,jmn ,kllm, + & rlonuo,rlatvo,rlonun,rlatvn, + & ktotal,iik,jjk,jk,ik,intersec,airen) + + implicit none + + + +c --------------------------------------------------------- +c Prepare l' interpolation des variables d'une grille LMDZ +c dans une autre grille LMDZ en conservant la quantite +c totale pour les variables intensives (/m2) : ex : Pression au sol +c +c (Pour chaque case autour d'un point scalaire de la nouvelle +c grille, on calcule la surface (en m2)en intersection avec chaque +c case de l'ancienne grille , pour la future interpolation) +c +c on calcule aussi l' aire dans la nouvelle grille +c +c +c Auteur: F.Forget 01/1995 +c ------- +c +c --------------------------------------------------------- +c Declarations: +c ============== +c +c ARGUMENTS +c """"""""" +c INPUT + integer imo, jmo ! dimensions ancienne grille + integer imn,jmn ! dimensions nouvelle grille + integer kllm ! taille du tableau des intersections + real rlonuo(imo+1) ! Latitude et + real rlatvo(jmo) ! longitude des + real rlonun(imn+1) ! bord des + real rlatvn(jmn) ! cases "scalaires" (input) + +c OUTPUT + integer ktotal ! nombre totale d'intersections reperees + integer iik(kllm), jjk(kllm),jk(kllm),ik(kllm) + real intersec(kllm) ! surface des intersections (m2) + real airen (imn+1,jmn+1) ! aire dans la nouvelle grille + + + + +c Autres variables +c """""""""""""""" + integer i,j,ii,jj,k + real a(imo+1),b(imo+1),c(jmo+1),d(jmo+1) + real an(imn+1),bn(imn+1),cn(jmn+1),dn(jmn+1) + real aa, bb,cc,dd + real pi + + pi = 2.*ASIN( 1. ) + + + +c On repere les frontieres des cases : +c =================================== +c +c Attention, on ruse avec des latitudes = 90 deg au pole. + + +c ANcienne grile +c """""""""""""" + a(1) = - rlonuo(imo+1) + b(1) = rlonuo(1) + do i=2,imo+1 + a(i) = rlonuo(i-1) + b(i) = rlonuo(i) + end do + + d(1) = pi/2 + do j=1,jmo + c(j) = rlatvo(j) + d(j+1) = rlatvo(j) + end do + c(jmo+1) = -pi/2 + + +c Nouvelle grille +c """"""""""""""" + an(1) = - rlonun(imn+1) + bn(1) = rlonun(1) + do i=2,imn+1 + an(i) = rlonun(i-1) + bn(i) = rlonun(i) + end do + + dn(1) = pi/2 + do j=1,jmn + cn(j) = rlatvn(j) + dn(j+1) = rlatvn(j) + end do + cn(jmn+1) = -pi/2 + +c Calcul de la surface des cases scalaires de la nouvelle grille +c ============================================================== + do ii=1,imn + 1 + do jj = 1,jmn+1 + airen(ii,jj) = (bn(ii)-an(ii))*(sin(dn(jj))-sin(cn(jj))) + end do + end do + +c Calcul de la surface des intersections +c ====================================== + +c boucle sur la nouvelle grille +c """""""""""""""""""""""""""" + ktotal = 0 + do jj = 1,jmn+1 + do j=1, jmo+1 + if((cn(jj).lt.d(j)).and.(dn(jj).gt.c(j)))then + do ii=1,imn + 1 + do i=1, imo +1 + if ( ((an(ii).lt.b(i)).and.(bn(ii).gt.a(i))) + & .or. ((an(ii).lt.b(i)-2*pi).and.(bn(ii).gt.a(i)-2*pi) + & .and.(b(i)-2*pi.lt.-pi) ) + & .or. ((an(ii).lt.b(i)+2*pi).and.(bn(ii).gt.a(i)+2*pi) + & .and.(a(i)+2*pi.gt.pi) ) + & )then + ktotal = ktotal +1 + iik(ktotal) =ii + jjk(ktotal) =jj + ik(ktotal) =i + jk(ktotal) =j + + dd = min(d(j), dn(jj)) + cc = cn(jj) + if (cc.lt. c(j))cc=c(j) + if((an(ii).lt.b(i)-2*pi).and. + & (bn(ii).gt.a(i)-2*pi)) then + bb = min(b(i)-2*pi,bn(ii)) + aa = an(ii) + if (aa.lt.a(i)-2*pi) aa=a(i)-2*pi + else if((an(ii).lt.b(i)+2*pi).and. + & (bn(ii).gt.a(i)+2*pi)) then + bb = min(b(i)+2*pi,bn(ii)) + aa = an(ii) + if (aa.lt.a(i)+2*pi) aa=a(i)+2*pi + else + bb = min(b(i),bn(ii)) + aa = an(ii) + if (aa.lt.a(i)) aa=a(i) + end if + intersec(ktotal)=(bb-aa)*(sin(dd)-sin(cc)) + end if + end do + end do + end if + end do + end do + + + +c TEST INFO +c DO k=1,ktotal +c ii = iik(k) +c jj = jjk(k) +c i = ik(k) +c j = jk(k) +c if ((ii.eq.10).and.(jj.eq.10).and.(i.eq.10).and.(j.eq.10))then +c if (jj.eq.2.and.(ii.eq.1))then +c write(*,*) '**************** jj=',jj,'ii=',ii +c write(*,*) 'i,j =',i,j +c write(*,*) 'an,bn,cn,dn', an(ii), bn(ii), cn(jj),dn(jj) +c write(*,*) 'a,b,c,d', a(i), b(i), c(j),d(j) +c write(*,*) 'intersec(k)',intersec(k) +c write(*,*) 'airen(ii,jj)=',airen(ii,jj) +c end if +c END DO + + return + end Index: /trunk/LMDZ.VENUS/libf/phyvenus/interp_horiz.F =================================================================== --- /trunk/LMDZ.VENUS/libf/phyvenus/interp_horiz.F (revision 1517) +++ /trunk/LMDZ.VENUS/libf/phyvenus/interp_horiz.F (revision 1517) @@ -0,0 +1,154 @@ +c +c $Id: interp_horiz.F 1403 2010-07-01 09:02:53Z fairhead $ +c + subroutine interp_horiz (varo,varn,imo,jmo,imn,jmn,lm, + & rlonuo,rlatvo,rlonun,rlatvn) + +c=========================================================== +c Interpolation Horizontales des variables d'une grille LMDZ +c (des points SCALAIRES au point SCALAIRES) +c dans une autre grille LMDZ en conservant la quantite +c totale pour les variables intensives (/m2) : ex : Pression au sol +c +c Francois Forget (01/1995) +c=========================================================== + + IMPLICIT NONE + +c Declarations: +c ============== +c +c ARGUMENTS +c """"""""" + + integer imo, jmo ! dimensions ancienne grille (input) + integer imn,jmn ! dimensions nouvelle grille (input) + + real rlonuo(imo+1) ! Latitude et + real rlatvo(jmo) ! longitude des + real rlonun(imn+1) ! bord des + real rlatvn(jmn) ! cases "scalaires" (input) + + integer lm ! dimension verticale (input) + real varo (imo+1, jmo+1,lm) ! var dans l'ancienne grille (input) + real varn (imn+1,jmn+1,lm) ! var dans la nouvelle grille (output) + +c Autres variables +c """""""""""""""" + real airetest(imn+1,jmn+1) + integer ii,jj,l + + real airen (imn+1,jmn+1) ! aire dans la nouvelle grille +c Info sur les ktotal intersection entre les cases new/old grille + integer kllm, k, ktotal + parameter (kllm = 400*200*10) + integer iik(kllm), jjk(kllm),jk(kllm),ik(kllm) + real intersec(kllm) + real R + real totn, tots + + logical firstcall, firsttest, aire_ok + save firsttest + data firsttest /.true./ + data aire_ok /.true./ + + + + + +c initialisation +c -------------- +c Si c'est le premier appel, on prepare l'interpolation +c en calculant pour chaque case autour d'un point scalaire de la +c nouvelle grille, la surface de intersection avec chaque +c case de l'ancienne grille. + + + call iniinterp_horiz (imo,jmo,imn,jmn ,kllm, + & rlonuo,rlatvo,rlonun,rlatvn, + & ktotal,iik,jjk,jk,ik,intersec,airen) + + do l=1,lm + do jj =1 , jmn+1 + do ii=1, imn+1 + varn(ii,jj,l) =0. + end do + end do + end do + +c Interpolation horizontale +c ------------------------- +c boucle sur toute les ktotal intersections entre les cases +c de l'ancienne et la nouvelle grille +c + PRINT *, 'ktotal 1 = ', ktotal + + do k=1,ktotal + do l=1,lm + varn(iik(k),jjk(k),l) = varn(iik(k),jjk(k),l) + & + varo(ik(k), jk(k),l)*intersec(k)/airen(iik(k),jjk(k)) + end do + end do + +c Une seule valeur au pole pour les variables ! : +c ----------------------------------------------- + do l=1, lm + totn =0. + tots =0. + do ii =1, imn+1 + totn = totn + varn(ii,1,l) + tots = tots + varn (ii,jmn+1,l) + end do + do ii =1, imn+1 + varn(ii,1,l) = totn/REAL(imn+1) + varn(ii,jmn+1,l) = tots/REAL(imn+1) + end do + end do + + +c--------------------------------------------------------------- +c TEST TEST TEST TEST TEST TEST TEST TEST TEST TEST +!! if (.not.(firsttest)) goto 99 +!! firsttest = .false. +!! ! write (*,*) 'INTERP. HORIZ. : TEST SUR LES AIRES:' +!! do jj =1 , jmn+1 +!! do ii=1, imn+1 +!! airetest(ii,jj) =0. +!! end do +!! end do +!! PRINT *, 'ktotal = ', ktotal +!! PRINT *, 'jmn+1 =', jmn+1, 'imn+1', imn+1 +!! +!! do k=1,ktotal +!! airetest(iik(k),jjk(k))= airetest(iik(k),jjk(k)) +intersec(k) +!! end DO +!! +!! +!! PRINT *, 'fin boucle' +!! do jj =1 , jmn+1 +!! do ii=1, imn+1 +!! r = airen(ii,jj)/airetest(ii,jj) +!! if ((r.gt.1.001).or.(r.lt.0.999)) then +!! ! write (*,*) '********** PROBLEME D'' AIRES !!!', +!! ! & ' DANS L''INTERPOLATION HORIZONTALE' +!! ! write(*,*)'ii,jj,airen,airetest', +!! ! & ii,jj,airen(ii,jj),airetest(ii,jj) +!! aire_ok = .false. +!! end if +!! end do +!! end do +!! ! if (aire_ok) write(*,*) 'INTERP. HORIZ. : AIRES OK' +!! 99 continue + +c FIN TEST FIN TEST FIN TEST FIN TEST FIN TEST FIN TEST FIN TEST +c--------------------------------------------------------------- + + + + + + + + + return + end