Skip to content

Commit

Permalink
use grid_type from file; remove --grid_type arg from fill_fraction, f…
Browse files Browse the repository at this point in the history 
…ix_nonadvective, check_nonadvective; require B grid for fix_nonadvective, check_nonadvective
  • Loading branch information
@aekiss
aekiss committed Oct 24, 2024
1 parent 806f4df commit c176f26
Show file tree
Hide file tree
Showing 3 changed files with 46 additions and 52 deletions.
25 changes: 9 additions & 16 deletions README.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -73,30 +73,25 @@ Options

```
usage: topogtools fill_fraction --input <input_file> --output <output_file>
--fraction <frac> [--grid_type <type>]
--fraction <frac>
```

Cells with a fraction of sea area smaller than <frac> will have their depth set
to zero. Can produce non-advective cells and/or new seas.

Options
* `--grid_type <type>` Arakawa type of horizontal grid ('B' or 'C'; default is 'B')

### check_nonadvective

```
usage: topogtools check_nonadvective --input <input_file>
[--grid_type <type>
--vgrid <vgrid> --vgrid_type <type>
[--vgrid <vgrid> --vgrid_type <type>
--potholes --coastal-cells]
```

Check for non-advective cells. There are two types of checks available: potholes
and non-advective coastal cells. Checking for non-advective coastal cells should
only be needed when using a B-grid.
Check topography for non-advective cells. There are two types of checks
available: potholes and non-advective coastal cells. B-grid connectivity rules
are assumed. Aborts if input_file is not on a B-grid.

Options
* `--grid_type <type>` Arakawa type of horizontal grid ('B' or 'C'; default is 'B')
* `--vgrid <vgrid>` vertical grid (default 'ocean_vgrid.nc')
* `--vgrid_type <type>` can be mom5 or mom6 (default 'mom5')
* `--potholes` check for potholes
Expand All @@ -106,17 +101,15 @@ Options

```
usage: topogtools fix_nonadvective --input <input_file> --output <output_file>
[--grid_type <type>
--vgrid <vgrid> --vgrid_type <type>
[--vgrid <vgrid> --vgrid_type <type>
--potholes --coastal-cells]
```

Fix non-advective cells. There are two types of fixes available: potholes and
non-advective coastal cells. Fixes to non-advective coastal cells should only be
needed when using a B-grid.
Fix non-advective cells. There are two types of checks available: potholes and
non-advective coastal cells. B-grid connectivity rules are assumed. Aborts if
input_file is not on a B-grid.

