module sst_module !version 19.0.4 !process sst data for ships !Written by Galina Chirokova for JHT project !This module processes SST data for SHIPS and includes: ! defifinition for 4 derived types (structures) ! - SST parameters for each SSt type ! type data_sst_params ! - all SST arrays used by iships ! type data_sst_set ! - SST related parameters for each model, ! including SHIPS, LGEM, Ric2o, Kaplan RII, Rozof RII, ! etclass, PrSEF, and AHI ! type model_sst_params ! - collection of all SST arrays for a single model ! type model_sst_set ! ! subroutines: ! init_data_sst_set - to initialise all SST data and parameters ! init_model_sst_set - to initialise all SST data and parameters ! select_best_sst - select best SSt from several SSt datasets and climo ! get_cooling_mpi - for input sst data get J Cione cooled and non-cooled mpi, ! and J Cione cooled and non-cooled vsst (mpi) ! get_model_sst - output for a given model final_sst and final_vsst to ! use base on reuested SST type and use/no_use of J Cione ! cooling ! get_sst_print_variables - prints SST data to char variables ! ! changelog: ! v19.0.1 - MD bug fixes: 200 format statements; ! missing outer loop in get_model_sst print section; ! add k to argument list for get_sst_print_variables; ! Fix bug in subroutine get_cooling_mpi; ! Fix bug so SST cooling not applied for missing SST ! v 19.0.2 - GC May 2019 ! - added loop to get_cooling_mpi to only run for ! non-missing SST values ! - cleaned comments ! v 19.0.3 - GC May 2019 ! - corrected loop to get_cooling_mpi to only run for ! non-missing SST values ! v 19.0.4 - MD Dec 2023 - legcay form of write statement changed implicit none private public get_cooling_mpi,select_best_sst,get_model_sst public get_sst_print_variables public data_sst_set public model_sst_set public data_sst_params public model_sst_params public init_data_sst_set public init_model_sst_set !declare SST structure type data_sst_params ! data lags in days ! read for SST from the last field in lsdiag integer :: isst_lag ! sst lag for the best sst - same as the lag of final sst used integer :: isst_best_lag end type data_sst_params type data_sst_set ! input sst for type (WSST, DSST, NSST, NSTA, DSTA) integer :: isst_inp ! input climo sst for type (WSST, DSST, NSST, NSTA, DSTA) integer :: isst_climo ! best sst for type (WSST, DSST, NSST, NSTA, DSTA) ! selected from data and climo integer :: isst_best ! ! final sst for each input type ! -- sst non-cooled real :: sst_final ! -- sst cooled real :: sstc_final ! -- vsst (rmpi adj by tr speed) non-cooled real :: vsst_final ! -- vsst (rmpi adj by tr speed) cooled real :: vsstc_final ! character variables for printing ! sst character(len=4) :: csst ! J Cione Colled SST character(len=4) :: csstc ! rmpi character(len=4) :: cvsst ! J Cione cooled rmpi character(len=4) :: cvsstc end type data_sst_set ! model sst - parameters - sst/vsst type , icione and maybe some other ! parameters grouped together type model_sst_params ! FLAGS --------------------------------- ! model sst type character(len=4) :: csst_type ! model vsst type character(len=4) :: cvsst_type ! ! J Cione model-specific flag for model sst integer :: icione_sst ! ! J Cione model-specific flag for model vsst integer :: icione_vsst end type model_sst_params ! model sst - these variables are used for each model included iniships: ! ships, lgem, ric2, rii_kaplen, rii_rozof, etclass, ahi, PrSEF type model_sst_set ! -- SST - VSST ------------------------- ! model sst - with or without J Cione cooling based on the flag ! use this for processing real :: sst_final ! model vsst - with or without J Cione cooling based on the flag ! use this for processing real :: vsst_final ! ! model sst - NOT cooled real :: sstn ! model sst - COOLED J CIone real :: sstc ! model vsst - NOT cooled real :: vsstn ! model vsst - COOLED J Cione real :: vsstc ! ! ----- PRINITING -------------------------- ! character variables for printing ! THES VARIABLES ARE ONLY USED FOR PRINITING SHIPS.TXT ! use sst and vsst for processing ! sst - NON cooled for th emodel character(len=4) :: csstn ! J Cione cooled SST for the model character(len=4) :: csstc ! rmpi NON cooled for teh model character(len=4) :: cvsstn ! J Cione COOLED rmpi for the model character(len=4) :: cvsstc end type model_sst_set contains subroutine init_data_sst_set(mft, var_sst_set,var_sst_params,& rmiss, imiss, cmiss) ! this subroutine initalizes all values in the model sst set ! to ! real - rmiss ! integer - imiss ! char - 'none' implicit none ! input ------------------- ! number of forecast times integer :: mft type(data_sst_set),dimension(0:mft) :: var_sst_set type(data_sst_params) :: var_sst_params real :: rmiss integer :: imiss character(len=4) :: cmiss ! local ------------------- ! iterator integer :: k ! sst lag in days for input sst - read from the last field in lsdiag file var_sst_params%isst_lag = imiss ! sst lag in days for the best sst - same as the lag of final sst used var_sst_params%isst_best_lag = imiss do k = 0,mft ! input sst for type (WSST, DSST, NSST, NSTA, DSTA) var_sst_set(k)%isst_inp = imiss ! input climo sst for type (WSST, DSST, NSST, NSTA, DSTA) var_sst_set(k)%isst_climo = imiss ! best sst for type (WSST, DSST, NSST, NSTA, DSTA) ! selected from data and climo var_sst_set(k)%isst_best = imiss ! ! final sst for each input type ! -- sst non-cooled var_sst_set(k)%sst_final = rmiss ! -- sst cooled var_sst_set(k)%sstc_final = rmiss ! -- vsst (rmpi adj by tr speed) non-cooled var_sst_set(k)%vsst_final = rmiss ! -- vsst (rmpi adj by tr speed) cooled var_sst_set(k)%vsstc_final = rmiss ! character variables for printing ! sst var_sst_set(k)%csst = cmiss ! J Cione Colled SST var_sst_set(k)%csstc = cmiss ! rmpi var_sst_set(k)%cvsst = cmiss ! J Cione cooled rmpi var_sst_set(k)%cvsstc = cmiss enddo return end subroutine init_data_sst_set subroutine init_model_sst_set(mft, temp_sst_set, temp_sst_params, & rmiss, imiss, cmiss) ! this subroutine initalizes all values in the model sst set ! to ! real - rmiss ! integer - imiss ! char - 'none' implicit none ! input --------------------- ! number of forecast times integer :: mft type(model_sst_set),dimension(0:mft) :: temp_sst_set type(model_sst_params) :: temp_sst_params real :: rmiss integer :: imiss character(len=4) :: cmiss ! local ----------------------- ! iterator integer :: k ! FLAGS --------------------------------- ! model sst type temp_sst_params%csst_type = cmiss ! model vsst type temp_sst_params%cvsst_type = cmiss ! ! J Cione model-specific flag for model sst temp_sst_params%icione_sst = imiss ! ! J Cione model-specific flag for model vsst temp_sst_params%icione_vsst = imiss ! do k = 0,mft ! -- SST - VSST ------------------------- ! model sst - with or without J Cione cooling based on the flag ! use this for processing temp_sst_set(k)%sst_final = rmiss ! model vsst - with or without J Cione cooling based on the flag ! use this for processing temp_sst_set(k)%vsst_final = rmiss ! ! model sst - NOT cooled temp_sst_set(k)%sstn = rmiss ! model sst - COOLED J CIone temp_sst_set(k)%sstc = rmiss ! model vsst - NOT cooled temp_sst_set(k)%vsstn = rmiss ! model vsst - COOLED J Cione temp_sst_set(k)%vsstc = rmiss ! ! ----- PRINITING -------------------------- ! character variables for printing ! THES VARIABLES ARE ONLY USED FOR PRINITING SHIPS.TXT ! use sst and vsst for processing ! sst - NON cooled for th emodel temp_sst_set(k)%csstn = cmiss ! J Cione cooled SST for the model temp_sst_set(k)%csstc = cmiss ! rmpi NON cooled for teh model temp_sst_set(k)%cvsstn = cmiss ! J Cione COOLED rmpi for the model temp_sst_set(k)%cvsstc = cmiss enddo return end subroutine init_model_sst_set subroutine get_sst_print_variables(k, mft, temp_sst_set) ! print sst, sstc, vsst. vsstc values ! to corresponding character variables that can be used for ! output implicit none !input ! number of forecast times integer :: mft ! iterator integer :: k ! input/output type(data_sst_set),dimension(0:mft) :: temp_sst_set ! local ! temp vsst var integer :: ivsst ! temp vsstc var integer :: ivsstc ! temp imiss integer,parameter :: imiss = 9999 ! initialize temporary integer variables ivsst = imiss ivsstc = imiss ! sst if (temp_sst_set(k)%sst_final .gt. 50.0) then temp_sst_set(k)%csst = ' N/A' else write(temp_sst_set(k)%csst,205) temp_sst_set(k)%sst_final 205 format(f4.1) endif !vsst !vsst - NON cooled if (temp_sst_set(k)%vsst_final .gt. 200.0) then temp_sst_set(k)%cvsst = ' N/A' else ivsst = ifix(0.49 + temp_sst_set(k)%vsst_final ) write( temp_sst_set(k)%cvsst,200) ivsst 200 format(i4) endif !vsstc - cooled if (temp_sst_set(k)%vsstc_final.gt. 200.0) then temp_sst_set(k)%cvsstc = ' N/A' else ivsstc = ifix(0.49 + temp_sst_set(k)%vsstc_final) write(temp_sst_set(k)%cvsstc,200) ivsstc endif return end subroutine get_sst_print_variables subroutine select_best_sst(mft, isst_max_lag_days ,& sst1_set, sst1_params, & sst2_set, sst2_params, & sst3_set, sst3_params,lulg ) ! select best available sst from sst1_set, sst2_set sst3_set and climo ! isst_inp from sst1_set SST used iv in valid range ! else isst_inp sst2_set is used if in valied range ! else isst_inp sst3_set is used if in valied range ! else isst_climo from sst1_set is used implicit none !input ---------------------- ! logical unit for log file integer :: lulg ! max allowed SST lag before using climo integer :: isst_max_lag_days type(data_sst_set),dimension(0:mft) :: sst1_set type(data_sst_set),dimension(0:mft) :: sst2_set type(data_sst_set),dimension(0:mft) :: sst3_set type(data_sst_params) :: sst1_params type(data_sst_params) :: sst2_params type(data_sst_params) :: sst3_params integer :: mft !output ----------------------- !selected best available sst integer :: ibest_sst_lag_days !local ----------------------- !iterator integer :: k ! set best SST equal to climo do k = 0,mft sst1_set(k)%isst_best = sst1_set(k)%isst_climo !write(lulg,*), "INFO: sst_module.f90: select_best_sst ",& ! sst1_set(k)%isst_inp, sst2_set(k)%isst_inp,sst3_set(k)%isst_inp,& ! sst1_params%isst_lag,sst2_params%isst_lag,sst3_params%isst_lag enddo do k = 0,mft ! if input SST1 is in valid range, and not too old, set best_sst to ! sst1 if ( ( sst1_set(k)%isst_inp .gt. 0 ) .and. ( sst1_set(k)%isst_inp .lt. 360) & .and. ( sst1_params%isst_lag .lt. isst_max_lag_days ) ) then sst1_set(k)%isst_best = sst1_set(k)%isst_inp sst1_params%isst_best_lag = sst1_params%isst_lag ! ELSE if input SST2 is in valid range, and not too old, set best_sst to ! sst1 else if ( ( sst2_set(k)%isst_inp .gt. 0 ) .and. (sst2_set(k)%isst_inp .lt. 360) & .and. ( sst2_params%isst_lag .lt. isst_max_lag_days ) ) then sst1_set(k)%isst_best = sst2_set(k)%isst_inp sst1_params%isst_best_lag = sst2_params%isst_lag ! ELSE if input SST3 is in valid range, and not too old, set best_sst to ! sst1 else if ( ( sst3_set(k)%isst_inp .gt. 0 ) .and. (sst3_set(k)%isst_inp .lt. 360) & .and. ( sst2_params%isst_lag .lt. isst_max_lag_days ) ) then sst1_set(k)%isst_best = sst3_set(k)%isst_inp sst1_params%isst_best_lag = sst3_params%isst_lag ! ELSE use climo; climo is provided as a 4th sst_det, and can be different ! from sst1, 2, and/or 3 else sst1_set(k)%isst_best = sst1_set(k)%isst_climo sst1_params%isst_best_lag = isst_max_lag_days endif enddo return end subroutine select_best_sst subroutine get_cooling_mpi(mft,rlat,ispeed,cmagt,& mpispd,ibasin,rmiss,ivmpi2,& mpi_sst_set ) implicit none !input------------START----------- integer :: mft integer :: ispeed real :: rlat(-2:mft) real :: cmagt(-2:mft) integer :: mpispd ! flag to adjust MPI for storm translational speed ! 0 to use mean speed adjusted mpi routine (mpical) ! 1 to use original mpi routine (mpicalo) and add ! fraction of storm speed from ! Schewerdt formula a =1.5*c**0.63 !icione cooling 0 - no, 1 - error (not used); 2 - yes !integer :: icione !basin number integer :: ibasin !missing values real :: rmiss !input SST mpi_sst_set%isst_inp - ! part of the declared structure mpi_sst_set type(data_sst_set),dimension(0:mft) :: mpi_sst_set !K Emanuel MPI - NOT used - placeholder integer,dimension(0:mft) :: ivmpi2 !input------------END----------- !output-----------START----------- ! these are output variables ! thes SHOULD be commented out ! since these are part of the mpi_set_sst structure that is ! already declared ! non cooled vsst ! !!mpi_sst_set(k)%vsst_final = rmpi ! ! J Cione cooled vsst !!mpi_sst_set(k)%vsstc_final = rmpic ! ! non cooled vsst !mpi_sst_set(k)%sst_final = sstt ! ! J Cione cooled vsst !mpi_sst_set(k)%sstc_final = ssttc ! !output-----------END----------- !local ---START--------------- !temporary lat real :: rlattm !temporary translation speed real :: ctm,ctemp,cadj !interator integer :: k !temporary point SST real :: sstt !temporary point cooled SST real :: ssttc !empirical MPI real :: rmpi !empirical cooled MPI real :: rmpic ! temp imiss integer,parameter :: imiss = 9999 !internal---END------------------ ! Calculate the mpi from the SST or the modified SST ! only do calculation for non-missing SST. the isst_best ! shoudl only be missing when there is no forecast track. ! do 10 k=0,mft if ( mpi_sst_set(k)%isst_best .ge. imiss ) then ! SST is missing, nothing to do. The return values ! weere already initialized by init_data_sst_set ! subroutine return else ! convert input sst to real sstt = 0.1*float(mpi_sst_set(k)%isst_best) ! !dropped support for icione = 1 rlattm = abs(rlat(k)) ctm = cmagt(k)/1.944 !upper limit 20 kt is needed since icione equation does not !work well for faster translation speeds if (ctm .gt. 20.0) ctm = 20.0 ! ! Calculate the SST cooling from J. Cione's equation ssttc = 1.1222793*sstt + 0.1425625*ctm - 0.0778590*rlattm & -3.3705640 ! make sure cooled SST is colder that SST if ( ssttc .gt. sstt ) ssttc = sstt !adjust MPI for storm translational speed if (mpispd .eq. 0) then call mpical(sstt ,rmpi, ibasin) call mpical(ssttc,rmpic,ibasin) else call mpicalo(sstt ,rmpi ,ibasin) call mpicalo(ssttc,rmpic,ibasin) ! if (cmagt(k) .lt. rmiss .and. cmagt(k) .gt. 0.0) then ctemp = cmagt(k) if (ctemp .gt. 40.0) ctemp=40.0 cadj = 1.5*(ctemp**0.63) rmpi = rmpi + cadj rmpic= rmpic+ cadj endif endif ! non cooled vsst mpi_sst_set(k)%vsst_final = rmpi ! J Cione cooled vsst mpi_sst_set(k)%vsstc_final = rmpic ! non cooled vsst mpi_sst_set(k)%sst_final = sstt ! J Cione cooled vsst mpi_sst_set(k)%sstc_final = ssttc ! ! Save appropriate MPI for LGE model ! if (icione .eq. 0) then ! vsstl(k) = rmpip ! else ! vsstl(k) = rmpipc ! endif ! vsstl(k) = float(ivmpi2(k)) ! if (vsstl(k) .lt. 20.0) vsstl(k) = 20.0 endif 10 continue return end subroutine get_cooling_mpi subroutine get_model_sst(mft,& temp_model_sst, & temp_model_params, & wsst_set, & dsst_set, & dsta_set, & nsst_set, & nsta_set, lulg ) implicit none !input -------------------------- ! different input sst data ! Reynolds Weekly type(data_sst_set),dimension(0:mft) :: wsst_set ! Reynolds Daily type(data_sst_set),dimension(0:mft) :: dsst_set ! Reynolds Daily spatially averaged type(data_sst_set),dimension(0:mft) :: dsta_set ! NCODA Daily type(data_sst_set),dimension(0:mft) :: nsst_set ! NCODA Daily spatially averaged type(data_sst_set),dimension(0:mft) :: nsta_set ! current model type(model_sst_set),dimension(0:mft) :: temp_model_sst type(model_sst_params) :: temp_model_params ! number of forecast times integer :: mft ! logical unit for log file integer :: lulg ! local ------------------------------------ ! iterator integer :: k ! NON cooled vsst integer - temporary integer :: ivsstn ! cooled vsst integer - temporary integer :: ivsstc !output ----------------------- ! for calulations !temp_model_sst%sst_final !temp_model_sst%vsst_final ! for printing ships output only !temp_model_sst%cvsstn !temp_model_sst%cvsstc ! define cooled and non_cooled sst and vsst do k = 0,mft ! WSST if (temp_model_params%csst_type .eq. 'wsst' ) then ! sst non cooled temp_model_sst(k)%sstn = wsst_set(k)%sst_final ! sstc cooled temp_model_sst(k)%sstc = wsst_set(k)%sstc_final ! vsst non cooled temp_model_sst(k)%vsstn = wsst_set(k)%vsst_final ! vsstc cooled temp_model_sst(k)%vsstc = wsst_set(k)%vsstc_final ! DSST elseif (temp_model_params%csst_type .eq. 'dsst' ) then ! sst non cooled temp_model_sst(k)%sstn = dsst_set(k)%sst_final ! sstc cooled temp_model_sst(k)%sstc = dsst_set(k)%sstc_final ! vsst non cooled temp_model_sst(k)%vsstn = dsst_set(k)%vsst_final ! vsstc cooled temp_model_sst(k)%vsstc = dsst_set(k)%vsstc_final ! DSTA elseif (temp_model_params%csst_type .eq. 'dsta' ) then ! sst non cooled temp_model_sst(k)%sstn = dsta_set(k)%sst_final ! sstc cooled temp_model_sst(k)%sstc = dsta_set(k)%sstc_final ! vsst non cooled temp_model_sst(k)%vsstn = dsta_set(k)%vsst_final ! vsstc cooled temp_model_sst(k)%vsstc = dsta_set(k)%vsstc_final ! NSST elseif (temp_model_params%csst_type .eq. 'nsst' ) then ! sst non cooled temp_model_sst(k)%sstn = nsst_set(k)%sst_final ! sstc cooled temp_model_sst(k)%sstc = nsst_set(k)%sstc_final ! vsst non cooled temp_model_sst(k)%vsstn = nsst_set(k)%vsst_final ! vsstc cooled temp_model_sst(k)%vsstc = nsst_set(k)%vsstc_final ! NSTA elseif (temp_model_params%csst_type .eq. 'nsta' ) then ! sst non cooled temp_model_sst(k)%sstn = nsta_set(k)%sst_final ! sstc cooled temp_model_sst(k)%sstc = nsta_set(k)%sstc_final ! vsst non cooled temp_model_sst(k)%vsstn = nsta_set(k)%vsst_final ! vsstc cooled temp_model_sst(k)%vsstc = nsta_set(k)%vsstc_final else write(lulg,*) "ERROR: sst_module.f90: get_model_sst ",& " Non supported SST type ", & temp_model_params%csst_type, & " Will STOP now " stop endif enddo ! select final sst and vsst based on icione flag do k = 0,mft if ( temp_model_params%icione_sst .eq. 2 ) then ! sst final temp_model_sst(k)%sst_final = temp_model_sst(k)%sstc ! vsst final temp_model_sst(k)%vsst_final = temp_model_sst(k)%vsstc else ! use non_cooled sst ! sst final temp_model_sst(k)%sst_final = temp_model_sst(k)%sstn ! vsst final temp_model_sst(k)%vsst_final = temp_model_sst(k)%vsstn endif enddo do k = 0,mft ! prepare variables for printing ! ! ! char SST - non cooled !temp_model_sst(k)%csstn ! sst - NON cooled if (temp_model_sst(k)%sstn .gt. 50.0) then temp_model_sst(k)%csstn = ' N/A' else write(temp_model_sst(k)%csstn,205) temp_model_sst(k)%sstn 205 format(f4.1) endif !vsst ! vsst - non cooled ! char VSST - non cooled !temp_model_sst(k)%cvsstn if (temp_model_sst(k)%vsstn .gt. 200.0) then temp_model_sst(k)%cvsstn = ' N/A' else ivsstn = ifix(0.49 + temp_model_sst(k)%vsstn ) write( temp_model_sst(k)%cvsstn,200) ivsstn 200 format(i4) endif !vsstc - cooled ! char VSSTC - cooled !temp_model_sst(k)%cvsstc if (temp_model_sst(k)%vsstc .gt. 200.0) then temp_model_sst(k)%cvsstc = ' N/A' else ivsstc = ifix(0.49 + temp_model_sst(k)%vsstc) write(temp_model_sst(k)%cvsstc,200) ivsstc endif enddo return end subroutine get_model_sst end module sst_module