Using interpolated images for spectrograms would be nice
Created by: jrsmith02
For visualizing short, quickly varying things, using spectrograms with high overlap and interpolation, and plotting them as images is quite useful. For scatter plotting, I have been using NonUniformImage (note: I thought this allows both lin and log y-axis, but perhaps that depends on python version because currently at LLO I'm only able to make lin).
Here is an example of scattering fringes at LLO:
And here is a complete, sloppy, python script to reproduce this result.
from gwpy.timeseries import TimeSeries
from gwpy.time import from_gps
from numpy import sqrt, array, arange, nonzero, argmax, log10
from pylab import savefig, specgram
import sys, string
import matplotlib.pyplot as plt
from matplotlib.image import NonUniformImage
import matplotlib.mlab as mlab
from os import makedirs, path
# read input argumentsi (update this to use argparse)
ifo = 'L1'
start_time = 1129118930
startutc = str(from_gps(start_time)) # get UTC time
dur = 30 # duration in seconds
end_time = int(start_time) + int(dur)
witness_base = 'LSC-SRCL_IN1_DQ'
witness_chan = ifo + ':' + witness_base
# load timeseries for witness channel
witness=TimeSeries.fetch(witness_chan, start_time, start_time+dur, verbose=True)
witness=witness.highpass(10,gpass=3) # highpass the witness data
# Calculate spectrogram
secsPerFFT = 1.0 # Hz
overlap = 0.9 # fractional overlap
Fs = witness.sample_rate.value
NFFT = int(round(Fs*secsPerFFT))
noverlap = int(round(overlap*NFFT))
Pxx, freq, t, im = specgram(witness.value,NFFT=NFFT,Fs=witness.sample_rate.value,noverlap=noverlap,scale_by_freq='magnitude',detrend=mlab.detrend_linear,window=mlab.window_hanning)
fig = plt.figure(figsize=(12,12))
ax1 = fig.add_subplot(111)
# Plot Spectrogram
im1 = NonUniformImage(ax1, interpolation='bilinear',extent=(min(t),max(t),10,55),cmap='jet')
im1.set_data(t,freq,20.0*log10(Pxx))
im1.set_clim(-160,-20)
ax1.images.append(im1)
cbar1 = fig.colorbar(im1)
plt.yscale('linear')
plt.ylabel('Frequency [Hz]',fontsize=24)
xlab = 'Time [seconds] from ' + str(startutc) + ' UTC'
plt.xlabel(xlab,fontsize=24)
plt.title(witness_chan,fontsize=24)
plt.xlim(0,max(t))
plt.ylim(1,25)
# save figure
dir = ifo + '-' + str(start_time) + '-' + str(dur) + '/'
if not path.exists(dir):
makedirs(dir)
filename = dir + ifo + '-' + witness_base + '-' + str(start_time) + '-' + str(dur) + '.png'
savefig(filename)
plt.close(fig)