Skip to content
Snippets Groups Projects
Commit bf045ac7 authored by Lionel GUEZ's avatar Lionel GUEZ
Browse files

Move initialization of shapefiles to a separate procedure,

init_shapefiles. This allows us to have variables for DBF field
numbers, as module variables. We avoid the burden of passing this
large number of variables as arguments to write_eddy. Safer to have
variables for DBF field numbers than relying on the order of those
fields.

Remove file_list. This was a duplication of information.

Use named constant for missing max speed value. It is used in two procedures.
parent a2da4d88
No related branches found
No related tags found
No related merge requests found
......@@ -11,13 +11,13 @@ src_test_get_1_outerm = good_contour.f test_get_1_outerm.f derived_types.f get_1
src_test_set_max_speed = test_set_max_speed.f derived_types.f set_max_speed.f good_contour.f max_speed_contour_ssh.f mean_speed.f spherical_polygon_area.f inside_4.f
src_test_get_snapshot = test_get_snapshot.f get_snapshot.f dispatch_snapshot.f write_eddy.f send_snapshot.f receive_snapshot.f local_extrema.f set_max_speed.f outermost_possible_level.f get_1_outerm.f max_speed_contour_ssh.f good_contour.f spherical_polygon_area.f mean_speed.f inside_4.f set_all_outerm.f derived_types.f
src_test_get_snapshot = test_get_snapshot.f get_snapshot.f dispatch_snapshot.f write_eddy.f send_snapshot.f receive_snapshot.f local_extrema.f set_max_speed.f outermost_possible_level.f get_1_outerm.f max_speed_contour_ssh.f good_contour.f spherical_polygon_area.f mean_speed.f inside_4.f set_all_outerm.f derived_types.f init_shapefiles.f
src_test_set_all_outerm = test_set_all_outerm.f derived_types.f set_all_outerm.f local_extrema.f get_1_outerm.f good_contour.f spherical_polygon_area.f
src_test_weight = test_weight.f weight.f derived_types.f
sources := $(shell cat ${makefile_dir}/file_list)
sources := $(sort ${src_test_local_extrema} ${src_test_get_1_outerm} ${src_test_set_max_speed} ${src_test_get_snapshot} ${src_test_set_all_outerm} ${src_test_weight}) test_good_contour.f test_inside_4.f test_max_speed_contour_ssh.f test_mean_speed.f
lib_list = contour_531 numer_rec_95 shapelib_03 netcdf95 geometry jumble netcdff fortrangis shp fortranc nr_util
......
......@@ -6,34 +6,24 @@ program test_get_snapshot
use derived_types, only: snapshot
use dispatch_snapshot_m, only: dispatch_snapshot
use get_snapshot_m, only: get_snapshot
use init_shapefiles_m, only: init_shapefiles
use jumble, only: new_unit
use netcdf, only: nf90_nowrite
use netcdf95, only: nf95_open, find_coord, nf95_inquire_dimension, &
nf95_get_var, nf95_close
use nr_util, only: pi, assert
use shapelib, only: shpt_point, shpt_polygon, shpfileobject, ftdouble, &
shpclose, ftinteger
use shapelib_03, only: shp_create_03, dbf_add_field_03
use shapelib, only: shpfileobject, shpclose
implicit none
type(snapshot) s
real, parameter:: deg_over_rad = pi / 180.
TYPE(shpfileobject) hshp_extremum
! shapefile extremum_$m. The fields in the DBF file are, in that
! order: ssh, date index, eddy index, interpolated, cyclone,
! sufficient amplitude, speed.
TYPE(shpfileobject) hshp_extremum ! shapefile extremum_$m
TYPE(shpfileobject) hshp_outermost ! shapefile outermost_contour_$m
TYPE(shpfileobject) hshp_max_speed ! shapefile x_speed_contour_$m
TYPE(shpfileobject) hshp_outermost
! shapefile outermost_contour_$m. The fields in the DBF file are,
! in that order: area, ssh, date index, eddy index, twice, radius4.
TYPE(shpfileobject) hshp_max_speed
! shapefile x_speed_contour_$m. The fields in the DBF file are, in
! that order: area, ssh, date index, eddy index.
integer ifield, unit_isolated, unit_number_eddies
integer unit_isolated, unit_number_eddies
real:: min_amp = 0.
! minimum amplitude of ssh, between extremum and outermost contour,
......@@ -71,46 +61,8 @@ program test_get_snapshot
step = &
[longitude(2) - longitude(1), latitude(2) - latitude(1)] * deg_over_rad)
call shp_create_03("extremum_1", shpt_point, hshp_extremum)
call dbf_add_field_03(ifield, hshp_extremum, 'ssh', ftdouble, nwidth = 13, &
ndecimals = 6)
call dbf_add_field_03(ifield, hshp_extremum, 'date_index', ftinteger, &
nwidth = 4, ndecimals = 0)
call dbf_add_field_03(ifield, hshp_extremum, 'eddy_index', ftinteger, &
nwidth = 5, ndecimals = 0)
call dbf_add_field_03(ifield, hshp_extremum, 'interpolat', ftinteger, &
nwidth = 1, ndecimals = 0)
call dbf_add_field_03(ifield, hshp_extremum, 'cyclone', ftinteger, &
nwidth = 1, ndecimals = 0)
call dbf_add_field_03(ifield, hshp_extremum, 'suff_amp', ftinteger, &
nwidth = 1, ndecimals = 0)
call dbf_add_field_03(ifield, hshp_extremum, 'speed', ftdouble, nwidth = 13, &
ndecimals = 6)
call shp_create_03("outermost_contour_1", shpt_polygon, hshp_outermost)
call dbf_add_field_03(ifield, hshp_outermost, 'area', ftdouble, nwidth = 20, &
ndecimals = 6)
call dbf_add_field_03(ifield, hshp_outermost, 'ssh', ftdouble, nwidth = 13, &
ndecimals = 6)
call dbf_add_field_03(ifield, hshp_outermost, 'date_index', ftinteger, &
nwidth = 4, ndecimals = 0)
call dbf_add_field_03(ifield, hshp_outermost, 'eddy_index', ftinteger, &
nwidth = 5, ndecimals = 0)
call dbf_add_field_03(ifield, hshp_outermost, 'twice', ftinteger, &
nwidth = 1, ndecimals = 0)
call dbf_add_field_03(ifield, hshp_outermost, 'radius4', ftinteger, &
nwidth = 2, ndecimals = 0)
call shp_create_03("max_speed_contour_1", shpt_polygon, hshp_max_speed)
call dbf_add_field_03(ifield, hshp_max_speed, 'area', ftdouble, nwidth = 20, &
ndecimals = 6)
call dbf_add_field_03(ifield, hshp_max_speed, 'ssh', ftdouble, nwidth = 13, &
ndecimals = 6)
call dbf_add_field_03(ifield, hshp_max_speed, 'date_index', ftinteger, &
nwidth = 4, ndecimals = 0)
call dbf_add_field_03(ifield, hshp_max_speed, 'eddy_index', ftinteger, &
nwidth = 5, ndecimals = 0)
call init_shapefiles(hshp_extremum, hshp_outermost, hshp_max_speed)
call new_unit(unit_isolated)
open(unit_isolated, file = "isolated_nodes_1.csv", status = "replace", &
action = "write")
......
dispatch_snapshot.o : write_eddy.o send_snapshot.o derived_types.o
get_1_outerm.o : spherical_polygon_area.o outermost_possible_level.o good_contour.o derived_types.o
get_snapshot.o : set_all_outerm.o set_max_speed.o receive_snapshot.o derived_types.o
receive_snapshot.o : derived_types.o
send_snapshot.o : derived_types.o
set_all_outerm.o : local_extrema.o get_1_outerm.o derived_types.o
set_max_speed.o : spherical_polygon_area.o mean_speed.o max_speed_contour_ssh.o inside_4.o good_contour.o derived_types.o
get_1_outerm.o : spherical_polygon_area.o outermost_possible_level.o good_contour.o derived_types.o
test_get_snapshot.o : get_snapshot.o dispatch_snapshot.o derived_types.o
test_good_contour.o : good_contour.o
test_inside_4.o : inside_4.o
test_get_1_outerm.o : get_1_outerm.o derived_types.o
test_get_snapshot.o : init_shapefiles.o get_snapshot.o dispatch_snapshot.o derived_types.o
test_local_extrema.o : local_extrema.o
test_max_speed_contour_ssh.o : max_speed_contour_ssh.o
test_mean_speed.o : mean_speed.o
test_set_all_outerm.o : set_all_outerm.o derived_types.o
test_set_max_speed.o : set_max_speed.o derived_types.o
test_get_1_outerm.o : get_1_outerm.o derived_types.o
write_eddy.o : derived_types.o
weight.o : derived_types.o
test_weight.o : weight.o derived_types.o
weight.o : derived_types.o
write_eddy.o : init_shapefiles.o derived_types.o
test_good_contour.o : good_contour.o
test_inside_4.o : inside_4.o
test_max_speed_contour_ssh.o : max_speed_contour_ssh.o
test_mean_speed.o : mean_speed.o
......@@ -59,6 +59,7 @@ module derived_types
end type snapshot
real, parameter:: missing_ssh = 1e4 ! flag for undefined contour
real, parameter:: missing_speed = 1e4 ! flag for missing max speed
contains
......
......@@ -15,19 +15,13 @@ contains
type(snapshot), intent(in):: s
TYPE(shpfileobject), intent(inout):: hshp_extremum
! shapefile extremum_$m. We assume that the fields in the DBF file
! are, in that order: ssh, date index, eddy index, interpolated,
! cyclone, sufficient amplitude, speed.
TYPE(shpfileobject), intent(inout):: hshp_extremum ! shapefile extremum_$m
TYPE(shpfileobject), intent(inout):: hshp_outermost
! shapefile outermost_contour_$m. We assume that the fields in the
! DBF file are, in that order: area, ssh, date index, eddy index,
! twice, radius4.
! shapefile outermost_contour_$m
TYPE(shpfileobject), intent(inout):: hshp_max_speed
! shapefile x_speed_contour_$m. We assume that the fields in the
! DBF file are, in that order: area, ssh, date index, eddy index.
! shapefile x_speed_contour_$m
integer, intent(in):: unit_isolated
! logical unit for file isolated_nodes_$m.csv
......
derived_types.f
dispatch_snapshot.f
get_snapshot.f
good_contour.f
inside_4.f
local_extrema.f
max_speed_contour_ssh.f
mean_speed.f
outermost_possible_level.f
receive_snapshot.f
send_snapshot.f
set_all_outerm.f
set_max_speed.f
get_1_outerm.f
spherical_polygon_area.f
test_get_snapshot.f
test_good_contour.f
test_inside_4.f
test_local_extrema.f
test_max_speed_contour_ssh.f
test_mean_speed.f
test_set_all_outerm.f
test_set_max_speed.f
test_get_1_outerm.f
write_eddy.f
weight.f
test_weight.f
......@@ -10,7 +10,7 @@ contains
use, intrinsic:: ieee_arithmetic, only: ieee_value, IEEE_QUIET_NAN
use contour_531, only: convert_to_ind, null_polyline
use derived_types, only: snapshot, null_ssh_contour
use derived_types, only: snapshot, null_ssh_contour, missing_speed
use netcdf, only: nf90_nowrite
use netcdf95, only: nf95_open, nf95_close, nf95_inq_varid, nf95_get_var, &
nf95_get_att
......@@ -103,7 +103,7 @@ contains
step)
else
s%list_vis(i)%speed_cont = null_ssh_contour()
s%list_vis(i)%max_speed = 1e4
s%list_vis(i)%max_speed = missing_speed
s%list_vis(i)%radius4 = 0
end if
end do
......
module init_shapefiles_m
implicit none
integer, protected:: ifield_extr_ssh, ifield_extr_date, &
ifield_extr_eddy_index, ifield_extr_interp, ifield_extr_cycl, &
ifield_extr_suff_amp, ifield_extr_speed
integer, protected:: ifield_out_area, ifield_out_ssh, ifield_out_date, &
ifield_out_eddy_index, ifield_out_twice, ifield_out_radius4
integer, protected:: ifield_max_speed_area, ifield_max_speed_ssh, &
ifield_max_speed_date, ifield_max_speed_eddy_index
contains
subroutine init_shapefiles(hshp_extremum, hshp_outermost, hshp_max_speed)
use shapelib, only: shpt_point, shpt_polygon, shpfileobject, ftdouble, &
ftinteger
use shapelib_03, only: shp_create_03, dbf_add_field_03
TYPE(shpfileobject), intent(out):: hshp_extremum ! shapefile extremum_$m
TYPE(shpfileobject), intent(out):: hshp_outermost
! shapefile outermost_contour_$m
TYPE(shpfileobject), intent(out):: hshp_max_speed
! shapefile x_speed_contour_$m
!---------------------------------------------------------------------
call shp_create_03("extremum_1", shpt_point, hshp_extremum)
call dbf_add_field_03(ifield_extr_ssh, hshp_extremum, 'ssh', ftdouble, &
nwidth = 13, ndecimals = 6)
call dbf_add_field_03(ifield_extr_date, hshp_extremum, 'date_index', &
ftinteger, nwidth = 4, ndecimals = 0)
call dbf_add_field_03(ifield_extr_eddy_index, hshp_extremum, 'eddy_index', &
ftinteger, nwidth = 5, ndecimals = 0)
call dbf_add_field_03(ifield_extr_interp, hshp_extremum, 'interpolat', &
ftinteger, nwidth = 1, ndecimals = 0)
call dbf_add_field_03(ifield_extr_cycl, hshp_extremum, 'cyclone', &
ftinteger, nwidth = 1, ndecimals = 0)
call dbf_add_field_03(ifield_extr_suff_amp, hshp_extremum, 'suff_amp', &
ftinteger, nwidth = 1, ndecimals = 0)
call dbf_add_field_03(ifield_extr_speed, hshp_extremum, 'speed', ftdouble, &
nwidth = 13, ndecimals = 6)
call shp_create_03("outermost_contour_1", shpt_polygon, hshp_outermost)
call dbf_add_field_03(ifield_out_area, hshp_outermost, 'area', ftdouble, &
nwidth = 20, ndecimals = 6)
call dbf_add_field_03(ifield_out_ssh, hshp_outermost, 'ssh', ftdouble, &
nwidth = 13, ndecimals = 6)
call dbf_add_field_03(ifield_out_date, hshp_outermost, 'date_index', &
ftinteger, nwidth = 4, ndecimals = 0)
call dbf_add_field_03(ifield_out_eddy_index, hshp_outermost, 'eddy_index', &
ftinteger, nwidth = 5, ndecimals = 0)
call dbf_add_field_03(ifield_out_twice, hshp_outermost, 'twice', &
ftinteger, nwidth = 1, ndecimals = 0)
call dbf_add_field_03(ifield_out_radius4, hshp_outermost, 'radius4', &
ftinteger, nwidth = 2, ndecimals = 0)
call shp_create_03("max_speed_contour_1", shpt_polygon, hshp_max_speed)
call dbf_add_field_03(ifield_max_speed_area, hshp_max_speed, 'area', &
ftdouble, nwidth = 20, ndecimals = 6)
call dbf_add_field_03(ifield_max_speed_ssh, hshp_max_speed, 'ssh', &
ftdouble, nwidth = 13, ndecimals = 6)
call dbf_add_field_03(ifield_max_speed_date, hshp_max_speed, 'date_index', &
ftinteger, nwidth = 4, ndecimals = 0)
call dbf_add_field_03(ifield_max_speed_eddy_index, hshp_max_speed, &
'eddy_index', ftinteger, nwidth = 5, ndecimals = 0)
end subroutine init_shapefiles
end module init_shapefiles_m
......@@ -7,53 +7,59 @@ contains
subroutine write_eddy(e, k, i, hshp_extremum, hshp_outermost, hshp_max_speed)
use, intrinsic:: ieee_arithmetic, only: ieee_is_nan
use derived_types, only: eddy
use derived_types, only: eddy, missing_speed
use init_shapefiles_m, only: ifield_extr_ssh, ifield_extr_date, &
ifield_extr_eddy_index, ifield_extr_interp, ifield_extr_cycl, &
ifield_extr_suff_amp, ifield_extr_speed, ifield_out_area, &
ifield_out_ssh, ifield_out_date, ifield_out_eddy_index, &
ifield_out_twice, ifield_out_radius4, ifield_max_speed_area, &
ifield_max_speed_ssh, ifield_max_speed_date, &
ifield_max_speed_eddy_index
use nr_util, only: pi
use shapelib, only: shpfileobject, shpt_polygon
use shapelib_03, only: shp_append_point_03, dbf_write_attribute_03, &
shp_append_object_03, shp_append_null_03
type(eddy), intent(in):: e
integer, intent(in):: k ! date index
integer, intent(in):: i ! eddy index
TYPE(shpfileobject), intent(inout):: hshp_extremum
! shapefile extremum_$m. We assume that the fields in the DBF file
! are, in that order: ssh, date index, eddy index, interpolated,
! cyclone, sufficient amplitude, speed.
TYPE(shpfileobject), intent(inout):: hshp_extremum ! shapefile extremum_$m
TYPE(shpfileobject), intent(inout):: hshp_outermost
! shapefile outermost_contour_$m. We assume that the fields in the
! DBF file are, in that order: area, ssh, date index, eddy index,
! twice, radius4.
! shapefile outermost_contour_$m
TYPE(shpfileobject), intent(inout):: hshp_max_speed
! shapefile x_speed_contour_$m. We assume that the fields in the
! DBF file are, in that order: area, ssh, date index, eddy index.
! shapefile x_speed_contour_$m
! Local:
integer ishape
real, parameter:: rad_over_deg = 180. / pi
!-------------------------------------------------------------
call shp_append_point_03(ishape, hshp_extremum, e%coord_extr * rad_over_deg)
call dbf_write_attribute_03(hshp_extremum, ishape, 0, e%ssh_extr)
call dbf_write_attribute_03(hshp_extremum, ishape, 1, k)
call dbf_write_attribute_03(hshp_extremum, ishape, 2, i)
call dbf_write_attribute_03(hshp_extremum, ishape, 3, &
call dbf_write_attribute_03(hshp_extremum, ishape, ifield_extr_ssh, &
e%ssh_extr)
call dbf_write_attribute_03(hshp_extremum, ishape, ifield_extr_date, k)
call dbf_write_attribute_03(hshp_extremum, ishape, ifield_extr_eddy_index, &
i)
call dbf_write_attribute_03(hshp_extremum, ishape, ifield_extr_interp, &
merge(1, 0, e%interpolated))
call dbf_write_attribute_03(hshp_extremum, ishape, 4, &
call dbf_write_attribute_03(hshp_extremum, ishape, ifield_extr_cycl, &
merge(1, 0, e%cyclone))
call dbf_write_attribute_03(hshp_extremum, ishape, 5, &
call dbf_write_attribute_03(hshp_extremum, ishape, ifield_extr_suff_amp, &
merge(1, 0, e%suff_amp))
if (ieee_is_nan(e%max_speed)) then
call dbf_write_attribute_03(hshp_extremum, ishape, 6, 1e4)
call dbf_write_attribute_03(hshp_extremum, ishape, ifield_extr_speed, &
missing_speed)
! (Cannot write NaN to dbf file.)
else
call dbf_write_attribute_03(hshp_extremum, ishape, 6, e%max_speed)
call dbf_write_attribute_03(hshp_extremum, ishape, ifield_extr_speed, &
e%max_speed)
end if
if (e%interpolated) then
......@@ -75,17 +81,26 @@ contains
end if
end if
call dbf_write_attribute_03(hshp_outermost, ishape, 0, e%out_cont%area)
call dbf_write_attribute_03(hshp_outermost, ishape, 1, e%out_cont%ssh)
call dbf_write_attribute_03(hshp_outermost, ishape, 2, k)
call dbf_write_attribute_03(hshp_outermost, ishape, 3, i)
call dbf_write_attribute_03(hshp_outermost, ishape, 4, merge(1, 0, e%twice))
call dbf_write_attribute_03(hshp_outermost, ishape, 5, e%radius4)
call dbf_write_attribute_03(hshp_max_speed, ishape, 0, e%speed_cont%area)
call dbf_write_attribute_03(hshp_max_speed, ishape, 1, e%speed_cont%ssh)
call dbf_write_attribute_03(hshp_max_speed, ishape, 2, k)
call dbf_write_attribute_03(hshp_max_speed, ishape, 3, i)
call dbf_write_attribute_03(hshp_outermost, ishape, ifield_out_area, &
e%out_cont%area)
call dbf_write_attribute_03(hshp_outermost, ishape, ifield_out_ssh, &
e%out_cont%ssh)
call dbf_write_attribute_03(hshp_outermost, ishape, ifield_out_date, k)
call dbf_write_attribute_03(hshp_outermost, ishape, ifield_out_eddy_index, &
i)
call dbf_write_attribute_03(hshp_outermost, ishape, ifield_out_twice, &
merge(1, 0, e%twice))
call dbf_write_attribute_03(hshp_outermost, ishape, ifield_out_radius4, &
e%radius4)
call dbf_write_attribute_03(hshp_max_speed, ishape, ifield_max_speed_area, &
e%speed_cont%area)
call dbf_write_attribute_03(hshp_max_speed, ishape, ifield_max_speed_ssh, &
e%speed_cont%ssh)
call dbf_write_attribute_03(hshp_max_speed, ishape, ifield_max_speed_date, &
k)
call dbf_write_attribute_03(hshp_max_speed, ishape, &
ifield_max_speed_eddy_index, i)
end subroutine write_eddy
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment