| [808] | 1 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 2 | ! MODULE optcld : |
|---|
| 3 | ! |
|---|
| 4 | ! regroupe les variables communes au routines/fonctions pour l'extrapolation/interpolation |
|---|
| 5 | ! des proprietes optiques des gouttes. |
|---|
| 6 | ! |
|---|
| 7 | ! Contient les routines : |
|---|
| 8 | ! INITIALISATION - lecture de la look-up table et initialisation des variables communes. |
|---|
| 9 | ! LOCATE - Recherche de valeur dans une table. |
|---|
| 10 | ! INTERPOLEMOI - interpolation d'une valeur a partir d'un jeu de donnees. |
|---|
| 11 | ! EXTRAPOLEMOI - extrapolation des valeurs en dehors de la look-up table. |
|---|
| 12 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 13 | MODULE optcld |
|---|
| 14 | IMPLICIT NONE |
|---|
| 15 | REAL,SAVE :: A_ex,A_ab,A_gg, B_ex,B_ab,B_gg |
|---|
| 16 | REAL,SAVE :: fmin_ex,fmin_ab,fmin_gg |
|---|
| 17 | REAL,SAVE :: r0cld |
|---|
| 18 | REAL,SAVE :: lmin,lmax |
|---|
| 19 | INTEGER,SAVE :: npts |
|---|
| 20 | REAL,ALLOCATABLE,SAVE :: tq_ex(:),tq_ab(:),tq_gg(:),tq_wln(:) |
|---|
| 21 | REAL,ALLOCATABLE,SAVE :: ltq_ex(:),ltq_ab(:),ltq_gg(:), & |
|---|
| 22 | ltq_wln(:) |
|---|
| 23 | REAL,SAVE :: frac_c(3),fhvi_c,fhir_c,lseuil_c |
|---|
| 24 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 25 | |
|---|
| 26 | CONTAINS |
|---|
| 27 | |
|---|
| 28 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 29 | ! SUBROUTINE iniqcld : |
|---|
| 30 | ! |
|---|
| 31 | ! Initialisation des variables commune et lecture de la look-up table. |
|---|
| 32 | ! |
|---|
| 33 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 34 | SUBROUTINE iniqcld() |
|---|
| 35 | IMPLICIT NONE |
|---|
| 36 | ! ---------- LOCAL |
|---|
| 37 | INTEGER :: i,iun |
|---|
| 38 | LOGICAL :: ok |
|---|
| 39 | CHARACTER*100 :: tmp |
|---|
| 40 | iun=-1 |
|---|
| 41 | ! rechercher une unite logique libre sur les 200 premieres. |
|---|
| 42 | DO i=1,200 |
|---|
| 43 | INQUIRE (UNIT=i, OPENED=ok) |
|---|
| 44 | IF (.not.ok) THEN |
|---|
| 45 | iun=i |
|---|
| 46 | exit |
|---|
| 47 | ENDIF |
|---|
| 48 | ENDDO |
|---|
| 49 | ! Une belle condition d'arret ==> aucun unite dispo. |
|---|
| 50 | IF (iun.eq.-1) THEN |
|---|
| 51 | PRINT*,"CATASTROPHE !" |
|---|
| 52 | PRINT*,"Impossible de trouver une unite logique libre", & |
|---|
| 53 | " sur les 200 premieres." |
|---|
| 54 | PRINT*,"Je ne peux pas lire optcld.table." |
|---|
| 55 | STOP "Je m'arrete (et salement en plus)..." |
|---|
| 56 | ENDIF |
|---|
| 57 | ! Lecture de la table !!! |
|---|
| 58 | OPEN(iun,file='optcld.table') |
|---|
| 59 | DO i=1,4 |
|---|
| 60 | READ(iun,*) tmp |
|---|
| 61 | ENDDO |
|---|
| 62 | READ(iun,*) npts |
|---|
| 63 | ! Petite pause dans la lecture pour allouer les tableaux |
|---|
| 64 | ALLOCATE(tq_ex(npts)) |
|---|
| 65 | ALLOCATE(tq_ab(npts)) |
|---|
| 66 | ALLOCATE(tq_gg(npts)) |
|---|
| 67 | ALLOCATE(tq_wln(npts)) |
|---|
| 68 | ALLOCATE(ltq_ex(npts)) |
|---|
| 69 | ALLOCATE(ltq_ab(npts)) |
|---|
| 70 | ALLOCATE(ltq_gg(npts)) |
|---|
| 71 | ALLOCATE(ltq_wln(npts)) |
|---|
| 72 | ! Reprise de la |
|---|
| 73 | READ(iun,*) tmp |
|---|
| 74 | READ(iun,'(ES14.7)') r0cld |
|---|
| 75 | READ(iun,*) tmp |
|---|
| 76 | READ(iun,'(2(ES14.7,2X))') lmin,lmax |
|---|
| 77 | DO i=1,3 |
|---|
| 78 | READ(iun,*) tmp |
|---|
| 79 | ENDDO |
|---|
| 80 | READ(iun,'(2(ES14.7,2X))') A_ex,B_ex |
|---|
| 81 | READ(iun,*) tmp |
|---|
| 82 | READ(iun,'(2(ES14.7,2X))') A_ab,B_ab |
|---|
| 83 | READ(iun,*) tmp |
|---|
| 84 | READ(iun,'(2(ES14.7,2X))') A_gg,B_gg |
|---|
| 85 | DO i=1,3 |
|---|
| 86 | READ(iun,*) tmp |
|---|
| 87 | ENDDO |
|---|
| 88 | READ(iun,'(ES14.7)') fmin_ex |
|---|
| 89 | READ(iun,*) tmp |
|---|
| 90 | READ(iun,'(ES14.7)') fmin_ab |
|---|
| 91 | READ(iun,*) tmp |
|---|
| 92 | READ(iun,'(ES14.7)') fmin_gg |
|---|
| 93 | DO i=1,3 |
|---|
| 94 | READ(iun,*) tmp |
|---|
| 95 | ENDDO |
|---|
| 96 | READ(iun,*) (frac_c(i),i=1,3) |
|---|
| 97 | DO i=1,2 |
|---|
| 98 | READ(iun,*) tmp |
|---|
| 99 | ENDDO |
|---|
| 100 | READ(iun,*) fhvi_c,fhir_c,lseuil_c |
|---|
| 101 | READ(iun,*) tmp |
|---|
| 102 | |
|---|
| 103 | DO i=1,npts |
|---|
| 104 | READ(iun,'(4(ES23.15,1X))') & |
|---|
| 105 | tq_wln(i),tq_ex(i),tq_ab(i),tq_gg(i) |
|---|
| 106 | ! ------------ on passe tout en log pour les interpolations |
|---|
| 107 | ltq_wln(i) = alog(tq_wln(i)) |
|---|
| 108 | ltq_ex(i) = alog(tq_ex(i)) |
|---|
| 109 | ltq_ab(i) = alog(tq_ab(i)) |
|---|
| 110 | ltq_gg(i) = alog(tq_gg(i)) |
|---|
| 111 | ENDDO |
|---|
| 112 | CLOSE(iun) |
|---|
| 113 | |
|---|
| 114 | WRITE(*,*) & |
|---|
| 115 | "LECTURE LOOK-UP optcld.table... TERMINEE :)" |
|---|
| 116 | |
|---|
| 117 | END SUBROUTINE iniqcld |
|---|
| 118 | |
|---|
| 119 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 120 | ! SUBROUTINE locate(xx,n,x,j) : |
|---|
| 121 | ! |
|---|
| 122 | ! Recherche de l'indice j de la valeur la plus proche par defaut de x dans le tableau xx de |
|---|
| 123 | ! dimension n |
|---|
| 124 | ! |
|---|
| 125 | ! ARGUMENTS D'ENTREE : |
|---|
| 126 | ! xx : tableau dans lequel rechercher l'indice. |
|---|
| 127 | ! x : valeur recherchee |
|---|
| 128 | ! n : dimension de xx |
|---|
| 129 | ! |
|---|
| 130 | ! ARGUMENT DE SORTIE : |
|---|
| 131 | ! j : indice de la valeur la plus proche PAR DEFAUT de x |
|---|
| 132 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 133 | SUBROUTINE locate(xx,n,x,j) |
|---|
| 134 | IMPLICIT NONE |
|---|
| 135 | ! ---------- INPUT |
|---|
| 136 | INTEGER,INTENT(in) :: n |
|---|
| 137 | REAL ,INTENT(in) :: x,xx(n) |
|---|
| 138 | INTEGER,INTENT(out) :: j |
|---|
| 139 | ! ---------- LOCAL |
|---|
| 140 | INTEGER jl,jm,ju |
|---|
| 141 | jl=0 |
|---|
| 142 | ju=n+1 |
|---|
| 143 | DO WHILE (ju-jl.gt.1) |
|---|
| 144 | jm=(ju+jl)/2 |
|---|
| 145 | IF (jm.eq.0) STOP "ALERTE jm=0 !!" |
|---|
| 146 | IF((xx(n).ge.xx(1)).eqv.(x.ge.xx(jm))) THEN |
|---|
| 147 | jl=jm |
|---|
| 148 | ELSE |
|---|
| 149 | ju=jm |
|---|
| 150 | ENDIF |
|---|
| 151 | ENDDO |
|---|
| 152 | IF (x.eq.xx(1))THEN |
|---|
| 153 | j=1 |
|---|
| 154 | ELSE IF (x.eq.xx(n)) THEN |
|---|
| 155 | j=n-1 |
|---|
| 156 | ELSE |
|---|
| 157 | j=jl |
|---|
| 158 | ENDIF |
|---|
| 159 | RETURN |
|---|
| 160 | END SUBROUTINE locate |
|---|
| 161 | |
|---|
| 162 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 163 | ! SUBROUTINE interpolemoi(ii,x,xin,yin,npts,yout,ver,loga) : |
|---|
| 164 | ! |
|---|
| 165 | ! Interpolation d'une valeur dans un jeu de donnees xin(npts),yin(npts). |
|---|
| 166 | ! |
|---|
| 167 | ! ARGUMENTS D'ENTREE : |
|---|
| 168 | ! ii : ind dans xin de la valeur la plus proche de x par defaut (locate est ton ami) |
|---|
| 169 | ! x : abscisse de la valeur a interpoler. |
|---|
| 170 | ! xin,yin : jeu de valeurs pour l'interpolation. |
|---|
| 171 | ! npts : nombre de points de xin,yin. |
|---|
| 172 | ! ver : type d'interpolation (0 = lineaire / 1 = quadratique) |
|---|
| 173 | ! loga : interpolation en espace log |
|---|
| 174 | ! |
|---|
| 175 | ! ARGUMENT DE SORTIE : |
|---|
| 176 | ! yout : valeur interpolee en x. ! |
|---|
| 177 | ! |
|---|
| 178 | ! NOTES : |
|---|
| 179 | ! - Si loga est utilisee alors xin et yin doivent alors representer les logarithmes du |
|---|
| 180 | ! jeu de donnees. |
|---|
| 181 | ! - Quelque soit la valeur de loga, yout est la valeur interpolee dans l'espace normal |
|---|
| 182 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 183 | SUBROUTINE interpolemoi(ii,x,xin,yin,npts,yout,ver,loga) |
|---|
| 184 | IMPLICIT NONE |
|---|
| 185 | ! ---------- INPUT |
|---|
| 186 | INTEGER,INTENT(in) :: npts |
|---|
| 187 | INTEGER,INTENT(inout) :: ii,ver ! ces 2 variables sont suceptible de changer |
|---|
| 188 | LOGICAL,INTENT(in) :: loga |
|---|
| 189 | REAL ,INTENT(in) :: x,xin(npts),yin(npts) |
|---|
| 190 | ! ---------- OUTPUT |
|---|
| 191 | REAL ,INTENT(out) :: yout |
|---|
| 192 | ! ---------- LOCAL |
|---|
| 193 | REAL :: myx,ytmp,denom |
|---|
| 194 | |
|---|
| 195 | ! ---------- on pourrait ameliorer la condition ici et tester le cas |
|---|
| 196 | ! x plus proche de x(i+1) et utiliser dans ce cas l'interpolation |
|---|
| 197 | ! quadratique... mais est ce vraiment nécessaire ??? |
|---|
| 198 | ! Si on est sur les 2 premiers ou les 2 derniers points : |
|---|
| 199 | ! ===> interpolation lineaire |
|---|
| 200 | IF (ii.eq.1.or.ii.eq.npts-1) ver = 0 |
|---|
| 201 | |
|---|
| 202 | ! Interpolation lineaire |
|---|
| 203 | IF (ver.eq.0) THEN |
|---|
| 204 | myx = x |
|---|
| 205 | IF (loga) myx=alog(myx) |
|---|
| 206 | denom = (xin(ii+1)-xin(ii)) |
|---|
| 207 | ytmp=((yin(ii+1)-yin(ii))*(myx-xin(ii)))/denom+yin(ii) |
|---|
| 208 | IF (loga) THEN |
|---|
| 209 | yout=exp(ytmp) |
|---|
| 210 | ELSE |
|---|
| 211 | yout=ytmp |
|---|
| 212 | ENDIF |
|---|
| 213 | ELSE |
|---|
| 214 | ! ------------ Recherche de l'indice le plus proche de la valeur. |
|---|
| 215 | ! Permet de choisir si l'on interpole avec : |
|---|
| 216 | ! i-1;i;i+1 OU i;i+1;i+2 |
|---|
| 217 | IF (x-xin(ii).gt.xin(ii+1)-x) ii = ii+1 |
|---|
| 218 | myx=x |
|---|
| 219 | IF (loga) myx=alog(myx) |
|---|
| 220 | ytmp = (myx-xin(ii))*(myx-xin(ii+1)) / & |
|---|
| 221 | ((xin(ii-1)-xin(ii)) * & |
|---|
| 222 | (xin(ii-1)-xin(ii+1))) * & |
|---|
| 223 | yin(ii-1) & |
|---|
| 224 | + & |
|---|
| 225 | (myx-xin(ii-1))*(myx-xin(ii+1)) / & |
|---|
| 226 | ((xin(ii)-xin(ii-1)) * & |
|---|
| 227 | (xin(ii)-xin(ii+1))) * & |
|---|
| 228 | yin(ii) & |
|---|
| 229 | + & |
|---|
| 230 | (myx-xin(ii-1))*(myx-xin(ii)) / & |
|---|
| 231 | ((xin(ii+1)-xin(ii-1)) * & |
|---|
| 232 | (xin(ii+1)-xin(ii))) * & |
|---|
| 233 | yin(ii) |
|---|
| 234 | IF (loga) THEN |
|---|
| 235 | yout=exp(ytmp) |
|---|
| 236 | ELSE |
|---|
| 237 | yout=ytmp |
|---|
| 238 | ENDIF |
|---|
| 239 | ENDIF |
|---|
| 240 | RETURN |
|---|
| 241 | END SUBROUTINE interpolemoi |
|---|
| 242 | |
|---|
| 243 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 244 | ! SUBROUTINE extrapolemoi(x,A,B,yout,loga) : |
|---|
| 245 | ! |
|---|
| 246 | ! Extrapolation lineaire d'une valeur a partir de coeff A et B issus d'un jeu de donnees |
|---|
| 247 | ! xin(npts), yin(npts). [ f(x) = A*x + B OU exp(A*log(x)+B) ] |
|---|
| 248 | ! |
|---|
| 249 | ! ARGUMENTS D'ENTREE : |
|---|
| 250 | ! x : abscisse de la valeur a extrapoler. |
|---|
| 251 | ! A : coeff directeur de la droite. |
|---|
| 252 | ! B : ordonnee a l'origine de la droite. |
|---|
| 253 | ! loga : interpolation en espace log |
|---|
| 254 | ! |
|---|
| 255 | ! ARGUMENT DE SORTIE : |
|---|
| 256 | ! yout : valeur extrapolee en x. |
|---|
| 257 | ! |
|---|
| 258 | ! NOTES : |
|---|
| 259 | ! - Si loga est utilisee alors A et B n'ont pas la meme signification (voir forme f(x)) |
|---|
| 260 | ! - Quelque soit la valeur de loga, yout est la valeur extrapolee dans l'espace normal |
|---|
| 261 | !%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|---|
| 262 | SUBROUTINE extrapolemoi(x,A,B,yout,loga) |
|---|
| 263 | IMPLICIT NONE |
|---|
| 264 | ! ---------- INPUT |
|---|
| 265 | REAL ,INTENT(in) :: x,A,B |
|---|
| 266 | LOGICAL,INTENT(in) :: loga |
|---|
| 267 | ! ---------- OUTPUT |
|---|
| 268 | REAL ,INTENT(out) :: yout |
|---|
| 269 | |
|---|
| 270 | IF (loga) THEN |
|---|
| 271 | yout = exp(A*alog(x)+B) |
|---|
| 272 | ELSE |
|---|
| 273 | yout = A*x+B |
|---|
| 274 | ENDIF |
|---|
| 275 | |
|---|
| 276 | END SUBROUTINE extrapolemoi |
|---|
| 277 | |
|---|
| 278 | END MODULE optcld |
|---|
