How to use py-accessors

First off, we'll be using pylab to plot up the solutions



In [7]:

    
%pylab inline









    



Populating the interactive namespace from numpy and matplotlib

Loading calibration solutions from a measurement set

The default method for accessing calibrations is from a measurement set format (which is a directory).



In [114]:

    
from askap.accessors.calibration import Source



In [112]:

    
ls -ld *.tab/









    



drwxr-xr-x@ 7 ban115  at-aipsusr  238 24 Oct 15:43 caldata1.tab//
drwxr-xr-x@ 7 ban115  at-aipsusr  238 24 Oct 15:21 calibdata.tab//



In [113]:

    
src = Source('calibdata.tab')

The most recent solution ID is available as a property:



In [15]:

    
src.most_recent_solution_id









    Out[15]:





0

If you're just intrested in the most recent solution you can get the accessor with



In [103]:

    
ac = src.most_recent_accessor

Single value access

Raw access to the bandpass, gain, leakage, and full Jones matrix is available with the get() commands, e.g.



In [18]:

    
ant = 3
beam= 0
chan = 4



In [19]:

    
ac.get_gain(ant, beam)









    Out[19]:





<askap.accessors._accessors.JonesJTerm at 0x1099b5ed0>

Aaah, but waht is this freaky JonesJTerm? - it just has 4 useful properties



In [20]:

    
dir(ac.get_gain(ant, beam))









    Out[20]:





['__class__',
 '__delattr__',
 '__dict__',
 '__doc__',
 '__format__',
 '__getattribute__',
 '__hash__',
 '__init__',
 '__module__',
 '__new__',
 '__reduce__',
 '__reduce_ex__',
 '__repr__',
 '__setattr__',
 '__sizeof__',
 '__str__',
 '__subclasshook__',
 '__weakref__',
 'g1',
 'g1IsValid',
 'g2',
 'g2IsValid']



In [21]:

    
g = ac.get_gain(ant, beam)
print 'Gain was', g.g1, g.g2, 'Valid flags:', g.g1IsValid, g.g2IsValid









    



Gain was (1+0j) (1+0j) Valid flags: False False

The bandpass is similar, but needing a channel number too



In [22]:

    
bp  = ac.get_bandpass(ant, beam, chan)



In [23]:

    
print 'Bandpass was', bp.g1, bp.g2, 'Valid flags:', bp.g1IsValid, bp.g2IsValid









    



Bandpass was (-0.45052626729-0.704844295979j) (0.12813809514+0.942421197891j) Valid flags: True True

Leakage is a different data type



In [26]:

    
lk =  ac.get_leakage(ant, beam)
dir(lk)









    Out[26]:





['__class__',
 '__delattr__',
 '__dict__',
 '__doc__',
 '__format__',
 '__getattribute__',
 '__hash__',
 '__init__',
 '__module__',
 '__new__',
 '__reduce__',
 '__reduce_ex__',
 '__repr__',
 '__setattr__',
 '__sizeof__',
 '__str__',
 '__subclasshook__',
 '__weakref__',
 'd12',
 'd12IsValid',
 'd21',
 'd21IsValid']



In [28]:

    
print 'Leakage was', lk.d12, lk.d21, 'Valid flags:', lk.d12IsValid, lk.d21IsValid









    



 Leakage was 0j 0j Valid flags: False False

Jones is a 2x2 numpy array, unmasked



In [105]:

    
ac.get_jones(ant, beam, chan)









    Out[105]:





array([[-0.44971299-0.70564002j,  0.05768348-0.12381425j],
       [ 0.06421472+0.00320966j, -0.93064499-0.0232438j ]])

Sliced access

If you want access to the whole data as one big numpy block, try this:



In [31]:

    
bp_g1, bp_g2 = ac.bandpass[:,:,:]



In [32]:

    
bp_g1.shape









    Out[32]:





(36, 1, 8)

The bandpass property is sliceable, and returns g1 and g2 as a tuple. They are each masked numpy arrays of size Nant x Nbeam x Nchan. If you feel the need, you can plot them up too (although this case is not especially interesting):



In [61]:

    
subplot(2,2,1)
plot(abs(bp_g1[:, 0, :]).T)
ylabel('Amplitude')
title('G1')
subplot(2,2,3)
plot(degrees(angle(bp_g1[:, 0, :]).T))
ylabel('Phase (deg)')

xlabel('channel')
subplot(2,2,2)
plot(abs(bp_g2[:, 0, :]).T)
title('G2')
subplot(2,2,4)
plot(degrees(angle(bp_g2[:, 0, :]).T))
xlabel('channel')









    Out[61]:





<matplotlib.text.Text at 0x10d5612d0>

Maybe un-transpose that puppy up a little



In [106]:

    
subplot(2,2,1)
plot(abs(bp_g1[:, 0, :]))
ylabel('Amplitude')
title('G1')
subplot(2,2,3)
plot(degrees(angle(bp_g1[:, 0, :])))
ylabel('Phase (deg)')

xlabel('Antenna')
subplot(2,2,2)
plot(abs(bp_g2[:, 0, :]))
title('G2')
subplot(2,2,4)
plot(degrees(angle(bp_g2[:, 0, :])))
xlabel('Antenna')









    Out[106]:





<matplotlib.text.Text at 0x11052de90>

Hmm, looks like the same solution for every channel - yeah, we saw that before!

You can get the full Jones matrix, which is Nant x Nbeam x Nchan x 2 x 2 and is masked in all the right places



In [41]:

    
jones = ac.jones[:,:,:]



In [42]:

    
type(jones)









    Out[42]:





numpy.ma.core.MaskedArray



In [43]:

    
jones.shape









    Out[43]:





(36, 1, 8, 2, 2)

WARNING: If you don't specify the slice size, it might return a very big array with most of it masked out, e.g.



In [45]:

    
d12, d21= ac.leakage[:,:]



In [46]:

    
d12.shape









    Out[46]:





(128, 128)

Currently there's no way to fix this, as the underlying library only tells you you're out of bounds by throwing an exception



In [47]:

    
g1, g2 = ac.gain[:,:]



In [48]:

    
g1.shape









    Out[48]:





(128, 128)

Gain does the same for this dataset

Another measurement set, this time with gain in it



In [49]:

    
ac2 = Source('caldata1.tab').most_recent_accessor



In [53]:

    
g1, g2 = ac2.gain[:,:]



In [54]:

    
g1.shape









    Out[54]:





(36, 1)

Aah, look, the gain is a sensible size in this example



In [60]:

    
plot(abs(g1))
plot(abs(g2))
ylabel('Amplitude')
xlabel('Antenna')
figure()
plot(degrees(angle(g1)))
plot(degrees(angle(g2)))
ylabel('Phase (deg)')
xlabel('Antenna')









    Out[60]:





<matplotlib.text.Text at 0x10d134b10>



In [69]:

    
d12, d21 = ac2.leakage[:,:]



In [72]:

    
d12.shape









    Out[72]:





(128, 128)

But leakage isn't....

Accessing calibration solutions stored in a parset

Parsets are just text files



In [92]:

    
ls -l *.in









    



lrwxr-xr-x  1 ban115  at-aipsusr  24 28 Oct 16:12 rndgains.in@ -> ../functests/rndgains.in

Here's what a parset-formatted gain solution looks like



In [101]:

    
with open('rndgains.in', 'rU') as f:
    lines = f.readlines()
print '*'*20, 'First 15 lines', '*'*20
print ''.join(lines[0:15])

print '\n\n', '*'*20, 'Last 15 lines', '*'*20
print ''.join(lines[-15:])









    



******************** First 15 lines ********************

# This is an automatically generated file with random complex gains
# for 36 antennae and 4 polarisation products
# 1 feeds will be simulated

gain.g11.0.0 = [0.352738,0.0]
gain.g22.0.0 = [0.199759,0.0]
leakage.d12.0.0 = [-0.0769765,-0.161752]
leakage.d21.0.0 = [-0.0598945,-0.222759]
gain.g11.1.0 = [-0.720074,-0.854905]
gain.g22.1.0 = [0.0643493,-0.717078]
leakage.d12.1.0 = [-0.0758581,-0.0900622]
leakage.d21.1.0 = [-0.0250298,0.27892]
gain.g11.2.0 = [0.728084,0.96271]
gain.g22.2.0 = [0.765411,0.340206]



******************** Last 15 lines ********************
gain.g22.32.0 = [0.954835,-0.507665]
leakage.d12.32.0 = [0.0428212,-0.255874]
leakage.d21.32.0 = [0.0718756,-0.0382147]
gain.g11.33.0 = [0.844343,-0.441219]
gain.g22.33.0 = [0.747404,-0.592294]
leakage.d12.33.0 = [-0.0419436,0.021918]
leakage.d21.33.0 = [-0.0363359,0.028795]
gain.g11.34.0 = [-0.900695,-0.0145212]
gain.g22.34.0 = [0.319942,1.22673]
leakage.d12.34.0 = [0.0991746,0.00159891]
leakage.d21.34.0 = [0.0675752,0.259098]
gain.g11.35.0 = [0.343066,-1.10439]
gain.g22.35.0 = [-0.65907,0.260133]
leakage.d12.35.0 = [0.0464101,-0.149402]
leakage.d21.35.0 = [0.271098,-0.107002]

If the data is contained in a parset, you can get it by setting type='parset' in the Source constructor:



In [73]:

    
src = Source('rndgains.in', type='parset')



In [107]:

    
ac3 = src.most_recent_accessor

Note that the leakage slice is still huge (128x128), even though there are only 36 antennas and 1 beam in the parset:



In [108]:

    
d12, d21 = ac3.leakage[:,:]



In [109]:

    
d12.shape









    Out[109]:





(128, 128)

Note how plotting abs(d12) gives only 36 valid values, but plotting angle(d12) gives 128 antennas - somehow the masking is not working. Oh well, caveat emptor.



In [110]:

    
subplot(2,2,1)
plot(abs(d12[:, 0]))
ylabel('Amplitude')
title('d12')
subplot(2,2,3)
plot(np.degrees(np.angle(d12[:, 0])))
ylabel('Phase (deg)')

xlabel('Antenna')
subplot(2,2,2)
plot(abs(d21[:, 0]))
title('d21')
subplot(2,2,4)
plot(degrees(angle(d21[:, 0])))
xlabel('Antenna')









    Out[110]:





<matplotlib.text.Text at 0x1107e7ed0>