diff --git a/examples/README b/examples/README index 691846b..6cea825 100644 --- a/examples/README +++ b/examples/README @@ -69,7 +69,10 @@ show_pull_histogram.py Displays pull histograms generated by the cspull script. show_response.py - Display one or several CTA response functions. + Displays one or several CTA response functions. + +show_sensitivity.py + Displays sensitivity information generated by the cssens script. show_ts_distribution.py Displays Test Statistic histograms generated using the cstsdist diff --git a/examples/show_sensitivity.py b/examples/show_sensitivity.py new file mode 100644 index 0000000..1cd538a --- /dev/null +++ b/examples/show_sensitivity.py @@ -0,0 +1,115 @@ +#! /usr/bin/env python +# ========================================================================== +# This script displays a differential sensitivity plot generated by cssens. +# +# Copyright (C) 2015 Juergen Knoedlseder +# +# This program is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program. If not, see . +# +# ========================================================================== +import sys +import math +import gammalib +try: + import matplotlib.pyplot as plt +except: + sys.exit("This script needs matplotlib") + + +# ===================== # +# Read sensitivity data # +# ===================== # +def read_sensitivity(filename): + """ + Read sensitivity information from CSV file. + """ + # Read filename + csv = gammalib.GCsv(filename,',') + + # Create dictionnary + sensitivity = {} + for column in range(csv.ncols()): + name = csv[0,column] + values = [] + for row in range(csv.nrows()-1): + values.append(float(csv[row+1,column])) + sensitivity[name] = values + + # Add linear energy values + if sensitivity.has_key("loge"): + name = 'energy' + values = [] + for value in sensitivity["loge"]: + values.append(math.pow(10.0, value)) + sensitivity[name] = values + + # Return + return sensitivity + + +# ===================== # +# Show sensitivity data # +# ===================== # +def show_sensitivity(sensitivity, title="Sensitivity"): + """ + Plot sensitivity. + """ + # Create figure + plt.figure() + + # Set plot attributes + plt.title(title) + plt.loglog() + plt.grid() + + # Show differential sensitivity + plt.plot(sensitivity['energy'], sensitivity['eflux'], 'ro-') + + # Set labels + plt.xlabel("Energy (TeV)") + plt.ylabel(r"Sensitivity (erg cm$^{-2}$ s$^{-1}$)") + + # Show plot + plt.show() + + # Return + return + + +# ======================== # +# Main routine entry point # +# ======================== # +if __name__ == '__main__': + """ + """ + # Print usage information + usage = "Usage: show_sensitivity [filename]" + if len(sys.argv) < 1: + sys.stdout.write(usage+"\n") + sys.exit() + + # Extract parameters + if len(sys.argv) == 1: + filename = 'sensitivity.dat' + else: + filename = sys.argv[1] + + # Set title + title = 'Differential sensitivity ('+filename+')' + + # Read sensitivity data + sensitivity = read_sensitivity(filename) + + # Show sensitivity data + show_sensitivity(sensitivity, title=title)