diff --git a/Overlap/CMakeLists.txt b/Overlap/CMakeLists.txt index 4751624a6c493624263c762ff9caff938b54f46d..9ca60079b328c4423bf00fc0d46ff2238727cb01 100644 --- a/Overlap/CMakeLists.txt +++ b/Overlap/CMakeLists.txt @@ -7,7 +7,7 @@ if(MPI_Fortran_HAVE_F08_MODULE) add_executable(eddy_graph eddy_graph.f90 get_snapshot.f90 overlap.f90 dispatch_snapshot.f90 recv_snapshot.f90 candidate_overlap.f90 spher_polygon_area.f90 send_snapshot.f90 unit_edge_m.f90 config_graph.F90 - read_grid.F90 read_snapshot.f90 read_eddy.f90) + read_grid.F90 read_snapshot.f90 read_eddy.f90 inters_lines.f90) target_link_libraries(eddy_graph PRIVATE EZMPI::ezmpi Jumble::jumble Shapelib_03::shapelib_03 Contour_531::contour_531 gpc_f Geometry::geometry MPI::MPI_Fortran Numer_Rec_95::numer_rec_95) @@ -21,7 +21,7 @@ endif() add_subdirectory(Tests) target_sources(test_overlap PRIVATE overlap.f90 spher_polygon_area.f90 candidate_overlap.f90 unit_edge_m.f90 config_graph.F90 read_grid.F90 - read_snapshot.f90 read_eddy.f90) + read_snapshot.f90 read_eddy.f90 inters_lines.f90) target_sources(test_spher_polygon_area PRIVATE spher_polygon_area.f90) target_sources(test_read_snapshot PRIVATE read_snapshot.f90 read_eddy.f90 read_grid.F90 config_graph.F90) diff --git a/Overlap/inters_lines.f90 b/Overlap/inters_lines.f90 new file mode 100644 index 0000000000000000000000000000000000000000..cf1ad92343c4d63247069e581e9db05efc7d240d --- /dev/null +++ b/Overlap/inters_lines.f90 @@ -0,0 +1,43 @@ +module inters_lines_m + + implicit none + +contains + + logical function inters_lines(polyline_1, polyline_2, min_inters) + + ! Do the two polylines intersect with a sufficient area of + ! intersection? + + ! Libraries: + use contour_531, only: polyline + use gpc_f, only: gpc_polygon_clip_f, GPC_INT, polygon + + use spher_polygon_area_m, only: spher_polygon_area + use spher_polyline_area_m, only: spher_polyline_area + + type(polyline), intent(in):: polyline_1, polyline_2 + real, intent(in):: min_inters ! minimum intersection between contours + + ! Local: + type(polygon) res_pol + + !------------------------------------------------------------------------- + + call gpc_polygon_clip_f(GPC_INT, polygon(nparts = 1, & + part = [polyline_1], hole = [.false.]), polygon(nparts = 1, & + part = [polyline_2], hole = [.false.]), res_pol) + + if (res_pol%nparts == 0) then + ! polyline_1 and polyline_2 do not overlap + inters_lines = .false. + else + ! polyline_1 and polyline_2 overlap + inters_lines = spher_polygon_area(res_pol) >= min_inters & + * min(abs(spher_polyline_area(polyline_1)), & + abs(spher_polyline_area(polyline_2))) + end if + + end function inters_lines + +end module inters_lines_m diff --git a/Overlap/overlap.f90 b/Overlap/overlap.f90 index 9bf4af49ca4dc929b2d95448b1b966e8b0b6aebd..56b003162add94959ee0a01066a65ac9dd28f260 100644 --- a/Overlap/overlap.f90 +++ b/Overlap/overlap.f90 @@ -13,15 +13,13 @@ contains ! Libraries: use contour_531, only: polyline - use gpc_f, only: gpc_polygon_clip_f, GPC_INT, polygon use jumble, only: twopi + use inters_lines_m, only: inters_lines use candidate_overlap_m, only: candidate_overlap use config_graph_m, only: dist_lim, min_inters_speed use derived_types, only: snapshot use read_grid_m, only: nlon, nlat, periodic - use spher_polygon_area_m, only: spher_polygon_area - use spher_polyline_area_m, only: spher_polyline_area use unit_edge_m, only: unit_edge type(snapshot), intent(inout):: flow(:) ! (max_delta + 1) @@ -41,7 +39,6 @@ contains integer i2 ! eddy index of successor integer l, n_select type(polyline) polyline_1, polyline_2 - type(polygon) res_pol integer, allocatable:: selection(:) ! identifying numbers of a selection of extrema @@ -99,20 +96,11 @@ contains + floor((e1%extr%coord(1) - e2%extr%coord(1)) / twopi & + 0.5) * twopi - call gpc_polygon_clip_f(GPC_INT, polygon(nparts = 1, & - part = [polyline_1], hole = [.false.]), polygon(nparts = 1, & - part = [polyline_2], hole = [.false.]), res_pol) - - if (res_pol%nparts /= 0) then - ! polyline_1 and polyline_2 overlap - if (spher_polygon_area(res_pol) >= min_inters_speed & - * min(abs(spher_polyline_area(polyline_1)), & - abs(spher_polyline_area(polyline_2)))) then - write(unit_edge, fmt = *) & - (k - delta) * e_overestim + i1, k * e_overestim + i2 - e1%delta_out = min(e1%delta_out, delta) - e2%delta_in = min(e2%delta_in, delta) - end if + if (inters_lines(polyline_1, polyline_2, min_inters_speed)) then + write(unit_edge, fmt = *) (k - delta) * e_overestim + i1, & + k * e_overestim + i2 + e1%delta_out = min(e1%delta_out, delta) + e2%delta_in = min(e2%delta_in, delta) end if end associate end DO