c23456789012345678901234567890123456789012345678901234567890123456789012 c Subroutine rapidge(perri,ishgc,id200,vsstc, + ipw19,sbtri,rirpcri,irhcn,icflx,vmaxri,sname, + atcfid,ismon,isday,isyr,istime,lush,probridisc) c c Subroutine rapidge.f c Last updated: February 9, 2024 c Author: John Kaplan c Purpose: This subroutine is used to compute the 2024 multi-lead time E. Pacific c SHIPS-RII that was initially adopted for operational implementation c by NHC in 12/2014. The multi-lead time SHIPS-RII provides estimates c of the probability of RI for the following lead times and RI thresholds: c 12-h/20-kt, 24-h/25-kt, 24-h/30-kt,24-h/35-kt, 24-h/40-kt, 36-h/45-kt, c and 48-h/55-kt, 72-h/65-kt based upon linear discriminant analysis. The c above SHIPS-RII probabilities are averaged together with the logistic c and Bayesian RI estimates for the same lead time and RI threshold in routine c rapidprint_72h.f to obtain a consensus RI probability. c c 2024 SHIPS multi-lead time RII notes: c c 1) Janauary 2024: c a) The RII was re-derived using the updated 1998-2023 E. Pacific SHIPS c data base (operational data was used for cases from 2000-2023 and reanalysis c data was used for cases from 1995-2000). c b) Only cases with potential values greater or equal to the RI threshold were used. c c) SSTs used to compute MPI were as follows: c a) Reynolds weekly SSTs from 1998-2014 c b) NCODA daily-averaged 1 day lagged SSTs (NSST) from 2015-2023 c d) SHGC is used in place of SHDC in this versiom of the RII c e) TPW19 averaged over each given lead time is used in this version of the SHIPS-RII c whereas in previous versions the TPW at T=0 h was used for all forecast lead times. c c 2) Nudging of the scaled RI predictors was used (i.e. the range of the c scaling predictors was expanded by 10%) when the RII was derived. c c 3) As noted above the SHIPS multi-lead time RII provides estimates c of RI probabilty for the following forecast lead-times/RI c thresholds: 12-h/20-kt, 24-h/25-kt, 24-h/30-kt, 24-h/35-kt, c 24-h/40-kt, 36-h/45-kt, 48-h/55-kt, 72-h/65-kt. c c 4) Probabilities are set to zero whenever any of the scaled c values is < 0. c c 5) The climatological heat content is used to run the SHIPS-RII c when the observed heat content is missing. c c 6) Notes: 1) as of July 29, 2016 the model was modified to produce a 72-h c RI probability of a 65-kt intensity change. c 2) GFS model-derived TPW is now used in place of satellite-derived c estimates in the real-time version of the RII models. c 7) In 2021 the SHIPS-RII was updated so that the developmental c probabilities from 5 bins rather than 4 were used to make c the operational RII forecasts c c Subroutine inputs: perri,ishgc,id200,Vsstc,ipw19, c sbtri,rirpcri,irhcn,icflx,vmaxri,sname,ismon, c isday,isyr,istime, and lush c c Subroutine output: probridisc c c Notes: This subroutine assumes that missing input values c are > 900. (e.g., perri=999. if missing) except for c the cflx or tpw parameters where values >9000. are assumed c to be missing. Also, missing output values including probabilities c are coded as 999. c c perri : 12 h intensity change observed for the preceding 12 h (kt) c c ishgc: Array containing the multi-layer vertical shear (kt) computed c from 0-500 km radius after the vortex has been removed (scaled *10) c c id200: Array containing 200 mb divergence (10**-7s-1) from 0-1000 km c c vsstc : Array containing MPI (kt) c Note that the E. Pacific MPI that is used when computing the potential c intensity uses the Whitney and Hobgood (2007) MPI algorithm and c does not use the Cione inner-core cooled SST. c However, the MPI used also include a correction for storm speed c based upon the algorithm developed by Schwerdt et al. (1979). c c ipw19: Array containing % area with 0-500 km radius upshear with precipitable water c <45 mm (scaled *10). Note this version of the RII uses the time-averaged Ipw19 c sbtri : std. dev. of the 50-200 km GOES channel 4 brightness c temperatures at t=0 h c rirpcri: second principle component computed from GOES channel 4 c brightness temperatures at time t= 0 c irhcn : Array containing Reynolds heat content (kj/cm2) (time-averaged) c Icflx : Array containing Boundary-layer latent heat flux (watts/m2) difference c predictor (q10layr-q10)*Vmx0 where q10 is obtained by bring air down c from 1000 mb and then allowing air to cool as it reaches c 95% RH in the core. q10layr is obtained using the same c methodology except using RH= 1/3* (RH_sfc + Rh_850-700 mb c + RH_700-500 mb) (time-averaged). c see JHT 2010 year 1 final report for more details c vmaxri: Maximum sustained wind at time t=0 (kt) c sname : Storm name (A10) c atcfid: ATCF storm ID (A8) c ismon : Month of forecast (I2) c Isday : Day of month of forecast (I2) c Isyr : Year of forecast (I2) c istime: Time of day of forecast (UTC) (I2) c lush : Output unit number of SHIPS log file c ratscd: The ratio of the probabilty of RI (prbscl) and the dependent sample mean c probability of RI that was computed using over-water cases from 1989-2020 (clrisc) c Rvar : Array containing the predictor magnitudes of each of the RII predictors. c sclvar: Array containing the scaled magnitudes of each RII predictor. c sclvrd: Array containing the scaled discriminant magnitudes of each RII predictor. c rmnval: Array containg the min. value minus 10% at which RI occured for each RII predictor c rmxval: Array containg the max. value plus 10% at which RI occured for each RII predictor c avgval: The mean value at which RI occurred for each RII predictor c sclav : The binned magnitudes of the discrminant RII values c prbri : The discrminant RII probabilites corresponding to each of the binned c RII values in the array sclav. c scaled: The case specific value of the discrminant RII . c prbrid: The case specific magnitudes of the discriminant RII probabilites c prscld: Array containing the RII probabilities computed for each case c for each RI threshold c Nindx = number of RI index thresholds c Mwrite= The index number used to print out the scaled and weighted values on c the log sheet (mwrite=3 is for the 24-h/ 30 kt RI threshold) c rvarsrt: array containg sorted RII variables c labsrt : Character array containing RII predictor labels c labdatsrt: Character array containing truncated predictor labels c rmnvalsrt: Array containg sorted min. value minus 10% at which RI occured for each RII predictor. c rmxvalsrt: Array containg sorted max. value plus 10% at which RI occured for each RII predictor. c dwgtsrt : Array containing sorted discriminant weights of each of the RII predictors c Sclvrdsrt: Array containing the total weighted discriminant magnitudes for each c RII predictor for each forecast period. c totwgt: Array containing weights of each of the predictors for each RI threshold c totalwgt: Array containing sum of the weights of each predictor for the 8 RI thresholds c prbridisc: Array containing the RII probabilities for the 5 lead times (12-h,24-h,36-h,48-h,72-h) c with row 1 corresponding to the 12-h lead time and row 5 corresponding to c the 72-h lead time. c Setzero: Flag for determining if RII probability is set to zero for a case if scaled c predictor value lies outside range of historical RI values for a given predictor. c Setzero is true if RI probability is set to zero (default); otherwise it is set false. c Probridisc: SHIPS-RII probabilities in a 8(row) X 1 (column) array for use in computing ensemble c RII that includes the Bayesian and Logistic RI models. c RII thresholds contained in the array are (from row 1 to row 8) as follows: c 12h/20kt,24h/25kt,24h/30kt,24h/35kt,24h/40kt,36h/45kt,48h/55kt,72h/65kt c Parameter (Mft=28) Parameter (Nlt=5) Parameter (mxh=100) Parameter (Nx=5) Parameter (Nindx=8) Parameter (Mxprob=4) Parameter (Nthrss=10) Real Rvar(nindx,Nthrss), sclvar(Nindx,Nthrss) Real Rvarsrt(Nindx,Nthrss) Real Sclvrd(Nindx,Nthrss) Real Ratsc(nindx),Ratscd(nindx),prscld(nindx) Real Rmnval(Nindx,Nthrss), Rmxval(Nindx,Nthrss) Real Avgval(Nindx,Nthrss) Real Sclav(Nindx,0:NX), Prbri(Nindx,0:Nx) Real Scaled(Nindx), Suscld(Nindx) Real sumwgt(nindx) Real dwgt(Nindx,Nthrss), Clrisc(Nindx) Real dwgtsrt(Nindx,Nthrss) Real totwgt(nindx,Nthrss) Real totalwgt(nthrss) Integer Ithrs(Nindx) Real Rmnvalsrt(Nindx,Nthrss), Rmxvalsrt(Nindx,Nthrss) Real sclvarsrt(Nindx,Nthrss), sclvrdsrt(nindx,nthrss) c c Integer ishgc(0:Mft), id200(0:Mft) Integer irhcn(0:Mft), icflx(0:Mxh), Ipw19(0:Mft) Real vsstc(0:Mft) Real shgcrilt(nlt), d200rilt(nlt) Real rhcnrilt(nlt), potrilt(nlt) Real cflxrilt(nlt), tpwrilt(nlt) Real prbridisc(nlt,mxprob), nprob(nlt) Real probridisc(nindx) c Character*6 Labdat(Nthrss) Character*8 atcfid Character*10 sname character*28 thlabs(nthrss) Character*6 Labdatsrt(nindx,nthrss), labdattmp Character*28 Labtmp Character*28 labsrt(nindx,nthrss) c Logical Flag(nindx), newper, newpot, operational Logical setzero Logical skip(nindx) Logical debug Logical sclflg(nindx,nthrss) c Data setzero/.true./ Data debug/.false./ Data nudgefactor/1.0/ Data sclmax /1.0/, Sclmin/0.0/ Data ithrs /20,25,30,35,40,45,55,65/ c c RII predictor ranges c c 12-h/20kt c c label: PER , D200 , VP10 , SHGC , PC02 , NOHC , VMAX , IR02 , PW19 , CFLX / DATA (rmnval(1, j), j = 1, nthrss) * / -22.0, -27.5, 27.9, 79.2,-188.1, 1.8, * 27.0, 27.0, 0.0,-105.6/ DATA (rmxval(1, j), j = 1, nthrss) * / 60.5, 174.9, 145.2, 369.6, 220.0, 107.8, 148.5, * 378.4, 469.7, 942.7/ DATA (avgval(1, j), j = 1, nthrss) * / 15.3, 65.7, 83.1, 165.3, -9.1, 28.8, 66.6, * 98.3, 15.1, 190.4/ c c 24-h/25kt,30kt,35kt,40kt c c 24-h/25kt c c label: PER , D200 , VP10 , SHGC , PC02 , NOHC , VMAX , IR02 , PW19 , CFLX / DATA (rmnval(2, j), j = 1, nthrss) * / -33.0, -33.0, 29.7, 84.6,-302.5, 0.0, 22.5, * 20.7, 0.0, -82.5/ DATA (rmxval(2, j), j = 1, nthrss) * / 49.5, 170.5, 148.5, 401.5, 220.0, 107.8, 137.5, * 378.4, 595.1, 800.8/ DATA (avgval(2, j), j = 1, nthrss) * / 11.2, 62.0, 92.1, 172.9, -15.3, 27.6, 56.4, * 114.1, 15.3, 157.7/ c c c 24-h/30kt c c label: PER , D200 , VP10 , SHGC , PC02 , NOHC , VMAX , IR02 , PW19 , CFLX / DATA (rmnval(3, j), j = 1, nthrss) * / -22.0, -33.0, 36.9, 84.6,-225.5, 2.7, 22.5, * 20.7, 0.0, -82.5/ DATA (rmxval(3, j), j = 1, nthrss) * / 44.0, 170.5, 148.5, 348.7, 220.0, 107.8, 132.0, * 378.4, 581.9, 800.8/ DATA (avgval(3, j), j = 1, nthrss) * / 12.3, 64.9, 92.3, 169.5, -13.3, 29.1, 57.4, * 110.6, 12.8, 161.6/ c c 24-h/35kt c label: PER , D200 , VP10 , SHGC , PC02 , NOHC , VMAX , IR02 , PW19 , CFLX / c c label: PER , D200 , VP10 , SHGC , PC02 , NOHC , VMAX , IR02 , PW19 , CFLX / DATA (rmnval(4, j), j = 1, nthrss) * / -16.5, -33.0, 40.5, 84.6,-189.2, 2.7, 27.0, * 20.7, 0.0, -74.8/ DATA (rmxval(4, j), j = 1, nthrss) * / 44.0, 170.5, 147.4, 348.7, 184.8, 107.8, 132.0, * 378.4, 383.9, 733.7/ DATA (avgval(4, j), j = 1, nthrss) * / 13.1, 67.1, 92.1, 166.8, -13.5, 29.8, 58.4, * 106.8, 11.1, 162.6/ c c c 24-h/40kt c c label: PER , D200 , VP10 , SHGC , PC02 , NOHC , VMAX , IR02 , PW19 , CFLX / DATA (rmnval(5, j), j = 1, nthrss) * / -11.0, -33.0, 42.3, 84.6,-189.2, 2.7, 31.5, * 20.7, 0.0, -74.8/ DATA (rmxval(5, j), j = 1, nthrss) * / 44.0, 167.2, 147.4, 348.7, 163.9, 107.8, 126.5, * 308.0, 367.4, 733.7/ DATA (avgval(5, j), j = 1, nthrss) * / 13.9, 68.8, 92.2, 163.0, -14.8, 31.5, 59.1, * 102.5, 10.6, 161.6/ c c 36-h/45kt c c label: PER , D200 , VP10 , SHGC , PC02 , NOHC , VMAX , IR02 , PW19 , CFLX / c DATA (rmnval(6, j), j = 1, nthrss) * / -11.0, -26.4, 51.3, 86.4,-302.5, 4.5, 22.5, * 20.7, 0.0, -71.5/ DATA (rmxval(6, j), j = 1, nthrss) * / 38.5, 158.4, 150.7, 334.4, 218.9, 106.7, 110.0, * 378.4, 353.1, 583.0/ DATA (avgval(6, j), j = 1, nthrss) * / 11.3, 67.5, 99.2, 168.7, -18.4, 30.9, 51.6, * 116.9, 9.7, 142.4/ c c 48-h/55kt c label: PER , D200 , VP10 , SHGC , PC02 , NOHC , VMAX , IR02 , PW19 , CFLX / c DATA (rmnval(7, j), j = 1, nthrss) * / -11.0, -11.0, 62.1, 98.1,-302.5, 6.3, 22.5, * 30.6, 0.0, -58.3/ DATA (rmxval(7, j), j = 1, nthrss) * / 33.0, 154.0, 151.8, 301.4, 184.8, 106.7, 88.0, * 402.6, 270.6, 453.2/ DATA (avgval(7, j), j = 1, nthrss) * / 9.6, 68.4, 105.0, 169.9, -35.0, 31.3, 45.9, * 127.5, 7.8, 123.1/ c c 72-h/65kt c c label: PER , D200 , VP10 , SHGC , PC02 , NOHC , VMAX , IR02 , PW19 , CFLX / DATA (rmnval(8, j), j = 1, nthrss) * / -16.5, -8.8, 59.4, 107.1,-240.9, 7.2, 27.0, * 30.6, 0.0, -45.1/ DATA (rmxval(8, j), j = 1, nthrss) * / 27.5, 154.0, 152.9, 279.4, 136.4, 99.0, 82.5, * 392.7, 190.3, 339.9/ DATA (avgval(8, j), j = 1, nthrss) * / 7.5, 64.3, 108.7, 178.2, -46.0, 31.0, 40.5, * 144.0, 10.7, 114.7/ c c RII probabilities c c 12-h/20kt RI probilities c DATA (sclav (1, j), j = 0, nx) * / 0.00, 3.36, 8.20, 9.01, 9.69, 10.77/ DATA (prbri (1, j), j = 0, nx) * / 0.00, 1.67, 15.84, 30.99, 49.07, 69.33/ c c c 24-h/25kt,30kt,35kt,40kt RI probabilities c c 24-h/25kt c DATA (sclav (2, j), j = 0, nx) * / 0.00, 5.98, 10.53, 11.21, 11.85, 12.82/ DATA (prbri (2, j), j = 0, nx) * / 0.00, 3.43, 25.04, 40.53, 54.69, 76.92/ c c 24-h/30kt c DATA (sclav (3, j), j = 0, nx) * / 0.00, 3.64, 8.87, 9.49, 10.12, 11.12/ DATA (prbri (3, j), j = 0, nx) * / 0.00, 2.24, 20.17, 38.06, 46.64, 76.13/ c c 24-h/35kt c DATA (sclav (4, j), j = 0, nx) * / 0.00, 3.04, 8.30, 8.96, 9.64, 10.64/ DATA (prbri (4, j), j = 0, nx) * / 0.00, 1.59, 17.72, 28.90, 43.28, 80.00/ c c 24-h/40kt c DATA (sclav (5, j), j = 0, nx) * / 0.00, 2.21, 7.68, 8.37, 9.08, 10.11/ DATA (prbri (5, j), j = 0, nx) * / 0.00, 1.06, 15.37, 23.83, 42.07, 75.64/ c c 36-h/45kt RI probabilities c DATA (sclav (6, j), j = 0, nx) * / 0.00, 3.26, 9.43, 10.05, 10.76, 11.89/ DATA (prbri (6, j), j = 0, nx) * / 0.00, 1.67, 25.62, 29.54, 46.11, 78.64/ c c c 48-h/55kt RI probabilities c DATA (sclav (7, j), j = 0, nx) * / 0.00, 2.90, 10.24, 10.99, 11.70, 12.83/ DATA (prbri (7, j), j = 0, nx) * / 0.00, 1.43, 21.14, 29.72, 55.26, 74.07/ c DATA (sclav (8, j), j = 0, nx) * / 0.00, 1.75, 7.59, 8.36, 9.06, 9.99/ DATA (prbri (8, j), j = 0, nx) * / 0.00, 1.01, 14.22, 20.75, 54.10, 53.57/ c c c 72-h/65kt RI probabilities c c c RII weights c c 12-h/20kt weights c c label: PER , D200, VP10, SHGC, PC02, NOHC, VMAX, IR02, PW19, CFLX/ DATA (dwgt(1, j), j = 1, nthrss) * / 8.53, 1.76, 2.33, 2.44, 0.11, 1.10, 2.11, * 1.96,-0.09,-2.48/ c c 24-h/25kt,30kt,35kt,40kt weight c c 24-h/25kt c label: PER , D200, VP10, SHGC, PC02, NOHC, VMAX, IR02, PW19, CFLX/ DATA (dwgt(2, j), j = 1, nthrss) * / 4.21, 2.69, 4.41, 3.39,-0.08, 1.29, 2.34, * 3.20,-0.06,-1.81/ c c 24-h/30kt c label: PER , D200, VP10, SHGC, PC02, NOHC, VMAX, IR02, MTPW, CFLX/ c c label: PER , D200, VP10, SHGC, PC02, NOHC, VMAX, IR02, PW19, CFLX/ DATA (dwgt(3, j), j = 1, nthrss) * / 4.41, 2.96, 4.15, 3.14,-0.19, 1.66, 2.40, * 2.07, 0.66,-3.61/ c c 24-h/35kt c c label: PER , D200, VP10, SHGC, PC02, NOHC, VMAX, IR02, PW19, CFLX/ DATA (dwgt(4, j), j = 1, nthrss) * / 4.33, 3.09, 4.29, 3.27,-0.12, 1.67, 2.35, * 1.69, 0.39,-3.74/ c c DATA (dwgt(4, j), j = 1, nthrss) c 24-h/40kt c label: PER , D200, VP10, SHGC, PC02, NOHC, VMAX, IR02, PW19, CFLX/ DATA (dwgt(5, j), j = 1, nthrss) * / 3.49, 2.91, 3.94, 3.70, 0.08, 2.89, 2.08, * 1.17, 0.37,-3.64/ c c 36-h/45kt weights c label: PER , D200, VP10, SHGC, PC02, NOHC, VMAX, IR02, PW19, CFLX/ DATA (dwgt(6, j), j = 1, nthrss) * / 3.42, 3.16, 5.33, 3.71,-0.32, 1.38, 2.50, * 1.69, 0.91,-3.29/ c c 48-h/55kt weights c label: PER , D200, VP10, SHGC, PC02, NOHC, VMAX, IR02, PW19, CFLX/ DATA (dwgt(7, j), j = 1, nthrss) * / 2.24, 3.29, 6.13, 3.76, 0.75, 0.79, 2.13, * 2.18, 1.18,-2.81/ c c 72-h/65kt weights c label: PER , D200, VP10, SHGC, PC02, NOHC, VMAX, IR02, PW19, CFLX/ DATA (dwgt(8, j), j = 1, nthrss) * / 0.98, 2.46, 7.70, 2.44, 0.32, 1.13, 0.61, * 1.23, 0.70,-1.97/ c c Data clrisc/6.3,12.5,8.6,6.2,4.2,6.7,5.9,4.7/ c Data ideltat /12/, itimeinc/6/ Data imiss/9999/ c c label: PER , D200, VP10, SHGC, PC02, NOHC, VMAX, IR02, MTPW, CFLX/ c Data newper/.false./, newpot/.false./, operational/.true./ Data Thlabs /'12 HR PERSISTENCE (KT)','D200 (10**7s-1)', + 'POT = MPI-VMAX (KT)', 'MULTI-LAYER SHEAR (KT)', + '2nd PC OF IR BR TEMP','HEAT CONTENT (KJ/CM2)', + 'MAXIMUM WIND (KT)', + 'STD DEV OF IR BR TEMP ', + '%area of TPW <45 mm upshear', + 'BL DRY-AIR FLUX (W/M2)'/ c labdat Data Labdat /'deltv6','d200','potint','shgc','irpc', + 'rhcn','vmax','btstd','tpw','cflx'/ c c Write(lush,*)'setzero=',setzero c c Compute mean RII values for each of the time dependent c RII predictors for each of the forecast lead c times (12-h, 24-h, 36-h, and 48-h) c c c Do 50 J=1,nlt If(j.lt.nlt)then Jmx=J Elseif(J.eq.nlt)then Jmx=J+1 Endif nft = (Jmx*ideltat)/itimeinc nshgc = 0 nd200 = 0 nrhcn = 0 npot = 0 ncflx = 0 ntpw = 0 sumshgc = 0. sumd200 = 0. sumrhcn = 0. sumcflx = 0. sumpot = 0. sumtpw = 0. c Do 25 I=0,nft If(Ishgc(I).lt.9999)then nshgc = nshgc + 1 sumshgc = sumshgc + float(Ishgc(I)) Endif If(Id200(I).lt.9999)then nd200 = nd200 + 1 sumd200 = sumd200 + Float(Id200(I)) Endif If(Irhcn(I).lt.9999)then nrhcn = nrhcn + 1 sumrhcn = sumrhcn + Float(Irhcn(I)) Endif If(Icflx(I).ne.imiss)then ncflx = ncflx + 1 sumcflx= sumcflx + Float(Icflx(I)) Endif If(vsstc(I).lt.200.and.vmaxri.lt.200)then npot = npot + 1 sumpot = sumpot + (vsstc(I)-vmaxri) Endif If(Ipw19(I).lt.9999)then ntpw = ntpw + 1 sumtpw= sumtpw + (Ipw19(I)) Endif 25 Continue nmin = nft/2 + 1 If(nshgc.ge.nmin)then shgcrilt(J) = sumshgc/float(nshgc) Else shgcrilt(J) = 9999. Endif If(nd200.ge.nmin)then d200rilt(J) = sumd200/float(nd200) Else d200rilt(J) = 9999. Endif If(ncflx.ge.nmin)then cflxrilt(J) = sumcflx/Float(ncflx) Else cflxrilt(J) = 9999. Endif If(npot.ge.nmin)then potrilt(J) = sumpot/float(npot) Else potrilt(J) = 9999. Endif If(nrhcn.ge.nmin)then rhcnrilt(J) = sumrhcn/float(nrhcn) Else rhcnrilt(J) = 9999. Endif If(ntpw.ge.nmin)then tpwrilt(J) = sumtpw/float(ntpw) Else tpwrilt(J) = 9999. Endif If(debug)then write(40,*)'perri=',perri,'sbtri=',sbtri write(40,*)'rirpcri=',rirpcri,'vmaxri=',vmaxri write(40,*)'nrhcn=',nrhcn,'nmin=',nmin,'sumrhcn=', + sumrhcn write(40,*)'j=',j,'shgcrilt(j)=',shgcrilt(J), + 'd200rilt=',d200rilt(J),'cflxrilt=', + cflxrilt(J),'potrilt=',potrilt(J), + 'rhcn=',rhcnrilt(j),'tpwrilt=',tpwrilt(j) Endif 50 continue c c Put appropriate RI predictor values c in array for each RII probability threshold. c Since nindx values 2-5 all correspond to the c 24-h lead time they are assigned the same c (i.e. 24-h) lead-time RI predictor values c nindx = 1 (12-h lead time) c nindx = 2-5 (24-h lead time) c nindx = 6 (36-h lead time) c nindx = 7 (48-h lead time) c nindx = 8 (72-h lead time) c c label: PER, D200 , VP10 , SHG C , PC02 , NOHC , VMAX , IR02 , MTPW , CFLX / c Do 250 M=1,nindx c Rvar(m,1) = Perri Rvar(m,5) = rirpcri*100. Rvar(m,7) = vmaxri Rvar(m,8) = sbtri*10. c If(M.eq.1)then Rvar(m,2) = d200rilt(1) Rvar(m,3) = potrilt(1) Rvar(m,4) = shgcrilt(1) Rvar(m,6) = rhcnrilt(1) Rvar(m,9) = tpwrilt(1) Rvar(m,10) = cflxrilt(1) Elseif(M.ge.2.and.M.le.5)then Rvar(m,2) = d200rilt(2) Rvar(m,3) = potrilt(2) Rvar(m,4) = shgcrilt(2) Rvar(m,6) = rhcnrilt(2) Rvar(m,9) = tpwrilt(2) Rvar(m,10) = cflxrilt(2) Elseif(M.eq.6)then Rvar(m,2) = d200rilt(3) Rvar(m,3) = potrilt(3) Rvar(m,4) = shgcrilt(3) Rvar(m,6) = rhcnrilt(3) Rvar(m,9) = tpwrilt(3) Rvar(m,10) = cflxrilt(3) Elseif(M.eq.7)then Rvar(m,2) = d200rilt(4) Rvar(m,3) = potrilt(4) Rvar(m,4) = shgcrilt(4) Rvar(m,6) = rhcnrilt(4) Rvar(m,9) = tpwrilt(4) Rvar(m,10) = cflxrilt(4) Elseif(M.eq.8)then Rvar(m,2) = d200rilt(5) Rvar(m,3) = potrilt(5) Rvar(m,4) = shgcrilt(5) Rvar(m,6) = rhcnrilt(5) Rvar(m,9) = tpwrilt(5) Rvar(m,10) = cflxrilt(5) Endif Do 200 N=1,nthrss If(Rvar(m,N).lt.900 .or. (Rvar(m,N).lt.9000 .and. + Labdat(N).eq.'cflx') .or. (Rvar(m,N).lt.9000 + .and. Labdat(N).eq.'tpw'))then If(Labdat(N).eq.'deltv6' .and. newper .or. + Labdat(N).eq.'potint' .and. newpot .or. + Labdat(N).eq.'vmax')then If(Rvar(m,N) .LE. Avgval(M,N))then Sclvar(M,N) = 1 - + ((Avgval(M,N)- Rvar(m,N))/ + (avgval(M,N) - rmnval(M,N))) Else Sclvar(M,N) = 1 - + ((Rvar(M,N) - Avgval(M,N))/ + (rmxval(M,N) - avgval(M,N))) endif Elseif(labdat(N).eq.'shgc'.or. + labdat(n).eq.'btstd' .or. + labdat(n).eq.'irpc' .or. + labdat(n).eq.'cflx' .or. + labdat(n).eq.'tpw')then Sclvar(M,N) = (Rmxval(M,N) - Rvar(M,N))/ + (Rmxval(M,N) - Rmnval(M,N)) Elseif(Labdat(N).eq.'d200'.or. (labdat(N).eq. + 'potint' .and. .not.newpot) + .or.labdat(N).eq.'rhl'.or. + labdat(N).eq.'pxcnt'.or. + (labdat(N).eq.'deltv6' .and. .not.newper) + .or.labdat(N).eq.'rhcn')then Sclvar(M,N)=(Rvar(M,N) - Rmnval(M,N))/ + (Rmxval(M,N) - Rmnval(M,N)) Endif Else Sclvar(M,N) = 999. Endif 200 Continue 250 Continue Do 340 Mx=1,nindx Flag(mx)= .false. skip(mx)= .false. 340 continue Do 350 M=1,nindx Scaled(M)= 0.0 Do 300 NN=1,nthrss If(Sclvar(M,NN).gt.900)then Flag(m) = .true. Elseif(Sclvar(M,NN).Gt.sclmax + .and. sclvar(m,nn).lt.900.)then Sclvar(M,NN) = Sclmax Elseif(sclvar(M,NN).lt.sclmin)then If(setzero)then Scaled(M)= 0. Sclvar(M,NN) = Sclmin skip (M) = .true. Else Sclvar(M,NN)= Sclmin Endif Endif If(.not.skip(M) .and. .not.flag(M))then Scaled(m) = Scaled(m) + + Sclvar(M,NN)*dwgt(M,NN) endif 300 Continue 350 Continue Do 540 M=1,nindx If(.not.flag(m))then Do 525 n=0,nx-1 If(Scaled(M) .Ge. Sclav(m,nx))then prscld(M) = prbri(M,nx-1) + ( (prbri(M,nx) - + prbri(M,nx-1))/ + (Sclav(M,nx) - Sclav(M,nx-1))* + (Scaled(M) - Sclav(M,nx-1)) ) Elseif(Scaled(M) .ge. sclav(m,n) .and. + scaled(M).lt.sclav(m,n+1))then Prscld(M) = Prbri(m,n) + + ((Prbri(m,n+1) - Prbri(m,n))/ + (sclav(m,n+1) - Sclav(m,n))*(Scaled(M)-Sclav(m,n)) ) Endif 525 Continue If(prscld(M).gt.100.)prscld(M)=100. ratscd(M) = prscld(M)/clrisc(M) Else Scaled(M) = 999. Prscld(M) = 999. ratscd(M) = 999. Endif 540 Continue c c Calculate relative weights of the discrminant predictors c Do 555 M=1,nindx Sumwgt(M) = 0.0 Do 550 K=1,nthrss sumwgt(M) = sumwgt(M) + dwgt(M,K) 550 Continue 555 Continue c Do 565 M=1,nindx suscld(M) = 0.0 Lcount = 0 Do 560 L=1,nthrss If(sclvar(M,L).Lt.sclmin + .and. .not.flag(m)) Sclvar(M,L)=sclmin If(sclvar(M,L).gt.sclmax + .and. .not.flag(M))Sclvar(M,L)=sclmax If(sclvar(M,L).lt.900)then Lcount = Lcount + 1 If(operational)then sclvrd(M,L) = Float(nthrss)*dwgt(m,L)/sumwgt(M) + *sclvar(M,L) else sclvrd(M,L) = sclvar(M,L)*dwgt(M,L) Endif suscld(M) = suscld(M) + sclvrd(M,L) Else sclvrd(M,L) = 999. Endif 560 continue c If(Lcount.lt.nthrss)suscld(M) = 999. If(scaled(M).eq.0)suscld(M)=0.0 565 Continue c write(lush,575) atcfid,sname,ismon,isday,isyr,istime 575 format(/,' ** 2024 E. PACIFIC RI INDEX',1x,a8, + 1x,a10,1x,i2.2,'/',i2.2,'/',i2.2,2x,i2.2,' UTC **') c c Check to insure that probri25 > probri30 >probri35 > probri40) c for the 24-h lead time and correct probabilities based on the c climatological rate of RI c n24bgn = 2 n24end = 5 Do 580 N=n24bgn,n24end-1 If(prscld(N).lt.prscld(N+1))then prscld(N) = prscld(N) + (clrisc(N)-clrisc(N+1)) prscld(N+1) = prscld(N+1) - (clrisc(N)-clrisc(N+1)) Endif 580 continue c c c Do 585 N=n24bgn,n24end-1 If(prscld(N+1) .gt. prscld(N)) then prscld(N+1)=prscld(N) Endif 585 Continue c c Adjustment factor applied so probabilities of higher thresholds are < c lower thresholds by a percentage called nudgefactor c Do 586 N=n24bgn,n24end-1 If(prscld(N+1).lt.900 .and. prscld(N).lt.900)then If(prscld(N+1) .ge. prscld(N))then prscld(N+1) = prscld(N) - nudgefactor If(prscld(N+1).lt.0.0)prscld(N+1)=0.0 Endif Endif 586 continue c Do 590 N=1,nindx If(.not.flag(N))then If(prscld(N).gt.100)prscld(N)=100. If(prscld(N).lt.0)prscld(N) =0.0 Endif probridisc(N) = prscld(N) 590 continue c c Initialize arrays for sorting purposes c Do 655 m=1,nindx Do 654 nn=1,nthrss totwgt(m,nn)=0.0 654 continue 655 Continue c Do 656 j=1,nthrss totalwgt(j)= 0.0 656 continue c Do 625 JJ=1,nthrss totalwgt(JJ) = 0. Do 600 II=1,nindx totalwgt(JJ) = totalwgt(JJ) + + abs(dwgt(II,JJ)) totwgt(II,JJ) = totalwgt(JJ) 600 continue 625 Continue Do 650 m=1,nindx Do 640 kk=1,nthrss rvarsrt(m,kk) = rvar(m,kk) labsrt(m,kk) = thlabs(kk) labdatsrt(m,kk) = labdat(kk) totwgt(m,kk) = totalwgt(kk) rmnvalsrt(m,kk) = rmnval(m,kk) rmxvalsrt(m,kk) = rmxval(m,kk) sclvarsrt(m,kk) = sclvar(m,kk) sclvrdsrt(m,kk) = sclvrd(m,kk) dwgtsrt(m,kk) = abs(dwgt(m,kk)) 640 continue 650 continue c c Sort output arrays c Do 675 m=1,nindx Do 670 J=1,nthrss-1 K = nthrss - J Do 660 L=1,K If(totwgt(m,L).lt.totwgt(m,L+1))then tottmp = totwgt(m,L) wgttmp = abs(dwgtsrt(M,L)) rmntmp = rmnvalsrt(M,L) rmxtmp = rmxvalsrt(M,L) scltmp = sclvarsrt(M,L) sclvrdtmp = sclvrdsrt(M,L) labtmp = labsrt(M,L) labdattmp = labdatsrt(M,L) vartmp = rvarsrt(M,L) c totwgt(m,L) = totwgt(m,L+1) dwgtsrt(m,L) = abs(dwgtsrt(m,L+1)) rmnvalsrt(m,L)= rmnvalsrt(m,L+1) rmxvalsrt(m,L)= rmxvalsrt(m,L+1) sclvarsrt(m,L)= sclvarsrt(M,L+1) sclvrdsrt(m,L)= sclvrdsrt(M,L+1) labsrt(M,L) = labsrt(M,L+1) rvarsrt(M,L) = rvarsrt(M,L+1) labdatsrt(M,L) = labdatsrt(M,L+1) c totwgt(m,L+1) = tottmp dwgtsrt(m,L+1) = wgttmp rmnvalsrt(m,L+1)= rmntmp rmxvalsrt(m,L+1)= rmxtmp sclvarsrt(m,L+1)= scltmp sclvrdsrt(m,L+1)= sclvrdtmp labsrt(m,L+1) = labtmp rvarsrt(m,L+1)= vartmp labdatsrt(m,L+1) = labdattmp endif 660 continue 670 continue 675 continue c Do M=1,nindx Do Npred=1,nthrss sclflg(M,npred) = .false. Enddo Enddo c c mwrite = 3 Do 680 NN=1,nlt Nprob(NN) = 0 680 Continue Do 800 M=1,nindx If(M.eq.1)then write(lush,691) ithrs(mwrite) 691 format(1x,'(SHIPS-RII PREDICTOR TABLE for ',i2,' KT OR MORE', + ' MAXIMUM WIND INCREASE IN NEXT 24-h)') write(lush,*) sumtotwgt = 0.0 Do NN=1,nthrss sumtotwgt = sumtotwgt + sclvrdsrt(Mwrite,NN) if (setzero. and. sclvarsrt(Mwrite,NN).eq.sclmin)then sclflg(Mwrite,NN) = .true. endif Enddo Do 699 NN=1,nthrss If(prscld(mwrite).le.100.)then If(prscld(Mwrite).gt.0. + .and. .not. sclflg(mwrite,nn))then sclvrdsrt(mwrite,nn)= + (sclvrdsrt(mwrite,nn)/sumtotwgt)*prscld(mwrite) Elseif(.not.sclflg(mwrite,nn).and. + prscld(mwrite).eq.0.)then sclvrdsrt(mwrite,nn)= 999.0 Elseif(sclflg(mwrite,nn))then sclvrdsrt(mwrite,nn) = sclmin Endif else mxthrsh= nthrss do mm=1,mxthrsh sclvrdsrt(mwrite,nn) = 999. enddo endif If(Labdatsrt(mwrite,NN).eq.'shgc' .or. + labdatsrt(mwrite,NN).eq.'btstd' .or. + labdatsrt(mwrite,NN).eq.'tpw')then rvarsrt(mwrite,NN) = rvarsrt(mwrite,NN)*0.1 rmxvalsrt(mwrite,NN) = rmxvalsrt(mwrite,NN)*0.1 rmnvalsrt(mwrite,NN) = rmnvalsrt(mwrite,NN)*0.1 elseif(labdatsrt(mwrite,NN).eq.'irpc')then rvarsrt(mwrite,NN)= rvarsrt(mwrite,NN)*0.01 rmxvalsrt(mwrite,NN) = rmxvalsrt(mwrite,NN)*0.01 rmnvalsrt(mwrite,NN) = rmnvalsrt(mwrite,NN)*0.01 endif If(Labdatsrt(mwrite,NN).eq.'shgc'.or. + Labdatsrt(mwrite,NN).eq.'btstd'.or. + Labdatsrt(mwrite,NN).eq.'tpw'.or. + Labdatsrt(mwrite,NN).eq.'irpc'.or. + Labdatsrt(mwrite,NN).eq.'cflx')then Write(lush,698)labsrt(mwrite,NN), + rvarsrt(mwrite,NN), + rmxvalsrt(Mwrite,NN),rmnvalsrt(Mwrite,NN), + sclvarsrt(Mwrite,NN),sclvrdsrt(Mwrite,NN) Else If(NN.eq.1)then Write(lush,697) 697 format(5x,'Predictor',' Value ', + ' RI Predictor Range ', + ' Scaled Value(0-1)', ' % Contribution') Endif Write(lush,698)labsrt(mwrite,NN) + ,rvarsrt(mwrite,NN), + rmnvalsrt(Mwrite,NN),rmxvalsrt(Mwrite,NN), + sclvarsrt(Mwrite,NN),sclvrdsrt(Mwrite,NN) Endif 698 Format(1x,A28,':',1x,F6.1,4x,F5.1,2x,'to' + 2x,F5.1,6x,F6.2,9x,F5.1) 699 Continue Endif iprobt = nint(prscld(m)) if (iprobt .gt. 100.0 .and. iprobt .lt. 900) iprobt = 100. if (iprobt .gt. 900.) then iprobt = 999. ratscd(M)= 999. Endif If(m.eq.1)ihr=12 If(m.ge.2.and.m.le.5)ihr=24 If(m.eq.6)ihr=36 If(m.eq.7)ihr=48 If(m.eq.8)ihr=72 If(m.eq.1)write(lush,*) write(lush,790) ithrs(M),ihr,iprobt,ratscd(m), + clrisc(m) 790 Format(1x,'SHIPS Prob RI for ',I2,'kt/',I3,'hr', + ' RI threshold=',i4,'% is ',F5.1, + ' times climatological mean (',F4.1,'%)') If(M.eq.1)then numlt = 1 indx = 1 nprob(numlt) = nprob(numlt) + 1 Elseif(M.ge.2.and.M.le.5)then numlt = 2 indx = M-1 nprob(numlt) = nprob(numlt) + 1 Elseif(M.eq.6)then numlt = 3 indx = 1 nprob(numlt) = nprob(numlt) + 1 Elseif(M.eq.7)then numlt = 4 indx = 1 nprob(numlt) = nprob(numlt) + 1 Elseif(M.eq.8)then numlt= 5 indx = 1 nprob(numlt) = nprob(numlt) + 1 Endif Prbridisc(numlt,indx) = iprobt 800 Continue Do 850 nn = 1,nlt If(nprob(nn).lt.mxprob)then mnprob = nprob(nn)+1 Do 825 numindx = mnprob,mxprob prbridisc(nn,numindx)= 999. 825 continue Endif c If(debug)then If(nn.eq.1)then write(40,*) endif write(40,*)'nn=',nn,'probridisc=', + (prbridisc(nn,num),num=1,mxprob) endif 850 Continue If(debug)then Do 950 NN=1,nindx write(40,*)probridisc(NN) 950 continue Endif return end c