Skip to content

Commit

Permalink
add use_all_stars option to pp_calibrate
Browse files Browse the repository at this point in the history
  • Loading branch information
@mommermi
mommermi committed Aug 19, 2019
1 parent 93079e0 commit 5bfc4fe
Show file tree
Hide file tree
Showing 3 changed files with 87 additions and 54 deletions.
110 changes: 63 additions & 47 deletions catalog.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -180,7 +180,7 @@ def add_fields(self, field_names, field_arrays, field_types=None):

def download_catalog(self, ra_deg, dec_deg, rad_deg,
max_sources, save_catalog=False,
max_mag=21):
max_mag=21, use_all_stars=False):
"""
download existing catalog from VIZIER server using self.catalogname
input: ra_deg, dec_deg, rad_deg, max_sources, (display_progress),
Expand Down Expand Up @@ -263,7 +263,8 @@ def download_catalog(self, ra_deg, dec_deg, rad_deg,
self.data['mag'] = self.data['rp1mag'] # use rmag for astrometry

# clip self.data to enforce magnitude error limits
self.data = self.data[self.data['e_rp1mag'] <= 0.03]
if not use_all_stars:
self.data = self.data[self.data['e_rp1mag'] <= 0.03]

# --------------------------------------------------------------------
# use astroquery vizier query for SkyMapper
Expand Down Expand Up @@ -318,7 +319,8 @@ def download_catalog(self, ra_deg, dec_deg, rad_deg,
# self.data.rename_column('e_zPSF', 'e_zmag')
# e_zPSF does not seem to exist...

self.data = self.data[self.data['e_rmag'] <= 0.03]
if not use_all_stars:
self.data = self.data[self.data['e_rmag'] <= 0.03]

# --------------------------------------------------------------------

Expand Down Expand Up @@ -430,11 +432,12 @@ def download_catalog(self, ra_deg, dec_deg, rad_deg,

# filter columns to only have really good detections
# see the Vizier webpage for a description of what the flags mean
Qflags = set('ABC') # only A, B, or C flagged detections
qmask = [True if not set(item).difference(Qflags) else False
for item in self.data['Qflg']]
# filter columns to only have really good detections
self.data = self.data[qmask]
if not use_all_stars:
Qflags = set('ABC') # only A, B, or C flagged detections
qmask = [True if not set(item).difference(Qflags) else False
for item in self.data['Qflg']]
# filter columns to only have really good detections
self.data = self.data[qmask]

# rename column names using PP conventions
self.data.rename_column('_2MASS', 'ident')
Expand Down Expand Up @@ -620,21 +623,22 @@ def download_catalog(self, ra_deg, dec_deg, rad_deg,
return 0

# apply some quality masks
try:
mask_primary = self.data['mode'] == 1
mask_clean = self.data['clean'] == 1
mask_star = self.data['type'] == 6
mask_bright = self.data['fiberMag_g'] < max_mag
mask = mask_primary & mask_clean & mask_star & mask_bright
except TypeError:
if self.display:
print('no data available from {:s}'.format(
self.catalogname))
logging.error('no data available from {:s}'.format(
self.catalogname))
return 0

self.data = self.data[mask]
if not use_all_stars:
try:
mask_primary = self.data['mode'] == 1
mask_clean = self.data['clean'] == 1
mask_star = self.data['type'] == 6
mask_bright = self.data['fiberMag_g'] < max_mag
mask = mask_primary & mask_clean & mask_star & mask_bright
except TypeError:
if self.display:
print('no data available from {:s}'.format(
self.catalogname))
logging.error('no data available from {:s}'.format(
self.catalogname))
return 0

self.data = self.data[mask]

# rename column names using PP conventions
self.data.rename_column('objID', 'ident')
Expand Down Expand Up @@ -970,7 +974,7 @@ def read_database(self, filename):
def lin_func(self, x, a, b):
return a*x + b

def transform_filters(self, targetfilter):
def transform_filters(self, targetfilter, use_all_stars=False):
"""
transform a given catalog into a different filter band; crop the
resulting catalog to only those sources that have transformed magnitudes
Expand Down Expand Up @@ -1007,13 +1011,16 @@ def transform_filters(self, targetfilter):
self['e_umag'].data])

# sort out sources that do not meet the C&G requirements
keep_idc = (mags[1]-mags[2] > 0.08) & (mags[1]-mags[2] < 0.5) & \
(mags[0]-mags[1] > 0.2) & (mags[0]-mags[1] < 1.4) & \
(mags[0] >= 14.5) & (mags[0] < 19.5) & \
(mags[1] >= 14.5) & (mags[1] < 19.5) & \
(mags[2] >= 14.5) & (mags[2] < 19.5)
filtered_mags = np.array([mags[i][keep_idc]
for i in range(len(mags))])
if not use_all_stars:
keep_idc = ((mags[1]-mags[2] > 0.08) & (mags[1]-mags[2] < 0.5) &
(mags[0]-mags[1] > 0.2) & (mags[0]-mags[1] < 1.4) &
(mags[0] >= 14.5) & (mags[0] < 19.5) &
(mags[1] >= 14.5) & (mags[1] < 19.5) &
(mags[2] >= 14.5) & (mags[2] < 19.5))
filtered_mags = np.array([mags[i][keep_idc]
for i in range(len(mags))])
else:
filtered_mags = mags

# ... derive a linear best fit and remove outliers (>3 sigma)
ri = np.array(filtered_mags[1]) - np.array(filtered_mags[2])
Expand Down Expand Up @@ -1103,7 +1110,10 @@ def transform_filters(self, targetfilter):
self['e_imag'].data])

# sort out sources that do not meet the C&G requirements
keep_idc = (mags[0]-mags[1] > 0.08) & (mags[0]-mags[1] < 0.5)
if not use_all_stars:
keep_idc = (mags[0]-mags[1] > 0.08) & (mags[0]-mags[1] < 0.5)
else:
keep_idc = [True]*len(mags[0])

# transformed magnitudes; uncertainties through Gaussian and C&G2008
nmags = np.array([np.empty(len(mags[0])),
Expand Down Expand Up @@ -1169,6 +1179,9 @@ def transform_filters(self, targetfilter):
if mags[0][idx] > 18 or mags[1][idx] > 17:
keep = False

if use_all_stars:
keep = True

if keep:
# 0.064 (sig: 0.035) is a systematic offset
# between UKIRT_Z and SDSS_Z
Expand Down Expand Up @@ -1374,9 +1387,10 @@ def transform_filters(self, targetfilter):
# transform magnitudes to Johnson-Cousins, Vega system
# using http://gea.esac.esa.int/archive/documentation/GDR2/Data_processing/chap_cu5pho/sec_cu5pho_calibr/ssec_cu5pho_PhotTransf.html

colormask = ((self.data['BPmag']-self.data['RPmag'] > -0.5) &
(self.data['BPmag']-self.data['RPmag'] < 2.75))
self.data = self.data[colormask]
if not use_all_stars:
colormask = ((self.data['BPmag']-self.data['RPmag'] > -0.5) &
(self.data['BPmag']-self.data['RPmag'] < 2.75))
self.data = self.data[colormask]

g = self.data['Gmag'].data
e_g = self.data['e_Gmag'].data
Expand Down Expand Up @@ -1425,17 +1439,18 @@ def transform_filters(self, targetfilter):
# transform magnitudes to Johnson-Cousins, Vega system
# using http://gea.esac.esa.int/archive/documentation/GDR2/Data_processing/chap_cu5pho/sec_cu5pho_calibr/ssec_cu5pho_PhotTransf.html

if targetfilter == 'g':
colormask = ((self.data['BPmag']-self.data['RPmag'] > -0.5) &
(self.data['BPmag']-self.data['RPmag'] < 2.0))
elif targetfilter == 'r':
colormask = ((self.data['BPmag']-self.data['RPmag'] > 0.2) &
(self.data['BPmag']-self.data['RPmag'] < 2.7))
elif targetfilter == 'i':
colormask = ((self.data['BPmag']-self.data['RPmag'] > 0) &
(self.data['BPmag']-self.data['RPmag'] < 4.5))
if not use_all_stars:
if targetfilter == 'g':
colormask = ((self.data['BPmag']-self.data['RPmag'] > -0.5) &
(self.data['BPmag']-self.data['RPmag'] < 2.0))
elif targetfilter == 'r':
colormask = ((self.data['BPmag']-self.data['RPmag'] > 0.2) &
(self.data['BPmag']-self.data['RPmag'] < 2.7))
elif targetfilter == 'i':
colormask = ((self.data['BPmag']-self.data['RPmag'] > 0) &
(self.data['BPmag']-self.data['RPmag'] < 4.5))

self.data = self.data[colormask]
self.data = self.data[colormask]

g = self.data['Gmag'].data
e_g = self.data['e_Gmag'].data
Expand Down Expand Up @@ -1500,8 +1515,9 @@ def match_with(self, catalog,
astropy.table functionality; how about astropy.coord matching?
"""

this_tree = spatial.KDTree(list(zip(self[match_keys_this_catalog[0]].data,
self[match_keys_this_catalog[1]].data)))
this_tree = spatial.KDTree(
list(zip(self[match_keys_this_catalog[0]].data,
self[match_keys_this_catalog[1]].data)))

# kd-tree matching
if tolerance is not None:
Expand Down
7 changes: 5 additions & 2 deletions doc/functions.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -252,7 +252,7 @@ the logical order:
the measured FWHMs.


.. function:: pp_calibrate ([-minstars int/float], [-catalog string], [-filter string], [-maxflag integer], [-instrumental], [-solar], images)
.. function:: pp_calibrate ([-minstars int/float], [-catalog string], [-filter string], [-maxflag integer], [-instrumental], [-solar], [-use_all_stars], images)

photometric calibration of each input frame in one specific filter

Expand Down Expand Up @@ -293,7 +293,10 @@ the logical order:
using the PANSTARRS, APASS, and SDSS catalogs. The
threshold of solar-like colors is defined by the
`_pp_conf.solcol` parameter; the default is the
actual color index +- 0.2 mag.
actual color index +- 0.2 mag.
:param use_all_stars: if used, no quality checks are performed on
calibration stars and all stars are used in the
calibration.
:param images: images to run `pp_calibrate` on


Expand Down
24 changes: 19 additions & 5 deletions pp_calibrate.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -59,7 +59,8 @@ def create_photometrycatalog(ra_deg, dec_deg, rad_deg, filtername,
preferred_catalogs,
min_sources=_pp_conf.min_sources_photometric_catalog,
max_sources=1e4, mag_accuracy=0.1,
solar=False, display=False):
solar=False, use_all_stars=False,
display=False):
"""create a photometric catalog of the field of view"""

for catalogname in preferred_catalogs:
Expand All @@ -68,6 +69,7 @@ def create_photometrycatalog(ra_deg, dec_deg, rad_deg, filtername,
# load catalog
n_sources = cat.download_catalog(ra_deg, dec_deg, rad_deg,
max_sources,
use_all_stars=use_all_stars,
save_catalog=True)

if display:
Expand All @@ -83,7 +85,7 @@ def create_photometrycatalog(ra_deg, dec_deg, rad_deg, filtername,
return None

# reject non-solar colors, if requested by user
if solar:
if solar and not use_all_stars:
sol_gr = 0.44 # g-r
sol_ri = 0.11 # r-i
n_rejected = 0
Expand All @@ -100,7 +102,8 @@ def create_photometrycatalog(ra_deg, dec_deg, rad_deg, filtername,
n_rejected += cat.reject_sources_with(
(cat['rmag']-cat['imag']) > sol_ri+_pp_conf.solcol)
elif 'PANSTARRS' in cat.catalogname:
cat.transform_filters('g') # derive Sloan griz
cat.transform_filters('g')
# derive Sloan griz
n_rejected += cat.reject_sources_with(
(cat['_gmag']-cat['_rmag']) < sol_gr-_pp_conf.solcol)
n_rejected += cat.reject_sources_with(
Expand Down Expand Up @@ -152,7 +155,8 @@ def create_photometrycatalog(ra_deg, dec_deg, rad_deg, filtername,
('GAIA' in catalogname and
filtername not in {'G', 'RP', 'BP'})):

n_transformed = cat.transform_filters(filtername)
n_transformed = cat.transform_filters(filtername,
use_all_stars=use_all_stars)
if n_transformed == 0:
raise ValueError(('unable to transform {:s} to {:s}'.format(
cat.catalogname, filtername) +
Expand Down Expand Up @@ -202,6 +206,7 @@ def create_photometrycatalog(ra_deg, dec_deg, rad_deg, filtername,


def derive_zeropoints(ref_cat, catalogs, filtername, minstars_external,
use_all_stars=False,
display=False, diagnostics=False):
"""derive zeropoint for a number of catalogs based on a reference catalog"""

Expand Down Expand Up @@ -470,6 +475,7 @@ def fchi2(zp): return np.sum([(zp-residuals)**2/residuals_sig])
def calibrate(filenames, minstars, manfilter, manualcatalog,
obsparam, maxflag=3,
magzp=None, solar=False,
use_all_stars=False,
display=False, diagnostics=False):
"""
wrapper for photometric calibration
Expand Down Expand Up @@ -555,6 +561,7 @@ def calibrate(filenames, minstars, manfilter, manualcatalog,
ref_cat = create_photometrycatalog(ra_deg, dec_deg, rad_deg,
filtername, preferred_catalogs,
max_sources=2e4, solar=solar,
use_all_stars=use_all_stars,
display=display)

if ref_cat == None:
Expand Down Expand Up @@ -634,7 +641,9 @@ def calibrate(filenames, minstars, manfilter, manualcatalog,

# match catalogs and derive magnitude zeropoint
zp_data = derive_zeropoints(ref_cat, catalogs, filtername,
minstars, display=display,
minstars,
use_all_stars=use_all_stars,
display=display,
diagnostics=diagnostics)

# zp_data content
Expand Down Expand Up @@ -680,6 +689,9 @@ def calibrate(filenames, minstars, manfilter, manualcatalog,
parser.add_argument('-solar',
help='restrict to solar-color stars',
action="store_true", default=False)
parser.add_argument('-use_all_stars',
help='ignore all quality checks and use all stars',
action="store_true", default=False)
parser.add_argument('images', help='images to process', nargs='+')
args = parser.parse_args()
minstars = float(args.minstars)
Expand All @@ -689,6 +701,7 @@ def calibrate(filenames, minstars, manfilter, manualcatalog,
instrumental = args.instrumental
man_magzp = args.magzp
solar = args.solar
use_all_stars = args.use_all_stars
filenames = args.images

# manfilter: None: instrumental magnitudes, False: no manfilter provided
Expand Down Expand Up @@ -720,4 +733,5 @@ def calibrate(filenames, minstars, manfilter, manualcatalog,
calibration = calibrate(filenames, minstars, manfilter,
manualcatalog, obsparam, maxflag=maxflag,
magzp=man_magzp, solar=solar,
use_all_stars=use_all_stars,
display=True, diagnostics=conf.diagnostics)

0 comments on commit 5bfc4fe

Please sign in to comment.