Context Navigation

extrapol.F @ 1322

Last change on this file since 1322 was 1299, checked in by Laurent Fairhead, 15 years ago

Nettoyage general pour se rapprocher des normes et éviter des erreurs a la
compilation:

tous les FLOAT() sont remplacés par des REAL()
tous les STOP dans phylmd sont remplacés par des appels à abort_gcm
le common control défini dans le fichier control.h est remplacé par le module control_mod pour éviter des messages sur l'alignement des variables dans les déclarations
des $Header$ remplacés par des $Id$ pour svn

Quelques remplacements à faire ont pu m'échapper

General cleanup of the code to try and adhere to norms and to prevent some
compilation errors:

all FLOAT() instructions have been replaced by REAL() instructions
all STOP instructions in phylmd have been replaced by calls to abort_gcm
the common block control defined in the control.h file has been replaced by the control_mod to prevent compilation warnings on the alignement of declared variables
$Header$ replaced by $Id$ for svn

Some changes which should have been made might have escaped me

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 5.0 KB

Rev	Line
[630]	1	!
[1299]	2	! $Id: extrapol.F 1299 2010-01-20 14:27:21Z fairhead $
[630]	3	!
	4	C
	5	C
	6	SUBROUTINE extrapol (pfild, kxlon, kylat, pmask,
	7	. norsud, ldper, knbor, pwork)
	8	IMPLICIT none
	9	c
	10	c OASIS routine (Adaptation: Laurent Li, le 14 mars 1997)
	11	c Fill up missed values by using the neighbor points
	12	c
	13	INTEGER kxlon, kylat ! longitude and latitude dimensions (Input)
	14	INTEGER knbor ! minimum neighbor number (Input)
	15	LOGICAL norsud ! True if field is from North to South (Input)
	16	LOGICAL ldper ! True if take into account the periodicity (Input)
	17	REAL pmask ! mask value (Input)
	18	REAL pfild(kxlon,kylat) ! field to be extrapolated (Input/Output)
	19	REAL pwork(kxlon,kylat) ! working space
	20	c
	21	REAL zwmsk
	22	INTEGER incre, idoit, i, j, k, inbor, ideb, ifin, ilon, jlat
	23	INTEGER ix(9), jy(9) ! index arrays for the neighbors coordinates
	24	REAL zmask(9)
	25	C
	26	C We search over the eight closest neighbors
	27	C
	28	C j+1 7 8 9
	29	C j 4 5 6 Current point 5 --> (i,j)
	30	C j-1 1 2 3
	31	C i-1 i i+1
	32	c
	33	c
	34	IF (norsud) THEN
	35	DO j = 1, kylat
	36	DO i = 1, kxlon
	37	pwork(i,j) = pfild(i,kylat-j+1)
	38	ENDDO
	39	ENDDO
	40	DO j = 1, kylat
	41	DO i = 1, kxlon
	42	pfild(i,j) = pwork(i,j)
	43	ENDDO
	44	ENDDO
	45	ENDIF
	46	c
	47	incre = 0
	48	c
	49	DO j = 1, kylat
	50	DO i = 1, kxlon
	51	pwork(i,j) = pfild(i,j)
	52	ENDDO
	53	ENDDO
	54	c
	55	C* To avoid problems in floating point tests
	56	zwmsk = pmask - 1.0
	57	c
	58	200 CONTINUE
	59	incre = incre + 1
	60	DO 99999 j = 1, kylat
	61	DO 99999 i = 1, kxlon
	62	IF (pfild(i,j).GT. zwmsk) THEN
	63	pwork(i,j) = pfild(i,j)
	64	inbor = 0
	65	ideb = 1
	66	ifin = 9
	67	C
	68	C* Fill up ix array
	69	ix(1) = MAX (1,i-1)
	70	ix(2) = i
	71	ix(3) = MIN (kxlon,i+1)
	72	ix(4) = MAX (1,i-1)
	73	ix(5) = i
	74	ix(6) = MIN (kxlon,i+1)
	75	ix(7) = MAX (1,i-1)
	76	ix(8) = i
	77	ix(9) = MIN (kxlon,i+1)
	78	C
	79	C* Fill up iy array
	80	jy(1) = MAX (1,j-1)
	81	jy(2) = MAX (1,j-1)
	82	jy(3) = MAX (1,j-1)
	83	jy(4) = j
	84	jy(5) = j
	85	jy(6) = j
	86	jy(7) = MIN (kylat,j+1)
	87	jy(8) = MIN (kylat,j+1)
	88	jy(9) = MIN (kylat,j+1)
	89	C
	90	C* Correct latitude bounds if southernmost or northernmost points
	91	IF (j .EQ. 1) ideb = 4
	92	IF (j .EQ. kylat) ifin = 6
	93	C
	94	C* Account for periodicity in longitude
	95	C
	96	IF (ldper) THEN
	97	IF (i .EQ. kxlon) THEN
	98	ix(3) = 1
	99	ix(6) = 1
	100	ix(9) = 1
	101	ELSE IF (i .EQ. 1) THEN
	102	ix(1) = kxlon
	103	ix(4) = kxlon
	104	ix(7) = kxlon
	105	ENDIF
	106	ELSE
	107	IF (i .EQ. 1) THEN
	108	ix(1) = i
	109	ix(2) = i + 1
	110	ix(3) = i
	111	ix(4) = i + 1
	112	ix(5) = i
	113	ix(6) = i + 1
	114	ENDIF
	115	IF (i .EQ. kxlon) THEN
	116	ix(1) = i -1
	117	ix(2) = i
	118	ix(3) = i - 1
	119	ix(4) = i
	120	ix(5) = i - 1
	121	ix(6) = i
	122	ENDIF
	123	C
	124	IF (i .EQ. 1 .OR. i .EQ. kxlon) THEN
	125	jy(1) = MAX (1,j-1)
	126	jy(2) = MAX (1,j-1)
	127	jy(3) = j
	128	jy(4) = j
	129	jy(5) = MIN (kylat,j+1)
	130	jy(6) = MIN (kylat,j+1)
	131	C
	132	ideb = 1
	133	ifin = 6
	134	IF (j .EQ. 1) ideb = 3
	135	IF (j .EQ. kylat) ifin = 4
	136	ENDIF
	137	ENDIF ! end for ldper test
	138	C
	139	C* Find unmasked neighbors
	140	C
	141	DO 230 k = ideb, ifin
	142	zmask(k) = 0.
	143	ilon = ix(k)
	144	jlat = jy(k)
	145	IF (pfild(ilon,jlat) .LT. zwmsk) THEN
	146	zmask(k) = 1.
	147	inbor = inbor + 1
	148	ENDIF
	149	230 CONTINUE
	150	C
	151	C* Not enough points around point P are unmasked; interpolation on P
	152	C will be done in a future call to extrap.
	153	C
	154	IF (inbor .GE. knbor) THEN
	155	pwork(i,j) = 0.
	156	DO k = ideb, ifin
	157	ilon = ix(k)
	158	jlat = jy(k)
	159	pwork(i,j) = pwork(i,j)
[1299]	160	$ + pfild(ilon,jlat) * zmask(k)/ REAL(inbor)
[630]	161	ENDDO
	162	ENDIF
	163	C
	164	ENDIF
	165	99999 CONTINUE
	166	C
	167	C* 3. Writing back unmasked field in pfild
	168	C ------------------------------------
	169	C
	170	C* pfild then contains:
	171	C - Values which were not masked
	172	C - Interpolated values from the inbor neighbors
	173	C - Values which are not yet interpolated
	174	C
	175	idoit = 0
	176	DO j = 1, kylat
	177	DO i = 1, kxlon
	178	IF (pwork(i,j) .GT. zwmsk) idoit = idoit + 1
	179	pfild(i,j) = pwork(i,j)
	180	ENDDO
	181	ENDDO
	182	c
	183	IF (idoit .ne. 0) GOTO 200
	184	ccc PRINT*, "Number of extrapolation steps incre =", incre
	185	c
	186	IF (norsud) THEN
	187	DO j = 1, kylat
	188	DO i = 1, kxlon
	189	pwork(i,j) = pfild(i,kylat-j+1)
	190	ENDDO
	191	ENDDO
	192	DO j = 1, kylat
	193	DO i = 1, kxlon
	194	pfild(i,j) = pwork(i,j)
	195	ENDDO
	196	ENDDO
	197	ENDIF
	198	c
	199	RETURN
	200	END

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format