Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

introduction of mask file for generating mask when appropriate #8

Merged
merged 4 commits into from
Mar 22, 2021
Merged

introduction of mask file for generating mask when appropriate #8

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
f59b5c3
changes for creating masks at runtime
Jan 28, 2021
66356cb
added ice_mesh_mod
Jan 29, 2021
6ff987c
implemented area correction factors as option
Feb 1, 2021
0f90ce1
more cleanup
Feb 1, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion cicecore/cicedynB/general/ice_init.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1604,7 +1604,7 @@ subroutine input_data
grid_type /= 'rectangular' .and. &
grid_type /= 'cpom_grid' .and. &
grid_type /= 'regional' .and. &
grid_type /= 'latlon' ) then
grid_type /= 'setmask' ) then
if (my_task == master_task) write(nu_diag,*) subname//' ERROR: unknown grid_type=',trim(grid_type)
abort_list = trim(abort_list)//":20"
endif
Expand Down
287 changes: 1 addition & 286 deletions cicecore/cicedynB/infrastructure/ice_grid.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -40,7 +40,7 @@ module ice_grid
private
public :: init_grid1, init_grid2, &
t2ugrid_vector, u2tgrid_vector, &
to_ugrid, to_tgrid, alloc_grid
to_ugrid, to_tgrid, alloc_grid, makemask

character (len=char_len_long), public :: &
grid_format , & ! file format ('bin'=binary or 'nc'=netcdf)
Expand Down Expand Up @@ -368,11 +368,6 @@ subroutine init_grid2
else
call popgrid ! read POP grid lengths directly
endif
#ifdef CESMCOUPLED
elseif (trim(grid_type) == 'latlon') then
call latlongrid ! lat lon grid for sequential CESM (CAM mode)
return
#endif
elseif (trim(grid_type) == 'cpom_grid') then
call cpomgrid ! cpom model orca1 type grid
else
Expand Down Expand Up @@ -882,286 +877,6 @@ subroutine popgrid_nc

end subroutine popgrid_nc

#ifdef CESMCOUPLED
!=======================================================================

! Read in kmt file that matches CAM lat-lon grid and has single column
! functionality
! author: Mariana Vertenstein
! 2007: Elizabeth Hunke upgraded to netcdf90 and cice ncdf calls

subroutine latlongrid

! use ice_boundary
use ice_domain_size
use ice_scam, only : scmlat, scmlon, single_column
use ice_constants, only: c0, c1, p5, p25, &
field_loc_center, field_type_scalar, radius
#ifdef USE_NETCDF
use netcdf
#endif

integer (kind=int_kind) :: &
i, j, iblk

integer (kind=int_kind) :: &
ni, nj, ncid, dimid, varid, ier

type (block) :: &
this_block ! block information for current block

integer (kind=int_kind) :: &
ilo,ihi,jlo,jhi ! beginning and end of physical domain

real (kind=dbl_kind) :: &
closelat, & ! Single-column latitude value
closelon, & ! Single-column longitude value
closelatidx, & ! Single-column latitude index to retrieve
closelonidx ! Single-column longitude index to retrieve

integer (kind=int_kind) :: &
start(2), & ! Start index to read in
count(2) ! Number of points to read in

integer (kind=int_kind) :: &
start3(3), & ! Start index to read in
count3(3) ! Number of points to read in

integer (kind=int_kind) :: &
status ! status flag

real (kind=dbl_kind), allocatable :: &
lats(:),lons(:),pos_lons(:), glob_grid(:,:) ! temporaries

real (kind=dbl_kind) :: &
pos_scmlon,& ! temporary
pi, &
puny, &
scamdata ! temporary

character(len=*), parameter :: subname = '(lonlatgrid)'

#ifdef USE_NETCDF
!-----------------------------------------------------------------
! - kmt file is actually clm fractional land file
! - Determine consistency of dimensions
! - Read in lon/lat centers in degrees from kmt file
! - Read in ocean from "kmt" file (1 for ocean, 0 for land)
!-----------------------------------------------------------------

call icepack_query_parameters(pi_out=pi, puny_out=puny)
call icepack_warnings_flush(nu_diag)
if (icepack_warnings_aborted()) call abort_ice(error_message=subname, &
file=__FILE__, line=__LINE__)

! Determine dimension of domain file and check for consistency

if (my_task == master_task) then
call ice_open_nc(kmt_file, ncid)

