subroutine get_sal_clm(mask_ij,xlats_ij,xlons_ij,ny,nx,iy,im,id,ih,dsearch,miss_fill,bmiss,sal_clm) ! ! Abstract: get Salinity climatology at the valid time (atime) and target resolution ! (nx,ny) ! ! input & output ! real, dimension(nx*ny), intent(in) :: xlats_ij ! latitudes of target grids (nx*ny) real, dimension(nx*ny), intent(in) :: xlons_ij ! latitudes of target grids (nx*ny) integer, dimension(nx*ny), intent(in) :: mask_ij ! mask at target grids (0 =water, 1 = land, 2 = sea ice) (nx*ny) real, dimension(nx,ny), intent(out) :: sal_clm ! Salinity climatology valid at atime (nx,ny) real, intent(in) :: dsearch,bmiss integer, intent(in) :: iy,im,id,ih,nx,ny,miss_fill ! local declare real, allocatable, dimension(:,:) :: sal_clm0 ! Salinity climatology at the valid time integer, allocatable, dimension(:,:) :: cmask ! surface mask of Salinity climatology (set to zero = water everywhere) real, allocatable, dimension(:) :: cxlats ! latitudes of Salinity Climatology real, allocatable, dimension(:) :: cxlons ! longitudes of Salinity Climatology real, dimension(nx*ny) :: sal_clm_ij ! Salinity climatology at target grids (nx*ny) integer :: nxc,nyc,mon1,mon2,sfcflag character (len=6), parameter :: fin_sal_clm='salclm' ! Salinity climatology file name ! ! get which two months used and their weights from atime ! call get_tim_wei(iy,im,id,ih,mon1,mon2,wei1,wei2) ! ! get the dimensions of the Salinity climatology & allocate the related arrays ! call get_dim_nc(fin_sal_clm,nyc,nxc) allocate( sal_clm0(nxc,nyc),cmask(nxc,nyc),cxlats(nyc),cxlons(nxc) ) ! ! get sal_clm at the analysis time from monthly climatology & cxlats, cxlons ! call get_sal_clm_ta(sal_clm0,cxlats,cxlons,nyc,nxc,mon1,mon2,wei1,wei2) ! ! get mask at 1 degree based on if salinity value (valid or not) ! do j = 1, nyc do i = 1, nxc if ( sal_clm0(i,j) < 5.0 .or. sal_clm0(i,j) > 50.0 ) then ! non-water cmask(i,j) = 1 else ! water (ocean, lake,river) cmask(i,j) = 0 endif enddo enddo ! ! get sal_clm (nx by ny lat/lon) valid at atime ! if ( nx == nxc .and. ny == nyc ) then sal_clm(:,:) = sal_clm0(:,:) write(*,'(a,2F9.3)') 'same dimensions, sal_clm, min : ',minval(sal_clm),maxval(sal_clm) else write(*,'(a,4I8)') 'different dimensions,nx,ny,nxc,nyc : ',nx,ny,nxc,nyc sfcflag=0 call lalo_to_tile(sal_clm0, cmask, cxlats, cxlons, nyc, nxc, & sal_clm_ij, mask_ij, xlats_ij,xlons_ij,ny, nx, & sfcflag,2500.0,miss_fill,32.0) write(*,'(a,2F9.3)') 'done with lalo_to_tile for sal_clm, min : ',minval(sal_clm_ij),maxval(sal_clm_ij) sal_clm(:,:) = reshape (sal_clm_ij, (/nx,ny/) ) endif end subroutine get_sal_clm subroutine read_woa05_sal_nc(filename,sal,xlats,xlons,nlat,nlon,itime) ! abstract: read woa05 salinity monthly climatology (netcdf) use netcdf implicit none ! This is the name of the data file we will read. character (len=6), intent(in) :: filename integer, intent(in) :: nlat,nlon integer, intent(in) :: itime real, dimension(nlat), intent(out) :: xlats real, dimension(nlon), intent(out) :: xlons real, dimension(nlon,nlat), intent(out) :: sal ! Local variables integer :: ncid,ntime integer, parameter :: ndims = 4 character (len = *), parameter :: lat_name = "LAT" character (len = *), parameter :: lon_name = "LON" character (len = *), parameter :: z_name = "DEPTH" character (len = *), parameter :: t_name = "TIME" character (len = *), parameter :: sal_name="SALT" integer :: no_fill,fill_value integer :: time_varid,lon_varid, lat_varid, z_varid, sal_varid ! The start and count arrays will tell the netCDF library where to read our data. integer, dimension(ndims) :: start, count character (len = *), parameter :: units = "units" character (len = *), parameter :: sal_units = "psu" ! PSU (Practical SalinitUnit). 1 PSU = 1g/kg character (len = *), parameter :: time_units = "day_of_year" character (len = *), parameter :: lat_units = "degrees_north" character (len = *), parameter :: lon_units = "degrees_east" character (len = *), parameter :: z_units = "downward" integer :: missv ! Loop indices integer :: i,j ! Open the file. call nc_check( nf90_open(filename, nf90_nowrite, ncid) ) ! Get the varids of time, latitude, longitude & depth coordinate variables. call nc_check( nf90_inq_varid(ncid, "TIME", time_varid) ) call nc_check( nf90_inq_varid(ncid, "LAT", lat_varid) ) call nc_check( nf90_inq_varid(ncid, "LON", lon_varid) ) call nc_check( nf90_inq_varid(ncid, "DEPTH", z_varid) ) ! Read the time, latitude and longitude data. ! call nc_check( nf90_get_var(ncid, time_varid, ntime) ) call nc_check( nf90_get_var(ncid, lat_varid, xlats) ) call nc_check( nf90_get_var(ncid, lon_varid, xlons) ) ! Get the varids of the sal netCDF variables. call nc_check( nf90_inq_varid(ncid, "SALT",sal_varid) ) ! Read 1 record of nlat*nlon values, starting at the beginning ! of the record (the (1, 1, 1, rec) element in the netCDF file). start = (/ 1, 1, 1, itime /) count = (/ nlon, nlat, 1, 1 /) write(*,*) 'WOA05 itime : ',itime ! Read the sal data from the file, one record at a time. call nc_check( nf90_get_var(ncid, sal_varid, sal, start, count) ) ! Close the file. This frees up any internal netCDF resources ! associated with the file. call nc_check( nf90_close(ncid) ) ! If we got this far, everything worked as expected. Yipee! print *,"*** SUCCESS reading file ", filename, "!" end subroutine read_woa05_sal_nc subroutine get_dim_nc(filename,nlat,nlon) ! abstract: get dimensions of sal array use netcdf implicit none character (len=6), intent(in) :: filename integer, intent(out) :: nlat,nlon ! Local variables integer :: ncid integer :: LatDimID,LonDimID ! Open the file. call nc_check( nf90_open(filename, nf90_nowrite, ncid) ) ! Get dimensions call nc_check( nf90_inq_dimid(ncid,"LAT",LatDimID) ) call nc_check( nf90_inq_dimid(ncid,"LON",LonDimID) ) call nc_check( nf90_inquire_dimension(ncid,LatDimID,len=nlat) ) call nc_check( nf90_inquire_dimension(ncid,LonDimID,len=nlon) ) write(*,'(a,1x,a6,2I8)') 'get_dim_nc, file, nlat, nlon : ',filename,nlat,nlon ! Close the file. This frees up any internal netCDF resources ! associated with the file. call nc_check( nf90_close(ncid) ) ! If we got this far, everything worked as expected. Yipee! print *,"*** SUCCESS get dimensions from nc file ", filename, "!" end subroutine get_dim_nc subroutine get_sal_clm_ta(sal_clm_ta,xlats,xlons,nlat,nlon,mon1,mon2,wei1,wei2) !$$$ ! Abstract: get Tf/SST climatology at analysis time ! Created by Xu Li, March, 2019 implicit none ! Input integer, intent(in) :: nlat,nlon,mon1,mon2 real, intent(in) :: wei1,wei2 ! Output real, dimension(nlon,nlat), intent(out) :: sal_clm_ta real, dimension(nlon), intent(inout) :: xlons real, dimension(nlat), intent(inout) :: xlats !input/output data file names character (len=6), parameter :: fin_sal_clm='salclm' ! Local declare real, dimension(nlon,nlat) :: sal_clm1,sal_clm2 integer :: i,j ! ! read in RTG SST climatology without bitmap (surface mask) for mon1 and mon2 ! call read_woa05_sal_nc(trim(fin_sal_clm),sal_clm1,xlats,xlons,nlat,nlon,mon1) call read_woa05_sal_nc(trim(fin_sal_clm),sal_clm2,xlats,xlons,nlat,nlon,mon2) ! ! sal_clim at the analysis time ! do j = 1, nlat do i = 1, nlon if ( abs(sal_clm1(i,j)) < 50.0 .and. abs(sal_clm2(i,j)) < 50.0 ) then sal_clm_ta(i,j) = wei1*sal_clm1(i,j)+wei2*sal_clm2(i,j) elseif ( abs(sal_clm1(i,j)) < 50.0 .and. abs(sal_clm2(i,j)) > 50.0 ) then sal_clm_ta(i,j) = sal_clm1(i,j) elseif ( abs(sal_clm1(i,j)) > 50.0 .and. abs(sal_clm2(i,j)) < 50.0 ) then sal_clm_ta(i,j) = sal_clm2(i,j) endif enddo enddo end subroutine get_sal_clm_ta