dyn3d/Inter_barxy/inter_barx.f

module inter_barx_m

  implicit none

contains

  function inter_barx(dlonid, fdat, rlonimod)

    ! From dyn3d/inter_barx.F, version 1.1.1.1, 2004/05/19 12:53:06

    ! Authors: Robert Sadourny, P. Le Van

    ! INTERPOLATION BARYCENTRIQUE BASEE SUR LES AIRES
    ! VERSION UNIDIMENSIONNELLE, EN LONGITUDE .

    use nr_util, only: assert_eq, pi

    REAL, intent(in):: dlonid(:) ! (idatmax)
    ! abscisses des interfaces des mailles donnees
    
    real, intent(in):: fdat(:) ! (idatmax) champ de donnees

    real, intent(in):: rlonimod(:) ! (imodmax)
    ! Abscisses des interfaces des mailles modele. L'indice 1
    ! correspond a l'interface mailLE 1 / maille 2. Les abscisses sont
    ! exprim\'ees en degres.

    real inter_barx(size(rlonimod)) ! champ du modele

    ! Local: 
    INTEGER idatmax, imodmax
    REAL xxid(size(dlonid)+1), xxd(size(dlonid)+1), fdd(size(dlonid)+1)
    REAL fxd(size(dlonid)+1), xchan(size(dlonid)+1), fdchan(size(dlonid)+1)
    REAL xxim(size(rlonimod))
    REAL x0, xim0, dx, dxm
    REAL chmax
    integer idat ! indice du champ de donnees, de 1 a idatmax
    INTEGER imod ! indice du champ du modele, de 1 a imodmax
    integer i, ichang, id0, id1, nid, idatmax1

    !-----------------------------------------------------

    idatmax = assert_eq(size(dlonid), size(fdat), "inter_barx idatmax")
    imodmax = size(rlonimod)

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

    chmax = maxval(xxim)
    IF(chmax < 6.50) THEN
       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) = dlonid(idat)
    ENDDO

    chmax = maxval(xxd(:idatmax))
    IF(chmax < 6.50) THEN
       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

    i = 2
    DO while (xxd(i) >= xxd(i-1) .and. i < idatmax)
       i = i + 1
    ENDDO
    IF (xxd(i) < xxd(i-1)) THEN
       ichang = i
       ! *** reorganisation des longitudes entre 0. et 360. degres ****
       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
    end IF

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

    id0 = 0
    id1 = 0

    DO idat = 1, idatmax
       IF (xxd(idatmax1- idat) < 360.) exit
       id1 = id1 + 1
    ENDDO

    DO idat = 1, idatmax
       IF (xxd(idat) > 0.) exit
       id0 = id0 + 1
    END DO

    IF(id1 /= 0) then
       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
    end IF

    IF(id0 /= 0) then
       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
    else
       DO idat = 1, idatmax
          xxid(idat) = xxd(idat)
          fxd (idat) = fdd(idat)
       ENDDO
    end IF
    xxid(idatmax1) = xxid(1) + 360.
    fxd (idatmax1) = fxd(1)

    ! initialisation du champ du modele

    inter_barx(:) = 0.

    ! iteration

    x0 = xim0
    dxm = 0.
    imod = 1
    idat = 1

    do while (imod <= imodmax)
       do while (xxim(imod) > xxid(idat))
          dx = xxid(idat) - x0
          dxm = dxm + dx
          inter_barx(imod) = inter_barx(imod) + dx * fxd(idat)
          x0 = xxid(idat)
          idat = idat + 1
       end do
       IF (xxim(imod) < xxid(idat)) THEN
          dx = xxim(imod) - x0
          dxm = dxm + dx
          inter_barx(imod) = (inter_barx(imod) + dx * fxd(idat)) / dxm
          x0 = xxim(imod)
          dxm = 0.
          imod = imod + 1
       ELSE
          dx = xxim(imod) - x0
          dxm = dxm + dx
          inter_barx(imod) = (inter_barx(imod) + dx * fxd(idat)) / dxm
          x0 = xxim(imod)
          dxm = 0.
          imod = imod + 1
          idat = idat + 1
       END IF
    end do

  END function inter_barx

end module inter_barx_m
1	guez	98	module inter_barx_m
2
3			implicit none
4
5			contains
6
7			function inter_barx(dlonid, fdat, rlonimod)
8
9	guez	225	! From dyn3d/inter_barx.F, version 1.1.1.1, 2004/05/19 12:53:06
10	guez	98
11	guez	122	! Authors: Robert Sadourny, P. Le Van
12	guez	98
13	guez	122	! INTERPOLATION BARYCENTRIQUE BASEE SUR LES AIRES
14			! VERSION UNIDIMENSIONNELLE, EN LONGITUDE .
15	guez	98
16	guez	122	use nr_util, only: assert_eq, pi
17	guez	98
18	guez	225	REAL, intent(in):: dlonid(:) ! (idatmax)
19			! abscisses des interfaces des mailles donnees
20
21			real, intent(in):: fdat(:) ! (idatmax) champ de donnees
22	guez	98
23	guez	225	real, intent(in):: rlonimod(:) ! (imodmax)
24			! Abscisses des interfaces des mailles modele. L'indice 1
25			! correspond a l'interface mailLE 1 / maille 2. Les abscisses sont
26			! exprim\'ees en degres.
27	guez	98
28	guez	225	real inter_barx(size(rlonimod)) ! champ du modele
29	guez	98
30	guez	225	! Local:
31	guez	98	INTEGER idatmax, imodmax
32			REAL xxid(size(dlonid)+1), xxd(size(dlonid)+1), fdd(size(dlonid)+1)
33	guez	122	REAL fxd(size(dlonid)+1), xchan(size(dlonid)+1), fdchan(size(dlonid)+1)
34			REAL xxim(size(rlonimod))
35	guez	98	REAL x0, xim0, dx, dxm
36	guez	208	REAL chmax
37	guez	225	integer idat ! indice du champ de donnees, de 1 a idatmax
38			INTEGER imod ! indice du champ du modele, de 1 a imodmax
39			integer i, ichang, id0, id1, nid, idatmax1
40	guez	98
41			!-----------------------------------------------------
42
43			idatmax = assert_eq(size(dlonid), size(fdat), "inter_barx idatmax")
44			imodmax = size(rlonimod)
45
46	guez	122	! REDEFINITION DE L'ORIGINE DES ABSCISSES
47			! A L'INTERFACE OUEST DE LA PREMIERE MAILLE DU MODELE
48	guez	98	DO imod = 1, imodmax
49			xxim(imod) = rlonimod(imod)
50			ENDDO
51
52	guez	208	chmax = maxval(xxim)
53	guez	122	IF(chmax < 6.50) THEN
54	guez	98	DO imod = 1, imodmax
55			xxim(imod) = xxim(imod) * 180./pi
56			ENDDO
57			ENDIF
58
59			xim0 = xxim(imodmax) - 360.
60
61			DO imod = 1, imodmax
62			xxim(imod) = xxim(imod) - xim0
63			ENDDO
64
65			idatmax1 = idatmax +1
66
67			DO idat = 1, idatmax
68			xxd(idat) = dlonid(idat)
69			ENDDO
70
71	guez	208	chmax = maxval(xxd(:idatmax))
72	guez	122	IF(chmax < 6.50) THEN
73	guez	98	DO idat = 1, idatmax
74			xxd(idat) = xxd(idat) * 180./pi
75			ENDDO
76			ENDIF
77
78			DO idat = 1, idatmax
79	guez	122	xxd(idat) = AMOD(xxd(idat) - xim0, 360.)
80	guez	98	fdd(idat) = fdat (idat)
81			ENDDO
82
83			i = 2
84			DO while (xxd(i) >= xxd(i-1) .and. i < idatmax)
85			i = i + 1
86			ENDDO
87			IF (xxd(i) < xxd(i-1)) THEN
88			ichang = i
89	guez	122	! * reorganisation des longitudes entre 0. et 360. degres **
90	guez	98	nid = idatmax - ichang +1
91			DO i = 1, nid
92	guez	122	xchan (i) = xxd(i+ichang -1)
93			fdchan(i) = fdd(i+ichang -1)
94	guez	98	ENDDO
95			DO i=1, ichang -1
96			xchan (i+ nid) = xxd(i)
97			fdchan(i+nid) = fdd(i)
98			ENDDO
99			DO i =1, idatmax
100			xxd(i) = xchan(i)
101			fdd(i) = fdchan(i)
102			ENDDO
103			end IF
104
105	guez	122	! translation des champs de donnees par rapport
106			! a la nouvelle origine, avec redondance de la
107			! maille a cheval sur les bords
108	guez	98
109			id0 = 0
110			id1 = 0
111
112			DO idat = 1, idatmax
113	guez	122	IF (xxd(idatmax1- idat) < 360.) exit
114	guez	98	id1 = id1 + 1
115			ENDDO
116
117			DO idat = 1, idatmax
118	guez	122	IF (xxd(idat) > 0.) exit
119	guez	98	id0 = id0 + 1
120			END DO
121
122	guez	122	IF(id1 /= 0) then
123	guez	98	DO idat = 1, id1
124			xxid(idat) = xxd(idatmax - id1 + idat) - 360.
125			fxd (idat) = fdd(idatmax - id1 + idat)
126			END DO
127			DO idat = 1, idatmax - id1
128			xxid(idat + id1) = xxd(idat)
129			fxd (idat + id1) = fdd(idat)
130			END DO
131			end IF
132
133			IF(id0 /= 0) then
134			DO idat = 1, idatmax - id0
135			xxid(idat) = xxd(idat + id0)
136			fxd (idat) = fdd(idat + id0)
137			END DO
138
139			DO idat = 1, id0
140	guez	122	xxid (idatmax - id0 + idat) = xxd(idat) + 360.
141			fxd (idatmax - id0 + idat) = fdd(idat)
142	guez	98	END DO
143			else
144			DO idat = 1, idatmax
145	guez	122	xxid(idat) = xxd(idat)
146			fxd (idat) = fdd(idat)
147	guez	98	ENDDO
148			end IF
149			xxid(idatmax1) = xxid(1) + 360.
150			fxd (idatmax1) = fxd(1)
151
152	guez	122	! initialisation du champ du modele
153	guez	98
154			inter_barx(:) = 0.
155
156			! iteration
157
158	guez	122	x0 = xim0
159			dxm = 0.
160	guez	98	imod = 1
161			idat = 1
162
163			do while (imod <= imodmax)
164	guez	122	do while (xxim(imod) > xxid(idat))
165			dx = xxid(idat) - x0
166			dxm = dxm + dx
167	guez	98	inter_barx(imod) = inter_barx(imod) + dx * fxd(idat)
168	guez	122	x0 = xxid(idat)
169	guez	98	idat = idat + 1
170			end do
171	guez	122	IF (xxim(imod) < xxid(idat)) THEN
172			dx = xxim(imod) - x0
173			dxm = dxm + dx
174	guez	98	inter_barx(imod) = (inter_barx(imod) + dx * fxd(idat)) / dxm
175	guez	122	x0 = xxim(imod)
176			dxm = 0.
177	guez	98	imod = imod + 1
178			ELSE
179	guez	122	dx = xxim(imod) - x0
180			dxm = dxm + dx
181	guez	98	inter_barx(imod) = (inter_barx(imod) + dx * fxd(idat)) / dxm
182	guez	122	x0 = xxim(imod)
183			dxm = 0.
184	guez	98	imod = imod + 1
185			idat = idat + 1
186			END IF
187			end do
188
189			END function inter_barx
190
191			end module inter_barx_m