Ignore:
Timestamp:
09/24/14 12:50:19 (10 years ago)
Author:
dneise
Message:
deleted private stuff, that accidentally found it's way into pyfact.py
File:
1 edited

Legend:

Unmodified
Added
Removed
  • fact/tools/pyscripts/new_pyfact/pyfact.py

    r17972 r17973  
    88import numpy as np
    99import ROOT
    10 import pylab
    11 
    12 from scipy.signal import firwin, iirdesign, lfilter
    1310
    1411########## BUILDING OF THE SHARED OBJECT FILES ###############################
     
    248245    f.next(3)
    249246    print "Event Id:", f.cols['EventNum']
    250 
    251     import matplotlib.pyplot as plt
    252     plt.ion()
    253     """
    254     for i in range(f.header['NPIX']):
    255         plt.cla()
    256         plt.title("Event %d"%(f.cols['EventNum']))
    257         plt.plot(f.cols['CalibData2D'][i], '.:', label='pixel %d'%i)
    258         plt.legend()
    259         answer = raw_input('anykey for next pixel; "q" to quit. :')
    260         if 'q' in answer.lower():
    261             break
    262     """
    263 
    264 
    265 
    266     from scipy import signal
    267     from scipy.signal import freqs, iirfilter
    268     #b, a = iirfilter(4, [1, 10], 1, 60, analog=True, ftype='cheby1')
    269     #b, a = signal.iirfilter(17, [50, 200], rs=60, btype='band',
    270     #                                analog=True, ftype='cheby2')
    271     #w, h = signal.freqs(b, a, 1000)
    272     #w, h = freqs(b, a, worN=np.logspace(-1, 2, 1000))
    273 
    274     b, a = iirfilter(N=4, # the order of the filter
    275             Wn = [0.15, 0.25], # critical frequencies
    276             rp = 1,       # maximum ripple in passband
    277             rs = 60,      # minimum attenuation in stopband
    278             btype='bandstop', # {'bandpass', 'lowpass', 'highpass', 'bandstop'}
    279             #analog=True,
    280             ftype='cheby1')
    281     #w, h = freqs(b, a, worN=np.logspace(-1, 5, 1000))
    282     w, h = freqs(b, a)
    283     plt.figure()
    284     plt.semilogx(w, 20 * np.log10(abs(h)))
    285     plt.xlabel('Frequency')
    286     plt.ylabel('Amplitude response [dB]')
    287     plt.grid()
    288     plt.show()
    289 
    290 
    291 
    292 
    293 
    294     fig1, (ax11,ax12) = plt.subplots(nrows=2)
    295     fig2, (ax21,ax22) = plt.subplots(nrows=2)
    296 
    297     flt = firwin(13, [0.15,0.25])
    298     t = np.arange(f.cols['CalibData2D'].shape[1])/2.*1e-9
    299     fig3, (ax31, ax32) = plt.subplots(nrows=2)
    300     ax31.plot(flt)
    301     #w, h = freqs(flt, [1.0])
    302     w, h = freqs(flt, [1.0], worN=np.logspace(-1, 5, 1000))
    303     ax32.semilogx(w, 20 * np.log10(abs(h)))
    304     ax32.set_xlabel('Frequency')
    305     ax32.set_ylabel('Amplitude response [dB]')
    306     ax32.grid()
    307 
    308 
    309 
    310     for i in range(f.header['NPIX']):
    311         data = f.cols['CalibData2D'][i]
    312         fil_data = lfilter(flt, [1.0], data)
    313 
    314         ax11.cla()
    315         ax12.cla()
    316         ax11.set_title("Event %d"%(f.cols['EventNum']))
    317         ax11.plot(t, data, '.:', label='pixel %d'%i)
    318         ax11.legend()
    319         ax11.set_ylim(-20,40)
    320         ax12.plot(t[:fil_data.size], fil_data, '.:', label='pixel %d'%i)
    321         ax12.legend()
    322         ax12.set_ylim(-20,40)
    323 
    324         """
    325         Pxx, freqs, t, plot = pylab.specgram( f.cols['CalibData2D'][i],
    326                 NFFT=32, Fs=2000000000,
    327                 detrend=pylab.detrend_none,
    328                 window=pylab.window_hanning,
    329                 noverlap=16)
    330         """
    331         ax21.cla()
    332         ax22.cla()
    333         dat, freqs, bins, im = ax21.specgram(data,
    334                 NFFT=64, Fs=2000000000,
    335                 detrend=pylab.detrend_none,
    336                 window=pylab.window_hanning,
    337                 noverlap=60,
    338                 interpolation='nearest')
    339         ax21.axis('tight')
    340         ax21.set_ylabel("frequency [Hz]")
    341         ax21.set_xlabel("time [s]")
    342 
    343         dat2, freqs2, bins2, im2 = ax22.specgram(data,
    344                 NFFT=64, Fs=2000000000,
    345                 detrend=pylab.detrend_none,
    346                 window=pylab.window_hanning,
    347                 noverlap=60,
    348                 interpolation='nearest')
    349         ax22.axis('tight')
    350         ax22.set_ylabel("frequency [Hz]")
    351         ax22.set_xlabel("time [s]")
    352 
    353         plt.show()
    354 
    355         answer = raw_input('anykey for next pixel; "q" to quit. :')
    356         if 'q' in answer.lower():
    357             break
Note: See TracChangeset for help on using the changeset viewer.