| 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 |
|---|