IPython solution:
hdulist = pyfits.open('[SOURCE]_filtered_gti.fits')
events=hdulist[1].data
hist(log10(events['ENERGY']),100)
yscale('log')
xlabel('log10(Energy/MeV)')
Note that we plot the histogram with both axes in log-scale in this example.
hist(events['TIME'],100)
xlabel('Time (s)')