001 /* 002 * $Id: StormAODT.java,v 1.1 2012/01/04 20:40:50 tommyj Exp $ 003 * 004 * This file is part of McIDAS-V 005 * 006 * Copyright 2007-2012 007 * Space Science and Engineering Center (SSEC) 008 * University of Wisconsin - Madison 009 * 1225 W. Dayton Street, Madison, WI 53706, USA 010 * https://www.ssec.wisc.edu/mcidas 011 * 012 * All Rights Reserved 013 * 014 * McIDAS-V is built on Unidata's IDV and SSEC's VisAD libraries, and 015 * some McIDAS-V source code is based on IDV and VisAD source code. 016 * 017 * McIDAS-V is free software; you can redistribute it and/or modify 018 * it under the terms of the GNU Lesser Public License as published by 019 * the Free Software Foundation; either version 3 of the License, or 020 * (at your option) any later version. 021 * 022 * McIDAS-V is distributed in the hope that it will be useful, 023 * but WITHOUT ANY WARRANTY; without even the implied warranty of 024 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 025 * GNU Lesser Public License for more details. 026 * 027 * You should have received a copy of the GNU Lesser Public License 028 * along with this program. If not, see http://www.gnu.org/licenses. 029 */ 030 031 package edu.wisc.ssec.mcidasv.data.cyclone; 032 033 import java.util.List; 034 035 import ucar.unidata.data.grid.GridUtil; 036 import ucar.unidata.util.IOUtil; 037 import ucar.unidata.util.StringUtil; 038 import visad.FlatField; 039 040 /** 041 * Created by IntelliJ IDEA. User: yuanho Date: Feb 20, 2009 Time: 3:09:14 PM To 042 * change this template use File | Settings | File Templates. 043 */ 044 045 public class StormAODT { 046 047 /** _more_ */ 048 StormAODTInfo.IRData odtcurrent_v72IR; 049 050 /** _more_ */ 051 StormAODTInfo.DataGrid areadata_v72; 052 053 /** _more_ */ 054 boolean lauto = false; 055 056 /** _more_ */ 057 int idomain_v72, ixdomain_v72, ifixtype_v72, rmwsizeman_v72; 058 059 /** _more_ */ 060 int oland_v72; 061 062 /** _more_ */ 063 boolean osearch_v72; 064 065 /** _more_ */ 066 int ostartstr_v72; 067 068 /** _more_ */ 069 float osstr_v72; 070 071 /** _more_ */ 072 static double c1 = 1.191066E-5; 073 074 /** _more_ */ 075 static double c2 = 1.438833; 076 077 public StormAODT() { 078 079 } 080 081 /** 082 * _more_ 083 * 084 * 085 * @param satgrid 086 * _more_ 087 * @param cenlat 088 * _more_ 089 * @param cenlon 090 * _more_ 091 * @param posm 092 * _more_ 093 * @param curdate 094 * _more_ 095 * @param cursat 096 * _more_ 097 * @param g_domain 098 * _more_ 099 * 100 * @return _more_ 101 */ 102 public StormAODTInfo.IRData aodtv72_drive(FlatField satgrid, float cenlat, 103 float cenlon, int posm, double curdate, int cursat, 104 String g_domain, int satId, int satChannel, boolean isTemperature) { 105 106 float ftmps, flats, flons, cenlon2; 107 108 int radius, irad, np, ii, jj, length; 109 int idomain = 0; 110 111 /* 112 * Set miscoptions flags in AODT 113 */ 114 115 int eyeSize = -99; 116 oland_v72 = 0; /* allow AODT operation over land */ 117 osearch_v72 = false; /* search for maximum curved band position */ 118 rmwsizeman_v72 = eyeSize; /* eye size parameter */ 119 odtcurrent_v72IR = new StormAODTInfo.IRData(); 120 odtcurrent_v72IR.domain = idomain_v72; 121 122 /* 123 * Set initial classification flag and value in AODT 124 */ 125 126 ostartstr_v72 = 0; /* user defined initial classification flag */ 127 osstr_v72 = 0.0f; /* starting initial classification value */ 128 129 /* 130 * Set image date/time info in AODT 131 */ 132 133 int iaodt = aodtv72_setIRimageinfo(curdate, cursat); 134 135 /* 136 * Get storm center lat/lon 137 */ 138 if (lauto == true) { 139 // aodtv72_runautomode( nauto, fauto, imagefile, &cenlat, &cenlon, 140 // &posm ); 141 } 142 143 /* 144 * Set center location in AODT positioning method (1=interpolation, 145 * 4=extrapolation, 0=error) 146 */ 147 // posm = 1; 148 aodtv72_setlocation(cenlat, cenlon, posm); 149 150 /* 151 * Set domain FLAG in AODT 152 */ 153 154 if (g_domain.equalsIgnoreCase("AUTO")) { 155 idomain = 0; 156 } 157 if (g_domain.equalsIgnoreCase("Atlantic")) { 158 idomain = 1; 159 } 160 if (g_domain.equalsIgnoreCase("Pacific")) { 161 idomain = 2; 162 } 163 if (g_domain.equalsIgnoreCase("Indian")) { 164 idomain = 2; 165 } 166 167 iaodt = aodtv72_setdomain(idomain); 168 169 /* 170 * Retrieve temperatures from image. This to be done in IDV 171 */ 172 173 GridUtil.Grid2D g2d = null; 174 float[][] temps = null; 175 float[][][] satimage = null; 176 float[][] lons = null; 177 float[][] lats = null; 178 int numx = 123; 179 int numy = 123; 180 181 try { 182 g2d = GridUtil.makeGrid2D(satgrid); 183 lons = g2d.getlons(); 184 lats = g2d.getlats(); 185 186 } catch (Exception re) { 187 } 188 189 /* now spatial subset numx by numy */ 190 GridUtil.Grid2D g2d1 = spatialSubset(g2d, cenlat, cenlon, numx, numy); 191 192 satimage = g2d1.getvalues(); 193 float[][] temp0 = satimage[0]; 194 int imsorc = satId, imtype = satChannel; 195 196 if (isTemperature) 197 temps = temp0; 198 else 199 temps = im_gvtota(numx, numy, temp0, imsorc, imtype); 200 201 /* 202 * Load the IR image information in AODT init areadata_v72 203 */ 204 205 aodtv72_loadIRimage(temps, g2d1.getlats(), g2d1.getlons(), numx, numy); 206 207 /* 208 * Set eye and cloud temperature values in AODT, return position for IR 209 * image data read 210 */ 211 212 StormAODTInfo.IRData tvIR = aodtv72_seteyecloudtemp( 213 StormAODTInfo.keyerM_v72, areadata_v72); 214 215 odtcurrent_v72IR.warmt = tvIR.warmt; 216 odtcurrent_v72IR.warmlatitude = tvIR.warmlatitude; 217 odtcurrent_v72IR.warmlongitude = tvIR.warmlongitude; 218 odtcurrent_v72IR.eyet = tvIR.eyet; 219 odtcurrent_v72IR.cwcloudt = tvIR.cwcloudt; 220 odtcurrent_v72IR.cwring = tvIR.cwring; 221 222 /* 223 * Determine scene type Set scene type 224 */ 225 226 float[] oscen = StormAODTSceneType.aodtv72_calcscene(odtcurrent_v72IR, 227 areadata_v72); 228 229 odtcurrent_v72IR.cloudt = oscen[0]; 230 odtcurrent_v72IR.cloudt2 = oscen[1]; 231 odtcurrent_v72IR.eyestdv = oscen[2]; 232 odtcurrent_v72IR.cloudsymave = oscen[3]; 233 odtcurrent_v72IR.eyefft = (int) oscen[4]; 234 odtcurrent_v72IR.cloudfft = (int) oscen[5]; 235 // { alst, Aaveext, Estdveye, Aavesym, eyecnt, rngcnt}; 236 float[] oscen1 = StormAODTSceneType.aodtv72_classify(odtcurrent_v72IR, 237 rmwsizeman_v72, areadata_v72, osstr_v72, osearch_v72); 238 239 odtcurrent_v72IR.eyescene = (int) oscen1[0]; 240 odtcurrent_v72IR.cloudscene = (int) oscen1[1]; 241 odtcurrent_v72IR.eyesceneold = -1; 242 odtcurrent_v72IR.cloudsceneold = -1; 243 odtcurrent_v72IR.eyecdosize = oscen1[2]; 244 odtcurrent_v72IR.ringcb = (int) oscen1[3]; 245 odtcurrent_v72IR.ringcbval = (int) oscen1[4]; 246 odtcurrent_v72IR.ringcbvalmax = (int) oscen1[5]; 247 odtcurrent_v72IR.ringcblatmax = oscen1[6]; 248 odtcurrent_v72IR.ringcblonmax = oscen1[7]; 249 odtcurrent_v72IR.rmw = oscen1[8]; 250 251 /* 252 * Determine intensity 253 */ 254 255 iaodt = aodtv72_calcintensity(idomain); 256 if (iaodt == 71) { 257 throw new IllegalStateException("center location is over land"); 258 } 259 260 /* 261 * Print out all diagnostic messages to screen 262 */ 263 return odtcurrent_v72IR; 264 } 265 266 public static class AODTResult { 267 } 268 269 /** 270 * _more_ 271 * 272 * @param g2d 273 * _more_ 274 * @param cenlat 275 * _more_ 276 * @param cenlon 277 * _more_ 278 * @param numx 279 * _more_ 280 * @param numy 281 * _more_ 282 * 283 * @return _more_ 284 */ 285 286 GridUtil.Grid2D spatialSubset(GridUtil.Grid2D g2d, float cenlat, 287 float cenlon, int numx, int numy) { 288 float[][] lats = g2d.getlats(); 289 float[][] lons = g2d.getlons(); 290 float[][][] values = g2d.getvalues(); 291 float[][] slats = new float[numx][numy]; 292 float[][] slons = new float[numx][numy]; 293 float[][][] svalues = new float[1][numx][numy]; 294 295 int ly = lats[0].length; 296 int ly0 = ly / 2; 297 int lx = lats.length; 298 int lx0 = lx / 2; 299 int ii = numx / 2, jj = numy / 2; 300 301 for (int j = 0; j < ly - 1; j++) { 302 if (Float.isNaN(lats[lx0][j])) 303 continue; 304 if ((lats[lx0][j] > cenlat) && (lats[lx0][j + 1] < cenlat)) { 305 jj = j; 306 } 307 } 308 for (int i = 0; i < lx - 1; i++) { 309 if (Float.isNaN(lons[i][ly0])) 310 continue; 311 if ((lons[i][ly0] < cenlon) && (lons[i + 1][ly0] > cenlon)) { 312 ii = i; 313 } 314 } 315 int startx = ii - (numx / 2 - 1); 316 int starty = jj - (numy / 2 - 1); 317 if (startx < 0) 318 startx = 0; 319 if (starty < 0) 320 starty = 0; 321 322 for (int i = 0; i < numx; i++) { 323 for (int j = 0; j < numy; j++) { 324 slats[i][j] = lats[i + startx][j + starty]; 325 slons[i][j] = lons[i + startx][j + starty]; 326 svalues[0][i][j] = values[0][i + startx][j + starty]; 327 } 328 } 329 330 return new GridUtil.Grid2D(slats, slons, svalues); 331 } 332 333 /** 334 * _more_ 335 * 336 * @param nx 337 * _more_ 338 * @param ny 339 * _more_ 340 * @param gv 341 * _more_ 342 * @param imsorc 343 * _more_ 344 * @param imtype 345 * _more_ 346 * 347 * @return _more_ 348 */ 349 float[][] im_gvtota(int nx, int ny, float[][] gv, int imsorc, int imtype) 350 351 /** 352 * im_gvtota 353 * 354 * This subroutine converts GVAR counts to actual temperatures based on the 355 * current image set in IM_SIMG. 356 * 357 * im_gvtota ( int *nvals, unsigned int *gv, float *ta, int *iret ) 358 * 359 * Input parameters: *nvals int Number of values to convert *gv int Array of 360 * GVAR count values 361 * 362 * Output parameters: *ta float Array of actual temperatures *iret int 363 * Return value = -1 - could not open table = -2 - could not find match 364 * 365 * 366 * Log: D.W.Plummer/NCEP 02/03 D.W.Plummer/NCEP 06/03 Add coeff G for 2nd 367 * order poly conv T. Piper/SAIC 07/06 Added tmpdbl to eliminate warning 368 */ 369 { 370 int ii, ip, chan, found, ier; 371 double Rad, Teff, tmpdbl; 372 float[][] ta = new float[nx][ny]; 373 int iret; 374 String fp = "/ucar/unidata/data/storm/ImgCoeffs.tbl"; 375 376 iret = 0; 377 378 for (ii = 0; ii < nx; ii++) { 379 for (int jj = 0; jj < ny; jj++) { 380 ta[ii][jj] = Float.NaN; 381 } 382 } 383 384 /* 385 * Read in coefficient table if necessary. 386 */ 387 String s = null; 388 try { 389 s = IOUtil.readContents(fp); 390 } catch (Exception re) { 391 } 392 393 int i = 0; 394 StormAODTInfo.ImgCoeffs[] ImageConvInfo = new StormAODTInfo.ImgCoeffs[50]; 395 for (String line : StringUtil.split(s, "\n", true, true)) { 396 if (line.startsWith("!")) { 397 continue; 398 } 399 List<String> stoks = StringUtil.split(line, " ", true, true); 400 401 ImageConvInfo[i] = new StormAODTInfo.ImgCoeffs(stoks); 402 ; 403 i++; 404 } 405 int nImgRecs = i; 406 found = 0; 407 ii = 0; 408 while ((ii < nImgRecs) && (found == 0)) { 409 410 tmpdbl = (double) (ImageConvInfo[ii].chan - 1) 411 * (ImageConvInfo[ii].chan - 1); 412 chan = G_NINT(tmpdbl); 413 414 if ((imsorc == ImageConvInfo[ii].sat_num) && (imtype == chan)) { 415 found = 1; 416 } else { 417 ii++; 418 } 419 420 } 421 422 if (found == 0) { 423 iret = -2; 424 return null; 425 } else { 426 427 ip = ii; 428 for (ii = 0; ii < nx; ii++) { 429 for (int jj = 0; jj < ny; jj++) { 430 431 /* 432 * Convert GVAR count (gv) to Scene Radiance 433 */ 434 Rad = ((double) gv[ii][jj] - ImageConvInfo[ip].scal_b) / 435 /* ------------------------------------- */ 436 ImageConvInfo[ip].scal_m; 437 438 Rad = Math.max(Rad, 0.0); 439 440 /* 441 * Convert Scene Radiance to Effective Temperature 442 */ 443 Teff = (c2 * ImageConvInfo[ip].conv_n) 444 / 445 /* 446 * -------------------------------------------------- 447 * ----- 448 */ 449 (Math.log(1.0 450 + (c1 * Math.pow(ImageConvInfo[ip].conv_n, 451 3.0)) / Rad)); 452 453 /* 454 * Convert Effective Temperature to Temperature 455 */ 456 ta[ii][jj] = (float) (ImageConvInfo[ip].conv_a 457 + ImageConvInfo[ip].conv_b * Teff + ImageConvInfo[ip].conv_g 458 * Teff * Teff); 459 } 460 461 } 462 463 } 464 465 return ta; 466 467 } 468 469 /** 470 * _more_ 471 * 472 * @param x 473 * _more_ 474 * 475 * @return _more_ 476 */ 477 public int G_NINT(double x) { 478 return (((x) < 0.0F) ? ((((x) - (float) ((int) (x))) <= -.5f) ? (int) ((x) - .5f) 479 : (int) (x)) 480 : ((((x) - (float) ((int) (x))) >= .5f) ? (int) ((x) + .5f) 481 : (int) (x))); 482 } 483 484 /** 485 * _more_ 486 * 487 * @param keyerM_v72 488 * _more_ 489 * @param areadata 490 * _more_ 491 * 492 * @return _more_ 493 */ 494 StormAODTInfo.IRData aodtv72_seteyecloudtemp(int keyerM_v72, 495 StormAODTInfo.DataGrid areadata) 496 /* 497 * Routine to search for, identify, and set the eye and cloud temperature 498 * values for the AODT library. Temperatures are set within AODT library. 499 * Inputs : none Outputs: none Return : -51 : eye, CWcloud, or warmest 500 * temperature <-100C or >+40C 0 : o.k. 501 */ 502 { 503 StormAODTInfo.IRData ird = StormAODTSceneType.aodtv72_gettemps( 504 keyerM_v72, areadata); 505 if (ird == null) { 506 throw new IllegalStateException( 507 "eye, CWcloud, or warmest temperature <-100C or >+40C"); 508 } 509 510 return ird; 511 512 // return iok; 513 } 514 515 /** 516 * _more_ 517 * 518 * @param temps 519 * _more_ 520 * @param lats 521 * _more_ 522 * @param lons 523 * _more_ 524 * @param numx 525 * _more_ 526 * @param numy 527 * _more_ 528 * 529 * @return _more_ 530 */ 531 public int aodtv72_loadIRimage(float[][] temps, float[][] lats, 532 float[][] lons, int numx, int numy) 533 /* 534 * Subroutine to load IR image data grid values (temperatures and positions) 535 * into data structure for AODT library Inputs : temperature, latitude, and 536 * longitude arrays centered on storm position location along with number of 537 * columns (x) and rows (y) in grid Outputs: none (areadata_v72 structure 538 * passed via global variable) Return : 0 : o.k. 539 */ 540 { 541 StormAODTInfo sinfo = new StormAODTInfo(); 542 /* allocate space for data */ 543 areadata_v72 = new StormAODTInfo.DataGrid(temps, lats, lons, numx, numy); 544 return 0; 545 } 546 547 /** 548 * _more_ 549 * 550 * @param indomain 551 * _more_ 552 * 553 * @return _more_ 554 */ 555 int aodtv72_setdomain(int indomain) 556 /* 557 * set current ocean domain variable within AODT library memory Inputs : 558 * domain flag value from input Outputs: none Return : -81 : error 559 * determining storm basin 560 */ 561 { 562 int domain; 563 float xlon; 564 565 /* obtain current storm center longitude */ 566 xlon = odtcurrent_v72IR.longitude; 567 if ((xlon < -180.0) || (xlon > 180.0)) { 568 return -81; 569 } 570 571 ixdomain_v72 = indomain; 572 /* determine oceanic domain */ 573 if (indomain == 0) { 574 /* automatically determined storm basin */ 575 if (xlon >= 0.0) { 576 domain = 0; /* atlantic and east pacific to 180W/dateline */ 577 } else { 578 domain = 1; /* west pacific and other regions */ 579 } 580 } else { 581 /* manually determined storm basin */ 582 domain = indomain - 1; 583 } 584 585 /* assign ocean domain flag value to AODT library variable */ 586 idomain_v72 = domain; 587 588 return 0; 589 } 590 591 /** 592 * _more_ 593 * 594 * @param ilat 595 * _more_ 596 * @param ilon 597 * _more_ 598 * @param ipos 599 * _more_ 600 * 601 * @return _more_ 602 */ 603 int aodtv72_setlocation(float ilat, float ilon, int ipos) 604 /* 605 * set current storm center location within from AODT library memory Inputs 606 * : AODT library current storm center latitude and longitude values and 607 * location positioning method : 1-forecast interpolation 2-laplacian 608 * technique 3-warm spot 4-extrapolation Outputs: none Return : -21 : 609 * invalid storm center position 21 : user selected storm center position 22 610 * : auto selected storm center position 611 */ 612 { 613 int iret; 614 615 /* assign current storm center latitude value to AODT library variable */ 616 odtcurrent_v72IR.latitude = ilat; 617 /* assign current storm center longitude value to AODT library variable */ 618 odtcurrent_v72IR.longitude = ilon; 619 /* assign current storm center positioning flag to AODT library variable */ 620 odtcurrent_v72IR.autopos = ipos; 621 if ((odtcurrent_v72IR.longitude < -180.) 622 || (odtcurrent_v72IR.longitude > 180.)) { 623 iret = -21; 624 } 625 if ((odtcurrent_v72IR.latitude < -90.) 626 || (odtcurrent_v72IR.latitude > 90.)) { 627 iret = -21; 628 } 629 630 iret = 21; /* user selected image location */ 631 if (ipos >= 1) { 632 iret = 22; 633 } 634 635 return iret; 636 } 637 638 /** 639 * _more_ 640 * 641 * @param datetime 642 * _more_ 643 * @param sat 644 * _more_ 645 * 646 * @return _more_ 647 */ 648 int aodtv72_setIRimageinfo(double datetime, int sat) 649 /* 650 * set IR image date/time within AODT library memory Inputs : AODT library 651 * IR image date/time/satellite information Outputs: none Return : 0 : o.k. 652 */ 653 { 654 /* assign IR image date to AODT library variable */ 655 odtcurrent_v72IR.date = datetime; 656 657 /* assign IR image satellite type to AODT library variable */ 658 odtcurrent_v72IR.sattype = sat; 659 660 return 0; 661 } 662 663 /** 664 * _more_ 665 * 666 * @param idomain 667 * _more_ 668 * 669 * @return _more_ 670 */ 671 public int aodtv72_calcintensity(int idomain) 672 /* 673 * Compute intensity values CI, Final T#, and Raw T#. Inputs : global 674 * structure odtcurrent_v72 containing current analysis Outputs : none 675 * Return : 71 : storm is over land 0 : o.k. 676 */ 677 { 678 int iok; 679 int iret; 680 int strength; 681 682 if ((odtcurrent_v72IR.land == 1)) { 683 iok = aodtv72_initcurrent(true, odtcurrent_v72IR); 684 iret = 71; 685 } else { 686 /* calculate current Raw T# value */ 687 odtcurrent_v72IR.Traw = aodtv72_Tnoraw(odtcurrent_v72IR, idomain); 688 odtcurrent_v72IR.TrawO = odtcurrent_v72IR.Traw; 689 /* check for spot analysis or full analysis using history file */ 690 /* if(hfile_v72==(char *)NULL) { */ 691 if (true) { 692 /* perform spot analysis (only Traw) */ 693 odtcurrent_v72IR.Tfinal = odtcurrent_v72IR.Traw; 694 odtcurrent_v72IR.Tfinal3 = odtcurrent_v72IR.Traw; 695 odtcurrent_v72IR.CI = odtcurrent_v72IR.Traw; 696 odtcurrent_v72IR.CIadjp = aodtv72_latbias(odtcurrent_v72IR.CI, 697 odtcurrent_v72IR.latitude, odtcurrent_v72IR.longitude, 698 odtcurrent_v72IR); 699 /* 700 * printf("%f %f %f %f\n",odtcurrent_v72IR.CI,odtcurrent_v72IR. 701 * latitude 702 * ,odtcurrent_v72->IR.longitude,odtcurrent_v72->IR.CIadjp); 703 */ 704 odtcurrent_v72IR.rule9 = 0; 705 /* odtcurrent_v72->IR.TIEraw=aodtv72_TIEmodel(); */ 706 /* odtcurrent_v72->IR.TIEavg=odtcurrent_v72->IR.TIEraw; */ 707 /* odtcurrent_v72->IR.TIEflag=aodtv72_tieflag(); */ 708 } else { 709 } 710 711 iret = 0; 712 } 713 714 return iret; 715 } 716 717 /** 718 * _more_ 719 * 720 * @param initval 721 * _more_ 722 * @param latitude 723 * _more_ 724 * @param longitude 725 * _more_ 726 * @param odtcurrent_v72IR 727 * _more_ 728 * 729 * @return _more_ 730 */ 731 float aodtv72_latbias(float initval, float latitude, float longitude, 732 StormAODTInfo.IRData odtcurrent_v72IR) 733 /* 734 * Apply Latitude Bias Adjustment to CI value Inputs : initval - initial CI 735 * value latitude - current latitude of storm Outputs : adjusted MSLP value 736 * as return value 737 */ 738 { 739 float initvalp; 740 741 float initvalpx; 742 float value; /* lat bias adjustement amount (0.00-1.00) */ 743 int sceneflag; /* 744 * contains lat bias adjustment flag 0=no adjustment 745 * 1=intermediate adjustment (6 hours) 2=full adjustment 746 */ 747 748 sceneflag = aodtv72_scenesearch(0); /* 749 * 0 means search for EIR based 750 * parameters... cdo, etc 751 */ 752 value = 1.0f; /* this value should be return from scenesearch() */ 753 /* printf("sceneflag=%d value=%f\n",sceneflag,value); */ 754 odtcurrent_v72IR.LBflag = sceneflag; 755 /* initvalp=aodtv72_getpwval(0,initval); TLO */ 756 initvalp = 0.0f; 757 if (sceneflag >= 2) { 758 /* EIR scene */ 759 if ((latitude >= 0.0) 760 && ((longitude >= -100.0) && (longitude <= -40.0))) { 761 /* do not make adjustment in N Indian Ocean */ 762 return initvalp; 763 } 764 /* apply bias adjustment to pressure */ 765 /* initvalp=-1.0*value*(-20.60822+(0.88463*A_ABS(latitude))); */ 766 initvalp = value * (7.325f - (0.302f * Math.abs(latitude))); 767 } 768 769 return initvalp; 770 } 771 772 /** 773 * _more_ 774 * 775 * @param type 776 * _more_ 777 * 778 * @return _more_ 779 */ 780 int aodtv72_scenesearch(int type) { 781 int curflag = 1, flag, eirflag; 782 float eirvalue, civalue, ciadjvalue, latitude; 783 double curtime, xtime, curtimem6, mergetimefirst, mergetimelast, firsttime = -9999.0; 784 785 float pvalue; 786 787 /* 788 * if(((odthistoryfirst_v72==0)&&(ostartstr_v72==TRUE))&&(hfile_v72!=(char 789 * *)NULL)) { 790 */ 791 792 if (true) { 793 flag = 2; 794 pvalue = 1.0f; 795 return flag; 796 } 797 798 return flag; 799 } 800 801 /** 802 * _more_ 803 * 804 * @param redo 805 * _more_ 806 * @param odtcurrent_v72IR 807 * _more_ 808 * 809 * @return _more_ 810 */ 811 int aodtv72_initcurrent(boolean redo, StormAODTInfo.IRData odtcurrent_v72IR) 812 /* 813 * initialize odtcurrent_v72 array or reset values for land interaction 814 * situations 815 */ 816 { 817 int b, bb; 818 // char comm[50]="\0"; 819 820 if (!redo) { 821 // odtcurrent_v72=(struct odtdata *)malloc(sizeof(struct odtdata)); 822 odtcurrent_v72IR.latitude = 999.99f; 823 odtcurrent_v72IR.longitude = 999.99f; 824 odtcurrent_v72IR.land = 0; 825 odtcurrent_v72IR.autopos = 0; 826 // strcpy(odtcurrent_v72IR.comment,comm); 827 // diagnostics_v72=(char *)calloc((size_t)50000,sizeof(char)); 828 // hfile_v72=(char *)calloc((size_t)200,sizeof(char)); 829 // fixfile_v72=(char *)calloc((size_t)200,sizeof(char)); 830 831 // b=sizeof(float); 832 // bb=sizeof(double); 833 // fcstlat_v72=(float *)calloc((size_t)5,b); 834 // fcstlon_v72=(float *)calloc((size_t)5,b); 835 // fcsttime_v72=(double *)calloc((size_t)5,bb); 836 } 837 838 odtcurrent_v72IR.Traw = 0.0f; 839 odtcurrent_v72IR.TrawO = 0.0f; 840 odtcurrent_v72IR.Tfinal = 0.0f; 841 odtcurrent_v72IR.Tfinal3 = 0.0f; 842 odtcurrent_v72IR.CI = 0.0f; 843 odtcurrent_v72IR.eyet = 99.99f; 844 odtcurrent_v72IR.warmt = 99.99f; 845 odtcurrent_v72IR.cloudt = 99.99f; 846 odtcurrent_v72IR.cloudt2 = 99.99f; 847 odtcurrent_v72IR.cwcloudt = 99.99f; 848 odtcurrent_v72IR.warmlatitude = 999.99f; 849 odtcurrent_v72IR.warmlongitude = 999.99f; 850 odtcurrent_v72IR.eyecdosize = 0.0f; 851 odtcurrent_v72IR.eyestdv = 0.0f; 852 odtcurrent_v72IR.cloudsymave = 0.0f; 853 odtcurrent_v72IR.eyescene = 0; 854 odtcurrent_v72IR.cloudscene = 0; 855 odtcurrent_v72IR.eyesceneold = -1; 856 odtcurrent_v72IR.cloudsceneold = -1; 857 odtcurrent_v72IR.rule9 = 0; 858 odtcurrent_v72IR.rule8 = 0; 859 odtcurrent_v72IR.LBflag = 0; 860 odtcurrent_v72IR.rapiddiss = 0; 861 odtcurrent_v72IR.eyefft = 0; 862 odtcurrent_v72IR.cloudfft = 0; 863 odtcurrent_v72IR.cwring = 0; 864 odtcurrent_v72IR.ringcb = 0; 865 odtcurrent_v72IR.ringcbval = 0; 866 odtcurrent_v72IR.ringcbvalmax = 0; 867 odtcurrent_v72IR.CIadjp = 0.0f; 868 odtcurrent_v72IR.rmw = -99.9f; 869 /* odtcurrent_v72->IR.TIEflag=0; */ 870 /* odtcurrent_v72->IR.TIEraw=0.0; */ 871 /* odtcurrent_v72->IR.TIEavg=0.0; */ 872 /* odtcurrent_v72->IR.sst=-99.9; */ 873 // if(!redo) odtcurrent_v72->nextrec=NULL; /* added by CDB */ 874 875 return 0; 876 } 877 878 /** 879 * Compute initial Raw T-Number value using original Dvorak rules 880 * 881 * @param odtcurrent 882 * @param idomain_v72 883 * @return return value is Raw T# 884 */ 885 886 float aodtv72_Tnoraw(StormAODTInfo.IRData odtcurrent, int idomain_v72) 887 /* 888 * Compute initial Raw T-Number value using original Dvorak rules Inputs : 889 * global structure odtcurrent_v72 containing current analysis Outputs : 890 * return value is Raw T# 891 * 892 * ODT SCENE/TEMPERATURE TABLE BD | WMG OW DG MG LG B W CMG CDG | TEMP |30.0 893 * 0.0 -30.0 -42.0 -54.0 -64.0 -70.0 -76.0 -80.0+| 894 * ---------------------------------------------------------------| Atl EYE 895 * | 3.5 4.0 4.5 4.5 5.0 5.5 6.0 6.5 7.0 | EMBC | 3.5 3.5 4.0 4.0 4.5 4.5 896 * 5.0 5.0 5.0 | CDO | 3.0 3.0 3.5 4.0 4.5 4.5 4.5 5.0 5.0 | 897 * ---------------------------------------------------------------| Pac EYE 898 * | 4.0 4.0 4.0 4.5 4.5 5.0 5.5 6.0 6.5 | EMBC | 3.5 3.5 4.0 4.0 4.5 4.5 899 * 5.0 5.0 5.0 | CDO | 3.0 3.5 3.5 4.0 4.5 4.5 4.5 4.5 5.0 | 900 * ---------------------------------------------------------------| Cat diff 901 * | 0 1 2 3 4 5 6 7 8 | add | 0.0 0.0 0.0 0.0 0.0-->0.5 0.5-->1.0 1.5 | 902 * (old) add |-0.5 -0.5 0.0 0.0-->0.5 0.5 0.5-->1.0 1.0 | (new) 903 * ---------------------------------------------------------------| 904 */ 905 { 906 907 double[][] eno = { 908 { 1.00, 2.00, 3.25, 4.00, 4.75, 5.50, 5.90, 6.50, 7.00, 7.50, 909 8.00 }, /* original plus adjusted > CDG+ */ 910 { 1.50, 2.25, 3.30, 3.85, 4.50, 5.00, 5.40, 5.75, 6.25, 6.50, 911 7.00 } }; /* adjusted based */ 912 double[][] cdo = { 913 { 2.00, 2.40, 3.25, 3.50, 3.75, 4.00, 4.10, 4.20, 4.30, 4.40, 914 4.70 }, 915 { 2.05, 2.40, 3.00, 3.20, 3.40, 3.55, 3.65, 3.75, 3.80, 3.90, 916 4.10 } }; 917 double[] curbnd = { 1.0, 1.5, 2.5, 3.0, 3.5, 4.0, 4.5 }; 918 double[] shrdst = { 0.0, 35.0, 50.0, 80.0, 110.0, 140.0 }; 919 double[] shrcat = { 3.5, 3.0, 2.5, 2.25, 2.0, 1.5 }; 920 921 double[][] diffchk = { 922 { 0.0, 0.5, 1.2, 1.7, 2.2, 2.7, 0.0, 0.0, 0.1, 0.5 }, /* 923 * shear 924 * scene 925 * types... 926 * original 927 * Rule 8 928 * rules 929 */ 930 { 0.0, 0.5, 1.7, 2.2, 2.7, 3.2, 0.0, 0.0, 0.1, 0.5 }, /* 931 * eye scene 932 * types... 933 * add 0.5 934 * to Rule 8 935 * rules 936 */ 937 { 0.0, 0.5, 0.7, 1.2, 1.7, 2.2, 0.0, 0.0, 0.1, 0.5 } }; /* 938 * other 939 * scene 940 * types 941 * ... 942 * subtract 943 * 0.5 944 * from 945 * Rule 946 * 8 947 * rules 948 */ 949 double eyeadjfacEYE[] = { 0.011, 0.015 }; /* 950 * modified wpac value to be 951 * closer to atlantic 952 */ 953 double symadjfacEYE[] = { -0.015, -0.015 }; 954 double dgraysizefacCLD[] = { 0.002, 0.001 }; 955 double symadjfacCLD[] = { -0.030, -0.015 }; 956 957 int diffchkcat; 958 int ixx, cloudcat, eyecat, diffcat, rp, xrp, rb; 959 float incval, lastci, lasttno, lastr9, lastraw; 960 float xpart, xparteye, xaddtno, eyeadj, spart, ddvor, dvorchart, ciadj; 961 float sdist, cloudtemp, eyetemp, fftcloud; 962 float t1val, t6val, t12val, t18val, t24val, delt1, delt6, delt12, delt18, delt24; 963 float t1valraw, t1valrawx, txvalmin, txvalmax; 964 double curtime, xtime, firsttime, firstlandtime; 965 double ttime1, ttime6, ttime12, ttime18, ttime24, t1valrawxtime; 966 StormAODTInfo.IRData odthistory, prevrec; 967 boolean oceancheck, adjustshear, firstland; 968 boolean t1found = false, t6found = false, t12found = false, t18found = false, t24found = false; 969 boolean first6hrs = false; 970 float symadj, dgraysizeadj, deltaT; 971 972 cloudtemp = odtcurrent.cloudt; 973 eyetemp = odtcurrent.eyet; 974 cloudcat = 0; 975 eyecat = 0; 976 lastci = 4.0f; 977 xpart = 0.0f; 978 979 for (ixx = 0; ixx < 10; ixx++) { 980 /* compute cloud category */ 981 if ((cloudtemp <= StormAODTInfo.ebd_v72[ixx]) 982 && (cloudtemp > StormAODTInfo.ebd_v72[ixx + 1])) { 983 cloudcat = ixx; 984 xpart = (float) (cloudtemp - StormAODTInfo.ebd_v72[cloudcat]) 985 / (float) (StormAODTInfo.ebd_v72[cloudcat + 1] - StormAODTInfo.ebd_v72[cloudcat]); 986 } 987 /* compute eye category for eye adjustment */ 988 if ((eyetemp <= StormAODTInfo.ebd_v72[ixx]) 989 && (eyetemp > StormAODTInfo.ebd_v72[ixx + 1])) { 990 eyecat = ixx; 991 } 992 /* eyetemp=Math.min(0.0,eyetemp); */ 993 } 994 if (odtcurrent.eyescene == 1) { 995 /* for pinhole eye, determine what storm should be seeing */ 996 /* 997 * eyetemp=pinhole(odtcurrent_v72->IR.latitude,odtcurrent_v72->IR.longitude 998 * ,eyetemp); 999 */ 1000 /* 1001 * eyetemp=(9.0-eyetemp)/2.0; / this matches DT used at NHC (jack 1002 * beven) 1003 */ 1004 eyetemp = (float) (eyetemp - 9.0) / 2.0f; /* 1005 * between +9C (beven) and 1006 * measured eye temp (turk) 1007 */ 1008 odtcurrent.eyet = eyetemp; 1009 } 1010 1011 /* category difference between eye and cloud region */ 1012 diffcat = Math.max(0, cloudcat - eyecat); 1013 1014 /* if scenetype is EYE */ 1015 rp = odtcurrent.ringcbval; 1016 rb = odtcurrent.ringcb; 1017 fftcloud = odtcurrent.cloudfft; 1018 1019 if (odtcurrent.cloudscene == 3) { 1020 /* CURVED BAND */ 1021 rp = Math.min(30, rp + 1); /* added 1 for testing */ 1022 xrp = rp / 5; 1023 incval = 0.1f; 1024 if (xrp == 1) { 1025 incval = 0.2f; 1026 } 1027 ddvor = (float) curbnd[xrp]; 1028 xaddtno = incval * (float) (rp - (xrp * 5)); 1029 /* 1030 * printf("rp=%d xrp=%d rb=%d ddvor=%f xaddtno=%f\n",rp,xrp,rb,ddvor 1031 * ,xaddtno); 1032 */ 1033 ddvor = ddvor + xaddtno; 1034 if (rb == 5) { 1035 ddvor = Math.min(4.0f, ddvor + 0.5f); 1036 } 1037 if (rb == 6) { 1038 ddvor = Math.min(4.5f, ddvor + 1.0f); 1039 } 1040 diffchkcat = 2; /* added for test - non-eye/shear cases */ 1041 } else if (odtcurrent.cloudscene == 4) { 1042 /* POSSIBLE SHEAR -- new definition from NHC */ 1043 ixx = 0; 1044 ddvor = 1.0f; 1045 sdist = odtcurrent.eyecdosize; /* shear distance */ 1046 while (ixx < 5) { 1047 if ((sdist >= shrdst[ixx]) && (sdist < shrdst[ixx + 1])) { 1048 spart = (float) ((sdist - shrdst[ixx]) / (shrdst[ixx + 1] - shrdst[ixx])); 1049 xaddtno = (float) ((spart * (shrcat[ixx + 1] - shrcat[ixx]))); 1050 ddvor = (float) (shrcat[ixx] + xaddtno); 1051 ixx = 5; 1052 } else { 1053 ixx++; 1054 } 1055 } 1056 diffchkcat = 0; /* added for test - shear cases */ 1057 } else { 1058 /* EYE or NO EYE */ 1059 if (odtcurrent.eyescene <= 2) { 1060 /* EYE */ 1061 xaddtno = (float) (xpart * (eno[idomain_v72][cloudcat + 1] - eno[idomain_v72][cloudcat])); 1062 /* 1063 * cloud category must be white (-70C) or below for full 1064 * adjustment; value will be merged in starting at black (-64C) 1065 * / if(cloudcat<5) { / gray shades / xparteye=0.00; } else 1066 * if(cloudcat==5) { / black / xparteye=xpart; } else { / white 1067 * and colder / xparteye=1.00; } 1068 */ 1069 eyeadj = (float) eyeadjfacEYE[idomain_v72] 1070 * (eyetemp - cloudtemp); 1071 /* symadj=-0.02*(odtcurrent_v72->IR.cloudsymave); */ 1072 symadj = (float) symadjfacEYE[idomain_v72] 1073 * (odtcurrent.cloudsymave); 1074 /* 1075 * printf("EYE : cloudsymave=%f symadj=%f\n",odtcurrent_v72->IR. 1076 * cloudsymave,symadj); 1077 */ 1078 ddvor = (float) eno[idomain_v72][cloudcat] + xaddtno + eyeadj 1079 + symadj; 1080 /* 1081 * printf("EYE : xaddtno=%f eyeadj=%f symadj=%f ddvor=%f\n",xaddtno 1082 * ,eyeadj,symadj,ddvor); 1083 */ 1084 ddvor = Math.min(ddvor, 9.0f); 1085 /* printf("ddvor=%f\n",ddvor); */ 1086 if (odtcurrent.eyescene == 2) { 1087 ddvor = Math.min(ddvor - 0.5f, 6.5f); /* LARGE EYE adjustment */ 1088 } 1089 /* 1090 * if(odtcurrent_v72->IR.eyescene==3) 1091 * ddvor=Math.min(ddvor-0.5,6.0); / LARGE RAGGED EYE adjustment 1092 */ 1093 diffchkcat = 1; /* added for test - eye cases */ 1094 /* printf("ddvor=%f\n",ddvor); */ 1095 } else { 1096 /* NO EYE */ 1097 /* CDO */ 1098 xaddtno = (float) (xpart * (cdo[idomain_v72][cloudcat + 1] - cdo[idomain_v72][cloudcat])); 1099 /* dgraysizeadj=0.002*odtcurrent_v72->IR.eyecdosize; */ 1100 dgraysizeadj = (float) dgraysizefacCLD[idomain_v72] 1101 * odtcurrent.eyecdosize; 1102 /* 1103 * printf("CDO : dgraysize=%f symadj=%f\n",odtcurrent_v72->IR.eyecdosize 1104 * ,dgraysizeadj); 1105 */ 1106 /* symadj=-0.03*(odtcurrent_v72->IR.cloudsymave); */ 1107 symadj = (float) symadjfacCLD[idomain_v72] 1108 * (odtcurrent.cloudsymave); 1109 /* 1110 * printf("CDO : cloudsymave=%f symadj=%f\n",odtcurrent_v72->IR. 1111 * cloudsymave,symadj); 1112 */ 1113 ddvor = (float) cdo[idomain_v72][cloudcat] + xaddtno 1114 + dgraysizeadj + symadj; 1115 ddvor = ddvor - 0.1f; /* bias adjustment */ 1116 /* 1117 * printf("CDO : xaddtno=%f dgraysizeadj=%f symadj=%f ddvor=%f\n" 1118 * ,xaddtno,dgraysizeadj,symadj,ddvor); 1119 */ 1120 ciadj = 0.0f; 1121 if (odtcurrent.cloudscene == 0) { /* CDO */ 1122 if (lastci >= 4.5) { 1123 ciadj = Math.max(0.0f, Math.min(1.0f, lastci - 4.5f)); 1124 } 1125 if (lastci <= 3.0) { 1126 ciadj = Math.min(0.0f, Math.max(-1.0f, lastci - 3.0f)); 1127 } 1128 /* printf("CDO : lastci=%f xaddtno=%f\n",lastci,ciadj); */ 1129 ddvor = ddvor + ciadj; 1130 } 1131 if (odtcurrent.cloudscene == 1) { /* EMBEDDED CENTER */ 1132 ciadj = Math.max(0.0f, Math.min(1.5f, lastci - 4.0f)); 1133 /* printf("EMBC : lastci=%f xaddtno=%f\n",lastci,ciadj); */ 1134 ddvor = ddvor + ciadj; /* changed from 0.5 */ 1135 } 1136 if (odtcurrent.cloudscene == 2) { /* IRREGULAR CDO (PT=3.5) */ 1137 ddvor = ddvor + 0.3f; /* additional IrrCDO bias adjustment */ 1138 ddvor = Math.min(3.5f, Math.max(2.5f, ddvor)); 1139 } 1140 diffchkcat = 2; /* added for test - non-eye/shear cases */ 1141 } 1142 } 1143 1144 dvorchart = ((float) (int) (ddvor * 10.0f)) / 10.0f; 1145 // odtcurrent_v72IR.TrawO=dvorchart; 1146 1147 return dvorchart; 1148 1149 } 1150 1151 }