source: dynamico_lmdz/simple_physics/phyparam/param/interpol.F @ 4187

Last change on this file since 4187 was 4176, checked in by dubos, 5 years ago

simple_physics : copy code from FH

File size: 2.5 KB
Line 
1      SUBROUTINE interpol(xinp,yinp,ninp,xout,yout,nout)
2      IMPLICIT NONE
3
4c=======================================================================
5c
6c  objet:
7c  interpolation lineaire
8c
9c  ATTENTION!:
10c  les tableaux x doivent etre croissants
11c
12c  Arguments:
13c  ----------
14c
15c  entree:
16c  -------
17c  ninp        INT  nombre de points des donnes
18c  xinp(ninp)  REAL abbcisses des donnees (doivent etre croissantes)
19c  yinp(ninp)  REAL valeurs des donnees
20c
21c  sorties:
22c  --------
23c  nout        INT  nombre de points interpolles
24c  xout(nout)  REAL abbcisses (doivent etre croissantes)
25c  yout(nout)  REAL valeurs interpollees
26c
27c=======================================================================
28
29c-----------------------------------------------------------------------
30c    0. Declarations :
31c    -----------------
32
33      INTEGER ninp,nout
34      INTEGER iinp,iout,iout0,iout1
35      REAL xinp(ninp),yinp(ninp)
36      REAL xout(nout),yout(nout)
37
38c-----------------------------------------------------------------------
39c   traitement des valeurs telles que xout<xinp(1):
40c   -----------------------------------------------
41
42      DO 100 iout=1,nout
43         IF(xout(iout).lt.xinp(1)) THEN
44              yout(iout)=yinp(1)
45c   si on veut extrapoler aux bords:
46c    1        +(xout(iout)-xinp(1))*(yinp(1)-yinp(2))/(xinp(1)-xinp(2))
47          ELSE
48              iout0=iout
49              GO TO 200
50          ENDIF
51100   CONTINUE
52
53      RETURN
54
55200   CONTINUE
56
57c-----------------------------------------------------------------------
58c   traitement du centre:
59c   ---------------------
60
61      iinp=1
62      DO 300 iout=iout0,nout
63          IF(xout(iout).ge.xinp(ninp)) THEN
64             iout1=iout
65             GO TO 400
66          ENDIF
67320       CONTINUE
68          IF(  (xout(iout).ge.xinp(iinp)  ) .AND.
69     1         (xout(iout).lt.xinp(iinp+1))        ) THEN
70               yout(iout)=yinp(iinp)+
71     1         (xout(iout)-xinp(iinp))*
72     2         (yinp(iinp+1)-yinp(iinp))/(xinp(iinp+1)-xinp(iinp))
73          ELSE
74               iinp=iinp+1
75               GO TO 320
76          ENDIF
77300   CONTINUE
78
79      RETURN
80
81400   CONTINUE
82
83c-----------------------------------------------------------------------
84c   traitement des points tels que xout>xinp(ninp):
85c   -----------------------------------------------
86
87      DO 500 iout=iout1,nout
88          yout(iout)=yinp(ninp)
89c   si on veut extrapoler aux bords:
90c    1    +(xout(iout)-xinp(ninp))*
91c    2    (yinp(ninp-1)-yinp(ninp))/(xinp(ninp-1)-xinp(ninp))
92500   CONTINUE
93
94      RETURN
95      END
Note: See TracBrowser for help on using the repository browser.