-
Notifications
You must be signed in to change notification settings - Fork 0
/
Interpol.py
70 lines (51 loc) · 1.88 KB
/
Interpol.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
# -*- coding: utf-8 -*-
#@author: astrolitterbox
import numpy as np
import pyfits
from scipy import ndimage
from scipy.interpolate import RectBivariateSpline
#from PIL import Image
import scipy as sp
import matplotlib.pyplot as plt
import math
#import cv
#from scipy.spatial import KDTree
imgFile = pyfits.open('fpC-005115-r2-0103.fit')
maskFile = pyfits.open('data/masks/IC3376_mask_r.fits')
imgFile.info()
img = imgFile[0].data
img = img - 1000 #soft bias subtraction, comment out if needed
head = imgFile[0].header
mask = maskFile[0].data
maskedImg = np.ma.array(img, mask = mask)
NANMask = maskedImg.filled(np.NaN)
badArrays, num_badArrays = sp.ndimage.label(mask)
print num_badArrays
data_slices = sp.ndimage.find_objects(badArrays)
#data_slices = np.ma.extras.notmasked_contiguous(maskedImg)
import inpaint
#def replace_nans(array, max_iter, tol,kernel_size=1,method='localmean'):
#data = np.array([[11,22.0,353],[1,np.NaN,343],[11,0,343]])
filled = inpaint.replace_nans(NANMask, 5, 0.5, 2, 'idw')
hdu = pyfits.PrimaryHDU(filled, header = head)
hdu.writeto('IC3376R_masked.fits')
#hdu = pyfits.PrimaryHDU(NANMask)
#hdu.writeto('NanMask.fits')
#invdisttree = idw.Invdisttree(maskedImg, maskedImg)
#print data_slices
#for slice in data_slices:
# dy, dx = slice
# coords = float(dy.start), float(dy.stop), float(dx.start), float(dx.stop)
# print 'slice', coords
# submask = cutoutMask[(dy.start):(dy.stop),(dx.start):(dx.stop)]
# subarray = cutout[(dy.start):(dy.stop+1),(dx.start):(dx.stop+1)]
#print np.where(submask == 1)
#print maskCoords, 'maskCoords'
#y = np.where(cutoutMask == 1)[0]
#x = np.where(cutoutMask == 1)[1]
#print y, x, 'xy'
#points = np.zeros((y.shape[0], 2))
#for i in range(0, y.shape[0]-1):
# points[i] = [y[i], x[i]]
# print points
#ret = sp.ndimage.interpolation.map_coordinates(maskedImg, points, output=None, order=3, mode='nearest', cval=0.0, prefilter=True)