PartMC 2.1.2
|
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_aero_data module. 00007 00008 !> The aero_data_t structure and associated subroutines. 00009 module pmc_aero_data 00010 00011 use pmc_spec_file 00012 use pmc_mpi 00013 use pmc_util 00014 use pmc_netcdf 00015 #ifdef PMC_USE_MPI 00016 use mpi 00017 #endif 00018 00019 integer, parameter :: AERO_NAME_LEN = 15 00020 integer, parameter :: AERO_SOURCE_NAME_LEN = 100 00021 00022 !> Aerosol material properties and associated data. 00023 !! 00024 !! The data in this structure is constant, as it represents physical 00025 !! quantities that cannot change over time. 00026 !! 00027 !! Each aerosol species is identified by an index <tt>i = 00028 !! 1,...,n_spec</tt>. Then \c name(i) is the name of that species, 00029 !! \c density(i) is its density, etc. The ordering of the species is 00030 !! arbitrary and should not be relied upon (currently it is the 00031 !! order in the species data file). The only exception is that it is 00032 !! possible to find out which species is water from the \c i_water 00033 !! variable. 00034 !! 00035 !! The names of the aerosol species are not important to PartMC, as 00036 !! only the material properties are used. The names are used for 00037 !! input and output, and also for communication with MOSAIC. For the 00038 !! MOSAIC interface to work correctly the species must be named the 00039 !! same, but without the \c _a suffix. 00040 type aero_data_t 00041 !> Number of species. 00042 integer :: n_spec 00043 !> Number of sources. 00044 integer :: n_source 00045 !> Water species number (0 if water is not a species). 00046 integer :: i_water 00047 !> Len n_spec, species names. 00048 character(len=AERO_NAME_LEN), pointer :: name(:) 00049 !> Length n_spec, mosaic_index(i) a positive integer giving the 00050 !> mosaic index of species i, or 0 if there is no match. 00051 integer, pointer :: mosaic_index(:) 00052 !> Len n_spec, densities (kg m^{-3}). 00053 real(kind=dp), pointer :: density(:) 00054 !> Len n_spec, number of ions in solute. 00055 integer, pointer :: num_ions(:) 00056 !> Len n_spec, molecular weights (kg mol^{-1}). 00057 real(kind=dp), pointer :: molec_weight(:) 00058 !> Len n_spec, kappas (1). 00059 real(kind=dp), pointer :: kappa(:) 00060 !> Len n_source, source names. 00061 character(len=AERO_SOURCE_NAME_LEN), pointer :: source_name(:) 00062 end type aero_data_t 00063 00064 contains 00065 00066 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00067 00068 !> Allocate storage for aero_data. 00069 subroutine aero_data_allocate(aero_data) 00070 00071 !> Aerosol data. 00072 type(aero_data_t), intent(out) :: aero_data 00073 00074 aero_data%n_spec = 0 00075 aero_data%n_source = 0 00076 allocate(aero_data%name(0)) 00077 allocate(aero_data%mosaic_index(0)) 00078 allocate(aero_data%density(0)) 00079 allocate(aero_data%num_ions(0)) 00080 allocate(aero_data%molec_weight(0)) 00081 allocate(aero_data%kappa(0)) 00082 allocate(aero_data%source_name(0)) 00083 aero_data%i_water = 0 00084 00085 end subroutine aero_data_allocate 00086 00087 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00088 00089 !> Allocate storage for aero_data parameters given the number of 00090 !> species. 00091 subroutine aero_data_allocate_size(aero_data, n_spec, n_source) 00092 00093 !> Aerosol data. 00094 type(aero_data_t), intent(out) :: aero_data 00095 !> Number of species. 00096 integer, intent(in) :: n_spec 00097 !> Number of sources. 00098 integer, intent(in) :: n_source 00099 00100 aero_data%n_spec = n_spec 00101 aero_data%n_source = n_source 00102 allocate(aero_data%name(n_spec)) 00103 allocate(aero_data%mosaic_index(n_spec)) 00104 allocate(aero_data%density(n_spec)) 00105 allocate(aero_data%num_ions(n_spec)) 00106 allocate(aero_data%molec_weight(n_spec)) 00107 allocate(aero_data%kappa(n_spec)) 00108 allocate(aero_data%source_name(n_spec)) 00109 aero_data%i_water = 0 00110 00111 end subroutine aero_data_allocate_size 00112 00113 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00114 00115 !> Frees all storage. 00116 subroutine aero_data_deallocate(aero_data) 00117 00118 !> Aerosol data. 00119 type(aero_data_t), intent(inout) :: aero_data 00120 00121 deallocate(aero_data%name) 00122 deallocate(aero_data%mosaic_index) 00123 deallocate(aero_data%density) 00124 deallocate(aero_data%num_ions) 00125 deallocate(aero_data%molec_weight) 00126 deallocate(aero_data%kappa) 00127 deallocate(aero_data%source_name) 00128 00129 end subroutine aero_data_deallocate 00130 00131 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00132 00133 !> Copy structure. 00134 subroutine aero_data_copy(aero_data_from, aero_data_to) 00135 00136 !> Source aerosol data. 00137 type(aero_data_t), intent(in) :: aero_data_from 00138 !> Destination aerosol data. 00139 type(aero_data_t), intent(inout) :: aero_data_to 00140 00141 call aero_data_deallocate(aero_data_to) 00142 call aero_data_allocate_size(aero_data_to, aero_data_from%n_spec, & 00143 aero_data_from%n_source) 00144 00145 aero_data_to%i_water = aero_data_from%i_water 00146 aero_data_to%name = aero_data_from%name 00147 aero_data_to%mosaic_index = aero_data_from%mosaic_index 00148 aero_data_to%density = aero_data_from%density 00149 aero_data_to%num_ions = aero_data_from%num_ions 00150 aero_data_to%molec_weight = aero_data_from%molec_weight 00151 aero_data_to%kappa = aero_data_from%kappa 00152 aero_data_to%source_name = aero_data_from%source_name 00153 00154 end subroutine aero_data_copy 00155 00156 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00157 00158 !> Returns the number of the species in aero_data with the given name, or 00159 !> returns 0 if there is no such species. 00160 integer function aero_data_spec_by_name(aero_data, name) 00161 00162 !> Aero_data data. 00163 type(aero_data_t), intent(in) :: aero_data 00164 !> Name of species to find. 00165 character(len=*), intent(in) :: name 00166 00167 integer i 00168 logical found 00169 00170 found = .false. 00171 do i = 1,aero_data%n_spec 00172 if (name == aero_data%name(i)) then 00173 found = .true. 00174 exit 00175 end if 00176 end do 00177 if (found) then 00178 aero_data_spec_by_name = i 00179 else 00180 aero_data_spec_by_name = 0 00181 end if 00182 00183 end function aero_data_spec_by_name 00184 00185 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00186 00187 !> Returns the number of the source in aero_data with the given name, or 00188 !> adds the source if it doesn't exist yet. 00189 integer function aero_data_source_by_name(aero_data, name) 00190 00191 !> Aero_data data. 00192 type(aero_data_t), intent(inout) :: aero_data 00193 !> Name of source to find. 00194 character(len=*), intent(in) :: name 00195 00196 integer :: i 00197 logical :: found 00198 character(len=AERO_SOURCE_NAME_LEN), pointer :: new_source_name(:) 00199 00200 found = .false. 00201 do i = 1,aero_data%n_source 00202 if (name == aero_data%source_name(i)) then 00203 found = .true. 00204 exit 00205 end if 00206 end do 00207 if (found) then 00208 aero_data_source_by_name = i 00209 else 00210 allocate(new_source_name(aero_data%n_source + 1)) 00211 do i = 1,aero_data%n_source 00212 new_source_name(i) = aero_data%source_name(i) 00213 end do 00214 new_source_name(i) = name(1:AERO_SOURCE_NAME_LEN) 00215 deallocate(aero_data%source_name) 00216 aero_data%source_name => new_source_name 00217 aero_data%n_source = aero_data%n_source + 1 00218 aero_data_source_by_name = aero_data%n_source 00219 end if 00220 00221 end function aero_data_source_by_name 00222 00223 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00224 00225 !> Fills in aero_data%i_water. 00226 subroutine aero_data_set_water_index(aero_data) 00227 00228 !> Aero_data data. 00229 type(aero_data_t), intent(inout) :: aero_data 00230 00231 integer :: i 00232 00233 do i = 1,aero_data%n_spec 00234 if (aero_data%name(i) == "H2O") then 00235 aero_data%i_water = i 00236 end if 00237 end do 00238 00239 end subroutine aero_data_set_water_index 00240 00241 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00242 00243 !> Fills in aero_data%mosaic_index. 00244 subroutine aero_data_set_mosaic_map(aero_data) 00245 00246 !> Aero_data data. 00247 type(aero_data_t), intent(inout) :: aero_data 00248 00249 integer, parameter :: n_mosaic_spec = 19 00250 character(AERO_NAME_LEN), parameter, dimension(n_mosaic_spec) :: 00251 mosaic_spec_name = [ 00252 "SO4 ", "NO3 ", "Cl ", "NH4 ", "MSA ", "ARO1 ", 00253 "ARO2 ", "ALK1 ", "OLE1 ", "API1 ", "API2 ", "LIM1 ", 00254 "LIM2 ", "CO3 ", "Na ", "Ca ", "OIN ", "OC ", 00255 "BC "] 00256 00257 integer :: i_spec, i_mosaic_spec, i 00258 00259 aero_data%mosaic_index = 0 00260 do i_spec = 1,aero_data%n_spec 00261 i_mosaic_spec = 0 00262 do i = 1,n_mosaic_spec 00263 if (aero_data%name(i_spec) == mosaic_spec_name(i)) then 00264 i_mosaic_spec = i 00265 end if 00266 end do 00267 if (i_mosaic_spec > 0) then 00268 aero_data%mosaic_index(i_spec) = i_mosaic_spec 00269 end if 00270 end do 00271 00272 end subroutine aero_data_set_mosaic_map 00273 00274 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00275 00276 !> Read aero_data specification from a spec file. 00277 subroutine spec_file_read_aero_data(file, aero_data) 00278 00279 !> Spec file to read data from. 00280 type(spec_file_t), intent(inout) :: file 00281 !> Aero_data data. 00282 type(aero_data_t), intent(inout) :: aero_data 00283 00284 integer :: n_species, species, i 00285 character(len=SPEC_LINE_MAX_VAR_LEN), pointer :: species_name(:) 00286 real(kind=dp), pointer :: species_data(:,:) 00287 00288 !> \page input_format_aero_data Input File Format: Aerosol Material Data 00289 !! 00290 !! A aerosol material data file must consist of one line per 00291 !! aerosol species, with each line having: 00292 !! - species name (string) 00293 !! - density (real, unit kg/m^3) 00294 !! - ions per fully dissociated molecule (integer) - used to 00295 !! compute kappa value if the corresponding kappa value is 00296 !! zero 00297 !! - molecular weight (real, unit kg/mol) 00298 !! - kappa hygroscopicity parameter (real, dimensionless) - if 00299 !! zero, then inferred from the ions value 00300 !! 00301 !! This specifies both which species are to be recognized as 00302 !! aerosol consituents, as well as their physical properties. For 00303 !! example, an aerosol material data file could contain: 00304 !! <pre> 00305 !! # species dens (kg/m^3) ions (1) molec wght (kg/mole) kappa (1) 00306 !! SO4 1800 0 96e-3 0.65 00307 !! NO3 1800 0 62e-3 0.65 00308 !! Cl 2200 1 35.5e-3 0 00309 !! NH4 1800 0 18e-3 0.65 00310 !! </pre> 00311 !! 00312 !! Note that it is an error to specify a non-zero number of ions 00313 !! and a non-zero kappa value for a species. If both values are 00314 !! zero then that species has zero hygroscopicity parameter. If 00315 !! exactly one of kappa or ions is non-zero then the non-zero 00316 !! value is used and the zero value is ignored. 00317 !! 00318 !! See also: 00319 !! - \ref spec_file_format --- the input file text format 00320 !! - \ref output_format_aero_data --- the corresponding output format 00321 00322 allocate(species_name(0)) 00323 allocate(species_data(0,0)) 00324 call spec_file_read_real_named_array(file, 0, species_name, species_data) 00325 00326 ! check the data size 00327 n_species = size(species_data, 1) 00328 if (.not. ((size(species_data, 2) == 4) .or. (n_species == 0))) then 00329 call die_msg(428926381, 'each line in ' // trim(file%name) & 00330 // ' should contain exactly 5 values') 00331 end if 00332 00333 ! allocate and copy over the data 00334 call aero_data_deallocate(aero_data) 00335 call aero_data_allocate_size(aero_data, n_species, 0) 00336 do i = 1,n_species 00337 aero_data%name(i) = species_name(i)(1:AERO_NAME_LEN) 00338 aero_data%density(i) = species_data(i,1) 00339 aero_data%num_ions(i) = nint(species_data(i,2)) 00340 aero_data%molec_weight(i) = species_data(i,3) 00341 aero_data%kappa(i) = species_data(i,4) 00342 call assert_msg(232362742, & 00343 (aero_data%num_ions(i) == 0) .or. (aero_data%kappa(i) == 0d0), & 00344 "ions and kappa both non-zero for species " & 00345 // trim(aero_data%name(i)) // " in " // trim(file%name)) 00346 if (species_name(i) == "H2O") then 00347 aero_data%i_water = i 00348 call warn_assert_msg(945800387, & 00349 aero_data%density(i) == const%water_density, & 00350 "input H2O density not equal to const%water_density (" & 00351 // trim(real_to_string(aero_data%density(i))) // " /= " & 00352 // trim(real_to_string(const%water_density)) // ")") 00353 call warn_assert_msg(945800387, & 00354 aero_data%molec_weight(i) == const%water_molec_weight, & 00355 "input H2O molec_weight not equal " & 00356 // "to const%water_molec_weight (" & 00357 // trim(real_to_string(aero_data%molec_weight(i))) // " /= " & 00358 // trim(real_to_string(const%water_molec_weight)) // ")") 00359 end if 00360 end do 00361 deallocate(species_name) 00362 deallocate(species_data) 00363 call aero_data_set_water_index(aero_data) 00364 call aero_data_set_mosaic_map(aero_data) 00365 00366 end subroutine spec_file_read_aero_data 00367 00368 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00369 00370 !> Read a list of species from the given file with the given name. 00371 subroutine spec_file_read_species_list(file, name, aero_data, species_list) 00372 00373 !> Spec file. 00374 type(spec_file_t), intent(inout) :: file 00375 !> Name of line. 00376 character(len=*), intent(in) :: name 00377 !> Aero_data data. 00378 type(aero_data_t), intent(in) :: aero_data 00379 !> List of species numbers. 00380 integer, pointer :: species_list(:) 00381 00382 type(spec_line_t) :: line 00383 integer :: i, spec 00384 00385 call spec_file_read_line_no_eof(file, line) 00386 call spec_file_check_line_name(file, line, name) 00387 allocate(species_list(size(line%data))) 00388 do i = 1,size(line%data) 00389 spec = aero_data_spec_by_name(aero_data, line%data(i)) 00390 if (spec == 0) then 00391 call spec_file_die_msg(964771833, file, & 00392 'unknown species: ' // trim(line%data(i))) 00393 end if 00394 species_list(i) = spec 00395 end do 00396 call spec_line_deallocate(line) 00397 00398 end subroutine spec_file_read_species_list 00399 00400 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00401 00402 !> Determines the number of bytes required to pack the given value. 00403 integer function pmc_mpi_pack_size_aero_data(val) 00404 00405 !> Value to pack. 00406 type(aero_data_t), intent(in) :: val 00407 00408 pmc_mpi_pack_size_aero_data = & 00409 pmc_mpi_pack_size_integer(val%n_spec) & 00410 + pmc_mpi_pack_size_integer(val%n_source) & 00411 + pmc_mpi_pack_size_integer(val%i_water) & 00412 + pmc_mpi_pack_size_string_array(val%name) & 00413 + pmc_mpi_pack_size_integer_array(val%mosaic_index) & 00414 + pmc_mpi_pack_size_real_array(val%density) & 00415 + pmc_mpi_pack_size_integer_array(val%num_ions) & 00416 + pmc_mpi_pack_size_real_array(val%molec_weight) & 00417 + pmc_mpi_pack_size_real_array(val%kappa) & 00418 + pmc_mpi_pack_size_string_array(val%source_name) 00419 00420 end function pmc_mpi_pack_size_aero_data 00421 00422 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00423 00424 !> Packs the given value into the buffer, advancing position. 00425 subroutine pmc_mpi_pack_aero_data(buffer, position, val) 00426 00427 !> Memory buffer. 00428 character, intent(inout) :: buffer(:) 00429 !> Current buffer position. 00430 integer, intent(inout) :: position 00431 !> Value to pack. 00432 type(aero_data_t), intent(in) :: val 00433 00434 #ifdef PMC_USE_MPI 00435 integer :: prev_position 00436 00437 prev_position = position 00438 call pmc_mpi_pack_integer(buffer, position, val%n_spec) 00439 call pmc_mpi_pack_integer(buffer, position, val%n_source) 00440 call pmc_mpi_pack_integer(buffer, position, val%i_water) 00441 call pmc_mpi_pack_string_array(buffer, position, val%name) 00442 call pmc_mpi_pack_integer_array(buffer, position, val%mosaic_index) 00443 call pmc_mpi_pack_real_array(buffer, position, val%density) 00444 call pmc_mpi_pack_integer_array(buffer, position, val%num_ions) 00445 call pmc_mpi_pack_real_array(buffer, position, val%molec_weight) 00446 call pmc_mpi_pack_real_array(buffer, position, val%kappa) 00447 call pmc_mpi_pack_string_array(buffer, position, val%source_name) 00448 call assert(183834856, & 00449 position - prev_position <= pmc_mpi_pack_size_aero_data(val)) 00450 #endif 00451 00452 end subroutine pmc_mpi_pack_aero_data 00453 00454 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00455 00456 !> Unpacks the given value from the buffer, advancing position. 00457 subroutine pmc_mpi_unpack_aero_data(buffer, position, val) 00458 00459 !> Memory buffer. 00460 character, intent(inout) :: buffer(:) 00461 !> Current buffer position. 00462 integer, intent(inout) :: position 00463 !> Value to pack. 00464 type(aero_data_t), intent(inout) :: val 00465 00466 #ifdef PMC_USE_MPI 00467 integer :: prev_position 00468 00469 prev_position = position 00470 call pmc_mpi_unpack_integer(buffer, position, val%n_spec) 00471 call pmc_mpi_unpack_integer(buffer, position, val%n_source) 00472 call pmc_mpi_unpack_integer(buffer, position, val%i_water) 00473 call pmc_mpi_unpack_string_array(buffer, position, val%name) 00474 call pmc_mpi_unpack_integer_array(buffer, position, val%mosaic_index) 00475 call pmc_mpi_unpack_real_array(buffer, position, val%density) 00476 call pmc_mpi_unpack_integer_array(buffer, position, val%num_ions) 00477 call pmc_mpi_unpack_real_array(buffer, position, val%molec_weight) 00478 call pmc_mpi_unpack_real_array(buffer, position, val%kappa) 00479 call pmc_mpi_unpack_string_array(buffer, position, val%source_name) 00480 call assert(188522823, & 00481 position - prev_position <= pmc_mpi_pack_size_aero_data(val)) 00482 #endif 00483 00484 end subroutine pmc_mpi_unpack_aero_data 00485 00486 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00487 00488 !> Write the aero species dimension to the given NetCDF file if it 00489 !> is not already present and in any case return the associated 00490 !> dimid. 00491 subroutine aero_data_netcdf_dim_aero_species(aero_data, ncid, & 00492 dimid_aero_species) 00493 00494 !> Aero_data structure. 00495 type(aero_data_t), intent(in) :: aero_data 00496 !> NetCDF file ID, in data mode. 00497 integer, intent(in) :: ncid 00498 !> Dimid of the species dimension. 00499 integer, intent(out) :: dimid_aero_species 00500 00501 integer :: status, i_spec 00502 integer :: varid_aero_species 00503 integer :: aero_species_centers(aero_data%n_spec) 00504 character(len=(AERO_NAME_LEN * aero_data%n_spec)) :: aero_species_names 00505 00506 ! try to get the dimension ID 00507 status = nf90_inq_dimid(ncid, "aero_species", dimid_aero_species) 00508 if (status == NF90_NOERR) return 00509 if (status /= NF90_EBADDIM) call pmc_nc_check(status) 00510 00511 ! dimension not defined, so define now define it 00512 call pmc_nc_check(nf90_redef(ncid)) 00513 00514 call pmc_nc_check(nf90_def_dim(ncid, "aero_species", & 00515 aero_data%n_spec, dimid_aero_species)) 00516 aero_species_names = "" 00517 do i_spec = 1,aero_data%n_spec 00518 aero_species_names((len_trim(aero_species_names) + 1):) & 00519 = trim(aero_data%name(i_spec)) 00520 if (i_spec < aero_data%n_spec) then 00521 aero_species_names((len_trim(aero_species_names) + 1):) = "," 00522 end if 00523 end do 00524 call pmc_nc_check(nf90_def_var(ncid, "aero_species", NF90_INT, & 00525 dimid_aero_species, varid_aero_species)) 00526 call pmc_nc_check(nf90_put_att(ncid, varid_aero_species, "names", & 00527 aero_species_names)) 00528 call pmc_nc_check(nf90_put_att(ncid, varid_aero_species, "description", & 00529 "dummy dimension variable (no useful value) - read species names " & 00530 // "as comma-separated values from the 'names' attribute")) 00531 00532 call pmc_nc_check(nf90_enddef(ncid)) 00533 00534 do i_spec = 1,aero_data%n_spec 00535 aero_species_centers(i_spec) = i_spec 00536 end do 00537 call pmc_nc_check(nf90_put_var(ncid, varid_aero_species, & 00538 aero_species_centers)) 00539 00540 end subroutine aero_data_netcdf_dim_aero_species 00541 00542 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00543 00544 !> Write the aero source dimension to the given NetCDF file if it 00545 !> is not already present and in any case return the associated 00546 !> dimid. 00547 subroutine aero_data_netcdf_dim_aero_source(aero_data, ncid, & 00548 dimid_aero_source) 00549 00550 !> Aero_data structure. 00551 type(aero_data_t), intent(in) :: aero_data 00552 !> NetCDF file ID, in data mode. 00553 integer, intent(in) :: ncid 00554 !> Dimid of the source dimension. 00555 integer, intent(out) :: dimid_aero_source 00556 00557 integer :: status, i_source 00558 integer :: varid_aero_source 00559 integer :: aero_source_centers(aero_data%n_source) 00560 character(len=(AERO_SOURCE_NAME_LEN * aero_data%n_source)) 00561 :: aero_source_names 00562 00563 ! try to get the dimension ID 00564 status = nf90_inq_dimid(ncid, "aero_source", dimid_aero_source) 00565 if (status == NF90_NOERR) return 00566 if (status /= NF90_EBADDIM) call pmc_nc_check(status) 00567 00568 ! dimension not defined, so define now define it 00569 call pmc_nc_check(nf90_redef(ncid)) 00570 00571 call pmc_nc_check(nf90_def_dim(ncid, "aero_source", & 00572 aero_data%n_source, dimid_aero_source)) 00573 aero_source_names = "" 00574 do i_source = 1,aero_data%n_source 00575 aero_source_names((len_trim(aero_source_names) + 1):) & 00576 = trim(aero_data%source_name(i_source)) 00577 if (i_source < aero_data%n_source) then 00578 aero_source_names((len_trim(aero_source_names) + 1):) = "," 00579 end if 00580 end do 00581 call pmc_nc_check(nf90_def_var(ncid, "aero_source", NF90_INT, & 00582 dimid_aero_source, varid_aero_source)) 00583 call pmc_nc_check(nf90_put_att(ncid, varid_aero_source, "names", & 00584 aero_source_names)) 00585 call pmc_nc_check(nf90_put_att(ncid, varid_aero_source, "description", & 00586 "dummy dimension variable (no useful value) - read source names " & 00587 // "as comma-separated values from the 'names' attribute")) 00588 00589 call pmc_nc_check(nf90_enddef(ncid)) 00590 00591 do i_source = 1,aero_data%n_source 00592 aero_source_centers(i_source) = i_source 00593 end do 00594 call pmc_nc_check(nf90_put_var(ncid, varid_aero_source, & 00595 aero_source_centers)) 00596 00597 end subroutine aero_data_netcdf_dim_aero_source 00598 00599 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00600 00601 !> Write full state. 00602 subroutine aero_data_output_netcdf(aero_data, ncid) 00603 00604 !> Aero_data to write. 00605 type(aero_data_t), intent(in) :: aero_data 00606 !> NetCDF file ID, in data mode. 00607 integer, intent(in) :: ncid 00608 00609 integer :: dimid_aero_species, dimid_aero_source 00610 00611 !> \page output_format_aero_data Output File Format: Aerosol Material Data 00612 !! 00613 !! The aerosol material data NetCDF dimensions are: 00614 !! - \b aero_species: number of aerosol species 00615 !! - \b aero_source: number of aerosol sources 00616 !! 00617 !! The aerosol material data NetCDF variables are: 00618 !! - \b aero_species (dim \c aero_species): dummy dimension variable 00619 !! (no useful value) - read species names as comma-separated values 00620 !! from the 'names' attribute 00621 !! - \b aero_source (dim \c aero_source): dummy dimension variable 00622 !! (no useful value) - read source names as comma-separated values 00623 !! from the 'names' attribute 00624 !! - \b aero_mosaic_index (dim \c aero_species): indices of species 00625 !! in MOSAIC 00626 !! - \b aero_density (unit kg/m^3, dim \c aero_species): densities 00627 !! of aerosol species 00628 !! - \b aero_num_ions (dim \c aero_species): number of ions produced 00629 !! when one molecule of each species fully dissociates in water 00630 !! - \b aero_molec_weight (unit kg/mol, dim \c aero_species): molecular 00631 !! weights of aerosol species 00632 !! - \b aero_kappa (unit kg/mol, dim \c aero_species): hygroscopicity 00633 !! parameters of aerosol species 00634 !! 00635 !! See also: 00636 !! - \ref input_format_aero_data --- the corresponding input format 00637 00638 call aero_data_netcdf_dim_aero_species(aero_data, ncid, & 00639 dimid_aero_species) 00640 call aero_data_netcdf_dim_aero_source(aero_data, ncid, & 00641 dimid_aero_source) 00642 00643 call pmc_nc_write_integer_1d(ncid, aero_data%mosaic_index, & 00644 "aero_mosaic_index", (/ dimid_aero_species /), & 00645 long_name="MOSAIC indices of aerosol species") 00646 call pmc_nc_write_real_1d(ncid, aero_data%density, & 00647 "aero_density", (/ dimid_aero_species /), unit="kg/m^3", & 00648 long_name="densities of aerosol species") 00649 call pmc_nc_write_integer_1d(ncid, aero_data%num_ions, & 00650 "aero_num_ions", (/ dimid_aero_species /), & 00651 long_name="number of ions after dissociation of aerosol species") 00652 call pmc_nc_write_real_1d(ncid, aero_data%molec_weight, & 00653 "aero_molec_weight", (/ dimid_aero_species /), unit="kg/mol", & 00654 long_name="molecular weights of aerosol species") 00655 call pmc_nc_write_real_1d(ncid, aero_data%kappa, & 00656 "aero_kappa", (/ dimid_aero_species /), unit="1", & 00657 long_name="hygroscopicity parameters (kappas) of aerosol species") 00658 00659 end subroutine aero_data_output_netcdf 00660 00661 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00662 00663 !> Read full state. 00664 subroutine aero_data_input_netcdf(aero_data, ncid) 00665 00666 !> Aero_data to read. 00667 type(aero_data_t), intent(inout) :: aero_data 00668 !> NetCDF file ID, in data mode. 00669 integer, intent(in) :: ncid 00670 00671 character(len=1000) :: name 00672 integer :: dimid_aero_species, n_spec, varid_aero_species, i_spec, i 00673 integer :: dimid_aero_source, n_source, varid_aero_source, i_source 00674 character(len=((AERO_NAME_LEN + 2) * 1000)) :: aero_species_names 00675 character(len=((AERO_SOURCE_NAME_LEN + 2) * 1000)) :: aero_source_names 00676 00677 call pmc_nc_check(nf90_inq_dimid(ncid, "aero_species", & 00678 dimid_aero_species)) 00679 call pmc_nc_check(nf90_Inquire_Dimension(ncid, & 00680 dimid_aero_species, name, n_spec)) 00681 call assert(739238793, n_spec < 1000) 00682 00683 call pmc_nc_check(nf90_inq_dimid(ncid, "aero_source", & 00684 dimid_aero_source)) 00685 call pmc_nc_check(nf90_Inquire_Dimension(ncid, & 00686 dimid_aero_source, name, n_source)) 00687 call assert(739238793, n_source < 1000) 00688 00689 call aero_data_deallocate(aero_data) 00690 call aero_data_allocate_size(aero_data, n_spec, n_source) 00691 00692 call pmc_nc_read_integer_1d(ncid, aero_data%mosaic_index, & 00693 "aero_mosaic_index") 00694 call pmc_nc_read_real_1d(ncid, aero_data%density, & 00695 "aero_density") 00696 call pmc_nc_read_integer_1d(ncid, aero_data%num_ions, & 00697 "aero_num_ions") 00698 call pmc_nc_read_real_1d(ncid, aero_data%molec_weight, & 00699 "aero_molec_weight") 00700 call pmc_nc_read_real_1d(ncid, aero_data%kappa, & 00701 "aero_kappa") 00702 00703 call pmc_nc_check(nf90_inq_varid(ncid, "aero_species", & 00704 varid_aero_species)) 00705 call pmc_nc_check(nf90_get_att(ncid, varid_aero_species, "names", & 00706 aero_species_names)) 00707 ! aero_species_names are comma-separated, so unpack them 00708 do i_spec = 1,aero_data%n_spec 00709 i = 1 00710 do while ((aero_species_names(i:i) /= " ") & 00711 .and. (aero_species_names(i:i) /= ",")) 00712 i = i + 1 00713 end do 00714 call assert(852937292, i > 1) 00715 aero_data%name(i_spec) = aero_species_names(1:(i-1)) 00716 aero_species_names = aero_species_names((i+1):) 00717 end do 00718 call assert(729138192, aero_species_names == "") 00719 00720 call pmc_nc_check(nf90_inq_varid(ncid, "aero_source", & 00721 varid_aero_source)) 00722 call pmc_nc_check(nf90_get_att(ncid, varid_aero_source, "names", & 00723 aero_source_names)) 00724 ! aero_source_names are comma-separated, so unpack them 00725 do i_source = 1,aero_data%n_source 00726 i = 1 00727 do while ((aero_source_names(i:i) /= " ") & 00728 .and. (aero_source_names(i:i) /= ",")) 00729 i = i + 1 00730 end do 00731 call assert(840982478, i > 1) 00732 aero_data%source_name(i_source) = aero_source_names(1:(i-1)) 00733 aero_source_names = aero_source_names((i+1):) 00734 end do 00735 call assert(377166446, aero_source_names == "") 00736 00737 call aero_data_set_water_index(aero_data) 00738 00739 end subroutine aero_data_input_netcdf 00740 00741 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 00742 00743 end module pmc_aero_data