source: LMDZ6/branches/Amaury_dev/libf/dyn3d/tetaleveli1j1.F90 @ 5220

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

  ! FH 2008/05/09 On elimine toutes les clefs physiques dans la dynamique
USE lmdz_dimensions, ONLY: iim, jjm, llm, ndm
  USE lmdz_paramet
  IMPLICIT NONE



  !================================================================

  ! 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)
  !


  INTEGER :: lt(ip1jmp1), lb(ip1jmp1)
  REAL :: ptop, pbot, aist(ip1jmp1), aisb(ip1jmp1)
  !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)) > &
              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>=pres .AND. pbot<=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))==0.OR. &
              pgcm(i, lt(i))==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)>pres) THEN
      ! Qpres(i)=1e33
      Qpres(i) = 1e+20
      !IM         PRINT*,'i,pgcm(i,1),pres =',i,pgcm(i,1),pres
    endif
  enddo

  !

END SUBROUTINE tetaleveli1j1