source: CPL/oasis3/trunk/src/mod/oasis3/src/mod_findpos.f @ 1677

Last change on this file since 1677 was 1677, checked in by aclsce, 12 years ago

Imported oasis3 (tag ipslcm5a) from cvs server to svn server (igcmg project).

File size: 5.7 KB
Line 
1      MODULE mod_findpos
2      CONTAINS
3      SUBROUTINE findpos ( px, py, kpts, pax, pay, ki, kj,
4     $                 kndx, kndy, plx, ply, ldnorth, ldsouth)
5C****
6C               *****************************
7C               * OASIS ROUTINE  -  LEVEL 3 *
8C               * -------------     ------- *
9C               *****************************
10C
11C**** *pos* - Computes position of points in a longitude-latitude grid
12C
13C     Purpose:
14C     -------
15C**   Interface:
16C     ---------
17C       *CALL*  *findpos* ( px, py, kpts, pax, pay, ki, kj,
18C     $                 kndx, kndy, plx, ply, ldnorth, ldsouth) 
19C
20C     Input:
21C     -----
22C                px      : longitudes of target grid
23C                py      : latitudes of target grid
24C                kpts    : dimension of target grid
25C                pax     : longitudes of source grid
26C                pay     : latitudes of source grid
27C                ki, kj  : dimension of source grid
28C
29C     Output:
30C     ------
31C                plx     : extended longitudes of source grid
32C                ply     : extended latitudes of source grid
33C                ldnorth : TRUE if north pole is a point of the grid
34C                ldsouth : TRUE if south pole is a point of the grid
35C
36C     Workspace:
37C     ---------
38C     Local variables
39C                ji, jj, jn, ilo, ihi, zxmin, zxmax, zymin, zymax
40C                iimin, iimax, ijmin, ijmax
41C
42C     Externals:
43C     ---------
44C     None
45C
46C     Reference:
47C     ---------
48C     See OASIS manual (1995)
49C
50C     History:
51C     -------
52C       Version   Programmer     Date      Description
53C       -------   ----------     ----      ----------- 
54C       2.0       O. Marti       96/07/15  Created
55C
56C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
57C
58C 
59      IMPLICIT none
60C
61C* ---------------------------- Argument declarations -------------------
62C
63      INTEGER, INTENT ( in) :: kpts, ki, kj
64      REAL, DIMENSION ( kpts), INTENT ( inout) :: px, py
65      INTEGER, DIMENSION ( kpts), INTENT ( out) :: kndx, kndy
66      REAL, DIMENSION ( ki), INTENT ( in) :: pax
67      REAL, DIMENSION ( kj), INTENT ( in) :: pay
68      REAL, DIMENSION ( -1: ki + 2), INTENT ( out) :: plx
69      REAL, DIMENSION ( -1: kj + 2), INTENT ( out) :: ply
70      LOGICAL, INTENT ( out) :: ldnorth, ldsouth
71C
72C* ---------------------------- Local declarations ----------------------
73C
74      INTEGER :: ji, jj, jn
75      INTEGER :: ilo, ihi
76      REAL    :: zxmin, zxmax, zymin, zymax
77      INTEGER :: iimin, iimax, ijmin, ijmax
78C
79C* ---------------------------- Poema verses ----------------------------
80C
81C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
82C
83C*    1. Initialization
84C        --------------
85C
86C
87C* Extension of source grid
88C
89C For longitudinal, by periodicity
90C
91      plx ( 1: ki) = pax ( 1: ki)
92      plx (  -1   ) = plx ( ki - 1) - 360.0
93      plx (   0   ) = plx ( ki    ) - 360.0
94      plx ( ki + 1) = plx (   1   ) + 360.0
95      plx ( ki + 2) = plx (   2   ) + 360.0
96C
97C  For latitudinal
98C
99      ply ( 1: kj) = pay ( 1: kj)
100C
101C* Point at the south pole ?
102C
103      IF ( abs ( ply ( 1) - ( -90.0)) .LT. 1.0e-1) THEN
104          ply ( -1) = -180.0 - ply ( 3)
105          ply (  0) = -180.0 - ply ( 2) 
106          ply (  1) = -90.0
107          ldsouth = .TRUE.
108      ELSE
109          ply ( -1) = -180.0 - ply ( 2)
110          ply (  0) = -180.0 - ply ( 1)
111          ldsouth = .FALSE.
112      ENDIF
113C
114C* Point at the north pole ?
115C
116      IF ( abs ( ply ( kj) - 90.0 ) .LT. 1.0E-1 ) THEN
117          ply ( kj    ) = 90.0
118          ply ( kj + 1) = 180.0 - ply ( kj - 1)
119          ply ( kj + 2) = 180.0 - ply ( kj - 2)
120          ldnorth = .TRUE.
121      ELSE
122          ply ( kj + 1) = 180.0 - ply ( kj)
123          ply ( kj + 2) = 180.0 - ply ( kj - 1)
124          ldnorth = .FALSE.
125      ENDIF
126C
127C Verif
128      zxmin = MINVAL ( plx) ; zxmax = MAXVAL ( plx)
129      zymin = MINVAL ( ply) ; zymax = MAXVAL ( ply)
130C
131C* Set TARGET grid within the limits of the source grid
132C
133      zxmin = plx ( 1)
134      zxmax = plx ( 1) + 360.0
135      DO jn = 1, kpts
136        IF ( px ( jn) .LT. zxmin) THEN
137            px ( jn) = px ( jn) + 360.0
138        ELSE IF ( px ( jn) .GT. zxmax) THEN
139            px ( jn) = px ( jn) - 360.0
140        ENDIF
141      END DO
142C
143C 2. Seek the POSITION of TARGET point by bisection
144C
145      iimin = ki + 10 ; iimax = -10 ; ijmin = kj + 10 ; ijmax = -10
146C
147      DO jn = 1, kpts
148C
149C* Longitude
150C
151        ilo = -1
152        ihi = ki + 1
153        ji  = ( ilo + ihi) * 0.5
154 1111   CONTINUE
155        IF ( ( ihi - ilo) .GT. 1 ) THEN
156            IF ( px ( jn) .LE. plx ( ji) ) THEN
157                ihi = ji
158            ELSE
159                ilo = ji
160            ENDIF
161            ji = ( ilo + ihi) * 0.5 
162            GO TO 1111
163        ENDIF
164        kndx ( jn) = ilo
165C
166C* Latitude
167C
168        ilo = -1
169        ihi = kj + 1
170        jj  = ( ilo + ihi) * 0.5
171 2222   CONTINUE
172        IF ( ( ihi - ilo) .GT. 1 ) THEN
173            IF ( py ( jn) .LE. ply ( jj) ) THEN
174                ihi = jj
175            ELSE
176                ilo = jj
177            ENDIF
178            jj = ( ilo + ihi) * 0.5
179            GO TO 2222
180        ENDIF
181        kndy ( jn) = ilo
182C
183        iimin = min ( iimin, kndx ( jn))
184        iimax = max ( iimax, kndx ( jn))
185        ijmin = min ( ijmin, kndy ( jn))
186        ijmax = max ( ijmax, kndy ( jn))
187      END DO
188C
189      IF ( iimin .LT.  0 ) 
190     $    CALL halte ('STOP in pos iimin incorrect')
191      IF ( iimax .LT.  0 ) 
192     $    CALL halte ('STOP in pos iimax incorrect')
193      IF ( ijmin .LT.  0 ) 
194     $    CALL halte ('STOP in pos ijmin incorrect')
195      IF ( ijmax .LT.  0 ) 
196     $    CALL halte ('STOP in pos ijmax incorrect')
197C
198C*    3. End of routine
199C        --------------
200      RETURN
201      END SUBROUTINE findpos
202      END MODULE mod_findpos
Note: See TracBrowser for help on using the repository browser.