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 KRG2RGY(OLDFLD, NEWFLD, NUMGG) 12C 13C----> 14C**** KRG2RGY 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 = KRG2RGY(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 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 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 S.Curic ECMWF November 2006 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" 77#include "grfixed.h" 78C 79C Parameters 80C 81 INTEGER JPACK, JUNPACK, JPREALB 82 PARAMETER (JPACK = 4000000) 83 PARAMETER (JUNPACK = 6000000) 84#ifdef REAL_8 85 PARAMETER (JPREALB = 8) 86#else 87 PARAMETER (JPREALB = 8) 88#endif 89C 90C Local variables 91C 92 INTEGER IFSEC0(JPGRIB_ISEC0),IFSEC1(JPGRIB_ISEC1), 93 X IFSEC2(JPGRIB_ISEC2),IFSEC3(JPGRIB_ISEC3),IFSEC4(JPGRIB_ISEC4) 94 INTEGER NFSEC2(JPGRIB_ISEC2) 95 INTEGER ISIZE, LOOP 96 LOGICAL LFIRST 97 CHARACTER*1 HTYPE 98#ifndef _CRAYFTN 99#if (!defined sgi) 100#ifdef POINTER_64 101 INTEGER*8 IZFSEC4 102#else 103 INTEGER IZFSEC4 104#endif 105#endif 106#endif 107C 108 REAL PLAT(JPMAXNG) 109 REAL ZFSEC2(5000), ZFSEC3(2), ZFSEC4(1) 110C 111 POINTER( IZFSEC4, ZFSEC4 ) 112C 113 INTEGER IRET, ILENB, IPUNP, IWORD 114C 115C Externals 116C 117#ifdef POINTER_64 118 INTEGER*8 JMALLOC 119#else 120 INTEGER JMALLOC 121#endif 122 INTEGER KINTRG 123 EXTERNAL KINTRG, JMALLOC 124C 125 DATA LFIRST/.TRUE./ 126 SAVE LFIRST, IZFSEC4 127C 128C -----------------------------------------------------------------| 129C* Section 1. Initialise 130C -----------------------------------------------------------------| 131C 132 100 CONTINUE 133C 134 KRG2RGY = 0 135C 136 IF( (NUMGG*2).GT.JPMAXNG ) THEN 137 CALL INTLOG(JP_ERROR,'KRG2RGY: Gaussian N too big = ',NUMGG) 138 CALL INTLOG(JP_ERROR,'KRG2RGY: Maximum handled = ',(JPMAXNG/2)) 139 KRG2RGY = -1 140 GOTO 900 141 ENDIF 142C 143C First time, dynamically allocate memory for the field values 144C 145 IF( LFIRST ) THEN 146 LFIRST = .FALSE. 147 ISIZE = JUNPACK * JPREALB 148 CALL INTLOG(JP_DEBUG,'KRG2RGY: Allocate memory = ',ISIZE) 149C 150 IZFSEC4 = JMALLOC(ISIZE) 151#ifdef hpR64 152 IZFSEC4 = IZFSEC4/(1024*1024*1024*4) 153#endif 154 IF( IZFSEC4.EQ.0 ) THEN 155 CALL INTLOG(JP_ERROR,'KRG2RGY: JMALLOC fail IZFSEC4',JPQUIET) 156 KRG2RGY = -1 157 GOTO 900 158 ENDIF 159 ENDIF 160C 161C -----------------------------------------------------------------| 162C* Section 2. Unpack the input GRIB 163C -----------------------------------------------------------------| 164C 165 200 CONTINUE 166C 167C Old field 168C 169 IRET = 1 170 IPUNP = JUNPACK 171 ILENB = JPACK 172 IFSEC3(1) = 0 173 IFSEC3(2) = INT(RMISSGV) 174 ZFSEC3(1) = 0.0 175 ZFSEC3(2) = RMISSGV 176 CALL GRIBEX (IFSEC0,IFSEC1,IFSEC2,ZFSEC2,IFSEC3,ZFSEC3,IFSEC4, 177 X ZFSEC4,IPUNP,OLDFLD,ILENB,IWORD,'D',IRET) 178 IF( IRET.GT.0 ) THEN 179 WRITE(*,*) 'KRG2RGY: Old field gribex return code = ',IRET 180 KRG2RGY = -1 181 RETURN 182 ENDIF 183C 184C -----------------------------------------------------------------| 185C* Section 3. Create new field 186C -----------------------------------------------------------------| 187C 188 300 CONTINUE 189C 190 DO LOOP = 1, JPMAXNG 191 NFSEC2(LOOP) = IFSEC2(LOOP) 192 ENDDO 193 HTYPE = 'R' 194 IF( LOCTAHEDRAL ) HTYPE = 'O' 195 CALL JGETGG( NUMGG, HTYPE, PLAT, NFSEC2(23), IRET) 196 IF( IRET.NE.0 ) THEN 197 WRITE(*,*) 'KRG2RGY: Problem calling JGETGG, status = ',IRET 198 KRG2RGY = -1 199 RETURN 200 ENDIF 201 NFSEC2(3) = NUMGG * 2 202 NFSEC2(4) = NINT(PLAT(1)*1000.0) 203 NFSEC2(7) = NINT(PLAT(NUMGG*2)*1000.0) 204c EMOS-199: NFSEC2(8) adjusted for reduced_gg/octahedral 205c NFSEC2(8) = NINT((360.0-(360.0/REAL(NUMGG*4)))*1000.0) 206 NFSEC2(8) = NINT((360.0-(360.0/REAL(NFSEC2(23+NUMGG))))*1000.0) 207 NFSEC2(10) = NUMGG 208C 209 IRET = KINTRG(IFSEC2, NFSEC2, ZFSEC4, NEWFLD) 210 IF( IRET.LE.0 ) THEN 211 WRITE(*,*) 'KRG2RGY: New field creation failed' 212 KRG2RGY = -1 213 RETURN 214 ENDIF 215 IF( LDOUBLE ) THEN 216 CALL INTLOG(JP_DEBUG, 217 X 'KRG2RGY: Set NIAREA - double interpolation',JPQUIET) 218c NIAREA(1) = NFSEC2(4)*100 219c NIAREA(2) = NFSEC2(5)*100 220c NIAREA(3) = NFSEC2(7)*100 221c NIAREA(4) = NFSEC2(8)*100 222 NIAREA(1) = 0 223 NIAREA(2) = 0 224 NIAREA(3) = 0 225 NIAREA(4) = 0 226 NIGAUSS = NUMGG 227 RIGAUSS(1:NINS) = 0 228 MILLEN (1:NINS) = 0 229 NINS = NUMGG*2 230 ENDIF 231 232C 233 KRG2RGY = IRET 234C -----------------------------------------------------------------| 235C* Section 9. Return 236C -----------------------------------------------------------------| 237C 238 900 CONTINUE 239C 240 RETURN 241 END 242