source: LMDZ4/branches/IPSL-CM4_IPCC_branch/libf/dyn3d/inter_barx.F @ 5342

!
! $Header$
!
       SUBROUTINE inter_barx ( idatmax,xidat,fdat,imodmax,ximod,fmod ) 

c      .... Auteurs :  Robert Sadourny ,  P. Le Van  .....
c
       IMPLICIT NONE
c    ----------------------------------------------------------
c        INTERPOLATION BARYCENTRIQUE BASEE SUR LES AIRES
c            VERSION UNIDIMENSIONNELLE  ,   EN  LONGITUDE .
c    ----------------------------------------------------------
c
c     idat : indice du champ de donnees, de 1 a idatmax
c     imod : indice du champ du modele,  de 1 a  imodmax
c     fdat(idat) : champ de donnees (entrees)
c     fmod(imod) : champ du modele (sorties)
c     xidat(idat): abscisses des interfaces des mailles donnees
c     ximod(imod): abscisses des interfaces des mailles modele
c      ( L'indice 1 correspond a l'interface mailLE 1 / maille 2)
c      ( Les abscisses sont exprimes en degres)


      INTEGER idatmax, imodmax
      REAL xidat(idatmax),fdat(idatmax),ximod(imodmax),fmod(imodmax)

c    ...  Variables locales ... 
    
      REAL xxid(idatmax+1), xxd(idatmax+1), fdd(idatmax+1)
      REAL  fxd(idatmax+1), xchan(idatmax+1), fdchan(idatmax+1) 
      REAL  xxim(imodmax)

      REAL x0,xim0,dx,dxm
      REAL chmin,chmax,pi

      INTEGER imod,idat,i,ichang,id0,id1,nid,idatmax1
      
      pi = 2. * ASIN(1.)

c  -----------------------------------------------------
c   REDEFINITION DE L'ORIGINE DES ABSCISSES
c    A L'INTERFACE OUEST DE LA PREMIERE MAILLE DU MODELE  
c  -----------------------------------------------------
      DO imod = 1,imodmax
       xxim(imod) = ximod(imod)
      ENDDO

      CALL minmax( imodmax,xxim,chmin,chmax)
       IF( chmax.LT.6.50 )   THEN
c        PRINT  3
c        PRINT *,'   conversion des longit. ximod (donnees en radians)'
c     , ,' en degres  .' 
c        PRINT  3
        DO imod = 1, imodmax
        xxim(imod) = xxim(imod) * 180./pi
        ENDDO  
       ENDIF

      xim0 = xxim(imodmax) - 360.

      DO imod = 1, imodmax
       xxim(imod) = xxim(imod) - xim0
      ENDDO

      idatmax1 = idatmax +1

      DO idat = 1, idatmax
       xxd(idat) = xidat(idat)
      ENDDO

      CALL minmax( idatmax,xxd,chmin,chmax)
       IF( chmax.LT.6.50 )  THEN
c        PRINT  3
c        PRINT *,'   conversion des longit. ximod (donnees en radians)'
c     , ,' en degres  .' 
c        PRINT  3
        DO idat = 1, idatmax
        xxd(idat) = xxd(idat) * 180./pi
        ENDDO  
       ENDIF

      DO idat = 1, idatmax
       xxd(idat) = AMOD( xxd(idat) - xim0, 360. )
       fdd(idat) = fdat (idat)
      ENDDO
c       PRINT *,' xxd redef. origine abscisses '
c       PRINT 2,(xxd(i),i=1,idatmax)

      DO i = 2, idatmax
        IF( ( xxd(i) - xxd(i-1)).LT.0. )  THEN
         ichang = i
         GO TO 5
        ENDIF
      ENDDO
      GO TO 6
c
c  ***  reorganisation  des longitudes entre 0. et 360. degres ****
c
 5    nid = idatmax - ichang +1
      DO i = 1, nid
        xchan (i) = xxd(i+ichang -1 )
        fdchan(i) = fdd(i+ichang -1 )
      ENDDO
       DO i=1,ichang -1
        xchan (i+ nid) = xxd(i)
        fdchan(i+nid) = fdd(i) 
       ENDDO
      DO i =1,idatmax
       xxd(i) = xchan(i)
       fdd(i) = fdchan(i)
      ENDDO

 6    continue


c  ------------------------------------------------
c    translation des champs de donnees par rapport
c    a la nouvelle origine, avec redondance de la
c       maille a cheval sur les bords
c -----------------------------------------------

      id0 = 0
      id1 = 0

      DO idat = 1, idatmax
       IF ( xxd( idatmax1- idat ).LT.360.)   GO TO 10
       id1 = id1 + 1
      ENDDO

 10   DO idat = 1, idatmax
       IF (xxd(idat).GT.0.) GO TO 20
       id0 = id0 + 1
      END DO

 20   IF( id1.EQ.0 ) GO TO 30
      DO idat = 1, id1
       xxid(idat) = xxd(idatmax - id1 + idat) - 360.
       fxd (idat) = fdd(idatmax - id1 + idat)     
      END DO
      DO idat = 1, idatmax - id1
       xxid(idat + id1) = xxd(idat)
       fxd (idat + id1) = fdd(idat)
      END DO

  30  IF(id0.EQ.0) GO TO 40
      DO idat = 1, idatmax - id0
       xxid(idat) = xxd(idat + id0)
       fxd (idat) = fdd(idat + id0)
      END DO

      DO idat = 1, id0
       xxid (idatmax - id0 + idat) =  xxd(idat) + 360.
       fxd  (idatmax - id0 + idat) =  fdd(idat)   
      END DO
      GO TO 50
 
 40   DO idat = 1, idatmax
       xxid(idat)  = xxd(idat)
       fxd (idat)  = fdd(idat)
      ENDDO

 50   xxid(idatmax1) = xxid(1) + 360.
      fxd (idatmax1) = fxd(1)

c  ------------------------------------
c   initialisation du champ du modele

      DO imod = 1, imodmax
       fmod(imod) = 0.
      END DO
 
c      PRINT *,' id0 id1 ',id0,id1
c      PRINT *,' xxim apres translation  '
c      PRINT 2,(xxim(i),i=1,imodmax)
c      PRINT *,' xxid apres translation '
c      PRINT 2,(xxid(i),i=1,idatmax)
c ---------------------------------------
c iteration

      x0   = xim0
      dxm  = 0.
      imod = 1
      idat = 1
 
 100  IF (xxim(imod).LT.xxid(idat)) THEN
       dx   = xxim(imod) - x0
       dxm  = dxm + dx
       fmod(imod) = (fmod(imod) + dx * fxd(idat)) / dxm
       x0   = xxim(imod)
       dxm  = 0.
       imod = imod + 1
       IF (imod.LE.imodmax) GO TO 100
  
      ELSE IF (xxim(imod).GT.xxid(idat)) THEN
       dx   = xxid(idat) - x0
       dxm  = dxm + dx
       fmod(imod) = fmod(imod) + dx * fxd(idat)
       x0   = xxid(idat)
       idat = idat + 1
       GO TO 100
 
      ELSE
       dx   = xxim(imod) - x0
       dxm  = dxm + dx
       fmod(imod) = (fmod(imod) + dx * fxd(idat)) / dxm
       x0   = xxim(imod)
       dxm  = 0.
       imod = imod + 1
       idat = idat + 1
       IF (imod.LE.imodmax) GO TO 100
      END IF
     

3      FORMAT(1x,70(1h-))
2      FORMAT(1x,8f8.2)

       RETURN
       END