photometry

This repo is intended to store all the archive codes and part of the results of (155140) 2005 UD observation at Seoul National University Astronomical Observatory 1-m telescope on 2018-10-12 & 2018-10-13, which are published in Ishiguro M. et al. 2019 (in prep).

NOTE

I decided not to upload the results and the explanation PDF files because the paper is not submitted yet. I will upload these before the final publication. (Currently they're at the local repository but I haven't yet add/commit/push'ed.)

Notebooks

Notebooks are made to give you how I reduced the data.

Note, however, that because I am not a well trained coder, it's messy. Especially because I don't know how to deal with the memory leak problem of matplotlib (or is it astropy?), I had to halt the program from time to time and re-run the notebook, as the memory gets full after few hundred images were analyzed. Thus, especially the outputs in Photometer is not what you'll get if you run everythin from the beginning. Please refer to phot_log_all.txt to see all the outputs (I copied-and-pasted outputs from each run). I know it's crazy but this was the best way I could find (facepalm). If you know a better way, please let me know. I'll be more than happy.

Notebook	Explanations	link
Reducer	preprocessing, CR rejection, and WCS implementation.	link
Photometer	ephemerides/PS1 query, centroid stars, Pillbox Aperture photometry, zeropoint fitting, target Circular Aperture photometry, standardization, and reduced magnitude calculation.	link
latex_generator	To automatically generate the LaTeX report file	link
light_curve_plotter	Simple plotting	link

Explanations on the files

• phot_targ_vis_inspected.csv: Which you may be interested in. See below.
• phot_targ.csv and vis_insp.csv: original photometry and visual inspection list.
• summary_reduced_all.csv: reduced files explanation for my purpose.
• ephem.csv: JPL HORIZONS ephemerides.
• SNUO1Mpy, ysfitsutilpy, and ysphotutilpy : The frozen snapshots of my personal packages. All these are publicly available via GitHub, but I recommend you to use the frozen in versions to completely reproduce my results. I used photutils v0.7. Due to some internal API issue with PillBoxAperture object, you cannot use photutils version lower than this...
• Reducer.ipynb, Reducer.html: preprocessing, CR rejection, and WCS implementation.
• shell scripts for astrometry (auto-generated by Reducer)
• Photometer.ipynb, Photometer.html, phot_log_all.txt: ephemerides/PS1 query, centroid stars, Pillbox Aperture photometry, zeropoint fitting, target Circular Aperture photometry, standardization, and reduced magnitude calculation. All the logs are saved in phot_log_all.txt (I did it in a rather manual copy-and-paste manner).
• latex_generator.ipynb: To automatically generate the LaTeX report file
• light_curve_plotter.ipynb: Simple plotting

The following files

1. 2005UD_photometry_Interpretation.ppdx/pdf

   • The PDF/PPTX files which explains how I did the photometry uwing pill-box aperture. There were some "failing" cases, and I tried to explain those as honestly as possible.
   • For written explanation, please see below.

2. report.pdf (1.1GB)

   • The pdf file which contains all the information of field stars, our target, and the zero-point fitting.

3. target_and_zeropoint_only.pdf (92MB)

   • Same as report but only the target and zeropoint images. You may check what "visual inspections" I did by comparing "visually rejected images" names and figures in this file.

phot_targ_vis_inspected.csv explained

phot_targ_vis_inspected.csv will be the one you are interested in. The columns include

column	explanation	comments
(all of the columns explained here)	The columns from JPL HORIZONS ephemerides.	Queried by astroquery.jplhorizons
jd_target	The lighttime corrected UT time on the object at the mid-exposure time.
trail	The trail length from ephemerides in pixel unit. Derived from the ephemerides, not from the image analysis.
xcentroid, ycentroid	Gaussian2D fitted centroid of the object.
fwhm, fwhm_arcsec	The seeing FWHM estimated from field stars in pixel and arcsec units*.
n_stars	number of stars used for the FWHM estimation and zero-point calculation.
isnear	if there is any PS1 object near the target (critical separation = 6*FWHM)	Using DAOGROUP algorithm implemented to photutils
xcenter, ycenter	just duplicate of xcentroid, ycentroid
aperture_sum, aperture_sum_err	The aperture sum and its error (including Poisson, read noise, and 1% uncertainty from flat error, but not the Poisson term from dark frames)
msky, ssky	the estimated modal value of the sky (using SExtractor estimator formula, which is a modified version of Doodson A 1917, Biometrika) and sample standard deviation from 3-sigma 5-iteration clippings
nrej, nsky	The number of rejected sky pixels from the above clipping, and the final number of sky pixels used for sky estimation
source_sum, source_sum_err	The sky-subtracted version of the aperture sum and its error (on top of aperture_sum_err, sky estimation error (i.e., 1-sigma CI of msky) and sky standard deviation noise are propagated as people do in IRAF)
mag, merr	The instrumental mag and its error (assuming 1st order Taylor expansion approximation using source_sum_err)
zp_fit, dzp	The zero point and its error estimated from n_stars field stars using weighted mean method.	math: weighted mean/error = minimum chi^2 fit and error.
mzp, zpStd	The zero point as zp_fit but using simple mean (no weighting) and its sample standard deviation. zpStd is not included in the dm_red below.
m_std	Thus standardized magnitude (its error is not shown because it's identical to dm_red)
m_red, dm_red	The reduced magnitude and its error (ontop of merr, dzp is propagated as 1st order approximation)	mzp and zpStd are not used (zp_fit and dzp are used)
seroius	0 = no problem visible; 1 = something happened but benevolent (faint star in the annulus but not aperture, etc); 2 = seriously affected (field object in aperture)
vis_insp_comment	A = Strange centroiding, B = PS1-uncatalogued object nearby, C = PS1-catalogued object nearby (isnear is 1)

Original Script in the Older Draft of the Paper

In the original draft of our paper, we had the following explanation on our photometry, but it was removed over time (as it was too much for this paper).

click to see

The photometric data were preprocessed in the standard manner for CCD data. The original object images were bias- and dark-subtracted and flat-fielded using the dome flat. The cosmic ray was then removed using the L. A. Cosmic algorithm \citep{2001PASP..113.1420V} implemented by \texttt{astroscrappy}\footnote{\url{https://github.com/astropy/astroscrappy} version 1.0.8 with a separable median filter, which is not available in the original \texttt{L. A. Cosmic} algorithm but accelerates the cosmic ray rejection calculation by more than an order of magnitude, with specifically tuned parameters.}. Then, the World Coordinate System (WCS) information was appended in each image header by the offline version of \texttt{astrometry.net} \citep{2010AJ....139.1782L}. After preprocessing, we queried the Pan-STARRS1 DR1 (hereafter DR1) catalog \citep{2016arXiv161205243F} $ r $ magnitude in the range between 10.0 and 15.2 mag, and toggled flags if there were DR1-cataloged objects near the target asteroid to avoid contamination of photometric signal of the asteroid. In addition, we discarded the extracted objects from the DR1 catalog if any pairs of stars were close to each. We only selected objects that were (1) not recognized as quasar, galaxy, or variables based on DR1 catalog's flags and (2) observed several times in the shorter wavelengths (at least three times for the $g$- and $r$-band and once for the $i$-band). Finally, we were left with a minimum of 5 to a maximum of more than 20 stars in each image. The magnitudes of these stars were used for the photometric calibration as explained below.

The aperture shape for each star is defined as a pill-box, a combination of a rectangle and two half-ellipses, similar to \texttt{TRIPPy} \citep{2016AJ....151..158F}. The height of the rectangle is where it faces the half-ellipses. The position angle of the aperture is obtained by fitting the 2-dimensional elliptical Gaussian to field stars with the initial guess from the ephemerides, and the sigma-clipped median of the angles of the field stars is used for the aperture position angles. After testing a number of combinations of the semi-major and minor axis lengths of the half-ellipses, we empirically determined the appropriate solution of half-circles with a radius $ 1.75F $, where $ F $ is the full-width at half-maximum (FWHM) of the point sources, to enclose a sufficient amount of the stellar signal even when the tracking accuracy of the telescope mount was not perfect. Therefore, the aperture was set as a combination of a rectangle with width $ L $, the expected trail length of the asteroid with respect to the field stars during the exposure time, retrieved from the ephemerides, and height $ 3.50 F $, with two half-circles with radius $ 1.75F $. The instrumental magnitude was then calculated by subtracting the sky value estimated from the locally defined pill-box annulus with inner and outer radii of $ 4F $ and $ 6F $, respectively, while keeping the same width ($L$) of the rectangle. For the asteroid, which is the target being tracked, we used a circular aperture of radius $2F$ and similarly selected the circular annulus with inner and outer radii of $ 4F $ and $ 6F $ to determine the local sky flux. Changing the sizes of the apertures affected the results only within much less than the estimated 1-$\sigma$ error bars.

The magnitudes in DR1 were converted into Johnson--Cousins $R_\mathrm{C}$ magnitudes by the transformation formula given in \cite{2012ApJ...750...99T}. By using the instrumental and catalog magnitudes, we determined the photometric zero point at each image. The instrumental magnitudes were then converted to the standardized magnitude using the photometric zero point. We ignored the color-term for the atmospheric extinction, which appeared to be negligible based on the stars of $ 0 \lesssim g-r \lesssim 1 $ from our analysis (zero point slope $ \lesssim 0.05 $ for $ g-r $).

% The below can be an item for discussion??? %This would not be true especially at the airmass $ X \gtrsim 2 $. However, the scatter due to the low signal-to-noise (S/N) ratios for both the asteroid and field stars are more dominant than the color-dependent term; it is visible in Figure \ref{fig:lightcurve} with large scatter and large error-bars (e.g., the rotational phase $\sim 0.7 $ in 2018 October 13 data, where the airmass was as high as $ \sim 2.6 $).

The observed $R_\mathrm{C}$ magnitudes were converted into reduced magnitudes (hypothetical magnitudes at the unit heliocentric distance of 1 au and the observer's distance of 1 au), which is given by,

$$ m_\mathrm{R}(1,1,\alpha)=R - 5\log(r_\mathrm{h} \Delta) ~, $$

\noindent where $r_\mathrm{h}$ and $\Delta$ are the heliocentric and the observer's distances in au at the epoch of our observation. Since our photometric data were taken at the opposition (i.e., $\alpha \sim 0^\circ$), we ignored the $\alpha$-dependency of the magnitude and derived the absolute magnitude $H_\mathrm{R} := m_\mathrm{R}(1,1,0)$ in $R_\mathrm{C}$-band. %Note that the absolute magnitude is formally defined in $V$-band but we regard $H_\mathrm{V} = H_\mathrm{R}$ To obtain the lightcurve, we further corrected the light time to take into account subtle rotation of the asteroid while light traveled to the observer's location.

Finally, we visually inspected each image with the locations of photometric apertures to check if our photometric results were affected by unexpected problems (such as close encounters with background objects that are not listed in the DR1 or imperfect centering of objects due to the low S/N, etc.). Of 810 images, 41 images were discarded due to such unexpected problems. In the period analysis, data points with a large scatter (less than 10 data points) were automatically rejected. The source codes, script files, output files and ancillary explanations for the detailed algorithm for the photometric analyses, except for the period analysis, are publicly available via GitHub service\footnote{\url{https://github.com/ysBach/IshiguroM_etal_155140_2005UD/tree/master/photometry}, including the frozen snapshots of the dependent packages developed for our purpose (\texttt{ysfitsutilpy}, \texttt{ysphotutilpy}, and \texttt{SNUO1Mpy}).}.

Data Analysis Strategy

1. Centroiding. I haven't yet found perfect algorithm to do centroiding, and I guess IRAF has some highly subjectively tuned internal algorithm for the centroiding which I haven't found. From short experience, I found elliptical Gaussian fitting for a centroid box of cbox = 8 * FWHM_initial_guess. In this code, I did put mask = data < data.min() + 5*ssky. For each iteration, new Gaussian fitting is done for cbox-sized box, and the centroid position is updated. If the distance to the new centroid from the old one is larger than 1/2 FWHM, move only 1/2 FWHM to that direction. Halt the iteration when shift is << 1 pixel, and iteration is done up to 10 times only.

2. FWHM. It's difficult to find the true FWHM from non-sidereal tracking mode, since stars are trailed. What I did is, after the centroiding to the field stars (PS1 catalogued stars) as described above, get the sigma_y (minor axis) values of the fittings, and get sigma-clipped average of that. And regard 2*sqrt(2)*sigma_y as the FWHM. In reality, trailed PSF (e.g., trailed Gaussian) is rather pill-box shaped and thus the Gaussian fit may be different from the real one. See VerešP+ 2012 AJ (ads) and Fraser+2016 (paper, TRIPPy github).

3. Target Aperture. For target, I used circular r_ap = 1.75*FWHM aperture. The sky annulus is defined as inner and outer radii by 4 and 6 FWHM, respectively. Normally it is recommended by classical IRAF to use r_ap = max(3, f_ap*FWHM) for f_ap=1.5 to 3 depending on the user's preference. The reason for minimum 3 pixel radius is because IRAF uses approximated aperture sum and if r_ap < 3, the non-random error in aperture photometry gets larger than about 0.1 mag (for 1 pixel radius, it gets up to 0.5 mag error). As python or IDL codes has no such problem, I didn't put lower limit on aperture. Also 1.75 is selected by my experience using the SNUO 1m telescope.

4. Star Aperture. Since stars are trailed, I developed a code to do pill-box aperture photometry. It consists of two half-ellipses (a=major, b=minor, theta=angle) with one rectangle (width=trail length, height=b of the ellipse, theta=theta of the ellipse). The theta is obtained from the sigma-clipped mean value of 2D Gaussian fits as I used in FWHM estimation above. The trail length is simply obtained from ephemerides, and a = b = 1.75 FWHM is used. If you see the aperture overplotted on stellar objects, you realize it is too large compared to visually acceptable aperture. I intended this because the tracking error sometimes make the stellar PSF be U- or V-shaped, so non-negligible stellar fluxes were lost for some images while target flux does not. Thus, sacrifising the photometric accuracy of stellar photometry, I enlarged aperture sizes.

5. Catalog Magnitude. The PS1 query was done by VizieR (only has PS1 DR1 data, though DR2 is already available). The conversion from PS1 filter to Johnson-Cousins R filter was adopted from Tonry et al. 2012. Only stars with PS1 r-mag between 10-15.2 were used with n_observation in PS1 >= 3(r), 3(g), 1(i). Extended objects were removed by one of the simplest method recommended by Pan-STARRS.

   For details which can be subjective, please see the Photometer in the zip file and the source codes in the package ysphotutilpy (ysphotutilpy.queryutil.PanSTARRS1 and its methods drop_for_diff_phot, select_filters, check_nearby, drop_star_groups. I tried to explain as detailed as possible in the source codes.)

6. Standardization. The standardization was done by comparing the R-instrumental mag (R_inst) and catalog mag (R_cat) of the field stars. The linearity curve (R_inst VS R_cat) and its residual, as well as the color-term plot (R_cat - R_inst VS g_PS1 - r_PS1) were used to visually inspect whether the stars were selected well.

7. Visual Inspection. It's virtually impossible to look at each image thoroughly. Thus, I let my source code to draw "report figure" for each star and target. Then automatically generate report LATEX files and compile them. While looking at the generated report PDF file, I found some images failed in correct centroiding or some coding bug, because (1) my code is not perfect yet, (2) objects were too faint to be fitted correctly, etc. This kind of "error-rate" is which I should resolve in the future, but for now I just ignored them (2 images had unexpected serious problems which I cannot find reason).

8. Data Goodness. For the light curve below, "good" means (1) visually not affected at all (even serious = 1 are removed) and (2) zpStd < 0.04 is used as just a heuristic criterion.

If you need, I can provide all the data (original, combined calibration frames, reduced object frames, etc) and intermediate products for your reference. All files are already uploaded to my personal web cloud.

And if you have any question or request, please don't hesitate to contact me.

Yoonsoo P. Bach, ysbach93@gmail.com