17 use camp_chem_spec_data
19 use camp_util,
only: string_t, split_string, split_char
37 character(len=GAS_NAME_LEN),
allocatable :: name(:)
41 integer,
allocatable :: mosaic_index(:)
43 integer :: i_camp_water = 0
52 subroutine gas_data_initialize(gas_data, camp_core)
57 class(camp_core_t),
intent(in) :: camp_core
59 type(chem_spec_data_t),
pointer :: chem_spec_data
61 type(string_t),
allocatable :: gas_spec_names(:)
62 type(property_t),
pointer :: property_set
63 character(len=:),
allocatable :: prop_name
68 camp_core%get_chem_spec_data(chem_spec_data), &
69 "No chemical species data in camp core.")
72 gas_spec_names = chem_spec_data%get_spec_names( &
73 spec_phase = chem_spec_gas_phase)
76 allocate(gas_data%name(
size(gas_spec_names)))
79 prop_name =
"is gas-phase water"
80 do i_spec = 1,
size(gas_spec_names)
81 gas_data%name(i_spec) = gas_spec_names(i_spec)%string
83 chem_spec_data%get_property_set( &
84 gas_spec_names(i_spec)%string, &
86 "Missing property set for gas species "// &
87 gas_spec_names(i_spec)%string)
88 if (property_set%get_logical(prop_name, bool_val))
then
89 call assert_msg(423633615, gas_data%i_camp_water == 0, &
90 "More than one gas-phase water species specified")
91 gas_data%i_camp_water = i_spec
95 call assert_msg(134440820, gas_data%i_camp_water /= 0, &
96 "No gas-phase water species specified.")
99 allocate(gas_data%mosaic_index(
size(gas_spec_names)))
100 gas_data%mosaic_index(:) = 0
102 end subroutine gas_data_initialize
113 if (
allocated(gas_data%name))
then
130 character(len=*),
intent(in) :: name
137 if (name == gas_data%name(i))
then
158 integer,
parameter :: n_mosaic_species = 77
159 character(len=GAS_NAME_LEN),
parameter,
dimension(n_mosaic_species) &
160 :: mosaic_species = [ &
161 "H2SO4 ",
"HNO3 ",
"HCl ", &
162 "NH3 ",
"NO ",
"NO2 ", &
163 "NO3 ",
"N2O5 ",
"HONO ", &
164 "HNO4 ",
"O3 ",
"O1D ", &
165 "O3P ",
"OH ",
"HO2 ", &
166 "H2O2 ",
"CO ",
"SO2 ", &
167 "CH4 ",
"C2H6 ",
"CH3O2 ", &
168 "ETHP ",
"HCHO ",
"CH3OH ", &
169 "ANOL ",
"CH3OOH ",
"ETHOOH ", &
170 "ALD2 ",
"HCOOH ",
"RCOOH ", &
171 "C2O3 ",
"PAN ",
"ARO1 ", &
172 "ARO2 ",
"ALK1 ",
"OLE1 ", &
173 "API1 ",
"API2 ",
"LIM1 ", &
174 "LIM2 ",
"PAR ",
"AONE ", &
175 "MGLY ",
"ETH ",
"OLET ", &
176 "OLEI ",
"TOL ",
"XYL ", &
177 "CRES ",
"TO2 ",
"CRO ", &
178 "OPEN ",
"ONIT ",
"ROOH ", &
179 "RO2 ",
"ANO2 ",
"NAP ", &
180 "XO2 ",
"XPAR ",
"ISOP ", &
181 "ISOPRD ",
"ISOPP ",
"ISOPN ", &
182 "ISOPO2 ",
"API ",
"LIM ", &
183 "DMS ",
"MSA ",
"DMSO ", &
184 "DMSO2 ",
"CH3SO2H ",
"CH3SCH2OO ", &
185 "CH3SO2 ",
"CH3SO3 ",
"CH3SO2OO ", &
186 "CH3SO2CH2OO ",
"SULFHOX "]
188 integer spec, mosaic_spec, i
190 gas_data%mosaic_index = 0
193 do i = 1,n_mosaic_species
194 if (gas_data%name(spec) == mosaic_species(i))
then
198 if (mosaic_spec > 0)
then
199 gas_data%mosaic_index(spec) = mosaic_spec
208 subroutine spec_file_read_gas_data(file, gas_data)
215 integer :: n_species, species, i
216 character(len=SPEC_LINE_MAX_VAR_LEN),
allocatable :: species_name(:)
217 real(kind=
dp),
allocatable :: species_data(:,:)
241 if (
size(species_data, 2) /= 0)
then
242 call die_msg(614290516,
'each line in ' // trim(file%name) &
243 //
' must only contain the species name')
247 n_species =
size(species_data, 1)
256 end subroutine spec_file_read_gas_data
278 character,
intent(inout) :: buffer(:)
280 integer,
intent(inout) :: position
285 integer :: prev_position
287 prev_position = position
302 character,
intent(inout) :: buffer(:)
304 integer,
intent(inout) :: position
309 integer :: prev_position
311 prev_position = position
331 integer,
intent(in) :: ncid
333 integer,
intent(out) :: dimid_gas_species
335 integer :: status, i_spec
336 integer :: varid_gas_species
337 integer :: gas_species_centers(gas_data_n_spec(gas_data))
338 character(len=((GAS_NAME_LEN + 2) * gas_data_n_spec(gas_data))) &
342 status = nf90_inq_dimid(ncid,
"gas_species", dimid_gas_species)
343 if (status == nf90_noerr)
return
350 gas_data_n_spec(gas_data), dimid_gas_species))
351 gas_species_names =
""
352 do i_spec = 1,gas_data_n_spec(gas_data)
353 gas_species_names((len_trim(gas_species_names) + 1):) &
354 = trim(gas_data%name(i_spec))
355 if (i_spec < gas_data_n_spec(gas_data))
then
356 gas_species_names((len_trim(gas_species_names) + 1):) =
","
359 call pmc_nc_check(nf90_def_var(ncid,
"gas_species", nf90_int, &
360 dimid_gas_species, varid_gas_species))
361 call pmc_nc_check(nf90_put_att(ncid, varid_gas_species,
"names", &
363 call pmc_nc_check(nf90_put_att(ncid, varid_gas_species,
"description", &
364 "dummy dimension variable (no useful value) - read species names " &
365 //
"as comma-separated values from the 'names' attribute"))
369 do i_spec = 1,gas_data_n_spec(gas_data)
370 gas_species_centers(i_spec) = i_spec
372 call pmc_nc_check(nf90_put_var(ncid, varid_gas_species, &
373 gas_species_centers))
380 subroutine gas_data_output_netcdf(gas_data, ncid)
385 integer,
intent(in) :: ncid
387 integer :: dimid_gas_species
408 "gas_mosaic_index", (/ dimid_gas_species /), &
409 long_name=
"MOSAIC indices of gas species")
411 end subroutine gas_data_output_netcdf
421 integer,
intent(in) :: ncid
423 character(len=1000) :: name
424 integer :: dimid_gas_species, n_spec, varid_gas_species, i_spec, i
425 character(len=:),
allocatable :: gas_species_names
427 call pmc_nc_check(nf90_inq_dimid(ncid,
"gas_species", dimid_gas_species))
428 call pmc_nc_check(nf90_inquire_dimension(ncid, dimid_gas_species, name, &
430 call assert(719237193, n_spec < 1000)
432 if (
allocated(gas_data%name))
deallocate(gas_data%name)
433 if (
allocated(gas_data%mosaic_index))
deallocate(gas_data%mosaic_index)
434 allocate(gas_data%name(n_spec))
435 allocate(gas_data%mosaic_index(n_spec))
440 call pmc_nc_check(nf90_inq_varid(ncid,
"gas_species", varid_gas_species))
441 allocate(
character(len=((GAS_NAME_LEN + 2) * 1000)) :: gas_species_names)
442 call pmc_nc_check(nf90_get_att(ncid, varid_gas_species,
"names", &
447 do while ((gas_species_names(i:i) /=
" ") &
448 .and. (gas_species_names(i:i) /=
","))
451 call assert(173021381, i > 1)
452 gas_data%name(i_spec) = gas_species_names(1:(i-1))
453 gas_species_names = gas_species_names((i+1):)
455 call assert(469721220, gas_species_names ==
"")