status = nf90_inq_dimid (ncid, 'ni', dimid)
status = nf90_inquire_dimension(ncid, dimid, len=ni)
status = nf90_inq_dimid (ncid, 'nj', dimid)
status = nf90_inquire_dimension(ncid, dimid, len=nj)
end if

! Determine start/count to read in for either single column or global lat-lon grid
! If single_column, then assume that only master_task is used since there is only one task

if (single_column) then
! Check for consistency
if (my_task == master_task) then
if ((nx_global /= 1).or. (ny_global /= 1)) then
write(nu_diag,*) 'Because you have selected the column model flag'
write(nu_diag,*) 'Please set nx_global=ny_global=1 in file'
write(nu_diag,*) 'ice_domain_size.F and recompile'
call abort_ice (subname//'ERROR: check nx_global, ny_global')
endif
end if

! Read in domain file for single column
allocate(lats(nj))
allocate(lons(ni))
allocate(pos_lons(ni))
allocate(glob_grid(ni,nj))

start3=(/1,1,1/)
count3=(/ni,nj,1/)
status = nf90_inq_varid(ncid, 'xc' , varid)
if (status /= nf90_noerr) call abort_ice (subname//' inq_varid xc')
status = nf90_get_var(ncid, varid, glob_grid, start3, count3)
if (status /= nf90_noerr) call abort_ice (subname//' get_var xc')
do i = 1,ni
lons(i) = glob_grid(i,1)
end do

status = nf90_inq_varid(ncid, 'yc' , varid)
if (status /= nf90_noerr) call abort_ice (subname//' inq_varid yc')
status = nf90_get_var(ncid, varid, glob_grid, start3, count3)
if (status /= nf90_noerr) call abort_ice (subname//' get_var yc')
do j = 1,nj
lats(j) = glob_grid(1,j)
end do

! convert lons array and scmlon to 0,360 and find index of value closest to 0
! and obtain single-column longitude/latitude indices to retrieve

pos_lons(:)= mod(lons(:) + 360._dbl_kind,360._dbl_kind)
pos_scmlon = mod(scmlon + 360._dbl_kind,360._dbl_kind)
start(1) = (MINLOC(abs(pos_lons-pos_scmlon),dim=1))
start(2) = (MINLOC(abs(lats -scmlat ),dim=1))

deallocate(lats)
deallocate(lons)
deallocate(pos_lons)
deallocate(glob_grid)

status = nf90_inq_varid(ncid, 'xc' , varid)
if (status /= nf90_noerr) call abort_ice (subname//' inq_varid xc')
status = nf90_get_var(ncid, varid, scamdata, start)
if (status /= nf90_noerr) call abort_ice (subname//' get_var xc')
TLON = scamdata
status = nf90_inq_varid(ncid, 'yc' , varid)
if (status /= nf90_noerr) call abort_ice (subname//' inq_varid yc')
status = nf90_get_var(ncid, varid, scamdata, start)
if (status /= nf90_noerr) call abort_ice (subname//' get_var yc')
TLAT = scamdata
status = nf90_inq_varid(ncid, 'area' , varid)
if (status /= nf90_noerr) call abort_ice (subname//' inq_varid area')
status = nf90_get_var(ncid, varid, scamdata, start)
if (status /= nf90_noerr) call abort_ice (subname//' get_var are')
tarea = scamdata
status = nf90_inq_varid(ncid, 'mask' , varid)
if (status /= nf90_noerr) call abort_ice (subname//' inq_varid mask')
status = nf90_get_var(ncid, varid, scamdata, start)
if (status /= nf90_noerr) call abort_ice (subname//' get_var mask')
hm = scamdata
status = nf90_inq_varid(ncid, 'frac' , varid)
if (status /= nf90_noerr) call abort_ice (subname//' inq_varid frac')
status = nf90_get_var(ncid, varid, scamdata, start)
if (status /= nf90_noerr) call abort_ice (subname//' get_var frac')
ocn_gridcell_frac = scamdata
else
! Check for consistency
if (my_task == master_task) then
if (nx_global /= ni .and. ny_global /= nj) then
write(nu_diag,*) 'latlongrid: ni,nj = ',ni,nj
write(nu_diag,*) 'latlongrid: nx_g,ny_g = ',nx_global, ny_global
call abort_ice (subname//'ERROR: ni,nj not equal to nx_global,ny_global')
end if
end if

! Read in domain file for global lat-lon grid
call ice_read_nc(ncid, 1, 'xc' , TLON , diag=.true.)
call ice_read_nc(ncid, 1, 'yc' , TLAT , diag=.true.)
call ice_read_nc(ncid, 1, 'area', tarea , diag=.true., &
field_loc=field_loc_center,field_type=field_type_scalar)
call ice_read_nc(ncid, 1, 'mask', hm , diag=.true.)
call ice_read_nc(ncid, 1, 'frac', ocn_gridcell_frac, diag=.true.)
end if

if (my_task == master_task) then
call ice_close_nc(ncid)
end if

!$OMP PARALLEL DO PRIVATE(iblk,this_block,ilo,ihi,jlo,jhi,i,j)
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi

do j = jlo, jhi
do i = ilo, ihi
! Convert from degrees to radians
TLON(i,j,iblk) = pi*TLON(i,j,iblk)/180._dbl_kind

! Convert from degrees to radians
TLAT(i,j,iblk) = pi*TLAT(i,j,iblk)/180._dbl_kind

! Convert from radians^2 to m^2
! (area in domain file is in radians^2 and tarea is in m^2)
tarea(i,j,iblk) = tarea(i,j,iblk) * (radius*radius)
end do
end do
end do
!$OMP END PARALLEL DO

!-----------------------------------------------------------------
! Calculate various geometric 2d arrays
! The U grid (velocity) is not used when run with sequential CAM
! because we only use thermodynamic sea ice. However, ULAT is used
! in the default initialization of CICE so we calculate it here as
! a "dummy" so that CICE will initialize with ice. If a no ice
! initialization is OK (or desired) this can be commented out and
! ULAT will remain 0 as specified above. ULAT is located at the
! NE corner of the grid cell, TLAT at the center, so here ULAT is
! hacked by adding half the latitudinal spacing (in radians) to
! TLAT.
!-----------------------------------------------------------------

!$OMP PARALLEL DO PRIVATE(iblk,this_block,ilo,ihi,jlo,jhi,i,j)
do iblk = 1, nblocks
this_block = get_block(blocks_ice(iblk),iblk)
ilo = this_block%ilo
ihi = this_block%ihi
jlo = this_block%jlo
jhi = this_block%jhi

do j = jlo, jhi
do i = ilo, ihi

if (ny_global == 1) then
uarea(i,j,iblk) = tarea(i,j, iblk)
else
uarea(i,j,iblk) = p25* &
(tarea(i,j, iblk) + tarea(i+1,j, iblk) &
+ tarea(i,j+1,iblk) + tarea(i+1,j+1,iblk))
endif
tarear(i,j,iblk) = c1/tarea(i,j,iblk)
uarear(i,j,iblk) = c1/uarea(i,j,iblk)
tinyarea(i,j,iblk) = puny*tarea(i,j,iblk)

if (single_column) then
ULAT (i,j,iblk) = TLAT(i,j,iblk)+(pi/nj)
else
if (ny_global == 1) then
ULAT (i,j,iblk) = TLAT(i,j,iblk)
else
ULAT (i,j,iblk) = TLAT(i,j,iblk)+(pi/ny_global)
endif
endif
ULON (i,j,iblk) = c0
ANGLE (i,j,iblk) = c0

ANGLET(i,j,iblk) = c0
HTN (i,j,iblk) = 1.e36_dbl_kind
HTE (i,j,iblk) = 1.e36_dbl_kind
dxt (i,j,iblk) = 1.e36_dbl_kind
dyt (i,j,iblk) = 1.e36_dbl_kind
dxu (i,j,iblk) = 1.e36_dbl_kind
dyu (i,j,iblk) = 1.e36_dbl_kind
dxhy (i,j,iblk) = 1.e36_dbl_kind
dyhx (i,j,iblk) = 1.e36_dbl_kind
cyp (i,j,iblk) = 1.e36_dbl_kind
cxp (i,j,iblk) = 1.e36_dbl_kind
cym (i,j,iblk) = 1.e36_dbl_kind
cxm (i,j,iblk) = 1.e36_dbl_kind
enddo
enddo
enddo
!$OMP END PARALLEL DO

call makemask
#else
call abort_ice(subname//'ERROR: USE_NETCDF cpp not defined', &
file=__FILE__, line=__LINE__)
#endif

end subroutine latlongrid
#endif

!=======================================================================

! Regular rectangular grid and mask
Expand Down
Loading