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 FIXAREA( ) 12C 13C----> 14C**** FIXAREA 15C 16C Purpose 17C ------- 18C 19C Fixup input/output field area definitions. 20C 21C 22C Interface 23C --------- 24C 25C IRET = FIXAREA( ) 26C 27C Input 28C ----- 29C 30C None. 31C 32C 33C Output 34C ------ 35C 36C None. 37C 38C 39C Method 40C ------ 41C 42C If default (0/0/0/0) selected for input, input area is set 43C - to global for spherical harmonics 44C - to global for lat/long grid 45C - to global for gaussian grid 46C 47C If default (0/0/0/0) selected for output, output area is set 48C - to same as input area for lat/long grid 49C - to same as input area for regular gaussian grid 50C - to global for reduced gaussian grid 51C 52C (Currently, subareas are not supported for reduced gaussian fields. 53C Should work OK; but need to correct setup of values in GRIB 54C product for number of points in each latitude row) 55C 56C Output area is adjusted to fit the given grid step. 57C 58C 59C Externals 60C --------- 61C 62C AREACHK - Match input/output field area definitions according to 63C grid specification. 64C GETENV - Get value of an environment variable 65C INTLOG - Logs output messages 66C 67C 68C Author 69C ------ 70C 71C J.D.Chambers ECMWF Jan 1995 72C 73C----< 74C 75 IMPLICIT NONE 76C 77C Function arguments 78C 79#include "parim.h" 80#include "nifld.common" 81#include "nofld.common" 82#include "grfixed.h" 83#include "intf.h" 84#include "current.h" 85C 86C Parameters 87C 88 INTEGER JPROUTINE 89 PARAMETER (JPROUTINE = 19200 ) 90C 91C Local variables 92C 93 INTEGER INORTH, ISOUTH, IEAST 94 INTEGER ITEMP, IRET 95 REAL EW, NS, NORTH, SOUTH, EAST, WEST, FACTOR3, FACTOR4 96 LOGICAL LDEFIN, LDEFOUT, LGLOBAL, LOVERDE 97 CHARACTER*20 OVERIDE 98 INTEGER IBLANK 99C 100C Externals 101C 102 INTEGER AREACHK 103C 104C ------------------------------------------------------------------ 105C* Section 1. Initialise 106C ------------------------------------------------------------------ 107C 108 100 CONTINUE 109 FIXAREA = 0 110 FACTOR3 = 0. 111 FACTOR4 = 0. 112 LDEFIN = .FALSE. 113 LGLOBAL = .FALSE. 114C 115C ------------------------------------------------------------------ 116C* Section 2. Fixup input area if default (0/0/0/0). 117C ------------------------------------------------------------------ 118C 119 200 CONTINUE 120C 121C Fixup input ocean area definition 122C 123 IF( NILOCAL.EQ.4 ) THEN 124 FACTOR3 = 0.0 125 FACTOR4 = 0.0 126C 127 IF( NIFORM.EQ.1 ) THEN 128 IF( ISEC1(60).EQ.1 ) FACTOR3 = 1.0 129 IF( ISEC1(60).EQ.2 ) THEN 130 IF( ISEC1(47).EQ.160 ) FACTOR3 = 1000.0 131 ENDIF 132 IF( ISEC1(60).EQ.3 ) FACTOR3 = 1000000.0 133 IF( ISEC1(60).EQ.4 ) FACTOR3 = 1000000.0 134C 135 IF( ISEC1(61).EQ.1 ) FACTOR4 = 1.0 136 IF( ISEC1(61).EQ.2 ) THEN 137 IF( ISEC1(47).EQ.160 ) FACTOR4 = 1000.0 138 ENDIF 139 IF( ISEC1(61).EQ.3 ) FACTOR4 = 1000000.0 140 IF( ISEC1(61).EQ.4 ) FACTOR4 = 1000000.0 141C 142 IF( (FACTOR3.EQ.0).OR.(FACTOR4.EQ.0) ) THEN 143 FIXAREA = 2 144 GOTO 900 145 ENDIF 146 NORTH = REAL(ISEC1(62))/FACTOR4 147 WEST = REAL(ISEC1(63))/FACTOR3 148 SOUTH = REAL(ISEC1(64))/FACTOR4 149 EAST = REAL(ISEC1(65))/FACTOR3 150 ELSE 151 IF( NIOCO3.EQ.1 ) FACTOR3 = 1.0 152 IF( NIOCO3.EQ.2 ) THEN 153 IF( NIVCDEF.EQ.160 ) FACTOR3 = 1000.0 154 ENDIF 155 IF( NIOCO3.EQ.3 ) FACTOR3 = 1000000.0 156 IF( NIOCO3.EQ.4 ) FACTOR3 = 1000000.0 157C 158 IF( NIOCO4.EQ.1 ) FACTOR4 = 1.0 159 IF( NIOCO4.EQ.2 ) THEN 160 IF( NIVCDEF.EQ.160 ) FACTOR4 = 1000.0 161 ENDIF 162 IF( NIOCO4.EQ.3 ) FACTOR4 = 1000000.0 163 IF( NIOCO4.EQ.4 ) FACTOR4 = 1000000.0 164C 165 IF( (FACTOR3.EQ.0).OR.(FACTOR4.EQ.0) ) THEN 166 FIXAREA = 2 167 GOTO 900 168 ENDIF 169 NORTH = REAL(NIOCO4F)/FACTOR4 170 WEST = REAL(NIOCO3F)/FACTOR3 171 SOUTH = REAL(NIOCO4L)/FACTOR4 172 EAST = REAL(NIOCO3L)/FACTOR3 173 ENDIF 174 NIAREA(1) = NINT( NORTH * PPMULT + 0.1) 175 NIAREA(2) = NINT( WEST * PPMULT + 0.1) 176 NIAREA(3) = NINT( SOUTH * PPMULT + 0.1) 177 NIAREA(4) = NINT( EAST * PPMULT + 0.1) 178 GOTO 300 179 ENDIF 180C 181C Sinisa bug fix for grid global fields 182C If input is lat/lon, check if it is global to within 183C a tolerance of 0.1 degrees 184 IF ( NIREPR .EQ. JPREGULAR ) THEN 185 IF ( NIAREA(4).GT.0 ) THEN 186 IEAST = NIAREA(4) 187 ELSE 188 IEAST = JP360 + NIAREA(4) 189 ENDIF 190 IF(NIAREA(1).EQ.JP90.AND.NIAREA(3).EQ.-JP90.AND. 191 X NIAREA(2).EQ.0.AND. 192 X IABS(JP360 - IEAST - NIGRID(1)).LT.1000) THEN 193 NIAREA(1) = 0 194 NIAREA(2) = 0 195 NIAREA(3) = 0 196 NIAREA(4) = 0 197csinisa 198cs ELSEIF(IABS(JP360 - IEAST - NIGRID(1)).LT.1000) THEN 199cs NOAREA(1) = NIAREA(1) 200cs NOAREA(2) = 0 201cs NOAREA(3) = NIAREA(3) 202cs NOAREA(4) = 36000000 203 ENDIF 204 ENDIF 205 206 LDEFIN = ( (NIAREA(1) .EQ. 0) .AND. (NIAREA(2) .EQ. 0) .AND. 207 X (NIAREA(3) .EQ. 0) .AND. (NIAREA(4) .EQ. 0) ) 208 209C 210 IF( LDEFIN ) THEN 211C 212 IF ( (NIREPR.EQ.JPSPHERE) .OR. (NIREPR.EQ.JPSPHROT) ) THEN 213C 214C Spectral input .. 215 EW = 0.0 216 NS = 0.0 217C 218 ELSE IF ( NIREPR .EQ. JPREGULAR ) THEN 219C 220C Regular lat/long grid .. 221 EW = FLOAT( NIGRID(1) ) / PPMULT 222 NS = FLOAT( NIGRID(2) ) / PPMULT 223C 224 ELSE 225C 226C Gaussian grid .. 227 EW = FLOAT( NIGAUSS ) 228 NS = 0.0 229 ENDIF 230C 231 NORTH = 0.0 232 WEST = 0.0 233 SOUTH = 0.0 234 EAST = 0.0 235 IRET = AREACHK( EW, NS, NORTH, WEST, SOUTH, EAST ) 236 IF( IRET.NE.0 ) THEN 237 FIXAREA = IRET 238 GOTO 900 239 ENDIF 240 NIAREA(1) = NINT( NORTH * PPMULT + 0.1) 241 NIAREA(2) = NINT( WEST * PPMULT + 0.1) 242 NIAREA(3) = NINT( SOUTH * PPMULT + 0.1) 243 NIAREA(4) = NINT( EAST * PPMULT + 0.1) 244C 245 ENDIF 246C 247C ------------------------------------------------------------------ 248C* Section 3. Fixup output area if default (0/0/0/0). 249C ------------------------------------------------------------------ 250C 251 300 CONTINUE 252C 253C See if environment variable has been specified to override 254C output area specification 255C 256 LOVERDE = .FALSE. 257 CALL GETENV('OVERRIDE_OUTPUT_AREA', OVERIDE) 258 IBLANK = INDEX(OVERIDE, ' ') 259 IF( IBLANK.GT.1 ) THEN 260 IF( OVERIDE(1:2).EQ.'ON' ) THEN 261 LOVERDE = .TRUE. 262 CALL INTLOG(JP_DEBUG, 263 X 'FIXAREA: OVERRIDE_OUTPUT_AREA is ON',JPQUIET) 264 ENDIF 265 ENDIF 266C 267 LDEFOUT = ( (NOAREA(1) .EQ. 0) .AND. (NOAREA(2) .EQ. 0) .AND. 268 X (NOAREA(3) .EQ. 0) .AND. (NOAREA(4) .EQ. 0) ) 269 X .OR. LOVERDE 270C 271C Fixup output ocean area definition 272C 273 IF( NILOCAL.EQ.4 ) THEN 274 IF( LDEFOUT ) THEN 275 IF( NIFORM.EQ.1 ) THEN 276 SOUTH = REAL(ISEC1(64))/FACTOR4 277 NORTH = REAL(ISEC1(62))/FACTOR4 278 WEST = REAL(ISEC1(63))/FACTOR3 279 EAST = REAL(ISEC1(65))/FACTOR3 280 281 EW = FLOAT( NOGRID(1) ) / PPMULT 282 NS = FLOAT( NOGRID(2) ) / PPMULT 283 284 IF( ISEC1(60).EQ.3 ) THEN 285 WEST = 0.0 286 EAST = 360.0 287 IRET = AREACHK( EW, NS, NORTH, WEST, SOUTH, EAST ) 288 IF ( IRET .NE. 0 ) THEN 289 FIXAREA = IRET 290 GOTO 900 291 ENDIF 292 ELSEIF( ISEC1(60).EQ.4 ) THEN 293 WEST = -90.0 294 EAST = 90.0 295 ENDIF 296 297 IF( ISEC1(61).EQ.3 ) THEN 298 SOUTH = 0.0 299 NORTH = 360.0 300 ELSEIF( ISEC1(61).EQ.4 ) THEN 301 SOUTH = -90.0 302 NORTH = 90.0 303 ENDIF 304 ELSE 305 NORTH = REAL(NIOCO4F)/FACTOR4 306 WEST = REAL(NIOCO3F)/FACTOR3 307 SOUTH = REAL(NIOCO4L)/FACTOR4 308 EAST = REAL(NIOCO3L)/FACTOR3 309 EW = FLOAT( NOGRID(1) ) / PPMULT 310 NS = FLOAT( NOGRID(2) ) / PPMULT 311 IF( NIOCO3.EQ.3 ) THEN 312 WEST = 0.0 313 EAST = 360.0 314 IRET = AREACHK( EW, NS, NORTH, WEST, SOUTH, EAST ) 315 IF ( IRET .NE. 0 ) THEN 316 FIXAREA = IRET 317 GOTO 900 318 ENDIF 319 ELSEIF( NIOCO3.EQ.4 ) THEN 320 WEST = -90.0 321 EAST = 90.0 322 ENDIF 323 324 IF( NIOCO4.EQ.3 ) THEN 325 SOUTH = 0.0 326 NORTH = 360.0 327 ELSEIF( NIOCO4.EQ.4 ) THEN 328 SOUTH = -90.0 329 NORTH = 90.0 330 ENDIF 331 ENDIF 332 333 NOAREA(1) = NINT( NORTH * PPMULT + 0.1) 334 NOAREA(2) = NINT( WEST * PPMULT + 0.1) 335 NOAREA(3) = NINT( SOUTH * PPMULT + 0.1) 336 NOAREA(4) = NINT( EAST * PPMULT + 0.1) 337 ELSE 338 NORTH = REAL(NOAREA(1)) / PPMULT 339 WEST = REAL(NOAREA(2)) / PPMULT 340 SOUTH = REAL(NOAREA(3)) / PPMULT 341 EAST = REAL(NOAREA(4)) / PPMULT 342 ENDIF 343 GOTO 900 344 ENDIF 345C 346 IF( LDEFOUT ) THEN 347C 348C If input is gaussian, check if it is global to within 349C a tolerance of 0.1 degrees. 350C 351 IF ( (NIREPR.EQ.JPGAUSSIAN) .OR. (NIREPR.EQ.JPQUASI) ) THEN 352 INORTH = NINT( RIGAUSS(1) * PPMULT ) 353 ISOUTH = -INORTH 354c EMOS-199: adjusted for reduced_gg/octahedral 355c IEAST = JP360 - (JP90/NIGAUSS) 356 IEAST = JP360 - (JP360/MILLEN(NIGAUSS)) 357 LGLOBAL = (NIAREA(2).EQ.0) 358 X .AND. (IABS(NIAREA(1)-INORTH).LE.1000) 359 X .AND. (IABS(NIAREA(3)-ISOUTH).LE.1000) 360 X .AND. (IABS(NIAREA(4)-IEAST ).LE.1000 .OR. 361 X NIREPR.EQ.JPGAUSSIAN .OR. .TRUE. ) 362C NOTE: a similar (but useless) test was done here before adding 363C .TRUE. be replicate the same (useless) behaviour 364 ENDIF 365C 366C Spectral output .. 367C 368C .. skip section fitting output area to grid 369C 370 IF ( (NOREPR.EQ.JPSPHERE).OR.(NOREPR.EQ.JPSPHROT) ) GOTO 900 371C 372C Regular lat/long grid .. 373C 374 IF ( (NOREPR.EQ.JPREGULAR).OR.(NOREPR.EQ.JPREGROT) ) THEN 375 IF( LDEFIN .OR. LGLOBAL.OR.(NIREPR.EQ.JPREDLL) ) THEN 376 EW = FLOAT( NOGRID(1) ) / PPMULT 377 NS = FLOAT( NOGRID(2) ) / PPMULT 378 NORTH = 90.0 379 WEST = 0.0 380 SOUTH = -NORTH 381 EAST = 360.0 382C 383C Use GRIB header values for reduced lat/long grids 384C (maybe 'Mediterranean' sub-area). 385C 386 IF( NIREPR.EQ.JPREDLL ) THEN 387 NORTH = REAL(NIAREA(1))/PPMULT 388 WEST = REAL(NIAREA(2))/PPMULT 389 SOUTH = REAL(NIAREA(3))/PPMULT 390 EAST = REAL(NIAREA(4))/PPMULT 391 ENDIF 392 IRET = AREACHK( EW, NS, NORTH, WEST, SOUTH, EAST ) 393 IF ( IRET .NE. 0 ) THEN 394 FIXAREA = IRET 395 GOTO 900 396 ENDIF 397 NOAREA(1) = NINT( NORTH * PPMULT + 0.1) 398 NOAREA(2) = NINT( WEST * PPMULT + 0.1) 399 NOAREA(3) = NINT( SOUTH * PPMULT + 0.1) 400 NOAREA(4) = NINT( EAST * PPMULT + 0.1) 401C 402 ELSE 403 NOAREA(1) = NIAREA(1) 404 NOAREA(2) = NIAREA(2) 405 NOAREA(3) = NIAREA(3) 406 NOAREA(4) = NIAREA(4) 407 ENDIF 408C 409 ENDIF 410C 411C Gaussian .. 412C 413 IF ( (NOREPR.EQ.JPGAUSSIAN) .OR. (NOREPR.EQ.JPQUASI) ) THEN 414C 415C Is output grid specification the same as the input? 416 IF( NIGAUSS.EQ.NOGAUSS) THEN 417 IF( (NOREPR.EQ.JPGAUSSIAN) .AND. LGLOBAL ) THEN 418 NOAREA(1) = JP90 419 NOAREA(2) = 0 420 NOAREA(3) = -JP90 421c EMOS-199: adjusted for reduced_gg/octahedral 422c NOAREA(4) = JP360 - (JP90/NOGAUSS) 423 NOAREA(4) = JP360 - (JP360/NOLPTS(NOGAUSS)) 424 ELSE 425 NOAREA(1) = NIAREA(1) 426 NOAREA(2) = NIAREA(2) 427 NOAREA(3) = NIAREA(3) 428 NOAREA(4) = NIAREA(4) 429 ENDIF 430C 431C Skip section fitting output area to grid 432 GOTO 900 433C 434 ELSE 435C 436C Different grid resolutions 437 EW = FLOAT( NOGAUSS ) 438 NS = 0.0 439 IF( (NOREPR.EQ.JPGAUSSIAN) .AND. LGLOBAL ) THEN 440 NORTH = 90.0 441 WEST = 0.0 442 SOUTH = -90.0 443c EMOS-199: adjusted for reduced_gg/octahedral 444c EAST = 360.0 - (360.0/(EW*4.0)) 445 EAST = 360.0 - (360.0/FLOAT(NOLPTS(NOGAUSS))) 446 ELSE 447 NORTH = FLOAT( NIAREA(1) ) / PPMULT 448 WEST = FLOAT( NIAREA(2) ) / PPMULT 449 SOUTH = FLOAT( NIAREA(3) ) / PPMULT 450 EAST = FLOAT( NIAREA(4) ) / PPMULT 451 ENDIF 452 IRET = AREACHK( EW, NS, NORTH, WEST, SOUTH, EAST ) 453 IF ( IRET .NE. 0 ) THEN 454 FIXAREA = IRET 455 GOTO 900 456 ENDIF 457 NOAREA(1) = NINT( NORTH * PPMULT + 0.1) 458 NOAREA(2) = NINT( WEST * PPMULT + 0.1) 459 NOAREA(3) = NINT( SOUTH * PPMULT + 0.1) 460 NOAREA(4) = NINT( EAST * PPMULT + 0.1) 461C 462C If regular lat/long input, check if output west longitude 463C reachs to full globe 464 ITEMP = NIAREA(4)+NIGRID(1) 465 IF( (NIREPR.EQ.JPREGULAR .AND. ITEMP.EQ.JP360) .OR. 466 X (NIREPR.EQ.JPSPHERE) .OR. 467 X (NIREPR.EQ.JPSPHROT) ) THEN 468c EMOS-199: adjusted for reduced_gg/octahedral 469c NOAREA(4) = JP360 - (JP90/NOGAUSS) 470 NOAREA(4) = JP360 - (JP360/NOLPTS(NOGAUSS)) 471 ENDIF 472 ENDIF 473 ENDIF 474C 475 ENDIF 476C 477C ------------------------------------------------------------------ 478C* Section 4. Now fixup output areas to correspond to the grid. 479C ------------------------------------------------------------------ 480C 481 400 CONTINUE 482C 483 IF ( (NOREPR.EQ.JPREGULAR).OR.(NOREPR.EQ.JPREGROT) ) THEN 484C 485C Regular lat/long grid .. 486 EW = FLOAT( NOGRID(1) ) / PPMULT 487 NS = FLOAT( NOGRID(2) ) / PPMULT 488C 489 ELSE 490C 491C Gaussian grid .. 492 EW = FLOAT( NOGAUSS ) 493 NS = 0.0 494 ENDIF 495C 496 NORTH = FLOAT( NOAREA(1) ) / PPMULT 497 WEST = FLOAT( NOAREA(2) ) / PPMULT 498 SOUTH = FLOAT( NOAREA(3) ) / PPMULT 499 EAST = FLOAT( NOAREA(4) ) / PPMULT 500 IRET = AREACHK( EW, NS, NORTH, WEST, SOUTH, EAST ) 501 IF ( IRET .NE. 0 ) THEN 502 FIXAREA = IRET 503 GOTO 900 504 ENDIF 505 NOAREA(1) = NINT( NORTH * PPMULT + 0.1) 506 NOAREA(2) = NINT( WEST * PPMULT + 0.1) 507 NOAREA(3) = NINT( SOUTH * PPMULT + 0.1) 508 NOAREA(4) = NINT( EAST * PPMULT + 0.1) 509 510C EMOS-199: [F|N] grids assume max(pl)=4*N, however O grids have 511C max(pl)=4*+16 and AREACHK doesn't know of the difference 512C The code below only checks if the adjusted EAST matches what 513C would happen on a [F|N] global grid, and corrects it for an O grid. 514 IF ( (NOREPR.EQ.JPQUASI) 515 X .AND. (HOGAUST.EQ.'O' 516 X .AND. NOAREA(2).EQ.0 517 X .AND. NOAREA(4).EQ.(JP360-JP90/NOGAUSS) )) 518 X NOAREA(4) = JP360 - (JP360/NOLPTS(NOGAUSS)) 519 520C 521C ------------------------------------------------------------------ 522C* Section 9. Closedown. 523C ------------------------------------------------------------------ 524C 525 900 CONTINUE 526C 527 RETURN 528 END 529