Commit 97cf14031ff7d4f02715cb8977282950efdbf815

Authored by Jürgen Knödlseder
1 parent c6aa9b8e

Add comments

Showing 1 changed file with 15 additions and 5 deletions
examples/make_spectrum.py
@@ -125,17 +125,21 @@ def plot_spectrum(spectrum, plotfile): @@ -125,17 +125,21 @@ def plot_spectrum(spectrum, plotfile):
125 ul_ed_engs = [] 125 ul_ed_engs = []
126 ul_eu_engs = [] 126 ul_eu_engs = []
127 ul_flux = [] 127 ul_flux = []
128 -  
129 - # Loop over rows of the file 128 +
  129 + # Get number of spectral points
130 nrows = table.nrows() 130 nrows = table.nrows()
  131 +
  132 + # Loop over rows of the file
131 for row in range(nrows): 133 for row in range(nrows):
132 134
133 - # Get TS 135 + # Get TS value, flux and flux error
134 ts = c_ts.real(row) 136 ts = c_ts.real(row)
135 flx = c_flux.real(row) 137 flx = c_flux.real(row)
136 e_flx = c_eflux.real(row) 138 e_flx = c_eflux.real(row)
137 139
138 - # Switch 140 + # If the TS > 9 and the flux error is smaller than the flux then
  141 + # take the fitted flux and flux error and append them to the
  142 + # arrays
139 if ts > 9.0 and e_flx < flx: 143 if ts > 9.0 and e_flx < flx:
140 144
141 # Add information 145 # Add information
@@ -145,7 +149,8 @@ def plot_spectrum(spectrum, plotfile): @@ -145,7 +149,8 @@ def plot_spectrum(spectrum, plotfile):
145 eu_engs.append(c_eu.real(row)) 149 eu_engs.append(c_eu.real(row))
146 e_flux.append(c_eflux.real(row)) 150 e_flux.append(c_eflux.real(row))
147 151
148 - # 152 + # ... otherwise take the upper limit and append it to the upper
  153 + # limit array
149 else: 154 else:
150 155
151 # Add information 156 # Add information
@@ -161,10 +166,15 @@ def plot_spectrum(spectrum, plotfile): @@ -161,10 +166,15 @@ def plot_spectrum(spectrum, plotfile):
161 # Plot the spectrum 166 # Plot the spectrum
162 plt.loglog() 167 plt.loglog()
163 plt.grid() 168 plt.grid()
  169 +
  170 + # Plot error bars for flux points and upper limits for upper limit
  171 + # points
164 plt.errorbar(energies, flux, yerr=e_flux, xerr=[ed_engs, eu_engs], 172 plt.errorbar(energies, flux, yerr=e_flux, xerr=[ed_engs, eu_engs],
165 fmt='ro') 173 fmt='ro')
166 plt.errorbar(ul_energies, ul_flux, xerr=[ul_ed_engs, ul_eu_engs], 174 plt.errorbar(ul_energies, ul_flux, xerr=[ul_ed_engs, ul_eu_engs],
167 yerr=1.0e-11, uplims=True, fmt='ro') 175 yerr=1.0e-11, uplims=True, fmt='ro')
  176 +
  177 + # Add plot labels
168 plt.xlabel('Energy (TeV)') 178 plt.xlabel('Energy (TeV)')
169 plt.ylabel(r'E$^{2}$dN/dE (erg cm$^{-2}$ s$^{-1}$)') 179 plt.ylabel(r'E$^{2}$dN/dE (erg cm$^{-2}$ s$^{-1}$)')
170 180