Philippe 02/05/2016: moved all LIBTOOLS files in LIBTOOLS directory
[MNH-git_open_source-lfs.git] / LIBTOOLS / tools / diachro / src / DIAPRO / interpxyz.f90
1 !     ######spl
2 MODULE MODI_INTERPXYZ
3 INTERFACE
4 !     #####################################################################
5       SUBROUTINE INTERPXYZ(PAX,PAY,PAZ,PCHAMP,        &
6                            PX,PY,PZ,                  &
7                            PXOR,PYOR,PDX,PDY,         &
8                            PZL,OTRAJ_GROUP,           &
9                            PRESX,PRESY,PRESZ,PRESCHAMP)
10 !     #####################################################################
11 !
12 !
13 ! entrees
14 !
15 REAL, DIMENSION(:,:,:),    INTENT(IN)     :: PAX,PAY,PAZ,PCHAMP
16                                                                  !
17                                                                  !
18                                                                  !
19 REAL,                      INTENT(INOUT)     :: PX,PY,PZ            !
20 REAL,                      INTENT(IN)     :: PXOR,PYOR,PDX,PDY   !
21 REAL, DIMENSION(:,:,:),    INTENT(IN)     :: PZL                 !
22 LOGICAL,                   INTENT(IN)     :: OTRAJ_GROUP
23 !
24 ! sorties
25 !
26 REAL,                      INTENT(OUT)    :: PRESX,PRESY,PRESZ,PRESCHAMP
27 !
28 !
29 END SUBROUTINE INTERPXYZ
30 !
31 END INTERFACE
32 !
33 END MODULE MODI_INTERPXYZ
34 !     ######spl
35       SUBROUTINE INTERPXYZ(PAX,PAY,PAZ,PCHAMP,        &
36                            PX,PY,PZ,                  &
37                            PXOR,PYOR,PDX,PDY,         &
38                            PZL,OTRAJ_GROUP,           &
39                            PRESX,PRESY,PRESZ,PRESCHAMP)
40 !     #####################################################################
41 !
42 !
43 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
44 ! BUT DE LA ROUTINE : interpoler les trois champs (3D) LG?M 
45 ! (ou trois champs 3D quelconques ecrits sur les points de masse) 
46 ! en un point M, de coordonnees cartesiennes (x,y,z) 
47 ! a priori non-situe sur un point de grille. 
48 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!                           
49 !
50 !
51 !
52 !!!!!!!!!!!!!!!!!!!!!!
53 ! Declarations
54 !!!!!!!!!!!!!!!!!!!!!!
55 !
56 IMPLICIT NONE
57 !
58 ! entrees
59 !
60 REAL, DIMENSION(:,:,:),    INTENT(IN)     :: PAX,PAY,PAZ,PCHAMP
61                                                                  !
62                                                                  !
63                                                                  !
64 REAL,                      INTENT(INOUT)     :: PX,PY,PZ            !
65 REAL,                      INTENT(IN)     :: PXOR,PYOR,PDX,PDY   !
66 REAL, DIMENSION(:,:,:),    INTENT(IN)     :: PZL                 !
67 LOGICAL,                   INTENT(IN)     :: OTRAJ_GROUP
68 !
69 ! sorties
70 !
71 REAL,                      INTENT(OUT)    :: PRESX,PRESY,PRESZ,PRESCHAMP
72 !
73 ! locales
74 !
75 INTEGER                              :: II,IJ,IK,JK             !
76 INTEGER                              :: IKU                     !
77 REAL                                 :: ZEPS1,ZEPS2,ZEPS3       !
78 REAL                                 :: ZXREL,ZYREL             !
79 REAL, DIMENSION(SIZE(PZL,3))         :: ZZLXY                   !
80 !
81 !
82 ! initialisations des variables locales
83 !
84 IKU=SIZE(PZL,3)
85 !
86 !
87 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
88 ! 1. Recherche de la maille contenant le point M(PX,PY,PZ) -> II,IJ,IK
89 !    Position de M au sein de la maille                    -> ZEPS1,ZEPS2,ZEPS3
90 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
91 !
92 !
93 ! 1.a partie horizontale
94 !
95 ZXREL=(PX-PXOR)/PDX+2
96 ZYREL=(PY-PYOR)/PDY+2
97 !
98 II=FLOOR(ZXREL)
99 IJ=FLOOR(ZYREL)
100 !
101 ZEPS1=ZXREL-REAL(II)
102 ZEPS2=ZYREL-REAL(IJ)
103 !
104 !
105 ! 1.b partie verticale
106 !
107 ! 1.b.1 altitude des niveaux du modele sur la verticale (PX,PY)
108 !
109 DO JK=1,IKU
110   ZZLXY(JK)=ZEPS2*(ZEPS1*(PZL(II+1,IJ+1,JK))+(1-ZEPS1)*(PZL(II,IJ+1,JK)))     &
111              + (1-ZEPS2)*(ZEPS1*(PZL(II+1,IJ,JK))+(1-ZEPS1)*(PZL(II,IJ,JK)))
112 ENDDO
113 !
114 IK=999
115 DO JK=2,IKU
116   IF (ZZLXY(JK).GE.PZ) THEN
117     IK=JK-1
118     EXIT 
119   ENDIF
120 ENDDO
121 !
122 IF (IK==1) THEN
123   print *,'la particule est sous le sol'
124   print *,' on la remonte a zs + dz/2 = ', ZZLXY(2)
125   PZ=ZZLXY(2)
126 ENDIF
127 !
128 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
129 !!!!!!!!!!!!!!Emergency exit!!!!!!!!!!!!!!
130 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
131 IF (IK==999) THEN
132    PRINT*,'PROBLEME AU POINT',II,IJ
133    PRINT*,'XREL, YREL, Z =',ZXREL,ZYREL,PZ
134    PRINT*,'ZZLXY(IKU)',ZZLXY(IKU)
135    STOP
136 END IF   
137 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
138 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
139 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
140 !
141 ZEPS3=(PZ-ZZLXY(IK))/(ZZLXY(IK+1)-ZZLXY(IK))
142 !
143 !------------------------------------------------------------------------------
144 !
145 !*    2. INTERPOLATION DES CHAMPS
146 !
147 PRESX=  ZEPS3 *                                                             &
148       (  ZEPS2*(ZEPS1*(PAX(II+1,IJ+1,IK+1))+(1-ZEPS1)*(PAX(II,IJ+1,IK+1)))  &
149        + (1-ZEPS2)*(ZEPS1*(PAX(II+1,IJ,IK+1))+(1-ZEPS1)*(PAX(II,IJ,IK+1)))  &
150       )                                                                     &    
151       + (1-ZEPS3) *                                                         &
152       (  ZEPS2*(ZEPS1*(PAX(II+1,IJ+1,IK))+(1-ZEPS1)*(PAX(II,IJ+1,IK)))      &
153        + (1-ZEPS2)*(ZEPS1*(PAX(II+1,IJ,IK))+(1-ZEPS1)*(PAX(II,IJ,IK)))      &
154       )
155 !
156 PRESY=  ZEPS3 *                                                             &
157       (  ZEPS2*(ZEPS1*(PAY(II+1,IJ+1,IK+1))+(1-ZEPS1)*(PAY(II,IJ+1,IK+1)))  &
158        + (1-ZEPS2)*(ZEPS1*(PAY(II+1,IJ,IK+1))+(1-ZEPS1)*(PAY(II,IJ,IK+1)))  &
159       )                                                                     &    
160       + (1-ZEPS3) *                                                         &
161       (  ZEPS2*(ZEPS1*(PAY(II+1,IJ+1,IK))+(1-ZEPS1)*(PAY(II,IJ+1,IK)))      &
162        + (1-ZEPS2)*(ZEPS1*(PAY(II+1,IJ,IK))+(1-ZEPS1)*(PAY(II,IJ,IK)))      &
163       )
164 !
165 PRESZ=  ZEPS3 *                                                             &
166       (  ZEPS2*(ZEPS1*(PAZ(II+1,IJ+1,IK+1))+(1-ZEPS1)*(PAZ(II,IJ+1,IK+1)))  &
167        + (1-ZEPS2)*(ZEPS1*(PAZ(II+1,IJ,IK+1))+(1-ZEPS1)*(PAZ(II,IJ,IK+1)))  &
168       )                                                                     &    
169       + (1-ZEPS3) *                                                         &
170       (  ZEPS2*(ZEPS1*(PAZ(II+1,IJ+1,IK))+(1-ZEPS1)*(PAZ(II,IJ+1,IK)))      &
171        + (1-ZEPS2)*(ZEPS1*(PAZ(II+1,IJ,IK))+(1-ZEPS1)*(PAZ(II,IJ,IK)))      &
172       )
173 IF (OTRAJ_GROUP) THEN
174   PRESCHAMP=  ZEPS3 *                                                         &
175         (  ZEPS2*(ZEPS1*(PCHAMP(II+1,IJ+1,IK+1))+(1-ZEPS1)*(PCHAMP(II,IJ+1,IK+1)))  &
176          + (1-ZEPS2)*(ZEPS1*(PCHAMP(II+1,IJ,IK+1))+(1-ZEPS1)*(PCHAMP(II,IJ,IK+1)))  &
177         )                                                                     &    
178         + (1-ZEPS3) *                                                         &
179         (  ZEPS2*(ZEPS1*(PCHAMP(II+1,IJ+1,IK))+(1-ZEPS1)*(PCHAMP(II,IJ+1,IK)))      &
180          + (1-ZEPS2)*(ZEPS1*(PCHAMP(II+1,IJ,IK))+(1-ZEPS1)*(PCHAMP(II,IJ,IK)))      &
181         )
182 ENDIF
183 !
184 !------------------------------------------------------------------------------
185 !
186 !
187 END SUBROUTINE INTERPXYZ