Subarray determination


In [1]:
%pylab qt


Populating the interactive namespace from numpy and matplotlib

In [2]:
rcParams['image.interpolation'] = 'nearest'
rcParams['image.cmap'] = 'hot'

In [3]:
import pyfits

reading in filenames that are fits


In [4]:
fitsfiles = !ls /Users/maye/data/irtf/bigdog/scan4/limb*.fits

In [5]:
# plot 1 test image
data = pyfits.getdata(fitsfiles[0])
imshow(data)
plt.axis('off')


Out[5]:
(-0.5, 1023.5, 649.5, -0.5)

The next plot determines the best top row for where the data starts. As soon as the rows become indistinguishable no refraction from the slit should be visible anymore.


In [10]:
for line in range(245,255):
    plot(data[line, :], label=line)
legend(loc='best')


Out[10]:
<matplotlib.legend.Legend at 0x10dc26b10>

We believe row 247 (the red line) still shows too many differences so we should include only above it. Now the same for the higher row number:


In [11]:
for line in range(390, 400):
    plot(data[line, :], label=line)
legend(loc='best')


Out[11]:
<matplotlib.legend.Legend at 0x10dc26ad0>

It looks like up to row 394 (0-indexing) the data is reliable. Creating an appropriate subframe reader:


In [6]:
def read_scan1(fname):
    data = pyfits.getdata(fname, 0)
    return data[248:394, :]

In [7]:
data = read_scan1(fitsfiles[0])

In [8]:
data.shape


Out[8]:
(146, 1024)

Limb master

Above shape is the shape of one spectral image. Now reading in all data into a 3D cube.


In [9]:
cube = []
for fitsfile in fitsfiles:
    # fish out the image number from the filename
    nr = int(fitsfile.split('-')[1][:5])
    cube.append(read_scan1(fitsfile))
    if nr > 104:
        break
# making a contiguous C-array out of the slow list container
cube = np.array(cube)
cube.shape


Out[9]:
(51, 146, 1024)

Above dimensions indicate that I have stacked 105 spectra with each 146x1024 pixels.


In [11]:
limb_master = cube.mean(axis=0)
limb_master.shape


Out[11]:
(146, 1024)

Again, the dimensions indicate the correct result that I averaged the first axis, and therefore have a 'moon_master' image:


In [12]:
imshow(limb_master, aspect='auto')
colorbar(orientation='horizontal')
savefig('limb_master.png', dpi=200)

In [18]:
pwd


Out[18]:
u'/Users/maye/Dropbox/src/diviner/notebooks'

In [13]:
from scipy.ndimage import median_filter

In [15]:
imshow(limb_master, aspect='auto')
colorbar(orientation='horizontal')
title('Limb master raw')
savefig('limb_master_original.png', dpi=200)

In [16]:
limb_master_filtered = median_filter(limb_master, 2)

In [17]:
imshow(limb_master_filtered, aspect='auto')
colorbar(orientation='horizontal')
title('Limb master median filter (size 2)')
savefig('limb_master_median_filtered.png', dpi=200)

Sky master

Starting to read from 107 (106 might be okay as well).


In [23]:
cube = []
for fitsfile in fitsfiles:
    # fish out the image number from the filename
    nr = int(fitsfile.split('-')[1][:5])
    if nr < 107:
        continue
    cube.append(read_scan1(fitsfile))
# making a contiguous C-array out of the slow list container
sky_cube = np.array(cube)
sky_cube.shape


Out[23]:
(33, 146, 1024)

In [24]:
sky_master = sky_cube.mean(axis=0)
sky_master.shape


Out[24]:
(146, 1024)

In [145]:
imshow(sky_master, aspect='auto', vmax=20000)
colorbar(orientation='horizontal')
title('sky master, raw, clipping at 20000')
savefig('sky_master_raw.png', dpi=200)

In [26]:
# do a minimal pepper/salt filter
sky_master_filtered = median_filter(sky_master, 3)

In [144]:
imshow(sky_master_filtered, aspect='auto', vmax=20000)
colorbar(orientation='horizontal')
title('sky master, median size 3 filter, clipping at 20000')
savefig('sky_master_median_filtered.png', dpi=200)

Taking an upper limit of 20000 counts to exclude hot pixels.


In [28]:
sky_master_filtered[sky_master_filtered > 20000]= np.nan

Regions

For comparison with Raquel's analysis, I also will split the spatial region into 4 bands, each with 36 pixels width. First, let's substract the sky.


In [29]:
master = moon_master_filtered - sky_master_filtered

In [18]:
rcParams['image.aspect'] = 'auto'

In [31]:
imshow(master)
colorbar(orientation='horizontal')
title('Background("sky")-subtracted master')
savefig('master_subtracted.png',dpi=200)

Now, having the master, let's cut out regions.


In [32]:
master.shape


Out[32]:
(146, 1024)

In [33]:
region1 = master[:36, :]
region2 = master[37:73, :]
region3 = master[74:110, :]
region4 = master[110:,:]

In [34]:
for region in [region1,region2,region3,region4]:
    print region.shape


(36, 1024)
(36, 1024)
(36, 1024)
(36, 1024)

In [35]:
spec1 = region1.mean(axis=0)
spec2 = region2.mean(axis=0)
spec3 = region3.mean(axis=0)
spec4 = region4.mean(axis=0)

In [36]:
for i,spec in enumerate([spec1,spec2,spec3,spec4]):
    plot(spec, label=i+1)
legend(loc='best')
title('master spectra for each region')
savefig('master_spectra_in_counts.png',dpi=200)

In [37]:
rcParams['lines.linewidth'] = 1.5

Regions before subtracting.


In [19]:
moon_regions = []
moon_regions.append(limb_master_filtered[:36, :])
moon_regions.append(limb_master_filtered[37:73, :])
moon_regions.append(limb_master_filtered[74:110, :])
moon_regions.append(limb_master_filtered[110:,:])

In [47]:
sky_regions = []
sky_regions.append(sky_master_filtered[:36, :])
sky_regions.append(sky_master_filtered[37:73, :])
sky_regions.append(sky_master_filtered[74:110, :])
sky_regions.append(sky_master_filtered[110:, :])

In [20]:
waves = pd.read_excel('lucy_spectral_cal.xlsx', parse_cols=0)

In [22]:
for i, moon_region in enumerate(moon_regions):
    plot(waves, moon_region.mean(axis=0), 
                 label='limb{}'.format(i+1))
title('Limb regions master spectra')
legend(loc='best')
xlabel('Wavelengths [micron]')
ylabel('Digital counts unscaled')
savefig('limb_separated_regions.png', dpi=200)

In [24]:
moon_masters = []
for moon in moon_regions:
    moon_masters.append(moon.mean(axis=0))

In [26]:
for t in [(3,4),(2,4),(1,4), (2,3), (1,3), (1,2)]:
    plot(waves, moon_masters[t[0]-1]/moon_masters[t[1]-1],
                 label='{}/{}'.format(*t))
title('Limb regional master spectra ratios')
legend(loc='best')
xlabel("Wavelengths [micron]")
savefig('individual_limb_ratios.png', dpi=200)

In [132]:
fig, axes = subplots(nrows=1)
for i in range(4):
    axes.plot(waves, moon_masters[i]/sky_masters[i], label=i)
axes.legend(loc='best')


Out[132]:
<matplotlib.legend.Legend at 0x11f7adf90>

In [137]:
for i in range(0,1000, 100):
    plot(moon_master[:, i:i+100].mean(axis=1), label=i)
legend(loc='best')


Out[137]:
<matplotlib.legend.Legend at 0x126a1b950>

In [ ]: