source: LMDZ6/trunk/libf/dyn3d/tetaleveli1j.F90 @ 5248

!================================================================
!================================================================
SUBROUTINE tetaleveli1j(ilon,ilev,lnew,pgcm,pres,Qgcm,Qpres)
  !================================================================
  !================================================================

  ! FH 2008/05/09 On elimine toutes les clefs physiques dans la dynamique
   ! USE dimphy
  IMPLICIT none

#include "dimensions.h"
  !cccc#include "dimphy.h"

  !================================================================
  !
  ! Interpoler des champs 3-D u, v et g du modele a un niveau de
  ! pression donnee (pres)
  !
  ! INPUT:  ilon ----- nombre de points
  !     ilev ----- nombre de couches
  !     lnew ----- true si on doit reinitialiser les poids
  !     pgcm ----- pressions modeles
  !     pres ----- pression vers laquelle on interpolle
  !     Qgcm ----- champ GCM
  !     Qpres ---- champ interpolle au niveau pres
  !
  !================================================================
  !
  !   arguments :
  !   -----------

  INTEGER :: ilon, ilev
  logical :: lnew

  REAL :: pgcm(ilon,ilev)
  REAL :: Qgcm(ilon,ilev)
  real :: pres
  REAL :: Qpres(ilon)

  !   local :
  !   -------

  !IM 211004
  ! INTEGER lt(klon), lb(klon)
  ! REAL ptop, pbot, aist(klon), aisb(klon)
  !
#include "paramet.h"
  !
  INTEGER :: lt(ip1jm), lb(ip1jm)
  REAL :: ptop, pbot, aist(ip1jm), aisb(ip1jm)
  !MI 211004
  save lt,lb,ptop,pbot,aist,aisb

  INTEGER :: i, k
  !
  ! PRINT*,'tetalevel pres=',pres
  !=====================================================================
  if (lnew) then
  !   on réinitialise les réindicages et les poids
  !=====================================================================


  ! Chercher les 2 couches les plus proches du niveau a obtenir
  !
  ! Eventuellement, faire l'extrapolation a partir des deux couches
  ! les plus basses ou les deux couches les plus hautes:
  DO i = 1, ilon
  !IM      IF ( ABS(pres-pgcm(i,ilev) ) .LT.
     IF ( ABS(pres-pgcm(i,ilev) ) .GT. &
           ABS(pres-pgcm(i,1)) ) THEN
        lt(i) = ilev     ! 2
        lb(i) = ilev-1   ! 1
     ELSE
        lt(i) = 2
        lb(i) = 1
     ENDIF
  !IM   PRINT*,'i, ABS(pres-pgcm),ABS(pres-pgcm)',
  !IM  .i, ABS(pres-pgcm(i,ilev)),ABS(pres-pgcm(i,1))
  END DO
  DO k = 1, ilev-1
     DO i = 1, ilon
        pbot = pgcm(i,k)
        ptop = pgcm(i,k+1)
  !IM         IF (ptop.LE.pres .AND. pbot.GE.pres) THEN
        IF (ptop.GE.pres .AND. pbot.LE.pres) THEN
           lt(i) = k+1
           lb(i) = k
        ENDIF
     END DO
  END DO
  !
  ! Interpolation lineaire:
  !
  DO i = 1, ilon
  ! interpolation en logarithme de pression:
  !
  ! ...   Modif . P. Le Van    ( 20/01/98) ....
  !   Modif Frédéric Hourdin (3/01/02)

    IF(pgcm(i,lb(i)).EQ.0.OR. &
          pgcm(i,lt(i)).EQ.0.) THEN
  !
    PRINT*,'i,lb,lt,2pgcm,pres',i,lb(i), &
          lt(i),pgcm(i,lb(i)),pgcm(i,lt(i)),pres
  !
    ENDIF
  !
    aist(i) = LOG( pgcm(i,lb(i))/ pres ) &
          / LOG( pgcm(i,lb(i))/ pgcm(i,lt(i)) )
    aisb(i) = LOG( pres / pgcm(i,lt(i)) ) &
          / LOG( pgcm(i,lb(i))/ pgcm(i,lt(i)))
  enddo


  endif ! lnew

  !======================================================================
  !    inteprollation
  !======================================================================

  do i=1,ilon
     Qpres(i)= Qgcm(i,lb(i))*aisb(i)+Qgcm(i,lt(i))*aist(i)
  !IM      PRINT*,'i,Qgcm,Qpres',i,Qgcm(i,lb(i)),aisb(i),
  !IM  $   Qgcm(i,lt(i)),aist(i),Qpres(i)
  enddo
  !
  ! Je mets les vents a zero quand je rencontre une montagne
  do i = 1, ilon
  !IM      if (pgcm(i,1).LT.pres) THEN
     if (pgcm(i,1).GT.pres) THEN
        ! Qpres(i)=1e33
        Qpres(i)=1e+20
  !IM         PRINT*,'i,pgcm(i,1),pres =',i,pgcm(i,1),pres
     endif
  enddo

  !
  RETURN
END SUBROUTINE tetaleveli1j