20 integer,
parameter :: GAS_NAME_LEN = 100
33 character(len=GAS_NAME_LEN),
pointer :: name(:)
37 integer,
pointer :: mosaic_index(:)
51 allocate(gas_data%name(0))
52 allocate(gas_data%mosaic_index(0))
64 integer,
intent(in) :: n_spec
66 gas_data%n_spec = n_spec
67 allocate(gas_data%name(n_spec))
68 allocate(gas_data%mosaic_index(n_spec))
80 deallocate(gas_data%name)
81 deallocate(gas_data%mosaic_index)
94 character(len=*),
intent(in) :: name
100 do i = 1,gas_data%n_spec
101 if (name == gas_data%name(i))
then
122 integer,
parameter :: n_mosaic_species = 77
123 character(len=GAS_NAME_LEN),
parameter,
dimension(n_mosaic_species) &
124 :: mosaic_species = [ &
125 "H2SO4 ",
"HNO3 ",
"HCl ", &
126 "NH3 ",
"NO ",
"NO2 ", &
127 "NO3 ",
"N2O5 ",
"HONO ", &
128 "HNO4 ",
"O3 ",
"O1D ", &
129 "O3P ",
"OH ",
"HO2 ", &
130 "H2O2 ",
"CO ",
"SO2 ", &
131 "CH4 ",
"C2H6 ",
"CH3O2 ", &
132 "ETHP ",
"HCHO ",
"CH3OH ", &
133 "ANOL ",
"CH3OOH ",
"ETHOOH ", &
134 "ALD2 ",
"HCOOH ",
"RCOOH ", &
135 "C2O3 ",
"PAN ",
"ARO1 ", &
136 "ARO2 ",
"ALK1 ",
"OLE1 ", &
137 "API1 ",
"API2 ",
"LIM1 ", &
138 "LIM2 ",
"PAR ",
"AONE ", &
139 "MGLY ",
"ETH ",
"OLET ", &
140 "OLEI ",
"TOL ",
"XYL ", &
141 "CRES ",
"TO2 ",
"CRO ", &
142 "OPEN ",
"ONIT ",
"ROOH ", &
143 "RO2 ",
"ANO2 ",
"NAP ", &
144 "XO2 ",
"XPAR ",
"ISOP ", &
145 "ISOPRD ",
"ISOPP ",
"ISOPN ", &
146 "ISOPO2 ",
"API ",
"LIM ", &
147 "DMS ",
"MSA ",
"DMSO ", &
148 "DMSO2 ",
"CH3SO2H ",
"CH3SCH2OO ", &
149 "CH3SO2 ",
"CH3SO3 ",
"CH3SO2OO ", &
150 "CH3SO2CH2OO ",
"SULFHOX "]
152 integer spec, mosaic_spec, i
154 gas_data%mosaic_index = 0
155 do spec = 1,gas_data%n_spec
157 do i = 1,n_mosaic_species
158 if (gas_data%name(spec) == mosaic_species(i))
then
162 if (mosaic_spec > 0)
then
163 gas_data%mosaic_index(spec) = mosaic_spec
172 subroutine spec_file_read_gas_data(file, gas_data)
179 integer :: n_species, species, i
180 character(len=SPEC_LINE_MAX_VAR_LEN),
pointer :: species_name(:)
181 real(kind=dp),
pointer :: species_data(:,:)
201 allocate(species_name(0))
202 allocate(species_data(0,0))
207 if (
size(species_data, 2) /= 0)
then
208 call
die_msg(614290516,
'each line in ' // trim(file%name) &
209 //
' must only contain the species name')
213 n_species =
size(species_data, 1)
217 gas_data%name(i) = species_name(i)(1:gas_name_len)
219 deallocate(species_name)
220 deallocate(species_data)
224 end subroutine spec_file_read_gas_data
247 character,
intent(inout) :: buffer(:)
249 integer,
intent(inout) :: position
254 integer :: prev_position
256 prev_position = position
272 character,
intent(inout) :: buffer(:)
274 integer,
intent(inout) :: position
279 integer :: prev_position
281 prev_position = position
302 integer,
intent(in) :: ncid
304 integer,
intent(out) :: dimid_gas_species
306 integer :: status, i_spec
307 integer :: varid_gas_species
308 integer :: gas_species_centers(gas_data%n_spec)
309 character(len=((GAS_NAME_LEN + 2) * gas_data%n_spec)) :: gas_species_names
312 status = nf90_inq_dimid(ncid,
"gas_species", dimid_gas_species)
313 if (status == nf90_noerr)
return
320 gas_data%n_spec, dimid_gas_species))
321 gas_species_names =
""
322 do i_spec = 1,gas_data%n_spec
323 gas_species_names((len_trim(gas_species_names) + 1):) &
324 = trim(gas_data%name(i_spec))
325 if (i_spec < gas_data%n_spec)
then
326 gas_species_names((len_trim(gas_species_names) + 1):) =
","
329 call
pmc_nc_check(nf90_def_var(ncid,
"gas_species", nf90_int, &
330 dimid_gas_species, varid_gas_species))
331 call
pmc_nc_check(nf90_put_att(ncid, varid_gas_species,
"names", &
333 call
pmc_nc_check(nf90_put_att(ncid, varid_gas_species,
"description", &
334 "dummy dimension variable (no useful value) - read species names " &
335 //
"as comma-separated values from the 'names' attribute"))
339 do i_spec = 1,gas_data%n_spec
340 gas_species_centers(i_spec) = i_spec
342 call
pmc_nc_check(nf90_put_var(ncid, varid_gas_species, &
343 gas_species_centers))
350 subroutine gas_data_output_netcdf(gas_data, ncid)
355 integer,
intent(in) :: ncid
357 integer :: dimid_gas_species
378 "gas_mosaic_index", (/ dimid_gas_species /), &
379 long_name=
"MOSAIC indices of gas species")
381 end subroutine gas_data_output_netcdf
391 integer,
intent(in) :: ncid
393 character(len=1000) :: name
394 integer :: dimid_gas_species, n_spec, varid_gas_species, i_spec, i
395 character(len=((GAS_NAME_LEN + 2) * 1000)) :: gas_species_names
397 call
pmc_nc_check(nf90_inq_dimid(ncid,
"gas_species", dimid_gas_species))
398 call
pmc_nc_check(nf90_inquire_dimension(ncid, dimid_gas_species, name, &
402 call
assert(719237193, n_spec < 1000)
407 call
pmc_nc_check(nf90_inq_varid(ncid,
"gas_species", varid_gas_species))
408 call
pmc_nc_check(nf90_get_att(ncid, varid_gas_species,
"names", &
411 do i_spec = 1,gas_data%n_spec
413 do while ((gas_species_names(i:i) /=
" ") &
414 .and. (gas_species_names(i:i) /=
","))
417 call
assert(173021381, i > 1)
418 gas_data%name(i_spec) = gas_species_names(1:(i-1))
419 gas_species_names = gas_species_names((i+1):)
421 call
assert(469721220, gas_species_names ==
"")
subroutine pmc_mpi_unpack_gas_data(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
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 pmc_mpi_pack_gas_data(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine gas_data_allocate_size(gas_data, n_spec)
Allocate storage for gas species with the given size.
subroutine gas_data_input_netcdf(gas_data, ncid)
Read full state.
subroutine pmc_nc_check(status)
Check the status of a NetCDF function call.
subroutine gas_data_set_mosaic_map(gas_data)
Fills in gas_data%mosaic_index.
subroutine pmc_mpi_pack_string_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
subroutine pmc_nc_read_integer_1d(ncid, var, name, must_be_present)
Read a simple integer array from a NetCDF file.
subroutine gas_data_deallocate(gas_data)
Free all storage.
subroutine pmc_mpi_pack_integer(buffer, position, val)
Packs the given value into the buffer, advancing position.
The gas_data_t structure and associated subroutines.
subroutine pmc_nc_write_integer_1d(ncid, var, name, dimids, dim_name, unit, long_name, standard_name, description)
Write a simple integer array to a NetCDF file.
subroutine spec_file_read_real_named_array(file, max_lines, names, vals)
Read an array of named lines with real data. All lines must have the same number of data elements...
Common utility subroutines.
subroutine pmc_mpi_unpack_integer(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
Wrapper functions for MPI.
subroutine pmc_mpi_unpack_integer_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
integer function pmc_mpi_pack_size_string_array(val)
Determines the number of bytes required to pack the given value.
integer function gas_data_spec_by_name(gas_data, name)
Returns the number of the species in gas with the given name, or returns 0 if there is no such specie...
Reading formatted text input.
subroutine gas_data_allocate(gas_data)
Allocate storage for gas species.
Wrapper functions for NetCDF. These all take a NetCDF ncid in data mode and return with it again in d...
subroutine pmc_mpi_unpack_string_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
integer function pmc_mpi_pack_size_integer_array(val)
Determines the number of bytes required to pack the given value.
subroutine gas_data_netcdf_dim_gas_species(gas_data, ncid, dimid_gas_species)
Write the gas species dimension to the given NetCDF file if it is not already present and in any case...
integer function pmc_mpi_pack_size_gas_data(val)
Determines the number of bytes required to pack the given value.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
subroutine pmc_mpi_pack_integer_array(buffer, position, val)
Packs the given value into the buffer, advancing position.