1C Copyright 1981-2016 ECMWF.
2C
3C This software is licensed under the terms of the Apache Licence
4C Version 2.0 which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
5C
6C In applying this licence, ECMWF does not waive the privileges and immunities
7C granted to it by virtue of its status as an intergovernmental organisation
8C nor does it submit to any jurisdiction.
9C
10
11      INTEGER FUNCTION KRG2RG(OLDFLD, NEWFLD, NUMGG)
12C
13C---->
14C**** KRG2RG
15C
16C     Purpose
17C     -------
18C
19C     Interpolate a field based on one reduced gaussian grid to a
20C     field based on a different reduced gaussian grid.
21C
22C
23C     Interface
24C     ---------
25C
26C     IRET = KRG2RG(OLDFLD, NEWFLD, NUMGG)
27C
28C     Input
29C     -----
30C
31C     OLDFLD - Old field in GRIB format
32C     NUMGG  - Number of the new gaussian field
33C
34C
35C     Output
36C     ------
37C
38C     NEWFLD - New field in GRIB format
39C
40C     Function returns:
41C       - the size in bytes of the new GRIB product if all is well
42C       - -1, otherwise.
43C
44C
45C     Method
46C     ------
47C
48C     Unpack input GRIBS.
49C     Create new field and pack it into GRIB format.
50C
51C
52C     Externals
53C     ---------
54C
55C     GRIBEX  - Decode and encode GRIB products.
56C     KINTRG  - Interpolate reduced gaussian field to reduced gaussian
57C     JMALLOC - Dynamically allocate memory
58C     INTLOG  - Log messages
59C
60C
61C     Author
62C     ------
63C
64C     J.D.Chambers     ECMWF     November 2000
65C
66C----<
67C
68      IMPLICIT NONE
69C
70C     Function arguments
71C
72      INTEGER NUMGG, OLDFLD(*), NEWFLD(*)
73C
74#include "parim.h"
75#include "nifld.common"
76#include "nofld.common"
77C
78C     Parameters
79C
80      INTEGER JPACK, JUNPACK, JPINTB, JPREALB
81      PARAMETER (JPACK   = 4000000)
82      PARAMETER (JUNPACK = 6000000)
83#ifdef INTEGER_8
84      PARAMETER (JPINTB = 8)
85#else
86      PARAMETER (JPINTB = 4)
87#endif
88#ifdef REAL_8
89      PARAMETER (JPREALB = 8)
90#else
91      PARAMETER (JPREALB = 8)
92#endif
93C
94C     Local variables
95C
96      INTEGER IFSEC0(JPGRIB_ISEC0),IFSEC1(JPGRIB_ISEC1),
97     X IFSEC2(JPGRIB_ISEC2),IFSEC3(JPGRIB_ISEC3),IFSEC4(JPGRIB_ISEC4)
98      INTEGER NFSEC2(JPGRIB_ISEC2),NFSEC3(JPGRIB_ISEC3)
99      INTEGER ISIZE, LOOP
100      LOGICAL LFIRST
101      CHARACTER*1 HTYPE
102#ifndef _CRAYFTN
103#if (!defined sgi)
104#ifdef POINTER_64
105      INTEGER*8 IZFSEC4, IPFSEC4
106#else
107      INTEGER   IZFSEC4, IPFSEC4
108#endif
109#endif
110#endif
111C
112      REAL PLAT(JPMAXNG)
113      REAL ZFSEC2(5000), ZFSEC3(2), ZFSEC4(1)
114      REAL PFSEC3(2), PFSEC4(1)
115C
116      POINTER( IZFSEC4, ZFSEC4 )
117      POINTER( IPFSEC4, PFSEC4 )
118C
119      INTEGER IRET, ILENB, IPUNP, IWORD
120C
121C     Externals
122C
123#ifdef POINTER_64
124      INTEGER*8 JMALLOC
125#else
126      INTEGER JMALLOC
127#endif
128      INTEGER KINTRG
129      EXTERNAL KINTRG, JMALLOC
130C
131      DATA LFIRST/.TRUE./
132      SAVE LFIRST, IZFSEC4, IPFSEC4
133C
134C     -----------------------------------------------------------------|
135C*    Section 1.   Initialise
136C     -----------------------------------------------------------------|
137C
138  100 CONTINUE
139C
140      KRG2RG = 0
141C
142      IF( (NUMGG*2).GT.JPMAXNG ) THEN
143        CALL INTLOG(JP_ERROR,'KRG2RG: Gaussian N too big = ',NUMGG)
144        CALL INTLOG(JP_ERROR,'KRG2RG: Maximum handled = ',(JPMAXNG/2))
145        KRG2RG = -1
146        GOTO 900
147      ENDIF
148C
149C     First time, dynamically allocate memory for the field values
150C
151      IF( LFIRST ) THEN
152        LFIRST = .FALSE.
153        ISIZE = JUNPACK * JPREALB
154        CALL INTLOG(JP_DEBUG,'KRG2RG: Allocate memory = ',ISIZE)
155C
156        IZFSEC4 = JMALLOC(ISIZE)
157#ifdef hpR64
158        IZFSEC4 = IZFSEC4/(1024*1024*1024*4)
159#endif
160        IF( IZFSEC4.EQ.0 ) THEN
161          CALL INTLOG(JP_ERROR,'KRG2RG: JMALLOC fail IZFSEC4',JPQUIET)
162          KRG2RG = -1
163          GOTO 900
164        ENDIF
165C
166        IPFSEC4 = JMALLOC(ISIZE)
167#ifdef hpR64
168        IPFSEC4 = IPFSEC4/(1024*1024*1024*4)
169#endif
170        IF( IPFSEC4.EQ.0 ) THEN
171          CALL INTLOG(JP_ERROR,'KRG2RG: JMALLOC fail IPFSEC4',JPQUIET)
172          KRG2RG = -1
173          GOTO 900
174        ENDIF
175      ENDIF
176C
177C     -----------------------------------------------------------------|
178C*    Section 2.   Unpack the input GRIB
179C     -----------------------------------------------------------------|
180C
181  200 CONTINUE
182C
183C     Old field
184C
185      IRET = 1
186      IPUNP = JUNPACK
187      ILENB = JPACK
188      IFSEC3(1) = 0
189      IFSEC3(2) = INT(RMISSGV)
190      ZFSEC3(1) = 0.0
191      ZFSEC3(2) = RMISSGV
192      CALL GRIBEX (IFSEC0,IFSEC1,IFSEC2,ZFSEC2,IFSEC3,ZFSEC3,IFSEC4,
193     X             ZFSEC4,IPUNP,OLDFLD,ILENB,IWORD,'D',IRET)
194      IF( IRET.GT.0 ) THEN
195        WRITE(*,*) 'KRG2RG: Old field gribex return code = ',IRET
196        KRG2RG = -1
197        RETURN
198      ENDIF
199C
200C     -----------------------------------------------------------------|
201C*    Section 3.   Create new field and put in GRIB format
202C     -----------------------------------------------------------------|
203C
204  300 CONTINUE
205C
206      DO LOOP = 1, JPMAXNG
207        NFSEC2(LOOP) = IFSEC2(LOOP)
208      ENDDO
209      HTYPE = 'R'
210      IF( LOCTAHEDRAL ) HTYPE = 'O'
211      CALL JGETGG( NUMGG, HTYPE, PLAT, NFSEC2(23), IRET)
212      IF( IRET.NE.0 ) THEN
213        WRITE(*,*) 'KRG2RG: Problem calling JGETGG, status = ',IRET
214        KRG2RG = -1
215        RETURN
216      ENDIF
217      NFSEC2(3)  = NUMGG * 2
218      NFSEC2(4)  = NINT(PLAT(1)*1000.0)
219      NFSEC2(7)  = NINT(PLAT(NUMGG*2)*1000.0)
220c     EMOS-199: adjusted for reduced_gg/octahedral
221c     NFSEC2(8)  = NINT((360.0-(360.0/REAL(NUMGG*4)))*1000.0)
222      NFSEC2(8)  = NINT((360.0-(360.0/REAL(NFSEC2(23+NUMGG))))*1000.0)
223      NFSEC2(10) = NUMGG
224C
225      IRET = KINTRG(IFSEC2, NFSEC2, ZFSEC4, PFSEC4)
226      IF( IRET.LE.0 ) THEN
227        WRITE(*,*) 'KRG2RG: New field creation failed'
228        KRG2RG = -1
229        RETURN
230      ENDIF
231      IFSEC4(1) = IRET
232C
233      IRET = 1
234      IPUNP = JUNPACK
235      ILENB = JPACK
236C
237C     If output 'accuracy' specified, use it for the new field
238C
239      IF( LNOACC ) THEN
240        IFSEC4(2) = NOACC
241      ENDIF
242C
243C     There may be 'missing' values in the new field if the old field
244C     has 'missing' values.
245C
246      IF( IFSEC1(5).EQ.192 ) THEN
247        IFSEC1(5) = 192
248        NFSEC3(1) = 0
249        NFSEC3(2) = INT(RMISSGV)
250        PFSEC3(1) = 0.0
251        PFSEC3(2) = RMISSGV
252      ENDIF
253C
254      CALL GRIBEX (IFSEC0,IFSEC1,NFSEC2,ZFSEC2,NFSEC3,PFSEC3,IFSEC4,
255     X             PFSEC4,IPUNP,NEWFLD,ILENB,IWORD,'C',IRET)
256      IF( IRET.NE.0 ) THEN
257        WRITE(*,*) 'KRG2RG: New field gribex return code = ',IRET
258        KRG2RG = -1
259        RETURN
260      ENDIF
261C
262      KRG2RG = IWORD * JPINTB
263C
264C     -----------------------------------------------------------------|
265C*    Section 9.   Return
266C     -----------------------------------------------------------------|
267C
268  900 CONTINUE
269C
270      RETURN
271      END
272