diff --git a/Inst_eddies/Tests/CMakeLists.txt b/Inst_eddies/Tests/CMakeLists.txt index 17a699e0de3856617d41b1de50cbe2d8e6909e09..3f04fb14beb6640794a34538d54326a5801e20e6 100644 --- a/Inst_eddies/Tests/CMakeLists.txt +++ b/Inst_eddies/Tests/CMakeLists.txt @@ -17,7 +17,10 @@ add_executable(test_set_all_outerm ${CMAKE_SOURCE_DIR}/derived_types.f90 set_all_outerm.f90 local_extrema.f90 get_1_outerm.f90 good_contour.f90 ${CMAKE_SOURCE_DIR}/spher_polyline_area.f90 nearby_extr.f90 - get_var.f90 ${CMAKE_CURRENT_LIST_DIR}/test_set_all_outerm.f90 config.f90) + get_var.f90 ${CMAKE_CURRENT_LIST_DIR}/test_set_all_outerm.f90 + config.f90 ${CMAKE_SOURCE_DIR}/shp_tr_create.f90 + ${CMAKE_SOURCE_DIR}/write_eddy.f90 + ${CMAKE_SOURCE_DIR}/shp_tr_close.f90) target_link_libraries(test_set_all_outerm geometry NetCDF95::netcdf95 numer_rec_95 shapelib_03 contour_531 jumble nr_util diff --git a/Inst_eddies/Tests/test_set_all_outerm.f90 b/Inst_eddies/Tests/test_set_all_outerm.f90 index 7227b95b2ab5e8fb0fbeb6818efaea85ff997c0c..326a92be908b78d9e4470fee22e80630b398919f 100644 --- a/Inst_eddies/Tests/test_set_all_outerm.f90 +++ b/Inst_eddies/Tests/test_set_all_outerm.f90 @@ -8,15 +8,14 @@ program test_set_all_outerm use netcdf95, only: nf95_open, nf95_close, nf95_get_var, find_coord, & nf95_inquire_dimension use nr_util, only: pi, assert, deg_to_rad, rad_to_deg, twopi - use shapelib, only: shpt_point, shpt_polygon, shpfileobject, ftdouble, & - shpclose, ftinteger - use shapelib_03, only: shp_create_03, dbf_add_field_03, shp_append_point_03, & - dbf_write_attribute_03, shp_append_null_03, shp_append_object_03 use config_m, only: config, max_radius_deg, min_radius - use derived_types, only: snapshot + use derived_types, only: snapshot, shp_tr, null_ssh_contour, missing_speed use get_var_m, only: get_var use set_all_outerm_m, only: set_all_outerm + use shp_tr_close_m, only: shp_tr_close + use shp_tr_create_m, only: shp_tr_create + use write_eddy_m, only: write_eddy implicit none @@ -35,11 +34,9 @@ program test_set_all_outerm ! (1 - max_radius(1):nlon + max_radius(1), nlat) if the grid is periodic ! in longitude, else (nlon, nlat). Sea-surface height, in m. - TYPE(shpfileobject) hshp_extremum ! shapefile extremum_$m - TYPE(shpfileobject) hshp_outermost ! shapefile outermost_contour_$m - integer ifield, ishape real step_deg(2), step(2) ! longitude and latitude steps, in degrees and rad logical periodic ! grid is periodic in longitude + TYPE(shp_tr) hshp !-------------------------------------------------------------- @@ -86,75 +83,22 @@ program test_set_all_outerm corner = [lon_min, lat_min] * deg_to_rad, & min_area = pi * (min_radius * 1e3)**2) - call shp_create_03("SHP_triplet/extremum", 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, 'valid', ftinteger, & - nwidth = 1, ndecimals = 0) - - call shp_create_03("SHP_triplet/outermost_contour", shpt_polygon, & - hshp_outermost) - call dbf_add_field_03(ifield, hshp_outermost, 'r_eq_area', ftdouble, & - nwidth = 10, ndecimals = 4) - 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) + do i = 1, s%number_vis_extr + s%list_vis(i)%speed_cont = null_ssh_contour() + s%list_vis(i)%max_speed = missing_speed + s%list_vis(i)%radius4 = 0 + end do + + call shp_tr_create(hshp, shp_tr_dir = "SHP_triplet") do i = 1, s%number_vis_extr - call shp_append_point_03(ishape, hshp_extremum, & - s%list_vis(i)%coord_extr * rad_to_deg) - call dbf_write_attribute_03(hshp_extremum, ishape, 0, & - s%list_vis(i)%ssh_extr) - call dbf_write_attribute_03(hshp_extremum, ishape, 1, 1) - call dbf_write_attribute_03(hshp_extremum, ishape, 2, i) - call dbf_write_attribute_03(hshp_extremum, ishape, 3, & - merge(1, 0, s%list_vis(i)%interpolated)) - call dbf_write_attribute_03(hshp_extremum, ishape, 4, & - merge(1, 0, s%list_vis(i)%cyclone)) - call dbf_write_attribute_03(hshp_extremum, ishape, 5, & - merge(1, 0, s%list_vis(i)%valid)) - - if (s%list_vis(i)%interpolated) then - call shp_append_null_03(ishape, hshp_outermost) - else - if (s%list_vis(i)%out_cont%n_points == 0) then - call shp_append_null_03(ishape, hshp_outermost) - else - call shp_append_object_03(ishape, hshp_outermost, shpt_polygon, & - s%list_vis(i)%out_cont%points * rad_to_deg) - end if - end if - - if (s%list_vis(i)%out_cont%area >= 0) then - call dbf_write_attribute_03(hshp_outermost, ishape, 0, & - sqrt(s%list_vis(i)%out_cont%area / 1e6 / pi)) - else - call dbf_write_attribute_03(hshp_outermost, ishape, 0, - 100.) - end if - - call dbf_write_attribute_03(hshp_outermost, ishape, 1, & - s%list_vis(i)%out_cont%ssh) - call dbf_write_attribute_03(hshp_outermost, ishape, 2, 1) - call dbf_write_attribute_03(hshp_outermost, ishape, 3, i) + call write_eddy(s%list_vis(i), hshp, k = 1, i = i) end do print *, "s%number_vis_extr = ", s%number_vis_extr - CALL shpclose(hshp_extremum) - print *, 'Created shapefile "SHP_triplet/extremum".' - CALL shpclose(hshp_outermost) - print *, 'Created shapefile "SHP_triplet/outermost_contour".' + CALL shp_tr_close(hshp) + print *, 'Created shapefiles in SHP_triplet.' print *, "Average number of points per outermost contour: ", & sum(s%list_vis%out_cont%n_points) / real(s%number_vis_extr)