source: LMDZ.3.3/branches/LF/libf/dyn3d/calc_uvtq.F

      SUBROUTINE  calc_uvtq (scalarq, plev,xtlon,ylat,champi,jlat ,
     ,  rlon,rlat, phislmd,pslmd,tslmd,pls,pk,pks,p,   champsor   ) 
c
c     Auteur :  P. Le Van  
cc
      IMPLICIT NONE
c
#include "dimensions.h"
#include "paramet.h"
#include "comconst.h"
#include "comvert.h"
#include "para_netcdf.h"

      INTEGER  jlat

      LOGICAl scalarq
      REAL*8  plev(levs)
      REAL    xtlon(lons),ylat(lats),champi(lons,lats,levs)
      REAL    phislmd(iip1,jlat),rlon(iip1),rlat(jlat) ,
     ,        pslmd(iip1,jlat),tslmd(iip1,jlat),pls(iip1,jlat,llm) ,
     ,        pk(iip1,jlat,llm),pks(iip1,jlat),p(iip1,jlat,llmp1) ,
     ,        champsor(iip1,jlat,llm),alpha(iip1,jjp1,llm),
     ,        beta(iip1,jjp1,llm),pkf(iip1,jjp1,llm)
c
c
      REAL champint(iip1,jjp1)
      LOGICAL invlev,invlon,invlat
      REAL  pmbar
      INTEGER ip180
      COMMON/invlog/invlev,invlon,invlat,ip180,pmbar

      REAL  ps,phisol,ts
      COMMON/deuxD/ps(lons,lats),phisol(lons,lats),ts(lons,lats)
  
      INTEGER l,i,j,lll,ind,iipjlat
c
      REAL uu(lons,lats,levs)
      REAL phisolmd(iip1,jjp1)
      SAVE phisolmd
      REAL ax(levs),ay(levs),yder(levs),bx,by
      REAL champhor(iip1,jjp1,levs)
      REAL prefkap,unskap

      EXTERNAL grille_m, gr_int_dyn, pression, exner_hyb, spline, splint 
c

      IF( scalarq )  GO TO 100

       CALL grille_m(lons, lats, xtlon, ylat, ps ,
     ,               iim, jlat, rlon, rlat, champint )


       CALL gr_int_dyn(champint, pslmd, iim, jlat )

       CALL grille_m(lons, lats, xtlon, ylat, phisol,
     ,             iim, jlat, rlon, rlat, champint  )
       CALL gr_int_dyn(champint, phisolmd, iim, jlat )

       CALL grille_m(lons, lats, xtlon, ylat, ts,
     ,             iim, jlat, rlon, rlat, champint)
       CALL gr_int_dyn(champint, tslmd, iim, jlat  )

       DO j = 1, jlat
       DO i = 1, iim
        pslmd(i,j)    = pslmd(i,j)*( 1.0+ (phisolmd(i,j)-phislmd(i,j))
     ,                             /287.0/tslmd(i,j))
       ENDDO
        pslmd(iip1,j) = pslmd(1,j)
       ENDDO

       iipjlat = iip1 * jlat
       CALL pression( iipjlat, ap, bp, pslmd, p )

       CALL exner_hyb(iipjlat,pslmd,p,alpha,beta,pks,pk,pkf )


c    ....  Calcul de pls , pression au milieu des couches ,en Pascals ....
c
      PRINT *,' Pref kappa unskap  ',preff,kappa

      prefkap =  preff  ** kappa
      unskap  =  1./ kappa

       DO l = 1, llm
        DO j = 1,jlat
         DO i = 1, iim
          pls(i,j,l) = preff * ( pk(i,j,l)/cpp) ** unskap
         ENDDO
          pls(iip1,j,l) = pls(1,j,l)
        ENDDO
       ENDDO

c
 100   CONTINUE
c
c
       DO l = 1, levs
        DO j = 1, lats
          DO i= 1, lons
           uu (i,j,l) = champi (i,j,l)
          ENDDO
        ENDDO
       ENDDO

        IF ( invlon )   THEN
c
c   .....    on tourne les longitudes  pour  avoir de - pi a pi ....
c

         DO l = 1, levs
          DO j = 1,lats

           DO i = ip180,lons
           ind  = i-ip180 +1
           uu (ind,j,l) = champi (i,j,l)
           ENDDO
   
           DO i= ind +1,lons
           uu (i,j,l)  = champi (i-ind,j,l)
           ENDDO

          ENDDO
        ENDDO

        ENDIF
c

c   *****   fin  de   IF(invlon)   ****
         
         DO l = 1, levs
           DO j = 1, lats
              DO i = 1, lons
               champi (i,j,l) = uu (i,j,l)
              ENDDO
           ENDDO
         ENDDO

         IF ( invlat )    THEN

          DO l = 1, levs
           DO j = 1, lats
              DO i = 1, lons
               champi (i,lats-j+1,l) = uu (i,j,l)
              ENDDO
           ENDDO
          ENDDO

         ENDIF

c   .....   Interpol.  horizontale  ......
c   **************************************
c
      DO l = 1, levs

      CALL grille_m(lons, lats, xtlon, ylat, champi(1,1,l),
     ,               iim, jlat, rlon, rlat, champint    )
      lll = l
      IF( invlev) lll = levs -l +1
      CALL gr_int_dyn(champint, champhor(1,1,lll), iim, jlat)
c
      ENDDO

c    ...  Interpolation  verticale ....   
c    **********************************
c
cc   Interpolation verticale par spline

c     pmbar = 100. ( si les donnees de pression sont en mb ) ou = 1.
c

      DO j=1,jlat
        DO  i=1,iim
          DO l=1,levs
            ax(l)= plev(levs-l+1) *  pmbar
            ay(l)= champhor(i,j,levs-l+1)
          ENDDO
          CALL spline (ax,ay,levs,1.e30,1.e30,yder)
          DO  l=1,llm
           bx = pls(i,j,llm-l+1)
            CALL splint (ax,ay,yder,levs,bx,by)
            champsor(i,j,llm-l+1)= by
          ENDDO
        ENDDO
        DO l = 1, llm
          champsor(iip1, j, l) = champsor(1, j, l)
        ENDDO
      ENDDO


      RETURN
      END