diff --git a/drawbicyc/include/cases/Cumulus_Congestus/plots.hpp b/drawbicyc/include/cases/Cumulus_Congestus/plots.hpp index 637a626..a2d4095 100644 --- a/drawbicyc/include/cases/Cumulus_Congestus/plots.hpp +++ b/drawbicyc/include/cases/Cumulus_Congestus/plots.hpp @@ -17,6 +17,14 @@ const std::vector series_Cumulus_Congestus({ //"RH_max", "cloud_top_height", "total_droplets_number", + "cl_meanr", + "cl_nc", + "cl_nr", + "nc", + "nr", + "com_mom0", + "com_mom1", + "com_mom2", "cloud_cover_rico", //"surf_flux_latent", //"surf_flux_sensible", @@ -34,20 +42,27 @@ const std::vector series_Cumulus_Congestus({ }); std::vector profs_Cumulus_Congestus({ - "00rtot", + "actrw_rw_cl_conc", "rliq", + "actrw_rw_cl", + "actrw_mom1", + "actrw_mom2", + "disp_r", + //"00rtot", + //"actrw_reff_cl", + //"rliq", //"thl", //"wvar", - "prflux", + //"prflux", //"clfrac", //"sd_conc", - "cl_nc", + //"cl_nc", //"cl_nc_up", //"w", //"u", //"v", - "coal_tele", - "base_prflux_vs_clhght", + //"coal_tele", + //"base_prflux_vs_clhght", //"non_gccn_rw_cl", //"gccn_rw_cl," //, "N_c", diff --git a/drawbicyc/include/gnuplot_profs_set_labels.hpp b/drawbicyc/include/gnuplot_profs_set_labels.hpp index f2cfd8e..b13abd3 100644 --- a/drawbicyc/include/gnuplot_profs_set_labels.hpp +++ b/drawbicyc/include/gnuplot_profs_set_labels.hpp @@ -13,7 +13,8 @@ void gnuplot_profs_set_labels(Gnuplot &gp, std::string plt, const bool normalize if (plt == "rliq") { - gp << "set title 'liquid water [g/kg]'\n"; + //gp << "set title 'liquid water [g/kg]'\n"; + gp << "set title 'liquid water in cloudy cells [g/kg]'\n"; } if (plt == "gccn_rw") { @@ -79,6 +80,10 @@ void gnuplot_profs_set_labels(Gnuplot &gp, std::string plt, const bool normalize { gp << "set title 'activated (rw>rc) droplets mean dry radius'\n"; } + if (plt == "actrw_reff_cl") + { + gp << "set title 'effective radius of activated droplets (r > rc) in cloudy cells'\n"; + } else if (plt == "rv") { gp << "set title 'rv [g/kg]'\n"; @@ -214,4 +219,31 @@ void gnuplot_profs_set_labels(Gnuplot &gp, std::string plt, const bool normalize { gp << "set title 'mean cloud droplet radius [um]'\n"; } + else if (plt == "actrw_mom1") + { + gp << "set title 'mean radius of actrw droplets [um]'\n"; + } + else if (plt == "actrw_mom2") + { + gp << "set title 'variance of actrw droplet size distribution [um^2] * 1e3'\n"; + } + else if (plt == "actrw_sd") + { + gp << "set title 'standard deviation of actrw droplet size distribution [um]'\n"; + } + else if (plt == "actrw_rw_cl_conc") + { + gp << "set title 'concentration of actrw droplet in cloudy cells [cm^-3]'\n"; + } + else if (plt == "actrw_rw_cl") + { + gp << "set title 'mean radius of actrw droplet in cloudy cells [um]'\n"; + } + else if (plt == "disp_r") + { + gp << "set title 'relative dispersion of r_{drop} [1]'\n"; + } + + + } diff --git a/drawbicyc/include/plot_prof.hpp b/drawbicyc/include/plot_prof.hpp index 2e7d086..f825827 100644 --- a/drawbicyc/include/plot_prof.hpp +++ b/drawbicyc/include/plot_prof.hpp @@ -28,7 +28,7 @@ void plot_profiles(Plotter_t plotter, Plots plots, std::string type, const bool } Gnuplot gp; string file = plotter.file + "_" + type + "_profiles_" + prof_start_s + "_" + prof_end_s+".svg"; - init_prof(gp, file, 3, 5); + init_prof(gp, file, 2, 4); string prof_file = plotter.file + "_" + type + "_profiles_" + prof_start_s + "_" + prof_end_s+".dat"; std::ofstream oprof_file(prof_file); @@ -69,6 +69,9 @@ void plot_profiles(Plotter_t plotter, Plots plots, std::string type, const bool typename Plotter_t::arr_t res(rhod.shape()); typename Plotter_t::arr_t res_tmp(rhod.shape()); typename Plotter_t::arr_t res_tmp2(rhod.shape()); + //typename Plotter_t::arr_t m1(rhod.shape()); + //typename Plotter_t::arr_t m2(rhod.shape()); + //typename Plotter_t::arr_t m0(rhod.shape()); blitz::Array res_prof_sum(n["z"]); // profile interpolate to the uniform grid summed over timesteps blitz::Array res_prof(n["z"]); // profile interpolate to the uniform grid blitz::Array res_prof_hlpr(n["z"]); // actual profile @@ -84,13 +87,44 @@ void plot_profiles(Plotter_t plotter, Plots plots, std::string type, const bool std::cout << at * n["outfreq"] << std::endl; res = 0; +// if (plt == "rliq") +// { +// // liquid water content +// res += plotter.h5load_ra_timestep(at * n["outfreq"]) * 1e3; // aerosol +// res += plotter.h5load_rc_timestep(at * n["outfreq"]) * 1e3; // cloud +// res += plotter.h5load_rr_timestep(at * n["outfreq"]) * 1e3; // rain +// res_prof_hlpr = plotter.horizontal_mean(res); // average in x +// } if (plt == "rliq") { - // liquid water content - // res += plotter.h5load_ra_timestep(at * n["outfreq"]) * 1e3; // aerosol + // liquid water content in cloudy cells res += plotter.h5load_rc_timestep(at * n["outfreq"]) * 1e3; // cloud res += plotter.h5load_rr_timestep(at * n["outfreq"]) * 1e3; // rain - res_prof_hlpr = plotter.horizontal_mean(res); // average in x + + typename Plotter_t::arr_t snap(plotter.h5load_rc_timestep(at * n["outfreq"])); + res_tmp = iscloudy_rc_rico(snap); + + res *= res_tmp; + prof_tmp = plotter.horizontal_sum(res_tmp); + res_prof_hlpr = where(prof_tmp > 0, plotter.horizontal_sum(res) / prof_tmp, 0); + } + if (plt == "actrw_rw_cl_conc") + { + // concentration of activated droplets (r > rc) in cloudy cells + { + typename Plotter_t::arr_t snap(plotter.h5load_rc_timestep(at * n["outfreq"])); + res_tmp = iscloudy_rc_rico(snap); + } + { + auto tmp = plotter.h5load_timestep("actrw_rw_mom0", at * n["outfreq"]); + typename Plotter_t::arr_t snap(tmp); + res_tmp2 = snap; + res_tmp2 *= rhod / 1e6; // per cm^3 + } + // cloudy only + prof_tmp = plotter.horizontal_sum(res_tmp); + res_tmp2 *= res_tmp; + res_prof_hlpr = where(prof_tmp > 0 , plotter.horizontal_sum(res_tmp2) / prof_tmp, 0); } if (plt == "gccn_conc") { @@ -930,6 +964,52 @@ void plot_profiles(Plotter_t plotter, Plots plots, std::string type, const bool res_prof_hlpr = plotter.horizontal_mean(res); // average in x } + if (plt == "actrw_mom1") + { + // mean radius of actrw droplets + res = plotter.h5load_timestep("actrw_rw_mom1", at * n["outfreq"]) * 1e6; + res_tmp = plotter.h5load_timestep("actrw_rw_mom0", at * n["outfreq"]); + res = where(res > 0 , res / res_tmp, res); + res_prof_hlpr = plotter.horizontal_mean(res); // average in x + } + else if (plt == "actrw_mom2") + { + // variance of actrw droplet size distribution + res = plotter.h5load_timestep("actrw_rw_mom2", at * n["outfreq"]) * 1e6 * 1e3; // convert to (um)^2 * 1e3 + res_tmp = plotter.h5load_timestep("actrw_rw_mom0", at * n["outfreq"]); + res = where(res > 0 , res / res_tmp, res); + res_prof_hlpr = plotter.horizontal_mean(res); // average in x + } + else if (plt == "disp_r") + { + // relative dispersion (std dev / mean) of droplet radius distribution averaged over cells away from walls + //typename Plotter_t::arr_t m0(plotter.nowall(typename Plotter_t::arr_t(plotter.h5load_timestep("cloud_rw_mom0", at * n["outfreq"])), distance_from_walls)); + typename Plotter_t::arr_t m0(plotter.h5load_timestep("actrw_rw_mom0", at * n["outfreq"])); + typename Plotter_t::arr_t m1(plotter.h5load_timestep("actrw_rw_mom1", at * n["outfreq"])); + typename Plotter_t::arr_t m2(plotter.h5load_timestep("actrw_rw_mom2", at * n["outfreq"])); + + //calculate mean radius, store in m1 + m1 = where(m0 > 0, + m1 / m0, 0.); + + // calculate stddev of radius, store in m2 + m2 = where(m0 > 0, + m2 / m0 - m1 / m0 * m1 / m0, 0.); + + // might be slightly negative due to numerical errors + m2 = where(m2 < 0, 0, m2); + m2 = sqrt(m2); // sqrt(variance) + + res = where(m1 > 0 , m2 / m1 , 0); + + typename Plotter_t::arr_t snap(plotter.h5load_rc_timestep(at * n["outfreq"])); + res_tmp = iscloudy_rc_rico(snap); + res *= res_tmp; + prof_tmp = plotter.horizontal_sum(res_tmp); + + res_prof_hlpr = where(prof_tmp > 0, plotter.horizontal_sum(res) / prof_tmp, 0); + } + // ====================================================================================== diff --git a/drawbicyc/include/plot_series.hpp b/drawbicyc/include/plot_series.hpp index 2167e77..d6c79cb 100644 --- a/drawbicyc/include/plot_series.hpp +++ b/drawbicyc/include/plot_series.hpp @@ -1350,12 +1350,12 @@ void plot_series(Plotter_t plotter, Plots plots, std::string type) try { // cloud fraction (cloudy if N_c > 20/cm^3) - typename Plotter_t::arr_t snap(plotter.h5load_timestep("cloud_rw_mom0", at * n["outfreq"])); + typename Plotter_t::arr_t snap(plotter.h5load_timestep("actrw_rw_mom0", at * n["outfreq"])); snap *= rhod; // b4 it was specific moment snap /= 1e6; // per cm^3 snap = iscloudy(snap); // cloudiness mask - typename Plotter_t::arr_t snap_m0(plotter.h5load_timestep("cloud_rw_mom0", at * n["outfreq"])); - typename Plotter_t::arr_t snap_m1(plotter.h5load_timestep("cloud_rw_mom1", at * n["outfreq"])*1e6); // in microns + typename Plotter_t::arr_t snap_m0(plotter.h5load_timestep("actrw_rw_mom0", at * n["outfreq"])); + typename Plotter_t::arr_t snap_m1(plotter.h5load_timestep("actrw_rw_mom1", at * n["outfreq"])*1e6); // in microns snap_m0 *= snap; snap_m1 *= snap; auto tot_m0 = blitz::sum(snap_m0); diff --git a/drawbicyc/include/plots.hpp b/drawbicyc/include/plots.hpp index 1fa2b4a..9f48fc0 100644 --- a/drawbicyc/include/plots.hpp +++ b/drawbicyc/include/plots.hpp @@ -60,7 +60,7 @@ class Plots fields.insert(fields.end(), fields_Lasher_Trapp.begin(), fields_Lasher_Trapp.end()); } else if(type == "cumulus_congestus") { - // profs.insert(profs.end(), profs_Cumulus_Congestus.begin(), profs_Cumulus_Congestus.end()); + profs.insert(profs.end(), profs_Cumulus_Congestus.begin(), profs_Cumulus_Congestus.end()); series.insert(series.end(), series_Cumulus_Congestus.begin(), series_Cumulus_Congestus.end()); // fields.insert(fields.end(), fields_Cumulus_Congestus.begin(), fields_Cumulus_Congestus.end()); }