1 | SUBROUTINE interpol(xinp,yinp,ninp,xout,yout,nout) |
---|
2 | IMPLICIT NONE |
---|
3 | |
---|
4 | c======================================================================= |
---|
5 | c |
---|
6 | c objet: |
---|
7 | c interpolation lineaire |
---|
8 | c |
---|
9 | c ATTENTION!: |
---|
10 | c les tableaux x doivent etre croissants |
---|
11 | c |
---|
12 | c Arguments: |
---|
13 | c ---------- |
---|
14 | c |
---|
15 | c entree: |
---|
16 | c ------- |
---|
17 | c ninp INT nombre de points des donnes |
---|
18 | c xinp(ninp) REAL abbcisses des donnees (doivent etre croissantes) |
---|
19 | c yinp(ninp) REAL valeurs des donnees |
---|
20 | c |
---|
21 | c sorties: |
---|
22 | c -------- |
---|
23 | c nout INT nombre de points interpolles |
---|
24 | c xout(nout) REAL abbcisses (doivent etre croissantes) |
---|
25 | c yout(nout) REAL valeurs interpollees |
---|
26 | c |
---|
27 | c======================================================================= |
---|
28 | |
---|
29 | c----------------------------------------------------------------------- |
---|
30 | c 0. Declarations : |
---|
31 | c ----------------- |
---|
32 | |
---|
33 | INTEGER ninp,nout |
---|
34 | INTEGER iinp,iout,iout0,iout1 |
---|
35 | REAL xinp(ninp),yinp(ninp) |
---|
36 | REAL xout(nout),yout(nout) |
---|
37 | |
---|
38 | c----------------------------------------------------------------------- |
---|
39 | c traitement des valeurs telles que xout<xinp(1): |
---|
40 | c ----------------------------------------------- |
---|
41 | |
---|
42 | DO 100 iout=1,nout |
---|
43 | IF(xout(iout).lt.xinp(1)) THEN |
---|
44 | yout(iout)=yinp(1) |
---|
45 | c si on veut extrapoler aux bords: |
---|
46 | c 1 +(xout(iout)-xinp(1))*(yinp(1)-yinp(2))/(xinp(1)-xinp(2)) |
---|
47 | ELSE |
---|
48 | iout0=iout |
---|
49 | GO TO 200 |
---|
50 | ENDIF |
---|
51 | 100 CONTINUE |
---|
52 | |
---|
53 | RETURN |
---|
54 | |
---|
55 | 200 CONTINUE |
---|
56 | |
---|
57 | c----------------------------------------------------------------------- |
---|
58 | c traitement du centre: |
---|
59 | c --------------------- |
---|
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 |
---|
67 | 320 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 |
---|
77 | 300 CONTINUE |
---|
78 | |
---|
79 | RETURN |
---|
80 | |
---|
81 | 400 CONTINUE |
---|
82 | |
---|
83 | c----------------------------------------------------------------------- |
---|
84 | c traitement des points tels que xout>xinp(ninp): |
---|
85 | c ----------------------------------------------- |
---|
86 | |
---|
87 | DO 500 iout=iout1,nout |
---|
88 | yout(iout)=yinp(ninp) |
---|
89 | c si on veut extrapoler aux bords: |
---|
90 | c 1 +(xout(iout)-xinp(ninp))* |
---|
91 | c 2 (yinp(ninp-1)-yinp(ninp))/(xinp(ninp-1)-xinp(ninp)) |
---|
92 | 500 CONTINUE |
---|
93 | |
---|
94 | RETURN |
---|
95 | END |
---|