Skip to content

Commit 255a97a

Browse files
magpowellmagpowell
committed
aerosols
1 parent 72843a0 commit 255a97a

File tree

4 files changed

+273
-93
lines changed

4 files changed

+273
-93
lines changed

‎include_tilt/tilt_utils.h‎

Lines changed: 30 additions & 24 deletions
Original file line numberDiff line numberDiff line change
@@ -62,17 +62,20 @@ void restore_bkg_profile(const int n_x, const int n_y,
6262
std::vector<Float>& var,
6363
std::vector<Float>& var_w_bkg);
6464

65-
void restore_bkg_profile_bundle(const int n_col_x, const int n_col_y, const int n_lay, const int n_lev,
65+
void restore_bkg_profile_bundle(const int n_col_x, const int n_col_y,
66+
const int n_lay, const int n_lev,
6667
const int n_lay_tot, const int n_lev_tot,
67-
const int n_z_in, const int n_zh_in, const int bkg_start_z, const int bkg_start_zh,
68+
const int n_z_in, const int n_zh_in,
69+
const int bkg_start_z, const int bkg_start_zh,
6870
Array<Float,2>* p_lay_copy, Array<Float,2>* t_lay_copy, Array<Float,2>* p_lev_copy, Array<Float,2>* t_lev_copy,
69-
Array<Float,2>* lwp_copy, Array<Float,2>* iwp_copy, Array<Float,2>* rel_copy, Array<Float,2>* dei_copy,
70-
Gas_concs& gas_concs_copy,
71+
Array<Float,2>* lwp_copy, Array<Float,2>* iwp_copy, Array<Float,2>* rel_copy, Array<Float,2>* dei_copy, Array<Float,2>* rh_copy,
72+
Gas_concs& gas_concs_copy, Aerosol_concs& aerosol_concs_copy,
7173
Array<Float,2>* p_lay, Array<Float,2>* t_lay, Array<Float,2>* p_lev, Array<Float,2>* t_lev,
72-
Array<Float,2>* lwp, Array<Float,2>* iwp, Array<Float,2>* rel, Array<Float,2>* dei,
73-
Gas_concs& gas_concs,
74-
std::vector<std::string> gas_names, bool switch_liq_cloud_optics, bool switch_ice_cloud_optics
75-
);
74+
Array<Float,2>* lwp, Array<Float,2>* iwp, Array<Float,2>* rel, Array<Float,2>* dei, Array<Float,2>* rh,
75+
Gas_concs& gas_concs, Aerosol_concs& aerosol_concs,
76+
std::vector<std::string> gas_names, std::vector<std::string> aerosol_names,
77+
bool switch_liq_cloud_optics, bool switch_ice_cloud_optics, bool switch_aerosol_optics
78+
);
7679

7780
void compress_columns_weighted_avg(const int n_x, const int n_y,
7881
const int n_out,
@@ -90,12 +93,17 @@ void tilt_fields(const int n_z_in, const int n_zh_in, const int n_col_x, const i
9093
const Array<Float,1> zh, const Array<Float,1> z,
9194
const Array<Float,1> zh_tilt, const Array<ijk,1> path,
9295
Array<Float,2>* p_lay_copy, Array<Float,2>* t_lay_copy, Array<Float,2>* p_lev_copy, Array<Float,2>* t_lev_copy,
93-
Gas_concs& gas_concs_copy, const std::vector<std::string> gas_names);
96+
Array<Float,2>* rh_copy,
97+
Gas_concs& gas_concs_copy, const std::vector<std::string> gas_names,
98+
Aerosol_concs& aerosol_concs_copy, const std::vector<std::string> aerosol_names, const bool switch_aerosol_optics
99+
);
94100

95101
void compress_fields(const int compress_lay_start_idx, const int n_col_x, const int n_col_y,
96102
const int n_z_in, const int n_zh_in, const int n_z_tilt,
97103
Array<Float,2>* p_lay_copy, Array<Float,2>* t_lay_copy, Array<Float,2>* p_lev_copy, Array<Float,2>* t_lev_copy,
98-
Gas_concs& gas_concs_copy, std::vector<std::string> gas_names);
104+
Array<Float,2>* rh_copy,
105+
Gas_concs& gas_concs_copy, std::vector<std::string> gas_names,
106+
Aerosol_concs& aerosol_concs_copy, std::vector<std::string> aerosol_names, const bool switch_aerosol_optics);
99107

100108
void create_tilted_columns(const int n_x, const int n_y, const int n_lay_in, const int n_lev_in,
101109
const std::vector<Float>& zh_tilted, const std::vector<ijk>& tilted_path,
@@ -112,17 +120,15 @@ void create_tilted_columns_levlay(const int n_x, const int n_y, const int n_lay_
112120
const std::vector<Float>& zh_tilted, const std::vector<ijk>& tilted_path,
113121
std::vector<Float>& var_lay, std::vector<Float>& var_lev);
114122

115-
void tica_tilt(
116-
const Float sza, const Float azi,
117-
const int n_col_x, const int n_col_y, const int n_col,
118-
const int n_lay, const int n_lev, const int n_z_in, const int n_zh_in ,
119-
Array<Float,1> xh, Array<Float,1> yh, Array<Float,1> zh, Array<Float,1> z,
120-
Array<Float,2> p_lay, Array<Float,2> t_lay, Array<Float,2> p_lev, Array<Float,2> t_lev,
121-
Array<Float,2> lwp, Array<Float,2> iwp, Array<Float,2> rel, Array<Float,2> dei,
122-
Gas_concs gas_concs,
123-
Array<Float,2>& p_lay_out, Array<Float,2>& t_lay_out, Array<Float,2>& p_lev_out, Array<Float,2>& t_lev_out,
124-
Array<Float,2>& lwp_out, Array<Float,2>& iwp_out, Array<Float,2>& rel_out, Array<Float,2>& dei_out,
125-
Gas_concs& gas_concs_out,
126-
std::vector<std::string> gas_names,
127-
bool switch_cloud_optics, bool switch_liq_cloud_optics, bool switch_ice_cloud_optics
128-
);
123+
void tica_tilt(const Float sza, const Float azi,
124+
const int n_col_x, const int n_col_y, const int n_col,
125+
const int n_lay, const int n_lev, const int n_z_in, const int n_zh_in ,
126+
Array<Float,1> xh, Array<Float,1> yh, Array<Float,1> zh, Array<Float,1> z,
127+
Array<Float,2> p_lay, Array<Float,2> t_lay, Array<Float,2> p_lev, Array<Float,2> t_lev,
128+
Array<Float,2> lwp, Array<Float,2> iwp, Array<Float,2> rel, Array<Float,2> dei, Array<Float,2> rh,
129+
Gas_concs gas_concs, Aerosol_concs aerosol_concs,
130+
Array<Float,2>& p_lay_out, Array<Float,2>& t_lay_out, Array<Float,2>& p_lev_out, Array<Float,2>& t_lev_out,
131+
Array<Float,2>& lwp_out, Array<Float,2>& iwp_out, Array<Float,2>& rel_out, Array<Float,2>& dei_out, Array<Float,2>& rh_out,
132+
Gas_concs& gas_concs_out, Aerosol_concs aerosol_concs_out,
133+
std::vector<std::string> gas_names, std::vector<std::string> aerosol_names,
134+
bool switch_cloud_optics, bool switch_liq_cloud_optics, bool switch_ice_cloud_optics, bool switch_aerosol_optics);

‎src_test/test_rt_tilt_prep_mc.cpp‎

Lines changed: 100 additions & 35 deletions
Original file line numberDiff line numberDiff line change
@@ -1,33 +1,68 @@
11
/*
2-
* This file is a stand-alone executable developed for the
3-
* testing of the C++ interface to the RTE+RRTMGP radiation code.
4-
*
5-
* It is free software: you can redistribute it and/or modify
6-
* it under the terms of the GNU General Public License as published by
7-
* the Free Software Foundation, either version 3 of the License, or
8-
* (at your option) any later version.
9-
*
10-
* This software is distributed in the hope that it will be useful,
11-
* but WITHOUT ANY WARRANTY; without even the implied warranty of
12-
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13-
* GNU General Public License for more details.
14-
*
15-
* You should have received a copy of the GNU General Public License
16-
* along with this software. If not, see <http://www.gnu.org/licenses/>.
17-
*/
18-
19-
#include <boost/algorithm/string.hpp>
20-
#include <chrono>
21-
#include <cmath>
22-
#include <iomanip>
23-
#include <cuda_profiler_api.h>
24-
25-
#include "Status.h"
26-
#include "Netcdf_interface.h"
27-
#include "Array.h"
28-
#include "Gas_concs.h"
29-
#include "tilt_utils.h"
30-
#include "types.h"
2+
* This file is a stand-alone executable developed for the
3+
* testing of the C++ interface to the RTE+RRTMGP radiation code.
4+
*
5+
* It is free software: you can redistribute it and/or modify
6+
* it under the terms of the GNU General Public License as published by
7+
* the Free Software Foundation, either version 3 of the License, or
8+
* (at your option) any later version.
9+
*
10+
* This software is distributed in the hope that it will be useful,
11+
* but WITHOUT ANY WARRANTY; without even the implied warranty of
12+
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13+
* GNU General Public License for more details.
14+
*
15+
* You should have received a copy of the GNU General Public License
16+
* along with this software. If not, see <http://www.gnu.org/licenses/>.
17+
*/
18+
19+
#include <boost/algorithm/string.hpp>
20+
#include <chrono>
21+
#include <cmath>
22+
#include <iomanip>
23+
#include <cuda_profiler_api.h>
24+
25+
#include "Status.h"
26+
#include "Netcdf_interface.h"
27+
#include "Array.h"
28+
#include "Gas_concs.h"
29+
#include "Aerosol_optics_rt.h"
30+
#include "tilt_utils.h"
31+
#include "types.h"
32+
33+
void read_and_set_aer(
34+
const std::string& aerosol_name, const int n_col_x, const int n_col_y, const int n_lay,
35+
const Netcdf_handle& input_nc, Aerosol_concs& aerosol_concs)
36+
{
37+
if (input_nc.variable_exists(aerosol_name))
38+
{
39+
std::map<std::string, int> dims = input_nc.get_variable_dimensions(aerosol_name);
40+
const int n_dims = dims.size();
41+
42+
if (n_dims == 1)
43+
{
44+
if (dims.at("lay") == n_lay)
45+
aerosol_concs.set_vmr(aerosol_name,
46+
Array<Float,1>(input_nc.get_variable<Float>(aerosol_name, {n_lay}), {n_lay}));
47+
else
48+
throw std::runtime_error("Illegal dimensions of \"" + aerosol_name + "\" in input");
49+
}
50+
else if (n_dims == 3)
51+
{
52+
if (dims.at("lay") == n_lay && dims.at("y") == n_col_y && dims.at("x") == n_col_x)
53+
aerosol_concs.set_vmr(aerosol_name,
54+
Array<Float,2>(input_nc.get_variable<Float>(aerosol_name, {n_lay, n_col_y, n_col_x}), {n_col_x * n_col_y, n_lay}));
55+
else
56+
throw std::runtime_error("Illegal dimensions of \"" + aerosol_name + "\" in input");
57+
}
58+
else
59+
throw std::runtime_error("Illegal dimensions of \"" + aerosol_name + "\" in input");
60+
}
61+
else
62+
{
63+
throw std::runtime_error("Aerosol type \"" + aerosol_name + "\" not available in input file.");
64+
}
65+
}
3166

3267
void read_and_set_vmr(
3368
const std::string& gas_name, const int n_col_x, const int n_col_y, const int n_lay,
@@ -164,6 +199,11 @@ void print_command_line_options(
164199
"cf4", "no2"
165200
};
166201

202+
std::vector<std::string> aerosol_names = {
203+
"aermr01", "aermr02", "aermr03", "aermr04", "aermr05", "aermr06", "aermr07",
204+
"aermr08", "aermr09", "aermr10","aermr11"
205+
};
206+
167207
Status::print_message("###### Starting Script ######");
168208

169209
////// FLOW CONTROL SWITCHES //////
@@ -172,6 +212,7 @@ void print_command_line_options(
172212
{"cloud-optics" , { false, "Enable cloud optics (both liquid and ice)."}},
173213
{"liq-cloud-optics" , { false, "liquid only cloud optics." }},
174214
{"ice-cloud-optics" , { false, "ice only cloud optics." }},
215+
{"aerosol-optics" , { false, "aerosol optics." }},
175216
{"tilt-sza" , { false, "tilt provided value of sza in input file. IN DEGREES. '--tilt-sza 50': use a sza of 50 degrees" }},
176217
{"tilt-azi" , { false, "tilt provided value of azi in input file. FROM POS Y, CLOCKWISE, IN DEGREES. '--tilt-azi 240': use of azi of 240 degrees" }}};
177218

@@ -185,6 +226,7 @@ void print_command_line_options(
185226
bool switch_cloud_optics = command_line_switches.at("cloud-optics" ).first;
186227
bool switch_liq_cloud_optics = command_line_switches.at("liq-cloud-optics" ).first;
187228
bool switch_ice_cloud_optics = command_line_switches.at("ice-cloud-optics" ).first;
229+
const bool switch_aerosol_optics = command_line_switches.at("aerosol-optics").first;
188230
const bool switch_tilt_sza = command_line_switches.at("tilt-sza" ).first;
189231
const bool switch_tilt_azi = command_line_switches.at("tilt-azi" ).first;
190232

@@ -276,6 +318,27 @@ void print_command_line_options(
276318
read_and_set_vmr("cf4" , n_col_x, n_col_y, n_lay, input_nc, gas_concs);
277319
read_and_set_vmr("no2" , n_col_x, n_col_y, n_lay, input_nc, gas_concs);
278320

321+
Array<Float,2> rh;
322+
Aerosol_concs aerosol_concs;
323+
324+
if (switch_aerosol_optics)
325+
{
326+
rh.set_dims({n_col, n_lay});
327+
rh = std::move(input_nc.get_variable<Float>("rh", {n_lay, n_col_y, n_col_x}));
328+
329+
read_and_set_aer("aermr01", n_col_x, n_col_y, n_lay, input_nc, aerosol_concs);
330+
read_and_set_aer("aermr02", n_col_x, n_col_y, n_lay, input_nc, aerosol_concs);
331+
read_and_set_aer("aermr03", n_col_x, n_col_y, n_lay, input_nc, aerosol_concs);
332+
read_and_set_aer("aermr04", n_col_x, n_col_y, n_lay, input_nc, aerosol_concs);
333+
read_and_set_aer("aermr05", n_col_x, n_col_y, n_lay, input_nc, aerosol_concs);
334+
read_and_set_aer("aermr06", n_col_x, n_col_y, n_lay, input_nc, aerosol_concs);
335+
read_and_set_aer("aermr07", n_col_x, n_col_y, n_lay, input_nc, aerosol_concs);
336+
read_and_set_aer("aermr08", n_col_x, n_col_y, n_lay, input_nc, aerosol_concs);
337+
read_and_set_aer("aermr09", n_col_x, n_col_y, n_lay, input_nc, aerosol_concs);
338+
read_and_set_aer("aermr10", n_col_x, n_col_y, n_lay, input_nc, aerosol_concs);
339+
read_and_set_aer("aermr11", n_col_x, n_col_y, n_lay, input_nc, aerosol_concs);
340+
}
341+
279342
for (const auto& gas_name : gas_names) {
280343
if (!gas_concs.exists(gas_name)) {
281344
continue;
@@ -320,6 +383,8 @@ void print_command_line_options(
320383
Array<Float,2> p_lay_out = p_lay;
321384
Array<Float,2> p_lev_out = p_lev;
322385
Gas_concs gas_concs_out = gas_concs;
386+
Aerosol_concs aerosol_concs_out = aerosol_concs;
387+
Array<Float,2> rh_out = rh;
323388

324389
Array<Float,2> lwp_out;
325390
lwp_out.set_dims({n_col, n_z_in});
@@ -336,13 +401,13 @@ void print_command_line_options(
336401
n_lay, n_lev, n_z_in, n_zh_in ,
337402
xh, yh, zh, z,
338403
p_lay, t_lay, p_lev, t_lev,
339-
lwp, iwp, rel, dei,
340-
gas_concs,
404+
lwp, iwp, rel, dei, rh,
405+
gas_concs, aerosol_concs,
341406
p_lay_out, t_lay_out, p_lev_out, t_lev_out,
342-
lwp_out, iwp_out, rel_out, dei_out,
343-
gas_concs_out,
344-
gas_names,
345-
switch_cloud_optics, switch_liq_cloud_optics, switch_ice_cloud_optics
407+
lwp_out, iwp_out, rel_out, dei_out, rh_out,
408+
gas_concs_out, aerosol_concs_out,
409+
gas_names, aerosol_names,
410+
switch_cloud_optics, switch_liq_cloud_optics, switch_ice_cloud_optics, switch_aerosol_optics
346411
);
347412

348413
std::vector<Float> z_out_compress = linspace(z.v()[0], z.v()[n_z_in - 1], n_z_in);

‎src_test/test_rte_rrtmgp_rt.cu‎

Lines changed: 28 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -352,6 +352,22 @@ void solve_radiation(int argc, char** argv)
352352
Array<Float,2> p_lev(input_nc.get_variable<Float>("p_lev", {n_lev, n_col_y, n_col_x}), {n_col, n_lev});
353353
Array<Float,2> t_lev(input_nc.get_variable<Float>("t_lev", {n_lev, n_col_y, n_col_x}), {n_col, n_lev});
354354

355+
356+
357+
if (input_nc.variable_exists("col_dry") & switch_tica)
358+
{
359+
std::string error = "col_dry is not supported in tica mode";
360+
throw std::runtime_error(error);
361+
362+
}
363+
364+
if (input_nc.variable_exists("tsi") & switch_tica)
365+
{
366+
std::string error = "tsi is overwritten in tica mode";
367+
throw std::runtime_error(error);
368+
369+
}
370+
355371
// Fetch the col_dry in case present.
356372
Array<Float,2> col_dry;
357373
if (input_nc.variable_exists("col_dry"))
@@ -459,21 +475,21 @@ void solve_radiation(int argc, char** argv)
459475

460476
if (switch_tica)
461477
{
462-
std::cout << "tica_azi: " << tica_azi << std::endl;
463-
std::cout << "tica_sza: " << tica_sza << std::endl;
464478
for (int icol=1; icol<=n_col; ++icol)
465479
{
466480
mu0({icol}) = 1.0;
467481
azi({icol}) = 0.0;
468482
}
469-
// what to do about col dry and aerosols?
470-
// what to do about tsi
471483

472484
std::vector<std::string> gas_names = {
473485
"h2o", "co2", "o3", "n2o", "co", "ch4", "o2", "n2", "ccl4", "cfc11",
474486
"cfc12", "cfc22", "hfc143a", "hfc125", "hfc23", "hfc32", "hfc134a",
475487
"cf4", "no2"
476488
};
489+
std::vector<std::string> aerosol_names = {
490+
"aermr01", "aermr02", "aermr03", "aermr04", "aermr05", "aermr06", "aermr07",
491+
"aermr08", "aermr09", "aermr10","aermr11"
492+
};
477493

478494
Array<Float,1> xh;
479495
Array<Float,1> yh;
@@ -495,6 +511,8 @@ void solve_radiation(int argc, char** argv)
495511
Array<Float,2> p_lay_out = p_lay;
496512
Array<Float,2> p_lev_out = p_lev;
497513
Gas_concs gas_concs_out = gas_concs;
514+
Aerosol_concs aerosol_concs_out = aerosol_concs;
515+
Array<Float,2> rh_out = rh;
498516

499517
Array<Float,2> lwp_out;
500518
lwp_out.set_dims({n_col, n_z_in});
@@ -511,13 +529,13 @@ void solve_radiation(int argc, char** argv)
511529
n_lay, n_lev, n_z_in, n_zh_in ,
512530
xh, yh, zh, z,
513531
p_lay, t_lay, p_lev, t_lev,
514-
lwp, iwp, rel, dei,
515-
gas_concs,
532+
lwp, iwp, rel, dei, rh,
533+
gas_concs, aerosol_concs,
516534
p_lay_out, t_lay_out, p_lev_out, t_lev_out,
517-
lwp_out, iwp_out, rel_out, dei_out,
518-
gas_concs_out,
519-
gas_names,
520-
switch_cloud_optics, switch_liq_cloud_optics, switch_ice_cloud_optics
535+
lwp_out, iwp_out, rel_out, dei_out, rh_out,
536+
gas_concs_out, aerosol_concs_out,
537+
gas_names, aerosol_names,
538+
switch_cloud_optics, switch_liq_cloud_optics, switch_ice_cloud_optics, switch_aerosol_optics
521539
);
522540

523541
lwp_out.expand_dims({n_col, n_lay});

0 commit comments

Comments
 (0)