PartMC  2.2.1
bin_grid.F90
Go to the documentation of this file.
00001 ! Copyright (C) 2005-2011 Nicole Riemer and Matthew West
00002 ! Licensed under the GNU General Public License version 2 or (at your
00003 ! option) any later version. See the file COPYING for details.
00004 
00005 !> \file
00006 !> The pmc_bin_grid module.
00007 
00008 !> The bin_grid_t structure and associated subroutines.
00009 module pmc_bin_grid
00010 
00011   use pmc_constants
00012   use pmc_util
00013   use pmc_spec_file
00014   use pmc_mpi
00015   use pmc_netcdf
00016 #ifdef PMC_USE_MPI
00017   use mpi
00018 #endif
00019 
00020   !> 1D grid of size bins.
00021   !!
00022   !! The grid of bins is logarithmically spaced in volume, an
00023   !! assumption that is quite heavily incorporated into the code. At
00024   !! some point in the future it would be nice to relax this
00025   !! assumption.
00026   type bin_grid_t
00027      !> Number of bins.
00028      integer :: n_bin
00029      !> Len n_bin, bin center radii (m^3).
00030      real(kind=dp), pointer :: center_radius(:)
00031      !> Len (n_bin + 1), bin edge radii (m^3).
00032      real(kind=dp), pointer :: edge_radius(:)
00033      !> Bin logarithmic width, equal to <tt>log(edge_radius(i+1)) -
00034      !> log(edge_radius(i))</tt> for any \c i (dimensionless).
00035      real(kind=dp) :: log_width
00036   end type bin_grid_t
00037 
00038 contains
00039   
00040 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00041 
00042   !> Allocates a bin_grid.
00043   subroutine bin_grid_allocate(bin_grid)
00044 
00045     !> Bin grid.
00046     type(bin_grid_t), intent(out) :: bin_grid
00047 
00048     bin_grid%n_bin = 0
00049     allocate(bin_grid%center_radius(0))
00050     allocate(bin_grid%edge_radius(0))
00051 
00052   end subroutine bin_grid_allocate
00053 
00054 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00055 
00056   !> Allocates a bin_grid of the given size.
00057   subroutine bin_grid_allocate_size(bin_grid, n_bin)
00058 
00059     !> Bin grid.
00060     type(bin_grid_t), intent(out) :: bin_grid
00061     !> Number of bins.
00062     integer, intent(in) :: n_bin
00063 
00064     bin_grid%n_bin = n_bin
00065     allocate(bin_grid%center_radius(n_bin))
00066     allocate(bin_grid%edge_radius(n_bin + 1))
00067 
00068   end subroutine bin_grid_allocate_size
00069 
00070 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00071 
00072   !> Frees all memory.
00073   subroutine bin_grid_deallocate(bin_grid)
00074 
00075     !> Bin_grid to free.
00076     type(bin_grid_t), intent(inout) :: bin_grid
00077 
00078     deallocate(bin_grid%center_radius)
00079     deallocate(bin_grid%edge_radius)
00080 
00081   end subroutine bin_grid_deallocate
00082 
00083 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00084 
00085   !> Copies a bin grid.
00086   subroutine bin_grid_copy(bin_grid_from, bin_grid_to)
00087 
00088     !> Bin_grid to copy from.
00089     type(bin_grid_t), intent(in) :: bin_grid_from
00090     !> Bin_grid to copy to.
00091     type(bin_grid_t), intent(inout) :: bin_grid_to
00092 
00093     if (bin_grid_from%n_bin /= bin_grid_to%n_bin) then
00094        call bin_grid_deallocate(bin_grid_to)
00095        call bin_grid_allocate_size(bin_grid_to, bin_grid_from%n_bin)
00096     end if
00097     bin_grid_to%center_radius = bin_grid_from%center_radius
00098     bin_grid_to%edge_radius = bin_grid_from%edge_radius
00099     bin_grid_to%log_width = bin_grid_from%log_width
00100 
00101   end subroutine bin_grid_copy
00102 
00103 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00104 
00105   !> Returns the size of a bin in a grid.
00106   real(kind=dp) function bin_grid_bin_size(bin_grid, i_bin)
00107 
00108     !> Bin grid.
00109     type(bin_grid_t), intent(in) :: bin_grid
00110     !> Bin number to return size of.
00111     integer, intent(in) :: i_bin
00112 
00113     bin_grid_bin_size = bin_grid%log_width
00114 
00115   end function bin_grid_bin_size
00116 
00117 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00118 
00119   !> Convert a concentration f(vol)d(vol) to f(ln(r))d(ln(r))
00120   !> where vol = 4/3 pi r^3.
00121   subroutine vol_to_lnr(r, f_vol, f_lnr)
00122     
00123     !> Radius (m).
00124     real(kind=dp), intent(in) :: r
00125     !> Concentration as a function of volume.
00126     real(kind=dp), intent(in) :: f_vol
00127     !> Concentration as a function of ln(r).
00128     real(kind=dp), intent(out) :: f_lnr
00129     
00130     f_lnr = f_vol * 4d0 * const%pi * r**3
00131     
00132   end subroutine vol_to_lnr
00133   
00134 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00135 
00136   !> Generates the bin grid given the range and number of bins.
00137   subroutine bin_grid_make(bin_grid, n_bin, r_min, r_max)
00138     
00139     !> New bin grid.
00140     type(bin_grid_t), intent(inout) :: bin_grid
00141     !> Number of bins.
00142     integer, intent(in) :: n_bin
00143     !> Minimum radius (m^3).
00144     real(kind=dp), intent(in) :: r_min
00145     !> Minimum radius (m^3).
00146     real(kind=dp), intent(in) :: r_max
00147 
00148     integer :: i_bin
00149 
00150     call assert_msg(538534122, n_bin > 0, &
00151          "bin_grid requires a positive n_bin, not: " &
00152          // trim(integer_to_string(n_bin)))
00153     call assert_msg(966541762, r_min > 0d0, &
00154          "bin_grid requires a positive r_min, not: " &
00155          // trim(real_to_string(r_min)))
00156     call assert_msg(711537859, r_min < r_max, &
00157          "bin_grid requires r_min < r_max, not: " &
00158          // trim(real_to_string(r_min)) // " and " &
00159          // trim(real_to_string(r_max)))
00160     call bin_grid_deallocate(bin_grid)
00161     call bin_grid_allocate_size(bin_grid, n_bin)
00162     call logspace(r_min, r_max, bin_grid%edge_radius)
00163     do i_bin = 1,n_bin
00164        bin_grid%center_radius(i_bin) &
00165             = exp(0.5d0 * log(bin_grid%edge_radius(i_bin)) &
00166             + 0.5d0 * log(bin_grid%edge_radius(i_bin + 1)))
00167     end do
00168     bin_grid%log_width = (log(r_max) - log(r_min)) / real(n_bin, kind=dp)
00169 
00170   end subroutine bin_grid_make
00171 
00172 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00173 
00174   !> Find the bin number that contains a given particle.
00175   !!
00176   !! This assumes logarithmically spaced bins. If a particle is below
00177   !! the smallest bin then its bin number is 0. If a particle is above
00178   !! the largest bin then its bin number is <tt>n_bin + 1</tt>.
00179   integer function bin_grid_particle_in_bin(bin_grid, radius)
00180 
00181     !> Bin_grid.
00182     type(bin_grid_t), intent(in) :: bin_grid
00183     !> Radius of particle.
00184     real(kind=dp), intent(in) :: radius
00185 
00186     call assert(448215689, bin_grid%n_bin >= 1)
00187     bin_grid_particle_in_bin = logspace_find(bin_grid%edge_radius(1), &
00188          bin_grid%edge_radius(bin_grid%n_bin + 1), bin_grid%n_bin + 1, &
00189          radius)
00190 
00191   end function bin_grid_particle_in_bin
00192   
00193 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00194 
00195   !> Make a histogram with of the given weighted data, scaled by the
00196   !> bin sizes.
00197   subroutine bin_grid_histogram_1d(x_bin_grid, x_data, weight_data, hist)
00198 
00199     !> x-axis bin grid.
00200     type(bin_grid_t), intent(in) :: x_bin_grid
00201     !> Data values on the x-axis.
00202     real(kind=dp), intent(in) :: x_data(:)
00203     !> Data value weights.
00204     real(kind=dp), intent(in) :: weight_data(size(x_data))
00205     !> Histogram to compute.
00206     real(kind=dp), intent(out) :: hist(x_bin_grid%n_bin)
00207 
00208     integer :: i_data, i_bin
00209 
00210     hist = 0d0
00211     do i_data = 1,size(x_data)
00212        i_bin = bin_grid_particle_in_bin(x_bin_grid, x_data(i_data))
00213        if ((i_bin >= 1) .and. (i_bin <= x_bin_grid%n_bin)) then
00214           hist(i_bin) = hist(i_bin) &
00215                + weight_data(i_data) / bin_grid_bin_size(x_bin_grid, i_bin)
00216        end if
00217     end do
00218 
00219   end subroutine bin_grid_histogram_1d
00220 
00221 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00222 
00223   !> Read the specification for a bin_grid from a spec file and
00224   !> generate it.
00225   subroutine spec_file_read_bin_grid(file, bin_grid)
00226 
00227     !> Spec file.
00228     type(spec_file_t), intent(inout) :: file
00229     !> Bin grid.
00230     type(bin_grid_t), intent(inout) :: bin_grid
00231 
00232     integer :: n_bin
00233     real(kind=dp) :: d_min, d_max
00234 
00235     !> \page input_format_bin_grid Input File Format: Diameter Axis Bin Grid
00236     !!
00237     !! The bin grid is logarithmic in diameter, consisting of
00238     !! \f$n_{\rm bin}\f$ bins with centers \f$c_i\f$ (\f$i =
00239     !! 1,\ldots,n_{\rm bin}\f$) and edges \f$e_i\f$ (\f$i =
00240     !! 1,\ldots,(n_{\rm bin} + 1)\f$) such that \f$e_{i+1}/e_i\f$ is a
00241     !! constant and \f$c_i/e_i = \sqrt{e_{i+1}/e_i}\f$. That is,
00242     !! \f$\ln(e_i)\f$ are uniformly spaced and \f$\ln(c_i)\f$ are the
00243     !! arithmetic centers.
00244     !!
00245     !! The diameter axis bin grid is specified by the parameters:
00246     !!   - \b n_bin (integer): The number of bins \f$n_{\rm bin}\f$.
00247     !!   - \b d_min (real, unit m): The left edge of the left-most bin,
00248     !!     \f$e_1\f$.
00249     !!   - \b d_max (real, unit m): The right edge of the right-most bin,
00250     !!     \f$e_{n_{\rm bin} + 1}\f$.
00251     !!
00252     !! See also:
00253     !!   - \ref spec_file_format --- the input file text format
00254     !!   - \ref output_format_bin_grid --- the corresponding output format
00255 
00256     call spec_file_read_integer(file, 'n_bin', n_bin)
00257     call spec_file_read_real(file, 'd_min', d_min)
00258     call spec_file_read_real(file, 'd_max', d_max)
00259     call bin_grid_make(bin_grid, n_bin, diam2rad(d_min), diam2rad(d_max))
00260 
00261   end subroutine spec_file_read_bin_grid
00262 
00263 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00264 
00265   !> Determines the number of bytes required to pack the given value.
00266   integer function pmc_mpi_pack_size_bin_grid(val)
00267 
00268     !> Value to pack.
00269     type(bin_grid_t), intent(in) :: val
00270 
00271     pmc_mpi_pack_size_bin_grid = &
00272          pmc_mpi_pack_size_integer(val%n_bin) &
00273          + pmc_mpi_pack_size_real_array(val%center_radius) &
00274          + pmc_mpi_pack_size_real_array(val%edge_radius) &
00275          + pmc_mpi_pack_size_real(val%log_width)
00276 
00277   end function pmc_mpi_pack_size_bin_grid
00278 
00279 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00280 
00281   !> Packs the given value into the buffer, advancing position.
00282   subroutine pmc_mpi_pack_bin_grid(buffer, position, val)
00283 
00284     !> Memory buffer.
00285     character, intent(inout) :: buffer(:)
00286     !> Current buffer position.
00287     integer, intent(inout) :: position
00288     !> Value to pack.
00289     type(bin_grid_t), intent(in) :: val
00290 
00291 #ifdef PMC_USE_MPI
00292     integer :: prev_position
00293 
00294     prev_position = position
00295     call pmc_mpi_pack_integer(buffer, position, val%n_bin)
00296     call pmc_mpi_pack_real_array(buffer, position, val%center_radius)
00297     call pmc_mpi_pack_real_array(buffer, position, val%edge_radius)
00298     call pmc_mpi_pack_real(buffer, position, val%log_width)
00299     call assert(385455586, &
00300          position - prev_position <= pmc_mpi_pack_size_bin_grid(val))
00301 #endif
00302 
00303   end subroutine pmc_mpi_pack_bin_grid
00304 
00305 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00306 
00307   !> Unpacks the given value from the buffer, advancing position.
00308   subroutine pmc_mpi_unpack_bin_grid(buffer, position, val)
00309 
00310     !> Memory buffer.
00311     character, intent(inout) :: buffer(:)
00312     !> Current buffer position.
00313     integer, intent(inout) :: position
00314     !> Value to pack.
00315     type(bin_grid_t), intent(inout) :: val
00316 
00317 #ifdef PMC_USE_MPI
00318     integer :: prev_position
00319 
00320     prev_position = position
00321     call pmc_mpi_unpack_integer(buffer, position, val%n_bin)
00322     call pmc_mpi_unpack_real_array(buffer, position, val%center_radius)
00323     call pmc_mpi_unpack_real_array(buffer, position, val%edge_radius)
00324     call pmc_mpi_unpack_real(buffer, position, val%log_width)
00325     call assert(741838730, &
00326          position - prev_position <= pmc_mpi_pack_size_bin_grid(val))
00327 #endif
00328 
00329   end subroutine pmc_mpi_unpack_bin_grid
00330 
00331 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00332 
00333   !> Check whether all processors have the same value.
00334   logical function pmc_mpi_allequal_bin_grid(val)
00335 
00336     !> Value to compare.
00337     type(bin_grid_t), intent(inout) :: val
00338 
00339 #ifdef PMC_USE_MPI
00340     if (.not. pmc_mpi_allequal_integer(val%n_bin)) then
00341        pmc_mpi_allequal_bin_grid = .false.
00342        return
00343     end if
00344 
00345     if (val%n_bin == 0) then
00346        pmc_mpi_allequal_bin_grid = .true.
00347        return
00348     end if
00349 
00350     if (pmc_mpi_allequal_real(val%edge_radius(1)) &
00351          .and. pmc_mpi_allequal_real(val%edge_radius(val%n_bin))) then
00352        pmc_mpi_allequal_bin_grid = .true.
00353     else
00354        pmc_mpi_allequal_bin_grid = .false.
00355     end if
00356 #else
00357     pmc_mpi_allequal_bin_grid = .true.
00358 #endif
00359 
00360   end function pmc_mpi_allequal_bin_grid
00361 
00362 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00363 
00364   !> Write the aero_diam dimension to the given NetCDF file if it is
00365   !> not already present and in any case return the associated dimid.
00366   subroutine bin_grid_netcdf_dim_aero_diam(bin_grid, ncid, &
00367        dimid_aero_diam)
00368 
00369     !> Bin_grid structure.
00370     type(bin_grid_t), intent(in) :: bin_grid
00371     !> NetCDF file ID, in data mode.
00372     integer, intent(in) :: ncid
00373     !> Dimid of the aero_diam dimension.
00374     integer, intent(out) :: dimid_aero_diam
00375 
00376     integer :: status, varid_aero_diam
00377     integer :: dimid_aero_diam_edges, varid_aero_diam_edges, 
00378          varid_aero_diam_widths
00379     real(kind=dp) :: aero_diam_centers(bin_grid%n_bin), 
00380          aero_diam_edges(bin_grid%n_bin + 1)
00381     real(kind=dp) :: aero_diam_widths(bin_grid%n_bin)
00382 
00383     status = nf90_inq_dimid(ncid, "aero_diam", dimid_aero_diam)
00384     if (status == NF90_NOERR) return
00385     if (status /= NF90_EBADDIM) call pmc_nc_check(status)
00386 
00387     ! dimension not defined, so define now define it
00388     call pmc_nc_check(nf90_redef(ncid))
00389 
00390     call pmc_nc_check(nf90_def_dim(ncid, "aero_diam", &
00391          bin_grid%n_bin, dimid_aero_diam))
00392     call pmc_nc_check(nf90_def_dim(ncid, "aero_diam_edges", &
00393          bin_grid%n_bin + 1, dimid_aero_diam_edges))
00394 
00395     call pmc_nc_check(nf90_enddef(ncid))
00396 
00397     aero_diam_centers = rad2diam(bin_grid%center_radius)
00398     aero_diam_widths = bin_grid%log_width
00399     aero_diam_edges = rad2diam(bin_grid%edge_radius)
00400 
00401     call pmc_nc_write_real_1d(ncid, aero_diam_centers, &
00402          "aero_diam", (/ dimid_aero_diam /), unit="m", &
00403          long_name="aerosol diameter axis bin centers", &
00404          description="logarithmically spaced centers of diameter axis grid, " &
00405          // "so that aero_diam(i) / aero_diam_edges(i) = " &
00406          // "0.5 * aero_diam_edges(i+1) / aero_diam_edges(i)")
00407     call pmc_nc_write_real_1d(ncid, aero_diam_edges, &
00408          "aero_diam_edges", (/ dimid_aero_diam_edges /), unit="m", &
00409          long_name="aerosol diameter axis bin edges", &
00410          description="logarithmically spaced edges of diameter axis grid, " &
00411          // "with one more edge than center")
00412     call pmc_nc_write_real_1d(ncid, aero_diam_widths, &
00413          "aero_diam_widths", (/ dimid_aero_diam /), unit="m", &
00414          long_name="aerosol diameter axis bin widths", &
00415          description="base-e logarithmic widths of diameter axis grid, " &
00416          // "so that aero_diam_widths(i)" &
00417          // "= ln(aero_diam_edges(i+1) / aero_diam_edges(i)) and " &
00418          // "all bins have the same width")
00419 
00420   end subroutine bin_grid_netcdf_dim_aero_diam
00421 
00422 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00423 
00424   !> Write full state.
00425   subroutine bin_grid_output_netcdf(bin_grid, ncid)
00426     
00427     !> bin_grid to write.
00428     type(bin_grid_t), intent(in) :: bin_grid
00429     !> NetCDF file ID, in data mode.
00430     integer, intent(in) :: ncid
00431 
00432     integer :: dimid_aero_diam
00433 
00434     !> \page output_format_bin_grid Output File Format: Bin Grid Data
00435     !!
00436     !! The bin grid data NetCDF dimensions are:
00437     !!   - \b aero_diam: number of bins (grid cells) on the diameter axis
00438     !!   - \b aero_diam_edges: number of bin edges (grid cell edges) on
00439     !!     the diameter axis --- always equal to <tt>aero_diam + 1</tt>
00440     !!
00441     !! The bin grid data NetCDF variables are:
00442     !!   - \b aero_diam (unit m, dim \c aero_diam): aerosol diameter axis
00443     !!     bin centers --- centered on a logarithmic scale from the edges, so
00444     !!     that <tt>aero_diam(i) / aero_diam_edges(i) =
00445     !!     sqrt(aero_diam_edges(i+1) / aero_diam_edges(i))</tt>
00446     !!   - \b aero_diam_edges (unit m, dim \c aero_diam_edges): aersol
00447     !!     diameter axis bin edges (there is one more edge than center)
00448     !!   - \b aero_diam_widths (dimensionless, dim \c aero_diam):
00449     !!     the base-e logarithmic bin widths --- <tt>aero_diam_widths(i)
00450     !!     = ln(aero_diam_edges(i+1) / aero_diam_edges(i))</tt>, so
00451     !!     all bins have the same width
00452     !!
00453     !! See also:
00454     !!   - \ref input_format_bin_grid --- the corresponding input format
00455 
00456     call bin_grid_netcdf_dim_aero_diam(bin_grid, ncid, &
00457          dimid_aero_diam)
00458 
00459     ! no need to write any more data as it's all contained in the
00460     ! dimension and associated variables
00461 
00462   end subroutine bin_grid_output_netcdf
00463 
00464 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00465 
00466   !> Read full state.
00467   subroutine bin_grid_input_netcdf(bin_grid, ncid)
00468     
00469     !> bin_grid to read.
00470     type(bin_grid_t), intent(inout) :: bin_grid
00471     !> NetCDF file ID, in data mode.
00472     integer, intent(in) :: ncid
00473 
00474     integer :: dimid_aero_diam, n_bin
00475     character(len=1000) :: name
00476     real(kind=dp), allocatable :: aero_diam_centers(:)
00477     real(kind=dp), allocatable :: aero_diam_edges(:)
00478     real(kind=dp), allocatable :: aero_diam_widths(:)
00479 
00480     call pmc_nc_check(nf90_inq_dimid(ncid, "aero_diam", dimid_aero_diam))
00481     call pmc_nc_check(nf90_Inquire_Dimension(ncid, dimid_aero_diam, name, &
00482          n_bin))
00483 
00484     call bin_grid_deallocate(bin_grid)
00485     call bin_grid_allocate_size(bin_grid, n_bin)
00486 
00487     allocate(aero_diam_centers(n_bin))
00488     allocate(aero_diam_edges(n_bin + 1))
00489     allocate(aero_diam_widths(n_bin))
00490 
00491     call pmc_nc_read_real_1d(ncid, aero_diam_centers, &
00492          "aero_diam")
00493     call pmc_nc_read_real_1d(ncid, aero_diam_edges, &
00494          "aero_diam_edges")
00495     call pmc_nc_read_real_1d(ncid, aero_diam_widths, &
00496          "aero_diam_widths")
00497 
00498     bin_grid%center_radius = diam2rad(aero_diam_centers)
00499     bin_grid%edge_radius = diam2rad(aero_diam_edges)
00500     bin_grid%log_width = aero_diam_widths(1)
00501 
00502     deallocate(aero_diam_centers)
00503     deallocate(aero_diam_edges)
00504     deallocate(aero_diam_widths)
00505 
00506   end subroutine bin_grid_input_netcdf
00507 
00508 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
00509   
00510 end module pmc_bin_grid