forked from deapplegate/wtgpipeline
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathadam_make_backmask_ims.py
executable file
·177 lines (163 loc) · 6.62 KB
/
adam_make_backmask_ims.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
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
#! /usr/bin/env python
#adam-does# makes BACKMASK directory with images used to mask reflections, cross-talk, areas of coherent offsets in astrometry, and maybe even asteroids
# mask in difference coadds, but place regions in SUPA*_[0-9]OCF.fits images, so regions are placed on individual chip exposures with backmask.py
# SUPAXXXXXXX-all.coadd.fits : difference images for each exposure ("coadd_cluster_SUPAXXXXXXX"-"coadd_cluster_all")
# SUPAXXXXXXX-all.coadd.smoothed.fits: smoothed difference images
#adam-call_example# adam_make_backmask_ims.py /nfs/slac/g/ki/ki18/anja/SUBARU/MACS0416-24/W-C-RC/SCIENCE/
#adam-comments# very useful
#BEFORE RUNNING CODE (1): set make_smoothed if you want smoothed images for all types of coadd differences (all, gabodsid, gabrot)
#BEFORE RUNNING CODE (2): set make_gabodsid and make_gabrot to determine the types of coadd differences you want to make (makes "all" always), can make gabodsid & gabrot as well
make_smoothed=0
make_gabodsid=0
make_gabrot=0
from numpy import *
import scipy.ndimage as ndimage
import matplotlib.pyplot as plt
import glob
import sys ; sys.path.append('/u/ki/awright/InstallingSoftware/pythons')
import imagetools
import astropy, astropy.io.fits as pyfits
import os
import re
input_dir=sys.argv[1]
input_dir=input_dir.replace('//','/')
if not input_dir.endswith('/'):
input_dir+='/'
if not os.path.exists(input_dir):
raise Exception("input_dir=sys.argv[1]:" +input_dir+ " isn't a directory")
dir_splits=input_dir.split('/')
try:
dir_splits.remove('')
dir_splits.remove('')
except ValueError:
pass
cluster=dir_splits[-3]
filter=dir_splits[-2]
dir_backmask=input_dir+'BACKMASK/'
if not os.path.exists(dir_backmask):
os.makedirs(dir_backmask)
supa_fls=glob.glob(input_dir+'coadd_'+cluster+'_SUPA*/coadd.fits')
regex_cluster=imagetools.getRegex(cluster)
if make_gabodsid==1:
gabodsid_fls=glob.glob(input_dir+'coadd_'+cluster+'_gabodsid[0-9][0-9][0-9][0-9]/coadd.fits')
gabodsid_dict={}
for fl in gabodsid_fls:
coadd_gabodsid_dir=fl.split('/')[-2]
fl_match=re.match('coadd_'+regex_cluster+'_gabodsid([0-9][0-9][0-9][0-9])',coadd_gabodsid_dir)
gabodsid=fl_match.group(1)
gabodsid_dict[gabodsid]=fl
if make_gabrot==1:
gabrot_fls=glob.glob(input_dir+'coadd_'+cluster+'_gab[0-9][0-9][0-9][0-9]-rot[0-1]/coadd.fits')
gabrot_dict={}
for fl in gabrot_fls:
coadd_gabrot_dir=fl.split('/')[-2]
fl_match=re.match('coadd_'+regex_cluster+'_gab([0-9][0-9][0-9][0-9])-rot([0-1])',coadd_gabrot_dir)
gab=fl_match.group(1)
rot=fl_match.group(2)
gabrot='gab'+gab+'-rot'+rot
gabrot_dict[gabrot]=fl
#SHNT: now make a dict that matches keys I'll get from the coadd_fl's I'm looping over
# to the matching gabrot_fl and gabodsid_fl
all_fls=glob.glob(input_dir+'coadd_'+cluster+'_all/coadd.fits')
all_fl=all_fls[0]
all_head=pyfits.open(all_fl)[0].header
all_img=imagetools.GetImage(all_fl)
all0s=all_img==0
smoothed_fls=[]
place_reg_fls=[]
orig_dir=input_dir.replace('/SCIENCE/','')+"_*/SCIENCE/"
for supa_fl in supa_fls:
supa_dir=os.path.dirname(supa_fl)
supa=supa_dir[-11:]
split_fl_str=orig_dir+supa+'_[0-9]*OCF.fits'
split_fls=glob.glob(split_fl_str)
split_fl=split_fls[0]
if len(split_fls)!=10:
raise Exception("split_fl: "+split_fl+" isn't returning 10 matches")
split_head=pyfits.open(split_fl)[0].header
gab=str(split_head['GABODSID'])
rot=str(split_head['ROTATION'])
gabrot='gab'+gab+'-rot'+rot
supa_img=imagetools.GetImage(supa_fl)
supa0s=supa_img==0
#make exposure-all images
im0s=logical_or(supa0s,all0s)
difference_img=supa_img-all_img
difference_img[im0s]=0
diffim=difference_img.copy()
diffim[im0s]=nan
hdu=pyfits.PrimaryHDU(asarray(diffim,dtype=float32))
hdu.header=all_head
difference_fl=dir_backmask+supa+'-all.coadd.fits'
hdu.writeto(difference_fl,clobber=True,output_verify='ignore')
# make smoothed image
if make_smoothed!=0:
smooth_img = ndimage.gaussian_filter(difference_img, sigma=3)
smooth_img[im0s]=nan
hdu=pyfits.PrimaryHDU(asarray(smooth_img,dtype=float32))
hdu.header=all_head
smoothed_fl=difference_fl.replace('coadd','coadd.smoothed')
hdu.writeto(smoothed_fl,clobber=True,output_verify='ignore')
print 'wrote 2 files: '+smoothed_fl+'\n\tand '+difference_fl
smoothed_fls.append(smoothed_fl)
place_reg_fls.append(split_fl_str)
else: print 'wrote file: '+difference_fl
#make gabodsid-all images
if make_gabodsid==1:
if gabodsid_dict.has_key(gab):
gabodsid_fl=gabodsid_dict[gab]
gabodsid_img=imagetools.GetImage(gabodsid_fl)
gabodsid0s=gabodsid_img==0
im0s=logical_or(supa0s,gabodsid0s)
difference_img=supa_img-gabodsid_img
difference_img[im0s]=0
diffim=difference_img.copy()
diffim[im0s]=nan
hdu=pyfits.PrimaryHDU(asarray(diffim,dtype=float32))
hdu.header=all_head
difference_fl=dir_backmask+supa+'-gabodsid'+gab+'.coadd.fits'
hdu.writeto(difference_fl,clobber=True,output_verify='ignore')
# make smoothed image
if make_smoothed!=0:
smooth_img = ndimage.gaussian_filter(difference_img, sigma=3)
smooth_img[im0s]=nan
hdu=pyfits.PrimaryHDU(asarray(smooth_img,dtype=float32))
hdu.header=all_head
smoothed_fl=difference_fl.replace('coadd','coadd.smoothed')
hdu.writeto(smoothed_fl,clobber=True,output_verify='ignore')
print 'wrote 2 files: '+smoothed_fl+'\n\tand '+difference_fl
smoothed_fls.append(smoothed_fl)
place_reg_fls.append(split_fl_str)
else: print 'wrote file: '+difference_fl
else: print "gabodsid "+gab+" not found!"
#make gabodsid_rot-all images
if make_gabrot==1:
if gabrot_dict.has_key(gabrot):
gabrot_fl=gabrot_dict[gabrot]
gabrot_img=imagetools.GetImage(gabrot_fl)
gabrot0s=gabrot_img==0
im0s=logical_or(supa0s,gabrot0s)
difference_img=supa_img-gabrot_img
difference_img[im0s]=0
diffim=difference_img.copy()
diffim[im0s]=nan
hdu=pyfits.PrimaryHDU(asarray(diffim,dtype=float32))
hdu.header=all_head
difference_fl=dir_backmask+supa+'-'+gabrot+'.coadd.fits'
hdu.writeto(difference_fl,clobber=True,output_verify='ignore')
# make smoothed image
if make_smoothed!=0:
smooth_img = ndimage.gaussian_filter(difference_img, sigma=3)
smooth_img[im0s]=nan
hdu=pyfits.PrimaryHDU(asarray(smooth_img,dtype=float32))
hdu.header=all_head
smoothed_fl=difference_fl.replace('coadd','coadd.smoothed')
hdu.writeto(smoothed_fl,clobber=True,output_verify='ignore')
print 'wrote 2 files: '+smoothed_fl+'\n\tand '+difference_fl
smoothed_fls.append(smoothed_fl)
place_reg_fls.append(split_fl_str)
else: print 'wrote file: '+difference_fl
else: print "gabrot "+gabrot+" not found!"
if make_smoothed!=0:
for place_reg_fl,smoothed_fl in zip(place_reg_fls,smoothed_fls):
print "ipython -i -- backmask.py %s %s" % (smoothed_fl,place_reg_fl)