diff --git a/Tests/Stdin/read_snapshot_nml.txt b/Tests/Stdin/read_snapshot_nml.txt new file mode 100644 index 0000000000000000000000000000000000000000..f4f7c4fee69eb4d70a80dedf31691122fcac3700 --- /dev/null +++ b/Tests/Stdin/read_snapshot_nml.txt @@ -0,0 +1,5 @@ + &MAIN_NML + CORNER = 16.12500 , -38.87500 , + NLON = 20, + NLAT = 20 + / diff --git a/Tests/read_snapshot.f b/Tests/read_snapshot.f index 3b0086239f42c09f847470deefe55ba1e85a0802..dfee05862c2d57414747063004acf3ab22427369 100644 --- a/Tests/read_snapshot.f +++ b/Tests/read_snapshot.f @@ -4,13 +4,14 @@ module read_snapshot_m contains - subroutine read_snapshot(s, hshp_extremum, hshp_outermost, hshp_max_speed) + subroutine read_snapshot(s, hshp_extremum, hshp_outermost, hshp_max_speed, & + corner, nlon, nlat) use, intrinsic:: ISO_C_BINDING use contour_531, only: convert_to_ind use derived_types, only: snapshot, eddy - use nr_util, only: pi + use nr_util, only: pi, assert, deg_to_rad use read_eddy_m, only: read_eddy use shapelib, only: shpfileobject, shpgetinfo @@ -23,12 +24,16 @@ contains TYPE(shpfileobject), intent(inout):: hshp_max_speed ! shapefile max_speed_contour_1 + real, intent(in):: corner(:) ! (2) longitude and latitude of the + ! corner of the whole grid, in rad + + integer, intent(in):: nlon, nlat + ! Local: integer ishape, k, i integer shapetype, dbffieldcount, dbfrecordcount real(c_double) minbound(4), maxbound(4) type(eddy) e - real, parameter:: deg_to_rad = pi / 180. !--------------------------------------------------------------------- @@ -43,17 +48,20 @@ contains end do allocate(s%ind_extr(2, s%number_vis_eddies)) - forall (i = 1:s%number_vis_eddies) s%ind_extr(:, i) & - = nint(convert_to_ind(s%list_vis(i)%coord_extr, & - corner = [0.125, -59.875] * deg_to_rad, & + = nint(convert_to_ind(s%list_vis(i)%coord_extr, corner, & step = [0.25, 0.25] * deg_to_rad)) - allocate(s%extr_map(120, 180)) + allocate(s%extr_map(nlon, nlat)) s%extr_map = 0 - forall (i = 1:s%number_vis_eddies) & - s%extr_map(s%ind_extr(1, i), s%ind_extr(2, i)) & - = merge(i, -i, s%list_vis(i)%suff_amp) + + do i = 1, s%number_vis_eddies + s%extr_map(s%ind_extr(1, i), s%ind_extr(2, i)) & + = merge(i, -i, s%list_vis(i)%suff_amp) + end do + + call assert(count(s%extr_map /= 0) == s%number_vis_eddies, & + "read_snapshot: one eddy per grid point") s%number_eddies = s%number_vis_eddies diff --git a/Tests/test_get_snapshot.f b/Tests/test_get_snapshot.f index 5ad99eec09f1af030186adaeb0189166eb2ed558..009742d2260082250809a0c69f1595fc0a1f2f38 100644 --- a/Tests/test_get_snapshot.f +++ b/Tests/test_get_snapshot.f @@ -11,14 +11,12 @@ program test_get_snapshot 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 nr_util, only: pi, assert, deg_to_rad use shapelib, only: shpfileobject, shpclose implicit none type(snapshot) s - real, parameter:: deg_over_rad = pi / 180. - TYPE(shpfileobject) hshp_extremum ! shapefile extremum_$m TYPE(shpfileobject) hshp_outermost ! shapefile outermost_contour_$m TYPE(shpfileobject) hshp_max_speed ! shapefile max_speed_contour_$m @@ -57,9 +55,9 @@ program test_get_snapshot call get_snapshot(s, min_amp, m = 1, n_proc = 1, k_end = 1, max_delta = 4, & nlon = nlon, nlat = nlat, k = 1, max_radius = [20, 20], & - corner = [longitude(1), latitude(1)] * deg_over_rad, & + corner = [longitude(1), latitude(1)] * deg_to_rad, & step = & - [longitude(2) - longitude(1), latitude(2) - latitude(1)] * deg_over_rad) + [longitude(2) - longitude(1), latitude(2) - latitude(1)] * deg_to_rad) call init_shapefiles(hshp_extremum, hshp_outermost, hshp_max_speed) diff --git a/Tests/test_local_extrema.f b/Tests/test_local_extrema.f index 5cc2a1af37fc2447f7a5efcdf1f9f8a841923d83..06a79f5bdb0e7ee3296bc50985a871123cb7561f 100644 --- a/Tests/test_local_extrema.f +++ b/Tests/test_local_extrema.f @@ -5,7 +5,7 @@ program test_local_extrema use netcdf, only: nf90_nowrite, nf90_clobber, NF90_FLOAT, NF90_INT use netcdf95, only: nf95_open, nf95_close, nf95_inq_varid, nf95_get_var, & nf95_put_var, nf95_create, nf95_def_dim, nf95_def_var, nf95_put_att, & - nf95_enddef, nf95_gw_var + nf95_enddef, nf95_gw_var, nf95_get_att use nr_util, only: assert implicit none @@ -28,6 +28,7 @@ program test_local_extrema ! level of innermost contour, for each extremum logical, allocatable:: cyclone(:) ! (n_extr) + real Fill_Value !--------------------------------------------------------------------- @@ -51,6 +52,8 @@ program test_local_extrema call nf95_inq_varid(ncid, "adt", varid) call nf95_get_var(ncid, varid, ssh) + call nf95_get_att(ncid, varid, "_FillValue", Fill_Value) + where (ssh == Fill_Value) ssh = huge(0.) call nf95_close(ncid) @@ -59,10 +62,13 @@ program test_local_extrema call new_unit(unit) open(unit, file = "test_local_extrema.csv", status = "replace", & action = "write") - write(unit, fmt = *) '"longitude index" "latitude index" cyclone ', & - '"longitude (degrees east)" ', & - '"latitude (degrees north)" "ssh (m)" "innermost_level (m)"' - ! (title line) + + ! Title lines: + write(unit, fmt = *) '"longitude index" "latitude index"' + write(unit, fmt = *) '"" "" "" "degrees east" "degrees north" "m" "m"' + write(unit, fmt = *) 'ilon ilat cyclone longitude latitude ssh ', & + 'innermost_level' + do i = 1, size(ind_extr, 2) write(unit, fmt = *) ind_extr(1, i), ind_extr(2, i), cyclone(i), & longitude(ind_extr(1, i)), latitude(ind_extr(2, i)), & diff --git a/Tests/test_read_snapshot.f b/Tests/test_read_snapshot.f index 1eeacf646d5d2a42ec4264573a89411b18b7ec19..a5e86e8514ecf3bc37c59e47d5030caec2ff5d9e 100644 --- a/Tests/test_read_snapshot.f +++ b/Tests/test_read_snapshot.f @@ -1,9 +1,12 @@ program test_read_snapshot + use, intrinsic:: ISO_FORTRAN_ENV + use derived_types, only: snapshot use dispatch_snapshot_m, only: dispatch_snapshot use init_shapefiles_m, only: init_shapefiles use jumble, only: new_unit + use nr_util, only: deg_to_rad use read_snapshot_m, only: read_snapshot use shapelib, only: shpfileobject, shpclose use shapelib_03, only: shp_open_03 @@ -14,22 +17,34 @@ program test_read_snapshot TYPE(shpfileobject) hshp_extremum ! shapefile extremum_1 TYPE(shpfileobject) hshp_outermost ! shapefile outermost_contour_1 TYPE(shpfileobject) hshp_max_speed ! shapefile max_speed_contour_1 - integer unit_isolated, unit_number_eddies + integer unit_isolated, unit_number_eddies, i + + real:: corner(2) = [0.125, - 59.875] + ! longitude and latitude of the corner of the whole grid, in degrees + + integer:: nlon = 120, nlat = 120 + namelist /main_nml/ corner, nlon, nlat !------------------------------------------------------------------------- + write(unit = error_unit, nml = main_nml) + write(unit = error_unit, fmt = *) "Enter namelist main_nml." + read(unit = *, nml = main_nml) + write(unit = *, nml = main_nml) + call shp_open_03("extremum_1_old", pszaccess = "rb", hshp = hshp_extremum) call shp_open_03("outermost_contour_1_old", pszaccess = "rb", & hshp = hshp_outermost) call shp_open_03("max_speed_contour_1_old", pszaccess = "rb", & hshp = hshp_max_speed) - call read_snapshot(s, hshp_extremum, hshp_outermost, hshp_max_speed) + call read_snapshot(s, hshp_extremum, hshp_outermost, hshp_max_speed, & + corner * deg_to_rad, nlon, nlat) CALL shpclose(hshp_extremum) CALL shpclose(hshp_outermost) CALL shpclose(hshp_max_speed) 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") @@ -57,4 +72,10 @@ program test_read_snapshot close(unit_number_eddies) print *, 'Created "number_eddies_1.csv".' + print *, "s%ind_extr:" + + do i = 1, s%number_vis_eddies + print *, s%ind_extr(:, i) + end do + end program test_read_snapshot diff --git a/Tests/tests.json b/Tests/tests.json index 15e2a0f67c463fa315a0672bf9e9ad39463090d3..792bf929abd2119d2f678444f9784389f6c46502 100644 --- a/Tests/tests.json +++ b/Tests/tests.json @@ -211,6 +211,13 @@ { "args": "$compil_prod_dir/test_read_snapshot", "required": ["$input_dir/Test_read_eddy/*"], - "title": "Read_snapshot" + "title": "Read_snapshot", + "stdin": "$stdin_dir/read_snapshot_nml.txt" + }, + { + "args" : ["$compil_prod_dir/test_local_extrema", + "$input_dir/2006_01/h_region_4.nc"], + "title" : "Local_extrema_missing", + "description": "test_local_extrema with input file containing missing values." } ] diff --git a/get_snapshot.f b/get_snapshot.f index 648b9afb5ab31463a2688b8e52b4dcb27f8c5ad3..835efd138817525ba603da2a8084e420de7fbf5f 100644 --- a/get_snapshot.f +++ b/get_snapshot.f @@ -36,7 +36,7 @@ contains ! eddy in longitude and latitude, in number of grid points real, intent(in):: corner(:) ! (2) longitude and latitude of the - ! corner of the global grid, in rad + ! corner of the whole grid, in rad real, intent(in):: step(:) ! (2) longitude and latitude steps, in rad ! Local: