#!/usr/bin/env python
# plot_flame_dist.py
# By Nathan C. Hearn
# September 24, 2008
#
# Plotting utility for flame_dist_analysis output
#
# Requires QuickFlash/utils to be in PYTHONPATH environment variable
try :
import psyco
except ImportError :
pass
data_start_keyword = "DistribData"
default_filename_extension = "png"
png_filename_extension = "png"
png_matplotlib_backend = "Agg"
svg_filename_extension = "svg"
svg_matplotlib_backend = "SVG"
ps_filename_extension = "ps"
ps_matplotlib_backend = "PS"
pdf_filename_extension = "pdf"
pdf_matplotlib_backend = "PDF"
mean_value_suffix = "_Mean_Value"
dispersion_suffix = "_Dispersion"
rms_suffix = "_RMS_Difference"
density_suffix = "_Density" ## Initially, a double-underscore was used
area_suffix = "_Area"
def read_file(filename) :
from plot_gen import read_keywords_table
## Read file data and attributes
keyword_values, table_dict \
= read_keywords_table(filename, table_start_line=data_start_keyword,
ignore_leading_comment_mark=True)
return keyword_values, table_dict
def build_plot(output_file, matplotlib_backend, plot_title, x_axis_title,
min_x, max_x, y_axis_title, pos_array_list, data_array_list,
legend_list) :
from plot_gen import Plot
num_lists = len(pos_array_list)
if (len(data_array_list) != num_lists) :
raise RuntimeError("List count mismatch")
if (len(legend_list) != num_lists) :
raise RuntimeError("List count mismatch")
p = Plot(plot_title, x_axis_title, y_axis_title)
data_bounds_set = False
min_y = None
max_y = None
for list_index in range(num_lists) :
data_array = data_array_list[list_index]
local_min_y = min(data_array)
local_max_y = max(data_array)
if (not data_bounds_set) :
min_y = local_min_y
max_y = local_max_y
data_bounds_set = True
else :
if (local_min_y < min_y) :
min_y = local_min_y
if (local_max_y > max_y) :
max_y = local_max_y
p.add_plot(pos_array_list[list_index], data_array,
legend_text=legend_list[list_index])
p.set_limits(min_x=min_x, max_x=max_x, min_y=min_y, max_y=max_y)
p.show_plot(output_file, matplotlib_backend=matplotlib_backend)
def scale_data(data_array, scale_factor) :
new_array = list()
for elem in data_array :
new_array.append(elem * scale_factor)
return new_array
def offset_data(data_array, offset_value) :
new_array = list()
for elem in data_array :
new_array.append(elem + offset_value)
return new_array
def main() :
from sys import exit, stderr
from os import path
from optparse import OptionParser
## Set up the command line parser
usage = str("%prog [ options ] output_prefix input_file1[,legend_1] ...")
version = str("%prog 0.1")
argparser = OptionParser(usage, version=version)
output_extension_str = str("output-extension")
argparser.add_option("--" + output_extension_str, type="string",
dest="output_extension",
help="Filename extension for output files")
argparser.set_defaults(output_format=None)
plot_title_str = str("plot-title")
argparser.add_option("--" + plot_title_str, type="string",
dest="plot_title",
help="Text for plot titles")
argparser.set_defaults(plot_title=None)
time_legend_str = str("time-legend")
argparser.add_option("--" + time_legend_str, action="store_true",
dest="time_legend",
help="Add time to plot legend")
argparser.set_defaults(time_legend=False)
center_flame_str = str("center-flame")
argparser.add_option("--" + center_flame_str, type="string",
dest="center_flame",
help="Center flame (choose mass or mass-fraction)")
argparser.set_defaults(center_flame=None)
rescale_pos_str = str("position-rescale")
argparser.add_option("--" + rescale_pos_str, type="string",
dest="position_rescale",
help="Comma-separated scale factor and unit name")
argparser.set_defaults(position_rescale=None)
opts, args = argparser.parse_args()
argc = len(args)
min_args = int(2)
if (argc < min_args) :
stderr.write("\n")
argparser.print_help(file=stderr)
stderr.write("\n")
exit(1)
## Get the required arguments
argptr = int(0)
user_output_prefix = str(args[argptr])
argptr += 1
file_list = list()
legend_list = list()
while (argptr < argc) :
file_legend_entry = str(args[argptr])
argptr += 1
file_legend_list = file_legend_entry.split(",")
num_words = len(file_legend_list)
if (num_words < 1) :
stderr.write("\nError: Unable to parse data file list\n\n")
exit(-1)
filename = file_legend_list[0]
legend_text = None
if (num_words < 2) :
legend_text_end = filename.rfind("_hdf5")
if (legend_text_end < 1) :
legend_text_end = len(filename)
legend_text = filename[0:legend_text_end]
elif (num_words == 2) :
legend_text = file_legend_list[1]
else :
stderr.write("\nError: Improper format for filename-legend pair [ "
+ file_legend_entry + " ]\n\n")
exit(-1)
file_list.append(filename)
legend_list.append(legend_text)
num_files = len(file_list)
## Read the optional arguments
output_extension = opts.output_extension
if (output_extension is None) :
output_extension = default_filename_extension
else :
if (output_extension == png_filename_extension) :
matplotlib_backend = png_matplotlib_backend
elif (output_extension == svg_filename_extension) :
matplotlib_backend = svg_matplotlib_backend
elif (output_extension == ps_filename_extension) :
matplotlib_backend = ps_matplotlib_backend
elif (output_extension == pdf_filename_extension) :
matplotlib_backend = pdf_matplotlib_backend
else :
stderr.write("\nError: Output extension [ " + output_extension
+ " ] not recognized\n\n")
exit(-1)
user_plot_title = opts.plot_title
plot_title_prefix = str()
if (user_plot_title is not None) :
plot_title_prefix = user_plot_title + " - "
add_time_legend = opts.time_legend
center_flame_method = opts.center_flame
recenter_flame = False
if (center_flame_method is not None) :
if (not ((center_flame_method == "mass")
or (center_flame_method == "mass-fraction"))) :
stderr.write("\nError: Flame centering method [ "
+ center_flame_method + " ] not recognized\n\n")
exit(-1)
recenter_flame = True
position_scale = 1.0
title_unit_suffix = " (cm)"
if (opts.position_rescale is not None) :
words = opts.position_rescale.split(",")
num_words = len(words)
if ((num_words < 1) or (num_words > 2)) :
stderr.write("\nError: Position scale requires comma-separated "
+ "scale factor and (optional) unit name\n\n")
exit(-1)
position_scale = float(words[0])
if (num_words > 1) :
title_unit_suffix = " (" + words[1] + ")"
else :
title_unit_suffix = str()
abs_pos_x_title = "Position" + title_unit_suffix
## Get other quantities
output_prefix = user_output_prefix
if (len(output_prefix) < 1) :
stderr.write("\nError: Improper output_prefix\n\n")
exit(-1)
if (path.split(output_prefix)[1] != "") :
if (output_prefix[-1] != "_") :
output_prefix += "_"
output_suffix = "." + output_extension
## Read the files
table_list = list()
sim_times_list = list()
flame_center_mass_list = list()
flame_center_mass_fraction_list = list()
area_fractions_present_list = list() ## booleans
area_min_flame_value_list = list()
area_min_flame_value_str_list = list()
area_max_flame_values_list = list()
area_max_flame_values_str_list = list()
flame_min_pos = None
flame_max_pos = None
bin_array_min_pos = None
bin_array_max_pos = None
user_dataset_list = None
rms_ref_value = None
time_data_filename = output_prefix + "time_data.csv"
time_data_file = open(time_data_filename, "w")
time_data_file.write("Filename" + "\t" + "Sim_Time" + "\t"
+ "Flame_Center_Mass_Fraction" + "\t"
+ "Flame_Center_Mass" + "\t" + "Min_Flame_Pos"
+ "\t" + "Max_Flame_Pos" + "\n")
for filename in file_list :
keyword_value, table = read_file(filename)
sim_time = float(keyword_value.get_first_value("SimTime"))
raw_flame_center_mass \
= float(keyword_value.get_first_value("FlameCenterMass"))
flame_center_mass \
= raw_flame_center_mass * position_scale
raw_flame_center_mass_fraction \
= float(keyword_value.get_first_value("FlameCenterMassFract"))
flame_center_mass_fraction \
= raw_flame_center_mass_fraction * position_scale
raw_local_flame_min_pos \
= float(keyword_value.get_first_value("MinFlamePos"))
local_flame_min_pos \
= raw_local_flame_min_pos * position_scale
raw_local_flame_max_pos \
= float(keyword_value.get_first_value("MaxFlamePos"))
local_flame_max_pos \
= raw_local_flame_max_pos * position_scale
time_data_file.write(filename + "\t" + str(sim_time) + "\t"
+ str(raw_flame_center_mass_fraction) + "\t"
+ str(raw_flame_center_mass) + "\t"
+ str(raw_local_flame_min_pos) + "\t"
+ str(raw_local_flame_max_pos) + "\n")
if (flame_min_pos is None) :
flame_min_pos = local_flame_min_pos
flame_max_pos = local_flame_max_pos
else :
if (local_flame_min_pos < flame_min_pos) :
flame_min_pos = local_flame_min_pos
if (local_flame_max_pos > flame_max_pos) :
flame_max_pos = local_flame_max_pos
local_bin_array_min_pos \
= float(keyword_value.get_first_value("BinArrayMinPos")) \
* position_scale
local_bin_array_max_pos \
= float(keyword_value.get_first_value("BinArrayMaxPos")) \
* position_scale
if (bin_array_min_pos is None) :
bin_array_min_pos = local_bin_array_min_pos
bin_array_max_pos = local_bin_array_max_pos
else :
if (local_bin_array_min_pos < bin_array_min_pos) :
bin_array_min_pos = local_bin_array_min_pos
if (local_bin_array_max_pos > bin_array_max_pos) :
bin_array_max_pos = local_bin_array_max_pos
distrib_datasets_index \
= keyword_value.find_first_pair("DistribDatasets")
if (distrib_datasets_index is not None) :
local_dataset_str \
= keyword_value[distrib_datasets_index].get_value()
local_dataset_list = local_dataset_str.split(",")
local_dataset_list.sort()
if (user_dataset_list is None) :
user_dataset_list = local_dataset_list
elif (local_dataset_list != user_dataset_list) :
stderr.write("\nError: Incompatible user dataset lists\n\n")
exit(-2)
elif (user_dataset_list is None) :
user_dataset_list = list()
elif (len(user_dataset_list) > 0) :
stderr.write("\nError: Incompatible user dataset lists\n\n")
rms_diff_index \
= keyword_value.find_first_pair("RootMeanSqDiffRefValue")
if (rms_diff_index is not None) :
local_ref_value = float(keyword_value[rms_diff_index].get_value())
if (rms_ref_value is None) :
rms_ref_value = local_ref_value
elif (local_ref_value != rms_ref_value) :
stderr.write("\nError: Incompatible RMS reference value\n\n")
table_list.append(table)
sim_times_list.append(sim_time)
flame_center_mass_list.append(flame_center_mass)
flame_center_mass_fraction_list.append(flame_center_mass_fraction)
area_flame_min_index \
= keyword_value.find_first_pair("AreaFlameMinvalue")
if (area_flame_min_index is not None) :
area_fractions_present_list.append(True)
min_flame_value_str = keyword_value[area_flame_min_index]
max_flame_values_str \
= keyword_value.get_first_value("AreaFlameMaxValues")
max_flame_values_str_list = max_flame_values_str.split(",")
max_flame_values_str_list.sort()
max_flame_values = list()
for value_string in max_flame_values_str_list :
max_flame_values.append(float(value_string))
area_min_flame_value_list.append(float(min_flame_value_str))
area_min_flame_value_str_list.append(min_flame_value_str)
area_max_flame_values_list.append(max_flame_values)
area_max_flame_values_str_list.append(max_flame_values_str)
else :
area_fractions_present_list.append(False)
area_min_flame_value_list.append(None)
area_min_flame_value_str_list.append(None)
area_max_flame_values_list.append(None)
area_max_flame_values_str_list.append(None)
if (add_time_legend) :
for file_index in range(num_files) :
sim_time = sim_times_list[file_index]
new_legend = legend_list[file_index] + " " \
+ ("%(time)3.2f" % { "time" : sim_time })
legend_list[file_index] = new_legend
## Get the position lists
pos_array_list = list()
pos_array_center_mass_list = list()
pos_array_center_massfract_list = list()
for file_index in range(num_files) :
orig_pos_array = scale_data(table_list[file_index]["Bin_Center"],
position_scale)
center_mass = flame_center_mass_list[file_index]
center_massfract = flame_center_mass_fraction_list[file_index]
center_mass_offset = -1.0 * center_mass
center_massfract_offset = -1.0 * center_massfract
pos_array_list.append(orig_pos_array)
pos_array_mass = offset_data(orig_pos_array, center_mass_offset)
pos_array_center_mass_list.append(pos_array_mass)
pos_array_massfract = offset_data(orig_pos_array,
center_massfract_offset)
pos_array_center_massfract_list.append(pos_array_massfract)
## Density distribution
dens_output_file = output_prefix + "density" + output_suffix
dens_data_list = list()
for data_table in table_list :
dens_data_list.append(data_table["Mean_Density"])
dens_plot_title = plot_title_prefix + "Mean Density"
dens_plot_y_title = "Density (g / cc)"
build_plot(dens_output_file, matplotlib_backend, dens_plot_title,
abs_pos_x_title, flame_min_pos, flame_max_pos,
dens_plot_y_title, pos_array_list, dens_data_list, legend_list)
## Burned mass fraction
burned_massfract_output_file = output_prefix + "burned_mass_fraction" \
+ output_suffix
burned_massfract_list = list()
for data_table in table_list :
burned_massfract_list.append(data_table["Burned_Mass_Fraction"])
burned_massfract_title = plot_title_prefix + "Burned Mass Fraction"
burned_massfract_y_title = "Mass Fraction"
build_plot(burned_massfract_output_file, matplotlib_backend,
burned_massfract_title, abs_pos_x_title, flame_min_pos,
flame_max_pos, burned_massfract_y_title, pos_array_list,
burned_massfract_list, legend_list)
## User dataset distribution
if (user_dataset_list is not None) :
for user_dataset in user_dataset_list :
## Mean value
mean_output_file = output_prefix + user_dataset + "_mean" \
+ output_suffix
mean_name = user_dataset + mean_value_suffix
mean_value_list = list()
for data_table in table_list :
mean_value_list.append(data_table[mean_name])
mean_title = plot_title_prefix + "Mean Value"
mean_y_title = user_dataset + " Mean Value"
build_plot(mean_output_file, matplotlib_backend, mean_title,
abs_pos_x_title, bin_array_min_pos, bin_array_max_pos,
mean_y_title, pos_array_list, mean_value_list,
legend_list)
## Dispersion
dispersion_output_file = output_prefix + user_dataset \
+ "_dispersion" + output_suffix
dispersion_name = user_dataset + dispersion_suffix
dispersion_list = list()
for data_table in table_list :
dispersion_list.append(data_table[dispersion_name])
dispersion_title = plot_title_prefix + "Dispersion"
dispersion_y_title = user_dataset + " Dispersion"
build_plot(dispersion_output_file, matplotlib_backend,
dispersion_title, abs_pos_x_title, bin_array_min_pos,
bin_array_max_pos, dispersion_y_title, pos_array_list,
dispersion_list, legend_list)
## RMS Value
rms_output_file = output_prefix + user_dataset + "_rms_diff" \
+ output_suffix
rms_name = user_dataset + rms_suffix
rms_list = list()
for data_table in table_list :
rms_list.append(data_table[rms_name])
rms_title = plot_title_prefix
rms_y_title = user_dataset
if (rms_ref_value is None) :
rms_title += "Population Standard Deviation"
rms_y_title += " Population Standard Deviation"
else :
rms_title += "RMS Difference from " + str(rms_ref_value)
rms_y_title += " RMS Difference from " + str(rms_ref_value)
build_plot(rms_output_file, matplotlib_backend,
rms_title, abs_pos_x_title, bin_array_min_pos,
bin_array_max_pos, rms_y_title, pos_array_list,
rms_list, legend_list)
## Flame area distribution
for file_index in range(num_files) :
if (area_fractions_present_list[file_index]) :
pass
if (__name__ == "__main__") :
main()