Options
* `--grid_type <type>` Arakawa type of horizontal grid ('B' or 'C'; default is 'B')
* `--vgrid <vgrid>` vertical grid (default 'ocean_vgrid.nc')
* `--vgrid_type <type>` can be mom5 or mom6 (default 'mom5')
* `--potholes` fix potholes
Expand Down
30 changes: 19 additions & 11 deletions src/topography.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -49,18 +49,9 @@ type(topography_t) function topography_constructor(filename, grid_type) result(t
character(len=*), intent(in) :: filename
character(len=1), intent(in), optional :: grid_type

character(len=1):: file_grid_type
integer(int32) :: ncid, depth_id, frac_id, geolon_id, geolat_id, dids(2), history_len ! NetCDF ids

if (present(grid_type)) then
if ( grid_type == 'B' .or. grid_type == 'C' ) then
topog%grid_type = grid_type
else
call handle_error(nf90_einval, .true., "grid_type must be B or C")
end if
else
topog%grid_type = 'B'
end if

write(output_unit,'(3a)') "Reading topography from file '", trim(filename), "'"

topog%original_file = filename
Expand All @@ -84,6 +75,19 @@ type(topography_t) function topography_constructor(filename, grid_type) result(t
call handle_error(nf90_get_att(ncid, depth_id, 'maximum_depth', topog%max_depth), isfatal=.false.)
call handle_error(nf90_get_att(ncid, depth_id, 'nonadvective_cells_removed', topog%nonadvective_cells_removed), isfatal=.false.)

if (present(grid_type)) then
topog%grid_type = grid_type ! grid_type arg overrides value in file
else
call handle_error(nf90_get_att(ncid, depth_id, 'grid_type', file_grid_type), isfatal=.false.)
if ( file_grid_type == 'B' .or. file_grid_type == 'C' ) then
topog%grid_type = file_grid_type
end if
end if
write(output_unit,*) " grid_type = ", topog%grid_type
if (.not.( topog%grid_type == 'B' .or. topog%grid_type == 'C' )) then
call handle_error(nf90_einval, .true., "grid_type must be B or C")
end if

! Get sea area fraction
call handle_error(nf90_inq_varid(ncid, 'sea_area_fraction', frac_id))
allocate(topog%frac(topog%nxt, topog%nyt))
Expand Down Expand Up @@ -541,6 +545,10 @@ subroutine topography_nonadvective(this, vgrid_file, vgrid_type, potholes, coast
integer(int32) :: im, ip, jm, jp
integer(int32) :: nseas

if ( this%grid_type /= 'B' ) then
call handle_error(nf90_einval, .true., "nonadvective: grid_type must be B")
end if

vgrid = vgrid_t(vgrid_file, vgrid_type)
write(output_unit,*) 'Zeta dimensions', 2*vgrid%nlevels + 1, vgrid%nlevels
allocate(zw(0:vgrid%nlevels))
Expand Down Expand Up @@ -651,7 +659,7 @@ subroutine topography_nonadvective(this, vgrid_file, vgrid_type, potholes, coast
if (fix .and. (coastal .or. potholes)) then
this%nonadvective_cells_removed = 'yes'
if (changes_made .and. this%lakes_removed == 'yes') then
! Check if new lakes were created new lakes
! Check if new lakes were created
call this%number_seas(number_of_seas = nseas, silent=.true.)
if (nseas > 1) then
write(output_unit,'(a)') "WARNING: new seas have been created. To fix, rerun deseas again."
Expand Down
43 changes: 18 additions & 25 deletions src/topogtools.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,7 @@ program topogtools
'usage: topogtools [--help] <command> [<args>] ', &
' ', &
'Collection of tools to edit and manipulate ocean model topographies. ', &
'See ''topogtools --help <command>'' to read about a specific subcommand. ', &
'See ''topogtools --help <command>'' to read about a specific subcommand. ', &
' ', &
'Available commands: ', &
' gen_topo - Generate a new topography file from a bathymetry ', &
Expand All @@ -46,7 +46,7 @@ program topogtools
' tripolar grid ', &
' --longitude-offset <value> offset (in degrees) between the central ', &
' longitude of the ocean horizontal grid and of ', &
' the bathymetry grid (default ''0.0'') ', &
' the bathymetry grid (default ''0.0'') ', &
'']
help_deseas = [character(len=80) :: &
'usage: topogtools deseas --input <input_file> --output <output_file> ', &
Expand All @@ -67,52 +67,45 @@ program topogtools
'Can produce non-advective cells. ', &
' ', &
'Options ', &
' --vgrid <vgrid> vertical grid (default ''ocean_vgrid.nc'') ', &
' --vgrid_type <type> can be ''mom5'' or ''mom6'' (default ''mom5'') ', &
' --vgrid <vgrid> vertical grid (default ''ocean_vgrid.nc'') ', &
' --vgrid_type <type> can be ''mom5'' or ''mom6'' (default ''mom5'') ', &
'']
help_fill_fraction = [character(len=80) :: &
'usage: topogtools fill_fraction --input <input_file> --output <output_file> ', &
' --fraction <frac> [--grid_type <type>] ', &
' --fraction <frac> ', &
' ', &
'Cells with a fraction of sea area smaller than <frac> will have their depth set ', &
'to zero. Can produce non-advective cells and/or new seas. ', &
' ', &
'Options ', &
' --grid_type <type> Arakawa type of horizontal grid (''B'' or ''C''; default ''B'')', &
'']
help_fix_nonadvective = [character(len=80) :: &
'usage: topogtools fix_nonadvective --input <input_file> --output <output_file> ', &
' [--grid_type <type> ', &
' --vgrid <vgrid> --vgrid_type <type> ', &
' [--vgrid <vgrid> --vgrid_type <type> ', &
' --potholes --coastal_cells] ', &
' ', &
'Fix non-advective cells. There are two types of fixes available: potholes and ', &
'non-advective coastal cells. Fixes to non-advective coastal cells should only be', &
'needed when using a B-grid. ', &
' ', &
'Options ', &
' --grid_type <type> Arakawa type of horizontal grid (''B'' or ''C''; default ''B'')', &
' --vgrid <vgrid> vertical grid (default ''ocean_vgrid.nc'') ', &
' --vgrid_type <type> can be ''mom5'' or ''mom6'' (default ''mom5'') ', &
' --vgrid <vgrid> vertical grid (default ''ocean_vgrid.nc'') ', &
' --vgrid_type <type> can be ''mom5'' or ''mom6'' (default ''mom5'') ', &
' --potholes fix potholes ', &
' --coastal-cells fix non-advective coastal cells ', &
'']
help_check_nonadvective = [character(len=80) :: &
'usage: topogtools check_nonadvective --input <input_file> ', &
' [--grid_type <type> ', &
' --vgrid <vgrid> --vgrid_type <type> ', &
' [--vgrid <vgrid> --vgrid_type <type> ', &
' --potholes --coastal_cells] ', &
' ', &
'Check for non-advective cells. There are two types of checks available: potholes', &
'and non-advective coastal cells. Checking for non-advective coastal cells should', &
'only be needed when using a B-grid. ', &
' ', &
'Options ', &
' --grid_type <type> Arakawa type of horizontal grid (''B'' or ''C''; default ''B'')', &
' --vgrid <vgrid> vertical grid (default ''ocean_vgrid.nc'') ', &
' --vgrid_type <type> can be ''mom5'' or ''mom6'' (default ''mom5'') ', &
' --vgrid <vgrid> vertical grid (default ''ocean_vgrid.nc'') ', &
' --vgrid_type <type> can be ''mom5'' or ''mom6'' (default ''mom5'') ', &
' --potholes check for potholes ', &
' --coastal-cells check for non-advective coastal cells ', &
' --coastal-cells check for non-advective coastal cells (for B grid) ', &
'']
help_mask = [character(len=80) :: &
'usage: topogtools mask --input <input_file> --output <output_file> ', &
Expand All @@ -132,12 +125,12 @@ program topogtools
call set_args('--input:i "unset" --output:o "unset" --vgrid "ocean_vgrid.nc" --vgrid_type "mom5" --level 0', &
help_min_max_depth, version_text)
case('fill_fraction')
call set_args('--input:i "unset" --output:o "unset" --fraction 0.0 --grid_type "B"', help_fill_fraction, version_text)
call set_args('--input:i "unset" --output:o "unset" --fraction 0.0', help_fill_fraction, version_text)
case ('fix_nonadvective')
call set_args('--input:i "unset" --output:o "unset" --grid_type "B" --vgrid "ocean_vgrid.nc" --vgrid_type "mom5" --potholes F &
call set_args('--input:i "unset" --output:o "unset" --vgrid "ocean_vgrid.nc" --vgrid_type "mom5" --potholes F &
&--coastal-cells F', help_fix_nonadvective, version_text)
case ('check_nonadvective')
call set_args('--input:i "unset" --grid_type "B" --vgrid "ocean_vgrid.nc" --vgrid_type "mom5" --potholes F --coastal-cells F', &
call set_args('--input:i "unset" --vgrid "ocean_vgrid.nc" --vgrid_type "mom5" --potholes F --coastal-cells F', &
help_check_nonadvective, version_text)
case ('mask')
call set_args('--input:i "unset" --output:o "unset"', help_mask, version_text)
Expand Down Expand Up @@ -202,19 +195,19 @@ program topogtools
write(error_unit,'(a)') "ERROR: sea area fraction must be larger than 0 and smaller than 1"
error stop
end if
topog = topography_t(file_in, grid_type=sget('grid_type'))
topog = topography_t(file_in)
call topog%fill_fraction(sea_area_fraction)
call topog%update_history(get_mycommand())
call topog%write(file_out)

case ('fix_nonadvective')
topog = topography_t(file_in, grid_type=sget('grid_type'))
topog = topography_t(file_in)
call topog%nonadvective(vgrid, sget('vgrid_type'), lget('potholes'), lget('coastal-cells'), fix=.true.)
call topog%update_history(get_mycommand())
call topog%write(file_out)

case ('check_nonadvective')
topog = topography_t(file_in, grid_type=sget('grid_type'))
topog = topography_t(file_in)
call topog%nonadvective(vgrid, sget('vgrid_type'), lget('potholes'), lget('coastal-cells'), fix=.false.)

case ('mask')
Expand Down

0 comments on commit c176f26

Please sign in to comment.