! subroutines for reading all data from EOHC file as a single array !writen by Galina Chirokova 08/2018 !v2019CIRA_1.0.1 ! ! v 1.0.1 GC - removed links to filter2 subroutine ! v 1.0.2 GC - changed variable names from n_data_arr ! to n_eohc_data_arr ! v 1.1.0 SS - updated make_lsdiag_mftline to 7day ! v 1.1.1 MD - Removed comma from legacy write statements (Jan 2024) ! module eohc_module use ncoda_climo2_module,only : filter_out_bad_climo implicit none private public :: get_eohc_var_names public :: basin_upper_to_lower public :: find_eohc_file public :: make_lsdiag_mftline public :: get_final_eohc_profile public :: get_basin_specific_params_from_upper contains subroutine get_final_eohc_profile(mft,ilon,ilat,& iyear4_lsdiag,imonth_lsdiag,iday_lsdiag,ifhour_lsdiag,& fn_eohc, lu_eohc, & eohc_profile,aakm,tlat_inp,tlon_inp,& n_eohc_data_arr_max,n_eohc_data_arr,var_names_arr,& var_names_lsdiag_KT_final,var_names_lsdiag_KT,& var_names_lsdiag_NT_final,var_names_lsdiag_NT,& var_scales_arr,var_scales_arr_final,nsta,& eohc_profile_arr,eohc_days_old_fheader,& var_names_arr_save,lu_log_info,lu_log_debug,ierr_adj,& err_read_eohc_master,ierr_int,ierr_nsta,& nshay_data_start_year,& nshay_data_end_year,& reset_nshay_to_missing,& develop,izp,& iyear4_eohc, imonth_eohc, iday_eohc,& ifhour_eohc) !this subroutine returns EOHC profile that is: ! - interpolated to provided lat,lon location ! - filtered to remove bad data based on NTMX and NDFR ! - the profile also includes added NSTA variables - ! averaged temperature implicit none ! input !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! number of forecast times integer :: mft ! eohc filename character(len=200) :: fn_eohc !logical unit number for eohc_file integer :: lu_eohc ! integer lat array for lsdiag integer:: ilat(-2:mft) ! integer lon array for lsdiag integer:: ilon(-2:mft) ! eohc profile at given lat,lon real:: eohc_profile(n_eohc_data_arr_max) ! Distance (km) around the center for area-averaged SST real :: aakm ! number of EOHC variables - MUST be 26 for now integer :: n_eohc_data_arr ! latitude of the point to interpolate real :: tlat_inp ! longitude of the point to interpolate real :: tlon_inp !25 max EOHC variables (26 used in climo since dsst is also used) integer :: n_eohc_data_arr_max ! array of variables names for EOHC character(len=4) :: var_names_arr(n_eohc_data_arr_max) ! saved array of variables names for EOHC ! thsi is used if there is no available recent EOHC file ! and the subroutine to get names is not called at current step character(len=4) :: var_names_arr_save(n_eohc_data_arr_max) ! array of LSDIAG names for 25 eohc variables character(len=4) :: var_names_lsdiag_NT(n_eohc_data_arr_max) ! array of experimental LSDIAG names for 25 eohc variables character(len=4) :: var_names_lsdiag_KT(n_eohc_data_arr_max) ! array of scales for variables to write to LSDIAG files real :: var_scales_arr(n_eohc_data_arr_max) ! array of scales for variables to write to LSDIAG files + ! NSTA scale added after NSST scale real :: var_scales_arr_final(n_eohc_data_arr_max+1) ! !4-digit year from lsdiag header integer :: iyear4_lsdiag !2-digit month from lsdiag header integer :: imonth_lsdiag !2-digit day from lsdiag header integer :: iday_lsdiag !forecast hour from lsdiag header integer :: ifhour_lsdiag ! the 1st and last year of N Shay data ! All NTMX, NDTX, NTFR, NDFR values will be set to missing for ! NShay data ! 1st year of NShay data integer :: nshay_data_start_year ! last year of NShay data integer :: nshay_data_end_year ! flag to reset NShay NTMX, NDTX, NTFR, NDFR to missing logical :: reset_nshay_to_missing ! flag for developmenal vs real-time data ! -- develop = .true. - code will stop on critical error ! = .false. - code will stop on critical error logical :: develop ! periodicity flag ! for AL and EP/CP set izp = 0 integer :: izp ! output ! area-averaged SST NSTA estinmated for a single point real :: nsta ! filtered eohc profile at given lat,lon ! this is eohc_profile_filteres +NSTA added after NSST real :: eohc_profile_final(n_eohc_data_arr_max+1) ! array of LSDIAG names for 25 eohc variables +NSTA added after NSST character(len=4) :: var_names_lsdiag_NT_final(n_eohc_data_arr_max+1) ! array of experimental LSDIAG names for 25 eohc variables+KSTA added after ! KSST character(len=4) :: var_names_lsdiag_KT_final(n_eohc_data_arr_max+1) ! var_names_arr + NSTA character(len=4) :: var_names_arr_final(n_eohc_data_arr_max+1) ! integer missing value ! eohc_profile_arr - contains eohc data for all 25 vars for all 20 ! forecst ! times for all 25 eohc variables + NSTA added after NSST real:: eohc_profile_arr(n_eohc_data_arr_max+1,-2:mft) !eohc_file age in days calculated from eohc header time integer :: eohc_days_old_fheader ! !errors ! = 0 if values were adjusted ! = 1 if d_floor and t_max needed for adjustment ! are not available integer :: ierr_adj ! master error read_eohc integer :: err_read_eohc_master ! llintsp error integer :: ierr_int ! error getting area averaged NSTA integer :: ierr_nsta ! ! data from EOHC file ! heaader date integer :: iyear4_eohc, imonth_eohc, iday_eohc ! heaader hour integer :: ifhour_eohc !local ! filtered eohc profile at given lat,lon real:: eohc_profile_filtered(n_eohc_data_arr_max) ! eohc profile with NSTA added after NSST real:: eohc_profile_nsta(n_eohc_data_arr_max+1) ! iterator integer :: ivar ! temporary variable to store eohc variable name character(len=4) :: cvar ! temporary store the index of NSST variable in th earray integer :: sst_var_num ! missing value real real,parameter :: rmiss = -999.9 ! missing value integer real,parameter :: imiss = 9999 ! current forecast time iterator integer :: iftime ! ------------------------------------------------------------------------------------ ! - populate array eohc_profile_arr(n_eohc_data_arr,-2:mft) ! that array contains all 25 EOHC variables for all forecast ! times ! for a single lsdiag case ! ------------------------------------------------------------------------------------ ! initialize eohc_profile to missing eohc_profile(:) = rmiss ! initialize eohc_profile_filtered to missing eohc_profile_filtered(:) = rmiss ! initialize eohc_profile_final to missing eohc_profile_final(:) = rmiss ! initialize eohc_profile_nsta to missing eohc_profile_nsta(:) = rmiss ! initialize eohc_profile_arr to missing eohc_profile_arr(:,:) = rmiss !------------------------------------------------------------------------------------- ! get interpolated and filtered EOHC variables for each non-missing ! lat,lon point for all forecast times from -2 to mft !if ( .false. ) then !initialize eohc age to missing eohc_days_old_fheader = imiss if ( fn_eohc(1:4) .ne. 'none' ) then do iftime = -2,mft !----------------- if (( ilat(iftime) .lt. 9000 ) .and. (ilon(iftime) .lt. 9000 )) then tlat_inp = 0.1*float(ilat(iftime)) tlon_inp = -0.1*float(ilon(iftime)) !if ( loncon .eq. 'W' ) then ! tlon_inp = -0.1*float(ilon(iftime)) !endif ! get interpolated and filtered EOHC data for a single point ! AND add to that profile area-averaged temperature using ! 5-point avarage with radius aakm in km call get_eohc_profile(fn_eohc, iyear4_lsdiag, lu_eohc,& eohc_profile,aakm,n_eohc_data_arr,tlat_inp,& tlon_inp,var_names_arr,n_eohc_data_arr_max,nsta,& eohc_profile_filtered,lu_log_info,lu_log_debug,ierr_adj,& err_read_eohc_master,ierr_int,ierr_nsta,& nshay_data_start_year,& nshay_data_end_year,& reset_nshay_to_missing,& develop,izp,& iyear4_eohc, imonth_eohc, iday_eohc,& ifhour_eohc) ! calulate eohc age from the header data just once ! at forecast time t = 0 if ( iftime == 0 ) then if (err_read_eohc_master == 0 ) then call tdiff(iyear4_lsdiag,imonth_lsdiag,iday_lsdiag,ifhour_lsdiag,& iyear4_eohc,imonth_eohc,iday_eohc,ifhour_eohc,eohc_days_old_fheader) eohc_days_old_fheader = eohc_days_old_fheader/24 else eohc_days_old_fheader = imiss endif endif else eohc_profile(:) = rmiss var_names_arr = var_names_arr_save endif !----------------- write(lu_log_debug,*) 'INFO: get_final_eohc_profile: Interpolated profile & &for tlon,tlat, iftime ',tlon_inp, tlat_inp, iftime !----------------- ! populate eohc_profile_arr with added NSTA call get_eohc_profile_with_NSTA(n_eohc_data_arr_max,n_eohc_data_arr,var_names_arr,& eohc_profile,eohc_profile_filtered,eohc_profile_final,& var_names_lsdiag_KT_final,var_names_lsdiag_KT,& var_names_lsdiag_NT_final,var_names_lsdiag_NT,& var_scales_arr,var_scales_arr_final,nsta,& eohc_profile_nsta,lu_log_info,lu_log_debug,ierr_adj) !eohc_profile_arr(:,iftime)) eohc_profile_arr(:,iftime) = eohc_profile_nsta enddo else ! eohc file is none eohc_days_old_fheader = imiss nsta = rmiss ! populate eohc_profile_arr with added NSTA - ! should be all missing in this case but get the correct ! variables order call get_eohc_profile_with_NSTA(n_eohc_data_arr_max,n_eohc_data_arr,var_names_arr,& eohc_profile,eohc_profile_filtered,eohc_profile_final,& var_names_lsdiag_KT_final,var_names_lsdiag_KT,& var_names_lsdiag_NT_final,var_names_lsdiag_NT,& var_scales_arr,var_scales_arr_final,nsta,& eohc_profile_nsta,lu_log_info,lu_log_debug,ierr_adj) eohc_profile_arr(:,iftime) = eohc_profile_nsta endif return end subroutine get_final_eohc_profile subroutine get_eohc_profile(fn_eohc,iyear4,lu_eohc,& eohc_profile,aakm,n_eohc_data_arr,tlat_inp,& tlon_inp,var_names_arr,n_eohc_data_arr_max,nsta,& eohc_profile_filtered,lu_log_info,lu_log_debug,ierr_adj,& err_read_eohc_master,ierr_int,ierr_nsta,& nshay_data_start_year,& nshay_data_end_year,& reset_nshay_to_missing,& develop,izp,& iyear4_eohc, imonth_eohc, iday_eohc,& ifhour_eohc) !this subroutine returns EOHC profile that is: ! - interpolated to provided lat,lon location ! - filtered to remove bad data based on NTMX and NDFR implicit none ! input -------------------------------------------- ! eohc filename character(len=200) :: fn_eohc !logical unit number for eohc_file integer :: lu_eohc !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! eohc profile at given lat,lon real:: eohc_profile(n_eohc_data_arr_max) ! Distance (km) around the center for area-averaged SST real :: aakm ! number of EOHC variables - MUST be 26 for now integer :: n_eohc_data_arr ! latitude of the point to interpolate real :: tlat_inp ! longitude of the point to interpolate real :: tlon_inp !25 max EOHC variables (26 used in climo since dsst is also used) integer :: n_eohc_data_arr_max ! array of variables names for EOHC character(len=4) :: var_names_arr(n_eohc_data_arr_max) ! 4-digit year of eohc data integer :: iyear4 ! the 1st and last year of N Shay data ! All NTMX, NDTX, NTFR, NDFR values will be set to missing for ! NShay data ! 1st year of NShay data integer :: nshay_data_start_year ! last year of NShay data integer :: nshay_data_end_year ! flag to reset NShay NTMX, NDTX, NTFR, NDFR to missing logical :: reset_nshay_to_missing ! flag for developmenal vs real-time data ! -- develop = .true. - code will stop on critical error ! = .false. - code will stop on critical error logical :: develop ! periodicity flag ! for AL and EP/CP set izp = 0 integer :: izp !output --------------------------------------------- ! filtered eohc profile at given lat,lon real:: eohc_profile_filtered(n_eohc_data_arr_max) ! area-averaged SST NSTA estinmated for a single point real :: nsta !errors ! = 0 if values were adjusted ! = 1 if d_floor and t_max needed for adjustment ! are not available integer :: ierr_adj ! master error read_eohc integer :: err_read_eohc_master ! llintsp error integer :: ierr_int ! error getting area averaged NSTA integer :: ierr_nsta ! ! data from EOHC file ! heaader date integer :: iyear4_eohc, imonth_eohc, iday_eohc ! heaader hour integer :: ifhour_eohc ! local ------------------------------------------- ! iterator integer :: ivar ! temporary variable to store eohc variable name character(len=4) :: cvar ! float missing values real,parameter :: rmiss = -999.9 !initialize nsta to missing nsta = rmiss ! get interpolated EOHC data for a single point ! AND calculate for that profile area-averaged temperature using ! 5-point avarage with radius aakm in km write(lu_log_debug,*) 'INFO: get_eohc_profile: & &input to get_interp_eohc_profile fn_eohc,tlat_inp, tlon_inp ',& fn_eohc,tlat_inp, tlon_inp call get_interp_eohc_profile(fn_eohc,lu_eohc,& eohc_profile,aakm,n_eohc_data_arr,tlat_inp,& tlon_inp,var_names_arr,n_eohc_data_arr_max,nsta,& err_read_eohc_master,ierr_int,ierr_nsta,& lu_log_info,lu_log_debug,izp,& iyear4_eohc, imonth_eohc, iday_eohc,& ifhour_eohc) !DEBUG do ivar = 1,n_eohc_data_arr cvar = var_names_arr(ivar) write(lu_log_debug,*) "INFO: get_eohc_profile: & &Input Profile to get_filtered_profile ",& cvar, eohc_profile(ivar) enddo !filter depth of the isoterms in the single profile to !remove unrealistic isotherms with t > max T in the profile ! and points below th eocean floor call get_filtered_profile(var_names_arr,iyear4, eohc_profile,& eohc_profile_filtered,n_eohc_data_arr,n_eohc_data_arr_max,& lu_log_info,lu_log_debug,ierr_adj,& nshay_data_start_year,& nshay_data_end_year,& reset_nshay_to_missing,& develop) !DEBUG do ivar = 1,n_eohc_data_arr cvar = var_names_arr(ivar) write(lu_log_debug,*) "INFO: get_eohc_profile: & &Output Profile from get_filtered_profile ",& cvar, eohc_profile(ivar) enddo return end subroutine get_eohc_profile subroutine get_eohc_profile_with_NSTA(n_eohc_data_arr_max,n_eohc_data_arr,var_names_arr,& eohc_profile,eohc_profile_filtered,eohc_profile_final,& var_names_lsdiag_KT_final,var_names_lsdiag_KT,& var_names_lsdiag_NT_final,var_names_lsdiag_NT,& var_scales_arr,var_scales_arr_final,nsta,& eohc_profile_arr_iftime,lu_log_info,lu_log_debug,ierr_adj) ! add to interpolated EOHC profile NSTA implicit none !input ! eohc profile at given lat,lon real:: eohc_profile(n_eohc_data_arr_max) !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! number of EOHC variables - MUST be 26 for now integer :: n_eohc_data_arr !25 max EOHC variables (26 used in climo since dsst is also used) integer :: n_eohc_data_arr_max ! array of variables names for EOHC character(len=4) :: var_names_arr(n_eohc_data_arr_max) ! area-averaged SST NSTA estinmated for a single point real :: nsta ! filtered eohc profile at given lat,lon real:: eohc_profile_filtered(n_eohc_data_arr_max) ! filtered eohc profile at given lat,lon ! this is eohc_profile_filteres +NSTA added after NSST real :: eohc_profile_final(n_eohc_data_arr_max+1) ! array of LSDIAG names for 25 eohc variables character(len=4) :: var_names_lsdiag_NT(n_eohc_data_arr_max) ! array of experimental LSDIAG names for 25 eohc variables character(len=4) :: var_names_lsdiag_KT(n_eohc_data_arr_max) ! array of scales for variables to write to LSDIAG files real :: var_scales_arr(n_eohc_data_arr_max) !output ! array of LSDIAG names for 25 eohc variables +NSTA added after NSST character(len=4) :: var_names_lsdiag_NT_final(n_eohc_data_arr_max+1) ! array of experimental LSDIAG names for 25 eohc variables+KSTA added after ! KSST character(len=4) :: var_names_lsdiag_KT_final(n_eohc_data_arr_max+1) ! var_names_arr + NSTA character(len=4) :: var_names_arr_final(n_eohc_data_arr_max+1) ! array of scales for variables to write to LSDIAG files + ! NSTA scale added after NSST scale real :: var_scales_arr_final(n_eohc_data_arr_max+1) ! eohc_profile_arr - contains eohc data for all 25 vars for all 20 ! forecast ! times for all 25 eohc variables + NSTA added after NSST real :: eohc_profile_arr_iftime(n_eohc_data_arr_max+1) !error ! = 0 if values were adjusted ! = 1 if d_floor and t_max needed for adjustment ! are not available integer :: ierr_adj ! local ! iterator integer :: ivar ! temporary variable to store eohc variable name character(len=4) :: cvar ! temporary store the index of NSST variable in th earray integer :: sst_var_num ! missing value real,parameter :: rmiss = -999.9 ! populate eohc_profile_arr with added NSTA !1st pass add variables before NSTA and NSTA eohc_profile_arr_iftime(:) = rmiss do ivar = 1,n_eohc_data_arr cvar = var_names_arr(ivar) !write(lu_log_debug,*) "INFO: get_eohc_profile_with_NSTA: & ! &Input Profile1 for adding NSTA ",cvar, eohc_profile(ivar) if ( trim(cvar) .ne. 'NSST' ) then eohc_profile_final(ivar) = eohc_profile_filtered(ivar) var_names_arr_final(ivar) = var_names_arr(ivar) var_names_lsdiag_KT_final(ivar) = var_names_lsdiag_KT(ivar) var_names_lsdiag_NT_final(ivar) = var_names_lsdiag_NT(ivar) var_scales_arr_final(ivar) = var_scales_arr(ivar) else if ( trim(cvar) == 'NSST' ) then eohc_profile_final(ivar) = eohc_profile_filtered(ivar) var_names_arr_final(ivar) = var_names_arr(ivar) var_names_lsdiag_KT_final(ivar) = var_names_lsdiag_KT(ivar) var_names_lsdiag_NT_final(ivar) = var_names_lsdiag_NT(ivar) var_scales_arr_final(ivar) = var_scales_arr(ivar) eohc_profile_final(ivar+1) = nsta var_names_arr_final(ivar+1) = 'NSTA' var_names_lsdiag_KT_final(ivar+1) = 'KSTA' var_names_lsdiag_NT_final(ivar+1) = 'NSTA' var_scales_arr_final(ivar+1) = 10.0 sst_var_num = ivar ! break loop after adding NSTA exit endif enddo !2nd pass add variables after NSTA do ivar = (sst_var_num+2),(n_eohc_data_arr+1) write(lu_log_debug,*) "INFO: get_eohc_profile_with_NSTA: & &Input Profile2 for adding NSTA ", ivar, eohc_profile(ivar-1) eohc_profile_final(ivar) = eohc_profile_filtered(ivar-1) var_names_arr_final(ivar) = var_names_arr(ivar-1) var_names_lsdiag_KT_final(ivar) = var_names_lsdiag_KT(ivar-1) var_names_lsdiag_NT_final(ivar) = var_names_lsdiag_NT(ivar-1) var_scales_arr_final(ivar) = var_scales_arr(ivar-1) enddo eohc_profile_arr_iftime = eohc_profile_final return end subroutine get_eohc_profile_with_NSTA subroutine get_eohc_var_names(n_eohc_data_arr_max,var_names_arr,lu_log_info,lu_log_debug) !this subroutine assigns names to all variables form single eohc file implicit none ! input integer :: n_eohc_data_arr_max !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! output character(len=4) :: var_names_arr(n_eohc_data_arr_max) !arrays of 2d lons a regular 0.25 deg grid var_names_arr(1) = 'NLON' !arrays of 2d lats a regular 0.25 deg grid var_names_arr(2) = 'NLAT' !sst var_names_arr(3) = 'NSST' ! max temeperature in profile var_names_arr(4) = 'NTMX' ! depth of the max temperature var_names_arr(5) = 'NDTX' ! depth of mixed layer var_names_arr(6) = 'NDML' ! depth of 32 deg C isotherm var_names_arr(7) = 'ND32' ! depth of 30 deg C isotherm var_names_arr(8) = 'ND30' ! depth of 28 deg C isotherm var_names_arr(9) = 'ND28' ! depth of 26 deg C isotherm var_names_arr(10) = 'ND26' ! depth of 24 deg C isotherm var_names_arr(11) = 'ND24' ! depth of 22 deg C isotherm var_names_arr(12) = 'ND22' ! depth of 20 deg C isotherm var_names_arr(13) = 'ND20' ! depth of 18 deg C isotherm var_names_arr(14) = 'ND18' ! depth of 16 deg C isotherm var_names_arr(15) = 'ND16' ! depth of the ocean var_names_arr(16) = 'NDFR' ! temperature of the lowest point in the profile var_names_arr(17) = 'NTFR' ! ohc relative to 26 deg C var_names_arr(18) = 'NOHC' ! ohc relative to 20 deg C var_names_arr(19) = 'NO20' ! depth-averaged temperature assuming constant mixing depth = 25 m var_names_arr(20) = 'NT05' ! depth-averaged temperature assuming constant mixing depth = 50 m var_names_arr(21) = 'NT10' ! depth-averaged temperature assuming constant mixing depth = 80 m var_names_arr(22) = 'NT16' ! depth-averaged temperature assuming constant mixing depth = 100 m var_names_arr(23) = 'NT20' ! depth-averaged temperature assuming constant mixing depth = 125 m var_names_arr(24) = 'NT25' ! depth-averaged temperature assuming constant mixing depth = 200 m var_names_arr(25) = 'NT40' return end subroutine get_eohc_var_names subroutine get_eohc_array(fn_eohc,lu_eohc,& n_lons_dat, n_lats_dat,& data_2dr_all_tmp,& iday,imonth,iyr,ihour,& lon_west, lon_east, dlon, lat_south, lat_north, dlat,& ihfcst,smin,smax,krecl,ohc_source,err_read_eohc_master,& lu_log_info,lu_log_debug) !this subroutine combines all EOHC data from a single file into asingle !array; that is needed ! - to avoid having to repeadly type all 25 variables ! - to reuse subroutine from ncoda_climo2 implicit none ! local --------------------------------- !data grid parameters integer,parameter :: n_lons_max = 1440 integer,parameter :: n_lats_max = 522 integer,parameter :: n_eohc_data_arr = 25 ! input ------------------------------------ !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! eohc filename character(len = 200) :: fn_eohc ! this assignment will only happen ONCE ! at the 1st call to subroutine character(len = 200) :: last_eohc_file = 'none' !logical unit number for eohc file integer :: lu_eohc !output ------------------------------------ ! array containing all data arrays from a single EOHC datafile real :: data_2dr_all_tmp(n_lons_max,n_lats_max,n_eohc_data_arr) !io variables ! actual number of lon points in data integer :: n_lons_dat ! actual number of lat points in data integer :: n_lats_dat ! master error read_eohc integer :: err_read_eohc_master !header data !analysis time !ihour = 0 !placeholder for model analysis time in hours integer :: iday,imonth,iyr,ihour real :: lon_west, lon_east, dlon, lat_south, lat_north, dlat !ihfcst = 0 !placeholder for model forecast time in hours integer :: ihfcst !ohc26 min/max values real :: smin,smax !krecl = nlats * nlons = total number of data lines integer :: krecl !source of ohc data (ie where the data originally came from) ! for EOHC files that should be "NCODA_fullTprof_v2017_2" character(len = 200) :: ohc_source ! local --------------------------------- !arrays of data !arrays of 2d lons a regular 0.25 deg grid real :: elon_2dr(n_lons_max,n_lats_max) !arrays of 2d lats a regular 0.25 deg grid real :: elat_2dr(n_lons_max,n_lats_max) ! sst real :: sst_2dr(n_lons_max,n_lats_max) ! max temeperature in profile, deg C real :: t_max_2dr(n_lons_max,n_lats_max) ! depth of the max temperature, m real :: d_tmax_2dr(n_lons_max,n_lats_max) ! ohc relative to 20 deg C real :: ohc20_2dr(n_lons_max,n_lats_max) ! ohc relative to 26 deg C real :: ohc26_2dr(n_lons_max,n_lats_max) ! depth of mixed layer, m real :: dml_2dr(n_lons_max,n_lats_max) ! depth of the ocean, m real :: d_floor_2dr(n_lons_max,n_lats_max) ! temperature of the lowest point in the profile, deg C real :: t_floor_2dr(n_lons_max,n_lats_max) ! depth of 26 degC - 32 degC isotherms, m real :: d16_2dr(n_lons_max,n_lats_max) real :: d18_2dr(n_lons_max,n_lats_max) real :: d20_2dr(n_lons_max,n_lats_max) real :: d22_2dr(n_lons_max,n_lats_max) real :: d24_2dr(n_lons_max,n_lats_max) real :: d26_2dr(n_lons_max,n_lats_max) real :: d28_2dr(n_lons_max,n_lats_max) real :: d30_2dr(n_lons_max,n_lats_max) real :: d32_2dr(n_lons_max,n_lats_max) ! depth-averaged temperature (m) assuming constant mixing depth =200m real :: td200_2dr(n_lons_max,n_lats_max) ! depth-averaged temperature (m) assuming constant mixing depth =125m real :: td125_2dr(n_lons_max,n_lats_max) ! depth-averaged temperature (m) assuming constant mixing depth =100m real :: td100_2dr(n_lons_max,n_lats_max) ! depth-averaged temperature (m) assuming constant mixing depth =80m real :: td080_2dr(n_lons_max,n_lats_max) ! depth-averaged temperature (m) assuming constant mixing depth =50m real :: td050_2dr(n_lons_max,n_lats_max) ! depth-averaged temperature (m) assuming constant mixing depth =25m real :: td025_2dr(n_lons_max,n_lats_max) ! array containing all data arrays from a single EOHC datafile ! this array is used to save data between calls to this subroutine real :: data_2dr_all_local(n_lons_max,n_lats_max,n_eohc_data_arr) ! float missing value used in EOHC data files real,parameter :: rmiss = -999.9 !save local variables from header integer :: iday_save integer :: imonth_save integer :: iyr_save integer :: ihour_save real :: lon_west_save real :: lon_east_save real :: lat_south_save real :: lat_north_save real :: dlon_save real :: dlat_save real :: ihfcst_save real :: smin_save real :: smax_save real :: krecl_save character(len = 200) :: ohc_source_save integer :: n_lons_dat_save integer :: n_lats_dat_save ! keep values until next call save data_2dr_all_local save iday_save,imonth_save,iyr_save,ihour_save save lon_west_save, lon_east_save, dlon_save save lat_south_save, lat_north_save, dlat_save save ihfcst_save save smin_save,smax_save,krecl_save,ohc_source_save save n_lons_dat_save save n_lats_dat_save save last_eohc_file if ( trim(last_eohc_file) == (fn_eohc) ) then write(lu_log_debug,*) "INFO: get_eohc_array: & &using saved arrays from the old eohc file ", fn_eohc if (err_read_eohc_master .ne. 0 ) return data_2dr_all_tmp = data_2dr_all_local iday = iday_save imonth = imonth_save iyr = iyr_save ihour = ihour_save lon_west = lon_west_save lon_east = lon_east_save lat_south = lat_south_save lat_north = lat_north_save dlon = dlon_save dlat = dlat_save ihfcst = ihfcst_save smin = smin_save smax = smax_save krecl = krecl_save ohc_source = ohc_source_save n_lons_dat = n_lons_dat_save n_lats_dat = n_lats_dat_save else ! initialize all data to missing data_2dr_all_tmp(:,:,:) = rmiss write(lu_log_info,*) "INFO: get_eohc_array: Will now read EOHC data from ", & &fn_eohc write(lu_log_debug,*) "INFO: get_eohc_array: Will now read EOHC data from ", & &fn_eohc, " with nlons, nlats ",n_lons_max, n_lats_max call read_eohc_data(fn_eohc,lu_eohc,n_lons_max, n_lats_max,& n_lons_dat, n_lats_dat, elon_2dr, elat_2dr,& sst_2dr, t_max_2dr, d_tmax_2dr, dml_2dr, & d32_2dr, d30_2dr, d28_2dr, d26_2dr, & d24_2dr, d22_2dr, d20_2dr, d18_2dr, d16_2dr,& d_floor_2dr, t_floor_2dr, ohc26_2dr, ohc20_2dr,& td025_2dr, td050_2dr, td080_2dr, td100_2dr,& td125_2dr, td200_2dr,& iday,imonth,iyr,ihour,& lon_west, lon_east, dlon, lat_south, lat_north, dlat,& ihfcst,smin,smax,krecl,ohc_source, err_read_eohc_master,& lu_log_info,lu_log_debug) ! save current filename into last open file last_eohc_file = fn_eohc write(lu_log_info,*) "INFO: get_eohc_array: Done reading EOHC data from ", fn_eohc if (err_read_eohc_master .ne. 0 ) then write(lu_log_debug,*) "ERROR: get_eohc_array: AFter reading EOHC data & &err_read_eohc_master = ", err_read_eohc_master return endif write(lu_log_info,*) "INFO: get_eohc_array: Will now write EOHC variables to 3d array " data_2dr_all_tmp(:,:,1) = elon_2dr data_2dr_all_tmp(:,:,2) = elat_2dr data_2dr_all_tmp(:,:,3) = sst_2dr data_2dr_all_tmp(:,:,4) = t_max_2dr data_2dr_all_tmp(:,:,5) = d_tmax_2dr data_2dr_all_tmp(:,:,6) = dml_2dr data_2dr_all_tmp(:,:,7) = d32_2dr data_2dr_all_tmp(:,:,8) = d30_2dr data_2dr_all_tmp(:,:,9) = d28_2dr data_2dr_all_tmp(:,:,10) = d26_2dr data_2dr_all_tmp(:,:,11) = d24_2dr data_2dr_all_tmp(:,:,12) = d22_2dr data_2dr_all_tmp(:,:,13) = d20_2dr data_2dr_all_tmp(:,:,14) = d18_2dr data_2dr_all_tmp(:,:,15) = d16_2dr data_2dr_all_tmp(:,:,16) = d_floor_2dr data_2dr_all_tmp(:,:,17) = t_floor_2dr data_2dr_all_tmp(:,:,18) = ohc26_2dr data_2dr_all_tmp(:,:,19) = ohc20_2dr data_2dr_all_tmp(:,:,20) = td025_2dr data_2dr_all_tmp(:,:,21) = td050_2dr data_2dr_all_tmp(:,:,22) = td080_2dr data_2dr_all_tmp(:,:,23) = td100_2dr data_2dr_all_tmp(:,:,24) = td125_2dr data_2dr_all_tmp(:,:,25) = td200_2dr write(lu_log_info,*) "INFO: get_eohc_array: Will now copy EOHC variables tempa vars to save " ! copy header data to temp variable to save till next run data_2dr_all_local = data_2dr_all_tmp iday_save = iday imonth_save = imonth iyr_save = iyr ihour_save = ihour lon_west_save = lon_west lon_east_save = lon_east lat_south_save = lat_south lat_north_save = lat_north dlon_save = dlon dlat_save = dlat ihfcst_save = ihfcst smin_save = smin smax_save = smax krecl_save = krecl ohc_source_save = ohc_source n_lons_dat_save = n_lons_dat n_lats_dat_save = n_lats_dat endif return end subroutine get_eohc_array subroutine read_eohc_data(fn_eohc,lu_eohc,n_lons_max,n_lats_max,& n_lons_dat,n_lats_dat,elon_2dd,elat_2dd,& sst_2dd, t_max_2dd, d_tmax_2dd,dml_2dd,& d32_2dd,d30_2dd,d28_2dd,d26_2dd,d24_2dd,d22_2dd,d20_2dd,d18_2dd,d16_2dd,& d_floor_2dd,t_floor_2dd,ohc26_2dd,ohc20_2dd, & td025_2dd,td050_2dd,td080_2dd,td100_2dd,td125_2dd,td200_2dd,& iday,imonth,iyear,ihour,& lon_west, lon_east, dlon, lat_south, lat_north, dlat,& ihfcst,smin,smax,krecl,ohc_source, & err_read_eohc_master,lu_log_info,lu_log_debug) !read data from a single eohc file implicit none !input ------------------------------------------------- !grid parameters !integer,parameter :: n_lons_max = 1440 !integer,parameter :: n_lats_max = 522 integer :: n_lons_max integer :: n_lats_max !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug character(len = 200) :: fn_eohc integer :: lu_eohc !output ----------------------------------------------- ! actual number of lons in the file integer :: n_lons_dat ! actual number of lats in the file integer :: n_lats_dat !arrays of 2d lons,lats on a regular 0.25 deg grid real :: elon_2dd(n_lons_max,n_lats_max) real :: elat_2dd(n_lons_max,n_lats_max) !arrays of derived variables on 2d regular 0.25 deg lat/lon grid ! sst real :: sst_2dd(n_lons_max,n_lats_max) ! max temperature profile real :: t_max_2dd(n_lons_max,n_lats_max) ! depth of th emax T in a profile real :: d_tmax_2dd(n_lons_max,n_lats_max) ! OHC relative to 20 deg C real :: ohc20_2dd(n_lons_max,n_lats_max) ! OHC relative to 26 deg C real :: ohc26_2dd(n_lons_max,n_lats_max) ! depth of mixed layer real :: dml_2dd(n_lons_max,n_lats_max) ! depth of the ocean floor real :: d_floor_2dd(n_lons_max,n_lats_max) ! temperature of the lowest level in th eprofile real :: t_floor_2dd(n_lons_max,n_lats_max) ! depth of isotherms d16 - d32 real :: d16_2dd(n_lons_max,n_lats_max) real :: d18_2dd(n_lons_max,n_lats_max) real :: d20_2dd(n_lons_max,n_lats_max) real :: d22_2dd(n_lons_max,n_lats_max) real :: d24_2dd(n_lons_max,n_lats_max) real :: d26_2dd(n_lons_max,n_lats_max) real :: d28_2dd(n_lons_max,n_lats_max) real :: d30_2dd(n_lons_max,n_lats_max) real :: d32_2dd(n_lons_max,n_lats_max) !depth-avaraged temperature, deg C assuming constant dmix ! for dmix = 200 m real :: td200_2dd(n_lons_max,n_lats_max) ! for dmix = 125 m real :: td125_2dd(n_lons_max,n_lats_max) ! for dmix = 100 m real :: td100_2dd(n_lons_max,n_lats_max) ! for dmix = 80 m real :: td080_2dd(n_lons_max,n_lats_max) ! for dmix = 50 m real :: td050_2dd(n_lons_max,n_lats_max) ! for dmix = 25 m real :: td025_2dd(n_lons_max,n_lats_max) !io errors ! master error read_eohc integer :: err_read_eohc_master !header data !analysis time !ihour = 0 !placeholder for model analysis time in hours integer :: iday,imonth,iyear,ihour real :: lon_west, lon_east, dlon, lat_south, lat_north, dlat !ihfcst = 0 !placeholder for model forecast time in hours integer :: ihfcst !ohc26 min/max values real :: smin,smax !krecl = nlats * nlons = total number of data lines integer :: krecl !source of ohc data character(len = 200) :: ohc_source ! local ------------------------------------------- !io errors ! error open eohc file integer :: err_open ! error read eohc file header integer :: err_read_header ! error read eohc file data integer :: err_read_data ! error incorrect recored length in eohc file integer :: err_rec_length !missing value used for all variables in non-filled !eohc data real,parameter :: rmiss = -999.9 ! integer missing value integer,parameter :: imiss = 9999 ! iterator integer :: ii ! iterator integer :: jj integer :: line_count !open/read errors - initialize to missing err_open = imiss err_read_header = imiss err_read_data = imiss err_rec_length = imiss err_read_eohc_master = imiss !set all variables to missing elon_2dd = rmiss elat_2dd = rmiss !SST sst_2dd = rmiss !OHC; cp = 4000 J * K-1; pho = 1025 kg*m-3 ohc20_2dd = rmiss ohc26_2dd = rmiss !depth of mixed layer; negative values indicate temperature inversion !dml defined as upper level with (SST - Tlev) >= 0.5 deg C dml_2dd = rmiss !depth of the ocean d_floor_2dd = rmiss !temperature at the ocean floor t_floor_2dd = rmiss !depth of isotherms d16_2dd = rmiss d18_2dd = rmiss d20_2dd = rmiss d22_2dd = rmiss d24_2dd = rmiss d26_2dd = rmiss d28_2dd = rmiss d30_2dd = rmiss d32_2dd = rmiss !depth-averaged temperature td200_2dd = rmiss td125_2dd = rmiss td100_2dd = rmiss td080_2dd = rmiss td050_2dd = rmiss td025_2dd = rmiss write(lu_log_info,*) "INFO: read_eohc_data: & &reading data from file fn_eohc ", trim(fn_eohc) !open data file open(unit = lu_eohc,file = trim(fn_eohc),status='old',& form='formatted', IOSTAT=err_open) 950 format(/,'CRITICAL ERROR: read_eohc_data: Error opening file fn_eohc : ',i5,2x,a100) if (err_open /= 0) then write(lu_log_info,*) 'ERROR: read_eohc_data: fn_eohc >>',trim(fn_eohc),'<<' write(lu_log_info,950) err_open, trim(fn_eohc) close(lu_eohc) err_read_eohc_master = err_open return endif !read header read (lu_eohc, '(i4.4,1x,2i2.2,1x,i2.2,6f9.3,2x,i3.3,2f9.3,1x,i10,1x,a25)',& IOSTAT=err_read_header) & iyear, imonth, iday, ihour, & lon_west, lon_east, dlon, lat_south, lat_north, dlat, ihfcst, & smin,smax, krecl, ohc_source write(lu_log_info,*) "INFO: read_eohc_data: EOHC header year, month, day, hour, err ",& &iyear, imonth, iday, ihour,err_read_header !error reading header if (err_read_header /= 0) then 951 format(/,'CRITICAL ERROR: read_eohc_data: Error reading header from file: ',a60,1x,i10) write(lu_log_info,951) fn_eohc, err_read_header close(lu_eohc) err_read_eohc_master = err_read_header return endif n_lats_dat = int( (lat_north-lat_south)/dlat ) +1 n_lons_dat = int( (lon_east-lon_west)/dlon ) +1 write(lu_log_debug,1001) n_lats_dat,n_lons_dat,n_lats_dat*n_lons_dat, krecl 1001 format(/,'INFO: read_eohc_data: Test dimensions: N_Lats = ',i10, ', N_Lons = ',i10, & ', N_Lats*N_Lons = ', i20, ', krecl = ',i10) write(lu_log_debug,*) " INFO: header n_lats_dat, n_lons_dat ", n_lats_dat, n_lons_dat if (n_lats_dat*n_lons_dat .ne. krecl ) then write(lu_log_info,1002) n_lats_dat,n_lons_dat,n_lats_dat*n_lons_dat, krecl 1002 format(/,'ERROR: read_eohc_data: lat/lon numbers do not match: N_Lats = ',i6, & ', N_Lons = ',i6, ', N_Lats*N_Lons = ', i10, ', krecl = ',i10) err_rec_length= 1 close(lu_eohc) err_read_eohc_master = err_rec_length return endif !read data write(lu_log_debug,*) 'INFO: read_eohc_data: & &EOHC an_lats_dat, n_lons_dat ', n_lats_dat,n_lons_dat write(lu_log_debug,*) 'INFO: read_eohc_data: & &EOHC lat_north,lat_south ',lat_north,lat_south write(lu_log_debug,*) 'INFO: read_eohc_data: & &EOHC lon_east,lon_west ', lon_east,lon_west write(lu_log_debug,*) 'INFO: read_eohc_data: EOHC dlat,dlon', dlat, dlon line_count = 2 do jj=1,n_lats_dat do ii=1,n_lons_dat !keep all lons 0 - 360 read(lu_eohc,'(f9.3,f9.3,23f10.3)',IOSTAT=err_read_data) & elon_2dd(ii,jj), elat_2dd(ii,jj), & sst_2dd(ii,jj), t_max_2dd(ii,jj), d_tmax_2dd(ii,jj), & dml_2dd(ii,jj), d32_2dd(ii,jj), d30_2dd(ii,jj),& d28_2dd(ii,jj), d26_2dd(ii,jj), d24_2dd(ii,jj), d22_2dd(ii,jj),& d20_2dd(ii,jj), d18_2dd(ii,jj), d16_2dd(ii,jj), d_floor_2dd(ii,jj),& t_floor_2dd(ii,jj), ohc26_2dd(ii,jj), ohc20_2dd(ii,jj), td025_2dd(ii,jj),& td050_2dd(ii,jj), td080_2dd(ii,jj), td100_2dd(ii,jj), td125_2dd(ii,jj),& td200_2dd(ii,jj) line_count = line_count +1 !if (line_count .gt. 0 ) then ! write(lu_log_debug,*) 'DEBUG: ',& ! & ii*jj,line_count,ii,jj,elon_2dd(ii,jj),elat_2dd(ii,jj),td200_2dd(ii,jj) !endif enddo enddo write(lu_log_info,*) 'INFO: read_eohc_data: error status for reading EOHC data ', err_read_data ! set err_read_eohc_master to err_read_data for normal output err_read_eohc_master = err_read_data !error reading data if (err_read_data > 0) then 952 format(/,'ERROR: read_eohc_data: Error reading data from file: ',a100, i10) write(lu_log_info,952) fn_eohc, err_read_data close(lu_eohc) err_read_eohc_master = err_read_data return else if (err_read_data < 0) then 953 format(/,'INFO: read_eohc_data: EOF reached for eohc file: ',a100, i10) write(lu_log_info,952) fn_eohc, err_read_data close(lu_eohc) ! this is normal EOF condition err_read_eohc_master = 0 return !else if (err_read_data == 0) then else continue endif !close data file close (unit=lu_eohc) write(lu_log_info,*) 'INFO: read_eohc_data: normal return, & &err_read_eohc_master = ', err_read_eohc_master return end subroutine read_eohc_data subroutine basin_upper_to_lower(basin_upper, basin_lower,lu_log_info,lu_log_debug) ! convert upper case basin to lowercase basin implicit none ! input ! upper_case basin name character(len=2) :: basin_upper !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug !output ! lower_case basin name character(len=2) :: basin_lower select case (basin_upper) case('AL') basin_lower = 'al' case('EP') basin_lower = 'ep' case('CP') basin_lower = 'cp' case('WP') basin_lower = 'wp' case('SH') basin_lower = 'sh' case('IO') basin_lower = 'io' end select return end subroutine basin_upper_to_lower subroutine get_basin_specific_params_from_upper( basin_upper,& izp,& lu_log_info,lu_log_debug) ! get basin name using in EOHC filename from upper case basin implicit none !input ----------------------- !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! upper case basin name character(len=2) :: basin_upper ! output ---------------------- ! EOHC file abbreviation for basin character(len=3) :: eohc_basin ! periodicity flag ! for AL and EP/CP set izp = 0 integer :: izp select case (basin_upper) case('AL') izp = 0 case('EP') izp = 0 case('CP') izp = 0 case('WP') izp = 1 case('SH') izp = 1 case('IO') izp = 1 end select return end subroutine get_basin_specific_params_from_upper subroutine get_eohc_basin_from_upper(basin_upper, eohc_basin,lu_log_info,lu_log_debug) ! get basin name using in EOHC filename from upper case basin implicit none !input !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! upper case basin name character(len=2) :: basin_upper ! output ! EOHC file abbreviation for basin character(len=3) :: eohc_basin select case (basin_upper) case('AL') eohc_basin = 'atl' case('EP') eohc_basin = 'pac' case('CP') eohc_basin = 'pac' case('WP') eohc_basin = 'glb' case('SH') eohc_basin = 'glb' case('IO') eohc_basin = 'glb' end select return end subroutine get_eohc_basin_from_upper subroutine get_eohc_file_list_name(basin_upper, eohc_file_list_name,lu_log_info,lu_log_debug) ! get basin name using in EOHC filename from upper case basin implicit none ! input !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! upper case basin name character(len=2) :: basin_upper ! output ! name of the file containing list of eohc files character(len=13) :: eohc_file_list_name select case (basin_upper) case('AL') eohc_file_list_name = 'eohc_atl.list' case('EP') eohc_file_list_name = 'eohc_pac.list' case('CP') eohc_file_list_name = 'eohc_pac.list' case('WP') eohc_file_list_name = 'eohc_glb.list' case('SH') eohc_file_list_name = 'eohc_glb.list' case('IO') eohc_file_list_name = 'eohc_glb.list' end select return end subroutine get_eohc_file_list_name subroutine get_eohc_fname(basin_upper, ihra, iyear4,& imonth,iday,fn_eohc_base,year_str, date_str,lu_log_info,lu_log_debug) ! construct EOHC filename from date/time and uppercase basin implicit none ! input ! upper case basin name character(len=2) :: basin_upper !basin specified in the eohc filename character(len=3) :: eohc_basin !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! 4-digit year integer, intent(in) :: iyear4 ! 2-digit month integer, intent(in) :: imonth ! 4-digit day integer, intent(in) :: iday ! number of hours to add to iyear4,imonth,iday integer :: ihra ! output ! eohc filename character(len = 21) :: fn_eohc_base ! local ! date string character(len = 8) :: date_str ! 4-digit year string character(len = 4) :: year_str ! 4-digit year for time + ihra hours integer :: iyear4_add ! 2-digit month for time + ihra hours integer :: imonth_add ! 2-digit day for time + ihra hours integer :: iday_add ! forecast hour: only needed for tadd subroutine integer,parameter :: ifhour = 0 ! forecast hour for time + ihra hours integer :: ifhour_add call get_eohc_basin_from_upper(basin_upper, eohc_basin,lu_log_info,lu_log_debug) write(lu_log_debug,*) "INFO: get_eohc_fname: eohc_basin ", eohc_basin if ( ihra .eq. 0 ) then iyear4_add = iyear4 imonth_add = imonth iday_add = iday ifhour_add = ifhour else call tadd(iyear4 ,imonth, iday ,ifhour ,ihra,& iyear4_add,imonth_add,iday_add,ifhour_add) endif write(date_str,'(i4.4,i2.2,i2.2)')iyear4_add,imonth_add,iday_add write(year_str,'(i4.4)')iyear4_add fn_eohc_base = 'eohc_'//eohc_basin//'_'//date_str//'.dat' write(lu_log_debug,*) "INFO: get_eohc_fname: fn_eohc_base ", fn_eohc_base return end subroutine get_eohc_fname subroutine read_eohc_list_from_file( basin_upper, max_num_eohc_files, & &eohc_fnames_list,ListOpenErr,ListReadErr,& &lu_log_info,lu_log_debug,lu_list) ! reads eohc data from list file implicit none ! input -------------------------------------- ! max possible number of eohc files integer :: max_num_eohc_files ! logical unit number for eohc_file_list_name integer :: lu_list !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! upper case 2-leter basin name character(len=2) :: basin_upper !output --------------------------------------- ! list of eohc files character(len=200) :: eohc_fnames_list(max_num_eohc_files) ! error flags ! error opening list file integer :: ListOpenErr ! error reading list file integer :: ListReadErr !local ---------------------------------- ! name of file containing list of eohc files for current basin character(len = 13) :: eohc_file_list_name integer :: list_idx ! integer missing value integer :: imiss = 9999 ! initialize errors to missing ListOpenErr = imiss ListReadErr = imiss !initialize eohc_fnames_list to none eohc_fnames_list(:) = 'none' ! get the name of eohc file list for current basin call get_eohc_file_list_name(basin_upper, eohc_file_list_name,lu_log_info,lu_log_debug) ! open EOHC filename list open(file=eohc_file_list_name,unit=lu_list,form='formatted',status='old',iostat=ListOpenErr) if (ListOpenErr /= 0 ) then write(lu_log_info,*) "ERROR: read_eohc_list_from_file: & &unable to open file ", eohc_file_list_name write(lu_log_info,*) "ERROR: read_eohc_list_from_file: & &ListOpenErr = ", ListOpenErr write(lu_log_info,*) "ERROR: read_eohc_list_from_file: & &will return eohc_file_list_name filled with none " endif ! read all eohc files from list do list_idx = 1, max_num_eohc_files read(lu_list,100,iostat=ListReadErr) eohc_fnames_list(list_idx) 100 format(a200) ! exit loop if end of file !if (ListReadErr .ne. 0 ) exit ! iostst = 0: read successful if (ListReadErr == 0 ) then continue ! iostat > 0 - read eror else if (ListReadErr > 0 ) then write(lu_log_info,*) 'ERROR: read_eohc_list_from_file: & &error while reading eohc list file. ListReadErr = ', ListReadErr ! exit loop if read error exit ! iostat < 0 - end of file reached else if (ListReadErr < 0 ) then write(lu_log_info,*) 'INFO: read_eohc_list_from_file: & &end of eohc list file is reached. ListReadErr = ', ListReadErr ! exit loop if end of file is reached ! this is normal outcome - so set error to zero exit else write(lu_log_info,*) 'WARNING: read_eohc_list_from_file: & &Should not be here. Something wrong with reading eohclist data' write(lu_log_info,*) 'WARNING: read_eohc_list_from_file: & &Should not be here. ListReadErr = ', ListReadErr ! exit loop if something unexpected exit endif enddo ! close list file ! close EOHC filename list close(unit=lu_list) return end subroutine read_eohc_list_from_file subroutine find_eohc_file(eohc_root,get_eohc_file_from_list,& basin_upper, iyear4, imonth, iday, & ieohc_start_year, min_days_old, max_days_old,& fn_eohc,n_days_old,lu_log_info,lu_log_debug,lu_list,& ListOpenErr,ListReadErr) !this subroutine finds the eohc file for provided iyear4, imonth, iday !if that file does not exist, it will find the most recent file but no !older than max_days_old !if no recent eohc file exists, then the returned eohc filename will be !'none' !subroutime with return the found eohc filename and ! the n_days_old = how old is the file relative to provided date/time implicit none ! input ------------------------------------- ! option to use EOHC file list to check if file exist: ! - true : list of files will be used to determine if file exists ! - false : inquire statemnt will be used to see if file exists in ! filesystem logical :: get_eohc_file_from_list ! max possible number of eohc files for a single basin ! ( 365*(2019-2005) + 10 ) ! must be defined here to be able to save eohc_fnames_list_local integer,parameter :: max_num_eohc_files = 10000 ! eohc root directory, ! ONLY needed if not using list character(len=200) :: eohc_root ! full path to eohc file, ! ONLY needed if not using list character(len=300) :: path_eohc !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug !logical unit number for file containing list of eohc files integer :: lu_list ! upper case basin name (AL, EP, WP, IO, or SH) character(len=2) :: basin_upper ! first yesr for which eohc data will be processed ! this is used to speedup the code by NOT searcging eohc data for the ! years for which these data are not available integer :: ieohc_start_year ! 4-digit year integer, intent(in) :: iyear4 ! month integer, intent(in) :: imonth ! day integer, intent(in) :: iday ! max allowed age of EOHC file in days integer, intent(in) :: max_days_old ! min allowed age of EOHC file in days integer, intent(in) :: min_days_old ! output ------------------------------------- ! full eohc filename with path character(len = 200) :: fn_eohc ! age of the found EOHC file integer :: n_days_old ! error flags for ! open eohc list file integer :: ListOpenErr ! read eohc list file integer :: ListReadErr ! local ----------------------------------- ! eohc filename without path character(len = 21) :: fn_eohc_base ! 1st eohc filename from list character(len=21) :: first_eohc_file ! temp eohc filename without path character(len = 21) :: temp_fn_eohc_base ! date string (i.e. '20170630') character(len = 8) :: date_str ! year string (i.e. '2017') character(len = 4) :: year_str !true/false if eohc file does/does not exist logical :: res ! local variablefor converting dates to find old EOHC files integer :: n_hours_old ! number of hours to subtract from current day if ! current EOHC file is not available integer :: ihra ! index at which short fn_eohc name starts integer :: basename_idx ! list of eohc files character(len=200) :: eohc_fnames_list(max_num_eohc_files) ! index for list of eohc_files integer :: list_idx ! variables to determine if saved eohc list available ! will be initialized to .false. only when the subroutine is called 1st ! time logical :: have_saved_eohc_list_atl = .false. logical :: have_saved_eohc_list_pac = .false. logical :: have_saved_eohc_list_glb = .false. ! save eohc file list and information about which list is saved save eohc_fnames_list save have_saved_eohc_list_atl save have_saved_eohc_list_pac save have_saved_eohc_list_glb ! initialize list errors ListOpenErr = 0 ListReadErr = 0 ! initialize eohc filename to none fn_eohc = 'none' ! print to verify saved list exists ! if ( have_saved_eohc_list_atl) then ! write(lu_log_debug,*) "INFO: find_eohc_file: have_saved_eohc_list_atl is True" ! else ! write(lu_log_debug,*) "INFO: find_eohc_file: have_saved_eohc_list_atl is False" ! endif ! if ( have_saved_eohc_list_pac) then ! write(lu_log_debug,*) "INFO: find_eohc_file: have_saved_eohc_list_pac is True" ! else ! write(lu_log_debug,*) "INFO: find_eohc_file: have_saved_eohc_list_pac is False" ! endif ! if ( have_saved_eohc_list_glb) then ! write(lu_log_debug,*) "INFO: find_eohc_file: have_saved_eohc_list_pac is True" ! else ! write(lu_log_debug,*) "INFO: find_eohc_file: have_saved_eohc_list_pac is False" ! endif ! read the file with the list of EOHC file names if not already ! saved ! !NOTE: the code is using single variable eohc_fnames_list to store !file lists for all basin. Thus, when data for the new basin are !read it makes sure that only file_list for current basin is "saved" ! that shoud work as long as basins processed in groups ie all ! Atlantic, then all east/central Pacific, then all global) ! if ( get_eohc_file_from_list ) then select case (basin_upper) ! Atlantic case ('AL') !if ( ( have_saved_eohc_list_atl ) ) then ! continue if ( ( .not. have_saved_eohc_list_atl ) ) then !else call read_eohc_list_from_file(basin_upper,max_num_eohc_files,& &eohc_fnames_list,ListOpenErr,ListReadErr,& &lu_log_info,lu_log_debug,lu_list) have_saved_eohc_list_atl = .true. have_saved_eohc_list_pac = .false. have_saved_eohc_list_glb = .false. endif ! east Pacific case ('EP') !if ( ( have_saved_eohc_list_pac ) ) then ! continue !else if ( .not. have_saved_eohc_list_pac ) then call read_eohc_list_from_file(basin_upper,max_num_eohc_files,& &eohc_fnames_list,ListOpenErr,ListReadErr,& &lu_log_info,lu_log_debug,lu_list) have_saved_eohc_list_atl = .false. have_saved_eohc_list_pac = .true. have_saved_eohc_list_glb = .false. endif ! central Pacific case ('CP') if ( .not. have_saved_eohc_list_pac ) then call read_eohc_list_from_file(basin_upper,max_num_eohc_files,& &eohc_fnames_list,ListOpenErr,ListReadErr,& &lu_log_info,lu_log_debug,lu_list) have_saved_eohc_list_atl = .false. have_saved_eohc_list_pac = .true. have_saved_eohc_list_glb = .false. endif ! ! west Pacific ! ! west Pacific ! case ('WP') ! if ( .not. have_saved_eohc_list_glb ) then ! call read_eohc_list_from_file(basin_upper,max_num_eohc_files,& ! &eohc_fnames_list,ListOpenErr,ListReadErr,& ! &lu_log_info,lu_log_debug,lu_list) ! have_saved_eohc_list_atl = .false. ! have_saved_eohc_list_pac = .false. ! have_saved_eohc_list_glb = .true. ! endif ! ! Indian Ocean ! case ('IO') ! if ( .not. have_saved_eohc_list_glb ) then ! call read_eohc_list_from_file(basin_upper,max_num_eohc_files,& ! &eohc_fnames_list,ListOpenErr,ListReadErr,& ! &lu_log_info,lu_log_debug,lu_list) ! have_saved_eohc_list_atl = .false. ! have_saved_eohc_list_pac = .false. ! have_saved_eohc_list_glb = .true. ! endif ! ! Southern Hemisphere ! case ('SH') ! if ( .not. have_saved_eohc_list_glb ) then ! call read_eohc_list_from_file(basin_upper,max_num_eohc_files,& ! &eohc_fnames_list,ListOpenErr,ListReadErr,& ! &lu_log_info,lu_log_debug,lu_list) ! have_saved_eohc_list_atl = .false. ! have_saved_eohc_list_pac = .false. ! have_saved_eohc_list_glb = .true. ! endif end select ! chack if the returned eohc file list does not have valid files ! and set all saved lists to .false. first_eohc_file = eohc_fnames_list(1) if ( first_eohc_file(1:4) == 'none' ) then have_saved_eohc_list_atl = .false. have_saved_eohc_list_pac = .false. have_saved_eohc_list_glb = .false. ! if list has no valid eohc files then return return endif endif !find most recent available eohc file !if the eohc file is not available then find the latest available eohc file !but no more that 1 week old ! initialize n_days_old to zero n_days_old = min_days_old ! search for most recent eohc_file, if not found ! then search for previous days file. If no eohc ! files withing max_days_old ara available, will return none do while (n_days_old .le. max_days_old) n_hours_old = n_days_old *24 ! !hour to add to current time ihra = (-1)*n_hours_old !get eohc basename call get_eohc_fname(basin_upper, ihra, iyear4,& imonth,iday,fn_eohc_base,year_str, date_str,lu_log_info,lu_log_debug) ! use list of all eohc files to determin if file exists ! find matching eohc file in the list if that filename exists in ! the list res = .false. if ( get_eohc_file_from_list ) then ! EOHC data start in 2005 if ( iyear4 >= ieohc_start_year ) then do list_idx = 1, max_num_eohc_files ! get index of the start of the file basename fn_eohc = eohc_fnames_list(list_idx) basename_idx = len(trim( fn_eohc ))-20 temp_fn_eohc_base=trim( fn_eohc(basename_idx:)) !write(lu_log_debug,*) "DEBUG : find_eohc_file: & ! &filename from EOHC file list : TmpBase ",& ! &temp_fn_eohc_base, " Base ", fn_eohc_base ! set flag to true if filename found if ( temp_fn_eohc_base == fn_eohc_base ) then res = .true. ! file found => exit the loop exit else res = .false. endif enddo endif else ! use inquire statemnt to determine if file exists ! check file existence on the filesystem path_eohc=trim(eohc_root)//year_str//"/"//date_str//"/" fn_eohc=trim(path_eohc)//trim(fn_eohc_base) !write(lu_log_debug,*) "INFO: find_eohc_file: & ! &EOHC filename for basin, iyear4, imonth, iday is",& ! &basin_upper,iyear4, imonth, iday, fn_eohc !check if found eohc file exists inquire (file = fn_eohc, exist = res ) endif ! exit loop if EOHC file exists or look for EOHC file from previous ! day but no older than max_days_old if ( res ) then write(lu_log_debug,*) 'INFO: find_eohc_file: & &EOHC file eixsts confirmed ', fn_eohc exit else fn_eohc = 'none' endif ! increment days old by one and go to beginning of the loop ! to search for next file n_days_old = n_days_old+1 enddo return end subroutine find_eohc_file subroutine get_interp_eohc_profile(fn_eohc,lu_eohc,eohc_profile,& aakm, n_eohc_data_arr,tlat_inp,& tlon_inp,var_names_arr,n_eohc_data_arr_max,nsta,err_read_eohc_master,& ierr_int,ierr_nsta,lu_log_info,lu_log_debug,izp,& iyear4,imonth,iday,ifhour) ! this subroutine obtains interpolated EOHC profile for ! a given lat, lon, and date implicit none ! parameters ------------------- ! maximum allowed numberof lats, lons in eohc array integer, parameter :: n_lons_max = 1440 integer, parameter :: n_lats_max = 522 !if no recent eohc data avialble integer, parameter :: max_days_old = 7 ! float missing value used in EOHC files real, parameter :: rmiss = -999.9 ! input --------------------------- ! number of data arrays to read from EOHC file integer,intent(in):: n_eohc_data_arr ! maximum number of data arrays avilable in EOHC file ! should be == 25 !25 max EOHC variables (26 used in climo since dsst is also used) integer,intent(in):: n_eohc_data_arr_max ! Distance (km) around the center for area-averaged SST real :: aakm ! input EOHC filename character(len=200) :: fn_eohc ! logical unit for eohc filename integer :: lu_eohc ! array of variable names in th esame order as data in data_2dr_all(:,:,:) character(len=4) :: var_names_arr(n_eohc_data_arr_max) ! input lon, lat real :: tlat_inp, tlon_inp !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! periodicity flag ! for AL and EP/CP set izp = 0 integer :: izp ! output ---------------------------------------------------------- ! EOHC profile (all data for a single point) , only includes selected ! variables real :: eohc_profile(n_eohc_data_arr_max) ! area averaged sst at the lat,lon of eohc_profile real :: nsta ! error flags ! master error read_eohc integer :: err_read_eohc_master ! llintsp error integer :: ierr_int ! error getting area averaged NSTA integer :: ierr_nsta ! ! data from EOHC file ! heaader date integer :: iyear4, imonth, iday ! heaader hour integer :: ifhour ! local --------------------------------------------------------- ! single 3d array containing all EOHC data from a single file real :: data_2dr_all(n_lons_max,n_lats_max,n_eohc_data_arr) ! temporary array to store single array ( i.e. NSST or NOHC, etc ) from EOHC data !real, allocatable :: data_tmp(:,:) real :: data_tmp(n_lons_max,n_lats_max) ! domain boundaries eohc header real :: flon_west, flon_east, fdlon, flat_south, flat_north,fdlat ! data from EOHC header integer :: fihfcst real :: fsmin real :: fsmax integer :: fkrecl ! EOHC source (NCODA) character(len = 200) :: fohc_source ! first 4 charachters of eohc_filename character(len = 4) :: fn_eoch_4 ! temp lat,lon real :: tlat, tlon ! actual number of lats, lons in eohc array integer :: n_lons_dat, n_lats_dat ! iterator integer :: idata ! iterator integer :: i ! ----------------------------- !if valid recent eohc file exists, then fn_eoch_4 = trim(fn_eohc(1:4)) if ( fn_eoch_4 .eq. 'none' ) then !no recent eohc data avialble; all eohc data are set to rmiss !climatology should be used write(lu_log_debug,*) "INFO: get_interp_eohc_profile: no eohc file available " continue else !set all eohc data to missing data_2dr_all(:,:,:) = rmiss data_tmp(:,:) = rmiss write(lu_log_debug,*) "INFO: get_interp_eohc_profile: & &>>>>>>>>>>>>>>>> Will now get EOHC array ", fn_eohc call get_eohc_array(fn_eohc,lu_eohc,& n_lons_dat, n_lats_dat,& data_2dr_all,& iday,imonth,iyear4,ifhour,& flon_west, flon_east, fdlon, flat_south, flat_north,fdlat,& fihfcst,fsmin,fsmax,fkrecl,fohc_source,err_read_eohc_master,& lu_log_info,lu_log_debug ) write(lu_log_debug,*) "INFO: get_interp_eohc_profile: <<<<<<<<<<<<<<<<< Just got EOHC array, in if loop " endif write(lu_log_debug,*) "INFO: get_interp_eohc_profile: <<<<<<<<<<<<<<<<< Just got EOHC array " if ( err_read_eohc_master .ne. 0 ) then ! set output arrays to missing eohc_profile(:) = rmiss nsta = rmiss write(lu_log_info,*) "INFO: get_interp_eohc_profile: error reading eohc file, " write(lu_log_info,*) "INFO: get_interp_eohc_profile: & &err_read_eohc_master = ", err_read_eohc_master write(lu_log_info,*) "INFO: get_interp_eohc_profile: & &Returning all missing values for eohc interpolated profile and NSTA" return endif !get interpolated values of the eohc data array for lat, lon write(lu_log_debug,*) "INFO: get_interp_eohc_profile: Will now make sure lon between 0 - 360 " ! make sure lon is between 0 - 360 before calling interpolation ! 0-360 is the convention used by llintohc if (tlon_inp < 0 ) then tlon = 360.+tlon_inp else if (tlon_inp .ge. 360 ) then tlon = tlon_inp - 360. else tlon = tlon_inp endif tlat = tlat_inp write(lu_log_debug,*) 'INFO: get_interp_eohc_profile: n_lons_dat, n_lats_dat ',& & n_lons_dat, n_lats_dat, tlon, tlat, n_eohc_data_arr,n_eohc_data_arr_max write(lu_log_debug,*) 'INFO: get_interp_eohc_profile: & &flat_south,flat_north,fdlat,flon_west,flon_east,fdlon ',& &flat_south,flat_north,fdlat,flon_west,flon_east,fdlon !!!!!!!!!!!!!!reference for llintsp call ! subroutine llintsp(f1,slon1,slat1,dlon1,dlat1,im1,jm1,i1,j1, ! + fsp,rlon,rlat,izp,ierr) !!!!!!!!!END !!!!!reference for llintsp call eohc_profile(:) = rmiss do idata = 1, n_eohc_data_arr_max data_tmp = data_2dr_all(:,:,idata) call llintsp(data_tmp,flon_west,flat_south,fdlon,fdlat,n_lons_max,n_lats_max,& n_lons_dat,n_lats_dat,eohc_profile(idata),tlon,tlat,izp,ierr_int) if (ierr_int .ne. 0 ) then write(lu_log_info,1003) var_names_arr(idata), ierr_int,tlon,tlat,iyear4,imonth,iday 1003 format(/,'WARNING: get_interp_eohc_profile: & &llintsp exited with errors:ierr_int,tlon,tlat,ityr,itmon,itday ',& &a60,i5,f9.3,f9.3,i5,i5,i5) endif enddo write(lu_log_debug,*) 'DEBUG: get_interp_eohc_profile: eohc_profile ', eohc_profile write(lu_log_debug,*) "INFO: get_interp_eohc_profile: Will now add NSTA " ! get NSTA - area averaged NSST ! initialize nsta to missing if ( ierr_int > 1 ) then eohc_profile(:) = rmiss nsta = rmiss write(lu_log_debug,*) 'INFO: get_interp_eohc_profile: ierr_int = ',ierr_int write(lu_log_debug,*) 'INFO: get_interp_eohc_profile: & &will set eohc_profile(idata) and nsta to missing and return ' return else nsta = rmiss do idata = 1,n_eohc_data_arr_max if ( var_names_arr(idata) == 'NSST' ) then data_tmp = data_2dr_all(:,:,idata) call estimate_area_averaged_sst(aakm,data_tmp,flon_west,flat_south,& fdlon,fdlat,n_lons_max,n_lats_max,& n_lons_dat,n_lats_dat,eohc_profile(idata),& tlon,tlat,izp,nsta,ierr_nsta,lu_log_info,lu_log_debug) ! break cycle after getting nsta exit endif enddo endif write(lu_log_debug,*) "INFO: get_interp_eohc_profile: Will now do normal return from get_interp_eohc_profile " return end subroutine get_interp_eohc_profile subroutine get_filtered_profile(var_names_arr, iyear4, eohc_profile,& eohc_profile_filtered,n_eohc_data_arr,& n_eohc_data_arr_max,lu_log_info,lu_log_debug,ierr_adj,& nshay_data_start_year,& nshay_data_end_year,& reset_nshay_to_missing,& develop) ! this subroutine uses filter_out_bad_climo subroutine from ncoda_climo2.f90 ! to filter out bas data from EOHC profile ! will work with only selected EOHC variables ! however, do to correct filtering d_floor, t_max, and dt_max variables ! msut be availabe as these are used tinfiltering rules implicit none ! input ------------------------------------- ! input var_names_arr for current eohc_profile character(len=4) :: var_names_arr(n_eohc_data_arr_max) ! current number of eohc variables integer :: n_eohc_data_arr ! max possible number of eohc variables (should be ==25)) integer :: n_eohc_data_arr_max ! input non-filtered eohc_profile real :: eohc_profile(n_eohc_data_arr_max) !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! 4-digit year of eohc data integer :: iyear4 ! ! the 1st and last year of N Shay data ! All NTMX, NDTX, NTFR, NDFR values will be set to missing for ! NShay data ! 1st year of NShay data integer :: nshay_data_start_year ! last year of NShay data integer :: nshay_data_end_year ! flag to reset NShay NTMX, NDTX, NTFR, NDFR to missing logical :: reset_nshay_to_missing ! flag for developmenal vs real-time data ! -- develop = .true. - code will stop on critical error ! = .false. - code will stop on critical error logical :: develop ! output ------------------------------------- ! output filtered eohc_profile real :: eohc_profile_filtered(n_eohc_data_arr_max) !error ! = 0 if values were adjusted ! = 1 if d_floor and t_max needed for adjustment ! are not available integer :: ierr_adj ! local ------------------------------------- ! TODO: make rmiss and imiss global variables ! real missing value real,parameter :: rmiss = -999.9 ! real missing value integer,parameter :: imiss = 9999 !local variables character(len=4) :: var_names_arr_local(n_eohc_data_arr_max) ! iterator integer :: ivar ! ocean depth - temp variable real :: ndfr ! max profile temperature - temp variable real :: ntmx ! profile SST - temp variable real :: nsst ! temp variable name character(len=4) :: cvar ! temp variable name character(len=3) :: cvar3 ! temp variable name starting with X character(len=4) :: cvarX ! initialize local variable names to actual variable names var_names_arr_local = var_names_arr ! replace actual variable names by climo var names - for compatibility ! with local variable names in filter_out_bad_climo ! this does not change variables names for the user do ivar = 1,n_eohc_data_arr cvar = var_names_arr(ivar) !write(lu_log_debug,*) "INFO: get_filtered_profile: & ! &Input Profile to filter_out_bad_climo ", cvar, eohc_profile(ivar) if (cvar (2:4) == 'SST' ) then ! special name for SST cvarX = 'XNST' else ! get variable names starting with 'X' for compatibility ! with names used for filtering inside ! filter_out_bad_climo cvar3 = cvar(2:4) cvarX = 'X'//cvar(2:4) var_names_arr_local(ivar)=cvarX endif enddo ! save temporary variables for ntmx and ndfr: these variables ! are used for filtering, thus it is convinient to save them in a ! temporary variable ! --also, nsst variable will be used to replace NSST in ! filtereg profile in case filetring does not work do ivar = 1,n_eohc_data_arr cvar = var_names_arr_local(ivar) select case (cvar) case ('XTMX') ntmx = eohc_profile(ivar) case ('XDFR') ndfr = eohc_profile(ivar) case ('XNST') nsst = eohc_profile(ivar) end select enddo ! initialize filtered profile to non-filtered values eohc_profile_filtered = eohc_profile ! initialize ierr_adj to large value ierr_adj = imiss ! filter profile if ocean depth amd max profile temperature ara ! available if ( ( ndfr .gt. rmiss+1 ) .and. ( ntmx .gt. rmiss+1 ) ) then write(lu_log_debug,*) "INFO: get_filtered_profile: Will now call filter_out_bad_climo " call filter_out_bad_climo(var_names_arr_local,& eohc_profile_filtered,n_eohc_data_arr,ierr_adj,& lu_log_info,lu_log_debug,develop) endif write(lu_log_debug,*) "INFO: get_filtered_profile: ierr_adji,ntmx,ndfr ",& &ierr_adj,ntmx,ndfr ! if if ocean depth amd max profile temperature are not available, ! or filtering subroutine returned errors, ierr_adj will be non-zero ! in that case save minimal profile as filtered profile ! set all depths and depth-averaged temperature to missing since ! these are likely to have non-physical values without filtering ! ! NOTE: this situation should never occur in practice if ndfr and ntmx ! are provided. Thus, this is very unlikely to occur for just a single ! case if ( ierr_adj /= 0 ) then write(lu_log_info,*) "WARNING: get_filtered_profile:& & Please provide non-missing& & ndfr and ntmx to filter out bad data. Cannot filter out & & bad data, ierr_adj = ", ierr_adj write(lu_log_info,*) "INFO: get_filtered_profile: Will use reduced profile& & with only SST and OHC." ! set to missing everything except: ! NLON,NLAT,NSST,NDML, ND20, ND26, NDML - take these from ! non-filtered profile eohc_profile_filtered = rmiss do ivar = 1,n_eohc_data_arr cvar = var_names_arr(ivar) select case (cvar) case ('NLON') eohc_profile_filtered(ivar) = eohc_profile(ivar) case ('NLAT') eohc_profile_filtered(ivar) = eohc_profile(ivar) case ('NSST') eohc_profile_filtered(ivar) = eohc_profile(ivar) case ('NSTA') eohc_profile_filtered(ivar) = eohc_profile(ivar) case ('NOHC') eohc_profile_filtered(ivar) = eohc_profile(ivar) case ('NO20') eohc_profile_filtered(ivar) = eohc_profile(ivar) case ('NDML') eohc_profile_filtered(ivar) = eohc_profile(ivar) case ('ND26') eohc_profile_filtered(ivar) = eohc_profile(ivar) case ('ND20') eohc_profile_filtered(ivar) = eohc_profile(ivar) end select enddo endif ! set extra variables to missing in old data if ( reset_nshay_to_missing ) then if ( iyear4 >= nshay_data_start_year ) then if ( iyear4 <= nshay_data_end_year ) then write(lu_log_info,*) "INFO: get_filtered_profile: Will set NTMX, NDTX, NDFR, NTRF& & to missing in 1995 - 2004 N Shay data" ! set to missing everything except: ! NLON,NLAT,NSST,NDML, ND20, ND26, NDML - take these from ! non-filtered profile do ivar = 1,n_eohc_data_arr cvar = var_names_arr(ivar) select case (cvar) case ('NTMX') eohc_profile_filtered(ivar) = rmiss case ('NDTX') eohc_profile_filtered(ivar) = rmiss case ('NDFR') eohc_profile_filtered(ivar) = rmiss case ('NTFR') eohc_profile_filtered(ivar) = rmiss end select enddo endif endif endif return end subroutine get_filtered_profile subroutine make_lsdiag_mftline(mft,var,var_name,line_var_name,& vscale,n_days_old,ilat,ilon,& climo,lu_log_info,lu_log_debug,& iyear4, & reset_nshay_to_missing, & nshay_data_start_year,& nshay_data_end_year) ! make full lsdiag lines for all eohc variables and NSTA implicit none !input ------------------------------------------ ! number of forecast times integer, intent(in):: mft ! data forlsdiag variables for all forecast times real, intent(in):: var(-2:mft) ! scale for converting lsdiag variable from float to integer real, intent(in):: vscale !set to ! true : variable age info is not recorded ! false : variable age info is recorded logical :: climo !logical unit number for log file integer :: lu_log_info !logical unit number for debug log file integer :: lu_log_debug ! lsdiag variable name character(len=4) :: var_name ! age of lsdiag variable in days integer :: n_days_old ! array of integer lats for all forecast times integer :: ilat(-2:mft) ! array of integer lons for all forecast times integer :: ilon(-2:mft) ! 4-digit year of eohc data integer :: iyear4 ! ! the 1st and last year of N Shay data ! All NTMX, NDTX, NTFR, NDFR values will be set to missing for ! NShay data ! 1st year of NShay data integer :: nshay_data_start_year ! last year of NShay data integer :: nshay_data_end_year ! flag to reset NShay NTMX, NDTX, NTFR, NDFR to missing logical :: reset_nshay_to_missing ! output -------------------------------------- ! lsdiag line character(len=165) :: line_var_name ! local ------------------------------------- ! integer verion of lsdiag variable integer :: ivar(-2:mft) !iterator integer :: i !iterator integer :: k ! eohc float missing value ! TODO: move definition of rmiss value to main program real,parameter :: rmiss = -999.9 ! lsdiag integer missing value integer,parameter :: imiss = 9999 ! initialize integer verion of lsdiag variable to missing ivar(:) = imiss ! convert valid values of lsdiag variable to integer do i = -2, mft ivar(i) = imiss if ( ( ilat(i) .ne. imiss ) .and. (ilon(i) .ne. imiss ) ) then if ( var(i) .gt. rmiss+1 ) then ivar(i) = nint((var(i))*vscale) endif endif enddo ! make lsdiag line if ( climo ) then ! do not add field for variable age (n_days_old) for climo write (line_var_name,400) (ivar(k),k=0,mft),var_name 400 format(10x,29(1x,i4),1x,a4) 401 format(10x,29(1x,i4),1x,a4,1x,i4) else if ( reset_nshay_to_missing ) then ! !commented out part will set n_days_old to missing ! ! will skip that for now ! if ( ( iyear4 >= nshay_data_start_year ) .and. & ! & ( iyear4 <= nshay_data_end_year ) ) then ! ! add field for variable age (n_days_old) for regular variables ! n_days_old = imiss ! endif ! add field for variable age (n_days_old) for regular variables write (line_var_name,401) (ivar(k),k=0,mft),var_name,n_days_old else ! add field for variable age (n_days_old) for regular variables write (line_var_name,401) (ivar(k),k=0,mft),var_name,n_days_old endif write(lu_log_debug,*) "INFO: make_lsdiag_mftline: & &new lsdiag lines : var_name ",' ', line_var_name return end subroutine make_lsdiag_mftline end module eohc_module