PartMC  2.3.0
extract_gas.F90
Go to the documentation of this file.
1 ! Copyright (C) 2009-2011 Matthew West
2 ! Licensed under the GNU General Public License version 2 or (at your
3 ! option) any later version. See the file COPYING for details.
4 
5 !> \file
6 !> The extract_gas program.
7 
8 !> Read NetCDF output files and write out the gas mixing ratios in text
9 !> format.
10 program extract_gas
11 
12  use pmc_gas_state
13  use pmc_output
14  use pmc_mpi
15  use getopt_m
16 
17  character(len=PMC_MAX_FILENAME_LEN) :: in_prefix, out_filename
18  character(len=PMC_MAX_FILENAME_LEN), allocatable :: filename_list(:)
19  character(len=1000) :: tmp_str
20  type(gas_data_t) :: gas_data
21  type(gas_state_t) :: gas_state
22  integer :: index, i_repeat, i_spec, out_unit
23  integer :: i_file, n_file
24  real(kind=dp) :: time, del_t
25  character(len=PMC_UUID_LEN) :: uuid, run_uuid
26  real(kind=dp), allocatable :: times(:), gas_mixing_ratios(:,:)
27  type(option_s) :: opts(2)
28 
29  call pmc_mpi_init()
30 
31  opts(1) = option_s("help", .false., 'h')
32  opts(2) = option_s("output", .true., 'o')
33 
34  out_filename = ""
35 
36  do
37  select case(getopt("ho:", opts))
38  case(char(0))
39  exit
40  case('h')
41  call print_help()
42  stop
43  case('o')
44  out_filename = optarg
45  case( '?' )
46  call print_help()
47  call die_msg(478715112, 'unknown option: ' // trim(optopt))
48  case default
49  call print_help()
50  call die_msg(935521190, 'unhandled option: ' // trim(optopt))
51  end select
52  end do
53 
54  if (optind /= command_argument_count()) then
55  call print_help()
56  call die_msg(744333329, 'expected exactly one non-option prefix argument')
57  end if
58 
59  call get_command_argument(optind, in_prefix)
60 
61  if (out_filename == "") then
62  out_filename = trim(in_prefix) // "_gas.txt"
63  end if
64 
65  call gas_data_allocate(gas_data)
66  call gas_state_allocate(gas_state)
67 
68  allocate(filename_list(0))
69  call input_filename_list(in_prefix, filename_list)
70  n_file = size(filename_list)
71  call assert_msg(579059629, n_file > 0, &
72  "no NetCDF files found with prefix: " // trim(in_prefix))
73 
74  call input_state(filename_list(1), index, time, del_t, i_repeat, uuid, &
75  gas_data=gas_data, gas_state=gas_state)
76  run_uuid = uuid
77 
78  allocate(times(n_file))
79  allocate(gas_mixing_ratios(n_file, gas_data%n_spec))
80 
81  do i_file = 1,n_file
82  call input_state(filename_list(i_file), index, time, del_t, i_repeat, &
83  uuid, gas_data=gas_data, gas_state=gas_state)
84 
85  call assert_msg(390171757, uuid == run_uuid, &
86  "UUID mismatch between " // trim(filename_list(1)) // " and " &
87  // trim(filename_list(i_file)))
88 
89  times(i_file) = time
90  gas_mixing_ratios(i_file, :) = gas_state%mix_rat
91  end do
92 
93  write(*,'(a,a)') "Output file: ", trim(out_filename)
94  write(*,'(a)') " Each row of output is one time."
95  write(*,'(a)') " The columns of output are:"
96  write(*,'(a)') " column 1: time (s)"
97  do i_spec = 1,gas_data%n_spec
98  write(*,'(a,i2,a,a,a)') " column ", i_spec + 1, ": gas ", &
99  trim(gas_data%name(i_spec)), " mixing ratio (ppb)"
100  end do
101 
102  call open_file_write(out_filename, out_unit)
103  do i_file = 1,n_file
104  write(out_unit, '(e30.15e3)', advance='no') times(i_file)
105  do i_spec = 1,gas_data%n_spec
106  write(out_unit, '(e30.15e3)', advance='no') &
107  gas_mixing_ratios(i_file, i_spec)
108  end do
109  write(out_unit, '(a)') ''
110  end do
111  call close_file(out_unit)
112 
113  deallocate(times)
114  deallocate(gas_mixing_ratios)
115  call gas_data_deallocate(gas_data)
116  call gas_state_deallocate(gas_state)
117 
118  call pmc_mpi_finalize()
119 
120 contains
121 
122  subroutine print_help()
123 
124  write(*,'(a)') 'Usage: extract_gas [options] <netcdf_prefix>'
125  write(*,'(a)') ''
126  write(*,'(a)') 'options are:'
127  write(*,'(a)') ' -h, --help Print this help message.'
128  write(*,'(a)') ' -o, --out <file> Output filename.'
129  write(*,'(a)') ''
130  write(*,'(a)') 'Examples:'
131  write(*,'(a)') ' extract_gas data_0001'
132  write(*,'(a)') ''
133 
134  end subroutine print_help
135 
136 end program extract_gas
subroutine input_state(filename, index, time, del_t, i_repeat, uuid, aero_data, aero_state, gas_data, gas_state, env_state)
Read the current state.
Definition: output.F90:509
subroutine die_msg(code, error_msg)
Error immediately.
Definition: util.F90:133
character function getopt(optstring, longopts)
Definition: getopt.F90:131
subroutine input_filename_list(prefix, filename_list)
Find all NetCDF (.nc) filenames that match the given prefix.
Definition: output.F90:580
subroutine close_file(unit)
Close a file and de-assign the unit.
Definition: util.F90:225
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
Definition: util.F90:76
subroutine pmc_mpi_finalize()
Shut down MPI.
Definition: mpi.F90:88
subroutine pmc_mpi_init()
Initialize MPI.
Definition: mpi.F90:55
subroutine gas_data_deallocate(gas_data)
Free all storage.
Definition: gas_data.F90:75
program extract_gas
Read NetCDF output files and write out the gas mixing ratios in text format.
Definition: extract_gas.F90:10
Wrapper functions for MPI.
Definition: mpi.F90:13
subroutine gas_state_deallocate(gas_state)
Free all storage.
Definition: gas_state.F90:63
The gas_state_t structure and associated subroutines.
Definition: gas_state.F90:9
subroutine gas_state_allocate(gas_state)
Allocate storage for gas species.
Definition: gas_state.F90:36
subroutine print_help()
subroutine gas_data_allocate(gas_data)
Allocate storage for gas species.
Definition: gas_data.F90:45
Current state of the gas mixing ratios in the system.
Definition: gas_state.F90:26
Write data in NetCDF format.
Definition: output.F90:68
subroutine open_file_write(filename, unit)
Open a file for writing with an automatically assigned unit and test that it succeeds. The file should be closed with close_file().
Definition: util.F90:205