source: CPL/oasis3/trunk/src/mod/oasis3/src/redglo.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: 11.1 KB
Line 
1      SUBROUTINE redglo 
2     $    (pzgg, pzgr, kredu, kinip, klon, klat, 
3     $     kmask, kland, ktronca, pmask, cdmsk)
4C*****
5C               *****************************
6C               * OASIS ROUTINE  -  LEVEL 3 *
7C               * -------------     ------- *
8C               *****************************
9C
10C**** *redglo* - Global gaussian grid transformation: reduced --> full grid
11C
12C     Purpose:
13C     -------
14C     Linear interpolation of field on reduced grid to full grid 
15C
16C     ***NOTE***
17C     ----------
18C     The flag cdmsk asks for extrapolation or not on the reduced
19C     gaussian grid. If cdmsk  = SEALAND sea values are extended
20C     to continental areas using the reduced grid sea-land mask.
21C     If cdmsk = LANDSEA, the opposite is performed. If cdmsk = 
22C     NOEXTRAP, THEN no extrapolation is performed.
23C
24C**   Interface:
25C     ---------
26C       *CALL*  *redglo (pzgg, pzgr, kredu, kinip, klon, klat,
27C                        kmask, kland, ktronca, pmask)
28C
29C     Input:
30C     -----
31C                pzgg    : field on reduced grid (real 2D)
32C                kredu   : number of points on reduced grid
33C                kinip   : number of longitudes per latitude value (integer 1D)
34C                klon    : number of longitudes of global grid
35C                klat    : number of latitudes of both grids
36C                pmask   : reduced grid mask value
37C                ktronca : truncature of gaussian grid
38C                cdmsk   : extrapolation flag
39C 
40C     Output:
41C     ------
42C                pzgg    : field on global grid  (real 2D)
43C                kmask   : mask on reduced grid (integer 1D)
44C
45C     Workspace:
46C     ---------
47C                pzgr    : work array to store field on reduced grid (real 1D)
48C                kland   : number of masked points per latitude circle on
49C                          reduced grid (integer 1D)
50C
51C     External:
52C     --------
53C     None
54C
55C     Reference:
56C     ---------
57C     See OASIS manual (1995)
58C
59C     History:
60C     -------
61C       Version   Programmer     Date      Description
62C       -------   ----------     ----      ----------- 
63C       1.0       L. Terray      94/01/01  created
64C       2.0       L. Terray      95/10/10  modified: new structure
65C       2.1       L. Terray      95/11/11  modified: choice of extrapolation
66C       2.3       S. Valcke      99/03/29  modified: error message for maskr
67C       2.3       S. Valcke      99/04/30  added: printing levels
68C       2.5       S. Valcke      01/03/23  added: nulrd opening
69C       2.5       S. Valcke      01/03/23  added: netCDF maskr files
70C
71C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
72C
73C* ---------------- Include files and USE of modules---------------------------
74C
75      USE mod_kinds_oasis
76      USE mod_parameter
77      USE mod_unitncdf
78      USE mod_unit
79      USE mod_gauss
80      USE mod_printing
81      USE mod_label
82      INCLUDE 'netcdf.inc'
83C
84C* ---------------------------- Argument declarations -------------------
85C
86      REAL (kind=ip_realwp_p) pzgg(klon,klat), pzgr(klon*klat)
87      INTEGER (kind=ip_intwp_p) kmask(kredu), kland(klat), kinip(klat)
88      CHARACTER*8 cdmsk
89C
90C* ---------------------------- Local declarations -------------------
91C
92      CHARACTER*8 clname
93      INTEGER (kind=ip_intwp_p) il_ncid, il_varid
94C
95C* ---------------------------- Poema verses ----------------------------
96C
97C %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
98C
99C*    1. Initialization
100C        --------------
101C
102      IF (nlogprt .GE. 2) THEN
103          WRITE (UNIT = nulou,FMT = *) ' '
104          WRITE (UNIT = nulou,FMT = *) ' '
105          WRITE (UNIT = nulou,FMT = *) 
106     $    '           ROUTINE redglo  -  Level 3'
107          WRITE (UNIT = nulou,FMT = *) 
108     $    '           **************     *******'
109          WRITE (UNIT = nulou,FMT = *) ' '
110          WRITE (UNIT = nulou,FMT = *) 
111     $    ' Go from reduced to global gaussian grid'
112          WRITE (UNIT = nulou,FMT = *) ' '
113          WRITE (UNIT = nulou,FMT = *) ' '
114      ENDIF
115      iflag = 0
116C
117C* Zeroes work array
118C
119      CALL szero(pzgr,klon*klat)
120C
121C* Zeroes some argument arrays
122C
123      CALL izero(kmask,kredu)
124      CALL izero(kland,klat)
125      zmask = pmask - 1.0
126C
127C* Open and READ in reduced gaussian grid maskr file
128C
129      IF (nlogprt .GE. 1) THEN
130          WRITE (UNIT = nulou,FMT = *) 
131     $    '   opening the reduced gaussian grid mask file '
132          WRITE (UNIT = nulou,FMT = *) ' '
133      ENDIF
134      WRITE(clname,'(''MSKRD'',I3.3)') ktronca
135      IF (lncdfgrd) THEN
136        istatus=NF_OPEN(crednam//'.nc', NF_NOWRITE, il_ncid)
137        IF (istatus .NE. NF_NOERR) THEN
138            WRITE (UNIT = nulou,FMT = *) 
139     $       ' ===>>>> : error in opening maskr.nc. We STOP!'
140            WRITE (UNIT = nulou,FMT = *) ' '
141            CALL HALTE('STOP in redglo')   
142        ENDIF
143C
144        istatus=NF_INQ_VARID(il_ncid, clname, il_varid)
145        istatus=NF_GET_VARA_INT (il_ncid, il_varid, 1, kredu, kmask)
146        IF (istatus .ne. NF_NOERR) THEN
147            WRITE (UNIT = nulou,FMT = *) 
148     $       ' ===>>>> : error in reading variable', clname
149            WRITE (UNIT = nulou,FMT = *) 
150     $         'in maskr.nc, the reduced grid mask file '
151            WRITE (UNIT = nulou,FMT = *) 
152     $            'Note that since version 2.3, the locator for the'
153            WRITE (UNIT = nulou,FMT = *) 
154     $            'reduced mask has to be MSKRDxxx WHERE xxx' 
155            WRITE (UNIT = nulou,FMT = *) 
156     $            'is half the number of latitude lines.'
157            WRITE (UNIT = nulou,FMT = *) ' '
158            CALL HALTE('STOP in redglo')   
159        ENDIF
160        CALL hdlerr(NF_CLOSE(il_ncid), 'redglo')
161      ELSE
162          iunit = nulrd
163          OPEN (UNIT = nulrd,FILE = crednam,STATUS='UNKNOWN',
164     $      FORM ='UNFORMATTED',ERR = 140,IOSTAT = iost)
165          IF (nlogprt .GE. 1) THEN
166              WRITE (UNIT = nulou,FMT = 1001) nulrd, crednam
167              WRITE (UNIT = nulou,FMT = *) ' '
168          ENDIF
169 140      CONTINUE
170C* Formats
171C
172 1001 FORMAT(10X,' open unit = ',I3,'    file ',A8,' ok')
173          IF (iost .ne. 0) THEN
174              WRITE (UNIT = nulou,FMT = *) 
175     $         ' ===>>>> : error in opening maskr, the red-grid file'
176              WRITE (UNIT = nulou,FMT = *) 
177     $         ' =======   =====                  ====='
178              WRITE (UNIT = nulou,FMT = *) 
179     $         ' logical unit ',nulrd,' error number = ',
180     $                     iost
181              WRITE (UNIT = nulou,FMT = *) 
182     $         ' We STOP. Verify the file ', crednam
183              WRITE (UNIT = nulou,FMT = *) ' '
184              CALL HALTE('STOP in redglo')   
185          ENDIF
186C
187          CALL locrint(clname, kmask, kredu, iunit, iflag)
188          IF (iflag .NE. 0) THEN
189              CALL prtout('WARNING: problem in reading unit =', 
190     $            iunit, 2)
191              CALL prcout('Error in reading field linked to string =',
192     $                clname, 1)
193              WRITE (UNIT = nulou,FMT = *) 
194     $            'Note that since version 2.3, the locator for the'
195              WRITE (UNIT = nulou,FMT = *) 
196     $            'reduced mask has to be MSKRDxxx WHERE xxx' 
197              WRITE (UNIT = nulou,FMT = *) 
198     $            'is half the number of latitude lines.'
199              CALL HALTE ('STOP in redglo')
200          ENDIF
201      ENDIF
202C* Assign reduced grid values to array pzgr
203C
204      DO 110 jj = 1, klat
205        DO 120 ji = 1, klon
206          pzgr(klon*(jj-1)+ji) = pzgg(ji,jj)
207 120    CONTINUE
208 110  CONTINUE
209C
210C* Test if extrapolation is required
211C
212      IF (cdmsk .EQ. 'SEALAND' .OR. cdmsk .EQ. 'LANDSEA') THEN
213C
214C
215C*    2. Check if there is any sea point and mask land points
216C        ----------------------------------------------------
217C
218C* First get number of masked points on each latitude circle
219C
220          indice = 0
221C* Extrapolation of land values from sea values, we set mask equal 1
222          IF (cdmsk .EQ. 'SEALAND') imask = 1
223C* Extrapolation of sea values from land values, we set mask equal 0
224          IF (cdmsk .EQ. 'LANDSEA') imask = 0
225          DO 210 jk = 1, klat
226            DO 220 ji = 1, kinip(jk)
227              IF (kmask(indice + ji) .EQ. imask) THEN
228                  kland(jk) = kland(jk) + 1
229              ENDIF
230 220        CONTINUE
231            indice = indice + kinip(jk)
232 210      CONTINUE
233C
234C* Then fill up continental or sea values unless all points are land or sea
235C
236          indice = 0
237          DO 230 jk = 1, klat
238            IF (kland(jk) .EQ. kinip(jk)) GO TO 250
239            DO 240 ji = 1, kinip(jk)
240              IF (kmask(indice + ji) .EQ. imask) 
241     $            pzgr(indice + ji) = pmask
242 240        CONTINUE
243 250        indice = indice + kinip(jk)
244 230      CONTINUE
245C
246C
247C*    3. Assign sea values to land points or land values to sea points
248C        -------------------------------------------------------------
249C
250          indice = 0
251          DO 310 jk = 1, klat
252            IF (kland(jk) .EQ. kinip(jk)) GOTO 315
253            DO 320 ji = 1, kinip(jk)
254              inow = indice + ji
255              IF (pzgr(inow) .GT. zmask) THEN
256                  iinf = indice + 1
257                  isup = indice + kinip(jk)
258                  ieast = inow + 1
259                  iwest = inow - 1
260                  icont = 0
261 325              CONTINUE
262                  IF (ieast .GT. isup) ieast = ieast - kinip(jk)
263                  IF (iwest .LT. iinf) iwest = iwest + kinip(jk)
264                  IF (pzgr(ieast) .LT. zmask .and.
265     $                pzgr(iwest) .LT. zmask) THEN
266                      pzgr(inow) = (pzgr(ieast) + pzgr(iwest))/2.
267                  ELSEIF (pzgr(ieast) .GT. zmask .and.
268     $                    pzgr(iwest) .LT. zmask) THEN
269                      pzgr(inow) = pzgr(iwest)
270                  ELSEIF (pzgr(ieast) .LT. zmask .and.
271     $                    pzgr(iwest) .GT. zmask) THEN
272                      pzgr(inow) = pzgr(ieast)
273                  ELSE
274                      icont = icont + 1
275                      ieast = ieast + 1
276                      iwest = iwest - 1
277                      if (icont .GT. kinip(jk)) 
278     $                    CALL HALTE ('STOP in redglo')
279                      GOTO 325
280                  ENDIF
281              ENDIF
282 320        CONTINUE
283 315        indice = indice + kinip(jk)
284 310      CONTINUE
285C* No extrapolation case
286        ELSE IF (cdmsk .EQ. 'NOEXTRAP') THEN
287          CONTINUE
288      ENDIF
289C
290C
291C*    2. Reduced to full linear interpolation 
292C        ------------------------------------
293C
294C* Interpolate from reduced to global
295C 
296      indice = 0
297      DO 270 jk = 1, klat
298        DO 280 ji = 1, klon
299          zxi = 1 + ((ji - 1) * kinip(jk)) / FLOAT(klon)
300          im = INT(zxi)
301          zdx = zxi - im
302          im = 1 + MOD(im + kinip(jk) - 1,kinip(jk))
303          ip = 1 + MOD(im,kinip(jk))
304          pzgg(ji,jk) = pzgr(indice + im) * (1.-zdx) 
305     $                + pzgr(indice + ip) * zdx
306 280    CONTINUE
307        indice = indice + kinip(jk)
308 270  CONTINUE
309C
310C
311C*    3. End of routine
312C        --------------
313C
314      IF (nlogprt .GE. 2) THEN
315          WRITE (UNIT = nulou,FMT = *) ' '
316          WRITE (UNIT = nulou,FMT = *) 
317     $    '          --------- End of routine redglo ---------'
318          CALL FLUSH (nulou)
319      ENDIF
320      RETURN
321      END
Note: See TracBrowser for help on using the repository browser.