Skip to content

Commit

Permalink
enable running broadband and image(XYZ) ray tracer variants both, ins…
Browse files Browse the repository at this point in the history 
…tead of one of the two
  • Loading branch information
@MennoVeerman
MennoVeerman committed Aug 19, 2024
1 parent 82df7cb commit 8966d91
Show file tree
Hide file tree
Showing 4 changed files with 61 additions and 50 deletions.
Binary file modified data/mie_lut_broadband.nc
View file
Binary file not shown.
13 changes: 10 additions & 3 deletions include_test/Radiation_solver_bw.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -145,8 +145,10 @@ class Radiation_solver_shortwave
Array<Float,3>& sw_bnd_flux_dn_dir, Array<Float,3>& sw_bnd_flux_net) const;

void load_mie_tables(
const std::string& file_name_mie,
const bool switch_broadband);
const std::string& file_name_mie_bb,
const std::string& file_name_mie_im,
const bool switch_broadband,
const bool switch_image);

#ifdef __CUDACC__
void solve_gpu(
Expand Down Expand Up @@ -244,11 +246,16 @@ class Radiation_solver_shortwave
std::unique_ptr<Optical_props_2str_rt> cloud_optical_props;
std::unique_ptr<Optical_props_2str_rt> aerosol_optical_props;

Array_gpu<Float,1> mie_phase_angs;
Array_gpu<Float,1> mie_phase_angs_bb;
Array_gpu<Float,3> mie_cdfs_bb;
Array_gpu<Float,4> mie_angs_bb;
Array_gpu<Float,4> mie_phase_bb;

Array_gpu<Float,1> mie_phase_angs;
Array_gpu<Float,3> mie_cdfs;
Array_gpu<Float,4> mie_angs;
Array_gpu<Float,4> mie_phase;

#endif
};
#endif
54 changes: 31 additions & 23 deletions src_test/Radiation_solver_bw.cu
View file
Original file line number Diff line number Diff line change
Expand Up @@ -733,31 +733,39 @@ Radiation_solver_shortwave::Radiation_solver_shortwave(
}

