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