31 character(len=GAS_NAME_LEN),
allocatable :: name(:)
35 integer,
allocatable :: mosaic_index(:)
48 if (
allocated(gas_data%name))
then 65 character(len=*),
intent(in) :: name
72 if (name == gas_data%name(i))
then 93 integer,
parameter :: n_mosaic_species = 77
94 character(len=GAS_NAME_LEN),
parameter,
dimension(n_mosaic_species) &
95 :: mosaic_species = [ &
96 "H2SO4 ",
"HNO3 ",
"HCl ", &
97 "NH3 ",
"NO ",
"NO2 ", &
98 "NO3 ",
"N2O5 ",
"HONO ", &
99 "HNO4 ",
"O3 ",
"O1D ", &
100 "O3P ",
"OH ",
"HO2 ", &
101 "H2O2 ",
"CO ",
"SO2 ", &
102 "CH4 ",
"C2H6 ",
"CH3O2 ", &
103 "ETHP ",
"HCHO ",
"CH3OH ", &
104 "ANOL ",
"CH3OOH ",
"ETHOOH ", &
105 "ALD2 ",
"HCOOH ",
"RCOOH ", &
106 "C2O3 ",
"PAN ",
"ARO1 ", &
107 "ARO2 ",
"ALK1 ",
"OLE1 ", &
108 "API1 ",
"API2 ",
"LIM1 ", &
109 "LIM2 ",
"PAR ",
"AONE ", &
110 "MGLY ",
"ETH ",
"OLET ", &
111 "OLEI ",
"TOL ",
"XYL ", &
112 "CRES ",
"TO2 ",
"CRO ", &
113 "OPEN ",
"ONIT ",
"ROOH ", &
114 "RO2 ",
"ANO2 ",
"NAP ", &
115 "XO2 ",
"XPAR ",
"ISOP ", &
116 "ISOPRD ",
"ISOPP ",
"ISOPN ", &
117 "ISOPO2 ",
"API ",
"LIM ", &
118 "DMS ",
"MSA ",
"DMSO ", &
119 "DMSO2 ",
"CH3SO2H ",
"CH3SCH2OO ", &
120 "CH3SO2 ",
"CH3SO3 ",
"CH3SO2OO ", &
121 "CH3SO2CH2OO ",
"SULFHOX "]
123 integer spec, mosaic_spec, i
125 gas_data%mosaic_index = 0
128 do i = 1,n_mosaic_species
129 if (gas_data%name(spec) == mosaic_species(i))
then 133 if (mosaic_spec > 0)
then 134 gas_data%mosaic_index(spec) = mosaic_spec
143 subroutine spec_file_read_gas_data(file, gas_data)
150 integer :: n_species, species, i
151 character(len=SPEC_LINE_MAX_VAR_LEN),
allocatable :: species_name(:)
152 real(kind=dp),
allocatable :: species_data(:,:)
176 if (
size(species_data, 2) /= 0)
then 177 call die_msg(614290516,
'each line in ' // trim(file%name) &
178 //
' must only contain the species name')
182 n_species =
size(species_data, 1)
191 end subroutine spec_file_read_gas_data
213 character,
intent(inout) :: buffer(:)
215 integer,
intent(inout) :: position
220 integer :: prev_position
222 prev_position = position
237 character,
intent(inout) :: buffer(:)
239 integer,
intent(inout) :: position
244 integer :: prev_position
246 prev_position = position
266 integer,
intent(in) :: ncid
268 integer,
intent(out) :: dimid_gas_species
270 integer :: status, i_spec
271 integer :: varid_gas_species
273 character(len=((GAS_NAME_LEN + 2) * gas_data_n_spec(gas_data))) &
277 status = nf90_inq_dimid(ncid,
"gas_species", dimid_gas_species)
278 if (status == nf90_noerr)
return 286 gas_species_names =
"" 288 gas_species_names((len_trim(gas_species_names) + 1):) &
289 = trim(gas_data%name(i_spec))
291 gas_species_names((len_trim(gas_species_names) + 1):) =
"," 294 call pmc_nc_check(nf90_def_var(ncid,
"gas_species", nf90_int, &
295 dimid_gas_species, varid_gas_species))
296 call pmc_nc_check(nf90_put_att(ncid, varid_gas_species,
"names", &
298 call pmc_nc_check(nf90_put_att(ncid, varid_gas_species,
"description", &
299 "dummy dimension variable (no useful value) - read species names " &
300 //
"as comma-separated values from the 'names' attribute"))
305 gas_species_centers(i_spec) = i_spec
307 call pmc_nc_check(nf90_put_var(ncid, varid_gas_species, &
308 gas_species_centers))
315 subroutine gas_data_output_netcdf(gas_data, ncid)
320 integer,
intent(in) :: ncid
322 integer :: dimid_gas_species
343 "gas_mosaic_index", (/ dimid_gas_species /), &
344 long_name=
"MOSAIC indices of gas species")
346 end subroutine gas_data_output_netcdf
356 integer,
intent(in) :: ncid
358 character(len=1000) :: name
359 integer :: dimid_gas_species, n_spec, varid_gas_species, i_spec, i
360 character(len=((GAS_NAME_LEN + 2) * 1000)) :: gas_species_names
362 call pmc_nc_check(nf90_inq_dimid(ncid,
"gas_species", dimid_gas_species))
363 call pmc_nc_check(nf90_inquire_dimension(ncid, dimid_gas_species, name, &
365 call assert(719237193, n_spec < 1000)
367 if (
allocated(gas_data%name))
deallocate(gas_data%name)
368 if (
allocated(gas_data%mosaic_index))
deallocate(gas_data%mosaic_index)
369 allocate(gas_data%name(n_spec))
370 allocate(gas_data%mosaic_index(n_spec))
375 call pmc_nc_check(nf90_inq_varid(ncid,
"gas_species", varid_gas_species))
376 call pmc_nc_check(nf90_get_att(ncid, varid_gas_species,
"names", &
381 do while ((gas_species_names(i:i) /=
" ") &
382 .and. (gas_species_names(i:i) /=
","))
385 call assert(173021381, i > 1)
386 gas_data%name(i_spec) = gas_species_names(1:(i-1))
387 gas_species_names = gas_species_names((i+1):)
389 call assert(469721220, gas_species_names ==
"")
subroutine pmc_mpi_pack_integer_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
An input file with extra data for printing messages.
integer function pmc_mpi_pack_size_integer_array(val)
Determines the number of bytes required to pack the given value.
Wrapper functions for NetCDF. These all take a NetCDF ncid in data mode and return with it again in d...
subroutine pmc_mpi_unpack_integer_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine ensure_integer_array_size(x, n, only_grow)
Allocate or reallocate the given array to ensure it is of the given size, preserving any data and/or ...
integer, parameter gas_name_len
Maximum length of the name of a gas.
subroutine pmc_mpi_pack_string_array(buffer, position, val)
Packs the given value into the buffer, advancing position.
The gas_data_t structure and associated subroutines.
subroutine assert(code, condition_ok)
Errors unless condition_ok is true.
subroutine pmc_mpi_unpack_gas_data(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
subroutine pmc_nc_check(status)
Check the status of a NetCDF function call.
subroutine gas_data_input_netcdf(gas_data, ncid)
Read full state.
subroutine die_msg(code, error_msg)
Error immediately.
subroutine gas_data_set_mosaic_map(gas_data)
Fills in gas_data%mosaic_index.
subroutine ensure_string_array_size(x, n)
Allocate or reallocate the given array to ensure it is of the given size, without preserving data...
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...
elemental integer function gas_data_n_spec(gas_data)
Return the number of gas species.
subroutine pmc_mpi_pack_gas_data(buffer, position, val)
Packs the given value into the buffer, advancing position.
Reading formatted text input.
subroutine pmc_mpi_unpack_string_array(buffer, position, val)
Unpacks the given value from the buffer, advancing position.
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...
integer function pmc_mpi_pack_size_gas_data(val)
Determines the number of bytes required to pack the given value.
subroutine pmc_nc_read_integer_1d(ncid, var, name, must_be_present)
Read a simple integer array from a NetCDF file.
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...
Common utility subroutines.
Wrapper functions for MPI.
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.
integer function pmc_mpi_pack_size_string_array(val)
Determines the number of bytes required to pack the given value.