In [1]:
from kbmodpy import kbmod as kb
from trajectoryFiltering import *
import numpy as np
import random as rd
import math
import matplotlib.pyplot as plt
import os
%matplotlib inline
In [2]:
path = '../../fraser/chip_7/'
In [3]:
files = os.listdir(path)
In [4]:
max_imgs = 5
In [5]:
files.sort()
files = [path+f for f in files]
files = files[:max_imgs]
files
Out[5]:
In [6]:
images = [kb.layered_image(f) for f in files]
Object Generation Parameters
In [7]:
p = kb.psf(1.5)
object_count = 10
x_range = (5,1750)
y_range = (5, 3650)
angle_range = (0.1, 0.5)
velocity_range = (200, 300)
flux_range = (800, 1500)
Search Parameters
In [8]:
min_observations = 25
angle_steps = 120
velocity_steps = 100
search_margin = 1.2
In [9]:
print(angle_range[1]-angle_range[0])
print(velocity_range[1]-velocity_range[0])
Generate Random trajectories within bounds
In [10]:
results_key = []
for _ in range(object_count):
traj = kb.trajectory()
traj.x = int(rd.uniform(*x_range))
traj.y = int(rd.uniform(*y_range))
ang = rd.uniform(*angle_range)
vel = rd.uniform(*velocity_range)
traj.x_v = vel*math.cos(ang)
traj.y_v = vel*math.sin(ang)
traj.flux = rd.uniform(*flux_range)
results_key.append(traj)
Add a real object to the list of trajectories
In [11]:
if path == '../../HITS/test_35/4,6tempExp/new_header/':
real_result = kb.trajectory()
real_result.flux = 5300
real_result.x = 3123
real_result.y = 3043
real_result.x_v = 2425
real_result.y_v = 1050
results_key.append(real_result)
results_key
Out[11]:
Test that clustering is not able to collapse together too many unique trajectories
In [12]:
len(cluster_trajectories(results_key, dbscan_args=dict(eps=0.007, min_samples=1))[1])
Out[12]:
Use the generated trajectories to add objects into the images
In [13]:
for t in results_key:
add_trajectory(images, t, p)
Sanity check
In [14]:
len(match_trajectories(results_key, results_key, 0.01, 1)[0])
Out[14]:
In [15]:
stack = kb.image_stack(images)
Mask out stars and bad pixels
In [16]:
flags = ~0 # mask pixels with any flags
flag_exceptions = [32,39] # unless it has one of these special combinations of flags
master_flags = int('100111', 2) # mask any pixels which have any of
# these flags in more than two images
In [17]:
stack.apply_mask_flags(flags, flag_exceptions)
In [18]:
stack.apply_master_mask(master_flags, 2)
In [19]:
images = [i.science() for i in stack.get_images()]
Calculate masked percentage of an image to estimate probablilty of placing an object under a mask
In [20]:
img = images[1]
percent_masked = img[np.where(img==-9999.99)].size/img.size
percent_masked
Out[20]:
In [21]:
np.shape(images[1])
Out[21]:
In [ ]:
fig = plt.figure(figsize=(12,12))
plt.imshow(images[0] , origin='lower', vmin=-50, vmax=200)#cmap=plt.cm.Greys_r,
plt.xlabel('X Pixels')
plt.ylabel('Y Pixels')
plt.colorbar()
In [ ]:
search = kb.stack_search(stack, p)
search.set_debug(True)
In [ ]:
search_ang_r = (angle_range[0]/search_margin,
angle_range[1]*search_margin)
search_vel_r = (velocity_range[0]/search_margin,
velocity_range[1]*search_margin)
#search.gpu(angle_steps,velocity_steps, *angle_range, *velocity_range, min_observations)
In [ ]:
#search.get_results(0,20)
In [ ]:
#search.region_search(2400, 1040, 50, 255, 3)
res = search.region_search(2000, 1000, 800, 60, 3)
res
In [ ]:
matched = []
for r in res:
if any(abs(r.ix-t.x)<=1 and abs(r.iy-t.y)<=1 for t in results_key ):
matched.append(r)
len(matched)
In [ ]:
matched
In [ ]:
pooled = search.get_psi_pooled()
In [ ]:
pooled_imgs = [np.array(im, copy=False) for im in pooled[1]]
In [ ]:
fig = plt.figure(figsize=(12,12))
plt.imshow(pooled_imgs[0] [2700:3600,2900:3800], origin='lower', vmin=-1, vmax=1)#cmap=plt.cm.Greys_r,
plt.xlabel('X Pixels')
plt.ylabel('Y Pixels')
plt.colorbar()
In [ ]:
fig = plt.figure(figsize=(12,12))
plt.imshow(pooled_imgs[1] [1350:1800,1450:1900], origin='lower', vmin=-1, vmax=1)#cmap=plt.cm.Greys_r,
plt.xlabel('X Pixels')
plt.ylabel('Y Pixels')
plt.colorbar()
In [ ]:
fig = plt.figure(figsize=(12,12))
plt.imshow(pooled_imgs[3] [300:500,300:500], origin='lower', vmin=-1, vmax=1)#cmap=plt.cm.Greys_r,
plt.xlabel('X Pixels')
plt.ylabel('Y Pixels')
plt.colorbar()
In [ ]: