source: CPL/oasis3/trunk/src/mod/oasis3/src/inigrd.f

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

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

File size: 16.7 KB
Line 
1      SUBROUTINE inigrd
2C****
3C               *****************************
4C               * OASIS ROUTINE  -  LEVEL 0 *
5C               * -------------     ------- *
6C               *****************************
7C
8C**** *inigrd*  - Initialize gcm's grids
9C
10C     Purpose:
11C     -------
12C     Get atmospheric and oceanic grids and masks as well as scale factors
13C     for interpolation purposes.
14C
15C     ***NOTE***
16C     ----------
17C     Grids, masks and scale factors must be ordered with the OASIS convention.
18C     They go from South to North and from Greenwhich to the east.
19C     EXCEPTION: array maskr (mask for reduced gaussian grid) must
20C     be ordered as in the atmospheric code (North to South in
21C     Arpege case).
22C
23C**   Interface:
24C     ---------
25C       *CALL*  *inigrd*
26C
27C     Input:
28C     -----
29C     None
30C
31C     Output:
32C     ------
33C     None
34C
35C     Workspace:
36C     ---------
37C     None
38C
39C     Externals:
40C     ---------
41C     locread, locrint
42C
43C     Reference:
44C     ---------
45C     See OASIS manual (1995)
46C
47C     History:
48C     -------
49C       Version   Programmer     Date      Description
50C       -------   ----------     ----      ----------- 
51C       1.0       L. Terray      94/01/01  created
52C       2.0       L. Terray      95/08/31  modified : new structure
53C       2.3       S. Valcke      99/04/30  added: printing levels
54C       2.5       S. Valcke      01/03/23  added: netCDF auxilary files
55C
56C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
57C
58C* ---------------------------- Include files ---------------------------
59C
60      USE mod_kinds_oasis
61      USE mod_parameter
62      USE mod_memory
63      USE mod_string
64      USE mod_unitncdf
65      USE mod_analysis
66      USE mod_label
67      USE mod_unit
68      USE mod_printing
69      INCLUDE 'netcdf.inc'
70C
71C* ---------------------------- Local declarations ----------------------
72C
73      CHARACTER*8 clwork
74      CHARACTER*8 clstrg,clstrg1,clstrg2,clstrg3,clstrg4
75      CHARACTER*8 clstrg5,clstrg6
76      INTEGER (kind=ip_intwp_p) il_varid
77      INTEGER (kind=ip_intwp_p) ist(2), icnt(2)
78C
79C* ---------------------------- Poema verses ----------------------------
80C
81C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
82C
83C
84C* This routine will be called only if one field (at least) goes through Oasis
85C
86      IF (ig_nfield .ne.0) THEN
87C
88C*    1. Initialization
89C        --------------
90C
91      IF (nlogprt .GE. 1) THEN
92          WRITE (UNIT = nulou,FMT = *) ' '
93          WRITE (UNIT = nulou,FMT = *) ' '
94          WRITE (UNIT = nulou,FMT = *)
95     $    '           ROUTINE inigrd  -  Level 0'
96          WRITE (UNIT = nulou,FMT = *)
97     $    '           **************     *******'
98          WRITE (UNIT = nulou,FMT = *) ' '
99          WRITE (UNIT = nulou,FMT = *)
100     $    ' get gcms grids, masks and surfaces'
101          WRITE (UNIT = nulou,FMT = *) ' '
102          WRITE (UNIT = nulou,FMT = *) ' '
103      ENDIF
104      ist(:)=0
105      icnt(:)=0
106C
107C* Initializa error codes
108C
109      infos = 0
110      iflag = 0
111C
112C
113C*    2. Get gcm's grids, masks and mesh surfaces
114C        ----------------------------------------
115C
116!$omp parallel do default (shared)
117!$omp+ private (jf,iadrold,isizold,iadrnew,isiznew)
118!$omp+ private (clwork,icount,clstrg,clstrg1)
119!$omp+ private (clstrg2,clstrg3,clstrg4,clstrg5)
120!$omp+ private (clstrg6)
121!$omp+ private (istatus,iunit1,inunit2,iunit3)
122!$omp+ private (iunit4,iunit5,iunit6,iunit7)
123!$omp+ private (iflag,ist,icnt,nc_file,n_intty)
124!$omp+ private (n_reaty,il_varid)
125!$omp+ reduction (+:infos)
126
127      DO 200 jf = 1, ig_nfield
128C
129C* Print
130C
131        IF (nlogprt .GE. 2) THEN
132            CALL prtout
133     $          ('Start reading source grid data for field nbr=',jf,2)
134        ENDIF
135C
136C* Assign local variables to array sizes and pointers
137C
138        iadrold = nadrold_grid(jf)
139        isizold = nsizold(jf)
140        iadrnew = nadrnew_grid(jf)
141        isiznew = nsiznew(jf)
142C
143C* Get the right locator string for initial grids files
144C
145        clwork = cficbf(jf)
146        icount = ilenstr(clwork,jpeight)
147        clstrg1 = clwork(1:icount)//cglonsuf
148C
149C* Read longitudes of initial grid associated to field jf
150C
151        iunit1 = nulgr
152        IF (lncdfgrd) THEN
153C Next BLOCK done only once per field
154            ist(1)=1 ; ist(2)=1
155            icnt(1)=nlonbf(jf) ; icnt(2)=nlatbf(jf)
156C
157            nc_file = nc_grdid
158            CALL hdlerr
159     $          (NF_INQ_VARID(nc_file, clstrg1, il_varid), 'inigrd')
160            CALL hdlerr
161     $          (NF_INQ_VARTYPE(nc_file, il_varid, n_reaty), 'inigrd')
162            IF (n_reaty .eq. NF_FLOAT) THEN
163                 IF ( ip_realwp_p .ne. ip_single_p ) THEN
164                    WRITE (nulou, *)
165     $             'Incoherence: Oasis is compiled in double precision'
166                    WRITE (nulou, *)
167     $   'but Grid auxiliary file contains single precision REAL'
168                    CALL HALTE('STOP in getfld')
169                 ENDIF
170               CALL hdlerr(NF_GET_VARA_REAL
171     $              (nc_file, il_varid, ist, icnt, 
172     $               xgrold(iadrold:iadrold+isizold-1)),'inigrd')
173            ELSE IF (n_reaty .eq. NF_DOUBLE) THEN
174                 IF ( ip_realwp_p .ne. ip_double_p ) THEN
175                    WRITE (nulou, *)
176     $             'Incoherence: Oasis is compiled in single precision'
177                    WRITE (nulou, *)
178     $   'but Grid auxiliary file contains double precision REAL'
179                    CALL HALTE('STOP in getfld')
180                 ENDIF               
181               CALL hdlerr(NF_GET_VARA_DOUBLE 
182     $              (nc_file, il_varid, ist, icnt, 
183     $               xgrold(iadrold:iadrold+isizold-1)),'inigrd')
184            ELSE
185                CALL prcout ('Problem with type of array =',clstrg1,1) 
186                CALL HALTE('STOP in inigrd')
187            ENDIF
188        ELSE
189            CALL locread (clstrg1,xgrold(iadrold),isizold,iunit1,iflag)
190            IF (iflag .NE. 0) THEN
191                CALL prtout('Problem in reading unit =',iunit1,2)
192                CALL prcout
193     $          ('Problem with array linked to string =',clstrg1,1)
194                infos = infos + iflag
195            ENDIF
196        ENDIF
197        IF (nlogprt .GE. 2) THEN
198            CALL prtout('Read source grid longitudes for field #',jf,1)
199        ENDIF
200C
201C* Read latitudes of initial grid associated to field jf
202C
203        clstrg2 = clwork(1:icount)//cglatsuf
204        IF (lncdfgrd) THEN
205            CALL hdlerr
206     $          (NF_INQ_VARID(nc_file, clstrg2, il_varid), 'inigrd')
207            IF (n_reaty .eq. NF_FLOAT) THEN
208                CALL hdlerr(NF_GET_VARA_REAL 
209     $              (nc_file, il_varid, ist, icnt, 
210     $              ygrold(iadrold:iadrold+isizold-1)), 'inigrd')
211            ELSE IF (n_reaty .eq. NF_DOUBLE) THEN
212                CALL hdlerr(NF_GET_VARA_DOUBLE 
213     $              (nc_file, il_varid, ist, icnt, 
214     $              ygrold(iadrold:iadrold+isizold-1)), 'inigrd')
215            ENDIF
216        ELSE
217            CALL locread(clstrg2,ygrold(iadrold),isizold,iunit1,iflag)
218            IF (iflag .NE. 0) THEN
219                CALL prtout('Problem in reading unit =',iunit1,2)
220                CALL prcout
221     $          ('Problem with array linked to string =',clstrg2,1)
222                infos = infos + iflag
223            ENDIF
224        ENDIF
225        IF (nlogprt .GE. 2) THEN
226            CALL prtout('Read source grid latitudes for field # ',jf,1)
227        ENDIF
228C
229C* Read initial mask associated to field jf
230C
231        clstrg3 = clwork(1:icount)//cmsksuf
232        iunit2 = nulma
233        IF (lncdfgrd) THEN
234            CALL hdlerr
235     $          (NF_INQ_VARID(nc_mskid, clstrg3, il_varid), 'inigrd')
236            CALL hdlerr
237     $          (NF_INQ_VARTYPE(nc_mskid, il_varid, n_intty), 'inigrd')
238            IF (n_intty .eq. NF_INT) THEN
239                CALL hdlerr(NF_GET_VARA_INT 
240     $              (nc_mskid, il_varid, ist, icnt, 
241     $              mskold(iadrold:iadrold+isizold-1)), 'inigrd')
242            ELSE
243                CALL prcout ('Problem with type of array =',clstrg3,1) 
244                CALL HALTE('STOP in inigrd')
245            ENDIF
246        ELSE
247            CALL locrint (clstrg3,mskold(iadrold),isizold,iunit2,iflag)
248            IF (iflag .NE. 0) THEN
249                CALL prtout('Problem in reading unit =',iunit2,2)
250                CALL prcout
251     $          ('Problem with array linked to string =',clstrg3,1)
252                infos = infos + iflag
253            ENDIF
254        ENDIF
255        IF (nlogprt .GE. 2) THEN
256            CALL prtout
257     $          ('Read source grid mask for field number =',jf,1)
258        ENDIF
259C
260C* Read initial mesh surface associated to field jf
261C
262        clstrg4 = clwork(1:icount)//csursuf
263        iunit3 = nulsu
264        IF (lncdfgrd) THEN
265            istatus = 
266     $          NF_INQ_VARID(nc_surid, clstrg4, il_varid)
267            IF (istatus .ne. NF_NOERR .and. lsurf(jf)) THEN
268                CALL prtout('Cannot read surface information for fld '
269     $              , jf, 1)
270                CALL HALTE ('STOP in inigrd')
271            ELSE IF (istatus .eq. NF_NOERR) THEN
272                IF (n_reaty .eq. NF_FLOAT) THEN
273                    istatus = NF_GET_VARA_REAL(nc_surid, 
274     $                   il_varid, ist, icnt, 
275     $                  surold(iadrold:iadrold+isizold-1))
276                    IF (istatus .ne. NF_NOERR .and. lsurf(jf)) THEN
277                        CALL prtout('Cannot read source grid surface
278     $ information for field ' , jf, 1)
279                        CALL HALTE ('STOP in inigrd')
280                    ENDIF
281                ELSE IF (n_reaty .eq. NF_DOUBLE) THEN
282                    istatus = NF_GET_VARA_DOUBLE(nc_surid, 
283     $                   il_varid, ist, icnt, 
284     $                  surold(iadrold:iadrold+isizold-1))
285                    IF (istatus .ne. NF_NOERR .and. lsurf(jf)) THEN
286                        CALL prtout('Cannot read source grid surface
287     $ information for field ' , jf, 1)
288                        CALL HALTE ('STOP in inigrd')
289                    ENDIF
290                ENDIF
291            ENDIF
292        ELSE
293            CALL locread (clstrg4, surold(iadrold), isizold,
294     $          iunit3,iflag)
295            IF (iflag .NE. 0 .and. lsurf(jf)) THEN
296                CALL prtout('Problem in reading unit =',iunit3,2)
297                CALL prcout('Problem with array linked to string ',
298     $              clstrg4,1)
299                infos = infos + iflag
300            ENDIF
301        ENDIF
302        IF (nlogprt .GE. 2) THEN
303            CALL prtout
304     $          ('Read source grid surface for field number ',jf,1)
305        ENDIF
306C
307C* Get the right locator string for final grids files
308C
309!$omp critical
310        clwork = cficaf(jf)
311        icount = ilenstr(clwork,jpeight)
312        clstrg5 = clwork(1:icount)//cglonsuf
313        IF (nlogprt .GE. 2) THEN
314            CALL prtout
315     $          ('Start reading target grid data for field nbr =',jf,2)
316        ENDIF
317C
318C* Read longitudes of final grid associated to field jf
319C
320        iunit4 = nulgr
321        IF (lncdfgrd) THEN
322C Next BLOCK done only once per field
323            ist(1)=1 ; ist(2)=1
324            icnt(1)=nlonaf(jf) ; icnt(2)=nlataf(jf)
325C
326            nc_file = nc_grdid
327            CALL hdlerr
328     $          (NF_INQ_VARID(nc_file, clstrg5, il_varid), 'inigrd')
329            IF (n_reaty .eq. NF_FLOAT) THEN
330                CALL hdlerr(NF_GET_VARA_REAL
331     $              (nc_file, il_varid, ist, icnt, 
332     $              xgrnew(iadrnew:iadrnew+isiznew-1)), 'inigrd')
333            ELSE IF (n_reaty .eq. NF_DOUBLE) THEN
334                CALL hdlerr(NF_GET_VARA_DOUBLE 
335     $              (nc_file, il_varid, ist, icnt, 
336     $              xgrnew(iadrnew:iadrnew+isiznew-1)), 'inigrd')
337            ENDIF
338        ELSE
339            CALL locread (clstrg5,xgrnew(iadrnew),isiznew,iunit4,iflag)
340            IF (iflag .NE. 0) THEN
341                CALL prtout('Problem in reading unit =',iunit4,2)
342                CALL prcout
343     $          ('Problem with array linked to string =',clstrg5,1)
344                infos = infos + iflag
345            ENDIF
346        ENDIF
347!$omp end critical
348        IF (nlogprt .GE. 2) THEN
349            CALL prtout
350     $          ('Read target grid longitudes for field number =',jf,1)
351        ENDIF
352C
353C* Read latitudes of final grid associated to field jf
354C
355!$omp critical
356        clstrg6 = clwork(1:icount)//cglatsuf
357        iunit5 = nulgr
358        IF (lncdfgrd) THEN
359            CALL hdlerr
360     $          (NF_INQ_VARID(nc_file, clstrg6, il_varid), 'inigrd')
361            IF (n_reaty .eq. NF_FLOAT) THEN
362                CALL hdlerr(NF_GET_VARA_REAL 
363     $              (nc_file, il_varid, ist, icnt, 
364     $              ygrnew(iadrnew:iadrnew+isiznew-1)), 'inigrd')
365            ELSE IF (n_reaty .eq. NF_DOUBLE) THEN
366                CALL hdlerr(NF_GET_VARA_DOUBLE 
367     $              (nc_file, il_varid, ist, icnt, 
368     $              ygrnew(iadrnew:iadrnew+isiznew-1)), 'inigrd')
369            ENDIF
370        ELSE
371            CALL locread (clstrg6,ygrnew(iadrnew),isiznew,iunit5,iflag)
372            IF (iflag .NE. 0) THEN
373                CALL prtout('Problem in reading unit =',nulgr,2)
374                CALL prcout
375     $          ('Problem with array linked to string =',clstrg6,1)
376                infos = infos + iflag
377            ENDIF
378        ENDIF
379!$omp end critical
380        IF (nlogprt .GE. 2) THEN
381            CALL prtout
382     $          ('Read target grid latitudes for field number =',jf,1)
383        ENDIF
384C
385C* Read final mask associated to field jf
386C
387!$omp critical
388        clstrg = clwork(1:icount)//cmsksuf
389        iunit6 = nulma
390        IF (lncdfgrd) THEN
391            CALL hdlerr
392     $          (NF_INQ_VARID(nc_mskid, clstrg, il_varid), 'inigrd')
393            IF (n_intty .eq. NF_INT) THEN
394                CALL hdlerr(NF_GET_VARA_INT 
395     $              (nc_mskid, il_varid, ist, icnt, 
396     $              msknew(iadrnew:iadrnew+isiznew-1)), 'inigrd')
397            ENDIF
398        ELSE
399            CALL locrint (clstrg,msknew(iadrnew),isiznew,iunit6,iflag)
400            IF (iflag .NE. 0) THEN
401                CALL prtout('Problem in reading unit =',iunit6,2)
402                CALL prcout
403     $          ('Problem with array linked to string =',clstrg,1)
404                infos = infos + iflag
405            ENDIF
406        ENDIF
407        IF (nlogprt .GE. 2) THEN
408            CALL prtout
409     $          ('Read target grid mask for field number =',jf,1)
410        ENDIF
411C
412!$omp end critical
413C* Read final mesh surface associated to field jf
414C
415!$omp critical
416        clstrg = clwork(1:icount)//csursuf
417        iunit7 = nulsu
418        IF (lncdfgrd) THEN
419            istatus = 
420     $          NF_INQ_VARID(nc_surid, clstrg, il_varid)
421            IF (istatus .ne. NF_NOERR .and. lsurf(jf)) THEN
422                CALL prtout('Cannot read surface information for fld '
423     $              , jf, 1)
424                CALL HALTE ('STOP in inigrd')
425            ELSE IF (istatus .eq. NF_NOERR) THEN
426                IF (n_reaty .eq. NF_FLOAT) THEN
427                    istatus = NF_GET_VARA_REAL(nc_surid, 
428     $              il_varid, ist, icnt, 
429     $                  surnew(iadrnew:iadrnew+isiznew-1))
430                    IF (istatus .ne. NF_NOERR .and. lsurf(jf)) THEN
431                        CALL prtout('Cannot read target grid surface
432     $ information for field ' , jf, 1)
433                        CALL HALTE ('STOP in inigrd')
434                    ENDIF
435                ELSE IF (n_reaty .eq. NF_DOUBLE) THEN
436                    istatus = NF_GET_VARA_DOUBLE(nc_surid, 
437     $               il_varid, ist, icnt, 
438     $                  surnew(iadrnew:iadrnew+isiznew-1))
439                    IF (istatus .ne. NF_NOERR .and. lsurf(jf)) THEN
440                        CALL prtout('Cannot read source grid surface
441     $ information for field ' , jf, 1)
442                        CALL HALTE ('STOP in inigrd')
443                    ENDIF
444                ENDIF
445            ENDIF
446        ELSE
447            CALL locread (clstrg, surnew(iadrnew), isiznew,
448     $          iunit7,iflag)
449            IF (iflag .NE. 0 .and. lsurf(jf)) THEN
450                CALL prtout('Problem in reading unit =',iunit7,2)
451                CALL prcout('Problem with array linked to string ',
452     $              clstrg,1)
453                infos = infos + iflag
454            ENDIF
455        ENDIF
456!$omp end critical
457        IF (nlogprt .GE. 2) THEN
458            CALL prtout
459     $          ('Read target grid surface for field number ',jf,1)
460        ENDIF
461 200  CONTINUE
462      IF (lncdfgrd) THEN
463          CALL hdlerr(NF_CLOSE(nc_grdid), 'inigrd')
464          CALL hdlerr(NF_CLOSE(nc_mskid), 'inigrd')
465          CALL hdlerr(NF_CLOSE(nc_surid), 'inigrd')
466      ENDIF
467C
468C* Finish up if problem in reading
469C
470      IF (infos .NE. 0) CALL HALTE('STOP in inigrd')
471C
472C
473C*    3. End of routine
474C        --------------
475C
476      IF (nlogprt .GE. 1) THEN
477          WRITE (UNIT = nulou,FMT = *) ' '
478          WRITE (UNIT = nulou,FMT = *) 
479     $    '          --------- End of routine inigrd ---------'
480          CALL FLUSH (nulou)
481      ENDIF
482      ENDIF
483      RETURN
484      END
485
486
Note: See TracBrowser for help on using the repository browser.