source: CPL/oasis3/trunk/src/mod/oasis3/src/inigrd.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: 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.