void Radiation_solver_shortwave::load_mie_tables(
const std::string& file_name_mie,
const bool switch_broadband)
const std::string& file_name_mie_bb,
const std::string& file_name_mie_im,
const bool switch_broadband,
const bool switch_image)
{
Netcdf_file mie_nc(file_name_mie, Netcdf_mode::Read);
const int n_bnd_sw = this->get_n_bnd_gpu();
const int n_re = mie_nc.get_dimension_size("r_eff");
const int n_mie = mie_nc.get_dimension_size("n_ang");

if (switch_broadband)
{
Array<Float,3> mie_cdf(mie_nc.get_variable<Float>("phase_cdf", {n_bnd_sw, 1, n_mie}), {n_mie, 1, n_bnd_sw});
Array<Float,4> mie_ang(mie_nc.get_variable<Float>("phase_cdf_angle", {n_bnd_sw, 1, n_re, n_mie}), {n_mie, n_re, 1, n_bnd_sw});

Array<Float,4> mie_phase(mie_nc.get_variable<Float>("phase", {n_bnd_sw, 1, n_re, n_mie}), {n_mie, n_re, 1, n_bnd_sw});
Netcdf_file mie_nc(file_name_mie_bb, Netcdf_mode::Read);
const int n_bnd_sw = this->get_n_bnd_gpu();
const int n_re = mie_nc.get_dimension_size("r_eff");
const int n_mie = mie_nc.get_dimension_size("n_ang");

Array<Float,3> mie_cdf(mie_nc.get_variable<Float>("phase_cdf", {n_bnd_sw, n_mie}), {n_mie, 1, n_bnd_sw});
Array<Float,4> mie_ang(mie_nc.get_variable<Float>("phase_cdf_angle", {n_bnd_sw, n_re, n_mie}), {n_mie, n_re, 1, n_bnd_sw});

Array<Float,4> mie_phase(mie_nc.get_variable<Float>("phase", {n_bnd_sw, n_re, n_mie}), {n_mie, n_re, 1, n_bnd_sw});
Array<Float,1> mie_phase_ang(mie_nc.get_variable<Float>("phase_angle", {n_mie}), {n_mie});

this->mie_cdfs = mie_cdf;
this->mie_angs = mie_ang;
this->mie_cdfs_bb = mie_cdf;
this->mie_angs_bb = mie_ang;

this->mie_phase = mie_phase;
this->mie_phase_angs = mie_phase_ang;
this->mie_phase_bb = mie_phase;
this->mie_phase_angs_bb = mie_phase_ang;
}
else

if (switch_image)
{
Netcdf_file mie_nc(file_name_mie_im, Netcdf_mode::Read);
const int n_bnd_sw = this->get_n_bnd_gpu();
const int n_re = mie_nc.get_dimension_size("r_eff");
const int n_mie = mie_nc.get_dimension_size("n_ang");
const int n_sub = mie_nc.get_dimension_size("sub_band");

Array<Float,3> mie_cdf(mie_nc.get_variable<Float>("phase_cdf", {n_bnd_sw, n_sub, n_mie}), {n_mie, n_sub, n_bnd_sw});
Array<Float,4> mie_ang(mie_nc.get_variable<Float>("phase_cdf_angle", {n_bnd_sw, n_sub, n_re, n_mie}), {n_mie, n_re, n_sub, n_bnd_sw});

Expand Down Expand Up @@ -1126,9 +1134,9 @@ void Radiation_solver_shortwave::solve_gpu_bb(
const int n_gpt = this->kdist_gpu->get_ngpt();
const int n_bnd = this->kdist_gpu->get_nband();

const int n_mie = (switch_cloud_mie) ? this->mie_angs.dim(1) : 0;
const int n_re = (switch_cloud_mie) ? this->mie_angs.dim(2) : 0;

const int n_mie = (switch_cloud_mie) ? this->mie_angs_bb.dim(1) : 0;
const int n_re = (switch_cloud_mie) ? this->mie_angs_bb.dim(2) : 0;
const int cam_nx = radiance.dim(1);
const int cam_ny = radiance.dim(2);
const Bool top_at_1 = p_lay({1, 1}) < p_lay({1, n_lay});
Expand Down Expand Up @@ -1285,9 +1293,9 @@ void Radiation_solver_shortwave::solve_gpu_bb(

if (switch_cloud_mie)
{
mie_cdfs_sub = mie_cdfs.subset({{ {1, n_mie}, {1,1}, {band, band} }});
mie_angs_sub = mie_angs.subset({{ {1, n_mie}, {1, n_re}, {1,1}, {band, band} }});
mie_phase_sub = mie_phase.subset({{ {1, n_mie}, {1, n_re}, {1,1}, {band, band} }});
mie_cdfs_sub = mie_cdfs_bb.subset({{ {1, n_mie}, {1,1}, {band, band} }});
mie_angs_sub = mie_angs_bb.subset({{ {1, n_mie}, {1, n_re}, {1,1}, {band, band} }});
mie_phase_sub = mie_phase_bb.subset({{ {1, n_mie}, {1, n_re}, {1,1}, {band, band} }});
}

raytracer.trace_rays_bb(
Expand All @@ -1298,7 +1306,7 @@ void Radiation_solver_shortwave::solve_gpu_bb(
mie_cdfs_sub,
mie_angs_sub,
mie_phase_sub,
mie_phase_angs,
mie_phase_angs_bb,
rel,
dynamic_cast<Optical_props_2str_rt&>(*optical_props).get_tau(),
dynamic_cast<Optical_props_2str_rt&>(*optical_props).get_ssa(),
Expand Down
44 changes: 20 additions & 24 deletions src_test/test_rte_rrtmgp_bw.cu
View file
Original file line number Diff line number Diff line change
Expand Up @@ -245,7 +245,8 @@ void solve_radiation(int argc, char** argv)
{"output-optical" , { false, "Enable output of optical properties." }},
{"output-bnd-fluxes", { false, "Enable output of band fluxes." }},
{"lu-albedo" , { false, "Compute spectral albedo from land use map" }},
{"broadband" , { false, "Compute broadband fluxes" }},
{"image" , { true, "Compute XYZ values to generate RGB images" }},
{"broadband" , { false, "Compute broadband radiances" }},
{"profiling" , { false, "Perform additional profiling run." }},
{"delta-cloud" , { false, "delta-scaling of cloud optical properties" }},
{"delta-aerosol" , { false, "delta-scaling of aerosol optical properties" }},
Expand All @@ -265,12 +266,13 @@ void solve_radiation(int argc, char** argv)
const bool switch_output_optical = command_line_options.at("output-optical" ).first;
const bool switch_output_bnd_fluxes = command_line_options.at("output-bnd-fluxes").first;
const bool switch_lu_albedo = command_line_options.at("lu-albedo" ).first;
const bool switch_image = command_line_options.at("image" ).first;
const bool switch_broadband = command_line_options.at("broadband" ).first;
const bool switch_profiling = command_line_options.at("profiling" ).first;
const bool switch_delta_cloud = command_line_options.at("delta-cloud" ).first;
const bool switch_delta_aerosol = command_line_options.at("delta-aerosol" ).first;
const bool switch_cloud_cam = command_line_options.at("cloud-cam" ).first;
const bool switch_raytracing = command_line_options.at("raytracing" ).first;
const bool switch_cloud_cam = command_line_options.at("cloud-cam" ).first;
const bool switch_raytracing = command_line_options.at("raytracing" ).first;

if (switch_longwave)
{
Expand Down Expand Up @@ -682,18 +684,16 @@ void solve_radiation(int argc, char** argv)
if (switch_broadband)
{
radiance.set_dims({camera.nx, camera.ny});

if (switch_cloud_mie)
rad_sw.load_mie_tables("mie_lut_broadband.nc", switch_broadband);
}
else
if (switch_image)
{
XYZ.set_dims({camera.nx, camera.ny, 3});

if (switch_cloud_mie)
rad_sw.load_mie_tables("mie_lut_visualisation.nc", switch_broadband);
}

if (switch_cloud_mie)
rad_sw.load_mie_tables("mie_lut_broadband.nc", "mie_lut_visualisation.nc", switch_broadband, switch_image);


Array_gpu<Float,2> liwp_cam;
Array_gpu<Float,2> tauc_cam;
Array_gpu<Float,2> dist_cam;
Expand Down Expand Up @@ -781,7 +781,7 @@ void solve_radiation(int argc, char** argv)
cudaEventDestroy(start);
cudaEventDestroy(stop);

Status::print_message("Duration shortwave solver: " + std::to_string(duration) + " (ms)");
Status::print_message("Duration shortwave solver (broadband version): " + std::to_string(duration) + " (ms)");
};

auto run_solver = [&](const bool tune_step)
Expand Down Expand Up @@ -856,7 +856,7 @@ void solve_radiation(int argc, char** argv)
cudaEventDestroy(start);
cudaEventDestroy(stop);

Status::print_message("Duration shortwave solver: " + std::to_string(duration) + " (ms)");
Status::print_message("Duration shortwave solver (image version): " + std::to_string(duration) + " (ms)");
};

if (switch_broadband)
Expand All @@ -871,7 +871,7 @@ void solve_radiation(int argc, char** argv)
cudaProfilerStop();
}
}
else
if (switch_image)
{
// tune step
run_solver(true);
Expand All @@ -893,30 +893,26 @@ void solve_radiation(int argc, char** argv)

if (switch_raytracing)
{
output_nc.add_dimension("gpt_sw", n_gpt_sw);
output_nc.add_dimension("band_sw", n_bnd_sw);

auto nc_sw_band_lims_wvn = output_nc.add_variable<Float>("sw_band_lims_wvn", {"band_sw", "pair"});
nc_sw_band_lims_wvn.insert(rad_sw.get_band_lims_wavenumber_gpu().v(), {0, 0});

if (switch_broadband)
{
Array<Float,2> radiance_cpu(radiance);
output_nc.add_dimension("gpt_sw", n_gpt_sw);
output_nc.add_dimension("band_sw", n_bnd_sw);

auto nc_sw_band_lims_wvn = output_nc.add_variable<Float>("sw_band_lims_wvn", {"band_sw", "pair"});
nc_sw_band_lims_wvn.insert(rad_sw.get_band_lims_wavenumber_gpu().v(), {0, 0});

auto nc_var = output_nc.add_variable<Float>("radiance", {"y","x"});
nc_var.insert(radiance_cpu.v(), {0, 0});
nc_var.add_attribute("long_name", "shortwave radiance");
nc_var.add_attribute("units", "W m-2 sr-1");
}
else
if (switch_image)
{
Array<Float,3> xyz_cpu(XYZ);
output_nc.add_dimension("gpt_sw", n_gpt_sw);
output_nc.add_dimension("band_sw", n_bnd_sw);
output_nc.add_dimension("n",3);

auto nc_sw_band_lims_wvn = output_nc.add_variable<Float>("sw_band_lims_wvn", {"band_sw", "pair"});
nc_sw_band_lims_wvn.insert(rad_sw.get_band_lims_wavenumber_gpu().v(), {0, 0});

auto nc_xyz = output_nc.add_variable<Float>("XYZ", {"n","y","x"});
nc_xyz.insert(xyz_cpu.v(), {0, 0, 0});

Expand Down

0 comments on commit 8966d91

Please sign in to comment.