21 integer,
parameter :: BIN_GRID_TYPE_INVALID = 0
23 integer,
parameter :: BIN_GRID_TYPE_LOG = 1
25 integer,
parameter :: BIN_GRID_TYPE_LINEAR = 2
39 real(kind=dp),
allocatable :: centers(:)
41 real(kind=dp),
allocatable :: edges(:)
43 real(kind=dp),
allocatable :: widths(:)
56 bin_grid%type = bin_grid_type_invalid
58 allocate(bin_grid%centers(0))
59 allocate(bin_grid%edges(0))
60 allocate(bin_grid%widths(0))
72 integer,
intent(in) :: n_bin
74 bin_grid%type = bin_grid_type_invalid
75 bin_grid%n_bin = n_bin
76 allocate(bin_grid%centers(n_bin))
77 allocate(bin_grid%edges(n_bin + 1))
78 allocate(bin_grid%widths(n_bin))
90 deallocate(bin_grid%centers)
91 deallocate(bin_grid%edges)
92 deallocate(bin_grid%widths)
104 type(bin_grid_t),
intent(inout) :: bin_grid_to
108 bin_grid_to%type = bin_grid_from%type
109 bin_grid_to%n_bin = bin_grid_from%n_bin
110 bin_grid_to%centers = bin_grid_from%centers
111 bin_grid_to%edges = bin_grid_from%edges
112 bin_grid_to%widths = bin_grid_from%widths
123 real(kind=dp),
intent(in) :: r
125 real(kind=dp),
intent(in) :: f_vol
127 real(kind=dp),
intent(out) :: f_lnr
129 f_lnr = f_vol * 4d0 * const%pi * r**3
141 integer,
intent(in) :: type
143 integer,
intent(in) :: n_bin
145 real(kind=dp),
intent(in) :: min
147 real(kind=dp),
intent(in) :: max
149 real(kind=dp) :: c1, c2
152 "bin_grid requires a non-negative n_bin, not: " &
155 if (
type == bin_grid_type_log)
then
157 "log bin_grid requires a positive min value, not: " &
161 "bin_grid requires min < max, not: " &
168 if (n_bin == 0)
return
169 if (
type == bin_grid_type_log)
then
170 call
logspace(min, max, bin_grid%edges)
174 call
logspace(c1, c2, bin_grid%centers)
175 bin_grid%widths = (log(max) - log(min)) /
real(n_bin, kind=dp)
176 elseif (bin_grid%type == bin_grid_type_linear)
then
177 call
linspace(min, max, bin_grid%edges)
180 call
linspace(c1, c2, bin_grid%centers)
181 bin_grid%widths = (max - min) /
real(n_bin, kind=dp)
183 call
die_msg(678302366,
"unknown bin_grid type: " &
201 real(kind=dp),
intent(in) :: val
203 call
assert(448215689, bin_grid%n_bin >= 0)
205 if (bin_grid%n_bin == 0)
return
206 if (bin_grid%type == bin_grid_type_log)
then
208 bin_grid%edges(bin_grid%n_bin + 1), bin_grid%n_bin + 1, val)
209 elseif (bin_grid%type == bin_grid_type_linear)
then
211 bin_grid%edges(bin_grid%n_bin + 1), bin_grid%n_bin + 1, val)
213 call
die_msg(348908641,
"unknown bin_grid type: " &
228 real(kind=dp),
intent(in) :: x_data(:)
230 real(kind=dp),
intent(in) :: weight_data(size(x_data))
235 integer :: i_data, x_bin
238 do i_data = 1,
size(x_data)
240 if ((x_bin >= 1) .and. (x_bin <= x_bin_grid%n_bin))
then
242 + weight_data(i_data) / x_bin_grid%widths(x_bin)
258 real(kind=dp),
intent(in) :: x_data(:)
262 real(kind=dp),
intent(in) :: y_data(size(x_data))
264 real(kind=dp),
intent(in) :: weight_data(size(x_data))
269 integer :: i_data, x_bin, y_bin
272 do i_data = 1,
size(x_data)
275 if ((x_bin >= 1) .and. (x_bin <= x_bin_grid%n_bin) &
276 .and. (y_bin >= 1) .and. (y_bin <= y_bin_grid%n_bin))
then
279 / x_bin_grid%widths(x_bin) / y_bin_grid%widths(y_bin)
288 subroutine spec_file_read_radius_bin_grid(file, bin_grid)
296 real(kind=dp) :: d_min, d_max
325 end subroutine spec_file_read_radius_bin_grid
350 character,
intent(inout) :: buffer(:)
352 integer,
intent(inout) :: position
357 integer :: prev_position
359 prev_position = position
377 character,
intent(inout) :: buffer(:)
379 integer,
intent(inout) :: position
384 integer :: prev_position
386 prev_position = position
417 if (val%n_bin == 0)
then
444 integer,
intent(in) :: ncid
446 character(len=*),
intent(in) :: dim_name
448 character(len=*),
intent(in) :: unit
450 integer,
intent(out) :: dimid
452 character(len=*),
intent(in),
optional :: long_name
454 real(kind=dp),
intent(in),
optional :: scale
456 integer :: status, varid, dimid_edges, varid_edges, varid_widths, i
457 real(kind=dp) :: centers(bin_grid%n_bin), edges(bin_grid%n_bin + 1)
458 real(kind=dp) :: widths(bin_grid%n_bin)
459 character(len=(len_trim(dim_name)+10)) :: dim_name_edges
460 character(len=255) :: use_long_name
462 status = nf90_inq_dimid(ncid, dim_name, dimid)
463 if (status == nf90_noerr)
return
468 dim_name_edges = trim(dim_name) //
"_edges"
469 if (present(long_name))
then
470 call
assert_msg(125084459, len_trim(long_name) <= len(use_long_name), &
471 "long_name is longer than " &
473 use_long_name = trim(long_name)
475 call
assert_msg(660927086, len_trim(dim_name) <= len(use_long_name), &
476 "dim_name is longer than " &
478 use_long_name = trim(dim_name)
482 call
pmc_nc_check(nf90_def_dim(ncid, dim_name, bin_grid%n_bin, dimid))
483 call
pmc_nc_check(nf90_def_dim(ncid, dim_name_edges, bin_grid%n_bin + 1, &
487 centers = bin_grid%centers
488 edges = bin_grid%edges
489 widths = bin_grid%widths
490 if (bin_grid%type == bin_grid_type_log)
then
491 if (present(scale))
then
492 centers = centers * scale
493 edges = edges * scale
496 unit=unit, long_name=(trim(use_long_name) //
" grid centers"), &
497 description=(
"logarithmically spaced centers of " &
498 // trim(use_long_name) //
" grid, so that " // trim(dim_name) &
499 //
"(i) is the geometric mean of " // trim(dim_name_edges) &
500 //
"(i) and " // trim(dim_name_edges) //
"(i + 1)"))
502 (/ dimid_edges /), unit=unit, &
503 long_name=(trim(use_long_name) //
" grid edges"), &
504 description=(
"logarithmically spaced edges of " &
505 // trim(use_long_name) //
" grid, with one more edge than center"))
507 (/ dimid /), unit=
"1", &
508 long_name=(trim(use_long_name) //
" grid widths"), &
509 description=(
"base-e logarithmic widths of " &
510 // trim(use_long_name) //
" grid, with " // trim(dim_name) &
511 //
"_widths(i) = ln(" // trim(dim_name_edges) //
"(i + 1) / " &
512 // trim(dim_name_edges) //
"(i))"))
513 elseif (bin_grid%type == bin_grid_type_linear)
then
514 if (present(scale))
then
515 centers = centers * scale
516 edges = edges * scale
517 widths = widths * scale
520 unit=unit, long_name=(trim(use_long_name) //
" grid centers"), &
521 description=(
"linearly spaced centers of " // trim(use_long_name) &
522 //
" grid, so that " // trim(dim_name) //
"(i) is the mean of " &
523 // trim(dim_name_edges) //
"(i) and " // trim(dim_name_edges) &
526 (/ dimid_edges /), unit=unit, &
527 long_name=(trim(use_long_name) //
" grid edges"), &
528 description=(
"linearly spaced edges of " &
529 // trim(use_long_name) //
" grid, with one more edge than center"))
531 (/ dimid /), unit=unit, &
532 long_name=(trim(use_long_name) //
" grid widths"), &
533 description=(
"widths of " // trim(use_long_name) &
534 //
" grid, with " // trim(dim_name) //
"_widths(i) = " &
535 // trim(dim_name_edges) //
"(i + 1) - " // trim(dim_name_edges) &
538 call
die_msg(942560572,
"unknown bin_grid type: " &
553 integer,
intent(in) :: ncid
555 character(len=*),
intent(in) :: dim_name
557 character(len=*),
intent(in) :: unit
559 character(len=*),
intent(in),
optional :: long_name
561 real(kind=dp),
intent(in),
optional :: scale
578 integer,
intent(in) :: ncid
580 character(len=*),
intent(in) :: dim_name
582 real(kind=dp),
intent(in),
optional :: scale
584 integer :: dimid, varid, n_bin, type
585 character(len=1000) :: name, description
586 real(kind=dp),
allocatable :: edges(:)
588 call
pmc_nc_check(nf90_inq_dimid(ncid, dim_name, dimid))
589 call
pmc_nc_check(nf90_inquire_dimension(ncid, dimid, name, n_bin))
590 call
pmc_nc_check(nf90_inq_varid(ncid, dim_name, varid))
591 call
pmc_nc_check(nf90_get_att(ncid, varid,
"description", description))
593 allocate(edges(n_bin + 1))
596 if (
starts_with(description,
"logarithmically"))
then
597 type = bin_grid_type_log
599 type = bin_grid_type_linear
601 call
die_msg(792158584,
"cannot identify grid type for NetCDF " &
602 //
"dimension: " // trim(dim_name))
605 if (present(scale))
then
606 call
bin_grid_make(bin_grid, type, n_bin, scale * edges(1), &
607 scale * edges(n_bin + 1))
609 call
bin_grid_make(bin_grid, type, n_bin, edges(1), edges(n_bin + 1))
real(kind=dp) function, dimension(x_bin_grid%n_bin, y_bin_grid%n_bin) bin_grid_histogram_2d(x_bin_grid, x_data, y_bin_grid, y_data, weight_data)
Make a 2D histogram with of the given weighted data, scaled by the bin sizes.
subroutine pmc_nc_write_real_1d(ncid, var, name, dimids, dim_name, unit, long_name, standard_name, description)
Write a simple real array to a NetCDF file.
integer function pmc_mpi_pack_size_real_array(val)
Determines the number of bytes required to pack the given value.
An input file with extra data for printing messages.
integer function pmc_mpi_pack_size_integer(val)
Determines the number of bytes required to pack the given value.
subroutine die_msg(code, error_msg)
Error immediately.
subroutine bin_grid_allocate_size(bin_grid, n_bin)
Allocates a bin_grid.
real(kind=dp) function, dimension(x_bin_grid%n_bin) bin_grid_histogram_1d(x_bin_grid, x_data, weight_data)
Make a histogram with of the given weighted data, scaled by the bin sizes.
integer function linspace_find(min_x, max_x, n, x)
Find the position of a real number in a 1D linear array.
subroutine bin_grid_output_netcdf(bin_grid, ncid, dim_name, unit, long_name, scale)
Write a bin grid to the given NetCDF file.
logical function starts_with(string, start_string)
Checks whether a string starts with a given other string.
subroutine pmc_mpi_unpack_real_array_alloc(buffer, position, val)
Unpacks the given value from the buffer to an allocated array, advancing position.
subroutine pmc_nc_check(status)
Check the status of a NetCDF function call.
subroutine bin_grid_make(bin_grid, type, n_bin, min, max)
Generates the bin grid given the range and number of bins.
integer function bin_grid_find(bin_grid, val)
Find the bin number that contains a given value.
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
subroutine logspace(min_x, max_x, x)
Makes a logarithmically spaced array of length n from min to max.
logical function pmc_mpi_allequal_real(val)
Returns whether all processors have the same value.
subroutine bin_grid_netcdf_dim(bin_grid, ncid, dim_name, unit, dimid, long_name, scale)
Write a bin grid dimension to the given NetCDF file if it is not already present and in any case retu...
real(kind=dp) function interp_linear_disc(x_1, x_n, n, i)
Linear interpolation over discrete indices.
subroutine spec_file_read_integer(file, name, var)
Read an integer from a spec file that must have the given name.
subroutine bin_grid_allocate(bin_grid)
Allocates a bin_grid.
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine bin_grid_input_netcdf(bin_grid, ncid, dim_name, scale)
Read full state.
real(kind=dp) elemental function diam2rad(d)
Convert diameter (m) to radius (m).
Common utility subroutines.
integer function pmc_mpi_pack_size_bin_grid(val)
Determines the number of bytes required to pack the given value.
logical function pmc_mpi_allequal_integer(val)
Returns whether all processors have the same value.
subroutine vol_to_lnr(r, f_vol, f_lnr)
Convert a concentration f(vol)d(vol) to f(ln(r))d(ln(r)) where vol = 4/3 pi r^3.
subroutine pmc_nc_read_real_1d(ncid, var, name, must_be_present)
Read a simple real array from a NetCDF file.
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
logical function pmc_mpi_allequal_bin_grid(val)
Check whether all processors have the same value.
Wrapper functions for MPI.
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
1D grid, either logarithmic or linear.
The bin_grid_t structure and associated subroutines.
Reading formatted text input.
Wrapper functions for NetCDF. These all take a NetCDF ncid in data mode and return with it again in d...
subroutine pmc_mpi_pack_real_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine pmc_mpi_unpack_bin_grid(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
character(len=pmc_util_convert_string_len) function real_to_string(val)
Convert a real to a string format.
integer function logspace_find(min_x, max_x, n, x)
Find the position of a real number in a 1D logarithmic array.
subroutine spec_file_read_real(file, name, var)
Read a real number from a spec file that must have the given name.
subroutine linspace(min_x, max_x, x)
Makes a linearly spaced array from min to max.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
subroutine bin_grid_copy(bin_grid_from, bin_grid_to)
Copies a bin grid.
subroutine pmc_mpi_pack_bin_grid(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine bin_grid_deallocate(bin_grid)
Frees all memory.