Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

fiasco.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: 17.2 KB

Line
1	SUBROUTINE fiasco (poutp, ptlon, ptlat, ptsurf,
2	$ ktmsk, ktlon, ktlat, cdtper, ktper,
3	$ pinpt, pslon, pslat, pssurf,
4	$ ksmsk, kslon, kslat, cdsper, ksper,
5	$ plonz, platz, psgrb, psgra, psfrb, psfra,
6	$ kfield, cdint, cdtyp, cddim)
7	C****
8	C *****************************
9	C * OASIS ROUTINE - LEVEL 3 *
10	C * ------------- ------- *
11	C *****************************
12	C
13	C**** fiasco - Interpolation driver
14	C
15	C Purpose:
16	C -------
17	C Interpolation from source grid to target grid
18	C Use Fast Scalar Interpolator (FSCINT) package from Yves Chartier
19	C or ANAIS software from Olivier Thual et al.
20	C
21	C** Interface:
22	C ---------
23	C CALL fiasco (poutp, ptlon, ptlat, ptsurf,
24	C ktmsk, ktlon, ktlat, cdtper, ktper,
25	C pinpt, pslon, pslat, pssurf,
26	C ksmsk, kslon, kslat, cdsper, ksper,
27	C plonz, platz, psgrb, psgra, psfrb, psfra,
28	C kfield, cdint, cdtyp, cddim)*
29	C
30	C Input:
31	C -----
32	C ptlon : longitudes of target grid (real 2D)
33	C ptlat : latitudes of target grid (real 2D)
34	C ptsurf : target grid square surfaces (real 2D)
35	C ktmsk : mask of target grid (integer 2D)
36	C ktlon : number of longitudes for target grid
37	C ktlat : number of latitudes for target grid
38	C cdtper : target grid periodicity
39	C ktper : number of overlapped points for target grid
40	C pinpt : input field on source grid (real 2D)
41	C pslon : longitudes of source grid (real 2D)
42	C pslat : latitudes of source grid (real 2D)
43	C pssurf : source grid square surfaces (real 2D)
44	C ksmsk : mask of source grid (integer 2D)
45	C kslon : number of longitudes for source grid
46	C kslat : number of latitudes for source grid
47	C cdsper : source grid periodicity
48	C ksper : number of overlapped points for source grid
49	C kfield : current field number
50	C cdint : type of interpolation to be performed
51	C cdtyp : type of source grid
52	C cddim : type of field (scalar or vector)
53	C
54	C Output:
55	C ------
56	C poutp : output field on target grid (real 2D)
57	C
58	C Workspace:
59	C ---------
60	C These are work arrays passed as arguments :
61	C plonz, platz, psgrb, psgra, psfrb, psfra
62	C
63	C Externals:
64	C ---------
65	C rgoptc, rgoptr, namset, nagset, cxgaig, igscint, namsst,
66	C nagsst, naflux.
67	C
68	C Reference:
69	C ---------
70	C See OASIS manual (1995)
71	C
72	C History:
73	C -------
74	C Version Programmer Date Description
75	C ------- ---------- ---- -----------
76	C 1.0 L. Terray 94/01/01 created
77	C 1.1 L. Terray 94/08/01 modified: interpolation differ
78	C for scalar or vector in fscint
79	C 2.0 L. Terray 95/12/15 modified: new structure
80	C 2.2 L. Terray 97/11/28 added: grid type A,B,L in fscint
81	C 2.3 L. Terray 99/03/01 corrected: calls to na(g-m)sst with
82	C inclusion of pointer ipdeb
83	C 2.3 S. Valcke 99/04/30 added: printing levels
84	C 2.3 L. Terray 99/08/11 modified: fscint extrapolation
85	C for periodic Z grids
86	C 2.3 L. Terray 99/09/15 changed periodicity variables
87	C
88	C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
89	C
90	C* ---------------- Include files and USE of modules---------------------------
91	C
92	USE mod_kinds_oasis
93	USE mod_anais
94	USE mod_parameter
95	USE mod_timestep
96	USE mod_unit
97	USE mod_printing
98	CAC
99	USE mod_intbi
100	USE mod_intlin
101	CAC
102	C
103	C* ---------------------------- Argument declarations -------------------
104	C
105	REAL (kind=ip_realwp_p) poutp(ktlon,ktlat),
106	$ ptlon(ktlon,ktlat), ptlat(ktlon,ktlat)
107	REAL (kind=ip_realwp_p) pinpt(kslon,kslat),
108	$ pslon(kslon,kslat), pslat(kslon,kslat)
109	REAL (kind=ip_realwp_p) pssurf(kslon,kslat),
110	$ ptsurf(ktlon,ktlat)
111	REAL (kind=ip_realwp_p) plonz(ktlon), platz(ktlat),
112	$ psgrb(kslon,kslat),
113	$ psgra(ktlon,ktlat), psfrb(kslon,kslat), psfra(ktlon,ktlat)
114	INTEGER (kind=ip_intwp_p) ktmsk(ktlon,ktlat), ksmsk(kslon,kslat)
115	CHARACTER*8 cdsper, cdtper
116	C
117	C* ---------------------------- Local declarations ----------------------
118	C
119	CHARACTER*1 cdtyp, clref, cltyp
120	CHARACTER*6 cddim, cldim
121	CHARACTER*8 cdint, clord, clflg, clind
122	LOGICAL llinit, llsym
123	C
124	C* ---------------------------- Poema verses ----------------------------
125	C
126	C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
127	C
128	C* 1. Initialization
129	C --------------
130	C
131	IF (nlogprt .GE. 2) THEN
132	WRITE (UNIT = nulou,FMT = *) ' '
133	WRITE (UNIT = nulou,FMT = *) ' '
134	WRITE (UNIT = nulou,FMT = *)
135	$ ' ROUTINE fiasco - Level 3'
136	WRITE (UNIT = nulou,FMT = *)
137	$ ' ************ *****'
138	WRITE (UNIT = nulou,FMT = *) ' '
139	WRITE (UNIT = nulou,FMT = *)
140	$ ' Interpolation package driving routine'
141	WRITE (UNIT = nulou,FMT = *) ' '
142	WRITE (UNIT = nulou,FMT = *) ' '
143	ENDIF
144	C
145	C* Assign local value to grid and field type variables to deal
146	C with FSCINT set-up
147	C
148	cltyp = cdtyp
149	cldim = cddim
150	C
151	C * Set up interpolation method
152	C
153	IF (TRIM(cdint) == 'BICUBIPS' .OR. TRIM(cdint) == 'BILINIPS') THEN
154	IF ( cltyp == 'Z' .AND. TRIM(cdsper) == 'P' ) THEN
155	plonz = pslon ( :, kslat/2)
156	platz = pslat ( 1, :)
157	clflg='IPSL'
158	ELSE
159	WRITE (UNIT = nulou,FMT = *) ' '
160	WRITE (UNIT = nulou,FMT = *)
161	$ 'BICUBIPS and BILINIPS works ONLY',
162	$ ' with Z grid and P periodicity'
163	WRITE (UNIT = nulou,FMT = *) ' '
164	CALL HALTE ('STOP in fiasco')
165	ENDIF
166	IF ( cldim == 'VECTOR' ) THEN
167	IF (TRIM(cdint) == 'BICUBIPS') THEN
168	WRITE(UNIT = nulou, FMT= *)
169	$ 'Using intbi, vector mode'
170	CALL intbi ( pinpt, plonz, platz, kslon, kslat
171	$ , poutp, ptlon, ptlat, ktlon*ktlat
172	$ , .TRUE.)
173	ENDIF
174	IF (TRIM(cdint) == 'BILINIPS') THEN
175	WRITE(UNIT = nulou, FMT= *)
176	$ 'Using intlin, vector mode'
177	CALL intlin ( pinpt, plonz, platz, kslon, kslat
178	$ , poutp, ptlon, ptlat, ktlon*ktlat
179	$ , .TRUE.)
180	ENDIF
181	ELSE
182	IF (TRIM(cdint) == 'BICUBIPS') THEN
183	WRITE(UNIT = nulou, FMT= *)
184	$ 'Using intbi, scalar mode'
185	CALL intbi ( pinpt, plonz, platz, kslon, kslat
186	$ , poutp, ptlon, ptlat, ktlon*ktlat
187	$ , .FALSE.)
188	ENDIF
189	IF (TRIM(cdint) == 'BILINIPS') THEN
190	WRITE(UNIT = nulou, FMT= *)
191	$ 'Using intlin, scalar mode'
192	CALL intlin ( pinpt, plonz, platz, kslon, kslat
193	$ , poutp, ptlon, ptlat, ktlon*ktlat
194	$ , .FALSE.)
195	ENDIF
196	ENDIF
197	RETURN
198	C
199	ELSE IF (TRIM(cdint) .EQ. 'BILINEAR') THEN
200	clflg = 'FSCINT'
201	CALL rgoptc('INTERP', 'LINEAIR', .TRUE.)
202	ELSEIF (TRIM(cdint) .EQ. 'BICUBIC') THEN
203	clflg = 'FSCINT'
204	CALL rgoptc('INTERP', 'CUBIQUE', .TRUE.)
205	ELSEIF (TRIM(cdint) .eq. 'NNEIBOR') THEN
206	clflg = 'FSCINT'
207	CALL rgoptc('INTERP', 'VOISIN', .TRUE.)
208	ELSEIF (TRIM(cdint) .EQ. 'SURFMESH') THEN
209	clflg = 'ANAISM'
210	llinit = linit(kfield)
211	ELSEIF (TRIM(cdint) .EQ. 'GAUSSIAN') THEN
212	clflg = 'ANAISG'
213	llinit = linit(kfield)
214	ELSE
215	WRITE (UNIT = nulou,FMT = ) ' WARNING*'
216	WRITE (UNIT = nulou,FMT = *)
217	$ ' ===>>>> : unknown kind of interpolator'
218	WRITE (UNIT = nulou,FMT = *)
219	$ ' ======= ============'
220	WRITE (UNIT = nulou,FMT = *)
221	$ ' --->>> CDINT = ',cdint
222	WRITE (UNIT = nulou,FMT = *) ' '
223	WRITE (UNIT = nulou,FMT = *)
224	$ ' We STOP!!! check variable syntax '
225	WRITE (UNIT = nulou,FMT = *) ' '
226	CALL HALTE ('STOP in fiasco')
227	ENDIF
228	IF (clflg .eq. 'FSCINT') THEN
229	CALL rgoptc('INTERP', clord, .FALSE.)
230	ELSEIF (clflg .EQ. 'ANAISM') THEN
231	clord = 'SURFMESH'
232	ELSEIF (clflg .EQ. 'ANAISG') THEN
233	clord = 'GAUSSIAN'
234	ELSE
235	clord ='UNKNOWN'
236	ENDIF
237	IF (nlogprt .GE. 2) THEN
238	CALL prcout
239	$ ('Current interpolator is clord = ', clord, 1)
240	ENDIF
241	C
242	C* set up value to be given in case of extrapolation in fscint
243	C
244	IF (clflg .eq. 'FSCINT') THEN
245	IF (cltyp .EQ. 'Z') THEN
246	C
247	C* No extrapolation for Z grids (global grid)
248	C
249	CALL rgoptc ('EXTRAP', 'OUI', .TRUE.)
250	ELSE
251	CALL rgoptc ('EXTRAP', 'VOISIN', .TRUE.)
252	ENDIF
253	CALL rgoptc ('EXTRAP', clind , .false.)
254	IF (nlogprt .GE. 2) THEN
255	CALL prcout
256	$ ('Current extrapolation in FSCINT is clind =', clind, 1)
257	ENDIF
258	ENDIF
259	C
260	C* set up grid descriptors for fscint interpolator
261	C
262	IF (clflg .eq. 'FSCINT') THEN
263	IF (cltyp .EQ. 'Z' .OR. cltyp .EQ. 'Y') THEN
264	C
265	C* here ig* descriptors are arbitrary (date of louis XVI's death)
266	C
267	ig1 = 21
268	ig2 = 0
269	ig3 = 1
270	ig4 = 1793
271	zdllg = 1.0
272	zdllt = 1.0
273	zlgin = 0.0
274	zltin = 0.0
275	clref = 'L'
276	CALL cxgaig (clref, ig1z, ig2z, ig3z, ig4z,
277	$ zltin, zlgin, zdllt, zdllg)
278	DO 110 ji = 1, kslon
279	plonz(ji) = pslon(ji,1) + 1.0
280	110 CONTINUE
281	DO 120 jj = 1, kslat
282	platz(jj) = pslat(1,jj) + 1.0
283	120 CONTINUE
284	ELSE IF (cltyp .EQ. 'G'.or. cltyp .EQ. 'A'
285	$ .OR. cltyp .EQ. 'B') THEN
286	ig1 = 0
287	ig2 = 0
288	ig3 = 0
289	ig4 = 0
290	llsym = .TRUE.
291	ELSE IF (cltyp .EQ. 'L') THEN
292	zdllt = pslat(1,2) - pslat(1,1)
293	zdllg = pslon(2,1) - pslon(1,1)
294	zlgin = pslon(1,1)
295	zltin = pslat(1,1)
296	CALL cxgaig (cltyp, ig1, ig2, ig3, ig4,
297	$ zltin, zlgin, zdllt, zdllg)
298	llsym = .TRUE.
299	ELSE
300	WRITE (UNIT = nulou,FMT = ) ' WARNING*'
301	WRITE (UNIT = nulou,FMT = *)
302	$ ' ===>>> : grid type not implemented yet'
303	WRITE (UNIT = nulou,FMT = *)
304	$ ' ====== ==== =========== '
305	WRITE (UNIT = nulou,FMT = *)
306	$ ' --->>> cltyp = ',
307	$ cltyp
308	WRITE (UNIT = nulou,FMT = *) ' '
309	WRITE (UNIT = nulou,FMT = *)
310	$ ' We STOP !!! Check source grid type'
311	WRITE (UNIT = nulou,FMT = *) ' '
312	CALL HALTE ('STOP in fiasco')
313	ENDIF
314	ELSEIF (clflg .eq. 'ANAISM') THEN
315	IF (llinit) THEN
316	C
317	C* Mask values for both models
318	C
319	ismsq = 1
320	itmsq = 1
321	C
322	C* If no intersection, give value imesh to nmesh array element
323	C
324	imesh = 0
325	iwlun = nulcc
326	C
327	C* Get pointers for ANAISM interpolator parameters
328	C
329	ipdeb = (naismfl(kfield)-1)ig_maxwoaig_maxgrd+1
330	isdeb = (naismfl(kfield)-1)*ig_maxgrd+1
331	C
332	C* Set up surfmesh interpolation
333	C
334	CALL namset (pslon, pslat, ksmsk, ismsq, kslon, kslat,
335	$ cdsper, ksper,
336	$ ptlon, ptlat, ktmsk, itmsq, ktlon, ktlat,
337	$ cdtper, ktper,
338	$ cltyp,
339	C
340	C* Specific ANAISM parameters
341	C
342	$ amint(ipdeb),
343	$ nmint(ipdeb),
344	$ naismvoi(kfield), psgrb, psgra,
345	$ niwtm(kfield), iwlun, imesh,
346	$ nmesh(isdeb),
347	$ naismfl(kfield))
348	C
349	C* Switch initialization flag
350	C
351	linit(kfield) = .FALSE.
352	ENDIF
353	ELSEIF (clflg .eq. 'ANAISG') THEN
354	IF (llinit) THEN
355	ismsq = 1
356	itmsq = 1
357	iwlun = nulgg
358	C
359	C* Get pointers for Anaisg interpolator parameters
360	C
361	ipdeb = (naisgfl(kfield)-1)ig_maxnoaig_maxgrd+1
362	CALL nagset (pslon, pslat, pssurf,
363	$ ksmsk, ismsq, kslon, kslat,
364	$ ptlon, ptlat, ptsurf,
365	$ ktmsk, itmsq, ktlon, ktlat,
366	$ cltyp,
367	C
368	C* Specific Anaisg parameters
369	C
370	$ agint(ipdeb),
371	$ ngint(ipdeb),
372	$ naisgvoi(kfield), psfrb, psfra,
373	$ varmul(kfield), niwtg(kfield), iwlun,
374	$ naisgfl(kfield))
375	C
376	C* Switch off initialization flag
377	C
378	linit(kfield) = .FALSE.
379	ENDIF
380	ELSE
381	WRITE (UNIT = nulou,FMT = ) ' WARNING*'
382	WRITE (UNIT = nulou,FMT = *)
383	$ ' ===>>> : unknown interpolation package'
384	WRITE (UNIT = nulou,FMT = *)
385	$ ' ====== ======= ======='
386	WRITE (UNIT = nulou,FMT = *)
387	$ ' --->>> clflg = ',
388	$ clflg
389	WRITE (UNIT = nulou,FMT = *) ' '
390	WRITE (UNIT = nulou,FMT = *)
391	$ ' We STOP !!! Check interpolation package'
392	WRITE (UNIT = nulou,FMT = *) ' '
393	CALL HALTE ('STOP in fiasco')
394	ENDIF
395	C
396	C
397	C* 2. Interpolation
398	C -------------
399	C
400	IF (clflg .eq. 'FSCINT') THEN
401	IF (cltyp .eq. 'Z' .OR. cltyp .EQ. 'Y') THEN
402	CALL igscint (poutp, ktlon, ktlat, ptlat, ptlon,
403	$ pinpt, kslon, kslat, cltyp, clref,
404	$ ig1z, ig2z, ig3z, ig4z, .true.,
405	$ plonz, platz, cldim)
406	ELSEIF (cltyp .eq. 'G'.or. cltyp .eq. 'A'
407	$ .OR. cltyp .eq. 'B'.or. cltyp .eq. 'L') THEN
408	CALL rgscint (poutp, ktlon, ktlat, ptlat, ptlon,
409	$ pinpt, kslon, kslat, cltyp,
410	$ ig1, ig2, ig3, ig4, llsym, cldim)
411	ELSE
412	WRITE (UNIT = nulou,FMT = ) ' WARNING*'
413	WRITE (UNIT = nulou,FMT = *)
414	$ ' ===>>> grid type not implemented yet'
415	WRITE (UNIT = nulou,FMT = *)
416	$ ' ====== ==== =========== '
417	WRITE (UNIT = nulou,FMT = *)
418	$ ' --->>> cltyp = ',
419	$ cltyp
420	WRITE (UNIT = nulou,FMT = *) ' '
421	WRITE (UNIT = nulou,FMT = *)
422	$ ' WE STOP !!! check source grid type'
423	WRITE (UNIT = nulou,FMT = *) ' '
424	CALL HALTE ('STOP in fiasco')
425	ENDIF
426	ELSEIF (clflg .eq. 'ANAISM') THEN
427	itmsq = 1
428	ipdeb = (naismfl(kfield)-1)ig_maxwoaig_maxgrd+1
429	CALL namsst (poutp, ktmsk, itmsq, ktlon, ktlat,
430	$ amint(ipdeb), nmint(ipdeb), naismvoi(kfield),
431	$ pinpt, kslon, kslat)
432	ELSEIF (clflg .eq. 'ANAISG') THEN
433	itmsq = 1
434	ipdeb = (naisgfl(kfield)-1)ig_maxnoaig_maxgrd+1
435	CALL nagsst (poutp, ktmsk, itmsq, ktlon, ktlat,
436	$ agint(ipdeb), ngint(ipdeb), naisgvoi(kfield),
437	$ pinpt, kslon, kslat)
438	ELSE
439	WRITE (UNIT = nulou,FMT = ) ' WARNING*'
440	WRITE (UNIT = nulou,FMT = *)
441	$ ' ===>>> : unknown interpolation package'
442	WRITE (UNIT = nulou,FMT = *)
443	$ ' ====== ======= ======='
444	WRITE (UNIT = nulou,FMT = *)
445	$ ' --->>> clflg = ',
446	$ clflg
447	WRITE (UNIT = nulou,FMT = *) ' '
448	WRITE (UNIT = nulou,FMT = *)
449	$ ' We STOP !!! Check interpolation package'
450	WRITE (UNIT = nulou,FMT = *) ' '
451	CALL HALTE ('STOP in fiasco')
452	ENDIF
453	C
454	C
455	C* 3. End of routine
456	C --------------
457	C
458	IF (nlogprt .GE. 2) THEN
459	WRITE (UNIT = nulou,FMT = *) ' '
460	WRITE (UNIT = nulou,FMT = *)
461	$ ' --------- End of routine fiasco ---------'
462	CALL FLUSH (nulou)
463	ENDIF
464	RETURN
465	END

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

Context Navigation

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

Download in other formats: